Truncated Nuclear Norm Minimization for Image Restoration Based On Iterative Support Detection

Recovering a large matrix from limited measurements is a challenging task arising in many real applications, such as image inpainting, compressive sensing and medical imaging, and this kind of problems are mostly formulated as low-rank matrix approximation problems. Due to the rank operator being non-convex and discontinuous, most of the recent theoretical studies use the nuclear norm as a convex relaxation and the low-rank matrix recovery problem is solved through minimization of the nuclear norm regularized problem. However, a major limitation of nuclear norm minimization is that all the singular values are simultaneously minimized and the rank may not be well approximated \cite{hu2012fast}. Correspondingly, in this paper, we propose a new multi-stage algorithm, which makes use of the concept of Truncated Nuclear Norm Regularization (TNNR) proposed in \citep{hu2012fast} and Iterative Support Detection (ISD) proposed in \citep{wang2010sparse} to overcome the above limitation. Besides matrix completion problems considered in \citep{hu2012fast}, the proposed method can be also extended to the general low-rank matrix recovery problems. Extensive experiments well validate the superiority of our new algorithms over other state-of-the-art methods.


Introduction
In many real applications such as machine learning [3][4][5], computer vision [6], and control [7], etc., we seek to recover an unknown (approximately) low-rank matrix from limited information. This problem can be naturally formulated as the following model where X ∈ R m×n is the decision variable, the linear map A: R m×n → R p (p < mn) and vector b ∈ R p are given. However, this is usually NP-hard due to the nonconvexity and discontinuous nature of the rank function. In paper [8], Fazel et.al firstly solved rank minimization problem by approximating the rank function using the nuclear norm (i.e. the sum of singular values of a matrix). Moreover, theoretical studies show that the nuclear norm is the tightest convex lower bound of the rank function of matrices [9]. Thus, an unknown (approximately) low-rank matrixX can be perfectly recovered by solving the optimization problem where X * = min(m,n) i=1 σ i (X) is the nuclear norm and σ i (X) is the i−th largest singular value of X, under some conditions on the linear transformation A.
As a special case, the problem (2) is reduced to the well-known matrix completion problem (3) [10,11], when A is a sampling (or projection/restriction) operator.
where M ∈ R m×n is the incomplete data matrix and Ω is the set of locations corresponding to the observed entries. To solve the kind of problems, we can refer to [10][11][12][13][14][15] for some breakthrough results. Nevertheless, they may obtain suboptimal performance in real applications because the nuclear norm may not be a good approximation to rank-operator, because all the non-zero singular values in rank-operator have the equal contribution, while the singular values in nuclear norm are treated differently by adding them together. Thus, to overcome the weakness of nuclear norm, Truncated Nuclear Norm Regularization (TNNR) was proposed for matrix completion, which only minimizes the smallest min(m, n)−r singular values [1]. The similar truncation idea was also proposed in our previous work [2]. Correspondingly, the problem can be formulated as where X r is defined as the sum of min(m, n) − r minimum singular values. In this way, ones can get a more accurate and robust approximation to the rank-operator on both synthetic and real visual data sets.
In this paper, we aim to extend the idea of TNNR from the special matrix completion problem to the general problem (2) and give the corresponding fast algorithm. More important, we will consider how to fast estimate r, which is usually unavailable in practice.
Throughout this paper, we use the following notation. We let < ·, · > be the standard inner product between two matrices in a finite dimensional Euclidean space, · be the 2-norm, and · F be the Frobenius norm for matrix variables. The projection operator under the Euclidean distance measure is denoted by P and the transpose of a real matrix by ⊺. Let X = U ΣV ⊺ is the singular value decomposition (SVD) for X, where Σ = diag(σ i ), 1 ≤ i ≤ min{m, n}, and σ 1 ≥ · · · ≥ σ min{m,n} .

1.A. Related Work
The low-rank optimization problem (2) has attracted more and more interests in developing customized algorithms, particularly for lager-scale cases. We now briefly review some influential approaches to these problems.
The convex problem (2) can be easily reformulated into the semi-definite programming (SDP) problems [16,17] to make use of the generic SDP solvers such as SDPT3 [18] and SeDuMi [19] which are based on the interior-point method. However, the interior-point approaches suffer from the limitation that they ineffectively handle large-scale problems which was mentioned in [9,20,21]. The problem (2) can also be solved through a projected subgradient approach in [9], whose major computation is concentrated on singular values decomposition. The method can be used to solve large-scale cases of (2). However, the convergence may be slow, especially when high accuracy is required. In [9,22], UV-parametrization (X = U V ⊺ ) based on matrix factorization is applied in general low-rank matrix reconstruction problems. Specifically, the low-rank matrix X is decomposed into the form U V ⊺ , where U ∈ R m×r and V ∈ R n×r are tall and thin matrices. The method reduces the dimensionality from mn to (m + n)r. However, if the rank and the size are large, the computation cost may also be very high. Moreover, the rank r is not known as a priori for most of applications, and it has to be estimated or dynamically adjusted, which might be difficult to realize. More recently, the augmented Lagrangian method (ALM) [23,24] and the alternating direction method of multipliers (ADMM) [25] are very efficient for some convex programming problems arising from various applications. In [26], ADMM is applied to solve (2) with AA * = I.
As an important special case of problem (2), the matrix completion problem (3) has been widely studied. Cai et al. [14,15] used the shrinkage operator to solve the nuclear norm effectively. The shrinkage operator applies a soft-thresholding rule to singular values, as the sparse operator of a matrix, though it can be applied widely in mangy other approaches. However, due to the above-mentioned limitation of nuclear norm, TNNR (4) is proposed to replace the nuclear norm. Since X r in (4) is non-convex, it can not be solved simply and effectively. So, how to change (4) into a convex func-tion is critical. Obviously, it is noted that where U ΣV ⊺ is the SVD of X, U = (u 1 , · · · , u m ) ∈ R m×m , Σ ∈ R m×n and V = (v 1 , · · · , v n ) ∈ R n×n . Then L r = (u 1 , · · · , u r ) T , R r = (v 1 , · · · , v r ) T and the optimization problem (4) can be rewritten as: While the problem (5) is still non-convex, they can get a local minima by an iterative procedure proposed in [1] and we will review the procedure in more details later.
A similar idea of truncation in the context of the sparse vectors, is also implemented on the sparse signals by our previous work in [2], which tries to adaptively learn the information of the nonzeros of the unknown true signal. Specifically, we present a sparse signal reconstruction method, Iterative Support Detection (ISD, for short), aiming to achieve fast reconstruction and a reduced requirement on the number of measurements compared to the classical l 1 minimization approach. ISD alternatively calls its two components: support detection and signal reconstruction. From an incorrect reconstruction, support detection identifies an index set I containing some elements of supp(x) = {i :x i = 0}, and signal reconstruction solves To obtain the reliable support detection from inexact reconstructions, ISD must take advantage of the features and prior information about the true signalx. In [2], the sparse or compressible signals with components having a fast decaying distribution of nonzeros, are considered.

1.B. Contributions and Paper Organization
Our first contribution is the estimation of r, on which TNNR heavily depends (in X r ). Hu et.al [1] seeks the best r by trying all the possible values and this leads high computational cost. In this paper, motivated by Wang et.al [2], we propose singular value estimation (SVE) method to obtain the best r, which can be considered as a special implementation of iterative support detection of [2] in case of matrices.
Our second contribution is to extend TNNR from matrix completion to the general low-rank cases. In [1], they have only considered the matrix completion problem.
The third contribution is based on the above two. In particular, a new efficient algorithmic framework is proposed for the low-rank matrix recovery problem. We name it LRISD, which iteratively calls its two components: SVE and solving the low-rank matrix reconstruction model based on TNNR.
The rest of this paper is organized as follows. In Section 2, computing framework of LRISD and theories of SVE are introduced. In Section 3, Section 4 and Section 5, we introduce three algorithms to solve the problems mentioned in subsection (2.A). Experimental results are presented in Section 6. Finally, some conclusions are made in Section 7.

Iterative Support Detection for Low-rank Problems
In this section, we first give the outline of the proposed algorithm LRISD and then elaborate the proposed SVE method which is a main component of LRISD.

2.A. Algorithm Outline
The main purpose of LRISD is to provide a better approximation to the model (1) than the common convex relaxation model (2). The key idea is to make use of the Truncated Nuclear Norm Regularization (TNNR) defined in (4) and its variant (5) [1]. While ones can passively try all the possible values of r which is the number of largest few singular values, we proposed to actively estimate the value of r.
In addition, we will consider the general low-rank recovery problems beyond the matrix completion problem. Specifically, we will solve three models equality-model (6), unconstrained-model (7) and inequality-model (8): where µ > 0 and δ > 0 are parameters reflecting the level of noise. The models (7) and (8) consider the case with noisy data. Here A is a linear mapping such as partial discrete cosine transformation (DCT), partial discrete walsh hadamard transformation (DWHT), discrete fourier transform (DFT). The general framework of LRISD, as an iterative procedure, starts from the initial r = 0, i.e. solving a plain nuclear norm minimization problem, and then estimates r based on the recovered result. Based on the estimated r, we solve a resulted TNNR model (6), or (7) or (8).
Using the new recovered result, we can update the r value and solve a new TNNR model (6), or (7) or (8). Our algorithm iteratively calls the r estimation and the solver of the TNNR model.
As for solving (6), (7) and (8), we following the idea of [1]. Specifically, A simple but efficient iterative procedure is adopted to decouple the L r , X and R r . We set the initial guess X 1 . In the l-th iteration, we first fix X l and compute L l r and R l r as described in (5), based on the SVD of X l . Then fix L l r and R l r to update X l+1 by solving the following problems, respectively: In [1], the authors have studied to solve the special case of matrix completion problems. For the general problems (9), (10) and (11), we will extend the current state of the art algorithms to solve them in Section 3, Section 4 and Section 5, respectively. In summary, the procedure of LRISD, as a new multistage algorithm, is summarized in the following Algorithm 1. By alternately running the SVE and solving the corresponding TNNR models, the iterative scheme will converge to a solution of a TNNR model, whose solution is expected to be better than that of the plain nuclear norm minimization model (2).
2. Repeat until r reaches a stable value.
Step 1. Estimate r via SVE on X re .
Step 2. Initialization: In the l − th iteration: a) Compute L l r and R l r of (9) ((10) or (11)) according to the current X l in the same way as (5) mentioned. b) Solve the model (9) ((10) or (11)) and obtain X l+1 . Goto a) until Step 3. Set X re = X l+1 .
3. Return the recovered matrix X re .
In the following content, we will explain the implementation of SVE of Step 1 in more details, and extend the existing algorithms for nuclear norm regularized models to TNNR based models (9)-(11) in the procedure b) of Step 2.

2.B. Singular Value Estimate
In this subsection, we mainly focus on Step 1 of LRISD, i.e. describing the process of SVE to estimate r, which is the number of largest few singular values. While it is feasible to find the best r via trying all possible r as done in [1], this procedure is not computationally efficient. Thus, we aim to quickly give an estimate of the best r.
As we have known, for (approximately) low-rank matrices or images, these singular values often have a feature that they all have a fast decaying distribution (as showed in Fig 1). To take advantage of this feature, we can extend our previous work ISD [2] from detecting the large components of sparse vectors to the large singular values of low-rank matrices. In particular, SVE is nothing but a specific implementation of support detection in cases of low-rank matrices, with the aim to acquire the estimation of the true r.  Now we present the process of SVE and the effectiveness of SVE. It is noted that as showed in the Algorithm 1, SVE is repeated several times until a stable estimate r is obtained. For each time, given the reference image X re , we can obtain the singular value vector S of X re by SVD. A natural way to find the positions of the true large singular values based on S, which is considered as an estimate of the singular value vector of the true matrixX, is based on thresholding due to the fast decaying property of the singular values. The choice of ǫ should be case-dependent. In the spirit of ISD, one can use the so called "last significant jump" rule to set the threshold value ǫ to detect the large singular values and minimized the false detections, if we assume that the components of S are sorted from large to small. The straightforward way to apply the "last significant jump" rule is to look for the largest i such that where τ is a prescribed value, and St is defined as absolute values of the first order difference of S. This amounts to sweeping the decreasing sequence {S i } and look for the last jump larger than τ . For example, the selected i = 4, then we set ǫ = S 4 , However, in this paper, unlike the original ISD paper [2], we propose to apply the "last significant jump" rule on absolute values of the first order difference of S, i.e., St, instead of S. Specifically, we look for the largest i such that where κ will be selected below, and Stt is defined as absolute values of the second order difference of S . This amounts to sweeping the decreasing sequence {St i } and look for the last jump larger than κ. For example, the selected i = 4, then we set ǫ = S 4 . We set the estimation rankr to be the cardinality of I, or a close number to it. Specifically, St is computed to obtain jump sizes which count on the change of two neighboring components of S. Then, to reflect the stability of these jumps, the difference of St need to be considered as we just do, because the few largest singular values jump actively, while the small singular values would not change much. The cut-off threshold κ is determined via certain heuristic methods in our experiments: synthetic and real visual data sets. Note that in subsection 6.F, we will present a reliable rule for determining threshold value κ. (9) and (11) In this section, we extend the existing ADMM method in [25] originally for the nuclear norm regularized model to solve (9) and (11) under common linear mapping A (AA * = I), and give closed-form solutions. The extended verison of ADMM is named as TNNR-ADMM, and the original ADMM for the corresponding nuclear norm regularized low-rank matrix recovery model is denoted as LR-ADMM. In addition, we can deduce that the resulting subproblems are simple enough to have closed-form solutions and can be easily achieved to high precision. We start this section with some preliminaries which are convenient for the presentation of algorithms later.

TNNR-ADMM For
When AA * = I, we present the following conclusions [26]: where (I + αA * A) −1 denotes the inverse operator of (I + αA * A) and α > 0. Definition 3.1( [26]): When A satisfies AA * = I, for δ ≥ 0 and Y ∈ R m×n , the projection of Y onto B δ is defined as In particular, when δ = 0, Then, we have the following conclusion: Definition 3.2: For the matrix X ∈ R m×n , X have the singular value decomposition as following: where (s) + = max{0, s}. Theorem 3.2( [15]): For each τ ≥ 0 and Y ∈ R m×n , we have the following conclusion: Definition 3.3: Denote A * be the adjoint operator of A satisfying the following condition:

3.A. Algorithmic Framework
The problems (9) and (11) can be easily reformulated into the following linear constrained convex problem: where δ ≥ 0. In particular, the above formulation is equivalent to (9) when δ = 0. The augmented Lagrangian function of (14) is: where Z ∈ R m×n is the Lagrange multiplier of the linear constraint, β > 0 is the penalty parameter for the violation of the linear constraint. The idea of ADMM is to decompose the minimization task in (15) into three easier and smaller subproblems such that the involved variables X and Y can be minimized separately in altering order. In particular, we apply ADMM to solve (15), and obtain the following iterative scheme: Ignoring constant terms and deriving the optimal conditions for the involved subproblems in (16), we can easily verify that the iterative scheme of the TNNR-ADMM approach for (9) and (11) is as follows.

3.B. The Analysis of Subproblems
According to the analysis above, the computation of each iteration of TNNR-ADMM approach for (9) and (11) is dominated by solving the subproblems (17) and (18). We now elaborate on the strategies for solving these subproblems based on abovementioned preliminaries.
First, the solution of (17) can be obtained explicitly via Theorem 3.2: which is the closed-form solution. Second, it is easy to obtain: Combining (20) and equipped with Definition 3.1, we give the final closed-form solution of the subproblem (18): When δ = 0, it is the particular case of (9) and can be expressed as: Therefor, when the TNNR-ADMM is applied to solve (9) and (11), the generated subproblems all have closedform solutions. Besides, some remarks are in order.
• The convergence of the iterative scheme is well studied in [30]. Here, we omit the convergence analysis. (10) In this section, we consider the model (10), which has attracted a lot of attention in certain multi-task learning problems [21,[31][32][33]. While TNNR-ADMM can be applied to solve this model, it is preferred for the noiseless problems. For the simple version of the model (10), i.e. the one based on the common nuclear norm regularization, many accelerated gradient techniques [34,35] based on [36] are proposed. Among them, an accelerated proximal gradient line search (APGL) method proposed by Beck et al. [35] has been extended to solve TNNR based matrix completion model in [1]. In this paper, we can extend APGL to solve the more general TNNR based the low-rank recovery problem (10).

4.A. TNNR-APGL with Noisy Data
For completeness, we give a short overview of the APGL method. The original model aims to solve the following problem: where G(X), T (X) meet these conditions: • G : R m×n → R is a continuous convex function, possibly nondifferentiable function.
• F : R m×n → R is a convex and differentiable function. In other words, it is continuously differentiable with Lipschitz continuous gradient L(F ) (L(F ) > 0 is the Lipschitz constant of ∇F ).
By linearizing F (X) at Y and adding a proximal term, APGL constructs an approximation of T (X). More specially, we have where τ > 0 is a proximal parameter and g = ∇F (Y ) is the gradient of F (X) at Y.

4.A.1. Algorithmic Framework
For convenience, we present the model (10) again and define F (X) and G(X) as follows Then, we can conclude that each iteration of the TNNR-APGL for solving model (10) requires solving the following subproblems.
During the above iterate scheme, we update τ k+1 and Y k+1 via the approaches mentioned in [20,35]. Then, based on (21), we can easily drive the TNNR-APGL algorithmic framework as follows.

4.A.2. The Analysis of Subproblems
Obviously, the computation of each iteration of the TNNR-APGL approach for (10) is dominated by the subproblem (22). According to Theorem 3.2, we get Then, the closed-form solution of (22) is given by By now, we have applied TNNR-APGL to solve the problem (10) and obtain closed-form solutions. In addition, the convergence of APGL is well studied in [35] and it has a convergence rate of O( 1 k 2 ). In our paper, we also omit the convergence analysis.

5.
TNNR-ADMMAP for (9) and (11) While the TNNR-ADMM is usually very efficient for solving the TNNR based models (9) and (11), its convergence could become slower with more constraints in [37]. Inspired by [1], the alternating direction method of multipliers with adaptive penalty (ADMMAP) is applied to reduce the constrained conditions, and adaptive penalty [38] is used to speed up the convergence. The resulted algorithm is named as "TNNR-ADMMAP", whose subproblems can also get closed-form solutions.

5.A. Algorithm Framework
Two kinds of constrains have been mentioned as before: Our goal is to transform (9) and (11) into the following form: where P and Q are linear mapping, x, y and c could be either vectors or matrices, and F and G are convex functions.
In order to solve problems easily, b was asked to be a vector formed by stacking the columns of matrices. On the contrary, if A is a linear mapping containing sampling process, we can put AY = b into a matrix form sample set. Correspondingly, we should flexibly change the form between matrices and vectors in the calculation process. Here, we just provide the idea and process of TNNR-ADMMAP. Now, we match the relevant function to get the following results: where P and Q : R m×n → R 2m×2n and C = C1 + C2. Denote B δ,2 = {ζ ∈ R p : ζ ≤ δ} and ξ ∈ R m×n that is the matrix form of ζ = AX − b. When δ = 0, it reflects the problem (9).
Then, the problems (9) and (11) can be equivalently transformed to min X X * − T r(L l r Y R l r ⊺ ), s.t P (X) + Q(Y ) = C.
(24) So the augmented Lagrangian function of (24) is: The Lagrangian form can be solved via linearized ADMM and a dynamic penalty parameter β is preferred in [38]. In particular, due to the special property of A (AA * = I), here, we use ADMMAP in order to handle the problem (25) easily. Similarly, we use the following adaptive updating rule on β [38]: where β max is an upper bound of {β k }. The value of ρ is defined as where ρ 0 ≥ 1 is a constant and ε is a proximal parameter. In summary, the iterative scheme of the TNNR-ADMMAP is as follows: Algorithm 5: TNNR-ADMMAP for (9) and (11) 1. Initialization: set X 1 = Data (the matrix form of b) ,Y 1 = X 1 , Z 1 = zeros(2m, 2n), and input β 0 , ε, ρ 0 .

5.B. The Analysis of Subproblems
Since the computation of each iteration of the TNNR-ADMMAP method is dominated by solving the subproblems (27) and (28), we now elaborate on the strategies for solving these subproblems. First, we compute X k+1 . Since the special form of P and Q, we can give the following solution by ignoring the constant term: Second, we concentrate on computing Y k+1 . Obviously, Y k+1 obeys the following rule: It can be solved as: (31) where Q * is the adjoint operator of Q which is mentioned in (13). Let where W ij ∈ R m×n , according to (13), we have Q(Y ), W = Y, Q * (W ) . More specially, Thus, the adjoint operator Q * is denoted as The left side in (31) can be shown as Then, we apply the linear mapping A (AA * = I) on both sides of (32), and we obtain Combining (32) and (33), we get Similarly, according to the property of Q * in (32), we can get the transformation for the right side in (31) Based on the above from (32) and (34), we achieve Some remarks are in order.
• The convergence of the iterative schemes of AD-MMAP is well studied in [38,39]. In our paper, we omit the convergence analysis.
Overall, when TNNR-ADMM and TNNR-APGL are applied to solve (9)-(11), the generated subproblems all have closed-form solutions. As mentioned before, TNNR-ADMMAP is used to speed up the convergence of (9) and (11) when there are too many constraints. When one problem can be solved simultaneously with the three algorithms, TNNR-ADMMAP is in general more efficient, in the case of matrix completion [1], and in our test problems.

Experiments and Results
In this section, we present numerical results to validate the effectiveness of SVE and LRISD. In summary, there are two parts certified by the following experiments. On one hand, we illustrate the effectiveness of SVE on both real visual and synthetic data sets. On the other hand, we also illustrate the effectiveness of LRISD which solves TNNR based low-rank matrix recovery problems on both synthetic and real visual data sets. Since the space is limited, we only discuss the model (9) using LRISD-ADMM in our experiments. If necessary, you can refer to [1] for extensive numerical results to understand that ADMMAP is much faster than APGL and ADMM without sacrificing the recovery accuracy, in cases of matrix completion. Similarly, for the low-rank matrix recovery, we have the same conclusion according to our experiments. Since the main aim of the paper is to present the effectiveness of SVE and LRISD, here, we omit the detailed explanation.
All experiments were performed under Windows 7 and Matlab v7.10 (R2010a), running on a HP laptop with an Intel Core 2 Duo CPU at 2.4 GHz and 2GB of memory.

6.A. Experiments and Implementation Details
We conduct the numerical experiments under the following four classes, where two representative linear mapping A: matrix completion and partial DCT, are used. The first two cases are to illustrate the effectiveness of SVE.
Here we compared our algorithm LRISD-ADMM with that proposed in [1], which we name as "TNNR-ADMM-TRY", on the matrix completion problem. The main difference between LRISD-ADMM and TNNR-ADMM-TRY is the way of determining the best r. The former is to estimate the best r via SVE while the latter one is to try all the possible r values and pick the one of the best performance.
The last two are to show the better recovery quality of LRISD-ADMM compared with the solution of the common nuclear norm regularized low-rank recovery models, for example, (2), whose corresponding algorithm is denoted as LR-ADMM as above.
(1) Compare LRISD-ADMM with TNNR-ADMM-TRY on matrix completion problems. These experiments are conducted on real visual data sets.
(2) Compare the real rank r withr which is estimated by SVE under different situations, where A is a two-dimensional partial DCT operator (AA * = I). These experiments are conducted on synthetic data sets.
(3) Compare LRISD-ADMM with LR-ADMM on the generic low-rank situations, where A is also a twodimensional partial DCT operator. These experiments are conducted on synthetic data sets under different problem settings.
(4) Compare LRISD-ADMM with LR-ADMM on the generic low-rank situations, where A is also a two-dimensional partial DCT operator. These experiments are conducted on real visual data sets.
In all synthetic data experiments, we generate the sample data as follows: b = AX * + ω, where ω is Gaussian white noise of mean zeros and standard deviation std. The MATLAB script for generating X * as follows X * = randn(m, r) * randn(r, n).
where r is a prefixed integer. Moreover, we generate the index set Ω in (5) randomly in matrix completion experiments. And, the partial DCT is also generated randomly.
In the implementation of all the experiments, we use the criterion X l+1 − X l 2 F / Data 2 F ≤ 10 −2 to terminate the iteration of Step 2 (in Algorithm 1) in LR-ADMM and LRISD-ADMM. In addition, we terminate the iteration of b) in Step 2 by the criterion: In our experiments, we set β = 0.001 empirically, which works quit well for the tested problems. The other parameters in TNNR-ADMM are set to their default values (we use the Matlab code provided online by the author of [1]). Besides, we use the PSNR (Peak Signal to Noise Ratio) to evaluate the quality of an image. As color images have three channels (red, green, and blue), PSNR is defined as 10 × log 10 ( 255 2 MSE ). Where M SE = SE/3T , SE = error 2 red + error 2 green + error 2 blue and T is the total number of missing pixels. For grayscale images, PSNR has the similar definition.

6.B. The Comparison between LRISD-ADMM and TNNR-ADMM-TRY on Matrix Completion
In this subsection, to evaluate the effectiveness of SVE, we compare the proposed LRISD-ADMM with TNNR-ADMM-TRY as well as LD-ADMM on matrix completion problems. As the better recovery quality of TNNR-ADMM-TRY than LR-ADMM on the matrix completion problem has been demonstrated in [1], we will show that the final estimated r of LRISD-ADMM via SVE is very close to the one of TNNR-ADMM-TRY.
We test three real clear images and present the input images, the masked images and the results calculated via three different algorithms: LR-ADMM, TNNR-ADMM-TRY and LRISD-ADMM. The recovery images are showed in Fig 2 and the numerical value comparison about time and PSNR are shown in Table 1. We can see that compared to LR-ADMM, both TNNR-ADMM-TRY and LRISD-ADMM achieve better recovery quality as expected. While the TNNR-ADMM-TRY achieves the best recovery quality as expected due to its trying every possible r value, its running time is extremely longer than LR-ADMM. Our proposed LRISD-ADMM can achieve almost the same recovery quality as TNNR-ADMM-TRY, but with a significant reductions of computation cost. In fact, if the best precision is expected, we use the estimatedr by LRISD-ADMM as a reference to search the best r around it, with the reasonably extra cost of computation. This will we can also find the best r, which is the same as TNNR-ADMM-TRY. Here, for convenience in notation, we name the process LRISD-ADMM-ADJUST in Table 1 and Fig 2. 6.C. The Effectiveness Analysis of SVE in Twodimensional Partial DCT: Synthetic Data In subsection 6.B, we showed that the estimatedr by LRISD-ADMM is very close to the best r by TNNR-ADMM-TRY, on matrix completion problems. Here we further confirm the effectiveness of the proposed SVE of LRISD-ADMM, by conducting experiments for the generic low-rank operator A on synthetic data, where A is a two-dimensional partial DCT (discrete cosine transform) operator. We compare the best r (true rank) with the estimatedr under different settings.
In the results below, r,r, and sr denote, the rank of the matrix X * , estimated rank, and sample ratio taken, respectively. We set the noise level std = 0.9, the sample ratios sr = 0.5 and choose κ = 10 in all of the tests. The reason of setting std = 0.9 is that we want to well illustrate the robustness of SVE to noise. Next, we compare r withr under different settings. For each scenario, we generated the model by 3 times and report the results. As shown in Fig 3, the proposed SVE performs the rationality and effectiveness to estimater in Step 1 in Algorithm 1. Even if there is much noise, this method is still valid, namely,r is (approximately) equivalent to the real rank. That is to say, the proposed SVE is pretty robust to the corruption of noise on sample data. In practice, we can achieve the ideal results in other different settings. To save space, we only illustrate the effectiveness of SVE using the above-mentioned situations.

6.D. The Comparison between LRISD-ADMM and LR-ADMM: Synthetic Data
In this subsection, we compare the proposed LRISD-ADMM with LR-ADMM on partial DCT data in general low rank matrix recovery cases. We will illustrate some numerical results to show the advantages of the proposed LRISD-ADMM in terms of better recovery quality.
We evaluate the recovery performance by the Relative Error as Reer = X re − X * F / X * F . We compare the reconstruction error under different conditions: different noise levels (std) and different sample ratios (sr) taken which are shown respectively in Fig 4(a) and Fig 4(b). In addition, we compare the recovery ranks which are obtained via the above two algorithms in Fig 4(c). For each scenario, we generated the model by 10 times and report the average results.
• We fix the matrix size to be m = n = 300, r = 15, sr = 0.5, and run LRISD-ADMM and LR-ADMM under different noise levels std. The results are shown in Fig 4(a) • We fix the matrix size to be m = n = 300, std = 0.5, and run LRISD-ADMM and LR-ADMM under different sr. The results are shown in Fig 4(b) • We set the matrix size to be m = n = 300, sr = 0.5, std = 0.5, and run LRISD-ADMM and LR-ADMM under different r. The results are shown in Fig 4(c) It is easy to see from Fig 4(a), as the noise level std increases, the total Reer becomes larger. Even so, LRISD-ADMM can achieve much better recovery performance than LR-ADMM. This is because the LRISD model better approximate the rank function than the nuclear norm. Thus, we illustrate that LRISD-ADMM is more robust to noise when it deals with low-rank matrices. With the increasing of sample ratio sr, the total Reer reduces in Fig 4(b). Generally, LRISD-ADMM does better than LR-ADMM, because LRISD-ADMM can approximately recover the rank of the matrix as showed in Fig 4(c).

6.E. The Comparison between LRISD-ADMM and LR-ADMM: Real Visual Data
In this subsection, we test three images: door, window and sea, and compare the recovery images by LRISD-ADMM and general LR-ADMM on the partial DCT operator. In all tests, we fix sr = 0.6, std = 0. The SVE process during different stages to obtainr is depicted in Fig 5. For three images, we set κ = 100, 125, 20 and generater = 8, 10, 7 in LRISD-ADMM. Moreover, Fig 6 shows the recovery results of the two algorithms.
As illustrated in Fig 5, it is easy to see SVE returns an stabler in merely three iterations. And, the estimatedr is a good estimate to the number of largest few singular values. From Fig 6, it can be seen that LRISD-ADMM outperforms the general LR-ADMM in terms of smaller PSNR. More important, using eyeballs, we can see the better fidelity of the recoveries of LRISD-ADMM to the true signals, in terms of better recovering sharp edges.

6.F. A note on κ
We note that the thresholding κ plays a critical rule for the efficiency of the proposed SVE. For real visual data, we can use κ = √ m * n/(3 * s), s = 0.5, 0.8, 1, 3, 5. For synthetic data, κ is denoted as κ = s * √ m * n/30, s = 1, 2, 3. The above heuristic, which works well in our experiments, is certainly not necessarily optimal; on the Table 1: We compare the PSNR and time between LR-ADMM, TNNR-ADMM-TRY, LRISD-ADMM and LRISD-ADMM-ADJUST. In TNNR-ADMM-TRY, they search the best r via testing all possible values of r. We use the estimated rankr as the best r in LRISD-ADMM. In LRISD-ADMM-ADJUST, we use the estimated rankr as a reference to search the best r.  images here. The first column are original images. The second column are masked images. The masked images is obtained by covering 50% pixels of the original image in our test. The third column depicts images recovered by LR-ADMM. The fourth column depicts images recovered by TNNR-ADMM-TRY and LRISD-ADMM-ADJUST. The fifth column depicts images recovered by LRISD-ADMM where we just use the estimatedr directly. Noticing the fourth column, we get the same image by applying two different methods TNNR-ADMM-TRY and LRISD-ADMM-ADJUST. The reason is the values of r calculated in the two methods are the same. But the procedure they find r is different. In TNNR-ADMM-TRY, they search the best r via testing all possible values (1)(2)(3)(4)(5)(6)(7)(8)(9)(10)(11)(12)(13)(14)(15)(16)(17)(18)(19)(20). In LRISD-ADMM-ADJUST, we use the estimated rankr as a reference to search around for the best r.
other hand, it has been observed that LRISD is not very sensitive to κ. Of course, the "last significant jump"  based thresholding is only one way for estimate the number of the true nonzero (or large) singular values, and one can try other available effective jump detection methods [2,40,41].

Conclusion
This paper introduces the singular values estimation (SVE) to estimate a appropriate r (in X r ) that the estimated rank is (approximately) equivalent to the best rank. In addition, we extend TNNR from matrix completion to the general low-rank cases (we call it LRISD). Both synthetic and real visual data sets are discussed. Notice that SVE is not limited to thresholding. Effective support detection guarantees the good performance of LRISD. Therefore future research includes studying specific signal classes and developing more effective support detection methods. We can find the observed estimatedr are 7, 9, 8. Compared to the best r, which are 8,10,7, estimatedr is approximately equivalent to the best r.