3.1 The general linear model

A linear model takes the form \[ \mathbf Y = X\boldsymbol \beta + \boldsymbol \varepsilon , \] where \(\mathbf Y\) is the \(n \times 1\) vector of responses, \(X\) is a known \(n \times p\) matrix (the design matrix), \(\boldsymbol \beta \) is a \(p \times 1\) vector of unknown parameters and \(\boldsymbol \varepsilon \) is the vector of random errors. The model may be a regression model, in which the columns of \(X\) are measured variables, or a model for a designed experiment, in which they are indicator (dummy) variables for the treatments. The assumptions are

(i): \(E(\varepsilon _i) = 0\), so that \(E(\mathbf Y) = X\boldsymbol \beta \);
(ii): \(\Var (\varepsilon _i) = \sigma ^2\), a common variance (homoscedastic errors);
(iii): \(\Cov (\varepsilon _i, \varepsilon _j) = 0\) for \(i \ne j\);
(iv): for tests and confidence intervals, the \(\varepsilon _i\) are normally distributed: \(\boldsymbol \varepsilon \sim N(\mathbf 0, \sigma ^2I)\).

Here \(X\) and \(\boldsymbol \beta \) are non-random while \(\mathbf Y\) and \(\boldsymbol \varepsilon \) are random, and \(\Var (\mathbf Y) = \Var (\boldsymbol \varepsilon ) = \sigma ^2I\).

Least squares

The method of least squares estimates \(\boldsymbol \beta \) by minimising \[ S(\boldsymbol \beta ) = \|\mathbf Y - X\boldsymbol \beta \|^2 = \sum _{i=1}^n\bigl (Y_i - E(Y_i)\bigr )^2 \] with respect to \(\boldsymbol \beta \).

Theorem 3.1.1. \(\hat {\boldsymbol \beta }\) minimises \(S(\boldsymbol \beta )\) if and only if it satisfies the normal equations \[ X^\top X\hat {\boldsymbol \beta } = X^\top \mathbf Y . \] Solutions always exist, and the fitted vector \(\hat {\mathbf Y} = X\hat {\boldsymbol \beta }\) is the same for every solution. If \(X\) has full column rank \(p\), the solution is unique: \(\hat {\boldsymbol \beta } = (X^\top X)^{-1}X^\top \mathbf Y\).

Proof. Let \(\hat {\boldsymbol \beta }\) satisfy the normal equations, so that \(X^\top (\mathbf Y - X\hat {\boldsymbol \beta }) = \mathbf 0\). For any \(\boldsymbol \beta \), \[ \|\mathbf Y - X\boldsymbol \beta \|^2 = \|\mathbf Y - X\hat {\boldsymbol \beta }\|^2 + \|X(\hat {\boldsymbol \beta } - \boldsymbol \beta )\|^2 + 2(\hat {\boldsymbol \beta } - \boldsymbol \beta )^\top X^\top (\mathbf Y - X\hat {\boldsymbol \beta }) , \] and the last term is zero. So \(S(\boldsymbol \beta ) \ge S(\hat {\boldsymbol \beta })\), with equality exactly when \(X\boldsymbol \beta = X\hat {\boldsymbol \beta }\); this shows that \(\hat {\boldsymbol \beta }\) is a minimiser, that every minimiser has the same fitted vector, and (taking \(\boldsymbol \beta \) to be another minimiser) that every minimiser satisfies the normal equations. Solutions exist because \(X^\top \mathbf Y\) lies in the column space of \(X^\top \), which equals the column space of \(X^\top X\). If \(X\) has rank \(p\), \(X^\top X\) is invertible. □

Geometrically, \(\hat {\mathbf Y}\) is the orthogonal projection of \(\mathbf Y\) onto the column space \(\mathcal C(X)\), and we write \(\hat {\mathbf Y} = H\mathbf Y\), where \(H\) is the projection matrix (or hat matrix) onto \(\mathcal C(X)\). When \(X\) has full rank, \(H = X(X^\top X)^{-1}X^\top \). In every case \(H\) is symmetric and idempotent (\(H^2 = H\)), and its rank, which equals its trace, is the rank of \(X\). The residual vector is \(\mathbf e = \mathbf Y - \hat {\mathbf Y} = (I - H)\mathbf Y\), and the residual (error) sum of squares is \[ SSE = \mathbf e^\top \mathbf e = \mathbf Y^\top (I - H)\mathbf Y = \sum _{i=1}^n(Y_i - \hat Y_i)^2 . \] When \(X\) has full rank, \(E(\hat {\boldsymbol \beta }) = \boldsymbol \beta \) and \(\Cov (\hat {\boldsymbol \beta }) = \sigma ^2(X^\top X)^{-1}\), and under normality \(\hat {\boldsymbol \beta } \sim N\bigl (\boldsymbol \beta , \sigma ^2(X^\top X)^{-1}\bigr )\).

Quadratic forms

Sums of squares are quadratic forms \(\mathbf Y^\top A\mathbf Y\) in the observations, and two results about them carry most of the theory.

Lemma 3.1.2. If \(E(\mathbf Y) = \boldsymbol \mu \) and \(\Var (\mathbf Y) = \Sigma \), then for any symmetric matrix \(A\), \[ E(\mathbf Y^\top A\mathbf Y) = \operatorname {tr}(A\Sigma ) + \boldsymbol \mu ^\top A\boldsymbol \mu . \]

Proof. \(\mathbf Y^\top A\mathbf Y\) is a scalar, so it equals its trace, and \(\operatorname {tr}(A\mathbf Y\mathbf Y^\top )\) by the cyclic property. Since \(E(\mathbf Y\mathbf Y^\top ) = \Sigma + \boldsymbol \mu \boldsymbol \mu ^\top \), \[ E(\mathbf Y^\top A\mathbf Y) = \operatorname {tr}\bigl (A(\Sigma + \boldsymbol \mu \boldsymbol \mu ^\top )\bigr ) = \operatorname {tr}(A\Sigma ) + \boldsymbol \mu ^\top A\boldsymbol \mu . \] □

Applied to \(SSE\) with \(\Sigma = \sigma ^2I\) and \(\boldsymbol \mu = X\boldsymbol \beta \), for which \((I - H)X\boldsymbol \beta = \mathbf 0\), the lemma gives \(E(SSE) = \sigma ^2\operatorname {tr}(I - H) = \sigma ^2(n - r)\), where \(r\) is the rank of \(X\). So \(MSE = SSE/(n - r)\) is an unbiased estimator of \(\sigma ^2\), and no normality is needed for this.

Definition 3.1.3. If \(Z_1, \dots , Z_k\) are independent with \(Z_i \sim N(\delta _i, 1)\), then \(\sum Z_i^2\) has the non-central chi-square distribution with \(k\) degrees of freedom and non-centrality parameter \(\lambda = \sum \delta _i^2\), written \(\chi ^2_k(\lambda )\). Its distribution depends on the \(\delta _i\) only through \(\lambda \), its mean is \(k + \lambda \), and \(\chi ^2_k(0) = \chi ^2_k\). If \(V_1 \sim \chi ^2_{k_1}(\lambda )\) and \(V_2 \sim \chi ^2_{k_2}\) are independent, \((V_1/k_1)/(V_2/k_2)\) has the non-central \(F\) distribution \(F_{k_1,k_2}(\lambda )\).

Theorem 3.1.4. Let \(\mathbf Y \sim N(\boldsymbol \mu , \sigma ^2I_n)\), and let \(A\) and \(B\) be symmetric idempotent matrices of ranks \(a\) and \(b\).

(i): \(\mathbf Y^\top A\mathbf Y/\sigma ^2 \sim \chi ^2_a(\lambda )\), with \(\lambda = \boldsymbol \mu ^\top A\boldsymbol \mu /\sigma ^2\).
(ii): If \(AB = 0\), then \(\mathbf Y^\top A\mathbf Y\) and \(\mathbf Y^\top B\mathbf Y\) are independent.

Proof. (i) A symmetric idempotent matrix has eigenvalues \(0\) and \(1\), so by the spectral theorem \(A = PP^\top \), where the \(n \times a\) matrix \(P\) has orthonormal columns, \(P^\top P = I_a\). Then \(\mathbf Y^\top A\mathbf Y = \|P^\top \mathbf Y\|^2\), and \(\mathbf Z = P^\top \mathbf Y/\sigma \sim N(P^\top \boldsymbol \mu /\sigma , I_a)\) has independent unit-variance normal components. Hence \(\mathbf Y^\top A\mathbf Y/\sigma ^2 = \|\mathbf Z\|^2 \sim \chi ^2_a(\lambda )\) with \(\lambda = \|P^\top \boldsymbol \mu \|^2/\sigma ^2 = \boldsymbol \mu ^\top A\boldsymbol \mu /\sigma ^2\).

(ii) Similarly write \(B = QQ^\top \). Since \(AB = 0\), \[ 0 = P^\top (PP^\top QQ^\top )Q = (P^\top P)(P^\top Q)(Q^\top Q) = P^\top Q . \] So the jointly normal vectors \(P^\top \mathbf Y\) and \(Q^\top \mathbf Y\) have covariance \(\sigma ^2P^\top Q = 0\) and are independent, and so are their squared lengths \(\mathbf Y^\top A\mathbf Y\) and \(\mathbf Y^\top B\mathbf Y\). □

[next] [front] [up]

Questions on this section

Stuck on something here? Ask below and it stays attached to this topic.