Least squares, projections and the SVD

Worked solutions

Problem 1

Here $a, b \in \mathbb{R}^n$ and $a \neq 0$.

1. $p = ta$ for some number $t$.$p$ lies on the line through $a$, and every point of that line is a multiple of $a$.
2. $a^\top(b - ta) = 0$.The closest point is the foot of the perpendicular: the error $b - p$ must be orthogonal to the line. If it were not, sliding $p$ along the line would shorten it.
3. $a^\top b - t\,a^\top a = 0$, so $t = \dfrac{a^\top b}{a^\top a}$.Distribute $a^\top$; $a^\top a = \|a\|^2 > 0$, so it can be divided by.
4. $p = \dfrac{a^\top b}{a^\top a}\,a = \dfrac{a\,(a^\top b)}{a^\top a} = \dfrac{aa^\top}{a^\top a}\,b$.$a^\top b$ is a number, so it can move to the other side of $a$; then regroup $a(a^\top b) = (aa^\top)b$.
5. $p = \dfrac{a^\top b}{a^\top a}a$, $P = \dfrac{aa^\top}{a^\top a}$$P$ is $n\times n$: a column times a row, divided by a number. Check: $a = (1, 0)^\top$ gives $P = \begin{bmatrix}1&0\\0&0\end{bmatrix}$, which keeps the first coordinate and zeroes the second.

Problem 2

1. $P^2 = \dfrac{aa^\top aa^\top}{(a^\top a)^2} = \dfrac{a\,(a^\top a)\,a^\top}{(a^\top a)^2} = \dfrac{aa^\top}{a^\top a} = P$.The middle $a^\top a$ is a number and cancels one factor of the denominator. Geometrically, projecting a point already on the line leaves it where it is.
2. $P^\top = \dfrac{(aa^\top)^\top}{a^\top a} = \dfrac{(a^\top)^\top a^\top}{a^\top a} = \dfrac{aa^\top}{a^\top a} = P$.$(XY)^\top = Y^\top X^\top$, and $(a^\top)^\top = a$.
3. Column $j$ of $P$ is $\dfrac{a_j}{a^\top a}\,a$.Column $j$ of $aa^\top$ is $a$ times the $j$-th entry of $a^\top$. Every column is a multiple of $a \neq 0$, so the column space is the line through $a$ and the rank is $1$.
4. $Pa = \dfrac{a\,(a^\top a)}{a^\top a} = a$.So $a$ is an eigenvector with eigenvalue $1$.
5. For $x \perp a$: $Px = \dfrac{a\,(a^\top x)}{a^\top a} = 0 = 0\cdot x$.$a^\top x = 0$ is what $x \perp a$ means. So every nonzero $x$ orthogonal to $a$ is an eigenvector with eigenvalue $0$.
6. The vectors orthogonal to $a$ form a subspace of dimension $n - 1$, so $\lambda = 0$ has $n - 1$ independent eigenvectors; with $a$ that makes $n$, so there are no other eigenvalues.They are the solutions of one nonzero equation $a^\top x = 0$ in $n$ unknowns. $n$ independent eigenvectors account for all $n$ eigenvalues, counted with multiplicity.
7. $P^2 = P$ and $P^\top = P$ by direct multiplication; rank $1$ (every column is a multiple of $a$); $Pa = a$ so $\lambda = 1$ with eigenvector $a$, and $Px = 0$ for every $x\perp a$, so $\lambda = 0$ with multiplicity $n-1$In $\mathbb{R}^2$ with $a = (1,1)^\top$: $P = \tfrac12\begin{bmatrix}1&1\\1&1\end{bmatrix}$, $P(1,1)^\top = (1,1)^\top$ and $P(1,-1)^\top = 0$. $P$ is symmetric, so its eigenvectors for $1$ and $0$ are orthogonal, as on the previous page.

Problem 3

1. $Pb = A\hat x$ for some $\hat x \in \mathbb{R}^k$.$Pb$ lies in $\operatorname{col}(A)$, and every point there is $A$ times some vector.
2. $b - A\hat x \perp \operatorname{col}(A)$, that is $A^\top(b - A\hat x) = 0$.As in Problem 1, the closest point is where the error is perpendicular to the whole subspace. Orthogonal to every column of $A$ means every entry of $A^\top(b - A\hat x)$, one column dotted with the error per entry, is $0$.
3. $A^\top A\hat x = A^\top b$, so $\hat x = (A^\top A)^{-1}A^\top b$.These are the normal equations, the same ones the matrix-calculus page reached by setting the gradient to $0$. $A^\top A$ is invertible because the columns of $A$ are independent.
4. $Pb = A\hat x = A(A^\top A)^{-1}A^\top b$ for every $b$, so $P = A(A^\top A)^{-1}A^\top$.Read the matrix off. With one column $a$, $A^\top A = a^\top a$ is a number and this is Problem 1's $P$.
5. $P^2 = A(A^\top A)^{-1}(A^\top A)(A^\top A)^{-1}A^\top = A(A^\top A)^{-1}A^\top = P$.The inner $(A^\top A)(A^\top A)^{-1}$ is $I$.
6. $P^\top = A\big((A^\top A)^{-1}\big)^\top A^\top = A(A^\top A)^{-1}A^\top = P$.Transpose reverses the product; the transpose of an inverse is the inverse of the transpose, and $(A^\top A)^\top = A^\top A$.
7. $P = A(A^\top A)^{-1}A^\top$; $A^\top(b - Pb) = A^\top b - A^\top b = 0$$A^\top Pb = (A^\top A)(A^\top A)^{-1}A^\top b = A^\top b$. Each entry of $A^\top(b - Pb)$ is one column of $A$ dotted with the residual, so all of them are $0$.

Problem 4

Each point $(x, y)$ asks for $c + m\,x = y$. The three equations, in the unknowns $c$ and $m$, are $A\,(c, m)^\top = b$ with $A = \begin{bmatrix}1&0\\1&1\\1&2\end{bmatrix}$ and $b = (1, 2, 2)^\top$.

1. Row $i$ of $A$ is $(1, x_i)$ and entry $i$ of $b$ is $y_i$.$c$ is multiplied by $1$ in every equation and $m$ by the point's $x$. Three equations, two unknowns, and the points are not on one line, so there is no exact solution.
2. $A^\top A = \begin{bmatrix}1+1+1&0+1+2\\0+1+2&0+1+4\end{bmatrix} = \begin{bmatrix}3&3\\3&5\end{bmatrix}$.Entry $(i,j)$ is column $i$ of $A$ dotted with column $j$.
3. $A^\top b = (1+2+2,\ 0+2+4)^\top = (5, 6)^\top$.Each column of $A$ dotted with $b$.
4. $3c + 3m = 5$ and $3c + 5m = 6$.The normal equations $A^\top A\,(c, m)^\top = A^\top b$ written out.
5. $2m = 1$, so $m = \tfrac12$; then $3c = 5 - \tfrac32 = \tfrac72$, so $c = \tfrac76$.Subtract the first equation from the second to remove $c$, then substitute back.
6. $A = \begin{bmatrix}1&0\\1&1\\1&2\end{bmatrix}$, $b = (1,2,2)^\top$; $A^\top A = \begin{bmatrix}3&3\\3&5\end{bmatrix}$, $A^\top b = (5,6)^\top$; $c = \tfrac76$, $m = \tfrac12$Check in the first equation: $3\cdot\tfrac76 + 3\cdot\tfrac12 = \tfrac72 + \tfrac32 = 5$. The fitted line $y = \tfrac76 + \tfrac12x$ passes at heights $\tfrac76$, $\tfrac53$, $\tfrac{13}6$ over $x = 0, 1, 2$.

Problem 5

1. $\hat x = (\tfrac76, \tfrac12)^\top$ and $A\hat x = (\tfrac76,\ \tfrac76 + \tfrac12,\ \tfrac76 + 1)^\top = (\tfrac76, \tfrac53, \tfrac{13}6)^\top$.Row $i$ of $A$ is $(1, x_i)$, so $A\hat x$ lists the line's heights at the three $x$ values.
2. $r = (1 - \tfrac76,\ 2 - \tfrac53,\ 2 - \tfrac{13}6)^\top = (-\tfrac16, \tfrac13, -\tfrac16)^\top$.Subtract entry by entry; each entry is how far the point is above the line.
3. First entry of $A^\top r$: $-\tfrac16 + \tfrac13 - \tfrac16 = 0$.The first column of $A$ is all ones, so this is the sum of the residuals: a least-squares line with an intercept always balances its errors.
4. Second entry: $0\cdot(-\tfrac16) + 1\cdot\tfrac13 + 2\cdot(-\tfrac16) = \tfrac13 - \tfrac13 = 0$.The second column of $A$ is the $x$ values, so the residuals are also uncorrelated with $x$.
5. $r = (-\tfrac16, \tfrac13, -\tfrac16)^\top$ and $A^\top r = (0, 0)^\top$This is Problem 3's orthogonality on actual numbers: $r$ is orthogonal to both columns, so to all of $\operatorname{col}(A)$. The squared error is $\|r\|^2 = \tfrac1{36} + \tfrac19 + \tfrac1{36} = \tfrac16$, and no other line does better.

Problem 6

1. $A^\top A = \begin{bmatrix}25&20\\20&25\end{bmatrix}$, with eigenvalues $45$ and $5$ and eigenvectors $(1,1)^\top$ and $(1,-1)^\top$.The previous page computed all of this (its Problem 10). $A = U\Sigma V^\top$ gives $A^\top A = V\Sigma U^\top U\Sigma V^\top = V\Sigma^2V^\top$, so the eigenvectors of $A^\top A$ are the $v_i$ and its eigenvalues are the $\sigma_i^2$.
2. $\sigma_1 = \sqrt{45} = 3\sqrt5$ and $\sigma_2 = \sqrt5$.Square roots of the eigenvalues, largest first.
3. $v_1 = \tfrac1{\sqrt2}(1,1)^\top$, $v_2 = \tfrac1{\sqrt2}(1,-1)^\top$, so $V = \tfrac1{\sqrt2}\begin{bmatrix}1&1\\1&-1\end{bmatrix}$.The columns of $V$ must be unit vectors, so divide each eigenvector by its length $\sqrt2$, in the order of the singular values.
4. $Av_1 = \tfrac1{\sqrt2}(3,\ 4 + 5)^\top = \tfrac1{\sqrt2}(3, 9)^\top$, and $u_1 = \dfrac{Av_1}{\sigma_1} = \dfrac{(3, 9)^\top}{\sqrt2\cdot3\sqrt5} = \tfrac1{\sqrt{10}}(1, 3)^\top$.$Av_i = \sigma_iu_i$ (Before you start), so $u_i = Av_i/\sigma_i$. Computing $U$ this way keeps its signs matched to $V$.
5. $Av_2 = \tfrac1{\sqrt2}(3,\ 4 - 5)^\top = \tfrac1{\sqrt2}(3, -1)^\top$, and $u_2 = \dfrac{(3, -1)^\top}{\sqrt2\cdot\sqrt5} = \tfrac1{\sqrt{10}}(3, -1)^\top$.Same rule with $\sigma_2 = \sqrt5$.
6. $u_1^\top u_2 = \tfrac1{10}(3 - 3) = 0$ and $u_1^\top u_1 = u_2^\top u_2 = \tfrac{1 + 9}{10} = 1$.They come out orthonormal automatically: $u_i^\top u_j = v_i^\top A^\top Av_j/(\sigma_i\sigma_j) = \sigma_j^2\,v_i^\top v_j/(\sigma_i\sigma_j)$, which is $0$ for $i \neq j$ and $1$ for $i = j$.
7. $U\Sigma = \tfrac1{\sqrt{10}}\begin{bmatrix}3\sqrt5&3\sqrt5\\9\sqrt5&-\sqrt5\end{bmatrix} = \tfrac1{\sqrt2}\begin{bmatrix}3&3\\9&-1\end{bmatrix}$.Multiplying by a diagonal matrix on the right scales the columns of $U$ by $3\sqrt5$ and $\sqrt5$; $\sqrt5/\sqrt{10} = 1/\sqrt2$.
8. $U\Sigma V^\top = \tfrac12\begin{bmatrix}3&3\\9&-1\end{bmatrix}\begin{bmatrix}1&1\\1&-1\end{bmatrix} = \tfrac12\begin{bmatrix}6&0\\8&10\end{bmatrix} = \begin{bmatrix}3&0\\4&5\end{bmatrix}$.$V^\top = V$ here, and the two factors of $\tfrac1{\sqrt2}$ make $\tfrac12$. Multiplying back is the check that nothing slipped.
9. $\sigma_1 = 3\sqrt5$, $\sigma_2 = \sqrt5$; $V = \tfrac1{\sqrt2}\begin{bmatrix}1&1\\1&-1\end{bmatrix}$; $u_i = Av_i/\sigma_i$ gives $U = \tfrac1{\sqrt{10}}\begin{bmatrix}1&3\\3&-1\end{bmatrix}$; check $U\Sigma V^\top = A$When the singular values are distinct and nonzero, as here, the SVD is unique up to signs: replacing $v_i$ by $-v_i$ flips $u_i$ too, and the product is unchanged.

Problem 7

1. $A = U\Sigma V^\top = \sigma_1u_1v_1^\top + \sigma_2u_2v_2^\top = \sigma_1u_1v_1^\top$.Multiply out column by row: $U\Sigma$ has columns $\sigma_iu_i$, and a product of matrices is the sum of each column of the first times the matching row of the second. The $\sigma_2$ term is $0$.
2. $Ax = \sigma_1u_1\,(v_1^\top x)$ for every $x \in \mathbb{R}^2$.$v_1^\top x$ is a number, so every output is a multiple of $u_1$.
3. $\operatorname{col}(A) = \operatorname{span}\{u_1\}$, so the rank is $1$.Step 2 puts every $Ax$ on the line through $u_1$, and $x = v_1$ gives $Av_1 = \sigma_1u_1 \neq 0$, so the whole line is reached. The rank is the dimension of the column space.
4. $Ax = 0$ exactly when $v_1^\top x = 0$, that is when $x$ is a multiple of $v_2$.$\sigma_1u_1 \neq 0$, so the product in step 2 is $0$ only when the number $v_1^\top x$ is. In $\mathbb{R}^2$ the vectors orthogonal to the unit vector $v_1$ are the multiples of $v_2$, since $V$ is orthogonal.
5. rank $1$; $\operatorname{col}(A) = \operatorname{span}\{u_1\}$; null space $=\operatorname{span}\{v_2\}$, because $Av_2 = \sigma_2u_2 = 0$In general the rank is the number of nonzero singular values, the $u_i$ that go with them span the column space, and the $v_i$ that go with zero singular values span the null space. The rank count agrees: $2$ columns $= 1$ (rank) $+ 1$ (null space dimension).

Problem 8

1. Independent columns make every $\sigma_i > 0$, so $\Sigma$ is invertible and $\Sigma^+ = \Sigma^{-1}$.$\sigma_i = \|Av_i\|$, because $Av_i = \sigma_iu_i$ with $\|u_i\| = 1$; and $Av_i \neq 0$ since $v_i \neq 0$ and only $x = 0$ gives $Ax = 0$ when the columns are independent.
2. $A^\top A = V\Sigma U^\top U\Sigma V^\top = V\Sigma^2V^\top$.Substitute $A = U\Sigma V^\top$ and $A^\top = V\Sigma U^\top$ ($\Sigma$ is diagonal, so $\Sigma^\top = \Sigma$); $U^\top U = I$ removes the middle.
3. $(A^\top A)^{-1} = V\Sigma^{-2}V^\top$.Multiply to check: $V\Sigma^2V^\top\,V\Sigma^{-2}V^\top = V\Sigma^2\Sigma^{-2}V^\top = VV^\top = I$, using $V^\top V = I$ and then $VV^\top = I$.
4. $(A^\top A)^{-1}A^\top = V\Sigma^{-2}V^\top\,V\Sigma U^\top = V\Sigma^{-2}\Sigma U^\top = V\Sigma^{-1}U^\top$.$V^\top V = I$ again, and diagonal matrices multiply entry by entry: $\sigma_i^{-2}\sigma_i = \sigma_i^{-1}$.
5. For Problem 4, $(A^\top A)^{-1} = \tfrac16\begin{bmatrix}5&-3\\-3&3\end{bmatrix}$.The $2\times2$ inverse of $\begin{bmatrix}3&3\\3&5\end{bmatrix}$: determinant $15 - 9 = 6$, swap the diagonal, negate the off-diagonal.
6. $A^+b = (A^\top A)^{-1}A^\top b = \tfrac16\begin{bmatrix}5&-3\\-3&3\end{bmatrix}\begin{bmatrix}5\\6\end{bmatrix} = \tfrac16(25 - 18,\ -15 + 18)^\top = (\tfrac76, \tfrac12)^\top$.Step 4 lets us use $A^\top b = (5, 6)^\top$ from Problem 4 instead of computing an SVD of that $A$.
7. $(A^\top A)^{-1}A^\top = V\Sigma^{-1}U^\top = A^+$; $A^+b = (\tfrac76, \tfrac12)^\top$The same $c$ and $m$ as Problem 4: when the columns are independent, $A^+b$ is the least-squares solution $\hat x$. $A^+$ is still defined when they are not, which is why robust solvers such as NumPy's lstsq go through the SVD; there $A^+b$ is the least-squares solution of smallest norm.

Problem 9

1. The Eckart–Young theorem: among all matrices of rank at most $1$, $A_1 = \sigma_1u_1v_1^\top$ is closest to $A$ in $\|\cdot\|_2$, and $\|A - A_1\|_2 = \sigma_2$.Quoted here, not proved (the spectral-norm form is Mirsky's). In $\|\cdot\|_2$ it is a closest matrix, not the only one; in the Frobenius norm it is the unique closest when $\sigma_1 > \sigma_2$. The same holds for rank $j$, with $A_j = \sum_{i \le j}\sigma_iu_iv_i^\top$ and error $\sigma_{j+1}$.
2. $A_1 = 3\sqrt5\cdot\tfrac1{\sqrt{10}}(1,3)^\top\cdot\tfrac1{\sqrt2}(1,1) = \tfrac{3\sqrt5}{\sqrt{20}}\begin{bmatrix}1&1\\3&3\end{bmatrix} = \tfrac32\begin{bmatrix}1&1\\3&3\end{bmatrix}$.$u_1$, $v_1$ and $\sigma_1$ from Problem 6; $v_1^\top = \tfrac1{\sqrt2}(1,1)$ is a row, and a column times a row gives the matrix of products. $\sqrt{10}\sqrt2 = \sqrt{20} = 2\sqrt5$.
3. $A - A_1 = \begin{bmatrix}3 - \tfrac32&0 - \tfrac32\\4 - \tfrac92&5 - \tfrac92\end{bmatrix} = \tfrac12\begin{bmatrix}3&-3\\-1&1\end{bmatrix}$.Subtract entry by entry.
4. $\tfrac12\begin{bmatrix}3&-3\\-1&1\end{bmatrix} = \sqrt5\cdot\tfrac1{\sqrt{10}}(3,-1)^\top\cdot\tfrac1{\sqrt2}(1,-1) = \sigma_2u_2v_2^\top$.$A = \sigma_1u_1v_1^\top + \sigma_2u_2v_2^\top$ (Problem 7, step 1), so the error is exactly the dropped term; $\sqrt5/\sqrt{20} = \tfrac12$.
5. $\|\sigma_2u_2v_2^\top x\| = \sigma_2\,\lvert v_2^\top x\rvert \le \sigma_2$ for unit $x$, with equality at $x = v_2$.$\|u_2\| = 1$, and $\lvert v_2^\top x\rvert \le \|v_2\|\,\|x\| = 1$ with equality when $x$ points along $v_2$; so the largest stretch is $\sigma_2$.
6. $A_1 = \sigma_1u_1v_1^\top$ and $\|A - A_1\|_2 = \sigma_2 = \sqrt5$Here $A_1 = \tfrac32\begin{bmatrix}1&1\\3&3\end{bmatrix}$. The error is the largest singular value thrown away, so a matrix whose singular values fall off fast is well approximated by few terms: the idea behind low-rank compression and PCA.

Problem 10

$P = A(A^\top A)^{-1}A^\top$ is the projection of Problem 3, with $P^2 = P$ and $P^\top = P$.

1. $(I-P)^2 = I - 2P + P^2 = I - 2P + P = I - P$, and $(I - P)^\top = I - P$.Expand, then $P^2 = P$; the transpose of $I$ is $I$ and $P^\top = P$. So $I - P$ is a projection matrix of the same kind as $P$.
2. $A^\top(I-P) = A^\top - A^\top A(A^\top A)^{-1}A^\top = A^\top - A^\top = 0$.So every $(I-P)b$ is orthogonal to every column of $A$: $I - P$ lands in $\operatorname{col}(A)^\perp$.
3. For $w \in \operatorname{col}(A)^\perp$: $Pw = A(A^\top A)^{-1}(A^\top w) = 0$, so $(I-P)w = w$.$A^\top w = 0$ is what orthogonal to every column means. $I - P$ leaves that subspace fixed, so with step 2 its range is exactly $\operatorname{col}(A)^\perp$. And $b - (I-P)b = Pb$ lies in $\operatorname{col}(A)$, orthogonal to that subspace, so $(I-P)b$ is the foot of the perpendicular from $b$, as in Problem 3; it is the residual $b - Pb$.
4. $b = Pb + (I-P)b$, so $\|b\|^2 = \|Pb\|^2 + 2\,(Pb)^\top(I-P)b + \|(I-P)b\|^2$.Expand $\|w + z\|^2 = (w+z)^\top(w+z) = \|w\|^2 + 2w^\top z + \|z\|^2$.
5. $(Pb)^\top(I-P)b = b^\top P^\top(I-P)b = b^\top(P - P^2)b = 0$.$P^\top = P$ and $P^2 = P$, so $P^\top(I - P) = P - P^2$ is the zero matrix. The cross term vanishes.
6. $I - P$ is the projection onto $\operatorname{col}(A)^\perp$, and $\|b\|^2 = \|Pb\|^2 + \|(I-P)b\|^2$Pythagoras: $b$ splits into perpendicular pieces. For Problem 4, $\|b\|^2 = 9$, $\|r\|^2 = \tfrac16$ (Problem 5) and $\|A\hat x\|^2 = \tfrac{49}{36} + \tfrac{100}{36} + \tfrac{169}{36} = \tfrac{53}6$; indeed $\tfrac{53}6 + \tfrac16 = 9$.