A Sketched Generalized Krylov Subspace Method for Large-Scale Regularization

Abstract: The generalized Krylov subspace (GKS) method is an effective projection technique for large-scale Tikhonov regularization with a general regularization matrix. As the subspace expands, however, two computational bottlenecks limit scalability: the thin QR factorizations of the tall projected matrices formed by the forward operator and the regularization matrix applied to the basis, and the full reorthogonalization of each new basis vector against all previous columns. We propose a sketched variant, named sGKS, that addresses both bottlenecks. The QR factorizations are performed on compressed matrices of much smaller row dimension, maintained incrementally via rank-one updates. Moreover, we observe that explicit reorthogonalization can be skipped entirely without compromising the quality of the approximation subspace, since no step of GKS relies intrinsically on the orthogonality of the basis. The resulting algorithm is independent of the choice of sketching operator and preserves the approximation quality of the original method: we show that, in the absence of sketching in the projected solve, sGKS produces iterates identical to those of standard GKS, and that the sketched projected solve delivers quasi-optimal residual norms controlled by the embedding quality. For more challenging problems where the loss of basis orthogonality becomes significant, we show that incorporating a small number of iterative refinement steps in the projected solve restores the spectral properties of the basis and recovers the full accuracy of the unsketched method. Numerical experiments on image deblurring, X-ray computerized tomography, seismic travel-time tomography, and dynamic computerized tomography demonstrate that sGKS matches the reconstruction quality of standard GKS while significantly reducing per-iteration costs and overall wall-clock time.

• The paper introduces sGKS, which uses randomized sketching to reduce costly QR factorizations and eliminate explicit basis reorthogonalization in large-scale regularization.
• The method achieves near-optimal residuals and reconstruction quality with significantly lower per-iteration computational cost compared to traditional GKS.
• Experimental results in image deblurring, CT, and seismic tomography validate sGKS’s scalability, effective error control, and minimal need for iterative refinement.

Sketched Generalized Krylov Subspace Method for Large-Scale Regularization

Introduction

The paper "A Sketched Generalized Krylov Subspace Method for Large-Scale Regularization" (2606.18073) introduces sGKS, a sketched algorithmic framework for efficiently solving large-scale, ill-posed least squares problems via Tikhonov regularization with general (non-identity) regularization operators. The motivation stems from the computational bottlenecks in the Generalized Krylov Subspace (GKS) method, specifically the high cost of repeated QR factorizations of projected matrices and frequent reorthogonalization of the basis, which inhibit the scalability of GKS for large and high-dimensional inverse problems. By leveraging randomized sketching, sGKS achieves a substantial reduction in both per-iteration cost and the overall computational burden, while carefully analyzing the preservation of solution quality and convergence properties.

Background: GKS and Computational Limitations

The GKS method projects Tikhonov-regularized least-squares problems

$$\min_{x \in \mathbb{R}^n} \|A x - b\|_2^2 + \lambda \|L x\|_2^2$$

onto a sequence of nested subspaces generated by iterated gradients of the Tikhonov functional, enabling flexible and powerful regularization with arbitrary $L$. Each iteration requires the solution of a reduced projected Tikhonov problem and the expansion of the search subspace. The primary computational cost arises in the incremental QR factorization of the projected forward operator and regularization operator, as well as in maintaining an orthonormal basis via full or double reorthogonalization---a cost that grows as $\mathcal{O}(m d_k + p d_k)$ per iteration for ambient dimensions $m, p, d_k$. As the dimension $d_k$ increases, these bottlenecks dominate total runtime, especially when $A$ and $L$ are large and structured (e.g., in tomography or image deblurring).

Randomized Sketching for Dimensionality Reduction

The sGKS method applies randomized subspace embeddings (sketching) to reduce the ambient dimension of all costly operations. A sketching operator $S \in \mathbb{R}^{s \times n}$, typically constructed as a subsampled randomized trigonometric transform (SRTT), projects vectors and matrices to lower-dimensional spaces while approximately preserving inner products up to user-defined accuracy $\varepsilon$. Thin QR factorizations of projected forward and regularization matrices are then performed on compressed versions, resulting in $\mathcal{O}(s d_k)$ cost per iteration for $L$0. This sketch-and-solve strategy is agnostic to the embedding used and generalizes prior work on randomized Golub-Kahan methods to the setting of arbitrary $L$1. The entire basis construction phase is also decoupled from the need for strict orthogonality; explicit reorthogonalization is bypassed as the subspace expansion does not rely on the orthonormality property in the GKS setting.

Formulation of sGKS and Theoretical Properties

At each iteration, sGKS:

1. Expands the basis by appending the Tikhonov gradient without explicit reorthogonalization.
2. Maintains sketches of the forward and regularization projected matrices for incremental QR updates, thus reducing memory and arithmetic complexity.
3. Solves the projected sketched Tikhonov-regularized least-squares problem at each step.

The authors provide a formal equivalence result: if sketching is not applied in the projected subproblem, sGKS iterates match GKS exactly, assuming initial subspace equality. When sketching is employed, sGKS produces solutions whose residuals satisfy quasi-optimality: the deviation in the objective from the unsketched iteration is controlled by the sketch accuracy parameter, with relative error guarantees derived from subspace embedding properties. Although basis orthogonality is not preserved and numerical rank can deteriorate for ill-conditioned projected problems, the loss of solution accuracy can be fully restored by minimal iterative refinement (e.g., one FOSSILS-based refinement step), as demonstrated in the dynamic tomography experiments.

Numerical Experiments

Extensive experiments are conducted on canonical large-scale inverse problems:

• Image Deblurring: For $L$2 images, sGKS achieves virtually identical reconstruction error to GKS for a wide range of sketch sizes, even at high noise, and matches the accuracy with significantly lower iteration costs.
• X-ray CT: On underdetermined $L$3-dimensional problems, sGKS reconstructs Shepp-Logan phantoms at the same error floor as GKS (RRE $L$4 vs. $L$5) but with wall-clock time and per-iteration cost reduced by a factor reflecting the $L$6 ratio.
• Seismic Travel-Time Tomography: On sparse, overdetermined systems, sGKS is compared against hybrid LSQR and its randomized variant, demonstrating superior accuracy and speedups, especially as $L$7.
• Dynamic CT: sGKS achieves a full $L$8 speedup over GKS (33s vs. 99s for 400 iterations) on a $L$9-dimensional inverse problem, and, when equipped with a single refinement, precisely overlays the GKS convergence trajectory, indicating perfect correction of sketch-induced bias.

The experiments reveal that sketch size determines the tradeoff between computational savings and error, but even modest sizes (small multiples of subspace dimension) suffice for near-optimality. The singular value decay in the computed basis and alignment of the principal subspaces are empirically confirmed, elucidating the connection between basis rank and the point of error stagnation.

Implications and Future Directions

sGKS demonstrates that randomized sketching---when combined with careful analysis of projection structure and basis expansion---enables fully scalable iterative regularization for extreme-scale inverse problems, even with general regularizers. Theoretically, this shows that explicit orthogonality can be forfeited without loss of convergence or accuracy so long as the subspace is consistently enriched and a small number of refinement steps can be deployed to recover any lost spectral properties. Practically, sGKS is naturally suited for parallelism and integration with distributed-memory solvers due to its reduced communication overhead. Future work includes the development of adaptive, iteration-wise selection of sketch size, extension to nonlinear/structured regularization such as TV or mixed norms, and high-performance, communication-avoiding implementations for 3D PDE-constrained inversion.

Conclusion

The sGKS method provides a principled, efficient, and generalizable framework for solving large-scale regularized inverse problems by removing the major algorithmic obstacles of GKS. Through randomized sketching and the elimination of superfluous orthogonality constraints, sGKS achieves optimal statistical and computational performance, supported by strong theoretical error bounds and empirical validation over a diverse suite of inverse problems. The findings establish sGKS as a viable and preferable alternative for practitioners confronting the scalability limits of projection-based regularization in large-dimensional settings (2606.18073).