Multisketched Least Squares Solver

• The multisketched solver innovatively composes fast CountSketch and Gaussian sketches to reduce large least-squares problems to manageable subproblems.
• It employs diverse strategies—including fresh sketches per iteration, adaptive sizing, and averaged estimators—to enhance convergence and computational efficiency.
• The approach addresses finite-precision challenges by integrating mixed precision and iterative refinements to maintain stability even in ill-conditioned settings.

A multisketched least squares solver is a randomized solver for overdetermined least-squares problems that uses more than one sketch, either by composing different sketching operators, drawing fresh sketches across iterations, averaging independent sketch-based estimators, or solving a sequence of sketched subproblems with increasing fidelity. In the most explicit recent formulation, multisketching composes a CountSketch with a Gaussian sketch, forming a sketch operator $S = S_2 S_1$, and then solves the reduced least-squares problem by QR on the compressed system (Higgins et al., 19 Aug 2025). Taken together, adjacent lines of work on iterative sketching, adaptive sketch sizing, sketch-based preconditioning, and multilevel or distributed averaging suggest that the term denotes a family of randomized least-squares solvers rather than a single canonical algorithm (Epperly, 2023, Chen et al., 8 Sep 2025).

1. Lineage within randomized least-squares

Modern multisketched solvers sit inside the broader development of randomized numerical linear algebra for problems of the form

$$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$

Three baseline paradigms are central: sketch-and-solve, sketch-and-precondition, and iterative sketching. Sketch-and-solve compresses the problem once and solves the smaller system directly. Sketch-and-precondition computes a QR factorization of a sketched matrix $SA=QR$ and uses the resulting triangular factor to precondition the original problem. Iterative sketching instead refines an iterate through repeated sketched corrections, and in one standard form uses updates satisfying

$$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$

with $SA=QR$ (Epperly, 2023).

A key precursor is the iterative Hessian sketch for constrained least squares, which showed that one-shot sketching of both $A$ and $y$ is sub-optimal for solution approximation, while repeated Hessian-only sketching can achieve geometric error contraction with sketch dimension proportional to the statistical complexity of the solution and only a logarithmic number of iterations (Pilanci et al., 2014). Later stability work altered the interpretation of fast randomized solvers: standard sketch-and-precondition methods in their most commonly used form were shown to be numerically unstable for ill-conditioned least-squares problems, whereas appropriately implemented iterative sketching was proved forward stable when $\kappa u \ll 1$ (Meier et al., 2023, Epperly, 2023). This suggests that multisketching emerged not only as a way to reduce arithmetic cost, but also as a mechanism for navigating the tension among dimension reduction, convergence rate, and finite-precision reliability.

2. Core architectural patterns

The literature uses several distinct multisketch patterns.

Pattern	Representative construction	Source
Composed sketches	CountSketch followed by Gaussian, with $S=S_2S_1$	(Higgins et al., 19 Aug 2025)
Fresh sketches per iteration	Independent sketches in iterative sketching or iterative Hessian sketch	(Epperly, 2023, Pilanci et al., 2014)
Increasing sketch sizes	Sequence $m_1<m_2<\cdots<m_K$ of sketched subproblems	(Chen et al., 8 Sep 2025)
Averaged independent sketches	$$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$0 in distributed LS	(Garg et al., 2024)
Multilevel corrections	Level corrections $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$1 across sketch sizes	(Haas et al., 25 May 2026)

In the explicit two-stage multisketch design, the first sketch is a CountSketch with sketching dimension $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$2, and the second is a Gaussian sketch with $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$3 rows. The resulting operator is

$$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$4

Given $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$5, the solver computes

$$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$6

then forms an economy QR factorization $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$7, and finally solves

$$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$8

The rationale is asymmetric: CountSketch is extremely fast to apply but needs embedding dimension $$x = \arg\min_{y \in \mathbb{R}^n}\|Ay-b\|, \qquad A \in \mathbb{R}^{m\times n},\; m\gg n.$$9 for subspace-embedding guarantees, whereas the Gaussian sketch preserves subspace geometry with high probability for $SA=QR$0 rows but is more expensive per application. Composing them quickly reduces a large matrix and then applies the denser sketch only to the already-compressed intermediate matrix (Higgins et al., 19 Aug 2025).

A different multisketch interpretation appears in sequential schemes. SLSE-FRS constructs and solves a sequence of sketch-and-solve subproblems with increasing sketch sizes $SA=QR$1, typically with $SA=QR$2, and uses each solution to initialize the next stage. All sketches are extracted from a single initial transformation $SA=QR$3, applied once, after which later stages operate on progressively larger subsamples of the transformed data (Chen et al., 8 Sep 2025). In this form, multisketching is not a product of operators but a schedule of sketches whose fidelity increases as the estimate improves.

3. Conditioning, stability, and mixed precision

Finite-precision behavior is a defining issue for multisketched least-squares solvers. Standard sketch-and-precondition algorithms such as Blendenpik- or LSRN-style pipelines were shown to become unstable because the preconditioner $SA=QR$4 is repeatedly applied inside the iterative method; while $SA=QR$5 may be well-conditioned, $SA=QR$6 inherits the ill-conditioning of $SA=QR$7, so triangular-solve errors accumulate over iterations. A modification called sketch-and-apply forms

$SA=QR$8

solves the unpreconditioned problem $SA=QR$9 with LSQR, and recovers $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$0; under $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$1 this yields backward error comparable to LSQR on a well-conditioned matrix (Meier et al., 2023).

Iterative sketching changes the stability picture. With a careful QR-based implementation of the sketched matrix, triangular solves, and perturbation analysis via Wedin’s theorem, iterative sketching was proved forward stable: when $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$2, the forward error and residual converge to the same accuracy as classic Householder QR after $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$3 iterations (Epperly, 2023). This is particularly relevant for multisketched solvers that use fresh sketches at every iteration, because it separates algorithmic speed from the instability mechanism identified for repeated right-preconditioner application.

Mixed precision introduces another layer. In a mixed-precision sketching preconditioner, one computes

$$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$4

in sketching precision $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$5 and an economic QR factorization $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$6 in $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$7, then uses $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$8 as a preconditioner for the original tall-skinny matrix. The finite-precision analysis gives explicit conditions ensuring that the computed $$x_{i+1}=x_i+d_i,\qquad R^\top R\, d_i = A^\top(b-Ax_i),$$9 remains a high-quality preconditioner. In particular, if $SA=QR$0 and $SA=QR$1, then the computed $SA=QR$2 is an effective preconditioner and $SA=QR$3 is $SA=QR$4; if $SA=QR$5, the low precision regularizes $SA=QR$6 by inflating its smallest singular values, and $SA=QR$7 ceases to be a quality preconditioner even though $SA=QR$8 is still typically smaller than $SA=QR$9 (Carson et al., 2024).

The same analysis carries into GMRES-based iterative refinement for least squares. Using the split-preconditioned augmented system

$A$0

the paper shows that, under the same conditioning assumptions, mixed-precision sketching plus split-preconditioned FGMRES-based iterative refinement can drive both the relative forward error and the residual to the level of the working precision unit roundoff (Carson et al., 2024).

4. Adaptive, averaged, and bias-reduced multisketching

One important extension of multisketching is aggregation. In the distributed and small-space setting, each machine computes its own sketch-based estimator $A$1, and a central server forms

$A$2

The key design goal is not merely low sketch error but low estimator bias. Using a sparse leverage score sparsified sketching matrix, or LESS embedding, with sketch size $A$3 and row sparsity $A$4, the method uses one or two passes over the data and $A$5 space to construct a nearly-unbiased least-squares estimator. The theoretical guarantee is that $A$6 and $A$7, so aggregating $A$8 independent estimators yields an $A$9-approximate least-squares solution (Garg et al., 2024).

Another form of multisketching is adaptive sketch growth. For regularized least-squares problems, adaptive sketching methods for iterative Hessian sketch and preconditioned conjugate gradient begin with a very small sketch, even an initial embedding dimension equal to $y$0, and increase the sketch size only when the observed progress is insufficient. The theory shows that the required sketch size scales at most with the effective dimension, that the adaptive method converges linearly, and that for common embeddings the critical sketch sizes are $y$1 for Gaussian or sub-Gaussian sketches and $y$2 for SRHT (Lacotte et al., 2020, Lacotte et al., 2021). In this setting, a multisketched solver is effectively a self-tuning schedule over sketch dimensions rather than a fixed sketch architecture.

Regularized least-squares with many values of $y$3 yields yet another variant. Sketch-to-precondition methods for Tikhonov regularization can reuse one sketch for a set of $y$4 regularization parameters, compute preconditioners from that shared sketch, and ensure that LSQR converges in $y$5 iterations for $y$6 accuracy. A second algorithm uses a sketch of size $y$7 when the statistical dimension is much smaller than $y$8, and avoids Gram-matrix computation, which the authors describe as more stable than existing algorithms in that setting (Meier et al., 2022). This broadens the meaning of multisketching to include one-sketch-many-solves workflows.

5. High-performance implementations and empirical behavior

The most detailed high-performance implementation is a GPU CountSketch designed specifically for multisketching. Rather than relying on sparse matrix–matrix multiplication libraries, the implementation uses a custom CUDA kernel in which each thread handles one row of $y$9 and atomically adds or subtracts that row into the target row of the sketch according to a random mapping and Rademacher sign. Row-major order is used for $\kappa u \ll 1$0 and intermediate results to optimize memory coalescing, and the CountSketched output is immediately fed to the Gaussian sketch in the multisketch pipeline (Higgins et al., 19 Aug 2025).

Empirically, this implementation reaches 50–60% of peak memory bandwidth, compared with roughly 20% for library SpMM. On large problems, multisketching is reported to be up to 77% faster than normal equations, specifically Cholesky on the Gram matrix $\kappa u \ll 1$1, for cases such as $\kappa u \ll 1$2 and $\kappa u \ll 1$3. Gaussian-only sketching is slower and can run out of memory for large $\kappa u \ll 1$4, while SRHT is less competitive on the GPU because of complex memory access and synchronization demands (Higgins et al., 19 Aug 2025).

The numerical trade-off is explicit. If $\kappa u \ll 1$5 is an $\kappa u \ll 1$6 subspace embedding, then the sketched least-squares solution $\kappa u \ll 1$7 satisfies

$\kappa u \ll 1$8

so the relative residual norm can increase by an $\kappa u \ll 1$9 multiplicative factor. At the same time, sketch-and-solve methods, including the multisketch pipeline, remain stable as long as $S=S_2S_1$0, the same regime as QR-based solvers, whereas normal equations are stable only for $S=S_2S_1$1 (Higgins et al., 19 Aug 2025).

Sequential and iterative multisketch schemes show a different empirical profile. SLSE-FRS, which integrates Sketch-and-Solve and Iterative-Sketching by solving sketched subproblems of increasing size and then refining on the full problem, is reported to achieve OLS-level error faster than the state-of-the-art methods compared in that work, converging in about half the time of IDS and about a third the time of PCG; 2 or 3 inner iterations per stage are sufficient in the experiments, and CountSketch offers lower initialization cost with similar convergence to SRHT (Chen et al., 8 Sep 2025). RCGLS extends iterative sketching into a conjugate-gradient framework, drawing a fresh sketch at each iteration, proving linear convergence in expectation, and developing an implementation that avoids the full-dimensional vector operations that make naive sketched CGLS impractical on sparse matrices (Zeng et al., 24 May 2026). These results indicate that practical multisketched solvers now range from one-shot approximate direct methods to fully iterative randomized Krylov-like algorithms.

6. Trade-offs, limitations, and open problems

Multisketching does not remove all fundamental limitations of sketched least squares. In small-space and streaming models, there are lower bounds showing that a single $S=S_2S_1$2-size sketch cannot achieve an $S=S_2S_1$3-approximation; the distributed bias-reduction strategy circumvents this only by using multiple independent estimators and averaging them (Garg et al., 2024). Likewise, solution quality is not a monotone function of adding more sketch levels. In the multilevel sketch-and-solve framework, the variance of correction terms decreases faster than that of simple SAS estimators, and with antithetic coupling the per-level variance decays as $S=S_2S_1$4 instead of $S=S_2S_1$5, yet the overall computational cost of the multilevel estimator is slightly higher than that of simple averaging, making a naive multilevel strategy unattractive for overdetermined least-squares problems (Haas et al., 25 May 2026).

Stability also remains unresolved at the strongest level. Forward stability has been established for iterative sketching under appropriate conditioning, but backward stability comparable to Householder QR is still not proved for any truly fast sketch-based least-squares solver (Epperly, 2023). Standard sketch-and-precondition remains unreliable for ill-conditioned problems because its instability arises from repeated preconditioner application rather than from poor sketch quality per se (Meier et al., 2023). Mixed precision improves throughput only in a safe conditioning regime: once $S=S_2S_1$6 exceeds $S=S_2S_1$7, preconditioner quality and iterative-refinement convergence degrade, and higher precision becomes necessary (Carson et al., 2024).

A plausible implication is that “multisketched least squares solver” names a design space rather than a settled method class. Within that space, the main axes are now clear: composition versus iteration, fixed versus adaptive sketch size, direct residual distortion versus iterative high-accuracy recovery, and throughput versus finite-precision robustness. The current literature supports multisketching as a practical route to faster least-squares computation, but it also shows that accuracy, stability, and cost depend delicately on embedding quality, conditioning, bias, and implementation details.