4.6 Lack of fit and pure error
The residual plots of Figure 8 suggested that a straight line has the wrong shape for the catalyst data. When the data contain replicates, repeated observations of \(y\) at the same value of \(x\), the suspicion can be tested formally.
The idea is this. The residual sum of squares measures two things at once: the random scatter of \(y\) about its true mean, and any systematic failure of the line to pass through the true means. Replicates measure the first on its own, because the scatter among observations at the same \(x\) can owe nothing to the shape of the model. If what is left over is much larger than replicate scatter would explain, the model has the wrong shape.
The decomposition
Suppose there are \(m\) distinct values \(x_1, \dots , x_m\) of the regressor, with \(n_i\) observations \(y_{i1}, \dots , y_{in_i}\) at \(x_i\), and \(n = \sum n_i\). Let \(\bar y_i\) be the mean of the observations at \(x_i\), and \(\hat y_i\) the fitted value there. The residual of \(y_{ij}\) splits as \[ y_{ij} - \hat y_i = (y_{ij} - \bar y_i) + (\bar y_i - \hat y_i) . \]
Theorem 4.6.1. \(SS_{\mathrm {Res}} = SS_{\mathrm {PE}} + SS_{\mathrm {LOF}}\), where \[ SS_{\mathrm {PE}} = \sum _{i=1}^m \sum _{j=1}^{n_i}(y_{ij} - \bar y_i)^2 , \qquad SS_{\mathrm {LOF}} = \sum _{i=1}^m n_i(\bar y_i - \hat y_i)^2 . \] \(SS_{\mathrm {PE}}\) is the pure-error sum of squares, with \(n - m\) degrees of freedom, and \(SS_{\mathrm {LOF}}\) the lack-of-fit sum of squares, with \(m - 2\).
Proof. Square the split residual and add over \(i\) and \(j\). The cross term is \[ 2\sum _{i=1}^m (\bar y_i - \hat y_i)\sum _{j=1}^{n_i}(y_{ij} - \bar y_i) = 0 , \] since deviations from a mean add to zero. The degrees of freedom are those of the distribution theory below: \(\sum _i (n_i - 1) = n - m\) for pure error, and \((n - 2) - (n - m) = m - 2\) for lack of fit. □
\(SS_{\mathrm {LOF}}\) is a weighted sum of squared distances between the fitted line and the means at each \(x\): it is small when the line passes close to the group means. At least one \(x\) must be replicated for there to be any pure error (\(n > m\)), and there must be more distinct \(x\) values than parameters (\(m > 2\)) for lack of fit to be measurable.
Theorem 4.6.2. Under A1–A4 with \(E(y_{ij}) = \mu _i\), where the means \(\mu _i\) need not lie on a line,
- (i)
- \(SS_{\mathrm {PE}}/\sigma ^2 \sim \chi ^2_{n - m}\), so \(MS_{\mathrm {PE}} = SS_{\mathrm {PE}}/(n - m)\) is unbiased for \(\sigma ^2\) whatever the true means;
- (ii)
- \(SS_{\mathrm {PE}}\) and \(SS_{\mathrm {LOF}}\) are independent;
- (iii)
- if the straight-line model is correct, \(\mu _i = \beta _0 + \beta _1 x_i\), then \(SS_{\mathrm {LOF}}/\sigma ^2 \sim \chi ^2_{m - 2}\), and \[ F_0 = \frac {SS_{\mathrm {LOF}}/(m - 2)}{SS_{\mathrm {PE}}/(n - m)} \sim F_{m - 2,\, n - m} . \]
Proof. Let \(\mathbf {Z}\) be the \(n \times m\) matrix of group indicators, with \(z_{(ij),k} = 1\) if \(k = i\) and \(0\) otherwise, and let \(\mathbf {P} = \mathbf {Z}(\mathbf {Z}^\top \mathbf {Z})^{-1}\mathbf {Z}^\top \) be its hat matrix. Regressing \(y\) on \(\mathbf {Z}\) fits a separate mean to each group, so \(\mathbf {Py}\) replaces each \(y_{ij}\) by \(\bar y_i\), and \(\mathbf {y}^\top (\mathbf {I} - \mathbf {P})\mathbf {y} = SS_{\mathrm {PE}}\). Like any hat matrix, \(\mathbf {P}\) is symmetric and idempotent with trace \(m\) (Theorem 3.4.2).
Both columns of \(\mathbf {X}\), the ones and the \(x\) values, are constant within groups, so they are combinations of the columns of \(\mathbf {Z}\); hence \(\mathbf {PX} = \mathbf {X}\) and, taking transposes and using symmetry, \(\mathbf {X}^\top \mathbf {P} = \mathbf {X}^\top \). Then \(\mathbf {PH} = \mathbf {H}\) and \(\mathbf {HP} = \mathbf {H}\), and so \(\mathbf {P} - \mathbf {H}\) is symmetric, with \[ (\mathbf {P} - \mathbf {H})^2 = \mathbf {P} - \mathbf {H} - \mathbf {H} + \mathbf {H} = \mathbf {P} - \mathbf {H}, \qquad (\mathbf {I} - \mathbf {P})(\mathbf {P} - \mathbf {H}) = \mathbf {0} . \] Its trace is \(m - 2\), and \(\mathbf {y}^\top (\mathbf {P} - \mathbf {H})\mathbf {y} = \mathbf {y}^\top (\mathbf {I} - \mathbf {H})\mathbf {y} - \mathbf {y}^\top (\mathbf {I} - \mathbf {P})\mathbf {y} = SS_{\mathrm {Res}} - SS_{\mathrm {PE}} = SS_{\mathrm {LOF}}\).
Write \(\mathbf {y} = \boldsymbol \mu + \boldsymbol \varepsilon \), where \(\boldsymbol \mu \) has entry \(\mu _i\) for each observation in group \(i\); it is a combination of the columns of \(\mathbf {Z}\), so \((\mathbf {I} - \mathbf {P})\boldsymbol \mu = \mathbf {0}\).
(i) \(SS_{\mathrm {PE}} = \boldsymbol \varepsilon ^\top (\mathbf {I} - \mathbf {P})\boldsymbol \varepsilon \), and Lemma 3.7.1(ii) applies, with \(\operatorname {tr}(\mathbf {I} - \mathbf {P}) = n - m\).
(ii) \(SS_{\mathrm {PE}}\) is a function of \((\mathbf {I} - \mathbf {P})\mathbf {y}\) and \(SS_{\mathrm {LOF}} = \bigl ((\mathbf {P} - \mathbf {H})\mathbf {y}\bigr )^\top (\mathbf {P} - \mathbf {H})\mathbf {y}\) of \((\mathbf {P} - \mathbf {H})\mathbf {y}\). These are jointly normal with covariance \(\sigma ^2(\mathbf {I} - \mathbf {P})(\mathbf {P} - \mathbf {H}) = \mathbf {0}\), so independent by (F1).
(iii) If \(\boldsymbol \mu = \mathbf {X}\boldsymbol \beta \) then \((\mathbf {P} - \mathbf {H})\boldsymbol \mu = \mathbf {X}\boldsymbol \beta - \mathbf {X}\boldsymbol \beta = \mathbf {0}\), so \(SS_{\mathrm {LOF}} = \boldsymbol \varepsilon ^\top (\mathbf {P} - \mathbf {H})\boldsymbol \varepsilon \) and Lemma 3.7.1(ii) gives \(\chi ^2_{m - 2}\). The \(F\) distribution follows from (i), (ii) and (F3). □
When the line is wrong, \((\mathbf {P} - \mathbf {H})\boldsymbol \mu \neq \mathbf {0}\), and a calculation like that of Lemma 3.1.4 gives \[ E(SS_{\mathrm {LOF}}) = (m - 2)\sigma ^2 + \boldsymbol \mu ^\top (\mathbf {P} - \mathbf {H})\boldsymbol \mu , \] where the second term is \(\sum _i n_i(\mu _i - \tilde \mu _i)^2 > 0\), \(\tilde \mu _i\) being the value at \(x_i\) of the least-squares line fitted to the true means. So \(MS_{\mathrm {LOF}}\) estimates \(\sigma ^2\) when the model is right and something larger when it is wrong, while \(MS_{\mathrm {PE}}\) estimates \(\sigma ^2\) in either case. A large \(F_0\) is evidence of lack of fit.
The lack-of-fit test
To test \(H_0\): the straight-line model is adequate, against \(H_1\): it is not, compute \[ F_0 = \frac {MS_{\mathrm {LOF}}}{MS_{\mathrm {PE}}} \] and reject \(H_0\) at level \(\alpha \) if \(F_0 > F_{\alpha , m - 2, n - m}\). The calculations go in an extension of the analysis of variance table:
| Source | SS | df | MS | \(F_0\) |
| Regression | \(SS_R\) | \(1\) | \(MS_R\) | |
| Residual | \(SS_{\mathrm {Res}}\) | \(n - 2\) | ||
| Lack of fit | \(SS_{\mathrm {LOF}}\) | \(m - 2\) | \(MS_{\mathrm {LOF}}\) | \(MS_{\mathrm {LOF}}/MS_{\mathrm {PE}}\) |
| Pure error | \(SS_{\mathrm {PE}}\) | \(n - m\) | \(MS_{\mathrm {PE}}\) | |
| Total | \(SS_T\) | \(n - 1\) |
If \(H_0\) is not rejected, there is no evidence against the model, and \(MS_{\mathrm {Res}}\), which pools both sources and has more degrees of freedom, is used to estimate \(\sigma ^2\) in the tests and intervals that follow. If \(H_0\) is rejected, the model should be revised before any inference is drawn from it.
Everything extends to multiple regression, where the replicates are repeated observations at the same combination of regressor values. With \(m\) distinct combinations and \(p\) parameters, \(SS_{\mathrm {LOF}}\) has \(m - p\) degrees of freedom; the proof is unchanged, with \(\operatorname {tr}(\mathbf {P} - \mathbf {H}) = m - p\).
Remark. Replicates must be genuine: separate runs of the process at the same settings, each subject to all the sources of variation that affect any other run. Measuring the same sample twice, or splitting one batch, captures only measurement error, which is smaller than the true run-to-run variation. \(MS_{\mathrm {PE}}\) is then too small, and the test finds lack of fit that is not there.
Example 4.6.3. Test the straight-line model for the catalyst data of Example 4.3.1 for lack of fit, at \(\alpha = 0.05\).
Solution. There are \(m = 5\) distinct concentrations, each with \(n_i = 2\) observations, \(n = 10\).
Pure error. Each pair differs by \(1\), so each pair contributes \((\tfrac 12)^2 + (\tfrac 12)^2 = 0.5\), and \(SS_{\mathrm {PE}} = 5 \times 0.5 = 2.5\) on \(n - m = 5\) degrees of freedom.
Lack of fit. From Example 4.3.1, \(SS_{\mathrm {Res}} = 15.1\), so \(SS_{\mathrm {LOF}} = 15.1 - 2.5 = 12.6\) on \(m - 2 = 3\) degrees of freedom. As a check, directly, with group means \(4.5, 7.5, 10.5, 10.5, 10.5\) and fitted values \(5.7, 7.2, 8.7, 10.2, 11.7\), \[ SS_{\mathrm {LOF}} = 2\bigl [(-1.2)^2 + 0.3^2 + 1.8^2 + 0.3^2 + (-1.2)^2\bigr ] = 2 \times 6.3 = 12.6 . \] The table. \(SS_R = \hat \beta _1 S_{xy} = 1.5 \times 30 = 45\).
| Source | SS | df | MS | \(F_0\) |
| Regression | 45.0 | 1 | 45.0 | |
| Residual | 15.1 | 8 | ||
| Lack of fit | 12.6 | 3 | 4.2 | 8.40 |
| Pure error | 2.5 | 5 | 0.5 | |
| Total | 60.1 | 9 |
\(F_{0.05,3,5} = 5.41\), and \(F_0 = 8.40 > 5.41\) (the \(P\)-value is \(0.021\)). We reject \(H_0\): the straight line does not adequately describe how yield depends on concentration.
This confirms what the plot of residuals against fitted values suggested. The means rise from \(4.5\) to \(10.5\) over the first three concentrations and then stay level, and a straight line cannot follow that. Notice too that the line’s \(MS_{\mathrm {Res}} = 1.89\) is nearly four times the pure-error estimate \(0.5\) of \(\sigma ^2\), because the lack of fit has inflated it. A curved model, for instance one with an \(x^2\) term, should be tried.
Questions on this section
Stuck on something here? Ask below and it stays attached to this topic.