portfolio_opt.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cvxpy as cp\n",
    "import numpy as np\n",
    "from scipy.linalg import sqrtm\n",
    "from scipy.stats import norm as Normal_Dist\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Formulating the portfolio optimization problem as a SOCP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "ones = np.ones(4)\n",
    "r_hat = np.array([0.12, 0.10, 0.07, 0.03])\n",
    "Sigma = np.array([\n",
    "    [0.0064, 0.0008, -0.0011, 0.0],\n",
    "    [0.0008, 0.0025, 0.0, 0.0],\n",
    "    [-0.0011, 0.0, 0.0004, 0.0],\n",
    "    [0.0, 0.0, 0.0, 0.0]\n",
    "])\n",
    "Sigma_root = sqrtm(Sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Here, we form the constraints and solve the problem within the function (so it is useful to use in a for-loop\n",
    "# later on)\n",
    "def optimization(epsilon, r_hat, Sigma_root):\n",
    "    x = cp.Variable(4)\n",
    "    ones = np.ones(4)\n",
    "    norm_inv_cdf_eps = Normal_Dist.ppf(epsilon)\n",
    "    q = -.01\n",
    "\n",
    "    constraints = [\n",
    "        ones*x == 1,\n",
    "        x <= .4,\n",
    "        x[3]<=.2,\n",
    "        x >= .05,\n",
    "        cp.norm(cp.matmul(Sigma_root,x))*(norm_inv_cdf_eps) >= (q-r_hat*x)\n",
    "    ]\n",
    "    obj = cp.Maximize(r_hat*x)\n",
    "    prob = cp.Problem(obj, constraints)\n",
    "    prob.solve()\n",
    "    return prob.value, x.value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal value: 0.09183017078042251\n",
      "Optimal soln: [0.29864447 0.29659824 0.35475729 0.05      ]\n"
     ]
    }
   ],
   "source": [
    "opt_val, opt_soln = optimization(1e-4, r_hat, Sigma_root)\n",
    "print(\"Optimal value:\", opt_val)\n",
    "print(\"Optimal soln:\", opt_soln)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Now lets solve the problem for a variety of values of ε from $10^{-6}$ to $10^{-1}$ and observe the optimal values of $\\hat{r}^{T}x$ versus ε."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# epsilons generated between 1e-6 and 1e-1 in log scale\n",
    "epsilons = np.logspace(-6,-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "rhat_vals = [] # optimal return values (i.e. the prob.value outputs)\n",
    "opt_x_vals = [] # optimal portfolio allocations\n",
    "#going over the values of epsilon\n",
    "for eps in epsilons:\n",
    "    rhat, optx = optimization(eps,r_hat, Sigma_root)\n",
    "    rhat_vals.append(rhat)\n",
    "    opt_x_vals.append(optx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA44AAAFECAYAAABoG2BvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3XeYXVW9//H3dya9955J7yQkMCR0ASmhBgQhYkEF8V7lWn8qV6+o4FVQLCBcFREUVAQEJNIhIfQ2CSGVJJM+aZOekDpl/f6YE++YG5IAOTlT3q/nyZOz9157zXcrT3I+WWuvFSklJEmSJEl6N3m5LkCSJEmSVLMZHCVJkiRJ+2RwlCRJkiTtk8FRkiRJkrRPBkdJkiRJ0j4ZHCVJkiRJ+2RwlCRJkiTtk8FRkiRJkrRPBkdJkiRJ0j41yHUBudShQ4fUu3fvXJchSZIkSTkxZcqUtSmljvtrV6+DY+/evSkqKsp1GZIkSZKUExGx5EDaOVVVkiRJkrRPBkdJkiRJ0j4ZHCVJkiRJ+2RwlCRJkiTtk8FRkiRJkrRPBkdJkiRJ0j4ZHCVJkiRJ+5TV4BgRYyNibkQUR8TVe7l+YkRMjYjyiLhoj2uXRcT8zK/Lqp0/MiJmZPq8OSIic75dRDydaf90RLTN5rNJkiRJUn2RteAYEfnArcCZwFDgYxExdI9mS4FPA3/Z4952wPeAMcBo4HvVguCvgc8BAzK/xmbOXw1MTCkNACZmjiVJkiRJH1CDLPY9GihOKS0EiIi/AuOA2bsbpJQWZ65V7nHvGcDTKaX1metPA2MjYjLQKqX0aub8XcD5wOOZvk/K3P9HYDLwrYP/WJIkSTXLlh1lLFm3LddlSNqHJg3z6d+pRa7LeN+yGRy7A8uqHZdQNYL4fu/tnvlVspfzAJ1TSiszn1cBnd9rwZIkSbVFSom3Sjbxl9eW8I+3VrK9rCLXJUnah5E92/D3Lx6X6zLet2wGx5xJKaWISHu7FhFXAlcCFBQUHNK6JEmSPqgtO8r4+7QV3PPaUmav3EyzRvmcP6obHxrYify8yHV5kt5F66YNc13CB5LN4Lgc6FntuEfm3IHee9Ie907OnO/xLn2ujoiuKaWVEdEVKN1bxyml24DbAAoLC/caLiVJkmqa6SUb+ctrS5nw1gq27apgaNdW/PD8wxg3shstm9TuL6SSar5sBsc3gAER0YeqcDceuPQA730S+FG1BXFOB/4zpbQ+IjZHxNHAa8CngF9l2kwALgOuz/z+8MF5DEmSpNx4Z2c5D09bzj2vL2Xm8s00bZjPeYd349IxBYzo0ZrM4vKSlHVZC44ppfKIuIqqEJgP3JFSmhUR1wJFKaUJEXEU8BDQFjg3In6QUhqWCYjXURU+Aa7dvVAO8AXgD0BTqhbFeTxz/nrgvoi4HFgCXJytZ5MkScqmmcs38efXljJh2nK27qpgcJeWXDduGONGdaeVo4uSciBSqr+zNQsLC1NRUVGuy5AkSWLrznImvLWCe15fyvSSTTRpmMe5I6pGF0f2bOPooqSsiIgpKaXC/bWrk4vjSJIk1RYpJe4vKuG6R2ezZUc5gzq35AfnDeP8Ud1r/WIakuoOg6MkSVKOlG7ZwbcfnMEzc0o5um87vnHGII4oaOvooqQax+AoSZKUA4/PWMm3H5rBtl0VXHPOUD59bG/y3E5DUg1lcJQkSTqENm0v4/sTZvHQm8sZ0aM1P7/4cPp3apnrsiRpnwyOkiRJh8iL89fyjb+9RemWnXzl1AF88eT+NMzPy3VZkrRfBkdJkqQs276rgusfn8MfX1lCv47NeegLxzKiR5tclyVJB8zgKEmSlEVvLt3A1+97i4Vrt/LZ4/rwzbGDaNIwP9dlSdJ7YnCUJEnKgl3llfxq0nxufbaYrq2b8pfPjeHYfh1yXZYkvS8GR0mSpINs3uotfPXeacxasZmLjuzBNecOpVUT92SUVHsZHCVJkg6SisrEHS8u4qdPzaVl4wb89pNHcsawLrkuS5I+MIOjJEnSQbBs/Ta+fv9bvL5oPacP7cyPPjKcDi0a57osSTooDI6SJEkfQEVl4g8vL+bGJ+fSIC+48aOHc+ER3YmIXJcmSQeNwVGSJOl9mrtqC996YDrTlm3k5EEd+eEFw+nepmmuy5Kkg87gKEmS9B7tLK/g1mcX8OvJxbRs0pCbxo/kvMO7Ocooqc4yOEqSJL0HU5Zs4OoHpjO/9B3OH9mNa84dRrvmjXJdliRllcFRkiTpAGzdWc5Pn5zLH19ZTNdWTbjzM0dx8qBOuS5Lkg4Jg6MkSdJ+PDdvDd9+cAYrNm3nk0f34ptjB9OisV+jJNUf/oknSZL0LjZs3cV1j8zmwTeX069jc+7//DEU9m6X67Ik6ZAzOEqSJO0hpcQ/pq/kBxNmsWl7Gf9xSn++eHJ/mjTMz3VpkpQTBkdJkqRqVm7azn89NJOJb5dyeI/W/OmKMQzp2irXZUlSThkcJUmSgMrKxJ9fX8oNj79NeWUl/3X2ED5zXB/y89xiQ5IMjpIkqd6buXwT//X3mUxbtpHj+rfnxxeMoKB9s1yXJUk1hsFRkiTVW5u2l/Hzp+Zy96tLaNe8ET/76OF85IjuRDjKKEnV5WWz84gYGxFzI6I4Iq7ey/XGEXFv5vprEdE7c75RRNwZETMi4q2IOClzvmVETKv2a21E/DJz7dMRsabatSuy+WySJKn2Sinx0JslfPhnz3H3q0v4xNG9mPj1k7jwyB6GRknai6yNOEZEPnArcBpQArwRERNSSrOrNbsc2JBS6h8R44EbgEuAzwGklIZHRCfg8Yg4KqW0BRhZ7WdMAR6s1t+9KaWrsvVMkiSp9pu3egvf/ftMXlu0nsN7tuHOTx/F8B6tc12WJNVo2ZyqOhooTiktBIiIvwLjgOrBcRzw/cznvwG3RNU/8w0FJgGklEojYiNQCLy++8aIGAh0Al7I4jNIkqQ6YuvOcm6eOJ/fv7iI5o0b8KMLhjP+qJ7kufiNJO1XNoNjd2BZteMSYMy7tUkplUfEJqA98BZwXkTcA/QEjsz8/nq1e8dTNcKYqp27MCJOBOYBX00pVf/5kiSpHkop8cTMVVz7yGxWbtrBxYU9+NbYwbRv0TjXpUlSrVFTF8e5AxgCFAFLgJeBij3ajAc+We34H8A9KaWdEfF54I/AKXt2HBFXAlcCFBQUHPzKJUlSjbF47VaumTCL5+etYXCXltxy6SiO7NUu12VJUq2TzeC4nKpRwt16ZM7trU1JRDQAWgPrMqOIX93dKCJepmoUcffx4UCDlNKU3edSSuuq9Xs78JO9FZVSug24DaCwsDDtrY0kSarddpRV8D+TF/Cb5xbQKD+Pa84ZyqeO6UWD/KyuCyhJdVY2g+MbwICI6ENVQBwPXLpHmwnAZcArwEXApJRSiohmQKSUtkbEaUD5HovqfAy4p3pHEdE1pbQyc3geMOegP5EkSarxnn27lO9NmMXS9ds47/BufOfsIXRu1STXZUlSrZa14Jh5Z/Eq4EkgH7gjpTQrIq4FilJKE4DfA3dHRDGwnqpwCVWL3jwZEZVUhc5P7tH9xcBZe5z7UkScB5Rn+vp0Fh5LkiTVUIvWbuVHj83h6dmr6dexOX+5YgzH9u+Q67IkqU6If11bpn4pLCxMRUVFuS5DkiR9AJu2lXHTxPnc9cpiGjfI44un9OeK4/vSqIHTUiVpfyJiSkqpcH/tauriOJIkSftUVlHJn19dwi8nzmfT9jIuKezJ104fSKeWTkuVpIPN4ChJkmqVlBLPzi3lvx+dw4I1Wzm2X3v+6+yhDO3WKtelSVKdZXCUJEm1xturNvPDR+bwYvFa+nZozu2fKuTDQzoREbkuTZLqNIOjJEmq8dZs2cnPn57HvW8spWWThlxzzlA+cXQv32OUpEPE4ChJkmqsHWUV3PHSIv7n2QXsKKvgsmN78+UPD6BNs0a5Lk2S6hWDoyRJqnFSSjw6YyXXP/42JRu2c+qQzvznWYPp17FFrkuTpHrJ4ChJkmqUacs2ct0js5myZAODu7Tkz1eM4Tj3Y5SknDI4SpKkGmHBmnf4+VPzeHTGSjq0aMwNFw7noiN7kp/nwjeSlGsGR0mSlFMrNm7npmfm87epJTRukMeXTunPlR/qR4vGfk2RpJrCP5ElSVJOrHtnJ7c+u4A/vboEgMuO6c0XTu5HhxaNc1yZJGlPBkdJknRIbdlRxu9eWMTvX1jI9rIKLjqyB18+dSDd2zTNdWmSpHdhcJQkSYfEjrIK7n5lCf8zuZgN28o4a3gXvnbaIPp3cqVUSarpDI6SJCmryioqub+ohJsnzmfV5h2cOLAj3zh9EMN7tM51aZKkA2RwlCRJWVFZmXhkxkp+8fQ8Fq3dyhEFbfjl+JEc3bd9rkuTJL1HBkdJknRQpZR4dm4pP31yHnNWbmZwl5bc/qlCPjykExFurSFJtZHBUZIkHRQpJV6Yv5abJs5nypINFLRrxi8vGcl5h3cjz70YJalWMzhKkqQPJKXE5LlruGnifKYt20jX1k344fmHcclRPWmYn5fr8iRJB4HBUZIkvS8pJZ6ZU8rNE+czY/kmurdpyo8uGM6FR3ancYP8XJcnSTqIDI6SJOk9qaxMPDV7FTdNLGbOys0UtGvGTy4cwQVHdHeEUZLqKIOjJEk6IBWVicdnruRXE4uZu3oLfTo052cfPZxxI7vRwMAoSXWawVGSJO1TRWXikekr+NWkYopL36F/pxbcNH4k54zoRr6L3khSvWBwlCRJe1VeUcnD01Zw67PFLFy7lUGdW3LLpaM487CuBkZJqmcMjpIk6V/sKq/k728u55Zni1m6fhtDurbiN584gtOHdnFbDUmqpwyOkiQJgC07yvjr68v4/YuLWLV5B8O7t+Z3nyrk1CGdiDAwSlJ9ltXgGBFjgZuAfOD2lNL1e1xvDNwFHAmsAy5JKS2OiEbAb4FCoBL4ckppcuaeyUBXYHumm9NTSqXv1lc2n0+SpLqgdPMO7nx5MX96dQlbdpRzTN/2/PjC4Zw0sKOBUZIEZDE4RkQ+cCtwGlACvBERE1JKs6s1uxzYkFLqHxHjgRuAS4DPAaSUhkdEJ+DxiDgqpVSZue/jKaWiPX7ku/UlSZL2YsGad/jd8wt5cOpyyisrOfOwrlx5Yl8O79km16VJkmqYbI44jgaKU0oLASLir8A4oHpwHAd8P/P5b8AtUfVPm0OBSQCZ0cSNVI0+vr6Pn7fXvlJK6aA8jSRJdcTUpRv47XMLeGr2ahrl53HxUT244vi+9O7QPNelSZJqqGwGx+7AsmrHJcCYd2uTUiqPiE1Ae+At4LyIuAfoSdX00578b3C8MyIqgAeAH2bC4bv1tbb6D4yIK4ErAQoKCg7Ok0qSVMNVViaenVvKb59byOuL19O6aUOuOrk/lx3bmw4tGue6PElSDVdTF8e5AxgCFAFLgJeBisy1j6eUlkdES6qC4yeperfxgKSUbgNuAygsLHQ0UpJUp+0qr+Thacu57fmFzC99h+5tmnLNOUO55KieNG9cU78GSJJqmmz+jbGcqlHC3Xpkzu2tTUlENABaA+syI4hf3d0oIl4G5gGklJZnft8SEX+hakrsXe/WVxaeS5KkGm/PFVIHd2nJLy8ZydkjutIwPy/X5UmSaplsBsc3gAER0YeqUDceuHSPNhOAy4BXgIuASSmlFBHNgEgpbY2I04DylNLsTCBsk1JaGxENgXOAZ/bVVxafT5KkGmfJuq388eUl3F+0jC07yzm2X3tuuGgEJw7o4AqpkqT3LWvBMfOe4VXAk1Rtx3FHSmlWRFwLFKWUJgC/B+6OiGJgPVXhEqAT8GREVFIVOj+ZOd84c75hps9ngN9lrr1bX5Ik1WkpJV5ZsI47XlrExLdLyY/g7BFdufz4Pozo4QqpkqQPLurzoFxhYWEqKtpzVw9JkmqHHWUVPPTmcv7w0mLmrt5C++aNuHRMAZ84uhedWzXJdXmSpFogIqaklAr318634iVJqmVWbtrOXa8s4Z7Xl7JxWxlDurbiJxeN4LzDu9GkYX6uy5Mk1UEGR0mSaoGUElOXbuCOlxbzxMxVpJQ4bWhnPntcH0b3aef7i5KkrDI4SpJUg+0qr+TRGSu486XFTC/ZRKsmDbj8+D588uhe9GzXLNflSZLqCYOjJEk1UOmWHdzz2jL+9NoS1mzZSb+Ozbnu/MP4yKju7r8oSTrk/JtHkqQaYvfqqH9+bSlPzlpFeWXipEEd+cxxfTihfwfy8pyOKknKDYOjJEk5tnHbLv42pYS/vLaUhWu30rppQz59bG8uHVNA344tcl2eJEkGR0mSciGlxLRlG/nTq0t5ZPoKdpZXckRBG3720cM5e0RXV0eVJNUoBkdJkg6hrTvLeXjaCv782hJmrdhM80b5XHRkDz4+phdDu7XKdXmSJO2VwVGSpENg7qot/OnVJTz05nLe2VnO4C4t+eH5h3H+qO60cLEbSVIN599UkiRlyc7yCh6fsYo/vbqEoiUbaNQgj3OGd+XjRxdwREFb916UJNUaBkdJkg6yt1dt5t43lvH3N5ezYVsZvds34ztnDeHCI3vQrnmjXJcnSdJ7ZnCUJOkg2LyjjAnTVnB/0TLeKtlEw/zg9KFdGD+6J8f1cysNSVLtZnCUJOl9Sinx2qL13PfGMh6buZIdZZUM7tKSa84Zyvmjuju6KEmqMwyOkiS9R6s27eCBqSXcV7SMJeu20bJxAy48ogcXF/ZkRI/WvrsoSapzDI6SJB2AXeWVTHp7Nfe+sYzn5q2hMsGYPu348ocHcOZhXWnayH0XJUl1l8FRkqR9KC7dwr1vLOPBqctZt3UXnVo25t8+1I+LC3vSu0PzXJcnSdIhYXCUJGkP67fu4h9vreDBN5fz1rKNNMgLPjykE5cc1ZMTB3SkQX5erkuUJOmQMjhKkkTVnouT5pTywNTlTJ5bSnllYnCXlnznrCGcP6o7HVs2znWJkiTljMFRklRvpZSYsmQDD765nEfeWsHmHeV0atmYzx7fhwtGdWdI11a5LlGSpBrB4ChJqneWrNvKg1OX8/dpy1mybhtNG+ZzxrDOfOSIHhzXvwP57rkoSdK/MDhKkuqFTdvKeGTGCh6cupwpSzYQAcf2a8+XThnAGYd1oUVj/0qUJOnd+LekJKnO2llewXNz1/DQm8uZOKeUXRWVDOzcgm+NHcz5o7rRtXXTXJcoSVKtYHCUJNUp5RWVvLJwHROmreCJWavYsqOcDi0a8Ymje/GRI7ozrFsrIpyKKknSe5HV4BgRY4GbgHzg9pTS9XtcbwzcBRwJrAMuSSktjohGwG+BQqAS+HJKaXJENAPuB/oBFcA/UkpXZ/r6NPBTYHmm+1tSSrdn8/kkSTVDZWVi6tINTHhrBY/NWMnad3bRsnEDTh/WhfNGduO4fu3dQkOSpA8ga8ExIvKBW4HTgBLgjYiYkFKaXa3Z5cCGlFL/iBgP3ABcAnwOIKU0PCI6AY9HxFGZe25MKT2bCZcTI+LMlNLjmWv3ppSuytYzSZJqjpQSs1Zs5h9vreCR6StZvnE7jRvkceqQzpx7eDdOGtSRJg3zc12mJEl1QjZHHEcDxSmlhQAR8VdgHFA9OI4Dvp/5/DfglqiaPzQUmASQUiqNiI1AYUrpdeDZzPldETEV6JHFZ5Ak1TDFpe/wj7dW8I/pK1i4ZisN8oITB3bkG2cM4tShnV3kRpKkLMjm367dgWXVjkuAMe/WJqVUHhGbgPbAW8B5EXEP0JOqqaw9gdd33xgRbYBzqZoKu9uFEXEiMA/4akqp+s/ffd+VwJUABQUFH+T5JEmHSMmGbTwyfSUTpq1g9srNRMDRfdrzuRP6MnZYF9o2b5TrEiVJqtNq6j/L3gEMAYqAJcDLVL3TCEBENADuAW7ePaIJ/AO4J6W0MyI+D/wROGXPjlNKtwG3ARQWFqZsPoQk6f0r2bCNx2es4rGZK3lz6UYARvZswzXnDOXsEV3p3KpJjiuUJKn+yGZwXE7VKOFuPfjfhWv2bFOSCYOtgXUppQR8dXejiHiZqlHE3W4D5qeUfrn7REppXbXrtwM/ORgPIUk6dJau28ZjM1fy+IyVvFWyCYBh3VrxjTMGce6IbhS0b5bjCiVJqp+yGRzfAAZERB+qAuJ44NI92kwALgNeAS4CJqWUUmb11EgpbY2I04Dy3YvqRMQPqQqYV1TvKCK6ppRWZg7PA+Zk6bkkSQfRorVbeWzGSh6fuZKZyzcDMKJHa64+czBnHtaFXu2b57hCSZKUteCYeWfxKuBJqrbjuCOlNCsirgWKUkoTgN8Dd0dEMbCeqnAJ0Al4MiIqqQqdnwSIiB7Ad4C3gamZfbh2b7vxpYg4DyjP9PXpbD2bJOmDKS59h8dnrOTRGSt5e9UWoGoa6nfOGsLYw7rQs50ji5Ik1SRRNSu0fiosLExFRUW5LkOS6ryUEvNL3+GxGSt5bMZK5q1+B4Aje7XlrOFdGXtYF7q3aZrjKiVJqn8iYkpKqXB/7Wrq4jiSpFqusjIxffkmnpq1iidnrWLBmq1EwFG92/H9c4cy9rCudGntAjeSJNUGBkdJ0kGzq7yS1xat46lZq3l69mpWbd5Bfl4wpk87Pn1sb84Y1oVOroYqSVKtY3CUJH0gW3eW89y8NTw5axWT3i5ly45ymjbM50MDO3L6sM6cMrgTbZq5z6IkSbWZwVGS9J6tfWcnE+es5slZq3mxeC27yitp26whY4d14fRhXTi+fweaNsrPdZmSJOkgOaDgGBH9gJKU0s6IOAkYAdyVUtqYzeIkSTXH0nXbeHLWKp6avYqiJRtICXq0bconxvTi9GGdKezVlgb5ebkuU5IkZcGBjjg+ABRGRH/gNuBh4C/AWdkqTJKUWxWViTeXbmDi26VMnLP6nyuhDunaii+dMoAzhnVhSNeWZLZGkiRJddiBBsfKzL6MFwC/Sin9KiLezGZhkqRDb/OOMl6Yt5aJc1bz7NxSNmwro0FeMLpPOy4u7MkZw9xjUZKk+uhAg2NZRHwMuAw4N3OuYXZKkiQdSkvWbWXinFImvr2a1xaup7wy0aZZQ04e1IlTBnfixIEdad3UP/IlSarPDjQ4fgb4N+C/U0qLIqIPcHf2ypIkZUt5RSVTl25k4pzVTHy7lOLSqimoAzq14PIT+nDqkM6M6tnG9xUlSdI/HVBwTCnNBr5U7XgRcEO2ipIkHVwbt+3i+flrmTRnNZPnrWHjtjIa5gdj+rTn42MK+PDgzhS0dwqqJEnauwNdVfU44PtAr8w9AaSUUt/slSZJer8qKxOzVmxm8txSJs9bw5tLN1CZoF3zRnx4cGc+PKQTJwzoQMsmTkGVJEn7d6BTVX8PfBWYAlRkrxxJ0vu1cdsuXpi/lmfnlvL8vLWsfWcnACN6tOaqUwZw0qCOHN6jDfl5roIqSZLemwMNjptSSo9ntRJJ0nvybqOKbZo15MQBHTlpUEdOHNiRDi0a57pUSZJUyx1ocHw2In4KPAjs3H0ypTQ1K1VJkvZqn6OKJ/fnQ4M6MbKno4qSJOngOtDgOCbze2G1cwk45eCWI0mqrryikrdKNvL8vLW8MH8N05Zt/Oeo4gkDOnKyo4qSJOkQ2G9wjIg84NcppfsOQT2SVO8tXbeN5+ev4YX5a3i5eB1bdpaTFzCiRxtHFSVJUk7sNzimlCoj4puAwVGSsmDzjjJeWbCOF+av4YX5a1mybhsA3ds05ZzDu3LCgI4c2689bZo1ynGlkiSpvjrQqarPRMT/A+4Ftu4+mVJan5WqJKkOq5p+uokX51dNP31z2UYqKhPNG+VzTL/2fPa4PpwwoAN9OjQnwlFFSZKUewcaHC/J/P7FaucS4D6OkrQfKSUWrt3Ky8Vreal4HS8vWMvmHeVEwIjurfn3D/XjhAEdGFXQlkYN8nJdriRJ0v9xQMExpdQn24VIUl2yevMOXipey4vFa3llwTpWbtoBQLfWTTjzsK6cMLADx/XrQNvmTj+VJEk13wEFx4j41N7Op5TuOrjlSFLttGlbGa8srBpNfKl4LQvWVM3qb9usIcf0a89V/auCYq/2zZx+KkmSap0Dnap6VLXPTYAPA1MBg6Okemn7rgqKlqz/59TTmcs3UZmgacN8RvdpxyVH9eTYfh0Y2rUVea5+KkmSarkDnar6H9WPI6IN8NesVCRJNVB5RSXTl2/i5cz006lLNrKropIGecGogjb8xykDOK5/B0b2bON7ipIkqc450BHHPW0F9vveY0SMBW4C8oHbU0rX73G9MVWjlkcC64BLUkqLI6IR8FugEKgEvpxSmpy550jgD0BT4LHMtRQR7aha9bU3sBi4OKW04X0+n6R6LqXEgjVb//me4qsLqvZTBBjStRWXHduLY/t3YHTvdjRv/H7/KJUkSaodDvQdx39QtYoqQB4wFLh/P/fkA7cCpwElwBsRMSGlNLtas8uBDSml/hExHriBqhVcPweQUhoeEZ2AxyPiqJRSJfDrzPXXqAqOY4HHgauBiSml6yPi6szxtw7k+SQJ/nVBm5eK17J6804Aerar2k/xuP4dOKZve9q3aJzjSiVJkg6tA/1n8hurfS4HlqSUSvZzz2igOKW0ECAi/gqMA6oHx3HA9zOf/wbcElWrRgwFJgGklEojYiNQGBHLgFYppVczfd4FnE9VcBwHnJTp64/AZAyOkvZh844yXlu4npcyQXF+6TtA1YI2x/bvwPGZBW0K2jfLcaWSJEm5daDB8ayU0r+EsIi4Yc9ze+gOLKt2XAKMebc2KaXyiNgEtAfeAs6LiHuAnlRNZe1J1bTV6oG1JNMHQOeU0srM51VA570VFRFXAlcCFBQU7KN8SXVNSomZyzfz9OxVvFi8lrdKNlFRmWjSMI/Rfdrz0cIeHNe/A0O6uKCNJElSdQcaHE/j/47enbmXcwfLHcAQoAhYArwMVBzozZl3HtO7XLsNuA2gsLBwr20k1R0VlYmpSzfwxMxVPDFzFcs3bicv4PCebfjCSf04rn8HRhW0oXGD/FyXKkmSVGPtMzhGxL8DXwD6RsT0apdaAi/tp+/lVI0S7tYjc25vbUoiogHQGliXUkrAV6vV8TIwD9iQ6Wdvfa6OiK4ppZUR0RW6fgPhAAAgAElEQVQo3U99kuqosopKXl24jidmruKp2atZs2UnjRrkceKADnzl1AGcOqQzbZs3ynWZkiRJtcb+Rhz/QtX7gz+marGZ3baklNbv5943gAER0YeqcDceuHSPNhOAy4BXgIuASZnRwmZApJS2RsRpQPnuRXUiYnNEHE3V4jifAn61R1/XZ35/eD/1SapDdpRV8ML8tTwxcxXPzFnNpu1lNGuUz8mDOzF2WBdOHtyJFq5+KkmS9L7s81tUSmkTsAn4WEQcDwxIKd0ZER0iok9KadE+7i2PiKuAJ6najuOOlNKsiLgWKEopTQB+D9wdEcXAeqrCJUAn4MmIqKQqdH6yWtdf4H+343g88wuqAuN9EXE5VdNbLz7g/xUk1Urv7Czn2bdLeWLWKp59u5Rtuypo3bQhpw7pzNjDunDCgA40aegUVEmSpA8qqmaF7qdRxPeo2lNxUEppYER0A+5PKR2X7QKzqbCwMBUVFeW6DEnvQWVl4slZq3hgagnPz1/LrvJKOrRozBnDqsLi0X3b0zA/L9dlSpIk1QoRMSWlVLi/dgc6b+sCYBQwFSCltCIiWn6A+iTpPUkp8dy8Nfz0ybnMWrGZ7m2a8okxvThzeBeOKGhLvqugSpIkZc2BBsdd1VcqjYjmWaxJkv7FlCXrueGJuby+aD092zXl5xcfzriR3Q2LkiRJh8iBBsf7IuK3QJuI+BzwWeD27JUlSTB7xWZufGouk94upWPLxlw3bhiXHFVAowZORZUkSTqUDig4ppRuzKxuuhkYBFyTUno6q5VJqrcWr93Kz5+ex4S3VtCqSQO+OXYQnz62N80auSqqJElSLhzwt7BMUHwaICLyIuLjKaU/Z60ySfXOqk07uGnifO4rWkaj/Dy+cFI/Pn9iP1o3a5jr0iRJkuq1fQbHiGgFfBHoTtU+iU9njv8f8BZgcJT0gW3YuotfP7eAP768mMqU+MSYAr54Sn86tWyS69IkSZLE/kcc7wY2AK8AVwDfBgI4P6U0Lcu1Sarj3tlZzh0vLuJ3zy/knV3lXDCqO189dSA92zXLdWmSJEmqZn/BsW9KaThARNwOrAQKUko7sl6ZpDqrrKKSu19Zwq3PFrNu6y5OH9qZ/3fGIAZ2dpcfSZKkmmh/wbFs94eUUkVElBgaJX0QS9Zt5Ut/ncZbyzZybL/2fOOMQYwqaJvrsiRJkrQP+wuOh0fE5sznAJpmjgNIKaVWWa1OUp2RUuKBqcv53sMzyc8Lbrl0FOeM6JbrsiRJknQA9hkcU0r5h6oQSXXXpu1lfOehGTwyfSWj+7Tjl5eMpFubprkuS5IkSQfITdEkZdXri9bz1XunsWrzDr5xxiD+7UP9yM+LXJclSZKk98DgKCkryisquXnifG55tpgebZvxt387xncZJUmSaimDo6SDbum6bXz53jd5c+lGLjyiBz8YN4wWjf3jRpIkqbbym5ykg+qhN0v47t9nEQG/+tgozj3cBXAkSZJqO4OjpINi844yvvv3mTw8bQVH9W7LLy4ZSY+2zXJdliRJkg4Cg6OkD2zKkvV8+a/TWLlpB187bSBfOKkfDfLzcl2WJEmSDhKDo6T3rbyikl9NKuZXk+bTvW1T7vv8MRzZywVwJEmS6hqDo6T3Zdn6bXzl3mlMWbKBj4zqzg/GDaNlk4a5LkuSJElZYHCU9J6klLivaBnXPTKHAG4aP5JxI7vnuixJkiRlkcFR0gFbvXkHVz8wnWfnruHovu346UWH07OdC+BIkiTVdQZHSfuVUmLCWyu45uFZ7Cir4HvnDuWyY3qTlxe5Lk2SJEmHgMFR0j6te2cn3314Jo/NWMWogjbc+NHD6dexRa7LkiRJ0iGU1fXyI2JsRMyNiOKIuHov1xtHxL2Z669FRO/M+YYR8ceImBERcyLiPzPnB0XEtGq/NkfEVzLXvh8Ry6tdOyubzybVB0/NWsUZv3yep2ev5ptjB3H/548xNEqSJNVDWRtxjIh84FbgNKAEeCMiJqSUZldrdjmwIaXUPyLGAzcAlwAfBRqnlIZHRDNgdkTck1KaC4ys1v9y4KFq/f0ipXRjtp5Jqi82bS/jB/+YxYNTlzO0ayv+dMUYBndpleuyJEmSlCPZnKo6GihOKS0EiIi/AuOA6sFxHPD9zOe/AbdERAAJaB4RDYCmwC5g8x79fxhYkFJakrUnkOqh5+et4Zt/m86ad3bypVP6c9UpA2jUIKuTEyRJklTDZfPbYHdgWbXjksy5vbZJKZUDm4D2VIXIrcBKYClwY0pp/R73jgfu2ePcVRExPSLuiAh3IZfeg607y/nOQzP41B2v06JJAx7892P52umDDI2SJEnK7juOH8BooALoBvQBvh4RfXdfjIhGwHnA/dXu+TXQj6qprCuBn+2t44i4MiKKIqJozZo1WSpfql1eW7iOM296gb+8vpTPndCHR/7jeA7v2SbXZUmSJKmGyOZU1eVAz2rHPTLn9tamJDMttTWwDrgUeCKlVAaURsRLQCGwMHPfmcDUlNLq3R1V/xwRvwMe2VtRKaXbgNsACgsL0/t+OqkO2FFWwY1PzuX3Ly2iZ9tm3HvlMYzu0y7XZUmSJKmGyeaI4xvAgIjokxkhHA9M2KPNBOCyzOeLgEkppUTV9NRTACKiOXA08Ha1+z7GHtNUI6JrtcMLgJkH6TmkOmnq0g2cffML3P7iIj4+poDHv3yCoVGSJEl7lbURx5RSeURcBTwJ5AN3pJRmRcS1QFFKaQLwe+DuiCgG1lMVLqFqNdY7I2IWEMCdKaXp8M8geRrw+T1+5E8iYiRVC+ss3st1SVS9y3jjU3P5w8uL6dKqCXd9djQnDuyY67IkSZJUg0XVAF/9VFhYmIqKinJdhnTIPD9vDf/54AyWb9zOp47pxTfHDqZF42zOWJckSVJNFhFTUkqF+2vnN0apHti4bRfXPTKHB6aW0Ldjc+7/t2M4qrfTUiVJknRgDI5SHZZS4rEZq/jehJls2FbGF0/ux3+cMoAmDfNzXZokSZJqEYOjVEet3ryD7/59Jk/NXs1h3Vvxx8+OZli31rkuS5IkSbWQwVGqY1JK3PvGMv77sTnsKq/kP88czOXH96FBfk3dtlWSJEk1ncFRqkMWr93Kfz44g1cWrmNMn3Zcf+EI+nRonuuyJEmSVMsZHKU6oLyikjteWsTPn55Hw7w8fnTBcMYf1ZO8vMh1aZIkSaoDDI5SLTdn5Wa+9cB0ppds4tQhnfnh+YfRpXWTXJclSZKkOsTgKNVSO8oquGVSMb95bgGtmzbklktHcfbwrkQ4yihJkqSDy+Ao1UKT55ZyzcOzWLp+Gx85ojvfPXsobZs3ynVZkiRJqqMMjlItsmrTDq57ZDaPzlhJ347N+csVYzi2f4dclyVJkqQ6zuAo1QLlFZXc9coSfv70PHZVVPL10wZy5Yf60rhBfq5LkyRJUj1gcJRquGnLNvKdh2Ywa8VmPjSwI9eOG0av9m6xIUmSpEPH4CjVUJu2l/HTJ9/mz68tpVPLxtx66RGcNbyLi99IkiTpkDM4SjVMSomHp63gh4/OYf3WnXz62N587bSBtGzSMNelSZIkqZ4yOEo1yII17/Ddv8/k5QXrOLxHa/7wmaM4rHvrXJclSZKkes7gKNUAO8oq+J/JC/jN5AU0bpjHdecfxqWjC8jPc1qqJEmScs/gKOXYc/PWcM3DM1mybhvnj+zGt88eQqeWTXJdliRJkvRPBkcpR5at38aPHpvD4zNX0bdDc/58xRiOc09GSZIk1UAGR+kQ21FWwW+eW8CvJy8gL8I9GSVJklTjGRylQySlxBMzV/HDR+ewfON2zhnRlW+fNYRubZrmujRJkiRpnwyO0iEwb/UWfvCPWbxUvI7BXVry1yuP5ui+7XNdliRJknRADI5SFm3aXsYvn5nHXa8soUXjBlw7bhiXji6gQX5erkuTJEmSDpjBUcqCysrE/VOW8ZMn5rJ+2y4uHV3A108fRLvmjXJdmiRJkvSeZTU4RsRY4CYgH7g9pXT9HtcbA3cBRwLrgEtSSosjoiFwO3BEpsa7Uko/ztyzGNgCVADlKaXCzPl2wL1Ab2AxcHFKaUM2n0/am6lLN/D9CbOYXrKJo3q35Y/njuaw7q1zXZYkSZL0vmUtOEZEPnArcBpQArwRERNSSrOrNbsc2JBS6h8R44EbgEuAjwKNU0rDI6IZMDsi7kkpLc7cd3JKae0eP/JqYGJK6fqIuDpz/K1sPZ+0p9ItO7jh8bk8MLWEzq0ac9P4kZx3eDciItelSZIkSR9INkccRwPFKaWFABHxV2AcUD04jgO+n/n8N+CWqPqWnYDmEdEAaArsAjbv5+eNA07KfP4jMBmDow6BXeWV/OHlRdw8sZhd5ZX8+0n9uOrk/jRv7ExwSZIk1Q3Z/GbbHVhW7bgEGPNubVJK5RGxCWhPVYgcB6wEmgFfTSmtz9yTgKciIgG/TSndljnfOaW0MvN5FdD5ID+P9C9SSjwzp5QfPzaHhWu38uHBnfivc4bSp0PzXJcmSZIkHVQ1dUhkNFXvMHYD2gIvRMQzmdHL41NKyyOiE/B0RLydUnq++s0ppZQJlv9HRFwJXAlQUFCQ1YdQ3TVz+Sb++9E5vLJwHf06NufOTx/FyYM75bosSZIkKSuyGRyXAz2rHffInNtbm5LMtNTWVC2ScynwREqpDCiNiJeAQmBhSmk5QEqpNCIeoipkPg+sjoiuKaWVEdEVKN1bUZkRytsACgsL9xoupXezatMObnyq6j3GNk0bcu24YXxsdAEN3V5DkiRJdVg2v+2+AQyIiD4R0QgYD0zYo80E4LLM54uASSmlBCwFTgGIiObA0cDbEdE8IlpWO386MHMvfV0GPJyVp1K9tG1XOb94eh4n3ziZCdNWcOUJfZn8jZP51DG9DY2SJEmq87I24ph5Z/Eq4EmqtuO4I6U0KyKuBYpSShOA3wN3R0QxsJ6qcAlVq7HeGRGzgADuTClNj4i+wEOZVSobAH9JKT2Rued64L6IuBxYAlycrWdT/VFRmXhgagk/e2ouqzfv5OwRXbl67GB6tmuW69IkSZKkQyaqBvjqp8LCwlRUVJTrMlRDvVy8lh8+OofZKzczsmcbvnvOEI7s1S7XZUmSJEkHTURMSSkV7q9dTV0cR8qZBWve4cePzeGZOaV0b9OUmz82inNHdHU/RkmSJNVbBkcpY/3WXdz0zDz+/NpSmjTM51tjB/OZ43rTpGF+rkuTJEmScsrgqHpvZ3kFd728hJsnzWfrznIuHVPAV04dSIcWjXNdmiRJklQjGBxVb1VUJh6etpyfPTWP5Ru3c9Kgjnz7rCEM7Nwy16VJkiRJNYrBUfVOSoln55bykyfm8vaqLRzWvRXXXzicEwZ0zHVpkiRJUo1kcFS9MnXpBq5//G1eX7SeXu2b8auPjeLs4V3Jy3PhG0mSJOndGBxVLxSXbuEnT8zlqdmr6dCiMdeNG8YlRxXQqEFerkuTJEmSajyDo+q0lZu288un53P/lGU0a9SAr582kM8e34fmjf1PX5IkSTpQfntWnbRx2y5+PXkBf3h5MSnBp4/twxdP7kd7V0qVJEmS3jODo+qUHWUV3PnSYn49uZgtO8u5YGR3vnraQHq2a5br0iRJkqRay+CoOqG8opK/TSnhl8/MZ9XmHZw8qCPfHDuYIV1b5bo0SZIkqdYzOKpWq6xMPDpjJb94Zh4L12xlVEEbbho/kjF92+e6NEmSJKnOMDiqVkop8dTs1fzi6Xm8vWoLAzu34DefOJIzhnUmwq01JEmSpIPJ4KhaJaXE5Hlr+MXT85hesok+HZpz0/iRnDOiG/nuxShJkiRlhcFRtcbLC9bys6fmMWXJBnq0bcpPLhrBR0Z1p0G+ezFKkiRJ2WRwVI1XtHg9P3tqHq8sXEeXVk344fmHcXFhTxo1MDBKkiRJh4LBUTXW9JKN/OypeTw3bw0dWjTmmnOGcumYApo0zM91aZIkSVK9YnBUjTNn5WZ+/vQ8np69mjbNGnL1mYP51DG9aNbI/1wlSZKkXPCbuGqM4tIt/OKZ+Tw6fSUtGzfga6cN5DPH9aZlk4a5Lk2SJEmq1wyOyrni0ne49dliHp62nCYN87nq5P587oS+tG5mYJQkSZJqAoOjcmbuqi38atJ8Hp2xkiYN8rnihL58/sS+tG/RONelSZIkSarG4KhDbvaKzfxq0nwen7mK5o3y+fcP9ePy4/sYGCVJkqQayuCoQ2ZGySZunjSfp2evpmXjBnzplP585rg+tG3eKNelSZIkSdoHg6Oybtqyjdw8cT6T3i6lVZMGfOXUAXzm2D6+wyhJkiTVElkNjhExFrgJyAduTyldv8f1xsBdwJHAOuCSlNLiiGgI3A4ckanxrpTSjyOiZ6Z9ZyABt6WUbsr09X3gc8CaTPffTik9ls3n075NWbKemyYW8/y8NbRp1pBvnDGITx7Ti1aukipJkiTVKlkLjhGRD9wKnAaUAG9ExISU0uxqzS4HNqSU+kfEeOAG4BLgo0DjlNLwiGgGzI6Ie4CdwNdTSlMjoiUwJSKertbnL1JKN2brmXRgXlu4jpsnzeel4nW0b96Iq88czCeO7kWLxg5wS5IkSbVRNr/JjwaKU0oLASLir8A4oHpwHAd8P/P5b8AtERFUjSY2j4gGQFNgF7A5pbQeWAmQUtoSEXOA7nv0qRxIKfHygnXcPHE+ry1aT4cWjfmvs4dw6ZgCmjUyMEqSJEm1WTa/0XcHllU7LgHGvFublFJ5RGwC2lMVIsdRFRKbAV/NhMZ/iojewCjgtWqnr4qITwFFVI1MbtizqIi4ErgSoKCg4H0+mnarrEw8PWc1/zN5AW8t20jnVo353rlD+djoApo0zM91eZIkSZIOgpo6FDQaqAC6AW2BFyLimWqjly2AB4CvpJQ2Z+75NXAdVaOV1wE/Az67Z8cppduA2wAKCwtTlp+jziqrqGTCtBX85rkFzC99h57tmvLD8w/joiN7GBglSZKkOiabwXE50LPacY/Mub21KclMS21N1SI5lwJPpJTKgNKIeAkoBBZmFs55APhzSunB3R2llFbv/hwRvwMeOfiPpO27KrivaBm3Pb+Q5Ru3M7hLS24aP5Kzh3elQX5ersuTJEmSlAXZDI5vAAMiog9VAXE8VYGwugnAZcArwEXApJRSioilwCnA3RHRHDga+GXm/cffA3NSSj+v3lFEdE0prcwcXgDMzNJz1Uubtpfxp1eXcMeLi1i3dRdH9mrLteOGccrgTlT93yJJkiSprspacMy8s3gV8CRV23HckVKaFRHXAkUppQlUhcC7I6IYWE9VuISq1VjvjIhZQAB3ppSmR8TxwCeBGRExLdN297YbP4mIkVRNVV0MfD5bz1aflG7ZwR0vLubPry5hy85yPjSwI184qR+j+7QzMEqSJEn1RKRUf1/zKywsTEVFRbkuo0Zatn4bv31+AfcVlVBWUclZw7vy7x/qx2HdW+e6NEmSJEkHSURMSSkV7q9dTV0cRzkyd9UWfj25mH9MX0lewIVH9ODzH+pHnw7Nc12aJEmSpBwxOIqUEq8tWs/vnl/IxLdLadYon88c25srTuhLl9ZNcl2eJEmSpBwzONZj5RWVPD5zFb97YSHTSzbRrnkjvnLqAC47pjdtmzfKdXmSJEmSagiDYz20dWc59xUt4/cvLqJkw3b6dGjOf19wGBce4R6MkiRJkv4vg2M9Urp5B394eTF/enUJm3eUU9irLdecM5RTh3QmL88VUiVJkiTtncGxHpi/egu3Pb+Qh6etoKyykrHDunDFCX05slfbXJcmSZIkqRYwONZRKSVeWbiO3z2/kGfnrqFJwzzGj+7JZ4/rQ29XSJUkSZL0Hhgc65jyikoem7mK3z2/kBnLN9G+eSO+dtrA/9/evcfaVZdpHP8+0GnF4qAc6wXKcIoQZ6iIhqbMwFw63k1IUQQBMRGGqOjoxGgmYmZiZiYxqDMOk4IJYoQiUbQ2hOAF/cNLkEuUVu4SxqZUKMHQMoUJ0kxr+84fezk5U9rdc9lr3873k+yctdf6rd969zlPdva71zp7894/PZYj/cAbSZIkSbNg4zgmntm5m3V3PcbaO7bw+NM7OW7JYi476yTe+fqj/cAbSZIkSXNi4zjiNm97lrV3bGH9xq08t2sPK5cdyT+tXs4b//hlfuCNJEmSpJ6wcRxBVcVtm7ZzzW2P8OOHt7Hw0ENY/bqjuOj0SZYfdcSgy5MkSZI0ZmwcR8jOXXu46Z7Huea2R/jVk8/y0sMX8rE3ncAFpx7LkhctGnR5kiRJksaUjeMIeOKZnVx/56/5+s8f5enndrP8qD/kC+eczBknv5JFC/z/RUmSJEntsnEcYnc/uoNrbt/CLfc/wd4q3nLiK7jo9ElWLjuSxP9flCRJktQfNo5DZveevdzywG+49vZHuPvRp3nRogVceNok7zttkmOOfOGgy5MkSZI0D9k4Dpm/u+FubnngN0xOvJB/Xr2cd52ylMMX+WeSJEmSNDh2JEPmwtMmOfuUpfz1q/06DUmSJEnDwcZxyJx63MSgS5AkSZKk/+eQQRcgSZIkSRpuNo6SJEmSpK5sHCVJkiRJXdk4SpIkSZK6snGUJEmSJHXVauOY5G1JHk6yKcml+9m+KMk3m+0/SzLZrP+DJNcluT/JQ0k+dbA5kyxr5tjUzLmwzccmSZIkSfNFa41jkkOBLwJvB04Ezk9y4j7DLgZ2VNXxwOXA55r15wCLquok4BTgg0kmDzLn54DLm7l2NHNLkiRJkuaozTOOK4FNVbW5qnYB3wDO3GfMmcB1zfJ64I1JAhSwOMkC4DBgF/DfB5qz2ecNzRw0c76jvYcmSZIkSfNHm43j0cBjU+5vbdbtd0xV/Q54Bpig0wD+FngCeBT4t6r6ry5zTgBPN3Mc6FgAJPlAkg1JNmzbtm32j06SJEmS5olh/XCclcAe4ChgGfCJJMf1YuKqurqqVlTViiVLlvRiSkmSJEkaawtanPtx4Jgp95c26/Y3ZmtzWeoRwFPAe4DvV9Vu4MkktwMr6Jxt3N+cTwEvTrKgOeu4v2M9z8aNG7cn+fWUVUfQOes5HQcbO5ftLwW2T7OOYTKT398wHWu2c81mv+nuM51x3caYr+E51lzmmum+/cpXt+3mq7/Hmm/5gtHMmPnq3XjztX/9ypj5Gs98HTutUVXVyo1OU7qZzhnDhcC9wPJ9xvwtcFWzfB6wrln+JHBts7wY+CXw2m5zAt8CzmuWrwI+PIuar+7V2LlsBza09Xdp8zaT398wHWu2c81mv+nuM51xB8mQ+RqSY81lrpnu2698ddtuvvp7rPmWr2bbyGXMfPVuvPlq/+/er+OYr9G7tXapanXO/H0E+AHwEJ2m8MEk/5JkdTPsK8BEkk3Ax4Hff73GF4HDkzwI3EWnibzvQHM2+3wS+Hgz10Qz90x9u4dj57p9FPXzMfXyWLOdazb7TXef6YzrNsZ8Dc+x5jLXTPftV75mcqxRYb56N958PZ/56t1487V//XpM5mt+5guANN2vhkiSDVW1YtB1aDyZL7XJfKltZkxtMl9q06jna1g/HGe+u3rQBWismS+1yXypbWZMbTJfatNI58szjpIkSZKkrjzjKEmSJEnqysZRkiRJktSVjaMkSZIkqSsbxxGT5JAkn0lyRZL3DboejZckq5L8NMlVSVYNuh6NnySLk2xIcsaga9F4SfInzXPX+iQfGnQ9Gj9J3pHky0m+meQtg65H4yXJcUm+kmT9oGs5EBvHPkpyTZInkzywz/q3JXk4yaYklx5o/8aZwFJgN7C1rVo1enqUrwKeBV6A+dIUPcoXdL5zd107VWpU9SJfVfVQVV0CvBs4vc16NXp6lLGbqur9wCXAuW3Wq9HSo3xtrqqL2610bvxU1T5K8pd0XpR/tape06w7FPhP4M10XqjfBZwPHApcts8Uf9PcdlTVl5Ksr6qz+1W/hluP8rW9qvYmeTnw71V1Qb/q13DrUb5OBibovDGxvaq+05/qNex6ka+qejLJauBDwPVV9fV+1a/h16uMNft9AfhaVf2iT+VryPU4X0P7+n7BoAuYT6rq1iST+6xeCWyqqs0ASb4BnFlVlwHPu5QryVZgV3N3T3vVatT0Il9T7AAWtVGnRlOPnr9WAYuBE4GdSb5XVXvbrFujoVfPX1V1M3Bzku8CNo76Pz16DgvwWeAWm0ZN1ePXYEPLxnHwjgYem3J/K3Bql/E3Alck+Qvg1jYL01iYUb6SnAW8FXgxcGW7pWkMzChfVfUPAEkupDm73Wp1GnUzff5aBZxF502v77VamcbFTF+DfRR4E3BEkuOr6qo2i9PIm+lz2ATwGeD1ST7VNJhDxcZxxFTVc8BQX/+s0VVVN9J5c0JqTVWtHXQNGj9V9RPgJwMuQ2OsqtYAawZdh8ZTVT1F5/9nh5YfjjN4jwPHTLm/tFkn9YL5UpvMl9pkvtQ2M6Y2jV2+bBwH7y7ghCTLkiwEzgNuHnBNGh/mS20yX2qT+VLbzJjaNHb5snHsoyQ3AHcCr06yNcnFVfU74CPAD4CHgHVV9eAg69RoMl9qk/lSm8yX2mbG1Kb5ki+/jkOSJEmS1JVnHCVJkiRJXdk4SpIkSZK6snGUJEmSJHVl4yhJkiRJ6srGUZIkSZLUlY2jJEmSJKkrG0dJkmYpyZ4k90y5XTqLOVYkWdMsX5jkyt5XKknS3CwYdAGSJI2wnVX1urlMUFUbgA09qkeSpFZ4xlGSpB5LsiXJ55Pcn+TnSY5v1p+T5IEk9ya5tVm3Ksl39jPHZJIfJbkvyQ+T/FGzfm2SNUnuSLI5ydn9fXSSpPnIxlGSpNk7bJ9LVc+dsu2ZqjoJuBL4j2bdp4G3VtXJwOqDzH0FcF1VvRb4GrBmyrZXAn8OnAF8thcPRJKkbrxUVZKk2et2qeoNU35e3izfDqxNsg648SBz/4ibu+wAAAD9SURBVBlwVrN8PfD5Kdtuqqq9wC+TvHzmZUuSNDOecZQkqR2173JVXQL8I3AMsDHJxCzn/p8py5nlHJIkTZuNoyRJ7Th3ys87AZK8qqp+VlWfBrbRaSAP5A7gvGb5AuCnbRUqSdLBeKmqJEmzd1iSe6bc/35V/f4rOV6S5D46ZwfPb9b9a5IT6Jwl/CFwL/BXB5j7o8C1Sf6eTpN5Uc+rlyRpmlJVBx8lSZKmLckWYEVVbR90LZIk9YKXqkqSJEmSuvKMoyRJkiSpK884SpIkSZK6snGUJEmSJHVl4yhJkiRJ6srGUZIkSZLUlY2jJEmSJKkrG0dJkiRJUlf/C0A0QjtoCMpCAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1080x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plotting returns\n",
    "plt.figure(figsize=(15,5))\n",
    "plt.plot(epsilons, rhat_vals)\n",
    "plt.xlabel('Epsilon')\n",
    "plt.ylabel('Returns')\n",
    "plt.xscale('log')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe that as the risk tolerance $\\epsilon$ decreases our optimal mean returns, $\\hat{r}^Tx^*$ decreases as well. Thus less risk implies less reward meaning there is a risk-reward tradeoff. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3sAAAFECAYAAACefieyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3XmUnPdd5/v3t6q7q1qrHW/Rai2WZVuy7DiSN7IYkmBnnKMEJkBiGwXIjQ13wh3gxjgcPM5MfDgDhGGIQ2YShwTHzECSy4VcD2PIcAi5SQbHkXNJDHbAaJwYy5uk1tL7UlW/+8dTvaq71d3qp6u7+v06p05VPfWr5/m2/FinP/ptkVJCkiRJktRcCo0uQJIkSZI0/wx7kiRJktSEDHuSJEmS1IQMe5IkSZLUhAx7kiRJktSEDHuSJEmS1IQMe5IkSZLUhAx7kiRJktSEDHuSJEmS1IRaGl3AbJ1//vlpy5YtjS5DkiRJkhriW9/61rGU0gVnarfkwt6WLVt44oknGl2GJEmSJDVERDw3k3YO45QkSZKkJmTYkyRJkqQmZNiTJEmSpCZk2JMkSZKkJmTYkyRJkqQmZNiTJEmSpCZk2JMkSZKkJpRb2IuIz0TEkYj4+yk+j4h4ICIORcSTEXFNXrVIkiRJ0nKTZ8/eQ8At03z+VmBH/XEn8J9zrEWSJEmSlpXcwl5K6avA8WmavB14OGW+AZwTEevyqkeSJEmSlpOWBl57A/D8mPeH68deakw5c/cLf/EL/PmhP6daqza6FEmSJEnz5NYdt/LRt3600WXMWSPD3oxFxJ1kQz3ZvHlzg6uZ3HMnn2OgOtDoMiRJkiTNk2pa2p05jVyN8wVg05j3G+vHTpNSejCltDeltPeCCy5YkOJm43du+R0uWnVRo8uQJEmSNI9+64d/q9ElnJVGhr1HgAP1VTmvB06llJbcEE5JkiRJWoxyG8YZEX8E3AScHxGHgQ8BrQAppU8AjwL/AjgE9AI/nVctkiRJkrTc5Bb2UkrvPsPnCfhXeV1fkiRJkpazRg7jlCRJkiTlxLAnSZIkSU3IsCdJkiRJTciwJ0mSJElNyLAnSZIkSU3IsCdJkiRJTSi3rRckNY+LirCnBDtbYUUBSgFtkT2f9prR16WAUiEoR1AuhP+6JEmSlpTWxw7A67/Q6DLmzLAnaUQp4Io22NOWhbs9JXhNucB5hdppbRNQI0hE/UiMfBJAkCiQ6i3TgtQvSZI0r1ZtbXQFZ8WwJy1TG1tOD3U7WmoU65mtWg9vRU4PepBFu6JBTpIkNbM9/67RFZwVw57U5M4pwK627LG7BFeV4KpSgbVjeusqFGmhOu57RUOcJEnSkmbYk5rEivoQzN0l2N2WPfaUi6wrjoa4KoX68MrxvXUTg54kSZKWPsOetMS0BVzWmoW6XfVQd1W5yMUto4Etm0sHxdN66yYfkilJkqTmY9iTFqkVAZe1weVtWY/d5W1wZanI1pbqyLy6GpAonBbqCg7BlCRJWvYMe1KDrS2MD3TZUMwim8f11EGNwmnDLQsjn0qSJEnjGfakBbKumPXU7RwT6q4sFbioOBrWhrcyOL2njtPm2UmSJEnTMexJ86gUcElrPdTVn68oBTvbglUxGtaqBMHpAa7gVgaSJEmaJ4Y9aQ4uKGaLpOxsywLdZW2wq1RkU3F0Ph1kWxoUU5WI8QHObQ0kSZKUN8OeNIXVBdjRmj0ubYNLW+GyUnBpa7CmML6XjkmGXkJ9S4M47bAkSZKUO8OelrW2gO2tWZC7tC0LdllvXZELx+xPl8j2qCum2hS9dPbUSZIkaXEx7KnptQVsa4FL2rL5dMOPy0pFNp427LJAgRqFCb10AbRQs5dOkiRJS4ZhT02hPWDbmCA3HOx2lopsKFQpjAlpVQpAmmLYpSteSpIkqTkY9rRknFPIhlxua82et7fCjja4tK3IuuL44Jb10KXTeugAigY6SZIkLQOGPS0aBWBjy/hAt60VdrQV2N4KawvjQ9p0gc4eOkmSJC13hj0tqDUF2NoCW1uzx3Cgu7RUZHOxRuuYxU9qQI1itqLlJAx0kiRJ0tQMe5pXpYAtY8Lc1tYs3G1vK7Btkt65bP4ck86fK8CkvXaSJEmSzsywp1lpBTa1ZoHu4tbRXrptrcH2tgIXTZg7VyOoUZiyd875c5IkSVI+DHsapy1gcwtsqQe6La1w8XCgayty0YSVLbP954oUqRKT9s5NPqdOkiRJUr4Me8vMqsh65Da3jD5vaYUtrcG21tN75s4U5rL95wxzkiRJ0mJj2GsiAby6CJvrvXHjAl1bgS0tsGbCnLnhRVAMc5IkSVJzMewtIasLWXDbNOaxuTV73tZWZMOE1Sxh7AIok8+NcxEUSZIkqTkZ9haJUowPcZtasoVQhnvlNrbA6hgf2IaHWBaoThnYXABFkiRJWp4MewugPWBDS7Zh+NhH1jNXYFMrnFc4PZRVKBCkKQObQywlSZIkTcWwNx96nuMHWvv5oTXjw9zFrVmP3DlTBrnpe97cNFySJEnSXBn25sN3fpU/PPfIyFuDnCRJkqRGKzS6gKZw+QegdOHI2xZqzpWTJEmS1FCGvflw7tVQLDe6CkmSJEkaYdiTJEmSpCZk2JMkSZKkJmTYkyRJkqQmZNiTJEmSpCaUa9iLiFsi4h8j4lBEfHCSzzdHxF9HxN9GxJMR8S/yrEeSJEmSlovcwl5EFIGPA28FrgDeHRFXTGh2L/CFlNJrgHcB/ymveiRJkiRpOcmzZ+9a4FBK6dmU0iDwOeDtE9okYE399VrgxRzrkSRJkqRloyXHc28Anh/z/jBw3YQ2/xb4HxHx88BK4M051iNJktQQCRhseRVdpUvoKl9Cd3k7XeVLGGg9r9GlSZrGmiO9vHbd0t1PO8+wNxPvBh5KKf2HiLgB+IOI2J1Sqo1tFBF3AncCbN68uQFlSpIknVmNFnpKW+iqh7ks1O2gq30HQ8U1ow1TjSCRoti4YiWd0bmFpf3/aJ5h7wVg05j3G+vHxnovcAtASumxiCgD5wNHxjZKKT0IPAiwd+/elFfBkiRJM9XTtpljq6+ns3xpFurad9LTtpEUo79eFVKVGgWIGP/lKOAvNNLid+UFa87caBHLM+wdBHZExFaykPcu4LYJbf4ZeBPwUERcDpSBoznWJEmSNCcDLa/i6OrXcWT16zmy9o30tm3MPpiml65mz52kBsot7KWUKhHxfuBLQBH4TErpqYj4MPBESukR4P8EPhURv0g2nP2nUkr+Q5ckSWq4SmEFx1Zdx5E1r+fomjdyqv1yACJVxwc7e+kkLVK5ztlLKT0KPDrh2H1jXj8N/ECeNUiSJM1EjRaOr7yGo2tex5E1b+D4ymuyIZnjlxJwnp2kJaPRC7RIkiQ1RCI41X4FR1ffyJE1b+DY6huoFlZASvVhmfUdqiLPnaokKT+GPUmStCwkCpxs38Wx1TeMPIZXyBw3NDOCRExzJklaGgx7kiSpKdUocnLFbo6tvpGjq26gY/X1VIqrgNPn3Tk0U1IzMuxJkqSmUKOFkyv3cHRV1mvXsfo6KoWVgOFO0vJk2JMkSUtSosCJFXs4uuZ1HF19Ix2rrqNaaAcMd5IEhj1JkrREJKC7tC3bCmH16zmy5g0Oy5SkaRj2JEnSotXfcn492L2eI2tvoq/11YDhTpJmwrAnSZIWjUqhnWOrrs/C3Zqb6Gy/DDDcSdJcGPYkSVLDJIITK67ilTU3ZRuZr9rrRuaSNE8Me5IkaUENFVZxZM0beHntm3n5nB9moOVVbmQuSTkw7EmSpNx1t13My+e8mZfWvoVjq28kRYsbmUtSzgx7kiRp3tVooWPVXl5e+xZeOudmusvbgPFz7xyaKUn5MuxJkqR5MVA8l1fW/iAvr30zr6x9E0PF1fW5d6M9dgY8SVo4hj1JkjQn2abmV3JkzU28vPaHOL7ytRCFCcMznXsnSY1i2JMkSTPW1/pqXlnzRo6suYkja25isGUtAJFqI8HO3jtJWhwMe5IkaUrVKHNs9XW8suaNvLL2TXSVdwCT7XtnD54kLTaGPUmSNCIBneWdHFlzE6+s/UGOrbqeWqHNuXeStAQZ9iRJWub6Wy7g6OrX8cqaN3Bk7Q/S33ohgHPvJGmJM+xJkrTMDBbP4ejq6zm6+nUcXfMGusqXAJMNzbT3TpKWMsOeJElNbqiwko5V19bD3es52X5FfdXM2riNzA13ktRcDHuSJDWZapQ4vvK1HFnzAxxd/XpOrLyaFC31eXeMWTXToZmS1MwMe5IkLXHVaOPEiqs5tvoGjq5+HR2r9tUXVUlAGp1vZ7iTpGXFsCdJ0hJTKaygY+Vejq2+jmOrbuDEymuycMfERVWCsStoSpKWF8OeJEmL3GDxHI6tupaOVddxbPUNnFxxZRboUiJI44ZjOu9OkjTMsCdJ0iLT13IhHauv49iq6zm2+kY623dmH0yYc0fEuAVWJEkay7AnSVIDJYKu8naOr9xHx6q9HFt9Iz2liwFOWy3TOXeSpNkw7EmStIAqhXZOrHgNHav20rFqH8dX7WOouAaYbJ87w50kae4Me5Ik5SQBfa0b6Fi1l+Or9tKx6lpOtV8xEujcxFySlCfDniRJ86QWrZxsvyILdiv30bH6evpbLwSGh2Qybiim4U6SlCfDniRJc5CAnraLObHyNRxf+RqOr3otp9p3j2yBUEhVag7JlCQ1kGFPkqQZGCyu5cTKqzm+8hqOr7yGEytfy2DLOcDkvXY1e+0kSQ1m2JMkaYLh4ZgnVl5T77XbS09pS/ZhSgQ1F1KRJC16hj1J0rJWo0hX+RJOrryKEyuu4sTKqznVvmtkOObERVSyve3stZMkLX6GPUnSspEo0FXezskVe0aD3YorqRZKQD3YUYAY3dvORVQkSUuVYU+S1JQSQXdpGydW7snC3cqrObliD9VCO+DqmJKk5mfYkyQtecNDMU+t2M3JFbs5uWIPJ1fuoVJYCUwV7JxnJ0lqboY9SdKSUol2Otsvq4e63ZxasYdT7ZdRGxmKabCTJAkMe5KkRWyweA4nV+ziVHu9x27lVXSVto0Eucnn2BnsJEkCw54kaRFIBD2liznVfjmn2ndxqv1yTq68ir629SNtTlsVE+fYSZI0HcOeJGlBDRVW0dl+OSdXXFF/3k1n++UjC6dMto8dGOwkSZotw54kKRejvXVXZI8VV3BqxW562zaOtJlsGKb72EmSND9yDXsRcQvwUaAI/F5K6dcnafPjwL8FEvCdlNJtedYkSZpfCRhouYDO9svobN/JqfbL6Gy/nM72y+ytkySpgXILexFRBD4OvAU4DByMiEdSSk+PabMD+BXgB1JKJyLiwrzqkSSdvcHiGrrKOznVvnMk0HWuuILB4tqRNvbWSZK0OOTZs3ctcCil9CxARHwOeDvw9Jg27wM+nlI6AZBSOpJjPZKkGaoU2ukqX0Jn+bKRHrvOFbvoa71opE0W6mLcFgdgb50kSYtFnmFvA/D8mPeHgesmtLkUICL+J9lQz3+bUvqLiSeKiDuBOwE2b96cS7GStBwNFVbSVd5BV/uldJYvpau8g84VV9DbtmG0UaoRJIdgSpK0xDR6gZYWYAdwE7AR+GpEXJlSOjm2UUrpQeBBgL1796aFLlKSlrps+GU90LXvyJ5XXE5f66tHG00R6ogC/sUrSdLSk2fYewHYNOb9xvqxsQ4Dj6eUhoDvRcQzZOHvYI51SVJTSkB/60VZT135kvpjB13tl9Pfev5Iu0i1rP3EzccNdZIkNZUzhr2IKAH/Etgytn1K6cNn+OpBYEdEbCULee8CJq60+UXg3cDvR8T5ZMM6n51p8ZK0HNWile7SlpFA112+hK7ypXS1X0KlsHKkXaQqEKeFutNCniRJakoz6dn7f4BTwLeAgZmeOKVUiYj3A18im4/3mZTSUxHxYeCJlNIj9c9+OCKeBqrA3Smljtn+EJLUjAaK59Jd3j6ml+4Sutp30tu2cdxQy0KqUKM4fvVLnFMnSdJyN5OwtzGldMtcTp5SehR4dMKx+8a8TsAv1R+StOxUo0x3eSvdpW10lbeNhLvu8iUMFdeMNpxqPh1Qi0ZPv5YkSYvRTH5D+Jv6oil/l3s1ktSEEgV62zbQVd5Od2k73eVt2ev2HfS1rhvXdqpeOufTSZKk2ZpJ2Hsd8FMR8T2yYZxB1im3J9fKJGkJSQR9revrvXRb6Clty16Xt9PTtplaoW2k7VRz6cBeOkmSNH9m8lvFW3OvQpKWgCzQraO7vJWe0ha6RwLdJZMEuhpMMezSuXSSJGkhnDHspZSei4irgNfXD30tpfSdfMuSpMaoUaSvbUO9d25LPdRdTHf7cKArjTaeZh6dK15KkqRGm8nWC/8aeB/wJ/VD/yUiHkwpfSzXyiQpJ5Vop6e0mZ7SVnpKF9NdujgLduXt9LatI40dSjlNoHMenSRJWsxmMozzvcB1KaUegIj4DeAxwLAnacbOPhQVqBZKVKOdaqE88qgVylSjPO5YtVAe0y57Hmg5j57yVnpK28ZtMA7Tz6Ez0EmSpKVqJmEvyPbAG1atH5O0DAwUX8WpFZfT2X45XeXtVAvt1KJELVqpFtpGXtcKZWqFErVoo1p/zh6t9UcLNGhoY6QaQW3yVS5xDp0kSWpOMwl7vw88HhF/Wn//DuDT+ZUkqRFq0UpXaTunVlzBqfYrONV+OZ0rrhzXCxapStT7uRIBRPY8SYBaTFIUSDiHTpIkLS8zWaDltyPiK2RbMAD8dErpb3OtSlJuEtDfciGdKy4fCXWnVuymq7ydFK31RrVsj5UJPXEpig5plCRJWiKmDHsRsSal1BkRrwK+X38Mf/aqlNLx/MuTdDaGCqvpbN+ZPco76Wy/jFMrdjPYcs5Im0KqnL63m/PUJEmSlrzpevb+EHgb8C3Gr60Q9ffbcqxL0ixUop2u9h0jga6zfSedK3bR13rRSJtItex/5Am9dW7iLUmS1Jym/C0vpfS2+vPWhStH0nSq0UZ3adtooGvfSWf7LnraNoyGuCm2CnDfN0mSpOVlJvvs/VVK6U1nOiZp/lQK7XSVL6GrfCmd5R10tV9KZ/vl9LRtGhPqEkHt9JUkHYIpSZIkpp+zVwZWAOdHxLmMbrewBtiwALVJTW+wuIau8o4s1LXvoKu8g872y+lrWz/aaMpQFyTcMkCSJEmTm65n7y7gF4D1ZPP2hsNeJ/C7OdclNY0E9LWuo7t8Sb23bnsW8Novn7CtQQ0mGX5pqJMkSdJcTDdn76PARyPi51NKH1vAmqQlqRplustb64HuErpL2+lqv3RkI/JhkapATLKtgXPqJEmSNH9mss/exyJiN3AFUB5z/OE8C5MWowT0t15Ed2lbFujK27Nw176T3tZ14+bTFahSo3jahuOn9dxJkiRJOZjJAi0fAm4iC3uPAm8Fvg4Y9tS0Botr6S5to7u8beQ5C3fbZtRLRwS1M//vJUmSJOVmJr+NvhO4CvjblNJPR8RFwH/Jtywpf5Vop6e8pd5LNxzqttNdvmTcpuMjvXST7EdnL50kSZIWq5mEvb6UUi0iKhGxBjgCbMq5LmleVArt9JS20F3aSndpa/a6vJWe8vZxG44DFFJl0mGX9tJJkiRpKZrJb7BPRMQ5wKfIVuXsBh7LtSppFoYKq7IgV95Cd2kLPcPBrrx93GqXkA27TMTo3LoxJuu5kyRJkpaqmSzQ8r/XX34iIv4CWJNSejLfsqRRCehvuZDe0sV0ly6mZ+Sxle7yNgZbzh3XfsoeOhx2KUmSpOVjJgu0/Ajw5ZTSqZTS9yPinIh4R0rpiwtQn5aJarTR27axHuK20FPaXH/eSk9pM9VCebTxNCtdgj10kiRJEsxsGOeHUkp/OvwmpXSyvkKnYU8zloCBlgvoKW2mt20zPaVNo4GuvJW+lovGDa2ccoNxcA6dJEmSNAMz+Y15sp2e/U1bpxksrhkJctlz/XVpK71tG6kWSuPaZ8MtC5POn3ODcUmSJOnszHSBlt8GPl5//6/IFmrRMjNUWEVv20Z6S5vobdtUH3a5iZ7SFnpLmxkqrh7XPtuDbup5cg63lCRJkvIzk9+2fx74N8Dn6+//kizwqclkPXP1IFfaQG/bJnraNtFb2kxvaRNDxTXj2keqESRqU4Q5F0ORJEmSGmcmq3H2AB+MiNXZ29Sdf1mab4mgv+UC+to20FvaSG/bBnrbNtLXtoGe0sX0tm2kUlw17jvTzpsjG2qZFqB2SZIkSbM3k9U4rwQeBl5Vf38MeE9K6e9zrk2zUI0SfW3rs2GW9SA3EuhKm+lrezW1aB33nTMNs3TenCRJkrR0zWQY5yeBX0op/TVARNwEPAjcmGNdGiNRoL/1wnqY20Bf63p629Zn70ub6GvbwEDLeRO+NLw9weQLoIDDLCVJkqRmNpOwt3I46AGklL4SEStzrGlZScBg8VX1ILeevrYN4173ljbR33IBacJiJmfqlXN7AkmSJGl5m0kaeDYi/g3wB/X3dwDP5ldS80gEAy3n0de2jr7W9fXnddlz23r62jbR13oRtULbhC/WKFCbctNwsFdOkiRJ0vRmEvZ+Bvh3wJ/U33+tfkx1x3oHOXnuj9N3bvuYILeRvtYLSBPmyc0kyBGFbPilJEmSJM3RTFbjPAH8HwtQy5L1/VO9/POrf8kgJ0mSJGnRmDLsRcR/g6lX1k8p7c+loiVo9wWrOdIzQH8Vg5wkSZKkRWG6nr3fWrAqlrhyS5GYqidPkiRJkhpgyrCXUvp/F7IQSZIkSdL8mW4Y598x+TDOAGoppatyq0qSJEmSdFamG8b5tkmOBbAJ+JV8ypEkSZIkzYfphnE+N/w6Il4D3Ab8GPA94P/OvzRJkiRJ0lxNuXRkRFwaER+KiH8APgb8MxAppR9MKf3uTE4eEbdExD9GxKGI+OA07f5lRKSI2Dvrn0CSJEmSdJrphnH+A9kG6m9LKR0CiIhfnOmJI6IIfBx4C3AYOBgRj6SUnp7QbjXwr4HHZ1m7JEmSJGkK020K96PAS8BfR8SnIuJNZHP2Zupa4FBK6dmU0iDwOeDtk7S7H/gNoH8W55YkSZIkTWO6OXtfBL4YESvJQtovABdGxH8G/jSl9D/OcO4NwPNj3h8GrhvbICKuATallP57RNw9lx9AkiRpMUspURnoY6DzJIOdp+jvOslg50kGuk5RHehrdHmSptG1YQPvec97Gl3GnE03jBOAlFIP8IfAH0bEuWSLtNwDnCnsTSsiCsBvAz81g7Z3AncCbN68+WwuK0mSlItatcpgdycDnScZ6DqZhbuuUwx0nqQ6NDjSLiIgglSrNbBaSTOxY8eORpdwVs4Y9sZKKZ0AHqw/zuQFsm0ahm2sHxu2GtgNfCUiAF4NPBIR+1NKT0y47sg19+7dO9nef5IkSQumOjRIz5GX6Dn6Ur3H7iQDPV2QRn9NKRQK1CYJdCmlce0kLV7XXntto0s4K7MKe7N0ENgREVvJQt67yLZvACCldAo4f/h9RHwF+MDEoCdJktRotWqVvo4jdL/yAt2vvEBfx5EstEUQU/TSTRb0JGkh5Rb2UkqViHg/8CWgCHwmpfRURHwYeCKl9Ehe15YkSTobKSX6T3bQ/cqLdL/yAr1HXqJWrQDZMMw03DOX0uhrSVpk8uzZI6X0KPDohGP3TdH2pjxrkSRJms5gdyfdL2c9dz1HXqQykC0UHoXCuJ47w52kpSLXsCdJkrRYDfX10HPkJbpfeZGeV15gsKcLmCTcORxT0hJl2JMkScvCYE83PUdfoufIi/QeeYmB7k4AIgqkZLiT1HwMe5IkqSkNdneNhLueIy9N3XOXDHeSmpNhT5IkLXkppXq4y4Jd79GXGOzpBhyWKWn5MuxJkqQlp9LfR//JDvpOHqfvxDF6j77EUG8PYLiTpGGGPUmStGilWo2BzpP0neyg/2QH/SePM3Cyg6H+vpE2hjtJmpxhT5IkLQpZb93xerDLQl1/54nR8BZBcPrWB4Y7SZqcYU+SJOUqpURtaJChvl4qfb0M9Y9/rvT1Mth9iqG+3pHvTOytq58Id7iTpJkz7EmSpDlJtRqVgT4q/X1ThLg+Kv3Z+1Stnvb9iCCiQK12+mf21knS2TPsSZKkESMBrh7UKv19VAb6sl65/izYVQfqnw0MTHqOKBRGzjXttVIipdODniRpfhj2JElqYiklqoMDVAf6R4JbZczran9/9jzQP32AiwIRMWkv3GnXtFdOkhYFw54kSUtIqtWoDg5QGejPAtrIcz3EDfRT7e8bPdbfd9qCJsOiUCCA2gzCWUo1pjiNJGmRMuxJktQgwwuXVAb6pwhwo0Eu653rm7LnDYbnwAW1WoIZLGWSajUXPJGkJmbYkyRpHtSq1SyQDfZTHRigMvw80F8/1k+lPpyyVg92lcEBpusuKxQKJGY+LDKbA2d8kyRlDHuSJNWllEjVStbLNjhAdWCgHuDqIW5wgMqYY7X6c2Wwn1qlMuV5h1edTKRZzWebyfBKSZKmYtiTJDWVlBK1ShbYqkPDwWxwTGgbPj44GtiGRj+fPowFhULMqrdtuCZXnZQkLTTDniRpUcl616qjAawezGqDg6MhbfjYSJvBcaHtTEMZs4VJor7oyGyGPab6fDhJkhY/w54kaV7VqtUshNUfIyFtaDicDY6+rr+vjX0/NHjGXrPhYZEAtVSbdt7bZFyYRJK0HBj2JElAfUn/oUFqQ0NUK/XnMSEsex6a8D5rN5ugBnOfwzZSq8MiJUk6I8OeJC1hqVbL5qcNDVKrDGWPoUGqlaExYS07PhLQ6m2GQ9rIZ9WZhaeRIZBzDGpgWJMkaSEY9iRpAWWLhwyHsqGRUHbasZHj9d60se1GAtoQterUK0COly0swlmGNHAIpCRJS4VhT5KmMNxrVqsOB61KFrIqY8JZpTI+rNWPDYe4NOZ4tf7ZTA1vkM2cFhI57adxYRFJkpYZw56kppVSotLXQ9/xYwz2do+EtZGQNuZ1qo4/Vq1USDMc1pgJolAPZ/WNrc92c2s3yJYkSWfDsCepKWTBrpe+E0fpO36MvuNH6TtxjEp/32ltIwr1HjPOekjjmApIteTwRkmStGgY9iQtSUP1HrvhUNd3/Oi4YBcKewYkAAATXElEQVSFwpQhLhsSuVCVSpIkNYZhT9KiV6sM0X3kJfo6jtB34hj9J44x1Nc78vlkwW5+euskSZKWLsOepEVpqL+Xrhf+mc4Xvk/3Ky+MzJ8z2EmSJM2MYU/SopBSYqDzJJ0vPEfXi8/Re+wV4PRwZ7CTJEmaGcOepIZJtRq9HUfofOH7dL3wHANdpwBGFk8ZbiNJkqTZM+xJWlC1SoXulw9nAe/Ff6Yy0A/17QqGud2AJEnS2TPsScpdrVrh1PPf49Tzz9Lz8mFq1SoRBVKq99oZ7iRJkuadYU9Sbga7Ozn+v77LiWf/kcpA/7j5dyNBT5IkSbkw7EmaV6lWo+vlwxz/p6foeun50z6TJEnSwjDsSZoXlf4+jj/7j5z4X99lsKeLKBQaXZIkSdKyZtiTNGcpJXo7XuH4Pz3NqeefzXru6itp2osnSZLUWIY9SbNWHRri5HP/xIlD36XvZIeLrUiSJC1Chj1JZ5RSItWqDHbVF1z53jPUKkMj++G52IokSdLiY9iTmlhKiVStUq0MUhsaojY0RHVokNrQYPZcGft+aPR4/XX2+QDVoaHRYZlj9sRzPzxJkqTFy7AnLTIpJWqVoewxNES1/lyr1ANZZUxoq9QD2nCbervhz6uVoRkOqwyiEARZ82l76gx4kiRpCWlra+Oaa65hzZo1s/7uM888MzKSqRHK5TIbN26ktbV1Tt837ElnKdVq9R6wodGQNjaojXlfq1Sy0FapjPScjX42eo6ZikKBIEikrJdtzkEskWoJY5wkSWo211xzDVu3bmXlypWzDm7r1q1rWNhLKdHR0cHhw4fZunXrnM5h2NOykWo1atVKPXANh7CpX1crQ6TKUBbM6oEsVSvjAt244Y0zEBH1vzDqAe0sV6xMtZoBTZIkaRpr1qyZU9BrtIjgvPPO4+jRo3M+R65hLyJuAT4KFIHfSyn9+oTPfwn434AKcBT4mZTSc3nWpMUrWwQk6yVLI6GsUg9oQ6e/rlRI1bEhbXy7NOZYtVIh1aqzqicLZtlecfMRzEZ+RodBSpIkLailFvSGnW3duYW9iCgCHwfeAhwGDkbEIymlp8c0+1tgb0qpNyJ+DvhN4CfyqklzNxzEUrUy2jtWrZAq499nAWxCSJvYbrjN8PF6m2qlMvthiPWesmy22XCYmp+VIbNzzS4gSpIkSTNxyy238I1vfIPXve51/Nmf/Vku18izZ+9a4FBK6VmAiPgc8HZgJOyllP56TPtvAHfkWE9TyoYmVkdDWLU6GqyqY4NXdVxQS9Xq+HA24btpQqCrVucQxKj/a8RIGEtnXvxj1n8AWU+ZfWWSJEmaiU996lPzer4PfehDc/re3XffTW9vL5/85CfntZ6x8gx7G4Dnx7w/DFw3Tfv3An8+2QcRcSdwJ8DmzZvnq755N7zM/UhYqlZHQ1V1OGCNhq6J70deVya2yZ5TtXJ6+zkPLRxefbHeI0YaCU7zaXjREMOYJEmSlqODBw/y3ve+l29+85tUq1WuvfZaPv/5z/OmN72Jr3zlK7lee1Es0BIRdwB7gTdO9nlK6UHgQYC9e/cuutzwyCOP8O3vfOes53RFFMjW7oiRgDQf88Qm5+qLkiRJUt727dvH/v37uffee+nr6+OOO+5g9+7dC3LtPMPeC8CmMe831o+NExFvBn4VeGNKaSDHenKzdetWnn76aQYGzq78lGpuYSZJkiQ1mfvuu499+/ZRLpd54IEHFuy6hRzPfRDYERFbI6INeBfwyNgGEfEa4JPA/pTSkRxrydWVV15JuVxudBmSJEmSFqGOjg66u7vp6uqiv79/wa6bW9hLKVWA9wNfAr4LfCGl9FREfDgi9tebfQRYBfxfEfHtiHhkitNJkiRJ0pJ01113cf/993P77bdzzz33LNh1c52zl1J6FHh0wrH7xrx+c57XlyRJkqRGevjhh2ltbeW2226jWq1y44038uUvf5kPfehD/MM//APd3d1s3LiRT3/609x8883zeu1FsUCLJEmSJC2E973vfTNuu27durPe2PzAgQMcOHAAgGKxyOOPPw7AD/3QD53VeWcizzl7kiRJkqQGMexJkiRJUhMy7EmSJElSEzLsSZIkSVITMuxJkiRJUhMy7EmSJElSEzLsSZIkSdIC+va3v80NN9zArl272LNnD5///OdzuY777EmSJElaNr7RNYt987pePmOTH925btY1rFixgocffpgdO3bw4osv8trXvpabb76Zc845Z9bnmo49e5IkSZKUk4MHD7Jnzx76+/vp6elh165dDA4OsmPHDgDWr1/PhRdeyNGjR+f92vbsSZIkSVJO9u3bx/79+7n33nvp6+vjjjvuYPfu3SOff/Ob32RwcJDt27fP+7UNe5IkSZKUo/vuu499+/ZRLpd54IEHRo6/9NJL/ORP/iSf/exnKRTmf9ClwzglSZIkKUcdHR10d3fT1dVFf38/AJ2dndx666382q/9Gtdff30u1zXsSZIkSVKO7rrrLu6//35uv/127rnnHgYHB/mRH/kRDhw4wDvf+c7cruswTkmSJEnKycMPP0xrayu33XYb1WqVG2+8kc997nN89atfpaOjg4ceegiAhx56iKuvvnper23YkyRJkrRsXL86zbjtunXriJjFVg2TOHDgAAcOHACgWCzy+OOPjxzPm8M4JUmSJKkJGfYkSZIkqQkZ9iRJkiSpCRn2JEmSJKkJGfYkSZIkqQkZ9iRJkiSpCRn2JEmSJGkBPffcc1xzzTVcffXV7Nq1i0984hO5XMd99iRJkiQtG+u/smF+T3jbzPftG7Zu3Toee+wxSqUS3d3d7N69m/3797N+/fp5Lc2ePUmSJEnKycGDB9mzZw/9/f309PSwa9cunnnmGUqlEgADAwPUarVcrm3PniRJkiTlZN++fezfv597772Xvr4+7rjjDnbv3s3zzz/PrbfeyqFDh/jIRz4y7716YNiTJEmSpFzdd9997Nu3j3K5zAMPPADApk2bePLJJ3nxxRd5xzvewTvf+U4uuuiieb2uwzglSZIkKUcdHR10d3fT1dVFf3//uM/Wr1/P7t27+drXvjbv1zXsSZIkSVKO7rrrLu6//35uv/127rnnHg4fPkxfXx8AJ06c4Otf/zo7d+6c9+s6jFOSJEmScvLwww/T2trKbbfdRrVa5cYbb+Spp57i7rvvJiJIKfGBD3yAK6+8ct6vbdiTJEmStGy8eNMLM267bt06IuKsrnfgwAEOHDgAQLFY5PHHHwfg5ptvPqvzzoTDOCVJkiSpCRn2JEmSJKkJGfYkSZIkqQkZ9iRJkiSpCRn2JEmSJKkJGfYkSZIkqQkZ9iRJkiRpAT333HNcc801XH311ezatYtPfOITuVzHffYkSZIkLRsbPrVhXs+XPpRm/Z1169bx2GOPUSqV6O7uZvfu3ezfv5/169fPa2327EmSJElSTg4ePMiePXvo7++np6eHXbt28cwzz1AqlQAYGBigVqvlcu1ce/Yi4hbgo0AR+L2U0q9P+LwEPAy8FugAfiKl9P08a5IkSZKkhbJv3z7279/PvffeS19fH3fccQe7d+/m+eef59Zbb+XQoUN85CMfmfdePcixZy8iisDHgbcCVwDvjogrJjR7L3AipXQJ8B+B38irHkmSJElqhPvuu4+//Mu/5IknnuCXf/mXAdi0aRNPPvkkhw4d4rOf/SyvvPLKvF83z2Gc1wKHUkrPppQGgc8Bb5/Q5u3AZ+uv/xh4U0REjjVJkiRJ0oLq6Oigu7ubrq4u+vv7x322fv16du/ezde+9rV5v26ewzg3AM+PeX8YuG6qNimlSkScAs4DjuVYVy527txJX19fo8uQJEmSNEZbWxvt7e0NreGuu+7i/vvv53vf+x733HMPH/zgBznvvPNob2/nxIkTfP3rX+cXf/EX5/26S2I1zoi4E7gTYPPmzQ2uZnJvfetbG12CJEmSpAm++93vcu655zbs+g8//DCtra3cdtttVKtVbrzxRp566inuvvtuIoKUEh/4wAe48sor5/3aeYa9F4BNY95vrB+brM3hiGgB1pIt1DJOSulB4EGAvXv3zn5tU0mSJEliblslnI0DBw5w4MABAIrFIo8//jgAN998c+7XznPO3kFgR0RsjYg24F3AIxPaPAK8p/76ncCXU0qGOUmSJEk6S7n17NXn4L0f+BLZ1gufSSk9FREfBp5IKT0CfBr4g4g4BBwnC4SSJEmSpLOU65y9lNKjwKMTjt035nU/8GN51iBJkiRJy1GewzglSZIkqeGW6kyxs63bsCdJkiSpaZXLZTo6OpZc4Esp0dHRQblcnvM5lsTWC5IkSZI0Fxs3buTw4cMcPXq00aXMWrlcZuPGjXP+vmFPkiRJUtNqbW1l69atjS6jIRzGKUmSJElNyLAnSZIkSU3IsCdJkiRJTSiW2qo0EXEUeG7MobXAqRl+/Uxtz+bz84FjM6xjMZnNn99iutZczzWX7830OzNpN10b76/Fc62zOddsv7tQ99d0n3t/Ley1ltv9BUvzHvP+mr/23l+TW6h7zPurOe+vi1NKF5yxVUppST+AB+er7dl8DjzR6D+LvP/8FtO15nquuXxvpt+ZSbsz3EPeX4vkWmdzrtl+d6Hur+k+9/5a2Gstt/ur/tmSu8e8v+avvfdX/v/dF+o63l9L79EMwzj/2zy2PdvPl6KF/Jnm81pzPddcvjfT78yk3XRtvL8Wz7XO5lyz/e5C3V+zudZS4f01f+29v07n/TV/7b2/JrdQP5P31/K8v4AlOIxzsYqIJ1JKextdh5qT95fy5P2lvHmPKU/eX8rTUr+/mqFnb7F4sNEFqKl5fylP3l/Km/eY8uT9pTwt6fvLnj1JkiRJakL27EmSJElSEzLsSZIkSVITMuxJkiRJUhMy7C2AiChExK9FxMci4j2NrkfNJSJuioivRcQnIuKmRtej5hMRKyPiiYh4W6NrUXOJiMvrf3f9cUT8XKPrUfOJiHdExKci4vMR8cONrkfNJSK2RcSnI+KPG13LVAx7ZxARn4mIIxHx9xOO3xIR/xgRhyLig2c4zduBjcAQcDivWrX0zNP9lYBuoIz3l8aYp/sL4B7gC/lUqaVqPu6vlNJ3U0o/C/w48AN51qulZ57usS+mlN4H/CzwE3nWq6Vlnu6vZ1NK78230rPjapxnEBFvIPtF+uGU0u76sSLwDPAWsl+uDwLvBorAv59wip+pP06klD4ZEX+cUnrnQtWvxW2e7q9jKaVaRFwE/HZK6faFql+L2zzdX1cB55H9Y8KxlNKfLUz1Wuzm4/5KKR2JiP3AzwF/kFL6w4WqX4vffN1j9e/9B+C/ppT+vwUqX4vcPN9fi/b3+5ZGF7DYpZS+GhFbJhy+FjiUUnoWICI+B7w9pfTvgdOGOUXEYWCw/raaX7Vaaubj/hrjBFDKo04tTfP099dNwErgCqAvIh5NKdXyrFtLw3z9/ZVSegR4JCL+O2DY04h5+jssgF8H/tygp7Hm+XewRcuwNzcbgOfHvD8MXDdN+z8BPhYRrwe+mmdhagqzur8i4keBm4FzgN/NtzQ1gVndXymlXwWIiJ+i3ouca3Va6mb799dNwI+S/UPVo7lWpmYx29/Bfh54M7A2Ii5JKX0iz+K05M3277DzgF8DXhMRv1IPhYuKYW8BpJR6gUU9nldLV0rpT8j+QUHKTUrpoUbXoOaTUvoK8JUGl6EmllJ6AHig0XWoOaWUOsjmgy5aLtAyNy8Am8a831g/Js0H7y/lyftLefL+Ut68x5Snpru/DHtzcxDYERFbI6INeBfwSINrUvPw/lKevL+UJ+8v5c17THlquvvLsHcGEfFHwGPAzog4HBHvTSlVgPcDXwK+C3whpfRUI+vU0uT9pTx5fylP3l/Km/eY8rRc7i+3XpAkSZKkJmTPniRJkiQ1IcOeJEmSJDUhw54kSZIkNSHDniRJkiQ1IcOeJEmSJDUhw54kSZIkNSHDniRpWYmIakR8e8zjg3M4x96IeKD++qci4nfnv1JJks5OS6MLkCRpgfWllK4+mxOklJ4AnpineiRJyoU9e5IkARHx/Yj4zYj4u4j4ZkRcUj/+YxHx9xHxnYj4av3YTRHxZ5OcY0tEfDkinoyIv4qIzfXjD0XEAxHxNxHxbES8c2F/OknScmTYkyQtN+0ThnH+xJjPTqWUrgR+F/id+rH7gJtTSlcB+89w7o8Bn00p7QH+K/DAmM/WAa8D3gb8+nz8IJIkTcdhnJKk5Wa6YZx/NOb5P9Zf/0/goYj4AvAnZzj3DcCP1l//AfCbYz77YkqpBjwdERfNvmxJkmbHnj1Jkkalia9TSj8L3AtsAr4VEefN8dwDY17HHM8hSdKMGfYkSRr1E2OeHwOIiO0ppcdTSvcBR8lC31T+BnhX/fXtwNfyKlSSpDNxGKckablpj4hvj3n/Fyml4e0Xzo2IJ8l64d5dP/aRiNhB1hv3V8B3gDdOce6fB34/Iu4mC4Y/Pe/VS5I0Q5FSOnMrSZKaXER8H9ibUjrW6FokSZoPDuOUJEmSpCZkz54kSZIkNSF79iRJkiSpCRn2JEmSJKkJGfYkSZIkqQkZ9iRJkiSpCRn2JEmSJKkJGfYkSZIkqQn9/5rxMDv7CmfgAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1080x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plotting allocations (Note: the sum of all the portfolio allocations is 1, so we do not plot that!)\n",
    "plt.figure(figsize=(15, 5))\n",
    "y1 = np.array([item[0] for item in opt_x_vals])\n",
    "plt.plot(epsilons, y1, c = 'gray')\n",
    "plt.fill_between(epsilons, 0, y1, facecolor = 'gray', label = 'x1')\n",
    "\n",
    "y2 = np.array([item[1] for item in opt_x_vals])\n",
    "plt.plot(epsilons, y1+y2, c= 'lightblue')\n",
    "plt.fill_between(epsilons, y1, y1+y2, facecolor = 'lightblue', label = 'x2')\n",
    "\n",
    "y3 = np.array([item[2] for item in opt_x_vals])\n",
    "plt.plot(epsilons, y1+y2+y3, c= 'orange')\n",
    "plt.fill_between(epsilons, y1+y2, y1+y2+y3, facecolor = 'orange', label = 'x3')\n",
    "\n",
    "\n",
    "\n",
    "y4= np.array([item[3] for item in opt_x_vals])\n",
    "plt.plot(epsilons, y1+y2+y3+y4, c= 'green')\n",
    "plt.fill_between(epsilons, y1+y2+y3, y1+y2+y3+y4, facecolor = 'green', label = 'x3')\n",
    "\n",
    "plt.legend(loc = 'best')\n",
    "plt.xlabel('Epsilon')\n",
    "plt.ylabel('Allocation')\n",
    "plt.xscale('log')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe that as our risk tolerance $\\epsilon$ goes down we allocate more weight on the risk-free asset ($x_4$ shown in green).\n",
    "For higher values of epsilon most of our weight is on assets with higher returns $x_1 (gray), x_2(blue)$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lets perform a Monte Carlo simulation and see how often we incur a loss."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Just regenerating the optimal x again with 10e-4 as our epsilon\n",
    "optimal_x_part_1 = optimization(1e-4, r_hat, Sigma_root)[1]\n",
    "optimal_x_part_1 = np.array(optimal_x_part_1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Generating 10000 random samples. Note: could have used np.random.randn as well\n",
    "returns = []\n",
    "np.random.seed(777)\n",
    "for i in range(10000):\n",
    "    rand_r = np.random.multivariate_normal(r_hat, Sigma)\n",
    "    returns.append(np.dot(rand_r, optimal_x_part_1))\n",
    "returns = np.array(returns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Frequency')"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEKCAYAAAAIO8L1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAF95JREFUeJzt3X20XWV94PHvjwAB3whIhNvchICmWuoLxIjMOHZUCsU3wBERq5CxTNOOdLVddjriy4zOrHFG1jhFqa02HaqBliJSKSnFcRBfuqYLhARSXosEBJJLIhERFBBI+M0f+7l4uOyTu++9Z59z7r3fz1pnnb2f/XJ/2ffk/s7zPHs/T2QmkiRNtMegA5AkDScThCSplglCklTLBCFJqmWCkCTVMkFIkmqZICRJtUwQkqRaJghJUq09Bx3ATBx44IG5fPnyQYchSbPKxo0bf5iZiyfbb1YniOXLl7Nhw4ZBhyFJs0pE3NNkP5uYJEm1TBCSpFomCElSLROEJKmWCUKSVMsEIUmqZYKQJNUyQUiSapkgJEm1TBBSn42MLiMinvXac+G+teUjo8sGHbLmqVk91IY0G20f28IhH7r8WeX3nP22ruXSIFiDkIbdgr2sWWggrEFILRgZXcb2sS29OdmuJ61ZaCBMEFILujUjgX/YNXvYxCRJqmWCkCTVMkFIkmqZIKQGuj274J1EmsvspJYa2N2zC9JcZQ1CklTLBCFJqmUTkzQT5SlnaS4yQUgz4VPOmsNsYpIk1TJBSJJqmSAkSbVaTRARcXdE3BQRmyJiQyk7ICKujIg7yvv+pTwi4tyI2BwRN0bEyjZjkyTtXj9qEG/MzCMyc1VZPwu4KjNXAFeVdYA3AyvKaw3w+T7EJknqYhBNTCcC68ryOuCkjvLzs3INsCgiRgYQnySJ9hNEAv83IjZGxJpSdlBmbivL24GDyvISoHOGla2lTJI0AG0/B/GvMnMsIl4EXBkR/9y5MTMzInIqJyyJZg3AsmUOlCZJbWm1BpGZY+X9fuBS4CjgB+NNR+X9/rL7GLC04/DRUjbxnGszc1Vmrlq8eHGb4UvSvNZagoiI50bE88eXgeOAm4H1wOqy22rgsrK8Hji93M10NPBQR1OUpInKMB8OQa62tNnEdBBwaRmnZk/gwsz8PxFxHXBxRJwB3AOcUva/AngLsBl4FHh/i7FJs5/DfKhlrSWIzLwLeFVN+QPAMTXlCZzZVjxSp5HRZWwf2/Ks8oOXLGXb1nsHEJE0fBysT/NS1wmAPv0OR2eVChOE1MlmG+lpjsUkSaplgpAk1TJBSJJqmSAkSbVMEJKkWiYISVItE4QkqZYJQpJUywQhSaplgpAk1TJBSJJqmSAkSbVMEJKkWiYISVItE4QkqZYJQpJUywQhSaplgpAk1TJBSJJqmSAkSbVMEJKkWiYISVItE4TmtJHRZUTEs16SJrfnoAOQ2rR9bAuHfOjyZ5Xfc/bbBhCNNLtYg5DmmgV71daaRkaXDToyzTLWIKS5ZteT1prUE63XICJiQUTcEBGXl/VDI+K7EbE5Ir4cEXuX8oVlfXPZvrzt2CRJ3fWjien3gNs61s8GzsnMlwAPAmeU8jOAB0v5OWU/Sb1i05OmqNUmpogYBd4KfBL4YFS3j7wJ+PWyyzrgE8DngRPLMsAlwOciIjIz24xRmjdsetIUtV2D+AzwH4GnyvoLgR9n5s6yvhVYUpaXAFsAyvaHyv6SpAFoLUFExNuA+zNzY4/PuyYiNkTEhh07dvTy1JKkDm3WIF4HnBARdwMXUTUtfRZYFBHjTVujwFhZHgOWApTt+wEPTDxpZq7NzFWZuWrx4sUthq/ZxAfipN5rrQ8iMz8MfBggIt4A/IfMfG9EfAU4mSpprAYuK4esL+tXl+3ftP9BTflAnNR7g3hQ7kNUHdabqfoYzivl5wEvLOUfBM4aQGySpKIvD8pl5reBb5flu4Cjavb5GfCufsQjqUO5/XWig5csZdvWewcQkIaFT1JL8523v6oLx2KSJNUyQUiSapkgJEm1TBCSpFomCElSLROEJKmWCUKSVMsEIUmqZYKQJNUyQUiSajVKEBHxirYDkSQNl6Y1iD+NiGsj4gMRsV+rEUmShkKjBJGZrwfeSzWhz8aIuDAijm01MknSQDXug8jMO4CPUc3n8K+BcyPinyPi37QVnCRpcJr2QbwyIs4BbqOaOvTtmflLZfmcFuOTJA1I0/kg/hj438BHMvOx8cLMvC8iPtZKZJKkgWqaIN4KPJaZuwAiYg9gn8x8NDMvaC06SdLANO2D+Aawb8f6c0qZJGmOapog9snMn46vlOXntBOS1N3I6DIi4lkvSb3XtInpkYhYmZnXA0TEq4HHJjlG6rntY1ucP1nqk6YJ4veBr0TEfUAABwPvbi0qSdLANUoQmXldRLwMeGkpuj0zn2wvLEnSoDWtQQC8BlhejlkZEWTm+a1EJUkauEYJIiIuAF4MbAJ2leIETBCSNEc1rUGsAg7PzGwzGEnS8Gh6m+vNVB3TkqR5omkN4kDg1oi4Fnh8vDAzT2glKknSwDVNEJ+Y6okjYh/gH4CF5edckpkfj4hDgYuAFwIbgdMy84mIWEjVp/Fq4AHg3Zl591R/riSpN5rOB/Ed4G5gr7J8HXD9JIc9DrwpM18FHAEcHxFHA2cD52TmS4AHgTPK/mcAD5byc8p+kqQBaTrc928ClwB/VoqWAH+7u2OyMj48x17llVRDhF9SytcBJ5XlE8s6Zfsx4RgK0uAs2Kt2WJOR0WWDjkx90rSJ6UzgKOC7UE0eFBEvmuygiFhA1Yz0EuBPgDuBH2fmzrLLVqpkQ3nfUs6/MyIeomqG+mHDGCX10q4nHdZknmt6F9PjmfnE+EpE7ElVG9itzNyVmUcAo1QJ5mXTirJDRKyJiA0RsWHHjh0zPZ0kqYumCeI7EfERYN8yF/VXgL9r+kMy88fAt4B/ASwqCQaqxDFWlseo5rweT0D7UXVWTzzX2sxclZmrFi9e3DQESdIUNU0QZwE7gJuA3wKuoJqfuquIWBwRi8ryvsCxVFOWfgs4uey2GrisLK8v65Tt3/TBPEkanKaD9T0F/Hl5NTUCrCv9EHsAF2fm5RFxK3BRRPw34AbgvLL/ecAFEbEZ+BFw6hR+liSpx5qOxfR9avocMvOwbsdk5o3AkTXld1H1R0ws/xnwribxSJLaN5WxmMbtQ/WH/IDehyNJGhZNH5R7oOM1lpmfAd7acmySpAFq2sS0smN1D6oaxVTmkpAkzTJN/8j/r47lnVTDbpzS82gkSUOj6V1Mb2w7EEnScGnaxPTB3W3PzD/qTTiSpGExlbuYXkP1MBvA24FrgTvaCEqSNHhNE8QosDIzfwIQEZ8A/j4z39dWYJKGVBnldaKDlyxl29Z7BxCQ2tI0QRwEPNGx/kQpkzTfOMrrvNE0QZwPXBsRl5b1k/j53A2SpDmo6V1Mn4yIrwGvL0Xvz8wb2gtLkjRoTUdzBXgO8HBmfhbYWuaWlloxMrqsdjYzSf3T9DbXj1PdyfRS4ItU04f+JfC69kLTfLZ9bIvt3NKANa1BvAM4AXgEIDPvA57fVlCSpMFrmiCeKJP3JEBEPLe9kCRJw6Bpgrg4Iv6MarrQ3wS+wdQmD5IkzTJN72L6dJmL+mGqfoj/nJlXthqZJGmgJk0QZcrQb5QB+0wKkjRPTNrElJm7gKciYr8+xCNJGhJNn6T+KXBTRFxJuZMJIDN/t5WoJEkD1zRBfLW8JEnzxG4TREQsy8x7M9Nxl9SKkdFlbB/bMugwJNWYrAbxt8BKgIj4m8x8Z/shaT7xiWlpeE3WSd05+M1hbQYiSRoukyWI7LIsSc9UJhKa+BoZXTboyDRNkzUxvSoiHqaqSexblinrmZkvaDU6SbOHEwnNObtNEJm5oF+BSJKGy1Tmg5AkzSMmCElSrdYSREQsjYhvRcStEXFLRPxeKT8gIq6MiDvK+/6lPCLi3IjYHBE3RsTKtmKTJE2uzRrETuAPMvNw4GjgzIg4HDgLuCozVwBXlXWANwMrymsN8PkWY5MkTaK1BJGZ2zLz+rL8E+A2YAlwIjD+ZPY64KSyfCJwflauoZp7YqSt+CRJu9eXPoiIWA4cCXwXOCgzt5VN24GDyvISoHPMha2lbOK51kTEhojYsGPHjtZiVm+NjC6rvUde0vBqOljftEXE84C/AX4/Mx/u/KOQmRkRU3oALzPXAmsBVq1a5cN7s4RDakizT6s1iIjYiyo5/FVmjo8G+4PxpqPyfn8pHwOWdhw+WsokSQPQ5l1MAZwH3JaZf9SxaT2wuiyvBi7rKD+93M10NPBQR1OUJKnP2mxieh1wGtVEQ5tK2UeATwEXR8QZwD3AKWXbFcBbgM3Ao8D7W4xNkjSJ1hJEZv4/njkabKdjavZP4My24pEkTY1PUkuSapkgJEm1TBCSpFomCEnt6jKRkJMJDb/WH5STNM91mUgIfFBy2FmDkCTVMkFIkmqZICRJtUwQkqRaJghJUi0ThCSplglCPeXEQNLc4XMQ6iknBtKUlIfoJjp4yVK2bb13AAGpkwlC0uB0eYjOLxTDwSYmSVItE4QkqZYJQpJUywQhSaplgpAk1TJBSJJqmSAkSbVMEJqybk9L+8S0eqbLLHTOQNdfPiinKev2tDT4gJN6xAfohoI1CElSLROEJKmWCUKSVMsEIUmqZYKQJNVqLUFExF9ExP0RcXNH2QERcWVE3FHe9y/lERHnRsTmiLgxIla2FZckqZk2axBfAo6fUHYWcFVmrgCuKusAbwZWlNca4PMtxiVJaqC1BJGZ/wD8aELxicC6srwOOKmj/PysXAMsioiRtmKTJE2u330QB2XmtrK8HTioLC8BtnTst7WUSZIGZGCd1JmZQE71uIhYExEbImLDjh07WohM47oNqSFpfuj3UBs/iIiRzNxWmpDuL+VjwNKO/UZL2bNk5lpgLcCqVaumnGDUXLchNRzuQJof+l2DWA+sLsurgcs6yk8vdzMdDTzU0RQlSRqANm9z/WvgauClEbE1Is4APgUcGxF3AL9a1gGuAO4CNgN/DnygrbgkzWKO8tpXrTUxZeZ7umw6pmbfBM5sKxZJc4SjvPaVT1JLkmqZICRJtUwQkqRaJghJs5+d161wylFJs5+d162wBiFJqmWCkCTVMkFIkmqZIOSgfJJq2UktB+WTVMsahCSplgliHrEpSdJU2MQ0j9iUJGkqrEFIkmqZICRJtUwQkqRaJog5yM5oqXAQvxmxk3oOsjNaKroN4vfpd9R+aTp4yVK2bb23H5HNCiYISfOPo782YhOTJKmWCWIWs69BUptsYprF7GuQ1CZrEJKkWiYISVItE8QsYF+D1CddnpvYc+G+8/J5CvsgZgH7GqQ+2c3tr/Px/6A1CElSLROEJKnWUCWIiDg+Im6PiM0Rcdag42lLtz6Fbu2ckuaGbv/3h7UvY2j6ICJiAfAnwLHAVuC6iFifmbcONrLpGxldxvaxLbXbbOeU5oDSqT1RtzGdZlt/4tAkCOAoYHNm3gUQERcBJwJ9TRDd/qgv2Hsfdj3xs9pjdrdtNn0YJE3RFAcDnI5uf5P6MbDgMCWIJUDnVdgKvLatH9arb/e722YikOapqQ4G2KUmMm5Qf18iM1v/IU1ExMnA8Zn578r6acBrM/N3Juy3BlhTVl8K3N6nEA8EftinnzVVwxwbGN9MDHNsYHwzMcjYDsnMxZPtNEw1iDFgacf6aCl7hsxcC6ztV1DjImJDZq7q989tYphjA+ObiWGODYxvJoY5tnHDdBfTdcCKiDg0IvYGTgXWDzgmSZq3hqYGkZk7I+J3gK8DC4C/yMxbBhyWJM1bQ5MgADLzCuCKQcfRRd+btaZgmGMD45uJYY4NjG8mhjk2YIg6qSVJw2WY+iAkSUNkXiaIyYb0iIiFEfHlsv27EbG8Y9uHS/ntEfFrTc/Zj/gi4tiI2BgRN5X3N3Uc8+1yzk3l9aI+x7Y8Ih7r+Plf6Djm1SXmzRFxbszgCaMZxPfejtg2RcRTEXFE2daTa9cwvl+JiOsjYme59btz2+qIuKO8VneU9+T6TTe2iDgiIq6OiFsi4saIeHfHti9FxPc7rt0R04ltJvGVbbs6YljfUX5o+RxsLp+LvfsdX0S8ccJn72cRcVLZ1rPrNy2ZOa9eVB3gdwKHAXsD/wQcPmGfDwBfKMunAl8uy4eX/RcCh5bzLGhyzj7FdyTwC2X55cBYxzHfBlYN8NotB27uct5rgaOBAL4GvLnf8U3Y5xXAnb28dlOIbznwSuB84OSO8gOAu8r7/mV5/15dvxnG9ovAirL8C8A2YFFZ/1LnvoO4dmXbT7uc92Lg1LL8BeDfDyK+Cb/nHwHP6eX1m+5rPtYgnh7SIzOfAMaH9Oh0IrCuLF8CHFO+lZ0IXJSZj2fm94HN5XxNztl6fJl5Q2beV8pvAfaNiIXTjKOnsXU7YUSMAC/IzGuy+h9xPnDSgON7Tzm21yaNLzPvzswbgacmHPtrwJWZ+aPMfBC4Eji+h9dv2rFl5vcy846yfB9wPzDpQ1j9iq+b8nt/E9XnAKrPRWufvYbxnQx8LTMfnWYcPTUfE0TdkB5Luu2TmTuBh4AX7ubYJufsR3yd3glcn5mPd5R9sVRT/9M0myFmGtuhEXFDRHwnIl7fsf/WSc7Zr/jGvRv46wllM712TeOb6rG9un49+QxHxFFU36Dv7Cj+ZGl6OmcGX1hmGt8+EbEhIq4Zb76h+r3/uHwOpnPOXsY37lSe/dnrxfWblvmYIOa8iPhl4GzgtzqK35uZrwBeX16n9TmsbcCyzDwS+CBwYUS8oM8xTCoiXgs8mpk3dxQP+trNCqU2cwHw/swc/5b8YeBlwGuomk8+NKDwDsnqqeVfBz4TES8eUBxdlev3CqpnwcYN9PrNxwTRZEiPp/eJiD2B/YAHdnNso2FC+hAfETEKXAqcnplPf4vLzLHy/hPgQqoqcd9iK81yD5QYNlJ9w/zFsv/oJOdsPb6O7c/6Bteja9c0vqke26vrN6PPcEn2fw98NDOvGS/PzG1ZeRz4IoO5dp2/w7uo+pSOpPq9Lyqfgymfs5fxFacAl2bmk+MFPbx+0zIfE0STIT3WA+N3iZwMfLO0764HTo3qTphDgRVUHYS9HCZk2vFFxCKq/6RnZeY/ju8cEXtGxIFleS/gbcDNTN1MYlsc1ZwfRMRhVNfurszcBjwcEUeXppvTgcumEduM4itx7UH1n/Tp/oceXrum8XXzdeC4iNg/IvYHjgO+3sPrN+3Yyv6XAudn5iUTto2U96Bq3+/7tSvXbGFZPhB4HXBr+b1/i+pzANXnos3P3mTew4QvJz28ftMzqN7xQb6AtwDfo/oW+9FS9l+BE8ryPsBXqDqhrwUO6zj2o+W42+m4W6TunP2OD/gY8AiwqeP1IuC5wEbgRqrO688CC/oc2zvLz94EXA+8veOcq6g++HcCn6M8wDmA3+0bgGsmnK9n165hfK+har9+hOob7i0dx/5GiXszVTNOT6/fdGMD3gc8OeFzd0TZ9k3gphLfXwLP6/e1A/5lieGfyvsZHec8rHwONpfPxcIB/W6XU9U49phwzp5dv+m8fJJaklRrPjYxSZIaMEFIkmqZICRJtUwQkqRaJghJUi0ThNQhfj7q580R8Xfl2ZLd7b8oIj7Qr/ikfjJBSM/0WGYekZkvpxpV88xJ9l9ENULslIw/NCgNMxOE1N3VdAy4FhF/GBHXlYHT/ksp/hTw4lLr+J8R8YaIuLzjmM9FxL8ty3dHxNkRcT3wrqjmmTg7Iq6NiO+ND2AYEb9cyjaVn7Wib/9iqcNQzUktDYvyDf8Y4LyyfhzV8CBHUc27sD4ifgU4C3h5Zo5PLvSGSU79QGauLPv+NrBnZh4VEW8BPg78KvDbwGcz86/KsA3WNjQQJgjpmfaNiE1UNYfbqOZdgGrso+OAG8r686gSxr1TPP+XJ6x/tbxvpBpuAaqay0fLwItfzTLXgtRvNjFJz/RYqQ0cQlVTGO+DCOB/lP6JIzLzJZl5Xs3xO3nm/6t9Jmx/ZML6+Hwduyhf2DLzQuAE4DHgiuiYOlbqJxOEVCOrGb1+F/iDMhz014HfiIjnAUTEkqjmpv4J8PyOQ+8BDi8j/i6iaqaakjLa7V2ZeS7V6KKvnNm/Rpoem5ikLjLzhoi4EXhPZl4QEb8EXF2NvMxPgfdl5p0R8Y8RcTPVVJF/GBEXU42++X1+3iQ1FacAp0XEk8B24L/35B8kTZGjuUqSatnEJEmqZYKQJNUyQUiSapkgJEm1TBCSpFomCElSLROEJKmWCUKSVOv/A1lW2E029kfxAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the histogram of random returns \n",
    "plt.hist(returns, 50, ec='black')\n",
    "plt.xlabel('Returns')\n",
    "plt.ylabel('Frequency')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Mean of the returns: 0.09194666595656943\n"
     ]
    }
   ],
   "source": [
    "print(\"Mean of the returns:\", returns.mean())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Percentage of loss: 0.0007\n"
     ]
    }
   ],
   "source": [
    "print(\"Percentage of loss:\", (returns < 0).sum() / returns.size)"
   ]
  }
 ],
 "metadata": {
  "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.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}