如何获取按测量日期拆分、跨基因型和土壤平均的emmeans？

我对leaf_area做了平方转换，拟合了包含交互项的线性混合模型：

fit <- lmer(leaf_area^2 ~genotype*soil_type + date_measurement + light + (1|repetition) + (1|y_position) + (1|x_position), data = dataset)

想获取各测量日期下，跨基因型和土壤类型平均的emmeans，执行了以下命令：

fit.emm <- emmeans(fit, ~ genotype*soil_type + date_measurement, type = "response")

但结果里date_measurement被自动平均成了27.4（实际测量日期是21、23、28、30、35 das）。即使去掉基因型和土壤类型的交互项，重新拟合模型并运行emmeans：

fit <- lmer(leaf_area^2 ~genotype + soil_type + date_measurement + light + (1|repetition) + (1|y_position) + (1|x_position), data = dataset)
fit.emm <- emmeans(fit, ~ genotype + soil_type + date_measurement, type = "response")

问题依然存在——date_measurement还是被平均了。

变量类型说明：

• date_measurement、light、x_position、y_position：数值型
• genotype、soil_type：因子型

解决方法 #

核心原因是date_measurement是数值型变量，emmeans默认会对数值型协变量取均值，而非按原始离散取值拆分。以下两种方法可以实现需求：

方法1：将date_measurement转换为因子型 #

在拟合模型前，把数值型日期转成因子，emmeans会自动按每个日期水平计算结果：

# 转换变量类型
dataset$date_measurement <- factor(dataset$date_measurement)

# 重新拟合含交互项的模型
fit <- lmer(leaf_area^2 ~genotype*soil_type + date_measurement + light + (1|repetition) + (1|y_position) + (1|x_position), data = dataset)

# 直接按日期计算跨基因型和土壤类型的平均emmeans
fit.emm <- emmeans(fit, ~ date_measurement, type = "response")

# 若需先保留基因型和土壤类型的分组结果，再汇总平均
fit.emm_full <- emmeans(fit, ~ genotype*soil_type*date_measurement, type = "response")
fit.emm_avg <- emmeans(fit.emm_full, ~ date_measurement)

方法2：用at参数指定日期取值 #

如果不想修改原始变量类型，可通过at参数指定要计算的具体日期，再对基因型和土壤类型取平均：

# 提取数据中实际存在的日期值
dates <- unique(dataset$date_measurement)

# 指定日期计算emmeans，自动跨基因型和土壤类型平均
fit.emm <- emmeans(fit, ~ date_measurement, 
                   at = list(date_measurement = dates),
                   type = "response")

关键提示 #

• 数值型协变量在emmeans中默认取均值计算边际均值，这就是你得到27.4这个平均日期结果的原因。
• 两种方法各有适用场景：方法1适合后续不需要用日期做数值分析的情况；方法2保留了日期的数值属性，适合后续还要做回归等分析的场景。

内容的提问来源于Stack Exchange，提问作者Elisa