Source-linked AI summary
Practical sketching algorithms for low-rank matrix approximation
Joel A. Tropp, Alp Yurtsever, Madeleine Udell, Volkan Cevher
TL;DR
The paper addresses low-rank approximation when only a random linear sketch of the input matrix is available, especially in single-pass or streaming settings. It develops practical sketch-based algorithms, including structured approximations, and reports error guarantees, parameter guidance, and empirical comparisons. With k = 2r + 1 and ℓ = 4r + 2, the expected Frobenius error is at most twice the optimal rank-r error, typically using rank ≈2r and storage about 6r(m + n).
Problem
The paper studies how to obtain accurate low-rank approximations from restricted sketch access, addressing practical questions about stability, parameters, structure preservation, and informative error bounds.
Method
The paper uses a simple random sketch with range and co-range summaries, then computes low-rank and structured approximations from those summaries.
Results
The expected Frobenius error is at most twice the optimal rank-r error for k = 2r + 1 and ℓ = 4r + 2, with typical approximation rank ≈2r and storage about 6r(m + n).
Takeaways & Limitations
The methods provide practical sketching algorithms with informative error bounds and parameter guidance for single-pass and streaming matrix-approximation settings.
Takeaways & Limitations
The algorithms are intended for single-pass or streaming settings; when multiple passes are possible, recommended methods are significantly more accurate.
Abstract
from arXiv · showhide
This paper describes a suite of algorithms for constructing low-rank approximations of an input matrix from a random linear image of the matrix, called a sketch. These methods can preserve structural properties of the input matrix, such as positive-semidefiniteness, and they can produce approximations with a user-specified rank. The algorithms are simple, accurate, numerically stable, and provably correct. Moreover, each method is accompanied by an informative error bound that allows users to select parameters a priori to achieve a given approximation quality. These claims are supported by numerical experiments with real and synthetic data.
1. Motivation.
The paper studies low-rank approximation when access to a matrix is restricted to a small random linear sketch, with emphasis on practical, reliable algorithms and informative guarantees.
- Motivation: The target is a low-rank approximation whose error is comparable to the best rank-r approximation.The desired approximation may have rank exceeding r when that improves accuracy.
- Motivation: A sketch is a random linear image of the matrix with dimension much smaller than its total dimension mn.The small sketch necessarily has a substantial null space, motivating randomized sketching maps.
- Motivation: Sketching is motivated by single-pass, streaming, and ordered-update settings where the full matrix cannot be revisited or retained.The sketch can be built during one scan, updated as columns arrive, or maintained across linear updates.
- Motivation: The framework uses range and co-range information to construct approximations from random test matrices and a stored sketch.The basic computation factors the range sketch, solves a least-squares problem using the co-range sketch, and forms a rank-k approximation.
- Motivation: For k = 2r + 1 and ℓ = 4r + 2, the expected Frobenius error is at most twice the optimal rank-r error.The resulting approximation typically has rank ≈2r, with storage about 6r(m + n).
- Motivation: The paper extends the sketching framework to structured approximations, including symmetry and positive semidefiniteness, using projections onto structured matrix families.It also develops informative error bounds, parameter guidance, pseudocode, and computational-cost analyses.
- Motivation: The methods target single-pass and streaming settings rather than all low-rank approximation problems.When multiple passes are available, the paper recommends methods whose error decreases exponentially with the number of passes.
- Motivation: The paper addresses practical questions that prior work often neglected, including numerical stability, parameter selection, achievable accuracy, and predictive error bounds.Experiments report that pursuing theoretical metrics can make newer algorithms perform worse in practice.
2. Background.
The background establishes notation for fields, matrix operations, norms, rank, singular values, expectations, and random matrix ensembles used throughout the paper.
- Notation: The scalar field F is either the real field R or the complex field C; the star and dagger denote conjugate transpose and pseudoinverse.The Frobenius norm is denoted by ∥·∥F.
- Notation: A matrix has rank r when its rank does not exceed r, and ⟦M⟧r denotes an optimal Frobenius-norm rank-r approximation.The optimal approximation need not be unique.
- Approximation error: The paper introduces tail energy to quantify the Frobenius error of a best rank-r approximation.The associated singular-value identity follows from the Eckart–Young theorem.
- Notation: Expectation notation distinguishes randomness over all variables from expectation over a specified random variable.Nonlinear functions bind before the expectation.
- Random matrices: Random test matrices summarize an input matrix through left and right multiplication, whose dimensions and distributions determine potential approximation accuracy.The paper introduces a central random-matrix ensemble and notation for treating real and complex cases together.
- Random matrices: Standard normal matrices are Gaussian matrices with independent mean-zero, variance-one entries in the real case.The complex version is formed from independent real standard normal matrices as G = G1 + iG2.
3. Sketching the Input Matrix.
The sketching scheme uses independent random test matrices to store range and co-range summaries, supporting low-rank approximation, linear updates, and practical tradeoffs among accuracy, cost, and stability.
- Design parameters: The target rank r controls the desired approximation quality, while the algorithm permits rank somewhat larger than r for improved accuracy.Choosing r itself requires prior knowledge about the input matrix and is outside the paper’s scope.
- Design parameters: The sketch stores a range summary of size k and a co-range summary of size ℓ, with r ≤ k ≤ ℓ and k ≤ n, ℓ ≤ m.Larger k and ℓ improve approximations but increase storage and arithmetic.
- Constructing the sketch: The basic sketch uses independent standard normal test matrices Ω ∈ F^n×k and Ψ ∈ F^ℓ×m.The paper analyzes Gaussian matrices in detail because they support informative error bounds.
- Constructing the sketch: The range sketch Y = AΩ captures the action of A, while the co-range sketch W = ΨA captures the action of A∗; both are necessary.These summaries are the core data retained for later approximation.
- Implementation: The sketch is implemented as an abstract data type whose private test and sketch matrices support initialization, linear updates, and queries for structured approximations.The update rule preserves the sketch under A ← θA + ηH.
- Implementation: For dense A, forming a Gaussian sketch costs Θ((k + ℓ)mn) flops, while for sparse A it costs Θ((k + ℓ) nnz(A)) flops.The stored sketch contains mk + ℓn floating-point numbers.
- Choosing distributions: Alternative test-matrix distributions trade randomness, storage, arithmetic, numerical stability, and error-bound quality differently.Orthonormal matrices can improve errors and stability at large dimensions, while ultra-sparse Rademacher matrices reduce costs but are less reliable.
- Choosing distributions: Except for ultra-sparse Rademacher matrices, the listed distributions often behave similarly to Gaussian matrices in practice.An exhaustive comparison of distributions is outside the paper’s scope.
4. Low-Rank Approximation from the Sketch.
The paper reconstructs a low-rank approximation from two random sketches using a numerically stable factorization, then derives error bounds and spectrum-dependent parameter choices.
- Parameter choices: Theorem 4.3 supports concrete guidance for selecting (k, ℓ) across applications and predicts performance for many input matrices.The paper presents these bounds as a key contribution and summarizes the recommendations in Table 1.
- Basic algorithm: The method factors the range sketch Y = QR, discards R, solves X := (ΨQ)†W using the co-range sketch, and returns  = QX.When ℓ≫k, ΨQ is well-conditioned, so X can be computed accurately through least squares.
- Basic algorithm: The formulation avoids the badly conditioned least-squares problem associated with an algebraically equivalent prior-work approximation.The paper distinguishes algebraic equivalence from numerical equivalence.
- Error analysis: Theorem 4.3 provides an accurate Frobenius-norm error bound whose minimum captures decay in the input matrix’s tail energy.The analysis identifies automatic exploitation of spectral decay as a source of empirical performance.
- Parameter choices: For k ≈2r, the method constructs a rank-k approximation with nearly the quality of a best rank-r approximation, balancing sketch size and approximation quality.The specific choice is k = 2r + α and ℓ = 2k + α.
- Parameter choices: Parameter recommendations depend on spectral behavior: k♮, ℓ♮ is a broad single choice, while k♯, ℓ♯ is unwise without sharp spectral decay.The paper also recommends using spectrum-aware choices when singular values are flat, slowly decaying, or rapidly decaying.
5. Low-Rank Approximations with Convex Structure.
The paper adds symmetry and positive-semidefiniteness by projecting an initial sketch-based approximation onto closed convex structure sets, while retaining error guarantees and controlled costs.
- Structure via projection: The basic approximation has no special properties beyond rank, so structured problems require reforming it to satisfy constraints such as symmetry or positive-semidefiniteness.The structured methods target approximations in the same family as the input matrix.
- Structure via projection: Projection onto a closed convex set containing A never increases Frobenius-norm error: ∥A − ΠC(Âin)∥F ≤ ∥A − Âin∥F.This fact underlies the structured approximation guarantees.
- Symmetric approximation: For conjugate-symmetric inputs, the method projects QX onto the symmetric set and represents the result as Âsym = USU*.The resulting approximation has rank 2k in the stated algorithmic form.
- Error and cost: Algorithm 5 uses Θ(kn) working storage and Θ(kℓn) arithmetic, while Algorithm 6 has similar costs.For Algorithm 6, both costs are explicitly Θ(kn) storage and Θ(kℓn) arithmetic.
- Error and cost: The symmetric and PSD methods preserve the low-rank error analysis because they combine Theorem 4.3 with the convex-projection error fact.The same parameter recommendations remain valid for the symmetric method, and the PSD method has an identical error bound.
- PSD approximation: For psd inputs, the method symmetrizes, computes an eigendecomposition, zeros negative eigenvalues, and returns Â+ = (UV)D+(UV)*.This produces a rank-(2k) positive-semidefinite approximation.
6. Fixed-Rank Approximations from the Sketch.
The fixed-rank procedure projects a sketch-based approximation onto rank-r matrices, preserving structural properties while retaining error guarantees tied to the best rank-r approximation.
- General fixed-rank approximation: A good initial approximation can be replaced by a fixed-rank matrix without significantly increasing the error.The projection preserves conjugate symmetry and positive semidefiniteness when those structures are imposed.
- General fixed-rank approximation: The Frobenius-norm error satisfies ∥A − JÂ_inK_r∥_F ≤ τ_{r+1}(A) + 2∥A − Â_in∥_F.The bound combines the tail energy beyond rank r with twice the initial approximation error.
- General fixed-rank approximation: The fixed-rank approximation can be computed as QJXK_r rather than truncating the full approximation QX.Because X is much smaller than the full m × n approximation, this formulation is more computationally efficient.
- General fixed-rank approximation: The fixed-rank method has storage Θ(k(m + n)) numbers and arithmetic cost Θ(kℓ(m + n)) flops.Its additional truncated-SVD and matrix–matrix multiplication work costs Θ(k^2n) flops.
- Error bounds: The method automatically exploits spectral decay while retaining the parameter recommendations used for the basic sketching method.The paper identifies this as the first theory providing these benefits for the fixed-rank approximation.
- Structured approximations: Structured variants produce fixed-rank conjugate-symmetric and positive-semidefinite approximations with error bounds identical to Corollary 6.4.Both variants use truncated eigenvalue decompositions; their storage is Θ(kn) and arithmetic cost is Θ(kℓn).
7. Computational Experiments.
The experiments evaluate sketch-based approximations across synthetic and real matrices, comparing methods, structured variants, and parameter choices under repeated randomized trials.
- Experimental protocol: 20 independent trials per input matrix produce each reported average relative error.The same matrix is retained while the complex standard normal test matrices are redrawn.
- Input matrices: The test suite includes matrices with low-rank, noisy, polynomial-decay, exponential-decay, and real-data spectra.The real-world Data matrix is complex, positive semidefinite, dimension n = 25,921, and exact rank 250.
- Comparisons: For rank-1 approximation of Data, Algorithm 7 achieves relative errors 3–6 orders of magnitude better than any competitor.Across matrices with good low-rank approximations, Algorithm 7 dominates alternatives when storage is adequate.
- Comparisons: For poorly approximated matrices, Boutsidis et al. perform best especially with small storage, while Algorithm 7 remains competitive and reliable overall.Algorithm 7 has the best all-around behavior, whereas the Boutsidis et al. method is strongest on poor low-rank approximations.
- Structured approximations: Structured Algorithms 8 and 9 always improve on Algorithm 7 and match or exceed the Boutsidis et al. method across all examples.The advantage is most significant for matrices with poor low-rank approximations.
- Parameter selection: Parameter choice (4.9) performs reasonably well universally, while matching choices (4.7), (4.9), and (4.10) to spectral decay nearly attains oracle performance.The same tuned recommendations also provide near-oracle performance for the structured approximations.
Appendix A. Analysis of the Low-Rank Approximation.
The appendix bases its analysis on classical expectation formulas for standard normal matrices, treating real and complex fields through the same underlying principles.
- The arguments require classical formulae for expectations of functions of a standard normal matrix.
- The real-case expectation results are taken from Halko et al. [19, Prop. A.1 and A.6].
- The complex case follows from the same principles, so its details are omitted.
A.1. Facts about Random Matrices.
This appendix subsection introduces the standard-normal matrix facts and notation used to support the paper’s randomized approximation analysis.
- Fact A.1 considers a standard normal matrix G ∈ F^{t×s} and conforming matrices B and C.
- The quantities α and β are specified by (2.2), while the function f is introduced in (2.3).
- The analysis also relies on randomized low-rank approximation results from Halko et al. [19, Sec. 10].Those results are stated using the notation established in Sections 3 and 4.
A.2. Results from Randomized Linear Algebra.
The subsection imports randomized linear-algebra bounds for constructing approximation subspaces and develops a second Frobenius-norm error analysis under explicit distributional assumptions.
- Fact A.2 applies to A ∈ F^{m×n}, a natural number ϱ with ϱ < k − α, and a standard-normal test matrix Ω.
- The matrix Q computed by (4.1) is the randomized subspace object whose behavior is covered by the imported result.
- The cited result follows from Halko et al. [19, Thm. 10.5], with Fact A.1 handling both real and complex cases.
- A second Frobenius-norm error bound for approximation (4.3) is established while distributional assumptions on test matrices are stated explicitly.
A.3. Proof of Theorem 4.3: Frobenius Error Bound.
The proof develops the approximation framework while numerical figures compare fixed-rank sketching methods by storage cost and relative error.
- Figure 7.2 compares Algorithm 7 with alternative literature methods for fixed-rank approximation as storage cost varies.The comparison uses oracle performance and relative error for rank r = 5 unless otherwise stated.
- Figure 7.3 compares unstructured, conjugate symmetric, and positive-semidefinite fixed-rank approximations under varying storage costs.All three methods are evaluated through oracle relative error for the same input-matrix settings.
- The proof’s main strategy is to establish that X approximates Q∗A, while QX approximates A through the relation A ≈ Q(Q∗A).
A.3.1. Decomposition of the Approximation Error.
The approximation error is decomposed by expressing the approximation as ˆA = QX and separating projection error from the error in estimating Q∗A.
- Lemma A.3 defines the approximation as ˆA = QX and decomposes its approximation error into two parts.
- The proof uses the Pythagorean theorem to justify the error decomposition.
- A matrix P with orthonormal columns is constructed to represent the complementary subspace needed for analyzing X ≈ Q∗A.
A.3.2. Approximating the Second Factor.
The proof represents the sketch using orthogonal components, assumes Gaussian test matrices for the expectation analysis, and derives a Frobenius-norm error bound.
- The complementary projector satisfies P P∗ = I − QQ∗, and Ψ1 = ΨP and Ψ2 = ΨQ separate the two subspaces.
- The sketch decomposes as W = Ψ1(P∗A) + Ψ2(Q∗A), enabling recovery of Q∗A by left-multiplying with Ψ† when Ψ2 has full column-rank.
- For a Gaussian test matrix Ψ independent of Ω, X is an unbiased estimator of Q∗A, and Ψ2 has full column-rank almost surely when ℓ ≥ k.
- The expectation and expected squared Frobenius-norm analyses use independence and the zero mean of Ψ1 to control the approximation error.
- The proof assumes Ω and Ψ are independent standard normal matrices and concludes with the bound F ≤ (1 + f(k, ℓ)) · (1 + f(ϱ, k)) · τ 2.