QR factorisation

A mirror can swing an arrow onto an axis without changing its length. Watch a mirror turn until the arrow lands.

Gram–Schmidt, written down

Orthonormal bases and Gram–Schmidt turned the columns of

into three perpendicular unit arrows, $\mathbf{q}_1=(1,1,0)/\sqrt2$, $\mathbf{q}_2=(1,-1,2)/\sqrt6$ and $\mathbf{q}_3=(-1,1,1)/\sqrt3$. On the way it worked out some numbers and then threw them away: the shadows it took off, and the lengths it divided by. Keep them.

Run Gram–Schmidt backwards. Each column of $\mathbf{A}$ is its leftover plus its shadows, and its leftover is a length times a $\mathbf{q}$. So

where $r_{pq}=\mathbf{q}_p\cdot\mathbf{a}_q$ is a shadow for $p<q$, and $r_{pp}$ is the length of the $p$-th leftover. Column $q$ uses only $\mathbf{q}_1$ to $\mathbf{q}_q$. In matrices,

with the $\mathbf{q}$‘s as the columns of $\mathbf{Q}$ and the numbers $r_{pq}$ in an upper triangle $\mathbf{R}$. This is the QR factorisation. On this page $\mathbf{Q}$ and $\mathbf{R}$ mean these two factors; the bold $\mathbf{R}$ of the signal processing pages, a correlation matrix, is a different thing.

For our $\mathbf{A}$:

Watch the matrices rebuild $\mathbf{A}$.

Looking ahead to The determinant as area: the diagonal of $\mathbf{R}$ multiplies to $1.4142\cdot1.2247\cdot1.1547=2.000$, which is the size of $\det\mathbf{A}=-2$. $\mathbf{Q}$ keeps volumes, so the triangle $\mathbf{R}$ carries all of the volume.

Least squares with QR

In Least squares we solved $\mathbf{A}^\top\mathbf{A}\hat{\mathbf{x}}=\mathbf{A}^\top\mathbf{y}$. Put $\mathbf{A}=\mathbf{Q}\mathbf{R}$ into it. Since $\mathbf{Q}^\top\mathbf{Q}=\mathbf{I}$, the left side is $\mathbf{R}^\top\mathbf{R}\hat{\mathbf{x}}$ and the right side is $\mathbf{R}^\top\mathbf{Q}^\top\mathbf{y}$. Cancel $\mathbf{R}^\top$:

$\mathbf{R}$ is a triangle, so this is one back substitution.

Try it on the line fit of that page, with columns $(1,1,1)$ and $(0,1,2)$ and heights $\mathbf{y}=(1,3,2)$. Gram–Schmidt gives $\mathbf{q}_1=(1,1,1)/\sqrt3$ with $r_{11}=\sqrt3$. The shadow of $(0,1,2)$ on $\mathbf{q}_1$ is $r_{12}=3/\sqrt3=\sqrt3$, which leaves $(-1,0,1)$, of length $r_{22}=\sqrt2$, so $\mathbf{q}_2=(-1,0,1)/\sqrt2$. Then $\mathbf{Q}^\top\mathbf{y}=(6/\sqrt3,\ 1/\sqrt2)=(3.4641,\ 0.7071)$, and back substitution gives

The same line, $y=1.5+0.5t$.

Why bother, when the normal equations gave the same answer? Because $\mathbf{A}^\top\mathbf{A}$ squares the difficulty. When two columns are nearly parallel, small rounding errors in $\mathbf{A}^\top\mathbf{A}$ become large errors in $\hat{\mathbf{x}}$. Here is an extreme case. Take a tiny number $\varepsilon=10^{-8}$ and

whose exact answer is $(1,1)$. The matrix $\mathbf{A}^\top\mathbf{A}$ has $1+\varepsilon^2$ on its diagonal, and a computer working to about 16 digits rounds $1+10^{-16}$ to exactly 1. It stores $\begin{bmatrix}1&1\\1&1\end{bmatrix}$, which has no inverse, and the normal equations fail. QR never forms $\mathbf{A}^\top\mathbf{A}$, and it returns $(1,1)$.

Mirrors that clear a column

Programs rarely run Gram–Schmidt itself; rounding errors slowly spoil the right angles. They use mirrors. The picture at the top of the page is one step.

A mirror is a reflection across a line, an orthogonal matrix from the Gram–Schmidt page, so it keeps lengths. Choose the mirror that sends the first column onto the first axis. The column $(3,4)$ has length 5, so the mirror must send it to $(5,0)$. The mirror line lies halfway between the two, along $(2,1)$. Its matrix, called a Householder reflection, is

where $\mathbf{v}$ is at a right angle to the mirror line. Here $\mathbf{v}^\top\mathbf{v}=20$, and $\mathbf{H}=\begin{bmatrix}0.6&0.8\\0.8&-0.6\end{bmatrix}$. Check: $\mathbf{H}\,(3,4)=(1.8+3.2,\ 2.4-2.4)=(5,0)$. This $\mathbf{H}$ is local to this page; signal processing pages use bold $\mathbf{H}$ for a filter’s matrix.

Apply the mirror to a whole matrix whose first column is $(3,4)$:

The first column is cleared below the top. The result is upper triangular, so it is $\mathbf{R}$. Reflecting twice brings you back, $\mathbf{H}^2=\mathbf{I}$, so $\mathbf{A}=\mathbf{H}\mathbf{R}$ with $\mathbf{Q}=\mathbf{H}$. Larger matrices use one mirror per column. The diagonal again carries the area: $5\cdot3.4=17$, and $\det\begin{bmatrix}3&1\\4&7\end{bmatrix}=21-4=17$.

Here $r_{22}=-3.4$ is negative. That is only a convention: a mirror that sends $(3,4)$ to $(-5,0)$ works equally well. Programs choose the sign that keeps $\mathbf{v}$ large, $\mathbf{v}=\mathbf{a}+\lVert\mathbf{a}\rVert\mathbf{e}_1$, here $(8,4)$. The other choice subtracts two nearly equal numbers when the column is close to the axis. For $\mathbf{a}=(1,0.01)$, $\lVert\mathbf{a}\rVert=1.00005$, and $1-1.00005=-0.00005$ keeps only the last few digits of the two numbers.

Worked example

1. From Gram–Schmidt. The $3\times3$ matrix above has $\mathbf{R}$ with diagonal $1.4142$, $1.2247$, $1.1547$, whose product is 2.
2. A line fit. $\mathbf{R}=\begin{bmatrix}1.7321&1.7321\\0&1.4142\end{bmatrix}$, $\mathbf{Q}^\top\mathbf{y}=(3.4641,0.7071)$, so $c_1=0.5$ and $c_0=1.5$.
3. Already orthonormal. If the columns of $\mathbf{A}$ are already orthonormal, $\mathbf{Q}=\mathbf{A}$ and $\mathbf{R}=\mathbf{I}$.
4. A mirror. $\mathbf{H}=\begin{bmatrix}0.6&0.8\\0.8&-0.6\end{bmatrix}$ sends $(3,4)$ to $(5,0)$, and $\mathbf{H}\mathbf{H}=\mathbf{I}$.
5. Determinant from R. The size of $\det\mathbf{A}$ is the product of the sizes of the diagonal of $\mathbf{R}$: 17 for the mirror example.
6. Predict. As two columns become nearly parallel, the second leftover of Gram–Schmidt shrinks, so $r_{22}$ falls towards 0. That small number is the warning sign that Conditioning measures.

Where you’ll meet this

Library routines for least squares, such as NumPy’s lstsq, factor $\mathbf{A}$ with orthogonal matrices instead of forming $\mathbf{A}^\top\mathbf{A}$. For a tall matrix, Householder QR costs about twice as much as the normal equations, and buys back the accuracy that squaring would lose. Statistics software fits regression lines this way. Some Kalman filters carry a triangular factor like $\mathbf{R}$ in place of a covariance matrix, for the same reason.

Reference card