tutorials/T10_Classification-Clustering.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# T10 Classifcation and clustering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cross-validating support vector machine classifier "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In today's tutorial, we will use the constructed data example from class to check how well a trained classifier separtes the data. \n",
    "\n",
    "The data was obtained through (hypothetical) extracellular recordings *in vivo* and post-experiment reconstruction of the recorded neurons. Firing rates (in spk/s), spike width is extracted from the extracellular recordings and the reconstruction yields the identity of the neuron, either an excitatory pyramidal cell or an inhibitory fast-spiking basket cell. \n",
    "\n",
    "First we will use the the **Support Vector Machine (SVM)** classifier (from the scikit-learn implementation) to learn the separation of the data in the two neuron types. \n",
    "\n",
    "Later, we will check the performance of the classifer through cross-validation. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1. Construction of the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "np.random.seed(1020203)\n",
    "spikeWidthPy = np.random.normal(350,50,200)\n",
    "spikeWidthFS = np.random.normal(190,40,100)\n",
    "\n",
    "spikeWidth = np.concatenate((spikeWidthPy,spikeWidthFS))\n",
    "\n",
    "firingRatePy = np.random.exponential(3,200)\n",
    "firingRateFS = np.random.exponential(10,100)\n",
    "firingRateFS +=3\n",
    "\n",
    "firingRate = np.concatenate((firingRatePy,firingRateFS))\n",
    "\n",
    "neuronType = np.zeros(300)\n",
    "neuronType[200:]=1\n",
    "\n",
    "Xdata = np.column_stack((spikeWidth,firingRate))\n",
    "Ydata = neuronType"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2. Data visualization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYkAAAJOCAYAAABC2PssAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdeZxU1Z3//9cHUXEBxA1aDKKCxAXE0G7EBRdcos5k1+gYTfKzJcmYzYlxi2PyTdRJMi4xKmaSjCbR6EQnmbjEKGqDcQElEsEFJQouNODKDgp8fn+ce6nb1XWrq7tuLd31fj4e9eiqe2/dOnW6+37u2c3dERERKaRPrRMgIiL1S0FCRERSKUiIiEgqBQkREUmlICEiIqn61joBWTKzp4EdgHm1TouISA8yAnjT3ffL39GrggSww8CBA4eOHTt2aK0TIiLSU8yaNYulS5cW3NfbgsS8sWPHDm1tba11OkREeowJEyYwderUgjUwapMQEZFUNQ8SZtZkZjeb2ZtmtsbMnjOzwxP7zcwuNbOFZrbazFrNbO9aplky9NIUSI76dw/bRKQu1DRImNk2wKOAAScAewLnAEsSh50HnBtt3z/a94CZ9a9uaiVzL02BWz4F910QgoN7eH7LpxQoROpErdskzgPa3P3ziW2vxE/MzIBvAFe4+53RtjMIgeJU4MYqplWyNuIoOPDLMP2G3LbpN4RtI46qXbpEZKNaB4mPA/eZ2e3AEcBC4BfAdR5mHtwVGALcH7/B3Veb2TRgPAoSPZsZHHd5eB4HigO/HLaZ1S5dIrJRrdskdgO+ArwMHAtcA1wBfDXaPyT6uTjvfYsT+0REpEJqXZLoAzzl7hdEr582s5GEIPGzxHH585lbgW0AtLW10dzcnPqBLS0ttLS0dD/Fkp24DSKuYoJciUKlCZG6UOsg0QY8l7fteeDr0fNF0c8hwGuJY3akY+kCgKamJjROooeY92AuQMTVThC2jTgaRh5du7SJCFD7IPEoMCpv2x7Aguj5K4RAMRF4EsDM+gGHAt+uUhqlUkYeDafdGRqp41LDcZcrQIjUkVq3SVwFHGRmF5nZCDP7DPA14DqAqPH6auB8M/ukme0D3ASsAG6tUZolSyOPbl+tZKYAIVJHalqScPcnzezjwGXAd4FXo5/XJw77EbAFIXAMAqYDx7j78ionV0Sk4dS6ugl3vwe4p8h+By6NHiIiUkW1rm4SEZE6piAhIiKpFCQagSbRE5FuUpDo7TSJnoiUoeYN11JhmkRPRMqgINHbaRI9ESmDqptERCSVShK9nSbRE5EyKEj0dppET0TKoCDR22kSPREpg4JEI8gPBppET0RKpIZrERFJpSAhIiKpFCRERCSVgoSIiKRSkOhpNFmfiFRRJkHCzJrN7Jtm9mMz+76ZfdbMts3i3JKgyfpEpMrK6gJrZmcS1qR+BZgJzAX6AYcA3zGzOcB33f3VMtMpoMn6RKTqyh0nsRXwUXdfXWinmY0FRhLWrpZyabI+EamysoKEu1/Xyf5Z5ZxfRERqK6s2iR+Z2QAz29TMHjSzt8zsX7I4tyTkT9YXVz3FbRQiIhnLalqOY9z9PDP7BPA68BngYeC3GZ1fQJP1iUjVZRUkNo1+fgz4nbu/Y6ojz54m6xORKssqSNxlZi8Aq4GvmNkOwJqMzi1JmqxPRKqorDYJM2sCcPfzgYOBZnf/AFgF/HP5yRMRkVoqtyTxKzMbBLQC9wF/BXD3lcDKMs8tIiI1Vm4X2OPNrB8wAfgE8BMze5UQMO7TIDoRkZ6t7DYJd19DFBQAzGxX4HjgZ2Y2xN0PKPczRESkNrIaJ3F8/NzdX3H364E/E6bnEBGRHiqrWWC/a2ZHxi/M7Dzgn9z9/YzOL9J7aCZf6UGyChL/BFxmZoea2Q+BA6NtIpKkmXylh8lknIS7v2Vm/wRMIcwG+2l3zRMh0oFm8pUeptypwpcDyWCwGbAb8Gkzc3cf0MXzXQj8ELjO3f812mbAvwMtwCBgOvBVd3+2nLSL1IRm8pUepqzqJnfv7+4DEo9+7r51vL0r5zKzg4CzgGfydp0HnAucA+wPLAEeMLP+5aRdREQ6l9W0HJjZJwm9mRx4xN3/2IX3DgRuAb4EXJLYbsA3gCvc/c5o2xmEQHEqcGNW6RepivyZfCFXolBpQupQVl1grwcmAbOBOcAkMyu61kSenwN3uPtDedt3BYYA98cbogWOpgHjy0q01Jdq9Piph15F+TP5Hnd5ro1i3oPVTYtICbLq3XQ4cKy7/7e7/zdhNtgJpbzRzM4CRgDfLbB7SPRzcd72xYl90tNVo8dPvfQqimfyjUsNcRvFaXdqokapS1lVN80FhgELotcfomPbQgdmNgq4DDi0kzEV+T2lrMA2ANra2mhubk49UUtLCy0tLZ0lTaqpGj1+6qlXkWbylR4kqyCxHfC8mc2IXu8PPG5mfwJw97QxEwcD2wNzEutPbAIcZmaTgL2jbUOA1xLv25GOpQsAmpqaaG1t7ebXkJqoRo8f9SoS6ZasgsQlnR9S0B+Bp/K2/TfwEqGE8SKwCJgIPAkQTSh4KPDtbn6miIiUKKsg8RSw2t03mNkewIeBP0drS6Ry9/eA95LbzGwl8I67z4leXw1cFC1q9CJwMbACuDWjtEutVaPHj3oViXRLVkFiGnBotLbEg4SgcTJwWgbn/hGwBXAducF0x7j78gzOLfWgGmt3a31wkW7JKkiYu68ysy8B17r7j8xsVndO5O4T8l47cGn0kN7mpSmh4TheuxvCBT3rtbu1PrhIt2TVBdbM7GBCyeGeaNsmGZ1bqqXa4wiS3VLjABF3S533YPYX75FHt69WUq8ikU5lVZL4OnAB8Ad3f9bMdgMezujcUg3xBTtZHRPX4VeqD389dUsVkYKymgV2GqFdIn79MvC1LM4tVVKLC7a6pYrUvXJngf05oQ1idoF9WxEar9e6+y3lfE6vEde/xxdA98pUq3SHLtgiUkC5JYnrCavSjSbM2fQm0A8YCQwAfkWYuE9qUZ1T79QtVaTulRUk3H0W8Fkz2xpoBpqA1cDz7j43g/T1HvVe/16LC7a6pYrUvazaJFYArVmcq9eq9+qcWlyw1S1VpO5ltp6E9HC1umBrsjuRupbVOAnpTH51Tlz1FE9dXQ960jiCelgbQqQBZBokoh5NUogWm8lOvawNIdIAMqluMrPxwC+ArYFhZrYvcLa7fyWL8/cKqn/PTr13AhDpRbJqk7gKOBaI14/4u5kdltG5e4/eUv9e6/Ee9d4JQKQXyay6yd1fy9u0PqtzSx1RVY9IQ8mqJPFaVOXkZrYZYUqO5zM6t9STeqjq0SA8karJKkhMAq4BhgKvA/cDao/ojeqhqkeD8ESqJqvqplHufpq7D3b3Hd39X4A9Mzq3SHtxJ4A4MMWBKzm9STW7yKo7rvRiWQWJa0vcJj1dvYz3KDamo5rtJmqjkV6u3FlgDwbGAzuY2bcSuwagRYd6p55Q1VPNdpN6aKMRqaBy2yQ2I4yN6Av0T2xfBny6zHNLPeoJ4z2q2W5SD200IhVU7iywU4GpZnaTuy/IKE1S73ryeI96WsNDpAfIqnfTKjP7MbA3YT0JANz9yIzOL1K6ZLsJwJDR0XMHB2ZMzm4ND3XHlV4uqyBxC3A7cCKhO+wZhAWIRKpvY7vJpFxQGDIapk8O+7NsL+gJbTQiZcgqSGzn7r80s68nqqCmZnRuka5JtptAuJuvVHtBT2ijESlDVkHig+hnm5mdACwEds7o3CJdF1+gq9Ettye30Yh0Iqsg8QMzGwicSxgfMQD4ZkbnrqxaT1YnlaP2ApGylR0kzGwTYKS73w0sBY4oO1XVEg+EStYnxxeVrBo2pXbUXiBStrKDhLuvN7N/IkwX3rNoIFTtVbIkp/YCkbJlVd30mJn9jNDDaWW80d3/ltH5K0MDoWqrGiU5tReIlCWrIDE++vn9xDYHNE5C0qkkJ1L3MgkS7t5z2iGS1LDZuUpWB6kkJ1L3MluZrkfKb9g87vLcne28B2udutrTDKciDS+r6qaeSQ2bxVW6OkglOZG6V9MgYWYXAJ8ERgFrgSeAC9x9TuIYA/4daAEGAdOBr7r7s5kkQg2b6SpdHaQuqiJ1L5MgYWafLLB5KTDb3ZcUeesE4HrgScAIDd9TzGwvd38nOuY8wiC9M4G5wCXAA2Y2yt2XZ5F+qZHeUJLTYEzp5bIqSXwJOBh4OHo9gVAq2MPMvu/uvyn0Jnc/NvnazE4nBJePAndFpYhvAFe4+53RMWcAS4BTgRszSr8UUo3qoJ5cktNgzLqzeNkaFry9il2225LBA/p1/gbpVFZBYgOwp7svBjCzwcANwIHANKBgkCigP6Ex/d3o9a7AEOD++AB3X21m0wjdbhsvSFTzzlXVQcWpC2/dWLl2Hd+47WkeeD5UXJjBiaObuOJTY9hq88Zuei1XVrk3PA4QkSXAHu7+jpl9kPamAq4BZgGPR6+HRD8X5x23GBjarZT2ZNW+c+0N1UGVpC68dWHl2nVMvHIqC5eu2bjNHe56pg2Aa0/9SK2S1itkFSQeMbO7gd9Hrz8FTDOzrYD3SjmBmV0JHAIc4u7r83bnT+VpBbYB0NbWRnNzc+rntLS00NLSUkqS6k8t7lyrWR2UVkoC1ftLqm/cPqtdgEi6e3YbFy9bo6qnMmQVJL5KCAwfJVzAfw3c6e5OCRP+mdlVwCnAEe7+cmLXoujnEOC1xPYd6Vi6AKCpqYnW1taupr9nqMWda7Wqt4qVkqA+6/3VhbfmFi9bw5TnCl4KgPArevWdVQoSZchqxLUDd0SPLjGzawgBYoK7v5C3+xVCoJhI6AGFmfUDDgW+XU6apQTVrN5KKyUdMCncdtRjvb/abGpuwdurClcpJAzbdsuqpKW3yrIL7H8Q7vAteri7D+jkfdcBpwMfB941s7gNYoW7r3B3N7OrgYvM7AXgReBiYAVwaxZp71EqeedaqMSAV696q1gpKRxQX/X+cX4lV8Cb96DabKpsl+22xCx9bamJe+6oUkSZsqpu+hFwkrs/38X3fSX6mT8HxveASxPn3gK4jtxgumMacoxEpe5ci5YY7sh9BtTHBbrW1PW1bgwe0I8TRzdtbKRO2mlgP64+Zb8apKp3ySpILO5GgMDdO73SRFVZl5ILGo2rUr2NijWI735UafNYldt2kVZKco+qmyYXLj3Ne7D6jdrq+lpXrvjUGCA0Usd/LkfvNZirTx6r7q8ZyCoHnzKz24E/EqbXAMDd/zej80usEr2N0qp6jr0M/nJh2DbqBBj4odz+Yy+DfzwUPjuLO+tipaQ4Pfnbt9wOHv5B9e/o1fW1rmy1eV+uPfUjXLxsDa++s4ph22ogXZayChIDgFXAMYltDihI9GTxhXvUCTD3HjhwUmhInn4DvPdq2BaXbMq9sy5WSoLC20ccBaveLv65PWXajJ6Szjo2eEA/BYcKyKp30xeyOI/USLEG8VPvCBevuERx4KRE0Phy7sJW6M5696M6fk6xC1+ppaTk9mJ39JVqO8i6A4HaOKSOlRUkzOw8d/+RmV1LgcFt7v61cs4vVdJZg3ifPl2vXnnvVbj107lgkQw0p90BVGGOpkq1HWTdgUBtHB1oDqb6UW5JIm6sfqrchDSkeqliKLdBPP/OOq6KGjI6bJt+AwwaDu/ODyWRl6bAjMnl3yWXckdfibaDrDsQqI1jo5Vr13H+nc/kGqE1B1PNlZXr7n6XmW0C7OPuGtzWFfVWxVCsqif1Yuxw3BXt2y5GHAW3fDoEiEWzYY/j4cU/hwAxaHgob86YnM1dci0Hs2XUgWDjHfO2WzA4o6T1ZOff+Uy77qyag6n2yg7N7r7ezMZlkZiG0pOqGPIvxi9NiUoIk6MG5KNz7RTNX8h9r1EnwNJXc+d5d34uQGRxl9zZHX05bQdTfwyHnhuq2gA2bIBH/hMOz+ZeqP0ds3NJ39/wxb738UHz2Wy6SZ+GnN5j8bI13D2743gH0BxMtZRV+e1pM/sTYYK/lfFGdYEtoidVMeRfjEceHXo5zZgcAsi8BxMN2VHQwEMQARi8T/i5OFpw8N0F2aYtKXlH392SxtQfh661z/8ftEwL235+WCgZQSaBInnHfHifZ/hi3/v41brjeHrpybk75gab3mPB26tSR05rDqba6ZPRebYF3gaOBE6KHidmdG6pByOPzgUvMzj+ilBSiNsc4nEV8cC79xLzMS5vCwEiDhYv3gu3nZY+lwKE0kpyv3u4eOdve2lK8TSfdmcu8MaBubPqvEPPzVWX/fywXIAYMjrsK6ZQuvPSmH/HPHXDvpzx/nf4/rrTuXvOIhYvX9sunYuXrWHGK++weFnhmU7rTXfTG0+xUYiZ5mCqFXWBrZWePoPovCmh9BBzh79cEEoPQ5vhjadgj4/Bkufgvflh4NsuHw3BYsjo8N60RvpC7TW3nRbes/BpOOWWsK2UNpzutB306RNKEMnSw5DRYVufIvdVJbYzFbpjnrphX6D9HfPKYRM4/9a/9ZhG3HIbnYtNsXHi6CaVImpEXWDL1d0eSl2pCqmXXlCx/LvjLbcLVU+xOECccgvcdz7MuDEMeptxY/i+x/wQ/nplevoLtdfEvaXm3hMuvFB/bTgltjMVm5Quecfc0xpxs0hvhyk2EoFGaqPc25Hnop+N2QW2nB5KpXajrLdeUBAC1IzJcMDZsOAxWDy74zGDhoVxETPyV5j1aHsUVA77t9z8S/G5RxyVe502TUi8rRKlrg0b2lcxQa7qqVhposR2plLumHtaI25W6dUUG/Wn3CBxMnA3sI27X5NBenqW/DvHeHxAfOdYzgjjuPSQ/Iz3XoVthnX9Djrrkkgc4HY/MlQxJYPE/mdDH8s1Wsd1+/HzePseHwuNwwufzk35EXePjXtKjTqhe+krVVq+LHw6FyDyG64z6uHU2R1zvTTiljqoLev0aoqN+lFukBhnZrsAXzSzXxMmYNzI3d8p8/z1rdCd45DR4Y4XQlXL9G4MGssvPRx7GSz4a64NoCt30JUoibw0JQoQF4bvN3ifXM+lVx8NgSy2aHYocbz3WmiwjjV/AQbtEtKRDB5xlVL888BJbFxLYsFfw/mynApj1Alw8m/D++N8OeJiGPM5+OfrcqWGs6bC/321eIBItjMVmxCRzu+YS62SqpRS2heSAaTW6a2VRhgZXm6QmAzcB+wGzKR9kPBoe2NZNDsEB0vcTRfrxVNIfgnFPXc33lVZj8dIXlzn3tM+QGyzS3geN04P+2goGSx4NHcMhPeOnBgecXpii2bnzg3hr+i4y3KlNAjTfMTBrZypMOLPufFQ2OWQkNbB+4QSDsAW2+QC6/0XwTO/g9GfTv+sUidETLw/7Y651o24xdoXrvjUmIIB5Li9h/DnOYs6nKtQenv6xbWRRoaXO+L6p8BPzewGd/9yRmnqOdr1UJoE8x8NVS/JevgDJnX9AlaohALt76qhtDvorMdjxAFv7j2w8/7w+pO5fWNOhWnRZy2aHXozJaubDpgUBc8bQr7FJa582wyDU3+fa/swC3fl8Tni7rjlToUx7szwPeLAlgx4e3ys64G1XTVckQkRS1SrRtzO2hfWrNvAA4l1peMAcvzegzlpTFPR9GZ9ca1VsOlpnQrKkVUX2MYLENCxh9KGDfCfI0NPHggXyKwaVoeMhmMT56rVQKvkV0kGCIAPjYMJF0HTvvDyw+0vsgdMCmMrYsm76w7tFnlBMC24dWUqjELtDw4MHp1rU4kDxIGTQl53p5G8lNlpS1SrRtzO2hemJAJE0n3PLeaJC47i4hP3Sk1vVhfXWt7J97ROBeXKajBdY0oO1gK4/8JcgIBw4fvLBV2vbsofQzHqhOhcF4b9pQwISztXXPV022khqCWPKzYwLTbi6NDGkO+As2GDQ+sPQ917fo/o3Y9qP6DtiItzjdOLZueqZuLqpuk3wD9SVsQrNa2xuIrsvuh34Q5/Ph9+95kQILbcLnfsltvBMZeFz84fFDevC5+ZocED+rH/8G1LuvBkMfCu2KA2KNDXPd6eaKAulN7OLq5dSXMcbOJfURxszr/zmZLP0V2lNNJXUrUHV/auyrNa2Ngb6YFcG8QBk8LPGZNz8xvF9e+lKDSGIr7Qx6WHUu+gC50rvoO//V+6NjAtlnYF+Ud0EZ1RoGfTyw/CHhNzgeLwb8NO+4W7+2QX2JETo21TcmmHUBUUL2eKhyq9UgNlsl0mvu2M07jNLvBeYpqQVW/Dr44JYz2g/e9yxuQwzXlnv8saDJTM8s568IB+NA3ox8KlHS9Cg/tvzpIVa7vVQF1uD6i4amnLzTap6Z18rRrpa1V6UpDITPSPn1+tMmMy7etoSpDlVNSFznXyb0OASBuYVqzL7LwpuWCYFLfDjPoYzI16MS2aHbVDEAXLie3THz8vtG3hrJCmPT4WekUtnhMu6MlBe6n3tHnMwveZfkPu/QdMgjdmwhtRldmBUTCYPjkXIKDjr66Uj6zB7LSlVOOUWn+/eNmaggECYPHytUzca3C7NolYZw3q3b24drg4Urw0M+u19xi05WYVa6eoVaeCWrWDKEhkpdDF+Pgrol483bgoZLmWdf77+vSJRkNf0LHOfN6DxbvMJv87DzgbXn0sV1qIu5PGAShO97GXdQwQSXFQiksUL00JPYxGnQD7nREG5k2f3P6O/8BJuaVN889TKLjtflT7ks2rjyYa1M8OU57Hpk8ObStts3IB8cBJuSqzzmS93kQnOqvG+dZbK7jy/hdT70Dzg8eCt4tXl5x+8C7069unyw3qxS6uB++2XYF3hO/2zdtn8dg/ctW4ncXpSb+d2emddrkN3ll1KuhK4K5V6UlBIktZXthrpbMus/EEfhsbdy9o3z33vgtyPXni90L7O+qk/C61cckDwuu598DI49u/J78RP3megsHtjpDueIDcotm550deEgL5xov5FSGgedS+khT3tiqlmquKfwudVeNc/Mc5PDrv7Xbb7nqmjXUbnL59rMOF7ssTdk/9LDMYNbh/txvU8y+uscf+8TYHXf5guwtt/p1zqfLbKSB3p71y7Tq+cdvTPPD8ko3fp9j4j7TvVW6ngq5WHXX2O65k6UlBohF1VmfeWc+c0+4Id9XxYLp4ZHRcejhgUunVLMmgNGR0rjoo2evopT+3f8+i2SH9x1+RS1Ox4ObR8zgAxUFt0WzAciWYuL1k9yPh9tNy50kO9jtgUv3MFRXZZbvideDJu/Ck/DEN8UV1xvz0MbDJ6pTujIpOXly/9runmf5K7rPyg1d3AkQh8Z32qvfX8ekbHuftle93+ExIH/9RrM6/uyPDu1p1VKyqDuDs38wEKtNOoSDRKJJVMXGd+c77h6qgZLfaLbcL8ymlnSO+W3/j6dyFfe49YVwDRHf+R5dezZI2JmTx7PbVQ9C+3WBG1CFgj4np50kGt1PvCA3r06NFj+L1Lmb+d67kE7//9n8J7SqjPgYDh7VvB0l+tzJl0cd/zhtLmbngXbbfejPeWvF+h/0D+/Vl6Zp1XUzX2tR9nx+/CzNeeSeTO9ZkgEgqNCCvHO50qLLKd9czbbQtXcNTC95t977utusUs3jZmtQAeNczbVx8Yseqo2JVdUmVaKdQkKg3afXq0HH7tJ+EC3op2x7+Qe5COOIo2PkAeH1GrofTbkeGOv/kfEqjTsjNFQUhoCTv1gcNz3VfjRceOuLi9r1/yqlmWTS7fUO4kwgEk7t2sY5HwLdrTLbcCOn8GWdHnQCf/U0YaR0bdULHdpBuyKKXypJla/j4dY+mNjDHuhogOvOZyU8A3Utz8gL79Kvvdv6GDBULELFkgEi6e3Yb33pzBVc+8GK7arLxu2/HN48eiWNsudkmrHp/fUnBo7M2n7ReXoXaQdJKFlm2UyhI1JNi9epQfH2FeVPgyV+FnkALnw4NyPMeDP398y+E770aAkRcCojPBdFFOTF30ql3hO1xldFxl7PxLvzd+eG4eGK+A7/csRTS2eSCyaqv/JIDwIAPwYQLYdU74TNGTsy1GySDj3vueySr0N6YCV/8S/vR0PHnx6WcEUe1z+e4uixugE8Gy9tOC/NOjTg6/Tt1IoteKqUEiErqSpoLBcVxwwZVI5mZcIeL/jinQ6B57B9vFww+J40pHjy33GyTop+3xaaF9+e3g7y94n0m/XZmapqzmgRSQaKepNWrb+xGmrK+QvIiPzjallwwp9A020NGh0nr7r+wfZfWgcPaTyMRX/hGHJVbNzopbgQuNKK4WNA74uLcNOHTb4Ch++e6o8bBYmgzPBl1rT01b3xCMsjEM7fGebLbkSG9L/45BMNfHQOHfSe8J388SKHv9OqjYRW8ufeE6re4Oi2eZHHuPazc7yyeHX1BuHN87Hsbz7l48CFFqyOy6KUy542lZQeInQYWHgfRVfnVI4WqYwoFxbS79izFf4ldHMra8TxWWkkkdtczbbzy1kqu+NQY9hk6sMP+Ve+vL/r+N95b3aFUkp+v8XTy1RivoSBRT4rVq4cDiq+vALDL+PaNswdOKtyzKG78fe2x3LZkw3GyS+ytnw6vdz8qt/ochOqmd+eH54X+UtOC3qgTQrXWqrej0dcXwcNRT6J4XMSg4bnxCnGwKtY999Q7csHtd58J3+Xd+eE8rz8ZvkOhOZSSJZkDJrXvGrvN8LCq3pDRG3tHrRt5PI8s2YIjnv4v5jz5GnOAL/a9j1X7ncV3Zgzi7jkPFq1CymJK7ZllXGD3GTqAHfpvzsMvvNntc+T71v/M4qrPjuX/3f1cu9LC0R/ekRP37bwevVL6GKwvN0IQuud2JUgAzFm4jBOv/Ss7DezHH7/6UXYscYZfaN+F98Bdt2WLzTahde6buXzdczCnH7QLo4b0r8p4DQWJ3mTI6AKL/ER/jRsnnEs02sZ36UNGh8n4Cg2Sy1/PYmOJZZ/cbK+LZucm4sufW6lQ0MsPbsn/lkHDOk4KGJ+zWA+muFQQT88ed3H90Pjc90w2UMdpzB/45t5x2dJEaemb73yGuxa3cUnfVXyx730A/Grdcfzi2eNYuCzX4JpWHZPFaN29durf6TFp5ryxrNvvTfPovLc5/ppHOvQaeuD5JRu7mqYZv3vXL8ClyiJATLTS/F8AACAASURBVNxrMCeMbuLxl9/u8uw6AAuXruH4ax7h3q8fCrCxNFCsETrZhfeJl9/psO+B5xbzwHOLw7/X3kM4fu/B3Pfc4opNAqkgUU/SuqbGw0w39syh4/oK7u174Wy88N/QfjGk+EL44n25UsCHxre/aiXbLwp1iY3v9uML/j8ezF1sSxkwVih4JEdm5x9b7H3Ji7477YZIJ8dDpE21Xmjg2y4fLXj8yvfXcfecwv/YC5etaf/ZkfwqpMED+qVOqX3EqB02NmoWuwv89WMLUvfVSjJAdMUPP7EPW27Wl1ffWcW/3vq3oj2rKqFYtdvg/ptvvCCX4+2V73PgZe3nITt+78EcMWoHHp7b/RKde+gJdtKYJp644KiKTQKpIFFPik3nAIXnXxp1QnRnfkHu+LjX0bCPhuqm6ZNzbQAQ7uLfnR+6wO40LjFVxdnR3ElHt58rKn9cwKBdojaCqPE2fk9+gCg2HiNtmvBYvB5FqXMe5U/bvuDR9uMhik2z3qERfXL7c0Rp2WrGZL67yWuA88W+f+FX644D2Fii+P6608kPFAWrkFJuSR964U0eeuHNdneDK9au2zhf0ar31xedt6gnOuXnT/DQuRN4Z8XaqgeIk8Y08d0T9+L/3f1cwbv6xcsrl54/P1te4Em6e3YbF5+4F/sP3zazcyYpSNSTYtM5QPvtp9yS6+4678HchW2LaJxDXJ1z6h3tewK9NKVjIIqn4x6RWAioUM+f/At9MiAU6upaLOglSzeFSkHxehTxoj2dpSWemyleBnXR7Lw2E8LiRfF7CpV4kund/aiQp3Fajvp3Vn6wni8+/V9AqGIKQSH4Yt/7mLphX6Zu2LfdKeMqpOTkdJ1dIOKqqpkL3u1wl1ts3qKeaPGytez3/ft5P4u6oS5qW7qGx19+i7MP371m7SZZiEdcH7v3kIqc37w7FW01YGZfAb4NNAHPAt9w90fyjmk9/PDDD29tba1BCmusK+tYd/XYcpY/LfRZ+eM25k2BWz4d9k+4qH2QO+0OIApAnaUlPn/c0H7sZfDIT2Dl27npNEYcxfzpf2Lqhn0Zt8ugjb1P4ov4yGVPMGjM8SG98TKt8bKj7lz/ixuZ/so7TN0whmT/mU/0f54/LN+rw9c/fu/B9N2kT0mT01XLgbtuy/y3Vlb0Tlmqq9yR1hMmTGDq1KlT3X1Ch3P3hCBhZicDvwW+Avw1+vkFYC93fzVxXOMGiUrqSlDp7jlfeiBcPeMR1GmfUUpaUo5ZMviQDuMLhgzYnDE7b8MDz3ds+Cv0z1aoz/9xew1mvcP9eXXXx+89GMwyH0EskuakMU3dGmndG4LEdOAZdz8rse0l4A53vyCxTUFCUo2//MGSxwZ09s+2ODGx2w9S6rQn7rljp717quGwkduz3r3dJH9Sv35/9kFcNeWlbvX6MoMnLjiqy43XxYJE3a9MZ2abAeOA+/N23Q+Mr36Kivv5z39e6yT0ClnnY1cHoHW2Ulq8+lp8bCH1ECAAhi1/rmLdTCVbJ41pYv9dt+PWsw5i+oVHcePp47o0+0wlVsar+yABbA9sAuS39i0GKtNSUwYFiWxknY9dHYBW6j9bscFx9eIPD8+o+zQ2mgH9OlZlHr/34HbjGwYP6MegLTfr0u+uEivj9aTeTflZVbANsK2tjebm5tSTtLS00NLSknHSpN6N26VrcwWV+s/W2ejZUj/rxNFNfH78Lnx28hOZN2z3Xfpa2WlMY8BNX9ifM/77yexP3ovdetZB7NB/c2a99h4AYz+0TcEqoq7+fVViZbyeECTeAtbTsdSwIx1LFzQ1NaE2Ccm3z9CBXZqvqNR/tmJTOE/cc0emvLAkdXT1Xf96CKs/WN9uANTRKUuDxn5/9kH8+vEFJXfZnLjnjjw9ZVFqGgcP2Lys8Qknjmni8FE7ctKYypy/J4qD/pPz32FRge++08B+G3vVddZttdjf104D+9G2bE3FRlrH6j5IuPv7ZjYTmAj8PrFrInBnbVIlPdEfv/rRLvVuKlWxpSzTVlc7cXRTwcnfrj55LBOvnFowmMX11fvvuh0Xn7iGb/3PrKKN0TsN7MfVp+zH4dcWT+OKteuYu2g5v3l8fmo7ypABm7PfsEHc9+yigvlUaMW5iXvuyNWn7FfyCnMH7bYt504cBRZmQl39wXq22HQT3nhvNWs+WM/dz7QxJfo95efLt47Zg7dWvM/Pp/6jrLagbbfclHdWfdBh+6Z9jA/Wb+h0evrxu2/HVSePZfCAfgWnc4/ncuqKzn53lRppHespvZtOBn5D6Pr6KDAJ+BKwt7svSBxX895Nzc3NPPXUUzX7/N6ikvk4542l/O3Vd/nIsPbjJMr9Zyt0ju6sHVHKEpvJc+dfgI1QIrn65LFstXnfdnnZ2fd85a0VXBxNi51/sY/Xwy72/pLyIErfhcd/mBeXrADSq1vSzh8Hkfx0pOX3pMN354f3Pt+uAf+kMU18a+Ie7dKw9eZ9U39f44/9OGddchWb9+3D5pv24U9/DwMeO/u9Fvp7644s/kbT9PgusLBxMN15hMF0c4Bvuvu0vGMUJHqJ3paP3fkHL/U9nV04u5OXlbggVfIiV+pndTVPk8cVysdqfqdK6hVBohRm9vrAgQOHjh07tmZpmDlzJuPGjavZ5/cWysfsKC+z0ZvzcdasWSxduvQNd985f19vCxJPAzsA82qYjCag504EUz+Uj9lRXmajN+fjCOBNd98vf0evChIiIpKtnjCYTkREakRBQkREUilIiIhIKgWJTpjZYWb2JzN7w8zczM7M229mdqmZLTSz1WbWamZ75x2zuZlda2ZvmdnK6HwdehH0ZmZ2gZk9aWbLzOxNM7vLzPbJO0Z52Qkz+6qZPRPl4zIze9zMTkjsVx52g5ldGP1//yyxTXmJgkQptiaMy/g6sLrA/vOAc4FzgP2BJcADZpZcrf5q4FPA54BDgQHA3Wa2SQXTXW8mANcTZu49ElgHTDGz5JqLysvOvQ58B/gI0Aw8BPzRzOIh4srDLjKzg4CzgGfydikvAdxdjxIfwArgzMRrI3SJuyixbQtgOXB29Hog8D5wWuKYDwEbgGNr/Z1qmJdbE+bkOkl5WXZevgOcrTzsVt4NBP5BuHFpBX4WbVdeRg+VJMqzK2HiwY1rXbj7amAaubUuxgGb5h3zGvA8dbgeRhX1J5Rk4zm8lZddZGabmNkphID7GMrD7vg5YfGyh/K2Ky8jdT/BX52Lp3AstNbF0MQx6wmz2eYfU3frYVTRNcAs4PHotfKyRGY2mpBv/Qil20+4+2wziy9MysMSmNlZhEFkpxfYrb/HiIJENkpa66Ibx/RKZnYlcAhwiLuvz9utvOzcXGAssA2hPvxmM5uQ2K887ISZjQIuAw519/eLHNrweanqpvLEK9wXW+tiEWFlve2LHNMwzOwqQiPfke7+cmKX8rJE7v6+u89z96c8rPE+C/gmysOuOJiQB3PMbJ2ZrQMOB74SPY+ni234vFSQKM8rhD+UifEGM+tH6OXwWLRpJvBB3jE7A3smjmkIZnYNcCohQLyQt1t52X19gM1RHnbFH4HRhBJZ/HgKuC16/iLKy6DWLef1/iA0CsZ/RKuAS6Lnw6L93wGWAZ8E9iH8kS0E+ifOcQPwBnA0sB/wMOHub5Naf78q5uN1UT4dSbg7ix9bJ45RXnaej1cQLlTDCRe5ywm9aY5XHpadt61EvZuUl4l8qXUC6v1B6N/vBR43RfsNuJTQXW4NMBXYJ+8c/YBrCUXYVcBdwIdq/d2qnI+F8tCBSxPHKC87z8ebgAXAWkK//SkkulsqD8vK2/wgobx01yywIiKSTm0SIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWChDQ0M5tkZp+PnreaWXOFPudeM9umwPZLzezfoudnmtlOiX3zzSx/yodC597PzH7RzXTdZmYju/NeaQwKEtLQ3H2yu/+6Cp/zMXd/r5PDzgR26uSYQi4kDOjqjhsIi+uIFKQgIb2KmW1lZveY2d/NbI6ZnRxtn29m/2FmM6LHiGj7xjv5xDn6mNnNZvaDaM2GH0dLrz5jZmcX+MzzzOxr0fOrzOyh6PlRZvbbxOdvHz2/yMzmmtkUYFS07dOEleZuMbNZZrZFdPpzzOxvZjbbzD5c4LP7A2Pc/e+Fvk+UB8PT8gV4BDjazDQjtBSkICG9zXHAQnff1933Ae5L7Fvm7gcAPyMsO1lIX+AW4EV3vxj4ErDU3fcnLGF5lpntmveeaYT5lCBc6Lc2s00J06E/kjzQzMYBpxDm+flkdE7c/Q7CBHOnuftYDwvcALzl7h8h3PG3C2aJz5uTmhs5BfPF3TcA84B9SziHNCAFCeltZhPujP/DzA5196WJfb9L/Dw45f03AnPc/YfR62OAz5vZLGA6sB2QX4c/ExgX3dWvJSwI1EwIHI/kHXso8Ad3X+Xuy4A/dfJ9/jfxGcML7G8C3uzkHFA8X5bQvWouaQAKEtKruPuLhGUlZwOXm9klyd0pz5MeA46IpoWGMMnbOdHd/Vh339Xd70++wd0/AOYDX4je/whwBLA7YSnLDsnswldaG/1cT+FFwlYTJpmLraP9/3W/KI3F8qVfdB6RDhQkpFeJegetcvffAj8BPpLYfXLi5+P57438ErgX+H1UT/8X4MtR9RFmtoeZbVXgfdMI1UHTCEFiEjDLO86gOQ34hJltEZU8TkrsW05Y+7srnicswRmbT/SdzewjhLWaO8uXPYBnu/i50iDUWCW9zWjgx2a2gbAgzJcT+zY3s+mEm6PPpZ3A3a80s4HAb4DTCNU8fzMzI1TtfLzA2x4BLgIed/eVZraGjlVNuPvfzOx2wpoDC/KOuQmYbGarSa8Oyz/fC2Y20Mz6u/ty4E5y1WNPEhbPgZR8MbPBwGp3byvl86TxaKpwaQhmNh9odvf8Ret7PDP7JrDc3bs8ViJ67zJ3/2X2KZPeQNVNIj3fDeTaLrrqPeDmDNMivYxKEiIikkolCRERSaUgISIiqRQkREQklYKEiIikUpAQEZFUChIiIpJKQUJERFIpSIiISCoFCRERSaUgISIiqRQkREQklYKEiIikUpAQEZFUChIiIpJKQUJERFIpSIiISCoFCRERSaUgISIiqRQkREQklYKEiIikUpAQEZFUChIiIpJKQUJERFL1rXUCsmRmTwM7APNqnRYRkR5kBPCmu++Xv6NXBQlgh4EDBw4dO3bs0FonRESkp5g1axZLly4tuK+3BYl5Y8eOHdra2lrrdIiI9BgTJkxg6tSpBWtg1CYhIiKpFCRERCSVgoSIiKRSkBARkVQKEiIikqq39W4SkV5s+Pn3ZH7O+VeckPk5exOVJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSqXeTiNRcJXotSTZUkhARkVQKEiIikkpBQkREUilIiIhIKgUJERFJpSAhIiKpFCRERCRVTYOEmV1gZk+a2TIze9PM7jKzffKOMTO71MwWmtlqM2s1s71rlWYRkUZS65LEBOB6YDxwJLAOmGJm2yaOOQ84FzgH2B9YAjxgZv2rm1QRkcZT0xHX7n5s8rWZnQ4sBT4K3GVmBnwDuMLd74yOOYMQKE4FbqxuikVEGkutSxL5+hPS9G70eldgCHB/fIC7rwamEUofIiJSQfU2d9M1wCzg8ej1kOjn4rzjFgNDC52gra2N5ubm1A9oaWmhpaWlzGSKiDSGugkSZnYlcAhwiLuvz9vt+YcX2AZAU1MTra2t2SdQRKQB1UV1k5ldBXwOONLdX07sWhT9HJL3lh3pWLoQEZGM1TxImNk1hEboI939hbzdrxACxcTE8f2AQ4HHqpZIEZEGVdPqJjO7Djgd+DjwrpnFJYYV7r7C3d3MrgYuMrMXgBeBi4EVwK01SbSISAOpdZvEV6KfD+Zt/x5wafT8R8AWwHXAIGA6cIy7L69GAkVEGlmtx0lYCcc4IWBcWun0iIhIezVvkxARkfqlICEiIqkUJEREJJWChIiIpFKQEBGRVLXuAisidWT4+feUdNz8K06ocEqkXqgkISIiqRQkREQklYKEiIikUpAQEZFUChIiIpJKQUJERFIpSIiISCoFCRERSaUgISIiqRQkREQklYKEiIik0txNIiIlaNR5rWpekjCzw8zsT2b2hpm5mZ2Zt/+maHvy8USNkisi0lBqHiSArYE5wNeB1SnHTAGaEo+PVSdpIiKNrebVTe5+L3AvhFJDymFr3X1R1RIlIiJAfZQkSnGImS0xsxfN7L/MbMdaJ0hEpBHUvCRRgvuA/wVeAYYDPwAeMrNx7r62lgkTEent6j5IuPttiZezzWwmsAA4gRA82mlra6O5uTn1fC0tLbS0tGSeThHpqNQeQVK/6j5I5HP3hWb2OjCy0P6mpiZaW1urmygRkV6qp7RJbGRm2wNDgbZap0VEpLereUnCzLYGRkQv+wDDzGws8E70uBS4kxAUhgOXA0uAP1Q7rSIijaYeShLNwNPRYwvge9Hz7wPrgdHA/wEvAjcDc4GD3X15TVIrItJAal6ScPdWwIoccmyVkiIiInlqHiSkZ2jUeWtEGl1FqpvMrI+ZDajEuUVEpHoyCxJmdquZDTCzrYDngLlm9u2szi8iItWXZUliL3dfBnycMBfTMOD0DM8vIiJVlmWQ2NTMNiUEif9z9w8Az/D8IiJSZVkGiRuB+cBWwDQz2wVYluH5RUSkyjLr3eTuPwV+mti0wMyOyOr8IiJSfWUHCTP7VieHXFnuZ4iISG1kUZLon8E5RESkDpUdJNz9e1kkRERE6k8W1U0/Lbbf3b9W7meIiEhtZFHdNDODc0iD0TQf2ahVPmoxocaRRXXTzcnXZraVu68s97wiIlJ7WU7LcbCZPQc8H73e18yuz+r8IiJSfVkOpruaMK332wDu/nfgsAzPLyIiVZbpLLDu/lrepvVZnl9ERKory/UkXjOz8YCb2WbA14iqnkREpGfKMkhMAq4BhgKvA/cDX8nw/CIimVNPreKyDBL/Cfyru78LYGaDom1fzPAzRESkirJskxgTBwiA6Pl+nb3JzA4zsz+Z2Rtm5mZ2Zt5+M7NLzWyhma02s1Yz2zvDdIuISIosg0SfqPQAgJltS2klla2BOcDXgdUF9p8HnAucA+wPLAEeMDPNGSUiUmFZVzc9ZmZ3EBYb+izww87e5O73Elayw8xuSu4zMwO+AVzh7ndG284gBIpTCWtYiIhIhWRWknD3XwOfAhYDbwKfdPfflHnaXYEhhEbw+HNWA9OA8WWeW0REOpFlSQJ3fw54LsNTDol+Ls7bvpjQi6qDtrY2mpubU0/Y0tJCS0tLNqmTDnpCTxHNGyVSukyDRAXlr5VtBbYB0NTURGtra8UTJCLSCDIdcV0Bi6KfQ/K270jH0oWIiGSs3oPEK4RAMTHeYGb9gEOBx2qVKBGRRlHz6iYz2xoYEb3sAwwzs7HAO+7+qpldDVxkZi8ALwIXAyuAW2uSYBGRBlLzIAE0Aw8nXn8vetwMnAn8CNgCuA4YBEwHjnH35dVNpohI46l5kHD3VkJDdNp+By6NHtJgatlbqjf1guoJvc6kPtV7m4SIiNSQgoSIiKRSkBARkVQKEiIikkpBQkREUtW8d5NIT9ebekFJ+Xrb34NKEiIikkpBQkREUilIiIhIKgUJERFJpSAhIiKp1LtJpEp6W68XaQwqSYiISCoFCRERSaUgISIiqRQkREQklYKEiIikUu+mBqcVy+qPfieNIevfc6V6xdV9ScLMLjUzz3ssqnW6REQaQU8pScwFJiRer69ROkREGkpPCRLr3F2lBxGRKqv76qbIbmb2hpm9Yma3mdlutU6QiEgj6AlBYjpwJnA8cBYwBHjMzLarZaJERBpB3Vc3ufufk6/N7AngZeAM4Mr849va2mhubk49X0tLCy0tLVkns2rU80VEqqnug0Q+d19hZs8CIwvtb2pqorW1tbqJEhHppXpCdVM7ZtYP+DDQVuu0iIj0dnUfJMzsJ2Z2uJntamYHAncAWwE31zhpIiK9Xk+obtoZ+B2wPfAm8ARwkLsvqGmqREQaQN0HCXc/pdZpEBFpVHUfJBqBeiyJSL2q+zYJERGpHQUJERFJpSAhIiKpFCRERCSVgoSIiKRS76ZuKLU3UqVWihIRqRaVJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWCRMaWz7qv1knoFZSP2VFeZqNR81FBImMr/t6Yf0hZUz5mR3mZjUbNR83dlJD1CnFacU5EeroeU5Iws6+Y2StmtsbMZprZobVOk4hIb9cjgoSZnQxcA1wG7Ac8BvzZzIbVNGEiIr1cjwgSwLeAm9z9v9z9eXc/B2gDvlzjdImI9Gp1HyTMbDNgHHB/3q77gfHVT5GISOMwd691Gooys52AN4DD3X1aYvslwGnuPiqx7fWBAwcOHTt2bLc+64mX3y43uby/aB6bDRlR9nkanfIxO8rLbNR7Ph6023bdfu+sWbNYunTpG+6+c/6+nhQkDnP3RxLb/x34nLt/OLHtaWA4xXttLSRUVVVKU4XP3yiUj9lRXmajN+fjCOBNd98vf0dP6AL7FrAeGJK3fUdgcXJDoS8oIiLdV/dtEu7+PjATmJi3ayKhl5OIiFRITyhJAFwJ/MbMZgCPApOAnYDJNU2ViEgv1yOChLvfbmbbARcT6gXnAB9z9wW1TZmISO9W99VNMXe/3t2Hu/vm7j4u2dOpkszsMDP7k5m9YWZuZmfm7Tczu9TMFprZajNrNbO9847Z3MyuNbO3zGxldL4OvQh6MzO7wMyeNLNlZvammd1lZvvkHaO87ISZfdXMnonycZmZPW5mJyT2Kw+7wcwujP6/f5bYprykBwWJGtqaUHL5OrC6wP7zgHOBc4D9gSXAA2bWP3HM1cCngM8BhwIDgLvNbJMKprveTACuJ4xtORJYB0wxs20TxygvO/c68B3gI0Az8BDwRzMbE+1XHnaRmR0EnAU8k7dLeQng7nqU+ABWAGcmXhuhS9xFiW1bAMuBs6PXA4H3CWM64mM+BGwAjq31d6phXm5N6LV2kvKy7Lx8BzhbeditvBsI/INw49IK/CzarryMHipJlGdXQtfcjaPB3X01MI3caPBxwKZ5x7wGPE9jjxjvTyjJvhu9Vl52kZltYmanEALuYygPu+PnwB3u/lDeduVlpEc0XNexeOzG4rzti4GhiWPWE8Z75B+TP/ajkVwDzAIej14rL0tkZqMJ+daPULr9hLvPNrP4wqQ8LIGZnUUYRHZ6gd36e4woSGQjf9i6FdiWr5RjeiUzuxI4BDjE3dfn7VZedm4uMBbYhlAffrOZTUjsVx52wsxGEWaVPtTDWKw0DZ+Xqm4qz6LoZ7HR4IuATYDtixzTMMzsKkIj35Hu/nJil/KyRO7+vrvPc/en3P0CQonsmygPu+JgQh7MMbN1ZrYOOBz4SvQ8nsit4fNSQaI8rxD+UDaOBjezfoReDvFo8JnAB3nH7AzsSYONGDeza4BTCQHihbzdysvu6wNsjvKwK/4IjCaUyOLHU8Bt0fMXUV4GtW45r/cHoVEw/iNaBVwSPR8W7f8OsAz4JLAP4Y9sIdA/cY4bCJMUHk1YNOlhwt3fJrX+flXMx+uifDqScHcWP7ZOHKO87DwfryBcqIYTLnKXE3rTHK88LDtvW4l6NykvE/lS6wTU+4PQv98LPG6K9htwKaG73BpgKrBP3jn6AdcSirCrgLuAD9X6u1U5HwvloQOXJo5RXnaejzcBC4C1hH77U0h0t1QelpW3+UFCeele/1OFi4hI7ahNQkREUilIiIhIKgUJERFJpSAhIiKpFCRERCSVgoSIiKRSkJCGZmaTzOzz0fNWM2uu0Ofca2bbFNh+qZn9W/T8TDPbKbFvvpnlT/lQ6Nz7mdkvupmu28xsZHfeK41BQUIamrtPdvdfV+FzPubu73Vy2JmEtdu76kLCgK7uuIGwuI5IQQoS0quY2VZmdo+Z/d3M5pjZydH2+Wb2H2Y2I3qMiLZvvJNPnKOPmd1sZj+I1mz4cbT06jNmdnaBzzzPzL4WPb/KzB6Knh9lZr9NfP720fOLzGyumU0BRkXbPk1Yae4WM5tlZltEpz/HzP5mZrPN7MMFPrs/MMbd/17o+0R5MDwtX4BHgKPNTDNCS0EKEtLbHAcsdPd93X0f4L7EvmXufgDwM8Kyk4X0BW4BXnT3i4EvAUvdfX/CEpZnmdmuee+ZRphPCcKFfmsz25QwHfojyQPNbBxwCmGen09G58Td7yBMMHeau4/1sMANwFvu/hHCHX+7YJb4vDmpuZFTMF/cfQMwD9i3hHNIA1KQkN5mNuHO+D/M7FB3X5rY97vEz4NT3n8jMMfdfxi9Pgb4vJnNAqYD2wH5dfgzgXHRXf1awoJAzYTA8UjesYcCf3D3Ve6+DPhTJ9/nfxOfMbzA/ibgzU7OAcXzZQndq+aSBqAgIb2Ku79IWFZyNnC5mV2S3J3yPOkx4IhoWmgIk7ydE93dj3X3Xd39/uQb3P0DYD7whej9jwBHALsTlrLskMwufKW10c/1FF4kbDVhkrnYOtr/X/eL0lgsX/pF5xHpQEFCepWod9Aqd/8t8BPgI4ndJyd+Pp7/3sgvgXuB30f19H8BvhxVH2Fme5jZVgXeN41QHTSNECQmAbO84wya04BPmNkWUcnjpMS+5YS1v7viecISnLH5RN/ZzD5CWKu5s3zZA3i2i58rDUKNVdLbjAZ+bGYbCAvCfDmxb3Mzm064Ofpc2gnc/UozGwj8BjiNUM3zNzMzQtXOxwu87RHgIuBxd19pZmvoWNWEu//NzG4nrDmwIO+Ym4DJZraa9Oqw/PO9YGYDzay/uy8H7iRXPfYkYfEcSMkXMxsMrHb3tlI+TxqPpgqXhmBm84Fmd89ftL7HM7Nv/Y9XowAAFURJREFUAsvdvctjJaL3LnP3X2afMukNVN0k0vPdQK7toqveA27OMC3Sy6gkISIiqVSSEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJFXfWicgS2b2NLADMK/WaRER6UFGAG+6+375O3pVkAB2GDhw4NCxY8cOrXVCRER6ilmzZrF06dKC+3pbkJg3duzYoa2trbVOh4hIjzFhwgSmTp1asAZGbRIiIpJKQUJERFIpSIiISCoFCRERSVXRIGFmh5nZn8zsDTNzMzszb7+Z2aVmttDMVptZq5ntnXfM5mZ2rZm9ZWYro/PtXMl0i4hIUOneTVsDc4BfR4985wHnAmcCc4FLgAfMbJS7L4+OuRr4Z+BzwNvAlcDdZjbO3ddnmdjh599T0nHzrzghy48VEalbFS1JuPu97n6hu98BbEjuMzMDvgFc4e53uvsc4AygP3BqdMxA4EvAt939AXf/G3A6MAY4upJpFxGR2rZJ7AoMAe6PN7j7amAaMD7aNA7YNO+Y14DnE8eIiEiF1DJIDIl+Ls7bvjixbwiwHniryDEiIlIh9TDi2vNeW4Ft+VKPaWtro7m5OfWNLS0ttLS0dCmBIiKNqpZBYlH0cwjwWmL7juRKF4uATYDtgTfzjplW6KRNTU1oWg4RkWzUsrrpFUIQmBhvMLN+wKHAY9GmmcAHecfsDOyZOEZERCqkoiUJM9uaMAUthIA0zMzGAu+4+6tmdjVwkZm9ALwIXAysAG4FcPelZvZL4MdmtoRcF9hngCmVTLuIiFS+uqkZeDjx+nvR42bC2IgfAVsA1wGDgOnAMYkxEgDfBNYBt0fHPgh8PusxEiIi0lFFg4S7txIamdP2O3Bp9Eg7Zg1wTvQQEZEq0txNIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJJWChIiIpFKQEBGRVAoSIiKSSkFCRERSKUiIiEgqBQkREUmlICEiIqkUJEREJFVNg4SZzTczL/C4J9p/U4F9T9QyzSIijaRvjT9/f2CTxOsmYCbwP4ltU4DTE6/fr0K6RESEGgcJd38z+drMvgQsA36f2LzW3RdVNWEiIgLUUZuEmRnwJeC37r4qsesQM1tiZi+a2X+Z2Y41SqKISMOpdXVT0kRgV+AXiW33Af8LvAIMB34APGRm49x9bdVTKCLSYOopSJwFPOnus+IN7n5bYv9sM5sJLABOIASPDtra2mhubk79kJaWFlpaWrJJsYhIL1cXQSKqQvpn4KvFjnP3hWb2OjAy7ZimpiZaW1uzTaCISIOqlzaJM4G1wG3FDjKz7YGhQFsV0iQi0vBqHiSiBuv/D7jN3Zcntm9tZj8xs4PNbLiZTQDuApYAf6hNakVEGks9VDdNIFQf/Uve9vXAaODzwDaE0sPDwGeTwURERCqn5kHC3R8GrMD21cCx1U+RiIjEal7dJCIi9UtBQkREUilIiIhIKgUJERFJpSAhIiKpFCRERCSVgoSIiKRSkBARkVQKEiIikkpBQkREUilIiIhIKgUJERFJ1eUgYWZ9zGxAJRIjIiL1paRZYM3sVmASYfrumcBAM7vS3X9cycTVq+Hn31PScfOvOKHCKRERqaxSSxJ7ufsy4OPAvcAw4PSKpUpEROpCqUFiUzPblBAk/s/dPwC8cskSEZF6UGqQuBGYD2wFTDOzXYBllUqUiIjUh5LaJNz9p8BPE5sWmNkRlUmSiIjUi6JBwsy+1cn7r8wwLSIiUmc6K0n0r0oqRESkLhUNEu7+vWolRERE6k9n1U0/Lbbf3b+WbXJERKSedFbdNLOSH25mlwL/nrd5sbsPifZbtL8FGARMB77q7s9WMl0iIhJ0Vt10c/K1mW3l7iszTsNcYELi9frE8/OAc4Ezo+MuAR4ws1HuvjzjdIiISJ6SxkmY2cFm9hzwfPR6XzO7PqM0rHP3RYnHm9FnGPAN4Ap3v9Pd5wBnEBrTT83os0VEpIhSB9NdDRwLvA3g7n8HDssoDbuZ2Rtm9oqZ3WZmu0XbdwWGAPfHB7r7amAaMD6jzxYRkSJKngXW3V/L27S+4IFdM51QlXQ8cBYhKDxmZttFzwEW571ncWKfiIhUUEkjroHXzGw84Ga2GfA1oqqncrj7n5OvzewJ4GVCtdIT8WF5b7MC2zZqa2ujubk59TNbWlpoaWnpVnpFRBpNqUFiEnANMBT4/9u7/yCryvuO4+8PSf3RgWKUIWx0EmJDpFYcKBurAdzViUVDmqnJtGrTTEwnIRYxxcQmWjVl2jSxHaJhqjbRjAPtJNEREytaqTF2MYJooEEliDQVTcGlkJbhR4NQ8Ns/nnPxcNmzexd299yz9/OaOcPe85x7zve7l93vnnOe8zybSZeA5gx0MBGxR9JPgQnAg9nqcUD+LGYsR55dHNLW1kZXV9dAh2Zm1pIavdz0NWBuRLw9IsYC1wALBjoYSScAE4FuYBOwFbiorn0GsHKgj21mZkdq9Ezi7IjYUXsRETskTTnWg0taACwFfk46Q7iZNNLs4ogISV8HbpS0AdgI3ATsAb5zrMc2M7O+NVokRkh6W61QSDq5H+/tzWnAd4ExwHbSfYhzI+LVrP1vgROBO3jzYbrf8TMSZmZDo9Ff9F8j9TpaQrpp/AfAXx/rwSPi8j7aA5ifLZXjaU7NrOoanU/iHyStBi4k9S76SESsH9TIzMysdA1fMsqKgguDmVkLafhhOjMzaz0uEmZmVshFwszMCrlImJlZIRcJMzMr5CJhZmaFXCTMzKyQi4SZmRVykTAzs0IuEmZmVshFwszMCrlImJlZIRcJMzMr5CJhZmaFXCTMzKyQi4SZmRVykTAzs0IuEmZmVshFwszMCrlImJlZoVKLhKQbJP1Y0i5J2yUtlXRW3TaLJEXdsqqsmM3MWknZZxKdwJ3A+4ELgQPA45JOrtvucaAtt3xwCGM0M2tZby3z4BExM/9a0seBncA0YGmuaV9EbB3K2MzMrPwziXqjSDHtqFs/XdI2SRsl3S1pbAmxmZm1nFLPJHqwEFgLPJ1btwz4HrAJGA98GXhC0tSI2DfkEZqZtZCmKRKSbgWmA9Mj4mBtfUTcm9vsBUlrgFeBWaTicZju7m7a29sLjzN79mxmz549YHGbmQ1nTVEkJN0GXA5cEBEv97ZtRLwmaTMwoaf2trY2urq6Bj5IM7MWVHqRkLSQVCA6I2JDA9uPAU4Fugc7NjOzVlf2cxJ3AJ8ErgB2SBqXLSOz9pGSFkg6T9J4SZ2kXk/bgO+XFriZWYsou3fTHFKPph+Szgxqy3VZ+0FgEvBPwEZgMfAScF5E7B7yaM3MWkzZz0moj/a9wMzetmk1469/pKHtXrll1iBHYmatoOwzCTMza2IuEmZmVshFwszMCpXeBdYav89gZjbUfCZhZmaFXCTMzKyQi4SZmRVykTAzs0IuEmZmVshFwszMCrlImJlZIT8nMUx5jCczGwg+kzAzs0IuEmZmVshFwszMCrlImJlZId+4bnG+wW1mvfGZhJmZFXKRMDOzQr7clLN77TJGTb647DCO2WDkUcZlqbvuuovZs2cP2P7KNFxycR7NZSjyUEQM6gGGkqSujo6Ojq6urqN6//FtE2j7xNcHNqgSdC+e1/R5NFJM2tvbWb169RBEM/iGSy7Oo7kMVB6dnZ0sX758eUR01rf5cpOZmRWqzOUmSXOAPwPagJ8C8yLiR+VGZUerkctX3Vt2DkEkZtabShQJSZcBC4E5wFPZv49KOjMifl5qcNbyXPBsOKtEkQA+ByyKiLuz19dIuhj4E+CG8sKyZtKKz3wMdM597a97y07GX//IsPoeWu+avkhIOg6YCiyoa3oMeP/QR2RDqdFfglU4dpm5DLSyCnJ/vodlFbKhLNy1ot2f/fVX0/dukvQOYAvQERFP5tZ/CfhYRJyRW7d59OjRp06ePPmojvXkymc4btx7jjXk0u3f+jPn0YtzTz+loe1WvfzfA3bMsj+Tgcq5v3k0etxG9ecz6e3Ya9asYerUqQMR0hEajXEgPpP853Es3+u1a9eyc+fOLRFxWn1blYrE+fkb1ZL+ArgiIibm1v0EGE/vZ0ivAd0FbW29tFWJ82g+wyUX59FcBiqP9wDbI2JKfUPTX24CfgEcBMbVrR8L/Fd+RU8JmpnZ0Wv65yQiYj+wBriorukiYOXQR2Rm1jqqcCYBcCvwj5KeBVYAVwHvAL5RalRmZsNcJYpERNwn6RTgJtI1uHXAByPi1XIjMzMb3pr+clNNRNwZEeMj4viImJrv6XSsJM2RtEnS65LWSJoxUPseLJLOl/SQpC2SQtKVde2SNF/Sa5L2SuqS9JslhVtI0g2Sfixpl6TtkpZKOqtum6bPRdLVkp7P8tgl6WlJs3LtTZ9DTyT9efb/6/bcuqbPJYsv6patufamzyFPUpukxdnPyOuS1kvqyLUPWj6VKRKDJfc091eAKaT7HI9KemepgfVtJOmM6k+BvT20fwH4PHAN8D5gG/ADSaOGLMLGdAJ3kp55uRA4ADwu6eTcNlXIZTPwReC3gHbgCeBBSWdn7VXI4TCSzgU+DTxf11SVXF4iXXmoLZNybVXJAUknkS6zC5gF/AYp7m25zQYvn4ho6QV4Bri7bt2/A18tO7Z+5LAHuDL3WqRucTfm1p0I7AY+U3a8feQyktSb7XeHQS7/A3ymijkAo4H/IBXuLuD2Kn0ewHxgXUFbJXLIxfYVYEUv7YOaT0ufSeSe5n6srqnqT3O/m9Rl+FBeEbEXeJLmz2sU6Qx3R/a6crlIeouky0kFbyUVzAG4C1gSEU/Ura9SLqdnl2M3SbpX0unZ+irlAPB7wDOS7pO0TdJaSXMlKWsf1HxaukgAY4C3UPe8Rfa6/rmMKqnFXsW8FgJrgaez15XJRdIkSXuAfaSed5dGxAtUKAcASZ8mPVx1cw/NVcnlGeBK4BLSJbNxwMqsA0xVcqg5nTSo6cvATNLPyC3A1Vn7oOZTid5NQ6D+sXP1sK6KKpWXpFuB6cD0iDhY11yFXF4CJgMnAR8FFkvqzLU3fQ6SziBd3pgR6RmlIk2dS0Q8mn8taRXpl+wngFW1zere1lQ55IwAVkdEbTDTn0iaQCoSt+e2G5R8Wv1MouGnuSum1oujMnlJug24ArgwIl7ONVUml4jYHxE/i4jaD/Ra4FoqlANwHukMe52kA5IOAB3AnOzr2kBCVcjlkIjYQ5qHZgLV+jwg3W9YX7fuRaDWuWZQ82npIhHD92nuTaT/OIfyknQCMIMmzEvSQuAPSQViQ11zpXKpMwI4nmrl8CCpF9Dk3LIauDf7eiPVyeWQLMaJpF+4Vfo8IPVsOqNu3XuB2nNig5tP2Xfuy16Ay4D9wKdIXcsWknoLvavs2PqIeyRv/hD/EvhS9vU7s/YvAruAjwBnkX7IXwNGlR17XR53ZHFeSPpLqLaMzG3T9LmQrhHPIA0wOQn4KvAGcElVcuglty6y3k1VyYU0tUAH6abubwMPZzG/qyo55HJ5H/B/wI2ke0W/D+wErh6Kz6T0b0AzLKSbQq+QbjiuIY04W3pcfcTcSbreWL8sytpF6gbYDbwOLAfOKjvuHvLoKYcA5ue2afpcgEWkv+z2kfqoPw7MrFIOveRWXySaPpfcL8n9pFGkHwDOrFIOdfnMAp7LYt0IfJZsFO/Bzqfphwo3M7PytPQ9CTMz652LhJmZFXKRMDOzQi4SZmZWyEXCzMwKuUiYmVkhFwmrPEmflfSipG9L+rCk6wu2K2wbxNg6JfV7JE5JUyR96yiPuaeXtmWSTi1omyvpk0dzTBu+/JyEVZ6kDaQnmzf1ss1bI+LAIB2/cN+S5gN7ImJBP/d5P/DliHjuKOLZExEje1h/IrA8Is4peN+vkuYtmNLfY9rw5VFgrdIkfYM0lPJDku4hzUPRHhFzJS0iTfwzBfg3SS/Ute0izSI3DvhCRCyRNII0smYHaUycEcA9EbGk7rhdpHFxpmXH3kiag/040iB4HyNN/HIVcFDSH5FmDdtAGka8NjjbvIhYUbfvUcDZtQKRTVO5MGsO4HzSPCh/mR3rDNLcAXMi4o3cfsYAS0nF5hHSU/pdWdstwIdJMwE+FhHXRcQvJb0i6ZyIeLaxT8CGOxcJq7SIuErSxcAFEfEL1c31TRoI7QMRcbCHtjbS0OQTgYeAJaSxb8aTxl8aSxpt856Cw58UER0Akt4GnBsRIelTpKLz+ayIHTqTkPQd4LaIeCqbIvdfSGOG5bWTpqatuY40Ts8KSSNJwy4AnAOcSRoOZFkW+5LsOG/PcropIn6QbX8JaUrVk4FLgYlZvCfljrWaNAaVi4QBLhI2/N0fR85NUfNg9pf3+uyXKqSicX+2fqukf+1l3/flvj4NuE9SG+lsoujS1weAM9+cVIxfkzQqInbntmkDtuderwBulfRt4HsRsTl7/7ORDasu6btZ7EuAXwF+SCosy3P7mUYqOG+QCs23JD1CGvyuZhupaJoBvnFtw9//9tK2L/e16v7t777/jjQI3iTSvNYnFLxnBHBeREzOllPrCgTA3vz7I+IW0ijFJwKrJNV+idffUKy9PkAaqHJmrSGbuvM/I815cYB0FvIAaWrMZbl9nJAd3wxwkTCr9xTwUUkjsrOLzgbfN5o02iik2c9qdpPm7a55DJhbeyFpcg/7epE0JHRtm1+PiBci4m9Il4NqReIcSe/O7qNclsUOqVj8MTAx15vrErJikF2yGh0R/wzMIw0xX/NeDr/UZS3ORcLscA8Am0m/KL9Jmit5ZwPvmw/cL+lHpBkPa5YCl2aT188gDfHcLul5SetJN7YPE2nipdHZDWyAeZLWSXqO9Fd+bWrOp0nzWKwjXd76fm4fB4HLgQskzQEu5s0zhlHAw5KeJw0pfW3u8NNIw5ybAe4Ca3YESSMjYo+kU0g3cKdFxNa+3jfAMVwL7I6IHp+VyObOvi4iPtTAvo4ndW1t72O7KcDnIuLjRxGyDVO+cW12pIezHj/HAX811AUi8/ekGciOWUTsI/WY6ssY4OaBOKYNHz6TMDOzQr4nYWZmhVwkzMyskIuEmZkVcpEwM7NCLhJmZlbIRcLMzAr9P4B2esgqwwWCAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x720 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "fig = plt.figure(figsize=(6,10))\n",
    "fig.subplots_adjust(hspace=0.4)\n",
    "ax = fig.add_subplot(311)\n",
    "#ax.plot(spikeWidth[:200],firingRate[:200],'o')\n",
    "ax.plot(spikeWidth[:200],firingRate[:200],'o')\n",
    "ax.plot(spikeWidth[200:],firingRate[200:],'x')\n",
    "ax.set_ylabel('firing rate (spk/s)')\n",
    "ax.set_xlabel('spike width (us)')\n",
    "\n",
    "ax1 = fig.add_subplot(312)\n",
    "ax1.hist(spikeWidth,bins=30)\n",
    "ax1.set_xlabel('spike width (us)')\n",
    "ax1.set_ylabel('cells')\n",
    "\n",
    "ax2 = fig.add_subplot(313)\n",
    "ax2.hist(firingRate,bins=30)\n",
    "ax2.set_xlabel('firing rate (spk/s)')\n",
    "ax2.set_ylabel('cells')\n",
    "\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 3. Implementation of the Support Vector Machine classification\n",
    "\n",
    "The `svm` package from the scikit-learn library (`sklearn`) is used to implement the support vector machine classifier below. \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([299.09546636,   6.41904046])"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.mean(Xdata,axis=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7d9cec179100>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAFOCAYAAACPJcp7AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOydd3gc5bW43yPZcu9NLmAZ2+BKAJsaNwgdbmKTSwqBUEPaDZeEBEgI/CDV4SYkkJAESGgh3OReAlxCMTXGJsbgBu4ONpa7bCRXybascn5/zIx2djW7Oyvtane1532efXb3m/bNrHTO9532iapiGIZhFB5F2e6AYRiGkR1MARiGYRQopgAMwzAKFFMAhmEYBYopAMMwjALFFIBhGEaBYgqgHSAi1SJyTLb74SEi3xORP7ify0RERaRDG137CyLySoLtc0XkugTbHxWRH2Wmd4lJ1rd8IMTznyEiW1t5jWS/4SwR2eL+X5wYsF1FZFRr+tBa0vEc0oEpgABEZIqILBCRfSKyW0T+KSIni8jpIlIjIj0CjlkmIv/hE3hLY7b3F5EjIlKe7v6qandV/TDd520pqvoTVc2KIFPVP6vquWH2FZGrROStTPcpW2RD0MU+/ywJ258D/+H+Xyxr42vnFaYAYhCRnsDzwK+BvsBQ4C6gVlXfBrYCn445ZgIwDvhvX3M3t93jMmBjBrtuFBBtNaPKJcQhjMwaDqzKcF/axfM3BdCcYwFU9b9VtUFVD6nqK6q63N3+GPDFmGO+CLygqlW+tj8BV8bs83i8iwaZSvxTXREZJSJvurOSShH5q2+/plGWa8K4X0ReEJEDIvKOiIz07XuuiKxzz/Nb95yBo3UROUVEFovIfhHZKSL3xPT1ehHZLiI7ROQm33F3isgTcc75aREpF5EJIlIkIreKyAYRqRKR/xGRvnGOe1NEPu1+nuJe/0L3+9ki8p77OWpULyLniMha935/A4jbPhb4PXC6ayrY67tcn3jPL6Bf/ysiFe7554nIeN+2ZL9FYN/iXOdOEXlKRJ4Qkf3AVe7v87aI7HV/g9+ISIm7/zz30Pfd+/us236xiLznHrNARI6Pd83WPv94fXC33SQiu9x+X53gunNF5Mci8k/gIOCZOkeKyLvus/s/EekrIp1EpBoodq+5IcR9TRHHXHSm+/0aEVkjIntE5GURGe7bV0Xk6yLyAfCBuGacePfi9ufnIrLZ/f/5vYh0CfG42wxTAM35F9AgIo+JyAUi0idm+5+AqSJyNIA4I5LLaC7cnwA+JyLFrrDpAbzTin79EHgF6AMMw5mhxOPzOLOWPsB64MduX/sDTwHfBfoB64AzEpznXuBeVe0JjAT+J2b7mcBo4FzgVhE5O9ENuP8cPwPOVtWVwA3ATGA6MATYA9wf5/A3gRnu52nAh+5x3vc3A67XH/gb8H2gP7AB+DiAqq4BvgK87ZoKevsODXx+cXgJ5xkMBJYCf47Znui3COxbAj6F8/v1dq/TAHzTPf504BPA19z7m+Ye8zH3/v4qIicBDwNfxvn9HwCeE5FOSa4LLXj+QX1wv5cCvXBm19cC9wf8n/m5Arge539ok9v2ReAanL+beuA+Va1V1e6+a8ZV3AAich7OrP3TqvoPEZkJfA+4BBgAzCd6Vg/O3+upODP+ZPfyM5wB5QnAKHefOxL1qc1RVXvFvICxwKM45p564DlgkG/7a8D33M/nAJVAR/d7GaBAB3e/84DZwG3A2UB5nGs2Hedrmwtc535+HHgQGBZwrAKj3M+PAn/wbbsQWOt+/iKOwPO2CbDFu0bAeefhCK/+cfo6xtd2N/BH9/OdwBMx+34bWO3vP7AG+ITv+2Cgzv8MfNs+ASx3P88BrgMWut/fBC5xP18FvOW734Ux97vV90yb9vXtE/f5hfi76e3ea6+Qv0XcvgWc+05gXpLr3wg8E/R34X7/HfDDmGPWAdND3FvKzz9OH2YAh4j+O98FnBbnunOBHwS0zfZ9HwccAYqDrhnn/+W7OMpkoq/9JeBa3/cinFnHcN9xZ4W5F/f3rAFG+radDmz0Hbs1zN9VJl82AwhAVdeo6lWqOgyYgDPK+JVvF78Z6ArgSVWtCzjV4zj/EJ/HmRG0hptx/qjeFZFVInJNgn0rfJ8PAt6oaAiOwAdAnb/ERJEI1+KMYNaKyCIRuThm+xbf503u+ePxHeB+VfVfbzjwjGuO2IujEBqAQQHHvw0cKyKDcEZUjwNHuSPpU3CUVSxB97slYL9Y4j2/KNzZ3WxxTFj7gXJ3U/8Q52pJ36K2i8ixIvK8a4LaD/wk5tqxDAdu8p63+8yPIvHv5tGS5x+PKlWt932P+4xdgp5L7N9eRxLfeyw3Av+jqit8bcOBe33PZjfO/9zQBH2Jdy8DgK7AEt/55rjtOYMpgCSo6lqckZzfofs0MNS1G15CfNv+34CLgA9VdVOcfTxq3PeuvrZSXz8qVPVLqjoEZwr/W0k9umIHjvkIcJxq/u+xqOoHqvp5HPPGz4CnRKSbb5ejfJ+PBrYnuPa5wPc9O7LLFuACVe3te3VW1W0BfTkILAH+E1ipqkeABcC3gA2qWhnnfpv66N6vv8+tLYV7GY5Z5mwcM0CZd6kQxybrWxCx/f0dsBYYrY6Z7ntJrr0F+HHM8+6qqrFmjuYXbtnzTxdBv1Ps314dzkw8LJcCM0XkRl/bFuDLMc+ni6ouSNKXICpxZgfjfefqpRETVU5gCiAGERnjOnWGud+PwhnBL/T2UdUaHFvsI8AmVV0cdC53v7NwpssJUdWPgG3A5e7I8hocu7vXr0u9PuHYyhVntJwKLwATRWSmOM7mr+NTMrGIyOUiMkBVGwHPSeq/5u0i0lUcx+fVwF+bnSTCKuB8HBvpJ9223wM/9hxtIjJARD6V4BxvAv9BxN48N+Z7LC8A40XkEvd+byD6fncCw8R1nLaAHkAtUIWjuH+SwrHJ+hb2+vuBahEZA3w1ZvtOIk5TgIeAr4jIqeLQTUQuEjesWRyn9aMJrpfq8w/qQ7q4XETGiUhX4AfAU6qayv/Ddhyz1g0i8jW37ffAd92/Z0Skl4hc2pLOuf8zDwG/FJGB7vmGun6HnMEUQHMO4Dh53hGRGhzBvxK4KWa/x3CmjHEjewBUdbGqJo1GcPkSjqmkChiPM8LyONntUzWOT+I/VTWlsFJ3lHYpjr2+Csd2uhhHiAVxPrDKvea9wOdU9bBv+5s4js3XgZ+ratwEIPf67wMXAw+JyAXuOZ8DXhGRAzjP+tQEp3gTR+jNi/M93v3Odu93NPBP3y5v4CimChFpyQj2cRzzwzYc/8bCxLun1LcwfBtnFnIAR9jEKuA7gcdcE8Rn3IHKl4Df4Awi1uOYKD2OStKHlJ5/UB+S31Jo/oQzM68AOuMo0JRQ1c04SuAWEblOVZ/Bmen+xTWprQQuaEUfb8F5xgvd870GHNeK86UdcR0SRgHiRjBtBb6gqv9I4bgynJyGjjH2TyNPcWdB7wPHx/FnGe0QmwEUGCJynoj0dkP/PJtx6JGr0T5R1SOqOtaEf2FhCqDwOB0n5rwS+Ddgpqoeym6XDMPIBmYCMgzDKFBsBmAYhlGg5E1BIxFZhpNEsT7bfTEMw8gjRgEfqWqz0th5owCAAb169Rp6wgknDE2+q2EYhgHw3nvvsW/fvsBt+aQA1p9wwglD586dm+1+GIZh5A0zZszgzTffDLScmA/AMAyjQDEFYBiGUaCYAjAMwyhQ8skHYBhGO6Kuro6tW7dy+PDh5DsboejcuTPDhg2jY8eOofY3BWAYRlbYunUrPXr0oKysDKcattEaVJWqqiq2bt3KiBEjQh1jJiDDMLLC4cOH6devnwn/NCEi9OvXL6UZlSkAwzCyhgn/9JLq8zQFYGQfVVjzd+c9TLthpIni4mJOOOEEJkyYwKWXXsrBgwfb7NqLFy/mhhuClzEoKyujsjLxEhVh9kmGKQAj+6x9Hv56Ocz5bkTYqzrf/3q5s90oeKpr6/nLu5uZ/dIa/vLuZqprW78URZcuXXjvvfdYuXIlJSUl/P73vw91XENDqovxNWfy5Mncd999rT5Pa8ioAhCRwSLymIh8JCKHRWS1iEz3bRcRuVNEtovIIRGZ6y3HZhQQYy6GU78K7/wuogTmfNf5fupXne1GQbOofDen/uQ1fvD8an7/5of84PnVnPqT11hUvjtt15g6dSrr16/n9ttv5957721qv+2227jvvvuYO3cuZ555JpdddhkTJ04EYObMmUyaNInx48fz4IMPNh3TvXt3brnlFiZNmsTZZ5/Nu+++y4wZMzjmmGN47rnnAJg7dy4XX+z8bVdVVXHuuedy4okn8uUvfxl/leZ410gHGVMAItIbZ3k5wVkYfSzwDWCXb7ebcZZa/AbOkoe7gFe9NUqNAkEEzv9pRAnc1Tsi/M//qbPdKFiqa+u56pF3qalt4OARZ+R98EgDNbUNbnvrZwL19fW89NJLTJw4kWuvvZbHHnsMgMbGRv7yl7/whS98AYB3332XH//4x6xevRqAhx9+mCVLlrB48WLuu+8+qqqqAKipqWHGjBksWbKEHj168P3vf59XX32VZ555hjvuuKPZ9e+66y6mTJnCsmXL+OQnP8nmzZubtsW7RjrIZBjozcAOVf2ir61pDVtxvBU3ArNV9W9u25U4SuAy4IEM9s3INTwl8M7vIm0m/A3g+fe3x3UDqcLzy7fz2ZOPbtG5Dx06xAknnAA4M4Brr72WkpIS+vXrx7Jly9i5cycnnngi/fr1A+CUU06JCrG87777eOaZZwDYsmULH3zwAf369aOkpITzzz8fgIkTJ9KpUyc6duzIxIkTKS8vb9aPefPm8fTTTwNw0UUX0adPn6TXSAeZVAAzgTki8lfgTGA78AfgfnXmNyOAUqBpIXFVPSQi84AzMAVQWHhmHz9zvmtKwKC8qqZp5B/LwSMNlFe23HHr+QBiue6663j00UepqKjgmmuuaWrv1q1b0+e5c+fy2muv8fbbb9O1a1dmzJjRFILZsWPHpoicoqIiOnXq1PS5vj54xhIUwZPoGukgkwrgGOBrwC+B2cAJwK/dbb/BEf4AO2OO2wkElnzesXEtk8eVxb3g9TOnc/0lM1rcYSNLqMLqZ2DjPBgxDcbNcr6/8zvYtcr5XmBK4J2Nu9lw9KXZ7kZGObZTPVXVtUn369+9E106FnGorrHZti4di+jfvSTUeeIRdOy0cy7ktu/fTn19Hb9+8BGqqmvZd+gIR+obm/bfUlFJtx69ONRYzPIly1m4cCH7Dh1p2u69HzxSj9RG32vs+U45/eM89Mhj3HTzd3ntlZfZs2cPu6trE16jUZXd1bVI5+j+19TW8+Q7ERPSrv3xn00mFUARsFhVvWHdMhEZDXwdRwF4xE7uJKANgMH9ezP3gVvS3lEjy1SsiBb+Is47OO19R8Hg47PbxzZm+JoH2r0CCMu54wZxz6v/CtwmIpw7flDar1lSUsKUadPp1as3xcXFgft84pxzeezhh5h22mRGjT6WSSef0uLrfefW27j+mi9y5pTTOGPKVIYddVTarxFExtYEFpFNwKuqep2v7Qrg96raTUSOwVmc/BRVXeTb5wWgUlWvjDnf3OknHTfdFEA7RNVRAqUTo0f68doLgIoXf8ob572W7W5klGM77WPEqGND7bts815u+MsyVJVDdY106ViEiHDf507kxKN7p71vjY2NnDXlNP74+JOMHDUq7efPJBvX/4t/1fZq+v6jr36WNcsWvqmqM2L3zeQM4J/AcTFtxwKb3M8bgQrgHGARgIh0BqYC38lgv4xcQyR4hB+v3Sg4Tjy6Ny/fOJVXVu1ky56DHNWnK+eOH0TXkvSLsHVr13DZpbO48OJP5Z3wT5VMKoBfAgtE5Dbgr8CJwA3A9wBUVUXkV8BtIrIW+BfwfaAaeDKD/TIMIw/pWtKBmSdmfkXY48aMZcmKtRm/Ti6QMQWgqotEZCbwE+B2YLP7/lvfbncDXYD7gT7AO8C5qnogU/0yDMMwHDJaDlpVXwBeSLBdgTvdl2EYhtGGWC0gw8g1Ni1gbcdx2e6FUQDYgjCGkWNUrHqT7e08AsjIDWwGYBhGwTKwV1dmnHEK008/mTOnnMa7C99u0XlOHH8sVSFLM69Y/j6vvjynRddJN6YADMPIbVQp+eCFwPUiAttToEuXLsxd8C5vvr2I2+/8IT+88/ZWdjY5K5e/z2uvmAIwDMNISsn6F+nxf1fT9R+3R60X0fUft9Pj/66mZP2LabnOgQP76d3bKcJWXV3NrIvP58wppzH11Em8+PzfAafK5+c+PZPpp5/MlFNO4pm//W/UOQ4dOsRnZv0bjz/yR2pqarjhq9dz9vSPc+bHT+XF5//OkSNHmP3jH/Ds355ixhmnNDu+rTEfgGEYOc2RURdy6KTr6bLUqYV/8Mwf0vUft9Nl6YMcOul6joy6sMXnPnToEDPOOIXDhw+zc2cFzzzvjMw7d+7M40/+Dz169qSqspLzPzGNCy66mDdee4XSwYP5y9+eBWD/vn1N56qpqeZLV1/BZz//BT572eX86M7bmTJ9Bvf97kH27d3LOTOmMP3Ms7j1tjt4b9lSfvaLX7XiqaQHUwCGYeQ2Ihw884cAdFn6YJMiOHTS9U57K8qEeCYggEXvLOTr11/LW+8uRVX50V138PY/36KoqIgd27eza9dOxo6bwP+77Vbuuv02zj3/Ak7/+JSmc13xuUv5jxu/xaWf/TwA/3jjdea8+AK/vc8R9LW1h9m2ZUuL+5oJzARkGEbu41MCHq0V/rGcfOpp7K6qorLyI576639TWVnJ6/PfZu6CdxkwcCC1hw8zavRoXp/3NmPHj+dHd97Of83+cdPxp5x2Oq+/+nLTal6qyiNP/IW5C95l7oJ3eX/Neo4dMyZt/U0HpgAMw8h9XJu/nyifQBr4YN06Ghoa6Nu3H/v372fAgAF07NiR+fPmssVdoWvHju106dqVz3zuMr5+wzdZ7ltL4Nbb7qBv335855vOQu9nfeJs/vDAb5sUwvL3nX279+hB9YHcKHZgCsAwjNzGFf6ezb/qpp1NPoHWKgHPBzDjjFO47qrL+c0Df6C4uJh//8zneG/pUj4x7Qye+utfGH2sU9dyzapVnDtjCjPOOIV7fv4zvnXzrVHn+/HPfk7t4cPc+f3vcdMt36Ouro5pp01myikn8dMf3gXAlKnTWbdurTmBDcMwklGy/sUm4e+Zffw+gfqjTufI6ItadO5d+4JXE+vXvz9z3nizWfvRw8s46+xzmrUvWxVZr+DXv3+o6fM9993fbN8+ffvy2pv/bEl3044pAMMwcpojoy7kwKcecaJ9PJu/qwTqjzq9VVFAhY4pAMMwchuR4BF+vHYjNOYDKARUYcfywEzKwHbDMAoCUwCFQMUKWPKws9C6L5OS1c847RUrsts/oyBRIFNL0hYqqhq8oHocTAEUAqUTnQXXN86LKIHVz0QWYi+dmO0eGgVIbWMxB/btNSWQJlSVA/v2UtsYvIh9EOYDKAREYNws5/PGec4LHOE/blbBLbhu5Abb6rrCrt10qvwI+wtsPYqjVLfVdQ19jCmAQsFTAp7wBxP+uUgBLQbTQBGb67pnuxsFjZmACgXP7OPH7xMwcoKKVW+y/az7st0No0CwGUAhEGvzHzcr8h1sJmAYBYopgEKgYkW08I/1CfQdBYOPz24fDcNoc0wBFAKlE2HSNc67L5OScbMc4W9RQIZRkJgCKAREgkf48doNwygIzAlsGIZRoJgCMAzDKFBMARiGYRQopgAMwzAKFFMAhmEYBYopAMMwjALFFIBhGEaBYgrAMAyjQDEFYBiGUaBkTAGIyJ0iojGvCt92cffZLiKHRGSuiIzPVH8MwzCMaDI9A1gHDPa9/EVnbgZuAr4BnAzsAl4VkR4Z7pNhGIZB5msB1atqRWyjiAhwIzBbVf/mtl2JowQuAx7IcL8MI/dY9qeCWQzGyA0yPQM4RkS2ichGEfmLiBzjto8ASoFXvB1V9RAwDzgjw30yjJykYsdWWwzGaFMyOQN4B7gKWAsMBL4PLHDt/KXuPjtjjtkJDI13wh2Ve5n8xbviXvD6mdO5/pIZLe+xYRhGAZExBaCqL/m/i8hC4EPgSmCht1vMYRLQ1sTg/r2Z+8At6eymYRhGwdJmYaCqWg2sAkYDnl+gNGa3gTSfFRiGYRgZoM0UgIh0BsYAO4CNOErgnJjtU4EFbdUnwzCMQiaTeQA/F5HpIjJCRE4FngK6AY+pqgK/Am4VkUtEZALwKFANPJmpPrULVGHHcuc9TLthGEYcMjkDGAb8N04uwNNALXCaqm5yt98N3APcDyzGyRM4V1UPZLBP+U/FCljyMKx+JiLsVZ3vSx52thuGYYQgk07gzyXZrsCd7ssIS+lEGDENNs5zvo+b5Qj/jfOcdlvg3TCMkNii8PmGiCP0wRH6niIYMc1pF8le3wzDyCusGFw+4lcCHib8DcNIEVMA+Yhn8/fj9wkYhmGEwExA+YYn/D2bv98HADYTMAwjNKYA8o2KFdHCP9Yn0HcUDD4+u300DCMvMAWQb5ROhEnXOO/eSN9TAn1HWRSQYRihMQWQb4gEj/DjtRuGYcTBnMCGYRgFiikAw8gFbDEYIwuYAjCMHKBix1ZWd5yQ7W4YBYYpAMPIEXpPvT7bXTAKDFMAhmEYBYopAMMwjALFFIBhGEaBYgrAMAyjQDEFYBiGUaCYAjAMwyhQTAEYhmEUKKYADMMwChRTAIZhGAWKKQDDMIwCxRSAYRhGgWIKwDAMo0AxBWAYhlGgmAIwDMMoUEwBGEa2UIUdy2Hp4zzd9TNR7cMqXne2G0YGMQVgGNmiYgUseZia7Wsjwl6Vk9bczbRlNzJs5xvZ7Z/R7jEFYBjZonQijJhGNw5yVv+9TcJ/zKYnWDv8crYOOivbPTTaOR2y3QHDKFhEYNwsajYuZsymJxiz6QkA1g6/nKVjb3a2G0YGSaoARGQyMBUYAhwCVgKvqeruDPfNMNo/IhyQ7nTTg01NJvyNtiKuCUhErhKRpcB3gS7AOmAXMAV4VUQeE5Gj26abhtFOUaWHVkc1nbTmbnMAG21CohlAN+DjqnooaKOInACMBjZnomOG0e5RhdXP0I2DTWYfzwcANhMwMk/cGYCq3h9P+Lvb31PV18NeSES+JyIqIr/xtYmI3Cki20XkkIjMFZHx4btvGHlMxQrYOI8aujYJ+6Vjb2bt8MsZs+kJiwIyMk7SKCARuVtEeopIRxF5XUQqReTyVC4iIqcBXwKWx2y6GbgJ+AZwMo6J6VUR6ZHK+Q0jLymdCJOu4YB0j4z0XSUw78RfWRSQkXHChIGeq6r7gYuBrcCxwHfCXkBEegF/Bq4F9vjaBbgRmK2qf1PVlcCVQA/gstB3YBj5iggMPr65mUeEraWfMPOPkXHCKICO7vuFwH+3IPrnQeApVY2dz44ASoFXvAbX5DQPOCPFaxiGYRgpEiYP4O8ishYnBPRrIjIAOBzm5CLyJWAUcEXA5lL3fWdM+05gaND5dlTuZfIX74p7vetnTuf6S2aE6ZqRi6g6dvHSidGj33jtuXoNw8gT4ioAERmsqjtU9VYR+RmwX1UbROQg8KlkJxaR44CfAFNV9UiCXWPj3SSgDYDB/Xsz94Fbkl3ayFfc0giMmAbjZjmC2I2UYeM8mHSNYzLJ9WsYRp6QyAT0sIgsFJHZwMdwBDOqWqOqFSHOfTrQH1gpIvUiUg9Mx5lF1ANV7n6lMccNpPmswCgE3NIIbJznCGS/YB4xzdmeD9cwjDwh7gxAVS8Qkc7ADGAW8HMR2QzMAeaoarL4/2eBxTFtjwAf4MwM/gVUAOcAiwDc600lBSez0Y5wSyMAjkDeOM/57B+t58M1DCNPSOgDUNXDuAIfQERGABcAvxGRUlU9JcGxe4G9/jYRqQF2uxE/iMivgNtcH8O/gO8D1cCTLb4jI7/xBLQnmCH9grktrmEYeUCYPIALvM+qulFVfwu8hFMSorXcDdwD3I8zWxiME3Z6IA3nNvIRzyTjxzPV5NM1DCMPCBMFdLuI1HphnCJyM3Cmqv4u1Yup6oyY7wrc6b6MQifWHj9uVuQ7pGeU3hbXSIVNC9ruWoYRQxgF8EngeRH5DnA+MMZtM4z04pZGiLLH++31fUelJwoo09dIgco18/hK30f4TPJdDSPtJFUAqlopIp8EXgOWAP/ujtwNI724pRGiYvE9Ad13VPqigDJ9jRSob1Q+c/JRbXpNw/BIlAdwgOh4/BLgGODfRURVtWemO2cUGF5phLDtuXoNw8gTEoWBWkE2wzCMdkyoJSFF5BKcqB8F5qvqsxntlWEYhpFxwoSB/hb4CrACZznIr4jI/ZnumGEYhpFZwswApgMTPMeviDyGowwMwzCMPCZMOeh1gH/t36NovrCLYRiGkWeEUQD9gDXuco1zgdXAABF5TkSey2jvjMJAFXYsb56JG6/dMIy0EMYEdEfGe2HkDtmol28lmg0jK4SZASzGifx5E9gB9AIWqOqbbpvRnvCEsb82jieMlzzsbE83VqLZMLJCmBnAPGCqiPQBXsdRCJ8FvpDJjhlZwi+MIbpWTqaEsZVoNoysEEYBiKoeFJFrgV+r6t0i8l6mO2ZkiWwJYyvRbBhtThgTkIjI6Tgj/hfctuLMdcnIOn4l4JFpYWwlmg2jzQmjAP4T+C7wjKquEpFjgH9ktltGVmlrYRxr87/ol819AoZhpJ0w1UDn4fgBvO8fAjdkslPtlmxE2KRKNurl51iJZsMoFOLOAETkQREJ9PiJSDcRuUZEzBGcCtmIsEmVeMLYG5FnKgpo0jXRysW7rle62TCMtJNoBvBbnNXAJuLUAPoI6AyMBnoCDwN/zngP2xPZiLBJlWzUyy/UEs1v3WOLwRhZJVE56PeAz4hId2Ayznq9h4A1qrqujfrXvsiHcMdCFcZZoGJ/LZ85zxaDMbJHGB9ANTA3810pECzc0TCMHCFMFJCRTizc0TCMHCHUgjBGmshGhI1hGEYcQisAEemmqjWZ7Ey7x8IdDcPIIcKsCHaGiKwG1rjfP+auEmakioU7GoaRQ4TxAfwSOA+oAlDV94FpmTs7Q/cAACAASURBVOxUu8WLpIk188RrNwzDyCChnMCquiWmqSEDfTGM1LCFZAyjVYRRAFtE5AxARaRERL6Naw4yjKySD5nVhpHDhFEAXwG+DgwFtgInAF/LZKcMIxS2kIxhtIowUUDHqWpUzR8R+Tjwz8x0yTBCkg+Z1YaRw4SZAfw6ZJthtD3ZWLvAMNoJcWcA7iIwZwADRORbvk09sQVhjFwhXma1KQHDSEqiGUAJ0B1HSfTwvfYD/57sxCLydRFZLiL73dfbInKRb7uIyJ0isl1EDonIXBEZ37rbMQoKW0jGMFpFomqgbwJvisijqrqpBefeCtwCfICjaK4EnhWRSaq6HLgZuAm4ClgH3AG8KiLHqeqBFlzPaCtyZWEby6w2jFYRxgl8UET+CxiPsx4AAKp6VqKDVPX/YppuE5GvAqeLyArgRmC2qv4NQESuBHYBlwEPhL8Fo83xwi/9gtc/Gp90TdsI3mysXWAY7YgwTuA/A2uBEcBdQDmwKJWLiEixiHwOx6S0wD1XKfCKt4+qHsJZevKMVM5tZIFcCb/M58xqdzEYw8gmYWYA/VT1jyLynz6z0JthTu6uJvY2zsyhGpilqivcxDKAnTGH7MTJNzByGQu/bDW2GIyRC4RRAHXu+w7XibsdGBby/OtwEsd6A58GHhORGb7tsV46CWhrYkflXiZ/8a64F7t+5nSuv2RG3O1GGrGFbQwj7wmjAH4kIr1wHLa/xgkD/WaYk6vqEWC9+3WxiJzsHvtjt60U8NcZGkjzWUETg/v3Zu4Dt4S5tJFpLPzSMPKehD4AESkGRqvqPlVdqapnquokVX2uFdfrBGwEKoBzfNfqDEzF8REYuUw+hV9awTjDiEtCBaCqDcAnW3JiEZktIlNFpExEJorIT4EZwJ9VVYFfAbeKyCUiMgF4FMdP8GRLrme0IfHCLz0lkEtF2NqyYJwpGyPPCGMCWiAivwH+CjStCKaqS5McVwo84b7vA5YDF6jqy+72u4EuwP1AH+Ad4FzLAcgD8in80h+xBNHLcKY7YilXwmMNIyRhFIAXsfMDX5sCyfIArkqyXYE73ZeRT3hhlmHbs0lbRiy1pbIxjDSQVAGo6plt0RHDyBhtFbFk4bFGnhFqRTDDyGviRSxlwiZv1UmNPMIUgNG+aeuIpbZUNobRSsL4AAwjf2nLgnGxysbvAwCbCRg5R1IFICKXBDTvA1ao6q70d8kw0khbRixZdVIjzwgzA7gWOB34h/t9BrAQOFZEfqCqf8pQ3wyj9Yg4Qj62THW89taQT+GxhkE4H0AjMFZVP62qnwbGAbXAqTj1/g0jt2mrZLB8rk5qFCRhFECZqvrr8+wCjlXV3UQKxRlG7pIr5asNI8cIYwKaLyLPA//rfv80ME9EugF7M9Yzw0gXFp9vGIGEmQF8HadOzwnAicDjwNdVtcaSxIy8IYfi8w+8/l/sK+rd5tc1jFjCZAIr8JT7Moz8JIfKV9fU1vO7Yx/klDa9qmE0J+kMwK3W+YGI7BOR/SJyQET2t0XnDCMt5GD56lNG9G3zaxpGLGF8AHcD/6aqazLdGcPICBafbxiBhFEAO034G3mNxecbRiBhFMBiEfkr8CxO/D8Aqvp0xnplGOkkn8pXG0YbEkYB9AQOAuf62hRoHwpANTgbNF67YRhGOyFMFNDVbdGRrGGrOBmGUaDEVQAicrOq3i0iv8YZ8UehqjdktGdtha3iZBhGgZJoBuA5fhe3RUeyhmWJGoZRoMRVAKr6dxEpBiao6nfasE9tT1stGWgEY34Yw8gKCRPBVLUBmNRGfcketopTdmmrap2GYUQRJgpomYg8h1MMrsZrbDdhoLaKU/YxP4xhZIUwCqAvUAWc5WtrP2GgliWafcwPYxhZwcJALUs0NzA/jGG0ORYGalmiiWkrB20OVes0jEIh0QxgtfvevsNAjcS0RaKc+WEMIyskUgCfBZ4HeqvqvW3UHyPXaAsHbQH5YQ68/l/Z7oJhNJFIAUwSkeHANSLyOBA1BHPXBDbaO23hoC0wP8zsY//HFoMxcoJECuD3wBzgGGAJ0QpA3XajEMi0g7aA/DA1tfW2GIyRM8RNBFPV+1R1LPCwqh6jqiN8LxP+hYQlyhlGuyRMGOhX26IjRo5iDlrDaLeESQQzCpkCctAaRqGRdFH4liIi3xWRRe5C8h+JyN9FZELMPiIid4rIdhE5JCJzRWR8pvpktADPQesf6XtKwHPcGs5Macfy5maxeO2GkQNkTAEAM4DfAmfglJGoB14TEb8H7GbgJuAbwMnALuBVEemRwX4ZqeA5YmPNPPHaCxUraGfkIRkzAanqef7vInIFsA/4OPB3ERHgRmC2qv7N3edKHCVwGfBApvpmGGnHCtrlJLV1DSzfto+q6lr6de/E8UN70aljcba7lTO0pQ+gB86MY4/7fQRQCrzi7aCqh0RkHs6swRSAkT9YQbuco7yyhsfeLkdVOdKglBQLL67YwZWnl1HWv1u2u5cTtKUCuBd4D3jb/V7qvu+M2W8nMDToBDsq9zL5i3fFvcD1M6dz/SUzWtfLXMcWT8ldrKBdzlBb18Bjb5dTW9/Y1HakQQHlsbfLufWCMXTqYDOBNlEAInIPMAWY4i4y4yfWOyYBbQAM7t+buQ/ckoEe5hG2iH3uYgXtcobl2/ahcRzvqsqKrfuYXGYJeRlXACLyS+BzwJmq+qFvU4X7Xgps8bUPpPmswPAwW3NuYvkSOUVVda074m/OkQalsvpIG/coN8moAhCRe3GE/wxVXRuzeSOOEjgHWOTu3xmYCrTvNYhbg9macxPLl8gpEgXdlhQL/buXtFlfcpmMKQARuR+4ApgJ7BERz+ZfrarVqqoi8ivgNhFZC/wL+D5QDTyZqX61C8zWnHsUWEG7XKa2roF3Pkxcq3LisF5t1JvcJpN5AF/Difx5Hdjhe33bt8/dwD3A/TjrDgwGzlXVAxnsV/7TnmvzxEucamyE1c8572H2b2ssXyJnWL5tH4nmAKcd088cwC4ZUwCqKnFed/r2UVW9U1UHq2pnVZ2uqisz1ad2Qayt+aJfRnwCmVICbZnlGi+havEf4MM3nHdLtDISkMj+72DK2COTMwAjE8SzNXtKIBOCsC2zXP1Obu96q5+BXauh51Dn3d9uzm8jhn7dO1FSHCzkzf4fjRWDyzeyYWtuy8ijRE7usTNhzbN56/yunDPbFoNpA44f2osXV+wgyAwkImb/92EzgHwjG7bm2FnGC99sPgvJxPX8jJsFRUXB7Xkg/AHqG9UWg2kDOnUs5srTy+jUoahpJlBSLHTqUOS2m/3fwxSAEY54QjkTwjeek7uxMfec31YFNCcp69+NWy8Yw8XHD2Ha6AFcfPwQbr1gjJWAiMFMQIZDshITgyY45hc/mchyTZRQVbUB9m/LrUQry8zOWTp1KLZs3ySYAjAckgmygeMcB2ymhW88J3dNpXP9geNyK9HKMrONPMYUgOGQSJDFCv9MCt94Tu7J18Ha52HMxbmVaGWZ2UYeYwrAcEgWfbNzZdtEHnnO7FiKimDcJ8Pv35ZYZraRp5gT2IiQKPom17Ncs+mMbc+Z2Ua7xhSAEaG1giybQjhbSzJmIzPbMNKEKQDDIZkga2xMLtyzuS5uvAziTDtjs5GZbRhpwnwAhkOycsbSwanFkyjcMZsRMcl8GJlaRc2qgBp5jCkAwyGZIBs0AbQ+sXDPdkRMPGdsJmP14zmhW+qcVmXYzjfYOuisZsoqsN0wWoGZgAyHZCUmvDIMycpBtGXGcCzxfBiDJmTHPNQChu18g2nLbuSkNXdHmdFOWnM305bdyLCdb2S3g0a7wmYARnjChDsGCeFFDzlx/EVF0fulcxH7ZEsyjp3pvOd4rP7WQWexdvjljNn0BABLx97MSWvuZsymJ1g7/HJnBmAYacIUgBGeZIueBwnhRQ85SWTzfwGjzoEhH4sct3EeTLoakNYrgjBLMuZDrL4IS8feDMCYTU80KYK1wy932nOtv0ZeYyYgIxxhwh2DhPDk65w6/ge2wbJHYdXTsMo9T9k0qFqfngghz4cRZI6adI1jBsqXWH2fEvAoZOFfW9fAovLdzFm5g0Xlu6mta8h2l9oNNgMwwhFmhO13JIMTHlo6EabcBIv+AB+thvL5zrbhUwF1vqfDDp/IGVs6MbF5KNdmAq7N389Ja+4uSCVQXlnDY2+Xo6ocaVBKioUXV+zgytPLrLJnGrAZQDbJp1LCyUbYngnHcyT7cwJE4JQvRZ9v0/yI8M+0AM6BWP3KObPZV9Q7+Y6u8Pds/k+ev7zJJxDlGC4AausaeOztcmrrG5uWeDzSoNTWN7rtNhNoLTYDyCb5VEo4aISdqIQ0CmVTnftQiLtId1uMvnMgVr++UVly3lNJ9xu2840m4e+N+P0+gV19J7O19BOZ7m5OsHzbPjSOwlNVVmzdZ+WeW4kpgGyS76WEkyow18Fb7nO89hgCB7ZHvq96GsZfklklkO5Y/QxQW9fA8m37mHNgLOuH/hAZdTGdfMpq6dibHeFfQFFAiRZ3P9KgVFYfaeMetT9MAWSTbCdOtZakCux4Z+DvVwAHtjszA08xlM+PPIdcvN/GxkgZan8Ya7z2FhBr515YPArZtC7azi1SMCN/D29x9yAlYIu7pwdTANkmn0sJJ1NgALs3RB/Tc2jk/vqNdKKAsrmgSzLWPu+UwKhc5zizi4oc4f/WL5zVySC4THVI/HZuD0fgKY+9Xc6tF4wp2DVsbXH3zGNO4GyT76WE42X+gnMf5fOcEf9JV8GAsY7Q9JaWLD0e+o2Ck67OXXPXmIsdpbV/myP0/cK/51BneysIY+cuVGxx98xjM4Bskix7NR9mAvEUWJ+RkfvqOxKWPOLE/Xfr77T3GQl7NkSc3bl6n0VFzsjfE/ovfstp7zk0MiNoBWbnToy3uPuKrfuorD5C/+4lTBzWy4R/mrAZQDbJgfDEVpEoOWz3BmdkP26WM9IfMc31BYjTvntDcmd3roTJekrATzLhn6Dvwypeb2r37NxBtAc7dzqSuLzF3c+fUMrksr4m/NOIzQCySQ6EJ7YKT4ENHOfU2on1CRR1iOQF+NvLY3wF8Ub/8aKMVrmmpUlXw+CPRfZPd30hD8/s4+etXyRWArF9d/vnxfjPO/FXbC39REp2bi9SqKq6ln7dO3H80F506pi7wtCSuHIfmwFkk2QVOHPVLOJROhGOOcup9bP6WUcAizjKYOA4x3lasdzZtyVVQuMt8uIpkKr1mV94Jtbmf+E9zX0CIfseVNQtrJ27vLKG2XPW8sLy7cz7oJIXlm9n9py1lFfWpO9e04glceUHNgMwWo4I9B7ufPaE8vhZzgh912rnuzew9UbufvyF5OKdPyjKqGwaTWUkvH0ylT+x9vmI8PdG/H6fwNrn4y9W7+t7KVAak+DlkczOnY+RQpbElR+YAmgNiTJhM2GKyEUGH+9E+ZTPd+P6Y5K+SidGm23A8QF4DmAIpwT8YbLjZ0W2ZTp/wovy8cf7e0rAywOIR0Df49Xz8ezcQeSjMDXndn5gJqDW0FZr4OaKMzQIESeTt2xqdLuX8bvmWccM5An/4VMdpeF3dq/5e/x7iBdlBG2z8ExRkTPCj7X1x2v3E9D3ltTzyUdh2t6d2+0FUwCtoa0WIs/mYuuhifln7zvSMdVsnOeEgMbu1sxXEHAPiaKMVj3jlJHwk0v5EzF9r5CBLS7qlo/C9PihvZA4ytiSuHKHjCoAEZkmIs+JyDYRURG5Kma7iMidIrJdRA6JyFwRGZ/JPqWV2LDNeMsktpa2UjQtIda847FpfnMh55mKvHtY86zjKxg4DgaOj57NqDozg6Aoo7JpkTIS8dYmyDYBIb5Lx97cpARSWdoxH4WpJXHlB5n2AXQHVgKPu69YbgZuAq4C1gF3AK+KyHGqeiDDfUsPQTZqT1h5tMQnEHtMa2oGZdJX4TfvlE11zEGrnnaE86b5MTtLpEKodw8DxzlKYMkfnfcR05znt+ZZZ5+ew5z2nSsjkVF9R0YyjOOtTZBKWYlMPJ84Ib4tKermCdPYkEoRaVNhmmoYqiVx5T4ZVQCq+iLwIoCIPOrfJs6Q5kZgtqr+zW27EtgFXAY8kMm+pY0gG7U/Rryl5Z2DYuDHzmxZzaBMlp32Btue8PeEcdUGxwcwYJyzFoB3rbJp0cdPvs4n7Ic671UbIpE3+7c2n+UMPj69+ROZeD4JKpC2pKhbtoVpS2P6Ezm3jeyTTR/ACKAUeMVrUNVDwDzgjGx1KiViTTEX/AK69I2OEfe2DxznLEsYllizT1AyUlhzRyZNSJ4w9pd03rnSEf4Dx8HJ1zU33fhZ86yj2EZMixRX878HzXLSnT+R6PmUTXO+t9QBr8r+F+9MmhEchmxlxFpMf3jybfnKbIaBlrrvO2PadwJDgw7YUbmXyV+8K+4Jr585nesvmZGWzoUi1s5bsQIO7YaSbs3rxvjNGGGIZ/bx4tG9UTMknwlkouy03zzivydvMZgTr4IiiTm3K+wGjIWTvxR9D7GzG4+xMzMfauvNrmoqmz+f3sfA0keiZzipzA4qVtCzcS/bh18eSZQLyAjOZcKEoU4c2iuvspQzQT5mPudCHkDsX5YEtAEwuH9v5j5wS+Z7FJZYO29sfXwPbySb6kg7yL/gmZZStXmnu+x0MrPJMWc50T1+5Vg+31FgH61xlKH/Hmoqg6+z+A8R30AmV02rWBFJXvMYOzOSvNbSpLPSidTQlTGbngCcPICgjOBcJlkY6oaPanhhxY4WC758K3ERRD4m60F2FUCF+14KbPG1D6T5rCA3ibXzeiNJz4bt0XNoc8dwGIL8C2uejXZ8hrV5x4unb6kSGDTBMfEELQYzcBwcdxFoffQI33P4eoLTP/LetTpSYsH/7rVnfNW0AAE3/+eR1csGjG3Z7EmEA9KdLUdfwphNTzQpgqCM4Fwl0cIsHYtg1fZ91DdGtqUi+PJx1BxEPibrQXZ9ABtxlMA5XoOIdAamAguy1alW4YU2+oU/RGrgpxKemCgG3rNTh7V5hzlXqngj5h5DokNgewyJmLv8IbIvfqv5SB6c/bxQUG+mNOWmiE/Aa/eUTSZCbVWbL1wDEeFfNtUxWflJ5dq+dX098kX4Q+IwVEXi3kay9Qzak28hH5P1IPN5AN1F5AQROcG91tHu96PVUZe/Am4VkUtEZALwKFANPJnJfmUMzyfQM8aF4Y1gU0nYSmep6EyUnfb+6f3r+/q/C+EKwHlmtMnXOe/jZkVMXP72ydclPk9r8MxTZVObZzSDE3bqLWLjkYridG3+flqSEZwtEsX0jx/Sk7oWCr5MLIaTLSdsPibrQeZNQJOBf/i+3+W+HsOJ/b8b6ALcD/QB3gHOzZscgFg8s4h/pNuaKKB0hTpmoux06fFOWYdmsf7A8CmORaWxsbngXPSQI8y9Egp+M1qsOc37XjoxvearZvfiPp9BE5pfZ8BYp+roprea/6aQvA+q9NDqqEJwng8A8mcmEC8MdfnWfazZsb9F6/ame9ScTXNSvi5fmdEZgKrOVVUJeF3lbldVvVNVB6tqZ1WdrqorM9mnjOKZM4JG2p5ZJCzpDHXMVNnpeMcd3O1Ezrz1i8jMwyujvGu149hNYfScNvNVotDNpozmGIX20RpH+McmnYWdPVWsoBsHo2z+Lc0Izia1dQ0s37qPyupa+nUvYeJQJwehNVnK6Rw1Z9uclK+Zz7kQBdR+yPcFXlKhYoUT0+8VffPoMQQ+Wg29hsG+rREHuOcb8ZRAxYpw0Tue+apsqmOKgejooT4jIxFYyRRZssglj7JpblnrpyMKoe/I6PP3HelcO0QU0B7pFT3Sb2FGcLZINrJuaZZyOkfNueCEzXayXkswBZBOEmR/piVUMZ20tvzBoAmOeeSjNdHtB7Y77SddAwvvdZSAlw/hlXnYuTK8MvSUKuoUlfOE97hZjgDevcFdHSxEOGhsmK7fnOPdiyf8RRylHTsjiA1BDREFVCudA2dfmY7/T0d4ZZjwxpYKvnSWuMgVJ2w2Mp9b8zubAihUWlv+oGJFtPD3MmY3zXfalz7iCH8/3nVSUYbe/o2NzcNO97jCP6x/xZuNqTYP6fQSzvxy2r/Wwe4NzvKTuVKELwnpsoeHHVm3VPCla9ScKFQ1l52wraW1v3M2w0CNbNLa8hCeoBw+NVLiQcT5Do4ZKDYaqjWVOj3/ihdR5YWD+rOswyAC/UZFt3mKac8GZ5bh2fVFnOzfEdMcJZCpaq9pJp328LYYWbe0xIU/4qe+Ic7SnOS2E7Y1pON3thlAodLa8hClvoJs4CgEvx29x5DoWj6pRM54+M1R6cqyVnWievysetq5gfIA5ZfuDOo2IJ328FwdWQeNfBXoWCwINLUBnDqiL/9Yuytvs4zjkY7f2WYAhUyiOP1kq5BBJIIo6DwHtqeWdxB0Pc9Mteghp33sTOgRkGMRNsvam+U0xfy7lUm95SzLApRfvAzqHI7hDzNqDxsvn4trEcQb+dY1KCLC+eNLmTZ6AKcd0w9EWPhhFfM+qOSF5duZPWct5ZU1bd7nTJCO2ZkpgEImkXBLZRWyoPP4F3GBiBLwzxr8BF1v0ISIiWfRQ7D6WTiQJMs6keLyFpgZMc0x7YyPUVqxkT6ZyKBuA5KFV4Iye85aXli+vZlgjFUMAFeeXkZJseCdslic82QivDGRYvK2PfnuZuob45h8VOlQXMSZxw3gnY27OZKHWcZhlXM6wmjNBFSoxAq3WDONV6I5WQ2eROfx1y2CxA5gv4mnpjKyToAXOvrRmuYRRxDxCXgF8cIUqRv7b86xsUpr94bovIimENRpTsSR9xzAveZIIGQIahuSKLwSYKErGD28qJ5HFmxscu34HYrnjSt1Z3rqnNLL8k4ziRyaQNS2eHgj31TMI7V1DSzdvIe1FU7+6XGlPZh0dJ+smIpSceqmI4zWFEChEq88BEQEahgfQZjzhIn68Y7zCsP5Q0fHfApeuimyrxeqGZRlnSjUc8S0aOEfpLSEyH14IaiqTlSTPwS170jHl1A+P30VSdNEovDKU0f0ZeGHVYHHxZZ08BTDc8ujy300NEJDY2Naq1wmCjd9dEE5kFjwe3gj38qQ5pHyyhoeWbAx6t4/2FXNy6squPqMEZT179Zm1UpTrSiajjBaUwCFSpiktTAO0HQmv4k4I39P+IMj/P95T8yOGukLOP3zLxmZTHHtWB5OaXkzFtXmSmX3hsiaxDkYDlrWvxvfPGc0L6/cyUfVtQzo3onzJgxiwfqqUII0DOlMsFq+bR8NjcH9qm9sJOx8o1GVQ0caWFtxIG5deU9J1NY18OiC8sBaRnUNyqMLNvKFU4fz5Lub26S8xPJt+2iMY9qK96wThdF6imt3TS3FPfoPDjqvKYBCJUzSWpBtf9UzzujX7wD2hOSO5RFF0JLkN6+aqh9v5N9jCIw+NyJ4PaEdpGiSKa5UlVYrIqYqXvwpb5z3WmrPIQ3EmhJ27T/M6h37OWVE37hRPamSzgSrnfsPR5WU9hOnOZD6RnhpVUXCfTzzyNJNe+L7EoCGRuVPCze1uNR1qnz4UTV1cbqT6FkH5WD4f/+9h+ro0GvQqKBjzQlsBBPPAVo+zzGHrHo6uXO4Nde7MGbUP+UmGHJCJC7fiyYKqmeULHIn3jHxMqO9qKdklU1zhETx4e98uDtt10lnGOihI5l3ynYoitTm2bH3MC+s2JFQuTSoM6MIoqXVSuNRW9fAqu37427vWEToZx30+yMSKOttBmAEE8+2rzhKoHw+IM5swCvH4JlDWrJco/96Xu0gP0v+6NTkT2ZeSubcjie0k2ZGX918zYA0VSRNt405kQO0QRs5tawfSzbviTJrgKBo3NLOwQh1DY3MWbmj1f3uUpL5sWijKt855zhKiouYPWdtqJlFvH384bTp+O2Wb9uX8M9ICR9ym+j3j8UUQD6Taj2foHZvhCs4yV1e+6AJTsTMmIujzSTjZ0E/zwE6L7LIuxdDDy1brtFfktlbK9hTBt6ykH6BG++8LXVKJ3Iel02LOHxbk9gWQDpLGO8/dISXV+3kg10H4pp4Ghph4cYqrjy9jH2H6qLsxjv2Huaxt8upb2gkkR7o6IYeKsrLqyrSYhsf1LMLHYocE04s8Wz5LWFdxQEUQgvIePjDadPx21VV1yZUvuOH9AxtbkqUHxCLmYDyhWSJUp4tM5E5JijWfsdyx6Sz5BGoWB45x5pnnTV9Y0ssiDg1ccZfEttB562ltXI8ob5zZbQALypyRv5hyy97iiQ2/DRRDoJ/H+86/rIPfUdGC/9US0LHIZ0lGxZuqGL2nHUs27KX6trExzUqPPnuZiYO6xVVfsFzHCcSHQKMG9wTwXGUtrTfsbHuxw3qTnFRsDhKl/BvVKisPsLO/YfT4ANxEsyCfrs/vvUh8z/YxfPLt/O/i7fw/PJtSRenSRTT37FYGDmge+ieJTpXLDYDyEWCRuqe8B44LrKgij9RavEfHEHplTAumxot7FSdV1lANIvH5oUwaGL0CHzQhGjnrneuWBt7+fxI5cxY52gqM5LSiXDS1dFF2TyB22ck4N5HvBF3ayqyxnMeQ0bKfKerZMP+Q0eahWomI97511VUJz4O2Lz7YNyInTD9Dpz1iHDeuFJeXl0R1d7o2uFTcQTHo6M7an93Y8v9IMUCHYqLOGVEX96JE07boPDSyuhlzTsUkXB2kCimvyjFjOtkeSBR5w59VqPtSJYV6y2oEltj/4Vv+soXBzhFlz7ijGb9o1yvBIJXDvnFbzW3xftnE0HOYa+kgkesSSSVrGIvgmjJI82zbWOLtaWbIMW28uloBejNxBobU/dzxJCuQmsvr9qZfKck5/dG5IvKdycVtnsO1sU1ESXrd6JZz8urK/jWOcdy8fFDmDZ6ABcfP4RTR/RNi/AHEHfUHi/aKAwjB3bn1gvGNNUbCkt9Iwlnrtd2NAAAFe9JREFUSOlcUCboXKgGxhfZDCAXCbJHBwl7iAhqf+x8U3VOom3oI6Y5o+DSidGj3PFuiWT/OTzhH2vOibWxA81GGrHO0WTJWbEj6DD7J4vaifVpJGqvWNHc99B3pKNsNs13Xt6KYN4+3tKfXmZxvJlNAtJVaO2j6tpQ+8U7f+yIvDUk63eyWc+6igNRs4dF5bvTErbasVg49Zj4SXBhKCkWJgxxYuwT/XaJqG9oZOmmPZw+sn+zbelcUMZ/rie7lFC5Y/MHQfuZAshFEsWdxwr7oIgZbYyYerxjB46LnHPxH6L3XxVQ18afiesdt2O5Iyj9tnTP5DR8CvQf7ThLY52j8e7Hv8oXRAvjPiMdvRK7vz+ZK9nqXrHJX0sfcdonXe34MWJLRMQqNk+RgrsmwIfRSnjgOPjwDdbthWe6fY6yno1cPPQw3T9wzjn/qC/zypGLqC7fHRgdkq4VsQZ078TWPYdC7etR36jUNTSy/9CRZtmnrUGhqd9BETKpznpSMWckQnDi+lOLcoo5h+83aWm/GhReXFHB4F5dAk1BQTH9Qc8RSBp95J1rTrcSdh6oDEyOMAWQqwTZo4OE/Vu/cIRS2VTwShpvessRyH4mX+e8r34mUld//zYnwarcd43hU5zj/deMJ3C3vx8xOfUf7QhVb3QdG3ETdD/9RkWv8uWZiry+nXR1dN/6jQo3qyhz1yTYOM/5/xwfwtcx5mLoUxY9o/CKxXl92L8t8j5iGosGXMraiue5YvfLDNxVwg/0C9SteIYriubxp8bz+OmHUzjYUEzJ8u2B9t90rYh13vhBLNuyN+52f4lkj0aFl1dV8NLKCtLnZoXGRqW2vqEpoig2QiZRIlrQ7CHoGbWUnfsOJ92nSODiiUOa+SJifxOvXw//c2PKJqUGDZ9IFuQved7193i/aWuij0wB5CpB9mhP2HsCc9FDEWHuN7l4SsDPmmed0Xasfd8vkAeMgdh8Ec+5HCRwPaFaNs0R/BDfORp0P1Xro53SY2dGK6bYuv1V6yMKJll27o7lkTLP5b4ZRNWGiK/Dv3+Qk9gLe/UrIZfq0Zdw1Yt9qKn/InUdhGs6zOEa5gDwcP35/KD+Cjw/TKLs0XRM+3t2KeGC8aWBGbCfPH4IJw7vzdJNe3hxRQUNvpleGGEqwKiB3Smvqgk1em5Q+K856ygqkqj9vWeQKBEtaNaz/9ARFpXvpl/3EopFGNK7C506FLFw427qGxpD+weONCgbkpSBFuCi4wdz2jH9OHF476S/SVn/bnz7vGP5+cv/Sl0JNDYmdZbHrw3U/N5amp1sTuBcJMjROnBcxPzgjconXxdp37nSFYozo8914T0RwV213jF/eOGVsZmtXQdEErouvCdi5/YEd2yYpLfv+FnRNvDYTNu4WcXzAV9k0ovfigj/A9sd2zs4272lGWMzeuNl53rLOcZyYHvw/ol+hwC2LX7OtWWLK+wj/KD+cqLDmLzTBWePtnRFLI/yyhreWLeLDjH/zaeN6MuJR/emU4diOhQXUdyC/3YFNlUdDF2LBxwlEE9Z1Dc0MrxfN0qKJamz0x/aun3vYbbsOcQ7G3fTu0sJ371gDB87qndaBZgInHR0HyD8b9KzcwnXfHxEtMM1BPWNsHN/Yt9NKgldzjkbefKdzUlDTv2YAshFgpKZYoU9RGLkPZu8anP7/ppnI6Wdvexdv73cz6b5kYSuoPj7RAI31fvxlEn5/Gg/ADg1f/yMn9W8BAQkL/kQK7bK5zvKJe7+PvxKy1va0vd+3O43+Lb8GWjkjg5/ijr0jg5PEGRWiRchE7b+exD+UWKsGX/hxt385KU1TedtuflEOe2Yfk0JYK2hEdjwUTWNCoN6dGb0wO6cN760aXF5j0Shrc8t386R+kY+efwQOsZqvVYwZVT/ZoI+9rfxZiT+38qbxXnRS588fjD/9rHBDEjixK+urUv4u6eS0AVOkt8Hu6pTWvjGTEC5SFCxMk8gx0aX+IuxefZ9L1fAb+IZOzNilgkql+A5c/2CK9acE0/gJlMCiYqv9RnZvMTCB68EXyOoL0ElH9yReVOIq9+Ec2C7MzMYf0nibF5PaXmzoBiz2bbuE7mmeg7DZSefKF7WZPa5o8OfuKaDYwrym4Eg2Mbd2kzgZKPEugZl3geVCc/RoUhoaNS4noAjDUpDo6Y0C0iEF9u/Ze8hOghsrKxh8+6DHDOge5MzM1lo68urKvj3SUfx+ZOP5k8LNzWZtrys4XhZxYn4+KjoyJzY36ZDkfBMozad2/9bDe7V2X1+SnFxkevwPsJH1fGjjt7fuo/V2/dR5zvX508+mn2H66iqruXA4Xo6FkvKjutYk1AiTAHkIi1JZkqlBEJQOeTxl0Sct/1GRztv/Qom1Ro7yfq9x1dHaOzMiJ+j51CnAJxfiY2blfge/PcLbn5CwD+Pt/JXohIR/tIUO1dGl8fuO4pefSfyhzmvcV3xC1E2f88cdE2HObzTOJaXG0/2PYZoG3eq9d+DSHWUGERxkXDOuEG8srqCoLXVS4qFgynMSlKhXgFV3t+6j5Xb9vHC8h1cdUZZ0tDW97fspWL/YSr2HY7yAyjOz3TKiL5sqjrItr3JHb+4x/hDUIN+G8/O7zVFLaTjzqpjHd4diiShf6Au5lyPvl0eUWLiPp8WEqZgnSmA9kIqJY5bUsM/XQu/JDpnxYqI8PdMXfGukWxWgUbX7B8703Gaf7TGXfnrY4nv2a+0/PfltncHPnbGBdywYBSv6WRA6FrciAh0GP8pblh9HPPlJIC4kT2JRu+ePXfC0F4JC4y1NB49tl+De3Xm9TW7aKC5BhARunYsTts6AvFoUGhoaOSRBRs5um/XpPtujyPcVWHBht0p2bdViTLPpWJ/r3OFt0e0wzv1Z+Yd0Rrh7/UjWTKhKYD2QiqzhpbMMNK58Eu8cyYZdQeavoLuYcjHIlFAfoV1yvWRWYs3ywlTIiIOJw+oZ9xFx/DxLTWUVxdT1r2Bi4+qpVsH5dMjnPaV779LzcQrA6NIEo3ePXvupqqaKDPDks17WOcuXTimtAfjh/TkRfHGjOE5qk8XTi7r29Sv8soagoxAHd21fz+qrk3bOgLJqGtQNnzU+oXbU7EAdYwxz6VjZqVoWpaN9ubV8U7lufRiCZNMaArACEdrauyEPWeSUXdKZEJhBdCtg/LZEc1Hot06KNPX/D9O6NybJWXfDDw2zOjdb2ZQ1Si79ge7qpmzqoILxg9mzqodoQVWSbFwclnfZuaOIFuzAIN7d6Zvt478PfOyP2vUNShl/bqyqHx3xP5eRNwFWsKeMx0uasXJT4gn5BGJWuPZw0v2i7emAVgUUEZ48Om52e5Cu6JFzzM2FDVZe4ZYcs5TUd/9USX1DY0pmRmCnJp1DcqLK3dw0tF9GD2we9LbOvDenGa+iETLMTYq/GPtLn75WmAlgXbFPa99wLPvbWPeB5WsdJ2zyTjw3py424TUZiHxKMKLUGpeJ+iqM0ZwVWzdHxcv2W/z7oMUde4RmFZuM4AM8OCzb3L9JTOy3Y12Q3t5nkERP+kodFbfqCzc6NTMSWYJqn5/Djd942tR5qhEyzHWNypvra9MW0G2XMfTx2Ejb6rfn0OPE84PPlea+lRUJJw5ZiBnjhnYLDkNdRT4pOF9qD5cz8pt+6KUzpEGpVGVjv2PDpzymgIwDB/VdcLzWztRfqCYsh4NXDyslu4dW/+vHDarszUkO1+xQO8uHZvCS70aM+t3JS4BXSjCP1e54rThTQrbnzkcO6AoLko44wicG+aEAhCRrwHfAQYDq4AbVXV+4qMMo3XsPFjEz1Z148P9xRzTs4HzBtfyrcU9UYWDDUV0LW7kh+9359Ep+zi5f13o83pKZJWvEFyqWZ2ZoEGdkE9IbwXQWIpIj+mjvZLMqQuOzb9IhCtOG87oQT2abQ8aUASF8EYumqNrAovIZ4F7ga8Bb7nvL4nIOFXdnNXOGe2Wx9d35o73vH8s4b09ytObOuMfKB1scP5nrnqrF+9eXEW3DskF5aLKjlz1Vi9XiZzdVAhuTGmPNomgScbumiPsP5zeCqB+ikUYP6Qnq7bvS7isZCET77F0KIIzRvYHJGlNqJQHFHHWA8gFJ/C3gEdV9SFVXaOq3wB2AF/Ncr+MdsrOg0Wu8BciAt//ORpVeH5Lp6Tnra4TrnqrFzX1RU3Kw1vsZOW2xAk5xUKTE69jsVAkmfNTv7xyZ8ZmI8VFMLxfVzMbtQABzhwzMFRNqBaEqQbunNUZgIiUAJOAn8dsegU4o+17ZBQCP1uVWsncgw1FlFcnL9D2/NZOLY77Pm98KZ07Fjc5+Cr2H2bBhsSLlxQLKY+yFWcBmUzNRo40KPsO1TO0T5eU1ygodBqV0MuAJgsh9v42SoqFIhFqP9oUuIyeZNMuKSJDgG3AdFWd52u/A/iCqh7na9vaq3uXoScce3QWepoaS9aWM2lMWba70W5I9/NcuacD1fXhh9fFAsO7NzCwc2KTyeaaYrYfbNmkuneXEvp2iyTt7D9cR1VNbVyFIkCvLh3Zf7geIGGst58jFevpO3wMB480hD4mFYpE6NethEacUWp750jFekpKRyXdz8nVk8BkOz+xfwfxaFRl8+6Dgb9hkQh9u5VQ36B0LBaqNv+LQzUHtqnqsGb9yhEFMM3v9BWR/wd8XlXH+NqWAWX8//bOP7ausozjn2+72S5r6SLDOX7oFucA3WS/IKJMYYIohkTQBJBIZsyyoZmGaCAwQvaHigQzJEw2jSSbQICISkCJP5Y527jKjy1jrRlM1G7KZGzOrSvrurV9/ON973Z2em97b1t2e+95PsnJPfe873nPe769Pc95fz3P4K2WPYTuo3IzlbFRj2phVPUcP/n9F9RMaJxS9AlmfT1v7GjF+gd1iFPbOHnquKYpM/IPuCV9WA9I6u89tPdvfcmoTaqprTv3Qx+nUL8U0Nd1YHfvwTd31za8+6zaiZOmanz9GUPdSl/XAevt3PeXuqkzL0EqPdbgUOS0AurOufBjhQYfSyrxeE+Xxr+rAcNGobxRpa/rALUNQ7+x9/ccOdB/tOu/45rOmlG4c8+s9+DenX0FonelqalvbIrTO4VUE/v57fj+3W39Rw8n+xxnAPvMbG66jHIPAu8H+oD3po6/BzjFHWC+yjuO4zjDp6zW1MyOAVuAq1JJVwGbT3+NHMdxskO5WwAAq4BHJb0I/BlYBpwNrC1rrRzHcaqcshsAM3tK0pnA3YS+3nbgGjPbVd6aOY7jVDdjYkDFzB42s2lmVmdm85MzgsYCkj4h6VlJb0gySYtT6ZK0UtIeSd2SNkn6cCpPnaSHJO2X9HYsb8CofLUj6U5JL0nqlLRP0nOSZqXyuJ5FIunrkrZHPTsltUr6XCLdtRwmku6K/++rE8eqSs8xYQAqgAZCy+SbQL7JzbcD3wKWAxcDbwF/kJRcw/1D4AvATcBC4Azg13onZmKMbS4HHias81gE9AIbJCWnUriexfNv4A5gHrAA2Ag8IynnP9u1HAaSPgosAbankqpLTzPzrYQN6AIWJ76LMEVxReLYBOAwsDR+bwKOEdY25PKcR3CZcnW576nMejYQZoJd63qOmqYHgKWu5bD1awL+TnhB2QSsjserTk9vAYyc6YRprCcimZtZN9DMydXM84HxqTz/AnbgK54bCS3R/8XvrucwkVQr6UaCUd2MazlcfgI8bWYbU8erTs+yDwJXAbk1DHtTx/cC5yTy9BHWPaTzpNdAZI0HgW1Aa/zuepaIpNkE/eoJLdTrzKxNUu6B41oWiaQlhIVTX86TXHW/TTcAo0d6SXUxgVpLD+ZaRUhaBVwGXGZm6VW2rmfxvAbMASYR+p7XS7o8ke5aFoGk84HvAQstrFEqRNXo6V1AIye3bHuw1cxvArXA5EHyZApJDxAGyRaZ2T8SSa5niZjZMTN73cxeNrM7CS2q23AtS+VSgg7tknol9QKfBL4W93Pe+apGTzcAI+efhD/6idXMkuoJo/+51cxbgOOpPOcCF5LBFc+SHgS+RHj4v5pKdj1HTg1Qh2tZKs8Aswmtqdz2MvBk3N9JtelZ7lHoStgIg2q5H8QR4J64/76YfgfQCVwPzCL8YPYAjYky1hAc310JzAX+SHhTqy33/Z1mLX8UtVpEeJPKbQ2JPK5n8Xp+n/AAmkZ4eN1LmHHyWddyVPTdRJwFVI16lr0ClbAR5q5bnm1dTBewkjBF7CjwJ2BWqox64CFCM/II8BxwXrnvrQxa5tPRgJWJPK5n8XquA3YBPYQ56RtITDd0LUesb9oAVJWeZXUH7TiO45QPHwNwHMfJKG4AHMdxMoobAMdxnIziBsBxHCejuAFwHMfJKG4AHMdxMoobAKdqkbRM0i1xf5OkBe/QdZ6XNCnP8ZWSvh33F0s6O5HWISntLiBf2XMl/XSY9XpS0geHc66TDdwAOFWLma01s5+dhutcY2YHh8i2mBDrulTuIiwqGg5rCAFMHCcvbgCcikHSREm/kfSKpHZJN8TjHZLuk/Ri3GbE4yfewBNl1EhaL+k70X/+/TFE5XZJS/Nc83ZJ34j7D0jaGPc/JemxxPUnx/0Vkl6TtAE4Px77IiFa1+OStkmaEItfLmmrpDZJF+S5diPwETN7Jd/9RA2mFdIFaAGulORef528uAFwKonPAHvM7CIzmwX8NpHWaWaXAKsJIfnyMQ54HNhpZncDXwUOmdnFhPB+SyRNT53TTPC1A+Eh3iBpPMGNdUsyo6T5wI0E/y/XxzIxs6cJTsVuNrM5FoKIAOw3s3mEN/VTDFXieu0F1ThJXl3MrB94HbioiDKcDOIGwKkk2ghvtPdJWmhmhxJpTyQ+Ly1w/o+BdjP7bvz+aeAWSduAF4AzgXSf+RZgfnwb7yEEXllAMAotqbwLgV+Z2REz6wSeHeJ+fpm4xrQ86VOBfUOUAYPr8hbD63pyMoAbAKdiMLOdhJB7bcC9ku5JJhfYT7IZuCK68IXg2Gt5fCufY2bTzez3yRPM7DjQAXwlnt8CXAF8gBDmb0A1S7ilnvjZR/7gTN0Ex2I5ejn1f7Y+1nEwXepjOY4zADcATsUQZ9EcMbPHgB8A8xLJNyQ+W9PnRh4Bngd+HvvFfwfcGrt0kDRT0sQ85zUTumiaCQZgGbDNBnpSbAaukzQhthiuTaQdJsQ/LoUdhPCEOTqI9yxpHiFG7VC6zAT+WuJ1nYzgg0NOJTEbuF9SPyHoxq2JtDpJLxBeam4qVICZrZLUBDwK3EzoetkqSYTuls/nOa0FWAG0mtnbko4ysPsHM9sq6SmC7/ddqTzrgLWSuincRZUu71VJTZIazeww8AtOdlm9RAhQAgV0kTQF6Daz/xRzPSd7uDtop+KR1AEsMLN0IO6KR9JtwGEzK3ktQDy308weGf2aOdWAdwE5zthmDSfHCkrlILB+FOviVBneAnAcx8ko3gJwHMfJKG4AHMdxMoobAMdxnIziBsBxHCejuAFwHMfJKG4AHMdxMsr/AcrEO7QSlMD3AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn import svm #1\n",
    "from matplotlib.colors import ListedColormap #2\n",
    "\n",
    "\n",
    "clf = svm.SVC(kernel='rbf') #3 creates instance of svm \n",
    "clf.fit(Xdata, Ydata) #4 trains the SVM on the data \n",
    "\n",
    "\n",
    "# figure production \n",
    "fig=plt.figure(figsize=(6,5)) #5\n",
    "markershape = ['o','x'] #6\n",
    "markercolor = ['C0','C1'] #7\n",
    "markerlabel = ['Pyramidal','Basket'] #8\n",
    "\n",
    "gridRes=0.1 #9 \n",
    "cmap = ListedColormap(markercolor) #10\n",
    "maxvalX = 450. #11\n",
    "minvalX = 50   #12\n",
    "maxvalY = 60   #13\n",
    "minvalY = 0    #14\n",
    "xgrid,ygrid=np.meshgrid(np.arange(minvalX,maxvalX,gridRes),np.arange(minvalY,maxvalY,gridRes)) #15\n",
    "Z = clf.predict(np.array([xgrid.flatten(),ygrid.flatten()]).T) #16 predicts class for the plot area \n",
    "Z = Z.reshape(xgrid.shape) #17\n",
    "plt.contourf(xgrid,ygrid,Z,alpha=0.4,cmap=cmap) #18\n",
    "for yind,yval in enumerate(set(Ydata)): #19\n",
    "    plt.scatter(x=Xdata[Ydata==yval,0], y=Xdata[Ydata==yval,1], c=markercolor[yind],marker=markershape[yind], label=markerlabel[yind],s=50) #20\n",
    "plt.xlabel('spike width (us)');plt.ylabel('firing rate (spk/s)')#21\n",
    "plt.title('SVM using spike width and rate, with '+clf.kernel+' kernel') #22\n",
    "plt.legend() #23\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4. Check the performance of the classifier \n",
    "\n",
    "We want to evaluate the classification error of the SVM. In other words, we want to assess how well the SVM classifies new data-points. We will do this through *cross-validation*. The idea is to train the classifier  on only 90 % of the data and test whether the remaining 10 % of the data are correctly classfied. \n",
    "\n",
    "Write a routine with fits 10 times the above SVM classifier to 90 % of the data, the test dataset. Use the remaining 10 %  - test dataset  - to evaluate the performance of the classifier. That means, calculate the percentage of erronous classification. Note, the split between training and test dataset should be different in each run.\n",
    "\n",
    "Use the function `clf.predict(xTest)` to calculate the predictions with the trained classifier from the test dataset.\n",
    "\n",
    "**Hint:** Split the dataset in train- and test subsets using the scikit-learn function `train_test_split()`. Check online which type of inputs this function takes and how to use it. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iteration number : 0\n",
      "iteration number : 1\n",
      "iteration number : 2\n",
      "iteration number : 3\n",
      "iteration number : 4\n",
      "iteration number : 5\n",
      "iteration number : 6\n",
      "iteration number : 7\n",
      "iteration number : 8\n",
      "iteration number : 9\n",
      "error rate :  0.03666666666666666\n"
     ]
    }
   ],
   "source": [
    "# your code goes here \n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "nIterations = 10\n",
    "\n",
    "allErrorRates = []\n",
    "clf = svm.SVC(kernel=j)\n",
    "for n in range(nIterations):\n",
    "    print('iteration number :',n)\n",
    "    xTrain, xTest, yTrain, yTest = train_test_split(Xdata,Ydata, test_size=0.1)\n",
    "    clf.fit(xTrain, yTrain)\n",
    "    #print('...fitting done')\n",
    "    yPredict = clf.predict(xTest)\n",
    "    correctClassifications = sum(yTest == yPredict)\n",
    "    errorRate = 1. - correctClassifications/len(yTest)\n",
    "    allErrorRates.append(errorRate)\n",
    "print('error rate : ', np.mean(allErrorRates))\n",
    "#print(sum(np.equal(yTest,yPredict)))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 5. Different kernel functions\n",
    "Test the performance for different kernel functions. Check the documentation for the scikit Support Vector Classification to see which kernel functions are available. Chose three and compare their performance by using the cross-validation routine programmed in step 4."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iteration number : 0\n",
      "iteration number : 1\n",
      "iteration number : 2\n",
      "iteration number : 3\n",
      "iteration number : 4\n",
      "iteration number : 5\n",
      "iteration number : 6\n",
      "iteration number : 7\n",
      "iteration number : 8\n",
      "iteration number : 9\n",
      "Kernel function : linear\n",
      "error rate :  0.009999999999999998\n",
      "iteration number : 0\n",
      "iteration number : 1\n",
      "iteration number : 2\n",
      "iteration number : 3\n",
      "iteration number : 4\n",
      "iteration number : 5\n",
      "iteration number : 6\n",
      "iteration number : 7\n",
      "iteration number : 8\n",
      "iteration number : 9\n",
      "Kernel function : poly\n",
      "error rate :  0.033333333333333326\n",
      "iteration number : 0\n",
      "iteration number : 1\n",
      "iteration number : 2\n",
      "iteration number : 3\n",
      "iteration number : 4\n",
      "iteration number : 5\n",
      "iteration number : 6\n",
      "iteration number : 7\n",
      "iteration number : 8\n",
      "iteration number : 9\n",
      "Kernel function : rbf\n",
      "error rate :  0.02666666666666666\n"
     ]
    }
   ],
   "source": [
    "# your code goes here \n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "nIterations = 10\n",
    "\n",
    "kernels = ['linear','poly','rbf']\n",
    "\n",
    "for j in kernels:\n",
    "    allErrorRates = []\n",
    "    clf = svm.SVC(kernel=j)\n",
    "    for n in range(nIterations):\n",
    "        print('iteration number :',n)\n",
    "        xTrain, xTest, yTrain, yTest = train_test_split(Xdata,Ydata, test_size=0.1)\n",
    "\n",
    "        clf.fit(xTrain, yTrain) \n",
    "        #print('...fitting done')\n",
    "\n",
    "        yPredict = clf.predict(xTest)\n",
    "        correctClassifications = sum(yTest == yPredict)\n",
    "        errorRate = 1. - correctClassifications/len(yTest)\n",
    "        allErrorRates.append(errorRate)\n",
    "    print('Kernel function :', j)\n",
    "    print('error rate : ', np.mean(allErrorRates))\n",
    "    #print(sum(np.equal(yTest,yPredict)))\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The end"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}