Canlapan_Inah_HW4.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4c38567e",
   "metadata": {},
   "source": [
    "# HW 4"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7d5b7969",
   "metadata": {},
   "source": [
    "## 2. Comparing Classifiers\n",
    "\n",
    "### 2.1 Divorce Classification "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "040482de",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.neighbors import KNeighborsClassifier \n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.naive_bayes import GaussianNB\n",
    "from sklearn.decomposition import PCA\n",
    "from mlxtend.plotting import plot_decision_regions\n",
    "from sklearn.model_selection import GridSearchCV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ef093ea5",
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('data/marriage.csv', header = None)\n",
    "data.rename(columns = {54: 'y'}, inplace = True)\n",
    "\n",
    "# Split test/train data\n",
    "X = data.iloc[:,0:54]\n",
    "y = data.iloc[:,54:55]\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=4264)\n",
    "\n",
    "X_train = np.array(X_train)\n",
    "X_test = np.array(X_test)\n",
    "y_train = np.array(y_train).reshape(len(y_train)).astype(int)\n",
    "y_test = np.array(y_test).reshape(len(y_test)).astype(int)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "15a69278",
   "metadata": {},
   "source": [
    "### Part A\n",
    "Report testing accuracy for each of the three classifiers. Comment on their performance: which performs the best and make a guess why they perform the best in this setting.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "e0a95a95",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 1.0\n",
      "Test Accuracy: 0.9411764705882353\n"
     ]
    }
   ],
   "source": [
    "# Logistic Regression - from heart.py\n",
    "lr_model = LogisticRegression(max_iter=200, solver='liblinear')\n",
    "lr_model.fit(X_train, y_train)\n",
    "lr_score_train = lr_model.score(X_train, y_train)\n",
    "lr_score_test = lr_model.score(X_test, y_test)\n",
    "\n",
    "print('Training Accuracy: ' + str(lr_score_train))\n",
    "print('Test Accuracy: ' + str(lr_score_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "1e60f1f9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 0.9779411764705882\n",
      "Test Accuracy: 0.9705882352941176\n"
     ]
    }
   ],
   "source": [
    "# KNN - from plot_classifier_comparison.py\n",
    "knn = KNeighborsClassifier(n_neighbors=5)\n",
    "knn.fit(X_train, y_train)\n",
    "knn_score_train = knn.score(X_train, y_train)\n",
    "knn_score_test = knn.score(X_test, y_test)\n",
    "\n",
    "print('Training Accuracy: ' + str(knn_score_train))\n",
    "print('Test Accuracy: ' + str(knn_score_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "5050ab49",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 0.9779411764705882\n",
      "Test Accuracy: 0.9705882352941176\n"
     ]
    }
   ],
   "source": [
    "# Naive Bayes - from plot_classifier_comparison.py\n",
    "naive_bayes = GaussianNB(var_smoothing=0.001)\n",
    "naive_bayes.fit(X_train, y_train)\n",
    "nb_score_train = naive_bayes.score(X_train, y_train)\n",
    "nb_score_test = naive_bayes.score(X_test, y_test)\n",
    "print('Training Accuracy: ' + str(nb_score_train))\n",
    "print('Test Accuracy: ' + str(nb_score_test))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e35033da",
   "metadata": {},
   "source": [
    "### Part B\n",
    "Now perform PCA to project the data into two-dimensional space. Build the classifiers (Naive Bayes, Logistic Regression, and KNN) using the two-dimensional PCA results. Plot the data points and decision boundary of each classifier in the two-dimensional space. Comment on the difference between the decision boundary for the three classifiers. Please clearly represent the data points with different labels using different colors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "4274eeaa",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Use PCA package to perform PCA\n",
    "pca = PCA(n_components=2)\n",
    "data_redux = pca.fit_transform(X)\n",
    "\n",
    "# Split into test/train sets \n",
    "X_train, X_test, y_train, y_test = train_test_split(data_redux, y, test_size=0.2, random_state=4264)\n",
    "\n",
    "X_train = np.array(X_train)\n",
    "X_test = np.array(X_test)\n",
    "y_train = np.array(y_train).reshape(len(y_train)).astype(int)\n",
    "y_test = np.array(y_test).reshape(len(y_test)).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "2fce14f7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 0.9779411764705882\n",
      "Test Accuracy: 0.9705882352941176\n"
     ]
    }
   ],
   "source": [
    "# Logistic Regression - from heart.py\n",
    "lr_model = LogisticRegression(max_iter=200, solver='liblinear')\n",
    "lr_model.fit(X_train, y_train)\n",
    "lr_score_train = lr_model.score(X_train, y_train)\n",
    "lr_score_test = lr_model.score(X_test, y_test)\n",
    "\n",
    "print('Training Accuracy: ' + str(lr_score_train))\n",
    "print('Test Accuracy: ' + str(lr_score_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "54c535e3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 0.9852941176470589\n",
      "Test Accuracy: 0.9705882352941176\n"
     ]
    }
   ],
   "source": [
    "# KNN - from plot_classifier_comparison.py\n",
    "knn = KNeighborsClassifier(n_neighbors=5)\n",
    "knn.fit(X_train, y_train)\n",
    "knn_score_train = knn.score(X_train, y_train)\n",
    "knn_score_test = knn.score(X_test, y_test)\n",
    "\n",
    "print('Training Accuracy: ' + str(knn_score_train))\n",
    "print('Test Accuracy: ' + str(knn_score_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "38ef7b60",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training Accuracy: 0.9779411764705882\n",
      "Test Accuracy: 0.9705882352941176\n"
     ]
    }
   ],
   "source": [
    "# Naive Bayes - from plot_classifier_comparison.py\n",
    "naive_bayes = GaussianNB()\n",
    "naive_bayes.fit(X_train, y_train)\n",
    "nb_score_train = naive_bayes.score(X_train, y_train)\n",
    "nb_score_test = naive_bayes.score(X_test, y_test)\n",
    "\n",
    "print('Training Accuracy: ' + str(nb_score_train))\n",
    "print('Test Accuracy: ' + str(nb_score_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "52bab359",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 1.0, 'Naive Bayes')"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlIAAADSCAYAAACB8C6CAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA0GklEQVR4nO3deXxcdb3/8dcnM9m3Jk2TbqGhC4VKEcuqKF4EBbSAgnIpq/7kAlVQUC+CIHoVxA0vXNlu4WJBhIIi0FtARVmqApfW0rLaFkqXtE3TtEmbNM3+/f1xTmCSZp05M9+zfJ6PxzxIZzJnPhPmPef7/Z7v+R4xxqCUUkoppUYvy3YBSimllFJBpQ0ppZRSSqkkaUNKKaWUUipJ2pBSSimllEqSNqSUUkoppZKkDSmllFJKqSRpQ2oUROROEfluEs/bT0RaRCSWjrr8SkSeEpELbNehok1EnhORC23XodRIiMh3RORu23WokQttQ0pE1ovICV5u0xhziTHmh6N9bWPMRmNMkTGmezSvJyJfFJFutxG2W0RWicjcZGq3wRhzsjHmXtt1KHv6Z0FEzhKRRhH5uIjUiIgRkSf6Ped+Efm++/O/uL9zW7/f+ZuIfDET78F9vYUi0uFmsfcWqY6RGhn3M79NRAoT7rtQRJ4byfONMT8yxnje8Hc7FG3uZ3eXiCwVkdlev04UhbYhFSIvGmOKgDHA7cAiERnj9YvoTkGlmzs6eRvwGWPM8wkPHS0ixwzx1D3A+SJSk876RuCnboeoKJmOkYqUOPB120UM4FJ3fzIWeA74td1ywiFyDSkRyRWRm0Vki3u7WURyEx6/UkS2uo9d6PaGp7uPLRSR692fK0RkiYg0ichOEfmriGSJyK+B/YD/dVv+Vyb0vOPuc8tF5FfuazSKyGPD1W2M6cH50BcCMxLey89FZKPbA7pTRPJH8V7uEJEnRWQPcJyITBSRR0Rku4i8KyJfS9jWkSKy3B0Z2yYiv3Dvz3NHEHa4f4tlIlLlPvbeIRX3b3OtiGwQkXoRuU9ESt3Hev8+F7jvpUFErkn6f7LyHRG5CLgJONEY80K/h38KXD/E05uAhcD3PKhjgoi8KiLfSnVbSg3hZ8C3Buv0isgtIrLJ/T79h4h8LOGx74vI/e7PfxCRS/s9d5WInO7+fKCIPO3ug1aLyJkjKc4Y0wUsAmYlbPdIEXnR/R7fKiK3ikiO+9htInJTvzr+V0Qud38e9b4jTCLXkAKuAY4GDgU+CBwJXAsgIicB3wBOAKYDHx9iO98EaoFxQBXwHcAYY84DNgKnuL3Wnw7w3F8DBcAHgErgP4cr2h0x+hLQCWxw7/4JcID7XqYDk4DrRvFezgZuAIqBF4D/BVa52zkeuFxETnR/9xbgFmNMCTANeNi9/wKgFKjG6eVcAuwd4LW+6N6OA6YCRcCt/X7no8BM97WvE5GDhviTqOCYD/wQON4Ys3yAx28DDpChD8XfAJwhIjOTLUKcEa3ngVuNMT9377vd3XEMdHu13ya+4u6w/iEiZyRbh4qE5TgjPoM12JfhfG+XAw8AvxWRvAF+7wFgXu8/RGQWMAV4QpxDh0+7v1Pp/t7tIvKB4YpzG0jnAC8l3N0NXAFUAB/G+R7+ivvYvcA8Eclyn1/hPv6ge18y+47QiGJD6hzgB8aYemPMduA/gPPcx84EfmWMecMY0+o+NphOYAIwxRjTaYz5qxnBhQtFZAJwMnCJMabRfe7zQzzlaBFpAtqAnwPnGmPqRUSAfwOuMMbsNMY0Az8CzhrFe3ncGPN3d7RrNjDOGPMDY0yHMWYdcFfC9jqB6SJSYYxpMca8lHD/WGC6MabbGPMPY8zuAV7rHOAXxph1xpgW4GrgLHFH6Vz/YYzZa4xZhRPKDw7xd1HB8UmcL+zXBnm8DaehNOiolDGmDrgT+EGSNczC2bF9zxizIGG7XzHGjBnkdkjC8/8LZyS4EvgusFCGPhyp1HXAZSIyrv8Dxpj7jTE7jDFdxpibgFycTmR/jwKHisgU99/nAL83xrQDc4H1xphfudtZATwCfH6Imv7L3Z+0AJeSsF9wv7tfcre1Hvhv3A64MeZlYBdOIwmc/cJzxphtwBEkt+8IjSg2pCby/ogO7s8TEx7blPBY4s/9/Qx4G/iTiKwTkatG+PrVwE5jTOMIf/8lY8wYoAxYDPQOAY/DGdX6R28PGviDez+M7L0k3jcFmJjYI8cZZatyH/8yzujXP93Dd72T3n8N/BFn7tYWEfmpiGQP8FoD/d3jCdsHqEv4uRVn1EoF3yU4n5273Q7AQO4CqkTklCG28xPgRBFJpoF9DrAZ+F0Sz8UYsyJhx/ck8Bvg9GS2paLBGPM6sATYZ98gIt8UkbfEmfTdhDOqXzHANpqBJ3i/UXIWzmcPnO/so/p9Z58DjB+irK+5+5M8nIbY70TkELemA8SZrlInIrtxOuaJNd0LnOv+fC7vz69Kdt8RGlFsSG3B+R/faz/3PoCtwOSEx6oH24gxptkY801jzFTgFOAbItLbWh9qZGoTUD7YsfMhXq8FZ5j1PBH5ENCAcwjtAwk96FJ3IuFI30tinZuAd/v1yIuNMZ92X3+tMWYeTo/8JzgBLHRH1P7DGDML+AhOOM8f4LUG+rt3AdtG83dQgVSP05P9GM4JE/swxnTi9I5/CAzY2DLG7ABudn9ntL6Pk5kHJOHECnHmFbYMcntjiO2ZwepUKsH3cI4cTOq9w50P9W2cowZlbsNmF4N/nh7EOaz2YSAfeNa9fxPwfL/v7CJjzPzhijLG9Bhj/oozGPAp9+47gH8CM9zDcN/pV9P9wGluR+Yg4LGEOka97xiuxiAJe0MqW5zJ0L23OM6H8loRGece570O5wMCzrHbL4nIQSJS4D42IBGZKyLT3R72bpzjy71n8WzDmQe0D2PMVuApnGPZZSKSLSLHjuTNuDuSu4Hr3MNxdwH/KSKVbk2TEo5Lj/i9uF4GdovIt0UkX0RiInKwiBzhbvtcERnnvm6T+5xuETlORGa7O6fdOMO4A53N9CBwhYjsLyJFOL2dh9xJjyrkjDFbgE8AJ4nIYHMCf41ziOOkITb1C5wG+3vz5+T9kxVqhnheJ/AFnJM1ft0718M4S5oUDXJ7b66JiHxeRIrEOWniUzg98sXDvnEVacaYt4GHgK8l3F2M04ncDsRF5DqgZIjNPInTCf0Bzndmj3v/Epy5hee5+5FsETlCRji31G2YzQJ6OwzFON/hLSJyIM7cxsT3Uoszt+vXwCPGmN65sEntO0ZSY1CEvSH1JM6oTe/t+zjzMJYDr+LM2Vjh3ocx5imcuRDP4rTUX3S30z7AtmcAf8Y51vwicLsx5jn3sRtxGmtNMvDZQefhfLH/E6e3fvko3tPNwKfd4dhvu3W+5A7F/hn3OPso3wvGOZX7FJwJkO/i9N7vxhlyBmfn9oaItOBMHjzLGNOGM4z8O5wAvoUzmfd+9nUPTgCXuttvAy4bxftWAWeM2YTTmPq8iNw4wOPdOD348iG2sRvnLL/E36nGOVS8eZjX78A5HFcJ3NPbmBqhr7vbb8I5rP9vCXlXaig/wGnA9/ojTmd6Dc7nto0hppG486F+j3Pi0AMJ9zfjjCadhTPiX4cz4pM7wGZ63do74orzfXytu68AZ2L82UAzTif9oQGefy/OfNr3lk1IYd8RGmKGnx8dWW7L/nUgN+gjJ2F6L0olEpFrge3GmP+2XYtSYeYePbkfqEkYGYs8bUj1IyKfw5ncV4jT+u4xxnzWalFJCtN7UUopZY97EtEiYJUxJtmzZ0Mp7If2knExzrHrd3CO4w47cc/HwvRelFJKWeAe0WjCWfLnZqvF+JCOSCmllFJKJUlHpJRSSimlkqQNKaWUUkqpJMWH/xXv3bV0XWiOJzbWbyXv5QV87bOH2y5FJesjl1lfWDFMmei1YfVrfHDHU3zh2GEv/aX8RjORtFcWfpfbLzzadhnKS1nZcPQlg2ZCR6RSVFY5gfU7W22XoZTvNG54kzlTK22XoVRGSSw2/C+pUNGGlAeMsd55U8p3iqqmsHrzSC8pqZRSwaQNKQ/sJcd2CUr5TnF5Fdt37x3+F5VSKsC0IeWB9p44uoyEUn3lFxVTtytUV4JQSql9WJlsPhDBUJrdQ14MnOsA+4sxhrZu2NWZhel3ke6swjKaW9soKcy3VJ0KoyBnAqCkrILVzZ0WKlNhFYRMVE/Znx6ELLRzHRW+aUiVZvcwpjCPHomDDwOCMeSZLtjTRlNn38mEpVMO5tV3Xuejh0y1VJwKoyBnwnnY0NWjOxPlnSBkomTsBPbESiju3mG7GpUhvjm0lxfDv+EAEKFH4uQNcEJGaeUk1m1rznxNKtSCnAnnYaGNvMzWpEItCJlo6+qiO66f+yjxTUNKRPwbjl4iAw4nT6iZwRu1enaS8laQM9HL+L1+FSiByATAAIe6VXj5piHlF8v/9gxfPuWjfOnTH+ahu385oueICF3daS5MKUuSyUSvdsmno7MrTZUplXmp5EGFkzakEnR3d3PbDd/h+tt/w4LHn+e5px5jwzurR/TcdonWUG5DUwtnXHUnO3btsV2KSqNUMgGQlZNHZ0R6GZqJ8Es1D1ETlUz4ZrL5aHz9/M+xa/fufe4vLSnhlvseTXq7q197hQn71TChegoAHz/5NF589o9MmTZz2Oe29kRrNdv7nniBxrpN3Lvk73zjnE/ZLify/JiJXsY4E8/9eJaVlzQT/pKOTHiRhyiJSiYC2ZDatXs3My66dZ/71y64NKXt7qivY9z4Se/9u6JqAqtffWVEz+2J5dLd3UMsFv5BvoamFpY8v4w7Tq9g/pJlXDD3GMaWFtouK9L8mAmAkuoDeePdtRz1gZqU6vA7zYT/pCMTI8mDSBZd3T1Jv0ZYRCkT4d/rj8JAi2qOtBddXH0g/9ywzeuSfOm+J15g7vQsZlbmMnd6Fvcu+bvtklSapJIJgJyicnbtafeyJF/STETDSPJQVDmFzdt3Zaok34pSJrQhlaCiagLb6za/9++GbVspr6wa0XOzsqMxF6S3l3H+HKdncf6cQpY8vyz0x8CjKpVMAIypqGRdfUs6SvMNzUR0jCQPEs+mpyfaI1JRy4Q2pBLMPPhQtmx4l7rajXR2dvD8U49z9L+cOKLnjp08g1Ubwr8AW28vo6LIOSpcURQPfW8jylLJBMDY8ZPZ0BDOL89emonoSDUPURG1TARyjlS6xOJxvvKdH3HNJfPo6e7mU587i5rpI5tEOGZcFRuWhbvnDfDcijVsqW/ngdfq+9w/cduaUE8mjKpUMgHOYY+OrnCvbq6ZiI5U8xAVUctEIBtSpSUlA04YLC0pSXnbRx57PEcee/yon5dfWMw7EbjQ/eKbUpu8rNLDj5kApyHVHgv30iCaCX9KVyZSyUNURC0TgWxIpXI6dzrt7Qnkn1OFgF8zAc7yB0plmp8zocLFszlSIhITkVdEZIlX2wyaTrIHPKtDRZNmwrGnJ9t2CconopCJWHYuezt0Nf8o8XKy+deBtzzcXuDklE9i2069eLF6T+QzAUCWjtSq94Q+E+XjJ7NxZ6vtMlQGedKQEpHJwGeAu73YXlBlxeORWAJBDU8z0ZeO1KqoZCIrK0sPZ0eMVyNSNwNXAoMuniEiF4nIchFZvnTxgx69rL+MmXIwK96ps12G8oeb0UwAUDhhOms21g//iyrsbkYzoUIo5YaUiMwF6o0x/xjq94wxC4wxhxtjDj/21HmpvqwvlVdN4p26fa/tpKJFM9FXPK+ANp0zEmmaCRVmXoxIHQOcKiLrgUXAJ0Tkfg+2m3G/+O4V/OvHD+biz/1LUs8fU1FFbWP4L4ehhqWZSFA6voY3axs9q0kFkmZChVbKDSljzNXGmMnGmBrgLOAZY8y5KVdmwSdPO5Pr73ggpW3sNbkeVaOCSjPRV0l5Jdt3t3lUkQoizYQKs0BfImZX4w5u+Nq57G7a6cn2Zh/+YYpLy1LahsRintSiVDL8mIm8gkK2N2tDStnhx0yocPG0IWWMec4YM9fLbQ7lmcd+Q8+WVfzlUf+MELd2ZekZSmnQ0NTCGVfdGbiLXmomIDe/gObOQPfZfEkzMTI2MqH7ADtsZSKw3267GnfwytO/4+bTJ/PK07/zrLeRKskrZs/eDttlhM59T7xAY92m0F700gt+zQRAd4/uWLymmRiejUxk5+bR0qGfdxtsZSKwDalnHvsNp0yHGVX5nDId3/TAx0w7lFVvb7ZdRqg0NLWw5Pll3HF6BUueXxa4Hnim+DUTAG2E+3p7maaZGBkbmRAEJLC71sCymYlA/t/u7WWcfVgpAGcfVuqbHvjYiVN4c/Mu22UMKKiHAu574gXmTs9iZmUuc6dnaQ98AH7OBIBk+fOrRjMRXn7PhF9pJkbPn99uw+jtZYwtcq7hNbYo25Pexo1XzueKc+dSu/4dzj1+Dn/4/ejPzBhbNYl36/15mZggHgro7WWcP6cQgPPnFGoPfAB+zgRAV/egazBapZkIL79nwq80E6MXyItgvfbyX/nr1jYefLW2z/1jtv+Vz33pa0lv9+qf3pFqacSzc2jp8t+fNXHYc/6SZVww9xjGlhbaLmtYvb2MiiLnb1pRFH+vt/GNcz5luTr/8HMmAGIlVTQ0tVAxpsiT7XlBMxFufs+EH2kmkuO/Pf4IXHfHb22XMDQfHh/vO+zZFpgv3edWrGFLfTsPvNb3EiMTt60JRP2Z4vdMxLJzfHcdSs1EuPk9E36kmUhOIBtSftfS468/a28v4+EziwFn2PPMh4PR21h806W2S1AeKK4+kJXvvMSEilLbpQCaCaX600wkz39DJyHQE8ujvaPTdhnvGWrYU6lMKCguo7lNM6GUX2kmkueboRNjDBgDIrZLGZwxI1porWjidN7dupMDp1RloKjh2R72VMkJUyZKyitY96p/TsLQTARTmDLhN5qJ5PmmIdXWDXmmix7i/gyJMWSZLtpGMM1jTPVMVr37hG8aUraHPUeioamFi398PwuuPs/3w8iZEqZMdLS1UZqfnf6aRkgzEUxBycSupiag0nY1o6KZSJ5vGlK7OrNgTxt5MRAfBsQYQ1u3W+cwikrKaFqvq5uP5kOfeMqt9n4cYcpEW2sL04p1UU7NRGqCkom3N20laA0pW8KQCd80pAxCU2cM/DONImmFJaVsadSLtI70Qx/UU27TLUyZaK7fxMTyAttlWKeZSE1QMtETwEN7toQhEzrZPA3yCorY1e7PBQgzZTTL9esqzeHX2lDLjEljbZdhlWZCqb7CkgltSKVJa0+u7RKsGumH3vaKtEplimZCqb7CkgltSKVJp2TT5bMFCDNlNB96PeU2GnJLx1G73Z/XoMwEzYRSfYUpE76ZIxU28ZJxNDa3Mq6s2HYpGTea5fr1lNtoKJs8g9W1r3PYzGrbpVihmVCqrzBlQhtSaTKmZjavvP0Snzoieg2p0Xzog3DKrUpd1X5TefWFJs62XYglmgml+gpTJlJuSIlINXAfMB7oARYYY25JdbtBl5NXyN7maB7a8/uHPt00E/uKx7PpjvCJTJoJzYTqK0yZ8GKOVBfwTWPMQcDRwFdFZJYH2w20sRMm82Ztk+0ylB2aiX6adtTT3R3NjoUCNBMqxFJuSBljthpjVrg/NwNvAZNS3W7QFY8pZ/tuXZQzijQTfa1e9hwr77mKi0882HYpyhLNhAozT8/aE5Ea4EPA/w3w2EUislxEli9d/KCXL+tb7Vm6knPURT0Tm95czs4XH2bhv3+WmdUVtstRPhD1TKjw8WyyuYgUAY8Alxtjdvd/3BizAFgAcNfSdZGYLbGnW1eXiLKoZ6Jhay2v/O4WHv/+mWRlaRaUZkKFkyffbiKSjROO3xhjfu/FNsPAxKK9KGeURT0TjfWbef6/r+GeKz6jjSgFRCsTXT3Q0xPtq1tEScrfcOJcOfJ/gLeMMb9IvaTwiOcX07i71XYZKsOinonWlmb+8Mtv8+tvnhzJddTUvqKWifxx+7Gpvsl2GSpDvOgqHgOcB3xCRFa6t097sN3AK548k7W1222XoTIvspkwxvDColtY8NUTtBGlEkUqE1nxbHp69MhkVKQ8R8oY8zdAPKgldKpqDmL5qhc5ctYU26WoDIpqJowxPP/QHXxsfAfTJ4+zXY7ykahmQkWDTl5Io5LyCmp3hufQXkNTC2dcdadvLhSp/GXXjnomtq7mq6cdabuUjNFMKNVXFDOhDak0ysnNY09nzHYZnrnviRdorNuU8QtFRjGYQbT80Tv52mlH4EyHiQbNhFJ9RTET2pBKs1aTbbsET/ReqfuO0ysGvUJ3utgKphq59W+9ygEFzVRXldkuJWM0E0r1FdVMaEMqzUxWOBpSvVfqnlmZ+94VujPBZjDVyBhjeOWRW7jmrI/aLiWjNBNK9RXVTHi2IKctN146j5aW5n3uLyoq5upb7a+M2xPPZ8/edgrzg7umVO+H9OEznbOwzp9TyJkPL+OCuccwtrQwra/dN5ht3Lvk7/tcGXykGppauPjH97Pg6vPSXrdNmc7Eqj8t4opTDyUvNxydhpHQTASL3/cTYRDlTAR+RKqlpZmpF/5yn9tAobEhv2IyW3fss4BvoPR+SCuKnHZ3RVE8I72N3mCeP8f5MJ8/pzCl3kZUDodkMhN7mnex47Vn+OScqZ5v2880E8Hi9/1EGEQ5E4EfkRopWz2SsVNns3zt44E+Hfy5FWvYUt/OA6/V97l/4rY1Sbf6R2KoYI72dROHfucvyUwvye+8yMSqJffwsy993OvSfE8zEU46cpW8KGciMg2p3h5Jf+vuviytr1teOYF3Vrak9TXSbfFNl1p5XS+D6eXQb1ikmonX//IQYzu2Ul11gNel+Z5mIpxs7SfCIMqZCEVDqm7TOrq7u/vc19hQz42XzrPeiygsGUPt3u7hf1Htw6tg2jx2b0s6M9HctIMXfvvfHFbeytVfjt5olE2aieT5eT+hkpeRTJSNGfK5gZ8jBdDd3U1uxX59btlF5b45/t3RrZcKsMnWsXub0pWJtS//mT/e8i3mH1XEVWce7VG1KtM0E/7bTyi7UslE4EekioqK2bToWrKLyvvcH8srADoA57h3Y0M9r946f5/fyc9AjR2Sl4FXUYOxdezelnRlYunCHzEx1shj155KXk52pBbeDBvNhMNP+wll15CZOO8zQz438A2pq299kGu+OHfI49otLc2MP+t6civ26/P4loWXQ176/wR7TQ7d3T3EYqEYAAwcW8fubfE6EzvqNvPXu67lK58+lM8c9UHAWTtKBZdm4n1+2U8ou1LJhNVPhxdnSIy0FxGLxWhv2NjndzpbdlJUMW3UdY9WTslYdu7ew7iy4rS/VpBEZQ2b0fBjJlb+9hf86vITqdTPb9ppJvblx0yozAlCJqw2pLw4Q2KgXkRbQy11D36HNoFrvjiXxoZ6cru7icVijK9+f72b9orKjEwyLKk5hNfXvcRxh+mOKFHieh1hPJyQDL9lwhjDuOKcfRpRelgvPTQT+/JbJlRmBSEToTjW1NuL6L31dHWQlV/MpHNuZOqFvyS7qJzciv32OWMjU6qq92fVxl1WXtuvbC/pH3ZeZaKnu5uCmJ51mgmaifTy+35C7SsomQjFgd/E3gPA5vVricez2fHH26hva6WjeSd76zdierrYvH7tPj2OdCstH8eG3Xsz9npBoGvYpJdXmdi64W0OrRmbqbIjTTORXn7fT6h9BSUToWhIJarbtI7uri46mndC806q/vV6Gp66he2LfwoY4vFsOlt20lpRSVFRZg61SVYWbV0ZealAiOIaNrb0rp0zUCYanrgJ0901ZCa6Ojooyo1Zqj46NBOZk2omVGYEKROeNKRE5CTgFiAG3G2M+bEX201Gd3c32eWTiBWVAZBdUc2E835OZ8MmsrJgUs0M1t19GTcsXJKxmkSEDsnJ2Ov5nZdL+vuVXzLRu3ZOZ0d7n0xUnPx18iv3o71h47CZ6NYz9NJOM5E5XmRCpV+QMpFyQ0pEYsBtwCeBWmCZiCw2xrw53HOLiooHnDA4mh5A/200NtQjecXE8kvo6Wgd8XbSbU9Ptu0SfCPsa9j4KRONDfVkF5XT1dWZVCYm1Eznj/fezWePmTWq56nR0UwMzm+ZUJkRpEx4MSJ1JPC2MWYdgIgsAk4Dhg2IF2dC9N/GNV+cS8GJV9Dd3c22x36C6XIWWzM9XcSyc1N+vWRJbjEtre0UFdirwS+8XMPGp6fG+iYTvWvn9B7O6M2E6emivWEjsdjQh+1y8wtoLphM4+5WykoKUq5NDczLy1z4MA8QokyozAhSJrw4a28SsCnh37XufdaMr57KpJoZxOPZZOfkkp2TSywetzpxsHjiNNbX7bD2+mGVeGqsj/g+E7F4nEk1M0aUiWnHnMbjL67JQJUqVT7NA4QsEyo4MpEJLxpSAy0os8+kChG5SESWi8jypYszsyZHLK+ALQsvZ8vCy6lbdC3r7r6MdXdfZmXyoMRy2NvemfHX9aOGphbOuOrOlE9l9fGpsb7PRGIehsvEpKkH8NBLG2hpbc9IjVHkRSZ8nAcIWSaG32gOe9s7vCsygoKUCS8O7dUC1Qn/ngxs6f9LxpgFwAKAu5auS2n26lAr3SYeC8+H95b2L6qYZnVRtappB7P81WUccdAUazX4hVcLrPn41FjfZKJ55/Y+80t6MzGaPIgIR5/9Lb53/+3cdNEnUylTDcKLTPg4D5DhTIx0HwHJZWI4ldNms+Kd3zFr/wmebC+KgpQJLxpSy4AZIrI/sBk4Czjbg+0OaqiVbv16lkVp+Tg2NuikxsQewvwlyZ/K6vNTY0OXiYn7H8CapRNZvWEbM6dUpbw99T4vMuHzPECGM2F7HxGPZ9PZo2e7JitomUj50J4xpgu4FPgj8BbwsDHmjVS3GzbZubm0dutlNfr2ELKSPm491KmxtoU1E0eeMZ+L73iW3Xt0cVkveZEJP+cBwpsJlR5By4Qn60gZY54EnvRiW2G210T7jD0vewh+PzU2jJkoKCrmk5f8kKsX3sRtXz3Jdjmh4FUm/J4HCGcmlPeCmInQrWzuZ51ZebR3dJKbE801pbxcYM3LJRTUyFVOrmFNyXReX7eVg6fq/I9UeZUJzYMKiyBmQhtSGZRTWsGuPW1URrQhFYResxre4addyKLHr+d6bUilTDOhVF9BzEQgG1JerHRrQzyvmKbmPVSW+bvOdNFec/pkMhN5BYXUteqihV7QTKRHUPcRKpiZCGRDyuYyBqkoq5nFqvVPcsB+lbZLUSGT6UzExk1l07ZGqqvKMvq6So1EUPcRKpi8WJBTjVDV5Bre2NRkuwyrvFqMU9m1/1En8cgLb9suI/A0D0r1FcRMaEMqg/IKimja2227DKt8fAkLNQpjqybyVn2b7TICT/OgVF9BzIQ2pDKsQ/Jsl2CNzy9hoUapM16IMQZjdOHBZGgelOorqJkI5Byp4Qx1eQDbx86jvMvx+SUsQi0dmag46MOc/L3fcOO5RzGrZnxkl/VIlubBLj/vJ6IqqJkIZUNqqMsD2LanJ8d2CVYE4BIWoZaOTMw48gSqph3CHS//maaHHuPLnz6MEw+bmkqZkaF5sM/P+4koCnIm9NBehuWVVbG1YZftMjLO75ewUMkpGVvJESefzdHzf8H/vNjAMys32C4pEDQP4WYifewhOUHORChHpPysaNIBrK19hQkVpbZLyaggLrKmRq6otIwTLv4+P//Rv3HswZOJx3WdqaFoHsJt56a3OW5itL7jUxXkTGhDKsMmTZvFsr/9hWMPtV1JZt3z3S9y8Y/vZ8HV5/l+mFaNnohzQe45p3+FC39xO8fMruHkw/dn8rgxdgvzqcU3XUpDU4tmIqTam3cwcZou/jkaQc6EHtrLsIKSUrbtit5p40E8pVWN3n4HfoiDzv0BW2b8K5fc/iw3PbrMdkm+pZlQqq+gZiKUI1J+vjxAPJ5NG9GacJ54Suv8JcGYPBg2mcxE2TjnGnynXXM3f1v4I/bsbacwP9fz1wkyzYR9ft5PRFGQMxHKhpTfT13t7knv9v02PBrUU1rDJFOZ6D3E17u21MwT5nHl/9zErV896b3HbNBMqP78vp9IN82Ed0LZkPJ6fZDmpp0s+tm/M+/Kn1NUmvq1xdI9IpU4PGr7gxjkU1rDJNOZ6G00jZu8P6/nV/P6u/XMnlo1+sI9oplQ/aUzE0GgmfBOKOdI9a4P0v82UGhGYtlTDxHf9hovP7nIk/qkoIzG3a2ebKs/v60MG+RTWsPEZiY+dvYVXP4/S60t+6GZUAPx+34inTQT3gplQ8pLzU07Wb30UW763CRWL32Ull2NKW8zv2IytdubUi9uAH2HR+1/EJ9bsYYHXmvn8Nvq37s98Fo7z61YY7UulbzRZiKencOJX/0RF/5qJT9/aCm7WvZmqFKHZkKlW/9MtLWmp6PsFc2Et1I6tCciPwNOATqAd4AvGWOaPKjLN5Y99RCnzIDplfmcMmMPLz+5iE/Mm5/SNiunHcKKdxYxe9pEj6p0+HF4dPFNl1p5XVs0EwMrq5zIqV+7kdUrX+Lsn9zJVf96DB87uDrttWom7ItiJp5Z9TJ85CO2yxqQZsJ7qY5IPQ0cbIw5BFgDXJ16Sf7R28uYN8dZWG3enFJPRqVi8TgdXd7POA/68GhIaCYGkRWLcdBhx3Dadfdx27ObuOS2p+ns6k7rRY81E74QuUzUrl5JY7M/R6U0E95LaUTKGPOnhH++BHw+tXL8pbeXMbbQuRjr2MJsTplByqNS5ZUTWfuX5I7DDyXIK8OGhWZieCLC8fNvYMNrL/HJ797BlafN5qQjZ5KV5f1MA82EfVHMxIcq4ZFn/sGcmekfdR0tzYT3vDxr7/8BDw32oIhcBFwEcO43r+fYU+el/IKDnXXRvHO7J+uDrH3l77xS38ZDr9b23U7d31NqSMXicfZ2e38JjaAPj4ZQRjMx1FlIXq2Z42Umpsw+mnE1s/jNHx7mgecf5+cXHs/4sSWj2sZwNBO+E4lM7G7pZLd5d/QFZ4BmwnvDNqRE5M/A+AEeusYY87j7O9cAXcBvBtuOMWYBsADgrqXrPBnLH+rq3TcsXJLy9i/+6f0pb2MwHaILFAaVXzOR7jyA95koKC7h45//Ms2Nn+Xf/uubLP7+F6yuN6WSo5noa+WfFvGtQzs82b7yv2EbUsaYE4Z6XEQuAOYCx5t0TnYImfasfNo7OsnNybZdiholzYT3SsormHrsGaxYU8thPjwcooammVBRltKkBBE5Cfg2cKoxxp8z63wqr6yS+sYW22Uoj2kmRq93BGrqYcexePkmy9Uor2kmVNilOrvzVqAYeFpEVorInR7UFAljambz2rvbbJehvKeZSFJBUTFb9uhhvRDSTKhQS/WsveleFRI1sVh2WpZAUHZpJlJUNoXN25uYNG6M7UqURzQTKuwCfa29IF+9e9LUmfzj9w/w2Y/arkSFRZDz0Gv/o07k0Zd+xaWnHG67FBUCYciE8r9AN6SCfPXu7NxcdrWFd0TKb1cWj4Ig56HXuIn7seIPmb2ETKZoJjIvDJkIqzDlQa+1Z1FXLM92CWmTeGVxpUajLWcse9vDd+q4ZkKp94UpD4EekRqJoRZks91bae3Jsfr66ZJ4ZfH5S+xew0n15ec89Co78GheenMFx30oPFNrNBP+FYRMhE3Y8hD6htRQC7LZZsT71c39oO+Vxdu4/XfPsnJtbSiGcIPOz3noVV41mbVvPcdxH7JdiXc0E/4VhEyETf883Lvk75z/mY8E9lCfHtqzKL9iIhvqdtouw1O9PY3z5zhBOH9OIY88/SL1mzeGYghXpV/l5BpWbgrPGmuaCaXeN1Aeljy/jDseeS6wh/q0IWVRYWUNG+oabZfhqf5XFi/Ny+L46k4+OiWbJc8vY8euPZYrVH4nIrRmj6G9o9N2KZ7QTCj1vv55qCiKc/JUeOAPL3DH6RWBzESoD+3deOk8GhvqefXWvhdTjeUVkG+ppkQTps5ixQvPcuyh02yX4pn+VxbfvacNuro4sLKduQfkce+Sv+sVxi1q2l63Tx4AupobLFQzuMqDj+WZlS9z8pEH2C4lZZoJf0tLJvQqOIPqnwdwMlEQ6+pzqC9ImQh1Q6qlpZnxZ11PbsV+fe7fsvByyLP/1otKy1jb2Ga7DE8lXlm8oamFM6+8hYfPLKaiKE5DSxdnPhz8iYVBZiSLiV+8eZ/7N9x6fuaLGcL+hxzJ0/c/HoqGlGbC39KRidaGWqrKw3OyhJcS8wB9MwHOob6gZSL0h/ZisRjtDRv73DpbdvpiQbZYPE6HCeeEcxh4CHfu9KxAHgMPiyzJ2icP7Q0byRJ/fRXE49ns7M4nbNe31Uz4TzoyEevpoCAvnGdley0MmbA/LJNm46un7nNfe0Wlb05r3Ut415IaaAgXYOK2NYEatg2T0rEVTKqZsc/97WMrLFQztMKaD/L6uq3MnjbRdime0Uz4T5AyEUZhyEToG1K+J+G9SGv/IVylRmPGEZ/g0T/8LFQNKc2EUn2FIRP+Gs+PooJydu4O1hkKSmVCYckYaveEt6OhlAqHUI9IBeGClQUVk9i0bRvlJcGYVKeCLQiZSNRTMpltO3dTVV5iuxQVUkHLhPKfUDek/DIPaihjpxzEaxve5IMzJtkuRUVAEDKRqKRmNv/csEwbUiptgpYJ5T96aM+ycROn8MbmXbbLUMqXJk2dyQtrd9guQymlBqUNKcuyc3Pp7LFdhVL+VFJWwfpdGhCllH950pASkW+JiBERPV80CXu7wjuhtqGphTOuujNwS/6nSjPhna7C8aE6IUMzoZlQfQU9Eyk3pESkGvgksDH1cqKpU3Jtl5A29z3xQmAvRJkszYS39jviRP73/962XYZnNBNK9RX0THgxIvWfwJVAuJYgzqD2WD5dXd22y/Bc71W+g3ohyhRoJjw0oWYGy9bvtl2GJzQTmgnVVxgykVJDSkROBTYbY1Z5VE8kFVRUs6Fup+0yPNe79L9zIcpgLfmfLM2E90SEllgJnSHobGgmoiFslzZKpzBkYtiGlIj8WUReH+B2GnANcN1IXkhELhKR5SKyfOliPd00UfHEGazdHK4zk3p7GefPcdbHOn9OYWB7G/1pJjJv3KxjeHblOttlpEQzMaLthCIT2T3ttksIhLBkYtiGlDHmBGPMwf1vwDpgf2CViKwHJgMrRGT8INtZYIw53Bhz+LGnzvPyPQTepGkHsmxduEakwnAhysFoJjJv6gc/zJ9e3Wq7jJRoJqKTidzwXoveU2HJRNILchpjXgMqe//thuRwY0yDB3VFSk5uLh0hOGyRKAwXohwtzUT6xLNz2NmVhzEGCej1KTUTmgnVV1gyEeqVzYMiFs9mV3swdw6DCcOFKJW/5E06iLdrtzOjunL4X/YhzYRSfYUlE54tyGmMqdFeRnJEhA5t04aOZsJb+eOqWb9NrwIQZJoJFUa6srlPtEue7RKU8rWpHziMv7xeZ7sMpZTqQxtSPhEvqaRuRzjWylEqHXJy82jozLFdhlJK9aENKZ/IKxtHY3Or7TKU8rdYjq7Ro5TyFW1I+URuUTlbdrTYLkMp3wvqWXtKqXDShpRPVE45gDdrm2yXoZSv9eSW0KQjt8rnjOmxXYLKIG1I+UTp2Ere3aE7CKWGkj9mHE0te22XodSgOtrbKMnWhlSUaEPKJ7Kysmjr1uVwlVIqyLq7OinI0e/yKNGGlI90ZukSCEoNp6dHJ5srpfxDG1I+0pVTTEurXuxSqcGMmfIBVryzzXYZSin1Hm1I+Uh2XhHtnV22y1DKt/ILi9nTrhlRSvmHNqR8pGBcNe9s0asnKKWUUkGhDSkfKZ1Qwzt1zbbLUMq3smIx2jq6bZehlFLv0YaUj0yYMp2VG3baLkMp36qYUM3qOl24VinlH9qQ8pFYPE5rl/4vUWowIkJWlmZE+deOus1MqSi0XYbKoLjtArx046XzaGnZ99BYUVExV9/6oIWKRq8nlmu7BBUiYciEUl5KdyZ2N2xlWlVxyttRwRGqhlRLSzNTL/zlPvevu/syC9Ukp6Nb18hR3glDJhIZY2jXKVIqBWHLhLJPx8h9po08vbq9UoMQEV24VinlKyk3pETkMhFZLSJviMhPvSgqykomTWf1xnrbZagUaCbSS0Rsl6BGSTOhwiylQ3sichxwGnCIMaZdRCq9KSu6CiurebfuRQ6cUmW7FJUEzYRSfWkmVNilOiI1H/ixMaYdwBijQykpmjR1FivX77BdhkqeZiLN9nZBT0+P7TLUyEUrEzo1I3JSnWx+APAxEbkBaAO+ZYxZlnpZySkqKh5wwmBRUXDOoMgrLOLdZr0ERoBpJtIsq2AM25taqCovsV2KGplIZWJ37T+ZMVsH3aJk2IaUiPwZGD/AQ9e4zy8DjgaOAB4WkalmgNnSInIRcBHAud+8nmNPnZdK3QMKw+ncWVlZmFioTqYMHc2EXQce93muXvhD7vnGXNulKJdm4n2drbsZWzolra+h/GXYPbYx5oTBHhOR+cDv3UC8LCI9QAWwfYDtLAAWANy1dJ2OfQ6hrUsn0/qZZsKussqJ9FR9gP/752aOOnCS7XIUmgkVbanOkXoM+ASAiBwA5AB61d0UdWbpopwB9hiaibQ76oyLueHhl3WuVDA8hmZChViqDal7gKki8jqwCLhgoOFaNTrZZRPZ2rDLdhkqOZqJDIjHs6n5+Jk88X9rbJeihqeZUKGW0mQcY0wHcK5HtShX4bhq1tWtYUJFqe1S1ChpJjKnbPx+7K1babsMNQzNhAo7Xdnch8onT2ftVh2RUkqpoDF6uDlytCHlQ/GcXFr2dtouQyml1Cjl0267BJVh2pDyofLKCaxraLVdhlK+tmdnPRXFemKG8pfsmJ51HTVWFiyqKM6x8bKBUl9VDYW6qFtUaCZGr7anmQ/Ong2FwV1cVA0uqJmYMGGifneHzTBrO0pYTp4QkYvcNUhCIUzvJ0zvJUjC9ncP0/sJ03sJkrD93cP0foL8XsJ0aO8i2wV4LEzvJ0zvJUjC9ncP0/sJ03sJkrD93cP0fgL7XsLUkFJKKaWUyihtSCmllFJKJSlMDalAHlsdQpjeT5jeS5CE7e8epvcTpvcSJGH7u4fp/QT2vYRmsrlSSimlVKaFaURKKaWUUiqjAt2QEpEviMgbItIjIof3e+xqEXlbRFaLyIm2ahwtETnJrfltEbnKdj2jJSL3iEi9e4HS3vvKReRpEVnr/rfMZo1hppnwH82EXZoJ/wlbJgLdkAJeB04HlibeKSKzgLOADwAnAbeLSCzz5Y2OW+NtwMnALGCe+16CZCHO3zzRVcBfjDEzgL+4/1bpoZnwn4VoJmzSTPjPQkKUiUA3pIwxbxljVg/w0GnAImNMuzHmXeBt4MjMVpeUI4G3jTHr3CumL8J5L4FhjFkK7Ox392nAve7P9wKfzWRNUaKZ8B/NhF2aCf8JWyYC3ZAawiRgU8K/a937/C6odQ+nyhizFcD9r14/IfOC+tkKat3D0UzYF9TPVlDrHk5gM2HlWnujISJ/BsYP8NA1xpjHB3vaAPcF4fTEoNatMkgzEYi6VQZpJgJRd2j5viFljDkhiafVAtUJ/54MbPGmorQKat3D2SYiE4wxW0VkAlBvu6Ag00wEou7haCY8pJkIRN3DCWwmwnpobzFwlojkisj+wAzgZcs1jcQyYIaI7C8iOTgTIRdbrskLi4EL3J8vAAbrIar00Uz4i2bCPs2EvwQ3E8aYwN6Az+G0ztuBbcAfEx67BngHWA2cbLvWUbynTwNr3NqvsV1PEvU/CGwFOt3/N18GxuKchbHW/W+57TrDetNM+O+mmbD+99dM+OwWtkzoyuZKKaWUUkkK66E9pZRSSqm004aUUkoppVSStCGllFJKKZUkbUgppZRSSiVJG1JKKaWUUknShpRSSimlVJK0IaWUUkoplSRtSCmllFJKJen/A+7lcuUFk0HKAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x216 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# http://rasbt.github.io/mlxtend/user_guide/plotting/plot_decision_regions/\n",
    "import matplotlib.gridspec as gridspec\n",
    "gs = gridspec.GridSpec(1, 3)\n",
    "fig = plt.figure(figsize=(10,3))\n",
    "\n",
    "ax = plt.subplot(gs[0, 0])\n",
    "fig = plot_decision_regions(X_test, y_test,  clf=lr_model, legend=2)\n",
    "plt.title('Logistic Regression')\n",
    "\n",
    "ax = plt.subplot(gs[0, 1])\n",
    "fig = plot_decision_regions(X_test, y_test,  clf=knn, legend=2)\n",
    "plt.title('KNN, k=5')\n",
    "\n",
    "ax = plt.subplot(gs[0, 2])\n",
    "fig = plot_decision_regions(X_test, y_test,  clf=naive_bayes, legend=2)\n",
    "plt.title('Naive Bayes')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "baf9b80f",
   "metadata": {},
   "source": [
    "### 2.2 Handwriteen Digits Classification\n",
    "\n",
    "(b) Comment on the performance of the classifier and give your explanation why some of them perform better than the others.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "3b3f7a75",
   "metadata": {},
   "outputs": [],
   "source": [
    "import scipy.io as spio\n",
    "import random\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.neural_network import MLPClassifier\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.metrics import confusion_matrix\n",
    "from sklearn.metrics import accuracy_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "05d22b38",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import Data\n",
    "mnist = spio.loadmat('data/mnist_10digits.mat',squeeze_me=True)\n",
    "\n",
    "X_train = mnist['xtrain'] # each row is an image\n",
    "X_test = mnist['xtest']\n",
    "y_train = mnist['ytrain'].ravel()\n",
    "y_test = mnist['ytest'].ravel()\n",
    "\n",
    "# Standardize X_train and X_test\n",
    "X_train = X_train/255\n",
    "X_test = X_test/255\n",
    "\n",
    "# Downsample\n",
    "seed = 4624\n",
    "index = random.sample(range(60000), 5000)\n",
    "\n",
    "X_train_down = X_train[index]\n",
    "y_train_down = y_train[index].ravel()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7e6c0ec1",
   "metadata": {},
   "source": [
    "### Part A\n",
    "Report confusion matrix, precision, recall, and F-1 score for each of the classifiers. For precision, recall, and F-1 score of each classifier, we will need to report these for each of the digits. So you can create a table for this. For this question, each of the 5 classifier, KNN, logistic regression, SVM, kernel SVM, and neural networks, accounts for 10 points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "a83d5690",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best n_neighbors: 1\n"
     ]
    }
   ],
   "source": [
    "# Tune KNN using GridSearchCV \n",
    "# https://medium.datadriveninvestor.com/k-nearest-neighbors-in-python-hyperparameters-tuning-716734bc557f\n",
    "\n",
    "n_neighbors = list(range(1,30))\n",
    "\n",
    "#Convert to dictionary\n",
    "hyperparameters = dict(n_neighbors=n_neighbors)\n",
    "\n",
    "#Create new KNN object\n",
    "knn_2 = KNeighborsClassifier()\n",
    "\n",
    "#Use GridSearch\n",
    "clf = GridSearchCV(knn_2, hyperparameters, cv=10)\n",
    "\n",
    "#Fit the model\n",
    "best_model = clf.fit(X_train_down,y_train_down)\n",
    "\n",
    "#Print The value of best Hyperparameters\n",
    "print('Best n_neighbors:', best_model.best_estimator_.get_params()['n_neighbors'])\n",
    "\n",
    "# Current results = 1; previous results = 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "004c28f4",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "KNN, k=3\n",
      "Training Accuracy: 0.9363833333333333\n",
      "Test Accuracy: 0.9366\n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "           0       0.95      0.99      0.97       980\n",
      "           1       0.90      1.00      0.94      1135\n",
      "           2       0.97      0.91      0.94      1032\n",
      "           3       0.92      0.93      0.92      1010\n",
      "           4       0.95      0.91      0.93       982\n",
      "           5       0.93      0.91      0.92       892\n",
      "           6       0.96      0.97      0.97       958\n",
      "           7       0.94      0.94      0.94      1028\n",
      "           8       0.97      0.87      0.92       974\n",
      "           9       0.91      0.93      0.92      1009\n",
      "\n",
      "    accuracy                           0.94     10000\n",
      "   macro avg       0.94      0.94      0.94     10000\n",
      "weighted avg       0.94      0.94      0.94     10000\n",
      "\n",
      "[[ 966    1    1    0    0    3    7    1    1    0]\n",
      " [   0 1130    1    1    1    1    1    0    0    0]\n",
      " [  16   32  937    9    5    1    3   21    8    0]\n",
      " [   3   11   10  936    1   26    0   10    8    5]\n",
      " [   1   22    0    0  896    0    9    1    2   51]\n",
      " [   6    8    1   33    3  815   11    3    3    9]\n",
      " [   8    5    0    1    6    5  933    0    0    0]\n",
      " [   0   35    5    1    4    0    0  964    0   19]\n",
      " [  15    9    9   27   11   28    2   12  847   14]\n",
      " [   6    8    2   10   20    0    2   19    0  942]]\n",
      "Logistic Regression\n",
      "Training Accuracy: 0.9392666666666667\n",
      "Test Accuracy: 0.9256\n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "           0       0.95      0.97      0.96       980\n",
      "           1       0.96      0.98      0.97      1135\n",
      "           2       0.93      0.90      0.91      1032\n",
      "           3       0.90      0.92      0.91      1010\n",
      "           4       0.94      0.94      0.94       982\n",
      "           5       0.90      0.87      0.88       892\n",
      "           6       0.94      0.95      0.95       958\n",
      "           7       0.93      0.92      0.93      1028\n",
      "           8       0.88      0.88      0.88       974\n",
      "           9       0.91      0.92      0.91      1009\n",
      "\n",
      "    accuracy                           0.93     10000\n",
      "   macro avg       0.92      0.92      0.92     10000\n",
      "weighted avg       0.93      0.93      0.93     10000\n",
      "\n",
      "[[ 955    0    2    4    1   10    4    3    1    0]\n",
      " [   0 1110    5    2    0    2    3    2   11    0]\n",
      " [   6    9  930   14   10    3   12   10   34    4]\n",
      " [   4    1   16  925    1   23    2   10   19    9]\n",
      " [   1    3    7    3  921    0    6    5    6   30]\n",
      " [   9    2    3   35   10  777   15    6   31    4]\n",
      " [   8    3    8    2    6   16  912    2    1    0]\n",
      " [   1    7   23    7    6    1    0  947    4   32]\n",
      " [   9   11    6   22    7   29   13   10  855   12]\n",
      " [   9    8    1    9   21    7    0   21    9  924]]\n",
      "SVM\n",
      "Training Accuracy: 0.9142\n",
      "Test Accuracy: 0.9121\n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "           0       0.94      0.98      0.96       980\n",
      "           1       0.95      0.98      0.97      1135\n",
      "           2       0.90      0.90      0.90      1032\n",
      "           3       0.87      0.90      0.89      1010\n",
      "           4       0.89      0.92      0.91       982\n",
      "           5       0.87      0.86      0.87       892\n",
      "           6       0.94      0.95      0.94       958\n",
      "           7       0.91      0.92      0.91      1028\n",
      "           8       0.91      0.84      0.88       974\n",
      "           9       0.91      0.86      0.88      1009\n",
      "\n",
      "    accuracy                           0.91     10000\n",
      "   macro avg       0.91      0.91      0.91     10000\n",
      "weighted avg       0.91      0.91      0.91     10000\n",
      "\n",
      "[[ 959    0    2    2    0   11    5    1    0    0]\n",
      " [   0 1116    4    5    0    4    1    3    2    0]\n",
      " [  14   13  931   17   16    2   15    9   14    1]\n",
      " [   6    3   19  908    2   31    2   11   23    5]\n",
      " [   2    4    8    1  908    1   13    8    1   36]\n",
      " [  12    3    5   58   10  766   11    4   18    5]\n",
      " [   8    4   13    0    9   14  906    0    4    0]\n",
      " [   2   13   23    7    7    1    0  942    3   30]\n",
      " [  10    8   15   29   17   43    9   10  820   13]\n",
      " [   5    9    9   11   46    6    0   45   13  865]]\n",
      "Kernel SVM\n",
      "Training Accuracy: 0.9526\n",
      "Test Accuracy: 0.9545\n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "           0       0.96      0.98      0.97       980\n",
      "           1       0.98      0.99      0.98      1135\n",
      "           2       0.96      0.94      0.95      1032\n",
      "           3       0.94      0.96      0.95      1010\n",
      "           4       0.93      0.96      0.95       982\n",
      "           5       0.95      0.95      0.95       892\n",
      "           6       0.96      0.96      0.96       958\n",
      "           7       0.96      0.94      0.95      1028\n",
      "           8       0.94      0.94      0.94       974\n",
      "           9       0.94      0.92      0.93      1009\n",
      "\n",
      "    accuracy                           0.95     10000\n",
      "   macro avg       0.95      0.95      0.95     10000\n",
      "weighted avg       0.95      0.95      0.95     10000\n",
      "\n",
      "[[ 964    0    2    0    0    8    4    1    1    0]\n",
      " [   0 1123    2    3    0    1    3    0    3    0]\n",
      " [   8    0  970   10    9    0    7    8   20    0]\n",
      " [   0    1    7  968    0   10    1   10   11    2]\n",
      " [   1    0    5    0  941    0    5    3    2   25]\n",
      " [   3    1    3   20    4  843    8    2    5    3]\n",
      " [   9    3    2    0    8   10  924    0    2    0]\n",
      " [   1   12   18    3    7    0    0  962    3   22]\n",
      " [   3    1    3   13    7   16    6    4  918    3]\n",
      " [  10    6    1    9   31    2    0    9    9  932]]\n",
      "Neural Networks\n",
      "Training Accuracy: 0.95985\n",
      "Test Accuracy: 0.9572\n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "           0       0.96      0.98      0.97       980\n",
      "           1       0.98      0.98      0.98      1135\n",
      "           2       0.96      0.95      0.95      1032\n",
      "           3       0.95      0.96      0.96      1010\n",
      "           4       0.95      0.96      0.96       982\n",
      "           5       0.97      0.93      0.95       892\n",
      "           6       0.95      0.96      0.96       958\n",
      "           7       0.97      0.93      0.95      1028\n",
      "           8       0.95      0.96      0.95       974\n",
      "           9       0.94      0.95      0.95      1009\n",
      "\n",
      "    accuracy                           0.96     10000\n",
      "   macro avg       0.96      0.96      0.96     10000\n",
      "weighted avg       0.96      0.96      0.96     10000\n",
      "\n",
      "[[ 963    1    0    0    0    4    6    1    3    2]\n",
      " [   0 1112    3    4    0    0    3    2   11    0]\n",
      " [  12    2  976    7    7    0   11    9    8    0]\n",
      " [   0    0   10  973    0    4    2    8   12    1]\n",
      " [   1    0    3    0  944    0    9    2    2   21]\n",
      " [   8    1    1   20    2  829   12    3    9    7]\n",
      " [   9    3    3    0    7    9  924    0    3    0]\n",
      " [   2    9   19    7    7    0    0  960    0   24]\n",
      " [   5    1    3    5    6    6    9    5  931    3]\n",
      " [   7    6    1   11   16    2    1    4    1  960]]\n"
     ]
    }
   ],
   "source": [
    "# from plot_classifier_comparison.py\n",
    "names = ['KNN, k=3', 'Logistic Regression', 'SVM', 'Kernel SVM', 'Neural Networks']\n",
    "\n",
    "classifiers = [KNeighborsClassifier(3), \n",
    "               LogisticRegression(solver='newton-cg'),\n",
    "               SVC(kernel=\"linear\"),\n",
    "               SVC(kernel='rbf'),\n",
    "               MLPClassifier(alpha=1, hidden_layer_sizes = (20, 10))\n",
    "              ]\n",
    "\n",
    "# Iterate over classifiers\n",
    "for i in range(len(names)):\n",
    "    clf_name = names[i]\n",
    "    clf = classifiers[i]\n",
    "       \n",
    "    print(clf_name)\n",
    "    \n",
    "    # Fit the model\n",
    "    if (clf_name == 'KNN, k=3') | (clf_name == 'SVM') | (clf_name == 'Kernel SVM'):\n",
    "        clf.fit(X_train_down, y_train_down)\n",
    "    else:\n",
    "        clf.fit(X_train, y_train)\n",
    "        \n",
    "    # Predict Test Records \n",
    "    y_pred_train = clf.predict(X_train)\n",
    "    y_pred_test = clf.predict(X_test)\n",
    "    \n",
    "    # Accuracy \n",
    "    acc_train = accuracy_score(y_train, y_pred_train)\n",
    "    acc_test = accuracy_score(y_test, y_pred_test)\n",
    "    print('Training Accuracy: ' + str(acc_train))\n",
    "    print('Test Accuracy: ' + str(acc_test))\n",
    "    \n",
    "    # Scores \n",
    "    # https://scikit-learn.org/stable/modules/generated/sklearn.metrics.classification_report.html\n",
    "    print(classification_report(y_test, y_pred_test))\n",
    "    \n",
    "    # Confusion Matrix\n",
    "    # https://scikit-learn.org/stable/modules/generated/sklearn.metrics.confusion_matrix.html\n",
    "    print(confusion_matrix(y_test, y_pred_test))\n"
   ]
  }
 ],
 "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.9.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}