如何检验随机效应是否显着?

机器算法验证 混合模式 lme4-nlme 随机效应模型 咕噜咕噜
2022-02-14 13:35:28

我试图了解何时使用随机效果以及何时不需要。我被告知一个经验法则是,如果你有 4 个或更多组/个人,我这样做(15 个单独的驼鹿)。其中一些驼鹿进行了 2 或 3 次试验,总共进行了 29 次试验。我想知道他们在高风险环境中的行为是否有所不同。所以,我想我会将个体设置为随机效应。然而,我现在被告知没有必要将个体作为随机效应包括在内,因为他们的反应没有太大的变化。我想不通的是如何测试在将个人设置为随机效应时是否确实存在某些问题。也许最初的问题是:我可以做哪些测试/诊断来确定个体是否是一个很好的解释变量,它是否应该是一个固定效应 - qq 图?直方图?散点图?我会在这些模式中寻找什么。

我将模型与个人一起作为随机效应运行并且没有,但随后我阅读了http://glmm.wikidot.com/faq他们指出:

不要将 lmer 模型与相应的 lm 拟合或 glmer/glm 进行比较;对数似然不相称(即,它们包含不同的附加项)

在这里,我假设这意味着您无法在具有或不具有随机效应的模型之间进行比较。但是我真的不知道我应该在它们之间进行比较。

在我的具有随机效应的模型中,我还试图查看输出以查看 RE 具有什么样的证据或意义

lmer(Velocity ~ D.CPC.min + FD.CPC + (1|ID), REML = FALSE, family = gaussian, data = tv)

Linear mixed model fit by maximum likelihood 
Formula: Velocity ~ D.CPC.min + FD.CPC + (1 | ID) 
   Data: tv 
    AIC    BIC logLik deviance REMLdev
 -13.92 -7.087  11.96   -23.92   15.39
Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.00000  0.00000 
 Residual             0.02566  0.16019 
Number of obs: 29, groups: ID, 15

Fixed effects:
              Estimate Std. Error t value
(Intercept)  3.287e-01  5.070e-02   6.483
D.CPC.min   -1.539e-03  3.546e-04  -4.341
FD.CPC       1.153e-04  1.789e-05   6.446

Correlation of Fixed Effects:
          (Intr) D.CPC.
D.CPC.min -0.010       
FD.CPC    -0.724 -0.437

您会看到我作为随机效应的个体 ID 的方差和 SD = 0。这怎么可能?0 是什么意思?是对的吗?那么我的朋友说“因为没有使用ID的变化,因为随机效应是不必要的”是正确的吗?那么,那我会用它作为固定效果吗?但是,变化如此之少的事实是否意味着它无论如何也不会告诉我们太多信息?

2个回答

估计值ID' 的方差 = 0,表明组间变异水平不足以保证在模型中加入随机效应;IE。你的模型是退化的。

当您正确识别自己时:很可能,是的;ID因为随机效应是不必要的。几乎没有什么东西可以用来测试这个假设:

  1. 您可以比较(REML = F始终使用)具有和不具有随机效应的模型之间的 AIC(或您最喜欢的 IC),并了解情况如何。
  2. 您将查看anova()这两个模型的输出。
  3. 您可以使用原始模型定义的后验分布进行参数引导。

请注意,选择 1 和 2 有一个问题:您正在检查参数空间边界上的某些东西,因此实际上它们在技术上并不合理。话虽如此,我认为您不会从它们那里得到错误的见解,而且很多人都在使用它们(例如,lme4 的开发人员之一 Douglas Bates 在他的书中使用它们,但清楚地说明了这个关于正在测试的参数值的警告在可能值集的边界上)。选择 3 是 3 中最乏味的,但实际上可以让您真正了解正在发生的事情。有些人也很想使用非参数引导程序,但我认为鉴于您从一开始就做出参数假设,您不妨使用它们。

我不确定我要建议的方法是否合理,所以如果我错了,那些对此主题了解更多的人会纠正我。

我的建议是在您的数据中创建一个额外的列,其常量值为 1:

IDconst <- factor(rep(1, each = length(tv$Velocity)))

然后,您可以创建一个使用此列作为随机效应的模型:

fm1 <- lmer(Velocity ~ D.CPC.min + FD.CPC + (1|IDconst), 
  REML = FALSE, family = gaussian, data = tv)

此时,您可以将(AIC)您的原始模型与随机效应ID(我们称之为fm0)与不考虑的新模型进行比较,ID因为IDconst对于您的所有数据都是相同的。

anova(fm0,fm1)

更新

user11852 要求举个例子,因为在他/她看来,上述方法甚至不会执行。相反,我可以证明该方法有效(至少在lme4_0.999999-0我目前使用的情况下)。

set.seed(101)
dataset <- expand.grid(id = factor(seq_len(10)), fac1 = factor(c("A", "B"),
  levels = c("A", "B")), trial = seq_len(10))
dataset$value <- rnorm(nrow(dataset), sd = 0.5) +
      with(dataset, rnorm(length(levels(id)), sd = 0.5)[id] +
      ifelse(fac1 == "B", 1.0, 0)) + rnorm(1,.5)
    dataset$idconst <- factor(rep(1, each = length(dataset$value)))

library(lme4)
fm0 <- lmer(value~fac1+(1|id), data = dataset)
fm1 <- lmer(value~fac1+(1|idconst), data = dataset)

anova(fm1,fm0)

输出:

  Data: dataset
  Models:
  fm1: value ~ fac1 + (1 | idconst)
  fm0: value ~ fac1 + (1 | id)

      Df    AIC    BIC  logLik  Chisq Chi Df Pr(>Chisq)
  fm1  4 370.72 383.92 -181.36                      
  fm0  4 309.79 322.98 -150.89 60.936      0  < 2.2e-16 ***

根据最后的测试,我们应该保持随机效应,因为fm0模型具有最低的 AIC 和 BIC。

更新 2

顺便说一句,NW Galwey 在第 213-214 页的“混合建模简介:超越回归和方差分析”中提出了同样的方法。