diff --git a/submissions/tingjun/Figure_1.png b/submissions/tingjun/Figure_1.png new file mode 100644 index 0000000..aced1b0 Binary files /dev/null and b/submissions/tingjun/Figure_1.png differ diff --git a/submissions/tingjun/Figure_2.png b/submissions/tingjun/Figure_2.png new file mode 100644 index 0000000..b613d5b Binary files /dev/null and b/submissions/tingjun/Figure_2.png differ diff --git a/submissions/tingjun/Figure_3.png b/submissions/tingjun/Figure_3.png new file mode 100644 index 0000000..d657484 Binary files /dev/null and b/submissions/tingjun/Figure_3.png differ diff --git a/submissions/tingjun/README.md b/submissions/tingjun/README.md new file mode 100644 index 0000000..591a503 --- /dev/null +++ b/submissions/tingjun/README.md @@ -0,0 +1,36 @@ +# Stochastic Volatility Model Submission + +## Model Overview + +The **Stochastic Volatility Model** describes the evolution of a security price $S$ and its instantaneous variance $V = \sigma^{2}$, where $\sigma$ is the volatility. The continuous-time stochastic processes are defined as: +- $dS = \mu S dt + \sigma S dw$ +- $dV = \phi V dt + \xi V dz$ + +In this specific implementation, we assume the **volatility of volatility** is zero ($\xi = 0$). The drift of variance, $\phi$, is driven by the absolute log-return of the asset: $\phi = b \times |\Delta \ln(S)|$, where $b$ is a free scaling parameter. + +## Parameters + +There are 3 free parameters in the model: +- Initial Annual Volatility: $\sigma = 0.15$ +- Annual Expected Return: $\mu = 0.02$ +- Volatility Scale Factor: $b = 0.0225$ + +## Simulation +The simulation consists of 500 samples, each spanning 10,000 days. Starting from day 3, the discrete-time updates are calculated as follows: + +- Variance Update: +$V_i = V_{i-1} \exp\left( b \cdot |\ln(S_{i-1} / S_{i-2})| \right)$ +- Price Update: +$S_i = S_{i-1} \exp\left( \mu - \sigma_i^{2}/2 + \sigma_i \epsilon \right)$ + +Note: +* $\sigma_i = \sqrt{V_i}$ +* $\epsilon\sim \mathcal{N}(0,1)$ +* Data generation logic can be found in the [price_generator.ipynb](price_generator.ipynb). + +## Results +The model achieves a high $R^{2}$ with the simulated paths. +* **Global R²**: **0.998** + + + \ No newline at end of file diff --git a/submissions/tingjun/dataset_part1.parquet b/submissions/tingjun/dataset_part1.parquet new file mode 100644 index 0000000..dbc1995 Binary files /dev/null and b/submissions/tingjun/dataset_part1.parquet differ diff --git a/submissions/tingjun/dataset_part2.parquet b/submissions/tingjun/dataset_part2.parquet new file mode 100644 index 0000000..4ae4e5d Binary files /dev/null and b/submissions/tingjun/dataset_part2.parquet differ diff --git a/submissions/tingjun/dataset_part3.parquet b/submissions/tingjun/dataset_part3.parquet new file mode 100644 index 0000000..7b2f3d2 Binary files /dev/null and b/submissions/tingjun/dataset_part3.parquet differ diff --git a/submissions/tingjun/price_generator.ipynb b/submissions/tingjun/price_generator.ipynb new file mode 100644 index 0000000..06fbef7 --- /dev/null +++ b/submissions/tingjun/price_generator.ipynb @@ -0,0 +1,351 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 182, + "id": "0f4d0758", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "End price after 5000000 days: 27.77\n" + ] + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAk0AAAHHCAYAAACiOWx7AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAXi9JREFUeJzt3QecE2X6wPFn6U1AQJp0UVGkCAiigiAIIhaU8287RUQ9PPBEThAsgKKHggWUoh4ioCLtsID0jrL0XqUKCEsTls4uS/6f59UJk2x2d7KbbNrv+/mEkMns5M1kyjPPWybO5XK5BAAAAOnKkf7bAAAAIGgCAABwiEwTAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAgAMETQAAAA4QNAEAADhA0IQsq1Spkjz11FMhWZN9+vSRuLg4CSe7d+82ZRo5cqREEqvc77//fkg+f/78+ebz9RmBoetT95Fwc/fdd8uzzz4b6mJEpWAej6dPny6FChWSw4cPS6wiaEKa1q9fL3/729+kYsWKki9fPrnyyivlzjvvlE8++SQq1trQoUNDHthooNK+fXu56qqrzDouXbq0NG7cWHr37h12ZXVygrYeOXLkkLJly0qLFi0iIgg6ffq09O3bV2rWrCkFChSQIkWKSKNGjWT06NESTnea0m3Avp7TeuiJM1z98ssvMnPmTHnllVdSBey+HmPHjk21jM2bN8tdd91lTuDFihWTJ554wueJ/OLFi9K/f3+pXLmy2b/09/3222+D/h2j1V133SVVq1aVfv36SazKFeoCIDwtXrxYmjZtKhUqVDBXhHoy37t3ryxZskQGDRokL7zwgnverVu3mpNkpNFApESJEiHLkm3fvl1uuukmyZ8/vzz99NPmRHfgwAFZtWqVvPfee/Lmm2+GTVmd0qD6ySefNIHGrl27TLnvuOMO+emnn6RVq1bp/q0Gi2fPnpU8efJIdjp48KA0a9bMnIgfeeQR6dy5s5w7d07+97//Sbt27WTq1KnyzTffSM6cOSXUdB199dVXHtOeeeYZqV+/vjz33HPuaRpMKF2fuXKF12F+wIABZn3rydfbo48+arJQdg0bNvR4vW/fPrMeNLD9z3/+I6dOnTLZUb3IW7Zsmcf289prr8m7775rjmG6r/3www/y2GOPmWBMf2v47x//+Ie8/PLL5vh02WWXxd4q1Bv2At7uvvtu1xVXXOE6duxYqvcOHjwYNiusd+/emgbI1N9Wr17ddfvttwe8TLt27TJl+vLLL9Od75///KcrV65crt27d2e4joNVVl/lHjBgQKb+Xv+2U6dOHtPWrVtnprdo0SLNvzt79qwrJSXFFSotW7Z05ciRw/XDDz+keu/ll1825X/33XeztUy6PnS9OFGwYEFXu3btXJFAt2vd5ocPH57pbe/555935c+f3/Xbb7+5p82aNcv8/Weffeaetm/fPlfu3Lk9tsmLFy+6GjVq5CpXrpzrwoULrmhUsWLFoG4PBw8edOXMmdP1xRdfuGJR5KUHkC127Ngh1atXl6JFi6Z6r2TJkunWoVtVCD///LP861//kiuuuMIsR69QkpKS5Pjx4yYbcfnll5tH9+7dPapA0mrb4rSt0JdffmmyG1rOvHnzyvXXXy/Dhg1LVeaNGzfKggUL3NUATZo0cb+vZezSpYuUL1/eLEOvijX7o+l+O51Pv7te9ep31MyETnO6jsuVK2eqP9NbxxmVdefOnfLQQw+ZagqtWrr55ptNZsebZk+0fcs111xjqirKlCkjDz74oClHWvR30QyGXr1PmjRJ/FWjRg2TIdOsk/231SqX119/3VT5aplPnDiR5u++dOlSk33QbaVgwYKmikWznXZbtmwxVcm6DvS71atXT3788ccMy6eZ0xkzZpjf8L777kv1vlZDXH311ea316xNcnKy+QytUvWm30E/W6/CLefPnzdVrbr96Hak25Nu7zrdTr+3Zrg0o6X7nc6r7UcC3abJagP466+/yt///nez3er++cYbb5jfWrPJ999/vxQuXNhklz/44INUy3T6nXzR7fLChQvSvHnzdKtK9TiRFs0A3nPPPSYLbtHl6XY9fvx49zTNKunv9c9//tNjfTz//PMmWxUfH59uWRMSEszvrPuofk/dX3Td6HHI/hmtW7c2VdE6j1azazVvSkqKx7J0f73hhhtk3bp1cvvtt5ttXtffxIkTzfu6bzdo0MBkna+99lqZPXu2x99bv5tu5//3f/9nfp/ixYvLiy++aPbrjDg9nul+WbduXZNB0s/Q/dd7XytZsqTZB/W7x6LwytsibOiJXA8qGzZsMDt7ZmgVnh54NY2rJ6fPP//cBBZa9acHPE2ta9WHpuv1MzSQCgQNkPTEoydBrZqYPHmyOXDqAaJTp05mnoEDB5ryaTWGpvBVqVKlzPOZM2fMge333383gZ6WVcvcs2dPU32mf6v0JKMHUQ0OO3bsKNddd5189913JnByuo714Dh37lwT5KUlvbJq1dItt9xiyqwBqh5IR40aZb67HpAfeOABM58exPVEM2fOHFMtoQfbkydPyqxZs8xvrAd7b/o3Wm04btw487305OCvY8eOmYd3VYyeWDQQ0wBDT7ZpVclp+bTcesLSMuv2pNVoU6ZMMa+VBpS33nqrCcB69OhhAis9ebZp08acYK114ItuGyqtbU+3H63O0W1Y2+LoyVmXpwHkZ5995lHu77//3nwXq9pHtzf9HXT70MBTtw+tQvroo49M0KLz2+l2oOXW4EkDzWC2S3r44YdNebTqSgOZt99+2wSD+p10W9QTqgZw+vtotZZWh2XmO3nT/Ui3UV8XCkrXc7du3UyAoCfvd955x7SLs+g+eejQIRMUe9MqSj2eWFavXm22BS2j93zW+7fddluaZW3btq3ZtnTf099CP1e3xz179rh/G72A0/2ya9eu5ll/w169epkAWo9rdrof6Las24de5OhxSv+v61kDGj2G6Lamf6cXABrAeld/acCkn63BvB5TP/74Y7NcbXuXFqfHM/1uWj2qVaf6+yvd13S7t/Y1i/42Gf3WUSvUqS6Ep5kzZ5oUrD4aNmzo6t69u2vGjBmupKSkDNPBWi2lm5ZWe2g63KLLiYuLc3Xs2NE9TVPkmiq3Vz3NmzfP/L0+Z1Tt5at67syZM6nKqGWpUqWKoyqvvn37miqPX3/91WN6jx49zPrYs2ePef3999+bz+7fv7/H99H0v5PquQ0bNphqBp23du3arhdffNEs8/Tp06nmTausXbp0MX+/aNEi97STJ0+6Kleu7KpUqZK72mvEiBFmvg8//DDVMqzfyF5Fkpyc7Hr44YdN+fR3d0L/tkOHDq7Dhw+7Dh065Fq6dKmrWbNmZvoHH3zg8dvqb+H9O3n/7rou9Xvo9uVdTWzfrvQzatSo4Tp37pzH+7fccovr6quvTrfMbdq0MZ/pqxraMmnSJDPPxx9/bF7r+tDXkydPTlWlbd/GvvrqK1PtZ/9t1Keffmr+/pdffvFYdzrvxo0bXf5Kr3pOl6v7iPf+8txzz6XaB3XftFdD6jrR39++bH++ky+33Xabq27duqmma1WbVuEOGzbM9eOPP7oGDhzoqlChgvmsKVOmuOdbvny5+ZzRo0enWka3bt3Me9Z20Lp161T7vNL9S+fT/Tkt+t2dVBf6Otb84x//cBUoUMBje9R9V5c3ZswY97QtW7a4f/clS5a4p1vbl6/j3H333Zeqil+nr127Ns3jsdPjmR5/Chcu7Kja8j//+Y/53HBqqpFdqJ5Dmg16NdOkV5Vr1641PVBatmxpruadVHuoDh06eAwHoOlnPY7rdIs2rtWrRq1iChRNcVsSExPlyJEj5kpLP0NfZ2TChAmm55RWB+nfWg/NMmj2ZeHChWY+varVTISm++3fx95IPj2aDVuzZo2pJtGUv6bBNTuiWaT//ve/jpahZdArZ/sVs17xahZAl7lp0yYzTTMumr3wVTbvIRu0akSvhDWbo8u3X+ln5IsvvjDVPZrC199br1L1KlyvpO00G2f/nXzRTIBW6+nfelcTW2X+448/zNW9XoFr5sz6rY4ePWq2123btpkr7LTo36j0GrRa72n2QGkmRtelZuAserWvV+qawbFvR5rlqFatmsd2ZGUV582b5/E5uo1qVXJ20Mbj3vug976p61yriuz7pr/fyZv+LrpfedPsh1aTarbl3nvvNZkN/f11W/r3v//tnk+rSJVWMXnTqlH7PPrsZD5fdNvULKJWFetvm958Fmv702OHZne0Ks1O90t743Ndt7qOdX3qvmKx/u/rmGhlyi3W/mzPsGX2eKZl0apR3Y4zcvlfv6EuJ9ZQPYc0aVpeqyH0JKqBk1bRaBpeU8d6ss/oAG9vc6C0/YTSenXv6ekdmPylJ2ptc6FBnx687DRossqRFj3RatsDPWD7oml69dtvv5lqI6unkv1g6JS2w9DeUHrw0gBHAxUNUDXo0W7S6bX9sMpgP+BarCoJfV+rPrXdkpbLSU8qTf1rj6Rp06Z5tJ1yQqsrtXpJgxoNNjQw1CoSb/rdMmK1tUqvelh7IOrJXtvk6COt30uD/fQCIj3h+Wq/Z71nn1fXoVbdjBkzxlTH6YlZ9xNtP2MPmnQ70uqNjLYjf9ZJoPjaNzWY0GDQe7oGOpn9Tr44HcLBajumVYjaBknbFllBiq/2U1bbHmsefXYyny/6m2oVlQZsehGj7QS1ak2rcbWK2KLVd9o2TwN3K6i2eF+gafm9L1B0/fo6Hipfx0RtX2en1erac9neziqzxzNtwqDVw9rLVfcXvVjSixEdZiCt3zAuzMbIyw4ETciQXnFpAKUPPcnrgUyvXrzHEvKWVhdtX9PtB9K0dkTvxpVpnWi1Tl6vhD/88ENzQNLy65WYBnzeDR990Xk006aNW33RdRBouk600aU+tIu1DvegbR0yCpqCQTM02ghZgzcNmqwrcyf0xOCkzBllmZyyfk9te6Pl9sVX13Z7cKltM/SkYrXb8abvKftFgmYMtP2PBpaaHdSTjW5ztWrV8iib/p66HfrifbIM1Dpxwtc+mNb+at83/f1O3rQ9kz8XSNbyNKOo25ZepChti+NNp2mgZWWXdF7NfGn57ccU62+18XZ6NMOpWS/dPjQLpkG5XlBogHTjjTeaxtWaHdQG02+99ZZ7rDUdMkTHoPI+1vhzPHQaXDoJWpwezzQ7rBfD+l11u9aHdqrRQFHbSdod++s39A6yYwFBE/xiNcD0ddAKFCv1690LTbMmGdGGvXp1qVWI9qtpX9UGaR1w9OCnmZaMTv7amFUbVuu89myTjlsV6HWcVlm1DL4+z6oasBrc6nfSXmiaDcmdO3e6n69X1VpNolfWWk2nGcZQjPVjNU7Xhupp/RZVqlQxz/qdMhNg6nfUE6E2pPUVNGmgrhkl3Sa1sblF59WTslbRadWonkitRvr28muGVoP4aLkiz+p30sBSq4qdsqqorCyJZkD0/ytWrEg1r47RVLt2bfdr/f/w4cNNZswe8Op+YL3v5PtqtkkfmrHRv9EehV9//bWputMsnGYZ7duO1VM0GLQM9oykZlo1KEqv04DT45nSC0wNFPWhy9Xsk14caMBov/jYtWuXCZjSyl5FM9o0wSfrCs2bVXfuTxWUv/REr1dfVl27RQdKzIh11WYvu6bJ9YrJm1Yb+RoeQFPSWrWnV1zedH7tMq20G7z+3z6cgZ5knY6YvmjRIhPEOFnHaZVVy6AnC3v3aW2XoD0V9UBqnSy0OknbHwwePDjVMnz9znqA1e7HmnHS0ZadZOgCrU6dOuYEob17vL+7VWa9OtZsmB7YfQXyGd3uQXse6nfV7UOrRr1pIKS9wvQq3Z4J0ioRrabWIF2rV3U7sFfNWduRtqfy1T5N29Po7xRpsvqdNIuqWQrv9jq+fif9nBEjRpju7VaGydqW9bfS3mUWvXjR30mDfHtVsQbT9uOGbjeffvqpCb70t0+LVut7d+XX4EOraK0qP1/HGm3K4OQ4lVlDhgzxeG0da9IbONbp8cxeDWtt47rulXc158qVK1MNOhoryDTBJ21gqAcO7V6tV4d6MNBuqnplrSdjX+PUBIrW6evBTw8IejWrBys9SDppL6H18NbVknav1SssPcDrydX7pKrdZjXg0e7WehWl82iDVu3yrJkqzULo+D06n54MtGu1duPX9gN6laWfodkH7eau0zRA0atOJ43NlbaZ0IOPjpVkHZw0ta9ZD61msDeeTqus+tl6Wwg9aOqQA/p3mkrXK0G9ordGatcUuy5XG2VrkKUNQ/U76ZAHejWpJxhvWu1kpee1CkIDk+ykZdfvrOtZr/B1m9OTp2bRtC2JdRLQE4lme7TaSEd+1uyTDsWgJwptC6OZkfToetHMia4D7fKt60ZPEvpbajZBgyHdJrzpdN1GtZpaP9u7a7sGm1ptp1k7vQjRbUWDai2/Ttfy++o6H86y+p102ArNWup2Zx/BXINSq2pdq810f9LtTbdR73GCXn31VdM8QKuwtcG47uPaTV9/A/txSavzdB/S9/TiRJsXaFWbXqxkNMK7BmBaFg04dL/WMmvGVbcrqzG3Bl2agdRODbrv6bFKA+hg3nZH92vtnKPtjHT71oyXbrP2amFvTo9n2jlAq0H1uKLrTjP7un3rvmfftg8dOmSqrL0bpceMUHffQ3iaNm2a6+mnn3ZVq1bNVahQIVeePHlcVatWdb3wwgupupmmNeSAdg+2s7rNapd0O/1b7RJrp/O0bdvWdN29/PLLTTde7aLvZMgB7bJcs2ZNV758+Uy3+/fee8/d5V671VsSEhJMt+TLLrvMvGfv0q/d9nv27Gm+s373EiVKmC7s77//vsewC0ePHnU98cQTpqtukSJFzP9Xr17taMgB7Z6toxXfcMMN5m919GLtZv3UU0+5duzY4TFvemXVef/2t7+5ihYtar5z/fr1Pbpp27tHv/baa6Ybv35W6dKlzd9Zn5XWqMxDhw4103V0bH9HBPdmDSswYcKENN/zHmri559/dt15553mu+t2or/tJ5984jGPfocnn3zSfCf9bldeeaXrnnvucU2cONHlhP7effr0MUM7aDd7/axbb73VNXLkSI/hDex0evny5U2Z3377bZ/z6Lai258uN2/evGZb1i73b775pisxMdGvdRfIIQec7INKtzMte2a+U1q027wOE2GnXfEbN25s7kKgI4br/vbAAw+4Vq5c6XMZeizQIQr0+KDb/eOPP272EW865IZ2j9djlO7HWuavv/46wzIeOXLE/B56/NP1ovtngwYNXOPHj0+1D998881mmylbtqx7aBbv7djXelRaLt2vvXlvD9bvtmnTJrPP6vap671z586pRo73NSK4k+OZ7iu6TkuWLGnm0WORHncPHDjgsaxhw4aZ9X7ixAlXLIrTf0IduAEAYoNmerRKVbNT3r3B4JuOCK4Df2o1ZqgbX994443m99OONbGINk0AgGyj1Z9aja69MxFZpk+fbhqj62jisYo2TQCAbKXd2RF57rrrLtOGLJaRaQIAAHCANk0AAAAOkGkCAABwgKAJAADAARqCB4iOmLx//34zYmy03DIBAIBo53K5zI25dWBVa0DgtBA0BYgGTBndrBIAAIQnvTWPjoaeHoKmANEMk7XS9ZYTAAAg/J04ccIkPazzeHoImgLEqpLTgImgCQCAyOKkaQ0NwQEAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAABwiaAAAAHCBoAgAAcICgCQAAwAGCJgAAAAcImmKYy+WSc8kpoS4GAAARgaAphv1r7Bqp9sZ02XXkdKiLAgBA2CNoimGT1+43z6MW7w51UQAACHsETQAAAA4QNAEAADhA0AQAAOAAQRMAAIADBE0AAAAOEDQBAAA4QNAEAADgAEETAACAAwRNAAAA4R40DRs2TGrWrCmFCxc2j4YNG8q0adPc7587d046deokxYsXl0KFCknbtm3l4MGDHsvYs2ePtG7dWgoUKCAlS5aUbt26yYULFzzmmT9/vtSpU0fy5s0rVatWlZEjR6Yqy5AhQ6RSpUqSL18+adCggSxbtiyI3xwAAESakAZN5cqVk3fffVdWrlwpK1askDvuuEPuv/9+2bhxo3n/pZdeksmTJ8uECRNkwYIFsn//fnnwwQfdf5+SkmICpqSkJFm8eLGMGjXKBES9evVyz7Nr1y4zT9OmTWXNmjXSpUsXeeaZZ2TGjBnuecaNGyddu3aV3r17y6pVq6RWrVrSsmVLOXToUDavEQAAELZcYebyyy93DR8+3HX8+HFX7ty5XRMmTHC/t3nzZpcWOT4+3ryeOnWqK0eOHK6EhAT3PMOGDXMVLlzYdf78efO6e/fururVq3t8xsMPP+xq2bKl+3X9+vVdnTp1cr9OSUlxlS1b1tWvXz/H5U5MTDRl0+dIUfGVKebR+4cNoS4KAAAh4c/5O2zaNGnWaOzYsXL69GlTTafZp+TkZGnevLl7nmrVqkmFChUkPj7evNbnGjVqSKlSpdzzaIboxIkT7myVzmNfhjWPtQzNUuln2efJkSOHeW3N48v58+fN59gfAAAgeoU8aFq/fr1pr6TtjTp27CjfffedXH/99ZKQkCB58uSRokWLesyvAZK+p/TZHjBZ71vvpTePBjlnz56VI0eOmIDN1zzWMnzp16+fFClSxP0oX758FtcEAAAIZyEPmq699lrT1mjp0qXy/PPPS7t27WTTpk0S7nr27CmJiYnux969e0NdJAAAEES5JMQ0m6Q92lTdunVl+fLlMmjQIHn44YdN1dnx48c9sk3ae6506dLm//rs3cvN6l1nn8e7x52+1t56+fPnl5w5c5qHr3msZfiimTF9AACA2BDyTJO3ixcvmvZCGkDlzp1b5syZ435v69atZogBbfOk9Fmr9+y93GbNmmUCIq3is+axL8Oax1qGBm36WfZ5tAz62poHAAAgV6iruFq1amUad588eVLGjBljxlTS4QC0nVCHDh3MUADFihUzgdALL7xgApmbb77Z/H2LFi1McPTEE09I//79TRuk119/3YztZGWBtJ3U4MGDpXv37vL000/L3LlzZfz48fLTTz+5y6GfodWC9erVk/r168vAgQNNg/T27duHbN0AAIDwEtKgSTNETz75pBw4cMAESTrQpQZMd955p3n/o48+Mj3ZdFBLzT5pr7ehQ4e6/16r1aZMmWLaQmkwVbBgQRP8vPXWW+55KleubAIkHfNJq/10bKjhw4ebZVm0KvDw4cNmfCcNvGrXri3Tp09P1TgcAADErjgddyDUhYgG2htPAz9tFK5ZsUhQqcef2banbqkkfe6rHuriAAAQ1ufvsGvTBAAAEI4ImgAAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAABwiaAAAAHCBoAgAAcICgCQAAwAGCJgAAAAcImgAAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAABwiaAAAAHCBoAgAAcICgCQAAwAGCJgAAAAcImgAAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAABwiaAAAAHCBoAgAAcICgCQAAwAGCJgAAAAcImgAAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAADCPWjq16+f3HTTTXLZZZdJyZIlpU2bNrJ161aPeZo0aSJxcXEej44dO3rMs2fPHmndurUUKFDALKdbt25y4cIFj3nmz58vderUkbx580rVqlVl5MiRqcozZMgQqVSpkuTLl08aNGggy5YtC9I3BwAAkSakQdOCBQukU6dOsmTJEpk1a5YkJydLixYt5PTp0x7zPfvss3LgwAH3o3///u73UlJSTMCUlJQkixcvllGjRpmAqFevXu55du3aZeZp2rSprFmzRrp06SLPPPOMzJgxwz3PuHHjpGvXrtK7d29ZtWqV1KpVS1q2bCmHDh3KprUBAADCWZzL5XJJmDh8+LDJFGkw1bhxY3emqXbt2jJw4ECffzNt2jS55557ZP/+/VKqVCkz7dNPP5VXXnnFLC9Pnjzm/z/99JNs2LDB/XePPPKIHD9+XKZPn25ea2ZJs16DBw82ry9evCjly5eXF154QXr06JFh2U+cOCFFihSRxMREKVy4sESCSj1+Ms9P3VJJ+txXPdTFAQAg2/lz/g6rNk1aYFWsWDGP6d98842UKFFCbrjhBunZs6ecOXPG/V58fLzUqFHDHTApzRDpSti4caN7nubNm3ssU+fR6UqzVCtXrvSYJ0eOHOa1NQ8AAIhtuSRMaGZHq81uvfVWExxZHnvsMalYsaKULVtW1q1bZ7JG2u5p0qRJ5v2EhASPgElZr/W99ObRwOrs2bNy7NgxU83na54tW7b4LO/58+fNw6LLAgAA0StsgiZt26TVZz///LPH9Oeee879f80olSlTRpo1ayY7duyQq666SkLZiP3NN98M2ecDAIDsFRbVc507d5YpU6bIvHnzpFy5cunOq22P1Pbt281z6dKl5eDBgx7zWK/1vfTm0brL/Pnzm6q/nDlz+pzHWoY3rSbU6kTrsXfvXr+/NwAAiBwhDZq0DboGTN99953MnTtXKleunOHfaO83pRkn1bBhQ1m/fr1HLzftiacB0fXXX++eZ86cOR7L0Xl0utLG4nXr1vWYR6sL9bU1jzcdukA/w/4AAADRK1eoq+TGjBkjP/zwgxmryWqDpK3YNQOkVXD6/t133y3Fixc3bZpeeukl07OuZs2aZl4dokCDoyeeeMIMRaDLeP31182yNbBROq6T9orr3r27PP300yZAGz9+vOlRZ9HhBtq1ayf16tWT+vXrm956OvRB+/btQ7R2AABAWHGFkH68r8eXX35p3t+zZ4+rcePGrmLFirny5s3rqlq1qqtbt26uxMREj+Xs3r3b1apVK1f+/PldJUqUcP373/92JScne8wzb948V+3atV158uRxValSxf0Zdp988omrQoUKZp769eu7lixZ4vi7aJm07N5lC2cVX5liHr1/2BDqogAAEBL+nL/DapymSMY4TQAARJ6IHacJAAAgXBE0AQAAOEDQBAAA4ABBEwAAgAMETQAAAA4QNAEAADhA0AQAAOAAQRMAAIADBE0AAAAOEDQBAAA4QNAEAADgAEETAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAAEETAADIqo37E2XGxoSYX5G5Yn4NAACAdLX++Gfz/GPnW6VmuaIxu7aongMAAI5sP3QqptcUQRMAAIADBE0AAAAOEDQBAAA4QNAEAADgAEETAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAgAMETQAAAA4QNAEAADhA0AQAAOAAQRMAAIADBE0AAAAOEDQBAAA4QNAEAADgAEETAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAgAMETQAAAA4QNAEAADhA0AQAABDuQVO/fv3kpptukssuu0xKliwpbdq0ka1bt3rMc+7cOenUqZMUL15cChUqJG3btpWDBw96zLNnzx5p3bq1FChQwCynW7ducuHCBY955s+fL3Xq1JG8efNK1apVZeTIkanKM2TIEKlUqZLky5dPGjRoIMuWLQvSNwcAAJEmpEHTggULTEC0ZMkSmTVrliQnJ0uLFi3k9OnT7nleeuklmTx5skyYMMHMv3//fnnwwQfd76ekpJiAKSkpSRYvXiyjRo0yAVGvXr3c8+zatcvM07RpU1mzZo106dJFnnnmGZkxY4Z7nnHjxknXrl2ld+/esmrVKqlVq5a0bNlSDh06lI1rBAAAhC1XGDl06JBLi7RgwQLz+vjx467cuXO7JkyY4J5n8+bNZp74+HjzeurUqa4cOXK4EhIS3PMMGzbMVbhwYdf58+fN6+7du7uqV6/u8VkPP/ywq2XLlu7X9evXd3Xq1Mn9OiUlxVW2bFlXv379HJU9MTHRlEufI0XFV6aYR+8fNoS6KACACDhf/G/lXle08ef8HVZtmhITE81zsWLFzPPKlStN9ql58+bueapVqyYVKlSQ+Ph481qfa9SoIaVKlXLPoxmiEydOyMaNG93z2JdhzWMtQ7NU+ln2eXLkyGFeW/N4O3/+vPkM+wMAAESvsAmaLl68aKrNbr31VrnhhhvMtISEBMmTJ48ULVrUY14NkPQ9ax57wGS9b72X3jwa6Jw9e1aOHDliqvl8zWMtw1d7rCJFirgf5cuXz/I6AAAgnLk0JxPDwiZo0rZNGzZskLFjx0ok6Nmzp8mMWY+9e/eGukgAACCIckkY6Ny5s0yZMkUWLlwo5cqVc08vXbq0qTo7fvy4R7ZJe8/pe9Y83r3crN519nm8e9zp68KFC0v+/PklZ86c5uFrHmsZ3rQXnj4AAEBsCGmmyeVymYDpu+++k7lz50rlypU93q9bt67kzp1b5syZ456mQxLoEAMNGzY0r/V5/fr1Hr3ctCeeBkTXX3+9ex77Mqx5rGVoFaB+ln0erS7U19Y8AAAgtuUKdZXcmDFj5IcffjBjNVnth7SNkGaA9LlDhw5mKABtHK6B0AsvvGACmZtvvtnMq0MUaHD0xBNPSP/+/c0yXn/9dbNsKxPUsWNHGTx4sHTv3l2efvppE6CNHz9efvrpJ3dZ9DPatWsn9erVk/r168vAgQPN0Aft27cP0doBAADhJKRB07Bhw8xzkyZNPKZ/+eWX8tRTT5n/f/TRR6Ynmw5qqT3WtNfb0KFD3fNqtZpW7T3//PMmmCpYsKAJft566y33PJrB0gBJx3waNGiQqQIcPny4WZbl4YcflsOHD5vxnTTwql27tkyfPj1V43AAABCb4nTcgVAXIhpoTzzNjGmjcM2IRYJKPf7MtD11SyXpc1/1UBcHABDm54sPHqolbeteansca+fvsOk9BwAAEM4ImgAAABwgaAIAAAh20HTu3Lms/DkAAED0Bk06flHfvn3lyiuvlEKFCsnOnTvN9DfeeEO++OKLYJQRAAAg8oKmt99+W0aOHGnGRNJBIS16vzjtxg8AABCN/A6aRo8eLZ9//rk8/vjjZowkS61atWTLli2BLh8AAEBkBk2///67VK1a1We1XXJycqDKBQAAENlBk96yZNGiRammT5w4UW688cZAlQsAAISZuDiJaX7fRkVvM6K3KdGMk2aXJk2aZG6iq9V2ejsTAACAaOR3pun++++XyZMny+zZs8193jSI2rx5s5l25513BqeUAAAAkXjD3kaNGsmsWbMCXxoAAIBoyTQtX75cli5dmmq6TluxYkWgygUAABDZQVOnTp1k7969qaZrGyd9DwAAIBr5HTRt2rRJ6tSpk2q69pzT9wAAAKKR30FT3rx55eDBg6mmHzhwQHLlylQTKQAAgOgLmlq0aCE9e/aUxMRE97Tjx4/Lq6++Su85AAAQtfxODb3//vvSuHFjqVixonswyzVr1kipUqXkq6++CkYZAQAAIi9ouvLKK2XdunXyzTffyNq1ayV//vzSvn17efTRRyV37tzBKSUAAECIZaoRkg5q+dxzzwW+NAAAAJEcNP3444/SqlUrk0nS/6fnvvvuC1TZAAAAIitoatOmjSQkJEjJkiXN/9MSFxcnKSkpgSwfYpTL5ZI+P26UcpcXkGcbVwl1cQAAcBY06Y15ff0fCJZ1+xJlVPxv5v8ETQCAiBtyIDk5WZo1aybbtm0LXokAETmTRMYSABDBQZO2adKecwAAIPa4XBLT/B7c8u9//7t88cUXwSkNAABAtAw5cOHCBRkxYoTMnj1b6tata4YfsPvwww8DWT4AAIDIDJo2bNjgvmHvr7/+mqr3HAAAQDTyO2iaN29ecEoCAAAQLUHTuHHjzOCWSUlJphddx44dg1cyAACASAyahg0bJp06dZKrr77a3G9u0qRJsmPHDhkwYEBwSwgAABBJvecGDx4svXv3lq1bt8qaNWtk1KhRMnTo0OCWDgAAINKCpp07d0q7du3crx977DHTk+7AgQPBKhsAAEDkBU3nz5/3GF4gR44ckidPHjl79mywyoYYRkdMAEBENwR/4403pECBAu7X2iD8nXfekSJFirinMU4TAACI6aCpcePGpj2T3S233GKq7SyM0wQAACTWg6b58+cHtyQAAADRdO85AACAWETQhLDEDXkAAOGGoAnicrlYCwAAZICgCQAAIBhBU3JycprvHTlyxN/FAQAARGfQ9Mgjj/iszjl48KA0adIkUOUCAACI7KBpz5498swzz3hMS0hIMAFTtWrV/FrWwoUL5d5775WyZcuaMZ6+//57j/efeuopM93+uOuuuzzm+eOPP+Txxx+XwoULS9GiRaVDhw5y6tQpj3nWrVsnjRo1knz58kn58uWlf//+qcoyYcIEU36dp0aNGjJ16lS/vgsCizG/AAARHzRpMLF48WLp2rWreb1//365/fbbTaAxfvx4v5Z1+vRpqVWrlgwZMiTNeTRI0vvbWY9vv/3W430NmDZu3CizZs2SKVOmmEDsueeec79/4sQJadGihVSsWFFWrlwpAwYMkD59+sjnn3/unke/z6OPPmoCrtWrV0ubNm3MY8OGDX59HwAAEL38uo2KuuKKK2TmzJly2223mdcaqNSpU0e++eYbcz86f7Rq1co80pM3b14pXbq0z/c2b94s06dPl+XLl0u9evXMtE8++UTuvvtuef/9900GS8ult3sZMWKEuVde9erVZc2aNeZ2L1ZwNWjQIBOcdevWzbzu27evCcIGDx4sn376qV/fCQAARKdM9Z7TKi4NKjQgqV+/vsn+5MyZM/Cl+2sk8pIlS8q1114rzz//vBw9etT9Xnx8vKmSswIm1bx5cxO8LV261D2P3gJGAyZLy5YtzS1hjh075p5H/85O59Hp6d3AWLNY9gcAAIjxTNPll1/us43JmTNnZPLkyVK8eHGPNkaBotmfBx98UCpXriw7duyQV1991WSmNJjRIE3bUmlAZZcrVy4pVqyYeU/ps/69XalSpdzv6XfTZ2uafR5rGb7069dP3nzzTYkGtB8CACBAQdPAgQMlFLSnnkXbTNWsWVOuuuoqk31q1qyZhFLPnj3d7bqUZpo0AwcAAGI4aGrXrp2EgypVqkiJEiVk+/btJmjStk6HDh3ymOfChQsm22W1g9JnHQ7Bznqd0TxptaWy2lrpAwAAxIZM9Z6bMWNGqunaOHzatGkSTPv27TNtmsqUKWNeN2zYUI4fP256xVnmzp0rFy9elAYNGrjn0R519kE5tT2WtpHSqjlrnjlz5nh8ls6j0wEAADIVNPXo0UNSUlJSTddARd/zh46npD3Z9KF27dpl/q9jQel72pttyZIlsnv3bhPU3H///VK1alXTSFtdd911pt3Ts88+K8uWLZNffvlFOnfubKr1tOeceuyxx0wjcB1OQIcmGDdunOktZ69ae/HFF00vvA8++EC2bNlihiRYsWKFWRYAAECmgqZt27bJ9ddfn2q6Dgyp1Wb+0MDkxhtvNA+lgYz+v1evXqahtw5Ked9998k111xjgp66devKokWLPKrFtAeffrZW1+lQAzoUgn0MpiJFipgsmAZk+vf//ve/zfLtYzndcsstMmbMGPN3Om7UxIkTzUCbN9xwA1sJAAB/iUvdJyym+D1OkwYhO3fulEqVKnlM14CpYMGCfi1LRxH3dUsWi69qQG/aU04DnvRoA3INttLz0EMPmQfCQ6zvmACAKMg0aRVZly5dzBAA9oBJMziaFQICgZgJAMKPK+08R0zwO2jS+7ZpRkmrxHT8I31o2yIdq0lH4QYAAIhGmaqe03u1ae+ytWvXSv78+U31l466DQAAEK38DpqsEaT1Jrj6AAAAiAWZuvfcggUL5N577zXd//WhbZkyamgNAAAQU0HT119/bW5uW6BAAfnXv/5lHlpFp13+M+rFBgAAEDPVc++8845pDP7SSy+5p2ng9OGHH0rfvn3NYJIAAAAS65kmHaNJq+a8aRWdDiAJAAAQjfwOmsqXL5/qPm1q9uzZ5j0AAIBo5Hf1nA5iqdVxeo84vf2I0nu+jRw50tzTDQAAIBr5HTQ9//zzUrp0aXNz2/Hjx5tpOril3ghXRwsHAoHbqAAAomKcpgceeMA8AAAAYoXfbZqqVKkiR48eTTX9+PHj5j0AAIBo5HfQtHv3bklJSUk1/fz58/L7778HqlwAAACRWT33448/uv8/Y8YMcw86iwZR2qOuUqVKgS8hAABAJAVNbdq0cd93rl27dh7v5c6d2wRM2jgcCIw4ViQAIDKDposXL5rnypUry/Lly6VEiRLBLBcAAEBk955j1G8AABCLHDcEj4+PlylTpnhMGz16tMk8lSxZUp577jnTGBwAACCmg6a33npLNm7c6H69fv166dChgzRv3lx69OghkydPln79+gWrnAgil8vF+gUAIFBBk942pVmzZu7XY8eOlQYNGsh///tf6dq1q3z88cfuEcKBrGJEcABAxAZNx44dk1KlSrlfL1iwQFq1auV+fdNNN8nevXsDX0IAAIBICpo0YLIagSclJcmqVavk5ptvdr9/8uRJM/QAEGjJKX/23ARCXY198SJV2UAscxw03X333abt0qJFi6Rnz55SoEABadSokfv9devWyVVXXRWsciKGtRq0KNRFQIzTgOnx4Uul9Sc/SwqBExCzHA850LdvX3nwwQfl9ttvl0KFCsmoUaMkT5487vdHjBghLVq0CFY5EcO2HzoV6iIAsnjHn/fc3HXklFQteRlrBIhBjoMmHcxy4cKFkpiYaIKmnDlzerw/YcIEMx0AACAa+T24pf2ec3bFihULRHkAAAAiu00TAABALCNoAgAAcICgCWEpLtQFANLAAPpA7CJoQliKY0hwAECYIWgCAACOuGJ8PRE0AQAAOEDQBAAA4ABBEwD4IdarJ4BYRtAEAADgAEET6KkGAIADBE0AAAAOEDQBAAA4QNAEAH5gRPDotWjbYdl37Eyoi4EwlivUBQAAINQWbz8iT3yxzPx/97utQ10chCkyTQhL3HsOQHZavvsYKxzhHTQtXLhQ7r33XilbtqzpwfX99997vO9yuaRXr15SpkwZyZ8/vzRv3ly2bdvmMc8ff/whjz/+uBQuXFiKFi0qHTp0kFOnTnnMs27dOmnUqJHky5dPypcvL/37909VlgkTJki1atXMPDVq1JCpU6cG6VsDAIBIFNKg6fTp01KrVi0ZMmSIz/c1uPn444/l008/laVLl0rBggWlZcuWcu7cOfc8GjBt3LhRZs2aJVOmTDGB2HPPPed+/8SJE9KiRQupWLGirFy5UgYMGCB9+vSRzz//3D3P4sWL5dFHHzUB1+rVq6VNmzbmsWHDhiCvAQAAIkecxLaQtmlq1aqVefiiWaaBAwfK66+/Lvfff7+ZNnr0aClVqpTJSD3yyCOyefNmmT59uixfvlzq1atn5vnkk0/k7rvvlvfff99ksL755htJSkqSESNGSJ48eaR69eqyZs0a+fDDD93B1aBBg+Suu+6Sbt26mdd9+/Y1QdjgwYNNwAYA7mMTY4IDMSts2zTt2rVLEhISTJWcpUiRItKgQQOJj483r/VZq+SsgEnp/Dly5DCZKWuexo0bm4DJotmqrVu3yrFjx9zz2D/Hmsf6HACxjR5zAMK695wGTEozS3b62npPn0uWLOnxfq5cuaRYsWIe81SuXDnVMqz3Lr/8cvOc3uf4cv78efOwVwMicOJiPQcMIFtxzEFEZ5rCXb9+/Uzmy3poA3MAQGTiOg0RHTSVLl3aPB88eNBjur623tPnQ4cOebx/4cIF06POPo+vZdg/I615rPd96dmzpyQmJrofe/fuzcK3BQAA4S5sgyatUtOgZc6cOR5VYNpWqWHDhua1Ph8/ftz0irPMnTtXLl68aNo+WfNoj7rk5GT3PNrI+9prrzVVc9Y89s+x5rE+x5e8efOaYQ7sDwDRj/ZNQOwKadCk4ylpTzZ9WI2/9f979uwx4zZ16dJF3n77bfnxxx9l/fr18uSTT5oecTocgLruuutMr7dnn31Wli1bJr/88ot07tzZ9KzT+dRjjz1mGoHrcAI6NMG4ceNMb7muXbu6y/Hiiy+aXngffPCBbNmyxQxJsGLFCrMsAACAkDcE18CkadOm7tdWINOuXTsZOXKkdO/e3YzlpEMDaEbptttuM8GNDkBp0SEFNLhp1qyZ6TXXtm1bM7aTRdsbzZw5Uzp16iR169aVEiVKmAEz7WM53XLLLTJmzBgzvMGrr74qV199tRnW4IYbbsi2dQFPcbQwAJCNaAiOsA+amjRpYsZjSotmm9566y3zSIv2lNOAJz01a9aURYsWpTvPQw89ZB6xKL3fAABigZ5vgIht04Tsk5RC0AQAQEYImiDfLtvDWgAcIjELxC6CJgBAzDty6tJgxUBaCJoAADHvy192x/w6QMYImgDAD7QXBmIXQRMAAIADBE0AAAAOEDQhLFEFgnBF7zkgdhE0AQAAR1wxvp4ImgAAABwgaAIAAHCAoAkAEDben7FV2o1YJhdSLoa6KEAqBE0A4AdXzLfqCK7B87bLgl8Py/yth4P8SYD/CJoAAGEnmUwTwhBBEwBkINZ7DAH4E0ETACDsEKgiHBE0AQAAOEDQBAB+YERwIHYRNCEscRsVAEC4IWgCAIQdMnoIRwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAgAMETQhLcRIX6iIAAOCBoAkAAMABgiYAAAAHCJoAwA+MH5Q9XNx9DmGIoAkAAMABgiaEpUi8jcqoxbul9w8bxEUqAgCiUq5QFwDw5fyFixG3Ynr/uNE8312jjDSoUjzUxQEiGtceCEdkmhCWZmxMkEh16vyFUBcBABAEBE0ISxcvukJdBMAnGigDsYugCWGJkAkAEG4ImoAAoy1G9Llw8VIbuz9OJ4W0LLEinC+ctLPHjsOn5EJK5LW9zKo4iW0ETQhL9EBDOLmQcukUTps1jF2+V5p9sEBe+HY1KyPGEDQBgB9ZD+6LmD3COaMxdP528zxtQ+R2WIml4WACiaAJYSnaq7jW7TsuR06dD3UxkAmxftKASFIEDomCwGCcJiDAMor31u49LvcP+cX8f/e7rVn/EYDq4hCscwlfB0/E7gVPXIxfNJBpArLZ4h1HWecRxj4CRoyfM4CYRtAEABmxB01ETUDMImhCWArn1HxWMThipCNqyg5UiYanuBjf/gmagADjYA8A0Smsg6Y+ffpIXFycx6NatWru98+dOyedOnWS4sWLS6FChaRt27Zy8OBBj2Xs2bNHWrduLQUKFJCSJUtKt27d5MIFz3uDzZ8/X+rUqSN58+aVqlWrysiRI7PtOyI2e88hcrODVM8hlsXFdqIpvIMmVb16dTlw4ID78fPPP7vfe+mll2Ty5MkyYcIEWbBggezfv18efPBB9/spKSkmYEpKSpLFixfLqFGjTEDUq1cv9zy7du0y8zRt2lTWrFkjXbp0kWeeeUZmzJiR7d8VsVGFFevp7UhEQ3AAETHkQK5cuaR06dKppicmJsoXX3whY8aMkTvuuMNM+/LLL+W6666TJUuWyM033ywzZ86UTZs2yezZs6VUqVJSu3Zt6du3r7zyyismi5UnTx759NNPpXLlyvLBBx+YZejfa2D20UcfScuWLbP9+yL6M03RHBBGK6pcAUREpmnbtm1StmxZqVKlijz++OOmuk2tXLlSkpOTpXnz5u55tequQoUKEh8fb17rc40aNUzAZNFA6MSJE7Jx40b3PPZlWPNYy0jL+fPnzXLsD0AREkU3bSYAIDaFddDUoEEDU502ffp0GTZsmKlKa9SokZw8eVISEhJMpqho0aIef6MBkr6n9NkeMFnvW++lN48GQWfPnk2zbP369ZMiRYq4H+XLlw/Y9wYQzrdRASLHxYuugN5YOC7GLxrCunquVatW7v/XrFnTBFEVK1aU8ePHS/78+UNatp49e0rXrl3drzXIInCCE8kXyEVFmou2+uIYP2cgwjwwbLEcPnFOFnRvKrlzhnWeJCJE1BrUrNI111wj27dvN+2ctIH38ePHPebR3nNWGyh99u5NZ73OaJ7ChQunG5hpTzudx/4AnJi8bj8rKtIQ5yJC6W2b9ieek18PngzI8uIktkVU0HTq1CnZsWOHlClTRurWrSu5c+eWOXPmuN/funWrafPUsGFD81qf169fL4cOHXLPM2vWLBPgXH/99e557Muw5rGWEQvC8co5khveZlT0s0kp2VUUxMj+AiB7hHXQ9PLLL5uhBHbv3m2GDHjggQckZ86c8uijj5p2RB06dDBVZPPmzTMNw9u3b2+CHe05p1q0aGGCoyeeeELWrl1rhhF4/fXXzdhOmilSHTt2lJ07d0r37t1ly5YtMnToUFP9p8MZIHQiN2RC9LdpImpC5AnUdhsX45t/WLdp2rdvnwmQjh49KldccYXcdtttZjgB/b/SYQFy5MhhBrXU3mza602DHosGWFOmTJHnn3/eBFMFCxaUdu3ayVtvveWeR4cb+Omnn0yQNGjQIClXrpwMHz6c4QYQNLF+0IlEtGlCpOO4EwNB09ixY9N9P1++fDJkyBDzSIs2HJ86dWq6y2nSpImsXr1aYhXncCB9EVxbHLFY54F1NjlF+k3dLK1qlJHa5T17nfsjLsbPGGFdPQdEJs6wAMLLoNnb5LOFO6XNkF9CXZSIRtCEsBx3g6tMhG2bpjDcX4CMbAtU77m42F7XBE0IS9F8q5EcsX7UiUD23pz8etm0zqP4GJBdIrkXcrgiaEJYYl+Hv34/flZenrBWNh84EdTtkZgXkcK+3R4+dT4gy4yT2EbQhLDcCc4lB27Y/+xGwBcaHb9aKRNX7pPWHy8KUQkQSJsPBKY6KZbZ80zJKWSdAoGgCWFpzhbPUdqjCZmK4NiS8GeG6SLnhqiQNxenp3CsnosLx6vsbMRWibDcCcjWIFy3R9qkZY+rriiUTZ8Uvbh+CDyCJoTluBvR3IAx/NY2/BrcktWFCBGcw2icxDKCJoTlPhDJIVMklx2+8ZsiEgWjB2JcGJ4vshNBUxRIueiS1XuOSdKFzDWejvF9INsxzk+EY4dBhIjihH3IEDRFgUFztskDQxfLS+PWhLooQFSeIMK5unj6hgOybt/xUBcD2XyhHCpxEtsImqLA5wt3mOef1h+QqBG+56gMZXR+jfWDTsSPCB5Gv+DG/YnS8etVct9gbo0RKz5bsEOq954uG35PzHDeMI71IxZBUxTI6kH8fCar9eAbIxlH+Yjg4RMzya4jpyVahdN6Dif9pm0x49i99t36ELVpipNYRtAEAH4I1SnjQOJZc9PVQyfPuaclJF76P2KLk3CITFPg5QrCMpHNojHwJ6sc+QNNlimcX4oUyC3RwPM2KqHZ4Z74YplsP3RKFvx6SCb981YzbcKKfRIt5m89JJ8v3BnqYkQMJwFRMI6jcRLbCJqAbFa7QlHZGcXVKmv2Hpc2Q36RAnlyyqa37oqKQNujTVOIzhoaMKlVe45HZVXwU18uD3URorpa+bJ8ueTkuQshLU80oHouCkRj5B8XxVeAN5QtItFs7pZD5vlMUkq2fm4we7jZB7cEwoGTgNk+BwFTYBA0RaFlu/6Q179fLyfPJUs4+nnbEdlx+M+rZkSfUHXPd8VgQM999mKXo+q5CIz1v1u9T7qOW5PpcQeDjaApCv3fZ/Hy9ZI98sHMXyXcaBfpv3+xVJp9sCDd+UK1r+/944w88cVSWbTtcIhKEPlClZUJ9MdOW39A/j58qRw+ed6rTZOEjXAePwrB9cfppIxnCsLmERfk7f+lcWtl0urfZdyKvRKOCJqi2O6j4dduZsuBk2F9Mug6fo0s2nbENLpF5pxNCs8rRH89/80q+Xn7EfnP1M1he8W+43D47ePR4vvVv8udHy6QnWGaFT/goOdkJLd5O+4kKAwBgibARrMKWRW5h6nAmLAyPK8QMyvxbHJEn3yQOV3GrZFth05Jj0kZj4cUrgIV7IdinLK4MMro2hE0RYFIGmzMaVFdUbwunX7E4h1HpPXHi2Tt3si6RUa0NTjNEReXpZOP3uLkyRHLzDAMiDznk1NC99kXUmT4op2y7aCzDH2wjqPhmmkNBYKmKOZ0Qy9VOK9kl3CP78KpfI/9d6ls3H9CHv3vklAXJablSLVN+LeR3D/kF1n462HzeyLyhPKidPiiXfL2T5vlzo8WZurvg9HMIS6bukKEazKAoCmKhePFQTjdtyutrEK4ye6u+3YXY6B7VnLKRfl0wY40b3qb1W3COm85argbwdXac7ccjMrtJXXQnH1W/XYsS38fsExTgJYTDQiaokBcVuuoszGQOXYmKazSwboeztnS76E8QIabD2dulVpvzpTdYToQp/52Wb2S1r8f8fMueXfaljRvepsjR+C2x2NRGjg1+2C+PD1yhUxcuS/qMg2ZDZo1gMzqbW6yutkF5TgaJzGNoAnZWiXVb+qWsFrjWmVS7Y3pcvTUnw3AA3GhnNGJPFKOOR/P3S4nz1+Q92dulXCjJ6S2wxabtkJZCZw6fr3S3AA1o5O3fRiFrOwvN/adJav2ZC17EEybD5wwDd/9deKvdmyzNx+UaJPZoOnliWvl5n5zZPLa/SEbviNQHRhC0Zs5LkwPlARNEWL/8bPy4axfPW7WmZUNPVQN+5JSwqs7evzOo+Z5xsaDHreqQHik5k+fv2CyQOv3JXpM//34WXM7ER0e4nwWBsGzfveM9r1AGrV4t4Sjlb/9Ia0GLZJb353r832tvtS2WemJwtq5TJ+8J6363Tz/e/zaTH92Vo/TwTjOx9m2l+dGrzBj28VSUw6Cpgihg+x9PGebdPxqZcCXHQ6b5nGvartAXSHpSdd72b5k9fPsVXzRJjmbR+bVkYB/XLvftJPRAVq1vdG9g39O80QW7ME0V+85HtDAMWeYXkJbt785dd5370etvtTM3r5j6Z0k019Tm/afkCN/ZXUjRVbbtGXlQjHLmaZADTngY1rbYfEyc9NB6TxmVWA+JNVnujwSAJqhTWvbzE4ETRHCusGr/WadgdoJrg/xvdC0S23tt2bJl7/sSnOeA4mZu9qv3nuGWbYGT8E8uERjtYRFD4xpmb7hgMlMBLLK6ZO52+Rf3642N/3V6iJfDtnG08qO7EYgqydyhmnDOafBwf7jaWe70/sttiaclLs/XiT13p4tkWRPkDIp3s4mpZhbTGnHhEAJxvhicV7bSbDWj70RvGbtHhy6WB4c6rvdYXYiaIoGcVk7CTS7rqSEknapVW9O3pTmPA37zZX5W/+8Es6M346eCWqmKBx73WWHjl+vMlVlHUYG7g71MzYmmGddblpX2vbsYUo6Z+rlu/+Q1QEI6OwfERfkoEnba+lo1Nk9EnUgGl6nlxlZtvsPiUS6HWYHzdjoLabes7WxC5fqOY/bCImn9Pa/rDiXfCl4/H7Nn1Wdvx4MffMJgqYo5nSHcXKovJBy0aTm+08PXUPuEb9kvi1IRueDeVkIyLyDpsw0pA32AHmBvHr1JZA317Qfg9PahBMSbZmmNA7aesPqhz6NlweGLg769w8krZrU0ajvyOD+jIHmNAGWXtbNn/On9izTNlRjlu5x/kdRbM5f1aOj438LYEPw4HMF8EMOnbiUxeSGvch26TUa93c7X7jtsGkEOnT+DgmUxDP+BRc5s3AhnFHQlNVMkf2Ec+hE+m02zmRj+ycNFur2nS1XvzZNfgvivQhPB3AsKa2myOgEPWnVvgxPLPbgNesH4MydGXxli/5nK7tmOD9fuMOjA0Koetc57R6f3rhh6QVU3utCL8C0+vXV73zfpiQ2c7ee23PWM00B6j1n2/7jvH6YlABGTXqvR0vTaqGtAUkLmaYo5jSV6SRecHJzyEDytR9mpS1IRkFRVqsm7MXNlUF013961rvwO02JHzh+zt148vYB8zP9edmZqbFXh6R1PLafWNI6aAdynJ/Mnhd8ZYuSUy4tbODsbfKfqVuk+YcLshzAa5VeVoxd7uyegQXz5krzvfQyI7ot2p11ePEQqpt3B4N25mkyYJ57iJOM7u0WrIbg2o7Un57Y6bkYwN9nje2WUYXypb2dhRJBU8wObunf/K99t0ECLc7Pre/3dBqg+mI/EGUUb2X19GqvIsqdM7i7Ve8fNkjdt2c5OujtTbenk3NNshBwZYXLx2+qAeMu24CbTgLILA8SGKRgTLtte8vs4rVKL9QuphNbT/+rrZpl2gbP174s3XlU6r4921RZhpoGOlkd6FWHjdl99Ix8vmine5pWT9obx9uPza4gtiN98otlmWzTFOf4N/dXrXJF0wgkwyfvSNAUozzSrSFKhKfX9dpXUJdWTyo7HUhu7LI9qRvvZrDTZbVDkz2tnCvIvaNGxf8mx88ky+jFl9o+2I1bvsccnAPZLiCzjWG1jYK/mYJrShVKMyB6fPhSuX3APLlgm77t4CmfY+HYf4Wsjxzue7lZ5Sves2eapm9IkGtemyYHbW09QsF+UZDerhTI3lp6o+pnRq0wt5/R3pShpsFbk/fnZ1iNqUNl6CM99mOuVk8e9TFS/LT1B2TZrqw1nk9vs9+SkLmbAHsLZPWcR5YtTAf9ImiC/LLj0gk/LXlzZbypaNsMHVPHqZmbEtK8rURmDr66k73w7WrpMWm9OVnb08YZVXlktRu4fSgIe3f4YErrK73yv/WmGmBLwomAVqtlZoyU+v+Z4/fgfvYTij1drxbvOCr7jp2Vk3+NQK20g4K9nZB9TKDA3Y7CFbBlZVS1Yd8UdcRyHeenwX/mOFqentC16iWQnRF0v9KhAhzN6/V1NItSqcdPsmhb+oNi+jL8510ZZt20k4PeBic7exrOSmd4Eb1Iuemd2eaRnowurDST+vw3wRn/KKvivNs0BSm4OXLq0rlhi4ML5uxC0BThtNGsdQuDzHbB/2GN79S33oDz9e/XmxNvRju5HqS1bYaO3uz0gP3SuLVS5+1ZHtOsq4vM7If2NhLaMNkzaEr/b7N+U9ZLn3XKj98jfsefI5JnRkYl1qvUzBzQ9hw9Iw/4GA/l/z6NT/NvihfMk+Z7k/xsaxOoTPySv0Z7D+erVl8X6VnZFh/+PN5UvbwycZ0EysGT5xxnJbwzelYj7yf8qAqyy6iqe9j8HfLWlE3Z2tPwje89myrYR8R2euxL73u5Ajj2UeBuo+IZGH48Z1tAlpve59gv0rLrQtQJgqYI12/an2Mc+cvJHdf1BpxfL9ljruRcfvR48of3SUPb6qjqZQv7fXNZHcjSokFeevXw3tKtcnC5pM+PG80VfFpKFMrr/n8OP7JWj/53iamCyEwWJ71eTGrhr0ckORPBwnNfrTCjYHvblM7VXiCbHNir3q4smj/Ty9lvGxDVvswst2nK0pK8lusjakpvPTu9GFqQwe1OMjteTkbfP9Cxqa9qK7sVu0N/H79G/eeZscBOnEt2vB9k1FlEh3gJhGBcK/zj65Xu6v9AO3r6UnB0Y4VL7ZvCCUFThLOP6eEtvRPODj/S2VPWHcgwY7Hg10vjHOXJQkPoY38NQ1C5REG/by5rp9U3/txkde8fabfZWbcvUUYu3u0ehNOXEoXyZLpNk44m/pkf1Zr26gtv9t/pxNlkSclEK81g3X9Pq1ImrNibYfsc++eXKZIv0yeS/Llz+fxNtApZq4yGzPPcZtITqBv2WjRI1lHqfe1Wej+9rDp3IXBDQHgP/JrRBUZ2Cpf2wToWWM0+Mx0P65E7g6AprWp1f9dvoHq22TNWSRm0k1yx+w/Hw1d409qKcB8wmKApilUtWcgc8EbH7zZVLnZLdv6RamfUcZP+NmyxfLXkN4+dc/3viRneFPUb2wB1GR0QnPRSSavK0Kmh87d7nJB8ZX/s7anSu7rP6BYsytaLPFPtozK6onbKfkDTY84Fe8F8tBHyJb2sjGYdfcv4O38481fpNnGduT2KU/aS2NsxOZE3dw6fv4lWIasBM7I+9ENm3dB7hsmM2k+Ouq/+ejAwjXOt3Vf3+67j15hbmKRHBwJ1GjSlV6Xs9DZPgbqHWEYNrrPbxJWp29b5kitHjkx9L++LV6127vXDBhm/wvdwEcFqb5TWra70guRvn8bLzf2ctcGz876nYdjebijUBUDmToy6M/i66vBusDpozjbp9cNGaf6RZ51/tdKXpTpRfrZwh6z47Zips79r4KWGn95ZH71K977ppmZj7FkfJ4FGWlraPtvbYq9G6x/N+lW+WZo626bly6hNk9Mbh9oDRl+3W9F1rr32Ln2W/zv7+AzGyOk+ca0ZYuCndQfSnc9+0tXbVmibDztfAYsGqZ/M2Sb7HfSQ0/Yjvm6AbH3l9G6O/NnCnX6P+bXSdv8prf5wSk8m9pGms9wmwt4N3JX6/nvaxiMzWRZ7O7xqb0yXFh8tlEB66stl5r5dejGUXgBTo89MR2VU78/0rJrJTDd8J9ua9z6r2653taO9rZVe9PWdskk2/H7pWJQRzXw+8nl8hsGFHvOccLoJWLcF8UXLoheuPt/z+oDJ6/ab2obuabRhcxI0acCsWVzdfl/7br3PDOz8rRlX9zbsN9d9QeJk//Qu223vzfN4HZ4hk0h4jh4FD97VEte8Pk0K5MkpV11xqWu25aQtWPll+xH3sPT2DIRurOUu96y604yEve3CVtvJt2iB3Kk+R8cV2f1ua5+/lJ5A9MA14qmbMvVLphfMPPbfpbKkZzMpXSSf6R2mQaEvmklzeSXHNJCzD8yXL3dOR59rry665d258l7bmnLn9aXc07wb4nvXTmp1lH52FR+/l9M2N+NX7JNfth/NsOu/fVwmPb5aN3pOz4tj15ghE75b87vM/XcTv9u4eB84fTmTlPXMgj89wiZ69abbceiUz/3F18nR13btOXzFn896c1WtWtWqW1Wv0uVyy1UlPKojs3pPxKyyfn/7ccHaJjWbccOVRWTdvrSzQ+/8tEnmbD7k86RnZW+1G7732FOVSxSSYul0DnByUZUnVw6PbU2zpO1GLEvzuPPY8CWycf8J+eLnXbKgWxMzsvr9ta405dTjnxY3l9fOqZlP9cYPGwJyL0qn1WF6kalDSWjPSH8GJvaubbdntvR4f2vVS9ufk6Dpzg8XyLZDp6RW+aJSp0JRd41Bp6ZVTeZHL8KaXltS/plBT77vVu/z2RGlfuViqaZrcNZm6C+mw8ysrrenmVEKp7GZ7AiawpwenKq+Ns1nI2CtNvM1v32HsV+J6YFAr2ZVqcKXGi6rSav3yYhffFe9ZGasn7lbDpleJXqVenXJQrLj8OmAnDitZT9av3yGVz+r9x5LdTL/5NEb5Z6aZcwOqQdlu7Tuvv7O1M0eDeifHb1CmlUrKf0erGG61HvLbwvG9HYfXf/qcr/stWZS8rJLbXTS4ytj4ytgsk5e+tvqzWk7j8l4PJtuE9bKgIdqpRpjaufh044CE029e5+4rOqEtEZ5vr7XpUb6Fi2zHlgbVCkmeXPl9Jmi99XeLTOZO38OwNqot3Z5z0aoeiNV97L+uga2T1M6blb5ywvII58vkY5NrpINtuxroGnWrYBX0O+Ebi/aG87aH2Z3vV3OnPf9m2lg9d9Fvo8JVV6dKlvfvsv9u9m1HZZ2L0t7G6CMqgfTCs61yrFC8QKppmvAZLFGvz+ffFHa1i1nhgC4LF8uWdS9qc9t4QevHp7WcdJfTgbstKQVMKXHO5izt3/TMcy898tpG9LPTGvAZI2JpQ9fmZ/hT9Zz1BPa2/99Fp+qPBrIaqBn1UxYgZ6vwEkviqtcUTDg2desinNF0xj1ATBkyBAZMGCAJCQkSK1ateSTTz6R+vXrZ/h3J06ckCJFikhiYqIULuxfz6/0aMbmnk9+djx//7/VTDNVq8GLtZMEwtreLUw6Nk/OODPgYlq0EW5Wey8Fmu7MOnhcMMZC0Z50mrVq17Cix3rp2+YGub92WdNg1JcVrzc3VUoP1StnDjjpNU63rOl1p3y+cKff9wQsXTifJPhokK0Nr51Un+mV6Y3li7qzLKrPvddLn8mbHH2+nnC1Ctga2bvR1SXk9muuSLexfXbSTG5avRM1KF7QvYnUf8f/dhvZadZLjeVO2wlH1/m1r3sGA5pd0Ayo07Y43m0mJ/yjodzY13PYkMzQINVJezvLiKfqScMqJeS6XhkHN4/WryDf/jXg7X21ykqX5ldL+WIFzLAM/g6HkVmdm1Y1+4xecGWWbnczX2pseuulZds7rczxVi9+NKuu95zMjLW9Wkitt9KusvVHwyrFZdjf60j+PDlTbX+Wn/51m7T+2Nl5Lq1MY1b4c/4maLIZN26cPPnkk/Lpp59KgwYNZODAgTJhwgTZunWrlCxZMiRBk9a5WylkwNuDN16ZbQd+AAi18f9o6LPaLysImjJJA6WbbrpJBg8ebF5fvHhRypcvLy+88IL06NEjYCvdHzX7zMjU4JUAAESj3QHONvlz/qb33F+SkpJk5cqV0rx580srJ0cO8zo+PnUd/fnz582Ktj+C4d8trg3KcgEAiEQXQ9jcg6DpL0eOHJGUlBQpVepSryilr7V9k7d+/fqZyNR6aEYqGNrdUknCxbWlPIcpCDZtAwCR5telXzUMANktj4P7kQbD5QVy+3XXhUCj91wm9ezZU7p27ep+rZmmYAVOwWj4BgBAJEm56Ar5oJcETX8pUaKE5MyZUw4e9LyDtb4uXbp0qhWXN29e8wAAAMEX6oBJUT33lzx58kjdunVlzpxL3Yi1Ibi+btiwYah+HwAAECbINNlodVu7du2kXr16ZmwmHXLg9OnT0r59+9D9QgAAICwQNNk8/PDDcvjwYenVq5dp/F27dm2ZPn16qsbhAAAg9jC4ZYAEa5wmAAAQPIzTBAAAEGA0BAcAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAABwiaAAAAHOA2KgHicrncI4sCAIDIYJ23rfN4egiaAuTkyZPmuXz58oFaJAAAyMbzuN4OLT3cey5ALl68KPv375fLLrtM4uLiJNBRsAZje/fu5b52QcR6zh6sZ9ZzNGF7jvx1rRkmDZjKli0rOXKk32qJTFOA6IouV66cBJNuJNwMOPhYz9mD9cx6jiZsz5G9rjPKMFloCA4AAOAAQRMAAIADBE0RIG/evNK7d2/zDNZzpGN7Zj1HE7bn2FrXNAQHAABwgEwTAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBU5gbMmSIVKpUSfLlyycNGjSQZcuWhbpIUWfhwoVy7733mtFgdTT377//PtRFikr9+vWTm266yYyaX7JkSWnTpo1s3bo11MWKOsOGDZOaNWu6BwBs2LChTJs2LdTFinrvvvuuOX506dIl1EWJKn369DHr1f6oVq1ayMpD0BTGxo0bJ127djVdLFetWiW1atWSli1byqFDh0JdtKhy+vRps241QEXwLFiwQDp16iRLliyRWbNmSXJysrRo0cKsfwSO3plAT+ArV66UFStWyB133CH333+/bNy4kdUcJMuXL5fPPvvMBKsIvOrVq8uBAwfcj59//llChSEHwphmlvTKfPDgwe772+l9d1544QXp0aNHqIsXlfQq5rvvvjNZEATX4cOHTcZJg6nGjRuzuoOoWLFiMmDAAOnQoQPrOcBOnTolderUkaFDh8rbb78ttWvXloEDB7KeA5hp0uz/mjVrJByQaQpTSUlJ5kqxefPmHve309fx8fEhLRsQCImJie4TOoIjJSVFxo4da7J5Wk2HwNPsaevWrT2O1Qisbdu2meYTVapUkccff1z27NkT4E9wjhv2hqkjR46YA16pUqU8puvrLVu2hKxcQCBo1lTbftx6661yww03sFIDbP369SZIOnfunBQqVMhkT6+//nrWc4BpQKpNJ7R6DsGrcRk5cqRce+21pmruzTfflEaNGsmGDRtM+8jsRtAEICRX53rQC2XbhGimJxitztBs3sSJE6Vdu3amGpTAKXD27t0rL774ommfpx11EBytWrVy/1/bjGkQVbFiRRk/fnxIqpsJmsJUiRIlJGfOnHLw4EGP6fq6dOnSISsXkFWdO3eWKVOmmF6L2mgZgZcnTx6pWrWq+X/dunVNJmTQoEGmsTICQ5tPaKccbc9k0doB3a61Her58+fNMRyBVbRoUbnmmmtk+/btEgq0aQrjg54e7ObMmeNRpaGvaZuASORyuUzApFVFc+fOlcqVK4e6SDFDjx16EkfgNGvWzFSDakbPetSrV8+0udH/EzAFr+H9jh07pEyZMhIKZJrCmA43oGl13RHr169vemRog8727duHumhRtxPar1p27dplDnraQLlChQohLVu0VcmNGTNGfvjhB9MWISEhwUwvUqSI5M+fP9TFixo9e/Y0VRq67Z48edKs8/nz58uMGTNCXbSootuwd3u8ggULSvHixWmnF0Avv/yyGUdPq+T2799vhuDRgPTRRx+VUCBoCmMPP/yw6Zbdq1cvc4LRrqzTp09P1TgcWaNj2TRt2tQjWFUasGoDRARu0EXVpEkTj+lffvmlPPXUU6zmANEqoyeffNI0mtWAVNuBaMB05513so4Rcfbt22cCpKNHj8oVV1wht912mxnrTf8fCozTBAAA4ABtmgAAABwgaAIAAHCAoAkAAMABgiYAAAAHCJoAAAAcIGgCAABwgKAJAADAAYImAAAQthYuXGhGBS9btqzExcXJ999/n6nbOL3//vvmvnV58+aVK6+8Ut555x2/l0PQBCAm6SjkegDWR+7cuc1I+zpq9ogRI8y92gCEh9OnT0utWrVkyJAhmV7Giy++KMOHDzeB05YtW+THH380tyfzFyOCA4jZoOngwYPmNi56d3r9v96mqF+/ftKoUSNzUM2ViztNAeEkLi7O3PS7TZs27ml6M+rXXntNvv32Wzl+/Li59997773nvmXT5s2bze2ENmzYINdee22WPp9ME4CYpWn60qVLm1R9nTp15NVXXzU3FJ42bZr7voMffvih1KhRw9yMtXz58vLPf/7T3OTZugIuXLiwTJw40WO5Wn2g8+sNcwEEV+fOnSU+Pl7Gjh0r69atk4ceekjuuusu2bZtm3l/8uTJUqVKFZkyZYpUrlxZKlWqJM8884z88ccffn8WQRMA2Nxxxx2mKmDSpEl/HiRz5JCPP/5YNm7cKKNGjZK5c+dK9+7dzXsaGD3yyCMmW2Wnr//2t7/JZZddxroFgmjPnj1mf5swYYLJEF911VXy8ssvmxv7Wvvlzp075bfffjPzjB492lwQrVy50uyj/iL3DABeqlWrZq5YVZcuXdzT9Qr17bfflo4dO8rQoUPNNL1iveWWW+TAgQNSpkwZOXTokEydOlVmz57NegWCbP369aZ6XRt422mVXfHixc3/tY2ivtaAyZrviy++kLp168rWrVv9qrIjaAIAHz1ttO2E0uBH2zlp49ETJ07IhQsX5Ny5c3LmzBkpUKCAaUxavXp1k4Xq0aOHfP3111KxYkVp3Lgx6xUIMq0qz5kzp8kc6bNdoUKFzLNezGj7RHtgdd1117kzVf4ETVTPAYAXbTiqbR92794t99xzj2lE+r///c8cmK0ePElJSe75NdtktYHSKoH27du7gy4AwXPjjTeaTJNmeKtWrerx0PaK6tZbbzUXOzt27HD/3a+//mqe9QLHHwRNAGCjbZY05d+2bVsTJGlq/4MPPpCbb77ZXKnu378/1fr6+9//btpMaNunTZs2Sbt27VinQACzSWvWrDEPtWvXLvN/zRLpPvn444/Lk08+adoh6nvLli0z2eGffvrJzN+8eXPT0ePpp5+W1atXm/36H//4hxlixLtaL0MuAIhB7dq1c911112uAwcOuPbt2+dauXKl65133nEVKlTIdc8997guXLjgWrNmjUsPkwMHDnTt2LHDNXr0aNeVV15pph07dsxjeY899pgrT548ZpkAAmfevHlmn/N+6D6skpKSXL169XJVqlTJlTt3bleZMmVcDzzwgGvdunXuZfz++++uBx980OzfpUqVcj311FOuo0eP+l0WxmkCELPjNGk7JKXtHS6//HLTa+6xxx4zmSLtNac++ugjGTBggBn/RdspWVe1x44dk6JFi3pkqJo1aybjx483XZ4BRB+CJgAIgK+++kpeeuklU32XJ08e1ikQheg9BwBZoL3odLiBd99917STIGACohcNwQEgC/r372/GddKeOj179mRdAlGM6jkAAAAHyDQBAAA4QNAEAADgAEETAACAAwRNAAAADhA0AQAAOEDQBAAA4ABBEwAAgAMETQAAAA4QNAEAAEjG/h9qeSvR/S9s5AAAAABJRU5ErkJggg==", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Calculated annualized volatility: 0.6397\n", + "Target annualized volatility: 0.1500\n" + ] + } + ], + "source": [ + "import numpy as np\n", + "\n", + "# Model parameters\n", + "sigma_annual = 0.15 # annualized volatility\n", + "mu_annual = 0.02 # average annual return\n", + "b = 0.0225 # parameter for V update\n", + "# Non-model parameters\n", + "trading_days = 252 # number of trading days in a year\n", + "N = 10000 # number of days per sample\n", + "num_samples = 500 # number of samples\n", + "S0 = 100 # initial stock price\n", + "\n", + "\n", + "all_prices = []\n", + "for sample in range(num_samples):\n", + " np.random.seed(42 + sample) # Different seed for each sample\n", + " \n", + " # Daily parameters\n", + " daily_return = mu_annual / trading_days\n", + " daily_vol = sigma_annual / np.sqrt(trading_days)\n", + " V = daily_vol ** 2 # initial variance (daily)\n", + "\n", + " # Simulate stock prices with adaptive volatility starting from day 3\n", + " prices = [S0]\n", + " for i in range(N):\n", + " if i >= 2:\n", + " # Compute previous return z\n", + " z = np.log(prices[-1]) - np.log(prices[-2])\n", + " V = V * np.exp(b * np.abs(z))\n", + " V = min(V, (10 * daily_vol)**2) # cap V\n", + " sigma_i = np.sqrt(V)\n", + " Z = np.random.normal(0, 1)\n", + " log_return = daily_return - 0.5 * sigma_i**2 + sigma_i * Z\n", + " new_price = prices[-1] * np.exp(log_return)\n", + " if not np.isfinite(new_price) or new_price <= 0:\n", + " new_price = 1e-10 # reset to small positive value to prevent NaN/inf\n", + " prices.append(new_price)\n", + " \n", + " all_prices.extend(prices)\n", + "\n", + "prices = all_prices # For compatibility with second cell\n", + "\n", + "# Output the end price of the last sample\n", + "print(f\"End price after {N * num_samples} days: {prices[-1]:.2f}\")\n", + "\n", + "# Plot the stock price over time\n", + "import matplotlib.pyplot as plt \n", + "plt.plot(range(len(prices)), prices)\n", + "plt.xlabel(\"Day\")\n", + "plt.ylabel(\"Stock Price\")\n", + "plt.title(f\"Simulated Stock Price Over Time ({num_samples} samples)\")\n", + "plt.show()\n", + "\n", + "# Calculate annualized volatility over all\n", + "daily_returns = [np.log(prices[i+1] / prices[i]) for i in range(len(prices)-1)]\n", + "annualized_volatility = np.std(daily_returns) * np.sqrt(252)\n", + "print(f\"Calculated annualized volatility: {annualized_volatility:.4f}\")\n", + "print(f\"Target annualized volatility: {sigma_annual:.4f}\")" + ] + }, + { + "cell_type": "code", + "execution_count": 183, + "id": "3ff1f563", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " → 3854789 clean windows\n" + ] + }, + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
tickerdateTsigmaz
0simulated450.0950280.155202
1simulated950.1102070.147910
2simulated1450.123517-0.289159
3simulated1950.087454-0.239317
4simulated2450.141630-0.043257
..................
3854784simulated49999291300.688426-0.906785
3854785simulated50000591300.7127720.588940
3854786simulated50001891300.7812630.377564
3854787simulated50003191300.802100-0.739435
3854788simulated50004491300.8124470.324417
\n", + "

3854789 rows × 5 columns

\n", + "
" + ], + "text/plain": [ + " ticker date T sigma z\n", + "0 simulated 4 5 0.095028 0.155202\n", + "1 simulated 9 5 0.110207 0.147910\n", + "2 simulated 14 5 0.123517 -0.289159\n", + "3 simulated 19 5 0.087454 -0.239317\n", + "4 simulated 24 5 0.141630 -0.043257\n", + "... ... ... ... ... ...\n", + "3854784 simulated 4999929 130 0.688426 -0.906785\n", + "3854785 simulated 5000059 130 0.712772 0.588940\n", + "3854786 simulated 5000189 130 0.781263 0.377564\n", + "3854787 simulated 5000319 130 0.802100 -0.739435\n", + "3854788 simulated 5000449 130 0.812447 0.324417\n", + "\n", + "[3854789 rows x 5 columns]" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Saved to dataset_part1.parquet, dataset_part2.parquet, dataset_part3.parquet\n" + ] + } + ], + "source": [ + "import pandas as pd\n", + "\n", + "HORIZONS = 5*(np.arange(26)+1) # does 1 to 26 weeks, can also do [5, 10, 20, 40, 80, 160]\n", + "\n", + "# Create price series from simulated prices\n", + "price = pd.Series(prices, index=range(len(prices)))\n", + "\n", + "ticker = \"simulated\"\n", + "\n", + "ret = np.log(price).diff().dropna().values\n", + "\n", + "rows = []\n", + "scale = np.sqrt(252)\n", + "\n", + "for T in HORIZONS:\n", + " i = 0\n", + " while i + T <= len(ret):\n", + " window = ret[i:i+T] # T points\n", + " if len(window) < T * 0.8:\n", + " break\n", + "\n", + " x = window.sum() # total price change over the period\n", + " sigma = np.std(window, ddof=0) * scale # use np.std to get std over period, ddof=1 means divisor is N-1\n", + " z_raw = x / np.sqrt(T / 252.0)\n", + "\n", + " # REJECT BAD WINDOWS\n", + " if not (np.isfinite(sigma) and sigma > 0 and np.isfinite(z_raw)):\n", + " i += T\n", + " continue\n", + "\n", + " rows.append({ # append row of data for this period\n", + " \"ticker\": ticker,\n", + " \"date\": i + T - 1,\n", + " \"T\": T,\n", + " \"z_raw\": float(z_raw),\n", + " \"sigma\": float(sigma)\n", + " })\n", + " i += T\n", + "\n", + "if not rows:\n", + " print(\" [no data]\")\n", + "else:\n", + " df = pd.DataFrame(rows)\n", + "\n", + " # CLEAN BEFORE DE-MEANING \n", + " df = df[np.isfinite(df['z_raw']) & np.isfinite(df['sigma']) & (df['sigma'] > 0)]\n", + "\n", + " # NOW de-mean safely, this step groups by ticker and T, and subtracts the group mean \n", + " df[\"z\"] = df.groupby([\"ticker\", \"T\"])[\"z_raw\"].transform(lambda g: g - g.mean())\n", + "\n", + " df = df.drop(columns=\"z_raw\")\n", + " df = df.dropna().reset_index(drop=True) # Final clean\n", + "\n", + " print(f\" → {len(df)} clean windows\")\n", + " display(df)\n", + "\n", + " # Split into 3 parts\n", + " n = len(df) // 3\n", + " part1 = df.iloc[:n]\n", + " part2 = df.iloc[n:2*n]\n", + " part3 = df.iloc[2*n:]\n", + "\n", + " # Save 3 small files\n", + " part1.to_parquet(\"dataset_part1.parquet\", compression=None)\n", + " part2.to_parquet(\"dataset_part2.parquet\", compression=None)\n", + " part3.to_parquet(\"dataset_part3.parquet\", compression=None)\n", + " print(\"Saved to dataset_part1.parquet, dataset_part2.parquet, dataset_part3.parquet\")" + ] + } + ], + "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.14.2" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}