如何判断lme4/nlme中线性混合模型是否过度确定？附案例疑问

模型代码与警告信息 #

library(nlme)
library(lme4)

m1 <- lmer( distance ~ age + (age|Subject), data = Orthodont )
m2 <- lmer( distance ~ age + I(age^2) + (age|Subject), data = Orthodont )
m3 <- lmer( distance ~ age + I(age^2) + I(age^3) + (age|Subject), data = Orthodont )

m1nlme <- lme(distance ~ age, random = ~ age|Subject, data = Orthodont)
m2nlme <- lme(distance ~ age + I(age^2), random = ~ age|Subject, data = Orthodont)
m3nlme <- lme(distance ~ age + I(age^2) + I(age^3), random = ~ age|Subject, data = Orthodont)
m4nlme <- lme(distance ~ age + I(age^2) + I(age^3), random = ~ age + I(age^2) + I(age^3)|Subject, data = Orthodont)

仅m3触发收敛警告：

In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,:Model failed to converge with max|grad| = 0.00762984 (tol = 0.002, component 1)

1. 收敛警告的含义与处理建议 #

• 这个警告表示模型的对数似然优化过程未达到预设的收敛标准：目标函数对参数的梯度最大绝对值（0.0076）超过了容忍阈值（0.002）。梯度趋近于0才代表找到似然函数的极值点，未收敛意味着当前得到的参数估计可能不稳定、不可靠。
• 不能直接忽略该警告，建议尝试以下操作：
  • 检查多项式项的共线性：三次项I(age^3)和低次项（age、age²）高度相关，会干扰参数估计的稳定性，可通过计算方差膨胀因子（VIF）验证。
  • 更换优化器或调整控制参数：比如给lmer加上control = lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 1e5))，增加迭代次数或换用更稳健的优化算法。
  • 对比m3nlme的结果：nlme包的lme和lme4的lmer用的优化算法不同，lme的结果可能更稳定，可用来交叉验证。

2. m2模型未报错的原因：参数计数错误与混合模型的估计逻辑 #

你的参数计数存在错误，且对混合模型的参数估计逻辑理解有误：

• 正确的参数计数：
  • 固定效应：截距、age系数、age²系数，共3个参数；
  • 随机效应协方差结构：(age|Subject)代表随机截距+随机age斜率的二元正态分布，需要估计3个参数（截距方差、age斜率方差、截距与斜率的协方差）；
  • 残差方差：1个参数；
  • 总计3+3+1=7个参数。
• 为什么m2没报错：
  混合模型的待估参数是群体层面的固定效应和随机效应分布参数，不是每个受试者的个体随机效应值（BLUPs）。BLUPs是基于估计出的协方差结构计算得到的，不属于待估参数范畴。
  虽然每个受试者只有4个观测值，但总样本量是27×4=108，足够支撑7个参数的估计。而m3因为加入了三次项，固定效应项间共线性大幅增强，加上随机结构的参数，导致优化过程难以收敛，才触发了警告。

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