Low-rank matrix completion by Riemannian optimization---extended version

Bart Vandereycken

Source-linked AI summary

Low-rank matrix completion by Riemannian optimization---extended version

Bart Vandereycken

arXiv:1209.3834v1Sep 18, 2012math.OCmath.DGmath.NA

TL;DR

Matrix completion seeks a low-rank matrix from sparse, possibly noisy entries, but standard convex methods scale poorly for very large datasets. This paper minimizes least-square fitting error directly on the fixed-rank matrix manifold using retraction-based nonlinear conjugate gradients, achieving competitive numerical performance with fast asymptotic convergence.

Problem
Matrix completion seeks the lowest-rank matrix agreeing with observed entries, with noisy settings requiring a relaxed misfit formulation; large-scale instances make off-the-shelf convex optimization poorly scalable.
Method
The paper directly minimizes least-square fitting error over the fixed-rank manifold Mk using retraction-based optimization and an efficient nonlinear conjugate-gradient implementation.
Results
The numerical experiments show the approach is very competitive with state-of-the-art matrix-completion solvers and achieves very fast asymptotic convergence factors without parameter tuning.
Takeaways & Limitations
Riemannian optimization can efficiently handle large-scale low-rank matrix completion without requiring vector-space structure for Mk.
Takeaways & Limitations
The proposed approach fixes the manifold rank, although an integrated manifold-based solution that selects rank is desirable.

Abstract

from arXiv · show

The matrix completion problem consists of finding or approximating a low-rank matrix based on a few samples of this matrix. We propose a new algorithm for matrix completion that minimizes the least-square distance on the sampling set over the Riemannian manifold of fixed-rank matrices. The algorithm is an adaptation of classical non-linear conjugate gradients, developed within the framework of retraction-based optimization on manifolds. We describe all the necessary objects from differential geometry necessary to perform optimization over this low-rank matrix manifold, seen as a submanifold embedded in the space of matrices. In particular, we describe how metric projection can be used as retraction and how vector transport lets us obtain the conjugate search directions. Finally, we prove convergence of a regularized version of our algorithm under the assumption that the restricted isometry property holds for incoherent matrices throughout the iterations. The numerical experiments indicate that our approach scales very well for large-scale problems and compares favorably with the state-of-the-art, while outperforming most existing solvers.

1 Introduction

The paper formulates noisy low-rank matrix completion as least-squares optimization over the fixed-rank manifold and solves it with a Riemannian nonlinear conjugate-gradient method. The approach targets large-scale problems where convex solvers scale poorly, while its broader comparison across geometries remains open.

Matrix completion seeks a lowest-rank matrix agreeing with observed entries, with a relaxed formulation allowing sampling-set misfit under noise.
Large-scale applications motivate specialized low-rank algorithms because off-the-shelf convex optimization scales poorly with matrix dimension.
1.1 The proposed method: optimization on manifolds: The proposed method directly optimizes the robust completion objective over matrices of prescribed rank k.
1.1 The proposed method: optimization on manifolds: LRGeomCG adapts nonlinear conjugate gradients to the fixed-rank manifold using Riemannian gradients, vector transport, retraction, and Armijo line search.
1.2 Relation to existing manifold-related methods: The paper reports that its manifold-based approach is competitive with state-of-the-art solvers, while comparisons among manifold geometries remain an open question.

2 Differential geometry for low-rank matrix manifolds

The paper equips the smooth manifold of rank-k matrices with an embedded Euclidean metric and the geometric operations needed for optimization. Metric projection supplies a local retraction, while tangent-space transport enables conjugate-gradient updates.

The fixed-rank set M_k is a smooth submanifold of R^m×n with dimension (m+n−k)k and explicitly characterized tangent spaces.
Using the Frobenius inner product on tangent matrices turns M_k into a Riemannian manifold.
The Riemannian gradient is obtained by orthogonally projecting the Euclidean gradient onto the current tangent space.
2.2 Metric projection as retraction: Metric projection of X+ξ onto M_k is chosen as the retraction and is computed from the SVD, but only locally because global projections may be nonunique or undefined.
2.4 Vector transport: Vector transport moves tangent vectors between nearby tangent spaces, allowing the PR+ conjugate direction to combine the current gradient with the transported previous direction.
The convergence analysis additionally enforces that the search direction remains sufficiently related to the gradient.

3 Implementation details

The implementation exploits compact low-rank factorizations, sparse sampling operations, and small dense matrix decompositions. These choices make retraction, gradient, transport, and conjugate-direction computations practical at large scale, although sparse operations dominate measured runtime.

Low-rank iterates and tangent vectors are stored using compact SVD-based and factored representations, reducing operations to the low-rank regime.
The Riemannian gradient is computed by projecting the sampling residual onto the tangent space, with cost 2(n + 2m)k^2 + 4|Ω|k flops.
Non-linear CG: PR+ conjugate directions transport the previous gradient and direction, update with Polak–Ribiére, and restart when the direction is nearly orthogonal to the gradient.
Initial guess for line search: An exact tangent-space least-squares step provides an effective line-search initial guess, so backtracking was almost never necessary in experiments.
Retraction: Retraction uses compact QR factorizations and an SVD of a small 2k-by-2k matrix, reducing the cost to 14(m+n)k^2 + C_SVDk^3 flops.
Computational cost: All experiments observed nbArmijo = 0, while sparse operations were slower than theoretical estimates because of poor data locality and were penalized by about C_sparse ≃ 5.
Computational cost: Compared with LMAFit, LRGeomCG is more costly per iteration but may compete through faster convergence; the implementation relies mainly on dense linear algebra and sparse matrix-vector products.

4 Convergence for a modified cost function

The convergence analysis modifies the cost function with regularization so iterates remain in a compact subset of the fixed-rank manifold. Under this modification, accumulation points are critical and satisfy tangent-space reconstruction conditions, while whole-space recovery requires additional assumptions.

Every accumulation point of the algorithm is a critical point, and the gradient norm converges to zero.This follows from the smooth retraction, Armijo-type line search, gradient-related directions, and compactness argument.
Reconstruction on the tangent space: Every accumulation point satisfies tangent-space agreement between the sampled iterate and the sampled target.The condition is expressed as PTX∗Mk PΩ(X∗) = PTX∗Mk PΩ(A).
Regularization is introduced because compactness of the iterates cannot generally be established for the original objective without modifying the algorithm or adding assumptions.
The modified algorithm keeps its iterates in a closed and bounded subset of the fixed-rank manifold, which is compact.The proof uses a level set of the regularized objective and bounds the singular values.
The regularization can be negligible when σ1(Xi) = O(1) and σk(Xi) ≫ εmach, but rank drops require transferring the problem to a lower-rank manifold in practice.The paper reports that rank drops did not occur in its experiments, so convergence of that modified procedure is outside its scope.
Reconstruction on the whole space: Tangent-space agreement does not by itself imply whole-space recovery because PΩ usually has a large null space; low-rank completion assumptions are needed.The discussion identifies restricted-isometry-type conditions as the relevant way to quantify whether observations suffice for completion.

5 Numerical experiments

The experiments compare LRGeomCG with LMAFit across size, rank, oversampling, noise, and initialization settings. LRGeomCG generally needs fewer iterations and is faster at high accuracy, while LMAFit can be faster at low precision; hybrid warm starts reduce LRGeomCG’s transient phase.

Influence of size and rank for fixed oversampling: LRGeomCG needs fewer iterations than LMAFit, but its first 20 iterations exhibit slower transient convergence.The transient behavior reflects the near-absence of nonlinear CG acceleration during the initial phase.
Influence of size and rank for fixed oversampling: At high accuracy LRGeomCG is faster, whereas LMAFit is faster when only low precision is required.LRGeomCG’s iterations cost about twice as much, so its asymptotic convergence advantage matters most at tighter tolerances.
Influence of size and rank for fixed oversampling: For difficult problems with OS = 3, LRGeomCG is about 50% faster than LMAFit.The relative advantage depends critically on the convergence factor, which is determined by oversampling.
Hybrid strategy: A hybrid using LMAFit for the first I iterations and then LRGeomCG reduces the 10^-12 solve time from 97 sec. or 147 sec. to about 80 sec. for I = 10, 20, 30, 40.For I = 20, the hybrid has almost no transient behavior and is always faster than either standalone solver; β accelerates almost immediately when I > 10.
Influence of oversampling: As OS approaches 2, both methods slow; LRGeomCG improves with larger OS, while LMAFit’s convergence factor stagnates.In timing, LRGeomCG is significantly faster for smaller and larger OS, whereas LMAFit is only slightly slower around OS = 7...9.
Influence of noise: With noise, both methods solve to the noise level, and iterations remain undisturbed until the final iterations except for the very noisy case ε = 1.Stagnation detection is effective for both methods, and LRGeomCG remains faster when higher accuracy is required.
Exponentially decaying singular values: Random initialization across increasing ranks gives unsatisfactory accuracy, whereas the homotopy strategy performs more favorably and lets LRGeomCG attain a smaller error.The experiment attributes the poor behavior to local optimizers being far from global ones.

6 Conclusions

The paper presents direct least-squares optimization over the fixed-rank manifold for large-scale matrix completion. Experiments show competitive performance, while the fixed-rank requirement remains a drawback.

The method directly minimizes matrix-completion least-squares error over Mk, the set of matrices with rank k.
Retraction-based Riemannian optimization enables nonlinear conjugate gradients despite Mk lacking vector-space structure.The embedded-manifold geometry supports an efficient implementation.
Numerical experiments show the approach is very competitive with state-of-the-art matrix-completion solvers.The method also achieves fast asymptotic convergence factors without parameter tuning.
The approach requires fixing the manifold rank, although an integrated method that determines rank within the manifold framework would be desirable.

A Proof of Proposition 2.3

The proof constructs a second-order retraction and uses it to derive a second-order model of the objective. This model yields the Riemannian Hessian.

A second-order retraction is introduced so Euclidean derivatives can be used to derive the Riemannian Hessian.
The proof therefore proceeds by constructing the retraction first and then deriving the objective’s second-order model.
A second-order Taylor expansion of the specific objective produces the model used to obtain the Riemannian Hessian.

A.1 An explicit second-order retraction on Mk

The appendix constructs an explicit second-order retraction on the fixed-rank manifold and establishes its smoothness, retraction properties, and second-order approximation behavior.

Second-order retractions approximate the exponential map and can be characterized by additional second-order conditions.
The construction maps tangent vectors through explicit matrix operations built from a rank-k factorization X = UΣV^T.
The resulting mapping is shown to be smooth in tangent-vector coordinates and locally smooth on the tangent bundle.
The construction is independent of the particular factorization matrices U, V, and Σ used to represent X.
The image has rank at most k, and continuity of rank ensures rank k sufficiently near the base point.
The appendix also proves that the mapping satisfies the second-order retraction requirements and notes that metric projection admits a corresponding second-order Taylor series.

A.2 Second-order model of f on Mk

The paper derives a second-order model for the least-squares objective on Mk and decomposes its Hessian into linear contributions, including a symmetric operator component.

The Taylor expansion of the explicit second-order retraction yields a second-order model for any smooth objective on Mk.
For matrix completion, the Riemannian gradient is the tangent-space projection of PΩ(X − A).
The Hessian is represented as the sum of two operators corresponding to separate contributions in the inner-product expression.
The first operator is linear, and imposing symmetry determines its coefficient as c = 1/2.
The assembled Hessian expression completes the proof of Proposition 2.3.