Quadratic Programming
In this chapter we turn the LP problem from Chapter 10 into a Quadratic Programming (QP) problem, and the first MIP model from Chapter 11 into a Mixed Integer Quadratic Programming (MIQP) problem. The chapter shows how to
- define quadratic objective functions,
- incrementally define and solve problems.
Chapter 7 shows how to formulate and solve the same examples with Mosel and in Chapter 17 the QP problem is input and solved directly with Xpress-Optimizer.
Problem description
The investor may also look at his portfolio selection problem from a different angle: instead of maximizing the estimated return and limiting the portion of high-risk investments he now wishes to minimize the risk whilst obtaining a certain target yield. He adopts the Markowitz idea of getting estimates of the variance/covariance matrix of estimated returns on the securities. (For example, hardware and software company worths tend to move together, but are oppositely correlated with the success of theatrical production, as people go to the theater more when they have become bored with playing with their new computers and computer games.) The return on theatrical productions are highly variable, whereas the treasury bill yield is certain.
Question 1: Which investment strategy should the investor adopt to minimize the variance subject to getting some specified minimum target yield?
Question 2: Which is the least variance investment strategy if the investor wants to choose at most four different securities (again subject to getting some specified minimum target yield)?
The first question leads us to a Quadratic Programming problem, that is, a Mathematical Programming problem with a quadratic objective function and linear constraints. The second question necessitates the introduction of discrete variables to count the number of securities, and so we obtain a Mixed Integer Quadratic Programming problem. The two cases will be discussed separately in the following two sections.
QP
To adapt the model developed in Chapter 2 to the new way of looking at the problem, we need to make the following changes:
- New objective function: mean variance instead of total return.
- The risk-related constraint disappears.
- Addition of a new constraint: target yield.
The new objective function is the mean variance of the portfolio, namely:
s,t 
SHARES VARst·fracs ·fract
where VARst is the variance/covariance matrix of all shares. This is a quadratic objective function (an objective function becomes quadratic either when a variable is squared, e.g., frac12, or when two variables are multiplied together, e.g., frac1 · frac2).
The target yield constraint can be written as follows:
s SHARES RETs·fracs 
TARGET
The limit on the North-American shares as well as the requirement to spend all the money, and the upper bounds on the fraction invested into every share are retained. We therefore obtain the following complete mathematical model formulation:
s,t SHARES VARst·fracs ·fracts NA fracs 0.5s SHARES fracs = 1s SHARES RETs·fracs TARGET
s 
SHARES: 0 
fracs 0.3
Implementation with BCL
The estimated returns and the variance/covariance matrix are given in the data file foliocppqp.dat:
! trs haw thr tel brw hgw car bnk sof elc 0.1 0 0 0 0 0 0 0 0 0 ! treasury 0 19 -2 4 1 1 1 0.5 10 5 ! hardware 0 -2 28 1 2 1 1 0 -2 -1 ! theater 0 4 1 22 0 1 2 0 3 4 ! telecom 0 1 2 0 4 -1.5 -2 -1 1 1 ! brewery 0 1 1 1 -1.5 3.5 2 0.5 1 1.5 ! highways 0 1 1 2 -2 2 5 0.5 1 2.5 ! cars 0 0.5 0 0 -1 0.5 0.5 1 0.5 0.5 ! bank 0 10 -2 3 1 1 1 0.5 25 8 ! software 0 5 -1 4 1 1.5 2.5 0.5 8 16 ! electronics
We may read this datafile with the function XPRBreadarrlinecb: all comments preceded by ! and also empty lines are skipped. We read an entire line at once indicating the format of an entry (`g') and the separator (any number of spaces or tabulations).
For the definition of the objective function we now use a quadratic expression (type XPRBquadExp). Since we now wish to minimize the problem, the optimization is started with the method minim (with empty string argument indicating the default algorithm).
#include <iostream>
#include "xprb_cpp.h"
using namespace std;
using namespace ::dashoptimization;
#define DATAFILE "foliocppqp.dat"
#define TARGET 9 // Target yield
#define MAXNUM 4 // Max. number of different assets
#define NSHARES 10 // Number of shares
#define NNA 4 // Number of North-American shares
double RET[] = {5,17,26,12,8,9,7,6,31,21}; // Estimated return in investment
int NA[] = {0,1,2,3}; // Shares issued in N.-America
double VAR[NSHARES][NSHARES]; // Variance/covariance matrix of
// estimated returns
int main(int argc, char **argv)
{
int s,t;
XPRBprob p("FolioQP"); // Initialize a new problem in BCL
XPRBlinExp Na,Return,Cap,Num;
XPRBquadExp Variance;
XPRBvar frac[NSHARES]; // Fraction of capital used per share
FILE *datafile;
// Read `VAR' data from file
datafile=fopen(DATAFILE,"r");
for(s=0;s<NSHARES;s++)
XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g ", VAR[s], NSHARES);
fclose(datafile);
// Create the decision variables
for(s=0;s<NSHARES;s++)
frac[s] = p.newVar(XPRBnewname("frac(%d)",s+1), XPRB_PL, 0, 0.3);
// Objective: mean variance
for(s=0;s<NSHARES;s++)
for(t=0;t<NSHARES;t++) Variance += VAR[s][t]*frac[s]*frac[t];
p.setObj(Variance); // Set the objective function
// Minimum amount of North-American values
for(s=0;s<NNA;s++) Na += frac[NA[s]];
p.newCtr(Na >= 0.5);
// Spend all the capital
for(s=0;s<NSHARES;s++) Cap += frac[s];
p.newCtr(Cap == 1);
// Target yield
for(s=0;s<NSHARES;s++) Return += RET[s]*frac[s];
p.newCtr(Return >= TARGET);
// Solve the problem
p.minim("");
// Solution printing
cout << "With a target of " << TARGET << " minimum variance is " <<
p.getObjVal() << endl;
for(s=0;s<NSHARES;s++)
cout << s << ": " << frac[s].getSol()*100 << "%" << endl;
return 0;
}
This program produces the following solution output (notice that the default algorithm for solving QP problems is Newton-Barrier, not the Simplex as in all previous examples):
Reading Problem FolioQP
Problem Statistics
3 ( 0 spare) rows
10 ( 0 spare) structural columns
24 ( 0 spare) non-zero elements
76 quadratic elements
Global Statistics
0 entities 0 sets 0 set members
Presolved problem has: 3 rows 10 cols 24 non-zeros
Barrier starts
Matrix ordering - Dense cols.: 9 NZ(L): 88 Flops: 504
Its P.inf D.inf U.inf Primal obj. Dual obj. Compl.
0 1.90e+01 1.00e+03 3.70e+00 8.7840000e+02 -3.8784000e+03 4.4e+04
1 2.70e-01 2.49e+00 5.24e-02 8.2924728e+00 -2.6937881e+03 3.2e+03
2 6.03e-04 2.98e-03 1.47e-04 5.2434336e+00 -9.6149827e+01 1.0e+02
3 3.01e-05 1.35e-04 6.66e-06 4.1180308e+00 -1.0473900e+01 1.5e+01
4 6.22e-07 5.11e-15 5.55e-17 1.6293994e+00 -2.4446102e+00 4.1e+00
5 9.58e-07 1.50e-15 5.55e-17 7.2626663e-01 3.0497504e-01 4.2e-01
6 2.67e-08 1.33e-15 5.55e-17 5.6812712e-01 5.2570789e-01 4.2e-02
7 2.75e-07 1.69e-15 2.78e-17 5.5898050e-01 5.5589613e-01 3.1e-03
8 3.23e-08 1.47e-15 5.55e-17 5.5758318e-01 5.5727441e-01 3.1e-04
9 9.99e-09 1.05e-15 5.55e-17 5.5740957e-01 5.5738509e-01 2.4e-05
10 1.77e-09 8.51e-16 5.55e-17 5.5739375e-01 5.5739329e-01 4.6e-07
11 3.95e-12 5.14e-16 5.55e-17 5.5739341e-01 5.5739341e-01 4.8e-10
Barrier method finished in 0 seconds
Uncrunching matrix
Its Obj Value S Ninf Nneg Sum Inf Time
0 .557393 B 0 0 .000000 0
Optimal solution found
With a target of 9 minimum variance is 0.557393
0: 30%
1: 7.15391%
2: 7.38246%
3: 5.46363%
4: 12.6554%
5: 5.91228%
6: 0.332458%
7: 30%
8: 1.09983%
9: 6.76156e-09%
MIQP
We now wish to express the fact that at most a given number MAXNUM of different assets may be selected into the portfolio, subject to all other constraints of the previous QP model. In Chapter 11 we have already seen how this can be done, namely by introducing an additional set of binary decision variables buys that are linked logically to the continuous variables:
s SHARES: fracs buys
Through this relation, a variable buys will be at 1 if a fraction fracs greater than 0 is selected into the portfolio. If, however, buys equals 0, then fracs must also be 0.
To limit the number of different shares in the portfolio, we then define the following constraint:
s SHARES buys MAXNUM
Implementation with BCL
We may modify the previous QP model or simply append the following lines to the program of the previous section, just after the solution printing: the problem is then solved once as a QP and once as a MIQP in a single program run.
XPRBvar buy[NSHARES]; // 1 if asset is in portfolio, 0 otherwise
// Create the decision variables
for(s=0;s<NSHARES;s++)
buy[s] = p.newVar(XPRBnewname("buy(%d)",s+1), XPRB_BV);
// Limit the total number of assets
for(s=0;s<NSHARES;s++) Num += buy[s];
p.newCtr(Num <= MAXNUM);
// Linking the variables
for(s=0;s<NSHARES;s++) p.newCtr(frac[s] <= buy[s]);
// Solve the problem
p.minim("g");
// Solution printing
cout << "With a target of " << TARGET << " and at most " << MAXNUM <<
" assets, minimum variance is " << p.getObjVal() << endl;
for(s=0;s<NSHARES;s++)
cout << s << ": " << frac[s].getSol()*100 << "% (" << buy[s].getSol()
<< ")" << endl;
When executing the MIQP model, we obtain the following solution output:
Reading Problem FolioQP
Problem Statistics
14 ( 514 spare) rows
20 ( 0 spare) structural columns
54 ( 5056 spare) non-zero elements
76 quadratic elements
Global Statistics
10 entities 0 sets 0 set members
Presolved problem has: 14 rows 20 cols 53 non-zeros
LP relaxation tightened
Barrier starts
Matrix ordering - Dense cols.: 9 NZ(L): 158 Flops: 698
Its P.inf D.inf U.inf Primal obj. Dual obj. Compl.
0 5.10e+01 1.00e+03 9.70e+00 5.4900000e+03 -1.8490000e+04 2.7e+05
1 1.73e+00 1.55e+01 3.30e-01 2.6470516e+01 -1.1284767e+04 1.9e+04
2 1.09e-02 3.97e-02 2.08e-03 3.2753756e+00 -1.7764326e+03 1.8e+03
3 1.08e-04 1.93e-04 1.05e-05 3.1700892e+00 -8.6725881e+00 1.2e+01
4 7.00e-06 3.11e-15 1.11e-16 1.0441961e+00 -7.4322539e-01 1.8e+00
5 1.07e-06 1.46e-15 1.11e-16 5.9697690e-01 4.6389066e-01 1.3e-01
6 2.83e-07 6.69e-16 1.11e-16 5.6159255e-01 5.4968539e-01 1.2e-02
7 6.08e-08 4.14e-15 1.11e-16 5.5787840e-01 5.5703242e-01 8.5e-04
8 1.63e-08 1.59e-14 1.11e-16 5.5745525e-01 5.5736598e-01 8.9e-05
9 4.94e-09 2.63e-13 1.11e-16 5.5739897e-01 5.5739241e-01 6.6e-06
10 8.47e-11 5.50e-13 5.55e-17 5.5739347e-01 5.5739341e-01 6.4e-08
Barrier method finished in 0 seconds
Optimal solution found
*** Heuristic solution found: 1.419000 Time: 0 ***
Generating cuts
Its Type BestSoln BestBound Sols Add Del Gap GInf Time
1 K 1.419000 .557393 1 3 0 154.58% 7 0
2 K 1.419000 .557393 1 3 2 154.58% 7 0
3 K 1.419000 .557393 1 5 3 154.58% 7 0
4 K 1.419000 .557393 1 6 5 154.58% 7 0
5 K 1.419000 .557393 1 8 6 154.58% 7 0
6 K 1.419000 .557393 1 9 8 154.58% 7 0
7 K 1.419000 .557393 1 11 9 154.58% 7 0
8 K 1.419000 .557393 1 12 11 154.58% 7 0
9 K 1.419000 .557393 1 14 12 154.58% 7 0
10 K 1.419000 .557393 1 15 14 154.58% 7 0
11 K 1.419000 .557393 1 17 15 154.58% 7 0
12 K 1.419000 .557393 1 18 17 154.58% 7 0
13 K 1.419000 .557393 1 20 18 154.58% 7 0
14 K 1.419000 .557393 1 21 20 154.58% 7 0
15 K 1.419000 .557393 1 23 21 154.58% 7 0
16 K 1.419000 .557393 1 24 23 154.58% 7 0
17 K 1.419000 .557393 1 26 24 154.58% 7 0
18 K 1.419000 .557393 1 27 26 154.58% 7 0
19 K 1.419000 .557393 1 29 27 154.58% 7 0
20 K 1.419000 .557393 1 30 59 154.58% 7 0
Cuts in the matrix : 1
Cut elements in the matrix : 10
Cuts in the cutpool : 60
Cut elements in the cutpool: 556
*** Heuristic solution found: 1.419000 Time: 0 ***
Branch Parent Solution Entity Value/Bound Active GInf Time
*** Solution found **** Time: 0
40 39 1.248762 5 0 0
*** Search completed *** Time: 0 Nodes: 55
Number of integer feasible solutions found is 2
Best integer solution found is 1.248762
Uncrunching matrix
With a target of 9 and at most 4 assets, minimum variance is 1.24876
0: 30% (1)
1: 20% (1)
2: 0% (0)
3: 0% (0)
4: 23.8095% (1)
5: 26.1905% (1)
6: 0% (0)
7: 0% (0)
8: 0% (0)
9: 0% (0)
The log of the Branch-and-Bound search tells us this time that 2 integer feasible solutions have been found (one of which by the MIP heuristics) and a total of 55 nodes have been enumerated to complete the search.With the additional constraint on the number of different assets the minimum variance is more than twice as large as in the QP problem.
If you have any comments or suggestions about these pages, please send mail to docs@dashoptimization.com.