{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import turicreate as gl\n", "import pandas as pd\n", "import numpy as np\n", "import csv\n", "import math\n", "import seaborn as sns; sns.set()\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline \n", "import scipy\n", "import random\n", "import dask.dataframe as dd\n", "from sklearn import decomposition\n", "from mpl_toolkits.mplot3d import Axes3D\n", "from sklearn.cluster import KMeans" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Single-cell RNA-seq analysis\n", "**In this problem, we analyze single-cell RNA-seq data and determine structure in this high-dimensional data. The data set consists of 272 cells, each row corresponds to the RNA-seq measurements of a particular gene. Each entry corresponds to the normalized transcript compatibility count (TCC) of an equivalence class. An equivalence class is a set of short RNA sequences. The TCC counts the number of reads of sequences which are compatible with each equivalence class for a given cell. The entries have been normalized so that each row in the matrix sums to 1.**" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "filename = \"Trapnell.csv\"\n", "# Trapnell = pd.read_csv(filename, nrows = 20).T\n", "# Trapnell = np.genfromtxt(filename,delimiter=',',dtype=float)\n", "with open(filename,'r') as dest_f:\n", " data_iter = csv.reader(dest_f,\n", " delimiter = \",\")\n", " data = [data for data in data_iter]\n", "Trapnell = np.asarray(data, dtype = float)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### (a) Determine cell clusters by applying k-means clustering to the data. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### 1. Using PCA to reduce demensions from 1,065,024 to 2" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "pca = decomposition.PCA(n_components=2)\n", "Trapnell_pca = pca.fit_transform(Trapnell)" ] }, { "cell_type": "code", "execution_count": 55, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZoAAAEeCAYAAACzJ9OtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XtYVNX+P/D3cAfRKBogMT2GggYi3shbkH1FlFAU0zxeMC0vp5tZh45m2jE9ah1PZnq0UtPUstQQwgQVUjPFGz9LU7yWqSgIKKIo11m/PzgzMTAzzAyz58b79Tw+j7P3nr0/ewPrs/daa68lE0IIEBERScTB0gEQEZF9Y6IhIiJJMdEQEZGkmGiIiEhSTDRERCQpJhoiIpIUEw1Z3O+//460tDS1ZdXV1di4cSPu3btnoajsT1JSEgYMGICQkBD07t0bly5dMvt1vnr1KoKCgur969ChA8LCwhATE4MPPvgAxcXFGr8vhMAPP/yAl19+GU8//TRCQkLQs2dPTJ48Gfv27Wvw+J999hmCgoLwxBNPoKKiwtSnR1o4WToAatrOnDmDZ599Fn/9618xaNAg1fI333wTaWlpGDJkiAWjsx8XL17ErFmz4OnpidGjR8PBwQEtW7a02HX29/fHsGHDVJ+FELh37x4OHTqENWvWYM+ePdi6dSuaNWum2qakpARvvfUW9uzZA29vb/Tp0wc+Pj7Iy8vDDz/8gH379uGFF17AW2+9pfW43333Hdzd3VFcXIxdu3YhNjZW0vOkGkw0ZFG3b99GZWVlveVFRUUWiMZ+5eTkQKFQYPTo0Zg+fbpquaWus7+/P1599dV6yxUKBSZPnoz9+/fjiy++wEsvvQSgJhG9/vrrOHDgAJ577jnMnDkT7u7uqu8VFhZiwoQJWLNmDVq1aoXRo0fX2/evv/6K8+fPY+rUqVizZg22bNnCRGMmrDojagKU1UQPPvighSPRzcHBARMnTgQA/Pjjj6rlSUlJOHDgAPr27Yu5c+eqJRkAePjhh7F06VLIZDJ8+umnGm9ekpOTAQDR0dHo2bMnDh8+jCtXrkh4NqTEREOoqqrC8uXLMXjwYHTu3Bnh4eF44YUXkJWVVW/bmzdvYsGCBXj66acRGhqK6OhoLFmyBKWlpWrbnTt3DomJiYiMjERISAi6du2KUaNGYefOnaptli1bhoSEBADA+vXrERQUhMOHDyMoKAhHjhwBAPTo0QPjxo1TfaeiogKffvopYmJi0KlTJ/Tq1QtvvvlmvQJj2bJlCAoKQlZWFkaMGIGQkBBER0fXi1NJedyUlBRs3rwZgwYNQqdOnTBw4ECkpKQAADIzMxEfH4/OnTsjOjoaX375Zb395Obm4t1330X//v3RqVMndOnSBfHx8di0aZPadklJSQgKCsL+/fuxfPlyPPnkk+jSpQuee+457NmzR+vPqq7s7Gy88sor6Nu3L0JCQtCjRw9MmDABhw4dUm3z9NNPY+bMmQCAhQsXIigoSHV9zH2d9eHr6wsAuHXrlmrZ1q1bAQBTp06FTCbT+L3HHnsMs2fPxuzZs1F3ZK2qqirs2LEDDz/8MDp27IiYmBgIIVT7JWk5/vOf//ynpYMgy5o7dy4+//xzBAQEICoqCq1bt8ZPP/2Eb7/9Ft27d0erVq0AAAUFBRg5ciT279+Pxx9/HP/3f/+H+/fvIzk5GSdOnEBsbCwcHBxw4sQJ/PWvf8WVK1fQr18/9O7dGw888AAOHTqE77//HiEhIWjbtq3q+GfOnEHnzp0xYsQIhIeHo0WLFsjNzcWdO3cwadIk9O7dGx07dkRlZSUmTZqELVu24NFHH8XAgQMhl8uxc+dObNu2DZGRkfD29gYAHDlyBEeOHMFPP/0EuVyOqKgo+Pj4IDo6WuM1yM3NxbZt23D16lVs374dkZGR6NSpE44dO4a0tDTcu3cP77//Prp3744ePXrgxIkT2LlzJzp27IjHHnsMQE1D97PPPotff/0Vffr0wZNPPglfX18cPXoUGRkZePDBBxEaGgqgpiorMzMTFy9exO7duzFgwAB06NAB2dnZ+Pbbb/HII4/g8ccf1/lzy8jIwOTJk3H79m30798fPXr0gIuLCw4cOIDU1FQ8/fTTkMvlEELA1dUVv//+O/r27YshQ4YgPDwc/v7+Zr/OJSUlWL9+Pfz9/REfH69xm59//hnff/89goKCMGzYMNy9exfvvfcePDw8MGfOHDg4aL8/Dg0NxWOPPQZHR0e15Xv37sWWLVsQHx+PyMhItGrVCuvWrcPvv/+O8ePH69wnmYCgJu3OnTuiQ4cOYsyYMWrLT5w4IQIDA8Wrr76qWpaYmCgCAwPF2rVr1badPXu2CAwMFDt37hRCCDFx4kTx+OOPiwsXLqht9/3334vAwEDxxhtvqJYdOnRIBAYGivnz56ttO3bsWBEYGChu376tWrZq1SoRGBgoPvjgA6FQKNRiDQ4OFsOHD1ct+/jjj0VgYKAYPny4qK6ubvA6KOPo2LGjOHnypGr5119/LQIDA0VgYKDYs2ePavnhw4dFYGCgmDZtWr3rcODAAbV9//LLLyIwMFA899xzqmXffvut6njHjx9XLb906ZLo3r276N69u9q5axIdHS3Cw8NFQUGB2vLPPvtMBAYGiv/85z/1jlf3Z2fu63zlyhURGBgoxo4dq3F9WVmZGDFihAgMDBQbNmwQQghx4cIFERgYKAYPHtzg/rV59dVXRWBgoMjOzlYte/nll0VgYKDIzMw0er+kH3YGaOIUCgWEELh27RquX7+ORx55BADQqVMnZGRkwM/PD0BNVcru3bvxl7/8Bc8//7zaPqZMmYIHH3wQcrkcAPD8889j+PDhCAgIUNvuiSeeAGB8A/TWrVvRvHlzvP7662rVJ8oqrtTUVJw/fx7t27dXrYuKijLobrVbt24ICQlRfe7atSsAoG3btnjqqadUyzt37gyg5klIaciQIQgNDUXv3r3V9hkaGgo3NzeN5x0TE4OwsDDV5zZt2mDMmDFYuXIl9u7dq7U3mEKhwJtvvgkXFxc8/PDDauts4Trn5uZi2bJlqs9CCBQVFWH//v3Izc1F165d8dxzzwGoeQoCoNYDzRAlJSXYs2cP/P390aVLF9Xy2NhY7N69G1u2bMHTTz9t1L5JP0w0TVyLFi0QExOD77//HlFRUejSpQsiIiLQr18/tGvXTrXd5cuXce/ePbVCUcnf31+tJ9OTTz4JoKaq7cyZM7h8+TJ+//13ZGdnA6h5R8ZQpaWl+P333yGXy/HJJ5/UW19YWAigpkqqdgHo7+9v0HHatGmj9lnZ6KysPlRydXUFALV3Mbp3747u3bujuLgYOTk5qvP++eefUV5ervG8w8PD6y1TVq+dOXNGa6JxcHBAVFQUgJpC+/z587h8+TIuXLiAw4cPA6hJRoYy13XOzc3F8uXL1c6nWbNmaNu2LUaNGoWEhAQ4OzsDALy8vAD8mXAMlZaWhoqKCsTExKglzn79+sHT0xM//vgjbty4AR8fH6P2Tw1joiG8//77CAkJQVJSkqrOffHixQgJCcH8+fPRsWNH3L59GwDg6enZ4P6uX7+OefPm4YcffoAQAg4ODvjLX/6Cbt264fTp00bFePfuXQA1yat2AVWXMk4lNzc3g45TtzeTkouLS4PfvX37NhYuXIjt27ejsrISMpkM/v7+6Nmzp9bzVjZ816Z8QlGeszZnz57F/PnzVQ36zs7OCAgIQEhICC5dulSvQVwf5rrO4eHh2LBhg17b+vn5wdnZGdeuXUNlZaUqAWmSl5cHT09Ptd9TZW+zVatWYdWqVRq/t23bNkyZMsWAMyBDMNEQnJ2dMXHiREycOBHXrl3DgQMHkJ6ejp9++glTpkxBZmamqtpCW2+ie/fuwcPDA0IITJ48GRcuXMCUKVPQv39/tG/fHm5ubigsLMSWLVuMitHDwwNAzVODpt5e1iAxMRH79u3DqFGjEBcXh8DAQFWBl5qaqvE7ZWVl9ZbduXMHgO6uyHfv3sXEiRNx584d/OMf/0Dv3r3x2GOPwcXFBb/88gu2b99u1DlY43V2d3dH9+7dkZWVhePHj2t8ClSaM2cODhw4gE8++QRPPvkkrly5gv/3//4ffH191ao+lUpLS7F9+3Zs3boVkydP1tqjjRqHiaaJu3LlCrZs2YIuXbqgX79+aNmyJUaMGIERI0Zg/PjxOHToEK5evYq2bdvC2dkZJ06cqLeP/Px8REREYOTIkRgzZgzOnTuH6Ohoteo0oObtdABqd9r6/mE3b94cLVu2xIULF1BWVlbvDjo5ORlXrlzBsGHD6lVzmUNJSQn27duHkJAQzJ07V23d1atXUV5ervEJ4+TJk/V6aB0/fhzAn1Vomhw6dAiFhYWqG4TaNF1nfVnrdR42bBiysrLwySefaE00Fy5cwMGDB+Hu7q5qi1E+zYwaNUr18mddJ0+exB9//IHDhw+jZ8+e0pxAE8c+fU2cm5sbVq1ahaVLl6q1N1RUVKCgoAAuLi6Qy+VwdXVFdHQ0Ll68WO+pRFmX36tXL1UVU92G6OLiYnzwwQcAat5pUHJyqrnXqfuCnbJ6pPbyYcOGobi4GIsXL1Zrf7hw4QLee+89rF27VlWfb27Ozs5wcHBASUmJ2nUsKyvDvHnzANQ/RwDYvHmzKjEANeO+bdiwAb6+vujbt6/W4ynbiOpe52vXrqmqvGpfZ11x143NGq/zkCFD0KVLFxw4cABz5sxBeXm52vpLly7h5ZdfRmVlJV5++WXVk+R3330HABg8eLDWfSuHwjH2aZsaxieaJk4ul2P8+PFYu3YtYmNjERkZCQcHB+zfvx8XL17ESy+9pPqjfeutt5CdnY133nkHO3fuRPv27XHy5EkcPXoU/fv3R0xMDBQKBUJDQ3Hs2DGMHj0aXbt2xa1bt5CRkYGKigq4u7urvYinbKNIS0uDh4cHhg0bhvbt26uWv/322+jTpw8SEhIwefJk/PTTT9iwYQOys7MRHh6OkpISpKen4/79+/j3v/+tVxuSFNzd3REVFYWdO3dixIgR6NOnD+7du4c9e/agsLAQDzzwAO7cuQOFQqHWO0smk2HkyJEYOHAghBDYtWsXysrK8MEHH6iSiSbdunWDv78/UlJScOvWLXTo0AHXr19HZmYmXF1dIZPJtA5MWZutXGeZTIYVK1bgxRdfxDfffIPdu3fjqaeewoMPPog//vgD+/btQ2VlJcaOHavqFXns2DFcvnwZXbp0waOPPqp138OGDcPHH3+M3bt34/bt23jggQfMdFZNB59oCImJifjnP/8JT09PbNu2DZs3b0azZs2waNEiTJs2TbWdr68vtmzZgueeew5nz57F+vXrce3aNfztb3/DkiVLANT0HlqxYgXi4+Nx9epVbNiwAceOHUNERAS+/fZb9OnTB5cuXcLly5cB1PRWUnaj/fLLL1VVc1OnTkXnzp1x4MABVVuBm5sb1q9fj1dffRXl5eX46quvsG/fPnTt2hXr16/XeddqDgsWLMD48eNx584dbNy4Efv370enTp2wadMmDB06FGVlZaoeYUpTp07FuHHjsGfPHuzcuROdO3fGxo0bERkZqfNYHh4eWLt2LQYMGIBTp05h48aNOH36NIYMGYLvvvsOHTp0wLFjxxp8Q9+WrvNDDz2ETZs2YcGCBWjXrh0OHjyI9evXIzs7G3379sWaNWswe/ZsVXWs8mmmoQFD/fz80Lt3b5SXl6u+Q6YlE8ZU5BJRoyQlJWHmzJmYOXNmvfeSiOwNn2iIiEhSTDRERCQpJhoiIpIU22iIiEhSfKIhIiJJMdEQEZGkmtwLm7dulUKhqF9b6O3tiaIi3YMYWitbjd1W4wYYuyXYatyA7cbu7W2aF3ObXKJRKITGRKNcZ6tsNXZbjRtg7JZgq3EDth17Y7HqjIiIJMVEQ0REkmKiISIiSTHREBGRpJhoiIhIUkw0REQkKYskmtTUVMTExGDAgAEa5yXPyclBfHw8oqOjMWvWLNVMgTdu3MDkyZMxdOhQjBo1ClevXjV36EREZCCzJ5r8/HwsWbIEX331FZKTk/HNN9/gwoULatskJiZizpw52LlzJ4QQ2Lx5M4CaGR779euH5ORkxMXFYfHixeYOn4iIDGT2RHPw4EH07NkTXl5e8PDwQHR0NNLT01Xrc3NzUVZWhrCwMABAfHw80tPTcfPmTZw5cwajRo0CAAwfPhyvv/66ucMnIiIDmX1kgBs3bkAul6s++/j4qKbv1bReLpcjPz8fV65cQcuWLbFo0SIcO3YMcrkcs2fPNvj4uoZUkMubG7w/a2Grsdtq3ABjtwRbjRuw7dgby+yJRqFQqOb0BgAhhNpnbeurqqpw+vRpvPrqq5g5cya2bNmCGTNmYMOGDQYdv6jorsahIOTy5igouGPEGVmercZuq3EDjN0SbDVuwHZjN1VyNHvVmZ+fHwoKClSfCwoK4OPjo3V9YWEhfHx8IJfL0axZM/Tr1w8AEBsbq/YkRERE1snsiaZ3797IysrCzZs3cf/+fezatQsRERGq9f7+/nB1dUV2djYAICUlBREREWjdujX8/Pywb98+AMCePXsQHBxs7vCJiMhAZk80vr6+mD59OhISEjB06FDExsYiNDQUkyZNwsmTJwEAixcvxsKFCzFw4EDcu3cPCQkJAIBly5Zh9erViI2Nxfr167FgwQJzh09ERAZqclM5s43Gethq3ABjtwRbjRuw3dhtto2GiIiaFiYaIiKSFBMNERFJiomGiIgkxURDRESSYqIhIiJJMdEQEZGkmGiIiEhSTDRERCQpJhoiIpIUEw0REUmKiYaIiCTFRENERJJioiEiIkkx0RARkaSYaIiISFJMNEREJCkmGiIikhQTDRERScrJ0gGQaWWdykPSvosoKimHdwtXxEcGoFewn6XDIqImjInGjmSdysMXaWdQUaUAABSVlOOLtDMAwGRDRBbDqjM7krTvoirJKFVUKZC076KFIiIiYqKxK0Ul5QYtJyIyByYaO+LdwtWg5URE5sBEY0fiIwPg4qT+I3VxckB8ZICFIiIiYmcAu6Js8GevMyKyJkw0dqZXsB8TCxFZFVadERGRpJhoiIhIUkw0REQkKSYaIiKSFBMNERFJiomGiIgkxURDRESSskiiSU1NRUxMDAYMGIAvv/yy3vqcnBzEx8cjOjoas2bNQlVVldr606dPIyQkxFzhEhFRI5g90eTn52PJkiX46quvkJycjG+++QYXLlxQ2yYxMRFz5szBzp07IYTA5s2bVevu37+PefPmobKy0tyhExGREcyeaA4ePIiePXvCy8sLHh4eiI6ORnp6ump9bm4uysrKEBYWBgCIj49XW79o0SKMHz/e3GETEZGRzD4EzY0bNyCXy1WffXx8cOLECa3r5XI58vPzAQCZmZkoKyvDwIEDjT6+t7en1nVyeXOj92tpthq7rcYNMHZLsNW4AduOvbHMnmgUCgVkMpnqsxBC7bO29QUFBVi5ciXWrVvXqOMXFd2FQiHqLZfLm6Og4E6j9m0plo7d2OmjLR13YzB287PVuAHbjd1UydHsVWd+fn4oKChQfS4oKICPj4/W9YWFhfDx8cHevXtRXFyMMWPGIC4uDgAQFxeHu3fvmi94qkc5fbRycjXl9NFZp/IsHBkRWQuzJ5revXsjKysLN2/exP3797Fr1y5ERESo1vv7+8PV1RXZ2dkAgJSUFERERGDEiBHIyMhASkoKUlJSVOs8PbVXhZH0OH00ETXE7InG19cX06dPR0JCAoYOHYrY2FiEhoZi0qRJOHnyJABg8eLFWLhwIQYOHIh79+4hISHB3GGSnjh9NBE1xCLz0QwePBiDBw9WW7Zq1SrV/zt06ICtW7fq3MfZs2cliY0M493CVWNS4fTRRKTEkQGoUTh9NBE1hDNsUqNw+mgiaggTDTUap48mIl1YdUZERJJioiEiIkkx0RARkaSYaIiISFJMNEREJCkmGiIikhQTDRERSYqJhoiIJMVEQ0REkmKiISIiSTHREBGRpJhoiIhIUkw0REQkqQYTzb1795Cbm1tv+fnz5yUJiIiI7IvORHPgwAFERkYiLi4OI0aMQH5+vmrdW2+9JXlwRERk+3Qmmv/85z/YuHEjDh8+jF69eiEhIQElJSUAACGEWQIkIiLbpjPRCCEQFBQER0dHvPHGG3jyyScxbdo0KBQKc8VHREQ2TmeicXBwwMWLF1WfZ86cCSEE3n33XSYbIiLSi85E8/rrr2P06NFIS0sDADg6OuLjjz/GuXPn2BmArF7WqTwkrjiAiYt+QOKKA8g6lWfpkIiaJCddK5988klkZmaisrJStaxFixbYtGmTKvkcOXIE4eHh0kZJZKCsU3n4Iu0MKqpqnryLSsrxRdoZAECvYD9LhkbU5DTYvdnT0xMPPvig+pccHPDMM88AABYuXChNZESNkLTvoirJKFVUKZC076KWbxCRVBr9wiZ7n5E1KiopN2g5EUmn0YlGJpOZIg4ik/Ju4WrQciKSDoegIbsUHxkAFyf1X28XJwfERwZYKCKipktnZwAiW6Vs8E/adxFFJeXwbuGK+MgAdgQgsoBGJxq20ZC16hXsx8RCZAUaXXX2xBNPmCIOIiKyU3o90fz2229YtWoViouL1Z5gPvnkE8ycOVOy4IiIyPbplWhmzJiB0NBQ9OjRg73MiIjIIHolmvv37+Odd96ROhYiIrJDeiWaNm3a4MaNG/Dx8ZE6HiIik8g6lcdeh1ZCr0SjUCgQGxuL4OBguLr++cLbJ598YtRBU1NTsXLlSlRVVWH8+PEYM2aM2vqcnBzMmjULpaWl6N69O+bOnQsnJydkZ2dj4cKFqKyshJeXFxYsWAB/f3+jYiAi+8Wx7qyLXokmKioKUVFRJjlgfn4+lixZgqSkJLi4uGDUqFF44okn0K5dO9U2iYmJmD9/PsLCwvD2229j8+bNGD16NBITE7FixQp06NABW7duxfz587Fy5UqTxEVE9kPXWHdMNOanV/fmYcOGqUZorqqqQteuXTFs2DCjDnjw4EH07NkTXl5e8PDwQHR0NNLT01Xrc3NzUVZWhrCwMABAfHw80tPTUVFRgWnTpqFDhw4AgKCgIFy/ft2oGIjIvnGsO+uiV6LZv38/hg8fjoyMDGRmZuLZZ59FRkaGUQe8ceMG5HK56rOPjw/y8/O1rpfL5cjPz4eLiwvi4uIA1FTlLV++HP379zcqBiKybxzrzrroVXW2dOlSbNy4UVW9df78eSQmJhpV0CsUCrUu0kIItc8Nra+oqMCMGTNQVVWFKVOmGHx8b29Prevk8uYG789a2Grstho3wNgtQd+4n48NxvItv6C8slq1zNXZEc/HBlvs3G31mpuCXommsrJSrQ2lffv2qK6u1vEN7fz8/HDs2DHV54KCArXebH5+figoKFB9LiwsVK0vLS3F3/72N3h5eWHlypVwdnY2+PhFRXehUNQfNkcub46CgjsG788a2Grstho3wNgtwZC4g1t7IWFgUL1eZ8GtvSxy7rZ8zU1Br0Tj5uaGkydPolOnTgCAkydPwt3d3agD9u7dG8uWLcPNmzfh7u6OXbt2Yd68ear1/v7+cHV1RXZ2Nrp164aUlBREREQAqOkk0KZNG8ydOxcODhx4moi041h31kOvRJOYmIipU6eiTZs2AIDff/8dS5cuNeqAvr6+mD59OhISElBZWYlnn30WoaGhmDRpEl577TV06tQJixcvxjvvvIO7d+8iODgYCQkJOH36NDIzM9GuXTtVRwQfHx+sWrXKqDiIiMg8ZELP4ZeLi4vxyy+/QKFQICwsrN70zraCVWfWw1bjBhi7Jdhq3IDtxm6WqrOUlBTExcVh7dq1assvXboEAJgwYYJJgiAiIvulM9H88ccfAIBz586ZJRgiIrI/OhPNa6+9BgBYuHChallFRQUKCwvRsmVLaSMjIiK7oFfXrd27d2PevHm4e/cuBg4ciLi4OHzxxRdSx0ZERHZAr0Tz6aefYuTIkdi1axfCwsKwZ88epKSkSB0bERHZAb0SjRACQUFBOHjwICIiIuDp6Qk9O6sREVETp1eicXBwwI4dO/DTTz+hT58+2LdvH2faJCIiveiVaP7xj39g8+bNmD59OuRyOVauXIlZs2ZJHRsREdkBvUYG6N69O9atW6f6/PXXX0sVDxER2Rm9Es3x48fx4Ycf4vbt22ptM6mpqZIFRtbBWqbDtZY4iMhweiWaOXPmID4+Ho8//jjbZpoQa5kO11riIOvHGxLrpFeicXJy4nAzTZC1TIdrLXGQdeMNifXSqzNA+/btcfbsWaljIStjLdPhWkscZN103ZCQZen1RHPlyhUMHz4cLVu2hKvrn1Ohso3GthhareDdwlVjYW7u6XCtJQ6ybrwhsV56JZrp06dLHQdJzJhqhfjIALXvAICLkwPiIwOkD9gK4yDrxhsS66VX1Vl4eDjc3Nzw22+/ISwsDM7OzggPD5c6NjIhY6oVegX7YfygDqo/VO8Wrhg/qIPZ67utJQ6ybvGRAXBxUi/SeENiHfR6oklKSsKaNWtQXl6OqKgovPTSS5g+fTpGjhwpdXxkIsZWK1jLdLjWEgdZL+XvB3udWR+9Es2GDRvwzTffYOzYsfD29kZSUhJefPFFJhoboq1aAaipVgP4B0q2jzck1kmvROPg4ABPT0/V50ceeQSOjo6SBUWmFx8ZgFWppzWu25RxDvfLqlD9v3dxi0rK8fn2mm2t5Y+W70cQ2S69Eo2XlxdycnJUL2t+9913eOCBByQNjExHWUhrc/d+Vb1l1QL4avfZeoW5JQp8vh9BhuBNifXRK9G8/fbbmDZtGi5fvoy+ffvC1dUVK1askDo2MoG6hbQhSsuqde7LXAW+Lb6wycLOMnhTYp30SjQBAQFISUnBpUuXUF1djbZt28LZ2Vnq2MgENBXSptxXYwr8vdlXsG77KRSVlMPFSYbKagEhAAcZEBnWEuOiOwCwvfcjWNhZji3elDQFeiWasrIyZGZmori4GABw9OhRAMCYMWOki4xMwhSFce2788YcQ9d+Kqr+HKxVIYA9x68BAMZFd2iwI4O1FSAs7CzH1m5Kmgq9Es3UqVNRUlKLhfcQAAAd4klEQVSCVq1aqZbJZDImGhvg6e6ksQ1G3+/qU/VW94U4TdVGAAyuwtv38zWMi+6gsyODtRTeWafy8NXus/WqG2tjYSc9vrRpnfRKNPn5+dixYwdHbrZBjZly++79KqzZfhoKHbuo+0Lchp1nVE8jwJ/VRs5OMoOr8JTH7RXspzXRWEPhnXUqD59vP63qtacNCzvpcRQJ66RXogkMDERhYSHkcrnU8ZCJ6brD1oeuJFO3kTvrVJ5aklGqqFKgwoiHKoda9zXWfKeatO9ig0mGhZ3xlE/IN0vK8VADHSv40qZ10ivRDBw4EIMGDUJgYCCcnP78yvr16yULjCyjmZujXsnJQVbzNKHsNt0r2M/ko+RGhrVU/T80wFtjEgsN8AZg2V5eDT1VmSMee+3lZkzHCr60aX30SjT//e9/MWXKFLRu3VrqeMjC9H0CUtR5uXNTxjmdbUEyGeDs6KBX9VndXmcAcOJikcZtT1wssngvL12dFbxbuOLfL/WR9PiWPn9dcTU2+bFjhX3QK9G4u7tj0qRJUsdCNqpaaH7pszYhAJlM6HxicnFy0DpYpq7eRJYujOIjAzS20TjIZCirqMLERT9I+pRh6fPXxFTJj73I7INeiaZ379748ssvERUVBRcXF9VyLy8vyQIj05DJagp5a1BeKVBRWY1Jgx9Hr2A/nLpcrHqPpqGCWFcbjaULI2XMtXuduTo7oqqqWvVZW0Frirt+S5+/JqZKftbcNkf60yvRrF27FhUVFZg3b55qmUwmQ05OjmSBUeNlncqDi5MM5ZVWkmkACPw5tM1T3R5FcOs/b1ayTuUhccUBjYWurt5E2t7NMWdhVLddoOY81J/c6ha0prjrzzqVBweZ5k4bliyMTZX8TNmLzF7bsWyBXonmxIkTUsdBJvZnIWY9SUZJU9VZQ4VuQ72JrK1Lqz4Frb53/doKSOU105RkLH3+pnoSqf1z16fXmTbW2o7VVOhMNCkpKYiLi8PatWs1rp8wYYIkQVHjmXLoGXPQp9DV1pvIGru06lPQ6kpGiSsOaHzRtXYBqe1n7CCDxSeGM+WTiPLnLpc3R0HBHaPiscZ2rKZEZ6K5fPkyAODcuXNmCYZMx5obS12d608xoW9Vi7a7+4a6tJq72kRTQQv82R0baLh96Yu0M3Bxrt9TT1lAavuuQlj+Lt3akr81tmM1JToTzZEjR5CQkAAhBGQymdpb5o0ZJSA1NRUrV65EVVUVxo8fX28om5ycHMyaNQulpaXo3r075s6dCycnJ1y7dg2JiYkoKipC27ZtsXjxYjRr1szoOOyZrkKsMfvT1h5gCCcNUxlpi7eZm6Oq3cbT3anevDn6VH8YU23S2MTUK9gPP524hpw/itWWHziZh3atvNAr2E9rMlKqedFV8zplXJZum9LFmt5nsfZrZe8cdK0cO3YsxowZA19fX3h4eCAhIQETJkzAQw89hDZt2hh1wPz8fCxZsgRfffUVkpOT8c033+DChQtq2yQmJmLOnDnYuXMnhBDYvHkzAGDu3LkYPXo00tPTERISwqkKdIiPDICTo2mGDGrm5qiaj72xSQbQ3EYTHxkATfcu5ZUKVQFx935VvS7Eyrv72pSdCiYu+gGJKw7gq91ntT4V1N0261SeKjEpj1tUUo5Vqafx6kf7VLORNiTrVF69JKMpXhdnnX+CWimTn4uT+vct3TZjarV/PhPn79L7+tfVFK6VNdP5Wx4dHY3o6GhcvnwZK1asQFRUFJ5++mksWbIEZ86cMeqABw8eRM+ePeHl5QUPDw9ER0cjPT1dtT43NxdlZWUICwsDAMTHxyM9PR2VlZU4evQooqOj1ZaTZr2C/dSGcGmM0rJqrEo9bbI2H+VTypA3U1SF+4WrxRq7YVc1NLYL1Ks/NCUJbe/tKJ9sam/7RdoZbMo4p/FcS8uq8UXaGb0KO12jJBSVlCPrVB7W7shp8P2jZm6OWgvIXsF+GD+og+qu3LuFq8XbZkyp7s+y4NZ9va9/XfZ+raydXr3Obt26hfLycri7uwMASktLcfv2baMOeOPGDbUx03x8fNR6tdVdL5fLkZ+fj1u3bsHT01M1BI5yOWmWdSrPKnucyVDzlFJapl64VzYiidWu/jCkE4SDDBqfdHR9X/lEMuSp9lq3yTqV12C15erU02jop+Pi5IDRUUEAtLd11K2e0tRFvO73n48NVutWbq1M3YBvSFUeu0Kbll6JJjY2FiNHjkRUVBSEEEhPT8fIkSONOqBCoVBr31G2/zS0vu52gHHtRN7enlrXyeXNDd6ftagbe/JPWRaKRDuZDPB0d8ade5VqyxvzpOTq7IjnY4NV539Tz3YpV2dHlFcaN+Co8hiafl/2Zl/B+vSzDe5DV5KRAXj4QXckDOqIp7o9CgA6E1vdYyvPq6ikHGt35EAIoPp/dZ5FJeVYvuUXvDKis2rf+tqbfQXr03JQeOu+Kj4A9ZYZul9ttP0si0rK8Y9Ps0x6rNo0Xcf16WfRorlbo45ny+VLY+mVaKZNm4bg4GAcOnQIADBjxgxERkYadUA/Pz8cO3ZM9bmgoAA+Pj5q6wsKClSfCwsL4ePjg4ceegh37txBdXU1HB0d631PX0VFd6HQ0NDQmK6TlqYp9oJb943al76DahpDCNRLMoZykMng7uqA0rJq1Z1mcGsv1fk/pKXR19PdCa7Ojmp3qLp6buny0P+eoDT9vqzbfsroBKa0ZsbTqv8b8jup6diaqh7LK6uxbvspg55q6naoKLh1H0s2HYeiVn1nwa37WLb5Z5TcKTPJ3b+2n6UUx6pN03Wse80MfeKx1fLFVMlRr0QDAP3790f//v0bfcDevXtj2bJluHnzJtzd3bFr1y61EQf8/f3h6uqK7OxsdOvWDSkpKYiIiICzszO6d++OHTt2YPDgwUhOTkZERESj47EXdX/xjZ3wTKokA/xZxaWp8HB11m8EAwGB0VFBWv+otb2/8df+gRq/o22eG0DzpHENNSA3tqdfMzcNXfL0ZMixDY1TUzWWQkOjminfTdGnV56+L7caoqGu0Hz503DGdXlpBF9fX0yfPh0JCQkYOnQoYmNjERoaikmTJuHkyZMAgMWLF2PhwoUYOHAg7t27h4SEBADAu+++i82bNyMmJgbHjh3D66+/bu7wrdLe7Cv1GrXvlxk3q6aU7t6v0vhH7OLkgISBHdGvS8sGOzAIAazZflqtl1htvYL90KfTnx0hHGRAgH8LJO27WK9nma4Ge+8Wrvh4WgQmDX7coAbkxnaXDe/oa/R3DTm2oXFKmcS0qduA39CxNHUEMabzgLbjKZfrajsizWSiMVMw2iB7rDr7x6dZGqvKZNDdFtCQxlYx6UPTlAAAMHHRD3p9v+6Iz/pMPe0oA2QOMq092nSNIq2k7fdFn+MraRrwVJ9ja6Pp2E6OMgiFUOsW7ursiISB2p8KNXnx/R/07touxdQIyg4Ouo6lzzb60HQda/9cdP1ufl6r2rM2Wy1fTFV1ZvYnGjK9Qi3tMY29gygqKVe9WCgVhah5ibHuXae+d9x17yT16XVWLbR3m25st1d97sKVx/FwrV9N1pg7Y01deCfEdMTEWPWnsldGdDb4/Ax5f0qK35f4yIB6I0rUrcY01dv/DXWFbuiJh+rTu42GrNfDD7prfKJxkNXcNVcb2anLQab/k0VjaKprj48M0KsLMKBekDT2ycsUd+LK81iz/bTWUZX//VIfrde2MefQ0HhwgHF31/qONNGvS0tJ2il6BfuhRXM3ndNKmPLtf11doU05jltTwURjBxIGdcSyzT9raKwFjG9abvxQM4aoPZCk6g9cz7q/2gVJY4beaUxDfG36jqpsbcOi6GpI19Ywr+wmbo53TepOK1GXuRKAtY3jZguYaOzAU90eRcmdMo130Hq8WG81lI23F64W48DJPL0mbKtbkDTUU0mX8koFsk7lNbrA0HdUZVMUjKZ6sVCfaRouXC3GnuPX1L4nhFBNZGcqxp6TOROANY3jZguYaOxEr2A/nV11peLdwhV371c1+t0RpYoqBfb9fE3n05RyYE9NBUndwsaQQUCrqoVJuubqO6qypoIxNMAbSfsuYlXq6QYLSlN2s9XnLfwTF4vqfc/UQ+3rOqchTzXcMM0EYJ2YaOyIqUds1sZBBrwQ++ddrKnbcRpKDJp6qdVWu7DJOpVnUAI2xfUzpEqsbqyGJA5TDtGiT0O6OYba13VO+oyOQNaJvc7siKYRaqWgEFCNfAyYvk2hoXdpNPVS06ZXsB883fW/nzLFuRg7UrCh72eYsuDXpyeVOXpbcd4Y+8REY0c0dcuUSu2X4UzZ2Ori5IDIsJY6E2ZFlQKrUk9rfGFTk7/2D6y3PydHGerOomCqhmNjRwo2tJA1ZcGvT3I0x1D77Dpsn1h1Zmfq1lGbolpr0uDHNb60qbzb/vdLfTQ2FOtLU5tLu1ZeDb4oqm+bhLZGYk3LTFW/b0xbgaG90Ew9XTKg+3qYo7GdXYftExMN6eTiJNPZ0UBZMI6L7qCWHOrOhqmJrjfUlQW1tre9lfRtk9Dn/RJLM7SQNXXBr09ylLqxnV2H7RMTjZ1r7GjM4/83FLw+d9ua5kap26PqxMUig+ZF0ae7sr3U3xtTyNpjLyt7PKemjonGzo2OCsLq7af1eielLk93p0a989FQgaHPG+p1C19N7Kn+noUs2SMmGjunLLQ2ZZwzaNoA5dD6dffTmCqNuk84+s70qCx8tQ12yPp7aghnzLQsjt78P7Y6uipgeOwNdRCQ4g9RU5IwZhRhaykwmtLvi7UwNu6GRmM2B1u+5qbAJ5omSFd7S91BJU1VsGt6R6S8strglwtZtUSGMuWLrWQcvkfTBOn7PsSGnWewKvV0oyeSUn7XkOVEpsLfPctjommC9HmhMOtUnsb3YoydL4Uv4pGl8HfP8lh11gRoq/7SVW2gK5kYcyeoqdeaq7MjG/JJcnwJ1PKYaOycroEaAe29yHQlE2Mnkqp7PH17nRE1Bl8CtTwmGjuh7alFW0Po59/nwNFBpnWkYF0jQRt7J1j3KcpWe+KQ7WEnEstiG40d2Jt9BV+kndHYaK8tWVQrhM6RgrWNBC3VVL1EZL/4RGMH1qflaE0ahs5Ro9yW1Q1EZCpMNHag8NZ9jcuLSsrRr0tLg0ZV1jV2GRGRMVh1ZgceftBd43LvFq4ap9/Vhj1xiEgKTDR2IGFQR60vYOqqNuvYxsvgybmIiAzFqjM78FS3R1Fyp0xrrzNNycbFSYbEv3a1QLREhlH2qLxZUo6H2FZok5ho7IS29hRtL6uNH9TBnOERGUXXe2C2OBhrU8WqsybAxfnPH3MzN0dWkZHN0DUgpr6UycoUY/aRcZho7JjyD6z2PDSlZdXYlHGOf2RkE0wxIKYpkhU1DhONHdP0BwYAd+9X8Y6ObIIpBsTk6M2Wx0Rjx3T9IfGOjmyBvlNa6MLRmy2PicaONfSHxDs6sna1p7SQwbhu+KZIVtQ47HVmxzT1OKuNd3RkC5Q9Ko0dhJXDKVkeE40dU/4hfbX7LErLqtXW8Y6OmhIOp2RZTDR2TvkHxvcIiMhSzJ5orl27hsTERBQVFaFt27ZYvHgxmjVrprZNRUUFZs2ahV9//RVubm5YvHgxAgICUFpairfffhu//fYbAGDq1Kl45plnzH0KNol3dERkKWbvDDB37lyMHj0a6enpCAkJwYoVK+pts2HDBri7uyMtLQ1vv/02Zs6cCQD47LPP0LJlS6SmpmLdunVYuHAhCgsLzX0KRERkALMmmsrKShw9ehTR0dEAgPj4eKSnp9fbbu/evRgyZAgAoEePHrh58yauXbuG8PBwjBs3DgDg7e0NLy8vJhoiIitn1qqzW7duwdPTE05ONYeVy+XIz8+vt92NGzcgl8tVn+VyOfLy8tCnTx/Vsh07dqCiogLt2rUzKAZvb0+t6+Ty5gbty5rYauy2GjfA2C3BVuMGbDv2xpIs0aSlpWHhwoVqy9q0aQOZTKa2rO5nABBCqC0XQsDB4c+Hr7S0NCxYsACrV69WJS19FRXdhUIh6i235fnrbTV2W40bYOyWYKtxA7Ybu6mSo2SJZtCgQRg0aJDassrKSjzxxBOorq6Go6MjCgoK4OPjU++7vr6+uHHjBlq3bg0AKCwsVG23YcMGrFmzBmvWrEFQUJBU4RMRkYmYtY3G2dkZ3bt3x44dOwAAycnJiIiIqLddZGQkUlJSAADHjh2Dq6srWrZsiYyMDKxbtw6bNm1ikiEishEyIUT9eiQJ5ebmYsaMGSgqKsIjjzyCDz/8EA888AA2bdqEGzduYNq0aSgvL8ecOXPw66+/wsXFBfPnz0dwcDCGDBmCmzdvwtvbW7W/+fPno1OnTnofn1Vn1sNW4wYYuyXYatyA7cZuqqozsycaS2OisR62GjfA2C3BVuMGbDd2UyUaDqpJRESSYqIhIiJJMdEQEZGkmGiIiEhSTDRERCQpJhoiIpIUEw0REUmKiYaIiCTFRENERJJioiEiIkkx0RARkaSYaIiISFJMNEREJCkmGiIikhQTDRERSYqJhoiIJMVEQ0REkmKiISIiSTHREBGRpJhoiIhIUkw0REQkKSYaIiKSFBMNERFJiomGiIgkxURDRESScrJ0AObm4CAzap21s9XYbTVugLFbgq3GDdh27I0lE0IISwdBRET2i1VnREQkKSYaIiKSFBMNERFJiomGiIgkxURDRESSYqIhIiJJMdEQEZGkmGiIiEhSTDRERCSpJpVorl27hjFjxmDgwIH429/+htLS0nrbVFRUIDExEYMGDcKwYcNw8eJFAEBpaSmmTZuGwYMHY/Dgwfj+++9tJnals2fP4plnnjFLvKmpqYiJicGAAQPw5Zdf1lufk5OD+Ph4REdHY9asWaiqqgKg33lKzdjYlT766CMsW7bMXOGqGBt3dnY2nn32WcTFxWH8+PHIzc01d+hGx37s2DHEx8dj8ODBmDp1Km7fvm0TcSudPn0aISEh5gpXjbGxb9u2DX379kVcXBzi4uKwZMmShg8mmpDJkyeL7du3CyGEWL58ufjggw/qbbN69Woxe/ZsIYQQR44cESNGjBBCCPHhhx+KRYsWCSGEKCwsFH369BEFBQVmirxxsQshxLZt20Tfvn1Fv379JI81Ly9P9OvXT9y6dUuUlpaKwYMHi/Pnz6tt88wzz4jjx48LIYSYOXOm+PLLL4UQ+p2ntcZeUlIiZs6cKUJDQ8XHH39sM3H369dP5OTkCCGE2LJli5g6darNxN6/f3/Vtv/+97/Ff/7zH5uIWwgh7t27J0aNGiUCAwPNFrNSY2J/7733RGpqqkHHazJPNJWVlTh69Ciio6MBAPHx8UhPT6+33d69ezFkyBAAQI8ePXDz5k1cu3YN4eHhGDduHADA29sbXl5eKCwstInY79y5g8zMTHz44YdmiffgwYPo2bMnvLy84OHhgejoaLV4c3NzUVZWhrCwMLXz0fc8rTF2AMjMzMRf/vIXTJgwwawxNybuiooKTJs2DR06dAAABAUF4fr16zYROwDs2LED7dq1Q2VlJfLz89GiRQubiBsAFi1ahPHjx5st3toaE/vJkyexbds2DB48GH//+9/1eopsMonm1q1b8PT0hJNTzYDVcrkc+fn59ba7ceMG5HK56rNcLkdeXh769OmDli1bAqj55a6oqEC7du1sIvbmzZtj2bJleOSRR8wSb904fHx81OLVFGd+fr7e5yklY2MHgKFDh2Ly5MlwdHQ0X8Ba4tI3bhcXF8TFxQEAFAoFli9fjv79+5svcA2xGXLNnZ2dcfbsWURGRuLw4cNmqxpubNyZmZkoKyvDwIEDzRZvbY2JXS6X46WXXsJ3332HRx55BO+9916Dx7PLaQLS0tKwcOFCtWVt2rSBTKY+THfdzwAghFBbLoSAg8Of+TgtLQ0LFizA6tWrVQWiKUkZu7koFIp6cdT+rG193e0AzecpJWNjt7TGxl1RUYEZM2agqqoKU6ZMMU/QesbW0PqgoCAcPHgQX3/9NaZPn46vv/7aquMuKCjAypUrsW7dOrPEqUljrvl///tf1fIXX3wRUVFRDR7PLp9oBg0ahB9//FHt3+eff447d+6guroaAFBQUAAfH5963/X19cWNGzdUnwsLC1XbbdiwAe+//z7WrFmjqmqwldjNyc/PDwUFBarPdeOtu14Z50MPPaTXeUrJ2NgtrTFxl5aW4sUXX0RVVRVWrlwJZ2dn8wWuITZ9Yy8vL0dGRoZq+ZAhQ3D27FnzBK0hLn3j3rt3L4qLizFmzBjV02RcXBzu3r1r9bHfuXNHLUEKIfR6grfLRKOJs7Mzunfvjh07dgAAkpOTERERUW+7yMhIpKSkAKjp0eLq6oqWLVsiIyMD69atw6ZNmxAUFGRTsZtb7969kZWVhZs3b+L+/fvYtWuXWrz+/v5wdXVFdnY2ACAlJQURERF6n6c1xm5pjYk7MTERbdq0wUcffQQXFxebid3JyQlz587Fr7/+CqCmNqBr165WH/eIESOQkZGBlJQU1d9rSkoKPD09rT52Dw8PrF69Gr/88gsAYOPGjXo90TSpXmdXr14VY8eOFYMGDRITJ04UxcXFQgghvvrqK/HRRx8JIYQoKysTb731loiJiRFDhw4Vv/76qxBCiMGDB4s+ffqIIUOGqP6dOHHCJmJXunLlill6nQkhxHfffSeeeeYZMWDAAPHZZ58JIYR48cUXVdcsJydHDB8+XERHR4s33nhDlJeX6zxPczI2dqWPP/7Y7L3OjI371KlTIjAwUMTExKh+r1988UWbiF0IIY4ePSqGDRsmhgwZIiZNmiSuX79uE3HXZoleZ0I07poPHTpUDBw4UEydOlWUlJQ0eCzOsElERJJqMlVnRERkGUw0REQkKSYaIiKSFBMNERFJiomGiIgkxURDRESSYqIhqiU9PR3jxo3D0qVLkZycbLE4zHH8pUuX6jVOFVFj2eVYZ0SNNW3aNLs9fl5eHhYsWIAff/wR8fHxkh2HSImJhpq8pUuXIjU1FV5eXmjTpg0AYMaMGWjfvj1eeOEFdOrUCRMmTMDBgwdx7949vPLKK0hPT8e5c+fg4+ODTz75BB4eHrh48SL+9a9/obi4GNXV1Rg3bhyeffZZHD58GEuWLMGjjz6K8+fPo6qqCnPnzkW3bt1w7NgxLFq0CAqFAgAwZcoUREdHqx3/2LFj+OCDD3D//n04Ozvj9ddfR0REBJKSkrB79244ODjgjz/+gJubG95//30EBAToPN+tW7ciPDwcAQEBZp8ojJomVp1Rk5aRkYFdu3YhOTkZX3/9tcaBDSsqKvDwww9j69atGDp0KN555x3MmjULO3bswN27d5GZmYmqqiq89tprePPNN5GUlISNGzfi888/x88//wwAOHHiBCZOnIjk5GTEx8erZiVctmwZJkyYgKSkJCxYsACHDh1SO/atW7fw2muvYdasWUhNTcX777+PxMREXLlyBQBw9OhRzJ49G9u3b0fnzp3x2WefNXjOr7zyCsaOHWuRkb2paeJvGjVpWVlZiIqKUs2DM3z4cI3bKSdja926NQIDA+Hr6wsHBwe0atUKt2/fxqVLl3D58mW8/fbbiIuLw9ixY1FWVobTp08DAFq2bImOHTsCAB5//HHVk8SgQYPw3nvv4c0338SpU6fwxhtvqB33xIkTaN26NTp37gwAaN++Pbp27YojR44AAIKDg+Hn51dvv0TWhFVn1OTVHu5P25DntYfO1zSMfnV1NZo3b64ajReoGVq9efPm+Pnnn+Hm5qZarpx7BwBGjRqFfv364cCBA9i/fz+WL1+uNtNhdXV1vfluhBCoqqqCs7Oz1v0SWRM+0VCTFhERgfT0dJSUlEChUKglCkO0bdsWbm5uqu9fv34dsbGxqiHstRk1ahRycnIQHx+PefPmoaSkRG0ekLCwMPz22284ceIEAOD8+fM4evQowsPDjYqTyBL4RENNWmRkJM6ePYvhw4ejRYsW6NChA27dumXwflxcXLBixQr861//wurVq1FVVYVp06ahW7duOHz4sNbv/f3vf8eCBQvw0UcfQSaT4ZVXXkGrVq1U6x966CEsXboU8+bNQ1lZGWQyGRYuXIi2bdvi+PHjRp0zkblxmgAiIpIUn2iI7MyCBQu0PkXNnDkTPXv2NHNE1NTxiYaIiCTFzgBERCQpJhoiIpIUEw0REUmKiYaIiCTFRENERJL6/7FYySLUvhY0AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# fig = plt.figure()\n", "# ax = fig.add_subplot(111, projection='3d')\n", "sns.set()\n", "_ = plt.scatter(Trapnell_pca[:,0], Trapnell_pca[:,1])\n", "_ = plt.xlabel(\"dimension_1\")\n", "_ = plt.ylabel(\"dimension_2\")\n", "_ = plt.title(\"scatter map after PCA\", size =20)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### 2. Determine which number of clusters is the best" ] }, { "cell_type": "code", "execution_count": 56, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAEeCAYAAAC+OaPqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3XtclGX+//EXw2E4eUAcUFAQNXBL1ApFMyFPIaJ4yL6Zhw7b1qpfS7dkU/uluWl+c00r85Cb5YbGZtsKWkCmRmmmQSeFVi0FFREGGFDOAzP37w9kaIQQFGYAP8/Hg4fOxTX3fd3X3Mx77uu+5r5tFEVREEIIIVqYytoNEEIIcWuQwBFCCGEREjhCCCEsQgJHCCGERUjgCCGEsAgJHCGEEBYhgWMlBw4c4M9//jPDhg2jf//+3HvvvcydO5cDBw5Yu2ltwoYNGwgICGD//v2mMoPBwI4dOygtLW2wXnsyatQogoKCrN2MFpORkcGjjz7KnXfeyV133cWWLVtuannHjh0jICCAVatWNVMLf196ejoJCQktvp62RALHCl5++WXmzZvHL7/8wujRo3n88ce55557+P7775k3bx4vvviitZvY6g0ZMoT58+fj5+dnKnvuued4+eWXqaqqsmLLRHN6/vnnOXr0KMOHD2fWrFltJlxPnjzJxIkT+f77763dlFbFztoNuNUcO3aMHTt2EBYWxrp167Czq30JioqKeOSRR9i1axehoaGMGTPGii1t3YKDgwkODjYry8/Pt1JrREtJS0ujW7duvPXWW9ZuSpNcvnyZyspKazej1ZEjHAtLSkoCYObMmWZhA9ChQweee+45AD7//HNLN02IVqeyspLOnTtbuxmimUjgWFjNp57Tp0/X+/ugoCBef/11HnvsMbNyg8HAe++9R2RkJIMGDSI0NJSoqCguXLhgVk+v17NlyxbGjx9P//79CQ4OZu7cuZw4ccKs3n/+8x8CAgJISEjgiSeeIDAwkJEjR5qWV1xczNq1axkzZgz9+/dnxIgRLF++vN6jiOjoaKZOnWoaZ58xY8Z1x65Pnz5NQEAAixcvNis/deoUAQEB3HfffWblRqOR4OBgZs2aBdQ9NxMQEMC3334LwODBg5k9e7bZ8ysqKnj99dcZNWoUgYGBhIeH88EHHzTYxhqzZ89m1KhRZGdn89xzzxEcHMzAgQOZOXMmx44dq1M3ICCAK1eumJVnZmYSEBDAvHnzTGU125CRkcGaNWu49957GThwINOnT+fEiRMYjUb+8Y9/MGrUKAYNGsS0adPqrO+3/fnoo48ycOBAhg8fzrJly+p9rRr7ui5evJiAgACOHz/O+PHjCQwMZPr06SiKQl5eHkuXLmXs2LEEBgZy7733EhUVxblz5xrVn43ZR2v6BqqHpwICAhg1atR1l71//35mz55NUFAQwcHBPPbYYyQnJzf4nKa8ZnD9/X3Dhg088sgjALz//vsEBASYvW5paWnMmzeP4OBgBgwYwKRJk4iJieHaq4zV/H1s2bKFoKAggoKC2L59OwCHDx/m0UcfZdiwYQwcOJCJEyfy9ttvo9frr9tH1iRDahY2fPhwoqOjefXVV8nIyGDChAkMGDAAW1tbABwdHQkPDzd7jqIo/PnPf+bQoUP07duXadOmUVBQQHx8PEePHuXf//43np6eVFRU8Pjjj/Pdd9/h7+/Pww8/TF5eHvv37+fQoUO8/vrrdYbpVq5ciYeHB7NnzyYzM5OePXtSVFTEjBkzOH36NMOGDeP+++8nMzOTXbt2cejQIf71r3/h4eEBwNatW3nttde44447mD59OpWVlSQmJrJw4UIqKiqYPHlyvf3g7++Pl5cX33zzjVn50aNHAbh06RIXLlygZ8+eABw/fpzCwkJCQ0PrXd78+fPZvXs3Fy9e5Mknn6R3795mv1+1ahVGo5Fx48ahUqn45JNPWLFiBVVVVaY3h4aUlJQwY8YMnJycmDx5Mnl5ecTHx/PEE08QHx+Pj4/PdZfxexYuXMjly5eJiIjg0qVLfPbZZ/zpT39i1KhRJCUlERYWRkVFBXv27GHOnDkkJibi6elpen55eTmzZs2iR48ezJw5kxMnTvDhhx9y7NgxPv74Y1xdXQGa9LrWmDt3LoGBgQwfPhxnZ2f0ej1PPvkkp0+fZuzYsYwbN47z58/z6aefcvjwYRISEho8ImnsPlpzju6tt96ia9euTJ8+nQ4dOjTYj2+//Tbr1q3D3d2dsLAw1Go1n3zyCY899hhvv/0299577w2/RjUas78PGTKEKVOmsHv3bgYOHMiIESPw9vYG4Msvv2T+/PnY29tz//3306VLFw4dOsRLL73Ezz//zMsvv2y2vkOHDvH5558zZcoU8vLyGDhwICkpKcyZMwc3NzfGjx+PWq3myJEjrFu3jnPnzvHKK6/c9Ha2GEVY3PLlyxV/f3/Tz1133aU8+eSTynvvvadcunSpTv2PPvpI8ff3V5555hmloqLCVL53717F399fefnllxVFUZS33npL8ff3VxYvXqxUVlaa6p04cUIZMGCAEhQUpBQVFSmKoigff/yx4u/vr4SEhCilpaVm63vppZcUf39/ZceOHWbl+/fvN7WjxpAhQ5QxY8aYre/SpUtK//79lalTpzbYDy+++KLi7++vnD171lQ2Z84cZdCgQYq/v7/yn//8x1T+5ptvKv7+/srp06fNHn/++eemOrNmzVL8/f2Vy5cv13leaGiokpubaypPS0tTAgIClIkTJzbYxt8ud+7cuYperzeVb968WfH391def/31BtugKIpy4cIF0zKubdvIkSPN6j/77LOm/SI7O9tUvmHDBsXf31/ZuXOnqWzkyJGKv7+/8vTTTysGg8FUvnbtWsXf319Zt26dqawpr+vzzz+v+Pv7K/Pnzzere/DgQcXf31954403zMrfeeedepd9rabso4qiKP7+/kpkZGSDy1QURTl79qxy++23K+PGjVO0Wq2pPCMjQxk0aJAyYcIERVEU5ejRo4q/v7+ycuVKU52mvGaN3d/rW09paakydOhQZejQocqFCxdM5QaDQXn66acVf39/JSkpyWzb/f39lQMHDpi1q6bu+fPnTWV6vV6ZNGmS8oc//EG5cuXKdfvLWmRIzQpeeukl3n77bUaMGIG9vT3FxcV8+eWXrF69mjFjxvDaa69hNBpN9T/99FMAli5dioODg6k8IiKCOXPmcNdddwGwe/dunJyceOGFF8zOD/Xv358ZM2Zw5coV9u3bZ9aW0NBQnJycTI+rqqqIjY3ltttuY+bMmWZ1R48ezV133cXnn39OcXExUH30pdPpSE9PN9Xr1q0bCQkJ1x2yqjlaqTnKMRqNpKSkMHnyZOzt7c2GQg4fPoy3tze33XZbg8v8PQ8++CBdu3Y1Pb799tvx9PSsMyTZkD/+8Y/Y29vXaX9GRsYNtanG1KlT6dixo+lxzesZERFhdiQzYMAAAC5evGj2fBsbG6KiolCpav+cn376aVxdXdm7dy/Q9Ne1RlhYmNnjmv3y559/pry83FQ+Y8YMkpKSmDFjRoPbeiP7aGMkJiZSVVXFvHnz0Gg0pnJfX1+ef/55HnjggWY5iX8z+/vBgwfR6XT86U9/okePHqZylUplOnf78ccfmz3H0dGxzlF9zWvw3Xffmcrs7e35xz/+wbFjx657JGhNMqRmJffddx/33XcfJSUlpKSk8M0333Dw4EHOnTvH1q1bMRqNREVFAdVj2F5eXmZvPlD9RvOXv/wFqB6bv3DhAnfddZdpCOW37r77bt59911OnjxpVl5zqF8jPT2d0tJSDAYDGzZsqLOciooKDAYDp06d4u677+ahhx5i69atREZGEhgYSEhICKGhoQQGBl63D4YNG4aDgwNHjx5lxowZpKamcuXKFe69915SU1NJSUkBqmf8nDhxgunTp193mb+nV69edco6d+5Mdnb2DS+jpp9vdtz82uG4mg8Av31TAlCr1fWuz8PDwzT0WMPBwYF+/fqRkpJCUVER2dnZTXpda1y7f9xzzz307NmTL774guHDh3PPPfcQEhLCfffdR/fu3RvczhvdRxuj5jmDBg2q87ub2W+udTP7e2pqKlB9Dqe+18DW1rbOtnfr1s003F7jwQcfZP/+/Tz//PNs3ryZESNGEBISwtChQ80+kLZGEjhW5uLiQmhoKKGhoTz//PP8+9//5sUXX2THjh3Mnz8fJycnrly5YvbpvD4lJSUAv/vppmZs/refSqH2TaxGzYnTs2fPNjgV9fLlywA8++yz+Pr68q9//Yvjx4/z008/sWHDBvz8/Fi+fDnDhg373WU4OzszePBgjh07htFo5OjRo6hUKgYPHsx3333Htm3b0Gq1pKSkYDAYfvf8TWNcu5034to/ZhsbG4A6J3ubytnZuVHr+z2/t2+4uLgAUFpa2uTXtYajo6PZYycnJ3bt2sXmzZtJSEhg37597Nu3D5VKxdixY/nb3/72u+dwbnQfbYya7asvyJrTzezvRUVFQO2IRX2u1/9QfWT9/vvvs23bNo4cOUJ0dDTR0dF07tyZ+fPn15kw05pI4FhQcXExU6dOxc/Pj7fffrvO721sbHjwwQdJTEzk8OHDZGdn4+fnh7Ozs+mP9VqlpaU4Ozub3ly0Wm299Wr+IK83xbRmOZMmTWLNmjXX3SYbGxumTZvGtGnTyM/P58iRI3z++efs27ePuXPncvDgQbp06fK7zw8JCeHrr7/mv//9L99++y39+vWjY8eOBAcHs23bNlJSUjh06BBOTk4MHTr0uu1pLa4NoRt5E22sa2dX1dBqtdjY2NCxY8cmv64N6dKlCy+88AJLly7l1KlTHDp0iLi4OD777DNUKhWvv/56vc9rrn20PjWhXVJSgpubm9nvysvLcXBwMBtyrE9jXrOb2d9r2rh9+/YGg6kxhgwZwpAhQygtLSUlJYWkpCR2797NypUr8fHxuakPZy1JzuFYkKurK0VFRRw5coS8vLwG66pUKtNYtL+/P1lZWeTm5tapN3nyZMLCwnB1daVHjx6kp6ej0+nq1Ks5H9K3b98G1+vn54eDgwNpaWn1fnLfvn07mzZtoqCggIKCAjZs2MDu3bsBcHd3Z+LEibz55ptMnTqVsrIyfv755wbXFxISAsCRI0f48ccfGTx4MFA9vGJnZ8e3337L119/zdChQ5vlKKWl1RyV/PbyOgDnz59vsXVmZWXVmdZcWFhIRkYGvXv3xsnJqUmva0OSk5NZuXIl58+fx8bGhn79+vHkk0/y0Ucf4ezsbBoGrU9z7aP18ff3B6pnM15r5cqVDBw48HfP1zX2NWvK/l5z9PtbNdO8a4bWfquwsJBVq1YRFxfX4HYC/POf/zSFurOzMyEhISxbtozly5cD5ud2WhsJHAubOXMmer2eZ555pt5PegcOHODIkSOMHTvWNDwQGRmJoiisXbsWg8FgqpuQkMC5c+dMn5amTJlCeXk5r7zyitnlXdLS0tixYwcdO3a87ncZ1Go148eP59dff+W9994z+92xY8dYs2YNH3/8MZ06dcLFxYX333+f9evXU1hYaFY3KysLAC8vrwbX17t3b3x8fIiJiaGoqIghQ4YA1W9Ot99+O5988gk5OTl1vpdTn5oT+tb8hnfNpXa++OILU1lFRQXbtm1rsXUaDAY2btxoeqwoCuvWraOsrIwHHngAaNrr2pDc3Fyio6N59913zcrz8vKoqKioc87nWs2xj9ZnwoQJqFQqtmzZYhaa58+fJyEhgZ49e9Y5z1Wjsa9ZU/b3mgkRv90Xa/6m33nnHbNJBwB///vfef/99xv1weTw4cNs2bKFH3/80ay8ZjLJ9f7mrEmG1Cxs7ty5nD59ms8++4z777+fe++9l169elFVVcVPP/3E999/T+/evXnppZdMz5k2bRr79u0jNjaWU6dOERwcTE5ODvv27cPb29s0ceDJJ5/k8OHD7N27l1OnTjF06FDy8/PZv38/iqKwfv36Ro1xP//88/zwww+8+uqrHDhwgAEDBpjWZ2dnxyuvvIJKpcLBwYFnnnmGlStXMmHCBMaOHYujoyPJycmcOHGCSZMm1fk+TH1CQkLYsWOH6fxNjSFDhpg+sTZmiKBmUsXSpUsZPnx4o75f09ymTZvGBx98wCuvvMJPP/2Em5sbBw4coEOHDr97ruZmubu7Exsbyy+//EL//v35/vvv+fHHHwkKCjLrg8a+rg0ZM2YMd955JzExMZw+fZpBgwZRXFzMZ599BlTPjmtIc+2j1+rTpw/z58/nzTffZNKkSYwcORJFUYiPj6eiooLVq1f/7nMb+5o1ZX+v2RcTEhJwdnZmypQp3HbbbaxcuZJFixYxZcoUxowZg4eHB99++y0nTpwgMDCQP/7xj9fd1qeffppjx47xyCOPMG7cODw9Pfn111/54osv6NOnD5GRkU3uP0uRIxwLs7W15c033+Stt95ixIgRnDhxgvfff5+PPvqIiooKnnvuOXbv3m02Dmxra8vmzZtZuHAh5eXl7Ny5k6NHjzJx4kRiYmJMn0rVajXbt2/nmWeeobKykpiYGI4ePcrIkSP58MMPG31tti5durBr1y7++Mc/kpOTQ3R0NCkpKYwaNYpdu3aZXcNs9uzZrF+/nh49ehAfH8/OnTvR6/UsWbKk0V9AqxlWCwgIMPuEXbOegICA686AApgzZw4DBw7k66+/ZufOnY1ad3Pr168fW7dupX///iQkJLBnzx6GDRvG9u3b68w2ai4ajYb33nuPyspKduzYwaVLl3jyySd55513zKZxN+V1/T0ODg68/fbbPPnkk+h0Onbu3EliYiIDBw4kOjr6ul+ubK59tD7/+7//y/r16+nevTtxcXHs3buXAQMGsGPHDgYOHPi7z2vKa9bY/d3b25uFCxdiY2PDzp07TR+cwsPD2bFjB0OHDuXQoUPs2LGDkpIS5s2bx/bt203nuRpSs03Dhw/n6NGjvPfee5w6dYpHHnmEnTt3ttgHm+Zgo9zsFBshhBCiEeQIRwghhEVI4AghhLAICRwhhBAWIYEjhBDCIiRwhBBCWIQEjhBCCIuQL35eVVBQgtHYdmeIu7u7kp9ffP2KtwDpC3PSH+akP2rdTF+oVDa4uV3/e0O/JYFzldGotOnAAdp8+5uT9IU56Q9z0h+1LNkXMqQmhBDCIiRwhBBCWIQEjhBCCIuQwBFCCGEREjhCCCEsQgJHCCFuUZa+WYBMixZCiDbGqChU6A2UVVRV/+gNlFdUUVpRRbneQGl5FeX6q48rrtbTX637m8cujva88tRQ1PYtc6+ma0ngCCGEBRkVheKySkrLq2oD45pQKK8wXA2P+kOjvMLA9Y5NbABHtS1OajucHOxwUtvh4mSPprMTjg52OKlt6ePTBXtbyw10SeAIIUQzMRoVLpfo0RWVU3ClgoKi6h9dUbnp/wVFFRiu82VLBzsVTmo7HNV2OKttcXSwo5OLc50Acap5XPPjUPtY7WCLysamwfVoNB3IzS1qzi5okASOEEI0QpXBSGFRBQXFV0PEFCjlV0OlgsvFeozXnBexs1XRpYMatw5q+vbohFsHNW6ualyc7K8Gh3loODrYYmfBow5LksARQtzy9JWG6iC5UveIRHf13ysl+jrPU9vb0qVjdZjc3ssNtw6OpnCp+XF1ssfmOkcatwqLB05WVhZRUVHk5+fj5+fH2rVrcXExvwCcXq/nhRdeIDU1FUdHR9auXUufPn1QFIU1a9bwxRdfoFKpePnll7n77rsBGD16NK6urqZlbNmyhe7du1t024QQrYvRqFBUqqewWM/lEj2XSyqoNEJm9hVTkBQUVVBcVlnnuc5qO9yuhomvpytuHRxx66D+TaA44qS2lTBpAosHzooVK5gxYwYRERFs3LiRTZs2ERUVZVYnOjoaJycnEhISSE5OZsmSJezatYvPPvuMM2fOEB8fz7lz5/jzn/9MfHw8RUVF2NvbExcXZ+nNEUJYmKIolOsNFBZXH3VcLtFz+TeBUvt/PUWleuqb+dvB2R63DmrcOzrS17uT2RFJzY+jgwwANTeL9mhlZSXJycls3LgRgKlTpzJr1qw6gZOUlMSCBQsAGDx4MDqdjqysLL788kvGjx+PSqXCz8+P7t2788MPP1BWVoaiKEyfPp2KigqeeuopwsPDLblpQoibVGUwXhMgFfWGyZUSPfoqY53n26ps6OTqQCcXB9w7OuLXvSOdXBzo7OpARxc1nVwd6OjiwG293LlcWGqFLRQWDZyCggJcXV2xs6terUajIScnp049rVaLRqMxPdZoNGRnZ6PVavHw8KhT7uTkxIgRI1i0aBF5eXnMnDkTf39/+vTp0/IbJYRoUJXBSN7lcnRXyusNk5qQqW9YC8DVyZ5OLtVh0bdHJzq7qOno4mAKl04uDnRyVePiaNeo4S0HC33nRNTVYoGTkJDA6tWrzcp8fX3r7BD17SCKopiVK4qCSqXCaDTWWz5mzBjGjBkDQI8ePRg7diyHDx9uUuC4u7tev1Irp9F0sHYTWg3pC3Mt3R+KoqC7Uk5WbgmZucVk5RZz8eq/2fmldaYBO9ipcOtYfU7Ep3vH6mGsq4/dOjjSuYOaLh0d6eSqxt6u+Wdsyf5Ry5J90WKBEx4eXmdYq7KykuDgYAwGA7a2tuTm5podsdTw9PREq9Xi4+MDQF5eHh4eHnTr1g2tVmuqV1P+xRdf0LVrVwIDA02/qzmKaqz8/OI2fVMmS8+nb82kL8w1Z3+UlleRU1BKtq6UHF31v9X/L6Oi0mCq52CnwsPNme5dnBnUtyvdujjTtVN1gHRyccDR4fon25XKKgoLqpql3b8l+0etm+kLlcqmyR/ULTqkZm9vT1BQEPHx8UycOJHY2FhCQkLq1AsNDSUuLo6goCBSUlJQq9V4eXkREhLCxx9/zIQJE8jMzCQjI4PAwED+/e9/8+GHH7Jp0yZ0Oh0HDx4kOjrakpsmRLtRWWUkt7DsmkApJbugzGxqsI0NdO3kiGcXZ/x7dqZbF2c8uzjTzc0Zt47q637pUNx6LD4NY/ny5SxevJjNmzfTvXt31q1bB0BMTAxarZYFCxYwe/Zsli1bRkREBA4ODqxZswaAcePGcfz4cSIjIwFYtWoVjo6OTJ8+nVOnTjFhwgSMRiOLFi3C29vb0psmRJthVBQKrlSQXXDtkUopeZfLzWZ2dXS2x7OLMwP7uJtCxbOLMx6dnVpkuEu0XzaKpS8X2krJkFr7IX1Rq7S8kjIjnDyT95thsDK0BaVmM73U9rZ4dnGqDhQ3Z7p1caabuzOebk44O9pbcQuan+wftdr1kJoQomVU6A1k5V89YZ9XwsXcEi7mlVBQVGGqo7KxQdO5egjs9l5u1cNfV386uzrIFxhFi5PAEaINqawycim/OkxqgiUzt5i8y+WmOna2Krzcnenn0xmvri78oXdXnOxs0HR2arfX6BJtgwSOEK1QlcGItqCMi3klXMwtvvpvCdqCMtPFIW1VNnh2ccave0fuHdAd764ueGtc0XR2xFZVGywyhCRaCwkcIazIqCjkFZaZhsBqAiZbV0qVoTpYbACNmxPeXV0I6udxNVhc6NbFWY5YRJsigSOEBSiKQkFRBZm5NUNhxWTmlXApvwR9Ze3Je/eOjnhrXAjs7Y63xgXvrq50c3e22B0ZhWhJEjhCtICCogpSTmm5WBMwecWUVdR+MbKTqwPeXV0IHeh9NVhc8OrqgpNa/iRF+yV7txDNKP9yOfHHznHopyyqDAoujnZ4a1wZekc3elwNFW+NK65O7WuqsRCNIYEjRDPIKyzj06PnOHz8EgDDA7sTHuyDh5uTTDcW4ioJHCFugraglE++Occ3qdnY2EDIQC/Ch/rQtZOTtZsmRKsjgSPEDcjWlfLJkQyOpuVga2vDfXd6Ex7sQ5eOjtZumhCtlgSOEE1wMa+ET49kcOy/OdjbqhgT1INxwT50dlVbu2lCtHoSOEI0Qqa2mL1HMkg5qcXB3pZxQ3wIG+JDRxcHazdNiDZDAkeIBpzPKWLv1xl8dzoXRwdbxg/z5f7BPengLEEjRFNJ4AhRj/RLV9j7dQY//pqHk9qOyOG9GBPUU6YzC3ETJHCE+I0zFy+z90gGx8/k4+Jox+QRfoy5u0e7u0S/ENYggSME8EtmIXu+ziAtXYerkz0PhPZm1F095Jv/QjQji1/5Lysri5kzZzJu3Djmzp1LSUlJnTp6vZ6oqCjCw8OZMmUKZ86cMfv9qVOniIiIMCt79913GTduHGFhYezbt69Ft0G0H6fOF/D3mB9YveN7zucU8eDIPqyZO4yIYb0kbIRoZhb/i1qxYgUzZswgIiKCjRs3smnTJqKioszqREdH4+TkREJCAsnJySxZsoRdu3YBEBsby2uvvYa9fe0Qx/Hjx9mzZw9xcXEUFxfz0EMPMWTIEDp37mzRbRNtg6Io/PdcAXu+zuD0hUI6uTgwfVRfQgd5o3aQi2QK0VIseoRTWVlJcnIyYWFhAEydOpXExMQ69ZKSkoiMjARg8ODB6HQ6srKyKCoq4sCBA6xbt86s/ldffcXYsWNRq9W4u7szZMgQkpKSWnx7RNuiKAqpZ/NZveN71v7rR7QFpcwYcxuvzhnG/UN8JGyEaGEWPcIpKCjA1dUVO7vq1Wo0GnJycurU02q1aDQa02ONRkN2djZeXl5s2LCBzMzMOvUDAwPr1BcCqoPm+Jl89nydQfqlK3TpqGb2/f7cO6A79nYSMkJYSosFTkJCAqtXrzYr8/X1rXMhw/oubKgoilm5oiioVL9/MGY0GuuUNVS/Pu7urk2q3xppNB2s3YRWQ6PpgNGocCwtmw/3n+JM5mU8ujgz/8GBjArywd7u1rpxmewb5qQ/almyL1oscMLDwwkPDzcrq6ysJDg4GIPBgK2tLbm5uXh4eNR5rqenJ1qtFh8fHwDy8vLqrVejW7du5Obmmh7n5ubi5+fXpPbm5xdjNCpNek5rIrcRruXu7spnX59lz9cZZOYW49HZicfH92PYHd2ws1VRWFB3okp7JvuGOemPWjfTFyqVTZM/qFv0Y569vT1BQUHEx8cD1RMAQkJC6tQLDQ0lLi4OgJSUFNRqNV5eXr+73JBZOsIdAAAgAElEQVSQEPbt20dZWRk6nY6jR48ybNiwltkI0WoZFYWUk1qefu0LNsWmUmkw8qcJf2DVU8GMGOAlt2MWwsosPktt+fLlLF68mM2bN9O9e3fTBICYmBi0Wi0LFixg9uzZLFu2jIiICBwcHFizZk2DyxwwYACRkZFMmzaNqqoqnnnmGTw9PS2xOaIVUBSFn87kE/vVWc5ri+nh4cpTkbczpJ8nKpXci0aI1sJGUZS2O47UjGRIre1RFIWfzxWw+6uznM26gqazI5Pv7U1EaF90+cXWbl6rcSvuGw2R/qhl6SE1+WabaJNOXyhk91dnOXWhkC4d1Tw6LoDhgd2xs1VhK0c1QrRKEjiiTUm/dIXdX50lNV1HJxcHZo71J2Sg1y0360yItkgCR7QJF7TFxB46yw+/5OHqZM//jOzLyLu8UdvL92iEaCskcESrdim/hLjD6Xz7Xy1OajumjPBjTFBPuc6ZEG2Q/NWKVklbWMbew+kcScvGwc6WCff4EjbEBxe5TYAQbZYEjmhVdFfK+eRIBoeOX0KlsuH+wT0JH+pLR7nDphBtngSOaBUuF1fw6dFzJP2QhaIohA7yImJYL9w6qK3dNCFEM5HAEVZVXFZJwtFzHPg+k6oqheGB3Zg4vBddOzlZu2lCiGYmgSOsorS8in3J59mXfIEKvYHgOzyZNNwPzy7O1m6aEKKFSOAIiyrXV3Hgu0wSj52npLyKuwM0TL7XD29N279atxCiYRI4wiL0lQaSfrjIp0fPUVRaycA+7kwe0RvfbnKZeCFuFRI4okVVGYwc+imLvUcyKCzWc3svN6aM6E0f707WbpoQwsIkcESLMBiNHEnNZu/XGeRdLqdvj048NfEO+vm6WbtpQggrkcARzcqoKHz73xziDqWTU1BGr24deCQsgDv8utR7d1chxK1DAkc0C0VR+OGXPHYfOsvF3BJ6aFx4emogg27rKkEjhAAkcEQzyCssI3rfaU6czadbF2fmTLqDoH4eqCRohBC/IYEjbpjBaOTz5ExiD5/FBhumj76N0Xd7Y6uSWwUIIeqyeOBkZWURFRVFfn4+fn5+rF27FhcXF7M6er2eF154gdTUVBwdHVm7di19+vQx/f7UqVM8++yzfPrpp6ay0aNH4+pa+12OLVu20L1795bfoFtU+qUr/DPxJOdzihnUtyszx/rj3snR2s0SQrRiFg+cFStWMGPGDCIiIti4cSObNm0iKirKrE50dDROTk4kJCSQnJzMkiVL2LVrFwCxsbG89tpr2NvXXjW4oKAAe3t74uLiLLott6Kyiip2HzrLge8y6ejiwLzJ/bk7QCPnaYQQ12XRsY/KykqSk5MJCwsDYOrUqSQmJtapl5SURGRkJACDBw9Gp9ORlZVFUVERBw4cYN26dWb1T5w4gaIoTJ8+nSlTppCQkNDyG3ML+vGXPF7cdowDKZncN8ibVX8aSlA/DwkbIUSjWPQIp6CgAFdXV+zsqler0WjIycmpU0+r1aLRaEyPNRoN2dnZeHl5sWHDBjIzM83q6/V6RowYwaJFi8jLy2PmzJn4+/ubDcNdj7t727+0ikbTMt/az79cxtbYExw5fgnfbh1Y8ugQ+vXq0iLrai4t1RdtlfSHOemPWpbsixYLnISEBFavXm1W5uvrW+fTcH2fjhVFMStXFAVVAyeix4wZw5gxYwDo0aMHY8eO5fDhw00KnPz8YoxGpdH1WxuNpgO5uUXNukyjovDlDxf595dnqDIoPBDam7AhPtjZqpp9Xc2pJfqiLZP+MCf9Uetm+kKlsmnyB/UWC5zw8HDCw8PNyiorKwkODsZgMGBra0tubi4eHh51nuvp6YlWq8XHxweAvLy8euvV+OKLL+jatSuBgYGmspqjKHFjMnOL+WfiSc5cvMIffN14ZFwAnm5yJWchxI2z6Dkce3t7goKCiI+PB6onAISEhNSpFxoaapoAkJKSglqtxsvL63eXe/HiRTZu3IjRaCQvL4+DBw9y3333tcg2tHf6SgMff3mGFe8lk6Mr408T/sCi6YMkbIQQN83ihwHLly9n8eLFbN68me7du5smAMTExKDValmwYAGzZ89m2bJlRERE4ODgwJo1axpc5vTp0zl16hQTJkzAaDSyaNEivL29LbE57Upaho7oxFNoC8sYHtiN/xnZlw5ya2chRDOxURSl7Z64aEa38jmcK6V6PjzwK9+kZePp5sQjYQH8oZVPCmiIjNGbk/4wJ/1Rq92cwxGtn6IofH0im11f/EpZRRUT7unFxHt8sbeztXbThBDtkATOLSpbV8r7iSc5eb6Qvj068WhYgNx1UwjRoiRwbjFVBiMJR8+x98g57O1UPBIWQMggL7nQphCixUng3EJ+ySzkn4mnyMorYXA/Dx4ecxudXdXWbpYQ4hYhgXMLKC2v5KOkM3z5YxbuHdUsmDaAgX27WrtZQohbjAROO6YoCskntcTs/4UrpXruH9yTySP8cHSQl10IYXnyztNO5V0uY8e+0xw/k49vtw4sfHAgvt3k+lFCCOuRwGlnDEYj+1My2X3o6k3RRvVldFAPuSmaEMLqJHDakYzsK/wz4RTncooY2MedWfcHyE3RhBCthgROO1Cur+IfcSfYe+is3BRNCNFqSeC0A2v/9SNns64w8k5vHgjtg7OjvKxCiNZH3pnaOG1hGWezrvBYxO2EBHazdnOEEOJ3yZnkNi4tXQdAcH8JGyFE6yaB08alns3HvaNaroMmhGj1JHDasCqDkZPnC7jDz10mCAghWj0JnDbsbNYVyioM9Pdru/euEULcOiweOFlZWcycOZNx48Yxd+5cSkpK6tTR6/VERUURHh7OlClTOHPmDAAlJSUsWLCAiRMnMnHiRD799FPTc959913GjRtHWFgY+/bts9j2WFNaug4bG/hDLzdrN0UIIa7L4oGzYsUKZsyYQWJiIv3792fTpk116kRHR+Pk5ERCQgJLly5lyZIlAGzduhUvLy/27t3L9u3bWb16NXl5eRw/fpw9e/YQFxfHBx98wJo1aygsLLT0pllcarqO3l4dcXG0t3ZThBDiuiwaOJWVlSQnJxMWFgbA1KlTSUxMrFMvKSmJyMhIAAYPHoxOpyMrK4shQ4Ywe/ZsANzd3encuTN5eXl89dVXjB07FrVajbu7O0OGDCEpKcli22UNxWWVZFy6wh1t+FbQQohbi0W/h1NQUICrqyt2dtWr1Wg05OTk1Kmn1WrRaDSmxxqNhuzsbIYPH24qi4+PR6/X07dvXz744AMCAwPr1G+Kpt6b29pO/ngRBRhxV080muqLctb8K6QvriX9YU76o5Yl+6LFAichIYHVq1eblfn6+taZTVXf7CpFUczKFUVB9ZuLTyYkJPDKK6/wzjvvYGdnh9ForLMMVRMvVpmfX4zRqDTpOdZ05KeLOKnt6OxkS25uERpNB3Jzi6zdrFZB+sKc9Ic56Y9aN9MXKpVNkz+ot1jghIeHEx4eblZWWVlJcHAwBoMBW1tbcnNz8fDwqPNcT09PtFotPj4+AOTl5ZnqRUdHs23bNrZt20ZAQAAA3bp1Izc31/T83Nxc/Pz8WmrTrE5RFNLSddzey02uAi2EaDMa9W5VUlLCihUrePTRRyksLGTZsmX1zi67Hnt7e4KCgoiPjwcgNjaWkJCQOvVCQ0OJi4sDICUlBbVajZeXF/v372f79u3ExMSYwgYgJCSEffv2UVZWhk6n4+jRowwbNqzJ7WsrsvJLKSiqkOnQQog2pVGBs3LlSjp27Eh+fj5qtZri4mKWLVt2Qytcvnw5u3btYvz48aSkpLBw4UIAYmJieOONNwCYPXs2er2eiIgIVq1axZo1awB48803qaioYM6cOUyaNIlJkyZx4sQJBgwYQGRkJNOmTePhhx/mmWeewdPT84ba1xbUXM7mDgkcIUQbYqMoynVPXEyePJnY2FjTv0ajkQkTJpiOVNqDtnQOZ92uH8krLOeVp4aaymRcupb0hTnpD3PSH7UsfQ6nUUc4156ANxgMTT4pL5pHZZWB0+cL5ehGCNHmNGrSwODBg/n73/9OeXk5hw4dYufOnQQHB7d020Q9TmdeRl9llPM3Qog2p1GHKYsWLcLZ2ZkOHTqwfv16AgIC+Otf/9rSbRP1SDurw1ZlQ4BPZ2s3RQghmqRRRzhvvvkmzz33HP/7v//b0u0R15GaruO2Hp1wdJB75wkh2pZGHeG098vEtBWFxRVk5hbTv7e7tZsihBBN1qiPyT169OCPf/wjd911Fy4uLqbyxx9/vMUaJuoyTYeW66cJIdqgRgVO587V5wsuXrzYoo0RDUtL19HR2Z6enm3rum9CCAGNDJyaa6JdvHiRqqoqfH19W7RRoi6jopCWoeMOvy6o5O6eQog2qFGBc+7cOebNm4dWq8VoNOLm5sbbb79Nnz59Wrp94qoLOcUUlVbKdGghRJvVqEkDf/vb3/jTn/5EcnIy3333HXPnzmXFihUt3TbxG6np+YCcvxFCtF2NCpz8/HymTJlievzAAw9QUFDQYo0SdaWl6+jp4UonV7W1myKEEDekUYFjMBjMbtms0+larEGirnJ9Fb9kXpbhNCFEm9aoczizZs3ioYceIjw8HBsbG+Lj43n00Udbum3iqpPnCzEYFbl+mhCiTWtU4Dz00EP4+vpy6NAhjEYjL730Uru+30xrk3ZWh4Oditt6yOVshBBtV6OG1HJyckhMTCQqKooHH3yQ6OhosztsipaVmqEjwMcNezu5QrcQou1q1DvY888/T+/evQHw9vZmyJAhLF26tEUbJqrlFZaRoyuV8zdCiDavUUNqBQUFPPLIIwCo1Woee+wxYmNjb2iFWVlZREVFkZ+fj5+fH2vXrjW7XA6AXq/nhRdeIDU1FUdHR9auXUufPn0oKSlh6dKlnD17FoA5c+YQEREBwOjRo3F1rf0G/pYtW+jevfsNtbE1SZW7ewoh2olGBY7BYCAnJ8d02+a8vDwacaPQeq1YsYIZM2YQERHBxo0b2bRpE1FRUWZ1oqOjcXJyIiEhgeTkZJYsWcKuXbvYunUrXl5evPHGG+Tn5zNp0iSCg4OxtbXF3t6euLi4G2pTa5aWrqNLRzXd3Z2t3RQhhLgpjQqcxx57jMmTJzNixAgAvvnmmxu6H05lZSXJycls3LgRgKlTpzJr1qw6gZOUlMSCBQuA6pu/6XQ6srKyGDJkCH5+fgC4u7vTuXNn8vLy0Gq1KIrC9OnTqaio4KmnniI8PLzJ7WttDEYjP58rYHA/DTZyORshRBt33cBRFIXJkyfTv39/9u/fj0ql4oknniAgIKDJKysoKMDV1RU7u+rVajQacnJy6tTTarVoNBrTY41GQ3Z2NsOHDzeVxcfHo9fr6du3L5mZmYwYMYJFixaRl5fHzJkz8ff3b/OX3knPKqKsooo7/OR2BEKItq/BwPn111956qmnePHFFxk2bBiffPIJNjY2xMTE8H//939mAXCthIQE00U/a/j6+tb5pF7fJ3dFUczKFUVBpaqd35CQkMArr7zCO++8g52dHWPGjGHMmDFA9a0Uxo4dy+HDh5sUOO7ure8KzPu+u4jKBkbc3ZMOzg7Xra/RdLBAq9oG6Qtz0h/mpD9qWbIvGgycNWvWsHDhQkaOHMnHH3+MjY0Nn376KTk5OfzlL39pMHDCw8PrDGtVVlYSHByMwWDA1taW3NxcPDw86jzX09MTrVaLj48PUH3OqKZedHQ027ZtY9u2baajrC+++IKuXbsSGBhYu2F2TbsjZn5+MUbjjZ2XainJaZfo1b0j5SUVlJdUNFhXo+lAbm6RhVrWuklfmJP+MCf9Uetm+kKlsmnyB/UGp0VfunSJyMhIAI4dO8bo0aNRqVR0796d4uLiJjfQ3t6eoKAg4uPjAYiNjSUkJKROvdDQUNMEgJSUFNRqNV5eXuzfv5/t27cTExNjNqR38eJFNm7ciNFoJC8vj4MHD3Lfffc1uX2tSUl5JWcvXZHp0EKIdqPBw4DfDmP98MMP/L//9/9MjysqGv7E/XuWL1/O4sWL2bx5M927d2fdunUAxMTEoNVqWbBgAbNnz2bZsmVERETg4ODAmjVrAHjzzTepqKhgzpw5puWtXLmS6dOnc+rUKSZMmIDRaGTRokV4e3vfUPtai58zClAUmQ4thGg/GgycTp06cfLkSYqLi8nNzWXw4MEAfP/996Yp0k3l7e1NdHR0nfKHH37Y9H+1Ws2rr75ap86ePXt+d7kvv/zyDbWntUpLz8dJbUdvr47WbooQQjSLBgPn2Wef5bHHHqO4uJhFixbh7OzMtm3b2LJli2lqs2h+iqKQmq7jdl83bFVyORshRPvQYOAMGjSIr776ivLycjp2rP6kfeedd/LRRx/Rq1cvS7TvlpStK0V3pYIJw2Q4TQjRflx3KpeDgwMODrVTcu+6664WbZCA1LNyORshRPsj4zWtUFqGDk83JzSdnazdFCGEaDYSOK1MZZWRk+cL6C9XFxBCtDMSOK3Mr5mF6CuNMpwmhGh3JHBamdR0HbYqG/r5yt09hRDtiwROK5OaruO2Hp1wdGjapXmEEKK1k8BpRS4XV3BBWyzDaUKIdkkCpxVJy6ieDi0TBoQQ7ZEETiuSlq6jg7M9PT1b360ShBDiZkngtBJGRSEtXccdvbqgkrt7CiHaIQmcViJTW8yV0ko5fyOEaLckcFqJ1HS5nI0Qon2TwGklUs/m00PjSmdXtbWbIoQQLUICpxWo0Bv4JfOy3N1TCNGuWTxwsrKymDlzJuPGjWPu3LmUlJTUqaPX64mKiiI8PJwpU6Zw5swZAEpKSnj66aeZOHEikydP5siRI6bnvPvuu4wbN46wsDD27dtnse1pDifPF2AwKtzRWwJHCNF+WTxwVqxYwYwZM0hMTKR///5s2rSpTp3o6GicnJxISEhg6dKlLFmyBID33nsPX19f9u7dy2uvvcZf//pXAI4fP86ePXuIi4vjgw8+YM2aNRQWFlp0u25GWroOBzsV/j06WbspQgjRYiwaOJWVlSQnJxMWFgbA1KlTSUxMrFMvKSmJyMhIAAYPHoxOpyMrK4v58+ezcOFCADIzM+nUqfoN+quvvmLs2LGo1Wrc3d0ZMmQISUlJltmoZpCarsPfpzP2drbWbooQQrQYiwZOQUEBrq6u2NlVXydMo9GQk5NTp55Wq0Wj0ZgeazQasrOzAbCzs+OJJ55g7ty5PP7446b6Hh4e9dZv7fIul5GtK6V/LxlOE0K0by12hciEhARWr15tVubr64vNNV9qvPYxgKIoZuWKoqBS1Wbjtm3buHjxItOnT+fOO+/EaDTWWcZv6zeGu7t1vt3//Zl8AEbc3RONpsNNLetmn9+eSF+Yk/4wJ/1Ry5J90WKBEx4eTnh4uFlZZWUlwcHBGAwGbG1tyc3NNTsyqeHp6YlWq8XHxweAvLw8PDw8+Pbbb+nVqxceHh54e3tz55138ssvv9CtWzdyc3NNz8/NzcXPz69J7c3PL8ZoVG5gS2/ON8ezcOugxlEFublFN7wcjabDTT2/PZG+MCf9YU76o9bN9IVKZdPkD+oWHVKzt7cnKCiI+Ph4AGJjYwkJCalTLzQ0lLi4OABSUlJQq9V4eXmRlJTE1q1bgephtNTUVAIDAwkJCWHfvn2UlZWh0+k4evQow4YNs9yG3SCD0ch/Mwq4w69LvUd6QgjRnlj8pivLly9n8eLFbN68me7du7Nu3ToAYmJi0Gq1LFiwgNmzZ7Ns2TIiIiJwcHBgzZo1AMybN48XXniBiRMnYmtry9KlS/H29sbb25vIyEimTZtGVVUVzzzzDJ6enpbetCZLv1REaUWVfP9GCHFLsFEUxfLjSK2QNYbU4g6ns+dwOm8sGIGrk/1NLUuGCWpJX5iT/jAn/VGrXQ+pCXOp6fn06t7xpsNGCCHaAgkcKyktr+Rs1hW5WKcQ4pYhgWMlP2cUoCjI+RshxC1DAsdKUtN1OKlt6e3V0dpNEUIIi5DAsQLl6t09+/m4YWcrL4EQ4tYg73ZWkK0rJf9KOf17u1u7KUIIYTESOFaQJnf3FELcgiRwrCA1XYeHmxMenZ2s3RQhhLAYCRwLq6wycvJ8gRzdCCFuORI4FvbrxcvoK40yHVoIccuRwLGwtHQdtiob+vm4WbspQghhURI4Fpaank8f7044qS1+3VQhhLAqCRwLulyi53xOsQynCSFuSRI4FvRzhkyHFkLcuiRwLCj1rA5XJ3t8u8ntbYUQtx4JHAtRFIW0DB2393JDJXf3FELcgix+5jorK4uoqCjy8/Px8/Nj7dq1uLi4mNXR6/W88MILpKam4ujoyNq1a+nTpw8lJSUsXryYjIwMbG1t+etf/8o999wDwOjRo3F1rb0Z0JYtW+jevbtFt60hF7TFXCnR099PLmcjhLg1WfwIZ8WKFcyYMYPExET69+/Ppk2b6tSJjo7GycmJhIQEli5dypIlSwB477338PX1Ze/evbz22mv89a9/BaCgoAB7e3vi4uJMP60pbADS5PyNEOIWZ9HAqaysJDk5mbCwMACmTp1KYmJinXpJSUlERkYCMHjwYHQ6HVlZWcyfP5+FCxcCkJmZSadOnQA4ceIEiqIwffp0pkyZQkJCgoW2qPFSz+rw1rjg1kFt7aYIIYRVWHRIraCgAFdXV+zsqler0WjIycmpU0+r1aLRaEyPNRoN2dnZeHl5YWdnxxNPPME333zD3/72N6B6CG7EiBEsWrSIvLw8Zs6cib+/P3369LHMhl1Hhd7AL5mFjL67h7WbIoQQVtNigZOQkMDq1avNynx9fbG55oT5tY+h+gT7b8sVRUGlqj0Y27ZtGxcvXmT69OnceeedjBkzhjFjxgDQo0cPxo4dy+HDh5sUOO7urtevdINS/ptDlUHhnkE90GhaboZaSy67rZG+MCf9YU76o5Yl+6LFAic8PJzw8HCzssrKSoKDgzEYDNja2pKbm4uHh0ed53p6eqLVavHx8QEgLy8PDw8Pvv32W3r16oWHhwfe3t7ceeed/PLLL5w/f56uXbsSGBhYu2F2Tdu0/PxijEblBrb0+r7+MRN7OxWeHRzIzS1qkXVoNB1abNltjfSFOekPc9IftW6mL1QqmyZ/ULfoORx7e3uCgoKIj48HIDY2lpCQkDr1QkNDiYuLAyAlJQW1Wo2XlxdJSUls3boVqB52S01NJTAwkIsXL7Jx40aMRiN5eXkcPHiQ++67z2LbdT1p6Tr8e3bGwd7W2k0RQgirsfgsteXLl7Nr1y7Gjx9PSkqKaRJATEwMb7zxBgCzZ89Gr9cTERHBqlWrWLNmDQDz5s0jNzeXiRMn8tRTT7F06VK8vb2ZPn06Go2GCRMmMGvWLBYtWoS3t7elN61euivlXMovlcvZCCFueTaKorTMOFIb01JDal/9lMX2hJO8/MQQvDUtd55IhglqSV+Yk/4wJ/1Rq10Pqd2KUtN1uHVQ49XV5fqVhRCiHZPAaUFGo8J/M3Tc0atLvbPxhBDiViKB04LSs69QUl4lVxcQQggkcFpU2lkdNsDtveTunkIIIYHTglIzdPh260AHZwdrN0UIIaxOAqeFlJZXcfbiFfr3luE0IYQACZwW899zBRgVRW5HIIQQV0ngtJC09HwcHWzp7dXR2k0RQohWQQKnBSiKQmq6jj/4umFnK10shBAggdMitAVl5F0ul+nQQgjxGxI4LSA1vfrunnL9NCGEqCWB0wLS0nVoOjvi4eZs7aYIIUSrIYHTzKoMRv57vkBmpwkhxDUkcJrZmYuXqdAbZDhNCCGuIYHTzFLTddiqbOjnK5ezEUKI35LAaWap6Tr6eHXESd1id+8WQog2yeKBk5WVxcyZMxk3bhxz586lpKSkTh29Xk9UVBTh4eFMmTKFM2fOmP2+qqqKhx56iP/85z+msr179zJ+/Hjuv/9+du7c2eLbUZ8rpXrOZxfJdGghhKiHxQNnxYoVzJgxg8TERPr378+mTZvq1ImOjsbJyYmEhASWLl3KkiVLzH6/ceNGMjIyTI9zcnJYv349H3zwAbGxsXz44Yf8+uuvLb0pdfycrkMB+veWCQNCCHEtiwZOZWUlycnJhIWFATB16lQSExPr1EtKSiIyMhKAwYMHo9PpyMrKAuD777/n5MmTjBw50lT/yJEjDB06lM6dO+Ps7ExYWFi9y21paek6XBzt8PXsYPF1CyFEa2fRwCkoKMDV1RU7u+rzGxqNhpycnDr1tFotGo3G9Fij0ZCdnU1xcTGrV6/m5ZdfbrC+h4dHvcttSYqikJqh4w6/LqhUcndPIYS4Voud2U5ISGD16tVmZb6+vnVutVzfrZcVRTErVxQFlUrFihUr+POf/0zXrl3N6huNxjr1m3pLZ3d31ybVv1bGpStcLtYzbIAXGo11jnCstd7WSPrCnPSHOemPWpbsixYLnPDwcMLDw83KKisrCQ4OxmAwYGtrS25uLh4eHnWe6+npiVarxcfHB4C8vDw0Gg3ffPMNp0+fZsOGDVy6dImjR49iZ2dHt27dSElJMT3/95bbkPz8YoxG5Qa2tNqh7y4A0NPdmdzcohtezo3SaDpYZb2tkfSFOekPc9IftW6mL1QqmyZ/ULfokJq9vT1BQUHEx8cDEBsbS0hISJ16oaGhxMXFAZCSkoJarcbb25vDhw8TFxdHXFwco0aN4plnniEyMpJ77rmHb775Bp1OR1lZGfv27at3uS0pLT0f764udOnoaNH1CiFEW2HxWWrLly9n165djB8/npSUFBYuXAhATEwMb7zxBgCzZ89Gr9cTERHBqlWrWLNmTYPL9PT05C9/+QuPPPIIkydPZsKECQwYMKDFt6VGRaWBUxcuy3RoIYRogI2iKDc+jtSO3MyQ2omz+azf9RPP/s9Aq02JlmGCWtIX5qQ/zEl/1GrXQ2rtVVq6DjtbFbf17GztpgghRKslgdMMUtN1BPTshNre1tpNEcVWs1IAAAuLSURBVEKIVksC5ybprpSTlVfCHXI7AiGEaJAEzk1Kk7t7CiFEo0jg3KS0DB2dXB3w1rhYuylCCNGqSeDcBKNRIS1dR/9eXZp8ZQMhhLjVSODchHM5RZSUV3FHbxlOE0KI65HAuQmpZ/OxAW7vJYEjhBDXI4FzEwqK9QT4dKajs4O1myKEEK2e3Af5Jswa60+VwWjtZgghRJsggXMTVCobHFTyZU8hhGgMGVITQghhERI4QgghLEICRwghhEVI4AghhLAICRwhhBAWIYEjhBDCImRa9FUqVdu/Flp72IbmIn1hTvrDnPRHrRvtixt5ntxiWgghhEXIkJoQQgiLkMARQghhERI4QgghLEICRwghhEVI4AghhLAICRwhhBAWIYEjhBDCIiRwhBBCWIQEjhBCCIuQwGnj3nrrLSIiIoiIiGDNmjXWbk6r8eqrr7J48WJrN8OqDh48yNSpUwkPD2flypXWbo7VxcXFmf5WXn31VWs3xyqKi4uZMGECmZmZABw5coSJEydy//33s379+hZfvwROG3bkyBEOHz7M7t27iY2NJS0tjc8//9zazbK6b775ht27d1u7GVZ14cIFli9fzqZNm9izZw8///wzX375pbWbZTVlZWWsWrWK6Oho4uLiSElJ4ciRI9ZulkX99NNPPPzww2RkZABQXl7O0qVL2bRpE/Hx8aSmprb4PiKB04ZpNBoWL16Mg4MD9vb29OnTh6ysLGs3y6oKCwtZv349c+bMsXZTrOrzzz9n/PjxdOvWDXt7e9avX8/AgQOt3SyrMRgMGI1GysrKqKqqoqqqCrVabe1mWdSuXbtYvnw5Hh4eABw/fhxfX1969uyJnZ0dEydOJDExsUXbIFeLbsNuu+020/8zMjJISEggJibGii2yvmXLlvGXv/yFS/+/vXuNieJqAzj+X1hwEUzFqtQEvPUiH2pAjQ1TaBFUhFBB0A+LrUG8X9ZYZOMqokTRFW/BVhJMDDFto4BEUaNoTYNa1qwltTGIgjG0odoikAYlRl2QpR+IU/fF16Zadgo+v4SEc3bOzDND2Ic5szynsVHrUDTV0NCAl5cXy5Yto7GxkSlTpvD5559rHZZm/Pz8WL16NXFxcfj4+DB58mQmTpyodVhutW3bNpd2c3Mzw4YNU9vDhw+nqampV2OQO5x+4NatWyxYsIC1a9cyevRorcPRTGlpKSNGjEBRFK1D0VxnZyd2ux2r1UpJSQnV1dWv9TRjXV0dR48e5fz581RWVuLh4UFhYaHWYWnK6XSi0/21xEBXV5dLuzdIwunjrly5wvz588nIyCApKUnrcDRVXl7OpUuXSExM5Msvv6SiogKr1ap1WJoYOnQoiqIwZMgQDAYD06ZNo7q6WuuwNGOz2VAUhTfffBNvb2+Sk5OpqqrSOixNvfXWW7S0tKjtlpYWdbqtt8iUWh/W2NjIypUrycvLk7/qgYMHD6rfHzt2jKqqKjIzMzWMSDtRUVFYLBba2trw9fWlsrKSqVOnah2WZoKDg9m1axcPHz7Ex8eHiooKxo8fr3VYmgoJCeGXX36hoaGBwMBATp06xezZs3v1mJJw+rDCwkIcDge5ublqn9FoJCUlRcOoxH9BSEgIixYtYu7cuXR0dBAeHt7rbyb/ZREREdy4cYPk5GS8vLwYP348S5Ys0TosTQ0YMIDc3FxWrVqFw+EgMjKS2NjYXj2mrPgphBDCLeQZjhBCCLeQhCOEEMItJOEIIYRwC0k4Qggh3EISjhBCCLeQhCP6vDt37jBu3DhKS0td+gsLC//VitHR0dFcu3btX9vfy3jw4AFGo5H4+HjOnTvX4/X6+npWrVrFzJkzSUhI4LPPPuPHH38Euq/ThAkTXvrYFy5c4Isvvnjp8ULI/+GIfsHDw4MdO3YwadIkxo4dq3U4vaa2tpY//vjjuVXBf/75Z1JTU9m+fTsfffQR0F05e9myZRQVFeHj4/NKx7527Rr3799/pX2I15skHNEvGAwG0tLSMJvNFBcX4+3t7fL6unXrePfdd1m4cGGPdnR0NJ988gmXL1/m/v37LFq0iJ9++onr16+j1+spKCggICAAgMOHD1NXV0d7eztpaWnMmTMH6F57pqCggI6ODgwGAxaLhQkTJrBv3z6uXr1Kc3Mz48aNY/ny5WzYsIH29na6urqYM2cOn376aY/z+e6778jPz8fpdOLr68v69evx8/MjMzOTpqYmEhMTKSkpwWAwqGMOHDjA7Nmz1WQDoCgKe/bscdkOYN++fbS2trJp06Ye7XPnzlFQUIBOp8PT05O1a9fi7e1NcXExnZ2dDBo0iPT0dEpLSykqKsLpdDJ48GA2btzI22+/zbp167h37x63b99mypQpREVFkZubi9PpBGDp0qXMmDHjVX/kog+ShCP6jeXLl2O328nLy8NisfyjsQ6HgyNHjlBeXk5GRgZlZWUEBwezcuVKysrK1OUOBgwYQFlZGU1NTSQlJRESEqKW///666/x9/fn1q1bpKWlqVNev/32G6dOnUKv15OZmUl0dDRLliyhpaUFq9VKSkoKHh5/zW7X19eTnZ1NcXExQUFB2O12VqxYwdmzZ9m6dSs5OTmcOHGixznU1NRgNpt79EdGRgKoi279nZ07d7J7925CQ0Ox2Wz88MMPmEwmjEYjra2tpKenU1VVxfHjxzl06BA+Pj7YbDZMJhNnzpwButdaOX36NACpqamkpaURHx9PXV0dJSUlknBeU5JwRL/h4eHBrl27mDVrFhEREf9obExMDABBQUEMHTqU4OBgAEaOHOkyjWQ0GgEICAggPDwcu92Op6cnzc3NzJ8/X91Op9Px66+/AhAaGope3/2rNn36dCwWC9XV1SiKQlZWlkuyAbh8+TJhYWEEBQUBqEU4a2pqXljNV6fTqXcRryI+Ph6TyURkZCTh4eEsXry4xzYXLlygoaFBvR4AbW1t3Lt3D4BJkyap/XFxcWzZsoWKigo+/PBD1qxZ88oxir5JPjQg+pURI0awefNmLBYLra2tar9Op+PZKk4dHR0u456dgvPy8vq/+382OTidTvR6PU6nE0VROHHihPp15MgRdb2igQMHqmOioqL49ttviYuLo7a2lpkzZ3L37l2XY/xv2XjoLh3/5MmTF557aGgoV69e7dGfn5/PyZMnXfpedD3S09M5fPgw77//PseOHXvulJ/T6SQxMVE937KyMo4ePcobb7zR45yNRiMnT54kPDwcm81GQkICDofjheci+idJOKLfiY2N5eOPP+arr75S+/z9/ampqQGgqanppUvTP11T5vfff8dut6MoCoqicOnSJerr6wG4ePEiCQkJPH78uMf4jIwMysvLiY+PJzs7Gz8/P/VO6ClFUbDZbNy+fRvofvDf2Nj4tyt2Lly4kNLSUmw2m9r3/fff880336h3bE/5+/tz/fp1urq6ePDgAefPnwfgyZMnREdH8+jRI1JSUsjOzubmzZu0t7fj6empJr2IiAhOnz5Nc3MzAEVFRaSmpj43LqPRSG1tLcnJyeTk5NDW1uZSFl+8PmRKTfRLWVlZXLlyRW3PmzcPs9nMjBkzCAwMJCws7KX263A4SEpKoqOjg6ysLMaMGQPAli1bWLNmDV1dXeoHDXx9fXuMX7FiBRs2bKCkpARPT0+mTZvG5MmTXbZ55513yM7OxmQy0dnZicFgYP/+/QwaNOiFsY0aNYr9+/ezd+9eduzYgdPpZMiQIRQUFPDee++5PMNJSEigsrKSmJgYAgIC+OCDD9TYMzMzMZvN6PV6dDodVqsVb29vwsLCMJvN5OTksHHjRhYvXsyCBQvQ6XT4+fmRn5//3Ck/s9mM1Wpl79696HQ6TCYTgYGBL3P5RR8n1aKFEEK4hUypCSGEcAtJOEIIIdxCEo4QQgi3kIQjhBDCLSThCCGEcAtJOEIIIdxCEo4QQgi3kIQjhBDCLf4EBz/Qh4kEwHcAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "Nc = range(1,11)\n", "kmeans_list = [KMeans(n_clusters=i) for i in Nc]\n", "score = [kmeans_list[i].fit(Trapnell_pca).score(Trapnell_pca) for i in range(len(kmeans_list))]\n", "_ = plt.plot(Nc,score)\n", "_ = plt.xlabel(\"Numbers of Clusters\")\n", "_ = plt.ylabel(\"Score\")\n", "_ = plt.title(\"Scores with numbers of clusters\", size = 20)\n", "plt.show() " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### 3. Clustering the cells and plot it by color" ] }, { "cell_type": "code", "execution_count": 57, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZoAAAEeCAYAAACzJ9OtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3Xl0VPX9+P/nXWZNAmFJgiwiognKqiIiIqgFg4giIEIFcd9XarHibvUj/VlaW2sFa/2hgiJIEcQqoIBIWQUVUFBBVgOEhJBt9nvv+/vHJAOTTMIkZBIC78c5Hs/cO/fe150M9zXvXRFCCCRJkiQpQdSGDkCSJEk6uclEI0mSJCWUTDSSJElSQslEI0mSJCWUTDSSJElSQslEI0mSJCWUTDRS3G666SaysrIoLi5O6HVM02TGjBl4vd6EXudUMnfuXK688kq6dOlCnz592LVrV71/zr/++itZWVmV/uvUqRM9evRg8ODBvPzyyxQWFsY8XgjB0qVLuf/++7niiivo0qULvXv35q677mL58uXHvP6//vUvsrKyuOiiiwgGg3V9e1I19IYOQJIqevTRR/nss8+49tprGzqUk8Ivv/zCk08+SXJyMjfeeCOqqtK6desG+5zbtGnDsGHDIq+FEHi9XtasWcNbb73FsmXLmDNnDklJSZH3FBcX89hjj7Fs2TJatGjBJZdcQnp6OgcOHGDp0qUsX76c22+/nccee6zK63788ce4XC4KCwtZvHgxQ4YMSeh9SkfIRCOdcA4dOtTQIZxUtm7dimVZ3HjjjYwfPz6yvaE+5zZt2vDggw9W2m5ZFnfddRcrVqzgnXfe4b777gPCieiRRx5h5cqVjBo1iokTJ+JyuSLH5efnc+utt/LWW2/Rtm1bbrzxxkrn/v7779m2bRv33HMPb731Fh9++KFMNPVIVp1J0kmuvJqoWbNmDRxJ9VRV5bbbbgPgq6++imyfO3cuK1eupG/fvjz//PNRSQagZcuW/P3vf0dRFN544w1CoVClc8+bNw+A7Oxsevfuzdq1a9m7d28C70Y6mkw0EgAFBQW89NJLXHHFFXTr1o3s7GxeeeUVPB5PlcesXbuWrKws/u///q/Svscff5ysrCy2bt0a2bZ7924efvhhLr/8crp06cIVV1zBc889R15eXuQ9WVlZrFu3DoALL7yQm266KbIvGAzyxhtvMHjwYLp27crFF1/Mo48+WumB8Y9//IOsrCxWr17NyJEj6dKlC9nZ2VXeS/l9zJ8/n9mzZ3PVVVfRtWtXBg0axPz58wFYsmQJw4cPp3v37mRnZ/Pee+9VOk9OTg7PPvssAwYMoGvXrpx33nkMHz6cmTNnRr1v7ty5ZGVlsWLFCl577TUuvfRSzjvvPEaNGsWyZcuq/Lwr2rBhAw888AB9+/alS5cuXHjhhdx6662sWbMm8p4rrriCiRMnAjBp0iSysrIin099f87xyMjIAODw4cORbXPmzAHgnnvuQVGUmMedeeaZPP300zz99NNUnFXLMAw+/fRTWrZsyTnnnMPgwYMRQkTOKyWerDqTyMvLY9SoUeTk5HDRRReRnZ3Nli1bmDp1Khs3buTf//43un58X5WCggJuueUWDh8+THZ2Nunp6fz000/MnDmTtWvX8vHHH2Oz2XjggQf46KOPyMnJ4c477+TMM88EIBQKceedd7JmzRq6devG2LFjOXToEJ999hn/+9//mD59OpmZmVHX/P3vf8+ZZ57JTTfdhMfjiarzj2XatGns3r2bq6++mt69e/PRRx/x2GOP8eOPPzJ9+nSys7Pp2bMnH3/8MX/84x/JyMhgwIABQLih+/rrr8fn8zFw4EBOO+00cnNzWbRoEc899xymaTJ27Nio673yyits376da665Bk3TWLRoEffeey//93//x4gRI6qN9YsvvuChhx6iefPmDBgwgKSkJLZt28ZXX33FunXrmDNnDueccw7jxo1j3bp1LFmyhL59+9KjRw969erVoJ9zdXbv3g0cSTilpaV8++23uN1uzjvvvGqPHTNmTMztX331FYcOHWLs2LEoisLAgQN59tlnmTt3Lg899BCaptU6XilOQjrlTZgwQWRmZopp06ZFbX/66adFZmamWLRokRBCiLFjx4rMzExRVFQkhBBizZo1IjMzU7z44ouVzvmHP/xBZGZmii1btgghhJg+fbrIzMwUc+bMiXrf888/LzIzM8WyZcsi2ypeRwgh3nzzTZGZmSlefvllYVlWZPumTZtE586dxYgRIyLbXn31VZGZmSlGjBghTNM85v2X38c555wjNm/eHNn+wQcfiMzMzErxrV27VmRmZoqHH3640me1cuXKqHNv3LhRZGZmilGjRkW2/ec//4lc79tvv41s37Vrl+jZs6fo2bNn1L3Hkp2dLXr16iXy8vKitv/rX/8SmZmZ4i9/+Uul61X8+9b357x3716RmZkpxo4dG3O/3+8XI0eOFJmZmWL69OlCCCG2b98uMjMzxTXXXHPM81flwQcfFJmZmWLDhg2Rbffff7/IzMwUS5YsqfV5pfjJEs0pLhgM8vnnn3PGGWdwyy23RO27++67adasGWlpacd9HcuyANi4cSPXXXdd5Ffk+PHjuffee495jTlz5pCSksIjjzwSVX1SXsW1YMECtm3bxtlnnx3ZN3DgQFQ1/trhCy64gC5dukRen3/++QB06NCByy67LLK9e/fuQLiqrNy1115Lt27d6NOnT9Q5u3XrhtPpjNnwPnjwYHr06BF53b59e8aMGcOUKVP48ssvq+wNZlkWjz76KHa7nZYtW0btu+iii4DaN/TXx+eck5PDP/7xj8hrIQSHDh1ixYoV5OTkcP755zNq1CiASFf62paSiouLWbZsGW3atIkqEQ0ZMoTPP/+cDz/8kCuuuKJW55biJxPNKW7Pnj14vd6oB165Nm3aRPVSOh7Z2dn885//ZNasWSxevJi+ffvSr18/+vfvf8wk4/F42LlzJ2lpaUydOrXS/vz8fCDcu+roB2CbNm1qFGP79u2jXpc3Ordt2zZqu8PhAIgai9GzZ0969uxJYWEhW7duZc+ePezcuZPvvvuOQCCAaZqVrterV69K27p16wbAjz/+WGWiUVWVgQMHAuGH9rZt29izZw/bt29n7dq1wJHEXhP19Tnn5OTw2muvRd1PUlISHTp0YPTo0YwbNw6bzQZAamoqQK3Hbn322WcEg0EGDx4clTgvv/xykpOT+eqrrzh48CDp6em1Or8UH5loTnFFRUUAJCcnJ/Q6GRkZzJkzhylTprBkyRIWLFjAggULsNlsDB8+nCeffDLyAK+otLQUCLclHf2Aqqj8Xso5nc4axVixN1M5u91+zGOLioqYNGkSn3zyCaFQCEVRaNOmDb1792bLli0xjylvhzhaeQml/J6r8tNPP/Hiiy9GGvRtNhsdO3akS5cu7Nq1q1KDeDzq63Pu1asX06dPj+u9rVq1wmazsW/fPkKhUCQBxXLgwAGSk5Ojvsvlvc3efPNN3nzzzZjHffTRR9x99901uAOppmSiOcWVV0lU1VPI6/Xidrtj7iv/hRjroebz+Spta9euHS+99BKmafL999+zYsUK5s6dy6xZs0hJSWHChAkxr1N+/Z49e8bs7XUimDBhAsuXL2f06NEMHTqUzMzMyANvwYIFMY/x+/2VtpWUlADVd0UuLS3ltttuo6SkhD/84Q/06dOHM888E7vdzsaNG/nkk09qdQ8n4ufscrno2bMnq1ev5ttvv41ZCiz3zDPPsHLlSqZOncqll17K3r17+eabb8jIyIiq+izn8Xj45JNPmDNnDnfddVeVPdqk4ycTzSmuQ4cO2Gw2Nm3aVGlfbm4u/fr144YbbuCFF16otL/812WsKUwqdoVdsmQJK1as4Pe//z3Jycl0796d7t27M2LECC677DI2bNhQZYwpKSm0bt2a7du34/f7K/2CnjdvHnv37mXYsGGVqrnqQ3FxMcuXL6dLly48//zzUft+/fVXAoFAzGS8efNmsrOzo7Z9++23wJEqtFjWrFlDfn4+t912W2TcSblffvkFiJ38j+VE/ZyHDRvG6tWrmTp1apWJZvv27axatQqXyxVpiykvzYwePToy+LOizZs3s3v3btauXUvv3r0TcwOSHEdzqnM4HGRnZ/PLL7/w4YcfRu0rr6e/+OKLYx7bvn17NE1jzZo1USWYL7/8kh9++CHqvTt27GDmzJmVxpSUN6i3bt06sq08gR098G7YsGEUFhYyefLkqPaH7du388c//pFp06ZF6vPrm81mQ1VViouLo9pt/H5/JEHHGkQ4e/bsSGIA2LlzJ9OnTycjI4O+fftWeb3yKsaKDf779u2LVHkZhhFX3BVjOxE/52uvvZbzzjuPlStX8swzzxAIBKL279q1i/vvv59QKMT9998fKUl+/PHHAFxzzTVVnrt8KpyK332pbskSjcRjjz3Ghg0beOqpp1i0aBFnn302mzdv5uuvv2bAgAEMHjw45nHlYzgWLVrEyJEj6d+/P3v37mXp0qVccMEFUaWUG264gdmzZzN58mTWrVtHVlYWhw4dYuHChbjdbu66667Ie8vbLp544gkuueQSxo0bx1133RUZx7FhwwZ69epFcXExCxcuxOfz8ec//znh7UxVcblcDBw4MPI5XHLJJXi9XpYtW0Z+fj5NmzalpKQEy7KiemcpisINN9zAoEGDEEKwePFi/H4/L7/8cpXtVRDuHdemTRvmz5/P4cOH6dSpE/v372fJkiU4HA4URalyYsqjNZbPWVEUXn/9de644w5mzZrF559/zmWXXUazZs3YvXs3y5cvJxQKMXbs2EjPyfXr17Nnzx7OO+882rVrV+W5hw0bxquvvsrnn39OUVERTZs2rae7OrXIEo1ERkYGH374IaNGjeKnn37i3XffZd++fdx777288sor1R770ksvcdNNN1FYWMj06dPJycnh1Vdf5corr4x6X9OmTZkxYwa//e1v2bVrF++88w5ffvkl/fr1Y/bs2XTq1Cny3nvuuYfu3buzcuXKSFuB0+nk3Xff5cEHHyQQCPD++++zfPlyzj//fN59991qf7XWh5deeombb76ZkpISZsyYwYoVK+jatSszZ87kuuuuw+/3R3qElbvnnnu46aabWLZsGYsWLaJ79+7MmDGD/v37V3stt9vNtGnTuPLKK/nhhx+YMWMGW7Zs4dprr+Xjjz+mU6dOrF+//pgj9BvT59y8eXNmzpzJSy+9xFlnncWqVat499132bBhA3379uWtt97i6aefjrSzlJdmjjVhaKtWrejTpw+BQCByjFT3FFGbylxJkmpt7ty5TJw4kYkTJ1YauyRJJyNZopEkSZISSiYaSZIkKaFkopEkSZISSrbRSJIkSQklSzSSJElSQslEI0mSJCXUKTdg8/BhD5ZVubawRYtkDh2qfiLDE1Vjjb2xxg0y9obQWOOGxht7ixZ1Mzj3lEs0liViJpryfY1VY429scYNMvaG0FjjhsYd+/GSVWeSJElSQslEI0mSJCWUTDSSJElSQslEI0mSJCWUTDSSJElSQslEI0mSJCVUgySaBQsWMHjwYK688sqYa5Nv3bqV4cOHk52dzZNPPhlZLfDgwYPcddddXHfddYwePZpff/21vkOXJEmSaqjeE01ubi6vvPIK77//PvPmzWPWrFls37496j0TJkzgmWeeYdGiRQghmD17NhBeCfLyyy9n3rx5DB06lMmTJ9d3+JIkSVIN1XuiWbVqFb179yY1NRW32012djYLFy6M7M/JycHv99OjRw8Ahg8fzsKFCykoKODHH39k9OjRAIwYMYJHHnmkvsOXJEmSaqjeZwY4ePAgaWlpkdfp6els2rSpyv1paWnk5uayd+9eWrduzZ/+9CfWr19PWloaTz/9dI2vX92UCmlpKTU+34miscbeWOMGGXtDaKxxQ+OO/XjVe6KxLCuyrjeAECLqdVX7DcNgy5YtPPjgg0ycOJEPP/yQxx9/nOnTp9fo+ocOlcacCiItLYW8vJJa3FHDa6yxN9a4QcbeEBpr3NB4Y6+r5FjvVWetWrUiLy8v8jovL4/09PQq9+fn55Oenk5aWhpJSUlcfvnlAAwZMiSqJCRJkiSdmOo90fTp04fVq1dTUFCAz+dj8eLF9OvXL7K/TZs2OBwONmzYAMD8+fPp168fp59+Oq1atWL58uUALFu2jM6dO9d3+JIkSVIN1XuiycjIYPz48YwbN47rrruOIUOG0K1bN+688042b94MwOTJk5k0aRKDBg3C6/Uybtw4AP7xj3/w73//myFDhvDuu+/y0ksv1Xf4kiRJUg2dcks5yzaaE0djjRtk7A2hscYNjTf2RttGI0mSJJ1aZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkomGkmSJCmhZKKRJEmSEkommpOIvnYNTYdcSYsOrWl2SU8cH81p6JAkSZLQGzoAqW7o69eResNQFJ8PAHXbzyQ/cj9KYSH+W+9o4OgkSTqVyRLNSSLpxeciSaac6vORNOmPYJoNE5QkSRIy0Zw09C3fx9yu+PwoBQX1HI0kSdIRMtGcJMx27WNuF7qGaNq0nqORJEk6Qiaak4T3D08gXK6obcLlxnfHPWC3N1BUkiRJMtGcNIJXXkXJ5L9jpqcjdB0rKRnvvQ/gnfh0Q4cmSdIpTvY6O4kERo4mcP0olNIShDsJNK2hQ5IkSZKJ5qSjKIiUJg0dhSRJUoSsOpMkSZISSiYaSZIkKaFkopEkSZISSiYaSZIkKaFkopEkSZISSiYaSZIkKaFkopEkSZISqkESzYIFCxg8eDBXXnkl7733XqX9W7duZfjw4WRnZ/Pkk09iGEbU/i1bttClS5f6CleSJEk6DvWeaHJzc3nllVd4//33mTdvHrNmzWL79u1R75kwYQLPPPMMixYtQgjB7NmzI/t8Ph8vvPACoVCovkOXJEmSaqHeE82qVavo3bs3qampuN1usrOzWbhwYWR/Tk4Ofr+fHj16ADB8+PCo/X/605+4+eab6ztsSZIkqZbqfQqagwcPkpaWFnmdnp7Opk2bqtyflpZGbm4uAEuWLMHv9zNo0KBaX79Fi+Qq96WlpdT6vA2twWK3LPjiC1i1Clq3hlGjoAbLEsjPvGE01tgba9zQuGM/XvWeaCzLQlGUyGshRNTrqvbn5eUxZcoU3n777eO6/qFDpViWqLQ9LS2FvLyS4zp3Q2mw2P1+Uq+/Fu2HzSgeD8Ltht9PoOijTzC69Tjm4fIzbxiNNfbGGjc03tjrKjnWe9VZq1atyMvLi7zOy8sjPT29yv35+fmkp6fz5ZdfUlhYyJgxYxg6dCgAQ4cOpbS0tP6Cl6K43pyCvmkjqseDAqheL2pJMU3uuBlE5WQuSdKpqd4TTZ8+fVi9ejUFBQX4fD4WL15Mv379IvvbtGmDw+Fgw4YNAMyfP59+/foxcuRIvvjiC+bPn8/8+fMj+5KTq64KkxLLOXsmit9XabuaewB1184GiEiSpBNRvSeajIwMxo8fz7hx47juuusYMmQI3bp1484772Tz5s0ATJ48mUmTJjFo0CC8Xi/jxo2r7zCleBxVxRlFiKr3SZJ0ylGEOLXqOGQbTd1xTXkN959eQPUdKdUIwOx4FodXbThmspGfecNorLE31rih8cbeaNtopJOH7/a7MHr2wnInlS0fnYRITaX43+/KEo0kSRFyhU2p9ux2iuZ8jG3NKvSv12JltCIwZCgkJTV0ZJIknUBkopGOj6IQuvgSQhdf0tCRSJJ0gpJVZ5IkSVJCyUQjSZIkJZRMNJIkSVJCyUQjSZIkJZRMNJIkSVJCyUQjSZIkJZRMNJIkSVJCyUQjSZIkJZRMNJIkSVJCyUQjSZIkJZRMNJIkSVJCyUQjSZIkJZRMNJIkSVJCHTPReL1ecnJyKm3ftm1bQgKSJEmSTi7VJpqVK1fSv39/hg4dysiRI8nNzY3se+yxxxIenCRJktT4VZto/vKXvzBjxgzWrl3LxRdfzLhx4yguLgbgFFsBWpIkSaqlahONEIKsrCw0TeN3v/sdl156KQ8//DCWZdVXfJIkSVIjV22iUVWVX375JfJ64sSJCCF49tlnZbKRGodAAKW0pKGjkKRTWrWJ5pFHHuHGG2/ks88+A0DTNF599VV+/vln2RlAOqEpRYWk3HEzLc9sTYvM9jTr3xv9m/UNHZYknZL06nZeeumlLFmyhFAoFNnWpEkTZs6cGUk+69ato1evXomNUpJqqOnoEeibN6KUfXf1rVtoOuIaDq9Yh9W2XQNHJ0mnlmN2b05OTqZZs2bRB6kqV199NQCTJk1KTGSSVEva5k3oW39ACQajtiuhEM5p/26gqCTp1HXcAzZl7zPpRKPt2onQtErblWAQ/aetDRCRJJ3ajjvRKIpSF3FIUp0xOneJVJkdzXI6CV14UQNEJEmnNjkFjXTSsc7sSHBANsLpimwTqgpJSfhvuqXhApOkU5RMNNJJqfiN/x/PI49itjoNq0kTAtcO4/DnXyGat2jo0CTplFNtr7N4yDYa6YRks+H73WP4fienSpKkhnbcJZqLLpJ13pIkSVLV4irR7NixgzfffJPCwsKoEszUqVOZOHFiwoKTJEmSGr+4Es3jjz9Ot27duPDCC2UvM0mSJKlG4ko0Pp+Pp556KtGxSJIk1QlLBPAEVhEwfgIUnPq5JDl6oyi2hg7tlBRXG0379u05ePBgomORJEk6bkJYHPbOxBfaiCW8WMKDN/QNh70fys5LDSSuEo1lWQwZMoTOnTvjcDgi26dOnVqriy5YsIApU6ZgGAY333wzY8aMidq/detWnnzySTweDz179uT5559H13U2bNjApEmTCIVCpKam8tJLL9GmTZtaxSBJ0skpYPyCaZUAR88wb2JahwiZe7HrpzdUaKesuBLNwIEDGThwYJ1cMDc3l1deeYW5c+dit9sZPXo0F110EWeddVbkPRMmTODFF1+kR48ePPHEE8yePZsbb7yRCRMm8Prrr9OpUyfmzJnDiy++yJQpU+okLkmSTg6GmQtUnhlCYGJYB7EjE019i6vqbNiwYZEZmg3D4Pzzz2fYsGG1uuCqVavo3bs3qampuN1usrOzWbhwYWR/Tk4Ofr+fHj16ADB8+HAWLlxIMBjk4YcfplOnTgBkZWWxf//+WsUgSdLJS1ObEus3tIKOqjSp/4Ck+BLNihUrGDFiBF988QVLlizh+uuv54svvqjVBQ8ePEhaWlrkdXp6Orm5uVXuT0tLIzc3F7vdztChQ4FwVd5rr73GgAEDahWDJEknL4ctC6VSolFQFBsOvWODxHSqi6vq7O9//zszZsyIVG9t27aNCRMm1OpBb1lWVBdpIUTU62PtDwaDPP744xiGwd13313j67dokVzlvrS0lBqf70TRWGNvrHGDjL0hxBt3aurt/HrwI/zBAyiA09GWtunDsNtSExtgNRrrZ14X4ko0oVAoqg3l7LPPxjTNWl2wVatWrF9/ZKXDvLw80tPTo/bn5eVFXufn50f2ezwe7r33XlJTU5kyZQo2W827Kh46VIplVe55kpaWQl5e41zyt7HG3ljjBhl7Q6hZ3A6aOEaTbPcDCqrioKgQoGHuuzF/5nUhrqozp9PJ5s2bI683b96My+Wq5oiq9enTh9WrV1NQUIDP52Px4sX069cvsr9NmzY4HA42bNgAwPz58yP7J0yYQPv27fnb3/6G3W6v1fUlSTp1qIoTVXEc+41SQsVVopkwYQL33HMP7du3B2Dnzp38/e9/r9UFMzIyGD9+POPGjSMUCnH99dfTrVs37rzzTh566CG6du3K5MmTeeqppygtLaVz586MGzeOLVu2sGTJEs4666xIR4T09HTefPPNWsUhSZIk1Q9FxDmCqbCwkI0bN2JZFj169Ki0vHNjIavOThyNNW6QsTeExho3NN7Y66rqrNoSzfz58xk6dCjTpk2L2r5r1y4Abr311joJQpIkSTp5VZtodu/eDcDPP/9cL8FIkiRJJ59qE81DDz0EwKRJkyLbgsEg+fn5tG7dOrGRSZIkSSeFuHqdff7557zwwguUlpYyaNAghg4dyjvvvJPo2CRJkqSTQFyJ5o033uCGG25g8eLF9OjRg2XLljF//vxExyZJkiSdBOJKNEIIsrKyWLVqFf369SM5OVlOty1JkiTFJa5Eo6oqn376Kf/73/+45JJLWL58uVxpU5IkSYpLXInmD3/4A7Nnz2b8+PGkpaUxZcoUnnzyyUTHJkmSJJ0E4poZoGfPnrz99tuR1x988EGi4pFOMEII9ovwAlKnKWqDlmTLa2tlYVqSGpe4Es23337LX//6V4qKiqLaZhYsWJCwwKSGt900eM5XRF5ZoklTVJ5zNeUsLa6vTZ3xF8L/nnCy/WMdYUDbS036veynaQfZTigdIYSJN/gt/tD3gIXDdg5J9p4oSs0n35XqVlxPjGeeeYbhw4dz7rnnyraZU4RXWDziLaSUIw/zX8u2zU5ujluJq9b1uAkB84e5ObxNxQqGv3u/rtD4z1Vuxqzz4JDrWEmES95FvnkEzRzAAMAbXEfQ+IVm7htR6un7KsUWV6LRdV1ON3OKWW4EMalcYjAQfGkEGGyr3ezdNbVvtUbxriNJBkBYCoYPfp5to+sdlZfslU49hrWfoLmP8iQTZmJahwkaO3DYzqrqUKkexJXmzz77bH766adExyLVgxJhsdM08B+je3q+ZRKIsT0A5FtWQmKLpXC7iohxOcOnkL9F/kqVwkLmAaDyGlmCUFkpR2pIcZVo9u7dy4gRI2jdujUOx5G1HWQbTeMRFIK/+EtYagSwoWAiuMnuZozdHbM6tLNmwwH4Kmx3lu2rL80zrZg/h3S3oGXX+kt40olNVZJR0BBU/E7oaIqsX21ocSWa8ePHJzoOKcFeC5SyzAgQAkJlVWIzgl7SFJVse+VqsPM0G5mazo+mESnZOICzNZ3z6zHRtLrIpPnZFvnHbUxzAAAgAElEQVRbjqo+UwW2JEHWSFltJoU59DMpUXQQ0d8JBRWnvVMDRSWVi6vuoVevXjidTnbs2EGPHj2w2Wz06tUr0bFJdSQoBAtDfoIVtvuBGQFvzGMUReHPrlRutSdxhqrRXtG4xe5msiu1XjuEKApc+x8vnUaHsCUJVLvgjGyD6xd5sSfXWxjSCU5RdJq5R6GpLQEN0FGVpqS6r0dV6qc9UapaXCWauXPn8tZbbxEIBBg4cCD33Xcf48eP54Ybbkh0fFId8AoRo1k/LAcLj2WyS1isDgVwKwpX2Jy0UjXsisJoh5vRDne9xluRPQUumxzgssmxWo0kKUxXm9MiaRymVQJYqEoT2Uv2BBFXiWb69OnMmjWL5ORkWrRowdy5c+XszY1IU0Uhhdj/4ATwoLeI33kLmRHy8WbQy1hPAQsDFVtnGpYZgNwNKgU/qchp9qTqaGoKmtpUJpkTSFwlGlVVSU4+Uk9x2mmnoWlawoKS6paiKNxkd/H3YOVqMgHsEGbUawP4U7CUbjYbrdUjX5ECy+KDoJd1ZpAWisoou5teuj3h8W//WOfL8U4EIExIbmNx9QyfHLApVRIy8yj1LyVk7UPBjsvWjSRHHxRFPq8aUlwlmtTUVLZu3Rr5hfDxxx/TtGnThAYm1Y3DlsVzviJei5FkjuXP3iNrnB+2LG73FjA35GOXZbLBDPGMr4j/1OK8NVHwo8rSB5wESxRCJQqGV6HoF5X5I9wxuz2fKA5tUfn8HiezLnezbLyDoh3y13WimVYRhd4PCFk5gEAQwBv6lmL/Zw0d2ikvrhLNE088wcMPP8yePXvo27cvDoeD119/PdGxScfJEIL7vYfJFVaMEQZhClTZfrNZHBn8NivopUSIqOFwfuDNgIerbS6ctaimECKcSALFEChW2PSGHcuArncE6Xh1OOLv37ZhVuhcJiwFz37Y8p5O55uMGGduWDmrNP77WxdmIBxrwY8q2+fZGPaJl5adT+Ds2Mh5g98gKn3TDQLGL5hWMZoquzk3lLgSTceOHZk/fz67du3CNE06dOiAzSbnDzrRrTaCHK4myUDVSQbCw98CQvDfoI95IR+xHukqCjstg3Pi6PJcIiy+M4KsNIL8uEvQ7Q4nrhw3IkDUWLt9K1206mUy/BMfnv0KwqycxIQJKyY6adnFS8Z5J9bD+6vHHBi+o2YyMBVCHsGqZx1cO+fEavs6mRhmLlQaRwMKOqZ1WCaaBhRXovH7/SxZsoTCwkIAvv76awDGjBmTuMik47bHMvEfx/EqcJ/3ML9WMUsAhKekaRZjHilLCBSIVLfOCHh4J+jFIFySueqmdtj36ggRuyR0YJ3Gto902g802btcx/BWfJ+CFRSsedHB0P80/MM7UAQ/vG0jZ5XG4Z9j1Ugr7F8n2wkSSdPSCFkHqJhsBAaa2qxhgpKAOBPNPffcQ3FxMW3bto1sUxRFJpoTXGtVrbbEciypKNUmGR3opOq0Uo88QLeYIf7mL2WbZeAErrW56KzpzAh6Ka8Ba/GtA/cBDbWKJFNu81s2rv3Qx3ev2yjcrkKlnnMKeZsa/uGdt0nlP1e7sY7R+9rRVHZeSKQk+wUEQlsqzA6gY9fOkKWZBhZXosnNzeXTTz+V3QUbm+N8rhVQ9fgbnfBUNH90HfkH/J7fw5uhI50DfMBHIR+LQkpUySrlFxtqMI7vkgK6C4b/18u0zsmIGHV3yac1bLWZEPDpOBdWoPr70V2C7ndXHDIrxcsSQaxYX4CjaGoqqe4bKPEvxbD2o2DDaetKsqNvPUUpVSWuRJOZmUl+fj5paWmJjkeqQ79Y8TeUq0B31cZGKxT5PVhVkrEBpysaphAsCfq52u7if6EA/w5V7oEWBIIVzuQ6EN96NuUPZmczOOu6ENvm2Igq1SiCHvcHMYOwcaqNre+FOxOcPSzE+Y8E62XmgOKdCt6DVSUZgeYEBGSNDtHjvsRMmWNaxZQGVhA0dqEodly2HrjtF5wUU+OHzAMU+xdjWgUcKlWw6x1JcQ5AVZwx32/TMmie9FuEEPKH8Qkkrn/xgwYN4qqrriIzMxNdP3LIu+++m7DApOO3pQaJxoKoJFOdEPCLMEHAtqDB50aAg5ZRbQFK5UjNuekSCFWgWBUfBEfO0KafQcchR3oI5H2nUbHqTNHCvdb+O8bF/rUapj+8f+MbdnYt1hm5xEvCp2U7xrOs87ggF4wP4WqRmGozy/JS4JmBIAAIhAjgCa7GsPJo6hqckGvGQwhB0NyNaeWjqc2wax1qnPhMq4RC74eIskpXAQSM7ZjeYpon3VjtsTLJnFjiSjT//Oc/ufvuuzn99NMTHY9Uh7ZW7Bd8DLWphAoAPx4jyZSfWyc8GHTfAA9d/toMtWKbhgKtLzHp+XCAtv2PROM5oFCyt/JDShgKP822EfIokSQDYAbC79+1UKfjNYnt/tzkDEFSK0Hpr2U3cBRFh7RuFrnrVVpfYiakhOUNbSx7EB/9FzAIGNsarEuvJfwc9s7GtIoAEwUNRXHTzD0aTU2K+zy+4HcxZmO2MKxDhMxcbFpGncYtJU5cicblcnHnnXcmOhapjtVXi0C8v9UtoA0KZ3dUOO13Xgr+loRV9oxUbXDBI0F6/q5y1Kpe9UWEGf6vopBHYf9aLeGJRlFg8HQf/xnkxgwcFaQCqgZfTQxX8QgD+k/2kzUyOh7DB7nfaOguQXoPi5rWdoXMHGKtw6KgYVj5DZJoSgMrMK0Cyn+6CCyEKKbE/wWp7qFxn8ewCoh1b6BgWkUy0TQicSWaPn368N577zFw4EDs9iNTjqSmpiYsMOn4+IQgFTjU0IEcxQI8CJ51NUUZL+Amha/fDhdrOg4xaJZlsW+1RtEuhRbnWqR3Dz+oXC0FLbtYHPxORRxV3aY5Be2uMNm1UMesUDrSXIKU0+uno0DLzha3/VjKd1Nt7Pufjvs0i52f2jC8SlRcXz7qJON8D6kdwwlp+3ydZY84I/WK9iaCITN9tDg3/rhVJYlYw24FVoOtwxII/UTl8rEgaO5ECCvuKjSbdhpBczdUGsFloWs1by82rVKC5h4UdBx6BxRFjgWsL3ElmmnTphEMBnnhhRci2xRFYevWrQkLTKo9r7C4x1tIcUMHEkMh4ZKWA0g7By78fbgE4y+AWZe5KdlT1iXbgozzTa5+34fugoFv+PjoWjfBErAMUFQ47SKT/i/7yfkqCcNLdBLSqdf1amxJcOGjIXg0xPb5OrsWVn6IWQb8+IGN3k8GObxdYelDzqiBnSEPzB/h4uZNnqi2pdwNKltn2jD9Ch2Hhmj/GxNFBU9gDQHjZyoX9zRsajq61jIxN3tMVZVxa9ZO5bJ3xRvcUDbav/xYHbvWHr2G42I8gXV4gqtRCHeTL0HQ1DUMu972mMdKxy+uRLNp06ZExyHVoXlBHwcsk0Q+Zo9u3K8pS4hwndNRvnzUSeF2FSt0ZPuBDRrr/mynzzNBmrQXjF3vYe9SjdIclfTzTdJ7hCMY9omXz+9ykb9FRVEgpZ3FwKl+nA00Ri9YqsSch00YCoGi8P1tnVF5ah1QCBTB/GEuzh0b4uxhBt++buebv9kj09ns+ETn9AEGV0zdiye4jlhVS3atPU1cV9X5fcXLrp9FwKhYqlGwae1q1CFAVVw0TxpT1qNuJ5pmx651Jcles7WwQuZ+PME1gBk1RU2Rbx4tk+9BUeLrBSnVXrWf8Pz58xk6dCjTpk2Luf/WW29NSFDS8VlhBBPaPnM8SUYFihEcvRSVGYRdi/WoJANg+hV+nGmjzzPhu9FscFpvk18WqPz6lY4QBhnnWTQ9Q3D9Yi/egwqWCcmnVf7lHCyBkr0qyW0tHAmuUWp7qYEwHZW2a47wom0A3jwVYcSYWseAA+t08r/X+G6KRdEOFfOoMTohr8LuL3R2Ls0n9eJY7Rc6Dv1MVKXy9etLiqM/ITMHIXxlHRVsKIqNJs6BNT6XpjahqetqANLSUsjLKznGEZX5Qj8Qu60HguZuHHrHGp9TqplqE82ePXsA+Pnnn+slGKlupNTx+AkFIonBguOa1kYBmleIT5hUORNzoFjhne5JaA5od5nBTx/aQISTk/ZXO+0HGFz5Lz+KCu70yglGWLDqOTvfv21H1cEKwTljQ/R9MYBazaQC+9ZobH3PhuGDs64z6HCVUe37j9bkdEHXu4J894/y9sxwG4oZhPwfVNr/xuSMKw12fKpjeCrPdgBgeJXwbAgxeukaXti7OIPUi+OLp76pqpsWSbcQMLYTMvPQteY49cyGaxMRFXvllW1GII4xCFSqG9UmmnXr1jFu3LjI4Cdx1IpTx9NPfcGCBUyZMgXDMLj55psrTWWzdetWnnzySTweDz179uT5559H13X27dvHhAkTOHToEB06dGDy5MkkJcXfXfJUMcLuZJMveFwJ4WgzXKl8YgT4LOSn8DinG7AAjxCkHvX90V2Q1t3i4DcVp5kRCBM8+8OJ6Ye3owdsGgbs/kJn20c6mSPCDwwhYP9ajQPrNNzpFiV7FX54x47pVyK/aX9834arueDCCUEsE/yHFBxNBVpZIWDdn+1895odww+IcAmidR+Dq2f44+4VVrKn/I3Kkf8L+Pr/c3DOaIMOgw1avGGR/70a1TU76rMKKShq5c9b0cHZpAnhJYsrPigFdv3M+IJMIEXRcdo64bR1qvU5hBCEzF8JGDtQFDtNQz0JDxeuGYctE7+xHSpVJlvYdTlkoz5U+89m7NixjBkzhoyMDNxuN+PGjePWW2+lefPmtG/fvlYXzM3N5ZVXXuH9999n3rx5zJo1i+3bt0e9Z8KECTzzzDMsWrQIIQSzZ88G4Pnnn+fGG29k4cKFdOnSRS5VUIWLdAcjbXWz/LIKPOUvYU7Id9xJBsAOfBLy8TtvIWNycpgX9BIUgtN/U/4QENH/jxrUWfmBbHgVtr4ffviYIfhklItPRrtYO8nOiolOvp4cPZMygOFT2PiGnS0zdKadm8T0C5J4KzOZlc/aKdqt8O2r9vAx4kjpYs8XOu90S2LP0mMXa0KlsOMTPWa8woJfV2gEi0FYoqxXmqDKhvIYOUjV4ZxRdtz2Cwj/VlQI/6U0kh2Xo6n1MCVCggkhKPb/l0LfPHyhDXiDa9n+6xR8wZp3QLJrZ2LX2nMkSSmATrKjP6riquZIqa5UW6LJzs4G4K233uKDDz5AVcN56bLLLmPUqFG1uuCqVavo3bt3pGt0dnY2Cxcu5IEHHgAgJycHv99Pjx49ABg+fDivvvoqI0eO5Ouvv+af//xnZPvYsWOZMGFCreI42R2KNbikFixgZx2dC8I9zmYEveHSlglbCbDYG6DLX9oR/VStQYm57Bm95R0b+9dqkcQS8lQTRzGseMIZVZr4fpqdgh9VYi/GGJ5qZuGtLq6d4yWtmkH3e77Uqh338+0/bSx9yBkeQ1ShBHf0a90t6PlogA2vOCI9mK0Q9H0hQPMsC7gEp55FwPgFFBWnnommNkUIC4GBUvZg9Rtb8AU3IjBw6lm0sPpVHfwJImjuJGDs5EgpxEIIi5LA5zhsNWuDUhSFpq5rCJq7CBjbUbDhsnWusou0aXnwBtcRNHahqkm47T1xnAClxMYsru4Whw8fJhAI4HKFs7/H46GoqKhWFzx48GDUnGnp6elRvdoq7k9LSyM3N5fDhw+TnJwcmQKnfLtUmSUES41jTCVcBY2qmk2PX3nx+egqvQCws8Sis6hRaonQ3YJOo8MPo63v2yqVXqqiOahUZWX6FfYu09GrKQwaPlj/VztdYyQaIeCrPzj48QNbtat/Hvq+fARq5VgVTaC7wgmlx31Bzn8wRLc7Quz9UscIQLv+RlRvOl1rGenGLIRJiX8pvtD3gIWmpKAqTQlZ+yivYvMED7Nj33aa2EfXaHljIQTe4AZ8oW8QIohNa0uS41KEKC0bGJpaq2lmquIP/Ujlqq7yfVtx23vU6HyKouDQO+DQO1T7PsvyUOCdjhB+wMI0D1PkO0Cy/RLcjgtqdE3piLgSzZAhQ7jhhhsYOHAgQggWLlzIDTfcUKsLWpYV1b5TcfK7qvbHmiSvNu1ELVpUXa2QlpZS4/OdKI6Ovcg0CZXm1/gcTkBVFLzi+KvIYklVVbxC4K9wfn+yhdBAOUa7rKKGH+aaA8wA2NzQ4QqFPve4UDXQ4nluKmBzAapSaZBn+RuMalenVijeES4pVPy+bPsMfv4QzGobx5QK/4/el9wKBv0NTu8Lya0chEccwWnjqjtn2K+5c/EbWylPKqYowhQVfxCaBEMF2FP3kprc9dgnLbMv7xO8oU0IEX74B80dBL07Cf80sVDQ0PUkOrS+HZt+/FV3gYMOAqWx9hiUBr7EUnZwesYoNC325Jq1deDQGoQIEN2v0sATWkW71n1Q1XAHD0uYlHp/Jhg6jNPRiiRnh2M+jxrz8+V4xZVoHn74YTp37syaNWsAePzxx+nfv3+tLtiqVSvWr18feZ2Xl0d6enrU/ry8vMjr/Px80tPTad68OSUlJZimiaZplY6L16FDpVhW5QdpbbtOnghatkzm8wMF7LRM2qkaPVS9ViWTAER1+KhraajsjVENpypwuGuAFt9WfGiU/+oPx6RocNlf/RilCv4ChTaXmpx2kcmhgvC7zxpp49DPldtk3OkWGReaHPpBo/k5Jj1/F2TNCw5+/aryRJ0Aqk2gqJQlogo/blRB8y4GYKv0fVk71UnIE6uxuvwzPdYPI0H6hSHS+gfwAb68Y7z9KJblpcizhXj+6kKEOFTwEyHfGXGe28Nhz3cxzi0oT2oCi5BRyK6cj0h1Xxd/4FVQzEzgByp3dgCw8Pr3sOPXDytdyzAP4Q19g2kVYdfa4bJ3r3Km51gKPT8T8zMUCgcO7sKmnYZplXDYOxOrbGlYBQ1NbU6qeySqYq98LI33+VJXyTHucu6AAQN46qmneOqpp2qdZCA8nc3q1aspKCjA5/OxePFi+vU7Umfcpk0bHA4HGzZsAMJjefr164fNZqNnz558+umnAMybNy/quFOVR1iMzMnhCV8RUwKlPOMr4jZvISm1KO0dnWLqeu5bG3COqmHGaLzQVWjirzyNSvSvfwUrpPDVY04UXdC6j0mrC82ocZ9dbgnRvJOJooUb1xVVhNeBuS+Iq4XgrOtCXPJsgPTuFr2eCBAep1c5HiukcMaVBmcNNVDt0fs1J1z4aOxRStVVlx2bABUufrJ2I6BMUYJCvFVhGqoa/wMkPOdYfMLTzBx/5atdb4fLdh5UeU8WQXM3ljiyumrA2EmB9z38oe8JmXvwBNdQ4HkH06qmsa4CTYn9uQjMsul+oNi/EEt4CFftWQhCGFY+nsDquK9zqqn3BSsyMjIYP34848aN47rrrmPIkCF069aNO++8k82bNwMwefJkJk2axKBBg/B6vYwbF643ePbZZ5k9ezaDBw9m/fr1PPLII/Ud/glnSsDD9mAQH+GvvQ/YJ0w8x1kysQNdVZ26GvYXAuYbgajVOjXC43OeczfhlrdDtLviSJKoiumD/z3p5NObXEw/P4nCX45kmuLdKgU/a2UPfAVhKRh+WPOSgy3vhrssz7o8iW+n2Pj8LleVXZU1lyCjp8nAf/m55I8BkttYaE7Bab0Nhn7krXIusqyRBro7VuwKsVO3OOo/cDUXBEtrl+I1NTVq1Ht1FEXFZesS97nDy9/FmzzqrkSc4ryUFkk3o1T5LVTKqrnCJfES/2LCJaDyGEws4SubFSA+R3ryHU3Fpp2GpjZBiFDZRKYV79PEb2yJ+zqnGkUksq7kBHSyVZ1dVZKPL8Y/7ljlg5rQgdvtblIUldcCpXU2JqciB3CxbucZZxPUsuLJW1lJBA7H8RtIETQ5QzBmjQdFgYW3O9n5Xz1qzrPYqqnKUgXOZoIbV3twVjNnbKzvi7Dgi/uc7FyoH9XOU3nNHc0RDsE0qNB9O3ztmzd5ImN6aqLEvxxfaCPR1U06qpKEJTwoKKDYOD1jBD5P/JNSeoPrKQ18Fee7VdJT6vYHYLF/Mf7QD1T8RiuKi5ZJd6MoKqZVxCHPO8SqalOVFFomxz/7vC+4mdLAciBckrFprWnqGoKquLBEkPzS14k1N4aCk7SU+2Kes7E+X+qq6kxO8tPIWdWkk+Oajwx4J+glBSVhSQbC7UJrjCBLjQADbOG69K53hNjwNzsidIyEIRS8uXDoe5WWXS1yv9biSDLlKicAVYf2Vxr0fSFQbZKp8owqDJjiJ3e9ysrnHOSur9zNWXOGF0NT7bDpX3asqJoyBTMUno6nNssbJDv6oSlN8IbWYwk/Nu00Uhz90bU0DOswCBNNbUGyuwk+T00eehrx9kfU1Zq3mx5Lkv1iAqHtZdPZlMegk+IYEOnlpmCnqp9WSg2n43HZu+K0nYNhFaAqLrSjqhlVxY6uZmBY+yscpeKwZdXoOqeSxr/W6ynuYt1eZS122+P485ZPNZNXh1UhVfEDf/WXsKasS3bP8UFsrviuq2gQLJvGJSnGHGdVHBVz2zljQlz1tp+UdrW/Z0WBZpll86nFOI2iwIWPBTE8ClawchyWAb5Dtas+UxQFt+M8WibfSXrKgzRzXx8ZK6KrzdC1ltX2jBIihC+0BU9gNQHjF0RZo5NDPzvOCHRSnJfVKvbqaGoKzZNupkXTi9DVDBx6Js3cN+C0HYlLVV3YtDZUfqTpZW09NaMoOjYtPSrJlGviGoSCkyMDQG1oSlOSHX1qfJ1ThSzRNHIPOJLZIoooNM2oiTQFsA8roeNi6pIXeNZXzEOOZC7ckkSwOM6HrQXpPcJ3eMH4AJ/f7arQ6yz2eJVKVEHrS47/kxICPr7ezaGtR0+nE844mhOu/LcPRxNo29/kx9kixlxn0LpPfHGEl0vegS+4OTwY03YOTr1TjcbHlDOtQgq8H5R1Xw6hYENVU2jm/i2amkyKI5uSwEKiy8g2VMWBJfzoakuSHf2waa1rfO3Y92YBVmRmZU1NIq3FQDSrd5XHNHUOptA3F8MqQEFFYOKyda5Re1Q8dLUZLZLvIBD6CdMqRNcycOgda/W5nypkomnkWqoan7ZrxzW797C/QrenRE8XGJ70pKphdTUXAKYGPJT+rgVVJ4dw4lA0gWYPr1qpl/Ve7TDI5OLnAqx50YEww4MeraPbhqsjIO87lbOPs2fuwW9UDm9TK5RWFFSboMd9Qc4YGE4iZ1xp0PJck/zvj8xkoLsFHa8xUFRB/g8qzTtZ1U7kWRr4smxwZvgvEDL34de2kOoaUeOBk8W+RQjho/zDEoQwrSI8gZWkOK/ArreDQMXK2BBCaKQl31tnE2ZaIkipfyn+smUGdDWNFOcAbFqrYx6rqm6aJ40lZOZhiRJ0NT1h0/Goih2XPf5xSKc6mWhOAi5VrXVbTG01ReGPriZ8YwR5J+Q79gFxMvxQsLWqh6Qgua1Fs7Msmpwh6HxziJado++8660hzh0TojRHweYWvNc7udqpaI6cWuH7t+1c/Gyw4lI5NVK4Q415vBVSKNpx5L5UHYbO9bFlho2f5ujo9nD70Nb3bHz4myQUFXSXYOBUP237VS7hGNZhfKFNRJdXDULmAYLmrhpNmSJEiJC1n5g9qUI/keK8An9oS4z9ABYBYwfOOmqfKPLNI2Tup/y+DOsghd4PaZ50MxBfw7RNSwNqvgKnlDiyjeYkcb5uq7c/phu4xZFEF81W54urWZqoYq4xAAXPfpW8TRoX/j5YKcmU0+zQtIPAnRFemTPe/neGr6wEdBxadLJijqfRXYL086IThuaArreHuP4zH9fM8fHd63YKf1ExfAohj4IvX+XTcS5K91XOXCFjL7FLfSECxo7ju4mjlV0iPG6kcsITmFii2qkU4maY+YTMA5WuIzDxBb+rk2tIDUMmmpPELfYkklDqpYjqBd4IlPK0r5gLtbpLcDZggMtJx2sMNHvs5CDM8KwA7/RIYunDDvyHqz/nGVeatL3MBOXYyabZ2VbUEsq10bKrRaueJprjqOupAt0tOOfGqtPy3qVa1IzR5SyDyOzUR1MUZ9myxBWpqErNZu5WFBs2tTWVE5eKUw9P82/X20Um6axwNDatbpZDNq3CKu7JwrBqPqWSdOKQieYkkaFqvJXUjGtsTs5QNS7S7NRV7fRA3VHpi+IHvjGDoCh01Y4/vanAAM3Ow85k+r/sJ+288od17AGQwlD4eY6NuVe7j1kKufwVP+40ERlMqTnLzluefJTwDAL9/lS7iUgrGjzDR9fbgziaWehuwZlXGYz83IujadXHeA+qxBpQbwUVPPsrl1zCk0PGKtGouGydaxxzE1c2quIuSyYKCjY0tTlJZT2p7FoHNC2N6Np2Gw79zLKqquOnaS2rGHSqocfRRiOduGQbzUkkXdV42HmkHvtpTyErrOOr3NKBszSdZUagUjuQH9hgBHnFlcozvmJWmWVLLhN/B4EMFO50JnOxZiOpvOU7BYYv8LFzkcai210VxpocYYXCVWm7v9DoMKjqnlopbQRj1nrYNs/GoR9UWnS2aNnFZONUO3mbwo3uFzwSJK1r3bR06S7o81yQPs/FP51Mq4tMYg2dtiWJmG00imIj1T2CQt+88PrPAAhSnIPQ1GoyWhU0tSktkm4nYGzHtIrRtZbYtTOOjFNRVJq5rscX3FQ2Al7DZe+GUz+3xteqiq6mYtc7EDR2caQri4KCjttWs9mapROLTDQnsasdLlb4ji/RDNQdNFFUbFTuxWYHmqoqqqLworspv1oGa40gThTWGEHWmcHIlDMVZypQAaeiMMmVyplVlIg6ZJtkXGCSu0GLOeYEIOSDgh+rTzQAtiQ4d0z0ZzFwSiKHotZM80yLjtcY7PivjuEN36vmFKR2tOgwOHaRzaa1omXSXRjWAYQwsGmtI92Ba6N8Vczq9rsd5+N2nF/raxxLU+fVeIJr8YU2IoSBXTudZGd/VLVuFvKTGoZMNCex8zU7OrXv5twCeNCZgp/yKGEAABUuSURBVEDwaoxJbVTgCv3IzLhtVZ229vBXKtvm5MOgl3khPz4huEi3c7FmZ4HhZ79lcq6mM75VOilF1VdXDZ7uY+lDTnYt1st+uEcnHJsbUjvWd5+7xPjNP/y0vVTnh3dsGH6FzBEhutwWqrbdSFHUOhu7ciJQFI1kRx85+PEkIxPNScyuKNxpT+KNoKfG3Z9twD/dzXAr4QkhX3Y15SlfUaRKTAGedTWhhRq7mU9XFH7rSOK3jqTItsOWhaYopCgK3TUbrex28qg+0TiawFVv+/Hmw6zLkvDlE5kfTNEEjqaCM7ITPWKofigqdBpt0Gn0yXE/JxLTKiFg/IwQBg79zCpX15QSQyaak9woh5uWisrrwVIOlzUCVJV07IBeNtn8U64UWh1VpdVVt/Gf5BZsNQ0sBOdqNmw1GHAyPeBhetBb9oVTcCowLegk3tYEd0sYudjLl793sHdZ+Cyn/8bgsskBtNhLgEgSAL7gj5QEFpW9svAE1+KydUvIdDlSbDLRnAJ+Y3fyG3u4ius5XxFfGpUbqZ3ATXY35+l2MlUdXVGwhGCDGWKnZdBO1eil2emq17z/77dGkPeCXoJQNk2OwCvg7v37ed/VLDJr87EktxYMed+PVdYcU92oeUkCsISfksBiKg5s9YU249DPxq63aajQTiky0ZxiLtIdrDGClWZktoDf2Jy0Knt6lwiLh7yF7LdMwutJKjRXFF5zN6NZFdVlVZkX8sWcAbr4/7V370FRnXkax7+n7yAQRtJeiJG4XsiMJqbUGFddXLNRQAdQ1Fpr1E2ZcYymLMlonPISZ0eToCY1XqITZ11MWWuMbk1WISgwGdyxzMREMZXEqBRjOZuJitwEBIGmb+/+gXZEQLl4Ghp+nyqr0t2n+zyHIjx9bu/r9ZLvdTO8jTevSMGI1mq4gq35G1sd7nwpGj+R+2h6mOdNVvoajNx9tMkGTLurZADW1Vbyf14PDhouJqhDcV152Vbf9jk1brUw5ZEBqO1Z0yGJLuVhzyMrWiJ7ND1EideDB+inGdgd/CMOO2v5s7ueYDRmWIJ43vTDnB0nXPWc8za9XNgLfOZ2opS673Dz9/pnk5ULHleTvRo3tHlvRoi2sJieoPmbfk2+UQ+E/qRourmrXg+/qbvJ914PGtBbM/B6UBjzrb2YYwmmWHmI0AyNiiOtvuVRKL20euB9n1izjWMuB9953Tho2JMxA2sjIgiul2+VQj8GzUaYLZYqxx/5YdpsA0HmkXLYzI+kaLoBpRT5HhffelxEaAYmmqxYNQ23UiyvraBCKd93uuvKS0ptJbNNQaS769DQ8KCIM9tYbg3BpGkUNzcWym2jDaZWn7y/w6Jp7AwO58/uej511xOuGUg02xgfFhaQ09uKwGIzR2M2Pka9+9Lty5sHYTI+2tmxehQpmgDnVoplxcV8UVt7+6Q9bOcWO4LDKfR6qL2rZHzvAQ657wzt3/DqH10OgtBYagvhMc3Id82UjQastLVvDnGzpjHVbGOq2fbghYV4yIyGEIItbZ9pUzwccjFAgDvmquOLurq7TtpDNYr1dTcpVZ4H3A75g3ogw1WHRykW23px7yzrJuBVSy/6P4QBNIUQPYsUTYA76qrH0cyVWzeUl0c0QytnYmngoqFwxpus/HtQGAM1Iyagv2ZgtTWUJKuMNyWEaDv5ehrgvPepkuvNXDl2P3bNQNDt/x5vsjI+5N79GiGEaDvZowlwU0w2bM2cnA/VDNz0tn6EMyuw3BbSpsuWhRCiNWSPJsDNtATxheamoL6eOhoKw4DGb4LC+N7rIcjtoK6Z91mApw1mrigPjxuM/JslmKdNMmiY6FrcnnJu1f8Fl+cqlY4QrIZR2MzD5QtRgJGiCXBWTeO/IiPJKrrBN24XdoOR501Wwg0GBhtM/KdWQ73yNhpIMxjYE/wjBsiJfdGFebyVVNR+iLo9Qp7T5cDJ/+JRVW2aRkAphdtbhNt7A5OhNyZDfykqP5O/NN2AUdMYZ7IyztT4nIpN09hsC2Ojo4qryosG/JPRwgpbKOFtHK9MCH+rcZ5BNZmr1U2t8yzBljEYtAfvgXuVk8ra/8HjLfOdzTQZehMePBuDJucg/UX+2nRjlzwuXq27SbFq2J9RwF88To4561Ayxpjo4lzuQpobPkbDgMdb0arPuFV/Ere35HZhNfxze8u45fjzQ80q7k+Kphvb7Kim9p7vhF5gr6uW/7jPMDNCdAVGQ3izzys8GLWQVn1GvSufxlMEAHhwuAvky5YfSdF0U7XKy3ctXN7sBQ676qhW3WMKZNE99bKOpenRfSNW0z9gMPRq7i1NqCYlc4f87vuTFE03ZUS778CXZjSutvE+GyH8yWyMJMwWj0HrBRjRNCM2UzRhtrhWf4bFGEVzQ8CajY/LBQF+JBcDdFNWTeM5o4VTHmezt3S6UPTR5HuG6Nps5qFYTUNQqha7PYIbN5qbQq9lobbnKa/58PYBZDcNk5WbCLX9iy55RfOkaLqxVbZQXqmtoPCeQ2QW4B9NFiJkqkoRADRNQ9N6YTCYodm5WltmNDxCRMhC6pwX8XiLMRrsBFlGYNBkcFd/kqLpxsINBg706s1/O+s46KzlFgojMOX2lABC9AQGzUYv66jOjtGjSdF0c5qmMdcazL9agqhGYUPDIsemhRB+5PeD9IWFhcybN4+4uDiWLl1KTU3Ty2ydTierVq0iPj6emTNncvnyZQBqampISUkhISGBhIQEjh075u/4AUvTNMI0g5SMEMLv/F40GzZs4Gc/+xk5OTmMGDGC9957r8ky+/fvJygoiOzsbNauXcuaNWsA2LNnD5GRkWRmZrJv3z42bdpEWVmZvzdBCCFEG/i1aFwuF3l5ecTGxgKQnJxMTk5Ok+VOnDhBYmIiAM8++yzl5eUUFhYyduxYFixYAEBERATh4eFSNEII0cX59RxNRUUFISEhmEwNq7Xb7RQXFzdZrqSkBLvd7ntst9spKipiwoQJvueysrJwOp0MGTKkTRkiIlo+CW63t2+a4q4gULMHam6Q7J0hUHNDYGfvKN2KJjs7m02bNjV6LioqqslNUs3dNKWUavS8UgrDXYNAZmdnk5qaSlpamq+0WuvGjVt4vU3vLLHbQyktrW7TZ3UVgZo9UHODZO8MgZobAjf7wypH3YomPj6e+Pj4Rs+5XC6ee+45PB4PRqOR0tJS+vTp0+S9ffv2paSkhIEDBwJQVlbmW27//v3s3buXvXv3Eh0drVd8IYQQD4lfz9GYzWbGjBlDVlYWAOnp6cTExDRZbtKkSWRkZABw9uxZrFYrkZGR5Obmsm/fPg4ePCglI4QQAUJTfh7C9Nq1a6xevZobN27Qv39/tm7dyiOPPMLBgwcpKSkhJSWF+vp6fv3rX3P+/HksFgtvvvkmw4cPJzExkfLyciIiInyf9+abb/LUU0+1ev1y6KzrCNTcINk7Q6DmhsDN/rAOnfm9aDqbFE3XEai5QbJ3hkDNDYGb/WEVjYyqKIQQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQlemzg7gbwaD1q7XurpAzR6ouUGyd4ZAzQ2Bnb2jNKWU6uwQQgghui85dCaEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCVz2qaAoLC5k3bx5xcXEsXbqUmpqaJss4nU5WrVpFfHw8M2fO5PLlywDU1NSQkpJCQkICCQkJHDt2LGCy31FQUMD06dP9kjczM5Np06YxdepUDhw40OT1/Px8kpOTiY2NZd26dbjdbqB126m39ma/Y/v27ezcudNfcX3am/vLL79k9uzZJCUl8eKLL3Lt2jV/R2939rNnz5KcnExCQgJLlizh5s2bAZH7josXLzJixAh/xW2kvdmPHDnCxIkTSUpKIikpiW3btj14ZaoHWbx4sTp69KhSSqldu3apt99+u8kyaWlpav369Uoppc6cOaPmzJmjlFJq69atavPmzUoppcrKytSECRNUaWmpn5J3LLtSSh05ckRNnDhRTZ48WfesRUVFavLkyaqiokLV1NSohIQEdenSpUbLTJ8+XX311VdKKaXWrFmjDhw4oJRq3XZ21exVVVVqzZo16umnn1bvvvtuwOSePHmyys/PV0op9Yc//EEtWbIkYLK/8MILvmXfeecd9dvf/jYgciulVG1trZo7d64aNmyY3zLf0ZHsGzduVJmZmW1aX4/Zo3G5XOTl5REbGwtAcnIyOTk5TZY7ceIEiYmJADz77LOUl5dTWFjI2LFjWbBgAQARERGEh4dTVlYWENmrq6s5fvw4W7du9UveU6dOMW7cOMLDwwkODiY2NrZR3mvXruFwOHjmmWcabU9rt7MrZgc4fvw4TzzxBAsXLvRr5o7kdjqdpKSk8OSTTwIQHR3N9evXAyI7QFZWFkOGDMHlclFcXExYWFhA5AbYvHkzL774ot/y3q0j2b/99luOHDlCQkICr732Wqv2IntM0VRUVBASEoLJ1DBgtd1up7i4uMlyJSUl2O1232O73U5RURETJkwgMjISaPjldjqdDBkyJCCyh4aGsnPnTvr37++XvPfm6NOnT6O8zeUsLi5u9Xbqqb3ZAWbMmMHixYsxGo3+C9xCrtbmtlgsJCUlAeD1etm1axcvvPCC/4I3k60tP3Oz2UxBQQGTJk3i9OnTfjs03NHcx48fx+FwEBcX57e8d+tIdrvdziuvvMLHH39M//792bhx4wPX1y2nCcjOzmbTpk2NnouKikLTGg/Tfe9jAKVUo+eVUhgMP/RxdnY2qamppKWl+f4gPkx6ZvcXr9fbJMfdj1t6/d7loPnt1FN7s3e2juZ2Op2sXr0at9vNyy+/7J/Qrcz2oNejo6M5deoUhw4d4pe//CWHDh3q0rlLS0vZvXs3+/bt80vO5nTkZ/673/3O9/yiRYuYMmXKA9fXLfdo4uPjOXnyZKN/77//PtXV1Xg8HgBKS0vp06dPk/f27duXkpIS3+OysjLfcvv372fLli3s3bvXd6ghULL7U79+/SgtLfU9vjfvva/fydm7d+9Wbaee2pu9s3Ukd01NDYsWLcLtdrN7927MZrP/gjeTrbXZ6+vryc3N9T2fmJhIQUGBf0I3k6u1uU+cOEFlZSXz5s3z7U0mJSVx69atLp+9urq6UUEqpVq1B98ti6Y5ZrOZMWPGkJWVBUB6ejoxMTFNlps0aRIZGRlAwxUtVquVyMhIcnNz2bdvHwcPHiQ6Ojqgsvvb+PHj+fzzzykvL6euro5PPvmkUd7HHnsMq9XKl19+CUBGRgYxMTGt3s6umL2zdST3qlWriIqKYvv27VgsloDJbjKZ2LBhA+fPnwcajgaMGjWqy+eeM2cOubm5ZGRk+P5/zcjIICQkpMtnDw4OJi0tjW+++QaADz74oFV7ND3qqrOrV6+q+fPnq/j4ePXSSy+pyspKpZRSH374odq+fbtSSimHw6F+9atfqWnTpqkZM2ao8+fPK6WUSkhIUBMmTFCJiYm+f+fOnQuI7HdcuXLFL1edKaXUxx9/rKZPn66mTp2q9uzZo5RSatGiRb6fWX5+vpo1a5aKjY1VK1asUPX19ffdTn9qb/Y73n33Xb9fddbe3BcuXFDDhg1T06ZN8/1eL1q0KCCyK6VUXl6emjlzpkpMTFS/+MUv1PXr1wMi990646ozpTr2M58xY4aKi4tTS5YsUVVVVQ9cl8ywKYQQQlc95tCZEEKIziFFI4QQQldSNEIIIXQlRSOEEEJXUjRCCCF0JUUjhBBCV1I0QtwlJyeHBQsWsGPHDtLT0zsthz/Wv2PHjlaNUyVER3XLsc6E6KiUlJRuu/6ioiJSU1M5efIkycnJuq1HiDukaESPt2PHDjIzMwkPDycqKgqA1atXM3ToUH7+85/z1FNPsXDhQk6dOkVtbS3Lli0jJyeHv/71r/Tp04ff//73BAcHc/nyZd566y0qKyvxeDwsWLCA2bNnc/r0abZt28bjjz/OpUuXcLvdbNiwgdGjR3P27Fk2b96M1+sF4OWXXyY2NrbR+s+ePcvbb79NXV0dZrOZV199lZiYGA4fPsyf/vQnDAYDf//737HZbGzZsoXBgwffd3s/+ugjxo4dy+DBg/0+UZjomeTQmejRcnNz+eSTT0hPT+fQoUPNDmzodDp59NFH+eijj5gxYwavv/4669atIysri1u3bnH8+HHcbjfLly9n5cqVHD58mA8++ID333+fr7/+GoBz587x0ksvkZ6eTnJysm9Wwp07d7Jw4UIOHz5MamoqX3zxRaN1V1RUsHz5ctatW0dmZiZbtmxh1apVXLlyBYC8vDzWr1/P0aNHGTlyJHv27HngNi9btoz58+d3ysjeomeS3zTRo33++edMmTLFNw/OrFmzml3uzmRsAwcOZNiwYfTt2xeDwcCAAQO4efMm3333Hd9//z1r164lKSmJ+fPn43A4uHjxIgCRkZH8+Mc/BuAnP/mJb08iPj6ejRs3snLlSi5cuMCKFSsarffcuXMMHDiQkSNHAjB06FBGjRrFmTNnABg+fDj9+vVr8rlCdCVy6Ez0eHcP99fSkOd3D53f3DD6Ho+H0NBQ32i80DC0emhoKF9//TU2m833/J25dwDmzp3L5MmT+eyzz/j000/ZtWtXo5kOPR5Pk/lulFK43W7MZnOLnytEVyJ7NKJHi4mJIScnh6qqKrxeb6OiaItBgwZhs9l8779+/To//elPfUPYt2Tu3Lnk5+eTnJzMG2+8QVVVVaN5QJ555hn+9re/ce7cOQAuXbpEXl4eY8eObVdOITqD7NGIHm3SpEkUFBQwa9YswsLCePLJJ6moqGjz51gsFt577z3eeust0tLScLvdpKSkMHr0aE6fPt3i+1577TVSU1PZvn07mqaxbNkyBgwY4Hu9d+/e7NixgzfeeAOHw4GmaWzatIlBgwbx1VdftWubhfA3mSZACCGErmSPRohuJjU1tcW9qDVr1jBu3Dg/JxI9nezRCCGE0JVcDCCEEEJXUjRCCCF0JUUjhBBCV1I0QgghdCVFI4QQQlf/D83vWhlF1Dy1AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "kmeans = KMeans(n_clusters=4)\n", "kmeans.fit(Trapnell_pca)\n", "_ = plt.scatter(Trapnell_pca[:,0],Trapnell_pca[:,1], c=kmeans.labels_, cmap='rainbow') \n", "_ = plt.xlabel(\"dimension_1\")\n", "_ = plt.ylabel(\"dimension_2\")\n", "_ = plt.title(\"cluster map after PCA\", size =20)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### (b) Which genes are good markers for each cluster? \n", "#### Using Logistic Regression to determine the good markers" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "Trapnell_df = pd.DataFrame(Trapnell);\n", "Tarpnell_sf = gl.SFrame(data=Trapnell)" ] }, { "cell_type": "code", "execution_count": 61, "metadata": {}, "outputs": [], "source": [ "Trapnell_label = kmeans.labels_" ] }, { "cell_type": "code", "execution_count": 62, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,\n", " intercept_scaling=1, max_iter=100, multi_class='warn',\n", " n_jobs=None, penalty='l2', random_state=None, solver='warn',\n", " tol=0.0001, verbose=0, warm_start=False)" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from sklearn.linear_model import LogisticRegression\n", "#Trapnell_sf = gl.SFrame(data=Trapnell_df)\n", "logreg = LogisticRegression(penalty=\"l2\")\n", "logreg.fit(Trapnell_df,Trapnell_label)" ] }, { "cell_type": "code", "execution_count": 63, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "cluster0 good maker gene's position: 13895 \t coef value: 0.294141443435\n", "cluster1 good maker gene's position: 94839 \t coef value: 0.168223717667\n", "cluster2 good maker gene's position: 13895 \t coef value: 0.626170928448\n", "cluster3 good maker gene's position: 668 \t coef value: 0.0395964722591\n" ] } ], "source": [ "cluster_0_largest_feature = np.argmax(logreg.coef_[0])\n", "cluster_1_largest_feature = np.argmax(logreg.coef_[1])\n", "cluster_2_largest_feature = np.argmax(logreg.coef_[2])\n", "cluster_3_largest_feature = np.argmax(logreg.coef_[3])\n", "print \"cluster0 good maker gene's position:\",cluster_0_largest_feature, \"\\t coef value:\", max(logreg.coef_[0])\n", "print \"cluster1 good maker gene's position:\",cluster_1_largest_feature, \"\\t coef value:\", max(logreg.coef_[1])\n", "print \"cluster2 good maker gene's position:\",cluster_2_largest_feature, \"\\t coef value:\", max(logreg.coef_[2])\n", "print \"cluster3 good maker gene's position:\",cluster_3_largest_feature, \"\\t coef value:\", max(logreg.coef_[3])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### c) Map the cells to a differentiation tree. " ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAekAAAFHCAYAAACF2LQJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XtwlNX9x/HP5rKBhFy5BsLGZEEoF4lJSEhoO+Vi1dppf0PV1tpWmbHt1GmndkQQtbajU622DmO9V62lVWsr+uvoWLwAtrYkITFchIFiLoQkEANCEnPfJPv8/kCeH5Hllmxynt19v2YYO5Nk86Gd6cfzfc5zjsuyLEsAAMBxokwHAAAAgVHSAAA4FCUNAIBDUdIAADgUJQ0AgENR0gAAOBQlDQCAQ1HSAAA4FCUNAIBDUdIAADgUJQ0AgENR0gAAOBQlDQCAQ1HSAAA4FCUNAIBDUdIAADgUJQ0AgENR0gAAOBQlDQCAQ1HSAAA4VIzpALhwvtY2Hdm8RV0H69Xf2amYhATFZ3o0eflSxSYnm44HAAgSl2VZlukQOD/tVdVqfPkVtWzfIblcsnw++2sut1uyLKXm5Srj6hVKnDnDYFIAQDBQ0iGiaeNbqntuvfw+n3S2/8lcLkW53bpo5Q1Kv/Ly0QsIAAg6xt0h4ERB/1H+Xt+5v9my5O/tVd1z6yWJogaAEMbGMYdrr6o+sYI+n4I+xcmibq+qHqFkAICRRkk7XOPLr5wYcQ+B3+dT44ZXg5wIADBaKGkH87W2ndgkNtRtA5allsrt6mtrC24wAMCooKQd7MjmLZLLNbwPcbnUvPnd4AQCAIwqStrBug7WD3rNaigsn09dB+uDlAgAMJooaQfr7+wM0ud0BOVzAACji5J2sJiEhKB8zp6qam3cuFEtLS1B+TwAwOjgPWkHi8/0yOV2D2/kHRMjTZygF154Qbfccos8Ho+Ki4tVXFyswsJCpaSkBC8wACCoOHHMwXytbXr/ph/K6usb8me4YmO18NmnFJucrL6+Pu3atUulpaUqKSlRZWWlsrKyBpV2UlJSEP8GAIDhoKQdbt99D+h4ecXQXsNyuZRWWKDPrV0d8Ms+n0+7du1SSUmJSkpKtGPHDs2YMUNFRUUqLi5WQUGBEhMTh/k3AAAMFSXtcO1V1dpz593y9/Ze8M9GxcVp3q/uOe/LNnp7e7Vz5067tHfu3KlZs2apuLhYRUVFKigoUEKQnpMDAM6Nkg4B9uUaF1DUUXFxw75ko6enRzt27LBL+4MPPtDs2bPt8fjChQsVHx8/5M8HAJwdJR0inHALVnd3tyorK+1n2nv27NHcuXPt8Xh+fr7Gjh0b1N8JAJGMkg4h7VXVatzwqloqtzviPumuri69//77KikpUWlpqfbu3av58+fb4/G8vDyNGTNmxHMAQLiipENQX1ubmje/q66D9erv7FBMwjjFZ3o0edkSxSYnG8vV2dlpl/bWrVu1f/9+LViwwC7t3NxcxcXFGcsHAKGGksaI6ejoUHl5uT0er6qqUk5OjoqKirR48WLl5OTI7XabjgkAjkVJY9S0t7dr27Zt9ni8pqZGubm59jPtnJwcxcbGmo4JAI5BScOYtra2QaVdV1envLw8ezx+ySWXUNoAIholDcdoaWlReXm5tm7dqpKSEjU2Nio/P99+5WvevHmKieEkWwCRg5KGYx0/flxlZWX2M+3Dhw9r4cKFdmnPnTtX0dHRpmMCwIihpBEyPv74Y5WVldnj8ebmZhUUFNilPWfOHEVFcbEbgPBBSSNkHT161C7skpISHTt2TIWFhXZpz549m9IGENIoaYSN5uZmu7BLSkrU2tqqoqIie/f4xRdfTGkDCCmUNMJWU1OTXdqlpaVqb2/XokWL7JX2zJkz5XK5TMcEgDOipBExDh06NKi0u7q67FV2cXGxvF4vpQ3AUShpRKyGhoZBz7T7+/sHjcezsrIobQBGUdKAJMuy7NI++Z62JHuVXVRUpMzMTEobwKiipIEALMtSXV3doI1o0dHR9rnjRUVF8ng8pmMCCHOUNHAeLMtSbW3toPF4XFycvcpevHixpk2bZjomgDBDSQNDYFmWampq7NF4aWmpEhIS7NIuLi7W1KlTTccEEOIoaSAILMvShx9+qNLSUm3dulVlZWVKSkoa9Ex7ypQppmMCCDGUNDAC/H6/9u/fP+iVr7S0tEHPtCdNmmQ6JgCHo6SBUeD3+7Vv3z67sLdt26YJEyYMGo9PmDDBdEwADkNJAwYMDAxo37599jPt8vJypaen24VdVFSktLQ00zEBGEZJAw4wMDCgPXv22M+0KyoqlJGRYZf2okWLlJqaajomgFFGSQMO1N/fr927d9vj8YqKCnk8HnsjWmFhoVJSUkzHBDDCKGkgBPT19WnXrl32RrTKykplZWUNKu2kpCTTMQEEGSUNhCCfz6ddu3bZp6Ht2LFDM2bMsMfjBQUFSkxMNB0TwDBR0kAY6O3t1c6dO+3S3rlzp2bNmmVvQisoKFBCQoLpmAAuECUNhKGenh7t2LHDfqa9a9cuzZ49W8XFxVq8eLHy8/MVHx9vOiaAc6CkgQjQ3d2tyspK+5n2nj17NHfuXHs8np+fr7Fjx5qOCeAzKGkgAnV3d6uiosJeae/du1fz58+3x+N5eXkaM2aM6ZhAxKOkAaizs1Pvv/++fZ/2/v37tWDBAru0c3NzFRcXZzomEHEoaQCn6ejoUHl5uT0er6qqUk5Ojn32eE5Ojtxut+mYQNijpAGcU3t7u7Zt22aPx2tqapSbm2s/087JyVFsbKzpmEDYoaQBXLC2trZBpV1XV6e8vDx7PH7JJZdQ2kAQUNIAhq2lpUXl5eX2hSGNjY3Kz8+3T0SbN2+eYmJiTMcEQg4lDSDojh8/rrKyMvuZ9uHDh7Vw4UK7tOfOnavo6GjTMQHHo6QBjLhjx47ZhV1aWqrm5mYVFBTYpT1nzhxFRUWZjgk4DiUNYNQdPXrULu2SkhIdO3ZMhYWFdmnPnj2b0gZESQNwgObm5kGl3draqqKiInv3+MUXX0xpIyJR0gAcp6mpadB4vL29XYsWLbJX2jNnzpTL5TIdExhxlDQAxzt06NCg0u7q6rJX2cXFxfJ6vZQ2whIlDSDkNDQ02IVdUlKi/v7+QePxrKwsShthgZIGENIsy7JL++Qfy7LsVXZRUZEyMzMpbYQkShpAWLEsS3V1dYM2okVHR9vnjhcVFcnj8ZiOCZwXShpAWLMsS7W1tYPG43FxcfYqe/HixZo2bZrpmEBAlDSAiGJZlmpqauwjTEtLS5WQkGCXdnFxsaZOnWo6JiCJkgYQ4SzL0ocffqjS0lJt3bpVZWVlSkpKGvRMe8qUKaZjIkJR0gBwCr/fr/379w965SstLW3QM+1JkyaZjokIQUkDwFn4/X7t27fPLuxt27ZpwoQJg8bjEyZMMB0TYYqSBoALMDAwoH379tnPtMvLy5Wenm4XdlFRkdLS0kzHRJigpAFgGAYGBrRnzx77mXZFRYUyMjLs0l60aJFSU1NNx0SIoqQBIIj6+/u1e/duezxeUVEhj8djb0QrLCxUSkqK6ZgIEZQ0AIygvr4+ffDBB/bBKpWVlcrKyhpU2klJSaZjwqEoaQAYRT6fT7t27bJLe8eOHZoxY4Y9Hi8oKFBiYqLpmHAIShoADOrt7dXOnTvt0t65c6dmzZplb0IrKChQQkKC6ZgwhJIGAAfp6enRjh077Gfau3bt0uzZs1VcXKzFixcrPz9f8fHxpmNilFDSAOBg3d3dqqystA9X2bNnj+bOnWuPx/Pz8zV27FjTMTFCKGkACCHd3d2qqKiwS3vv3r2aP3++PR7Py8vTmDFjTMdEkFDSABDCOjs79f7779vPtP/73/9qwYIFdmnn5uYqLi7OdEwMESUNAGGko6ND5eXl9kq7qqpKOTk59tnjOTk5crvdpmPiPFHSABDG2tvbtW3bNnsjWk1NjXJzc+1n2jk5OYqNjTUdE2dASQNABGlraxtU2nV1dcrLy7PH45dccgml7SCUNABEsJaWFpWXl9sXhjQ2Nio/P98+EW3evHmKiYkxHTNiUdIAANvx48dVVlZmP9M+fPiwFi5caJf23LlzFR0dbTpmxKCkAQBndOzYMbuwS0tL1dzcrIKCAru058yZo6ioKNMxwxYlDQA4b0ePHh1U2h9//LEKCwvt0p49ezalHUSUNABgyJqbm+3SLikpUWtrq4qKiuzd4xdffDGlPQyUNAAgaJqamgattNvb27Vo0SJ7pT1z5ky5XC7TMUMGJQ0AGDGHDh0aVNpdXV32Kru4uFher5fSPgtKGgAwahoaGuzCLikpUX9//6DxeFZW1oiXtq+1TUc2b1HXwXr1d3YqJiFB8ZkeTV6+VLHJySP6uy8UJQ0AMMKyLLu0T/6xLMteZRcVFSkzMzNopd1eVa3Gl19Ry/Ydkssly+ezv+ZyuyXLUmperjKuXqHEmTOC8juHi5IGADiCZVmqq6sbtBEtOjraPne8qKhIHo9nSJ/dtPEt1T23Xn6fTzpb7blcinK7ddHKG5R+5eVD/JsEDyUNAHAky7JUW1s7qLTj4uLsVfbixYs1bdq0c37OiYL+o/y9vnN+70lRcXGOKGpKGgAQEizLUk1NjbZu3arS0lKVlpYqPj7eLu3i4mJNnTp10M+0V1Vrz513y9/be8G/LyouTvN+dY/R0TclDQAISZZlqaqqSiUlJdq6davKysqUlJQ06Jl2yx/W63h5xdlH3GficimtsECfW7s6+OHPNwIlDQAIB36/X/v377fH43vKy/Xri2Yodhgbz1yxsVr47FPGdn1T0gCAsNSw4VXV/+WvUn//kD/D5XbLc903lbHif4KY7PxxVhsAICx11zcMq6AlyfL51HWwPkiJLhwlDQAIS/2dnUH6nI6gfM5QUNIAgLAUk5AQpM8ZF5TPGQpKGgAQluIzPSdOEhsGl9ut+MyhHaASDJQ0ACAsTVq2dGivXp3KsjR52ZLgBBoCShoAEJbcKclKzb1UGuorWC6XUvNyjV66QUkDAMJWxjXfUNQQR95Rbrcyrl4R5EQXmMHobwcAYAQlzpyhi1becMHPpk+e3W36NixKGgAQ1tKvvFxlY+M0EOU69+jb5XLM5RqSFGM6AAAAI+mdd97RhuoP9b3fP62jr7+hlsrt3CcNAIBp7e3tWrJkidatW6cvfOELkqS+tjY1b35XXQfr1d/ZoZiEcYrP9GjysiVGN4kFQkkDAMLW2rVr5fP59NBDD5mOMiSMuwEAYam8vFxvvfWWtmzZYjrKkLFxDAAQdnp6erRq1Srde++9SklJMR1nyChpAEDYefjhh3XxxRfrqquuMh1lWBh3AwDCyt69e/X888/rnXfeMR1l2FhJAwDCRn9/v1atWqW1a9dqypQppuMMGyUNAAgbzzzzjOLj43XdddeZjhIUjLsBAGHh4MGDevTRR/X666/LNdRLNRyGlTQAIORZlqU1a9bo5ptvVlZWluk4QUNJAwBC3t/+9je1trbqBz/4gekoQcWJYwCAkHbkyBEtX75cL774oubNm2c6TlCxkgYAhLS77rpL1113XdgVtMTGMQBACHvzzTe1d+9ePfzww6ajjAjG3QCAkNTW1qalS5fq0UcfVVFRkek4I4KSBgCEpNWrV0uSHnzwQcNJRg7jbgBAyCktLdXmzZv17rvvmo4yotg4BgAIKd3d3Vq1apXuu+8+JSUlmY4zoihpAEBIWbdunebNm6fLL7/cdJQRx7gbABAydu/erZdeekmbN282HWVUsJIGAISEkzdc3XnnnZo4caLpOKOCkgYAhISnnnpKqampuvbaa01HGTWMuwEAjldbW6snnnhCb7zxRtjccHU+WEkDABzNsiytXr1aP/nJT5SZmWk6zqiipAEAjvbiiy+qu7tbN910k+koo44TxwAAjvXRRx/psssu01//+lfNmTPHdJxRx0oaAOBIlmXpzjvv1He/+92ILGiJjWMAAId64403VF1drccff9x0FGMYdwMAHKelpUXLli3Tk08+qYKCAtNxjKGkAQCOc+uttyouLk733Xef6ShGMe4GADjKv//9b7333nvasmWL6SjGsXEMAOAY3d3dWrNmje6//34lJiaajmMcJQ0AcIzf/OY3uvTSS7V8+XLTURyBcTcAwBF27typV155hTH3KVhJAwCM6+vr06pVq3T33Xdr/PjxpuM4BiUNADDuiSee0OTJk7VixQrTURyFcTcAwKjq6mr9/ve/15tvvhlRN1ydD1bSAABj/H6/brvtNv3sZz9TRkaG6TiOQ0kDAIz585//rP7+ft14442mozgSJ44BAIw4fPiwLr/8cm3YsEGzZs0yHceRWEkDAEadZVlau3atVq5cSUGfBRvHAACj7rXXXlNDQ4Oefvpp01EcjXE3AGBUHT9+XMuWLdMzzzyjvLw803EcjZIGAIyqn/70p0pOTtY999xjOorjMe4GAIyaf/7znyorK+Poz/PExjEAwKjo7OzU7bffrgceeEAJCQmm44QEShoAMCoeeOABFRQU6Etf+pLpKCGDcTcAYMRVVlbq9ddf1+bNm01HCSmspAEAI8rn8+m2227TL3/5S6WlpZmOE1IoaQDAiHrsscc0ffp0fe1rXzMdJeQw7gYAjJgPP/xQf/jDH7jhaohYSQMARsTAwIBWrVqlW2+9VdOmTTMdJyRR0gCAEbF+/XpFR0fre9/7nukoIYsTxwAAQdfY2KgrrrhCf//73zVjxgzTcUIWK2kAQFBZlqXbb79d3//+9ynoYaKkAQBB9eqrr+qjjz7SzTffbDpKyGPcDQAImmPHjmnp0qVav369cnJyTMcJeaykAQBB84tf/EIrVqygoIOE96QBAEGxadMmbd++XZs2bTIdJWww7gYADFtHR4eWLl2q3/72t/riF79oOk7YYNwNABi2+++/X5///Ocp6CBj3A0AGJaKigpt3LhRW7ZsMR0l7LCSBgAMWU9Pj1atWqV77rlHKSkppuOEHUoaADBkjzzyiGbMmKGrrrrKdJSwxLgbADAk+/bt05/+9Ce9/fbb3HA1QtjdDQC4YAMDA/r617+ub33rW/rOd75jOk7YYtwNALhgzz77rMaMGaNvf/vbpqOENVbSAIALUl9fr6985St67bXXlJ2dbTpOWGMlDQA4b5ZlafXq1frRj35EQY8CShoAcN5efvlltbS06Ic//KHpKBGBcTcA4LwcPXpUy5cv1/PPP6/58+ebjhMRWEkDAM7Lz3/+c1177bUU9CjiPWkAwDm9/fbb2r17t9atW2c6SkRh3A0AOKtPPvlES5cu1e9+9zsVFxebjhNRKGkAwFmtWbNGlmXpwQcfNB0l4jDuBgCcUWlpqTZt2sQNV4awcQwAEFBPT49uu+02/epXv1JycrLpOBGJkgYABLRu3TrNmTNHV1xxhekoEYtxNwDgNHv27NFf/vIXbdq0yXSUiMZKGgAwSH9/v1atWqU77rhDkyZNMh0nolHSAIBBnn76aSUnJ+ub3/ym6SgRj3E3AMB24MABPfbYY3rjjTfkcrlMx4l4rKQBAJL+/4arH//4x8rMzDQdB6KkAQCfeumll9TZ2ambbrrJdBR8ihPHAABqbm7WZZddppdeeklz5swxHQefYiUNANBdd92l66+/noJ2GDaOAUCE+8c//qH9+/frkUceMR0Fn8G4GwAiWGtrq5YtW6bHH39chYWFpuPgMyhpAIhgq1atUmxsrO6//37TURAA424AiFD/+c9/9K9//YsbrhyMjWMAEIG6u7u1Zs0a3XfffUpMTDQdB2dASQNABHrooYe0YMECXXbZZaaj4CwYdwNAhPnggw/08ssva/Pmzaaj4BxYSQNABOnr69Ott96qu+66SxMmTDAdB+dASQNABHnyySc1ceJEXX311aaj4Dww7gaACFFdXa2nnnpKGzdu5IarEMF70gAQAfx+v6655hpdeeWVXKARQhh3A0AEeOGFF9Tb26uVK1eajoILwEoaAMJcU1OTvvzlL2vDhg2aNWuW6Ti4AKykASCMWZalO+64QzfeeCMFHYLYOAYAYez111/XgQMH9OSTT5qOgiFg3A0AYer48eNatmyZnn76aeXn55uOgyGgpAEgTN1yyy1KTEzUvffeazoKhohxNwCEoffee0+lpaXccBXi2DgGAGGmq6tLa9as0a9//WslJCSYjoNhYCUNACHI19qmI5u3qOtgvfo7OxWTkKD4TI8mL1+qB9etU35+vpYsWWI6JoaJZ9IAEELaq6rV+PIratm+Q3K5ZPl89tdcbrcs/4B2trXpf359n6bl5xlMimCgpAEgRDRtfEt1z62X3+eTzvJ/3Zak6Lg4XbTyBqVfefnoBUTQMe4GgBBwoqD/KH+v75zf65Lk7+1V3XPrJYmiDmFsHAMAh2uvqj6xgj6Pgj7VyaJur6oeoWQYaZQ0ADhc48uvnBhxD4Hf51PjhleDnAijhZIGAAfztbad2CQ21O1DlqWWyu3qa2sLbjCMCkoaABzsyOYtkss1vA9xudS8+d3gBMKooqQBwMG6DtYPes1qKCyfT10H64OUCKOJ3d0A4DB+v1+HDh1STU2Nev77XyUG4TP7OzuC8CkYbZQ0ABjS1tammpqaQX8OHDigAwcOKCUlRV6vV98YEx+Uko5JGBeET8Foo6QBYAT19fWpvr5eNTU1qq2tHVTI3d3dys7Oltfrldfr1Ve/+lV5vV5lZWXZZ243vvK/qn/pb8MaebvcbsVneoL1V8Io4sQxABgmy7J07Nix01bFNTU1OnTokKZMmWKX78lC9nq9mjJlilzn2BTma23T+zf9UFZf35DzuWJjtfDZpxSbnDzkz4AZrKQB4Dz19PTowIEDdgGfXBnX1tbK5XINWhVfe+218nq9yszM1JgxY4b8O90pyUrNvVTHyyuG9hqWy6XUvFwKOkRR0gBwCsuy1NTUdFoR19TU6MiRI5o+fbpdxEVFRbr++uvl9XqVlpZ2zlXxUGVc8w217twlf2/vBf9slNutjKtXjEAqjAbG3QAiUkdHx2nPiE9u3Bo3bpy8Xq+ys7MHrY49Ho9iYsysbezLNS6gqKO4ZCPkUdIAwtbAwIAaGhoCjqdbW1tPe0Z8spiTkpJMRw/ofG/BksulKLebgg4DlDSAkNfS0hJwPF1fX6/x48efVsJer1dTp05VVFTonefUXlWtxg2vqqVye8D7pGVZSs3LVcbVK5Q4c4bBpAgGShpASPD5fDp48GDAHdT9/f0Bx9PZ2dkaO3as6egjoq+tTc2b31XXwXr1d3YoJmGc4jM9mrxsCZvEwgglDcAxLMvSkSNHAr5T3NTUpPT09NPG016vVxMnThyxTVuASZQ0gFHX3d09qIRPfVYcGxsbcDydmZkpt9ttOjowqihpACPC7/fr8OHDp42ma2trdezYMXk8ntOK2Ov1KjU11XR0wDEoaQDD8sknnwQcT9fV1Sk5OXlQAZ/8k5GRoejoaNPRAcejpAGcU39/v33+9GcLubOz87TNWif/OW4clzoAw0FJA5B0YtPW8ePHA46nGxoaNGnSpIDj6fT0dDZtASOEkgYiTE9PT8BXmWpra2VZVsDxdGZmZti+ygQ4GSUNhCHLsvTRRx+dVsI1NTVqbm5WRkbGaStir9er8ePHsyoGHISSBkJYZ2fnaRu2amtrVVtbq7FjxwZ8lcnj8Sg2NtZ0dADngZIGHG5gYECHDh0KeNLWyfOnP7sizs7OVjKnTgEhj5IGHKK1tTXgePrgwYNKS0sLOJ6eNm1aSJ4/DeD8UNLAKOrr67M3bX12TN3T0xNwPJ2dna34+HjT0QEYQEkDQWZZlj7++OOA4+nDhw8rPT094HvFkydPZtMWgEEoaWCIuru7deDAgYCbtqKiogJeBJGZmam4uDjT0QGECEoaOAu/36+mpqaA4+mjR48OOn/61FVxWlqa6egAwgAlDUhqb28/rYRramp04MABJSUlBRxPT58+XTExMaajAwhjlDQiRn9/vxoaGgLuoG5vb1dWVtZp4+ns7GwlJiaajg4gQlHSCDsnz5/+7Mq4oaFBEyZMCDieTk9P51UmAI5DSSMk9fb2Bjx/uqamRn6//7RXmLxer7Kysjh/GkBIoaThWJZlqbm5OeDu6aamJk2dOjXgDuoJEybwKhOAsEBJw7iuri57NH3qiLq2tlZxcXEBnxN7PB653W7T0QFgRFHSGBV+v/+M50+3tLTooosuOm087fV6lZKSYjo6ABhDSSOo2traAo6nDxw4oJSUlIDj6WnTpik6Otp0dABwHEoaF6yvr0/19fUBd1B3d3efdgnEyU1bCQkJpqMDQEihpBGQZVk6duxYwPH0oUOHNGXKFGVnZ59WyFOmTGHTFgAECSUd4Xp6euzzpz+7acvlcgVcFWdmZmrMmDGmowNA2KOkA/C1tunI5i3qOliv/s5OxSQkKD7To8nLlyo2Odl0vAtmWZZ9/vRny/jIkSOaPn36Gc+fZlUMAOZQ0qdor6pW48uvqGX7DsnlkuXz2V9zud2SZSk1L1cZV69Q4swZBpMG1tHRccbzp8eNGxdw97TH4+H8aQBwKEr6U00b31Ldc+vl9/mks/1X4nIpyu3WRStvUPqVl49ewE8NDAzY509/djzd2tp6xvOnk5KSRj0rAGB4KGmdLOg/yt/rO/c3fyoqLm5Ei7qlpSXgeLq+vl7jx48POJ6eOnUq508DQBiJ+JJur6rWnjvvlr+394J/NiouTvN+dc+QR98+n++M50/39fWdVsIn/zPnTwNAZIj4kt533wM6Xl5x9hH3mbhcSiss0OfWrj7jt1iWpSNHjgR8p7ipqUnp6ekBD/iYOHEim7YAIMJF9I4hX2vbiU1iQ/33FMtSS+V29bW1qd/tPu0Z8cl/xsbGDirgwsJC+1Umzp8GAJxJRK+kG1/5X9W/9LdBu7gvVJ+kf3zSqlfq6+TxeAKOqFNTU4MXGgAQMSJ6Jd11sH5YBS1JsZKuW7pM99+xhvOnAQBBFdFbgfs7O4PyOXEuUdAAgKCL6JKOCdKFDzEJ44LyOQAAnCqiSzo+03PiJLFhcLndis/0BCkRAAD/L6JLetKypUPf2X2SZWnysiXdQdGNAAABX0lEQVTBCQQAwCkiuqTdKclKzb1UGur7yC6XUvNyQ/LSDQCA80V0SUtSxjXfUNQQR95Rbrcyrl4R5EQAAJwQ8SWdOHOGLlp5g6Li4i7o506e3e3E27AAAOEhot+TPunkJRmhcAsWACByRPSJY5/VXlWtxg2vqqVye0jeJw0ACC+UdAB9bW1q3vyuug7Wq7+zQzEJ4xSf6dHkZUvYJAYAGDWUNAAADhXxG8cAAHAqShoAAIeipAEAcChKGgAAh6KkAQBwKEoaAACHoqQBAHAoShoAAIeipAEAcChKGgAAh6KkAQBwKEoaAACHoqQBAHAoShoAAIeipAEAcChKGgAAh6KkAQBwKEoaAACHoqQBAHAoShoAAIeipAEAcChKGgAAh6KkAQBwKEoaAACHoqQBAHAoShoAAIeipAEAcChKGgAAh6KkAQBwKEoaAACHoqQBAHAoShoAAIeipAEAcKj/A6vuxus5y1iHAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import networkx as nx\n", "centers = kmeans.cluster_centers_\n", "distance = scipy.spatial.distance.cdist(centers,centers)\n", "G = nx.from_numpy_matrix(distance)\n", "T = nx.minimum_spanning_tree(G)\n", "nx.draw(T)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.1" } }, "nbformat": 4, "nbformat_minor": 2 }