Skip to content

Instantly share code, notes, and snippets.

@andlima
Last active April 19, 2019 02:56
Show Gist options
  • Save andlima/804218eeef71b993946b4c135638a4cb to your computer and use it in GitHub Desktop.
Save andlima/804218eeef71b993946b4c135638a4cb to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Populating the interactive namespace from numpy and matplotlib\n"
]
}
],
"source": [
"%pylab inline\n",
"\n",
"import seaborn as sns\n",
"import pandas as pd\n",
"import scipy.stats"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Expectation–maximization algorithm\n",
"\n",
"## Definition\n",
"\n",
"In statistics, an expectation–maximization (EM) algorithm is an iterative method to find maximum likelihood or maximum a posteriori (MAP) estimates of parameters in statistical models, where the model depends on unobserved latent variables. The EM iteration alternates between performing an expectation (E) step, which creates a function for the expectation of the log-likelihood evaluated using the current estimate for the parameters, and a maximization (M) step, which computes parameters maximizing the expected log-likelihood found on the E step. These parameter-estimates are then used to determine the distribution of the latent variables in the next E step.\n",
"\n",
"Source: https://en.wikipedia.org/wiki/Expectation%E2%80%93maximization_algorithm\n",
"\n",
"## This notebook\n",
"\n",
"This notebook implements a version of EM for a dataset with **a single categorical latent variable**, which can be used for clustering. It's important to notice that EM can be used for more than one latent variable and also for numerical latent variables."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Randomized dataset"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"def randomized_dataset(n_observations, n_features, n_groups, seed=0):\n",
" '''Generate a random dataset.\n",
"\n",
" - Each observation belongs to a group (a cluster)\n",
" - The group is identified as a categorical variable\n",
" - Each group is a multivariate gaussian distribution with dimension ``n_features``\n",
" - Each observation has features that come from the respective group's distribution\n",
"\n",
" Args:\n",
" n_observations: number of observations in the dataset\n",
" n_groups: number of groups (latent categorical variables)\n",
" n_features: number of features\n",
"\n",
" Returns:\n",
" X (array[n_observations, n_features]): random dataset\n",
" labels (array[n_observations]): identifiers for the groups\n",
" probs (array[n_groups]): probability a random observation belongs to a certain group\n",
" mu (array[n_groups, n_features]): for each group, the average value of each feature\n",
" sd (array[n_groups, n_features]): for each group, the standard deviation of each feature\n",
" '''\n",
"\n",
" random.seed(seed)\n",
"\n",
" probs = exp(randn(n_groups))\n",
" probs = probs / sum(probs)\n",
" assert np.allclose(sum(probs), 1.0)\n",
"\n",
" mu = random.randint(128, 256, size=(n_groups, n_features))\n",
" sd = random.randint(2, 8, size=(n_groups, n_features))\n",
"\n",
" labels = random.choice(range(n_groups), size=n_observations, p=probs)\n",
"\n",
" X = mu[labels] + sd[labels] * random.randn(n_observations, n_features)\n",
" scale_factor = np.max(X, axis=0)\n",
"\n",
" return (X, labels, probs, mu, sd)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"# Create a dataset with 2 features, so we can view the clusters\n",
"\n",
"X, labels, probs, mu, sd = randomized_dataset(n_observations=1024, n_features=2, n_groups=4, seed=37)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Plot observations"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEFCAYAAADzHRw3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzsnXmYVNWZ/z91a+2luruquwFZFAH7sCMIoghKXDHGPUZHYzLJjKPzJI4+mRkTjf6SmZiZ0UkmEzOTPcYxoyGD0SxqEhdEQTSANEsjHFYBWaSX6uq99t8ftVDdfWuv6u6qOp88eeDu53Tj9773Pe9iCIVCKBQKhaL40UZ7AAqFQqHID0rQFQqFokRQgq5QKBQlghJ0hUKhKBGUoCsUCkWJoARdoVAoSgQl6IoRRwhxjxBiuxDifSHELiHEL4QQZyY5/2UhxOwU9/xnIcRnchjTOiHEJwt1/pBra4UQa7O5VqFIhmm0B6AoL4QQ3wIWAJ+QUh4VQmjAp4F3hBBLpZQfDr1GSvnxVPeVUv6//I+2YDiA80d7EIrSQwm6YsQQQkwG7gGmSCldAFLKIPC0EOI84EHgC0KID4A/A/OBh4DvAJ+UUm4RQnwF+CugG3gLuEFKOVUI8RTQIqX8lhBiAPg34ErgDOBxKeUPhBBVwA+Ac4D6yD1ul1LKJGOeAPwQmAkEgR9KKZ+IOz418tzqoduRa58GGiKnvySlfAT4OVAhhNgGnAc0Ad+NjMkIPCGlfFIIsTKyvxeoBlYAP42MPwi8B9wd+RkqFMrlohhRlgK7o2I+hNeA5XHbLVLKWVLKF6I7hBBXAX8JLCEshPYEz7ECbVLKZcAnge8IIWzA1UCnlPJCKWUTsBn4Yooxfx/YK6WcCVwI/I0QYkaKa6LcBRyUUi4iLMbnCCFqgc8B/VLKcwED8BzwFSnlecAlwD8IIS6I3GMu8BdSyvnAdYA9ct2SyPFpaY5FUQYoC10x0pgT7LcC8XUo1uuc83FgjZSyE0AI8d/AZQnu99vIn1sj966SUj4nhDgohLgXmAGsBN5JMd7LgQcApJRuwgKLECLFZQD8EXg5sj7wGmHRdgshHHHnNAHTgSfj7lkBLAR2A0ellIcj+zcA/yKEWAe8CvynlHJ/OgNRlAfKQleMJO8StlIn6Bz7GLAxbrtH5xw/YYs2SiDJs/oBpJTRl4RBCPG3wM+APuBZ4JdD7qeHn7gXjRBimhCiJu54aMg9LNG/SCk3A2cDPwamApsirqV4jIBbSnlu9P/ABYTdMhD3c5BSHiL8IvpXoAZ4TQhxbYrxK8oIJeiKEUNKeQx4AvilEGJSdL8Q4nPAzcBjKW7xEnBzxG0BYV96JtXlrgKeklL+DJDAtYQFNRmvEXaREHnu64R92FE6AUtcFM5fRA8IIf4NeERK+RvgPmAXYQvfDxiFEIbIOPqFEJ+OXDMFaCHsUhpE5IX0c+AVKeWXgT8Bi9KevaLkUYKuGFGklA8C/wv8VgjRIoTYR9itcWGcayHRtWuBnxCOiNkC1BK2ttPlW8DdQogdhF06WwlbvMn4IjArcs3bwL9KKd+LG5ObsEvmD0KIzUS+DCL8J3CuEKIF2AIcAlYDJ4BNhAXeDlwP/HXkGa8Qfgm8rTOWpwm/gN4XQrwXmf8TOucpyhSDKp+rKBaEEIuBZdEoEyHEl4ClUspbR3dkCsXYQC2KKoqJvcCXhRB/Q9jVcgT4m9EdkkIxdlAWukKhUJQIyoeuUCgUJYISdIVCoSgRRtWH3traPSL+HoejEpcrk2CI4kPNsTQohzlCecyzkHNsbLTr5k+UhYVuMqUKNS5+1BxLg3KYI5THPEdjjmUh6AqFQlEOKEFXKBSKEkEJukKhUJQIStAVCoWiRFCCrlAoFCWCEvQ08Aa8tPa14w14R3soCoVCkRBVyyUJgWCA5/e/xI7WXbg8nTisdcxvnMNNM67BqJV+2JVCoSgulKAn4fn9L7Huww2x7Q6PK7Z9S9N1ozUshUKh0EW5XBLgDXjZ0bpL99jOtl3K/aJQKMYcStAT4PZ04/J06h7rGOjE7eke4REpFApFcpSgJ6DWasdhrdM95rTVUWtN1HBeoVAoRgcl6AmwGC3Mb5yje2xewxwsRovuMYVCoRgt1KJoEm6acQ0Q9pl3DHTitNUxr2FObL9CoVCMJZSgJ8GoGbml6Tqun74Kt6ebWqtdWeYKhWLMogQ9DSxGC42V9aM9DIVCoUiK8qErFApFiaAEXaFQKEoEJegKhUJRIihBVygUihIh6aKoEMIMPAlMBazAo8AR4HtAAPAAn5FSfiSEeAK4CIimUF4vpXQXaNwKhaLM8fkC9PV4qay2YDarYnmQOsrl00C7lPJOIUQ90AwcAu6VUm4TQtwNfBn4ErAIuEpK2VbQESsUirImGAyyce0BDu1to6fLQ3WNlbObGlh26XQ0rbydDqkEfQ3wXNy2H7hNSnki7voBIYQGnAP8WAgxHviZlPLJvI9WoVCUPRvXHmDnlmOx7Z4uT2x7+eXnjNawxgSGUCiU8iQhhB34HfATKeWzkX3LgJ8BFwMDwH3AfwBG4A3g81LKHcnu6/cHQiaT+lRKhcfvxTXgxmGrxWpSiU2K8sXn9fP9x9fhdvUPO1bnqOBvH1iJ2VIW6TUGvZ0pZy6EmAK8AHw/TsxvBb4KXCOlbBVCGIHvSin7IsfXAguApILucvVlNINsaWy009pafNURM2mwUaxzzAQ1x9Ih23m6Xf26Yg7g7uzn8Acd1Doqch1eXijk77KxUb84YKpF0fHAK8AXpZSvR/Z9GrgbWCml7Iic2gSsFkIsIhw5sxz4n/wMvXxRDTYUisFUVluorrHS0+UZdqzabqWyury/YFOtIDwEOIBHhBDrhBDrCUe42IHnI/v+SUq5G3gGeBd4E3haSqnfHUKRFqrBhkIxHLPZyNlNDbrHpjY1lH20S1ILXUp5H2HfeEqklI8Dj+djUIr0Gmyo+jKKcmTZpdMB+GBvGz3dHqrtVqZGolzKnbJYPShGog02OjyuYcdUgw1FOaNpGssvP4ell0xTcehDKO+gzTGMarChUCTHbDZS66hQYh6HstDHMKrBhkKhyAQl6GMY1WBDoVBkghL0IkA12FAoFOmgfOh5wBvw0trXnlUoYS7XKhQKRTzKQs+BTDI583mtQqFQ6KEEPQcSZXIGQwFuFTdmdS2oLFCFQpEdyuWSJckyOTcc+zOr97xAIBjI+FqVBaooVny+AG5XPz6f/r97ReFRFroO3oA3ZVRJskzOIEHWH38nFqWSybUqC1RRbKj65GMHJehxZOLXTpbJGWVn2y6un74Ki9Ey6CWhskAVpUSp1icvxo5IStDjyMSvHc3kjD9/KB0DnXQMuFl/7J1hL4m5DbN569jbw64phizQoMeD3+3GVFuLZrUW/DrF2MXnC3Bor36Tsg/2trH0kmlFI4ZRivmLQwl6hFR+7ailHc9NM64hGAqw4difCRIcdp3TVsfao2/y9vFNsX3Rl8Qlk5axcvLyQVmgc+pnsmLShXgD3jEp6qFAgNY1q+lp3oq/owOT00n1wkU03nIbBmPi/2izvU4x9unr8eqWsgXo6fbQ1+MdM/XJ06WYvziUoEfIxq9t1IzcKm4kFIL1x98Zdp3FaOWd41t077mzbTePXPD3XD99FR0Dnaw7+jYtbXtYf+zdhK6e0bZwW9espvO1V2Pb/vb22Pa42+7I+3WKsc9Yrk+ejcuk2L84ylrQ8+XXvqXpOoyakZ1tu2gfcFFnqaHCXMmJ3pMJr+nwuHhm96+5+uzLWHf07UEvhKGunkJbuOksAgc9Hnqat+oe62lupuHGT+q+ZJJd596wnvrrb8RYUZn94BWjSrQ+ebxFG2Uk6pPriXYuLpNi/+IoS0FPtPg5yyl4+8S7w84f6tceKoBGzRhzv2xv24Xb00WXtyflOLacambLqWa0BNGjUVdP55o1BbFwM1kE9rvd+Ds6dO/jd3XEvhyGfkEkuy40MMCpXz7DGZ+/K+s5ZMJof+GUKiNVnzxevI1GQ0LRzsVlMpa/ONKhLAU90eKnVQv/sjQMBAnhtNYxv3FurLphMgF8fv9LvHXstJUd0vGpJ0LP/w5hV09nd3tWlnE6ZLIIbKqtxeR04m9vH3YfY50D1yt/pHfn9tgXRNW8BdRddgXG6mpMDif+juHXAfTt2UPQ4ymowCoffmEpdH1yPYvbajPRfqo3dk5UtAf6vBw/6ta9z6E0XCaj/cWRK2Un6MkWPz3BcEJPkBAAs+tnxoTNG/CyWr7An0++Fzs/KoCBYIBd7XvyPlanrY6Kfh+u9nbdFt9Ry9gyblzG9850EVizWqleuGjQl0IUY1UV7nVrT4+rvR33urW4163FVF+PIYlYBzpdWc8hXZQPf2SI1ifPN3oWdyK3yL73WxPep6crPZdJMXdEKjtBT7b4OZTNHzVzw/SrefHQq2xvbUl43c62Xbi9+e/uPa9hDn/q2MJZVRq1vcOteJPDiam2NuP7Bj0eOk5+QHevS/dfQPwicLybovGW24Dwl4Hf1YHJ4aRq/nx6d2xP+KyYRW80QmB4BqHJ4USrqMB76tQgV0i+3CPZ+v4VY4Nki5SZYjCAxZrawi7mjkhJBV0IYQaeBKYCVuBR4AjhRtEBwAN8Rkr5kRDiLuBuwA88KqV8sYDjzpp0EoKieAIefrX3N2z+qDnpeZ3eLmqtNbg9XTmNTUMjSJB6m4N5DXP4xNlX8M0/f4fAJCsL9/YPO79iwfyMxCjQ18ep1c/Sv2c3flcHd1YZ2T/RzPpF1YS0098ATlsdNaZKTq1+RtdN0XDjJwf5zN3r3kj5bIPZTEhH0LXKSo48+vXTrpoF52LAQM/25mHPzYZ0fP+F/Dood3JNzkm2SJkpoRB4PQHSXYMv1BdHIUlloX8aaJdS3imEqAeagUPAvVLKbUKIu4EvCyEeB/4OWAzYgA1CiFellPn5TeSRdBKC4pGu/SnPqbc5mOk4h7dPbEp6ntVopdJow+XV9/Etn7SUS6dcTK3VjtFg5Jk9z+HydrJ+UTUA0455qO4N0lOlcXCSjY9duyqtOUR9yO4N6wkNDMT223v8LNzrB+CtxaejeOY1zMH9/PNJ3RRREUzmWx80Bo8X+7KL6JcyZt1rlZV4jx4Z9Az32tcHXTfouffek9Z840k2vmy/cBSpyVdyTrJFykyptlvSWtQsxgzRKKkEfQ3wXNy2H7hNSnki7voB4Hzg7YiAe4QQ+4H5wOZkN3c4KjGZRuYH1th4WrDurr+Nyu1mtny4g9a+9ojHfDhWo4WuNFwpS888lyunX5xS0D0BDxZj+Edu1SxomgGP30tDpZPFk+dz54KbMWpGPH4vP9r8vzF/fUgz8NZiOxvPraaqP0BvhRFnTQN3nTkFq+n0P9D4OcZz8CdP6vq+o8w45uXdc6Gupp7Fk+dz69lXsPOnD+ie27d9G86/+RzG2JeBnd4LL+DEiy8lnbu1sZ45938BAG+HC2NVJTv+Xv8ZevTv2E7A40k4x8QkHl/jhUsZP7khw/sVnsznOPb4429adCNNKiosrLphLpD+PGcvmMim9YeG7R8/sYaBPi/uzgGdq3Tuc+4kJk6sS3g8GAjyyu/fR7acxN3ZT21dBWLuBK68djaaMbsM0ZH+XSYVdCllD4AQwk5Y2B+OirkQYhnwReBi4Cog3uzsBlKaPi5XX3ajzpDGRjutrYOF+ZrJV3PFGZfh9nTz+pE3WX98eLjiknGL2Nn+Pm5vYlfK0gnnsWrilQT6Alg1S2xhNRHd3vDKvCfohSAsGbeQ22fdjMVoobW1m+f2vcifT2zWvY/fZMBtD//KZjtn0eXyEPZ66c8Rwj7k1neGzy0ee2+Av9szjvGf/Tw9f/gD2//zHwh06q8XeFtb2fEv/07DtTdgmTAhvFh67U1Ytu8cZG0Pu66rmz0//jmOVdfgPXECzWbF05q+b9TT1sbAyY/o7PFl7FevvvYm6ga8g3z/1QsXUn3tTbo/s9Ek0e+xmPD5Ary//bjusd3bjzP//MlMnFiX9jwXXjiF/n6v7iJlIBDirT/tY2/L8JwPs0XD7wvGzl944ZSkz9zw2r5BLyG3q59N6w/R3+/NKkO0kL/LRC+KlIuiQogpwAvA96WUz0b23Qp8FbhGStkqhOgC4p9gB9JbeRxFoq3dbmm6HqNmYkdrCx2ezkHhioa9Bt0sUACn1cFt4kaMmpFAIEAooa2fmC2ntrFkwrmcVTOFP3ywVre+y1CWTjgv7UbRyXzI8fh2tPDhA/+ou3A5lL5tzRzZ1ozBZqN22XLqb7iJYF9v0mtCAwN0vvbq6S8FgwE0La3nARgsFt7/xjfxtrWnFXY4dFF13G13DPL9q4XQwpFOck4maJrG0kumMWvBGRCCGkdFzBWiafCxjzdhtRmHCf6SFWcz0OdLy3VS7BmiUVItio4HXgG+KKV8PbLv04QXP1dKKaNKsQn4phDCRnjxdBbQUrBR55lkzZhvabqOg12HOdYz3OKY3xhOOAoEA6yWL+AN+jJ+dogQ39/xcwyR/6Ui/iWSDun6uIG0xTVKaGCAzrWvERgYSOulMfjiUEbPCw0M4I34/5OFHSaLOVeMDPlMzknHF58sKsVqTS+Qr9gzRKOkmu1DgAN4RAjxCGAE5gKHgeeFEABvSim/JoR4AlhPuGnGV6WU6Tm2xhB6zZiNmpEvL76X5/b9jh1tu3B7unHa6sJRKNOuoLWvndePvDUoPj0bQpH/JcPkD7G46ixM/lD4N6FDNIvVbrBAuxsDUDV/Ae431upfkAf6du9Cq6sj6EodOaSHyeHE3+nCWFuHrUngO/Yh3mMfpryup3nroLDDoMfDR888TffG0185UfHvk5JgX69KLBoB8pmck0nWZy5RKcWeIRrFEApl7ibIF62t3SPy8Hz5sqJiWW2u5MVDr8YyRrNxtWSCIRhixdZeZp8MYe0e0BUkZ30lP3pnNTs/amHOxqPMPjSAxR/CABisNowOB/6TJ5I/KAdMjY34WxMndSTD2NBAoKsLvJl3arIvu4jxn/4sbc+voWfr1oQZqXrUXX7FmEssKgUfOpy2rPX83pqmpTVPny/A6p9s0hVZe42VW+86P+nLIdNolaE+9CjzFk8aiz503c/5skssyoWoBb9m7+/SDnvMByu29gyKQ9dzOfxi+69Z9+EGLt7SzcJ9g2PWQ56BsJhbrOAtTCSpv7U1YfJQKgJt2SeOdG98G8/Ro0kXZBOhEosKRz6Sc7J1g2QTMunzBZizcCLBQJAjBzpyzhD1+QJ0tPXi8wVG1PeuBJ3E1Qb19idLmS8EVUEzsxMUbYwKkt9kYOsH23B2+pl2dHgC0mkK/CWRIHmo0KTjntFDJRYVnlzdIGaLhs+rkyVt1hK6QTJx0+iJ/1kz6pl33mSqa6wZi/Gg+0VeCiPZHKOsBT1Rsa3rp63itwf/qFuEK5PSAfng3um309/9Ld1jflcHvo4OTr72Mqu27MPeG0y+rOr1YrBaCXkKY6XHJy2NKMH0C6HFoxKLRp6oG6SuNrXI+30BQhn+ajONVtET/11bj6NphqzcLKPdHKMsBT1qeb9+5C3dOuT7Og8OimqJ7g+GQoRCwYL7zKPU2xw0TjiTY04ngQSZjp2vv4rnzfWpg/4BrbaOylmz6Xl3Y/4HWwg0DUIhTE4nFTOa6N8v8bfrRNNoWkJR1xxOjDYrvhPD1w+qFy5U7pY8kI6veqglXOuo4MzpTl3LNXrugT2t+P36v1efN6jrcsnETZPvUMWxEPpYVoI+1CJPFCaoF6IIsOnkFgYCI1fNYE79bH5/9A0s4wPM0lnrS1UYayhBd2fxiDlgv2g5BAL07d5N96Z30Wz64muZNFnfh65phPp68XW6MNhsAIQ8HkzOeqoXLlShjDmSia96qOXqdvUntFyHnquHvcaKxWrE7eof9CLJJFol36GKYyH0sawEfWj970wt7ZESc4tmYemE8zAQYt2Hb2OYb2bAXxGr5eKvqWT8+RdRu+KSgoYjjjaeQwfxfnjaPx7sj7h0rFbw+WIZnw033ULb82uG1aohGIy5l6L77cuWM/6OO5VlngfSdS8ks1z37DjBkhVTsVrNKc+Nx2Iz8dxT7w17kWQSMpnvUMWxEPo4tltY55GRXszMBqtmwaJZ8AV9tLTt5t0htVz+95p6fnGtk+dvmIT95pvY8pufjfKIC0u8mA/C46Fy/gKmPPgw4267A81iCUerVKYuo9cv81+3vhxJ5V7w+U4vjiezXH3eIBte3Z/WuQBV1Vbqx1XRfqo3dl70RbJx7QGCwSChUAiT5bS0mS0ac8+bOCxaJSr+emTTzCLf98uGshH0kV7MzAQDBs5tmIsn6MUb9BIihMvbiWfIF0G0lkurv4tfv/8Clr2HR2nEo0/ftmY+ePgrnFr9DKFAAL/bTSCNxKZoZIsiNzJJ749arok4ftgVewEkO7eyyszl181koF8/I/uDvW1seHUfLe8dxx8XGePzBjEYDLpRJssunc68xZOw11gxGMKunHmLJ2XdzCLf98uUsnG5ZFIHHcKlbis0G50+N3ZzNR6/F28o88SXdHBYazncnX7oncNax4fH9rJYp+lFORGtDwPQcOMn0ypxoCJb8kMm7gWz2cjEMx26BbQAenu8Mf9yMpeJzxfgt88mXjPq7vJwMMNFyXw3s4i/n81iZsDrG9E49LKx0KN10NPFiJFOX9iS6/b1FEzMAc5xTKfTk77VOMNxNidMvXRXlc2vLyk9zeEGJNULF6U8V0W25IdM3QvLr5gxyA0Sz9AXwFAr1xy5Ti8efSj9vfrWe6qiYNF4+XyJr9lsxNlQNeIFvcrGQgdiFQp3tu2iYyAcXz7DcTZWzcL7HZKOgU6ctjoG/B56/YUv7asZNFZMvIBrp13FPtdB3a8Hm9FKpakSlyc8trn1s5CuA/hNBg4m6GRUbvjb2/C5OmJRK91btxIYUgIgWhlSRbbkj0x6b1qtJmbNPyOtxcp4K7fL1c/Lz+3El2aGs8EQrvs2lHwuSo7lBhhlJejRqoqfOPsK1uz7HXs7DrD5ZDMOax0zHeewcPx8nNY6vrFJP5EnEZpBI5hpBgQQDAUxGDQqzBUJuyjNqRd8qukG+v0eaq12Xtj/Eif7PgKIdDIKMfvgABb/6WtS12wsPTpff5Xxd3yGcbfdQdDjpWv9m4OOhwYGQDOoYlx5JFN3xbJLpxMKhZA7T8asbZNFIxQKEQwGh/m4zWYjJrMxo25FiUpT5WNRMtOSAj6vf1hYZaEpK0GP8uKhVwdVR+zwuNh4chMbT26i2lyd1j0cllrmNszmoonn86Md/4PLm92C647WFq6fvorrp61iX+dBTvScIBgXTvneqR00n2ph+cSlXHP2Fbx74vS4w31ADVj9OjcuM3p37CBwYx9tzz9H19vrdc9RtVsKQ7rp/ZqmYTAYBrlO/N4gLe8dx2DQz8zMtAWdOfKC8PtCsW0xb0JeFiXTDdOMCv+RAx24Xf1Zt9/LhrJzwqYKX+zx9aR1n7NqpnDx5Atwe7rpTNAjNB06PJ2sli/w/P4XOdZzfJCYRwkS5K3j7/AfzT/EG9fJyOQPMe3YmGvbOir4XR2cWv0s7nVrE2aNqgiX0SWTUMcoZrORs2bU61yR4BneYEzMo9uJIlwyIZOxR4Xf7Qq7Q+PDKgtN2Vno+Qpf3NbWwra2/PTwSLeW+kd9pwZtV/UHsJd5pEsMs5n+Pe8nPUVFuIwumWZSxkrw7gsLadQ/Xl1jZeo5YZE/vK+dnm4PVXYrAwO+QeGKUQ7uaeW8ZWdRUZm9Dz3dsY92+n/ZCXqm4Ytjmd4KI91VGrV6om61QoGKcI1JPB78KearIlxGl0wzKYe6OKL+8bOmO1lxRRMAF6wML1D6fQH+78ktus/t7fGy5sktTJvZmLXbI92xj3b6f9m5XDINXxzLRCNd9DDUO0d4NGMAQ4J/zppG7ccuVREuo0wmoY7JLN0jBzpiLo6o/7661pYwLBLCop6L2yPdsSdLjBqJ9P+yE3QIhy+umHghWpFNv8JoG7Zv/aJqmpsqcFdpBAB3lcaembVoA4WLmx+zJIg0qr14JePv+IyKcMkzPl8At6tf1/ediAtWTqN+XFXs3WswQP24Ki5YOW3QeZk2mt68/pCuu2UoiXz16ZBOFuhop/+nahJtBp4EphJu/vyolPJ3kWPfAaSU8oeR7SeAi4Boz6XrpZRjcgXKqBm57MyL2XD83dEeSkb0B4bXG4/Wedl4bjVV/QF6K4xcWruQwE9fHYURji1M9Q2qqmIByKYjUJSNr++n/VRvbDsUgvZTvby77uCgSJFM3DPpFvSC3Nwe6YZpRgX+6IEO3J39g+LzCx3DnsqH/mmgXUp5pxCiHmgWQrwDPA00Af8ed+4i4CopZfb9xEaQWqudOmsNrgwyNMciRjRqrbX0+F3UWWtZNGkeH5+2iqPO91KmwZcyxjoHkx94EAIBQn6/ss7zSDZNHILBIBte28/72/T72g5dMMykamKqgl7x5MPtkSpMMyr8dTdXcPiDDiqrLRiNhqxfgpmQStDXAM/FbfuBauDrwNXRnUIIDTgH+LEQYjzwMynlk6ke7nBUYjKNzH9ojY322N8DwQC/2P4H+vzFn2UZDAa4+8NxeLa2Emg/hKWxm4GlARrOX8zJP/xptIc3agQ6XRx7/Jv4OlxYGxtwnn8+Z3/+s0Uh7PH/VscaPq+fIwd0mowQtkjrbq7AbBkuK3/8TQu7tur3GYCw5WyzmHE2VOHz+unu8vDxG+ZRUWFhb8tJ3J391NZV0DR3AldeOxvNeFoE62orqHVUxMIEkzFrwUQmTqxLY6ZhomOx11h155WKGU3h9oZ//E2L7kuwosLCqhvmZnzfRCQdoZSyB0AIYScs7A9LKQ8Bh4QQV8edWgV8D/gPwAi8IYTYIqXckez+Llf26fWJ+oDqMbT79kg3eS4kK97rpm/f6UQab2sbJ158CYOK5sAX6W7kOdXKiRdfon/AG2uqnYigx4Pf7cZUWzsqETGF7BSfKz5fgI+OdSUUTndnP4c/6Bhmvfp8Ad7fnljMIWw59w14efOXcpAVO/FMB9fdcS5+XzDmpmjv6B3qf/A1AAAgAElEQVR2/ZnTnbrWvNmi4fcFY26PhRdOSevnm4tbKUr0d5ls/ru3H2f++ZMzdr8keumnfOUIIaYALwDfl1I+m+C0PuC7Usq+yDVrgQVAUkHPhkR9QG+acQ1GLfUPpRjqoqeDyR+iptvP7IP6fTwL1Te0mEmWKRoKBGhds5qe5q34OzowOZ1UL1xE4y23FYVVX0iGilum9VLScYlYbCY2vXWQlvdOC19Pl4e9LSc5uPcUs+afkTTbM1FdmSUrpjLQ58/YZ53P3qAjGcqYalF0PPAK8EUp5etJTm0CVgshFhGOnFkO/E9eRjiEoV2Hov0+AW5pum7QuVErvsZ/+j/gsVwXPRUmf4jq3gALZB9nn/CmbgqtGEQ0U9QybtywY61rVsdK8QL429tj26ms+lInUTz4UKK+7aELf+mk77ef6qWrU9/y93uDKcU02YJltBtSuuQ7OWgkOxmlstAfAhzAI0KIRyL7rpZSDvrJSyl3CyGeAd4FfMDTUsq8m8HJrOudbbu4fvoqLEbLMCu+vsLB9NqzueWc64oyscgQDLFiaw/TjnmoUSKeNYkyRYMeDz3NW3WvKff6L8nEzWCAEGCPWMMXrJzGhtf26bopEi1wDnpWirDDdMQ03boyyUjHoq6stqQdrZLJAm+upPKh3wfcl+DY14dsPw48nreR6ZDMuu4Y6MTt6aaxsn6YFd/W30FbfwfbWlu48IwlzG2YzVvH3i7kUPPKiq09qkxuHkiUKep3u/F36C/0JbPqy4Fk4hYKwbW3LWD8pBrMZiMbXtuX0E0RdYkc2NOatC55Mkaq0XIyi7qq2sL2zUc5vL89I996JqWGc6GoUv+TWddOWx21VntSK94T8LDuww1cMmkZKycvZ+upHXR5uwo97JxQBbgywGDA5KynasECDBjo2b4Nv6sj1kw6UTy6qbY2Ybejcq//kkzc7DXWmJinclMsWTE1spW4MbvZoiW10keq0XIyi9paYR4UrZOubz3fnZESUVSCHk3b14tQmdcwB4vRQmtfe0ofeUv7bh5e+ve4Pd00tyZuaTUWKLYCXJO/8lVO/s/P8Z9IHtWQb0xOJxP/7ktYGhtjVnjDzbekFbGiWa1UL1w0yIcepdzrv6TrLkjlptjw6n72tnyU9Fli3gQMBgN7dpzQFfZM3RO5JPHoWdRnTndyOEHIZrq+9Xy4hJJRVIIO4bT9YCjEppNbGIg0UbYarYRCQfp9/fiCXuostUnrk3cMdNLW384H7rHfZDlpAa6xhtGIbcqZg2KEddG0hCVus8Xf20vXhjcHWeGa1Zq2qyR6XU9zc1pWfTmRjrsgqZvCbuXY4cRrVkPdFktWTGXDq/s5fthFb483Y/dEPkIO9Szqvh4vu5r1E6NGyh2UiqITdKNmRDMYYmIOYVfKm8c28u7J9/AGvCnj0p22OkIhQ051zEeKZK3mgoDPxNhpcBEI4jl+DO/xBNa5wcAZ9/09rU//HH9H9lmsVYvOo3dXy+Bqkh4Pna+9SigYZPztd2Z8T4PRyLjb7qDhxk+Oahz6WCQdd0EySz5Zg2iAj39yHvXjTjeWsVrNXPaJWVlb2PkMOYy3qEcyWiVbiqs6FckjXTwBDyFCeCJib0wwvXkNc3BW1GLRRv8XkA56Bbi2nWPjF9c4+NkNDeFjlYYk3smRwVhXh6/Dldj6DoUIdnfhd+l/tkYx1NWFQyj00DQab7094bVdG98mmEMMftSqV2I+nFSNlBMVr1p+xYwkFQgtYCBhc4tMGzdn00QjXUa78FY6FJ2FnkkceY2lhul1Z/NBz2Ha+1w4bXXMawgnIT2//yU8weJYbNQrwOU3nRa8txbbaZlRwR0vd4xqSGOg00Xb6mcSn6BpVDQJDFZruMdnAmoWLqJvz258J4Z/3lomTSbY35+w1ntoYABvayu2yZMzHr8iN5JZ8omsd4/Hz//9bEvMLbJkxdkM9PkyssrjLflCJ/GMVLRKthSdoGcSR97pdfOJaVcyfdJEDhw7HisTMNazRU0GE/7QcD+K32TAbdf/lXVVJ/a1h8h/42itqhqD1UpgiOskmfVtmTQZU3Xynq2WSZPp3bE9HEZoNEIgYlFpGpVnnsnEBx7Cd+pU0nuoOP3RRW/hb6gQmi1GvJ5AbPEz6hbZveMEfm8wLb+3nq/8rOnOpCGHubpFRipaJVuKzuWSSYOKaCij1WShsbI+5lt3e7rHdGKRnpinvCZZs4vonzZbLLTPMuXMcFejLNFsNqZ8+SGMDof+CUbjabeJpmGZciZnPvgwfrc7qXXuPfZhOHwwFIqJedWi85j27e+y8LvfRrNYMDc2hueig8Fmw9zYmPW8FIUhKoS33nU+N3/2vITn+YcIfLKGFFFfeVS8e7o87Go+gdWmb/R4PH7+/OZBgnlYkM/GHTQSFJ2FDuFIFwhnh3YMdGLRLLruk7n1s3QXSGutdiyaZVDD5VJg/aJqDKEQ8/YPYNRxqGuVlUz6ysOx0L6Bw4c58o2vZfUsv6sD36lTBDoTuL+CQSb941cgEMQ6eTIme7iYkKm2FlN9vX5p3wTRL57Dh9Esp3+PmtVK7bLldK59bdi5tcuWK//3GMZsNrJjyzG8nvR82YnCAZP5yj39PuYsmsjelpODwh99aZQQKHaKzkKHcKTLLU3X8fDSv+drFzzAkvELdc8LJfn4NpTgh3lIM7BtZiWGBKujAZcLzWKJCZ7J6QyLaCLM5oSWsMnhDAu1M0Gru1CIkz/9Eb3bt2KsrIztjsZ865LAcopma8bTeOtfUHf5FZjq68NfHfX11F1+BY23/kXi+ShGHZ8vwPEj6X8d63UnguRx7709XuacOxGrTb+GS66Lo2OZorTQo1iMFmqtdt7vkLrHd7W/jzdw9bD9bk93SuvcgIHQqMeNZE6yuPUBuw2DPezDjlYXTBYPbqyqovrcRbjXrR12rHrhQkx2e8KEHIBAR4dugSu9mO+q+fPDvvM0szVVmGFxkkkzCkgcDpgqhBADo9qsebQoakGH9Oq7TKJ+0P5kC6uaQaNCs9EbyL5W+2hQX+HgU+fcSPOpHRyc9IZu3PruCSCPvsYtTdfRumY13RuT17MJuN3UXXYFBpMxYcLNaXHemrBD0tACV4nE+JTxmYyzNTNJHlKMPulUXownUThgsrj3s86pp2XrhwnvOVZixgtB0Qt6OvVdosQ3xUhUQiAYChadmAMsmbSA3R17kR376VwUtsKnHfNQ3Rukp0rj4CQr6xdV42zbxbVTPpawumA8Jmc9ZqczqSUcFeea5Zdw5OsP694nUYGreDEOBQIQDGGw2WKLpgabjZplF6lszRIimRDXj6vCO+BPOxwwUQhhKBTi/ebEiUxjJWa8EBS9oKdT3yUQDLBm7+8GNcWY2zCbSyYto6V9Nx0DndjN1fT5+7OKMBlNnFZHOOrHEDr9M0gSt94x0EnnqeMJqwvGE28Zp7KELY2NCRc70ylw1bpm9bBFztDAAAZNK/sGE6XGskunU1FhYff243R3e6iqsjK1qZ7ll88gEAjR5eoHA9TUVSRN1dcLIQRY/ZNNCa8xmTWWrDg773MaKxS9oMPwqJf4BCJvwMsPNz/Pmx++Gzu/w+PirWNvc/Gki5hZdw4bT26myzc2234lY0HDXG6feRNd3m5+tHN4PxG9uHWnrY66cRPpS1BdEMDorMe+aFFGlnEuBa5UPfLyQtM0rrx2Nr09Axza105vj4cP9rYRCoXX6D/Yl1lp2vi4d7erP6k7x+8LMtDnw2otCekbRknMKhr1cv30VTGXitFg5Pn9L7G9tSWhjz2+wFcxYrdU89jm72UUUz+vYQ62isSLmfZlyxl/x51ZCWi2Ba5UPfLy45Xfvz+o0FVvj5f3mwfXAMqmBksqH321PffkorFMSQh6FIsxnEAE6TWCLmYx19DYcPzdtM4LERr01QLJxTdbF0e2kSeqHnl54fMFkEmKdQ0lk7ZvyXz0AGeLxpL1n0OJCXqUsZ7anw+CpJftdtHEpVx25sWxsgdRChn2l2nkiapHXj74fAE+OtaF25V+B65MwwyXXTqdUCiE3Hk6schs0RDzJoyZmiuFIlWTaDPwJDAVsAKPSil/Fzn2HUBKKX8Y2b4LuBvwR857sYDjTkq6BbysmgVPiWWLRqm3OWIWuVFLbJGMlbA/VY+8tBlad8VgSNxseiiZhhlqmsaKK5q4YOX0cOPpENSMwTT9QpDKQv800C6lvFMIUQ80CyHeAZ4GmoB/BxBCTAD+DlgM2IANQohXpZSj4tNIVcDLaa1jfuNcgsEAbx1/Z4RHlx80g0YwNNxKd1od3DP/czRWOlPWhR9LqESh0SWX7j7pMLRGebpiDtmHGZrNRuobkxeDKzVSCfoa4Lm4bT9QDXwdiE/BPB94OyLgHiHEfmA+sDnZzR2OSkymwrw1LzjrXF7e+8aw/edPOpe/WXwHNbZqnty6uiDPHgn0xBzC8z53WjHXqbDDZP2a0wCNjfaEx0qFkZxjMBDkld+/j2w5ibuzn9q6CsTcCVx57ezUnafSxOf1cyRB6zY9LFYTPq+f2roKmvI8lpFmpP+9JhV0KWUPgBDCTljYH5ZSHgIOCSHiBb0GiC+20Q2kXMlyuQqXwLNq4pX09fnY2bYL10AnZqMFQrDp2Db2tx1mbsNMdrbuLtjzRwKb0UqlqRKXp5PGSieznbNYNfFKWluLLwQzHRob7SU7tygjPccNr+0bZDm7Xf1sWn+I/n5v3gpYuV39SX3mFZVm+vt82GvCiUFLVkxloM8f+1po7+jNyzhGmkL+LhO9KFIuigohpgAvAN+XUj6b4LQuIP4JdiC9LhQFIj6U8TeHX+TND4bGoafvatEwYNRM+IK+Qgw1azwBL19a9AUsRjPTJ02ky1W8UTuKkSdVd590I0tSkSyU0F5j5ea/PA+vJzDI3WO16hfWUiQn1aLoeOAV4ItSyteTnLoJ+KYQwkZ48XQW0JK3UebIrlN7dfenW4ArSIjgGBNzCCcJRX3lVpMFUIKuSJ9Cd/eJkiyUcGpTAxWVFioqdS5UZEwqC/0hwAE8IoR4JLLvainloO8nKeVJIcQTwHrCJXm/KqVM3MVgBHF7umnr0/ffFWM1xXjm1M+MJVIpFJkykk2Ph9Zdqa2rYMp0Z8mHEY40hlAmy815prW1u+AP7/H28OCGRxPGbRdT6GI0SchhraPCXEGft49OrxuHtY4LzjqXVROvTBqiWOwoH3r+GepDjzJv8aSCNIGIRtOcNdVJpzv9WPRipMA+dN2GDsW5dJwB/X5P0iQch61uBEeTG0GC3HvuXcypn8mxnuO4vJ2ECNHhcfHy3jd4fv9Loz1ERZGx7NLpzFs8CXuNFYMh7NOet3hSwSznWOs2S0nmNI46Jf9TrbXaqbc5aB/Qj0nv6O8sGitdw8CWk9vY49qne3xn2y6un76qqOLPFaPLWG96rMiMkhd0o8GIzWyFBB59b8hLsbjSg4TYeDJxadBoQ49aqz32pxJ3RTrEVyxUFC8lL+jP73+JY93pFwIaDZZOOI99roNpV03U0HTdSA5rLWuPvkVL255Y3ff5janT/xUKRWlQ0j70YinStfmjZmY7m9I+P9GaQIW5kreOvUOHxxXzra/7cIPyrSsUZUJJC3q6RbpGm2AoSJ/fw4qJF6Kl8StxWuu4eNKF1NscGDBQb3Nw1YyL6ffpRw3sbNuFNzD21wgUCkVulLTLJVWRrnQwMDIu9gOdB/j83DtYn0axsPmNc7ml6bpBPVKN1UFe2b9e9/yobz1aK16hUJQmJW2hR/uN5oLdPDJJO25fN99t/nHScwwYmFg1gYsmLsUb8MYaeliMFhy2WhxW/RDMoc2y9fAGvLT2tStLXqEoYkraQodwv9HKSjN/PrIt1m/UZqrgWM/x1BcDXb5uLJoFr05YY6L92ZKqaUWIEMd7T/LNTd+ONYeOLnhaTambZesRCAZ4fv9Lgxpoq4VUhaI4KXlBN2pG/nLhp7jijMuG9RuNNpV2WGvp9ffjSdCSLpFoL52wiE0nm/EER76GSnTBE+CWpuuA5M2yE/H8/pcGvQT07qtQKIqDkhf0KPH9RoFhTaV/e+CPSXuQWjQLVeZKOj3uIUJpSMvvHY9ZM2MAvEEfddYaOj1dWc7qdDIR6DfLThaHniwKSCUpKRTFR0n70FMR74O+acY1LJ1wXsJzvUEvgVCA8ycs4sEl93NL03UYNSMrp1yU8XN9QR/eoI+lE87jwSX347Q6sp5DdMEznvh5JSNZFJDefRUKxdimrAU9HqNm5DZxIw5L4touXd5u/nzyPV48dLqZsdNWh9WYXau0/Z0Hc164dVhrs662GI0C0iOdhVSFQjG2UIIeh8VoYcG4uSnP23ZqBz3entM7sqxYGbWCb5pxDSsnL4/FlZsN6Rf3P8cxPWu3SLKXSbKFVIVCMTYpGx/6UOJjuC1GS2z7E2dfAYRFu9Or79vu9HbxL5v/k3Mb5tLnH8i6sFfUCh7q+642V/LioVfZ2baL9gFXwkYcVqOVW87JbeEym4VUhUIxNik7QR8apldnqaXSUkm/r39Q2N4DS/6Or7/zeMIIF7enizePbcxpLEOt4PiF23iBX3v0Ld2WeReesYQKc24FlTJdSFUoFGOXshP0oWF6Lm8nLu/phcFo2F4wFMCAbg35vGAzWvnEtCuSnmMxWqi12rlk8kWEQvB+x56CWdFDo4AUCkXxUVaCnkmxrh1tuwoaX+4JeOnx9lFh0rew9RJ+5jbM5JLJy3HaapUVrVAohpFS0IUQZuBJYCrhBtCPAu8DTxEuc9ICfEFKGRRC/A6oB3xAv5Ty6sIMOzsyKdbl9nRTY7bT5StM6F6qKBK9hJ+3jr2DZjCqhB+FQqFLOlEunwbapZQrgKuB/wL+A3g4ss8AXB85dwawXEq5cqyJOSQP0xuKxWhhIEHmaD6YWz8roZWdKuFH1VtRKBR6pCPoa4BH4rb9wHnAm5HtPwCXCyHGA3XA74UQG4QQn8jrSPNAJjHfnoBn2IKoEQ2r0RorWXvxxAuZWHVGVr72UJJrVMKPQqHIhpQuFyllD4AQwg48BzwMfEtKGY2j6wZqAQvwbeC7gBN4WwixSUp5KtG9HY5KTKaRKQDV2Bh2b9xdfxuhzX7e/ODdhOdGwxiH4qis49+ueJA+fz8OWy2/3Pkbjh8/kdV49rj2UOP4FFbTcCu9xm+lodJJa1/78HlUOpk+aaLuddE5ljJqjqVDOcxzpOeY1qKoEGIK8ALwfSnls0KIx+MO24FO4CTwQymlHzglhGgGBJBQ0F2uvqwHngmNjXZaW09btddOuZp3j+gX1aq11OBOEH/e3ufiw4/aaKysp62nm3cPb8t6TK19HRw4djxhZMkc5yzW9Q2vLTPbOYsulwcYPPahcyxF1BxLh3KYZyHnmOhFkdLlEnGlvAJ8WUr5ZGR3sxBiZeTvVwPrgcuB/4tcUw3MBXbnNOoC8eKhVxNGsCxonJOwtkr8Qqbb051T44xUi6JDs0frbQ5WTl6uEn4UCkVC0rHQHwIcwCNCiKgv/T7gCSGEhbBoPyelDAghrhJCvAsEgYeklG0FGXUOJFtwtBmtXDd9FZrBmLKueK3VntSaB7AYLHhD+guYqVLrVcKPQqHIlHR86PcRFvChXKJz7v35GFQhSbbgGI0NTycd3mK0ML9hTtLSuXpiXm9zxO41tPyAQqFQ5EJZJRZB8j6jiWqrJBLcW5qu42DX4bS7H9VZa3hg8b1UmCpSdglSnYQUCkWmlF21xUwqDKaqK27UjHx58b1cPOlC6qw1GDBQZ61J+Gy3p5t+vyeWNNThcREiFCs38Pz+l2LnpnOOQqFQxFN2Fjrkt8KgUTNyq7iRG2dcg9vTTYXJymObv5fwC6DCZE3ZJQhQnYQUCkXGlKWgF2LBMb64VbJmzf1+T1pJQ6nOUYW0FArFUMpS0KMUqsJgsi+AQCiQ0ocPJDwnlw5FCoWitClrQS8Uyb4AjBiTWvDR8xKd0+vv57cH/qgWRxUKxTCUoBeQRF8A6fjwo39/98TmQUXCPAFPTOhvabouFvpY48+ur6lCoSgdlKCPAun48I2akeunr2J7a4tu1ccdrbsIBAPsat+Dy9NJQ6WTOc5ZynJXKMoYJeijSCofvtvTTafHrXusw+MalNTU2tceq/2i6qUrFOVJ2cWhFxPJ6rdrCX51ql66QlG+KEEfwyRLggoS1N2v6qUrFOWLEvQxjl7VxYsnXYjDom+5p6riqFAoShflQx/jJFpATacipEKhKC+UoBcJQxdQh4Y+NlY6mR2JclEoFOWJEvQiZajlPn3SxEgnI4VCUa4oH3qRE7Xc9XqMKhSK8kIJep7xBry09rWr0EGFQjHiKJdLnlANKRQKxWiTUtCFEGbgSWAqYAUeBd4HngJCQAvwBSllUAjxNeAawA/cL6XcVJhh5wePL4C7x0OF1US/x09tdbgeirvHQ221Fas5fSGONqSIEm1IASpzU6FQjAzpWOifBtqllHcKIeqBZmAb8LCUcp0Q4ofA9UKIw4T7jC4FpgC/BpYUaNw5EQgG+dXa/TTvbaW9y4NmgGAILCYDGAx4fUGcdgszz3Jy+xXnUGk1x8RfT+iTNZ5WDSkUCsVIkY6grwGei9v2A+cBb0a2/wBcCUjgFSllCDgihDAJIRqllK35HHA++NXa/by25cPYdjAU/tPrDxH+6ICObi8bW06yde8pGusq6e334ur24qyxsrCpkVsvnYFRCy9BJGs8nW1DCtVAWqFQZEpKQZdS9gAIIeyEhf1h4FsR4QboBmqBGqA97tLo/oSC7nBUYjIVzr884PXj6vIw4PXT2BjOnnT3eGje25bBPYIcPdUT227v8vDalg8JYeCem+djs5io8VtpqHTS2tc+7PrGSifTJ01MOwolEAzwi+2/ZvOHO2jr66Ch0smSyfO5c8HNKX3x0TmWMmqOpUM5zHOk55jWoqgQYgrwAvB9KeWzQojH4w7bgU6gK/L3ofsT4nL1ZTbaNIl3qXR0eWh0VDBvmpMQsFW20tmTewTK61uOsm3vqZi1Psc5K1btMJ7ZzlmR+PD0YsTX7P3dIF98a187L+99g74+X1JffGOjndbW0q7houZYOpTDPAs5x0QvinQWRccDrwBflFK+HtndLIRYKaVcB1wNvAHsBx4XQnwLmAxoUsr0TeE8MtSlcsrVz+vvHcv7c6LWOsCtl+beeFr54hUKRS6kY6E/BDiAR4QQj0T23Qc8IYSwALuB56SUASHEeuAdwvHtXyjEgFPh8QVo3juybvvmvW3cfMn0nBtPF8IXr1Ao8key4IixQDo+9PsIC/hQLtE59+vA13MeVQ64ezx0dI1sCryrewB3j4dxjsqcGk9H65+naiCtUChGlqFuXL3giGwIBoN8+9v/xv79+zCbzXzlK48wefKUrO9XcpmitdVWnDW599e0mtP/0TjstlgMey4kq3+uqigqFKNH1I3b3uUhxGl366/W7s/pvuvXr8Pr9fKjH/2ce+65l//6r+/kdL+SE3Sr2cjCpsac71NpNWJJU9QXNjXk7fNLr/75ysnLVRVFhWKUSObGbd7bhscXyPreO3ZsY+nSCwGYO3cee/bszvpeUKKp/7deOgMI/7Bd3QPU19pwdQ3g12/yo4urx5f0uAFw1thY2NQQe14+SKeBtEKhGDmSuXHj3a3Z0NvbS1VVdWxb0zT8fj8mU3bSXJKCbtQ0br+8iWuXTeXAh2427v6I1s6BvN2/vsbKfZ+cT6OjsmALI7n44hUKRf6IunHbdUQ9V3drVVUVfX2nw7dDoVDWYg4lKuiBYJBfvr6PjTtPMODNwCxPk4VNjUwepxYoFYpyIOrGjQ+FjpKru3XevAW8/fZ6LrvsClpadjJtWm5f+yUp6L9au5+1eYo7NxsN2KssdHZ7cNjz72JRKBRjn6Fu3HxpwcUXf4zNm//MPfd8nlAoxEMPfS2n+5WcoHt8AbbKU3m7ny8QorPLw+JZjXxm1Uwqrea83VuhUBQHUTfuzZdMz2scuqZp/OM/PpSHEUbul7c7jRHcPR46uvPbXCIIbNrdymPPNBMI5t+Fo1AoigOr2ci4Aq6d5UrJCXpttRWnvTBRIUdP9fCLVySnXH05hSopFApFISg5l4vVbGSRGKe7gJEP1m87wVvbTlCfp0wxhUKhyBclqUS3XjqDC+aOL8i9ozWD85UpplAoFPmiJAXdqGl89qqZ1OehBEAqcs0UUygUinxRkoIOYdfLgnMaCv6caKaYQqEofbwBL6197XgD+Q28yBcl50OPxzACz8hXYS6FQjF2CQQDPL//JXa07sLl6cRhrWN+Y7jfQapOYqnYtauFH/zgCf7rv36c8zhLStDjaxUDbNtX+P4alTYTJuNIvDoUCsVo8fz+lwZ1EuvwuGLbyTqJpeKZZ/6HP/3pZWy2ipzHCCUi6Hq1ipvOrNOtvZBvjp7q4Vdr93P75U3A2C+Ar1AoMqOQncQmTZrMN7/573zjG/8vlyHGKAlBH9pyrr3LwzstH43Y85v3tnHDirP5zfpDeS+Ar1AoRpdCdhJbufIyTpw4nsvwBlH0gj4aLeeG4uoe4NlX97Gx5WRsX3y/0aj1rlAoio9i6iSWlqALIZYCj0kpVwohFgE/JNzGfhtwn5QyKIT4HVAP+IB+KeXVhRp0PKPRcm4oDruVPYc7dI9F+40q94tCUZxEO4nF+9CjjLVOYil9AUKIB4CfArbIrh8D90spVwBu4PbI/hnAcinlypESc8hfy7lcmHmmA1eC+jEqrFGhKH6KpZNYOhb6AeAm4BeR7clSyo2Rv78NXC+EeBWoA34vhKgD/k1K+WLeR6tDslrFhcZiMnDR/IncfMl09hxxFaQAvkKhGH0K2UnsjDMm8uMfP5WXexlCoVDKk4QQU4HVUsoLhBAbgQellG8KIb4P2IGHgE8B3wWchIX+IimT12hRR3MAAAvCSURBVLH1+wMhkyl3V0QgEOTJ3+/i3ZYTtHb2k8aU8sY4RwUXzD2DYCjEixsODTt+3Ypp3HXDvJEbkEKhKAd0Y6WzEXRBWLgDwGagFngAsEgpeyPn/x/wPSnl+mT3bW3tzqv0enwBjrf28I2n38vnbdPisvMmYTAYdAvgj0SUS2OjndbW7oI/ZzRRcywdymGehZxjY6NdV9CziXK5Bvi8lPK4EOJ7wB+Ay4EvAtcIIaqBuUBu7auzZMA7OnVVtu1r59G7lua9AL5CoVCkSzaCvg94WQjRB7whpXwZQAhxlRDiXcL9IB6SUhY+TZOwVd7RNcBr733Ijv1tI5JMpEd89+9sO4ArFApFLqQl6FLKD4ALIn//PfB7nXPuz+vIUhCfHTpaIh6PWvxUKBSjTdEmFg3NDh1tcu3+rVAoxj5Bjwe/242pthbNOvYMuKIU9LGQHRolvnORQqEoTUKBAK1rVtPTvBV/Rwcmp5PqhYtovOU2DMbsDDm/38+//us/ceLECXw+L5/97F+xfPklOY2zKAV9LGSHAiw6p4HPXj0Te+XYyRRTKBT5p3XNajpfezW27W9vj22Pu+2OrO75pz+9TE1NHY888g3c7k4+97k7chb0oqwaNdrZoTaLEZvFSPO+Nv75qc08+9peAsHgqI1HoVAUjqDHQ0/zVt1jPc3NBD3ZGZcf+9jl3HXXPbFtozF3+7ooBT2aHTrSaBpMcFQw4A0w4A0QQvUWVShKHb/bjb9Dv1aT39WB3+3O6r6VlZVUVlbR19fLww9/mbvu+ttchgkUqaBDuBH05Ysnx/qGapEwe0sBnUjBIJx09eseU71FFYrSxFRbi8np1D/mcGKqrc363h99dJJ7772Hq676OFdeuSrr+0QpWkE3ahq3X97E/OnhOsTBSM6p1z8641FFuBSK0kSzWqleuEj3WPXChVlHu3R0tPOlL32Rv/3be/nEJ67PZYgxinJRNIrHF2DHgfbRHgag4tAVilKm8ZbbgLDP3O/qwORwUr1wYWx/Njz99M/p7u7mqad+ylNP/RSAb3/7CaxWW4orE1PUgj5Wol1AxaErFKWMwWhk3G130HDjJ/MWh37//f/A/ff/Q55GGKaoBT0a7TLSmaJTxlXTN+AfVoRLoVCUNprVimXcuNEeRkKKWtBHuha6o9rCeTPHceulM/AHQqoIl0KhGFMUtaADMcs4WrbWYjYWrOKipp2uWGk1G1URLoVCMaYoekGPRrtEy9ZWV1r4zfqDNO9to71rIOX1NotGKAQeX+rEINX4WaFQjGWKNmxxKFGLudJq4vbLm3j0rqVcNHdCyusa6yrTEvN4VMy5QqEYi5SMoOux54gr4bHaKjMXn3sGvf36zZ2ToWLOFYryxOcL4Hb14xujBl3Ru1wSkSqksavXx879Hbh6Mhd0FXOuUJQXwWCQjWsPcGhvGz1dHqprrJzd1MCyS6ej5dBiMhAI8Nhjj3L06GE0zchDD32NSZMmZ32/krXQUxXwCgGuLK1sFXOuUJQXG9ceYOeWY/REjMSeLg87txxj49oDOd337bfDbZd/8IMn+au/upvvfe8/crpfyQp6Pgt4aYZwi+36GhuXL56sYs4VijLC5wtwaK9+R80P9rbl5H65+OKVPPDAV4FwXReHoz7re0GaLhchxFLgMSnlSiHEIuCHgAfYBtwnpQwKIb5GuIG0H7hfSrkpp5Hlgajw7jjQTqurn1CSc+uqLXQmcL9csnASVy2ZomLOFYoypK/HG7PMh9LT7aGvx0utoyLr+5tMJh599Gu89dY6Hn30sazvA2lY6EKIB4CfAtECAz8mLNgrADdwe0TkLwGWArcB/53TqPJENKTxvx+4lPNnj094nmaA+dOd/PNfnc/HFk2ivsaGZjhtkd9++TmMc1QqMVcoypDKagvVCdy31XYrldW5N7h5+OF/4pe//DWPPfYo/f36FV3TIR0L/QBwE/CLyPZkKeXGyN/fBq4HnMArUsoQcEQIYRJCNEopx0afOGD/h50JjwVD8Nb2k1jMJu68UuD5WEBlgSoUCgDMZiNnNzWwc8uxYcemNjVgzkEj/vjHl2htPcWdd34Om82Gpmk5LbKmFHQp5a+FEFPjdh0UQlwipXwTuBaoAmqA+LKH3UAtkFTQHY5KTKbCC+aJtl46ulMvgO440M7dN1fQaDGR/Trz6NHYaB/tIRQcNcfSoZjmef2nzqWiwsLelpO4O/upraugae4Errx2NpoxsQCnmuPNN1/Hgw8+yP3334Pf7+fhh7/K5MkNWY8zm7DFzwHfjbhiNhP2pXcB8SO3A4lN4gguV18Wj88cR20FTnvqIl5tnf0c+KC9KFP6GxvttLZ2j/YwCoqaY+lQjPM876KzmH/+ZPp6vFRWWzCbjbR39CY8P905Pvzwo4O207km0YsiG9v+GuDzUsprgHrgVcKul6uEEJoQ4kxAk1LqLwuPAjaLKa2IFxVfrlAokmE2G6l1VOTkZikk2Vjo+4CXhRB9wBtSypcBhBDrgXcIvyS+kL8h5of4Il6Jaryo+HKFQlHMGEKhZMF8haW1tXtEHh7/6ePxBejoGuC1LUfZcaBjWE1zYw4LEqNJMX7CZoqaY+lQDvMs5BwbG+0Gvf0lm/qfCKvZyBn1Vdx51Uw8PhXNolAoSoeyE/R4VE1zhUJRShSnf0GhUCgUw1CCrlAoFCWCEnSFQqEoEZSgKxQKRYkwqmGLCoVCocgfykJXKBSKEkEJukKhUJQIStAVCoWiRFCCrlAoFCWCEnSFQqEoEZSgKxQKRYmgBF2hUChKhJIpziWEWAo8JqVcKYQYB/wEcABG4DNSygNCiLuAuwE/8KiU8sXRG3F2DJnnucAPCc9nL/DXUspgsc5TCGEGngSmAlbgUeB94CkgBLQAX4jM8Wv/v727CdGqiuM4/pmScGNBLWxhMJv42wstKorypVkIYS6CNolJLapFZGW4EMpBCKKQ3pyKNrYpioHMhS0CKVQsCje2ifhTLXLRC6VmBlmN1eLcB56GO4Q2zNO9ne/qOffexfnyP/w459yHe5TDVmaUQ8sPj6LPZ8scjkfxIs4oJ4DdnZnfdbWOtHtm5t7m3gY8lJk3Ne1Oes5Ry4+NMHt6MUNvjsPbhcXNpR14IzNXYxuWR8SleBgrcCueiohOHU/U4rkdT2TmSmVAreu450Ycy8xVWIuX8By2NdfGcHtEXItbcCPW4+UR9fdcaHPcqQTcBPZga8frSLunZhJyr1JLHfdscxxp9vQi0PEl7hhqr8CyiHgPd+EAbsCHmflrZp7EF7hmoTv6L5nteQQXR8SYco7r77rt+RYmh9ozuA4Hm/a7WIOV2JeZf2bmUSyKiH8+Y/C/QZvj+sz8pGkvwmndriMtnhFxCZ7G5qHrXfZsq+VIs6cXgZ6ZbythNmAcJzJzjbKc3YoLcXLomVO4aKH6OB+0eH6OKXyGpcrg6axnZv6cmaciYgl2KzOcscwcfJ9i4NIrx8z8BiLiZmzC8zrsSKvnJF7Fo4rLgM56zjFex40we3oR6C0cw97m9zu4Hj8ps9gBS/DjAvdrvtmJVZm5HK/hWR33jIjLsB+vZ+ab+GPo9sClb44i4k7lfci6zPxexx35u6cy+bgcr2AaV0bECzru2VLLkWZPb16KzuID3KYMpNX4FIfxZEQsVvabr1BesnWZ48pgga+V5V5nPSNiKfZhU2a+31w+EhETmXlA2afcryxZd0TEM1iG8zLzh1H0+Wxpc4yIjcoLs4nMPN482tk6Mmctr2rujWM6Mzc3+8ud9JzDcaTZ09dA34JdEfGAstTZkJknImIKh5SVyeOZeXqUnZwH7sN0RMzgN9yfmd922PMx5d8BkxEx2Jt8BFMRcYGytbQ7M89ExCF8pDg+OJLenhuzHc/H1fgKeyICDmbm9g7XkfZars3MX4Yf6uF4vccIs6d+PrdSqVR6Ql/30CuVSuV/Rw30SqVS6Qk10CuVSqUn1ECvVCqVnlADvVKpVHpCDfRKpVLpCTXQK5VKpSf8BWvrClWyzu3EAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"def plot_observations(X, labels, title):\n",
" '''Plot observations coloring according to the given labels.'''\n",
"\n",
" df = pd.DataFrame(X, columns=[f'x{k+1}' for k in range(X.shape[1])])\n",
" df['label'] = labels.astype(int)\n",
" for label in set(labels):\n",
" x1, x2 = df.query('label == @label')[['x1', 'x2']].values.T\n",
" scatter(x1, x2, label=label)\n",
" plt.title(title)\n",
" legend()\n",
" show()\n",
"\n",
"plot_observations(X, labels, title='Original clusters')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Algorithm"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### E and M steps"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"def e_step(X, probs, mu, sd):\n",
" '''Perform the Expectation step (E-step) for the EM-algorithm.\n",
"\n",
" Args:\n",
" X (array[n_observations, n_features]): the feature matrix\n",
" probs (array[n_groups]): the current best estimation for the parameter\n",
" mu (array[n_groups, n_features]): the current best estimation for the parameter\n",
" sd (array[n_groups, n_features]): the current best estimation for the parameter\n",
"\n",
" Returns:\n",
" T (array[n_groups, n_observations]): the probability an observation belongs to a group\n",
" '''\n",
"\n",
" # The probability density function for the multivariate gaussian\n",
" pdf = scipy.stats.multivariate_normal.pdf\n",
"\n",
" T = array([probs_ * pdf(X, mean=mu_, cov=sd_)\n",
" for (probs_, mu_, sd_) in zip(probs, mu, sd)])\n",
" T /= sum(T, axis=0)\n",
"\n",
" assert np.allclose(T.sum(axis=0), 1.0)\n",
"\n",
" return T\n",
"\n",
"\n",
"def m_step(T, X, probs, mu, sd):\n",
" '''Perform the Maximization step (M-step) for the EM-algorithm.\n",
" \n",
" Update ``probs``, ``mu`` and ``sd`` inplace based on ``T`` and ``X``.\n",
"\n",
" Args:\n",
" T (array[n_groups, n_observations]): the probability an observation belongs to a group\n",
" X (array[n_observations, n_features]): the feature matrix\n",
" probs (array[n_groups]): to be updated to the next best estimation\n",
" mu (array[n_groups, n_features]): to be updated to the next best estimation\n",
" sd (array[n_groups, n_features]): to be updated to the next best estimation\n",
" '''\n",
"\n",
" for j in range(len(probs)):\n",
" probs[j] = mean(T[j], axis=0)\n",
" mu[j] = sum(T[j].reshape(-1, 1) * X, axis=0) / sum(T[j])\n",
" sd[j] = sum(T[j].reshape(-1, 1) * (X - mu[j]) * (X - mu[j]), axis=0) / sum(T[j])\n",
"\n",
" assert np.allclose(sum(probs), 1.0)\n",
" return (probs, mu, sd)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Initialization, Iteration and Post-processing"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"# Initialize the variables for parameter estimation\n",
"# probs_hat: 1/n_groups for each group \n",
"probs_hat = ones(probs.shape) / len(probs)\n",
"mu_hat = random.random(mu.shape)\n",
"sd_hat = ones(sd.shape)\n",
"assert np.allclose(sum(probs_hat), 1.0)\n",
"\n",
"# Scale X to avoid numerical problems\n",
"scale_factor = np.max(X, axis=0)\n",
"XX = X / scale_factor\n",
"assert np.allclose(XX.max(), 1.0)\n",
"\n",
"# Iterate alternating between E-step and M-step\n",
"for iteration in range(1000):\n",
" T = e_step(XX, probs_hat, mu_hat, sd_hat)\n",
" m_step(T, XX, probs_hat, mu_hat, sd_hat)\n",
"\n",
"# Post-processing: scaling back mu_hat and sd_hat\n",
"mu_hat = mu_hat * scale_factor\n",
"sd_hat = (sd_hat * (scale_factor ** 2)) ** 0.5"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Results"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"* Original parameter values (ground truth)\n"
]
},
{
"data": {
"text/plain": [
"array([0.20602828, 0.42700046, 0.30769962, 0.05927164])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"array([[163, 195],\n",
" [170, 213],\n",
" [191, 216],\n",
" [245, 215]])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"array([[2, 2],\n",
" [5, 4],\n",
" [6, 2],\n",
" [6, 4]])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"\n",
"* Estimated parameter values\n"
]
},
{
"data": {
"text/plain": [
"array([0.05371094, 0.32170901, 0.20373055, 0.4208495 ])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"array([[245.59807655, 215.79282986],\n",
" [191.2382189 , 216.08071049],\n",
" [163.09387915, 195.16652376],\n",
" [170.57574173, 212.59164509]])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"array([[6.81597025, 3.41986827],\n",
" [6.09132375, 1.94990394],\n",
" [1.99913102, 1.9359254 ],\n",
" [4.93429354, 4.15610579]])"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"\n",
"Note: the order of the groups can be different\n"
]
}
],
"source": [
"print('* Original parameter values (ground truth)')\n",
"display(probs)\n",
"display(mu)\n",
"display(sd)\n",
"\n",
"print()\n",
"print()\n",
"print('* Estimated parameter values')\n",
"display(probs_hat)\n",
"display(mu_hat)\n",
"display(sd_hat)\n",
"\n",
"print()\n",
"print()\n",
"print('Note: the order of the groups can be different')"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAEJCAYAAAC5Tb0qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl8FEX+//HXTEKuSeIBCisIImitq0lAUUBBwINdT/zigch9xKCoKIqggrocooDKiiAe3CDKjfexqyByiIhARChE5fD+gUqOmSQmmd8fE0IUDImZzDTM++mjH6a7Zqqr2vaTT6q7q11+vx8REXEed7gbICIih6YALSLiUArQIiIOpQAtIuJQCtAiIg6lAC0i4lAK0CIiQWaMaW6MWXaI7VcZYz42xqw2xqQfrh4FaBGRIDLG3Au8AMT9YXsN4EmgPdAGuNkYU6e8uhSgRUSC60ug4yG2nwFst9b+Yq0tAD4EWpdXUXQ1NK7KCrL26vFGOUizlEOd85FpXeaicDfBMWKSa7qqWkdqgzYVjjmbdi4vd3/W2oXGmFMOUZQM7Cuzng0cU15dyqBFREIjC0gqs54E/FreFxyZQYuIhJLLVeUkvCK2AKcZY44HcoALgXHlfUEBWkQinstVfYMJxpibgERr7XPGmIHA2wRGL6Zaa78t77sK0CIS8dwEN4O21u4AWpT8/GKZ7a8Cr1a0HgVoEYl4IRriqDQFaBGJeO5qHOKoCgVoEYl4Ts2gnflrQ0RElEGLiES5osLdhENSgBaRiOfUIQ4FaBGJeG6HBmiNQYuIOJQyaBGJeC6H5qoK0CIS8aLcCtAiIo7kCvKj3sHizF8bIiKiDFpERI96i4g4lO6DFhFxKKfeB60ALSIRz6kXCRWgRSTiaQxaRMShnDoG7cxfG9WouLiY4aPH0KV3Or0y+rNr9ze/K1+weCmduvemS690lq9YCcBWu43OPfqQ3n8AXp8PgOemTmfDpsyQtz+YIvlYpDQ5gykvjQfg5AZ1mb5gAtPnT2DoyIGl/7P2G9CDOUsnM3PRRM5K+/tBdbS5+HxefOVZZi2exLU3XglA7TonMHPRRKbNe4oTa9cC4Ir/u5R/XXVRiHpWdZF4XrhdrgovoRRxGfR7yz4gP7+AOVOfZ2PmZ4wd/xQTHh8DwJ49e5nz8nxenjmV/IICuvftR8vm57L41dd46P7BrP1kPavXrCUtNYVvv/ueJqkpYe5N1UTqseiV0ZkrO7bH5w0EkkHD+vP0uCmsW7OBoaMG0q59K77/9geatWhClw79qHPSiTwxeQQ3XZ1RWkd0dBSDHuxP56sy8PnymLlwIsv+t4r2V7Zj2uS5uFwu2l/Zjnmzl9L2kgsY1P/hMPW28iLxvHDqGHS1Z9DGGEdl6es3bqTV+c0BSEs5i8+3bC0ty9z8OU3TUomJiSEpMZH6J9dj2xfbSYhPwJeXh8/nIz4+jmenTCO9V49wdSFoIvVY7N71LXdlDC1dPyPldNat2QDAh8s+okWrc2jaLJVVH3wMwA/f/URUdBTHHX9M6XcaNm7A7h3fkp2VQ+FvhXz68SbOPjcVb66P+IQ44hPi8Hnz6N73BuZMWxDaDlZRJJ4XLperwksoVUvwNMacaoxZYoz5BvjKGLPLGPO6Meb06thfZeTmekn0JJauu91RFBYWApCTm0tioqe0zJOQQHZOLjd1uo658xbw674sah5/PAnx8WyxlhGPjuX95StC3odgidRj8d83P6CwsKh0vez/dN5cL4lJHjxJCeRk5x7YnuMlMenAsUpM9JBdpjw310dSsoc3lv6X5hecw7ktmrDmw3WcfEpd3C43Q0cNpOONV1Rzz4IjEs8Lpw5xVFd2+wIw2lpbz1p7irW2PjACmFZN+6swjyeBXK+3dL3YX0x0dGCkJ9HjwZt7oCzX6yU5KZETatVizKjhDB44gCkzZ9OnZzfmLVzCsCGDmDv/yMqOytKxCPAXF5f+nOBJIDsrh9xsLx5PwoHtiYHt++Xk5OJJPFDu8cSTnZWDz+vjoUGP8fDgsXTrewPPT5hF39u68siw8bRu15L4+LjQdKoKIvG8cFXin1CqrgAdZ639qOwGa+2aatpXpTRNS2XFytUAbMz8jNMaNSotSznzH3yyYSP5+flk5+Tw1dc7aNzo1NLyFStXk5ZyFslJSRQUFADg8+WFtgNBpGMRsHXzdpq1aAJAq7bNWb92E5+uy+T8Nuficrmoc9KJuF0ufv1lX+l3vt6+k/qn1CP5mCSia0RzTvM0Nn6yubS88ekNyc/L55td3xEbG4Pf7ycqyk2NmBoh719lReJ54Xa5K7yEUnVdJNxojJkKvAXsA5KAy4FN1bS/Cru4bRtWf/QxXXvfjB8/Ix58gBlz5lK/Xj3atWlNl07X0yP9For9fu64NYPY2FgAioqKWLBkKeNGjwSgZfPz6NIrnbYXtgpnd6pExyJg3MiJPPToIGrE1OCr7Tt5943lFBcXs37tJmYvnoTL7eaRYYE7Pi7vcAnxCfEsnPsq40ZMZPKscbjdLhbPe4OfftxTWmff27oyauiTALyy8G1mL57E5kxL1r7ssPSxMnReOIfL7/cHvVJjjAu4BmgFJANZwEpgsbX2sDssyNob/EbJEa9ZSsdwN8Ex1mUuCncTHCMmuWaVxx2uO6dXhWPOgk+mhWyco1oy6JIgvLhkERFxtCg9SSgi4kxOnSzJmb82REREGbSIiFPn4lCAFpGI59QhDgVoEYl4Tp2LQwFaRCKeMmgREYfSGLSIiEMpgxYRcSiNQYuIOFSwMuiS+e8nAWlAPtDXWru9TPk9QGegGHjEWlvu09Z6UEVEIl4QJ+y/hsBsni2BIcDj+wuMMccCdwAtgfbA+MNVpgAtIhEviBP2tyIwi+f+KZablSnLBXYCnpKl+KBv/7Fdf6k3IiJHkSBm0MkEpljer8gYU3YoeTfwObAeeOpwlSlAi0jEC+IbVbIIzH+/n9taW1jy82XA34CGQH3gGmPMeeVVpgAtIhI8Kwm8nARjTAsgs0zZL4APyLfW5gG/AseWV5nu4hCRiOcO3l12i4FLjTGrABfQyxgzENhurX3FGHMJsMYYUwx8CLxbXmUK0CIS8aLcwRlMsNYWA/3+sHlrmfKHgIcqWp8CtIhEPKc+6q0xaBERh1IGLSIRz61HvUVEnMmpQxwK0CIS8TSbnYiIQzk0PitAO13xb7+FuwmOsS5zUbib4BjnpV0f7iY4xoav36tyHcqgRUQcSvNBi4g4lC4Siog4lIY4REQcyqHxWQFaREQZtIiIQ+kioYiIQymDFhFxKIfGZ81mJyLiVMqgRSTiBWvC/mBTgBaRiOfUIQ4FaBGJeE69SOjMvF5ERJRBi4joPmgREYfSZEkiIg4V5XZmgNYYtIiIQymDFpGIpyEOERGHcugIhwK0iIgyaBERh3JofFaAFhHRk4QOUVxczPDRY+jSO51eGf3Ztfub35UvWLyUTt1706VXOstXrARgq91G5x59SO8/AK/PB8BzU6ezYVNmyNsfbJs2f06f/gMO2v762+/SqWdfbup9M/MWLQHA6/XS97Y76ZZ+C9u2fwnA+o2bmDrrxZC2uTpE+nlRI6YGo8c/wMxFT/PMzDHUP6UuzS84mzlLn2Hmoqfpf3fvg76TmOTh6WmjmfLyeCbPGkvNWscBcM0NlzFz0dPcP/zAeTV6/AN4EhNC1p/KclXin1CKuAD93rIPyM8vYM7U57nztlsYO/6p0rI9e/Yy5+X5zHphMpMnPMn4ic9QUFDA4ldf46H7B9P6gpasXrOWPXt/5tvvvqdJakoYe1J102a/yL9HjyG/oOCgsieensRzTz3BjGcnMnPuPLKyslm1dh1tW1/A/ffcxeJXX8fv9/PivAV0ueHaMLQ+uCL9vOh44xV4vT66d7yNRx+ewJB/38Fd92Uw7O5H6d7xNpo1T6Oxafi771x97T/Zbr+mT6c7eef1ZfTI6ATAlR3b0+Pa2zmhTi2SkhNp3a456z/OJDfHG46uVYjLVfEllCIuQK/fuJFW5zcHIC3lLD7fsrW0LHPz5zRNSyUmJoakxETqn1yPbV9sJyE+AV9eHj6fj/j4OJ6dMo30Xj3C1YWgObluXZ4YPfKQZac1akR2Ti75BQX4/X5wQUJ8PD5fHr68POLj43jjnf9yUZvWxMbGhrjlwRfp50Wjxg34cNlaAHZ+tZuGjeqzdfN2ko9NIrpGNDGxMRQXFf/uO1/Yr0nwBLJiT2IChb8VAZDnyyM2NoYa0dH4/X46XH8Zi156LbQdqiS3y1XhJaTtCuneHCA310uiJ7F03e2OorCwEICc3FwSEz2lZZ6EBLJzcrmp03XMnbeAX/dlUfP440mIj2eLtYx4dCzvL18R8j4EyyXt2hAdHXXIssanNqRzr3Q6dunBhRe0JDkpiRbnnsPen39m/uKlXNvhKt7/YAWnN27M8MfGMW32kT3MEennhd2ynQsvagFASpMzOLFOLb78YgcTpjzC4nen8+P3/4+vv9z1u+/s+3UfLVs3Y+E7U+lxcyeWzHsDgBcmzmH0U0P539sruLzDxSyd/yY9M27k/hF30uDUk0PetyNZtQRoY8z7xphVf1hWG2NWVcf+KsPjSSDXe+BPrWJ/MdHRgWuliR4P3twDZbleL8lJiZxQqxZjRg1n8MABTJk5mz49uzFv4RKGDRnE3PkLQt6H6rZt+5esWLWGNxa+xJsLX+bnX37hnffex+12M2TgAEY/PIy33v0fN11/Lc9Pn8ntGX35/sef2LFrd7ib/pdF+nmxZN6b5OZ4eWHuE7S5pCXf7PqOXhk3cm373lzVtiu7dnxD9/QbfvedjDt6MP3Zl7i2fW9u6X4v4555GIAN6z7jrpuH8c7ry2h6biq7dn7HCbVrMemJqWTc3i0MvTs8t9tV4SWk7aqmeocAiUA3oHPJcmPJv8OqaVoqK1auBmBj5mec1qhRaVnKmf/gkw0byc/PJzsnh6++3kHjRqeWlq9YuZq0lLNITkqioGTc1ufLC20HQiDR4yE2Noa42FiioqI4/rjjyMrKLi3f+/Mv7Ni1m7ObpJGXl0+UOwqXC3wlF8qORJF+XpyZ+nc+XZdJ384Dee/tD/ly2068Xh9eb+C/6Z6ffiY5OfF338nal01Odi4AP+/55Xd/ZQD0vuUmpk+eS1xcLMVFRfj9EO+JD02HKsnlclV4CaVquc3OWvuRMWYWkGqtXVwd+/irLm7bhtUffUzX3jfjx8+IBx9gxpy51K9Xj3ZtWtOl0/X0SL+FYr+fO27NKB1fLSoqYsGSpYwrGbNt2fw8uvRKp+2FrcLZnaB645138Xp9XHfN1Vx3zdX06HcbNWrU4OS6J9HhistKP/f89Jmk9wxkQjd0vIZ+d93D32rXxpzWOFxNr7JIPy927fiGWwf2onv6DWRn5fDvweM4q8kZPDNzDPn5BWRn5fDgPY8B8MzMMdze534mPTGNhx69hxu6dSA6Oorh9z1eWt9JdWuTlJyI3fIlLpeLOifV5ulpo5n4+NRwdbFcwUqMjTFuYBKQBuQDfa2128uUXwY8VLK6HuhvrfX/WX0uv/9Py8KmIGuv8xoVJsW//RbuJjiGu0aNcDfBMc5Luz7cTXCMDV+/V+XwOrX7mArHnN4z7/3T/RljOgJXW2t7GmNaAPdZazuUlCUBq4G21to9xph7gWnW2v/3Z/VF3EVCEZE/CuIQRyvgLQBr7RqgWZmy84FM4HFjzArgx/KCMyhAi4gQ5XZVeDmMZGBfmfUiY8z+oeRaQDtgMHAZcKcx5vTyKtOj3iIS8YJ48S8LSCqz7rbWFpb8vBf42Fr7A4Ax5gOgCbDtzypTBi0iEjwrgcsBSsagyz73/wlwljGmVklW3QL4vLzKlEGLSMQL4t1zi4FLS575cAG9jDEDge3W2leMMfcBb5d8dp619rPyKlOAFpGIF6whDmttMdDvD5u3lil/CXipovUpQItIxHPobKMK0CIiTp0PWgFaRCKeQ+OzArSIiN5JKCLiUA6NzwrQIiJOzaD1oIqIiEMpgxaRiFeBOTbCotwAbYx5H/jTafistRcFvUUiIiHm0BGOw2bQD4eiESIi4eTUMehyA7S1dvn+n40xFwApwFSghbX2g2pum4hIRKvQRUJjzABgJDCQwFR6zxpj7qnOhomIhIrLVfEllCp6F0dP4J9ArrV2L3Au0Lu6GiUiEkpH+lu9i6y1BWXW84CiamiPiEjIOfWt3hUN0MuNMeMAjzHmGuAV4H/V1ywREalogB4EfAFsBLoBbwAagxaRo4JTx6Ar9KCKtbbYGLMY+B4oANaWec+WiMgRzam32VX0Lo7rgQ1AdyAD2GCM+Vd1NkxEJFSO6AwaGAqcY639HsAY04DAOPRb1dUwCTi/WZdwN8Ex1mycF+4mOMZHn1b4rUlSAU591LuiY9C/AT/sX7HW7gQ0xCEiUo0ONxdH95IfvwZeNcbMIBCYOxO4YCgicsRz6hj04YY42pX8O6dkubxkPZfAK8VFRI54Do3Ph52Lo9eflRlj4oPfHBGR0HM5dAy6QhcJjTFXEZiLI5FA5hwFJAAnVF/TRERC44jMoMt4EkgH7gZGAdcAnupqlIhIKDl1DLqid3H8aq19H1gDHGOtHQxosn4ROSo49T7oigZonzHmdGAL0NYYEwPEVF+zRERC50ifLOkBAmPQrxLInH8EllRXo0REQsmpGXRl3knoIvDkYC6wCzi7epsmIhLZ9E5CEYl4LndFBxNCq8LvJBQROVo59CaOCt9mJyJy1HLqgyrOzOtFREQZtIiIhjhERBzKqU8SKkCLSMRzO3QMWgFaRCRIjDFuYBKQBuQDfa212w/xmdeBpdbayeXVp4uEIhLxgvgk4TVAnLW2JTAEePwQnxkJHF+RdilAi0jEC+JcHK0oeVertXYN0KxsoTHmOqAYeLMi7VKAFhFxV2IpXzKwr8x6kTEmGsAYcxZwE/BgZZoVUYqLixk+egxdeqfTK6M/u3Z/87vyBYuX0ql7b7r0Smf5ipUAbLXb6NyjD+n9B+D1+QB4bup0NmzKDHn7q6pGTA1Gjb+fGYsmMHHmY5x8Sl3Ou+BsZi+dxIxFE7j17oNfopOY5GH8CyN5/qUnmL5wAqlN/wFAhxsuY8aiCQwZfkfpZ0eNvx9PYkLI+hMskX5e/NGmzzbT+5bbD9r+xtv/5abeN9Ot7y2MeHQcxcXFeL1e+tw6gK59+rHti8Bw6/oNm5g6c06om/2XBTGDzgKSyqy7rbX7X7DdHagLvAf0BAYaY/5VXmURd5HwvWUfkJ9fwJypz7Mx8zPGjn+KCY+PAWDPnr3MeXk+L8+cSn5BAd379qNl83NZ/OprPHT/YNZ+sp7Va9aSlprCt999T5PUlDD3pvI63ng5Pq+PHh1vp8Gp9Rjy79s5ruaxPHDnaL7evpMp88bT2DRku/269Dtd+1zH2pXreXHaIhqcWo9H/jOULlf148qOl9Lz2jt4/Nl/k5ScSNo5Z/Lpx5nk5njD2MO/JtLPi7KmzprDa2++Q3xc3O+25+Xl8/Szz7PwxRnEx8Vx79CHWf7hKoqKimjb+gKand2ERa+8zuCBdzDn5fk88vDQMPWg8oJ4l91K4CpgnjGmBVD629pae+/+n40xDwM/WGvfKq+ykGXQxpjYUO2rPOs3bqTV+c0BSEs5i8+3bC0ty9z8OU3TUomJiSEpMZH6J9dj2xfbSYhPwJeXh8/nIz4+jmenTCO9V49wdaFKGjZuwMplawHY+dU3nNKoPnbzdo45NonoGtHExsZQVFT8u+/MnrKAhS++BkBUVBQF+QUA5PnyiY2NITo6Gr/fT4fr/8Xil14PbYeCJNLPi7JOrluXJx8dedD2mJgazHr+mdLAXVRURGxMDAnx8SXHIY/4+DjeePtdLmp7IbGxjvhfvkKCmEEvBvKMMasIvInqLmPMQGPM1X+lXUHPoEveX/g08BvwgLX25ZKiN3HAW1hyc70kehJL193uKAoLC4mOjiYnN5fExANv8vIkJJCdk8tNna5j7JNPUbNmTWoefzwJ8fFssZZps+bQqmUL2rVpHY6u/CXbtnxJ64ta8v47K0lpcgYn1qnFl1/s5D9TRvHrL1l8sfUrdny563ffycnOBaBmreMY+eR9jBsxCYApE+fwyFMP8N7bH3JZh4tZOv8temTcSO2/ncCL0xay86tvDtq/U0X6eVHWpRe15dvvvj9ou9vtpmbNwM0HL85bgNfno2Xzc/H7/XywchXzFi3h9n7pPDFhErek92b46LHUq3cSvbt1CXEPKi9YGbS1thjo94fNWw/xuYcrUl91ZNAPAE2B5kCGMWZ/SuGIO8E9ngRyvQf+BC/2FxMdHfg9lejx4M09UJbr9ZKclMgJtWoxZtRwBg8cwJSZs+nTsxvzFi5h2JBBzJ2/IOR9qIql894kNyeX5+Y+zoWXtOSbXd/TM6MT17XvQ4e23di941u6pV9/0Pcam4ZMnjOOp8dOZf1HmwDYsO4zBt78IO++voym56awe+e3nFC7Js88MY3027uFumtVEunnRUUVFxcz7j8TWb12HU88OhKXy4Xb7WbI3Xfy6PAHefOd/3LTDdfx3NQZ3H5LOj/88BM7du06fMXh5tAZ+6sjQBdYa3+21u4FOgC3GWPacWDi/7BqmpbKipWrAdiY+RmnNWpUWpZy5j/4ZMNG8vPzyc7J4auvd9C40aml5StWriYt5SySk5IoKAj8me/z5YW2A1V0Zurf+XTdZ9zc+W7ee/tDvty2A6/Xh88buMi156e9JCcn/e47DRs34LGJD3L/naNYtXztQXX2uqUz0ye/RFxcHEVFxfj9kOCJD0l/giXSz4uKGv7oWAoKCvjPmEcOGqPe+/Mv7Ny1m3OappGXl0+UOwqX6+g9FqFQHRcJdxhjngCGWWuzjTEdgbeBY6thX5V2cds2rP7oY7r2vhk/fkY8+AAz5sylfr16tGvTmi6drqdH+i0U+/3ccWtG6ThaUVERC5YsZdzowNhcy+bn0aVXOm0vbBXO7lTarh3fcMvAnnRPv4HsrBz+PXgcKU3+zsSZj1GQX0B2Vg4P3RO4ODZx5mMM6PMAt9/bh9jYGAY92B8IDHkMvDlwp9Df6tYmKTmRbVu+xOVyUeekE3lq2iNMenxa2Pr4V0T6eVGe199+F5/Xxz/OMCx+5XXObpJK3/4DAOjS6XoubnshAM9Pm0F6r+4AdLruGvoNuJs6tU/EnNY4bG2vKHeUI/7AP4jL7w9uYltyz19XYJ611luyrTZwn7X2zorUUZC11xHZthO0SLsh3E1wjDUb54W7CY7hLy4KdxMcI/bYE6scXTMnvVjhmJNy600hi+ZBz6BL7vmb/odtPwIVCs4iIqHm0MnsIu9BFRGRI0XEPagiInIQh6bQCtAiEvGc+k5CBWgRiXhODdAagxYRcShl0CIS8Rw6BK0ALSLi1CEOBWgRiXh6q7eIiFM5Mz4rQIuIKIMWEXEoBWgREady6A3HCtAiEvGcmkE79PeGiIgogxaRiKf7oEVEHEoBWkTEqTQGLSIilaEMWkQinkMTaAVoERGn3manAC0iEc8V5czRXme2SkRElEGLiGg2u0poltIx3E1wjHWZi8LdBHEglzsq3E04qmgMWkTEofSgioiIQ7nczrwc58xWiYiIMmgREV0kFBFxKI1Bi4g4le7iEBFxpmDdZmeMcQOTgDQgH+hrrd1epvwu4MaS1Testf8urz5dJBQRcbsqvpTvGiDOWtsSGAI8vr/AGHMq0AU4H2gJtDfGpJbbrCp1SkTkKOByuSq8HEYr4C0Aa+0aoFmZst3Av6y1RdbaYqAGkFdeZRriEBEJ3hB0MrCvzHqRMSbaWltorf0N2GOMcQFjgU+ttdvKq0wBWkQiXhAf9c4Cksqsu621hftXjDFxwFQgG7j1cJVpiENEJHhWApcDGGNaAJn7C0oy56XARmtthrW26HCVKYMWEQnefdCLgUuNMasIDJz0MsYMBLYDUUAbINYYc1nJ5++z1q7+s8oUoEUk4gVrLo6Si3/9/rB5a5mf4ypTnwK0iEQ8p043qjFoERGHUgYtIqK5OEREnMmpQxwK0CIS8fRWbxERqRRl0CIiDh3iOOoz6JQmZzDlpfEAnNygLtMXTGD6/AkMHTmwdNyp34AezFk6mZmLJnJW2t8PqqPNxefz4ivPMmvxJK698UoAatc5gZmLJjJt3lOcWLsWAFf836X866qLQtSzqisuLmb46DF06Z1Or4z+7Nr9ze/KFyxeSqfuvenSK53lK1YCsNVuo3OPPqT3H4DX5wPguanT2bAp86D6jyQ6FgdE4rEI4mRJQXVUZ9C9MjpzZcf2+LyBE2bQsP48PW4K69ZsYOiogbRr34rvv/2BZi2a0KVDP+qcdCJPTB7BTVdnlNYRHR3FoAf70/mqDHy+PGYunMiy/62i/ZXtmDZ5Li6Xi/ZXtmPe7KW0veQCBvV/OEy9rbz3ln1Afn4Bc6Y+z8bMzxg7/ikmPD4GgD179jLn5fm8PHMq+QUFdO/bj5bNz2Xxq6/x0P2DWfvJelavWUtaagrffvc9TVJTwtybqtGxOCAij0UkZ9DGmHhjTEwo9lXW7l3fclfG0NL1M1JOZ92aDQB8uOwjWrQ6h6bNUln1wccA/PDdT0RFR3Hc8ceUfqdh4wbs3vEt2Vk5FP5WyKcfb+Lsc1Px5vqIT4gjPiEOnzeP7n1vYM60BaHtYBWt37iRVuc3ByAt5Sw+33LggafMzZ/TNC2VmJgYkhITqX9yPbZ9sZ2E+AR8eXn4fD7i4+N4dso00nv1CFcXgkbH4oBIPBYut6vCSyhVS4A2xjQ0xiwxxkw2xlwCbAG2GGOurI79/Zn/vvkBhYUH5iMp++eJN9dLYpIHT1ICOdm5B7bneElMSixdT0z0kF2mPDfXR1KyhzeW/pfmF5zDuS2asObDdZx8Sl3cLjdDRw2k441XVHPPgiM310ui50Bf3e4oCgsDE2/l5OaSmOgpLfMkJJCdk8tNna5j7rwF/Lovi5rHH09CfDxbrGXEo2N5f/mKkPchWHQsDojIY+FyVXwJoerKoKcBTwKrgQXAeUBT4L5q2l+F+IuLS39O8CQuC2l5AAAHaklEQVSQnZVDbrYXjyfhwPbEwPb9cnJy8SQeKPd44snOysHn9fHQoMd4ePBYuvW9gecnzKLvbV15ZNh4WrdrSXx8pR65DwuPJ4Fcr7d0vdhfTHR0YNQr0ePBm3ugLNfrJTkpkRNq1WLMqOEMHjiAKTNn06dnN+YtXMKwIYOYO//I+guiLB2LAyLyWERYgI621i631s4Allhrf7LWZgGFh/tiddq6eTvNWjQBoFXb5qxfu4lP12Vyfptzcblc1DnpRNwuF7/+cmC+7a+376T+KfVIPiaJ6BrRnNM8jY2fbC4tb3x6Q/Lz8vlm13fExsbg9/uJinJTI6ZGyPtXWU3TUlmxMjCR1sbMzzitUaPSspQz/8EnGzaSn59Pdk4OX329g8aNTi0tX7FyNWkpZ5GclERBQQEAPl+5L4dwNB2LAyLxWETaRUJrjHkBuNla2xPAGDME+KGa9lch40ZO5KFHB1EjpgZfbd/Ju28sp7i4mPVrNzF78SRcbjePDAvc8XF5h0uIT4hn4dxXGTdiIpNnjcPtdrF43hv89OOe0jr73taVUUOfBOCVhW8ze/EkNmdasvZlh6WPlXFx2zas/uhjuva+GT9+Rjz4ADPmzKV+vXq0a9OaLp2up0f6LRT7/dxxawaxsbEAFBUVsWDJUsaNHglAy+bn0aVXOm0vbBXO7lSJjsUBEXksHPqot8vv9we90pI3215lrV1aZltXYJG11vvn3wxIbdAm+I06Qq3LXBTuJog4WkxyzSpH1182r69wzDnuzLNDFs2rJYMumRN16R+2za6OfYmIVJXL5cxHQo7q+6BFRCoiWBP2B5sCtIiIQ8egnflrQ0RElEGLiGg+aBERp1KAFhFxJldUVLibcEgagxYRcShl0CIiGuIQEXEmXSQUEXEqPUkoIuJMoZ6Iv6IUoEVENMQhIuJMGoMWEXEqjUGLiDiUQ8egnflrQ0RElEGLiGgMWkTEoVxuZ87FoQAtIuLQi4TObJWIiCiDFhEJ1pOExhg3MAlIA/KBvtba7WXK04EMoBAYaa19rbz6lEGLiLhcFV/Kdw0QZ61tCQwBHt9fYIypA9wBXAD8ExhtjIktrzIFaBGJeC53VIWXw2gFvAVgrV0DNCtTdh6w0lqbb63dB2wHUsurzJFDHJt2LnfmPS8iclSKSa4ZrJiTDOwrs15kjIm21hYeoiwbOKa8ypRBi4gETxaQVGbdXRKcD1WWBPxaXmUK0CIiwbMSuBzAGNMCyCxTthZobYyJM8YcA5wBfFZeZS6/319dDRURiShl7uJIBVxALwIBe7u19pWSuzhuJpAcP2KtXVhefQrQIiIOpSEOERGHUoAWEXEoR95mF06HexIoEhljmgOPWWvbhrst4WKMqQFMBU4BYgk8BfZKWBsVJsaYKOB5wABFQC9r7ZfhbdXRSRn0wf70SaBIZIy5F3gBiAt3W8KsK7DXWtsauAx4OsztCaerAKy1FwAPAk+EtzlHLwXog5X3JFAk+hLoGO5GOMB8YFiZ9cI/++DRzlq7hMCdCAANgB/D2JyjmgL0wQ75JFC4GhNuJbcB/RbudoSbtTbHWpttjEkCFgBDw92mcLLWFhpjZgATCBwPqQYK0Acr70kgiWDGmJOB94FZ1toXw92ecLPW9gBOB543xnjC3Z6jkQL0wcp7EkgilDGmNvAOMNhaOzXc7QknY0w3Y8x9JateoJjAxUIJsoj9070ci4FLjTGrOPAkkMj9wHHAMGPM/rHoy6y1vjC2KVwWAdOMMR8ANYA7rbV5YW7TUUlPEoqIOJSGOEREHEoBWkTEoRSgRUQcSgFaRMShFKBFRBxKAVpExKEUoCXojDFtjTHLyimfbozpGaz6/mKdDY0xUyr6eZFwUICWSNUAaBTuRoiUR08SSrUxxrQBRgEJwLHAXdbapSXFVxpjbgdigBHW2nkl8wyPBdoCUcB0a+2TldjlFcaYW4HawChr7XPGmLrAlJL9n1RS54PAU8CpxpiJ1tr+Ve6sSDVQBi3V6XYCLzw4G+gLjCxTlgA0B/4J/McYUwdIByj5/HlAB2NM60rsL66kzisI/GIA6AzMtda2AFKAO40xtYA7gHUKzuJkCtBSnboCZ5XMXXE3kFimbIa1ttBa+x2wmkBgvQS42hizAfgIqEcgqFbUUmutH9gM1AKw1o4Ddhlj7gH+QyBj18xrckRQgJbqtIJAJvwJgYzWVaas7BSubgJzTkcB91prm1hrmwAtCLxmqqIKAUqCNADGmMcJZMs7CWTwe/7QDhHH0hi0VJfjCVyIa03g3Y6PEgjA+3U2xiwC6hN4a01fAhft0o0xrxJ479+HQL8qtuNSoJ+1dpUx5gqgbkk7CtH5Lw6nDFqqy88ELs5tBrYQeAlCQpmJ3XMIZNavARnW2j3AZOAL4FNgHTDNWrusiu0YDcwyxnwG3FZSb8OSNh1rjJlVxfpFqo2mGxURcSj9iSdHDGPMWAJDFn+0zlrbN9TtEaluyqBFRBxKY9AiIg6lAC0i4lAK0CIiDqUALSLiUArQIiIO9f8B5/n3UpNYzWkAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 2 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# The predicted label is the one with the largest probability\n",
"labels_hat = array(T).argmax(axis=0)\n",
"\n",
"df = pd.DataFrame(list(zip(labels, labels_hat)), columns=['label', 'label_hat'])\n",
"count_df = df.assign(Q=1).groupby(['label', 'label_hat'])['Q'].count()\n",
"unstacked = count_df.unstack().fillna(0.)\n",
"unstacked = unstacked / unstacked.sum()\n",
"sns.heatmap(unstacked, annot=True, fmt='.1%');"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEFCAYAAADzHRw3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzsnXmYVNWZ/z91a+2luruquwFZFAH7sCMIoghKXDHGPUZHYzLJjKPzJI4+mRkTjf6SmZiZ0UkmEzOTPcYxoyGD0SxqEhdEQTSANEsjHFYBWaSX6uq99t8ftVDdfWuv6u6qOp88eeDu53Tj9773Pe9iCIVCKBQKhaL40UZ7AAqFQqHID0rQFQqFokRQgq5QKBQlghJ0hUKhKBGUoCsUCkWJoARdoVAoSgQl6IoRRwhxjxBiuxDifSHELiHEL4QQZyY5/2UhxOwU9/xnIcRnchjTOiHEJwt1/pBra4UQa7O5VqFIhmm0B6AoL4QQ3wIWAJ+QUh4VQmjAp4F3hBBLpZQfDr1GSvnxVPeVUv6//I+2YDiA80d7EIrSQwm6YsQQQkwG7gGmSCldAFLKIPC0EOI84EHgC0KID4A/A/OBh4DvAJ+UUm4RQnwF+CugG3gLuEFKOVUI8RTQIqX8lhBiAPg34ErgDOBxKeUPhBBVwA+Ac4D6yD1ul1LKJGOeAPwQmAkEgR9KKZ+IOz418tzqoduRa58GGiKnvySlfAT4OVAhhNgGnAc0Ad+NjMkIPCGlfFIIsTKyvxeoBlYAP42MPwi8B9wd+RkqFMrlohhRlgK7o2I+hNeA5XHbLVLKWVLKF6I7hBBXAX8JLCEshPYEz7ECbVLKZcAnge8IIWzA1UCnlPJCKWUTsBn4Yooxfx/YK6WcCVwI/I0QYkaKa6LcBRyUUi4iLMbnCCFqgc8B/VLKcwED8BzwFSnlecAlwD8IIS6I3GMu8BdSyvnAdYA9ct2SyPFpaY5FUQYoC10x0pgT7LcC8XUo1uuc83FgjZSyE0AI8d/AZQnu99vIn1sj966SUj4nhDgohLgXmAGsBN5JMd7LgQcApJRuwgKLECLFZQD8EXg5sj7wGmHRdgshHHHnNAHTgSfj7lkBLAR2A0ellIcj+zcA/yKEWAe8CvynlHJ/OgNRlAfKQleMJO8StlIn6Bz7GLAxbrtH5xw/YYs2SiDJs/oBpJTRl4RBCPG3wM+APuBZ4JdD7qeHn7gXjRBimhCiJu54aMg9LNG/SCk3A2cDPwamApsirqV4jIBbSnlu9P/ABYTdMhD3c5BSHiL8IvpXoAZ4TQhxbYrxK8oIJeiKEUNKeQx4AvilEGJSdL8Q4nPAzcBjKW7xEnBzxG0BYV96JtXlrgKeklL+DJDAtYQFNRmvEXaREHnu64R92FE6AUtcFM5fRA8IIf4NeERK+RvgPmAXYQvfDxiFEIbIOPqFEJ+OXDMFaCHsUhpE5IX0c+AVKeWXgT8Bi9KevaLkUYKuGFGklA8C/wv8VgjRIoTYR9itcWGcayHRtWuBnxCOiNkC1BK2ttPlW8DdQogdhF06WwlbvMn4IjArcs3bwL9KKd+LG5ObsEvmD0KIzUS+DCL8J3CuEKIF2AIcAlYDJ4BNhAXeDlwP/HXkGa8Qfgm8rTOWpwm/gN4XQrwXmf8TOucpyhSDKp+rKBaEEIuBZdEoEyHEl4ClUspbR3dkCsXYQC2KKoqJvcCXhRB/Q9jVcgT4m9EdkkIxdlAWukKhUJQIyoeuUCgUJYISdIVCoSgRRtWH3traPSL+HoejEpcrk2CI4kPNsTQohzlCecyzkHNsbLTr5k+UhYVuMqUKNS5+1BxLg3KYI5THPEdjjmUh6AqFQlEOKEFXKBSKEkEJukKhUJQIStAVCoWiRFCCrlAoFCWCEvQ08Aa8tPa14w14R3soCoVCkRBVyyUJgWCA5/e/xI7WXbg8nTisdcxvnMNNM67BqJV+2JVCoSgulKAn4fn9L7Huww2x7Q6PK7Z9S9N1ozUshUKh0EW5XBLgDXjZ0bpL99jOtl3K/aJQKMYcStAT4PZ04/J06h7rGOjE7eke4REpFApFcpSgJ6DWasdhrdM95rTVUWtN1HBeoVAoRgcl6AmwGC3Mb5yje2xewxwsRovuMYVCoRgt1KJoEm6acQ0Q9pl3DHTitNUxr2FObL9CoVCMJZSgJ8GoGbml6Tqun74Kt6ebWqtdWeYKhWLMogQ9DSxGC42V9aM9DIVCoUiK8qErFApFiaAEXaFQKEoEJegKhUJRIihBVygUihIh6aKoEMIMPAlMBazAo8AR4HtAAPAAn5FSfiSEeAK4CIimUF4vpXQXaNwKhaLM8fkC9PV4qay2YDarYnmQOsrl00C7lPJOIUQ90AwcAu6VUm4TQtwNfBn4ErAIuEpK2VbQESsUirImGAyyce0BDu1to6fLQ3WNlbObGlh26XQ0rbydDqkEfQ3wXNy2H7hNSnki7voBIYQGnAP8WAgxHviZlPLJvI9WoVCUPRvXHmDnlmOx7Z4uT2x7+eXnjNawxgSGUCiU8iQhhB34HfATKeWzkX3LgJ8BFwMDwH3AfwBG4A3g81LKHcnu6/cHQiaT+lRKhcfvxTXgxmGrxWpSiU2K8sXn9fP9x9fhdvUPO1bnqOBvH1iJ2VIW6TUGvZ0pZy6EmAK8AHw/TsxvBb4KXCOlbBVCGIHvSin7IsfXAguApILucvVlNINsaWy009pafNURM2mwUaxzzAQ1x9Ih23m6Xf26Yg7g7uzn8Acd1Doqch1eXijk77KxUb84YKpF0fHAK8AXpZSvR/Z9GrgbWCml7Iic2gSsFkIsIhw5sxz4n/wMvXxRDTYUisFUVluorrHS0+UZdqzabqWyury/YFOtIDwEOIBHhBDrhBDrCUe42IHnI/v+SUq5G3gGeBd4E3haSqnfHUKRFqrBhkIxHLPZyNlNDbrHpjY1lH20S1ILXUp5H2HfeEqklI8Dj+djUIr0Gmyo+jKKcmTZpdMB+GBvGz3dHqrtVqZGolzKnbJYPShGog02OjyuYcdUgw1FOaNpGssvP4ell0xTcehDKO+gzTGMarChUCTHbDZS66hQYh6HstDHMKrBhkKhyAQl6GMY1WBDoVBkghL0IkA12FAoFOmgfOh5wBvw0trXnlUoYS7XKhQKRTzKQs+BTDI583mtQqFQ6KEEPQcSZXIGQwFuFTdmdS2oLFCFQpEdyuWSJckyOTcc+zOr97xAIBjI+FqVBaooVny+AG5XPz6f/r97ReFRFroO3oA3ZVRJskzOIEHWH38nFqWSybUqC1RRbKj65GMHJehxZOLXTpbJGWVn2y6un74Ki9Ey6CWhskAVpUSp1icvxo5IStDjyMSvHc3kjD9/KB0DnXQMuFl/7J1hL4m5DbN569jbw64phizQoMeD3+3GVFuLZrUW/DrF2MXnC3Bor36Tsg/2trH0kmlFI4ZRivmLQwl6hFR+7ailHc9NM64hGAqw4difCRIcdp3TVsfao2/y9vFNsX3Rl8Qlk5axcvLyQVmgc+pnsmLShXgD3jEp6qFAgNY1q+lp3oq/owOT00n1wkU03nIbBmPi/2izvU4x9unr8eqWsgXo6fbQ1+MdM/XJ06WYvziUoEfIxq9t1IzcKm4kFIL1x98Zdp3FaOWd41t077mzbTePXPD3XD99FR0Dnaw7+jYtbXtYf+zdhK6e0bZwW9espvO1V2Pb/vb22Pa42+7I+3WKsc9Yrk+ejcuk2L84ylrQ8+XXvqXpOoyakZ1tu2gfcFFnqaHCXMmJ3pMJr+nwuHhm96+5+uzLWHf07UEvhKGunkJbuOksAgc9Hnqat+oe62lupuHGT+q+ZJJd596wnvrrb8RYUZn94BWjSrQ+ebxFG2Uk6pPriXYuLpNi/+IoS0FPtPg5yyl4+8S7w84f6tceKoBGzRhzv2xv24Xb00WXtyflOLacambLqWa0BNGjUVdP55o1BbFwM1kE9rvd+Ds6dO/jd3XEvhyGfkEkuy40MMCpXz7DGZ+/K+s5ZMJof+GUKiNVnzxevI1GQ0LRzsVlMpa/ONKhLAU90eKnVQv/sjQMBAnhtNYxv3FurLphMgF8fv9LvHXstJUd0vGpJ0LP/w5hV09nd3tWlnE6ZLIIbKqtxeR04m9vH3YfY50D1yt/pHfn9tgXRNW8BdRddgXG6mpMDif+juHXAfTt2UPQ4ymowCoffmEpdH1yPYvbajPRfqo3dk5UtAf6vBw/6ta9z6E0XCaj/cWRK2Un6MkWPz3BcEJPkBAAs+tnxoTNG/CyWr7An0++Fzs/KoCBYIBd7XvyPlanrY6Kfh+u9nbdFt9Ry9gyblzG9850EVizWqleuGjQl0IUY1UV7nVrT4+rvR33urW4163FVF+PIYlYBzpdWc8hXZQPf2SI1ifPN3oWdyK3yL73WxPep6crPZdJMXdEKjtBT7b4OZTNHzVzw/SrefHQq2xvbUl43c62Xbi9+e/uPa9hDn/q2MJZVRq1vcOteJPDiam2NuP7Bj0eOk5+QHevS/dfQPwicLybovGW24Dwl4Hf1YHJ4aRq/nx6d2xP+KyYRW80QmB4BqHJ4USrqMB76tQgV0i+3CPZ+v4VY4Nki5SZYjCAxZrawi7mjkhJBV0IYQaeBKYCVuBR4AjhRtEBwAN8Rkr5kRDiLuBuwA88KqV8sYDjzpp0EoKieAIefrX3N2z+qDnpeZ3eLmqtNbg9XTmNTUMjSJB6m4N5DXP4xNlX8M0/f4fAJCsL9/YPO79iwfyMxCjQ18ep1c/Sv2c3flcHd1YZ2T/RzPpF1YS0098ATlsdNaZKTq1+RtdN0XDjJwf5zN3r3kj5bIPZTEhH0LXKSo48+vXTrpoF52LAQM/25mHPzYZ0fP+F/Dood3JNzkm2SJkpoRB4PQHSXYMv1BdHIUlloX8aaJdS3imEqAeagUPAvVLKbUKIu4EvCyEeB/4OWAzYgA1CiFellPn5TeSRdBKC4pGu/SnPqbc5mOk4h7dPbEp6ntVopdJow+XV9/Etn7SUS6dcTK3VjtFg5Jk9z+HydrJ+UTUA0455qO4N0lOlcXCSjY9duyqtOUR9yO4N6wkNDMT223v8LNzrB+CtxaejeOY1zMH9/PNJ3RRREUzmWx80Bo8X+7KL6JcyZt1rlZV4jx4Z9Az32tcHXTfouffek9Z840k2vmy/cBSpyVdyTrJFykyptlvSWtQsxgzRKKkEfQ3wXNy2H7hNSnki7voB4Hzg7YiAe4QQ+4H5wOZkN3c4KjGZRuYH1th4WrDurr+Nyu1mtny4g9a+9ojHfDhWo4WuNFwpS888lyunX5xS0D0BDxZj+Edu1SxomgGP30tDpZPFk+dz54KbMWpGPH4vP9r8vzF/fUgz8NZiOxvPraaqP0BvhRFnTQN3nTkFq+n0P9D4OcZz8CdP6vq+o8w45uXdc6Gupp7Fk+dz69lXsPOnD+ie27d9G86/+RzG2JeBnd4LL+DEiy8lnbu1sZ45938BAG+HC2NVJTv+Xv8ZevTv2E7A40k4x8QkHl/jhUsZP7khw/sVnsznOPb4429adCNNKiosrLphLpD+PGcvmMim9YeG7R8/sYaBPi/uzgGdq3Tuc+4kJk6sS3g8GAjyyu/fR7acxN3ZT21dBWLuBK68djaaMbsM0ZH+XSYVdCllD4AQwk5Y2B+OirkQYhnwReBi4Cog3uzsBlKaPi5XX3ajzpDGRjutrYOF+ZrJV3PFGZfh9nTz+pE3WX98eLjiknGL2Nn+Pm5vYlfK0gnnsWrilQT6Alg1S2xhNRHd3vDKvCfohSAsGbeQ22fdjMVoobW1m+f2vcifT2zWvY/fZMBtD//KZjtn0eXyEPZ66c8Rwj7k1neGzy0ee2+Av9szjvGf/Tw9f/gD2//zHwh06q8XeFtb2fEv/07DtTdgmTAhvFh67U1Ytu8cZG0Pu66rmz0//jmOVdfgPXECzWbF05q+b9TT1sbAyY/o7PFl7FevvvYm6ga8g3z/1QsXUn3tTbo/s9Ek0e+xmPD5Ary//bjusd3bjzP//MlMnFiX9jwXXjiF/n6v7iJlIBDirT/tY2/L8JwPs0XD7wvGzl944ZSkz9zw2r5BLyG3q59N6w/R3+/NKkO0kL/LRC+KlIuiQogpwAvA96WUz0b23Qp8FbhGStkqhOgC4p9gB9JbeRxFoq3dbmm6HqNmYkdrCx2ezkHhioa9Bt0sUACn1cFt4kaMmpFAIEAooa2fmC2ntrFkwrmcVTOFP3ywVre+y1CWTjgv7UbRyXzI8fh2tPDhA/+ou3A5lL5tzRzZ1ozBZqN22XLqb7iJYF9v0mtCAwN0vvbq6S8FgwE0La3nARgsFt7/xjfxtrWnFXY4dFF13G13DPL9q4XQwpFOck4maJrG0kumMWvBGRCCGkdFzBWiafCxjzdhtRmHCf6SFWcz0OdLy3VS7BmiUVItio4HXgG+KKV8PbLv04QXP1dKKaNKsQn4phDCRnjxdBbQUrBR55lkzZhvabqOg12HOdYz3OKY3xhOOAoEA6yWL+AN+jJ+dogQ39/xcwyR/6Ui/iWSDun6uIG0xTVKaGCAzrWvERgYSOulMfjiUEbPCw0M4I34/5OFHSaLOVeMDPlMzknHF58sKsVqTS+Qr9gzRKOkmu1DgAN4RAjxCGAE5gKHgeeFEABvSim/JoR4AlhPuGnGV6WU6Tm2xhB6zZiNmpEvL76X5/b9jh1tu3B7unHa6sJRKNOuoLWvndePvDUoPj0bQpH/JcPkD7G46ixM/lD4N6FDNIvVbrBAuxsDUDV/Ae431upfkAf6du9Cq6sj6EodOaSHyeHE3+nCWFuHrUngO/Yh3mMfpryup3nroLDDoMfDR888TffG0185UfHvk5JgX69KLBoB8pmck0nWZy5RKcWeIRrFEApl7ibIF62t3SPy8Hz5sqJiWW2u5MVDr8YyRrNxtWSCIRhixdZeZp8MYe0e0BUkZ30lP3pnNTs/amHOxqPMPjSAxR/CABisNowOB/6TJ5I/KAdMjY34WxMndSTD2NBAoKsLvJl3arIvu4jxn/4sbc+voWfr1oQZqXrUXX7FmEssKgUfOpy2rPX83pqmpTVPny/A6p9s0hVZe42VW+86P+nLIdNolaE+9CjzFk8aiz503c/5skssyoWoBb9m7+/SDnvMByu29gyKQ9dzOfxi+69Z9+EGLt7SzcJ9g2PWQ56BsJhbrOAtTCSpv7U1YfJQKgJt2SeOdG98G8/Ro0kXZBOhEosKRz6Sc7J1g2QTMunzBZizcCLBQJAjBzpyzhD1+QJ0tPXi8wVG1PeuBJ3E1Qb19idLmS8EVUEzsxMUbYwKkt9kYOsH23B2+pl2dHgC0mkK/CWRIHmo0KTjntFDJRYVnlzdIGaLhs+rkyVt1hK6QTJx0+iJ/1kz6pl33mSqa6wZi/Gg+0VeCiPZHKOsBT1Rsa3rp63itwf/qFuEK5PSAfng3um309/9Ld1jflcHvo4OTr72Mqu27MPeG0y+rOr1YrBaCXkKY6XHJy2NKMH0C6HFoxKLRp6oG6SuNrXI+30BQhn+ajONVtET/11bj6NphqzcLKPdHKMsBT1qeb9+5C3dOuT7Og8OimqJ7g+GQoRCwYL7zKPU2xw0TjiTY04ngQSZjp2vv4rnzfWpg/4BrbaOylmz6Xl3Y/4HWwg0DUIhTE4nFTOa6N8v8bfrRNNoWkJR1xxOjDYrvhPD1w+qFy5U7pY8kI6veqglXOuo4MzpTl3LNXrugT2t+P36v1efN6jrcsnETZPvUMWxEPpYVoI+1CJPFCaoF6IIsOnkFgYCI1fNYE79bH5/9A0s4wPM0lnrS1UYayhBd2fxiDlgv2g5BAL07d5N96Z30Wz64muZNFnfh65phPp68XW6MNhsAIQ8HkzOeqoXLlShjDmSia96qOXqdvUntFyHnquHvcaKxWrE7eof9CLJJFol36GKYyH0sawEfWj970wt7ZESc4tmYemE8zAQYt2Hb2OYb2bAXxGr5eKvqWT8+RdRu+KSgoYjjjaeQwfxfnjaPx7sj7h0rFbw+WIZnw033ULb82uG1aohGIy5l6L77cuWM/6OO5VlngfSdS8ks1z37DjBkhVTsVrNKc+Nx2Iz8dxT7w17kWQSMpnvUMWxEPo4tltY55GRXszMBqtmwaJZ8AV9tLTt5t0htVz+95p6fnGtk+dvmIT95pvY8pufjfKIC0u8mA/C46Fy/gKmPPgw4267A81iCUerVKYuo9cv81+3vhxJ5V7w+U4vjiezXH3eIBte3Z/WuQBV1Vbqx1XRfqo3dl70RbJx7QGCwSChUAiT5bS0mS0ac8+bOCxaJSr+emTTzCLf98uGshH0kV7MzAQDBs5tmIsn6MUb9BIihMvbiWfIF0G0lkurv4tfv/8Clr2HR2nEo0/ftmY+ePgrnFr9DKFAAL/bTSCNxKZoZIsiNzJJ749arok4ftgVewEkO7eyyszl181koF8/I/uDvW1seHUfLe8dxx8XGePzBjEYDLpRJssunc68xZOw11gxGMKunHmLJ2XdzCLf98uUsnG5ZFIHHcKlbis0G50+N3ZzNR6/F28o88SXdHBYazncnX7oncNax4fH9rJYp+lFORGtDwPQcOMn0ypxoCJb8kMm7gWz2cjEMx26BbQAenu8Mf9yMpeJzxfgt88mXjPq7vJwMMNFyXw3s4i/n81iZsDrG9E49LKx0KN10NPFiJFOX9iS6/b1FEzMAc5xTKfTk77VOMNxNidMvXRXlc2vLyk9zeEGJNULF6U8V0W25IdM3QvLr5gxyA0Sz9AXwFAr1xy5Ti8efSj9vfrWe6qiYNF4+XyJr9lsxNlQNeIFvcrGQgdiFQp3tu2iYyAcXz7DcTZWzcL7HZKOgU6ctjoG/B56/YUv7asZNFZMvIBrp13FPtdB3a8Hm9FKpakSlyc8trn1s5CuA/hNBg4m6GRUbvjb2/C5OmJRK91btxIYUgIgWhlSRbbkj0x6b1qtJmbNPyOtxcp4K7fL1c/Lz+3El2aGs8EQrvs2lHwuSo7lBhhlJejRqoqfOPsK1uz7HXs7DrD5ZDMOax0zHeewcPx8nNY6vrFJP5EnEZpBI5hpBgQQDAUxGDQqzBUJuyjNqRd8qukG+v0eaq12Xtj/Eif7PgKIdDIKMfvgABb/6WtS12wsPTpff5Xxd3yGcbfdQdDjpWv9m4OOhwYGQDOoYlx5JFN3xbJLpxMKhZA7T8asbZNFIxQKEQwGh/m4zWYjJrMxo25FiUpT5WNRMtOSAj6vf1hYZaEpK0GP8uKhVwdVR+zwuNh4chMbT26i2lyd1j0cllrmNszmoonn86Md/4PLm92C647WFq6fvorrp61iX+dBTvScIBgXTvneqR00n2ph+cSlXHP2Fbx74vS4w31ADVj9OjcuM3p37CBwYx9tzz9H19vrdc9RtVsKQ7rp/ZqmYTAYBrlO/N4gLe8dx2DQz8zMtAWdOfKC8PtCsW0xb0JeFiXTDdOMCv+RAx24Xf1Zt9/LhrJzwqYKX+zx9aR1n7NqpnDx5Atwe7rpTNAjNB06PJ2sli/w/P4XOdZzfJCYRwkS5K3j7/AfzT/EG9fJyOQPMe3YmGvbOir4XR2cWv0s7nVrE2aNqgiX0SWTUMcoZrORs2bU61yR4BneYEzMo9uJIlwyIZOxR4Xf7Qq7Q+PDKgtN2Vno+Qpf3NbWwra2/PTwSLeW+kd9pwZtV/UHsJd5pEsMs5n+Pe8nPUVFuIwumWZSxkrw7gsLadQ/Xl1jZeo5YZE/vK+dnm4PVXYrAwO+QeGKUQ7uaeW8ZWdRUZm9Dz3dsY92+n/ZCXqm4Ytjmd4KI91VGrV6om61QoGKcI1JPB78KearIlxGl0wzKYe6OKL+8bOmO1lxRRMAF6wML1D6fQH+78ktus/t7fGy5sktTJvZmLXbI92xj3b6f9m5XDINXxzLRCNd9DDUO0d4NGMAQ4J/zppG7ccuVREuo0wmoY7JLN0jBzpiLo6o/7661pYwLBLCop6L2yPdsSdLjBqJ9P+yE3QIhy+umHghWpFNv8JoG7Zv/aJqmpsqcFdpBAB3lcaembVoA4WLmx+zJIg0qr14JePv+IyKcMkzPl8At6tf1/ediAtWTqN+XFXs3WswQP24Ki5YOW3QeZk2mt68/pCuu2UoiXz16ZBOFuhop/+nahJtBp4EphJu/vyolPJ3kWPfAaSU8oeR7SeAi4Boz6XrpZRjcgXKqBm57MyL2XD83dEeSkb0B4bXG4/Wedl4bjVV/QF6K4xcWruQwE9fHYURji1M9Q2qqmIByKYjUJSNr++n/VRvbDsUgvZTvby77uCgSJFM3DPpFvSC3Nwe6YZpRgX+6IEO3J39g+LzCx3DnsqH/mmgXUp5pxCiHmgWQrwDPA00Af8ed+4i4CopZfb9xEaQWqudOmsNrgwyNMciRjRqrbX0+F3UWWtZNGkeH5+2iqPO91KmwZcyxjoHkx94EAIBQn6/ss7zSDZNHILBIBte28/72/T72g5dMMykamKqgl7x5MPtkSpMMyr8dTdXcPiDDiqrLRiNhqxfgpmQStDXAM/FbfuBauDrwNXRnUIIDTgH+LEQYjzwMynlk6ke7nBUYjKNzH9ojY322N8DwQC/2P4H+vzFn2UZDAa4+8NxeLa2Emg/hKWxm4GlARrOX8zJP/xptIc3agQ6XRx7/Jv4OlxYGxtwnn8+Z3/+s0Uh7PH/VscaPq+fIwd0mowQtkjrbq7AbBkuK3/8TQu7tur3GYCw5WyzmHE2VOHz+unu8vDxG+ZRUWFhb8tJ3J391NZV0DR3AldeOxvNeFoE62orqHVUxMIEkzFrwUQmTqxLY6ZhomOx11h155WKGU3h9oZ//E2L7kuwosLCqhvmZnzfRCQdoZSyB0AIYScs7A9LKQ8Bh4QQV8edWgV8D/gPwAi8IYTYIqXckez+Llf26fWJ+oDqMbT79kg3eS4kK97rpm/f6UQab2sbJ158CYOK5sAX6W7kOdXKiRdfon/AG2uqnYigx4Pf7cZUWzsqETGF7BSfKz5fgI+OdSUUTndnP4c/6Bhmvfp8Ad7fnljMIWw59w14efOXcpAVO/FMB9fdcS5+XzDmpmjv6B3qf/A1AAAgAElEQVR2/ZnTnbrWvNmi4fcFY26PhRdOSevnm4tbKUr0d5ls/ru3H2f++ZMzdr8keumnfOUIIaYALwDfl1I+m+C0PuC7Usq+yDVrgQVAUkHPhkR9QG+acQ1GLfUPpRjqoqeDyR+iptvP7IP6fTwL1Te0mEmWKRoKBGhds5qe5q34OzowOZ1UL1xE4y23FYVVX0iGilum9VLScYlYbCY2vXWQlvdOC19Pl4e9LSc5uPcUs+afkTTbM1FdmSUrpjLQ58/YZ53P3qAjGcqYalF0PPAK8EUp5etJTm0CVgshFhGOnFkO/E9eRjiEoV2Hov0+AW5pum7QuVErvsZ/+j/gsVwXPRUmf4jq3gALZB9nn/CmbgqtGEQ0U9QybtywY61rVsdK8QL429tj26ms+lInUTz4UKK+7aELf+mk77ef6qWrU9/y93uDKcU02YJltBtSuuQ7OWgkOxmlstAfAhzAI0KIRyL7rpZSDvrJSyl3CyGeAd4FfMDTUsq8m8HJrOudbbu4fvoqLEbLMCu+vsLB9NqzueWc64oyscgQDLFiaw/TjnmoUSKeNYkyRYMeDz3NW3WvKff6L8nEzWCAEGCPWMMXrJzGhtf26bopEi1wDnpWirDDdMQ03boyyUjHoq6stqQdrZLJAm+upPKh3wfcl+DY14dsPw48nreR6ZDMuu4Y6MTt6aaxsn6YFd/W30FbfwfbWlu48IwlzG2YzVvH3i7kUPPKiq09qkxuHkiUKep3u/F36C/0JbPqy4Fk4hYKwbW3LWD8pBrMZiMbXtuX0E0RdYkc2NOatC55Mkaq0XIyi7qq2sL2zUc5vL89I996JqWGc6GoUv+TWddOWx21VntSK94T8LDuww1cMmkZKycvZ+upHXR5uwo97JxQBbgywGDA5KynasECDBjo2b4Nv6sj1kw6UTy6qbY2Ybejcq//kkzc7DXWmJinclMsWTE1spW4MbvZoiW10keq0XIyi9paYR4UrZOubz3fnZESUVSCHk3b14tQmdcwB4vRQmtfe0ofeUv7bh5e+ve4Pd00tyZuaTUWKLYCXJO/8lVO/s/P8Z9IHtWQb0xOJxP/7ktYGhtjVnjDzbekFbGiWa1UL1w0yIcepdzrv6TrLkjlptjw6n72tnyU9Fli3gQMBgN7dpzQFfZM3RO5JPHoWdRnTndyOEHIZrq+9Xy4hJJRVIIO4bT9YCjEppNbGIg0UbYarYRCQfp9/fiCXuostUnrk3cMdNLW384H7rHfZDlpAa6xhtGIbcqZg2KEddG0hCVus8Xf20vXhjcHWeGa1Zq2qyR6XU9zc1pWfTmRjrsgqZvCbuXY4cRrVkPdFktWTGXDq/s5fthFb483Y/dEPkIO9Szqvh4vu5r1E6NGyh2UiqITdKNmRDMYYmIOYVfKm8c28u7J9/AGvCnj0p22OkIhQ051zEeKZK3mgoDPxNhpcBEI4jl+DO/xBNa5wcAZ9/09rU//HH9H9lmsVYvOo3dXy+Bqkh4Pna+9SigYZPztd2Z8T4PRyLjb7qDhxk+Oahz6WCQdd0EySz5Zg2iAj39yHvXjTjeWsVrNXPaJWVlb2PkMOYy3qEcyWiVbiqs6FckjXTwBDyFCeCJib0wwvXkNc3BW1GLRRv8XkA56Bbi2nWPjF9c4+NkNDeFjlYYk3smRwVhXh6/Dldj6DoUIdnfhd+l/tkYx1NWFQyj00DQab7094bVdG98mmEMMftSqV2I+nFSNlBMVr1p+xYwkFQgtYCBhc4tMGzdn00QjXUa78FY6FJ2FnkkceY2lhul1Z/NBz2Ha+1w4bXXMawgnIT2//yU8weJYbNQrwOU3nRa8txbbaZlRwR0vd4xqSGOg00Xb6mcSn6BpVDQJDFZruMdnAmoWLqJvz258J4Z/3lomTSbY35+w1ntoYABvayu2yZMzHr8iN5JZ8omsd4/Hz//9bEvMLbJkxdkM9PkyssrjLflCJ/GMVLRKthSdoGcSR97pdfOJaVcyfdJEDhw7HisTMNazRU0GE/7QcD+K32TAbdf/lXVVJ/a1h8h/42itqhqD1UpgiOskmfVtmTQZU3Xynq2WSZPp3bE9HEZoNEIgYlFpGpVnnsnEBx7Cd+pU0nuoOP3RRW/hb6gQmi1GvJ5AbPEz6hbZveMEfm8wLb+3nq/8rOnOpCGHubpFRipaJVuKzuWSSYOKaCij1WShsbI+5lt3e7rHdGKRnpinvCZZs4vonzZbLLTPMuXMcFejLNFsNqZ8+SGMDof+CUbjabeJpmGZciZnPvgwfrc7qXXuPfZhOHwwFIqJedWi85j27e+y8LvfRrNYMDc2hueig8Fmw9zYmPW8FIUhKoS33nU+N3/2vITn+YcIfLKGFFFfeVS8e7o87Go+gdWmb/R4PH7+/OZBgnlYkM/GHTQSFJ2FDuFIFwhnh3YMdGLRLLruk7n1s3QXSGutdiyaZVDD5VJg/aJqDKEQ8/YPYNRxqGuVlUz6ysOx0L6Bw4c58o2vZfUsv6sD36lTBDoTuL+CQSb941cgEMQ6eTIme7iYkKm2FlN9vX5p3wTRL57Dh9Esp3+PmtVK7bLldK59bdi5tcuWK//3GMZsNrJjyzG8nvR82YnCAZP5yj39PuYsmsjelpODwh99aZQQKHaKzkKHcKTLLU3X8fDSv+drFzzAkvELdc8LJfn4NpTgh3lIM7BtZiWGBKujAZcLzWKJCZ7J6QyLaCLM5oSWsMnhDAu1M0Gru1CIkz/9Eb3bt2KsrIztjsZ865LAcopma8bTeOtfUHf5FZjq68NfHfX11F1+BY23/kXi+ShGHZ8vwPEj6X8d63UnguRx7709XuacOxGrTb+GS66Lo2OZorTQo1iMFmqtdt7vkLrHd7W/jzdw9bD9bk93SuvcgIHQqMeNZE6yuPUBuw2DPezDjlYXTBYPbqyqovrcRbjXrR12rHrhQkx2e8KEHIBAR4dugSu9mO+q+fPDvvM0szVVmGFxkkkzCkgcDpgqhBADo9qsebQoakGH9Oq7TKJ+0P5kC6uaQaNCs9EbyL5W+2hQX+HgU+fcSPOpHRyc9IZu3PruCSCPvsYtTdfRumY13RuT17MJuN3UXXYFBpMxYcLNaXHemrBD0tACV4nE+JTxmYyzNTNJHlKMPulUXownUThgsrj3s86pp2XrhwnvOVZixgtB0Qt6OvVdosQ3xUhUQiAYChadmAMsmbSA3R17kR376VwUtsKnHfNQ3Rukp0rj4CQr6xdV42zbxbVTPpawumA8Jmc9ZqczqSUcFeea5Zdw5OsP694nUYGreDEOBQIQDGGw2WKLpgabjZplF6lszRIimRDXj6vCO+BPOxwwUQhhKBTi/ebEiUxjJWa8EBS9oKdT3yUQDLBm7+8GNcWY2zCbSyYto6V9Nx0DndjN1fT5+7OKMBlNnFZHOOrHEDr9M0gSt94x0EnnqeMJqwvGE28Zp7KELY2NCRc70ylw1bpm9bBFztDAAAZNK/sGE6XGskunU1FhYff243R3e6iqsjK1qZ7ll88gEAjR5eoHA9TUVSRN1dcLIQRY/ZNNCa8xmTWWrDg773MaKxS9oMPwqJf4BCJvwMsPNz/Pmx++Gzu/w+PirWNvc/Gki5hZdw4bT26myzc2234lY0HDXG6feRNd3m5+tHN4PxG9uHWnrY66cRPpS1BdEMDorMe+aFFGlnEuBa5UPfLyQtM0rrx2Nr09Axza105vj4cP9rYRCoXX6D/Yl1lp2vi4d7erP6k7x+8LMtDnw2otCekbRknMKhr1cv30VTGXitFg5Pn9L7G9tSWhjz2+wFcxYrdU89jm72UUUz+vYQ62isSLmfZlyxl/x51ZCWi2Ba5UPfLy45Xfvz+o0FVvj5f3mwfXAMqmBksqH321PffkorFMSQh6FIsxnEAE6TWCLmYx19DYcPzdtM4LERr01QLJxTdbF0e2kSeqHnl54fMFkEmKdQ0lk7ZvyXz0AGeLxpL1n0OJCXqUsZ7anw+CpJftdtHEpVx25sWxsgdRChn2l2nkiapHXj74fAE+OtaF25V+B65MwwyXXTqdUCiE3Hk6schs0RDzJoyZmiuFIlWTaDPwJDAVsAKPSil/Fzn2HUBKKX8Y2b4LuBvwR857sYDjTkq6BbysmgVPiWWLRqm3OWIWuVFLbJGMlbA/VY+8tBlad8VgSNxseiiZhhlqmsaKK5q4YOX0cOPpENSMwTT9QpDKQv800C6lvFMIUQ80CyHeAZ4GmoB/BxBCTAD+DlgM2IANQohXpZSj4tNIVcDLaa1jfuNcgsEAbx1/Z4RHlx80g0YwNNxKd1od3DP/czRWOlPWhR9LqESh0SWX7j7pMLRGebpiDtmHGZrNRuobkxeDKzVSCfoa4Lm4bT9QDXwdiE/BPB94OyLgHiHEfmA+sDnZzR2OSkymwrw1LzjrXF7e+8aw/edPOpe/WXwHNbZqnty6uiDPHgn0xBzC8z53WjHXqbDDZP2a0wCNjfaEx0qFkZxjMBDkld+/j2w5ibuzn9q6CsTcCVx57ezUnafSxOf1cyRB6zY9LFYTPq+f2roKmvI8lpFmpP+9JhV0KWUPgBDCTljYH5ZSHgIOCSHiBb0GiC+20Q2kXMlyuQqXwLNq4pX09fnY2bYL10AnZqMFQrDp2Db2tx1mbsNMdrbuLtjzRwKb0UqlqRKXp5PGSieznbNYNfFKWluLLwQzHRob7SU7tygjPccNr+0bZDm7Xf1sWn+I/n5v3gpYuV39SX3mFZVm+vt82GvCiUFLVkxloM8f+1po7+jNyzhGmkL+LhO9KFIuigohpgAvAN+XUj6b4LQuIP4JdiC9LhQFIj6U8TeHX+TND4bGoafvatEwYNRM+IK+Qgw1azwBL19a9AUsRjPTJ02ky1W8UTuKkSdVd590I0tSkSyU0F5j5ea/PA+vJzDI3WO16hfWUiQn1aLoeOAV4ItSyteTnLoJ+KYQwkZ48XQW0JK3UebIrlN7dfenW4ArSIjgGBNzCCcJRX3lVpMFUIKuSJ9Cd/eJkiyUcGpTAxWVFioqdS5UZEwqC/0hwAE8IoR4JLLvainloO8nKeVJIcQTwHrCJXm/KqVM3MVgBHF7umnr0/ffFWM1xXjm1M+MJVIpFJkykk2Ph9Zdqa2rYMp0Z8mHEY40hlAmy815prW1u+AP7/H28OCGRxPGbRdT6GI0SchhraPCXEGft49OrxuHtY4LzjqXVROvTBqiWOwoH3r+GepDjzJv8aSCNIGIRtOcNdVJpzv9WPRipMA+dN2GDsW5dJwB/X5P0iQch61uBEeTG0GC3HvuXcypn8mxnuO4vJ2ECNHhcfHy3jd4fv9Loz1ERZGx7NLpzFs8CXuNFYMh7NOet3hSwSznWOs2S0nmNI46Jf9TrbXaqbc5aB/Qj0nv6O8sGitdw8CWk9vY49qne3xn2y6un76qqOLPFaPLWG96rMiMkhd0o8GIzWyFBB59b8hLsbjSg4TYeDJxadBoQ49aqz32pxJ3RTrEVyxUFC8lL+jP73+JY93pFwIaDZZOOI99roNpV03U0HTdSA5rLWuPvkVL255Y3ff5janT/xUKRWlQ0j70YinStfmjZmY7m9I+P9GaQIW5kreOvUOHxxXzra/7cIPyrSsUZUJJC3q6RbpGm2AoSJ/fw4qJF6Kl8StxWuu4eNKF1NscGDBQb3Nw1YyL6ffpRw3sbNuFNzD21wgUCkVulLTLJVWRrnQwMDIu9gOdB/j83DtYn0axsPmNc7ml6bpBPVKN1UFe2b9e9/yobz1aK16hUJQmJW2hR/uN5oLdPDJJO25fN99t/nHScwwYmFg1gYsmLsUb8MYaeliMFhy2WhxW/RDMoc2y9fAGvLT2tStLXqEoYkraQodwv9HKSjN/PrIt1m/UZqrgWM/x1BcDXb5uLJoFr05YY6L92ZKqaUWIEMd7T/LNTd+ONYeOLnhaTambZesRCAZ4fv9Lgxpoq4VUhaI4KXlBN2pG/nLhp7jijMuG9RuNNpV2WGvp9ffjSdCSLpFoL52wiE0nm/EER76GSnTBE+CWpuuA5M2yE/H8/pcGvQT07qtQKIqDkhf0KPH9RoFhTaV/e+CPSXuQWjQLVeZKOj3uIUJpSMvvHY9ZM2MAvEEfddYaOj1dWc7qdDIR6DfLThaHniwKSCUpKRTFR0n70FMR74O+acY1LJ1wXsJzvUEvgVCA8ycs4sEl93NL03UYNSMrp1yU8XN9QR/eoI+lE87jwSX347Q6sp5DdMEznvh5JSNZFJDefRUKxdimrAU9HqNm5DZxIw5L4touXd5u/nzyPV48dLqZsdNWh9WYXau0/Z0Hc164dVhrs662GI0C0iOdhVSFQjG2UIIeh8VoYcG4uSnP23ZqBz3entM7sqxYGbWCb5pxDSsnL4/FlZsN6Rf3P8cxPWu3SLKXSbKFVIVCMTYpGx/6UOJjuC1GS2z7E2dfAYRFu9Or79vu9HbxL5v/k3Mb5tLnH8i6sFfUCh7q+642V/LioVfZ2baL9gFXwkYcVqOVW87JbeEym4VUhUIxNik7QR8apldnqaXSUkm/r39Q2N4DS/6Or7/zeMIIF7enizePbcxpLEOt4PiF23iBX3v0Ld2WeReesYQKc24FlTJdSFUoFGOXshP0oWF6Lm8nLu/phcFo2F4wFMCAbg35vGAzWvnEtCuSnmMxWqi12rlk8kWEQvB+x56CWdFDo4AUCkXxUVaCnkmxrh1tuwoaX+4JeOnx9lFh0rew9RJ+5jbM5JLJy3HaapUVrVAohpFS0IUQZuBJYCrhBtCPAu8DTxEuc9ICfEFKGRRC/A6oB3xAv5Ty6sIMOzsyKdbl9nRTY7bT5StM6F6qKBK9hJ+3jr2DZjCqhB+FQqFLOlEunwbapZQrgKuB/wL+A3g4ss8AXB85dwawXEq5cqyJOSQP0xuKxWhhIEHmaD6YWz8roZWdKuFH1VtRKBR6pCPoa4BH4rb9wHnAm5HtPwCXCyHGA3XA74UQG4QQn8jrSPNAJjHfnoBn2IKoEQ2r0RorWXvxxAuZWHVGVr72UJJrVMKPQqHIhpQuFyllD4AQwg48BzwMfEtKGY2j6wZqAQvwbeC7gBN4WwixSUp5KtG9HY5KTKaRKQDV2Bh2b9xdfxuhzX7e/ODdhOdGwxiH4qis49+ueJA+fz8OWy2/3Pkbjh8/kdV49rj2UOP4FFbTcCu9xm+lodJJa1/78HlUOpk+aaLuddE5ljJqjqVDOcxzpOeY1qKoEGIK8ALwfSnls0KIx+MO24FO4CTwQymlHzglhGgGBJBQ0F2uvqwHngmNjXZaW09btddOuZp3j+gX1aq11OBOEH/e3ufiw4/aaKysp62nm3cPb8t6TK19HRw4djxhZMkc5yzW9Q2vLTPbOYsulwcYPPahcyxF1BxLh3KYZyHnmOhFkdLlEnGlvAJ8WUr5ZGR3sxBiZeTvVwPrgcuB/4tcUw3MBXbnNOoC8eKhVxNGsCxonJOwtkr8Qqbb051T44xUi6JDs0frbQ5WTl6uEn4UCkVC0rHQHwIcwCNCiKgv/T7gCSGEhbBoPyelDAghrhJCvAsEgYeklG0FGXUOJFtwtBmtXDd9FZrBmLKueK3VntSaB7AYLHhD+guYqVLrVcKPQqHIlHR86PcRFvChXKJz7v35GFQhSbbgGI0NTycd3mK0ML9hTtLSuXpiXm9zxO41tPyAQqFQ5EJZJRZB8j6jiWqrJBLcW5qu42DX4bS7H9VZa3hg8b1UmCpSdglSnYQUCkWmlF21xUwqDKaqK27UjHx58b1cPOlC6qw1GDBQZ61J+Gy3p5t+vyeWNNThcREiFCs38Pz+l2LnpnOOQqFQxFN2Fjrkt8KgUTNyq7iRG2dcg9vTTYXJymObv5fwC6DCZE3ZJQhQnYQUCkXGlKWgF2LBMb64VbJmzf1+T1pJQ6nOUYW0FArFUMpS0KMUqsJgsi+AQCiQ0ocPJDwnlw5FCoWitClrQS8Uyb4AjBiTWvDR8xKd0+vv57cH/qgWRxUKxTCUoBeQRF8A6fjwo39/98TmQUXCPAFPTOhvabouFvpY48+ur6lCoSgdlKCPAun48I2akeunr2J7a4tu1ccdrbsIBAPsat+Dy9NJQ6WTOc5ZynJXKMoYJeijSCofvtvTTafHrXusw+MalNTU2tceq/2i6qUrFOVJ2cWhFxPJ6rdrCX51ql66QlG+KEEfwyRLggoS1N2v6qUrFOWLEvQxjl7VxYsnXYjDom+5p6riqFAoShflQx/jJFpATacipEKhKC+UoBcJQxdQh4Y+NlY6mR2JclEoFOWJEvQiZajlPn3SxEgnI4VCUa4oH3qRE7Xc9XqMKhSK8kIJep7xBry09rWr0EGFQjHiKJdLnlANKRQKxWiTUtCFEGbgSWAqYAUeBd4HngJCQAvwBSllUAjxNeAawA/cL6XcVJhh5wePL4C7x0OF1US/x09tdbgeirvHQ221Fas5fSGONqSIEm1IASpzU6FQjAzpWOifBtqllHcKIeqBZmAb8LCUcp0Q4ofA9UKIw4T7jC4FpgC/BpYUaNw5EQgG+dXa/TTvbaW9y4NmgGAILCYDGAx4fUGcdgszz3Jy+xXnUGk1x8RfT+iTNZ5WDSkUCsVIkY6grwGei9v2A+cBb0a2/wBcCUjgFSllCDgihDAJIRqllK35HHA++NXa/by25cPYdjAU/tPrDxH+6ICObi8bW06yde8pGusq6e334ur24qyxsrCpkVsvnYFRCy9BJGs8nW1DCtVAWqFQZEpKQZdS9gAIIeyEhf1h4FsR4QboBmqBGqA97tLo/oSC7nBUYjIVzr884PXj6vIw4PXT2BjOnnT3eGje25bBPYIcPdUT227v8vDalg8JYeCem+djs5io8VtpqHTS2tc+7PrGSifTJ01MOwolEAzwi+2/ZvOHO2jr66Ch0smSyfO5c8HNKX3x0TmWMmqOpUM5zHOk55jWoqgQYgrwAvB9KeWzQojH4w7bgU6gK/L3ofsT4nL1ZTbaNIl3qXR0eWh0VDBvmpMQsFW20tmTewTK61uOsm3vqZi1Psc5K1btMJ7ZzlmR+PD0YsTX7P3dIF98a187L+99g74+X1JffGOjndbW0q7houZYOpTDPAs5x0QvinQWRccDrwBflFK+HtndLIRYKaVcB1wNvAHsBx4XQnwLmAxoUsr0TeE8MtSlcsrVz+vvHcv7c6LWOsCtl+beeFr54hUKRS6kY6E/BDiAR4QQj0T23Qc8IYSwALuB56SUASHEeuAdwvHtXyjEgFPh8QVo3juybvvmvW3cfMn0nBtPF8IXr1Ao8key4IixQDo+9PsIC/hQLtE59+vA13MeVQ64ezx0dI1sCryrewB3j4dxjsqcGk9H65+naiCtUChGlqFuXL3giGwIBoN8+9v/xv79+zCbzXzlK48wefKUrO9XcpmitdVWnDW599e0mtP/0TjstlgMey4kq3+uqigqFKNH1I3b3uUhxGl366/W7s/pvuvXr8Pr9fKjH/2ce+65l//6r+/kdL+SE3Sr2cjCpsac71NpNWJJU9QXNjXk7fNLr/75ysnLVRVFhWKUSObGbd7bhscXyPreO3ZsY+nSCwGYO3cee/bszvpeUKKp/7deOgMI/7Bd3QPU19pwdQ3g12/yo4urx5f0uAFw1thY2NQQe14+SKeBtEKhGDmSuXHj3a3Z0NvbS1VVdWxb0zT8fj8mU3bSXJKCbtQ0br+8iWuXTeXAh2427v6I1s6BvN2/vsbKfZ+cT6OjsmALI7n44hUKRf6IunHbdUQ9V3drVVUVfX2nw7dDoVDWYg4lKuiBYJBfvr6PjTtPMODNwCxPk4VNjUwepxYoFYpyIOrGjQ+FjpKru3XevAW8/fZ6LrvsClpadjJtWm5f+yUp6L9au5+1eYo7NxsN2KssdHZ7cNjz72JRKBRjn6Fu3HxpwcUXf4zNm//MPfd8nlAoxEMPfS2n+5WcoHt8AbbKU3m7ny8QorPLw+JZjXxm1Uwqrea83VuhUBQHUTfuzZdMz2scuqZp/OM/PpSHEUbul7c7jRHcPR46uvPbXCIIbNrdymPPNBMI5t+Fo1AoigOr2ci4Aq6d5UrJCXpttRWnvTBRIUdP9fCLVySnXH05hSopFApFISg5l4vVbGSRGKe7gJEP1m87wVvbTlCfp0wxhUKhyBclqUS3XjqDC+aOL8i9ozWD85UpplAoFPmiJAXdqGl89qqZ1OehBEAqcs0UUygUinxRkoIOYdfLgnMaCv6caKaYQqEofbwBL6197XgD+Q28yBcl50OPxzACz8hXYS6FQjF2CQQDPL//JXa07sLl6cRhrWN+Y7jfQapOYqnYtauFH/zgCf7rv36c8zhLStDjaxUDbNtX+P4alTYTJuNIvDoUCsVo8fz+lwZ1EuvwuGLbyTqJpeKZZ/6HP/3pZWy2ipzHCCUi6Hq1ipvOrNOtvZBvjp7q4Vdr93P75U3A2C+Ar1AoMqOQncQmTZrMN7/573zjG/8vlyHGKAlBH9pyrr3LwzstH43Y85v3tnHDirP5zfpDeS+Ar1AoRpdCdhJbufIyTpw4nsvwBlH0gj4aLeeG4uoe4NlX97Gx5WRsX3y/0aj1rlAoio9i6iSWlqALIZYCj0kpVwohFgE/JNzGfhtwn5QyKIT4HVAP+IB+KeXVhRp0PKPRcm4oDruVPYc7dI9F+40q94tCUZxEO4nF+9CjjLVOYil9AUKIB4CfArbIrh8D90spVwBu4PbI/hnAcinlypESc8hfy7lcmHmmA1eC+jEqrFGhKH6KpZNYOhb6AeAm4BeR7clSyo2Rv78NXC+EeBWoA34vhKgD/k1K+WLeR6tDslrFhcZiMnDR/IncfMl09hxxFaQAvkKhGH0K2UnsjDMm8uMfP5WXexlCoVDKk4QQU4HVUsoLhBAbgQellG8KIb4P2IGHgE8B3wWchIX+IimT12hRR3MAAAvCSURBVLH1+wMhkyl3V0QgEOTJ3+/i3ZYTtHb2k8aU8sY4RwUXzD2DYCjEixsODTt+3Ypp3HXDvJEbkEKhKAd0Y6WzEXRBWLgDwGagFngAsEgpeyPn/x/wPSnl+mT3bW3tzqv0enwBjrf28I2n38vnbdPisvMmYTAYdAvgj0SUS2OjndbW7oI/ZzRRcywdymGehZxjY6NdV9CziXK5Bvi8lPK4EOJ7wB+Ay4EvAtcIIaqBuUBu7auzZMA7OnVVtu1r59G7lua9AL5CoVCkSzaCvg94WQjRB7whpXwZQAhxlRDiXcL9IB6SUhY+TZOwVd7RNcBr733Ijv1tI5JMpEd89+9sO4ArFApFLqQl6FLKD4ALIn//PfB7nXPuz+vIUhCfHTpaIh6PWvxUKBSjTdEmFg3NDh1tcu3+rVAoxj5Bjwe/242pthbNOvYMuKIU9LGQHRolvnORQqEoTUKBAK1rVtPTvBV/Rwcmp5PqhYtovOU2DMbsDDm/38+//us/ceLECXw+L5/97F+xfPklOY2zKAV9LGSHAiw6p4HPXj0Te+XYyRRTKBT5p3XNajpfezW27W9vj22Pu+2OrO75pz+9TE1NHY888g3c7k4+97k7chb0oqwaNdrZoTaLEZvFSPO+Nv75qc08+9peAsHgqI1HoVAUjqDHQ0/zVt1jPc3NBD3ZGZcf+9jl3HXXPbFtozF3+7ooBT2aHTrSaBpMcFQw4A0w4A0QQvUWVShKHb/bjb9Dv1aT39WB3+3O6r6VlZVUVlbR19fLww9/mbvu+ttchgkUqaBDuBH05Ysnx/qGapEwe0sBnUjBIJx09eseU71FFYrSxFRbi8np1D/mcGKqrc363h99dJJ7772Hq676OFdeuSrr+0QpWkE3ahq3X97E/OnhOsTBSM6p1z8641FFuBSK0kSzWqleuEj3WPXChVlHu3R0tPOlL32Rv/3be/nEJ67PZYgxinJRNIrHF2DHgfbRHgag4tAVilKm8ZbbgLDP3O/qwORwUr1wYWx/Njz99M/p7u7mqad+ylNP/RSAb3/7CaxWW4orE1PUgj5Wol1AxaErFKWMwWhk3G130HDjJ/MWh37//f/A/ff/Q55GGKaoBT0a7TLSmaJTxlXTN+AfVoRLoVCUNprVimXcuNEeRkKKWtBHuha6o9rCeTPHceulM/AHQqoIl0KhGFMUtaADMcs4WrbWYjYWrOKipp2uWGk1G1URLoVCMaYoekGPRrtEy9ZWV1r4zfqDNO9to71rIOX1NotGKAQeX+rEINX4WaFQjGWKNmxxKFGLudJq4vbLm3j0rqVcNHdCyusa6yrTEvN4VMy5QqEYi5SMoOux54gr4bHaKjMXn3sGvf36zZ2ToWLOFYryxOcL4Hb14xujBl3Ru1wSkSqksavXx879Hbh6Mhd0FXOuUJQXwWCQjWsPcGhvGz1dHqprrJzd1MCyS6ej5dBiMhAI8Nhjj3L06GE0zchDD32NSZMmZ32/krXQUxXwCgGuLK1sFXOuUJQXG9ceYOeWY/REjMSeLg87txxj49oDOd337bfDbZd/8IMn+au/upvvfe8/crpfyQp6Pgt4aYZwi+36GhuXL56sYs4VijLC5wtwaK9+R80P9rbl5H65+OKVPPDAV4FwXReHoz7re0GaLhchxFLgMSnlSiHEIuCHgAfYBtwnpQwKIb5GuIG0H7hfSrkpp5Hlgajw7jjQTqurn1CSc+uqLXQmcL9csnASVy2ZomLOFYoypK/HG7PMh9LT7aGvx0utoyLr+5tMJh599Gu89dY6Hn30sazvA2lY6EKIB4CfAtECAz8mLNgrADdwe0TkLwGWArcB/53TqPJENKTxvx+4lPNnj094nmaA+dOd/PNfnc/HFk2ivsaGZjhtkd9++TmMc1QqMVcoypDKagvVCdy31XYrldW5N7h5+OF/4pe//DWPPfYo/f36FV3TIR0L/QBwE/CLyPZkKeXGyN/fBq4HnMArUsoQcEQIYRJCNEopx0afOGD/h50JjwVD8Nb2k1jMJu68UuD5WEBlgSoUCgDMZiNnNzWwc8uxYcemNjVgzkEj/vjHl2htPcWdd34Om82Gpmk5LbKmFHQp5a+FEFPjdh0UQlwipXwTuBaoAmqA+LKH3UAtkFTQHY5KTKbCC+aJtl46ulMvgO440M7dN1fQaDGR/Trz6NHYaB/tIRQcNcfSoZjmef2nzqWiwsLelpO4O/upraugae4Errx2NpoxsQCnmuPNN1/Hgw8+yP3334Pf7+fhh7/K5MkNWY8zm7DFzwHfjbhiNhP2pXcB8SO3A4lN4gguV18Wj88cR20FTnvqIl5tnf0c+KC9KFP6GxvttLZ2j/YwCoqaY+lQjPM876KzmH/+ZPp6vFRWWzCbjbR39CY8P905Pvzwo4O207km0YsiG9v+GuDzUsprgHrgVcKul6uEEJoQ4kxAk1LqLwuPAjaLKa2IFxVfrlAokmE2G6l1VOTkZikk2Vjo+4CXhRB9wBtSypcBhBDrgXcIvyS+kL8h5of4Il6Jaryo+HKFQlHMGEKhZMF8haW1tXtEHh7/6ePxBejoGuC1LUfZcaBjWE1zYw4LEqNJMX7CZoqaY+lQDvMs5BwbG+0Gvf0lm/qfCKvZyBn1Vdx51Uw8PhXNolAoSoeyE/R4VE1zhUJRShSnf0GhUCgUw1CCrlAoFCWCEnSFQqEoEZSgKxQKRYkwqmGLCoVCocgfykJXKBSKEkEJukKhUJQIStAVCoWiRFCCrlAoFCWCEnSFQqEoEZSgKxQKRYmgBF2hUChKhJIpziWEWAo8JqVcKYQYB/wEcABG4DNSygNCiLuAuwE/8KiU8sXRG3F2DJnnucAPCc9nL/DXUspgsc5TCGEGngSmAlbgUeB94CkgBLQAX4jM8Wv/v727CdGqiuM4/pmScGNBLWxhMJv42wstKorypVkIYS6CNolJLapFZGW4EMpBCKKQ3pyKNrYpioHMhS0CKVQsCje2ifhTLXLRC6VmBlmN1eLcB56GO4Q2zNO9ne/qOffexfnyP/w459yHe5TDVmaUQ8sPj6LPZ8scjkfxIs4oJ4DdnZnfdbWOtHtm5t7m3gY8lJk3Ne1Oes5Ry4+NMHt6MUNvjsPbhcXNpR14IzNXYxuWR8SleBgrcCueiohOHU/U4rkdT2TmSmVAreu450Ycy8xVWIuX8By2NdfGcHtEXItbcCPW4+UR9fdcaHPcqQTcBPZga8frSLunZhJyr1JLHfdscxxp9vQi0PEl7hhqr8CyiHgPd+EAbsCHmflrZp7EF7hmoTv6L5nteQQXR8SYco7r77rt+RYmh9ozuA4Hm/a7WIOV2JeZf2bmUSyKiH8+Y/C/QZvj+sz8pGkvwmndriMtnhFxCZ7G5qHrXfZsq+VIs6cXgZ6ZbythNmAcJzJzjbKc3YoLcXLomVO4aKH6OB+0eH6OKXyGpcrg6axnZv6cmaciYgl2KzOcscwcfJ9i4NIrx8z8BiLiZmzC8zrsSKvnJF7Fo4rLgM56zjFex40we3oR6C0cw97m9zu4Hj8ps9gBS/DjAvdrvtmJVZm5HK/hWR33jIjLsB+vZ+ab+GPo9sClb44i4k7lfci6zPxexx35u6cy+bgcr2AaV0bECzru2VLLkWZPb16KzuID3KYMpNX4FIfxZEQsVvabr1BesnWZ48pgga+V5V5nPSNiKfZhU2a+31w+EhETmXlA2afcryxZd0TEM1iG8zLzh1H0+Wxpc4yIjcoLs4nMPN482tk6Mmctr2rujWM6Mzc3+8ud9JzDcaTZ09dA34JdEfGAstTZkJknImIKh5SVyeOZeXqUnZwH7sN0RMzgN9yfmd922PMx5d8BkxEx2Jt8BFMRcYGytbQ7M89ExCF8pDg+OJLenhuzHc/H1fgKeyICDmbm9g7XkfZars3MX4Yf6uF4vccIs6d+PrdSqVR6Ql/30CuVSuV/Rw30SqVS6Qk10CuVSqUn1ECvVCqVnlADvVKpVHpCDfRKpVLpCTXQK5VKpSf8BWvrClWyzu3EAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEFCAYAAADzHRw3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzsnXmcXFWZ97+191LV3dVLAtnI2gfIRhYIS0JCJoiIAiIIA7ihiI4wOL7vjA6C4gw6wow6gCIjGnGBAUHwRUHZwpIFCCEhG+R09kC23qqX6qX294+q26nuvrf26u6qOt/Phw9969577rldnd997nOexRSJRFAoFApF4WMe7QkoFAqFIjcoQVcoFIoiQQm6QqFQFAlK0BUKhaJIUIKuUCgURYISdIVCoSgSrKM9AUXxIoSYCuwFtsd9bALulVKuznLsvwBPSikfFkK8C6yQUnYYHFsNPC2lXJnmNa4EbpZSrkjx+BXAT6WUc9K5Ttz5XwLsUsoHMjlfoVCCrsg3fVLKM7QNIcREYIcQYpOUclsuLhA/vgFu4KxcXCvPLAV2jPYkFIWLEnTFiCKlPCyE2A00CiEWAl8EKoFOKeUFQogvAv9A1B3YRtRC3iWEmAD8BpgAHATGaWMKISJAg5SyVQjxr8DngCCwG/g88GugPGbJLwIagXuBOsAC3Ke9MQgh/g24Lnbt3Ub3IYS4Afg/QAhojV0zfv/DwA4p5X8N3RZCfBX4CuAH+oGbAAFcClwohOiTUv5MCPFt4FOx38UB4B+klEeEEK8C7cCpwM+Bw8DtQDg2n3+WUr6e+JtQFCPKh64YUYQQ5wAzgbdiH80m6i65QAixnKgwLpNSLgDuAZ6OHfcz4E0p5WzgH4mK2dCxLyUq4OfE3B77gZuBL3DiTcEEPAl8S0q5CFgO/F8hxNlCiMuICugZwLlAtcE9zAfuBj4qpZwHPAN8O8X7twD/HTv3TOAXwFIp5dOxcX4SE/PPAnOBs2Lzfg74ZdxQHinl6VLK+4H/JCr2i4E7gBWpzEVRfCgLXZFvNMsYon9vrcB1UsoPhBAA26SUXbH9lxAV+w2xfQBuIUQtsAr4vwBSyj1CiDU611oFPCGl9MSO+wYM+PI1GoEZwOq4a5QDC4DTgaeklN2x81YTfXgM5e+A56WUH8Su89+x41ck+2VIKUNCiCdi9/gs8DzwqM6hHyfqJtoUm6cFqIjbvzbu58eAp2PjvUj0QagoQZSgK/LNIB+6Dt64ny3A76SU3wQQQpiJulg8QISoda0R1BkrGDuO2Pk1QM2QYyxE3Tvxfv3xQCdRSzfZNfSuUw6cMuSYofO1az9IKa8XQswh+gD6FvAZ4NM687xbSvnz2DUcRNcCNAZ+b1LKb8cePhcSfUP5PxTGmoEixyiXi2Is8Tzw90KIk2PbXwFejv38N+DLAEKIKcAFOue/BFwhhKiKbd8JfIOoAFuEECZAAn1CiOtjY00muhC5CPgrcJUQoib2MPmMwTxfAVbFzfMmhlvFLcDi2DUmEHXtIISoF0J8ALTFLPvbgTNj5wQBW9zv4ktx9/JvwO+GTkQIYRVCHAAqpJQPEl1/mBd7AChKDGWhK8YMUsoXhBB3Ay8KIcJAF3CFlDIihPga8GshxPvAh8C7Ouc/J4Q4HVgfc1PsBG4EeoGNse1lwGXAvUKIfyEqoHdIKdcDCCHmApuIvhVsBRp0rrNdCPHPwN9i1zkK3EDUnaNxP/CIEEISXdBcEzu3VQhxF/CyEKKPqIjfGDvnr8CPY2PeDUwE3owt+h4ian0PnUtQCPF14FEhRIDowugNUkqf4S9aUbSYVPlchUKhKA6Uy0WhUCiKBCXoCoVCUSQk9KELIWzAamAq4ADuAvYQjZ01EfUx3hILxbqR6OJQELhLSvmXPM5boVAoFENIZqFfT3Q1fhlwMfBT4AfAbVLK84jGxV4qhDiJaLzuecBFwH+oVXaFQqEYWZJFuTxBNKtOIwh8KmaR24GTgONEY17Xx1bWfUKIPcA84O1EgweDoYjVasl48gqFQlGimPQ+TCjoUkovgBDCRVTYb4+J+SlEY347icb1To/9rNGNQdp0PB5Pb0ozzyUNDS5aWrpH/LojibrH4kDdY3GQj3tsaHDpfp50UTSWePEK0Qy+RwGklAellLOAB4EfE40Xjr+CC9AtZapQKBSK/JBQ0GMp0S8A34yrRveMEGJW7JBuookMG4FlQoiyWO3p01BlQBUKhWJESeZDv41o/Yg7hBB3xD77NvCwEMJPNAPvS1LKY0KI+4gWDDID35ZS9udr0gqFQqEYTjIf+q3ArTq7ztM59iHgoRzNS6FQKBRpohKLFAqFokhQgq5QKBRFghJ0hUKhKBKUoKdJIBCi09NHIBAa7akoFArFIFQ99BQJh8NsWLOX/U2teLt8OKscTGus59yVMzCb1XNRoVCMPkrQU2TDmr1s33R4YNvb5RvYXrpqltFpCoVCMWIo0zIFAoEQ+5tadfcdaGpV7heFQjEmUIKeAr1eP94u/Y5e3m4fvV7/CM9IoVAohqMEPQUqnHacVfrVgJ0uBxVOu+4+hUKhGEmUoKeAzWZhWmO97r6pjfXYbKoEsEKhGH3UomiKnLtyBhD1mXu7fThdDqbGolwUCoViLKAEPUXMZjNLV81iyfLp9Hr9VDjtyjJXKBRjCiXoaWKzWah2l4/2NBQKhWIYyoeuUCgURYISdIVCoSgSlKArFApFkaAEXaFQKIqEhIuiQggbsBqYCjiAu4BDwP1ACPABn5VSHo+1oDuPaJ9RgMuklJ15mrdCoVAA4AuE6PT6qHY6cJR45FmyKJfrgTYp5WeEEHXAFmA/cIuU8l0hxE3AN4FvAAuBi6SU+kVPFAqFIoeEwmEeX7OHLU0ttHf5qK1ysKCxgatXzsRSohVQkwn6E8CTcdtB4Bop5dG48/uFEGZgFvALIcR44FdSytU5n61CoVDEeHzNHl7a9OHAdluXb2D72lWNozWtUcUUiUSSHiSEcAHPAA9JKR+NfXYu8CvgfKCfaDPpHwMW4BXgBinltkTjBoOhiNVa2q9I6RLwB+nu8uGqcmCzqzQCRWnS7w/ytXvW0OzpG7ZvnLucn/3LSsqK+9+HSe/DpHcshJgMPA08ECfmVwPfBi6RUrYIISzAvVLK3tj+NcB8IKGgezy9ad1BLmhocNHS0p38wDFGOg02CvUe00HdY3GQ6T02e3pp0RFzgNaOPvYeaGOcuyLb6eWEfHyPDQ0u3c+TLYqOB14AbpZSvhz77HrgJmCFlLI9dmgj8JgQYiHRyJmlwG9yM3UFqAYbCkU81U4HtVUO2nTKWrtdZVQ79aujFjvJVg5uA9zAHUKIV4UQa4lGuLiAp2KffU9K+T7wCPAm8BrwWynlznxOvJRQDTYUisE4bBYWNDbo7lvQWF+y0S4JLXQp5a1EfeNJkVLeA9yTi0kpBpNKgw1VX0ZRaly9ciYAW5pa8XT343aVsaCxfuDzUqSoVw2KBa3Bhp6oqwYbilLFYjZz7apGPrV8hopDj1GawZoFhmqwoVAY47BZGOeuKHkxB2WhFwyqwYZCoUiGEvQCQTXYUCgUyVCCXmCoBhsKhcII5UPPMYFAiE5PX0ahhNmcq1AoFMpCzxHpZHLm8lyFQqHQUIKeI7LJ5FRZoAqFIhco8y8JqbhBEmVy7tp2FJ8vkNG5KgtUUYj4AiGaPb341N/uiKMsdAPScYMkyuQM+MOse3EPf/fx03T3qyxQRbGg6pOPPuq3bIDmBtHEVnODbFizd9ixWianEUcOegZZ2vFWf6JzCyEL1B/y09Lbhj/kH5HzFGMXrT55W5ePCCfqkz++Zs9oTy0rCumNQ1noOiRzgyxZPn1QDLjNZmHCFDdNO47pntPj9dPr9eOqduha/VNn1bHjnSPDzhvLWaChcIin9jzLtpadeHwduB01zGuYzRUzL8FiNp5zpucpxja+QIgtTS26+7Y0tfKp5TMKLpOzEN84xuasRplU3CBDWXrhTKx2/V+nZmmve2mPrtUPMHfxRFxVDkwmcFU5mL1wArMXTDT0oY92iONTe57l1Q/X0e7zECFCu8/Dqx+u46k9z+blPMXYptPro93g34ynu59Or/6+kSBTC7sQ3ziUhR4jEAgNZGBmUgzL4bBy2ryTB0WraJwyq4431uzhvXePDtsHcGB3G1d+fhGnzTuZcCTM++8e5eCeNnZuPjLMd5/vEEd/yE+nr5tqhwu7Rd/d4w/52daiXx15e+tOLpvxUd1zMz1PMfYZi/XJs7GwC/WNo+QF3Uggp8xw896W4S6URG6QofVWKp12JpziJhwO646l4e3y8fufv0EwEMFmNxPwhwftiw9hzFeIYzqukE5fNx5fh+447f0dAw+EoQ+GVM5rqKjL+B4Uo4dWnzy+x6fGSNQn9wVCwyouZtNzNJU3jrHSESmekhd0I4G02qJPcJMJIhFwuuxMEw3DimHFW/Y2m4Wlq2Zx5rJprHtxD0cOeWjacRyTbve/wQQD0d6u8WIez4GmVhade0pavv100FwhGporBOCqxksHHVvtcOF21NDu8wwbx+2oZs0Hr7OjddewB4Mm7r7Q8H8odoudaod+W61ck8pbiCJ9Rqo+ebx4Wy0mXSv88mXTsrKwx+IbRyqUtKAnWvwMBqLCqvXQnjKzfsACDgRCeLv62b7pQw7ubR/m+nh77f5BC6Qp9OFOirfbR1tzT15CHNN1hdgtduY1zB70ANAot1Xw+uE3Bra1B0NfsI8rZl4CBr+LFJ55WaMWZPNLvuuT67lQKspsfNDsHThGs8LbO/t1xRhSs7BH+40jU0pa0BMtfg5l985jnL1iGm+vPTDgnolHs+zDoTAH97YbjJI5lU47e3Y1D7wxDCWbEMdMXCFXzLwEiAp+e38HtWU1nF57KjvbdumO89axd3i/vQlfWP/33R/yGbpqckU6byGKzNHqk+caPReKkWhv3q1vqAHYbZaULOxC7IiUrEm0DVgNTAUcwF3AIaJ9RUOAD/islPK4EOJGos2jg8BdUsq/5HHeOSHR4udQAv4wa1/Yze6dzQmPO7C7jR6dKJhscZTbeN9gURUyD3H0h/z4QwFqHNW6ol5bVjPgChnqqriq8VIum/HRQUK87sibhtfq8ht3Pq911Oi6aj4+7UK8gd6sBV4tyBY2iRYp80UhdkRKZqFfD7RJKT8jhKgDtgD7gVuklO8KIW4CvimEuAf4R2AxUAasE0K8KKUcvVilFNA6AelFpuhx+MBwn/FQerx+Kp32rEXdZjcTDIRxuhxMnu6m6T39B4nJBKefcXLajS76gn080fQMTZ69dPg6sZv1xWxu/WwsJgtPND2j66qwW+wD1nsi33oyjFw1bxx9G3/In7V7RC3Iji56i5bpkGiRMu25+ENpLWrm640jHyQT9CeAJ+O2g8A1UkrNVLQC/cBZwPqYgPuEEHuAecDbOZ5vzomPTOlO8AdjtZrp7TGuyaLhqnIwabqb9981jmoBsNpMWK1m+vtCWK1mTGYGBHxqYz1nLptGf28Au8PCmmclQYPF0ggw/6wpKYcsan7kN468Pcj9of1sN9sJhAPUltUwtz5qIT+y60neOvbOwLFGropEvvWh1Nir6PR3J3XVaAuo2bpHEj1s4t9CFLklV8k5iRYp06W2KrVFzWwfQqNBQkGXUnoBhBAuosJ+uybmQohzgZuB84GLgM64U7uB6mQXd7srsFpH/hfV0DD4H+8n/34hAX+Q7i4fb7y6h3feODTsnHlnTmLP+810dfQnHPu0+RNYfO4pSQU9GIgQjCU6BINhGsY7ufIzC6mpq8RmtxIOhXn+mZ1sffsD/D7jhIjq6jJOmVqLzT74qxx6jxoPb/lDQsENhoOcN+VMvrDgKv74/nP8cNO9tPbprwlsa9vB9a7LqCpzDnz2JfdVHPAe4FDHEcLoP4TGVdTxgwu/RYevEyLR5dB1zxu7auJ5r/19qtxXJbxHI84+5Qyea3pl2OdLppzBxJPGpnWe7j2ONR7603bd0MGKcjs3Xj4XSP0ez5s/kWfW7hv2+fQJVXT3+mlJ8m/zxDgTmDShxnB/KBRm9Z938uaOo7R09NFQU87Zc07mhk/MxmLJLNdjpL7HpIuiQojJwNPAA1LKR2OfXQ18G7hEStkihOgC4mfsAvTfb+PweHozmnQ2NDS4aGkx9uUuXjYVfyDEftmCt9s/EK545rKpeLt9CQW9cc54Fpwzmd6+9N0tLce9/O+vN7LsQkHduEre2XAwJVfQ+EnVdHT2DfrM6B79IT9vHnw34Xhhwqw7tJH97R9y2Du8HEE87X0d3PrsnSwcN4+rGi/FYo66Zg50DI8MiMdmsvPbTX9ie8tOPP4OamzVhuGMQ2npbWfv4SPMmDiBvYePJPWtx/v9PzrhI/T2BgYt5M6tn81HJ3wk4d/EaJHsb3Ws4wuEWL9V/294/dYjXHzWZCZNqEn5Hj9xzhS6vf1s2d1Kp9dPbdWJRcpgKMLvn5es1ym/UWa34A+EBhY1P3HOlITXfPSlpkEPoWZPH8+s3Udvnz9p/Loe+fgejR4QyRZFxwMvADdLKV+OfXY90cXPFVJKzXTbCHxfCFFGdPH0NGBHbqY+siTq3bn0wpnsb2rRjRV3Vjk4/6JGzGYzXUMENlU62vr582NboxspxPHZ7GaWXpj6insiP/JQkom5Rk+wh7VH3mBf10H+aeFNhguPg8buOcrhnhMLvB2BzgRHD8ZutvPSwdf42bbdtPa2G/rWE4Uoxi/kqoXQ/JFKcs6kFMfSXDfb9rbR6fVT43Qwb2bdgOvGYobPf+xUysusw6JSLl82HW+vPyXXSaFmiGoks9BvA9zAHUKIOwALMAc4CDwlhAB4TUr5XSHEfcBaovVhvi2lTO39Z4yi17vT4bBxqlF6/8y6gdj0vTIHq/EpxK5Pa2xI6DuPVjRsx2SKUGV3EQj7qbFX4/GnJurpcNh7hN/sfDyjBdF43I5qOnxd2Mw2/OHhbzq+sI91R0+4Z/R86/6Qn8fk0wn9/vkMj1REyWVyztCQRY/XxyubD2Mxmwas5kRRKRWO1CK0CzVDVCOZD/1W4NZUBpJSPgQ8lItJjWWGpfe7HDjKrBzY3crOzalZtbnAajPRtOM4Rw51DKvlEgqHeFz+P946umlY3LfDnL8Mt+1t72HGRDiVp5EBHb4uIkR0xTwRbxx9m4unruSvB9awtWWH4ZvItpadhMIhdrYNz2RVyUW5JVfJOelazUOjUtJZ3CzUDFGNkk4syoShLpmtGw+xc4txfHi+0EoF6NVy+d3WP/L64fW65xkl9uSKbMQcIJLh+b6Qj59sfpBjvYnzBNp9HtYeGR4eCSq5KB/kIjknU6s53QgbTfjnzajjlS3DjbOxnCGqoQQ9RYbWbLHZLFQ47XnJCjXCajdD5ERZgni0Wi4Rc4iNH2wdsTmNJZp7k7u6zJh1o29UclF+yEVyjrPChsNuod8/PNorkdWcanEuPeGfPM5JT1+ADq8vqwxRXyDE0dYeQoHQiDwMlKAPYahwJypXm07pgGz52JVzcFaV8YfVm3T3a7Vc/I5ewzDDYieVtwOjUEqVXJRfsknO+eNr+3TFHIyt5nTcNEYlBS5YOJGLzpyc0UNo0EOi20eta2SaY5S8oGsCXlZhHVSnpdJpZ2pjPSYTg7oJaS6OSCRCJJyDqlspUOl0MOEUd/Rnl52e7uH+ZafLQVmFhb8deN3QCo3HbrIRiAQzdnGMJGZMRIgmAJ1WK9h4dDP+yPDfQaL7dttrmNtwGttb3tddFFbJRdmTjq9aO9ZVbVxMLhQO8+iLTbz2rv7aVJndwuXLpuvuS9VNk0j4t+1p49MXzMzIss6mdG82lKygD7W8h9Yh7/H62bn5CCaDh6ncfsyw1G2umdpYh8ViYsOavfh9QYNj6vnTgWcHpc8nwh9JnvU6Vmgor+fL8z6H01ZBp89Lf9DHpuYtw4472XmSYbhlJBLBF/Jzet2prD86PIlpbv1s5W7JkHR81cOPLWPezDquXTVr2LGPr9mj68vW8AdCeHv9WMymYQ+SVBc38xHVMpqhjyUr6EProBuJc8RAs0dKzOvGVbJ01cxh89Ww2c2IuSfxwYQdbDj81ojMaaQ53tfCr3Y8Qltf+8Cirtlkxmay4gv7qStzM7d+NpdN/yj/b9/f2N66k7b+weGTHYFO3jr2Dg6zg4nOCfQFevH4OgeSi7TqkYr0SccaHX5sP69sPsyeDzv5zucXD4h6KsW4apwOnn/7A7btaR32IEk1wiYfUS2jGfpYkoKeqA76aGMyRx8ilU47U2bUMv+syfj6g4bzdZTZODpxF+uOpWaZFypHegZHEoUjYXwRP/VltXxj0T9Q7agCopEqF09dyQ82/jed/q5h4/jCPg57j3D+xHNYOfl8FYeeJelYo4mO/aDZy6MvNvGZi04FUivGVVlu45XNJ4yc+AfJ1StnEo5EKLOb6Y8ZX2V2C+fNPWnQ4mY+6p6PZuhjSQp6LhYzzWYI58FI/9TnFmG1Wtj+zocc3NPG+1uPUeG06zamhuhi6J7Du6MpXyVIa387d264h3MnnjUQS94X9CUs1Quws20Xn4xVi1RkTjrWaDKR3rK7lU+vjEaDJBJFE7D4tAb2fqifYbylqZVQKDzMXdPvD2EymYa5dnJd93w0m2OUpKCnUwcdouGCdpuV3l4/LpcDe5mVtuaenM+r0umgpraCt17bNyhJyUjMASpcNtoY2TrRYw1/xD8sCzRZGV8V1ZIb0rFGq50OapwOPF79f3edXv/AAyCRKNptZja932K4nN/W1c/mNHzY+ah7PlrNMUpS0NOtgw7Q2+PHUW6l/mQnxw8Pf5XPBVMbo+KSjjto0vQaNlpGxp8/1tnavGMgljxZGV8V1ZIb0rFGHTYLZzTWD3KTxDO0rO1QUbTborHoPp08jKF0GpS6TuTDzmXd8/iHhMVuI+QPjEgcev4CIsc4566cwdzFE3FVOTCZosW1GueM59T546mMtXKz2aO/Hq0Wua8vyH7ZRq839xEiZeXWlGLbK53R+bqqHMxZNIH1VX/N+VwKFY+/g8fk04TCIa6YeQnnTzwPi8GfuIpqyR1Xr5zJqsWTqKsqw2yCuqoyVi2epGuNXrtqFpPHOXVGGf4A0ETxrhuXcOcXzqSyLHX702xQ3C6XPmxfIESzpxdfwLi8tcNm4eT6yhHLMC1JCx1OpPCfuWwq617cw+GDHpp2HMdZ5WDKjFpOmVnP68/LEYtm6e8L8uar+1iyfLqhO6jS5eCqLyzC7wtR4bTz1L5naD7Uij1QQcDWT0RZ6rx17B3KreVc1XgpZpOJkE5c+kTnBBXVkkPScVlYzGa+8/nFPPKiZMOO4/gD2oKlmXAkQigcHubjdtgs2G2WtDoWGaWI5MKHnW5JgX5/kGZP74g0yihZQdd4e+0BmnYcH9j2dvl4f+sx3t+auEGFHjNPq2f/njZCgcySdbT0fSN3UE+3j8ce2sh1Xz2bEEEOvNXLrNbzsfnLCNj76XIf59iUXWAa+8lC+WR7604+csoK3jii3zCrP9hHKBLCUqoryXkiVZeFxWzGarEMiDlAvz/MmncOYzaZdBNv0u1YVGY3E4kw4J4ps1s4d0iES6akW1Jg2942Wjx9GXdrSoeSdblA7sIXy8otzF5wMouXTstYzAG6u6Lp+4vOm0pNrf4/jP6+IL9/4A3WvtxE1ZGJ2P0VmDBj91dQf3waJx06NePrFwvt/R38Qf7JsBCZtiCqGB2ShTrquTAcNgvzZtanfI1+f3iQr73fH8KsE+GSLunMXRP+Zk8fEU4I/+Nr9mQ1h0SUtKDnqhZLf1+Ig3vb2bT+EM6q7Pxzf3tqB797YAMd7cbdnHz9IQ7v1g/ZcnnGYwqV9NdKjb2Kg90fGO53O6rVgugokkqoYzyhcJhHX2pi6+6okGr+8boqB3+3aCIrF02M8987KLPr//1v2tVMd292zdtTnXsmD61cUNIul3TDFxPh7fKx573jhqUCUqW9JbVwyL4e/RIAdn8ZtkAZfsuJB4LDbMeXZn3xQsZudSSsvDjLPUMtiI4i6SbeDHVxaP7xeTPquO5CAcBVK6K1YfyBEN9dre9q6/D6uXP12yw6NXO3x2iWFEiFkjbltPDFXGJUKiDXOMr1/b9+ez8B2+BmUdWOpP26iwpfwEeVXd8Cd1gcXDVL1T0fTbRQRz2GLlomLJ61t33A0tX89w3uCmoTvCV7vNm5PVKduyb8euQzW7SkBR2i4YuzF5yMKYUenmOJqbP0H0Td7uODol0cZgeBUOEU4soFHYFO3bR/gHNOPpNym3GFP0X6pBK+N5QrV0xn8jgnmpFsNsHkcU6uXDG4emK67hmrxURFmS3p9bNxe6QSppnOQyuXJGsSbQNWA1OJNn++S0r5TGzfTwAppXwwtn0fcB6grTZdJqVMvfvvKGE2m5l/1pRR6TqUDXL78UHb5vIQzVWHolEucZwxbg4bj20eyamNSbQCXipcMXekG74Xz+Mv7+GDZu/AdjgSrefy5Kv7BkWKZOKeiR/XiGzcHqmGaWoCv21vG60dfYOyRdMpNZwOyXzo1wNtUsrPCCHqgC1CiDeA3wKNwH/GHbsQuEhKOTarXiWgwmnH6XLg7R6ZZhU5YUgwTbXThfucOgLtNXj6O6hyuJhXP5tLZ3yU3Z59WTdvLmRq7FV8fcFXCEXCKlwxh2RS8zu6wLnbsMb50NT8dDJRU6nQqJELt0eyME1N+G/6VDl7D7RR7XRgtZgyfgimQjJBfwJ4Mm47CDiBO4GLtQ+FEGZgFvALIcR44FdSytXJLu52V2C1jvw/roaGE/7VcCjMC39+D5+vsN0SnpZeFjXPoWyqmY0fbMXj6+S9domzsoxFE+fw4r61oz3FUaPD38VP3v05nr5O6itqOXPSPD4z/1MF0RQ6/m91LNHvD7Jtb5vuvm1727jpU+WU2YfLy0N/2m6Y+g9Ry9lit9FQX0m/P4iny8eNl8+lotzOmzuO0trRR31NOWfPOZkbPjEbi+WECB5t7aE9RaPsvPkTmDShJqVjgYG5uKscuveVjNmN44Ho/es9BCvK7dx4+dy0xx1KwplJKb0AQggXUWGajTisAAAgAElEQVS/XUq5H9gvhLg47tBK4H7gx0Tr/r0ihNgkpdyWaHyPxzg0Lx2Gto1LREODi5aWEzHI617anVZNl7HM5s372eFbN+BDb+1r57mmV3CYx3an8pGgvS/apailt43nml6htzeQtCm0P+Sn09c9aiV2h/6tjhV8gRD7DnfS4unT3d/a0cfeA23DrFdfIMT6rYn/rbldZfj7/dz7v7sGWbGnTnHzresW4g+EBtwU7e2DI8JCgRC1Ln33jNkUfamtjbk9PnHOlJR+t9m4lTS07zHR/a/feoSLz5qcsvvF6EGf9FEjhJgMPA08IKV81OCwXuBeKWVv7Jw1wHwgoaBnS6J+n+YUftljuS56JoT7TcNCFgHDBJtSJlFT6FA4xFN7nmVby048vg7cjhrmNcweKM9bqsSLW1uXLyqSOnl0Ru6MVGqcV5RZefr1vbz8zuA65+t3HOOdpmaWzptgmO2ZyD2z/IwJXHTWlLR91rlsJTcSoYwJVS/mPnkB+GYSF0ojsE4IYYktpC4F8r4Sp3Xx0eLItX6fG9bsHXZsIBCi09NHwH8ifnskmzyPBBGg9ug07H2VJZ9clIxE2aJP7XmWVz9cR7vPQ4QI7T4Pr364jqf2PDvCsxxbaOKmWcDJ6qUMjX5JFMqn8UGzl/Xb9ctu9PvDSUMOjSJQrr2wcaAsb6rkOjloJEIZk1notwFu4A4hxB2xzy6WUg5615JSvi+EeAR4EwgAv5VS7sx6dglIZF1rNVFsNsswK76qpoyTJlWz9MKZOU0sGguYMVPfcgp1LVMI2PtUbReMG0cblc/1h/xsa9H/001k1Rc7icRtqDvjyhXTefSlJl03hZEFHU+/P7FQJurLmcva5qlY1NVOR8rXGYnGF8l86LcCtxrsu3PI9j3APVnPKEUSWdfe7mhNlGp3+bBenF0d/XR19LO/qYVT553M1Fl17HjHuBFtIWLCNFDbBeDYKe+P8oxGD6PG0Ublczt93Xh8HbpjlXJTjETiFonA/73mDKZPrMZhs/DoS02GbgrNXbJpVzMdCRq3JCIV90QuapsnCplM1M80kW89340vCjb1P5F17XQ5qHDaE1rxAX+Y7ZsOM2fRBOYunsjeXS0JOwMVKi7PeI5PkiVXWrfGXsUZ4+YNahzd3t+RtCl0om5HpdwUI5G41VaVDYh5MjfF5cuiRobJsN9QtDJiIis9XfdEpjHfiSzqRP1ME/nW89EdKZ6CFfREXYemNtZjs1no9PQldacc3N3G1TeeRa/Xx95dxbNAqqFX22UkMGNibv1ctrbmdV1clxpHFf965tdx2qONFK5qvJTLZnw0pYiVRN2OSrkpRqrugmRuikdf3M2GHYlLU5879yTMJhPrth3VFfZU3RO5iFDRs6jnzag1DNlM5A6KJ5fdkeIpWEGHaNp+JBJBbj820IjCZjcTiUTw+QIEgyEqXXZ6uo0tb2+3j66OPo7lqa3caKNX22UkCBPBYcnvn5dR0bFyawXl1sHp/XaLPWVXiWa9p2rVlwqpuAsSZ3Y62HWw3XD8uiGCe/myaTz64m52HfTQ4fWl7Z7IRYSKnkXd6fXx6hZ9N20+C2+lQkELutlsxmQyDeoqFPCH2fHOkQGRtxmU0tRwuhwQIaHoFzLxtV3sZhv+8MgkULntNTR5hkcbabisTkwmE12BzOOsfWE/J1WM41hv86DPj/Yc4+5N9/PNxbdkFGZoMVvSsupLhVTcBYks+VOnuA2tcxNw65XzmDTuhEurwmHjSx8/PSOXSTLXTypWdDzxFnW65QhGkoKObUvmI4//v9Eb1tTGepzVZViTCH8hEiEcDT+IRCuPlVnKRuzaM9zT6PAbl/LpDnpTEnMTxlXT3PZqfCH9B/Fh7xGe3P1M8okmQLPqlZgPRhM3I0E0Ch38+wsbE4TtOcBkMmxukW7IYbpFvdJhtApvpUJBW+jpxJGXV9qZMLmG5qPddHX04XQ5mBpLQtqwZu9AI+hiwoSZ+uZpgIn2kw7gDfVis9oIRHJvpduwESAwECa4t2N/wjrsbns1mEyGESUaJ1WO52iPvlXXWDuTt469Y3juttadfHLmJUqQR5hElryR9e7tD/CdX20ccLtcvmw63l5/2la5dr18W9H5jlbJlIIW9HTiyHu9fs5cNo1TptZy8ED7QJmAQCDEPplaQZ9CpbZ5CnXNUwjY++mpbeHw5PdyGps+0TmBaVVTWHfkzYGY72RCPX/cXMKREK8ffkN3vwkTJ1eOpy843P9fZnGwYvrZrDppJbvadxuWyu30dZdsmOFYQG/hL14I27v6MZshFGagv6jm51637Qg+fzilhUyjxc8zZtUPyjjVmDejNmsrOt/RKplS0H6GdBpUaKGMNruVanf5QM2XXq+/aP3nGmbMA31H3cdO4aRDp+IwOzBhwu2oodJamdX4fYFedrbt0t3nMDuwm09YyGUWB+dPPI8rZl7C8klLDceMEOFIz7FhD4bF4xbwH0vv4IaF11BuK2de/WzDMUo5zHCsognhXTcuYcnp4wgZvBj3+8Mp9+GMz2CNPycCrFo8iVpX1BrXWtdt29vGoy81EQpn/1aeiTsonxS0hQ7RSBeIZod6u31YbGZd98kps+p0C3fZHZboikwJJVO6POMJzGjmpkWfo6GiFovJwl0bf5ywbVsi2hNY4/6wn6/O+zxOWyVWs52GitoBF0htWTW1DrduzLdRhuf+rgODtq9qvJR9XQfTSh5SjA3kodTbJRgtZCZa/Ny6u427blxCKBzhlc2HB0oVZFOPZaxT0BY6RCNdlq6axdU3nsXff3kJjaePS+t8vy9UUmIO0dj0bm8/dosNu8VOKBIiYLC4GI9R1cZaRw1uh34pUhMmHtj2a3654xE2HN2IxXTiH6QW862HnpjD8BosFrOFby6+hfMnnkONowoTJurK3KyYtLTkwwzHMp1eHx1pLEwaLWQmW/xs8fSybY9+4EQ+mzWPFgVvoWvYbBYqnHYO7dNv5HBwdxtnrxj+5RVbPZdUCJtDuJxlVDtchMIhHpNP4/EltpZMmDhj3BzdRch5DXMAdJNxNGHWClwBg8rW6sV8n157KjvbdqWcrWkxW7hafJJPzrxEhRkWCIkWLfUwWshMtvgZimB4jdGOGc8HRSPokFp9l6Ekyjh1lFnBBL6+4LB9YxlnlYOPfPJ03t96jPffHd5azxK2Me34Gdgtdp5oeiZhpIhGbVkNVzVeSrm1PGHCjbbPhEnXyh5a4Moo5vuxXU+z9sjwBdNEbpR0kocUo0uieHU9jMIBE40zf1Ydv/zLe4ZjjnbMeD4oKkFPpb6LRnxTjHg/fHe3j8pKB45yC+0tI5sunytOn3cyu3c2c9DgVRPA2lJFb3+fYWXBocytn025tTxhwo22b3/nIe579xe64xgVuNLEOBQO8UTTMwOLrJovvdZRw7yGOcqNUkRcvXImFeV21m89Qnt3P64KG2fMqMdqM7N1d1vK4YBGIYTBUJjDLT2G582bWTdmFjNzRVEJeir1XcKhMOte2q3bFGPJ8ul0dfTR0+3jpT8XXoVCV1U0tj4CSbsw9XT7ON7uSRpeqCekiSxhu8XOtOophoudySJPtFrkGpqVP6f+tKQdhhSFhcVs5oZPzKbb28+W3a10eP3sPNDOgsYGvvfFs2jv7AOTiYaa8oS1V/RCCAFu+583E15/1aJJOb2fsUBRCToMj3qJTyAKBEI884etg8ROa4oRiUQgAu+9e4QcRDONOJdcPZeG8S56e/z87Y87kh7vdDkYX+s2rCzottfw1fk3DIpKSZVMC1wlqkW+s20X/pBf+caLjNV/3skrcXVRtAgUeaiD3v5AWkW14uPemz29CRdda5x2aqtGLnN6pCg6QdeiXpYsnz7gUrFYTIOaXOgRX+Cr0LDaTBzY3cZrf21KeXF3amM9FWXlhsI7f9wcJrpOynhOmRS4UrXISwtfIMSbO4av8UC0c5FGJmGGyRZdF8wa3RT9fFF0gq5hs1modkcr7qXSCLpQxRwgHIqwc3NqTTriXUyQv8qCmRS4UrXIS4tOr4+WDv1G03qkU1Qr0WLp5HFOrr2wuOLPNRIKeqw/6GpgKuAA7pJSPhPb9xNASikfjG3fCNwEBGPH/SWP806ZYmsErUeqLqLGOeM5/6LGQQlW+a4smE7kiapFXjr4AiH8wTB11WW0dqRW3jndMMNBZQa6+6mpdHBGYz3XrpqVcj30QiOZhX490Cal/IwQog7YIoR4A/gt0cbQ/wkghDgJ+EdgMVBGtGH0i1LKUQ/uTrWAl9VmJhgoXCs9Edpi6bkrZ2A2+EMeKyF/qhb56JJpd59UGVp3xW5P/RrphhmO1Xor+SSZoD8BPBm3HQScwJ3AxXGfnwWsjwm4TwixB5gHvJ27qWZGssShSqeN6aeO4/ChDtqbjUOcxjJms76V7qxy8LEr51IVV7tmrKNqkY8OuejukwpDm074kjSEjifT0rT56g40FknWJNoLIIRwERX226WU+4H9Qoh4Qa8C4lMNu4HqZBd3uyuwWvMvNKfPn8DGtft191lsVhx2KwFfYSUPxWPkcjl9/gROnX3yyE4mh0wk8RtDQ0Px+9RH6h4f+tN23e4+FeV2brx8bk6u0e8PGrZu06PcYcXnD1JfU87Zc07mhk/MxmIpTFfJSH2PSRdFhRCTgaeBB6SUjxoc1gXEz9gFJA5wBjyekUncWXDOZPr6/NHEoSGWepenj00bDo7IPPKJzW7GUWajp9tHdU05k2fUsuCcybS0ZN4RaCzT0OAq2nvTGKl79AVCrN+qHzSwfusRLj5rck5cFc2eXlo8xougNU47XT3+gcSgy5dNw9sbGHCVtLcX5ht0Pr5HowdEskXR8cALwM1SypcTHLoR+L4Qoozo4ulpQPJg6BFCC2VcdO4pPPnwO0VZtyUYCPPJ6+ditVk4ZWotHZ2pRw8oSptUuvvkwmWRKJSwrqqM73x+MX2+4CBfd4XDlvV1S4lk7y+3AW7gDiHEq7H/yoceJKU8BtwHrAXWAN+WUo58Z+Ik+H0hvN3FJ+YQTRSqcpdHa73bizYaVZEHNKHVI5f1TpK1bnNV2MdUbfFCJJkP/VbgVoN9dw7Zfgh4KGczywMVTjt2uyVaMncIVpsJMBVspMuUmXUDiVQKRTokitnOdY/MoXVX6mvKmTejbtRbtxULypSLEQ5FOHX+yby3JbUEndGmwmmjrydApcuBo8zKgd2t7Nx8BGeVg9PnT2DBOZMNQxQViqGMVI/MoaGEM6bW0a3cgzmjpAS91+vXtc4hGikS8IeYvXACctsxgsGxbalHIjBr9nisNhPvbTnRRNnb5WPj2v309flZumrWKM5QUUiMdMy2FkpYZrdS3EvbI0tJmXAVTjtVNcYFefbLFg7ubhvzYg7Q1xOgacdxdu9s1t1/oKmVQCBEIBCi09NHoMg6syjyw1jrkalIj5Ky0C0WE/YEC4bBYLjgFk2NatB0d/l4/fndHDnkGVYmWLliFIripKQEfcOavbTGVXErZmx2M007BrtitAJlyhWjUBQnJWOqFVKRriq3A5MptWNt9vS+Qs0Vo1Aoio+SEfRUi3SNBfz9ISKR1I4Vc09i7uKJuKqiDwFXlYN5iycZumKMeqsqFIrCp2RcLsmKdKWC1WoekQXT/r4gFouJUMhY1Z1VDiZMcXPW+dNxOKyDGno01DvZ19SSUm9VhUJRPJSMha71G82GYDCMZQQegVZbYjG3Ws309wVo2nGMP/zqbda9tBuLxRTNErVZsNmthveq9VY1QkXFKBSFS8lY6BDtN1pebuf9rUcG+o2eMita0e/g7ja83T7KK2z09gQMxwgZFGWsbaikPUGHcT1MZojoGPzBQGJ/S/xbgtFiZ6LeqnqEw+FBbfpUVIxCUXiUlKCbzWY+evkc5p01acA9oVmrZ68I0ev1Y3dYkhbwslrNlFVEKxtqQnn2iun88Tfv0N6SWgVJkxkuv+4Mtr59mGMfdNDbE8BZ5aC/P0Awg3Z4B5paWbJ8+sD96PVWTWSZb1izV7d5NqioGIWiUChJ00vrNxovcNpn5RX2pK6ZYDBMOBRm1uzxXPXFxSxdNQur1cJHLp+d8hwiYXj6d+9y/HAn00QD19x4Fh+7cm5GYg7RuHO9xU69ex1KogggFRWjUBQOJSnoyTh35QxmL5yQMHSwN5ap+fbaAwOfOavK0g4j7On2s3PzEXZuOUyVuxynQdW7ZNjs5owXOxNFAKmoGIWicFCCroPZbOb8jzRy+hnJu/3s29VCX29U8IKBkK5PPBUOxCzkbBduM0GLANJDRcUoFIVDSfnQjQgEQoP8zNr2khXTMVvM7NvVQo+Bldrj9fOH1Ztii6n+jMMaNUt46GKmxWpKukgK0QYXvV4/1e5h5eqTokUAxfvQNZJFxSgUirFDSQv60MiOSpedsnIbvv7goEiPT31+EY8++JahWPd6/Vm7JTRLeOhiZlmFlbfXHtBtn6d3fqakGxWjUCjGHiUt6EMjO3q6/fR0nxBmLdIjHI7k3Tk11BLWFjOBAYHv8vSx5a1DuhUWs7Wk042KUSgUY4+SFfR0arscaGrLOPokFWx2M2cum2a4PxwO89Zr+wbeJLSF14A/jKsqt5Z0/INEoVAUFkkFXQhhA1YDU4k2gL4LeA94GIgQbQb9NSllWAjxDFAHBIA+KeXF+Zl29qRT26Wnx0d5ZbRDUD4IBsL09wZwOPS/jqFvElqdlsY5J3H+RbOUJa1QKIDUHAnXA21SymXAxcBPgR8Dt8c+MwGXxY6dCSyVUq4Yy2IOiSM7hmKzmQn4DVJEc0Ai/3eiN4mjhzx5m5NCoSg8UhH0J4A74raDwCLgtdj2X4FVQojxQA3wZyHEOiHEx3M60xyTTm2XgD88LNLEbI66SkymaKGs2nGVsUbT6ZPI/61ixBUKRaokdblIKb0AQggX8CRwO/BfUkpN4bqBasAO/Ai4F6gF1gshNkop9XukAW53BVbryLsLGhpcAFz26TMwYWKbTrdzDZvdrFuK1lVdzo1fX0Z/fxBXlYOXn9vFxrX7056LzW7mY5fPwVGub6HXVJdT7S6n0zO8kW51TTmnTK3FptOFSbvHYkbdY3Gg7jF3pLQoKoSYDDwNPCClfFQIcU/cbhfQARwDHpRSBoFmIcQWQACGgu7xpFb3JJc0NLhoaTnRlnbxslN4b+sRgoHhol3htBtawF0dfRw50km1u5yWVi/vbT2S0XwC/jAffNCRcCFyyoxa3RjxyTNq6dDpmD70HosRdY/FgbrHzMfUI6nLJeZKeQH4ppRydezjLUKIFbGfLwbWAquAP8TOcQJzgPezmvUI8PbaA7piDtGszVQyKLNpnuGqSh4/fu7KGcOaWMxdPFHFiCsUikGkYqHfBriBO4QQmi/9VuA+IYSdqGg/KaUMCSEuEkK8CYSB26SUY7rnW6IFR5vdzJLl0zGbTUkzKCuc9oTWPIDJhG4XolTix1WMuEKhSIVUfOi3EhXwoSzXOfbruZjUSJHIstZCCVPJoLTZLEybVcfOLUcNrzVUzG12M6fOO3lgnKHlBxQKhSJdSjaxCBK3pTNKxTcS3KUXzuLY4S7amlNrcuEos7Fk+XQA1r20O2FjCdV8QqFQpEJJq0Gi0EWjVHwj69lsNnPl5xcxe+EEKp0OMBH9vwE9sZBDLWlIe6ho5QY2rNk7cGwqxygUCkVJW+iQ26JUWtndcy5I3v3I6XJgd1gSNpbQLPhkxygXjUKhACXoeVlwjK+Hkqgsrd8XSilpKNkxqvaKQqEAJegD5KsoVaI3gFAoktSHDxgeUxmLrlEoFApQgp53Er0BmM2JLXjtOKNjfL4gb722Ty2OKhQKQAn6iGH0BpCKD1/7ede2o4PKEAT84QGhj39gKBSK0kQJ+iiTig/fbI4mOe1vaiXgH+562bXtKPtkCz3dfpxVDk6fP4EF50xWVrtCUWKof/FjhGRhkYmSoAL+8ECnJW+Xj41r96uQRoWiBFGCXiCkU78doi6cQCCUxxkpFIqxhhL0AiGd+u2gaqUrFKWI8qEXEEMXUCtdDvr7A7r9ThN1QVIoFMWJEvQCQm8B9a3X9iUNe1QoFKWBEvQCJD4EUi/s8bRYlItCoSgtlKAXOHpW+4QJNUXfBUahUAxHCXqRkK/SBQqFonBQUS55JBAI0enpU+GDCoViRFAWeh5QDSkUCsVokFTQhRA2YDUwFXAAdwHvAQ8DEWAH8DUpZVgI8V3gEiAIfF1KuTE/0849YZ+PYGcn5vJywn19mMvLCXV2EgGs1dWE+/qwVldjdiRP7tEaUmhoDSkAlq6ala9bUCgUJU4qFvr1QJuU8jNCiDpgC/AucLuU8lUhxIPAZUKIg0T7jC4BJgN/BM7M07xzRiQUouWJx/Bu2Uywrc24mzNgravDuWAhdZd+kpDXqyvwiRpPq4YUCoUin6Qi6E8AT8ZtB4FFwGux7b8CHwEk8IKUMgIcEkJYhRANUsoWo4Hd7gqs1pEXt4YG18DP+x5aTcdLL57YaSDmAMG2NjpeepGu9WsJ9/twNNRTe9ZZTLvhc5gs0ftob+3B223ckKLMbqO2vjKt+Qb8Qbq7fLiqHNjsqXnJ4u+xWFH3WByoe8wdSdVBSukFEEK4iAr77cB/xYQboBuoBqqAtrhTtc8NBd3j6c1s1hmguVROmjmJ9q5oSnywu5vm9RvSH6uvHwBfcwtH//Is3vYOxl/3WcwOB4FACKfLuGlFvz+Qckhhpr74hgZX0YctqnssDtQ9Zj6mHimZe0KIycDTwANSykeFEPfE7XYBHUBX7Oehn48qg1wq7e0cbainfM58IkTwbn6HUEf2U+zesJ6+XbtwLlxIw1XXpNS0IhWUL16hUKRDKoui44EXgJullC/HPt4ihFghpXwVuBh4BdgD3COE+C9gEmCWUuo7k0eQliceG+RS8TW34FvzUs6vE2xvG7jOuZ/+eyC7xtPKF69QKNIlFQv9NsAN3CGEuCP22a3AfUIIO/A+8KSUMiSEWAu8QTS+/Wv5mHA6hH0+vFs2j+g1vVu2UP/JK7NuPJ2o/rlqDq1QjA6+QIhOr49qpwPHGDSoUvGh30pUwIeyXOfYO4E7s55Vjgh2dhJsb8/dgBYLhBInCQU97QQ7O7GPG5dV9qZW/zxZA2mFQpF/QuEwj6/Zw5amFtq7fNRWOVjQ2MDVK2diySK3JBwO86Mf/ZA9e3Zjs9n41rfuYNKkzOswFXWWi7W6GmttbVZjWOrrmXr3j3Cde15SMQewumuxVldndU1IXP9cVVJUKEaWx9fs4aVNH9LW5SMCtHX5eGnThzy+Zk9W465d+yp+v5//+Z9f85Wv3MJPf/qTrMYrakE3Oxw4FyzMaoxQayuHfngX3s3vpHS8c8GClJKPUuHclTOYu3girioHJhO4qhzMXTwxLV+8QqHIDl8gxJYm/WC9LU2t+LIo7bFt27ssWXIOAHPmzGXXrvczHgtKIPW/4aprAAYSh0xlDiL9+r5pI8IeT+IDTCastXU4FywYuF4uSKWBtEKhyC+dXh/tButZnu5+Or0+xrkrMhq7p6eHykrnwLbZbCYYDGK1ZibNRS/oGpFYwlC6Yp4Ma10dE275J+wNDTmzzIeiKikqFKNHtdNBbZWDNh1Rd7vKqHZm/u++srKS3t4T+TiRSCRjMYcid7nAibDFUC4XR+NwLlhI2aRJeRNzhUIxujhsFhY0NujuW9BYn1W0y9y583nzzfUA7NixnenTZ2Y8FhS5hZ7PsEWtrksuXSwKhWJscvXKqNBuaWrF092P21XGgsb6gc8z5fzzL+Dtt9/iK1+5gUgkwm23fTer8Ypa0IOdndGCW3nAVFZOw1XXDNRwUSgUxYvFbObaVY18avmMnMahm81m/vmfb8vBDGPj5WykMYi1uhpLTU1exg4c/pDjv/8t/uZmwr7c+uUVCsXYxGGzMM5dMSaTiqDILXQtbLHzlTV5Gb9r3et0rX1tkPtFWewKhWK0KGoLHWDcNddhmzgpP4PHIme0srotTzyWn+soFApFChS9oJssFiZ8dWTKyni3bFHuF4VCMWoUvaAD2Ny1WLIsAZAKWh0XhUKhGA1KQtDNDgeWuGysfJGrOi4KhWJs4g/5aeltwx/yj/ZUdCnaRVGtQ5EmsOHenrxf01xRgSmLLC+FQjE2CYVDPLXnWba17MTj68DtqGFew2yumHkJFnN2gRA7d+7g5z+/j5/+9BdZz7Po1GdohyJrbS3ls0Te4tHj8X9wiJYnHmPcNdcBgx8qKpNUoShcntrzLK9+uG5gu93nGdi+qvHSjMd95JHf8Pzzz1FWlpvSHkXnctFS/YNtbRCJEGxro/vN9PuGZop3yxZCfb00P/YIB75zGwe+/U0OfOc2mh97hEgK5XcVCsXYwh/ys61lp+6+7a07s3K/TJw4ie9//z8zPn8oRSXoo9GhaChBTzvN//vIsIeKCmtUKAqTTl83Hp9+7+H2/g46fZk3gF6x4u+yKsY1lFSbRC8B7pZSrhBCLAQeBHzAu8CtUsqwEOIZoA4IAH1SyotzNssUyXmHogywut30GtQ01trTKfeLQlE4VDtcuB01tPuGl9GuLauh2uEahVnpk9RCF0L8C/BLoCz20S+Ar0splwGdwLWxz2cCS6WUK0ZDzCE3HYqypVycRsigfroKa1QoCg+7xc68htm6++bWz8ZuGTvtIFNxuewFrojbniSl1JzS64GlQojxQA3wZyHEOiHEx3M8z5TIRYeijLHbqV6xknHXXGv4UFFhjQpFYXLFzEtYMWkpdWVuTJioK3OzYtJSrph5yWhPbRAmrfFDIoQQU4HHpJRnCyE2AP8qpXxNCPEA4AJuAz4N3AvUEhX686SUzYnGDQZDEas1t7VPIqEQ+1f/hvaNb+NrbYVwOKfjJ8IxroHas84iEg5z7Lm/Dtt/8scvYfqNN4zYfBQKRREkRL0AAAshSURBVG7xBf14+jtxl1XjsI6qZW7S/TADQRdEhTsEvA1UA/8C2KWUPbHj/wDcL6Vcm2jclpbu5BfPkLDPR6/cxZH7smu6mgk1K1eB2YR3yxaCnnas7tqB9nQjUbyrocFFS0vmCzWFgLrH4kDdY8Zj6gp6JsurlwA3SCmPCCHuB/4KrAJuBi4RQjiBOUB23U5zQN6eFknwbn2Xqf/2feo/eaWKQ1coFCNGJoK+G3hOCNELvCKlfA5ACHGREOJNIAzcJqVszeE8UyLs8xFob6fj5Rfp2b51RJKJ9NAWP+3jxmEfN25U5qBQKEqPlARdSnkAODv285+BP+sc8/WcziwNBmWHjpKIx6MWPxUKxWhQFKn/WnboWMG5YIFysSgUihGn4AV9LGSHaljcblyLFqvG0QpFkTLW6zMVvKCPhexQALOrisnfuh17Xd1oT0WhUOQYvaJ/2badDAaD/Md/fI+jR48SCPj53Oe+yNKly7OaZ8ELupYdOtq+83B3Fx/e8wPVW1ShKEKGunW1+kzAQHXVdHn++eeoqqrhjjv+nc7ODr7wheuyFvSCL841qtmhQ1BFuBSK4iORWzebtpMXXLCKG2/8ysC2xZK9fV3wgg7QcNU11Ky6EKvm7jDHbstmG5X5qN6iCkXxkMitm019poqKCioqKunt7eH227/JjTd+NZtpAkUi6CaLhXHXXEflvPnRD7R0/0BgVOajinApFMVDoqJ/2YYoHz9+jFtu+QoXXfQxPvKRj2Y8jkZRCDpEX4t6tm0d7WkAKg5doSgmErl1swlRbm9v4xvfuJmvfvUWPv7xy7KZ4gAFvyiqMVaiXUDFoSsUxYYWiqxXnylTfvvbX9Pd3c3DD/+Shx/+JQA/+tF9OBxlSc40JqXiXPkil8W5wj4fB75z24hHu9gnTyHc2zsqRbiMUAWPigN1j2OPTOLQx3pxrjGJ9lo0YhmjZjP2iZOY8q+3R9vMjeFkA4VCkRvMDseYrs9UND50iI92qQezGVNZ5q8uSQmH8X9wiNannhj4kpWYKxSK0aRoLHQ4Ee2ila21OJ20PfN01O/Vlrz4o6msjEgkAmmEHKo+oQqFYqxQVBa6hmYxWyoqGHfNdUz9t+/jOndp0vNsDePSEnNQIYoKhWLsUJSCrkefNO63Ya6upmrZckI93rTHVSGKCoVirFASgp4spDHc2UXPjm2EMgh7VCGKCkXpEAiE6PT0EQiERnsquhSVD92I5AW8IoQ8ntQHNJmw1tZlHYeqUCgKg3A4zIY1e9nf1Iq3y4ezysG0xnrOXTkDszlzuzgUCnH33XfxwQcHMZst3Hbbd5k4cVLG45WEoOcypNFaV8eEW/4Je0ODsswVihJhw5q9bN90eGDb2+Ub2F66albG465fvxaAn/98NZs3b+L++3/MD3/444zHS0nQhRBLgLullCuEEAuBBwEf8C5wq5QyLIT4LtEG0kHg61LKjRnPKg9olnTftq34mpsTHmupcRPq0LfYnQsWUjYp8yeoQqEoLAKBEPub9KPkDjS1smT5dGy2zBIJzz9/BefGAjaOHz+G251dP4Wk7wpCiH8BfgloQd2/ICrYy4BO4NqYyC8HlgDXAD/LalZ5QAtpnPejuzEnWsQ0m6mcN48p37uL6hUrB2LarXX11Ky6ULlYFIoSo9frx9ulH/3m7fbR6/VnNb7VauWuu77LT37yn1xwwd9lN1YKx+wFrgB+F9ueJKXcEPt5PXAZUAu8IKWMAIeEEFYhRIOUsiWr2eWBUE8v4a4u4wPCYbpefw2z3c746z875ltOKRSK/FLhtOOscuiKutPloMJpz/oat9/+PdraWvnylz/P73//BOXl5RmNk1TQpZR/FEJMjftonxBiuZTyNeATQCVQBcSvOHYD1UBCQXe7K7BaR7bmSchnx9FQj6858bOmb9tWar/8BSwOF0yqH6HZ5Y6GBtdoTyHvqHssDgrhHk+fP4GNa/cP+/y0+ROYMKEm6flG9/inP/2J48ePc9NNN1FebsJqtTB+fDWODI3HTBZFvwDcG3PFvE3Ul94FxM/YBXQkG8jj6c3g8tnR0OCifN4Z+JIskPpaWzm258MxXbfBiEIreJQJ6h6Lg0K5xwXnTKavz8+Bpla83T6cLgdTG+tZcM7kpPNPdI8LF57LD37wPT796WsIBoPcfPM/0dXlBxK7cYweEJkI+iXADVLKI0KI+4G/AseBe4QQ/wVMAsxSyuS59qPEiVKYmw1DGVXCkEKh0DCbzSxdNYsly6fT6/VT4bRnvBAaT3l5Of/+7z/MwQyjZCLou4HnhBC9wCtSyucAhBBrgTeILrR+LWczzAPxNV+OP/JbujesH3aMShhSKBRDsdksVLsz82+PBEVTDz1Vhr7+REIhWp54TLdw/WjWNM+GQnmNzQZ1j8WBuseMxyzueuiZMrRCo4pmUSgUhUrJC7rGWC9cr1AoFMkoieJcCoVCUQooQVcoFIoiQQm6QqFQFAlK0BUKhaJIUIKuUCgURcKoxqErFAqFIncoC12hUCiKBCXoCoVCUSQoQVcoFIoiQQm6QqFQFAlK0BUKhaJIUIKuUCgURYISdIVCoSgSirLaohBiCXC3lHKFEGIc8BDgBizAZ6WUe4UQNwI3AUHgLinlX0Zvxukz5B7PAB4kei9NwJeklOFCvkchhA1YDUwFHMBdwHvAw0AE2AF8LXaf3yXaSSsIfF1KuXE05pwuBvd4CLgfCBFt7/hZKeXxQv0u9e5RSvlMbN+1wC1SynNi20Vzj8CbjILuFJ2FHut1+kugLPbRPcAjUsrzgduBU4UQJwH/CJwHXAT8hxCiYIqg69zjd4F/k1IuJfoHdUmh3yNwPdAmpVwGXAz8FPgxcHvsMxNwmRBiIbAcWAJcA/xslOabCXr3eC9RkVsBPAV8s8C/S717JGaEfJHo90gR3uOo6E7RCTqwF7gibvs8YJIQ/7+9s3nRKYzD8DVMshkL0myo2ej2lQ1ZCFnY4A8glA0LmTKyUJislESYjQ0bSlOGBQs10cw0pGzGRrpjwwKFGWOUr/GxeM5bb+O1MKM5nWd+1+o95zyL5+r3dPd8vJ2ju8AuoB9YCzyw/dX2KPAcWDXdHZ0CEx2HgPmSmkgf6P5O9R2vA5111+PAamCguL4DbAbWA722f9l+CTRLWjitPZ08jRx32H5cXDcDX6h2Lf9wlLQAOAV01N3PypGScie7QLd9gxRoNdqAEdubScvZI8A8YLSuzRhQmS9CN3B8BnQBT4FW0uCpuuMn22OSWoAe0iynyXbtXRU1n8p6NnK0/RpA0jqgHThHXo6dwGXgEMmjRk6Oxykpd7IL9Aa8B24Vv28Da4CPpJlsjRbgwzT3639yAdhgeylwBThLBo6SFgN9wFXb14CfdY9rPpX2bOCIpO2kM5Fttt+SkSNp8rEEuAh0A8slnScjx6KOpeROloeiE7gPbCUNpo3AE+ARcFLSXNKe8zLSIVtVGSYNFoBXpOVepR0ltQK9QLvte8XtIUmbbPeT9ir7SMvW05LOAIuAWbbfldHnf6WRo6TdpEOzTbaHi6aVreVf6riieNYGdNvuKPaXc3IsJXdmQqAfBi5J2k9a7uy0PSKpCxgkrVKO2f5SZienyF6gW9I48A3YZ/tNxR2Pkv4h0Cmptj95EOiSNIe0vdRj+4ekQeAhyfNAKb2dHBMdZwMrgRfATUkAA7ZPVLiWjeq4xfbn+kYVH6+NHPdQQu7E63ODIAgyYSbsoQdBEMwIItCDIAgyIQI9CIIgEyLQgyAIMiECPQiCIBMi0IMgCDIhAj0IgiATfgPJlhTWCzI98gAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plot_observations(X, labels, title='Original clusters')\n",
"plot_observations(X, labels_hat, title='Predicted clusters')"
]
}
],
"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.6.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment