验证中心极限定理

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45

library(moments)
options(digits = 3) # 保留小数点后3位
options(scipen = 200) # 200个数字以内都不使用科学计数法

# 原始发布为0-1分布数据(概率为0.5, 期望为0.5, 方差为0.25), 可设置为任何分布
X = rbinom(10000, size = 1, prob = 0.5)
X_mean = mean(X) # 原始分布的均值
X_var= var(X) # 原始分布的方差



m=10000 # 抽取的次数
par(mfrow=c(2,3)) # 两行三列绘图

list_result = list()


for(n in c(5,10,50,100,500,1000)){ # 每次抽样的样本数
  
  s=c()
  
  for(j in 1:m){ # 抽样m次
    sample=sample(X,n) # 抽取n个样本
    sample_mean = mean(sample) # 求样本的均值
    s[j]=sample_mean
  }
  
  s_mean = mean(s) # 样本均值的均值
  s_var = var(s) # 样本均值的方差
  # 绘制样本均值的直方图, 将数据划分为 50 个区间
  hist(s,breaks = 50,main = paste('n = ',n,sep='')) 
  
  result =c(
    n,
    X_mean,
    s_mean,
    X_var,
    s_var,
    X_var/n,
    skewness(s),# 偏度
    kurtosis(s)) # 峰度
  
  cat(result)
  cat('\n')
}