¿Cómo puedo probar si un efecto aleatorio es significativo?

34

Estoy tratando de entender cuándo usar un efecto aleatorio y cuándo es innecesario. Me han dicho que una regla de oro es si tienes 4 o más grupos / individuos que yo hago (15 alces individuales). Algunos de esos alces fueron experimentados 2 o 3 veces para un total de 29 ensayos. Quiero saber si se comportan de manera diferente cuando están en paisajes de mayor riesgo que no. Entonces, pensé que establecería al individuo como un efecto aleatorio. Sin embargo, ahora me dicen que no hay necesidad de incluir al individuo como un efecto aleatorio porque no hay mucha variación en su respuesta. Lo que no puedo entender es cómo probar si realmente hay algo que se tiene en cuenta al configurar el individuo como un efecto aleatorio. Tal vez una pregunta inicial es: ¿Qué prueba / diagnóstico puedo hacer para determinar si Individual es una buena variable explicativa y debería ser un efecto fijo - qq parcelas? histogramas? ¿gráfico de dispersión? ¿Y qué buscaría en esos patrones?

Ejecuté el modelo con el individuo como un efecto aleatorio y sin él, pero luego leí http://glmm.wikidot.com/faq donde dicen:

no compare los modelos lmer con los ajustes lm correspondientes, ni glmer / glm; las probabilidades de registro no son proporcionales (es decir, incluyen diferentes términos aditivos)

Y aquí supongo que esto significa que no se puede comparar entre un modelo con efecto aleatorio o sin él. Pero de todos modos no sabría realmente qué debo comparar entre ellos.

En mi modelo con efecto aleatorio, también estaba tratando de ver la salida para ver qué tipo de evidencia o importancia tiene el 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

Usted ve que mi varianza y SD de la ID individual como efecto aleatorio = 0. ¿Cómo es eso posible? ¿Qué significa 0? ¿Está bien? Entonces, ¿mi amigo que dijo "ya que no hay variación usando ID como efecto aleatorio es innecesario" es correcto? Entonces, ¿lo usaría como un efecto fijo? ¿Pero el hecho de que haya tan poca variación no significa que no nos va a decir mucho de todos modos?

Kerry
fuente
Para obtener una variación exacta de 0 de un efecto aleatorio, consulte stats.stackexchange.com/questions/115090 .
ameba dice Reinstate Monica

Respuestas:

21

La estimación, IDla varianza = 0, indica que el nivel de variabilidad entre grupos no es suficiente para garantizar la incorporación de efectos aleatorios en el modelo; es decir. Su modelo es degenerado.

Como te identificas correctamente: muy probablemente sí; IDcomo efecto aleatorio es innecesario Se te ocurren pocas cosas para probar esta suposición:

  1. Puede comparar (usando REML = Fsiempre) el AIC (o su IC favorito en general) entre un modelo con y sin efectos aleatorios y ver cómo funciona.
  2. Observaría la anova()salida de los dos modelos.
  3. Puede hacer un arranque paramétrico utilizando la distribución posterior definida por su modelo original.

Tenga en cuenta que las opciones 1 y 2 tienen un problema: está buscando algo que esté en los límites del espacio de parámetros, por lo que en realidad no son técnicamente sólidos. Habiendo dicho eso, no creo que obtengas ideas erróneas de ellos y mucha gente los use (por ejemplo, Douglas Bates, uno de los desarrolladores de lme4, los usa en su libro, pero establece claramente esta advertencia sobre los valores de los parámetros que se están probando en el límite del conjunto de valores posibles). La opción 3 es la más tediosa de las 3, pero en realidad te da la mejor idea de lo que está sucediendo. Algunas personas sienten la tentación de usar bootstrap no paramétrico también, pero creo que dado el hecho de que está haciendo suposiciones paramétricas para comenzar, también podría usarlas.

usεr11852 dice Reinstate Monic
fuente
66
El paquete RLRsim es una forma realmente conveniente de probar efectos aleatorios usando pruebas de razón de probabilidad basadas en simulación.
atrichornis
@atrichornis: +1. Paquete interesante; No lo sabía. Acabo de echar un vistazo a su código, bastante sencillo, podría decir. Desearía que lo incorporen (o algo así) lme4especialmente ahora que mcmcsamp()está roto y las personas solo tienen sus propias implementaciones de arranque ad-hoc para obtener algunos valores p decentes, etc.
usεr11852 dice Reinstate Monic el
Es cierto, los modelos mixtos no son sencillos en R. Muchas aproximaciones y soluciones alternativas ... ¿Aunque reúno SAS, etc., solo paso por alto algunas de las mismas incertidumbres? Ben Bolker es coautor de ambos paquetes, puede tener sus razones para no incorporarlo. ¡Probablemente el tiempo!
atrichornis
44
La rutina de carga en el límite del espacio de parámetros tiene su propio conjunto de problemas y problemas que conducen a inconsistencias . El bootstrap no es una panacea, y no debe arrojarse a la bolsa ligeramente, suponiendo que lo resuelva todo.
StasK
2
Eche un vistazo, el argumento es muy sutil. Por lo que puedo recordar, se reduce al hecho de que está haciendo el arranque desde una distribución que es diferente de la nula; y dadas las distribuciones no estándar obtenidas en el límite, se violan las condiciones de regularidad y la distribución de arranque no converge al objetivo. Creo que la rutina de arranque no paramétrica todavía se puede construir aquí eliminando los medios grupales de los residuos. Sin embargo, con la violación de la independencia de las observaciones entre grupos, puede surgir otra capa de complicaciones.
StasK
3

No estoy seguro de que el enfoque que voy a sugerir sea razonable, por lo que aquellos que saben más sobre este tema me corrigen si me equivoco.

Mi propuesta es crear una columna adicional en sus datos que tenga un valor constante de 1:

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

Luego, puede crear un modelo que use esta columna como su efecto aleatorio:

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

En este punto, podría comparar (AIC) su modelo original con el efecto aleatorio ID(llamémoslo fm0) con el nuevo modelo que no tiene en cuenta, IDya que IDconstes el mismo para todos sus datos.

anova(fm0,fm1)

Actualizar

user11852 estaba pidiendo un ejemplo, porque en su opinión, el enfoque anterior ni siquiera se ejecutará. Por el contrario, puedo demostrar que el enfoque funciona (al menos con el lme4_0.999999-0que estoy usando actualmente).

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)

Salida:

  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 ***

Según esta última prueba, debemos mantener el efecto aleatorio ya que el fm0modelo tiene el AIC más bajo y el BIC.

Actualización 2

Por cierto, este mismo enfoque es propuesto por NW Galwey en 'Introducción al modelado mixto: más allá de la regresión y el análisis de varianza' en las páginas 213-214.

VLC
fuente
¿Has probado tu idea? Por favor, demuéstrame que estoy equivocado, pero creo que tu idea ni siquiera se ejecutará. Si IDconstes igual para todos sus datos, entonces no tiene ninguna agrupación. Necesita un factor de agrupación para tener al menos un nivel muestreado y la forma en que configura el modelo no tiene ninguno. Tal vez podría creer la razón de usar una "agrupación aleatoria", pero ese es un juego de pelota diferente en conjunto. Pruebe su enfoque con algunos datos ficticios. Creo firmemente que con su configuración propuesta lmer()no se ejecutará. (Yo uso lme4_0.99999911-1)
usεr11852 dice Reinstate Monic
@ user11852 Por favor, vea mi actualización y háganos saber si este enfoque también funciona lme4_0.99999911-1.
VLC
No, no funcionará. Y dije, " no debería " porque conceptualmente no tienes un modelo mixto para empezar. (Lo siento, podría sonar agresivo en mi comentario anterior). Supongo que lo que quieres intentar es configurar elZmatriz para ser un vector unitario, pero a menos que encuentre una manera de hacerlo explícitamente en R (es decir, escriba su propia función de desviación) no tiene suerte.
usεr11852 dice Reinstate Monic el
3
Sí, hice lo que me sugieres; No funcionará / computará. Error in lFormula(formula = value ~ fac1 + (1 | idconst), data = dataset) : grouping factors must have at least 1 sampled level. Y como dije, conceptualmente está mal. No se trata de engañar al software para que entregue algunos números, se trata de si lo que usted dice es razonable. No tiene un segundo modelo mixto para comparar si en ese modelo el efecto aleatorio es, por construcción, una constante. También podría excluirlo e intentar un modelo lineal en su lugar.
usεr11852 dice Reinstate Monic el
1
Actualizar concertación definiendo una variable aleatoria de un solo grupo en lme4. Esto se puede hacer si se establece la opción: control=lmerControl(check.nlev.gtr.1="ignore"). Ben Bolker lo menciona aquí: github.com/lme4/lme4/issues/411 .
Robin Beaumont
1

Me gustaría responder la pregunta más "inicial".

Si sospecha que existe algún tipo de heterogeneidad en la variación entre cualquiera de las variables dependientes debido a algunos factores, debe continuar y graficar los datos utilizando diagramas de dispersión y de caja. Algunos patrones comunes para verificar, pongo esta lista a continuación de varias fuentes en la web.

Patrones de heterocedasticidad

Además, trace su variable dependiente por grupos de factor / tratamiento para ver si hay una varianza constante. De lo contrario, es posible que desee explorar efectos aleatorios o regresiones ponderadas. Por ej. Este cuadro a continuación es un ejemplo de una variación en forma de embudo en mis grupos de tratamiento. Así que elijo ir a los efectos aleatorios y probar los efectos en pendientes e intersecciones.

Diagrama de caja para verificar la heterocedasticidad

A partir de aquí, las respuestas anteriores responden a su pregunta principal. También hay pruebas que verifican la heterocedasticidad, una de ellas está aquí: https://dergipark.org.tr/download/article-file/94971 . Pero no estoy seguro de si existen pruebas para detectar la heterocedasticidad a nivel de grupo.

Arun
fuente
Utilice solo el campo "Su respuesta" para proporcionar respuestas a la pregunta del OP. CV es un sitio estricto de preguntas y respuestas, no un foro de discusión. La última parte en negrita de tu publicación es una nueva pregunta, no una respuesta a esta pregunta. Si tiene una nueva pregunta, haga clic en el gris ASK QUESTIONen la parte superior y pregúntela allí. Como eres nuevo aquí, es posible que quieras hacer nuestro recorrido , que tiene información para nuevos usuarios.
gung - Restablece a Monica