使用循环找到P（-1.5

Question

我从R中的标准正态分布创建了一个大小为n = 20的iid随机样本。我有一个m = 10 ^ 5场景的蒙特卡罗模拟。我试图找到样本的过度峰度介于-1.5和1.5之间的概率。到目前为止，我有这个代码：

set.seed(1234)
n<-10
m<-10^5

sample.normal<- rep(0,m)
for(i in (1:m)){
x<-rnorm(n,0,1)
sample.normal[i]<-kurtosis(x,na.rm=FALSE)-3
}

#finding p(-1.5<y2<1.5)
y<-sample.normal
Z<-ecdf(y)
Z(1.5)-Z(-1.5)

但是，现在我想将样本大小更改为n = 10和n = 200之间的值范围，以显示样本大小如何增加，更多值介于-1.5和1.5之间。所以我想在10,20,30,40 ...... 200处计算Z（1.5）-Z（-1.5）（记住m = 10 ^ 5）（我最终表明样本量更大） ，过度峰度倾向于标准法线的零点。如何使用循环或不同方法做到这一点？任何建议，将不胜感激。

Answer 1

如果要将代码自动化为n的其他值，则以下函数简化了模拟，然后在sapply的n = 10, 20, ..., 200循环中调用它。

请注意，其中一个参数是rdist，默认为rnorm。如果你想将参数传递给RNG函数rdist使用点...参数。

library(e1071)

set.seed(1234)

simKurtosis <- function(n, m = 10^5, rdist = rnorm, ...){
  y <- replicate(m, kurtosis(rdist(n, ...)) - 3)
  Z <- ecdf(y)
  Z(1.5) - Z(-1.5)
}

ex_kurt <- sapply(seq(10, 200, by = 10), simKurtosis)