QR Factorization via Householder Reflectors

Linear Algebra · Hard · Free problem

Let $A \in \mathbb{R}^{n \times p}$ be a full-rank matrix with $n \geq p$.

  1. Construct the first Householder reflector $H_1 = I - 2\mathbf{u}\mathbf{u}^T$ (where $\mathbf{u}$ is a unit vector) that zeroes out entries 2 through $n$ of the first column of $A$.
  1. Prove that $H$ is both symmetric and orthogonal.
  1. Derive the flop count for computing the full QR factorization via Householder reflectors, showing it is $O(np^2)$.
  1. Explain the numerical advantages of Householder QR over the classical Gram-Schmidt process.

Hints

  1. To zero entries below the diagonal in column $\mathbf{a}_1$, you need a reflection that maps $\mathbf{a}_1$ to $\alpha \mathbf{e}_1$ where $\alpha = \pm\|\mathbf{a}_1\|$. The reflection hyperplane is perpendicular to $\mathbf{a}_1 - \alpha \mathbf{e}_1$.
  2. For symmetry and orthogonality: expand $(I - 2\mathbf{u}\mathbf{u}^T)^2$ and use $\mathbf{u}^T\mathbf{u} = 1$.
  3. At step $k$, applying the reflector to the trailing submatrix costs $O((n-k)(p-k))$ flops. Sum over $k = 1, \ldots, p$ and extract the leading term in $n$ and $p$.

Worked Solution

How to Think About It: QR factorization writes $A = QR$ where $Q$ is orthogonal and $R$ is upper triangular. The Householder approach builds $R$ one column at a time: at each step, you apply a reflection that zeros out everything below the diagonal in the current column. Each reflection is an orthogonal transformation, so the product of reflections gives you $Q$. The beauty is that reflections are numerically stable -- they preserve vector norms exactly (in exact arithmetic) -- unlike Gram-Schmidt, which can lose orthogonality when columns are nearly dependent.

Formal Solution:

Part 1: Constructing the first Householder reflector.

Let $\mathbf{a}_1$ be the first column of $A$. We want $H_1 \mathbf{a}_1 = \alpha \mathbf{e}_1$ where $\mathbf{e}_1 = (1, 0, \ldots, 0)^T$ and $\alpha = \pm\|\mathbf{a}_1\|_2$.

Define the vector:

$\mathbf{v} = \mathbf{a}_1 - \alpha\, \mathbf{e}_1, \quad \alpha = -\text{sign}(a_{11})\|\mathbf{a}_1\|_2$

The sign choice $\alpha = -\text{sign}(a_{11})\|\mathbf{a}_1\|$ avoids catastrophic cancellation when $a_{11}$ and $\|\mathbf{a}_1\|$ have the same sign.

Set the unit vector:

$\mathbf{u} = \frac{\mathbf{v}}{\|\mathbf{v}\|_2}$

Then the Householder reflector is:

$H_1 = I - 2\mathbf{u}\mathbf{u}^T$

Verification: $H_1 \mathbf{a}_1 = \mathbf{a}_1 - 2\mathbf{u}(\mathbf{u}^T\mathbf{a}_1)$. Since $\mathbf{v} = \mathbf{a}_1 - \alpha\mathbf{e}_1$, we have $\mathbf{u}^T\mathbf{a}_1 = \frac{\mathbf{v}^T\mathbf{a}_1}{\|\mathbf{v}\|}$, and $\mathbf{v}^T\mathbf{a}_1 = \|\mathbf{a}_1\|^2 - \alpha a_{11}$. Working through the algebra:

$\|\mathbf{v}\|^2 = \|\mathbf{a}_1\|^2 - 2\alpha a_{11} + \alpha^2 = 2(\alpha^2 - \alpha a_{11})$

So $\mathbf{u}^T\mathbf{a}_1 = \frac{\alpha^2 - \alpha a_{11}}{\|\mathbf{v}\|} = \frac{\|\mathbf{v}\|^2}{2\|\mathbf{v}\|} = \frac{\|\mathbf{v}\|}{2}$.

Therefore $H_1\mathbf{a}_1 = \mathbf{a}_1 - 2\mathbf{u}\cdot\frac{\|\mathbf{v}\|}{2} = \mathbf{a}_1 - \mathbf{v} = \alpha\mathbf{e}_1$. The entries below position 1 are zeroed.

Part 2: Proving $H$ is symmetric and orthogonal.

Symmetry: $H^T = (I - 2\mathbf{u}\mathbf{u}^T)^T = I^T - 2(\mathbf{u}\mathbf{u}^T)^T = I - 2\mathbf{u}\mathbf{u}^T = H$.

Orthogonality: $H^T H = H H = (I - 2\mathbf{u}\mathbf{u}^T)(I - 2\mathbf{u}\mathbf{u}^T)$

$= I - 4\mathbf{u}\mathbf{u}^T + 4\mathbf{u}(\mathbf{u}^T\mathbf{u})\mathbf{u}^T = I - 4\mathbf{u}\mathbf{u}^T + 4\mathbf{u}\mathbf{u}^T = I$

using $\mathbf{u}^T\mathbf{u} = 1$ (since $\mathbf{u}$ is a unit vector).

Since $H^T = H$ and $H^T H = I$, we also have $H^2 = I$, so $H$ is an involution (its own inverse). Geometrically, it is a reflection through the hyperplane perpendicular to $\mathbf{u}$. $\blacksquare$

Part 3: Flop count for full QR.

The full QR factorization applies $p$ successive Householder reflectors. At step $k$ ($k = 1, \ldots, p$), the reflector $H_k$ acts on the trailing $(n - k + 1) \times (p - k + 1)$ submatrix.

At step $k$: - Computing $\mathbf{v}_k$: $O(n - k)$ flops. - Applying $H_k$ to the remaining $p - k$ columns: for each column, compute $\mathbf{u}_k^T \mathbf{c}_j$ (an inner product costing

(n-k)$ flops) and then the rank-1 update $\mathbf{c}_j - 2(\mathbf{u}_k^T\mathbf{c}_j)\mathbf{u}_k$ (another (n-k)$ flops). Total per column: $4(n-k)$ flops. Over $p - k$ columns: $4(n-k)(p-k)$ flops.

Summing over all $p$ steps:

$\text{Total flops} \approx \sum_{k=1}^{p} 4(n - k)(p - k) = 4\sum_{j=0}^{p-1}(n - 1 - j)(p - 1 - j)$

For $n \gg p$, the dominant term is:

$\approx 4\sum_{j=0}^{p-1} n(p - j) \approx 4n \cdot \frac{p^2}{2} = 2np^2$

So the total flop count is $\mathbf{2np^2 + O(np)}$ -- this is $O(np^2)$.

For the square case $n = p$, this gives n^3/3$ flops (accounting for the quadratic terms), which matches the standard textbook result.

Part 4: Numerical advantages over Gram-Schmidt.

  1. Backward stability: Householder QR is backward stable -- the computed $\hat{Q}$ and $\hat{R}$ satisfy $A + \delta A = \hat{Q}\hat{R}$ where $\|\delta A\| = O(\epsilon_{\text{mach}}) \|A\|$. Classical Gram-Schmidt (CGS) is NOT backward stable; it can produce a $Q$ that is far from orthogonal when columns of $A$ are nearly linearly dependent.
  1. Orthogonality preservation: Householder reflections are orthogonal transformations applied to the matrix. In floating point, $\hat{Q}$ from Householder satisfies $\|\hat{Q}^T\hat{Q} - I\| = O(\epsilon_{\text{mach}})$. With CGS, $\|\hat{Q}^T\hat{Q} - I\|$ can be $O(\kappa(A) \cdot \epsilon_{\text{mach}})$ where $\kappa(A)$ is the condition number -- disastrous for ill-conditioned matrices.
  1. Modified Gram-Schmidt (MGS) improves on CGS and is equivalent to Householder in some metrics, but Householder remains the gold standard for libraries like LAPACK because it also has superior cache performance (column-oriented memory access pattern).
  1. No explicit $Q$ needed: In many applications (least squares), you only need $R$ and the ability to apply $Q^T$ to a vector. Householder stores each reflector as a single vector $\mathbf{u}_k$, using $O(n)$ storage per step rather than forming the full $Q$ matrix.

Answer: The Householder reflector $H = I - 2\mathbf{u}\mathbf{u}^T$ with $\mathbf{u} = (\mathbf{a}_1 + \text{sign}(a_{11})\|\mathbf{a}_1\|\mathbf{e}_1) / \|\cdots\|$ zeros the subdiagonal entries of the first column. It is symmetric ($H^T = H$) and orthogonal ($H^2 = I$) because $\mathbf{u}$ is a unit vector. The full QR decomposition costs np^2$ flops. Householder is preferred over Gram-Schmidt because it is backward stable and preserves orthogonality to machine precision regardless of the condition number of $A$.

Intuition

Householder QR is the workhorse of numerical linear algebra. The idea is elegant: instead of orthogonalizing columns one at a time (Gram-Schmidt), you apply a sequence of reflections that each zero out an entire subdiagonal segment. Each reflection is a simple rank-1 perturbation of the identity, cheap to store and apply. The product of orthogonal matrices is orthogonal, so the accumulated transformations automatically give you an orthogonal $Q$ without any re-orthogonalization step.

The numerical stability advantage is the real reason this method dominates. In finance, when you run PCA on a returns matrix or solve a least-squares regression with correlated predictors, the matrix can be ill-conditioned. Gram-Schmidt will silently produce garbage orthogonal factors in this regime, while Householder keeps working. This is why every serious numerical library (LAPACK, NumPy, MATLAB) uses Householder QR as the default.

Open the full interactive solver →