对 lmer 和 glmer 模型使用 aictab()
using aictab() for lmer and glmer models
我正在尝试使用 aictab()
输出进行 AICc 模型选择。这个帖子是相似的() but does not apply to my problem because it used glm.nb
models and I am using lmer
models and also glmer(family=binomial)
. The aictab
documentation (https://www.rdocumentation.org/packages/AICcmodavg/versions/2.2-1/topics/aictab)说函数可以处理lm
、lme
和glm
模型:那包括lmer
和glmer
吗?
我的代码在几天前确实有效,但最近出现故障并返回此错误代码
Error in aictab.default() : Function not yet defined for this object class
m1 <- lmer(y ~ x + (1|z), data=df)
m2 <- lmer(y ~ w + (1|z), data=df)
m3 <- lmer(y ~ v + (1|z), data=df)
cand.set <- list(m1, m2, m3)
names <- list("x", "w", "v")
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
tl;dr 您加载了 lmerTest
包,因此您的模型具有不同的 class,这令人困惑 aictab()
。您可以确保在拟合模型时加载了 lme4
和 而不是 lmerTest
,或者使用 bbmle::AICctab()
,这似乎更稳健。 (可能值得就此联系 AICcmodavg
的包维护者......)
设置 lme4
:
library(lme4)
library(AICcmodavg)
ss <- transform(sleepstudy, junk1=rnorm(nrow(sleepstudy)),
junk2=rnorm(nrow(sleepstudy)))
m1 <- lmer(Reaction ~ Days + (1|Subject), data=ss, REML=FALSE)
m2 <- update(m1, . ~ . - Days + junk1)
m3 <- update(m1, . ~ . - Days + junk2)
cand.set <- list(m1, m2, m3)
names <- c("Days", "junk1", "junk2") ## this should be a vector, not a list ...
这很好用:
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
## K AICc Delta_AICc AICcWt Cum.Wt LL
## Days 4 1802.31 0.00 1 1 -897.04
## junk2 4 1918.47 116.16 0 1 -955.12
## junk1 4 1918.71 116.40 0 1 -955.24
现在加载 lmerTest
并改装模型(我们可以这样做,例如 m1 <- as(m1, "lmerModLmerTest")
但改装很容易)。
library(lmerTest)
m1 <- lmer(Reaction ~ Days + (1|Subject), data=ss, REML=FALSE)
m2 <- update(m1, . ~ . - Days + junk1)
m3 <- update(m1, . ~ . - Days + junk2)
cand.set <- list(m1, m2, m3)
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
Error in aictab.default(cand.set, modnames = names, second.ord = TRUE, :
Function not yet defined for this object class
bbmle::AICctab()
函数更健壮一点,因为它仅依赖于为 class 定义的 logLik
方法(默认情况下它给出 table只有 delta-AIC 和 df)
library(bbmle)
AICctab(m1, m2, m3, mnames=names, base=TRUE, weights=TRUE, logLik=TRUE)
## logLik AICc dLogLik dAICc df weight
## Days -897.0 1802.3 58.2 0.0 4 1
## junk2 -955.1 1918.5 0.1 116.2 4 <0.001
## junk1 -955.2 1918.7 0.0 116.4 4 <0.001
我正在尝试使用 aictab()
输出进行 AICc 模型选择。这个帖子是相似的(glm.nb
models and I am using lmer
models and also glmer(family=binomial)
. The aictab
documentation (https://www.rdocumentation.org/packages/AICcmodavg/versions/2.2-1/topics/aictab)说函数可以处理lm
、lme
和glm
模型:那包括lmer
和glmer
吗?
我的代码在几天前确实有效,但最近出现故障并返回此错误代码
Error in aictab.default() : Function not yet defined for this object class
m1 <- lmer(y ~ x + (1|z), data=df)
m2 <- lmer(y ~ w + (1|z), data=df)
m3 <- lmer(y ~ v + (1|z), data=df)
cand.set <- list(m1, m2, m3)
names <- list("x", "w", "v")
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
tl;dr 您加载了 lmerTest
包,因此您的模型具有不同的 class,这令人困惑 aictab()
。您可以确保在拟合模型时加载了 lme4
和 而不是 lmerTest
,或者使用 bbmle::AICctab()
,这似乎更稳健。 (可能值得就此联系 AICcmodavg
的包维护者......)
设置 lme4
:
library(lme4)
library(AICcmodavg)
ss <- transform(sleepstudy, junk1=rnorm(nrow(sleepstudy)),
junk2=rnorm(nrow(sleepstudy)))
m1 <- lmer(Reaction ~ Days + (1|Subject), data=ss, REML=FALSE)
m2 <- update(m1, . ~ . - Days + junk1)
m3 <- update(m1, . ~ . - Days + junk2)
cand.set <- list(m1, m2, m3)
names <- c("Days", "junk1", "junk2") ## this should be a vector, not a list ...
这很好用:
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
## K AICc Delta_AICc AICcWt Cum.Wt LL
## Days 4 1802.31 0.00 1 1 -897.04
## junk2 4 1918.47 116.16 0 1 -955.12
## junk1 4 1918.71 116.40 0 1 -955.24
现在加载 lmerTest
并改装模型(我们可以这样做,例如 m1 <- as(m1, "lmerModLmerTest")
但改装很容易)。
library(lmerTest)
m1 <- lmer(Reaction ~ Days + (1|Subject), data=ss, REML=FALSE)
m2 <- update(m1, . ~ . - Days + junk1)
m3 <- update(m1, . ~ . - Days + junk2)
cand.set <- list(m1, m2, m3)
aictab(cand.set, modnames=names, second.ord=TRUE, nobs=NULL, sort=TRUE)
Error in aictab.default(cand.set, modnames = names, second.ord = TRUE, : Function not yet defined for this object class
bbmle::AICctab()
函数更健壮一点,因为它仅依赖于为 class 定义的 logLik
方法(默认情况下它给出 table只有 delta-AIC 和 df)
library(bbmle)
AICctab(m1, m2, m3, mnames=names, base=TRUE, weights=TRUE, logLik=TRUE)
## logLik AICc dLogLik dAICc df weight
## Days -897.0 1802.3 58.2 0.0 4 1
## junk2 -955.1 1918.5 0.1 116.2 4 <0.001
## junk1 -955.2 1918.7 0.0 116.4 4 <0.001