我想制作一个“广义线性模型的固定/随机模型”(family =“二项式”),因为我有一个数据库,其中观察来自总体并且有分组结构。然后我使用函数glmer
来自lme4
包,我也读到我可以使用 glmmPQL
库中的函数 MASS
(遥远,2006)。
当我想使用 Hausman 测试(Greene,2012)证明随机模型与固定模型的使用合理时,我的问题就出现了,我没有找到允许我执行类似于 phtest
的特定函数。测试包中的特色功能plm
.
如何证明使用随机模型的合理性?
最佳答案
这是对 plm::phtest
函数的简单调整。我已经对我实际更改的唯一代码行进行了评论。 使用风险自负,如果可能的话,根据 plm::phtest
的结果进行测试。我只是调整了代码,没有考虑它是否真的做了正确的事情!
phtest_glmer <- function (glmerMod, glmMod, ...) { ## changed function call
coef.wi <- coef(glmMod)
coef.re <- fixef(glmerMod) ## changed coef() to fixef() for glmer
vcov.wi <- vcov(glmMod)
vcov.re <- vcov(glmerMod)
names.wi <- names(coef.wi)
names.re <- names(coef.re)
coef.h <- names.re[names.re %in% names.wi]
dbeta <- coef.wi[coef.h] - coef.re[coef.h]
df <- length(dbeta)
dvcov <- vcov.re[coef.h, coef.h] - vcov.wi[coef.h, coef.h]
stat <- abs(t(dbeta) %*% as.matrix(solve(dvcov)) %*% dbeta) ## added as.matrix()
pval <- pchisq(stat, df = df, lower.tail = FALSE)
names(stat) <- "chisq"
parameter <- df
names(parameter) <- "df"
alternative <- "one model is inconsistent"
res <- list(statistic = stat, p.value = pval, parameter = parameter,
method = "Hausman Test", alternative = alternative,
data.name=deparse(getCall(glmerMod)$data)) ## changed
class(res) <- "htest"
return(res)
}
示例:
library(lme4)
gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
data = cbpp, family = binomial)
gm0 <- glm(cbind(incidence, size - incidence) ~ period + herd,
data = cbpp, family = binomial)
phtest_glmer(gm1,gm0)
## Hausman Test
## data: cbpp
## chisq = 10.2747, df = 4, p-value = 0.03605
## alternative hypothesis: one model is inconsistent
这似乎也适用于 lme
模型:
library("nlme")
fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
fm0 <- lm(distance ~ age*Subject, data = Orthodont)
phtest_glmer(fm1,fm0)
## Hausman Test
## data: Orthodont
## chisq = 0, df = 2, p-value = 1
## alternative hypothesis: one model is inconsistent
关于r - Hausman 对 lme4 中的 "glmer"的规范测试,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/23630214/