在R中进行分层自助抽样，使用>25层。

Question

在R中进行分层自助抽样，使用>25层。

rrandomgroup-bystatistics-bootstrapsample-size

3

我有大约25个不同组的数据。为了看看如果我有不同的样本量，每个组的方差会如何变化，我正在尝试进行分层自助法。例如，在样本量为5的情况下，它应该为每个组产生1000个包含5个重新抽样点的集合。我希望在可能的范围内收集最小的样本大小，每组为5到30。

我遇到的问题是我必须对每个组进行子集，并在单独的组上运行自助法，然后将R输出复制并粘贴到Excel中。（我对R和编码相当新手）。这需要太长时间。我需要自动化自助法以识别组，并以某种方式将1000组的统计信息保存到数据框中。这有意义吗？

以下是我目前的代码：...

#sample data
set.seed(1234)
df <- data.frame(g.name = as.factor(sample(c(LETTERS),100, replace = T)),
            C.H = as.numeric(sample(c(1:9),100, replace=T)))

#subset data by group... here only a three examples
Agroup=subset(df,C.H=='A')
Bgroup=subset(df,C.H=='B')
Cgroup=subset(df,C.H=='C')

#Bootstrap selecting a sample size of "i", "B" number of times. i.e. I am 
selecting sample sizes from 5 to 30, 1000 times each. I then apply var() to 
the sample, and take the multiple variances(or the variance of the 
variances). C.H is the measurement ranging from 1 to 9.  

B=1000
cult.var=(NULL)
for (i in 5:30){
 boot.samples=matrix(sample(Agroup$C.H,size=B*i, 
replace=TRUE),B,i)
  cult.var[i]=var(apply(boot.samples,1,var))
}
print(cult.var)

这个方法可行，但是需要大量复制和粘贴。我认为我需要使用for循环按组进行引导抽样，或者想出其他办法。我确实找到了一种方法可以单独进行分层抽样而不进行引导抽样。所以也许我可以想出如何重复执行1000次...

使用boot()函数的示例并不适用于我的情况。我已经尝试过一些调整，但基本没有效果。我不确定如何编写函数，这也可能是我无法解决问题的原因。

- andemexoax

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Matt Tyers · Accepted Answer

这是我的一点尝试...

# generating data
set.seed(1234)
df <- data.frame(g.name = as.factor(sample(c(LETTERS),100, replace = T)),
                 C.H = as.numeric(sample(c(1:9),100, replace=T)))

boot.samples <- with(df, tapply(C.H, g.name, function(x) lapply(5:30, function(i) replicate(1000, sample(x,size=i,replace=T)))))

str(boot.samples$A)
## List of 26
##  $ : num [1:5, 1:1000] 7 7 3 7 7 7 3 3 2 7 ...
##  $ : num [1:6, 1:1000] 7 2 2 2 3 7 7 2 2 7 ...
##  $ : num [1:7, 1:1000] 2 3 2 7 2 3 7 2 3 3 ...
##  $ : num [1:8, 1:1000] 7 7 3 3 3 7 2 7 7 3 ...
##  $ : num [1:9, 1:1000] 2 2 2 7 2 7 3 3 3 7 ...
## ...and so on

variances <- lapply(boot.samples, function(y) sapply(y, function(x) apply(x, 2, var)))
    str(variances)
## List of 26
##  $ A: num [1:1000, 1:26] 3.2 5.8 6.2 3.2 0.3 4.8 5 5.8 6.7 3.2 ...
##  $ B: num [1:1000, 1:26] 3.2 0.8 4.7 5.3 5.3 5.3 1.2 4.7 4.2 3.8 ...
##  $ C: num [1:1000, 1:26] 9 4.8 2.7 9.8 8.3 9.8 10.2 10.2 9 12.3 ...
##  $ D: num [1:1000, 1:26] 8.3 7.5 9.8 3.8 3.5 3.5 5.7 3.7 6.7 3.2 ...
## ...and so on

variancesvariances <- lapply(variances, function(x) apply(x, 2, var))
str(variancesvariances)
## List of 26
##  $ A: num [1:26] 3.15 2.27 1.53 1.3 1.03 ...
##  $ B: num [1:26] 4.32 3.54 2.83 2.46 2.09 ...
##  $ C: num [1:26] 13.06 10.08 8.46 6.98 5.59 ...
##  $ D: num [1:26] 4.9 3.7 3.02 2.39 2.07 ...
## ...and so on

按照宣传，看起来会随着样本量的增加而下降...让我们做一张漂亮的图片。

cols <- rainbow(26)
plot(NA, xlim=c(1,26), ylim=c(0,max(unlist(variancesvariances))))
for(i in 1:26) {
  lines(variancesvariances[[i]], col=cols[i])
  text(1,variancesvariances[[i]][1],names(variancesvariances)[i],col=cols[i])
}

请注意，可以使用as.data.frame(variancesvariances)将其转换为数据框。

这次我理解正确了吗？