Gamma GLM Likelihood For Independent Observations

Likelihood for Independent Gamma Responses

A gamma generalized linear model specifies a positive mean for each observation and a dispersion parameter that controls variability around that mean. For independent responses, the joint likelihood is the product of their gamma densities. Writing it out carefully matters: the dispersion appears not only in the exponential term, but also in the gamma normalization and the power of the response. Those terms contribute to the dispersion score.

Use the gamma parameterization with mean $\mu_i$ and common dispersion $\phi>0$, so that the shape is $1/\phi$, the scale is $\phi\mu_i$, and $\operatorname{Var}(Y_i)=\phi\mu_i^2$. The observations $y_1,\ldots,y_n$ are assumed positive. Let $x_i$ be a fixed $p$-vector of predictors, let $\beta$ be a $p$-vector of regression coefficients, and set $\eta_i=x_i^{\mathsf T}\beta$. A link $g$ relates the mean to the linear predictor by $g(\mu_i)=\eta_i$. Write $h=g^{-1}$, so $\mu_i=h(\eta_i)$, and assume the inverse link is differentiable at each predictor considered and gives $\mu_i>0$.

Definition: Under the independent gamma GLM described above, the likelihood is the product of the individual gamma densities, with each observation having mean $\mu_i=h(x_i^{\mathsf T}\beta)$ and all observations sharing dispersion $\phi$. The coefficient and dispersion scores are the partial derivatives of the log-likelihood with respect to $\beta$ and $\phi$, respectively.

For $y_i>0$, the gamma density in this parameterization is

f(y_i;\mu_i,\phi) = \frac{y_i^{1/\phi-1}\exp\!\left(-y_i/(\phi\mu_i)\right)} {\Gamma(1/\phi)(\phi\mu_i)^{1/\phi}}, \qquad y_i>0.

Independence gives the likelihood by multiplication. Taking logarithms yields a sum, which is the convenient form for differentiation.

L(\beta,\phi) = \prod_{i=1}^n \frac{y_i^{1/\phi-1}\exp\!\left(-y_i/(\phi\mu_i)\right)} {\Gamma(1/\phi)(\phi\mu_i)^{1/\phi}},

\ell(\beta,\phi) = \sum_{i=1}^n \left[ \left(\frac{1}{\phi}-1\right)\log y_i -\frac{y_i}{\phi\mu_i} -\log\Gamma\!\left(\frac{1}{\phi}\right) -\frac{1}{\phi}\log(\phi\mu_i) \right].

Here $\Gamma$ is the gamma function. Define the digamma function by $\psi(a)=\Gamma'(a)/\Gamma(a)$ for $a>0$. All score formulas below concern interior parameter values $\phi>0$ and means in the domain of the inverse link.

Regression-Coefficient Score for a General Link

The mean contribution to one observation’s log-likelihood has a particularly simple derivative. Holding $\phi$ fixed, the terms depending on $\mu_i$ are $-y_i/(\phi\mu_i)-(1/\phi)\log\mu_i$. Differentiate them with respect to $\mu_i$, then use the chain rule through $\eta_i=x_i^{\mathsf T}\beta$.

Result: If $h$ is differentiable at every $\eta_i$, the coefficient score is $$ U_\beta(\beta,\phi) = \frac{\partial\ell}{\partial\beta} = \frac{1}{\phi} \sum_{i=1}^n x_i\,h'(\eta_i)\frac{y_i-\mu_i}{\mu_i^2}. $$ If $g$ is differentiable and $g'(\mu_i)\ne 0$, then $h'(\eta_i)=1/g'(\mu_i)$, giving the equivalent form $$ U_\beta(\beta,\phi) = \frac{1}{\phi} \sum_{i=1}^n x_i\,\frac{y_i-\mu_i}{\mu_i^2g'(\mu_i)}. $$

To verify the formula, differentiate the log-likelihood first with respect to an individual mean. The gamma normalization terms that depend on $\mu_i$ are included in this derivative.

\frac{\partial\ell}{\partial\mu_i} = \frac{y_i}{\phi\mu_i^2} -\frac{1}{\phi\mu_i} = \frac{y_i-\mu_i}{\phi\mu_i^2}.

Since $\partial\mu_i/\partial\beta=h'(\eta_i)x_i$, the chain rule gives

\frac{\partial\ell}{\partial\beta} = \sum_{i=1}^n \frac{\partial\ell}{\partial\mu_i} \frac{\partial\mu_i}{\partial\beta} = \frac{1}{\phi}\sum_{i=1}^n x_i h'(\eta_i)\frac{y_i-\mu_i}{\mu_i^2}. \qquad\square

The expression uses the derivative of the inverse link, not the derivative of the link itself. The reciprocal form follows from differentiating $g(h(\eta))=\eta$: when $g'(\mu_i)\ne0$, $g'(\mu_i)h'(\eta_i)=1$. No particular link has been assumed.

Worked Example: Coefficient Score With a Nonlinear Link

Consider one observation with $x_1=1$, $\eta_1=\beta$, and the differentiable link $g(\mu)=\mu^2$ on $\mu>0$. At $\beta=4$, the fitted mean is $\mu_1=2$. Let $y_1=3$ and $\phi=1/4$. Calculate the coefficient score.

The inverse link is $h(\eta)=\sqrt{\eta}$ for $\eta>0$, so $h'(4)=1/(2\sqrt{4})=1/4$. Substitution into the general score gives

U_\beta = \frac{1}{1/4}(1)\left(\frac14\right) \frac{3-2}{2^2} = 4\cdot\frac14\cdot\frac14 = \frac14.

As a check, $\partial\ell/\partial\mu_1=(3-2)/((1/4)2^2)=1$, and $d\mu_1/d\beta=h'(4)=1/4$. Their product is $1/4$, as above.

Dispersion Score

The dispersion derivative must account for every occurrence of $\phi$, including those inside $\log(\phi\mu_i)$ and $\Gamma(1/\phi)$. Differentiating one summand of the log-likelihood gives

Result: For $\phi>0$, the dispersion score is $$ U_\phi(\beta,\phi) = \frac{\partial\ell}{\partial\phi} = \frac{1}{\phi^2} \sum_{i=1}^n \left[ \frac{y_i}{\mu_i}-1 +\log\!\left(\frac{\phi\mu_i}{y_i}\right) +\psi\!\left(\frac{1}{\phi}\right) \right]. $$

For the derivative, treat the means as fixed when differentiating with respect to $\phi$. In a GLM, $\mu_i$ depends on $\beta$, not on $\phi$. The four terms in the summand differentiate as follows:

\begin{aligned} \frac{\partial}{\partial\phi} \left[\left(\frac{1}{\phi}-1\right)\log y_i\right] &=-\frac{\log y_i}{\phi^2},\\ \frac{\partial}{\partial\phi}\left[-\frac{y_i}{\phi\mu_i}\right] &=\frac{y_i}{\phi^2\mu_i},\\ \frac{\partial}{\partial\phi}\left[-\log\Gamma\!\left(\frac{1}{\phi}\right)\right] &=\frac{\psi(1/\phi)}{\phi^2},\\ \frac{\partial}{\partial\phi}\left[-\frac{1}{\phi}\log(\phi\mu_i)\right] &=\frac{\log(\phi\mu_i)-1}{\phi^2}. \end{aligned}

Adding these contributions and combining the two logarithms produces the stated score. In particular, omitting the derivative of the gamma normalization would remove the digamma term and would not give the derivative of this likelihood.

Worked Example: Evaluating the Dispersion Score

Suppose there are two independent observations with fitted means $\mu_1=2$, $\mu_2=3$, observed values $y_1=2$, $y_2=6$, and common dispersion $\phi=1/2$. Evaluate $U_\phi$.

Here $1/\phi=2$ and $\psi(2)=1-\gamma$, where $\gamma\approx0.577216$ is Euler’s constant, so $\psi(2)\approx0.422784$. The bracket for the first observation is

\frac{2}{2}-1+\log\!\left(\frac{(1/2)2}{2}\right)+\psi(2) = 0+\log(1/2)+0.422784 \approx -0.270363.

For the second observation the bracket is

\frac{6}{3}-1+\log\!\left(\frac{(1/2)3}{6}\right)+\psi(2) = 1+\log(1/4)+0.422784 \approx 0.036490.

Since $1/\phi^2=4$, the total score is

U_\phi = 4(-0.270363+0.036490) \approx -0.935492.

Equivalently, the sum of the two brackets is $1+2\psi(2)-\log(8)$, which gives the same value. This is the score evaluated at the specified parameters; it is not, by itself, a dispersion estimate.

Why the Scores Have Mean Zero

At the true parameter values, each score has expectation zero under the model. This is a useful check on the formulas, though it does not imply that a score computed from a particular finite sample will equal zero.

Result: For fixed $\mu_i>0$ and $\phi>0$, if $Y_i$ has the gamma distribution specified above, then the expected coefficient and dispersion score contributions for observation $i$ are both zero.

For the coefficient contribution, the mean parameterization gives $E[Y_i]=\mu_i$. Therefore

E\!\left[ \frac{1}{\phi}x_i h'(\eta_i) \frac{Y_i-\mu_i}{\mu_i^2} \right] = \frac{x_i h'(\eta_i)}{\phi\mu_i^2} \left(E[Y_i]-\mu_i\right) =0.

For the dispersion contribution, put $k=1/\phi$ and $Z_i=Y_i/(\phi\mu_i)$. Then $Z_i$ has a gamma distribution with shape $k$ and scale 1. The identity $E[\log Z_i]=\psi(k)$ follows by differentiating the gamma integral

\Gamma(k)=\int_0^\infty z^{k-1}e^{-z}\,dz.

Indeed, for $k>0$, differentiation under this integral is justified locally: on any compact interval of positive $k$-values, the differentiated integrand is bounded near zero by an integrable constant times $z^{a-1}|\log z|$ for some $a>0$, and at infinity by an integrable exponential tail times a polynomial and $\log z$. Thus

\Gamma'(k) = \int_0^\infty z^{k-1}e^{-z}\log z\,dz, \qquad E[\log Z_i]=\frac{\Gamma'(k)}{\Gamma(k)}=\psi(k).

Because $Y_i=\phi\mu_i Z_i$, it follows that $E[Y_i/\mu_i]=\phi E[Z_i]=1$ and $E[\log(\phi\mu_i/Y_i)]=-\psi(1/\phi)$. Substitution into the dispersion bracket gives

E\!\left[ \frac{Y_i}{\mu_i}-1+ \log\!\left(\frac{\phi\mu_i}{Y_i}\right)+ \psi\!\left(\frac{1}{\phi}\right) \right] = 1-1-\psi\!\left(\frac{1}{\phi}\right)+\psi\!\left(\frac{1}{\phi}\right) =0. \qquad\square

Common Pitfall and Practical Checks

A frequent error is to treat the gamma shape $1/\phi$ as if it were the dispersion $\phi$, or to differentiate only the exponential part of the density. The shape and dispersion are reciprocals in this parameterization, and the normalization depends on dispersion too. Retaining the $\log\Gamma(1/\phi)$ term is essential for the correct dispersion score.

A second error is to multiply by $g'(\mu_i)$ in the coefficient score. The chain rule uses $d\mu_i/d\eta_i=h'(\eta_i)$; if expressed using the link derivative, this is $1/g'(\mu_i)$, provided that derivative is nonzero. The two worked calculations above illustrate how to check this chain rule and how to assemble the dispersion terms.

These scores are the derivatives of the full independent-observation log-likelihood. In likelihood-based fitting, stationary interior estimates satisfy score equations, but evaluating the equations, solving them, and assessing dispersion estimates are distinct tasks. The formulas here apply when all observed responses are positive, the means are positive, and the parameters are interior; they do not assign a gamma density to a zero response.

Check Your Understanding

Use the likelihood and score formulas to answer the following questions.

Why does independence make the joint log-likelihood a sum, and which part of the gamma density contributes the digamma term to the dispersion score?
For a differentiable inverse link $h$, identify the chain-rule factor in the coefficient score. How is it written in terms of $g'$ when $g'(\mu_i)\ne0$?
If $y_i=\mu_i$, what is the contribution of observation $i$ to the coefficient score? Does that force its dispersion-score contribution to be zero?
Explain why the expected dispersion-score contribution is zero even though its value for an individual observed response need not be zero.
What condition on the responses is required for the gamma likelihood written here to apply?

Gamma GLM Likelihood For Independent Observations

What You'll Learn