{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "6b0b3c2f-3aaa-43fe-a8d0-2b696df5c1a6",
   "metadata": {
    "id": "6b0b3c2f-3aaa-43fe-a8d0-2b696df5c1a6",
    "tags": []
   },
   "source": [
    "```{index} single: application; portfolio\n",
    "```\n",
    "```{index} single: application; investment\n",
    "```\n",
    "```{index} single: solver; ipopt\n",
    "```\n",
    "\n",
    "# Markowitz portfolio optimization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5da22c67-5c34-4c3a-90a4-61222899e855",
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "executionInfo": {
     "elapsed": 19868,
     "status": "ok",
     "timestamp": 1686079839753,
     "user": {
      "displayName": "Marcos Domínguez Velad",
      "userId": "18131612151504710302"
     },
     "user_tz": -120
    },
    "id": "5da22c67-5c34-4c3a-90a4-61222899e855",
    "outputId": "d02c0ec3-0992-411e-857c-30e70d7c2768",
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using default Community Edition License for Colab. Get yours at: https://ampl.com/ce\n",
      "Licensed to AMPL Community Edition License for the AMPL Model Colaboratory (https://colab.ampl.com).\n"
     ]
    }
   ],
   "source": [
    "# Install AMPL and solvers\n",
    "%pip install -q amplpy pandas\n",
    "\n",
    "SOLVER_QO = \"ipopt\"\n",
    "\n",
    "from amplpy import AMPL, ampl_notebook\n",
    "\n",
    "ampl = ampl_notebook(\n",
    "    modules=[\"coin\"],  # modules to install\n",
    "    license_uuid=\"default\",  # license to use\n",
    ")  # instantiate AMPL object and register magics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "b13edf26",
   "metadata": {
    "executionInfo": {
     "elapsed": 24,
     "status": "ok",
     "timestamp": 1686079839755,
     "user": {
      "displayName": "Marcos Domínguez Velad",
      "userId": "18131612151504710302"
     },
     "user_tz": -120
    },
    "id": "b13edf26",
    "tags": []
   },
   "outputs": [],
   "source": [
    "from IPython.display import Markdown, HTML\n",
    "import numpy as np\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a82bd88d-7b16-4c82-b0b7-5dbd0d81b71c",
   "metadata": {
    "id": "a82bd88d-7b16-4c82-b0b7-5dbd0d81b71c",
    "tags": []
   },
   "source": [
    "## Problem description and model formulation\n",
    "\n",
    "A canonical stochastic optimization problem is the so-called portfolio selection problem, also known as **Markowitz portfolio optimization**. Assume an investor has an initial capital $C$ that she wants to invest in $n$ possible risky assets, each of them with an unknown return rate $r_i$, $i=1,\\dots,n$, or in another risk-free asset with guaranteed return rate $R$. Let $x$ be the vector whose $i$-th component $x_i$ describes the amount invested in asset $i$ and $\\tilde{x}$ the amount invested in the risk-free asset. We consider a stochastic model where the return of the $n$ risky assets is then a random vector $r$ with known expected values $\\mu = \\mathbb E r $ and covariance\n",
    "\n",
    "$$\n",
    "    \\Sigma = \\mathbb{E} [ (r-\\mu)(r-\\mu)^\\top].\n",
    "$$\n",
    "\n",
    "The return of the investment $y = R \\tilde{x} + r^\\top x$ then also becomes a random variable with mean\n",
    "\n",
    "$$\n",
    "    \\mathbb{E} y = R \\tilde{x} + \\mathbb{E} r^\\top x = R \\tilde{x} + \\mu^\\top x\n",
    "$$\n",
    "\n",
    "and variance\n",
    "\n",
    "$$\n",
    "    \\mathrm{Var}(y) = \\mathbb{E}(y-\\mathbb{E}y)^2 = x^\\top \\Sigma x.\n",
    "$$\n",
    "\n",
    "The variance of the return of the investment is one possible way to quantify the risk of the investment $x$.\n",
    "\n",
    "The problem the investor is facing is how to select a portfolio that achieves a good compromise between risk and expected return. More specifically, one could try to maximize the expected return $\\mathbb{E} y$ subject to an upper bound on the tolerable risk, obtaining the following optimization problem:\n",
    "\n",
    "$$\n",
    "\\begin{align*}\n",
    "    \\max \\quad  & R \\tilde{x} + \\mu^\\top x \\\\\n",
    "    \\text{s.t.}\\quad\n",
    "    & \\sum_{i=1}^n x_i + \\tilde{x}  = C  \\\\\n",
    "    & x^\\top \\Sigma x \\leq \\gamma^2 \\\\\n",
    "    & \\tilde{x} \\geq 0 \\\\\n",
    "    & x_i \\geq 0 & \\forall \\, i=1,\\dots,n.\n",
    "\\end{align*}\n",
    "$$\n",
    "\n",
    "The first constraint describes the fact that the total amount invested must be equal to the initial capital. The second constraint ensures that the variance of the chosen portfolio is upper bounded by a parameter $\\gamma^2$, which captures the risk the investor is willing to undertake. The last nonnegativity constraint excludes the possibility of short-selling.\n",
    "\n",
    "One can easily show that the quadratic constraint $x^\\top \\Sigma x \\leq \\gamma^2$ is convex thanks to the fact that $\\Sigma$ is positive semidefinite, being a covariance matrix. The Markowitz optimization problem is thus convex. Let us implement it in AMPL."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0e194e8d",
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 76
    },
    "executionInfo": {
     "elapsed": 22,
     "status": "ok",
     "timestamp": 1686079839756,
     "user": {
      "displayName": "Marcos Domínguez Velad",
      "userId": "18131612151504710302"
     },
     "user_tz": -120
    },
    "id": "0e194e8d",
    "outputId": "7ab2b05d-8d94-44dd-dc67-dc4eabfa1045",
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/markdown": [
       "**Solution:** $\\tilde x = 0.158$, $x_1 = 0.561$,  $x_2 = 0.142$,  $x_3 = 0.139$"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/markdown": [
       "**Maximizes objective value to:** $1.17$"
      ],
      "text/plain": [
       "<IPython.core.display.Markdown object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Specify the initial capital, the risk threshold, and the guaranteed return rate.\n",
    "C = 1\n",
    "gamma = 1\n",
    "R = 1.01\n",
    "\n",
    "# Specify the number of assets, their expected return, and their covariance matrix.\n",
    "n = 3\n",
    "mu = np.array([1.2, 1.1, 1.3])\n",
    "Sigma = np.matrix([[1.5, 0.5, 2], [0.5, 2, 0], [2, 0, 5]])\n",
    "Sigma_df = pd.DataFrame(Sigma, index=range(n), columns=range(n))\n",
    "\n",
    "# If you want to change the covariance matrix Sigma, ensure you input a semi-definite positive one.\n",
    "# The easiest way to generate a random covariance matrix is first generating a random m x m matrix A\n",
    "# and then taking the matrix A^T A (which is always semi-definite positive)\n",
    "# m = 3\n",
    "# A = np.random.rand(m, m)\n",
    "# Sigma = A.T @ A\n",
    "\n",
    "\n",
    "def markowitz(gamma, mu, Sigma):\n",
    "    m = AMPL()\n",
    "    m.eval(\n",
    "        \"\"\"\n",
    "      reset;\n",
    "      # Assets\n",
    "      set A;\n",
    "      \n",
    "      param Sigma{A, A};\n",
    "      param mu{A};\n",
    "      param C;\n",
    "      param R;\n",
    "      param gamma;\n",
    "\n",
    "      var xtilde >= 0;\n",
    "      var x{A} >= 0;\n",
    "\n",
    "      maximize objective:\n",
    "        sum{i in A} mu[i]*x[i] + R*xtilde;\n",
    "      \n",
    "      subject to bounded_variance:\n",
    "        sum{i in A, j in A} x[i] * Sigma[i,j] * x[j] <= gamma**2;\n",
    "\n",
    "      subject to total_assets:\n",
    "        sum{i in A} x[i] + xtilde = C;\n",
    "    \"\"\"\n",
    "    )\n",
    "\n",
    "    m.set[\"A\"] = list(range(n))\n",
    "\n",
    "    m.param[\"C\"] = C\n",
    "    m.param[\"R\"] = R\n",
    "    m.param[\"mu\"] = mu\n",
    "    m.param[\"gamma\"] = gamma\n",
    "    m.param[\"Sigma\"] = Sigma_df\n",
    "\n",
    "    m.option[\"solver\"] = SOLVER_QO\n",
    "    m.solve(verbose=False)\n",
    "\n",
    "    return m\n",
    "\n",
    "\n",
    "model = markowitz(gamma, mu, Sigma)\n",
    "\n",
    "x = model.var[\"x\"]\n",
    "display(\n",
    "    Markdown(\n",
    "        f\"**Solution:** $\\\\tilde x = {model.var['xtilde'].value():.3f}$, $x_1 = {x[0].value():.3f}$,  $x_2 = {x[1].value():.3f}$,  $x_3 = {x[2].value():.3f}$\"\n",
    "    )\n",
    ")\n",
    "display(\n",
    "    Markdown(\n",
    "        f\"**Maximizes objective value to:** ${model.obj['objective'].value():.2f}$\"\n",
    "    )\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "dWzByezMVYJb",
   "metadata": {
    "executionInfo": {
     "elapsed": 19,
     "status": "ok",
     "timestamp": 1686079839757,
     "user": {
      "displayName": "Marcos Domínguez Velad",
      "userId": "18131612151504710302"
     },
     "user_tz": -120
    },
    "id": "dWzByezMVYJb"
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.9"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}