Uso de R para GLM con distribución Gamma

Question

Actualmente tengo un problema para entender la sintaxis de R para ajustar un GLM utilizando la distribución Gamma.

Tengo un conjunto de datos, donde cada fila contiene 3 co-variables ( $X_1, X_2, X_3$ ), una variable de respuesta ( $Y$ ), y un parámetro de forma ( $K$ ). Quiero modelar la escala de la distribución Gamma como una función lineal de las 3 covariables, pero no entiendo cómo establecer la forma de la distribución a $K$ para cada fila de datos.

Una situación que me parece análoga es que para una distribución binomial, el MLG requiere que el número de ensayos ( $N$ ) se conoce para cada entrada de datos.

Answer 1

Utilicé el forma.gamma función de MASA descrito por Balajari (2013) para estimar después el parámetro de forma y ajustar las estimaciones de los coeficientes y las predicciones en el MLG. Te aconsejo que leas la conferencia ya que, en mi opinión, es muy clara e interesante en lo que respecta al uso de la distribución gamma en los MLG.

    glmGamma <- glm(response ~ x1, family = Gamma(link = 
                    "identity")
    library(MASS)
    myshape <- gamma.shape(glmGamma)
    gampred <- predict(glmGamma , type = "response", se = TRUE, 
                      dispersion = 1/myshape$alpha) 
    summary(glmGamma, dispersion = 1/myshape$alpha)

Answer 2

El MLG gamma habitual contiene el supuesto de que el parámetro de forma es constante, del mismo modo que el modelo lineal normal asume una varianza constante.

En lenguaje GLM, el parámetro de dispersión, $\phi$ en $\text{Var}(Y_i)=\phi\text{V}(\mu_i)$ suele ser constante.

En términos más generales, tiene $a(\phi)$ pero eso no ayuda.

Tal vez sería posible utilizar un GLM Gamma ponderado para incorporar este efecto de un parámetro de forma especificado, pero aún no he investigado esta posibilidad (si funciona es probablemente la forma más fácil de hacerlo, pero no estoy nada seguro de que lo haga).

Si tuvieras un GLM doble podrías estimación ese parámetro en función de las covariables... y si el software de doble glm te permitiera especificar un offset en el término de varianza podrías hacerlo. Parece que la función dglm en el paquete dglm permiten especificar un desplazamiento. No sé si le permitirá especificar un modelo de varianza como (digamos) ~ offset(<something>) + 0 sin embargo.

Otra alternativa sería maximizar la probabilidad directamente.

---

> y <- rgamma(100,10,.1)

> summary(glm(y~1,family=Gamma))

Call:
glm(formula = y ~ 1, family = Gamma)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-0.93768  -0.25371  -0.05188   0.16078   0.81347  

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 0.0103660  0.0003486   29.74   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for Gamma family taken to be 0.1130783) 

    Null deviance: 11.223  on 99  degrees of freedom
Residual deviance: 11.223  on 99  degrees of freedom
AIC: 973.56

Number of Fisher Scoring iterations: 5

La línea donde dice:

   (Dispersion parameter for Gamma family taken to be 0.1130783)

es el que quieres.

Que $\hat\phi$ está relacionado con el parámetro de forma de la Gamma.