[Paper Review] Faster Least Squares Approximation

Petros Drineas, Michael W. Mahoney+2 more|arXiv (Cornell University)|Oct 7, 2007

Sparse and Compressive Sensing Techniques19 references16 citations

TL;DR

This paper presents two randomized algorithms that accelerate least squares approximation by preprocessing data with the Randomized Hadamard Transform, followed by either uniform sampling of constraints or sparse random projection. The algorithms achieve relative-error approximation to the optimal least squares solution in O(nd ln d) time, significantly faster than the standard O(nd²) time for large n ≫ d.

ABSTRACT

Least squares approximation is a technique to find an approximate solution to a system of linear equations that has no exact solution. In a typical setting, one lets $n$ be the number of constraints and $d$ be the number of variables, with $n \gg d$. Then, existing exact methods find a solution vector in $O(nd^2)$ time. We present two randomized algorithms that provide very accurate relative-error approximations to the optimal value and the solution vector of a least squares approximation problem more rapidly than existing exact algorithms. Both of our algorithms preprocess the data with the Randomized Hadamard Transform. One then uniformly randomly samples constraints and solves the smaller problem on those constraints, and the other performs a sparse random projection and solves the smaller problem on those projected coordinates. In both cases, solving the smaller problem provides relative-error approximations, and, if $n$ is sufficiently larger than $d$, the approximate solution can be computed in $O(nd \log d)$ time.

Motivation & Objective

• Address the computational bottleneck of exact least squares solutions in large-scale overconstrained problems where n ≫ d.
• Develop faster algorithms that maintain relative-error approximation guarantees for the optimal solution vector and residual norm.
• Reduce the O(nd²) time complexity of traditional methods to O(nd ln d) while preserving accuracy.
• Provide theoretical runtime and error bounds for randomized algorithms using structured random projections and sampling.

Proposed method

• Apply the Randomized Hadamard Transform (RHT) to preprocess the input matrix A and vector b, reducing coherence and enabling efficient sampling.
• Use uniform random sampling of r constraints from the RHT-transformed system to solve a smaller least squares problem.
• Apply a sparse random projection matrix T with sparsity parameter q to project the transformed data into a lower-dimensional space.
• Solve the reduced least squares problem on the sampled or projected data to obtain an approximate solution vector.
• Leverage properties of the SVD and normal equations to relate the approximate solution to the true optimal solution.
• Use concentration inequalities and matrix tail bounds (e.g., Markov’s inequality) to bound the expected error in the residual and solution vector.

Experimental results

Research questions

• RQ1 Can randomized sampling and projection techniques be used to accelerate least squares approximation while maintaining relative-error guarantees?
• RQ2 What is the minimum number of sampled or projected constraints needed to achieve (1+ε)-relative error in the residual norm?
• RQ3 How does the condition number of A and the alignment of b with the column space of A affect the error in the approximate solution?
• RQ4 Can the Randomized Hadamard Transform be used to preprocess data so that uniform sampling or sparse projection yields accurate low-rank approximations?
• RQ5 What is the theoretical runtime of such algorithms, and how does it scale with n and d?

Key findings

• The proposed algorithms achieve a (1+ε)-relative error approximation to the optimal residual norm ‖Ax−b‖₂ with probability at least 0.8.
• The solution vector error is bounded by ‖x_opt − x̃_opt‖₂ ≤ √ε (κ(A)√(γ⁻²−1)) ‖x_opt‖₂, where γ is the fraction of b's norm in the column space of A.
• The algorithms run in O(nd ln d) time when n is sufficiently larger than d and less than e^d, offering a significant speedup over the standard O(nd²) time.
• The failure probability can be reduced to 1−δ by repeating the algorithm O(log(1/δ)) times, with only a constant-factor increase in runtime.
• The sparse random projection method uses a matrix T with expected O(nkq) non-zero entries, and the time complexity is bounded by O(nd) + O(nd log(nkq)) + O(kd²).
• Theoretical guarantees are established using concentration bounds and matrix perturbation analysis, particularly through Lemma 9 and Lemma 10 on expected error in projected inner products.