15.097 Student Project: DC Programming

DC Programming: The Optimization Method You Never Knew You Had To Know May 18, 2012 1 1.1 Introduction What is DC Programming? Recall that in the lectures on Support Vector Machines and Kernels, we repeat­ edly relied on the use of convex optimization to ensure that solutions existed and could be computed. As we shall see later in the lecture, there are many cases in which the assumption that the objective function and constraints are convex (or quasi-convex) are invalid, and so the methods on convex functions that we have developed in class prove insufficient. To deal with these problems, we develop a theory of optimization for a superclass of convex functions, called DC - Difference of Convex - functions. We now define such functions formally. Definition 1.1. Let f be a real valued function mapping Rn to R. Then f is a DC function if there exist convex functions, g, h : Rn → R such that f can be decomposed as the difference between g and h: f (x) = g(x) − h(x) ∀x ∈ Rn In the remainder of this lecture, we will discuss the solutions to the following - the DC Programming Problem (DCP): minimize n x∈R f0 (x) subject to fi (x) ≤ 0, i = 1, . . . , m. 1 (1) where fi : Rn → R is a differentiable DC function for i = 0, . . . , m. 1.2 Some Intuition About DC Functions Before we continue with the discussion of the solution to (1), we develop some in­ tuition regarding DC functions. Recall that a function f : Rn → R is convex if for every x1 , x2 ∈ Rn and every α ∈ [0, 1] , f (αx1 +(1−α)x2 ) ≤ αf (x1 )+(1−α)f (x2 ). In particular, as you may recall, if f is a twice-differential function, then it is con­ vex if and only if its Hessian matrix is positive-semidefinite. To get a sense of what DC functions can look like, we will look at some common convex functions and the DC functions they can form. Example 1.2. Consider the convex functions f1 (x) = (a) f1 (x) = 1 x (b) f2 (x) = x2 1 x and f2 (x) = x2 . (c) f = 1 x − x2 Example 1.3. Consider the convex functions f1 (x) = abs(x) and f2 (x) = −log(x). (d) f1 (x) = abs(x) (e) f2 (x) = −log(x) (f) f = abs(x) + log(x) Notice that while in these examples, the minimum is easy to find by inspection in the convex functions, it is less clear in the resulting DC function. Clearly, then, having DC functions as part of an optimization problem adds a level of complexity to the problem that we did not encounter in our dealings with convex functions. Fortunately, as we shall soon see, this complexity is not unsurpassable. 2 1.3 1.3.1 How Extensive Are These Functions? Hartman Theorem 1.4. The three following formulations of a DC program are equivalent: 1. sup{f (x) : x ∈ C}, f, C convex 2. inf{g(x) − h(x) : x ∈ Rn }, g, h convex 3. inf{g(x) − h(x) : x ∈ C, f1 (x) − f2 (x) ≤ 0}, g, h, f1 , f2 , C all convex. Proof. • We will show how to go from formulation (1) to formulation (2). Define an indicator function to be: � 0 if x ∈ C IC (x) = (2) ∞ otherwise Then sup{f (x) : x ∈ C} = inf{IC (x) − f (x) : x ∈ Rn }. • We will show how to go from formulation (3) to formulation (1). We have that inf{g(x) − h(x) : x ∈ C, f1 (x) − f2 (x) ≤ 0}, g, h, f1 , f2 , C all convex. Then we can write: αt = inf{g(x) + tmax{f1 (x), f2 (x)} − h(x) − tf2 (x) : x ∈ C} for some value of t such that α = αt' for all tf > t. It can be shown that such a t always exists. • Finally, it is clear that (2) is a special case of (3). Thus, we have shown conversions (1) → (2) → (3) → (1), indicating that the three formulations are equivalent. Next we demonstrate just how large the class of DC functions is. Theorem 1.5 (Hartman). A function f is locally DC if there exists an ε-ball on which it is DC. Every function that is locally DC, is DC. Proposition 1.6. Let fi be DC functions for i = 1, . . . , m. Then the following are also DC: � 1. i λi fi (x), for λi ∈ R 2. maxi fi (x) 3 3. mini fi (x) 4. i fi (x) 5. fi , twice-continuously differentiable 6. If f is DC and g is convex, then the composition (g ◦ f ) is DC. 7. Every continuous function on a convex set, C is the limit of a sequence of uniformly converging DC functions. 2 2.1 Optimality Conditions Duality Before we can discuss the conditions for global and local optimality in the canonical DC Programming problem, we need to introduce some notions regarding duality in the DCP, and develop some intuition of what this gives us. To begin with this, we introduce conjugate functions and use them to demonstrate the relationship between the DCP and its dual. Definition 2.1. Let g : Rn → R. Then the conjugate function of g(x) is g ∗ (y) = sup{xT y − g(x) : x ∈ Rn } To understand what why conjugate functions are so important, we need just one more definition. Definition 2.2. The epigraph of a function, g : Rn → R is the set of points lying on or on top of its graph: epi(g) = {(x, t) ∈ Rn × R : g(x) ≤ t} Note that g is convex if and only if epi(g) is a convex set. Having defined the epigraph, we can now give a geometric interpretation of the conjugate function: The conjugate function g ∗ ‘encloses’ the convex hull of epi(g) with g’s supporting hyperplanes. In particular, we can see that when f is differ­ entiable, df = argsup{xT y − g ∗ (y) : y ∈ Rn } = y(x) dx 4 That is, y is a dual variable for x, and so can be interpreted as (approximately) the gradient (or slope in R2 ) of f at x. This should be a familiar result from the convex analysis lecture. Convex Set C Supporting Hyperplanes Image by MIT OpenCourseWare. Courtesy of Dimitri Bertsekas. Used with permission. (g) A Convex fn f (x) and the affine function (h) A convex set being enclosed by supporting hyper­ meeting f at x. The conjugate, f ∗ (y) of f (x) planes. If C is the epigraph of a function f ∗ , then the is the point of intersection. intersections of C with the hyperplanes is the set of values that f ∗ takes. Theorem 2.3. Let g : Rn → R such that g(x) is lower semi-continuous and convex on Rn . Then g(x) = sup{xT y − g ∗ (y) : y ∈ Rn } where g ∗ (y) is the conjugate of g(x). We provide this theorem without proof, and omit further discussion of lower semicontinuity, but we can safely assume that the functions that we will deal with satisfy this. Note that this condition implies that g ∗∗ = g, that is, the conjugate of g conjugate is g, which means that the definitions of g and g ∗ are symmetric. (What does this mean for minimizing g?) We now demonstrate the relationship between DCP and its dual problem. First note the form of the conjugate function f ∗ for f ≡ g − h: (f (x))∗ = ((g − h)(x))∗ = sup{xT y − (g − h)(x) : x ∈ Rn } = h∗ (y) − g ∗ (y) Let α be the optimum value to DCP. α = inf{g(x) − h(x) : x ∈ X} = inf{g(x) − sup{xT y − h∗ (y) : y ∈ Y } : x ∈ X} = inf{inf{g(x) − xT y + h∗ (y) : x ∈ X} : y ∈ Y } = inf{h∗ (y) − g ∗ (y) : y ∈ Y } 5 Thus, we have that the optimal value to DCP is the same as the optimal value for its dual! Now that is symmetry! This means that we can solve either the primal or the dual problem and obtain the solution to both - the algorithm that we will employ to solve the DCP will crucially rely on this fact. 2.2 Global Optimality Conditions Definition 2.4. Define an ε-subgradient of g at x0 to be ∂ε g(x0 ) = {y ∈ Rn : g(x) − g(x0 ) ≥ (x − x0 )T y − ε ∀x ∈ Rn } Define a differential of g at x0 to be ∂g(x0 ) = ∂ε g(x0 ) ε>0 Given these two definitions, we have the following conditions for global optimality: Theorem 2.5 (Generalized Kuhn-Tucker). . Let x∗ be an optimal solution to the (primal) DCP. Then ∂h(x∗ ) ⊂ ∂g(x∗ ). Let y ∗ be an optimal solution to the dual DCP. Then ∂g ∗ (y ∗ ) ⊂ ∂h∗ (y ∗ ) Proof. This condition essentially follows from the equivalence of the primal and dual optima. We showed before that if α is the optimum value of the DCP, then α = inf{g(x) − h(x) : x ∈ X} = inf{h∗ (y) − g ∗ (y) : y ∈ Y } Then if α is finite, we must have that dom g ⊂ dom h and dom h∗ ⊂ dom g ∗ where dom g = {x ∈ Rn : g(x) < ∞}, the domain of g. That is, h (respectively, g ∗ ) is finite whenever g (respectively, h∗ ) is finite. Note that we require this inclusion because we are minimizing the objective function, and so if there existed an x ∈ R∗ such that g(x) < ∞, h(x) = ∞, then g(x)−h(x) would be minimized at x, yielding an objective value of −∞. Note also that this statement is not an in and only if statement, as we work under the convention that ∞ − ∞ = ∞. Thus, we have that if x∗ is an optimum to the primal DCP, then x∗ ∈ dom g, and by weak duality, g(x∗ ) − h(x∗ ) ≤ h∗ (y) − g ∗ (y), ∀y ∈ dom h∗ 6 and so (for x∗ ∈ dom h) if x∗ ∈ ∂h(x∗ ), then xT y ≥ h(x∗ ) + h∗ (y) ≥ g(x∗ ) + g ∗ (y) where the first inequality is by definition of ∂h(x∗ ) and the second inequality follows from the weak duality inequality presented. We can also think of this in terms of the interpretations of the dual problem as well. We have that if g, h are differentiable, and so ∂h(x∗ ) = ∅ = ∂g(x∗ ), then ∂h(x∗ ) is just the set of gradients of h at x∗ , and by equality of the primal and dual optima, it is the set of y ∗ that optimize the dual problem. Thus, this optimality condition in terms of subdifferentials is analogous to the one we discussed in terms of domains. Corollary 2.6. Let P and D be the solution sets of the primal and dual problems of the DCP, respectively. Then: x∗ ∈ P if and only if ∂ε h(x∗ ) ⊂ ∂ε g(x∗ ) ∀ε > 0. y ∗ ∈ D if and only if ∂ε g ∗ (y ∗ ) ⊂ ∂ε h∗ (y ∗ ) ∀ε > 0. Theorem 2.7. Let P and D be the solution sets of the primal and dual problems of the DCP, respectively. Then: {∂h(x) : x ∈ P} ⊂ D ⊂ dom h∗ and {∂g ∗ (y) : y ∈ D} ⊂ P ⊂ dom g Note that this theorem implies that solving the primal DCP implies solving the dual DCP. 2.3 Local Optimality Conditions We would like to construct an algorithm to find global optimal solutions based on the conditions discussed in the previous section. However, finding an algorithm that does this efficiently in general is an open problem, and most approaches are combinatorial, rather than convex-based, and so rely heavily on the formulation of a given problem, and are often inefficient. Thus, we present local optimality conditions, which (unlike the global optimality conditions) can be used to create a convex-based approach to local optimization. We present these theorems with­ out proof as, although they are crucially important to solving DC programming 7 problems, their proofs do not add much more insight. Thus, we refer further in­ vestigation to either Hurst and Thoai or Tao and An. Theorem 2.8 (Sufficient Local Optimality Condition 1). Let x∗ be a point that admits a neighborhood U (x) such that ∂h(x) ∩ ∂g(x∗ ) 6= ∅ ∀x ∈ U (x) ∩ dom g Then x∗ is a local minimizer of g − h. Theorem 2.9 (Sufficient Local Optimality Condition 2: Strict Local Optimal­ ity). Let int(S) refer to the interior of set S. Then if x∗ ∈ int(dom h) and ∂h(x∗ ) ⊂ int(∂g(x∗ )), then x∗ is a strict local minimizer of g − h. Theorem 2.10 (DC Duality Transportation of a Local Minimizer). Let x∗ ∈ dom ∂h be a local minimizer of g − h and let y ∗ ∈ ∂h(x∗ ). Then if g ∗ is differ­ entiable at y ∗ , y ∗ is a local minimizer of h∗ − g ∗ . More generally, if y ∗ satisfies Theorem 2.9, then y ∗ is a local minimizer of h∗ − g ∗ . 3 Algorithms As discussed earlier, the conditions for global optimality in DC programs do not wield efficient general algorithms. Thus, while there are a number of popular tech­ niques - among them, branch-and-bound and cutting planes algorithms, we omit discussion of them and instead focus on the convex-based approach to local opti­ mization. In fact, although there has not been an analytic result to justify this, according to the DC Programming literature, the local optimization approach of­ ten yields the global optimum, and a number of regularization and starting-point choosing methods exist to assist with incorporating the following local optimiza­ tion algorithm to find the global optimum in different cases. 3.1 DCA-Convex Approach to Local Optimization We now present an algorithm to find local optima for a general DC program. First, we offer the algorithm in raw form, then we will explain each step in an iteration, and finally, we will state a few results regarding the effectiveness and efficiency of the algorithm. 8 3.1.1 DCA 1. Choose x0 ∈ dom g 2. for k ∈ N do: 3. choose yk ∈ ∂h(xk ) 4. choose xk+1 ∈ ∂g ∗ (yk ) 5. if min{|(xk+1 − xl )i |, | 6. then return xk+1 7. (xk+1 −xl )i |} (xk )i ≤ δ: end if 8. end for 3.1.2 DCA Explanation and Intuition Let us go over each step of the DCA algorithm in more detail. The overarching method of the algorithm is to create two sequences of variables, {xk }k , {yk }k so that {xk } converges to a local optimum of the primal problem, x∗ , and {yk } converges to the local optimum of the dual problem, y ∗ . The key idea is to manipulate the symmetry of the primal and dual problem in order to follow a variation on the typical sub gradient-descent method used in convex optimization. Let us now consider each step individually. • Choose x0 ∈ dom g: Since we are utilizing a descent approach, the convergence of the algorithm is independent of the starting point of the sequences that the algorithm creates. Thus, we can instantiate the algorithm with an arbitrary choice of x0 so long as it is feasible. • Choose yk ∈ ∂h(xk ): We have that ∂h(xk ) = arg min{h∗ (y) − g ∗ (yk−1 ) − xTk (y − yk−1 ) : y ∈ Rn }. Moreover, since this is a minimization over y, we hold yk−1 , xk constant, and so we have that ∂h(xk ) = arg max{xTk y − h∗ (y) : y ∈ Rn }. Computing this, however is just an exercise in convex optimization, since by the local optimality conditions, xk is (approximately) a subgradient of h∗ − g ∗ , and so serves a role similar to that in a typical subgradient descent algorithm, which we solve quickly and efficiently. Since we maximize given xk (which improves together with yk−1 ), we guarantee that (h∗ − g ∗ )(yk − yk−1 ) ≤ 0 ∀k ∈ N, and 9 due to the symmetry of duality, that yk converges to a critical point of (h∗ − g ∗ ), i.e., a local minimizer. • Choose xk+1 ∈ ∂g ∗ (yk ): Given the symmetry of the primal and dual problem, this is entirely symmetric to the step finding yk . Thus, we choose xk+1 ∈ arg max{xT yk − g(x) : x ∈ Rn }. −xl )i • If min{|(xk+1 − xl )i |, | (xk+1 |} ≤ δ: then return xk+1 : (xk )i Although we can guarantee convergence in the infinite limit of k, complete convergence may take a long time, and so we approximate an optimal solu­ tion within a predetermined bound, δ. Once the solution (the change in xk or yk ) is small enough, we terminate the algorithm and return the optimal value for xk+1 - recall that solving for xk is equivalent to solving for yk . Figure 1: An example for the Proximal Point subgradient descent method. This has been shown to be equivalent to a regularized version of DCA, and so offers valuable intuition into how DCA works. Courtesy of Dimitri Bertsekas. Used with permission. 3.1.3 Well-Definition and Convergence Here we present a few results regarding the effectiveness and efficiency of the DCA algorithm. We will give the results without proof, and direct anyone interested in further delving into this matter to Thoai. Lemma 3.1. The sequences {xk } and {yk } are well defined if and only if dom ∂g ⊂ dom ∂h dom ∂h∗ ⊂ dom ∂g ∗ Lemma 3.2. Let h be a lower semi-continuous function on Rn and {xk } be a sequence of elements in Rn such that (i) xk → x∗ ; (ii) There exists a bounded 10 sequence {yk } such that yk ∈ ∂h(xk ); (iii) ∂h(x∗ ) 6= ∅. Then lim h(xk ) = h(x∗ ) k→∞ 4 Applications to Machine Learning The literature has references to many uses of DC Programming in Operations Research, Machine Learning, and Economics. One interesting use of DC Programming is discussed in the 2006-paper “A DCProgramming Algorithm for Kernel Selection”. In this paper, the authors discuss a greedy algorithm to learn a kernel from a convex hull of basic kernels. While this approach had been popularized before, it was limited to a finite set of basic kernels. The authors comment that the limitation was due to the non-convexity of a critical maximization involved in conducting the learning, but find that the optimization problem can be formulated as a DC Program. In particular, the objective function used to weigh basic kernels is DC as the limit of DC functions. Another interesting use of DC Programming is discussed in the 2008-paper “A DC programming approach for feature selection in support vector machines learn­ ing”. Here, DC Programming is employed in an SVM algorithm that attempts to choose optimally representative features in data while constructing an SVM classifier simultaneously. The authors equate this problem to minimizing a zeronorm function over step-k feature vectors. Using the DC-decomposition displayed below, the authors employ the DCA algorithm to find local minima and applied it to ten datasets - some of which were particularly sparse - to find that the DCA algorithm created consistently good classifiers that often had the highest correct­ ness rate among the tested classifiers (including standard SVMs, for example). Despite this, DCA consistently used less features than the standard SVM and other classifiers, and as a result, was more efficient and used less CPU capacity than many of the other commonplace classifiers. Thus, all in all, DCA proved a greatly attractive algorithm for classifying data, and particularly excellent for very large and sparse datasets. 11 5 Conclusion As we can see, DC Programming is quite young. Although the oldest result presented in this lecture dates back to the 1950s (Hartman’s), many of the results and algorithms discussed here were developed in the 1990s, and their applications are still very new to the scientific community. However, as can be seen from the examples that we presented, DC Programming has tremendous potential to expand and expedite many of the algorithms and techniques that are central to Machine Learning, as well as other fields. Thus, it is likely that the near future will bring many more algorithms inspired by DCA and DCA-related approaches as well as the various combinatoric global approaches that are used to solve DC Programming problems, but were not discussed here. References [1] H. A. L. T. H. M. L. V. V. Nguyen and T. P. Dinh. A dc programming approach for feature selection in support vector machines learning. Advances in Data Analysis and Classification, 2(3):259–278, 2008. [2] A. A. R. H. C. A. M. M. Pontil. A dc-programming algorithm for kernel selec­ tion. Proceedings of the 23rd International Conference on Machine Learning, 2006. [3] P. D. TAO and L. T. H. AN. Convex analysis approach to d. c. programming: Theory, algorithms and applications. ACTA MATHEMATICA VIETNAM­ ICA, 22(1):289–355, 1997. [4] R. H. N. Thoai. Dc programming: An overview. Journal of Optimization Theory and Application, 193(1):1–43, October 1999. 12