Canonical Exponential-Family Regression Models

Why the Natural Parameter Matters

In regression, observations may have different distributions because their predictors differ. The exponential-family representation puts those distributions into a common form, so that the contribution of each observation to the likelihood can be differentiated systematically. When the natural parameter is a linear function of the regression coefficients, the resulting score and Fisher information have a particularly direct form.

This is a theoretical derivation, not a recipe for fitting a model. The key point is to track how the likelihood depends on the coefficient vector through each observation's natural parameter. A canonical linear predictor does not generally make the response mean linear in the predictors.

Definition: Let \(Y_1,\ldots,Y_n\) be independent responses with densities or mass functions \(f_i(y_i;\theta_i,\phi)=\exp\{[y_i\theta_i-b(\theta_i)]/a(\phi)+c_i(y_i,\phi)\}\). Here \(\theta_i\) is the natural parameter, \(a(\phi)>0\) is a scale factor, and \(\phi\) is treated as known. A canonical exponential-family regression model specifies \(\theta_i=x_i^T\beta\), where \(x_i\) is a fixed \(p\)-dimensional vector and \(\beta\) is the unknown \(p\)-dimensional coefficient vector. The permissible values of each \(\theta_i\) lie in the open natural-parameter space of the family.

The vectors \(x_i\) are regarded as fixed. We assume the response supports do not depend on \(\beta\), the natural-parameter space is open, and the density or mass function is regular enough to differentiate its normalizing integral or sum. In particular, the relevant derivatives can be passed under that integral or sum. These are the standard regular exponential-family conditions needed for the moment and score calculations below.

Moments From the Exponential-Family Form

The function \(b\) determines the first two moments. To see why, fix one observation and temporarily write its response as \(Y\) and its natural parameter as \(\theta\). Normalization says that integrating or summing the density over its support gives one. Differentiating that identity with respect to \(\theta\) gives

0=E\!\left[\frac{Y-b'(\theta)}{a(\phi)}\right], \qquad\text{so}\qquad E[Y]=b'(\theta).

Differentiating the normalization identity a second time gives

0=E\!\left[ \frac{(Y-b'(\theta))^2}{a(\phi)^2} -\frac{b''(\theta)}{a(\phi)} \right].

Since the first derivative established \(E[Y]=b'(\theta)\), this second identity implies \(\operatorname{Var}(Y)=a(\phi)b''(\theta)\). Applying these facts separately to the independent responses gives

E[Y_i]=b'(\theta_i), \qquad \operatorname{Var}(Y_i)=a(\phi)b''(\theta_i).

Thus the mean is generally a function of the natural parameter, not the natural parameter itself. The coefficient vector enters the mean through \(b'(x_i^T\beta)\).

Worked Example: Bernoulli Responses

For a Bernoulli response with success probability \(p\), the canonical representation has \(a(\phi)=1\), \(\theta=\log\{p/(1-p)\}\), and \(b(\theta)=\log(1+e^\theta)\). Differentiating gives \(b'(\theta)=e^\theta/(1+e^\theta)=p\) and \(b''(\theta)=e^\theta/(1+e^\theta)^2=p(1-p)\). The general moment formulas therefore return \(E[Y]=p\) and \(\operatorname{Var}(Y)=p(1-p)\).

For example, at \(\theta=\log 2\), \(p=2/3\), and the variance is \((2/3)(1/3)=2/9\). The predictor is linear in the log odds, not in the success probability.

Deriving the Regression Score and Information

Write \(X\) for the \(n\)-by-\(p\) matrix whose \(i\)th row is \(x_i^T\). Independence makes the log-likelihood the sum of the individual log-likelihoods. With \(\theta_i=x_i^T\beta\), its terms that depend on \(\beta\) are

\ell(\beta)= \sum_{i=1}^n \left\{ \frac{y_i x_i^T\beta-b(x_i^T\beta)}{a(\phi)} +c_i(y_i,\phi) \right\}.

Theorem: Under the regularity conditions above, the score vector for \(\beta\) is \(U(\beta)=\nabla_\beta\ell(\beta)=\sum_{i=1}^n x_i\{y_i-b'(\theta_i)\}/a(\phi)\). The observed information and Fisher information are both \(I(\beta)=\sum_{i=1}^n b''(\theta_i)x_ix_i^T/a(\phi)\).

To derive the score, apply the chain rule to each summand. Since \(\nabla_\beta\theta_i=x_i\), differentiating \(y_i\theta_i-b(\theta_i)\) gives \(x_i\{y_i-b'(\theta_i)\}\). For the second derivative, \(x_i\) is fixed and \(\nabla_\beta b'(\theta_i)=b''(\theta_i)x_i\). Therefore

\nabla_\beta\ell(\beta) = \sum_{i=1}^n\frac{x_i\{y_i-b'(\theta_i)\}}{a(\phi)}, \qquad -\nabla_\beta^2\ell(\beta) = \sum_{i=1}^n\frac{b''(\theta_i)}{a(\phi)}x_ix_i^T.

The negative Hessian is the observed information: in this canonical model it does not depend on the observed responses. Taking its expectation therefore leaves it unchanged, giving the stated Fisher information.

There is also a direct covariance check. The \(i\)th score contribution is \(x_i\{Y_i-b'(\theta_i)\}/a(\phi)\), with expectation zero. Independence makes the covariance between contributions from distinct observations zero. Using \(\operatorname{Var}(Y_i)=a(\phi)b''(\theta_i)\), the covariance of contribution \(i\) is \(b''(\theta_i)x_ix_i^T/a(\phi)\). Summing yields

E[U(\beta)]=0, \qquad \operatorname{Cov}(U(\beta)) = \sum_{i=1}^n\frac{b''(\theta_i)}{a(\phi)}x_ix_i^T =I(\beta).

This proves the score has mean zero and verifies the information calculation from its covariance. In matrix form, with \(W\) diagonal and \(W_{ii}=b''(\theta_i)/a(\phi)\), the result is \(I(\beta)=X^TWX\).

Identifiability consequence: If \(b''(\theta_i)>0\) for every \(i\), then \(I(\beta)\) is positive definite exactly when \(X\) has full column rank. Indeed, for every vector \(v\), \(v^TI(\beta)v=\sum_i [b''(\theta_i)/a(\phi)](x_i^Tv)^2\). This is positive for every nonzero \(v\) precisely when no nonzero \(v\) is orthogonal to every row of \(X\).

Worked Score and Information Calculations

Worked Example: Bernoulli Regression With an Intercept and Slope

Take three independent Bernoulli responses with \(x_1=(1,0)^T\), \(x_2=(1,1)^T\), and \(x_3=(1,2)^T\). Let \(\beta=(0,\log 2)^T\), and suppose the observed responses are \(y=(1,0,1)\). Their natural parameters are \(0,\log 2,\log 4\), so their success probabilities are \(1/2,2/3,4/5\).

For Bernoulli responses, \(a(\phi)=1\) and \(b'(\theta_i)=p_i\). The score is

U(\beta) =\sum_i x_i(y_i-p_i) = \begin{pmatrix}1\\0\end{pmatrix}\frac12 +\begin{pmatrix}1\\1\end{pmatrix}\left(-\frac23\right) +\begin{pmatrix}1\\2\end{pmatrix}\frac15 = \begin{pmatrix}1/30\\-4/15\end{pmatrix}.

The information weights are \(p_i(1-p_i)=1/4,2/9,4/25\). Thus

I(\beta) =\sum_i p_i(1-p_i)x_ix_i^T = \begin{pmatrix} 569/900 & 122/225\\ 122/225 & 194/225 \end{pmatrix}.

For instance, the lower-right entry is \(0^2(1/4)+1^2(2/9)+2^2(4/25)=2/9+16/25=194/225\). Each entry is determined by predictor vectors and model probabilities, not by the realized response values.

Worked Example: Poisson Regression With an Intercept

Let \(Y_i\) be independent Poisson responses with natural parameters \(\theta_i=x_i\beta\), using scalar predictors \(x=(1,1,2)\). For the Poisson family, \(a(\phi)=1\), \(b(\theta)=e^\theta\), and the mean is \(\mu_i=e^{\theta_i}\). Set \(\beta=\log 2\), giving means \((2,2,4)\). For observed responses \(y=(1,3,5)\), the score is

U(\beta)=\sum_i x_i(y_i-\mu_i) =(1)(1-2)+(1)(3-2)+(2)(5-4)=-1+1+2=2.

Because \(b''(\theta_i)=e^{\theta_i}=\mu_i\), the scalar Fisher information is

I(\beta)=\sum_i \mu_i x_i^2 =(2)(1^2)+(2)(1^2)+(4)(2^2)=2+2+16=20.

The score uses the observed-minus-mean residuals, while the information weights each squared predictor by the corresponding model mean.

Worked Example: Gaussian Responses and the Natural Parameter

Suppose independent responses have normal distributions with known variance \(\sigma^2=4\), and set their natural parameters to \(\theta_i=x_i\beta\) for \(x=(1,2)\) and scalar \(\beta=1/4\). In this representation, \(\theta_i=\mu_i/\sigma^2\), so the means are \(\mu_i=\sigma^2\theta_i=4x_i\beta=(1,2)\). For observed values \(y=(2,1)\), the score and information are

U(\beta)=\sum_i\frac{x_i(y_i-\mu_i)}{\sigma^2} =\frac{(1)(1)+(2)(-1)}{4}=-\frac14, \qquad I(\beta)=\sum_i\frac{x_i^2}{\sigma^2} =\frac{1+4}{4}=\frac54.

Here \(b''(\theta)=\sigma^2\), so the general information weight \(b''(\theta_i)/a(\phi)\) equals \(4\) when \(a(\phi)=1\), consistent with the formula above after expressing the Gaussian score through its mean.

What Canonical Does—and Does Not—Mean

The canonical condition is \(\theta_i=x_i^T\beta\): it is the natural parameter that is linear in the predictors. Since \(E[Y_i]=b'(\theta_i)\), the response mean is \(b'(x_i^T\beta)\), which need not be linear. In the Bernoulli case this gives a logistic mean; in the Poisson case it gives an exponential mean. For a normal response with known variance, the natural parameter is the mean divided by that variance.

A common mistake is to see a linear predictor and conclude that the expected response is linear in the predictors. That reasoning silently replaces \(b'(\theta_i)\) by \(\theta_i\), which is valid only for particular forms of \(b\), not for exponential-family models in general. The correction is to keep the two steps separate: first set \(\theta_i=x_i^T\beta\), then obtain the mean as \(b'(\theta_i)\).

The score and information formulas matter because they connect the distributional form of each response to the geometry of coefficient-based inference. They show exactly how each observation contributes: its score is its predictor vector multiplied by a centered response, and its information is a weighted outer product of that predictor vector. The weights come from the curvature \(b''\) of the family’s normalizing function. These identities provide the theoretical structure for likelihood inference without requiring a particular applied fitting procedure.

Check Your Understanding

Use the natural-parameter representation and the derived formulas to answer the following.

For a canonical model, express \(E[Y_i]\) and \(\operatorname{Var}(Y_i)\) in terms of \(b\), \(a(\phi)\), and \(\theta_i\).
Why does the observed information equal the Fisher information in the canonical model described here?
For Bernoulli responses, what are \(b'(\theta)\) and \(b''(\theta)\), and what do they represent?
State the condition on the design matrix that makes the information matrix positive definite when all \(b''(\theta_i)\) are positive.
Why does a linear natural parameter not generally imply a linear response mean?