如何在emmeans中为纵向研究设置有效自定义配对对比？

我开展了一项纵向研究：第-3天设置两种处理（MS only、HF only），之后每组个体在第0天至第2天被进一步分为四种处理（MS:Control、HF:Control、MS:WT、HF:WT）。两类处理完全不重叠——第-3天的处理在0-2天不存在，0-2天的处理在-3天也不存在。我想设置自定义配对对比，移除所有实际不存在的无关对比，但手动指定均值矩阵或使用by变量都无法实现。当前拟合模型及emmeans输出如下：

lm.hf <- lm(y ~ day*xtrt + sex, data=df)
emm <- emmeans(lm.hf, ~day*xtrt)
emm

day xtrt       emmean   SE  df lower.CL upper.CL
 -3  MS only     18.54 1.06 132    16.44    20.64
 0   MS only    nonEst   NA  NA       NA       NA
 2   MS only    nonEst   NA  NA       NA       NA
 -3  HF only     15.92 1.04 132    13.87    17.97
 0   HF only    nonEst   NA  NA       NA       NA
 2   HF only    nonEst   NA  NA       NA       NA
 -3  MS:Control nonEst   NA  NA       NA       NA
 0   MS:Control  13.63 1.20 132    11.26    16.00
 2   MS:Control  10.65 1.20 132     8.29    13.02
 -3  HF:Control nonEst   NA  NA       NA       NA
 0   HF:Control  14.81 1.34 132    12.16    17.46
 2   HF:Control  13.45 1.47 132    10.55    16.35
 -3  MS:WT      nonEst   NA  NA       NA       NA
 0   MS:WT        5.66 1.16 132     3.37     7.95
 2   MS:WT        9.66 1.34 132     7.00    12.32
 -3  HF:WT      nonEst   NA  NA       NA       NA
 0   HF:WT        9.49 1.34 132     6.84    12.14
 2   HF:WT       13.90 1.34 132    11.25    16.54

Results are averaged over the levels of: sex 
Results are given on the sqrt (not the response) scale. 
Confidence level used: 0.95 

若需更清晰了解处理设计，可参考附加图表。

核心问题是处理分组与时间点完全交叉但存在大量无数据组合，导致emmeans生成了无意义的nonEst均值。解决思路是先筛选有效均值子集，再针对性设置对比。

步骤1：筛选有效emmeans #

先移除nonEst行，只保留有实际估计值的处理-时间组合：

# 筛选emmean不为NA的有效均值
emm_valid <- emm[!is.na(emm@emmeans$emmean)]
# 查看筛选结果
emm_valid

这一步会得到仅包含实际存在的组合的均值集合，比如-3 MS only、-3 HF only、0 MS:Control等。

步骤2：手动设置自定义配对对比 #

基于筛选后的emm_valid，可根据研究需求构建对比矩阵。示例如下：

# 自定义对比矩阵，长度需与emm_valid行数一致
contrast_matrix <- list(
  # -3天两种处理的对比
  "MS_only_vs_HF_only" = c(1, -1, 0, 0, 0, 0, 0, 0, 0),
  # 0天MS组内Control与WT的对比
  "MS_Control_vs_MS_WT_0d" = c(0, 0, 1, 0, 0, -1, 0, 0, 0),
  # MS:Control在0天与2天的时间差异对比
  "MS_Control_0d_vs_2d" = c(0, 0, 1, 0, -1, 0, 0, 0, 0),
  # 0天HF组内Control与WT的对比
  "HF_Control_vs_HF_WT_0d" = c(0, 0, 0, 1, 0, 0, -1, 0, 0)
)

# 应用自定义对比
contrasts(emm_valid, method = contrast_matrix)

注意：对比矩阵的元素顺序需与emm_valid的均值行顺序对应，1和-1代表待对比的两组，0代表不参与该对比的组。

步骤3：按时间自动分组生成对比 #

如果需要按时间点分组做组内对比，可结合by参数拆分有效均值后自动生成配对对比：

# 按day分组，在每个时间点内做处理间配对对比
emm_by_day <- emmeans(emm_valid, ~ xtrt | day)
# 生成每组内的两两配对对比
contrast(emm_by_day, method = "pairwise")

此方法会自动在每个时间点内生成有意义的对比，避免跨时间点的无效组合。

内容的提问来源于stack exchange，提问作者ASAP Scramz