PaperlayerView exact revision on arXiv

Source-linked AI summary

Improving CUR Matrix Decomposition and the Nyström Approximation via Adaptive Sampling

Shusen Wang, Zhihua Zhang

arXiv:1303.4207v7cs.LGmath.NA

TL;DR

Low-rank approximation methods face challenges with structural preservation, subset selection, and storing dense factors for large matrices. This paper develops adaptive-sampling error bounds and corresponding CUR and Nyström algorithms, reporting improved approximation guarantees and sample complexity.

  • Problem

    Selecting effective column subsets is difficult, while column-based decompositions can require infeasibly dense coefficient matrices for large sparse data.

  • Method

    The paper establishes a general adaptive-sampling error bound and applies it to develop CUR and modified Nyström approximation algorithms.

  • Results

    The proposed methods significantly improve on classical counterparts, while the Nyström analysis gives a constant-factor Frobenius-norm relative-error bound.

  • Takeaways & Limitations

    The CUR algorithm achieves relative-error guarantees using c = 2kϵ^-1(1 + o(1)) columns and r = cϵ^-1(1+ϵ) rows, fewer than the cited prior bounds.

  • Takeaways & Limitations

    The paper’s lower-bound discussion does not apply when the intersection matrix is not W†.

Abstract

from arXiv · show

The CUR matrix decomposition and the Nyström approximation are two important low-rank matrix approximation techniques. The Nyström method approximates a symmetric positive semidefinite matrix in terms of a small number of its columns, while CUR approximates an arbitrary data matrix by a small number of its columns and rows. Thus, CUR decomposition can be regarded as an extension of the Nyström approximation. In this paper we establish a more general error bound for the adaptive column/row sampling algorithm, based on which we propose more accurate CUR and Nyström algorithms with expected relative-error bounds. The proposed CUR and Nyström algorithms also have low time complexity and can avoid maintaining the whole data matrix in RAM. In addition, we give theoretical analysis for the lower error bounds of the standard Nyström method and the ensemble Nyström method. The main theoretical results established in this paper are novel, and our analysis makes no special assumption on the data matrices.

1. Introduction

The paper motivates representing large matrices with actual columns and rows rather than SVD bases, whose limitations include poor interpretability, sparsity loss, and infeasible storage. It develops an adaptive-sampling framework and improved CUR and Nyström methods with relative-error guarantees, lower sampling or storage demands, and new lower-bound analyses.

  • Motivation: SVD-based approximations can lose sparsity and produce basis vectors with limited concrete meaning, motivating representations built from actual matrix columns and rows.These limitations make decomposition and interpretation challenging for large sparse matrices.
  • Column selection: Selecting the best subset of columns is hard because there are (n choose c) possible choices.Approximate polynomial-time methods therefore seek multiplicative error bounds, including relative-error bounds when the factor is 1+ϵ.
  • CUR motivation: CUR avoids the dense coefficient matrix X that can make CX decompositions infeasible to compute and store when matrix dimensions are near equal.Existing CUR methods can nevertheless require many sampled rows and columns; subspace sampling uses O(kϵ−2 log k) columns and O(kϵ−4 log2 k) rows for relative-error bounds w.h.p.
  • Technical contribution: Theorem 5 extends adaptive sampling analysis from projection onto one space to simultaneous projection onto column and row spaces.The bound supports relative-error CUR and Nyström guarantees regardless of which relative-error column-selection algorithm is used.
  • Proposed algorithms: The proposed CUR algorithm requires much fewer columns and rows than subspace sampling while achieving a relative-error bound.The paper reports that the method outperforms subspace sampling theoretically and empirically.
  • Proposed algorithms: Sampling 2k columns enables the modified Nyström method to achieve a relative-error bound in expectation for target rank k.The paper also establishes lower error bounds for standard and ensemble Nyström methods and compares them with the modified method’s upper bound.

2. Notation

This section establishes notation for matrices, low-rank approximations, norms, leverage scores, pseudoinverses, projections, and computational costs. It also defines the SVD-based quantities used throughout the paper.

  • Matrix notation: The notation defines identity matrices, vectors of ones, zero vectors or matrices, matrix rows, columns, and contiguous column submatrices.For A ∈ R^(m×n), a(i) denotes row i, a_j column j, and A_i:j the columns from i through j.
  • Low-rank approximation: For rank ρ = rank(A) and k ≤ ρ, the truncated SVD components define A_k, the best rank-k approximation, and σ_i(A), the i-th largest singular value.For SPSD matrices, the SVD coincides with the eigenvalue decomposition and U_A = V_A.
  • Matrix norms: The section defines Frobenius, spectral, and nuclear norms, and uses ||·||_ξ to denote any of these three norms.The spectral norm equals σ_A,1, while the nuclear norm is the sum of the singular values.
  • Leverage scores: SVD-based leverage scores quantify the influence of columns and rows in rank-k approximation, with higher scores indicating greater influence.Column scores are defined from the right singular vectors and row scores from U_A,k.
  • Pseudoinverses, projections, and costs: The Moore-Penrose inverse supports projections onto the column space of X and row space of Y, while full and truncated SVDs cost O(mn^2) and O(mnk) flops, respectively.Computing A† also costs O(mn^2) flops for an m×n general matrix with m ≥ n; matrix multiplication can be parallelized.

3. Previous Work

Previous work established adaptive sampling and near-optimal column selection methods that underpin the paper’s CUR and Nyström algorithms. It also developed deterministic and randomized CUR methods and the first known relative-error algorithm for standard Nyström approximation.

  • Adaptive sampling: Adaptive sampling selects additional columns according to the residual after an initial subset, with guaranteed error reduction; the paper generalizes its bound in Theorem 5.The earlier adaptive-sampling lemma is presented as a special case of the paper’s more general result.
  • CUR algorithms: Prior CUR work included Stewart’s deterministic SCRA algorithm, which was effective but time-expensive and lacked a data-independent error bound.SCRA selects columns and rows through truncated pivoted QR decompositions.
  • CUR algorithms: Drineas et al. introduced a two-stage randomized CUR algorithm with a high-probability relative-error bound based on leverage-score sampling of columns and rows.Its running time is dominated by the truncated SVD of A, O(mnk).
  • Nyström approximation: Gittens and Mahoney established the first and only known relative-error algorithm for the standard Nyström method, using subspace sampling and scaled columns.Their result concerns relative error in nuclear norm for SPSD matrices.

4. Main Results

The paper establishes a general adaptive-sampling error bound that jointly supports CUR and modified Nyström approximations. It also derives lower bounds for conventional Nyström methods and extends expected relative-error guarantees to high-probability bounds.

  • Adaptive sampling bound: Theorem 5 generalizes adaptive sampling by bounding projection error from selected columns and rows simultaneously.This joint guarantee addresses the combined column-and-row selection required by CUR and Nyström approximation.
  • Adaptive sampling bound: Corollary 7 applies the adaptive bound to any relative-error column-selection algorithm, yielding expected relative-error guarantees for CUR and modified Nyström methods.For CUR, additional rows are adaptively sampled from the residual; for modified Nyström, additional columns are sampled adaptively.
  • CUR approximation: The adaptive CUR algorithm combines near-optimal column selection with adaptive row sampling to obtain a tighter theoretical bound than existing algorithms.It also avoids loading the entire data matrix into RAM and can achieve linear-time scaling in appropriate single- or multi-processor settings.
  • Modified Nyström approximation: The modified Nyström method provides a constant-factor Frobenius-norm relative-error bound for symmetric matrices without requiring positive semidefiniteness.Using the less efficient optimal column-selection algorithm further reduces the required column count.
  • Lower bounds: For conventional Nyström methods, there exist SPSD matrices whose spectral, Frobenius, and nuclear-norm approximation errors obey the lower bounds in Table 3 under any sampling strategy.The spectral- and Frobenius-norm lower bounds can depend on matrix dimension m, so they are not constant-factor bounds.
  • High-probability guarantees: Expected relative-error bounds extend to high-probability guarantees through Markov’s inequality and repeated sampling with selection of the best sample.The failure-related quantity decays exponentially with the number of repetitions t.

5. Empirical Analysis

Empirical results show that adaptive sampling improves approximation accuracy over competing methods in CUR and Nyström settings, with efficiency depending on sample size and matrix structure. The modified Nyström intersection matrix consistently reduces error but requires additional computation.

  • CUR matrix decomposition: Adaptive CUR sampling achieves much lower approximation error than subspace sampling across all tested cases.It outperforms SCRA on Farm Ads and Gisette, underperforms SCRA on Enron, and is comparable on Dexter.
  • CUR matrix decomposition: Subspace and adaptive CUR sampling are much more efficient than SCRA, especially when c and r are large.Adaptive sampling is comparable to subspace sampling for small c and r but becomes less efficient as they grow.
  • Nyström approximation: Adaptive and subspace Nyström sampling significantly outperform uniform sampling when leverage-score heterogeneity is large.Their advantage over uniform sampling diminishes when leverage scores are more homogeneous.
  • Nyström approximation: For reasonably small c, specifically c < 10k, adaptive and subspace sampling outperform uniform sampling by selecting better basis columns.The experiments identify this effect across the tested Nyström settings.
  • Nyström approximation: The modified Nyström method always yields much lower error than the standard method, but its intersection matrix costs more time to compute.The standard method uses U = W† with O(c3) cost, whereas the modified method uses U = C†A(C†)T with O(mc2)+TMultiply(m2c) cost.

6. Conclusion · Appendix A. The Dual Set Sparsification Algorithm

The paper introduces a more general relative-error bound for adaptive sampling and uses it to develop improved CUR and Nyström algorithms. Appendix A presents the deterministic dual set sparsification algorithm and its implementation and runtime analysis.

  • 6. Conclusion: The adaptive sampling analysis yields novel CUR and Nyström algorithms with significant improvement over classical counterparts.The conclusion identifies the relative-error bound as the basis for these algorithms.
  • 6. Conclusion: The relative-error CUR algorithm selects c = 2kϵ−1(1 + o(1)) columns and r = cϵ−1(1+ϵ) rows from the original matrix.
  • Appendix A. The Dual Set Sparsification Algorithm: Appendix A includes the deterministic dual set sparsification algorithm as an important step in the near-optimal column selection algorithm.The algorithm is attributed to Boutsidis et al. (2011) and is presented in Algorithm 3 with bounds in Lemma 14.
  • Appendix A. The Dual Set Sparsification Algorithm: Given r with k < r < n, Algorithm 3 deterministically computes nonnegative weights with at most r nonzero entries.The inputs include columns of an arbitrary matrix X and a decomposition of the identity.
  • Appendix A. The Dual Set Sparsification Algorithm: Each iteration performs an eigenvalue decomposition of the SPSD matrix Aτ, enabling efficient computation of the required matrix expression and φ(L, Aτ).The eigenvalue decomposition also provides the eigenvalues used to compute φ(L, Aτ) directly.
  • Appendix A. The Dual Set Sparsification Algorithm: The algorithm runs in r iterations, with total running time O(rnk2) + TMultiply(nl).This combines O(k3) eigenvalue decompositions, O(k2) comparisons across n items, and computation of ∥xi∥2 for each xi.
  • Appendix A. The Dual Set Sparsification Algorithm: The near-optimal column selection algorithm combines randomized SVD via random projection, dual set sparsification, and adaptive sampling.The appendix separately identifies the costs of these three steps and gives their total complexity.

Appendix B. Proofs of the Adaptive Sampling Bounds

Appendix B provides proofs for four main results on adaptive sampling bounds, organized across four subsections.

  • Appendix B. Proofs of the Adaptive Sampling Bounds: Proofs of Theorems 5, 8, and 10 and Corollary 7 appear in Appendices B.1, B.3, B.4, and B.2, respectively.The appendix is divided into four corresponding subsections.

B.1 The Proof of Theorem 5

Theorem 5 is proved through its equivalent formulation, Theorem 15, which establishes an adaptive sampling bound by constructing sampled vectors within the span of the selected columns. The proof relies on a lemma showing that top right singular vectors of AR†R are annihilated by A − AR†R.

  • Theorem 15: Theorem 5 is equivalently expressed as Theorem 15, so the proof establishes Theorem 15 under the paper’s column-space convention.The passage explicitly states that the authors prove Theorem 15 instead of Theorem 5.
  • Theorem 15: Theorem 15 analyzes adaptive sampling after selecting c1 columns, forming residual B = A − C1C†, and sampling c2 additional columns independently with probabilities p_i.The resulting matrix C contains all c1 + c2 selected columns of A, and the theorem’s expectation is taken over C2.
  • Proof construction: The proof constructs random vectors as linear functions of sampled columns and shows that they lie in span(C1) + span(C2) = span(C).Consequently, the column space of the constructed random matrix W is contained in span(C).
  • Proof completion: The argument completes the bound by placing the column space of F inside span(W) and using orthogonality between A(I − R†R) and (I − CC†)AR†R.The subsequent equalities invoke Lemma 16 and the preceding construction.
  • Lemma 16: Lemma 16 shows that each top right singular vector v_j of AR†R satisfies (A − AR†R)v_j = 0.The proof derives this from the row-space projection AR†R and the orthogonality of v_j to the complementary right-singular subspace.

B.2 The Proof of Corollary 7 · B.3 The Proof of Theorem 8 · B.4 The Proof of Theorem 10

These sections derive CUR and Nyström error bounds from projection properties, Lemma 2, and Corollary 7, then account for the computational costs of the corresponding algorithms. The proofs also establish the auxiliary projection inequality used in the CUR analysis.

  • B.2 The Proof of Corollary 7: The CUR error bound follows after verifying that rank(CC†A) = rank(C) = ρ ≤c, so Theorem 5 applies.The proof uses that C is formed from columns of A and its column space is contained in A's column space.
  • B.2 The Proof of Corollary 7: For symmetric A, Theorem 15 is applied to CT1 to derive the Nyström approximation error bound.The argument treats CT1 as consisting of rows of A and takes expectation with respect to C2.
  • B.2 The Proof of Corollary 7: Lemma 17 bounds projection terms by exploiting nested column spaces and the orthogonality of PC and ¯PC.The proof invokes the monotonicity of projection and complementary-projection norms when span(M) ⊂ span(N).
  • B.3 The Proof of Theorem 8: Theorem 8's error bound follows directly from Lemma 2 and Corollary 7.Its total complexity combines near-optimal column selection, adaptive sampling, pseudoinverse computation, and multiplication of C†AR†.
  • B.4 The Proof of Theorem 10: Theorem 10's error bound likewise follows immediately from Lemma 2 and Corollary 7.The proof selects c1 = O(kϵ−1) columns for C1 and c2 = O(kϵ−2) columns for C2 before constructing the intersection matrix.
  • B.4 The Proof of Theorem 10: Theorem 10's complexity accounts for adaptive sampling, the O(mc2) term, multiplication TMultiply(m2c), and construction of U = C†A(C†)T.These components are combined into the stated total time complexity.

Appendix C. Proofs of the Lower Error Bounds · C.1 Construction of the Adversarial Cases · C.1.1 The Adversarial Case for the Spectral Norm Bound

Appendix C constructs adversarial cases for proving lower error bounds of standard and ensemble Nyström methods. The spectral-norm case uses a positive definite matrix B whose singular-value structure is analyzed through its Frobenius and spectral norms.

  • Appendix C. Proofs of the Lower Error Bounds: Appendix C. Proofs of the Lower Error Bounds develops adversarial cases for the lower bounds of the standard and ensemble Nyström methods.The appendix states that Theorems 20, 21, 22, 24, and 25 are used to prove Theorem 12.
  • C.1 Construction of the Adversarial Cases: C.1 Construction of the Adversarial Cases considers separate constructions for spectral-norm, Frobenius-norm, and nuclear-norm bounds.The passage explicitly distinguishes the spectral norm case from the Frobenius and nuclear norm cases.
  • C.1.1 The Adversarial Case for the Spectral Norm Bound: C.1.1 The Adversarial Case for the Spectral Norm Bound constructs an m×m positive definite matrix B.The matrix construction is introduced before its properties are established in Lemma 18.
  • C.1.1 The Adversarial Case for the Spectral Norm Bound: The constructed matrix uses a parameter α ∈[0, 1), and xT Bx > 0 for every nonzero x ∈Rm.These conditions establish the stated positivity property of B.
  • C.1.1 The Adversarial Case for the Spectral Norm Bound: The proof analyzes the best rank-k approximation Bk through the squared Frobenius norm of B and the singular values of B.Because B is SPSD, the proof does not distinguish between its singular values and eigenvalues.
  • C.1.1 The Adversarial Case for the Spectral Norm Bound: The spectral norm of B is its largest singular value, and its maximum is attained at x = 1 √m1m, which is the top singular vector.The proof identifies the spectral norm as σ1 = λ1 and then studies the projection onto the subspace orthogonal to the top singular vector.
  • C.1.1 The Adversarial Case for the Spectral Norm Bound: For every j > 1, the proof characterizes the j-th eigenvalue and eigenvector of B as the corresponding singular value and singular vector.The argument concludes after establishing these properties for the remaining spectral components.

C.1.2 The Adversarial Case for The Frobenius Norm and Nuclear Norm Bounds

This section constructs a block-diagonal adversarial matrix from repeated copies of a structured matrix B and characterizes its singular values and best rank-k approximation error. The result provides the basis for Frobenius- and nuclear-norm lower-bound analysis.

  • Adversarial construction: The adversarial construction uses a p × p matrix B with unit diagonal and off-diagonal entries α, repeated across k blocks in an m × m block-diagonal matrix A, where m = kp.A is defined as A = BlkDiag(B, ..., B) with k blocks.
  • Singular-value characterization: Lemma 19 identifies Ak as the best rank-k approximation to the constructed matrix A.The lemma applies directly to the matrix A defined in (10).
  • Singular-value characterization: The top k singular values of A are equal, while the remaining m − k singular values are also equal, as specified by Lemma 19.The passage gives σ1(A) = · · · = σk(A) = 1 + pα −α and σk+1(A) = · · · = σm(A) = 1 −α.
  • Approximation error: The residual of the best rank-k approximation has squared Frobenius norm ∥A − Ak∥_F^2 = (1 −α)(m −k).This residual expression is stated in Lemma 19 for the constructed adversarial matrix.

C.2 Lower Bounds of the Standard Nystr¨om Method

This section establishes lower bounds for the standard Nyström method on structured SPSD matrices. The bounds include a family with repeated block structure whose Frobenius and nuclear norm lower bounds improve on the basic construction by a factor of k.

  • Theorem 20: Theorem 20 lower-bounds the standard Nyström approximation error for an m × m matrix with unit diagonal and off-diagonal entries α ∈ [0, 1).The constructed matrix is SPSD for the stated range of α.
  • Theorem 20: Theorem 20 derives lower bounds for the residual in the Frobenius, spectral, and nuclear norms.The proof computes the residual’s Frobenius, spectral, and nuclear norms from its singular values.
  • Theorem 21: Theorem 21 applies the construction to an m × m SPSD block-diagonal matrix and strengthens the Frobenius and nuclear norm lower bounds by a factor of k.Here, k < m is an arbitrary positive integer, and the bound follows by distributing sampled columns across the k blocks.
  • Theorem 22: Theorem 22 shows that there exists an m × m SPSD matrix for which the standard Nyström method has lower bounds in the considered norms.The construction takes α → 1, using the matrices from the preceding lower-bound results.

C.3 Lower Bounds of the Ensemble Nystr¨om Method

This section establishes lower bounds for the ensemble Nyström method under non-overlapping sampling. The bounds apply to structured SPSD matrices and show that its relative-error ratio can remain bounded below by a constructed hard instance.

  • Lemma 23: For an SPSD matrix B with unit diagonal and constant off-diagonal α ∈ [0, 1), the ensemble Nyström approximation error has a lower bound.The result assumes t non-overlapping samples, each containing c columns.
  • Theorem 24: For the structured matrix A, the ensemble Nyström method’s approximation error is likewise lower bounded when t samples of c non-overlapping columns are selected.The proof uses the block-diagonal structure of A and applies the bound from Lemma 23.
  • Theorem 25: There exists an m×m SPSD matrix A for which the ensemble Nyström method’s relative-error ratio is lower bounded.This conclusion follows from Theorem 24 and Lemma 19 by taking α → 1.
Loading 1303.4207v7…