Financial Portfolio Optimization with amplpy

Description: Financial Portfolio Optimization with amplpy and amplpyfinance

Tags: amplpy, amplpyfinance, finance

Notebook author: Filipe Brandão <fdabrandao@gmail.com>

Model author: Harry Markowitz

References:

1. amplpy: https://amplpy.readthedocs.io/

2. amplpyfinance: https://amplpyfinance.readthedocs.io/

3. Cornuejols, G., and Tütüncü, R. (2018). Optimization Methods in Finance (2nd edition): Bond Dedication example. Cambridge University Press.

# Install dependencies
%pip install -q amplpy pandas matplotlib yfinance amplpyfinance

# Google Colab & Kaggle integration
from amplpy import AMPL, ampl_notebook

ampl = ampl_notebook(
    modules=["gurobi"],  # modules to install
    license_uuid="default",  # license to use
)  # instantiate AMPL object and register magics

Downloading data and callculate mu and S

from pypfopt import expected_returns, risk_models
import yfinance as yf

tickers = [
    "MSFT",
    "AMZN",
    "KO",
    "MA",
    "COST",
    "LUV",
    "XOM",
    "PFE",
    "JPM",
    "UNH",
    "ACN",
    "DIS",
    "GILD",
    "F",
    "TSLA",
]

ohlc = yf.download(tickers, period="max")
prices = ohlc["Adj Close"].dropna(how="all")
mu = expected_returns.capm_return(prices)
S = risk_models.CovarianceShrinkage(prices).ledoit_wolf()

[*********************100%***********************]  15 of 15 completed

Minimize volatility

Using EfficientFrontierWithAMPL.min_volatility:

from amplpyfinance import EfficientFrontierWithAMPL

ef = EfficientFrontierWithAMPL(None, S, weight_bounds=(0, 1), solver="gurobi")
ef.min_volatility()
_, sigma1, _ = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.01442725331
11 barrier iterations
No basis.
Annual volatility: 12.0%

Using amplpy directly:

from amplpy import AMPL
import pandas as pd

ampl = AMPL()
ampl.eval(
    r"""
    set A ordered;
    param S{A, A};
    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;
    minimize portfolio_variance:
        sum {i in A, j in A} w[i] * S[i, j] * w[j];
    s.t. portfolio_weights:
        sum {i in A} w[i] = 1;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
print(f"Annual volatility: {sigma*100:.1f}%")

Gurobi 9.5.1: optimal solution; objective 0.01442725331
11 barrier iterations
No basis.
Annual volatility: 12.0%

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Maximize return for a target risk

Using EfficientFrontierWithAMPL.efficient_risk:

from amplpyfinance import EfficientFrontierWithAMPL

ef = EfficientFrontierWithAMPL(mu, S, solver="gurobi")
ef.efficient_risk(target_volatility=0.15)
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.2821809144
10 barrier iterations
No basis.
No dual variables returned.
Expected annual return: 28.2%
Annual volatility: 15.0%
Sharpe Ratio: 1.75

Using amplpy directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param target_volatility;
    param market_neutral default 0;
    set A ordered;
    param S{A, A};
    param mu{A} default 0;
    
    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;
    maximize portfolio_return:
        sum {i in A} mu[i] * w[i];
    s.t. portfolio_variance:
        sum {i in A, j in A} w[i] * S[i, j] * w[j] <= target_volatility^2;
    s.t. portfolio_weights:
        sum {i in A} w[i] = if market_neutral then 0 else 1;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.param["target_volatility"] = 0.15
ampl.param["market_neutral"] = False
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.2821809144
10 barrier iterations
No basis.
No dual variables returned.
Expected annual return: 28.2%
Annual volatility: 15.0%
Sharpe Ratio: 1.75

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Minimizing volatility for a given target return

Using EfficientFrontierWithAMPL.efficient_return:

from amplpyfinance import EfficientFrontierWithAMPL

ef = EfficientFrontierWithAMPL(mu, S, weight_bounds=(None, None), solver="gurobi")
ef.efficient_return(target_return=0.07, market_neutral=True)
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.7924129059
1 simplex iterations
Gurobi 9.5.1: optimal solution; objective 0.008882889913
9 barrier iterations
No basis.
Expected annual return: 7.0%
Annual volatility: 9.4%
Sharpe Ratio: 0.53

Using amplply directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param target_return;
    param market_neutral default 0;

    set A ordered;
    param S{A, A};
    param mu{A} default 0;

    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;

    minimize portfolio_variance:
        sum {i in A, j in A} w[i] * S[i, j] * w[j];
    s.t. portfolio__return:
        sum {i in A} mu[i] * w[i] >= target_return;
    s.t. portfolio_weights:
        sum {i in A} w[i] = if market_neutral then 0 else 1;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.param["target_return"] = 0.07
ampl.param["market_neutral"] = True
ampl.param["lb"] = -1
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.008882889913
9 barrier iterations
No basis.
Expected annual return: 7.0%
Annual volatility: 9.4%
Sharpe Ratio: 0.53

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Maximize the Sharpe Ratio

Using EfficientFrontierWithAMPL.max_sharpe:

from amplpyfinance import EfficientFrontierWithAMPL

ef = EfficientFrontierWithAMPL(mu, S, solver="gurobi")
ef.max_sharpe()
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.3017535499
13 barrier iterations
No basis.
Expected annual return: 25.5%
Annual volatility: 12.9%
Sharpe Ratio: 1.82

Using amplpy directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param risk_free_rate default 0.02;

    set A ordered;
    param S{A, A};
    param mu{A} default 0;

    var k >= 0;
    var z{i in A} >= 0;  # scaled weights
    var w{i in A} = z[i] / k;

    minimize portfolio_sharpe:
        sum {i in A, j in A} z[i] * S[i, j] * z[j];
    s.t. muz:
        sum {i in A} (mu[i] - risk_free_rate) * z[i] = 1;
    s.t. portfolio_weights:
        sum {i in A}  z[i] = k;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.3017535499
13 barrier iterations
No basis.
Expected annual return: 25.5%
Annual volatility: 12.9%
Sharpe Ratio: 1.82

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Maximize quadratic utility

Quadratic utility: $\max_w w^T \mu - \frac \delta 2 w^T \Sigma w$

Using EfficientFrontierWithAMPL.max_quadratic_utility:

ef = EfficientFrontierWithAMPL(mu, S, weight_bounds=(0, 1), solver="gurobi")
ef.max_quadratic_utility(risk_aversion=2, market_neutral=False)
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.2710034744
11 barrier iterations
No basis.
Expected annual return: 31.2%
Annual volatility: 20.1%
Sharpe Ratio: 1.45

Using amplpy directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param risk_aversion default 1;
    param market_neutral default 0;

    set A ordered;
    param S{A, A};
    param mu{A} default 0;
    
    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;

    maximize quadratic_utility:
        sum {i in A} mu[i] * w[i]
        - 0.5 * risk_aversion * sum {i in A, j in A} w[i] * S[i, j] * w[j];
    s.t. portfolio_weights:
        sum {i in A} w[i] = if market_neutral then 0 else 1;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.param["risk_aversion"] = 2
ampl.param["market_neutral"] = False
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.2710034744
11 barrier iterations
No basis.
Expected annual return: 31.2%
Annual volatility: 20.1%
Sharpe Ratio: 1.45

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Maximize return for a target risk with sector constraints

Using EfficientFrontierWithAMPL.add_sector_constraints:

sector_mapper = {
    "MSFT": "Tech",
    "AMZN": "Consumer Discretionary",
    "KO": "Consumer Staples",
    "MA": "Financial Services",
    "COST": "Consumer Staples",
    "LUV": "Aerospace",
    "XOM": "Energy",
    "PFE": "Healthcare",
    "JPM": "Financial Services",
    "UNH": "Healthcare",
    "ACN": "Misc",
    "DIS": "Media",
    "GILD": "Healthcare",
    "F": "Auto",
    "TSLA": "Auto",
}

sector_lower = {
    "Consumer Staples": 0.1,  # at least 10% to staples
    "Tech": 0.05  # at least 5% to tech
    # For all other sectors, it will be assumed there is no lower bound
}

sector_upper = {"Tech": 0.2, "Aerospace": 0.1, "Energy": 0.1, "Auto": 0.15}

ef = EfficientFrontierWithAMPL(mu, S)
ef.add_sector_constraints(sector_mapper, sector_lower, sector_upper)
ef.efficient_risk(target_volatility=0.15)
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.2768198377
9 barrier iterations
No basis.
No dual variables returned.
Expected annual return: 27.7%
Annual volatility: 15.0%
Sharpe Ratio: 1.71

Using amplpy directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param target_volatility;
    param market_neutral default 0;

    set A ordered;
    param S{A, A};
    param mu{A} default 0;

    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;
    
    maximize portfolio_return:
        sum {i in A} mu[i] * w[i];
    s.t. portfolio_variance:
        sum {i in A, j in A} w[i] * S[i, j] * w[j] <= target_volatility^2;
    s.t. portfolio_weights:
        sum {i in A} w[i] = 1;

    set SECTORS default {};
    set SECTOR_MEMBERS{SECTORS};
    param sector_lower{SECTORS} default -Infinity;
    param sector_upper{SECTORS} default Infinity;
    s.t. sector_constraints_lower{s in SECTORS: sector_lower[s] != -Infinity}:
        sum {i in SECTOR_MEMBERS[s]} w[i] >= sector_lower[s];
    s.t. sector_constraints_upper{s in SECTORS: sector_upper[s] != Infinity}:
        sum {i in SECTOR_MEMBERS[s]} w[i] <= sector_upper[s];
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.param["target_volatility"] = 0.15
sectors = set(sector_mapper.values())
ampl.set["SECTORS"] = sectors
for sector in sectors:
    ampl.set["SECTOR_MEMBERS"][sector] = [
        ticker for ticker, s in sector_mapper.items() if s == sector
    ]
ampl.param["sector_lower"] = sector_lower
ampl.param["sector_upper"] = sector_upper
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.2768198377
9 barrier iterations
No basis.
No dual variables returned.
Expected annual return: 27.7%
Annual volatility: 15.0%
Sharpe Ratio: 1.71

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

Maximize return for a target risk with cardinality constraints

Using EfficientFrontierWithAMPL:

ef = EfficientFrontierWithAMPL(mu, S)
ef.ampl.param["card_ub"] = 3
ef.efficient_risk(target_volatility=0.15)
mu1, sigma1, sharpe1 = ef.portfolio_performance(verbose=True)

Gurobi 9.5.1: optimal solution; objective 0.2665926133
926 simplex iterations
204 branch-and-cut nodes
No dual variables returned.
Expected annual return: 26.7%
Annual volatility: 15.0%
Sharpe Ratio: 1.64

Using amplpy directly:

from amplpy import AMPL

ampl = AMPL()
ampl.eval(
    r"""
    param target_volatility;
    param market_neutral default 0;

    set A ordered;
    param S{A, A};
    param mu{A} default 0;
    
    param lb default 0;
    param ub default 1;
    var w{A} >= lb <= ub;
    
    maximize portfolio_return:
        sum {i in A} mu[i] * w[i];
    s.t. portfolio_variance:
        sum {i in A, j in A} w[i] * S[i, j] * w[j] <= target_volatility^2;
    s.t. portfolio_weights:
        sum {i in A} w[i] = if market_neutral then 0 else 1;
    
    var y{A} binary;
    s.t. w_lower{i in A}:
        lb * y[i] <= w[i];
    s.t. w_upper{i in A}:
        w[i] <= ub * y[i];

    param card_ub default Infinity;
    s.t. card_limit:
        sum {i in A} y[i] <= card_ub;
    """
)
ampl.set["A"] = ef.tickers
ampl.param["S"] = pd.DataFrame(ef.cov_matrix, index=ef.tickers, columns=ef.tickers)
ampl.param["mu"] = ef.expected_returns
ampl.param["card_ub"] = 3
ampl.param["target_volatility"] = 0.15
ampl.option["solver"] = "gurobi"
ampl.solve()
sigma2 = ampl.get_value("sqrt(sum {i in A, j in A} w[i] * S[i, j] * w[j])")
mu2 = ampl.get_value("sum {i in A} mu[i] * w[i]")
risk_free_rate = 0.02
sharpe2 = (mu2 - risk_free_rate) / sigma2
print(f"Expected annual return: {mu2*100:.1f}%")
print(f"Annual volatility: {sigma2*100:.1f}%")
print(f"Sharpe Ratio: {sharpe2:.2f}")

Gurobi 9.5.1: optimal solution; objective 0.2665926133
926 simplex iterations
204 branch-and-cut nodes
No dual variables returned.
Expected annual return: 26.7%
Annual volatility: 15.0%
Sharpe Ratio: 1.64

pd.Series(ef.clean_weights()).plot.barh()

<AxesSubplot:>

previous

Employee Scheduling Optimization

next

Google Hashcode 2022