¿Cómo evaluar el ajuste de un GLMM binomial equipado con lme4 (> 1.0)?

19

Tengo un GLMM con una distribución binomial y una función de enlace logit y tengo la sensación de que un aspecto importante de los datos no está bien representado en el modelo.

Para probar esto, me gustaría saber si los datos están bien descritos por una función lineal en la escala logit. Por lo tanto, me gustaría saber si los residuos se comportan bien. Sin embargo, no puedo averiguar en qué trama de residuos trazar y cómo interpretar la trama.

Tenga en cuenta que estoy usando la nueva versión de lme4 ( la versión de desarrollo de GitHub ):

packageVersion("lme4")
## [1] ‘1.1.0’

Mi pregunta es: ¿Cómo inspecciono e interpreto los residuos de un modelo mixto lineal generalizado binomial con una función de enlace logit?

Los siguientes datos representan solo el 17% de mis datos reales, pero el ajuste ya lleva alrededor de 30 segundos en mi máquina, así que lo dejo así:

require(lme4)
options(contrasts=c('contr.sum', 'contr.poly'))

dat <- read.table("http://pastebin.com/raw.php?i=vRy66Bif")
dat$V1 <- factor(dat$V1)

m1 <- glmer(true ~ distance*(consequent+direction+dist)^2 + (direction+dist|V1), dat, family = binomial)

La trama más simple ( ?plot.merMod) produce lo siguiente:

plot(m1)

ingrese la descripción de la imagen aquí

¿Esto ya me dice algo?

binomial residuals glmm lme4-nlme Henrik
fuente
1
Me podría encontrar tiempo para volver y tomar una grieta en esto, pero creo que el general, la respuesta es que es difícil de hacer mucho con los residuos de los modelos binarios. Mi principal descubrimiento tan lejos de hacer zoom sobre un poco en el terreno hay muy arriba, y añadiendo una línea suavizada (utilizando type=c("p","smooth")en plot.merMod, o en movimiento a ggplotsi desea que los intervalos de confianza) es que parece que hay un patrón pequeño, pero significativo, que se podría arreglarse adoptando una función de enlace diferente. Eso es todo hasta ahora ...
Ben Bolker
@BenBolker Gracias. ¿Y no puedes simplemente publicar esto y el enlace a freakonomics como respuesta a la pregunta? Entonces al menos obtendrías los 150 puntos.
Henrik
3
Encontré que este hilo de CV, stats.stackexchange.com/questions/63566/… , es muy útil. El puesto se explica cómo crear un gráfico de desviaciones binned en R.
Nova
@Henrik ¿Podría explicarme cómo funciona el modelo true ~ distance*(consequent+direction+dist)^2 + (direction+dist|V1)? Will la estimación dar modelo de interacción entre distance*consequent, distance*direction, distance*disty la pendiente de directiony dist que varía con V1? ¿Qué (consequent+direction+dist)^2denota el cuadrado en ?
ABC
@Henrik Ejecuté tu código y muestra el Warning message: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model failed to converge with max|grad| = 0.123941 (tol = 0.001, component 1). Por qué ?
ABC

Respuestas:

18

Respuesta corta ya que no tengo tiempo para mejorar: este es un problema desafiante; Los datos binarios casi siempre requieren algún tipo de agrupamiento o suavizado para evaluar la bondad del ajuste. Fue algo útil usar fortify.lmerMod(desde lme4, experimental) junto con, ggplot2y particularmente, geom_smooth()dibujar esencialmente la misma gráfica residual-ajustada que tiene arriba, pero con intervalos de confianza (también reduje un poco los límites y para acercar ( -5,5) región). Eso sugirió alguna variación sistemática que podría mejorarse ajustando la función de enlace. (También intenté trazar los residuos contra los otros predictores, pero no fue demasiado útil).

Intenté ajustar el modelo con todas las interacciones de 3 vías, pero no fue una gran mejora ni en la desviación ni en la forma de la curva residual suavizada.

(logístico(X))λλ

## uses (fragile) internal C calls for speed; could use plogis(),
##  qlogis() for readability and stability instead
logitpower <- function(lambda) {
    L <- list(linkfun=function(mu)
              .Call(stats:::C_logit_link,mu^(1/lambda),PACKAGE="stats"),
              linkinv=function(eta)
              .Call(stats:::C_logit_linkinv,eta,PACKAGE="stats")^lambda,
              mu.eta=function(eta) {
                  mu <-  .Call(stats:::C_logit_linkinv,eta,PACKAGE="stats")
                  mu.eta <-  .Call(stats:::C_logit_mu_eta,eta,PACKAGE="stats")
                  lambda*mu^(lambda-1)*mu.eta
              },
              valideta = function(eta) TRUE ,
              name=paste0("logit-power(",lambda,")"))
    class(L) <- "link-glm"
    L
}

λ

Ver también: http://freakonometrics.hypotheses.org/8210

Ben Bolker
fuente
3

Este es un tema muy común en los cursos de bioestadística / epidemiología, y no hay muy buenas soluciones para él, básicamente debido a la naturaleza del modelo. A menudo, la solución ha sido evitar diagnósticos detallados utilizando los residuos.

Ben ya escribió que los diagnósticos a menudo requieren binning o suavizado. La agrupación de residuos está (o estaba) disponible en el brazo del paquete R, consulte, por ejemplo, este hilo . Además, hay algunos trabajos realizados que utilizan probabilidades predichas; Una posibilidad es el diagrama de separación que se ha discutido anteriormente en este hilo . Esos podrían o no ayudar directamente en su caso, pero podrían ayudar a la interpretación.

JTT
fuente
-1

Puede usar AIC en lugar de gráficos residuales para verificar el ajuste del modelo. Comando en R: AIC (modelo1) le dará un número ... así que debe comparar esto con otro modelo (con más predictores, por ejemplo) - AIC (modelo2), que dará otro número. Compare las dos salidas y querrá el modelo con el valor AIC más bajo.

Por cierto, cosas como AIC y el índice de probabilidad de registro ya se enumeran cuando obtiene el resumen de su modelo glmer, y ambos le brindarán información útil sobre el ajuste del modelo. Desea un número negativo grande para la razón de probabilidad logarítmica para rechazar la hipótesis nula.

usuario108972
fuente
3
Esto sería más útil si OP intentara comparar modelos competitivos, pero no parece que eso sea lo que están tratando de hacer, y AIC no se puede usar para evaluar el ajuste absoluto del modelo.
Patrick Coulombe