如何解决线性混合模型中畸形拟合(Singular fit)的问题

由Kaizong Ye,Weilong Zhang撰写

假设我们有一个模型

mod <- Y ~ X*Condition + (X*Condition|subject)

#Y = logit变量
#X =连续变量
#Condition=值A和B,伪编码; 重复设计
#因此所有参与者都要同时符合这两个条件
#subject=不同主题的随机效果

summary(model)
Random effects:
 Groups  Name             Variance Std.Dev. Corr             
 subject (Intercept)      0.85052  0.9222                    
         X                0.08427  0.2903   -1.00            
         ConditionB       0.54367  0.7373   -0.37  0.37      
         X:ConditionB     0.14812  0.3849    0.26 -0.26 -0.56
Number of obs: 39401, groups:  subject, 219

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)    
(Intercept)       2.49686    0.06909   36.14  < 2e-16 ***
X                -1.03854    0.03812  -27.24  < 2e-16 ***
ConditionB       -0.19707    0.06382   -3.09  0.00202 ** 
X:ConditionB      0.22809    0.05356    4.26 2.06e-05 ***

在这里,我们观察到奇异拟合,因为截距和x随机效应之间的相关性是-1。处理该模型的一种方法是删除高阶随机效应(例如X:ConditionB),并查看在测试奇异性时是否有区别另一种是使用贝叶斯方法,例如blme软件包以避免奇异性。

什么是首选方法,为什么?

我之所以这样问,是因为使用第一个或第二个会导致不同的结果-在第一种情况下,我将删除X:ConditionB随机效应,并且无法估计X和X:ConditionB随机效应之间的相关性。另一方面,使用blme允许我保留X:ConditionB并估计给定的相关性。

通常有以下几种解决办法:

1.处理此模型的一种方法是删除高阶随机效应(高Variance ),并查看在测试奇异性时是否有所不同。

当获得奇异拟合时,这通常表明模型过度拟合-也就是说,随机效应结构太复杂而无法由数据支持,这自然导致建议删除随机效应中最复杂的部分结构(通常是随机斜率)。这种方法的好处在于,它生成了一个更加简约的模型,但并没有过度拟合。

2.另一种是使用贝叶斯方法,例如blme软件包以避免奇异性。

如果希望使用最大随机效应结构来拟合模型,并且lme4获得奇异拟合,那么在贝叶斯框架中拟合相同的模型可能很好地通过检查迹线图以及各种参数的好坏来告知lme4为什么会出现问题估计收敛。采用贝叶斯方法的优点是,这样做可能会发现原始模型的问题。 (数据不支持最大随机效应结构的原因),或者可能揭示lme4无法拟合模型的原因。

简而言之,以上两种方法都有其优点。

3.与其他线性模型一样,固定效应中的共线性可能导致奇异拟合。

那将需要通过删除条款来修改模型。但是,在lmer中,当估计随机效应方差非常接近零并且(非常宽松地)数据不足以拖动时,也可以在非常简单的模型中触发该警告(或“边界(奇异)拟合”警告)。估计远离零起始值。

两种方法的正式答案大致相似。删除估计为零的字词。但是有时候,可以忽略不计的方差是合理的,但是希望将其保留在模型中。例如,因为正在故意寻找可能存在较小差异的区间,或者可能正在进行多个类似实验,并且希望始终提取所有差异,可以通过lmerControl禁止显示这些警告,可以将其设置为不使用相关测试。例如,可以包括

control=glmerControl(check.conv.singular = .makeCC(action = "ignore",  tol = 1e-4))

这保留了默认公差(makeCC需要),但是抑制了奇异拟合测试。 (默认值为action =“ warning”,它将运行测试并发出警告)。

4.尝试重新调整预测变量x,例如,

lmerfit <- glmer(y~ I(x/10) +

                  (1 + I(x/10) | s), data = d)

并且还可以使用nlme :: lme(),即

 lmefit <- lme(y~ x, data = d,

               random = ~ x| s)

要么

 lmefit <- lme(y~ I(x/10), data = d,

               random = ~ I(x/10) | s)

并且还可能将优化算法更改为lme()中的最优化;

或者

library(GLMMadaptive)



model4 <- mixed_model(fixed = y ~ x, random = ~ 1 | group,

               data = data,

               family = binomial(link="logit"))



library(glmmTMB)



model4 <- glmmTMB(y ~ x + (1 | group),

                data = data,

                family = binomial(link = "logit"))


可下载资源

关于作者

Kaizong Ye拓端研究室(TRL)的研究员。

本文借鉴了作者最近为《R语言数据分析挖掘必知必会 》课堂做的准备。

​非常感谢您阅读本文,如需帮助请联系我们!

 
QQ在线咨询
售前咨询热线
15121130882
售后咨询热线
0571-63341498

关注有关新文章的微信公众号


永远不要错过任何见解。当新文章发表时,我们会通过微信公众号向您推送。

技术干货

最新洞察

This will close in 0 seconds