obtener grados de libertad de lmer

11

He ajustado un modelo de lmer con lo siguiente (aunque con salida impresa):

Random effects:
 Groups        Name        Std.Dev.
 day:sample (Intercept)    0.09
 sample        (Intercept) 0.42
 Residual                  0.023 

Realmente me gustaría construir un intervalo de confianza para cada efecto usando la siguiente fórmula:

(n1)s2χα/2,n12,(n1)s2χ1α/2,n12

¿Hay alguna manera de obtener convenientemente los grados de libertad?

usuario1357015
fuente
1
Creo que quieres comprobar lmerTest . Hay una serie de aproximaciones para aproximar el df en un modelo de efectos mixtos para los efectos fijos (por ejemplo , Satterthwaite , Kenward-Roger, etc.). Para su caso, me parece que usted complica demasiado su vida. Asumes que cada efecto es gaussiano. Simplemente use la desviación estándar para obtener el intervalo de confianza de su elección.
usεr11852
3
@ usεr11852 En un modelo de efectos mixtos, se supone que cada efecto es gaussiano, pero el parámetro es la varianza de la distribución gaussiana, no la media. Por lo tanto, la distribución de su estimador será muy asimétrica, y el intervalo de confianza de desviaciones estándar normales ± ~ 2 no será apropiado.
Karl Ove Hufthammer
1
@KarlOveHufthammer: Tienes razón al señalar esto; Veo lo que tú (y probablemente el OP) quieres decir. Pensé que estaba preocupado por los medios y / o la realización de los efectos aleatorios cuando mencionó los grados de libertad.
usεr11852
los grados de libertad son "problemáticos" para los modelos mixtos, consulte: stat.ethz.ch/pipermail/r-help/2006-May/094765.html y stats.stackexchange.com/questions/84268/…
Tim

Respuestas:

17

En cambio, simplemente crearía intervalos de confianza de probabilidad de perfil . Son confiables y muy fáciles de calcular con el paquete 'lme4'. Ejemplo:

> library(lme4)
> fm = lmer(Reaction ~ Days + (Days | Subject),
            data=sleepstudy)
> summary(fm)
[]
Random effects:
 Groups   Name        Variance Std.Dev. Corr
 Subject  (Intercept) 612.09   24.740       
          Days         35.07    5.922   0.07
 Residual             654.94   25.592       

Ahora puede calcular los intervalos de confianza de probabilidad de perfil con la confint()función:

> confint(fm, oldNames=FALSE)
Computing profile confidence intervals ...
                               2.5 %  97.5 %
sd_(Intercept)|Subject        14.381  37.716
cor_Days.(Intercept)|Subject  -0.482   0.685
sd_Days|Subject                3.801   8.753
sigma                         22.898  28.858
(Intercept)                  237.681 265.130
Days                           7.359  13.576

También puede usar la rutina de arranque paramétrica para calcular los intervalos de confianza. Aquí está la sintaxis R (usando el parmargumento para restringir para qué parámetros queremos intervalos de confianza):

> confint(fm, method="boot", nsim=1000, parm=1:3)
Computing bootstrap confidence intervals ...
                              2.5 % 97.5 %
sd_(Intercept)|Subject       11.886 35.390
cor_Days.(Intercept)|Subject -0.504  0.929
sd_Days|Subject               3.347  8.283

Los resultados naturalmente variarán un poco para cada ejecución. Puede aumentar nsimpara disminuir esta variación, pero esto también aumentará el tiempo que lleva estimar los intervalos de confianza.

Karl Ove Hufthammer
fuente
1
Buena respuesta (+1). También mencionaría el hecho de que también se pueden obtener CIs de arranque paramétrico en este caso. Este hilo contiene una discusión muy interesante sobre el tema.
usεr11852
@ usεr11852 Gracias por la sugerencia. Ahora he agregado un ejemplo usando el bootstrap paramétrico.
Karl Ove Hufthammer
6

Los grados de libertad para los modelos mixtos son "problemáticos". Para leer más sobre él, puede consultar el lmer, los valores p y todo lo publicado por Douglas Bates. También las preguntas frecuentes sobre modelos mixtos de r-sig resumen los motivos por los que es molesto:

  • En general, no está claro que la distribución nula de la razón calculada de sumas de cuadrados sea realmente una distribución F, para cualquier elección de grados de libertad del denominador. Si bien esto es cierto para casos especiales que corresponden a diseños experimentales clásicos (anidados, parcelas divididas, bloques aleatorios, etc.), aparentemente no es cierto para diseños más complejos (desequilibrados, GLMM, correlación temporal o espacial, etc.).
  • Para cada receta simple de grados de libertad que se ha sugerido (rastro de la matriz del sombrero, etc.) parece haber al menos un contraejemplo bastante simple donde la receta falla gravemente.
  • Otros esquemas de aproximación df que se han sugerido (Satterthwaite, Kenward-Roger, etc.) aparentemente serían bastante difíciles de implementar en lme4 / nlme,
    (...)
  • Debido a que los autores principales de lme4 no están convencidos de la utilidad del enfoque general de las pruebas con referencia a una distribución nula aproximada, y debido a la sobrecarga de cualquier otra persona que esté cavando en el código para habilitar la funcionalidad relevante (como un parche o un complemento -on), es poco probable que esta situación cambie en el futuro.

Las preguntas frecuentes también ofrecen algunas alternativas

  • use MASS :: glmmPQL (usa reglas antiguas de nlme aproximadamente equivalentes a las reglas 'internas-externas' de SAS) para GLMM, o (n) lme para LMM
  • Adivine el denominador df de las reglas estándar (para diseños estándar) y aplíquelos a las pruebas t o F
  • Ejecute el modelo en lme (si es posible) y use el denominador df que se informa allí (que sigue una simple regla 'interna-externa' que debería corresponder a la respuesta canónica para diseños simples / ortogonales), aplicada a las pruebas t o F. Para la especificación explícita de las reglas que usa lme, consulte la página 91 de Pinheiro y Bates: esta página está disponible en Google Books
  • utilizar SAS, Genstat (AS-REML), Stata?
  • Suponga el denominador infinito df (es decir, la prueba de Z / chi-cuadrado en lugar de t / F) si el número de grupos es grande (> 45? Se han planteado varias reglas generales para determinar qué tan grande es "aproximadamente infinito", incluyendo [en Angrist y Pischke's "Econometría principalmente inofensiva", 42 (en homenaje a Douglas Adams)

Pero si está interesado en los intervalos de confianza, hay mejores enfoques, por ejemplo, basados ​​en bootstrap como lo sugiere Karl Ove Hufthammer en su respuesta, o los propuestos en las Preguntas frecuentes.

Tim
fuente
"Adivina el denominador df de las reglas estándar (para diseños estándar) y aplícalos a las pruebas t o F"; Realmente me gustaría si alguien pudiera dar más detalles al respecto. Por ejemplo, para un diseño anidado (p. Ej. Pacientes versus controles, varias muestras por sujeto; con ID de sujeto como efecto aleatorio), ¿cómo obtenemos los grados de libertad para tal diseño?
Arnaud A