Prueba y pruebas de múltiples efectos aleatorios

8

Tengo curiosidad acerca de cómo el paquete lmerTest en R, específicamente la función "rand", maneja las pruebas de efectos aleatorios. Considere el ejemplo del pdf lmerTest en CRAN que usa el conjunto de datos "zanahorias" incorporado:

#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)

El modelo especifica dos variaciones aleatorias (la intersección y "sens2"), ambas anidadas en "Consumidor", y la covarianza entre la intersección y "sens2". A continuación se muestra la salida (no proporcionada en el pdf) para los componentes aleatorios de la ejecución de lmer:

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

Lo que se espera dada la especificación del modelo. La salida de la función rand sigue:

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

Dada la tabla de efectos aleatorios, creo que lmerTest está evaluando la pendiente aleatoria para "sens2", pero también podría ser la covarianza entre la pendiente y la intersección. La prueba para la intercepción aleatoria no está incluida. Calculé otro modelo con solo la intersección aleatoria (sin pendiente aleatoria o covarianza), y obtuve lo siguiente de la declaración "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

Aquí se proporciona la prueba de la varianza aleatoria asociada con la intercepción. Entonces, ¿alguien sabe lo que está probando la prueba de los componentes de varianza aleatoria del primer modelo? ¿Hay alguna forma que no pueda ver en la documentación para probar los tres componentes aleatorios? Debo mencionar que la página para la prueba de rand en inside-R.org tiene la siguiente descripción confusa (que no veo en el pdf sobre CRAN):

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

¿Es posible que la descripción de "Valores" lo tenga al revés y que solo se reporten efectos significativos? Ejecuté el procedimiento de "paso" y no estaba claro si los tres componentes de varianza aleatoria se consideraron en la ejecución.

Cualquier idea sobre el asunto es muy apreciada.

Joe

EDITAR: La trama se complica un poco. Se me ocurrió comprobar una estructura de covarianza "diagonal" (sin covarianza entre la intersección aleatoria y la pendiente) usando lo siguiente (basado en esta excelente publicación ):

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

La salida de lmer para las variaciones aleatorias, usando VarCorr, es la siguiente:

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

Que omite correctamente la covarianza (correlación) entre la pendiente aleatoria y la intercepción. La ejecución de la función "rand" desde lmerTest produce el siguiente resultado:

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

Por lo tanto, probará los dos componentes de varianza para este modelo. Pero la pregunta permanece con respecto al primer modelo con la covarianza aleatoria. ¿Qué es la prueba lmerTest?

BBJonz
fuente

Respuestas:

2

La documentación sobre la lmerTest::rand()función es definitivamente breve.

Por lo que he podido reunir, creo que pone a prueba la hipótesis de que la variación del efecto aleatorio (es decir, una intercepción variable (1 | Consumer) ) es significativa versus la nula de que no hay variación entre el nivel de grupo, , donde para es el indicador de grupo. (Estoy siguiendo la notación de Gelman & Hill (2007), ver cap. 12).H0:σα2=0αj[i]N(μα,σα2)j=1,,J


No soy un experto, así que el código me confundió un poco. Específicamente, no tengo claro qué hace la elimRandEffsfunción, pero supongo que está convirtiendo en un término fijo (que se agrupa) y luego comparando esto con el modelo original. Esperemos que alguien con mayor conocimiento pueda aclarar esto.αj[i]α

En el lado teórico, randdebe estar realizando algo como la prueba propuesta en Lee y Braun 2012 . Sin embargo, a diferencia de su extensión para probar efectos aleatorios a la vez (sección 3.2), la salida parece probar solo un efecto aleatorio a la vez. Un resumen más simple de la misma idea se encuentra en las diapositivas 10-12 que se encuentran aquí .0<rqrand


Por lo tanto, su intuición de que "lmerTest está evaluando la pendiente aleatoria para 'sens2' [y] también podría ser la covarianza entre la pendiente y la intercepción" es correcta en esas randpruebas si las variaciones de efecto aleatorio son significativamente diferentes de cero.

Pero es incorrecto decir que "la prueba para la intercepción aleatoria no está incluida". El RE en su primera especificación:

 (1 + sens2 | Consumer) 

supone una correlación distinta de cero entre la intersección y la pendiente, lo que significa que varían juntas y, por lo tanto, rand()prueba esa especificación con un modelo sin RE (es decir, reduciendo a la regresión clásica).

La segunda especificación

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

da dos líneas de salida porque los efectos aleatorios son aditivamente separables. Aquí randprueba (en la primera fila de salida) un modelo con una intersección agrupada / fija con una pendiente aleatoria según su especificación. En la segunda fila, la prueba es contra una pendiente agrupada con una intercepción aleatoria. Entonces, al igual que la stepfunción, randprueba los RE independientes uno por uno.

Todavía estoy desconcertado por inside-R.org, tenga en cuenta que

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

Eso no está en la documentación del paquete, así que no sé de dónde vino ni dónde se encontraría dicha prueba en la salida.

EDITAR

Creo que estoy equivocado sobre el modelo nulo en un modelo de pendiente / intersección correlacionado como en la primera especificación. La stepdocumentación dice:

en la parte aleatoria si existe correlación entre pendiente e intersección, el modelo simplificado contendrá solo una intersección. Es decir, si la parte aleatoria del modelo inicial es (1 + c | f), este modelo se compara con (1 | f) mediante el uso de LRT.

Me imagino que el mismo principio está en juego rand.

Frijoles Tony
fuente
1

La documentación de lmerTest describe rand()como ceder

"... un vector de estadísticas de Chi cuadrado y los correspondientes valores p de pruebas de razón de probabilidad"

Por lo tanto, creo que es una prueba de razón de probabilidad. Es decir, simplemente, la comparación de un modelo con un efecto aleatorio dado con ese mismo modelo sin el efecto aleatorio.

Con respecto al ejemplo, rand()compara el modelo con una pendiente aleatoria de sens2 dentro de Consumer con un modelo con la intercepción aleatoria de Consumer .

Calcule los dos modelos:

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

Examine la salida de rand(m):

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

Realice una prueba de razón de probabilidad comparando el modelo m con el modelo 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 *

De hecho, el anova() Chisq es ligeramente menor que el de rand(m), pero por lo demás, la salida es esencialmente idéntica. Quizás mi interpretación es inexacta, pero siempre supuse que así era como la rand()función generaba su salida.

JMB
fuente