---
title: Applications to Linear Models
module: Orthogonality and Least Squares
moduleNumber: 6
lessonNumber: 5
order: 605
summary: >
  Curve fitting is a least-squares problem in statistical notation. The
  least-squares line, polynomial fits, and multiple regression all reduce to
  X beta = y with a design matrix X built from the data, solved by the same
  normal equations.
topics: [Orthogonality and Least Squares]
sources:
  - book: Lay
    ref: "§6.6 Applications to Linear Models"
---

Fitting a formula to data is the most common use of least squares. A relationship
$y = \beta_0 + \beta_1 x$ is posited, the data rarely fall exactly on any such
line, and the parameters $\beta_0, \beta_1$ are chosen to make the fit as good as
possible. Written in matrix form the problem is
[least squares](/linear-algebra/orthogonality-least-squares/least-squares-problems)
with statistical notation: $X\boldsymbol{\beta} = \mathbf{y}$, where $X$ is the
**design matrix**, $\boldsymbol{\beta}$ the **parameter vector**, and $\mathbf{y}$
the **observation vector**.

## The least-squares line

Given data points $(x_1, y_1), \dots, (x_n, y_n)$, the line $y = \beta_0 +
\beta_1 x$ assigns to each $x_j$ a **predicted** value $\beta_0 + \beta_1 x_j$.
The gap between the observed $y_j$ and the predicted value is the **residual**.

> **Definition (Least-squares line).** The **least-squares line** $y = \beta_0 +
> \beta_1 x$ is the line minimizing the sum of the squares of the residuals,
> $\sum_{j=1}^n (y_j - \beta_0 - \beta_1 x_j)^2$. It is also called the **line of
> regression of $y$ on $x$.**

$$
% caption: A least-squares line minimizes the total squared vertical residual;
% each residual is the gap between the observed data point and the point directly
% above or below it on the line.
\begin{tikzpicture}[>=stealth, scale=1.0, font=\small]
\definecolor{acc}{HTML}{4A6FA5}
\draw[->, black] (-0.3,0) -- (6.0,0) node[right] {$x$};
\draw[->, black] (0,-0.3) -- (0,3.4) node[above] {$y$};
\draw[acc, thick] (0.3,0.6) -- (5.3,2.85) node[above, acc] {f\/itted line};
% data points and residuals
\foreach \x/\y/\ly in {1.0/1.6/0.94, 2.2/1.2/1.46, 3.4/2.6/1.99, 4.6/2.1/2.51} {
  \fill[black] (\x,\y) circle (1.9pt);
  \draw[black, thick] (\x,\y) -- (\x,\ly);
}
\node[black] at (3.7,2.35) {residual};
\end{tikzpicture}
$$

If the points were exactly on the line, the parameters would satisfy $n$ equations
$\beta_0 + \beta_1 x_j = y_j$, i.e. $X\boldsymbol{\beta} = \mathbf{y}$ with
$$
X = \begin{bmatrix} 1 & x_1 \\ 1 & x_2 \\ \vdots & \vdots \\ 1 & x_n \end{bmatrix},
\qquad
\boldsymbol{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \end{bmatrix},
\qquad
\mathbf{y} = \begin{bmatrix} y_1 \\ y_2 \\ \vdots \\ y_n \end{bmatrix}.
$$
The squared distance $\lVert X\boldsymbol{\beta} - \mathbf{y} \rVert^2$ is the
sum of squared residuals, so the least-squares $\boldsymbol{\beta}$ solves the
normal equations $X^\top X\boldsymbol{\beta} = X^\top \mathbf{y}$.

> **Worked example.** Fit the least-squares line to four points:
>
> | $x_j$ | $2$ | $5$ | $7$ | $8$ |
> | --- | --- | --- | --- | --- |
> | $y_j$ | $1$ | $2$ | $3$ | $3$ |
>
> Build the design matrix and observation vector
> $$
> X = \begin{bmatrix} 1 & 2 \\ 1 & 5 \\ 1 & 7 \\ 1 & 8 \end{bmatrix},
> \qquad
> \mathbf{y} = \begin{bmatrix} 1 \\ 2 \\ 3 \\ 3 \end{bmatrix},
> $$
> and compute
> $$
> X^\top X = \begin{bmatrix} 4 & 22 \\ 22 & 142 \end{bmatrix},
> \qquad
> X^\top \mathbf{y} = \begin{bmatrix} 9 \\ 57 \end{bmatrix}.
> $$
> Inverting the $2 \times 2$ matrix,
> $$
> \boldsymbol{\beta}
>   = \frac{1}{84}\begin{bmatrix} 142 & -22 \\ -22 & 4 \end{bmatrix}
>     \begin{bmatrix} 9 \\ 57 \end{bmatrix}
>   = \begin{bmatrix} 2/7 \\ 5/14 \end{bmatrix},
> $$
> so the least-squares line is $y = \tfrac{2}{7} + \tfrac{5}{14}x$.

A common preprocessing step subtracts the mean $\bar{x}$ from every $x_j$, forming
$x^* = x - \bar{x}$. In this **mean-deviation form** the two columns of $X$ become
orthogonal, $X^\top X$ is diagonal, and the normal equations decouple into two
independent one-line solves.

> **Worked example.** Center the four points at $\bar{x} = \tfrac{22}{4} = 5.5$,
> giving deviations $(-3.5, -0.5, 1.5, 2.5)$. The centered column is orthogonal to
> the constant column, so each weight is a single ratio:
> $$
> \beta_0 = \bar{y} = \tfrac{9}{4},
> \qquad
> \beta_1 = \frac{x^*\cdot\mathbf{y}}{x^*\cdot x^*} = \frac{7.5}{21} = \frac{5}{14}.
> $$
> The fitted line $y = \tfrac{9}{4} + \tfrac{5}{14}(x - 5.5)$ matches the earlier
> $y = \tfrac{2}{7} + \tfrac{5}{14}x$, now with a diagonal $X^\top X$.

## The general linear model

Statisticians write the fit with an explicit **residual vector**
$\boldsymbol{\epsilon} = \mathbf{y} - X\boldsymbol{\beta}$, so
$$
\mathbf{y} = X\boldsymbol{\beta} + \boldsymbol{\epsilon}.
$$
Any equation of this form is a **linear model**, and the least-squares fit
minimizes $\lVert \boldsymbol{\epsilon} \rVert$. The word "linear" refers to the
_parameters_ $\beta_i$, not the shape of the fitted curve. A model is linear as
long as it has the form
$$
y = \beta_0 f_0(x) + \beta_1 f_1(x) + \cdots + \beta_k f_k(x)
$$
for known functions $f_0, \dots, f_k$. Changing the $f_i$ changes the columns of
$X$ but leaves the normal equations unchanged.

### Polynomial fits

If the data curve rather than run straight, fit a polynomial. For a parabola
$y = \beta_0 + \beta_1 x + \beta_2 x^2$, each point contributes an equation
$y_j = \beta_0 + \beta_1 x_j + \beta_2 x_j^2 + \epsilon_j$, and the design matrix
gains a squared column:
$$
X = \begin{bmatrix}
  1 & x_1 & x_1^2 \\
  1 & x_2 & x_2^2 \\
  \vdots & \vdots & \vdots \\
  1 & x_n & x_n^2
\end{bmatrix},
\qquad
\boldsymbol{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \beta_2 \end{bmatrix}.
$$
A cubic $y = \beta_0 + \beta_1 x + \beta_2 x^2 + \beta_3 x^3$ adds a fourth column
of cubes. The columns are the sampled basis functions $1, x, x^2, x^3$; the model
stays linear in $\boldsymbol{\beta}$ because those functions are known once the
$x_j$ are fixed.

$$
% caption: The design matrix for a polynomial fit: each row samples the basis
% functions 1, x, x-squared at one data point, so column k holds the k-th power
% of every observation.
\begin{tikzpicture}[>=stealth, scale=1.0, font=\footnotesize]
\definecolor{acc}{HTML}{4A6FA5}
\draw[black, thick] (0,0) rectangle (4.2,2.6);
\foreach \x in {1.4,2.8} \draw[acc!55] (\x,0) -- (\x,2.6);
\node[acc] at (0.7,2.95) {$1$};
\node[acc] at (2.1,2.95) {$x$};
\node[acc] at (3.5,2.95) {$x^2$};
\node[black] at (0.7,2.2) {$1$}; \node[black] at (2.1,2.2) {$x_1$}; \node[black] at (3.5,2.2) {$x_1^2$};
\node[black] at (0.7,1.5) {$1$}; \node[black] at (2.1,1.5) {$x_2$}; \node[black] at (3.5,1.5) {$x_2^2$};
\node[black] at (0.7,0.8) {$\vdots$}; \node[black] at (2.1,0.8) {$\vdots$}; \node[black] at (3.5,0.8) {$\vdots$};
\node[black] at (0.7,0.3) {$1$}; \node[black] at (2.1,0.3) {$x_n$}; \node[black] at (3.5,0.3) {$x_n^2$};
\node at (4.6,1.3) {$=X$};
\end{tikzpicture}
$$

## Multiple regression

With two predictors $u$ and $v$ and a response $y$, a plane
$y = \beta_0 + \beta_1 u + \beta_2 v$ fits data $(u_1, v_1, y_1), \dots, (u_n, v_n,
y_n)$. Each observation gives $y_j = \beta_0 + \beta_1 u_j + \beta_2 v_j +
\epsilon_j$, and the model is $\mathbf{y} = X\boldsymbol{\beta} +
\boldsymbol{\epsilon}$ with
$$
X = \begin{bmatrix}
  1 & u_1 & v_1 \\
  1 & u_2 & v_2 \\
  \vdots & \vdots & \vdots \\
  1 & u_n & v_n
\end{bmatrix},
\qquad
\boldsymbol{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \beta_2 \end{bmatrix}.
$$
The least-squares solution is the **least-squares plane**. Cross terms and higher
powers keep the model linear: an equation like $y = \beta_0 + \beta_1 u + \beta_2 v +
\beta_3 u^2 + \beta_4 uv + \beta_5 v^2$ (a **trend surface**) is still linear in
$\boldsymbol{\beta}$, with columns given by the sampled functions $1, u, v, u^2,
uv, v^2$. Whatever the predictors, the fit obeys the same normal equations, and
whenever $X^\top X$ is invertible,
$$
\hat{\boldsymbol{\beta}} = (X^\top X)^{-1} X^\top \mathbf{y}.
$$

$$
% caption: A least-squares plane fits scattered points in space: the vertical
% residuals from each data point to the plane are minimized in the sum-of-squares
% sense.
\begin{tikzpicture}[>=stealth, scale=1.0, font=\small]
\definecolor{acc}{HTML}{4A6FA5}
% plane
\fill[acc!10] (-0.4,0.2) -- (4.4,-0.4) -- (5.6,1.4) -- (0.8,2.0) -- cycle;
\draw[acc, thick] (-0.4,0.2) -- (4.4,-0.4) -- (5.6,1.4) -- (0.8,2.0) -- cycle;
\node[acc] at (5.0,1.55) {least-squares plane};
% points above/below with residual stems
\foreach \x/\y/\fy in {1.4/1.4/0.95, 2.6/0.35/0.75, 3.6/1.7/0.72, 2.0/1.95/1.15} {
  \draw[black, thick] (\x,\fy) -- (\x,\y);
  \fill[black] (\x,\y) circle (1.9pt);
}
\end{tikzpicture}
$$

### Splitting the sum of squares

The projection geometry gives a bookkeeping identity used throughout regression.
The residual vector $\boldsymbol{\epsilon} = \mathbf{y} - X\hat{\boldsymbol{\beta}}$
is orthogonal to $\operatorname{Col} X$, and $X\hat{\boldsymbol{\beta}}$ lies in
it, so the Pythagorean theorem splits the total variation of $\mathbf{y}$:
$$
\underbrace{\lVert \mathbf{y} \rVert^2}_{\text{SS(T)}}
  = \underbrace{\lVert X\hat{\boldsymbol{\beta}} \rVert^2}_{\text{SS(R)}}
  + \underbrace{\lVert \mathbf{y} - X\hat{\boldsymbol{\beta}} \rVert^2}_{\text{SS(E)}}.
$$
The total sum of squares (SS(T)) divides into the part explained by the fit
(SS(R), the regression sum of squares) and the part left over (SS(E), the error
sum of squares). This is the orthogonal decomposition of $\mathbf{y}$ into
$\operatorname{Col} X$ and its complement, read as squared lengths, and it
underlies the analysis of variance.

$$
% caption: The Pythagorean split of the observation vector: y decomposes into the
% fitted part in Col X and the residual orthogonal to it, so the total sum of
% squares equals the regression plus error sums of squares.
\begin{tikzpicture}[>=stealth, scale=1.0, font=\small]
\definecolor{acc}{HTML}{4A6FA5}
\draw[acc, thick] (-0.5,0) -- (5.2,0);
\node[acc] at (4.8,-0.35) {$\operatorname{Col} X$};
\coordinate (o) at (0.4,0);
\coordinate (f) at (3.4,0);
\coordinate (y) at (3.4,2.2);
\draw[->, black, very thick] (o) -- (y) node[above right] {$\mathbf{y}$};
\draw[->, acc, very thick] (o) -- (f);
\draw[black, thick] (f) -- (y);
\draw[black] (3.4,0.24) -- (3.16,0.24) -- (3.16,0);
\node[black] at (1.5,1.4) {SS(T)};
\node[acc] at (1.9,-0.35) {SS(R)};
\node[black] at (3.75,1.1) {SS(E)};
\fill[black] (o) circle (1.6pt);
\end{tikzpicture}
$$

## The design-matrix pattern

The examples differ only in how the columns of $X$ are formed. A line samples
$\{1, x\}$; a parabola samples $\{1, x, x^2\}$; a plane samples $\{1, u, v\}$.
Once $X$ is built, the parameter vector solves the same normal equations, and
linear algebra handles all of them at once.

| Model | Fitted form | Columns of $X$ |
| --- | --- | --- |
| Least-squares line | $\beta_0 + \beta_1 x$ | $1,\ x$ |
| Parabola | $\beta_0 + \beta_1 x + \beta_2 x^2$ | $1,\ x,\ x^2$ |
| Cubic | $\beta_0 + \beta_1 x + \beta_2 x^2 + \beta_3 x^3$ | $1,\ x,\ x^2,\ x^3$ |
| Multiple regression | $\beta_0 + \beta_1 u + \beta_2 v$ | $1,\ u,\ v$ |
| Trend surface | $\beta_0 + \beta_1 u + \cdots + \beta_5 v^2$ | $1, u, v, u^2, uv, v^2$ |

Building a design matrix and solving the normal equations covers straight-line
regression, polynomial curve fitting, and multiple regression with one procedure.
Replacing the dot product by an abstract
[inner product](/linear-algebra/orthogonality-least-squares/inner-product-spaces)
extends the same least-squares idea to approximating functions, not just data
points.[^lay-66]

[^lay-66]: Lay, §6.6 — Applications to Linear Models: the least-squares line and regression coefficients, the general linear model $\mathbf{y} = X\boldsymbol{\beta} + \boldsymbol{\epsilon}$, polynomial curve fitting, and multiple regression with the design matrix $X$ and normal equations $X^\top X\boldsymbol{\beta} = X^\top\mathbf{y}$.