QP as LP: cutting planes

In (1) I described a simple linearization scheme for a QP model we we can solve it as a straight LP. For a simple (convex) quadratic function \(f(x)\), instead of solving:

\[ \min\> f(x)=x^2\]

we solve

\[\begin{align} \min\>&z\\                               &z \ge a_i x + b_i\end{align}\]

In this post I do things slightly different: instead of adding the linear inequalities ahead of time, we add them one at the time based on the previously found optimal point. This approach is called a cutting plane technique (2).

Example: portfolio model

We consider again the simple portfolio model:

\[\bbox[lightcyan,10px,border:3px solid darkblue]{ \begin{align} \min \>& \sum_t w_t^2\\                               & w_t = \sum_i r’_{i,t} x_i\\                               & \sum_i \mu_i x_i \ge R\\                               & \sum_i x_i = 1\\                               &x_i\ge 0\end{align}} \]

The model is linearized as an LP model:

\[\bbox[lightcyan,10px,border:3px solid darkblue] {\begin{align}\min \>& \sum_t z_t\\                               &z_t \ge a_{t,k} w_t + b_{t,k} & \text{Added $K\times T$ cuts}\\                               & w_t = \sum_i r’_{i,t} x_i\\                               & \sum_i \mu_i x_i \ge R\\                               & \sum_i x_i = 1\\                               &x_i, z_t\ge 0\end{align}} \]

Initially we start without cuts. Later on, during the Cutting Plane algorithm, we will add linear cuts in each iteration. The algorithm would look like:

1. \(k:=1\)
2. Solve model LP, let \(w^*_t\) be the optimal values for \(w_t\).
3. if \(k=\)MAXIT: STOP
4. Add the constraint \(z_t \ge a_{t,k} w_t + b_{t,k}\) where \(a_{t,k}=2 w^*_t\) and \(b_{t,k}=-(w^*_t)^2\). Note that we add one cut for each \(w_t\) here (our dataset has \(T=717\) time periods). 
5. \(k:=k+1\)
6. go to step 2

Here we stop simply when a certain number of iterations MAXIT has been exceeded. That can be refined by stopping when the objective does not seem to change much anymore. Another optimization could be to only add cuts that are different from the ones added before (for some \(t\) we may converge quicker than for others).

The algorithm converges very fast:

----    118 PARAMETER report2             obj(LP)         w'w iter1               1.02907546 iter2   0.21577087  0.46879016 iter3   0.33210537  0.48432099 iter4   0.37143835  0.41397603 iter5   0.40990988  0.41362293 iter6   0.41109694  0.41222051 iter7   0.41183602  0.41212331 iter8   0.41192720  0.41204267 iter9   0.41196991  0.41204112 iter10  0.41197620  0.41203782

Note that the optimal QP solution has an objective: 0.41202816.

This is pretty good performance especially because the dataset is not very small (it has 717 time periods and 83 stocks).

Here is a picture of the cuts introduced for the first element \(z_1=w_1^2\):

A combined approach

We can even combine the two methods:

1. Start with a coarse-grained (i.e. cheap) initial set of cuts. In (1) we used \(n=10\) inequalities per \(w_t\). For this experiment I reduced this to \(n=5\).
2. Then apply our cutting plane algorithm. Instead of MAXIT=10 iterations we now do 5 iterations.

This yields even faster convergence:

----    142 PARAMETER report3            obj(LP)         w'w iter1  0.32921509  0.41401219 iter2  0.40812966  0.41471624 iter3  0.41017228  0.41211351 iter4  0.41188869  0.41208158 iter5  0.41195624  0.41203588

References

1. QP as LP: piecewise linear functions, http://yetanothermathprogrammingconsultant.blogspot.com/2017/04/qp-as-lp-piecewise-linear-functions.html
2. J. E. Kelley, Jr. "The Cutting-Plane Method for Solving Convex Programs." J. Soc. Indust. Appl. Math. Vol. 8, No. 4, pp. 703-712, December, 1960.
3. Cardinality Constrained Portfolio Optimization: MIQP without MIQP Solver, http://yetanothermathprogrammingconsultant.blogspot.com/2016/02/cardinality-constrained-portfolio.html. Here this cutting plane method is applied to a MIQP model (not strictly “allowed” as this is no longer convex, but useful as a heuristic).