1 Definition and basic properties
1.1 The QR decomposition \(A=QR\)
For a given matrix \(A \in \mathbb{R}^{m\times n}\) (or \(A \in \mathbb{C}^{m\times n}\)), a QR factorization expresses \(A\) as \[ A = QR, \] where \(Q\) has orthonormal columns (or is unitary in the complex case) and \(R\) is an upper triangular matrix. When \(m \ge n\), \(Q\) is commonly taken as \(m\times n\) with orthonormal columns, while \(R\) is \(n\times n\).
The factorization separates geometric information (captured by \(Q\)) from scaling and coupling among columns (captured by \(R\)). This separation is the reason QR is widely used in numerical linear algebra: it enables stable transforms of problems into forms that are easy to solve.
1.2 Orthogonality/unitarity and its implications
In the real case, “orthogonal” means \(Q^TQ = I\) when \(Q\) has orthonormal columns. In the complex case, the analogous condition is \(Q^*Q = I\), where \(Q^*\) denotes the conjugate transpose. These identities imply that multiplying by \(Q\) or \(Q^*\) preserves Euclidean norms and inner products: \[
| \|Qv\|_2 = \|v\|_2,\qquad \langle Qv, Qw\rangle = \langle v, w\rangle. |
|---|
\] Consequently, QR transformations tend to control numerical errors better than factorizations that do not rely on norm-preserving steps.
1.3 Triangular structure of \(R\)
The matrix \(R\) is upper triangular, meaning \(r_{ij}=0\) for \(i>j\). This structure matters computationally because solving systems involving \(R\) reduces to back substitution. For least-squares problems, the triangular form makes it possible to solve for unknown coefficients while the remaining rows encode residual information.
In exact arithmetic, the triangularity is guaranteed by the way \(Q\) is constructed. In floating-point arithmetic, triangular entries that are theoretically zero may appear as small numerical artifacts; algorithms aim to keep these artifacts small.
1.4 Rank, pivoting, and reduced vs. full QR
The rank of \(A\) influences how QR behaves and how solutions are interpreted. If \(A\) is rank-deficient, some diagonal entries of \(R\) may be zero (or numerically tiny), which affects solvability and can require special handling (e.g., truncated or pivoted QR).
Two common variants differ in how many columns of \(Q\) are included:
- Reduced QR: \(Q\) is \(m\times n\) (when \(m\ge n\)), and \(R\) is \(n\times n\).
- Full QR: \(Q\) is \(m\times m\) and \(R\) is \(m\times n\), with additional structure in the extra rows of \(R\).
Column pivoting (often called QR with pivoting) reorders columns during factorization to improve numerical rank detection. The result has the form \[ A P = QR, \] where \(P\) is a permutation matrix. This variant is used when columns are not well conditioned relative to each other.
1.5 Geometric interpretation in vector spaces
Geometrically, QR factorization provides an orthonormal basis for the column space of \(A\). The columns of \(Q\) span the same subspace as the columns of \(A\) (up to rank issues), while \(R\) records how the original columns of \(A\) project onto that basis. In least-squares settings, this means the approximation \(Ax\) is expressed via projections onto an orthonormal basis, and the residual corresponds to the component orthogonal to the column space.
This viewpoint clarifies why QR is effective: projecting onto an orthonormal basis isolates “best approximation” structure from the algebraic complexity of the original matrix.
2 Computing QR factorization
2.1 Classical Gram–Schmidt method
The classical Gram–Schmidt procedure constructs orthonormal vectors by orthogonalizing each new column of \(A\) against the previously computed orthonormal set. If the columns of \(A\) are \(a_1,\dots,a_n\), Gram–Schmidt produces vectors \(q_1,\dots,q_n\) and coefficients that assemble into \(R\).
In exact arithmetic, classical Gram–Schmidt yields a valid QR factorization. In floating-point arithmetic, it can suffer from loss of orthogonality when columns are nearly linearly dependent, because computed orthogonal components may already contain rounding error.
2.1.1 Numerical stability considerations
The loss of orthogonality is especially pronounced when the matrix columns have strong correlations. As orthogonality degrades, the resulting \(Q\) may deviate from truly orthonormal behavior, which can lead to larger residuals in least-squares solutions compared with more stable approaches.
This motivates improved variants such as modified Gram–Schmidt or reflector-based methods.
2.1.2 Modified Gram–Schmidt
Modified Gram–Schmidt rearranges the order of operations so that orthogonalization is performed incrementally in a way that tends to preserve orthogonality better in finite precision. While it still relies on inner products and subtraction (and therefore can accumulate rounding error), it often provides significantly improved behavior over the classical algorithm.
In practice, modified Gram–Schmidt is commonly used when the implementation is simpler than Householder transformations, particularly in some medium-scale contexts.
2.2 Householder reflections
Householder methods build \(Q\) using a sequence of orthogonal reflectors. A Householder reflector has the form \[ H = I - 2\frac{vv^T}{v^Tv} \] in the real case (and the conjugate-transpose analog in the complex case). Each reflector is chosen to map a selected vector to a multiple of a coordinate vector, systematically introducing zeros below the diagonal of \(R\).
2.2.1 Construction of Householder vectors
To zero out entries in a column, an algorithm constructs a vector \(v\) such that applying \(H\) to a target subvector aligns it with a desired direction. Careful scaling is used to avoid cancellation and overflow. The choice of sign for the component used in building \(v\) is typically selected to enhance numerical robustness.
This construction converts a general matrix into upper triangular form while maintaining orthogonality by design.
2.2.2 Accumulating reflectors to form \(Q\)
Rather than forming \(Q\) explicitly, implementations often store the reflector vectors and their scalar factors. When needed, \(Q\) or \(Q^*\) is applied to a vector by sequentially applying the reflectors. This avoids storing a full \(m\times n\) orthonormal matrix and improves both memory usage and computational efficiency.
When \(Q\) is eventually required (for example, to compute certain diagnostics), it can be explicitly assembled from the stored reflectors.
2.3 Givens rotations
Givens rotations eliminate selected elements using plane rotations. A Givens rotation \(G(i,j,\theta)\) acts nontrivially only on coordinates \(i\) and \(j\) and rotates a 2D subvector to introduce a zero in a chosen position. For QR, a sequence of rotations targets entries below the diagonal, one at a time.
2.3.1 Zeroing specific entries systematically
The algorithm selects a target element \(a_{k\ell}\) below the diagonal and computes a rotation angle \(\theta\) so that the rotation makes this entry zero while preserving other entries in a controlled way. Repeating this process over the matrix yields an upper triangular \(R\).
This “localized” elimination is well-suited to situations where only certain entries need alteration.
2.3.2 Sparse and incremental update advantages
Because each rotation affects only two rows, Givens rotations can be advantageous for sparse matrices: they introduce structured fill-in and can be integrated into incremental updates. In some contexts—such as online or streaming updates—rotations can be applied as new rows arrive without recomputing the entire factorization.
In contrast to global reflector operations, the granularity of Givens updates can preserve sparsity more effectively depending on the matrix pattern.
2.4 Algorithmic complexity and operation counts
For dense matrices with \(m\ge n\), the dominant cost of QR factorization scales on the order of \(O(mn^2)\). The constant factors differ among algorithms: Householder reflectors typically have favorable performance and stability properties for dense arithmetic; Gram–Schmidt variants may have higher sensitivity; Givens rotations can increase operation counts but may reduce memory movement for sparse or structured problems.
Operation counts also depend on whether a reduced or full QR is computed, and on whether \(Q\) is formed explicitly or applied implicitly.
3 QR for least-squares problems
3.1 From \(Ax \approx b\) to a triangular solve
Given an overdetermined system \(Ax \approx b\) with \(A\in\mathbb{R}^{m\times n}\) and \(m\ge n\), the QR factorization enables a transformation: \[ A = QR \quad \Rightarrow \quad QRx \approx b. \] Multiplying by \(Q^T\) yields \[ Rx \approx Q^T b. \] With reduced QR and \(R\) upper triangular, the approximation becomes a system that can be solved by back substitution on the leading triangular part, while the remaining components of \(Q^Tb\) correspond to the residual norm.
This reduces a general least-squares task to a simple triangular solve plus residual evaluation.
3.2 Deriving the normal equations-free approach
Standard least-squares approaches often derive the normal equations \(A^TAx=A^Tb\). However, forming \(A^TA\) squares the condition number and can amplify rounding error. QR avoids this by working directly with orthogonal transformations:
- \(Q\) preserves vector norms and inner products.
- The least-squares solution is obtained from triangular structure in \(R\) rather than solving a potentially ill-conditioned normal equation system.
Thus, QR provides a numerically safer route to least-squares coefficients.
3.3 Overdetermined vs. underdetermined systems
| - Overdetermined (\(m>n\)): the system has more equations than unknowns, so one typically seeks the vector \(x\) minimizing \(\|Ax-b\|_2\). QR naturally supports this by mapping the problem to a triangular system involving \(R\). |
|---|
| - Underdetermined (\(m<n\)): the system has fewer equations than unknowns. One common objective is the minimum-norm solution among all \(x\) satisfying \(Ax=b\) (or among those minimizing \(\|Ax-b\|\)). QR can be adapted using QR factorizations of \(A^T\) or by using alternative decompositions, with the key idea still relying on orthogonal transformations. |
The precise form depends on which dimension is larger, but QR remains useful due to its geometric stability properties.
3.4 Handling rank deficiency
When columns of \(A\) are linearly dependent or nearly so, the least-squares minimizer may not be unique. Rank deficiency manifests as small or zero diagonal entries in \(R\). Practical strategies include:
- Truncated QR: use only the leading part of \(R\) corresponding to the detected rank.
- QR with column pivoting: reorder columns to expose the numerical rank more reliably, then compute a least-squares solution consistent with that rank.
These approaches help avoid division by extremely small values and provide more meaningful solutions.
3.5 Error analysis and conditioning insights
In floating-point arithmetic, orthogonal transformations are particularly important because they tend to introduce errors in a controlled manner. The conditioning of the least-squares problem is governed by the singular values of \(A\), but QR helps ensure that the computed solution is close to what would be obtained in exact arithmetic, up to modest perturbations.
A typical takeaway is that QR-based least-squares solvers generally achieve good backward stability: the computed solution corresponds to solving a slightly perturbed problem, with the perturbation size tied to machine precision and the magnitude of \(A\).
4 Variants and enhancements
4.1 Column pivoting in QR (QR with pivoting)
Pivoting modifies QR by permuting columns of \(A\) to improve stability and reveal effective rank. The factorization is written as \[ AP = QR, \] where \(P\) is a permutation matrix. During factorization, the algorithm chooses the next column based on estimates of which remaining columns are likely to contribute significantly to the span.
4.1.1 Estimating numerical rank
Numerical rank is assessed using the magnitudes of diagonal entries of \(R\) (or related criteria). If these entries decay below a threshold relative to the largest diagonal element, the remaining part of \(R\) is treated as representing near-null directions. This helps distinguish true dependencies from benign scaling effects.
Pivoting is particularly relevant for problems where columns have heterogeneous norms or where near-collinearity is present.
4.2 Economic/reduced QR factorization
Reduced QR computes only the portion of \(Q\) needed to represent the column space of \(A\). This lowers memory usage and can reduce runtime, especially when \(m\) is large. The trade-off is that \(Q\) is not a full orthogonal matrix, but for least-squares tasks and subspace computations the reduced form is typically sufficient.
4.3 QR factorization for special matrix types
QR methods adapt to matrix structure and storage patterns.
4.3.1 Dense matrices
For dense data, Householder-based QR is common because it achieves strong accuracy and efficiency under typical dense linear algebra cost models. The implementation often leverages blocked algorithms to improve cache performance.
4.3.2 Sparse matrices
Sparse matrices require careful attention to fill-in and memory traffic. Givens rotations can offer fine-grained updates, while sparse-aware strategies may still use reflector-like approaches depending on sparsity structure. In all cases, the goal is to preserve sparsity as much as possible while maintaining numerical quality.
4.3.3 Structured matrices
Structured matrices (e.g., banded or having exploitable patterns) can allow QR variants that reduce computational burden. By restricting operations to relevant bands or leveraging symmetry-like properties when applicable, algorithms can lower both time and storage relative to generic QR.
4.4 Iterative refinement and backward error perspectives
Even after QR produces an initial least-squares solution, iterative refinement can improve accuracy. The idea is to:
- Compute a solution \(x\).
- Form a residual \(r=b-Ax\).
- Solve a correction equation (often using the same QR factors).
- Update \(x\) and repeat until improvement saturates.
From a backward error perspective, QR’s use of orthogonal transformations typically yields favorable error behavior, making it a good foundation for refinement steps.
5 Implementation details in numerical computing
5.1 Storing \(Q\) implicitly (reflector/rotation form)
Explicitly storing a full \(Q\) can be expensive. Many implementations store only:
- Householder reflector vectors and scaling factors, or
- the parameters defining each Givens rotation.
When computing \(Q^Tb\) or applying \(Q\) to another matrix, the library applies these stored transformations in sequence. This preserves the numerical properties of the factorization while using less memory.
5.2 Forming \(R\) and extracting solution vectors
The upper triangular matrix \(R\) is produced during the elimination process. For Householder QR, the reflectors are applied to transform \(A\) so that the leading part becomes triangular; the resulting triangular entries can be extracted directly.
In least-squares, after computing \(y = Q^Tb\), the solution is obtained by solving \(Rx = y_{1:n}\) via back substitution. If rank deficiency is considered, the solver may instead apply a truncation or pivot-based selection before back substitution.
5.3 Practical considerations for scaling and normalization
Numerical behavior depends on how vectors are scaled during reflector or rotation construction. Many algorithms include safeguards:
- avoid overflow/underflow when computing norms,
- handle very small magnitudes safely,
- choose signs in reflector formulas to reduce cancellation.
Proper scaling helps keep the computed orthogonal transformations accurate and ensures that diagonal entries of \(R\) reflect the true problem rather than arithmetic artifacts.
5.4 Verification, residual norms, and diagnostics
After solving, the quality of the result is often assessed using residual norms: \[
| \|b-Ax\|_2. |
|---|
\] Additionally, in stable QR implementations, it is common to estimate whether orthogonality constraints are satisfied (e.g., checking that \(Q^TQ\) is close to identity in explicit constructions, or using internal diagnostics in implicit settings).
Diagnostics may also include checking the size pattern of \(R\)’s diagonal to identify potential rank issues.
5.5 Typical pitfalls and how to avoid them
Common pitfalls include:
- Forming normal equations instead of using QR, which can degrade accuracy.
- Using classical Gram–Schmidt in sensitive problems without mitigation, due to orthogonality loss.
- Ignoring rank deficiency, leading to unstable divisions during back substitution.
- Improper thresholding for truncation in pivoted QR, which can produce inconsistent solutions across scales.
Robust implementations address these through modified Gram–Schmidt, Householder reflectors, pivoting strategies, and careful rank thresholds.
6 Applications in applied mathematics
6.1 Solving linear systems and least squares in practice
Beyond pure least squares, QR is used wherever triangular solves after an orthogonal transformation are beneficial. It provides reliable solutions for overdetermined systems and supports robust residual evaluation, making it a standard tool in scientific computing.
6.2 Data fitting and regression
| In regression models, one often minimizes an objective proportional to \(\|Ax-b\|_2\). QR-based solvers avoid the numerical drawbacks of normal equations and deliver stable coefficient estimates. This is especially valuable when predictors are correlated or when scaling varies across features. |
|---|
6.3 Signal processing and orthogonal projections
Many signal processing tasks rely on orthogonal projections, filtering, and subspace computations. Because QR produces an orthonormal basis for the column space of \(A\), it helps represent projected components cleanly and compute residual parts corresponding to modeling error or noise.
6.4 Numerical methods for optimization subproblems
Optimization algorithms frequently require solving least-squares or linearized subproblems, such as in Gauss–Newton or trust-region methods. QR serves as a dependable mechanism for computing search directions, handling constraints via transformations, and stabilizing computations when the Jacobian matrix is ill-conditioned.
6.5 Integrating QR into broader algorithms (overview)
In larger numerical workflows, QR is often used as a building block:
- as a preprocessing step for solving iterative methods more stably,
- as a core routine in eigenvalue-related computations (where QR appears in related orthogonal transformations),
- as part of matrix rank estimation pipelines.
Because QR is orthogonality-based and algorithmically modular, it integrates well into modern linear algebra libraries and higher-level computational tools.