Source-linked AI summary

LSRN: A Parallel Iterative Solver for Strongly Over- or Under-Determined Systems

Xiangrui Meng, Michael A. Saunders, Michael W. Mahoney

arXiv:1109.5981v2cs.DScs.MSmath.NA

TL;DR

The paper addresses accurate min-length least squares for strongly over- or under-determined, possibly rank-deficient systems. It introduces LSRN, which uses random normal projection preconditioning with iterative solvers, and reports predictable performance, strong comparisons with competing solvers, and parallel scalability.

  • Problem

    The paper seeks an accurate solver for strongly over- or under-determined least squares systems when A may be rank-deficient.

  • Method

    LSRN uses random normal projections to construct a preconditioner, then applies LSQR or the Chebyshev semi-iterative method to compute the min-length solution.

  • Results

    LSRN’s preconditioned system is provably well-conditioned, making iteration counts predictable; it is faster than competing solvers on reported large dense and sparse problems.

  • Takeaways & Limitations

    LSRN supports rank-deficient systems and parallel settings while automatically benefiting from sparse matrices and fast linear operators.

  • Takeaways & Limitations

    Random normal projection costs O(mn^2) on dense matrices, whereas some fast randomized transforms cost roughly O(mn log m).

Abstract

from arXiv · show

We describe a parallel iterative least squares solver named \texttt{LSRN} that is based on random normal projection. \texttt{LSRN} computes the min-length solution to $\min_{x \in \mathbb{R}^n} \|A x - b\|_2$, where $A \in \mathbb{R}^{m \times n}$ with $m \gg n$ or $m \ll n$, and where $A$ may be rank-deficient. Tikhonov regularization may also be included. Since $A$ is only involved in matrix-matrix and matrix-vector multiplications, it can be a dense or sparse matrix or a linear operator, and \texttt{LSRN} automatically speeds up when $A$ is sparse or a fast linear operator. The preconditioning phase consists of a random normal projection, which is embarrassingly parallel, and a singular value decomposition of size $\lceil γ\min(m,n) \rceil \times \min(m,n)$, where $γ$ is moderately larger than 1, e.g., $γ= 2$. We prove that the preconditioned system is well-conditioned, with a strong concentration result on the extreme singular values, and hence that the number of iterations is fully predictable when we apply LSQR or the Chebyshev semi-iterative method. As we demonstrate, the Chebyshev method is particularly efficient for solving large problems on clusters with high communication cost. Numerical results demonstrate that on a shared-memory machine, \texttt{LSRN} outperforms LAPACK's DGELSD on large dense problems, and MATLAB's backslash (SuiteSparseQR) on sparse problems. Further experiments demonstrate that \texttt{LSRN} scales well on an Amazon Elastic Compute Cloud cluster.

1. Introduction.

The paper targets accurate min-length least squares solutions for strongly over- or under-determined, possibly rank-deficient systems. It introduces LSRN, a parallel randomized solver with predictable iterative performance and strong empirical results.

  • Problem: LSRN addresses strongly over- or under-determined least squares problems without assuming that A has full rank.The target is the unique min-length solution, including cases with infinitely many minimizers.
  • Method: LSRN uses random normal projections to build a well-conditioned preconditioner for LSQR or the Chebyshev semi-iterative method.The preconditioning process is embarrassingly parallel and speeds up automatically for sparse matrices and fast linear operators.
  • Method: The Chebyshev semi-iterative method is preferred on clusters with high communication cost.The paper presents it as the more suitable iterative method in that setting.
  • Results: LSRN has predictable performance because its preconditioned system is provably extremely well-conditioned.The paper connects this conditioning guarantee to a predictable number of iterations.
  • Results: On large dense systems, LSRN is faster than LAPACK’s DGELSD for strongly over-determined problems and much faster for strongly under-determined problems.On sparse systems, it runs significantly faster than competing solvers in both regimes.
  • Scope: The paper evaluates LSRN’s handling of Tikhonov regularization and its scalability on clusters.These topics are covered in the paper’s regularization, empirical evaluation, and scalability sections.

2. Least squares solvers.

The paper reviews direct and iterative least squares solvers and randomized preconditioning methods. It motivates LSRN by combining rank-deficiency support with spectrum-independent conditioning guarantees and parallel advantages.

  • Direct methods: The SVD computes the min-length solution accurately and robustly under rank-deficiency, while QR is less expensive but slightly less robust for rank determination.Complete orthogonal factorization provides another direct approach, and normal equations are cheaper but less accurate on ill-conditioned problems.
  • Iterative methods: Krylov methods started from the zero vector converge to the min-length solution; LSQR is preferred over normal-equation CG because it is numerically more stable.The CS method and LSMR are also identified as applicable Krylov subspace methods.
  • Iterative methods: The iteration count of CG-like methods is generally hard to predict because convergence depends on the condition number of A^T A.The associated computational bound can be as difficult as solving the least squares problem itself and may require reorthogonalization.
  • Randomized methods: Earlier randomized methods reduce problem size through transforms and sampling, but their accuracy depends on sample size and some methods assume full rank.The reviewed approaches include randomized Hadamard, SRFT, and Blendenpik-based methods.
  • Trade-offs: Random normal projection costs O(mn^2) on dense matrices, compared with roughly O(mn log m) for some fast randomized transforms.The paper identifies this as a potential issue for some applications, while citing sparse and parallel settings as advantages.
  • LSRN contributions: LSRN uses random normal projections whose strong probability bounds depend only on s/r; for s = 4r, its preconditioned condition number is less than 3 with high probability.This conditioning result is independent of the spectrum of A and supports predictable iteration counts.
  • LSRN contributions: The predictable conditioning enables the CS method, which avoids vector inner products and is particularly suitable for clusters with high communication cost.The paper notes that accurate singular-value bounds are otherwise difficult to obtain for efficient CS convergence.

3. Preconditioning for linear least squares.

The paper develops preconditioning conditions and a random-normal-projection preconditioner that preserves the min-length least-squares solution while making the iterative system well-conditioned and predictable.

  • Preconditioning conditions: Nonsingular left and right preconditioners preserve consistency with the original least-squares system, while range conditions ensure the same min-length solution.The sufficient conditions are range(M)=range(A) for left preconditioning and range(N)=range(A^T) for right preconditioning.
  • LSRN construction: LSRN applies a Gaussian random projection, computes a reduced SVD, forms a preconditioner, and solves the preconditioned system iteratively.For over-determined systems, it uses G∈R^s×m with s=⌈γn⌉, while under-determined systems use G∈R^n×s with s=⌈γm⌉.
  • Complexity and implementation: LSRN supports dense, sparse, and operator representations of A, with embarrassingly parallel projection stages and near-linear shared-memory scalability.For γ=2 and ε=10^-14, the analysis gives approximately 2n+200 matrix-vector multiplications; MPI scalability remains good when the SVD is not dominant.
  • Theoretical properties: The preconditioned system has the same min-length solution as the original system almost surely.This consistency result is stated explicitly for Algorithm 1 and follows from the random projection's full-rank property with probability 1.
  • Theoretical properties: Its spectrum depends only on a random matrix of size s×r, independent of A's spectrum, enabling strong singular-value concentration and predictable Krylov iterations.The concentration result also supports efficient use of the Chebyshev semi-iterative method.
  • Complexity and implementation: Although LSRN can require more floating-point operations than a direct SVD, matrix-matrix multiplication, sparsity, and operator access can make it faster in practice.The paper reports faster performance than DGELSD on large dense strongly over-determined systems and significantly faster performance than competing solvers on sparse systems.
  • Extensions: LSRN can also be extended to certain types of Tikhonov regularization.The regularized problem is transformed into an equivalent least-squares problem suitable for the solver.

5. Tikhonov regularization.

The paper reformulates Tikhonov-regularized least squares as ordinary least-squares problems, applying LSRN directly when over-determined and using a transformed under-determined formulation when nearly square.

  • Over-determined systems: Tikhonov regularization is equivalent to an ordinary least-squares problem of size (m+n)×n.The regularization matrix W controls the penalty term and is often chosen as λI_n.
  • Over-determined systems: When m≫n, LSRN can directly solve the equivalent over-determined formulation.This regime retains the strong over-determined structure needed by the randomized method.
  • Under-determined systems: When m≪n, the equivalent formulation is nearly square, so random sampling and random projection methods like LSRN do not perform well.The paper therefore introduces an alternative equivalent formulation for this regime.
  • Under-determined systems: The alternative formulation is under-determined of size m×(m+n), allowing LSRN to compute its min-length solution.The original solution is recovered as x*=W^-1z*, assuming W^-1 is easy to apply.
  • Cross-validation: For cross-validation over W, the random projection on A is reused, while each candidate W requires a projection on W, an SVD, and predictable iterations.This separates the reusable cost from the marginal cost of evaluating each regularization choice.

6. Numerical experiments.

The experiments evaluate LSRN’s conditioning, parameter choice, solver comparisons, sparse-data behavior, and scalability across synthetic and real-world least-squares problems. LSRN delivers predictable iteration behavior, benefits from rank deficiency and sparsity, and scales effectively in parallel settings.

  • 6.2. κ(AN) and number of iterations.: The number of LSQR iterations depends mainly on r/s rather than the original condition number, making convergence predictable from γ.The theoretical upper bound closely matches observed iteration counts, while κ(AN) remains small when s is sufficiently larger than r.
  • 6.2. κ(AN) and number of iterations.: The oversampling-factor experiment identifies γ in [1.8, 2.2] as optimal for the tested 10^5 × 10^3 problem.Increasing γ raises randn, mult, and svd costs linearly, while iteration time decreases.
  • 6.4. Dense least squares.: For large dense full-rank problems, LSRN becomes faster than DGELS and DGELSD, while Blendenpik remains fastest in the reported tests.LAPACK performance decreases substantially on under-determined systems, whereas LSRN’s performance changes little.
  • 6.4. Dense least squares.: Rank deficiency reduces LSRN’s running time because it takes advantage of the lower rank and requires fewer iterations.On a 10^6 × 10^3 problem, the reported time falls from 100 to 93 seconds when the rank is reduced to 800.
  • 6.5. Sparse least squares.: On sparse problems, LSRN benefits from faster matrix operations and leads competing solvers on sufficiently large over-determined and under-determined systems.For the 10^6 × 10^3 sparse over-determined problem, LSRN takes 23 seconds; MATLAB’s backslash takes 55 seconds on the same problem.
  • 6.6. Real-world problems.: On real-world sparse problems, LSRN offers predictable performance and the min-length solution, whereas MATLAB’s performance depends heavily on sparsity patterns and may not return that solution.The comparison includes structured matrices and the unstructured tnimg data, with limitations for competing solvers on rank-deficient or under-determined cases.
  • 6. Numerical experiments.: LSRN’s preconditioning is embarrassingly parallel and automatically benefits from sparse matrices and fast linear operators.The implementation relies on matrix-vector and matrix-matrix multiplications, so sparse operations reduce runtime without changing the algorithm.

7. Conclusion.

The paper acknowledges contributions from collaborators and references assistance related to MPI barrier behavior during multithreaded SVD.

  • The authors acknowledge learning about Coakley et al.'s least-squares algorithm after completing the manuscript's initial version.
  • They thank Mark Tygert for directing them to that reference.
  • They also credit Lisandro Dalcin with an MPI Barrier implementation intended to prevent idle processes from interrupting multithreaded SVD too frequently.
Loading 1109.5981v2…