Linear least squares

Linear least squares#

Linear least squares refers to the problem of minimizing $χ^{2}$ with a model that can be written as a linear combination of basis functions $f_{j} (x)$ as

m (x) = \sum_{j = 1}^{M} a_{j} f_{j} (x)

where $a_{j}$ are the model parameters. For example if we wanted to fit a polynomial we would choose $f_{j} (x) = x^{j - 1}$ and minimize $χ^{2}$ to find the polynomial coefficients ${a_{j}}$ .

We’ll follow the analysis of Numerical Recipes with some small changes in notation. Write out $χ^{2}$ and minimize it:

χ^{2} = \sum_{i = 1}^{N} {[\frac{y_{i} - \sum_{j = 1}^{M} a_{j} f_{j} (x_{i})}{σ_{i}}]}^{2}

\frac{\partial χ^{2}}{\partial a_{k}} = 0 \Rightarrow \sum_{i = 1}^{N} \frac{1}{σ_{i}^{2}} [y_{i} - \sum_{j = 1}^{M} a_{j} f_{j} (x_{i})] f_{k} (x_{i}) = 0

\sum_{i = 1}^{N} \frac{y_{i}}{σ_{i}} \frac{f_{k} (x_{i})}{σ_{i}} - \sum_{i = 1}^{N} \sum_{j = 1}^{M} a_{k} \frac{f_{j} (x_{i})}{σ_{i}} \frac{f_{k} (x_{i})}{σ_{i}} = 0

At this point it is useful to define the design matrix

A_{i j} = \frac{f_{j} (x_{i})}{σ_{i}}

and also the vector $\vec{d}$ which is the data weighted by the errors at each point:

d_{i} = \frac{y_{i}}{σ_{i}} .

Then we have

\sum_{i = 1}^{N} A_{i k} b_{i} - \sum_{i = 1}^{N} \sum_{j = 1}^{M} a_{k} A_{i j} A_{i k} = 0

\sum_{i = 1}^{N} (A^{T})_{k i} b_{i} - \sum_{j = 1}^{M} \sum_{i = 1}^{N} (A^{T})_{j i} A_{i k} a_{k} = 0

So in matrix form, we have

A^{T} d = A^{T} A a .

These are known as the normal equations for the linear least squares problem.

We’ve assumed here that the noise is uncorrelated between measurements (the noise on one measurement does not depend on the noise on a previous measurement). The normal equations can be generalized to include correlated noise, but we’ll leave this for later and for now focus on how we can solve these equations.

Exercise:

Show that $χ^{2}$ can be written in our matrix notation as

χ^{2} = (d - A a)^{T} (d - A a) .