R lmerTest 和多重随机效应检验

机器算法验证 lme4-nlme 混合模式 似然比
2022-03-31 12:14:34

我很好奇 R 中的 lmerTest 包，特别是“rand”函数，如何处理随机效应的测试。考虑使用内置“胡萝卜”数据集的 CRAN 上lmerTest pdf的示例：

#import lme4 package and lmerTest package
  library(lmerTest)
#lmer model with correlation between intercept and slopes
#in the random part
  m <- lmer(Preference ~ sens2+Homesize+(1+sens2|Consumer), data=carrots)
# table with p-values for the random effects
  rand(m)

该模型指定了两个随机方差（截距和“sens2”），都嵌套在“Consumer”中，以及截距和“sens2”之间的协方差。来自 lmer 运行的随机分量的输出（未在 pdf 中提供）如下：

Random effects:
Groups   Name        Variance Std.Dev. Corr
Consumer (Intercept) 0.195168 0.44178      
         sens2       0.002779 0.05271  0.18
Residual             1.070441 1.03462      
Number of obs: 1233, groups:  Consumer, 103

鉴于模型规范，这是预期的。rand 函数的输出如下：

Analysis of Random effects Table:
                 Chi.sq Chi.DF p.value  
sens2:Consumer   6.99      2    0.03 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

鉴于随机效应表，我认为 lmerTest 正在评估“sens2”的随机斜率，但它也可能是斜率和截距之间的协方差。不包括随机截距测试。我估计了另一个只有随机截距（没有随机斜率或协方差）的模型，并从“rand”语句中得到以下信息：

Analysis of Random effects Table:
           Chi.sq Chi.DF p.value    
Consumer   79.6      1  <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

此处提供了与截距相关的随机方差的检验。那么，有谁知道第一个模型中随机方差分量的测试是什么？有没有一种我无法从文档中看到的方法来测试所有三个随机组件？我应该提到inside-R.org的 rand 测试页面有以下令人困惑的描述（我在 CRAN 的 pdf 中没有看到）：

Values
Produces a data frame with tests for the random terms being non-significant.

Note
If the effect has random slopes, then first the correlations between itercept [sic] and slopes are checked for significance

“价值观”描述是否有可能倒退并且只报告显着影响？我运行了“步骤”程序，但不清楚运行中是否考虑了所有三个随机方差分量。

非常感谢您对此事的任何见解。

乔

编辑：情节变厚了一点。我突然想到通过使用以下方法检查“对角线”协方差结构（随机截距和斜率之间没有协方差）（基于这个优秀的帖子）：

m2 <- lmer(Preference ~ sens2+Homesize+(1|Consumer)+(0+sens2|Consumer), data=carrots)

使用 VarCorr 的随机方差的 lmer 输出如下：

Groups     Name        Std.Dev.
Consumer   (Intercept) 0.441807
Consumer.1 sens2       0.052719
Residual               1.034618

它正确地忽略了随机斜率和截距之间的协方差（相关性）。从 lmerTest 运行“rand”函数会产生以下输出：

Analysis of Random effects Table:
                 Chi.sq Chi.DF p.value    
Consumer         84.4      1  <2e-16 ***
sens2:Consumer    6.3      1    0.01 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

所以它将测试这个模型的两个方差分量。但问题仍然是关于具有随机协方差的第一个模型。什么是 lmerTest 测试？

2个回答

该lmerTest::rand()函数的文档绝对简洁。

从我能够收集到的内容来看，我认为它检验了随机效应变化（即不同的截距 (1 | Consumer) ）是显着的假设，而组级变化，其中 for是组指示符。（我遵循 Gelman & Hill (2007) 的符号，见第 12 章）。H0:σα2=0αj[i]∼N(μα,σα2)j=1,…,J

我不是专家，所以代码让我有点困惑。具体来说，我不清楚该elimRandEffs函数的作用，但我猜它正在将转换为固定（即合并）术语，然后将其与原始模型进行比较。希望有更多知识的人可以澄清这一点。αj[i]α

在理论上，rand必须执行类似于Lee 和 Braun 2012中提出的测试。然而，与他们对一次测试随机效应的扩展不同（第 3.2 节），输出似乎一次只测试一个随机效应。在此处找到的幻灯片 10-12 中对同一想法进行了更简单的总结。0<r≤qrand

rand因此，如果随机效应方差与零显着不同，您认为“lmerTest 正在评估 'sens2' [和] 的随机斜率也可能是斜率和截距之间的协方差”的直觉是正确的。

但是说“不包括随机截距的测试”是不正确的。您的第一个规范中的 RE：

 (1 + sens2 | Consumer)

假设截距和斜率之间的相关性不为零，这意味着它们一起变化，因此rand()针对无 RE 模型测试该规范（即，减少到经典回归）。

第二规范

 (1  | Consumer) + (0 + sens2 | Consumer)

给出两行输出，因为随机效应是加法可分的. 此处rand测试（在第一个输出行中）具有池化/固定截距的模型，该模型具有随机斜率，符合您的规范。在第 2 行中，测试是针对具有随机截距的合并斜率的。因此，就像step函数一样，rand一次测试一个独立的 RE。

我仍然对 inside-R.org 感到困惑，请注意

  Note
  If the effect has random slopes, then first the correlations between itercept [sic] and slopes are checked for significance

那不在包文档中，所以我不知道它来自哪里，也不知道在输出中会在哪里找到这样的测试。

编辑

我认为我对第一个规范中的相关斜率/截距模型中的空模型是错误的。step文档说：

在随机部分，如果斜率和截距之间存在相关性，则简化模型将只包含一个截距。也就是说，如果初始模型的随机部分是(1+c|f)，那么这个模型通过使用LRT与(1|f)进行比较。

我想同样的原则也适用于rand.

lmerTest文档描述rand()为屈服

“......卡方统计的向量和似然比检验的相应 p 值。”

因此，我认为这是一个似然比检验。也就是说，简单地说，将具有给定随机效应的模型与没有随机效应的同一模型进行比较。

关于示例，rand()将具有随机斜率的模型与Consumer内的sens2与具有随机截距的模型进行比较Consumer。

计算两个模型：

m <- lmer(Preference ~ sens2+Homesize+(1+sens2|Consumer), data=carrots, REML = FALSE)
m2 <- lmer(Preference ~ sens2+Homesize+(1|Consumer), data=carrots, REML = FALSE)

检查输出rand(m)：

rand(m)
Analysis of Random effects Table:
           Chi.sq Chi.DF p.value  
sens2:Consumer   6.99      2    0.03 *

执行似然比测试比较模型m和模型m2：

anova(m, m2, test = "Chi")
   Df  AIC    BIC    logLik  deviance  Chisq   Chi Df Pr(>Chisq)  
m   5 3751.4 3777.0 -1870.7   3741.4                           
m2  7 3748.7 3784.5 -1867.4   3734.7   6.6989    2     0.0351 *

实际上，anova() Chisq略小于rand(m)，但除此之外，输出基本相同。也许我的解释不准确，但我一直认为这就是rand()函数生成输出的方式。

其它你可能感兴趣的问题
上一篇通过移除方差接近于零的预测变量来证明特征选择 下一篇学习贝叶斯定理时“啊哈”时刻的资源