```python
from ipywidgets import interactive, interact

import numpy as np
import pandas as pd
import matplotlib as mpl
import matplotlib.pyplot as plt
import cvxpy as cvx
import ipyvolume.pylab as p3

import bokeh.plotting as bp
from bokeh.models import HoverTool

bp.output_notebook()

import itertools
""
import holoviews as hv
hv.extension('bokeh')

%matplotlib inline

plt.rc('figure', figsize=(14, 8))
```


```python

```


```python
def surface3d(f, xs, ys=None, show=True, logcolor=False, key=None):
    if ys is None:
        ys = xs
    xs, ys = np.meshgrid(xs, ys)
    zs = f(xs, ys)
    
    if logcolor:
        znorm = np.log1p(zs - zs.min())
    else:
        znorm = zs
    znorm /= znorm.max()
    color = mpl.cm.plasma(znorm)[..., 0:3]
    
    fig = p3.figure(key=key)
    p3.plot_surface(xs, ys, zs, color=color)
    
    if show:
        p3.show()
        

def isosurface3d(f, xs, ys=None, zs=None, show=True, key=None):
    if ys is None:
        ys = xs
    if zs is none:
        zs = xs
        
    xs, ys, zs = np.meshgrid(xs, ys, zs)
    density = f(xs, ys, zs)
        
    fig = p3.figure(key=key)
    
    if show:
        p3.show()

def line2d(f, xs, fig=None, ax=None, **kwargs):
    xs = np.array(xs)
    if fig is None:
        fig = plt.gcf()
    if ax is None:
        ax = plt.gca()
    ax.plot(xs, f(xs), **kwargs)
    return fig, ax

def scatter2d(f, xs, fig=None, ax=None):
    xs = np.array(xs)
    if fig is None:
        fig = plt.gcf()
    if ax is None:
        ax = plt.gca()
    ax.scatter(xs, f(xs), **kwargs)
    return fig, ax  
        
def scatter3d(xs, ys, zs, show=True, marker='sphere', color='red', key=None):
    p3.scatter(xs, ys, zs, marker=marker, color=color)

def plot_contours(f, x0, x1, levels=[1, 2, 3, 4], title=""):
    xx, yy = np.meshgrid(np.linspace(-2, 2, 50), np.linspace(-2, 2, 50))
    z = f(xx, yy)
    
    line = interpolate(x0, x1)
    curve = (hv.Curve((line[:, 0], line[:, 1]), label='x0 -> x1')
             .options(color='red', line_width=3, line_dash='dashed'))
    
    scatter = (hv.Scatter(([x0[0], x1[0]], [x0[1], x1[1]]), extents=(-2, -2, 2, 2))
               .options(color='red', size=10))
    
    contour = (hv.operation.contours(hv.Image(z, bounds=(-2, -2, 2, 2)), levels=levels)
                .options(line_width=3))
    
    return (hv.Overlay([contour, scatter, curve],  label=title)
            .options(height=450, width=700, legend_position='top_left', show_title=True))

def plot_convex1d(f, x0, x1, title=""):
    # Prepare data
    theta = np.linspace(0, 1, 50)
    xs = x0 * theta + (x1 * (1 - theta))
    ys = f(xs) 
    convex_bounds = f(x0) * theta + f(x1) * (1 - theta) 
    data = bp.ColumnDataSource({'x': xs, 'y': ys, 'theta': theta, 'bound': convex_bounds})
    
    # Plot figure
    hover = HoverTool(mode='vline', names=['f'], 
                      tooltips=[("t", "@theta"),
                                ("x", "@x"),
                                ("f(x)", "@y"), 
                                ('Convex Bound', "@bound")])
    p = bp.figure(plot_width=700, plot_height=450, tools=[hover], title=title)
    p.line('x', 'y', 
           source=data, line_width=5, name='f', legend='f(x0 * t + x1 * (1 - t))')
    p.line('x', 'bound', 
           source=data, line_width=5, color='red', line_dash='dotted', legend='f(x0) * t + f(x1) * (1 - t)')                       
    p.circle([x0, x1], [f(x0), f(x1)], color='red', size=20)
    p.legend.location = 'top_left'
    bp.show(p)
    return p
    
def interpolate(x0, x1, samples=50):
    x0 = np.array(x0)
    x1 = np.array(x1)
    theta = np.linspace(0, 1, samples)[:, None]
    return theta * x0 + (1 - theta) * x1
```

<center>
  <h1>Convex Optimization for Finance</h1>
  <h3>Scott Sanderson</h3>
  <h3>Email: ssanderson@quantopian.com, GitHub: [@ssanderson](https://github.com/ssanderson), Twitter: [@scottbsanderson](https://twitter.com/scottbsanderson)</h3>
</center>

# Goals

- Situate Convex Optimization within Broader Landscape
- Build **Geometric** Intuition for Convex Functions and Convex Sets
- Show Applications of Convex Optimization to Finance
- Provide Resources for Further Study


# Outline

- **Optimization**
- **Convex** Optimization
- Convex Optimization **for Finance**

# Optimization

Mathematical optimization is a family of techniques for finding a "best" value from a set of choices.

A well-posed optimization problem has at least two parts:

- A set $S$ of possible choices. Often called the **optimization domain**.

- An scalar-valued **objective function** $f$ defined on $S$, to be **minimized** or **maximized**:
  - **Loss/Cost Function** (minimization)
  - **Utility Function** (maximization)

We often also add a third component, a set of **constraints** on "valid" elements of $S$.

We write optimization problems like this:

\begin{align}
\underset{x \in S}{\text{maximize}}&& f(x)&\\
\text{subject to}&& c_i(x)& = 1&&\\
\end{align}

**In plain English:** *Find the value of $x$ in $S$ that maximizes $f(x)$ while satisfying all constraints $c_i$.*

## Examples


```python
def rosenbrock(x, y):
    """The Rosenbrock function for a = 1, b = 7.5.
    
    See https://en.wikipedia.org/wiki/Rosenbrock_function.
    """
    a = 1
    b = 5 
    return (a - x) ** 2 + b * (y - x ** 2) ** 2

surface3d(rosenbrock, np.linspace(-1.5, 1.5), np.linspace(-0.5, 1.5), logcolor=True)
```


```python
from scipy.optimize import minimize

minimize(lambda p: rosenbrock(*p), x0=[0.0, 0.0], method='nelder-mead')
```


```python

```


```python
def minimize_traced(f, *args, **kwargs):
    """Call scipy.minimize, and store the points it asks us to evaluate.
    """
    trace = []
    
    def traced_f(point):
        trace.append(tuple(point))
        return f(*point)

    result = minimize(traced_f, *args, **kwargs)
    return result, trace
    
result, trace = minimize_traced(rosenbrock, x0=[-1., 1.], method='nelder-mead')
trace[:10]
```


```python
def show_optimization_trace(x0, y0, method):
    func = rosenbrock; initial_point = [x0, y0]

    result, trace = minimize_traced(func, x0=initial_point, method=method)
    
    xs, ys = zip(*trace)
    zs = list(itertools.starmap(func, trace))
    
    surface3d(func, np.linspace(min(xs) * 1.1, max(xs) * 1.1), np.linspace(min(ys) * 1.1, max(ys) * 1.1), logcolor=True)
    scatter3d(np.array(xs), np.array(ys), np.array(zs), color=np.linspace(0.1, 1, num= 2 * len(xs)))
```


```python
interact(show_optimization_trace, x0=-1.4, y0=1.0, method=['nelder-mead', 'powell']);
```


```python
from scipy.optimize import show_options

show_options() 
```

## Convex Optimization

General purpose optimization is **hard**.

In the most extreme cases, it's **impossible**.

Under modest assumptions, it's *merely* computationally infeasible.


### Challenges

- Number of candidates grows exponentially in the dimensionality of the search space.
- Search algorithms can get "trapped" in local minima and/or saddle points.
- Complex feasible regions are hard to navigate.

We can do much better if we can assume that our objective/constraints are "well-behaved".

**Convex Functions** and **Convex Constraints** are a particularly important class of well-behaved inputs.

## Convex Functions

A function $f(x)$ is **convex** if:

$$f(x_0t + x_1(1 - t)) \leq f(x_0)t  + f(x_1)(1 - t)$$

where

$$t \in [0, 1]$$

**Intuition:**

- Left-hand side is the point "t-percent" of the way between $x_0$ and $x_1$.
- Right-hand side is the point "t-percent" of the way between $f(x_0)$ and $f(x_1)$.
- Inequality says that $f$ is changing faster in the output space than the input space.


```python
plot_convex1d(lambda x: x ** 2 + 2 * x - 1, -2, 2, "Convex Function");
```


```python
plot_convex1d(lambda x: x ** 3 - 2 * x ** 2 - x, -5, 5, "Non-Convex Function");
```

## Multivariate Convexity

Convexity generalizes naturally to arbitrary dimensions. 

The same inequality applies for $x \in \mathbb{R}^n$:

$$f(x_0t  + x_1(1 - t)) \leq f(x_0)t  + f(x_1)(1 - t) $$


```python
def convex_func(x, y): 
    return x ** 2 + 3 * y ** 2

plot_contours(convex_func, (-1, 0), (0, np.sqrt(3) / 3), np.linspace(0.5, 2, 4))
```


```python
def nonconvex_func(x, y): 
    return x ** 2 + 4 * y ** 2 - 3 * np.sin(x ** 2 + y ** 2)

plot_contours(nonconvex_func, (-0.5, 0.666), (0.5, 0.666), levels=np.linspace(0.1, 1, 4))
```

## Convex Sets

Closely-related to **convex functions** is the concept of **convex sets**:

Let $C$ be a subset of $\mathbb{R}^n$, and let $x_0, x_1 \in C$.

$C$ is convex iff:

\begin{align}
t x_0 + (1 - t)x_1 \in& C \\
\forall t \in& [0, 1]
\end{align}

**English:** If two points are in a convex set, every point on the line between them is also in the set.

## Useful Facts about Convexity

- The intersection of two convex sets is also convex.

- Level surfaces of convex **functions** are the **boundaries** of **convex sets**.

- In 2 dimensions, a convex set doesn't intersect with **lines** tangent to its boundary
- In 3 dimensions, a convex set doesn't intersect with **planes** tangent to its boundary.

- In general, an N-dimensional convex set doesn't intersect with its **supporting hyperplanes**.

- A scalar function is convex if it has positive second derivative.
- In general, a twice-differentiable function is convex iff its Hessian is positive semi-definite.

## Convex Optimization

Convexity is a very nice property for optimization.

If the objective function $f$ is convex, and the **feasible region** satisfying all constraints is a convex set, then:

- Local extrema are guaranteed to be global extrema.
- There are efficient (polynomial-time) algorithms to find extrema.

## CVX and Disciplined Convex Programming

Most low-level solvers require you to specify a problem in "canonical form" as a collection of matrices. 

This can be hard to get right, and hard to debug if you get it wrong.

Modelling languages solve this problem by translating from high-level expressions to low-level solver representations.

## CVX and Disciplined Convex Programming

[CVX](http://cvxr.com/cvx/) (MATLAB) and its siblings [CVXPY](https://github.com/cvxgrp/cvxpy) (Python) and [Convex.jl](https://github.com/JuliaOpt/Convex.jl) (Julia) implement a system called [Disciplined Convex Programming](http://dcp.stanford.edu/home).

DCP allows you to build expressions out of a set of primitives with known curvature and sign. The modeling language translates from expressions into the canonical form of a target backend.

## Convex Optimization for Portfolio Construction


```python
import cvxpy as cvx 

expected_returns = np.array([0.001, 0.002, 0.003, 0.004, 0.005])

def maximize_returns(expected_returns):
    """Construct the portfolio that maximizes expected portfolio returns given
    expected security-level returns.
    """
    weights = cvx.Variable(len(expected_returns))
    objective = cvx.Maximize(weights.T * expected_returns)
    problem = cvx.Problem(objective, [cvx.sum_entries(cvx.abs(weights)) <= 1])
    problem.solve()
    return weights.value.round(3).ravel()
```


```python
maximize_returns(expected_returns)
```


```python
def markowitz_portfolio(means, cov, risk_aversion):
    """Generate the optimal fully-invested portfolio for a given risk/returns tradeoff.
    """
    weights = cvx.Variable(len(means))
    
    expected_return = weights.T * means
    expected_vol = cvx.quad_form(weights, cov)
    
    utility = expected_return - risk_aversion * expected_vol
    objective = cvx.Maximize(utility)
    
    constraints = [
        cvx.sum_entries(weights) == 1,  # fully-invested
        weights >= 0,                   # long-only
    ]
    
    problem = cvx.Problem(objective, constraints)
    problem.solve()

    return np.array(weights.value.flat).round(4), expected_return.value, expected_vol.value
```


```python
expected_rets = np.array([0.001, 0.002, 0.003, 0.004, 0.005])
cov = np.array([[0.02, 0. ,  0.  , 0.  , 0.  ],
                [0.  , 0.02, 0.  , 0.  , 0.  ],
                [0.  , 0.  , 0.02, 0.  , 0.  ],
                [0.  , 0.  , 0.  , 0.02, 0.0 ],
                [0.  , 0.  , 0.  , 0. , 0.02]])

weights, rets, var = markowitz_portfolio(expected_rets, cov, 0.2)
print("Weights:", weights); print("Expected Return:", rets); print("Expected Variance:", var)
```


```python
import bokeh.plotting as bp

def plot_efficient_frontier(means, cov, min_aversion=0, max_aversion=1, nsamples=25):
    """Plot the expected return and variance for different levels of risk aversion.
    """
    samples = np.linspace(min_aversion, max_aversion, nsamples)
    
    portfolios = []
    returns = []
    stddevs = []
    for aversion in samples:
        portfolio, ret, var = markowitz_portfolio(means, cov, aversion)
        portfolios.append(portfolio)
        returns.append(ret)
        stddevs.append(var)
        
    fig = bp.figure()
    fig.scatter(stddevs, returns)
    fig.xaxis.axis_label = 'Risk'
    fig.yaxis.axis_label = 'Return'
    return fig
```


```python
bp.show(plot_efficient_frontier(expected_rets, cov))
```

## Example: Constraining Historical Downside Risk


```python
alphas = np.linspace(-2, 2, 8)
rets = np.random.standard_t(df=10, size=(500, 8)) * 0.1
pd.DataFrame(rets).describe()
```


```python
def maximize_alpha_constrain_downside(alphas, returns, percentile, max_loss):
    weights = cvx.Variable(returns.shape[1])

    # Number of worst-case return periods to sample.
    nsamples = round(returns.shape[0] * percentile)

    portfolio_rets = returns * weights
    avg_worst_day = cvx.sum_smallest(portfolio_rets, nsamples) / nsamples
    
    objective = cvx.Maximize(weights.T * alphas)
    constraints = [avg_worst_day >= max_loss]
    
    problem = cvx.Problem(objective, constraints)
    problem.solve()
  
    return weights.value.round(4).ravel()
    
result = maximize_alpha_constrain_downside(alphas, rets, percentile=0.05, max_loss=-0.05)
print("Portfolio:", result)

portfolio_rets = rets.dot(result)
worst_days = portfolio_rets[portfolio_rets <= np.percentile(portfolio_rets, 5)]
print("Average Bad Day:", worst_days.mean())
```

# Review

- Python has lots of great tools for optimization.
  - `scipy.optimize`
  - `cvxpy`
  - ... lots of other tools not covered.

- Convex functions "curve upward" when traversing lines between points.
- Convex sets always contain lines between points.

- Convex optimization can be used to **efficiently** solve a wide array of interesting problems.
- High-level modeling tools like `cvxpy` or `quantopian.optimize` mean that you don't need to be an expert to productively use convex solvers.

## Questions?

- **GitHub:** https://github.com/ssanderson
- **Twitter:** @scottbsanderson
- **Email:** scott@quantopian.com

03935050d0f9aa52c27b6b010f1d2eb0c39dd31b
ipynb
Jupyter Notebook
notebooks/Main.ipynb
ssanderson/convex-opt-webinar
1644755190388de9d4294cbbc70c62e5890e70fd
Apache-2.0
2018-06-28T16:00:41.000Z
2022-03-11T07:42:34.000Z
notebooks/Main.ipynb
ssanderson/convex-opt-webinar
1644755190388de9d4294cbbc70c62e5890e70fd
Apache-2.0
2020-03-24T15:46:42.000Z
2022-03-11T23:21:59.000Z
notebooks/Main.ipynb
ssanderson/convex-opt-webinar
1644755190388de9d4294cbbc70c62e5890e70fd
Apache-2.0
2018-07-14T16:44:58.000Z
2022-03-21T20:21:51.000Z
Qwen/Qwen-72B
1. YES
2. YES