[R] GLM Gamma Family logLik formula?

[Ben Bolker]

jseabold wrote:
> 
> Hello all,
> 
> I was wondering if someone can enlighten me as to the difference
> between the logLik in R vis-a-vis Stata for a GLM model with the gamma
> family.
> 
> Stata calculates the loglikelihood of the model as (in R notation)
> some equivalent function of
> 
> -1/scale *
> sum(Y/mu+log(mu)+(scale-1)*log(Y)+log(scale)+scale*lgamma(1/scale))
> 
> where scale (or dispersion) = 1, Y = the response variable, and mu is
> the fitted values given by the fitted model.
> 
> R seems to use a very similar approach, but scale is set equal to the
> calculated dispersion for the gamma model.  However, when I calculate
> the logLik by hand this way the answer differs slightly (about .5)
> from the logLik(glm.m1).
> 
> I haven't been able to figure out why looking at the help.  If anyone
> has any ideas, the insight would be much appreciated.
> 
> 

This is not a full answer, but the beginning of an answer about how to find
out the answer.

stats:::logLik.glm

shows that R computes the log-likelihood (a bit oddly) by
extracting the AIC component of a model, dividing by 2, and subtracting
it from the number of parameters (since AIC = -2 (logLik) + 2 p).

Gamma()$aic in turn gives the formula for the AIC of a gamma GLM:

function (y, n, mu, wt, dev) 
{
    n <- sum(wt)
    disp <- dev/n
    -2 * sum(dgamma(y, 1/disp, scale = mu * disp, log = TRUE) * 
        wt) + 2
}

and ?dgamma gives the 

 The Gamma distribution with parameters 'shape' = a and 'scale' = s
     has density

               f(x)= 1/(s^a Gamma(a)) x^(a-1) e^-(x/s)



By the way, I find it deeply confusing that Stata refers to what R 
calls the "shape" parameter (a, the thing that appears inside the 
Gamma or lgamma() function) as "shape" <http://www.stata.com/help.cgi?glm>.  
Oh well.

  good luck

    Ben Bolker

-- 
View this message in context: http://www.nabble.com/GLM-Gamma-Family-logLik-formula--tp24504030p24508346.html
Sent from the R help mailing list archive at Nabble.com.