基于nlme的分段增长曲线模型自相关问题处理及疑问

1. 判断lag2自相关是否构成问题 #

nlme的ACF()绘图自带模拟置信带（图中虚线），如果lag2的自相关点明显超出这个范围，就说明存在有统计学意义的自相关。另外可以手动验证：

• 提取模型的标准化条件残差（混合模型不能用普通残差），按受试者分组计算组内自相关，再取平均值：

  library(dplyr)
  res_df <- mydata %>%
    mutate(resid = resid(model, type = "pearson")) %>%
    group_by(subject) %>%
    summarise(acf_list = list(acf(resid, plot = F, lag.max = 3)$acf[-1])) %>%
    unnest(acf_list) %>%
    group_by(lag = row_number()) %>%
    summarise(mean_acf = mean(acf_list))
  

• 计算独立残差的95%置信区间：2/sqrt(7)（约0.755），如果lag2的平均自相关绝对值超过这个值，说明自相关确实存在。
• 更严谨的方法是做似然比检验：对比原模型和加入lag2自相关结构的模型，看拟合度是否显著提升。

2. 处理lag2自相关的方法 #

如果corAR1只能解决lag1问题，试试以下两种方案：

• AR(2)自相关结构：直接指定corARMA(p=2)，同时处理lag1和lag2的自相关：

  model_ar2 <- lme(score ~ time1 + time2 + time1 * group,
                   random = ~ 1 + time1 + time2 + time1 * group | subject,
                   data = mydata, 
                   correlation = corARMA(form = ~ timenum | subject, p=2),
                   control = lmeControl(opt = "optim", maxIter = 100, msMaxIter = 100))
  

• 无结构协方差矩阵：用corSymm拟合任意滞后的自相关，虽然参数多，但对小样本（每个受试者7个观测）可以尝试：

  model_unstr <- lme(score ~ time1 + time2 + time1 * group,
                     random = ~ 1 + time1 + time2 + time1 * group | subject,
                     data = mydata, 
                     correlation = corSymm(form = ~ timenum | subject),
                     control = lmeControl(opt = "optim", maxIter = 100, msMaxIter = 100))
  

• 用anova(model, model_ar2, model_unstr)比较模型的AIC/BIC，选择拟合最优的那个。

3. 未处理自相关的实际影响 #

• 核心问题是标准误估计偏倚：正自相关会低估标准误，增加I类错误（假阳性）；负自相关会高估标准误，增加II类错误（假阴性）。
• 你关注的time1*group效应本身不显著，且用两组t4-t1差值对比结果一致，说明这个效应的结论是稳健的。如果是负自相关导致标准误高估，理论上会增加II类错误，但你的交叉验证结果也没发现显著效应，所以这个风险很低。
• 但未处理自相关会影响置信区间的准确性，如果你需要报告效应的置信区间，还是建议修正自相关结构，保证推断的可靠性。

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