R中带有aov对象的独立t检验答案

【问题标题】：Independent t-test with an aov-object in RR中带有aov对象的独立t检验
【发布时间】：2020-05-21 10:37:55
【问题描述】：

是否有一种简单的方法可以使用 aov 对象计算独立 t 检验，类似于 pairs(emmeans(aov.object,~variable) 如何计算成对比较？

【问题讨论】：

标签： r t-test

【解决方案1】：

我正在尝试辨别问题的真正含义。这个答案是基于这个猜测：

是否可以使用emmeans() 进行所有独立的t 测试，比较多个样本，一次两个？

池 t 检验：否（相似但不完全）

这个问题的答案是“不完全是”。我将用 datasets 包中的InsectSprays 数据进行说明。我们可能适合这个aov 模型：

> fm1 <- aov(count ~ spray, data = InsectSprays)

然后使用emmeans比较喷雾：

> library("emmeans")
> pairs(emmeans(fm1, "spray"), adjust = "none")

 contrast estimate  SE df t.ratio p.value
 A - B      -0.833 1.6 66 -0.520  0.6045 
 A - C      12.417 1.6 66  7.755  <.0001 
 A - D       9.583 1.6 66  5.985  <.0001 
 A - E      11.000 1.6 66  6.870  <.0001 
 A - F      -2.167 1.6 66 -1.353  0.1806 
 B - C      13.250 1.6 66  8.276  <.0001 
 B - D      10.417 1.6 66  6.506  <.0001 
 B - E      11.833 1.6 66  7.391  <.0001 
 B - F      -1.333 1.6 66 -0.833  0.4080 
 C - D      -2.833 1.6 66 -1.770  0.0814 
 C - E      -1.417 1.6 66 -0.885  0.3795 
 C - F     -14.583 1.6 66 -9.108  <.0001 
 D - E       1.417 1.6 66  0.885  0.3795 
 D - F     -11.750 1.6 66 -7.339  <.0001 
 E - F     -13.167 1.6 66 -8.223  <.0001

这确实会执行 t 测试，一次比较两种方法（添加 adjust = "none" 表示没有多重性调整，因此每个测试都按原样进行，不参考其他测试） .如果您对每对样本分别运行t.test()，则所比较的均值正是要比较的均值。然而，这些结果与t.test() 的结果并不完全相同，因为fm1 中的误差方差是通过将所有 6 个样本汇集在一起来估计的，而t.test() 将使用可变性仅在特定测试中使用的两个样本中。我将在A - B 比较中明确显示这一点：

> t.test(count ~ spray, data = InsectSprays, subset = (spray %in% c("A","B")), var.equal = TRUE)

    Two Sample t-test

data:  count by spray
t = -0.45352, df = 22, p-value = 0.6546
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -4.643994  2.977327
sample estimates:
mean in group A mean in group B 
       14.50000        15.33333

注意t的比例略有不同，有22个d.f。而不是 66，因为只合并了两个样本的方差。

韦尔奇 t 检验：是

在方差不完全相等的不同假设下，可以使用类似的方法准确地获得所有t检验。这需要一个没有像aov() 那样指定等方差的模型。下面是使用gls() 模型的相同示例：

> library("nlme")
> fm2 <- gls(count ~ spray, data = InsectSprays, weights = varIdent(form = ~ 1 | spray))

> pairs(emmeans(fm2, "spray", df.method = "boot"), adjust = "none")
 contrast estimate    SE   df t.ratio p.value
 A - B      -0.833 1.837 21.8 -0.454  0.6547 
 A - C      12.417 1.477 14.7  8.407  <.0001 
 A - D       9.583 1.542 16.7  6.214  <.0001 
 A - E      11.000 1.451 13.9  7.580  <.0001 
 A - F      -2.167 2.252 20.5 -0.962  0.3473 
 B - C      13.250 1.358 15.5  9.754  <.0001 
 B - D      10.417 1.429 17.8  7.289  <.0001 
 B - E      11.833 1.330 14.5  8.894  <.0001 
 B - F      -1.333 2.177 19.5 -0.613  0.5472 
 C - D      -2.833 0.920 20.9 -3.078  0.0057 
 C - E      -1.417 0.758 21.6 -1.868  0.0754 
 C - F     -14.583 1.882 13.2 -7.748  <.0001 
 D - E       1.417 0.879 19.6  1.612  0.1229 
 D - F     -11.750 1.934 14.5 -6.076  <.0001 
 E - F     -13.167 1.862 12.7 -7.071  <.0001 

Degrees-of-freedom method: satterthwaite

请注意，A - B 比较的测试与 t.test() 的测试相同，默认情况下将方差视为不等：

> t.test(count ~ spray, data = InsectSprays, subset = (spray %in% c("A","B")))

    Welch Two Sample t-test

data:  count by spray
t = -0.45352, df = 21.784, p-value = 0.6547
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -4.646182  2.979515
sample estimates:
mean in group A mean in group B 
       14.50000        15.33333

【讨论】：

这非常有用。谢谢！