{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Segreg Example All Models" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from matplotlib import pyplot as plt\n", "import matplotlib\n", "import numpy as np\n", "from numpy.random import default_rng\n", "\n", "from segreg import analysis\n", "from segreg.model import OLSRegressionEstimator\n", "from segreg.model import OneBkptSegRegEstimator\n", "from segreg.model import TwoBkptSegRegEstimator" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "matplotlib.rcParams['legend.numpoints'] = 1\n", "matplotlib.rcParams['figure.figsize'] = (10,5)\n", "matplotlib.rcParams['axes.grid'] = True" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Create a Dataset Based on a Two-Bkpt Function" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAEvCAYAAAB2a9QGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA/gklEQVR4nO3dd3xUZdr/8c9NSGihl4QSadJCSFBCEVtQRCyICrGs6yqusu7q6vooCuoqViy4ru6662ML+9v1EQ0oIlYsWeyUXULvoIkh9DZASLt/f6QIyQQmmXLOZL7v14sXM2fOmXNN7szMlesux1hrEREREZHAauB0ACIiIiL1kZIsERERkSBQkiUiIiISBEqyRERERIJASZaIiIhIECjJEhEREQmChk4H4E27du1st27dgnqOgwcP0qxZs6CeQ2pHbeJOahf3UZu4k9rFfULVJkuWLNlprW1fdbsrk6xu3bqxePHioJ4jKyuLtLS0oJ5Dakdt4k5qF/dRm7iT2sV9QtUmxpgfvG1Xd6GIiIhIECjJEhEREQkCJVkiIiIiQeDKMVkiIiLeFBUVkZubS0FBgdOhVNOyZUtWr17tdBhylEC3SePGjenSpQvR0dE+7a8kS0REwkZubi7NmzenW7duGGOcDucYBw4coHnz5k6HIUcJZJtYa9m1axe5ubl0797dp2PUXSgiImGjoKCAtm3bui7BkvrPGEPbtm1rVUVVkiUiImFFCZY4pba/e0qyRERE6mjq1KlMnz79uPvMmTOHVatWhSii4JgxYwZ5eXmV92+88caAvabnn3+efv36cc011wTk+QC2bNnC//3f/1XeX7x4MbfddlvAnt9XGpMl4mIej4fMzEzy8/OJj48nPT2d2NhYp8MSkVqYM2cOF198MYmJiU6HclwlJSVERUV5fWzGjBkkJSXRqVMnAF555ZWAnfdvf/sbH374oc/jnHxRkWSNGTMGgNTUVFJTUwP2/L5SJUvExTIzM8nJyaGwsJCcnBwyMzOdDkkk4j322GP06dOHkSNHsnbt2srtM2bMYPDgwaSkpDBu3DgOHTrEN998w9y5c5k0aRIDBw5k48aNvPzyy9X28+bVV1+ld+/epKWlcdNNN3HrrbcCsGPHDsaNG8fgwYMZPHgwX3/9NVBWVbvhhhtIS0ujR48ePP/885XP9a9//YshQ4YwcOBAfvOb31BSUgJAbGwsDzzwAEOHDuXbb7/l4YcfZvDgwSQlJTFx4kSstcyaNYvFixdzzTXXMHDgQA4fPkxaWlrllVliY2O57777SElJYdiwYWzbtg2AjRs3MmzYMAYPHswDDzzg9Q/Em2++mU2bNnHJJZfw7LPPVqsMJiUlsWXLFrZs2UK/fv246aab6N+/P6NGjeLw4cMAbNiwgZEjR5KSksKpp57Kxo0bmTx5Ml9++SWnn346zz77LFlZWVx88cUA7N69m0svvZTk5GSGDRvGsmXLTvjzqyslWSIulp+fj7UWKJvZkp+f73BEIpFtyZIlzJw5k//+97+8/fbbLFq0qPKxMWPGsGjRIrKzs+nXrx+vvvoqw4cP55JLLuHpp59m6dKl9OzZk8svv7zaflXl5eXxyCOP8N133zF//nzWrFlT+djtt9/OHXfcwaJFi5g9ezY33nhj5WNr1qzh448/ZuHChTz00EMUFRWxevVq3nzzTb7++muWLl1KVFQUr7/+OlB2bb+kpCS+//57zjjjDG699VYWLVrEihUrOHz4MPPmzWP8+PGkpqby+uuvs3TpUpo0aXJMrAcPHmTYsGFkZ2dz1lln8fLLL1fGefvtt7No0aLKClhVL774Ip06deKLL77gjjvuOO7Pfv369dxyyy2sXLmSVq1aMXv2bACuueYabrnlFrKzs/nmm2/o2LEjTzzxBGeeeSZff/11ted98MEHOeWUU1i2bBmPP/44v/rVr4778/OHugtFXCw+Pp6cnBystRhjiI+PdzokEdd46L2VrMrbH9DnTOzUggfH9K/x8S+//JLLLruMpk2bAnDJJZdUPrZ69WquvfZa9u7di8fj4fzzz/f6HCtWrOD+++8/7n4LFy7k7LPPpk2bNgCkp6ezbt06AD799NNjxkPt37+fAwcOAHDRRRfRqFEjGjVqRIcOHdi2bRufffYZS5YsYfDgwQAcPnyYDh06ABAVFcW4ceMqn+uLL77gqaee4tChQ+zevZv+/ftXdrnVJCYmprJKNGjQIObPnw/At99+y5w5cwD4xS9+wV133XXc5zmR7t27M3DgwMrzbNmyhQMHDvDTTz9x2WWXAWXrWJ3IV199VZmgnXPOOezatYt9+/YB3n9+Xbp0qXPMSrJEXCw9Pb3amCwRcVZNM8x++9vf8u6775KSksKMGTPIysryut/111/PnDlzjtmvpKSEQYMGAWWJ2ymnnFLj+UtLS/n222+rVZQAGjVqVHk7KiqK4uJirLVcd911TJs2rdr+jRs3rhyHVVBQwO9+9zsWL15MQkICU6dO9Wm5gujo6MqfScU566phw4aUlpZW3j/6/FVf2+HDhysr/bXh7ZiK+L39/PyhJEvExWJjY5kwYYLTYYi40vEqTsFy1llncf311zN58mSKi4t57733+M1vfgOULXzZsWNHioqKeP311+ncuTMAzZs3r6w01bRfVFQUS5curdznp59+4o477mDPnj00b96c2bNnM2DAAABGjRrFX//6VyZNmgTA0qVLKys83px77rmMHTuWO+64gw4dOrB7924OHDhA165dj9mvIqFp164dHo+HWbNmMX78eK+vwRfDhg1j9uzZXHnllcycOdOnY7p168a8efMA+M9//sPmzZuPu3+LFi3o0qULc+bM4dJLL+XIkSOUlJQcN96zzjqL119/nT/+8Y9kZWXRrl07WrRoUavX5iuNyRIREfHRqaeeypVXXsnAgQMZN24cZ555ZuVj999/P0OHDuW8886jb9++lduvuuoqnn76aU455RQ2btzII4884nW/o3Xu3Jl7772XoUOHMnLkSBITE2nZsiVQtuTB4sWLSU5OJjExkRdffPG4MScmJvLoo48yatQokpOTOe+889i6dWu1/Vq1asVNN93EgAEDuPTSSyu7F6Gs+nbzzTdXDnz3xZ///Gf+9Kc/MWTIELZu3VoZ//GMGzeO3bt3M3DgQP7+97/Tu3fvEx7zz3/+k+eff57k5GSGDx9Ofn4+ycnJNGzYkOHDh/Pss88es//UqVMrf36TJ0/mH//4h0+vpy5MXUptwZaammorZi0ES1ZWFmlpaUE9h9SO2sRdKpaPaNKkCYcPH9byES4Sye+V1atX069fP6fD8CrQl9XxeDzExsZSXFzMZZddxg033FA59igcHDp0iCZNmmCMYebMmbzxxhu8++67IY0hGJc68vY7aIxZYq2ttkaEKlki4lXF8hHWWi0fIeKAqVOnMnDgQJKSkujevTuXXnqp0yHVypIlSxg4cCDJycn87W9/45lnnnE6pJDTmCwR8UrLR4g460QrybvdmWeeSXZ2ttNhOEqVLBHxKj4+vnLGjZaPEBGpPSVZIuJVeno6CQkJGGNISEjQ8hHiGm4cSyyRoba/eydMsowxCcaYL4wxq40xK40xt5dvb2OMmW+MWV/+f+sajt9ijFlujFlqjAnuaHYRCZiK5SPi4+OZMGGCBr2LKzRu3Jhdu3Yp0ZKQs9aya9cunxY8reDLmKxi4E5r7X+MMc2BJcaY+cD1wGfW2ieMMZOBycA9NTzHCGvtTp+jEhER8aJLly7k5uayY8cOp0OppqCgoFZfwBJ8gW6Txo0b12oF+BMmWdbarcDW8tsHjDGrgc7AWCCtfLd/AFnUnGSJiIj4LTo6mu7duzsdhldZWVnHXaldQs/pNqnVmCxjTDfgFOB7IK48AatIxDrUcJgFPjHGLDHGTPQjVhEREZGw4fNipMaYWODfwGPW2reNMXutta2OenyPtbbauCxjTCdrbZ4xpgMwH/i9tXaBl/0mAhMB4uLiBvm6BH9dVSzyJu6hNnEntYv7qE3cSe3iPqFqkxEjRnhdjNSndbKMMdHAbOB1a+3b5Zu3GWM6Wmu3GmM6Atu9HWutzSv/f7sx5h1gCFAtybLWvgS8BGUrvgd7NeNIXjHZrdQm7qR2cR+1iTupXdzH6TbxZXahAV4FVltr/3TUQ3OB68pvXwdUWyvfGNOsfLA8xphmwChghb9Bi4iIiLidL2OyTgeuBc4pX4ZhqTHmQuAJ4DxjzHrgvPL7GGM6GWM+KD82DvjKGJMNLATet9Z+FPBXISIiIuIyvswu/AowNTx8rpf984ALy29vAlL8CVBEpL6ruBh3fn4+8fHxuhi3SD2hFd9FRBxWcTHuwsJCXYxbpB7RBaJFRBwW6ItxqzIm4g6qZImIOCzQF+NWZcx5Ho+HjIwMpk2bRkZGBh6Px+mQxAFKskREHFZxMe6YmJiAXIw70JUxqT0lugLqLhQRcVzFxbgDJT4+npycHKy1AamMSe0p0RVQkiUiEta8jb9KT0+vtk1CS4mugJIsEZGwVtEtZa2t7JaaMGFCQCtjUntKdAWUZImIVArHWXnqlnKnQHcBS3jSwHdxJc3MESeE42DlQM9MFJHAUZIlrhSOX3YS/vypCjn1h0GgZyaKSOCou1BcSV0g4gR/BivXNDYq2NQtJeJeqmSJK6kLRALJ1yqTP1Uh/WEgIlWpkiWuFA4zc8JxkHSk8rXK5E9VSFP2RaQqJVniSuHQBeJU95DUXiiqTOHwh4GIhJaSLJE6UvdQ+AhFlcnbHwaqdopENo3JEqkjjRsLH07NwKvvs2S11IrI8amSJVJH6h4KH051P9f3aqe6zMtUVCybNGlCRkaGKpZSSUmWSB2Fw7gxcVZ9Hwxf35NIX1Ukm7169YroZFOqU5IljguHvwLdNLbGTbHI8flT7QyHdq7vSaSvlGxKTTQmSxznrcvBbdw0tsZNsQRDfRrnU1HtnDJlChMmTKhVkhQO7ezPWLf61M4anyk1USVLHBcOfwW6KUY3xRIMbhvn41RFKRza2Z8uc3/a2W1VvoqKpTFGlzaSY6iSJY4Lh78C3RSjm2IJBrclF05VlNTONXNbla8i2YyPj691xVLqNyVZ4riKLgc3/xXopovwuimWYHBbcuFU0qd2rpnbEnGRmqi7UBxX8VdgVlYWaWlpTofjlZtmEroplpr4053jtqUxnBrcHQ7t7I2vbe9PO2vAvYQLJVki4hdvX6r+jLdxW3LhtqTP7UJxnUi1iYQLJVki4hdvX6r1qTvHbUmf24Wi7dUmEi40JktE/OLtS9Vt46okdNT2Ij9TJUsimtumgocjb+Nj1J0THP4u3BuK3/f6vgCrSG0oyZKI5rY1mcKRty9VdecEh7+XbwnF77tTa2f5SomchNIJkyxjTALw/4B4oBR4yVr7nDGmDfAm0A3YAlxhrd3j5fjRwHNAFPCKtfaJgEUv4qf6NHbIH/588SihCh1/f1/d/vseivj0h5WEki9jsoqBO621/YBhwC3GmERgMvCZtbYX8Fn5/WMYY6KAF4ALgETg6vJjRVwhFONHwuHyIW5b3FG88/f31e3jpUIRn9sTTalfTljJstZuBbaW3z5gjFkNdAbGAmnlu/0DyALuqXL4EGCDtXYTgDFmZvlxqwIQu9QjRSWl/Hd7Mfuz84J2joKCAr7/7jv27t1Lq1atGDpsGM0Sz2bf/vJtLVvRN3EYc/2Iwds5vv/uO3bt8mBtM9Zu9rDi1bmc7bL1wL768TDFxa0q7//w42Halv8cVm0Nbrv4wtvPtXHjxo7G5ISK39dle6MpbHlyrX9fA/37HmihiG9n067sOryrcgxh26ZtA3KOKGNoWGwDEKHUJ6Yio/dpZ2O6AQuAJOBHa22rox7bY61tXWX/8cBoa+2N5fevBYZaa2/18twTgYkAcXFxg2bOnFnrF1MbHo9H/fAu8s9VR/jsx2KnwxARqbO+rSz3DG1WWY0T54Xqu37EiBFLrLWpVbf7PPDdGBMLzAb+YK3d7+MvkbedvGZ11tqXgJcAUlNTbbBX/nbz6uKR5psNO/nso+85J6Eh96afHpDnfOftt9mavxVrwRjoGN+RHTt3UlRUVLlPdHQ0EydODMj5Krz00kvVztG+XbtqsVx2+eUBPa+/Dh06xMcffcSOnTtp364d548eTdOmTQFYuHAhQ4YMcTQ+bz/XQLddOHFDm0QCb58jNb1356/axpMfrSG/WU+uHnJSiCOVmjj9Xe9TkmWMiaYswXrdWvt2+eZtxpiO1tqtxpiOwHYvh+YCCUfd7wK4pzYtjvMcKWbSrGV0b9eMK/taTu4QmL84ivf8REtTWJnmF+/5iX6dj11qIKFzh4Cdr0K/zq2rnSM9/fIqg8ovd2EVNZbkm3/l9ZHc2AYB/znVlrefq9MxOckNbRIJvH2O1PRz79GuGXMXruPReas4s1c7urRuGsJIxa1OOPDdlJWsXgVWW2v/dNRDc4Hrym9fB7zr5fBFQC9jTHdjTAxwVflxIgA8/sFq8vYdZnp6Mo2iAldi9zaANhQX3PV2jorZd1OmTGHChAkuTLDcr75fLFncqTYD8Rs0MNyQ1AiAu2cto7RU47PEt0rW6cC1wHJjzNLybfcCTwBvGWN+DfwIpAMYYzpRtlTDhdbaYmPMrcDHlC3h8Jq1dmWAX4OEqQXrdvB/3//IxLN6MKhrG7I2B+65nVq7ScsZ1MzbMhGAT0tH6OcqdeHvmli1XVi1fdMG3HdRIve+s5zXF/7ItcO6+vsSJMz5MrvwK7yPrQI418v+ecCFR93/APigrgFK/bS/oIjJs5fRs30z/ue83gF/fn0pu4+39YkArVkkARHoC5VD3T5Hrh6SwIcrtjLtg9Wc3as9J7VVt2Ek07ULxRGPzVtN/v4CnrliII2jo5wOR0LA2/pEWrPIncJhbbeqvK315sTvlzGGJ8clE2UMd83KVrdhhFOSJSH3xdrtvLk4h5vP7snAhFZOh1NNOH7BhANv41vcvjhmpArHxWnddKHyTq2a8McxiSzcvJsZ32wJyTnFnZRkSUjtO1TWTdgnrjm3j+zldDheheMXTDjwNnhdA9rdKRwrjE5NdqlJ+qAujOjTnqc+XsPmnQdDdl5xF10gWkLqoXkr2ekp5NXrBtOooTu7CcPxCyYc1DS+RWOw3Cc+/tjlTsKhwui2C5UbY3hiXDLn/enfTMrM5s3fnEZUAy1SGmmUZEnIzF+1jbf/8xO3nduLpM4tnQ6nRuH4BSM/83dGmdR+Vp0buHGyS1yLxjw0tj93vJlNxtebufHMHk6HJCGmJEtCYs/BQu59Zzn9Orbg1hEnOx3OcYXjF0wkqimZ8ndGmbgzYQlXlw7szAfL83nq47Wk9YnsRXQjkcZkSUg8OHclew4W8kx6CjEN3f1rp8VDw0NNY+fU3StuYozhscuSaBoTxZ2Z2RSXlDodkoSQu7/tpF74cPlW5mbncdu5vUjs1CLgz6/ZgJGppmRKMxbFbTo0b8zDY5PIztnLy18GcNVlcT0lWRJUuzxHuH/OCpI6t+C3aT2Dcg7NBoxMNSVTmrEobjQmuSMXDojn2fnrWLftgNPhSIhoTJYE1QPvrmR/QRH/lz6M6Kjg5PTqHopMNY2d03ii+qW+TGQwxvDI2CS+37SAO9/K5u3fDa/xM7G+vGZRkiVBNG9ZHu8v38qk8/vQJ7550M5TX2YD6oO1dpRMRQa3TWTw533aNrYRj16axG9f/w8vZm3k9+d6XyvQba9Z6k7dhRIUOw4c4Y9zVpCS0IrfnBXcacv1pXtI3Z4i1bmtUu3v+/SCAR0Zk9KJ5z9fz6q8/V73cdtrlrpTJUsCzlrLve8s52BhCc+kJ9MwSN2EFepLRcPfD1ZVwqQ+ClWl2tf3TyASoIcv6c+3G3dxV2Y2c245vdqM6/pSnRdVsiQI3l2ax/xV27hrVG9O7hC8bsL6xt9ZcaqESX0Uqkq1r++fQMxebd0shscvS2LV1v288MWGao/Xl+q8qJIlAbZtfwEPzl3JqSe14tdnaHXj2vB3EdRQdDGoWha5nGr7UFWqfX3/+PM+rfozHDNgAC98sYHzEuOOuQpGfanOi5IsCSBrLfe+vZyCohKmp6foOl3HUdMXlj8frKHoYtCAXP+Fa6Ja39ve1/ePP+/Tqj/Dfp0M3zfrwp1vZTP396e79nquUnfqLpSAmf2fn/hszXbuHt2XHu3d/6XhpGB07YWii0EDcv0Xrt269b3tnXj/7NuxlSfGDWDttgM89+n6gJ9PnKdKlgTE1n2Heei9lQzu1poJw7s5HY7rBeMLKxRdDBqQ679wTVbqe9s79f45p28c6YO68OK/N3J+/3hSEloFNQYJLVWyxG/WWu6ZvZziEsv09BQaqJvwhML10i8akOs/tX3kqulneP/FicS1aMydmdkUFJU4HKUEkipZ4rc3F+WwYN0OHh7bn65tmzkdTljwd5C7UzQgt4w/46rU9pGrpp9hyybRPDEumeteW8izn65jygX9HIhOgkFJlvgld88hHn1/Naf1aMsvh3Z1OpywoS+s8ObPIHC3tX0oBuKH62D/UDq7d3uuHpLAyws2MSoxnkFdWzsdkgSAugulzsq6CZdhreWp8cnqJpSIEa7jqrwJxUD8cB3sH2r3XZRIx5ZNuCszm8OF6jasD5RkSZ29/v2PfL1hF/de1I+ENk2dDkekVjweDxkZGUybNo2MjAw8Ho/Px4bruCpvQpEw1qekNJhiGzXkqfHJbN55kOmfrHU6HAkAJVlSJzm7D/H4B6s54+R2/GLISU6HIzXwJ5Go7/yprtSnQeChSBjrU1IabKef3I5rh3Xlta83s3DzbqfDET8pyZJaKy21TJqVTQNjeHJ8cuWHp7iPumlq5k91pWJc1ZQpU5gwYUJYjy8KRcJYn5LSUJh8QV+6tG7CpFnZHCosdjoc8YMGvkut/b9vt/Ddpt08OW4AnVs1cTocOQ5109Ssvq/75KtQDMR322B/t2vWqCFPj0/hqpe+46mP1jL1kv5OhyR1pEpWBAhkl9GWnQd54qM1pPVpzxWpCQGMUoJB3TQ1U3VF3GxYj7ZMOL0bM77ZwrcbdzkdjtSRKlkRIFDXHCsptdyVmU10VAOeuNyZbkJNBa+dcF2TKRRUXRG3u/v8vnyxZjuTZmXz8R/OolkjfWWHmxNWsowxrxljthtjVhy1LcUY860xZrkx5j1jTIsajt1Svs9SY8ziQAYuvgtUl1HG15tZ/MMepo7pT3zLxoEM0WcaY1Q79WnskEikaRITxfT0FH7ae5hpH652OhypA1+6C2cAo6tsewWYbK0dALwDTDrO8SOstQOttal1C1H8FYguo407PDz98VpG9uvA5ad2DnSIPtMYIxGJJKnd2nDjGd3513c/8tX6nU6HI7V0wiTLWrsAqDqPtA+woPz2fGBcgOOSAPJ37ElFN2Hj6Cgev2yAo7MJNcZIRCLNnaP60LN9M+6ZvYwDBUVOhyO1UNeB7yuAS8pvpwM1jYC2wCfGmCXGmIl1PJf4yd8uo5e/3MR/f9zLw2P706GFM92EFTRYWUQiTePosm7DrfsO8+g8dRuGE1PR9XLcnYzpBsyz1iaV3+8LPA+0BeYCt1lr23o5rpO1Ns8Y04Gyitfvyytj3s4xEZgIEBcXN2jmzJl1e0U+8ng8Gp/ig588pTz49WFSOkRx68BGflexSktL2bNnD0VFRURHR9O6dWsaNCjL9dUm7qR2cR+1iTsFu10y1xby/uYi/mdQI5LbaxC8L0L1XhkxYsQSb8Oi6pRkVXmsN/Ava+2QEzzHVMBjrZ1+ovOlpqbaxYuDO04+KyuLtLS0oJ4j3BWXlHL5378hd89hPrnjLNrFNvL7OTMyMo5ZmyghIaFyhpfaJDRqO0NT7eI+ahN3Cna7HCkuYcxfvmL/4WI+vuMsWjaJDtq56otQvVeMMV6TrDp1F5ZXpjDGNADuB170sk8zY0zzitvAKMq6GSVM/O+CTSzL3ccjY5MCkmCBBq67gWZoioSnRg3Lug13eI7wyLxVTocjPvBlCYc3gG+BPsaYXGPMr4GrjTHrgDVAHpBRvm8nY8wH5YfGAV8ZY7KBhcD71tqPgvEiJPDW5O/nz5+u46LkjlyU3DFgz6uB685ToisSvpK7tOJ3aT2ZtSSXz1ZvczocOYETdupaa6+u4aHnvOybB1xYfnsTkOJXdOKIopJS7nwrm5ZNonlkbLUeYr9ocUzn6XIyIuHn6G7+9nHx9O7QnclvL2f+Ha1p1TTG6fCkBho5J9X87YuNrMzbz4u/HESbZoF982qVbecp0RUJP0dfuSMvN4e09g15bWdbps5dyZ+vOsXp8KQGSrLkGCvz9vGXz9dz6cBOjE5ShaM+UqIrEn6qdvObvbnces5Q/vzpekYnddTntUspyZJKhcVl3YRtmsXoqu8iIiHg62xfb938vxxxMvNXbeO+d5YzuFtr2gZogpIETl0XI5Uw5/F4yMjIYNq0aWRkZODxePjL5+tZk3+AaZcPqOzj97afiIgEhq+zfb0txBwd1YBnrkhhf0ERD8xdGeLIxReqZEWoo/v3c3JyeO6f7/BaTlvGD+rCuf3iatwvMzPTa1dTbddeEhER32f71tTN3ze+BX8Y2ZunP17LBUl5XJzcqdo++nx2jipZEeroN3ZxKczKaUr72Eb88eLEGvc73geA1l4SEam9QCxr85uzepDcpSV/nLOCHQeOVHtcn8/OUZIVoY5+Y/+3uDN7ShvzxLgB1VYQ9vUDQGsviYjUXiCux9owqgHPpKdw8EgJ989ZTtUruejz2TlKsiJUxRt7T4NWrCiOY9zAeNL6dKhxvxN9AGiRURGR2qvoBpwyZQoTJkyoczder7jm/M+o3ny8chtzs/OOeUyfz85RkhWhYmNjufqXv2JZ0xQ6tmzC1EuTa9zPlw+AQPw1JiIidXfTmT045aRWPPDuSrbvL6jcrs9n52jgewR75pO1bNpxkH/9eijNG/t3oVGtvSQiElreBrRPT0/hwue+ZMrby3nlulSMMfp8dpAqWRFq8ZbdvPLVZn4x9CTO6NXO6XBERKSWvA1o79k+lknn9+GzNdt5+z8/OR1ixFOSFYEOFRZzV2Y2nVs14d4L+zkdjoiI1EFNA9onnN6d1K6tmfreSvL3FRzvKSTIlGRFoKc+WsuWXYd4enwKsY3UYywiEo5qGtAe1cAwPT2FopJSJr+9rNpsQwkdJVkR5rtNu5jxzRauH96N03q2dTocERGpo+MNaO/WrhmTR/cla+0OMhfnOhhlZFMZI4IcPFLMpFnZdG3blLtH93E6HBER8cOJBrT/6rRufLgin0fmreL0Xu3o3KpJCKMTUCUrojzx4Rpy9xzm6fEpNI1Rfi0iUp81aGB4enwKJdYyeba6DZ2gJCtCfL1hJ//87gduOL07Q7q3cTocEREJgZPaNmXKhf34cv1O3liY43Q4EUdJVgQ4UFDE3bOW0aNdMyadr25CEZFIcs2Qkzj95LY89v4qcnYfcjqciKIkKwI8/sFqtu47zNPpKTSOjnI6HBERCaEGDQxPjkvGGMPds5ZRWup/t6HH4yEjI4Np06aRkZGBx+MJQKT1j5Kseu7f63bwxsIcbjqrB4O6tg7pufUmFBFxhy6tm3L/Rf34dtMu/vX9D34/n7eFUKU6JVn12L7DRdwzaxknd4jljpG9Q35+vQlFRNzjysEJnNW7PdM+WMMPuw769Vw1LYQqx1KSVY89Om8VOzxHeMahbkK9CUVE3MMYw5PjBtAwyjAp079uw5oWQpVjKcmqpz5fs43MJbncfHYPUhJaORKD3oQiIu7SsWUTHhzTn4VbdpPxzZY6P8/xFkKVn2mxpHpo36EiJs9eTt/45tx2bi/H4khPT692hXgREXHWuFM78+HyrTz10RpG9GlPj/axtX6OEy2EKmWUZNVDD723kt0HC3nt+sE0auh7N6HH46mWFMXG1v7NV0FvQhER9zHG8PjlAxj17ALuyswm8+bhRDUwTodVL6m7sJ75ZGU+b//3J3434mSSOres1bEaqC4iEhniWjTmoUv6858f9/LqV5ucDqfeUpJVj+w5WMi976wgsWMLbh1xcq2P10B1EZHIMXZgJ0YlxjH9k3Vs2H7A6XDqJSVZ9cgDc1ey73Ah09NTiGlY+6bVQHURkchhjOGxywbQLCaKO9/Kprik1OmQ6h0lWfXEB8u38l52Hred04vETi3q9ByaLSIiElnaN2/EI5cmkZ27j/9doG7DQDvhwHdjzGvAxcB2a21S+bYU4EUgFtgCXGOt3e/l2NHAc0AU8Iq19onAhS4VdnqOcP+cFQzo3JKb03rW+Xk0UF1EJPJcnNyJD5fn8+z8dRxY+w12z08BmfwkvlWyZgCjq2x7BZhsrR0AvANMqnqQMSYKeAG4AEgErjbGJPoVbQSr6RI11lr+OGcFnoJinrkiheio4BcndbkcEZH65eGx/YkxJWT+0JiCI0Wa/BQgJ/xGttYuAHZX2dwHWFB+ez4wzsuhQ4AN1tpN1tpCYCYw1o9YI1pNM//eW7aVD1fkc8d5vekd19zRWEREJDy1jW3E6Y1+ZJdtxrLieE1+CpC6lj1WAJeU304HErzs0xnIOep+bvk2qQNvM/+2HyjggXdXMDChFTed2d3RWEREJLwNT2hKj6hdLC3uyG7bVJOfAsBUfFkedydjugHzjhqT1Rd4HmgLzAVus9a2rXJMOnC+tfbG8vvXAkOstb+v4RwTgYkAcXFxg2bOnFnX1+QTj8cTVn3Nu3btorCwsPJ+dHQMr//YjOU7S3h4eBM6xYZuDkPVWGJiYmjbtu1xjvBNuLVJpFC7uI/axJ3CvV1KS0vJ3bGHp5fH0Dwapp7elJhaLGjtRqFqkxEjRiyx1qZW3V6nFd+ttWuAUQDGmN7ARV52y+XYClcXIO84z/kS8BJAamqqTUtLq0toPsvKyiLY5wikqqux020I/128hvsu7McvzurhaCxjxowJyC9xuLVJpFC7uI/axJ3qS7t0HLCNm/7fYlbYBP4nrbfT4fjF6TapU5JljOlgrd1ujGkA3E/ZTMOqFgG9jDHdgZ+Aq4Bf1DnSCHf0zL9t+ws470//ZlDX1txwRui6Cb3FIiIi9ct5iXFcfmpnXvhiA6MS42p99RD52Qn7mIwxbwDfAn2MMbnGmF9TNlNwHbCGsupURvm+nYwxHwBYa4uBW4GPgdXAW9balcF5GZHDWsuUt5dTWFLK9PQUXW9KREQC7sGL+9MuNoY738rmSHGJT8do5nl1vswuvNpa29FaG22t7WKtfdVa+5y1tnf5v8m2fGCXtTbPWnvhUcd+UL5PT2vtY8F8IZFi1pJcPl+znXtG96V7u2ZOhyMiIvVQy6bRPHF5Mmu3HeD5z9b7dIxmnlenFd/DyNZ9h3n4vVUM6d6G607r5nQ4IiJSj43o24ErUrvw96yNLM3Ze8L9NfO8OiVZIeJvGdVayz2zl1NiLdPHp9BA3YQiIhJk91+cSHyLxtz51lIKio7fbajr31anJCtE/C2jvrkohwXrdjD5gr6c1LZpkKIUERH5WYvG0Tw5PpmNOw7y7Px1x91X17+trk6zC6X2/Cmj5u45xKPvr2Z4z7b8cmjXYIUoIiJSzZm92vOLoSfx0pebGNU/jkFd23jdTzPPq1MlK0TqWkYtLbXcPWsZ1lqeGp+sbkIREQm5ey/sR6eWTbgrcxmHC32bbShKskKmrmXU1xf+yDcbd3H/xYl0aa1uQhERCb3YRg15enwym3ceZPona50OJ2youzBE6lJG/XHXIaZ9sJoze7XjqsHeLg8pIiLys6pX5EhPTw/YZWWGn9yOX53Wlde+3sz5/eMZ0t17t6H8TJUslyottUyalU2UMTw5Lrmyq1FERKQmwV6r6p7RfUlo3ZS7MrM5VFgc0Oeuj5RkudQ/vt3C95t388cxiXRq1cTpcEREJAwEe62qZuXdhjl7DvHkh2sC+tz1kZIsF9q88yBPfrSGEX3akz6oi9PhiIhImAjFWlVDe7Tl+uHd+Me3P/DNxp0Bf/76REmWy5SUWiZlZhMT1YBpl6ubUEREfBeqtaruPr/s0m53z1qG54i6DWuige8uk/H1Zhb/sIc/XZFCfMvGTocjIiJhJFRrVTWJiWJ6ejLjX/yWx95fzbTLBwT9nOFIlSwX2bDdw1Mfr2VkvzguO6Wz0+GIiIjUaFDXNtx0Zg/eWPgjC9btcDocV1KS5RLFJaXcmZlN05goHr88Sd2EIiLiev9zXm96tm/G5NnL2F9Q5HQ4rqMkyyVe/nIz2Tl7eXhsEh2aq5tQRETcr3F0FNPTU8jfX8Bj81Y7HY7rKMlygXXbDvDs/HVckBTPmOSOTocjIiLis1NOas3NZ/fkzcU5fLFmu9PhuIqSLIcVlZRy51vZxDZuyCOXqptQRETCz+0je9EnrjmT317GvkPqNqygJMthL2ZtZPlP+3j00iTaxTZyOhwREZFaa9SwrNtwp6eQh+atdDoc19ASDg5albef5z9fz5iUTlw4IDDdhMG8bpWIiEhNBnRpyS1pPXn+8w1ckNSR8xLjnA7JcapkOaSwuJS7MrNp2SSGhy/pH7DnDfZ1q0RERGpy6zm96NexBfe+s5w9BwudDsdxSrIc8sIXG1i1dT+PX5ZE62YxdXoOj8dDRkYG06ZNIyMjA4/HE/TrVomIiNQkpmEDnklPYc/BQqa+p25DJVkOWPHTPl74YgOXndKZUf3rfl0pb1WrUFy3SkREpCaJnVpw27m9eHdpHh+t2Op0OI5SkhViR4pLuCszmzbNYpg6xr9uQm9Vq1Bdt0pERKQmv03rSVLnFtz3zgp2eY44HY5jNPA9xP7y2QbW5B/gtetTadk02q/nio+PJycnB2ttZdUqVNetEhERqUl0VAOeSR/ImL98xQPvruSFa051OiRHqJIVQtk5e/n7vzeSPqgL5/T1f9aFqlYiIuJWfeKbc/vIXry/fCvzluU5HY4jVMkKkYKism7CDs0bcf/FiQF5TlWtRETEzX5zVg8+WbWNP85ZwdDubWnfPLLWg1QlK0Se/XQd67d7eGJcMi2b+NdNKCIiEg4aRjXgmfRkDhaWcN87yyvHEUcKVbJCYMkPe3h5wSauHpLA2b3bn3B/LSgqIiL1xckdmnPXqN48/sEa5mbnMXZgZ6dDChlVsoKsoKiESZnZdGzZhHsv7OfTMVpQVERE6pNfn9GDU09qxQPvrmT7/gKnwwmZEyZZxpjXjDHbjTErjto20BjznTFmqTFmsTFmSA3HbjHGLK/YL5CBh4vpH69l086DPDU+meaNfesm1IKiIiLiFt4Wvq6tqAaG6ekpFBSVcG8EdRv6UsmaAYyusu0p4CFr7UDggfL7NRlhrR1orU2tU4RhbNGW3bz69WauHdaV009u5/NxWlBURETcIlC9Kz3ax3L36L58uno7s//zU4CjdKcTJlnW2gXA7qqbgRblt1sCkTk38zgOFRZzV2Y2XVo3YfIFfWt1rJZmEBERtwhk78qE4d0Y0q0ND723kq37DgcqRNeq68D3PwAfG2OmU5aoDa9hPwt8YoyxwP9aa1+q4/nCisfj4Ya/fsgPe5tyXcIebFEBNPJ94LqWZhARkWCqzQQrbwtf11WDBoan05MZ/ecvmTx7OTMmDK7suamPjC/9osaYbsA8a21S+f3ngX9ba2cbY64AJlprR3o5rpO1Ns8Y0wGYD/y+vDLm7RwTgYkAcXFxg2bOnFnX1+QTj8cTtBl732/ew9/XxjC83RHGdD5CTEwMbdu2Dcq56pNgtonUndrFfdQm7hRO7bJr1y4KCwsr7x/ve6q0tJQ9e/ZQVFREdHQ0rVu3pkED/+bNffpDEf9aXciEpBjO7hK8ZY1C1SYjRoxY4m1YVF2TrH1AK2utNWUp6D5rbYsTPMdUwGOtnX6i86WmptrFi4M7Tj4rK4u0tLSAP+/BI8UMm/oeUMrYRquINqXExMQwZcqUgJ+rvglWm4h/1C7uozZxp3Bql2nTplVLskL5PVVaarnmle9Z/tM+Pr7jLDq3ahKU84SqTYwxXpOsuqaiecDZ5bfPAdZ7OWEzY0zzitvAKGBF1f3CXdVZFw/PXc4BG81ZMVuINqUauC4iIq7j9ASrBg0MT41PxlrLPbOW1dvZhr4s4fAG8C3QxxiTa4z5NXAT8IwxJht4nPJuPmNMJ2PMB+WHxgFfle+zEHjfWvtRMF6Ek46edfHdlr28uSSP64YlkNqtjQaui4iIK7lhglVCm6bce1E/vtqwk9e//zHk5w+FEw58t9ZeXcNDg7zsmwdcWH57E5DiV3RhoGLWRaFtwFeFXWnVoIApFyXROLrev3QREQlTTk6wOnrQfVxcPKd178PjH6zm7N7tSWjT1JGYgkUrvvupouS6qCiBQzaGcV0O0Tg6yumwREREXOnoHqDc3ByGmvU0MIZJs7IpLT1+t2EgFkYNJSVZfkpPT+dIm56sK2nP8DaHuONXlzkdkoiIiGtVXXfr8K48/nhxP77btJt/fvdD5X7eEqpwu+ycLhDtp5KoRny2P55eHRry6u9Hq4olIiJSztt6XN7W3boiNYEPV+TzxIdrOLt3e7q1a1aZUFlrKxOqcLvsnCpZfnpk3ip2eI7wzBUpSrBERESO4q3y5G3Q/cGDBzmlaDUlRYVc+8In7Nt/wGtC5fSsyNpSJYvarXx7tM9Wb2PWklx+f87JJHdpFfxARUREwoi3RMnboPuMjAz25ecwNPoAXx7uzqSX3+c0LxWv9PT0at/XbqYkC7yWJE8062LvoUImv72cvvHN+f05vUIUqYiISPjw9ZI8FclYz6hdbClpzWc7WvD78efA5+9XK4CE02XnlGRRt4tfTp27kj0HC8m4fjAxDdXrKiIiUpWvlaeKZAwspzf6gTlHBvDA++uZdfP1RDUI32sbKsmi9he//GhFPnOW5vGHkb1I6twyRFGKiIiEF18rT0cnYyfFx/Nwv/7cPWc1L3+5iZvP7hmCSINDSRa+Z9oAuw8Wcv+c5fTv1IJbRpwcwihFRETqp6rJmLWWz9fv4U+frOPcvh3oFdfcwejqTkkWvmfaHo+Ha5/7mN2eGMa328qRw6cQHSZXXBcREQkXxhgevSyJUc8u4M7MbN7+7XAaRoXf0Jzwi9hBU1+by8oDjRnYcCuF27e4fhE0ERGRcNUuthGPjE1iWe4+/nfBJqfDqRMlWT7a6TnCvK3NaGcOMqDh1rBYBE1ERCScXZTckYuSO/LnT9exJn+/0+HUmpIsH1hruf+dFRQRxZmNttDAhMciaCIiIuHukbFJtGwSzZ1vZVNUUup0OLWiJMsHc7Pz+GhlPreP6EFy1/bHrFIrIiIiwdOmWQyPXjqAlXn7eeGLDU6HUysa+H4C2/cX8MC7KznlpFbcMrIfUQ0SnQ5JREQkooxOiufSgZ346+cbGNkvLmyWT1IlqwYej4fXXssg/al3OFhQyMMX9QrrBdFERETC2dRL+tO6WQx3ZWZTWBwe3YZKsmqQmZnJF1sO8kNxC05tmMviz993OiQREZGI1appDNMuG8Ca/AP85fP1TofjEyVZNdiYt5PvCxPo0OAA/aK2aSahiIiIw0YmxjHu1C78LWsjy3L3Oh3OCSnJ8sJay3elPSmhAWdGbyGqgWYSioiIuMEDYxJpFxvDnW9lc6S4xOlwjivikiyPx0NGRgb5+flkZGTg8Xiq7ZO5OJcNBxtxXoeDtGtsNZNQRETEJVo2ieaJccms3+7h2fnu7jaMuNmFmZmZ5OTk0KtXL3JycsjMzDzmkjo/7T3MI/NWMbR7G/5204U00GB3ERERVxnRpwNXpibw0oKNjOofx6kntXY6JK8irpKVn5+PtRag2qrt1lomz15GibU8PT4l6AlWRVVt2rRpNVbVREREpLr7Lu5HfIvG3JWZTUGRO7sNIy7Jio+Px5iy5Knqqu1vLMzhy/U7uffCfpzUtmnQY6moqhUWFlZW1UREROTEWjSO5snxyWzacZBnPlnrdDheRVySlZ6eTkJCAsaYY8Za5ew+xGPvr+KMk9txzdCTQhLL8apqIiIicnxn9mrPNUNP4pWvNrN4y26nw6km4pKs2NhYJkyYQHx8PBMmTCA2NpbSUss9s5dhjOGJcQMqK13BdryqmoiIiJzYlAv70blVE+7KzOZQYbHT4Rwj4pIsb17//ge+2biL+y7qR5fWwe8mrFBRVdO1EEVEROomtlFDnhqfzJZdh3jqI3d1G0bc7MKqfth1kMc/WMPZvdtz1eCEkJ67oqomIiIidTe8ZzuuO60rM77Zwvn94zmtZ1unQwIivJJVWmqZlLmMhlGh7SYUERGRwLrngr50bduUu2dnc/CIO7oNT5hkGWNeM8ZsN8asOGrbQGPMd8aYpcaYxcaYITUcO9oYs9YYs8EYMzmQgQfCjG+2sHDLbh4c05+OLZs4HY6IiIjUUdOYhjw9PoXcPYd54sM1TocD+FbJmgGMrrLtKeAha+1A4IHy+8cwxkQBLwAXAInA1caYRH+CDaT8g6U89fEazu3bgXGndnY6HBEREfHTkO5tuOH07vzzux/4esNOp8M5cZJlrV0AVJ0XaYEW5bdbAnleDh0CbLDWbrLWFgIzgbF+xBowJaWWV5YfoVHDKKZdrm5CERGR+mLS+X3o0a4Zd89axuFi62gsdR2T9QfgaWNMDjAdmOJln85AzlH3c8u3Oe7VrzaxYW8pD4/tT4cWjZ0OR0RERAKkcXQU069IYeu+w8xcU+hoLHWdXfhb4A5r7WxjzBXAq8DIKvt4Kw/VmFIaYyYCEwHi4uLIysqqY2jHV2ot//qmgJS2lhZ71pGV5e6LS0YSj8cTtHaXulO7uI/axJ3ULu4yuls0X+UWMu+TL4iNcabHylSsOH7cnYzpBsyz1iaV398HtLLWWlPW17bPWtuiyjGnAVOtteeX358CYK2ddqLzpaam2sWLF9f2tfisoKiET79YwMWjRgTtHFJ7WVlZpKWlOR2GVKF2cR+1iTupXdwllN/1xpgl1trUqtvr2l2YB5xdfvscwFs5aBHQyxjT3RgTA1wFzK3j+QKqcXSUY1mtiIiIBJ8bvutP2F1ojHkDSAPaGWNygQeBm4DnjDENgQLKu/mMMZ2AV6y1F1pri40xtwIfA1HAa9balcF5GSIiIhKJPB4PmZmZ5OfnEx8fT3p6OrGxsU6HBfiQZFlrr67hoUFe9s0DLjzq/gfAB3WOTkREROQ4MjMzycnJwVpLTk4OmZmZrrmaSkSv+C4iIiLhLT8/n4rx5dZa8vPzHY7oZ0qyREREJGzFx8dXrndpjCE+Pt7hiH6mJEtERETCVnp6OgkJCcTExJCQkEB6errTIVWq6zpZIiIiIo6LjY11zRisqlTJEhEREQkCJVkiIiIiQaAkS0RERCQIlGSJiIiIBIGSLBEREZEgUJIlIiIiEgRKskRERESCQEmWiIiISBAoyRIREREJAiVZIiIiIkGgJEtEREQkCJRkiYiIiASBkiwRERGRIFCSJSIiIhIESrJEREREgkBJloiIiEgQKMkSERERCQIlWSIiIiJBoCRLREREJAiUZImIiIgEgZIsERERkSBQkiUiIiISBEqyRERERIJASZaIiIhIECjJEhEREQmChifawRjzGnAxsN1am1S+7U2gT/kurYC91tqBXo7dAhwASoBia21qQKIWERERcbkTJlnADOCvwP+r2GCtvbLitjHmGWDfcY4fYa3dWdcARURERMLRCZMsa+0CY0w3b48ZYwxwBXBOgOMSERERCWv+jsk6E9hmrV1fw+MW+MQYs8QYM9HPc4mIiIiEDWOtPfFOZZWseRVjso7a/ndgg7X2mRqO62StzTPGdADmA7+31i6oYd+JwESAuLi4QTNnzqzVC6ktj8dDbGxsUM8htaM2cSe1i/uoTdxJ7eI+oWqTESNGLPE27tyXMVleGWMaApcDg2rax1qbV/7/dmPMO8AQwGuSZa19CXgJIDU11aalpdU1NJ9kZWUR7HNI7ahN3Ent4j5qE3dSu7iP023iT3fhSGCNtTbX24PGmGbGmOYVt4FRwAo/ziciIiISNk6YZBlj3gC+BfoYY3KNMb8uf+gq4I0q+3YyxnxQfjcO+MoYkw0sBN631n4UuNBFRERE3MuX2YVX17D9ei/b8oALy29vAlL8jE9EREQkLGnFdxEREalXPB4PGRkZ5Ofnk5GRgcfjcSQOJVkiIiJSr2RmZpKTk4O1lpycHDIzMx2JQ0mWiIiI1Cv5+flULFFlrSU/P9+ROJRkiYiISL0SHx9P2UVpwBhDfHy8I3EoyRIREZF6JT09nYSEBIwxJCQkkJ6e7kgcdV6MVERERMSNYmNjmTBhQlgvRioiIiIiNVCSJSIiIhIESrJEREREgkBJloiIiEgQKMkSERERCQIlWSIiIiJBoCRLREREJAiUZImIiIgEgZIsERERkSAwFRdQdBNjzA7ghyCfph2wM8jnkNpRm7iT2sV91CbupHZxn1C1SVdrbfuqG12ZZIWCMWaxtTbV6TjkZ2oTd1K7uI/axJ3ULu7jdJuou1BEREQkCJRkiYiIiARBJCdZLzkdgFSjNnEntYv7qE3cSe3iPo62ScSOyRIREREJpkiuZImIiIgETcQlWcaY0caYtcaYDcaYyU7HE6mMMQnGmC+MMauNMSuNMbeXb29jjJlvjFlf/n9rp2ONNMaYKGPMf40x88rvq00cZoxpZYyZZYxZU/6eOU3t4ixjzB3ln10rjDFvGGMaq01CzxjzmjFmuzFmxVHbamwHY8yU8u//tcaY84MdX0QlWcaYKOAF4AIgEbjaGJPobFQRqxi401rbDxgG3FLeFpOBz6y1vYDPyu9LaN0OrD7qvtrEec8BH1lr+wIplLWP2sUhxpjOwG1AqrU2CYgCrkJt4oQZwOgq27y2Q/l3zFVA//Jj/laeFwRNRCVZwBBgg7V2k7W2EJgJjHU4pohkrd1qrf1P+e0DlH1pdKasPf5Rvts/gEsdCTBCGWO6ABcBrxy1WW3iIGNMC+As4FUAa22htXYvahenNQSaGGMaAk2BPNQmIWetXQDsrrK5pnYYC8y01h6x1m4GNlCWFwRNpCVZnYGco+7nlm8TBxljugGnAN8DcdbarVCWiAEdHAwtEv0ZuBsoPWqb2sRZPYAdQEZ5N+4rxphmqF0cY639CZgO/AhsBfZZaz9BbeIWNbVDyHOASEuyjJdtml7pIGNMLDAb+IO1dr/T8UQyY8zFwHZr7RKnY5FjNAROBf5urT0FOIi6oRxVPsZnLNAd6AQ0M8b80tmoxAchzwEiLcnKBRKOut+FshKvOMAYE01ZgvW6tfbt8s3bjDEdyx/vCGx3Kr4IdDpwiTFmC2Vd6ecYY/6F2sRpuUCutfb78vuzKEu61C7OGQlsttbusNYWAW8Dw1GbuEVN7RDyHCDSkqxFQC9jTHdjTAxlA+DmOhxTRDLGGMrGmKy21v7pqIfmAteV374OeDfUsUUqa+0Ua20Xa203yt4bn1trf4naxFHW2nwgxxjTp3zTucAq1C5O+hEYZoxpWv5Zdi5l40rVJu5QUzvMBa4yxjQyxnQHegELgxlIxC1Gaoy5kLJxJ1HAa9bax5yNKDIZY84AvgSW8/P4n3spG5f1FnASZR9k6dbaqoMaJciMMWnAXdbai40xbVGbOMoYM5CyyQgxwCZgAmV/JKtdHGKMeQi4krKZ0v8FbgRiUZuElDHmDSANaAdsAx4E5lBDOxhj7gNuoKzd/mCt/TCo8UVakiUiIiISCpHWXSgiIiISEkqyRERERIJASZaIiIhIECjJEhEREQkCJVkiIiIiQaAkS0RERCQIlGSJiIiIBIGSLBEREZEg+P/rv5d4XrdwrgAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "two_bkpt_truth_params = [20, 20, 80, 20, 0.1, -0.1]\n", "truth_func = TwoBkptSegRegEstimator().get_func_for_params(two_bkpt_truth_params)\n", "\n", "indep = np.arange(101)\n", "seed = 65434789\n", "rng = default_rng(seed)\n", "dep = truth_func(indep) + 0.3 * rng.standard_normal(len(indep))\n", "\n", "plt.scatter(indep, dep, s=15, color=\"gray\");\n", "plt.plot(indep, truth_func(indep));\n", "plt.legend([\"data-generating function\"]);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Estimate All Models" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "ols_estimator = OLSRegressionEstimator()\n", "one_bkpt_estimator = OneBkptSegRegEstimator()\n", "two_bkpt_estimator = TwoBkptSegRegEstimator()\n", "\n", "estimators = [ols_estimator, \n", " one_bkpt_estimator, \n", " two_bkpt_estimator]\n", "\n", "for estimator in estimators:\n", " estimator.fit(indep, dep)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Plot the Fits" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlkAAAEvCAYAAAB2a9QGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAABdkUlEQVR4nO3dd1zU9R/A8dfn2IoiTlQwR87ce2W4V7aUzKyc2dDSLPtltmxowyxTGy60nNEwS9PMIjW3Zu6GpoKKA+ep7M/vj++BjEMOuAW8nz16yH3vO953H7jv+z5Taa0RQgghhBD2ZXJ1AEIIIYQQhZEkWUIIIYQQDiBJlhBCCCGEA0iSJYQQQgjhAJJkCSGEEEI4gCRZQgghhBAO4OnqAKwpW7asrlq1qkOvcfXqVYoXL+7Qa4jckTJxT1Iu7kfKxD1JubgfZ5XJzp07z2mty2Xe7pZJVtWqVdmxY4dDrxEZGUloaKhDryFyR8rEPUm5uB8pE/ck5eJ+nFUmSqlj1rZLc6EQQgghhANIkiWEEEII4QCSZAkhhBBCOIBb9skSQgghhPtJTEwkOjqauLg4V4dik4CAAA4ePGi38/n6+hIcHIyXl5dN+0uSJYQQQgibREdHU6JECapWrYpSytXh5OjKlSuUKFHCLufSWhMbG0t0dDTVqlWz6RhpLhRCCCGETeLi4ihTpkyBSLDsTSlFmTJlclWLJ0mWEEIIIWxWFBOsVLl97ZJkCSGEEKLA8/f3d3UIWUifLCHcmNlsJiIigpiYGIKCgggLC3PLDxIhhBBZSU2WEG4sIiKCqKgoEhISiIqKIiIiwtUhCSGEy02dOpX69etTv359PvzwwwzPnTp1ig4dOtC4cWNatWrFhg0bXBMkUpMlhFuLiYlBaw0YI1tiYmJcHJEQQrjWzp07CQ8PZ+vWrWitadWqFXfccUfa84sXL6Z79+5MmDCBixcv4uHh4bJYJckSwo0FBQURFRWF1hqlFEFBQa4OSQghAJj4/X4OnLxs13PWq1SSV/vcdtN9Nm7cyL333pu28PN9992XobaqRYsWDB06lMTERLp27Uq7du3sGmNuSHOhEG4sLCyMkJAQvL29CQkJISwszNUhCSGES6XW7menQ4cOrF+/nsqVKzNixAg+//xzJ0WWldRkCeHG/P39GTJkiKvDEEKILHKqcXKUDh06MHjwYF544QW01nz77bd88cUXac8fO3aMypUr8+ijjxIbG8uuXbt45JFHXBKrJFlCCCGEKDCaNm3K4MGDadmyJQDDhw+nSZMmac9HRkby3nvv4eXlhZ+fH4sWLXJVqJJkCSGsS50+ws/Pj/DwcJk+QgjhNsaOHcvYsWMzbDObzQAMGjSIQYMGAfZdVicvpE+WEMKq1OkjtNYyfYQQQuSBJFlCCKtk+gghhMgfSbKEEFYFBQWlrdMl00cIIUTuSZIlhLAqdfoIpZRMHyGEEHmQY5KllApRSv2qlDqolNqvlBpt2V5aKbVWKfWP5d/AbI4/qpTaq5TarZTaYe8XIIRwjNTpI4KCghgyZIh0end3ixZB1apgMhn/unBElRDCYEtNVhLwrNa6LtAaGKmUqge8AKzTWtcE1lkeZ6ej1rqx1rp5viMWQgiR0aJFMGIEHDsGHtr4d8QISbSEcLEckyyt9Smt9S7Lz1eAg0Bl4G5ggWW3BcA9DopRCCHEzbz0IoQkwMPFYHwJaO8N167BhAmujkwIuxs6dCjly5enfv36adtCQ0PZsSNrY9mkSZOYMmWKzeeOjIxk06ZNdokTctknSylVFWgCbAUqaK1PgZGIAeWzOUwDPymldiqlRuQjViGEEOldvwC/fwR3nYcHikEZE/yXDJ194R5fOHHc1REKYXeDBw9m9erVDjm3vZMsmycjVUr5A18DY7TWl1NHHdmgndb6pFKqPLBWKXVIa73eyvlHACMAKlSoQGRkpK3nzxOz2ezwa4jckTJxT1Iu7kedPcjJzz6mwulf8UhJIOW6J6a1ZjiUZHyt7eANHX1JKe/L5p++I9E7wNUhFwlF4W8lICCAK1euuDSGJk2acOzYMVJSUtJiSU5OZt68eYwcOZIrV64wc+ZMmjdvjtaa+Ph4rly5wvz581mxYgWLFi2ib9++NGjQgJ07d6btX65cOT755BM8PDz4/PPPee+992jbtm2W68fFxdlczjYlWUopL4wEa5HW+hvL5tNKqYpa61NKqYrAGWvHaq1PWv49o5T6FmgJZEmytNazgFkAzZs316GhoTa9gLyKjIzE0dcQuSNl4p6kXNxESjL8tQq2fgZHN4CnLzTqD60ew7RuDywaATrJ2Hd9AlzxwnS3F+0OvAQPfgnl67o2/iKgKPytHDx40KUzqKfy9/fHZDKlxeLh4UFiYiJbt25l/fr1PPnkk+zbtw+lFD4+PixYsIC1a9fyww8/4OPjk+3+TzzxBP7+/jz33HPZXtvX1zfDMj43k2OSpYwqq7nAQa311HRPrQAGAW9b/v3OyrHFAZPW+orl527A6zZFJoQQAq6dh12fw/a5cOk4BIRwuPoj1Og3EYqVNvYZ2MD4d8IEOH4cqlSBZ9+CO+rA0gEwpyuEhUPNrq57HaLw+fEFiNlr33MGNYCeb+fp0AEDBgDGAtKXL1/m4sWLAHzxxRcEBwezfPlyvLy8ctzfnmzpk9UOeBjoZJmGYbdSqhdGctVVKfUP0NXyGKVUJaXUKsuxFYCNSqk/gW3ASq21YxpShRCiMInZByuegqn14OdXoVQVuP9zeHo3UVX63kiwUg0cCEePQkqK8e/AgRDcDB79BUpXhcX3w5ZPwTKLvxCFTeZuTKmP69evz9GjR4mOjrZpf3vKsSZLa70RyO7Kna3sfxLoZfn5CNAoPwEKIURhl7oY9+lTJ2kZcJYOvofwjNoMnn7Q8H5oOQKC6ud8ImsCgmHIavhmBKz+H5z7C3q+Cx5eOR8rxM3kscbJUZYtW0bHjh3ZuHEjAQEBBAQYfRGbNGnCE088wV133cWaNWuoVKlStvuXKFGCy5cv2y0mmzu+CyGEcIwVS+cTcmI19+rdlDp7hSsegZToMhGaPpK1xiovfPyh/0JY9xr8Pg3OH4Gw+eBndQ5pIdzagAEDiIyM5Ny5cwQHBzNx4kQAAgMDadu2LZcvX2bevHkZjmnfvj1Tpkyhd+/erF27Ntv9+/TpQ79+/fjuu++YPn06t99+e75ilSRLCCFcJWYvbP2MsOgleJHEf4Swmo4c9ajDC+3H5Pm0qTVjMTExBAUFERYWZszY3/V1KFsbvh9t9NN6cBmUqWG/1yOEEyxZsiTLtmHDhlnd98UXX0zrHN+9e3e6d++e9lzfvn2ZPHlyhv1r1arFnj177BarJFlCCOFMyUnw10pjlOCx38HTjyP+Lfj1ak1OU9ZYK7JipXxdIiIigqioKLTWREVFERERwZAhQ4wnmwyEwKqw7CGY0xnu/wKq5e/busgq20RXFCmyQLQQQjjD1VjYMBWmNYIvH4GLUUbN0tgDVH78K3yqNMXb29sui3HHxMSgLR3ctdbExMRk3KFqO3h0HRQvD1/cY4xeFHaVmugmJCSkJbrCfURGRtK8ueNX+pOaLCGEcKRTe2DbZ7D3K0iKg2odoNe7UKsHmDwA8IcbNU12EBQUlFaTpZQiKCgo606lq8PwtRAx2BjFeO5v6DIxLSaRPzkmuqJIkCRLCCHsLTkJDv1gNAke3wRexaDRAGOUYIV6dr2UtWapsLCwLNus8g2AByNg9QuwaTrEHob7Zhsd5UW+2JToikJPkiwhhLCXq7Gwa74xcejlE8bcVt3ehCYPOWwkX3b9r2yuGfPwhN5ToFxt+PF/MK87DFgKpUIcEm9RYXOiKwo1SbKEEMIiz52VT+0xaq32RkByPFS7A3pNgVrdHd78ZrdmqZaPGk2IEYNhdicYsASCHd9npbDy9/e3axOwKJik47twS2azmfDwcCZPnkx4eDhms9nVIYkiIFedlZMTYf+3MK8HfHY77P/GGLn35BYYtALq9HJK/6agoKC0marz3Sx1a2cY/jN4F4PwXkY/MiHcTFRUFB07dqRu3brcdtttTJs2DYDQ0FB27NiRZf/XXnuNKVOm2Hz+yMhINm3aZJdYJckSbklG5ghXsKlW6Oo5WD8FPmxo1PpcOQXd3sL82HbCzzZk8tzlTv1iEBYWRkhIiN1GJlKuNgz/BSo3g6+Hwa+TZCke4VY8PT15//33OXjwIFu2bGHmzJkcOHDAbue3Z5IlzYXCLcnIHOEKN+2sfHK30SS472ujSbB6R7hzKtTsBiYPIsLDs5+byoEc0ixVvAw8shx+eAZ+ewfO/QP3fAxefva9jhB5ULFiRSpWrAhAiRIlqFu3LidOnABg4cKFPP3002mzuNetWzfDsbNnz+abb77hm2++oWfPnjRu3Jht27al7V++fHk+/fRTPDw8WLhwYb5nfZckS7glGZkj7MnWvlZZOivfd4+RVG2dBVFbwKs4NH3YGCVYrnaGYwvdFwNPH7h7JpStBT+/BhePwQOLoYT8LQr3cfToUf744w9atWoFwNWrV9m0aRPr169n6NChbN68OW3fGTNm8NNPP7F8+XJ8fHys7r9v3z4ef/xx/P39ee655/IdnyRZwi0VhJE5MqNzwXHTGdDTSasVMp+FnfNhThujOTCwGnSfbPS58g2weo1C+cVAKWg/BsrWhK+HWzrEL4WKDV0dmXAD72x7h0PnD9n1nHVK1+F/Lf9n075ms5m+ffvy4YcfUrJkScBY1xCgQ4cOXL58mYsXLwLwxRdfEBwczPLly/HyurE4enb724skWcItFYSRObbeuIXr2VzLdPKPdE2CCVCjE/SZBrd2BdPNu7AWhC8GeVanNwxdDUsGGB39+84xOvYL4SKJiYn07duXgQMHct9996VtTx0Ekvlx/fr12b17N9HR0VSrVi3H/e1Fkiwh8qjQNQ8VYjetZUpOhAPfGclV9Dbw9oemgyxNgrVsvoa1LwaFqrazYiN49Bcj0Vr6oLEkUNunjNouUSTZWuNkb1prhg0bRt26dRk7dmyG55YtW0bHjh3ZuHEjAQEBBAQYNc9NmjThiSee4K677mLNmjVUqlQp2/1LlCjB5cuX7RKrJFlC5FGhbB4qpKzWMpnPGE2C2+eCOcaYI6rH29D4wWybBHOr0NV2lgiCwSth+ROw9mUST+1n8eXmnDx9ruAnkaLA+P333/niiy9o0KABjRs3BmDSpEkABAYG0rZt27SO7Om1b9+eKVOm0Lt3b9auXZvt/n369KFfv35899130vFdCFcp1M1DhUyGWqYTu+CnZ415rZIToEZnuOsjm5oEc6tQ1nZ6F4N+4RBZE6/173EHG/mSPkRFJRb8JDKPUmss/fz8CA8Pl2TTwdq3b5/2d5Ver15Zm7CvXLnCa6+9lva4e/fudO/ePe1x3759mTx5coZjatWqxZ49e+wSqyRZQuRRQeg3JiySEuDgCtj6KURvN5oEmw02mgTL1nTYZQttbafJBJ1eYsWmQ/RKWsVwFrM45V5iYrxyPrYQSq2xrFmzZuGosRR2I0mWcLmC8C3QnfrWuFMsbu/KaaNJcMdcMJ+G0jWgxzuWJsGSDr98fmo7C0I5x1buzOfH/blff8dwlvBbwFBXh+QShbLGsgiIjIx0+DUkyRIuVxC+BbpT3xp3isUR7JJcRO80aq32fwspiUZTYKvHjdGCdm4SvJn81HYWhHI2kkhYcKoUD+rldDs7C7bVMNZBzEFBSCJtlVpjCXZY2kgUKpJkCZcrCN8C3SlGd4rFEfKcXCQlwIHlRnJ1Yid4l4DmQy1NgrfmOR5XJQMFoZwzJJHxE+CrYbDqOWOG+O6TwCP7W0x+kkh3S9BSayyVUvZZ2sjNpTZ/F0XW+oLdjKxdKFzOrgvcOog7xehOsThCrpOLKzHw62T4sD588yjEXYKe78HYA9Dr3XwlWOC6dTQLXDn7lIABS6D1SNj2GSy+3yiLbOQniXS3tU1Tk82goCCGDBlSYGvkbOHr60tsbGyuk43CQGtNbGwsvr6+Nh8jNVnC5QrCt0B3GknoTrE4gs2dxaN3WJoElzu0SdBVNUoFspxNHtBjkjG/2MpnYW43Y4b40tWy7JqfQQEFoZavsAoODiY6OpqzZ8+6OhSbxMXF5Sopyomvry/BwcE27y9JlnC51G+BkZGRhIaGujocq9xpJKE7xZKd/DTn3DS5SIo3kqqtn8LJXUaTYIvhRh+gMjUc8lpcNUKwIJSzNWazmYg9Gm+PMPrFfof37E6oBxbDLW0y7JefJLLQjtosALy8vDLMmO7uIiMjadKkicuuL0mWECJfrCVU+elvYzW5uHwKdobDjnlw9SyUqQm9pkCjB4xmKgcqkDVKLnSj7Csym/48nPgDAQv6GHORNX4wbb/8JJFSJqKgkCRLCJEv1hIquzTnaG3MabX1M6NDe0oy1OwGrUZAdeeNEiyoNUqukr7sYwlknmkgz4TsMmaJP/c3dHol32UnZSIKCkmyhBD5Yi2hyldzTlI87PvG6Dx98g/wKQktHnVok6Cwn8xlX6piVXjoJWPU4cYPjJGH980C7+KuDlUIh5MkSxRp7jYUvCCyllDlqTnn8imjOXBnuNEkWLaWpUlwAPhImUD+J+51xu+71bL38II7P4SyteGnCTCvBzy4DEpWcnp8QjiTJFmicFu0CCZMgOPHoUoVeOstGDgw7emCMOGju7N2U7W5OUdriNpm1Fod+M5oEqzVHVo9BtU7QhGdiyc7+Z241xm/79mWvVLQ5kmjNvKroTCrozHlQ+WmTo1PEjnhTDkmWUqpEOBzIAhIAWZpracppUoDy4CqwFHgfq31BSvH9wCmAR7AHK3123aLXoibWbQIRoyAa9eMx8eOGY8hLdGSoeCG/Nx48tQ/JjHOWKB562dwajf4BBjTL7QYBqWr5/4FFBH5/X11i9/3Wt1h2E+w+AEI7wX3fgq33eO0+OSLlXAmW2qykoBntda7lFIlgJ1KqbXAYGCd1vptpdQLwAvA/9IfqJTyAGYCXYFoYLtSaoXW+oA9X4QQVk2YQHzidZZ1K8Pl4h6YNKgUjWnNREwNrqGUIiooikuXLoEGhSKwVCCLDy7GQ3mglMKkTJiUCUW6n5XCRLqflQkTJqv7K6WIj4tn/W/ruXj+ImVKl6FTp04UL1Y8yzlye42bXdfaOW7GaTeeyydh+1xjPcFr56BcHeg9FRr2lyZBG+R3+Ra3mfqgwm3w6C+w9EGIGASxL8HtzzklPrdINEWRkWOSpbU+BZyy/HxFKXUQqAzcDYRadlsARJIpyQJaAv9qrY8AKKWWWo6TJEs4nD5+nJdHVObHNqVQKRptSpdo/PGR8a8JCMx4XOS2SPsHYwLKGj8ui1xm//PbEsJNkrJ44lHBCoXxHpm0iYURCzEpEwlxCfh97Wdbwpk5MbRc0xR/BXX5BKarZ1FaYwoqj6lUC1Txcpiu7MG0ad9NE9bk5GSO/neUuGtxFC9WnFtr3IqPt0+213V0wprT/raew9ZrpCb9Xe7qwvcrvieOOMqHlKfXPb24lnjN5sTaraY+8C8Hg76HFU/BL2/C2b8Ju3cyEd9+79D4HJbILVpE62efhTNnrHZNEEWTys3U+EqpqsB6oD5wXGtdKt1zF7TWgZn27wf00FoPtzx+GGiltR5l5dwjgBEAFSpUaLZ06dJcv5jcMJvN0g7vZuxZJikpKeyZO4y53fwZHRHD8JXn0IBWcKVcWVZ8MBUvby9Kly6NRpNCClrrGz+j0dr4GTC2a8v29P9lOubixYskJCakPe/p5UliUmKGY1FQKrBUhmtkPl/q9YAMz2fZls3xaXHrjK8n83+pcV2Lu0ZycnLasSaTCR9fH7TWJCQl4OHpYfUcWd6TDNdLxivhIl7x51ApcaQoD+K9ShHvHUCyydOIUVt/fZnPn5ySnOVayqSsllnm/4oyZfkvNVnL7rFCpSVo6R+n/y81KbT2fPrzpZ0jF+f3vxZNCfNRkr1KcqlUA7TJx+r1bhpfNq8p83NoiLsWR0pKCp4enkatcrqENv1/aLh29RopySl4enriX9wfD5NHlvckcNcf3PJlBH7X4ml4+DomDck+Pvz13HOc6dLF1b8GRZqz7vUdO3bcqbVunnm7zR3flVL+wNfAGK31ZRsXh7S2k9VPPa31LGAWQPPmzbWjZ/5259nFiyp7lsnz4c/zYzd/7tp4iWErzwHGL2Oipxc/t+/AkX+O4O3tzQPjH7DL9VJNnjyZhISEtMfe3t5ZvjmHhIQwpKd79QG5WZ+sXJfLpROwI7VJMNZoEmw5xpg4NI/D9q29r+PHj8/xOK3TJYKWBDBFZ/w//fb0+6foFJJ1cloimPp8sk6+cT6dYjyXcuOYDOdOTRStXCM1sU/WyTeSxXTHJackZ9yW7tr//PMP1WtUN/azcnyWa6c/r2X/tNcHVq9r7f3KsK8lntTXYcv7k/pc5vf+QsnipPiEkHz1LJh3k1y8DCkmj7TzWbtulmuki9FmSUB8LvaNy+a5YGCssdTKg2tjGb/oFB7x8dRbuJB6b76Zu5iEXbn6Xm9TkqWU8sJIsBZprb+xbD6tlKqotT6llKoInLFyaDQQku5xMHAyPwELkZNdp3exRq2h3PVyNLxQk0sBlwm4dImrZcqw5o472NewocP6e9htOgMny/fkjlrD8S2WUYIrQKdA7V7GKMFqHfI9SjCvTTypTW2FTeTpSEJvC3V1GPZ3YhcsGQAx/0C/eUYn+VxKn3CmJYI6axKcXRKaPvn75LNPSExMTKsR9fTyZPDQwVkT5dvbo9F81z6QJZ1Lc+emizT477oxqlkUabaMLlTAXOCg1npquqdWAIOAty3/fmfl8O1ATaVUNeAE8ADwoJX9hLCLqMtRjP51NAEE0PZcWw418OKvho3TFp6+HBGBtwOTnXxNZ1AQJcbBvq+MtQRj9oJvALR+wpg4NLCq3S5TEBJVYQeVmxod4pc8YPzf7U1o/WSukvS0JkJlzCrvhVeew6ldtnbGWuiKIdQrUy/rjvFl4dgx6hyPY1N9fyYOqcSSiYfxCq6S52uLwsGWmqx2wMPAXqXUbsu2FzGSqy+VUsOA40AYgFKqEsZUDb201klKqVHAGowpHOZprffb+TUIAcCl+Es8ue5JNJpZPWax5cctTk92CnVCld6laGOU4K4FlibBusZkkw3vz7ZJ0FqTJGDT1BFF5n0VEFAZhq6Gb0bAmheNpXh6TTEmNM2l/M6JZXNy/9ZbMGIE/teu8eLCk4x5+hYW9a7I4PvfynXMonCxZXThRqz3rQLobGX/k0CvdI9XAavyGqAQtkhMTuTZyGeJNkczp9sc6laoS90hdV0dVuGiNRzbZDQJHvwB0DeaBKvenmNtg7VpIgCZs0hk5V0c7v8CfnndWIrn/BEIWwDFSmd7iL0XKodcJPeWUYRxzz5L5z/O0PFgIh/fF0TXu0OpbPPVRGEkM76LAk9rzZtb32RrzFYmtZ9EswrNXB1SoWJKjoddn8PWWXB6L/iWgrajoPkwCLzF5vNkNz+RzFnkftxiVnSTCbq8BmVqwvejYU4XePBLKHur1d0dtlC5rQYOZEvlyoSGhvLi1RjuXn43b215i5mdZ+Y4T50ovJyzjL0QDhS+P5xv/vmGxxo+Rp8affJ9PrPZTHh4OJMnTyY8PByz2WyHKAugi1Gw9lXabB5mzGWkU6DPRzD2IHR9PVcJFhid11NvNqmd161tE66XmrAkJCRkqHV0iSYDYdAKuH4B5nSGI79Z3S27hcpd8fsVVDyIp5o8xYYTG1hzbI1TrinckyRZokD7+djPfLDzA3pW7cnIxiPtck63usE4m9ZwdCMsewimNYRNH3GxVH0YvBKe+B2aDQLvYnk6dVhYGCEhIXh7e6cNRLC2Tbie282Kfktbo0N8iSBYeJ8xPUgm1hIqV/5+DagzgHpl6vHOtne4nHDZadcV7kWaC0WBte/cPsZvGE/Dcg15o/0bdquSd7sbjDMkXIO9EcZagmf2g18gtH0aWgxj/+4jhFZtn+9LZNe/RfpguR+3WX4nvdLVjDUPvxpqNB+e/Ru6vQEmY4oOdxvZ62Hy4NU2rzJg5QCm7ZzGy21edkkcwrUkyRIF0inzKZ765SnK+JVhWsdp+Hj42O3cbnmDcZSLUbB9jjFK8PoFqFAf7poO9fulq7E64tIQc8st+hMVcG47ZYZvAAxYZow63DITzh+GvnPAp4RbjkCtV6YeA+sO5IsDX9CnRh8al2/s6pCEk0mSJQqcq4lXGfXLKOKS4pjddTZl/cra9fxue4Oxl9QmwW2fwaGVxrY6d0Krx41mmQLSSTe7ZMppi10XYu6YsKTx8IRe70K5WrDqeZjbHR5cCqXcc06qUY1HsfbYWiZunsiXfb7Ey5T3ebtEwSNJlihQklKSGPfbOA5fPMzHnT/m1kDrI43yw61vMPmRcA32LINtszM1CQ6HUiE5H+9mskumimRzb1HUYjiUrg5fDobZneCBJRDSwtVRZVHMqxgTWk3gqV+eYsH+BQxvMNzVIQknkiRLFChTdkxhw4kNvNTqJdpWbgtI81COLhyzNAl+DnEXoUIDuGsGNOgHXn6uji7PskumilRzb1FXoxMM/xkW3w/ze8M9Hxu/124mNCSULlW68Omfn9L9lu6ElCx4X2pE3sjoQlFgLD64mEUHF/FwvYfpX6d/2vYiPRowO1obQ92XDoSPGsPmmVA9FIb8CI9vgKYPF+gEC6yPJgProxhFIVauljHyMLg5fD0Mfp0EKblcKNoJXmj5Ap4mT97c+mbalwNR+ElNligQNkRv4J3t7xAaHMqzzZ7N8Jw0D6WTcBX2fAnbZsGZA1CsDLR/BpoPhYBgV0dnV9n1nSu0zb1FlE011cVKw8PLYeUz8Ns7xlI8d3+c5+lGHKFC8Qo83eRpJm+bzI///Uiv6r2y3Vdq5wsPSbKE2/v7wt+MWz+OWoG1eKfDO3hYhmynKizNQ/n6YL1wNF2T4CUIamjcZOr3BS9fh8btKpJMFQ02D2Tw9DaawcvWhrWvGM3kA5YYc2vZUX7+TvvX7s8PR37gne3v0K5yOwJ8AqzuJ4M3Cg9pLhRu7dz1c4xaN4rinsWZ3mk6xbyyfjMtLM1DuW721BqORMKSB2FaY9j8sdFHZegaeGy9MVN2IU2wRNGRq5pqpaDd0/DAYjj7l9Eh/tSfdo0nP90TUufOuhR/iQ92fpDtflI7X3hITZZwW9eTrvP0L09zMf4i83vMJ6i49W+khaVGw+YP1oSr8OdSo0nw7CGjSfD2sVyt9wBfrv6dmEW/EhR0UJoYRKGQp5rqOr1g6GpYMgDm9YD7ZkPdO296iK01VPlNgGqXrs0j9R4hfH84fWr0sbrWamGpnRdSkyXcVIpOYcLGCew7t4+3b3+bemXquTokh8txnbXz/8GaCTC1LqwcC54+RpPgMweg8yt8ufp3GQAgCp0811RXbGh0iC9f11gmauMHRu1vNmytobLHeoiPN3qcSsUr8frm10lMTszyfGGpnRdSkyXc1PQ/prP22Fqea/4cnap0cnU4TmG1I3dqk+DWz+Dv1cYSInXvMiYODWmZYeJQZzQxSIfcostVZZ+vmuoSFYx1N5c/CT+/Buf+gTs/NPpvZWLr309+JitO/x62DGrJctNywveHM6LhiAz7FZbaeSFJlnBD3/7zLXP2ziGsVhiP1HvE1eE4RHY3rLQP1ngz7FkKW2fBub+gWFno8JwxSrBkJavndEYTg3TIzb+CmqgW2LL38oN+86BsLfjtbaNGuP9CKF4mw262/v3kJwFK/x56RXlRq0otPvvzM7pX7c4tJW/J0zmFe5PmQuFWtp3axuubX6dNxTaMbzXebos+u5tsmybOH4HVL8LUerDyWeMGcc+n8Mx+6PRStgkWOKeJQTrk5l9BndetQJe9UtBxPPSdCyd2wpxORsf4dFzx93Pb6dvw8fDhjS1vyNxZhZTUZAm3cTrxNB9FfsQtJW9hSuiUQr3GV4YP25QU/E5uhsWr4e81RpNgvbuh5WNZmgRvxhlNDNIhN/8KarJSKMq+QT8odQssfRDmdIWwcLi1M+Cav59q5asxptkY3tjyBj8c+YE+Nfo49PrC+aQmS7iFC3EX+PTMp3iaPJnReQYlvUu6OiSHCgoKwptEmuvdPMkCHkhcZnzD7jAOxuwzmjeqtHK7xZqlQ27+2aPjtCsUmrIPaWF0iC8VAovCjCZ5J7H2Hvar1Y9G5Rrx3vb3uBh30WmxCOeQmizhcgnJCYz5dQwXky4S3jWc4BKFa2byLGIPM7DMXji+CG/iOOsdQlznV/FtNsAYMejGpEOuIT/9qvLTcdqVClXZlwoxpnj4ejj8OM6YIb7H2+Dh2Ftidu/hK21eof/3/Xl/5/u80e4Nh8YgnEuSLOFSWmte2/Qau87sYnDZwTQu39jVITlGSgoc+cUYJfjPWrxNHlD/Xmj5GOWCm7tdjZW4ufx0Ane3ZMUZHfHdsrO/Twlj0tKfX4VN0+H8YQibD77WZ2F3pFqBtRh02yDm7pvLXTXuokVQC6fHIBxDmguFS83aM4vvj3zPqMajaFY866R8BV78FaM5YmYLWNgXTu6GO543mgT7zjGaLiTBKnAKar8qa5zREd9tO/ubPKDbm3DXdPhvvdFP6/wRl4TyWKPHqOxfmdc3v05CcoJLYhD2J0mWcJkf//uRGbtn0Kd6nyzzxBR4sYfhx//B+3WN5gjfAGPW6Wf2QccXoWRFV0dY5JnNZsLDw5k8eTLh4eGYzWabjy2o/aqscUbC6PZJadNHjAWmr56B2Z3h2Canh+Dn6cfLrV/m6OWjzN071+nXF44hSZZwid1ndvPSxpdoWr4pr7V9rXBM1ZCSAv+shYX9YHpT2D4XaveE4euMjrYN73d6n6v8JBKFXX5qVwpNJ3CckzAWiKS02u3G32qxMrDgLvhjkdNDaFe5HT2r9WT23tkcueSaGjVhX5JkCaeLvhLN6F9HU6F4BT7s+CHeHllnXy5Q4i4bfa1mNIdF/SBmD4SON+a26jsbgpu7LDS3baZxA/mpXUntVzV+/HiGDBni+v5F+eCMhLHAJKVlasDwtVC1HXz3JKx91fjy5ETPt3geX09f3tgsc2cVBtLxXTjV5YTLjFw3ksSURGZ2nkmgb6CrQ8q7c/8aizTvXgQJZghuYSRX9e62umyHK7h9M40LFYp5n+zAGR3x3a2z/035BcLAr+DH5+H3DyH2X7j3M/BxTiJd1q8sY5uNZeLmiXx3+DvuufUep1xXOIYkWUWAu4zsSUxJ5LnI5zh++Tifdf2MagHVnB5DvqWkwOF1sPVT+PdnMHlB/b7QagRUdr+O+5JIZK+gTqUgnMDDC3pPhbK1Yc14CO8BA5ZBQGWnXP6+mvex4vAKpuyYQofgDpT2Le2U6wr7kySrCHCHNce01kzeOpnNpzbzetvXaVmxZZ7O47KEMe4y7F5s1FydPwz+QdBxAjQbDP7lHX/9PJJEInsFqnbFgbTWaA0pWqMh7WcsP2uyPp96jMay3bIfWbZBSkrGc6WkHp/+XOmum/6aZD6XTt2eep4czkX6WFNfR/r9bpw7c3wAKT53Ur5lAC13PkvSxx3Y2GI6FwIaZDyX5Ri05q9jiRzbdBStLech3bl0xmNSX5+2sl1rzS0pj/Bn/DiGfPcKbQNGZXiPrb0n1l+78V5pK8dYK6ss74u1cifrtUh/3QyvI93P1q6Vad8s58pSRvrGe56S7v3LUn43fkdITmRXaP7+RvIjxyRLKTUPuBM4o7Wub9nWCPgU8AeOAgO11petHHsUuAIkA0laa9d1TinC3KHJ6IsDXxDxdwTD6g/j3pr35vk8uU0YM/+Bpv9jzPohYvk35cbz6vy/+Oyag8++pZgSr5JQsTnXen7C9Zq90SZvUhI1+vw1y7UynislwwdcpptV2gestRtPxhtU+pixdi7jMtl+WAY170EFjOd/PXwZzeVsP+AynMty3b+OJxK95Vg2N8OcbtK2fcCmpHt9mV93hnNleZ8znivtxpvuXOlveDduLsaHdPobQOabSEq2sd44JuN7pS3XME6W4TXf7FxZrp3pJpX5XEB8fDxeG9ZmvG6K9fK7kQRYv67ISUlqqleZ5/Uet294hGcTn2BlSuvsdz+43y5XVQp8ynXgSJlfOPRPLfT1mqDApEChUApMSqEs+yqVdRtp+2U8hrRj0u9v/Eu6bSbLObE8b7Ick3qutOMyx5J6LRMoTDeet5wr9ef011XptqU9n/5aqdfPsG/68xiPTSbjRaTGeubUSbuUR16pnDrWKaU6AGbg83RJ1nbgOa31b0qpoUA1rfXLVo49CjTXWp/LTVDNmzfXO3bsyM0hubJ6Xwxb/9hLzVq1M3zIom/c7LL7FpL1QzDjB3t2H5aZv6WkP1f6m0hKpmNI/6GY5UaR8YPd2rcQjeboseNcv3497QPV18+PysHBWb45WL/xWL9JWruJZPctJM57L5dLzMY7viHFLw8BbcqQaKT+HBcfj5e3d4ZvIZm/pcbFxWO5X1goPDw9s401LxQp3GH6kyEea7jDYw/x2pMfUtowP6k7e3X1vJ1UZJD+Azb1BkDmD1jLz2T6YM1yYyHdB6xlW9o10m3PfK70Nx7Sb0t3DTJ9sGcba6abQWpMN7vhGTcg6687843HlOn8p06dpHKlSlZfd8bzpH/dVs6V9roy3ays3MRuvNb07zmYTDfeq8yv2+p7nqv3Sll5HVnLL0tyYMp4rtT3J+P7mPXn1PjSH5N6Hs/rsZRdORSfU9u53OZ5rrZ6BqVMGd7HzZs20a5d27TjUmNOf66MiUv6csn4ew0QlxRH3xV9UUrx9V1f4+Ph3itCuKPIyEhCQ0Mdfh2l1E5rFUk51mRprdcrpapm2lwbWG/5eS2wBsiSZLmrD3/+m0MxCbB/r8OuYe0mkuUPm+w+hG7yAZtu+83+QNN/MOlipYlPvEhSYhLeXp54lQzk/NWErB+W6T5Yb3wwmDJ844Gs588QU/rtCsz6GH8mz6eEqkqTEqPwDPBJdwPL+GF8KuYUlSuV58YHXNbXffDgAcxXrgAaBZQoUYKGDRpY/xZkKYis71/69/bGa/FOMlP71PfUi15KwPUorvqUY1fwk/wb0o94nzL0V4oHrNysrN1E0r9X2ZZfuvcy9Yab3Y0n403yxvWt3Yys/a5ZvYmmv36mD/70vwubN2+iXdu2Vt8zhUKZcr5+hnJJvZDIs8jIWEJDG7g6jCKmOAxbCSuepuTmdylp/g/umgFevml7lPRRlPG3XyLk6+nLy21e5tGfHmX2ntmMajLKbucWzpHXPln7gLuA74AwICSb/TTwk1JKA59prZ23EudNfD6sJRt/30TbNm0zfAtK/403a7Jh/SaSuapVbiI3nL56mgdXPk95U2kW95pLuWLlbrp/ZOR5QkMb3nQfc6dbLH2MTlv6GPXKf5+ss39ZRgkugcSrENIKWk2keN27aOrhRdP8nb3AC/Q1Ub6kb847ClHYefrAvZ9C2Zrwyxtw4aixNI8D+2W2rtiaPtX7MHffXHpV60X1UlKbXpDk2FwIYKnJ+iFdc2Ed4COgDLACeFprXcbKcZW01ieVUuUxarye0lqvz7yfZd8RwAiAChUqNFu6dGneXpGNzGZzgZ7bxt3Fp8Tz4ekPOZt4lmeCnqGytzEqJyUlhQsXLpCYmIiXlxeBgYGYjEZ055aJTqFM7E4qn/iB0hd2k6I8OVO+A9HBvTGXuNU5MRQQ8rfifqRMXK/s2U3UPfgBiV4B7G3wElf9qzqsXK4kX+HNk29S0asiT1d4GpOSKS5t5ay/lY4dO+atudAarfUhoBuAUqoW0Dub/U5a/j2jlPoWaMmNZsbM+84CZoHRJ8vRbajOaqctipJTkhkTOYaTiSeZ3mk6HYI7pD0XHh6eYUqBkJCQtI7rTimT6xeNea22zYYL/0GJStDpJUxNBxPkX46iMMFBbkdoyt+K+5EycQehcLInHksG0GLPBOg3j8iT/g4rl6R/knhl0ytcDL7IfTXvc8g1CiNX/63kKR221EyhlDIBL2GMNMy8T3GlVInUnzGSsn15D1UUFB/s/IDIqEj+1+J/GRIscOFIxzOH4IexMLUerHkR/CtAv3AYswc6jAP/mzdlFiYyC7wQdlKpibFkVplbYckDBEd957Ahm/fceg/NKzTn/R3vE3s91iHXEPaXY5KllFoCbAZqK6WilVLDgAFKqb+BQ8BJINyybyWl1CrLoRWAjUqpP4FtwEqt9WpHvAjhPiL+jmDBgQU8WOdBHqz7YJbnnbqGWUoyHFoFn98NH7eCPxbCbffCY+th2Bqof58x6WAR4w5TeghRaJSsBEN+hDp3cuvhefD9aEhKsPtllFK83OZlridd593t79r9/MIxbBldOCCbp6ZZ2fck0Mvy8xGgUb6iEwXKppObeGvLW9xe+XbGtRhndR+nTI55/YKRUG2bDRePWZoEXzYmDi1e1v7XK2BkFngh7My7GIQt4Nj8R7ll1wI4fwTu/xyK2W+mdrPZzG/f/Ebtq7VZ9d8qulXuRucane12fuEYMuO7sIvDFw/zbOSzVC9VnffueA9Pk/VfLYfOsn3mEGz7DP5cConXoEob6DoR6txZJGussiOzwAvhACYT/1V/mFuadYUVT8GcLvDgMmMkoh2kNvPXpCZHfY/y6oZXaXtLW/w8/exyfuEYkmSJfIu9HsvIdSPx9fRlZqeZFPcq7ryLpyTD36th62fw32/g4QMNwoy1BCtKRao1spyMEA7U6AEIrApLB8KczkaNVvXQfJ82tZnfAw+axjZlfdB6Zu2Zxeimo/N9buE4kmSJfIlLiuPpX58m9nos4T3Cqehf0TkXTmsSnAUXj0PJytD5FWg6SJoEhRCuVaW10SF+cX9Y2Bd6TYHm1r/Y2DraN30zf4X4Ctymb2P+vvn0rNaTWoG1HP2KRB7JZBtFlNlsJjw8nMmTJxMeHo7ZbM71fik6hZd/f5k9Z/cw+fbJ1C9b3/GBnz5gdCydWg9+egkCQiBsAYzeA7c/KwmWEMI9BN4Cw36C6h3hhzGwerxR856JraN9w8LCCAkJwdvbm5CQEKb0mYK/tz+vb36dFJ3i4Bcj8kpqsoooWxdavtl+H+/+mNVHV/NMs2doXaY14eHhNs+9lCspyfDXj0Z/q//Wg6cvNOgHrR6HIFlaRAjhpnxLwoClxhfCLR9D7GHoO8fYbmHraF9rzfzjWoxjwsYJfPX3V9xf+/5sw8jt3HjCfqQmq4iy9Q87u/1WHF7BZ3s+476a9zHktiGOmXvp2nn4fRp81BiWDYTYI9DlNRh7EO6eKQmWEML9eXhCz7eh91T492eY1x0uHEt7Oj/T2vSp3odWQa34cOeHnL12Ntv9ZG4815Ekq4iy9Q/b2n47Ynbw6qZXaRXUipdavYRSyr5zL53eDyueNpoE174CAVWMzqOj/4T2z9h1WLQQQjhFi2Hw0Ndw+QTM7gTHtwJZmwFzM9pXKcVLrV8iPjmed7a/k+1+Mjee60iSVUTZ+oedeb/WPVszJnIMwf7BvB/6Pl6WqRHyPclochJlz26G+XfCJ21hzzJoGAaP/w5DVkK9u41vhEIIUVDV6AjD1xnNhQvuhD1fpjUDjh8/niFDhuS6Ga9qQFUebfgoa46uYX201VXrnDsJtMhA7lpFlK3D+NPvdyn+Eg+tegiFYmbnmQT4BKTtl+e5l66dh12fw/Y51L8UZXRk7zIRmj4iNVZCiMKnbE0j0fryEfjmUTj3N4S+CKa813kMrT+UVf+t4q0tb9H87uYU8yqW4XmZG891JMkSNklMTuSZyGc4YT7B7G6zqVKySobncz330un9xtxWe76EpOtQ9Xb2BQ+k/n3jpMZKCFG4FSsND30DK8fC+veMROueT42Z43MhfYf2FkEt+NL0JZ/++Sljm4/NsJ/Mjec6cjcTOdJa8/qW19kes51J7SfRrEKzvJ0oOQn+WmUkV8c2gqcfNLwfWj0GFW7jXGSkJFhCiKLB0xvumg7lasNPLxvz/T2wBEraPtdg+tHfKkpRv0p9Pj/wOb2r96Z26doODF7YSvpkiRzN3TeX5f8u57GGj9GnRp/cn+Daedj4AUxrBF8+bHyYdH0dxh6Auz6CCrfZP2ghhHB3SkHbp2DAEjj7t9Eh/uRumw/P3KG9dkxtAnwCmLh5IslW5uQSzidJlripn47+xLRd0+hZrScjG4/M3cGn9sB3I2FqXfj5NShTHfovgtG7od1o6XMlhBAAtXvCsDWgTBDeEw5+b9NhmTu031LhFsa1GMfec3v58u8vHRmxsJEkWSJbe8/u5cWNL9KoXCPeaPdG2h/zTSUnwf7lEN4LPrsd9n4NjQbAE5th0PdQ904weTg8diGEKFCCGhhL8ZSvB8segg1TwVJLlR1ro8R7V+tNm4ptmLZrGmeunXFS8CI70gFGWHXKfIqnfnmKsn5l+ajTR/h4+Nz8gKuxsGs+bJ8Hl6OhVBXo9iY0eQj8Ap0SsxBCFGglKsDgH4wWgHUTjQ7xfaaBp/XP3+w6tL/c+mXuXXEvb297m6mhUx0dtbgJSbJEFuYEMyN/GUlCcgJzu8+ltO9NmvVO/QlbZ8HeCEiOh2p3QK/3oFZ3qbESQojc8vKDvnOhbC2InAwXjkL/hblalzWkZAiPN3qcabumERkVSWhIqKOiFTmQJEtkkJSSxLj14zhy8QifdPmEGqVqZN0pOREO/WCMEjy+GbyKQZOB0HIElK/r/KCFEKIwUQpCXzDm1Fr+pNEh/sEvoXwdm08xqN4gVh5ZyaStk2gZ1DLL3FnCOaRPlsjg3e3vsvHERia0nkCbSm0yPnn1HKyfAh82hIjBcOUUdHvLGCV45weSYAkhhD3V7wuDV0FSHMztCv/8bPOhXh5evNrmVU5dPcXM3TMdGKS4GanJEmkWHVzEkkNLGHzbYMJqpZsR+ORu2DYL9n5lNAlW7wh3ToWa3aRJUAghHCm4mdEhfvEDsDgMerwDrUbYdGjj8o0JqxXGwoML6V29N/XK1LNbWOknQk2dRT63SwIVBVKTJQBYH72ed7e/S6eQToxpOsZoEtz3NcztDrPuMEYMNnkIntwKjyw3hhznkGCZzWbCw8OZPHky4eHhmM1mZ7wUIYQoXAKCYehqqNUTfhwHK581PqNtMKbZGAJ9Au0+d1bqRKgJCQlERUURERFht3MXJpJkCf46/xfjfhtH7cDaTG76LB4bpsKHDeCroWCOge6TLE2CU3PVJ0D+CIUQwk58/I0O8O1Gw/Y5sCgMrl/M8bCS3iV5oeULHIg9wJJDS+wWTuaJUGNiYux27sJEkqwi7uy1s4z6ZRT+Hj7MSChOsenN4Nc3oVwdo6PlU39Am5HgVyrX55Y/QiGEsCOTyVgt4+6ZcHSj0U8r9nCOh3Wv2p32ldsz/Y/pxFy1z+dw5olQg4KC7HLewkaSrCLsevxlnlr5EJfMp5l55ADlD62BpoNg5HajSbBW93ytDC9/hEII4QBNHjI+o6+ehTmdjYTrJpRSTGg1gRSdwuStk+0SgrWJUEVW0vG9KDKfIWVHOC/+/TkHvBUfXfWgTqc3oPGD4Btgt8uEhYVl6RgphBDCDqq2t3SI7w+f3wN9PjSSr2wElwjmicZP8MHOD1h3fB2dq3TO1+WzmwhVZCRJVlFyYpcxt9X+b5hW0o+fSwXwv6p3E3r762AyGaNFwsPtNlpE/giFEMKBSleHYWuNKXW+Gwln/4Iur2U7KOnheg+nzZ3VumJrinsVd2q4RZE0FxZ2SQnG1AtzusDsjnDoB76p14V5pQLoX7s/Azu8kdYkKB3VhRCigPErBQO/ghbDYdNHxrqH8dZHcnuZvHilzSucvXaW6X9Md26cRZQkWYXVldMQ+Q58WB++HgbXYqHH22wduJA3rh6iXaV2vNDyhQyLPktHdSGEKIA8PKH3+9DzPfh7NczrARejrO7aqFwj7q99P4sPLmbfuX1ODrTokSSrsIneCd+MgA9ug8hJENTQ+JYzaidH6nbnmU0TqBpQlffueA9PU8bWYumoLoQQBVirEfBgBFw8ZizFE73D6m6jm46mrF9ZJm6eSFJKkpODLFokySoMkhJgz5cwuzPM6QSHVkHzoTBqJzz0FdTsyoWES4z8eSReJi9mdJ5BCe8SWU4jo0WEEKKAq9nF6Kfl5QfzexvdRTIp4V2CF1q+wKHzh1h0cJELgiw6cuz4rpSaB9wJnNFa17dsawR8CvgDR4GBWuvLVo7tAUwDPIA5Wuu37Re64Mpp2BkOO+aB+TSUudWoLm70APiWTNstITmB0b+O5uz1s8ztPpfK/pWtnk46qgshRCFQvo4x8nDZQ0Z3kdh/4Y7/GQtPW3S9pSt3BN/BzN0z6XpLVyr5V5KlchzAltGF84EZwOfpts0BntNa/6aUGgqMA15Of5BSygOYCXQFooHtSqkVWusD9gi8qEn/y18/MI4egcfx+ut7SEmEW7tCq8ehRqcs81pprXll0yv8ceYP3rvjPRqVa2TXWOQPUQgh3FDxsvDId/D9aIicDOf+NiYx9fIDjC4hL7Z6kXu+u4e3tr7FjE4z0gY/aa3TBj/JF+/8ybG5UGu9HjifaXNtYL3l57VAXyuHtgT+1Vof0VonAEuBu/MRa5H29ZdLCDj+Ew/Hz6dPzIfoQyuN0SRP7bI0CXaxOnHop3s+ZeWRlTzd5Gl6VO1hl1hkFKIQQhQAnj5wzyfQ+VVjLdr5vY0WEItK/pUY2Xgk66PX8/Pxn2XwkwPktU/WPuAuy89hQIiVfSoD6Yc3RFu2idy4EgO/Tqbv8Ve5T6/Cl3hW0YnpXk9Cz7ehTI1sD115ZCUf7/6Yu2rcxfAGw+0WkvwhCiFEAaEU3D7WWPfwzEGjQ3zM3rSnB9YdSJ3SdZi8dTKBQYEy+MnOVOrN8qY7KVUV+CFdn6w6wEdAGWAF8LTWukymY8KA7lrr4ZbHDwMttdZPZXONEcAIgAoVKjRbunRpXl+TTcxms/s2cWlNyct/U/nED5Q7uwmlk4kp0YBDAaHEFKsDyoS3tzdlypTJ9hRH4o4w/fR0qvpUZWSFkXgq+807GxsbS0JCQtrjnGKxlVuXSREm5eJ+pEzck7uXi/+VwzTY+xaeSVc5UO9ZYsu2BOBY/DHej3mf9v7t6UxnEhMT8fLyIjAwEFM+llZzB84qk44dO+7UWjfPvD1Pd16t9SGgG4BSqhbQ28pu0WSs4QoGTt7knLOAWQDNmzfXoaGheQnNZpGRkTj6GrmWFA/7v4Wtn8LJP8CnpDEkt8VwSvhU4GhEBDHHjhMUFESfPn2y/cWJuhLFKytfoXLJyoT3DKeUbym7hpm5T9bNYskNtywTIeXihqRM3JP7l0so3NELljxAg32ToNsb0GYUKMWpbadYfHAxj/d6nIblGro6ULtxdZnkKclSSpXXWp9RSpmAlzBGGma2HaiplKoGnAAeAB7Mc6SF2eVTxgjBneHGgp9la0GvKdBoAPgYyYs/2NQB8XLCZUauG0kKKczoNMPuCRbIKEQhhCiwSlaEIT/C8sfhp5eMpXh6T2VU41GsPbaWiZsnsvTOpXiZvFwdaaFgyxQOS4BQoKxSKhp4FfBXSo207PINEG7ZtxLGVA29tNZJSqlRwBqMKRzmaa33O+A1FExaQ9Q22PYZHPgOUpKhVg+j5qp6xwxDbW2VmJLI2MixRF2JYlbXWVQNqGr/uIUQQhRs3sWg33z49S3YMAXO/4d//y94sdWLjPl1DAsPLGRI/dx/kZaR51nlmGRprQdk89Q0K/ueBHqle7wKWJXn6AqjxDjY/42xUPOp3eATAC0fg5bDjcU+80hrzVtb3mLrqa282e5NWgS1sF/MQgghCheTCTq/bLScrBgFczrT+cEv6RjSkU/+/IRuVbtlO6didmQKiKwKdo+2guTySVj3hrHczfInIPE69J4KYw9Aj0n5SrAAFuxfwNf/fM2jDR7l7ltlpgwhhBA2aNQfBv0AcZdhTmdeDApFoXhry1vYMjAuPRl5npUkWY6kNRzfAhGD0R/UR294n3/iSrG6/JOYB62DFsPS+lzlx7pj65i6cyrdbunGqCaj7BC4EEKIIqNKK2OG+BKVCIoYzlNlWrDhxAbWHFuTq9PI+rdZSZLlCIlx8Mci+KwDzOsO//7C/hIdmK6GsljfxbZzfkR8lXU9qbzYH7ufFza8QIOyDXir/VuYlBSpEEKIXAq8BYb9BLd2ZsCWL6jnGcA7W9/hckKWFfOyJevfZmW/yZMEXDoBO+bCzvlwLRbK1YE7P4CG/fn+/Y9IwJhbyl7VqDFXY3hq3VME+gYyrdM0fD19831OIYQQRZRvSRiwFI+fXuLVXbMZULki07a9x8vt37DpcBl5npUkWfmV2iS47TM4sAJ0CtTuBa0eg2od0kYJBgUFpXUItEc16tXEq4xaN4prSdf4oucXlPUra49XI4QQoigzeUCPydQrW5OBm9/gi8PL6VOxHY1r2GdZtqJGkqy8SoyDfV8ZowRj9oBvALR50lhPMLBqlt3DwsKyDG3Nq+SUZJ5f/zz/XvyXmZ1nUjOwZj5eiBBCCJFJ86GMCqjM2g1jmfjrs3zp4Y9X1faujqrAkSQrty5Fw/a5sGuBpUmwLtz5ITS8H7yLZ3uYPatRp+yYwvro9bzU6iXaVW5nl3MKIYQo+Ow5V1Wxmt2ZkDCBp3ZMZsF3DzE89G1o9ICdIy7cJMmyhdZwfLOx3M3BHwB9o0mw6u15mjg0r5YeWsrCgwt5qO5D9K/T32nXFUII4f7sPVdV6G0P0uXU73yqN9D9+5GEnPsbOr5kzLMlciRJ1s0kXoe9libB03vBtxS0GWlpErzF6eH8fuJ33t72NncE38FzzZ9z+vWFEEK4N0fMVfVCm1fYvPxu3qx2G59ueB917h+499Obtt4Ig6Si1lyMgrWvwtS6xky4OgX6TIOxB40FNV2QYP1z4R+e/e1Zbi11K+92eBcPk4fTYxBCCOHeHDFXVYXiFXi66dNsSrrAqjZD4OD3EN7TmGRb3JQkWam0hqMbYdnDMK0hbPoIqrY3ZsJ94ndoNthY78kFzl0/x8h1IynmWYwZnWdQzMs1cQghhHBvjpqrqn/t/jQo24B3L/7BpbBwiD0MszvByT/scv7CSpoLE67B3gjYNgtO7wO/QGj7tDEbe6kqro6OuKQ4nv7laS7GXyS8RzhBxWUGXSGEENY5aq4qD5MHr7Z5lf4/9OeDS3/y2rCfYHF/mNcT7vsM6slybtYU2Zosn7gzsPYV+KAefP+0sbHPR0aTYNeJbpFgpegUJmycwL5z+5h8+2RuK3Obq0MSQghRRNUuXZtH6j3C1/98zU7ijKV4gurDl4/A+ilGi5DIoOglWVpDxBBab3kMNk03RgcOXgWPb4Rmg8DLz9URppnxxwx+OvYTY5uNpXOVzq4ORwghRBH3eKPHqVS8Eq9vfp0Ev1JGl5oGYfDLG/Dt45AU7+oQ3UrRS7KUguJlOV7lXhi9B/p/AVXbOXUaBlss/3c5s/fOpm/Nvgy6bZCrwxFCCCEo5lWMCa0ncOTSEebtmwdevnDfbGNahz1LYUEfMJ91dZhuo+glWQC93uO/6o9AqRBXR2LV9pjtTNw8kVYVWzGh9YS0kSJCCCGEq3UI7kD3qt2ZvWc2Ry8dNSop7hgH/cLh1J8wpxOcPuDqMN1C0Uyy3NjRS0cZ8+sYQkqEMDV0Kl4mL1eHJIQQQmTwvxb/w8fDhze3vJk2Lxf174Mhq4wmw7nd4J+1rg3SDUiS5UYuxl1k1C+j8FAezOw8k5LeJV0dkhBCCJFFuWLlGNNsDFtjtvLDkR9uPFG5GTz6K5SuCovvhy2fFukO8ZJkuYnE5ETGRI7hlPkUH3X6iJASeWvKNJvNhIeHM3nyZMLDwzGbzXaOVAghhIB+tfrRqFwj3tv+HhfiLtx4IqAyDFltLD+3+n+wciwkJ7ouUBeSJMsNaK15bfNr7Dy9kzfavUHj8o3zfK7UdasSEhLS1q0SQggh7M2kTLzS5hWuJFxh6s6pGZ/08Yf7v4B2Y2DHPFjYF65fsHqewkySLDcwZ+8cVhxewZONnqRX9V42H2et1soR61YJIYQQ1tQKrMWg2wax/N/lbI/ZnvFJk8mYd/Luj+HYJpjT1ZgpvgiRJMvF1hxdw0d/fESvar14vNHjuTrWWq2VI9atEkIIIbLzWKPHCPYPNubOSk7IukOTgfDId3At1liK578Nzg/SRSTJcqE9Z/cwYeMEmpRvwuvtXs/1VA3Waq0ctW6VEEIIYY2fpx8vt36Zo5ePMmfvHOs7VW0Hj64D/wrwxT2w63Onxugqsnahi5wwn+CpX56inF85Puz4IT4ePrk+R1BQEFFRUWit02qtHLVulRBCCJGdtpXb0qtaL+bsnUOPaj2oHlA9606lq8PwtRAxGFY8BWf/gq6vg8nD6fE6i9RkucCVhCuMWjeKxOREZnaZSWnf0nk6j9RaCSGEcBfjWozD19OXNza/cWPurMx8A+DBCGjxKGyeAUsHQvwV5wbqRFKT5WRJKUmM+20cRy8d5ZOun1jP9m0ktVZCCCHcRVm/soxtNpaJmyey/N/l3FvzXus7enhC7ylQrjb8+D+Y1wMGLHXbVVjyQ2qynEhrzdvb3ub3k7/zUuuXaF2xtatDEkIIIezmvpr30bR8U97f+T7n487ffOeWj8LACLh43OgQH7X95vsXQJJkOdGig4tY9tcyhtw2hL61+ma7n0woKoQQoiBKnTvrauJV3t/xfs4H3NoZhq0F72Iwvzfs/crxQTqRJFlOEhkVybvb36VLlS6MaTbmpvvKhKJCCCEKqhqlajC0/lBWHF7BllNbcj6gfB0Y/ouxJM/Xw+DXyYVmKZ4ckyyl1Dyl1Bml1L502xorpbYopXYrpXYopVpmc+xRpdTe1P3sGXhBcuj8IZ5f/zz1ytRj0u2TMKmbv+0yoagQQgh3kZfWlUcbPEqVElV4Y/MbxCfH53yR4mWMubQaPwS/vQ1fDYXE63aI3rVsqcmaD/TItO1dYKLWujHwiuVxdjpqrRtrrZvnKcIC7sy1M4xcN5KS3iWZ3mk6fp5+OR4jE4oKIYRwF3lpXfH19OXlNi9z/MpxZu2ZZduFPL3h7hnQZSLs/9ZoPrxyOp/Ru1aOSZbWej2QufeaBkpafg4ATto5rkLhWuI1Rq0bhTnBzMzOMylXrJxNx8nUDEIIIdxFXltXWldsTZ/qfZi3bx6HL9q4nI5S0H4M9F8IZw4aHeJj9uYxctfL6xQOY4A1SqkpGIla22z208BPSikNfKa1tjGdLdjMZjPLIpbxVfxXnPA9wXvt3qN26do2Hy9TMwghhHAks9lMREQEMTExBAUFERYWhr+/v9V9rU18bavnWjzH+hPreX3z64T3CM+xu0yaunfC0NWwZADM7Q5950Ad29f2dRcq2wnD0u+kVFXgB611fcvjj4DftNZfK6XuB0ZorbtYOa6S1vqkUqo8sBZ4ylIzZu0aI4ARABUqVGi2dOnSvL4mm5jN5mx/ofIrNjaWH6/9yNaUrXT16Eo7v3aUKVPGIdcqTBxZJiLvpFzcj5SJeypI5RIbG0tCwo11Br29vbO9T6WkpHDhwgUSExPx8vIiMDAQk8n2cXNbzFtYFLuIAaUH0LZEdnUy1nnHn6f+vkmUuPIvR6oPIirkHqO2y0bOKpOOHTvutNYtKq9J1iWglNZaK6Pz0CWtdckczvEaYNZaT8npes2bN9c7dji2n3xkZCShoaEOOfewacPYVmobNS7XoMn5Jnh7ezN+/HiHXKswcWSZiLyTcnE/UibuqSCVy+TJk7MkWY66T2mtGbpmKH9d+IsV96ygrF/Z3J0g8Tosf8Lop9XkIej9gdF/ywbOKhOllNUkK69TOJwE7rD83An4x8oFiyulSqT+DHQD9mXer6DLPOril8O/sD1gO0HXg2h0vpF0XBdCCOF2nDnASinFy21eJi4pjve2v5f7E3j5Qd95cMf/4I+F8MW9cC2HiU7dhC1TOCwBNgO1lVLRSqlhwKPA+0qpP4FJWJr5lFKVlFKrLIdWADZa9tkGrNRar3bEi3Cl9KMu9sfs5/kNz1M9oDp9ffri6+0rHdeFEEK4HWcPsKoeUJ3hDYaz6r9V/H7i99yfwGSCji/CfXMgervRIf7s3/YP1M5y7PiutR6QzVPNrOx7Euhl+fkI0Chf0RUAqaMu4kxxbCy3EZWi+KTrJ1T0r+jq0IQQQgirXDHAaniD4fz434+8vul1+pr7EhsTm2On+ywahkHgLbD0QZjTBe6fDzU6OTTu/JAZ3/MpKCiIFJXCpvKbiPOI4z51nyRYQgghRCbeHt680uYVTl47ydora/O0qonZbCb85/3MTOzH+ZRi6IX9YNtsB0adP5Jk5VO/fv3YX2U/533P05OePBX2lKtDEkIIIdxSi6AWVL9Wnb9K/sUlr0s3nXfL2kzzqV10ziX6MSupH1G+dWDVc7BqHCQnOfnV5EySrHz6/PDn/KX+YnTT0UweMrnADN8VQgghHM1aotTdtzteKV7sLLMTFNl2urc203z6iVHj8WJxcm9oMwq2zYLF90PcJWe+vBxJkpUP3x/+nk///JR7br2HYfWHuTocIYQQwq1YS5QeCXuETqoT533PcyHkAmFhYVaTMWszzWceFVmhYiXo/hb0+Qj++w3mdoPz/7nyJWcgSRZ5W/xy1+ldvLrpVVoEteCV1q+kFboQQgghDNYSJX9/f94d/C6tglqxzWcb103XrSZj1qaZyHZUZLNB8PC3cCUG5nSGY5tc9ZIzkCSL3C9+GXU5itG/jqayf2U+CP0ALw8vJ0UqhBBCFBzZzceVOndWfHI872x/x2oyZi2hSh0VOX78eIYMGZKxi061DvDoL+AXCAvugt2Lnf56M8vr2oWFSm4Wv7wUf4kn1z2JRjOz80wCfAKcFaYQQghRoISFhWVZIzHVLSVvYUTDEczYPYN7gu7BK8orw/qIeZpmokwNGP4zfPkILH+CalX6QocOxjxbLiBJFrYvfpmYnMjYyLFEm6OZ3XU2VUpWcXKkQgghRMGRU6I0pP4QVv23ii2JWwgLCcswd1ae+QXCQ9/Aqucov/9HiLsIxUrn/Xz5IM2F2DbzrdaaN7e+ybaYbUxsO5HmQVmWKBJCCCFELqTOnRVzLQZzA7P1ZsC88PCCOz9kZ7MpLkuwQGqyANtmvg3fH843/3xDK92K/V/uJzYoNnez1AohhBAii2YVmtG3Zl8+P/A5var3ok7pOvY5sVIkeZW0z7nySGqybPDzsZ/5YOcH1Na1CTkekqdZaoUQQghh3TPNniHAJ4DXN79Ockqyq8OxG0mycrD/3H7GbxhPw3INaXSqERj943PsIC+EEEII2wT4BDCuxTj2ntvLsr+WuTocu5Ek6yZirsYw6pdRlPErw7SO0wgOCrY6FFUIIYQQ+dO7Wm/aVGzDR398xOmrp10djl1IkpWNq4lXGbluJHFJcczoNIOyfmVt6iAvhBBCiNxTSvFy65dJSkni7W1vuzocu5CO71YkpSQx7rdxHL54mI87f8ytgbcCtnWQF0IIIUTehJQM4fFGjzNt1zR+Pf4rHat0dHVI+SI1WVZM2TGFDSc20DG5I799/pvNS+0IIYQQIn8G1RvEraVuZdK2SVxLvObqcPJFkqxMFh9czKKDi2immxEYFSgjCYUQQggn8vLw4tU2rxJzNYYZu2e4Opx8kSQrnfXR63ln+zuEhoRS81RNm5faEUIIIYT9NC7fmLBaYSw6uIj9sftdHU6eSZJl8feFvxn32zhqB9bmndvfoVJQJRlJKIQQQrjImGZjCPQJZOKmiSSlJLk6nDwpckmW2WwmPDycmJiYtL5W566fY+S6kfh7+TO903SKeRWTkYRCCCGEC5X0LskLLV/g4PmDLD201NXh5EmRG10YERFBVFQUNWvWJCoqisURi/m13K9cir/Egh4LqFC8AiAjCYUQQghX6161O98d/o7pf0ynyy1dCCpesFqVilxNVkxMTFpfqxSdwvL45eyP3c87t79D3TJ1nRpLaq3a5MmTZQSjEEIIkYlSigmtJpCiU5i0dZKrw8m1IpdkBQUFpfW12h+4nyi/KJ5r/pxL5uJIrVWTEYxCCCGEdcElgnmy8ZP8GvUr646vc3U4uVLkkqzUvlZ7UvZwKOAQ91S7h4frPeySWNLXqskIRiGEEMK6h+o9RK3AWkzaOglzQsFp9SlySZa/vz+39biNVUmraFOxDa+0fyWtZsvZ0teqyQhGIYQQwjovkzF31tlrZwvU3FlFLslKSknitc2vUd6rPFNCp+Bl8nJZLDKCUQghhLBNw3IN6V+7P4sPLmbfuX2uDscmRW50oafJk487f8yObTso6V3SpbHICEYhhBDCdk83fZpfjv/CxM0TWdJ7CZ4m905jilxNFkDVgKqU9Srr6jCEEEIIkQslvEvwQqsXOHT+EIsOLnJ1ODnKMclSSs1TSp1RSu1Lt62xUmqLUmq3UmqHUqplNsf2UEr9pZT6Vyn1gj0DF0IIIUTR06VKF+4IvoOZu2dy0nzS1eHclC01WfOBHpm2vQtM1Fo3Bl6xPM5AKeUBzAR6AvWAAUqpevkJVgghhBBFW+rcWQBvbX0rbZS+O8oxydJarwfOZ94MpHZoCgCspZItgX+11ke01gnAUuDufMQqhBBCCEFF/4qMbDyS9dHrWXtsravDyVZe+2SNAd5TSkUBU4DxVvapDESlexxt2SaEEEIIkS8D6w6kbum6vL3tba4kXHF1OFYpW6rZlFJVgR+01vUtjz8CftNaf62Uuh8YobXukumYMKC71nq45fHDQEut9VPZXGMEMAKgQoUKzZYudexikGazGX9/f4deQ+SOlIl7knJxP1Im7knKxfmOxx9nSswU2vu35/4y92d53lll0rFjx51a6+aZt+d17OMgYLTl5whgjpV9ooGQdI+Dsd6sCIDWehYwC6B58+Y6NDQ0j6HZJjIyEkdfQ+SOlIl7knJxP1Im7knKxTVObTvFooOLeKzDYzQq1yjDc64uk7w2F54E7rD83An4x8o+24GaSqlqSilv4AFgRR6vJ4QQQgiRxagmoyhfrDwTN08kMSXR1eFkYMsUDkuAzUBtpVS0UmoY8CjwvlLqT2ASlmY+pVQlpdQqAK11EjAKWAMcBL7UWu93zMsQQgghRFGk4zWtzK3458I/jF4wGrPZfdY2zLG5UGs9IJunmlnZ9yTQK93jVcCqPEcnhBBCCHETERER+ET5UKlsJTb5bWJOxBzGDBnj6rCAIjrjuxBCCCEKh5iYGLTWND7fGIViVcIqt5k7S5IsIYQQQhRYQUFBKKUollyM+hfrc8r3FGuOrnF1WIAkWUIIIYQowMLCwggJCcHb25vQUqHUKVWHt7e9zaX4S64OLc9TOAghhBBCuJy/vz9DhgxJe9wmtg0DVg5g2q5pdKCDCyOTmiwhhBBCFCL1ytTjoboPEfF3BEfijrg0FkmyhBBCCFGojGw8kqDiQSw9v5TEZNfNnSVJlhBCCCEKlWJexZjQagLVfaq7dIJS6ZMlhBBCiEInNCQUDhsJl6tITZYQQgghhANIkiWEEEII4QCSZAkhhBBCOIAkWUIIIYQQDiBJlhBCCCGEA0iSJYQQQgjhAJJkCSGEEEI4gCRZQgghhBAOIEmWEEIIIYQDSJIlhBBCCOEAkmQJIYQQQjiAJFlCCCGEEA4gSZYQQgghhANIkiWEEEII4QCSZAkhhBBCOIAkWUIIIYQQDiBJlhBCCCGEA0iSJYQQQgjhAJJkCSGEEEI4gCRZQgghhBAOIEmWEEIIIYQDeOa0g1JqHnAncEZrXd+ybRlQ27JLKeCi1rqxlWOPAleAZCBJa93cLlELIYQQQri5HJMsYD4wA/g8dYPWun/qz0qp94FLNzm+o9b6XF4DFEIIIYQoiHJMsrTW65VSVa09p5RSwP1AJzvHJYQQQghRoOW3T9btwGmt9T/ZPK+Bn5RSO5VSI/J5LSGEEEKIAkNprXPeyajJ+iG1T1a67Z8A/2qt38/muEpa65NKqfLAWuAprfX6bPYdAYwAqFChQrOlS5fm6oXkltlsxt/f36HXELkjZeKepFzcj5SJe5JycT/OKpOOHTvutNbv3JY+WVYppTyB+4Bm2e2jtT5p+feMUupboCVgNcnSWs8CZgE0b95ch4aG5jU0m0RGRuLoa4jckTJxT1Iu7kfKxD1JubgfV5dJfpoLuwCHtNbR1p5UShVXSpVI/RnoBuzLx/WEEEIIIQqMHJMspdQSYDNQWykVrZQaZnnqAWBJpn0rKaVWWR5WADYqpf4EtgErtdar7Re6EEIIIYT7smV04YBstg+2su0k0Mvy8xGgUT7jE0IIIYQokGTGdyGEEEIUKmazmfDwcGJiYggPD8dsNrskDkmyhBBCCFGoREREEBUVhdaaqKgoIiIiXBKHJFlCCCGEKFRiYmJInaJKa01MTIxL4pAkSwghhBCFSlBQEMaiNKCUIigoyCVxSJIlhBBCiEIlLCyMkJAQlFKEhIQQFhbmkjjyPBmpEEIIIYQ78vf3Z8iQIQV6MlIhhBBCCJENSbKEEEIIIRxAkiwhhBBCCAeQJEsIIYQQwgEkyRJCCCGEcABJsoQQQgghHECSLCGEEEIIB5AkSwghhBDCASTJEkIIIYRwAJW6gKI7UUqdBY45+DJlgXMOvobIHSkT9yTl4n6kTNyTlIv7cVaZ3KK1Lpd5o1smWc6glNqhtW7u6jjEDVIm7knKxf1ImbgnKRf34+oykeZCIYQQQggHkCRLCCGEEMIBinKSNcvVAYgspEzck5SL+5EycU9SLu7HpWVSZPtkCSGEEEI4UlGuyRJCCCGEcJgil2QppXoopf5SSv2rlHrB1fEUVUqpEKXUr0qpg0qp/Uqp0ZbtpZVSa5VS/1j+DXR1rEWNUspDKfWHUuoHy2MpExdTSpVSSn2llDpk+ZtpI+XiWkqpZyyfXfuUUkuUUr5SJs6nlJqnlDqjlNqXblu25aCUGm+5//+llOru6PiKVJKllPIAZgI9gXrAAKVUPddGVWQlAc9qresCrYGRlrJ4AVinta4JrLM8Fs41GjiY7rGUietNA1ZrresAjTDKR8rFRZRSlYGngeZa6/qAB/AAUiauMB/okWmb1XKw3GMeAG6zHPOxJS9wmCKVZAEtgX+11ke01gnAUuBuF8dUJGmtT2mtd1l+voJx06iMUR4LLLstAO5xSYBFlFIqGOgNzEm3WcrEhZRSJYEOwFwArXWC1voiUi6u5gn4KaU8gWLASaRMnE5rvR44n2lzduVwN7BUax2vtf4P+BcjL3CYopZkVQai0j2OtmwTLqSUqgo0AbYCFbTWp8BIxIDyLgytKPoQeB5ISbdNysS1qgNngXBLM+4cpVRxpFxcRmt9ApgCHAdOAZe01j8hZeIusisHp+cARS3JUla2yfBKF1JK+QNfA2O01pddHU9RppS6Ezijtd7p6lhEBp5AU+ATrXUT4CrSDOVSlj4+dwPVgEpAcaXUQ66NStjA6TlAUUuyooGQdI+DMap4hQsopbwwEqxFWutvLJtPK6UqWp6vCJxxVXxFUDvgLqXUUYym9E5KqYVImbhaNBCttd5qefwVRtIl5eI6XYD/tNZntdaJwDdAW6RM3EV25eD0HKCoJVnbgZpKqWpKKW+MDnArXBxTkaSUUhh9TA5qraeme2oFMMjy8yDgO2fHVlRprcdrrYO11lUx/jZ+0Vo/hJSJS2mtY4AopVRty6bOwAGkXFzpONBaKVXM8lnWGaNfqZSJe8iuHFYADyilfJRS1YCawDZHBlLkJiNVSvXC6HfiAczTWr/l2oiKJqVUe2ADsJcb/X9exOiX9SVQBeODLExrnblTo3AwpVQo8JzW+k6lVBmkTFxKKdUYYzCCN3AEGILxJVnKxUWUUhOB/hgjpf8AhgP+SJk4lVJqCRAKlAVOA68Cy8mmHJRSE4ChGOU2Rmv9o0PjK2pJlhBCCCGEMxS15kIhhBBCCKeQJEsIIYQQwgEkyRJCCCGEcABJsoQQQgghHECSLCGEEEIIB5AkSwghhBDCASTJEkIIIYRwAEmyhBBCCCEc4P9lQB7/TdO6gwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "func_arr = [estimator.model_function for estimator in estimators]\n", "u = one_bkpt_estimator.params[0]\n", "u1, u2 = two_bkpt_estimator.params[[0,2]]\n", "extra_pts_arr = [None, [u], [u1, u2]]\n", "\n", "analysis.plot_models(func_arr=func_arr,\n", " indep=indep,\n", " dep=dep,\n", " extra_pts_arr=extra_pts_arr,\n", " mark_extra_pts=True,\n", " scatter_size=15,\n", " scatter_color=\"gray\",\n", " marker=\"o\",\n", " legend=[\"ols\", \"1bkpt\", \"2bkpt\"]);" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }