使用 Tukey 方法调整时，如何计算 emmeans 的置信区间？

Question

我很难恢复到 emmeans 包给出的“CL”值。

我的第一个目标是进行一种功率研究：我正在研究植物，我想看看我可以在多大程度上改进下个季节的实验设计，使我的基因型之间有更大的差异。我的基因型被重复，这被用作阻断因子。我想研究重复次数（假设均值和残差 SE 保持不变）。最后我想在预算问题和统计能力之间找到折衷。

我正在考虑查看置信区间 (CI)，考虑到如果两种基因型的 CI 至少涵盖一个共同值，则这两种基因型没有显着差异。所以我想看看在玩重复次数时我能减少多少 CI。

我对使用 Tukey 方法进行多重比较的 CI 计算的猜测是：

μi ± (σ/(2*√(ni)))*q

with μi being the mean for the i-th genotype

ni the number of observations for the i-th genotype,

σ the residual standard error and

q the studentized range statistics with t (the number of genotypes) and dfE (the degrees of freedom of the residual) as argument for α = 0.05

相当于：

μi ± (SEi*q)/2

with SEi being the standard error of the i-th genotype

这是一个包含 6 个基因型和 3 个重复以及代码 I 运行 :

的小型示例数据集

library(emmeans)

# import DF
df <- structure(list(Geno = structure(c(6L, 3L, 2L, 1L, 4L, 5L, 2L, 
                                     4L, 5L, 1L, 3L, 6L, 2L, 3L, 1L, 5L, 6L, 4L),
                                   .Label = c("A", "B", "C", "D", "E", "F"),
                                   class = "factor"), variable = c(2.279628571, 3.925157143, 3.089, 2.26, 2.503,
                                                                   2.495114286, 2.867166667, 3.069238095, 3.884285714, 3.409595238, 
                                                                   3.710714286, 1.763142857, 2.865285714, 4.219214286, 3.263452381, 
                                                                   3.359428571, 2.335285714, 2.443), 
                 Rep = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L), 
                                 .Label = c("1", "2", "3"), class = "factor")), .Names = c("Geno", "variable", "Rep"),
            row.names = c(NA, 18L), class = "data.frame")

# Number of observation
N <- nrow(df)

# Number of treatments
t <- nlevels(df$Geno)

# Define the model
model <- lm(data=df, variable ~ Geno + Rep)

# Compute means and confidence intervals
E_means <- as.data.frame(emmeans(model, pairwise~Geno)$emmean)

# Extract Standard Errors
se <- E_means[,"SE"]

# Studentized range statistic
q <- qtukey(0.95, nmeans=t, df=E_means[,"df"])

expected_CI <- se*q/2
emmeans_CI <- ((E_means[,"upper.CL"] - E_means[,"lower.CL"])/2)

print(expected_CI)
print(emmeans_CI)

请注意，大部分时间我都在处理不平衡的数据（所以我对每种基因型都有不同的 SE）。我想首先说明一个简单（=平衡）的情况。

expected_CI 和 emmeans_CI 对于我到目前为止所做的每项测试总是不同的（差别不大，但仍然不同）所以我想我做的计算方式与emmeans 包。所以这是我的问题：它是如何在 emmeans 包中完成的？

如有任何帮助，我们将不胜感激！

Answer 1

您似乎将 EMM 与 EMM 的差异混为一谈。此外，您使用的公式仅适用于 balanced 单向设计。

以你的例子为例，如果你愿意尝试：

emm <- emmeans(model, pairwise ~ Geno)

然后做：

confint(emm$emmeans, adjust = "tukey")
confint(emm$contrasts, adjust = "tukey")

您会注意到一些事情。首先，EMM 的 CI 本身不使用指定的 Tukey 调整（它用 Sidak 代替），因为该调整仅适用于成对比较。其次，这两个摘要具有不同数量的结果（除非恰好有 3 个处理）和不同的标准误差。例如，有 4 个处理，有 4 个 EMM，但有 6 个成对比较；对于不平衡的设计，可能有 6 个不同的 SE 与这些比较相关。

如果你想了解成对比较的 Tukey 调整，你需要使用适合成对比较的统计数据——第二个总结。 Tukey 调整后的 CI 使用估计差异加减 sqrt(0.5*qtukey(.95, nmeans, df)) * SE，其中 SE 来自成对差异的结果。

使用 Tukey 方法调整时，如何计算 emmeans 的置信区间？

How are computed the confidence intervals of emmeans when adjusted with Tukey method?

r

confidence-interval

emmeans