为什么转换 emmeans 与 data.frame 的对比没有报告正确的 p 值？ 可重现的示例：

这可以很容易地完成，但您必须做的是获得基本输出，然后插入正确的 P 值。为了说明这一点，我将展示一个不同的例子，其中一个因素有两个以上的水平。

require(emmeans)
#> Loading required package: emmeans

warp.lm = lm(breaks ~ wool * tension,data = warpbreaks)
(cons = pairs(emmeans(warp.lm,~ wool * tension),simple = "each"))
#> $`simple contrasts for wool`
#> tension = L:
#>  contrast estimate   SE df t.ratio p.value
#>  A - B       16.33 5.16 48  3.167  0.0027 
#> 
#> tension = M:
#>  contrast estimate   SE df t.ratio p.value
#>  A - B       -4.78 5.16 48 -0.926  0.3589 
#> 
#> tension = H:
#>  contrast estimate   SE df t.ratio p.value
#>  A - B        5.78 5.16 48  1.120  0.2682 
#> 
#> 
#> $`simple contrasts for tension`
#> wool = A:
#>  contrast estimate   SE df t.ratio p.value
#>  L - M      20.556 5.16 48  3.986  0.0007 
#>  L - H      20.000 5.16 48  3.878  0.0009 
#>  M - H      -0.556 5.16 48 -0.108  0.9936 
#> 
#> wool = B:
#>  contrast estimate   SE df t.ratio p.value
#>  L - M      -0.556 5.16 48 -0.108  0.9936 
#>  L - H       9.444 5.16 48  1.831  0.1704 
#>  M - H      10.000 5.16 48  1.939  0.1389 
#> 
#> P value adjustment: tukey method for comparing a family of 3 estimates

# get the estimates,etc. into a data frame:
df = as.data.frame(cons)

# get the Tukey-adjusted P values:
pv = unlist(lapply(unlist(cons),function(X) as.data.frame(X)$p.value))

# replace the p values and display
df$p.value = pv
df
#>   tension wool contrast   estimate       SE df    t.ratio      p.value
#> 1       L    .    A - b 16.3333333 5.157299 48  3.1670322 0.0026768025
#> 2       M    .    A - B -4.7777778 5.157299 48 -0.9264108 0.3588672592
#> 3       H    .    A - B  5.7777778 5.157299 48  1.1203107 0.2681556374
#> 4       .    A    L - M 20.5555556 5.157299 48  3.9857208 0.0006572745
#> 5       .    A    L - H 20.0000000 5.157299 48  3.8779987 0.0009185485
#> 6       .    A    M - H -0.5555556 5.157299 48 -0.1077222 0.9936237722
#> 7       .    B    L - M -0.5555556 5.157299 48 -0.1077222 0.9936237722
#> 8       .    B    L - H  9.4444444 5.157299 48  1.8312771 0.1703517915
#> 9       .    B    M - H 10.0000000 5.157299 48  1.9389993 0.1388570254

^{由 reprex package (v1.0.0) 于 2021 年 3 月 15 日创建}

带有 combine = TRUE 的方法不适用于除 adjust = "none" 之外的任何东西，因为家庭规模是所有对比组合的规模。此外，Tukey 方法只能应用于单组成对比较。两组或更多组配对比较组合不构成一组配对比较，因此无法使用 Tukey 方法进行调整。

如果目标是向其他人展示结果，我仍然不建议这样做；因为查看这个数据框会让人非常不清楚 P 值调整是如何进行的，以及对哪些家庭进行了调整。在这个例子中，我们有六个比较系列； cons 的原始注释显示清楚地说明了这一点，而 df 的列表没有。

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

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

通常，当 simple 是一个列表或“每个”时，返回值是一个 emm_list 对象，每个条目对应于 简单的。但是，在 combine = TRUE 的情况下，所有元素都被合并 使用单个 emmGrid 对象中的一组对比 rbind.emmGrid.. 在这种情况下，adjust 参数设置调整 组合对比集的方法。

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

@H_386_7@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,+               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