Worst-case risk analysis

Covariance uncertainty

In this example we do worst-case risk analysis using CVXPY. Our setting is a single period Markowitz portfolio allocation problem. We have a fixed portfolio allocation $w \in {\bf R}^n$. The return covariance $\Sigma$ is not known, but we believe $\Sigma \in \mathcal S$. Here $\mathcal S$ is a convex set of possible covariance matrices. The risk is $w^T \Sigma w$, a linear function of $\Sigma$.

We can compute the worst (maximum) risk, over all possible covariance matrices by solving the convex optimization problem

$$ \begin{array}{ll} \mbox{maximize} & w^T\Sigma w \\ \mbox{subject to} & \Sigma \in \mathcal S, \quad \Sigma \succeq 0, \end{array} $$

with variable $\Sigma$.

If the worst-case risk is not too bad, you can worry less. If not, you'll confront your worst nightmare

Example

In the following code we solve the portfolio allocation problem

$$ \begin{array}{ll} \mbox{minimize} & w^T\Sigma_\mathrm{nom} w \\ \mbox{subject to} & {\bf 1}^Tw = 1, \quad \mu^Tw \geq 0.1, \quad \|w\|_1 \leq 2, \end{array} $$

and then compute the worst-case risk under the assumption that $\mathcal S = \left\{ \Sigma^\mathrm{nom} + \Delta \,:\, |\Delta_{ii}| =0, \; |\Delta_{ij}| \leq 0.2 \right\}$.

We might expect that $|\Delta_{ij}| = 0.2$ for all $i \neq j$. This does not happen however because of the constraint that $\Sigma^\mathrm{nom} + \Delta$ is positive semidefinite.

# Generate data for worst-case risk analysis.
import numpy as np
np.random.seed(2)
n = 5
mu = np.abs(np.random.randn(n, 1))/15
Sigma = np.random.uniform(-.15, .8, size=(n, n))
Sigma_nom = Sigma.T.dot(Sigma)
print "Sigma_nom ="
print np.round(Sigma_nom, decimals=2)

Sigma_nom =
[[ 0.58  0.2   0.57 -0.02  0.43]
 [ 0.2   0.36  0.24  0.    0.38]
 [ 0.57  0.24  0.57 -0.01  0.47]
 [-0.02  0.   -0.01  0.05  0.08]
 [ 0.43  0.38  0.47  0.08  0.92]]

# Form and solve portfolio optimization problem.
# Here we minimize risk while requiring a 0.1 return.
from cvxpy import *
w = Variable(n)
ret = mu.T*w 
risk = quad_form(w, Sigma_nom)
prob = Problem(Minimize(risk), 
               [sum_entries(w) == 1, 
                ret >= 0.1,
                norm(w, 1) <= 2])
prob.solve()
print "w ="
print np.round(w.value, decimals=2)

w =
[[-0.01]
 [ 0.13]
 [ 0.18]
 [ 0.88]
 [-0.18]]

# Form and solve worst-case risk analysis problem.
Sigma = Semidef(n)
Delta = Symmetric(n)
risk = quad_form(w.value, Sigma)
prob = Problem(Maximize(risk), 
               [Sigma == Sigma_nom + Delta, 
                diag(Delta) == 0, 
                abs(Delta) <= 0.2])
prob.solve()
print "standard deviation =", sqrt(quad_form(w.value, Sigma_nom)).value
print "worst-case standard deviation =", sqrt(risk).value
print "worst-case Delta ="
print np.round(Delta.value, decimals=2)

standard deviation = 0.168895026598
worst-case standard deviation = 0.422021584411
worst-case Delta =
[[ 0.    0.04 -0.2  -0.    0.2 ]
 [ 0.04 -0.    0.2   0.09 -0.2 ]
 [-0.2   0.2  -0.    0.12 -0.2 ]
 [-0.    0.09  0.12 -0.   -0.18]
 [ 0.2  -0.2  -0.2  -0.18 -0.  ]]