A Singular Value Thresholding Algorithm for Matrix Completion

1.1 Motivation

There is a rapidly growing interest in the recovery of an unknown low-rank or approximately lowrank matrix from very limited information. This problem occurs in many areas of engineering and applied science such as machine learning [1, 3, 4], control [42] and computer vision, see [48]. As a motivating example, consider the problem of recovering a data matrix from a sampling of its entries. This routinely comes up whenever one collects partially filled out surveys, and one would like to infer the many missing entries. In the area of recommender systems, users submit ratings on a subset of entries in a database, and the vendor provides recommendations based on the user’s preferences. Because users only rate a few items, one would like to infer their preference for unrated items; this is the famous Netflix problem [2]. Recovering a rectangular matrix from a sampling of its entries is known as the matrix completion problem. The issue is of course that this problem is extraordinarily ill posed since with fewer samples than entries, we have infinitely many completions. Therefore, it is apparently impossible to identify which of these candidate solutions is indeed the “correct” one without some additional information. In many instances, however, the matrix we wish to recover has low rank or approximately low rank. For instance, the Netflix data matrix of all user-ratings may be approximately low-rank because it is commonly believed that only a few factors contribute to anyone’s taste or preference. In computer vision, inferring scene geometry and camera motion from a sequence of images is a wellstudied problem known as the structure-from-motion problem. This is an ill-conditioned problem for objects may be distant with respect to their size, or especially for “missing data” which occur because of occlusion or tracking failures. However, when properly stacked and indexed, these images form a matrix which has very low rank (e.g. rank 3 under orthography) [21, 48]. Other examples of low-rank matrix fitting abound; e.g. in control (system identification), machine learning (multiclass learning) and so on. Having said this, the premise that the unknown has (approximately) low rank radically changes the problem, making the search for solutions feasible since the lowest-rank solution now tends to be the right one. In a recent paper [13], Candès and Recht showed that matrix completion is not as ill-posed as people thought. Indeed, they proved that most low-rank matrices can be recovered exactly from most sets of sampled entries even though these sets have surprisingly small cardinality, and more importantly, they proved that this can be done by solving a simple convex optimization problem. To state their results, suppose to simplify that the unknown matrix M ∈ R^n×n is square, and that

one has available m sampled entries {M_ij: (i, j) ∈ Ω} where Ω is a random subset of cardinality m. Then [13] proves that most matrices M of rank r can be perfectly recovered by solving the optimization problem

provided that the number of samples obeys

for some positive numerical constant C.¹ In (1.1), the functional X_∗ is the nuclear norm of the ‖ ‖ matrix M, which is the sum of its singular values. The optimization problem (1.1) is convex and can be recast as a semidefinite program [29,30]. In some sense, this is the tightest convex relaxation of the NP-hard rank minimization problem

since the nuclear ball X: X_∗ 1 is the convex hull of the set of rank-one matrices with {‖ ‖ ≤} spectral norm bounded by one. Another interpretation of Candès and Recht’s result is that under suitable conditions, the rank minimization program (1.3) and the convex program (1.1) are formally equivalent in the sense that they have exactly the same unique solution.

1.2 Algorithm outline

Because minimizing the nuclear norm both provably recovers the lowest-rank matrix subject to constraints (see [45] for related results) and gives generally good empirical results in a variety of situations, it is understandably of great interest to develop numerical methods for solving (1.1). In [13], this optimization problem was solved using one of the most advanced semidefinite programming solvers, namely, SDPT3 [47]. This solver and others like SeDuMi are based on interior-point methods, and are problematic when the size of the matrix is large because they need to solve huge systems of linear equations to compute the Newton direction. In fact, SDPT3 can only handle n n matrices with n 100. Presumably, one could resort to iterative solvers such as the method × ≤ of conjugate gradients to solve for the Newton step but this is problematic as well since it is well known that the condition number of the Newton system increases rapidly as one gets closer to the solution. In addition, none of these general purpose solvers use the fact that the solution may have low rank. We refer the reader to [40] for some recent progress on interior-point methods concerning some special nuclear norm-minimization problems. This paper develops the singular value thresholding algorithm for approximately solving the nuclear norm minimization problem (1.1) and by extension, problems of the form

where is a linear operator acting on the space of n₁ n₂ matrices and b R^m. This algorithm is A × ∈ a simple first-order method, and is especially well suited for problems of very large sizes in which the solution has low rank. We sketch this algorithm in the special matrix completion setting and let_Ω be the orthogonal projector onto the span of matrices vanishing outside of Ω so that the

(i, j)th component of_Ω(X) is equal to X_ij if (i, j) ∈ Ω and zero otherwise. Our problem may be expressed as

with optimization variable X ∈ R^n1×n2. Fix τ > 0 and a sequence {δ_k}_k≥1 of scalar step sizes. Then starting with Y ⁰ = 0 ∈ R^n1×n2, the algorithm inductively defines

(1.6)

until a stopping criterion is reached. In (1.6), shrink(Y, τ) is a nonlinear function which applies a soft-thresholding rule at level τ to the singular values of the input matrix, see Section 2 for details. The key property here is that for large values of τ, the sequence X^k converges to a solution which {} very nearly minimizes (1.5). Hence, at each step, one only needs to compute at most one singular value decomposition and perform a few elementary matrix additions. Two important remarks are in order:

1. Sparsity. For each k 0, Y ^k vanishes outside of Ω and is, therefore, sparse, a fact which can be used to evaluate the shrink function rapidly. 2. Low-rank property. The matrices X^k turn out to have low rank, and hence the algorithm has minimum storage requirement since we only need to keep principal factors in memory.

Our numerical experiments demonstrate that the proposed algorithm can solve problems, in Matlab, involving matrices of size 30, 000×30, 000 having close to a billion unknowns in 17 minutes on a standard desktop computer with a 1.86 GHz CPU (dual core with Matlab’s multithreading option enabled) and 3 GB of memory. As a consequence, the singular value thresholding algorithm may become a rather powerful computational tool for large scale matrix completion.

1.3 General formulation

The singular value thresholding algorithm can be adapted to deal with other types of convex constraints. For instance, it may address problems of the form

where each f_i is a Lipschitz convex function (note that one can handle linear equality constraints by considering pairs of affine functionals). In the simpler case where the f_i’s are affine functionals, the general algorithm goes through a sequence of iterations which greatly resemble (1.6). This is useful because this enables the development of numerical algorithms which are effective for recovering matrices from a small subset of sampled entries possibly contaminated with noise.

1.4 Contents and notations

The rest of the paper is organized as follows. In Section 2, we derive the singular value threshold-ing (SVT) algorithm for the matrix completion problem, and recasts it in terms of a well-known Lagrange multiplier algorithm. In Section 3, we extend the SVT algorithm and formulate a general iteration which is applicable to general convex constraints. In Section 4, we establish the convergence results for the iterations given in Sections 2 and 3. We demonstrate the performance and effectiveness of the algorithm through numerical examples in Section 5, and review additional implementation details. Finally, we conclude the paper with a short discussion in Section 6. Before continuing, we provide here a brief summary of the notations used throughout the paper. Matrices are bold capital, vectors are bold lowercase and scalars or entries are not bold. For instance, X is a matrix and X_ij its (i, j)th entry. Likewise, x is a vector and x_i its ith component. The nuclear norm of a matrix is denoted by X_∗, the Frobenius norm by X_F ‖ ‖ ‖ ‖ and the spectral norm by X₂; note that these are respectively the 1-norm, the 2-norm and the ‖ ‖ sup-norm of the vector of singular values. The adjoint of a matrix X is X^∗ and similarly for vectors. The notation diag(x), where x is a vector, stands for the diagonal matrix with x_i as {} diagonal elements. We denote by X, Y = trace(X^∗Y) the standard inner product between two

matrices (X ²_F = X, X). The Cauchy-Schwarz inequality gives X, Y X_F Y_F and it is

well known that we also have X, Y X_∗ Y₂ (the spectral and nuclear norms are dual from

one another), see e.g. [13,45].

2.1 The singular value shrinkage operator

Consider the singular value decomposition (SVD) of a matrix X ∈ R^n1×n2 of rank r

where U and V are respectively n₁ × r and n₂ × r matrices with orthonormal columns, and the

singular values σ_i are positive (unless specified otherwise, we will always assume that the SVD of a matrix is given in the reduced form above). For each τ 0, we introduce the soft-thresholding operator_τ defined as follows:

where t₊ is the positive part of t, namely, t₊ = max(0, t). In words, this operator simply applies a soft-thresholding rule to the singular values of X, effectively shrinking these towards zero. This is the reason why we will also refer to this transformation as the singular value shrinkage operator. Even though the SVD may not be unique, it is easy to see that the singular value shrinkage operator is well defined and we do not elaborate further on this issue. In some sense, this shrinkage operator is a straightforward extension of the soft-thresholding rule for scalars and vectors. In particular, note that if many of the singular values of X are below the threshold τ, the rank of_τ (X) may be considerably lower than that of X, just like the soft-thresholding rule applied to vectors leads to sparser outputs whenever some entries of the input are below threshold. The singular value thresholding operator is the proximity operator associated with the nuclear norm. Details about the proximity operator can be found in e.g. [35]. Theorem 2.1 For each τ 0 and Y ∈ R^n1×n2, the singular value shrinkage operator (2.2) obeys

Proof. Since the function h₀(X):= τ X_∗ + ¹₂ X Y ²_F is strictly convex, it is easy to see that

there exists a unique minimizer, and we thus need to prove that it is equal to_τ (Y). To do this, recall the definition of a subgradient of a convex function f: R^n1×n2 R. We say that Z is a subgradient of f at X₀, denoted Z ∈ ∂f(X₀), if

for all X. Now X̂ minimizes h₀ if and only if 0 is a subgradient of the functional h₀ at the point

where ∂ X̂_∗ is the set of subgradients of the nuclear norm. Let X ∈ R^n1×n2 be an arbitrary

matrix and UΣV ^∗ be its SVD. It is known [13,37,49] that

Set X̂:=_τ (Y) for short. In order to show that X̂ obeys (2.5), decompose the SVD of Y as

where U₀, V₀ (resp. U₁, V₁) are the singular vectors associated with singular values greater than τ (resp. smaller than or equal to τ). With these notations, we have

By definition, U₀^∗W = 0, W V₀ = 0 and since the diagonal elements of Σ₁ have magnitudes bounded by τ, we also have W₂ 1. Hence Y X̂ τ∂ X̂_∗, which concludes the proof. ‖ ‖ ≤ − ∈ ‖ ‖

2.2 Shrinkage iterations

We are now in the position to introduce the singular value thresholding algorithm. Fix τ > 0 and a sequence δ_k of positive step sizes. Starting with Y₀, inductively define for k = 1, 2,..., {}

(2.7)

until a stopping criterion is reached (we postpone the discussion this stopping criterion and of the choice of step sizes). This shrinkage iteration is very simple to implement. At each step, we only need to compute an SVD and perform elementary matrix operations. With the help of a standard numerical linear algebra package, the whole algorithm can be coded in just a few lines. Before addressing further computational issues, we would like to make explicit the relationship between this iteration and the original problem (1.1). In Section 4, we will show that the sequence {^k} X converges to the unique solution of an optimization problem closely related to (1.1), namely,

Furthermore, it is intuitive that the solution to this modified problem converges to that of (1.5) as τ as shown in Section 3. Thus by selecting a large value of the parameter τ, the sequence of → ∞ iterates converges to a matrix which nearly minimizes (1.1). As mentioned earlier, there are two crucial properties which make this algorithm ideally suited for matrix completion. • {^k} Low-rank property. A remarkable empirical fact is that the matrices in the sequence X have low rank (provided, of course, that the solution to (2.8) has low rank). We use the word “empirical” because all of our numerical experiments have produced low-rank sequences but we cannot rigorously prove that this is true in general. The reason for this phenomenon is, however, simple: because we are interested in large values of τ (as to better approximate the solution to (1.1)), the thresholding step happens to ‘kill’ most of the small singular values and produces a low-rank output. In fact, our numerical results show that the rank of X^k is nondecreasing with k, and the maximum rank is reached in the last steps of the algorithm, see Section 5.

Thus, when the rank of the solution is substantially smaller than either dimension of the matrix, the storage requirement is low since we could store each X^k in its SVD form (note that we only need to keep the current iterate and may discard earlier values).

Sparsity. Another important property of the SVT algorithm is that the iteration matrix Y ^k is sparse. Since Y ⁰ = 0, we have by induction that Y ^k vanishes outside of Ω. The fewer entries available, the sparser Y ^k. Because the sparsity pattern Ω is fixed throughout, one can then apply sparse matrix techniques to save storage. Also, if Ω = m, the computational cost | | of updating Y ^k is of order m. Moreover, we can call subroutines supporting sparse matrix computations, which can further reduce computational costs.

One such subroutine is the SVD. However, note that we do not need to compute the entire SVD of Y ^k to apply the singular value thresholding operator. Only the part corresponding to singular values greater than τ is needed. Hence, a good strategy is to apply the iterative Lanczos algorithm to compute the first few singular values and singular vectors. Because Y ^k is sparse, Y ^k can be applied to arbitrary vectors rapidly, and this procedure offers a considerable speedup over naive methods.

2.3 Relation with other works

Our algorithm is inspired by recent work in the area of ℓ₁ minimization, and especially by the work on linearized Bregman iterations for compressed sensing, see [9–11,23,44,51] for linearized Bregman iterations and [14–17,26] for some information about the field of compressed sensing. In this line of work, linearized Bregman iterations are used to find the solution to an underdetermined system of linear equations with minimum ℓ₁ norm. In fact, Theorem 2.1 asserts that the singular value thresholding algorithm can be formulated as a linearized Bregman iteration. Bregman iterations were first introduced in [43] as a convenient tool for solving computational problems in the imaging sciences, and a later paper [51] showed that they were useful for solving ℓ₁-norm minimization problems in the area of compressed sensing. Linearized Bregman iterations were proposed in [23] to improve performance of plain Bregman iterations, see also [51]. Additional details together with a technique for improving the speed of convergence called kicking are described in [44]. On the practical side, the paper [11] applied Bregman iterations to solve a deblurring problem while on the theoretical side, the references [9, 10] gave a rigorous analysis of the convergence of such iterations. New developments keep on coming out at a rapid pace and recently, [32] introduced a new iteration, the split Bregman iteration, to extend Bregman-type iterations (such as linearized Bregman iterations) to problems involving the minimization of ℓ₁-like functionals such as totalvariation norms, Besov norms, and so forth. When applied to ℓ₁-minimization problems, linearized Bregman iterations are sequences of soft-thresholding rules operating on vectors. Iterative soft-thresholding algorithms in connection with ℓ₁ or total-variation minimization have quite a bit of history in signal and image processing and we would like to mention the works [12, 39] for total-variation minimization, [24, 25, 31] for ℓ₁ minimization, and [5, 7, 8, 19, 20, 27, 28, 46] for some recent applications in the area of image inpainting and image restoration. Just as iterative soft-thresholding methods are designed to find sparse solutions, our iterative singular value thresholding scheme is designed to find a sparse vector of singular values. In classical problems arising in the areas of compressed sensing, and signal or image processing, the sparsity is expressed in a known transformed domain and soft-thresholding is applied to transformed coefficients. In contrast, the shrinkage operator_τ is adaptive. The SVT not

only discovers a sparse singular vector but also the bases in which we have a sparse representation. In this sense, the SVT algorithm is an extension of earlier iterative soft-thresholding schemes. Finally, we would like to contrast the SVT iteration (2.7) with the popular iterative soft-thresholding algorithm used in many papers in imaging processing and perhaps best known under the name of Proximal Forward-Backward Splitting method (PFBS), see [8,22,24,31,33] for example. The constrained minimization problem (1.5) may be relaxed into

for some λ > 0. Theorem 2.1 asserts that_λ is the proximity operator of λ X_∗ and Proposition

3.1(iii) in [22] gives that the solution to this unconstrained problem is characterized by the fixed point equation X =_λδ(X + δP_Ω(M X)) for each δ > 0. One can then apply a simplified

version of the PFBS method (see (3.6) in [22]) to obtain iterations of the form

Introducing an intermediate matrix Y^k−k, this algorithm may be expressed as

(2.10)

The difference with (2.7) may seem subtle at first—replacing X^k in (2.10) with Y ^k−1 and setting δ_k = δ gives (2.7) with τ = λδ—but has enormous consequences as this gives entirely different algorithms. First, they have different limits: while (2.7) converges to the solution of the constrained minimization (2.8), (2.10) converges to the solution of (2.9) provided that the sequence of step sizes is appropriately selected. Second, selecting a large λ (or a large value of τ = λδ) in (2.10) gives a low-rank sequence of iterates and a limit with small nuclear norm. The limit, however, does not fit the data and this is why one has to choose a small or moderate value of λ (or of τ = λδ). However, when λ is not sufficiently large, the X^k may not have low rank even though the solution has low rank (and one may need to compute many singular vectors), and Y ^k is not sufficiently sparse to make the algorithm computationally attractive. Moreover, the limit does not necessary have a small nuclear norm. These are reasons why (2.10) is not suitable for matrix completion.

2.4 Interpretation as a Lagrange multiplier method

In this section, we recast the SVT algorithm as a type of Lagrange multiplier algorithm known as Uzawa’s algorithm. An important consequence is that this will allow us to extend the SVT algorithm to other problems involving the minimization of the nuclear norm under convex constraints, see Section 3. Further, another contribution of this paper is that this framework actually recasts linear Bregman iterations as a very special form of Uzawa’s algorithm, hence providing fresh and clear insights about these iterations. In what follows, we set f_τ (X) = τ X_∗ +₂¹ X ²_F for some fixed τ > 0 and recall that we wish ‖ ‖ ‖ ‖

The Lagrangian for this problem is given by

where Y ∈ R^n1×n2. Strong duality holds and X^⋆ and Y ^⋆ are primal-dual optimal if (X^⋆, Y ^⋆) is a saddlepoint of the Lagrangian (X, Y), i.e. a pair obeying sup i_Xnf (X, Y) = (X^⋆, Y ^⋆) = i_Xnf sup (X, Y). (2.11) Y Y L L L (The function g₀(Y) = inf_X (X, Y) is called the dual function.) Uzawa’s algorithm approaches the problem of finding a saddlepoint with an iterative procedure. From Y₀ = 0, say, inductively define

(2.12) where δ_k≥1 is a sequence of positive step sizes. Uzawa’s algorithm is, in fact, a subgradient {} method applied to the dual problem, where each step moves the current iterate in the direction of the gradient or of a subgradient. Indeed, observe that

where X̃ is the minimizer of the Lagrangian for that value of Y so that a gradient descent update for Y is of the form

It remains to compute the minimizer of the Lagrangian (2.12), and note that

However, we know that the minimizer is given by_τ (_Ω(Y)) and since Y ^k =∑_Ω(Y ^k) for all k 0, D P P ≥ Uzawa’s algorithm takes the form

which is exactly the update (2.7). This point of view brings to bear many different mathematical tools for proving the convergence of the singular value thresholding iterations. For an early use of Uzawa’s algorithm minimizing an ℓ₁-like functional, the total-variation norm, under linear inequality constraints, see [12].

3.1 Linear equality constraints

Set the objective functional f_τ (X) = τ X_∗ +₂¹ X ²_F for some fixed τ > 0, and consider the ‖ ‖ ‖ ‖ following optimization problem:

where is a linear transformation mapping n₁ n₂ matrices into R^m (^∗ is the adjoint of). This A × A A more general formulation is considered in [13] and [45] as an extension of the matrix completion problem. Then the Lagrangian for this problem is of the form

where X ∈ R^n1×n2 and y ∈ R^m, and starting with y⁰ = 0, Uzawa’s iteration is given by

(3.3)

The iteration (3.3) is of course the same as (2.7) in the case where is a sampling operator extracting m entries with indices in Ω out of an n₁ × n₂ matrix. To verify this claim, observe that in this situation, =_Ω, and let M be any matrix obeying (M) = b. Then defining A^∗A P A Y ^k = ^∗(y^k) and substituting this expression in (3.3) gives (2.7).

3.2 General convex constraints

One can also adapt the algorithm to handle general convex constraints. Suppose we wish to

minimize f_τ (X) defined as before over a convex set X ∈. To simplify, we will assume that this convex set is given by

where the f_i’s are convex functionals (note that one can handle linear equality constraints by considering pairs of affine functionals). The problem of interest is then of the form

Just as before, it is intuitive that as τ, the solution to this problem converges to a minimizer → ∞ of the nuclear norm under the same constraints (1.7) as shown in Theorem 3.1 at the end of this section.

Put (X):= (f₁(X),..., f_m(X)) for short. Then the Lagrangian for (3.4) is equal to

where X R^n1×n2 and y R^m is now a vector with nonnegative components denoted, as usual, ∈ ∈ by y 0. One can apply Uzawa’s method just as before with the only modification that we will use a subgradient method with projection to maximize the dual function since we need to make sure that the successive updates y^k belong to the nonnegative orthant. This gives

Above, x₊ is of course the vector with entries equal to max(x_i, 0). When is an affine mapping of the form b (X) so that one solves − A

this simplifies to

(3.6)

and thus the extension to linear inequality constraints is straightforward.

3.3 Example

An interesting example concerns the extension of the Dantzig selector [18] to matrix problems.

Suppose we have available linear measurements about a matrix M of interest

where z ∈ R^m is a noise vector. Then under these circumstances, one might want to find the matrix which minimizes the nuclear norm among all matrices which are consistent with the data b.

Inspired by the work on the Dantzig selector which was originally developed for estimating sparse parameter vectors from noisy data, one could approach this problem by solving

where E is an array of tolerances, which is adjusted to fit the noise statistics [18]. Above, vec(A) vec(B), for any two matrices A and B, means componentwise inequalities; that is, A_ij B_ij for all indices i, j. We use this notation as not to confuse the reader with the positive semidefinite ordering.

In the case of the matrix completion problem where extracts sampled entries indexed by Ω, one can always see the data vector as the sampled entries of some matrix B obeying_Ω(B) = ^∗(b); P A the constraint is then natural for it may be expressed as

If z is white noise with standard deviation σ, one may want to use a multiple of σ for E_ij. In words, we are looking for a matrix with minimum nuclear norm under the constraint that all of its sampled entries do not deviate too much from what has been observed. Let Y₊ R^n1×n2 (resp. Y₋ R^n1×n2) be the Lagrange multiplier associated with the compo-∈ ∈ nentwise linear inequality constraints vec(^∗(r)) vec(E) (resp. vec(^∗(r)) vec(E)). Then A ≤ − A ≤ starting with Y_±⁰ = 0, the SVT iteration for this problem is of the form

where again [·]₊ is applied componentwise.

We conclude by noting that in the matrix completion problem where =_Ω and one A^∗A P observes_Ω(B), one can check that this iteration simplifies to

(3.10)

Again, this is easy to implement and whenever the solution has low rank, the iterates X^k have low rank as well.

3.4 When the proximal problem gets close

We now show that minimizing the proximal objective f_τ (X) = τ X_∗ +₂¹ X ²_F is the same as ‖ ‖ ‖ ‖ minimizing the nuclear norm in the limit of large τ’s. The theorem below is general and covers the special case of linear equality constraints as in (2.8).

Theorem 3.1 Let X_τ^⋆ be the solution to (3.4) and X_∞ be the minimum Frobenius-norm solution

Assume that the f_i(X)’s, 1 i m, are convex and lower semi-continuous. Then

Proof. It follows from the definition of X_τ^⋆ and X_∞ that

Summing these two inequalities gives

which implies that X_τ^{⋆ 2}_F is bounded uniformly in τ. Thus, we would prove the theorem if we

could establish that any convergent subsequence {X{_τ^⋆_k}_k≥}₁ must converge to X_∞.

Consider an arbitrary converging subsequence {X^k_τ^⋆_k} and set X_c:= lim_k→∞ X_τ^⋆_k. Since for

each 1 i m, f_i(X_τ^⋆_k) 0 and f_i is lower semi-continuous, X_c obeys

Furthermore, since X_τ^{⋆ 2}_F is bounded, (3.13) yields

An immediate consequence is lim_τ→∞ X_τ^⋆_∗ = X_∞∗ and, therefore, X_c∗ = X_∞∗. This

shows that X_c is a solution to (1.1). Now it follows from the definition of X_∞ that X_cF

X, while we also have X X because of (3.14). We conclude that X =

X and thus X_c = X_∞ since X_∞ is unique.

4.1 Convergence for matrix completion

We begin by recording a lemma which establishes the strong convexity of the objective f_τ.

Lemma 4.1 Let Z ∈ ∂f_τ (X) and Z^′ ∈ ∂f_τ (X^′). Then

Proof. An element Z of ∂f_τ (X) is of the form Z = τZ₀ + X, where Z₀ ∂ X_∗, and similarly ∈ ‖ ‖ for Z^′. This gives 〈 − ^′ − ^′〉 〈₀ −₀^′ − ^′〉 ‖ − ^′‖²_F

and it thus suffices to show that the first term of the right-hand side is nonnegative. From (2.6), we have that any subgradient of the nuclear norm at X obeys Z₀₂ 1 and Z₀, X = X_∗. In

particular, this gives

Whence, 〈₀ −₀^′ − ^′〉 〈₀ 〉 〈₀^{′ ′}〉 − 〈₀ ^′〉 − 〈₀^′ 〉

which proves the lemma. This lemma is key in showing that the SVT algorithm (2.7) converges. Theorem 4.2 Suppose that the sequence of step sizes obeys 0 < inf δ_k sup δ_k < 2. Then the sequence X^k obtained via (2.7) converges to the unique solution of (2.8). {} Proof. Let (X^⋆, Y ^⋆) be primal-dual optimal for the problem (2.8). The optimality conditions give

for some Z^k ∂f_τ (X^k) and some Z^⋆ ∂f_τ (X^⋆). We then deduce that ∈ ∈

and, therefore, it follows from Lemma 4.1 that

We continue and observe that because_ΩX^⋆ =_ΩM,

Therefore, setting r_k =_Ω(Y ^k Y ^⋆)_F,

since for any matrix X,_Ω(X)_F X_F. Under our assumptions about the size of δ_k, we have

2δ_k δ_k² β for all k 1 and some β > 0 and thus

Two properties follow from this: