1 Introduction to QR decomposition
1.1 Definition and basic notation
QR decomposition is a matrix factorization that represents a given matrix \(A\) as the product of two factors, \[ A = QR, \] where \(Q\) is orthogonal (real case) or unitary (complex case) and \(R\) is upper triangular (or has an upper-triangular form consistent with the matrix dimensions). For \(A\in \mathbb{R}^{m\times n}\), common conventions include the “full” form with \(Q\in \mathbb{R}^{m\times m}\) and \(R\in \mathbb{R}^{m\times n}\), and the “thin” (economy-size) form with \(Q\in \mathbb{R}^{m\times n}\) and \(R\in \mathbb{R}^{n\times n}\).
1.2 Geometric interpretation via orthogonal bases
Geometrically, the columns of \(Q\) provide an orthonormal basis for the relevant subspace associated with \(A\). In the full column-rank (tall) case, the column space of \(A\) is spanned by a subset of columns of \(Q\). The triangular factor \(R\) then stores how the original columns of \(A\) are expressed in that orthonormal basis. This viewpoint is particularly useful because orthonormal bases preserve lengths and angles, simplifying analysis and computation.
1.3 Relationship to orthogonal projections
Orthogonal projections onto column spaces are central in least squares and approximation. QR decomposition links projections to triangular systems: when \(A=QR\), solving projection-related problems often reduces to operations involving \(R\) (and occasionally \(Q^\top\)). The orthogonality of \(Q\) ensures that projecting onto the column space spanned by \(Q\) can be computed reliably through transpose (or conjugate transpose in complex settings) without explicitly forming projection matrices.
2 Mathematical foundations
2.1 Orthogonality and unitary matrices
A real matrix \(Q\) is orthogonal if \(Q^\top Q = I\). In the complex case, unitarity replaces transpose with conjugate transpose, \(Q^*Q=I\). Orthogonality implies that \(Q\) preserves inner products, and therefore distances and norms. These properties underpin numerical stability: operations can be performed using orthogonal/unitary transformations without amplifying errors dramatically.
2.2 Rank, column space, and null space
The rank of \(A\) is reflected in \(R\). In many practical QR variants, the diagonal entries of \(R\) indicate whether the corresponding directions in the column space are present or nearly absent. When \(A\) is rank-deficient, some diagonal elements become zero (exact arithmetic) or small (floating point), and the remaining nonzero structure identifies the effective column space. The null space structure is correspondingly related to the orthogonal complement of the column space.
2.3 Existence and uniqueness considerations
A QR factorization exists for any real or complex matrix under mild conditions. Uniqueness is limited: even when \(R\) has strictly positive diagonal entries, uniqueness can be enforced only up to sign/phase conventions depending on the chosen convention for \(Q\) and \(R\). In general, many decompositions may satisfy \(A=QR\), especially when \(A\) is rank-deficient.
2.4 Variants: full, thin, and economy-size QR
“Full QR” uses a square orthogonal/unitary \(Q\) of size \(m\times m\) to match the full range of transformations. “Thin QR” or “economy-size QR” reduces computation by forming only the columns of \(Q\) that are needed for subsequent tasks, typically yielding a smaller \(Q\) and a square \(R\) for \(m\ge n\). The choice affects memory footprint and runtime, while the numerical quality depends more strongly on the chosen algorithm (e.g., Householder vs. Gram–Schmidt) and on conditioning.
3 Computing QR decomposition
3.1 Gram–Schmidt orthogonalization
Gram–Schmidt constructs orthonormal vectors sequentially by subtracting projections onto previously computed directions.
3.1.1 Classical Gram–Schmidt
Classical Gram–Schmidt performs projection subtraction using the original coefficients computed from the current residual relative to already formed vectors. In exact arithmetic it produces orthonormal columns, but in floating point arithmetic it can suffer from loss of orthogonality, especially when columns of \(A\) are nearly linearly dependent.
3.1.2 Modified Gram–Schmidt
Modified Gram–Schmidt reorganizes the projection/subtraction steps to reduce the accumulation of rounding errors. It typically offers improved orthogonality compared with the classical form, though stability can still degrade in very ill-conditioned settings or for large problems without additional safeguards.
3.2 Householder reflections
Householder-based QR uses orthogonal reflections to systematically introduce zeros below the diagonal of \(R\). Each reflection is designed to transform a vector so that a chosen component aligns with an axis, allowing multiple eliminations while preserving norms.
3.2.1 Constructing Householder vectors
A Householder reflection has the form \(H = I - 2vv^\top/(v^\top v)\) (real case), where \(v\) is chosen so that \(Hx\) aligns with a coordinate direction for the vector \(x\) being transformed. In QR, successive Householder matrices are applied to submatrices, driving the lower triangular part toward zero and producing an upper-triangular \(R\). The overall orthogonal factor \(Q\) can be represented implicitly as a product of reflections.
3.2.2 Numerical stability considerations
Householder transformations are norm-preserving and tend to be backward stable in common floating point models. As a result, they usually maintain orthogonality of the computed \(Q\) more robustly than Gram–Schmidt for difficult matrices, making them the default choice in many numerical libraries.
3.3 Givens rotations
Givens rotations are plane rotations that zero out individual elements by acting on two coordinates at a time. They are especially convenient when working with sparse matrices or when updating a factorization after small changes.
3.3.1 Targeted zeroing strategies
Each Givens rotation is selected to annihilate a particular entry in \(A\) while affecting only a small subset of coordinates. By applying a sequence of such rotations, one systematically constructs \(R\) with zeros below the diagonal. This “elementwise” control can be advantageous for structured sparsity.
3.3.2 Sparse/structured applications
Because Givens rotations can be applied locally, they can preserve sparsity patterns more effectively than dense reflections in certain formats. That makes them useful in large-scale problems where fill-in would be costly, such as some iterative solvers and incremental update contexts.
3.4 Choosing an algorithm in practice
The choice depends on the problem structure and required robustness.
3.4.1 Accuracy vs. computational cost
Householder QR typically offers strong accuracy with moderate overhead. Modified Gram–Schmidt can be faster in some settings (especially with hardware and memory considerations) but may require reorthogonalization for high accuracy. Givens rotations can be competitive when sparsity or structured updates dominate cost.
4 QR for least squares and solving linear systems
4.1 Least-squares formulation
| Least squares seeks \(x\) minimizing \(\|Ax-b\|_2\). With QR, the problem is transformed into a structured system: if \(A=QR\) (thin QR for \(m\ge n\)), then |
|---|
\[
| \|Ax-b\|_2 = \|QRx-b\|_2 = \|Rx-Q^\top b\|_2 |
|---|
\] because orthogonality of \(Q\) preserves the Euclidean norm. The minimizer can be extracted by solving an upper-triangular system involving \(R\).
4.2 Solving overdetermined systems
For overdetermined systems with full column rank (more equations than unknowns), the least-squares solution corresponds to the exact solution of the triangular system \[ Rx = Q^\top b \] restricted to the leading rows that define the triangular structure. The remaining components of \(Q^\top b\) relate to the residual orthogonal to the column space of \(A\).
4.3 Rank-deficient and truncated QR
When columns are linearly dependent or nearly so, \(R\) may contain small or zero diagonal entries, complicating direct back substitution. Common strategies include:
- computing a truncated solution based on an estimated numerical rank;
- using pivoting variants to identify stable directions;
- applying regularization when appropriate.
These approaches aim to produce a meaningful minimizer without amplifying noise in nearly null directions.
4.4 Handling constraints and projections
4.4.1 Computing the residual norm from \(R\)
| Residuals can be evaluated efficiently using the orthogonal transform. If \(A=QR\) with thin \(Q\), then the residual norm \(\|Ax-b\|_2\) equals the norm of the components of \(Q^\top b\) corresponding to the orthogonal complement of the column space. In many implementations, this reduces to combining the trailing part of transformed data with the triangular system solution, avoiding explicit computation of the projection matrix. |
|---|
5 Numerical properties and conditioning
5.1 Backward/forward error perspectives
Numerical analysis of QR often focuses on backward stability: the computed factors correspond to an exact QR decomposition of a slightly perturbed matrix \(\tilde{A}\). This perspective is meaningful because it relates observed errors to small data perturbations rather than to intrinsic algorithmic failure. Forward error (deviation in the solution itself) then depends on the conditioning of the underlying problem.
5.2 Stability comparisons with normal equations
Normal equations solve \(A^\top A x = A^\top b\), squaring the condition number in the process. QR avoids this amplification by using orthogonal transformations directly on \(A\) and \(b\). As a result, QR-based least squares typically yields more reliable solutions in floating point arithmetic, particularly when \(A\) is ill-conditioned.
5.3 Conditioning effects on \(Q\) and \(R\)
Even with stable algorithms, the conditioning of \(A\) influences the interpretation of \(R\). Small singular values manifest as small diagonal entries or near-linear dependencies among columns, affecting the sensitivity of the computed solution to perturbations. While orthogonal transformations themselves are well-behaved, the final accuracy is limited by how well-posed the problem is.
5.4 Reorthogonalization and loss of orthogonality
In finite precision, computed vectors in \(Q\) may drift from exact orthogonality. For Gram–Schmidt variants, reorthogonalization (repeating the orthogonalization step) can restore orthogonality and improve accuracy. For Householder and Givens methods, orthogonality loss is typically smaller, though it can still occur due to rounding and truncation choices in practical implementations.
6 Applications across applied mathematics
6.1 Regression and data fitting
In statistics and machine learning contexts, QR decomposition provides a numerically stable way to compute least-squares regressions, including linear models and polynomial fits. Compared with solving via normal equations, QR tends to deliver improved numerical behavior when predictors are correlated or when the design matrix has a wide dynamic range.
6.2 Signal processing and orthogonal transforms
Orthogonality is valuable in signal representations because it preserves energy (norms) and simplifies interpretation of components. QR-related orthogonalization steps appear in algorithms that construct orthonormal bases, perform subspace tracking, or stabilize computations involving projections and orthogonal matching.
6.3 Numerical methods for differential equations
Many discretizations of differential equations lead to linear systems solved repeatedly or to time-stepping schemes requiring robust factorization. QR methods can support stable projections onto reduced spaces or facilitate solving least-squares problems arising in methods such as collocation and certain iterative updates.
6.4 Model order reduction and orthogonalization steps
In reduced-order modeling, it is common to build low-dimensional subspaces from simulation data. QR decomposition is frequently used to orthonormalize these snapshots, yielding a basis with improved numerical properties. The triangular factor can then help map between coordinate representations in the reduced space and the original system.
7 Rank, pivoting, and enhanced QR variants
7.1 Column pivoting
Column pivoting reorders columns of \(A\) during the factorization to improve numerical robustness. The goal is to choose, at each step, a column that is likely to contribute more to the rank structure, preventing growth of roundoff effects in later triangular solves.
7.2 QR with column permutation
Pivoting introduces a permutation matrix \(P\) such that \[ A P = QR. \] Equivalently, \(A = QRP^\top\). In least squares, the permutation modifies how unknown coefficients are ordered and can yield a more reliable determination of which components are supported by the data.
7.3 Estimating numerical rank using \(R\)
A practical use of pivoted QR is determining numerical rank: diagonal entries of \(R\) are compared against a tolerance derived from their magnitude and the precision level. Entries below the threshold are treated as effectively zero, producing a truncated factorization suited for rank-deficient least squares.
7.4 Pivoted QR for robust least squares
Pivoted QR is often preferred when \(A\) may be ill-conditioned or nearly rank-deficient. It helps avoid unstable back substitution by ensuring that the largest available contributions appear earlier in the triangular system. This leads to improved behavior of truncated solutions and more faithful residual minimization under uncertainty.
8 QR algorithm connections (eigenvalues and Schur forms)
8.1 From QR decomposition to QR iteration
For eigenvalue computation, QR decomposition serves as the inner engine of QR iteration applied to a square matrix \(A\). A typical step factorizes \(A_k = Q_k R_k\) and then forms \[ A_{k+1} = R_k Q_k. \] Under suitable conditions, this process converges such that \(A_k\) approaches an upper triangular or quasi-triangular form related to eigenvalues.
8.2 Shifts and accelerated convergence
Shifts modify the QR iteration by subtracting a scalar multiple of the identity before factorization and adding it back after. This improves convergence speed, especially when eigenvalues are clustered or when a targeted region of the spectrum is of interest. The shift can be selected via heuristics based on entries near the bottom of the current iterate.
8.3 Relation to Hessenberg reduction
To make QR iteration efficient, the initial matrix is often reduced to upper Hessenberg form, which preserves eigenvalues while introducing a sparse-like structure. Because QR steps become cheaper on Hessenberg matrices, the combination of Hessenberg reduction and shifted QR iteration is a standard route to eigenvalue and Schur form computations.
8.4 Practical implementation notes
In practice, numerical software balances stability, convergence control, and performance. Choices include the use of implicit shifts, handling of special cases (e.g., small subdiagonal entries), and termination criteria that relate to matrix norms or residuals of the Schur decomposition.
9 Implementation details and best practices
9.1 Data types, tolerances, and stopping criteria
Correct implementation depends on numeric precision (single vs. double) and on consistent tolerances. In rank estimation and iterative contexts, tolerances are often scaled by norms of \(A\) or by machine precision. Stopping criteria for iterative schemes (including QR iteration for eigenvalues) typically use thresholds on subdiagonal elements or changes in the iterate.
9.2 Complexity and memory considerations
The computational cost of QR depends on matrix dimensions. For dense matrices, Householder QR has a cost on the order of \(O(mn^2)\) when \(m\ge n\). Memory usage varies: storing \(Q\) explicitly can be expensive, so many implementations store Householder vectors compactly and apply \(Q\) or \(Q^\top\) implicitly when needed.
9.3 Verification and reconstruction checks
| A standard validation is to reconstruct \(A\) (or a transformed form) from computed factors and measure \(\|A-QR\|\) relative to \(\|A\|\). For least squares, verifying residual norms and comparing predicted residuals from triangular factors can detect dimension mismatches, incorrect transposes, or errors in permutation handling. |
|---|
9.4 Interfacing with libraries and APIs
In software libraries, QR interfaces typically return either explicit \(Q\) and \(R\) or a compact representation plus metadata describing the form. Implementations also differ in whether they compute full, thin, or economy-size factors, and in whether they apply column pivoting. Careful attention to returned shapes and documented conventions is essential when integrating QR routines into larger workflows.
10 Common pitfalls and troubleshooting
10.1 Loss of orthogonality in naive methods
Using the unmodified classical Gram–Schmidt procedure can lead to noticeable loss of orthogonality in \(Q\), producing inaccurate residuals or unstable least-squares solutions. Symptoms include a reconstructed \(A\) that deviates unexpectedly, or solutions that differ significantly under minor perturbations.
10.2 Interpreting signs/phase ambiguities in \(Q\)
Even with correct computation, the factors are not unique: columns of \(Q\) may differ by sign (or complex phase) while compensating changes occur in \(R\). When comparing results across implementations, this ambiguity can appear as “differences” in \(Q\) even though the product \(QR\) matches to numerical precision.
10.3 Interpreting “upper triangular” in rectangular cases
For rectangular matrices, “upper triangular” refers to a specific structural shape: the nonzero entries of \(R\) lie on and above a diagonal-like boundary appropriate to the factorization dimensions. Misinterpreting which sub-block is triangular (especially when using thin vs. full QR) can cause errors in back substitution or in downstream residual calculations.
10.4 Debugging incorrect shapes and dimensions
Many implementation bugs arise from mismatched matrix dimensions: confusing \(m\times n\) with \(n\times m\), using the wrong transpose convention, or treating a compact Householder representation as an explicit \(Q\). Debugging often involves checking the reported dimensions of returned factors, verifying that \(Q^\top Q\) (or \(Q^*Q\)) is close to identity, and confirming that \(QR\) reproduces \(A\) within tolerance.