A Sketched Generalized Krylov Subspace Method for Large-Scale Regularization

Technical Summary: A Sketched Generalized Krylov Subspace Method for Large-Scale Regularization

Problem Statement
The paper addresses large-scale ill-posed least-squares problems of the form $\min_{x} \|Ax - b\|_2^2$, where the matrix $A$ is severely ill-conditioned and $b$ is noise-contaminated. To stabilize the solution, Tikhonov regularization is employed, leading to the penalized problem:
$$\min_{x \in \mathbb{R}^n} \left\{ \|Ax - b\|_2^2 + \lambda \|Lx\|_2^2 \right\},$$
where $L$ is a general regularization matrix encoding prior information. While efficient methods exist for standard-form Tikhonov ($L=I$), handling general $L$ is computationally challenging. The Generalized Krylov Subspace (GKS) method [28] was proposed to handle general $L$ by expanding a search space using the gradient of the Tikhonov functional. However, as the subspace dimension grows, GKS faces two prohibitive computational bottlenecks:

1. Thin QR Factorizations: Updating the QR factorizations of the projected matrices $AV_k$ and $LV_k$ at every iteration.
2. Reorthogonalization: Explicitly reorthogonalizing new basis vectors against all previous columns to maintain numerical stability, which requires full inner product computations.

These operations scale unfavorably with iteration count and problem size, limiting the scalability of GKS for large-scale inverse problems.

Methodology: Sketched GKS (sGKS)
The authors propose sGKS, a randomized variant of the GKS method that utilizes sketching (randomized dimensionality reduction) to accelerate the most expensive steps while preserving approximation quality. The methodology rests on two primary observations:

1. Elimination of Explicit Reorthogonalization: The authors observe that the GKS algorithm does not intrinsically rely on the orthogonality of the basis vectors for its core steps (unlike classical Krylov methods for linear systems that rely on Arnoldi/Lanczos relations). Consequently, sGKS skips the explicit reorthogonalization of the basis entirely. Instead of maintaining an orthonormal basis $V_k$, sGKS maintains a non-orthonormal basis $W_k$ formed by appending the raw gradient vectors. This removes the $O(d_k^2 n)$ cost of reorthogonalization (where $d_k$ is the subspace dimension).

2. Sketch-and-Solve for Projected Problems: To reduce the cost of solving the projected Tikhonov problem and updating QR factorizations, sGKS employs sketching operators $S_A$ and $S_L$.

   • The projected least-squares problem is solved on compressed matrices $S_A A W_k$ and $S_L L W_k$ rather than the full $A W_k$ and $L W_k$.
   • The QR factorizations are performed on these sketched matrices of size $s \times d_k$ (where $s \ll m, p$) and updated incrementally via rank-one extensions.
   • The sketching operators used are Oblivious Subspace Embeddings (specifically Subsampled Randomized Trigonometric Transforms, SRTT), which depend only on the subspace dimension, not the subspace itself.

Key Contributions

• Algorithmic Design: The development of sGKS, which integrates randomized sketching into the GKS framework for general-form Tikhonov regularization.
• Theoretical Analysis:
  • Equivalence: The authors prove that if sketching is omitted in the projected solve, sGKS produces iterates mathematically identical to standard GKS (Proposition 1).
  • Quasi-Optimality: They demonstrate that when sketching is applied, the method yields quasi-optimal residual norms, with the error controlled by the embedding quality ($\epsilon$) of the sketching operators (Proposition 2).
  • Iterative Refinement: For cases where the loss of basis orthogonality becomes significant (leading to numerical rank deficiency), the paper shows that incorporating a small number of iterative refinement (IR) steps in the projected solve restores the spectral properties of the basis and recovers the full accuracy of the unsketched method.
• Computational Efficiency: The method reduces per-iteration costs by replacing expensive full-rank QR updates and reorthogonalization with cheaper operations on compressed matrices.

Numerical Results
The paper evaluates sGKS on four large-scale inverse problems:

1. Image Deblurring: sGKS with various sketch dimensions ($s$) achieves reconstruction quality comparable to standard GKS. Larger $s$ yields convergence trajectories nearly indistinguishable from GKS, while smaller $s$ results in slightly higher final errors but faster convergence.
2. X-ray CT Reconstruction: On a problem with $n=250,000$ unknowns and a large regularization operator ($p \approx 500,000$), sGKS matches the reconstruction accuracy of GKS (RRE $\approx 0.183$) while significantly reducing wall-clock time. The speedup is driven by eliminating reorthogonalization and reducing QR update costs from $O(p)$ to $O(s)$.
3. Seismic Travel-Time Tomography: On a sparse, overdetermined system, sGKS outperforms hybrid LSQR methods in both accuracy and speed. While the speedup over standard GKS is modest here (as gradient computation dominates), sGKS still provides a 10% reduction in time while maintaining superior reconstruction quality compared to LSQR variants.
4. Dynamic CT: In a space-time formulation with a massive regularization operator, the QR and reorthogonalization costs dominate the standard GKS runtime. Here, sGKS achieves a 3$\times$ speedup over standard GKS. The inclusion of iterative refinement (IR) allows sGKS to recover the exact convergence trajectory and optimal regularization parameter of the unsketched GKS, whereas sGKS without IR stagnates at a slightly higher error.

Significance and Claims
The paper claims that sGKS is a scalable alternative to standard GKS for large-scale regularization with general matrices. Its significance lies in:

• Scalability: It removes the computational bottlenecks (reorthogonalization and full-rank QR updates) that previously limited the applicability of GKS to very large problems.
• Accuracy Preservation: It maintains the reconstruction quality of the original method, with theoretical guarantees on residual norms and empirical evidence showing that iterative refinement can recover full accuracy even when basis orthogonality is lost.
• Generality: The approach is independent of the specific choice of sketching operator and applies to general regularization matrices $L$, not just standard-form problems.

The authors conclude that while the method introduces a sketching parameter ($s$) to tune the trade-off between speed and accuracy, it offers a robust framework for solving large-scale inverse problems where standard projection methods are computationally prohibitive. Future work suggested includes adaptive sketch dimension selection and extensions to nonlinear inverse problems.