r - 为什么将 emmeans 与 data.frame 对比而不报告正确的 p 值？

Question

我正在运行的对比度的 p 值未正确转换为 data.frame。为什么会这样，我该如何解决？

emmeans 的控制台输出：

> pairs(emmeans(lmer.mod, ~ Status*Stim*Treatment), simple = "each")
$`simple contrasts for Status`
Stim = 1, Treatment = None:
 contrast               estimate     SE   df t.ratio p.value
 Control - Subclinical  -0.24213 0.0571 57.5 -4.241  0.0002 
 Control - Clinical     -0.16275 0.0571 57.5 -2.851  0.0164 
 Subclinical - Clinical  0.07938 0.0571 57.5  1.390  0.3526 

emmeans 的 data.frame 的控制台输出：

> mod.EMM <- pairs(emmeans(lmer.mod, ~ Status*Stim*Treatment), simple = "each")
> as.data.frame(mod.EMM)
   Stim Treatment      Status               contrast     estimate         SE       df      t.ratio      p.value
1    1      None           .  Control - Subclinical -0.242125000 0.05709000 57.46544  -4.24111052 3.680551e-03
2    1      None           .     Control - Clinical -0.162750000 0.05709000 57.46544  -2.85076195 2.721389e-01
3    1      None           . Subclinical - Clinical  0.079375000 0.05709000 57.46544   1.39034857 1.000000e+00

可重现的例子：

model1 <- lm(uptake ~ Type + Treatment + conc + Type*Treatment, data=CO2)

library(emmeans)
pairs(emmeans(model1, ~ Type*Treatment), simple="each")
# $`simple contrasts for Type`
# Treatment = nonchilled:
#   contrast             estimate   SE df t.ratio p.value
# Quebec - Mississippi     9.38 1.85 79 5.068   <.0001 
# 
# Treatment = chilled:
#   contrast             estimate   SE df t.ratio p.value
# Quebec - Mississippi    15.94 1.85 79 8.610   <.0001 
# 
# 
# $`simple contrasts for Treatment`
# Type = Quebec:
#   contrast             estimate   SE df t.ratio p.value
# nonchilled - chilled     3.58 1.85 79 1.934   0.0566 
# 
# Type = Mississippi:
#   contrast             estimate   SE df t.ratio p.value
# nonchilled - chilled    10.14 1.85 79 5.477   <.0001 

as.data.frame(pairs(emmeans(model1, ~ Type*Treatment), simple="each"))
#    Treatment        Type             contrast  estimate       SE df  t.ratio      p.value
# 1 nonchilled           . Quebec - Mississippi  9.380952 1.851185 79 5.067538 1.036140e-05
# 2    chilled           . Quebec - Mississippi 15.938095 1.851185 79 8.609670 2.252161e-12
# 3          .      Quebec nonchilled - chilled  3.580952 1.851185 79 1.934410 2.265719e-01
# 4          . Mississippi nonchilled - chilled 10.138095 1.851185 79 5.476542 1.995066e-06

model1 <- lm(uptake ~ Type + Treatment + conc + Type*Treatment, data=CO2)
pairs(emmeans(model1, ~ Type*Treatment), simple="each")
# $`simple contrasts for Type`
# Treatment = nonchilled:
#   contrast             estimate   SE df t.ratio p.value
# Quebec - Mississippi     9.38 1.85 79 5.068   <.0001 
# 
# Treatment = chilled:
#   contrast             estimate   SE df t.ratio p.value
# Quebec - Mississippi    15.94 1.85 79 8.610   <.0001 
# 
# 
# $`simple contrasts for Treatment`
# Type = Quebec:
#   contrast             estimate   SE df t.ratio p.value
# nonchilled - chilled     3.58 1.85 79 1.934   0.0566 
# 
# Type = Mississippi:
#   contrast             estimate   SE df t.ratio p.value
# nonchilled - chilled    10.14 1.85 79 5.477   <.0001 

as.data.frame(pairs(emmeans(model1, ~ Type*Treatment), simple="each"))
#    Treatment        Type             contrast  estimate       SE df  t.ratio      p.value
# 1 nonchilled           . Quebec - Mississippi  9.380952 1.851185 79 5.067538 1.036140e-05
# 2    chilled           . Quebec - Mississippi 15.938095 1.851185 79 8.609670 2.252161e-12
# 3          .      Quebec nonchilled - chilled  3.580952 1.851185 79 1.934410 2.265719e-01
# 4          . Mississippi nonchilled - chilled 10.138095 1.851185 79 5.476542 1.995066e-06

从外部帮助更新：

“看起来pairs()的结果本身不是一个可以变成数据框的emmGrid对象，而是一个包含两个emmGrid对象的列表。如果你从列表中按位置提取这些对象中的任何一个，使用[[] ]，像这样，

pairs(emmeans(model1, ~ Type*Treatment), simple = "each")[[2]]

然后你可以 data.frame() 每个结果，它会是正确的。你最终会得到两个不同的数据框来保存涉及两个不同变量的对比，但是这些数据框中的每一个都有正确的 p 值。”

我希望有人对这个问题有更好的解决方法，这样我就可以将所有对比合并到一个 data.frame 中。

标签： rdataframeemmeans

Answer 1

您看到的不同 p 值反映了未调整的 p 值与为多重比较调整的 p 值。

?emmeans::pairs文档告诉我们：

通常，当 simple 是一个列表或“每个”时，返回值是一个 emm_list 对象，其中每个条目都与 simple 的条目相对应。但是，如果 combine = TRUE，则使用 rbind.emmGrid 将元素全部组合成单个 emmGrid 对象中的一组对比。在这种情况下，adjust 参数设置组合对比集的调整方法。

因此，通过您的可重现示例，您可以将所有简单的主要效果组合到一个数据帧中，并将combine参数设置为TRUE. adjust您可以通过设置参数在未调整和调整后的 p 值之间进行选择。

model1 <- lm(uptake ~ Type + Treatment + conc + Type*Treatment, data=CO2)

> pairs(emmeans(model1, ~ Type*Treatment), simple = "each", combine = TRUE,
+               adjust = "none")
 Treatment  Type        contrast             estimate   SE df t.ratio p.value
 nonchilled .           Quebec - Mississippi     9.38 1.85 79 5.068   <.0001 
 chilled    .           Quebec - Mississippi    15.94 1.85 79 8.610   <.0001 
 .          Quebec      nonchilled - chilled     3.58 1.85 79 1.934   0.0566 
 .          Mississippi nonchilled - chilled    10.14 1.85 79 5.477   <.0001 

这是带有 Bonferroni 调整的一个：

> pairs(emmeans(model1, ~ Type*Treatment), simple = "each", combine = TRUE,
+               adjust = "bonferroni")
 Treatment  Type        contrast             estimate   SE df t.ratio p.value
 nonchilled .           Quebec - Mississippi     9.38 1.85 79 5.068   <.0001 
 chilled    .           Quebec - Mississippi    15.94 1.85 79 8.610   <.0001 
 .          Quebec      nonchilled - chilled     3.58 1.85 79 1.934   0.2266 
 .          Mississippi nonchilled - chilled    10.14 1.85 79 5.477   <.0001 

P value adjustment: bonferroni method for 4 tests