Efficient Inexact Proximal Gradient Algorithm for Nonconvex Problems

1 Introduction

In regularized risk minimization, we consider optimization problems of the form

where f is the loss, and g is the regularizer. Typically, f is smooth and convex (e.g., square and logistic losses), and g is convex but may not be differentiable (e.g., ℓ1 and nuclear norm regularizers). The proximal gradient (PG) algorithm [Parikh and Boyd, 2014], together with its accelerated variant (APG) [Beck and Teboulle, 2009b; Nesterov, 2013], have been popularly used for solving this convex problem. Its crux is the proximal step proxg (·) = argminx 12 ‖x-·‖22 +ηg(x), which can often be easily computed in closed-form.

While convex regularizers are easy to use, the resultant predictors may be biased [Zhang, 2010]. Recently, there is growing interest in the use of nonconvex regularizers, such as the log-sum-penalty [Cand`es et al., 2008] and capped ℓ1-norm [Zhang, 2010] regularizers. It has been shown that these often lead to sparser and more accurate models [Gong et al., 2013; Lu et al., 2014; Zhong and Kwok, 2014; Yao et al., 2015]. However, the associated proximal steps become more difficult to compute analytically, and cheap closed-form solutions exist only for some simple nonconvex regularizers [Gong et al., 2013]. This is further aggravated by the fact that the state-of-the-art PG algorithm for nonconvex optimization, namely the nonmonotone accelerated proximal gradient (nmAPG) algorithm [Li and Lin, 2015], needs more than one proximal steps in each iteration.

When the optimization objective is convex, one can reduce the computational complexity of the proximal step by only computing it inexactly (i.e., approximately). Significant speedup has been observed in practice, and the resultant inexact PG algorithm has the same convergence guarantee as the exact algorithm under mild conditions [Schmidt et al., 2011]. However, on nonconvex problems, the use of inexact proximal steps has not been explored. Moreover, convergence of nmAPG hinges on the use of exact proximal steps.

In this paper, we propose a new PG algorithm for nonconvex problems. Unlike nmAPG, it performs only one proximal step in each iteration. Moreover, the proximal step can be inexact. The algorithm is guaranteed to converge to a critical point of the nonconvex objective. Experimental results on nonconvex total variation models and nonconvex low-rank matrix learning show that the proposed algorithm is much faster than nmAPG and other state-of-the-art, while still producing solutions of comparable quality.

The rest of the paper is organized as follows. Section 2 provides a brief review on the PG algorithm and its accelerated variant. The proposed algorithm is described in Section 3, and its convergence analysis studied in Section 4. Experimental results are presented in Section 5, and the last section gives some concluding remarks.

2 RelatedWork

In this paper, we assume that f in (1) is L-Lipschitz smooth (i.e., ‖∇f(x) - ∇f(y)‖2 ≤ L‖x - y‖2), and g is proper, lower semi-continuous. Besides, F = f +g in (1) is bounded from below, and lim‖x‖2→∞F(x) = ∞. Moreover, both f and g can be nonconvex.

First, we consider the case where f and g in (1) are convex. At iteration k, the accelerated proximal gradient (APG) algorithm generates xk+1 as

where θk = k-1 k+2 and η is the stepsize [Beck and Teboulle,

2009b; Nesterov, 2013]. When θk = 0, APG reduces to the plain PG algorithm.

On nonconvex problems, yk can be a bad extrapolation and the iterations in (2), (3) may not converge [Beck and Teboulle, 2009a]. Recently, a number of PG extensions have been proposed to alleviate this problem. The iPiano [Ochs et al., 2014], NIPS [Ghadimi and Lan, 2015], and UAG [Ghadimi and Lan, 2015] algorithms allow f to be nonconvex, but still requires g to be convex. The GD algorithm [Attouch et al., 2013] also allows g to be nonconvex, but does not support acceleration.

The current state-of-the-art is the nonmonotone APG (nmAPG) algorithm 1 [Li and Lin, 2015], shown in Algorithm 1. It allows both f and g to be nonconvex, and also uses acceleration. To guarantee convergence, nmAPG ensures that the objective is sufficiently reduced in each iteration:

where vk+1 = proxηg (xk - η∇f(xk)) and δ > 0 is a constant. A second proximal step has to be performed (step 8) if a variant of (4) is not met (step 5).

Algorithm 1 Nonmonotone APG (nmAPG).

Require: choose η ∈ (0, 1/L), a positive constant δ, ∆1 = F(x1), q1 = 1, and ν ∈ (0, 1);

• 1: x0 = x1 = xa1 = 0 and t0 = t1 = 1;
• 2: for k = 1, . . . ,K do

• 4: xak+1 = proxηg (yk - η∇f(yk));
• 5: if F(xak+1 ) ≤ ∆k - δ2 ‖xak+1 - yk‖22 then
• 6: xk+1 = xak+1 ;
• 7: else

• 10: end if
• 11: qk+1 = νqk +1;

• 13: ∆k+1 = 1 qk+1 (νqk∆k + F (xk+1));
• 14: end for
• 15: return xK+1.

3 Efficient APG for Nonconvex Problems

The proposed algorithm is shown in Algorithm 2. Following [Schmidt et al., 2011], we use a simpler acceleration scheme in step 3. Efficiency of the algorithm comes from two key ideas: reducing the number of proximal steps to one in each iteration (Section 3.1); and the use of inexact proximal steps (Section 3.2). Besides, we also allow nonmonotone update on the objective (and so F(yk) may be larger than F(xk)). This helps to jump from narrow curved valley and improve convergence [Grippo and Sciandrone, 2002;

1 A less efficient monotone APG (mAPG) algorithm is also proposed in [Li and Lin, 2015].

Wright et al., 2009; Gong et al., 2013]. Note that when the proximal step is inexact, nmAPG does not guarantee convergence as its Lemma 2 no longer holds.

Algorithm 2 Noconvex inexact APG (niAPG) algorithm.

Require: choose η ∈ (0, 1L ) and δ ∈ (0, 1 η -L); 1: x0 = x1 = 0; 2: for k = 1, . . . ,K do 3: yk = xk + k-1 k+2 (xk - xk-1); 4: ∆k =maxt=max(1,k-q),...,k F(xt); 5: if F(yk) ≤ ∆k then 6: vk = yk; 7: else 8: vk = xk; 9: end if 10: zk = vk - η∇f(vk); 11: xk+1 = proxηg (zk); // possibly inexact 12: end for 13: return xK+1.

3.1 Using Only One Proximal Step

Recall that on extending APG to nonconvex problems, the key is to ensure a sufficient decrease of the objective in each iteration. Let ρ = 1/η. For the standard PG algorithm with exact proximal steps, the decrease in F can be bounded as follows.

In nmAPG (Algorithm 1), there is always a sufficient decrease after performing the (non-accelerated) proximal descent from xk to xpk+1 (Proposition 3.1), but not necessarily the case for the accelerated descent from xk to xak+1 (which is generated by a possibly bad extrapolation yk). Hence, nmAPG needs to perform extra checking at step 5. If the condition fails, xpk+1 is used instead of xak+1 in step 9.

As the main problem is on yk, we propose to check F(yk) (step 5 in Algorithm 2) before the proximal step (step 11), instead of checking after the proximal steps. Though this change is simple, the main difficulty is how to guarantee convergence while simultaneously maintaining acceleration and using only one proximal step. As will be seen in Section 4.1, existing proofs do not hold even with exact proximal steps.

The following shows that a similar sufficient decrease condition can still be guaranteed after this modification.

In step 4, setting q = 0 is the most straightforward. The F(yk) value is then checked w.r.t. the most recent F(xk). The use of a larger q is inspired from the Barzilai-Borwein scheme for unconstrained smooth minimization [Grippo and Sciandrone, 2002]. This allows yk to occasionally increase the objective, while ensuring F(yk) to be smaller than the largest objective value from the last q iterations. In the experiments, q is set to 5 as in [Wright et al., 2009; Gong et al., 2013].

3.2 Inexact Proximal Step

Proposition 3.2 requires exact proximal step, which can be expensive. Inexact proximal steps are much cheaper, but the inexactness has to be carefully controlled to ensure convergence. We will propose two such schemes depending on whether g is convex. Note that f is not required to be convex.

Convex g. As g is convex, the optimization problem associated with the proximal step is also convex. Let hηg(x) ≡ 12 ‖x-zk‖22 +ηg(x) be the objective in the proximal step. For any zk, the dual of the proximal step at zk can be obtained as

where g∗ is the convex conjugate of g. In an inexact proximal step, the obtained dual variable ˜wk only approximately maximizes (5). The duality gap εk ≡ hηg(xk+1)-Dηg( ˜wk), where xk+1 = zk - η ˜wk, upper-bounds the approximation error of the inexact proximal step ϵk ≡ hηg(xk+1) - hηg (proxηg (zk)). To ensure the inexactness ϵk to be smaller than a given threshold τk, we can control the duality gap as εk ≤ τk [Schmidt et al., 2011].

The following shows that xk+1 satisfies a similar sufficient decrease condition as in Proposition 3.2. Note that this cannot be derived from [Schmidt et al., 2011], which relies on the convexity of f .

Nonconvex g. When g is nonconvex, the GD algorithm [Attouch et al., 2013] allows inexact proximal steps. However, it does not support acceleration, and nonmonotone update. Thus, its convergence proof cannot be used here.

As g is nonconvex, it is difficult to derive the dual of the corresponding proximal step, and the optimal duality gap may also be nonzero. Thus, we monitor the progress of F instead. Inspired by Proposition 3.1, we require xk+1 from an inexact proximal step to satisfy the following weaker condition:

where δ ∈ (0, ρ - L). This condition has also been used in the GD algorithm. However, it requires checking an extra condition which is impractical. 2

We could have also used condition (6) when g is convex. However, Proposition 3.3 offers more precise control, as it can recover (6) by setting εk = ρ-L-δ 2ρ ‖xk+1 - vk‖22 (note that δ < ρ - L). Besides, the duality gap εk is readily produced by primal-dual algorithms, and is often less expensive to compute than F .

2 Specifically, the condition is: ∃ wk+1 ∈ ∂g(xk+1) such that ‖wk+1 + ∇f(vk)‖22 ≤ b‖xk+1 - vk‖22 for some constant b > 0. However, the subdifferential ∂g(xk) is in general difficult to compute.

4 Convergence Analysis

Definition 4.1 ([Attouch et al., 2013]). The Frechet subdifferential of F at x is

The limiting subdifferential (or simply subdifferential) of F at x is ∂F(x) = {u : ∃xk → x, F(xk) → F(x), uk ∈ ˆ ∂F(xk)→u, as k →∞}.

Definition 4.2 ([Attouch et al., 2013]). x is a critical point of F if 0 ∈ ∇f(x) +∂g(x).

4.1 Exact Proximal Step

In this section, we will show that Algorithm 2 (where both f and g can be nonconvex) converges with a O(1/K) rate. This is the best known rate for nonconvex problems with firstorder methods [Nesterov, 2004]. A similar O(1/K) rate for ‖G(vk)‖22 is recently established for APG with nonconvex f but only convex g [Ghadimi and Lan, 2015]. Note also that no convergence rate has been proved for nmAPG [Li and Lin, 2015] and GD [Attouch et al., 2013] in this case. Besides, their proof techniques cannot be used here as their nonmonotone updates are different.

Theorem 4.1. The sequence {xk} generated from Algorithm 2 (with exact proximal step) have at least one limit point, and all limit points are critical points of (1).

Let G(v) = v-proxηg (v-η∇f(v)), the proximal mapping at v [Parikh and Boyd, 2014]. The following Lemma suggests that ‖G(v)‖22 can be used to measure how far v is from optimality [Ghadimi and Lan, 2015].

Lemma 4.2 ([Gong et al., 2013; Attouch et al., 2013]). v is a critical point of (1) if and only if G(v) = 0.

The following Proposition shows that the proposed Algorithm 2 converges with a O(1/K) rate.

Proposition 4.3. Let φ(k) = argmint=max(k-q,1),...,k ‖xt+1 - vt‖22 . (i) limk→∞‖G(vφ(k))‖22 = 0; and (ii) mink=1,...,K ‖G(vφ(k))‖22 ≤ 2(q+1)c1 (ρ-L)K , where c1 = maxt=1,...,q+1 F (xt) - inf F .

4.2 Inexact Proximal Step

Convex g. As in [Schmidt et al., 2011], we assume that the duality gap εk decays as O(1/k1+ς ) for some ς > 0. Let c ≡ ∑∞k=1 εk. Note that c <∞.

Theorem 4.4. The sequence {xk} generated from Algorithm 2 have at least one limit point, and all limit points are critical points of (1).

Proposition 4.5. Let ek ≡ xk+1 - proxηg (xk - η∇f(xk)), the difference between the inexact and exact proximal step solutions at iteration k. We have ‖ek‖22 ≤ 2εk.

Note that the proof techniques in [Schmidt et al., 2011] cannot be used here as f is not required to be convex. As in Proposition 4.3, we also use ‖G(vφ(k))‖22 to measure how far vφ(k) is from optimality.

When all εk's are zero, Proposition 4.6 reduces to Proposition 4.3. In general, the bound of mink=1,...,K ‖G(vφ(k))‖22 in Proposition 4.6 is larger due to the inexact proximal step. Nonconvex g. With inexact proximal steps, nmAPG no longer guarantees convergence, and its proof cannot be easily extended. On the other hand, GD allows inexact proximal steps but uses a different approach to control inexactness. Moreover, it does not support acceleration.

The following shows that Algorithm 2 generates a bounded sequence, and Corollary 4.8 shows that the limit points are critical points.

Theorem 4.7. The sequence {xk} generated from Algorithm 2 has at least one limit point.

Corollary 4.8. Let {xkj } be a subsequence of {xk} with limkj→∞xkj = x∗. If (i) xk+1 ̸= vk unless vk = proxηg (vk - η∇f(vk)), and (ii) limkj→∞ F(xkj ) = F(x∗), then x∗ is a critical point of (1).

Assumption (i), together with Lemma 4.2, ensures that the sufficient decrease condition in (6) will not be trivially satisfied by xk+1 = vk, unless vk is a critical point. Assumption (ii) follows from Definition 4.1, as the subdifferential is defined by a limiting process.

Proposition 4.9. [Attouch et al., 2013] Assumption (ii) is satisfied when (i) the proximal step is exact; or (ii) g is continuous or is the indicator function of a compact set.

When the proximal step is exact or when g is convex, G(·) has been used to measure the distance from optimality. However, this is inappropriate when g is nonconvex and the proximal step is inexact, as the inexactness can no longer be directly controlled. Instead, we will measure optimality via ak ≡ ‖xk+1 - vk‖22 .

When the proximal step is exact, xk+1 = proxηg (vk - η∇f(vk)), and ak = ‖G(vk)‖22 . Proposition 4.10 then reduces to Proposition 4.3 (but with a looser bound).

5 Experiments

In this section, we perform experiments when g is convex (Section 5.1) and nonconvex (Section 5.2).

5.1 Image Inpainting

The total variation (TV) model [Beck and Teboulle, 2009a] has been popularly used in image processing. Let y ∈ Rd be the vectorized input image and x ∈ Rd be the recovered one. We consider the TV model with nonconvex log-sum-penalty regularizer [Cand`es et al., 2008].

where M ∈ {0, 1}d is a mask such that Mij = 1 indicates that the corresponding pixel is observed, Dh and Dv are the horizontal and vertical partial derivative operators, ⊙ is the elementwise multiplication, and κ(α) = log(1 + |α|).

As suggested in [Yao and Kwok, 2016], (7) can be transformed as the minimization of f (x)+λ TV(x), where f (x) = 12 ‖M⊙(x-y)‖22 -λ[ TV(x)+∑di=1 κ([Dhx]i)+κ([Dvx]i)] is nonconvex but smooth, and TV(x) = ‖Dhx‖1 + ‖Dvx‖1 is the standard (convex) TV regularizer. Thus, we only need to handle the proximal step of the TV regularizer, which will be computed numerically by solving its dual using L-BFGS [Beck and Teboulle, 2009a].

The following solvers on the transformed problems are compared: (i) GDPAN [Zhong and Kwok, 2014], which performs gradient descent with the proximal average; (ii) nmAPG; (iii) the proposed niAPG, in which inexactness of the proximal step is controlled by decaying the duality gap εk at a rate of O(1/k1.5 ); and (iv) the exact version of niAPG(exact), which simulates an exact proximal step with a small duality gap (10-4 ). We do not compare with the GD algorithm [Attouch et al., 2013], as its inexactness condition is difficult to check and it does not use acceleration.

As a further baseline, we compare with the convex TV model: minx 1 2 ‖M ⊙ (x - y)‖22 + λ(‖Dvx‖1 + ‖Dhx‖1), which is solved using ADMM [Boyd et al., 2011]. We do not compare with CCCP [Yuille and Rangarajan, 2002], which is slow in practice [Yao and Kwok, 2016].

Experiments are performed on the "Lena" image. We normalize the pixel values to [0, 1], and then add Gaussian noise from N(0, 0.05). 50% of the pixels are randomly sampled as observed. For performance evaluation, we report the CPU time and root-mean-squared error (RMSE) on the whole image. The experiment is repeated five times. Results are shown in Table 1. 3 Figure 1 plots convergence of the objective. 4 As can be seen, the nonconvex TV model has better RMSE than the convex one. Among the nonconvex models, niAPG is much faster, as it only requires a small number of cheap inexact proximal steps (Table 5).

3 In all the tables, the boldface indicates the best and comparable results (according to the pairwise t-test with 95% confidence).

4 Because of the lack of space, the plot for λ = 0.01 is not shown.

5.2 Matrix Completion

In this section, we consider matrix completion with a nonconvex low-rank regularizer. As shown in [Lu et al., 2014; Yao et al., 2015; Yao and Kwok, 2015], it gives better performance than nuclear-norm based and factorization approaches. The optimization problem can be formulated as

where Oijs are the observations, Ωij = 1 if Oi,j is observed, and 0 otherwise, σi(X) is the ith leading singular value of X, and r is the desired rank. The associated proximal step can be solved with rank-r SVD [Lu et al., 2014].

Synthetic Data. The observed m×m matrix is generated as O = UV + G, where the entries of U ∈ Rm×k , V ∈ Rk×m (with k = 5) are sampled i.i.d. from the normal distribution N(0, 1), and entries of G sampled from N(0, 0.1). A total of ‖Ω‖1 = 2mk log(m) random entries in O are observed. Half of them are used for training, and the rest as validation set.

In the proposed niAPG algorithm, its proximal step is approximated by using power method [Halko et al., 2011], and inexactness of the proximal step is monitored by condition (6). Its variant niAPG(exact) has exact proximal steps computed by the Lancoz algorithm [Larsen, 1998]. They are compared with the following solvers on the nonconvex model (8): (i) Iterative reweighted nuclear norm (IRNN) [Lu et al., 2014] algorithm; (ii) Fast nonconvex low-rank learning (FaNCL) algorithm [Yao et al., 2015], using the power method to approximate the proximal step; (iii) nmAPG, in which the proximal step is exactly computed by the Lancoz algorithm. We also compare with other matrix completion algorithms, including the well-known (convex) nuclear-norm-regularized algorithms: (i) active subspace selection [Hsieh and Olsen, 2014] and (ii) ALT-Impute [Hastie et al., 2015]. We also compare with state-of-the-art factorization models (where the rank is tuned by the validation set): (i) R1MP [Wang et al., 2015]; and (ii) state-of-the-art gradient descent based AltGrad [Zhao et al., 2015]. We do not compare with the Frank-Wolfe algorithm [Zhang et al., 2012], which has been shown to be slower [Hsieh and Olsen, 2014].

Testing is performed on the non-observed entries (denoted ¯Ω). Three measures are used for performance evaluation: (i) normalized mean squared error (NMSE): √‖P¯Ω(X -UV )‖2F /√‖P¯Ω(UV )‖2F ; (ii) rank of X; and (iii) training time. Each experiment is repeated five times.

Table 2 shows the performance. Convergence for algorithms solving (8) is shown 5 in Figure 2. As has also been observed in [Lu et al., 2014; Yao et al., 2015], nonconvex regularization yields lower NMSE than nuclear-normregularized and factorization models. Again, niAPG is the fastest. Its speed is comparable with AltGrad, but is more accurate. Table 4 compares the numbers of proximal steps. As can be seen, both nmAPG(exact) and nmAPG require significantly fewer proximal steps.

Recommender Systems. We first consider the MovieLens data sets (Table 6), which contain ratings of different users on movies or musics. We follow the setup in [Wang et al., 2015; Yao et al., 2015; Yao and Kwok, 2015], and use 50% of the observed ratings for training, 25% for validation and the rest for testing. For performance evaluation, we use the root mean squared error on the test set ¯Ω: RMSE = √‖P¯Ω(O - X)‖2F /√‖¯Ω‖1, rank of the recovered matrix X, and CPU time. The experiment is repeated five times.

Table 3 shows the recovery performance. IRNN is not compared as it is too slow. Again, the nonconvex model consistently outperforms nuclear-norm-regularized and factorization models. R1MP is the fastest, but its recovery performance is poor. Figure 3(a) shows the convergence, and Table 4 compares the numbers of proximal steps. niAPG(exact) is faster than nmAPG due to the use of fewer proximal steps. niAPG is even faster with the use of inexact proximal steps.

5 Because of the lack of space, the plot for m = 500 is not shown.

Finally, we perform experiments on the large Netflix and Yahoo data sets (Table 6). We randomly use 50% of the observed ratings for training, 25% for validation and the rest for testing. Each experiment is repeated five times. We do not compare with nuclear-norm-regularized methods as they yield higher rank and RMSE than others. Table 7 shows the recovery performance, and Figures 3(b) show the convergence. Again, niAPG is the fastest and most accurate.

6 Conclusion

In this paper, we proposed an efficient accelerated proximal gradient algorithm for nonconvex problems. Compared with the state-of-the-art [Li and Lin, 2015], the proximal step can be inexact and the number of proximal steps required is significantly reduced, while still ensuring convergence to a critical point. Experiments on image denoising and matrix completion problems show that the proposed algorithm has comparable (or even better) prediction performance as the state-of-the-art, but is much faster.

Acknowledgments

This research project is partially funded by Microsoft Research Asia and the Research Grants Council of the Hong Kong Special Administrative Region (Grant 614513). The first author would thank helpful discussion and suggestions from Lu Hou and Yue Wang.