如何解决lme函数,用于二分变量的混合模型:在级别0,块1的反演中的奇异性
我正在考虑一些协变量来调整固定效应模型。关于模型的规范,这些协变量中的两个嵌套且具有固定的影响。看到下面的错误正在发生。
library(nlme)
library(lme4)
dados$VarCat=as.factor(dados$VarCat)
dados$VarX5=as.factor(dados$VarX5)
dados$VarX6=as.factor(dados$VarX6)
modelANew <- lme(log(Resp)~log(VarX1)+log(VarX2)+(VarX3)+(VarX4)+VarX5/VarX6,random = ~1|VarCat,dados,method="REML")
Error in MEEM(object,conLin,control$niterEM) :
Singularity in backsolve at level 0,block 1
变量X6是二分变量。在我看来,这似乎干扰了模型的收敛或估计。我该怎么解决?
解决方法
您的数据不平衡,导致固定效应模型的秩不足(或者,如果愿意,则为多重共线性)。当包含X5/X6
时,表示要估算X5
和X6
的所有组合的效果。但是:
with(dd,table(VarX6,VarX5))
VarX5
VarX6 A B H IND Q S T
0 2 9 94 155 0 1 15
1 0 0 0 0 8 0 0
只有VarX5=Q
处于VarX6=1
级别,而从未处于VarX6=0
级别。这意味着VarX6
变量及其与VarX5
的交互是冗余信息。
如注释中所指出,如果您在lme4::lmer()
中运行它,则会自动为您删除多余的列,并显示一条消息:
library(lme4)
m2 <- lmer(log(Resp)~log(VarX1)+log(VarX2)+(VarX3)+(VarX4)+
VarX5/VarX6 + (1|VarCat),dd,REML=TRUE)
固定效应模型矩阵的秩不足,因此会降低7列/系数
您可以通过attr(getME(m2,"X"),"col.dropped")
找出它删除了哪些列。
或者,如果您将其拟合到lm()
中(我知道您想拟合混合模型,但这是一个很好的诊断方法),您会发现它并没有抱怨,但会自动设置所有NA
的冗余系数:
m3 <- lm(log(Resp)~log(VarX1)+log(VarX2)+(VarX3)+(VarX4)+
VarX5/VarX6,data=dd)
coef(m3)
(Intercept) log(VarX1) log(VarX2) VarX3 VarX4
0.46921538 0.79476848 -0.45769296 1.85386835 -2.78321092
VarX5B VarX5H VarX5IND VarX5Q VarX5S
-0.04677216 0.21896140 0.24584351 -2.00226719 0.32677006
VarX5T VarX5A:VarX61 VarX5B:VarX61 VarX5H:VarX61 VarX5IND:VarX61
0.17474369 NA NA NA NA
VarX5Q:VarX61 VarX5S:VarX61 VarX5T:VarX61
NA NA NA
这个问题与Singularity in backsolve at level 0,block 1 in LME model非常相似。当您的设计如此不平衡时,“要做些什么”就不是一个简单答案的问题。
- 您可以自己从模型中删除条款(例如,在这种情况下,您无法真正估计关于
VarX6
的任何信息,因为它与VarX5
完全多余,因此请在其中替换VarX5/VarX6
VarX5
来建立您的模型。 - 您可以使用诸如
lmer
之类的功能来自动为您删除条款
您无法所做的实际上是估算VarX5/VarX6
-您的设计只是不包含该信息。有点像说“我想估计汽车颜色对速度的影响,但我只测量了红色汽车”。
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。