Efficient Generalized Conditional Gradient with Gradient Sliding for Composite Optimization

1 Introduction

This paper studies unconstrained composite optimization problems of the form

where X is a general vector space equipped with inner product (e.g. Rd equipped with ℓ2 norm). l(x) is a smooth convex loss function, i.e. it is continuously differentiable with L-Lipschitz continuous gradient:

r(x) stands for the regularizer which is a nonsmooth closed proper convex function. Problem (1) is of vital importance in machine learning because many sparse estimation problems fit into this model. For example, in convex relaxed low rank matrix completion problem [Cand`es and Recht, 2009], l(x) is the Frobenius norm between the observed matrix and the low rank estimation matrix, and r(x) is the trace norm of the low rank estimation matrix. In regularized empirical risk minimization problem, l(x) can be logistic loss or least square loss, and r(x) can vary from simple ones like ℓ1 norm [Tibshirani, 1996] to very complex form such as graph-guided lasso [Kim and Xing, 2009] and group lasso [Yuan and Lin, 2006] for inducing structured sparsity.

Many different lines of methods exist for solving such sparse optimization problems, among them is the popular proximal gradient (PG) based approach ([Beck and Teboulle, 2009]; [Nesterov, 2013]). This kind of methods can achieve the optimal rate of convergence under certain problem settings, hence it enjoys low iteration complexity. The periteration cost mainly comes from gradient evaluation and a proximal map (PM) related to the type of the regularizer. On the one hand, due to the optimal iteration complexity, the number of gradient evaluations is optimal for PGmethod. On the other hand, the proximal map itself is a quadratic optimization problem composed by r(x) and a quadratic term. For particular regularizers, it admits efficient evaluation. For example, the proximal map related to the lasso regularizer is simply the soft-thresholding. However, in some cases, evaluating the proximal map can be quite computational demanding. For instance, to solve the proximal map related the trace norm regularizer, it requires a full singular value decomposition (SVD) in each iteration [Cand`es and Recht, 2009]. As a result, the high per-iteration cost raised by the proximal map becomes the bottleneck of the PG method for those problems.

To address the high per-iteration cost raised by the PM, generalized conditional gradient method (GCG) ([Jaggi, 2013]; [Yu et al., 2014]; [Jaggi, 2011]; [Clarkson, 2010]), has been receiving increasing research interest. It only requires to evaluate a linear operator (LO) in each iteration, which is intuitively much easier to solve than the quadratic subproblem of PM. In fact, this intuition is elaborated by many structured sparse regularizer, e.g. LO requiring spectral norm versus PM requiring full SVD for trace norm regularized problem, in which the former is much more computational efficient [Yu et al., 2014]. As a result, although the iteration complexity of GCG based methods are inferior to their PG counterpart [Dudik et al., 2012], some studies have found that the low per-iteration cost can sometimes compensate for the extra iterations leading to less overall execution time than PG. Nevertheless, the inferior iteration complexity leads to extra gradient evaluations. Although the per-iteration cost for LO can be small enough to afford excess iterations, the increased demand for gradient evaluation can raise an inevitable trade-off, especially for large scale problems. In addition, the convergence results become weaker when approximation techniques are introduced to solve LO. As a result, these approximation techniques also bring about an increased count of gradient evaluations, which can counteract the efficiency gained from the accelerated LO evaluations.

In this paper, we therefore propose a novel algorithm called Generalized Conditional Gradient with Gradient Sliding (GCG-GS) and its refined variant for gauge regularized problems. We first extend a recent optimization scheme called gradient sliding to general unconstrained composite convex optimization problems. Instead of evaluating gradient on each iteration, we skip it from time to time, which can be viewed as many LO evaluations sharing the same gradient value. For gauge regularized problem ([Jaggi, 2013]; [Yu et al., 2014]), where efficient approximation techniques exist for handling r(x), we present an improved variant of the general GCG-GS algorithm to incorporate these techniques to further accelerate the algorithm. As a result, our algorithm has optimal count of gradient evaluations as their PG counterpart, and more importantly, it allows efficient approximation techniques to be used without increasing the optimal count of gradient evaluation. Experiments of CUR-like matrix factorization problem with group lasso penalty on four real-world datasets have demonstrated the efficiency of the proposed method.

2 Background

2.1 Generalized Conditional Gradient

Generalized conditional gradient method (GCG) is suitable for Problem (1). It is a generalization of the classic conditional gradient method (CG) [Frank and Wolfe, 1956], which solves constraint problem on a subset of X , e.g. convex compact set [Jaggi, 2013]. Also, to keep the nonsmooth part intact, GCG defines the following alternative duality gap at iteration k ([Yu et al., 2014]; [Bredies et al., 2005]),

instead of taking the sub-differential of nonsmooth r(x) as some CG methods proposed [Jaggi, 2013]. The duality gap is essential for constructing the linear operator because it is related to the optimal necessary condition and is always an upper bound approximation to the primal gap F (x) - F (x∗ ), as described in the following lemma [Yu et al., 2014]:

Lemma 1. For Problem 1 with any x ∈ X, the duality gap G(x) ≥ 0 and G(x) = 0 iff x satisfies the necessary optimal condition 1 . Also, the duality gap is always an upper bound of the primal gap that G(x) ≥ F (x)-F (x∗), where x∗ denotes the global optimal point.

1 We say that x satisfies the necessary optimal condition if 0 ∈ ∇l(x) + ∂r(x).

Hence, GCG solves Problem (1) by minimizing the duality gap on each iteration, which amounts to evaluating a linear operator (LO), i.e.

Note that additional assumptions should be made to F(x), otherwise the above linear subproblemmay diverge. To avoid introducing additional complexity, we assume that the solution sequences dk and xk (also at and ut to be introduced later) are finite whose maximum distance between each other is upper bounded by a positive constant Ds. This assumption is the same as Assumption 3 in [Yu et al., 2014], where more sophisticated equivalent assumptions are also discussed (see Proposition 3 in [Yu et al., 2014]).

Then, the next step xk+1 can be obtained by

where the step size αk can be set to deterministic sequence for example of order O( 1k ) or choosing by optimizing the following problem,

GCG method has O( 1K ) convergence rate for Problem (1), where K is the total number of iterations. Or equivalently, it needs O( 1ϵ ) iterations to find an ϵ accurate solution. 2 Apparently, this incurs additional count of gradient evaluations than their PG counterpart based on this one gradient evaluation per-iteration scheme.

2.2 Approximation Techniques for Accelerating Linear Operator Evaluation

The linear operator in Eq.(4) admits low per-iteration cost that allows them to afford excess count of iterations, which is the major motivation for adopting GCG method rather than PG method for some tasks. This evaluation can be further accelerated by various approximation techniques. Specifically, we consider the case when r(x) is the so-called generalized gauge function [Yu et al., 2014] defined as

where h(·) is a convex increasing function and κ(x) is a gauge function (convex, positively homogeneous, also see [Friedlander et al., 2014] for detail) satisfying

The convex compact subset C is the unit ball of the gauge κ(x), namely

This formal definition of r(x) is a little complicated. Indeed, one can simply take it as a generalized norm function.

The key idea to accelerate LO evaluation is to circumvent direct computation of κ(x) by approximation techniques.

2 The ϵ accurate solution refers to a solution x having the primal gap no larger than ϵ.

Following [Yu et al., 2014], we convert Eq.(4) to a constrained form, it gives

To avoid estimating ζ, it derives dk by direction and scalar separately. For estimating the direction ak, it follows

We leave the derivation of scalar later when used.

The approximation for accelerating the evaluation of LO lies in two parts. The first is to allow dk or ak to be calculated approximately, i.e. ak is relaxed to satisfy:

This approximation can further reduce the per-iteration cost of LO evaluation, and thus the reduction of overall running time is observed practically [Jaggi, 2013]. However, with such approximation, the iteration complexity guarantee is weakened by a factor associated with the degree of approximation allowed [Jaggi, 2011]. As a result, the count of gradient evaluations is increased which entails a trade-off between the time increase of gradient evaluation and the decrease of LO evaluation, especially for large scale problems whose gradient are expensive to compute.

The other approximation is related to the constraint set. When C is a convex hull of atomic domain A, Eq.(11) can be equivalently solved on A:

where κo (-∇l(xk)) = maxa∈A〈a,-∇l(xk)〉 is called the polar operator [Zhang et al., 2013]. Hence, the evaluation of LO is converted to the evaluation of the polar operator, which is more efficient to deal with than κ(x) itself and is the actual form adopted practically. A simple example of this polar operator is, when κ(x) is a norm, we can immediately find the associated κo (g) is its dual norm.

2.3 Conditional Gradient Sliding Algorithm (CGS)

Paper [Lan and Zhou, 2014] has proposed a gradient sliding technique for constraint smooth objective function to reduce gradient evaluations of CG algorithm. Under CG framework, although the requirement for solving linear operator is still O( 1ϵ ), gradient evaluations can be surprisingly reduced to match the iteration complexity of PG counterpart under the same problem settings.

The CGS separates into outer and inner iteration. On each outer iteration, two additional sequences are maintained. It can be seen as a variant to Nesterov's optimal gradient method [Nesterov, 2013], with the modification of the step (16) calling a CG subroutine rather than a gradient descent or proximal mapping procedure for PG, as shown in the following:

The inner loop, namely the CG subroutine, applies the classic CG algorithm to optimize the following subproblem:

According to CGmethod, on each inner iteration, it optimizes the duality gap:

where ut is the solution sequence of the inner loop subroutine 3 , and D is the constraint set of the problem [Lan and Zhou, 2014] considers. Here, we have slightly extend the notation of the duality gap to incorporate the additional variable v. The subprocess returns the latest ut once the duality gap is less than ηk. That is, it returns utk when

As a result, it can be viewed as if many LO evaluations can share the same gradient to maintain the same convergence rate, instead of updating the gradient for each LO evaluation.

3 The Proposed Algorithm

In this section, we propose our novel algorithm for Problem (1), called Generalized Conditional Gradient with Gradient Sliding (GCG-GS). We first present the GCG-GS for general regularizer. Our algorithm is related to [Lan and Zhou, 2014], but it can suit to more general composite optimization problems with the unconstrained domain. Although for some problems it can be equivalently transformed between regularization and constraint form by Lagrangian duality, the regularization form we consider here allows additional heuristic local optimization, which is hardly known for constraint form. More importantly, we will propose a refined GCG-GS algorithm, which admits various approximation techniques for further accelerating the LO evaluation. When the approximation techniques are involved, neither [Lan and Zhou, 2014] nor our general GCG-GS are applicable. This issue is mainly because the stopping criteria of Eq.(20) is either no longer computable, or computationally expensive to obtain. Our refined GCG-GS will handle this issue.

3.1 General GCG-GS

Algorithm Description:

To involve the gradient sliding scheme, our algorithm also separates into outer and inner loops. In the outer loop, we evaluate the gradient of the smooth part gk = ∇l(zk+1). Then the subroutine is called. As for the subroutine, we consider Φ(v) composed of φ(v) = 〈v, g〉 + β2 ||v - u||2 and nonsmooth r(v), i.e. Φ(v) = φ(v) + r(v). We use a different definition of G(ut, v) according to Eq.(3), which is to be optimized as follows:

Note that G(ut, vt) is the duality gap of Φ(v) at ut. Thus, the subroutine actually solves the minimization problem of Φ(v)

3 We use k to represent sequences related to outer loop and t for inner loop in this paper.

by GCG algorithm until certain duality gap is obtained. In addition, αt can be optimally chosen by solving

(22)

Finally, the next variable is obtained by

For clarity, we summarize the general GCG-GS in Algorithm (1) and Algorithm (2).

Algorithm 1 General-GCG-GS

Require: x0,K, γk, βk, ηk

• 1: for k = 0, 1, ...,K -1 do
• 2: zk+1 = (1-γk)yk +γkxk;
• 3: gk =∇l(zk+1);
• 4: xk+1 = GCG(gk, xk, βk, ηk);
• 5: yk+1 = (1-γk)yk +γkxk+1;
• 6: end for

Ensure: yK

Algorithm 2 GCG: General-GCG-GS subroutine

Require: g, u, β, η

• 1: let φ(x) = 〈g, x〉 + β2 ||x-u||2 , u0 = u;
• 2: for t = 0, 1, ..., do
• 3: G(ut, v) = r(ut) - r(v) - 〈v -ut,∇φ(ut)〉;
• 4: vt = argmaxv∈X G(ut, v);
• 5: if G(ut, vt) ≤ η then
• 6: break;
• 7: end if

• 9: ut+1 = (1-αt)ut +αtvt;
• 10: end for
• 11: Return: u+ = ut;

Convergence Analysis:

For the general GCG-GS algorithm, we have the following convergence guarantee. We first introduce the following notation:

Theorem 1. Under finite solution sequence assumptions, apply GCG-GS to Problem (1),

a) for any x, the output yK satisfies,

b) consider the inner loop, namely the LO evaluations, for a particular stage k, denote β = βk, the duality gap satisfies,

where Ds is the upper bound of the solution sequence.

Corollary 1. With the sequences setting as βk = 2L k+1 , γk = 2 k+2 , ηk = 2LD0 K(k+1) and denoting D0 = ||x0 - x||2 , we have: a)

For finding an ϵ solution, we get:

b) the total number of inner LO evaluations is

For finding an ϵ solution, we have:

Discussion:

There are various assignment of sequence, please see [Lan and Zhou, 2014]. For the particular sequence we adopt here, it is apparent that our number of gradient evaluations for finding an ϵ solution is O(√ 1ϵ ), which is the same as those optimal complexity of PG methods for Problem (1). In addition, in terms of the total number of LO evaluations TK, the proposed method maintains the same order of complexity as those plain GCG methods which is also optimal for GCG.

However, the above algorithm is only conceptual in some sense. Note that both the Subproblems (21) and (22) can be difficult to solve for some r(x). Also, even they were solvable, we often prefer avoiding directly computing r(x) by considering more efficient substitution such as polar operator in Section 2.2. Inspired by this, we will propose the refined GCG-GS algorithm in the next subsection to allow more efficient inner loop execution.

3.2 Refined GCG-GS for Gauge Regularized Problem

In this subsection, we follow the assumption as in Section 2.2, where the regularizer is a generalized gauge function defined by Eq.(7). Note that most practically used sparsity inducing and rank minimization regularizers can be seen as generalized gauge function.

Efficient Approximation Techniques:

Essentially, we apply the efficient approximation techniques introduced in Section 2.2 to minimize Φ(v). To efficiently minimize Eq.(21), we first convert it to constraint form:

Then we update vt by solving the direction and scalar separately, namely αtvt ≈ θtat (θt denotes the approximate scalar). The direction is updated by

Incorporating the approximation Eq.(12) and solving it on atomic domain A, we can obtain at more efficiently by

The scalar, denoted as θt, would be originally chosen as

where αt is a deterministic sequence to be specified in our next theorem. Again, to avoid direct evaluation of κ(ut), an upper substitution ρt is used as in [Yu et al., 2014]. This is achieved by choosing ρ0 ≥ κ(u0) and the update scheme ρt+1 = (1 - αt)ρt + θt. Then ρt ≥ κ(ut) can be held iteratively, see [Yu et al., 2014]. Thus, by alternatively using h(ρt) provided that h(·) is increasing convex and κ(·) is convex, θt can be obtained by

Finally, an additional local heuristic optimization can be adopted to further improve the practical performance, which is another motivation for using the regularized form rather than constrained form. Denoting such re-optimization by Improve, we adopt the following conceptual requirement for it, which is Relaxed assumption according to [Yu et al., 2014]:

Weighted Average as Return Value:

An important issue with the above approximation is that the duality gap Eq.(21) is either no longer computable, or even when we can compute it, it is unreasonable for us to directly evaluate it because we do all the above approximations to avoid computing κ(v) directly. As a result, the stopping criteria in general GCG-GS algorithm (also CGS algorithm) cannot be used. Furthermore, the choice of return value becomes a problem because the previous bound on duality gap only guarantees the minimum one. Again, as we cannot directly compute the duality gap, it also becomes unknown that on which particular ut the duality gap is small enough.

To solve the above stopping criteria problem, we propose a simple alternative by estimating a maximum iteration count m of the inner LO evaluation loop. As shown in our convergence analysis, different outer loops can share the same m.

As for the choice of return value, instead of returning a particular ut, we propose using the weighted average of ut as the returned value. We show such averaged ¯um can guarantee G(¯um, x) to be smaller than the desired η as long as the proper approximated m is used. In detail, the return value to the outer loop is

Intuitively, variables of later iterations gain more weights. This intuition is compatible to the analysis in [Jaggi, 2011], where the one achieving the smallest duality gap lies in the last third iterations. This average scheme is by observing the special construct of φ, namely, the corresponding G(ut, x) is convex in ut for arbitrary yet fixed x. Note that it also allows an online update ¯ut+1 = (1- 2 t+2 )¯ut + 2 t+2 ut, which is exactly what has been shown in the algorithm. We note that [Lacoste-Julien et al., 2013] has used the same weighted average as an update option in block coordinate conditional gradient method. To summarize together, Algorithm (3) and Algorithm (4) show the implementation details of the RefinedGCG-GS.

Algorithm 3 Refined-GCG-GS

Require: x0,m, βk, γk

• 1: for k = 0, 1, ...,K -1 do
• 2: zk+1 = (1-γk)yk +γkxk;
• 3: gk =∇l(zk+1);
• 4: xk+1 = Re-GCG(gk, xk, βk,m);
• 5: yk+1 = (1-γk)yk +γkxk+1;
• 6: end for

Ensure: yK

Algorithm 4 Re-GCG:Refined-GCG-GS subroutine

Require: Input from outer loop: g, u, β,m; Sequence αt

• 1: let φ(x) = 〈g, x〉 + β2 ||x-u||22 , u0 = u;
• 2: for t = 0, 1, ...,m-1 do
• 3: choose at satisfy 〈at,∇φ(ut)〉 ≤ ϵt -κo (-∇φ(ut));
• 4: θt = argminθ φ((1-αt)ut +θat)+(1-αt)h(ρt)+ αth( θ αt ));
• 5: ˜ut+1 = (1-αt)ut + θtat;
• 6: ˜ ρt+1 = (1 - αt)ρt + θt;
• 7: (ut+1, ρt+1) = Improve(˜ut+1, ˜ ρt+1, φ, r);
• 8: ¯ut+1 = (1- 2 t+2 )¯ut + 2 t+2 ut+1;
• 9: end for
• 10: Return: u+ = ¯um;

Convergence Analysis: 4

Theorem 2. Let the sequence settings be βk = 2L k+1 , γk = 2 k+2 , αt = 2 t+2 and the inner loop count m = ⌈ 6K(Ds+δκ(x)) D0 ⌉, where Ds is an upper bound on the distance of the solution path, D0 is the distance between x0 and x, δ is the constant satisfies ϵt ≤ δβkαt 2 . Then the following results hold for the above algorithm for any x in X . a) After K outer loops, the output yK satisfies

4 We provide the proof in an anonymous dropbox folder for interested reader. https://www.dropbox.com/sh/xhfmrxr0xss1mxl/AAC00hytB5tncUssuKRMyxca?dl=0

for finding an ϵ solution, and the number of FO evaluation requires

b) The total number of LO evaluation TK for finding an ϵ solution requires

Discussion:

The outer loop complexity in this subsection is exactly the same as the one in the previous subsection, despite all the approximation we make to efficiently evaluate LO. As a result, the count of gradient evaluation keeps unchanged. As a sharp comparison, the convergence rate of [Jaggi, 2011] is degenerated by a factor of 2 when dk is evaluated approximately. Apparently, the approximations made to the LO lead to the increasing count of evaluations of gradient.

We point out that, by properly restarting GCG-GS, this algorithm can also obtain optimal count of gradient evaluation for strongly convex problem [Lan and Zhou, 2014]. In this paper, we only discuss the problem under convex assumption due to space limitation. In fact, the extension to strongly convex case is straightforward in some sense.

4 Experiment

In this section, we demonstrate the efficiency of the proposed algorithm by a CUR-like matrix factorization task ([Mahoney and Drineas, 2009]; [Mairal et al., 2011]) regularized by group lasso penalty. This experiment was conducted by using MATLAB on a laptop computer of Intel Core i7 2.7GHz processor with 8 GB RAM.

Experimental Setup:

We consider the following CUR-like matrix factorization problem [Mairal et al., 2011].

where D is the input data matrix, Xi: and X:j denote the row vectors and column vectors, respectively, and || · ||∞ is the max norm. We set λ = 5 × 10-4 in our experiment. We utilized the following four real datasets as used in [Yu et al., 2014]: SRBCT, Brain Tumor 2, 9 Tumor and Leukemia 5 , which are of sizes 83 × 2308, 50 × 10367, 60 × 5762, and 72 × 11225, respectively. For comparison, we utilized GCG TUM algorithm in [Yu et al., 2014] 6 . For our GCGGS algorithm, we implemented the outer loop routine, which then called the same polar operator of GCG TUM for inner loop subroutine. Hence, the improved performance is gained purely from the gradient sliding scheme. We set our inner loop estimation m to 3 for all four datasets. Other input sequences were assigned exactly as the theoretical part. Note that we did not compare with PG based methods because they have already been shown to be less efficient than GCG TUM in [Yu et al., 2014].

5 Download from http://www.gems-system.org.

6 Download from http://users.cecs.anu.edu.au/∼xzhang/GCG.

Experimental Result:

Figure 1 shows the experimental results, where three sets of plots are drawn: objective function value versus total running time, versus polar operator evaluation time only and versus gradient evaluation time only. We sampled every 30 iterations for GCG TUM and every 10 outer loop iterations for GCG-GS. In general, our algorithm is much faster than GCG TUM algorithm in terms of convergence speed, as illustrated by the curve of objective function value versus total running time. Also, our algorithm requires much less time on gradient evaluation to achieve certain decrease of objective function value on all four datasets. Finally, the time requirements for polar operator evaluation of our algorithm are similar to GCG TUM on Brain Tumor 2 and SRBCT, superior than GCG TUM on Leucamia and inferior than GCG TUM on 9 Tumor.

Choice of Inner Iteration Count m:

In this subsection, we study the effect of different estimates of m. We run the GCG-GS algorithm with 400 outer loops on the four datasets, while vary the maximum inner iteration from 2 to 7. Figure 2 plots the objective function value versus number of iterations. With the different inner iteration number, the algorithm actually converges at similar outer loop count. To be specific, all lines begin to converge around 50 numbers of outer iteration. Although the algorithm converges to different objective function value with different m, the difference is below 0.005. In addition, we observed that m = 3 always yields relatively superior performance. Hence, m is not hard to tune practically.

5 Conclusion

In this paper, we have proposed a new algorithm under GCG framework. Our algorithm has optimal count of gradient evaluations as those PG method, which is an order superior than plain GCG methods. Also, it admits the incorporation of efficient approximation techniques for accelerating the evaluation of linear operator that CGS lacks. Meanwhile, this count of gradient requirement remains unchanged. Experiment on a CUR-like matrix factorization task with group lasso penalty on four real datasets have demonstrated the efficiency of the proposed method.

Acknowledgments

This work was supported by the Faculty Research Grant of Hong Kong Baptist University under Projects: FRG2/1415/075 and FRG1/14-15/041, and the National Science Foundation of China under Grant: 61272366.