{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "

CSCI-UA 9473 Introduction to Machine Learning

\n", "\n", "

Assignment 1: Gradient descent\n", "

\n", "\n", "\n", "\n", "__Given date:__ February 9\n", "\n", "__Due date:__ February 21\n", "\n", "__Total:__ 15pts" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Question 1. (5pts) Local vs global minimas and gradient descent\n", "\n", "We consider the following function.\n", "\n", "\\begin{align}\n", "F(x_1, x_2) = 3(1-x_1)^2\\exp(-(x_1^2) - (x_2+1)^2)\\\\ \n", " - 10(x_1/5 - x_1^3 - x_2^5)\\exp(-x_1^2-x_2^2)\\\\\n", " - (1/3)\\exp(-(x_1+1)^2 - x_2^2)\n", "\\end{align}\n", "\n", "The surface plot of this function is given below together with its contour plot. The function has a single global minimum located near $(0.23, -1.62)$ and shown in red in the contour plot.\n", "\n", "We want to implement gradient descent iterations on that function. Starting from a random initial point $(x_1, x_2)$, code the following updates \n", "\n", "\\begin{align}\n", "x_1^{(k+1)} = x_1^{(k)} - \\eta * \\text{grad}_{x_1} F(x_1, x_2)\\\\\n", "x_2^{(k+1)} = x_2^{(k)} - \\eta * \\text{grad}_{x_2} F(x_1, x_2)\n", "\\end{align}\n", "\n", "where $\\text{grad}_{x_i}$ represents the gradient of $F(x_1, x_2)$ with respect to $x_i$. Choose a sufficiently small learning rate and plot the iterates (in white) on the contour plot. Repeat your experiments for various initial iterates. " ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPkAAADyCAYAAACPiGNTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAACyoUlEQVR4nOz9eZRkaXrWCf6+u9i+m7m5u/kSHntExpaxZmZlSSUQSA0twfRANzAsgjkNNLSYZrpnRi1Q0z0DYhmWmT4ItegWg+gWINECLQhtpaWUmZWZlbFk7LuH7+bm7rbvd/3mj2tm4eHhHu4RGRGVleXPOXUyyuz6tWvX7vO97/cuzyuklOxiF7v44kL5Zl/ALnaxi1eLXZLvYhdfcOySfBe7+IJjl+S72MUXHLsk38UuvuDYJfkudvEFh7bN+7v5tV3s4tVDvMqT71ryXeziC45dku9iF19w7JJ8F7v4gmOX5LvYxRccuyTfxS6+4Ngl+S528QXHLsl3sYsvOHZJvotdfMGxS/Jd7OILjl2S72IXX3DsknwXu/iCY5fku9jFFxy7JN/FLr7g2CX5LnbxBccuyXexiy84dkn+TYCUEtM0sW2bXUnsXbxqbCcasYuXDNd1MU2Tbrc7eE1VVXRdR9M0VFVFiFeqIbCLbzOIbSzJrpl5SZBSYts2tm0jhMCyrMHrUkpc1x2Q2zAMotEoPp9vl/TfHnilP/CuJX8N6Lvn64nchxACIQSKogyOnZ6eZmpqilAoBOxa+l18NuyS/BXDtm0WFxdxHIexsTGEEAPrvRlZ+6RXVRVVVQdWvtPpDI7XNG3wv13S72I77JL8FWG9e+667sBN78OyLHRd35agm1l6x3GwbXtwjKZpA0uvKMou6XfxBHZJ/grgui6WZQ3c87717r/38OFDVlZWAIhEIiSTSZLJJMFgEOCJ4zeif74+NpJeCPGEpd8l/S52Sf4S0SdcP6jWt7590nY6HW7cuEEqleL8+fMAtFotKpUK9+/fHwTcOp0OpmkSDoe3/czNSG/bNpcuXeLUqVO7pN/FLslfFqSUWJaF4zhPEU8IQaPRYHl5maNHj5JKpTBNEykl0WiUaDTK5OQkruvSaDS4d+8e09PTuK5LLBYbWHqfz7ftdfQ/23GcwZ7etu3BwrNL+m8/7JL8JaCf++4H09aTxnVdlpaWaLfbvPXWW88kqqIoxONxIpEIk5OThEIh6vU6lUqFpaUlHMchHo+TTCZJJBLour7ttW1m6S3Leor0uq6jquou6b+A2CX5Z8DG3HffPe+j1Wpx48YNIpHIIO+9E/Tde0VRSCQSJBIJ9u7di+M41Go1KpUK8/PzSCmfIL2mbf9z9iP367/DRtL3g3iapj21SOziWw+7JH9BSClZWVkhFArh8/meIkI+n2d2dpZjx45hmiaVSuW5z78RqqqSSqVIpVKAl57rk352dhYhBIlEgmQyueNy2c1Ib5omhmGQz+cZHh4mFAoNLP0u6b/1sEvyF0A/uLawsMCePXvw+/2D92zb5s6dO7iuy4ULF9A0jWKx+Fw16jslkaZppNNp0uk04KXlqtUqxWKRdrvNlStXBqSPx+NPeRpbfXaf9JVKhWw2OyB931vZ6N7v4vONXZI/B7ZzzxuNBjdu3GBycnJQ+ALPTolthuc9vg9d1xkaGmJoaIharcbx48epVqusrq7y4MEDdF0fBPGi0eiOCKooyoD0/WsyTRPTNAfv75L+841dku8Qm+W+FUXBdV2klCwsLLC0tMTJkyeJRCJP/O2Lkvazwufzkc1myWazgFcTX6lUyOfzNBoN/H7/YD8fjUaf8iA2XnP//V3Sf2thl+TbYGPue/2etN9ocvXqVXw+HxcuXHhif9vH67Lk28Hv9zMyMsLIyAgAnU6HarXK4uIijUaDYDA4sPT9HP2ztg6bkb6/p++Tvp/KC4fDu6T/JmGX5M/As3Lf4FnGu3fvcujQoQFxNsM3y5Jvh2AwSDAYZHR0dFCs0w/itVotTNNkeXmZdDpNKBTacQluH1JK1tbW6HQ6TE5OAp6lX99ss0v6V49dkm+BZ+W+pZTMzMxQrVY5cODAMwkOnx9Lvt1nhkIhQqEQY2NjSCm5fPkyUkoePXpEu93etAR3u3MCTzTbSCkxDAPDMJBSPuHa91N2u3i52CX5Bqx3zzcLrhmGwY0bN4hGo+RyuR1XoX0eLfmz0P/uY2Nj7NmzByklzWbzqRLcfvQ+EAhsep713XabLZau6+4KaLxi7JJ8HaSUNBoNFhcX2bdv31MPWLFY5N69exw6dIihoSEePny4I/JuRvLtmlA+L4vCeoJuLMFtNpuUy2Xu3r2LaZqbluBu1VLbP+cu6V89dkneQ98971eVbSxNffjwIbVajbNnzw6sVj+6vh02klZKSaFQQNM0EonEpsG6zwOetdAoikIsFiMWiwHePdqsBNd13YH4xXbYivSdTofZ2VkmJibw+/27pH9OfNuTfGPuW1XVJ4jb6XS4fv06mUyGc+fOPdV48ryW3DRNbty4MaiSe/ToEZqmPZG//jxa8u2wVQnu/Pw8lUqF1dXVFyrB7W8bGo0GQoinBDR2Lf32+LYm+Va57z7BVlZWePjwIW+88QbJZPKpv39eklerVW7dusXBgwdJJpODzzVNk3K5PMhf92MCfr+fcDj8TXtwP8tC0y/BbbVaaJo2KNDZrAQ3Ho9v6824rjuIxq8X0NhVzdke35Ykf1buW1EUHMfh9u3bdLtdzp8/v2VwbafuOnjNKnfv3uX06dOEQqHBZ4NXtNLPX0spuX//PsAglfW8Ue2Xic9Kkv6efKsS3FKpxPT0NKqqPrMEd7O9/VaqOf2gKeyq5sC3Ick3uucbf/R2u02j0WB0dJSjR49uWwyynbWzLIvbt2/jOA5f+tKXts0L97vAYrEYmUxm06j28/aYvyhexpbBdd1Nv/P6ElzwtjH9EtyHDx8+tYV5VgCvj61I39e3F0JgGMZgu/DtQvpvK5Jv5p6vx9LSEnNzcwQCAfbs2bPt+friDFuhVqtx8+ZNJicnWVlZ2XHhx/rFY7Oo9voAl+u6z73XfR68LEu+HbYrwe10OszPzz8Rt9jJtW8M5N28eZOzZ88O3v92END4tiD5drlv27a5ffs2ABcuXOCTTz7Z0XnX7983ft78/DzLy8ucPn0aTdMoFAo7vt5neQibBbiq1epgr9t/vy8g+c2uKNspyTdiYwnuN77xDXw+35YluDslPTxZhvvtoJrzhSd5vzT1ypUrvPnmm0/9aPV6nZs3b7Jnzx7Gxsae69ybkdGyLG7duoXP5+P8+fOoqoplWa8sWq6q6lN73b6Vv3TpEj6fj1QqRTKZJBKJvPaH9kVJvhGKojA6OrplCW4oFBqQficluPDto5rzhSb5+tLUdrv91A86Pz9PPp/ftHNsJ9hI8v6CsXfvXkZHR7c87nnP+zzQdZ1sNsvc3Bznz5+n2+0OlGSazSahUGhA+mAw+Mof2pdF8vXYrAS33W5TqVQ2LcENBAI7Jv1OVXO+lUj/hST5Zn3f638M0zS5desWfr9/y86xnWCzVtNTp049pbL6zcx7BwKBJyxgnwwPHz6k2+0SjUYHZFgvfvGy8CpIvhFCCMLhMOFwmPHx8U2DlX3S7zQb0j/vTlRzPu9SWV84km8cSbTxplcqFW7fvs2BAwcYHh7+TJ/VD7xdv34dVVWfq9X0Wem3V7UobCTD+tLU27dvY9v2IIi3WV3Ai+B1kHwjnlWCaxgGn3zyyQtlKNaTvlwuMzw8/IRqTt/Sf9566b9QJN8YRNnonk9PT1MsFjlz5sxLyTd3Oh3y+TxHjhwhl8ttedxmpP08VLStL02dmppaJxS5wtLiHVotl+np6R0XrGyGbwbJN6L/PUOhEKVSidOnT29agtvP0+9EBbffJ7+VgMYP/uAP8iM/8iMcOXLklX63neALQfLtct+GYdDpdLAsi/Pnz3/mFVZKOUi3ZTKZZxL8RfDNcu9VUWbI/18yNjyNMgJdK0xbfi9za/8p09PTg9x1KpXacRrrZZD8Zd2LfrZhswxFn/QLCwu4rjsg/FZpyX4FXh8bI/flcvm1Fy5thW95km+X++53jgWDQfbv379jgm/1cPbTbUIIjh49yurq6rbn+mZbsp3Ase+gd/40oOIIFUc66FqdOD/HwfEGHPhRTMseWL9Go0EgEBiQfquI9ssg+ctKBW4kZh+qqj6xRdmoggsMCL++oehZ36tfqfh5wLcsybcaSdSH67o8ePCARqPBuXPnuHHjxo6DLn1LuvFHbDabXL9+ncnJScbHx6lWq5/JymxFgNdtyV23gtH+ATo08RFAxcRwXYIiTEh1cZ3fwWp9L3rwl54ov+2nsWZmZp4ov02lUoNOvZdlyV/GQrnTxWKzEtxarUa5XObRo0eoqjrIWsRisU0Xjn6E//OAb0mS91MbN27c4NixY5uWpt64cYOhoSHOnj07iLDvlOSbFbn03fMTJ04QjUYHxz1PtLaPVqvFtWvXnnALk8nkS69W2wlc16Lc/PMEZA0pdFThabMhJF3a2A4IYaDJOkb3LxCK/LT39iZprH5Eu99fHo/H6XQ6T0xgfbFrfHmW/EXOo+s6mUyGTCYDePvuS5cuPVGCm0gkBtsYRVGwLOu5So6FEIeBn1330j7gb0gp/7/rjvku4BeBmd5L/05K+f/a7tzfciRfn/uu1+tPEbxQKDA9Pf1U59jzkrzv2jmOw507d3AcZ6Cj3seLWNxCocCjR484evQoPp+PRqNBpVJhbm4OIQTJZBLHcV7bfq7U+r9iOTcIquBIG4SCjY4iXEBSdcOEpENAMVHlLTqdv0Mw+MNPnedZ5bcPHjx4ouvseXvopZTfVJJvhM/nQ9d1Dh8+DDwuwV1eXubevXv8g3/wD3Ach0uXLnHmzJkdfVcp5T3gTQAhhAosAT+/yaHvSym/73mu91uG5NtpnjuOw7179zAMgwsXLjwVIX0ekvf7lpvNJjdu3GB8fJzx8fFNu6B2SnIpJXfv3qXdbnP+/HmEENi2/cRElH612sLCAqurq5TL5W33vJ8FbfMOpv1b+JQw0EJFA1wMGcAnDABsaVHHQHc1DKHTMH6GIe0/w6fvf+a5+8GtcDjMgQMH0HV94PLOzMygKMrAg4nFYs8kXz/e8lnRj4h/Vmz8zTeW4P7tv/23+ZN/8k/yT/7JP2Fubo7f/u3fft7r/25gWko595kvlm8Rkm+X++6TcWxsbMvOsee15IVCgaWlJY4fPz5QP3nRc3a7XdrtNqOjoxw+fHhA8I3oV6v1685TqdSgiqvT6RCJRAaLwsvoPlvp/HX8wgVpU3Z1dBHFdF0M2cGnguWCKroIoCTDuG4NKYIU23+VXPzf7+gz+lZ44z6333VWKBS4f//+QAN+s/Lbb7a7vhHbLRb79+/H5/PxUz/1Uy/6EX8c+NdbvPeOEOIakAf+b1LKW9ud7HNP8n5wbSvV1MXFRebn559JRtg5IR3HodFoIKV8yj3fiJ1Y8nK5zJ07d/D7/ezbt2/bz++fF3hqz9toNCiXy9y8eRPHcQb7wBeRkFpo/QtsZxFFhFlzXFwUfEoLAZiuil+JYMl2z20HU5o4UiEgOtiyQN34dWL+7932c7YKmm3sOusH8frlt+FweED6l2XJXxfJbdt+YY9BCOED/hDw9J4IrgB7pJRNIcQfBH4BOLjdOT+3JN/OPbdtm263S7lc3paMsDOSt1otrl+/jq7rHDp0aNtzPovkUkpmZ2dZXV3l7NmzXLly5am/fRY2E37cWLhSrVafiPj2rXw/MLgVHNek0P3nhHFYdWKElRoQQNDClaAIlyXbJYwg2ntWu9LGkgECtHFli1rnH+6I5DslaF8DPpfLIaWk1WoNym9brRZCCAqFwmcqv32ZJH/WeVqt1lOlzc+BPwBckVKubHxDSllf9+9fEUL8uBAiI6UsPuuEn0uSb5f77jeCaJrG8ePHd/TDbdRu24h+QOzYsWPMz8/v6Dq3Wjgsy+LmzZsEAoEXKr7ZiYewsfusH/zpt2IahsHS0tKgK2s9Hjb/MbZsUyaAJU3CgOmCXwWXIEJ4aqnLVhohaoRFB0uqKMKh7foAhTXHIdD5DwwF/+NnXueLpL+EEEQiESKRCBMTE1SrVRYWFuh2u0+U3/a9mJ1mJbYj506xVb69j8+YI/8TbOGqCyFGgBUppRRCXAAUoLTdCT9XJH+WLFP//bm5OZaXlzl16hS3bt3a8Q/Xl3XaCNd1uXv3LoZhcP78eXRdf2EVVng89HDfvn3bDl14mVgf/JFS8o1vfGOgMtvtdonFYp6Vj4com79N040RUFr0H1WJd28c/Gh4JHcUmyU7wJjmQ+kF4upuEBcVmy4z7Z98JSTf7ByBQICpqakN5bdP6sWlUqkt89awPTl3iu3c9Re15EKIEPD7gb+47rX/AkBK+RPAHwX+khDCBjrAH5c7iPx+bkgupaRSqWAYBslk8qkHwzRNbt68STAY5K233hqUJz5vWmw92u02169fZ2Rk5ImAXT+6vh02kryfS3/R1tWtzvsif6+qKhMTE0xMTAzSWeVymcvVf4wWcZCKBgoEeguk1tt7O1JBE9BxNTThEX/ZTpJQPdELQ4IhwwSUOoZbpNz9mFTg7S2v5VVUvG2c097Xi1tbW3tCOmpj+a3rui+lFmE7kjebzRf6/aWUbSC94bWfWPfvHwN+7HnP+7kgeT/33Wg0aDabgx+vj606x/p57J1gI8n7SqzHjh0jkUg889jtzuk4Dnfv3sW27R3FB143+umsaDTK1dVHWFIhJbxgptpzzf09krccC78CbdeHT/FeqzoCTfiIKCZt14eFRgBQhM7D1v/EhVdM8u3OsVEvbuPWJRgMkkql6Ha7L6WddjvvsR84/Lzgm/o0bpRl0jTtCdJK6c3h2qpz7EUsueu63Lt3j06ns6US6/O4647jcPHiRUZHR5mcnHwpUeBXVdb6aeNnaboGfgU0pQ0EMKRLxYoRVxXCShW/4lW8WVLDR//fDiUnRlgUcVBRRZuW6yMkWnTcWermAjHfxDO/z2fB8wbMNm5dOp0O5XKZcrlMsVikVCo9ISjxvNjOkn+eSlrhm0jyzXLf6y1zt9vlxo0bJBKJLYNXz2vJ2+02Fy9eJJvNcuTIkS0fvp2SrFgs0m63uXDhwlPewGb4ZrddPuq8j4vEdUEA8+Y4hmxiS4Wi6+DKCFGlw0FlDZeeiysFQoCNQsnJIUQbAMON4AqXihvno9o/5XuH/tYru+7PEhVfX37b7XYH/ePrxzut76HfSZvpTgJv3/aWfKuJof0I+NraGvfv3+fIkSOD6PFmeB5L3mw2yefzvPnmm9sKImx3XiklDx8+pFqtEgqFdkTw58GrsOS3Gh9StprENRdHCqaNPVRtSVIX6ESAGq7UqUuFO90QGa0OGHTcMP11adr0kdM7KEKiCoOiHUCKLkXrLmVjgZR/a2v+WfCyFkfHcdA0bVB+u2fPHlzXHQTxFhYWkFJuW37bP89W+Dx1oMFrJvl2uW/w9t/dbpdz585tu3/aiSXvd6OVSiVyudyOFE+eRXLTNLl+/TrxeJxz587x0UcfbXs+2Lqz7XXhSv2rNF1BApOinUYRFgreQtK0TeI6IAJAizaCB90sx0JzWDIEogOAKXWq9hApfZWmHaPj+ghrLUxH4euVf873Df93r+T7vcqKt/XlteDVX2ysP1gfxOtnaZ71bLZarafiSt9MvDaSb5f7brfb3Lx5E2DQObYdtrO43W6Xa9euMTQ0xIEDB6hWqzu61q0s6foxR/1KrZ1i/Tn7hTKNRoN0Ov2UBNHLtuT57jx5c42AorNkpFCUIGCh9wJrYhBZV9EEmI6CJQR326MM+yRK76eQEuquJC4FTRnARBJ0BVJYlK3bfPDJ+6SjmYFQ5MvSjHtZUfGdLBaapj3Vcba++cTv9w/c9a0W7VarxcTEq/FqXgSvnOTb5b7hcefYoUOHBnnPneBZJO+7/EePHiWVSlEqlV443dYXaszn84MxR8+LPnFt2+bGjRsDgcVqtToYktC3GC/Svvos/OraL2C6HRxX0pF+ksIGAQHFW0j8vf+60iNA2/WhKdAVflYtlxE/vUo476EuWklMaXupNpkgpFQwFQX7wEPG9Tef0IzrdruUSqXPNL31dfeTr4fP52N4eHiQ1el0Oty9e5disUg+nx+U36ZSqUFg+EXcdSHELNAAHMCWUp7b8L4A/kfgDwJt4M9KKa9sPM9meKUk306WqZ96Mk2TCxcuoCgK09PTOz7/Zu56vwCkXq8/4fI/byS+f96+EoyiKAMd9ReBEIJms8nt27fZu3fvQAQwmUyyd+9ebNseTP8slUoIIQaa6Z9FOrljd5jrLhHWHEpWEL9iIUQHKRUU0UVKiU94zTJtx+g1pmhoivfakhkiprW8+9K7hqYbQREdBNBxJJoIs2ZmaDm3+a6J6BOlt5988smWpbc7/U6fpwaV9eW30Wh0UH7bV4X99NNPmZmZ4Z133nmR0/+eZ5So/gG8OvWDwFvA/9T777Z4ZSTvB9c++eQTLly48MzOsYmJiYGlex4rtrGKrdvtcv36dVKp1FMu//OS3LbtQS37xMQE4+PjO76uzdAvdT116hTRaPSpa+lP/hwaGqJcLrO6uooQ4qmKtZ1GgPv496u/TNcxkYTxKV60XBUSKYNAm66rE1G9VJlf9YjtyP4sMYEiYK47SkSp4NM8i2+4UYTwEVJr+BSDuW4CVXFo21W+Xv1Nvpz8/YC3CGuaxsGDXg9Ff3prP3+9Uw3419VY8rzn2Vh+20/R/vIv/zJ/5+/8Hf7u3/27/PRP//Sg7/wz4g8D/2uvwu1jIURCCDEqpVze7g9fOsk35r4dx3mqNHVpaYn5+fknVFbg+fOp/ekkAKVSibt3724ZkX/efvJarUY+n9+2u207SCl58OAB3W6Xt956a8dunKZpjI2NMTY29kTF2sLCAuBpjqXT6UEwaKvPvlK/j0Bj1RBkA03c3la/ablENXClHzBxpCCg9KS0hHe+ru1HCIGBQ8uMMqF5/RFNR2JKQVARtJwADhEU6iB9XKp9NCD5Rmyc3tputymXy08sZH3Xd/1C9s101zfDVsUwiqLw9ttvk06n+fEf/3GmpqaeZ0GWwG8IISTwT6WU//OG98eAhXX/f7H32usleV+WqU/sjT+MbdvcunULRVFeSmVY35I/ePCAarX6zIj8dg0qfbiuy/LyMq1Wi7fffvu5rOZGrI/Ex+PxHX/fjYG39eqi8Fhcoh8M6ld0rd8XAvza2odU7BquVNAVgRASXaiYrsqqEaLuJGlbgoAaYSJgEfGV+1cAgOP6QPWIX7VDZJwuQdXEkhaaIqlacVquhisMHFcg0Gi6azxo3eNg+NnWS4jHGvDrS2/7lWr9VFYqlXppjSUvc7HYLk8ejUafN3bzrpQyL4TIAl8VQtyVUr637v3NLnxH0dmXRnIpJYZhbNr3Dd6Ez1u3bjE1NfXSJIxd12VxcZFcLse5c+ee+QPuxJL33f1gMMjw8PBnIni/U+7AgQNks1muXLmyacS8bVr8+EeXWKo1ODM2ysnRLGP+Z39uX1wim80+YRH7+8J4PI5lWXx17TKqAMN10WW/fjvG7XoEW6qEhEXbVnGUILdbSUbsKLnAPP6eRTdcBa33LKsKLHZiTAYbqL1AXd0WSOGiCLDcEH61jUTnN4q/ui3JN2KjTHI/lVUsFllZWaFerzM0NEQqlfpMM91eh8JMu93ett13I6SU+d5/V4UQPw9cANaTfBFYH7IfxxOO2BYvjeR9Ym+8iVJKZmZmWFlZ2XSE0Iui7+bFYrHBnu9Z2KoLbf357ty5MxDDLxaf2aL7BDZaiHw+z9zc3BPfd7O02Hy1xk98fIXry6tkIyF+7sZt3puZI6wquJbF77ckZ8dHycW2fmA2s4i1Wo2Pi3codNYQikbYZ+Jg07Hi5Ds+/JpJSFEBi4CqARbS1cmbJvnOJEdinkXvSpcI4LgKqiJxUWnYwyBqABhSx3VDhJQ6lgQ/gqIRpiGWqJo1Er74ju/hRqxPZbmuSyaTwbKsJ0Ql+t7Li5SmfhZs5/Y/b8WbECIMKFLKRu/f3wNsFGj8JeAHhRA/gxdwq+1kPw4v2V3faC1N06TT6dDpdAbR853gWW5Vv569VCpx5MgRSqVt22k3vbb151sv7hAIBCiXyy8k39yvi++3rT5L9PH9R3P89OXrPCpVGE/EKDSaTMRjqIrCxcVlRsIBfvbabX7p9j1GolFifh9fmhrndG6UkG9rS98v7viGOoeuBmg6NmFhIB2NmZafYC+4ZmOgAWovR163LKIqtB0/t+tZjkZXCPT25k3Lh+gZrqV2kHSwjqpIhIxi9YJ0PtWg0Elj4tKRKj+/8lX+7Pgf2dE93A6u6xIIBMhkMoOZbq1W66nS1H4Q71U3CG3n9m9XEbcJhoGf751TA/6VlPLXNrSZ/gpe+uwhXgrtz+305K/sbvQtYzAY5MCBAzsmeD8tttlN6u9xY7EY586do9FofKZW063EHbaaO74Z+uQ1DGNQeLNZXXz/OCklP3XxKr/zcIbVepN0OESh0WIqmcCwba7mVzgylKLaahPWNXRV41fuPuRQJsX9tTL/QlxnKhlnJBrl7T1jHMlmBqmtPu7W8yy4TYThEvWpOK6gacVwMdBVF9sRhHy9whwsBF4O3IOOLeBBI0fSV0YDDEcloHpeUMNxsDoxxsI1qpaNrrl07QiakFRNhYjfBqnxaf0Bf9wyX4p7vJFU66PafWXYvkjketXbfn/5y55J9qzv9CJFTFLKR8CpTV5f32Yqgf/yuU/OK4quT09PUy6XOXv2LHfu3MG27R0LD25F8n676aFDhwYthc/ToLLRkj5L3OF5022VSoV79+49s9ZeCEHLNPmHH17i4VqJlmHiAi3TZF8qQaXTYaZU5Y3hIZrdLhFdxZGSS/OLnB0bwXElrnTxqxr/9uY9Dg+luLy0jOO6HBxKsy+V4O3JMUZjUX5i+ncRikQo0HW62E4KW+ngwweY2I730NuOiq9n2TXV+77SUUGzaEuVRivNIX0VKVXoiUoI4WKj0zR9aKoDCMqGj6Cq4tcspKuiCUHLafG/5X+ds+L5Zr5vhu3c442lqf3A5HqRyH4A73WUFm+2bf1m4qWS3DAMrl69SiKR4Ny5cyiK8lxEhKeJ29/Tr62tPdVu+rxpsT62E3fYqWhE34I/ePBg2yGKhXaXf/pbH1A1LFRFsFRvMhINMxaLstJsUe10GE/E6No2yYCfQr3BfHmFtyZydHvFRCD4xkKeM2MjqEJQNwyCmsa/vXGH/akkvzM9R0u2mU/kUYWOL9TBL4KsdiWJECg9cvsUH2DgOjqoNq4EX4/kaj99ZimgKiw1hkn6bVxMpKvjxZskFTOLT/dSan7NptgVhH0gZBhb1BFC59PGXc7Izx5k7Su+7hTrA5PAoNW0X7cRjUY/U+ntdtb6VbQJfxa8VJLn83n27dv3hDX7LCQ3TZMbN24QDoc3bTd93nNLKbl169a24g47WTz6Qxds2+bs2bPPJPh707P81O2HWEIhFgpwZbHA0WyGmN9Hvt5ECImmquiqQjIYYLpYodzucjo3jOE46IpKtdvlYanCl/aMY9g2bcsmqGlcXV7hzdFhdFWl1GqzHM9jOzYCgS5htaOg9m5bx7EIAF3LIeAHDR3oYFoqmu49mDY2PkCVfhwMqlLFbgnSMWgaOn29qGJHkHSDBP0dXCdK15aEfW26bhefAOnE6EqL9+UDznFus9uyY3xWtdZgMMjIyAiFQoEzZ84MVG/7pbfrp9hsVyyz3bWYpvmZsjKvAi+V5Pv27XuKdBuFILZDn2A7aQZ5Hkvebrdpt9uMj49vK+6w3Z680+lw7do1crncYDb1ZpBS8pPfuMKHsws0LYugP8CnSwXOjecQAhZqdVKhIHOVGnuScaI+H/eLZZCSkXAQCQRUjZmKR/q+VXddF0e6fJov8vbkGKbj0DRNZMimrhgEpIauWTSaAfSAi+xNPvL1LtPf6zhpGTZxH1iOhqb3egt65axN0yWogevCqqMTNCJ0bNBV2TuXQtVUCfig1HHway6u40NTTVrdIdrSwZKST5ylz+wiv4wilv45xBaqt329uL7r39/Pb7zub7VecnjJJN/sh3xea6soCktLS9Tr9W2bQXZ67n6zit/vZ8+ePdse/yx3vV9Z1x/DVCqVtsx//+hvvke+1sByHFqWTb5d5Z0943Qth/lihYRfY6ZU4dBQGp+qcnt1jeFImHrXIKAIIj6dWytFFAFHshlMx8WnqRTabfL1JufHR2lbFrqq4NdUrjrLBDWJIyQBF9q2D50uWo+YpqVRLaXo2qC2Q2i2wHY1/IoLWDiOgq553zum+7GwUNwAimKT74SJYAEGUoKrOqgKdDpxUL2Gl3JbIxXUaTkCV7FxHZ2G6/Az8x/zJ/a8UC038PJ04raaaLpx6EO5XCafz3P37l1CodATDSg7EXH8PPWSw2voQlNVdccD7yzLolgsEolEdpRy264lc724w/nz57l06dKOrmMzD6GfaltbW3uism7T/Helxo/+5nuYjoNA8rBUBtfl9MgQDcNkpVJDE1DsmuyNRqjX6yy2u0zFo1Q6XUYiYRzT5NZKkaFwCFu6aIog7NNYqDewHZd96SR10yQVDNCxbR5Za2ghgW1o+MIWXUMftIhqikujGqOtAMJF90sQIFUfKw6EO36GNYGLAz7vtzJcEwWwTAV84AhJx0qiu6sI1we9NtU1QyHo8wJzfr9NtR3E0QxcS8OvKNiqwy8VrvPHJ99+YaK+DEu+06q5zUpv+/rv3W6XUCiEZVlYlrWpW/55tOQvN7ewCXZqbWu1Gp988gmxWIxcLrejH2S7vdHly5cBOHfu3HONFdrortu2zbVr1+h0Ok+Vzm5cED6YmedHfvW36Fo2uqJwaXGZVDDIgWScumFSqtXoOA5oGkeGh7BVlVXTZl86RdW0SOsqrWaTh9U6Y5EQDdMkFQoS8umstTqkAn6atoWmCkJ+ja7r0BYGFb2DKsCHhs/RMVXYG4+iINgT3oPPFwVFkvFHoJcu6+/Dx5MZHrWjHI0fGXyPaNC7//vjXl+1lNDExOwM4ViP72VQCeGYXtGLcALY0vO8dHw4qoHjupTNGn/vyi/TarVeOCj1unXi+p8ZDocZHx/n5MmTnDt3btAKfOPGDS5dusTDhw8pl8uDZ/xFSC6EmBBC/I4Q4o4Q4pYQ4r/a5JjvEkLUhBBXe//7Gzs9/yt31zVNwzTNLf9GSsn8/DzLy8ucPn2a1dXV53LvN8Oz9vM7cf3WE7ffiTY5OcnY2NPpoPWW/J994wq/Oz0LUiIUwUdzi5zKDRPSdQqVKsK2WTMtxhNxcrEo+VqDutElGQrStm32ZlLUuwar7RqjkRDLzRYjAR+tRpO6ECQiQWaqDQ5kktS6BvFQgI5jsuTUiAf87E3G0H0ufnQ+rs4S8AveEPv5sLTAG/Es1XadlC9Eza6T0EKUbU/xpWF3caRkoaEx7D9CU86zZnsqrlovapf1J1hwWqy6Bt+Z3Mv11l3AW8QbtIiZYWxbw1Q6+FwNR+1gWwrCCuJTdT5uzfPd9+4ibYdoNDoQy3hdQaqX4Q0oikI4HCYej3P48OEnSm+np6dZXl7mq1/96kAc5Tk+zwb+GynlFSFEFLgshPiqlPL2huOee6IpvCZLvpW7blkW165do9Vqcf78eUKh0HPv4dejv2DcvXuX06dPP0Xw51FhdV2X1dVVrl69yrFjxzYleP/YjmXx13/lt3j/0RyulNQMk0+XCrw7NUFA01iq1lAchyXD5PDwEKOxCNOlCrZ0UISCX9MYjUZYqTdZqTdJBANUTZOjo1mCkQjBcIhQwEeh3mBIV6g2G8R8KqoGwbjGoUwGv65xrbyCX2i0e/cvrIRYa3v3PqR6ZAqrnhVO+R7vG9eMJgCOhPdXSvidfQQU77hGj+xx7XFp7XSzy4HgFCElQN1tgYB0MIWjuSAkqh0CBIqRpC4dqpZNzZL8suZtdXK53GDxvHz5Mo8ePaJarb50sYz1eBXtqv3S20OHDnH+/HkuXLhAOBzm2rVrnDp1in/+z//5js4ppVzuC0BIKRvAHbwOs5eC17In34y0/QaOvXv3Mjo6+sTxhmE89+fsRNyh34m2k55iwzCYn5/fUra5j1LX5B/9+u8ihQISZktVDMfhy3snaZkmC6UKEU1lzbI5PpzBr+l8ulRgKpWg2GwzEg0T9vl4UCrjU1SEInCkSy4UomPaaIrAEiA0jeFkAkUV+G0TIVzKnSq26YBfYcIfZc2sM1trkY34SQs/nxSqHE8NQQckfaknz+sIqd6WI6GHKRmeRW87nsfVMDW61jDRQJGi6eXCvXQbKCjku1UWOy7fMbSHle4sAD7Fz97QOHdas1h6l73+Sa51i+AK/FJHCsHXV+e5W1nljfQI8XicvXv3PlW4EggEBjXpL6LAsxVedi/5ZhgZGeHs2bPkcjn+2l/7a7Tb7ec+vxBiCjgNfGOTt597oim8BpJvpqW+sLDA0tLSpg0rO20JXY+dijvsxJJblsX169eRUm6rNXdxYYkfu3YHn+4j6Fe5vFggFQpwKjfs5ayrVUKajqtqHIgHaHZNriytcnYyR63TZSwWRdcUPpxbYH8qRcsyyYRDRH0axXqDsKKw1u0Q9/s8PTVVYLg2mk+n5ZqsuQ5+RWVI99E0OuxTo9xr1FFVPyk7wrJTG4w/ajieRe4Tuf+90utIXjI9BRhNKMy2WowaCSbSCvPdFQzXO8+IP8Fs22tQKXd1sr4kq2YF2xFMt9ZI6BF8qorr+nBxOeRPM2NVsRxwbYe/fvW3+Pnv/pODe7ixo65fuNIPdPU76mzb/kw16a9LeKIfXe/v558HQogI8G+BvyrXDTfs4YUmmsJrSqH13fV+P7mqqly4cGHTm7Vdt9hG2LbN1atXOXHixLbiDtude32p67Py3wA/++lNfvXuA2zHJajDN+aXeGN4iHQwyEyxgmEaBDSddDRCxO/jZn6Fumnx9p5xTMchoKo0TJN7SyVOjAxhu5JcMIpfVZmpVIgogo8Wl3hjOMP9cpmj2Qym47BmeJVxPldlXyDJlcoyDdtgLBUiQgTNahEXIRqO5w2ttqsIJCvdGppQKJuea951vbx4sGfRA4pG2fQsj9lbCHX8XF5p8e7oXlaMBgAJ3ZNvBljpNug6Cll/gvvNIobrMBUaRRWC65U8I/4Y81aNgKpzLDZEpWMxW63zh37hX/HnDp7mK/umSK0rIhLrNNLHx8cHNelra2tcvXoVIcTAym+Ww34WXudE0+cV+QQQQuh4BP+XUsp/t/F9+YITTeE1RtcbjQaffPIJQ0NDHD9+fMvVcKd78v6gQsuyOHPmzI7UW55lyZeXl7lx4wYnT5585qBCy3H4O7/1Pr927yFNw8CRkjtrZd6eHCMR8HN5YQkci1AgwFgyQUjX+cbcEh3b5nwuS9e2Wao2aJkWs6UqX5ocw69qtAyDrmXz9blFkoEgD2pN9qeSzFRqjEajXMwvo2sqI6EwHy4uevpq0uJsZpRSt4NP6lg4TITi3C43kapAE4Kqa3BEH4ZmCKMaYWZVkHJGCRJEEwqip0Uw5H98/+qWZ/UDQsOSLh8t1wjKJFE1MFCWiaoBFjsV1swmMSWL0nuUSmaLrgkOkqjqpyttTsQmmK61WGh2OJLMokUk/+7hHX74N36bv//+h3z1wSO6m8Rt+oUpfr+fc+fOcfLkSUKhEPl8nosXL3Ljxg2WlpbodDrb/vYvUyfuZRfD9EQa/xlwR0r5j7Y4ZqR3HOI5JprCa9qTt1otbt68uaNBgDsheV/cIZPJPFeX0WYkd12X+/fvDyahPMslLLc7/OhvvEep3aZhGDQNE8O0ODWSod41uLNS5EAsTCAUQgJSuny8UGBfKkHar3O3WEECo7Eofk3jzdEs99fKlNodTo+NcG+1xFsTYxiWxZF4lHg4RCoUwJGSU9ksHy0scmZ8lLdHc3x9eZF3JsdoWzb7Yyku5VeYykYY9yd5VKtTlyYToTiKq2JJPwWzxuFoivvNEnP1Fg9qTYb9CWzdu9cxNQhUQUKh61lt2RMeyQUTXC6tkfKF2RPw4hNjgRQV0xuC2LEdsnqGJWOVjB7nUaNGWPVzv13khDbK11eWCKo6h+JpTNshEQwRDfpRWgrTlQrXCiv86v1p9iTivD0xxvmJ3KCzbn3KTdf1gXLqZmIZfSWZzcYZv0xLvp3m+gsUw7wL/GnghhDiau+1vwZMwmebaAqv2F23bXvQ7/uVr3xlR4GP7Ui+XtwhnU4/V1R2s373a9eukUqlOH369DPdv3urRf7H9z6m0u54RKg3GYvHCAcD3Fkr0zAtTqQT+IJBHpYqHMikmClVODcxiiYUPp5fIqRpnMyNUDO6FFstOpZN27T48t5xWobNcCTEbLnCarPNvniYj+eXeHNsmE+XVzg2kuH8yCgPihX2DMU5kRrik/k8+4YTDPsjzIoacS1Ay7YIahq2dMmKBB9WFzk/5AU2Iz7v4RyLxHnUqFIwTMxVhclwmk7Tc9XTepiV3h69I73fIa57f9eyLD4oFDgSHyOo9Nx8ofGgsUbXtTmTHOdmpYAhHU7GR+naXRZaXUzX4WhsmLv1VXR00noEC5dUIExAqEjho2p2KSw1eVgq829u3mZPLM7vPTDF0aH0pr/LZmIZ/aEI/fLU9cqwr2ui6YvMQZNSfsDm8k7rj3mhiabwCi15s9kc5JdbrdaOI5tbkXwzcQd4vj38epLXajVu3rz5ROvqVvjqvWn+96u3WKjWGIlGeNQrRQ3oGl+fWUAXgvOjQ5RMmxsLy5ybyGE7LlGfD9uRfJrPcySVIKRrfDy3SMSnsyeRIBBVMWyHj2fz2K7LidEs5XaXdyZzFGt1zu5N0tRWODsZ5+ZiheO5IepmF8VNEFR1NFUlrgVo2xYT0Shzq02GEgEmwlGqjTbdkHdfOo63/3Z7C39c9+5dUNEoWwalKnzHcI60sAjLx49EsZdH75NsMpTkQbPIrWqJsWCcN0JT6KrLpcqSd35X5Ug0x7X6Ao9aZcbVJMtWkTPJHNdreUYCMUb8MZAKhuXSkAYhVUf3KSiKQFUEy+0WSb+P22smV5YLxH0+hlwb39gKJ0ayWy7E60kNTyvDAgNL/1mGPuxkT/55q3h76SQXQrC4uMjc3NxAjXV+fn7Hf78ZybcSd+gf/7yWfHFxkYWFhW1r46WU/OTHl/ng0QIP14oczGaodLrsy6QQwNdnFhgO6IzFYlxbLWE7Lu9OTVDtdrm9vMaZ8VEM22YqlUAV8OlykcPZNEORMF9/NE/EpxP2+TiSSRMP+rm4kEcIwQPlIXJPnU7EQiiSjhTEwj4WzDrv5A7x9fklDuVSHEmmuTS/zFg2yp5QnJnyIkbb5Ewix32jRszygm+LvWj4mtHXT/fu32gwxrRVAWC20aJsaBwazzJTbhBXA5Rt7+8L9SoAwZ48zLA/wlKnxlKnxunEGCejEzxsFbhfL1K3Dc4mJ3Cl5EGtQlTRud8okvKF2BvM8FEhz5AvTMoXJqhquAo4ioOuqShCJSOCOK5kplLjQCqJJV1+bXmV2Q8vkg0FSQaDfGX/Hk7nRtHVrcm2sTz13r17uK77ROdZ37V/ntTat33tupSSmzdvfqY53RtJ+yxxh/7xzyMcMTs7i67rW0b3++g6Dn/rq7/L1aUVVhpNJlMJXCkJaBrNrsFsucrhWBhd93FlpUQmHOT4SJaPZhYwHZcv75tktdni/mqJcxM5Gq02Q6EAUZ+frz+c5+BQiqFIiA9nF2mbFmpVQdck8TdXkXoXIaBr6WiKg0938EUsDLvBtHsTXUtiOA4xzU9I14lqPtqmTS4aYTQapmGYaALmWlVywQjLZp2wqpNvewHaeo/88Z77riJYbjexpctazWVSGSUdVCk3lhBSUO61sbWMXmGMq7ECRFU/16sFbOnybmaKuXaZum3Qth1UVOq2wYFAghWnzf7AML+9NM/JxDCGYzPTLBMQOlElQFjoxFQVoYBUJAFVZW8qjuk6zFfrTEZCBDWVi4vLhH0+HlWq/P/kNfYnExwbzfLu1DjxZ+i8CeGNxe5b+n7nWX/oQ/+9dDpNOBx+5rZtu8Bbs9l8bhHHV42XviefmJh4rukYG7He/d5O3KF//E4sebfbJZ/Pk06nOX78+DOvb6FS46cfLjLT6hLUdcbiUYK6TqXdpmM7SNthKhLCUlUelKq8MZQiEgjwuw/nyIRDHB/N8rUHswB8Zf8ebi6vUmp1OJSKIyT4NRVdCD6aXmQqHmE4HuWT0jzJU2Us6aILcFxAkdgoSFtiuyoIsDWD+KEqlUca/oTKZCJOudQmELEZCYW5urBCJhlixBfAF/IzGYjhlz5GgmFsxWahW2Gx5Vn2fmQ9F4qy0PJc2rppMN2o8qYcYthKMpWJs2zWWe7WKUlvcVD8OlgwqgS4Z9cJC41rlTymdDifGOd+vULDMtgXSJA3m5yMTfD+yjzHk8Pk2w06lsXheAbHlbQtB8OxaTgGmlQQUtA0LVpdk/FolMl4lHK9SbVrkItFcKWX4TBth68+nOF6YZVfu/sQXVU5PpLlO/ZOcCDz9LDB9YG3jZ1nhmEMpKP6JO0vCBsLoV6FUuurxkt31xOJxJaCiTshfr8WfCfiDrAzkveDdUNDQySTyWdex4czC/zUJ5/y6VqFoyNDveuBpmFQ6RikdRVHVak7EsV1OJyKY1oO14ornJvIgYSv3Z/l6HAGv6by8aNFfJrKqeEMM6UaD4o1JuMRrI5nFW0puTSbJ3q+jiEVpKsQtiW2raD2dNg6XT+q7rV2uraK8Hexx1eZXXY5mx5hzqjS0myGA2Hi/gDFTos3fDGur9SIZgPMVKuk9SCX1laYCseYiCe53ihQ7/UUpPwhFloNhBQstT2y265kvt0h1YrxoNxiPJyg4XQZ9iUpO0WGwxYNpcS438+QTFI1TQpGg4VqlSE9QMsyqdsWk2qci8Ul9oQTFDp1TMflZHoEx5J0DIuY6kdTBUgFVQokoEqB6hd0HYeZSpUxf4Co38+DUpk98Ri6qnJrZY2pZIJ0KMi9Ypm1ZgvDtvlwbgFNCA4NpTk5Osxbk2NE/L5nktPv9zM6OjoQiWw2m5RKJW7evInrugPXPh6P74jkL7NS72XgtUw17RNxJ3uf5xF3gGe76/1a9kKhwNmzZ1lZWXnmgvCvLl/nf792i+lihaOJKD5VZaXZYjgSodxqk1IVXKFSNQz2ppMsVxt0LRcFSIeCSFeyUm+xN52g0u6iSjg6nOHm8iqNjsloJIRdb+JDMF1ucmp8GAOD1ckajl/i2gJFg5btQ8FF7VWr2VLFtFQifhMUh3rbh1QFoWGDWstgIhkn7NdZ67SYiMcYdSOUW17Ry3zTs9o107PCqWCIi4sr7IknBoE4vWfhRkNh8i1v314xvUVI7fWrDgXDNO08bshGETZCqtjSRQqDOdlE0fwMqXFSeow7tVX26lFqXZNFWoRVHdN1cF3JmVSO9/ILHIllkMDdSom47mcoEMYnVCzHoWGYCBdy/hh7E3E6nS7XV1Y5nEkT0jU+mltiTyLOcCTMw1KZaqfDufGcJ55heP3uv3DrPvdLZX7h1j10VWFSU/jOUITT22i2CyEG88unpqaeakLpdDrk83kymQyhUOi5RSW2+Mz/CG+YoQr8pJTy7254/4WHHcJrInm/tHW7L/+84g6wtSV3HGdQXdcP1m11bNu0+CcffINfvn2frmlzYXKMar3O/bUSh7MZWp0urmURjMWodg0mEnFmS1Wy4TCqdLm1Wub8RA5HSpaqdRKhAMeyQ3z4aIHlWpPvOjjFe/dncYHj2RRN27sGw7GYHp1D61lsVwgUJFJ6+3FVSBRFIhVAgmEqWI6G42joPoem0+a+ssQZfS/Xl9fo+m3OBcPcqhRpYXAwkeZBo0RM9/Oo7gXYBl/fgaXVNmdyOVq2F33PBDyS+xWVQrdHdqMLuKypD/FrJlIKHKmgqS6a1BkJJChaNXK+KK5UCEo4reRo2DaGbeKXKppQ0E2HBH5+Z2mek4khFttNcOFoOkOp08F2HKI+nYjfh44nZdVxLJYaDeKqyrFsb1FYK7E3lWAoHGKxVqdrOxwdzuAiaZk2bdOmaVmcHBnCdF1urK4ylYpzv97gRqNN5uY9ArrOW5M5LkyMkQk/2+puHGX88ccfo2kaMzMzA9e879q/SAyqZ6D+CfD78QYoXBRC/NKGDrQXHnYIryi6vhH90tatGj1eVNyhf+6NlrzdbnPt2rWnatlVVX2q7XWuXOUfv/cNfvvhI0aiUfanUzRNk3yzzUQyQbfTZaFS48hIFtuVmI5Dud1hbzLBbLlKqdXh/MgQd5bXaJkWZyZH6XQtPny4wOnJUUzb5r07s4xFAzQch+Vah6Zhcnp8hBupaaTPAa3fPOJdk+0oCEXQMnT8mgPCe8+wNRypIFSJYysoisQJmTzorLInFScU0HlU9aLSlWad8VCUuOYn7Nf5oLCAAGYbnmVPB0LM1RpcmV/j3OgwaX8QvRc9n4zEeVAv4xcKc80iyVQLdAthK2iag3QF4/4RWm6bJaPI8egk0tVAaqgIfKpAN3z4nTim7RDwBamaBveqZU7FM9yplUgqOllfkG7HIIoPV8LNtTXSgRApf5COY2M4NnuTcYq1BqvtNhGfj/F4DMeV1LoGAsiEgqiK4pX89ryQ48MZGobF3dUibwxnaBgmUV0lqGn85sNZkgE/q80WP3/zHqlAgPFkjLcmxzmd87Tytnve+jPqpJSDGXWLi4v80A/9EIZh8P777/P222/vSMPgk08+AXjYk2WmNzzhDwPrSf7Cww7hNZS1wrNd6q3EHXYqLrDROheLRT799FPeeOONp5pVNh773vQs/+0vf5XfuD/NqbERRmMRCvUmpmWTCvixDIP5Wp0D2SEUoXCnsMZwJMJYLMrlhWUsx+XtsREuL62RCAbIxSNeV5gLYZ/OYqmGaruoAhYbXfbFYkgpMW2Hucgqjt/CtFQsQ8e2xED73O2pqgpFoW2sG/znqEhHRQhwbQVVkwgkFbVORTRZq7VRENi2ZK1t8o18nsuFAsv1JlktzDvZMaq9CHn/NoyFI1xeWqFdcYkKPyPBCLFe1H00ECCcaKFp0DV0XKmgOGHeiEwy113FlTDsS9AwHWqWyVynyO3GMkWzjqvaCBws4TLbqtJ2LE6ks9xsFDmdGWUynmau1QFX8qBaYq1RYzQQpNRpoQnB/mSSlUbLKwNudxiNRtBVlZVWC01VvCyHZREN+FGE4GGxjCIEx7IZyu0uNwurnBzN0rUdYn4fmhB8tJBnNBrh4FAaw7apdzpcK6zya/em+ZdXbvBXfuHX+X/8h9/iF27eZa5S2/bZE0IMuunOnj3Lz/zMzyCl5Gd+5mf4gR/4gR09v0tLS7D5IMP12GrY4Y7wWt31jdhK3OFZAxY2om+d+5NVyuXyloMP+5F7V0r+14tX+RefXKXS7vClqQlMx+F+scyBdIqWYdCxHQxNMJ5MIqXk0kKed6cmWG20uLG0whsjGWK6n0tzS7wxlOTOWoVzk2MU6y1mS1UOZlO4ps3NfIlDI2mk63BtscRoIooYcVnRKkjb+34G4JgqQdXFFzKRCPqz7GxbQ7MtVE1imwo2Cn7NReguTlfDMFWQsCxqZJQQCSfIfK3KVDDIkmsxEY1xv+KNPRoPxTifynG5vMxc3XuIs6Ew+UYT03a4tFCgbVlMjibYo8Zo+YtoCqiqi20LVNUl4fNTMg0yeoyiVedEZD936nk0RUUIhYlwgqV2g4Cw8Cs6EREgE4rTNi0UVN4d2sPllWWygTC5SIx71Qonh4apdbtEhM6EK7haKDAZDIOQdEyTg4k4Lcvi7lqJUyNZVATXl1c5PJQmoGp8PL9IOhzi6FCG5UaT6WKZM+OjdG0bv+aNGL5dbbI/nSATDtE0TLqWTbljEA/4GY/HKLXa3C+WOT02ws9cvc2v3JsmHQzwPYf383v2Tw2eoWcZn0gkQjwe58d//Me3fW63Od/GF1942CG8Aku+EzHHlyXu0D+234lmWRZnz57dsqJJURRqXYO/9Rvv8WPvfwMQnJsco2WaPCyW2JtMYJgmpUaDiN9H0OfDtB1my1W+NDXBvZUij4oVLuzJoTiCuVKNbCTEnbUK37V/ikuPlig22uzNJBCGhWHahHw6C6UaIUVHFTBvVVgeriNVzwXvp7FwFTqORrelD34UKb154obhLQaOEKCAtAVIaHd0hCqQUiAdhaJioMR0DsfjzHRaOKaNz3IZCgTxK16d+OV8gfPJUZz+/e09KhPROO3eGOiH1QrLrGJqBq6EdkcjqoXQUEnqMeY7JSxXcCK8l4vVOSbDQ9TtDiP+GAoqU/5hIkrEm2+uQNc16Lg2S+060/UKb6Qy1C2Dqtnl/PAo0pHEND8L7TYlx+ZAJk3FsdkXibLUaNI2ujwoljg5lMZxXC4uLnMwkyLs8/HxwiJD4RBHhtIs1RvMliq8OTaC5bgIwHZcriwV2BsJkgmHKbU6mLZDqdMlGwkzEY+xWKtzd63E6bEROpZNNOBDAXRV5Tv3Tu7oOYQXq3breZvbDTJ84WGH8Brd9fXtpjdu3KBerw/UYDY7fqcFLqZpsri4yMjICEeOHHlmyeH9UoV/+Ml1fu7aLQ4NZZhKxal3DRpdg1QohG2ZrNUbZOMxXKDRNfGpKm9kh7i+WMCwHS7sybFYrKMqCtXeA3Mik+Jrd2c5MzlK0zDx2w5N0yFfazORjHN0OMONxTVyQyHsAzaucBGqF2Dr663RmzratXWk7X0H11JACGxXwe6ooAjPVTc0rK7miTG6oEkFXRMcSKXJJkLoER8HIlHeGMpyu16j1bU4EohS7Ylx2JZLVPqZisRZ6pV8ZoLe7zAWiVKxG+gRB8sSoCqENB91t8Oo34tLALRtk5JpEVB0bjfyHI9NIF0Nw1K43Vhhvl3GRdKVNhJI+4NMRGLE/T4apsmJ1BBjwSj3SxVAcL24ykgv8h3XA6QCIZbaHRKhAI6icDCVpNU1uFZYYX80jOo6fDK/xGgvjrJUb7DSbHJ0ZAgpwXFdWqZn/c+OjxJQVfK1BpqiUO50mYhHyYZD3F0rMlepcW58FMvxBDNtx2UyEef/+T1fQV33PG3X5NIfxPg8OH/+PMBBIcReIYQP+ON4ww3X45eAPyM8vM1zDDuE17wnb7VaXLx4kVQq9VLaTVdWVpiZmSGdTj+hLrMZ/v2te/zIb37AzXKVd6bGCegaq80WYZ8Pw3HwuQ5d08JVVXyaTs20GI6ECPt8XFsqEA34eXNshCuzy4wlYp7LnsuS9oeYLtaYyiS4MrfMyXSU+6sNHBdS4SA+VUU64KgOd1JlpC5Rel9NOGLgiA38H1fQbutIF4Tj3R8hBO32462LrUhvAVAAQ8HVXUJ+DV0VlI0OpnSwcOm4Nm9kMuxPprjfaHE6k0UBZssVCs0mTtNkTA+jIjAcbxFOhvxoSQOJQPdJhKvgal6wMqQGuFFf4XhkD3tDo9xvrjIVHPau0dVYaDWZaZUQCPaEEyxZdfJmg5rToWg2aVgGSOjYFnfK3hbCr2vM1WucHh5BFwp7IjE+LRQI6hoVo0suHKZimCiaxophcH48RzQS5l65SjboJ4rLbKlErd1hPBZDUxQsx6HUk62+MOFp3C+2O8QD/l52JEY8EODKUoFSu8O58VEcKbEdh67tcGI0y1/77i8/5ZW+ipLW3pb0B4Ffx5N9+jdSyltCiP+iP/AQb9jhI7xhh/8L8Jef6zOe64p2gK3c9b54/fHjx7ft/d6O5FJKHjx4QKPR4NixYywvb72oGbbNT3x4iZ+6+CkBTeN4yov4LtcbpMMhSq0WIeni0wPMlmu8MTKEdFyE9DzZO4VVjuWGUSR89GiRN8eHuTyX5ysH9nB1dplU2Ku1LjfbHE1GuLXc4NTUCNfmCpycHGFxtUq51cY86SIQuC0BKjgdgeaC0B2kC67i7cQVW8EWArOte8TvP1OOBq7tEdsROIaKCLkoQiAFHE5kMG0HC68RxZSSoNDwawq1jjc3/urKKt81OcnXer0E2XCUaytFMj4fNA0OhGLMyVV0PwhX4HZ1zg+PU3da5Lsl7jdXAbjbXOVgyCsxvlnP807qAO+tznE0Osydxgon4znuNlY46ssihYaJp2OnodLtSTlLabHUrrM3kkRFYbXRIqjr3CyucmwoS7nV5szICJeXlhkJhWkaBpPxOJbjMlOtc2AojSYUKp0OKg5hRcExulS7HYqGhRCCs2MjtCyLe2slcgE/DdMkF4sQ1DU+nFvAp6qcH89huy6NroHpuHzv4X38xbfPbvosvarBClLKX8Ej8vrXXsqwQ3gNlrwviFipVDh//vyOxB12Eo0XQnDmzBl8Pt/WQhD1Bv/1L/46P/mNK+xLJzmUSdG1HepdAwkYloXPtjEVlYeVOidzwwRUjZsra6RCAWqdDkeGhzAtm08XlnlzYoRys8OX907y/v05prJJ5so1hsN+cgE/dVMQCfi4s7jGkdwQtVqbyXSc5gEH05Y4ikD299CuimlqWA0dxVIG6TPZ+yqmpWGv+1rSErjt3gNmKAjX++kc3eV4bJhyw+LOcpWQ9PaTLWlSNtrM1GtUDYOD6RSZQJAb+TXOZUdRpGCl7XWZjcQT3Kk0WbCqtOliGyqOFNiqZK5d41atxIHQJKN+b6Dg/tAIl8t5DodH8SkqM406GV+YO40VjsVGWDOaHA6MMdNp8WltDVu6WNKlYVtYjktA0Un5ggSFj5lKFdOxqVpdGobJ/lQKv1DwqRrTlSphn0bCpw9+s+V6g/3JBAKodQ1Cug9VVRmKxYhGIyz3tlh7IyEWikWu5VfYn4jTdVzSwSB+VeXrs4vEAwHezI1gOA7lVoeObfN/PHFkS4LDzjrQPm/NKfCKSd7tdrl06RI+n4/h4eEdy+8+S/zx0qVLTE5OcvDgQYQQW3ahXVrI8wP/+ud579EcX9ozQSIQoNLpogooNJokdB2fY7PQMSi2O5ybyNEyTD6eW+Tk6AgBRaXU9iSgbhfWeGffBLblEPP7uTSzxImJEW4urvBGJobZdQmGwixXG4ylYigCQopGtd3l6/YibhBQenlwh8d3XYDrKphNdRAAG3wTIXC7jx0tqYBj9fbqrsDVJarlkT4uwtyvVHgjnaVpWTS6NlF0JmMx9icTJEMBFCk4mEhR6XS5slTg3dwYHdNz0RUhsIWFHbNxXQVFCIQmGVGC5LteFD7fqHG3WuZ4ZJLZZhUXyVyrxonwJLOtGkndq9c2HAfD0Fgx2pTsDseiaQrdBvfrRdqOiVAkqirQdIVY0IeiCpaaTQ4nU2SDQYKKyu1SCZ+q0LEtDkRj3CyWSQaDmLZFMhjEsG1apkU84KfeNUiHgijAteVVhsIhjg4PUXNhrtXlzdEshm2jS0m72eCj+SXGomEOD6VomCalVpuu4/Bnz57iT5058czncjt3/UX25K8Dryy6Xi6XuXz5Mvv37yeXyz2XOONmJM/n89y6dYtTp049EY3frJ/8Z67c4D//2V+k1uny7t4Juo7NYq1OxO+j0DGYjATRheRWpUE8GOD02ChLtQZ3VtZ4d+8EILm+VmJP0iPssZEs+WqDUqtDqd1mMpvAsWyOJqLcLdRJhv1cmytwZm+OpXKd81NjXJlbZHXSwIkDwqvHBgYBNrmO7NJWkQ3v4RnkRRyQhuotCjZI4XkBoq0ge394IJJmMhTnG0vLxH1+8u0Gd6trWK6DISUrjTaVVhc/KtKFTxbznB4dJqiqzBRrKA68PZqj2umipiTSESg+yZ6YNywh5feCcSO+KAtmHRdJsdIg5XodX46UVA0HIQV3G6ucio/TMiAXjLHQrnEkmOJ2q0hA1TiTzpHyBdFQ8SsqIUUn6Q+wN5Ygouks1VsYrsudUomDqSQKgkOpNDdLFYbCIWJ+nWK7iwRMxyGs+1hrthkKh1CFwqWlZfamEkwl49wrlpipVLkwkUMKgabrCFXlXr3FwWSctE9nfq3EaqWKYdn8hfOn+E9OPB4ssRW+FdtM4RXNJ5+ZmXlC3KFarT73ZNP+ouC6Lvfu3cMwDM6fP/9U7nx9uq1tmvytr77Hv799n4OZFEPhMOVWh7ZpoioKLcNk2KdhuJKHlRqnx0a8+vRGE9txOJ4bpm1ZLFTrnBzJgBTcWl7l6EgW2ZXsyyRodEwcy2GlWsN0BMcnstyYX+XUnhE6hsWhTIqvLcxiHRCYwkZBIA1A7W3ye9F0YQsG+gwSXEtFM9wB8YUtkEJAS0XR101zaau4vQxhvuNZwDmnzYloloulPG+mhzFsh4ZwwLFwHMmq2Sai6ryZG+HW8ioXcjk+zRdoWzaW7TDTLqGkPdnmsNDJ+KJ0TZey4bDfP+rl0jtNEOALhrjTWONwMAW2y43aCsfDGW5311BcjXK3xkKrwWQozqNOjUPBFH4lxO8uznM4liGg6iy1mxiGQ0TzkfYF8esaa+06EZ/OqeFhLMsb4FjvmmiqwmQkwpWlFaZSCXRVoWnbtDsdhiIhFCG4vLTMkaE0Q+Gw971MiwsTORxX0ugaqIrCQk9iK6CpTJer6D4/fin5Px3Zx6RjcPHiRRKJBOl0mng8vimZX4UqzOvASye5YRhYlvWEuMPzzEODx9bZMAyuXbvG0NAQR44ceWYOfr5S5b/6hV9jrlLl3SlPtGC6XGY4EqFpWoxEwmB0meuatOwOX9rrFcAUWy2SwSCLtQa6onjNGgLaves9MTaMKhTajlcTXW42aXZMJjNJam2Detfg1GgC07DpmBbXfBWMIQkuqI6C6wNhe/luTImmgq266/xyz7gLIXBqGsR7hO6tia6pInAek99RUBxwVUnNMrA73vGPalV0oeBISQgNXfhJx+MgBdIF03JQHDg7OsKHc4tkQiFO50a4Wy4iYhJHcTkSz6CqKh8Xl5gIRFkw6mA0qJgmewLDaKrDvcYaAMtWl5FADMwmd9sVDvgSfGN1iUOBJPdsA5/QmNAjSEfjYjnP6fQI+XaDeqPGoXgaxw8rjSY+S2EkGCGeyZKvNQgqGndKRY6kMzS7JnujUW6slRiJRUiHQjxYKxH26YxEwnQtm/vFCidGskT9Pj6eX0RXVM720mGlZgtN0yi12hyMRdAVhWv5VbKRELqq8pffOcd37PPy4P0e82KxyMOHD/H7/YMe82AwOBi4sV0KLZf77PPYXzZeOskDgQCHDh164rUXmVFer9eZn58faLltBUVRuF6s8F9/8m8J+3S+tGeCUrvD/dUiJ8dGWG20GI+EsU2DO7UmYV3j7EQOw+51LAnBdKXK0WwGAVxZWGYqnSAoBHXDYLnZZCIe5+BwimqlTkATJLMpFFeQDgVZqza5X2oihn2UAh0c3dvjuipId91GWwUhBY4tEE2Q6sZNOOAoKC0XNywRtrcPl0J4bnvwsQdweCjNndYaw8Ew7Y53X8tGh/PZEaaLVfZE4xRMg/nlVaKaj5juJ6r7WGm0aRsmBzMpoj4fl/JLtMIWfl3jjWSay9UCF7JeKjKp+1k0YTKUYK7Xf/5WaoKJoMt8p8q+cIbFZpOwqtNyLOKBBLRa3O1WmfJHsDoGiqVw1y5yNJ7mYaOMLlSOp4cQrtc3vi+SJN9ocqtZZG88gYNksdHgeDYLrkRTFBYbTXy6yr5kgo8XlkiHQuRiUVYaLVYaLU6NDqEpKl+fXSAdDrE3maBr2RTbbVShYLsOB9NJyvUGFxfyHBhKE1BV/uI7Z3h78sm+hvU95pvpvyuK8sy4Urvd/lzuyV9LCu15ZpRLKalUKhSLRd566y2C63S5Nzv2f/n4Mv/zvTkODw+RDAW4ubyK6TgcHs5gWDZRBeqdNjONNsdHsnTbLRzXRUGwUKsT0nWOj2RpGiYLlRrHclnqXQNHQMzvo2zatC0Do9MmEPCRisRRETi2y2ylyoreopVxcWQH6QNFgqMCrkRu8OwUS+DqeL3TbYntc1F6+20AxQHRVnHDthdl71vvthiQfCgUptsj9ngwztXVVaYycWabNfzo1LsmC7JO07Y4nslSbhuUjQ5d2yYbDeMEPHWa2UoVf1SDgOBIJo0rvNVmpTcyqWR5wo5D/siA5HPNGk3bZCKYYK3TodBtcjI5ghQOH60ucS6T43JliZA/xFK9hemYRBWdQquB5VgcDsW4vVJifyxJ1TBYa7U5ksygKV2mq2VOpodpmxZt02S50SQXieFTFPZFo1zJF4j5/RxIJbizWqJtWZzOjdC1LS4uLrE3GWcoHKbW9bwrRQjCPp1EIMBKs8lMvcmp0WH8msqff+s0Z8aeXVcRDAYHjSh9/fe5ubmBQmw6nSaVSg0GKYBnyV/UXRdC/H3g+wETmAb+nJSyuslxs0CDXrRGSnluu3O/ktr1jeN8d+quO47D7du36Xa7jIyMPJPgLdPkR37lt7lfLLE3GiKka7w3PcdkMk7Y58OnKDQbTRbbXdq2zZf3TmDY3lyuhKJweSHP4WyaTDhMud2ha1qEAz7yjQb70yka7TbVrkE2Gub+yhr7Mimajovu2lQMg1jQTygRIGsIFqo1cME1QUiJ5gikCejg6BLZW/wHt0SCRKA2FVDXLQaOQGoCtSVw1y+WjkAxwfV593K+UufARGrg0qe1MGt6myuLBY4l01yvFDmVGWKl02Z/JEHXZ9MyvCF8VcMgqKuMZMNEVR9hn8bXVxYZjgUZD0XJalESoRBL7ToTagofOmFVZyQQ5WGvXXVPKMGq7XV8Xa8UeDvlubzXyytMhuOYlsv+aIpLpTzH40M8aJZ5O72H312a52g4zoNqiaCiMRmJMF+vcjyZJab6ubXmWfT75TInh4cRLoQ1jdlaHU1ReCOb4Wp+BceVnM2NstZqc2+txInhLCGfRqHuLVCq4nlZflVlulxhrdXmeCpO1O/jz791muMjzzf8oK//3mg0GB4eJpVKUS6XmZ+fHyjJNBoNGo3GZ9mTfxX4YSmlLYT4e8APAz+0xbG/ZydDFQbX/6JX9DzYSS16p9Ph4sWLxONx9u3b98xGgLlKlf/8Z3+JYqvNaDRCsWvw8fwS5yZyvXSKpNFo8KDRJODTe3PJLAr1BvRTYnvGSYVCLNXqaALyTe8BmUolaJsmlu0gXYfr+QLj6SSfrq4RCfr4eClPyK/z8eIiPk1hulFlIh7y3HTA1QSuAqgCqQjUtkDr9Ajbu9uK1fN4bIHoPv39RFPpB+HB9Sy9r6WiILzvAAS6UKx71vbm6ionE8N0bYe5So2QqmFLGPFF+Ggpjyug3O0Q0FVyiTAiIGi6FiqCqmUwEYkyFoyR88W5uLqMqigUzS4tx+b95QVU20cuEH/8e0oVXWqEVZ3j8WE+XVshF4hiug7j/gQPq1WulgqM6gHu1Uucjef43aV53khmmO22SAfDZMMRAkJlWPHx8dISzXYHn6ZQaDa4kMvR6JgUmi0M10FTVY5mM9xcWcNyXc6OjTBfrXNvrcRbEzkCusZ0qYJP00AIspEwuqJwo7BGpdPl7fERVEXhB7907rkJvh79wFtfSebYsWNcuHCB8fFxLl++zO/8zu/wV/7KX+G//+//eyqVynOdW0r5G1LKviX8GK8+/aXgtZB8O3WXUqnElStXOHz4MJOTk8907997NMcP//JvEtA0ArrK5cVlKobFl/dOIiWs1uq0222mG21O5kY4OpxhsVan0u5gOi6m63JuIkfbtrhZWCEVDnK/VOFoNkMuFqVpmCgIHlSqrHa6ZOMxHpYrnB3Psdxo8qU949S6Bu+Mj2PaDueGRwmgcCw95D1AAlSLQeRcl4rXqNEQA4st7Mf3Rekqj/Nm/ZQaCqrZc+Ft7zif0BmPxHB7BbD5WodyyyO55boYVe/fddvmZGaY5XKTm6trvJHJUOx0iIWC1G2ThrRpuTZB3Ss2eVgvsyeY4NP86kC6uT8madjneVJVw+BGsciZ+Dgxzc/t6hoLrTp7ginahk3HsfErPob8Ia6srnImPYotXZCCo+E0d8ol4j4/hXYLv6oxEYlRaRlIReNOq8GJ4SyL3TYjuh/NdlgslVms10gFA4RUjUwoyKNyFcO2OT82ynSpwnytxrtTE1iuy53VNUaingX1FnnBteVVfJrKO3vGsByXP3N0PweHto7t7ASbVbwJIYjFYvzlv/yXOXnyJD/5kz/JiRMndtRL/gz8n4Ff3eI9CfyGEOKyEOIv7ORkr8Vd3wp9LfW1tbUn2kM3K3Dx5JGv8OHsArbr4kPy8VyefakEQdeLkhfrTW8KqKLyztQ4Xdvh0kKew9kMhXqT0VgEp93m7lqRtmkxmYxTanU4Mz6K5TjMVauMRiN8spQn4/cR9OuE/T7OpRK0LYvxaIy5Sp1iu0XbNlmqNziXy3G71iBpWLgSzg2NslJrsGB4yidWL8Cm2gLRlFgbNP6EK9CaEjsEKI+JrXYETkB6CwbelJLTgSRLeG5yJhQh6Q9waW0ZBZgrNYn4VJquQ6fVZSwc5lap5O1tw1F8moquKeSbDQJ+Dcty2ZOIowpBxTBQBMw0q6hCMN3w6spbPdIfjqe5WytRMjp8ZXgPH6x6ZbGmLYlrIaDKo0aFr2T38LXGPFdLK0yFEwRsiYZOxehyPDXErfIa3zmyhw/mFziWyXJjbZWzwyMs1ut8aXSMD+YXeXN4mKuFFU6lUtSbLQzLxjFMXODEcJaHpTKFRosv7xn3ptasFjk5msWwHYK65gVP8wXGYlH2puK0TIu/+vZpfObzT8rdiJ0Uw4yNjXHu3Obb5N/3+34fhULhqddv3br1h6WUvwgghPjrePPK/+UWH/OulDIvhMgCXxVC3JVSvves634t/eSbwbZtbt68OZhztT41sbHApW1a/L9/+wMelSusNb0a57urRd7ZM4EQcG0hTypog4ChuNd88LBYotBocWEiR8M0mUzG0RSFj6oNRqMRxhMx/JpKMhhgoVpjsVrn+GiWm4UVTmbSqLpKu2uQCAZ4f2aOXDxKvWPQME3Oj+domiZ7xhM8rHjjqEb8Ae5Ua3y6VGAyGuNAOMl8pToguXAATUFtu48tt/TSZ0pXQdEfL2qqA9IRCEcSQKPT+wO5bhJuwudnrewRfiwUZrnV5szwKBdLefK1Nh3LYiIQJKBqXF9Z5chQmjuFIt+xZ4LpSoW82eSt0TEuFfIoftgbT7BiNnkrPU7Xseg02yhCJ6CoBNXHEeWZWo1jsSw36yvoqFxeXeZEJkvdNvhgaZE3kkPcrq0R1wPcqxQxZYvDyTQ3y2v83tEpfmdunnPDo1wqLHN2eBTXkUxFErw3v8DZkVGuFJb5zolJbq8UCek6uJKo7kMH7q+u0bRs3sxmWGu2eVAqc3ZsFNuVWK5LRNcHfeaZcIhat8vf+4Pfjeh2qPckrj4LPqtS62/+5m9u9Vaf4D8AfB/w3VuNQJJS5nv/XRVC/DxwAXgmyV+Lu95H/7r73WjZbJajR48+lXtcn3Kbr9T4G7/223yysESja2L0Xj8zNkqx1eLDWW/AgapqTKXTxAJ+vj4zT61jcH4ihwtee2Grw0dzi+yLhhiJRwjpOgFN46O5RVYaTd7IZqjXG5zIDiE1lWvLqygCvj67wEgkQjYU5vBQmgvjOS4u5rmzWqTeMag0DM6P5Kj2erEnggHmq3VmSlXennws3qH0fjLFEo+98/7eXBHozcdbGmmDUATBlsKBXGbwut1xiPd75R1JvtLgYDTOcMTrB7iWL3BhJEela9B1XPx6gLDqtZ1qhsF+f5jLS8sMR8KcTY8wXamwN5HAdl0mggnctoLRdblSWEPVAtwqV/BZPvxo+ITCVDjBXLPGtdIq5xJj3C4XkcDDapWMFsaRkkKryXAgTKVlMOmP4CJZa3d4e2ic9+cWOJxMc2llmbPDI+CAhsIny8ucGR3hfqnIu7lxriwV0FSFoK4RVTWCuo+yZSNUjZMjWUqdLg9KZQ5GwzTbbQr1OmFd42GpwrHhITLhINVOh3/0/d/DUCT8UmeTPytP3ul0nhksfhZ6Yo4/BPwhKdcv508cExZCRPv/Br4HuLnduV8byfsu+OrqKlevXuXYsWNbFg70Sf71mXn+u1/7bX7z/iMy4RBd22YyEWc4EuHjuUVmShWOxCOoqkbY78dwbL4+s8CBTIqTOa/yq9k1WazWWW40+cq+STQhWG00EQI+mltkMhHjaCZFpyele7tU5kZhlTdHvb8/PzbKcMSrpFpuNLm1vMbRTIYvT05wa2UN2/UKW1odm5NDWZLRxw04j5aK7PcHEP20Gp5FD3a8266sMy5K+zHJtV7UbTTi6Zn10epYHIp5muKlijcoISB9OL1OFseVBB0VrRcDSQQC3Fgp4aBwq9IkGArStm3MZpvZ1RKr7RYxXWN/IMVqs03bslnuTUvRdM/Jy0VifLycZ1iLMRZ4/N0sW3I66XWijYeiLNdbRDTPNT8cybDUaHC/3eDNVJaJSIxS02uGWW42mIp5EXjLcbhUWObcyCildocjqQyfLhcI+3UyoSAaXnHPcqNBxKezL5lgtdlhrdPlTG6UUCRCodUhqmoslCqM+HWEY1PvdPmx/8MfIBH0ym9f10RTKeVn+ZwfA6J4LvhVIcRPAAghckKIfofaMPCBEOIa8AnwH6SUv7bdiV/ZnnwjFEXhwYMHNJtNzp8//8zAhKIofHWhwK99cpOmYXJmfBTTcRiKhBAIPpidJxsKkvFp6MEgxVqdVqdLpdvly1OTWK7Dw2KZiUSM+WqNiUSMkWiEuWqNfKPFweEh5io13tkzRqfd5UGxzGQqycNyhXQoxPnxHPlanYV6k0AgyJ3VNU6NDhPSdJarDQr1JoVagzcyGQzX5s7KmkeQSoO9yQQAmWCAlUYX1RCcHhviWt7LeGi2wEaid/Fq0nuesOKArw1mCAKaSst1WCo12B987CpXmh2a3S4o0HE88t/NFzk4/HiYwFK5wan0MFeLBaaLXoQ3HQ4xW6vh82lcyOa4vrbCkZEMslmlUm0z226h6IKgolCwvP15X9010iv+WG416JgOp5OjPGgUuVsp0bFtzmVHaRgm+VaTI4k0Tcfg46U8b2ZGuFxcpmoYhITOw2qFN7PD3F5bJRNLM9/rjjsx5DWQZAJBPl0uMBzxIu8rtSaj0Qhz7Q770yl8qsZspYor4Y3erPZHpQpTqQSOlIyHgriui5CSv3RwgjvXrw3KVG3bfmkTTbci+U7nCmwFKeWBLV7P40kx0xN7PPW8534tltyyLBqNBo7jcPbs2WcSvG1a/M2vvsfPziyhKQpnxkeREgzLptjq8Gl+mWOpOCOhAFVXYrmSqmkR1HXOT4yxVK/zyfwSe9NJKp0ux0eyDIVDfDAzT71jMBoOogrB/lSSR2slrq6VOJAd6okFDLMvleRr03NMl6rsi0cJ6hrjiRim7fDx/CKT8ThHh9JUe0GftC/IVNJrwdybTHB5YZkTqRRptSfZ5EpC+An2au5H417xhNJ63CoOXgOLr+39HJHo45E/STz33KeqlBptCpUWb6TSlFteICmoawR7bkLM72euVOPqQoEvjY4PRGcWajX2JhM0mwatroHhOPh1jaCts9TtMhmNYEqX8UgYXQj2+SP4pSdM9ahWBeBYKkux0+bqygpvpcfp9OoeWqZNtDdA8W61xN5QEst1ubxa4EggimNLlltNRsMRrq+u8qXRSS4XCmiKSjIQwKco+IXK5eUCU6kEI+EI08UKuViU26tFJsJBdEXlzuoaIV1nbzJOvWsyW6kykYyjqQpISduyGI5F+fH/9Ps5f+YMZ86cIZVKUSwWmZ+fZ2lpacezzLfCdm7/ZyX6q8IrJ3mj0eDixYuEw+FthyUsVGr8yZ/+OX7u+m0ORMMczKSodrq0TS+S7bgOx+JRuhJuFivk4jFq3S6T0TDZcJCvz8xTqDd7C4MkEQhSM7p8MOu58HtScXyqCq7Dp0vLrLa7fMfeSTRFwZEuriv5YHaeyUSMdybGWGl3ubpUIOn3CiuGI2EyoRAfzy5xJjfCG9kMlxeWuV8oci43yoNVLwh3e6VMLPg4ANNoGRxLe3vrSNgjhFQFqXXSV/6AiiMlQ8EgxYa3JYsH/DyaKaMrCknfY4secx5r2I3FotyeW2NfMsl4/LE7bRgOTsfh/PAoAVfDZwserVUJB/wkAwFWSk1CPt0jvIBDiRQjoQRuV5AIRCk2TfapEVK2yoQvjGF4ewuB4O5ambPpUU+fzhFcLazwZnqYk6ksHywucr5XGosrSOoBqoaBI13OD43yu/PznB0ZZbnZZH8sQbnd4XKhwOmRYWK63xv7nEry6XKBUyPDWK7k0+UCe1IJhqNhFmp1yu02E4k4QV2j2u5g2A5v5kb4W9/7XQNJ5X6Z6qFDhxgbG2N0dHTQ7HTx4kUePHhAuVx+ru7IZ5F4p+rC3wy8EpL3b8Ty8jI3btzg5MmTRKPRZ5a2fji7wB/7336OUrvDV/btQSK5mi8QDwYoNFtMJWKkFcH9eovFepPT454Y4VA4hCUlHy3kmUjEODE67M3YMkyW6nVmS1W+Y9+kd556EwHcXCkxEgnznfv2kK81uLFcIOrzU6g3eWtijOFIhA/mFmlbNl+anKDUbnMtv8JELIbAk2qaXivjRyUa8ONIibAlmR4Rx+IxZgsV0sEAuqowu1bh+kyBiXiMpd5eWgBjvf27X1W9GWsSxvAP9uFxn49G22AqFCairust77rEegG4vnyy2pWE1nXolVptOpYNEpYqdWIhb3HpujaTgRjz1Tq6cNGEIBNO8ChfYaZSxXElxa5n7ZKxONPtNk3DRZouCVVnbyDMcqvJlUKBd7MTPKpWcKXkftEbYAhwqbDMmUiC+80m98pljqeHGA/FmK/VyQaDXC4s85WxCS4tLbPSbHFyOIuK4rWOhsLcXlvjQi5HudNhutHi3NgoIU3jRmEVXVUYT8RQhWC55sls/b6D+/i/f+WdJ/TY1sN1XYLBIBMTE7z55ptPWPlLly5x/fp1lpaW6HY3qUza4tneiG63+8JBt1eNV7Ind12Xu3fv0ul0BrPMtittPTM2yt/+g9/NtXyBn79xh3qny/HcMF3LYiTkp9Pp8LDeIheLDgbfdUyLatdgrdniwtgIiqYyU64xHA0zV62Ti0c5NjLEfKXGfKXGG5kUS7U6J0eyRAJ+3p+eIxYMcCiTwaepJIKe6P5sucrJkSFk1+CDmXlURfBd+/ZwcX6JtmXzxvAQqoRrSytMJuMk/D5uLa0CgslkjKFwmKurdd5IDpGNSu7Pe/vxtBog3/VIPhQNc2t6hYkpj+iFmve6NDSvwUVK/L0Hyq3bBFKPXXhNKBxKpLi0sjxwyWdXqrwd9wpCwj6d+bJXb95XYG2aJmFdx2hahPw9bXshOBRKcS2/wlA4RL7ZZCgUGsg1z/Vc9QOpFJcKy4Q0nbF0gpneEIOFtSITgRAFo8uxdJarhRUORmIstpssdS3ODI9wZaWA43hTYFZaLdLBoJcTn1tgfypJqdNBlwq241LtzYc7N5rj7lqRcqfLiXiUpmlyZ7XEoUyKmN9PqeUFCf2ayh85cZQ/dfrZYg8bA28bm1Ha7TalUom7d+9iWRbJZHLQcrrTvfzncS55H6+E5CsrK+i6zuHDhwcr33adaAFd47sOTPFdB6b4v3zHW/zr3/gtVvwhPno4w1qnS6lr8PbkGEIIlmp1YsEA87U6w5Ewb2ZTrLXazFTrnJvIsdJocmxkiICu8f70HEFd52AyhmGapIN+6obB1cIqh4ZSJAJBWpZJqdWhZZkENJ23JsYotdo8rNTZl4ozGo3xuw/nGI6EGY/HUFyodjxFkoVyjaPxCCsSLNfBh8L8ihe0up1f45196zqdXMGBTJKHxQpDoRBl2WJYD2GoLgU8kgdVncOZNHfWioje/SqUOhxKPnbtm/UW1aaJGhB0uo9D9A9nixwZTiOE4O5yEUXAbLmKqgjqXW/LcHFpmQOZCAFVxbW8eoS2ZXEom2Zltc1EPEanYrEnFudmyWsrXemROhMI8vXFPCeGsiiK4PrqCgD7QhFmikVs12Wm3uTd3BjvLS6y1m7z9liOh8UKlW6Xc6OjNE2TWytrnBrxil7eyo3RMi1urq4xFo2yJxrlo/lF/JrGV6YmubaYp1ZtcG5sFAHcXS2RCASIB/z8J8eP8J+dfGPb53G76HooFCIUCjExMYHjOFQqFVZXV3nw4AGBQOCJBWErfF5VYeAVueu5XI79+/c/4do87xzxqWiIdwIK/+B7voN/8af+KD/83d8BQnBxIU82GqHYbPNGNsN4PMbVlRL5RpNzE6O4UpIMBmhbFu8/mmdPMs5kOIgrQWoaa12TfKPFl/d6LvztFU+csJ+eGwoHub/mVYqdzKZRbIcPZxc4mo4jpYviCnRFYbFaRxOCw7Ewd0tNjvVqopOBIFO9CLsAKuUOWq+STZWCQC9IFui51g9m1wivc8VVBHrXM8/BXlmp60JUeWzJLVtSrneZCgRZLnlWNxUOUml2WZyrMBwKIcAb82TZnBkdxarb2I63/7RUhZOpLNOrFZKRYO96VM5lRnENl3bDIih1ItLHu9kJhgIh/IrKUMh7iG+uraGjcjDhRfXj4QiOC6OBIHsjET5YWORIOOJ9/2aHPfE4CrDSaBFVfZiOw9XCCr9nzxS3V4vcXF3jTG6EoWCIjxbzjMWinBsd4YOZeZq2w7uT47QMi4uLywxHw6TDQf7YqWM7Irh3/3aeQlNVlUwmw+HDhzl//jwHDhwYeKatVouHDx9SqVSe2st/FksuhPgfhBBLvdTZVSHEH9ziuP9ICHFPCPFQCPHf7vT8r63iTdO0HQtHVCoVWq3WE1JPk8kT/IkzJ6h2Orz/aJ73pmdZabZ4f2aekXCQkUgY25E0jC5ty6bc7vDO5BjVWo1lw2I8EWeuUiXj9zEci/KwWKbQaHJ+IkfHshmOhkFKbuS9AXmpYJD5ao2OZXNqbJhOxyCIoFSvU+xaHElGkbbERMOvqVxbLPDW1Bi351boWg77hpMIBR7Ml3hz/yiX88uUq20Wi3X2TMaxTe8h6ZoOmvl4MTQ6Bo/my0TGdNaviU798b2rtz3rrZo+ml0vSBdVBTXAtF1KK01yWpipSBSzZiJNl2qrS7ypoCmCpB4aEN6UXnDuwVKZrm1jqa6XQqtWaFsWhm1ztbDCZCyGjkJAVdkTi3N5eRlNCN7NTfDh4iJSgKpqTIZjTNca3Gs2eTs7wpWVFUwp2R+JIBBcyhcYi0bZm4jztUdzJAMBzuVGKDU7LNYbnMmN4FdUvj6/SMrvYzIa5u5aiVK7w+nRYTRV4fuOHOT73nhSs+BZeNE8uRBiYOVzuRxXrlwhkUgMrHwwGBxY+Zcg/fT/kVL+g2dci8r2gxE3xWsthtmJJV9cXOTevXskEolNSwQTwSDff+wwf/8PfS//7I/9Yf7pH/0+3p0YZanRwnId1podYn4/Z0ayTK+ucafW5EAm7Q2+H8oQ9+t8vLhMsdXinakJpPSCVKqicK9Y5sRolkwoxKf5ArqqsicZ9zrBFAUDQduFI4kozY7l1d6XquxNRFGEwOhahHw+XClRXIiq3t53drlCNhpmseS55Ck9QKX+uKhJdB5bhXqtjePCoWSaRvtxvfXsfJn96SRBXaPW9PauigMncp7ueSbxuEtsudxgrdJirVijWG9j9vflwPmxHNdnCiiqQBUCzVEo17pUOl0OZJJeL/5Qhkq3S1jXuVP04gnZSIRPFvNEFB/DgbBX4ONK5taKHIwlSAcCTMUSfLK4xNmRUY6m01wsrDAcibIvkQChka+3OBAJ43ccPp5f4tRwlpCuUWp1cIGj2QwKcHExz4F4lIlIhJulGvWuwbt7vG3PHzl+hO89uBfLsjBN0xt7tU2E/GVUvLmuO5hw2rfy+/btw7ZtfvEXf5G/9Jf+EtevX+drX/salvXZS2g3wQV6gxGllCbQH4y4LV5bMcx2JO+7RH3pqJs3b267KOiqypf2TjKpK/y5k0fpBsP8zsNZbi8u8bW5JVQheGtyDKenMmLYNjfWKkzGokymE6w0Wqw2W4zGIjS6BofSnvv58fwSR4czxPx+7q2WGI6GaZkWYZ/OsKZxt9LgSDbDraVV3hjJcH+1xIlkjJtzq+xLRykCpUabmOJFm6vtLl+amOCjFW9vO7NUxqc/vvWPZkuMD4dZrLfo9J6PTsWgWPMWgrBfp1U3mcSPFXFZ7gXVogE/xXwLn6oOqt6SoSDVXnVZ0/TuX75SI6irpHw+zLbnEVS7Bm+NjfHR/CJn9owy06iiaQonM1lSgSDhrI9owEep02Gl3eRej+wj4bAXNEsm8Ls2t+tN6HQ5lR0ekG26XGYiHGciFmOuVuPcyCggWKg3UHx+VFcyJgR319YY9vkJqSoL7TZDoSDL9SYHYhE0VeNGqcJINML+VIKleoO/+h1v85379uC6LlLKJwjeL3jp/2/js/VZi2E2LhRCCMLh8BOp4V//9V/n3/ybf8PCwgJ/+k//6ef9iB8UQvwZ4BLw30gpN/aqbjb0cEfji1+ru74VaddruR09enQgtbzTPbyiKFiWxcFMCre0xjvhKX74e38PlxaX+dW7D5gpVdEUhUKjyelsGhSFT+aXGIlGyEbCXtTZdlhrtii1O5yfyKGrCh/OLnJ4KI2uKKQDAbqdNjerLc6Mj3J1scDxsSyqUDieHeJhoUw2FuZRqcHeVBjNlSyv1lAFOBLsroOmelHk0XiUmN9PudkhFvTTWOmyVwY9K13yrHS13CEdDbJcbZKOhFhZM3n0YJWpE1mW8UiuC4W1UpOTx0ep9DTUhxNhqs0Ofk1lpd7Grwmqhs0bo2keLpQZGwrhUwVBR9LueJ5C27Y4PzrKg0IZ23G5r5YxbIdcPEq+3uDs2Ch1wyCc0Gn0Rj8v1+oEdB9nR0ZZajRYabZYabY4kc3iV7ziFk0R/N49U3w8v0THtjk5nEVTFK7kC/hVlbO5MfL1BnfKNd5IJqjXGt7wymCA6VqDk6NZAprG9eVV/t5//Pt4u2fN12sHgkfiPuEdxxk8N4qi7EibbSfYzhuwbZszZ87wQz+0uc7DVh1oP/qjPwrevPG/ide69DeBf4jXbroeLzz08LWRfKsUWq1W4+bNmxw+fHgw6L1//E4LFVRVxbIsrly5QiqVGog+fn88xvcfO0zbtPhwdp5vzC3yOw8esdLqcDo3gov0Rgk7Lo/KVZLBAG/2dOGmSxWOjwwRDwT4ZG6JyXCAkmFzdjKHcOGdPePcXlxjKByiY3j9oIlwgNV6i0gozPRcEdN2OTAa50Gpxp2HBQ7k4txdrRDz+8kvem2dUU3BAFYKLXJ7oiyUvAU8FQkSDwVYrjaJ+n2sAIbpEHXWjUvqWe/bt5fZe8C7d6FeGWom4mex3mYsk8IRkojix3YkXaFycnSEy3PLTA1HSPo0quUmbihArWvw5sQIV5YLHBse4tbqGkFN42G5QrXb5fToCIV6k4PhCPFImEuFVSr5AhdyOVwpKbXaBFSNWytrnMuNYtgO///23js6rvLO/3/dO71pNNKMumXZci+SjDEYs5QYDARwi1kMJCG7hCSb70JISLJpm0BO2oaz+SaQk4XsN7Ah+UGya5sSICGwJKGZaizLcpGbuqWZkTSj6eWW3x+jO0iyukc2Mnqfo3NkGN37zMx9P8+nvj+vNLeRbzFTU1JMS18QbzTKgkIXZXY7r7W0owKXzK0kEI9zoi/I4sJ8oimJOVYTkVCYo4kUP9pwKRcOavYZjsGnt6IoQ35UVSWdzrhWp0P20x2sMEYHGqqqerXfBUH4f8CzI7xsykMPzyjJR9JSb21tZdWqVacMPhxJT300aEMPV6xYccqEVACr0cCVi6q5clE1tyysyjShBCP875HjVOQ7Oej1s9BdSJ7JyJ72k6Rlhb+rmkNckni9pZ1yq5kCZx6GRJq0LHMyGMamNzLXk88Jbx9WvYFFZW72tXaztrqCE229LC0vYl9rN5GEwtJSD0eO+omGMqdgf2+QvmCMsgoLeVYrPcSIxdMstzpoJ0Nym9GAry2IYMgUy2gQ4iqCkJGSisYy1yuwW/EeD1Jd5co0zABWY+arLbZZaTjShaPSgMWkx6TTMSDnhs5kYoHNwTudXTisRmw6ESkWY4XLhU3UU1dUjMtqpjcRx2400NwXIJJKo9psNPp6WV1SgigKvNORedZWl5WgKCpJSaI7HEFEoLa0mBN9AbpCYURB4LzSYlTg1ZZ2yvIcLCwsYE/HSSKpNOeVeBB0eiJSFKvVQkcwxA+vuJgCKcmbb76Jw+HA7XZTWFg4qqDiYMLLskxzczM2my2jhivLSJKUtRQnQ/jxmlNisRgFBQWj/v+xIAhC6aABhlsZubPsHQYGIwKdZAYj3jKR658Vn1xVVY4cOUIsFhtRS33468eC3+/PDj0cieDDodeJLC90cVXNCr5y2UU0dvt48cgJ2gL9vHy8lXyziTWVxXQGQ5zoC7LMlYczz8E7bSdZ7HETSiZZXOxGkmT0iNRVlCArKumUzMXVlRxv7yEWT9Hu78eo19EdCFOZnyl46eqLUl3lor8745sbZZFk7P1aajnyfsDGIIr09kRZWOPOkhIgFkyyrNTDgZN++gbUYAocVo73ROk+EqB6npP5BXmU5OdjVI0kI2mSKRl/JEZ5QR7h/hSxAj1z3U6CfTEsHj0Oo4E8sw2daKCxt5d5TgcHerqw6XXIAiRkmZWeQpREioUlpUSlNMkBqSlFVaktKUaSZY73BgkmEiwsdOG2Wnm74+RAG68HvShyxN+DarfT3BdgdXkpRlHk5ROtGASBteUlhNISTd5uVpYWk5Ikfr71oywr9mSfmVAolK1FF0URt9uN2+3OkngwNP3/ZDLJypUrs2b7YJN+uFk/FuknMljhNPLk9wmCUEfG/G4BPgeZDjTgV6qqXjug/aYNRtQBj6iqemAiF5+2k3y4OoyWQkun0+zbl+kQqqurG7VMcCJDD1tbW/H5fCxduhSfzzfumlRVxWg0cuzYMQA8Hg8rS4tZWZqJUO/v8vJmSwf/XX+ASCLBSrcLvcnEvpNelpUU4TQZqe8N4rJY6OoPk2cwYdXp0asiiqRytDtD8CWVRdQf76Jufgn+cIxATzS7Bqsi0h7J+MJdvhjzKgpgQO2lvzuA3aQnkpRQ0hlmO1QD6dT7n0MknMAg67GZDIQDA00qpszXqMgqrc0hUpKCMaXjxMk+yitcWC0GvP0Rqt1zON7RR1Ivs7iokLf7OilMW1lRUsSbrZ0sqygiz2SiON+J3WrBaTHRHgiiSmka/b2ogDuVosnfy7IiN4VWC2+2d3LI62eOMy8jxuGwk5IV3mzrzIws8rjZ1+2lP5Hk/PLM9NACqwWBTICzwmpmjsvFfn8PkVSaC+eUYxBF7vjIxSwpet99EwQBp9OJ0+mkurqaZDJJT08Px48fJxaLkZ+fj9vtpqCgAFEUOXbsGOl0mmXLlmWfMe2U1+v1WVNeIz1kiKyRfTjhJ6IKM5ZgxFhQVXXEKN3gDrSBf58yGHEiOKPmeiqV4p133mHBggXjnrpjkVxRFA4cOIAgCJx//vlEo9Fx/XftCy0oKKCmpoaenh4OHDiALMu43W48noxG26KCfGp0Mr06I7u9fZkTvMiNgMCJviDlLicWo55CqxW70YAeEVVSiUlpSgocGApFDrX5cdksdPSEKDSaON7Vx/w5BZzwBpCiSRwWY2YSi6ySp3u/2UQvmql2GdjX6SfYl0m3NTd5cZW///D0BaIku2VWXDSHPYHOzHsbyHkX5Jnp688E7ryBCKIo0BUIM7c0n/kWA+FYkvLCPAqdVvpiCZaVFXHU10uqUOa8ihKCySTz8p282dLBIk8hr3W1YxJF7CY9HquNalcevdEoRlFATiR4zddDnslIXVkpLYEAJ7w91JYW0xdLsLq8FBGBV5rbMOp0XFpViS8apcnXy3kVpXSHIixw2Miz2Xijoyv7mpQs88VL17LAPbbpazKZhkgmB4NB/H4/x48fJ51OYzabWbFixaiHyODgncFgyJJdi9prv+t0uqypPxOnp8AZJLnf7ycajbJu3boJfRhaMG04tEh8cXFxNnUxnv+ufXFaF5HVaqWyspLKykrS6TQ9PT00NzcTDodJp9NUVVWxprKSa0WRaCrFC00neP7wMdw2K2lZoScaJyXLWAV9Rr7JIGKQRMyCjlAkgSffSqnTQSKcwjBQiqAbCDoaBRPVJXbqm7sodtkJ+iLZdQb6ojiUTAWalM48nKmUgmWAxHazgeSAX0+/lPXN+waaXtz5Dvr6E7jyLPTG4pR5nLSHQjgtJg40dWN1W6gsctLY6qW0JA+zQU+5y4GiQFyW6QlHcdoKOK+iBKNexxLFiVkn4EvLmUIhX+a0Pa+shFg6zbJCA6ok8cqJVgAuKC8hrUIkmQIV9nZ1s7KkCJfZzOst7ciqykWV5cRTadLJJO78fN7r8lFd6KI8z4EvEuVH113BvALXuM/HYIiiSEFBAS6Xi6amJiRJwuFwcOjQIdLpNAUFBbjd7jFr0YcH74an6FKpVFZ1eKRrfOhq1+F9c11VVY4dO0YoFMJms014t9PpdKd0BYVCIfbv339KJH4syWdVVbNR/ZG+HIPBQGlpKQaDgaNHj7JgwQLC4TBvvfUWdrsdj8fDxqUL2LpyCSd6Azy1/zDHevroCkVIyAp9yThOkxmLRY+QBqNOj02AxgMnybeZMeh0CCq0ecOUFztoae6lyJM5md0OK8cbuqlYmE8gkiDUlSDUn6BimYu+lgxxbRYjumjmBHHZLXh9GZInglEWeZw0+fpJy5kNwTwQbHO7bPTG4uTnmQmkEigxBYvJkJlfnpKZX1pIKJ4gaZAwGvToRZGUqlDlcfFeezfLyzy819nN0gIn9f4gdqMBj8NGhdlModWCNxLhWE+AVeUlNPb1s6zYTZ5ex94uH0lFYWl+HpF4nCWeQow6kVda2sgzm1hVWszJ/jDHegPUlhbhGxDRNOn0HO3p45c3XE+lyzn8K5oQVFXl8OHD6HQ6li9fjiAIzJ07F0mS6Ovro6uri8OHD2Oz2bK+/Gi6BsNTdKFQiO7ubpYtWzaqLx+NRqdsrk83pvUklySJhoYGbDYb5513Hm+88caE/3Z4Cs3r9XL8+HHq6upO2TFHMu0181xRFARBGLuPvb0dr9c7RNBCVTPa7T6fj9bWVgwGAx6Ph/+zdhU6g4G/HG3m7daTvNnSTp7FTH8qictkwuk00y/HqZ7jRi8IHDjuZV6pkxPd/czLd/Ku3E9Xdz9zq12YhMzDVGS1YdTrONmV2dQqrHZ65CAArjwL7U29VNQU4LRZ8A7kyHWCSPfRIPmlRnoHAnDaZmY2ZyLPBqOOSoudnlCU4iIHelGkvS9EVXE+FpOdiJTCojMQSMcpzrPTn0qytqqcI91+ajwFxFVYVV6CzWjAH42RliQauryEkynOryhDURUWegqx6EXe7PRi0uu4dG4lvZEoJ3oDzLVZORaNscSVj8tuZXdLO2lFZU15CRKZfnlUlZZAkF/duJHyQf3wk4Gqqhw8eBCj0ciCBQuGfNd6vZ6ioiKKiopQVZVIJEJPTw/79u0DoLCwELfbjcPhGPEZiUajHDhwgJqaGux2+6iFOH19fR8+cz0ajVJfX09VVRWlpWOPpBkJmgmuqionTpwgEAiwZs2aEVMnw0/yiRJcURSOHDmSLWQYfNJretp5eXksWLCAeDyO3+/P+vEL3W7WrVtF4IJa/na0hWcam9DrdaQFBbNNT280QufJfpx2M3q9AaNBR6ArjICKioDLbCYUGNA9aw9QPDc/m/RM9MSzprjDYsIHFBktqIPWF40kSMZlqg0Ojhv7iScl+vozpn8sFsOo16GLqXS0B4iZVVYXO+nwBplb6KQ3HMNsNqA36VBElTyziUg6jcNkxB8MUlWQB3oDqVQKnSjwWnM7pXl2UopCkd1GXVkJvkiUo/5easqKafL1UVNShNNs4t32k8TSEmvmlJFISyy0WBEViTfaTmIUBVa5XURTKVqCIZYWe4gkUzyyfVNWN32yUFWVAwcOYDabT2mKGg5BEHA4HDgcDubNm0cqlaK3t5fW1lYikQhOpzMbvNPr9USjURoaGli5cmWWwCMV4uzevZvjx4/nRGJqOjBtJO/q6mL58uXk5Q3dnScqkaP55Pv27cNsNrN69eoJReInSnBJkti/fz9Op3NIS+xosFgsI/rx0WiUC1wuPrrpCo6FYvyhsYnWeB9eKYyr2EqZJY/G492sripjX307C6o9HG3rwXsyRCyUIXmgL8bcivcDTRa9nvllBRzv7Mso2QCtjV7m1pYMWv9AZ1tKpVy0EHcLhAckocwGI8XoaW724co3IaVT9AdiFLsdJBWZMlceolkkLSioIvijMVxWE519fZm8fTKNrydAtbuAd9u7WOAuwGYyohtwwd5u6yQpy1wwp4yULDO/sACjKPJqc3smgDavknAyRXt/iOpCF/tO9jDXbqGysJD9Xj/BRJKFjsyU2W+vO58C09QGEWgBWJvNxvz58yf990ajkdLS0qxqTH9/f/Z7FUWRWCzGsmXLxjyh33vvPf7lX/6FN998E4/HM6X3Md2YNpIvXLjwFBNaM8En0iwgSRI+n4/FixdTUTH2xBiNoMMDbKMRNx6P09DQwNy5cykpKRnxNWNB8+O1h0PrPyYY5GOlTgJ5eprnVfLo7ga6InGcBUZSA4UrxoHd3uO0IdqtHG7OpP6MyvtrFWRwCANfzYBgYyIhIQ3Usut1In19md9NRj1Hj3hxOEwU5llw6XV0HupF1AlEZInyShf2tJET3gClsg2bw4RoVEkIMjp9Ro3FZTPT4u/B48zjkL8XvSgyr7CAYDzJRVUVKKpKRzBEntnE8d4AVqOBC0sq6E8mON7Tx+IiN/u6vCwtclNkt3HI68cXjXF+RSk94Qjz7VYKnE5ebe1ABNbNrcBuMnLnhauQoxEaGhpQFCWb5Rg8RHA0KIpCY2Nj9lQ+XWjzzlwuF7FYjL1791JeXk5HRwfHjx/H5XLhdrtxuVzZE7u+vp4777yTJ554grlz5572GqYLZ3S4glbaOh7Jg8Eghw4dwm63j0vwwRgrwKahv7+fgwcPsmzZMpzOqQV5BkMUxeyESy1vX1VWhjMQQKmuZMe+VsKmNG+mexCqdNT39uC06DGpAor8vovRcciLXi8iSQqpeIqu1gCmQh2p+PsZhrQviV4vUuC0EohkAnMMSEUVuuy0tPdRUGAlnkhTVeWmv7sPh82CrztNVbGLtnQYuwDWeApJUUGvw2Q0cLjbS5W7gL1dPjw2K6V5dkRRpNzpoCcS42hPHytKPXgjMZYWu7GbTDR0ewnGE5xfUUYomaSmtBiTTsfLJ1oRyBA5FIuRSKZwFWTGDlc4HcwvcKGi8q0rLsFpMYPHTVVVVdZ01qyjwXnv4c+LoihZK6yqquq0v8PBiMViNDQ0ZCXLgKyQhN/v58iRI+zduxefz8czzzzD008/TXV1dU7XkGuccZKPV8WmlbrW1NRw9OjRCV1XlmX0ej379u2jqKgIj8czYuS0u7ub1tZW6urqcqrHpaoqTU1NyLI8ZBpMTU2cZt+fefuEFyMqskUgYVGQHQKvJXzokirmYh0FKT29x8IsOL+Cw80+Qn0x4tEUC2vmEOjKDDgUBPB3hFi+dg4xWSYwoCSTGBBYtNky+faCAhs9kThWmxHFKrK3x0vSopBQFcxWPdWlTlKygl0U8YVDdPb0UmA209DtZ6nbhdVszmZFGrt8RNNp1lRm5rmX5dkx6nS81tyGKAhcXDWHhCSRlCQcJiO7Wzsoz7NT7S6gvbeP1v4INaXFtPeHWV1egtVgIJpO87PNV+MwmYZ8hsNN52AwmC12MZlMeDyebES8oaGBgoICKisrc/YdwvsW3rJly4ZEyjUhCbfbjaqqBINBHn/8cZxOJ5/61Kd47LHHpuQunCkI46hMTlmCUqsTHoyGhgbmzZs3YqphcKnrypUZza49e/Zw4YWjd9MN979jsRh+vx+/348gCHg8HoqKijCbzTQ3N9Pf38/KlStHLKM9nfe5f/9+HA4H8+fPP8XMTKUlPv7tx4jIEilJJinLqFYRISojmwSMDj2STiUhyVj0OtJpFV1URkgouJxmIj1JxJiKy24i1hHBlm/BvbSQYwd9qCpYPQb6k2ncxQ66IlHM+XpCSQmDRUQB0qqKTi9gUATmVxRiMxkx6/Wc6Okjnkhis5qJpSWqXXnEk0l6Y3EsOh2t0QQWg566shIi6TQd/SHm5OfxXkdXRm3H5aSzP0xzb4Da8hLagyHm5udhNRh4s7UDWVWpKy9BGdgw9KKIIIj8fOs1WEapOx8Ng7/XcDicrXrLy8vLmQRyPB5n3759LF26dEwL78iRI9x666089thjrFy5kv7+fqxW66i19BPEtOo4n1GSHzx4kNLSUlyuocUOkiSxb9++bCRbqzN+6623uOiii0Ze2DgBtmQyid/vx+fzEQqFsFgsLFmyJKcPRjKZpKGhgfLy8lGnwQC8uLuJB3a+RiiWAqsOfUwBs0BSVcAsgghKWkHRC2ASkFUF2SwgKJlRyKKkIogCSDKKTkQnCEiqiqiCioqqExBEASGlIBpEJFlFJwrZqS26NLjzrcwvLkAFDnX7cOh0JARwWSyUOuwE4gmO9fSxtMTDIZ+fMpsVqwgnYwl6EinqSorwx+KUOGxYDQZ2t3YAmRSbpKgoqoJeyLSRWnU6aitKSaQlOvpDFNkzkk33XbcBs2FqG6wsy+zbt4/CwkLMZnOW8Hl5eXg8nmxEfCpIJBLU19ePS/Dm5mZuvvlmHn30UVatWjWle42CaSX5WTfXY7EY+/btY968eUOCYKIojqplPZEAm8lkoqioiO7ububOnYvFYqG1tZVoNEpBQQEej4f8/Pwppz0ikQiNjY0sXLhwXJG/DesW86e/HORIqhchDf2CjKiK6JIDwxcEgYSkgqxgdhqJR2QQVXRmHaJRIGwATCAZRQQFVBFERUAVMr9nPhRAJyDLoEdA0qvolMyABJNBR3VJIX3RGE3+XqrzbMQRmeOwYTMaaPT66YvFOX9OGcFEgmVFHkwGPa+3tCMKUFtUSDwRRy+nicRi7AmEcJpMLC/xEEun6ewP47HbOODzUu10UO5ycbinl1AiwbwCF8tLPHz18nVZTfTJQpZl6uvrKSkpobw803JaXFyMqqpDIuIGgyEbvJuoO6YRfMmSJWMSvK2tjVtuuYWHH3441wSfdkzbSa4oyillqceOHcPhcFBcnGkI0WRwV6xYMeIHvHv3btatWzd0QYMIPl5/7/79+1mwYMGQ6jhFUejr68Pv9xMMBnE4HBQVFVFYWDhhiaC+vj6ampqG5E/Hg78nzGf/9X/oT6dAyOyucbOAAYG0DsS0gmwUsQkiUUFBF5eR7DqQVQQVHAYD0f4E6XwdslFAEFUkm5AZtaRTEQdGHKMXUPQqupgKBhGdDIvmuGkO9RNNpVlakIdoNGHS6TDpdexuyyjonD+nDElRCMeTWI0G6rsyMs2LPIWEk0k6giFKHHaO9fZR5bBhADpiCZKyzNyCfKLxBE6TCavVzDvtXdgMBhZ5CllRUsRdl16IfoqbqWbllZWVjVtvEY/H6enpwe/3k0qlsoUuo5WzDiZ4fn7+qNft7Ozk7//+7/mP//iPU57HHGFmmusjkbylpQWj0UhZWRltbW10dXVRV1eXnUs+HMNJrlUZjVfB1tvby9GjR1mxYsWYJNTaF30+H729vZjNZoqKisYseezq6qK9vZ2amhrMZvOIrxkN//eXL/GXN46is+iJJFOZOnqTgaQikUzLKBYdJlVF0okoSRmH2UAklsJo0hMfsLl0cQkRkHQidr2ekEVB1UFaryI5BFRBxZAQUFAREMizGOkjhVEUWORyYrRaCETj2EwG9nf7mZOfR2W+k0A8QWewnzJnHu3BfhZ5CrEZM6OAJVmhYqBdNs9sQlEU9nR2YzcamGO3kkqlSEgyot5AeyTK/IJ8Shx2qgry+crl6xCn6B5JkkR9fT0VFRWTTnXKskxvby89PT309/dnS5S1XvRkMsnevXtZvHjxKe7jYHR3d3PDDTfw05/+lMsuu2xK72MCOHdI3t7ejqIoRKNRJEli+fLlY56eGsknU6La0dFBV1cXNTU1o24eoyEajeLz+ejp6UEQhGyk3mKxZEQbW1oIBoNTDt6lJZnbv/j/4QtGSelVrEYDMUlGkBXUgZle6ATMeh0xVcGKQESvopLJr0tpGX0aZJ2KYhQxKKD2J5FcJsSEjIhIWi+jegwkLCCokLLAHLsVt8VMFPCGIhTnOfBHoiwoLMCoE6k/6UUUBcryMqWvFr2elKzQ0O3FZTFT5cqchP3xBJKSKUNd7CnEZTHT0dtHOC3htFgIJ5OUWEygKFS5nHzlsouG5JUn9Vml09TX11NZWZm1/KYKrUTZ7/fT25sZZZVIJFi4cCElJSWjPk8+n49t27Zx3333ccUVV5zWGsbBzCS5qqqkBvTANLS1tdHS0sKcOXOoqqoaNwC2e/du1q5dmzXRtYaA0e539OhREonEuJvHRJBIJLIRXU0+yGKxsGLFitO69s6n9/DI799Ej4jepqNfldHFJGSrHrMqkBBVjJJCwqQbiElkzPU8UUdIVTNjmhISilFAr4jIqJiMAqmEAnoRXVRGb9QhyirLVpZSvNxCV1Li1Y5ubEYjbpsFs8GAQRSIpNI0+Xsptlspc+YhCgL9sQQJWaajP8SKEg8Ok4nucIRwIolJr0dSFOYV5iPLCns6unCaTRTYrDgGqtZCiSR/N28On1y2EL/fP+IpOh7S6TR79+6lqqpqQkIgk0EqlWLPnj0UFhYSi8VIJBIjFrr09PSwbds2vve973HNNdfkdA0j4NwgeSQSYc+ePeTn51NbO7Hpq2+++SarVq3Kkns0gmtpLLvdPm798mSh+YTaw5mLwN0ddz6GPxQlFEli1AvEVRnZbsAgiqQHdOf0Bh1KNI0uISOb9ZlEuQiKUUSISwiSimrVg6JilwUiOhUrAklJxSKImA06PvfxpSxePI/y8nLi6TR7Orp46Wgzh7w9BOJxeqIxassyRSzt/SEkJXNvg07HnHwHibRM/cluSux2rCYDeWYTqqrij0TpDEWocNgpdjqQFIX+eJy4JLNt5VI+u3Z19r0OP0V1Oh0ej2fU4FgqlaK+vp558+blvEw0lUqxd+9eFixYkA2WaoUuPT09BAIBOjs7aW1t5emnn+aee+5h48aNOV3DKJj5JPf7/Rw9epTKykoikQhLliyZ0N83NDQgSRLFxcW43e4RT4FEIkFDQwMVFRVjprGmAu3alZWVWZ9QC9z5fD76+/uzKZzJBO7q69v49refQjIJmBRIAUgKRruRRDQJOgG71UhEVlAlBcWSua7VqMvov4eSiKpKyq7HEEqDXoeaksGixySKpJMSS+c72Hrt/CyhBp9SqqrS2O3nteZW9nR08W5HF3Pz89DrdORbzMiKTEcwTE8szmJPIXaTkZQkEU1LBGNxIskUK4oKMRiNNPl7sRkNpGWFz154HjfUjj3VZLiFpG2YTqcze4JXV1cPCZbmAiMRfDhUVeW9997j61//OoFAgMLCQv7t3/6Niy++OKdrGQEzk+SQ+UJbWlro6emhtraWWCxGZ2cny5cvH/PvBrfxxWIxfD4ffr8fg8GQ9ZNNJhOhUIgDBw6wZMmSMYMnU0E4HKaxsXHMa2spHO2UslgsQyqzRkM0GuUrdz1OIJAklJZAr0OUFBS9iJCSUcx6TDqBxICRYDbqiCoyhkgawWJE1r6zeBosBlBVdHEZdCJ2qx4lKfOTn2yjoqI4W44ZCARGNZtPhsLsbm7nxaMnCMTjdIUiCAjUlBaRkCQOeXsoctgIJZLYBSh0OOhNJmkN9FPlcmI1Gtleu4zNK8bfvAdDa9HUzPpUKpV15U7X3RoMzTqYP3/+mJtHOBzmhhtu4I477mD79u3Z2Mx4KdIcYGaSXFEU9uzZg16vZ8mSJYiiSDgcprm5mZqampFvNk6ALR6P4/P58Pl8pNNpJElixYoVU1bJHA1adH7lypUTVvtQVZVoNJo9pUYzS7XaeWdeKfd+4xniiTS6PBOKqiCpIKgqik5EBUxGkaSqYhNF4v1JBL0Oq0lHWMz8t0QwgewwIsYlBEHEIoCSkrnwwiq+9p0tp6xPM5t7enqy/fEej2dIliCcTPJ6czvvtHfyWnM73kiUxZ5CRCCdSGAwm2gPRYgmU5w/pwydIHDd0oVsXL54yp+3FukuKysjmUzS19eH0WjM5rwnm8UYDM06GI/g0WiUG2+8kdtuu20qgxFOFzOT5KqqcvLkySEfbDwe5/DhwyMWE0w0gq4JOPb09ODxeOjt7UWSJNxuN0VFRSMqd04GJ0+epLOzc0rR+cEYbJZq69Pr9Zw8eZLa2losFgv/9/vP8OqrxyAto7MbicsyqkkPkoJq1KGPJxH0etSkgmIQQSeiygpGiw45LoMKFpuBSFLGIAoo0TROh5n//M1tWK1jr13rj/f7/ciynFW7HdwBlpZl9nZ08ZejJ3ix6RgFdjudoTBFdhsVTge9sQT/sKaOqxdPvUFDy1UPT2XFYrFszltbn8fjGVXcYSRoBB/Pv4/H42zfvp2bb76ZT3/601N+L6eBmUlyyJhJg6+fSqXYt28fa9asGXqTCbaIKorCoUOHEAQhax0A2f5un89HPB7PPrCTKWHVxCnC4TArV67MqbmYTqc5evQoPp8Po9GYXZ/RYOHOTz5MKJrMDEBMpRHNekSjnlRSynSYWTJmtarIKFYTQjINsTSCI3O6WWWIqTJGow5k+Mjli7jjqx+d9Po0QkWjUVwuV9aPj0QiHDhwgJUrV9IVT/K3460c9PrY09nN966+nMuqq6b8uWj14uMVo6TTaXp7e/H7/UQikTE71Ab/jSZaMhbBE4kEt9xyC1u2bOFzn/tcToO22vUvvfRSkskkkiRxww038N3vfnf4y84dksuyzDvvvMPatWvfv8EECZ5Op2loaMDtdmcFHEeCVgTh8/kIh8O4XC6KiorGjIQrisLBgwfR6/UTEpCYDDT971AolG28GeyHvvZcC++80YkkgyLJCDoRg15Ey0uIegFJEFAlGYPDjNwbA1VFyTMhhBOIoh4lLWGwGHDlWfiP330Wo2nq1cpaf7xm1qfTaebPn09ZWdkQPz6eTk+60WQwJtoQMtL6NGXWQCCA2WzOxkE0y0uSJPbu3cvcuXPHTMGlUik++clPsmHDBu68886cExzed+PsdjvpdJq/+7u/4/777x/CAWZy7fpw7fXhMk2DNa/HSkVpPb7z58+fkJSzpumlPbBer5empiby8vIoKioacgIM3jxy3fiviQsC1NbWZh8izRdWVZUlS5Zw8MaHCYaTYDGCrJAGUBUw6DHpRCRFRZAV9JEUysDnJPbHQDeg5QYYUhKXXr74tAgO7/fHC4JAIBBgyZIl9Pf389577w3x4y2n4Sdr/QojKQdNZH0FBQXZOIwWB9m/fz+qqlJQUIDf7x83x55Op7ntttu47LLLpo3gkOGAVnWZTqdJp9PTdq9R1zCdJ3k6nT5FRXX37t1cdNFFE65gCwQCHD58eEoPxGBokXCfz0dfXx9Wq5X8/HxOnjzJvHnzTruqajhkWR6iXDLWe3z31Sbu+8oORIMOSRSRBRGDqJISdRhESIeiCOgwWQwkEBEBNRRFdVgzJj4ChU4L//mnu3LyAGlTaerq6oZkCSbix48HTTdtxYoVOVc31RRdDAYDsiwPcTsGHyKSJHH77bdTW1vLN7/5zWknnSzLrF69mmPHjvHP//zP/PjHPx7+kpl7ko+GiRL85MmTdHR0sGrVqtOKsEJmR83Pzyc/Px9VVbOnu8FgoLOzk3Q6nfGTx0h9TRSaNl1JScmElG3Ov2Qxl169gpf/UI9OLyInJUSbETGlIMsKoggYB1RiTAYMsoQkiJgFlZgKVrOej/+fy3PysHZ3d9Pe3s6qVatOqUsYrnM3WMllNEINRiQSYf/+/dNCcEmSOHjwINXV1ZSUlAxxO44cOZINyLrdbu655x6WLl16RggOGeuyvr6eYDDI1q1baWxsZMWKFdN+Xw3Tbq4PhqY40tXVhcfjGbX+W1VVjh8/TjQaZfXq1TkNgkGmZLGlpYU1a9ZgtVqzufh9+/YNEZuYinpMIpFg3759zJ8/f1IVW//83c0c2dNMZ1sAQZJJR1KYrEZSiooqSWDM5MNNyURWxFGIJzFJCuVVLi7fVDfptQ5HZ2cn3d3drFq1atzafIPBQElJyYiEGikfrxF8Mp17E4XWa15eXp4tWtLcjsLCwqwU8x/+8Ad+/OMfk06nqauro7u7e0pKwlNFfn4+l19+Oc8///wZJfm0muuSJA1RUZVlmUgkQnd3d7Z4ROv60h4GWZY5cOAAFovlFA3tXKCjo4Pu7m5qampGPLUHi01MNjWnPchLly4dM1o8Go43tvP17Q9hMBtJJmWsdiOxhIyqKIhWI0ogDIqKkGdDr8hI8TQ2u55v//Yfmb+o8rQ2w7a2Nnp7e6mpqTmt64xUxpqXl4fP5xtRM/90ofWaj9eKqigKd999N3a7nS9/+cv86U9/Yu3atdNONq2IKz8/n3g8zlVXXcXXvvY1rr/++sEvm7nRdY3kI00x0aKOXq+Xnp4ejEYjBQUFdHd3U1FRkRUHyBW0SS7xeHzCDSyTSc0FAoFsj/npPMj/+a87+N+d76LoDKiAYDCgJJIYdCqa0I6iKIg6HSgKn/7OtSy4oITe3l6sVuspm+ZEMDj6n2vtcL/fz8GDB7FaraiqOqV892jQTvCSkpIxS5oVReHrX/86AA888MAZ1UdvaGjgU5/6VLaK88Ybb+Q73/nO8JfNbJKnUqkJ+d/aw2AwGLKqLkVFRadVkKJBk+81m80sXLhwSg/XWKk5v9+fFZ883dgBwDe2/JT2pi6ioUSm60yvB1Qwm9HrBFL9ERzuPK6/7RJuvOtq4P1NU2uV1bIMY1WMaRtfMplk2bJlOX/4+/v7OXToEDU1NVit1qwf7/P5JuzHjwaN4MXFxWMeCIqicM899xAOh3nooYc+qAMQZi7JH330UebPn09dXd2YJ6c2jVI7BbUorjaO+HR8ZC0IVlxczJw5c6b8XgZjsNa63+9HVVUWLVpEUVFRTuIHsizzr1vvp+m9FlRZQRxwK0SDDiWVBkXlzgc+zuUfu2DUawyPhA93OzSFWS2Nl2u3KBgMcvjw4Wx133AM9uMDgQA2my2r0DOeFTJRgquqyve//326urp4+OGHcx7bgYxGwq233kp3dzeiKPLZz36Wu+66a7KXmbkkf/LJJ3n88cdpampi/fr1bN68mTVr1gwx2dvb2/H7/axcuXJUH1mrV5dlGY/HQ3FxMVarddz7a/n16urqnLctDg4OVlZW4vf7s6k5zWQ+HVXYcDjMA3f/Fy17vPR2BVElGQEVg9XM33/xGm740sSr2tLpdJbw8XicgoKC7BTOqVo2Y0FzXerq6iZk2WiBMU2hZ6x2VEVR2LdvHx6PZ8zMhaqq3HfffRw7doxHH300pwq9g9HV1UVXVxfnnXce4XCY1atX89RTT7Fs2djdeMMwc0muIR6P8/zzz7Nz50727dvHZZddxnXXXcezzz7LTTfddMocstGQSqWyJ3wqlcLtdlNcXDxiUCzXQxQGQyuv1ev1LFq0aMgEF+1h1eIMY+nAjwbtFNQsmxcefYV3n6+ntLqY7V/bhNUxdc14rdxT8xHz8/Ozaqe5OOn6+vo4cuTIhAk+EgZrtaXT6Wyjis1myxYujWWVqarK/fffz969e3n88cdPVy55Uti8eTN33HEHGzZsmMyfzXySD0YymeSpp57iK1/5CkVFRaxatYqPfexjXHzxxZP6MkYKihUXF+NwOLLFHDU1NTkdogAZU7GhoQGXyzXu9I7BbbLD5aRGgzZQoLa2Nif+/UhrLygoYO7cudlBAVqBkM1my5aIToUYvb29HDt2bEzdvslisB/f09OT1bcfzY9XVZUHH3yQ1157jf/5n//JSd3DRNHS0sKll15KY2PjZAu3zi2SA3znO9+htraWjRs38te//pVdu3bx+uuvc8EFF7BlyxYuu+yySX05sixnCR8IBBAEgaVLl2bLM3MFrcGmoqJi0vlVze0Y3JWmVYtp6OrqoqOjg9ra2pw/nJrCTXFx8Yhm7nArZLRW1NGgbU6rVq3K+dq1sUj5+fnY7fYhfvzgTUlVVR5++GFeeOEFdu3albONZiKIRCJcdtllfOtb3+JjH/vYZP/83CP5SJAkiVdffZUdO3bw8ssvs2rVKrZs2cL69esn7NcdOXKEZDJJSUkJfr+fUCg0oQaViUDz7yeisz4ehvvIWsFGOBymtrY25/6jZqJPZnMaHrjTCD+SazRaGWwuMJjgg3sLhvvxv/nNb0gmk7S3t/PCCy/k3AoaC+l0muuvv56rr76au+++eyqX+HCQfDBkWWb37t3s3LmTv/zlLyxbtowtW7awYcOGEQNuWp24zWYbovE2uEGlv78fp9OZbVCZDOE1BZrTrZ8fCVo5Zn9/P3q9PmebkgZNFeV0RBFTqVTWUkokEtlct9PpxO/309LSMmIZ7OlCS33m5eWN6xr94he/4L//+7/Jz88nEonw0ksv5bzwZiSoqsqnPvUpCgoK+NnPfjbVy3z4SD4YiqLwzjvvsGPHDl588UUWLFjApk2buOaaa3A4HITDYQ4fPkxZWdm46ZRgMIjX6yUQCEx4qIKmEqPlenMJLY2lKApLly5FVdVsai4YDI7YNTcZaCW2Y+maTRZavYCWTVAUhcWLF1NUVJTTHPRkRhPv2LGDRx55hOeeew673U4kEsl56exoeO2117jkkkuGFBL98Ic/5Nprr53MZT7cJB8MRVGor69n586d/OlPf6KwsJDm5mZ+85vfTGp0zfChCqOlvabTR1YUJVu+O5LC7OCuud7e3mweeaKpOa1fe7zhAVOF1shSVVVFX19fNnA30Vz3WFBVlcbGRux2+7gEf+qpp3jwwQd59tlnc55FOYOYJflIeP311/n0pz/N+vXrefvtt3G73WzZsoXrrrtuUqeW5ttp5bXaFJV4PJ4t9cy1jzw8yj3RNU40Nae1c06HewGZza+zs5O6urrsZ3O6gTsNqqpy4MABrFbruOOAn3vuOX7605/y3HPPTctGdtttt/Hss89SVFREY2Njzq8/CLMkHwmvv/46c+fOpaKiImv27ty5k2eeeYa8vDw2bdrExo0b8Xg8k4qwRyIRDh48SCwWy5rLuWpBhfcr8MrKyqYsIT1Sak4b0aypzE5Htxdk2n+18VZjuRCa6KZWETi44m40qKrKwYMHMZvNVFePrRv3wgsv8KMf/Yg//vGP06am+sorr2C327n11ltnSf5BglaJtmvXLp5++mmMRiObNm1i8+bNY47Egff9QKvVSnV19ZAHVRTFIfPOp4JkMpkdHJCrySCaYKTP5yOZTJJOp1m+fHnOdcsh08Hn8/mora2dVIxgpMDd8EYfjeAmk2ncARl//etfuffee3nuuedyPmFlOFpaWrj++utnSf5BhaqqtLW1sWvXLp566ikURWHjxo1s2bKFioqKIQ+SlkcuKioasZoqkUhky2tVVc2W10602EZLwS1atCjnEtLw/qTV0tJSgsEgiUQiWymWi5ns7e3t9PT0nHYr6vBGH63izuv1YjQax20vfvXVV/nmN7/Jc889N+khiFPBLMlnEDSxil27dvHkk08Sj8e57rrr2Lx5MzqdjnfeeYdLL710QjJQqVQqS3hJkrIn/GimqGZCT5eP7Pf7OXHixJBKs8EFQpFIJJuac7lckyZ8W1sbfX191NTU5DyCHggEsvUNBQUFYwYX33jjDb7yla/w7LPP5rwVeTTMknwGw+fz8eSTT/Lb3/6WI0eOcMMNN/CZz3xmSC36RKAVtni9XpLJZJbwmu5ZMBjMtltOR97W6/XS1tY2ZgZgpNFOWhR8PNK2trZmJ7nmuk1TE7rU6XQsWLCASCQyZPjD4Gk57777Ll/4whf4wx/+QGVlZU7XMRZmST7D0dnZyXXXXcfPf/5zjh49yq5du+ju7ubqq69m69atLF26dFIPtiRJ2dMzFothsViIRqOsWrUq5zX0kAmCnTx5ckiUezxoqTmv1zsk7TXS6dnc3Ew4HGbFihXTQvCmpiYEQRhxY43FYvj9fl566SUefPBBIpEIjzzyyHSPED4FsySf4dCKTwb7yMFgkGeeeYYnnniC5uZmNmzYwJYtW6itrZ3Ug97Z2UlLSwt2u51YLJY1RfPz83NST6+16E42CDYYY6XmOjo6iEajLF++fFoIfuTIEVRVHVfnvrGxkc9//vNs27aNt956iyuuuIIvfOELOV3PaLj55pv529/+Rk9PD8XFxXz3u9+drgkrsyQ/WwiHwzz33HPs2rWLpqYmrrjiCjZv3sz5558/5oPf1taWDVLp9fqsuez1egmFQuTn52f946kQqKWlhf7+/pyb0JqyTHt7O6qqUlVVRXFxcU7rwLU58lql3FgEP3ToEP/4j//I73//+8n2Z880zJL8g4DBPfENDQ1cdtllbN68mbVr12ZPUm3UUjQaHdXE1SaAeL3eIaWrE/GPtfSgplM3HSfssWPHSKVSzJ8/P+t6aA0q4+W5J3p9SZLGVaM5cuQIt956K4899lh28sw5jJlN8m9/+9s8/fTTiKJIUVERv/71r3M+R/xMI5FI8OKLL7Jz50727NnDunXr2LRpE8899xzbt29nzZo1EzLJh/vHdrs96x8PN8E1E1eWZZYuXZpzNRfthJUk6ZTrD89za4UtkxFj1AieTqfHXX9zczO33HILv/71rydVrjwZPP/889x1113Issztt9+eFXo8S5jZJA+FQtm00QMPPMDBgwd56KGHTveyHxikUilefPFFvvjFL2K1WjnvvPPYunUrl1566aSq5IbX02ty1R6PB51Ox6FDh9DpdJOO/k/03hP1kSVJyua5I5HIhGINmgWiCUaOdf22tja2b9/Or371q1MGY+YKsiyzaNEiXnzxRSoqKlizZg2/+93vzqZLMLMnqAzOC0ej0TM+B2q6YTQas+Nv7rjjDl555RV27NjBN7/5zWxP/Ec+8pFx/VpBEHA6nTidThYsWJCVq96zZw/JZJK8vLxxCTIVaGksURQnNOxRr9dTXFxMcXFxNtbQ1dXF4cOHcTqd2aEKg12JEydOTIjgnZ2d3HzzzTz44IPTRnCAt99+mwULFmRr42+66Saefvrpc9bvPyNjkr71rW/xm9/8BqfTyV//+tczccszii984QvZh3f9+vWsX78eWZZ5/fXX2bVrF/feey/Lly9ny5YtXHnlleO2rGpD8iwWS1b+2WAwUF9fj16vH5I/Ph2oqprVqpuKoKMoirjdbtxu9xApqWPHjmGz2SguLiYcDpNIJFi+fPmY1+/u7mb79u3cf//9rFu37rTe13jo7OwcUtVYUVHBW2+9Na33PJvIibl+5ZVX0t3dfcp//8EPfsDmzZuz//7Rj35EIpEYaT7zOQ1FUXj77bfZuXNntid+y5YtXH311aM2kWiywx6PZ8gDqdXT+3y+U5pTJoPJ1IpPFprKzZEjR7Klq2N1zfl8PrZt28Z99913RvLgO3bs4M9//jO/+tWvAPjtb3/L22+/zc9//vNpv/comNk++WC0trZy3XXXTXdhwQcaWk/8jh07eP7556msrGTTpk1ce+212X5orVOtvLx8TLmm4XLVGuHHsxS0dk6LxcL8+fOnxYUaXEgzuNFHk1vWNqaenh62bdvG9773Pa655pqcr2MkvPHGG9x77738+c9/BjKHD8A3vvGNM3L/ETCzSX706FEWLlwIwM9//nNefvlldu7cOalrfPWrX+WZZ57BaDRSXV3Nf/3Xf01p1tgHDZo4wo4dO/jjH/+Ix+Phyiuv5I9//CO//OUvJ5WF0OSqvV4v6XR6SHntYGiddna7fdx+7amipaWFUCg0YhpRa/Q5duwYX//611FVlTvvvJN/+qd/mpa1jARJkli0aBEvvfQS5eXlrFmzhscff5zly5efsTUMw8wm+bZt22hqakIURebOnctDDz006eaCF154gfXr16PX6/na174GMNKM5xkNVVV59dVXueWWW6isrMRut7Np0yauv/76SffEa3LVXq93SMrLZrNx4MCBCUkqTRUTrXXv7+/nxhtvpLa2lvb2djweT9Z8PhP44x//yBe/+EVkWea2227jW9/61hm79wiY2STPNZ588kl27tzJY489draXknPcd999rF27lksuuWRIT7zJZGLjxo0T6okfDi3l1d3dTV9fHw6HgwULFuB0OnNupre1tREIBMYleDgc5oYbbuCOO+5g+/btQGaTO9cyL5PALMkHY+PGjWzfvp1PfOITZ3spZwSDe+KffPJJAK6//voRe+JHgzZayOVyYbPZ8Hq9Q4Y2TqX9dDgm2o4ajUa58cYbue222/jkJz95Wvc8h/DhIPlEIvQ/+MEPePfdd3niiSc+lLv+4J74J554gkQiwfXXX8/mzZuZN2/eiJ+Jpic3fLTQSHLVxcXFU6qn1wQlxmviicfj3HjjjXz84x/ntttum9Q9poIdO3Zw7733cujQId5++23OP//8ab/nFPHhIPl4ePTRR3nooYd46aWXci6NPBOhqmq2J/6JJ54gGAxy7bXXsnnz5mxVnJaGKyoqGnc4oCYFrclVFxcXT0gKuqOjI9sNNxbBE4kEt9xyC1u2bOFzn/vcGdmkDx06hCiKfO5zn+Pf//3fZ0k+Cj4QJH/++ee5++67efnll09rOukM2tknjd7eXp566imeeOIJvF4v69ev59VXX+VnP/vZpBo8RpKCLi4uHrGefqKab6lUik984hNcddVV3HnnnWfcCrv88stnST4GPhAkX7BgAclkMqvKuXbt2inVv8+gnf200N7eztVXX43H4yEUCmV74icr36QVtWj95mazOUt4n8+H1+sdl+DpdJp/+Id/4OKLL+bLX/7yWXGzPuwkPyNlraeLY8eO5eQ6S5cuzcl1Puior6/nnnvuYfv27dme+J/+9KccOXIk2xO/evXqcQkvCAJ5eXnk5eVl5Zl8Ph9vvfUWkiRRXV2NLMujklySJG6//XbWrFkzbQSfaLXlhxkzguSzmBw2btyY/d3hcHDTTTdx0003EYvF+NOf/sQvf/lLGhsbsz3xF1544YTUZex2O+FwGIvFwqJFi+jt7aW+vh6dTpetthssJPn5z3+eZcuW8Y1vfGPaTvD//d//nZbrnks450g+u7OPDqvVyrZt29i2bVu2J/63v/0tX/rSl1i3bh1bt25l3bp1o+rFdXV1ZTXldDodDoeDqqqqbBXb/v37CQQC7N69m/b2dubOncu99977ocyEfJAwI3zyXGMG+GhnFKlUir/85S/s2rWLN954gwsvvJAtW7ZwySWXZBtKuru76ejoGFc00ufzcdddd7F//35KSkq4/fbbz0i6bCQ8+eST3Hnnnfj9fvLz86mrq8vWq3/AMBt4yzWmQvIPmJLItEGSpGxP/Kuvvsp5552XbRm97777xiS4oijZz+WBBx4gGAzS0tLCeeedd6aWP1MxvaaOqqpj/ZxTeOKJJ9Ty8nLVaDSqRUVF6lVXXTWhv5MkSZ0/f756/PhxNZlMqjU1NeqBAwemebVnH5Ikqd/73vfUiooKta6uTr3pppvU3/3ud6rf71ej0eiQn3A4rH7pS19SP/OZz6iyLJ/tpc80jMfD0/rJrRLgBxxbt26lo6ODZDKJ1+udsOk2WEnEaDRmlUTOdagDwpSNjY3s2bOHu+66i3feeYcrrriCW2+9lSeeeIJIJIKqqnz/+9+nr6+PBx98MOcCkyPhq1/9KkuWLKGmpoatW7cSDAan/Z4zFR8qkk8VIymJdHZ2nsUVnRno9XoeeeQRnE4noiiydu1afvKTn7B3716+8Y1v0NjYyNVXX80FF1xAU1MTDz/88GnNSZsMNmzYQGNjY3a+nNYTPotTcc5F16cD6ghxiw9zxFgURVavXs3q1av54Q9/yLPPPsv69evPGMEBrrrqquzva9eunbRGwYcJsyf5BFBRUUF7e3v23x0dHTNeVjpXEEWRTZs2Tcss9InikUce4aMf/ehZu/8HHbMknwDWrFnD0aNHaW5uJpVK8fvf/55NmzZN+jq33XYbRUVFrFixYhpWee7hyiuvZMWKFaf8DI6H/OAHP0Cv1/Pxj3/8LK70g40PZQptKsiFksgrr7yC3W7n1ltv/VDr3OUK51Bn4mye/FzCGZqSec4jV52JHxBMK8lnzfVZzEjccccdhMNhNmzYQF1d3RkVgpxpmI2uz2JGIlediR8GzJ7ks5jFOY5Zks9AtLe385GPfISlS5eyfPly7r///rO9pFl8gDFL8jOIm2++mYsuuoimpiYqKip4+OGHp3QdvV7PT37yEw4dOsSbb77JL37xCw4ePJjj1Z4dfPvb36ampoa6ujquuuoqTp48ebaXNOMxG10/B7B582buuOMONmzYcLaXcto410ddj4LZ6PosRkdLSwt79+7lwgsvPNtLyQnO9VHXZwOz0fUZjEgkwrZt2/jZz342hBwzHef6qOszjVlzfYYinU5z/fXXc/XVV3P33Xef7eVMCrOjrk/BWa14m8UHEELGhn0U6FNV9YuneS0z8ApgImPZ7VRV9Z7TXmQOIAjCXOA5VVVni/1PA7M++czExcAngfWCINQP/Fw7xWslgfWqqtYCdcA1giCszdE6Jw1BEBYO+ucm4PDZWsu5glmffAZCVdXXyJGJp2ZMucjAPw0DP2fTvPs3QRAWAwrQCszWq54mZs31WSAIgg7YAywAfqGq6tfO8pJmkUPMmuuzQFVVWVXVOqACuEAQhFkf+BzCLMlnkYWqqkHgb8A1Z3cls8glZkn+IYcgCB5BEPIHfrcAVzIb7DqnMBt4m0Up8OiAXy4C/6Oq6rNneU2zyCFmA2+zmMU5jllzfRazOMcxS/JZzOIcxyzJZzGLcxyzJJ/FLM5xzJJ8FrM4xzFL8lnM4hzHLMlnMYtzHLMkn8UsznH8/5YLmzfioBJ7AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "from mpl_toolkits.mplot3d import Axes3D \n", "\n", "import matplotlib.pyplot as plt\n", "from matplotlib import cm\n", "from matplotlib.ticker import LinearLocator, FormatStrFormatter\n", "import numpy as np\n", "\n", "\n", "fig = plt.figure()\n", "ax = fig.gca(projection='3d')\n", "\n", "# Make data.\n", "x = np.linspace(-3, 3, 100)\n", "y = np.linspace(-3, 3, 100)\n", "x1, x2 = np.meshgrid(x, y)\n", "F = 3*(1-x1)**2 * np.exp(-(x1**2) - (x2+1)**2)\\\n", " - 10*(np.true_divide(x1,5) - x1**3 - x2**5)*np.exp(-x1**2 - x2**2)\\\n", " - np.true_divide(1,3)*np.exp(-(x1+1)**2 - x2**2)\n", "\n", "# Plot the surface.\n", "surf = ax.plot_surface(x1, x2, F, linewidth=0, alpha=1, cmap = 'viridis')\n", "plt.show()\n", "\n" ] }, { "cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbgAAAEoCAYAAAAqrOTwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAZF0lEQVR4nO3df6jdd33H8dc7N2mbNP0hszMmjbM40ZZOlJX+s7FfdjVzonOboBtOcBD8Q1ZhgtpAnR2VDYds6P4wYJlKpwi1OLS6tnTDCasapdpqqnTKaBK11GHbNJVyk/f+uPfkntycc8+P7+f35/mAQG5yc873fO85n2fe3/M955i7CwCA1mzLvQEAAMRA4AAATSJwAIAmETgAQJMIHACgSQQOANCkwYEzs4vM7Otm9m0z+66ZfSDEhgEAMIQNfR2cmZmki939pJntkPRVSTe5+wMhNhAAgGVsH3oBvlbIk+tf7lj/xavHAQBZBXkOzsxWzOxBSY9LutfdvxbicgEAWNbgCU6S3P20pFea2eWS7jKza9394fHvMbODkg5Kkl1wwa/veMEvh7hqAEBnnnvs2BPufsWs7xv8HNx5F2j2fknPuPs/TPueC1+03/e9+11BrxcA0Icf3fTub7r7dbO+L8RZlFesT24ys52SbpD0yNDLBQBgiBCHKF8o6RNmtqK1YH7W3b8Q4HIBAFhaiLMovyPpVQG2BQCAYHgnEwBAkwgcAKBJBA4A0CQCBwBoEoEDADSJwAEAmkTgAABNInAAgCYROABAkwgcAKBJBA4A0CQCBwBoEoEDADSJwAEAmkTgAABNInAAgCYROABAkwZ/ojeA6bbvPZXsulZP7Ep2XUANCBywpJTxmses7SGA6A2BA+ZQWsyWMe02ED60isABm7QQs0VMur1EDy0gcOheb0Gbx+Z9QvBQIwKH7hC0xRE81IjAoQtELazx/UnsUCoCh2YRtTSIHUpF4NAcwpbPaN8TOpSAwKEJRK0sTHUoAYFD1Qhb+ZjqkAuBQ5UIW30IHVIjcKhKbWG7es9Pk13X0Z+8INl1DUHokAqBQzVKjlvKkC26DaWGj9AhNgKH4pUWthJitohJ21tS9AgdYiFwKFYpYastaPPYfJtKCB6hQ2gEDkXKHbcWo7aV8dubO3bb954icgiCwKE4ueLWW9SmKSF2THMIYXDgzGy/pE9K2iPpjKTD7v5PQy8X/SFs5Rntm5yhI3JYVogJblXSX7v7t8zsEknfNLN73f17AS4bncgRN8I2v5xTHdMcljU4cO7+Y0k/Xv/902Z2VNI+SQQOc0kdN8I2TK6pjmkOiwr6HJyZvVjSqyR9bcLfHZR0UJJWnve8kFeLiqWMG2ELK0foiBwWESxwZrZb0p2S3uXuT23+e3c/LOmwJO3cs99DXS/qlSpuucP2e89/JOrl3//Ey6Ne/ixX7/kpkUORggTOzHZoLW53uPvn5vk3u45tO/v7U1eeCbEZqEircYsds3mvM3X0Uk9zRA7zCHEWpUn6uKSj7v7hZS5jPHYSwWtdirilDFuOqM2yeZtSBS/lNEfkMEuICe43JL1V0kNm9uD6n93s7ncve4FMdxgiRdxKjNpWxrc3duxSTnNEDlsJcRblVyVZgG2ZiNi1Jfb0FjtutYVtktFtSBE6IoecqnonEw5lYisx49ZC2DZLEToih5yqCtxmTHd1iTm9xYpbi2HbLHboiBxy2Tb7W+qw69i2s7/QF+IWRszbm/ulGuhT1RPcNEx25cn96QCL6C1s42JOcykmOaY4jGt+3Bmf7Jju2hN6Mug5buNi7YcUk1xN/5lCXN2t+MQO0xC3c7E/ULuuV3liVzee14kvRuSY4pBKk8/BLYOXIMRTw2LDtAK0h9FlCqY7YA1THGrFBDcHzspsW87p7bW75/vYxLtPXhN5S4D2MJ4siLMyF8dp2+d77e7vzR23Zb4/NKY41IgVeiBil0fNJ5gMCVXOyAG1YVUOiNilk/IDNkMbcriRQ5WLY4rrF8/BRcLzdgCQF4FLgNgBQHoELrFeX2+3emIXh4rGjA41zvOcWgmHJVN9IjgQEoHLjOlueUd/8oIgJ5vc/8TLs71UoIR49YA3Ye4TZ0MUpPWTVFhgAKTU5kragFZjFzpyoc6m5BDcdLH2Tc1nwqIOba2ejWo1dqUhcudjn6BmrJiVaSF2pU5xEgv6OPYFalfvKomqY0fkynX/Ey+Pvg84PIkU6lsZMVGNsSNy5Ulxu4kbUuFlAg2q6aUHoV8fF+qlA9LGYt/DZ8WlCjpxQ0oErnE1xK7kyEnthi71lErckBqB60jJsSs9clI7octx+JW4IQcC16lR7EoK3eg5uVChGy2qsUIn1RG73M8nEjfkQuA6V+JUV8M0N7I5HiUEL3fQRggbciNwOKuk2MWInBT/g1InxSVG9EqJ2CSEDaUgcJiohNiFPmQppQvduJJjFBJhQ2kIHGbK/XxdjI/ayRG6FhE1lIzAYW45p7oY05xE6JZB1FALAoel5IodoUuPoKFWBA6D5TiEGesTwscX8x5j13LM+NDT/gQJnJndLul1kh5392tDXCbqk3qqizXNjbQcu5ZDBoyEmuD+RdJHJX0y0OWhcimnutihk84PQg3BI2LoXZDAuftXzOzFIS4LbWktdCNbxSNF/IgXMBvPwSGJlIcvU4ZuEuJTLp6H60uyDw8zs4NmdsTMjpx+9plUV4sCpfrcutUTu1jMgI4lC5y7H3b369z9upWdF6e6WhQsdeiIHaR8kz3Sq+fjn9GslJ9ETuggEbleBFlVzOzTkv5b0svM7JiZ/WWIy0VfRqFjqkMKRK59oc6ifEuIywFGWj37EmXhpJO2cRYlipYjdCO9BW+ehb7FfTK6TYSuPQQOVcj1dmAjtSzssRfpRS6/ln02wjTXHgKHquw6ti3Lx/ZMWvhSLeC1Lro1TsRMc20hcKhO7s+nG2ERXExNEzHTXBt4mQCqlfLlBQirhjNYt+89VXyIsTVWB1SP0NWr9MhJ5U+bmI5VAc0gcnWqJXKErj6sCGgKkUNMRK4unGSC8+w+vtjJGyf3lRWVUk5CwfxifUJ7DJxpWQ8Ch4WDNs+/LyF6uV5SgD5wpmX5CFynhkZtkcvPGTsih5iIXNkIXGdih22r68wVOiJXvloOT07CIcty5T+OhCR2Hz+TJW6lbAMnnyC2miPdKh71Hcgdts2IHFpF5MrCIcqGlRa2cbkPWwKx8LxcOVhdGlVy3Mal3k6muPK0OPW0eJtqxKO9QbXEbaS27QXmQeTyI3CNIRZAOYhcXgSuITXHreZtx/J6CEAPt7FUBK4RBAIAzsVZlA1oJW67j5/hrMpApp1Mwwve8+DMyjxYTYCGzPpsvJLOIu3t0F1vt7cE5dzbsZRWprdUmGDKihwQE/d0AEiEKS4tAlcxprfFML2Vo+eFvufbnhqBAwA0icBViultMb1Mb73cztoxxaVB4IDGlB45FnekQuDQvNIX/Bi2us097g/0icChGLzIGz1hko2PdzKpEM+/za/naWX8to9e+9bz/kB/CByKEGN6YzHfwL5AjzgmhCaxoJeJw3JIicAhO557Q68IflysLGgO0xsAKVDgzOyAmX3fzB41s/eGuExM19LEE/q2EDcAI4NXFzNbkfTPkv5A0jWS3mJm1wy9XLSPuAGIKcQKc72kR939h+7+nKTPSHpDgMtFw1qaQgGUKcQqs0/SY2NfH1v/s3OY2UEzO2JmR04/+0yAq0WteEkAgBRCrDQ24c/8vD9wP+zu17n7dSs7Lw5wtX2rdQIibgBSCfFC72OS9o99faWkEwEuF41pIW5bnda9emJXwi0BMEuIwH1D0kvN7CpJxyW9WdKfBbhczHBy37Zq3rar1rgt8jql8e8ldkB+gwPn7qtm9k5J/y5pRdLt7v7dwVuGZtQWtxAvvh1dBqED8gnyXpTufreku0NcFhZT+hRXU9xivKsEoQPyqfNMBZyjxBNOTu7b1n3cNl8+b8sEpFXeyoillBK5WGGT4sQtdXh6jxyT7LnYH3GVsSoiiNyRi3n9oeOWc6LqPXJAKnweXGNGkUn1vFyKqMaIW27b957if+9AZASuUbFDl2paDBm3EsI2jsgBcRG4xoUMXepDoC3HbYTI9Yufe3wErhOT4rRV9HI/n9dD3EZ6i9zqiV3F/0zQBgLXsdwRm6SnsI3rLXK942edRnkrHLrVa9x6xAKPFAgcikDc6t1uLIa4p0PgkB1x6xMLPWLjOThkQ9jOx3NxbeNnmxYTHLIgbpBY8BEXgUNyxG1rLd6mrfQSuV5uZ0k4RIlkWnzLLWAexC0PJjgkQdywlZYD0PJtKx0THKKK9RE3aM8oBC39fIlbXkxwiIa4YRmtRKGV21EzAofgTl15hrhhkNrjUPv2t4LAIagYYZOIW49WT+yqMhQ1bnOreA4OQcQKm9Rf3Fggz1XLpw/wcysPExwGI26IrfRpruRt6xkTHJYWM2xSn3Fjodza+P4p4f7Bz6tsBA4Lix02pLXr2PQDOSX/rHO9rICo1YPAYSGpFrwS/neeWqqFc6ugzfreEoO3eb/FuO8QtToROMwl5cJG3OJYJGyzLqPE0I1M25fz3q+IWTsIHLZU8kLWipgLaoiozbrcWu4jhKs/BA5T5Vi4epveYi26scI267pqiR36QOBwHhap+FoI21bXz30IJSBwOCv3otTL9BYjbrnDthmhQwkIHFiEEukhbJtx+BI5EbiOseCk0WPYJmGqQ2oErkMsMOmEjluNYduM0CEVAtcRFpQ0mNjmQ+gQ26DAmdmbJP2NpKslXe/uR0JsFMJiAUmDaW05PE+HWIZOcA9L+mNJHwuwLQiMxSINwhYOUx1CGhQ4dz8qSWYWZmsQRK2LQy2f+yVxGDI2QocQeA6uIS0sBiVHrvYXZ+8+vtz94+S+fOHl8CWGmBk4M7tP0p4Jf3XI3T8/7xWZ2UFJByVpx6XPm3sDMVtrD/ySIldj1JYN2SKXlyN6THVY1MzAufsNIa7I3Q9LOixJO/fs9xCX2bPWH+S5Ihf7DXljhC100Ba9ztSxI3SYF4coK9PTgzrFB1rW8jE1m+WI2jS5YsfhS8wy9GUCb5T0EUlXSPqimT3o7q8JsmU4R88P4EkRWiR6uT4mJXTYSoraNLlj1/PjBOcbehblXZLuCrQt2IQH63SlfrZXj1GbZrTthA65cIiyMDww6xQybDVHbZKcoZN4TPWMwBWCB2F9iNpicoROYqrrGYHLiAdcnUKFrYeoTULokAqBS4wHV51qnNYu+dGzC/+bp6/aGWFLJssdOonHY+sIXAI8iOpUywkjy4Rs3stKEbzdx89ke7cUprq2EbhIeMDUq/RpLWTQFr2uWMHLNc2NMNW1icAFwoOibqVPaymjtpXRdsQMXc73vpSY6lpC4AbgAVC30t9hpJSoTRIzdLmnuRFCVz8CtwDu6PUjamHFDl3uyEkcvqwZgZuBO3QbCFtcl/zo2aYjN8JUVxcCtwl33HbE+kiaUGGrPWqbxZrmSoucxFRXi+4Dx52zPaV/1lprYdssxjRXYuRGmOrK1VXguAO2K/anYhO2xfQWOYmprkRNB447WftqCJvUV9xGeozcCFNdGZoKHHemPsSOmlRn2Lb/z/GZ37P6kn0JtmRDrJNPasFUl1e1gePO0p8UYZPqids8QZv1b1IHL4RaprjNmOrSqyZw3Cn6lCpqUh1hWyZqsy4vduR6PlQ5CVNdOkUGjh86UoZNKj9uocM26bJrnOZqx1QXV/bA8YPFuNRhk/qO2+briRU5pritEbo4sgTuzAX8IHGuHGGTyo7bImHb7qd1y5NfliR98NIbdfNT90iSbr3sgFZtJfi2IQ4OX4aVfYJD33KFTSo7bou65ckv69eeOyFJ+tTPPqXtfvrsn99y+R/m3DQsialuOAKHLHKGTSo/bsselrxIp6X1uP1Ci09uKU46wWKY6pbXxgFsVGPXsW3NxC2WZeL2wUtvPO9Q5Kqt6LbLblzocmqLW+k/y9BKePzUhD2FZEp4YLb6CQA3P3XP2cOSI9v9tA49eU+mLUJMhG4+HKJEdDwQ0/mFVrRqK+fFbh61TW/g8OUsrDyIprT/ZdZyOGuZ0Nx62QE9dMFePXTBXr31l9569ve3XnYg2nUuItbbdbXyMoEQSnu8lYAJDlG0/kAr6fCktPZ82/jZkoucOcnk1hamug0EDkGVGrZapreR1Zfsi/4C71Rh6/nNlnPr/aUGBA7BlBq3Wo0CFDp0TGz96TV0BA6D9Ri2p6/ameww5XiQlo1drqjFnN54/m1xvR2+JHAYpMe45VTT9MWhybL1MNWxOmFpvceNBXy62PuG6S2cls++bPNWIbpWHxAYjvDXqcXQtXVrkERrD4IhWMzPlWJ/ML3F1VLo2rgVSKaVO35IRG5tHxC3trQQukFbb2YfMrNHzOw7ZnaXmV0eaLtQoJrv7LEXxp4j1/Nt70HNoRu61fdKutbdXyHpB5LeN3yTUKJa7+AppZpiSpH69jK95VVj6Aa9TMDdx9+q/AFJfzpsc1Ci2u7U05zcty3JO5qMFv3S3s4rlBwRJ27lqOm1dCHvNW+X9KVpf2lmB83siJkdOX3ymYBXC5SptWku14RK3MpV+n9+Z05wZnafpD0T/uqQu39+/XsOSVqVdMe0y3H3w5IOS9KFL9rvS20tkiv9DryoVFPcyHgQapzockeauJWv5BeMzwycu9+w1d+b2dskvU7Sq92dcDWktbiNpI7cSC2HLnNHTSJsNSoxdIOegzOzA5LeI+m33f1UmE0C4ssVOen8gJQQvBKiNkLc6lZS6Ia+F+VHJV0o6V4zk6QH3P0dg7cK2bU6vY3LGblxqYNXUsw2I27tKCF0Q8+i/NVQGwLkUErkxpUcoFgIW7tyho57Fc7Tw/Q27uS+bSywGbHv+5BjXeGeBaxjoU2L/1j0J/WLxfk8OJyjt+lts9GCW9phy5YQNaQ6bMk9DZiARTg8JjZsFvs/1ExwwBTjizET3fKIGrYSc5ojcDir98OTW+HQ5WKIGha169i24JEjcMACCN3WCBuGCD3NEThgCRy+3EDUEFqoaY7AAQP1FjuChhRCTHMEDgho8+LfQvAIGnIaMs0ROCCiGoNH0FCaZSNH4ICEpsUjdfiIGGqzzCFLAgcUgOAA4fGoAgBUY5HX6xI4AECTCBwAoEkEDgDQJAIHAGgSgQMANInA4azYHz4IACkROABAkwgcAKBJBA4A0CQCh3PwPByAVhA4AECTCBwAoEkEDufhMCWAFhA4AECTCBwmYooDUDsCBwBoEoHDVExxAGpG4AAATSJw2BJTHIBaETgAQJMGBc7M/tbMvmNmD5rZPWa2N9SGoRxMcQBqNHSC+5C7v8LdXynpC5JuGb5JKBGRA1CbQYFz96fGvrxYkg/bHAAAwtg+9ALM7DZJfyHpSUm/O3iLUKxTV57RrmM8bQugDua+9dBlZvdJ2jPhrw65++fHvu99ki5y9/dPuZyDkg6uf3mtpIeX2uL2PF/SE7k3ohDsizXshw3siw3siw0vc/dLZn3TzMDNy8x+RdIX3f3aOb73iLtfF+SKK8e+2MC+WMN+2MC+2MC+2DDvvhh6FuVLx758vaRHhlweAAChDH0O7u/M7GWSzkj6X0nvGL5JAAAMNyhw7v4nS/7Tw0OutzHsiw3sizXshw3siw3siw1z7Ytgz8EBAFASzvkGADQpW+B4m68NZvYhM3tkfX/cZWaX596mHMzsTWb2XTM7Y2Zdni1mZgfM7Ptm9qiZvTf39uRiZreb2eNm1v3Licxsv5n9h5kdXX983JR7m3Ixs4vM7Otm9u31ffGBLb8/1yFKM7t09E4oZvZXkq5x9y5PUjGzGyXd7+6rZvb3kuTu78m8WcmZ2dVaO2HpY5Le7e5HMm9SUma2IukHkn5f0jFJ35D0Fnf/XtYNy8DMfkvSSUmfnOelRy0zsxdKeqG7f8vMLpH0TUl/1On9wiRd7O4nzWyHpK9KusndH5j0/dkmON7ma4O73+Puq+tfPiDpypzbk4u7H3X37+fejoyul/Sou//Q3Z+T9BlJb8i8TVm4+1ck/V/u7SiBu//Y3b+1/vunJR2VtC/vVuXha06uf7lj/dfUdmR9Ds7MbjOzxyT9uXij5pG3S/pS7o1AFvskPTb29TF1upBhMjN7saRXSfpa5k3JxsxWzOxBSY9Lutfdp+6LqIEzs/vM7OEJv94gSe5+yN33S7pD0jtjbktus/bF+vcckrSqtf3RpHn2Q8dswp91e2QD5zKz3ZLulPSuTUfAuuLup9c/weZKSdeb2dRD2IPfbHnGhtww57f+q6QvSpr4PpYtmLUvzOxtkl4n6dXe8Gs3FrhP9OiYpP1jX18p6USmbUFB1p9vulPSHe7+udzbUwJ3/7mZ/aekA5ry3sY5z6Lkbb7WmdkBSe+R9Hp3P5V7e5DNNyS91MyuMrMLJL1Z0r9l3iZktn5ixcclHXX3D+fenpzM7IrRWeZmtlPSDdqiHTnPorxT0jlv8+Xux7NsTGZm9qikCyX9bP2PHujxjFIze6Okj0i6QtLPJT3o7q/JulGJmdlrJf2jpBVJt7v7bXm3KA8z+7Sk39HaO+j/VNL73f3jWTcqEzP7TUn/Jekhra2XknSzu9+db6vyMLNXSPqE1h4f2yR91t1vnfr9DR8NAwB0jHcyAQA0icABAJpE4AAATSJwAIAmETgAQJMIHACgSQQOANAkAgcAaNL/A75Kg3LFxMvdAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig1, ax = plt.subplots(constrained_layout=True)\n", "contour = ax.contourf(x1, x2, F,cmap = 'viridis')\n", "plt.scatter(0.23, -1.62,c='r',marker='X')\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# put your solution here\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Question 2. (5pts) Regression through the normal equations \n", "\n", "We consider the simple regression problem below, similar to the one discussed in class. \n", "Find the model that minimizes the sum of squares loss\n", "\n", "\\begin{align}\n", "\\ell(\\boldsymbol \\beta) = \\frac{1}{2N}\\sum_{i=1}^N (t_{\\text{noisy}}^{(i)} - (\\beta_0 +\\beta_1 x^{(i)}))^2\n", "\\end{align}\n", "\n", "__using the Normal Equations__. To do this:\n", "\n", "- Start by building the matrix $\\tilde{\\boldsymbol X}$ with\n", "\\begin{align}\n", "\\tilde{\\boldsymbol X} = \\left[\\begin{array}{cc}\n", "1 & x^{(1)} \\\\\n", "\\vdots & \\vdots \\\\\n", "1 & x^{(N)} \\end{array}\\right]\n", "\\end{align}\n", "- Then compute the matrix $\\tilde{\\boldsymbol X}^T\\tilde{\\boldsymbol X}$ and the vector $\\tilde{\\boldsymbol X}^T\\boldsymbol t$ where $\\boldsymbol t = [t_{\\text{noisy}}^{(1)}, \\ldots , t^{(N)}_{\\text{noisy}}]^T$\n", "\n", "- Finally solve the equations through\n", "\n", "\\begin{align}\n", "\\boldsymbol \\beta_{\\text{OLS}} = \\left(\\tilde{\\boldsymbol X}^T\\tilde{\\boldsymbol X}\\right)^{-1}(\\tilde{\\boldsymbol X}^T\\boldsymbol t) \n", "\\end{align}\n", "\n", "using the function __np.linalg.inv__ from the linear algebra package. Plot the result in green on top of the plot below and compare with the true (blue) (unknown) model. \n" ] }, { "cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAg+UlEQVR4nO3deXhU5f3+8feTECCBsO/LJOwBEkGIIFhXVBARFL/8ao0WrTXa1pZa+22FIIsK1qVWbV3Kj7pUU4tVIogU933DoEA2CFtI2MMSErIn83z/mNgCBgkyM2eW+3VdXmHOTE4+I83d43PO3MdYaxERkcAV4fQAIiLy3RTUIiIBTkEtIhLgFNQiIgFOQS0iEuCa+WKnnTp1svHx8b7YtYhISFqzZs1+a23nxp7zSVDHx8eTmZnpi12LiIQkY8z2Ez2npQ8RkQDXpKA2xtxujMkxxmQbY140xrT09WAiIuJx0qA2xvQEfgUkW2sTgUjgGl8PJiIiHk1d+mgGRBtjmgExwC7fjSQiIkc7aVBba3cCDwGFwG7gsLX2zeNfZ4xJNcZkGmMyi4uLvT+piEiYasrSR3tgCtAH6AG0MsZcd/zrrLWLrLXJ1trkzp0bvcJERCQ0padDfDxERHi+pqd7dfdNWfq4GNhmrS221tYCS4GxXp1CRCRYpadDaips3w7Wer6mpno1rJsS1IXA2caYGGOMAcYBeV6bQEQkmKWlQUXFsdsqKjzbvaQpa9RfAC8DXwFZDd+zyGsTiIgEs8LCU9v+PTTpk4nW2rnAXK/9VBGRUOFyeZY7GtvuJfpkoojI6ViwAGJijt0WE+PZ7iUKahGR05GSAosWQVwcGOP5umiRZ7uX+KSUSUQkrKSkeDWYj6cjahERL3hv4z4efivfJ/tWUIuInIbismp++eLX3PjMl6zM2k1FTZ3Xf4aWPkREvge327Iks4j7VuZRVevm9osHcusFfWnRLNLrP0tBLSJyijbvK2Pm0iy+LDjE6D4dWDg1iX6dW/vs5ymoRUSaqKq2nife38KT728mpnkzHrj6DKYl98LzoW3fUVCLiDTBZ1sOkJaRxdb95Vw5vAezJw2hU+sWfvnZCmoRke9wqLyGhSvz+NeaHfTuEM1zPxnF+QP92xCqoBYRaYS1lmVrd3HPilxKKmu59fx+zBg3gOjm3j9ZeDIKahGR42w/UM7sV7P5aNN+hvVux/NXJTGkRxvH5lFQi4g0qK13s/ijbTzydj5RkRHMnzyU686OIzLCtycLT0ZBLSICfF14iJlLs9iwp4zxQ7syb/JQureNdnosQEEtImGurKqWB9/YyPOfb6drbEv+ev1Ixg/t5vRYx1BQi0jYWpW9h3nLc9hbVsX0MfHccelAYltGOT3Wt6jrQ0SC1/e8qezuw5Xc/PdMbn1hDe1iolj6s7HMmzw0IEMadEQtIsHqm5vKfnO/wm9uKgsnrBytd1ue/6yAB9/YSL213HlZAjf9oA9RkYF9zGqstV7faXJyss3MzPT6fkVE/iM+vvFbYMXFQUHBtzbn7iplZkYW64pKOG9gZ+6dkoirY8y3v98hxpg11trkxp7TEbWIBKcm3lS2sqaeR97OZ/HH22gXHcWj1wxn8rAePu/n8CYFtYgEpybcVPaD/GJmv5pF0cFKfpjcm5kTE2gX09yPQ3qHglpEgtOCBceuUcN/bipbXFbNPStyWb5uF307t2JJ6tmM7tvRuVlPk4JaRILTNycM09I8yx0uF+57F/DSgHO57+EPqKypZ8a4Afz8wn4+KfP3J51MFJGQsHnfEWZlZLF620FG9enAwquS6N/Fd2X+3qaTiSISsqrr6nnivS08+f4WoptHcv/VSUwb2ZsIh/s5vElBLSJB6/OtB5iVkcXW4nImD+vBXZOG0DnWP2X+/qSgFpGgU1JRw30rN7Aks4he7aN59sazuGBQF6fH8hkFtYgEDWsty9d5yvwPVdRyy3l9mXHxAGKah3aUhfa7E5GQUXSwgrRXs/kwv5hhvdry3E9GMbRHW6fH8gsFtYgEtNp6N09/vI0/vZ1PpDHMvWIIPx4T73iZvz8pqEUkYK0tKmHm0izydpdyyZCuzJ88lB7tAqPM358U1CIScI5U1/HQGxt57rMCusS24KnrRjIhMbDK/P1JQS0iAeXNnD3MXZ7DntIqrj87jt+OH0SbAO2J9hcFtYgEhD2Hq5i7PJs3cvaS0C2Wx1NGMMLV3umxAkKTgtoY0w5YDCQCFviJtfYzH84lImGi3m1J/2I7D6zaSG29m99NGMTN5/YN+DJ/f2rqEfWjwCpr7f8YY5oDgdO2LSJBK293KTOXZrG2qIRzB3Ti3isTievYyumxAs5Jg9oY0wY4D7gBwFpbA9T4diwRCWWVNfU8+s4mFn+0lTbRUTzyw+FMGR5cZf7+1JQj6r5AMfCMMWYYsAaYYa0t9+lkIhKSPtpUTFpGNoUHK5g2shezJg6mfavgK/P3p6YsAjUDRgBPWmvPBMqBO49/kTEm1RiTaYzJLC4u9vKYIhLsDhyp5vYla7n+b6uJjDD84+bRPDhtmEK6CZpyRL0D2GGt/aLh8cs0EtTW2kXAIvD0UXttQhEJatZa/rVmBwtX5lFeXcevLurPzy/sT8uo4C7z96eTBrW1do8xpsgYM8hauxEYB+T6fjQRCXZbiz1l/p9vPUhyXHvum5rEgK6xTo8VdJp61ccvgfSGKz62Ajf6biQRCXY1dW6e+mALf3lvMy2aRbDwqiSuOSu0yvz9qUlBba1dCzR6ixgRkaN9WXCQmUuz2LzvCJPO6M6cK4bQJbal02MFNX0yUUS84nBFLX9YtYEXVxfSs100z9xwFhcmhG6Zvz8pqEXktFhrWbF+N/Nfy+VgeTU3n9uH2y8ZGPJl/v6kf5Mi8r0VHaxgzrJs3ttYTFLPtjx741kk9gyPMn9/UlCLyCmrq3fzzCcFPPxWPsbAXZOGMH1MHM3Uz+ETCmoROSXrd3jK/HN2lTIuoQt3X5lIzzAs8/cn/d+fiDRJeXUdd7+Wy5WPf0JxWTVPuMpZPG8aPTu0gvh4SE93esSQpSNqETmpt3P3MmdZNrtLq0gZ7eJ3petp87NUqKjwvGD7dkhN9fw5JcW5QUOUglpETmhvaRXzlufw7+w9DOzampevHcPIuA4Qf8V/Q/obFRWQlqag9gEFtYh8i9ttSV9dyAP/3kB1vZv/He8p82/erGG1tLCw8W880XY5LQpqETnGxj1lzFy6nq8KSzinf0cWXJlEfKfjyvxdLs9yx/FcLv8MGWYU1CICQFVtPY+9s4lFH24ltmUz/jhtGFNH9Gy8zH/BAs+a9NHLHzExnu3idQpqEeGTzfuZlZHF9gMVXD2iF2mXD6bDd/VEf7MOnZbmWe5wuTwhrfVpn1BQi4SxA0eqWfB6Hku/3kl8xxjSfzqac/p3ato3p6QomP1EQS0Shqy1vPLVTha8nktZVR23Xdif2y5SmX+gUlCLhJlt+8tJy8ji0y0HGNlQ5j9QZf4BTUEtEiZq6tz89YMt/LmhzP/eKxO5dpRLZf5BQEEtEgYyG8r8N+07wuVJ3Zl7xRC6tFGZf7BQUIuEsMOVtdy/agP/+MJT5v+36cmMG9zV6bHkFCmoRUKQtZaVWXuY91oOB45Uc9MP+vCbSwbSqoV+5YOR/tZEQsyOQxXMWZbDuxv2kdizDU9PP4ukXirzD2YKapEQUVfv5tlPC/jjm/kAzL58MDeMjVeZfwhQUIuEgKwdh5mZsZ7snaVclNCFu6cMpVf7GKfHEi9RUIsEsfLqOh5+K59nPtlGx9YtePzaEUxM6tZ4P4cELQW1SJB6J28vc5blsLOkkmtHu/j9hATaRkc5PZb4gIJaJMjsK61i/mu5vJ61mwFdWvPyrWNIju/g9FjiQwpqkSDhdlv+sbqQ+1dtoLrOzR2XDOSW8/v9t8xfQpaCWiQI5O8tY+bSLNZsP8SYvh1ZcFUifTu3dnos8RMFtUgAq6qt5y/vbuavH26hdYtmPDRtGFefqMxfQpaCWiRAfdpQ5l9woIKpZ/Yk7fLBdGzdwumxxAEKapEAc7C8hgWv5/HKVzuI6xjDCzeN5gcDmljmLyFJQS0SIKy1ZHy9k3tWeMr8f3FhP3550QCV+YuCWiQQFOwvJ+3VLD7ZfIAzXe24b2oSCd3aOD2WBAgFtYiDaurc/P+PtvLYO5toHhnBPVOGkjI6TmX+cgwFtYhD1mw/xKylWWzcW8Zlid2YN3koXVXmL41QUIv4WWlVLQ+s2kD6F4V0b9OSxT9O5uIhKvOXE2tyUBtjIoFMYKe1dpLvRhIJTdZaVmXvYe7yHPYfqeaGsfHccekgWqvMX07iVP4XMgPIA3SGQ+QU7SqpZM6ybN7O28eQ7m1YPD2ZM3q1c3osCRJNKgkwxvQCLgcW+3YckdBS77Y8/fE2Ln74Az7ZfIBZExNYfts5pxfS6ekQHw8REZ6v6elemlYCVVOPqB8BfgfEnugFxphUIBXA5XKd9mAiwS5752FmZWSxfsdhLhjUmXumJNK7w2mW+aenQ2oqVFR4Hm/f7nkMkJJyevuWgGWstd/9AmMmAROttT83xlwA/PZka9TJyck2MzPTa0OKBJOKmjr+9FY+T39SQPuYKOZeMZRJZ3T3Tj9HfLwnnI8XFwcFBae/f3GMMWaNtTa5seeackR9DjDZGDMRaAm0Mca8YK29zptDioSC9zbuY3ZGNjtLKvnRqN7cOWEwbWO8WOZfWHhq2yUknDSorbUzgZkARx1RK6RFjrKvrIq7X8tlxfrd9O/SmpduGcOoPj4o83e5Gj+i1nJjSNN1QSKnwe22LMks4r6VeVTVurn94oHcekFfWjTzUT/HggXHrlEDxMR4tkvIOqWgtta+D7zvk0lEgszmfZ4y/y8LDjG6TwcWTk2in6/L/L85YZiW5lnucLk8Ia0TiSFNR9QiTZWeDmlpVO3czRPjf8qTSZcR07I5D1x9BtOSe/mvzD8lRcEcZnSzNZGmaLgs7nN3LBNveIzHEidyed7HvBO/n/93Vm/dcUV8SkfUIk1waP5CFp5/M/864xJ6l+zhuZfmcP62r2BtHEy/1unxJMQpqEW+g7WWZWt3cc/4WZREx3Lr5/9ixif/JLqu2vMCXRYnfqCgFjmBwgMVpL2axUeb9jOsqoTnl9zFkOJtx75Il8WJHyioRY5TW+9m8UfbePSdfJpFRDB/8lCu21pC5D/2HvtCXRYnfqKgFjnK14WHmLk0iw17yhg/tCvzJg+le9toGBsPBl0WJ45QUIsAZVW1PPjGRp7/fDtdY1vy1+tHMn5ot2NfpMvixCEKagl7q7L3MG95DnvLqpg+Jp47Lh1IbEsv9nOInCYFtYSt3YcrmbMsh7dy95LQLZYnrxvBma72To8l8i0Kagk79W7L858V8OAbG6m3ljsvS+CmH/QhKlKf/5LApKCWsJK7q5SZGVmsKyrhvIGduXdKIq6Op1nmL+JjCmoJC5U19TzyTj6LP9pGu+goHr1mOJOH9dBHvyUoKKgl5H2QX8zsV7MoOljJD5N7M3NiAu1imjs9lkiTKaglZBWXVXPPilyWr9tF386tWJJ6NqP7dnR6LJFTpqCWkGOt5aXMIhau3EBlTT0zxg3g5xf2812Zv4iPKaglpGzed4RZGVms3naQUX06sPCqJPp38XGZv4iPKaglJFTX1fPk+1t44r0tRDeP5P6rk5g2sjcRETpZKMFPQS1B74utB5iVkcWW4nImD+vBXZOG0Dm2hdNjiXiNglqCVklFDfet3MCSzCJ6d4jm2RvP4oJBXZweS8TrFNQSdKy1LF+3i3tW5HKoopZbzu/Lr8cNJLq5ThZKaFJQS1ApOlhB2qvZfJhfzLBebXnuJ6MY2qOt02OJ+JSCWoJCbb2bpz/exp/ezifSGOZdMYTrx8QTqZOFEgYU1BLw1hWVcOfSLPJ2l3LJkK7MnzyUHu2inR5LxG8U1BKwjlTX8dAbG3nuswK6xLbgqetGMiGx28m/USTEKKglIL2Zs4e5y3PYU1rFj8+O47fjB6nMX8KWCngloOw5XMUtz2eS+vwa2kZH8crPxjL/yDpiEwZARATEx0N6utNjiviVjqglINS7LelfbOeBVRuprXfz+wkJ/PTcPkT980VITYWKCs8Lt2/3PAbdv1DChrHWen2nycnJNjMz0+v7ldCUt7uUmUuzWFtUwrkDOnHvlYnEdWzleTI+3hPOx4uLg4ICf44p4lPGmDXW2uTGntMRtTimsqaeR9/ZxOKPttImOopHfjicKcOPK/MvLGz8m0+0XSQEKajFER9tKiYtI5vCgxVMG9mLWRMH075VI2X+LlfjR9Qul++HFAkQCmrxqwNHqrn39Twyvt5J306tePHmsxnT7zvK/BcsOHaNGiAmxrNdJEwoqMUvrLX8a80OFq7Mo7y6jl9d1J+fX9ifllEn6ef45oRhWppnucPl8oS0TiRKGFFQi89tLfaU+X++9SBnxbdn4VVJDOga2/QdpKQomCWsKajFZ6rr6nnq/a08/t5mWkRFcN/UJH6YrDJ/kVN10qA2xvQG/g50A9zAImvto74eTILblwUHmbk0i837jjDpjO7MuWIIXWJbOj2WSFBqyhF1HXCHtfYrY0wssMYY85a1NtfHs0kQOlxRyx9WbeDF1YX0bBfNMzeexYUq8xc5LScNamvtbmB3w5/LjDF5QE9AQS3/Ya1lxfrdzH8tl4Pl1dx8bh9uv2QgMc21uiZyuk7pt8gYEw+cCXzRyHOpQCqAS9e4hpWigxXctSyb9zcWk9SzLc/eeBaJPVXmL+ItTQ5qY0xr4BXg19ba0uOft9YuAhaB5yPkXptQAlZdvZtnPing4bfyMQbmTBrC9LEq8xfxtiYFtTEmCk9Ip1trl/p2JAkG63eUMHNpFjm7ShmX0IW7r0ykp8r8RXyiKVd9GOBvQJ619mHfjySBrLy6jj++mc+zn26jU+sWPJEygssSux3bzyEiXtWUI+pzgOuBLGPM2oZts6y1K302lQSkt3P3MmdZNrtLq0gZ7eJ3ExJoozJ/EZ9rylUfHwM6XApje0urmLc8h39n72Fg19a8fO0YRsZ1cHoskbCha6fkhNxuS/rqQh749waq69387/hB3HxuX5o3042BRPxJQS2N2rinjJlL1/NVYQlj+3VkwVVJ9OnUyumxRMKSglqOUVVbz5/f3cRfP9hKbMtm/HHaMKaO6KmThSIOUlDLf3yyeT9pGVkUHKjg6hG9SLt8MB0aK/MXEb9SUAsHjlSz4PU8ln69k/iOMaT/dDTn9O/k9Fgi0kBnhcKYtZaX1+zg4oc/YPm6Xdx2YX9W/fq8/4Z0errn5rIREZ6v6elOjisStnREHaa27S8nLSOLT7ccYGRce+6bmsTAo8v809OPvQXW9u2ex6ASfxE/M9Z6v5YjOTnZZmZmen2/cvpq6tws+nALj727mRbNIvj9hASuHeX6dpl/fHzjN5WNi4OCAn+MKhJWjDFrrLXJjT2nI+owkllwkFkZWeTvPcLlSd2Ze8UQurQ5QZl/YeGpbRcRn1FQh4HDlbU8sGoD6V94yvz/Nj2ZcYO7fvc3uVyNH1GrwlbE7xTUIcxay8qsPcx7LYcDR6q56Qd9+M0lA2nVogl/7QsWHLtGDRAT49kuIn6loA5ROw5VMGdZDu9u2EdizzY8Pf0sknqdQpn/NycM09I8yx0ulyekdSJRxO8U1CGmrt7Ns596yvythdmXD+aGsfE0i/weV2KmpCiYRQKAgjqEZO88zJ1L15O9s5SLErpw95Sh9Gof4/RYInKaFNQhoLy6jj+9lc/Tn2yjY+sWPH7tCCYmqcxfJFQoqIPcuxv2cterOewsqeTa0S5+PyGBttEq8xcJJQrqILWvtIr5r+XyetZuBnRpzcu3jiE5XmX+IqFIQR1k3G7LP1YXcv+qDVTXufntpQNJPa+fyvxFQph+uwPRCcqQ8veWMe2vnzH71WwSe7Rl1Yxzue2iAQppkRCnI+pA00gZUtXPfsHje1ry1MFoWrdoxkPThnG1yvxFwoaCOtCkpR3zacBPXWeQNv4XbCtuydQze5B2+WA6tm7h4IAi4m8K6kDTUHp0qGUsCy76CS8nXULcoV28sGQ2P/jD1w4PJyJOUFAHGOtykdG6L/de9FNKW7TiF58u4ZefLaFlz+5OjyYiDlFQB5CC/eWk3fQIn5RHMWJnHvet+guD9m9XGZJImFNQB4DaejeLPtzKY+9sonlkNPf0KCXln48TcaDQU9SvMiSRsKagdtia7YeYtTSLjXvLuCyxG/MmD6Vrm5bwq2lOjyYiAUJB7ZDSqv+W+Xdv05LFP07m4iEnKfMXkbCkoPYzay2rsvcwd3kO+49Uc8PYeO64dBCtm1LmLyJhSengR7tKKpmzLJu38/YxpHsbFk9P5oxe7ZweS0QCnILaD+rdluc+LeChNzdiLaRNHMyN53zPMn8RCTsKah/L3nmYWRlZrN9xmAsGdeaeKYn07qAyfxFpOgW1j1TUfFPmX0D7mOb8+UdnMumM7urnEJFTpqD2gfc27mN2RjY7Syr50aje3DlhMG1jVOYvIt+PgtqL9pVVcfdruaxYv5v+XVrz0i1jGNVHZf4icnoU1F7gdluWZBZx38o8qmrd3H7xQG69oC8tmkU6PZqIhIAmBbUxZgLwKBAJLLbW/sGnUwWRzfvKmLk0iy8LDjG6TwcWTk2iX+fWTo8lIiHkpEFtjIkEHgcuAXYAXxpjlltrc309XCCrqq3nife38OT7m4lp3owH/ucMpo3spZOFIuJ1TTmiHgVsttZuBTDG/BOYAoRtUH+25QBpGVls3V/OlcN7MHvSEDqpzF9EfKQpQd0TKDrq8Q5g9PEvMsakAqkALpfLK8MFmkPlNSxcmce/1uygd4do/v6TUZw3sLPTY4lIiGtKUDf23/L2WxusXQQsAkhOTv7W88HMWsuytbu4Z0UuJZW13Hp+P2aMG0B0c50sFBHfa0pQ7wB6H/W4F7DLN+MEnu0Hypn9ajYfbdrP8N7teGFqEoO7t3F6LBEJI00pm/gSGGCM6WOMaQ5cAyz37VgOSk+H+Hhqm0Xx5ISbufSh9/i6sIT5k4fyys/GKqRFxO9OekRtra0zxtwGvIHn8rynrbU5Pp/MCenpkJrK1217MfPHf2JDlz6M3/wF8y5PoPvYeKenE5EwZaz1/nJycnKyzczM9Pp+fa2sfwIP9rmQ50dMpOuRg8x/6ynGb/rcczusggKnxxOREGaMWWOtTW7sOX0yscGq7D3MGz+Tva07MH3NCu746Hliayo9TxYWOjuciIS1sA/q3YcrmbMsh7dy95JQV8WTz/+WM3fnH/uiEL3cUESCQ9gGdb3b8vxnBTz4xkbqreXOyxK4qaiEqBd2HPvCmBjPXcBFRBwSlkGdu6uUmRlZrCsq4byBnbl3SiKujjFAP89V42lpnuUOl8sT0ikpTo8sImEsrIK6sqaeR97OZ/HH22gXHcWj1wxn8rAex/ZzpKQomEUkoIRNUH+QX8zsV7MoOljJD5N7M3NiAu1imjs9lojISYV8UBeXVXPPilyWr9tF386tWJJ6NqP7dnR6LBGRJgvZoHa7LS9lFnHfvzdQWVPPjHED+PmF/VTmLyJBJySDevO+I8zKyGL1toOM6tOBhVcl0b+LyvxFJDiFVFBX19XzxHtbePL9LUQ3j+T+q5OYNrI3EREq8xeR4BUyQf351gPMyshia3E5k4f14K5JQ+gcqzJ/EQl+TWnP84+G1joiIjxf09Ob9G0lFTX8/uX1XLPoc2rq3Dx741k89qMzFdIiEjIC44i6obWOigrP4+3bPY/hhNc0W2tZvs5T5n+oopZbzuvLjIsHENM8MN6SiIi3BEZ7Xny8J5yPd4LWuqKDFaS9ms2H+cUM69WWhVOTGNqj7feeV0TEaYHfnneidrrjttfWu/nbx9t45O18Io1h3hVDuH5MPJE6WSgiISwwgtrlavyI+qjWurVFJdz5yno27CnjkiFdmT95KD3aRftxSBERZwRGUC9YcOwaNfynta6sqpY/vpnPc58V0CW2BU9dN5IJid2cm1VExM8CI6i/OWF4XGvdG8PHMffhD9lbVsWPz47jt+MHEdsyytlZRUT8LDBOJh5n9+FK5i7L4c3cvSR0i2Xh1CRGuNp7cUIRkcAS+CcTG9S7LS98vp0H39hIbb2b309I4Kfn9iEqMnAu9xYR8beACerDFbVMf2Y1a4tKOHdAJ+69MpG4jq2cHktExHEBE9RtopsR3zGGG8bGM2X4cWX+IiJhLGCC2hjDI9ec6fQYIiIBR4u/IiIBTkEtIhLgFNQiIgFOQS0iEuAU1CIiAU5BLSIS4BTUIiIBTkEtIhLgfFLKZIwpBhopmG6STsB+L44TDPSeQ1+4vV/Qez5Vcdbazo094ZOgPh3GmMwTNUiFKr3n0Bdu7xf0nr1JSx8iIgFOQS0iEuACMagXOT2AA/SeQ1+4vV/Qe/aagFujFhGRYwXiEbWIiBxFQS0iEuACJqiNMROMMRuNMZuNMXc6PY8/GGOeNsbsM8ZkOz2LPxhjehtj3jPG5BljcowxM5yeydeMMS2NMauNMesa3vN8p2fyF2NMpDHma2PMCqdn8QdjTIExJssYs9YY8/3v7t3YvgNhjdoYEwnkA5cAO4AvgR9Za3MdHczHjDHnAUeAv1trE52ex9eMMd2B7tbar4wxscAa4MpQ/ns2nnvKtbLWHjHGRAEfAzOstZ87PJrPGWN+AyQDbay1k5yex9eMMQVAsrXW6x/yCZQj6lHAZmvtVmttDfBPYIrDM/mctfZD4KDTc/iLtXa3tfarhj+XAXlAT2en8i3rcaThYVTDP84fHfmYMaYXcDmw2OlZQkGgBHVPoOioxzsI8V/gcGeMiQfOBL5weBSfa1gCWAvsA96y1ob8ewYeAX4HuB2ew58s8KYxZo0xJtWbOw6UoG7sluMhf9QRrowxrYFXgF9ba0udnsfXrLX11trhQC9glDEmpJe5jDGTgH3W2jVOz+Jn51hrRwCXAb9oWNr0ikAJ6h1A76Me9wJ2OTSL+FDDOu0rQLq1dqnT8/iTtbYEeB+Y4OwkPncOMLlhzfafwEXGmBecHcn3rLW7Gr7uAzLwLOl6RaAE9ZfAAGNMH2NMc+AaYLnDM4mXNZxY+xuQZ6192Ol5/MEY09kY067hz9HAxcAGR4fyMWvtTGttL2ttPJ7f5Xettdc5PJZPGWNaNZwgxxjTCrgU8NrVXAER1NbaOuA24A08J5hestbmODuV7xljXgQ+AwYZY3YYY25yeiYfOwe4Hs8R1tqGfyY6PZSPdQfeM8asx3NA8pa1NiwuVwszXYGPjTHrgNXA69baVd7aeUBcniciIicWEEfUIiJyYgpqEZEAp6AWEQlwCmoRkQCnoBYRCXAKahGRAKegFhEJcP8HqxLLjrOnn/8AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "x = np.linspace(0,5,10)\n", "\n", "noise = np.random.normal(0,.3,len(x))\n", "beta_true = np.random.normal(0,1,2)\n", "\n", "t = beta_true[0] + beta_true[1]*x\n", "\n", "tnoisy = t+noise\n", "\n", "plt.scatter(x, tnoisy, c='r')\n", "plt.plot(x, t)\n", "plt.show()\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# put your code here\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Question 3 [5pts]. Successive orthogonalization\n", "\n", "We once again consider the dataset from Question 2. This time we will focus on an extended model containing not only the original features but also order 2 and 3 monomials in those features. \n", "\n", "We want to rely on _successive orthogonalization_ to derive the regression vector $\\mathbf{\\beta}$. \n", "\n", "Let $\\tilde{\\mathbf{X}}$ denote the feature matrix, with first column being the vector $[1,1,...,1]$, second column being the original feature vector $x$, third column containing the degree two monomials, and so on.\n", "\n", "Starting with $\\mathbf{z}_0 = \\mathbf{c}_0$ (the first column of $\\tilde{\\mathbf{X}}$), successive orthogonalization of the columns $\\mathbf{c}_j$ of $\\tilde{\\mathbf{X}}$ can be encoded as follows\n", "\n", "1. Initialize $\\mathbf{z}_0 = \\mathbf{c}_0$\n", "\n", "\n", "\n", "2. For j=1,2,..., D\n", "\n", " Compute the coefficient $\\hat{\\gamma}_{\\ell,j}$ for $\\ell=0,..., j-1$ as\n", " \n", " $$\\hat{\\gamma}_{\\ell, j} = \\frac{\\langle \\mathbf{z}_\\ell, \\mathbf{c}_j\\rangle }{\\langle \\mathbf{z}_\\ell, \\mathbf{z}_\\ell\\rangle}$$\n", " \n", " And define the new column $\\mathbf{z}_j$ as \n", " \n", " $$\\mathbf{z}_j = \\mathbf{c}_j - \\sum_{k=0}^{j-1}\\hat{\\gamma}_{kj}\\mathbf{z}_k\\quad \\quad (*)$$\n", "\n", "\n", "Once we have all the $\\mathbf{z}_j's$ we can compute the projection of the (noisy) target vector $\\mathbf{t}_{\\varepsilon}$ onto each $\\mathbf{z}_j$ as \n", "$\\langle \\mathbf{z}_j, \\mathbf{t}_{\\varepsilon}\\rangle/\\langle\\mathbf{z}_j,\\mathbf{z}_j \\rangle$\n", "\n", "Since the $\\mathbf{z}_{D+1}$ is the only $\\mathbf{z}_j$ that contains $\\mathbf{x}_{D+1}$ this gives us the coefficient $\\hat{\\beta}_{D+1}$. I.e. \n", "\n", "$$\\hat{\\beta}_{D+1} = \\frac{\\langle \\mathbf{z}_{D+1}, \\mathbf{t}_{\\varepsilon}\\rangle }{\\langle\\mathbf{z}_{D+1},\\mathbf{z}_{D+1} \\rangle} $$\n", "\n", "The remaining $\\hat{\\beta}_j$ can be derived by using the recursion (*) on the $\\mathbf{z}_j$, projecting the target vector $\\mathbf{t}_{\\varepsilon}$ on $\\mathbf{z}_j$ and expressing $\\mathbf{z}_j$ from the recursion as a linear combination of the $\\mathbf{c}_j$\n", "\n", "Implement a function that realizes successive orthogonalization and returns the vector of regression weights $\\mathbf{\\beta}$. " ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAg9UlEQVR4nO3de3iU9Z3+8fc3IZAECJHzKcOIQDAkQDEClaqrFksREaitLlHXujV1u/7WtnatEIsngtZaT62tRWsVm9Xqz3AQEBYPeKCoBFaYHCDEkITzOSQhJCSZ7/4xsYsYYAIzeSYz9+u6uJJ5MnlyzwW58+U5fGKstYiISOiKcjqAiIicnopaRCTEqahFREKcilpEJMSpqEVEQlyHYOy0Z8+e1u12B2PXIiJhaf369Qestb1a+lhQitrtdpOXlxeMXYuIhCVjTPmpPqZDHyIiIU5FLSIS4lTUIiIhTkUtIhLiVNQiIiFORS0icq5ycsDthqgo39ucnIDuPiiX54mIRIycHMjMhNpa3+Pyct9jgIyMgHwJrahFRM5FVtb/lfSXamt92wNERS0ici4qKlq3/SyoqEVEzoXL1brtZ0FFLSJyLrKzIT6erT2SeD3t275t8fG+7QGik4kiIueg7gc38oe9sfxxdwcSj1Ux5Wg58Q/dH7ATiaCiFhE5a2u/OEjWQg+lB2KZdlF/7puSQvzvbw7411FRi4i0UmXtceYtL+L1vB0kdY9jwW1juWxYixNKA0JFLSLiJ2stSzbu4qG3Cqk81sAdl1/AXVcNJa5jdFC/ropaRMQPFQdryVrk4aOtBxiVlMgr09NI6Z/QJl9bRS0ichoNTV7+/PE2nnqnmGhjeHDqCG4aP4joKNNmGVTUIiKnsHF7JffmeijaXcXVKX148LoR9OsW1+Y5VNQiIiepqW/k8ZVbeHltGb27duK5my5iUmpfx/KoqEVETrCqcC9zFuezp6qOm8cP4hffSSYhNsbRTCpqERFgz5E6HlhSwIqCPST36cqzGWMY4zrP6ViAilpEIpzXa8n5tJxfr9hCQ5OXeyYlc/ulg4mJDp0JGypqEYlYW/ZUMyt3ExsqKvnWkJ5kT09lUI/OTsf6GhW1iEScuoYmfvfeVv70QSkJcTE8ecMopo0egDFtd8lda6ioRSSirCk5QNZCD2UHa7n+ooHMnnwh3Tt3dDrWaamoRSQiHDp6nLnLCsndsBN3j3j+60fjuGRIT6dj+UVFLSJhzVpL7oadzF1WSHVdI3deMYQ7rxxCbExw53MEkopaRMJW2YGjZC3ysKbkIGNciTwyYyTJfbs6HavVVNQiEnYamrzM/7CUZ97dSsfoKOZOS2XmWBdRbTifI5BU1CISVtaXH2Z2rocte6uZnNaX+68dQZ+EWKdjnRO/i9oYEw3kATuttVOCF0lEpPWq6hr4zYot/PXTcvolxPLCLel8O6WP07ECojUr6ruAIqBtBrCKiPjBWsvKgj3cv6SA/dX13HqJm7uvTqZLp/A5YODXKzHGDASuAbKBnwc1kYiIn3YfOcacxQWsKtxLSr8Enr8lnZEDE52OFXD+/sh5CrgHaH+nS0Uk7DR5La+sLeM3K7fgtTB78nBum3A+HUJoPkcgnbGojTFTgH3W2vXGmH86zfMygUwAl8sVqHwiIl9RuKuKWQs9bNxeyeXDejF3WipJ3eOdjhVU/qyoJwBTjTGTgVggwRjzV2vtTSc+yVo7H5gPkJ6ebgOeVEQi2rHjTTz1bjEvfLSN8+JjePrG0Uwd1T9k53ME0hmL2lo7C5gF0Lyi/sXJJS0iEkwfFO/nvkUeth86xg3pScyaPJzE+NCezxFI4XNaVETCzoGaeh5eWsjiz3cxuFdn/pY5nnGDezgdq821qqittauB1UFJIiLSzFrLG+t3kL2siGPHm7jrqqH85IoL6NSh/cznCCStqEUkpJTur2H2Qg+flB5irLs782akMqR3ZF9wFp7XsohIu3O80csz725l0tMfUbirikdnpPFa5vjTl3RODrjdEBXle5uT01Zx25RW1CLiuHVlh5iV66FkXw1TRvZjzrUp9O56hvkcOTmQmQm1tb7H5eW+xwAZGcEN3Ma0ohaRsxOA1eyRYw3MyvXw/efWcux4E3/54cX8fuaYM5c0QFbW/5X0l2prfdvDjFbUItJ657iatdayzLObB98q5GBNPbdfej4/mziM+I6tqKSKitZtb8dU1CLSeqdbzZ6hqHccrmXO4gLe27yPtAHd+MutF5M6oFvrM7hcvh8QLW0PMypqEWm9s1jNNjZ5eenvZTyxqhiAX01J4V++Oejs53NkZ391VQ8QH+/bHmZU1CLSeq1czebvPMK9uZvI31nFlcN78/C0VAYkxp1bhi9X7llZvh8QLpevpMPsRCKoqEXkbPi5mj1a38iTq4p5cc02enTpxLMzxzA5rW/g5nNkZIRlMZ9MRS0irefHavb9zfu4b1E+OyuPMXOci19OGk63uBiHArdvKmoROTunWM3uq67jobcKWbppN0N7d+GNO77Jxe7uDgQMHypqEQkIr9fyt7ztPLK8iLpGL3dPHMaPL7+Ajh10u8a5UlGLyDkr2VfNrFwP68oOM35wd+ZNT2Nwry5OxwobKmoROWt1DU38YfUX/HF1CZ07deCx60fy/YsGRsQw/7akohaRs/JJ6UFmL/RQuv8o00b3574pKfTs0snpWGFJRS0irVJZe5x5y4t4PW8HSd3jWHDbWC4b1svpWGFNRS0ifrHWsmTjLh5eWsjh2gbuuPwC7rpqKHEdI3OYf1tSUYvIGW0/VEvWonw+LN7PqKREFtyWRkr/BKdjRQwVtYicUkOTlxc/3saT7xQTbQwPTh3BTeMHER2lk4VtSUUtIi3auL2Se3M9FO2uYmJKHx66bgT9up3jfA45KypqEfmKmvpGHl+5hZfXltG7ayeeu+kiJqX2dTpWRFNRi8g/rCrcy5zF+eypquPm8YP4xXeSSYjVfA6nqahFhL1VdTywpIC38/eQ3Kcrz2aMYYzrPKdjSTMVtUgE83otOZ9V8Njbmzne5OWeScncfulgYs52mL8EhYpaJEJt2VPNrNxNbKio5FtDepI9PZVBPTo7HUtaoKIWiTB1DU387r2t/OmDUhLiYnjyhlFMGz1A8zlCmIpaJIKsKTlA1kIPZQdr+d6YgWRdcyHdO3d0OpacgYpaJAIcOnqcucsKyd2wE3ePeHJ+NI4JQ3o6HUv8pKIWCWPWWnI37GTuskKq6xq584oh3HnlEGJjNJ+jPVFRi4SpsgNHyVrkYU3JQca4EnlkxkiS+3Z1OpacBRW1SJg53ujl+Y9KeebdrXSMjmLutFRmjnURpfkc7ZaKWiSMrC8/zOxcD1v2VjM5rS/3XzuCPgmxTseSc6SiFgkDVXUNPLZiMzmfVtAvIZYXbknn2yl9nI4lAXLGojbGxAIfAp2an///rbX3BzuYiJyZtZYV+Xu4f0kBB2rqufUSN3dfnUyXTlqDhRN//jbrgSuttTXGmBjgY2PM29baT4KcTUROY1flMeYsLuCdor1c2C+B529JZ1RSotOxJAjOWNTWWgvUND+Maf5jgxlKRE6tyWtZsLaMx1duoclaZk8ezm0TzqeD5nOELb/+f2SMiQbWA0OAZ621n7bwnEwgE8DlcgUyo4g0K9h1hNm5HjbuOMJlw3qRPS2VpO7xTseSIPOrqK21TcBoY0wisNAYk2qtzT/pOfOB+QDp6elacYsE0LHjTTz1TjEvfLyN8+JjePrG0Uwd1V/zOSJEq/6vZK2tBFYDk4IRRkS+7oPi/Vz91Af86cNSrh8zkHeS9nPdtAmY6GhwuyEnx+mIEmT+XPXRC2iw1lYaY+KAbwO/DnoykQh3oKaeh5cWsvjzXQzu1Zm/ZY5n3NoV8JNMqK31Pam8HDIzfe9nZDgXVoLKn0Mf/YCXm49TRwGvW2uXBjeWSOSy1vJ63nbmLd9M7fFG7rpqKD+54gI6dYiGK7P+r6S/VFsLWVkq6jDmz1Ufm4BvtEEWkYj3xf4aZud6+HTbIca6uzNvRipDep8wn6OiouVPPNV2CQu6Kl4kBNQ3NvHc6lKefb+E2JgoHpmRxg3pSV+fz+Fy+Q53nExXWoU1FbWIw9aVHWJWroeSfTVMGdmPOdem0LvrKeZzZGf7jkmfePgjPt63XcKWilrEIUdqG3h0xWZe/ayCAYlx/OWHF3NFcu/Tf9KXx6GzsnyHO1wuX0nr+HRYU1GLtDFrLUs37ebBtwo5dLSe2y89n59NHEZ8Rz+/HTMyVMwRRkUt0oZ2HK7lV4vyeX/LftIGdOOlH15M6oBuTseSEKeiFmkDjU1eXvp7Gb/972KMgfuuuZBbL3FrPof4RUUtEmT5O49wb+4m8ndWceXw3jx03QgGnqf5HOI//TgX8VdOju+W7agov27dPlrfyNylhUz9/cfsrarn2Zlj+PO/pKukpdW0ohbxR07OVy+LO8Ot2+9t3suvFhWws/IYM8e5+OWk4XSLi2nDwBJOjG/cdGClp6fbvLy8gO9XxDFud8s3mgwaBGVl/3i4r7qOB98qZNmm3Qzt3YV5M9K42N29zWJK+2WMWW+tTW/pY1pRi/jjDLdue72W19Zt55G3i6hv8HL3xGH8+PIL6NhBRxfl3KmoRfxxmlu3t+6tZlauh7zyw4wf3J1509MY3KtL22eUsKUf9yL+yM723ap9grqu3XjiP37L5Gc+omR/DY9dP5JXbx+vkpaA04paxB8n3bq9Nv0qsib9P0r3RTNtdD/um5JCzy6dnM0oYUtFLeKvjAwqp3+fecuLeD1vB0ld4lgwLY3LhvVyOpmEORW1iB+stSzZuIuH3iqk8lgDd1x+AXddNZS4jtFOR5MIoKIWOYOKg7VkLfLw0dYDjBrYjVf+dRwp/ROcjiURREUtcgoNTV7+/PE2nnqnmGhjeODaFG7+ppvok4f5iwSZilqkBRu3V3Jvroei3VVMTOnDg1NH0D8xzulYEqFU1CInqKlv5PGVW3h5bRm9u3biuZsuYlJqX6djSYRTUYs0W1W4lzmL89lTVcfN4wfxi+8kkxCr+RziPN3wIqGvlVPrWmvPkTrueGU9ty/IIyE2hjf/7RIeui5VJS0hQytqCW2tnFrXGl6vJefTcn69YgsNTV7umZTM7ZcOJkbD/CXEaHqehDY/p9a11uY9VczK9fA/FZV8a0hPsqenMqhH57Pen8i50vQ8ab/OMLWuteoamnjm3a3M/7CUrrEdeOIHo5j+jQEYo0vuJHSpqCW0nWZqXWutKTnA7IUeyg/W8r0xA8m65kK6d+4YgJAiwaWDcRLaWphaR3y8b7ufDh09zs9f/5yMFz7FADk/GsdvfzBKJS3thlbUEtpOmlqHy+UraT9OJFpryd2wk7nLCqmua+TOK4Zw55VDiI3RfA5pX1TUEvoyMlp9hUfZgaNkLfKwpuQgY1yJPDJjJMl9uwYpoEhwqaglrBxv9PL8R6U88+5WOkZHMXdaKjPHuojSfA5px1TUEjbWlx9mdq6HLXur+W5qXx6YOoI+CbFOxxI5Zypqafeq6hp4bMVmcj6toG9CLM/fks7ElD5OxxIJGBW1tFvWWlYW7GHO4gL219Rz6yVu7r46mS6d9M9awssZ/0UbY5KABUBfwAvMt9Y+HexgIqezq/IYcxYX8E7RXi7sl8Dzt6QzKinR6VgiQeHP0qMRuNtau8EY0xVYb4xZZa0tDHI2ka9p8loWrC3j8ZVbaLKW2ZOHc9uE8+mg+RwSxs5Y1Nba3cDu5verjTFFwABARS1tqnBXFbNyN7FxxxEuH9aLudNSSeoef+ZPFGnnWnUwzxjjBr4BfNrCxzKBTADXWdzeK3Iqx4438dQ7xbzw8TbOi4/h6RtHM3VUf83nkIjhd1EbY7oAbwI/tdZWnfxxa+18YD74pucFLKFEtA+K93PfIg/bDx3jhvQkZk0eTmK8bv2WyOJXURtjYvCVdI61Nje4kUTgQE09Dy8tZPHnuxjcqzOvZY5n/OAeTscScYQ/V30Y4M9AkbX2ieBHkkhmreWNvB1kLy+i9ngjd101lJ9ccQGdOmg+h0Quf1bUE4CbAY8x5vPmbbOttcuDlkoi0hf7a8ha6OGT0kOMdXdn3oxUhvTWfA4Rf676+BjQWRsJmvrGJp5bXcqz75cQGxPFIzPSuCE9SfM5RJrpFi5x1LqyQ8zK9VCyr4YpI/sx59oUenfVfA6RE6moxRFHaht4dMVmXv2sggGJcfzl1ou5Ynhvp2OJhCQVtbQpay3LPLt5YEkhh47Wc/ul5/OzicOI76h/iiKnou8OaTM7DtcyZ3EB723eR9qAbrz0w4tJHdDN6VgiIU9FLUHX2OTlpb+X8dv/LsYYuO+aC7n1Erfmc4j4SUUtQZW/8wj35m4if2cVVw7vzUPXjWDgeZrPIdIaKmoJiqP1jTy5qpgX12yjR5dOPDtzDJPT+mo+h8hZUFFLwL23eS+/WlTAzspjzBzn4peThtMtLsbpWCLtlg4Syqnl5IDbDVFRvrc5Oad9+r7qOv79vzZw20t5xHWM5o07vsm86WkqaZFzpBW1tCwnBzIzobbW97i83PcYICPjK0/1ei2vrdvOo28XUdfg5ecTh/HjywdrPodIgBhrAz+RND093ebl5QV8v9KG3G5fOZ9s0CAoK/vHw5J91czK9bCu7DDjB3dn3vQ0Bvfq0mYxRcKFMWa9tTa9pY9pRS0tq6g47fa6hib+sPoL/ri6hM6dOvDY9SP5/kUDdbJQJAhU1NIyl6vlFbXLxdovDpK10EPpgaNMG92f+6ak0LNLp7bPKBIhdDJRWpadDfFfvd658rxe3PNvT/DPz39Cg9fLy7eN5akbv6GSFgkyrailZV+eMMzKwlZUsOTSGTz0rX+h8kg0P778fH561TDiOupkoUhbUFHLqWVksP2708lalM+HxfsZ1T+RV6ankdI/welkIhFFRS0tamjy8uLH23jynWKijeGBa1O4+ZtuojXMX6TNqajlazZur+TeXA9Fu6uYmNKHB6eOoH9inNOxRCKWilr+oaa+kcdXbuHltWX07tqJ5266iEmpfZ2OJRLxVNQCwKrCvcxZnM+eqjpuHj+IX3wnmYRY3fotEgpU1BFub1UdDywp4O38PST36crvZ47hokHnOR1LRE6goo5QXq8l57MKHnt7M8ebvPznd5LJvGwwMRrmLxJyVNQRaMueamblbmJDRSUThvQge1oa7p6dnY4lIqegoo4gdQ1N/O69rfzpg1K6xnbgiR+MYvo3Bmg+h0iIU1FHiDUlB8ha6KHsYC3fGzOQrGsupHvnjk7HEhE/qKjD3KGjx5m7rJDcDTtx94gn50fjmDCkp9OxRKQVVNRhylrLwv/ZycNLC6mua+TOK4Zw55VDiI3RfA6R9kZFHYbKDhwla5GHNSUHGeNK5JEZI0nu29XpWCJyllTUYaShycv8D0t55t2tdIyO4uFpqWSMdRGl+Rwi7ZqKOkysLz/M7FwPW/ZW893UvjwwdQR9EmKdjiUiAaCibueq6hr4zYot/PXTcvomxPL8LelMTOnjdCwRCSAVdTu2In8P9y/JZ191Pbde4ubuq5Pp0kl/pSLhRt/V7dDuI8eYs7iAVYV7ubBfAvNvTmdUUqLTsUQkSM5Y1MaYF4EpwD5rbWrwI8mpNHktr6wt4zcrt9BkLbMnD+e2CefTQfM5RMKaPyvql4DfAwuCG0VOp3BXFbMWeti4vZLLhvUie1oqSd3jz/yJItLunbGorbUfGmPcbZBFWnDseBNPvVvMCx9tIzEuhqdvHM3UUf01n0MkggTsGLUxJhPIBHC5XIHabUT7sHg/WYs8bD90jBvSk5g1eTiJ8ZrPIRJpAlbU1tr5wHyA9PR0G6j9RqIDNfXMXVrIos93MbhXZ17LHM/4wT2cjiUiDtFVHyHEWssb63cwb3kRR+sb+Y+rhvKTf7pA8zlEIpyKOkSU7q9h9kIPn5QeYqy7O/NmpDKkt+ZziAic8bouY8yrwFog2Rizwxjzr8GPFTmON3p55t2tTHr6Iwp3VfHIjDRe61LKkLFpEBUFbjfk5DgdU0Qc5M9VH//cFkEi0bqyQ8zK9VCyr4YpI/sx59oUei95E36cCbW1vieVl0Nmpu/9jAznwoqIY4y1gT/vl56ebvPy8gK+33Bx5FgDj769mVc/q2BAYhxzp6VyxfDevg+63b5yPtmgQVBW1pYxRaQNGWPWW2vTW/qYjlG3IWstyz17eOCtAg7W1HP7pefzs4nDiO94wl9DRUXLn3yq7SIS9lTUbWTH4VrmLC7gvc37SBvQjb/cejGpA7p9/YkuV8sral2bLhKxVNRB1tjk5aW/l/HEqmIA7rvmQm69xH3q+RzZ2b5j0l8eowaIj/dtF5GIpKIOovydR7g3dxP5O6u4cnhvHrpuBAPPO8N8ji9PGGZl+Q53uFy+ktaJRJGIpaIOgqP1jTy5qpgX12yjR5dOPDtzDJPT+vo/nyMjQ8UsIv+gog6w9zfv475F+eysPMbMcS5+OWk43eJinI4lIu2YijpA9lXX8eBbhSzbtJuhvbvwxh3f5GJ3d6djiUgYUFGfI6/X8tq67Tz6dhF1jV7unjiMH19+AR07aJi/iASGivoclOyrZlauh3Vlhxk/uDvzpqcxuFcXp2OJSJhRUZ+FuoYm/rD6C/64uoTOnTrwm+tHcv1FAzXMX0SCQkXdSp+UHmT2Qg+l+48ybXR/7puSQs8unZyOJSJhTEXtp8ra48xbXsTreTtI6h7HgtvGctmwXk7HEpEIoKI+A2stSzbu4uGlhRyubeCOyy/grquGEtdRw/xFpG2oqE9j+6Fashbl82HxfkYlJbLgtjRS+ic4HUtEIoyKugUNTV5e/HgbT75TTLQxPDh1BDeNH0R0lE4WikjbU1GfZOP2Su7N9VC0u4qrU/rw4HUj6NctzulYIhLBVNTNauobeXzlFl5eW0bvrp147qaLmJTa1+lYIiIqaoBVhXuZszifPVV13Dx+EP/5nWS6xmo+h4iEhogu6r1VdTywpIC38/eQ3Kcrz2aMYYzrPKdjiYh8RUQWtddryfmsgsfe3szxJi/3TErm9ksHE3OqYf4iIg6KuKLesqeaWbmb2FBRybeG9CR7eiqDenR2OpaIyClFTFHXNTTxu/e28qcPSkmIi+HJG0YxbfQAzecQkZAXEUW9puQAWQs9lB2s5fqLBjJ78oV079zR6VgiIn4J66I+dPQ4c5cVkrthJ+4e8fzXj8ZxyZCeTscSEWmVsDx7Zq3lzfU7uOq3q1ny+S7uvGIIK356mX8lnZMDbjdERfne5uQEO66IyGmF3Yq67MBRshZ5WFNykDGuRB6ZMZLkvl39++ScHMjMhNpa3+Pyct9j0C+bFRHHGGttwHeanp5u8/LyAr7f02lo8jL/w1KeeXcrHaOj+OV3hzNzrIuo1szncLt95XyyQYOgrCxQUUVEvsYYs95am97Sx8JiRb2+/DCzcz1s2VvN5LS+3H/tCPokxLZ+RxUVrdsuItIG2nVRV9U18JsVW/jrp+X0TYjl+VvSmZjS5+x36HK1vKJ2uc5+nyIi56hdFrW1lpUFe7h/SQH7quu59RI3d1+dTJdO5/hysrO/eowaID7et11ExCHtrqh3VR5jzuIC3inay4X9Eph/czqjkhIDs/MvTxhmZfkOd7hcvpLWiUQRcVC7Keomr2XB2jIeX7mFJmuZPXk4t004nw6Bns+RkaFiFpGQ4ldRG2MmAU8D0cAL1tpHg5rqJIW7qpiVu4mNO45w+bBezJ2WSlL3+LaMICLimDMWtTEmGngWmAjsANYZY5ZYawuDHe7Y8SaeereYFz7aRmJcDE/fOJqpo/prPoeIRBR/VtRjgRJrbSmAMeY14DogqEX9QfF+7lvkYfuhY9yQnsSsycNJjNd8DhGJPP4U9QBg+wmPdwDjTn6SMSYTyARwncPlbAdq6nl4aSGLP9/F4F6deS1zPOMH9zjr/YmItHf+FHVLxxm+djujtXY+MB98dya2Noi1ljfydpC9vIja443cddVQfnLFBXTqEN3aXYmIhBV/inoHkHTC44HArkAHqXrlVX69IYrk/eXM8+QyJOU/oMOwQH8ZEZF2x5+iXgcMNcacD+wEbgRmBjRFTg7d/i2TNzt2w1W5hyishiGJiDQ740XI1tpG4E5gJVAEvG6tLQhoiqwsqK3FXbnbV9LguzswKyugX0ZEpD3y6zpqa+1yYHnQUmgYkojIKYXGLw441VUiGoYkIhIiRZ2d7Rt+dCINQxIRAUKlqDMyYP5834B+Y3xv58/XiUQREUJpKJOGIYmItCg0VtQiInJKKmoRkRCnohYRCXEqahGREKeiFhEJccbaVg+6O/NOjdkPtPDrvP3SEzgQwDjtgV5z+Iu01wt6za01yFrbq6UPBKWoz4UxJs9am+50jrak1xz+Iu31gl5zIOnQh4hIiFNRi4iEuFAs6vlOB3CAXnP4i7TXC3rNARNyx6hFROSrQnFFLSIiJ1BRi4iEuJApamPMJGPMFmNMiTHmXqfztAVjzIvGmH3GmHyns7QFY0ySMeZ9Y0yRMabAGHOX05mCzRgTa4z5zBizsfk1P+h0prZijIk2xvyPMWap01nagjGmzBjjMcZ8bozJC+i+Q+EYtTEmGigGJuL7refrgH+21hY6GizIjDGXATXAAmttqtN5gs0Y0w/oZ63dYIzpCqwHpoXz37MxxgCdrbU1xpgY4GPgLmvtJw5HCzpjzM+BdCDBWjvF6TzBZowpA9KttQG/ySdUVtRjgRJrbam19jjwGnCdw5mCzlr7IXDI6RxtxVq721q7ofn9any/LHmAs6mCy/rUND+Maf7j/OooyIwxA4FrgBeczhIOQqWoBwDbT3i8gzD/Bo50xhg38A3gU4ejBF3zIYDPgX3AKmtt2L9m4CngHsDrcI62ZIH/NsasN8ZkBnLHoVLUpoVtYb/qiFTGmC7Am8BPrbVVTucJNmttk7V2NDAQGGuMCevDXMaYKcA+a+16p7O0sQnW2jHAd4F/bz60GRChUtQ7gKQTHg8EdjmURYKo+Tjtm0COtTbX6TxtyVpbCawGJjmbJOgmAFObj9m+BlxpjPmrs5GCz1q7q/ntPmAhvkO6AREqRb0OGGqMOd8Y0xG4EVjicCYJsOYTa38Giqy1Tzidpy0YY3oZYxKb348Dvg1sdjRUkFlrZ1lrB1pr3fi+l9+z1t7kcKygMsZ0bj5BjjGmM3A1ELCruUKiqK21jcCdwEp8J5het9YWOJsq+IwxrwJrgWRjzA5jzL86nSnIJgA341thfd78Z7LToYKsH/C+MWYTvgXJKmttRFyuFmH6AB8bYzYCnwHLrLUrArXzkLg8T0RETi0kVtQiInJqKmoRkRCnohYRCXEqahGREKeiFhEJcSpqEZEQp6IWEQlx/wttQ6eUfiGSqAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "x = np.linspace(0,5,10)\n", "\n", "noise = np.random.normal(0,.3,len(x))\n", "beta_true = np.random.normal(0,1,2)\n", "\n", "t = beta_true[0] + beta_true[1]*x\n", "\n", "tnoisy = t+noise\n", "\n", "plt.scatter(x, tnoisy, c='r')\n", "plt.plot(x, t)\n", "plt.show()" ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }