equity-index-options-pricing-example.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "simple-prior",
   "metadata": {},
   "source": [
    "# Pricing equity index options\n",
    "In this notebook I price equity index options using volatility tracking, riskless yield curve construction, and options pricing code from this package. I ensure that all market and economic data required for pricing is obtained from [FRED](https://fred.stlouisfed.org) and [Yahoo-finance](https://finance.yahoo.com) using [pandas-datareader](https://pydata.github.io/pandas-datareader/). Since the lifetime of an option is usually measured using trading\n",
    "days rather than calendar days, I make use of [the pandas-market-calendars library](https://pandas-market-calendars.readthedocs.io/en/latest/index.html) to get correct trading calendars.\n",
    "\n",
    "This notebook can also be used as a tutorial to using the `pricing.curves`, `pricing.options`, `volatility.parameter_estimators`, and `volatility.volatility_trackers` packages in this repository."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "shaped-miniature",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Uncomment if running in Google Colaboratory, otherwise the import of the curves module in the cell below will fail\n",
    "#!git clone -l -s https://github.com/ilchen/options-pricing.git cloned-repo\n",
    "#%cd cloned-repo\n",
    "\n",
    "# Install the latest version of pandas-datareader and yfinance\n",
    "# !pip install pandas-datareader -U\n",
    "# !pip install yfinance -U\n",
    "# !pip install pandas-market-calendars -U"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "conscious-nebraska",
   "metadata": {},
   "source": [
    "## Constructing a riskless rate curve for USD\n",
    "Since I already covered it in a separate notebook, I will simply import it. The resultant YieldCurve object is assigned to the curve variable.\n",
    "\n",
    "I prepared two alternative rislless rate/yield curves &mdash; one based on US Treasury Yields (you can use it for pricing options maturing up to 30 years in the future) and another based on AMERIBOR rates (you can use it for pricing options maturing up to 2 years in the future). Since AMERIBOR term structure is no longer published by the Federal Reserve Bank of St. Louis, it cannot be used anymore unfortunately..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "square-algeria",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABJEAAAIzCAYAAACwSohtAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAAsTAAALEwEAmpwYAACTB0lEQVR4nOzdd3ydV2H/8c/R3pKt6W3LtrwSOztxpp0NZFCgjJIAbSFQVoHSUtoyCqVQOn6FUsouhB2gQFjZMdnDGU7ivWcsWd6SrH1+f9xrozi2ZXndR9Ln/Xrdl66e59x7v/fKBOvrc84TYoxIkiRJkiRJR5KV6QCSJEmSJElKPkskSZIkSZIk9csSSZIkSZIkSf2yRJIkSZIkSVK/LJEkSZIkSZLUL0skSZIkSZIk9csSSZIkDUkhhK+EED52lGPXhRCuPMy5eSGETSc2nSRJ0uBjiSRJUsKEEGIIYcpBxz4ZQvhen+//LoSwNoTQEkLYFEL48WGeq6XPrTeEsK/P928+2e/lZAohvD2EsDSEkN/nWGUIoSmEcG2M8V0xxk+fghwhhPD+EMILIYTW9M/jJyGE00/2a58qIYT8EMI3QwjrQwh7QwjPhhBecdCYK0IIy0IIbSGE+0MIE/qc+7cQwsr0Y5eFEN5ymNd5S/rP/9v7yXOk11p80J/77hDCr47wXK8PITySfq4Fhzj/tRDC8vT/ft52pFySJA11lkiSJA0yIYS3AjcDV8YYS4BzgHsPNTbGWLL/BmwAru9z7Pt9njPnVGTvTwgh+2jHxhi/AWwGPt7n8H8Cv40x3nGCox3JF4C/BN4PjAQagF8ArxroEyXl53AIOcBG4DKgHPgH4LYQwkSAEEIV8H/Ax0h9BguBvsVmK3B9+rFvBb4QQriw7wuEEEYAfwcsPlKQ/l4rxjirz5/50nTunxzhKXeQ+nPzucOcXwS8G3j6SLkkSRoOLJEkSRp8zgXujDGuBogxbo0xfm0gT7B/iVYI4SMhhK3A/4YQskIIfxtCWB1C2B5CuC2EMLLPY34SQtgaQtgdQngghDCrz7lXhhCWpGeabA4hfDh9/G0hhIcOeu0DM61CCN8OIfxPCOG3IYRW4EMhhMa+ZVII4TUhhEWHeStvB94dQjgjhHANcAXwwT7P/U99nue69AyaXemZJ7MP89kUph+7M4SwhNTnfbjPcSrwHuBNMcb7YowdMca2GOP3Y4yfS49Z0HdmzcGfSfrzeE8IYSWwMv15/NtBr/PLEMKH0vdHhxB+FkLYFlKz0d5/hHzlIYRb02PXhxD+IYSQ1TdHepbQzvRzveJQzxNjbI0xfjLGuC7G2Btj/DWwFjg7PeQ1wOIY409ijO3AJ4E5IYTp6cd/Isa4LP3Yx4EHgbkHvcxngS8CzYd7P0fzWge5FKgCfna4J4sx3hNjvA3Ycpjz/x1jvBdo7yeXJElDniWSJEmDz2PAW0IIfx1COGcgs3cOUkdqJscE4BbgfcCrSc02GQ3sBP67z/jfAVOBGlKzMr7f59w3gXfGGEuB04D7BpDjT4DPkJo18l/AduDqPudvBm491ANjjOtIzUT6FvAV4N0xxp0HjwshnJke806gEvgqcHvosxSuj08Ak9O3a0jNnDmcK4BNMcYnjjDmaLwaOB+YCfwQeEMIIaSzjyD1efwoXQD9itTsmDHp1/9AukA7lP8iNfunntTP9S3An/Y5fz6wnFTR8nngm/tf90hCCLWkZlztnzU0K50JSJVOwOr08YMfW0iqmFvc59h5pGbUfaW/1x7Ia5H62f0sPUaSJB0nSyRJkgaZGOP3SBU+1wC/B5pCCB85hqfqBT6Rnj2zD3gX8Pcxxk0xxg5SMzxet3+JVYzxWzHGvX3OzQkhlKefqwuYGUIoizHujDEOZOnPL2OMD6dnqbQD3wFuAkjPhLoG+MERHv+l9Os/G2P8xWHG3AJ8Ncb4eIyxJ8b4HaADuOAQY18PfCbGuCPGuJHU7JjDqQRePML5o/XZ9OvtIzVLJwKXpM+9Dng0xriFVPlSHWP8VIyxM8a4Bvg68MaDnzBdLr4R+Gj657YO+HdSpdx+62OMX48x9pD63EcBtUcKGkLIJVUgfifGuCx9uATYfdDQ3aSKwYN9hVQJdGefnF8G3htj7D3Saw/ktUIIRaQ+u28fxXNKkqSjYIkkSVLy9AC5Bx3LJVWUAJBeLnUlUEGq/Pn0EWajHM62dGmz3wTg5+nlXruApekstSGE7BDC59JL3fYA69KPqUp/fS3wSmB9COH3IYSDlyodycaDvv8ecH0IoZhUofNgjPGwRU2MMaazHmkvnQnAX+1/b+n3N47UjKuDjT4o0/ojPO92UsXL8Trweun38yPgTelDf8IfZn1NAEYf9D7+jkMXP1Wk/tz0zb+e1Aym/bb2ed229N2Sw4VMz4T6LtAJvLfPqRag7KDhZcDegx7/r6Rmqr0+/T4htd/QczHGxw7xeuNDn02yB/JapJa97SBVtO5/vq/0eb6/O9z7lCRJh2aJJElS8mwAJh50bBKHKDNijF0xxp8Az5H65Xwg4kHfbwReEWOs6HMriDFuJlVk3AhcSWp51P58IZ3jyRjjjaSWuv0CuC19vhUo2v8CIYS6/nKkX+9RUiXAzaRKi+O1kdTsor7vrSjG+MNDjH2RVMG03/gjPO+9wNgQwjlHGPOSz4DUMsKDHfyz+CGpWWATSC0527+nz0Zg7UHvozTG+MpDPGczqeJxQp9j40ltRj5g6WVu3yRVWL02xtjV5/RiYE6fscWklgP2XbL2j8ArgKtjjHv6PPYK4I/S+21tBS4E/j2E8KUY44aDNoc/qtdKeytwa5+yivQV+/Y/3z8fy+cgSdJwZokkSVLy/Bj4hxDC2JDa7PpKUle2+ikc2BD5VSGE0vT5V5DaD+bx43zdrwCfSRcXhBCqQwg3ps+Vklr+tZ1UIXLgF/AQQl4I4c0hhPJ0sbCH1FI5SC1bmpXe+LqA1DK4o3Er8DfA6aSuxHW8vg68K4Rwfkgp3v8ZHmLsbcBHQwgjQghjSS0dPKQY40pSS7F+GFKbleeFEApCCG8MIfxtetizwGtCCEUhtaH4n/cXNsb4DKkS6BukNlHflT71BLA3pDZEL0zPEDsthPCyzb/TS9RuI/UzLU3/XD9EaqbXsfgfYAapK/ztO+jcz4HTQgivTf+cP05qdtEygBDCR0kVkVfGGLcf9Ni3pZ/3jPRtIfCPwN8fJscRXyv9emOB+aSW6B1R+jMsIHUFuqz0zy+3z/m89PkA5KbP+3doSdKw5P8BSpKUPJ8CHgEeIrW59eeBN8cYX0if30NqCdMGYFf6/F/EGB96+VMNyBeA24G7Qgh7SW3gfX763K2kZkJtBpakz/V1M7AuvdTtXcCbAWKMK9Lv5x5gZfo9HY2fk15e12eZ1TGLMS4E3kFq/6SdwCpS5cWh/COp97oWuIv+Z0K9P/28/03q57Ea+CNSG2AD/D9Sy78aSZUa33/5UxzSD0jN/DqwH1S6GLqOVNmylj8UTeWHeDykCrBWYA2pz/4HpDYYH5B0AfXO9Otu7bMkbP/PeRupJY2fIfX5ns9L92n6Z1KzoFYdvJwsxrgrpq4wuDXGuJXUZ7Unxnjwvkf7P4P+XgtSfx4fjekrGPbjZmAfqZLskvT9r/c5f1f62IXA19L3Lz2K55UkacgJfWb4SpIkJUYIYTWpK77dk+kskiRJciaSJElKoBDCa0ntE3RfprNIkiQpJSfTASRJkvoKISwAZgI3H+Ul3yVJknQKuJxNkiRJkiRJ/XI5myRJkiRJkvpliSRJkiRJkqR+Ddo9kSoqKuKUKVMyHeMlWltbKS4uznSMl0liriRmgmTmSmImSGauJGaCZOZKYiZIZq4kZoJk5kpiJkhmriRmAnMNRBIzQTJzJTETJDNXEjNBMnMlMRMkM1cSM0EycyUxE5z6XE899VRzjLH6kCdjjIPy1tDQEJPm/vvvz3SEQ0piriRmijGZuZKYKcZk5kpiphiTmSuJmWJMZq4kZooxmbmSmCnGZOZKYqYYzTUQScwUYzJzJTFTjMnMlcRMMSYzVxIzxZjMXEnMFGMycyUxU4ynPhewMB6mi3E5myRJkiRJkvpliSRJkiRJkqR+WSJJkiRJkiSpX5ZIkiRJkiRJ6pclkiRJkiRJkvpliSRJkiRJkqR+WSJJkiRJkiSpX5ZIkiRJkiRJ6pclkiRJkiRJkvpliSRJkiRJkqR+WSJJkiRJkiSpX5ZIkiRJkiRJ6pclkiRJkiRJkvpliSRJkiRJkqR+WSJJkiRJkiSpX5ZIkiRJkiRJ6pclkiRJkiRJkvpliSRJkiRJkqR+WSJJkiRJkiSpX0O2RLp3aSOf+OULmY4hSZIkSZI0JAzZEumrv1/Ddx5dz77OnkxHkSRJkiRJGvSGZIm0s7WThet3ALBue2uG00iSJEmSJA1+Q7JEun95E70xdX+9JZIkSZIkSdJxG5Il0j1LGxlZnAfA2ua2DKeRJEmSJEka/IZcidTR3cPvl2/j2tPqqCrJY12zM5EkSZIkSZKO15ArkR5bs4PWzh6umlHLxMpi90SSJEmSJEk6AYZciXTPkkYKc7OZO7mSCZZIkiRJkiRJJ8SQKpFijNyztJFLG6ooyM1mUlURjXs6aOvsznQ0SZIkSZKkQW1IlUiLt+zhxd3tXDmjFoCJVcUArN/u5tqSJEmSJEnH46hKpBDCuhDC8yGEZ0MIC48w7twQQncI4XV9jr01hLAyfXtr+lh+COGOEMILIYR39xn7tRDCWcf6Zu5Z2kgIcPn0GgAmVqZKJDfXliRJkiRJOj45Axg7P8bYfLiTIYRs4F+Au/ocGwl8AjgHiMBTIYTbgUuAh4B/Bh4GvhxCmANkxxifHvC7SLtnaSNnjx9BZUk+8IeZSGvdF0mSJEmSJOm4nMjlbO8DfgY09Tl2DXB3jHFHjHEncDdwLdAFFAG5QEiP/TTwsWN98Rd37+OFzXu4cmbtgWMl+TlUleQ7E0mSJEmSJOk4HW2JFIG7QghPhRBuOfhkCGEM8EfA/xx0agywsc/3m9LH7gYmAo8BXwwh3AA8HWPcMrD4f3DP0lR3tX8/pP0mVhaxzj2RJEmSJEmSjkuIMfY/KIQxMcbNIYQaUgXQ+2KMD/Q5/xPg32OMj4UQvg38Osb40xDCh4GCGOM/pcd9DNgXY/y3Po/NBe4EbgT+ERgP3BpjvP0QOW4BbgGorq4++7bbbjtw7t8XttPU1svnLikkhHDg+Dee7+CF5h7+c37R0X8qx6ilpYWSkpKT/joDlcRcScwEycyVxEyQzFxJzATJzJXETJDMXEnMBMnMlcRMkMxcScwE5hqIJGaCZOZKYiZIZq4kZoJk5kpiJkhmriRmgmTmSmImOPW55s+f/1SM8ZxDnowxDugGfBL48EHH1gLr0rcWUkvaXg28Cfhqn3FfBd500GP/EngbcBWp/ZOygd/3l6OhoSHut7e9K079u9/GT/9qcTzYl+5bGSd85Nexpb3rZedOtPvvv/+kv8axSGKuJGaKMZm5kpgpxmTmSmKmGJOZK4mZYkxmriRmijGZuZKYKcZk5kpiphjNNRBJzBRjMnMlMVOMycyVxEwxJjNXEjPFmMxcScwUYzJzJTFTjKc+F7AwHqaL6Xc5WwihOIRQuv8+cDXwwkFF1KQY48QY40Tgp8C7Y4y/IDXD6OoQwogQwoj0Y+/s89wjgOuAW0ntkdRLaulcYX+5+npwxTY6e3pfsh/SfhMqUzOQ1rukTZIkSZIk6ZgdzZ5ItcBDIYRFwBPAb2KMd4QQ3hVCeNeRHhhj3EFqw+wn07dPpY/t93HgMzHGXlLl0iXA88B3B/Im7l7aSHlhLudMGPGycxMrU1doW+cV2iRJkiRJko5ZTn8DYoxrgDmHOP6Vw4x/20Hffwv41mHGfrDP/XZSM5UGpLunl/uXNXH59Bpysl/eiU2sSpVIa71CmyRJkiRJ0jE72quzJdbTG3axs63rZVdl268kP4fq0nzWOxNJkiRJkiTpmA36EumepY3kZgcubag67JiJlUWsa3ZPJEmSJEmSpGM1+EukJY1cUF9JaUHuYcdMrCxmrTORJEmSJEmSjtmgLpFWb2thTXMrVx3iqmx9TawqZtveDlo7uk9RMkmSJEmSpKFlUJdI9yxpBOCKw+yHtJ9XaJMkSZIkSTo+g7tEWtrIzFFljKkoPOK4iVVFAO6LJEmSJEmSdIwGbYnUE+Gp9Tu5sp+lbOBMJEmSJEmSpOM1aEukfd2R3ghX9bOUDaA4P4fq0nzWNVsiSZIkSZIkHYtBWyK1dUFtWT6njSk7qvGTKoudiSRJkiRJknSMBm2JtK87cuWMWkIIRzV+YlURa90TSZIkSZIk6ZgM2hIpwlHth7TfxKpimls6aOnoPnmhJEmSJEmShqhBWyIFYG595VGPP7C5tvsiSZIkSZIkDdigLZEKcwIFudlHPd4rtEmSJEmSJB27QVsileYd3V5I+02sKgJg/Xb3RZIkSZIkSRqoQVsiFeYMbHxRXg41pfmsdTmbJEmSJEnSgA3aEulYTKwqdk8kSZIkSZKkYzCsSqRJlcWsczmbJEmSJEnSgA2rEmlCVRHNLR3sbe/KdBRJkiRJkqRBZViVSJPSV2hzc21JkiRJkqSBGVYl0sSqVInk5tqSJEmSJEkDM6xKpAmVRQCs326JJEmSJEmSNBDDqkQqysuhtiyftc0uZ5MkSZIkSRqIYVUiAUysLGadM5EkSZIkSZIGZNiVSJOqil3OJkmSJEmSNEDDrkSaUFlMc0sne9u7Mh1FkiRJkiRp0Bh2JdKkqtTm2uvcF0mSJEmSJOmoDbsSaWJVMYD7IkmSJEmSJA3AsCuRJoxMl0jNlkiSJEmSJElHa9iVSIV52dSVFbDWmUiSJEmSJElHbdiVSAATq4pYv909kSRJkiRJko7W8CyRKotdziZJkiRJkjQAw7NEqipme2sne9q7Mh1FkiRJkiRpUBieJVKlm2tLkiRJkiQNxPAskaqKAFjnvkiSJEmSJElHZViWSBNGOhNJkiRJkiRpIIZliVSYl82o8gJLJEmSJEmSpKM0LEskSF+hbbslkiRJkiRJ0tEYviVSVZF7IkmSJEmSJB2l4VsiVRazo7WT3fu6Mh1FkiRJkiQp8YZviVSV2lx7vUvaJEmSJEmS+jV8S6TKVIm01s21JUmSJEmS+jVsS6QJlUUArGt2XyRJkiRJkqT+DNsSqSA3m9HlBS5nkyRJkiRJOgrDtkQCmFBZzFpLJEmSJEmSpH4N6xJpYlUx69wTSZIkSZIkqV/DukSaVFXEzrYudrd1ZTqKJEmSJElSog3rEmlC+gpt61zSJkmSJEmSdETDukSqr0qVSGtd0iZJkiRJknREw7pEGl9ZRFaANZZIkiRJkiRJRzSsS6T8nGzGjihizbaWTEeRJEmSJElKtGFdIgHUVxezZpszkSRJkiRJko5k2JdIk6qKWdvcSm9vzHQUSZIkSZKkxBr2JVJ9dQn7unpo3Nue6SiSJEmSJEmJNexLpMnpK7S5pE2SJEmSJOnwhn2JVF9dAuDm2pIkSZIkSUcw7Euk2rJ8ivKyWe1MJEmSJEmSpMMa9iVSCIFJVcWsabZEkiRJkiRJOpxhXyJBaknb2maXs0mSJEmSJB2OJRJQX1XMpp37aO/qyXQUSZIkSZKkRDqqEimEsC6E8HwI4dkQwsJDnH9zCOG59JhHQghz+py7NoSwPISwKoTwt32Ofz/9mH/uc+wfQgivPs73NGD11cXECOu3t53ql5YkSZIkSRoUcgYwdn6Msfkw59YCl8UYd4YQXgF8DTg/hJAN/DdwFbAJeDKEcHv6dffFGGeHEO4OIZQDRcD5McZ/OuZ3c4zqq/5whbZpdaWn+uUlSZIkSZIS74QsZ4sxPhJj3Jn+9jFgbPr+ecCqGOOaGGMn8CPgRqALKAwhZAG5QA/wKeATJyLPQE2qLgZwc21JkiRJkqTDONoSKQJ3hRCeCiHc0s/YPwd+l74/BtjY59wmYEyMcSmwDXga+BUwBciKMT591MlPoJL8HGrL8lmzzRJJkiRJkiTpUEKMsf9BIYyJMW4OIdQAdwPvizE+cIhx84EvAxfHGLeHEF4HXBtjfHv6/M2klqy996DH/Qp4J/CnwBzg7hjj1w/x/LcAtwBUV1effdtttw3s3R7Bvzyxj84e+NjcwmN+jpaWFkpKSk5YphMlibmSmAmSmSuJmSCZuZKYCZKZK4mZIJm5kpgJkpkriZkgmbmSmAnMNRBJzATJzJXETJDMXEnMBMnMlcRMkMxcScwEycyVxExw6nPNnz//qRjjOYc8GWMc0A34JPDhQxyfDawGGvocmwvc2ef7jwIfPehxN6afswH4VvrYnUDRkXI0NDTEE+mj//dcnP3JO2Nvb+8xP8f9999/4gKdQEnMlcRMMSYzVxIzxZjMXEnMFGMycyUxU4zJzJXETDEmM1cSM8WYzFxJzBSjuQYiiZliTGauJGaKMZm5kpgpxmTmSmKmGJOZK4mZYkxmriRmivHU5wIWxsN0Mf0uZwshFIcQSvffB64GXjhozHjg/4CbY4wr+px6EpgaQpgUQsgD3gjc3udxucAHgM8DhaSWzQFkA3n9ZTuR6quK2b2vix2tnafyZSVJkiRJkgaFo7k6Wy3w8xDC/vE/iDHeEUJ4F0CM8SvAx4FK4Mvpcd0xxnNijN0hhPeSmlmUTWqm0eI+z/0e4DsxxrYQwnNAUQjheeC3McZdJ+YtHp3J1ampYWubW6ksyT+VLy1JkiRJkpR4/ZZIMcY1pPYpOvj4V/rcfzvw9sM8/rfAbw9z7j/73I/Am/pNfJLU779C27ZWzpk4MlMxJEmSJEmSEulor8425I2pKCQ3O7C6uSXTUSRJkiRJkhLHEiktJzuLCZXFrNnWmukokiRJkiRJiWOJ1Ed9VTFrtjkTSZIkSZIk6WCWSH3UV5ewYUcb3T29mY4iSZIkSZKUKJZIfdRXF9PVE9m0c1+mo0iSJEmSJCWKJVIf9VXpK7S5ubYkSZIkSdJLWCL1UV9dAuDm2pIkSZIkSQexROpjZHEeFUW5rLZEkiRJkiRJeglLpIPUVxWz1uVskiRJkiRJL2GJdJBJVSUuZ5MkSZIkSTqIJdJB6quLadrbwd72rkxHkSRJkiRJSgxLpINMrk5doW1ts7ORJEmSJEmS9rNEOsj+K7RZIkmSJEmSJP2BJdJBJlQWEQJeoU2SJEmSJKkPS6SD5OdkM3ZEIWu2eYU2SZIkSZKk/SyRDqHeK7RJkiRJkiS9hCXSIdRXF7O2uZXe3pjpKJIkSZIkSYlgiXQI9dUl7OvqoXFve6ajSJIkSZIkJYIl0iHUVxUDuKRNkiRJkiQpzRLpEOqr95dIbq4tSZIkSZIElkiHVFdWQFFeNqudiSRJkiRJkgRYIh1SCIFJVcWsabZEkiRJkiRJAkukw6qvLmFts8vZJEmSJEmSwBLpsCZVFbNp5z7au3oyHUWSJEmSJCnjLJEOY3J1MTHC+u1tmY4iSZIkSZKUcZZIh1FfVQJ4hTZJkiRJkiSwRDqsSdXFAG6uLUmSJEmShCXSYZXk51BTms+abZZIkiRJkiRJlkhHUF9dzBqv0CZJkiRJkmSJdCT11SWs2dZKjDHTUSRJkiRJkjLKEukI6quK2b2vix2tnZmOIkmSJEmSlFGWSEdQn95ce62ba0uSJEmSpGHOEukI6qtKANxcW5IkSZIkDXuWSEcwdkQhudmB1W6uLUmSJEmShjlLpCPIyc5iQmWxM5EkSZIkSdKwZ4nUj/qqYtZscyaSJEmSJEka3iyR+jGpupgNO9ro7unNdBRJkiRJkqSMsUTqx+SqErp6Ipt27st0FEmSJEmSpIyxROpHfXUxAGvcXFuSJEmSJA1jlkj9qK8uAXBzbUmSJEmSNKxZIvVjRFEu5YW5rGm2RJIkSZIkScOXJVI/QghMqSlhVaPL2SRJkiRJ0vBliXQUGmpLWNG0lxhjpqNIkiRJkiRlhCXSUWioLWVXWxfb9nZkOookSZIkSVJGWCIdhWm1pQAsb9yb4SSSJEmSJEmZYYl0FBrq0iXSVkskSZIkSZI0PFkiHYWqknwqi/NY6ebakiRJkiRpmLJEOkoNtaUuZ5MkSZIkScOWJdJRmlZXysrGvfT2eoU2SZIkSZI0/FgiHaWG2lJaO3vYvGtfpqNIkiRJkiSdcpZIR2laXQkAK1zSJkmSJEmShiFLpKM0pSZ9hTZLJEmSJEmSNAxZIh2l8sJcRpUXsGKrJZIkSZIkSRp+LJEGIHWFtpZMx5AkSZIkSTrlLJEGYFpdKau3tdDd05vpKJIkSZIkSaeUJdIANNSW0tndy/odbZmOIkmSJEmSdEpZIg3AtNrU5truiyRJkiRJkoYbS6QBmFJTQgheoU2SJEmSJA0/lkgDUJiXzYSRRaywRJIkSZIkScPMUZVIIYR1IYTnQwjPhhAWHuJ8CCF8MYSwKoTwXAjhrD7n3hpCWJm+vTV9LD+EcEcI4YUQwrv7jP1a38cmUUNtKctdziZJkiRJkoaZgcxEmh9jPCPGeM4hzr0CmJq+3QL8D0AIYSTwCeB84DzgEyGEEcA1wEPAbODm9Ng5QHaM8eljfC+nxLS6UtZtb6OjuyfTUSRJkiRJkk6ZE7Wc7Ubg1pjyGFARQhhFqiy6O8a4I8a4E7gbuBboAoqAXCCkn+PTwMdOUJ6TZmptKT29kTXbWjMdRZIkSZIk6ZQ52hIpAneFEJ4KIdxyiPNjgI19vt+UPna443cDE4HHgC+GEG4Ano4xbhlY/FPvwBXa3BdJkiRJkiQNIyHG2P+gEMbEGDeHEGpIFUDvizE+0Of8r4HPxRgfSn9/L/ARYB5QEGP8p/TxjwH7Yoz/1uexucCdpGYz/SMwntSsptsPkeMWUsvlqK6uPvu22247pjd9PLp7I++8u41XTMrldQ15LznX0tJCSUnJKc/UnyTmSmImSGauJGaCZOZKYiZIZq4kZoJk5kpiJkhmriRmgmTmSmImMNdAJDETJDNXEjNBMnMlMRMkM1cSM0EycyUxEyQzVxIzwanPNX/+/KcOs5URxBgHdAM+CXz4oGNfBd7U5/vlwCjgTcBXDzcufewvgbcBV5HaPykb+H1/ORoaGmKmXPUfC+Kff/uJlx2///77T32Yo5DEXEnMFGMycyUxU4zJzJXETDEmM1cSM8WYzFxJzBRjMnMlMVOMycyVxEwxmmsgkpgpxmTmSmKmGJOZK4mZYkxmriRmijGZuZKYKcZk5kpiphhPfS5gYTxMF9PvcrYQQnEIoXT/feBq4IWDht0OvCV9lbYLgN0xxhdJzTC6OoQwIr2h9tXpY/ufewRwHXArqT2SekktnSvsL1cmNdSWstzlbJIkSZIkaRg5mj2RaoGHQgiLgCeA38QY7wghvCuE8K70mN8Ca4BVwNeBdwPEGHeQ2jD7yfTtU+lj+30c+EyMsZdUuXQJ8Dzw3eN+ZyfRtNpSNu7YR2tHd6ajSJIkSZIknRI5/Q2IMa4B5hzi+Ff63I/Aew7z+G8B3zrMuQ/2ud9OaqZS4jXUpTbXXtnUwhnjKjIbRpIkSZIk6RQ42quzqY8DV2jb6pI2SZIkSZI0PFgiHYNxI4soyM1ihfsiSZIkSZKkYcIS6RhkZwWm1JS4ubYkSZIkSRo2LJGOUUNtqTORJEmSJEnSsGGJdIym1ZbSuKeDXW2dmY4iSZIkSZJ00lkiHaP9V2hb0diS4SSSJEmSJEknnyXSMdp/hTb3RZIkSZIkScOBJdIxGlVeQGl+Diu2WiJJkiRJkqShzxLpGIUQaKgrdSaSJEmSJEkaFiyRjkNDbSkrG/cSY8x0FEmSJEmSpJPKEuk4TKstYWdbF9taOjIdRZIkSZIk6aSyRDoOB67QttUrtEmSJEmSpKHNEuk4NHiFNkmSJEmSNExYIh2HqpJ8KovzvEKbJEmSJEka8iyRjlNDrVdokyRJkiRJQ58l0nGaVpe6Qltvr1dokyRJkiRJQ5cl0nFqqC2ltbOHzbv2ZTqKJEmSJEnSSWOJdJym1ZUAsLLJJW2SJEmSJGnoskQ6TlP3X6Fta0uGk0iSJEmSJJ08lkjHqawgl9HlBaxwc21JkiRJkjSEWSKdAA11pSzfaokkSZIkSZKGLkukE6ChtpRV21ro8QptkiRJkiRpiLJEOgEaakvp7O6lqc0SSZIkSZIkDU2WSCfAtPTm2ptaejOcRJIkSZIk6eSwRDoBptSUEAJssUSSJEmSJElDlCXSCVCYl82EkUXORJIkSZIkSUOWJdIJ0lBbyua9lkiSJEmSJGloskQ6QabVlbK1LdLe1ZPpKJIkSZIkSSecJdIJMnNUGb0Rlm/dm+kokiRJkiRJJ5wl0gkye1wFAIs27cpoDkmSJEmSpJPBEukEGV1eQFkeLNq4O9NRJEmSJEmSTjhLpBMkhMCk8myecyaSJEmSJEkagiyRTqBJ5Vms2tZCS0d3pqNIkiRJkiSdUJZIJ9Ck8ixihBc2u6RNkiRJkiQNLZZIJ9Ck8mwAl7RJkiRJkqQhxxLpBCrLC4wdUciiTc5EkiRJkiRJQ4sl0gk2Z2wFizbuynQMSZIkSZKkE8oS6QSbPbacTTv3sb2lI9NRJEmSJEmSThhLpBNs9tgKAJ5zc21JkiRJkjSEWCKdYKePLScEeG6jJZIkSZIkSRo6LJFOsJL8HCZXl3iFNkmSJEmSNKRYIp0Es8eWs2jTbmKMmY4iSZIkSZJ0QlginQRnjKuguaWDLbvbMx1FkiRJkiTphLBEOgkObK69cVdGc0iSJEmSJJ0olkgnwYxRpeRmBxZtcnNtSZIkSZI0NFginQT5OdlMrytzc21JkiRJkjRkWCKdJLPHlvP8pt309rq5tiRJkiRJGvwskU6SOWMr2NvRzdrtrZmOIkmSJEmSdNwskU6S2ePKAVjk5tqSJEmSJGkIsEQ6SaZUl1CYm81zbq4tSZIkSZKGAEukkyQnO4vTx5SzyM21JUmSJEnSEGCJdBLNHlvOki176OrpzXQUSZIkSZKk42KJdBLNHldBR3cvy7fuzXQUSZIkSZKk42KJdBLNGZvaXNt9kSRJkiRJ0mBniXQSjR9ZREVRLs+5L5IkSZIkSRrkLJFOohACp48p59mNuzIdRZIkSZIk6bhYIp1kc8ZWsLKphX2dPZmOIkmSJEmSdMwskU6y2WPL6emNLN7ivkiSJEmSJGnwskQ6yc4YVwHAIjfXliRJkiRJg5gl0klWU1ZAXVmBm2tLkiRJkqRB7ahLpBBCdgjhmRDCrw9xbnwI4f70+edCCK/sc+6jIYRVIYTlIYRr0seqQwgPhRBeCCG8us/YX4YQRh/ne0qc2WPLec6ZSJIkSZIkaRAbyEykvwSWHubcPwC3xRjPBN4IfBkghDAz/f0s4FrgyyGEbOBNwFeA84APpMdeDzwTY9wy8LeRbHPGVbC2uZXdbV2ZjiJJkiRJknRMjqpECiGMBV4FfOMwQyJQlr5fDuwvgm4EfhRj7IgxrgVWkSqOuoAiIB/oCSHkkCqTPn8M7yHxZo8tB+C5zbsyG0SSJEmSJOkYhRhj/4NC+CnwWaAU+HCM8bqDzo8C7gJGAMXAlTHGp0IIXwIeizF+Lz3um8DvgLuBHwC1wEdIzVTaE2P8dj85bgFuAaiurj77tttuO/p3egq0tLRQUlLysuOtXZH33NvGa6fmcv3kvMTkyqQkZoJk5kpiJkhmriRmgmTmSmImSGauJGaCZOZKYiZIZq4kZgJzDUQSM0EycyUxEyQzVxIzQTJzJTETJDNXEjNBMnMlMROc+lzz589/KsZ4ziFPxhiPeAOuA76cvj8P+PUhxnwI+Kv0/bnAElKznL4E3NRn3DeB1x302BGkSqUS4OvAT4G5/eVqaGiISXP//fcf9txln78vvuM7T566MH0cKVemJDFTjMnMlcRMMSYzVxIzxZjMXEnMFGMycyUxU4zJzJXETDEmM1cSM8VoroFIYqYYk5kriZliTGauJGaKMZm5kpgpxmTmSmKmGJOZK4mZYjz1uYCF8TBdzNEsZ7sIuCGEsA74EXB5COF7B435c+C2dCn1KFAAVAGbgXF9xo1NH+vrY8BnSO2T9BDwVuCTR5FrUJk9tsLNtSVJkiRJ0qDVb4kUY/xojHFsjHEiqU2y74sx3nTQsA3AFQAhhBmkSqRtwO3AG0MI+SGEScBU4In9DwohTAXGxhgXkNojqZfU/kqFx/m+EmfOuAq27mmnaU97pqNIkiRJkiQN2ECuzvYSIYRPhRBuSH/7V8A7QgiLgB8Cb0vPglpMaobSEuAO4D0xxp4+T/MZ4O/T938I/AXwJPCFY82VVHPSm2svcjaSJEmSJEkahHIGMjg9Y2hB+v7H+xxfQmrZ26Ee8xlSZdGhzr2+z/0m4MKB5BlMZo0uJzsrsGjjLq6aWZvpOJIkSZIkSQNyzDORNDCFedlMrSlh0aZdmY4iSZIkSZI0YJZIp9CcsRU8v3n3/qvSSZIkSZIkDRqWSKfQ7HHl7GrrYsOOtkxHkSRJkiRJGhBLpFNoztgKAJ7duCujOSRJkiRJkgbKEukUmlZXSlFeNgvX7cx0FEmSJEmSpAGxRDqFcrOzOGfiSB5dsz3TUSRJkiRJkgbEEukUm1tfyaqmFpr2tmc6iiRJkiRJ0lGzRDrF5k6uBOCxNTsynESSJEmSJOnoWSKdYqeNLqMkP4dHV7ukTZIkSZIkDR6WSKdYTnYW500ayWPuiyRJkiRJkgYRS6QMmFtfydrmVrbudl8kSZIkSZI0OFgiZcD+fZEeXdOc4SSSJEmSJElHxxIpA2aMKqOswH2RJEmSJEnS4GGJlAHZWYHz6yt51H2RJEmSJEnSIGGJlCFz6yvZuGMfm3a2ZTqKJEmSJElSvyyRMuTAvkguaZMkSZIkSYOAJVKGTKstZURRrkvaJEmSJEnSoGCJlCFZWYEL6it5bPV2YoyZjiNJkiRJknRElkgZNHdyJVt2t7Nhh/siSZIkSZKkZLNEyqC59e6LJEmSJEmSBgdLpAyaUlNCVUm++yJJkiRJkqTEs0TKoBACF9SP5FH3RZIkSZIkSQlniZRhcydX0rS3gzXNrZmOIkmSJEmSdFiWSBnmvkiSJEmSJGkwsETKsElVxdSWuS+SJEmSJElKNkukDAshMLe+ksfXuC+SJEmSJElKLkukBJg7uZLmlk5WNrVkOookSZIkSdIhWSIlwNz6KsB9kSRJkiRJUnJZIiXAuJGFjKkotESSJEmSJEmJZYmUACEELqiv5LG12+ntdV8kSZIkSZKUPJZICTF3ciW72rpYtnVvpqNIkiRJkiS9jCVSQsydXAnAo2tc0iZJkiRJkpLHEikhxlQUMn5kkfsiSZIkSZKkRLJESpALJ1fy+Nrt9LgvkiRJkiRJShhLpASZO7mSve3dLNmyJ9NRJEmSJEmSXsISKUHm1u/fF6k5w0kkSZIkSZJeyhIpQWrKCqivLnZfJEmSJEmSlDiWSAkzt76SJ9ftpLunN9NRJEmSJEmSDrBESpi5kytp6ejm+c27Mx1FkiRJkiTpAEukhLngwL5ILmmTJEmSJEnJYYmUMFUl+TTUlrgvkiRJkiRJShRLpAS6ZGo1j6/ZQWtHd6ajSJIkSZIkAZZIiXTljFo6e3p5cGVzpqNIkiRJkiQBlkiJdM7EEZQX5nLP0sZMR5EkSZIkSQIskRIpNzuL+dOquW9ZEz29MdNxJEmSJEmSLJGS6sqZtexo7eSZDTszHUWSJEmSJMkSKakubagmNztwt0vaJEmSJElSAlgiJVRZQS4X1FdyzxJLJEmSJEmSlHmWSAl25YxaVm9rZc22lkxHkSRJkiRJw5wlUoJdMaMGgHuXNmU4iSRJkiRJGu4skRJs7IgiZowqc18kSZIkSZKUcZZICXfVjBoWrtvBztbOTEeRJEmSJEnDmCVSwl05s5beCPcvd0mbJEmSJEnKHEukhDttdDm1Zfnc45I2SZIkSZKUQZZICZeVFbh8ei2/X76Nju6eTMeRJEmSJEnDlCXSIHDVzBpaO3t4bM2OTEeRJEmSJEnDlCXSIHDh5CoKc7O5Z4lL2iRJkiRJUmZYIg0CBbnZXDK1inuWNhJjzHQcSZIkSZI0DFkiDRJXzqzlxd3tLN6yJ9NRJEmSJEnSMHTUJVIIITuE8EwI4deHOf/6EMKSEMLiEMIP+hx/awhhZfr21vSx/BDCHSGEF0II7+4z9mshhLOO5w0NVZdPryEEvEqbJEmSJEnKiIHMRPpLYOmhToQQpgIfBS6KMc4CPpA+PhL4BHA+cB7wiRDCCOAa4CFgNnBzeuwcIDvG+PQxvZMhrqokn7PGj7BEkiRJkiRJGXFUJVIIYSzwKuAbhxnyDuC/Y4w7AWKMTenj1wB3xxh3pM/dDVwLdAFFQC4Q0mM/DXzsWN7EcHHljFpe2LyHF3fvy3QUSZIkSZI0zISj2ag5hPBT4LNAKfDhGON1B53/BbACuAjIBj4ZY7wjhPBhoCDG+E/pcR8D9gH/CdwKzAD+FWgBzooxfrKfHLcAtwBUV1effdtttx3t+zwlWlpaKCkpOWnPv6Wll797aB9vmZnH5eNzE5PrWCQxEyQzVxIzQTJzJTETJDNXEjNBMnMlMRMkM1cSM0EycyUxE5hrIJKYCZKZK4mZIJm5kpgJkpkriZkgmbmSmAmSmSuJmeDU55o/f/5TMcZzDnkyxnjEG3Ad8OX0/XnArw8x5tfAz0nNLJoEbAQqgA8D/9Bn3MdIlVB9H5sL3EeqoPoP4KfADf3lamhoiElz//33n9Tn7+3tjZd9/r74lm8+PqDHnexcxyKJmWJMZq4kZooxmbmSmCnGZOZKYqYYk5kriZliTGauJGaKMZm5kpgpRnMNRBIzxZjMXEnMFGMycyUxU4zJzJXETDEmM1cSM8WYzFxJzBTjqc8FLIyH6WKOZjnbRcANIYR1wI+Ay0MI3ztozCbg9hhjV4xxLalZSVOBzcC4PuPGpo/19W5Ss5IuAHYDbwD+6ihyDTshBK6cUcujq7fT0tGd6TiSJEmSJGkY6bdEijF+NMY4NsY4EXgjcF+M8aaDhv2C1CwlQghVQAOwBrgTuDqEMCK9ofbV6WOkx44gNdPpVlJ7JPUCESg8rnc1hF05s5bOnl4eXLEt01EkSZIkSdIwMpCrs71ECOFTIYQb0t/eCWwPISwB7gf+Osa4Pca4g9SG2U+mb59KH9vv48BnYoy96ee4BHge+O6x5hrqzpkwgvLCXO72Km2SJEmSJOkUyhnI4BjjAmBB+v7H+xyPwIfSt4Mf8y3gW4d5vg/2ud9OaqaSjiAnO4vLp9dw/7Imenoj2Vmh/wdJkiRJkiQdp2OeiaTMuXJGLTvbunh6w85MR5EkSZIkScOEJdIgdGlDFbnZgXuWuKRNkiRJkiSdGpZIg1BpQS4X1Fe6L5IkSZIkSTplBrQnkpLjqpm1fPyXi1nVtJcpNaWZjiNJkiRJB8QY2dfVQ2tHD60d3bR2drOvs4eO7l46unvo6Oqls6eXjq709929B27dPb30xEhvb6SnF3p6U98fuN8LW7d28KumRWQFCAECgawsgEAIkBUgKwRys7PIzc4iLyeLvOyDv88iNydQkJNNUX4ORXnZFOZmU5SXTVFeDoV5qfu52c69kPazRBqkrplVxyduX8yvFr3IB6+yRJIkSZJ04nT39LKzrYudbZ2s2NlDz9JGdu/rYve+Lvbs6059bd//fRctHd20daYLo45u2rp6iPHYXjsnK5CVFcgO4Q/3999C6mt7ew/r2rYTY6Q3QiQSI/RGgNSx3hjp7ol09vTS2d17zJ9FbnagOD+HsoJcygtzKStM3S8ryKW8KJeyghzKClPnNjb3UL1lN1Ul+YwoyiMvxwJKQ4sl0iBVW1bA3PpKbl+0hQ9cOZUQvEqbJEmSpEOLMbKnvZtte9tp2tPBtpYOtu3tYEdr58tvbZ3saut66RM8vvAl3xbnZR8oTsoKcqkrK6AoP4fivGyK01/7fl+Ul5rpk5+TRX5u+mv6fl52Fvm5qe/zsrOO6nebBQsWMG/evAG9/+7eSFdPL13d6WIpXS61d/XQ1tnDvs6eAzOm2jp7aOvsTn9NlWN72lOF2Z72bpr2tBwo0tq7XlpQ/dvChw7cLyvIobIkn8riPEYW51FZkk91aT61ZfnUlRVQm75VFueR5ZW3NQhYIg1iN54xmo/87Hme37yb2WMrMh1HkiRJUgbs6+zhxd37eHF3Oy/ubmfr7n007umgaW87TXtTZdG2vR10HGI2Tk5WYERx3oGSY8boMiqL8xhRlEdlSerr+pVLueT8s9OzcHIpLcgZdEu8QgjkppezkXdin7uju4e97d3sauvivocfZ/zUmWxv7WR7S6qUa25JlXXrt7fx9IadbG/tfNksrZysQE1pPrXlBdSWFlBXXsDYEYWMHVHE2BGFjBtRRHlR7okNLh0DS6RB7NpZo/jYLxbzy2e3WCJJkiRJQ1Bvb2RbSwebdraxccc+Nu/ax5ZdfyiMXty97+WzhoCKolyqS/KpKcvn3IkjqSlNzYCpLs2nprTgwP2ygpx+Z/4s2LmCOeMqTtI7HPzyc7LJL8mmqiSfTSOymXfaqCOO7+7pZVtLB417Oti6u52mve1s3d1O454OGve0s3pbCw+taqalo/sljystyGHsiCLG9SmXJlYVMamqhLEjCgddsafByRJpECsvyuWyadX8+rkt/N0rZ5Dt9EdJkiRp0NnV1sm67W089mI3i+9fxaad+9i0s41NO/exeec+OnteOoNoZHEedWUFjKko4OwJFYwqL2RUecGBr3XlBRTkZmfo3ag/OdlZ6Z9VIYw79JgYI3v2dbNxZ9uBPwsbd6S+rtveykOrmmnr7PnDc2YFxo8sYlJVcepWnfpaX1VCPNbNqaRDsEQa5G48YzR3L2nk8bXbuXByVabjSJIkSTpIjKnZROu3t6Vvrazb/7W5lT3tfWacLFpOZXEeY0cUMnN0GVfPqn3JkqYxFYUU5lkQDXUhBMqLcikvKue0MeUvOx9jZEdrqnxc29zK2uYW1ja3smZbKw+vbn7JPk0F2TB9ycNMqy1lam0J0+pKmVZbSnVpvnvrasAskQa5K6bXUpyXze3PbrFEkiRJkjJoX2cPa5tbWb2thTXbWlnT3MLqbS2s3dZKa59ZI9lZgTEVhUyoLOLGM8YwobKIiZXFNK5ezB9dcylFef6apiMLIaQ27C7J5+wJI15yrrc3snVPe6pUam7l908voy03m3uXNfLjhRsPjKsoyqWhppSGuhKm1ZYyc3Q5M0aV+udPR+SfjkGuMC+bq2fV8dvnX+Qfb5xFfo7/KiFJkiSdTDtbO1nZ1MKKxr2samo5UBpt3rXvwJgQYHR5IfXVxZxzzkgmVRUfKIvGHGb/mgVNS/0FXsctKyswuqKQ0RWFXDSlinHta5k37wIAmls6WNG4lxVb97KiqYUVW/fyy2e3sDc9Gy4rQH11CaeNLmPW6HJmjUl9LS90U2+l+F+oIeCGM0bz82c288CKZq6aWZvpOJIkSdKQsLO1kxWNe1nZ1MLKxr2saGxhZVMLzS0dB8YU52VTX13CORNH8IbqcdRXp/ahmVRV7LIzJU5VST5VJfkvWcUSY2TL7naWbNnDC5t3s3jLbh5fu4NfPLvlwJhxIws5bXQ5p48t58xxI5g9tpzifOuE4cif+hBw8ZQqRhbn8ctnN1siSZIkSQPU0d3DqqYWlr24l+WNe1n64h6Wbd3Ltr0vLYum1JYyf1o1Dem9ZabWljK6vMB9ZTSohZBaXjmmovAlv09ub+lg8ZY9vLBlN4s3p77+7oWtQGrGUkNtKWeOH8GZ4ys4c1wFk6tLyPJiT0OeJdIQkJudxatOH8VPntpIS0c3JTbCkiRJ0svEGGna28GSLXv4zZpO/u/FZ1i2dQ+rt7XS05u6glVeThZTa0q4dGo10+pSRVGDZZGGocqSfC5tqObShuoDx3a2dvLspl08s2EXz2zYyW+e28IPn9gAQGlBDmeMSxVK502q5KwJFS7PHIL8iQ4RN5wxmu8+tp67l2zlj84cm+k4kiRJUkb19EbWNreyeMtulry4hyVbUrftrZ0Hxoyp2Mn0ulKumlnL9LoyZowqZWJlMTmH2K9IEowozmP+tBrmT6sBUpt4r2lu5ZkNO3l2Y6pc+tL9q+i9bxU5WYHTxpRz/qSRnDdpJOdMGEl5kXsrDXaWSEPE2eNHMKaikNuf3WKJJEmSpGGls7uXFY17eX7z7vSeLntYtnXPgcuc52YHGmpLuXx6DbNGlzFjVBnNq5/jVVfNz3ByaXDLygpMqSlhSk0Jf3zOOAD2tnfx9IZdPLF2O0+s3cH/PryOrz6whhBgWm1pulSqpKczZji9joUl0hCRlRW4fs5ovv7gGra3dFBZkp/pSJIkSdIJ193Ty8qmFp7ftJvnNu/i+U27Wbp1L53dqcKotCCHmaPKeNN545k1upyZo8qYUlNCXs5LZxct2ODSNOlkKC3I5bKGai5LL4Nr7+ph0cZdPLF2B0+s28FPntrEdx5dD8BXlj7IJVOruGhKFedOHOlm9IOAJdIQcuMZo/nK71fz2xe2cvMFEzIdR5IkSTouMUbWbW/j2Y07WbRxN89t2sXiLXvo2F8Y5ecwa0wZb7twIqePKWf22HLGjyxy7yIpQQpyszm/vpLz6ysB6Orp5fnNu/ne3U+ypSeHbz28lq8+sIa87CzOnjCCi9Ol0uljysl2o+7EsUQaQqbXlTK1poTbn91siSRJkqRBZ3tLB4s27eLZDbt4dtNuFm3cxe59XQAU5WVz2uhybrpgArPHlnP6mHImVhZ7NShpkMnNzuKs8SPYMzmPefPm0tbZzZPrdvLwqmYeXNnMv965nH+9czllBTlcPLWKedNqmDetmprSgkxHF5ZIQ0oIgRvPGM2/3bWCzbv2MaaiMNORJEmSpEPq7O5lyYt7eHr9Tp7ZuItnN+5k4459wB8uH/7K0+uYM7aCM8ZXMLWm1FkJ0hBUlJfzkuVvzS0dPLJ6Ow+t3MbvV2zjt89vBWD22HLmTavh8uk1zB5TboGcIZZIQ8wNc8bwb3et4FeLtvCuyyZnOo4kSZIEQNPedp5en7os+NMbdvLcpt0HlqWNKi/gzPEV3HT+BM4YV8HpY8u9NLg0TFWV5HPDnNHcMGc0MUaWvLiHBcu3cd+yJr5030q+eO9KKovzuGxaNZdPr+GSqdWUF3rVt1PF/zIPMeMrizhjXAW/fNYSSZIkSZnR0xtZvnUvT63fwW8WtfOxJ+47MMsoLzuLWWPKuOmCCZw9YQRnjR9BXbnLVCS9XAiBWaPLmTW6nPfMn8LO1k4eWJkqlO5b1sT/Pb2ZnKzABfWVXDOrlqtn1VFb5n9PTiZLpCHoxjNG84+/WsLKxr1MrS3NdBxJkiQNcfs6e1i0aRcL1+3gyXU7eXr9TvZ2dANQkR+YO7Wct1wwkbMmjGDW6DIKcr0Ck6SBG1Gcx41njOHGM8bQ0xt5duNO7lnaxJ2Lt/KxXy7mY79czBnjKrj2tDqumVXHpKriTEceciyRhqBXzR7Fp3+9hNsXbeGvrp6W6TiSJEkaYna0dvLkuh0HSqPFW3bT1RMBmFZbyg1njObciSM5e8IIVi16nPnzz85wYklDTXZW4OwJIzl7wkg+cu10VjXt5Y4XtnLn4kY+97tlfO53y2ioLeGaWalCadboMq/ceAJYIg1BNaUFXDSlil8+u4UPXdWQ6TiSJEka5Br3tPP42h08sXY7T6zdwYrGFgDycrI4Y2wF77iknnMmjuDs8SMpL3rp3iSr/aVN0ikwpaaU915eynsvn8rmXfu4a/FW7ly8lf++fxX/dd8qxo8s4vo5o7h+zmim1ZZaKB0jS6Qh6vo5o/mbnz7Hsxt3ZTqKJEmSBpEYI5t27uOxNanC6Il1O1i/vQ2Akvwczp4wghvPGMN5k0Yye2w5+TkuTZOULGMqCvnTiybxpxdNYntLB/csbeTXz73IV36/hv++fzVTa0q4fs5orps9ivrqkkzHHVQskYaoa0+r4x9+8QK3L9rCZW6LJEmSpMPYXxo9uno7j61J3bbsbgegoiiX8yaO5OYLJnD+pEpmjColJzsrw4kl6ehVluTzhnPH84Zzx9Pc0sHvXtjKrxZt4f/ds4L/uHsFp40p4/rZo3nV7FGZjjooWCINUWUFuVw+rYZfLXqRSy7yX4ckSZL0Bxt3tKULox08tmY7m3elrpxWWZzHBfWVvKt+JOdPqmRqTQlZWS75kDQ0VJXkc/MFE7j5ggm8uHsfv3nuRX713It89nfL+OzvljG1IosthRt41exRlBfm9v+Ew5Al0hB2wxmjuWPxVpZuL+DyTIeRJElSxjTuaeeR1c08smo7j67ZzqadqdJoRFEuF9RX8s7L6rmgPlUauU+IpOFgVHkhb7+knrdfUs+G7W386rktfP/hFfzdz5/nk79azNUza3ntWWO5ZGqVMzD7sEQawi6fXkNpfg6PbOnmPZkOI0mSpFNmV1snj63ZzsOrtvPI6mZWb2sFUsvTzp80krdfPIm5k6ucaSRJwPjKIt4zfwoz2Ujl1DP52VObuH3RFn793ItUl+bz6jNG85qzxjJjVFmmo2acJdIQVpCbzXVzRvOzhRvY095FWYHT8SRJkoaits5unli7gx8t6+Tfnn+QxVv2ECMU5WVz3qSRvOHccVw4uYqZo8osjSTpMEIIzB5bweyxFfz9q2Zy//ImfvbUJv734XV8/cG1zBxVxmvPHssfnTmGkcV5mY6bEZZIQ9wbzh3HD5/YwK8XvcifnD8+03EkSZJ0AnT39PL85t08tLKZh1Y18/SGnXT1RHICnD2xlA9c0cBFUyqZPbaCvByXYUjSQOXlZHHNrDqumVXHjtZObn92Mz97ejOf/vUS/uV3y7j2tDredN54LqgfOayWAVsiDXFzxpYztiTw44UbLZEkSZIGqRgj67a38dDKbTy0qplHVm9nb3s3ADNHlfFnF03iwilVtG98gWuumJvhtJI0tIwszuNtF03ibRdNYtnWPfzoiY3839OpJW+Tqop547njeO3ZY6kqyc901JPOEmmICyFwydhcfrhsF8u37mVaXWmmI0mSJOko7Grr5OFV23lw5TYeXNl84ApqYyoKeeVpo7h4ahUXTq6kss8vLQu2DJ9/DZekTJheV8Ynb5jF375iOr99/kV++MQGPvu7ZfzbXcu5elYdf3LeeObWVw7ZpcOWSMPAhaNz+OnKLn785EY+fv3MTMeRJEnSIXT19PLsxl08uGIbD6xs5rlNu+iNUFqQw4WTK3nXvMlcPKWKiZVFw2rphCQlUUFuNq85ayyvOWssKxv38sMnNvJ/z2ziN8+9yPiRRfzJ+eN547njqCgaWnsnWSINA6V5gatm1vLzZzbxkVdMIz8nO9ORJEmSBKzf3soDK5t5cMU2Hl29nb0d3WQFOGNcBe+7fCqXNlQxZ2yFl5eWpASbWlvKx6+fyd9cO407F2/l+49v4HO/W8b/u3sFrz5jDG+9cCIzRw+NK7tZIg0Trz9nHL99fiv3Lm3ilaePynQcSZKkYamts5vH1+xgwfImfr9iG+u2twGpJWrXzRnNpVOruHByFeVFXlVXkgabgtxsbjxjDDeeMYalL+7h1kfX8/NnNvHjhRs5b+JI3nrhRK6eVUvuIP6HAUukYeKSqdWMLi/gx09utESSJEk6RWKMrGxq4ffLt/H7Fdt4Yt0OOrt7KcjN4sLJVbztwolc2lDNpKpil6hJ0hAyY1QZn33N6fzttdO5beFGbn1sHe/5wdPUlRVw0wXjeeN54wflRtyWSMNEdlbgdWeP5b/uX8WWXfsYXVGY6UiSJElD0t72Lh5e1cyCdHH04u52AKbWlPCWCyZw2bRqzp04koJctxiQpKGuvCiXd1xaz59dPIn7lzXxnUfX8W93reCL967iujmjeMcl9cwYNXiWulkiDSOvO3scX7xvFT99ahPvv2JqpuNIkiQNCftnG92/rIn7lzexcN1Ounsjpfk5XDSlivdfUc2lDdWM8R/xJGnYys4KXDmzlitn1rKqqYVbH13HT5/axP89vZlLplZxy6X1XDylKvGzUi2RhpHxlUVcOLmS2xZu5L3zpwzZSw5KkiSdbO3dkbuXNHL/8iYWLGtiS3q20fS6Ut5+ST3zplVz9oQRg3rfC0nSyTGlpoRP3XgaH7qqge8/voFvP7KOm7/5BNPrSrnl0nqumz2avJxk/v+HJdIw84Zzx/GXP3qWx9Zs58IpVZmOI0mSNGisbW7lvmVN3L+sicdWt9EdF1Kcl81FU6p43xVTuayh2i0DJElHraIoj/fMn8LbL5nEL5/dwtcfWMOHblvE5+9Yzp9dPJE3njeesoJkXWjBEmmYuWZWHWUFOfx44UZLJEmSpCPo7O7libU7UsXR8ibWNrcCqX9BvnJCDm+58mzOmTgysf9aLEkaHPJzsnn9OeN43Vlj+f2KbXztgTX882+X8cV7V/Gm88YxM7s30xEPsEQaZgpys3n1mWP40ZMb+VRbl5ePlSRJ6qNpTzv3L2/ivmVNPLSymdbOHvJysphbX8nbLpzI5dNrGDeyiAULFvgPcpKkEyorKzB/eg3zp9fw/KbdfP3BNXzr4XWEGHly3/P8xWWTGTeyKKMZLZGGodefM45bH13PLxdt5i1zJ2Y6jiRJUsb09kZe2LKbe5emiqPnN+8GYFR5ATeeOYbLp9Vw4ZRKivL8a7Mk6dQ5fWw5X3zTmfz1NdP4+A8f5KcLN/HjJzdy4xmjefe8KUypKclILv/fcBg6bUw5s0aX8eMnN1oiSZKkYaets5uHVjZz37JUcdS0t4MQ4KzxI/jra6ZxxYwaptWWJv4KOZKkoW/cyCLeOiufz950AV97YA0/eGI9P39mM688bRTvnj+ZWaPLT2keS6Rh6vXnjOMTty/mhc27OW3Mqf1DJ0mSdKpt3rWP+5Y2cu+yJh5ZvZ3O7l5K83O4tKGay6fXMG9aNZUl+ZmOKUnSIdWVF/Dx62fy7vmT+dZDa7n10fX85vkXuWJ6De+5fApnjR9xSnJYIg1Trz5jDJ/57VJ+snCjJZIkSRpyensjz2/ezb1LG7l7aRNLX9wDwITKIm46fwJXzKjhXDfFliQNMlUl+fzNtdN556WT+c6j6/jWw2t5zZcf4eIpVXzwqqmcPWHkSX19S6Rhqrwol2tn1fGLZ7fw0VfOoCA3O9ORJEmSjkt7Vw8Pr2rmnqWN3Ls0tUwtK8DZE0bw0VdM54oZtUyuLnaZmiRp0CsvyuX9V0zlzy+exPceW8/XHljDa//nUS5tqOaDV07lzJM0M8kSaRh7w7njuH3RFu5cvJUbzxiT6TiSJEkDtm1vB/cubeSepU08tGob7V29lOTncGlDFVfOqGXetBpGFudlOqYkSSdFcX4O77xsMjfPncCtj67nq79fzR99+REun17DB69s4PSxJ3blkSXSMDa3vpJxIwu5beFGSyRJkjQoxBhZ1dTCXUsauWdpI89u3EWMMKaikDecM44rZ9Zy/qRKl6lJkoaVorwc3nXZZG66YALfeWQdX3tgDdd/6SGunFHLB6+aesI24LZEGsaysgJ/fPY4/uPuFWzc0ca4kUWZjiRJkvQyPb2Rx9Zs554ljdy9tJH129sAmD22nA9d2cCVM2uZXufV1CRJKsnP4T3zp/CWuRP434fX8Y0H1/CqLzZy7aw6PnDVVKbXlR3X81siDXOvPXss/++eFfzkqU186KqGTMeRJEkCoLWjmwdWbOPuJY3c+UIbrV2PkZedxYVTKnnHJfVcOaOWuvKCTMeUJCmRSgtSeya99cKJfOuhtXzrobXcsXgrrz5jNH919bRjnkRiiTTMjako5JKp1fxk4Ubef/kUcrKd+i1JkjKjaW879y5t4u4ljTy0qpnO7l4qinKZU53DzfNnc0lDNSX5/vVVkqSjVV6YywevauBPL5rIVx9Yw7ceWstvn9/KzXMn8N75UxgxwH0D/X9hcdP547nlu09xz9JGrj1tVKbjSJKkYWT1thbuWtzI3Uu28kx6f6NxIwu56fwJXDWzlnMnjuChBx9g3un+HUWSpGNVUZTHR66dzlvmTuD/3b2C/314Lbc9uZF3zZvMn100icK8o7tiuyWSuGJGLWNHFPK/D6+zRJIkSSdVb2/kmY07uWtJI3cvaWTNtlYATh9TzgevbOAq9zeSJOmkGVVeyOdfN4e3X1LP5+9Yxr/euZzvPrqeD141ldeeNbbf1UmWSCI7K/DWuRP5zG+XsmTLHmaOPr6NtiRJkvpq7+rh4VXN3LW4kXuXNdLc0klOVmDu5EreduFErpxRy+iKwkzHlCRp2GioLeUbbz2XJ9bu4LO/W8pHfvY833hwLX9z7fQjPs4SSQC8/pzUVdq+88g6/uV1szMdR5IkDXI7Wzu5d1kTdy/ZygMrmtnX1UNJfg7zplVz1cxa5k2robwwN9MxJUka1s6bNJL/+4sLuXPxVj5/x3LecevCI44/6hIphJANLAQ2xxivO8yY1wI/Bc6NMS5MH/so8OdAD/D+GOOdIYRq4OdABfAPMcZfpMf+EviLGOOWo82lE6O8KJfXnDWGnz61iY+8YjojB7i5liRJ0sYdbdy1pJG7Fm/lyXU76I1QV1bAa88ew9Uz67igvpK8HC/iIUlSkoQQuPa0UVwxo5YfP7mRm//l8GMHMhPpL4GlwCHXOoUQStNjHu9zbCbwRmAWMBq4J4TQALwJ+Arwf8BvgV+EEK4HnrFAypy3XTiR7z++gR89uYF3z5uS6TiSJCnhYowsb9zLnS80cteSrSzesgeAhtoS3j1vClfPquX0MeXubyRJ0iCQm53FTRdM4OYjjDmqEimEMBZ4FfAZ4EOHGfZp4F+Av+5z7EbgRzHGDmBtCGEVcB7QBRQB+UBPCCEH+ABw/dHk0ckxtbaUi6dU8d1H13PLJfX9bqglSZKGn97eyNMbUhtj37l4K+u3txECnDmugo++YjpXz6pjUlVxpmNKkqSTIMQY+x8Uwk+BzwKlwIcPXs4WQjgL+PsY42tDCAvSYxaGEL4EPBZj/F563DeB3wF3Az8AaoGPkJqptCfG+O1+ctwC3AJQXV199m233TaAt3rytbS0UFJSkukYLzOQXM80dfOFpzt4zxn5nFt38rbMGgqf1amSxEyQzFxJzATJzJXETJDMXEnMBMnMlcRMkMxcScwEh8/V3RtZsr2Hpxt7eLqphz2dkewAMyuzOasmmzNrsqkoOHn/+JTEzyuJmSCZuZKYCZKZK4mZIJm5kpgJkpkriZkgmbmSmAlOfa758+c/FWM855AnY4xHvAHXAV9O358H/Pqg81nAAmBi+vsFwDnp+18Cbuoz9pvA6w56/AhSpVIJ8HVSeyrN7S9XQ0NDTJr7778/0xEOaSC5unt648X/cm/84/955OQFikPjszpVkpgpxmTmSmKmGJOZK4mZYkxmriRmijGZuZKYKcZk5kpiphhfmqulvSv+etGW+L4fPB1P+/gdccJHfh1nfux38d3feyr+4plNcfe+zozkSookZooxmbmSmCnGZOZKYqYYk5kriZliTGauJGaKMZm5kpgpxlOfC1gYD9PFHM1Uk4uAG0IIrwQKgLIQwvdijDelz5cCpwEL0uvd64DbQwg3AJuBcX2ea2z6WF8fI7VM7k3AQ+kS6f+Aa44im06w7KzAW+dO5J9+s5QXNu/mtDHlmY4kSZJOkb2dkdsWbuSuxVt5YGUznd29jCzO4xWn13HNrDoumlJFQW52pmNKkqQM6bdEijF+FPgoQAhhHqmlajf1Ob8bqNr//UHL2fYBPwgh/AepjbWnAk/0GTsVGBtjXBBCmAO0AxEoPO53pmP2x+eM49/vWsF3HlnHv/7xnEzHkSRJJ9GWXfu4a/FW7li8lcfXtBF5jjEVhbz5/PFcM6uOcyaMcJ9ESZIEDOzqbC8RQvgUqSlOtx9uTIxxcQjhNmAJ0A28J8bY02fIZ4C/T9//IfAL4G+Bjx9rLh2/8sJcXnv2GG5buIm/fcV0KkvyMx1JkiSdQKu3tXDHC1u5a/FWFm3aDcDUmhKuq8/lna86n1mjy7yimiRJepkBlUgxxgWk9jwixnjIoifGOO+g7z9Dqiw61NjX97nfBFw4kDw6ed46dyLfe2wDP3pyI++ZPyXTcSRJ0nGIMfLC5j3cmZ5xtKqpBYA54yr4m2uncc2sOiZXl7BgwQKXskuSpMM6eZff0qA2tbaUS6ZW8d1H13PLpfXkOo1dkqRBpac3snDdDu5c3Midi7eyedc+srMC500cyc0XTODqWbWMKncHAUmSdPQskXRYb7twIn/+nYXcuXgr180enek4kiSpH53dvTyyupk7F2/lrsWNbG/tJC8ni0umVPGXV07lyhm1jCzOy3RMSZI0SFki6bDmT6thQmUR3354nSWSJEkJ1dbZze+Xb+OOxVu5b2kTezu6Kc7LZv70Gq49rY5502ooyfevfJIk6fj5NwodVlZW4C1zJ/LpXy/h+U27OX2seyRIkpQEu9u6uGdpI3cs3soDK7bR0d3LiKJcXnF6HdeeVseFk6soyM3OdExJkjTEWCLpiP74nLH8+13L+fYj6/j318/JdBxJkoatpj3t3LUktb/Ro6u3090bqSsr4I3njuPa00Zx7sQR5LiHoSRJOokskXREZQW5vO7ssfzoiY189JXTqSrJz3QkSZKGjQ3b2w5cUe3pDTuJESZVFfOOS+u5ZlYds8eUk5UVMh1TkiQNE5ZI6tdb5k7k1kfX88PHN/C+K6ZmOo4kSUNWjJEVjS2p4uiFrSx5cQ8AM0eV8cErG7j2tDqm1pQQgsWRJEk69SyR1K8pNSVcMrWK7z62nnfNm0yuU+UlSTphensjizbt4s7FqaVqa5tbCQHOHj+Cv3/lDK49rY5xI4syHVOSJMkSSUfnTy+ayJ99eyG/WrSF15w1NtNxJEka1Hp6I4+sbubOF7Zy5+JGtu5pJycrMHdyJX9+8SSunllLTVlBpmNKkiS9hCWSjsq8hhqm15XypftXceMZY8h2/wVJkgakvauHh1Y2c+firfzuuTZauh6nIDeLS6dW8zenTeOK6bWUF+VmOqYkSdJhWSLpqGRlBd53+VTe84On+c3zL3LDnNGZjiRJUuLtbe/i/uXbuPOFrSxY3kRrZw+l+TmcVpXN266Yw6UN1RTl+dcxSZI0OPi3Fh21V6Q38/zSfSu57vRRXg1GkqRDaG7p4J4lqf2NHl61nc6eXqpK8rnhjDFce1odc+sreeShB5h32qhMR5UkSRoQSyQdtayswHsvn8Jf/uhZ7ly8lVec7l9+JUkC2LSz7cDG2AvX7aA3wtgRhbxl7gSuOa2Os8aPcCm4JEka9CyRNCDXzR7NF+5ZyRfuXck1s+qcjSRJGpZijKxsauGOF7Zy5+KtLN6yB4BptaW8d/4UrjmtjpmjygjB/5+UJElDhyWSBiQ7PRvpQ7ct4p6ljVw9qy7TkSRJOiV6eyPPbtrFnYu3ctfiRtY2twJw1vgKPvqK6Vwzq46JVcUZTilJknTyWCJpwG6YM5ov3LuSL963kqtm1vqvrJKkIauzu5fH127nzsVbuXtJI417OsjJCsydXMmfXTyJq2fWUltWkOmYkiRJp4QlkgYsJzuL98yfwt/89DnuX97E5dNrMx1JkqQTpq2zm98v38adi7dy77Im9rZ3U5ibzaUNVVwzq44rptdSXpSb6ZiSJEmnnCWSjskfnTmGL967ki/eu4r502qcjSRJGtR2tHZyz9JG7lq8lQdXNtPR3UtFUS7XzKrj6pm1XDK1msK87EzHlCRJyihLJB2T3Ows3j1vCn/38+d5cGUzlzZUZzqSJEkDsnFHG3cvaeSuJVt5Ym3qimqjywt403njuXpWLedNHElOdlamY0qSJCWGJZKO2WvPHsOX7ktdqe2SqVXORpIkJVqMkaUv7uWuJamNsZe8+Icrqr173hSumVXHaWO8opokSdLhWCLpmOXnZPMX8ybzsV8u5tHV27lwSlWmI0mS9BLdPb0s29HDA79awl1LtrJp5z5CgHMmjODvXzmDq2bWekU1SZKko2SJpOPyx+eM40v3r+IL9660RJIkJcK+zh4eWLmNuxY3ct+yRna2dZGXs56Lp1Tx3vlTuGJGLdWl+ZmOKUmSNOhYIum4FORm867LJvOPv1rC42u2c359ZaYjSZKGoeaWDu5b2sRdSxp5cOU2Orp7KSvI4fLpNYxmO+/+o3mU5PvXHkmSpOPh36Z03N503nj++/7V/Nd9qyyRJEmnzLrmVu5aspW7lzSycP1OYoQxFYWpjbFn1nLupJHkZmexYMECCyRJkqQTwL9R6bgV5Gbzzkvr+cxvl/LU+h2cPWFkpiNJkoag3t7IMxt3cc/SRu5e0siqphYAZo4q4/2XT+WqmbXMGu3G2JIkSSeLJZJOiDdfMJ7/+f1qvnjvKr7zZ+dlOo4kaYjY19nDQ6uauWdJI/cua6S5pZOcrMD59SN58/njuXJGLeNGFmU6piRJ0rBgiaQToigvh3dcUs+/3LGMZzfu4oxxFZmOJEkapLbt7eD+Zan9jR5atY32rl5K83O4bFo1V82sZV5DDeVFuZmOKUmSNOxYIumEuXnuBL76wGr+/a7l3Ppn57mcQJJ0VGKMrGhs4Z6ljdyztJFnN+4iRhhdXsAbzhnHlTNrOX9SJXk5WZmOKkmSNKxZIumEKcnP4f2XT+VTv17CPUubuGpmbaYjSZISqrO7lyfW7jhQHG3auQ+A2WPL+cAVDVwxo8b9jSRJkhLGEkkn1M1zJ/DDJzbw6V8v4ZKpVRTkZmc6kiQpIXa2dvL7Fdu4Z2kjv1++jb0d3eTnZHHRlCrePW8KV8yoobasINMxJUmSdBiWSDqhcrOz+MT1s7jpm4/zzYfW8p75UzIdSZKUITFGVjW1cO+yJu5b2sTC9TvojVBVks8rTx/FlTNruWhKJUV5/nVEkiRpMPBvbTrhLp5axbWz6vjSfat4zVljGFVemOlIkqRTpKO7hyfW7uDepU3ct6yJDTvaAJg1uoz3zp/C5TNqmT2mnKwsl6lJkiQNNpZIOin+/lUzuG95E5/73TK+8MYzMx1HknQSNe1tZ8Hybdz2TDvvue9uWjt7yM/J4uIpVbzzsnoun17jPyhIkiQNAZZIOinGjSziXZfW88X7VnHTBRM4d+LITEeSJJ0gvb2R5zfv5r5lTdy/vInnNu0GYER+4MYzx3HF9BounFxFYZ774kmSJA0llkg6af5i3hR++tQmPvHLxfzqfReT7dIFSRq09rZ38eDKZu5b1sSC5dtobukgBDhzXAUfvrqBy6fX0rj8KebPPz3TUSVJknSSWCLppCnMy+ajr5zB+374DD9+ciN/cv74TEeSJB2lGCPLG/eyYPk2FixvYuG6nXT3RsoKcrhsWg2XT6/msoYaRhbnHXhM0wr/sUCSJGkos0TSSXXd7FF897H1/Oudy3jV6aMoL8rNdCRJ0mHsbe/i4VXNLFi+jd+v2MaLu9sBmF5Xyp9fMokrptdy1vgKcrKzMpxUkiRJmWCJpJMqhMAnr5/Fdf/1IP/vnhV88oZZmY4kSUqLMbL0xb0sWJFaovb0+tRso9L8HC6eWsUHrkzNNqorL8h0VEmSJCWAJZJOupmjy3jz+RP47mPredN545lWV5rpSJI0bO1q6+ShVc38Pj3bqGlvBwAzR5Vxy6X1zJtWw5njK8h1tpEkSZIOYomkU+JDVzVw+6It/OOvFvP9t5+f6TiSNGz0pK+kliqNmnh24y56I5QX5nLJ1CoubahmXkM1NWXONpIkSdKRWSLplBhRnMeHr27gY79czB0vbKUw04EkaQhr3NPOgyubeWDFNh5cuY2dbV2EAHPGVvC+y6dy2bRq5oyt8KqZkiRJGhBLJJ0ybzpvPN9/fAP/9JulfOJcf3GRpBOlvauHJ9bu4MGV23hgRTPLG/cCUFWSx/zpNVzWUM0lU6tfciU1SZIkaaAskXTK5GRn8ckbZvHGrz3G79bmcnWmA0nSIBVjZNPeXr7+wBoeWLmNJ9buoKO7l7zsLM6dNII/Oms6l06tZnpdKVnONpIkSdIJYomkU+qC+kpeNXsUv37hRd7duJeptW6yLUlHo3FPOw+vauahlc08tKo5vSH2UqbWlPDm8ydwaUMV50+qpDAvO9NRJUmSNERZIumU+8R1M3lg2Yu874fP8Iv3XERBrr/wSNLBWjq6eWLtdh5c2czDq5pZ0dgCwMjiPC6aUkV1TzPvuP5iRpW7y5wkSZJODUsknXI1ZQW84/R8/uOpvfzzb5fyqRtPy3QkScq4rp5eFm3cxcOrtvPwqmae3rCT7t5Ifk4W500ayWvPGsvFU6uYUVdGVlZgwYIFFkiSJEk6pSyRlBGzq3N4+8Wj+cZDa7loShXXzKrLdCRJOqV6eyNLXtzDo6u38/DqZp5Yu4O2zh5CgNNGl/OOS+u5ZEoVZ00Y4YxNSZIkJYIlkjLmr6+dxmNrt/M3P32O08eUM7rCf1GXNHTFGFnT3Mojq7fzyKpmHl2znV1tXQBMri7mtWeN5aIplZw/qZIRXkVNkiRJCWSJpIzJz8nmv950Ftd98UE+8KNn+eEtF5DtVYQkDRExRtZvb+OxNdt5dM12HluzncY9HQCMLi/gyhm1XDi5kgsnV1FXXpDhtJIkSVL/LJGUUZOqivn0q0/jQ7ct4r/uW8kHrmzIdCRJOiYxRjbu2PeS0ujF3e0AVJXkc0H9SOZOruSiyVVMqCwiBEtzSZIkDS6WSMq415w1lodWNvPFe1cyt76S8+srMx1Jkvq1f6bRE2t38Nja7Ty+Zgebd+0DoLI4jwvqK7lgciVz60cyubrE0kiSJEmDniWSEuFTrz6Npzfs5AM/fpbfvv8S9wORlDi9vZFV21p4fO0OHl+znSfW7qBpb2p52sjiPM6fNJJ3XlbPBfWVTK2xNJIkSdLQY4mkRCjJz+G/3nQWr/mfh/nIz57jqzef7S9gkjKqu6eXJS/u4c51Xfzg1oU8uW4HO9MbYdeVFTB3ciXnTRrJ+ZOcaSRJkqThwRJJiXH62HI+cu10/uk3S/neY+u5ee7ETEeSNIy0dHTzzIadPLluJ0+t38EzG3bR1tkDwITKvVw5ozZdGlUybmShpZEkSZKGHUskJcqfXTSJh1Y18+nfLOWciSOZMaos05EkDVFbd7ezcP0OFq7bycL1O1iyZQ+9EbICzBhVxuvPGcc5E0fQtWUZf3Tt/EzHlSRJkjLOEkmJkpUV+Lc/nsMrvvAg7/vhM9z+3osoyvOPqaTj09XTy5Ite3hq/U6e3rCTp9fvZEv6ymmFudmcMa6C986fwjkTR3Lm+ApKC3IPPHbBjhWZii1JkiQlir+dK3GqSvL5zzecwU3ffJz3//AZ/uems8nNzsp0LEmDSHNLB0+v38lTG3byzPpdLNq0i47uXgBGlRdw1oQR/Pn4EZw9YQSzRpf53xhJkiTpKFgiKZEumlLFp26Yxcd+uZi//ski/uP1Z5CV5f4jkl5uX2cPi7fs5tmNu3hm4y4WbdzFpp37AMjNDswaXc6bz5/A2RNGcNaECkaVF2Y4sSRJkjQ4WSIpsW6eO5E97d38653LKc7P4Z9efZob2UrDXE9vZM22lgNl0bMbd7Fs6156eiMAYyoKOWNcBW+ZO4Gzxo/gtDHlFORmZzi1JEmSNDQcdYkUQsgGFgKbY4zXHXTuQ8DbgW5gG/BnMcb16XNvBf4hPfSfYozfCSHkA78ExgJfjjF+OT32a8BXYoxPH9/b0lDxnvlT2NvezVd+v5rSglz+9hXTMx1J0inS2xtZv6ON5zbt4vlNu3lu824Wb95Na/qKaaX5OcwZV8G7LqvnjHEjmDOunJrSggynliRJkoaugcxE+ktgKXCoy2U9A5wTY2wLIfwF8HngDSGEkcAngHOACDwVQrgduAR4CPhn4GHgyyGEOUC2BZIO9pFrp9HS0ZUuknJ4z/wpmY4k6QSLMbJhRxtPbO3m0d8t5flNu3l+8272tncDkJ+TxczRZbzu7LGcNqacM8dXUF9V4jJXSZIk6RQ6qhIphDAWeBXwGeBDB5+PMd7f59vHgJvS968B7o4x7kg/z93AtcAuoAjIBfb/BvBp4F0Dfgca8kIIfOqG02hJL20rLcjhLXMnZjqWpGPU1dPLqqYWFm/Zw+Itu1m8ZQ9Lt+xhb0eqMMrLXsf0UaXcMGc0s8eWc/qYCqbWlrj5tSRJkpRhIcbY/6AQfgp8FigFPnzwcraDxn4J2Bpj/KcQwoeBghjjP6XPfQzYB/wncCswA/hXoAU4K8b4yX5y3ALcAlBdXX32bbfd1m/2U6mlpYWSkpJMx3iZJOY6lkzdvZH/fraDZ5p6eMfpeVw0Jrf/B52CXCdbEjNBMnMlMRMkM9epyrSvO7Jpby8b9vayfk8vG/b0sqmll/SF0sjLgnGlWUwoy2J8WRY1uR1MrSkmN0EzjJL484Nk5kpiJkhmriRmAnMNRBIzQTJzJTETJDNXEjNBMnMlMRMkM1cSM0EycyUxE5z6XPPnz38qxnjOoc71OxMphHAd0BRjfCqEMK+fsTeRWrp22ZHGxRi7gT9JPyYXuBO4MYTwH8B44NYY4+2HeNzXgK8BTJs2Lc6bd8Q4p9yCBQtIWiZIZq5jzXTxJT382bef5FuLd3D2nNO59rS6ROQ6mZKYCZKZK4mZIJm5TnSm/fsXLXtxD0u37mXpi3tYtnUPG3fsOzCmvDCXWaNHcNWcMmaNLmfW6DLqq0vI7lMYDYfP6kRJYq4kZoJk5kpiJjDXQCQxEyQzVxIzQTJzJTETJDNXEjNBMnMlMRMkM1cSM0Gych3NcraLgBtCCK8ECoCyEML3Yow39R0UQrgS+HvgshhjR/rwZmBen2FjgQUHPf+7Sc1KugDYDbwBuA94WYkkFeRm8/W3nMObv/E47//hM3zzbedwydTqTMeSho0YI9taOljZ2MLyrXtZ0biXZVv3snzrXvZ1pTa8zgowqaqY2WMreOO545leV8r0UWWMLi/wCouSJEnSINZviRRj/CjwUYD0TKQPH6JAOhP4KnBtjLGpz6k7gX8OIYxIf3/1/udKP24EcB2pvZOuB3pJbcBdeGxvR8NBcX4O3/7Tc3nj1x7jlluf4ntvP4+zJ4zMdCxpyNnV1smKxhZWNKbKov2l0c62rgNjRhbn0VBbwhvPG8eMujKmjyqlobaUgtzsDCaXJEmSdDIM5OpsLxFC+BSwML3s7F+BEuAn6X9l3hBjvCHGuCOE8GngyfTDPrV/k+20jwOfiTH2hhDuBN4DPA985VhzaXioKMrj1j8/jzd89THe/I3H+dxrZvPqM8dkOpY06OyfWbSqqYVVTS2sbEx/bWqhuaXjwLiS/Bwaaku49rQ6GmpLmVZbSkNdKVUl+RlML0mSJOlUGlCJFGNcQHo5Wozx432OX3mEx3wL+NZhzn2wz/12UjOVpKNSU1rAT941l/d8/2k+8ONnWbRpF3/3yhlewUk6hJ7eyKadbazZ1srqbS088EIHX1r6CCubWti97w8zi0rzc5hcU8L8adVMqSmhoS5VGI1yKZokSZI07B3zTCQpCapK8vne28/ns79dxrceXsuSLXv40p+cRXWpsyM0PLV0dLM2XRQduDW1snZ7K537L4kGlOTCzDGBV80exdSaEqbUlDC1ppTasnzLIkmSJEmHZImkQS83O4uPXz+T2WPL+dv/e47r/+shvnLz2ZwxriLT0aSToqunlw072li7rZW1za2saW5lbXMLa7a10rT3D0vQsrMC40cWMbm6mHnTqqmvLmZydQn11SU89+QjzJs3N4PvQpIkSdJgY4mkIePVZ45ham0J7/zuU7z+K4/y6VfP4g3njs90LOmY9PRGNu/cx7rtralbcxvrtqdKow072ujpjQfGjizOY1JVMZc1VDOpupj6qlRZNL6yiPwcN7iWJEmSdGJYImlImTW6nF+992Le/6Nn+MjPnmfRpt184vqZ/iKtROru6WXLrnbWbW9l/fZW1qaLonXbW9m4o42unj8URYW52UyoLGLmqDJedfooJlUVHyiMKoryMvguJEmSJA0XlkgackYU5/HtPz2Pf7trOf+zYDVLX9zDV246m9qygkxH0zDU1dNLY2sv9y9vYn1zK+u2t7F+eyvrt7excedLi6KC3CwmVhbTUFPK1TPrmFRVxITKYiZVFVNT6l5FkiRJkjLLEklDUnZW4CPXTmf2mHL+6ieLeNUXH+JvXzGd15w5hqwsfxHXidXe1cOGHW2sSy81W5cuidZvb2Pzrn2ppWcPPglAcV42EyqLmT6qlGtOq2NiZaoomlhZ7KbWkiRJkhLNEklD2itOH8XkmhL++ieL+PBPFvGth9byD6+awYVTqjIdTYPM7rYu1u9IlUP7C6P1O9rYsL2NrXvaXzK2vDCXiZVFzBlXwQ1zRtPevIFrLzqbCZXFVJXkWRRJkiRJGpQskTTkNdSW8vN3X8SvntvC5+9Yzp9843GumF7DR185nSk1pZmOp4SIMdK0t+Ml5VDqa2oJ2u59XS8ZX1Oaz4TKIi6aUpWaTVRVzISRRUyoLHrZHkULFrzIORNHnsq3I0mSJEknnCWShoWsrMCNZ4zhmll1/O/D6/jy/au45j8f5E3njeMDVzZQVZKf6Yg6Bbp6etm8c9+Bcmj99jbWbW9jw47UMrT2rt4DY7OzAmMqCplQWcT1c0YxYWQx4ytTJdH4kUUU5fmfT0mSJEnDi78FaVgpyM3mL+ZN5vXnjOUL967k+49v4BfPbOHd8yfzZxdNynQ8nQBHtT9RWkFuFhNGFjOhsphLp1anCqLKYiZWFjG6opDc7KwMvhNJkiRJShZLJA1LlSX5fOrG03jL3Il87nfL+Pwdy/n+YxuYP6qHM9o6vWR6wrV0dB+4wtm67a2sb059XbGljZ133PGSsWUFOUysKj6wP9H4yiImVhYzobLIK55JkiRJ0gBYImlYm1JTwjfeeg6PrG7mc79bxveW7ua2z9zLlTNreN3ZY7l0ajU5zkbJiF1tnX8oidIzidan9ydqbul4ydiqknwmVhYxszKbC2bVv6QoshCUJEmSpBPDEkkCLpxcxS/fcxG3/uo+1oU6fvnsFn77/FaqSvJ59Rmjee3ZY5kxqizTMYeU3t7URtYbdrQdmFW0vs/9gzeyrisrYGJVEVfOqGFCesnZhMrUPkUl+an/lC1YsIB586Zm4u1IkiRJ0pBniSSlhRCYUJbNW+fN4qOvmMGC5U389KlNfPuRdXzjobXMGl3Ga88ayw1njHYj7qPU0tHNhu1tbNzZxsYdbWxI3zbuaGPjzn10dh9+I+vUTKLiAxtZF+RmZ/CdSJIkSZIskaRDyMvJ4upZdVw9q44drZ3c/uxmfvr0Jj716yV86tdLqK8u5sxxIzhrQgVnjhvBtLpSsrOG1946vTHStKedzbv2sXnXPrbs2seWXe1s2pm+v3sfu9peOpuoND+H8ZVFTK0p5YoZtYwbWcS4EYVMrCxmzAg3spYkSZKkJLNEkvoxsjiPt100ibddNIllW/dw79ImntmwiwXLm/jZ05sAKMrLZs7YigOl0hnjKwbtbKUYI3s7umna007jng4a97TTtDf9Nf194952Xty5j+47733JY0vzcxgzopDRFYWcNaGCMRWpWUTjRxYxbmQh5YW5bmQtSZIkSYOUJZI0ANPryphel9obKcbIhh1tPLNhF09v2MkzG3bxld+vOXAJ+ZL8HOrKC6grKzj01/ICRhTlnZQZTDFGOrp7aevsoa2zm7bOHlo7utnb3s3Otk52tnayo62LXW2d7GzrYmdrZ5/jnbR39b7sOUvyc6gpy6emNJ8zx43g9PJu5s5uYHRFqjQaM6KQsoLcE/5eJEmSJEnJYIkkHaMQQnrPnmJefeYYAPZ19vD85t08t2kXm3buY+vudrbuaeehlc007W0n3S+9RG52ID8nm4LcLPJzssnPzaKgz9e9e9r55urH6emN9MZIby/0xP33Iz0x0tMLHV09tKYLo7bOngNl1uHzQ3lhLiOK8hhRlEtdWQEzRpUxoiiX2rICqkvzqS0roLasgJrSfIrzX/qfiwULFjBv7sQT9XFKkiRJkhLOEkk6gQrzsjlv0kjOmzTyZee6e3ppbulk6552tu7ex4u729m9r4uO7l7au3r+8LWrl47uHtrTX1u7Ijnt3WRnBbJDICsLckIW2VmBrBDSX6EgN5vivByK8lNfC/OyKc7Lpig/58Dx0vwcKoryGFmcR3lh7rDbx0mSJEmSdOwskaRTJCc768AyNsZVHPXjUpetv+jkBZMkSZIk6Sh4KSRJkiRJkiT1yxJJkiRJkiRJ/bJEkiRJkiRJUr8skSRJkiRJktQvSyRJkiRJkiT1yxJJkiRJkiRJ/bJEkiRJkiRJUr8skSRJkiRJktQvSyRJkiRJkiT1yxJJkiRJkiRJ/bJEkiRJkiRJUr8skSRJkiRJktQvSyRJkiRJkiT1yxJJkiRJkiRJ/bJEkiRJkiRJUr8skSRJkiRJktQvSyRJkiRJkiT1yxJJkiRJkiRJ/bJEkiRJkiRJUr8skSRJkiRJktSvEGPMdIZjEkLYCyzPdI6DVAHNmQ5xCEnMlcRMkMxcScwEycyVxEyQzFxJzATJzJXETJDMXEnMBMnMlcRMYK6BSGImSGauJGaCZOZKYiZIZq4kZoJk5kpiJkhmriRmglOfa0KMsfpQJ3JOYYgTbXmM8ZxMh+grhLAwaZkgmbmSmAmSmSuJmSCZuZKYCZKZK4mZIJm5kpgJkpkriZkgmbmSmAnMNRBJzATJzJXETJDMXEnMBMnMlcRMkMxcScwEycyVxEyQrFwuZ5MkSZIkSVK/LJEkSZIkSZLUr8FcIn0t0wEOIYmZIJm5kpgJkpkriZkgmbmSmAmSmSuJmSCZuZKYCZKZK4mZIJm5kpgJzDUQScwEycyVxEyQzFxJzATJzJXETJDMXEnMBMnMlcRMkKBcg3ZjbUmSJEmSJJ06g3kmkiRJkiRJkk4RSyRJkiRJkiT1yxJpgEIIifzMkpgriZmSys9KkiRJkpR0ifzFNYQwLYQwItM5+gohnBFCqIsx9mY6S19JzJXETOCfq6EkhBAyneFQkppLRyepPz9zDX5+VpKkUyGp/3+TxH8wT2ImSO7PsK/EfXAhhKuB24FL099n/ENMZ/oVcFP6+0R8bknMlcRM4J+rgQghXBhCeFcI4bIQQlGm8wCEEOaFEP4uhPDmEEJlTMgVARKcK4k/wyRmSurPz1xHnylxf67Az2ogQghnhhBuCCFMCCGUpo9l9P8Pk5gpnWF6COHiEEJ+CCE303kgmZkgmT/DJGZKai7/XA0o08UhhHeGEC4JIdTEGGOmM6VzXRFC+GgI4U0hhPExxt5M50pipnSuRP4MjyRR4dK/VH8WWA+8ASDTf/FKZ/occBdwVjpTb6ZLiCTmSmKmPrn8c3V0uW4AvgqcDXwIqM9kHoAQwvXAF4D/396Zx+s1Xnv8u3IyIycyGIMQFYIaaug11dQkDUqJkEQGaigRtAgx1dDqgEtruJ1R89yqqS1Kq2q6qLGKVFGkiJgjwbp/rLWT7W1yzvsmeZ+93PP8Pp/zOWcPsr/Wevbaz7ieXsCXgI1L16os71G5IvowIlNU/2Wu+pnClSvItmpEznURMAI4BPi2iCxfZcU+IpNz7QxcDUwCzgQmi0ifqniiMkFMH0ZkisqVy1VDTMOB84F1gNHAhSKyRoBytS1wHtAF2AD4nYisW7GtwjE5V0gftitVDfEDbA08B6zvx3cD+1TMtDnwNLCRH98HHB/AVuG4IjI5Ry5X9XMtDdwArOfH5wETgRWBXhUx9QAuAzb145OAKcDawDJ+rlPmCu3DiExR/Ze5PsXlKtuqYa7uwMWlb+EWwL3AJcCAzPQJrt7Ar4F1/XgM8L/AqUD/zBTbhxGZonLlctUw18nAJP+7D3Ak8CCwhp+Tiri+AZxcOp4M/Lvk1yrqMuGYIvuwvZ9IvVv/AvZQ1Yf9+FJgEFQ6+vqMMz3gxycDy4lI74p4CkXkisgEuVw1otlYA2gzEemLjbTsApwBHF/RKNDHQD9guIgMAMYCmwJHAOeJyIpaTT6pqFwRfRiRKar/Mlf9iliuINuqESmwPLARgKreBTwMzAD2F5GWzDRXc7CZbZ9xrkuxwagewA5QSZ0mIhPE9GFEpqhcuVzVoZIN3gRWcaYZqnoa1tY5TUT6qfdCVMD1AtBanFfVs7H2zqUiMiDldzAiUw1XKB/Wq8o7kUSkVUSWVtWnVfXe0qU7gLEiMryCF2A5n544XVUfLF16HNgEGJ6SJzJXRCbnyuWqfq5uItJTVd/FltiNB34FXKCquwBnY42iQQmZWkWkt6p+ABwGDAV+ClylqrsBJwDT8Y955grpw4hMUf2XuepnCleunCvbqn6usq1OBUaKyMkicobzXAusgTXaOiyTc5V9eD6wq4gcIiLfcZaHsSWTpKrTRGRyrnA+jMgUlSuXq8ZUssEVwJ4ickjp8kXY92allEw1XH8EthWRqWAdJqp6DnArsF5HZ6rhCuXDelV1cquRwFXAzSKyj4isW1xT1cewl3WsJJyhISK7OdN1IjJFLHdNwfQP4HvAgSKyciqmqFwRmZwrl6vGuK4Arve/78KWAP4BeNb5/oyNACVpaJT8d4uIHADMUNXNgZuwmWWo6gvYmuYVUzAF54row4hMUf2XuepnCleunCvbqn6usq32AWYChwOdgddVdTdV/QOwJIkqzxGZnKvswy8DD2Gj04OxZRejVfVCoIck2nk2IpNzhfNhRKaoXLlcNcS0rYgcVRyr6ovAzsDeInKon5uOxfZ1UjA515dE5IclrunONVlEppY6TLqQrt4Xjsm5QvqwUXWu6sEisiI2MjcBm7b4FWB1sdkjf/Tb7sXWwvbDXtxmM/UFjgP2waZTfhHrCe+jqpf7bXdiSa8GAc83mykqV0Qm58rlqn6u1bFpnF/FpuhOxJKtXoQ1fg4QkQ+ceTBmt2Yzlf3XivlvoIhcCdwCnC2WZO7fWILY05rNFJwrog8jMkX1X+aqnylcuXKubKv6uWpttQu2ZOUqVT2mdN84rEH2dkdk8ueVfbiC/74DuFpVbyrdN8Gvf9QRmfx54XwYkSkqVy5XDTFtB1wD3C82O2oqgKo+JCITsU64lYEPsFmvJzWbybk2A34OvCQiF6vqXs41TUS2AK52e3bGdsk+pyMyOVdIHy6MKutEwhKUvQs8oqoficjLwChgmIi8qqpPqupfReQ+LMdACrUAbwH/UNWZIvI6sD3wBRH5t6rerqqvisjdwLRETFG5IjJBLleNqBUbSbkHQEReAPbDOrN+iyUh3x+YBUzwGVPNVq3/XsL8tzs2SnUmtlTkDWBvVX06AVNkrog+jMgU1X+Zq35FLFeQbdWIFmSrkSLynqr+TUS+AhwMjFPVGR2UCf7Th89jPhwlIj/D6hTjgKOBUar6Vgdlgpg+jMgUlSuXq/rVG1uafB1wrIh8V1WPBvC2zWZYaoxlgJGq+kwCJoAlgBOBC4GLReRSVR3jXNPEdh/bABsoP1NVn+qgTBDXh41LK8zqDfwPcAzQ3Y8H+7nxFTL9APgZsIQfr4AFrm9UbKtwXBGZcrlqiEmwnB2jga5+bmNspPpLftwdWDKS/4CuQLdo5aoKrog+jMgU1X+Z69NfrrKtFrut+pB4p6OgTPPz4UbAzcBOftwfWK0jMwX3YTimiFy5XDXM1LNkowuB75WuVbK7mD+7V8kmVwKX1V7LTLF92OhPJTmRROZmI/8VtqRojIj0UOsFvBYYLyI9EzMVtjgX6/U+SkSWUNWXsJG7nSXhOtzIXBGZnCuXq8b1v1jjYksR6aaq92MBbV8/nqWq76QAacd/12H+W1JVZ6slOEyiqFwlhfFhRKao/stcC6Uw5QqyrRpRnd/npdR2pnmxozLVqNaHD2AJhyeKSHdVfVVVU85eDscU0YcRmSJzuXK5qlOq+p7/+TC2BGs5ETlKbCnUwSJSySoj9RliajOyDgI6i8iPnOtoEememeZyhfRho0raiVQ0qNW72rAEjw8BQ4CpItIFWBp4j0RLjUpMxfOexQJED+BHItIPy7r/IenW4c71SxQuEVkyGlOJrXarzwjlSiCercp8/h5egOWFGoHliQIbNZ/tbMnUTlzojfkvKVNErqJsRfJh+R0MxFTP96Y3uVwhIt2CcoUr62VFslWhT6mtiu/znI7MVGd574b5MAlXnbE9KZNzNdKe6LDlqoHYnporXFkvK6Kt5sP4IdYJcQxmszOBO/x8pVLV11R1d2AbjOtyVZ3V7OcuyFZVMrXFFdmH9UjmvSdNfIjI2sArqvq6iHRS1Y9FpLOqfigiy2FrEncE1sXWMO6nn9wCvRlMWwF/V9VXSuda1Na9DsCmvk3AAkcf4MBmMznD9liW/1uKwl01l4hsg63PPBb42P0XwVaDgXeA2ar6qp+rulz1AWaVepkRka6qOrtiW/UEPiqPhpe4+mEJAzfFpuq2AhNV9aFmczlHp6KzrWr/OUN3rJzPLp2LwNULmKOq75fOdVHVOVX5UERWA15T1beKymAApk2Aaar6WulcBP9ti1U87y0qqEG4tgHWVtVzSnG96ji6Mpak9B1VnePnKi1XztADi6Pl2FDE0apsFTa2O0uY+J5je0NM4WK7c4WL7zm2N8SUY3v9TH2x9sS7pXOF//oA76nqLLFd404ChqnqE81kquFrUdWP/O+ibb885tu3RWQX4L+xJYmPN5mlN9YWLLe9KmVyhuWwvFrvlGJDwVW5DxdZ2vx1f0OwXUkuA5b1cy3+ezssc/oyfjwA6J2AaSjwD2DT0rmiQ21bbDvHlf24Fc9jk4BrGPAi8MXSuU4lWyXncqaXsaC/Zg3TNhXaaifgAWzJwHFAz1K52r6icrUz8BtsxtE4YKPStUr858/a1Zluxz6Cg0rXvoj1fLf68RCgfwKmrYB9SsedSn9vW4X//Fk7YtvK3g5s7edagnBdD/y+iA+l93D7iny4oceG7wB9amxVFdNQ4FFgvdK5IrZX8r3xZw3HYvuw0rlK41XJXjOwXUAG1JSrSso78GUs8fPvsAZE5eXKn7Wzx4Zr3G6rlK5V9c0JF9v9WeHiOzm2N8IULrb7s8LFd3Jsb4Qpx/b6mXZ1O92BJRkvt1e3xXbaXM6PxwDrJCxXU0vH5dj+Bex7tKofjwAGJ2DaCbgB2+V6DLBUKS5UwuTP2gH73tyI7c7dwifbz5X4cLH+PyYwYmfgciyx8JWlYLYstq3syKT/w9Yp8lfg837creTUpYB7KmASLMnlNcBufq7Vf/p7wbsH2D0x107Ag1jm+snA1XgiTqCv+2/XlEz+7MHAY9iIwIZYR1Jfv9bLuVLbag2scjMEq0CfhnWcbgl0cf/tVoGtVgX+jq013wvLzfRNt11/4L4KbDUU26noTuDwmmvLO1MVthruPtwSmISt0R8cgGtHfw8/j1Um/sK8BO0rV+FDf3ZvLK/Xd7CdJvpXyYRVUB8ENvTjltK1qr43gs0s/T3zEhov4d+drlXFq5pyNQSY6j4skppWUt6B9UuxfUesEr1UDdOoCmy1nseGzwJfwRocpwNrYt/pKspWuNjuXOHie47tDXOFiu3+7FDxPcf2hplybK+faQXgKax9MxRb6vQjrKO0pzNVEa+2wiaF/A04vXS+BVtZkbxNiA2WPIYlqN4d2zxiU79WZTt1B2wZ5MZYx9XtwNJ+rV8V5aoZP01N3CQiLVjnCFhv6meBk0Xkcj8/UlVfKK2JTaEvAj1U9R4R6Q+cCvQSkT9ilZ7hatuwJ2Py58wSkX8C94jlH/oV9rJ+iI1SfUFVP0jF5dO+twemqOqzInIX8DnsY/2O2tLEnVX1lcT+A6sgv6iqj7qtBgHnisg0rHI4wvlSchVMTwCIJcs+CNsO9FVgZ1WdXoGtejnX/cD9IvI3bLrwLlh5H6+2TWlKrtWB7wO3YNtbHq6qZ/i1WcAeqvqP8jKIZsvL+zDgFFX9E/AnEVkBC/5PYVPVR6nqcyltJbb+fhvgGI9ZL2INxhNF5LfONVZVn07M1QkbIJiNzeochCXsvQNYDthLVf+euFyNBFZU1Qd9ecg3PcZfBzxJBd8bVVURmYWN6v/Wp1xfhI0OP4tVnHdQ1dcS+29pbKTzSFV9Qmxp8CisMgjwPhWUd2AV4AmP7f/EZpn+UETuA16qiKnM9QjwiJercdg3+lpgtNqWwSm5IsZ2CBbfPbYPJ15s70zM2N5CvNgOweJ74NjemxzbF4YrUmzvCjyvvlxORJ7B6qe7AYp1IL0oMjf3aiquAVhak2uBK0XkDFU9XG1J2wwR2UVVX/b6oTaby5+zJdah9QDwgIgMwjqT7vV2YFKmkjYGjlfV+31JWyvwPY+jL2Btwlcq8OFiVVMTa6vqR2q7f9wMfKCqJ2GdR1cBXTzgd0ppPFU9AviziNwP/BrrKbwWG2EZAbyZmqmkztg28KdguxIcjtlqJJa5PWXjZxa2/fytfuphbJTl5NI9r/jv1Lb6M7CUiNyKfawvxta4vgRsDcypwIcPAu+JyHF+vAE2SjwLGKiq0yG9rVT1r1iZPtiPH8DK/UrYlqXJGxmqeh62Reqj2K4EnxeRKX7tDcyPpOpA8mfNAs7GKoNFXHwLaxChqm+q6nP+d0pbzQGOUtVbRGQpbLbiPcAjWINjtdSNDOf6WC0nxc1YZflc5s0KfNsbGalj+3jgbm9M/wZ4Hrgfq1BvmLoDqcT1IZbQ/8fAWVhM/67zbQXMrsBWbwAHq+ptfnwd1kA8zY9nVlHesdG5/iJyBdYwvB77Pq8ErOWNjCq+z48AH4rIOD9ezvnWw5b7TIM0tvKO5SK2z4gS20Wkq3Och8WDyuO7WP6QWc5zc5TYLpYD5kOsAylEbPdOraLeHia2l8rVeKwDsPL4XrLVh1hy5RCxXSxHzUzgkCixvfAfln6iT5TYXuJ6HPNThNje15/1HBbbT/fjadjMrenASmodSJ3UlYCrn3NcClylqq8DXwXWE5GzSre+7/d93GwuEenr35CzgevEhX1flind+kYqJucqbHWiqt4glrPwOmxJ26+AdYAtgOkpfdgsLfZOJBHZXER2FZE9S6d7AOuKyObAf2EjVKNFZPlEFYmCaS8AVd0b+/jcqKrnqeoV2Au6FTblM9Xsh4JrjHMdAvwTmwZ3q9q2kXdgs5FSvQBz/afzEqZ18Wd/DVhaRIY2m6MNrtHOsiXW2fY7Vf2+qt4H3IqNKnROXK72VEuI+UNgUxG5HthAVSdjnUuji97mFBKRrUVkgoh8zU+dDwwSkVEAaiPX9wCTSr5NxXSAM8z0CthfgPOAjUVkbxHZHRhVquyn4JooIvuo6jRvUBRl5zEsIR4iskfKcl/i2s/thKq+jc0uOEpVL8KWG+xRgQ/3K53ug00jHogtyfgLsHXpA5+K6SAAVf0K1gj7s6r+QFXPxhpDe3qDMlUFteA60E8dh80sGAJc4Q3+Il61JPzmFOVqX/UElDJvK9lJWCV/rRQsNUwTnOkVYDy2BP5uVT1VVX+DlfUviEjPxLaaICJf9Ur9rcCuIlJMl98fG9AYm4LHmXYAThIbbQW4BFi9ythe4jpRLBkzqvpW1fHdmU4RkQGq+rSqvh0ktu+AzcgfqJ4MPUBsn8tUOt0Xm30+kApie4nrRBH5DICqjqTi+F6yVfEOnoI1nquO7TsA3xaRVdQG7+d2OFNdbC/8t4qq/gvYlxixveBaGevw+yPVx/ZhwKlFJwS2lHQJETkCQFWfxTriRotI94S2GoaVqz7OMdN/PwfsD6wjIieIyEjggFKZazbTqSLSiiX/f7PUGfMUvmueiIzGvjctbfxzi5vr22IzvQHwetbuqnqCqt6ELX/dnIR9Dc3UYv2Qi8gI4CdYT9vxInKaX7oKm5lxLfANYCKWl6jpjq1hmioiZwKo6kHYiEGhvthW601/AebDdax4jzNwBjaN8ufe8TAU+4g3fbu/+fjv+2CzIZzlPWxG0vrNZmmH63TvVb4TeEdEjvZbV8fW8qcuVyeIyLdU9U5V3QkLrLv6rUsBMxNW5kdgo7/dsfK+P3AbtoZ5cxE53G99H9vZrumdW/NhOrS45pXou4ETgKOwmXgPJeqAKLi6YVPjD625ZTbwkVgn78lYB2/TVcN1bA3XM6UOSQXepMkzSmuYujvTYX7pCmyd/iXAIVhOFiXdO1gwTSlVtkZh/izsJNjMg1Qdk2Wuo0VkMrbk4XrMLmf7rWth8aqpy8rnw9UNOK4oV2o7vQjwurNskYKnhqk7Fkcnq+oLaqPnM/zdA1u61YJN80/N9U0R2Rf4JdaoOAJbMgZmrzdTMLm2xPKJbC8iK2IDck9is32O8HuSxfYFcBWN6yK+/4UK4rsz7QCMKDO55lBBbK/hGlpjq6dKnWvJYnuJaUdn+oyfuxxbbncpiWP7fLi2EZE1YW58P7Zkq6TxvcS0nYisipWdX/u1c/x30the4hoBDCvKlfruZ1iC7aSxvcS0IzBcRFZV1ec8tr8qPqhP4the4toJs1dfVf0pNvvvMCwvEiSM7f69+RZwmc7bdfBJrM28hoj80M8ticWuVJ0iZa4ZNdda1GZIjQAOAH4K3FQqc81murzoPKq55WPgYxGZiH137ismRCTiukxtxm1Z00v10b7OmKSvoenSxZRcCfgM1ku6hR8PxAJrP6whPQXL6wMW9LsurmcvBNOyeOZ2Pz/J71u32UxtcF3PvOTQ3bEOih9jveFN52rDVv1rbDUMmIYtbZMKuYqM9kP9+HZshOqzFfpv2Zr7DsM6S1OXq238eCxwIPN21fsSNp3yNuAJbLZUFUwHAJ+rue8AbERoSIW2+gSXl/V3sU6utSrm2qjmvn39vqbv5tAG0/pYResHzEssKqTZNXJBTJvU3Lcf1cT2MtfXirjk1y/COt8eSBGv2rFX7Xu4Bza9vwdNju/t+RCbLn8tNtMgSWxvg+vA+byDR2Dxfe0UXCWb3IbNLj3K378VsE6tq0kY29vgOgLffZR5Oy8lje/tMHXG8j0mje1tcPXGkuUWG7wki+1tMC2JpXj4EfN2jEsS29vh6lNjq6TxfT5MU5ynFUsYfX7q2N6GrWrfwWSxfQFMR7r/emH5Qq/CkiAni+1t+LA3n9xkKVlsx5bzvYbveIYtp/s8sJ0fD/bYfhO2RDhJbF8A1yZYztnyfbsBz1Voq08wYbMn3/bYvmYwWx2E5e1NFq+a/VNsgbfI8tGLDVT1Sp861htbA7ivqj4mIl1VdbZPz21qT2U7TDc40+M+krEEcBJwvqo+WiHXjcBXVfXx0n09sKD2bkVMZVvNXWsuIr1U9a1mM9XJ1QX7KK4HPKOqL1fEdCO2pfETha18ZPgWVX2s2UzOVWw9+piILItVGG7CdvF5QFW/XuKfobauuQqmG7GOt8fVZgTiM1tuS/gOtsvl5y8DJpffywBca2AVnx+ksFcbTKsBd6nqVL+vRROM+LTDNJB5dloBG4k6N0C5WhX4q9oS1yL3wZyEcbSu99DvXVr/cyQtJdOqWFLMI8WW16yL2e75ZjO1w7UKloC1iFknA9eo5SZqNlPxTRmIdW7/CTgUa9h3wXKgvO2x4fUUsb0dLsUai5NU9Q0R+Tq2RD9FvGqLqRVrQLZisydTxvb2bDUZG7Q7EjirYlsJNlp+gqr+W2ypWNNnw9fBVdjqYKzz5njgnIpthbPsr6rviy1Hmp0wttf1Dvq9qWJ7W7bqhLW5umIz+R9OGNvrjVenAFcniu2t2BLujbG0KvtgAwHDgd+q6iS/rz9WrlLNjloQ1zCsrv41v28c1r54MgKTiHTDOt2OSVjvq4drWeBo4BepuJJoUXuhsNGKLlii7OJc0Tl1MZYcEGD9RX3WYmRa1f8uRoc7BePaMKCt1kvF1GC5SsbVgK1SjwIXXD38uAs2Q2uUH/fEEpBPDMTUA3gGmBDMVgXXOL/WGoxrNFYJ6xGAqShXyXzYgJ3G+nG3IFw9nWvvVLZq0F7jAzH1xGa7jg1qqwkJmVbCGlvFtu+rYPlgwGZHvQNciG/BHozrAqBfQKalsU6SXsG4ehXlLAjThSn916CtWv1c0+N7A0xR38E+wZh+iXXSR7RV34RMKztTdz+eBMwEDvPjnlh+nwmJbVUP19+AcQGZ9vLj7oG4nsJ2+INS2/H/y88irSEWS052E5Y48eJizTLz1gH3AXr62tervCe1qaqTaQlnusaZNBjXFQFtdXUKpga4yuVqmdJ60yqZCltdmYKphutc4AIRWVNtpt+f1GdKqSV2uwTfpSAI0/tYroWZKZgWgusdVZ2jCUZ9GuR6Xy0f2PsBmIpyNbOZLA0yFXZ6F+bmZInA9Z5zzWjr36qAq7BXqtHNem11MVaxT6KIMcuZbsbyrJwvIkNU9Z/ADWI5mg7FRvRbsOShqfJF1cvVGUsO3XSuBpn2wpb4NH2mSINc48Vm6qeI7fWWqz0kQaLcBrk6A+PcVk2N7w0yjQpqq9EpuBpg6gSMDRqv9kwYr4rvzYUisrqqngvsqKpnlb6D15Co3t4g17VUU2doj+ltmLvbchSua7B8wmiiVVhJtTA9T9gozkrY+sytsRxDh2PJQ9cu3fcLbFrZXTR5vWREpqhcEZmickVkaoPrCOdax+8p1nePBh4CVu9oTJkrM2Wu6rkiMkXlaodpEDZLcjawi9+/MbBixbaqhCsiU1SuiExRuSIyReWKyBSVawFMR2Lb0g/xe4rvzRhsI6OqvoOVckVkisxV1c+iGLIFS/68IvOW9BwK/AsY7MdnYNO/UyW3CscUlSsiU1SuiEx1cK3hx2OwRG5JknRGZMpcmSlzVc8VkSkq1wKYvo4lMO3HvOU8SZbiR+aKyBSVKyJTVK6ITFG5IjJF5VoA02F88nvzZdIn1w/HFZEpMlcVPwtjvNWxHtu+2C4EU2quT8HWUrdgu0GtlsCh4ZiickVkisoVkakBrguwacOfAwZ0RKbMlZkyV/VcEZmictXBdAy2o1BP50q1y1I4rohMUbkiMkXlisgUlSsiU1SuOpiK700LtpnECkFslZwrIlNkrip/GjXgjth2jHdi60u/jPXoTi3dMxD4abL/gYBMUbkiMkXlisjUINdPOjJT5spMmat6rohMUbkaYPoxiRpjUbkiMkXlisgUlSsiU1SuiExRuRpgyt/BgEyRuar+KZICtysR2Qw4DRijqg+JyE+ATYDNgHvEtjq/HNgC2EBE+qhqU5OIRmSKyhWRKSpXRKaF4NpQEmznGpEpc2WmzFU9V0SmqFwNMn0O6E2ChKsRuSIyReWKyBSVKyJTVK6ITFG5In5vonJFZIrMFUIN9MJtRmmLcKA/cKP/vRqWWPg8LG/Auil6wCIyReWKyBSVKyJTVK6ITJkrM2Wu6rkiMkXlisgUlSsiU1SuiExRuSIyReWKyBSVKyJTVK6ITJG5Ivw0YsQWoFfp7wHYTiXL+7lVsC0SWxM6NhxTVK6ITFG5IjJF5YrIlLkyU+aqnisiU1SuiExRuSIyReWKyBSVKyJTVK6ITFG5IjJF5YrIFJkrwk8n6pSqfqSqb/mhADOBGar6sojshSUq66Kqb9b7by6qIjJF5YrIFJUrIlNUrohMmSszZa7quSIyReWKyBSVKyJTVK6ITFG5IjJF5YrIFJUrIlNUrohMkbkiqNiabuH+Y5ELgJeBodhUr0cXE9dCKyITxOSKyAQxuSIyQUyuiEyQuRpRZqpfmat+RWSCmFwRmSAmV0QmiMkVkQlickVkgphcEZkgJldEJojJFZEJ4nKl1kJ1IomIAF2AJ/33dqr69GJm+9QzQUyuiEwQkysiE8TkisgEmSszNUeZ69PNBDG5IjJBTK6ITBCTKyITxOSKyAQxuSIyQUyuiEwQkysiE8TlqkqLOhNpInC/qj6+2IgWURGZICZXRCaIyRWRCWJyRWSCzNWIMlP9ylz1KyITxOSKyAQxuSIyQUyuiEwQkysiE8TkisgEMbkiMkFMrohMEJcrtRa1E0l0Uf6BJigiE8TkisgEMbkiMkFMrohMkLkaUWaqX5mrfkVkgphcEZkgJldEJojJFZEJYnJFZIKYXBGZICZXRCaIyRWRCeJypdYidSJlZWVlZWVlZWVlZWVlZWVlZXUM1b07W1ZWVlZWVlZWVlZWVlZWVlZWx1XuRMrKysrKysrKysrKysrKysrKale5EykrKysrKysrKysrKysrKysrq13lTqSsrKysrKysrKysrKysrKysrHaVO5GysrKysrKysrKysrKysrKystpV7kTKysrKysrKysrKysrKysrKympX/wcwbKZ3fLD/AQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1440x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Uncomment the below line if you want to price with a yield curve based on US Treasury Yields\n",
    "%run riskless-yield-construction-example.ipynb\n",
    "\n",
    "# Uncomment the below line if you want to price with a yield curve based on AMERIBOR rates\n",
    "#%run riskless-rate-curve-construction-ameribor.ipynb"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "recovered-omaha",
   "metadata": {},
   "source": [
    "## Estimating and forecasting volatility of the index, options on which we want to price\n",
    "Now that the riskless curve has been constructed, the next thing to do is obtain an accurate volatility estimate for the index whose options we want to price. There are various methods for doing it. The most common are the Exponentially Weighted Moving Average (EWMA) and GARCH(1, 1) methods. You can read about them on the Internet or delve into John C. Hull's [Risk Management and Financial Institutions](http://www-2.rotman.utoronto.ca/~hull/riskman/index.html) or [Options, Futures, and Other Derivatives](http://www-2.rotman.utoronto.ca/~hull/ofod/index.html). For the purposes of pricing options GARCH(1, 1) is preferred because it supports volatility forecasting for future maturities by incorporating mean reversion (and volatility of equities lends itself to mean reversion). In this notebook I therefore use GARCH(1, 1).\n",
    "\n",
    "### GARCH(1, 1)\n",
    "In order to use GARCH(1, 1) one first needs to estimate its parameters based on the historical data for the equity whose options you want to price. GARCH(1, 1) has three parameters &mdash; ω, α, and β. I use the maximum likelihood method to estimate them, which is an optimization task. The `volatility.parameter_estimators.GARCHParameterEstimator` class comes to rescure. _Implementation detail_: In that class I use [SciPy](https://docs.scipy.org/doc/scipy/tutorial/optimize.html) to optimize for maximum likelihood. Alternatively I could've used Google's TensorFlow 2.x.\n",
    "\n",
    "Constraint optimization for three variables is not fast. If you are not content with the speed of `volatility.parameter_estimators.GARCHParameterEstimator`, you can use its `volatility.parameter_estimators.GARCHVarianceTargetingParameterEstimator` counterpart. This latter class is \n",
    "faster because it uses the so-called _variance targeting_ method whereby it sets ω based on the sample variance of price changes. Then it optimises for only two variables instead of three as GARCHParameterEstimator does. It's marginally less accurate.\n",
    "\n",
    "*NB* When estimating ω, α, and β it's vital to use _adjusted closing_ prices, i.e. prices that have been adjusted for dividends and stock splits. Failing to do that will result in too high an estimate. Yahoo Finance provides such prices in the `Adj Close` column.\n",
    "\n",
    "#### Estimating GARCH(1, 1) parameters for S&P 500 (^GSPC)\n",
    "Since volatility tends to change with time, it's recommended to estimate it based on the most recent 90 or 180 days, or on the same number of days as the duration of the option being priced. In this notebook I use date from the past two years."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "wrapped-lease",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Objective function: 4190.47699 after 103 iterations\n",
      "Optimal values for GARCH(1, 1) parameters:\n",
      "\tω=0.000000830666, α=0.09031, β=0.90553\n"
     ]
    }
   ],
   "source": [
    "TICKER = '^GSPC'\n",
    "\n",
    "# I'll use price changes over the past two years to estimate GARCH(1, 1) ω, α, and β parameters\n",
    "start = BDay(1).rollback(date.today() - relativedelta(years=+2))\n",
    "data = yfin.Ticker(TICKER).history(start=start, actions=False, auto_adjust=False)\n",
    "asset_prices = data['Adj Close']\n",
    "\n",
    "vol_estimator = parameter_estimators.GARCHParameterEstimator(asset_prices)\n",
    "print('Optimal values for GARCH(1, 1) parameters:\\n\\tω=%.12f, α=%.5f, β=%.5f'\n",
    "      % (vol_estimator.omega, vol_estimator.alpha, vol_estimator.beta))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "technical-wales",
   "metadata": {},
   "source": [
    "##### Optional: Estimating with a faster estimator &mdash; GARCHVarianceTargetingParameterEstimator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "timely-morning",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Objective function: 4190.09183 after 34 iterations\n",
      "Optimal values for GARCH(1, 1) parameters:\n",
      "\tω=0.000001003151, α=0.08521, β=0.90459\n"
     ]
    }
   ],
   "source": [
    "var_target_vol_estimator = parameter_estimators.GARCHVarianceTargetingParameterEstimator(asset_prices)\n",
    "print('Optimal values for GARCH(1, 1) parameters:\\n\\tω=%.12f, α=%.5f, β=%.5f'\n",
    "      % (var_target_vol_estimator.omega, var_target_vol_estimator.alpha, var_target_vol_estimator.beta))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "considerable-portfolio",
   "metadata": {},
   "source": [
    "#### Forecasting volatility for S&P 500\n",
    "Now that we have obtined optimal values for GARCH(1, 1) parameters ω, α, and β, we can proceed to forecast S&P 500 volatility. I crated class `volatility.parameter_estimators.GARCHVolatilityTracker` to assist in it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "twenty-criticism",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Volatility of ^GSPC for term 0.6822 years: 0.14558\n"
     ]
    }
   ],
   "source": [
    "vol_tracker = volatility_trackers.GARCHVolatilityTracker(vol_estimator.omega, vol_estimator.alpha,\n",
    "                                                         vol_estimator.beta, asset_prices)\n",
    "\n",
    "# The 3rd Friday of next March, i.e. the 3rd Friday of next March\n",
    "from dateutil.relativedelta import FR\n",
    "maturity_date = curve.date + relativedelta(years=0 if curve.date.month < 3 else 1, month=3, day=1, weekday=FR(3))\n",
    "vol = vol_tracker.get_annual_term_volatility_forecast(curve.to_years(maturity_date))\n",
    "print('Volatility of %s for term %.4f years: %.5f' % (TICKER, curve.to_years(maturity_date), vol))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "under-lodging",
   "metadata": {},
   "source": [
    "#### Optional: Visualizing volatility term structure of S&P 500 index for 200 business days into the future\n",
    "In contrast to the `equity-options-pricing-example.ipynb` notebook where I visualized historical, and future volatilities, in this notebook I'll visualize the _volatility term structure_ for S&P 500, i.e. forecast average volatility for a certain future term expressed in years. Given the mean reverting nature of volatility, volatility for the longer terms trends towards the mean volatility of a given asset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "warming-tender",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABJ8AAAJpCAYAAAADlpNPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAAsTAAALEwEAmpwYAABRhUlEQVR4nO3debhkZ1Uv/u8KCRAIJEG4kbmZEQlDCINe/ZE4MIVJUBEViYigOKAGtb1XZBSjggNOiFxEHIh6UQw0EiMkIApCEgIJF5ApMhpAIBDmJu/vj9odKtV1Ttc+p98+Ved8Ps9znq7atWqdVW/Xrt61+t37rdZaAAAAAKCHw7a6AAAAAAC2L80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIA2EJVdbuqurCqPltVP73V9Sy7qtpVVa2qDt/g8/9XVb1gXq6q+seqevTBrBcA0HwCgB2nqi6f+rmiqr4wdf8HlqC+U6vq9R3zP7Wq/qJX/g34hSTntNau01p77uyDVfWNVfVPVfXJqvp0VZ1fVQ+Yevx6VfXyqrqsqj5SVb8wJ8clU3/Pl1bVi6rqqHnFVNW5VfXFqffEu2Ye//6q+s+q+lxVvayqrjdTy98Pj/1nVX3/pkZmk6rqpKr60PS21tqzWmuPnRffWrt/a+3Phud2fR8CwE6i+QQAO0xr7ah9P0k+kORBU9v+cpEcG511crBU1dW28Hcf7Nd+8yRvX+fxlyc5O8nXJ/kfSX46yWemHv/5JNdMcsMk35jkX9fI86Dh7/yEJCcm+eV1fudPTr0nbrdvY1V9Y5I/TvKoJMcl+XySP5x63h8k+fLw2A8k+aPhOQDADqb5BAAkSarqsKraXVXvrar/rqq/2TerZer0pB+pqg8kec0wM+Rfq+q3hxk576uqbx62f7CqPrbeKUxD3PuG083eX1U/UFXfkOR5Sb5pmHXz6SH2RVX1R1X1yqr6XJKThxk6j53J9/qp+99YVWcPM4YuHU63ul+S/5XkEUP+tw6xl1TVd0w998rZUfNe+7D9MVX1jqr6VFWdVVU3X+e1Priq3j6M07nD60xVvSbJyUl+f6jntjPPu36SWyT5k9bal4eff22tTc/I+UqSj7XWPt9a+1Rrba3mU5KktfbhJP+Y5I7rxa3hB5K8vLX2utba5UmenORhVXWdqrp2kocneXJr7fKhxjMzaVTNjseNhplY07Om7lpVn6iqI4b34i8Ps6c+VlUvrqqj5xVUVT88/D18dng/PX7Yfu3hdd5oahbXjWqdmW/73lPz3odVdffhfXS1qfiH7XsPAQBr03wCAPb5qSQPTXLvJDdK8qlMZrJMu3eSb0hy3+H+PZO8LcnXJfmrJGckuXuSWyf5wUyaKvud3jU0Bp6b5P6ttesk+eYkF7bW3pHkx5K8YZh1c8zU074/ya8muU6SdU+HqqrrJPnnJK8aXsutk7y6tfaqJM9K8tdD/juvPyTzX3tVPSSTJtbDktwgyb8keckatdx2eOxnhthXJnl5VV29tfZtw3P3zTT6j5mn/3eS9yT5i6p6aFUdN+dXvDnJI6vqRxZ5EVV10yQPSPKWdcJ+bWgE/WtVnTS1/RuTXNlsaa29N5OZTrcdfvbOvIa3Ds+5itbaR5K8IZNm1T7fn+T/tta+kuTU4efkJLdMclSS31+j1o8leWCS6yb54SS/XVUntNY+l+T+ST4yNYvrI+u85un69nsfttbenMnfx32mQh+V5MWL5ASAnUzzCQDY58eS/O/W2odaa19K8tQk3z1zmtlTW2ufa619Ybj//tban7bWvprkr5PcNMnTW2tfaq39UyaNiVuv8fuuSHLHqjqytfbR1tp6p54lyT8Ms36uaK198QCxD0zyX62157TWvtha+2xr7d8P8JwDmX7tP5bk11pr72it7c2koXWXNWY/PSLJntba2UNj5dlJjsyk4bau1lrLpAFzSZLnJPloVb2uqm6TJFV16yTPT3JSkt1V9Zhh+zWq6sszs4VeNswke32S1w41z/OLmTR8bjzkfnlV3Wp47Kgkl83EX5ZJQ/CoXPV0wOnH5vmrJI8c6q0k3zdsSyYzrH6rtfa+YYbVLyX5vppzymNrbU9r7b1t4rVJ/inJt67xOzfrzzJpqmaYtXXfqZoBgDVoPgEA+9w8yd8Ppxh9Osk7knw1k+v37PPBmedcOnX7C0nSWpvdtt/Mp2FWyiMyaeJ8tKr2VNXtD1Df7O9ez02TvHdE/CKmf//Nk/zu1Fh9Mkll0rCZdaMk/7nvTmvtiiHXvNj9DM3An2yt3Wr4vZ/L12bb/EiSM1trr8tkRs7ThwbUvZK8tbU23Sh66DCD5+attSdMNRBnf9+/D826Lw0X3/7XTGZKJcnlmcwwmnbdJJ89wGPzvDST09pumOT/y6QZ+S/DY1cZs+H24bnqezFJUlX3r6o3DqdXfnqo9fpr/M7N+oskDxpm7n1vkn9prX200+8CgG1D8wkA2OeDmZwGd8zUzzWHawTt0w7WL2utndVa+85MLpT9ziR/coDfMbv9c0muNXX/66dufzCT2TuL5DlQrnnP+2CSx8+M1ZGttX+b87yPZNI0SnLlLJ+bJvnwnNh1tdY+mMmpkPuu13R4kiOGx96f5H5Jfj3JC4Y/D4aWSWMtmVwY/cpTFavqlkmukeQ/hp/D983KGtw5a1xMvbX2qUxmKT0ik1PuzhhmeiUzY5bkZkn25qrNzlTVNTJpYj07yXHDaZqvnKp3M+/X/Z477AtvyOR0y0cl+fNN5AeAHUPzCQDY53lJfnXfqWNVdYPh2kYHXVUdV1UPGWaQfCmTWTNXDA9fmuQmVXX1A6S5MJOLXV9rOP1s+ppHr0hyw6r6meEUtOtU1T2n8u+qqsNmcn3fcLHrE5N89wF+9/OS/FINK7lV1dFV9T1rxP5NklOq6tur6ogkpw2veV6j6iqq6tiqelpV3Xq4CPf1kzwmyRuHkL/L5OLpDx0uhP2ZTK6zdKtMVqIbpaqOqar7VtU1q+rwqvqBTGYlvWoI+ctMZv586/B39/QkfzfMlPrcUM/Tq+raVfU/kzwk6zdo/irJD2Uy3tOnr70kyc9W1S2Ga4btu07X3pnnXz2T5tfHk+ytqvvnqtdkujTJ1611sfIDWOt9+OIkv5Dk+ExeLwBwAJpPAMA+v5vJ6mT/VFWfzaTBcc/1n7JhhyX5uUxmuHwyk4t5//jw2GsymS3zX1X1iXVy/HYm15S6NJNr8fzlvgdaa59N8p1JHpTkv5K8O5NrJyXJ3w5//ndVXTDcfnImDZtPJXlaDnAdn9ba32cys+iMqvpMkoszubj1vNh3ZXKdoN9L8omhpge11r683u8YfDnJrkwunr7v93wpk4txp7X2hkxmDT0lk+srvS7JuZk0c15SVXdd4HdMOyLJMzNp5nwiw0Xo911EfLgu149lMtYfy+R6Tk+Yev4TMrme1ccyaSD9+AGu5XVmkttkcn2u6VXjXphJ0+p1Sd6f5ItDLVcx/D3/dCYNvk9lMhZnTj3+zqGO9w2nSN5ooVGYWOt9+PcZTlFtrY1u8AHATlRfm90MAAAcSFW9N5PTLv95q2sBgFVg5hMAACyoqh6eyfWgXrPVtQDAqthvuVoAAGB/VXVukjskedSwaiEAsACn3QEAAADQjdPuAAAAAOhG8wkAAACAbnbcNZ+uf/3rt127dl1l2+c+97lc+9rXXjjHmPjtHLssdaxa7LLUsQyxy1LHMsQuSx3LELssdaxa7LLUsQyxy1LHMsQuSx3LELssdaxa7LLUsQyxy1LHMsQuSx3LELssdaxa7LLUsQyxy1LHwYg9//zzP9Fau8HcJ7XWdtTP3e52tzbrnHPO2W/besbEb+fYZalj1WKXpY5liF2WOpYhdlnqWIbYZalj1WKXpY5liF2WOpYhdlnqWIbYZalj1WKXpY5liF2WOpYhdlnqWIbYZalj1WKXpY5liF2WOg5GbJLz2hq9GKfdAQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Ry+1QUAAAAAsLx27d5z5e3Tjt+bU4f7l5x+ykLP13wCAAAA2GE221AaQ/MJAAAAYBs4lA2lMTSfAAAAAJbUsjaUxtB8AgAAADiE5jWUVqmZNJbmEwAAAMAm7WsorerspJ40nwAAAABmbIfT3ZaF5hMAAACwI2gobQ3NJwAAAGBlaSgtP80nAAAAYKloKG0vmk8AAABAdztthTe+RvMJAAAA2BArvLEIzScAAADgShpKHGyHbXUBAAAAAGxfZj4BAADANubi3Ww1zScAAABYMRpKrBLNJwAAAFgCGkpsV5pPAAAA0Mm8hpJmEjuN5hMAAACMoKEE41jtDgAAAIBuzHwCAABgx9s3m8m1luDg03wCAABg23Hxblgemk8AAACsBA0lWE2u+QQAAABAN2Y+AQAAsGXMZoLtT/MJAACAg0pDCZjmtDsAAAAAujHzCQAAgAOaN5vJTCZgEd1mPlXVTavqnKr6f1X19qp64rD9elV1dlW9e/jz2DWe/+gh5t1V9ehh2zWq6lVVdXFVPWEq9vlVdUKv1wIAALAd7dq9J7t278lFH77sytsAB1vP0+72JjmttXaHJPdK8hNVdYcku5O8urV2mySvHu5fRVVdL8lTktwzyT2SPGVoUt03yeuT3CnJo4bYOye5Wmvtgo6vBQAAAIAN6HbaXWvto0k+Otz+bFW9I8mNkzwkyUlD2J8lOTfJL848/b5Jzm6tfTJJqursJPdL8ukk10pyRJIaYp+R5Mc6vQwAAICV4ULfwDKq1lr/X1K1K8nrktwxyQdaa8cM2yvJp/bdn4p/UpJrttaeOdx/cpIvJPmdJC9O8g1JfjPJ5UlOaK099QC//3FJHpckxx133N3OOOOMqzx++eWX56ijjlr49YyJ386xy1LHqsUuSx3LELssdSxD7LLUsQyxy1LHqsUuSx3LELssdSxD7LLUsQyxy1LHqsUuSx3LELuVdVz04cuuvH3ckcmlX5jcPv7GR28qdjp+TOx0/JjYg1HzdhmLZR636fhlG7e1ci9D7HT8Th2Lk08++fzW2onzXkP35lNVHZXktUl+tbX2d1X16elmU1V9qrV27Mxz5jafWmvPnoo5IslZmcykelqSmyV5cWvtzPXqOfHEE9t55513lW3nnntuTjrppIVf05j47Ry7LHWsWuyy1LEMsctSxzLELksdyxC7LHWsWuyy1LEMsctSxzLELksdyxC7LHWsWuyy1LEMsVtZx+xspudcNDmBZd5spjGx0/FjYqfjx8QejJq3y1gs87hNxy/buK2Vexlip+N36lhU1ZrNp66r3Q0Nopcm+cvW2t8Nmy+tqhu21j5aVTdM8rE5T/1wvnZqXpLcJJPT86Y9IZNZUPdKclmSRyR5TZJ1m08AAABbzelxwE7Srfk0nFL3f5K8o7X2W1MPnZnk0UlOH/78hzlPPyvJs6ZWwrtPkl+ayn1skgdmcm2oByW5IklLcuRBfhkAAAAL0VACmK/nanf/M5MV6b6tqi4cfh6QSdPpO6vq3Um+Y7ifqjqxql6QJMOFxp+R5M3Dz9P3XXx88CuZnMZ3RSaNqm9NclGSP+/4egAAAAAYqedqd6/P11akm/Xtc+LPS/LYqfsvTPLCNXL/7NTtL2YyMwoAAOCgmjebyUwmgHF6znwCAAAAYIfresFxAACAZTO9IpXrMgH0Z+YTAAAAAN2Y+QQAAKw0q8wBLDfNJwAAYOloKAFsH067AwAAAKAbM58AAIBDwmwmgJ3JzCcAAAAAujHzCQAA2DCzmQA4EDOfAAAAAOjGzCcAAOAq5s1mMpMJgI0y8wkAAACAbsx8AgCAHWDfbCbXZQLgUDPzCQAAAIBuzHwCAIAVZTYTAKvAzCcAAAAAujHzCQAAlsS8VeYSs5kAWG1mPgEAAADQjZlPAADQkdlMAOx0Zj4BAAAA0I3mEwAAAADdOO0OAABGciodACzOzCcAAAAAujHzCQAAMn82k5lMALB5Zj4BAAAA0I2ZTwAAbFtmMwHA1jPzCQAAAIBuNJ8AAAAA6MZpdwAArJR9p9LtO40uiVPpAGCJmfkEAAAAQDdmPgEAsKXmXRQ8MZsJALYLM58AAAAA6EbzCQAAAIBunHYHAMBB51Q6AGAfM58AAAAA6MbMJwAAFmI2EwCwEWY+AQAAANCN5hMAAAAA3TjtDgBgB3MqHQDQm5lPAAAAAHSj+QQAAABAN067AwDYZuadSuc0OgBgq5j5BAAAAEA3mk8AAAAAdOO0OwCAFbDvVDor0gEAq8bMJwAAAAC6MfMJAGALzLsoeGI2EwCw/Zj5BAAAAEA3mk8AAAAAdOO0OwCAg8SpdAAA+zPzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuXPMJAGAdruMEALA5Zj4BAAAA0I3mEwAAAADdOO0OANhxnEoHAHDomPkEAAAAQDeaTwAAAAB047Q7AGBbmHcqndPoAAC2nplPAAAAAHSj+QQAAABAN5pPAAAAAHTjmk8AwNLadx2nfddwSuI6TgAAK8bMJwAAAAC60XwCAAAAoBun3QEAh8y+0+gSp9IBAOwUZj4BAAAA0I3mEwAAAADdaD4BAAAA0I1rPgEAm+I6TgAArMfMJwAAAAC60XwCAAAAoBun3QEA+3EqHQAAB4uZTwAAAAB0o/kEAAAAQDeaTwAAAAB045pPALBDuI4TAABbwcwnAAAAALrRfAIAAACgG80nAAAAALpxzScAWGHzruPkGk4AACwTM58AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuXHAcAJbMvouI77uAeBIXEQcAYGWZ+QQAAABAN5pPAAAAAHTTrflUVS+sqo9V1cVT2+5SVW+sqgur6ryquscaz310Vb17+Hn0sO0aVfWqqrq4qp4wFfv8qjqh1+sAAAAAYON6znx6UZL7zWz7jSRPa63dJcmvDPevoqqul+QpSe6Z5B5JnlJVxya5b5LXJ7lTkkcNsXdOcrXW2gV9XgIAHBy7du/Jrt17ctGHL7vyNgAA7ATdLjjeWntdVe2a3ZzkusPto5N8ZM5T75vk7NbaJ5Okqs7OpIn16STXSnJEkhpin5Hkxw5q4QAAAAAcNId6tbufSXJWVT07k1lX3zwn5sZJPjh1/0PDtr/NZMbTG5P8ZlU9OMkFrbV5DSwAAAAAlkC11voln8x8ekVr7Y7D/ecmeW1r7aVV9b1JHtda+46Z5zwpyTVba88c7j85yRdaa8+eijkiyVlJHpLkaUluluTFrbUz16jjcUkelyTHHXfc3c4444yrPH755ZfnqKOOWvh1jYnfzrHLUseqxS5LHcsQuyx1LEPsstSxDLHLUsfBjr3ow5clSY47Mrn0C5Ntx9/46IMWOx0/Jnat3MsQOx1vLIzFvPhDMRbLPG7T8cs2bmvlXobY6XhjYSzmxfts8dmykdjp+J06FieffPL5rbUT572GQ918uizJMa21VlWV5LLW2nVnnvPIJCe11h4/3P/jJOe21l4yFfPEJJcl+XAms6eemeQ1rbV7H6imE088sZ133nlX2XbuuefmpJNOWvh1jYnfzrHLUseqxS5LHcsQuyx1LEPsstSxDLHLUseBYqev2XTa8XvznIsmk4kvOf2UdeN7xU7Hj4ldK/cyxE7HGwtjMS/+UIzFMo/bdPyyjdtauZchdjreWBiLefE+W3y2bCR2On6njkVVrdl86nnB8Xk+kmRfg+jbkrx7TsxZSe5TVccOFxq/z7AtSTJse2CSF2dyDagrMrmW1JEd6wYAAABgA7pd86mqXpLkpCTXr6oPZbKC3Y8m+d2qOjzJFzOcCldVJyb5sdbaY1trn6yqZyR585Dq6fsuPj74lSS/2lq7oqrOSvITSS5K8rxerwUAAACAjem52t0j13jobnNiz0vy2Kn7L0zywjXy/uzU7S9mMjMKAAAAgCV0qE+7AwAAAGAH6TbzCQBWzeyFFE8d7q918UcAAODAzHwCAAAAoBvNJwAAAAC60XwCAAAAoBvNJwAAAAC6ccFxALY1FxEHAICtZeYTAAAAAN1oPgEAAADQjeYTAAAAAN245hMAK2fedZxcwwkAAJaTmU8AAAAAdKP5BAAAAEA3mk8AAAAAdKP5BAAAAEA3mk8AAAAAdGO1OwCWghXsAABgezLzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6MYFxwHoZt9FxPddQDyJi4gDAMAOY+YTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjQuOA7CwfRcQT1xEHAAAWIyZTwAAAAB0o/kEAAAAQDeaTwAAAAB0o/kEAAAAQDeaTwAAAAB0Y7U7gB3OCnYAAEBPZj4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3V7gC2ISvYAQAAy8LMJwAAAAC60XwCAAAAoBvNJwAAAAC60XwCAAAAoBvNJwAAAAC6sdodwIqwgh0AALCKzHwCAAAAoBvNJwAAAAC60XwCAAAAoBvNJwAAAAC60XwCAAAAoBur3QFsoXkr2Fm9DgAA2E7MfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALo5fKsLANhudu3ekyQ57fi9OXW4fcnpp2xlSQAAAFvGzCcAAAAAutF8AgAAAKAbzScAAAAAutF8AgAAAKAbzScAAAAAutF8AgAAAKCbw7e6AIBlt2v3nitvn3b83pw63L/k9FO2qiQAAICVYeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQzeFbXQDAVti1e8+Vt087fm9OHe5fcvopW1USAADAtmTmEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0M3hW10AwMGya/eeK2+fdvzenDrcv+T0U7aqJAAAgB3vgDOfqurrDkUhAAAAAGw/i5x298aq+tuqekBVVfeKAAAAANg2Fmk+3TbJ85M8Ksm7q+pZVXXbvmUBAAAAsB0csPnUJs5urT0yyY8meXSSN1XVa6vqm7pXCAAAAMDKOuAFx4drPv1gJjOfLk3yU0nOTHKXJH+b5BYd6wMAAABghS2y2t0bkvx5koe21j40tf28qnpen7IAAAAA2A4WuebTL7fWnjHdeKqq70mS1tqvd6sMAAAAgJW3yMyn3Un+ZmbbL2Vyyh1AV7t277ny9mnH782pw/1LTj9lq0oCAABghDWbT1V1/yQPSHLjqnru1EPXTbK3d2EAAAAArL71Zj59JMl5SR6c5Pyp7Z9N8rM9iwIAAABge1iz+dRae2uSt1bVX7bWzHQCAAAAYLQ1LzheVfuu8/SWqnrb7M+BElfVC6vqY1V18cz2n6qqd1bV26vqN9Z47v2q6l1V9Z6q2j21/S+H3/+sqW2/XFUPPVA9AAAAABx6651298ThzwduMPeLkvx+khfv21BVJyd5SJI7t9a+VFX/Y/ZJVXW1JH+Q5DuTfCjJm6vqzKHWL7TW7lRVZ1fV0UmuleSerbVnbrBGAAAAADpa77S7jw5//udGErfWXldVu2Y2/3iS01trXxpiPjbnqfdI8p7W2vuSpKrOyKRh9bIkR1bVYUmOSPLVJE9P8pSN1AcAAABAf+uddvfZqvrMnJ/PVtVnNvj7bpvkW6vq36vqtVV19zkxN07ywan7H0py49baO5J8PMkFSV6e5NZJDmutXbDBWgAAAADorFpr/ZJPZj69orV2x+H+xUnOSfLTSe6e5K+T3LJNFVFV353kfq21xw73H5XJqXU/OZP75Uken+SHk9w5ydmttT9Zo47HJXlckhx33HF3O+OMM67y+OWXX56jjjpq4dc1Jn47xy5LHasWuyx1bGXsRR++7Mrbxx2ZXPqF5PgbH71wbDI/fhlip+PHxE7HG4uNjcUyj9t0/LKN21q5lyF2Ot5YGIt58T5bfLZsJHY63lgYi3nxPlt8tmwkdjp+p47FySeffH5r7cR5r2HN5lNVXbe19pmqut68x1trn5z7xKvm2JWrNp9eleTXW2vnDPffm+RerbWPTz3nm5I8tbV23+H+Lw2/79emYh6S5K5J/irJ7tbaY6rqrCTf1Vr7/Ho1nXjiie288867yrZzzz03J5100oFezobit3PsstSxarHLUsdWxu7avefK26cdvzfPuejwXHL6KQvHJpkbvwyx0/FjYqfjjcXGxmKZx206ftnGba3cyxA7HW8sjMW8eJ8tPls2EjsdbyyMxbx4ny0+WzYSOx2/U8eiqtZsPq13wfG/yuRi4+cnaUlq6rGW5JbrPHctL0tycpJzquq2Sa6e5BMzMW9OcpuqukWSDyf5viTfv+/Bqjoiyc8kOSXJbYZakuRqQ751m08AAAAAHDrrXXD8gcOft9hI4qp6SZKTkly/qj6UyYXBX5jkhcPpd19O8ujWWquqGyV5QWvtAa21vVX1k0nOyqSh9MLW2tunUv9Ekj9rrX2+qt6W5FpVdVGSV7bWPr2RWgEAAADoY72ZT0mSqnp1a+3bD7RtVmvtkWs89INzYj+S5AFT91+Z5JVr5P2dqdstyVq/BwAAAIAttmbzqaqumeRamcxcOjZfO+3uupmsSAcAAAAA61pv5tPjM7m20o0yue7TvubTZ5L8ft+yAAAAANgO1rvm0+8m+d2q+qnW2u8dwpoAAAAA2CYOeM2n1trvVdUdk9whyTWntr+4Z2HAapleVvTU4fZaS4UCAACwcyxywfGnZLJq3R0yuQj4/ZO8PonmEwAAAADrOmyBmO9O8u1J/qu19sNJ7pzk6K5VAQAAALAtLNJ8+kJr7Yoke6vqukk+luSmfcsCAAAAYDs44Gl3Sc6rqmOS/Ekmq95dnuQNPYsCAAAAYHtY5ILjTxhuPq+qXpXkuq21t/UtCwAAAIDtYM3mU1WdsN5jrbUL+pQEAAAAwHax3syn56zzWEvybQe5FgAAAAC2mTWbT621kw9lIQAAAABsPwe85lNVHZHkx5P8f8Omc5P8cWvtKx3rArbYrt17rrx92vF7c+pw/5LTT9mqkgAAAFhBi6x290dJjkjyh8P9Rw3bHturKAAAAAC2h0WaT3dvrd156v5rquqtvQoCAAAAYPs4bIGYr1bVrfbdqapbJvlqv5IAAAAA2C4Wmfn080nOqar3JakkN0/yw12rAgAAAGBbWLP5VFWvTPJXSV6W5DZJbjc89K7W2pf6lwYAAADAqlvvtLs/TnJKkvcl+fMkt0ryTo0nAAAAABa1ZvOptfYPrbVHJtmV5KVJfijJB6rqT6vqOw9RfQAAAACssANecLy19vnW2l+31r4ryX2S3CXJq3oXBgAAAMDqO+AFx6vquCTfm+T7ktwwyd8kObVvWUAPu3bvufL2acfvzanD/UtOP2WrSgIAAGCbW++C4z+a5JGZXGj8pUl+vrX2b4eqMAAAAABW33ozn74pya8leXVr7YpDVA8AAAAA28iazafW2mMOZSEAAAAAbD8HvOA4AAAAAGyU5hMAAAAA3ax3wfHrrffE1tonD345AAAAAGwn611w/PwkLUnNeawluWWXigAAAADYNta74PgtDmUhAAAAAGw/6818ulJVHZvkNkmuuW9ba+11vYoCAAAAYHs4YPOpqh6b5IlJbpLkwiT3SvKGJN/WtTJgIbt277ny9mnH782pw/1LTj9lq0oCAACAKy2y2t0Tk9w9yX+21k5Octckn+5ZFAAAAADbwyLNpy+21r6YJFV1jdbaO5Pcrm9ZAAAAAGwHi1zz6UNVdUySlyU5u6o+leQ/exYFAAAAwPZwwOZTa+27hptPrapzkhyd5FVdqwIAAABgW1jkguM3m7r7/uHPr0/ygS4VAQAAALBtLHLa3Z4kLUkluWaSWyR5V5Jv7FgXAAAAANvAIqfdHT99v6pOSPKEbhUBAAAAsG0sstrdVbTWLkhyzw61AAAAALDNLHLNp5+buntYkhOSfKRbRQAAAABsG4tc8+k6U7f3ZnINqJf2KQdIkl2791x5+7Tj9+bU4f4lp5+yVSUBAADAhixyzaenHYpCAAAAANh+Fjnt7rZJnpRk13R8a+3b+pUFAAAAwHawyGl3f5vkeUlekOSrfcsBAAAAYDtZpPm0t7X2R90rAQAAAGDbOWyBmJdX1ROq6oZVdb19P90rAwAAAGDlLTLz6dHDnz8/ta0lueXBLwcAAACA7WSR1e5ucSgKAQAAAGD7WWTmU6rqm7P/ancv7lQTAAAAANvEAZtPVfXnSW6V5MJ8bbW7lkTzCQAAAIB1LTLz6cQkd2ittd7FwHa2a/eeK2+fdvzenLp7Ty45/ZQtrAgAAAD6W2S1u4uTfH3vQgAAAADYfhaZ+XT9JP+vqt6U5Ev7NrbWHtytKgAAAAC2hUWaT0/tXQQAAAAA29MBm0+ttddO36+qb0nyyCSvnf8MAAAAAJhYZOZTququSb4/yfckeX+Sl/YsCgAAAIDtYc3mU1XdNpMZTo9M8okkf52kWmsnH6LaAAAAAFhx6818emeSf0nywNbae5Kkqn72kFQFAAAAwLZw2DqPPSzJR5OcU1V/UlXfnqQOTVkAAAAAbAdrNp9aay9rrX1fktsnOSfJzyT5H1X1R1V1n0NUHwAAAAArbL2ZT0mS1trnWmt/1Vp7UJKbJHlLkl/sXhkAAAAAK2+h1e72aa19Ksnzhx/Y8Xbt3pMkOe34vTl1uH3J6adsZUkAAACwVA448wkAAAAANkrzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6ObwrS4Als2u3XuSJKcdvzenDrcvOf2UrSwJAAAAVpaZTwAAAAB0o/kEAAAAQDeaTwAAAAB0o/kEAAAAQDeaTwAAAAB0o/kEAAAAQDeaTwAAAAB0o/kEAAAAQDeaTwAAAAB0o/kEAAAAQDeaTwAAAAB0c/hWFwC97dq958rbpx2/N6cO9y85/ZStKgkAAAB2DDOfAAAAAOhG8wkAAACAbjSfAAAAAOimW/Opql5YVR+rqovnPHZaVbWquv4az310Vb17+Hn0sO0aVfWqqrq4qp4wFfv8qjqh1+sAAAAAYON6znx6UZL7zW6sqpsmuU+SD8x7UlVdL8lTktwzyT2SPKWqjk1y3ySvT3KnJI8aYu+c5GqttQs61A8AAADAJnVrPrXWXpfkk3Me+u0kv5CkrfHU+yY5u7X2ydbap5KcnUkT6ytJrpXkiCQ1xD4jyZMPZt0AAAAAHDyH9JpPVfWQJB9urb11nbAbJ/ng1P0PDdvOTrIryRuTPLeqHpzkgtbaRzqVCwAAAMAmVWtrTUA6CMmrdiV5RWvtjlV1rSTnJLlPa+2yqrokyYmttU/MPOdJSa7ZWnvmcP/JSb7QWnv2VMwRSc5K8pAkT0tysyQvbq2duUYdj0vyuCQ57rjj7nbGGWdc5fHLL788Rx111MKva0z8do5dljoOFHvRhy+78vZxRyaXfmFy+/gbH71ufK/Y6fgxsWvlXobY6XhjYSzmxR+KsVjmcZuOX7ZxWyv3MsROxxsLYzEv3meLz5aNxE7HGwtjMS/eZ4vPlo3ETsfv1LE4+eSTz2+tnTjvNRzK5tPxSV6d5PPDwzdJ8pEk92it/dfUcx6Z5KTW2uOH+3+c5NzW2kumYp6Y5LIkH07yzUmemeQ1rbV7H6imE088sZ133nlX2XbuuefmpJNOWvh1jYnfzrHLUseBYnft3nPl7dOO35vnXHR4kuSS009ZN75X7HT8mNi1ci9D7HS8sTAW8+IPxVgs87hNxy/buK2Vexlip+ONhbGYF++zxWfLRmKn442FsZgX77PFZ8tGYqfjd+pYVNWazadDdtpda+2i1tr/aK3taq3tyuR0uhOmG0+Ds5Lcp6qOHS40fp9hW5Jk2PbAJC/O5BpQV2Ry/agjD8HLAAAAAGCEbs2nqnpJkjckuV1VfaiqfmSd2BOr6gVJ0lr7ZCYXEn/z8PP0Yds+v5LkV1trV2TSlPrWJBcl+fM+rwQAAACAjTq8V+LW2iMP8PiuqdvnJXns1P0XJnnhGs/72anbX8xkZhQ7zOyUv1OH+2tNUwQAAAC2xiFd7Q4AAACAnUXzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6ObwrS4A9tm1e8+Vt087fm9OHe5fcvopW1USAAAAsElmPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQzeFbXQDb267de668fdrxe3Pq7j255PRTtrAiAAAA4FAy8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOjm8K0ugNWza/eeK2+fdvzenLp7Ty45/ZQtrAgAAABYVmY+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANCN5hMAAAAA3Wg+AQAAANDN4VtdAMth1+49SZLTjt+bU4fbl5x+ylaWBAAAAGwDZj4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdHL7VBdDHrt17rrx92vF7c+pw/5LTT9mqkgAAAIAdyMwnAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgG80nAAAAALrRfAIAAACgm27Np6p6YVV9rKountr2m1X1zqp6W1X9fVUds8Zz71dV76qq91TV7qntfzk891lT2365qh7a63UAAAAAsHE9Zz69KMn9ZradneSOrbU7JfmPJL80+6SqulqSP0hy/yR3SPLIqrpDVd0pyReG5969qo6uqhsmuWdr7WX9XgYAAAAAG9Wt+dRae12ST85s+6fW2t7h7huT3GTOU++R5D2ttfe11r6c5IwkD0nylSRHVtVhSY5I8tUkT0/ylE4vAQAAAIBNOnwLf/djkvz1nO03TvLBqfsfymR20zuq6uNJLkjy50luneSw1toF3StdErt277ny9mnH782pw/1LTj9lq0oCAAAAWFe11volr9qV5BWttTvObP/fSU5M8rA2U0BVfXeS+7XWHjvcf1QmzaefnIl7eZLHJ/nhJHdOcnZr7U/WqONxSR6XJMcdd9zdzjjjjKs8fvnll+eoo45a+HWNiT+YsRd9+LIrbx93ZHLpFya3j7/x0ZuKnY4fEzsd3yt2rTqWbSwOxbitlXsZYqfjjYWxmBfvs8Vny0Zip+ONhbGYF++zxWfLRmKn442FsZgX77PFZ8tGYqfjd+pYnHzyyee31k6c9xoOefOpqk7NpGn07a21z895zjcleWpr7b7D/V9Kktbar03FPCTJXZP8VZLdrbXHVNVZSb5rXs5pJ554YjvvvPOusu3cc8/NSSedtPDrGhN/MGNnZz4956LJxLV5M5/GxE7Hj4mdju8Vu1YdyzYWh2Lc1sq9DLHT8cbCWMyL99nis2UjsdPxxsJYzIv32eKzZSOx0/HGwljMi/fZ4rNlI7HT8Tt1LKpqzeZTzwuO76eq7pfkF5I8eJ0m0ZuT3KaqblFVV0/yfUnOnMpxRJKfSfIbSY5Msq97drUkV+9UOgAAAAAb0K35VFUvSfKGJLerqg9V1Y8k+f0k10lydlVdWFXPG2JvVFWvTJLhguQ/meSsJO9I8jettbdPpf6JJH82NK/eluRaVXVRkvNba5/u9XoAAAAAGK/bBcdba4+cs/n/rBH7kSQPmLr/yiSvXCP2d6ZutyTzfg8AAAAAS+CQnnYHAAAAwM6i+QQAAABAN5pPAAAAAHSj+QQAAABAN5pPAAAAAHSj+QQAAABAN5pPAAAAAHSj+QQAAABAN5pPAAAAAHSj+QQAAABAN5pPAAAAAHSj+QQAAABAN5pPAAAAAHRz+FYXsNPt2r3nytunHb83pw73Lzn9lK0qCQAAAOCgMfMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG40nwAAAADoRvMJAAAAgG4O3+oCtqNdu/dcefu04/fm1OH+JaefslUlAQAAAGwJM58AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6EbzCQAAAIBuNJ8AAAAA6ObwrS5gVezavefK26cdvzen7t6TS04/ZQsrAgAAAFh+Zj4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0I3mEwAAAADdaD4BAAAA0M3hW13AVtq1e0+S5LTj9+bU4fYlp5+ylSUBAAAAbCtmPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN1oPgEAAADQjeYTAAAAAN10az5V1Qur6mNVdfHUtutV1dlV9e7hz2PXeO6jh5h3V9Wjh23XqKpXVdXFVfWEqdjnV9UJvV4HAAAAABvXc+bTi5Lcb2bb7iSvbq3dJsmrh/tXUVXXS/KUJPdMco8kTxmaVPdN8vokd0ryqCH2zkmu1lq7oNNrAAAAAGATujWfWmuvS/LJmc0PSfJnw+0/S/LQOU+9b5KzW2ufbK19KsnZmTSxvpLkWkmOSFJD7DOSPPngVg4AAADAwVKttX7Jq3YleUVr7Y7D/U+31o4ZbleST+27P/WcJyW5ZmvtmcP9Jyf5QpLfSfLiJN+Q5DeTXJ7khNbaUxeo43FJHjfcvV2Sd82EXD/JJ0a8tDHx2zl2WepYtdhlqWMZYpeljmWIXZY6liF2WepYtdhlqWMZYpeljmWIXZY6liF2WepYtdhlqWMZYpeljmWIXZY6liF2WepYtdhlqWMZYpeljoMRe/PW2g3mPqO11u0nya4kF0/d//TM45+a85wnJfnlqftPTvKkmZgjkrwmyXWS/FaS/5vkwZuo87xe8ds5dlnqWLXYZaljGWKXpY5liF2WOpYhdlnqWLXYZaljGWKXpY5liF2WOpYhdlnqWLXYZaljGWKXpY5liF2WOpYhdlnqWLXYZaljGWKXpY6er6+1dshXu7u0qm6YJMOfH5sT8+EkN526f5Nh27QnZDIL6l5JLkvyiCSnHfRqAQAAANiUQ918OjPJo4fbj07yD3Nizkpyn6o6drjQ+H2GbUmSYdsDM2k+XSvJFUlakiM71g0AAADABnRrPlXVS5K8IcntqupDVfUjSU5P8p1V9e4k3zHcT1WdWFUvSJLW2iczuZD4m4efpw/b9vmVJL/aWrsik6bUtya5KMmfb6Lc53eM386xy1LHqsUuSx3LELssdSxD7LLUsQyxy1LHqsUuSx3LELssdSxD7LLUsQyxy1LHqsUuSx3LELssdSxD7LLUsQyxy1LHqsUuSx3LELssdfR8fX0vOA4AAADAznaoT7sDAAAAYAfRfAIAAACgG80nAAAAALrRfAIAAACgmx3ffKqqo6rqhKo65iDnvUVVPayqbn8w8/bMPXYsquqEHnnHGDMWVXXdqrpbVR17MPOOMTZvVV1/RO6Fxrmqjq2q647I2+29DGxOVR037PcnVNVxI563of26qo46wOPfUlU/V1X3GZN3gd87Km9VXW9E7rGfyw/ulHehmseMRVXduqoeXlV32GzezfwbvkDuw6duHzWsgnzA8VhkjKvqBlV116q604Hev4vWOxN7vZHvt5U67ux8PLtSx529jjnH5h6j13Fnr2PO4TmL7Nej9rvhOVv+vliGvENu+94Gco+xjfe9DR13Xqm1tqN+kvzh1O1vSfKBJOck+WCSB4zIc9HM/ZdN3X5Ikvcn+dMk70py6pzn3zTJGUn+Jcn/SnLEvFxjcye5fZJ/TLInya2SvCjJp5O8Kck3bHQskpww83O3JB9KctckJ2wi72Ombt8kyauHev8tyW3njNuYsfiLJNcfbt93qOOfk/xnku/ZRN5PJnlBkm/PsGLkOu+TMXnvPzz++mFc357kvcM4f/tm3stJbpTkxUkuS/LVIfYDSZ46/d7byHt5Kva4qffHcQvuR7dI8rAktx+5Hx91gMe/JcnPJbnPmLwL/N7ReZNcr8dYJHlwp7yL1rvwWCS5dZKHJ7nDZnMnOabX31+Sw6ffY0lOXGQ8FhnjJDcY9us7Hej9u0jNSe6S5I1J3pHJ59o/J3nnsO2EOXleNnV74f16Tp4PzNx/09TtH01yYZKnJPnXJLtnYo8f6vtgJsvzHjsvzwby/s9hHN6e5J5Jzs7ks/ODSb5pM2Mx/L1O/zw8yX/tu7+JvL88dfsOSf5jiL8kyT03MRbn5Gv/7j1qyPuCJBcl+amN5h1i9g7vsx/JAfbDkTWfmuS/h1rvn+R9mRwLfDDJIzcxxncY6n1Pki8n+fch/kVJjt7kWNwsk2O4jyd59/A7PjZs2zUTu+XHnel0zDnEbNvjzpFjvPAx5wZyb/lx58gxXviYc0zNGbHfbfd9LyP2O/uefS+b3PeGx++SEceda47TooHb5SfJBVO3z9k3WElumeS8mdjZg87pg8+Pz8S+Zer2vyW5xXD7+kneOqeOs5P82PAX+XvDc75uNtfY3Elel+RBSR457HTfl6SGba/exFhcMfzuc6Z+vjD8+ZpN5J2O/Zskj8tkRt53zda7gbG4aCZ21zqxY/K+K8lPZnJQ+uEkv5vkXmu838bkvTDJNyT5pkwOxO81bP+G6XHa4Di/JslJU+/r305y7STPTPL8Tb6X75IFP4yyjb8ADzELfwkeMxbZxl+AN5B7234BHlPzsP2ec55/r8zfTxferzNpds37OS3JJ9fJ++YkNxhuXzv7H7C/Psn9khyT5EmZ7Ce3ms2zgbxvymS//qYkn0jyLcP2E5L86ybH4itJXpHkhcPf8Z8m+ezw5ws3kXf683tPkvsPt++R5N82MRYXz8TuO664VpK3bTTvsP2iJA9M8peZ7Cv/kMkxxpEHGOMD1XzRMEa3SPKZqffEcQeo+UBj/MYkt5sa1z+b2q/+7ybH4g1JHpHkalPbrjaMxxvX+bs+0L/VXY470+mYc9i+bY87R47xwsecG8i95cedI8d44WPOMTVnxH633fe9jNjv7Hv2vWxy35uqeeHjzrV+FgraTj8zf4Hnr/XYcP8rmXxB+NM5P59dJ+/sF9i3zKnjwpn7P5jhQPwAb7p1c8+8id5zgNc3ZiwenuS1GQ6Qh23vPwhjPB07Oybzxm3MWLw9yXWH269Pctj0Y5vIOx17syS/kOSCTL6sPusg5f3geu+XDYzz7AfI+VO337nRmvfVlgU/jDLuQ26lvgDvG68s+CV45Fhs2y/AG8i9bb8Aj6k5ybtnnzv12HvmbBvzWfTFJM/IpOk1+/Ppmdi3Jjk2yddl/4Og2byzY3NyJv97fa/M+cwakXd6zN6x1uve4FjcPZMm5I9PbXv/GuO+0c/7NV/PRsYiyY2H2+ckueZw+2rZ/9+9hfPOqfnIJN+b5O8y2Q//ahM1Xzh1+yMzj83ue2PGeL3Pu9n3ydixWG//e/c6v3dLjjvT6Zhzdlu22XHnyDFe+JhzA7m3/Lhz5BgvfMw5puYx+90Gal6pfW/MfreBsbDvzY/dsfvecH/UcedaP1eeY7+D3L6q3pZJd3hXVR3bWvtUVR2W5OozsW9L8uzW2sWzSarqO2Y23bmqPjPkvUZV3bC19tGqunomB32zjqiqa7bWvpgkrbW/qKr/SnJWJl8yNpp7+v5vzTw2+/oWHovW2kur6qwkz6iqx2TSAGhzXteovEluUlXPHWJvUFVHtNa+sm+M5uQeMxZPS3JOVf1BJh3rv62qMzP5wvOqTeStqXH5QJLfSPIbw/mxj9hE3k9X1eOTXDfJp6rqZzP5n4HvSHL5nLEYM84fr6ofzOTLyMMymdmSqqrsf+23se/la7fW/n12Y2vtjVU1+16efs8c3lp7/xD7iaq6Yib2WUl+M5OZLrNmaz5sOLf7sEymxX58yPu5qpp9/nVaa/v+/p9dVecneVVVPSr7v6fH5E0mU1ovSpKq+nhr7fVD/AVVdeRM7Jix+OYkpyd5c2vtj4b8J7XWfnhODWPyTrtRa+0fh9g3zal3zFh8papu3Fr7cCbv3c8N27+U+e+hUblba69I8oqhxgdl0nz6g6o6q7X2/RvM+9XW2ieSfKKqLm+tvXeIvXSym1zFmDE+srX2ruHxN1XV84bbf1JVP7eJsfjHqtqTyfTqDw7bbprkh7L/51sybr++IJPZXefPJqmqx85sOjrJ+UPeNpX3qEx9Vk49/+jW2mXDazqnqh6e5KVJrreJvNOfB78089jsZ2EyYixaa2+uqu9M8lNVdU6SX8za//aNGeNbDv8eVSb/Dl6rtfb54bHZf/vGjMXPJvmnqnppJgfkrxn+7f6WTL5EbTRvctV/+76Qyb9Pf1NVRyd56CZyf6Cqfi3JdZK8s6qek0lT6zuSfHQmdswYv7eqnpzJ/wQ/LJP/KElVHZH9/w0ZOxbnV9UfJvmzXHX/e3QmDcBpy3Dc2euYM9nex529jjnH5l6G485ex5xjah6z342tedX2vTH7XWLf22hu+97XjD3unG/RLtV2+Uly85mfqw/br5/9T1/51iQ3WyPPiQv+vmMy/5oTP5vk3nO23zXJ2RvNneTxmXM9kUyuufI7Gx2LOTWek+RjB2GMHz3zc+yw/esz01He4DjfOsmvJ/n7JC9P8kdJ7ruZvEl+6yC8D+flvWmSPx5q/PrhPXJxJjNS5p2/PWacb5bJh9rFmZwbfcNh+9clefgmx/i5Q42PyKRJ8s3D7T1Jfn8m9quZzCr5bCanIO2r4+rZ/3+4/y3J3daoZfZ/CS7J5H8h3j/8uS/vUdn/fxDemv2v93GnTGZg/PdG8+7LPXX7oTOPXTxzf+GxGLYfluSJmex790jyvjXGZswYfzrJmcO+8fEk11qn3jFjfFImX3yfnuT3h7/Lp2Qy/fxJc2oek/sta7zuo5M8ehN5z0zya0O9r0nynExOo3xKkrM2McZ/l+TJQ67nZJiplsmB1rs2ORb3T/K84e/v5cPtha9jsdZ+neR2GU6bnBO/6PXcrpVhRtjUtu/PnOnqmXw+/ckm8j54+r07tf1WSX5hM2Mx8/iNMvkcnbvvjRzje8/8HLVvfJP8xEbHYmpf+PFMptv/XiYNs4Wvq7dO3v323bE/a/z9XTeTpuHu4X3+8Exmev7Bvvf/Bsf4mEy+KLwiya9m8h8P+8Zn7mkTI8bi6sMYvyqTWZMXZXL9lSckucZM7Jh/q7scd6bTMeewfccdd64zFps65lzn72/LjztHjvGmjznn1Txmv9tAzSu1743Z78aOxZw67Xub+NkO+97U9k0fd9aQaEeqYZWE1tont7qWVTN0UK/TWvvMgvEntNYu6FzWQVOT1QFuk8mXjE9tdT2Lqqrrt8nsjYOV79hMZoSs+/dcVffP5No3Nx42fTjJma21Vy74e47J5MP2DVPbbpdJM2i/11NVx7XWLl0g77Uy+bL8/qlt35/J3+sbZ2JvluTJrbUf3UjeYfuDk/xz+9oMhn3bb5XJh/5vLJD7mMyMxczjN0ryO5kcDN3yQPnWy1tV954JO7+1dnlNVq/47tbaHyyQd62xODqTZsNtkxyeyUUU/6G19s4RNc/7+3tSa+3Zi+YYkfe6SX4ik//d+/1MLlr5w5lcS+GZrbXZGRjz8h6T/cf4mEwuMnqHTBqfp7fWPjuMzzfMvg/H1LwRw9/tlfvpIvvRJn7XUa21ebM2l05VXa/HsUBVPbi1dmaHvAe93qq6dZI7Z3JK2v/bZK5jWmufPiiF7Z/78Nba3uH2UZlcePd9B2M8quoGmVwI96tDzoP6/nXcuXGOO5fPwT7mHHIudNy5gbz2vQ2y7y2fVdr39rPZbt6q/WTkSgnr5HngiNj9rg0zJndGrJQwMu+olRJG5D0hC66SMHaMM3J1vBE1j1opYUTehVdJ2MBYjFodb0TNo1ZK2ODrGb0y3iZ+18Iri231TxZcaW4DeRdaGW+J6h21Ot46eY7p+He1oZXxRuTf0Op4w3MfN2fbXbP2wgB33UzudWI/sNG8w+teaGGADYzFqNXxRtT8sOx/odq5CwOMzLvwwgAj856TEQsDjMy98MIAI/OemgUXBhiZd9TCAAfIPftv6lIfd86p96Acc66Re6WOO9PpmHPYtlLHnel0zDls2/Rx5xp5d+S+t8ZY2Pe+ts2+97VtB+U7X8YcG250YFb1JyNXSlgnz9Nm7s8edE4ffH58ZI2zuRdeKWFk3tdlxEoJI/IuvErCmLzDtlGr442oedRKCSPyLrxKwgbG4sKMWB1vRM2jVko4QO7ZA/wd+QV4jdy+AH9t20H5Ejyn5pX6Ajw1tpv6Epzk8XO2XTjv7yljVymZyZ0RCwOMzLvwwgAbGItRq+ONqHnhhQFG5l14YYCReUctDDAy98ILA2wg70ILA4zMO2phgAPknv03damPO+fkPSjHnGvkXqnjznQ65tz3Xp66vfTHnel0zDls2/Rx5xp5d+S+t8ZY2Pe+ts2+97VtB+U7X+Yca60ZO3ZAVv0n41dKuH0m10t47vDzi5l//Z2FV0nYQO4LZ+6vuVLCyLxvmbq9yEoJi+ZdeJWEMXln65ozLm/ZRM1jV0pYNO/CqyRsciwWWR1v0ZpHrZRwgP1t9gD/wuzAL8Br5PYF+GvbDsqX4Dk1r9QX4GHbmOXhb5/J/7AdNbP9fnPyjl0db6HcGbEy3si8C6+Mt4GxeMvU7UVWx1u05oVXxhuZd+GV8UbmfUsWXBlvkzWvuzLeyLwXTt1ed2W8Tb7f1lwZb2r7PZLcfbh9h0z+zdnv2hcH2PfmHXcumnfUceeIvBfO3F/3mHNk7rdM3V7kuHPRvGOPOxfNO/aYc6G8w+Njjzs3UvMix50bybvIMeeYsdjwcWeSF6/z2Kh9b9HcY/e9EXlH73sjxuItU7cPuO+NqHnUvjci76h9b+RYjNr3NljzAfe9DeY94L43ciw2s+99y7Bf32fRv+/WduZqdwuvlFBVv5hJh/iMTL4oJpOpfy+pqjNaa6dPhb8ti6+SMDb3wisljMy78EoJY/K2EaskjKw3GbFSwsjcC6+UMDLvwqskbGAsFl4pYWTusSslrOfLM/fHrIw3Ju+YlfHG5B2zMt7Y3GNWxhuTd8zKeGPyTjvQynhj845dHW/R3GNWxhuTd8zKeGPyJguujldVP53JNarekeT/VNUTW2v/MDz8rOy/8sjCq5SMzL3wynhja64FV8bbwFhMfyasuzremNxtxMp4I2teeGW8kXnHrIw3NvfCK+ONzLvwyngj845ZGS9V9ZRMZj0eXlVnZzJ79Zwku6vqrq21X50KH3PcOSbvwsedI/OOWZ1rbO4xx50L5x153Dmm3jHHnGPyJuOOO8fkHnPcOSbvmGPOsWOx0HHnMD5X+VVJTq7JtRXTWnvwzONj9r0xucfse2Pyjvm+N3Ysxux7C+ceue+NqXnMvjd2LMbse2Nyj9n3xuQds++NHYuFv/NV1Ztaa/cYbv9oJv++/n2Sp9TkOl+z31PnG9Op2g4/GbdCyX9kzvmOQ453z2wbtUrCyNwLr5QwMu+YFUoWzjunxnOy9ioJo/JmxEoJG8i90EoJI8d44VUSNlDvmJUSxtR8UFZKGJ7zgZn7z82CK+ONzLvwyngj8y68Mt5Gck/dfujMYxdvNO+w7bAssDLeyHo/nQVXxttAvSdlxOp4I2p+yxpxR2dmZbyRec/MgivjbWAsFlodL5N/u/atkLYryXlJnniA173QKiVjcmfEyngj8y68Mt7YsciI1fE2Ms7DY+uujDdyLO6dBVfG28BYHJ0FV8YbWfOY/XZM3oVXxhuZ95iMWBlvyH21TGZnfiZf+1/0I7P/ipdjjjvH5F34uHNk3lGrc43MPea4c+G8c+pc77hzTL1jjjlH15vFjzvH1DzmuHNM3jHHnKPGIgsed2bynx1/kckxw72HPz863J73nh2z7y2cO+P2vTF5x3zfGzsWY/a9UblH7HtjxmLMvje63iy+742pecy+NybvmH1v7Pti4e98uersuTcnucFw+9oZc62zRQN34k8m16K5+ZztN8+cJbKXIfcy5s2k63pdY7z93m9Djret8XNRki/NiV/0C/DCeTPuC/CYvKOWhh+Ze8wX4FFjPPW8A30BHlPvvWd+1vsCPLreLPgleGTNY74Aj8k75gvw2P3jmCzwJTgzU8OHOl6Vyf9oXrjJfbpL7lXLu4o1G4tDNhZvmXd7uL+Zmlcq77LWnPWPO5eu3lUc4y3Ke1gmX7zPTnKXYdvo/1w7lLlXLe9mcx9g39tRY7Gd8g553prk2EwaU+fNPPaWhfMcjGK2y0/2v7r8/TK56Os/ZnKx4edncuDynsy5jsSieXvm3ql5V7Hm7TAWw7ZLM7k44s1nfnZl5pocI1/DSuVdxZpXLe8q1twx72syHFhMbTs8k9Pqvjoy1+wF1Q9K7lXLu4o1H6q8q1hz57H49wz/eZCrXi/k6Iy7AOzsv9UrlXcVa172vKtY86Eci2HbTZL8bSazjxdeVOZAeQ9G7lXLu4o1G4utGYtMTsl7XyaL4LwvX5sldVRGNJV34jWf1nP3TP63OUnSWntVVd02k9NWbjxs/nAm11L56kbz9sy9g/OuYs0rPxaDV2QyC+bC2eCqOnfRxFX1uNba81c4b8/cOzJvz9yrlndO7h/KzDXOWmt7k/xQVf3xmLyZuk7BQc69anl75l71vD1zr1reJPn/WmtfGnJeMbX9iExOEVnU7L+pq5a3Z+6dmrdn7lXLOy93WmsfSvI9VXVKJqf1bcS89/LByL1qeXvmXrW8PXOvWt79crfWdq0Rd0UmqxAupIaO1Y4yXPjrIbnql/AzW2vvWMa8PXOvWt6euVctb8/cPWte53c+vrU29ovitsvbM/eq5e2Ze9XyHozcw3594yT/3lq7fGr7/drXLqy/7fOuYs2rlncVa+48FvdI0trkovR3yGQ28jtba6/cSXlXsWZj0T9v79wzv+fFrbUfOpg5VzFvz9yrlrdn7lXL2zt3kp0386nGryi2pXlXsWZj0T/vqtZ8AOutrraT8vbMvWp5e+Zetbybyl3jV4TblnlXseZVy7uKNXcei6dk3Mpf2zLvKtZsLPrn7Vzz2JW/tmXenrlXLW/P3KuWt3fuNbUNnAO4yj/Z4KptW5V3FWs2FsZig793Q+ckb7e8q1izsVj+scgGV23bbnlXseZVy7uKNR+CsRi9att2y7uKNRuLlR+LDa3att3yrmLNxmK1x2K9nx038ymT8xJvlOQ/Z7bfcHhs2fL2zL1qeXvmXrW8PXN3q7mq3rbWQ5msmrYj8vbMvWp5e+Zetbydcx/WhlOJWmuXVNVJSf5vVd18yL1T8vbMLW//3KuWN0n2tsm1FT9fVe9trX1m+D1fqKrN/Ju6anlXsWZj0T9vz9wnJnlikv+d5OdbaxdW1Rdaa6/dZL2rlrdn7lXL2zP3quXtnXuundh8+pkkr66qdyf54LDtZkluneQnlzBvz9yrlrdn7lXL2zN3r7zJ5MvzfZN8amZ7Jfm3HZS3Z+5Vy9sz96rl7Zn70qq6SxsufN5au7yqHpjkhUmO30F5V7HmVcu7ijX3HIsvV9W1WmufT3K3fRur6uhs7j90Vi1vz9yrlrdn7lXL2y13m1zs/Ler6m+HPy/NQfjuu2p5e+Zetbw9c69a3t6517Ljmk/t4K4o1j3vKtZsLPrnXdWacxBX/lrxvD1zr1renrlXLW/P3AdzRbhVztszt7z9c69a3uTgrvy1ynl75l61vD1zr1re3rnTDs7KXyuft2fuVcvbM/eq5e2de9aOXO0OAAAAgEPjsK0uAAAAAIDtS/MJAAAAgG523DWfAADWUlVfl+TVw92vT/LVJB8f7t+jtfblQ1DD05O8rrX2z71/FwDAoeCaTwAAc1TVU5Nc3lp79gKxhw8Xqd4WttvrAQC2ltPuAADWUVV3q6rXVtX5VXVWVd1w2H5uVf1OVZ2X5InD/d+uqvOq6h1Vdfeq+ruqendVPXN4zrWrak9VvbWqLq6qR8z5fS+qqu8ebl9SVU+rqguq6qKquv2c+NdV1V2m7r++qu48/K4XVtWbquotVfWQ4fFdVfUvQ84Lquqbh+0nDdvPTPL/FqkVAGARTrsDAFhbJfm9JA9prX18aMD8apLHDI9fvbV2YpJU1YOSfLm1dmJVPTHJPyS5W5JPJnlvVf12kpOSfKS1dsrwnKMXqOETrbUTquoJSZ6U5LEzj/+fJKcm+Zmqum2Sa7bW3lpVz0rymtbaY6rqmCRvqqp/TvKxJN/ZWvtiVd0myUuSnDjkOiHJHVtr76+qh2+gVgCA/Zj5BACwtmskuWOSs6vqwiS/nOQmU4//9Uz8mcOfFyV5e2vto621LyV5X5KbDtu/s6p+vaq+tbV22QI1/N3w5/lJds15/G+TPLCqjsikKfaiYft9kuwe6j43yTWT3CzJEUn+pKouGp57h6lcb2qtvX/qNYytFQBgP2Y+AQCsrTJpIn3TGo9/bub+l4Y/r5i6ve/+4a21/6iqE5I8IMkzq+rVrbWnH6CGfXm+mjnHbq21z1fV2UkekuR7M5ltta/2h7fW3nWVFzS5ltWlSe6cyX9EfnHe69lgrQAA+zHzCQBgbV9KcoOq+qYkqaojquobN5qsqm6U5POttb9I8puZnOZ2MLwgyXOTvLm19qlh21lJfqqqavjddx22H53ko621K5I8KsnVDnGtAMAOY+YTAMDarkjy3UmeO1zz6PAkv5Pk7RvMd3yS36yqK5J8JcmPH4wiW2vnV9Vnkvzp1OZnZFLr26rqsCTvT/LAJH+Y5KVV9UNJXpX9Z291rRUA2HmqtbbVNQAAsAnDLKVzk9x+mNEEALA0nHYHALDChhlM/57kf2s8AQDLyMwnAAAAALox8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOhG8wkAAACAbjSfAAAAAOjm/wfk5THBOO8h2wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1440x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from math import floor, ceil\n",
    "\n",
    "terms = [x / 100 for x in range(5, 505, 5)]\n",
    "term_vols = [vol_tracker.get_annual_term_volatility_forecast(term) * 100 for term in terms]\n",
    "series = pd.Series(term_vols, index=terms)\n",
    "\n",
    "ax = series.plot(figsize=(20,10), grid=True, kind='bar', title='Term structure of S&P 500 volatility')\n",
    "ax.set_ylim(bottom=floor(series.min()), top=ceil(series.max()))\n",
    "ax.yaxis.set_major_formatter(mtick.PercentFormatter())\n",
    "ax.set_ylabel('Annual Volatility')\n",
    "_ = ax.set_xlabel('Terms in years')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "renewable-botswana",
   "metadata": {},
   "source": [
    "## Pricing using Black-Scholes-Merton Pricer\n",
    "The Black-Scholes-Merton model can be used to price European options, including equity indices. With enough customizations it can also be used to price American call options (even those paying dividends or having a known\n",
    "dividend yield). Exchange and over-the-counter traded equity index options tend to be European.\n",
    "\n",
    "Let's put it to use to price an S&P 500 call and put options with a strike of 3900. Some of the stocks comprising S&P 500 pay dividends. It's essential to take into account the expected dividend yield of S&P 500 to achieve accurate pricing. Where do you get S&P 500's dividend yield? I calculate it based on forward dividend yields of individual S&P 500 companies in [this notebook](https://github.com/ilchen/US_Economic_Data_Analysis/blob/main/Stock_Market.ipynb). By running that notebook (please clone the corresponding repository), you will be able to obtain the forward dividend yield for the last business day. The class I use to calculate it [is here](https://github.com/ilchen/US_Economic_Data_Analysis/blob/main/stock_market/metrics.py#L13-L217). Alternatively, at the start of every month Aswath Damodaran @AswathDamodaran calculates the so-called equity risk premium (ERP) for S&P 500 and publishes it for download as a spreadsheet from [his website at STERN](https://pages.stern.nyu.edu/~adamodar/). `Spreadsheet to compute current ERP for current month`. From the date on S&P 500 dividends in tab `Buyback & Dividend computation`, you can\n",
    "calculate the current dividend yield as well as it's estimate for the next months.\n",
    "\n",
    "Market makers tend to calculate the life time of options using trading days rather than calendar days. The `pricing.options.BlackScholesMertonPricer` class supports both ways of pricing options. It's determined by the `holidays` parameter of its constructor, passing an array of `datetime64[D]` values implies calculating the life of\n",
    "the option using the number trading days to its expiry and dividing by 252. Not passing a value for this parameter\n",
    "will imply calculating the life of he option using the number of calendar days to its expiry and dividing by a value between 365 and 366. For pricing options on SP500, whose underlying futures are traded allmost on all non-weekend days, whether to use a custom trading calendar or not doesn't make a big difference on the pricing. Therefore in this notebook I use calendar days for calculating the lifetime of options. If you want to price with it, please make use of the holidays calendar created in the below cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "conscious-payday",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Taking the trading calendar into account when pricing\n",
    "import pandas_market_calendars as mcal\n",
    "holidays = mcal.get_calendar('NYSE').holidays().holidays"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "alone-deployment",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1219.21, σ: 0.1456, Δ: 0.973, Γ: 0.000, ν: 198.128\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 4.39, σ: 0.1456, Δ: -0.017, Γ: 0.000, ν: 198.128\n"
     ]
    }
   ],
   "source": [
    "strike = 4500.\n",
    "\n",
    "# Expected S&P 500 dividend yield\n",
    "q = 0.013891537607483252 # Updated for 9th April 2024\n",
    "\n",
    "cur_price = asset_prices[-1]\n",
    "pricer = options.BlackScholesMertonPricer(maturity_date, vol_tracker, strike, curve, cur_price,\n",
    "                                          ticker=TICKER, q=q)\n",
    "print(pricer)\n",
    "\n",
    "put_pricer = options.BlackScholesMertonPricer(maturity_date, vol_tracker, strike, curve, cur_price, is_call=False,\n",
    "                                              ticker=TICKER, q=q)\n",
    "print(put_pricer)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "knowing-mumbai",
   "metadata": {},
   "source": [
    "Sometimes you might want to price using _an implied volatility_ for the expiry date, implied volatilities are volatilities implied by option prices of market makers. You can easily use implied volatilities with \n",
    "the `BlackScholesMertonPricer` class, just pass the volatility value as the second parameter instead of\n",
    "a volatility tracker (I also use a custom trading calendar in the below cell):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "favorite-hamilton",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1245.56, σ: 0.2100, Δ: 0.925, Γ: 0.000, ν: 593.230\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 30.13, σ: 0.2100, Δ: -0.066, Γ: 0.000, ν: 593.230\n"
     ]
    }
   ],
   "source": [
    "impl_vol = .21\n",
    "\n",
    "pricer = options.BlackScholesMertonPricer(maturity_date, impl_vol, strike, curve, cur_price,\n",
    "                                          ticker=TICKER, q=q, holidays=holidays)\n",
    "print(pricer)\n",
    "\n",
    "put_pricer = options.BlackScholesMertonPricer(maturity_date, impl_vol, strike, curve, cur_price, is_call=False,\n",
    "                                              ticker=TICKER, q=q, holidays=holidays)\n",
    "print(put_pricer)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "subtle-question",
   "metadata": {},
   "source": [
    "#### Put-Call Parity\n",
    "An important theoretical property of the price of a European option is the so-called put-call parity, which is expressed as:\n",
    "$$c - p = S_0\\cdot e^{-q\\cdot T} - K\\cdot e^{-r\\cdot T}$$\n",
    "where:\n",
    "* c is the price of a European call option\n",
    "* p is the price of the corresponding European put option\n",
    "* S<sub>0</sub> is the present price of the underlying stock\n",
    "* q is the expected dividend yield of the index\n",
    "* K is the strike price\n",
    "* T is the liefteime of the option expressed in years\n",
    "* r is the riskless rate for maturity T\n",
    "\n",
    "Let's put it to the test:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "interesting-province",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Put-Call parity: 1215.43 = 1215.43\n"
     ]
    }
   ],
   "source": [
    "from math import exp\n",
    "price_difference = pricer.get_price() - put_pricer.get_price()\n",
    "\n",
    "right_hand_side_expr = pricer.s0 * exp(-pricer.q * pricer.T) - pricer.strike * exp(-pricer.r * pricer.T)\n",
    "\n",
    "assert np.allclose(price_difference, right_hand_side_expr, .1)\n",
    "print(f'Put-Call parity: {price_difference:.2f} = {right_hand_side_expr:.2f}')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "annoying-acrylic",
   "metadata": {},
   "source": [
    "## Pricing using Binomial Tree Pricer\n",
    "The Binomial Tree model can be used to price any European or American option, including on index some of whose constituent stocks pay dividends.\n",
    "\n",
    "Let's put it to use to price an S&P 500 call and put option with a strike of 3900 and compare resulting prices with those obtained from the BlackScholesMertonPricer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "surrounded-editor",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1244.56, σ: 0.2100, Δ: 0.926, Γ: 0.000, ν: 517.340\n",
      "Price difference between Black-Scholes-Merton and Binomial Treee pricers for calls: 0.998806\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 29.73, σ: 0.2100, Δ: -0.065, Γ: 0.000, ν: 517.340\n",
      "Price difference between Black-Scholes-Merton and Binomial Treee pricers for puts: 0.401021\n"
     ]
    }
   ],
   "source": [
    "bin_tree_pricer = options.BinomialTreePricer(maturity_date, impl_vol, strike, curve, cur_price,\n",
    "                                             ticker=TICKER, q=q)\n",
    "print(bin_tree_pricer)\n",
    "\n",
    "print('Price difference between Black-Scholes-Merton and Binomial Treee pricers for calls: %.6f'\n",
    "      % abs(pricer.get_price() - bin_tree_pricer.get_price()))\n",
    "bin_tree_put_pricer = options.BinomialTreePricer(maturity_date, impl_vol, strike, curve, cur_price, is_call=False,\n",
    "                                                 ticker=TICKER, q=q)\n",
    "print(bin_tree_put_pricer)\n",
    "\n",
    "print('Price difference between Black-Scholes-Merton and Binomial Treee pricers for puts: %.6f'\n",
    "      % abs(put_pricer.get_price() - bin_tree_put_pricer.get_price()))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "regular-headset",
   "metadata": {},
   "source": [
    "Now pricing with our volatility forecast:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "conservative-jacket",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1219.21, σ: 0.1456, Δ: 0.974, Γ: 0.000, ν: 222.560\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 4.37, σ: 0.1456, Δ: -0.017, Γ: 0.000, ν: 222.560\n"
     ]
    }
   ],
   "source": [
    "bin_tree_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, curve, cur_price,\n",
    "                                             ticker=TICKER, q=q)\n",
    "print(bin_tree_pricer)\n",
    "\n",
    "bin_tree_put_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, curve, cur_price, is_call=False,\n",
    "                                                 ticker=TICKER, q=q)\n",
    "print(bin_tree_put_pricer)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "exposed-silver",
   "metadata": {},
   "source": [
    "## Pricing with your own discount curve\n",
    "You might want to prefer to price using your own discount curve. The `YieldCurve` class of the `pricing.curves` package makes it fairly easy to do."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "sustainable-ambassador",
   "metadata": {},
   "source": [
    "The `parallel_shift` method of the `YieldCurve` class allows for pricing using a parallel shift to the yield curve. Let's price the same options but with a curve points shifted down by 200bps. As expected, this makes calls cheaper\n",
    "and puts more expensive."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "entertaining-tragedy",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1162.20, σ: 0.1456, Δ: 0.969, Γ: 0.000, ν: 279.509\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 5.89, σ: 0.1456, Δ: -0.022, Γ: 0.000, ν: 279.509\n"
     ]
    }
   ],
   "source": [
    "bin_tree_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, curve.parallel_shift(-200),\n",
    "                                             cur_price, ticker=TICKER, q=q)\n",
    "print(bin_tree_pricer)\n",
    "\n",
    "bin_tree_put_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, curve.parallel_shift(-200),\n",
    "                                                 cur_price, is_call=False, ticker=TICKER, q=q)\n",
    "print(bin_tree_put_pricer)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "excited-holmes",
   "metadata": {},
   "source": [
    "Let's now price with a flat discount curve, which I'll construct by reusing the offsets from earlier in this notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "micro-indicator",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "^GSPC EUROPEAN call option with strike 4500.0 and maturity 2025-03-21, price: 1129.71, σ: 0.1456, Δ: 0.966, Γ: 0.000, ν: 316.353\n",
      "^GSPC EUROPEAN put option with strike 4500.0 and maturity 2025-03-21, price: 6.94, σ: 0.1456, Δ: -0.025, Γ: 0.000, ν: 316.353\n"
     ]
    }
   ],
   "source": [
    "# Reusing 'relativedelta' offsets defined in Cell 3, which are for 0, 1 month, 3 months, 6 months, 1 year, 2 years,\n",
    "# 3 years, 5 years, 7 years, 10 years, 20 years, 30 years.\n",
    "\n",
    "# Defining a flat curve with a discount rate of 2% and yearly compounding frequency.\n",
    "custom_curve = curves.YieldCurve(today, offsets, [.02]*len(offsets), compounding_freq=1)\n",
    "                                 \n",
    "bin_tree_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, custom_curve,\n",
    "                                             cur_price, ticker=TICKER, q=q)\n",
    "print(bin_tree_pricer)\n",
    "\n",
    "bin_tree_put_pricer = options.BinomialTreePricer(maturity_date, vol_tracker, strike, custom_curve,\n",
    "                                                 cur_price, is_call=False, ticker=TICKER, q=q)\n",
    "print(bin_tree_put_pricer)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "greenhouse-retreat",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.9.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}