6

在下面的代码中,我使用 bootstrapping 来计算 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 的重采样,以便引导重采样模仿产生原始数据?

4

3 回答 3

5

EDIT2:Hack 被删除,因为它是一个错误的解决方案。相反,必须使用 boot 函数的参数 strata :

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,其中第二个样本的所有值都高于第一个样本的最大值?

上面的代码可以很好地得到置信区间的偏向估计,但是关于差异的显着性测试应该通过对完整数据集的排列来完成。

于 2010-09-01T07:30:32.443 回答
1

虽然在某些情况下实际的土壤床可以被认为是一个分层变量,但这不是其中之一。你只有一个操作,在植物群之间。因此,您的零假设是它们确实来自完全相同的人群。在这种情况下,将这些项目视为来自一组 11 个样本是引导的正确方法。

如果您有两个地块,并且在每个地块中以平衡的方式在不同季节尝试不同的肥料,那么这些地块将是统计样本,您希​​望将它们视为此类。但这里不是这样。

于 2010-09-01T16:34:45.143 回答
1

跟随 John,我认为使用 bootstrap 来测试这两个不同总体的总和是否显着不同的适当方法如下:

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)


b_x <- boot(x, sum, R = 10000)
b_y <- boot(y, sum, R = 10000)

z<-(b_x$t0-b_y$t0)/sqrt(var(b_x$t[,1])+var(b_y$t[,1]))
pnorm(z)

所以我们可以清楚地拒绝他们是同一种群的空值。我可能错过了一个自由度调整,我不确定自举在这方面是如何工作的,但这样的调整不会彻底改变你的结果。

于 2013-06-04T02:26:54.510 回答