Weighted Least-Squares Geometry In GLM Fitting

Why Fisher Scoring Has a Least-Squares Geometry

A generalized linear model is not generally fitted by minimizing ordinary squared residuals. Yet a Fisher-scoring iteration for a GLM can be computed by solving a weighted least-squares problem. The connection is local: at the current fitted means, the model supplies a working response and weights, and the next coefficient vector is the weighted least-squares fit to that response.

This representation is useful because it exposes the geometry of an iteration. After multiplying each design row and working response by the square root of its weight, the update is an ordinary least-squares projection in a transformed space. The weights and working response are recalculated as the fitted means change, so this geometry describes one step at a time rather than a fixed least-squares problem for the original data.

Definition: For independent GLM observations with linear predictors $\eta_i=x_i^{\mathsf T}\beta+o_i$, fitted means $\mu_i=g^{-1}(\eta_i)$, and positive prior weights $q_i$, define the working response and working weight at the current parameter value by $$ z_i=\eta_i+(y_i-\mu_i)\frac{d\eta_i}{d\mu_i}, \qquad W_i=\frac{q_i}{\phi V(\mu_i)} \left(\frac{d\mu_i}{d\eta_i}\right)^2. $$ Here $o_i$ is a known offset, $g$ is a differentiable link, $\operatorname{Var}(Y_i)=\phi V(\mu_i)/q_i$, and $W_i>0$ is evaluated at the current fit.

From the Score to a Weighted Least-Squares Step

Assume the observations belong to a regular one-parameter exponential-dispersion family, with differentiable positive variance function $V(\mu)$, and that the current means are interior points of the mean space. For this family, the derivative of the $i$th log-likelihood contribution with respect to $\mu_i$ is $q_i(y_i-\mu_i)/\{\phi V(\mu_i)\}$. Applying the chain rule to $\mu_i=g^{-1}(\eta_i)$ gives the score for $\beta$:

U(\beta)=\sum_{i=1}^n x_i \frac{q_i}{\phi V(\mu_i)} \frac{d\mu_i}{d\eta_i}(y_i-\mu_i).

The Fisher information is the expected outer product of the score contributions. Independence makes the cross-observation terms vanish, and the variance of $Y_i$ is $\phi V(\mu_i)/q_i$. Thus the information matrix is

\mathcal I(\beta)=X^{\mathsf T}WX, \qquad W=\operatorname{diag}(W_1,\ldots,W_n).

Fisher scoring adds an increment $\delta$ satisfying $\mathcal I(\beta)\delta=U(\beta)$. If $X^{\mathsf T}WX$ is nonsingular, this specifies a unique increment. The next result shows that the same update is obtained by fitting the working response.

Theorem: Under the conditions above, suppose $X^{\mathsf T}WX$ is nonsingular at the current parameter value. The Fisher-scoring update $\beta_{\mathrm{new}}=\beta+\delta$ is the unique minimizer of $$ \sum_{i=1}^n W_i\left(z_i-o_i-x_i^{\mathsf T}\beta_{\mathrm{new}}\right)^2. $$

Proof. By the definitions of $z_i$ and $W_i$,

W_i(z_i-\eta_i) = \frac{q_i}{\phi V(\mu_i)} \left(\frac{d\mu_i}{d\eta_i}\right)^2 (y_i-\mu_i)\frac{d\eta_i}{d\mu_i} = \frac{q_i}{\phi V(\mu_i)} \frac{d\mu_i}{d\eta_i}(y_i-\mu_i).

Consequently $U(\beta)=X^{\mathsf T}W(z-o-X\beta)$, while $\mathcal I(\beta)=X^{\mathsf T}WX$. The scoring equation for $\delta$ becomes

X^{\mathsf T}WX\delta=X^{\mathsf T}W(z-o-X\beta),

or, since $\beta_{\mathrm{new}}=\beta+\delta$,

X^{\mathsf T}W(z-o-X\beta_{\mathrm{new}})=0.

These are the weighted least-squares normal equations for the stated objective. Its quadratic part has matrix $X^{\mathsf T}WX$, which is positive definite under the nonsingularity assumption, so the solution is the unique minimizer. This is exactly the Fisher-scoring update. $\square$

The Projection View

Let $W^{1/2}$ be the diagonal matrix with entries $\sqrt{W_i}$, and define the transformed design matrix $A=W^{1/2}X$ and transformed working response $b=W^{1/2}(z-o)$. The weighted least-squares objective is

\sum_{i=1}^n W_i(z_i-o_i-x_i^{\mathsf T}\beta)^2 = \|b-A\beta\|^2.

Thus the fitted transformed response $A\beta_{\mathrm{new}}$ is the ordinary Euclidean orthogonal projection of $b$ onto the column space of $A$. The residual vector $b-A\beta_{\mathrm{new}}$ is orthogonal to every column of $A$, which is another way to express the normal equations.

Geometric result: For any real matrix $A$, the least-squares fitted vector is the orthogonal projection of $b$ onto $\mathcal C(A)$, the column space of $A$. If $A$ has full column rank, its coefficient vector is unique. If it is rank deficient, the projected fitted vector is still unique, but more than one coefficient vector may produce it.

Proof. The column space $\mathcal C(A)$ is a finite-dimensional subspace of Euclidean space, so each $b$ has a unique orthogonal decomposition $b=p+r$, where $p\in\mathcal C(A)$ and $r$ is orthogonal to every vector in $\mathcal C(A)$. Since $p\in\mathcal C(A)$, there is at least one coefficient vector $\beta_0$ such that $p=A\beta_0$. For any other $\beta$, $A\beta-p$ belongs to $\mathcal C(A)$, so it is orthogonal to $r$. Hence

\|b-A\beta\|^2 = \|r-(A\beta-p)\|^2 = \|r\|^2+\|A\beta-p\|^2 \geq \|r\|^2.

Equality holds exactly when $A\beta=p$. Therefore $p$ is the unique fitted vector minimizing the squared distance, although its coefficient representation need not be unique. If $A$ has full column rank, $A\beta=p$ has at most one solution, so the coefficient vector is unique as well. $\square$

In a GLM step, this projection is induced by the current weights: the rows of $X$ are rescaled by $\sqrt{W_i}$, and the working response, after removing the offset, is rescaled in the same way. Changing the fitted means generally changes $W$, $z$, and therefore the projection space and target vector for the next step.

Worked Updates in Three GLMs

Worked Example: A Logistic Fisher-Scoring Step

Consider two grouped binomial observations, with trial counts $q=(10,20)$ and observed proportions $p=(0.4,0.7)$. Fit an intercept-only logistic model, so $X$ is a two-row column of ones. At the current value $\beta=0$, both fitted probabilities are $\mu_i=0.5$. For the binomial proportion, $V(\mu)=\mu(1-\mu)$, $\phi=1$, and the logit link has $d\eta/d\mu=1/\{\mu(1-\mu)\}$.

The working weights are $W_i=q_i\mu_i(1-\mu_i)$, giving $W_1=10(0.25)=2.5$ and $W_2=20(0.25)=5$. The working responses are

z_1=0+\frac{0.4-0.5}{0.25}=-0.4, \qquad z_2=0+\frac{0.7-0.5}{0.25}=0.8.

There is no offset. The weighted least-squares estimate of the intercept is the weighted mean:

\beta_{\mathrm{new}} = \frac{2.5(-0.4)+5(0.8)}{2.5+5} = \frac{-1+4}{7.5} =0.4.

Checking the scoring equation directly, the score at zero is $10(0.4-0.5)+20(0.7-0.5)=3$, and the information is $10(0.25)+20(0.25)=7.5$. The increment is $3/7.5=0.4$, agreeing with the working-response calculation.

Worked Example: A Poisson Step with a Binary Covariate

Suppose three independent Poisson counts are $y=(1,4,2)$, with unit prior weights and log link. Use the design rows $(1,0),(1,0),(1,1)$, no offset, and current coefficient vector $\beta=(0,0)^{\mathsf T}$. All fitted means are $1$. For a Poisson model, $V(\mu)=\mu$, and $d\eta/d\mu=1/\mu$, so $W_i=\mu_i=1$.

The working response is $z_i=0+(y_i-1)/1$, giving $z=(0,3,1)$. All weights equal one, so the least-squares fit is an ordinary line fit to these three working responses. The first two rows have covariate value zero, so the fitted intercept is their average, $(0+3)/2=1.5$. The final row has covariate value one, and its fitted value must be $1$, so the slope is $1-1.5=-0.5$. Therefore $\beta_{\mathrm{new}}=(1.5,-0.5)^{\mathsf T}$.

For a check, $X^{\mathsf T}X=\begin{pmatrix}3&1\\1&1\end{pmatrix}$ and $X^{\mathsf T}z=(4,1)^{\mathsf T}$. The equations $3\beta_0+\beta_1=4$ and $\beta_0+\beta_1=1$ give $\beta_0=1.5$, $\beta_1=-0.5$, as above.

Worked Example: A Weighted Gaussian Projection

Take a Gaussian identity-link model with known $\phi=1$, prior weights $q=(1,4,1)$, design rows $(1,0),(1,1),(1,2)$, response $y=(1,3,2)$, and no offset. For this model $V(\mu)=1$, $d\mu/d\eta=1$, and the working weights are $W=q$. The working response is $z=y$, regardless of the current coefficient vector.

The weighted normal equations are

\begin{pmatrix}6&6\\6&8\end{pmatrix} \begin{pmatrix}\beta_0\\\beta_1\end{pmatrix} = \begin{pmatrix}15\\16\end{pmatrix}.

Indeed, the left matrix is $X^{\mathsf T}WX$, and the right side is $X^{\mathsf T}Wy$. Subtracting the first equation $6\beta_0+6\beta_1=15$ from the second gives $2\beta_1=1$, so $\beta_1=0.5$ and $\beta_0=2$. The fitted values are $(2,2.5,3)$, giving residuals $(-1,0.5,-1)$. Their weighted residuals are $(-1,2,-1)$; the inner products with the intercept and slope columns are $-1+2-1=0$ and $0+2-2=0$. This verifies the projection’s orthogonality directly.

What the Geometry Does—and Does Not—Mean

The weighted least-squares view makes the computation of a Fisher-scoring step concrete and shows how the current fit scales different observations in transformed space. It does not turn the GLM likelihood into a single fixed least-squares objective. For a non-Gaussian model, both the working responses and weights typically depend on the current parameter values; a new scoring iteration changes them.

The geometry also clarifies the role of the design. The transformed design matrix determines the subspace onto which the transformed working response is projected. A response value affects the target vector through $z$, while the current means and variance function determine the weights that reshape the design. These are distinct ingredients of the iteration.

Common pitfall: It is incorrect to say that a GLM Fisher-scoring step is ordinary least squares on the observed responses. The step fits the working responses, not generally $y_i$, and uses weights that depend on the current fitted means. For example, in the Poisson calculation above, the working response was $y_i-1$, not $y_i$. The weighted least-squares problem is a representation of one scoring step, not a replacement for the original GLM likelihood.

Check Your Understanding

Use the definitions and derivation above to answer the following.

For a log-link model, write the working response in terms of $y_i$, $\mu_i$, and the current linear predictor.
Why does the Fisher-scoring information matrix have the form $X^{\mathsf T}WX$? State the role of independence.
What condition makes the coefficient vector in a weighted least-squares step unique? What remains unique if that condition fails?
In the logistic example, explain why the working weights were $2.5$ and $5$, rather than both being $0.25$.
Explain why the transformed projection generally changes between GLM iterations.