Partial Regression With Nuisance Predictors

Why Residualize Before Estimating a Coefficient Block?

A coefficient for a target predictor block describes its association with several outcomes after accounting for nuisance predictors. One way to obtain that adjusted association is to fit all predictors and outcomes together. The Frisch–Waugh–Lovell result gives an equivalent route: first remove from the outcomes and target predictors everything explained by the nuisance predictors, then regress the residualized outcomes on the residualized target predictors.

This equivalence is useful for understanding adjustment, checking a computation, and separating the role of nuisance covariates from the target block. The residualization is an algebraic property of ordinary least squares. It does not require a Gaussian error distribution or a particular covariance structure among the outcomes.

Model and Residualization

Let \(Y\) be an \(n\times q\) matrix of responses, \(Z\) an \(n\times r\) matrix of nuisance predictors, and \(X\) an \(n\times p\) matrix of target predictors. Write the multivariate linear model as

Y=Z\Gamma+X\Theta+U,

where \(\Gamma\) is the \(r\times q\) nuisance coefficient matrix, \(\Theta\) is the \(p\times q\) target coefficient matrix, and \(U\) is the error matrix. Assume for now that \(Z\) has full column rank and that the combined design matrix \([Z\ X]\) has full column rank. Ordinary least squares minimizes the sum of squared entries of \(Y-Z\Gamma-X\Theta\), equivalently the trace of its residual sum-of-squares-and-cross-products matrix.

Definition: The projection matrix onto the column space of \(Z\) and its residual-maker matrix are

P_Z=Z(Z'Z)^{-1}Z', \qquad M_Z=I_n-P_Z.

For any \(n\)-row matrix \(A\), \(M_ZA\) is the matrix of residuals from projecting each column of \(A\) onto the columns of \(Z\). In particular, the residualized outcomes are \(M_ZY\), and the residualized target predictors are \(M_ZX\).

The residual-maker is symmetric and idempotent: \(M_Z'=M_Z\) and \(M_Z^2=M_Z\). Also \(Z'M_Z=0\). Thus every column of \(M_ZA\) is orthogonal to every column of \(Z\). If \(Z\) contains an intercept, residualizing against \(Z\) removes both the intercept component and the components explained by the other nuisance predictors. A separate intercept should not then be added to the residual-on-residual regression.

The Multivariate Frisch–Waugh–Lovell Result

Theorem: Under the full-rank assumptions above, the ordinary least-squares target coefficient matrix in the full regression is

\widehat{\Theta} =(X'M_ZX)^{-1}X'M_ZY.

This is exactly the coefficient matrix obtained by regressing \(M_ZY\) on \(M_ZX\), with no additional intercept. The corresponding nuisance coefficient matrix is

\widehat{\Gamma} =(Z'Z)^{-1}Z'(Y-X\widehat{\Theta}).

The inverse in the theorem exists under the stated assumptions. Indeed, for any \(a\in\mathbb R^p\),

a'X'M_ZXa =(Xa)'M_Z(Xa) =\|M_ZXa\|^2.

This quantity is zero exactly when \(Xa\) lies in the column space of \(Z\). If a nonzero \(a\) had this property, a nonzero linear combination of the columns of \(X\) would be explained by \(Z\), contradicting full column rank of \([Z\ X]\). Hence \(X'M_ZX\) is positive definite.

To prove the theorem, fix a candidate target coefficient matrix \(\Theta\). The least-squares choice of \(\Gamma\) for that fixed \(\Theta\) is obtained by projecting \(Y-X\Theta\) onto the columns of \(Z\):

\widehat{\Gamma}(\Theta) =(Z'Z)^{-1}Z'(Y-X\Theta).

Substituting this choice leaves the residual matrix

Y-Z\widehat{\Gamma}(\Theta)-X\Theta =M_Z(Y-X\Theta).

Therefore, minimizing the full regression criterion over both coefficient matrices is equivalent, after minimizing over \(\Gamma\), to minimizing

\|M_ZY-M_ZX\Theta\|_F^2

over \(\Theta\), where \(\|\cdot\|_F^2\) is the sum of squared matrix entries. Differentiating this quadratic criterion with respect to \(\Theta\), or setting its least-squares normal equations to zero, gives

X'M_Z(M_ZY-M_ZX\widehat{\Theta})=0.

Since \(M_Z'=M_Z\) and \(M_Z^2=M_Z\), this becomes \(X'M_ZY-X'M_ZX\widehat{\Theta}=0\). The invertibility just established yields the stated formula for \(\widehat{\Theta}\). Substitution into \(\widehat{\Gamma}(\Theta)\) gives the formula for \(\widehat{\Gamma}\). \(\square\)

Worked Examples

Worked Example: Adjusting for an Intercept

There are three observations, one target predictor \(x=(0,1,2)'\), and two outcomes with columns \(y_1=(1,2,4)'\) and \(y_2=(3,0,3)'\). Let \(Z=\mathbf 1_3\), so the only nuisance predictor is an intercept. Calculate the target coefficient row using residualization.

The predictor mean is \(1\), so

M_Zx=x-\bar{x}\mathbf 1_3=(-1,0,1)'.

The outcome means are \(7/3\) and \(2\). Thus

M_ZY= \begin{pmatrix} -4/3&1\\ -1/3&-2\\ 5/3&1 \end{pmatrix}.

The residualized predictor has squared length \(2\). Its inner products with the residualized outcome columns are

(-1,0,1) \begin{pmatrix} -4/3&1\\ -1/3&-2\\ 5/3&1 \end{pmatrix} =(3,0).

Consequently, the target coefficient matrix is the \(1\times2\) row \((3/2,0)\). The adjusted slope is \(1.5\) for the first outcome and \(0\) for the second. The intercepts can be recovered from the outcome means and predictor mean: they are \(7/3-(3/2)(1)=5/6\) and \(2-0(1)=2\), respectively.

Worked Example: A Nuisance Covariate Changes the Slope

Let \(Z\) contain an intercept and \(z=(-1,-1,1,1)'\). Take the target predictor \(x=(0,2,2,4)'\), and define \(v=(-1,1,-1,1)'\) and \(e=(1,-1,-1,1)'\). The two outcomes are

y_1=5\mathbf 1+2z+3v+e=(1,5,3,11)', \qquad y_2=-\mathbf 1+z-2v+2e=(2,-6,0,0)'.

First residualize the target predictor. Its projection onto the intercept is \(2\mathbf 1\); its projection onto \(z\) is \(z\), since \(z'z=4\) and \(z'(x-2\mathbf 1)=4\). Therefore

M_Zx=x-2\mathbf 1-z=(-1,1,-1,1)'=v.

The vectors \(v\) and \(e\) are each orthogonal to both columns of \(Z\), and \(v'e=0\), while \(v'v=4\). It follows that \(M_Zy_1=3v+e\) and \(M_Zy_2=-2v+2e\). Their inner products with \(v\) are \(12\) and \(-8\), so the adjusted target slopes are \(12/4=3\) and \(-8/4=-2\).

For comparison, fitting each outcome on an intercept and \(x\) alone gives slopes based on \(x-2\mathbf 1=(-2,0,0,2)'\). The inner products with \(y_1\) and \(y_2\) are \(20\) and \(-4\); the denominator is \(8\). The unadjusted slopes are therefore \(20/8=2.5\) and \(-4/8=-0.5\), not the adjusted slopes. The difference occurs because the unadjusted fit does not remove the association of \(x\) with the nuisance predictor \(z\).

Worked Example: Two Target Predictors and Two Outcomes

Let \(Z\) again contain an intercept and \(z=(-1,-1,1,1)'\). Let the two target predictor columns be \(a=(1,-1,1,-1)'\) and \(b=(1,-1,-1,1)'\), and set \(X=[a\ b]\). Both columns are orthogonal to \(Z\), to each other, and have squared length \(4\). Consider the data

Y= \begin{pmatrix} 6&3\\ 0&-1\\ 8&-5\\ 6&3 \end{pmatrix}.

These outcomes can be written as \(Y=\mathbf 1(5,0)+z(2,-1)+a(2,-1)+b(1,3)\). The nuisance part is removed by \(M_Z\), leaving \(M_ZY=a(2,-1)+b(1,3)\). Also \(M_ZX=X\) and \(X'X=4I_2\). Directly,

X'M_ZY= \begin{pmatrix} 8&-4\\ 4&12 \end{pmatrix}, \qquad \widehat{\Theta} =(4I_2)^{-1} \begin{pmatrix} 8&-4\\ 4&12 \end{pmatrix} = \begin{pmatrix} 2&-1\\ 1&3 \end{pmatrix}.

Thus the first target predictor has coefficients \(2\) and \(-1\) for the two outcomes, while the second has coefficients \(1\) and \(3\). The result estimates the whole coefficient block at once; the same residualized design is used for every outcome.

What the Result Does—and Does Not—Say

The Frisch–Waugh–Lovell result says that the target block from the joint ordinary least-squares fit can be recovered after removing the nuisance-predictor projections from both sides of the target relationship. Residualizing only \(Y\), or only \(X\), generally does not produce the same coefficient. Residualizing both is essential because the target predictors may themselves be associated with the nuisance predictors.

The coefficient matrix is an adjusted linear association: it describes the part of the outcome variation associated with target predictors after their linear overlap with \(Z\) is removed. This algebra alone does not establish a causal effect. Nor does it say that residualized OLS is interchangeable with every alternative fitting procedure. In particular, the result concerns the ordinary least-squares coefficient matrix; weighting observations or using a different criterion requires its own derivation.

Key takeaway: Fit the target block in the full model, or first residualize both \(Y\) and \(X\) against the nuisance design \(Z\) and fit one on the other. Under full column rank, both routes give \((X'M_ZX)^{-1}X'M_ZY\).

Common Mistakes

Residualizing only the outcomes. This leaves the target predictors’ nuisance-related component in the regression. Use \(M_ZY\) and \(M_ZX\).
Adding an intercept after residualizing against an intercept. The residualized columns already have zero sample means. Adding another intercept changes the stipulated residual-on-residual fit.
Interpreting the result as a new model with different target coefficients. The residual regression is a computational representation of the target block in the original joint model; it gives the same block, not a different estimand.
Assuming the identity applies to any estimator. The derivation uses the ordinary least-squares criterion and its projection geometry. A weighted or generalized least-squares fit must use the corresponding geometry.

Check Your Understanding

Use the model and residualization definitions above to answer the following.

In \(Y=Z\Gamma+X\Theta+U\), which two matrices must be residualized to recover \(\widehat{\Theta}\) by partial regression?
Why does residualizing \(Y\) against \(Z\) remove a nuisance component from every outcome at once?
If \(Z\) includes an intercept, why is an additional intercept ordinarily omitted from the regression of \(M_ZY\) on \(M_ZX\)?
State the full-rank condition that ensures \(X'M_ZX\) is invertible, and explain its connection to overlap between the target and nuisance predictor spaces.
In the second worked example, why does the unadjusted slope for the first outcome differ from the adjusted slope?