{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<!--NAVIGATION-->\n",
    "\n",
    "\n",
    "<a href=\"https://colab.research.google.com/github/saskeli/x/blob/master/bayes.ipynb\"><img align=\"left\" src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open in Colab\" title=\"Open and Execute in Google Colaboratory\"></a>\n",
    "\n",
    "|                                       -                                       |                                       -                                       |                                       -                                       |\n",
    "|-------------------------------------------------------------------------------|-------------------------------------------------------------------------------|-------------------------------------------------------------------------------|\n",
    "|  [Exercise 1 (blob classification)](<#Exercise-1-(blob-classification&#41;>)  | [Exercise 2 (plant classification)](<#Exercise-2-(plant-classification&#41;>) |  [Exercise 3 (word classification)](<#Exercise-3-(word-classification&#41;>)  |\n",
    "|       [Exercise 4 (spam detection)](<#Exercise-4-(spam-detection&#41;>)       |                                                                               |                                                                               |\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## ML: Naive Bayes classification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Classification* is one form of supervised learning. The aim is to annotate all data points with a label. Those points that have the same label belong to the same class. There can be two or more labels. For example, a lifeform can be classified (coarsely) with labels animal, plant, fungi, archaea, bacteria, protozoa, and chromista. The data points are observed to have certain features that can be used to predict their labels. For example, if it is has feathers, then it is most likely an animal.\n",
    "\n",
    "In supervised learning an algorithm is first given a training set of data points with their features and labels. Then the algorithm learns from these features and labels a (probabilistic) model, which can afterwards be used to predict the labels of previously unseen data.\n",
    "\n",
    "*Naive Bayes classification* is a fast and simple to understand classification method. Its speed is due to some simplifications we make about the underlying probability distributions, namely, the assumption about the independence of features. Yet, it can be quite powerful, especially when there are enough features in the data.\n",
    "\n",
    "Suppose we have for each label L a probability distribution. This distribution gives probability for each possible combination of features (a feature vector):\n",
    "\n",
    "$$P(features | L).$$\n",
    "\n",
    "The main idea in Bayesian classification is to reverse the direction of dependence: we want to predict the label based on the features:\n",
    "\n",
    "$$P(L | features)$$\n",
    "\n",
    "This is possible by [the Bayes theorem](https://en.wikipedia.org/wiki/Bayes%27_theorem):\n",
    "\n",
    "$$P(L | features) = \\frac{P(features | L)P(L)}{P(features)}.$$\n",
    "\n",
    "Let's assume we have to labels L1 and L2, and their associated distributions: $P(features | L1)$ and $P(features | L2)$. If we have a data point with \"features\", whose label we don't know, we can try to predict it using the ratio of posterior probabilities:\n",
    "\n",
    "$$\\frac{P(L1 | features)}{P(L2 | features)} = \\frac{P(features | L1)P(L1)}{P(features | L2)P(L2)}.$$\n",
    "\n",
    "If the ratio is greater than one, we label our data point with label L1, and if not, we give it label L2.\n",
    "The prior probabilities P(L1) and P(L2) of labels can be easily found out from the input data, as for each data point we also have its label. Same goes for the probabilities of features conditioned on the label."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first demonstrate naive Bayes classification using Gaussian distributions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.550760Z",
     "start_time": "2020-06-24T19:29:17.727236Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.721987Z",
     "start_time": "2020-06-24T19:29:18.551884Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJztnXeYFFXWxt87sad7hgyCMyAgBhAUBcEAJjCBioqorOIii6xiThgworu6BtQVEytmFHURcRVxRdHFdUWCoJJUEElKEpgc+3x/vPQ3qapzdZg+v+epZ6aruu493Qyn7j3RiAgURVGUpk9avAVQFEVRYoMqfEVRlBRBFb6iKEqKoApfURQlRVCFryiKkiKowlcURUkRVOEriqKkCKrwFUVRUgRV+IqiKClCRrwFqEubNm2kc+fO8RZDURQlqViyZMkOEWkb6H0JpfA7d+6MxYsXx1sMRVGUpMIY80sw71OTjqIoSoqgCl9RFCVFUIWvKIqSIiSUDV9RFCVeVFVVYdOmTSgvL4+3KLa4XC4UFBQgMzMzrPtV4SspSWUlMGsWMH8+0LIlMGoU0KNHvKVS4smmTZuQl5eHzp07wxgTb3EaISLYuXMnNm3ahC5duoQ1hip8JeXYsgU45hjg99+BoiIgIwN44gngqquAhx6Kt3RKvCgvL09YZQ8Axhi0bt0a27dvD3sMteErKccFFwCbNlHZA0B1NVBWBjz9NDBnTnxlU+JLoip7H5HKpwpfSSl++QVYvBioqWl8raQEeOSR2MukKLFCFb6SUmzcCGRn219fvz5moihKI+bOnYuDDjoI3bp1w4MPPhj18VXhKylF166AvyCM7t1jJ4ui1KWmpgZXXnklPvzwQ6xcuRJvvPEGVq5cGdU5VOErKcW++wKDBgFZWY2vud3ALbfEXiYlSZk+HejcGUhL48/p0yMa7uuvv0a3bt3QtWtXZGVl4cILL8Ts2bOjIqoPVfhKyvHaa8ChhwK5uYzQcbsBlwt48EHguOPiLZ2SFEyfDowbR6eQCH+OGxeR0t+8eTM6duz4/68LCgqwefPmaEj7/zgalmmM6QjgFQD7ABAAU0XkCSfnVJRAtGwJfP018MUXwJdfAs2aAcOHA+3axVsyJWmYOBEoLa1/rrSU5y+6KD4yBYHTcfjVAG4UkaXGmDwAS4wxH4tIdA1TihIixgADB/JQlJDZsCG080GQn5+PjRs3/v/rTZs2IT8/P+zxrHDUpCMiv4rI0r2/FwFYBSC6n0BRFCXWdOoU2vkgOPLII/Hjjz/i559/RmVlJWbMmIGzzjor7PGsiJkN3xjTGcDhABbGak5FURRH+Mtf6Pypi9vN82GSkZGBKVOm4NRTT0X37t1x/vnn45BDDolQ0AZzRHU0G4wxuQBmArhORAobXBsHYBwAdIrg6agoihIzfHb6iRNpxunUico+Qvv9kCFDMGTIkCgIaI3jCt8Ykwkq++ki8k7D6yIyFcBUAOjbt684LY+iKEpUuOiihHbQWuGoScew8MM0AKtEZLKTcymKoij+cdqGfyyAUQBOMsYs23s4t19REp7KSuDVV4HjjweOPBKYNAmIoPifoigh4KhJR0S+AJDY5eeUmFFeDpx4IvDddyxUBgDff8/SxF99BRxwQHzlU5SmjmbaKjHjySeB5ctrlT3Ah8Du3cDo0XETS1FSBlX4Ssx49lnWnW+I1wssWQL89lvsZVKUVEIVvhIz9uyxv5aVxZW+oqQyY8aMQbt27dCzZ09HxleFr8SMPn3sr4mw4GCiIBogrMSB0aNHY+7cuY6NrwpfiRn33NM4ORHgueuvZ8XKeFJTw45XHToA6elA+/bscWvVHUtRolwdGQBw3HHHoVWrVpEPZIMqfCVmHH008PLLQIsWrFDZrBmV/LhxfBjEm1GjgLvvpi9BBNi6Fbj3XmDkyHhLpiQaDlRHjgmq8JWYct55wLZtwLvvAm+8AWzeDDz2GFdJ8WT5cmD2bOuKtx98ACxdGh+5lMTEX3XkRCYmtXQUpS6ZmYzHTyRmzwYqKqyvlZfzAXXEEbGVSUlcHKiOHBNU4SsKGBpq56gViY8df9cu4J13+LN/f2DAANbxV+JPp04041idT2TUpKMoAIYOtXcau91AlMuSB+Sll9h/99prgdtvB04/nTuMnTtjK4dijQPVkQEAI0eOxNFHH401a9agoKAA06ZNi2zABugKX2lSVFUBCxYwm/eoo4C2bYO778gjgRNOAObPr58clpPDlXW/fo6Ia8mSJcD48TQl+aiqAlasAM4/H/jkk9jJoljjUHVkvPHGG5EL5wdd4StNhnffZV/ac84BLr6Y/wnHjw/eHDNrFnDddYweyswE8vKAq68G3nsvtqaURx6x9idUVbEH788/x04WxZ6LLgLWr6c5cP365KiUrCt8pUnw9df8D9cwcuLll4HcXMbTByIrC/jrX4H77weKi3lfPKKHvvuOSsSK7Gzgxx+BLl1iK5PSNNAVvtIkuO8+6zo9paXAU081fhA0pLyc4aHdu3NncPnlwKpVzsgaCH/KvKoKKCiInSyphiR4inWk8ukKX4kpu3ezFHJ2NnDssVxVR4NFi+yjbDIygJ9+Ag491Pp6eTkwcCBt5L6HxptvMlTz/fdjH0J6/fXAp582fkilpbGEdI8esZWnqVNczB1gQYELRUU70axZa+TnG+TlxVuy+ogIdu7cCVcEKemq8JWYIALcdhtr39dV8s8+G51M1pYtmRlrRWUl4C9b/cUXgZUr6+8QvF4q3EsuoVMuljb8k06i0p88mSv66mqalzwehmkq0aOsjBngP/0E5OQU4J57NqFbt+345RegTRvrUiDxxOVyoSCCLZ5JpC1M3759ZfHixfEWQwmToiJgxgzaoLt0YamCNm147YEHaBtvuGp1u4E5c9gBKxKefBK49dbG4xvDom2LFtnfe/jhwLJl1tdycxn107t3ZPKFw4oVwCuvsCPY8cczQicnJ/ZyNGWefhq4+WZrk1/r1iyzkeHwstjr5c43Nzf8Ha8xZomI9A34RhFx9ABwGoA1AH4CcKu/9/bp00eU5OSrr0SaNRPxeEQAkZwcEbdb5L33RCorRZo353mr4/jjI5+/ooLj+OYHRFwukZYtRVav9n9vt272sjVvLrJgQfhybdkisny5SHFx+GMoztGnj/2/fbNmIl9+6dzcNTUiDz0k0rq1SFaWSHa2yMUXi+zcGfpYABZLEPrY6Sbm6QCeAnA6gB4ARhpj1ALZxKioAE47DSgsrO1mVVbGVdOFF3L1XF1tf3806tRkZQHz5gHPP0+TSL9+NCGtWQMcdJD/e089lWGYVlRWAocd1vj8nDnAcccB+fn0Rbz7bn0fwvr19At07cqfbdsyiaqqKuyPCIC7qEmTWJ2xbVvWJvr228jGTGXq5jqEcz0SrruORQN37uTfWUUF8NZbNDE5Nm8wT4VwDwBHA/iozuvbANxm935d4Scnb70lkptrvUpyuUTuvJMrGLuVVH5+fOX/+WeRvLzGcrndIrfd1vj999zDa3Xf6/GI3Hwzr+/eLdKunUhaWuPxLrnEXo4lS3i9f3+RMWNEvv22/vWiIpGDD+Z36hszLY3jfvZZ1L6OlGLiRK6srf4u3W5+507w66/28+bmirzySmjjIRFW+ADyAWys83rT3nP/jzFmnDFmsTFm8fbt2x0WR3GC9evtVyTl5cDGjVwFW8W05+SwrGw86dyZGbbdu1OevDw6SG+8kX6HumzYADz4YGObb0kJMGUK8MMPdAIXFzeOpS8t5Qpu06bGMjz1FDN6X3sNWLiQ+QP9+wPPPcdzjzzCHcIvv9T/rn3O5Usv1aYt4XD11fy3bvi36XYDEybQru4En3xiv6ssLmaUmBPEPUpHRKYCmArQaRtncZQw6NaNirKoqPG1nBzgkEOAu+6imaW4uFZheTzAgQcCN93Ec2+8Abz6KpXYhRfS6evxxOYz9OnDSJ21a2maOugg6wiNt9+2T4qqquJnWLDAPu4/MxP44gt+Ph8bN9Z+Bz5qamgWu/xyfgeVlTSL2Sn1bduA1av50FKCZ599GCY8dix/ZmQwZPiOOxgp5RTp6f4jv+weBpHitMLfDKBjndcFe88pTYgzzmDhMSuFn5YGjB7NaJ3VqxmGOWsWHwSXXsrs2IoK4JhjuFPw+QAWLWJs9Ndf10b6hIsII15+/x3o2dN/iOb++/sfq7jY3g5fXQ18/LH/lXZlZeMHyfTp9g8RoPY78Ud6unXimRKYAw4APv+ctvSiIia2OR2Zc8op9n9HubksDeIIwdh9wj3AB8o6AF0AZAFYDuAQu/erDT95Wb5cpE2bWlu+x8Pf580LfO+111rbMzMzRUaNikyupUtFDjiA8jRvTvv32LGM6gmHTz+191cEe+zZ0/jzRzIeQB9EeXlk35USW+6/v7EvyOUS6ddPpKoqtLEQpA3f0eeYiFQbY64C8BGAdAAviMgKJ+dU4sOhh9I2PWsW7dgdOwIjRgRnA33xRftiYTNmAC+8EN6Ka/Nmxq833HlMn04TyBFHcMU9ZAh3GMEkV51wAs1Q330XXsRNVhawYwcLtAFMFvv9d56vrAx9PIA7hjvvpCnCaSorgX/8g/Hru3czouSOO+KTp5DsTJzIHeU999CU2Lw58Oc/8/t0bIcRzFMhVoeu8FMPr1fEGP+r1xdeCG/sW26xj4QARNLTObfHwxj+0tLgxt21S6RXr/BW4llZIjt2cJwXXuCKzp+MVkdaGsfJzRVp0UJk8mR+j05TWSkycGD9VakvSmjOHOfnV+xBgkTpKIpfjOGK2R/jx3MVHCrz5tm3LQToGBWhjXzhQmZcBkOLFvRLhOpYS0vjirh1a/oUrrySjtqGMqal0cdht+No25b5BYsWcady/fXB7U4qK+kT+eYb/z4DO2bMYM5EXYe0L0roj3+MT1cwJTRU4StxZ9Ik/1vYtLTwwtRatgz+veXlNC0Fa1Y5/fTQtt0uF53FL7zA13//u/1cmZlMIHv6aSp+X7q9x8Nt/5w5DCU9+ODgHzpPPcUHxcknM2GsQwcWhwuFqVPtHcjl5YxyURIbVfhK3Dn/fK587Sgt5Uo2VP7859DCOr3e4HcS3buz7aG/2jbp6WxT2KcP2xSuWsXMW4AhoHYr4qwsPuQuv5w7gVtu4Qr6b39jHkCozdRfeokx5YWFPIqL+X2OHAn85z/Bj2MVheUjLY3jKglOMHafWB1qw09dXn3VPvolL0/k3XdDH7O6WmTo0Pr1dfwdbndo0TtVVSJ3310/89V3GMOooI0bre8dN44+BCs5PB6RhQtD/7xWeL0iHTrYf+YBA4Ify59PxOUKrwaMEh2gNnwlmTjvPEabNLRFp6XRZj50aOhjpqfTbPHcc+xv260bf1qtyl0uYMyY0KoVZmQwwqK4mKWMfeV0s7KYNfvVV/WblRQXM3O3c2cmcIlFvH5aGu858shQP60127f737X4qyLakGuusW707nYDl13mP79BSQy0PLKSMKxZwxDJbduo+L1ehnfOmRO9ln5eL7MqZ8xgWGVNDRVWv36cJ4LeEqipAbZsoRmpofIrKWGphLVrazNqjaHSz8ykLLm5DNdcsKDW9BMpRUV0EtuFkLZoAezaFfx433/PZLkff6Tc1dV0Pj/wAB+wSnwItjxy3EsrKIqPgw5iI4ovv2TWbbduVMTRbD6SlkbH6Y03splIVRUdsEcdFfk86el8QFkxdSqwbl398gm+tZYx/Jxnnklbe7S6gAGsC3T88eyg1TAyJyuLDV5CoWdPYPlyNlLftYv/ZrEqf6FEjq7wFSUG9OxJB6wdxnCnMWAA8N570VX6P/3EB0pJSW1kkMsFtG8PLFmippimQLArfLXhK0oMCFTnxpcP8J//sCF7NOnWjaaYa65hZufBBzMzd9kyVfaphq7wFSUGXHUVzTrBlGNo0YKO1lj20VWSG13hK0oCcdNNwTuECwvDr6ujKP5Qha9EnaIiNgMZPJilk99+23+Lw1Sgc2dG3xx+eGDF36JFeDb80lLgs884jz4wFCs0SkeJKlu20EG4a1dtzZXPPgMef5yRIrGo6JhIlJTwgffzz7Slf/EFv5tHH2VvgIa2fbcbuOGG4M05JSXsafv228w38JV7MIblG0KNwlGaOMFkZ8Xq0Ezb5Kamhj1XrTIxc3JEHngg3hLGlgULRJo1q80gzs1l9u1XXzEL+IIL+L1kZDDr1u0WGTEiuFroXq/IXXfxHqtMX1/m8Ny5zn9OJf4gyExbddoqUWPcONZKt6OggO386rJiBbM9W7UCTj216ewACgv5ea3qzzRvzp2Q2826+u+/z/NnnAH06hXc+JMmsbaOXStFH337hpZNG4i1a1lts00bYOBATbZKFDTxSokpmzez8bY/du+u/b2wEDj7bJYf8PX3NIYZsKef7qyswVJdHbj3qB0zZtgXR6upAf75T5pbevUKXsn7KC8HHn44sLIHGHoZKZ9/ztLRixdz7+Dr++rxADNnMndASQ4cc9oaYx42xqw2xnxrjJlljGnh1FxK/Jk3L3Cp3kMOqf39gguYUVtWxhozRUV8CJx3HjtmxYLCQjqXhw5luYBPPqFCe/tt1ujPyqrtvbt9e2hjr1ljr5CLi1maIFzWrGHGcDAE03HMH//6Fx/AixbVZgZXV9N3sG0brzXctSmJi5NROh8D6CkihwL4AcBtDs6lxJm0NP8r4fR0FhoD6MD87DPr5iSVlcBjj9U/V1TEhKRFi+qXBxDhOFOnAh9+GHwk0PbtwBVXsF7+Ndewhs7rrwPDhtEEMno0FbIIZZw+ncXMCguDGx+gg7Zhs3IfHk9ktXLy8oKL58/KYlnlcPF6WaLZX9JYZSXw5JPhz6HEmGAM/ZEeAM4BMD3Q+9Rpm7xs3eq/Vd+f/1z73vfeo/PS7r1HHMH31dSI3H47HZs+5+c++9ARuXatSNeuPOd2s4Ryu3YiS5b4l3PzZo4RqK2ildP50UeD/z527bIvy5yXJ1JUFPp3XJcePfzL63aLdOtGOcJl5crgSksfc0xkn0WJHCRYeeQxAD60umCMGWeMWWyMWbw91H2zkjC0a0c7b8NCWpmZXDVPmVJ7rn17/+3w9t2XP++7j+GcZWW1jTu2bgXOPRc49lgWWCsupumkqIgmhpNO8t+oY8IENhEPNVahrIwr/WBp0YI1cTye2pW+x8PV+QcfRG5qefVVjlHXjJaRwRj/3r1ZvfKbbyhHuNTUBOe/aN8+/DmUGBPMU8HuADAPwPcWx7A675kIYBb2lnHwd+gKP7nxekVeeklk//3Z3LplS5Fbb23cHNzrFenSxXq1mJUl8tFHImVl9g1RjPHfPOTpp+3lC7VhuNXOIxR27aI8110n8uyzInv2hD6GHevWiVx+OXc6vXqJ/P3vwTdiD4bqapE2bfx/Jx6PyLx50ZtTCQ8EucKPKEpHRAb7u26MGQ3gDACD9gqlNGGMoc04kN3YF41z9NGNS/Z6vcD8+Vzl260uRex3CCUlrABpd1+4Gag5OcAf/hD6fS1a0F/gBF26AM8848zYAP0ukyezVaSVHT8nB7j4Yu6qlOTAySid0wBMAHCWiAQRQKakEitXWpcYqK6mGcfr9e+YtHsYZGXZ16RPSwMOO8y/XJmZjeXKyGDc+dix/u9tiowaRfNR1678HnwdyE45hRE8zzyjRd6SCSdt+FMA5AH42BizzBjzrINzKUnGm2/ahy1mZgKrV7NRuBU5ObUlBBqSns4wSjsuuMA+pDErC3j+eeDFF1lG2BjGm48cyRj05s3tx23KDB/Omvo7dnAHtWsX8NFHwKBBquyTDccSr0Skm1NjK8lPIEWRlsZErgEDmGjkM8W43XwQ3HknFXFVVW04Zloaq1J26mQ95mOP0RHc0IyUkUGzywMP1DqdL7yQc/pWtamOMan7wGtK6J+yEhf+8Af71nhVVay02asXSy9cdRWVfJ8+LDr2xRdsB3jMMfUfHCK8PnNm4zG3bQNuu816V5GZSb9DQ3mysuKj7H/7DfjrX5kMNmkSs5gVJSoE49mN1aFROqlDebnIoYc2jppxu0Xuuy/w/bNn28eIN2vG8evy3HMc2+r9aWmMokkE5s6tXxAtO5s5AO+8E2/JlEQGCRaHryj1yM7mSv2KKxibDtAx+OyzwB13BL7/H/+gPdmOTz+t/7qkxD4T1+utX+cnXhQW0l5eWlrb7LyighEyF10UenkHRWmIKnwlbuTl0a5eWEilu3Yto0KCYc8e+2siTMiqy8CB9o7e3FxGncSbt97yfz2UxC9FsUIVvpIQhBrtMXQoo3WsqKpijH9d+vbl0bD8si/kcvjw0OZ3gk2b7HctZWXMLFaUSFCFryQlY8dS4Td8UOTksOJmQQFf79nDzk9nnsnonQEDqPSbNWO8/XHHsWpnOC0Fo0337rXmrYbk5gKHHhpbeZSmhzZAUZKWH35gaOaqVYy0qaxkDP7jj1OBr1nDmjtlZbSLG8OwzuHDgauvBvLzgQ4d4v0paqmoYNKYla2+RQvuAOwim5TURhugKE2eAw9kGYX165kUdMAB9WPFhw8Hfv+9tlCaCE0mM2dyF9A34H+P2JKdTWfz4MF8QFVU8Fx2NjB3rip7JXJU4StJT+fOPOqyciXr7lttYEtKgCeeoJkn0ejZkyv5uXOZ3dq5M/0VgZrLKEowqMJXmiS//upfSSZyMlNGBvvbpjolJWwFuW4dS12cd559UxklOFThK02S7t2tO2oBzJ494ojYyqOExn//CwwZwnDd4mI6ra+5hjufo46Kt3TJi0bpKE2Sffdlv9WGYZgAo3Nuvjn2MinBUVREZe9regPw5549wGmn+U+4U/yjCl9psviKr+Xk0BSQl0fH54svsiuUkpi8+aZ9v4OaGjaZV8JDTTpKkyUvD5g3D/j+e2DhQkbwDBmiduBEZ80a+1V8cTEbzCvhoQpfafL07MlDSQ72358PZavKph4Pay4p4aEmHUVpArzxBh3VWVlsKj5pkr3TOtEZOdK+LHVaGpvYKOHhuMI3xtxojBFjTBun51KUVGTSJOCyy9glrKoK2LoVePBB4OST7W3hiUzz5sDs2VzN+8xvHg8jdf71L/5UwsNRk44xpiOAUwBscHIeRUlVtm9npy5fOWUfZWXAN98AH36YnDH9J53EBLTXX6fN/qCDuPLXrluR4bQN/zGwkflsh+dRlJRk7lz7ss/FxSypnIwKH2D9oPHj4y1F08Ixk44xZhiAzSKy3Kk5FCXVqa62Lh/ho6oq8Bjl5cDUqUxo6t2b/YK3bo2ejEriENEK3xgzD0B7i0sTAdwOmnMCjTEOwDgA6GTXfVpRHGbPHmDDBjo827aNtzTBM3iwvZ0+N5flCPxRWsrmMHVDIVevBp56CvjqKxaoU5oOEa3wRWSwiPRseABYB6ALgOXGmPUACgAsNcY0ejiIyFQR6Ssifdsm0/80BQBbA06YALRrRwUzaBDryycLpaXA6NFU9AMGsDzx0KHJ006wY0d2CWuYW5CVxZ4A557r//7HHmOhubpx7xUV/HcdMyb68irxJSb18Pcq/b4issPf+7QefnJRVMSaNBs31g8BdLuZLZkMtmPfA6qu0zMzE9hvPyrCZKhS6fUCf/0r8OijdNYaA4wYwYqgLVv6v3e//bizsSI7m9fatYu+zEp0CbYevsbhK2HzzDOsOtkw3ru0lB2pvN74yBUsS5fSbNEwwqWqCvjtN+Ddd+MjV6ikpbHx+44djGzZvRt45ZXAyh7gQ9uOzEz/vYOV5CMmCl9EOgda3SvJx8svc0VpRUkJsDzB3fWff06npxXFxYyASSbS09mf16pgnB3+msCkpXEHoDQddIWvhE1lpf21tDT/1xOBnBz7kMa0NPv+sk2Je++1bgbvdtM3kwi9fpXooQpfCZuzzrJXCF4vcNhhsZUnVM4+297s5HIBF18cW3niwdFHA6+9RvNPXl5tc/drrgFuvz3e0inRRhW+EjY33MCVoDH1z7vdwD33UHEkMu3bU86GES4eD+u1JFrPW6c491xg2zZgzhz2+/3tN2bvNvx3VZIfVfhK2OTn0+k5YABX+jk5jGF/9FHgxhvjLV1w3HIL2+gNGEDZe/cGnnsOmDYt3pLFlowMfgeDB2v5gqZMTMIyg0XDMpOXXbvoqO3Qgc5DRVFiR7BhmVoPX4kKLVsGFwaoKEr8UJOOoihKiqAKX1EUJUVQha8oipIiqMJXFEVJEVThK4qipAiq8BVFUVIEVfiKoigpgip8RVGUFEEVvqIoSoqgCl9RFCVFcFThG2OuNsasNsasMMY85ORcShPkm2+AN94AFiwAEqjmk6IkK47V0jHGnAhgGIDDRKTCGKOdMZOV7dvZO69jR7ZUcppff2Un8TVrWIlNBGjVCvjgA6BnT+fnV5QmipMr/CsAPCgiFQAgItscnEtxgp07gTPPpKI/4QSgoAAYPtzZRqci7Cz+7bdsjltUxH6DGzYAxx/vvwmroih+cVLhHwhgoDFmoTHmc2PMkQ7OpUSbmhoq2H//m13KCwv58/33qZCdMrF89hmwcSPnb0hFBdszKYoSFhEpfGPMPGPM9xbHMNBc1ArAUQBuBvCWMY176BhjxhljFhtjFm/fvj0ScZRo8tFHwC+/NG5MW1lJU8vnnzsz77Jl9s1wS0qA//3PmXkVJQWISOGLyGAR6WlxzAawCcA7Qr4G4AXQyAAsIlNFpK+I9G3btm0k4ijR5NNPaUqxoqTEOYXftq19o9zMTGDffZ2ZV1FSACdNOu8COBEAjDEHAsgCsMPB+ZRokpvLvndWZGay8asTnH22vbkoIwMYM8aZeRUlBXBS4b8AoKsx5nsAMwD8URKpn6Lin5EjqditSEsDRoxwZt7cXOD119lZ3Dd/Whpf33cfcOCBzsyrKCmAY2GZIlIJ4GKnxlcc5qCDgPHjgWefpQnHh8fDDuX77efc3GedxSidJ59kLP7++wNXXQUccYRzcypKCqBNzBV7RIB33wX+9jdg/XqgWzfg1luBM86It2SKotRBm5grkWMMcM45PBRFSXq0lk447N4NPPoocNJJdDL+61+A1xtvqRRFUfyiK/xQWbsWOOoo2rXLynhu3jxmos5C372pAAAdo0lEQVSezVIAiqIoCYiu8EPloouA33+vVfYAlf/8+cBLL8VNLEVRlECowg+FjRuB5cutzTelpcATTzgz7+bNwOWXMympVSvg4ouBn35yZi5FUZosatIJhe3bmQVaXm59fYcDeWUbNwKHH86CZdXVPDdjBvDee8CXX2r1SEVRgkZX+KFwwAFAVZX1NWOciRO/7TY6iX3KHmBhsaIi4Oqroz+foihNFlX4oZCXB1x6KZCT0/haTg5wxx3Rn/Odd6wrRwLAF1/QlJToeL0sxnbLLcD99wM//BBviRQlJVGTTqg89hhX3DNn1hb5EgGmTmX0TrSpu7JviDH2O45EYfduRjCtXctibBkZwF/+wszZ++8HZs3iQy0zkw7x005jKQVFUaKOZtqGy4YNtKF7PMDgwdar/mhw4omsEW/FgQeyVHEic+657FTVsOSx2w20bEnfhK8qp8cDHHYYw1yd+j4VpQkSbKatKvxEZ+FCJng1NN3k5ABvv81WgInKzp3skmXn5DbGujKmxwM89BBwxRV8jz/Ky7nb+u47duYaOZKRTIqSQgSr8HXvnOj0788uUwceCLhcVPSdOgHTpyeWsq+oANatownHx6ZN9rXtAfsyyCUlwIQJwAMP+J/Tp+Qvv5z1fiZM4APm3XdDlx+gr+GZZ1iszeXiz2ef1SxqpcmgK/xkYvNmOnA7dgy88o0V1dXAXXexsqUIXw8eDDz/PJV9fr79Cj8QbjewbZt17f3qair3rVut71uzhtdDYcwY4K23GlcHPf984IUXQhtLUWKIrvCbIvn5XN0nirIHgHHjmHBWXExFWVHBiJz+/blKHjrU/yrfHxkZNGlZ8e9/20co1dQA//hHaHN9/z3zG+oqe4CvZ8wAVqwIbTxFSUBU4UeLPXuAKVO4Spw0iQlTTZ1Nm4A33miseKurab9//XVg2jSgVy82NgEYjeNyAQMGcCUeCLsmLOvW2UcoVVQAK1cG/zkARgrZ9dKtqqKfQFGSHA3LjAZffw2cfDIVXWkpkJ0NPPgg8NxzwKhR8ZbOOebPt2+DWFJCJTp2LLBoEd/72Wc0kYwYwQYql10GvPqqfehpWpp9qGvXrnwYWJmLsrOBHj1C+yyVlfa2+pqaxA9/VZQgcGyFb4zpbYz5yhizzBiz2BjTz6m54kpVFTBkCFBYWLvSrahgcbU//xn4+ef4yuckLpd/85JvBW8MI40mTWLyVdeurCr6wgusTdS+feMHh9tNv4DdCv+UU+x3COnpfJiEwumn24/ndjM/QFGSHCdNOg8BuFdEegO4a+/rpseHH9qbAsKxJScTp55qvzrPzQVGjw48Ro8ewOrVwHXXAa1bc3Xevz9LTV/sp0NmRgbw8ce8x2cuysnhMX166A7bY44B+vblQ6wuLhfPH3NMaOMpSgLipMIXAM32/t4cwBYH54ofGzbYb/crKxOjqqUI8M9/0m7epQswbBjw1VeRj9usGcMhG66M3W6gXz+umoOheXPg4YdZfK68nLINHhz4vl696Ed45hm2Xnz4Yb4+++zQP4sxfHj/6U+UPyeHP8eM4flEcpQrSpg4FpZpjOkO4CMABnywHCMiv1i8bxyAcQDQqVOnPr/80ugtic28ecwmLSpqfM3lYn2diRNjL5cPEZo36kagGEOF9tRT/lfhS5bQ7p6TwzaHHTpYv+/DD4F776WjtHVrlk245hp7c0yiU1HBngetWnHHoSgJTkwybY0x8wC0t7g0EcAgAJ+LyExjzPkAxomI32VbUsbhe71cNW/c2DiRyO1mNMk++8RHNoDlH045pXG4IUBF/ttvXKnXpbSUjcoXLqTJJj2dn+3OO4Hbb3de5sWLmWn73Xd07l5/Pc1HiqJYEqzCh4g4cgDYg9oHigFQGOiePn36SFLy448iBQUieXkiGRkiubk8Pv443pKJjB0rYowIVXb9IzdXZPr0xvdccomIy9X4/W63yPvvOyvvSy9xnrS02nk9HpGbbnJ2XkVJYgAsliD0spM2/C0Ajt/7+0kAfnRwrvjSrRuwfj1j0v/6V6bjb90anB3aaXbtsi9h4KurX5c9e5htahXuWFrKz+cUu3ezfk5paf0QyZIS4OmngW+/dW5uRUkBnIzDvwzAE8aYDADl2Gunb7KkpzOrNJHq2wDMD5g719qkAzSOPtm40T6+HfBfy76oiA+Y9u2Dz67dsAF45RXOW1lpXxq5ooI9gydPDm5cRVEa4ZjCF5EvAPRxanwlSC66CLj7buYF1F01Z2YCxx7LSJe6tG9vH2YKsLxDQ7ZuZQGzOXMYLpmeTrv7XXfxdzv+8Q86d71ezpmV5T/Edds2+7EURQmIllZIFsrK/DdDsSM3F/jf/4A2beqf93qtFWybNjRFWUXYeDzATTfVP1daWlvRs7KSr4uKgEceoXnGjlWrgGuv5U7CJ4O/B01uLjBokP11RVECogo/0Xn7bfoI8vIY9TNiBGPNQ2H+/MYmnZoanrcKGX3xRaBzZ84J0MzidrNq5EUX1X/v668zfr7hw6i0lGUTNm+2lunpp4MvV5CWRoV/wQXBvT8QIsyA3rDB3r+hKE0QVfiJzPPPM05+7draei6zZgF9+lDJBssDD1jb8MvK6GD+/Xfax//yF47fogWrQ774ImP4r70W+M9/WAqhYQLSrFn2/oGMDD5UrPjxR/81dDIzmZDldrML1pdfBldsLRAzZ7K8dM+ewEEHsTH9vHmRj6soSYAWT0tUKitpPmlYibKmhpE0f/87a9MEw4YN/ucpKKCSLSuj2cbjoaIePpyHP3ztCa0wxt55e9hhnMPKjJOTQ0du8+b0GRx8sH8ZguVf/2Ixu7Ky2nNr1zLz+OOPtXyC0uTRFX6ismyZffXGigqGTgaLv8SvykoqwJISzldURCfsySfz4eKPTz7xX6Khqso+YeqKK6wdumlpVPLnnEObfbSUPcAHaF1l76O0FLjttujNoygJiir8RCVQ7ZZQarvccIO1OcQXUdMQEe4iApk6xo+3d7Smp9OU1Ly59fXOnVnuwe3mjiI9nT6Djh3ZQCXatWsKC5n1bMeXX0Z3PkVJQNSkk6gcfrh9rXmXC/jDH4If6+qrgQULqEh94Zm5ubST79plfU9VFWPu7Vbov/7q31TkcrECpj/OOovjzJzJXUWvXixD7C+UM1wC1fVJ1ro/ihICTX+Fv3YtcOGFVHBuN5XM99/HW6rAZGTQTt9wZZ6RwQJlV15pfd+yZaxg+eijtZU609NZLXPePDpgx4+nkr3iCnsbe3q6/xLDgaJb7B5WDWnWDLj0Ula7HDrUGWUP0C8wcKD1ziE9HTjvPGfmVZREIpj6C7E6ol5L56efRJo3r1+XxRjWZlm6NPD9K1aIXHGFyHHHiYwbJ/Ldd9GVLxg++EDk0EMpt8slMnq0yG+/NX5fRYXImWeyDk1GhkhWFt9/xRUiXq/12OvXi+TkWNfZAURathR54AGRmprG93q9It26Wd+XliYyalR0Pv+WLSITJ4r07y9y2mkis2ZZyxMMK1aINGsmkp5eK2tGhkibNiIbNkRHXkWJAwiylk7clXzdI+oK/7zz6iv7useAAf7vffllKsOMDL4/PZ2vp02LrozBEkjJ3XyzvfJOTxcZPpwPwIa8/jrvy862vtftFhk/3nrOuXMbz2kMlarVXKHyzTccq65sHo/IsGEi1dXhjbl2LR+abdqItGsnctVVIps3Ry6rosSR1FL4Xq/Ijh0ixcX1z9spMd/KrqTEerytW62rRQI8v2VLeHI6RWUlK1/afVbfqrtFC5F16xrfv2WLyIgRtQ+3hkd2Nt9TUyPy2WessPnNN7z3k09EDj+c46eni5x+usiqVZF/Jq9X5IADrOXxeKyrfCpKihKswk9+G/5rrwGdOgH77gu0bEmn39q1vBYorNDu+ptv+o8Sef318GR1il27Apdd8HoZqXLXXY2vdejAEgd2Y6SnA9Om8Xs+80zWzTn2WLb+O+QQYOlShjaWl7OeTjRCKVesALbYNEkrKQGmTIl8DkVJMZJb4T/7LBuFb9rE8MCqKibQ9OvH6I8TTrC/t0eP2tIBDdm61TpeG6BS++23yOQuKqLiLSjgQ+r009lsJFzS0wM/3AAq/dmzG58XqX1IWuEri7x5M2UvKuK55csZry/CzlDBOmqDYedO/+OFkmmsKAqAZFb4lZWM7GiYier1Mvtz8mTrfqsAIzYee8x+7COOsH8Y5OVxZRsuJSXAUUex/+rmzawBP3cucNJJwAcfhD5eYSEzRO2StBpi9b4XXvBf9hiwrntTXc3Y9nBj2EVYsmHkSD6c77qLD2qApQ8qKqzvS0/XrFhFCYdg7D6xOkKy4S9axA5Tdjbr/ffn+/73P5G+fUUyM3kcckjgTlSVlexg1dDha4xI+/aMiAmXxx+3d67us0/oESh3323vb7Cy4593Xv37vV6R/Pzg7rfzaUyZEvr34PXSYerx1Hbkys6mL+J//+N7Lr3U+rtyuxlxoyiKiARvw0/exCtfn1V/1wGuphctYuao10sTSiAyM7nyHDKEjTnS0njvvvvSRh1scw8rpk2zNxeVlADffMPiaMHy4ov2zUoa4vEA991X/1xxMU1Y4ZKZCbRtG/p98+dT9rqF1yoqeJx7Ls10zzzD1zNn0mQkwvlef50muViwYwfw7ruUc+BA7v4UJUmJSOEbY0YAuAdAdwD9RGRxnWu3AfgTgBoA14jIR5HM1YhDD6W5xqp4l8sFXHxx/XN2Kf52dOkCrFzJh8XatXzdv3/kKf/+lHNaWvDK20dDk1bD8XzyHncc8MQTjR2qLhcfjuHU2geohM88M/T7nnrKvspmcTHwxReUefp0NjRfvJhJWgMHBucrEOEYb71F/8bZZ7POv11HLSuefBKYMKH2+0lPB448krX/c3ODH0dREoVgtgF2B6joDwLwGYC+dc73ALAcQDaALgDWAkgPNF7IYZnvvMPtfd3tfmamyH77iezaFfK2KCbccAOTouzCDe1CRe0YPtw+16BFC5Hdu0XKyvyPMXKkfUimlVkI4Gdwu0XmzAnvezj6aPs5mjUT+ec/wxtXhCa5IUPqm4tyc0WOOir47/fTTxv/bfnMThdcEL5siuIAiEVYpoisEpE1FpeGAZghIhUi8jOAnwD0i2QuS845h47Oo4/mVr9ZM9ZvX7yYNd0TkRtuoNO44U7B7WbFxlBrvt99N8driNsN3HMPdzYul/8xHn+c5iqrcerSqhXwxz8Cxx/P+jzff88Io3A45hh701hlJXdw4fLww7VNX3xmv+JimssmTAhujAcesN49VVQw0mnnzvDlU5R4EcxTIdCBxiv8KQAurvN6GoDzAo0T9UzbRGXlSq42s7O5Cm3ZUmTyZPsSCIH4/HORrl25Is3L4wr50UdDG6+wUOSJJ7grsHPO3nlnePJZsX49P7vVCvrkkyMbe5997HcPbrdIVVXgMfw5sps3Z9CAoiQIiJbT1hgzD0B7i0sTRcQiqDs0jDHjAIwDgE6dOkU6XHLQvTv7zG7bxpj2/faLLIb9uONYKO2HH+gQ7tEjdMdyXh4bio8YwTyG33+vXeF6PJT51lvDl7Eh++3H3dnw4bUllqurufJ/+236TWbNok/jpJO4iwvWf7J9u/216mqu9gPtAAsK7NszVlQwWU1RkoyAWkZEBocx7mYAHeu8Lth7zmr8qQCmAkDfvn39hN00Qdq14xENjGHLvkjp0IHO6pdfpuLNymKbxREjIotOsuL445nENn8+o2EOP5yfYcIEZtJ6vVTQDz7Ia3Pn8uETiIIC+9LNOTn2ORZ1ufFGVvFs6FjOyKDzPj8/8BiKkmgEsw0IdKCxSecQ1HfaroMTTlul6TFjhrWz1OViXH4wPP209Rhud/BmKa9X5LLLeI/PUZ2bK9KpkxZbUxIOBGnSMSLhL6qNMecAeBJAWwC7ASwTkVP3XpsIYAyAagDXiciHgcbr27evLF68ONDblKZM794s2WCFy0VzTaCQSBHgqquYQQxwp5CezvDR6dNDM58tXsz+unv2AKecwrr52dnB368oMcAYs0REApYAiEjhRxtV+ApatKBytSI3l4XaDjgguLHWrWPj8poaRhN17x49ORUlgQhW4SdvLZ1EYf58rvw6dmQFyZkz/WcAK/7xl7VbXe2/IXtDunZlh68bblBlryjQnraR8fjjwMSJtdEsmzbRHDF/vjPle/fsoUli5Upg//2BUaOANm2iP0+8+O9/WcrCiqws5l00axZbmRSlCaEmnXDZto2hhValENxuKq/evaM333//S7OE18vIEV/y1ltvsRdsMlFWRgVet3+tCL9PO4W/3358mIZaIkNRUgA16TjNrFn2dVkqKujoixZlZVTqRUW1YYJlZdxZnH++/7jzRGLGDJpZ8vL4wLrwwtpyyEuXspGLHfn5quwVJUJU4YdLcbF9wbGaGnvHYzi88459gxOvF3j11ejN5RRPPw386U/Azz/zs1RV0d/Rpw8V/e7d9Vf8Ddm9O3ayKkoTRRV+uPir2piby05Q0eKXX+yrYpaXB25eEiuqq7nzGTuWztKFC2mqKS+3blZTXU1l/8wzNH/ZNTzJzAROPNF5+RWliaMKP1yOPJLZnw1jsjMygNatWdM9WnTrZp9h6nazr2wwlJUBn3/OssG+cgbR4vffgV69gEsuYc3/KVNYEmHUKCp+u7II5eU09bRuzcJsVsXjsrMZaaMoSkSowg8XY4APP2Qyj8vF6BGXCxgwgHVyolmGYNgw+/HS0hrX/rfi8cdZxuGss+gPaNeO5ROixdixjHv39SfwerminzWL35M/fA+DJ5/kAyI7m7skt5uN0z/6iLZ/RVEiQsMyIyEvj/Vmtm+nbTo/35kaK9nZwLx5wKBBXJmXlHDFbwwTiwJ18Xr55frhoz7Gj2fc+5Ahkcm3ezc7gVntGkpL2THKjpyc2gdWRgYbzWRl0cbv9TIRK9SS0YqiWKJhmclERQWV5w8/cOV73nmBi4kFCnfs3Zt14iPhhx/ofLXqPgZQaT/6KGvo133oZGayWNu33zIC5/77revQ5+XxPZ07RyanojRRgg3L1BV+MpGdDVxwQWj3FBayIqUd335b+3tJCfDmm8Dq1UzsuvDC4EIh8/Pto4gAjjVmDO30EycycczlAkaOZCXM5s2p5B980No5XV4OPPKIM8lsipJCqMJv6lh116qLb4ewcCFw6qmMnCkpoRnlxhvZ3WnQIP9zeDw0y7z6auNENLcbuOMO/j5sGA+vt3EOw3ff2YdlVlXRZKQoSkSo07ap88039rV9srKoqMvKgNNOY+6AL7GrtJS/DxvmPyHKxxNPsBGL281xc3K4I5kwgQ3E62KVsJaTwweBHWrHV5SI0RV+U6a4mKv2qirr6/vuS7v5O+/YJ5GJcOV+zTX+58rJYTTNsmXAp5/SZHP22ZwjGHr1oq3fyg+Qk8MoIEVRIkIVflPmrbfsFXlWFpOhWrViOGXDzk4+SkuBNVZ96m3o3Tu8GkLGAC+9xLDR8vLa1b7LRR/AuHGhj6koSj3UpNOUWbXKXpFXVta2Adx/f/+JXQcf7Ix8DRk0CPjyS+4M2rZlVM4ddwBffaUmHUWJAhEpfGPMCGPMCmOM1xjTt875k40xS4wx3+39eVLkoioh07WrvaL0eGrDHM89175MhDHBJXZFi8MOY42dbduY2zBxYnB9bBVFCUikK/zvAZwL4D8Nzu8AcKaI9ALwRwBJUN2rCTJypH2EjjG1IZ4uFxuEN29e2z7Q4+Hx3nuBE7sURUkKIrLhi8gqADANlIqI1M3kWQEgxxiTLSI21bEUR2jRgqUNzjmHztfSUq74jWECV91mIv37A5s3M3N4zRruDi64QBuOKEoTIhZO2+EAlqqyjxMnn8xqm6++yozYgw+miaZVq8bv9XiA0aOdkWPzZhZJ27GDD5czzgitmbiiKBET8H+cMWYegPYWlyaKyOwA9x4C4G8ATvHznnEAxgFAp06dAomjhEPr1sB118Vv/qeeAm66ibuMigqWSmjdGliwACgoiJ9cipJiBFT4IjI4nIGNMQUAZgG4RETW+hl/KoCpAGvphDOXksAsWsTkq7oZuEVFNC+dfTagtZMUJWY4EpZpjGkB4AMAt4rIf52YQ0kSJk+27vtbU8Ow0ZUrQxuvooLZw2vW2GcQK4piSaRhmecYYzYBOBrAB8aYj/ZeugpANwB3GWOW7T3aRSirkoysXm1fMiEzk0lfwSACPPQQ4/NPOAE44gg2hlmwIGqiKkpTJyKFLyKzRKRARLJFZB8ROXXv+ftFxCMivesc26IjspJUdO9u3+y9qir4xiaTJwP33ktzUGEhTULr1gGnnw6sWBE9eRWlCaOZtoqzXH894/wbkp7OiKEePQKPUVkJ3HefdenksjJg0qTI5VSUFEAVvuIsRx5JU4zLVdv/Ny8P6NiRpZeD4ccf7c1CXi/w2WdREVVRmjoaCK04z5VXMiIn3Dh8j8e+CBygdXYUJUhU4SuxIT+fDVXCoXNn9rq1iuhxuYA//Sm08UQYLrplC01KBx4YnlyKkmSowleSg5dfBk48kTZ7XzvFnBw+DEJJKluxgiWYt22jM7mqCujXjwXbWrd2RHRFSRTUhq8kB337Mv7+j39kdu6BBwJ33w18/XVtwbdA7NnDrlzr1rHRSmEhHyBffsmOXxrXrzRxdIWvJA/dugHTpoV//yuvWCeBVVUxCezrr+lfUJQmiq7wldTh88+tQzsBRvssWRJbeRQlxqjCV1KHDh3sk8AyMtSGrzR5VOErqcPYsdZJYABX+GecEVt5FCXGqMJXUofDDmPmr68JDMCVvdsNTJ+urRSVJo86bZXU4v77gVNPBf7+d2D9ehZhu/762DVqV5Q4ogpfST0GDuShKCmGmnQURVFSBFX4iqIoKYIqfEVRlBRBFb6iKEqKoApfURQlRTCSQAWjjDHbAfwSbznCpA2AHfEWIkL0MyQGyf4Zkl1+IPk+w34i0jbQmxJK4SczxpjFItI33nJEgn6GxCDZP0Oyyw80jc9ghZp0FEVRUgRV+IqiKCmCKvzoMTXeAkQB/QyJQbJ/hmSXH2gan6ERasNXFEVJEXSFryiKkiKowncAY8yNxhgxxrSJtyyhYIx52Biz2hjzrTFmljGmRbxlChZjzGnGmDXGmJ+MMbfGW55QMcZ0NMbMN8asNMasMMZcG2+ZwsUYk26M+cYY8368ZQkHY0wLY8w/9/5fWGWMOTreMkULVfhRxhjTEcApADbEW5Yw+BhATxE5FMAPAG6LszxBYYxJB/AUgNMB9AAw0hjTI75ShUw1gBtFpAeAowBcmYSfwce1AFbFW4gIeALAXBE5GMBhSO7PUg9V+NHnMQATACSdc0RE/i0i1XtffgWgIJ7yhEA/AD+JyDoRqQQwA8CwOMsUEiLyq4gs3ft7Eahk8uMrVegYYwoADAXwfLxlCQdjTHMAxwGYBgAiUikiu+MrVfRQhR9FjDHDAGwWkeXxliUKjAHwYbyFCJJ8ABvrvN6EJFSWPowxnQEcDmBhfCUJi8fBBY833oKESRcA2wG8uNcs9bwxpsm0QtMGKCFijJkHoL3FpYkAbgfNOQmLP/lFZPbe90wETQzTYymbAhhjcgHMBHCdiBTGW55QMMacAWCbiCwxxpwQb3nCJAPAEQCuFpGFxpgnANwK4M74ihUdVOGHiIgMtjpvjOkFrg6WG/ZLLQCw1BjTT0R+i6GIfrGT34cxZjSAMwAMkuSJ2d0MoGOd1wV7zyUVxphMUNlPF5F34i1PGBwL4CxjzBAALgDNjDGvicjFcZYrFDYB2CQivt3VP0GF3yTQOHyHMMasB9BXRJKmAJMx5jQAkwEcLyLb4y1PsBhjMkAn8yBQ0S8C8AcRWRFXwULAcJXwMoDfReS6eMsTKXtX+DeJyBnxliVUjDELAIwVkTXGmHsAeETk5jiLFRV0ha/UZQqAbAAf792lfCUil8dXpMCISLUx5ioAHwFIB/BCMin7vRwLYBSA74wxy/aeu11E5sRRplTlagDTjTFZANYBuDTO8kQNXeEriqKkCBqloyiKkiKowlcURUkRVOEriqKkCKrwFUVRUgRV+IqiKCmCKnxFUZQUQRW+oihKiqAKX1EUJUX4Py1x31zeBvHLAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f575d6587b8>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.datasets import make_blobs\n",
    "X,y = make_blobs(100, 2, centers=2, random_state=2, cluster_std=1.5)\n",
    "colors=np.array([\"red\", \"blue\"])\n",
    "plt.scatter(X[:, 0], X[:, 1], c=colors[y], s=50)\n",
    "for label, c in enumerate(colors):\n",
    "    plt.scatter([], [], c=c, label=str(label))\n",
    "plt.legend();\n",
    "#plt.colorbar();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.728361Z",
     "start_time": "2020-06-24T19:29:18.723459Z"
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.naive_bayes import GaussianNB\n",
    "from sklearn.naive_bayes import MultinomialNB\n",
    "model = GaussianNB()\n",
    "#model = MultinomialNB()\n",
    "model.fit(X, y);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Naive Bayes algorithm fitted two 2-dimensional Gaussian distribution to the data. The means and the variances define these distributions completely."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.733086Z",
     "start_time": "2020-06-24T19:29:18.729646Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Means: [[-1.64939095 -9.36891451]\n",
      " [ 1.29327924 -1.24101221]]\n",
      "Standard deviations: [[ 2.06097005  2.47716872]\n",
      " [ 3.33164807  2.22401384]]\n"
     ]
    }
   ],
   "source": [
    "print(\"Means:\", model.theta_)\n",
    "print(\"Standard deviations:\", model.sigma_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's plot these distributions. First we define a helper function to draw an ellipse that gives the standard deviation in each direction from the origo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.740310Z",
     "start_time": "2020-06-24T19:29:18.734352Z"
    }
   },
   "outputs": [],
   "source": [
    "def plot_ellipse(ax, mu, sigma, color=\"k\", label=None):\n",
    "    \"\"\"\n",
    "    Based on\n",
    "    http://stackoverflow.com/questions/17952171/not-sure-how-to-fit-data-with-a-gaussian-python.\n",
    "    \"\"\"\n",
    "    from matplotlib.patches import Ellipse\n",
    "    # Compute eigenvalues and associated eigenvectors\n",
    "    vals, vecs = np.linalg.eigh(sigma)\n",
    "\n",
    "    # Compute \"tilt\" of ellipse using first eigenvector\n",
    "    x, y = vecs[:, 0]\n",
    "    theta = np.degrees(np.arctan2(y, x))\n",
    "\n",
    "    # Eigenvalues give length of ellipse along each eigenvector\n",
    "    w, h = 2 * np.sqrt(vals)\n",
    "\n",
    "    ax.tick_params(axis='both', which='major', labelsize=20)\n",
    "    ellipse = Ellipse(mu, w, h, theta, color=color, label=label)  # color=\"k\")\n",
    "    ellipse.set_clip_box(ax.bbox)\n",
    "    ellipse.set_alpha(0.2)\n",
    "    ax.add_artist(ellipse)\n",
    "    return ellipse"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we do the actual plotting:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.822876Z",
     "start_time": "2020-06-24T19:29:18.741616Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEJCAYAAAC+I6F6AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAGspJREFUeJzt3XuQnfV93/HP9+zu2ZukRbeVZBB3cYvDxWxjDAyOzYQ6dm906vFMa8dm6sie2AFTOm1ajy9px9PMhFKwk9bBOJBCOmnadHDT2gltAjYE00QEAYaYi0BcBOiGJCSt9nbOt398z8NZrXZXu/vbc57nnPN+zTzznNvufnW0+3zO7/Y85u4CAGCpSnkXAABobQQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhS6CAxs51m5nNsb+VdHwBA6s67gAU4JOn2WR4/0uxCAAAnsiKvbDeznZLk7mfmWwkAYC6F7toCABRfK3Rt9ZrZJyWdLumopKck/cjdK/mWBQCQWqNr64xZnnpZ0g3u/sPmVgQAmKnoQfI1SQ9LekbSYUlnS/qipK2SxiR9wN2fnONrt9Zep8HBwcsvuOCCptQMAO3g8ccf3+fu6xfy2kIHyVzM7FZJt0i6392vP9nrR0ZGfNu2bY0vDADahJk97u4jC3ltqw62f7u2vybXKgAALRske2v7wVyrAAC0bJBcUdu/lGsVAIDiBomZXWhmJ7Q4zOxMSb9Vu3tfM2sCAJyoyOtIPiHpFjP7kaRXFLO2zpH0MUl9kr4v6db8ygMASMUOkgclnS/pMklXKcZDDkp6RNK9ku71VpxyBgBtprBBUltsyIJDACi4wo6RAABaA0ECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkhQ8SMzvNzH7XzN4ws3Ez22lmt5vZ6rxrAwBI3XkXMB8zO0fSo5KGJX1P0k8l/ZykmyR9xMyucvf9OZYIAB2v0EEi6T8qQuRGd/9W9qCZ3SbpZknfkPT5nGoDWlK1KlUqsc12u1qV3OO12T67bRa3p+/NpFJJ6uqKbbbb2evRnsyn/6YUSK018qKknZLOcffqtOdWSnpTkkkadvej832vkZER37ZtWwOrBfJVqUgTE9LkZOzHx6XR0diOHYvHsy0LhKX+6WdfO/17zHU7u18qSd3dUrks9fZKAwOx9fdLPT3xeLkctwmdYjCzx919ZCGvLXKL5EO1/QPTQ0SS3P2wmf2FpOskXSHpz5pdHNBsk5MRCseOSYcOSUePxu3RUWlq6vjXmsWBO9t6euIAXirlc6CuVuvbsWPS4cP1UJsZOn19ETCDg9LKlbHv74/HS4Uf1e1MRQ6S82v75+d4/gVFkJwnggRtJDvYHjsmHTkiHTggHTwYrYxMT089IFatiu6jIiuV6iFQLs/9Ovd66+roUen11+tBUypJK1ZIq1dLQ0P1Fk1vb+Prx/yKHCRDtf2hOZ7PHj9ltifNbKukrZJ0+umnL29lwDIaH4/AePttac+euD394Jl1Ba1alW+dzTC9JdXff/xz7vFevfmm9Oqr9feou1tas0basCHeo8FBWi7NVuQgSeLud0q6U4oxkpzLAd5VqUjvvCPt2ye99VZ88naP1kV/v7R2LeMEs8m6vfr6jn+8Uomusr17430slaR166SNG6VTTokQRmMVOUiyFsfQHM9njx9sQi1AkqmpCI433oh9pRKfpAcHpfXr866utXV1RZfXihVxv1qtB0u1GuH8nvdEsHRCqy4PRQ6S52r78+Z4fkttP9cYCpC7w4elXbuiK6ZSiU/Hq1fT9dJI2VhKFiyTk9Irr0g7dkSQnH12tFh6evKts50UOUgerO2vM7PSLNN/r5I0KumxPIoD5jI1FZ+GX345Zlf19EQXS9EHxNtVT0+MoUjS2Ji0fXuEzebN0mmnxcwwuhLTFDZI3H2HmT2gmJn1BUnfmvb0r0salPQ7J1tDAjRLpRKzjJ57LrpUVqyQhofzrgrTZWMslUoM2r/ySrQQL7qIbq8UhQ2Sml9RnCLlm2Z2raS/kfR+xRqT5yV9OcfagHcdOCA980x0Za1ZE+MfKK6urmglSjFL7pFHosvrnHPo8lqKQvfUuvsOSSOS7lEEyC2SzpF0h6QrOM8W8lapSM8+Kz36aMwYGh4mRFrNihUx4eGVV6SHH44PBVicwv/Ku/trkm7Iuw5gpslJ6cknYzxkwwb62VtZNmV4bEx67DHpkktiphcWpvBBAhRRtSo99VQsImQcpH309UW31/btsRB07dq8K2oNhe7aAorqtddiFToHmvbT0xOnYHniieNPS4O5ESTAIlWr0osv1qeUov309cU07j178q6kNRAkwCIdORLjIwyqt7eVK2MxKU6OIAEWqVo9+WvQ+kqlmJWHkyNIgEXq7z/xOhpoP6OjsVgRJ0eQAIvU2xsnADw01wUO0PKyU9YzBXhhCBJgCS64ILq4mNXTnvbtk848s776HfMjSIAl6O+X3ve+uK7IsWN5V4Pl4h4LTNetk86b67zjOAFBAizR2rXSBz4QrZIDBxgzaXUTEzHd99RTpcsuY1beYhAkQIKhIemqq2J1+549tE5akbu0f39cqfKyy6T3vpdT/i8WmQsk6uuTLr44rm3x9NPS7t0RMDMvCYtiqVaja3JiQjrjDOncc6VyOe+qWhNBAiyTNWukq6+OIHnppWihlMuxsI1PuMUxNhan+zeL8N+8mWuRpCJIgGXU1RVTRjdtioPV66/Hebmq1QiU/v68K+xM01sfK1dG99XwMC2Q5UKQAA1gFp9yL7pI2rIlppNmrRSzCJSBAa7d3kgTEzHuMTkZAX/aaTGQvmoVp/xfbgQJ0GA9PdFC2bQpVksfOhSXed23Lz4pm0Wo9PdzgEsxNRXBka3tGRiIsY9166IVwiysxuGtBZpoYCC2TZviPE5HjkSw7N4dM4eyU6+UyzFYXy7TapnN1FSMdYyN1c991tcXFxhbvz6uejgwkG+NnYQgAXLS1RWzu4aGpNNPj2A5ejRaLQcPxnbgQP1AaRanZ+nt7Zy+/UoluqjGx6OLKlur09cXq87POCNCY3CQ8ac8ESRAQXR1Rf/9qlVxLi8pDpxjY7E+5ejResC8886JCyC7u6Mbrbu7vhWZe4TD1FTspwdFJgvb4eHY9/fH1tOTT82YXcF/1YDOlg3M9/fH9OLNm+PxajU+qU9MxAF4YiJaMqOj9VbNbOcBM4uDtVl0mU3furrqj2djNfON2WQHffeop1KJ25XK8ffn+7cNDNRbE4OD0dLo6YkWV7kcYci4UfERJEALKpXioDvfosdqtR4y0w/u2ZY9N71VMDERt7Ovn77PmB2/lUpxwM9aClmrqFyO21lITd+y59AeCBKgTZVK9TEVoJGYDwIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhSyCAxszPNzOfZ/iDvGgEAoeiX2n1S0v2zPP6TZhcCAJhd0YNku7t/Pe8iAABzK2TXFgCgdRS9RfIeM/ucpLWS9kv6sbs/lXNNAIBpih4kv1Db3mVmD0n6tLu/mktFAIDjFLVra1TSv5V0uaTVte2Dkh6U9POS/szMBuf7Bma21cy2mdm2vXv3NrhcAOhcDQsSM9t5kim8M7f7sq919z3u/lV3/2t3P1jbfiTpOkn/T9K5kj4738939zvdfcTdR9avX9+ofyYAdLxGdm3tkDS2iNe/cbIXuPuUmd0l6f2SrpF0xxJrAwAsk4YFibtf26BvnfVTzdu1BQBojqKOkcznitr+pVyrAABIKmiQmNn7zOyE2szsWkk31+7eN/N5AEDzFXX6722StpjZo5Jerz12saQP125/xd0fzaUyAMBxihok90q6XtLfkvSLknok7Zb0h5J+y90fzrE2AMA0hQwSd/+upO/mXQcA4OQKOUYCAGgdBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSdOddADAn99jM4n62B1AoBAkay12amJDGxqSpqdgmJ6Xx8Xh85n5qqv51ZrGfziy2nh6ptze2cvn4293dsWWPd/NrDjQSf2FYHlkYjI1Jx45Jhw5Jhw/HVq0e35owk7q6pFIp9tmBv7c3Hpuv5ZG1UqpVqVKRRkfjZ1QqsVWrJ35Nb6+0cqW0alXs+/piy34egCQECRbPPQ7gR45I+/dLu3dHiGRKpWgNlMvSmjXLe7DOWiSl0sJbGlNT9XCbmqq3ckol6ZRTpOFhaWhIWrEiWjoAFoUgwcm5R2gcPizt2SPt21fvgiqXpYGB+LRfVFmLZ2Dg+MfdI2Cee64eLitXRrCsWRO3y+Xm1wu0GIIEs3OX3nknguO116LFUSpFl9CqVdEl1erMIlymB8zYmPTKK9KOHXF/7VrptNMiWHp786kTKDiCBMcbHY2uqp0746Da3R3BMTSUd2XNkY2fSPUuvCefjPvDw9Lpp0urV7dHkALLhCBBOHhQevFFae/eaHkMDRW7u6oZzKTBwdjco2vvr/4qxlHOPlvavJkxFUAECQ4fll54QXrrrThgDg/nXVExmcVg/IoVMT70wgvSyy9L558vbdpECwUdjSDpVKOj0ksvSa++Gl05GzbkXVHr6O6W1q2L9TBPPx2hcuGFEcJMJ0YHIkg60a5d0k9+Ep+ih4dZMb5UPT3x/o2NSU88EQPzl1zCoDw6Dh+fOol7jINs3x7rJ1avJkSWQ19fBMrhw9Jjj8VUaaCDECSdolKJVsjzz8dBj9OGLL9TTon9j38sHTiQby1AExEkncA9WiG7dtGP32grVsS6lMceI0zQMTiidILdu2Nbv56urGbo64tAefrpaAkCbY4gaXeTk9Izz9S7XdAcAwMxM+611/KuBGg4gqTdvf12nJmXmUTNt3p1TLGeeSp8oM0QJO1u926pvz/vKjpTd3eco2x0NO9KgIYiSNrd2Bin8ciTWXQvAm2s4UFiZj1mdpOZ3W1m281swszczD67gK/9tJn9pZkdMbNDZvaQmf2dRtfcVnp66qd8R/O5c/oUtL1mtEgGJd0u6TOSNkp6ayFfZGa3SrpH0iZJ35F0n6SflfTHZvbFRhTaljZsiGtuoPkqlejeGhzMuxKgoZoRJKOSPirpPe6+UdLvnuwLzOxKSbdI2iHpYne/2d2/IOlySW9LutXMzmxYxe1kzZroXqFV0nwHDsRp51m3gzbX8N9wd59w9x+4+5uL+LLP1/bfcPd3V3W5+05Jvy2pV9INy1dlG+vrky64IGZvoXnGx+unmwfaXFE/Kn24tv+TWZ77wYzX4GQ2b45rixw8mHclnWFqKlojP/MzTHRARyhckJjZoKRTJR2ZoxXzQm1/XvOqanGlkjQyEqut9+3Lu5r2NjYm7d8vXXopp+ZHxyhckEjKrul6aI7ns8fnXaptZlvNbJuZbdu7d++yFdeyenulyy+XNm6MtSXVat4VtZ8jR6SjR6UrrpBOPTXvaoCmWVCQmNnO2pTdhW73Nbrwk3H3O919xN1H1q9fn3c5xdDdLV18sbRlS1xSl9OdL49KJVohknTllTHBAeggCz2X+A5JY4v4vm8soZZM1uIYmuP57HE6/JfCLIJk3Trppz+N1smqVax+X4pqNcZC3KVzz40ZWoyJoAMtKEjc/dpGFzLtZx01s12STjWzTbOMk2yp7Z9vVk1tafXq6ILZv1969llpzx5paIhzci2Eu3ToUJzD7KyzpDPPjNlxQIcq4hiJJP15bf+RWZ77xRmvwVKZRcvk6qtjcHh8PALlnXcYQ5nN+HhMVti7N07Jf801MbWaEEGHK+pl8r4t6VOSvmxm92drSWqLEL8gaVzS3blV125KJWnTpphldPCg9Prr0ptvRpgMDMTK7E69jsnkZARrpRLvxYUXxsXBCA/gXU0JEjP7NUkX1O5eWtvfYGZX124/4u53Za9390fN7DZJ/0zSU2b23yWVJX1C0hpJv1pbnIjlVCrFQPGaNfVFjLt2xSdwd6lcjoNpuZx3pY1TrcYpZY4di9t9fdF9NTwcY0mdGqjAPJrVIvmIpA/OeOzK2pa5a/qT7n6LmT2taIFslVSV9NeSftPd/1cDa4UUYbFxY2wTEzEmsH9/hMqh2nyIrq4Ilt7e1j3ATk1FaIzV5pKUStHdd9ZZMWa0cmXr/tuAJjHvgIvujIyM+LZt2/Iuo31MTMTU4QMHYkzl0LQlP6VShFC2FeUgPDUVdU9MHH9a93I5gmPdumhxDAxwbixAkpk97u4jC3ltUcdIUGTlcr0L7Jxz4iA9OhqD0UePSocPx7jC/v3HXx3QLFoxs21LOXi7x9hFpRLdUJVK1FKtxjb9Z5fLERRZF1VfX30rStgBLYogQbru7jg4SzGbKeMe4TI+Hl1Hk5Oxn5ioP561EmY7O7HZiUE0837W8unrO/52b+/xt7v5VQcahb8uNI5Z/VP/0FzrS2uyFsX0loR7bFmLwax+u7ubC0YBBUGQoBhKJcYmgBbFXy4AIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkBAkAIAlBAgBIQpAAAJIQJACAJAQJACAJQQIASEKQAACSECQAgCQECQAgCUECAEhCkAAAkhAkAIAkDQ8SM+sxs5vM7G4z225mE2bmZvbZeb7mM7XXzLV9vtF1AwAWprsJP2NQ0u2127slvSVp8wK/9nuSts/y+LZlqAsAsAyaESSjkj4qabu7v2lmX5f0tQV+7f3ufk+jCgMApGt4kLj7hKQfNPrnAADy0YwWSYpLzexLkvok7ZL0oLu/nnNNAIBpih4kN824XzGzuyR9yd3H8igIAHC8ogbJy5J+VdIDkl6XNCTpakn/TtLnJK2S9I/n+wZmtlXS1trdI2b2XMOqXZh1kvblXENR8F7U8V7U8V7UFeG9OGOhLzR3P/mLzHYu5ptK+n13/+Qc3+vrisH2X3b3uxbxPWVmmyU9KWm1pEvd/cnFfH2ezGybu4/kXUcR8F7U8V7U8V7Utdp7sdAWyQ5Ji+lKemMJtZyUu79mZt+X9E8kXaMIFQBAjhYUJO5+baMLWYS9tf1grlUAACS15ilS3l/bv5RrFYt3Z94FFAjvRR3vRR3vRV1LvReFDBIzO6Fv0MxKZvavJH1AMQj1J00vLIG7t9QvRiPxXtTxXtTxXtS12nuxoMH25B9i9muSLqjdvVTSJZIelfRC7bFHpg+8m5lL+oliDGSXYtbWVZLeq1gpf727P9DwwgEAJ9WsIHlI0gfnecnvuftnpr3+NyX9nKQtktZIqkp6VdL/lXSbu7datxYAtK2mBAnmV1tk+U9rd7e4+4t51tMsZrZF0j+U9LcVHxo2SDog6TFJt7v7gzmW1zBmdpqkfyPpI5LWSnpT0v2Sft3dD+RZW7OY2VpJ10v6mKSflXSqpAlJT0u6W9Ld7l7Nr8L8mdknJd1bu7vo5RLNRJDkzMz+rqT/KemIpBXqrCD5A0mfkPSspEckvS3pfEl/T1KXpJvc/Zv5Vbj8zOwcRbfusOLs1j9VtL4/JOk5SVe5+/78KmyO2qUg/pMiRB9U9DhsUHywGJL0R5I+7h16gKqtmXta8XewQgQJ5mJm6xW/LA9J2qjo/uukIPmMpCfd/YkZj39Q0v+R5JLOdPc3cyivIczsTyVdJ+lGd//WtMdvk3SzpN9x97a/3o6ZfVgxhf9/T295mNlGSX+puNTEP3L3P8qpxNyYmSl+/8+S9D8k/XMVPEgKOWurg2QzM76QaxU5cfd7ZoZI7fEfKsK1LOnKZtfVKLXWyHWSdkr67RlPf03SUUmfMrO2XyPl7n/u7n88s/vK3d+S9O3a3Z9vemHFcKOkD0u6QfE7UXgESU5qn8b/gaTPdUJXxhJM1vZTuVaxvD5U2z8wywH0sKS/kDQg6YpmF1Yw7fh/vyBmdqGk35B0h7v/KO96FoogyYGZnSHpDkn3ufv38q6naGrvz7WKqd4t88e0AOfX9s/P8Xw2Hf68JtRSSGbWLemXandbaq1Yqtq//V7FeNG/zrmcRSnq2X/blpmVJP2eYnD9xpzLKRwz65X0+5J6Jf2LNpvFNFTbH5rj+ezxU5pQS1H9hmK92Pfd/U/zLqbJvirpMklXu/uxvItZDFokS2BmO83MF7HdN+3Lb1YMqv9yOxwkE9+Lmd+rS/GJ7CpJ/1XSrc36dyB/ZnajpFsUM9k+lXM5TWVm71e0Qv69u/8473oWixbJ0izpbMhmdp6kbyjmyH+/EYXlYFnODF0LkfskfVzSH0r6ZBtO/cxaHENzPJ89frAJtRSKmX1R0d37rKRr3f3tnEtqmlqX1n9WdHl+JedyloQgWYKEsyFfpOiyucHMbpjjNS/E7D9d7+73L/HnNM1ynBnazHoU3Vkfl/RfJP2Su1dSv28BZRdXm2sMZEttP9cYSluqXU77PyhOi3Stu+/JuaRmW6H678RY7e9/pu+Y2XcUg/BfalplC0SQNNdOSd+d47mPKdaS/DdJ79Re2/bMrKxogfx9xaeyG9p4RXO2Uv86MyvNWD+xUtGlN6pY2d8RzOxfKsZFtkv6BXfP+6qAeRjX3MeF9ynGTR5RfBApZLcXCxILYtr5yDppQWKvYsHVRxV/SFvbOEQksSBxOjP7iuJUMY9Luq6TurMWKuWKss1EiwR5+rYiRPYpzvL81Vma9Q+5+0NNrquRfkVxipRvmtm1kv5GcY2dDym6tL6cY21NY2afVoRIRdLDkm6c5f9+p7vf0+TSsAQECfJ0Vm2/TjH1cS4PNb6U5nD3HbXr7WQnbfyo4nxTd6iDTtqo+v99l6S5+vx/KOmeplSDJHRtAQCSsI4EAJCEIAEAJCFIAABJCBIAQBKCBACQhCABACQhSAAASQgSAEASggQAkIQgAQAk+f8qAxX9NbEwnQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f575aaa1160>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.xlim(-5, 5)\n",
    "plt.ylim(-15, 5)\n",
    "plot_ellipse(plt.gca(), model.theta_[0], np.identity(2)*model.sigma_[0], color=\"red\")\n",
    "plot_ellipse(plt.gca(), model.theta_[1], np.identity(2)*model.sigma_[1], color=\"blue\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Accuracy score* gives a measure about how well we managed to predict the labels. The maximum value is 1.0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.837418Z",
     "start_time": "2020-06-24T19:29:18.823858Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy score is 1.0\n"
     ]
    }
   ],
   "source": [
    "from sklearn.metrics import accuracy_score\n",
    "y_fitted = model.predict(X)\n",
    "acc=accuracy_score(y,y_fitted)\n",
    "print(\"Accuracy score is\", acc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The score was the best possible, which is not a surprise, since we tried to predict the data we had already seen! Later we will split our data into two parts: one for learning the model and the other for testing its predictive skills."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Another example"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's generate some more data using multivariate normal distributions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.944447Z",
     "start_time": "2020-06-24T19:29:18.839523Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJztnX+wXGd537/Pru61sr4moJUNxLL22jWh2DQQEBQXypgogKNk7LaTYaxee+QfU+Fr8KiTdDKApgmTGc3QhDaoY4yqYhHhvWOSplAzVLYxNNPQZLC5pjYYgbFBki3VqaWrECwJVbL09I+zh3vu2fc9533P7939fmbeuXd33z3nPWd3n+d9n+d5n0dUFYQQQiaPVt0DIIQQUg9UAIQQMqFQARBCyIRCBUAIIRMKFQAhhEwoVACEEDKhUAEQQsiEQgVACCETChUAIYRMKKvqHkASa9eu1dnZ2bqHQQghI8Pjjz9+TFUvdunbaAUwOzuLxcXFuodBCCEjg4gccu1LExAhhEwoVACEEDKhUAEQQsiE0mgfACGE1MHZs2dx+PBhnD59uu6hWFm9ejXWrVuHqampzMegAiCEkBiHDx/GRRddhNnZWYhI3cMZQlWxtLSEw4cP4/LLL898HJqACCEkxunTp9Htdhsp/AFARNDtdnOvUKgACCHEQFOFf0gR46MCcGBhAZidBVqt4O/CQt0jIoSQ/EyEAsgjwBcWgK1bgUOHANXg79atVAKEkPJ56KGH8PrXvx5XXnklPvGJTxR+fC8FICJ7RORFEXkq8twaEXlERJ4Z/H2V5b1bBn2eEZEteQfuSl4Bvn07cOrUyudOnQqeJ4SQsjh37hw+9KEP4cEHH8T+/ftx//33Y//+/YWew3cF8KcAros99xEAX1fV1wH4+uDxCkRkDYA/APCPAbwdwB/YFEXR5BXgzz3n9zxAkxEhE0cJP/rHHnsMV155Ja644gpMT0/jxhtvxAMPPJD7uFG8FICq/hWA47GnbwCwd/D/XgD/zPDW9wN4RFWPq+rfAXgEw4qkFLII8JCFheDzNLF+vf09ZZiMqFQIaSgl/eiPHDmCyy677OeP161bhyNHjuQd7QqK8AG8WlVfGPz/twBebehzKYDnI48PD54rHZugtj0fEn6m584Nv9bpADt2mN9XhsmIfghCGswI24kLdQKrqgLQPMcQka0isigii0ePHs09ph07AoEdJUmAh5g+UwBot4Hdu4G5OfP78qw4fMYyIt8vQsafMn70AC699FI8//zyvPnw4cO49NJi581FKID/KyKvBYDB3xcNfY4AuCzyeN3guSFUdbeqblDVDRdf7JTSOpG5uUBg93qASPA3SYCH2D678+eT35t1xZFlLDm/X4SQIijjRw/gbW97G5555hkcOHAAZ86cwRe+8AVcf/31uY4ZpwgF8GUAYVTPFgAmL8XDAN4nIq8aOH/fN3iuVEK7+c03B4/vuw84eDBd+APZP9OsK44yxkIIqYAyfvQAVq1ahbvvvhvvf//78YY3vAEf+MAHcPXVV+c65hCq6twA3A/gBQBnEdjxbwfQRRD98wyArwFYM+i7AcBnI++9DcCzg3ary/ne+ta3alb6fdVORzWwmget0wmeL/v9/b5qr6cqEvx1PWcZYyGE+LN//36/NxT9o3fENE4Ai+oq01071tHyKIBeb6XADFuv536Mmj7Txo+FkHHHWwHURF4FMLbZQIuwm8/NuZmLqqBJYyGEjAdjmwqCdnNCCElmbBVASX4ZQggZG8ZWAWQN/xxVuFOYEOLL2CoAIBD2Bw8Gsfuu4Z+jiGmn8M03B4qPyoAQYmOsFcCkYNoprIP92EwbQchoctttt+GSSy7BG9/4xtLOQQUwBqRFNjFtBCGjxy233IKHHnqo1HNQAYwBLpFNTBtBSHmU4YN797vfjTVr1uQ/UAJUAGOAKeIpTlPCX+msJuPGKGfrpQKoiDIFXzTiCQicv1GaEv46yj8UQmyMcrZeKoAKqELwhRFPqkHSuyaGv47yD4UQG6OcrZcKoAKqFnyu4a9Vm2NG+YdCiI1RzjpABVABTRR8dZhjRvmHQoiNsrIObN68Gddccw2efvpprFu3Dvfee2++AxqgAqiApgi+6Ix/y5bqzTF5fih0HpOmUlbWgfvvvx8vvPACzp49i8OHD+P2228vZsARqAAqoAl5ie68E7jppuUZv6nWMRCsSnyErU/frD8UOo9J0xnZrAOueaPraHnqATSNOvP59/vm2gim1u26F5+pqlBNEbUdCPFhUuoBcAVQEXXOEFzNOuEqxdU0tG1bNjOSrzmniT4UMv5omE+loRQxPiqAEcZVkKYJyqg55vhxt2MsLABLS/7nW1gAbr11pTnn1luTlUBTfChkcli9ejWWlpYaqwRUFUtLS1i9enWu40hTLxAANmzYoIuLi3UPo5GEdvHoDLzTMdvUZ2cDQWui2wWOHUvv2+sFKxeXY8b7Rlm71qw44uOI4nOthBRB6Hg9ffp03UOxsnr1aqxbtw5TU1MrnheRx1V1g9NBXG1FdbRx8gEUjY9dvN9XnZ4e7ttuB69F/RPdrurUVLpdX8TuR0jyAST5H5JgTWRC3AB9AOOPj118bg7YsyeYZYd0u8DevcH/0QibpaXAJNTtJkfq2Mwv3e5y3yJDN0c2yoKQBpO7KLyIvB7An0WeugLA76vqpyJ9rgXwAIADg6e+qKp/mPfck8z69WYTjE0whwJz27ZAyC8tBf8Dw47cM2eAmRm7SQYIQlhNZpmdO4P/42abMHRzZgY4cWL4eFHlRAiphtwrAFV9WlXfrKpvBvBWAKcAfMnQ9RthPwr//LjuLQhn4SLBPoCo/T1UBCbSHMdJMf0LC/aNZhdcAExPr3x+enpZcRBCqqNoE9BGAD9SVYt7kBSFy6aq6AYqX9avTzfhmMwy4TltG82OHw/MUdFx79lDkw4hdVC0ArgRwP2W164RkSdF5EERubrg8xIDpiR0LnQ6wKZNybtvbcoh7Zxr1gTCfseOQMk891zwHu7qJaQGXL3FaQ3ANIBjAF5teO0VAGYG/28C8EzCcbYCWASwuH79+uJd5GOCyy7cpEgdW2u1gkgg2+thBI7t3C47jV13EJcd+cPIIjKOwCMKqEgFcAOArzr2PQhgbVo/hoHacQkDtfXJ08JQUZtwT1M6obBNG3vZaSaqSmNBSNX4KIAiTUCbYTH/iMhrRII6VSLydgSmJ4v7kbjgEgbqUirSlzVr7I7jpaVAlKa9P23sSU7korKVsjgNIQWEgQKAiFwI4L0APhh57g4AUNVdAH4bwLyIvAzgZwBuHGgqkhGXMNDQsbpli90p60OnA+TdGPmTn9iVSKsVOIZF7IqkqPw/zC9ESEFOYFU9qapdVf37yHO7BsIfqnq3ql6tqm9S1Xeo6t8Ucd5JxjUMdG4uiNJJo9dLjsUPo4xOnvQfa5Rz5+z5hkIllTQ1KCr/D/MLEcJkcCOLT279NWuSjxXm7tm506xU+v1id9/GBXzL8VtYZA2FJtRoIKR2XJ0FdTQ6gd2xRbTY8gDZHJ9pkTFJEUJltjBvket157lnhIwyqCMKqIxGBeBGUkRLUiRQXFG4CMM0hVJW8ylIMz8/fC0U9mRSoAKYMJLCKpPCMUN8QyKjwrTdrk4JhHsI0hRbPBR1ampYaYXXR8VAxg0fBcB6AGNAqxWINRO2iJpozn7XGgC+5y6LVauAl1/Of5xuF/jZz1hngIwXPvUA6AQeA2yRK+22WTiLrHR25gmJTHMwl0ERwh8IQlHz7gUoMuU1IVVDBTAG2CJabLH/qitnuHWFRHa79qijYNugP1nfF8V1L0A02Z7qcL4kQpoOFcAYYAsJ7fXM/ePP5wmJtMX0pyEShJ3aQlmzKp9f+7VhJTA1NZyCutOx73twPTd3E5ORx9VZUEejEzgfPknXouGdobPVhTz5htLG7htt1O0OX69IEBVkcvaa7k/oQHZxCNvyHkUd7IRUDRgFRELSolzyJkUzvd+lhYnfbILZV7F0OvY9CqY6yfH7Y4oeSrsPPnWZCakKKoCqGeFYQpsQ63bdL6nf9wsHjYZgxpWHKWTT1trt5eykaRvU0q4lizBnRlHSRKgAqmTEpYBrzYC0S3I9TtS8lDddtUi2FYjpWrKac0ZY95MxxUcBcB9AXvIE0TcA2/BNJF2S63E6nSA76b592UpVRglzCLkku4sTv5ak8fd6gUOcewPIKMB9AFUy4nmFfWoGJF2S63FOnQJ27cov/IFA8GcR/sDwtSSNv4zwTu4fIE2ACiAvDcgrnEeYmEJIs4RHhsdJSikdkrTonJpKf38RtFor79PcHHDNNfb+RYZ3cv8AaQyutqI6Gn0A9Zw+zzHz2vVbrXzv9/UhzM8vX7PLe9rt5fdkhdFDpExAJ3DFzM8vh8HYJERJ3sKyhEnW4WYpRN+E5qt40pRA0v3j/gFSJlQAZRP9dXe79lST0f4lrRKaJkzKKETvM6Ov6lzt9srr9vlKcAVAysRHAdAH4MuddwI337xswF1aAs6cWdknbjAuMWfA+jUnzM/XVNpwx45i8vH4MjMTfBxVce7csr8lbtNP+0ps2mQ+pu15QsqCYaA+LCwEwt/lnoksh6jYciZH+2Qcz8KtX8PWs3fjFC78+dOd6Zexe8+q2sIW61AAddHpBNfrUis5/LhHPHKYNByGgZbF9u3u08zoFLyoSKF4uM+2bZg7+6fYjX+FHg5CcB49HMTui36n1ph1WxK6Cy8cDrVctar88ZTJqVNuwh9Y/rjLiBxmWCnJhKutKK0BOAjguwCegMEGBUAA/EcAzwL4DoC3pB2zcT6ArNtmi/AB+Gx5rdmb2O8HKR3iw5qeNpdrjPrQx7VNTy9/3La0Fd1u9vtt+mr4JPUj44NJ/tpa0QpgbcLrmwA8OFAE7wDwaNoxG6cAkjyc3a65InvUM2jqE+/nm6ymod7EJCFnu9Q6Hchlt6mpZeWX9BXKQtIxRygrCSmIpiqA/wRgc+Tx0wBem3TMxikAW/5gW9inay5ml35FJe2piCzD7fdHN4zUpaVdW9aFW9pxGzAfIBVSlwI4AODbAB4HsNXw+lcAvCvy+OsANhj6bQWwCGBx/fr1Zd6nbLgGyLvG+uXt55O2s0J8Z/NRc1DdgrpIoe57D8q413HFwgR2401dCuDSwd9LADwJ4N2x150UQLQ1bgXgg2uAvmu/Ecs6midL5/z8eK8Ekq69jHsdVSwj9jUiGahFAaw4KPBxAP8m9tzom4B8KHoFoFre1C3puDnOGa805jML9q0x0PSWpNCK+Cht9zou3LkJbfypXAEAuBDARZH//wbAdbE+vxlzAj+WdtyRVgBF+wBs58irEJLOX8B0MYtjNzz8uKwCwuwgSbeyKN2edpym7RwnxVOHArhiYPZ5EsD3AGwfPH8HgDsG/wuATwP40SBcNNH8o2UrgCoMoa7nyDKWotbySVPCAqaLWYR4eBk+ysOlKlidLVrqMl7+0mXmbvsK+H5tuAIYf2o3ARXVSlMAtlqEtjDNJpKlAK6JpClhAdPFJCHuYhZxKRk5Pb380Zny8DRtJREmnjON1fWj9AlIC/tnrX1MRgsqgDRcppZN/lUk5S72Xcsn3QubEd5DybimWbadJjqMeNH4UODHN53FdfmoRRa5fJS2jy1UAtGVgcn8FCqBUZjrED+oANJwnRI2dV2cpsB8ftW+4ToZFGPSYiVJkLmc1tWkMTNTjdAuutk+yqSv8Kh/vUk+qADScDUuN9Uz5vILjzpyTWv/aJ6ANIN7u53LNObrZ7ZdXnj66Ow26Ra45OAfhWZSfkXsmm7q15vkgwogDddZb1OnSK6//m7X/TqTtuEWICl8Ik3ThKGLsogLzlFPMxH/KuYxrTX9603yQQXggm9RlyaRZZfVCEkKm7C2uSSSlEC4eOl2ze9vt7PtV6hTEYQrp1Wrsh/HlK+Qu4PHAyqALPj8AprwayliCphFUlR0aSaTUdlC1ZTBtIktdHKbXrvgAvv7bJY87g4eL6gAyqQpv5a0zGmdjuqFF2aTMDlt/kVdXlzHZtmy4NpGZQWQtYmY9yIkpeJuwNeAZIAKoEyaspMmLa/w/HxykHmatAiJSt54PGackldGvs5kn8utW0A3vXFFMDpQAZRJlXvpkwRqktRSzTcljibksUnVCgvfu96S6GvjlEeoKY1O49HARwGwJKQvRZV3TCNeafzQoeBxWOvPdr6wHqOtvqBIcLz5efPrU1NBZXfAXMw+pMLC967MzQU1dc+fB/buHS4/2URGqX5ynpKVpKG4aoo6WiNXAFX5ANJMTWnjsL0/Wj+g2x3eIRXdH+BTwaSClZHvrY+uCLLMeMNjj7t/gCuA8QI0AZVMFVFALgI1yT5vCmtpt4efm54efs41G1tUIlTgG8lzCt+cRPF9cnUL36pap5OeuZQ0GyqAccBV2tmmxXmcwFFF4uMDMGVp83UWx/dnRJL6CM5bBXgaaUI8TZ/XLZjLaO226saN5mtvQqQzyQYVwDjgau/w3TXl2kKpaivPtXHj8HjDNJdha7WGx5sjlKcnh6zC24U8CVRHfSexrdlm9lkUAJVGM6ACGBdcflFlxTCGwjspW1t0PEmG8rAiimquYP4+NmsHJ5wEmO12utboid/2qjdfV9lcF5VpeyOzFLyh0igeKoBRx+dXUdYKIPwFu0oOl+MVUOy3j83ak0MqOO8tMOIO3aidP9rHViqiqNvatBY3oWXxtSS9J8lKSV9D8VABjDJZQl1M/TduLEY6JEm8qORwPVZRthRP53Jei1q0TU8PW7tGuXW7K+9BloCupPf4zlEYbZQPKoBRJsv0y7S/P2n2XpTZKC1Hg6kVZUvxDC91va2utyYaTTvq/oG4qybpemyrrqT76/t1Y5rqfFABjDJFxNMn/YLDtXfYJ+tU1raGT2qm1JtZlVF82pr1tuK8861L+jiaqgRMJTTTbmeajo5+hWzbScJ+SfMDrgDKwUcBcCdw0yhip3HSls3du4F77gl2+3Y6wbZZX3q94Dj79tl3Cpu49tpgN/PS0vJzqv7nB4DTp5NfX1gAZmeBVguYncX6NSeM3dbL88u7q7F8W9KIfxyu76uas2eBiy4Cut3kXcdLSz+/VQCALVuCxyZOnQJ27VrepL60BJyI3d5uN/iKzM2Z702nE3wVTM+HG9FJBbhqijraRK4Asu40dkmE47Jxq4jdv6Y2M1P8lloPv0h/6hbt4OTK24oT2sfmobSX8a0IrqUifK1hVbZw1u7Sd2oqX62B+FctyfnOKKDiQZUmIACXAfhLAPsBfA/ANkOfawH8PYAnBu33XY49kQpA1f9X4WJXj0utrKYXFyUS/srL9pTazECWcfWxWXs4oIJz2sOBQPin3afYxxHbm+ZlDwfqjSKq0nEdzhNsX01TBBYphqoVwGsBvGXw/0UAfgjgqlifawF8xffYE6sAfEkystokVZap6vT0svRLM/z6VGGP5yfymX6WsTfCYoR2XZwlnT7sW5cSqKqFtzDNHUUlUDyVKoChAwIPAHhv7DkqgDLJ4jj2jcZptdJLZkWnda7HjaeL8FVMJmGd1w5juW+ukUQ2S5fIsj4e9wRz4b6/NF1Mh2/x+CiAQp3AIjIL4FcBPGp4+RoReVJEHhSRqxOOsVVEFkVk8ejRo0UOb3yxOYhVA69exMn5c+bmAi9drxd4B9tt+/E7HeBVrwo8iknMzATH9eHMGeCmm4IxrF0beBZ9MDi8Fzb1MYsDaOEcZnEAC9jsd8z4/Rw4lJ87ZHaYu6ZJDsXeoUPASy8FmbfjNCU9tM0BDLiNcffu4BhJxwGYYrp2XDVFWgMwA+BxAP/C8NorAMwM/t8E4BmXY3IF4IhL3F7aWjvNbuFiVvHdGFZEiwXk9+e/MWymCZ29LseL5iUKp/yDa+/hQPostt9XwTnn4ZuKrDXBmZy02JuZyZdn0HrvSCGgahMQgCkADwP4Hcf+BwGsTetHBeBBmtRI+6Wl2TdcJFL0HFnrEdtaq2XONhqTVNaEcThgP3bcV2JRqKm5iAbvsykKW0vb6F3EtokiW5gGyseMZXJAZ/EBMGoonUoVAAAB8HkAn0ro8xoAMvj/7QCeCx8nNSqADGTdSOaSzct3d1BREif0LcR//YZz2GbfgnPJEi1KTNlFo4e6eFG7eDGIJOrFBNDgfSZFkdZcZsLz881IQeGzaSza0oS3y+vMHZRO1QrgXQAUwHewHOa5CcAdAO4Y9PkwghDRJwF8E8A/cTk2FUAGsqSSCPH5hcbjIU27gouIeQx3D5swKDurmSZpBZCQD8I468cJ7XfvShxPXGmkmUzS/PVFL6jytOnp4cwjaSkf2u30r17a/IM7h92o3ARUVqMCyEBd06QyjdameH+LRDAL7JPan7rFfvyZGWskklWhdF9yvwe9XmYLXdoMOxSsdUYVRb9eSf1s1+eSeyjJJ8HcQSuhAph06jCUlm2cjpKSVtppw1eSFItIHKtJySR0HJRvUrpp08floldVy731Li3NVRRa8XxyFob3OU25cQWwEioAUj1Jv/wisn+GuOYzyNJCZ3BE4nitAFSdKp/0u3dpr/uSU6qJNL0aWsia4Bzudu2VSNvtbP6LNAVIH8AwVACkepJmv0kGXFfJEp6jYknn5QPIcG/yBl9FE7uOWwu/Pkl9KPyHoQIgdso0DyUdO2sdgFZr2T5SUyIdq0nJ9d4lSPm0oC3bbWu1XHfbnq9dkGdp0a+PzQTkmRF8YqACIGbKcBC7KpSmhbIU0dIC+NNCY0RSVwBJmTRDEnPy4ahqv++VmqnO1m4P+0P6ffMWEM7+zVABEDN5QkSjGHbJrvhlmn7Bk1JRPbw/Ltc7iA5Kspy5FrG3nUJwLnGlEV5CU3MUxTdm+yxcJ3XTGBUAMVNEtTFfYd7pNE+qFN2yZFqNSHKboPLR17Zb3MMBp5VG2lekzpYlymeSN41RARAzRawAyvQ4ukqfOpPq28YdGuRV3VJgOkgiH33d76t2xFLwJmWlUcbHW6QiyRLnX9RidxShAiBmipgWFT1FDI8XCnWXimTz835KoKppbdrU3SSBEuwUvkKsP/8N7cmhlc5qh5VG2lfE9+Mseo7gcsvizxWpTEYNKgBiJ69htMhfdxg47iNxrrrKT6BHvaZXXVWsZDJJl34/ea9CVAKlKORM+roAw7cp44frrU77mmTZFhIurmwJ6EyF721fEa4AqABIHkxSKfy1dbvpRWOiTbX8IPaqHdEmaWSTQA5TfKs8zyHofd/qeuumplbuS4gL4VB5RRdwLpvDZmaybSKznX/coQIg5ZIW7+8SVhIKuSrMM+FegrLP49rCVYmrkd8lh4KjdMuyqvDR0fFLis4Nwo/AFDjmM2/waYwCogIgdZImccZ1G2taCxMA2aRW0v3LYd/I4hzNq6NdzD7dbvG+/Ukw95jwUQCFloQkZIh46cleL3gclo7csSMoOTlpnD0LnD49XF+x0wnuScj27cCpUyv7qJqP6VBf0dYl6a22iqOuLC0NX0Kc48eBvXuB6Wn/409NDb8vfhuJGSoAUj5zc8DBg8D588HfaN1gk4KYFE6eXCnMRYAtW1beH5+iuQ6S2tYl6a07dpRfqzg8v0232Wi1gM99Dtizxz7HIAm4LhXqaDQBTShNstfnbVnsGtFCwbZ7kdHDmWiRS/DtZL18l32A4fmzWgMnzcafBugDICONKfnLqLWpKXsiG59m8pBGy2/GFYajEhiS87ZQn4HDOqtwtkX6hvorOuQi4gFMuZImDSoAMvrYgr5HpUWzleUt5tvt2msM+Ib0OOadWJEBVQ7p/MbvG0/lkt8vXi7apqeKigeYlHBPG1QAZDxoUoTQ9LTqxo1+72m1isuA6iMxrVuFExRGWh1kOWkU4q462iUip8itGpMaAaSqXgpAgv7NZMOGDbq4uFj3MEhduHoeW63AwezLzAxw4oT/++qi2wV27lz2brZagbyLI2K+H7OzwKFDw8+HjvfBa7M4gEOYNXY7eHDlc7YhuA4pzsJCEPhkGmb8eEnndT3fOCIij6vqBpe+hUQBich1IvK0iDwrIh8xvH6BiPzZ4PVHRWS2iPOSMWVhIRBWaUxNAf2+mwQycfJktvfVxdISsHVrcH8Ae+jOmjXm55NiQCPhuM/BfFzT211DRF37hQFjqsFH2+2a+6kmzw9araDNzi7fLjJMbgUgIm0AnwbwGwCuArBZRK6KdbsdwN+p6pUA/gTAv8t7XjKmLCwEQi5tCtjtBvF/c3PZAtV7PbugrJKW50/w1KlgigwEQntqarjP0hKwdu2w5EuKAQ3DcbtdrIdZUZje7rKNwzUmP9T7oeAGgGPH7IJe1a4gzp0LXj90aKXOJDFcbUW2BuAaAA9HHn8UwEdjfR4GcM3g/1UAjgGB+Smp0Qcwgbja/dN2yya1MGlb3U7m0JCexdeRVi8xvM7oeWw+gJhxv7/x3uHU0o7Fz8JMFVmKt9jcEy5V08Lz2aJuJ8kngCqdwAB+G8BnI49vBnB3rM9TANZFHv8IwNq0Y1MBTCCusYCmfDkuAj1MJ51URquKFk2dmSX+MYwycn2vrbSWJa9Qf/4bQTec1177ee3jX5YabJ8k5H2CnYqoeTTqjLQCALAVwCKAxfXr15d5n0gTybICiBKdMobTwXic/Px8/aWvonH8WVuv53cM0z0rSvLmJE1wu2YwneRCMCFVKwCagEhxuJhzbIVxXSREUq7+qtqFFxYT7xiuZrKumlSTJW+F0jRND7malNKsXID3vrlc1FGXuGoFsArAjwFcDmAawJMAro71+RCAXYP/bwTw5y7HpgKYUHyNyq4zVR+TSZ2t1XJbIczM+F2PzwrAV5EU8JHbBHfevW5JNYeKXtDEi+mYNnGXrQQqVQDB+bAJwA8Hpp3tg+f+EMD1g/9XA/gvAJ4F8BiAK1yOSwVAnEgqP+VSm6CJrQgTkYvkyZKqoiR7imm2XMQiJO02FnU5rrEIZZujKlcAZTUqAOLEKMzqq27R1BNJCXL6fb9qLBXnWUj6aIs4BuC+oEkz57jq67Id0j4KgOmgyeiTN2F90bTbdY9g5TbYl14Ctm0z74zavj2oTZBGTXmWbR+loHZDAAAQQ0lEQVStiHtsf9p2D5evT3R7iqp5f4Fr5u5GfV1dNUUdjSsA4kSeJDKhI7Uo53A0CVzTTU55E/lU4OFMctu45hdKsnC5LmhcTFEuK4Cx9AGU1agAiDMuu4FCoWcTWr7J3kyt1UqOoBm1FpauLCITaUZsQ3MxpSR9DD46K82MZEt7PT1tvn1lQgVAJpukX2vaL9AUQhKdKUerm6e1UfVNRCu5x6fPrttzCyTPqYraGOaiz6ens+2CLhoqADLZJEUFFYFPOOmoKYHQJJZ0H0MJZzlGv3uX9maOBbUEcED73btyScI8i42i9FRTInxcoAIgk00V5gkfn0GWspB1tlCKZVBexloCOKH9qVtyK4Es+YaK/Cq4pG1qQsoJKgBCqtiC6WIKSpktN7Zl3DfRwwHz4XCgsOmxr1Av+quQ5leoGyoAQrLiIy1cAsyL3tBVZZuaUl21yus9gnPmW4FzhU2P6873Y9s6EQ0AqxMfBcB9AISEuAR7R0kL6FYF9u4FNm1KT5rfRM6eBV5+2est1loCeK6wAPikujZVMDcXlKKI1iLodoE9eyrdIlEIVACEhGzfHhRciRItwBJnx470spWnTgH79gUbqHq9oH8TNoqVxA58DB2srLTWwUnswHa3qjAOJNW1qYq5uaBYTTj/P3Zs9IQ/QAVAyDK+U8u5OeCOO9KVQFjd7ODBYIfu3r2juSJIo93GHO7HbmxFDwchOI8eDmI3tmJu/hcLk5CmKmSuVcdIDFdbUR2NPgBSKVmNyy7hIXEvZdZKYE1tUQN4RTuE646398U25qKvBXQCExLB9ReWN2YwLbGaSZE0PV2Eqa1ePbw5rkLhP4oUme46DSoAQkKqjhlMq9EbZxRDRKem3FNL20JjJkxR2BZ7ZdQw9lEAEvRvJhs2bNDFxcW6h0FGmdnZZRt8lF4vsMkXTasV/IZNmM5pG1+3Cxw/bj9W3ZiuZe1aYGlpuG+3G3hJQ8Joq6jDvdOpPNNolSR9LUyIrEzo6oOIPK6qG5zGle0UhIwIVccMJuUvNnkpbR7NnTv9hf/GjYFgrgKT0jIJf9PzvtFWY4Dta2ELCKsqookKgIw3VccMmgS6SBAtZJrdzs2tDBGN5tz3FebPPhvMyvv94TFMT/sdy4VofYG05PwuifOrCuSvAZue37q15ogmV1tRHY0+AJKbCtMWrzhnEfZt3zoHUR9DvDhtmc5mFz9G9J7XvZXXgumWFemiYBQQFQCpg1F2OEbDRaNlHpMEaFpl8rpadHw2pVzTZ5WmayuuhJkLHwVAJzAhTcfkNI0TOlGB9L510usFpp6wTuPx44E5LrR5xMcems/uuafUYdl88VHKihsoGjqBCRknTE7TOKHfwKVvnYR5lpaWgJ/+FLjvvkCq2sauCuza5V4AOCMu7oe0PgsLgSIxlV5uKrkUgIj8sYj8QES+IyJfEpFXWvodFJHvisgTIsIpPSEuhBLFZWoaOphHyZF69mxQrD7ENnbV0iOEXGICkvr45hFsCnlXAI8AeKOq/gqAHwL4aELf96jqm12XJoRMNFGJkkQ8ZCQ0rYwK0RDRJAlbsmIzRelESYvMGdXI1lwKQFW/qqphvthvAliXf0iEECdTTrc7vHnq9Glz37SEdb4UfTwgWcKWHBgfj8btdoMWj8y1MaqRrUX6AG4D8KDlNQXwVRF5XES2FnhOQsaTJMnR6wWx/vEcxAsLwMmT9vcVKbSLCh6Jj2nVquE+U1OBcijZyD43t5yw9dixoJ0/v+yiSKIJKaqzYLjbKxGRrwF4jeGl7ar6wKDPdgAvA7B9Iu9S1SMicgmAR0TkB6r6V5bzbQWwFQDWN/3uEVIWa9bY0yrYQlGS7A3hbynNpJSX6elAYCcpoihRRbJ9u7kAzblzwF//dZBGO1wVhUZ2oBHpI3bsMGe3aHqK6txhoCJyC4APAtioqqnhByLycQAnVPWTaX0ZBkomFte8OlGSEs70+8HfMkNEWy3gF34hOH6nA/zsZ8EUutVKTmwTjtk3YQ7QqNjMhYVAh9miXKvSU5WFgYrIdQB+D8D1NuEvIheKyEXh/wDeB+CpPOclZOw5ftzvecBub+h2A+mTZuien89nJjp/Ppj5qwZ/Q6GfJPyjdRWzrPgbZGQPTUj33RfovqWl5kcE5fUB3A3gIgRmnSdEZBcAiMgvici+QZ9XA/hfIvIkgMcA/HdVfSjneQkZb7IYlZMSy4UkGbrf+c7kGXioNIrKKzQ1FYzNNdzVxJo1jQvAH6mIINctw3U0poIgE0vWHEZZUim4VCeL5unJU80sPjbffEfx1moN1yAos96DA7bUSKZyEGUA5gIiZAyoIi/O/LxbMjfTuX2VgCnZW1llMU3nSirLVeB9rjvXHRUAISSdft9N+IuYhaPP7N02Ky+rIprIsAK1ZUSNjyFn5rc6EtBG8VEAzAVEyCiTx/69bVsgn9II5Vjcm5nkVDbtpAKC6CaRoK1dW97O5TVrhnMz2ArWxO9BToN9UomHpsFsoISMKnlKKy4sADfdlNyn3Q5i8ON0u8DMTBCB4xrjuLAA3HprkP8nyqpVgRM6a/1DIHBKnzmz/LjTCcJRbQLfhTw1GWuG2UAJmQTyhJsk9REJ9g3YBODSkn/Ws+3bh4U/YN745cvtty/XVmy3gS1bksNlo9jCXidkEyoVACGjSp4ENEl9wvKVrkLw1KlA6CYpgaTz5ZlpX3BBkC46XKmcOwfce69dsHe7K20zd9zhXpOxYeGmheDqLKij0QlMSAJ5wk1s7+12l/v4hmgmeTqTon3a7XIcwa7jc4m2qtuz6wEYBUTIBJBHKLm+1zWSJk359Pvm0pTT00EYpmksGzcWFyXUbucT1jYF1m43rtQoFQAhk0KevQJZN40lrQqSdjv1+ysVSLdrr4xuUgp5Wt5dWL6F72vERwEwCogQ4sfCQmDzN0UIFZWcLWtqCBt5x+U6ngYkp2MUECGkPObmgtTMrs7TLORJ8taKibUixpVWMiykQcnpXKACIIT447rbKWvkTFJm0+g55+dXPu73gc9/vvhdWPHrDcNObeMelYghV1tRHY0+AEJGmDKc1AXn7clM0rXVHDEEOoEJIbWTJUw16gzudoOW5Biu0/Fqc6LXnA3ORwHQBEQIKQffjWphaotwl/HSUlBZ5b77Asfqvn35E+0XaZqJ1laIFg4eoQrxVACEkHLwLWpjS22xbVtyFE4oWOPC/c47hx/HE8SVUaprhCrEUwEQQsrBFDkjEghe0+zbNkMOcw/ZWL9+ePVw6BDwmc+sfLxrVzWlumyV2RpYIZ4KgBDih6sZJRo5AwTCXwf7jkyz7ywz5FCwmlYPccJzxynaNDNC+aC5EYwQ4k7WFNQ2E05045Tp2En0esupqFstu4B3OU7Nm7eKhBvBCCHlkDUFtYtj1DRz7nbN7wuFdqh0XFcP8SyhDTXNVAUVACHEnawRLq6O0Xhkzc6dbvZ0l526nU6Q/nkETDNVkUsBiMjHReSIiDwxaJss/a4TkadF5FkR+UiecxJCaiRrhEtWx6irPd3UL75LePdu4J57zKGbE0ouH4CIfBzACVX9ZEKfNoAfAngvgMMAvgVgs6ruTzs+fQCENIy8ZSi3b/crJUm88fEBrCp7MADeDuBZVf0xAIjIFwDcACBVARBCGkYosLMI8rk5CvyGUYQP4MMi8h0R2SMirzK8fimA5yOPDw+eMyIiW0VkUUQWjx49WsDwCCGFYtsBS0aOVAUgIl8TkacM7QYAnwHwDwC8GcALAP593gGp6m5V3aCqGy6++OK8hyOEEGIh1QSkqr/uciAR+c8AvmJ46QiAyyKP1w2eI4QQUiN5o4BeG3n4zwE8Zej2LQCvE5HLRWQawI0AvpznvIQQQvKT1wn8RyLyZgAK4CCADwKAiPwSgM+q6iZVfVlEPgzgYQBtAHtU9Xs5z0sIISQnuRSAqt5sef7/ANgUebwPwL485yKEEFIs3AlMCCETChUAIYRMKFQAhBBSFg0vDl/FTmBCCJk84mkzwhoIQGM2z3EFQAghZZA1dXaFUAEQQkgZjEBxeCoAQggpg6TU2Q3xDVABEEJIGdhqIGzaNFzAPl4fuSKoAAghpAxsxWz27WuMb4BF4QkhpEpsBexFghTbOWFReEIIaSpZy2qWABUAIYRUSdb6yCVABUAIIVXiWui+ArgTmBBCqqYh9ZG5AiCEkAmFCoAQQiYUKgBCCJlQqAAIIWRCoQIghJAJpdE7gUXkKIBDdY/DwFoAx+oeRAY47uoZ1bFz3NVS5Lh7qnqxS8dGK4CmIiKLrlutmwTHXT2jOnaOu1rqGjdNQIQQMqFQARBCyIRCBZCN3XUPICMcd/WM6tg57mqpZdz0ARBCyITCFQAhhEwoVAA5EJG7ROQHIvI9Efmjusfjg4j8roioiKyteywuiMgfD+71d0TkSyLyyrrHlISIXCciT4vIsyLykbrH44KIXCYifyki+wff6W11j8kHEWmLyP8Wka/UPRZXROSVIvIXg+/290XkmirPTwWQERF5D4AbALxJVa8G8Mmah+SMiFwG4H0Anqt7LB48AuCNqvorAH4I4KM1j8eKiLQBfBrAbwC4CsBmEbmq3lE58TKA31XVqwC8A8CHRmTcIdsAfL/uQXiyE8BDqvoPAbwJFY+fCiA78wA+oar/DwBU9cWax+PDnwD4PQAj4wBS1a+q6suDh98EsK7O8aTwdgDPquqPVfUMgC8gmCw0GlV9QVW/Pfj/JQTC6NJ6R+WGiKwD8JsAPlv3WFwRkV8E8G4A9wKAqp5R1Z9UOQYqgOz8MoB/KiKPisj/FJG31T0gF0TkBgBHVPXJuseSg9sAPFj3IBK4FMDzkceHMSKCNEREZgH8KoBH6x2JM59CMKnJX1S3Oi4HcBTA5wamq8+KyIVVDoAFYRIQka8BeI3hpe0I7t0aBEvltwH4cxG5QhsQVpUy7o8hMP80jqRxq+oDgz7bEZgqFqoc2yQhIjMA/iuAf62qP617PGmIyG8BeFFVHxeRa+sejwerALwFwF2q+qiI7ATwEQD/tsoBEAuq+uu210RkHsAXBwL/MRE5jyCfx9GqxmfDNm4R+UcIZh1PiggQmFG+LSJvV9W/rXCIRpLuNwCIyC0AfgvAxiYo2gSOALgs8njd4LnGIyJTCIT/gqp+se7xOPJOANeLyCYAqwG8QkT6qnpTzeNK4zCAw6oarrL+AoECqAyagLLz3wC8BwBE5JcBTKPhSahU9buqeomqzqrqLIIv4FuaIPzTEJHrECzxr1fVU3WPJ4VvAXidiFwuItMAbgTw5ZrHlIoEs4J7AXxfVf9D3eNxRVU/qqrrBt/pGwH8jxEQ/hj87p4XkdcPntoIYH+VY+AKIDt7AOwRkacAnAGwpeGz0lHnbgAXAHhksHr5pqreUe+QzKjqyyLyYQAPA2gD2KOq36t5WC68E8DNAL4rIk8MnvuYqu6rcUzjzl0AFgYThR8DuLXKk3MnMCGETCg0ARFCyIRCBUAIIRMKFQAhhEwoVACEEDKhUAEQQsiEQgVACCETChUAIYRMKFQAhBAyofx/OryD63lkArkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f575a5c9828>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "cov=np.array([[ 4.68, -4.32],\n",
    " [-4.32,  4.68]])\n",
    "mean1 = [0,0]\n",
    "mean2 = [0,4]\n",
    "n=500\n",
    "x1 = np.random.multivariate_normal(mean1, cov, n).T\n",
    "x2 = np.random.multivariate_normal(mean2, cov, n).T\n",
    "X=np.vstack([x1.T,x2.T])\n",
    "y=np.hstack([[0]*n, [1]*n]).T\n",
    "plt.scatter(X[:n,0], X[:n,1], color=\"red\", label=0)\n",
    "plt.scatter(X[n:,0], X[n:,1], color=\"blue\", label=1)\n",
    "plt.legend();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The two clusters seem to be quite separate. Let's try naive Bayesian classification on this data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.948424Z",
     "start_time": "2020-06-24T19:29:18.945733Z"
    }
   },
   "outputs": [],
   "source": [
    "model = GaussianNB()\n",
    "#model = MultinomialNB()\n",
    "model.fit(X, y);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:18.952919Z",
     "start_time": "2020-06-24T19:29:18.949630Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Means: [[-0.11860779  0.12869494]\n",
      " [ 0.03765647  3.98121459]]\n",
      "Standard deviations: [[ 5.18579166  4.69249675]\n",
      " [ 5.0005392   4.84702521]]\n"
     ]
    }
   ],
   "source": [
    "print(\"Means:\", model.theta_)\n",
    "print(\"Standard deviations:\", model.sigma_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:19.115073Z",
     "start_time": "2020-06-24T19:29:18.954208Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy score is 0.822\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJzsnXd4FNX6x7+zvSX0JqE3RcRCEEWvXgVRQFGxY0dFEHtX9AfW6/WqF1FQsaKiYkOvoljAAjaaoBQRQZDQOyE92ff3xzfjtpnd2WST3ZDzeZ7zJJmdOXN2dnPec96qiQgUCoVCUf+wpXsACoVCoUgPSgAoFApFPUUJAIVCoainKAGgUCgU9RQlABQKhaKeogSAQqFQ1FOUAFAoFIp6ihIACoVCUU9RAkChUCjqKY50DyAeTZs2lfbt26d7GAqFQlFnWLhw4XYRaWbl3IwWAO3bt8eCBQvSPQyFQqGoM2iats7quUoFpFAoFPUUJQAUCoWinqIEgEKhUNRTMtoGoFAoFOmgrKwMeXl5KC4uTvdQTPF4PMjJyYHT6axyH0oAxKG4GCgsBBo1AjQt3aNRKBS1RV5eHrKystC+fXtoGfjPLyLYsWMH8vLy0KFDhyr3s9+qgFatAoYNA3JygN69genTrV+7bx9w0UVAw4ZAq1ZA587ArFk1N1aFQpFZFBcXo0mTJhk5+QOApmlo0qRJtXco+6UAWL0ayM0Fpk0DNmwAFizghP7f/1q7/qyzgHffBUpKgNJSYM0aYMgQYNmymh23QqHIHDJ18tdJxfiSEgCapr2kadpWTdOWhh1rrGnaF5qmrar82cjk2ksrz1mladql1R14PO6/HygoAILB0LHCQuD//o9qnXisWQN8+y0n/3BKSoDHHze/7osvgD59uGvo3Rv4/POqj1+hUChqg2R3AK8AOCXq2J0AZolIFwCzKv+OQNO0xgDGAugD4EgAY80ERSqYOxeoqIg9rmmc4OOxdi3gcsUer6gAfvvN+JpPPgHOOAOYNw/Ys4c7jjPOAD76KOmhG7JlC9/T5s2p6U+hUNQNZs6ciW7duqFz58545JFHUt5/UgJARL4FsDPq8OkAplT+PgXAGQaXngzgCxHZKSK7AHyBWEGSMtq1Mz5eVga0aBH/2u3bgb17Y4+7XMAxxxhfc8st3GGEU1QE3Hpr4rHGo7wcuOwyvp9TTwXatwcuuYTvQ6FQ7N9UVFRg9OjR+PTTT7F8+XK8+eabWL58eUrvkQobQAsR2VT5+2YARlNsawDrw/7OqzxWI9x9N+DzRR7zeKjHb9LE/LqdO4HLLzd+zeMBbrrJ+LXffzc+vmoVIJJ4vGbcdx/wzjtUP+3Zw5/vvguMHVv1PhUKRQ0wdSpXaDYbf06dWu0u582bh86dO6Njx45wuVw4//zz8eGHH1a733BSagQWEQFQjSkP0DRthKZpCzRNW7Bt27Yq9dG/PzBpEtC4MeD3A243MHQo8PLL8a977z1jd0+bDRg1CjjgAOPrzHYVLVpUz3306aeNdxaTJlW9T4VCkWKmTgVGjADWreOKb906/l1NIbBhwwa0adPm779zcnKwYcOG6o42glQIgC2aprUCgMqfWw3O2QCgTdjfOZXHYhCRySKSKyK5zZpZSmgXwaZNwOTJNPYuXgz8+iuwdSs/i+hdQTT5+cbqFREgXqzFPffE9u3zAWPGJD38mPGYHa/OzkKhUKSQMWNiV2qFhdWfAGqBVAiA/wHQvXouBWC0R/kMwABN0xpVGn8HVB5LKS+8AHTsSFXNzTcDXbsCX34JZGdbu/6UUwC7Pfa41wsMHmx+3ahR9Dxq2JC7jQYNqL4ZPbpq70MnN9f4eK9eKjBNocgY/vorueMWad26NdavD2nO8/Ly0Lp1ijXnImK5AXgTwCYAZaAe/woATUDvn1UAvgTQuPLcXAAvhF07HMAfle1yK/fr1auXWGXdOhGPR4Rr41DzeETWrrXcjVx3nYjfH7re7xe56CKRYDDxteXlItu382cqmD+f97fbORa7nX//9FNq+lcoFMYsX77c+snt2sVOPACPV4OysjLp0KGDrFmzRkpKSqRnz56ydOnShOMEsEAszulJpYIQkQtMXupncO4CAFeG/f0SgJeSuV8yvPeesVokGORrN99srZ8nnwROOw2YMoVeOBddxNW/lRW33R7fyJwsubnAokXAo48CP/8MHHYYcPvtQLduqbuHQqGoJg89RJ1/uBrI5+PxauBwOPD000/j5JNPRkVFBYYPH46DDz64moONukdKe0sjZWWRgV86Ism5TWoacNJJbJlA165UbSkUigzlwgv5c8wYqn3atuXkrx+vBoMGDcKgQYOq3Y8Z+00qiNNPN9bfOxx8TaFQKGqMCy9kFGkwyJ8pmPxrg/1GAHTrBtx5Jw22djub1wvccQdw4IHpHp1CoVBkHvuNCghggNSZZzIJHACcey5w6KHpHVNNs3UrcO+9wIcfMlitf3+qHxs1YtRwp07pHqFCochU9isBAAA9e7LVB/Lz6RK6ZUvIzvHii/zpdAL/+Q/tB8OGpW+MCoUic9lvVED1kVdfZfoKIyN3WRmjhq+6ivUNFAqFIholAOowc+bEBiBG43AAX31VO+NRKBR1CyUA6jBduzLyOBFG6a3TweLFjGO4+Wbghx/SPRqFIrMZPnw4mjdvjh49etTYPZQAqMNcfXX8HEUAE9mdcELtjCceDz8M9O3Lojrjx9NYfcMN6R6VQpG5XHbZZZg5c2aN3kMJgBrk229ZQ6BhQ+Dww4GPP05t/61bM9fRgQdylW+zsXm9QCDA9sEH6d8B/Pkn8MADtEkEgwzOKyykgXrhwvSOTaFIBTWQDRrHHXccGjduXP2O4rDfeQFlCl99xRQSRUX8e/Fi4LzzOOldYJZQowr06QOsWAFs20Z10PbtLEeZnc36B4FA6u5VVWbMME6lUVwMTJ9OTyaFoq6iZ4PW7XF6Nmgg8+PB1A6ghrj99tDkr1NYCNx2W82kcm7WjJN+x47AyJF0/TSb/OfMAY46irUSunZNzWolHi6XsQCw2xm7oFDUZepwNmglAGqKpUuNj2/enLgwfU3y3XdMe/3TT/ySrlrF1crEiTV3zzPOMM7T5HAA55+f+PqyMib0u/tuxjkot1ZFJlFD2aBrBSUAaoicHOPjgYA1z51UUVoK3Hgjo4NtNia5M1qt3HMPC98D/FlSYq3/deuYsTTe+c2bM2ZBt034/Vz5P/EE0Llz/P537WJg32WXAf/6Fw3HHTqYl+FUKGqbtm2TO55JKAFQQ4wda1wl7PbbORHXFkcdxRTXRUVUPUWrpXSKi4H16+lZFAhwrD17At9/b3z+li00cB94IL2MmjWLX3LznHOAvDzuNCZMYL6skSMTj//ee4E1a0Kr/oICYMcOCgSFIhN46CHj//VqZoPGBRdcgKOPPhorV65ETk4OXtTD/FOJ1cIB6WjJFITJRCZOFGncWMTtFsnKErnvPmuFZVLFBx8Y16kwaoGAyKBBsUV1/H6RlStj++7dW8TpjDzX5xOZOzf+mHbtEpk6VeS111g8JxHNmhmP1+kU2bu3as9FoUhEUgVhROT111n/RdP48/XXa2RYMVS3IIzaAdQg11zDZG15eUzZ8H//l7pSjuXlXBXHMyhb1ev7fFxRz54da58oKaHvfjgrVwK//BKbgqKoiGodM955B2jVCrjiCrZWrRLXOoi3W1JlMRWZQh3NBq0EQE1jtwNNm9LgmQpKSlhrOCuLGT+7dAG++ML4XF2nb4bNxn5uv50uo0a2ifLyWIP2Sy8Z6/xFgA0bjO+1eTOzkxYXs5WWUoBcdRUN02ZcdFHsuOx2qp8ywcVVoajLKAFQx7j0Uurai4s5Oa9eTS+bn3+OPffCC81X0C4XsGwZsHs3cNddNMYaTepOJ3DkkaG/y8uByZON+7TZ6GFkxPvvm1dmu/RS4+MAMG4c0KMHJ3unkwKrZUuW7FQoahKpCX/tFJKK8SkBUIfYsoWRvdGG3KIi4JFHYs+/8ELWEY5WlTgc9AzKzmb9Y7+fO4mGDSP98jWNnjvh9ZQ3buTq3Yzrr+fP4mLGGyxaxJ1BcbH5juTPP82T2gUCwLx5fN8PP0zht2ZN3fCwUNRdPB4PduzYkbFCQESwY8cOeKoZSKMigesQ69Zxgo5eqYswGjgat5tePK+8Qt38tm0MFLvvPmDAAE76eXmhiXnrVk74zZvTvnDccbyuTZtQn02aGPv0A8AhhwCNG1PXf8UVFCDBIK95+mnz95VIl2+zAf36sSkUtUFOTg7y8vKwbdu2dA/FFI/Hgxwzf3OLVFsAaJrWDcC0sEMdAfyfiIwPO+efAD4E8GflofdF5P7q3ru+0aWLsZrGbo9U04RTUcFUDGvXUn2zcCH1/SeeSHfK8FW5PrE/8ghw+eXG/fn93Fm88UbkTsTno7rmt99Cun6dffsoEA4/PFZVpWnU50e70SkU6cTpdKJDhw7pHkaNU20VkIisFJHDROQwAL0AFAKYbnDqHP08NflXjUaNgFGjYidLr5f1kMPZtInlMQMB4KOPqLbRJ/hgEJg1iz710RQUJA6ymjiRqSY8HgqEBg3oKdSvH4WLUaRzcTEjedu3D2Uw9XoTxw8oFIqaI9U2gH4AVovIuhT3q6jksccYEdu6NQXBySdTzRMeUVtWxtTLH39s7iaqe9RHEwhQZXPOOZys+/dnVlOdYJC6/SOOYCK3X36hamnkSF6zerXx/crLuXtZtYqlKo8+mu2llxjZq1Aoap9U2wDOB/CmyWtHa5q2BMBGALeKyLIU37te8PvvXDFv386/N22K9fT53/+o3ikvT65vp5Or+XHjQqmb161j8ZZXXwWOP552gfXr2bfDwUjg2bNpHP7mG3P7QGEh0KIFPYjuuos7kooK5iQ69VTgzTcjbQHBINVJHg/tFqlGhILzu+/oVTR0qHIrVdRDrEaMJWoAXAC2A2hh8Fo2gEDl74MArIrTzwgACwAsaNu2bfXC5PYzCgpEmjZltKG+htc0Rhvn54fOe+ghEZvNehSw3lq3jo0E1tsBB4icdVZs9K/bLTJqlMisWebX6m3wYONz/H6RTz8Njf+rr0RateJxr1fk4IONo5GrSmmpyMCB7N/hYBR0w4YiP/+cunsoFOkCaYoEHghgkYhsMRAye0VkX+XvnwBwaprW1EQgTRaRXBHJbdasWQqHV/d57z3q0sNVNyJcTb/zTuhYjx7UzSfLpk3mmUq3baMrZrQvf0kJ8PrrTH2bKMvp998bVzArKAiNf8MGYNAgjqWggDuR5cu5+zCLI0iWyZO5WykoCEVU794NnHVWzaTqVigylVQKgAtgov7RNK2lpnGDr2nakZX33ZHCe9cL1q839pcvKOBrOoMGUd2SLGbqG4DqHjM//uJi8/TX4WRlGR/Xq5ht2wace25snIMI3+NnnyW+hxVefNH4OW7erLKMKuoXKREAmqb5AZwE4P2wYyM1TdPzPZ4NYGmlDWACgPMrtyqKJOjd29hd0u+PdAN1OJiJ0Gr6iUTZSb1e4B//iH+9lRz9u3YZr7AdDu4iWrY0zz5aXs5dQSowE3R63IJCUV9IiQAQkQIRaSIie8KOPSsiz1b+/rSIHCwih4rIUSJi8m+uiEe/fgy2Cg/+83iAgw5iYFc4LVsmrrbldALHHgsceqj5OW43/f53xNmvxYsMDic/n7sATeO9fT7+rKgA9uyJP/lWVNCzKRVcdpmxIG3ShEZthaK+oFJB1CFsNvrv33UXXSfbt2cit6lTgblzmSpCp3Nn89z/AHMBjR5Nl84JE2InRKeT9Ya3bAGefz5+X8ns5TZt4vnBIIVL+/aJk9YBnPwPPtj6feIxahSQmxvy+vH5KJjefltlGFXUM6xai9PR6no9gJqmuFjk7LPpWdOgAT1yhg8XKS8XOfNMEbvd2BsnEBDp0kVk585QXy++yD4CAfZzyimRr48ZQ4+ZZD2L4rVwb6Z4zeUyzv2/e7fI9dezZkDz5iI33xy/RkB4LYaKCpHPPhO5916RSZMi36tCUZdBEl5AaZ/k4zUlAOJz3XV0k4wuynLffeaTtdPJYixFRexj3TqR/v15vsMh0qePyPz5sffauVOkQwdOxuH9VcXdNFkh8eyzseMpK6N7qNsdKSjatRPp2ZNupMOGiaxeLTJlikhODs/JyRF55ZVQP1u3inz+uUiS9T8UiowlGQGg8fzMJDc3VxYsWJDuYWQkIjT+Gqlmmjalzt7oo/X7QwbbkhIGWW3ZElLD2O1Mz/Dnn7E2hPx8FnCZOZPnvPOOdf1/ddA06udfeAH48kvgrbfoxVNcHN9uYLdT1QXE5i2aNAn49VcmqfN46GLasyejp5s0qdn3o1DUJJqmLRSRXEvnKgFQNykrow7d6OOz20N69ujj551HmwEATJsGXHllrAdPIAA8+2z8qkZLl1KPbrV4fKqw263ZDBLRuDHHHp4PyelkkryZM6vfv0KRLpIRAMoIXEdxOukRZEQwaLwydjgiSzb+8YexP/y+fXwtHp061f7kD6Rm8gdYojM6GV5ZGfD116E0G4n47ju633o8QE4O8NRTKpBMUbdQAqAOM2kS1Rl2O/92OqniMXP/bNYsMkCsZ09jd8hAIL5rKAAsWZI4fsAMp5PFZ8JLPfr9LP+4ejVwwglV6zcZ9GcWjcNBl9RELFhA19v58ykIN2xgRtZx41I6TIWiRlECoA5zzDHM73/ppQwSu+oqTkxGk5vNBhx1VOSxQYOAdu1CenKAv+fkMEFbPNasYYBYVWjdmiqkW29ljYNDD2U66SlTaJM47bSq9dugQaxQcjpjBaLPRyFjlJYiELCWnXTs2Fj7S2Ehs7XGc5lVKDIKq9bidDTlBVQ1nniC3kDhnjSBgMiyZZHn5efTY6h5c3oA+f0iV19tzSVy6dJYDySr7ZJL4vedl8exWO3P4WCiuUAg1rX0oYdE3npLpGNHusV26MC/N24UadEilJxO09jP8cfTLTQ8uZ4RuleRkYvtqlWJn59CUVMgTcngFBnCTTcxz/4hh9CjZdAg6qu7dw+dU1BAI+4jj7AUZHk5p7C2bVl4JhEHH8xaAeFqHCu43aHaBStXsrBMhw40vs6aRTVK587WA7J8PqBVK76HfftidfDvvcfV/gMPAO++CyxbRkN4q1bchdx1FyuV2e3cPXzzDXDbbXxWW7ea39csYriign0rFHUCq5IiHS2TdwBffinSqxdXwV27irz5ZrpHlBwTJkTuEvTm8fC1AQMYTBaepjmakhKRww5LbvUfCIhs2MDdSCAQGUfgdMbGGZi1nBymoe7cOTIWwKh5PCJZWSLZ2Wzffhv5Prp3j73G6eRuyIw5c2Kfn8/HYDSFIp1ABYLVLLNmGQdgPf98ukdmnZNOMp4s7fbICdXvF7njDvN+Bg60Pvm3aSPy/fe8bsgQ65HARq1lS07+ZtHO8Vp2tkhhIcexbZu50MnOFnnySePAOBGRTz5hRLWm8dxx4xiFrVCkk2QEgIoDqAK9e9PYGk2zZgyqqgv5ZC69lBk4rWS/dLuZJrlt29jXnnmGxlwjd9JoWremkXfBAnrOpCvzZnY2K5ydfjqwdy8/N6OANk2jUdzhoBpp+nTjDKvl5VQhJfu5L1/O3EiHHaaCzxSpQ8UB1DC//WZ8fPduRsvWJhs3smTjrl3JXXfttbHeMWYTmNMJfPWV8WuXXkodvhWPoA0bmHxOLzeZLkRCMQDZ2cA//2k8sYuEgsVmz6awM8LhSG7y37aN8QO9e7MITU4OC+pk8FpMsZ+iBEAVMFoJAzRI1lZd2aIi1rHt1Ak45RTggAO4Erc6ifTuzTQIfj8nQb+fxl8jF1JNMzcM+3ys63v33dwpVGf3U9W4gmQpLqYLrc6UKcYuoeEUFrKSWCo47zxg8WL2uWcPx/Pkk8xGqlDUJkoAVIEHH4wNoPL56MESPYnt28ecOa+9FpmuubqMHg18+iknj717+fOZZ5jCwSqXX87V6IwZwI8/0lMoPCZAx+kETj7ZvB+/n7sfXXteVWprV1BWxtiDlSt5T5fLmu/+0qWMmzArWmOFzZt5fXR5y4IC4L//rXq/CkWVsGosSEfLVCOwiMjdd4cKpNts9G0PTzcswiyTfj89UAIBeqOMH1/9excXmxdg79Spen2/8QbHnJ3NcbdsKbJoUeLr2revukE3Xc3hsO51FO0htGCB8edSXs7vwSuviOTminTrxlTau3fznBUrzGMcqvvZKRQiyRmB0z7Jx2uZJADWrqVb4EEHiRx5ZKzroc/HtMM6+fnG/+g+n8iSJdUby44d5hNXw4bV61tEpKBA5IsvRObOte7VcvDB8SfN6nj8WGlV8QaqThs0KPTef/xR5NBDuRBwuznph7uIut30FiooYBrrxo2N+xwxovqfnUKRjABQKiALvPgide3PPQesWAHMmxebCK2wELjjDv4rA1SrGOm0S0rogVIdGjWi50o0mgYcd1z1+gaozurfn3pys5w50YwaFd8QXJOeUfEK1tcUs2cD337LlBj9+zM3UjDIz3flykivqJISGutffz3kUWTE2rW1MnSF4m+UAEjAxInA1Vdbm2C2bw95lxQXh4RBOBUVsVkoEyFC/fO33/LaFSuYzTIarxf497+T6ztVjBwJHH98eu6dDo+i4mJg4EDWF7ZiPygoAL74gr+b2RC++or9Vgd9P6FQWEEJgDiUlDBVgNXVpd8fMg4PGED/cKNzzjrL+hj++otZO/v0YZK05s2Bc881nihatEhfUXO7nX7yfr/xa0ZeNkaul1UhXS6lhYWsxWzl++F0sv4xEH8BEG0ctsrGjcCZZ/I+ej3n776rWl+K+kPKBICmaWs1TftV07TFmqbFhElpZIKmaX9omvaLpmlHpOreNcXq1dbP9fnohmmzcQW2bRt95L3ekCrI7+ck3q8f/9bzz8+aZTyhi3CVuWIFJ5u9e/lz2TLjVd769bHFXWoTj4fZMKPVPcEg33P79nwGgQBVWDNnArfckjpBkA6srrY1jbukYNB8ku/QgcXpk6WkhJle//c/CqOKCqop//EPLjaMFiIKBZD6HcAJInKYGEehDQTQpbKNAGASVpM5NG8ev+Shy8WWlUUX0LvvZrDToYcCffuydGFFBX8fMYIr5Dfe4GQwZw7QsiWjUYcO5b0+/jiy/19+Adats74DcTjMawHUFqWlsa6kIhRy7dpRXVJYyER0Bx5IgbF3L20sZmha3YiujoemUdj362euMurWrWp9v/8+AwGjd0Ii/E4p91KFKVatxYkagLUAmsZ5/TkAF4T9vRJAq3h9ZoIX0FlnGbtcHnWUyJ49TJ1cVhY6/7DDYj1SfD6Rzz4LnbNnD91CjTyE8vJC582aJdKggbHHSHQxdo8nfvKy2mLIEOPxalrkmO12kQMOEFm4UGTXLpGzz65dL55UtGQ9m+z2+NcMHly1Z37vvfHv27597DXbt4u8/z6TGoZ/fxV1H6TJC0gAfK5p2kJN00YYvN4awPqwv/Mqj2U0U6awOIrbHVJfPPoo0y9kZ9MjR1dhrFzJnDnRK/bCQmD8+NDf06cb36uiIlSvF+Aq2WgH4vFQx+t2swiKxwMMHhx5j3TRqZOxvl8kcoVaUUG99bHHcidUWFj3VEFG1dTiUVFhrjJK1jYUzsEHG9tedKLzNI0fz/QTl11Gu0Hr1vRiUtQ/UikAjhWRI0BVz2hN06rkkKhp2ghN0xZomrZg27ZtKRxe1fD7Gcm7cSOrb+3YwXzxRuzcaT6Jhb+VXbuM9cAlJexfJzsbeOihyInG6+U/7GefAXl5wCefsH7vu++mXv1TXEzvo3h58aMZPTpxWoVwior4vr/+GujaNekhppXqeuzoeDxUhw0bVrXrzzwTaNrU+DWHg2pGnR9+YN4hPYI8P5+fr5nTgmL/JmUCQEQ2VP7cCmA6gCOjTtkAoE3Y3zmVx6L7mSwiuSKS28zI2T1NNG5MHa1RqgSdww4z1td7PMCQIaG/TzrJOEbA72fxlnBuugn46CP+Ex91FEsRLlrEc3fupF69dRX3UTt2AMOHU9A0aEB3Vz2p3MSJNNT27cvcR6efbi3RXadONEYmm92ysJD2jhtvtB57kG6s2mY0zfw9aRrw8MP0Jkq2uI6Oy8V8TLpzgY7XS9vS/feHjk2ebGyDKCqim7GinmFVVxSvAfADyAr7/XsAp0SdMxjApwA0AEcBmJeo30ywASTL5MnU5eu6Xq+XZQh37Yo8b8SIyEhhv1/k9NNj00kY8dlnTNHg91P336ePyF9/JTfO0lLm09fTWQCMLu7eXeTDD2OLnbjdHJ9VDjmkajr1/HyR335Lv34/la1rV5HXXot9ph6PyDnn8Hm9/75I794sdHPRRSKrV1t/1gUFLM4jwu/BmDG0XT35JO1N4Zx+uvEYs7NFPvjA+j0VmQtqOxUEgI4AllS2ZQDGVB4fCWBk5e8agIkAVgP4FUBuon5TKQB27xa58koaX71e/oOEG1xTyXffiZx7rsixx4r8+9+x/4QinOg//JD/kIMHs6KYlbQLf/wRO5HY7ZzMrQgPnXfeYa6f6IkgK8s8rYPbzQIqVmjUqGoC4Mwzef0RR1i7xuViczjSP9GbtX79+J7efFOkWTN+/9xukYsvZlqPs8+ONJDbbDT+//ln/Ge8YEHI6cDlEjnvvNiFRjSvvmqcosTjsVYLWpH51LoAqKmWKgEQDIocfnhk/hy7nf+Mxx/P440aMcFbaWlKblljDBtmPMkEArGlDuPxf/9n3I/dbj55BwIiy5db6/+EE8wnRLc7cucR3nw+kV9+YX4dny/kUWW3x1Yrc7m4Ym7alL/Xdj4gq81m43fN4+HPCy6g585zz4n07Gl8jcMhctVV5s93/fpYTzJN40LAzKtnyxbuRA46KCQEbDY+5wkTrH93FJmNEgBRfPWVsdtldPN6uXLPVP7803ySy8oSeest63299pr5M9E0Y3fF7OyQqiERP/0Uv1ZvTo7xcZ9P5KabRA48kH87nSKtWtHFdflykfvvF2nXjseuusq4rrEuMDQtc3cGmpY4E2nXrubP9+67za/PzRV5+GG6f7ZoQXXj2LGh2shZWXxuxx0nMnwjbbuoAAAgAElEQVQ4PyvF/oMSAFE89ZR5+uTo5vGIrFuXktumnNtui/X/D5/wpk2zrgYqLOQkGm/VHC4EfD6uWJPhscfM+z/kEOPPxOczzrQ6alRs/w89ZCxkfD7GI4waJTJjhkjz5rU/waeidepkXo/YTJcP8DsSvds1EugNGzKFtWL/IhkBUC9yAXXrZt010e1m6oVMZPly87w3FRX06DnhBLr46VWmDj2U1cIaNKBnzrBhfH+rV9PzpF8/8yhbm401fPV6uCOMojviMHKksWuq309X2miXWZstlFEznMJC4KWXmGqjeXPGYgwdyucRfS5Aj5Y2bYBJk+hVVduZQlPFn38ywd7w4Zyywzn6aPPvdDAYGT9iFn8QDDJCW1GPsSop0tFStQOoqKDe00zvHL0DWLMmJbdNOQ8/TDVVovHfd59I377xz/X5qOufNo16Y6NzsrOpj68O770XMnpqGu978cXcqXz3HdUU+uu5uSJt2xqPxeGIXO3bbFRlmBVX8Xho8MzLM+8zepdTnVaTqia/X+SjjyKfazx1oNWWlSXy9tvV+3wVmQeUCiiW7dup33c6+Y9z7LHGbnlDhsTvJz/fuidMqtm+nUZEMzWQ3lq0sGbz0FuXLsZ9er0i+/ZVf9x5eSKPP07D848/RqqpgkEK3A0b+He0R0y85vUynUS8SdzrpSAzmyzHj6cRuaYm71Q1tztk9/jrL6qAqlLNLPr7vn176LP44gs6RbRrR0P1b79V/7NX1D5KAMQhGOSOQIRudEceyQnE6xW57jqRoiLj67Zto7um7nbYrZvI99+nfHgJ+esvkTZt4v9jG7l3Jtt8Pk7Ytc2SJbGCWX/mRuMcPNi6wDBqZ57JzzGewTqTmqZx51adyd9u5/f9qadCz/2ZZyJtMjYbFxHLltX+d0BRPZQASBK9jqsZwSD9raNVSIFA7RuMg8GQh4xRczpFevWq2sTQogX9zzt3FnnxRcZOLF5c+/7hP/4ocvTRnORatRK5+WZjbx+3W+See+KreBK1Ro0YgGVFPZhJzUzl5HCY74icTpHTTuPut00bPtvLLxc5+WTz+xx4IOtEmy2MFJmHEgAp5qefjHXNLpfInXfW7lg++cRct69pkT7eyTavl/eoqBC55RauCPWoZk2j+umRR0I7qHA2baKX0HPPiWzcmPr3rQuE8PFmZVF19OKL5u6giZrTWb0dRF1qbdqIHHNM8s8qEKCqLdloc0V6UAIgxbz1lrlaRY9cTRXBIOMWRo0SueGGWDdAs0AwQOSKK+h6mchQbNYaNOA9Hn/cfJLw+URuvz1yTC+9FBIWPh9/T9ZlNBF79ohceGEo4Cs3V+Tnn0OvP/lk8qoRp9O6e3B9b3a7yMCBqf1MFTVDMgKgXriBVpcjjjDOlOjzMZ1xqhBhQrZTTwWefRZ46im6AT7wQOgcs6RiWVlMOLdhg3GyL4eDyeTMEsc5ncAll/D3xx+PTSGsU1jIcellDdevZwbQ4mK+VljI32+4IbVFznX3UBG6re7bFznG669nxtVNm4Cbb7aWqtnvj5/cryax2eiebJQUMBOpqGBN46qWrFRkJnXk65deunThpOz1ho45HMyiOXx46u7z44+sGFZQwIkuGOQk9/DDocn0kkuMc78Hg0D//sAxx9BPPhqHg5Pm9u3Gk15ZGbB5M7BnT2RKaiMcDk78AKtRiRiP59134/eTDIMHM8toWRmF8W+/MYXxn3+GzrHbmbX1sccopLp3Z+yDWYru3bvj1+f1eoETT0zdewgnGGT9iHTVM64qdb0ymyISJQAs8sYbTMXcrh1zr19yCdMyN2yYunt88IHxylvTmPcfYODW5ZdzcnI4QitIj4e7htNOiy3K4nIxMGjpUq6i9bKN0buJDz8ETjmFhWjiUV7OgiIAJ2SjSSxe7dtk+eUXYPHi2OI4paWc6KPRNArmZcu4IwoX3EYYve7zMUf/O++kb5eQDuz2UGH56OMDBtS9oj2K+CgBYBGHA7jjDq7Et20DXnwRaNUqtffw+YxVPDZbaJLSNE5677wTipwFuGq/915Gy4aXUdA0TtjRk3RpaWyEbGkp8OuvVEP5/carPZ8PuOaa0C5jyBDj86ILkejs2we89hrr1P78c+RrmzczennsWKqVWrTg5HvWWcb3KCtLHLXtdvNeZhOX283o6NxcRg/37Ml6DRMmAN9/zx3FDTfEv0ddRNP4bDt25PPu1Yvv/5JLgPnzeTwQ4HcsK4vf9cmT0z1qRcqxaixIR8sUI3Bt8fvvxgZcr5dpg8O5/HLj4Cans+pGYIAeH2+9RX/8Aw6IfX3QoEgvoF27jO/XqFGst9D8+TQ0BwJ04fT5GHBUUUHvJt2AbHWsHo/Igw9ae7YXXGDsHpmVxdTY8Zg6teqeValoZsn5qtP8fn5uvXpFBoPplJUxDqRdOwbKnXuuyMqViZ/zTz8xYv2ZZ2K/syLMtqtqENcsUF5AdZfnnuPEFgiw+Xwi//tf7Hk9ehj/Y6ciJfIZZ4hMnGjsa65pLFavEy/RXrduobFXVBgLFL+fXkTZ2YnHFV1UvmlT61HZK1YYjzM7m4nxwpk/X+TGG0WuvZYptnfsqJ5QzeTmcomcemrs83r++UhPMD31hll0cEUFhazfz8/G5+Pvs2fz9bVrRfr352sOBwP4dHfhggKmpJg8OblCOApjlACooyxcyJVWjx4sIjJhgnExGRFWkjLyX0/FStHjiV/QJTc3NI5rronfl90u8q9/MeraLD3FgQdad98MBDi2YcOsB+GVldF1VY/21TS+x8aNmY8onHHjOHnZbKEI8RNPZLCYXpMg3ZN2TQiBvXtDz6C0lJlCo8+z2czTpb/9tvEuqVEjBhS2aBH5fXU4mAdq7lzuCvUU1R4PY1CSKW6kiEQJgDqIrgLRJ3Cnk6tTsxXXokWxvvper0jr1jU/YTRrFhrHf/+b+Hy7XeTrr1OTosLrFVm6NLlne+21xiUuv/gi8rw//oivgrLb+bkkk2epLjSbLXLlvWaNucqrdWvjZzxwoPH52dkid91l/Mz02gTRx/1+kZkzk/uMFSGSEQDKCJwBiACjRtEDSITHyspYhP2uu4yvOfxwukV27UrDsddLH/9t28zvE+41lIh46bN1L6HFi4GZMxP3VVFBA3NVi56HU14OfPaZ9fMLCmiwj/auKikBHnww8tiMGfH7qqjg5yJCA6nVFOOZTjDI1N06TZqYp9DetAk47zx6V4Wjf2+N+OsvGv+jKSyM9ewC+Jm98ELicSuqjxIAGcCOHfzHikYE+Ppr8+v69aMveX4+8Pnn7MPoHwqgV8+DD/Kft0mTxGMyc+F0uxmYNm0a0Lcv72uFPXuA9u2tnZtoXMm4Im7ebB48t3Jl5N8ejzUBabcDb73F53nJJYy/yFRXUat++3PmAKtW8ffsbOCcc4zdY4NB4L33gN69OUl37Mj3/uuvxgLRZgMGDjSOTTFyRdYxC0RUpBirW4V0tPqiAioqMs9G2bmztT7uv988p43LxRq9paU0uCVKoWzWjjhC5IcfWBayQQPr19ntrMqVqpw7J58c/1kUFoq88orI6NEiTzxhrs6w2UReeCF03dat1oy9fj+9pHQ2bqSuO9WeOqloXbtaq5fsdIoMGEBD7X/+Q+NtPDuQyxU/gZ7bzef05Zf8vnTsGHm+yyXSvbtxyhG/n55XiqoBZQOoe1x2Waz+2eejO108Skv5z3rttcaTl8NBQ215Oc9/6SXjCdFm4/Vmk3QgwGtFaKxORp9vVOaxOk3TWJfBiC1b6Lqo65x1w6LZ/T0e5hTKy+P1egGbQMDYHqBprJ8QbaRcvVrk/POteTPVdvP5KLD79o0vhHUB5vFUX5j5fKFnKkLhevbZfLZZWaxTvGsXCxJ5vSHhEAhQwCtX0apTqwIAQBsAXwFYDmAZgBsMzvkngD0AFle2/7PSd30SAAUFLPLh8XAS8XhEbr01vjfEN99wlZadbW6Y9Psj0zlff73xeU5n/BWd1yvy9NPsY9Uq85VybWXWnDPH+Jlcfrnx+zDyatEnPYeDz7tHD5Fff+XE9Prr3EVcfz2FR3Y2J67mzUUeeIBC1ygr6g8/xPdostnSk33U6eQio7buFwiEqpgVF9Nrze0OuRbb7SJXXknvoz/+ELn3Xi5UZswwfq4K69S2AGgF4IjK37MA/A6ge9Q5/wTwcbJ91ycBoJOXx0kkUQ7+PXvMJ/1AgBNW48b0vgln8mTjHYAVFcHJJ7Oc4/vvc0sfPZG53RRaVn3mvV6u1jWNk2bfvtYDwTp1Etm8Ofa5mKktrARTaRqfWfTuYsMGqiQOPTTk3x4IUIWxdWvkuaNHm9/H7WaQ1OefZ6a6KJXN4QipcW6+2fhztdlEjjpKuXymmmQEQLWNwCKySUQWVf6eD2AFAJOck4pEtG7NrJ2NGsU/b/p04+NOJ3DhhfRo2bKF2UTDueACGoSjjZ3xCqfr5372GdM4DB0K5OXx39hmo4HP4wFuuYWGYaNspEbj7NIFWLOGBj+9kHtxceJrASaB69ePYwjHzECsTzvxEKF30DvvRB4/4ACmrVi5kmMtKKBXy6pVwIgRkefGS+7WtClw0EF83/t7UrXycn4+AFNIGH2uwSCNxz/9VLtjU4RIqReQpmntARwOwOgjPVrTtCWapn2qadrBqbxvfWTPHmNPnbIy4O23gY8+AvbujX09EGDW0ZNOogeGpsX3fLHZ6K0Rfa/S0pAA6NIFWLcOeOgh4KabrI1fhLl+Lr4YGDSI+Wi+/97atQAnj3Xr+F7C0SedRJhNwAUFoUynACf6p59mbqDoSaysjII23PPq/PPN77llC3DZZRQCTZtaG2dd5rjjKNSNXEDDWb68dsajMMDqViFRAxAAsBDAUIPXsgEEKn8fBGBVnH5GAFgAYEHbtm1rapdU51m2LL6qxe2mesUskliEutZ4ffh8Iocckni773aH8slUNWeO01YuDntFUtcEAvRcOfhgBiJVpxhOeJ+ffsp2xhnsL55O326PTCXx3nvW7uPxsF9dNWK38/dnnqFaJJNrFDscVa/AFt38/vTU1t6fQW17AQFwAvgMwM0Wz18LoGmi8+qjDSAZRo2KP+F6vYzUjUe8ieaZZ5hCwaz+rN58PkaPioj06RM6bkO5XIfxshJdZANayrMYIS2wybgP5Mv9GCMuFBu+bqYzDx9bVfTq4XYMj4fjv+Yaa4JM03i+zptvJl9hzOEQycmhbSDc7jNvHutQV6f4e001r5dCt7rV1BwOphVJ1gYQDDJS+LLLWAXv22+Tu35/p1YFAAANwKsAxsc5pyUArfL3IwH8pf8drykBEJ9gkJ4WRx1lPkmfckr8Ps45J/Zauz2UIGz58sSr6latRCrKgyILF8oPj34rAFfyU3Cx7ENoqVgCh+ThAMnGbsN+jsSP8hEGSXusFgdKxYZy0TS6MLZqFendY7enxpDq89Gg3LEjYynmzzd7v5G7E4+HnkV6WoqyMpEmTao2BrdbZMiQ2M+mvJyeSL178/137ChySOcCOQyLpCm2GvYVwF7JwTrRkNxuKpmmaYx9eOwxkZ49rX8OTZqE4kdcLnpsRe9Q582js0GrViLHH8/yqNHf+UsuCQloTeNnWNu1uTOZ2hYAxwIQAL8g5OY5CMBIACMrz7kWdBFdAuBHAH2t9K0EgDXmzTP37Bk5Mv61GzeyWLjuURQI8J8vvAD4yy+HfOPD/9n1rI9fvL6ZM4HfL5KdLTPsg+VwLJBCxC4Ri+GUW/FvwwmiHf4UASQIyF4EpBhOeeiCX2X9evr3jxrFILYuXWK9fRpjm2Rjl2G/DpTIXPSVfPhlN7LlL7SWt3G2HIZF0qNH5PN44onIXdEFeF3Woq0IIKvQUY7F12KzBWXcuMg0yr//Xr2U0W4336MZ8+cWy2THSCmCW3YhWwrhkZdxqThQKoCIFwXyGoZJEVxyLcaLD/ui7hEMa5H3bo31MhEjZSW6yNc4TgZiBidXlEszbJGZOEnOwbSIax+4P8gv3ssvy20X5onfH9tvdDv+eL6XYNB41T93bqxqyecT+eCD0DnfzQ1Ka99OcaIk4jyvl5+BopYFQE02JQCsEQxySx7tyunz0VaQiKIikddeo7pnypTY9MgiXKl98AFz5993HyOLr7yycgV83HEx24hSm1sKYLx1+AQnxworlMpwPB97fteukQMpLBS5/nrJt2dLMVzyCQbIY7hRKqDJE7hBfMiPnBhsRfJv2x0x/QYBqYAmBY0OEHnkEYarCiOD/X6qrx7EnVKCyKCCffDJY50nxTyfbdvM1Wm6gI0XZ5GVxRgEQ558Usbj+ojdFMfilZea3S7LlonsG3i2FGtu2YBW4kaR6X2cKBEbysSJYgGCcgDyZBsaSwnsEe/xWjzJcWGPvInzJB9+eQqj/xYml9le4YPy+yXo88vUjvdI7uFlpvEWgMh118X5EublybBDfomZ2AFmDRURkU8+kT3+VlIGm5TCIVNxvnhQKAB3ZBMmJPqm1w+UAKiHbNgQ8qMPBJix06iOQMrZtMl05gtqxhFPy3Bg5Qq1onJSKpbG2C7r0Cb2fE2LrCzSv3+EYryicjLXJ/X7cI8EsFd82Cd+5Mv/Ydzfr5s2r/fv/BK7dom0922Rlegs5TDWbZQ3aWa4hB08OFZn7/cz2nXNGpEbbjBOw9EVv8ln9oES9HipJxkz5m+BJK+9Jru8rWQLmhmOJR8Bbh0qP4MPcZrpTkhvh2GRPIOrpQd+kfdwulQYnLQXAfGgUDRUyH24VwSQQnikPdYIwF1DxDUul8gll8ibbxrvhLzeIL+PK1dS35WVxdSi99zDBYTHI3uRJbvQQC7EqxHX2mwiJd/Nj1nhBAH5Dn3+fs4vvlgL3/c6gBIA9Zi8POrt9dQPNc7Klea6D5fLVEE8B8fIqfhQemCJXI/xkgeDajEAl83z5/Ner70WfyKvbMVwyV/IkWIkYUH1+USmTxfZuVN2du9rOvn/3YqKYh7Frl2sHRAezT1unEiwoJAK8yOOkLVHnCkNfSVit1Nl0hIbZRcaSEX4c/J6Rc46i5126SIF8MbsREICUGNu8MrcHD+htwSwN+7Qe2KxCLiDqDB5n7uRLYdhkWRhj0zH6X8LhYsxRQCRA5AXe53NJsWTp0jbNkFxOEIqIQdKpItjtZRPfSs2V4ZeeCHs2D745Ch8//eh7GyR4LH/MBxnEJBO+F18vsTBk/UFJQAUtUd5OfMjGE3+I0eyEojVSdio6T6m27fXfPhsPB1NeNM0+nua8OefLDSzZ48wWVOvXhGW5bXeA+WKLt9Ihw4iU3LuknKnwQ7K42GCoUql+FIcZDiWzS16iuzb9/ezCQLSHUvFVmkbiG5eFMijuCXheyyCS9pijXTFCimrVA8VwCsn41Nxo0huxaPG1zocsqnV4XKO72Nxo0g8KJQL8ZpsR2PLlvsKaDINZwvAtz92rJjm8ggCMsb5b/nkk1r6vtcBlABQGLNlC2tOTpgQ8ttMBR9/zP9UfYvu9VLxvXUrJ8qqOo273SL/+AdngAEDqtZHTQmNs8+29mzefts4Z4fXy2o//fsb95+dHXLxAmQnGkgBvFJaORmXwSZFdp+Uzf6WiZHCdE8b0EqOwveVLrUhw28AeyUXP5naZsIn1QU4XK7AZNmBRn8f2wu/HIKf5TR8IHtQs1Vx5tn6iNfL1CLl5UJbkMlYi15UqUPDUQJAEcvbb3PS0dNjJlNR3QrLlnHFP2AA8wnv3s3jFRWhNJDVmRQyLXlOdnZsoqVwNmxgYeHGjY2v9/lYePeWW4x3HrqA+Pzzvw3sQUDy4ZNdaCC7Dz9e5JdfeK8XXjAUsn8hR2bjODnnwCUyPGemvIcz/17NJ9t0AVAGe2KbSjVb0OWSnZ16SWFOF3obzJ7NxEJG5zudzKSo+BslABSRmFU193qZCzlZioqY4nHsWPpj2mxUc7zwgsj48RQ2uo585kwWNUj3hF0Tzeczrk+5fj2NufGiuLKyaKVfuzZ2h+DxiJx0Epe+ubmR/Wgaw7PDjdA//hj/Xm43jRFmQjT8uNudnnSl0WMJH5PPx+/UJZeExqZpnPw//rhK/xL7M0oAKCKZMsVYDWGzidx2m/V+Fi5kVRi7ndcapRB1Ojm5NWtGo22qcgZkYrPZRC69NPY5XX11/PBpTaPdpLSU5y9YwIle0zgBX3UV9foffWT8ufn9zAWuEwyKtG1rfj+3m8UKjASAy8UYjgYN2Mcjj2TebgtgcEowyJ3m449TlRkeiGHGt9+KHH44v6tNmzLkej/PN60EgCKSF180rwJz443W+li/PrkqMHqIZqonAk1L7wrVqI0eHZrMRUQ6dDA/1+NhJNvy5bHPuLQ0cnIaM8a4D4dD5NFHI6+dNcv8nj4fd2ZGn0d0BFVpqXVjeG02l8vahB/OwoXGkWU335xcP3WMZASAqglcHxg82Djfs8cDnHuutT4mTmSuZKuIWCvs6nYzpaiV/Mh2O3D77cxpbVZMNh1MnAgMGRL6O17R5ZdfZl7pgw4KHdu5M1SQNzw1a5s2gM8X24fHA+TkRB474QTmETfCbgeuvjqUCxzg8/b5gBtuYDpXHacTOP30qle893j40+dLbc5rux3Iykrumvvvj81NXlgITJrEQtoKWJIS6WpqB5BCnn02VHvPZuNK6PrrrV8/eHDqV3TNmyeX9rJly+TOdzjYbLbaWdWOG8cV9PDhxq9rGs/Ryc8XGTqU7ykQoKvjK6+EXt+9O7b4sqZRlWEUrl1WRq8pfYfkcvEznzmTrweDIl98wfGNGGFeVm3HDla/CQTi1wmNbh4P83XceiszCTZtan5usmqm886LHOOmTcwG17gxo+vGjmXN0s6d+d0+8kgeN+orK8vYdrOfAKUCUhiyahXrGY4ZEwqussoDD1Q//aPe7HYGOiWT6rJly+RUPy1b0kunvFxk0qTqeyFZaU4nM8vF0//fcUfomZ52WqxA8/moztFZvFjkwANDnls9e4qsWEFV0Zw5zM8RXZZs4UKRf/2L73vbttjPcvt26tGvvlrk1VcNg9okGGT/L7/M78pFFyVOCwswu6DOtdfGfsaaxrwlLVpwIrYqCBo1ourmgANCfYZf63DE2qQMgsz+FlR79yb3/a9DKAGgSD3bttGzJXwSdjg46en6fn2S9XrN/7FtNk5AZobpVDVNoxB44IHMMUR7PDSMB4NcwZoJ1AEDYp//X38xwiwYpM6+XTtOoHrI8f33W/scf/6Z1+ifVSBAoRWebiOa556z/gy9Xi4wcnJCnjq6kPP5OOY33wwVb64pW1H49y38b5+POTn2Y5QAUNQMq1eLnHkmDcrNmkXmrCkv52r01ltpcLz0UmO1y2GH8fwPP0zOqFydlqjgcW02p5NG4nffjU2LoLcDD4x87l9+yWOaxmffsGGsgPX7Q6qeePToEXs/l4urdTN69kzuPUbvFFwukaOPDhVTNooI1osxpNrAHwiEclY3bkxBWWt5UtKDEgCK9LN1K9M46t5HXi8nvMWL+XpxsXmQVH1oTZsae2Y5HFTN6CxYYH2FrBdxMGPbNnO1W8uW5td17Fj995tI1WO308//1VdTWw5Nn0PqUeX5ZASA8gJS1AzNmrHY61NPsXL6Aw8Aq1cDhx7K191u4PPPgebN6d2RnU3PkUzy7qlJioqAHj0ij2kaizbffXfo2IMPxnqymLFrV/zXHQ7z11wu89eGDo1fONoKIvFfr6gAVqxg4WSjYteJ0LRYryOfj987/XVFDHG+EQpFNfF6gcsvZzOiVy9g40ZWgy8sBI49FujcGdi8uXbHmQ5KSoCff449fsghdBVt3Bjo25dCNNHkCfBZn3NO/HMaNgSOPhqYOzfSLdjrBa680vy6u+4CJkwASksTj6OqOJ1Abi4wbhwQDCZ//cCBfHbPPMPvUvPmwOOP87jCHKtbhXQ0pQKqhzzySOYFelWl2WxUcVU1qlZ31TUqIBB+DsDzundngr9Eqo7166may8qiWs7vZ0K64mJe+/TTjAj2+VjCa8ECXmdkO7Cq3gk3BJtdn50tMmOGuV0kXnM4GCl89NFUn02dGrJN1UOQhApIr9ObkeTm5sqCBQvSPQxFbVJcDPTvDyxYkFzgWabRsCEwYwZw0knWAuLM8Hi4Wg9Xi/h8VMt4PMDatcD69cCaNVTTtG8PTJkC9Olj3mdFBce2fj1w5JFA7948PmYMMH587HgPPZTt7bf5+YTTtCnVVmvXmt/P7aaKb9++kDrL6WSz24Hjjwf+8x+e06lT7D2Sxefjc/j+e6BBg+r1VQfRNG2hiORaOtmqpEhHUzuAekpFBb1e0r2KD29utzU/+PBV78yZIv/9b3LxDkYtJ0fkmGM4hlat6MNfUcHWtWvsuAIBZiM1YsWKULCY0yly7rl0Ac3Pjx8rYbNxdR5u1M/K4qrdSpCd3y9yzTUM0Dr2WLrDGuXkOeWU1BiB3e79PuWDGUhiB6BsAIrMY8sWGowzBX2V+t13QHm5tWvKyoBTTknN/fPyuMr++GPgxBNDBtlvvgE2bYodU1kZ8PzzwNixkcd37KBdYfduTpPBIDB9OlNTvPoqV+RmBudgkDuy4cNpn2jViqklZs3iTiSR4baggH3/9FP88956C7jwQuDLL2m4repuoKQEePNN2gEUpigBoMgs9u0Dnn3WmidIp07AX39VzWvE7bauYqqooMdSOvntNxo0mzenIDj8cGDdOk7k0ZSUAL//Hnv85Zc5oYZfU1ZGYZuXl9jIW1ICzJ/PptO4sbXxa1p8TyOdBg34/rZuBbZtA154AXjySeP3mYhgEHjjDeZm6tcvvhdUPSUlbqCapp2iadpKTdP+0DTtToPX3ZqmTat8/SdN09qn4r6K/Yzly4F27YB//9v8HE3jSvWxx4B336UHS7I0alT1ZGfppLycXlP9+nEi79XLOMmf3Q60bBnrTbN0qTdkRb0AACAASURBVPEKX4SeVxddlHiSdrsj//7nP7k7SeRm6fUCl14a/xyA73HWLO62WrQA/vtfYPt2oG3bxNeGY7fzuquvpndUmzZ0M1VEYlVXZNYA2AGsBtARgAvAEgDdo865BsCzlb+fD2Calb6VDaCeoUdsxtPt2u1MiSAisnFj8vpiTWOUciqDjarTquLxFAiEahKblZR0OJhs7z//YWUxEZGnnjIOKvP7RebNYyK76683v6/Px0CtaJYvp2dRIBBKTeFyhZLJeTyVhX0TsGRJKEeQ3s9//sPXhgwx/z7oNgg9VbjHE/tcNY0R2PUgIAy1GQkM4GgAn4X9fReAu6LO+QzA0ZW/OwBsB+iBFK8pAVCP2LTJ2qTsdjPJmc7gwclN5g0bUnBU1tpNW9MnpC+/rFounOHD+f4vvDC+0NQ0TsIXXCCycycn2PDUGB4Po2XPPZeTeP/+zKoZHv2raTzvvPPM0ygEg0wa99VXzFS6axczm06aZK3+dHk5DdxGQmfOHPYb/ZzsdhZ72bVLZOVKkWnTRB57jMV1jJ6F31+1Cnh1jNoWAGcDeCHs74sBPB11zlIAOWF/rwbQNFHfSgDUIzZvtj6Rh6dU3ruX+YkSedrY7UxxvHmzyOTJIq1bp1cAjB7NSfPzz6uWqdRuZ0pjq+U2fT76+K9fz4ydfj8zbF50EX8PXzH7fCLvvMOkcxMncsfUvHkotuHxx1O/kv7mG+PcUJrGMYpw/D4fdwc+H6vTGXk79e1r/Ayys0Xmzk3tuDOQZARAxllFNE0bAWAEALRNVu+nqLu0aAF07w4sXsx/VzP8fuC440J/Z2UB77/PoirLlgH/+x8wbx790w85hEZlux049VRgwwbaGNIdX6BpQMeOwMEHV10vXVEBPP00DeF//JH4fL0QyujR9OfXGTSIHjrR515/PZ/XsmXAO++EYgN27gTuvZf2hVtvrdrYjcjPN7YjiIRSXIwcCZx9Ng3czZoBBx5o3Nf55/N7ZBR/0aYN8MEHLMCzdSttCxdeaN2YXR1KS2m7ycrKnNQUViWFWYNSASlSxW+/MUmanjFSrz0cvjIdMCB29fnVVyIDB7JY+k03UcUTzauvZk6tW58vNbUVTjmFunWr53foEPtcGjUyPtfhEHn/fUYYG73esGFqa+vu3m28E/L7WdL01lv5u8PBVNgffWTeV1ER1Vp6unGHg30PGRJb9N7j4c5Aj3hOFfn53D0NHcr003rhH4eDJUFnz07t/cJALauAHADWAOiAkBH44KhzRiPSCPy2lb6VAKiHFBQwSOjhhxlINXky9fV9+rCqWXjtXRFODuG6YZeLQkRXDezeTd23UebNdDUrgmjoUE6y8c455hjr93S7RW67LfZ5xwu48/vNx+pwpL6oyqRJ/Cz1CdrvZ3qHq64yru1rVtFMhKkgpk4VGTZM5JZbRJ54Iv53oEuX1Km1duygsNXHbPQMfb4aq0pWqwKA98MgAL+Duv0xlcfuBzCk8ncPgHcA/AFgHoCOVvpVAkARl5IS49wxTqfI+eez9kB4Wch0T/zJtkRjThSZrBt7fT7aCnbujH2GU6ZUzQjdrFnNeNTMn88Jf+hQkddfZwprs93SSSdZ7/foo+O/H69XZN066/0FgyzBacSttya2Z9ntLIpTA9S6AKippgSAIi6//mpeVKYuTvipbi4XhaHTSQOqEcGgeQ1js+bzsUpYbbBihfln7HJZ34UcdFD89+TxiOTlJe6ntDSkjtI0qh2//TbynE6drD3HPn2Sfx4WSEYAqHoAirpL06bm0avJpBQ+6qjUjEfHZssMI19pKSN9y8qA666jMbxdO+CMM4DTTwfuu48BYF99Fb8fm40GdbebRufJk1njQWfvXraaoG1b42A3gEFj48ZZ66djR/PXNA3o1g1o3TpxP1dcAUycSMO5CPDrr0z58euvoXMCgcT9OJ1MzZ1urEqKdDS1A1Ak5OSTrSdb02vUulzcIQwdyjTI5eWsd5yKVXe7dqx69vXXmRNsFq/pxnaz1202ulX+619Uzbz6Kp+ZCEuEHnNMaJdx9NEiq1al/jO+6y7z8bVokfj6eInunE66uK5cmbifLVuMP1ObjbYGncmT49sbNE2kQQPWea4BoFRAinrDzp0i/fqFvDk8HvMJzekUefttkVmzGDwUTnFxbHBYVbyGXC5OEn371l7N45pqep7+cK+lQECkWzfGU7RoEalqs9loGygoSO1nvHu3+WfavHni6z//3LzOQM+eIYGWiB9/5MRtJihbtRI57TSqhIYP5zPLyuIza9aMtR2ys+mNZEXgVJFkBEDGxQEoFEnRqBEzR65fT3XGTz8Bd9xh7APerJl51Sy3G/jhB/qPT55Mf/czzqCa6euvgT//BDp0YH7/E080VzHpKqkFC6qWpC5TcDioJlq0iO9dZ98+5v6/8ko+4/DnEAwCe/YAgwcD11wDDBkCfPst8OmnTMh28cXJ5/QBmCCuTx/gxx8j7+d2A8OGJb7e7+c0bUTXrrH5jczo3Nk8O2kwyMysH33E7+P48cA99zBxXuvWzMKaCWrBaKxKinQ0tQNQJM3s2SH/7+jVv55XprqMHGltBW231w01kNHKf+RIkbVrzVUnDRvGN7RHeyg5HNxJfPBB1Z75qlVc7eufbVYWq5Tt2ZP42ooKrr6jx+j3MxVHMowaZc1ryu9P/U7IIlBGYEW95fjjaegMz/apF1sfPjw195g0ibuDRPj9wDHHZObKLx52OwvTezzmO538/PiG9ugaBeXl3DFcfHHVcvx37sydyKRJNF5PncrdxZNPslrZcccxO6wYrPRtNu5Cmjdn1bGsLK76b7uNmVWT4amnuLJPVGnMbuduMtOxKinS0dQOQFEltm8XOftsrvrtdlbA0jNiJqKkhJHE0QFn0Xz/feIcPh6PyMUXp39FX5Xm9TLKuE+f1LrUZmcnv+o2Yt8+Bm+Fxwj4/XTRNKO0lMGFb75pHC2eDD/9FN/Q6/NZ/86lGKgdgKJe06RJKH9NYSFXit26xb9GBHjoIer8O3ViHw8/bLyiBKgLT1RgpHFjIDfXmltgplFUxDq9RxzB1WwqSUUthilTmKsofDdRUMAV+qZN5vc9+WTmCmrVqnr3792b3xObwRRqt9O2kOg7lwEoAaDYf3E4rFWhAoAJEzjh5+dz8svPp0B46inj8486KnEFrc2bWdCkfXuqU+oahYXAK6+Y++FXBZeLBtHq8umnxoZ+t5vG4ppG01gl7vjjQwJS0/g5d+/OqmZ1ACUAFAqAk3/0hFJYCPzrX8bnZ2UlrjcbDLK+8T33AGPGcGJo0cJ41ZipFBUlF1QXz6PG5WK21lSUZmzTxnhnEgzyGdcGLVoAs2dT0P/yC71/5s3j71aCyjKAOvRNVChqCBHWnzVi61bz60aPpmE53oRWVETDZ+PGTK08Z451t8O6hsvFes533hmr5unenc84VdGv11wTu7uz2VgKs7YjbPXU4yeeyJ+JEKFa8t//Bl57LTIdd3k5XZkXLUpO8FYVq8aCdDRlBFbUGmZZMbt3j39dMCjy3/8mztzp87HsogizVOqFTdJt7E1Vc7lYOEanuJhRw48+WnNVuN59l+mss7JotD70ULquZjIlJSInnhhKbR0I8D38+isD1ho3DgWPtW5dpTTVUJHACkUlhYUizz/PqmHXXSeybJnxeZ98Euvf7fXSa8TqfeJV59K0UBlHEXqxfPop/+HTPXlXpeXkUGg6nYxyfeCBUHbMTz9lSgy9lOTVV/P51ASlpSILF9ZMCoqa4LHHjL3HOnUyji9o2DDpeAIlABQKEeaA6d495K5nt/Of7N13jc//6ivmtmnaVOTYY/l3MkyfHj9F8+mnx17ToYPxuW639RxH6Wher/Hq1Chvj8slcuqpkecVF7O0ZM+erOs7aZJ5euX9CbMCOw6HcdBgVpbIG28kdYtkBICyASj2XyZOZPCQrmOtqKBh98orjdM0/POfwNy51FXPmcO/k2HfPnP9vt0OnHVW7PFzzjG+pmVLa/rkdFFUxBKN4Xz/PfDoo7HnlpbSQKqnlAgGgf79aRz/5Rfg559ZXvLMMznt7c+Y6fVFjEuVlpUB27fX2HCUAFDsv7z9NieqaCoqaiZK85hjzF0m27UDzjsv9vhdd/E1v59/u92MG3j99eSNgCeeSINkbUUeL1gQOcbnnjN//3Y7sHo1f//sMz7/8M+msJBpqX/6qebGmwlcdhng9cYeb9bMOF5E05JfiCSBEgCK/RezcP2KipoJzurQgf/g+mQO0DOldWuuco1iEho2BJYs4W7l0kuZguG334Bjj2WKg2SCsPLz6ZL4zDOAzxf5WiqCr4zw+TjOTz9l4jMzSkpCRdznzOFuKZrSUu7A9mduuAE4/PDQ98/n4/f0o48YXBb+ufn93DXW5E7Qqq4oHU3ZABTV4r33YsP1bbbEnj3VIRikl0/fvtRvP/wwbRFVYfly4/Gb6eW7dQtd++ST1B/bbOyjZ8/0eh3Z7TSCB4MiTz1lbAgNBFgGMp2Ulorccw/rQ7hcIkceKXLeeawd8cor1lNHx6Oigs4F997Lymq7d/N4SQn/PuYYkX/+k9+jioqku4cyAisUwsnmlltCedmzskTatmUhk7rCd9+xtnGiCdZmE7nmGl4zfTonWK/X2JDsdKZHCLhczNa6dKmxx0ujRvR4KS0VmTaN3kMPPmitVGOquOAC8xxPfr9Ir14iRUW1N54qkIwA0Hh+ZpKbmysLFixI9zAUdZ0NG2igbN4c+Mc/6lYkLgBMnw5cdJFx6gOdBg1YljA7GzjggPjnppNAgMFONhttAE4nA+latQLef59ZP485Bli1isZ7m43t+eepXqtJ1q9nDp942Up9PkaARxvAMwhN0xaKSK6Vc+vYf4JCUQVat6a3zfHH173JHwBefDH+hO5wAPfey/QIM2emPnlbKtm3jxNsYSHX1ZoGPPEEJ/yePZneecWKkOdWMEiBcfnltJPUJL/9ljhKu7AQmDYt/jnvvgv06sXv3QUX8L1lKNX6b9A07T+apv2madovmqZN1zStocl5azVN+1XTtMWapqklvUJhhd27adBdvjz+eV4vV84ADdwZvKuPoaQEeOmlkOfSG28Yu0MCwC23mKfsSAVdupjfO5yGhtMcefxxGvMXLQI2bqQnWm4usGZN6saZQqq7HPoCQA8R6QngdwB3xTn3BBE5zOrWRKGo1yxZwiyit94aWZLRCIcDOOUU/n7yyVUruJJOli0L/R5vBe5wMANnTdG+PTBgQPzMrX4/8xAZUVQEjB0buVsLBrmbefDBlA41VVRLAIjI5yKil/75EUBO9YekUChwwQWsrxtP9eP10q7xxRehiXPhwtSmbzYjlbEG4avuK68079tmq/lEetOmcQy6r77TSb1/djYFw223sS60EX/8YaxirKjIWPfWVBaFHw7ATDkmAD7XNE0APCcik1N4X4Vi/2LDBvNVfyDA4LEOHYCcHObW13X+774LXHihuQrIZktdhslUqpn0bKoiwIcfUgAY9R8Msg7DffcxrXTz5sCNN3LXkyo8HtaAmDAh9Ky++QbYuZMOBPFSTbdqZV4jon371I0xhSQUAJqmfQmgpcFLY0Tkw8pzxgAoBzDVpJtjRWSDpmnNAXyhadpvIvKtyf1GABgBAG3btrXwFhSK/Qy73XyiLi4G7rgj1tBbXg6MGGE+AdntXD3XtHeQzcbJ2+lMXDBHRxcA8+YxSMzovdvtwLBhDKLKzw/tGubMAcaNo6oslWha6BmfeKK1a5o2BU47jcVgwtVwXi+F1Pz5tAdkUo1oq/6iZg3AZQB+AOCzeP44ALdaOVfFASjqLW3amPvSz5gRe/4vvzCQyswH3+tl8jUj/3srQWZWmt/PJHr9+jGQqmFDa336/XwP48cbJ0RL1LzeUDBVuikoEBk2jO/D5+N7c7kYhBcIMEvqihU1OgTUVjI4TdNOAXA7gCEiYri00DTNr2lalv47gAEAllbnvgrFfs8JJxgf17RQTp1wGjbkLsAIm43eRDffDMyYwb7btOE1Ph8Nm34/UxHcdlvVx1xQQF33rFnAjh30YrKictJTdrRqVbWUFS4XvW4yAZ8PmDqVhYSmTaOIKi0F9u6lC+xffzERXm3YaSxQXS+gpwFkgWqdxZqmPQsAmqYdoGnaJ5XntAAwV9O0JQDmAZghIjOreV+FYv/mpJMicwrpOJ3MvRNNmzbAYYfFVifzeICXX6ZrIsDEYrNncyLauZMG5AkTOGnPmgXEU7tqGtUcqcwr5HJRdfX118yDZJQjKBFlZTTSPvMMVS0XXQT88EPqxlgVsrOZ3yfaI0uEwiBTjMJWtwrpaEoFpKi3FBWxVkB42ga3m7lpgkHjazZuFDnkEKodsrOZAuPOO83P15k7l/mR7Pb4qpZ//IPnf/wx+69KbiFNY0oOp5PjPO44kUWL4qumErWsLJEePUJ9aBp/nzjR/D3Pny9y661s8+dX7TNKxOmnG483O1vk/fdr5p6SnAoo7ZN8vKYEgKJes3UrE6g1bMgiNTfdxEpi8QgGOaHOmMHrE7FkSfwiNuEtvDBJcbHIrFnx7Q5G7aCDmKTviSdE5szheK+4IrHwidccDuOcR16vyN69se/57rspIGw2Np+Px2bPFrnqKuYgmjs3uc/KiJdeik3mB1AwW/lsqogSAAqFwhpmxmajlfvAgSz3GM5NN3FCs2ok1oVIUVEoS2rfvlWf/OO17Gxm1Lz6au6mjjqKWVKNkr3Z7aHjmsax3n579Z5tcbHIEUdE7m58Pia4q0GSEQAqGZxCUZcpKKDR8eefgR49qP82q4MQzYcfAmeckdz9fD7gzjuZewhg9OvgwSzkYrPR4FleTtuD7qrpdtOlcuxYJnQ77zzq/AHmuj/kEOCtt4yrtIXjdCY+J3qsTiefkW4gdzppgLVinPZ6aVzW6xhUheJi2mDeeQdo1IhRxP36Vb0/CySTDE4JAIWirrJhAz139u7lJOfzsf34I9CpU/xrv/uOE1G83Dcul7Evv8tFQTNnDktX3n47E58tXcpsmoceCuTlMVBL0+gRk5PDibdFi1hDr9/PcZh5MVnB4Yi8XtOArCxOwFbjEaJxu4GHH6b3VB0iGQGQykhghUJRm9x4IydX3aWwsJAT3tVXswZvPG66Kf7kf9BB9PiZMyf2tbIyYMoU3nfVKqafeOCByIlST04HcCIGgHPPNfbyKSqqfpbWAw6gZ5OO30/h9PvvVe/Tbo+trLafUQdz4yoUCgD06Y/2Jw8GqV5J5Ge+ZIn5a243Szx26WI8MYtE9l9YyEydp50GrFtn3u9HHxkfDwarl6LCZmPmzXCKi5mbxwink+9Rj39wu41dW0VYkjH62CefAEOHAgMHsnZzdXYuaUYJAIWirhLt869jsyVON9C0qfFxu53F3tu1A0aPjp8ZM5oZM6iS2r3b+PV4+nuz2AIrtQ30mgHhlJcbCxWXi7EQ69cz589TT/H3l1+mzj8ri83nY2rqZs0ir7/pJu5kpk9n7YWRI4FTT2Ut5tmzjYP0Mhmr1uJ0NOUFpFDEYdSo2NQJTqfIuecmvnb8+Fjfe59PZOzYyPOmTKGrZ3Y2PWMS+ev7fHTxNGLAAONrnE66uhr1NXasdTdVK55Mgwebp43YtUvkrbdE3n7b2H101SpjjyeHg++hQQN6EvXvb3x9LQHlBqpQ1AP+v737D7K6rvc4/nyzwHoWQjTAu4qBkih7u0a2kmZSjvcaoiNjY+PvuXipe7MrYuNkiTXdmhrLX+ikZQQ0kzmhA15sGrxe05tMU/6gIgUENUJBU5dMZXdxt+O+++N9zt2zy/fs7uHsOd+zntdjZgfO93zPft/sLJ/39/v59X7rLffW1migM5n4s6XFfc+ewT/b0xPFz5uaej//hS+4Z7P7n9vZ6b5hQ6wZeOSRwZPARRclX3P79mjoC+f8jxkT33PTppiqOW5cxNPcHFM4Dzts+GoYZzLl/bzvvLN4veDCr8bGKCSfklISgAaBRUaq97wndtD81a96Z+CcdtrQBlTNYuD22mtjxs7hh8dW00kymdgKOe+734Urr+wt29j/3OOPT/4+M2fC9u2wfHkMHM+eHV0o+S2W//jHKMuYzUJLS0wPfe21aFaHQ0tLeZ+fOLF4t1uhri5Yty7GRmp8EFnTQEWkdF1dUff22Wd7+/bNYg3Cc88VH2MYqueei+mk+/bt/96B1DXIZGIQupw5+J2dkSjffHPwcxsbY2yh/xhCFagovIhUVmNjrCU4//wYWB01Kp4Sfv3rvo1/NhszZT75yRgsXbduaHf0XV3FB4CnTYt99ceNi5lKLS3x9zFj4qloypTY4G7WrDje2hp79Je7AKupKQZ+J02Kzd4mTIgYk564mpvjvPvvjyI206bFhnw1VhtYTwAiUp78NM7+3SPu0eg/+mhvd9G4cVHYZfkgRQF7euJu+9VX+x7PZGJF8Re/GN1f+/bBRz4SXUpPPBELzhYsKG32Uqmy2Uh03d2xAO6UU+Lf193dW3hn3bqodXzddb1FeBoaoptt06aKVgjTQjARqZ5Ro5Lvgh9+GDZs6DtW0NERTwRXXVW8T/7tt2NriFmzoq6AWXQzjR8fC8zmzo3tr/fujfd6emDlyliLUA2jR0cMeVu2wK23xljMscfGgrgZM+Dcc/tWYHvnnd4C8StWVCfWQSgBiEhlPPhg8spf91ipnJQA2tvh5JOjJnJHR+9TxRlnwKJFMH9+3D3/5S99P7dwYQwqz5w5eFzZbHTNbNgQieTSSweu9TuY5mb4znf6Htu8OTkpZrNx3RqhBCAilTFpUnSH9N9yYvRoOPTQ5M/ccUfMBsoP/uYXeG3cCD/+Mdx88/5FVvLnrVwZBePXroWnnooniNNPj7vztWuj+2nRIrjrrlgl3N4eXUVf/3okq49+dPj+7YcdVnwPohqqda4EICKVcfHF8I1v7H88m40N7F5+Ofr5C61enTzzp709Gs6GhuT3//a32AvouOPi6aC9PRr8t9+OO/H8TKUvfanvVhb5ZHLhhbBz5/AVbJ88OZ5W1q/vmwDzu6nWCM0CEpHSPfYYfOxj0S9/zDGxOVz/CSVTp8K99/bOmMlk4rhZ7CA6YwYsW9b3M8XWInR3x1dS4w/R2O/eHTuk5rudOjqioS/cgiKbTd4nac+emHo6nO66C845p3ffoYkT4fbboyZwjdAsIBEpzcaN8PGP9x3gbGqKu/2kgdiurqg9fN55+3cHZTKRTPKLx1avhs98JnmRWTFNTXDCCbF3f2FMpchkYOvWyszO+etfI8FMnz689ZSL0DoAEamcr3xl/4a2szP60pP6vRsbo8snqfHr7o5ZQXnnnw+XXRZ98/mN2Yqtvm1oiM3n8kXtD3RLabN4GqnU1MxDDomnpCo0/qVSAhCR0mzalHz8nXdiV8wkXV3JC8B6evomE7PYamLbNrjzzth18+qrI4n0N358TL1ctCgWo+UXpQ1k1Kjojmlq6k0ykyfDmjUDf+5dSoPAIlKa979//wVaEA18sa0P5s1L7ntvaoquof6mTYsviOmdP/0ptLX1Fo856CD43vf6Nvg33RSlKXfujCeLxsZo4M3g9dcj2Zx4YnQzvfoq/OY3MYXz7LMHTxzvUmUlADP7L+CzQFvu0FJ3X59w3jzgNqABWOHu3y7nuiKSoq99LWoJ9x8DuPzy3oHe/pqb4frrYenSaJx7euIzn/50jCcM5L3vjQI2P/hBzKp53/tgyZLY4qHQxIlx3kMP9W6Od+aZ0VX0wgtxvSlT4twjjohxgzpX1iBwLgG0u/tNA5zTADwL/AuwG3gSuNDdtw72/TUILFKj1qyJ4iivvBJ341deGYPAgxVwefrpmB2zb1/c+c+dO3xTLwWova0g5gDPu/sOADNbDSwABk0AIlKjzjsvyiV2dMRd/1Aqd0Fs8XzDDZWNTYZsOAaBrzCzp8xslZkdkvD+EcCugte7c8cSmdm/m9lGM9vY1tZW7DQRSZtZDMQOtfGXmjNoAjCzX5jZ5oSvBcD3gRnAbODPwM3lBuTuy9291d1bJ6ewl7aISL0YtAvI3Ye0bM3Mfgj8POGtl4AjC15PzR0TEZEUldUFZGbNBS/PBTYnnPYkcIyZHWVmY4ELgJ+Vc10RESlfuYPAN5jZbMCBncB/AJjZ4cR0z/nunjWzK4AHiWmgq9x9S5nXFRGRMpWVANz90iLHXwbmF7xeD+y3PkBERNKjrSBEROqUEoCISDVs2wa//CW88Ubakfw/JQARkUpqa4M5c+DDH44tNJqb4VvfSjsqQAlARKSyPvWp2EG1sxPefDOqkF1/fdQlTpkSgIhIpezaFQV0CquSQWyhccst6cRUQAlARKRSXn+9eCGYtjbYsSN2UW1thYULoypZFakegIhIpcyalVypbOzYGBeYPTt2Rs1mo5tozRp44AE49dSqhKcnABGRShk7NgrBNzX1bnvd2AiTJsGf/gR790bjD1Ewp6MDPv/5qoWnJwARkUq65JKoorZsGbz4YlRHW7wYpk5NPn/r1iihmVQGc5gpAYiIVNpJJ8E99/Q9NmFCdP/019hYtQLy6gISEUnDkiXRNVQok4ki90njBhWgBCAikoZrroGLL46SmgcfHH+edRbceGPVQlAXkIhIGhoaYPly+OY3Yft2OProKFZfRUoAIiJpmjIlvlKgLiARkTqlBCAiUqeUAERE6pQSgIhInVICEBGpU+buacdQlJm1AS+kHUeBScCetIM4QIo9HYo9HSM19uGIe5q7Tx7KiTWdAGqNmW1099a04zgQij0dij0dIzX2asetLiARkTqlBCAiUqeUAEqzPO0AyqDY06HY0zFSY69q3BoDEBGpU3oCEBGpU0oAB8DMFpvZNjPbYmY3pB1PqczsajNzM5uUdixDZWY35n7mT5nZf5vZxLRjGoiZzTOzinogXwAAAxhJREFU7Wb2vJl9Oe14hsrMjjSz/zOzrbnf7yVpx1QqM2sws9+b2c/TjqUUZjbRzNbkfs+fMbOTK31NJYASmdlpwALgg+7+j8BNKYdUEjM7EjgDeDHtWEr0EPABdz8eeBa4NuV4ijKzBuAO4EygBbjQzFrSjWrIssDV7t4CnAT85wiKPW8J8EzaQRyA24D/cffjgA9ShX+DEkDpLge+7e5dAO7+WsrxlGoZcA0wogZ/3P1/3T1XPZvHgCIFVWvCHOB5d9/h7t3AauKmoea5+5/d/Xe5v+8lGqHqblJfBjObCpwFrEg7llKY2cHAXGAlgLt3u/sblb6uEkDpZgKnmtnjZvaomZ2YdkBDZWYLgJfc/Q9px1KmfwMeSDuIARwB7Cp4vZsR1Ijmmdl04EPA4+lGUpJbiRucnrQDKdFRQBvwo1z31QozG1fpi6ogTAIz+wXwDwlvXUf8zA4lHo9PBO41s6O9RqZTDRL7UqL7pyYNFLu735875zqim+LuasZWb8xsPLAWuMrd30o7nqEws7OB19z9t2b2ibTjKdFo4ARgsbs/bma3AV8Gvlrpi0o/7v7Pxd4zs8uB+3IN/hNm1kPs39FWrfgGUix2M/sn4i7jD2YG0YXyOzOb4+6vVDHEogb6uQOY2ULgbOD0Wkm4RbwEHFnwemru2IhgZmOIxv9ud78v7XhKcApwjpnNBw4CJpjZT9z9kpTjGordwG53zz9trSESQEWpC6h064DTAMxsJjCWEbDplLs/7e5T3H26u08nfuFOqJXGfzBmNo94tD/H3TvTjmcQTwLHmNlRZjYWuAD4WcoxDYnF3cFK4Bl3vyXteErh7te6+9Tc7/cFwCMjpPEn9/9wl5kdmzt0OrC10tfVE0DpVgGrzGwz0A38a43fjb5b3A40Ag/lnmAec/fPpRtSMnfPmtkVwINAA7DK3bekHNZQnQJcCjxtZptyx5a6+/oUY6oXi4G7czcNO4DLKn1BrQQWEalT6gISEalTSgAiInVKCUBEpE4pAYiI1CklABGROqUEICJSp5QARETqlBKAiEid+juWObICM4YEiwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f575a2d24e0>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "y_fitted = model.predict(X)\n",
    "colors=np.array([\"red\", \"blue\"])\n",
    "plt.scatter(X[:,0], X[:,1], color=colors[y_fitted])\n",
    "plt.scatter([], [], color=\"red\", label=\"0\")\n",
    "plt.scatter([], [], color=\"blue\", label=\"1\")\n",
    "from sklearn.metrics import accuracy_score\n",
    "acc=accuracy_score(y,y_fitted)\n",
    "plt.legend()\n",
    "print(\"Accuracy score is\", acc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Even thought the score is quite good, we can see from the plot that the algorithm didn't have good models for the data. We can plot the models the algorithm used:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:19.194395Z",
     "start_time": "2020-06-24T19:29:19.116088Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZsAAAEJCAYAAABCNoqwAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAHh9JREFUeJzt3X2QJHWd5/H3t7v6aXp6umd6ehiYBofnWR8ApVFkPBCJ9YF1WY11T7xzFZBF99BDxYvgTneFWI2Ni/UUV7x1WVa4gA3ZDTZCYw3h9BRERcRGHpZYYWaZGaBnmKGn57G7p5+qv/fHL4tpmn6onq5fZlbV5xWRZFdmZdW3k5r69C/zl780d0dERCSmhqwLEBGR2qewERGR6BQ2IiISncJGRESiU9iIiEh0ChsREYlOYSMiItGlFjZm9gEz+4aZ/czMDpmZm9ldi2xzgZn9wMz2mdkRM3vSzD5tZo1p1S0iIstXSPG9vgCcDQwDA8CmhZ5sZn8A/DMwBvwjsA/4feBrwGbgj2IWKyIilWNpjSBgZhcTQubfgYuA+4F/cPcPz/HcVcnzOoHN7t6fLG8FfgK8FfiQu9+dSvEiIrIsqR1Gc/f73X2rl5duHwB6gLtLQZO8xhihhQTwpxHKFBGRCPLaQeAdyfy+OdY9CIwCF5hZS3oliYjIscpr2JyZzLfMXuHuU8B2wvmmU9IsSkREjk2aHQSWojOZH5xnfWl513wvYGbXANcAtLe3n7tp04L9EUREZIZHH310r7v3VOr18ho2y+butwK3AvT19Xl/f/8iW4iISImZPVfJ18vrYbRSy6VznvWl5QdSqEVERJYpr2HzTDI/Y/YKMysAJwNTwLY0ixIRkWOT17D5STJ/9xzrLgRWAA+5+3h6JYmIyLHK6zmbe4D/CVxuZt+YdVHnl5Ln/E1WxYmIlGNycpKBgQHGxsayLmVera2t9Pb20tTUFPV9UgsbM3sf8L7k4fpk/lYzuyP5ea+7fw7A3Q+Z2Z8QQucBM7ubMFzNZYRu0fcQhrAREcmtgYEBOjo62LhxI2aWdTmv4u4MDQ0xMDDAySefHPW90mzZnAN8dNayUzh6rcxzwOdKK9z9u2Z2EfB54A+BVsIQNp8F/rrMkQhERDIzNjaW26ABMDO6u7sZHByM/l6phY273wjcuMRtfgFcGqMeEZE05DVoStKqL68dBEREpIbktYOAiEjtefhhOFDBywO7uuD88xd92n333cd1111HsVjk6quv5oYbbqhcDWVS2IiIpOXAAeip2AgwUMa5lmKxyLXXXsuPfvQjent7Oe+887jssst47WtfW7k6yqDDaCIiNeyRRx7htNNO45RTTqG5uZnLL7+c733ve6nXobAREalhO3fu5MQTT3z5cW9vLzt37ky9DoWNiIhEp7AREalhGzZs4IUXXnj58cDAABs2bEi9DoWNiEgNO++889i6dSvbt29nYmKCu+++m8suuyz1OtQbTUQkLV1dZfUgW9LrLaJQKHDLLbfwrne9i2KxyFVXXcXrXve6ytVQJoWNiEhayrgmJoZLL72USy/NdjAWHUYTEZHoFDYiIhKdwkZERKJT2IiISHQKGxERiU5hIyIi0anrs4hISrK4w8BVV13F97//fdatW8dTTz1VuTdfIoWNiEhKMrjDAFdccQWf/OQn+chHPlK5Nz4GOowmIlLDLrzwQtasWZN1GQobERGJT2EjIiLRKWxERCQ6hY2IiESn3mgiIinJ4A4DfOhDH+KBBx5g79699Pb2ctNNN/Gxj32sckWUSWEjUmGTkzA+Hualn0dG4MiRME1OgvvRCcAMGhqgUIC2tjC1t0NLCzQ1ham5OTyW6pXFHQa+853vpP+mc1DYiCxDsRiCZHQUhoZg794QKDOZhRApTS0tYVlpXSlwSq83PByux5iagunpV75WczN0d4dp5coQSE1NcX9HkUpQ2IgsgXsIg8FB2L0bDh0Ky8xCiLS1hRA4VqUwms/UVAii3btDEJmFwFm3Do47DlatCi0kkbxR2IiUYWQkBMzzz4efGxvDl/zatUdbKWkoFEKYzQy0iQl44QXYti20fHp7Yf36EDxp1iZzc3csx/8jfGbTOiKFjcg8pqZCC2LHDjh8OLQYOjpCKyJPmpvDBKHm558PwdPaCiedBBs2hJ8lfa2trQwNDdHd3Z3LwHF3hoaGaE3hA5LrsDGzHcBr5lm9x93Xp1iO1InJSdi5E7ZuDedQ8hgw8ykUoDQyyeQkPPts+D02boTXvCYc5pP09Pb2MjAwwGAlu6BVWGtrK729vdHfJ9dhkzgI3DzH8uG0C5HaVizCiy/C00+Hn7u6wpd3tWpqCh0JpqfDYbbt2+G000LolFpCEldTUxMnn3xy1mXkQjX8Uzrg7jdmXYTUtsOH4fHHw8n/1atrq4dXQ0No7RSLIXCeew7OPruyow+LLEb9VqTu7doFv/hFaAGsW1dbQTNTY2No6axYAb/+dTi8NrtrtUgs1dCyaTGzDwMnASPAk8CD7l7MtiypdsUiPPNM+Gu/u7t2Q2a2lpbQqtm6FQ4ehNe/Xh0IJL5qCJv1wJ2zlm03syvd/adZFCTVb3ISHn00XLNy3HH110W4oSH83gcOwC9/CW9+c+jKLRJL3g+j3Q5cQgicduANwN8CG4F7zezs+TY0s2vMrN/M+vPcE0TS5w7/9m/hgsyenvoLmpm6usLv/9hjIYBFYrG0LuipJDP7CnA98F13f/9iz+/r6/P+/v74hUlV2LYt9Dg77risK8mPoaEQvOecU9/hK0eZ2aPu3lep18t7y2Y+30rmF2ZahVSdwUH47W/Dlf9yVHd36Pa9bVvWlUitqtawKR0X01FmKdv0NDz1VOja3NiYdTX509MDW7bA2FjWlUgtqtawKQ3Urb/DpGwHDoQvUg3TP7eGhjDt2ZN1JVKLchs2ZvY7ZvaqlouZbQRuSR7elWZNUt1eeEFBs5hVq0JX8Co8lSs5l+euzx8ErjezB4HngMPAqcDvAa3AD4CvZFeeVJOJiXBOors760ryrbk5tAAPHizvLpAi5cpz2NwPnAm8EdhMOD9zAPg54bqbO70au9JJJsbHw1z3elmcmc7bSOXlNmySCzZ10aZURLGoLr3lMgu3KhCpJP2dJ3WhoUHnIZZCvfWk0hQ2UhcKBYVNuaanq/vWCpJPChupC21t4SZoOhexsOnp0Ars7My6Eqk1ChupC2bhbpWHD2ddSb4dOhRuI62bq0mlKWykbqxbF0JH93CZ38QEpHCHYKlDChupG83NcOKJsG9f1pXkU+naGh1CkxgUNlJXzjgjXCV/8GDWleTL2Fho8WnUZ4lFYSN1pVAIX6ju6ixQMjUVwvfcc0NHCpEYFDZSd9rawhfroUMKnKkp2LsX3vCGMBq2SCwKG6lLXV3wlrfAkSNhLLB6NDoabpp29tnhXJZITAobqVtr1sDb3hauv3nppfrqpbZvXxjCZ/Nm9T6TdChspK61tkJfH5x6agicI0eyriiuiQnYvTvcqfSCC9TzTNKjQSmk7jU0hF5q3d3hltEvvRS+hGvp3jdTU7B/PzQ1hQ4SJ5ygXmeSLoWNSKK7O/y1v2cPPP106KG1alVo/SyLezhmFR5AaYw2S/5jFhIvwrf/1FQ4J2UWAvXEE0PgiKRNYSMyQ0MDHH98GG1gcBC2bAktnUIhBM+cA1ROTYVubZOTYRofh/ExGBsPP09MLD4KqFm46rSlBVpboKU1pFyhENKhra3soZinp8OwPOPjYdNNm0JLRiEjWVLYiMyhsRHWr4fjjgtdpHfvhoEBmBgr0lwcY2XTGIWRQ6H5Mzp6dMNSK6WxEQqNIUDa2hZvtZRaP1OTcGAMivtCariHyQxWrICuTuhImlutrS8H0PQ0DA+HzGtsDOFy/PGh151uFyB5oLARWYAZdLaO07lqH6ev2sWhAwd4caiJnXtbmWooQFMLze3dtLU4zU3LuIeBWWjFFAow17ki99BqemkQdr3IVNEYm2xkrK0LX7MW61jJ+pOa2bAhBIxuESB5o4+kyFzGx0P/4IGBcNWjGQ0rVtD1mk66NhpnTsPIEWPkyBRDB52hgwUODjcChplTaISmglNodAoFp/EY+n1OT8NU0ZgqwuRUA5NTBdzDJf7NrU73uinWNg+xkudotyKFiW6Y7IXiGijUUO8GqQkKG5ES93A2fdu2cKKmdOiqp+dVh8EaGqCjfZqO9mnWr50CxpmcgpEjjQyPNjB8pIEj4w0cGWvg0HAjU0n/ADNwN472Ejiq9Bal0zsNBitap2lvnaatZZoVreH92tumaWkubd8MrAkbjY7CE0+ExT09cMopoZmjbmeSAwobEfdwKf2WLeEcTFvbnAGzmKYCdHUU6eoovmpdsQiTU8bklFGctpffthQsVuqUZk5TIUxLOhRmBu3tYXIPPQQefjj04T7zzHAFq0JHMqSwkfq2fz8880w4ZNbREbqhRdDYCI2NTmtLCvemNgu/S0cHjIzAr34VwmbTptDSEcmAwkbqU7EIW7eGQ2bt7aHbWS0qtXaGh+Ghh+C008JwCeqiJilT2Ej9GRmBJ58Mh8x6esIJmFq3cmU4//Tss0dH31yxIuuqpI7Uwb8ykRn27IGf/zxckFIvQVPS0BAOEx45EvbBnj1ZVyR1pI7+pUnd27kT+vvDUACrVmVdTXY6O8P5nEcfhV27sq5G6oTCRurD/v3h0NnateGq/nrX3BwGg3v88fq9oY+kSmEjtW9k5GiLRpfWH1UohFZOf/8rh9wRiUBhI7VtehoeeyyMQrns4ZtrUGl8tcceq6+7x0nqch82ZtZrZt82s11mNm5mO8zsZjPTHdNlcfv2hQscOzqyriS/Vq0Ko43u3591JVLDch02ZnYq8ChwJfAI8DVgG3Ad8Esz686wPKkGO3aoi2852trCvhKJJNdhA/xvYB3wX939fe5+g7u/gxA6ZwJfzrQ6ybfR0XBTmpUrs64k/1aurI/7Yktmchs2SavmncAO4JuzVn8RGAH+2MzaUy5NqsX+/RoPrFylwdnUM00iyW3YABcn8x+6+yvOXLr7YeAXwArg/LQLkyoxPq7eZ0vR0BD2mUgEeQ6bM5P5lnnWb03mZ6RQi1SjYlEtm6VoaFCPNIkmz2HTmcwPzrO+tHzOYWzN7Boz6zez/sHBwYoXJ1WgpSUEjpSnWNQFrxJNnsNmWdz9Vnfvc/e+np6erMuRLLS3K2yWYnpaPfckmjyHTanl0jnP+tJyndGUua1eHS5anJzMupL8m5gI+0r3u5FI8hw2zyTz+c7JnJ7M5zunI/WuoQFOPjlcsCgLO3Qo3Ea6nkbBllTl+ZN1fzJ/p5m9ok4z6wA2A6PAw2kXJlXkuOPC4aGpqawrya+pqXAr6Vq9gZzkQm7Dxt2fBX4IbASunbX6JqAduNPdR1IuTapJWxu89rWwd2/4QpVXmp4O++b1r9fYcRJV3i9C+C/AQ8Bfm9klwG+BtxCuwdkCfD7D2qRanHRSGPl5+3b99T7b4GC4TXRvb9aVSI3LbcsGXm7d9AF3EELmeuBU4OvA+e4+lF11UlXOPDPcpXJIH5mXDQ2F8D1Dl6pJfHlv2eDuLxAG4hQ5dg0NcNZZ4WZhe/aEm6g1NmZdVTaKxXDobN26sE/UKUBSoE+Z1I/mZujrC62cvXthbCzritI3NhZ+902b4E1vCvf5EUlB7ls2IhXV0BDOUaxeHW4YduRIuLak1oe1cQ+DbDY0wFvfGn5/kRQpbKQ+rVkDmzfD00/Diy+GXmu1eoO1Q4dCqPb2hlZdS0vWFUkdUthI/WpthXPOCRd+btkS7ufS2hpCp9pbOu4hZMbHoacHzj033JFTJCMKG5HOTjjvvHCYaft22L07HG7q6Ki+VsD4eLgN9vQ0HH98CNLO+UZ8EkmPwkakpKsL3vjGcE3O4CA8/zwcPJj/4Bkbg+HhEDDt7eFQWU9P+FkkJxQ2IrO1t4dp48bwJb5379HgMQu92trashuOf2IinIOZmAiHyzo6Qu+ytWsVMJJbChuRhaxcGaaNG0OL5/DhcDHkvn3hsFvp3E5ra+hGXChU7u6gU1NhmpgI0/R0eL/2dli/PnRy6OhQwEhVUNiIlKvU4lm/PjyenITR0RBCQ0MhiEZGQjDMpRRMpXlprLb5xmxrbg4htnYtdHcffX/d6lqqkD61IseqqSmcfO/shBNOOLp8ejoE0cypNLLy9PTRcDEL54PMQoA0Nb1y0pX9UkMUNiKV1tAQOhPktUOBSAb0p5OIiESnsBERkegUNiIiEp3CRkREolPYiIhIdAobERGJTmEjIiLRKWxERCQ6hY2IiESnsBERkegUNiIiEp3CRkREolPYiIhIdAobERGJTmEjIiLRKWxERCQ6hY2IiESnsBERkehyGTZmttHMfIHp7qxrFBGR8hWyLmARTwDfnWP5U2kXIiIixy7vYfO4u9+YdREiIrI8uTyMJiIitSXvLZsTzOzjQDcwBPzS3Z/MuCYREVmivIfN7ybTy8zsAeCj7v58JhWJiMiS5fUw2ijwF8C5wOpkugi4H3g78GMza1/oBczsGjPrN7P+wcHByOWKiMhCooWNme1YpPvy7Omu0rbu/pK7/7m7/8bdDyTTg8A7gV8BpwFXL/T+7n6ru/e5e19PT0+sX1NERMoQ8zDas8DYEp6/a7EnuPuUmd0GvAW4EPj6MdYmIiIpihY27n5JpJcuHRNb8DCaiIjkR17P2Szk/GS+LdMqRESkbLkMGzN7k5m9qjYzuwT4TPLwrtnrRUQkn/La9fmrwOlm9hAwkCw7C3hH8vOfuftDmVQmIiJLltewuRN4P3Ae8B6gCdgD/BNwi7v/LMPaRERkiXIZNu7+98DfZ12HiIhURi7P2YiISG1R2IiISHQKGxERiU5hIyIi0SlsREQkOoWNiIhEp7AREZHoFDYiIhKdwkZERKJT2IiISHQKGxERiU5hIyIi0SlsREQkOoWNiIhEp7AREZHoFDYiIhKdwkZERKJT2IiISHQKGxERiU5hIyIi0SlsREQkOoWNiIhEp7AREZHoFDYiIhKdwkZERKJT2IiISHQKGxERiU5hIyIi0UUPGzNrMrPrzOx2M3vczCbMzM3s6jK2/aiZPWJmw2Z20MweMLP3xq5ZREQqK42WTTtwM3AFsB7YXc5GZvYV4A7geODvgLuANwD/YmafjFGoiIjEkUbYjAKXAie4+3rg24ttYGYXANcDzwJnuftn3P1a4FxgH/AVM9sYrWIREamo6GHj7hPufq+7v7iEzT6RzL/s7vtnvNYO4JtAC3Bl5aoUEZGY8tpB4B3J/L451t076zkiIpJzuQsbM2sHNgDD87SGtibzM9KrSkREliN3YQN0JvOD86wvLe9a6EXM7Boz6zez/sHBwYoVJyIiS1dW2JjZjqS7crnTXbELX4y73+rufe7e19PTk3U5IiJ1rVDm854FxpbwuruOoZaSUsulc571peUHlvEeIiKSorLCxt0viV3IjPcaMbOdwAYzO36O8zanJ/MtadUkIiLLk8dzNgA/SebvnmPde2Y9R0REci6vYfOtZP55M1tdWphcyHktMA7cnn5ZIiJyLMo9Z7MsZnYDsCl5eE4yv9LM3pb8/HN3v630fHd/yMy+CnwWeNLM7gGagQ8Ca4BPJRd4iohIFUglbAiHwy6ateyCZCq5beZKd7/ezP6V0JK5BpgGfgP8lbt/P2KtIiJSYamEjbu//Ri3u4MwGKeIiFSxvJ6zERGRGqKwERGR6BQ2IiISncJGRESiU9iIiEh0ChsREYlOYSMiItEpbEREJDqFjYiIRKewERGR6BQ2IiISncJGRESiU9iIiEh0ChsREYlOYSMiItEpbEREJDqFjYiIRKewERGR6BQ2IiISncJGRESiU9iIiEh0ChsREYlOYSMiItEpbEREJDqFjYiIRKewERGR6BQ2IiISncJGRESiU9iIiEh00cPGzJrM7Dozu93MHjezCTNzM7t6gW2uSJ4z3/SJ2HWLiEjlFFJ4j3bg5uTnPcBu4MQyt/0e8Pgcy/srUJeIiKQkjbAZBS4FHnf3F83sRuCLZW77XXe/I1ZhIiKSjuhh4+4TwL2x30dERPIrjZbNcpxjZp8GWoGdwP3uPpBxTSIiskR5D5vrZj0umtltwKfdfSyLgkREZOnyGjbbgU8BPwQGgE7gbcBfAh8HVgH/aaEXMLNrgGuSh+Nm9lS0auvLWmBv1kXUEO3PytL+rJwzK/li5u6LP8lsB/CaJbzuP7j7h+d5rRsJHQT+xN1vW8JrYmYnAk8Aq4Fz3P2JMrfrd/e+pbyXzE37srK0PytL+7NyKr0vy23ZPAss5bDVrmOoZVHu/oKZ/QD4z8CFhOAREZGcKyts3P2S2IUswWAyb8+0ChERKVs1DlfzlmS+bQnb3BqjkDqlfVlZ2p+Vpf1ZORXdl7kMGzN71XFCM2sws/8OvJVwAvC+cl/P3fUBrBDty8rS/qws7c/KqfS+LKuDwLLfxOwGYFPy8BzgbOAhYGuy7OczOwuYmQNPEc7J7CT0RtsMvJ4wIsH73f2H0QsXEZGKSCtsHgAuWuAp/8fdr5jx/L8C3gycDqwBpoHngf8HfNXdl3IITUREMpbKYTR3f7u72wLTFbOe/9/c/SJ3P8HdW919hbtvcvdPlhM0xzLS9IxtP2pmj5jZsJkdNLMHzOy9y/j1a5aZbVxkdO67s64xj8ys18y+bWa7zGzczHaY2c1mtjrr2qpNsu/m+/ztzrq+PDKzD5jZN8zsZ2Z2KNlXdy2yzQVm9gMz22dmR8zsSTP7tJk1lvu+eb2oc7mOaaRpM/sKcD3hQtK/A5qBy4F/MbNPufstccqtek8A351juS6kncXMTiUcQl5HGNX8aUIr/jrg3Wa22d2HMiyxGh3k6L/3mYbTLqRKfIFwKmOY8F23aaEnm9kfAP9MuPzlH4F9wO8DXyOc3vijst7V3WtuIoTEe4Djk8c3Ag5cvcA2FyTP+Xdg9YzlG4GhZEdvzPp3y9OU7BsH7si6lmqZgP+b7LNPzVr+1WT5t7KusZomYAewI+s6qmkCLiacojDg7cnn7q55nrsKeAkYB/pmLG8l/NHkwOXlvG8ue6Mtl7tPuPu97v7iEjYr3ZDty+6+f8Zr7QC+CbQAV1auSqk3SavmnYQvyG/OWv1FYAT4YzPTNWQSjbvf7+5bPUmNRXwA6AHudveX7yPmYWzKLyQP/7Sc963JsDlG70jmc3WpvnfWc+SVTjCzj5vZ/0jmZ2VdUE5dnMx/6O7TM1e4+2HgF8AK4Py0C6tyLWb24eTzd52ZXbyUcwmyoIW+Fx8k9A6+wMxaFnuhWj1nsyTJX5IbgOF5WkOlLtpnpFdVVfndZHpZ0gPxo+7+fCYV5VNpYMMt86zfSmj5nAH8OJWKasN64M5Zy7ab2ZXu/tMsCqoh835m3X3KzLYDrwNOAX670AupZRN0JvOD86wvLe9KoZZqMgr8BXAuYXDU1YQu7vcTjgX/WIeEXkGfs8q7HbiEEDjtwBuAvyWcT7zXzM7OrrSaULHPbG7DZpEujXNNC3bdk7ktZz+7+0vu/ufu/ht3P5BMDxL+Ov8VcBqwaHdzkWPl7je5+0/cfY+7j7r7U+7+CUKHizZC5yDJgTwfRktzpOlSOnfOs760/MAy3iOvKr6fk+b1bYRx7C4Evn6MtdWaev6cpe1bhMsYLsy6kCpXsc9sbsPGUxxp2t1HzGwnsMHMjp/jvM3pyXy+Y+1VK+J+1ujcr/ZMMp/v3F/Nfs4yoM9fZTwD9BE+s4/OXGFmBeBkYIoyBkbO7WG0DPwkmb97jnXvmfUcWVypR5WGFjrq/mT+TjN7xb89M+sgXCA3CjycdmE1SJ+/yljoe/FCQu/Jh9x9fLEXUtgc9a1k/vmZw4aY2UbgWsJFTbenX1Z+mdmbZn9pJssvAT6TPNS5tIS7P0u41flGwmdqppsIf4Xf6e4jKZdWlczsd+bqgJL8my2N9qHP3/LcQxhl//KZo/GbWSvwpeTh35TzQqkMxJmFpY40nWzzv4DPEoZwuIcwEsEHgW7CFd8armaGpHvz6YT9OpAsPoujffP/zN2/NMemdWuO4Wp+Szi3dTHh8NkFruFqymLhFvPXE673eA44DJwK/B7hCvcfEEaIn8iqxjwys/cB70sergfeRWgB/ixZttfdPzfr+fcQzu3eTRiu5jJCt+h7gP9Y1gWiWQ+dEHFIhgcIQynMN90xz3ZXAL8mXM19GPgp8N6sf588TsDHgO8TrogfJrT+nieMn/Qfsq4vrxNhnL7bgReBCcIX5c3MGCZJU1n78SLgO4Tx5Q4Ak4RzNT8CPkLyx7SmV+23Gxf5btwxxzabCeG9HzgC/Cvh6EVjue9bsy0bERHJD52zERGR6BQ2IiISncJGRESiU9iIiEh0ChsREYlOYSMiItEpbEREJDqFjYiIRKewERGR6BQ2IiIS3f8Hl10sAEaSL1cAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f575a31b390>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.xlim(-10, 10)\n",
    "plt.ylim(-15, 10)\n",
    "e1=plot_ellipse(plt.gca(), model.theta_[0], np.identity(2)*model.sigma_[0], color=\"red\", label=\"0\")\n",
    "e2=plot_ellipse(plt.gca(), model.theta_[1], np.identity(2)*model.sigma_[1], color=\"blue\", label=\"1\")\n",
    "plt.legend([e1, e2], [\"0\", \"1\"]);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The problem with naive Bayesian classification is that it tries to model the data using Gaussian distributions, which are aligned along the x and y axes. With this example data we would have needed Gaussian distributions which are \"tilted\"."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Text classification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We next try to classify a set of messages that were posted on a public forum. The messages were divided into groups by the topics. So, we have a data set ready for classification testing. Let's first load this data using scikit-learn and print the message categories."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:19.444349Z",
     "start_time": "2020-06-24T19:29:19.195328Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['alt.atheism',\n",
       " 'comp.graphics',\n",
       " 'comp.os.ms-windows.misc',\n",
       " 'comp.sys.ibm.pc.hardware',\n",
       " 'comp.sys.mac.hardware',\n",
       " 'comp.windows.x',\n",
       " 'misc.forsale',\n",
       " 'rec.autos',\n",
       " 'rec.motorcycles',\n",
       " 'rec.sport.baseball',\n",
       " 'rec.sport.hockey',\n",
       " 'sci.crypt',\n",
       " 'sci.electronics',\n",
       " 'sci.med',\n",
       " 'sci.space',\n",
       " 'soc.religion.christian',\n",
       " 'talk.politics.guns',\n",
       " 'talk.politics.mideast',\n",
       " 'talk.politics.misc',\n",
       " 'talk.religion.misc']"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.datasets import fetch_20newsgroups\n",
    "data = fetch_20newsgroups()\n",
    "data.target_names"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We concentrate on four message categories only. The tool `fetch_20newsgroups` allows us to easily split the data into training and testing data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:19.778817Z",
     "start_time": "2020-06-24T19:29:19.445242Z"
    }
   },
   "outputs": [],
   "source": [
    "categories = ['comp.graphics', 'rec.autos', 'sci.electronics', 'sci.crypt']\n",
    "train = fetch_20newsgroups(subset='train', categories=categories)\n",
    "test = fetch_20newsgroups(subset='test', categories=categories)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see what we got:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:19.783690Z",
     "start_time": "2020-06-24T19:29:19.780035Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training data: Data: <class 'list'> 2364 Target: <class 'numpy.ndarray'> 2364\n",
      "Test data: Data: <class 'list'> 1574 Target <class 'list'> 1574\n"
     ]
    }
   ],
   "source": [
    "print(\"Training data:\", \"Data:\", str(type(train.data)), len(train.data), \"Target:\", str(type(train.target)), len(train.target))\n",
    "print(\"Test data:\", \"Data:\", str(type(test.data)), len(test.data), \"Target\", str(type(test.data)), len(test.target))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use as features the frequencies of each word in the dataset. That is, there are as many features as there are distinct words in the dataset. We denote the number of features by $f$. As the features are now counts, it is sensible to use multinomial distribution instead of Gaussian. \n",
    "\n",
    "Let's try to model these messages using multinomial distributions. Each message category has its own distribution. A multinomial distribution has $f$ non-negative parameters $\\theta_1,\\ldots , \\theta_f$, which sum up to one. For example, the parameter $\\theta_3$ might tell the the probability of the word \"board\" appearing in a message of the category this distribution is describing.\n",
    "\n",
    "In scikit-learn there is a class `CountVectorizer` that converts messages in form of text strings to feature vectors. We can integrate this conversion with the model we are using (multinomial naive Bayes), so that the conversion happens automatically as part of the `fit` method. We achive this integration using the `make_pipeline` tool."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:20.367332Z",
     "start_time": "2020-06-24T19:29:19.784456Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy score is 0.920584498094\n"
     ]
    }
   ],
   "source": [
    "#from sklearn.feature_extraction.text import TfidfVectorizer  # an alternative feature extractor\n",
    "from sklearn.feature_extraction.text import CountVectorizer\n",
    "from sklearn.naive_bayes import MultinomialNB\n",
    "from sklearn.pipeline import make_pipeline\n",
    "\n",
    "#model = make_pipeline(TfidfVectorizer(), MultinomialNB())\n",
    "model = make_pipeline(CountVectorizer(), MultinomialNB())\n",
    "model.fit(train.data, train.target)\n",
    "labels_fitted = model.predict(test.data)\n",
    "print(\"Accuracy score is\", accuracy_score(labels_fitted, test.target))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The classifier seem to work quite well! Notice that now we used separate data for testing the model.\n",
    "\n",
    "Let's have a closer look at the resulting feature vectors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:20.734692Z",
     "start_time": "2020-06-24T19:29:20.368401Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Type of feature matrix: <class 'scipy.sparse.csr.csr_matrix'>\n",
      "  (0, 20579)\t1\n",
      "  (0, 19220)\t1\n",
      "  (0, 29697)\t1\n",
      "  (0, 6320)\t1\n",
      "  (0, 25926)\t1\n",
      "  (0, 34222)\t1\n",
      "  (0, 31398)\t1\n",
      "  (0, 17883)\t1\n",
      "  (0, 16809)\t1\n",
      "  (0, 34425)\t1\n",
      "  (0, 23460)\t1\n",
      "  (0, 21787)\t1\n",
      "  (0, 11068)\t1\n",
      "  (0, 29494)\t1\n",
      "  (0, 29505)\t1\n",
      "  (0, 18436)\t1\n",
      "  (0, 24025)\t1\n",
      "  (0, 25336)\t1\n",
      "  (0, 12577)\t1\n",
      "  (0, 27517)\t1\n",
      "  (0, 30641)\t1\n",
      "  (0, 5980)\t1\n",
      "  (0, 29104)\t1\n",
      "  (0, 27521)\t1\n",
      "  (0, 11100)\t1\n",
      "  :\t:\n",
      "  (0, 17310)\t1\n",
      "  (0, 25400)\t1\n",
      "  (0, 23118)\t1\n",
      "  (0, 31686)\t6\n",
      "  (0, 27158)\t1\n",
      "  (0, 18085)\t1\n",
      "  (0, 12580)\t1\n",
      "  (0, 2100)\t1\n",
      "  (0, 20381)\t1\n",
      "  (0, 32729)\t1\n",
      "  (0, 23854)\t2\n",
      "  (0, 11079)\t1\n",
      "  (0, 15109)\t2\n",
      "  (0, 20509)\t1\n",
      "  (0, 23858)\t1\n",
      "  (0, 26624)\t1\n",
      "  (0, 30377)\t1\n",
      "  (0, 16034)\t1\n",
      "  (0, 19099)\t1\n",
      "  (0, 13317)\t6\n",
      "  (0, 34790)\t6\n",
      "  (0, 9553)\t4\n",
      "  (0, 21852)\t5\n",
      "  (0, 18962)\t3\n",
      "  (0, 15373)\t1\n"
     ]
    }
   ],
   "source": [
    "vec=CountVectorizer()\n",
    "features=vec.fit_transform(train.data)\n",
    "print(\"Type of feature matrix:\", type(features))\n",
    "print(features[0,:])        # print the features of the first sample point"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The feature matrix is stored in sparse format, that is, only the nonzero counts are stored. How many words were in the first message?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-06-24T19:29:20.739668Z",
     "start_time": "2020-06-24T19:29:20.735790Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of words: 177\n",
      "Word 'it' appears in the first message 2 times.\n",
      "\n",
      "From: jgfoot@minerva.cis.yale.edu (Josh A. Goldfoot)\n",
      "Subject: Re: Organized Lobbying for Cryptography\n",
      "Organization: Yale University\n",
      "Lines: 21\n",
      "Distribution: inet\n",
      "Reply-To: jgfoot@minerva.cis.yale.edu\n",
      "NNTP-Posting-Host: minerva.cis.yale.edu\n",
      "X-Newsreader: TIN [version 1.1 Minerva PL9]\n",
      "\n",
      "Shaun P. Hughes (sphughes@sfsuvax1.sfsu.edu) wrote:\n",
      ": In article <1r3jgbINN35i@eli.CS.YALE.EDU> jgfoot@minerva.cis.yale.edu writes:\n",
      "[deletion]\n",
      ": >Perhaps these encryption-only types would defend the digitized porn if it\n",
      ": >was posted encrypted?\n",
      ": >\n",
      ": >These issues are not as seperable as you maintain.\n",
      ": >\n",
      "\n",
      ": Now why would anyone \"post\" anything encrypted? Encryption is only of \n",
      ": use between persons who know how to decrypt the data.\n",
      "\n",
      ": And why should I care what other people look at? \n",
      "\n",
      "I was responding to another person (Tarl Neustaedter) who held that the\n",
      "EFF wasn't the best organization to fight for crytography rights since the\n",
      "EFF also supports the right to distribute pornography over the internet,\n",
      "something some Crypto people might object to. In other words, he's\n",
      "implying that there are people who will protect any speech, just  as long\n",
      "as it is encrypted.\n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(\"Number of words:\", features[0,:].sum())\n",
    "col = vec.vocabulary_[\"it\"]   # Get the column of 'it' word in the feature matrix\n",
    "print(f\"Word 'it' appears in the first message {features[0, col]} times.\")\n",
    "print()\n",
    "print(train.data[0])   # Let's print the corresponding message as well\n",
    "#print(vec.get_feature_names())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### <div class=\"alert alert-info\">Exercise 1 (blob classification)</div>\n",
    "\n",
    "Write function `blob_classification` that gets feature matrix X and label vector y as parameters. It should then return the accuracy score of the prediction. Do the prediction using `GaussianNB`, and use `train_test_split` function from `sklearn` to split the dataset in to two parts: one for training and one for testing. Give parameter `random_state=0` to the splitting function so that the result is deterministic. Use training set size of 75% of the whole data.\n",
    "<hr/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### <div class=\"alert alert-info\">Exercise 2 (plant classification)</div>\n",
    "\n",
    "Write function `plant_classification` that does the following:\n",
    "\n",
    "* loads the iris dataset using sklearn (`sklearn.datasets.load_iris`)\n",
    "* splits the data into training and testing part using the `train_test_split` function so that the training set size is 80% of the whole data (give the call also the `random_state=0` argument to make the result deterministic)\n",
    "* use Gaussian naive Bayes to fit the training data\n",
    "* predict labels of the test data\n",
    "* the function should return the accuracy score of the prediction performance (`sklearn.metrics.accuracy_score`)\n",
    "<hr/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### <div class=\"alert alert-info\">Exercise 3 (word classification)</div>\n",
    "\n",
    "This exercise can give four points at maximum!\n",
    "\n",
    "In this exercise we create a model that tries to label previously unseen words to be either Finnish or English.\n",
    "\n",
    "Part 1.\n",
    "\n",
    "Write function `get_features` that gets a one dimensional np.array, containing words, as parameter. It should return a feature matrix of shape (n, 29), where n is the number of elements of the input array. There should be one feature for each of the letters in the following alphabet: \"abcdefghijklmnopqrstuvwxyzäö-\". The values should be the number of times the corresponding character appears in the word.\n",
    "\n",
    "Part 2.\n",
    "\n",
    "Write function `contains_valid_chars` that takes a string as a parameter and returns the truth value of whether all the characters in the string belong to the alphabet or not.\n",
    "\n",
    "Part 3.\n",
    "\n",
    "Write function `get_features_and_labels` that returns the tuple (X, y) of the feature matrix and the target vector. Use the labels 0 and 1 for Finnish and English, respectively. Use the supplied functions load_finnish() and load_english() to get the lists of words. Filter the lists in the following ways:\n",
    "\n",
    "* Convert the Finnish words to lowercase, and then filter out those words that contain characters that don't belong to the alphabet.\n",
    "* For the English words first filter out those words that begin with an uppercase letter to get rid of proper nouns. Then proceed as with the Finnish words.\n",
    "\n",
    "Use get_features function you made earlier to form the feature matrix.\n",
    "\n",
    "Part 4.\n",
    "\n",
    "We have earlier seen examples where we split the data into learning part and testing part. This way we can test whether the model can really be used to predict unseen data. However, it can be that we had bad luck and the split produced very biased learning and test datas. To counter this, we can perform the split several times and take as the final result the average from the different splits. This is called [cross validation](<https://en.wikipedia.org/wiki/Cross-validation_(statistics)>).\n",
    "\n",
    "Create `word_classification` function that does the following:\n",
    "\n",
    "Use the function `get_features_and_labels` you made earlier to get the feature matrix and the labels. Use multinomial naive Bayes to do the classification. Get the accuracy scores using the `sklearn.model_selection.cross_val_score` function; use 5-fold cross validation. The function should return a list of five accuracy scores.\n",
    "\n",
    "The cv parameter of `cross_val_score` can be either an integer, which specifies the number of folds, or it can be a *cross-validation generator* that generates the (train set,test set) pairs. What happens if you pass the following cross-validation generator to `cross_val_score` as a parameter: `sklearn.model_selection.KFold(n_splits=5, shuffle=True, random_state=0)`.\n",
    "\n",
    "Why the difference?\n",
    "<hr/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### <div class=\"alert alert-info\">Exercise 4 (spam detection)</div>\n",
    "\n",
    "This exercise gives two points if solved correctly!\n",
    "\n",
    "In the `src` folder there are two files: `ham.txt.gz` and `spam.txt.gz`. The files are preprocessed versions of the files from https://spamassassin.apache.org/old/publiccorpus/. There is one email per line. The file `ham.txt.gz` contains emails that are non-spam, and, conversely, emails in file `spam.txt` are spam. The email headers have been removed, except for the subject line, and non-ascii characters have been deleted.\n",
    "\n",
    "Write function `spam_detection` that does the following:\n",
    "\n",
    "* Read the lines from these files into arrays. Use function `open` from `gzip` module, since the files are compressed. From each file take only `fraction` of lines from the start of the file, where `fraction` is a parameter to `spam_detection`, and should be in the range `[0.0, 1.0]`.\n",
    "* forms the combined feature matrix using `CountVectorizer` class' `fit_transform` method. The feature matrix should first have the rows for the `ham` dataset and then the rows for the `spam` dataset. One row in the feature matrix corresponds to one email.\n",
    "* use labels 0 for ham and 1 for spam\n",
    "* divide that feature matrix and the target label into training and test sets, using `train_test_split`. Use 75% of the data for training. Pass the random_state parameter from `spam_detection` to `train_test_split`.\n",
    "* train a `MultinomialNB` model, and use it to predict the labels for the test set\n",
    "\n",
    "The function should return a triple consisting of\n",
    "\n",
    "* accuracy score of the prediction\n",
    "* size of test sample\n",
    "* number of misclassified sample points\n",
    "\n",
    "Note. The tests use the `fraction` parameter with value 0.1 to ease to load on the TMC server. If full data were used and the solution did something non-optimal, it could use huge amounts of memory, causing the solution to fail.\n",
    "<hr/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<!--NAVIGATION-->\n",
    "\n",
    "\n",
    "<a href=\"https://colab.research.google.com/github/saskeli/x/blob/master/bayes.ipynb\"><img align=\"left\" src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open in Colab\" title=\"Open and Execute in Google Colaboratory\"></a>\n"
   ]
  }
 ],
 "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.9"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}