Degrees-of-Freedom Correction for Residual Variance

Why Correct the Residual Sum of Squares?

In Gaussian regression, residuals measure the differences between observed responses and fitted means. Their squared sum is a natural measure of unexplained variation, but it tends to be too small as an estimate of the response variance: the fitted coefficients were chosen using the same data, so fitting absorbs some variation. The degrees-of-freedom correction accounts for this effect.

The key fact is geometric. With $n$ observations and $p$ linearly independent fitted coefficients, the residual vector lies in an $n-p$-dimensional subspace. Under the Gaussian model, its squared length has expectation $\sigma^2(n-p)$. Dividing the residual sum of squares by $n-p$, rather than by $n$, therefore gives an unbiased estimate of $\sigma^2$.

Model and Variance Estimator

Let $Y=(Y_1,\ldots,Y_n)^\mathsf{T}$ be the response vector, and let $X$ be a fixed $n$-by-$p$ design matrix of rank $p$, where $p<n$. The Gaussian identity-link model is

Y=X\beta+\varepsilon,\qquad \varepsilon\sim N_n(0,\sigma^2 I_n).

Here $\beta$ is the vector of $p$ regression coefficients, and $\sigma^2>0$ is the common response variance. Equivalently, conditional on the design, the responses are independent with $Y_i\sim N(x_i^\mathsf{T}\beta,\sigma^2)$. The least-squares estimate is $\hat\beta=(X^\mathsf{T}X)^{-1}X^\mathsf{T}Y$, which is well-defined because $X$ has full column rank. The fitted response vector is $X\hat\beta$, and the residual vector is $e=Y-X\hat\beta$.

Definition: The residual sum of squares is $\operatorname{RSS}=e^\mathsf{T}e$. The degrees-of-freedom-corrected residual variance estimator is $$ \hat\sigma^2_{\mathrm{res}}=\frac{\operatorname{RSS}}{n-p}. $$ The number $n-p$ is the residual degrees of freedom: the dimension remaining after fitting the $p$-dimensional mean model.

The count $p$ is the rank of the design matrix, not necessarily the number of variable names in a model description. For example, if a column is an exact linear combination of other columns, it adds no independent fitted direction. The full-rank assumption ensures that the $p$ columns represent $p$ distinct directions.

Why the Residual Space Has Dimension n minus p

Define the projection matrix onto the column space of $X$ by $H=X(X^\mathsf{T}X)^{-1}X^\mathsf{T}$. Then $X\hat\beta=HY$, so $e=(I_n-H)Y$. The matrix $H$ is a symmetric idempotent projection of rank $p$; $I_n-H$ projects onto the orthogonal complement of the model space and has rank $n-p$.

Because $X\beta$ belongs to the column space of $X$, $(I_n-H)X\beta=0$. Thus the residual vector is determined only by the projected error:

e=(I_n-H)\varepsilon,\qquad \operatorname{RSS}=\varepsilon^\mathsf{T}(I_n-H)\varepsilon.

This identity also explains why estimating the mean does not leave $n$ unconstrained residual components. The residuals must be orthogonal to every column of $X$, imposing $p$ independent linear constraints.

Distribution and Unbiasedness

Theorem: Suppose $Y=X\beta+\varepsilon$, with fixed $X$ of rank $p<n$, and $\varepsilon\sim N_n(0,\sigma^2 I_n)$ for $\sigma^2>0$. Then $$ \frac{\operatorname{RSS}}{\sigma^2}\sim\chi^2_{n-p}, \qquad E[\operatorname{RSS}]=(n-p)\sigma^2, \qquad E[\hat\sigma^2_{\mathrm{res}}]=\sigma^2. $$

Proof. The matrix $I_n-H$ is symmetric and idempotent, and its rank is $n-p$. By the spectral theorem, there is an orthogonal matrix $U$ such that $U^\mathsf{T}(I_n-H)U$ is diagonal with $n-p$ entries equal to 1 and $p$ entries equal to 0. Set $Z=U^\mathsf{T}\varepsilon/\sigma$. Since $\varepsilon/\sigma\sim N_n(0,I_n)$ and $U$ is orthogonal, the components of $Z$ are independent standard normal random variables. Therefore

\frac{\operatorname{RSS}}{\sigma^2} =\frac{\varepsilon^\mathsf{T}(I_n-H)\varepsilon}{\sigma^2} =\sum_{j=1}^{n-p}Z_j^2 \sim\chi^2_{n-p}.

A chi-square random variable with $n-p$ degrees of freedom has expectation $n-p$. Multiplying the displayed random variable by $\sigma^2$ gives $E[\operatorname{RSS}]=(n-p)\sigma^2$. Dividing both sides by $n-p$ gives $E[\hat\sigma^2_{\mathrm{res}}]=\sigma^2$. $\square$

The exact chi-square statement uses Gaussian errors, but the expectation calculation itself needs less. If the errors have mean zero and covariance matrix $\sigma^2 I_n$, then

E[\operatorname{RSS}] =E[\varepsilon^\mathsf{T}(I_n-H)\varepsilon] =\operatorname{tr}\{(I_n-H)E[\varepsilon\varepsilon^\mathsf{T}]\} =\sigma^2\operatorname{tr}(I_n-H) =(n-p)\sigma^2.

The trace identity follows by expanding the quadratic form into entries and using linearity of expectation. Normality is what supplies the chi-square distribution, not the unbiasedness calculation under these mean and covariance conditions.

Worked Example: Correcting an Intercept-Only Estimate

Suppose four independent observations from an intercept-only Gaussian model are $2,4,6,8$. Estimate the common mean by the sample average, then compare the residual variance estimate using $n-p$ with the maximum-likelihood estimate using $n$.

There is one fitted coefficient, the intercept, so $n=4$, $p=1$, and the residual degrees of freedom are $4-1=3$. The fitted mean is

\bar Y=\frac{2+4+6+8}{4}=5.

The residuals are $-3,-1,1,3$, and their squared sum is

\operatorname{RSS}=(-3)^2+(-1)^2+1^2+3^2=9+1+1+9=20.

The degrees-of-freedom-corrected estimate is $20/3\approx6.67$. The Gaussian maximum-likelihood estimate of the variance is $20/4=5$. These are different estimates from this particular data set; the unbiasedness claim is about repeated samples from the stated model. Across such samples, the residual sum of squares has expectation $3\sigma^2$, so its division by 3 has expectation $\sigma^2$.

Worked Example: Checking Residual Degrees of Freedom in Regression

Consider the model with an intercept and one slope, using the five design values $x=(-2,-1,0,1,2)$ and responses $y=(0,-1,5,3,8)$. Find the fitted line and the corrected residual variance estimate.

The proposed fitted line is $\hat y=3+2x$. Its fitted values are $(-1,1,3,5,7)$, giving residuals $e=(1,-2,2,-2,1)$. To verify these are least-squares residuals, check that they are orthogonal to both design columns:

\sum_{i=1}^5 e_i=1-2+2-2+1=0,\qquad \sum_{i=1}^5 x_i e_i=(-2)(1)+(-1)(-2)+0(2)+1(-2)+2(1)=0.

The residual vector is therefore orthogonal to the intercept and slope columns, so the displayed line is the least-squares fit. There are $p=2$ independent fitted coefficients and $n=5$ observations. Hence

\operatorname{RSS}=1^2+(-2)^2+2^2+(-2)^2+1^2=14,\qquad \hat\sigma^2_{\mathrm{res}}=\frac{14}{5-2}=\frac{14}{3}.

The divisor is 3, not 5 and not 4: estimating the intercept and slope uses two fitted directions. The corrected estimate is unbiased under the full-rank Gaussian model, although this observed value need not equal the true variance.

Extension to Known Precision Weights

Now suppose observations are independent Gaussian responses with a common dispersion $\sigma^2$ and known positive precision weights $w_i$:

Y_i\sim N(x_i^\mathsf{T}\beta,\sigma^2/w_i),\qquad w_i>0.

The weights describe relative precision: larger $w_i$ means smaller response variance. Let $W=\operatorname{diag}(w_1,\ldots,w_n)$, and assume $X$ has rank $p<n$. Weighted least squares minimizes $(Y-Xb)^\mathsf{T}W(Y-Xb)$. Write its minimized value as

\operatorname{WRSS}=(Y-X\hat\beta)^\mathsf{T}W(Y-X\hat\beta).

Theorem: Under the independent Gaussian model with known positive weights above and full-rank $X$, $$ \frac{\operatorname{WRSS}}{\sigma^2}\sim\chi^2_{n-p}, \qquad E[\operatorname{WRSS}]=(n-p)\sigma^2. $$ Thus $\operatorname{WRSS}/(n-p)$ is unbiased for the common dispersion $\sigma^2$.

Proof. Let $W^{1/2}$ be the diagonal matrix with entries $\sqrt{w_i}$, and define $Y^*=W^{1/2}Y$, $X^*=W^{1/2}X$, and $\varepsilon^*=W^{1/2}\varepsilon$. Since the weights are positive, $X^*$ also has rank $p$. The model becomes $Y^*=X^*\beta+\varepsilon^*$, and the components of $\varepsilon^*$ are independent $N(0,\sigma^2)$. Moreover,

\operatorname{WRSS} =(Y^*-X^*\hat\beta)^\mathsf{T}(Y^*-X^*\hat\beta).

The weighted least-squares fit is ordinary least squares for the transformed response and design. Applying the projection argument to $X^*$, the residual sum of squares divided by $\sigma^2$ has a chi-square distribution with $n-p$ degrees of freedom. Its expectation is therefore $(n-p)\sigma^2$, proving the result. $\square$

Worked Example: Correcting a Weighted Residual Sum

Three independent responses are $0,2,4$, with known precision weights $1,2,1$, in a weighted intercept-only model. Estimate the common mean and dispersion.

There is one coefficient, so $n-p=3-1=2$. The weighted fitted mean is

\hat\mu=\frac{1(0)+2(2)+1(4)}{1+2+1}=\frac{8}{4}=2.

The residuals are $-2,0,2$. Weighting their squared values gives

\operatorname{WRSS}=1(-2)^2+2(0)^2+1(2)^2=4+0+4=8.

The unbiased dispersion estimate is $8/2=4$. The weights affect the fitted mean and the residual sum, but because they are known positive precision weights, they do not change the residual degrees-of-freedom count: it remains $n-p$.

What the Correction Does—and Does Not—Say

The Gaussian likelihood estimates dispersion by dividing the residual sum of squares by $n$ (or the weighted residual sum by $n$ in the known-weight model). The corrected estimate instead divides by $n-p$. The former is the maximum-likelihood estimate; the latter is unbiased under the conditions established above. Unbiasedness alone does not make one estimate preferable for every purpose, and it does not assert that any single observed estimate equals the true dispersion.

The same logic applies conditionally when the design is random: if, conditional on $X$, the stated Gaussian model and full-rank condition hold, then the conditional expectation of the corrected estimate is $\sigma^2$. Averaging over $X$ then gives the same unconditional expectation, provided the model has the same $\sigma^2$ for the designs considered.

The correction should not be applied mechanically if the weights are not known precision weights, the response variances do not follow the stated model, or the fitted mean is not a rank-$p$ linear subspace as assumed. The formula depends on the model and fitting procedure, not merely on counting the rows in a data table.

Common Mistakes

Using $n$ because there are $n$ residuals. The residuals satisfy $p$ independent constraints from the fitted model. The expected residual sum of squares is $(n-p)\sigma^2$, not $n\sigma^2$.
Counting parameters without checking rank. The correct $p$ is the rank of the design matrix. Redundant columns do not add independent fitted directions.
Assuming the corrected estimate is the maximum-likelihood estimate. In the Gaussian model, the likelihood estimate uses divisor $n$; the degrees-of-freedom correction gives an unbiased estimate instead.
Changing the degrees of freedom just because weights appear. With known positive precision weights, the whitened model still has $n$ observations and a rank-$p$ mean space, so the divisor remains $n-p$.

Key takeaway: For a full-rank Gaussian identity-link model, the residual sum of squares has expectation $(n-p)\sigma^2$. With known positive precision weights, the weighted residual sum has the same expectation. Dividing by the residual degrees of freedom gives an unbiased estimate of dispersion.

Check Your Understanding

Use the assumptions of the tutorial when answering.

A Gaussian regression has 12 observations and a design matrix of rank 4. What is the residual degrees of freedom, and what divisor gives an unbiased residual variance estimate?
Why does the Gaussian result give a chi-square distribution for the residual sum of squares after scaling by $\sigma^2$?
A model description lists five columns, but one is an exact linear combination of the other four. How many independent fitted directions are there?
With known positive precision weights, why does whitening the responses preserve the $n-p$ residual degrees of freedom?
In what sense do the maximum-likelihood and degrees-of-freedom-corrected estimates differ?