Likelihood and Quasi-Likelihood for Tweedie Models

Why distinguish likelihood from quasi-likelihood?

A Tweedie regression specifies more than a pattern for the mean: when a full Tweedie distribution is assumed, it gives a probability model for each response. Quasi-likelihood uses a narrower set of information, namely the mean and variance functions. This distinction matters because the two approaches can yield the same regression-coefficient estimating equation while supporting different conclusions about dispersion, the power index, and likelihood-based comparisons.

Write the independent responses as $Y_1,\ldots,Y_n$, with $E[Y_i]=\mu_i$ and $\operatorname{Var}(Y_i)=\phi\mu_i^p$. Here $\phi>0$ is the dispersion parameter and $p$ is fixed for now. In a log-link regression, $\log\mu_i=\mathbf{x}_i^\mathsf{T}\boldsymbol{\beta}$, where $\mathbf{x}_i$ is the covariate vector for observation $i$. The same reasoning below applies to a differentiable link.

Two ways to use a Tweedie model

Definition: For a specified variance function $V(\mu)>0$, a quasi-likelihood $Q(y;\mu)$ is a function of the observed value $y$ and proposed mean $\mu$ satisfying $$ \frac{\partial Q(y;\mu)}{\partial\mu}=\frac{y-\mu}{V(\mu)}. $$ For a Tweedie variance function, $V(\mu)=\mu^p$. The scaled quasi-likelihood contribution is $Q(y;\mu)/\phi$. A quasi-likelihood is not, just by this definition, a log probability density or mass function.

A full Tweedie model, by contrast, specifies a normalized distribution for each response. In exponential-dispersion form its density or mass function can be written

f(y;\theta,\phi)=a(y,\phi)\exp\left\{\frac{y\theta-b(\theta)}{\phi}\right\},

with $b'(\theta)=\mu$ and $b''(\theta)=V(\mu)$. For a fixed power $p$, the Tweedie family supplies the corresponding carrier $a(y,\phi)$ as well as the mean–variance relationship. The carrier is essential to the normalized likelihood, even though it does not depend on the regression coefficients when $p$ and $\phi$ are fixed.

Worked Example: Computing a Log-Link Quasi-Score

Take $p=3/2$, $\phi=0.4$, and two independent observations with $y_1=0$ and $y_2=6$. Use an intercept-only log-link regression with known offsets $\log(1)$ and $\log(4)$. At $\beta=0$, the fitted means are $\mu_1=1$ and $\mu_2=4$.

For the log link, $\partial\mu_i/\partial\beta=\mu_i$. The quasi-score contribution is therefore $$ \frac{\partial\mu_i}{\partial\beta}\frac{y_i-\mu_i}{\phi\mu_i^p} =\frac{y_i-\mu_i}{\phi\mu_i^{p-1}}. $$

The first contribution is $(0-1)/(0.4\cdot1^{1/2})=-2.5$. The second is $(6-4)/(0.4\cdot4^{1/2})=2/(0.8)=2.5$. Their sum is zero, so $\beta=0$ solves the quasi-score equation for these data.

For a full Tweedie model with this same fixed $p$ and $\phi$, the regression-coefficient likelihood score is the same equation, as proved below. This equality does not say that the quasi-likelihood is the full log-likelihood.

When the coefficient equations agree

Theorem: Suppose the responses are independent and follow a regular exponential-dispersion model with density or mass function $f(y;\theta,\phi)=a(y,\phi)\exp\{(y\theta-b(\theta))/\phi\}$. Assume the support does not depend on the regression coefficients, and that differentiation of the log-likelihood is valid. For fixed $\phi$ and fixed variance function $V(\mu)=b''(\theta)$, the likelihood score for $\boldsymbol{\beta}$ equals the quasi-score $$ \mathbf{U}(\boldsymbol{\beta})=\sum_{i=1}^n \frac{1}{\phi}\frac{\partial\mu_i}{\partial\boldsymbol{\beta}} \frac{y_i-\mu_i}{V(\mu_i)}. $$

Proof. The log-likelihood is $\ell=\sum_i\log f(y_i;\theta_i,\phi)$. Since the carrier has no dependence on $\boldsymbol{\beta}$, differentiating the $i$th term gives

\frac{\partial\ell_i}{\partial\boldsymbol{\beta}} =\frac{y_i-b'(\theta_i)}{\phi} \frac{\partial\theta_i}{\partial\boldsymbol{\beta}} =\frac{y_i-\mu_i}{\phi} \frac{\partial\theta_i}{\partial\mu_i} \frac{\partial\mu_i}{\partial\boldsymbol{\beta}}.

Because $d\mu/d\theta=b''(\theta)=V(\mu)>0$, the inverse-function rule gives $\partial\theta_i/\partial\mu_i=1/V(\mu_i)$. Substitution and summation yield the stated score. The quasi-score follows by differentiating $\sum_i Q(y_i;\mu_i)/\phi$ with respect to $\boldsymbol{\beta}$. Thus the two score equations are identical. $\square$

In particular, for fixed $p$, Tweedie regression coefficients obtained by solving this equation agree with the full-likelihood coefficient estimates whenever the same solution is used. The factor $1/\phi$ does not change the roots, so this coefficient equation is also unchanged by the value of a common dispersion parameter. The result is about the estimating equation; it does not make every objective function or inferential procedure equivalent.

Worked Example: The Gamma Case and Quasi-Deviance

At $p=2$, the Tweedie model is Gamma. Take $\phi=1$ and two observations $y_1=1,y_2=3$. In an intercept-only model, the log-link score equation is $\sum_i(y_i/\mu-1)=0$, so the fitted mean is $\hat\mu=(1+3)/2=2$. The full Gamma likelihood gives the same fitted mean for fixed $\phi$.

For $V(t)=t^2$, integrate the defining quasi-likelihood derivative from $y$ to $\mu$: $$ Q(y;\mu)=\int_y^\mu\frac{y-t}{t^2}\,dt =1-\frac{y}{\mu}+\log\left(\frac{y}{\mu}\right),\qquad y>0. $$ The unit deviance contribution is $d(y,\mu)=-2Q(y;\mu)=2\{y/\mu-1-\log(y/\mu)\}$.

For $y_1=1,\mu=2$, this is $2(0.5-1-\log 0.5)=2(-0.5+0.693147)=0.386294$. For $y_2=3,\mu=2$, it is $2(1.5-1-\log 1.5)=2(0.5-0.405465)=0.189070$. The total is $0.575364$, to six decimal places.

This deviance is determined by the mean and variance function. It is not itself the full Gamma log-likelihood: the latter also includes the normalizing terms involving $y$ and $\phi$.

What the moment information can support

Result: Suppose $Y_i$ are independent, $E[Y_i]=\mu_i$, and $\operatorname{Var}(Y_i)=\phi V(\mu_i)$. If the derivatives $\partial\mu_i/\partial\boldsymbol{\beta}$ exist and the relevant variances are finite, the quasi-score has expectation zero and covariance $$ \operatorname{Cov}(\mathbf{U})= \frac{1}{\phi}\sum_{i=1}^n \frac{1}{V(\mu_i)} \frac{\partial\mu_i}{\partial\boldsymbol{\beta}} \left(\frac{\partial\mu_i}{\partial\boldsymbol{\beta}}\right)^\mathsf{T}. $$

Proof. Let $\mathbf{h}_i=\partial\mu_i/\partial\boldsymbol{\beta}$. The $i$th score contribution is $\mathbf{h}_i(Y_i-\mu_i)/\{\phi V(\mu_i)\}$. Its expectation is zero because $E[Y_i-\mu_i]=0$; consequently $E[\mathbf{U}]=\mathbf{0}$. Independence makes the covariance of the sum equal to the sum of the individual covariances. The $i$th covariance is

\frac{\mathbf{h}_i\mathbf{h}_i^\mathsf{T}}{\phi^2V(\mu_i)^2} \operatorname{Var}(Y_i) = \frac{\mathbf{h}_i\mathbf{h}_i^\mathsf{T}}{\phi V(\mu_i)}.

Summing proves the formula. This result needs only the stated first two moments and independence, not the full Tweedie density. $\square$

Under additional regularity and large-sample conditions, this covariance information supports approximate inference for regression coefficients. If the variance function is not trusted, a sandwich covariance estimate can instead use the empirical variability of the score contributions. Neither approach turns a quasi-score into a likelihood, and neither provides a likelihood for an individual response.

Worked Example: A Moment Dispersion Estimate Need Not Be the MLE

Use the Gamma data $1,3$, with an intercept-only fitted mean $\hat\mu=2$. A common Pearson moment estimate of dispersion divides the sum of squared Pearson residuals by the residual degrees of freedom $n-1=1$: $$ \hat\phi_{\mathrm{P}}= \frac{\sum_i (y_i-\hat\mu)^2/\hat\mu^2}{n-1} =\frac{1/4+1/4}{1}=0.5. $$ This calculation uses the mean–variance relationship and a degrees-of-freedom convention.

The full Gamma likelihood gives a different estimating equation. Set $k=1/\phi$, the Gamma shape. For observations with mean $\mu$, the density is $$ f(y;\mu,k)=\frac{k^k}{\Gamma(k)\mu^k}y^{k-1}\exp(-ky/\mu). $$ Profiling the mean at $\hat\mu=2$, the likelihood score equation for $k$ reduces to $$ \log k-\psi(k)=\log(2)-\frac{\log(1)+\log(3)}{2}, $$ where $\psi(k)=\frac{d}{dk}\log\Gamma(k)$ is the digamma function. The right side is $0.143841$. At the Pearson estimate $\phi=0.5$, we have $k=2$, and the left side is $\log 2-\psi(2)=0.270363$, not $0.143841$. Thus the Pearson moment estimate is not the full-likelihood estimate in this sample. The likelihood estimate instead solves the displayed equation.

Which inferential tools need a full distribution?

The distinction is especially important when the power index or dispersion is unknown. With $p$ fixed, the coefficient score can be obtained from the mean and variance functions, and the coefficient roots agree with the full-likelihood roots under the theorem’s assumptions. But choosing $p$ changes the variance function and, in the full Tweedie model, changes the response distribution and its carrier. Quasi-likelihood alone does not supply a normalized likelihood score for $p$. A full-likelihood profile estimate of $p$ therefore uses distributional information beyond the coefficient quasi-score.

Likewise, quasi-likelihood can be paired with moment-based dispersion estimation, but it does not by itself specify a likelihood estimate of $\phi$. The Gamma example shows that the Pearson and likelihood estimates can differ. Full-likelihood inference for dispersion uses the density’s normalizing terms as well as its mean–variance relationship.

Likelihood-ratio tests compare maximized full likelihoods under specified models. AIC and other likelihood information criteria also require actual likelihood values, including the normalizing terms; substituting a quasi-likelihood or quasi-deviance does not produce the usual likelihood criterion. Quasi-likelihood-based tests and criteria can be constructed with their own assumptions and calibration, but they should be identified as quasi-likelihood methods, not described as ordinary likelihood-ratio tests or AIC.

Common Pitfall

Common mistake: “The coefficient score equations match, so quasi-likelihood and full likelihood are interchangeable for every purpose.”

Correction: Matching coefficient equations establishes a specific equivalence for coefficient estimation under the stated model conditions. It does not provide the normalized response distribution, make dispersion estimates identical, or justify likelihood-ratio tests and AIC. Always check whether the calculation needs only the mean and variance, or the full density and its likelihood.

A practical rule is to begin by naming the target. For regression coefficients at fixed $p$, the quasi-score may be sufficient. For a probability, a likelihood-based estimate of dispersion or power, or a likelihood comparison, the full Tweedie distribution is needed. For uncertainty estimates, state whether the variance model is being assumed correct or whether a robust sandwich calculation is being used.

Check Your Understanding

Use the distinction between moment information and a full response distribution to answer the following.

Why do the full-likelihood and quasi-likelihood coefficient score equations agree for a fixed Tweedie power?
Which information in a full exponential-dispersion likelihood is absent from a quasi-likelihood defined only through its derivative?
For Gamma responses, calculate the unit deviance contribution when $y=2$ and $\mu=1$.
Why does a Pearson estimate of dispersion not have to equal the full-likelihood estimate?
Can ordinary AIC be calculated from a quasi-likelihood alone? Explain what is missing.