Examples

Computing least-squares solution - non-orthogonal case

Computing least-squares solution: non-orthogonal case

Problem: Consider the matrix equation $\left(\begin{array}{cc}4& 0\\ 0& 2\\ 1& 1\end{array}\right)\mathbf{x}=\left(\begin{array}{c}2\\ 0\\ 6\end{array}\right).$ There is no exact solution; find the best approximate “least-squares” solution and the resulting error.

Solution: Write $A$ and $\mathbf{b}$ for the matrix and vector. Next calculate $A^{\mathsf{T}}A=\left(\begin{array}{cc}17& 1\\ 1& 5\end{array}\right)$, and $A^{\mathsf{T}}\mathbf{b}=\left(\begin{array}{c}14\\ 6\end{array}\right)$. Now solve the matrix equation

$$A^{\mathsf{T}}A=\left(\begin{array}{cc}17& 1\\ 1& 5\end{array}\right)\left(\begin{array}{c}{x}_1\\ x_2\end{array}\right)=\left(\begin{array}{c}14\\ 6\end{array}\right).$$

We use row-reduction on the augmented matrix:

$$\left(\begin{array}{rrr}17& 1& 14\\ 1& 5& 6\end{array}\right)\sim \left(\begin{array}{rrr}17& 1& 14\\ 0& 84/17& 88/17\end{array}\right)\sim \left(\begin{array}{rrr}17& 1& 14\\ 0& 84& 88\end{array}\right).$$

Therefore $x_2=88/84=22/21$ and $x_1=(14-x_2)/17=16/21$.

To compute the error, first find $A\mathbf{x}$:

$$A\mathbf{x}=\left(\begin{array}{cc}4& 0\\ 0& 2\\ 1& 1\end{array}\right)\left(\begin{array}{c}16/21\\ 22/21\end{array}\right)=\left(\begin{array}{c}64/21\\ 44/21\\ 38/21\end{array}\right)$$

and then subtract from $\mathbf{b}$ and compute the norm:

$$|\mathbf{b}-A\mathbf{x}|=\left|\left(\begin{array}{c}2\\ 0\\ 6\end{array}\right)-\left(\begin{array}{c}64/21\\ 44/21\\ 38/21\end{array}\right)\right|=\left|\left(\begin{array}{c}-22/21\\ -44/21\\ -88/21\end{array}\right)\right|=\frac{1}{21}\sqrt{22^2+44^2+88^2}=\frac{\sqrt{10164}}{21}\approx \mathrm{4.80.}$$

Computing least-squares solution - orthogonal case

Computing least-squares solution: orthogonal case

Problem: Find the least-squares solution for the following matrix equation:

$$A\mathbf{z}=\mathbf{b}:\left(\begin{array}{rr}1& -3\\ 1& 1\\ 1& 0\\ 1& 4\\ 1& -2\end{array}\right)\mathbf{z}=\left(\begin{array}{r}4\\ -1\\ -1\\ -3\\ 1\end{array}\right).$$

Notice that the second column of the matrix and the vector represent the $x$ and $y$ values, respectively, of the data from Exercise 07-02.

Solution: The matrix columns are orthogonal. In this situation it is faster to compute the projection directly:

$${\mathrm{proj}}_{⟨{\mathbf{a}}_1,{\mathbf{a}}_2⟩}(\mathbf{b})=\frac{\mathbf{b}\cdot {\mathbf{a}}_1}{{\mathbf{a}}_1\cdot {\mathbf{a}}_1}{\mathbf{a}}_1+\frac{\mathbf{b}\cdot {\mathbf{a}}_2}{{\mathbf{a}}_2\cdot {\mathbf{a}}_2}{\mathbf{a}}_2=\frac{0}{5}{\mathbf{a}}_1+\frac{-27}{30}{\mathbf{a}}_2=\left(\begin{array}{c}27/10\\ -9/10\\ 0\\ -36/10\\ 18/10\end{array}\right).$$

Furthermore, we have $0/5$ and $-27/30$ for the coefficients of $\mathbf{z}$ that give the best approximate solution, since these are the weights on the columns of $A$ that yield the projection. Therefore $\mathbf{z}=\left(\begin{array}{c}0\\ -9/10\end{array}\right)$ is our solution.

Finally, the error is given by:

$$\left|\left(\begin{array}{r}4\\ -1\\ -1\\ -3\\ 1\end{array}\right)-\left(\begin{array}{c}27/10\\ -9/10\\ 0\\ -36/10\\ 18/10\end{array}\right)\right|=\left|\left(\begin{array}{c}-57/10\\ 19/10\\ 0\\ 76/10\\ -38/10\end{array}\right)\right|=\frac{\sqrt{10830}}{10}\approx \mathrm{10.41.}$$

Finding a best fit cubic

Finding a best fit cubic

Suppose we have some data like this:

Input:	1	2	3	4	5	6
Output:	5	9	8	4	6	8

Let’s fit a cubic polynomial $f(x)=a_0+a_{1}x+a_2{x}^2+a_3{x}^3$ to this data. We use the objects:

$$X=\left(\begin{array}{rrrr}1& 1& 1& 1\\ 1& 2& 4& 8\\ 1& 3& 9& 27\\ 1& 4& 16& 64\\ 1& 5& 25& 125\\ 1& 6& 36& 216\end{array}\right),\mathbf{y}=\left(\begin{array}{c}5\\ 9\\ 8\\ 4\\ 6\\ 8\end{array}\right),\mathbf{a}=\left(\begin{array}{c}{a}_0\\ a_1\\ a_2\\ a_3\end{array}\right).$$

First, compute:

$$X^{\mathsf{T}}X=\left(\begin{array}{cccc}6& 21& 91& 441\\ 21& 91& 441& 2275\\ 91& 441& 2275& 12201\\ 441& 2275& 12201& 67171\end{array}\right),A^{\mathsf{T}}\mathbf{y}=\left(\begin{array}{c}40\\ 141\\ 615\\ 3027\end{array}\right).$$

Row reduction of the augmented matrix $\left(\begin{array}{rr}{X}^{\mathsf{T}}X& A^{\mathsf{T}}\mathbf{y}\end{array}\right)$ yields:

$$\left(\begin{array}{rr}{X}^{\mathsf{T}}X& A^{\mathsf{T}}\mathbf{y}\end{array}\right)\sim \left(\begin{array}{rrrrr}254016& 0& 0& 0& -1354752\\ 0& 254016& 0& 0& 3828048\\ 0& 0& 254016& 0& -1275120\\ 0& 0& 0& 254016& 122304\end{array}\right)⇝\mathbf{a}=\left(\begin{array}{c}-16/3\\ 11393/756\\ -1265/252\\ 13/27\end{array}\right)$$

Therefore our best fit cubic polynomial is:

$$y=\frac{-16}{3}+\frac{11393}{756}x+\frac{-1265}{252}{x}^2+\frac{13}{27}{x}^3.$$

Using Penrose quasi-inverse to solve LLS

Using Penrose quasi-inverse to solve LLS

In an exercise on computing an SVD by hand we found the following SVD:

$$A=\left(\begin{array}{rr}-3& 1\\ 6& -2\\ 6& -2\end{array}\right)=\left(\begin{array}{rcc}1/3& 2/\sqrt{5}& 2/\sqrt{45}\\ -2/3& 1/\sqrt{5}& -4/\sqrt{45}\\ -2/3& 0& 5/\sqrt{45}\end{array}\right)\left(\begin{array}{rc}3\sqrt{10}& 0\\ 0& 0\\ 0& 0\end{array}\right)\left(\begin{array}{rc}-3/\sqrt{10}& 1/\sqrt{10}\\ 1/\sqrt{10}& 3/\sqrt{10}\end{array}\right)=U\Sigma V^{\mathsf{T}}.$$

Problem: Use the SVD above in order to (a) find the Penrose quasi-inverse $A^{+}$, and (b) find a vector $\mathbf{x}$ that minimizes $|A\mathbf{x}-\mathbf{b}|$ for $\mathbf{b}=\left(\begin{array}{c}2\\ 2\\ 4\end{array}\right)$.

Solution: The quasi-inverse is formed by transposing $\Sigma$ and inverting the entries. We obtain for (a):

$${\Sigma }^{+}=\left(\begin{array}{ccc}\frac{1}{3\sqrt{10}}& 0& 0\\ 0& 0& 0\end{array}\right),$$ $$\begin{array}{rl}{A}^{+}=V{\Sigma }^{+}{U}^{\mathsf{T}}& =\left(\begin{array}{rc}-3/\sqrt{10}& 1/\sqrt{10}\\ 1/\sqrt{10}& 3/\sqrt{10}\end{array}\right)\left(\begin{array}{ccc}\frac{1}{3\sqrt{10}}& 0& 0\\ 0& 0& 0\end{array}\right)\left(\begin{array}{rcc}1/3& -2/3& -2/3\\ 2/\sqrt{5}& 1/\sqrt{5}& 0\\ 2/\sqrt{45}& -4/\sqrt{45}& 5/\sqrt{45}\end{array}\right)\\ & \\ & =\left(\begin{array}{rcc}-1/30& 1/15& 1/15\\ 1/90& -1/45& -1/45\end{array}\right).\end{array}$$

Then (b) to solve the LLS, we simply compute $A^{+}\mathbf{b}$ and obtain $\left(\begin{array}{c}1/3\\ -1/9\end{array}\right)$.

Notice that after calculating $A^{+}$ once, we can solve any LLS problem very efficiently. That is a nice feature.