dispersión en summary.glm()

Question

He realizado un glm.nb por

glm1<-glm.nb(x~factor(group))

siendo el grupo una variable categórica y x una variable métrica. Cuando intento obtener el resumen de los resultados, obtengo resultados ligeramente diferentes, dependiendo de si utilizo summary() o summary.glm . summary(glm1) me da

    ...
Coefficients:
                    Estimate Std. Error z value Pr(>|z|)  
    (Intercept)       0.1044     0.1519   0.687   0.4921  
    factor(gruppe)2   0.1580     0.2117   0.746   0.4555  
    factor(gruppe)3   0.3531     0.2085   1.693   0.0904 .
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

    (Dispersion parameter for Negative Binomial(0.7109) family taken to be 1)

mientras que summary.glm(glm1) me da

    ...
Coefficients:
                    Estimate Std. Error t value Pr(>|t|)  
    (Intercept)       0.1044     0.1481   0.705   0.4817  
    factor(gruppe)2   0.1580     0.2065   0.765   0.4447  
    factor(gruppe)3   0.3531     0.2033   1.737   0.0835 .
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

    (Dispersion parameter for Negative Binomial(0.7109) family taken to be 0.9509067)

Entiendo el significado del parámetro de dispersión, pero no el de la línea

(Dispersion parameter for Negative Binomial(0.7109) family taken to be 0.9509067) .

En el manual dice, que sería la dispersión estimada, pero parece ser una mala estimación, ya que 0,95 no se acerca a 0,7109, ¿o la dispersión estimada es algo diferente al parámetro de dispersión estimado? Supongo que, tengo que establecer la dispersión en el summary.nb(x, dispersion=) a algo, pero no estoy seguro, si tengo que poner la dispersión a 1 (que dará el mismo resultado que summary() o si debo insertar una estimación del parámetro de dispersión, que en este caso conduce a summary.nb(glm1, dispersion=0.7109) ¿o algo más? O estoy bien con sólo usar el summary(glm1) ?

Answer 1

En primer lugar, no debe utilizar summary.glm en un objeto de la clase "negbin" . Si se observa el código de la función summary.glm justo en la parte superior verás el cómputo del dispersion . Tenga en cuenta que summary.glm sólo conoce los modelos que pueden ser ajustados por glm y, por lo tanto, se trata especialmente de las familias binomial y Poisson, en las que el parámetro de dispersión $\phi$ es se supone que sea igual a 1. Para otros modelos distintos a estos, $\phi$ se calcula a partir del objeto modelo, pero tenga en cuenta que esto se basa en la suposición de que esto es apropiado para una familia que es no binomial o Poisson. El family para el modelo ajustado por glm.nb es "Negative Binomial(theta)" . Por lo tanto, cuando se utiliza summary.glm en el modelo ajustado por glm.nb el en código

if (is.null(dispersion)) 
    dispersion <- if (object$family$family %in% c("poisson", 
        "binomial")) 
        1
    else if (df.r > 0) {
        est.disp <- TRUE
        if (any(object$weights == 0)) 
                warning("observations with zero weight not used for calculating dispersion")
            sum((object$weights * object$residuals^2)[object$weights > 
            0])/df.r
    }

la prueba de "poisson" o "binomial" falla y entonces calcula $\phi$ donde en realidad se supone que es igual a 1 por defecto para esta familia (según la definición de summary.negbin .

No hay ningún problema con esto, simplemente es más sencillo llamar al método correcto y suministrar un valor diferente para $\phi$ a través de un argumento dispersion .

En segundo lugar, usted no entiende la salida. Cuando veas

Negative Binomial(0.7109)

como he aludido anteriormente, la cifra citada entre paréntesis es $\hat{\theta}$ el parámetro de la distribución Binomial Negativa. Este valor es el estimado durante el ajuste. No es $\phi$ , el parámetro de dispersión, y por lo tanto los dos números no deben ser necesariamente iguales; son sólo dos números.

Como la dispersión calculada $\phi$ (siguiendo el código que cito más arriba) es bastante cercano a uno (~0,95), la suposición de que $\phi = 1$ utilizado para los errores estándar no es tan malo en summary.negbin . Por supuesto, puedes hacer

summary(glm1, dispersion = 0.9509)

y obtener la salida adicional que el negbin te da el método, más el valor calculado y no asumido de $\phi$ .

Answer 2

De Venables & Ripley (2002), Estadística Aplicada Moderna con S El modelo "Theta" define una distribución gamma con forma $\theta$ y escala $\frac{1}{\theta}$ Por lo tanto, la media es $1$ y la variación $\frac{1}{\theta}$ . Sea $E$ sea una variable aleatoria con esta distribución; una respuesta $Y$ se distribuye condicionalmente en $E$ como Poisson con media $\mu E$ , donde $\mu$ es una función de los predictores y coeficientes según su elección de enlace. Marginalmente, su distribución es binomial negativa, con función de masa

$$ f(y)=\frac{\Gamma(\theta +y)}{\Gamma(\theta) y!}\cdot\frac{\mu^y \theta^\theta}{(\mu+\theta)^{\theta+y}} $$

expectativa

$$\operatorname{E}Y=\mu$$

y la variación

$$\operatorname{Var} Y = \mu +\frac{\mu^2}{\theta}$$

Como señala @Momo, el parámetro de dispersión es otra cosa totalmente distinta, que dejarías variar para hacer una estimación de cuasi-verosimilitud. Para el modelo binomial negativo y el modelo de Poisson (verdadero), se fija correctamente en un valor de uno.