是否可以将大小不等的样本传递给 R 中的函数启动答案

【问题标题】：Is it possible to pass samples of unequal size to function boot in R是否可以将大小不等的样本传递给 R 中的函数启动
【发布时间】：2013-08-15 15:11:44
【问题描述】：

我目前正在R 中编写有关引导的教程。我选择了boot 包中的函数boot。我得到了 Efron/Tibshirani (1993) 的《Bootstrap 简介》一书，并复制了他们的一些示例。

在这些示例中，它们经常根据不同的样本计算统计数据。例如，他们有一个例子，他们有 16 只老鼠的样本。其中7只小鼠接受了旨在延长试验手术后存活时间的治疗。其余9只小鼠未接受治疗。收集每只老鼠的存活天数（数值如下）。

现在，我想使用引导方法来确定均值的差异是否显着。但是，如果我正确理解了boot 的帮助页面，我就不能只将两个样本大小不相等的不同样本传递给函数。我的解决方法如下：

#Load package boot
library(boot)
#Read in the survival time in days for each mouse
treatment <- c(94, 197, 16, 38, 99, 141, 23)
control   <- c(52, 104, 146, 10, 51, 30, 40, 27, 46)
#Call boot twice(!)
b1 <- boot(data = treatment,
           statistic = function(x, i) {mean(x[i])},
           R = 10000)
b2 <- boot(data = control,
           statistic = function(x, i) {mean(x[i])},
           R = 10000)
#Compute difference of mean manually
mean_diff <- b1$t -b2$t

在我看来，这个解决方案有点小题大做。我感兴趣的统计数据现在保存在向量mean_diff 中，但我不再获得boot 包的所有强大功能。我无法在mean_diff 等上拨打boot.ci 等。

所以我的问题基本上是我的 hack 是否是使用 R 中的 boot 包和比较两个不同样本的统计数据进行引导的唯一方法。还是有其他方法？

我考虑过传入一个包含 16 行和一个附加列“Group”的 data.frame：

df <- data.frame(survival=c(treatment, control), 
                 group=c(rep(1, length(treatment)), rep(2, length(control))))
head(df)
  survival group
1       94     1
2      197     1
3       16     1
4       38     1
5       99     1
6      141     1

但是，现在我必须告诉boot，它必须始终从前 7 行中抽取 7 个观测值，从后 9 行中抽取 9 个观测值，并将它们视为单独的样本。我不知道该怎么做。

【问题讨论】：

你为什么不用t.test?????
@Dwin 我知道我可以运行t.test(df$survival ~ df$group) 来替代boot。但是，这不是我的问题（实际上我的教程中有这一部分）。问题是关于我想将引导程序应用于比较两个样本的统计数据的一般情况。均值检验的差异只是一个例子。或者你有没有想到结合t.test 和boot 的东西？在这种情况下，如果您能分享该解决方案，那就太好了，因为我不太明白怎么做。
我在想，如果你能正确地对抽样进行分层，你可以使用 t.test(...)$t 作为启动统计数据。

标签： r

【解决方案1】：

我从来没有真正弄清楚 boot 的最大优势是什么，因为手动编写引导程序代码非常容易。例如，您可以使用replicate 尝试以下操作：

myboot1 <- function(){
    booty <- tapply(df$survival,df$group,FUN=function(x) sample(x,length(x),TRUE))
    sapply(booty,mean)
}
out1 <- replicate(1000,myboot1())

从中你可以很容易地得到一堆有用的统计数据：

rowMeans(out1) # group means
diff(rowMeans(out1)) # difference
mean(out1[1,]-out1[2,]) # another way of getting difference
apply(out1,1,quantile,c(0.025,0.975)) # treatment-group CIs
quantile(out1[1,]-out1[2,],c(0.025,0.975)) # CI for the difference

【讨论】：

正如我所写，boot.ci 非常简洁。此外，我的“hack”看起来并不比你的复杂，所以boot 并不是在使用不同样本的情况下使事情变得更复杂（即使我的 hack 是最好的解决方案）。但是，当然，你有一个观点。还有其他运行引导程序的解决方案（比如你的好那个），它们并不比使用boot 复杂。
实际上，感谢您的回答，我找到了一种将其与boot 结合的方法，所以+1。谢谢！

【解决方案2】：

这是?boot.return中的一个例子：

diff.means <- function(d, f)
{    n <- nrow(d)
     gp1 <- 1:table(as.numeric(d$series))[1]
     m1 <- sum(d[gp1,1] * f[gp1])/sum(f[gp1])
     m2 <- sum(d[-gp1,1] * f[-gp1])/sum(f[-gp1])
     ss1 <- sum(d[gp1,1]^2 * f[gp1]) - (m1 *  m1 * sum(f[gp1]))
     ss2 <- sum(d[-gp1,1]^2 * f[-gp1]) - (m2 *  m2 * sum(f[-gp1]))
     c(m1 - m2, (ss1 + ss2)/(sum(f) - 2))
}
grav1 <- gravity[as.numeric(gravity[,2]) >= 7,]
boot(grav1, diff.means, R = 999, stype = "f", strata = grav1[,2])

可以参考戴维森和欣克利的Section3.2。

【讨论】：

我认为这个例子只适用于大小相等的样本，但让它更通用应该很简单。我认为只需要一个gp2 <- seq(from=table(as.numeric(d$series))[1]+1,by=1, length.out=table(as.numeric(d$series))[2]) 并相应地调整其余部分。那谢谢啦。然而，这个例子让我想知道f 实际上是什么。我认为它类似于sample(1:26, 26, TRUE)，但这个例子似乎表明它是一个布尔向量。现在我有点困惑...试图查看源代码，但这对于周五早上来说看起来太复杂了....
它可以是数字索引、逻辑或权重的向量。您使用整数向量的减号从数据框或矩阵中删除项目。请注意，如果这些是合乎逻辑的，那么“-”就不是正确的否定方法。

【解决方案3】：

再想一想，我意识到我实际上可以将 Thomas 的回答与 boot 结合起来。这是一个解决方案：

b <- boot(data=df, 
           statistic = function(x, i) {
             booty <- tapply(x$survival,x$group,FUN=function(x) sample(x,length(x),TRUE))
             diff(sapply(booty,mean))*-1
           },
           R=10000)

诀窍是你提供给参数statistic 的函数必须接受一个参数 i 作为索引，但是你在你的函数中完全忽略了这个参数。相反，您自己进行采样。当然，这不是最有效的（因为boot 也必须进行采样），但我想在大多数情况下这应该不是什么大问题。

【讨论】：