3

I need to solve $Ax=b$ in lots of ways using QR decomposition.

$$A = \begin{bmatrix} 1 & 1 \\ -1 & 1 \\ 1 & 2 \end{bmatrix}, b = \begin{bmatrix} 1 \\ 0 \\ 1 \end{bmatrix}$$

This is an overdetermined system. That is, it has more equations than needed for a unique solution.

I need to find $\min ||Ax-b||$. How should I solve it using QR?

I know that QR can be used to reduce the problem to $$\Vert Ax - b \Vert = \Vert QRx - b \Vert = \Vert Rx - Q^{-1}b \Vert.$$

but what do I do after this?

Poperton
  • 6,594

2 Answers2

8

The most straightforward way I know is to pass through the normal equations:

$$A^T A x = A^T b$$

and substitute in the $QR$ decomposition of $A$ (with the convention $Q \in \mathbb{R}^{m \times n},R \in \mathbb{R}^{n \times n}$). Thus you get

$$R^T Q^T Q R x = R^T Q^T b.$$

But $Q^T Q=I_n$. (Note that in this convention $Q$ isn't an orthogonal matrix, so $Q Q^T \neq I_m$, but this doesn't matter here.) Thus:

$$R^T R x = R^T Q^T b.$$

If $A$ has linearly independent columns (as is usually the case with overdetermined systems), then $R^T$ is injective, so by multiplying both sides by the left inverse of $R^T$ you get

$$Rx=Q^T b.$$

This system is now easy to solve numerically.

For numerical purposes it's important that the removal of $Q^T Q$ and $R^T$ from the problem is done analytically, and in particular $A^T A$ is never constructed numerically.

Ian
  • 104,572
  • 1
    Is there any reason to make this so convoluted? From $Ax=b$ you have $QRx=b$, multiply by $Q^T$ on the left. – Martin Argerami Apr 12 '19 at 15:45
  • @MartinArgerami Because actually the least squares solution usually does not satisfy $Ax=b$. This simple perspective only shows you that this approach gives you a solution when a solution exists. Now you could argue directly that multiplying both sides by $Q^T$ furnishes an equation whose solution is the least squares solution. (Such an argument would resemble the usual geometric argument for deriving the normal equations.) This would make a good alternative answer to mine. – Ian Apr 12 '19 at 15:46
2

Note that $Rx$ has the form $$Rx = \begin{bmatrix} y_1 \\ y_2 \\ 0\end{bmatrix} $$ , so if $$ Q^{-1}b = \begin{bmatrix} z_1 \\ z_2 \\ z_3\end{bmatrix}$$ then $|| Rx - Q^{-1}b||$ will be minimal for $y_1 = z_1$, $y_2=z_2$. This set of equation is no longer overdetermined.

Using matrix notation, if tou write $R = \begin{bmatrix} R_1 \\ 0\end{bmatrix}$ and intoduce $P=\begin{bmatrix}1 & 0 & 0 \\ 0 & 1& 0\end{bmatrix}$, then you have $$ R_1x = PQ^{-1}b$$ $$ x = (R_1)^{-1}PQ^{-1}b$$

Adam Latosiński
  • 10,723