Quadratic Programming
In this chapter we turn the LP problem from Chapter 15 into a Quadratic Programming (QP) problem, showing how to
- transform a QP model into matrix format,
- input and solve QP problems with Xpress Optimizer.
Chapter 7 shows how to formulate and solve this example with Mosel and in Chapter 12 the same is done with BCL.
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. Which investment strategy should the investor adopt to minimize the variance subject to getting a minimum target yield of 9?
QP model
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:
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:
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 ∈ NA fracs ≥ 0.5
∑s ∈ SHARES fracs = 1
∑s ∈ SHARES RETs·fracs ≥ 9
∀ s ∈ SHARES: 0 ≤ fracs ≤ 0.3
Matrix representation
For the problem input into Xpress Optimizer, the mathematical model is transformed into the following constraint matrix (Table QP matrix).
| frac1 | frac2 | frac3 | frac4 | frac5 | frac6 | frac7 | frac8 | frac9 | frac10 | ||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | Oper. | RHS | ||
| MinNA | 0 | 10 | 13 | 16 | 19 | ≥ | 0.5 | ||||||
| Allfrac | 1 | 11 | 14 | 17 | 110 | 112 | 114 | 116 | 118 | 120 | 122 | = | 1 |
| Yield | 2 | 52 | 175 | 268 | 1211 | 813 | 915 | 717 | 619 | 3121 | 2123 | ≥ | 9 |
| ↑ | ↑ | ||||||||||||
| rowidx | matval | ||||||||||||
| colbeg | 0 | 3 | 6 | 9 | 12 | 14 | 16 | 18 | 20 | 22 | 24 | ||
As in the previous chapters, the superscripts for the matrix coefficients indicate the order of the entries in the arrays rowidx and matval, the first three entries of which are highlighted (printed in italics).
The coefficients of the quadratic objective function are given by the following variance/covariance matrix (Table Variance/covariance matrix).
| frac1 | frac2 | frac3 | frac4 | frac5 | frac6 | frac7 | frac8 | frac9 | frac10 | ||
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | ||
| frac1 | 0 | 0.1 | |||||||||
| frac2 | 1 | 19 | -2 | 4 | 1 | 1 | 1 | 0.5 | 10 | 5 | |
| frac3 | 2 | -2 | 28 | 1 | 2 | 1 | 1 | -2 | -1 | ||
| frac4 | 3 | 4 | 1 | 22 | 1 | 2 | 3 | 4 | |||
| frac5 | 4 | 1 | 2 | 4 | -1.5 | -2 | -1 | 1 | 1 | ||
| frac6 | 5 | 1 | 1 | 1 | -1.5 | 3.5 | 2 | 0.5 | 1 | 1.5 | |
| frac7 | 6 | 1 | 1 | 2 | -2 | 2 | 5 | 0.5 | 1 | 2.5 | |
| frac8 | 7 | 0.5 | -1 | 0.5 | 0.5 | 1 | 0.5 | 0.5 | |||
| frac9 | 8 | 10 | -2 | 3 | 1 | 1 | 1 | 0.5 | 25 | 8 | |
| frac10 | 9 | 5 | -1 | 4 | 1 | 1.5 | 2.5 | 0.5 | 8 | 16 |
Implementation with Xpress Optimizer
The following program folioqp.c loads the QP problem into Xpress Optimizer and solves it. Notice that the quadratic part of the objective function must be specified in triangular form, that is, either the lower or the upper triangle of the original matrix. Here we have chosen the upper triangle, which means that instead of 4·frac2·frac4 + 4·frac4·frac2 we only specify the sum of these terms, 8·frac2·frac4. Due to the input conventions of the Optimizer the values of the main diagonal also need to be multiplied with 2. As with the matrix coefficients, only quadratic terms with non-zero coefficients are specified to the Optimizer (hence the spaces in the array definitions below).
#include <stdio.h>
#include <stdlib.h>
#include "xprs.h"
int main(int argc, char **argv)
{
XPRSprob prob;
int s, status;
double objval, *sol;
/* Problem parameters */
int ncol = 10;
int nrow = 3;
int nqt = 43;
/* Row data */
char rowtype[] = {'G','E','G'};
double rhs[] = {0.5, 1, 9};
/* Column data */
double obj[] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
double lb[] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
double ub[] = {0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3,0.3};
/* Matrix coefficient data */
int colbeg[] = {0, 3, 6, 9, 12, 14, 16, 18, 20, 22, 24};
int rowidx[] = {0,1,2,0,1, 2,0,1, 2,0,1, 2,1,2,1,2,1,2,1,2,1, 2,1,2};
double matval[] = {1,1,5,1,1,17,1,1,26,1,1,12,1,8,1,9,1,7,1,6,1,31,1,21};
/* QP problem data */
int qcol1[] = {0,
1,1,1,1,1,1,1,1,1,
2,2,2,2,2, 2,2,
3, 3,3, 3,3,
4,4,4,4,4,4,
5,5,5,5,5,
6,6,6,6,
7,7,7,
8,8,
9};
int qcol2[] = {0,
1,2,3,4,5,6,7,8,9,
2,3,4,5,6, 8,9,
3, 5,6, 8,9,
4,5,6,7,8,9,
5,6,7,8,9,
6,7,8,9,
7,8,9,
8,9,
9};
double qval[] = {0.1,
19,-2, 4,1, 1, 1,0.5, 10, 5,
28, 1,2, 1, 1, -2, -1,
22, 1, 2, 3, 4,
4,-1.5,-2, -1, 1, 1,
3.5, 2,0.5, 1,1.5,
5,0.5, 1,2.5,
1,0.5,0.5,
25, 8,
16};
for(s=0;s<nqt;s++) qval[s]*=2;
XPRSinit(NULL); /* Initialize Xpress Optimizer */
XPRScreateprob(&prob); /* Create a new problem */
/* Load the problem matrix */
XPRSloadqp(prob, "FolioQP", ncol, nrow, rowtype, rhs, NULL,
obj, colbeg, NULL, rowidx, matval, lb, ub,
nqt, qcol1, qcol2, qval);
XPRSchgobjsense(prob, XPRS_OBJ_MINIMIZE); /* Set sense to maximization */
XPRSlpoptimize(prob, ""); /* Solve the problem */
XPRSgetintattrib(prob, XPRS_LPSTATUS, &status); /* Get solution status */
if(status == XPRS_LP_OPTIMAL)
{
XPRSgetdblattrib(prob, XPRS_LPOBJVAL, &objval); /* Get objective value */
printf("Minimum variance: %g\n", objval);
sol = (double *)malloc(ncol*sizeof(double));
XPRSgetlpsol(prob, sol, NULL, NULL, NULL); /* Get primal solution */
for(s=0;s<ncol;s++) printf("%d: %g%%\n", s, sol[s]*100);
}
XPRSdestroyprob(prob); /* Delete the problem */
XPRSfree(); /* Terminate Xpress */
return 0;
}A QP problem is loaded into the Optimizer with the function XPRSloadqp. This function takes the same arguments as function XPRSloadlp with four additional arguments at the end for the quadratic part of the objective function: the number of quadratic terms, nqt, the column numbers of the variables in every quadratic term (qcol1 and qcol2) and their coefficient, qval.
If we wish to load a Mixed Integer Quadratic Programming (MIQP) problem, then we need to use the function XPRSloadqglobal that takes the same arguments as XPRSloadqp plus the nine arguments for MIP problems introduced with function XPRSloadglobal in the previous chapter.
As opposed to the previous examples we now minimize the objective function. Notice that for solving and solution access we use the same functions as for LP problems. When solving MIQP problems, correspondingly we need to use the MIP solving and solution functions presented in Chapter 16.
Executing this program produces the following output:
Minimum variance: 0.557394 0: 30% 1: 7.15401% 2: 7.38237% 3: 5.46362% 4: 12.6561% 5: 5.91283% 6: 0.333491% 7: 29.9979% 8: 1.0997% 9: 6.97039e-06%
All but the last share are selected into the portfolio (the value printed for 9 is so close to 0 that Xpress Optimizer interprets it as 0 with its default tolerance settings).