NLA MT25, Least-squares
Flashcards
Via QR factorisation
See in particular ∆least-squares-with-qr.
@State a result that characterises the backward stability of finding the least squares solution to an overdetermined linear system $Ax \approx b$ via the QR decomposition.
The QR-based least-squares algorithm (i.e. compute the thin QR $A = QR$ via Householder reflectors, then solve $Rx = Q^\top b$) is backward stable. The computed $\hat x$ is the exact solution of the perturbed problem:
\[\hat x = \text{argmin} _ x \|(A + \Delta A)x - (b + \Delta b)\| _ 2\]for some $\Delta A$, $\Delta b$ with $\|\Delta A\| = O(\epsilon \|A\|)$ and $\|\Delta b\| = O(\epsilon \|b\|) = O(\epsilon)$.
Via normal equations
@Prove that
\[\hat x = \text{argmin} _ x \vert \vert Ax - b\| _ 2 \iff (A^\top A) \hat x = A^\top b\]
so that in particular we may solve the least-squares problem by solving the normal equations.
We minimise the squared objective:
\[\begin{aligned} \vert \vert Ax-b \vert \vert _ 2^2 &= (Ax - b)^\top (Ax -b) \\ &= (x^\top A^\top - b^\top)(Ax - b) \\ &= x^\top A^\top A x - x^\top A^\top b - b^\top Ax + b^\top b \\ \end{aligned}\]Since $x^\top A^\top b = b^\top A x$ (each is a $1 \times 1$ scalar, hence equal to its own transpose), the two cross terms combine to $-2\, b^\top A x$. Differentiating with respect to $x$ and setting to $0$, the minimum occurs when
\[\begin{aligned} &2A^\top A \hat x - 2A^\top b = 0 \\ \iff& A^\top A\hat x = A^\top b \end{aligned}\]What are the pros and cons of solving least squares problems $\min _ x \vert \vert Ax - b \vert \vert _ 2$ using the normal equations $A^\top A x = A^\top b$ and the Cholesky decomposition of $A^\top A$ (which exists, since $A^\top A \succ 0$ when $\text{rank}(A) = n$).
- Pros: Fast compared to e.g. using the QR decomposition.
- Cons: The forward error is $O(u \kappa _ 2(A)^2)$ rather than $O(u \kappa _ 2(A))$ for QR-LS. Forming $A^\top A$ explicitly already squares the condition number, even though the Cholesky solve of the resulting system is itself backward stable.
Via the SVD
@Prove that the full-rank least-squares problem $\min _ x \|Ax - b\| _ 2$ (with $A \in \mathbb R^{m \times n}$, $m \ge n$, $\mathrm{rank}(A) = n$) is solved via the SVD by $\hat x = V \hat\Sigma^{-1} U _ 1^\top b$, and find the minimal residual norm.
Setup: take the full SVD $A = U \Sigma V^\top$ with $U = [U _ 1 \mid U _ 2] \in \mathbb R^{m \times m}$ orthogonal ($U _ 1$ the first $n$ columns), $V \in \mathbb R^{n \times n}$ orthogonal, and $\Sigma = \begin{bmatrix} \hat\Sigma \\ 0 \end{bmatrix}$ where $\hat\Sigma = \mathrm{diag}(\sigma _ 1, \ldots, \sigma _ n) \succ 0$ (invertible since $\mathrm{rank}(A) = n$). Write $U^\top b = \begin{bmatrix} \hat s \\ \tilde s \end{bmatrix}$ with $\hat s = U _ 1^\top b \in \mathbb R^n$ and $\tilde s = U _ 2^\top b \in \mathbb R^{m-n}$.
Minimise the squared residual:
\[\begin{aligned} \min _ x \|Ax - b\| _ 2^2 &= \min _ x \|U \Sigma V^\top x - b\| _ 2^2 && (\star 1) \\ &= \min _ x \|\Sigma V^\top x - U^\top b\| _ 2^2 && (\star 2) \\ &= \min _ {y} \left\| \begin{bmatrix} \hat\Sigma \\ 0 \end{bmatrix} y - \begin{bmatrix} \hat s \\ \tilde s \end{bmatrix} \right\| _ 2^2 && (\star 3) \\ &= \min _ {y} \left( \|\hat\Sigma y - \hat s\| _ 2^2 + \|\tilde s\| _ 2^2 \right) && (\star 4) \\ &= \|\tilde s\| _ 2^2 && (\star 5) \end{aligned}\]attained at $\hat\Sigma y = \hat s$, i.e. $y = \hat\Sigma^{-1} \hat s$. Undoing the substitution,
\[\hat x = V y = V \hat\Sigma^{-1} \hat s = V \hat\Sigma^{-1} U _ 1^\top b,\]with minimal residual $\|A \hat x - b\| _ 2 = \|\tilde s\| _ 2 = \|U _ 2^\top b\| _ 2$.
The steps:
- $(\star 1)$ Substitute the full SVD $A = U \Sigma V^\top$.
- $(\star 2)$ $U$ is orthogonal, so the $2$-norm is invariant under $U^\top$: $\|U \Sigma V^\top x - b\| _ 2 = \|U^\top(U \Sigma V^\top x - b)\| _ 2 = \|\Sigma V^\top x - U^\top b\| _ 2$.
- $(\star 3)$ Change variables $y = V^\top x$; as $V$ is orthogonal this is a bijection of $\mathbb R^n$, so the min over $x$ equals the min over $y$. Insert the block forms of $\Sigma$ and $U^\top b$.
- $(\star 4)$ Since $\Sigma y = \begin{bmatrix} \hat\Sigma y \\ 0 \end{bmatrix}$, the squared norm splits into its top $n$ and bottom $m - n$ rows; the bottom rows give $\|0 - \tilde s\| _ 2^2 = \|\tilde s\| _ 2^2$, free of $y$.
- $(\star 5)$ The first term is $\ge 0$ and is killed by $y = \hat\Sigma^{-1} \hat s$ (valid as $\hat\Sigma$ is invertible); the second is constant.
Remark: this is exactly $\hat x = A^\dagger b$, since $A^\dagger = V \hat\Sigma^{-1} U _ 1^\top$ for full-column-rank $A$ (∆pseudoinverse-least-squares). The residual $\|\tilde s\| _ 2 = \|U _ 2^\top b\| _ 2$ is the component of $b$ orthogonal to $\mathrm{range}(A)$.
Bite-sized
Least-squares setup: given $A \in \mathbb R^{m \times n}$ with $m \ge n$ (typically $m \gg n$) and $b \in \mathbb R^m$, find $x$ minimising $\ \vert Ax - b\ \vert _ 2$. The system is overdetermined, so $Ax = b$ generally has no exact solution.
The QR-based LS solution: with thin QR $A = QR$, $\hat x = $ $R^{-1} Q^\top b$. The cost is one Householder QR ($\approx 2mn^2 - \tfrac{2}{3} n^3$) plus a triangular solve ($O(n^2)$).
Geometric interpretation of the least-squares solution.
$\hat x$ is the unique vector such that the residual $A\hat x - b$ is orthogonal to $\mathrm{range}(A)$ (the column space of $A$). Equivalently, $A \hat x$ is the orthogonal projection of $b$ onto $\mathrm{range}(A)$. The orthogonality condition $A^\top (A \hat x - b) = 0$ is exactly the normal equations.
Solving LS via the SVD: with full SVD $A = U \begin{bmatrix} \hat\Sigma \\ 0 \end{bmatrix} V^\top$ ($\hat\Sigma$ the $n \times n$ diagonal block of singular values), and writing $U^\top b = \begin{bmatrix} \hat s \\ \tilde s \end{bmatrix}$, the solution is $\hat x = $ $V \hat\Sigma^{-1} \hat s$. The residual norm is $\ \vert \tilde s\ \vert _ 2$.
The least-squares solution can be written via the pseudoinverse: $\hat x = $ $A^\dagger b$. For tall full-rank $A$, $A^\dagger = (A^\top A)^{-1} A^\top$.