R中的中心极限定理

Question

我希望模拟中心极限定理以进行演示，但我不确定如何在 R 中进行。我想创建 10,000 个样本大小为 n 的样本（可以是数字或参数），从我将选择的分布（均匀，指数等......）。然后我想在一个图（使用 par 和 mfrow 命令）中绘制原始分布（直方图）、所有样本的均值分布、均值的 QQ 图，以及第四张图（有四个，2X2 )，我不确定要绘制什么。你能帮我开始用 R 编程吗？我想一旦我有了模拟数据，我应该没问题。谢谢。

我最初的尝试如下，太简单了，我什至不确定是否正确。

r = 10000;
n = 20;

M = matrix(0,n,r);
Xbar = rep(0,r);

for (i in 1:r)
{
  M[,i] = runif(n,0,1);
}

for (i in 1:r)
{
  Xbar[i] = mean(M[,i]);
}

hist(Xbar);

Answer 1

CLT 声明给定 i.i.d.来自具有均值和方差的分布的样本，样本均值（作为随机变量）具有随着样本数量n 的增加而收敛到高斯分布的分布。在这里，我假设您要生成包含n 样本的r 样本集，每个样本集以创建样本均值的r 样本。一些代码如下：

set.seed(123) ## set the seed for reproducibility
r <- 10000
n <- 200      ## I use 200 instead of 20 to enhance convergence to Gaussian

## this function computes the r samples of the sample mean from the 
## r*n original samples
sample.means <- function(samps, r, n) {
  rowMeans(matrix(samps,nrow=r,ncol=n))
}

为了生成绘图，我们使用来自here 的ggplot2 和Aaron 的qqplot.data 函数。我们还使用gridExtra在一帧中绘制多个图。

library(ggplot2)
library(gridExtra)
qqplot.data <- function (vec) {
  # following four lines from base R's qqline()
  y <- quantile(vec[!is.na(vec)], c(0.25, 0.75))
  x <- qnorm(c(0.25, 0.75))
  slope <- diff(y)/diff(x)
  int <- y[1L] - slope * x[1L]

  d <- data.frame(resids = vec)

  ggplot(d, aes(sample = resids)) + stat_qq() + geom_abline(slope = slope, intercept = int, colour="red") + ggtitle("Q-Q plot")  
}

generate.plots <- function(samps, samp.means) {
  p1 <- qplot(samps, geom="histogram", bins=30, main="Sample Histogram")
  p2 <- qplot(samp.means, geom="histogram", bins=30, main="Sample Mean Histogram")
  p3 <- qqplot.data(samp.means)
  grid.arrange(p1,p2,p3,ncol=2)
}

然后我们可以将这些函数与 uniform 分布一起使用：

samps <- runif(r*n)  ## uniform distribution [0,1]
# compute sample means
samp.means <- sample.means(samps, r, n))
# generate plots
generate.plots(samps, samp.means)

我们得到：

或，泊松分布，均值 = 3：

samps <- rpois(r*n,lambda=3)
# compute sample means
samp.means <- sample.means(samps, r, n))
# generate plots
generate.plots(samps, samp.means)

我们得到：

或，指数分布，均值 = 1/1：

samps <- rexp(r*n,rate=1)
# compute sample means
samp.means <- sample.means(samps, r, n))
# generate plots
generate.plots(samps, samp.means)

我们得到：

请注意，样本均值直方图的均值看起来都像 Gaussians，其均值与原始生成分布的均值非常相似，无论是均匀分布、泊松分布还是指数分布，正如 CLT 预测的那样（也它的方差将是原始生成分布的方差的 1/(n=200)。

Answer 2

也许这可以帮助您入门。我对正态分布进行了硬编码，只显示了您建议的两个图：随机选择的样本的直方图和所有样本均值的直方图。

我想我的主要建议是使用列表而不是矩阵来存储样本。

r <- 10000
my.n <- 20

simulation <- list()

for (i in 1:r) {
  simulation[[i]] <- rnorm(my.n)
}

sample.means <- sapply(simulation, mean)

selected.sample <- runif(1, min = 1, max = r)

dev.off()
par(mfrow = c(1, 2))
hist(simulation[[selected.sample]])
hist(sample.means)