boot.ci() 为引导样本返回奇怪的置信区间

Question

我有两组正实数。

> dput(group1)
c(2.10753, 2.57251, 2.61687, 4.62551, 7.13166, 6347.63, 4.22139, 
10.7373, 2.11568, 2.71866, 4.09376, 10.9046, 109807, 5.87156, 
3.17082, 3.4703, 2.47262, 9.24319, 34.6945, 5.72567, 12.0134, 
108.33, 6.60707, 6.24304, 3.59048, 10.3174, 48.0265, 5.32097, 
3.77157, 6.67401, 22.633, 34.8186, 21.5315, 9.42882, 7.10627, ...)

> dput(group2)
c(4.88474, 65.4318, 128.101, 24.1271, 5.44262, 54.8987, 2.85175, 
14.1089, 172.23, 66.8563, 6.74067, 2.19603, 2.12985, 4.12735,
16.401, 3.22688, 15.6943, 4.32861, 36.4752, 7.33769, 75.855, 
62.7653, 35.1786, 3.71099, 29.0186, 34.4472, 19.1061, 2.75174, ...)

第 1 组包含约 1000 个值，第 2 组约 30,000 个。我对两组之间的中位数比率感兴趣，并使用以下 R 函数为 2000 个引导样本中的每一个计算该比率（对于 boot() 命令，请参见下面的函数输出）：

medianRatio <- function(x, i, noGroup1, noGroup2) {
    all <- x[i]
    currGroup1 <- all[1:noGroup1]
    currGroup2 <- all[c(noGroup1 + 1):length(all)]
    ratio <- median(currGroup1) / median(currGroup2)
    return(ratio)
}

boot() 函数的调用如下所示

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = c(group1, group2), statistic = medianRatio, R = noBs, 
    noGroup1 = length(group1), noGroup2 = length(group2))


Bootstrap Statistics :
    original      bias    std. error
t1*  1.08847 -0.08597889  0.05451763`

得到的 bootstrap 样本分布的平均值为 1.002，sd 为 0.054（直方图的目测证实了 1 附近的正态分布）。另外：

范围(Group1_Group2.BS$t) [1] 0.823311 1.198469

然而，当我在引导对象上运行 boot.ci() 时，报告的置信区间为

BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 2000 bootstrap replicates

CALL : 
boot.ci(boot.out = Group1_Group2.BS, type = "norm")

Intervals : 
Level      Normal        
95%   ( 1.068,  1.281 )  
Calculations and Intervals on Original Scale

我不明白这里发生了什么，因为报告的置信区间甚至没有涵盖引导样本的（对称）分布的平均值。我错过了什么？

Answer 1

了解如何计算正态置信区间可能有助于理清思路。我在this question的回答中找到了一个很好的解释。

bootstrap 正态置信区间是围绕统计量的观察值构建的，并带有偏差校正，以解决 bootstrap 统计量的中间值与观察值之间的差异。 CI 的计算方法是使用来自原始样本的感兴趣统计量的估计值，减去偏差，并加上 1.96 倍的引导标准误差。如果您使用来自 boot 对象的值“手动”执行此操作，您会看到您获得的值与来自 boot.ci 的值相同。

1.08847 - -0.08597889 - 1.96*0.05451763
[1] 1.067594
1.08847 - -0.08597889 + 1.96*0.05451763
[1] 1.281303

您的偏差大到足以使百分位 CI 和正常 CI 之间的差异变得明显。当我考虑 bootstrap 置信区间时，我希望使用 bootstrap 分布的平均值减去偏差，而不是 观察到的统计量减去偏差。我没有花足够的时间考虑这个问题，也没有带上我的引导笔记，但你可以更仔细地检查一下这个问题。

Answer 2

您正在错误地执行引导程序。您在函数 medianRatio 中假设 x[i] 的第一部分将保存第 1 组的值和第 2 组的第二部分。这是不正确的，因为 i 只是一个从 1 替换为 noGroup1+noGroup2 的样本。您需要在引导命令中使用strata 选项。我认为以下将起作用。 stype="f" 选项指定 boot 将生成一个长度为 noGroup1+noGroup2 的向量，该向量表示为 bootstrap 样本选择了每个观察值的次数（0,1,2,...）。此代码基于启动文档中的示例。

medRatio <- function(x, f, noGroup1){
 gp1 <- 1:noGroup1
 m1 <- median(rep(x[gp1],f[gp1]))
 m2 <- median(rep(x[-gp1],f[-gp1]))
 return(m1/m2)
}
n1 <- length(group1)
n2 <- length(group2)
boot(c(group1,group2),medRatio,R=2000,stype = "f",strata = rep(1:2,c(n1,n2)),noGroup1=n1)