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.