10.2 Orthogonal polynomials
The centred quadratic still had \(\mathbf {1}\) and \(\mathbf {u}^2\) correlated, and adding a cubic term would change every lower coefficient. Both problems disappear if the powers are replaced by polynomials that are orthogonal over the observed \(x\) values. Let \[ y_i = \alpha _0P_0(x_i) + \alpha _1P_1(x_i) + \dots + \alpha _kP_k(x_i) + \varepsilon _i , \] where \(P_j\) is a polynomial of degree exactly \(j\), \(P_0 = 1\), and \[ \sum _{i=1}^n P_j(x_i)P_l(x_i) = 0 \qquad (j \ne l) . \] The \(P_j\) can be built by the Gram–Schmidt process: \(P_j\) is \(x^j\) minus its projection on \(P_0, \dots , P_{j-1}\). Since each \(P_j\) has degree exactly \(j\), the columns \(P_0, \dots , P_k\) span the same space as \(1, x, \dots , x^k\), so by the argument of Theorem 10.1.2 the fitted values are those of the ordinary polynomial.
Theorem 10.2.1. If the columns \(\mathbf {p}_0, \dots , \mathbf {p}_k\) of \(\mathbf {X}\) are mutually orthogonal, then
- (i)
- \(\hat \alpha _j = \mathbf {p}_j^\top \mathbf {y}/\mathbf {p}_j^\top \mathbf {p}_j\), which does not depend on which other columns are in the model;
- (ii)
- \(\operatorname {Var}(\hat \alpha _j) = \sigma ^2/\mathbf {p}_j^\top \mathbf {p}_j\), and the \(\hat \alpha _j\) are uncorrelated;
- (iii)
- the regression sum of squares splits as \(SS_R = \sum _{j \ge 1}SS(\alpha _j)\) with \(SS(\alpha _j) = (\mathbf {p}_j^\top \mathbf {y})^2/\mathbf {p}_j^\top \mathbf {p}_j\), and \(SS(\alpha _j)\) is also the extra sum of squares for adding \(\mathbf {p}_j\) to a model containing any of the other columns.
Proof. \(\mathbf {X}^\top \mathbf {X} = \operatorname {diag}(\mathbf {p}_0^\top \mathbf {p}_0, \dots , \mathbf {p}_k^\top \mathbf {p}_k)\), so \(\hat {\boldsymbol \alpha } = (\mathbf {X}^\top \mathbf {X})^{-1}\mathbf {X}^\top \mathbf {y}\) has the components in (i), each involving only its own column, and \(\operatorname {Var}(\hat {\boldsymbol \alpha }) = \sigma ^2(\mathbf {X}^\top \mathbf {X})^{-1}\) is diagonal, which is (ii). For (iii), with \(\mathbf {p}_0 = \mathbf {1}\), \[ SS_R = \hat {\boldsymbol \alpha }^\top \mathbf {X}^\top \mathbf {y} - n\bar y^2 = \sum _{j=0}^k \hat \alpha _j\mathbf {p}_j^\top \mathbf {y} - \frac {(\mathbf {1}^\top \mathbf {y})^2}{n} = \sum _{j=1}^k \frac {(\mathbf {p}_j^\top \mathbf {y})^2}{\mathbf {p}_j^\top \mathbf {p}_j} , \] because the \(j = 0\) term is \((\mathbf {1}^\top \mathbf {y})^2/n\). The same calculation for any subset of the columns gives the sum of the same terms over the subset, so adding \(\mathbf {p}_j\) always adds exactly \(SS(\alpha _j)\). □
So with orthogonal polynomials, the order can be raised one term at a time without refitting, and the sequence of tests in consideration (2) above reads straight off one table.
Equally spaced values
When \(x_i = a + (i - 1)d\), \(i = 1, \dots , n\), write \(u_i = (x_i - \bar x)/d\), which takes the values \(-(n - 1)/2, \dots , (n - 1)/2\) in steps of one. The first two orthogonal polynomials are \[ P_1 = u , \qquad P_2 = u^2 - \frac {n^2 - 1}{12} . \] To check this: the \(u_i\) are symmetric about \(0\), so \(\sum u_i = \sum u_i^3 = 0\) (Theorem 10.1.3), which gives \(\sum P_1 = 0\) and \(\sum P_1P_2 = \sum u^3 - \frac {n^2 - 1}{12}\sum u = 0\). For \(\sum P_2 = 0\) we need \(\sum u_i^2 = n(n^2 - 1)/12\). With \(u_i = i - (n + 1)/2\), \begin {align*} \sum u_i^2 &= \sum i^2 - n\Bigl (\frac {n + 1}{2}\Bigr )^2 = \frac {n(n + 1)(2n + 1)}{6} - \frac {n(n + 1)^2}{4} \\ &= \frac {n(n + 1)\bigl [2(2n + 1) - 3(n + 1)\bigr ]}{12} = \frac {n(n + 1)(n - 1)}{12} . \end {align*}
Higher polynomials are found the same way. The polynomials are usually scaled by constants \(\lambda _j\) so that their values are small integers. For \(n = 7\), \(u = -3, \dots , 3\), \(P_2 = u^2 - 4\) and \(P_3 = (u^3 - 7u)/6\):
| \(u\) | \(-3\) | \(-2\) | \(-1\) | \(0\) | \(1\) | \(2\) | \(3\) | \(\sum P_j^2\) |
| \(P_1\) | \(-3\) | \(-2\) | \(-1\) | \(0\) | \(1\) | \(2\) | \(3\) | 28 |
| \(P_2\) | \(5\) | \(0\) | \(-3\) | \(-4\) | \(-3\) | \(0\) | \(5\) | 84 |
| \(P_3\) | \(-1\) | \(1\) | \(1\) | \(0\) | \(-1\) | \(-1\) | \(1\) | 6 |
Each row sums to zero, and the products of any two rows sum to zero, as can be checked directly.
Example 10.2.2. Fit polynomials up to the cubic to the irrigation data of Example 10.1.4, using orthogonal polynomials, and decide the order.
Solution. With \(y = 8, 16, 23, 26, 27, 29, 26\) and the table above, \begin {align*} \textstyle \sum P_1y &= -24 - 32 - 23 + 0 + 27 + 58 + 78 = 84 , \\ \textstyle \sum P_2y &= 40 + 0 - 69 - 104 - 81 + 0 + 130 = -84 , \\ \textstyle \sum P_3y &= -8 + 16 + 23 + 0 - 27 - 29 + 26 = 1 . \end {align*}
By Theorem 10.2.1, \(\hat \alpha _0 = \bar y = 22.143\), \(\hat \alpha _1 = 84/28 = 3\), \(\hat \alpha _2 = -84/84 = -1\) and \(\hat \alpha _3 = 1/6 = 0.167\), and these values do not change as terms are added or removed. The sums of squares are \begin {gather*} SS(\alpha _1) = \frac {84^2}{28} = 252, \qquad SS(\alpha _2) = \frac {84^2}{84} = 84, \\ SS(\alpha _3) = \frac {1^2}{6} = 0.167 . \end {gather*} From \(SS_T = 338.857\), the residual for the cubic model is \(338.857 - 252 - 84 - 0.167 = 2.690\) on \(3\) degrees of freedom.
| Source | SS | df | MS | \(F_0\) |
| Linear | 252.000 | 1 | 252.000 | 281.0 |
| Quadratic | 84.000 | 1 | 84.000 | 93.7 |
| Cubic | 0.167 | 1 | 0.167 | 0.19 |
| Residual | 2.690 | 3 | 0.897 | |
| Total | 338.857 | 6 |
The cubic term is far from significant (\(F_{0.05,1,3} = 10.13\)), so the quadratic is enough. Refitting the quadratic only moves the cubic sum of squares into the residual, \(2.690 + 0.167 = 2.857\) on \(4\) degrees of freedom, and the quadratic test becomes \(F_0 = 84/0.714 = 117.6\), which is \(t_0^2 = 10.84^2\) from Example 10.1.4. The fitted quadratic is \[ \hat y = 22.143 + 3P_1 - P_2 = 22.143 + 3u - (u^2 - 4) = 26.143 + 3u - u^2 , \] the same as before.
Questions on this section
Stuck on something here? Ask below and it stays attached to this topic.