机器算法验证 - lmer() 中的“模型无法收敛”警告 - 吾爱随笔录

lmer() 中的“模型无法收敛”警告

机器算法验证 r 混合模式 lme4-nlme

2022-02-05 02:57:18

使用以下数据集，我想看看响应（效果）是否会随着站点、季节、持续时间及其交互而发生变化。一些关于统计的在线论坛建议我继续使用线性混合效应模型，但问题是由于每个站内的重复是随机的，我几乎没有机会在连续的季节从完全相同的地点收集样本（例如，后季风 s1 的 repl-1 可能与季风不同）。它不同于临床试验（采用受试者内设计），您在不同季节重复测量同一受试者。但是，考虑到站点和季节是随机因素，我运行了以下命令并收到了警告消息：

Warning messages:
1: In checkConv(attr(opt, "derivs"), optpar,ctrl=controlpar,ctrl=controlcheckConv, 
: unable to evaluate scaled gradient
2: In checkConv(attr(opt, "derivs"), optpar,ctrl=controlpar,ctrl=controlcheckConv, 
: Model failed to converge: degenerate Hessian with 1 negative eigenvalues

谁能帮我解决这个问题？代码如下：

library(lme4)
read.table(textConnection("duration season  sites   effect
                          4d    mon s1  7305.91
                          4d    mon s2  856.297
                          4d    mon s3  649.93
                          4d    mon s1  10121.62
                          4d    mon s2  5137.85
                          4d    mon s3  3059.89
                          4d    mon s1  5384.3
                          4d    mon s2  5014.66
                          4d    mon s3  3378.15
                          4d    post    s1  6475.53
                          4d    post    s2  2923.15
                          4d    post    s3  554.05
                          4d    post    s1  7590.8
                          4d    post    s2  3888.01
                          4d    post    s3  600.07
                          4d    post    s1  6717.63
                          4d    post    s2  1542.93
                          4d    post    s3  1001.4
                          4d    pre s1  9290.84
                          4d    pre s2  2199.05
                          4d    pre s3  1149.99
                          4d    pre s1  5864.29
                          4d    pre s2  4847.92
                          4d    pre s3  4172.71
                          4d    pre s1  8419.88
                          4d    pre s2  685.18
                          4d    pre s3  4133.15
                          7d    mon s1  11129.86
                          7d    mon s2  1492.36
                          7d    mon s3  1375
                          7d    mon s1  10927.16
                          7d    mon s2  8131.14
                          7d    mon s3  9610.08
                          7d    mon s1  13732.55
                          7d    mon s2  13314.01
                          7d    mon s3  4075.65
                          7d    post    s1  11770.79
                          7d    post    s2  4254.88
                          7d    post    s3  753.2
                          7d    post    s1  11324.95
                          7d    post    s2  5133.76
                          7d    post    s3  2156.2
                          7d    post    s1  12103.76
                          7d    post    s2  3143.72
                          7d    post    s3  2603.23
                          7d    pre s1  13928.88
                          7d    pre s2  3208.28
                          7d    pre s3  8015.04
                          7d    pre s1  11851.47
                          7d    pre s2  6815.31
                          7d    pre s3  8478.77
                          7d    pre s1  13600.48
                          7d    pre s2  1219.46
                          7d    pre s3  6987.5
                          "),header=T)->dat1


m1 = lmer(effect ~ duration + (1+duration|sites) +(1+duration|season),
          data=dat1, REML=FALSE)

3个回答

“解决”您遇到的未收到有关收敛失败的警告的问题相当简单：您不使用默认的BOBYQA优化器，而是选择使用早期版本中默认使用的Nelder-Mead优化例程。1.0.x或者您安装软件包optimx，以便您可以直接使用L-BFGS-B例程或（与ver. 之前的版本nlminb相同）。例如：lme41

m1 = lmer(effect~duration+(1+duration|sites)+(1+duration|season), 
          data = dat1, REML = FALSE, 
          control = lmerControl(optimizer ="Nelder_Mead")
library(optimx)
m1 = lmer(effect~duration+(1+duration|sites)+(1+duration|season), 
          data = dat1, REML = FALSE, 
          control = lmerControl(
                           optimizer ='optimx', optCtrl=list(method='L-BFGS-B')))
m1 = lmer(effect~duration+(1+duration|sites)+(1+duration|season), 
          data = dat1, REML = FALSE, 
          control = lmerControl(
                           optimizer ='optimx', optCtrl=list(method='nlminb')))

一切正常（没有警告）。有趣的问题是：

为什么你一开始就收到这些警告和
为什么当你使用REML = TRUE你没有警告。

简而言之，1. 您收到这些警告是因为您duration将因子定义为固定效应和随机斜率sites以及season. 该模型有效地超出了估计斜率和您定义的截距之间的相关性的自由度。如果您使用稍微简单的模型，例如：

m1 = lmer(effect~duration+ (1+duration|sites) + (0+duration|season) + (1|season),
          data=dat1, REML = FALSE)

您不会遇到收敛问题。该模型将有效地估计每个不相关的随机截距和随机斜率season。

此外，2. 当您定义时，REML = FALSE您使用的是估计的最大似然而不是受限的最大似然。REML 估计试图“排除”固定效应的影响 $X$ 在开始寻找最佳随机效应方差结构之前（有关此问题的更多详细信息，请参见线程“什么是“限制最大似然”以及何时使用？ ”）。在计算上，这个过程基本上是通过将原始 LME 模型方程的两个部分相乘来完成的 $y = X\beta + Z\gamma + \epsilon$ 通过矩阵 $K$ 这样 $KX = 0$ ，即你改变了原来的 $y$ 到 $Ky$ 以及 $Z$ 到 $KZ$ . 我强烈怀疑这会影响设计矩阵的条件数 $Z$ 因此，可以帮助您摆脱最初发现自己的数字难题。

最后一点是，我不确定season从一开始就用作随机效果是否有意义。毕竟只有这么多的季节，所以你不妨把它们当作固定效果。

问题是统计的而不是技术的。实际上，我使用了随机效应模型而不是固定效应模型。我认为，没有一个因素应该被视为随机因素，因为我们需要至少 5 或 6 个水平或重复才能将一个因素视为随机效应（参见这里对于随机效应因子推荐的最小组数是多少？）。

上述数据集仅包含一式三份的样本/站点/季节，这对于随机效应模型来说是不够的。在数据集中，4 天和 7 天的持续时间属于在同一时间运行的两个独立的平行实验。因此，按持续时间（4 天和 7 天）对数据集进行吐痰，并以季节和地点为因素对每个持续时间执行 2 路方差分析，这足以在此处对效果（响应变量）进行建模。模型应该如下：

lm(day_4_effect~sites*season, data=dat1)

lm(day_7_effect~sites*season, data=dat1)

感谢 Bodo Winter ( http://www.bodowinter.com/tutorial/bw_LME_tutorial2.pdf ) 和 @usεr11852 帮助我解决了这个问题。

我在较大样本中建议的优化器更改取得了成功，但根据我的经验，此修复并不总是适用于相对较小的样本。Matuschek 等人，2017 年（Journal of Memory and Language）指出，当数据无法支持最大或接近最大模型时，模型收敛失败也可能是由于随机效应结构的过度规范造成的。在这种情况下，您需要简化随机效应结构，Matuschek et al., 2017 对此提出了一些建议。至少假设我已经正确解释了他们的工作，混合效果模型令人困惑。

其它你可能感兴趣的问题

上一篇为什么不修剪你的神经网络？下一篇饼图的问题