Column Subset Selection, Matrix Factorization, and Eigenvalue Optimization

1. Introduction

• Rank-revealing QR seeks well-conditioned columns spanning the matrix’s numerical range, while broader column subset selection asks for submatrices with distinguished properties.Source: §1 · p. 2
• Stable rank is an analytic surrogate for algebraic rank and is tightly connected to how many strongly independent columns can be extracted.Source: §1 · p. 2
• Bourgain–Tzafriri guarantees a large, well-conditioned submatrix for every standardized matrix, with stable rank interpolating between extreme cases.Source: §1 · pp. 2–3

1.1. Contributions.

• Randomized polynomial-time algorithms produce the column sets guaranteed by both Kashin–Tzafriri and Bourgain–Tzafriri.The paper identifies the lack of constructive proofs as the gap and notes that no known column-selection algorithm guarantees a constant-order condition number.Source: §1.1 · p. 3
• The work also develops algorithms for Pietsch and Grothendieck factorizations, approximation algorithms for two NP-hard matrix norms, and a connection to the maxcut semidefinite program.Source: §1.1 · p. 3

1.2. Overview.

• The paper develops the Kashin–Tzafriri algorithm first, then an analogous Bourgain–Tzafriri algorithm, alongside stable-rank discussion and computational appendices.Source: §1.2 · p. 3

2. The Kashin–Tzafriri Theorem

• The proof first samples columns and then uses a factorization to identify and remove redundant columns that inflate spectral norm.The factorization exposes a subset with small norm, but the proof itself is not constructive.Source: §2 · p. 3, §2.1 · p. 4

2.1. Intuitions.

• Random sampling need not reduce spectral norm: duplicated columns in the double-identity example can leave the sampled submatrix’s norm unchanged.Source: §2.1 · p. 3

2.2. The (∞, 2) operator norm.

• Exact calculation of the (∞, 2) operator norm is computationally difficult, although the paper states it can be approximated in polynomial time up to a small relative error.For a class of positive semidefinite matrices, estimating this norm within an absolute tolerance is NP-hard.Source: §2.2 · p. 4

2.3. Pietsch factorization.

• Pietsch factorization writes B=TD with nonnegative diagonal D normalized by trace(D^2)=1 and bounds ∥T∥ by KP∥B∥∞→2, but its standard proof offers no algorithmic insight.Source: §2.3 · p. 4

2.4. Proof of Kashin–Tzafriri.

• |τ| ≥ st. rank(A) and ∥Aτ∥ < 15; the proof obtains this bound by sampling σ and then selecting τ from Aσ using the factorization.Once the Pietsch factorization is known, the column selection is straightforward, making its computation the missing algorithmic step.Source: §2.4 · p. 5

3.1. Pietsch and eigenvalues.

• Theorem 3.1 relates an eigenvalue condition on D to whether B has a factorization B = TD with ∥T∥≤α; if no D satisfies the condition, no such factorization exists.Source: §3.1 · p. 5

3.2. Factorization via optimization.

• A convex program minimizes λmax(B∗B − α²F) over nonnegative diagonal trace-one F; a nonpositive objective yields B = TD with ∥T∥≤α.Because the maximum eigenvalue is convex, it can be minimized in polynomial time.Source: §3.2 · p. 6

3.3. Other formulations.

• An alternative formulation avoids choosing α in advance, and its optimal value is α⋆², with each optimizer having trace α⋆².Source: §3.3 · p. 6
• A minimizer of the alternative formulation yields an optimal Pietsch factorization by setting D = (E⋆/trace(E⋆))^1/2 and T = BD†.Source: §3.3 · p. 6
• The dual of the alternative formulation is the maxcut semidefinite program, connecting Pietsch factorization with graph partitioning.Source: §3.3 · p. 6
• The alternative program overestimates the (∞,2) norm by at most a multiplicative factor KP², providing an approximation algorithm for that norm.Source: §3.3 · p. 7

3.4. Algorithmic aspects.

• Entropic Mirror Descent is described for the convex program; although not theoretically the most efficient approach, preliminary experiments suggest performance rivals more sophisticated techniques.For an m × s matrix B, forming B∗B costs at most O(s²m), and solving the dual program costs at most eO(s^3.5), suppressing log-like factors.Source: §3.4 · p. 7

• With probability at least 4/5, Algorithm 1 returns τ with |τ| at least st.rank(A)/2 and ∥Aτ∥≤15.Source: §4 · p. 7
• The computational cost is eO(|τ|²m + |τ|^3.5), and is sublinear in matrix size when st.rank(A)=o(n^1/3.5).Better methods for solving (3.2) would strengthen this bound.Source: §4 · p. 7
• Repeated Norm-Reduce calls bound the probability of failure before reaching stable rank by less than 0.2.Source: §4 · p. 7

5.1. Preliminary results.

• For a uniformly random subset σ of cardinality s drawn from {1, 2, . . . , n}, the bound ∥Hσ×σ∥∞→1 ≤ s/8 holds with probability at least 1/9.Source: §5.1 · p. 8
• Grothendieck factorization writes G = D1 T D2 with nonnegative diagonal factors and relates ∥G∥∞→1 to ∥T∥.For Hermitian G, the factors can be taken equal.Source: §5.1 · p. 8

5.2. Proof of Bourgain–Tzafriri.

• The proof bounds ∥Hτ×τ∥ by 0.5, a condition that Proposition 5.1 connects to a condition-number bound for Aτ.Here H is the hollow Gram matrix of the standardized matrix A.Source: §5.1 · p. 8, §5.2 · p. 9

5.3. Grothendieck factorization via convex optimization.

• For Hermitian G, convex program (5.2) characterizes whether a factorization G = D T D exists with ∥T∥ ≤ α, and an O(1)-value feasible point suffices for the construction.The proof suggests α = s/4; the optimization need not be solved to optimality.Source: §5.3 · pp. 9–10
• The optimization approach also yields approximation algorithms for the (∞, 1) norm, whose exact computation is NP-hard in general.For symmetric matrices, the least feasible α approximates the norm within Grothendieck’s factor KG.Source: §5.3 · p. 10

5.4. An algorithm for Bourgain–Tzafriri.

• For an m × n standardized matrix A, Algorithm 2 produces a set τ = τ⋆ of column indices with probability at least 3/4, at computational cost eO(|τ|^2 m + |τ|^3.5).Source: §5.4 · p. 10

6. Future Directions

• Vershynin’s result guarantees a subset of size at least (1−ε)·st.rank(A) with κ(Aτ) ≤ C(ε), while Szarek gives matrices where selecting st.rank(A) columns forces κ(Aτ) = ω(1).Together, the results indicate that stable rank describes how large a well-conditioned subset can generally be.Source: §6 · p. 10
• The authors suggest simplifying Vershynin’s iterative construction into a column-selection algorithm that might support rank-revealing factorizations.They describe the potential impact on numerical linear algebra as significant.Source: §6 · p. 10

Appendix A. Random Reduction of Norms

• Appendix A studies how random selection reduces the (∞, 2) and (∞, 1) norms, including for structured matrices used in the two theorem proofs..Source: Appendix A · p. 12

A.1. Random Coordinate Models.

• The appendix contrasts fixed-size coordinate sampling, which has dependent selections, with independent Bernoulli sampling, then relates the models for monotonic matrix norms.The Bernoulli model has average cardinality δn; Proposition A.1 is identified as the Poissonization result.Source: §A.1 · p. 12

A.2. Reduction of the (∞, 2) norm.

• The subsection presents a general result for the (∞, 2) norm of uniformly sampled columns, with constants described as especially good.The basic argument is attributed to Bourgain and Tzafriri [BT91, Thm. 1.1], with modern proofs noted as simpler.Source: §A.2 · p. 12
• Corollary A.3 considers a standardized matrix and sets s ≤ ⌈2 st.rank(A)⌉ with δ = s/n as a key step in the Kashin–Tzafriri proof.The supplied passage gives the setup but not the corollary’s resulting bound.Source: §A.2 · p. 12
• The corollary’s proof uses standardization to relate stable rank to n/∥A∥² and applies Poissonization.The proof then bounds terms using the unit norms of A’s columns.Source: §A.2 · p. 13

A.3. Proof of Theorem A.2.

• The proof works with the (2, 1) norm, which is dual to the (∞, 2) norm, to use specialized methods.Source: §A.3 · p. 13
• It centers the selector terms, introduces an independent copy, and uses Jensen’s inequality and random Rademacher signs to symmetrize the expectation.The signs are independent and uniformly distributed over ±1.Source: §A.3 · p. 13
• A Rademacher comparison removes absolute values from the inner product, after which the remaining expectation bounds complete the estimate.Source: §A.3 · pp. 13–14

A.4. Reduction of the (∞, 1) norm.

• Rudelson and Vershynin’s Theorem A.4 gives detailed information about the (∞, 1) norm of a random principal submatrix, in a setting connected to graph cuts.The theorem is cited as being established using Theorem A.2 methods plus a decoupling argument [BT91, Prop. 1.9].Source: §A.4 · p. 14
• Corollary A.5 concerns an n-column standardized matrix with hollow Gram matrix H and a sample size s ≤ ⌈c · st. rank(A)⌉, with δ = s/n.The authors rely on this corollary in their proof of the Bourgain–Tzafriri theorem.Source: §A.4 · p. 14

Appendix B. Entropic Mirror Descent

• The appendix introduces Entropic Mirror Descent as a practical method for the convex minimizations over the probability simplex required by the Kashin–Tzafriri and Bourgain–Tzafriri algorithms.It proceeds to explain the method and its application to the objective functions in the paper.Source: Appendix B · p. 15

B.1. Convex analysis.

• The subdifferential contains vectors satisfying the stated inequalities, and its elements, called subgradients, describe directions and rates of ascent at a point.For a differentiable function, the gradient is the unique subgradient.Source: §B.1 · p. 15
• When J is differentiable at f, its gradient is the unique subgradient.This is the differentiable case of the section’s description of subgradients.Source: §B.1 · p. 15
• The section defines a function’s Lipschitz constant relative to a norm and invokes a dual-norm relation for that characterization [Roc70, Thm. 24.7].The dual norm is denoted |||·|||∗.Source: §B.1 · p. 15

B.2. Interior subgradient methods.

• Projected subgradient methods can increase the objective, leave the constraint set, and zigzag after projection back onto the set.The negative subgradient need not be a descent direction, so these methods are typically nonmonotone.Source: §B.2 · p. 15
• Interior subgradient methods use a constraint-tailored divergence to address some problematic behavior of projected subgradient methods.Each iteration computes a subgradient to linearize the objective, then minimizes that linearization penalized by divergence from the current iterate.Source: §B.2 · p. 15
• Entropic Mirror Descent reweights the current iterate by exp(−βθ) and normalizes the result to obtain the next iterate.The algorithm is presented as returning an approximate minimizer after a specified number of iterations.Source: §B.2 · p. 16
• The divergence penalty keeps successive iterates close, keeps them from approaching the constraint boundary, and with a careful β choice guarantees average progress toward the optimum set.Closeness matters because the linearization is useful only locally.Source: §B.2 · p. 16

B.3. Optimization on the probability simplex.

• Entropic Mirror Descent is designed for convex minimization over the probability simplex, using relative entropy as its divergence and a closed-form update.The update reweights simplex coordinates by exp(−βθ_j) and normalizes them.Source: §B.3 · p. 16
• The EMD efficiency estimate bounds the approximate-minimizer error using the objective’s ℓ1 Lipschitz constant and iteration count; an iteration-dependent step size adds a logarithmic factor.The cited estimate is stated for convex J on the simplex.Source: §B.3 · p. 16

B.4. Pietsch factorization via EMD.

• The Pietsch factorization problem is recast as simplex optimization, where a normalized maximal eigenvector of B∗B − α2 diag(f) supplies a subgradient through its componentwise squared magnitudes.The objective and subgradient can be evaluated from an eigenvector calculation; for sparse B, the matrix need not be formed explicitly.Source: §B.4 · p. 17
• For Pietsch factorization, L ≤ α2 and EMD ostensibly needs eO(α4) iterations for constant precision, though far fewer suffice in practice.The EMD application closely resembles multiplicative weights for maxcut, differing substantially in step sizes and final-solution construction [Kal07, Ch. 6].Source: §B.4 · p. 17

B.5. Grothendieck factorization via EMD. 

• For Grothendieck factorization, EMD uses a maximal eigenvector to form a subgradient; its Lipschitz constant satisfies L ≤ α, giving an apparent eO(α2) iteration count.The eigenvector calculations can be streamlined by exploiting the matrix structure.Source: §B.5 · p. 18