Lab5.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-success\">\n",
    "<h1>Feature Selection Lab</h1>\n",
    "<h2>ASI Lab - 01/06/2017</h2>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>1. Aims</h2>\n",
    "<p>\n",
    "<ul>\n",
    "<li>To pick some features from the leukaemia data.\n",
    "</ul>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.cm as cm\n",
    "import seaborn as sns\n",
    "import scipy\n",
    "import scipy.io as sio\n",
    "import time as time\n",
    "%matplotlib inline\n",
    "\n",
    "plt.rcParams['figure.figsize'] = (20.0, 10.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>2. Feature Selection</h2>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Like last week, this lab is quite open-ended. On the course webpage, you can download a (slightly reduced) version of the leukaemia data. It consists of some training data $X$ and associated labels $t$ as well as independent test data $X_test$ and $t_test$. Your task is to select some features and classify the data – in previous labs we have provided code for classiﬁcation (SVM, KNN, Bayes) and you may use whichever you like."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-warning\">\n",
    "Let's load the data from the <code>labdata.mat</code> file.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 94,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data = sio.loadmat('./Lab5/labdata.mat')\n",
    "X = data['X']\n",
    "X_test = data['X_test']\n",
    "t = data['t']\n",
    "t_test = data['t_test']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>2.1 Feature Selection</h3>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You are free to use any of the feature selection schemes we covered in class – scoring the features, PCA, or clustering them (you should have K-means code from last week to do this). In each case, you should use the test set to investigate how performance varies as you change the number."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h4>2.1.1 Scoring</h4>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use the score described in the lectures (only on the training data) to get a set of features (corre-sponding to a subset of the original ones). You should make sure you extract this subset from the test data too!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def compute_s(X, t):\n",
    "    for i in range(X.size):\n",
    "        train_x = np.delete(X, i)\n",
    "        train_t = np.delete(t, i)\n",
    "        pred = makePredictions(train_x, train_t, test_x, k)\n",
    "    X_0 = np.array([X[_] for _ in range(len(X)) if t[_]==0])\n",
    "    X_1 = np.array([X[_] for _ in range(len(X)) if t[_]==1])\n",
    "    s = np.abs(X_1.mean() - X_0.mean())/(X_0.var() - X_1.var())\n",
    "    return s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 118,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "38"
      ]
     },
     "execution_count": 118,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X[:,1].size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_best_features(X):\n",
    "    s = []\n",
    "    for x in range(X.shape[1]):\n",
    "        s.append((compute_s(X[:,x], t), x))\n",
    "    s.sort(key=lambda x : x[0], reverse=True)\n",
    "    return np.array([x[1] for x in s])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {},
   "outputs": [],
   "source": [
    "best_features_scoring = get_best_features(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 116,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlcAAAJACAYAAACg6EOnAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8nFWh//HvTCZ7s7VJ2nRLup7uG6ULrez7Jgh4EUVE\nNper96rXKyrX5cpVL4jX5acoXriCqGyCIKuAYEsXW1qWFtrTpvvetEmbNHsy8/sjaUibSZfJyTzz\nTD7v18vXK3POZOZbDzP99jnPMxOIRCICAACAG0GvAwAAACQTyhUAAIBDlCsAAACHKFcAAAAOUa4A\nAAAcolwBAAA4FIrll4wxQUm/lDRVUqOkm6215Z3mvyLpOklhSd+31j7lICsAAEDCi/XI1RWSMqy1\ncyXdLumewxPGmHxJ/yJprqTzJf2kpyEBAAD8IqYjV5LmS3pRkqy1S40xMzvN1UraIim7/X/hE3nA\nioqapPg004KCLFVV1XkdAw6xpsmF9Uw+rGly8ct6FhXlBLqbi7Vc5Uo62Ol2qzEmZK1tab+9TdL7\nklIk/eBEHrCgIEuhUEqMcRJLUVGO1xHgGGuaXFjP5MOaJhe/r2es5apaUuc/ebBTsbpIUomkEe23\nXzLGLLLWLjvWA/qhpZ6IoqIcVVTUeB0DDrGmyYX1TD6saXLxy3oeqwDGes7VIkkXS5IxZo6kVZ3m\nqiTVS2q01jZIOiApP8bnAQAA8JVYj1w9Jek8Y8xiSQFJNxpjviyp3Fr7jDHmXElLjTFhSW9IetlN\nXAAAgMQWU7my1oYlfeao4bWd5r8t6ds9yAUAAOBLfIgoAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJc\nAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIA\nAHCIcgUAAOAQ5QoAAMAhyhUAAIBDIa8DAMCx1G3bqx1P/l2BlBQNu/YcpRfmeR0JAI6JcgUgYdm7\n/6jN9z+v5qoaSdKmX/9Fo77wEY289TKPkwFA99gWBJCQ9rzypsp//mRHsZKkxr1VWnf3H3Vw1UYP\nkwHAsVGuACSkXc8sUqSxuct4S3Wdtj32mgeJAODEUK4AJKTW2obu5+q6nwMAr1GuACSknAll3c7l\nTxsTvyAAcJIoVwAS0sjbLlf+9K4lasD8yRp27dkeJAKAE8PVggASUqhfpk59+Jta/z+P68Bb6xUI\nBlUwa5zGfuVaBVN56wKQuHiHApCw0gvzNem/bvE6BgCcFLYFAQAAHKJcAQAAOMS2IAAgYYSbmrX5\ngRdUuWyNAikBDZg3RaXXn6dASorX0YATRrkCACSEcHOL3vzUD7X31RUdY7ueWazKJas1/d4vKxBk\nswX+wH+pAICEsOWhF48oVoftfGaRdv1lsQeJgNhQrgAACaFq2droE+GI9i18N75hgB6gXAEAEkPK\nMf5KOtYckGD4rxUAkBCKz54uBQJdxgNpqSq5ZI4HiYDYUK4AAAlhyFVnaui1ZyuQ+sGVgYH0VJXd\neJGKTp/mYTLg5HC1IAAgIQQCAU39n3/W4EtP055X31QgGNSgS+aq8LRJXkcDTgrlCgCQMAKBgIrP\nPUXF557idRQgZmwLAgAAOES5AgAAcIhtQQC+FYlEtPFXz2j3C/9Q0/6Dyi4dqOHXX6BBF832OhqA\nPoxyBcC31vzng9r4q6elcESSVFu+Q5VvWk1p/pwGXz7P43QA+iq2BQH4UlNVjXY8uaCjWB3WcrBW\nW377gkepAIByBcCn9i18V427K6POHSrfoUhra5wTAUAbyhUAX8oaPlCB9NSoc6m5WVKQtzcA3uDd\nB4Av5U8brf6zxkedKzxjugJRvkYFAOKBcgXAtyb/8DYVzB7f8aW+KdkZKrl8nsZ/6waPkyHR1G7Z\nI/ujR2Xv/qNq1m3zOg6SHFcLAvCtfqOH6LSnv6+9r6xQ7aZd6j93ovInj/Q6FhLM+p8+oY33Pq3m\nqhpJ0qb7/qLST16o8f/xSY+TIVlRrgD4WiAQ0MDzZnodAwmqcoVV+U+fUGttQ8dYS3WdNt73jApO\nHadBF87yMB2SFduCAICktfNPC44oVodFmlq06/klHiRCX0C5AgAkrdb6xm7nwg1NcUyCvoRyBQBI\nWgUzx3U7lzdlVByToC+hXAEAktbQfzpLhWdO6zJeMHu8Rnz6Eg8SoS/ghHYAQNIKhlJ06m+/rvU/\n+5MOLFurSCSs/OljNfpfr1ZKVrrX8ZCkKFcAgKSWkpmucV+7zusY6EPYFgQAAHCIcgUAAOAQ24IA\n0EOtjc3a8tCLqnrTKpga0sDzT1XJZafx/YZAH0W5AoAeaG1o0vLr/0v7FrzTMbbjyQUqW/KeJv3g\nVg+TAfAK24IA0AMb7336iGIlSWoNa+sfXlHlm2u9CQXAUxy5AoAeqOqmQIUbmrT7uaXqf4wPsfRa\nS22DNj/wnGrWbVdqXraGf/w85Y4v9ToW4HuUKwDoiWD351Ul8jlX9Tv3afkN31f1uxs7xnY88XeN\n//anNPxj53iYDPA/tgUBoAcGzJkYdTwlO0MlV34ozmlOnP3vPxxRrCSpuapG5T994pjfxwfg+ChX\nANADI269TMUXnHrEWDA9VWU3XaL8ySM9SnV8B1asizpet2mXdj6zKM5pgOTCtiAA9EAwNaRT/+92\n7XhqofYvWq1gWqpKLj9NhfMmex3tmCKRSPdzreE4JgGST0zlyhgTlPRLSVMlNUq62Vpb3mn+Iknf\nlhSQtELS56213b+SAcDHAikpGnr1mRp69ZleRzlh+dPHqLZ8R5fxzNJBGnJF4m5nAn4Q67bgFZIy\nrLVzJd0u6Z7DE8aYHEl3S7rUWjtb0mZJhT3MCQBwaOxXP6aciWVHjIXysjXq81fwhcZAD8W6LThf\n0ouSZK1daoyZ2WnuNEmrJN1jjBkp6X+ttRU9iwkAcCm7dKDmPnWnNt33rGo37FAoN0vDP3au8qeP\n8Toa4HuxlqtcSQc73W41xoSstS1qO0p1lqRpkg5JWmiMWWKtjX72ZLuCgiyFQikxxkksRUU5XkeA\nY6xpcmE92xXlaMhdt3idwgnWNLn4fT1jLVfVkjr/yYPtxUqS9ktabq3dLUnGmAVqK1rHLFdVVXUx\nRkksRUU5qqio8ToGHGJNkwvrmXxY0+Til/U8VgGM9ZyrRZIuliRjzBy1bQMetlLSJGNMoTEmJGmO\npPdjfB4AAABfifXI1VOSzjPGLFbbFYE3GmO+LKncWvuMMebrkl5qv+9j1trVDrICAAAkvJjKlbU2\nLOkzRw2v7TT/iKRHepALAADAl/iEdgAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ\n5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoV\nAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAA\nAIcoVwAAAA5RrgAAAByiXAEAADgU8joAgORQuWyNtj/+ulpq6tRv7DCNvPUyhfpleh0Lx9F88JAa\nKw4qc0ihUjLTvY4DJAXKFYAe2/jrZ7Tu7j+qpaa+Y2z380t16sN3KHNQfw+ToTutdY1adfuvVfHa\nSjXuPaCsskEquXyexn3jEwoEAl7HA3yNbUEAPdJcXasNv3r6iGIlSdWrNmr9jx7xKBWO552v/D9t\nf/Rvatx7QJJUt3m3Nvz8Sa3/8eMeJwP8j3IFoEe2P/66GnfujzpXtWJdnNPgRNTv3KeKv73VdSIS\n0a7nFisSicQ/FJBEKFcAeiSQ0v3bCLtLial69SY1HzgUda5h936F65vinAhILpQrAD0y9OozlTmk\nKOpc/ikmzmlwInInj1Bq/5yoc5klhQpmpsU5EZBcKFcAeiTUL1OjvvgRhfKyjxjPnzFW5t8/5lEq\nHEtmSaGKzzml60QwoMGXz+OEdqCHuFoQQI+VfeoiFZxitO2RV9VcXafc8aUqu/EiLu1PYFPu+bxS\nMtK09/DVgqWDNOTK+Rr1xau8jgb4HuUKgBN5k0cqb/JIr2PgBKWkp2rKjz6nltoGNVfVKL04X8G0\nVK9jAUmBcgUAfVgoO0Oh7AyvYwBJhXOuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAA\ngEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACH\nKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGu\nAAAAHKJcAQAAOBSK5ZeMMUFJv5Q0VVKjpJutteVR7vOcpKettb/qaVAAAAA/iPXI1RWSMqy1cyXd\nLumeKPe5U1JBrMEAAAD8KNZyNV/Si5JkrV0qaWbnSWPM1ZLCh+8DAADQV8S0LSgpV9LBTrdbjTEh\na22LMWaSpOskXS3pWyf6gAUFWQqFUmKMk1iKinK8jgDHWNPkwnomH9Y0ufh9PWMtV9WSOv/Jg9ba\nlvafPylpiKS/SSqT1GSM2WytPeZRrKqquhijJJaiohxVVNR4HQMOsabJhfVMPqxpcvHLeh6rAMZa\nrhZJukzSY8aYOZJWHZ6w1v774Z+NMd+RtPt4xQoAACBZxFqunpJ0njFmsaSApBuNMV+WVG6tfcZZ\nOgDwuZbaeu16bqlSc7M08LyZCqQkx+kPALoXU7my1oYlfeao4bVR7vedWB4fAJJB+S+e1JYHXlD9\n9gpJUu6kERp3x/UqPmuGx8kA9CY+RBQAesHuF/6hdXc90lGsJKl69Sat+tqv1Vxd62EyAL2NcgUA\nvWDHnxcq3NDUZbx+yx5t/r8XPEgEIF4oVwDQC5qrDnU/V1kdxyQA4o1yBQC9IKtsYLdzORNHxDEJ\ngHijXAFALxh5y2XKHN61YPWfO1FDrzrdg0QA4iXWj2IAABxDvzFDNePXX9GG//eUDr67QcH0VPWf\nM0Hj7/gkH8cAJDnKFQD0koIZYzXzga95HQNAnLEtCAAA4BDlCgAAwCHKFQAAgEOccwXAd/YtfEf7\nFr+ntLx+Gn79eQplZ3odCQA6UK4A+Ea4qVkrb7tHe15ZoUhTsyRp82+f18Q7b9bAc2d6nA4A2rAt\nCMA37N2PaPfzSzuKlSTVbdqtNd99UOFOYwDgJcoVAN/Y/8aqqOOH1m3TjicXxDkNAERHuQLgG611\njd3ONR3ju/wAIJ4oVwB8I2diadTx1IIclVwyJ85pACA6yhUA3xj1uSuVVXrU9/UFAxp69RnKivI9\nfgDgBa4WBOAbeZNGaOZD39Sm+57RofXbFeqXpYEXnKrSGy70OhoAdKBcAfCV3HHDNfXH/+x1DADo\nFtuCAAAADlGuAAAAHKJcAQAAOES5AgBIklobmvike8ABTmgHgD6uctkalf/0CR18Z4MCqSkqmDlO\n4+64Xtmlg7yOBvgS5QoA+rBDm3bprc/9WPXbKjrGdj2zSLUbd2resz9USma6h+kAf2JbEAD6sM2/\n+csRxeqw6tWbtOWhlzxIBPgf5QoA+rD67V2L1WG1m3bFMQmQPChXANCHpQ3IjWkOQPcoVwDQhw29\n9lyF8vt1Gc8cVqyyT1/sQSLA/yhXANCHDZg9XhO/e6P6jRveNpASVN6MMZry488rfUCet+EAn+Jq\nQQDo44Zde46GXHWGKpetUTA9TQWnjFUgEPA6FuBblCsAgIKpIRXOm+x1DCApsC0IAADgEOUKAADA\nIcoVAF8LNzWrpbbe6xgA0IFzrgD4UlNltd77jwdUufQ9tdQ1Knf8cI245VINumiO19EA9HGUKwC+\nE4lEtOLmu7R/0eqOsf2LVqt67VaFcrM5MRuAp9gWBOA7u55bov1L3usy3ry/Wlsf+qsHiQDgA5Qr\nAL5T8/4WKRyJOle3fU+c0wDAkShXAHwnc0hRt3PphflxTAIAXVGuAPjO0I+eqdxJI7uMBzPSNPjK\n0z1IBAAfoFwB8J1gakjTfvZFDZg/WYH0tutyskcNkfn6xzXkivkepwPQ13G1IABfyp1Yprl/+p6q\n125R0/5qFZxilJKR5nUsAKBcAfC33HGlXkcAgCOwLQgAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcA\nAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAA\nHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIdCsfySMSYo6ZeSpkpqlHSztba8\n0/yXJF3bfvN5a+13exoUAADAD2I9cnWFpAxr7VxJt0u65/CEMWakpI9LOk3SHEnnG2Om9DQoAACA\nH8RaruZLelGSrLVLJc3sNLdN0oXW2lZrbURSqqSGHqUEAADwiZi2BSXlSjrY6XarMSZkrW2x1jZL\n2meMCUi6W9Jb1tp1x3vAgoIshUIpMcZJLEVFOV5HgGOsaXJhPZMPa5pc/L6esZarakmd/+RBa23L\n4RvGmAxJD0iqkfS5E3nAqqq6GKMklqKiHFVU1HgdAw6xpsmF9Uw+rGly8ct6HqsAxrotuEjSxZJk\njJkjadXhifYjVk9Lesdae5u1tjXG5wAAAPCdWI9cPSXpPGPMYkkBSTcaY74sqVxSiqQzJKUbYy5q\nv//XrbVLepwWAAAgwcVUrqy1YUmfOWp4baefM2JOBAAA4GN8iCgAAIBDlCsAAACHKFcAAAAOUa4A\nAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoV4IFIa6vCLXwzFAAko1i//gZADGq37Nba7z+s\nA8vXKhyOKH/qKI3516uVP32s19EAAI5QroA4aa1v1Iqb7lL1qo0dY3t27Vf12q2a+8R/KmtYsYfp\nAACusC0IxMnm375wRLE6rH7zbm267y8eJAIA9AbKFRAntRt3djtXt21vHJMAAHoT5QqIk7T+ud3P\nFeTEMQkAoDdRroA4Kfv0xcoYUthlPJTXT8M+do4HiQAAvYFyBcRJxsD+mnzXZ5Q7ZaQUCEiS+plh\nmvCdG9R/1niP0wEAXOFqQSCOBp47U8Vnz9D+Je8p3NCkwtOnKpjKyxAAkgnv6kCcBYJBFc6b7HUM\nAEAvYVsQAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAA\nwCHKFQAAgEN8/Q0AdBJuata6ex7V/jdWqaW+UbkTyzTq81cqd1yp19EA+ATlCgDaRSIRrbztHu1+\nfmnHWM17m1X1ptWs392hfqOHeJgOgF+wLQgA7fYteEd7Xnmzy3jdxl3aeO/THiQC4EeUKwBot3/x\nakWaWqLOHVq/Pc5pAPgV5QoA2qXm9et2LiUnM45JAPgZ5QoA2g2//nxlDh/YdSIY0MDzZsY/EABf\nolwBQLvUnCxN/N5Nyho1uGMslJetsk9fotIbLvQwGQA/4WpBAOhk0IWzVHTmNG1//HU1V9eq5JI5\nyi4r8ToWAB+hXAHAUVIy0lR6/flexwDgU2wLAgAAOES5AgAklUgkopbaekXCYa+joI9iWxAAkBQi\nkYg2/OIp7fzzG6rfsU/pxfkquXiOxn71WgWCHEtA/FCuAABJofynf5K96w9Sa9sRq+bKaq1fu1Wt\n9Y2a8J0bPU6HvoQqjz4rEolo38J3tP7Hj2nrH19RuDn6J3MDSHzhllbtfHJBR7HqbNdfFqultt6D\nVOirOHKFPqm1rlErbv2RKv7+VsfXnWz6zbOa8qPPqWDGWI/TAThZTZXVqtu2J+pc/fYKHdqwU/lT\nRsU5FfoqjlyhT1pz54Pa+/LyI75Hrua9zXrvjv9VJBLxMBmAWKTmZiutKD/qXNqAXGUNLYpzIvRl\nlCv0SfsWroo6fuDtcu1b+G6c0wDoqZSMNBWfc0rUuaKzpiutf26cE6EvY1sQfVJrd+dftIbVsGt/\nfMMAcGLid29UuL5Re15+U037Dio1v5+KzpquyXd91uto6GMoV+iTcsaXqn7Hvi7jGYMGaOAFszxI\nBKCngmmpmvqTL6hhb5WqV21UzrjhyhzCdiDij21B9Ekjbr1MacUFRw6GUjT02rOVlt/Pm1AAnMgo\nLlDxOadQrOAZjlyhTyo6Y5pm/u9XtfmBF1S3ebdS+/fToIvnqPT6C7yOBgDwOcoV+qz+syeo/+wJ\nXscAACQZtgUBAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGu\nAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwKFQrL9ojAlK+qWkqZIaJd1srS3vNH+LpNsktUi6\n01r7bA+zAgAAJLyeHLm6QlKGtXaupNsl3XN4whgzSNIXJc2TdIGkHxhj0nsSFAAAwA96Uq7mS3pR\nkqy1SyXN7DQ3S9Iia22jtfagpHJJU3rwXADQKyKRiFobmxWJRLyOAiBJxLwtKClX0sFOt1uNMSFr\nbUuUuRpJeT14LgBwbuOvntGOJxeofkeF0ovzNeiSuRr75Y8qEOR0VACx60m5qpaU0+l2sL1YRZvL\nkXTgWA9WUJClUCilB3ESR1FRzvHvBF9hTZNLUVGO3r37Ea258yFFmtvetpr2HVTN+1uUFm7VrLs/\n43FCnCxeo8nF7+vZk3K1SNJlkh4zxsyRtKrT3DJJ/2WMyZCULmm8pNXHerCqqroeREkcRUU5qqio\n8ToGHGJNk0tRUY727j6gdQ++3FGsOtv46Osa9vmPKJSd6UE6xILXaHLxy3oeqwD2pFw9Jek8Y8xi\nSQFJNxpjviyp3Fr7jDHmZ5IWqu28rm9aaxt68FwA4ExTZY3qtu6JOle/vUI1dpsKZoyNcyr0hnBz\ni9b/5AntX7RKrQ1Nyps0QqO+cJWySwd6HQ1JLOZyZa0NSzr62PnaTvO/kfSbWB8fAHpLKDdbaQPy\nVH+ovstcakE/ZQ4r9iAVXItEIlp524+0+7mlHWMH31qvymVrNfuRbylzcKGH6ZDMOGsTQJ+Tkp6q\n4nOmR50rPGOaMory45woduGWVm1+8EWt/Ow9Wvm5H2v7468pEg57HSshVLy6Unv+urzL+CG7VRt/\n+WcPEqGv6Mm2IAD41sTvflqttY3a8/JyNVfWKCU3S8WnT9WUH33O62gnLNzSqhU336U9L/yjY2zn\nnxao4u/vaNrP/0WBQMDDdN7bt3i1Is2tUedq1m6Ncxr0JZQrAH1SMC1V0372RdXv2qcDK8uVO7FU\n2WUlXsc6KVsf/usRxeqwHU8u0KCL56jk4jkepEocoeyMbudSsvhca/QetgUB9GmZJYUquWSO74qV\nJFUufT/6RGtYFX9bGd8wCaj0kxcofVD/rhPBgIrPm9l1HHCEcgUAfnWMbT8+CFVKL8rXuDuuV+aQ\noo6xUE6Wyj51kYZ/4nwPkyHZsS0IAD5VeMZU7XxyQZfxQGqKBl5wqgeJEs+wa87SoAtna9sfXlFL\nXYNKLpmrnLHDvI6FJEe5AgCfGvbRs7R/wTva8dRCKdz23YiB1BQN//j5Kjp7hsfpEkdqTpZG3na5\n1zHQh1CuAMCnAsGgpv3iSxp40WxVvP62gilBDTx/lorOmdHnrxQEvES5AgAfCwQCGnzZPA2+bJ7X\nUQC044xHAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAA\nAIcoVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAO\nUa4AAAAcCnkdAACOFg5HtHzxFq17b6/C4YhKR/TX/HNGKpSa4nU0ADguyhWAhBKJRPTob1dqxZJt\nHWOrVu6SfX+vbvrCHAoWgITHtiCAhLJm1W699Y9tXcbXr6nQwlc3epAIAE4O5QpAQrHv7VU4HH1u\ny8bK+IYBgBhQrgAklGCg+7elYDAQxyQAEBvKFYCEMmXmYKV2c17VmHFFcU4DACePcgUgoYwYPUDz\nzxmp1NQP3p4CQWn6rKGafXqZd8EA4ARxtSCAhHPJVRM1YcpAvfPmToXDEY2dUKSJ00oUCLAtCCDx\nUa4AJKQRYwo1Ykyh1zEA4KSxLQgAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwB\nAAA4RLkCAABwiA8RBQD4UvXaLdr2h1cUam5R2uihGv6J85WSnup1LIByBQDwny2/e0lr7/ydmg8c\n6hjb+eeFOvWhbyqtIMfDZADbggAAn2k5VK/ynzxxRLGSpKpla7Xurj96lAr4AOUKAOAr2x57TfXb\nK6LOVb65Ns5pgK4oVwAAfwmHu52KtHY/B8QL5QoA4CtDrjpD6YP6R50rmDE2zmmQCFrrG7XloZe0\n8VdPq2H3fq/jcEI7AP+LRCLa/MDz2v3sEjXuO6DMYQM1/LpzVXLpXK+joRekFeRoxK2Xav2PHlVr\nXWPHeO6kkRrz5Y96mAxe2PnMIq39/sOq27RLklT+sz9p+CfO17hvfMKzTJQrAL5n//sPKv/5k1JL\nqyTp0Lrtqlq+Rq0NTRp69Rkep0NvGP35jyh/ymjteHKBgk1NSisr0YhbL1NaXj+voyGOGisO6L3/\nuF+Nuys7xpr2V2vDvX9Wv3HDNfQjp3uSi3IFwNdaDtVrx+OvdxSrjvHqOm196CUNuep0BQIBb8Kh\nVxV+aIoKPzRFRUU5qqio8ToOPLD5wZeOKFaHRZpatOf5pZ6VK865AuBrlcvWdHvlWE35NoXrm+Kc\nCEC8tNbUdjvXXN39XG+jXMGphvpmvf/OLu3cftDrKOgjMocUKpiZFnUuNbefgukcoAeSVf4pRurm\nyHS/MUPjnOYDvOvAiUgkouefel8r/7FdByvrFUoNasToAbrqE1NVWMw5EOg9OWa4+s+ZqH2vvdVl\nrvBDUxRISfEgFYB4KLl0rorOnKaKo17/2WOGauRnPuxRKo5cwZE3/rZRr7+0Xgcr6yVJLc1hrV9T\noUd/+5aP7zG4AAAXDElEQVQikYjH6ZDsJv/wNvWfM0FKaXtLC2amaeBFszXxu5/2OBmA3hQIBjXz\ngds14rMfVt60McqZUKqh/3S2Tv2/rytrWLFnuThyBSfeXblTkSif3bdlw36tXb1H4ycPin8o9BnZ\nZYM098//pb2vrtCh8h3qP2s8n3cE9BEpWema+J0bvY5xBMoVnKitaYw6Hg5Le3cf0vjJcQ6EPicQ\nCGjguTM18NyZXkcB0MexLQgn+hdmRR1PTUvRyDED4pwGAADvUK7gxKz5ZcrI7HogdNykYg0rK/Ag\nEQAA3mBbEE5MmTFYrc1hLfn7Ju3dXaOMzFSNnVCsy66Z5HU0AADiinIFZ6bPHqrps4equblVKSlB\nBYN8KjYAoO+JqVwZYzIlPSypWFKNpBustRVH3eduSfPbn+M+a+1vepgVPpGayucKAQD6rljPufqs\npFXW2g9JekjSHZ0njTFnSRptrZ2rtoL1NWMMJ94AAICkF2u5mi/pxfafX5B07lHzSyQd/vS+iKQU\nSc0xPhcAAIBvHHdb0Bhzk6QvHTW8R9LhL4+rkZTXedJa2yCpwRiTKulBtW0LHup5XAAAgMQWiOWr\nSYwxT0r6obV2mTEmT9Iia+2ko+5TIOkJSa9ba793vMdsaWmNhEKcqwMAAHyh26u2Yr1acJGkiyUt\nk3SRpIWdJ9tPeH9V0j3W2t+fyANWVdXFGCWxFBXlqKKixusYcIg1TS7JsJ6tDU3a+vuXVb+9Qtmj\nh2jYR89SMLXvXvydDGuKD/hlPYuKcrqdi/XVeK+kB40xb0hqknSdJBlj7lLb0ap5kkZKusUYc0v7\n79xord0U4/MBACQdXL1Jb3/hJ6p5f0vH2Nbf/VWn3Pdvyho+0MNkAA6LaVuwN1RU1CRGkB7yS+PG\niWNNk4vf13PpNd/WvgXvdBkfdOlczbz/ax4k8p7f1xRH8st6FhXldLstyNffAIBP1G7Zrcpla6LO\nVf7jfbUcqo9zIgDRUK4AwCdaDtUr3NAUda61oUnhRj7xBkgElCsA8IncccOVO7Es6lzepJFK7d/9\nCbYA4odyBQA+EUhJUdnNlyiUm3XEeFphnkbedpkCAb7PE0gEfffaXQDwoeHXnafMkiJte/RVNe49\noIwhAzT8+gs0YNZ4r6MBaEe5AgCfKTprmorOmuZ1DADdYFsQAADAIcoVAACAQ5QrAAAAhyhXAAAA\nDlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA7xxc0AEsKy\nNzbrzcXbVLW/Tv3y0jV5+mCddeEYBQIBr6MBwEmhXAHw3OLXNuqZx1erpTksSaqqrNf2zQdUX9es\nS66a6HE6ADg5bAsC8FQ4HNHyRVs7itVhkYj09vLtamxo8SgZAMSGcgXAUw31zarYeyjqXNX+eu3Y\neiDOiQCgZyhXADyVlh5SVnZa1Ln0zJAKBmTFOREA9AzlCoCnQqGgzMTiqHOjTRHlCoDvcEI7AM9d\n/k+TVV/XrLWr96ihvkWpqSkaZQp1zSeneR0NAE4a5QpAr4lEIlr46ka99/Yu1dU2qbAoW6edPUJj\nxh15pCo1NUWfuPVU7d1do43r9mvI8DwNKyvwKDUA9AzlCkCvefaJ97Tg5XJFIm23d22v1sbySl17\n4wyNnzywy/2LB+WoeFBOnFMCgFuccwWgV9RUN2rl0m0dxeqw2ppGLXptgzehACAOKFcAesWad3er\nprox6tzuHTVxTgMA8UO5AtArCgZkKtDNO0xGJmckAEhevMPBuYMH6rVs0RaFWyKaPL1Eg4fnex0J\nHhg9rkilI/pr84bKLnNjxkf/6AUASAaUKzj1xt826OVnrWprmiRJC17ZoJmnDdcV107mC3j7mEAg\noCuvm6LHH3pb27e0fcp6amqKxk8u1iUfmeBxOgDoPZSrONq+9YBWLtmmpuZWlY3qrxmzhykYTJ7C\nsXd3jf769FrV1TV3jDU2tGjx6xs1rDRfM08b7mE6eGHI8Hx98Rtn6J03d+hAZb1Gjhmg0lH9vY4F\nAL2KchUnr7+0Xq88a9XQ/iW0S/++WW8v26FPfW6WQqkpHqdzY/miLUcUq8MiYen9d3dTrvqoYDCg\n6bOGeh0DAOKGE9rjoGp/nV57cX1HsTps7eo9evWFdR6lcq+5KXyMudY4JgEAwDuUqzhYvnirag81\nRZ3bvH5/nNP0npFjC9XdaVWDh+XFNwwAAB6hXMVBJBzpdq71GHN+M3lGiSZMHdRlfFhZvk4/f7QH\niQAAiD/OuYqDKacM0cJXNnTZFpSk4Un0/WmBQEDX3zZLr724ThvX7VdLS1jDSvN19sVjlZ2d5nU8\nAADignIVByVDczX79DItfHWDwq0fHKkqHVmgsy4e62Ey90KhoM67dJzXMQAA8AzlKk4uu2aSho8o\n0Oq3d6m5qVVDhuXp9PNGKz2DJQAAIJnwN3scTZ05RFNnDvE6BgAA6EWc0A4AAOAQ5QoAAMAhyhUA\nAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiA8RBYAEV/3eZm39/ctqPlCjrBGDNeLW\nS5WW18/rWAC6QbkCgAS27dHX9P53HlBzZU3H2O7nluiUB76mfiMHe5gMQHfYFgSABNXa2Kzynz1x\nRLGSpJo1W7Tu7kc8SgXgeChXAJCg9ry4TLXlO6LOHVi5Ls5pAJwoyhUAAIBDlCsASFADL5yl7DFD\no87ln2LinAbAieKEdgDoZP+S1dr+6GtqqqpRVukgjbj1MmUNLfIkS0p6qkZ/8Sq9/+3/U3Nldcd4\nzoRSjf3qtZ5kAnB8lCsAaLf5wRe19s6H1FJd1zG25+XlmnHfvyl/8ihPMg376FnKmzRSW//wVzVX\nHVL2yMEacetlSs3J8iQPgOOjXAGApNaGJm361dNHFCtJqtu4S+U/eUIz7/+aR8mk3AmlmnTnLZ49\nP4CTwzlXACBpz6srVLtxV9S5g2+tVyQSiXMiAH5FuQIAtZ3f1J1AKCWOSQD4HeUKACQVnzVduZNG\nRJ0rmDVegUAgzokA+BXlCgAkBVJSNParH1P64AFHjOdNG61x3/iER6kA+BEntANAu0EXzlL+9NHa\n/MDzaqqqUb8xQ1V6/QVKyUjzOhoAH6FcAUAnGQP7a9zXOVIFIHYxlStjTKakhyUVS6qRdIO1tiLK\n/bIkLZZ0u7X2xZ4EdSEcjuhAVb0yMkLKyuZfogAAwL1Yj1x9VtIqa+13jDHXSrpD0r9Eud8vJCXE\n9cv/eGOzFr+2Sbt3VCs9I6SRplBXXjtFeQWZXkcDAABJJNYT2udLOnwk6gVJ5x59B2PMv6ntqNU7\nMT6HM6tW7tDTj6zSjq0H1doaUV1ts1av3KXf/+ZNPrsGAAA4ddwjV8aYmyR96ajhPZIOtv9cIynv\nqN85R9IYa+1txph5JxKkoCBLoV76LJl3V+xSU2Nrl/FN5fu1deMBzZwz3OnzFRXlOH08eI81TS6s\nZ/JhTZOL39fzuOXKWnu/pPs7jxljnpR0+E+eI+nAUb92k6RSY8zrksZJmmGM2W2tfbu756mqqutu\nqsf27qmJOh6JSOvW7FHpqAJnz1VUlKOKiujPB39iTZML65l8WNPk4pf1PFYBjPWcq0WSLpa0TNJF\nkhZ2nrTWXnf4Z2PMbyU9cqxi1dvy8jO1Y8vBLuOBgDSwJNeDRAAAIFnFes7VvZImGmPekHSrpO9K\nkjHmLmPMLFfhXDll7jClpXfdciwbNUCTppd4kAgAACSrmI5cWWvrJF0TZfzfo4x9KpbncGnqKUNU\nV9OkJQs2a9f2g8rISNVIM0BXXjtFwSBfaQEAANzpMx8iOvfMEZp9epkq99UqIzNV/XLSvY4EAACS\nUJ8pV5IUDAZUWNzP6xgAACCJ8cXNAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHK\nFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQLga42NLaqpblQkEvE6CgBI6mNf3AwgeRyq\nadTTj67ShrUVamxo0cDBOTrtrJGaOXe419EA9HGUKwC+E4lE9Pv7lmv92n0dY1s3HdDe3e8qIzNV\nk6aVeJgOQF/HtiAA31m7eo82rNvXZbyhvkXL39jiQSIA+ADlCoDv7Nh6UOFw9LkDlfXxDQMAR2Fb\nEIDvFA/K6XYuJy8jjkmAriKtrdr2+OuqWr5WKZnpGnr1mcqfNtrrWIgjyhUA35k0vUSlowq0ZUPV\nEeOh1KBmzB7qUSpAam1o0vJP/UD7XnurY2z7H1/R6H+9RqO/cJWHyRBPbAsC8J1gMKCP3zRT4ycP\nVHp6iiSpuCRHF10xQTPmDPM4Hfqy8p89cUSxkqSWQw3a8Is/q27rHo9SId44cgXAl/oXZeumL85V\nVWWdamuaNGhIrkIh/r0Ib1X+Y23U8eaqGm1//HWN/co/xTkRvEC5AuBrBf2zVNA/y+sYgCQp0t2V\nFpIird3PIbnwzzwAABzp7sT1UE6WSj48L85p4BXKFQAAjoz+0jUqOHXcEWOB1BQNv/585Rq+PaCv\nYFsQAABH0nKzNeex72rjb/6i6lUbFcxI06CL56jk4jleR0McUa4AAHAoJStdY/7laq9jwENsCwIA\nADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABw\niHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDl\nCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMChkNcB\nACCZtLaG9Y+FW7RhXYUCCmjUuELNnl+mYDDgdTQAcUK5ApBU9i9Zre2Pva7mg4eUPXKwRn7mw0ov\nzIvLc4fDEf3u18u1+q1dHWNvL9+h9Wsq9IlbTqVgAX0E5QpA0th0/3OyP3hYLTX1HWN7/rpcMx+4\nXf1GD+n151++aMsRxeqwd9/cqbembtcpc4b1egYA3uOcKwBJoaW2XhvvffqIYiVJh+w2rf/J43HJ\nsHHd/m7nytdUxCUDAO9RrgAkhZ1Pv6H6bXujzh1YuS4uGQLH2PUL8G4L9BlsCwInqamxRYtf36Rd\n26uVnhnSqacN17CyAq9j9XnB1NRu5wIpKXHJMG7yQK1Yuk2RyJHjwaA0YUpJXDIA8B7lCjgJh2oa\n9cDPl2rrpqqOsZVLt+nCK8Zr/tmjPEyGksvnaf3/PKbaDTu7zBXMNHHJMHXmEK17f6/eXLxV4XDb\nWDAloFnzSjVx2qC4ZADgPcoVcBJeenrNEcVKkhrqW/TaS+WaOXe4MjK7P3qC3pWSnqqx/3at3vvW\n/WqqONgxnjd9jMzXPx6XDIFAQNd8cromTR+sNe/uliRNmFqicZOKFTjWniGApEK5Ak7Clo1VUccP\nVtZrxZKtmsfRK08N+cjpyp8xRlseekktB2vVb+wwld5woVIy0uKWIRAIaMKUQZowhSNVQF8VU7ky\nxmRKelhSsaQaSTdYayuOus+nJH1WUoqkp6213+tZVCCxRY5/F8RBdlmJJnzrU17HANCHxXr9ymcl\nrbLWfkjSQ5Lu6DxpjBnVfp8zJc2SlGaMYb8Evjd8RPQT1/MKMjRjNp9hBACIvVzNl/Ri+88vSDr3\nqPlzJb0p6UFJf5e0yFrbHONzAQnjvMuMhpbmHzGWnhHSGeeNVlZ2/LaeAACJKxA5+prhoxhjbpL0\npaOG90j6Z2vtGmNMUNJWa+3QTr/zTUnXSjpNUqakNyTNstYe6O55WlpaI6FQfC6XBnqiob5ZLz+3\nVtu3HFBmVqrmnTlSY8YXex0LABBf3V6lctxzrqy190u6v/OYMeZJSTntN3MkHV2a9kt63VpbI6nG\nGLNG0lhJy7p7nqqquuNF8YWiohxVVNR4HQMORVvTOWeUHXGbNfcPXqPJhzVNLn5Zz6KinG7nYr1a\ncJGki9VWli6StDDK/OeNMRlqO6F9gqTyGJ8LAADAN2ItV/dKetAY84akJknXSZIx5i5JT1hrlxlj\n7ldbyQpI+p61ttJFYAAAgER23HOu4qWioiYxgvSQXw5n4sSxpsmF9Uw+rGly8ct6FhXldHvOFV8l\nCgAA4BDlCgAAwCHKFQAAgEOUKwAAAIf44mb4SiQS0ZLXN2n127tUV9ukooE5mn/2SJWO6u91NAAA\nJFGu4DPPPvGeFrxSrki47fb2LQe1Yd0+XXfzTI02hd6GAwBAbAvCR6oPNmjF0q0dxapj/ECDFr7C\nZ9QCABID5Qq+sfrtXTpU3RR1bveOxP9MFABA30C5gm/k5mV0O5eWzg43ACAxUK7gGxOmDNKw0vyo\nc2PGc74VACAxUK7gG8FgQJdfO1mDhuR2jIVSg5o4bZAuvnKCh8kAAPgAeynwlRGjB+hf7zhTK5Zs\nVfXBRo0cM0CjuEoQAJBAKFfwnVAoqNkfKvM6BgAAUbEtCAAA4BDlCgAAwCHKFQAAgEOUKwAAAIco\nVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4A\nAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHApEIhGv\nMwAAACQNjlwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAh0JeB0gWxpgrJV1jrb0u\nytwtkm6T1CLpTmvts/HOhxNjjMmU9LCkYkk1km6w1lYcdZ+nJRVKapZUb629KO5BcVzGmKCkX0qa\nKqlR0s3W2vJO87wufeQE1vOnkuar7XUrSR+21h6Me1CcFGPMbEn/ba0986jxyyR9S22vzwestb/x\nIF7MOHLlQPuL+geK8v+nMWaQpC9KmifpAkk/MMakxzchTsJnJa2y1n5I0kOS7ohynzGS5ltrz6RY\nJbQrJGVYa+dKul3SPYcneF36Urfr2e4USRe0vy7PpFglPmPMv0v6X0kZR42nSvofSedLOkPSrcaY\ngfFPGDvKlRuL1faXcjSzJC2y1ja2v9jLJU2JWzKcrPmSXmz/+QVJ53aebH+B50v6izHmDWPMpXHO\nhxPXsZbW2qWSZnaa43XpP92uZ/tRrTGS7jPGLDLGfNqbiDhJGyR9JMr4eEnl1toqa22TpDcknR7X\nZD3EtuBJMMbcJOlLRw3faK191BhzZje/liup87+gaiTl9UI8nKRu1nOPPlivaGuVprZ/Mf9UUn9J\ni4wxy6y1e3szK2Jy9Guv1RgTsta2RJnjdZn4jrWe2ZJ+LunHklIkvWaMedNa+64HOXGCrLV/MsaU\nRZny/euTcnUSrLX3S7r/JH+tWlJOp9s5kg44C4WYRVtPY8yT+mC9oq3Vbkm/an9D32uMeUuSkUS5\nSjxHv/aC7esWbY7XZeI71nrWSfqptbZOkowxf1PbuVmUK3/y/euTbcHet0zSh4wxGcaYPLUd7lzt\ncSZ0b5Gki9t/vkjSwqPmz5X0uCQZY/pJmiRpTdzS4WR0rKUxZo6kVZ3meF36z7HWc6zajiKntJ+v\nM1/SyvhHhCNrJI0xxvQ3xqSpbUtwiceZTgpHrnqJMebLatszfsYY8zO1/SUdlPRNa22Dt+lwDPdK\netAY84akJknXSZIx5i5JT1hrXzDGXGCMWSopLOkb1tp93sXFMTwl6TxjzGJJAUk38rr0teOt5+8k\nLVXbVbwPWWvf8zArYmCMuU5SP2vtfe1r+5LaXp8PWGt3eJvu5AQikYjXGQAAAJIG24IAAAAOUa4A\nAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAh/4/6GXVz4e4aqUAAAAASUVO\nRK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x113f303c8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,10))\n",
    "plt.scatter(X[:,best_features_scoring[0]], X[:,best_features_scoring[1]], c=t, cmap=cm.Spectral)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 98,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([196, 200,  22, 142, 280, 136,  29, 250, 223, 162, 193,  68, 203,\n",
       "        83, 271, 278,  82, 295,  18, 256,  16, 124,  70, 176, 218, 106,\n",
       "        48, 231,  54,  41, 224,  71, 119,  10,  47, 253, 121,  12,  67,\n",
       "       144, 246, 104, 166,  33, 157, 174, 254, 125, 198, 158,  80, 285,\n",
       "       190,  26,  98,  39, 227, 248,  20, 143, 141,  21,  30, 237, 214,\n",
       "       135,  17, 262,   7,  96, 178, 270, 115,  77, 129, 257,  45, 187,\n",
       "       131, 184,  55,  31, 123, 260,  85,  90, 282, 139, 132, 137,  11,\n",
       "       168,   8, 105,  32,   1, 299,  36, 274,  52, 258, 284, 242, 163,\n",
       "       116,  28, 213,  64, 281, 221, 122, 251,  84, 273, 111, 127, 240,\n",
       "       108, 206, 244,  72, 171,  60, 112, 265, 298, 180, 199, 110, 255,\n",
       "        65, 117,  97,  74, 294, 147, 287, 272, 165, 177, 153, 266, 148,\n",
       "        44,  40,  75, 149, 118,  94, 150, 235, 230,  23, 261, 283, 155,\n",
       "       263,  62, 145,  79,  35, 233, 211,  66, 297, 202, 102,  89, 159,\n",
       "       207, 293, 183, 243,  76,  51,  58, 170, 188,  87, 107, 114, 126,\n",
       "        53, 291, 217, 182, 152, 173, 269, 130, 189, 275, 232, 181, 279,\n",
       "        34, 245, 222, 201, 185,  63, 204,  88,  38, 276, 169, 128, 219,\n",
       "       264,   9, 292, 209, 138, 154,  92, 100,  69, 205,  81, 216,  14,\n",
       "       249, 228,  61, 267, 247,  46, 133, 241, 186,  56, 286, 101,  86,\n",
       "       194, 195, 134,   3, 113,  59,  95, 151, 192, 239, 164, 156,  57,\n",
       "        99,  37, 191, 289, 146,  43,  73,  93,  42,  25,  24, 120, 103,\n",
       "        78, 288,   6, 220, 290, 208, 210, 226, 175, 161, 109, 234, 277,\n",
       "       172,  27, 259, 252, 197, 160,   2, 238,   5, 140, 268, 225,   0,\n",
       "       212, 167,  15,  50,  91,  49, 229, 215, 179, 236,  19,  13,   4, 296])"
      ]
     },
     "execution_count": 98,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "best_features_scoring"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h4>2.1.2 Clustering</h4>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can use the K-means code you wrote last week. Remember to pass the data the correct way around and only cluster the training data. Transforming the test data is a little more complex –your kmeans should give you assignments of features to clusters. Use these alongside the test data to compute the cluster means for the test data (you do not need to re-assign the clusterings)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def euclidean_distance(p1, p2):\n",
    "    return np.linalg.norm(np.array(p1)-np.array(p2))\n",
    "\n",
    "def find_closest_centroid(datapoint, centroids):\n",
    "    # find the index of the closest centroid of the given data point.\n",
    "    distances = [euclidean_distance(datapoint, _) for _ in centroids]\n",
    "    return np.argmin(distances)\n",
    "\n",
    "def randomize_centroids(data, k):\n",
    "    centroids = data\n",
    "    np.random.shuffle(centroids)\n",
    "    return centroids[:k]\n",
    "\n",
    "MAX_ITERATIONS = 500\n",
    "\n",
    "def check_converge(centroids, old_centroids, num_iterations, threshold=0):\n",
    "    # if it reaches an iteration budget\n",
    "    if num_iterations > MAX_ITERATIONS:\n",
    "        return True\n",
    "    # check if the centroids don't move (or very slightly)\n",
    "    dist = [euclidean_distance(curr, old) for curr, old in zip(centroids, old_centroids)]\n",
    "    if all(dist[x]<=threshold for x in range(len(dist))):\n",
    "        return True\n",
    "    return False\n",
    "\n",
    "def update_centroids(centroids, cluster):\n",
    "    for i, c in enumerate(cluster):\n",
    "        centroids[i] = sum(c)/len(c)\n",
    "    return centroids"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 100,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# data : set of data points\n",
    "# k : number of clusters\n",
    "# centroids: initial list of centroids\n",
    "def kmeans(data, k=2, centroids=None):\n",
    "    \n",
    "    # randomize the centroids if they are not given\n",
    "    if not centroids:\n",
    "        centroids = randomize_centroids(data, k)\n",
    "\n",
    "    old_centroids = centroids[:]\n",
    "\n",
    "    iterations = 0\n",
    "    while True:\n",
    "        iterations += 1\n",
    "\n",
    "        # init empty clusters\n",
    "        clusters = [[] for i in range(k)]\n",
    "\n",
    "        # assign each data point to the closest centroid\n",
    "        for point in data:\n",
    "            # find the closest center of each data point\n",
    "            centroid_idx = find_closest_centroid(point, centroids)\n",
    "            \n",
    "            # assign datapoint to the closest cluster\n",
    "            clusters[centroid_idx].append(point)\n",
    "        \n",
    "        # keep the current position of centroids before changing them\n",
    "        old_centroids = centroids[:]\n",
    "        \n",
    "        # update centroids\n",
    "        centroids = update_centroids(clusters, clusters)\n",
    "        \n",
    "        # if the stop criteria are met, stop the algorithm\n",
    "        if check_converge(centroids, old_centroids, iterations):\n",
    "            break\n",
    "    return np.array(centroids), iterations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "new_data = kmeans(data=X.T)[0].T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlcAAAI+CAYAAACRwC8WAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VPW9//H3JJOdbCRhDzt8QUBEEZBFpW61irVetVZ7\ntbZaq73dN3uvt8uvdr92vWoXtWpta70udamC1hUBi4oIqHwh7DshZCN7JvP7IyEGMgl08p05mTmv\n5+PRxyNzvpM5709C4zvnnDkJhMNhAQAAwI0UrwMAAAAkE8oVAACAQ5QrAAAAhyhXAAAADlGuAAAA\nHKJcAQAAOBT0OsBh5eW1SX9PiMLCbFVW1nsdwzN+nt/Ps0vM7+f5/Ty75O/5k332kpLcQE9rHLmK\no2Aw1esInvLz/H6eXWJ+P8/v59klf8/v59kpVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ\n5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoV\nAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4FDQ6wD9UaipRbsefkmN+ytV\nfPp0DTzFeB0JAAAkCMrVUSqWr9Pam3+nQ3a7JKnsFw9r8HmzNOP2LyoljS8XAADoHacFu2hrDemd\n/767s1hJUltjs/Y8/qo2/PQvHiYDAACJgnLVxZ6nlqtm3ZaIa+WvrIlzGgAAkIgoV1007a/sca31\nUEMckwAAgERFuepiyAdnK60wN+Ja7qSRcU4DAAASEeWqi+yRgzX8IwukQOCI7ZnDijT2hkUepQIA\nAImEt78dZcoPrlfWqMHa/9wbaq2tV874ERp7/QUqmDHR62gAACABUK6OEggENO4zH9a4z3zY6ygA\nACABcVoQAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4FBUNxE1\nxqRIukPSdElNkq6z1pZ1WT9f0rclBSS9Kemz1tpw3+MCAAD0b9EeubpYUqa19jRJN0u67fCCMSZX\n0k8lXWitnS1pq6TiPuYEAABICNGWq/mSFkuStfY1STO7rM2VtFbSbcaYpZL2WWvL+5QSAAAgQUT7\ntwXzJFV3eRwyxgStta1qP0q1UNJJkg5JWmqMWWGt3dDbCxYWZisYTI0yTuIoKcn1OoKn/Dy/n2eX\nmN/P8/t5dsnf8/t19mjLVY2krl+xlI5iJUkVkl631u6VJGPMK2ovWr2Wq8rK+iijJI6SklyVl9d6\nHcMzfp7fz7NLzO/n+f08u+Tv+ZN99t6KY7SnBZdJ+pAkGWPmqP004GGrJE01xhQbY4KS5kh6N8r9\nAAAAJJRoj1w9JukcY8xytb8j8FpjzJcllVlrnzDGfFPSko7nPmStXecgKwAAQL8XVbmy1rZJ+sxR\nm9d3WX9Q0oN9yJU0wuGwatZtUXNlrQovmOV1HAAAEGPRHrnCcahas0nvfvseVb5uFW5p1XsTRmj4\nx8/VuM9c5HU0AAAQI5SrGAk1tejtL/xate9u7dxWs3GnDv3oAWUNL9KwRfO8CwcAAGKGP38TIzsf\nfP6IYnVYW0Ozdj+6NP6BAABAXFCuYqRxT0WPa00HqntcAwAAiY1yFSO5k0dJgUDEtazSQXFOAwAA\n4oVyFSNDF83VwDkndNueXpKvUVef50EiAAAQD5SrGAmkpOiUu76m4ZedqawRJUobmKuhH5ih6T//\nnIoilC4AAJAceLdgDGUUF2jG/35RbS2tCreENGRUcVL/KQAAAEC5iouUtKCUxpcaAAA/8OVpwXA4\nrLbWkNcxAABAEvLV4ZS2llat//4ftf/5VWqpqlXO2GEqveoclV6+0OtoAAAgSfiqXK35yu3a+dcX\nOx837a9S9drNCqSkaMSlZ3iYDAAAJAvfnBas27pHexev7LY9VNeoHX9+zoNEAAAgGfmmXB14da1a\nq+sirtVv2xfnNAAAIFn5plzlTRqpQEZaxLX0orw4pwEAAMnKN+WqcOYkFc2OfPPOwefNinMaAACQ\nrHxTriRp+s8/q5KFM5SSlS5JyhhcqNHXXaAJX7rM42QAACBZ+OrdglkjBmn2g99W9TtbVbd5l4pO\nm6qM4nyvYwEAgCTiq3J1WP6U0cqfMtrrGAAAIAn56rQgAABArFGuAAAAHKJcAQAAOES5AgAAcIhy\nBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADjky78tmGzKl67Rzv97US0V\ntcoePVhjPr1IOaOGeB0LAABfolwluK33PK313/+jWg81dG7b/483dfJvvqyCGRM9TAYAgD9xWjCB\nhRqatPl3Tx5RrCSpfutelf3qUY9SAQDgb5SrBLbvhVWq37In4lr122UKh8NxTgQAAChXCSyYldHj\nWiCNM74AAHiBcpXASs6YrrxpYyOuFZ46SYFAIM6JAAAA5SqBBVJTZW6+UlnDS47YXnDyRE265WqP\nUgEA4G+cO0pwg8+eqYJnJ2jrH55Wc0Wtcs0IjbzqHKWkp3kdDQAAX6JcJYGM4nyZr33M6xgAAECU\nq7gIt7Vp33Nv6EDNIWXPPVHZw4u9jpTw9jz9mnY+9IIadx9U5rAilV5xloZ8cJbXsQAAoFzFWtXq\njVp78+9UvbpMCoeVNjBXwy5eoKk/uJ4LzqO07Y9L9O63/6BQXaOk9ttOVLy6Vif8v09q5JVne5wO\nAOB3vr6gvWpNmexPH1TZ/z6q5spa568fDoW09mu/UfVbG6WOe061HKzVtj88o013Pu58f34QDoW0\n7d7FncXqsNbaem27d7HCbW0eJQMAoJ0vj1yFw2Gt+eod2v3IKwo1NEmSttz1d036z4+r9PKFzvaz\n+4llql6zKVIA7V/yusbfdLGzfflF/fb9qnl3a8S1mve2qmHXAWWXDopvKAAAuvDlkaut9zytHQ88\n11msJKlpT4XWf/+PajpQ7Ww/DbsqelxrqT7kbD9+EszNVnBAduS1nCwFcyOvAQAQL74sV/ufXxVx\ne9Peg9r2xyXO9lN8xnSlZmdGXMsZN8zZfvwkozhfRXOnRFwrmjtV6QUD4pwIAIAj+bJchRqbelxr\nq+957V9VMG2sBkd4B1tGSYFGf/JDzvbjN1NuvU6FsyZJh98QkBJQ4ezJmnLrp7wNBgCAfHrNVd7k\nUTq4bF237YH0NBWfeZLTfZ30q88ra3ixDry8Wm31TcoaO0xjrr9AxfOmOd2Pn2SXDtLcx3+gPU+t\n0KGNOzVg4ggNveA0BVJ8+bsCAKCf8WW5GvfZj+jAq2t1aP32I7YPvWCO89KTkhbU5Fuulm65WiUl\nuSovd/+uRD8KpKRo2EXzvI4BAEA3vixXWcOKNftPt6jsfx9TzTtblJKRrpLTp2ss794DAAB95Mty\nJUlZIwZp2o9u8DoGAABIMlykAgAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEA\nADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcCkbzScaYFEl3SJou\nqUnSddbasgjP+bukx621v+lrUAAAgEQQ7ZGriyVlWmtPk3SzpNsiPOdWSYXRBgMAAEhE0Zar+ZIW\nS5K19jVJM7suGmMuldR2+DkAAAB+EdVpQUl5kqq7PA4ZY4LW2lZjzFRJV0q6VNK3jvcFCwuzFQym\nRhkncZSU5HodwVN+nt/Ps0vM7+f5/Ty75O/5/Tp7tOWqRlLXr1iKtba14+OrJQ2X9IKk0ZKajTFb\nrbW9HsWqrKyPMkriKCnJVXl5rdcxPOPn+f08u8T8fp7fz7NL/p4/2WfvrThGW66WSVok6SFjzBxJ\naw8vWGu/fvhjY8x3JO09VrECAABIFtGWq8cknWOMWS4pIOlaY8yXJZVZa59wlg4AACDBRFWurLVt\nkj5z1Ob1EZ73nWheHwAAIFFxE1EAAACHKFcAAAAOUa4AAAAcivaCdvQTu59cpt1/W6aWqlrljBmq\n0dddqLxJI72OBQCAb1GuEtjGXzysjT97SG1NzZKkilfXqvyl1Tr5t19R4SnG43QAAPgTpwUTVEtt\nvbbdv7izWB3WsGO/Nt3+N49SAQAAylWC2vv3FWrcdSDiWvXaTXFOAwAADuO0YIKoWPGOttz1lA5t\n3Km0vGxllg7u8bmpmelxTAYAALqiXCWAihXvaNUN/6OmfZXvb3zdKq1wgFoqD3V7/sDZJ8QxHQAA\n6IrTgglgy11PHVmsOrQ2tihzePH7GwIBFc2fpknfuiaO6QAAQFccuUoAhzbujLg93NCkkddcppTU\nFDUfrFX+tDEa9uH5CqTQmQEA8ArlKgGk5WX3uJYzarCGX7wgjmkAAEBvOMSRAEoWnhxxe960sRq2\naG6c0wAAgN5QrhLAhC9eqhEfO0vBLkew8qaN1dQffVqB1FQPkwEAgKNxWjABBFJTddIvPqfamy7W\n/udXKXPIQA27aC7FCgCAfohylUByJ5Yqd2Kp1zEAAEAvOC0IAADgEOUKAADAIcoVAACAQ5QrAAAA\nhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5R\nrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhylU/EW5r8zoCAABwIOh1AD8LNbVo/a33q/zl\n1WqpqdOA8SM06hPna9iFp3kdDQAARIly5aG3P/9L7f7bq52Pm/YcVPWazUrNCGrwOad6mAwAAESL\n04IeqVq7Wfuee6Pb9tbqQ9p2/xIPEgEAABcoVx45uGytQnWNEdfqt+6NcxoAAOAK5cojOWOGSSmB\niGtphXlxTgMAAFzhmiuPDDp3pgpPMap8ff2RCykBDfngLG9CHUNrXaO2P/CsmitqVDBzogafc6oC\ngcgFEQAAv6JcORQOh1X99iaFW0MqmDFegdTUHp8bCAQ0/Wf/oTXf+I0q31ivcHOrMgYVatglp2vs\njR+OY+rjU770ba39xm9Vv2l3+4bUFJWcMV0z775ZqdkZ3oYDAKAfoVw5sv/Ft7ThJ39R1dtlUltY\neVNHa9xnL9Hwjyzo8XMGTByh0x79nqretKrbtk8lC2coY2D/OyUYDoX03nfufb9YSVKoTeUvvKX3\nbr1fU39wvXfhAADoZyhXDjTurdCar96uxp0HOrfVrN2id/7r9xowfrjyp43t8XMDgYAKZ05S4cxJ\n8YgalT3P/FM167ZEXKtYvi7OaQAA6N+4oN2BLXc/fUSxOqy5okbb//isB4ncaj5Y0+NaqLE5jkkA\nAOj/KFcONFf0XD6aDlTHMUlsDL1wrjIGFUZcy5syOr5hAADo5yhXDmSXDupxLWtESRyTxEbGwDyV\nXnW2AmlHXqCfWVqicTdd7FEqAAD6J665cmDM9Rdq19+W6tD67Udszx45SGOuv9CjVG5Nuvkq5YwZ\npr1Pr1BL1SFljxmqsdcv4sgVAABHoVw5EByQpZN/91XZH/5JVW9atYXaVHjSeI3/0uW9HtVKNKUf\nXajSjy70OgYAAP0a5cqRPDNSp977TYUamhRuCyuYk+l1JAAA4AHKlWOpWdxQEwAAP+OCdgAAAIco\nVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIctWDQ5t2af8Lq9RSU+d1FAAAkEC4FcNRGvZUaM1X79DB\n5e8oVN+ozGHFGnbxPE3+1icUCAS8jgcAAPo5jlwd5e0v/Erl/3hTofpGSVLj7gPafOcT2nT7ox4n\nAwAAiYBy1cXB19er4rV3uy+Ew9r79D/jHwgAACQcTgt2UfPOVoWbWiKuNZVX9fh5tRt2as+Ty5SS\nkabSK89WxsC8WEUEAAD9HOWqi6K5UxXMzVZrbX23teyRg7ttC4fDevdb92jHg8+rtab9c7b8/ilN\n/MpHNerq82KeFwAA9D+cFuwid+IIDTpnZrftqdkZGnHFB7pt3/XIy9py9987i5UkNe09KPujP6lu\n276YZgUAAP0T5eoo03/xOY365IeUM3aY0oryVDBzoqZ8/zqVXraw23P3LV4phdq6bW+uqNGOPz0b\nj7gAAKCfieq0oDEmRdIdkqZLapJ0nbW2rMv6lyRd0fHwaWvtd/saNF5SM9I07YefVjgUUqixRanZ\nGT3egqG1oanH12mt73kNAAAkr2iPXF0sKdNae5qkmyXddnjBGDNW0lWS5kqaI+lcY8yJfQ0ab4HU\nVAVzMnu9t1WeKY28kBLQwNmTY5QMAAD0Z9GWq/mSFkuStfY1SV0vVNoh6YPW2pC1NiwpTVJjn1L2\nU2Nvulh5U8d22z7o7FM09MK5HiQCAABei/bdgnmSqrs8DhljgtbaVmtti6QDxpiApJ9Kestau6Gv\nQfujjOICnfrnW7TpFw+ras0mpaQHVXTaFE344mXczR0AAJ+KtlzVSMrt8jjFWtt6+IExJlPSPZJq\nJd10PC9YWJitYDA1yjgeKsnVyLu+cvxPL8k99pOSmJ/n9/PsEvP7eX4/zy75e36/zh5tuVomaZGk\nh4wxcyStPbzQccTqcUkvWGt/fLwvWFnZ/d5SyaakJFfl5bVex/CMn+f38+wS8/t5fj/PLvl7/mSf\nvbfiGG25ekzSOcaY5ZICkq41xnxZUpmkVElnSMowxpzf8fxvWmtXRLkvAACAhBFVubLWtkn6zFGb\n13f5ODPqRPCFcDisbfct1v7nVynU0KhcM1LjPvsRZQ0r9joaAAB9wp+/gSfWfeO32nb/EikcliRV\nLF2rA6+u1awHblF26SCP0wEAED3u0I64q1m3RbseebmzWB12aP12bfrfRz1KBQCAG5QrxN3eJSvV\neqgh4lr1ui1xTgMAgFucFuxnGvYd1IYf/VmVb1iF29pUMGOCJnzloxowZqjX0ZxJzcroeS0zPY5J\nAABwj3LVj4QamvTGNT9U9VsbO7fVle1S9drNmvvYrUofmOdhOndKrzxbW37/lBp3H+i2VjRvmgeJ\nAABwh9OC/ci2+5YcUawOO7R+uzb/9gkPEsVGesEAmW9epYzBhZ3bAulBDblorsZ//t88TAYAQN9x\n5KofqV2/vce1us174pgk9kovX6iSM0/S9geeU6i+QcXzp6v4jOn82SAAQMKjXPUjaQU5Pa4F87Lj\nmCQ+MgcVauKXL/c6BgAATnFasB8p/fi5Si/qfl1VMDdLpZcv9CARAAD4V3Hkqh/JHT9ck797rcp+\n9lDnacCs0kEae+OHNXD2CVG9ZnNNnbb89kkdsjsUzM3S8EtOV/GCE13GBgAAXVCu+pnSyxZq2KJ5\n2v34qwqHQhr24QUK5kT314Qa9lTo9X//vmrWbu7ctvtvSzXhyx/V+M9d4ioyAADognLVD6Vmpqv0\nox/o8+tsvO2vRxQrSQrVN2nL755U6VVnKyNJbu0AAEB/wjVXSawqwm0dJKlpf6V2PfxynNMAAOAP\nlKskFkjp+dsbSOVbDwBALPBf2CRWeKqJuD1zRLFG8O5DAABignKVxCZ+/WMqnDXpiG3B/ByN//yl\nSstNvvtmAQDQH3BBexJLL8jVnIe/p+0PPKuadVsUHJCl0ivPVt7kUV5HAwAgaVGuklxqRprGfOoC\nr2MAAOAbnBYEAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5\nAgAAcIhylSRaqg/p0ObdCjW1eB0FAABf48/fJLjWugat/fpvdOCVNWoqr1LO2GEafskCTfjKRxUI\nBLyOBwCA71CuEtzqL/xae59c3vm4btMubfjZQ0rNztS4my72MBkAAP7EacEEVlu2SwdeWt19IdSm\n3Y+/Gv9AAACAcpXIqleXqbW2PuJa496DCre1xTkRAADgtKBDe595Tdv/+Kzqtu9XRlGehpw/R2Nu\nWBSza58KT5mg1LxshWq6F6zMoUUKpNCdAQCIN8qVI7v+tlRrv3anWjuKTt1G6eDK99RUXqXJ/311\nTPaZM2aYBp05Q3ueWHbkQjBVwy+eH5N9AgCA3nFow5Ft9y3uLFad2sLa+ejLaqmpi9l+T/rl5zXi\nig8oc8hABYKpGjCxVOZrH9OYGy6K2T4BAEDPOHLlQFtzi+rKdkVca9pdoQNL12roBXNisu/U7Ayd\n9MvPq7WuQc0VNcocWqSUNL6tAAB4hf8KOxBICyqYm62m/VXd1lIy0pQ9clDMMwRzshTMyYr5fgAA\nQO8oVw4EAgEVnzFddZt2d1sbOPsE5U8b2+d9hMNhbfnNE9r91HI17atU1ogSjbj0TI38+Dl9fu3j\n2XdzRY1SMtKUlpsd8/0BAJDIKFeOnPDta9VcXq39z69SqL5RSk1R4UyjqT++wcnrb/jJX7Txlw9L\nofbbKzTs2K/KVRsUamrWmE9d4GQfkexdvFKbf/O4qtdtUWpGmgbOmqzJ3/mEckYNidk+AQBIZJQr\nR1Iz03XKXV9X1dtlOrjiHeWMG65BZ5/i5DYMofom7Xr0lc5idVi4qUU7/vK8Rl97fkxuu3DwDas1\nX71dzeXV7Tlqpb1Pv6aG3Qc076kfcW0XAAAR8F9Hxwqmj1fB9PFOX7Nmww7Vb90bca1uyx41V9Qo\no6TA6T4lafv9SzqLVVfVq8u048HnNerfz3O+TwAAEh23YkgAmUMHKpiXE3EtvTBXwRhdB9Wwu7zH\ntUjXlwEAAMpVQsgaPFDFC06MuFZ8+olKzUyPyX4zSwp7zlQa+3dAAgCQiChXCWLaT25QyVknK6Wj\nSAUHZGrIhadp6q3Xx2yfIz52loL53Y+Y5U4epZFXxf5digAAJCKuuUoQGcUFmv3nb6nyrQ2qXrNZ\nA2dNVt7kUTHdZ8np03XCdz+prXf/XTVrNyslM12Fp07SCd+5NmZHywAASHSUqwRTOGOiCmdMjNv+\nRn7sLJVefqZq7A4FB2QpZ+TguO0bAIBERLnCMQVSU5V/wmivYwAAkBC45goAAMAhyhUAAIBDlCsA\nAACHKFcAAAAOcUE7YqrpQLV2P7ZUKdnpKrzhQq/jAAAQc5QrxMyG2x7UtnuXqGl/pSRp6x2Pa/zX\nrtDwixd4nAwAgNjhtCBiYs/TK7Txl490FitJqi3bpXf/+2417q3wMBkAALFFuUJM7H5imcJNLd22\nN+2v0rZ7F3uQCACA+KBcISZaaxt6XGuprY9jEgAA4otyhZgYMG54j2v508fHMQkAAPFFuUJMjP3M\nRRowsbTb9qL50zTi3073IBEAAPHBuwURE1nDijXzDzer7FePqPrtMgXSgxp+xnSN/OJlCqSmeh0P\nAICY8WW5CofDCgQCXsdIegPGD9dJv/p85+OSklyVl9d6mAgAgNjzTblqa27Re99/QOUvvqWW6kMa\nMH64Rl3zQQ27aJ7X0QAAQBLxTbla/cVfa/cjr3Q+btp7UNVrNimQmqqhF8zxMBkAAEgmvrigvWb9\nNu1f8nq37a019dr+wBIPEgEAgGTliyNXFUvXqvVQ5Psu1W/dG+c0gBt1W3Zr+5+fV1tjs4rmT9Pg\nc0/lWkIA6AeiLlfGmBRJd0iaLqlJ0nXW2rIu69dLukFSq6RbrbVP9TFr1HLGDZeCqVJrqNta+sA8\nDxIBfbPl7r9rw08fVEtlbefjIefP1sm//apSgrwbEwC81JfTghdLyrTWnibpZkm3HV4wxgyR9HlJ\n8ySdJ+mHxpiMvgTti5KFJ2ngTNN9IRDQ4PNmxT8Q0AcNew5ow8/+2lmsJEmhNu19aoU23/6Yd8EA\nAJL6Vq7mS1osSdba1yTN7LI2S9Iya22TtbZaUpmkE/uwrz4JBAI68WefVdH8aQqkp0mSMgYVasz1\nF2rc5y7xKhYQlR1/fl4tB2oirh1Yti7OaQAAR+vLNVd5kqq7PA4ZY4LW2tYIa7WS8vuwrz4bMG64\nTnvke6pctVF1W3erZMF0ZZQUeBkJiEpbS2sva93/WDYAIL76Uq5qJOV2eZzSUawireVKqurtxQoL\nsxWMw7UiJeedLOnkmO+nx/2X5B77SUnMz/O7mj1w5Qe09bdPqrW+sdva0Dkn9NuvcX/NFS9+nt/P\ns0v+nt+vs/elXC2TtEjSQ8aYOZLWdllbKen7xphMSRmSJkvq9XxFZWV9H6IkBr/fodzP8zudfeRQ\nDb/sTG374xKpLdy5OX/6eA371AX98mvs5++95O/5/Ty75O/5k3323opjX8rVY5LOMcYslxSQdK0x\n5suSyqy1TxhjfiVpqdqv6/ova233X7MBRGXqj29Q/knjtP8fbyrU2KK8KaM17qaLlV7oz98SAaA/\nibpcWWvbJH3mqM3ru6z/XtLvo319AD0LBAIaeeU5GnnlOV5HAQAcxRd3aAcAAIgXyhUAAIBDlCsA\nAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAA\nDlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByi\nXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkC\nAABwiHIFAADgUNDrAMC/as/fV2jPk8vVWtugARNHaNyNH1bGoEKvYwEAIIlyhQRjf/oXlf36UYWb\nWiRJ+//xhva/sEqn3vefyhk9xON0AABwWhAJpHHfQW27d3FnsTrs0PrtKvvVIx6lAgDgSJQrJIxd\njy1V84HqiGvVb5fFOQ0AAJFRrpAwUjPTe1wLBDnDDQDoHyhXSBgjLluorJGDIq4NnD0pzmkAAIiM\ncoWEEczJ1MSvfUzpJflHbC86fbrMN67yKBUAAEfiXAoSSunlC1U0d6q2379ELYcaVDhjvIZfcroC\nqaleRwMAQBLlKiGFw2Ht+OsLKn9+ldqaW5Q/fbzG3nCRgjmZXkeLi+wRJZr0nx/3OgYAABFRrhLQ\n2q/eqe1/ek4KhyVJ+xav1IGXVuvUP9+itAHZHqcDAMDfuOYqwVS89o52/N+LncXqsIP/fFeb73jc\no1QAAODgqMRSAAAXb0lEQVQwylWC2ffs691uonlY1aoNcU4DAACORrlKMIGUni/cDqRxUTcAAF6j\nXCWYEZcvVGpe5OuqiuZOi3MaAABwNMpVgsmdOELjb/qIggOy3t+YmqIhi+ZqzPUXehcMAABI4t2C\nCWnCly5TydmnaPcjLyvU3Kri+dM05PzZCgQCXkcDAMD3KFcJqmDaWBVMG+t1DAAAcBROCwIAADhE\nuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIF\nAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADgWj+SRjTJakByQNklQr6RprbflRz/mppPkd+/id\ntfb3fcwKAADQ70V75OpGSWuttQsk3S/plq6LxpiFksZba09Te8H6hjGmsE9JAQAAEkC05Wq+pMUd\nHz8j6eyj1ldI+mTHx2FJqZJaotwXAABAwjjmaUFjzKckfemozfskVXd8XCspv+uitbZRUqMxJk3S\nfWo/LXio73EBAAD6t0A4HP6XP8kY86ikH1lrVxpj8iUts9ZOPeo5hZIelvSStfZ7x3rN1tZQOBhM\n/ZezAAAAeCDQ00JUF7RLWibpQ5JWSjpf0tKuix0XvD8v6TZr7Z+O5wUrK+ujjJI4SkpyVV5e63UM\nz/h5fj/PLjG/n+f38+ySv+dP9tlLSnJ7XIu2XN0p6T5jzKuSmiVdKUnGmJ+o/WjVPEljJV1vjLm+\n43OutdZuiXJ/AAAACSGqcmWtrZd0WYTtX+/4cKWkn/chFwAAQELiJqIAjtBa16DWukavYwBAwor2\ntCCAJHPw9fXa+POHVLW6TIGUgApPMTI3X6m8yaO9jgYACYVyBUD1O8u1+qafqX77/s5t+xavVN3m\nPZr39I+VlpvtYToASCycFgSgLXc9dUSxOuzQhh3a8vunPEgEAImLcgVAjTvLe1xriFC6AAA9o1wB\nUHpxflRrAIDuKFcANOrqD0YsUZkjSjTmug95kAgAEhflCoDyThilqT/6tPJOHCcFAlJqigpOmagT\n/+cmZQ4p8joeACQU3i0IQJI0bNE8Db3gNFW9vUmB1FTlTxujQKDHP50FAOgB5QpAp0BKigpnTPA6\nBgAkNE4LAgAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEA\nADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABw\niHIFAADgEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDl\nCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIAAHCIcgUAAOAQ5QoAAMAhyhUA\nAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADAIcoVAACAQ5QrAAAA\nhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAA9MnmjRVa/uJm7dtT43WUfiHodQAA\nAJCYqg826K/3rtKmjQcUag0rMzOoSdMG64prT/Y6mqc4cgUAAKLy8AOrteG9coVaw5KkxsZWrX59\nl558aJ3HybwV1ZErY0yWpAckDZJUK+kaa215hOdlS1ou6WZr7eK+BAUAAP1H+b5D2rThQMQ1++5+\ntYXa4pyo/4j2yNWNktZaaxdIul/SLT0873ZJ4Sj3AQAA+qmK8jo1N4UirjXUNau5OfKaH0RbruZL\nOnwk6hlJZx/9BGPMV9V+1OrtKPcBAAD6qdHjBqpgYFbEteLBA5SR6d/Luo85uTHmU5K+dNTmfZKq\nOz6ulZR/1OecJWmCtfYGY8w8F0EBAED/kZmVppNOHa6Xny1TuMs5qrT0VM1eMEqBQMC7cB4LhMP/\n+lk7Y8yjkn5krV1pjMmXtMxaO7XL+p8ljZLUImmSpP2SrrbWru7pNVtbQ+FgMPVfzgIA6F1TU6tW\nvLJF4baw5pw+RllZaV5HQpIIh8N6+rF39eZr21VT3aiSQTmaf9Z4zTtzrNfR4qHH9hhtufqKpFxr\n7XeMMVdIOsNae2MPz71X0oPHuqC9vLw26a/NKinJVXl5rdcxPOPn+f08u8T8Xs7/2tKtevGZDaoo\nr5ckFRZl6Yxzx2v+B8bFZf987/07f7LPXlKS22O5ivaaqzslTTHGvCrp05K+K0nGmJ8YY2ZF+ZoA\nAIf27KzR04+801msJKmyokHPPPaetpRVeJgMSG5RXW1mra2XdFmE7V+PsO0T0ewDANA3K1/dqvq6\nlm7bmxpb9eaKHRozvsiDVEDy4yaiAJCkGhq6F6vjWQPQN5QrAEhSQ4bl9bg2eMiAOCYB/IVyBQBJ\nau6ZY1Q6uqDb9qGleZp/VnwuaAf8yL93+AKAJJeeEdS1/zFbSx5fr22bDiosaeSYQp2zaJKyc9K9\njgckLcoVACSxvPwsXXb1DK9jAL7CaUEAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADg\nEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHK\nFQAAgEOUKwD4F7Q0h/Tu23u0teygwuGw13EA9ENBrwMAQKJ4+dkyLX9piyrK65SSIo0cM1CLLpuq\nUeMGeh0NQD/CkSsAOA5rVu3WksffU0V5nSSprU3auumgHrr/LbU0hzxOB6A/oVwBwHFY/c8dao5Q\novbtrtU/l26NfyAA/RblCgCOw6FDzT2uVVc1xjEJgP6OcgUAx2FgcXaPa8NK8+KYBEB/R7kCgOMw\n98wxyivI6LZ97MQiTZ85woNEAPor3i0IAMdh5JiBuuITp+jl58q0e0e10tJTNXZCkS68bKpSUgJe\nx0tobW1hvfr8Jm1cX65Qa1gjRhVo4QfHKys73etoQFQoVwBwnCZOGaSJUwaptbVNKSkBSpUD4XBY\nf77rDa1+fVfntg3v7tcmW67rvniasrIoWEg8nBYEgH9RMJhCsXJk/bp9WvPmrm7bt22u1MtLyjxI\nBPQd5QoA4JkN7+xXW1vktV3bq+MbBnCEcgUA8EwwLbXntSD/iUJi4l8uAMAzM+eWKjM78uW/E6cM\ninMawA3KFQAgplpbQqo62KDWlu53uB88NE9nnT9RWdlpnduCaSk6dd5IzTl9dBxTAu7wbkEAQEy0\ntYX1lz+8odeXb1d1ZYMKBmZp6oyh+tAlU454Q8DCD07U1JOG6o0VOxQKtWnytCEaZ4o9TA70DeUK\nABATT/7fOi39x6bOxwf21+mljncAXnjp1COeWzIkV+d/5ISYZQmF2rR65U5VVzVq4uRBGjG6IGb7\nAihXAADnmptate6t3RHX1q7ao/M+PFlpvVzM7tK2TQf1yJ/e1u4d7e8+fC7daspJQ/Sxa09RKhfN\nIwb4VwUAcK66slGVFQ0R16oO1qu2Oj5/7LqtLazH/vJ+sZKkluaQVq/cpSVPvBeXDPAfyhUAwLm8\ngkzlF2b2uJabF3nNtXdW79HObZHvl2Xf3R+XDPAfyhUAwLmMzKBOOHFIxLUp04cqLT0+pwSrq3o+\nQtbU0BqXDPAfyhUAICYuvuJELTxvovIL2o9S5Rdmat7CMVp0+dRjfKY7U6YPUXZO5L9POGR4Xtxy\nwF+4oB0AEBOpwRR94sbZOuuCCao82KDCoixlZqYd+xMdKizK1oxZw7XsxS1HbM/Nz9D8s8bFNQv8\ng3IFAIipzKw0DR0e31LV1cUfO1GFRdl6d81eNdS3aNDgAZp/1liNmcC9tBAblCsAQFILBAI687wJ\nOvO8CV5HgU9wzRUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADgEOUKAADA\nIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHOIPNwMAgH6t7lCznn1yvbZvPihJGjmmUOcsmqQBuRke\nJ4uMcgUAAPqtluaQ7vn1Cm3bXNm5bcfWKu3YWqXPfGWe0jP6X5XhtCAA4AhtbWHt3FqpXTuqFQ6H\nvY4Dn1v24uYjitVh27dUatmLmz1IdGz9r+4BADyz5s1deuGZjdq5rUqBFKl0dKHOvWiSJk0Z7HU0\n+NSenTU9r+3oec1LHLkCAEiS9u6q0aN/els7t1VJksJt0vbNlXr4/tWqrmrwOB38KiMztce19Mz+\neYyIcgUAkCSteHmLDtU2d9tedbBBy17on6dfkPxOnj1SGRFKVHpGqk6ZM8KDRMdGuQIASJJqa5qi\nWgNiafT4gTp3kVFefmbntrz8DJ27aJLGTCj2MFnPojqeZozJkvSApEGSaiVdY60tP+o5n5B0o6RU\nSY9ba7/Xt6gAgFgqGJjV41phUXYckwBHOuPcCZo5d6TeXLFD4XBYp5w2st/ehkGK/sjVjZLWWmsX\nSLpf0i1dF40x4zqec6akWZLSjTFpfcgJAIixBWeN08Di7iVq0NABWnDWOA8SAe/LGZCh088ZrzPO\nndCvi5UUfbmaL2lxx8fPSDr7qPWzJb0h6T5JL0taZq1tiXJfAIA4KCzK1pXXnaJJUwcpJydNObnp\nOmH6EH38+lOVlc3vx8DxOuZpQWPMpyR96ajN+yRVd3xcKyn/qPViSadLmispS9KrxphZ1tqqvsUF\nAMTS6HFFuu4Lc9XY2KKAAhEvJAbQu2P+v8Zae7eku7tuM8Y8Kim342GupKNLU4Wkl6y1tZJqjTHv\nSZooaWVP+ykszFYw2PPbLZNFSUnusZ+UxPw8v59nl5jfz/P7eXbJ3/P7dfZofyVZJulDai9L50ta\nGmH9s8aYTLVf0H6CpLLeXrCysj7KKImjpCRX5eW1XsfwjJ/n9/PsEvP7eX4/zy75e/5kn7234hht\nubpT0n3GmFclNUu6UpKMMT+R9LC1dqUx5m61l6yApO9Zaw9GuS8AAICEEVW5stbWS7oswvavd/n4\nF5J+EX00AACAxMNNRAEAAByiXAEAADhEuQIAAHCIG5gAAACn9u2p0QtPb1RtTYPGTyrRCdOHKBAI\neB0rbihXAADAmZeWbNQ//m7V2NAqSXr1xc2advIwXXXdTKWm+uOEmT+mBAAAMXdgf51eeGZDZ7GS\npHCbtOaN3Vr6j00eJosvyhUAAHDi9eXbVV8X+U8Jb7IH4pzGO5QrAADgRFuorce11l7Wkg3lCgAA\nOHHC9CFKS4tcLUpHFcQ5jXcoVwAAwIkx44t08pzSbttHjinUGedN8CCRN3i3IAAAcObSfz9Jo8YU\navPGg6o71KThI/N1xnnjlZWV7nW0uKFcAQAAZwKBgGYtGK0LLpmm8vJar+N4gtOCAAAADlGuAAAA\nHKJcAQAAOES5AgAAcIhyBQAA4BDlCgAAwCHKFQAAgEOUKwAAAIcoVwAAAA5RrgAAAByiXAEAADhE\nuQIAAHCIcgUAAOAQ5QoAAMChoNcBAAAAotXWFtaqf+7U3t01yi/I1JwFo5WWnuppJsoVAABISLU1\njbrvzpXaWnawc9trr2zV5dfM0KixAz3LxWlBAACQkJ58aN0RxUqS9u2u1ZP/t07hcNijVJQrAACQ\ngEKtbdq8sSLi2vYtldq1rTrOid5HuQIAAAmnrS2s1pZQ5LVQWPX1TXFO9D7KFQAASDhp6akaVpof\nca1kyACNnVAc50Tvo1wBAICEdPq545Wbn3HEtvSMVM09Y7SCad69Y5B3CwIAkKRCoTa98MwGbXyv\nXM1NrRo6PF9nnDdeQ4bleR3NiUlTBuvaz87R8pc2q7KiQTkD0nXKnFJNOWmop7koVwAAJKkH71ml\nt1bu7Hy8c1u1NpdV6JP/MVuDhyZHwRo5plAjx5zidYwjcFoQAIAktKWsQuve2t1te8X+Or387CYP\nEvkH5QoAgCS08d39amlpi7i2b3dNnNP4C+UKAIAklJWd3uNaRiZXBcUS5QoAgCQ0a/4oDSzOjrg2\naergOKfxF8oVAABJKCMzqIsun6qikpwjts2eP0rzzxrnYbLkx3FBAACS1NQZwzTxhEFauXy7Ghta\nNOXEoRo6IjneJdifUa4AAEhi6RlBzV841usYvsJpQQAAAIcoVwAAAA5RrgAAAByiXAEAADhEuQIA\nAHCIcgUAAOAQ5QoAAMAhyhUAAIBDlCsAAACHKFcAAAAOUa4AAAAcolwBAAA4RLkCAABwiHIFAADg\nEOUKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOES5AgAAcIhyBQAA4FAgHA57nQEA\nACBpcOQKAADAIcoVAACAQ5QrAAAAhyhXAAAADlGuAAAAHKJcAQAAOBT0OkAyM8ZkSXpA0iBJtZKu\nsdaWR3hetqTlkm621i6Ob8rYOZ75jTE/lTRf7f8Wf2et/X3cgzpkjEmRdIek6ZKaJF1nrS3rsn69\npBsktUq61Vr7lCdBY+A4Zv+SpCs6Hj5trf1u/FPGzrHm7/Kcv0t63Fr7m/injI3j+N6fL+nbkgKS\n3pT0WWtt0twH6Djm/4qkKyW1SfqBtfYxT4LGkDFmtqQfW2vPPGr7IknfUvvPvHsS/Wf88eLIVWzd\nKGmttXaBpPsl3dLD826XlDQ/aLrodX5jzEJJ4621p6m9YH3DGFMY/5hOXSwps2OmmyXddnjBGDNE\n0uclzZN0nqQfGmMyPEkZG73NPlbSVZLmSpoj6VxjzImepIydHufv4lZJif5vPJLevve5kn4q6UJr\n7WxJWyUVexEyhnqbv0DSFySdJulcSb/wJGEMGWO+LukuSZlHbU+T9HO1z32GpE8bYwbHP2H8Ua5i\na76kw0einpF09tFPMMZ8Ve1Hrd6OY654Odb8KyR9suPjsKRUSS3xiRYznTNba1+TNLPL2ixJy6y1\nTdbaakllkpKpYPQ2+w5JH7TWhjqOWKRJaox/xJjqbX4ZYy5V+5GLpDk63UVvs8+VtFbSbcaYpZL2\nRTqCn+B6m79O0jZJOR3/a4t7utjbJOmSCNsnSyqz1lZaa5slvSrp9Lgm8winBR0xxnxK0peO2rxP\nUnXHx7WS8o/6nLMkTbDW3mCMmRf7lLETzfzW2kZJjR2/3dyn9tOCh2KdNcby9P7MkhQyxgStta0R\n1rp9TRJcj7Nba1skHTDGBNR+FOMta+0GT1LGTo/zG2Omqv200KVqP0WSbHr7d18saaGkkyQdkrTU\nGLMiyb7/vc0vtf9y8a7af4H8YbzDxZq19hFjzOgIS8n+M69HlCtHrLV3S7q76zZjzKOScjse5kqq\nOurTPiVplDHmJUmTJJ1sjNlrrV0d47jORTm/Ok4DPizpJWttMvzQqdH7M0tSSpcfsEevRfyaJLDe\nZpcxJlPSPWr/AXtTnLPFQ2/zXy1puKQXJI2W1GyM2ZpE11j2NnuFpNettXslyRjzitqLVjKVq97m\nP1/SUEljOh4vMcYss9aujGdAjyT7z7weUa5ia5mkD0laqfb/gy3tumitvfLwx8aYeyU9mIjFqhe9\nzt9xwfvzkm6z1v4p/vFiYpmkRZIeMsbMUfvpkMNWSvp+R8nIUPsh83XxjxgzPc7eccTqcUkvWGt/\n7FG+WOtxfmvt1w9/bIz5jqS9SVSspN7/3a+SNNUYU6z2/7DOkZRsFzX3Nn+lpAZJTdbasDGmSlKB\nBxm98J6kCcaYgWo/anm6pP/xNlJ8UK5i605J9xljXpXUrPbTAjLG/ETSwz74zaXX+dV+YfdYSdd3\nvItOkq611m7xIqwjj0k6xxizXO3vjLrWGPNltV938IQx5ldqL5kpkv6r49RosuhxdrWfDjlDUkbH\nO8ck6ZvW2hXeRI2JXr/33kaLuWP9u/+mpCUdz33IWptMv1RIx57/bEmvGWPa1H7d0XMeZo05Y8yV\nkgZYa3/X8XVYovafefdYa3d5my4+AuFwMr5JDQAAwBu8WxAAAMAhyhUAAIBDlCsAAACHKFcAAAAO\nUa4AAAAcolwBAAA4RLkCAABwiHIFAADg0P8H7Blp4AJXayQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10b883630>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10,10))\n",
    "plt.scatter(new_data[:,0], new_data[:,1], c=t, cmap=cm.Spectral)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "def plot_kmeans(data):\n",
    "    fig, axes = plt.subplots(2,2)\n",
    "    for k, ax in zip([2,3,4,5], axes.flatten()):\n",
    "        t0 = time.time()\n",
    "        centroids, iterations = kmeans(data, k=k)\n",
    "        t1 = time.time()\n",
    "        labels = [find_closest_centroid(x, centroids) for x in data]\n",
    "\n",
    "        #PLOTTING\n",
    "        #plt.title(\"K-means with k=\" + str(k) + \" time=\" + str(t1-t0) + \" iterations=\" + str(iterations))\n",
    "        ax.set_title(\"K-means with k=\" + str(k) + \" time=\" + str(t1-t0) + \" iterations=\" + str(iterations))\n",
    "\n",
    "        ax.scatter(data[:,0], data[:,1], c=labels, cmap=cm.Spectral)\n",
    "\n",
    "        ax.scatter(centroids[:,0], centroids[:,1], c=np.arange(k), cmap=cm.Spectral, s=300)\n",
    "        ax.axis([-6,6,-6,6])\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "plt.rcParams['figure.figsize'] = (20.0, 20.0)\n",
    "plot_kmeans(X.T)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}