Aún más sobre la presunta sobredispersión en modelos de Poisson

[Esta entrada continúa el ciclo al que he dedicado esta y esta otra entradas durante los últimos días.]

Las dos entradas anteriores de la serie se resumen en que:

  • el modelo de Poisson no recoge todas las fuentes de error que pueden existir en los datos y que
  • las soluciones al uso (como, p.e., usar modelos quasi-Poisson) son puros remiendos.

Si el error en el modelo de Poisson entra (también) en el término lineal, podemos modelar ese error explícitamente. Podría haber implementado la solución INLA o Stan del problema, pero me conformaré con la lme4. Primero, generaré los datos (igual que en las entradas anteriores) y añadiré una variable categórica que identifique cada registro:

n <- 1000
sigma <- .5
x <- rep(-5:5, each = n)

x_real <- -1 + .5 * x + rnorm(length(x), 0, sigma)
y <- sapply(
  x_real, 
  function(lambda) rpois(1, exp(lambda)))

datos <- data.frame(
    x = x,
    y = y,
    id = factor(1:length(x))
)

Como se aprecia, he añadido un error normal con \sigma = .5 en el término lineal. Y ahora,

library(lme4)
modelo_glmer <- glmer(
    y ~ x + (1 | id), 
    data = datos, 
    family = "poisson")
summary(modelo_glmer)

da

Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: y ~ x + (1 | id)
   Data: datos

     AIC      BIC   logLik deviance df.resid 
 23420.8  23442.8 -11707.4  23414.8    10997 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5140 -0.4656 -0.2288  0.1937  8.3973 

Random effects:
 Groups Name        Variance Std.Dev.
 id     (Intercept) 0.2763   0.5257  
Number of obs: 11000, groups:  id, 11000

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) -0.98043    0.02124  -46.16   <2e-16 ***
x            0.48867    0.00556   87.89   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
  (Intr)
x -0.780

Como se puede ver, no solo se recuperan (aproximadamente) los coeficientes originales, sino que tenemos estimaciones bastante precisas (0.52 vs 0.5) de la desviación estándar del error del término lineal.

Y aunque pensaba que esta iba a ser la última de las entradas de las serie, creo que la voy a cerrar otro día con una adicional sobre los intervalos de predicción.

4 comentarios sobre “Aún más sobre la presunta sobredispersión en modelos de Poisson

  1. Jose Luis Cañadas Reche 22 julio, 2020 15:40

    ¿Así de simple? Un efecto aleatorio aunque no haya grupos y no se repita ningún valor de la variable categórica?.
    Pensándolo bien es justo eso, modelar la varianza de la y en vez de dejarla en el término del error.

  2. Carlos J. Gil Bellosta 22 julio, 2020 17:56

    ¿Por qué tiene que haber grupos? Mira la definición matemática del modelo. Y si lo quieres ver en términos de medidas repetidas, piensa que es un caso en el que hay medidas que se repiten una única vez por sujeto.

  3. Jose Luis Cañadas Reche 22 julio, 2020 19:25

    Si, si matemáticamente es así, pero siempre lo había pensado en hay al menos 2 observaciones para alguna de los niveles.

Comenta

Your email address will not be published.

Puedes usar estas etiquetas y atributos de HTML: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>

Este sitio usa Akismet para reducir el spam. Aprende cómo se procesan los datos de tus comentarios.