使用R的interactions包探究nlme模型三向交互时遇错误求助

错误原因 #

interactions包的核心函数（包括sim_slopes和interact_plot）主要为lme4包的lmer/glmer模型设计，对nlme包的lme模型兼容性极差，两个错误的具体原因：

1. sim_slopes报错：函数在提取模型结果的列名时，lme的summary输出结构和它预期的（如lmer）不一致，导致关键变量tcol长度为0，触发条件判断错误。
2. interact_plot报错：函数调用了family()函数，但lme对象没有对应的方法（lme用corStruct/weights设置误差结构，而非glm类模型的family参数），因此抛出无适用方法的错误。

解决方案 #

方案1：改用lme4包构建模型（推荐） #

把lme模型换成lme4的lmer，interactions包可以完美支持：

# 加载包
library(lme4)
library(interactions)

# 构建lmer模型
mod1_lmer <- lmer(MS ~ age + sex + edu + GDST1*Time + HLI*Time + GDST1*HLI*Time + (1|ID), 
                  data = NuAge_long, na.action = na.omit)

# 运行sim_slopes
sim_slopes(model = mod1_lmer, pred = Time, modx = GDST1, mod2 = HLI, data = NuAge_long)

# 绘制交互图
interact_plot(model = mod1_lmer, pred = Time, modx = GDST1, mod2 = HLI, data = NuAge_long)

方案2：坚持使用nlme的lme模型，手动实现分析与绘图 #

如果不想更换模型包，可以用emmeans做简单斜率分析，用ggplot2手动绘制交互图：

（1）三向交互的简单斜率检验

library(emmeans)

# 取GDST1和HLI的10%、50%、90%分位数作为分组点
gdst_quant <- quantile(NuAge_long$GDST1, c(0.1, 0.5, 0.9), na.rm = TRUE)
hli_quant <- quantile(NuAge_long$HLI, c(0.1, 0.5, 0.9), na.rm = TRUE)

# 计算不同GDST1、HLI水平下Time的斜率及统计检验
emm_obj <- emmeans(mod1, ~ Time | GDST1 * HLI, 
                   at = list(GDST1 = gdst_quant, HLI = hli_quant))
slope_results <- pairs(emm_obj, reverse = TRUE)
print(slope_results)

（2）绘制三向交互图

library(ggplot2)

# 生成预测数据（控制协变量取均值/中位数）
pred_data <- expand.grid(
  Time = seq(min(NuAge_long$Time, na.rm = TRUE), max(NuAge_long$Time, na.rm = TRUE), length.out = 20),
  GDST1 = gdst_quant,
  HLI = hli_quant,
  age = mean(NuAge_long$age, na.rm = TRUE),
  sex = median(NuAge_long$sex, na.rm = TRUE),
  edu = mean(NuAge_long$edu, na.rm = TRUE)
)

# 生成预测值（level=0表示只取固定效应预测）
pred_data$MS_pred <- predict(mod1, newdata = pred_data, level = 0)

# 绘图
ggplot(pred_data, aes(x = Time, y = MS_pred, color = factor(GDST1), linetype = factor(HLI))) +
  geom_line(linewidth = 1) +
  labs(
    x = "时间",
    y = "认知能力预测值(MS)",
    color = "压力(GDST1)",
    linetype = "生活方式(HLI)"
  ) +
  scale_color_discrete(labels = c("低", "中", "高")) +
  scale_linetype_discrete(labels = c("低", "中", "高")) +
  theme_minimal()

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