---
title: Generalized Nyström Decomposition
url: https://www.emergentmind.com/topics/generalized-nystrom-decomposition
type: topic
---

# Generalized Nyström Decomposition

The generalized Nyström decomposition is a broad family of algorithms for scalable low-rank approximation of matrices, kernels, and operators. It generalizes the classical Nyström method, originally restricted to symmetric positive semidefinite (SPSD) kernel matrices, to encompass arbitrary rectangular, nonsymmetric, and even infinite-dimensional settings. Modern variants achieve near-optimal accuracy, provable error bounds, rigorous numerical stability, and efficient computational complexity for tasks such as approximate SVD, EVD, principal component analysis, compressed kernel evaluation, and matrix/tensor factorization. Central constructs include sketching matrices, sampled index subsets, and core matrices whose spectral decompositions are leveraged for interpolation.

## 1. Mathematical Formulation

The generalized Nyström decomposition approximates a given matrix $A \in \mathbb{R}^{m \times n}$ (not necessarily symmetric or square) by identifying index sets or sketching matrices that extract submatrices and project $A$ into lower-dimensional spaces. The foundation is the construction of a core matrix $A_M$ from sampled rows $I \subset \{1,\ldots,m\}$ and columns $J \subset \{1,\ldots,n\}$, both of size $s$, yielding $A_M = A(I,J) \in \mathbb{R}^{s \times s}$. Auxiliary matrices $C = A(:,J)$ and $R = A(I,:)$ are also defined [1305.0203].

The generalized rank-$s$ approximation has the canonical form:
\[
\widehat{A} = C \, A_M^{+} \, R
\]
where $A_M^{+}$ is the Moore–Penrose pseudoinverse of the core. For kernel matrices, this reduces to the classical Nyström formula $K \approx C \, W^{+} \, C^\top$ with $W = K_M = K(I,J)$ [1303.4207, 1111.6926]. For tensor data, multilinear generalizations invoke mode-wise sketches and core tensors in the Tucker format [2309.02877].

Generalized extensions allow for arbitrary sketching matrices $X \in \mathbb{R}^{n \times r}$, $Y \in \mathbb{R}^{m \times (r+\ell)}$ (Gaussian, SRHT, etc.), yielding the factors
\[
C = A X, \qquad R = Y^\top A, \qquad W = Y^\top A X
\]
and forming $\widehat{A} = C W^{+} R$ [2009.11392, 2410.21829].

## 2. Algorithmic Workflow and Core Variants

Generalized Nyström algorithms proceed as follows [1305.0203, 2410.21829, 2009.11392]:
- Select or generate sketching matrices (or coordinate index sets) for the column and row spaces—these may be chosen randomly (e.g., Gaussian, SRHT), adaptively, or via leverage scores [1605.07583].
- Compute the sampled factors $C$, $R$, and the core $A_M$/$W$.
- Compute the spectral (SVD/EVD) decomposition of the core.
- Interpolate singular vectors or eigenvectors back to the full space using explicit formulas, e.g.,
  \[
  \begin{align*}
  \text{Left singular}: & \quad \bar{U} = [U; F_M H \Lambda^{-1}] \\
  \text{Right singular}: & \quad \bar{W} = [H; B_M^\top U \Lambda^{-1}]
  \end{align*}
  \]
  for SVD [1305.0203].
- Form the low-rank approximation using the interpolated vectors and singular/eigenvalues.

Variants include:
- Generalized column-only sketch (GN-c): single column sketch followed by two QR decompositions and a pseudo-inverse [2410.21829].
- Recursive sampling with ridge leverage scores for kernel matrices, yielding linear scaling in $n$ and spectral-norm guarantees [1605.07583].
- Pivoted sampling using RRQR for well-conditioned cores [1305.0203].
- Greedy nuclear-score-based selection for sparse/structured cases [2407.01698].

## 3. Error Bounds and Theoretical Guarantees

Generalized Nyström methods enjoy explicit, tight error bounds in operator, Frobenius, and nuclear norms. For rectangular matrices $A$ and target rank $s$, under reasonable conditioning assumptions (e.g., nonsingular core), the spectral-norm error is bounded by quantities involving the singular values of $A$ and those of $A_M$ [1305.0203]:
\[
\|A - \widehat{A}\|_2 \leq \left(\frac{\sigma_{s+1}(A)}{\sigma_s(A_M)}\right) \left( \frac{\sigma_1(A)^2}{\sigma_s(A_{lg})} + 2\sigma_1(A) + \sigma_{s+1}(A) \right)
\]
where $A_{lg}$ is the best rank-$s$ SVD approximation.

Adaptive sampling schemes (with leverage scores, randomized low-rank preprocess, and residual sampling) attain relative-error bounds:
\[
E[\|A - \widehat{A}\|_F^2] \leq (1+\epsilon) \|A - A_k\|_F^2
\]
for the best rank-$k$ approximation $A_k$ [1303.4207, 1605.07583]. Recursive ridge-leverage sampling guarantees $\|K - \widetilde{K}\|_2 \leq \lambda$ with sample size $s = O(d_{\text{eff}}^\lambda \log(d_{\text{eff}}^\lambda/\delta))$ [1605.07583].

QR-based rank reduction on the core attains the optimal nuclear-norm error among all fixed-rank approximants generated from the landmarks selected [1708.03218]. For infinite-dimensional trace-class operators, tight expectation and tail bounds are derived in Hilbert–Schmidt, operator, and trace norms [2404.00960].

## 4. Numerical Stability and Computational Complexity

Numerical stability is achieved by using QR/SVD-based pseudo-inverse computations on the core, often with explicit $\epsilon$-thresholding to avoid amplification of floating-point error from ill-conditioned cores [2009.11392, 2309.02877, 2610.11493]. The error incurred from truncating small singular values is rigorously bounded.

Computational cost is minimized:
- Dense matrices: $O(mn \log n + r^3)$ for the generalized algorithm with SRHT/structured sketches [2009.11392].
- Sparse matrices: $O(\text{nnz}(A)r)$ for coordinate sampling.
- The entirety of the computation requires at most two passes over $A$ and limited memory proportional to $r$, enabling streaming and out-of-core calculation.
- For kernel matrices: $O(n s)$ kernel evaluations, $O(ns^2)$ additional runtime [1605.07583].
- For adaptive sampling: $O(m k^2 \epsilon^{-4} + T(m^2k \epsilon^{-2}))$ [1303.4207].

## 5. Extensions: Inductive Maps, Tensors, and Infinite-Dimensional Settings

The decomposition is extended to inductive kernel methods, enabling out-of-sample extension and rapid evaluation for new points [1206.4619]. After learning the dictionary $S^*$, arbitrary new data $x \in \mathbb{R}^d$ is mapped via $e(x) = [k(x,z_1),...,k(x,z_m)]^\top$ and embedded inductively as $\phi(x) = S^{*1/2}e(x)$, facilitating O(m) evaluation complexity.

The multilinear Nyström decomposition generalizes the sketch/interpolate paradigm to tensors in the Tucker format. For each mode, sketch matrices $X_k$, $Y_k$ are used, and the core tensor and mode-wise factors are computed in a single pass. Error bounds match the randomized HOSVD up to moderate constants, yielding efficient compression for $d$-way arrays [2309.02877].

Infinite-dimensional analogs for operators in separable Hilbert spaces use Gaussian process (GP) samples as sketch vectors and construct rank-constrained operator approximants $\widehat{\mathcal{A}} = \mathcal{A} \Omega G^+ (\mathcal{A}\Omega)^*$, with expectation/tail bounds matching SVD lower bounds up to dimension-dependent quality factors [2404.00960].

## 6. Adaptive Sampling, Structure-Preserving, and CUR Generalization

Generalized Nyström extends as well to CUR decompositions, selecting rows and columns adaptively based on leverage scores, residuals, or nuclear objectives [1303.4207, 2407.01698]. Greedy selection via nuclear scores maximizes trace reduction and yields exponentially tight bounds relative to determinantal point process (DPP) sampling:
\[
1 - \frac{L_K(G_k)}{D_s(K)} < \exp\left(\frac{k}{1+s}\right)
\]
For subset selection and Laplacian reduction, structure-preserving variants handle null modes and sparsity efficiently.

## 7. Applications, Empirical Results, and Error Estimation

Generalized Nyström is foundational for large-scale kernel learning, streaming PCA, adaptive beamforming, image denoising, compressed graph Laplacian computation, and surrogate model construction in scientific computing [1111.6926, 2309.02877, 2407.01698]. Empirical studies consistently demonstrate near-optimal accuracy, rapid runtimes, stability under ill-conditioning, and strong structure preservation.

To assess approximation accuracy without extra data passes, fast leave-one-out error estimators (LPO, LTO, LRO) are available. These express cross-validated error directly in terms of the already computed factors and core, and enable adaptive rank selection efficiently [2601.11493].

---

In summary, the generalized Nyström decomposition synthesizes randomized and adaptive sampling, core compression, spectral interpolation, and rigorous error control to deliver versatile, efficient, and mathematically principled low-rank approximation algorithms for matrices, kernels, tensors, and operators across symmetric and nonsymmetric, finite and infinite-dimensional contexts [1305.0203, 2009.11392, 2410.21829, 2404.00960].

Source: https://www.emergentmind.com/topics/generalized-nystrom-decomposition