在下面的代码中,我使用自举来计算CI和零值假设下的p值,即施用于番茄植物的两种不同肥料对植物产量没有影响(另一种选择是"改良"肥料更好).第一个随机样品(x)来自使用标准肥料的植物,而"改良"的样品已经用于第二个样品(y)来自的植物中.
x <- c(11.4,25.3,29.9,16.5,21.1)
y <- c(23.7,26.6,28.5,14.2,17.9,24.3)
total <- c(x,y)
library(boot)
diff <- function(x,i) mean(x[i[6:11]]) - mean(x[i[1:5]])
b <- boot(total, diff, R = 10000)
ci <- boot.ci(b)
p.value <- sum(b$t>=b$t0)/b$R
我不喜欢上面的代码是重新取样完成,好像只有一个11个样本的样本(将前5个分离为样本x,其余为样本y).您能否告诉我如何修改此代码,以便从第一个样本中替换大小为5的重新采样,并从第二个样本中单独重新采样大小为6,以便引导程序重新采样将模拟生成"第二个样本"的"单独样本"设计原始数据？
EDIT2:Hack被删除,因为它是一个错误的解决方案.相反,必须使用引导函数的参数层次:
total <- c(x,y)
id <- as.factor(c(rep("x",length(x)),rep("y",length(y))))
b <- boot(total, diff, strata=id, R = 10000)
...
请注意,您不会接近对您的p.value的正确估计:
x <- c(1.4,2.3,2.9,1.5,1.1)
y <- c(23.7,26.6,28.5,14.2,17.9,24.3)
total <- c(x,y)
b <- boot(total, diff, strata=id, R = 10000)
ci <- boot.ci(b)
p.value <- sum(b$t>=b$t0)/b$R
> p.value
[1] 0.5162
你如何解释两个样本的p值为0.51,其中第二个样本的所有值都高于第一个的最高值？
上面的代码可以很好地得到置信区间的偏差估计,但是关于差异的显着性检验应该通过整个数据集的排列来完成.