¿Hay alguna prueba para determinar si la sobredispersión de GLM es significativa?

44

Estoy creando GLM de Poisson en R. Para verificar la dispersión excesiva, estoy mirando la relación de desviación residual a grados de libertad proporcionada por summary(model.name).

¿Existe un valor de corte o prueba para que esta relación se considere "significativa"? Sé que si es> 1, los datos se dispersan en exceso, pero si tengo relaciones relativamente cercanas a 1 [por ejemplo, una relación de 1.7 (desviación residual = 25.48, df = 15) y otra de 1.3 (rd = 324, df = 253)], ¿debería cambiar a quasipoisson / binomial negativo? Encontré aquí esta prueba de importancia: 1-pchisq (desviación residual, df), pero solo he visto eso una vez, lo que me pone nervioso. También leí (no puedo encontrar la fuente) que una relación <1.5 es generalmente segura. Opiniones?

kto
fuente

Respuestas:

45

En el paquete R AER encontrará la función dispersiontest, que implementa una Prueba de sobredispersión de Cameron y Trivedi (1990).

Sigue una idea simple: en un modelo de Poisson, la media es y la varianza es también. Son iguales. La prueba simplemente prueba esta suposición como una hipótesis nula contra una alternativa donde donde la constante significa baja dispersión y significa sobredispersión. La función Es una función monotónica (a menudo lineal o cuadrática; la primera es la predeterminada). La prueba resultante es equivalente a probar vs. y el estadístico de prueba utilizado es un estadística que es asintóticamente normal normal bajo nulo.E(Y)=μVar(Y)=μVar(Y)=μ+cf(μ)c<0c>0f(.)H0:c=0H1:c0t

Ejemplo:

R> library(AER)
R> data(RecreationDemand)
R> rd <- glm(trips ~ ., data = RecreationDemand, family = poisson)
R> dispersiontest(rd,trafo=1)

Overdispersion test

data:  rd
z = 2.4116, p-value = 0.007941
alternative hypothesis: true dispersion is greater than 0
sample estimates:
dispersion 
    5.5658 

Aquí vemos claramente que hay evidencia de sobredispersión (se estima que c es 5.57) que habla bastante en contra del supuesto de equidispersión (es decir, c = 0).

Tenga en cuenta que si no lo usa trafo=1, en realidad hará una prueba de vs. con que, por supuesto, tiene el mismo resultado que la otra prueba aparte de que la estadística de prueba se desplaza por una. Sin embargo, la razón de esto es que este último corresponde a la parametrización común en un modelo cuasi-Poisson. H0:c=1H1:c1c=c+1

Momo
fuente
1
Tuve que usar glm(trips ~ 1, data = data, family = poisson)(es decir, en 1lugar de .mis datos), pero genial, gracias
Phil
12

Una alternativa es la odTestde la psclbiblioteca que compara las proporciones de log-verosimilitud de una regresión binomial negativa con la restricción de una regresión de Poisson . Se obtiene el siguiente resultado:μ=Var

>library(pscl)

>odTest(NegBinModel) 

Likelihood ratio test of H0: Poisson, as restricted NB model:
n.b., the distribution of the test-statistic under H0 is non-standard
e.g., see help(odTest) for details/references

Critical value of test statistic at the alpha= 0.05 level: 2.7055 
Chi-Square Test Statistic =  52863.4998 p-value = < 2.2e-16

Aquí se rechaza el nulo de la restricción de Poisson a favor de mi regresión binomial negativa NegBinModel. ¿Por qué? Porque la estadística de prueba 52863.4998excede 2.7055con a p-value of < 2.2e-16.

La ventaja de esto AER dispersiontestes que el objeto devuelto de la clase "htest" es más fácil de formatear (por ejemplo, convertir a LaTeX) que el 'odTest' sin clase.

Luke Singham
fuente
5

Otra alternativa es usar la P__dispfunción del msmepaquete. La P__dispfunción se puede usar para calcular las estadísticas de dispersión Pearson y Pearson después de ajustar el modelo con o .χ2glmglm.nb

Travesura_Monkey
fuente
2

Otra opción sería usar una prueba de razón de probabilidad para demostrar que un GLM de cuasipoisson con sobredispersión es significativamente mejor que un GLM de poisson regular sin sobredispersión:

fit = glm(count ~ treatment,family="poisson",data=data) 
fit.overdisp = glm(count ~ treatment,family="quasipoisson",data=data) 
summary(fit.overdisp)$dispersion # dispersion coefficient
pchisq(summary(fit.overdisp)$dispersion * fit$df.residual, fit$df.residual, lower = F) # significance for overdispersion
Tom Wenseleers
fuente