PaperlayerView exact revision on arXiv

Source-linked AI summary

Relative-Error CUR Matrix Decompositions

Petros Drineas, Michael W. Mahoney, S. Muthukrishnan

arXiv:0708.3696v1cs.DS

TL;DR

The paper addresses the difficulty of interpreting low-rank approximations whose components are linear combinations of matrix rows and columns. It develops randomized column-based and CUR decompositions using actual columns and rows, with subspace sampling. The algorithms achieve relative-error guarantees with polynomial sample sizes and runtime on the order of computing the top k singular vectors.

  • Problem

    Truncated-SVD components can be difficult to interpret, motivating low-rank approximations expressed using a small number of actual columns and rows.

  • Method

    The paper develops randomized column-based and CUR algorithms that select actual columns and rows using subspace sampling based on singular-vector subspace information.

  • Results

    The algorithms achieve relative-error guarantees with probability at least 1 − δ, use polynomially many sampled columns and rows, and run in time O(SVD(A, k)).

  • Takeaways & Limitations

    Actual columns and rows can serve in place of eigencolumns and eigenrows while providing improved interpretability in terms of the original data.

  • Takeaways & Limitations

    The column-selection result uses c > k and ε error, and does not include reference to regularization or conditioning.

Abstract

from arXiv · show

Many data analysis applications deal with large matrices and involve approximating the matrix using a small number of ``components.'' Typically, these components are linear combinations of the rows and columns of the matrix, and are thus difficult to interpret in terms of the original features of the input data. In this paper, we propose and study matrix approximations that are explicitly expressed in terms of a small number of columns and/or rows of the data matrix, and thereby more amenable to interpretation in terms of the original data. Our main algorithmic results are two randomized algorithms which take as input an $m \times n$ matrix $A$ and a rank parameter $k$. In our first algorithm, $C$ is chosen, and we let $A'=CC^+A$, where $C^+$ is the Moore-Penrose generalized inverse of $C$. In our second algorithm $C$, $U$, $R$ are chosen, and we let $A'=CUR$. ($C$ and $R$ are matrices that consist of actual columns and rows, respectively, of $A$, and $U$ is a generalized inverse of their intersection.) For each algorithm, we show that with probability at least $1-δ$: $$ ||A-A'||_F \leq (1+ε) ||A-A_k||_F, $$ where $A_k$ is the ``best'' rank-$k$ approximation provided by truncating the singular value decomposition (SVD) of $A$. The number of columns of $C$ and rows of $R$ is a low-degree polynomial in $k$, $1/ε$, and $\log(1/δ)$. Our two algorithms are the first polynomial time algorithms for such low-rank matrix approximations that come with relative-error guarantees; previously, in some cases, it was not even known whether such matrix decompositions exist. Both of our algorithms are simple, they take time of the order needed to approximately compute the top $k$ singular vectors of $A$, and they use a novel, intuitive sampling method called ``subspace sampling.''

1 Introduction

The paper develops interpretable low-rank approximations using actual columns and rows of a data matrix rather than difficult-to-interpret SVD basis vectors. Its randomized CUR algorithms provide relative-error guarantees with polynomial sample sizes and runtime comparable to computing top singular vectors.

  • 1 Introduction: Large data matrices are often compressed with truncated SVD, but its linear-combination basis vectors can be difficult to interpret in terms of the underlying data.The paper motivates using actual matrix columns and rows as interpretable alternatives.
  • 1 Introduction: CUR decompositions represent a matrix as CUR, where C and R contain actual columns and rows and U is constructed so CUR is close to A.The selected columns and rows can replace eigencolumns and eigenrows while retaining improved interpretability.
  • 1.1 Applications: DNA microarray and SNP analyses illustrate the need for bases consisting of actual genes, conditions, individuals, or SNPs rather than linear combinations of thousands of features.The paper presents CUR as a direct formulation of selecting representative SNPs and individuals.
  • 1.2 Our Main Results: The CUR result samples columns and rows with counts polynomial in k, 1/ε, and log(1/δ), while U is a weighted Moore-Penrose inverse of the C–R intersection.The algorithms run in time O(SVD(A, k)).
  • 1.3 Summary of Main Technical Result: Subspace sampling chooses columns according to row norms of the top k right singular vectors, enabling relative-error guarantees through a sampled subspace.The technical result uses random sampling to solve a generalized ℓ2 regression problem to accuracy 1 ± ε.

2 Review of Linear Algebra

This section reviews matrix norms, the SVD and best rank-k approximation, generalized inverses, projections, and sampling-matrix notation used by the paper.

  • Matrix norms: The section defines Frobenius and spectral norms and records inequalities relating them and their behavior under matrix multiplication.These norm relations support later approximation-error analyses.
  • Singular value decomposition: The SVD writes A as UΣV^T, with singular values ordered nonincreasingly and singular vectors given by the columns of U and V.The singular values and vectors provide the decomposition underlying the paper’s low-rank approximation framework.
  • Best rank-k approximation: A_k is the rank-k approximation minimizing spectral and Frobenius-norm distance among rank-k approximations to A.The paper denotes the time required to compute this best rank-k approximation by O(SVD(A, k)).
  • Generalized inverses and projections: The Moore-Penrose generalized inverse A^+ is uniquely defined by four conditions and yields the column-space projection P_A = AA^+.Weighted generalized inverses are also introduced for positive diagonal row and column weight matrices.
  • Sampling matrices: The section introduces the sampling-matrix formalism used to represent sampled and rescaled columns or rows of A.Sampling matrices select entries, while diagonal rescaling matrices compensate for the sampling probabilities; these operations produce C and R.

3 Relationship with Previous Related Work

The paper’s relative-error column and CUR methods extend prior structural and additive-error approaches, providing polynomial-time guarantees where earlier work lacked algorithms or simultaneous row-column results.

  • Numerical linear algebra: Numerical linear algebra methods used QR, quasi-Gram-Schmidt, maximum-volume selection, and rank-revealing factorizations to construct related column- and row-based decompositions.These approaches supplied structural or empirical properties, but the cited work did not provide the paper’s relative-error guarantees.
  • Alternative approximation methods: Random-projection and interpolative-decomposition methods offered alternative approximations, but their cited bounds were weaker than the paper’s relative-error bounds.The comparison specifically notes that a bound involving √m(k + 20)||A − A_k||_2 can exceed ||A||_F.
  • Additive-error methods: Prior randomized and sampling-based work mainly provided additive-error guarantees for column-based or CUR approximations.Examples include bounds of the form ||A − CX||_F ≤ ||A − A_k||_F + ε||A||_F and additive-error CUR analyses.
  • Theoretical computer science: Earlier work established existence of O(k^2/ε^2) columns for relative-error CX decompositions but offered no algorithm beyond exhaustive search.Those results also did not apply to selecting columns and rows simultaneously.
  • Theoretical computer science: The paper presents the first CUR matrix decomposition algorithm with relative error, when existence of such a representation was previously unknown.This distinguishes the contribution from earlier relative-error CX results.
  • Connections to other work: The paper’s algorithms build on related sampling and isotropic-position ideas while also addressing an open question about core sets for approximate matrix computation.The selected columns and rows can be viewed as a core set or as chosen variables and features.

4 Our Main Column-Based Matrix Approximation Algorithm

The column-based algorithm samples and rescales a small set of columns using subspace-informed probabilities, then projects A onto their span. It achieves relative Frobenius-error guarantees comparable to the best rank-k SVD approximation.

  • Algorithm: Algorithm 1 samples columns of A according to a carefully constructed nonuniform probability distribution and returns a rescaled column matrix C.The Exactly(c) and Expected(c) variants select exactly c columns or at most c columns in expectation.
  • Algorithm: O(SVD(A, k)) time suffices to compute the sampling probabilities used by Algorithm 1.The probabilities depend on approximately computed top-k right singular vectors.
  • Guarantees: 3200k^2/ε^2 columns yield the theorem's projection guarantee with probability at least 0.7 under Exactly(c), while Expected(c) uses O(k log k/ε^2) columns in expectation.The constant-probability guarantee can be amplified by independent repetitions.
  • Approximation: The approximation is A′ = CC^+A, where C^+ gives the least-squares coefficient matrix for expressing A using the selected columns.For a fixed C, C^+A is the minimum-length coefficient solution.
  • Guarantees: The rank of CC^+A may exceed k in the stated theorem, but a rank-at-most-k CX approximation within the same relative error can be obtained.The paper identifies this as a presentational choice that can be removed.
  • Sampling rationale: Subspace sampling uses row norms of the top-k right singular vectors rather than column norms of A to capture the relevant subspace.The sampling challenge is ensuring the sampled top-k right-singular-vector matrix has full rank.

5 Our Main Column-Row-Based Matrix Approximation Algorithm

The column-row algorithm first selects columns, then samples corresponding rows to form a CUR decomposition. Its guarantees preserve relative-error quality while using only sampled rows and columns.

  • Algorithm: Algorithm 2 takes A and a selected column matrix C, samples rows nonuniformly, and returns row matrix R and intersection matrix W.Both Exactly(c) and Expected(c) row-sampling variants are supported.
  • Implementation: The stated running time is O(mn), with improvements possible for sparse input matrices.The bound treats the sample counts as constants independent of m and n.
  • Guarantees: 3200c^2/ε^2 sampled rows give the theorem's guarantee with probability at least 0.7 under Exactly(c), while Expected(c) uses O(c log c/ε^2) rows in expectation.The constant-probability result can be amplified through independent repetitions.
  • Analysis: The row-sampling analysis reduces to approximating a least-squares problem by sampling rows from C and A.The sampled intersection determines U through the corresponding rescaled rows.
  • CUR construction: The CUR representation uses sampled columns C, sampled rows R, and a generalized inverse U of their intersection.The matrix U encodes inverse information from the intersection of C and R.
  • Combined result: Combining the column and row algorithms yields a CUR approximation whose failure probability can be reduced to at most δ.The construction repeats each stage ln(2/δ) times and selects the best outcome.

6 An Approximation Algorithm for Generalized ℓ2 Regression

The paper develops generalized ℓ2 regression as the technical foundation for its CX and CUR algorithms. The central questions concern selecting columns and rows that preserve relative-error residuals.

  • Motivation: Given columns C of A, the paper asks how to choose C so that ||A − CC^+A||_F is within relative error ε of ||A − A_k||_F.This connects column selection to the best rank-k SVD approximation.
  • Motivation: The paper also asks how to choose rows R and matrix U so that ||A − CUR||_F remains within relative error ε of ||A − CC^+A||_F.This separates the column-selection and row-selection approximation tasks.
  • Technical foundation: Generalized ℓ2 regression is introduced as the main technical result from which the column-based and column-row-based algorithms follow.The section presents a sampling algorithm, its approximation theorem, and the subspace-sampling probabilities.

6.1 Description of the Algorithm and Theorem

Algorithm 3 approximates generalized ℓ2 regression by sampling and rescaling rows of both the input matrix and target matrix. Theorem 5 gives probabilistic approximation bounds under rank and alignment conditions.

  • Algorithm: Algorithm 3 samples r corresponding rows of A and B, rescales them, and solves the resulting sampled regression problem.It returns an approximate residual value and coefficient matrix.
  • Theorem: Theorem 5 considers rank-no-greater-than-k A, target matrix B, optimal residual Z, and minimum-length solution X_opt = A^+B.These quantities define the original generalized regression problem.
  • Guarantees: 3200k^2/(βε^2) sampled rows give the stated bounds with probability at least 0.7 under Exactly(c), while Expected(c) uses O(k log k/(βε^2)) rows in expectation.The Expected(c) result applies under the same assumptions.
  • Guarantees: Equation (21) guarantees that substituting the sampled minimum-length solution into the original residual gives a good approximation to the original regression problem.Equations (22) and (23) provide additional bounds involving conditioning and target alignment.
  • Conditions: The relative-error guarantee requires A_k to be well-conditioned and a constant fraction of B's weight to lie in the column space of A.If B lies completely in that column space, then Z = 0, γ = 1, and the sampled solution is exact.

6.2 Discussion of the Method of “Subspace Sampling”

Subspace sampling chooses columns or rows using probabilities tied to the singular-vector subspace rather than solely to matrix entry magnitudes. This separates subspace information from singular-value scale and supports the paper’s relative-error analysis.

  • Sampling principle: Subspace sampling probabilities track the Euclidean norms of rows of the relevant left singular vectors.The probabilities are required to be close to, or not much smaller than, these row norms.
  • Sampling principle: Unlike norm-based sampling from A, subspace sampling separates information about singular-vector subspaces from information about singular-value magnitudes.Norms of rows or columns of A mix subspace information with the scale information in ΣA.
  • Interpretation: The row norms of UA,k represent how A disperses its singular-value information across the ambient constraint space.Sampling therefore biases selection toward the part of the space receiving that information.
  • Relation to regression: The method was first used in preliminary ℓ2 regression results and is contrasted with a more complex probability form for generalized regression.The alternative probabilities directly bias sampling toward where the residual matrix B has substantial weight.

6.3 Proof of Theorem 5

The proof of Theorem 5 combines sampling lemmas controlling subspace preservation, approximate matrix multiplication, and sampled residual norms. With sufficiently many subspace-sampled rows, these bounds yield the theorem’s claims, while Expected(c) sampling improves one sample requirement.

  • Sampling lemmas: r ≥400k^2/βǫ^2 gives the first lemma with probability at least 0.9 under subspace-sampling probabilities.The lemma concerns the deviation between the pseudoinverse of the sampled subspace and its transpose.
  • Sampling lemmas: Subspace sampling preserves the rank of the sampled singular-vector subspace when all singular values of SUA,k are positive.The proof concludes rank(SUA,k) = rank(UA,k) = rank(Ak).
  • Combining bounds: The proof conditions on simultaneous success of three lemmas, which occurs with probability at least 0.7 when r = 400k^2/βǫ^2.The probability follows from combining three events, each holding with probability at least 0.9.
  • Combining bounds: Equations (21) and (22) follow by combining the lemma bounds and setting ǫ′ = ǫ/8.The proof uses the assumed sample size r = 400k^2/βǫ^2.
  • Assumption and limitation: The bound for the residual-fitting term can be weak when B has most of its weight outside the column space of Ak.A stronger bound requires assuming that a constant fraction of B’s weight lies in the subspace spanned by Ak.
  • Expected(c) variant: For Expected(c) sampling, the first lemma holds with r = O(k log k/βǫ^2) rows rather than the quadratic-in-k requirement used for Exactly(c).This improvement uses a direct spectral-norm bound instead of passing through a Frobenius-norm bound.

7 Empirical Evaluation

The empirical evaluation tests CX and CUR sampling on SNP, recommendation-system, and term-document matrices using normalized Frobenius reconstruction errors. Across these domains, modest oversampling generally produces reconstruction close to or better than the best rank-k approximation baseline.

  • Experimental setup: Three datasets evaluate the methods: DNA SNP HapMap, Jester recommendation-system, and Reuters term-document matrices.The domains were chosen as representative modern scientific and internet data applications with approximately low-rank structure.
  • Experimental setup: The evaluation reports Θ1, Θ2, and Θ3 as Frobenius-error ratios normalized by ||A−Ak||F.Θ1 evaluates CX, Θ2 enforces rank no greater than k, and Θ3 evaluates CUR.
  • DNA SNP HapMap: For European HOXB data with k = 10, Θ1 = 1.36 at c = 10, Θ1 < 1.0 when c ≳17, and Θ1 = 0.62 at c = 30.Similar trends are reported for Θ2, Θ3, and the 17q25 region.
  • Jester recommendation system: For Jester column sampling with k = 15, Θ1 = 1.14 at c = 15, Θ1 ≤1 when c ≳29, and Θ1 = 0.99 at c = 30.Comparable behavior occurs when sampling rows, while CUR is noisier but remains near the baseline with modest oversampling.
  • Reuters term-document data: For Reuters, the worst reported point is Θ3 = 1.272 for k = 100 and c = 100, while Θ3 < 1.1 for c ≥300.The curves decrease more slowly with oversampling than those for the preceding datasets.

8 Conclusion

The paper presents randomized CX and CUR algorithms that express low-rank approximations through sampled columns and rows while providing relative-error guarantees. It closes by identifying open questions about determinism, stability, norms, optimal selection, and extensions.

  • Conclusion: The algorithms provide relative-error guarantees for low-rank approximations expressed using a small number of actual columns and/or rows.Their analysis uses subspace sampling and approximates a generalized ℓ2 regression problem by random sampling.
  • Open problems: Deterministic approximation algorithms and hardness results for optimally choosing columns or rows remain open problems.The paper also asks about deterministic optimization of the maximum volume of the induced parallelepiped.
  • Open problems: The paper asks whether post-sampling conditions can certify a 1+ǫ approximation and whether comparable results hold with regularization or conditioning.These questions target verification and extensions of the presented formulations.
  • Open problems: The authors identify heuristic variants, numerical stability, and generalization to other matrix norms as unresolved questions.These questions concern practical applicability and the scope of the theoretical guarantees.

A Approximating Matrix Multiplication

This section presents two complementary procedures for randomly sampling and rescaling columns and/or rows. It also uses these procedures to approximate products of two matrices.

  • Two complementary procedures randomly sample and rescale columns and/or rows from an input matrix.
  • The sampled columns and rows are used to construct approximations expressed through selected matrix elements.
  • The procedures support approximate multiplication by randomly sampling columns and rows from the input matrices.

A.1 Sampling Columns and Rows from Matrices

The paper gives two randomized methods for constructing a rescaled matrix from sampled columns, with analogous modifications for sampling rows. One samples exactly c columns, while the other uses inclusion probabilities that need not produce a fixed sample size.

  • Both algorithms take an m × n matrix A and a probability distribution {p_i} and construct a rescaled matrix from sampled columns.
  • Exactly(c): Exactly(c) selects exactly c columns in independent trials, choosing column i with probability p_i each time.
  • Expected(c): Expected(c) uses inclusion probabilities min{1, cp_i}, so the exact number of returned columns is unknown before execution.

A.2 Approximate Matrix Multiplication Algorithms

Algorithm 6 approximates AB by sampling corresponding column-row pairs and rescaling them, using either fixed-size or expected-size sampling. Its guarantees depend on norm-based, potentially nonuniform probabilities and include Frobenius- and spectral-norm bounds.

  • Algorithm 6 samples corresponding columns of A and rows of B, returning small sampled and rescaled matrices C and R.
  • Exactly(c) samples exactly c column-row pairs, whereas Expected(c) includes each pair at most once and chooses at most c pairs in expectation.
  • The algorithm forms C = ASD and R = DST B from the sampled columns of A and corresponding rows of B.
  • The approximation CR = AS D D S^T B ≈ AB uses nonuniform probabilities based on products of Euclidean norms of corresponding columns and rows.
  • When β = 1, one probability family is optimal for approximating AB by CR in the stated sense.
  • Quality guarantees: Theorems 6 and 7 provide quality-of-approximation results for Exactly(c) and Expected(c), including a new Frobenius-norm bound and a prior spectral-norm bound.
Loading 0708.3696v1…