Optimization Overview

In this post, I am summarizing all I know for optimization. It is a hard, primary problem once you define your model with some kind of objective function integrating a bunch of parameters and you start to wonder how you could learn the values of those parameters  to minimize (or maximize) the objective function.

Unconstrained Optimization

You wish to learn parameters to minimize the objective function and those parameters don’t have any constraints. If this is a convex problem, then there must be a global minimum and you can find it by finding all points whose first-order derivatives are zero. If the objective function is not convex, don’t panic. Gradient Descent (with lots of its variants such as momentum gradient descent, stochastic gradient descent, AdaGra) definitely help you to reach local minimum. People find that when the number of variables is large you will reach a not too bad local minimum. (missing reference here). Another family of approach is called quasi-Newton method which is my favorite. They have comparably fast speed as Newton Method does but its computation is way more economic. I implemented one unconstrained LBFGS algorithm in Python and some of its details are revealed here: https://czxttkl.com/?p=1257 

Constrained Optimization

When it comes to constrained optimization, there exist tons of approaches which can be applied in different conditions. Here is a summary of techniques you can equip: (from https://ipvs.informatik.uni-stuttgart.de/mlr/marc/teaching/13-Optimization/03-constrainedOpt.pdf):

Screenshot from 2016-02-22 12:53:14

First of all, let’s get familiar with the foundations of constrained optimization.

Suppose we are faced a constrained optimization problem:

\min\limits_{x \in \mathbb{R}^n} f(x)

s.t. c_i(x)=0, i \in \xi and g_i(x) \geq 0, i \in \iota

At a feasible point x, the inequality constraint i \in \iota is said to be active if g_i(x)=0 and inactive if the strict inequality g_i(x) > 0 is satisfied.

The active set \mathcal{A}(x) at any feasible x consists of the equality constraint indices from \xi together with the indices of the inequality constraints i for which g_i(x)=0, i.e., \mathcal{A}(x)=\{i \in \xi\} \cup \{i \in \iota|g_i(x)=0\}

LICQ: Given the point x and the active set \mathcal{A}(x), we say that the linear independence constraint qualification (LICQ) holds if the set of the active constraint gradient \{\nabla g_i(x), i \in \mathcal{A}(x)\} is linearly independent.

KTT (Karush-Kuhn-Tucker) conditions

KTT conditions are necessary conditions for a optimum solution in nonlinear programming. (Some of the constraints or the objective function is non-linear.) In other words, if a solution x^* is a local minimum of a function f which we aims to minimize under certain equality constraints \{c_i(x) | i \in \xi \} and inequality constraints \{g_i(x) | i \in \iota\}, then x^*, f(x^*), \{h_i(x^*)\} and \{c_i(x^*)\} must satisfy KTT conditions. The reverse, i.e., if a solution x satisfies KTT conditions then it is a local minimum, is not always true.

More formally to state KTT conditions:

Assume x^* is a local solution which achieves local minimum of f(x). If f(x), \{c_i(x) | i \in \xi\} and \{g_i(x) | i \in \iota \} are continuously differentiable, and LICQ holds at x^*, then there is a Lagrange multiplier vector \lambda^*, with components \lambda_i^*, i \in \xi \cup \iota, such that following conditions are satisfied at (x^*, \lambda^*):

  • Define: \mathcal{L}(x, \lambda) = f(x) - \sum\limits_{i \in \xi} \lambda_i c_i(x) - \sum\limits_{i \in \iota} \lambda_i g_i(x)
  • \nabla \mathcal{L}(x^*, \lambda^*)=\nabla f(x^*) - \sum\limits_{i \in \xi} \lambda_i^* \nabla c_i(x^*) - \sum\limits_{i \in \iota | g_i(x^*)=0} \lambda_i^* \nabla g_i(x^*)= 0
  • c_i(x^*) = 0, i \in \xi
  • g_i(x^*) \geq 0, i \in \iota
  • \lambda_i^* \geq 0, i \in \iota
  • \lambda_i^* c_i(x^*)=0, i \in \xi
  • \lambda_i^* g_i(x^*)=0, i \in \iota

In \nabla \mathcal{L}(x^*, \lambda^*), we ignore to subtract \sum\limits_{i \in \xi \cup \iota \backslash \mathcal{A}(x^*)} \lambda_i^* \nabla c_i(x^*) because KKT conditions indicate that\lambda_i=0 for i \in \xi \cup \iota \backslash \mathcal{A}(x^*), i.e., for those inactive inequality constraints g_i(x^*) > 0.

Note that, we require LICQ holds at x^* so that \lambda^* is unique. LICQ is called constraint qualification. There are other constraint qualifications that must hold to guarantee KKT conditions hold if x is a local minimum. (See https://www.math.washington.edu/~burke/crs/516/notes/cq_lec.pdf) For example, if x^* is a local minimum and Abadie constraint qualification (ACQ) holds at x^*, then KTT conditions still hold. Here ACQ is a weaker constraint qualification than LICQ because LICQ => ACQ.

Here is my rough understanding of KKT conditions: We want to let argmin_x \mathcal{L}(x, \lambda) = argmin_x f(x) subject to the same constraints. \mathcal{L}(x, \lambda) consists of f(x) and a bunch of -\lambda_i g_i(x) and -\lambda_i c_i(x). To make sure \mathcal{L}(x, \lambda) and f(x) have the same local solution x^*, you must do the following checks:

  • For a local minimum x^* of f(x), if \lambda_i^* g_i(x^*) > 0, then argmin_x \mathcal{L}(x, \lambda) = argmin_x f(x) is not guaranteed unless \lambda_i^* g_i(x^*) = 0. Put it in other words, if you want \mathcal{L}(x, \lambda) to have the same optimum as f(x), then \lambda_i g_i(x), the term added to\mathcal{L}(x, \lambda) besides f(x), should not influence f(x) and \mathcal{L}(x, \lambda) to reach their local minimum at x^* at the same time.
  • For equality constraints and active inequality constraints,  \lambda_i^* c_i(x^*) and \lambda_i^* g_i(x^*) are zeros which can be safely added to \mathcal{L}(x^*, \lambda^*) without affecting x^* being the local optimum for both \mathcal{L}(x, \lambda) and f(x).

Next, you should note that for the local minimum x^*,\nabla f(x^*) - \sum\limits_{i \in \xi} \lambda_i^* \nabla c_i(x^*) - \sum\limits_{i \in \iota | g_i(x^*)=0} \lambda_i^* \nabla g_i(x^*)= 0. This means, at x^*, the force to minimize f(x) is pulled by the force to ensure x^* to stay in feasible areas. This can be illustrated in the following example:

Screenshot from 2016-02-23 21:58:35

As you can see, at x^*=(1, 0), the force to minimize f(x) is the direction -\nabla f(x^*) = (-1, -\frac{1}{2}) (marked as purple arrow) . The force to pull x^* to stay the feasible region (the square in the middle) is \nabla c_1(x^*) and \nabla c_2(x^*) (marked as red arrows). \lambda^*=(\frac{3}{4}, \frac{1}{4}, 0, 0) makes sure \nabla f(x^*), \nabla c_1(x^*) and \nabla c_2(x^*) are linearly dependent at this point so that x* can’t be dragged any further by any force.
Untitled Diagram

Although knowing a feasible point x and a \lambda satisfying KTT conditions is generally not sufficient to conclude x is the solution for the constrained optimization problem, under certain conditions x indeed is the solution. For example, if you know additional information like Second Order Sufficient Conditions, then KTT conditions become sufficient to say x is a local solution. For another example, if Slater’s condition holds for the constrained problem, then (x, \lambda) \text{ satisfy KTT conditions } \Rightarrow x \text{ is a local solution}. (Ref: http://math.stackexchange.com/questions/379543/kkt-and-slaters-condition, https://www.cs.cmu.edu/~ggordon/10725-F12/slides/16-kkt.pdf p11). In the famous SVM model, the objective function and its constraints also satisfies Slater’s condition such that KTT conditions are sufficient and necessary. See https://czxttkl.com/?p=3114 for how we do optimization for SVM using KTT.

Actually, Slater’s condition is originally used to prove that the strong duality holds for the constrained problem. Let’s illustrate it by introducing Duality Theory.

Duality Theory

An optimization problem can usually be constructed as a primal problem, for which there exists a corresponding dual problem. The dual problem sometimes has magics such that the computation of the dual problem is easier, or the solution of the dual problem can hint bounds of the primal problem.

Let’s say we face the same optimization problem as follows.

\min\limits_{x \in \mathbb{R}^n} f(x)

s.t. c_i(x)=0, i \in \xi and g_i(x) \geq 0, i \in \iota

This problem can be converted into a primal problem incorporating all the constraints c_i(x)=0, i \in \xi and g_i(x) \geq 0, i \in \iota:

\min\limits_x p(x) = \min\limits_x \max\limits_{\lambda} \mathcal{L}(x, \lambda) = \min\limits_x \max\limits_{\lambda} f(x) - \sum\limits_{i \in \xi} \lambda_i c_i(x) - \sum\limits_{i \in \iota} g_i(x) subject to \lambda \geq 0

Why does the primal problem look like this? This is because the minimum of p(x) is the same as the minimum of the original constrained optimization problem. For those x s.t. c_i(x) \neq 0 \text{ or } g_i(x) <0,  \max\limits_{\lambda} \mathcal{L}(x, \lambda) = \infty therefore any x violating the constraints is definitely not the solution to the primal problem. For those x satisfying the constraints, \max\limits_{\lambda} \mathcal{L}(x, \lambda) =f(x). Therefore, it is equivalent to \min\limits_x p(x) as to do the original constrained problem.

Now we construct the dual problem of the primal problem. The dual problem is:

\max\limits_{\lambda \in \mathbb{R}^{|\xi \cup \iota|}} q(\lambda) = \max\limits_{\lambda} \inf\limits_x \mathcal{L}(x, \lambda) subject to \lambda \geq 0. (Here \inf h(x) is infimum function denoting the lower bound of h(x). For example, \inf\limits_x e^{-x} = 0. We don’t use \min\limits_x h(x) to denote the lower bound of h(x) because h(x) can be infinitely close to but never reach the lower bound for any x.)

Now we have Weak Duality which goes like this:

\max\limits_{\lambda} q(\lambda) = \max\limits_{\lambda} \inf\limits_x \mathcal{L}(x, \lambda) \leq \min\limits_x p(x) = \min\limits_x \max\limits_{\lambda} \mathcal{L} (x, \lambda) for any \lambda \geq 0 and any feasible x

This suggests that the optimum of the dual problem has some non-negative duality gap smaller than the optimum of the primal problem. If we can close the duality gap between the dual and the primal, we can solve the primal problem by solving the dual problem which sometimes is much easier. The duality gap is zero when Strong Duality holds for the optimization problem.

Ref: http://jackvalmadre.tumblr.com/post/15798388985

Here for the duality theory for linear programming specifically: https://www.math.washington.edu/~burke/crs/407/notes/section4.pdf

Lagrange Multiplier

Lagrange Multiplier is used for finding minimum/maximum of a function under equality constraints. Lagrange Multiplier utilizes an important fact that when a function achieves a minimum/maximum under certain equality constraints, the gradient of the function and the gradients of the constraints have the same direction. What Lagrange Multiplier essentially does is to add constraints, each multiplied by a Lagrange Multiplier, to the objective function to form an “auxiliary objective function”. By taking the derivative of the auxiliary objective function and set it to zero, we can find all stationary points of the auxiliary objective function. Among those stationary points, there must be at least one that achieves the minimum/maximum. A good tutorial of Lagrange Multiplier can be found here: http://www.slimy.com/~steuard/teaching/tutorials/Lagrange.html

Lagrange Multiplier method can be modified to solve constrained problems with inequality constraints. The multiplier and optimal solution can be found by using KKT conditions: http://www.csc.kth.se/utbildning/kth/kurser/DD3364/Lectures/KKT.pdf. One example of using Lagrange Multiplier to solve problems with inequality constraints can be found here: http://users.wpi.edu/~pwdavis/Courses/MA1024B10/1024_Lagrange_multipliers.pdf.

A family of constrained problem solvers converts constrained problems into unconstrained problems. In the family of such solvers, one kind of algorithms are called inexact method: the solution to the transformed unconstrained problem is close but not exactly the same to that to the original constrained problem. The other kind of algorithm is called exact method in which the solution does not change after the transformation. In the following introduced algorithms, Quadratic Penalty Method belongs to inexact method. \ell_1 Penalty Method and Augmented Lagrangian Method belong to exact method.

Quadratic Penalty Method


Screenshot from 2016-06-29 18:05:06

quamethod\ell_1 Penalty Method

When only equality constraints are present, Q(x;\mu) is smooth. But when inequality constraints are also present, Q(x;\mu) is not smooth any more. Moreover, \mu sometimes goes to be very large in iterations, likely to cause the Hessian matrix of Q(x;\mu), \nabla_{xx}^2 Q(x;\mu), become ill-conditioned.

\ell_1 Penalty Method

Screenshot from 2016-06-29 18:09:49

Screenshot from 2016-06-29 18:08:41

\ell_1 Penalty Method is exact but the non-smoothness of \phi_1(x;\mu) makes it difficult to find minimizer of it. Refer to the textbook to some hack ways to solve it.

Augmented Lagrangian Method

It reduces the possibility of ill-conditioning of Quadratic Penalty Method and preserves the smoothness as opposed to \ell_1 Penalty Method.

Screenshot from 2016-06-29 18:20:21

https://en.wikipedia.org/wiki/Augmented_Lagrangian_method

Stochastic Gradient Descent

Stochastic Gradient Descent has been popular among large-scale machine learning. So far the most canonical text about stochastic gradient descent I’ve found is: Bottou, L. (2012). Stochastic gradient descent tricks. In Neural Networks: Tricks of the Trade (pp. 421-436). Springer Berlin Heidelberg.

The core idea is that “use stochastic gradient descent when training time is the bottleneck“. What does this mean? Roughly understanding, gradient descent should result in higher training accuracy than its stochastic counterpart because the former uses the full dataset to calculate a less noisy weight update. However in an iteration, gradient descent will take a long time to traverse all samples to calculate weight update. During that, the weight update of stochastic gradient descent has happened many times and the objective function it optimizes has been improved many times! Therefore, the paper points out that “SGD needs less time than the other algorithms to reach a predefined expected risk”. Please understand Table 2 in the paper carefully.

Projected Gradient Descent

This post highly summarizes the difference between ordinary gradient descent and projected gradient descent: http://math.stackexchange.com/questions/571068/what-is-the-difference-between-projected-gradient-descent-and-ordinary-gradient

This paper showcases how to use projected gradient descent for non-negative matrix factorization:

Lin, C. J. (2007). Projected gradient methods for nonnegative matrix factorization. Neural computation, 19(10), 2756-2779.

This paper also points out that if all the constraints are just non-negativity constraint, then we can set parameters P=P'^2, then we can just optimize the objective w.r.t P', where P' \in \mathbb{R}.

However, I have not seen projected stochastic gradient descent been applied on any problem. I will keep searching why it is rarely used. My guess is that  stochastic-based updates are already very noisy and unstable. The projection performed after each stochastic update will make parameter updates even more unstable such that the optimization cannot find the correct path towards the optimal parameterization.

Quadratic Programming

Quadratic Programming is an optimization problem with a quadratic objective function and linear constraints. If constraints only consist of equalities, it looks like the following:

\min_x q(x) = \frac{1}{2} x^T G x + x^T c \\ \text{subject to } Ax=b, A \in \mathbb{R}^{m \times n}, m < n

Now we introduce an iterative way to solve the quadratic programming with equality constraints, which is covered in Chapter 16.3 of Numerical Optimization by Nocedal & Wright.

First, the optimal solution x^* \in \mathbb{R}^n can be represented as x^*=Yx_Y+Zx_Z, where x_Z \in \mathbb{R}^{n-m}, x_Y \in \mathbb{R}^{m}, Y \in \mathbb{R}^{n \times m}, Z \in \mathbb{R}^{n \times (n-m)} and AZ=0. Note that in x^*=Yx_Y+Zx_Z, Y, Z and x_Y are some constants determined by elimination algorithms (introduced soon). x_Z is the free variables that \in \mathbb{R}^{n-m}. Why can x be represented by x_Z? This is because A has rank m (m<n) therefore we can use n-m variables in x_Z to represent n variables in x. Since AZ=0, Ax=b indicates that AYx_Y=b. Therefore, Yx_Y is a solution satisfying Ax=b. Therefore, x^*=Yx_Y+Zx_Z can be seen as a sum of a particular solution to the constraint Ax=b (the first term) plus a transformation Z applied on x_Z. Since AZ=0 holds, any choice of x_Z will still satisfy the constraint Ax=b. The only difference of choices of x_Z remains to be whether x_Z minimizes the objective function.

To determine Y, Z and x_Y, we use methods called simple elimination or general elimination:

1. simple elimination says that: find a permutation matrix P such that AP= [B|N]. Here the columns of B \in \mathbb{R}^{m \times m} are a subset of columns of A that are linearly independent. N \in \mathbb{R}^{m \times (n-m)}. Now that if we set:Screenshot from 2016-06-13 14:25:26

, x can be represented as x = Yb + Zx_Z.

2. general elimination uses QR-factorization on A^T\Pi (\Pi is a permutation matrix to obtain Q_1, Q_2 and R:

Screenshot from 2016-06-13 14:21:25

Then, set x to Screenshot from 2016-06-13 14:24:24

(I am unclear of when to use simple elimination or general elimination as the book compares the two algorithms vaguely to me. )

Now let’s go back solving the quadratic programming problem once obtaining x^*=Yx_Y+Zx_Z. Plugging it to \min_z q(x)=\frac{1}{2} x^T G x + x^T c with all constants removed, we obtain min_{X_z} \frac{1}{2}x_Z^T Z^T G Z x_Z + x_Z^T (Z^T G Y x_Y + Z^T c). The textbook gives an iterative method to find desired x_Z using Conjugate Gradient method, with the assumption that Z^T G Z is positive definite:

Screenshot from 2016-06-13 15:00:54Screenshot from 2016-06-13 15:01:39

\ell_1 regularization

Please note that \ell_1 regularization refers to the problem of the form: \min\limits_x f(x) = L(x) + \lambda \left\Vert x\right\Vert_1. Such forms of problems often have solutions x^* with great sparsity. That is why \ell_1 regularization is often used to achieve feature selection. This kind of problem is different to another kind of problem (also introduced in this post): use \ell_1 as nonsmooth exact penalty function to solve constrained problems.

This paper gives a full overview of methods to solve \ell_1 regularization problems.  I am here giving a simple approach:

When given a \ell_1 regularization problem \min\limits_{x \in \mathbb{R}^n} f(x) = L(x) + \lambda \left\Vert x\right\Vert_1, we can write an equivalent constrained optimization problem: \min\limits_{x^+,x^-} L(x^+ - x^-) + \lambda(x^+ + x^-)\quad s.t. \quad x^+ = \max(0,x), x^- = -\min(0,x), x^+ \geq 0, x^- \geq 0

To explain why we can write such equivalent constrained optimization problem, let’s first check that x^+ - x^- = x and x^+ + x^- = |x|. Therefore L(x^+ - x^-)=L(x) and \lambda(x^+ + x^-) = \lambda\left\Vert x\right\Vert_1. Depending on the strength of \lambda, x^+ and x^- will be constrained to not fluctuate too radically, resulting the constrained \ell_1 norm of original x. Eventually, the \ell_1 regularization problem becomes a box constraint optimization problem (albeit the box is single side bounded, because you only require x^+ and x^- to be non-negative.)

There is an online implementation of such method using L-BFGS-B (box constraint of L-BFGS): https://gist.github.com/vene/fab06038a00309569cdd

Coordinate Descent

https://www.cs.cmu.edu/~ggordon/10725-F12/slides/25-coord-desc.pdf

Some Topics to Explore

Non convex global optimization: http://mathoverflow.net/questions/32533/is-all-non-convex-optimization-heuristic

Leave a comment

Your email address will not be published. Required fields are marked *