getting_started_notebooks/intro_tutorials/10_Introduction_to_Clustering.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"introduction\"></a>\n",
    "## Introduction to Clustering\n",
    "#### By Paul Hendricks\n",
    "-------\n",
    "\n",
    "Clustering is an important technique for helping data scientists partition data, especially when that data doesn't have labels or annotations associated with it. Since these data often don't have labels, clustering is often described as an unsupervised learning technique. While there are many different algorithms that partition data into unique clusters, we will show in this notebook how in certain cases the DBSCAN algorithm can do a better job of clustering than traditional algorithms such as K-Means or Agglomerative Clustering. \n",
    "\n",
    "We will also show how to cluster data with K-Means and DBSCAN in NVIDIA RAPIDS – an open-source data analytics and machine learning acceleration platform that leverages GPUs to accelerate computations. We will see that porting this example from CPU to GPU is trivial and that we can experience massive performance gains by doing so.\n",
    "\n",
    "In this notebook, we will describe K-Means and DBSCAN - several popular clustering algorithms - and show how to use the GPU accelerated implementation of these algorithms in RAPIDS.\n",
    "\n",
    "**Table of Contents**\n",
    "\n",
    "* [Introduction to Clustering](#introduction)\n",
    "* [Setup](#setup)\n",
    "* [Generating Data](#generate)\n",
    "* [Overview of K-Means and Agglomerative Clustering](#overview)\n",
    "* [Clustering using DBSCAN](#dbscan)\n",
    "* [Accelerating Clustering with RAPIDS](#accelerating)\n",
    "* [Clustering using GPU Accelerated K-Means](#gpukmeans)\n",
    "* [Benchmarking: Comparing GPU and CPU](#benchmarking)\n",
    "* [K-Means](#benchmarkingkmeans)\n",
    "  * [GPU](#benchmarkingkmeansgpu)\n",
    "  * [CPU](#benchmarkingkmeanscpu)\n",
    "* [DBSCAN](#benchmarkingdbscan)\n",
    "  * [GPU](#benchmarkingdbscangpu)\n",
    "  * [CPU](#benchmarkingdbscancpu)\n",
    "* [Conclusion](#conclusion)\n",
    "\n",
    "Before going any further, let's make sure we have access to `matplotlib`, a popular Python library for visualizing data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Matplotlib Version: 3.1.2\n"
     ]
    }
   ],
   "source": [
    "import os\n",
    "\n",
    "try:\n",
    "    import matplotlib; print('Matplotlib Version:', matplotlib.__version__)\n",
    "except ModuleNotFoundError:\n",
    "    os.system('conda install -y matplotlib')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"setup\"></a>\n",
    "## Setup\n",
    "\n",
    "This notebook was tested using the following Docker containers:\n",
    "\n",
    "* `rapidsai/rapidsai-nightly:0.8-cuda10.0-devel-ubuntu18.04-gcc7-py3.7` from [DockerHub - rapidsai/rapidsai-nightly](https://hub.docker.com/r/rapidsai/rapidsai-nightly)\n",
    "\n",
    "This notebook was run on the NVIDIA Tesla V100 GPU. Please be aware that your system may be different and you may need to modify the code or install packages to run the below examples. \n",
    "\n",
    "If you think you have found a bug or an error, please file an issue here: https://github.com/rapidsai/notebooks/issues\n",
    "\n",
    "Before we begin, let's check out our hardware setup by running the `nvidia-smi` command."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wed Jan 22 22:15:57 2020       \n",
      "+-----------------------------------------------------------------------------+\n",
      "| NVIDIA-SMI 418.39       Driver Version: 418.39       CUDA Version: 10.1     |\n",
      "|-------------------------------+----------------------+----------------------+\n",
      "| GPU  Name        Persistence-M| Bus-Id        Disp.A | Volatile Uncorr. ECC |\n",
      "| Fan  Temp  Perf  Pwr:Usage/Cap|         Memory-Usage | GPU-Util  Compute M. |\n",
      "|===============================+======================+======================|\n",
      "|   0  Quadro GV100        Off  | 00000000:15:00.0 Off |                  Off |\n",
      "| 29%   40C    P2    26W / 250W |  16565MiB / 32478MiB |      0%      Default |\n",
      "+-------------------------------+----------------------+----------------------+\n",
      "|   1  Quadro GV100        Off  | 00000000:2D:00.0  On |                  Off |\n",
      "| 33%   46C    P0    29W / 250W |    284MiB / 32470MiB |     22%      Default |\n",
      "+-------------------------------+----------------------+----------------------+\n",
      "                                                                               \n",
      "+-----------------------------------------------------------------------------+\n",
      "| Processes:                                                       GPU Memory |\n",
      "|  GPU       PID   Type   Process name                             Usage      |\n",
      "|=============================================================================|\n",
      "+-----------------------------------------------------------------------------+\n"
     ]
    }
   ],
   "source": [
    "!nvidia-smi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's see what CUDA version we have:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/bin/sh: 1: nvcc: not found\n"
     ]
    }
   ],
   "source": [
    "!nvcc --version"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's load some helper functions from `matplotlib` and configure the Jupyter Notebook for visualization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from matplotlib.colors import ListedColormap\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"generate\"></a>\n",
    "## Generating Data\n",
    "\n",
    "We'll generate some fake data using the `make_moons` function from the `sklearn.datasets` module. This function generates data points from two equations, each describing a half circle with a unique center. Since each data point is generated by one of these two equations, the cluster each data point belongs to is clear. The ideal clustering algorithm will identify two clusters and associate each data point with the equation that generated it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Scikit-Learn Version: 0.22\n",
      "(100, 2)\n"
     ]
    }
   ],
   "source": [
    "import sklearn; print('Scikit-Learn Version:', sklearn.__version__)\n",
    "from sklearn.datasets import make_moons\n",
    "\n",
    "\n",
    "X, y = make_moons(n_samples=int(1e2), noise=0.05, random_state=0)\n",
    "print(X.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's visualize our data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deZxU1Zn/8c8jjYHY3dCNCEKDaCKucUVcEwjGuCSBTMbMGBOJAUMIiaJAgOj8jMkkDrgwkbiFAHGbGP1lmZC4L1F/Go1iwAVZNSodjCA00k1E6fj8/qjbbVFdVV3LrapbVd/361Wv7qq6fe+5fbvrueec55xj7o6IiEjU7FbqAoiIiCSjACUiIpGkACUiIpGkACUiIpGkACUiIpGkACUiIpGkACVFY2ZtZrZfmvdvNLP/k+cxRptZcz77iCozu8fMvlrqcnQwsxVmNjrN+4+Y2XkFLsMwM3Mzqynmz0pxKEBVOTN71czeCYJHx+PaQhzL3Wvd/ZU070929/8sxLE7BB9I24Pz3GxmD5nZv2fx83kFwLgPxY7f9atmNjuTn3X309395gyPU/Dg4O6HuPsjwfEuM7Pbct1XJd9YSO505yAAn3P3B0tZADPr4e7/LNLhDnf3dWa2J3A6cK2ZHeju3y/S8QH6unu7mR0PPGRmy9393iIeXyTyVIOSlMzsXDN7wsz+28y2mtkrZnZC8Pp6M9sY3+RkZjcFzXQPmFmrmT1qZvvEve9m9tG4bW8ws7vNbDvwyeC1H8ZtP87MlpvZNjN72cxOC17/mpmtDI7xipl9I5fzc/e33P1W4JvAd82sX7r9m9kewD3AoLga0CAzG2lmTwa/ozfM7Foz2z3DMjwJrAAODY5xgpk9Y2ZvB19PiPt9dNaKgmvwuJldZWYtZvZXMzs9eO9HwMeJBd62oDwWXMeNwb6fN7NDE8tjZp80sxfinj9oZk/HPX/czD4ffP+qmX0quC4XA/8eHO+5uF3uE/wNtZrZ/cFNQVbM7DNmtiz4O1hvZpcl2WyCmW0Ifv/T4352NzObHfz9bDazO82sMcVxzg2ud2vw+/xytmWVkLm7HlX8AF4FPpXivXOBduBrQA/gh8DrwHXAh4BPA61AbbD9TcHzTwTvXwM8Hrc/Bz4at+3bwInEbpR6Ba/9MHh/ZPD+KcH7g4EDg/c+A3wEMGAU8A/gqOC90UBzmvPtLEPcaz2D8zw9l/0DRwPHEWuRGAasBC5McfxhQRlqgv2fGOz/ZKARaAHOCd7/UvC8X/CzjwDnxV2bncDXg2vzTWADYInbBs9PBZ4F+gbHPQjYO0n5egHvAHsGZfh7sN86oHfwXkd5XiX42wEuA25L2NcjwMvA8OBnHwHmpPi9pLxuwXsfC/4ODgPeBD6f8Pu8Hdgj2G5TXLkuBJ4Cmoj9Tf4UuD3JtdgD2AYcELy3N3BIqf8/q/2hGpQA/G9w99/x+Hrce3919597rPntDmAI8AN3f9fd7wfeAz4at/1d7v6Yu78LXAIcb2ZDUhz3d+7+hLu/7+47Et6bCCx29weC9//m7qsA3P0ud3/ZYx4F7idWY8iJu+8E3iIWILLev7s/6+5PuXu7u79K7ENwVDeHfQvYAiwEZrv7Q8QC41p3vzXY1+3AKuBzKfbxmrv/LLg2NxP7UB2QYtudxILMgcSC2Ep3fyPJuewAlhK7yRgBPA88TiyQHheUb3M35xbv5+6+xt3fAe4EjsjiZzvK9Ii7vxD8HTxPLBgl/n6/7+7b3f0F4OfEgjvAN4BL3L05+Ju8DDjTkidGvA8cama93f0Nd1+RbVklXOqDEojdjabqg3oz7vt3ANw98bXauOfrO75x9zYz2wIMin892bZJDAHuTvZG0JT1PWJ35rsBHwZeSLZtJsysJ9CfWMDIev9mNhyYR+wD/cPE/q+e7eawe7p7e8Jrg4DXEl57jVjtMZm/d3zj7v8wM9j1WhD3/sMWS365DhhqZr8FZrj7tiSbP0pQowm+byEWEN4Nnmfj73Hf/yNV+dIxs2OBOcSaQXcnVhP6vwmbxf8tvUasJgWwD/BbM3s/7v1/khDI3X27xZJlZgCLzOwJYHrHTZGUhmpQErbO2pKZ1RKrlWxIsW26qfTXE2tm24WZfQj4NXAVMMDd+xILZJZrgYFxxJr4ns5g/8nKfAOxms7+7l5PrD8ml/JsIPaBGm8o8Lcc9tWlnO4+392PBg4hFny/k+JnOwLUJ4LvHyUWoEaROkAVclmEXwBLgCHu3ge4ka6/3/ha+lA++JtbT6zptm/co5e7d/mduvt97n4KsZroKuBnYZ+IZEcBSsJ2hpmdFCQJ/CfwZ3dPV1NKZRHwNTM7OejoHmxmB/LBHfQmoD2o7Xw6l4KaWWPQEX4dMDdouupu/28C/cysT9xrdcT6L9qCMn4zl/IQC4TDzexsM6sJ7ugPBv6Qw77eBDrHnJnZMWZ2bFBb3A7sIFaTSOZPwAHE+gGfDpq69gGOBR5Lc7xhZpbXZ4qZ9Up4GLHf7xZ332FmI4Gzk/zo/zGzD5vZIcT6TO8IXr8R+JEFyTpm1t/MxiU57gAzGxskwrwLtJH69yNFogAlAL+3XcdB/TaPff2CWPPYFmLJAzllQrn708Q+aP6bWLLEo8A+7t4KXECsP6OF2IfVkix3/5yZtQHrgPOAi9z90uC4afcfNPncDrwS9NcNItYsdDaxBJGf8cGHY7bnvBn4LDAd2AzMBD7r7m/lsLtriPW1tJjZfKA+KFsLsSawzcRqicnKsR34C7DC3d8LXn6SWJ/XxhTH62hy22xmf8mhvBBrynwn4fERYArwAzNrBS4ldm0SPUrsej4EXBX0j0Ls97AEuD/4+aeIBdpEuxH7vW8g9rc7KjiulFBHxo9I3szsJmKZWP9R6rKISPlTDUpERCJJAUpERCJJTXwiIhJJqkGJiEgkRXqg7p577unDhg0rdTFERKSAnn322bfcvX/i65EOUMOGDWPp0qWlLoaIiBSQmSXOoAKoiU9ERCJKAUpERCJJAUpERCIp0n1QIiJRs3PnTpqbm9mxI3GFGOlOr169aGpqomfPnhltrwAlIpKF5uZm6urqGDZsGMESJ5IBd2fz5s00Nzez7777ZvQzauITEcnCjh076Nevn4JTlsyMfv36ZVXzVIASEcmSglNusv29KUBJxWtra2PNmjW0tbWVuigikgUFKImkMIJKe3s706bPYHDTEE4+9XQGNw1h2vQZtLcnrrQuUv4uu+wyrroq6RJfaW3dupXrr78+7+Nfe+21fPSjH8XMeOutXJYw60oBSiIlzKAyc9ZsHn16KfN+/0euufdx5v3+jzz69FJmzppdgJKLpBblWnwuAcrdef/993d57cQTT+TBBx9kn332Ca1sClASKWEFlba2NhYtWsSUy39Mw14DAGjYawBTLv8xixcv7vygiPIHh5S/QtXib7nlFg477DAOP/xwzjnnnC7vjx49unOauLfeeouOOU1XrFjByJEjOeKIIzjssMNYu3Yts2fP5uWXX+aII47gO9/5DgBXXnklxxxzDIcddhjf+973AHj11Vc56KCDmDJlCkcddRTr16/f5ZhHHnkkYc+dqgAlkZFpUMnEhg0bqG9o7NxPh4a9BlDft4HXX39dzX9ScIWoxa9YsYIf/ehHPPzwwzz33HNcc801Gf/sjTfeyNSpU1m+fDlLly6lqamJOXPm8JGPfITly5dz5ZVXcv/997N27Vqefvppli9fzrPPPstjjz0GwOrVqxk/fjzLli0LtaaUigKUREZ3QWXDhg0Z72vQoEFsa9lCy8Y3d3m9ZeObbNvawvU33KjmPymoMG+44j388MOceeaZ7LnnngA0NjZm/LPHH388l19+OXPnzuW1116jd+/eXba5//77uf/++znyyCM56qijWLVqFWvXrgVgn3324bjjjsup3LlQgJLI6C6oDBo0KON91dbWMnHiRK6/+MLO/bVsfJPrL76Qc845h1tvvTX0D458qbmxsoR5wxXP3btN166pqensI4ofd3T22WezZMkSevfuzamnnsrDDz+cdP/f/e53Wb58OcuXL2fdunVMnDgRgD322COnMucqlABlZovNbKOZvZjifTOz+Wa2zsyeN7OjwjiuVJba2lrGjx/PtbOndgkqEyZMoLa2Nqv9XTF3DqNGjmD62DFceNpJTB87hlEjR/DNyd8oyAdHrpRtWJnCvOGKd/LJJ3PnnXeyefNmALZs2dJlm2HDhvHss88C8Ktf/arz9VdeeYX99tuPCy64gLFjx/L8889TV1dHa2tr5zannnrqLjdqf/vb39i4cWNOZc1XWDWom4DT0rx/OrB/8JgE3BDScaVMdFc76PiQvvnmm3l1zUqmnHIckz95DNPGfpJRI0dwxdw5WR+zpqaGeVdfRfP613nwvntoXv86866+iqFDh+b1wRF2TUfZhpUpXS0+lxuuDocccgiXXHIJo0aN4vDDD2fatGldtpkxYwY33HADJ5xwwi4p33fccQeHHnooRxxxBKtWrWL8+PH069ePE088kUMPPZTvfOc7fPrTn+bss8/m+OOP52Mf+xhnnnnmLgEslfnz59PU1ERzczOHHXYY5513Xk7ntwt3D+UBDANeTPHeT4EvxT1fDezd3T6PPvpol/K2c+dOv2jadK/v09ebhu3n9X36+kXTpvvOnTt32e6iadP9qJNG+cLHlvmvV23wa+99wg8++hj/1rfPL0i5Eo+38LFlftRJo/yiadPzPpdstLa2en2fvp3l6HgsfGyZ9+nb4K2trTnvWwrjpZdeynjbjr+ZPn0bfMiw/bxP34a8/2bKXbLfH7DUk8SAYvVBDQbicxKbg9e6MLNJZrbUzJZu2rSpKIWTwsmkdpCsM3nvYfsy7b8XcNttt2VVU8m0dpOq+S9dTa0QNZ1C9VNINKSqxdfUaJ7uTBQrQCXr0fNkG7r7Ancf4e4j+vfvskS9FEFYTViZZjGF8SGdbT9Oth8chcrIKlQ/hURLbW0tw4cPz7lZr1oVK0A1A0PinjcBujWMmLA76zMNPGF8SOdau8n0g6NQNZ1C9VOIVIJiBaglwPggm+844G13f6NIx5Y04mtLYTdhZRp48v2QLlTtJpdzyUUuzY0i1SCsNPPbgSeBA8ys2cwmmtlkM5scbHI38AqwDvgZMCWM40ruEmtLg5qa+OmCBXzjB1eF9iGfTeDJ50O6ELWbxGbOQtZ01E8hkkKyzImoPJTFVzjJstgOPHqkj50weZdssl+v2uBDhu3nq1evzuk42WYxtba2+urVq7PKXgszEy5dpp4yssQ9uyw+6SqbLD6LvRdNI0aM8I4JDyU8bW1tDG4awrzf/3GXWkfLxjeZ+tnR/PSPS+kdjBhv2fgm08eOoXn963nVEtra2tiwYQODBg0qSL/KtOkzePTppZ3NfB21m1EjRzDv6syXIMhkP/HnAuxyXoU+Tym9lStXctBBB5W6GF1cdtll1NbWMmPGjKx+buvWrfziF79gypT8Gra+/OUvs3TpUnr27MnIkSP56U9/Ss+ePbtsl+z3Z2bPuvuIxG011VEVStck9qHevXl11Qog3M76QmcxhdGPk2lfVm1tLfvttx+Xfu+yuISSJo459jgGDW7SbBBSVsJabuPLX/4yq1at4oUXXuCdd95h4cKFeZdNAaoKDRo0iM0b30za4b/97be5Ysq5ZddZH0Y/TjZ9WV0TSh5h+/twwhnjSjIbhObxi6bG+nrMrMujsb4+r/1GcbmNM844o/P8Ro4cSXNzc17nCApQVct22435sy7YpcN//qwLqOnZk9WrVpZtZ30+NbVMM/VS1bSmz7uRJ+5Zwjvbtxdt8lnN4xdtLa2tOHR5tGQwdVAqUV9uY+fOndx6662cdlq62e8yowBVhTZs2ECfhkb6D25i6mdH861TT+TCz32SYQcdQkNjP7Zt21aVgwozzdRLV9Oq7dtAy8a/dz4v9GwQmsev+kR9uY0pU6bwiU98go9//ONZnFVyClBVpr29nZ9cex1vbfw7K/78J9ydg0ccy/UPPMnYc79B69tbq3r2gkz6stLVtNq2ttCw18DO54WcDaIY478kejzCy218//vfZ9OmTcybNy/b00pKAarKzJw1mz8te44bHniK6x54kvl3PcrmNzbwP/P+S7MXkFlfVqqa1tXTJnPi6WPpvcceRZkNQvP4VaeoLrexcOFC7rvvPm6//XZ22y2c0FI+nQuSt4477vj08oa9BnD+nGuYcsrxfP3rXy+LhIhi6OjLSuWKuXOYOWs208eOob5vA29v3cLw4Qfw5D1LeOmpx9m2tYUJEyYU9PcZX5NLHC6gefwqV/xyGz169ODII4/kpptu2mWbGTNm8G//9m/ceuutjBkzpvP1O+64g9tuu42ePXsycOBALr30UhobGzuX2zj99NO58sorWblyJccffzwQ+1+47bbb6NGjR9pyTZ48mX322afz577whS9w6aWX5nWuGgdVRdasWcPJp57ONfc+3uW9qaedyEP33Zv2Q1m6Shz3VOxxUGGN/5LMZTMOqrG+PmlCRENdHVu2bQu7aGUhm3FQqkFVkXR33K1bi9f3VEmDWRNrWt3VvMKWWJMrRs1NMletQSgs6oOqIqWeOVsp0eHTPH5SyRSgKlziAM5SzpytlOjC0XpDxRXlrpEoy/b3pj6oCtXe3s7MWbNZtGgR9Q2NbGvZwsSJE7li7hxqamqK3syWbv6/MOb6EymWv/71r9TV1dGvX79u073lA+7O5s2baW1tZd99993lPfVBVZn42kp85/nMWbOZd/VVRe8rySQlWgkaUg6amppobm5m06ZNpS5K2enVqxdNTU0Zb68aVAWKYm0limWSmEpKWpHypNnMq0gUB3CWOkFDulLSikSdAlQFKuTy5PnQ0ubRoqQVibpQmvjM7DTgGqAHsNDd5yS83we4DRhKrN/rKnf/eXf7VRNfcpk0yUR5AKealEpPTa4SJQVr4jOzHsB1wOnAwcCXzOzghM2+Bbzk7ocDo4GrzWz3fI9dbbJpkolybUUp0aUXxWZgkURhNPGNBNa5+yvu/h7wS2BcwjYO1FksJ7MW2AKooTtL2TTJaACnpBPVZmCReGEEqMFA/NKKzcFr8a4FDgI2AC8AU939fSRjuS6toNqKJKOkFSkHYQSoZCPVEju2TgWWA4OAI4BrzSzpmsdmNsnMlprZUo0z+ICaZKpLMZZwj3IzsAiEE6CagSFxz5uI1ZTifQ34jcesA/4KHJhsZ+6+wN1HuPuI/v37h1C86MrmQ0hNMtWhmKnfagaWqAsjQD0D7G9m+waJD2cBSxK2eR04GcDMBgAHAK+EcOyylMuHkJpkqkM+qd+51rrUDCyR5e55P4AzgDXAy8AlwWuTgcnB94OA+4n1P70IfCWT/R599NFeiS6aNt2POmmUL3xsmf961QZf+NgyP+qkUX7RtOlpf27nzp1+0bTp3qdvgw8Ztp/36dvgF02b7jt37ixSyaWQWltbvb5P386/i47HwseWeZ++Dd7a2pr05zr+Lur79PWmYft5fZ+++ruQsgIs9SQxQFMdFVkY4080jqgypVtQ8sLTTuLB++5JOl9hlMe8iWRCUx1FRBjJDmqSqUy59DPmmt0pUg4UoIos/kPone3b2fDXl3ln+3YlO0hO/YzK7pRKpnSdIqutreVrX/sal57zBbZu2Ux9Yz+2bdlM38Z+nHvuuaoVVblsl3CPv+FJbDLWDY+UO9WgSsGgccBA5t/1KNfd9wTz73qUxgEDk48ok6qSbep3trWuYoyvEgmLkiSKTJN0Stg6Vk9evHhxl1pXR2DrboVlkVJSkkREhNVnoDth6ZBJrUtLa0g5UoAqsnxnhNAic5JKquxOZfpJuVKAKrJ8Z4TQnbBkS5l+Uq7UB1UCmfQZJKP+K8mF/m4k6tQHFSG5TtKpO2HJheZxlHKl9J0S6ugzyJTGvEiush1fJRIFauIrM5p3TfKheRwlilI18SlAlZlc+69ERKJKAarC6E5YcqG/G4kiJUlUGM1oLtnobvycBn5LFKlNSKQKxI+fi++7nDFzJrvZbpoCSSJJTXwiFS7dOKjzTz+J4Ycdyflzf6KkGykZNfGJVKl04+c+1HsPvvjtGZoCSSIplABlZqeZ2WozW2dmSefcMbPRZrbczFaY2aNhHFdEupdu/sd/tG1j2IGH7PK6Bn5LVOQdoMysB3AdcDpwMPAlMzs4YZu+wPXAWHc/BPhivscVkcykmkniuu9OpaamJzu271pT0sBviYowalAjgXXu/oq7vwf8EhiXsM3ZwG/c/XUAd98YwnHLgrKjJAqumDuHUSNHMH3sGC487SSmjx3D6GOP4evnnacpkCSywkjTGQysj3veDBybsM1woKeZPQLUAde4+y3JdmZmk4BJAEOHDg2heMUTP8akV69eWiBOIqNj/scffP+yXcZBdQz81hRIEkV5Z/GZ2ReBU939vOD5OcBIdz8/bptrgRHAyUBv4EngM+6+Jt2+yyWLL9lqpcMPOIB/9vwQ3/qva5QdJZGnAbxSSoXM4msGhsQ9bwISe1ebgXvdfbu7vwU8BhwewrEjIdkaTW3t77PXkGHKjpKykMvAbzVfS6GFEaCeAfY3s33NbHfgLGBJwja/Az5uZjVm9mFiTYArQzh2yaVarXTGjxfwxD1LeGf79s5tlR0llUCrOkux5B2g3L0d+DZwH7Ggc6e7rzCzyWY2OdhmJXAv8DzwNLDQ3V/M99hRkG6MyR71fWjZ+PfO15QdJZVAqzpLsWgmiTylG6X/rVNP4OrfPsjew/ZVH5RUBK3OK4WgmSQKJN1qpR/72GFcctZnO9N6R40coewoKWta1VmKSfnOIUi3WumOHTuUHSUVQ6s6SzGpiS9EStWVaqBVnSVsqZr4VIMKUUeqrkglS9ZicM4553DexAm0tbXp5kxCoz4oEclKx6wUzetf5967fs9XvvIVbrnlFk79zOeUci6hUoASkZzU1taycNFinlz+vFLOpSDUByUiOVHKuYRFaeYiEiqlnEuhKUCJSE7SLYSolHMJgwKUlI3G+nrMrMujsb6+1EWrWOkmhE03SF3rSUkYFKCkbLS0tuLQ5dHS2lrSclWiTCeETbYQomZMkbAoSULKhpnhQCPQkuT9hro6tmzbVuRSVaZsB+NqkLrkI1WShAJUN/SPFx0dAcqI1Zy6vA9E+e+5XCg7TxI11tcnbakI66ZQWXxZStbEcf4FU3nppZe6tMdr4TapJMrOk0Slal5XgEoh2Zo3Dz7xJCeOGt3ZHr9jxw4t3CYVR9l5EhUKUEmkWyW3fedOLr/jLh59eikfHzVaC7cVUUNdHVbqQlQBZedJNgqZRasAlUS6Jo7avg24v8+E/7ic555b3iWITbn8xyxevFjNfQWwZds29TEVibLzJFOFbOYLJUCZ2WlmttrM1plZyuqDmR1jZv80szPDOG6hpGviaNvaQsNeA3F/n9r6vmqnL4LE8U8QS4joGXzteDTU1ZWukBUmfkLYB++7h+b1rzPv6quoqdECCFI8eQcoM+sBXAecDhwMfMnMDk6x3VzgvnyPWWipmjh+MnsqY/71LHrvsQdmu9G2bava6YsgVQdtO7GsvY6HUszD17GEjJr1qltH83rio6HAxw2jBjUSWOfur7j7e8AvgXFJtjsf+DWwMYRjFlx8E8fkTx7Dt049gYFDh/GVaRfTsvFNFv/wYg4//Ai104tIxeu4+Uu8SdxS4OOGEaAGA+vjnjcHr3Uys8HAvwA3drczM5tkZkvNbOmmTZtCKF5uOpo4Vq18iZsWLmD8OeN56t7fM/2zozvb4//fo4+onV5EpEDCaFBOlliV2JP9Y2CWu/+zow8hFXdfACyA2EDdEMqXk/b2dmbOms2iRYuob2hkW8sWxo8fzzcnf4OhQ4d21pDmXX0VP/j+ZRrMK1VDg9erU0NdHZZisG6hhFGDagaGxD1vAhIzBEYAvzSzV4EzgevN7PMhHLtgko2D+tOy51i4aHGXf0q104cvPjFCoiHT+fmkMnVk0SY+Ctn3m/dUR2ZWA6wBTgb+BjwDnO3uK1JsfxPwB3f/VXf7LtVUR5rqpfQ6pjUCzb0XFdnOzyfRV+gpjDJVsKmO3L0d+Dax7LyVwJ3uvsLMJpvZ5Hz3Xwqa6iVatvBBpyygrL0SSDV4XeP+ylvUVwgIZVCDu98N3J3wWtKECHc/N4xjFlL8OKjEGpRSyKUaZXLTNnz48BKVTiqVZpJIQlO9iOxK8/NVnygsEKoAlYKmehH5gG7aqk8Umv+0HlQ3lFJbGlHpvJUPdAy9WLx4MfV9G9i2tYUJEyZwxdw5mgKpTMUnI+3yevC1WOuuacFCEQmFbtoqR7obwY4aVKJiBqiqu+3RP5dIfjrG/Un5S9caEYUxiFXTB6VBhiKF09bWxrJly1i2bJlSziU0VROgks0MocUFRfLT3t7OhdOmsdeAgYz59Gkcd/wJDNh7EBdeNI329nba2tpYs2aNglYZSjmDeRGXtamKPijNDCFSGNOmz+CBx//EhVdd3zm7xPxZF7BpQzODB+zFmtWrO+eynDhxohIqJKmCzSRRDjQzROWIwtgMiWlra2PhooWdwQli/1MXzJ1Py6Y32fbuTrVYSF6qIkBpkGHliMLYDInZsGED9X0bkt749fpwLeNnXqppkSKsHG72qiJAaZChSPgGDRrEtq0tSW/8/tG2jWEHHrLL62qxiJZyuNmrigAF6WeGUEeuSPZqa2s5b+J5/HjGlF1u/ObPuoDdbDd2bN/1/0ktFpKtqkiSiBc/DqpXr15dFiVUR260pRv5HuW/5UrV3t7OjJkzWbBgAR/qvQf/aNtGz567c+CBB+If6q2lOSIsSv9LmkkiCa1vU36i9E8lH2hra2Pt2rUA7L///p03f5oWKZpSziBBbO01BagMFDJAKfW8PGmOvvKimVuiqds5+CISoKqmDyqRUs/LUymWnZbcdUyLpOBUPmogMll9VRuglHouItLVTqKT1RdKgDKz08xstZmtM7MuI/HM7Mtm9nzw+JOZHR7GcfOh1HMRkWjLu7fSzHoA1wGnAM3AM2a2xN1fitvsr8Aod28xs9OBBcCx+R47X5UicC4AABAWSURBVFfMncPMWbOZPnZMl45ctZ2LiJRWGDWokcA6d3/F3d8DfgmMi9/A3f/k7i3B06eAphCOm5F0Y5xqamqYd/VVNK9/nQfvu4fm9a93Bi3Nei4ilSrlRLAlLVVXYQSowcD6uOfNwWupTATuSfWmmU0ys6VmtnTTpk05Fyqb5TXiO3I167mIVLpkyUYAW0pcrkRhBKhkq1olzVE0s08SC1CzUu3M3Re4+wh3H9G/f/+cC5VLoGlra2PRokWd46JAc4iJSHWIwvIaicIIUM3AkLjnTUCXHG0zOwxYCIxz980hHDelXAONUs9FpFpFcQhHGAHqGWB/M9vXzHYHzgKWxG9gZkOB3wDnuPuaEI6ZVq6BRqnnIiLRkXeAcvd24NvAfcBK4E53X2Fmk81scrDZpUA/4HozW25mhZu/iNwDjVLPy0eqpQJ2T/JaqQcbikhuKnaqo1zn2Wtvb9ccYmUg7Zx8qV6P8N+6SDWrurn48g00GgcVbQpQIpWj6gJUBwWayqQAJVI5UgWoim+z6hjjJCIi5aVqJ4sVEZFoU4CSspRqUGHPJK+VerBhNUs31ZhIdxSgpCylGlT4XpLXSj3YsBplM9WYSCoKUCISOs1pKWFQgBKRULW1tbFw4UK+8M2L6LVHLHNWc1pKLhSgRCQ07e3tXDB1Ku+++y7XXnwR3/jkCG6+4gf8s71dc1pK1hSgRCQ0M2fNZvnqdVz/wJNcd98TXPOHR3ht1UvcNu9yzWkpWav4cVAiUhwdqwjM+/0fd1lF4Pw51zD1s6N5beWL3c5pqYH1Ek81KBEJRbpVBHru/iGOPHB/rpg7J+nPKutPklGAEpFQpFtF4P329/jJ/Pkp58FU1p8kowAlIqFIt1zNxAkTUzbZaSVrSUUBSkRCc8XcOYwaOYLpY8dw4WknMX3sGEaNHJGyaQ+0krWkVtVJEuqQFQlXTU0N866+ih98/7KM/7fimwbjg5Sy/qQqa1DqkBUprI5VBDK58dNK1pJKKAHKzE4zs9Vmts7MuvRqWsz84P3nzeyoMI6bK3XIikRLLk2DEtNYX4+ZdXk01teH+jOlkPeChWbWA1gDnAI0A88AX3L3l+K2OQM4HzgDOBa4xt2P7W7fYSxYmKitrY3BTUN2GasBsTu2iz43mjWrVjFw4MBQjynF1VhfT0tra+fzGiBZ3bihrk6TyEaMmt2zl3bxzhSf77n8TCGlWrAwjBrUSGCdu7/i7u8BvwTGJWwzDrjFY54C+prZ3iEcO2vpOmRrdu/FR/cfrua+MtfS2opD56M97vv4R3wQk2jIpmlQKl8YAWowsD7ueXPwWrbbAGBmk8xsqZkt3bRpUwjF21W6sRo733uXub+6V819IiIREEaAsiSvJdYRM9km9qL7Ancf4e4j+vfvn3fhEqXqkP3J7KmM+dez2HvYvhp/ISISAWGkmTcDQ+KeNwGJAxcy2aZorpg7h5mzZnPR50ZTs3svdr73LmP+9Sy+Mu1iYNfxF8OHDy9VMUVEqloYNahngP3NbF8z2x04C1iSsM0SYHyQzXcc8La7vxHCsXPSMVZjzapVvPfOP5hzx118deal9AimYdH4CxEpFw11dRh0edRAyiy9VD/TUFdXilNIKe8A5e7twLeB+4CVwJ3uvsLMJpvZ5GCzu4FXgHXAz4Ap+R43DAMHDmTSpEks/uHFGn9RQRL/+Wro+o8YxX9GkVxs2bYNd9/lAbCT1IlByX7G3SOX1Zp3mnkhFSLNPFF7ezszZ81m8eLF1PdtYNvWFiZMmMAVc+eknNhSRCTKopZG3p1UaeZVH6A6aPyFiFSKSglQqiIEOsZfiIhINFTlXHxSecpl6hYRyZwClFSExNkjNFuEVLNyydLrjpr4REQqTNSy8XKlGpRIAjUXikSDalAiCTqaCxOZmgtFiko1KBERiSQFKKkIldIpLCIfUBOfVIRK6RQWkQ+oBiUiIpGkACWSQM2FItGgACWSoFxmepbqU21DIBSgRALV9s8v5afaZkxRkoRIQOOfRKJFNSgREYkkBSiRLKgZUKR48gpQZtZoZg+Y2drga0OSbYaY2R/NbKWZrTCzqfkcU6SUqq0PQMpfOd9U5VuDmg085O77Aw8FzxO1A9Pd/SDgOOBbZnZwnscVEak6uQyBKOebqnwD1Djg5uD7m4HPJ27g7m+4+1+C71uBlcDgPI8rEjqNf5Ioa6yvTxpUGurqKnYIRL5ZfAPc/Q2IBSIz2yvdxmY2DDgS+HOexxUJXaX+k0tlqMYs024DlJk9CAxM8tYl2RzIzGqBXwMXunvKTwIzmwRMAhg6dGg2hxARkQrSbYBy90+les/M3jSzvYPa097AxhTb9SQWnP7H3X/TzfEWAAsARowYkeyGQaRkaog1+yV7XSRMyZr0DGgAtpSkRMWXbx/UEuCrwfdfBX6XuIGZGbAIWOnu8/I8nkhJ1aXoj0r1ukiuUiY3ZLmfcu5bzTdAzQFOMbO1wCnBc8xskJndHWxzInAOMMbMlgePM/I8rkhJpJqnDyjbVF6pbOU8t2ReLRPuvhk4OcnrG4Azgu8fJ3mriEjFqMYObCmd+A/UcqgJ5UpN5yIiZaaj1l7pNNWRiIhEkgKUiEgElXNyQ1jUxCciEkHlkMRQaKpBiYRAd7sCsHuSTE4zY3dTnlguFKBEQpBtKm85zzAtqe2k67glD16X7KmJT6QIks0KoLR0kfQUoESKIHGclBp8RLqnJj4REYkkBSgREYkkNfGJiISkJ8mbb3sWuyAVQgFKpAQaSP5BprT08vZelUxBVCwKUCJF0FBXlzRDr5KX6xbJlwKUSBEoCIlkT0kSIiIFogHZ+VGAEhGhMMEk5aq4GpCdETXxiYigRSejSDUoERGJpLwClJk1mtkDZrY2+NqQZtseZrbMzP6QzzFFqk22TU/q95BKkW8NajbwkLvvDzwUPE9lKrAyz+OJVJ1s+zHU7yGVIt8ANQ64Ofj+ZuDzyTYysybgM8DCPI8nIlI2tE5YfvJNkhjg7m8AuPsbZrZXiu1+DMwEur0qZjYJmAQwdOjQPIsnIpKZdIOpc6Xxb/nptgZlZg+a2YtJHuMyOYCZfRbY6O7PZrK9uy9w9xHuPqJ///6Z/IhI2VJ/UWaK8XvKdtHJYpWrmnVbg3L3T6V6z8zeNLO9g9rT3sDGJJudCIw1szOAXkC9md3m7l/JudQiFUKpzZmJ6u8pquWqFPn2QS0Bvhp8/1Xgd4kbuPt33b3J3YcBZwEPKziJZC7bfgz1e0ilyLcPag5wp5lNBF4HvghgZoOAhe5+Rp77F6l62fZjqN9DKkVeAcrdNwMnJ3l9A9AlOLn7I8Aj+RxTRESqg2aSEKlS6uCXqFOAEimhXPqLwgosiQN6G+Jej1rAimq/WlTLVSk0WaxICeXSX1SozLEWiGxGWvzvqbG+vnNWjI5gCqVZ/FH9fYWlACUiZUWp3dVDTXwiIhJJClAiUjBKxJB8KECJVKnEDv5CKPbM6gp8lUUBSqTMhJU5ljj3XJgZaR01JxL21Zj1nrLTEfhUc6sMSpIQKTOFyhwLc78pExlC2HfKWceJZSKmPb4SKcqKalAikpNS1VI6Amlis+GWdGXlg+CoGlX5UA1KRHJSTrWUKI/xktRUgxKRotJMC5Ip1aBEpKjck9VlslOI1W8lehSgRCR0hQ4g3SV0pDq+lBcFKBEJXannqEs8fkfKu5QX9UGJSE7KaSbvciqrfEA1KBHJSalrSdkop7LKB/KqQZlZo5k9YGZrg68NKbbra2a/MrNVZrbSzI7P57giIlL58m3imw085O77Aw8Fz5O5BrjX3Q8EDgdW5nlcERGpcPkGqHHAzcH3NwOfT9zAzOqBTwCLANz9PXffmudxRUSkwuUboAa4+xsAwde9kmyzH7AJ+LmZLTOzhWa2R6odmtkkM1tqZks3bdqUZ/FERKRcdRugzOxBM3sxyWNchseoAY4CbnD3I4HtpG4KxN0XuPsIdx/Rv3//DA8hIiKVptssPnf/VKr3zOxNM9vb3d8ws72BjUk2awaa3f3PwfNfkSZAiYiIQP5NfEuArwbffxX4XeIG7v53YL2ZHRC8dDLwUp7HFRGRCpdvgJoDnGJma4FTgueY2SAzuztuu/OB/zGz54EjgMvzPK6IiFQ4C2PixkIxs03AawXa/Z7AWwXad6lV8rlBZZ9fJZ8bVPb5VfK5QWHPbx9375J0EOkAVUhmttTdR5S6HIVQyecGlX1+lXxuUNnnV8nnBqU5P83FJyIikaQAJSIikVTNAWpBqQtQQJV8blDZ51fJ5waVfX6VfG5QgvOr2j4oERGJtmquQYmISIQpQImISCRVTYAysy+a2Qoze9/MUqZKmtlpZrbazNaZWVlMyZTFulyvmtkLZrbczJYWu5zZ6O46WMz84P3nzeyoUpQzVxmc32gzezu4VsvN7NJSlDMXZrbYzDaa2Ysp3i/ba5fBuZXzdRtiZn8M1uxbYWZTk2xT3Gvn7lXxAA4CDgAeAUak2KYH8DKxGdh3B54DDi512TM4tyuA2cH3s4G5KbZ7Fdiz1OXN4Hy6vQ7AGcA9xFbuPg74c6nLHfL5jQb+UOqy5nh+nyA2QfSLKd4v52vX3bmV83XbGzgq+L4OWFPq/7uqqUG5+0p3X93NZiOBde7+iru/B/yS2JpXUdftulxlJpPrMA64xWOeAvoGExaXg3L9O8uIuz8GbEmzSdleuwzOrWy5+xvu/pfg+1ZiC8sOTtisqNeuagJUhgYD6+OeN9P1AkVRJutyAThwv5k9a2aTila67GVyHcr1WkHmZT/ezJ4zs3vM7JDiFK0oyvnaZaLsr5uZDQOOBP6c8FZRr123y22UEzN7EBiY5K1L3L3LTOvJdpHktUjk4ac7tyx2c6K7bzCzvYAHzGxVcEcYNZlch8heqwxkUva/EJufrM3MzgD+F9i/4CUrjnK+dt0p++tmZrXAr4EL3X1b4ttJfqRg166iApSnWbsqQ83AkLjnTcCGPPcZinTnluG6XLj7huDrRjP7LbGmpigGqEyuQ2SvVQa6LXv8B4O7321m15vZnu5eCZORlvO1S6vcr5uZ9SQWnP7H3X+TZJOiXjs18e3qGWB/M9vXzHYHziK25lXUdbsul5ntYWZ1Hd8DnwaSZiJFQCbXYQkwPsgqOg54u6OZswx0e35mNtDMLPh+JLH/1c1FL2lhlPO1S6ucr1tQ7kXASnefl2Kzol67iqpBpWNm/wL8BOgP3GVmy939VDMbBCx09zPcvd3Mvg3cRyzTarG7ryhhsTM1B7jTzCYCrwNfhNi6XATnBgwAfhv879QAv3D3e0tU3rRSXQczmxy8fyNwN7GMonXAP4Cvlaq82crw/M4Evmlm7cA7wFkepFFFnZndTiybbU8zawa+B/SE8r92GZxb2V434ETgHOAFM1sevHYxMBRKc+001ZGIiESSmvhERCSSFKBERCSSFKBERCSSFKBERCSSFKBERCSSFKBERCSSFKBERCSS/j9qlJ44VKWSBwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# create figure\n",
    "figure = plt.figure()\n",
    "axis = figure.add_subplot(111)\n",
    "\n",
    "axis.scatter(X[y == 0, 0], X[y == 0, 1], \n",
    "             edgecolor='black',\n",
    "             c='lightblue', marker='o', s=40, label='cluster 1')\n",
    "axis.scatter(X[y == 1, 0], X[y == 1, 1], \n",
    "             edgecolor='black',\n",
    "             c='red', marker='s', s=40, label='cluster 2')\n",
    "axis.set_title('Empirical Data Points with Labels')\n",
    "\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"overview\"></a>\n",
    "## Overview of K-Means and Agglomerative Clustering\n",
    "\n",
    "There exist several algorithms for partitioning data into partitions, two of the more common of which are called K-Means and Agglomerative Clustering.\n",
    "\n",
    "The K-Means algorithm approaches the clustering problem by partitioning a set of data points into disjoint clusters, where each cluster is described by the mean of the samples in the cluster. The mean of the samples in a particular cluster is called a centroid; the K-Means algorithm finds the centroids and associates data points with centroids in such a way as to minimize the within-cluster sum-of-squares.\n",
    "\n",
    "For more information on the K-Means algorithm and its implementatin in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#k-means\n",
    "\n",
    "In the code cell below, we instantiate the `KMeans` algorithm from the `sklearn.cluster` module and apply it to our data using the `fit_predict` method. We see that `KMeans` identifies two centroids; one located at about (-0.23, 0.56) and the other located at (1.17, -0.05)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.23011109  0.56790752]\n",
      " [ 1.17408322 -0.05027964]]\n"
     ]
    }
   ],
   "source": [
    "from sklearn.cluster import KMeans\n",
    "\n",
    "\n",
    "km = KMeans(n_clusters=2, random_state=0)\n",
    "y_km = km.fit_predict(X)\n",
    "print(km.cluster_centers_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Agglomerative Clustering algorithm behaves a little bit differently and does not identify clusters using centroids. Instead, it recursively merges the pair of clusters that minimally increases a given linkage distance. Put another way, the Agglomerative Clustering algorithm identifies the two data points that are \"closest\" out of all the data samples. It then takes those two data points and identifies a third data point that is \"closest\" to those two data points. The algorithm continues in this fashion for each data point; finding the next data point that is \"closest\" to the preceeding cluster of data points, where the definition of \"closest\" depends on the distance metric chosen.\n",
    "\n",
    "For more information on the Agglomerative Clustering algorithm and its implementatin in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#hierarchical-clustering\n",
    "\n",
    "Below, we instantiate the `AgglomerativeClustering` algorithm from the `sklearn.cluster` module and apply it to our data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cluster import AgglomerativeClustering\n",
    "\n",
    "\n",
    "ac = AgglomerativeClustering(n_clusters=2,\n",
    "                             affinity='euclidean',\n",
    "                             linkage='complete')\n",
    "y_ac = ac.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can visualize the results of both algorithms applied to the data. Visually, we see that neither algorithm ideally clusters our data. The ideal algorithm for this unique set of data would recognize that both sets of samples are generated from two different equations describing two different half circles."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAADQCAYAAAAK/RswAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO2de3wU9dX/30cCgrmQhHu4ReoVLSAiDwgVxSqgLfZC+1hbtQIPReQOAmp/ilYtKKLgDXnAemu1Pl5aa1XQolIvVbCiBQVRRIhBwiWQBEUSOb8/djZuNrOb3exmL7Pn/XrNK7sz35k5M8l8cub7Pd9zRFUxDMMwDMPwEkck2wDDMAzDMIx4Yw6OYRiGYRiewxwcwzAMwzA8hzk4hmEYhmF4DnNwDMMwDMPwHObgGIZhGIbhOczBMVIaEdkqIt9PATt+KSIrk22HYSQSEfm1iLyWbDsaQ7KeWRE5U0RKEn1eN0TkeRG5NNl2JAtzcNKA4H/yInKhiJSLyBCXtmeKiIrIU0HrezvrX0mAyREjInkicoeIbBORKhH52PneNo7nmCsij8RyDFX9o6qeGy+bDCPeiMgrji4cmWxbEo2IFDv6luVf15TPrIj0F5HnRGSfiOwVkbdF5LI4nyPmlztVHaGqD8bLpnTDHJw0w/HG7wbOV9VXQzTbBZwuIm0C1l0KfNTU9kWDiLQA/gGcBAwH8oDTgT1A/ySaVodA0TSMVEREioHvAQqMTKoxTYCINEu2DX5EZCCwCngVOAZoA1wOjEimXYGIj4z//57xNyCdEJFxwG3AMFV9I0zTQ8BfgAud/ZoBPwf+GHS8E0TkRecNZJOI/Dxg2/ki8q6IVIjIdhGZG7DN/7Z0qdPzsltErgnY3l9E1jr77hSRhSHsvAToBvxYVT9Q1cOqWqaqv1PV51yu/wERuTHge52uYBGZLSKfi0ilcz1ni8hw4Grgv50eovectq1FZLmI7HD2udEvok63/OsicruI7AXmBnfVO9c/XkQ2O2/Nd4uI+O+3iNzm3JdPRWRi8NulYcSZS4B/AQ/ge5mpRUTaiMjfnOdxjfO3Hvi3fK7zvOwXkXtE5FURGet2EhE53TnGfufn6QHbXnGO/YbzrP3NOfcfA85dHNA+nP48ICL3Or0kB4CzwmkSsNr5uc8598DAZ1ZElojIgqBr+auITHc+F4nIkyKyy3lmJ4e517cCD6rqfFXdrT7eUdWfuzV2nv1jgq7tRudzWxF5Vr7tCfqniBwhIg/j08a/Odczy2k/wLm/+0TkPRE5M+j+3yQirwNfAj2cdWOd7b8WkddEZIGjWZ+KyIiA/Y8WkdWOfr7kaFpMPd9JR1VtSfEF2Ao8CewEejfQ9kygBF9PyFvOuvOAFcBY4BVnXTawHbgMyAL6AruBkwKO8118TnAv59w/crYV43tT/F+gFdAb+Bo40dn+JnCx8zkHGBDC1sfwCUVD1/595/MDwI3B1+p8Pt65nqIAG7/jfJ4LPBJ03L8A9zn3oT3wNvAbZ9uvgRpgknNvWjnrXgvYX4FngXx8QrQLGO5sGw98AHQBCoCXnPZZyf5bssWbC/AxMAE4FagGOgRse8xZjgJ6Os/Ja862tkAF8BPnb32Ks/9YZ/uvA9oWAuXAxU7bXzjf2zjbX3Hs+A7Q2nkGPgK+77R/CPiD07Yh/XkA2A8MwqdBLYlMk7ICrjvQ9jOc84nzvQD4CihyjvcOcC3QAugBbMH3Ihl8n48CvgHOCvO7OBNHl5zvChwT8P0BHB0Dfg8sAZo7y/cCbNyKo33O9874erfPc2w+x/neLuD+b8PXI57lHO+VoN9lNfA/QDN8vU6lAed7E1jg3IPB+P4uHgl1nemwWA9O+nAOvje0/0TSWH09PIUicjy+t7uHgpr8ANiqqn9Q1RpV/Tc+J2qUs/8rqvof9fWqvA88CgTH/Fyvql+p6nvAe/gcHfA9RMeISFtVrVLVf4Uwsw2wI5LriYBvgCOBniLSXFW3quonbg1FpAO+7uSpqnpAVcuA23F6vBxKVfVO5958FeKc81R1n6puA14G+jjrfw4sUtUSVS0H5sXh+gzDFREZDHQHHlfVd4BPgIucbc2AnwLXqeqXqvoBEBiTcR6wQVWfUtUaYDHwRYhTnQ9sVtWHnefiUWAj8MOANn9Q1U9UdT/wPPCJqr7kHPv/gFOcdmH1x+Gvqvq6o0EHI9SkUPwTn6PxPef7KOBNVS0FTsPnJNygqodUdQu+l7cLXY5TgM+5iJduVQOdgO6qWq2q/1TH23DhV8Bzqvqccw9eBNbi+x36eUBVNzj3tNrlGJ+p6v+q6jf4/g46AR1EpBu++3Ctcw9eA56J0zUmDXNw0ofxwHHAMv9QSAQ8DEwEzgKeDtrWHfgvp6tzn4jsA34JdAQQkf8SkZedLtv9zvmDA38DhfBLfL01AGMcWzc63dI/CGHfHnwPWMyo6sfAVHy9NWUi8piIFIVo3h3f282OgGu/D19Pjp/tEZw21PUXBe0fybEMo7FcCqxU1d3O9z/x7TBVO3xv86H+Huv8rTr/XEPNACoCPgta9xm+ngU/OwM+f+Xy3f+MhNUfFzsj1SRXnOt6DF+vE/gcQP+QfXegKMiWq4EOLocqBw4TJ93CN9z1MbBSRLaIyJwwbbsDPwuyc3CQLQ1pTa1mqeqXzsccfL/bvQHrIjlWymMOTvpQBpyN7w3kngj3eRhft/VzQX+44PvjfVVV8wOWHFW93Nn+J3wefFdVbY2vGzUix0pVN6vqL/A5DPOBJ0Qk26XpS8CwENvcOICvi9hPoBiiqn9SVf/brDrnxvkcyHZ8Q2ptA649T1VPCjxchDa5sQPf8JSfrjEcyzBCIiKt8PUYDhGRL0TkC2Aa0FtEeuMbOq0h9N9jnb9V5+UpsG0gpfierUC6AZ83wvSG9AfqP4PhNCmS5/VRYJSIdAf+C1+Pkd+WT4NsyVXV84IP4Ojom/h6xSLlS0LolqpWquoMVe2BrydsuoicHeKatgMPB9mZraqBPcSN1a0d+Hr8A+1Me90yByeNcLpThwLDReT2CNp/iq8L9xqXzc8Cx4nIxSLS3FlOE5ETne25+Dz6gyLSH6fLOxJE5Fci0k5VDwP7nNXfuDR9GN9D+6QTcHiEE5R4tYjUExdgHXCeiBSKSEd8PTb+cx4vIkPFN0X2IL63Rf85dwLF4swqUNUdwErgNvFNUz9CRL4jLtPuG8njwBQR6Swi+cDsOB3XMIL5Eb6/8574hkj7ACfiG5K5xBmKeApfoPxRInICviFrP38HvisiPxJfEPwVBL04BPAcPs24SESyROS/nfM+2wi7G9IfN8Jp0i58PSs9Qu2squ867ZYBK1TVr01vAxXim6TQSnyTBE4WkdNCHGoW8GsRuVKcmariS8PxWIj264CLnOMOJ2BYTUR+ICLHOI5lBb7fZaBuBV7PI8APRWSYc6yW4ptoEcohjRhV/QzfcNdcEWkhvpliP2xgt5THHJw0Q1W343NyRonI7yNo/5rjGAWvrwTOxTfOXIqv63I+vjgW8PX83CAilfiC7x6PwszhwAYRqQIWAReq6kEXG77GF4C4EXgR3wP+Nr5u57dcjvswvlifrfgclD8HbDsSX6zLbuda2uPrZgbf2D/AHhH5t/P5EnzBdB/g63Z+gvh1O/+vY9/7wLv4/jHU4O7kGUYsXIov7mWbqn7hX4C7gF86TstEfEG/X+B7hh7F14OJM6z1M+AWfEPGPfH9o/s6+ESqugdf7MwMp+0s4AcBQ2MRE4H+uBFSk5yelZuA153hmwEhjvEoPs35U8C+3+D7Z94H+BSfhizDd8/cbH8DnwYPBbaIb6blUnzPuRtTnOP7h+H+ErDtWHw92VX4eobuUdVXnG2/B37rXM9MR/svwKdru/C9HF5J/P6P/xIYiO93eyM+fa33d5BO+KOnDcNoIpypmEtUNbh73zASjojMBzqqar0Mt04vZwnwS1V9OeHGGSmDiPwZ2Kiq1yXblsZiPTiGEWecbu7znG78zsB11A/yNoyE4Az/9hIf/fFNAng6YPswEcl3hnevxhfXEmrmo+FRnCHC7zhD9sPx9Rb9paH9UhlLPGYY8UeA6/F18X6FL87h2qRaZGQyufiGZorwTVa4DfhrwPaB+IZs/EO2PwqTGsHwLh3xxWu1wdeLd7kTt5S22BCVYRiGYRiew4aoDMMwDMPwHGk5RNW2bVstLi5OthmGYYThnXfe2a2q7ZJtRyyY1hhG6hNKa9LSwSkuLmbt2rXJNsMwjDCISHDW27TDtMYwUp9QWmNDVIZhGIZheA5zcAzDMAzD8Bzm4BiGYRiG4TnSMgbHSAxVVVWUlpZSVFRETk5OwzsYhmFESabpTHV1NSUlJRw8WK96jdEALVu2pEuXLjRv3jyi9ubgZBiRiElNTQ2zZs9h+fLl5BUUUlG+lzFjxnDL/HlkZdmfjGEY4TGdCU1JSQm5ubkUFxfjq7FpRIKqsmfPHkpKSjj66KMj2seGqDKEmpoaps+YSecuXTl72Ag6d+nK9Bkzqampqdd21uw5vPr2Whb+7WUWvfAaC//2Mq++vZZp02fw0UcfUVVVlYQrMAwj1YmHzsyaPYeqqirPas3Bgwdp06aNOTdRIiK0adMmqp4vc3AyhHBiEkhVVRXLly9nws13UNC+AwAF7Tsw4eY7WLr0Ps46Z1hY0TIMI3OJh84suW8JRZ27NOggpTPm3DSOaO+bOTgZQDgxuf/+++u8JZWWlpJXUFjbzk9B+w4UtO/IVfc9HFK0DMPIXOKlM61ycpl974NhHSTDiARzcDIAn5gUuIpJXn4BpaWlteuKioqoKN9LednOOm3Ly3ZyoGI/Be07hhStxuLl7mjDyBTipTNff/UVxSecVLtvvLTGdCY8c+fOZcGCBVHvt2/fPu65556Yz3/XXXdxzDHHICLs3r075uNBnBwcEblfRMpEZH2I7SIii0XkYxF5X0T6BmwbLiKbnG3mpsdI8ENcU1PD4jvvYueOUlcxqdhXTlFRUe26nJwcxowZwz1XT61tX162kzvnTGHoTy+kVXY24C5a0QpINOP1hmE6kzo0lc4smPYbzh71i1qdgfpak4k6k8rOWWMcHFXl8OHDddYNGjSIl156ie7du8fNtnj14DwADA+zfQRwrLOMA+4FEJFmwN3O9p7AL0SkZ5xsyihCPcQzr5zFm+veZ+iP/5s750ypIyb3XD2V0aNH15vlcMv8eQzp348ZI4cyZfggJpwzgI7divnV9Ktr2wSKVmMFJNLxesNweADTmaTSVDozdfhgpo88i22bPuT8X42p086vNe3bt884nWkq5+yhhx6iV69e9O7dm4svvrje9jPPPLO2RMnu3bvx12PbsGED/fv3p0+fPvTq1YvNmzczZ84cPvnkE/r06cOVV14JwK233sppp51Gr169uO666wDYunUrJ554IhMmTKBv375s3769zjlPOeUU4l73TVXjsgDFwPoQ2+4DfhHwfRPQCRgIrAhYfxVwVUPnOvXUUzUTqKys1E2bNmllZWWDbadNn6F9Bw/RZavf1Sc3luqy1e/qKYPO0KNycnXZ6nf18fXbdOTo8ZrTOl87dO2uLY5sqRMnTdbq6uqQ5/R/vmLipHrH7jt4iE6bPiPkuQO3h7q2vNb5tfv4l2Wr39XW+QURXbOR2gBrNU76opp4nVHTmnrES2cCz7ljx47ac4fTEq/ozAcffBBx28Zcc0OsX79ejzvuON21a5eqqu7Zs0dVVa+77jq99dZbVVV1yJAhumbNGlVV3bVrl3bv3l1VVSdOnKiPPPKIqqp+/fXX+uWXX+qnn36qJ510Uu3xV6xYof/zP/+jhw8f1m+++UbPP/98ffXVV/XTTz9VEdE333wzrH3du3evtc0Nt/sXSmsSFYPTGQh010qcdaHW10NExonIWhFZu2vXriYzNNlUVVXxwQcfMHHS5Ii99lDBfVf8fhE1NdW0zM6hWVYWl866liWr1vDbpY/QrkNHJk28ojbfhNubwrXXzaVHjx7ccfvCOm9aM0YOZUj/ftwyf15UgYWBhAsyDB76MowIiVlnwLQmlNbEQ2egvtYcf8KJLLlvKS1btqzXq+PXmmv/328zTmcaq60NsWrVKkaNGkXbtm0BKCwsjHjfgQMHcvPNNzN//nw+++wzWrVqVa/NypUrWblyJaeccgp9+/Zl48aNbN68GYDu3bszYMCARtndGBLl4LjN7dIw6+uvVF2qqv1UtV+7dvWqoqc9gQ/9md8/l2XLljFwxEgWPvtKg12q4R7io3Ly2LpxQ+26VtnZtMrOoapif50x8XDduFlZWSy8bQEl27fx0ornKdm+jYW3LeDgwYO8/vrr5LTOj1pAwgUZBo/XG0aExKwzYFoTSmvioTMQvdbccP1c3nrrrYzTmaZyzlS1wenWWVlZtTEygXlnLrroIp555hlatWrFsGHDWLVqlevxr7rqKtatW8e6dev4+OOPGTPGN+yYHRBblQgS5eCUAF0DvncBSsOszzgCH/olL6/h7pVv8MW2rTyy8OYGvfbwMxIO8H93LQg7Jh7pm0JOTg7HHXccLVu2rBXIseMvZ9fOHSydO4dvAt76yst2sn/fXg4cOOBqc6ggw1Dj9YYRAaYzETBr9hxeeWtN1FoTq85AdFrTo0cPrr1uboM6U7GvnLy8PNcg3HTWmaZyzs4++2wef/xx9uzZA8DevXvrtSkuLuadd94B4Iknnqhdv2XLFnr06MHkyZMZOXIk77//Prm5uVRWVta2GTZsWJ3f5eeff05ZWVmjbI2VRDk4zwCXOLMcBgD7VXUHsAY4VkSOFpEWwIVO24zC99Avq/fQT5q3iH888ShfHTgQ1msP9xCPGzeOoQP/y3V4yU+0bwp138Be594X/8W2jzex/Mbf1p779hmXU1PzDSN/Mipkt3eo7uhA2wwjCkxnGmDfvn0sWbKEK36/qJ7WrHrysbBaE6vOQHRaE4nO3H3VFLp378Zxx58QcogtXXWmqZyzk046iWuuuYYhQ4bQu3dvpk+fXq/NzJkzuffeezn99NPrTNn+85//zMknn0yfPn3YuHEjl1xyCW3atGHQoEGcfPLJXHnllZx77rlcdNFFDBw4kO9+97uMGjWqjgMUisWLF9OlSxdKSkro1asXY8eObdT11cEtMCfaBXgU2AFU43tbGgOMB8Y72wXfLIZPgP8A/QL2PQ/4yNl2TSTn81rg36ZNm7Rdp851guD8S9tORXrn8/9sMCiuurpap02foa3zC7RrcQ9tnV+g06bPqA3uCxdEGE0gXri2LY48Utt26KTNjzxSc1rna6ucXB05erzet2pN2MC4SAIcowm4NlID4hxknGidUQ9qza8vu0wLO3R01ZqO3Y9uUGti0Rn/9ki0JrzOtNSibsXaKjtbWx6VrQXtO+hRuXkNak2q6Ew0QcYN3e9MJJog47jOcEjU4jXR2bFjhzY/8siQTsNtT78YceR8Yx/QSKP1N23apF2Ke7gKZLtORdqj58l61wuv1x6j9+ln6MjR4xs9a8H/gOe1ztcuxT00r3V+xj/g6UK8HZxkLF7SmsrKSs3Na61H5ea5ak12Xmu964XXI9KaWByBSLQmnM50KT5af/yTn2ivAYPqHCMWrUmkzkTj4PixF7xvScVZVEYY9u7dS4uWLVkwdVydrsjFc6bQomUrrrvkpxF3qfrjZKLtvoy0GzfcuPD+vXuYvvA+OhX7Kr0Gdn23zM5pVGBcOuewMIxUwj881KXHsfW0ZsHUcWQ1a8Y1F/4gIq1prM5AZFrTUPzJSy/9g8m33OU6zNYYrUl1nYnlfmcy4nN+0ot+/fqpPwmRF5g0eQovvvYGxcf35PXnnyE7rzX7dpfRpcex7Ny+lc2bNtGxY8eE2FJVVUVpaSlFRUUhH6bpM2by6ttra2OGyst2ctecKWz96EOWv/ZevfZXDBvExJtv59YrLqNk+7aIH9Kqqio6d+nKwr+9XGfMvrxsJzNGDo3qWEbiEZF3VLVfsu2IBS9pTVVVFR06FfGdk3tTVNyjjtYccUQzXn15FSeeeGLCnqmGtMZNZ+65eiq9j/sOL656hUUvvFZvn8ZoTaJ15sMPP+TEE0+M2/EyDbf7F0prsoJXGImlqqqKhx56qPbhuvjK/0d52ReIHMGMH3+fsWPGJsy5gW/fFMJxy/x5zJo9hxkjh5KXX0DFvnIuvvhiPln/HuVlO+uJRGX5Xv7vrgVRB8ZFEpDYkK2GYXyLHj7MlFvurKc1V/703EY5N5G8EIWiIa1x05nRo0dz7f/7Ld2Lj46b1pjOeBcbokoywQ9Xq+xsio7+Dp2KjyYvv5AJl49v1HGbsnaJW66KOxcvYuzYsfUi/m+b9htQZejA/4p61kI657AwjFSjtLSUNu07uGpNm3YdohrSSUR9p1D5t/Lz811nFzVWa0xnvIs5OEmkqqqKAwcOsL98j+vDdeirL+nWrVtUx0xkYbngcWG32jIDep9MybbPWHjbgjrZTCM9frrmsDCMVMH/spOXlxfyH3lllP/IExmz4hZ/Ek+tMZ3xMG6Rx6m+pPvMhuCI/aNycvXk/gPjUm+kKWqXREs8I/5tmmT6gs2iSipuM4P69f+vmPUhleo7xUtrEqkzjZlFlQgCa1FFQ3l5ud59990xn/+iiy7S4447Tk866SS97LLL9NChQ67toplFZUHGSSA4eG536edcf9nPKd9dRtt2HWrHmm+ZPy+qNxEvB+XGMtZvJAcLMk4ubkG6d181hWbVX7P5o4/qxLVEozUfffQRZw8b4RrkO3X4YF5a8XzaxqxEojOxalGqBhnPnTuXnJwcZs6cGdV+W7du5Qc/+AHr16+PeB+/A3LEEd8OIj333HOMGDEC8JWEOOOMM7j88svr7RtNkLENUSUYt1TlbYs6c8PDT9E8qzl/feqJ2rHmaId00rmwXEOEmyZZmJeHiNRbCvPykmCpYSSfcIUxN3/0ERs//KBOXEs0WuPlmJVwOpPI4X8/TaVtDz30EL169aJ3795cfPHF9bafeeaZ+B373bt3U1xcDMCGDRvo378/ffr0oVevXmzevJk5c+bwySef0KdPH6688koAbr31Vk477TR69erFddddB/gcoRNPPJEJEybQt29ftm/fXuec5513Xu319e/fn5KSkpiuEWwWVUII9PjDOSGtCwrJzs5udA9FoPAE9+Cku/CEo7yy0rVyokSQHtwwvEKkOpOXX0BFRUWje1kCY1aCp3B7OWYlMO4o8JpnzZ7DwtsWNMk5m0LbNmzYwE033cTrr79O27ZtXWtRhWLJkiVMmTKFX/7ylxw6dIhvvvmGefPmsX79etatWwf4qolv3ryZt99+G1Vl5MiRrF69mm7durFp0yb+8Ic/cM8994Q8R3V1NQ8//DCLFi1q9DX6sR6cJsTN47/zrrvZv9c9qDhWJ8SC5Qwj80i0zkD61ndqLJEWCU0HVq1axahRo2jbti0AhYWFEe87cOBAbr75ZubPn89nn31Gq1at6rVZuXIlK1eu5JRTTqFv375s3LiRzZs3A9C9e3cGDBgQ9hwTJkzgjDPO4Hvf+14UV+WOOThNiNtMgzfefY/jTzihyZyQTBMew8h0kqEzoaZwRzusni54afhfVRGRsG2ysrI4fPgwAAcPHqxdf9FFF/HMM8/QqlUrhg0bxqpVq1yPf9VVV7Fu3TrWrVvHxx9/zJgxYwDIzs4Oe97rr7+eXbt2sXDhwmgvyxVzcJqIcB7/Rx9tYmCfXk3ihKSz8DRl7h7D8CLJ0hk/6VhCoDE646W4o7PPPpvHH3+cPXv2ALgOURUXF/POO+8A8MQTT9Su37JlCz169GDy5MmMHDmS999/n9zc3DrVwocNG1anV+vzzz+nrKysQbuWLVvGihUrePTRR+sEH8eCOTgxEuphCRtrk1/I5EkTm9QJSSfhSUbwnmGkE6mqM+lELDrjpeH/k046iWuuuYYhQ4bQu3dvpk+fXq/NzJkzuffeezn99NPZvXt37fo///nPnHzyyfTp04eNGzdyySWX0KZNGwYNGsTJJ5/MlVdeybnnnstFF13EwIED+e53v8uoUaPqOEChGD9+PDt37mTgwIH06dOHG264IfaLdZs7Hu0CDAc2AR8Dc1y2Xwmsc5b1wDdAobNtK/AfZ1tEeTNSITdFQ9VnUylXRKoTae6eUDkvCnJzFai3FOTmJvIyjCAifZ6jWTJNa0xn4kesOhOvXDnR5MExbatPNHlw4iE4zYBPgB5AC+A9oGeY9j8EVgV83wq0jeacyRYd1cgellRIupfqRCLQDYm8kZrE28HJRK0xnYkP8dSZWJMLpmqiv3QhGgcnHkNU/YGPVXWLqh4CHgMuCNP+F8CjcThv0vCPe4/+7c18daCKrw4ccI2ot4DfhokkeC+RaeGNlCajtCYwvqZldg6ln35Cy+wc05lGEE+dSafh/0wnHg5OZyAwY0+Js64eInIUvi7mJwNWK7BSRN4RkXGhTiIi40RkrYis3bVrVxzMbjzbtm3jiObNmf3z87hp/CX85qx+PHjLDeQVtqkTUZ/OAb+JoqHgvby8PM9MzzRiJqO0prS0lNz8Ap554D5+c1a/Wq155oH7yG2dbzoTBfHQmXhOgvB1OhjREu19i4eD4zbfLJQVPwReV9XAsO1BqtoXGAFcISJnuO2oqktVtZ+q9mvXrl1sFsfIvUvuo2O3YhY9+wp3r3idRc++wmcbP2D5jb9Nu4j6ZNNQ8F5FRYVnpmcaMZNRWlNUVET5nt1s/XBDHa3Z+uEGyvfuqaczOTk5tUn+zPGvSyw6k5tfwKTJk+M2CaJly5bs2bPHnJwoUVX27NlDy5YtI94nHi5+CdA14HsXINR/nQsJ6jJW1VLnZ5mIPI2vG3p1HOxqEqqqqnjooYfq1HsqaN+BSfMWccWw0xk7Zmxt12VNTQ2zZs9h+fLl5BUUUlG+lzFjxkRdY8rr3DJ/HrNmz2HGyKH16uMcPHjQNTvzpX2PoerLLzn++OPrHKsgN5e9FRWJvgQjMWSU1gDo4cNMnr+4jtZMnr+YSSMG12lnWtMwjdGZ8rKd7Nm1k3c3bo5bBuMuXbpQUlJCskci0pGWLVvSpQiABgoAACAASURBVEuXiNvHXGxTRLKAj4Czgc+BNcBFqrohqF1r4FOgq6oecNZlA0eoaqXz+UXgBlV9Idw5k1kAL1yhufFnncYrL62kZ8+egHuxu3uunsqQ/v2aLLV3OhOqiJ3bfRx7xinuKcyx7t9UId7FNk1rvmXK8MH8I6CopWlN5ESjM3dfNYVN773L4udWe66AsZdosmKbqloDTARWAB8Cj6vqBhEZLyLjA5r+GFjpFxyHDsBrIvIe8Dbw94YEJ9mEG8s99NWXdOvWDfBWau9EESp4zy2I0sg8TGt8lJftpDJgKNy0Jjqi0Zk+xx9DQZu2NkSepsQl0Z+qPqeqx6nqd1T1JmfdElVdEtDmAVW9MGi/Lara21lO8u+bykSa8MlLqb2TjVsQpZGZmNaY1jQVbjqzeNEiKveVeyKDcSZimYwbQSTTMr2U2jtVsOmZRqZhWpN4AnXGSxmMM5GYY3CSQTLHxQMJNZbrx8bFmw4RsRicFCfeMTjJwLTGqKmpYeq06Tz44IO0LiikqmJ/bXCyBXCnBqG0xhycJsQ/s+H++++vF7VvD0ZstBCh2mV9c+BQGv5NexFzcBKHaU3TUHd2WgH79u7h15f+mt/dcD1lZWUhHU4jsZiDk0QaevsyoidcD04gNm08eZiDk3hMa+KLW8/YHTMnsOWD9RS2bWfT8VOEJptFZTSMxY4klsCqdOURVLE1DK9gWhM/Qs1Om7rgHuSII5j35AorG5PimINjGIZhGEGEm52WW1BIedkXNh0/xTEHJwriWYvEMAwjFKY1ySfc7LSqfeUUtO8I2HT8VMYcnAaoqqrigw8+YOKk+NUiMQzDCMa0JrUINUV88ezJDP3phbTKzq5dZ9PxUxNzcEJQU1PD9Bkz6dylK2d+/1yWLVvGwBEjWfjsKzbumgIU5OYiUG8pSKpVhhE9pjWpS3AeokkjBrO7tITzfzUGCJ0Tx3rgUgNzcEIwa/YcXn17LQv/9jJLXl7D3Svf4IttW3lk4c027poC7K2oQFVrl4LcXADKCXJ4nPWFeXmISL2lMC8vWZdgGIBpTSoTnN24tKSE//7pT5j1k3NdEy8GOqvWA5d8bJq4C1VVVXTu0rVOxXDweetTf3gWS1atoVV2NlOHD+algIJ3RupiiQETj00TbxjTmvQkmoKdlnCx6bFp4lEQLno+J7+A8rIvbNzVMIyYMa1JT9ym41vR09TDHBwXGoqeFznCapEYhhEzpjXewYqeph7m4LgQKnp+wdRxNG/enGsu/EG9gneGYRjRYlrjHcI5q/v37bUeuCQQFwdHRIaLyCYR+VhE6oX7i8iZIrJfRNY5y7WR7pss3Kr4fn/QQF579RVKtm9j4W0LLDW3R7Bg4/TBtMZIVUI6q9N+Q03NN1x73VwLNk40gTNRGrMAzYBPgB5AC+A9oGdQmzOBZxuzr9ty6qmnaqKorKzUTZs2aWVlZcLOacSfgtzcwAoOtUsBqDqfjfgCrNUY9UWj1AvTGiOZVFdX67TpM7RVdrYWduio2XmtdeTo8XrfqjXad/AQnTZ9RrJN9CShtCYePTj9gY9VdYuqHgIeAy5IwL4JwWq7eAP/tHKo6+HsTaZRRrSY1hgpTVZWFjdcP5esZllMX7iE+15ey6WzrqVtUWcLNk4C8XBwOgPbA76XOOuCGSgi74nI8yJyUpT7IiLjRGStiKzdtWtXHMw2wBJSGWmFaU2a4nWdCby+0tJSWhe24cRT+9dmOwYLNk4G8XBwxGVdcGKRfwPdVbU3cCfwlyj29a1UXaqq/VS1X7t27RptrOEjMCHVWecMo2OnIiZOmuyJMWJL6udZTGvSjODEd0Wdu3DZ6NHs27cv2abFBbfEfovvvIuK8j2uwcY23T+xxMPBKQG6BnzvAtRxUVW1QlWrnM/PAc1FpG0k+yYCr79duDFr9hxeeWsNA0eMZP++cnILCn0p4gcNTnsnp7yysn6wjbPeiJwUdBRNa9IMf5bmW59+kb5Dz+UbVZ557nmKOnfxRIbfwCzUi154jYV/e5k3173PccefUC/Y2Kb7u9OkOuMWmBPNAmQBW4Cj+TZ476SgNh35Nmtyf2AbvjeqBvd1W+IV+OcPCMtrna9dintoXut8nThpsm7YsEErKys9G/RXWVmpea3zddiFl2jv08/QZavf1Sc3luqy1e/qCX1P04mTJifbxJjACRwOXoDQwca5uck2O+UIdx8j3D/eQcae0Zrc1q31stGjtby83PM6s2z1uzpy9Ph6WnPKoDPSOug28Pqe3Fhauyxb/W7t/5LW+QXatbiHts4v0GnTZ2h1dXWyzU45YtUZ5xiuWhPz3ENVrRGRicAKfDMV7lfVDSIy3tm+BBgFXC4iNcBXwIWOUa77xmpTpAR63/602gumjuPhP/2J6q+/Rg8fpk37DlSU72XMmDHcMn+eJ6ZrlpaWclRuHqv/9hS3PrmiTtbNmXcsZdoPz+T3N9/kyTeNvRUVyTbBaCRe05pbp46jqEtXUPWszuQVFIIIL/75kXpac8XvFzFj5FBuuH5uWmpNuMR+rQsKmTTxCn5/802uJR2MxJCxtajC1YCZfN4ZHPPdPkyev9hz9URqamqYNn0G9y29j5y8fA59fZCzR/2CX02/mmaOqE4ZPoh/rHghbeveWN2p+BDrfbRaVD5Cac3SuXP47KMPmXnHUs/pDMC+ffvo3LUbCrRsdZSr1qRzja1w/0NmjBxKyfZt5tREQDz02mpRBRHK+26ZnUNNdXWtcwPeqicya/Yc3nj3Pe598V8s++e7LHr2FT7b+AGPLLwZ8D2clfv2WSCcYcQJN6356sABXnvur7XODXhLZwBu+N2N9Oh5Mnc+tzqk1qRz0G2oxH4Wa5M6ZKyDEyqt9taNGzgqN8+T9URCFYObNG8Rq558jB1bP/XEw1mQm4tAvaUgNzepdhmZiZvWlJd9QW7rfE/qDHyrNVMX3ONprXHLQm2lNVKHjHVwQnnfjy6+lYNfVnlyil+4MePmLY5kzs9GeOLh9Cf1C14s/iY6sqjvJPqjdY3I8WvN3VdNqdUVkSPYu6vMkzoDfq0p8LzWZGVlsfC2BZRs38ZLK56nZPs2brh+Llu2bPFEL1wiaEqdyVgHB3ze9/dO68ukEYMZM7g3E84ZyJb179H1mBNYPHuy57od27dvz56yna6iWnPoIJs/2uS5ujcpONU5bajBZbqZs96IjptvupFm1V9z+TkDGPO9Pky7YCh5hW08qTM1NTUsvvMudu4ozRitycnJoUePHlx73dw6OXG8MBW+qWlKncloBycrK4sj5AiO63UKMxf9L/e/8R/ufP6fNMtqxub332XSiMFM8VC34w2/u5GCdu3rieodMycwdsxYOnbsmGQL44/lxGmYUE6gN/71pAZXX/NbDrdoye1/XcXvHnqS2/+6irYdO3lSZ2bNnsOb695n6I//mzvnTKnnwHlVa9xy4rz69lpmzU6Zuq5JJRk6k7GzqCB8FPz0kWex6cMPqaio8MQUP/+13vr0i/z9keWsevIxcvILqCzfC6qUbPuM/Pz8ZJsZd2xGVcOEu0duFOTmRjTcZ7OofGSiziz828vkFbbhkYU3s+rJx8jOa0152U7GjRvH7Qtv80zPjR+bUdUwTaUzzrFtFlUwYfMY5BdSUVHhmeJ3/mttW9SZS2ddy5JVa7hmyUPc9/JaCtu0paysLNkmGilOQW6uxTI1gkzUmYL2HWiWlVWrNb9d+gjtOnRk0sQrPOfcQPjfsReCxhNJPHUmox2cUDOpvBLoF0jwtbbKzqbo6O9w8ECV567ViB82rBc7mawz4NOaVtk5VFXs99S1BpJJv+OmoKl0JmMcHLcaMJmUxyCTrtUwkkmw1mTSs5dJ1xpIpl53yuM2nTbVl2jqw7jVmwqsCeLfngk1QzLpWv1Y7amGIVwtmEbWhnGOG9daVMlY4qU1mfTsZdK1BuJ23VdMnKRr1qzRf//7356rNRYtTaUzzrFdtcbzQcbTZ8zk1bfX1ia3C5UOvaqqKmNqhgReK5Ax1w3QQoRql/XNgUNp+CzEg8K8PNdu4QJgb8D3aAOzMy3IOBKtyVSdycnJyZhrr6qqYtu2bdx977384f4/cORR2XxZWUFW8+ZcesklTLxiAt26dfP0PXCjqXQGQmuNpx0ci2wPTU1NDbNmz2H58uXkFRR6rtBfKCKN5I8mgt9rxGvmWSY5OKY1oclErZk+YyYvvvZGbSbn8rKdLJxxOVvWv09ufgFff3nA8/egIeI5wzUjZ1FZZHtoLGdDfTIxoNYtN0WozKJW6iI0pjWhyTStqaqqYtnyZfXKVEy/7V6aNW/O7c++4vl7EEyydCYuDo6IDBeRTSLysYjU+42JyC9F5H1neUNEegds2yoi/xGRdSISe3KbAMJFtu8r30NehmazDVWTykuF/ozIcEuE6B/CCx7PToUeLdOa9CITtaa0tJS8fPcyFbkFhZSXfeH5exBMsnQmZgdHRJoBdwMjgJ7AL0SkZ1CzT4EhqtoL+B2wNGj7WaraJ97d2aEi2xdMHUez5i04/oQTMzKVtr1tGumIaU36kYlaU1RURMW+cldnt2pfOQXtfVmcvXwPUoV49OD0Bz5W1S2qegh4DLggsIGqvqGq5c7XfwFd4nDeiAis9jr+rNO4YtjpdD/uRJa+8k7GdRP6sZwNBlBbjyu4i7gwiTY1gGlNmpGJWpOTk8PYMWO5Y+aEOs7u4tmTGfrTC2mVnV27zqv3IJhkFe6Nh4PTGdge8L3EWReKMcDzAd8VWCki74jIuFA7icg4EVkrImt37doVsXH+aq8bP/yALysruO3plxg3dx7NsrIyrpvQTybnbGhOiHHfZBqVJELW6XK2p2BRUtOaNCNTteaW+fM4e9AAp5BzHyacM5Ad27Zy/q/GAJlxDwIJV1CzKQsgx8OBcisl4RoCLSJn4ROdwQGrB6lqqYi0B14UkY2qurreAVWX4nQ39+vXL+qpXxUVFRS0bUen4qPrrA/sJjzuuOOiPWzacsv8ecyaPYcZI4eSl19Axb5yRo8enfaF/hoieCq4f+piOXX/kC2g1vcQS2oFW5vWpCGZqDVZWVncsXAhN95wA5s3b6ampoaHHn6EWT85N2PuQaT4H7Cm0Jp4ODglQNeA712AeoOKItILWAaMUNU9/vWqWur8LBORp/F1Q9cTnVgJ7CoNnsa5f99eDhw4QFVVVUZ40/Dt2+YN18/NiNwUoUiFwNlUJEV7tNJeazIx4DiTtSYnJ4dTTjkFgNNOO43f33xTxt2DZBKPIao1wLEicrSItAAuBJ4JbCAi3YCngItV9aOA9dkikuv/DJwLrI+DTbX406YDrl2ld8ycQE3NN4z8ySg6d+macYGAOTk5nin0Z8SPvQ03SQYpqzWB5Rks4Ngd0xq7B4kmZgdHVWuAicAK4EPgcVXdICLjRWS80+xaoA1wT9AUzQ7AayLyHvA28HdVfSFWm8CXXGr6jJl07tKVs4eNoHOXrhzWw3yvX19mjBzK1OGDmTRiMHt3fsHtz6zKiPwMhpHOpKLWuOnM9BkzufmmGy3g2DCSjGczGYdLm37D9XPZvHkzQ848i9uffSXjM49mSgp1oy6RpE5vTFZRP5mQybih8gxffPEFxxx7HPOfeKFOTE4m6gyY1mQqydIaT2Yybii5FEB2djatC9tkVH6GYEK9fe7bt69e5XXDe+ytqKhNrhUYVO0PuLbsxeGJJIldJAHHmYBpTWaTLK3xpIMTSXKpTMzPEExwCvVbn36RPz/5FEWdu9QRoUyLFchEAgUo1bIXpyqmM5FjWmP4SaTWeNLBiURUMjU/gx+3t8+/P7KctkVduPOF1ywmyTAawHQmMkxrjGThSQcnUlEJzDw6dfhgZowcypD+/TIiN0Hw2+dXBw7wjyceZfL8xRlTM8YwYsF0JjJMa4xk4dk67ZEkl8rk/AzBuTrKy74gt4Hu9kxKTmYYkWA60zCmNd9iQdaJxZM9OPCtqJRs38ZLK56nZPs2Ft62gKys+j5dJuYmCH77LGjfkYq9ezIqVqAwL682TXjgkmLlCYwUxnSmYUxrQgdZW8xR0+JZB8dPpopKJAR2nV/102EI1CsQ5+VYgZC1mFKrPIGRBpjOhCfTtWbqtOmseHU1N//57xZzlEjcoplTfTn11FPViB+VlZW6adMmLS8v12nTZ2jr/ALtWtxDW+cX6LTpM7S6ujrZJsadgtxcN99GC5yfRuwAazUF9CKWxbQmvmSa1lRXV+vESZO1+ZFHaseu3TU7r7WOHD1eH1+/TZetfldb5xdoZWVlss1Me0JpjWcT/RmNJxPGiUXEtUqjv5pjOj4XqUYmJPozYsPrWjN9xkxe/tcaJs5bVJsI8s45U+h+Qk8unXUtU4cP5qUVz3s25ihRZFSiPyM2rLvdMIxE4GWt8U+P9zs34AuknjRvEauefIwdWz+tE3MUWM/MiA/m4BhGCCwI2TCMxhIuEWR2XmvuucYXc9SyZUsLQG4izMExjCD8KcMtCNkwjMYSLhFkedlOzho4oDbNQGCWZwtAjh8Z5+BYN6DREFaewIgHpjWZTahEkHfNmcK4ceO4687FHDx4sMF6ZkbjyRgHx/IQREegOHtRqAtyc2uLvAUuXiouGckQW7oNw4nIcBHZJCIfi0i9V1zxsdjZ/r6I9I1033hhWhM5XtcZtyzWZw04jdsX3gZEVs8s1UlpnXGbWhXtAgwHNgEfA3Nctguw2Nn+PtA30n3dlsZM3Zw2fYb2HTxEl61+V5/cWKrLVr+rfU4/QydOmhz1sbxMdXW1Tps+Q/Na52vn7kfrUTm52uqobO1S3EPzWud7cipnKABVl4U0mUYeif1NeY3EeZo40Az4BOgBtADeA3oGtTkPeN7RnAHAW5Hu67bES2tOGXSGTps+ozG30ZNkms74p8cHTwmvrKzUvNb5tX8r/iWdppAnW2ec47tqTcw9OCLSDLgbGAH0BH4hIj2Dmo0AjnWWccC9UewbM27F3grad2DivEUsXXofEydNtrcrh8Dx4FPPHsYx3+1jBfGMVKE/8LGqblHVQ8BjwAVBbS4AHnJ0719Avoh0inDfmAmlNVf8fhFL7lvCvn374n3KtCTTdCbUbDErxtq0xGOIKuVFJ1w3YEH7jrz85r889TA1lkBxbpmdk/EF8TJhGCvN6AxsD/he4qyLpE0k+8ZMOK1plZPLtOnT433KtMN0pi6ZXoy1KYmHg5MQ0RGRcSKyVkTW7tq1KyoDw0WzH6jYz4SbMudhCkegOEdSEM/r7K2ocB228EIQsn8MHKAwybZEgbisC87IGKpNJPv6DhCj1uwvd6+z9PVXX/HUU0+bzpjO1CGaembpRrJ1Jh4OTkJER1WXqmo/Ve3Xrl27qAz0dwPeNWdKnW7AO+dMYehPL6RT8dEZ8zCFo27V345UhnAKvVoQL5OoM+09ybZEQQnQNeB7FyD4oQ3VJpJ9gdi1ZtRPf8qCab+ppzVnj/oFrQsKTWdMZ1zxYtLDZOtMPBychIhOrNwyfx6DT+3DhHMGMOGcgUz94Vl0P6Env5p+dUY+TG4EjgcfPFDF2aN+weLZk21sOE0JOcTm0jZNhuHWAMeKyNEi0gK4EHgmqM0zwCXObKoBwH5V3RHhvnFh4W23sW3Th0z5wZlcMWxQrdac/6sxpjOYzniNlNYZty74aBYgC9gCHM23sxNOCmpzPnVnNrwd6b5uSywF8K6YOEl7nnqa3vXC67XR6n0HD7EZDg7l5eV62ejRmpefr126H63ZuXnaKts3u8GrBfEyBRI8K4wmKLaJb5bUR/hmRF3jrBsPjHc+C76JC58A/wH6hdu3oaWxWjNt+gztPXCw3vjHv+gj72w2nQnCdMa7JFpnnHO6ak1GiY7qt9MTvV7FNloCp212Ke6huXmt9deXXabl5eUhpzga6YUXHJxEL43VGtMZd0xnvE8qOTgZW03c61Vso2X6jJm8+vba2umt/m7iIf37sfC2BSH3s/uYPoSroN4UOmDVxO35CMZ0xvskWmecc1o18UC8GNDVWELl7gg3VdOytaYfNu098ZjOfIvpTGaQSjqTsQ6O8S2NSRduBeLSDy9PezdSH9OZzCCVdMYcHCNsniC3WR+NeRMzDCOzMZ0xEo05OEbU6cK9UCDOMIzEYjpjJJr0T5VoxIVb5s9j1uw5zBg5lLz8Air2lTN69GjXdOF1E3V9Kz6WT8gwjHCYzhiJJGNnURnuRDpbYfqMmax68y1+NnEmxSecxMEDVRHNhjCalsK8PMorK+utL8jNTfgYuM2iMkJhOpPepJLOQGitsR4cow7+WR/hqKmp4bAe5qP33mXBlHF8WVlB8xYtGDtmjBWISzLllZXuUzRdxMgwkoXpTHqTLjpjDo4RNbNmz+Gfa/7NnS+8VpvL4u6rpnDEEUd4okCcYRjJJ5zOHDx4kM2bNwNw7LHH2jR8wxUbojKioqqqis5durLwby/XGRe/8IQiql3aNwcOpeHfWLqSjCRbYWyxISqjUYTSmd2lnzPtgqHoN99w5FE5fFlZQVbz5vzP2LEsuPUWe8FKEA3pTKKHsCzRnxEXQs1sqKZu5Vj/4ub0GMmjMC8PEam3FOblJds0w6gllM78/ZHldDv2RO584XWWv7aOe158k2O+24f/e/ovlhsnhfAPYQUvbk5PU2IOjhEVoXJZGOlBqgiPYYTDTWe+OnCAl/7vT8y84746eXEmz1/M/vK9LF++3HLjGHUwB8eIilC5LIzUIIv6KdIFC7Yz0gs3ndm6cQMtjzrKNS9ObkEhR+XkWG6cBJEuOpNq9hhpgFsuCyO5hBzzBvbiEx/DSCeCdWZ/+V6qqw+558XZuwdRJc+GWpuUdNOZmHpwRKRQRF4Ukc3OzwKXNl1F5GUR+VBENojIlIBtc0XkcxFZ5yznxWKPkRiysrJYeNsCSrZv46UVz1OyfVu9NoV8+8dusR5NT8ihp6RaFT9MazKPYJ35vGQ748dfzh0zJ9T26uwu/Zzf/urH1Bw6xFG5eRx/wolWjLMJSTediXWIag7wD1U9FviH8z2YGmCGqp4IDACuEJGeAdtvV9U+zvJcjPYYCSSwUnJz6nZVluMedGyxHkYjMa3JUAJ15pb58zh70AAmjRjMmMF9mHz+EPLbtuOeF99kyctrrBhnipAqFcVjdXAuAB50Pj8I/Ci4garuUNV/O58rgQ+BzjGe10gxDgVVjjVSh0BhSRXhaQSmNQZZWVncsXAhZV98wbN/eYrmWVnMvGOpFeNMAQJ1JFUqisfq4HRQ1R3gExegfbjGIlIMnAK8FbB6ooi8LyL3u3U7B+w7TkTWisjaXbt2xWi2YWQOgcKSKsLTCExrjFpycnJo164d+W3aWjHOFCEVdaRBB0dEXhKR9S7LBdGcSERygCeBqarqvwv3At8B+gA7gNtC7a+qS1W1n6r2a9euXTSnNoy0o4VLrhoRoYWkWhhf/DCtMaIhVMoKK8YZOV7XmQZnUanq90NtE5GdItJJVXeISCegLES75vgE54+q+lTAsXcGtPlf4NlojDfSk1Qr1JaK+BMnBhNKdgpyc13rwKTB0FMtpjVGNAROJZ9w8x215RzuuXoqo0ePJicnx7SmAbyuM7FOE38GuBSY5/z8a3ADERFgOfChqi4M2tbJ3+0M/BhYH6M9RooQ6kHIIn0KtaUTGSDWpjVGPdxSVlQfOMC/X3uV2xf6OulMa+JHuulMTLWoRKQN8DjQDdgG/ExV94pIEbBMVc8TkcHAP4H/AIedXa9W1edE5GF8XcYKbAV+EyBCIbH6MOmJv36JEPqtwQKUfaRSTanGEs9aVKY1RjiqqqooLS2lqKiI3Nzc2mfHtCY8XtAZCK01MfXgqOoe4GyX9aXAec7n1wjR46WqF8dyfsMwMgPTGiMc/qnkhhGIlWowDMMwDMNzWKkGw0hB/IkT3dYbhmHEA6/rjPXgGAmnAPdCbakaiZ8M/IkTg+9JNVjJC8OIENOa8HhdZ8zBMRKGP4tucN2SgtzclEwSFW8K8/Jcc06EE5GQtV9sFohhuBKYrbs8aH0aJbaMiWi1xqs6Y0NURsLwuqg0hE2PN4ymJxKd8Xp+HNMaH9aDY6QU4d48GtMDYhiGEUxDPRamNd7AenCMlKKhNw97KzEMo6mxHhBvYD04RlII9YaU6h53U77ZuR0boDDmIxtGZmI6E9mxwZs6k+q/Z8OjhHxDSrgl0dGUb3bh7knwfbFZIIbRMKYzURwb7+mMOTiGkSBiKVSXTmnTDcNILuG0JtzMKK/pjDk4hpEgvDA7wzC8glsvTrr3WPgJpzX+IalMwGJwjJQjVGKuwPwWgUtzsBkPhmFEhT8nTuDidwxCaU0W9bXGdCZ1MQfHSDncMmuWV1ZSXllZJ1mXf6nm22meBUH7mAgZhhEOt6BbN60BarUmUTpTSF0HC7yRYThRxOTgiEihiLwoIpudnwUh2m0Vkf+IyDoRWRvt/ob3CPWG5HdsGptZs9xln3hm5CzIza2t3xIsPLGKTkP3JJMxrTEaQyTPVGO0JhE648/EHO/zZJLOxNqDMwf4h6oeC/zD+R6Ks1S1j6r2a+T+hofYW1ERtos4kUQzJXNvRQU1NI24pdI9SUFMa4yoSaVnKlqdaaqA31S6J01NrA7OBcCDzucHgR8leH/DiImQKdtpnMPS3GKBmgrTGiNtibfOgMUdRkKsDk4HVd0B4PxsH6KdAitF5B0RGdeI/Q2jSQjZPd3I47n17qR7wboUIWFaIyLjRGStiKzdtWtXnMw3Mpl46wwuxzKtqU+D08RF5CWgo8uma6I4zyBVLRWR9sCLIrJRVVdHsT+OWI0D6NatWzS7Gh4nVM6HYEK9RcWbQuoLl4h4ppBfU5EqWqOqS4GlAP369fNWYhCj0aSazoBpTUM06OCo6vdDbRORnSLSSVV3iEgnoCzEMUqdn2Ui8jTQH1gNX18pwgAABehJREFURLS/s6+JTgYRTVK8wAe5MC8vbIKr4D+cpsgI4Q8MDMbq2IQnVbTGyCwi1Zpk6ExWiH39/7hNa8IT6xDVM8ClzudLgb8GNxCRbBHJ9X8GzgXWR7q/kZk0NhAuEQF0IWchxO0MhgumNUaT0BjNSFSgbqgJDTVxPYt3idXBmQecIyKbgXOc74hIkYg857TpALwmIu8BbwN/V9UXwu1vGKlAqGmTfhELFp29iTIsMzGtMTyJF6dnpwoxlWpQ1T3A2S7rS4HznM9bgN7R7G8YiaKA0Cnbw72NherWzsLerpoC0xojnYm3zvhzcRnhsVpURkazF59QRJtzwmq9GIYRKU2hM2Ba0xBWqsHIGBKVwTOTMoUahlGXRD7/pjXhsR4cI2NI1LRJm55pGJlLIp9/05rwWA+OYRiGYRiewxwcwzAMwzA8hzk4hmEYhmF4DnNwDMMwDMPwHObgGIZhGIbhOSTaefmpgIjsAj5r5O5tgd1xNKexmB11MTvq4gU7uqtqu3gak2hi0Bov/P7iidlRF7OjLrHa4ao1aengxIKIrFXVfmaH2WF2pI8d6Uaq3Dezw+zIZDtsiMowDMMwDM9hDo5hGIZhGJ4jEx2cpck2wMHsqIvZURezI71JlftmdtTF7KiLp+3IuBgcwzAMwzC8Tyb24BiGYRiG4XHMwTEMwzAMw3N43sERkZ+JyAYROSwiIaehichwEdkkIh+LyJwmsKNQRF4Ukc3Oz4IQ7baKyH9EZJ2IrI3j+cNen/hY7Gx/X0T6xuvcUdpxpojsd65/nYhc2wQ23C8iZSKyPsT2RN2Lhuxo8nvhnKeriLwsIh86z8oUlzYJuSfpiulM7XFNZ+qex7Tm23MkXmdU1dMLcCJwPPAK0C9Em2bAJ0APoAXwHtAzznbcAsxxPs8B5odotxVoG+dzN3h9wHnA84AAA4C3muB3EYkdZwLPNvHfxBlAX2B9iO1Nfi8itKPJ74Vznk5AX+dzLvBRMv4+0nkxnTGdCWGLac2350i4zni+B0dVP1TVTQ006w98rKpbVPUQ8BhwQZxNuQB40Pn8IPCjOB8/HJFc3wXAQ+rjX0C+iHRKgh1NjqquBvaGaZKIexGJHQlBVXeo6r+dz5XAh0DnoGYJuSfpiukMYDpTD9OaOjYkXGc87+BESGdge8D3Eurf+FjpoKo7wPeLBtqHaKfAShF5R0TGxenckVxfIu5BpOcYKCLvicjzInJSnG2IhETci0hJ6L0QkWLgFOCtoE2pdE/SFdMZ05lgUum5Stj9SJTOZDV2x1RCRF4COrpsukZV/xrJIVzWRT1/PpwdURxmkKqWikh74EUR2eh437EQyfXF5R7EwY5/46srUiUi5wF/AY6Nsx0NkYh7EQkJvRcikgM8CUxV1YrgzS67ZFSOCdOZhk1zWWc6E55Uea4Sdj8SqTOecHBU9fsxHqIE6BrwvQtQGk87RGSniHRS1R1Ol1tZiGOUOj/LRORpfN2tsQpPJNcXl3sQqx2Bf/Cq+pyI3CMibVU1kQXhEnEvGiSR90JEmuMTnT+q6lMuTVLiniQT05kGMZ2JnpR4rhJ1PxKtMzZE5WMNcKyIHC0iLYALgWfifI5ngEudz5cC9d74RCRbRHL9n4FzAdeo9yiJ5PqeAS5xotgHAPv9Xd1xpEE7RKSjiIjzuT++v9E9cbajIRJxLxokUffCOcdy4ENVXRiiWUrckzTHdMZ0JpiUeK4ScT+SojOxRCinwwL8GJ9X+DWwE1jhrC8Cngtodx6+qO5P8HU5x9uONsA/gM3Oz8JgO/BF/b/nLBviaYfb9QHjgfHOZwHudrb/hxAzQRJgx0Tn2t8D/gWc3gQ2PArsAKqdv40xSboXDdnR5PfCOc9gfN3A7wPrnOW8ZNyTdF1MZ0JfX6bqjHMe05pvbUi4zlipBsMwDMMwPIcNURmGYRiG4TnMwTEMwzAMw3OYg2MYhmEYhucwB8cwDMMwDM9hDo5hGIZhGJ7DHBzDMAzDMDyHOTiGYRiGYXiO/w88twVGikWnFgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 576x216 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 3))\n",
    "\n",
    "\n",
    "ax1.scatter(X[y_km == 0, 0], X[y_km == 0, 1],\n",
    "            edgecolor='black',\n",
    "            c='lightblue', marker='o', s=40, label='cluster 1')\n",
    "ax1.scatter(X[y_km == 1, 0], X[y_km == 1, 1],\n",
    "            edgecolor='black',\n",
    "            c='red', marker='s', s=40, label='cluster 2')\n",
    "ax1.set_title('K Means Clustering')\n",
    "\n",
    "\n",
    "ax2.scatter(X[y_ac == 0, 0], X[y_ac == 0, 1], c='lightblue',\n",
    "            edgecolor='black',\n",
    "            marker='o', s=40, label='cluster 1')\n",
    "ax2.scatter(X[y_ac == 1, 0], X[y_ac == 1, 1], c='red',\n",
    "            edgecolor='black',\n",
    "            marker='s', s=40, label='cluster 2')\n",
    "ax2.set_title('Agglomerative Clustering')\n",
    "\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"dbscan\"></a>\n",
    "## Clustering using DBSCAN\n",
    "\n",
    "Unlike K Means or Agglomerative Clustering, DBSCAN is a density-based approach to spatial clustering. It views clusters as areas of high density separated by areas of low density. This approach has several advantages; whereas K Means focuses on finding centroids and assoicating data points with that centroid in a spherical manner, the DBSCAN algorithm can identify clusters of any convex shape. Additionally, DBSCAN is robust to areas of low density. In the above visualization, we see that Agglomerative Clustering ignores the low density space space between the interleaving circles and instead focuses on finding a clustering hierarchy that minimizes the Euclidean distance. While minimizing Euclidean distance is important for some clustering problems, it is visually apparent to a human that following the density trail of points results in the ideal clustering. \n",
    "\n",
    "For more information on the DBSCAN algorithm and its implementation in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#dbscan"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cluster import DBSCAN\n",
    "\n",
    "\n",
    "db = DBSCAN(eps=0.2, min_samples=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's fit our model to the data and generate predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_db = db.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, let's visualize the model applied to our data. We see that the DBSCAN algorithm correctly identifies which half-circle each data point is generated from."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3de5hU1Znv8e+r4GmH7pZuAbm26BMwXoJIOiiigWAMSiaQyZhziI6aQEIIXkAwQPQ88XKiB4xy1Bh0GCFjzEyMx2SOJFHxFnWSaAQTYmy5xhsdEFpAoTUEWt7zR1cz1dVV1XXZVbWr9u/zPPXQVbW69trsrnprrfWutczdERERCZvDSl0BERGRZBSgREQklBSgREQklBSgREQklBSgREQklHqUugLp9OnTx4cOHVrqaoiISAG99NJL77h738THQx2ghg4dypo1a0pdDRERKSAzezPZ4+riExGRUFKAEhGRUFKAEhGRUAr1GJSISNgcOHCA5uZm9u3bV+qqlJ2qqioGDx5Mz549MyqvACUikoXm5mZqamoYOnQoZlbq6pQNd2fnzp00Nzdz3HHHZfQ76uITEcnCvn37OProoxWcsmRmHH300Vm1PBWgRESypOCUm2z/3xSgREQklBSgJJJaW1vZuHEjra2tpa6KSCCuv/56br311qx/791332Xp0qV5H/+uu+7iIx/5CGbGO++8k/frgQKUlJEggkpbWxtz513NoMFDOGfi+QwaPIS5866mra0twJqKdBbmL0S5BCh35+DBg50eGzt2LE8++STHHntsYHVTgJLQCzKozF+wkGdfXMOSn/+KOx77NUt+/iuefXEN8xcsLEDNJeoK9YXohz/8ISNGjODUU0/l4osv7vL8+PHjDy0T984779CxpmlTUxOjR49m5MiRjBgxgk2bNrFw4UL+/Oc/M3LkSL75zW8C8N3vfpdPfOITjBgxguuuuw6AN954gxNPPJFZs2YxatQotmzZ0umYp512GkGvnao0cwm9+KBS1+8Ydu/YztJr5jB/wUKW3JZ5l0ZrayvLly8/9DoAdf2OYdbNtzNv8gRuvOF6qquraW1tZevWrQwcOJDq6upCnZZEQFB/u/Gampq46aab+M1vfkOfPn3YtWtXxr97zz33MHv2bC666CL279/Phx9+yKJFi3jllVdYu3YtAI8//jibNm3ixRdfxN2ZPHkyzz33HA0NDWzYsIEf/OAHgXQJZkItKAm1jqAy6+bbuwSVFStWZNVlsnXrVmrr6g+9Toe6fsdQ27uOt956S91/Epgg/3bjPf3001xwwQX06dMHgPr6+ox/d8yYMdx8880sXryYN998kyOPPLJLmccff5zHH3+c0047jVGjRrF+/Xo2bdoEwLHHHssZZ5yRU71zoQAlodZdUNm6dWvGrzVw4ED27N7F7h3bOz2+e8d29ry7m6V336PuPwlMkH+78dy923TtHj16HBojip93dOGFF7Jy5UqOPPJIJk6cyNNPP5309b/1rW+xdu1a1q5dy+bNm5k+fToAvXr1yqnOuVKAklDrLqgMHDgw49eqrq5m+vTpLL1mzqHX6+hyufjii7n//vsD/7Yr0RXk3268c845hwcffJCdO3cCJO3iGzp0KC+99BIADz300KHHX3vtNY4//niuvPJKJk+ezMsvv0xNTQ179+49VGbixImd/ub/8pe/sGPHjpzqmi8FKAm1dEFl2rRpWY8R3bJ4EeNGNzJv8gTmnHcW8yZPYNzoRr4x8+sF+bYblDBngUlyQf/tdjj55JO59tprGTduHKeeeipz587tUubqq6/m7rvv5swzz+yU8v2Tn/yEU045hZEjR7J+/XouueQSjj76aMaOHcspp5zCN7/5TT7zmc9w4YUXMmbMGD72sY9xwQUXdApgqdx5550MHjyY5uZmRowYwVe/+tWczq8Td8/7BqwAdgCvpHjegDuBzcDLwKhMXvfjH/+4S2Xau3evb9iwwffu3dttuaamJr/8iiv9qN51PmTo8X5U7zq/au48P3DgQGDH37t3r9ce1dvvfe4P/tP1Ww/d7n3uD35U77pu61koBw4c8KvmzvPao3r74KHHe+1RvfM+d8nPq6++mnHZjusX5N9uuUv2/wes8WSxI9mD2d6ATwKj0gSoScCjsUB1BvC7TF5XAaryZPqBm6zcZZdf4U1NTQULFlfNneejzhp3KEjd+9wffNRZ4/yqufO6/d1MA24x6ySFkU2A6lCov49yVPQA1f76DE0ToP4Z+FLc/Q3AgO5eUwGq8mT6gVuKD+Zcvu0WsoUT1lZd1OUSoOS/ZBOgijUGNQiIn9XVHHusCzObYWZrzGxNS0tLUSonqQU59pFp2m3Q6bmZnkOPHj1YctutNG95iydXPUrzlrdYctut9OiRerpgISf+FioLTKRcFCtAJcuJ9GQF3X2Zuze6e2Pfvn0LXC1JpRAz4DP9wA3qgznXc6iurmb48OHdDmIXap5Lh0JlgYmUi2IFqGZgSNz9wYC+/oVYIVoGmX7gBvXBXOhljQrdwilUFphIuShWgFoJXGLtzgDec/dtRTq2ZKijK+ztt98uSMsg0w/cID6YC926geK0cFKlxd+yeFHery0SdoEEKDP7MfA8cIKZNZvZdDObaWYzY0UeAV6jPc38X4BZQRxXgpHYFTbshBM4rGdPauuP7lQuiJZBph+4+X4wF6J1kziWVYwWTi7jYhJNpd5u46KLLuKEE07glFNOYdq0aRw4cCDv1wwsi68QN2XxFUeyjLmPjvqET5x6ScGyx7KZB5VLem6QGXDpMvU0zyV6wprFd9111/l3v/vdrH/v9ddf95NPPjmr3zl48KB/+OGHnR775S9/6QcPHvSDBw/61KlTfenSpUl/N4xZfBJSqbrCrr59Gc88/H/Z9sbrQPAtg0wTETItl+z3gmrdpBvLStXC2bdv36HWllaBiK762lrMrMutvrY2r9cN43YbkyZNOnR+o0ePprm5Oa9zBNSCirqmpibvO2BQp1ZGx63vgEFeXVNbti2DIFo32bbEEltbR/5dL/+76hofdOxxWgWiQmTTggLap5sm3No/enPzyiuv+PDhw72lpcXd3Xfu3OnunVtQ48aN89WrV7u7e0tLix977LHu7n755Zf7j370I3d3/9vf/uYffPBBlxbUqlWr/Gtf+9qhVtJnP/tZf/bZZ/311193M/Pnn38+bf3279/vp512mj/33HNJn8+mBaWO7Ihbevc9hwb648drdu/Yzv6/fsCmjRvYs2dPWe6N1NG6ufGG63Pe3ymTsazhw4cfejzZ/j93LriSoSeezOQvfz3vvYByof2tKku+223cdNNNNDc384UvfIFhw4Z1KRO/3Qa0//1s2rSJhoaGjLbbmDVrFp/85Cc5++yzszir5NTFF2Gtra3cf//9jP/8F/newtmdusJunTODiy++mP79++fUxRYmuXYTQnaZeqm6S69cfCdP//QBqnpVF3V1dG1vX5ncw7vdxg033EBLSwtLlizJ9rSSUoCKsI7WwYVXfYs+Awcx++/Hc9nEscz53KfY/tYbfGPm10tdxZLLZiwrXWuruncdu3e8XdRVILS9fWUK63Yb9957L6tWreLHP/4xhx0WTGhRgIqwfv36seudFr7x6dNpWv0C7s5Jjadz078/zMG2NhoaGkpdxVDINOU9XWur9d3d1PXrX7RVIIoxD0xKI6zbbcycOZPt27czZswYRo4cyY033pj3uVr7+FQ4NTY2ekcmigRv7ryreeLXv2XOrUs7jZe8s7WZ//GPXyjqOEk5yGQsZ+68q3n2xTWHAkOyMahxoxsL/n+7ceNGzpl4Pnc89usuz8057yyeXPVop7Ezydy6des48cQTMypbX1vL7iQf7nU1NezasyfoqpWFZP9/ZvaSuzcmllWAiqjW1lYGDR5yaDC/w+4d27ni/LPY2txM7969S1jD8tTW1sb8BQtZsWIFNb3r2NmyncMOO5y6+qPZ+967TJs2jVsWLyr4RNt013fe5Ak0b3mrrMcVSymbACVdZROg1MUXUenGS/r0PaZkWzyXu/h5UU+tepQdb7/N21v/wlOPP1bUVSC0jp9UAqWZR1T8eEniN+xir5RdiWnQHZmDHUrRnXbL4kXMX7CQeZMnUNu7jj3v7j7UgpP8ZJJJJ11l22OnFlREheEbttKgC0vr+BVGVVUVO3fuzPrDNurcnZ07d1JVVZXx72gMKsLix0sSv2EX40MsWUJBsZIIRHJ14MABmpubO80vksxUVVUxePBgevbs2elxJUlEXLputFJ0sWkQX0Q6KEkiojLpRstnpYVcaTtzEemOAlSFC+tqAtrOXES6owBVwcK8mkAYkjSke9oqREpJAaqChb0bTduZh5cyLCUMFKAqWNi70ZQGHV5h7RqWaFGAKnPpumDKpRutFEkaklqYu4YlWgIJUGZ2npltMLPNZtblK5aZHWVmPzezP5pZk5l9JYjjRlmmXTDqRpNshb1rWKIj7wBlZocD3wfOB04CvmRmJyUUuwx41d1PBcYDt5nZEfkeO8oy7YJRN5pkK+xdwxIdQbSgRgOb3f01d98PPABMSSjjQI21L15VDewCNNqao1y6YNSNJpkql65hqXxBBKhBwJa4+82xx+LdBZwIbAX+BMx294PJXszMZpjZGjNb09LSEkD1Ko+6YKKtGKnf6hqWMAgiQCVb0jdx/aSJwFpgIDASuMvMapO9mLsvc/dGd2/s27dvANWrPOqCiaZipn6ra1jCIIgA1QwMibs/mPaWUryvAD/zdpuB14GPBnDsipDtN2J1wURTKVK/1TUspRREgFoNDDOz42KJD1OBlQll3gLOATCzY4ATgNcCOHZZy+cbsbpgoiXf1G+tCCHlKO8A5e5twOXAKmAd8KC7N5nZTDObGSv2v4AzzexPwFPAAnd/J99jl7t8vhGrCyZach131IoQUs603UaJaLsJyUaufy/ac0vKgbbbCBll4kk2chl31IoQUu7UH1QiHZl42954HfeD1PXrz5G9eikTT1K6ZfEi5i9YyLzJE7rsgJxMJl+Chg8fXoyqi+REAapEqqqqGH7CCcz7/Kep69uPve+9y1mTprBjyxvKxJOkOsYdb7zh+ox2QI6fjpDYLagvQVIO1MVXIvMXLOTgEVV8//Hf8v0nnueOXzzDmxvXcfiBvykTT9LKNPU7225BZfpJ2ChJogSUICHF0tbWxvwFC1mxYkWXbsGOjM+OMsuXL6e2rp49u3cxffr0TmVECklJEiESZIKEvvVKOplMR9DeTxJWClAlEMRSRZrfItlI1S2oTD8JMwWoEghiqSJ965UgaLqDhJnGoEokk7GBVDSGJUHR35KEgcagQiafpYr0rVeCooWHJcyUolNiHWMD2dD8FglSthOARYpFXXxlSmusSdBaW1szmgAsErRUXXwKUGUqnzEsEZEwUYCqUPrWKyLlTkkSFUo7nkpQUk361mRwKRUFKJGISzXpe9++fZoMLiWlwQqRiIuf9B2fcHP2uPEcPKKqy+PzFyxUIo4UhcagRCIs3UTdWeeewZKHn2bA0OM6Pa4JvBI0jUGJSBfpJn33qu2N+8Euj2syuBRLIAHKzM4zsw1mttnMki4GZ2bjzWytmTWZ2bNBHFdE8pNu4eL397yL2WFdHtdkcCmWvAOUmR0OfB84HzgJ+JKZnZRQpjewFJjs7icDX8z3uCKSv3RLHY04dSQrvnONlkCSkgmiBTUa2Ozur7n7fuABYEpCmQuBn7n7WwDuviOA45YNpelKmN2yeBHjRjcyb/IE5px3FvMmT2Dc6Eb+89lnkj6uJZCkWPJOkjCzC4Dz3P2rsfsXA6e7++VxZW4HegInAzXAHe7+wxSvNwOYAdDQ0PDxN998M6/6FVv8xNmqqirtVCplI9Wkb00Gl0JLlSQRxKekJXksMer1AD4OnAMcCTxvZi+4+8Yuv+i+DFgG7Vl8AdSvKJJtmz38hBP4sOd/U5qulIVUCxfnsqCxSBCC6OJrBobE3R8MJKb4NAOPufv77v4O8BxwagDHDo1kGwi2th2k35Ch2qlUIkld25KvIALUamCYmR1nZkcAU4GVCWUeBs42sx5m9nfA6cC6AI4dCqm2zb769mX85tGV/PX99w+VVZquVLpUK1NoBQrJVt4Byt3bgMuBVbQHnQfdvcnMZprZzFiZdcBjwMvAi8C97v5KvscOi/RzSY5i9463Dz2mNF2pdMl6E559cQ3zFySdgSKSklaSCEC62fiXTTyT2/7jSQYMPU57NknF0xbykgutJFFA6eaSfOxjI7h26t8rTVciIV1vgrq2JVvKdQ5Ium2z9+3bpzRdiYT4lSkSW1Dq2pZsqYsvYJozIlE3d97VPPvimkNJQ+ralu5oR10RKYqOOYErVqzo0pugCeqSjAKUiBRVfG8CoJ4FSUlJEiJSVNXV1Rx//PF8+7rrNSdKcqIAJSIFozlRkg918YlIQWhOlGRKXXwiUlSaEyX5UoASkYJIt1uv5kRJJhSgRCRn6VYsT7fCinbllUwoQElZq6+txcy63Opra0tdtYqW6YrlqXbr1XJfkgklSUhZM7Muu2NC+y6aYf7bLnfZrhahFVYkHU3UzZHeWOHWEaDqgd1Jnq+rqWHXnj1FrlVlU3aeBE1ZfFnSpmvlZTfgSW679+4tZbUqkrLzpFgUoFJINcHwqrnzkg4Ka3triQpl50mxKEAlkWoL91k3386yZf/Mp86deKhFtW/fPrW0JFKUnSfFoqWFk0jXhVHXrz/f+uf7ObJXNUuvmcPZ48Zz8IiqQ/3xHW/U+QsWamuBIqirqcHUjVd06fY/EwmKkiSSSDcIPOdzn+Kep1dzZK9ebHvjda6aMoG7n3hBg8VFVF9bm3RsqQ7YlfCYsvkKS0lEEoSCJkmY2XlmtsHMNptZylUgzewTZvahmV0QxHELJVUXxvcWzmbCP07lyF69AHA/SHVtbw0WF9nuvXuTJ0TQHpDib3U1NaWqZiRUV1czfPhwBScpiLwDlJkdDnwfOB84CfiSmZ2UotxiYFW+xyyG+AmGs88by6xzz6B/w1D+ae41h8qYHUbrnnc1WBwi7t7pphRzkfIVRAtqNLDZ3V9z9/3AA8CUJOWuAH4K7AjgmAXXo0cPltx2K81b3uKpVY/xta/NoKX5Tfbs2gm0B6EV37mGU08dqcFiEZECCCJJYhCwJe5+M3B6fAEzGwT8AzAB+ES6FzOzGcAMgIaGhgCql5/q6moGDhzIrG/MhLvv6TIofPNN3+Gaa/+nBotFRAIWRICyJI8ljkrfDixw9w/NkhWP+0X3ZcAyaE+SCKB+OWtra2P+goUsX76c2rp69uzexSWXXMI3Zn6dhoaGQy2kJbfdyo03XK/BYhGUOCHBCaKLrxkYEnd/MJCYHdAIPGBmbwAXAEvN7PMBHLugkk3W/e0f/si9y1d0eeNpsLjwOhaG7UHXZAglRJSeVl+RoOWdZm5mPYCNwDnAX4DVwIXu3pSi/L8Cv3D3h7p77VKuxaf1xsJHC8OGW7YLyIp0KFiaubu3AZfTnp23DnjQ3ZvMbKaZzcz39UtF642JZC7d6isrVqzQEmCSk0BWknD3R4BHEh67J0XZLwdxzEKLX28ssQWlFHKRzjL5Qjd8+PAS1U7KldbiS0HrjYlkTgvISiEoQKWh3UBFMqMvdFIIWosvA0qbDYeUa/BpU8JQ6JiWsWLFii5zAnv00LrUkpp21BWRotAXOslWqgClrzUiEqiOOYEi+YrcGJR2vhUpPr3vJBeRCVCa5S5SfJ3fd+cxYNAgrrhytt53kpHIBKhkyxY9++Ia5i9IuX2ViOSp8/vuN9z+i2d58jfPM2bsWYeClFpXkkokkiS0bJFI8aV731028UymfWUaVVVVnRZjnj59urL+IqigO+qGnZYtqmwdi8gm3upra0tdtUhrf9/VJX3f9e7Tj3+971955ner1ashKUUiQGmWe2VLuQV8kjlTUjwDBw7k3V07k77vWt97l7YDbVz2v+/Q2n2SUiQClGa5ixRfdXU1X770y9w6Z0an9933Fs5m1Ccn0Kv2KPVqSFqR6ei9ZfEi5i9YmHTnW00sFCmM/7PkNsaMPYvLJp5J7z79eH/Pe4w9fzLb3vgzf/vgfS3GLGlFIkkiXnwwqqqq6rJjrgZpy4/2iQq3trY25lw1l/vuu4+j6upp3fMe06ZN46Af5D9X/177R4mWOkpGG6xVBgWo8pDYU6G1+6SDAlQCpZ5XDi0iW97UxS6RTjNPRqnnlWPXnj24e5ebglN56Fi7T8FJEkU2QCn1XEQk3CIboJR6LiISboEEKDM7z8w2mNlmM+syDdzMLjKzl2O335rZqUEcN1/d7ZirNcJEREon7yQJMzsc2AicCzQDq4EvufurcWXOBNa5+24zOx+43t1P7+61i7VhYarsIqWfi4gUXiGTJEYDm939NXffDzwATIkv4O6/dffdsbsvAIMDOG5GMmkFJQ7SauVzEZHSCyJADQK2xN1vjj2WynTg0VRPmtkMM1tjZmtaWlpyrlSu+z+1trayfPnyQ3OjQGuEiYiUQhABypI8lrTf0Mw+RXuAWpDqxdx9mbs3untj3759c65Urq0gpZ+LiIRDEAGqGRgSd38w0OVT3MxGAPcCU9x9ZwDHTSmfVpDSz0VEwiGIALUaGGZmx5nZEcBUYGV8ATNrAH4GXOzuGwM4Zlr5tIKUfi4iEg55p6S5e5uZXQ6sAg4HVrh7k5nNjD1/D/Bt4GhgqZkBtCXL2AhKfCsol5WS0618LuGUbrkjSL43lJZCEgm3il2LL4iFYLVGWPlIt2AsJB8U1WKyIuEQucVitVJytChAiZSvyAWoDmoFRYMClEj5ShWgKr4p0TEJV0REyktkF4sVEZFwU4CSilBXU4NBl1tdTU3a50QkvCq+i0+iQeni5UljxJKOWlAiUnS5rpUp0aIAJSJFN+equax69jlu/skvtWOApKQAJSJF09bWxhVXzmbZvyxj144dLPjvk7jvlhuprT9aOwZIFwpQIlI08xcs5NcvreXuJ17g+088zx2/eIY317/Kj5bcrB0DpAsFKBEpio5dBi5fdEenXQauWHQHT//0Aba98bp2DJBOFKBEpCjS7TLQq/Yoll6b2Y4BmeySLZVBAUpEiiLdXmu7d2znU2POSLtjgDL/okcBSkSKItVea3ctnM2MGTO463t3pl3IOdddsqV8VfxisSISHrnuMtDa2sqgwUNY8vNfddnjbd7kCTRveUsTfctYZBeLFZHw6NGjB0tuu5Ubb7g+qxUkMtklW4tCV55Id/FpsFWkNDp2Gci01ZNu/EqZf5UrkgFKg60i5SXV+NXSazLL/JPyFEiAMrPzzGyDmW02sy4jltbuztjzL5vZqCCOmysNtoqUn1sWL2Lc6EbmTZ7AnPPOYt7kCYwb3Zg28086q6+txcy63Oprawv6u7nKO0nCzA4HNgLnAs3AauBL7v5qXJlJwBXAJOB04A53P7271y5EkkS6wdarPjeejevX079//0CPKeFTX1vL7r17uzxeV1OjldFDTiug5y7dztPdxYJ8fjeDeiVNkgiiBTUa2Ozur7n7fuABYEpCmSnAD73dC0BvMxsQwLGzlm6wtccRVXxk2HB191WgxG9/u/fuxaHLLVnQknDJdvxKylcQAWoQsCXufnPssWzLFEW6wdYD+//G4oceU3dfBUoMSCISfkEEKEvyWOJnQCZl2guazTCzNWa2pqWlJe/KJUo12Pq9hbOZ8I9TGTD0OK2qLCISAkEEqGZgSNz9wUDicsSZlAHA3Ze5e6O7N/bt2zeA6nXVMdh61efG89WzT2PO5z7FsR89iX+aew2AVlUWEQmBIALUamCYmR1nZkcAU4GVCWVWApfEsvnOAN5z920BHDsnHZMFN65fz/6/fsCin/ySS+d/m8NjM9k1t0JEKlFdTQ0GXW49odsMvVS/W1dTU7D65h2g3L0NuBxYBawDHnT3JjObaWYzY8UeAV4DNgP/AszK97hB6N+/PzNmzGDFd67R3IqIqaPrG63QbzaRUtu1Zw/u3uV2gK4JQ4lJQ6l+t5BZr5Ffiy/XtcGkvCitXCS1QqaQZ3j8pGnmkQ9QHTS3QkSiKqwBSk2EmI65FSIiEg6RXItPKl8plmURkWApQElF0koRIpkrRYZeJtTFJyIScWFNFFILSkREQkkBSkREQkkBSiQDSroQKT4FKKlIQQ/6KulCpPiUJCEVKayDviKSObWgREQklBSgREQklBSgREQklBSgRDIQ1pn2IpVMAUokA6XYC0ckH5UwNUIBSqQblfBGl+iphKkRSjMX6UbHGz2RldEbXaQcqQUlIiKhpAAlkid1AYoURl4ByszqzewJM9sU+7cuSZkhZvYrM1tnZk1mNjufY4qETbq+fgUpKTdh+sKVbwtqIfCUuw8DnordT9QGzHP3E4EzgMvM7KQ8jytSFsppQFoqS65TI8KUXJFvgJoC3Bf7+T7g84kF3H2bu/8+9vNeYB0wKM/jihSN5kCJlEa+WXzHuPs2aA9EZtYvXWEzGwqcBvwuz+OKFI3mOkk5qoTs024DlJk9CfRP8tS12RzIzKqBnwJz3D3lO97MZgAzABoaGrI5hIiIVJBuA5S7fzrVc2a23cwGxFpPA4AdKcr1pD04/Zu7/6yb4y0DlgE0NjYm+wIgEip1NTVJv5XWAbuLXx2JqPra2i7jREb73+GuktQof/mOQa0ELo39fCnwcGIBMzNgObDO3ZfkeTyR0Nm1Z0/S8ajdaJxKiidlckOWrxOmMdd8A9Qi4Fwz2wScG7uPmQ00s0diZcYCFwMTzGxt7DYpz+OKhIrW6pNKEaa/5bySJNx9J3BOkse3ApNiP/+a9gAsEjnJul2g/duogpcUS/wHcDm16rUWn0gBVUImlZQ/9/IcztdSRyIiEkoKUCIiFSBMyQ1BURefiEgFqMQxTbWgREQklBSgRAqoErtdJBhHJFkx3Mw4wpT03EEBSqSAcplTEqbtDqRwDtB1Uq3HHpd2GoMSKaFU86SUmi6iACVSUsnmSamDR6SduvhERCSUFKBERCSU1MUnIlICPUnenduz2BUJMQUokZCpI/kHl7u4z3oAAAgdSURBVFLTK8v+Ml0fr5gUoERKKOVmh1rtXEQBSqSUFIREUlOShIiIhJIClIhICGgFka4UoEREslSIYNIxaTvxlmylkajQGJSISJa0U3Jx5NWCMrN6M3vCzDbF/q1LU/ZwM/uDmf0in2OKiEg05NvFtxB4yt2HAU/F7qcyG1iX5/FEJCbXbiaNdUi5yDdATQHui/18H/D5ZIXMbDDwWeDePI8nIjG5jllorEPKRb4B6hh33wYQ+7dfinK3A/OBg3keT0SkImlzy666TZIwsyeB/kmeujaTA5jZ3wM73P0lMxufQfkZwAyAhoaGTA4hIlJU6VYAyZUmbXfVbYBy90+nes7MtpvZAHffZmYDgB1Jio0FJpvZJKAKqDWzH7n7P6U43jJgGUBjY6MWq5JISrWRoZZASq9Y/2/ZvpauZ27y7eJbCVwa+/lS4OHEAu7+LXcf7O5DganA06mCk4i00zhRbsL6/xbWeoVdvgFqEXCumW0Czo3dx8wGmtkj+VZORFLLdcxCYx1SLvKaqOvuO4Fzkjy+FZiU5PFngGfyOaaItMu1a0hdSlIutNSRiIiEkgKUiIiEkgKUSAiVepyoXFebKPX/W7nVK+y0WKxICOU6ThRUOnOyxVDrY4+bdd6QPkyp0mGpR6Kw1ivs1IISqSCFTGfeneR1w5oqXa4tQOlMLSgRqTjaDqMyqAUlIiKhpAAlIiWl7jhJRQFKRLpIlnVWKMVeBig+8Ck4hpsClEgFCSqdedeePbh7p1vQqdIdwYGE16vP6dUyFx/4tEZeuClJQqSCFDKdOejXTpnIEMBrp9wOg/ZsRCkPakGJSEGUsvusI5gmtox2pasvnbsz1d1XempBiUhBlFuqd8c8r0RhrW8UqAUlIqGjZYAEFKBEJITcPe8xr0ySOuLLSPioi09ESiJlIkNALadMAlx8mcQ1BqX0FKBEpCS0gKp0R118IlIQ5bbFRLnVNwrUghKRgii3FlK51TcK8mpBmVm9mT1hZpti/9alKNfbzB4ys/Vmts7MxuRzXBERqXz5dvEtBJ5y92HAU7H7ydwBPObuHwVOBdbleVwREalw+QaoKcB9sZ/vAz6fWMDMaoFPAssB3H2/u7+b53FFRKTC5RugjnH3bQCxf/slKXM80AL8wMz+YGb3mlmvPI8rIiIVrtsAZWZPmtkrSW5TMjxGD2AUcLe7nwa8T+quQMxshpmtMbM1LS0tGR5CREQqTbdZfO7+6VTPmdl2Mxvg7tvMbACwI0mxZqDZ3X8Xu/8QaQKUuy8DlgE0NjYmWxpLREQiIN8uvpXApbGfLwUeTizg7m8DW8zshNhD5wCv5nlcERGpcPkGqEXAuWa2CTg3dh8zG2hmj8SVuwL4NzN7GRgJ3JzncUVEpMLlNVHX3XfS3iJKfHwrMCnu/lqgMZ9jiYhItJh7eId5zKwFeLNAL98HeKdArx0mUTjPKJwjROM8o3COEI3zzOYcj3X3vokPhjpAFZKZrXH3im/VReE8o3COEI3zjMI5QjTOM4hz1GKxIiISSgpQIiISSlEOUMtKXYEiicJ5RuEcIRrnGYVzhGicZ97nGNkxKBERCbcot6BERCTEFKBERCSUIhOgzOyLZtZkZgfNLGXqo5mdZ2YbzGyzmaVcMzCssthE8g0z+5OZrTWzNcWuZy66uzbW7s7Y8y+b2ahS1DMfGZzjeDN7L3bd1prZt0tRz3yY2Qoz22Fmr6R4vuyvI2R0npVwLYeY2a9iG9E2mdnsJGVyv57uHokbcCJwAvAM0JiizOHAn2nfIuQI4I/ASaWue5bneQuwMPbzQmBxinJvAH1KXd8szqvba0P76iWPAgacAfyu1PUuwDmOB35R6rrmeZ6fpH2Hg1dSPF/W1zGL86yEazkAGBX7uQbYGOT7MjItKHdf5+4buik2Gtjs7q+5+37gAdo3ZSwn3W4iWaYyuTZTgB96uxeA3rFV9stFJfz9dcvdnwN2pSlS7tcRyOg8y567b3P338d+3kv7bumDEorlfD0jE6AyNAjYEne/ma7/2WGXySaSAA48bmYvmdmMotUud5lcm3K/fpnWf4yZ/dHMHjWzk4tTtaIq9+uYjYq5lmY2FDgN+F3CUzlfz7wWiw0bM3sS6J/kqWvdvctWIMleIsljocvDT3eeWbzMWHffamb9gCfMbH3sG19YZXJtyuL6pZFJ/X9P+7plrWY2Cfh/wLCC16y4yv06ZqpirqWZVQM/Bea4+57Ep5P8SkbXs6IClKfZXDFDzcCQuPuDga15vmbg0p1nhptI4u0rzuPuO8zsP2jvXgpzgMrk2pTF9Uuj2/rHv/nd/REzW2pmfdy9khYeLffrmJFKuZZm1pP24PRv7v6zJEVyvp7q4utsNTDMzI4zsyOAqbRvylhOut1E0sx6mVlNx8/AZ4CkmUYhksm1WQlcEssaOgN4r6O7s0x0e45m1t/MLPbzaNrfwzuLXtPCKvfrmJFKuJax+i8H1rn7khTFcr6eFdWCSsfM/gH4HtAX+KWZrXX3iWY2ELjX3Se5e5uZXQ6soj2jaoW7N5Ww2rlYBDxoZtOBt4AvQvsmksTOEzgG+I/Ye6MH8O/u/liJ6puRVNfGzGbGnr8HeIT2jKHNwAfAV0pV31xkeI4XAN8wszbgr8BUj6VKlQsz+zHtGWx9zKwZuA7oCZVxHTtkcJ5lfy2BscDFwJ/MbG3ssWuABsj/emqpIxERCSV18YmISCgpQImISCgpQImISCgpQImISCgpQImISCgpQImISCgpQImISCj9f424ovFy/xZgAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[y_db == 0, 0], X[y_db == 0, 1],\n",
    "            c='lightblue', marker='o', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 1')\n",
    "plt.scatter(X[y_db == 1, 0], X[y_db == 1, 1],\n",
    "            c='red', marker='s', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 2')\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"accelerating\"></a>\n",
    "## Accelerating Clustering with RAPIDS\n",
    "\n",
    "So how can we accelerate our code using RAPIDS? First, we cast our data to a `pandas.DataFrame` and use that to create a `cudf.DataFrame`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cuDF Version: 0.12.0b+1777.g162097c\n",
      "Pandas Version: 0.24.2\n"
     ]
    }
   ],
   "source": [
    "import cudf; print('cuDF Version:', cudf.__version__)\n",
    "import pandas as pd; print('Pandas Version:', pd.__version__)\n",
    "\n",
    "\n",
    "X_df = pd.DataFrame({'fea%d'%i: X[:, i] for i in range(X.shape[1])})\n",
    "X_gdf = cudf.DataFrame.from_pandas(X_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"gpukmeans\"></a>\n",
    "## Clustering using GPU Accelerated K-Means\n",
    "\n",
    "Next, we load the `KMeans` class from the `cuml` package and instantiate it in the same way we did with the `sklearn.cluster.KMeans` class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cuml import KMeans as KMeans_GPU\n",
    "\n",
    "\n",
    "km_gpu = KMeans_GPU(n_clusters=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `KMeans` class from `cuml` implements the same API as the `sklearn` version; we can use the `fit` and `fit_predict` methods to fit our `KMeans` model to the data and generate predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_km_gpu = km_gpu.fit_predict(X_gdf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, let's visualize our results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dfXxU5Z338c9Po4tLEjLIg4aAYFesD8WHRhSxVaEWtC3utvZ++VCrJd5IaRUklES7d7ftrixQSMFWRG6g1e1uXe+2e1etis9arS2NW6oij/VWSGMJQmIYKpXR6/4jk3SSzEzm4czMOTPf9+s1LzInJ3OukxPmN9e5ftfvMuccIiIifnNEoRsgIiISjwKUiIj4kgKUiIj4kgKUiIj4kgKUiIj4UlmhG5DMsGHD3NixYwvdDBERyaGXXnrpbefc8L7bfR2gxo4dS3Nzc6GbISIiOWRmb8bbrlt8IiLiSwpQIiLiSwpQIiLiS74egxIR8ZvDhw/T0tLCoUOHCt2UwBk0aBA1NTUcddRRKe2vACUikoaWlhYqKioYO3YsZlbo5gSGc459+/bR0tLCuHHjUvoZ3eITEUnDoUOHOPbYYxWc0mRmHHvssWn1PBWgRETSpOCUmXR/bwpQIiLiSwpQIiJF4Jvf/CbLli1L++c6OjpYtWpV1sf//ve/z9/93d9hZrz99ttZvx4oQImI5Fw4HGb79u2Ew+FCN6WfTAKUc44PPvig17bJkyfzxBNPcMIJJ3jWNgUoKWpDKysxs36PoZWVhW6alIBIJML8+gWMqhnN1GmXMqpmNPPrFxCJRLJ63XvvvZcJEyZwxhlncO211/b7/kUXXdRTJu7tt9+mu6bp5s2bmThxImeeeSYTJkxgx44dNDY28oc//IEzzzyTr33tawB85zvf4ZxzzmHChAn80z/9EwBvvPEGp5xyCnPmzOHss89m9+7dvY551lln4XXtVKWZi+8Mrayk/cCBfttDFRXs7+xM67XaDxzAxdlucV5fxGsLGxp5dmMzTQ8+TWjESNrb9rDqtnksbGikaXn6t+OgK8jcfvvtvPDCCwwbNoz9+/en/LOrV69m7ty5XHPNNbz33nu8//77LF68mFdffZVNmzYB8Nhjj7Fjxw42btyIc44ZM2bw3HPPMWbMGLZt28YPfvADT24JpkI9KPGd7qDS9xEvaIn4VTgcZt26dcxZtILQiJEAhEaMZM6iFaxfvz7j231PPfUUV1xxBcOGDQNg6NChKf/spEmTWLRoEUuWLOHNN9/kmGOO6bfPY489xmOPPcZZZ53F2WefzdatW9mxYwcAJ5xwAuedd15G7c6EApSULN3+k1xqbW2lMjS0Jzh1C40YSWVViNbW1oxe1zk3YLp2WVlZzxhR7Lyjq6++mgceeIBjjjmGadOm8dRTT8V9/VtvvZVNmzaxadMmdu7cSV1dHQCDBw/OqM2ZUoCSkqWemuRSdXU1ne37aW/b02t7e9seOjvaqa6uzuh1p06dyv3338++ffsA4t7iGzt2LC+99BIAP/nJT3q2v/7665x44oncfPPNzJgxg5dffpmKigoOxPzNT5s2rVcP749//CNtbW0ZtTVbClAiIjlQXl5OXV0dq26b1xOkusegZs6cSXl5eUave9ppp/H1r3+dCy+8kDPOOIP58+f322fBggXcddddnH/++b1Svv/zP/+T008/nTPPPJOtW7fyxS9+kWOPPZbJkydz+umn87WvfY1PfvKTXH311UyaNImPfOQjXHHFFb0CWCJ33HEHNTU1tLS0MGHCBG644YaMzi+WORdvCNkfamtrnRYsLD1mFj+xga7bD+lIlnCRMIEig+NI6diyZQunnHJKSvtGIhEWNjSyfv16KqtCdHa0M3PmTJYuWUxZWWnmqMX7/ZnZS8652r77etKDMrP1ZtZmZq8m+L6Z2R1mttPMXjazs704rhSnUEUFBv0eoYqKtF9rf2cnzrl+j3SzAXNN42HFqaysjKbly2jZvYsnNjxCy+5dNC1fVrLBKV1e3eL7ITA9yfcvBU6KPmYBd3l0XAmAVN98u/fr2+MJVVT4Mqh4SeNhxa28vJzx48dnfFuvVHkSoJxzzwHJkvEvB+51XX4NVJnZ8V4cW/wv1TfffL9JZ9JTU09HJH/y1c8cBcROO26Jbnur745mNouuXhZjxozJS+OkNGXSI9PEX5H8yVcWX7yk/bij0M65Nc65Wudc7fDhw3PcLOkr6D2EoLdfRP4qXz2oFmB0zPMaILNZapJTQe8hBL39IvJX+epBPQB8MZrNdx7wjnOu3+09kVLlZeailKZCL7dxzTXXcPLJJ3P66aczc+ZMDh8+nPVrepVm/mPgReBkM2sxszozm21ms6O7PAy8DuwE/jcwx4vjSnbi3Q4DSL2yV2pSffMt5TfpoKTDS/HxarmNa665hq1bt/LKK6/w7rvvsnbt2qzb5lUW31XOueOdc0c552qcc+ucc6udc6uj33fOua845z7knPuIc06zb30gYdacx8dJ9c3XL2/SycaxSjmISvpyNSbqx+U2Lrvssp7zmzhxIi0tLVmdI2i5DZF+ko1jqcKEpCMXY6J+X27j8OHD/Nu//RsrV67M+By7KUBJXH3TLoPSQwhVVMT9z5+L9icqo1QGxC5Hl8k6ViKJZLvcxu23305LSwuf/exnOemkk/rtE7vcBnQtG7Jjxw7GjBmT0nIbc+bM4eMf/zgf+9jH0jir+BSgJK6g9hTyGQiS1vKLfa4MQvFQtsttnHvuufziF79g2rRprF27lhNPPLHf6996663ceOONvba/8cYbAy638a1vfYu9e/dy9913p3NKCamauUjAae5XafHrchtr165lw4YN/PjHP+aII7wJLQpQJayM/oP9hrrVQaM6fqXFr8ttzJ49mz179jBp0iTOPPNMvv3tb2d9rlpuowQlXIKCroKKpb7cRLIlOmJvISZdFqTv8xz+Pr1cnkQGls5yG6n+LZWSdJbb0IflEpRs7ETyO44lxU1/S9nRLT7pRXN6UpdoTlTfW6f6fYpkRj0o6UW3hFKnT8elK5VMOukv3fcX9aBEAk7VLfJr0KBB7Nu3Tx/m0uScY9++fQwaNCjln1EPSvJKg8be0+8tv2pqamhpaWHv3r2FbkrgDBo0iJqampT3V4AqQfmsttCXlsOQoDvqqKMYN25coZtREhSgSpA+cYtIEGgMqoipwoCIBJl6UEVMt9NEJMjUgxIRwuEw27dv76m/JuIHClCSV0qJ9pdIJML8+gWMqhnN1GmXMqpmNPPrFxCJRAb+YZEcU4CSvPLLqrnSZWFDI89ubKbpwadZ+ejzND34NM9ubGZhQ2OhmybiTYAys+lmts3MdppZv79sMxtiZg+a2e/NbLOZfcmL44pI5sLhMOvWrWPOohWERowEIDRiJHMWrei13IJIoWQdoMzsSOBO4FLgVOAqMzu1z25fAV5zzp0BXAQsN7Ojsz12qUo1O0+30ySZ1tZWKkNDe4JTt9CIkVRWhWhtbS1Qy0S6eNGDmgjsdM697px7D7gPuLzPPg6osK7iVeV0reqgm9wZSnX9H91Ok2Sqq6vpbN9Pe9ueXtvb2/bQ2dFOdXV1gVom0sWLADUK2B3zvCW6Ldb3gVOAVuAVYK5z7gMPji0iGSovL6euro5Vt83rCVLtbXtYdds8Zs6cSXl5eYFbKKXOi3lQ8Ur69p1+Mw3YBEwBPgQ8bma/dM71+yhvZrOAWQBjxozxoHkiksjSJYtZ2NBI/YwpVFaF6OxoZ+bMmSxdsrjQTRPJfkVdM5sEfNM5Ny36/FYA59y/xuzzC2Cxc+6X0edPAY3OuY3JXlsr6sanFVSLWyEK6obDYVpbW6murlbPSfIu0Yq6Xtzi+y1wkpmNiyY+XAk80GefXcDUaENGAicDr3twbJGik+oYo5fKy8sZP368gpP4StYByjkXAb4KbAC2APc75zab2Wwzmx3d7Z+B883sFeBJoME593a2xw6ybOrkKTtPREpB1rf4cqmYb/HpNp0kksnfhtbZkiBLdItPxWJFioAKA0sxUqkjERHxJQUoEZ/RGKNIF93iK5Ay4k8g0wURjRmJdFEPqkAi9E8jdqj+k+Se1n6SoFCAEikCqdwW1NpPEjQKUCJFIJXCwFr7SYJG86AKJJt5UJrzIukKh8OMqhlN04NP91peo71tD/UzptCye5eqSEjB5LLUkeRZIUrhSLBp7ScJIgWoAlEqseST1n6SIFKAKhAtJij5pLWfJIg07UakRGjtJwkaJUkEkArNSja09pP4jYrFFpFQRUXcIqAav5JUdK/9JOJ3ClABpHEqESkFSpIQKTGJSh2pBJL4jQKUSIlIVOro0KFDKoEkvqRbfCIlIrbUUWjEyJ40849deBEfHD2o3/aFDY00LV9W6GZLCVMWn0gJSFbqaM4l59H086c4fuy4XttVAknyJaeljsxsupltM7OdZha38qSZXWRmm8xss5k968VxRSQ1yUodDa6swrkP+m1XCSQptKwDlJkdCdwJXAqcClxlZqf22acKWAXMcM6dBnw+2+OKSOqSlTo62NmB2RH9tqsEkhSaFz2oicBO59zrzrn3gPuAy/vsczXwM+fcLgDnXJsHxxWRFCUrdTThjDNZ/y+3qQSS+I4XAWoUsDvmeUt0W6zxQMjMnjGzl8zsix4cV0TSsHTJYi6cWEv9jCnMm34B9TOmcOHEWn757DNxt6sEkhRa1kkSZvZ5YJpz7obo82uBic65m2L2+T5QC0wFjgFeBD7lnNse5/VmAbMAxowZ89E333wzq/YVktZtEj9KVOpIJZCkUHJZ6qgFGB3zvAboO7LaArztnDsIHDSz54AzgH4Byjm3BlgDXVl8HrQvLxIGI2B/n23xyhSJ5EuiUkcqgSR+48Utvt8CJ5nZODM7GrgSeKDPPj8HPmZmZWb2t8C5wBYPju0bCRcRLGirRHJL1Sckl7IOUM65CPBVYANdQed+59xmM5ttZrOj+2wBHgVeBjYCa51zr2Z7bBEpjERVKVR9QrykiboeSboERrxtPv69iwxkfv0Cnt3YzJxFK3pVn7hwYq2qT0jaEo1BKUB5RAFKSkWyqhSqPiGZyGklCUnO+jy0bpMEWbKqFKo+IV5SgPJIqKKiXyDqDkbOuV4PpZhLkCWrSqHqE+IlBSiP7O/s7BeIFIykGCWrSqHqE+IlLbchImlbumQxCxsaqZ8xhcqqEJ0d7Vx77bXcUDeTcDisICWeUA9KRNJWVlZG0/JltOzexaO/eJAvfOEL3HvvvUz71GeUci6eUYASkYyVl5ezdt16Xtz0Mk0PPs3KR5+n6cGneXZjMwsb4q68I5IypZmLSMaUci5eUJq5iHhOKeeSSwpQIpIxpZxLLilAiUjGlHIuuaQAJYEytLISM+v3GFpZWeimFa2BKpYnWghRCx5KtpQkIYGStOahj/+WgygSibCwoZF169ZRGRpKZ/t+6urqWLpkMWVl/adQasFDyZSSJKToDOWvJaUA9aY8trChkWc3NqecPt694KGCU3HL5xpgClASWO30XyDSQdyVjSU94XCYdevW9SynAV2ZeXMWrWD9+vVaoLAEFWINMAWoJDTeIaVK6ePSV7o9ai8oQCWRcBn3Awd6gtXRcQKYgpgEndLHJVahetQKUBnqDlaH0W2mfIpd1kRyR+njEqtQPWpVM5dAiV2+xExhKpfiVSyfOXOm0sdLUGyPum9Jq7f37mHEiBE5Oa4nPSgzm25m28xsp5klvCFpZueY2ftmdoUXx5XSEW88EOAotGJxrsRWLH9iwyO07N5F0/JlcVPMpbh196hXLJjTq0d9R8PNhIaN4Nv//C85OW7W86DM7EhgO3AJ0AL8FrjKOfdanP0eBw4B651zPxnotQs9DyrpnJs4X/fbR/NyPKP5TyKF1dHRwajRY7AjjqAiNJRwRztTPncln/pCHQs/+8msCgMnmgflxUehicBO59zr0QPdB1wOvNZnv5uAnwLneHDMvAhVVGBxxpJCBWiLiEghtbW1MXTYcBb/dAPtbX8iNOI4jhk8GKBnHGr8+PGeHtOLW3yjgN0xz1ui23qY2SjgH4DVHhwvb2KXcY+9bdTOX28nxbvFpNtMIlJsusehDh0MUz3uQz3BKZeZnV4EqHgj1X3vuawAGpxz7w/4YmazzKzZzJr37t3rQfO8ERusYh/vxdnmnOs1mC9SbPJZTUD8oRCZnV4EqBZgdMzzGqBvzmEtcJ+ZvQFcAawys7+P92LOuTXOuVrnXO3w4cM9aF7mNFFXpLdCVBMQ/8h3YWAvkiTK6EqSmAr8ka4kiaudc5sT7P9D4KHAJ0loYD4vhlZW0n7gAGVAvLfAUEWFeqt5NL9+Ac9ubO6ZsNn9CfrCibU0LV9W6OZJllIt+Ot1YeCcFYt1zkWArwIbgC3A/c65zWY228xmZ/v6Utq6q3n0nRAN6FZqnqk+X/FKt2ecr8LAnkxocM49DDzcZ1vchAjn3PVeHFNE8iuVagJeZ3FJfsTW2YvtGS9saCxoz1iljkQkJarPV5xS7RkXIjFGAUpEUqL6fMVpoJ7xrl27CpYYowCVRGxhUs1xEtHy7sVooJ7xqrtW532ZjW5a8l18rTuLry9l7xWWlncvLomyM88/6wzuvffenrGpbu1te6ifMSWr8kaxtOS7BFKiCdIKToWl5d2LS6Ke8Zdn31jQhStVllhEpMR1V67/9re+2atnHA6HEy6zkY/EmJLqQakyhEjuqPxR8PXtGRc6MaakAlSyJdxFJDO9J3lO5/hRo7jp5rkqf1QkCpkYU1JJEipdJOK9eAPsy+bNorzsCF584XkOHTqkhIoikMvEmERJEgpQKECJZCocDjOqZnTcLK+vTDufU089jT/s3EllaCid7fupq6tj6ZLFWpXXhwqZmaksPikaGkv0j65JnqG4WV7lQ0K88+5fCjJ/RlLn5wr1ClASOBpL9I/q6mo69u+LP8mzfR/zm1arsKzPxdbh89sHiZIKUKoMIeKt8vJyrr/uepbNm9Ury+uOhpv528EVHD92XK/98zV/RlLj9wr1JRWg4k36DFVU0H7ggG4XiWTou03LKS87gq9MO585l0xi3mcuZvio0bx7MKzCsj6XSoX6QiqpABWPbheJZKesrIwXX3ieG+pu4GBHO0OGVNH8xCNMOOMMFZb1uerqat7Zv48tL23k3YMHe7b75YNESWXxxaPMvuDRNfOv2EywQYMGsbChkfXr11NZFaKzo52ZM2cqi88nIpEICxsaWb16NcdUVHDoz3/mE5+/mk99oY67v7Egr6skK808Ab3ZBY8KyAaLCsv6U7z5a8tvuZFd27fypeuv58uzb2TMmDF5uWYKUAkoQIlIqUk2f+2m6RdQVlbGkKHH5m3umuZBiYgIkDw54piKChruuscXKeeeBCgzm25m28xsp5n1OxMzu8bMXo4+fmVmZ3hxXC8o9VxESk2yRQr/8u67jP3waUDhU86zDlBmdiRwJ3ApcCpwlZmd2me3/wdc6JybAPwzsCbb43pF6w2JSKlJVKV82S03MvWKqzhm8OCefQuZcu7FTcWJwE7n3OsAZnYfcDnwWvcOzrlfxez/a6DGg+PmjAbhRaTYLV2ymIUNjdTPmEJlVYh3OvYTibzPLd+5s9d+hUw59+IW3yhgd8zzlui2ROqARxJ908xmmVmzmTXv3bvXg+alT3OjRKTYdS9S2LJ7F09seIQ/7t7NjbNmcfc3Fvhm7poXPSiLsy1u+puZXUxXgLog0Ys559YQvQVYW1ubVRqdekIiIsl1L1II/XtVsXPXCsGLANUCjI55XgP0u1lpZhOAtcClzrl9Hhx3QN09oX5tUU9IRKSfREu/F6w9HrzGb4GTzGwc8EfgSuDq2B3MbAzwM+Ba59x2D44pIiI5EturKqSsA5RzLmJmXwU2AEcC651zm81sdvT7q4FvAMcCq8wMIBJvUpaIiEg3T6YGO+ceBh7us211zNc3ADd4cax8CFVUxL0NqLlR/pJojPEo4HCc/TX2KBIsqtgYh97EgiHhGCPxs3Q09igSLEVd6khVIkREgquoe1DqCYmIBFdR96BERCS4FKBERMSXFKAksBKNMR4VZ5vGHgsnHA6zffv2glTDlmBTgJLASlSJ/r0421ShPv8ikQjz6xcwqmY0U6ddyqia0cyvX0AkEil00yQgFKBEJCcWNjTy7MZmmh582heL30nwKECJiOfC4TBr167ls1++hUGDu2q5FXrxOwkeBSgR8VQkEuHmuXP5y1/+wvdvu4UbL67lnqXf5v1IpKCL30nwKECJiKcWNjSyadtOVj3+IndueIGVDz3Dm1tf40dNiwq6+J0ET1FP1BWR/AqHw6xbt46mB58mNGIk0HVr76bFK5n76Yt4c8urBVv8ToJHPSgR8UxrayuVoaE9walbaMRIjjr6bzjrwycNuPid0tKlmwKUiHimurqazvb9PUuGd2tv28MHkff43h13UFYW/8aN0tKlLwUoEfFMeXk5dXV1rLptXk+Qam/bw6rb5lE3sy7prT2lpUtf5ly8hQn8oba21jU3Nxe6GSKShkgkwsKGRtavX09lVYjOjnZmzpzJ0iWLE/aewuEwo2pG9xq7gq7gVj9jCi27d2ncqoiZ2UvxFrFVkoSIeKqsrIym5cv49re+SWtrK9XV1QMGl2RjV91p6X5YglzyS7f4RCQnysvLGT9+fEo9n2RjV0pLT12xJZiUbIAaWlmJmfV7DK2sLHTTREpOsrErpaUPLNsEE78GNk8ClJlNN7NtZrbTzPqNaFqXO6Lff9nMzvbiuNnoXi6876Ndy4KLFMTSJYu5cGIt9TOmMG/6BdTPmMKFE2sHTEuXzBNM/J45mXWShJkdCWwHLgFagN8CVznnXovZ5zLgJuAy4FxgpXPu3IFeO5dJEmZGvDM3wM+JI5KaoZWV/T5slAHx/tuFKipU6dxHwuFwymNXkl2Cyfz6BTy7sZk5i1YQGjGyp9d64cRampYvy9cpJEyS8KIHNRHY6Zx73Tn3HnAfcHmffS4H7nVdfg1UmdnxHhw7J3S7L/ji9ZAj9O8xq9fsP+mMXUlqCSbxdFf96A5O3T/jp4K+XgSoUcDumOct0W3p7gOAmc0ys2Yza967d68HzUuf3rhEJCgyTTDJNLDlkxcByuJs63uPLJV9ujY6t8Y5V+ucqx0+fHjWjRMRKWaZJpgEIXPSiwDVAoyOeV4D9A29qeyTV4mWCw8VslEiIhlIlmCSKEMvEJmT8ZbGTudB19jz68A44Gjg98Bpffb5FPAIXTHgPGBjKq/90Y9+1OUD4FycR9evR4Io3jXVdZZid+DAAbdt2zZ34MABd/jwYXfL/HpXOaTK1Yw90VUOqXK3zK93hw8f7tm/e58hVSE3euyJbkhVqN8++QA0uzgxwJNSR9EsvRXAkcB659ztZjY7GgBXm5kB3wemA38GvuScGzA9L1+ljpTRV3yUxSelLp0MvUJnTibK4lMtPuK/mYHeuEQkmIJW2zCXaeaBt7+zM+4tRgUnEQmiIGTopUIBSkSkyAQhQy8VClBSNFRfUaRLIDL0UqDlNqRodFeP6Ms04VpK0NIli1nY0Ej9jCn91uUKCiVJSNFQNqZIf4XO0EuFFiwUESlB3bUNg0hjUCJxaDxLpPDUgxKJQ+NZIoWnHpQUjYT1FSsqCtouEcmMelBSNDSxWqS4qAclIiK+pAAlIiK+pAAlEofGs8TvEq3zVEwUoETiUAFh8atIJML8+gWMqhnN1GmXMqpmNPPrFxCJxFtMJtgUoEREAmRhQyPPbmym6cGnWfno8zQ9+DTPbmxmYUNjoZvmOZU6EulD64OJXwVtnadUaT0okRR1T9Lt+4gXtETyqVjWeUqVApRIhkphkFr8pVjWeUqVApRImnoPUk/n+FGjuOnmuUU5SC3+ku06T0H7UJVVgDKzoWb2uJntiP4birPPaDN72sy2mNlmM5ubzTFFCq33IPULrHjoWZ544UUmTb5AQUpybumSxVw4sZb5My7m5mnnM3/GxVw4sTbpOk9BzfzLtgfVCDzpnDsJeDL6vK8IUO+cOwU4D/iKmZ2a5XFFCiIcDrNu3TrmLFrRMw4QGjGSBSvW8MorLzPvlvkFbqGUCveBIxKJ4D4YONEtqJl/2Qaoy4F7ol/fA/x93x2cc2855/47+vUBYAswKsvjiuRMskm6XYPUobiD1FXDRnDPPfcE5vaJBFN3sPnuQ8+w6smNfPehZ5IGm0QfquYsWsH69et9/feabYAa6Zx7C7oCETAi2c5mNhY4C/hNkn1mmVmzmTXv3bs3y+aJpC/ZJN3q6mo69u+LO0h9sPMdhoSGFl0mlfhHJsEmyJl/AwYoM3vCzF6N87g8nQOZWTnwU2Cecy7hZBLn3BrnXK1zrnb48OHpHEIk58rLy7n+uutZNm9Wr0Hq7zXOZfKlMwh3vlN0mVTiD+FwmBdeeIHyIVVpBZsgZ/4NGKCcc59wzp0e5/FzYI+ZHQ8Q/bct3muY2VF0Bad/d879zMsTEMm37zYtp7zsCL4y7XzmXDKJeZ+5mOPGjKVt9xspZVKJpCM2weGG2V9m7563WPPNRt6PSXBIFmyyzfwrpGxv8T0AXBf9+jrg5313MDMD1gFbnHNNWR5PpODKysp48YXnuaHuBg52tDNkSBW/fvRBLjr3HJYuWRy4VF7xt75Zo3c9/mt27dzGun/5RyC1YNOd+Vc/Ywrzpl9A/YwpA2b++UFWpY7M7FjgfmAMsAv4vHNuv5lVA2udc5eZ2QXAL4FXgA+iP3qbc+7hgV5fpY7E78LhMK2trVRXVzNo0CAWNjSybt06KkND6WzfT11dHUuXLKasTGuDSvqSlTaac8kkho88jnDnO8ycOTOlv7PYv1c/9ZwSlTrK6n+Nc24fMDXO9lbgsujXz9OVBCVSdMrLyxk/fjwA8+sX9HzSDY0Y2fPJdmFDI03LlxW4pRJEyRIcRhx/PGtX38XkyZNTDjaxf69BoEoSIh4Iciqv+FeyBIcDHR1pBacgUoAS8UCQU3nFv4Kc4OAF3RgX8UDsJ92+YwV+T+UVf1u6ZDELGxqpnzGFyqoQnR3tPWNOxU7rQYl4pHsMqsgaLkoAAAy5SURBVPs2X/cn3Qsn1moMqkR5mZTg1wQHL2g9KJEcC2oqr3gvEonw1Ztu5rjjq7n4kmmeFGftTnAotuCUjHpQIh5L95NuMX8yLkWRSIRJky/glZdfJjR8BAfe6eCCyy6nbfcbXHTuOepNx6EelEieDPRJt3sib0dHRyCXQJDkbplfTzjyAXc+9ivufPxFVj70DH/a9QYjRo9VRmeaFKBE8qTfmjyjx/D4878K3BIIklg4HOaH9/yQBSvW9JpucNPilbzwyAOUVw5RRmcaFKBE8iS2ZM3in24AM+YtW6V5U0WktbWVqqHHxp1uMLhyCO+071dGZxoUoETyoO9E3va2P1GZ4I1M86aCq2u6QXvcibUdb7dx3XXXaZwxDQpQInnQdyJvaMRxHAjoEgiSWKKJtcvmzeIjH5nAiu+qXnY6FKBE8qBvyZpjBg9m6hVXcUfDzSVZIaCYxU43mDt9Mrd85iI+MXkSL77wvIoGp0lp5iJ50nci79utf+RbX/oftL/dxrDhI3tVCNAbWfBp+kDqEqWZK0CJ5EkkEmFhQyPr16/vVbLmG//rH2lra9MbmZQsBSgRn9An69Klax+fJuqK+EQplqwJIi9XRu43B06TslOiACUiEiMXwaT3su2alJ0q3eITEYnhdVX6ZMu218+YQsvuXSXfm9YtPpEikMltJy9vVRW7XKyMrMUsM5dVgDKzoWb2uJntiP4bSrLvkWb2OzN7KJtjipSiTG47adwjfbkIJsmWbdek7OSy7UE1Ak86504Cnow+T2QusCXL44mUpEzGMDTukb5cBJNSX7Y9G1mNQZnZNuAi59xbZnY88Ixz7uQ4+9UA9wC3A/Odc59O5fU1BiWS2RiGxj0yl4uVkRPNgdOk7C65GoMa6Zx7CyD674gE+60AFgIfDPSCZjbLzJrNrHnv3r1ZNk8kOBKNFWVy26nYxz1yOa6W6crIydpUVlZG0/JltOzexRMbHqFl9y6ali9TcBrAgAHKzJ4ws1fjPC5P5QBm9mmgzTn3Uir7O+fWOOdqnXO1w4cPT+VHRAJtoLGiTG47Feu4Rz7G1dINJum0SXPg0jNggHLOfcI5d3qcx8+BPdFbe0T/bYvzEpOBGWb2BnAfMMXMfuThOYgE2kBjRZmMYRTruEc+x9VSDSYa68udbMegvgPsc84tNrNGYKhzbmGS/S8CFmgMSqRLqmNFmYxhFNu4hx/H1fzYpiDK1RjUYuASM9sBXBJ9jplVm9nDWb62SNFLdawokzGMYhv38OO4mh/bVEyy+kt1zu0DpsbZ3gpcFmf7M8Az2RxTpJjEjhX1/QQeb6yo+7ZTOjL5GT9K93dVqm0qJqokIVJAfhkrCkK1Cb/8rvzepqLinPPt46Mf/agTKXaHDx92t8yvd0OqQm702BPdkKqQu2V+vTt8+PCAP3vgwAG3bds2d+DAgayOXTmkytWMPdFVDqnqOXa2r50L2fyuSqlNQQM0uzgxQMViRXwinbWCuhMg1q1bR2VoKJ3t+6mrq0s7ASLepNQ7b53LkYf/wvZt27J67Vzq/l1VVlbS2dnpi/WVtNZT5lQsVsTn0pkj40Vqc6LCqCNGjyUc+cDXadODBg1i9d1rOPnDp/imzqDmOHlPPSiRgPEqtXn79u1MnXYpKx99vmfbuwcPcuPFtax86Blfp03nohyRFI56UCJFwqvU5njVJtrb/kTFkKqcpU17kYyRiyUxxJ8UoEQCxqsyRvEy0MyOoH1vm+clkrwsUaS5R6VDAUokYLxMbe5bGPXrV36aj0yY4GnadDgc5n/OmsVTL/7Gk3GtZAH67b17GDGif83qIKTRSxzxUvv88lCauUh8Xqc2x6aUe/Xa3a9TMWSIqxo2wv1tRaWbMXO2u//VXe6nW1vd2ud+54ZUhTJKY79lfr07feIkt/a53/W81oRJH3PVJ4xzt8yv79eGeGn04h8ozVyk+OQytTmd1463b7xEhu81zuWED5/KdQu/AcC86RfwxIZH0q500dHRwajRY7AjjqAiNJRwRztTPncln/pCHQs/+8meZI64bWi4ibM+fBLfu+MOXyR8SOIkCQUoEclYovlY3/hf/8gJY8fFzTSc95mLWf3Ubzl0MJxxZmB3BuLin26gve1PhEYcxzGDBwN/DXrV1dW9sh3fj0T4UdMinvg//8HRfzOI9w//hRvqbvDV/K5SlShA6aqISMZi52PFpnvfMn9+wkSG8qoQb2zdzM/u+m7G41rd41CHDoapHvehnu2xyRx9kyl+1LSIN7e+xh2/eLZXWxc2NCo13aeUJCEiGUmW7v2zn/0X77Tvi5vIsH/PWyydc31Kq9QmkkqiSGwyxbsHD/LkT37MTYtXKjU9QBSgRCQjydK9h4SGcsXnPtcvgNx561yuufpq/tiyO+ulPwZamj02iL2xdTMVSk0PHAUoEcnIQPOxmpYv7xdALjr3HNbcfbcnyQmprHfVHcSWzrme/X96y/P5XZJbSpIQkYylUnLID0VUw+EwN8+dy6ZtO/nKv65UeSSfURafiHguSMvKB6mtpUYBSkRyxg+9pFQFqa2lQmnmIpIzQVpWPkhtLXVZJUmY2VAze9zMdkT/DSXYr8rMfmJmW81si5lNyua4IiJS/LLN4msEnnTOnQQ8GX0ez0rgUefch4EzgC1ZHldERIpctgHqcuCe6Nf3AH/fdwczqwQ+DqwDcM6955zryPK4IiJS5LINUCOdc28BRP/tX+ceTgT2Aj8ws9+Z2VozG5zoBc1slpk1m1nz3r17s2yeiIgE1YABysyeMLNX4zwuT/EYZcDZwF3OubOAgyS+FYhzbo1zrtY5Vzt8+PAUDyEiIsVmwCw+59wnEn3PzPaY2fHOubfM7HigLc5uLUCLc+430ec/IUmAEhERgexv8T0AXBf9+jrg5313cM79CdhtZidHN00FXsvyuCIiUuSyDVCLgUvMbAdwSfQ5ZlZtZg/H7HcT8O9m9jJwJrAoy+OKiEiR83UlCTPbC7yZw0MMA97O4esXUjGfGxT3+RXzuUFxn18xnxvk7vxOcM71SzrwdYDKNTNrjldeoxgU87lBcZ9fMZ8bFPf5FfO5Qf7PT8ttiIiILylAiYiIL5V6gFpT6AbkUDGfGxT3+RXzuUFxn18xnxvk+fxKegxKRET8q9R7UCIi4lMKUCIi4kslFaDM7PNmttnMPjCzhKmSZjbdzLaZ2U4zC0RZpjTW5nrDzF4xs01m5vvlige6Ftbljuj3XzazswvRzkykcG4Xmdk70Wu1ycy+UYh2ZsLM1ptZm5m9muD7gb1ukNL5BfnajTazp6Nr9202s7lx9snP9XPOlcwDOAU4GXgGqE2wz5HAH+iqwn408Hvg1EK3PYVzWwo0Rr9uBJYk2O8NYFih25viOQ14LYDLgEcAA84DflPodnt4bhcBDxW6rRme38fpKhL9aoLvB/K6pXF+Qb52xwNnR7+uALYX6v9dSfWgnHNbnHPbBthtIrDTOfe6c+494D661r3yuwHX5gqgVK7F5cC9rsuvgapo4WK/C+rfWUqcc88B+5PsEtTrBqR0foHlnHvLOfff0a8P0LXA7Kg+u+Xl+pVUgErRKGB3zPMW+l8cP0plbS4ABzxmZi+Z2ay8tS4zqVyLoF6vVNs9ycx+b2aPmNlp+WlaXgT1uqUj8NfOzMYCZwG/6fOtvFy/AZfbCBozewI4Ls63vu6c61dtPd5LxNnmi1z8ZOeWxstMds61mtkI4HEz2xr9NOhHqVwL316vAaTS7v+mq0ZZ2MwuA/4vcFLOW5YfQb1uqQr8tTOzcuCnwDznXGffb8f5Ec+vX9EFKJdk/aoUtQCjY57XAK1ZvqYnkp1bimtz4Zxrjf7bZmb/RdetJr8GqFSuhW+v1wAGbHfsm4Jz7mEzW2Vmw5xzxVCMNKjXLSVBv3ZmdhRdwenfnXM/i7NLXq6fbvH191vgJDMbZ2ZHA1fSte6V3w24NpeZDTaziu6vgU8CcbOQfCKVa/EA8MVoVtF5wDvdtzp9bsBzM7PjzMyiX0+k6//rvry3NDeCet1SEuRrF233OmCLc64pwW55uX5F14NKxsz+AfgeMBz4hZltcs5NM7NqYK1z7jLnXMTMvgpsoCvTar1zbnMBm52qxcD9ZlYH7AI+D11rcxE9N2Ak8F/R/zdlwH845x4tUHsHlOhamNns6PdXAw/TlVG0E/gz8KVCtTcdKZ7bFcCXzSwCvAtc6aIpVH5nZj+mK5NtmJm1AP8EHAXBvm7dUji/wF47YDJwLfCKmW2KbrsNGAP5vX4qdSQiIr6kW3wiIuJLClAiIuJLClAiIuJLClAiIuJLClAiIuJLClAiIuJLClAiIuJL/x/VRmPXWXcyTgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[y_km_gpu == 0, 0], X[y_km_gpu == 0, 1],\n",
    "            c='lightblue', marker='o', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 1')\n",
    "plt.scatter(X[y_km_gpu == 1, 0], X[y_km_gpu == 1, 1],\n",
    "            c='red', marker='s', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 2')\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"gpudbscan\"></a>\n",
    "## Clustering using GPU Accelerated DBSCAN"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we load the `DBSCAN` class from the `cuml` package and instantiate it in the same way we did with the `sklearn.cluster.DBSCAN` class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cuML Version: 0.12.0a+736.gd722488\n"
     ]
    }
   ],
   "source": [
    "import cuml; print('cuML Version:', cuml.__version__)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cuml import DBSCAN as DBSCAN_GPU\n",
    "\n",
    "db_gpu = DBSCAN_GPU(eps=0.2, min_samples=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `DBSCAN` class from `cuml` implements the same API as the `sklearn` version; we can use the `fit` and `fit_predict` methods to fit our `DBSCAN` model to the data and generate predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_db_gpu = db_gpu.fit_predict(X_gdf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, let's visualize our results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3de5hU1Znv8e+r4GmH7pZuAbm26BMwXoJIOiiigWAMSiaQyZhziI6aQEIIXkAwQPQ88XKiB4xy1Bh0GCFjzEyMx2SOJFHxFnWSaAQTYmy5xhsdEFpAoTUEWt7zR1cz1dVV1XXZVbWr9u/zPPXQVbW69trsrnprrfWutczdERERCZvDSl0BERGRZBSgREQklBSgREQklBSgREQklBSgREQklHqUugLp9OnTx4cOHVrqaoiISAG99NJL77h738THQx2ghg4dypo1a0pdDRERKSAzezPZ4+riExGRUFKAEhGRUFKAEhGRUAr1GJSISNgcOHCA5uZm9u3bV+qqlJ2qqioGDx5Mz549MyqvACUikoXm5mZqamoYOnQoZlbq6pQNd2fnzp00Nzdz3HHHZfQ76uITEcnCvn37OProoxWcsmRmHH300Vm1PBWgRESypOCUm2z/3xSgREQklBSgJJJaW1vZuHEjra2tpa6KSCCuv/56br311qx/791332Xp0qV5H/+uu+7iIx/5CGbGO++8k/frgQKUlJEggkpbWxtz513NoMFDOGfi+QwaPIS5866mra0twJqKdBbmL0S5BCh35+DBg50eGzt2LE8++STHHntsYHVTgJLQCzKozF+wkGdfXMOSn/+KOx77NUt+/iuefXEN8xcsLEDNJeoK9YXohz/8ISNGjODUU0/l4osv7vL8+PHjDy0T984779CxpmlTUxOjR49m5MiRjBgxgk2bNrFw4UL+/Oc/M3LkSL75zW8C8N3vfpdPfOITjBgxguuuuw6AN954gxNPPJFZs2YxatQotmzZ0umYp512GkGvnao0cwm9+KBS1+8Ydu/YztJr5jB/wUKW3JZ5l0ZrayvLly8/9DoAdf2OYdbNtzNv8gRuvOF6qquraW1tZevWrQwcOJDq6upCnZZEQFB/u/Gampq46aab+M1vfkOfPn3YtWtXxr97zz33MHv2bC666CL279/Phx9+yKJFi3jllVdYu3YtAI8//jibNm3ixRdfxN2ZPHkyzz33HA0NDWzYsIEf/OAHgXQJZkItKAm1jqAy6+bbuwSVFStWZNVlsnXrVmrr6g+9Toe6fsdQ27uOt956S91/Epgg/3bjPf3001xwwQX06dMHgPr6+ox/d8yYMdx8880sXryYN998kyOPPLJLmccff5zHH3+c0047jVGjRrF+/Xo2bdoEwLHHHssZZ5yRU71zoQAlodZdUNm6dWvGrzVw4ED27N7F7h3bOz2+e8d29ry7m6V336PuPwlMkH+78dy923TtHj16HBojip93dOGFF7Jy5UqOPPJIJk6cyNNPP5309b/1rW+xdu1a1q5dy+bNm5k+fToAvXr1yqnOuVKAklDrLqgMHDgw49eqrq5m+vTpLL1mzqHX6+hyufjii7n//vsD/7Yr0RXk3268c845hwcffJCdO3cCJO3iGzp0KC+99BIADz300KHHX3vtNY4//niuvPJKJk+ezMsvv0xNTQ179+49VGbixImd/ub/8pe/sGPHjpzqmi8FKAm1dEFl2rRpWY8R3bJ4EeNGNzJv8gTmnHcW8yZPYNzoRr4x8+sF+bYblDBngUlyQf/tdjj55JO59tprGTduHKeeeipz587tUubqq6/m7rvv5swzz+yU8v2Tn/yEU045hZEjR7J+/XouueQSjj76aMaOHcspp5zCN7/5TT7zmc9w4YUXMmbMGD72sY9xwQUXdApgqdx5550MHjyY5uZmRowYwVe/+tWczq8Td8/7BqwAdgCvpHjegDuBzcDLwKhMXvfjH/+4S2Xau3evb9iwwffu3dttuaamJr/8iiv9qN51PmTo8X5U7zq/au48P3DgQGDH37t3r9ce1dvvfe4P/tP1Ww/d7n3uD35U77pu61koBw4c8KvmzvPao3r74KHHe+1RvfM+d8nPq6++mnHZjusX5N9uuUv2/wes8WSxI9mD2d6ATwKj0gSoScCjsUB1BvC7TF5XAaryZPqBm6zcZZdf4U1NTQULFlfNneejzhp3KEjd+9wffNRZ4/yqufO6/d1MA24x6ySFkU2A6lCov49yVPQA1f76DE0ToP4Z+FLc/Q3AgO5eUwGq8mT6gVuKD+Zcvu0WsoUT1lZd1OUSoOS/ZBOgijUGNQiIn9XVHHusCzObYWZrzGxNS0tLUSonqQU59pFp2m3Q6bmZnkOPHj1YctutNG95iydXPUrzlrdYctut9OiRerpgISf+FioLTKRcFCtAJcuJ9GQF3X2Zuze6e2Pfvn0LXC1JpRAz4DP9wA3qgznXc6iurmb48OHdDmIXap5Lh0JlgYmUi2IFqGZgSNz9wYC+/oVYIVoGmX7gBvXBXOhljQrdwilUFphIuShWgFoJXGLtzgDec/dtRTq2ZKijK+ztt98uSMsg0w/cID6YC926geK0cFKlxd+yeFHery0SdoEEKDP7MfA8cIKZNZvZdDObaWYzY0UeAV6jPc38X4BZQRxXgpHYFTbshBM4rGdPauuP7lQuiJZBph+4+X4wF6J1kziWVYwWTi7jYhJNpd5u46KLLuKEE07glFNOYdq0aRw4cCDv1wwsi68QN2XxFUeyjLmPjvqET5x6ScGyx7KZB5VLem6QGXDpMvU0zyV6wprFd9111/l3v/vdrH/v9ddf95NPPjmr3zl48KB/+OGHnR775S9/6QcPHvSDBw/61KlTfenSpUl/N4xZfBJSqbrCrr59Gc88/H/Z9sbrQPAtg0wTETItl+z3gmrdpBvLStXC2bdv36HWllaBiK762lrMrMutvrY2r9cN43YbkyZNOnR+o0ePprm5Oa9zBNSCirqmpibvO2BQp1ZGx63vgEFeXVNbti2DIFo32bbEEltbR/5dL/+76hofdOxxWgWiQmTTggLap5sm3No/enPzyiuv+PDhw72lpcXd3Xfu3OnunVtQ48aN89WrV7u7e0tLix977LHu7n755Zf7j370I3d3/9vf/uYffPBBlxbUqlWr/Gtf+9qhVtJnP/tZf/bZZ/311193M/Pnn38+bf3279/vp512mj/33HNJn8+mBaWO7Ihbevc9hwb648drdu/Yzv6/fsCmjRvYs2dPWe6N1NG6ufGG63Pe3ymTsazhw4cfejzZ/j93LriSoSeezOQvfz3vvYByof2tKku+223cdNNNNDc384UvfIFhw4Z1KRO/3Qa0//1s2rSJhoaGjLbbmDVrFp/85Cc5++yzszir5NTFF2Gtra3cf//9jP/8F/newtmdusJunTODiy++mP79++fUxRYmuXYTQnaZeqm6S69cfCdP//QBqnpVF3V1dG1vX5ncw7vdxg033EBLSwtLlizJ9rSSUoCKsI7WwYVXfYs+Awcx++/Hc9nEscz53KfY/tYbfGPm10tdxZLLZiwrXWuruncdu3e8XdRVILS9fWUK63Yb9957L6tWreLHP/4xhx0WTGhRgIqwfv36seudFr7x6dNpWv0C7s5Jjadz078/zMG2NhoaGkpdxVDINOU9XWur9d3d1PXrX7RVIIoxD0xKI6zbbcycOZPt27czZswYRo4cyY033pj3uVr7+FQ4NTY2ekcmigRv7ryreeLXv2XOrUs7jZe8s7WZ//GPXyjqOEk5yGQsZ+68q3n2xTWHAkOyMahxoxsL/n+7ceNGzpl4Pnc89usuz8057yyeXPVop7Ezydy6des48cQTMypbX1vL7iQf7nU1NezasyfoqpWFZP9/ZvaSuzcmllWAiqjW1lYGDR5yaDC/w+4d27ni/LPY2txM7969S1jD8tTW1sb8BQtZsWIFNb3r2NmyncMOO5y6+qPZ+967TJs2jVsWLyr4RNt013fe5Ak0b3mrrMcVSymbACVdZROg1MUXUenGS/r0PaZkWzyXu/h5UU+tepQdb7/N21v/wlOPP1bUVSC0jp9UAqWZR1T8eEniN+xir5RdiWnQHZmDHUrRnXbL4kXMX7CQeZMnUNu7jj3v7j7UgpP8ZJJJJ11l22OnFlREheEbttKgC0vr+BVGVVUVO3fuzPrDNurcnZ07d1JVVZXx72gMKsLix0sSv2EX40MsWUJBsZIIRHJ14MABmpubO80vksxUVVUxePBgevbs2elxJUlEXLputFJ0sWkQX0Q6KEkiojLpRstnpYVcaTtzEemOAlSFC+tqAtrOXES6owBVwcK8mkAYkjSke9oqREpJAaqChb0bTduZh5cyLCUMFKAqWNi70ZQGHV5h7RqWaFGAKnPpumDKpRutFEkaklqYu4YlWgIJUGZ2npltMLPNZtblK5aZHWVmPzezP5pZk5l9JYjjRlmmXTDqRpNshb1rWKIj7wBlZocD3wfOB04CvmRmJyUUuwx41d1PBcYDt5nZEfkeO8oy7YJRN5pkK+xdwxIdQbSgRgOb3f01d98PPABMSSjjQI21L15VDewCNNqao1y6YNSNJpkql65hqXxBBKhBwJa4+82xx+LdBZwIbAX+BMx294PJXszMZpjZGjNb09LSEkD1Ko+6YKKtGKnf6hqWMAgiQCVb0jdx/aSJwFpgIDASuMvMapO9mLsvc/dGd2/s27dvANWrPOqCiaZipn6ra1jCIIgA1QwMibs/mPaWUryvAD/zdpuB14GPBnDsipDtN2J1wURTKVK/1TUspRREgFoNDDOz42KJD1OBlQll3gLOATCzY4ATgNcCOHZZy+cbsbpgoiXf1G+tCCHlKO8A5e5twOXAKmAd8KC7N5nZTDObGSv2v4AzzexPwFPAAnd/J99jl7t8vhGrCyZach131IoQUs603UaJaLsJyUaufy/ac0vKgbbbCBll4kk2chl31IoQUu7UH1QiHZl42954HfeD1PXrz5G9eikTT1K6ZfEi5i9YyLzJE7rsgJxMJl+Chg8fXoyqi+REAapEqqqqGH7CCcz7/Kep69uPve+9y1mTprBjyxvKxJOkOsYdb7zh+ox2QI6fjpDYLagvQVIO1MVXIvMXLOTgEVV8//Hf8v0nnueOXzzDmxvXcfiBvykTT9LKNPU7225BZfpJ2ChJogSUICHF0tbWxvwFC1mxYkWXbsGOjM+OMsuXL6e2rp49u3cxffr0TmVECklJEiESZIKEvvVKOplMR9DeTxJWClAlEMRSRZrfItlI1S2oTD8JMwWoEghiqSJ965UgaLqDhJnGoEokk7GBVDSGJUHR35KEgcagQiafpYr0rVeCooWHJcyUolNiHWMD2dD8FglSthOARYpFXXxlSmusSdBaW1szmgAsErRUXXwKUGUqnzEsEZEwUYCqUPrWKyLlTkkSFUo7nkpQUk361mRwKRUFKJGISzXpe9++fZoMLiWlwQqRiIuf9B2fcHP2uPEcPKKqy+PzFyxUIo4UhcagRCIs3UTdWeeewZKHn2bA0OM6Pa4JvBI0jUGJSBfpJn33qu2N+8Euj2syuBRLIAHKzM4zsw1mttnMki4GZ2bjzWytmTWZ2bNBHFdE8pNu4eL397yL2WFdHtdkcCmWvAOUmR0OfB84HzgJ+JKZnZRQpjewFJjs7icDX8z3uCKSv3RLHY04dSQrvnONlkCSkgmiBTUa2Ozur7n7fuABYEpCmQuBn7n7WwDuviOA45YNpelKmN2yeBHjRjcyb/IE5px3FvMmT2Dc6Eb+89lnkj6uJZCkWPJOkjCzC4Dz3P2rsfsXA6e7++VxZW4HegInAzXAHe7+wxSvNwOYAdDQ0PDxN998M6/6FVv8xNmqqirtVCplI9Wkb00Gl0JLlSQRxKekJXksMer1AD4OnAMcCTxvZi+4+8Yuv+i+DFgG7Vl8AdSvKJJtmz38hBP4sOd/U5qulIVUCxfnsqCxSBCC6OJrBobE3R8MJKb4NAOPufv77v4O8BxwagDHDo1kGwi2th2k35Ch2qlUIkld25KvIALUamCYmR1nZkcAU4GVCWUeBs42sx5m9nfA6cC6AI4dCqm2zb769mX85tGV/PX99w+VVZquVLpUK1NoBQrJVt4Byt3bgMuBVbQHnQfdvcnMZprZzFiZdcBjwMvAi8C97v5KvscOi/RzSY5i9463Dz2mNF2pdMl6E559cQ3zFySdgSKSklaSCEC62fiXTTyT2/7jSQYMPU57NknF0xbykgutJFFA6eaSfOxjI7h26t8rTVciIV1vgrq2JVvKdQ5Ium2z9+3bpzRdiYT4lSkSW1Dq2pZsqYsvYJozIlE3d97VPPvimkNJQ+ralu5oR10RKYqOOYErVqzo0pugCeqSjAKUiBRVfG8CoJ4FSUlJEiJSVNXV1Rx//PF8+7rrNSdKcqIAJSIFozlRkg918YlIQWhOlGRKXXwiUlSaEyX5UoASkYJIt1uv5kRJJhSgRCRn6VYsT7fCinbllUwoQElZq6+txcy63Opra0tdtYqW6YrlqXbr1XJfkgklSUhZM7Muu2NC+y6aYf7bLnfZrhahFVYkHU3UzZHeWOHWEaDqgd1Jnq+rqWHXnj1FrlVlU3aeBE1ZfFnSpmvlZTfgSW679+4tZbUqkrLzpFgUoFJINcHwqrnzkg4Ka3triQpl50mxKEAlkWoL91k3386yZf/Mp86deKhFtW/fPrW0JFKUnSfFoqWFk0jXhVHXrz/f+uf7ObJXNUuvmcPZ48Zz8IiqQ/3xHW/U+QsWamuBIqirqcHUjVd06fY/EwmKkiSSSDcIPOdzn+Kep1dzZK9ebHvjda6aMoG7n3hBg8VFVF9bm3RsqQ7YlfCYsvkKS0lEEoSCJkmY2XlmtsHMNptZylUgzewTZvahmV0QxHELJVUXxvcWzmbCP07lyF69AHA/SHVtbw0WF9nuvXuTJ0TQHpDib3U1NaWqZiRUV1czfPhwBScpiLwDlJkdDnwfOB84CfiSmZ2UotxiYFW+xyyG+AmGs88by6xzz6B/w1D+ae41h8qYHUbrnnc1WBwi7t7pphRzkfIVRAtqNLDZ3V9z9/3AA8CUJOWuAH4K7AjgmAXXo0cPltx2K81b3uKpVY/xta/NoKX5Tfbs2gm0B6EV37mGU08dqcFiEZECCCJJYhCwJe5+M3B6fAEzGwT8AzAB+ES6FzOzGcAMgIaGhgCql5/q6moGDhzIrG/MhLvv6TIofPNN3+Gaa/+nBotFRAIWRICyJI8ljkrfDixw9w/NkhWP+0X3ZcAyaE+SCKB+OWtra2P+goUsX76c2rp69uzexSWXXMI3Zn6dhoaGQy2kJbfdyo03XK/BYhGUOCHBCaKLrxkYEnd/MJCYHdAIPGBmbwAXAEvN7PMBHLugkk3W/e0f/si9y1d0eeNpsLjwOhaG7UHXZAglRJSeVl+RoOWdZm5mPYCNwDnAX4DVwIXu3pSi/L8Cv3D3h7p77VKuxaf1xsJHC8OGW7YLyIp0KFiaubu3AZfTnp23DnjQ3ZvMbKaZzcz39UtF642JZC7d6isrVqzQEmCSk0BWknD3R4BHEh67J0XZLwdxzEKLX28ssQWlFHKRzjL5Qjd8+PAS1U7KldbiS0HrjYlkTgvISiEoQKWh3UBFMqMvdFIIWosvA0qbDYeUa/BpU8JQ6JiWsWLFii5zAnv00LrUkpp21BWRotAXOslWqgClrzUiEqiOOYEi+YrcGJR2vhUpPr3vJBeRCVCa5S5SfJ3fd+cxYNAgrrhytt53kpHIBKhkyxY9++Ia5i9IuX2ViOSp8/vuN9z+i2d58jfPM2bsWYeClFpXkkokkiS0bJFI8aV731028UymfWUaVVVVnRZjnj59urL+IqigO+qGnZYtqmwdi8gm3upra0tdtUhrf9/VJX3f9e7Tj3+971955ner1ashKUUiQGmWe2VLuQV8kjlTUjwDBw7k3V07k77vWt97l7YDbVz2v+/Q2n2SUiQClGa5ixRfdXU1X770y9w6Z0an9933Fs5m1Ccn0Kv2KPVqSFqR6ei9ZfEi5i9YmHTnW00sFCmM/7PkNsaMPYvLJp5J7z79eH/Pe4w9fzLb3vgzf/vgfS3GLGlFIkkiXnwwqqqq6rJjrgZpy4/2iQq3trY25lw1l/vuu4+j6upp3fMe06ZN46Af5D9X/177R4mWOkpGG6xVBgWo8pDYU6G1+6SDAlQCpZ5XDi0iW97UxS6RTjNPRqnnlWPXnj24e5ebglN56Fi7T8FJEkU2QCn1XEQk3CIboJR6LiISboEEKDM7z8w2mNlmM+syDdzMLjKzl2O335rZqUEcN1/d7ZirNcJEREon7yQJMzsc2AicCzQDq4EvufurcWXOBNa5+24zOx+43t1P7+61i7VhYarsIqWfi4gUXiGTJEYDm939NXffDzwATIkv4O6/dffdsbsvAIMDOG5GMmkFJQ7SauVzEZHSCyJADQK2xN1vjj2WynTg0VRPmtkMM1tjZmtaWlpyrlSu+z+1trayfPnyQ3OjQGuEiYiUQhABypI8lrTf0Mw+RXuAWpDqxdx9mbs3untj3759c65Urq0gpZ+LiIRDEAGqGRgSd38w0OVT3MxGAPcCU9x9ZwDHTSmfVpDSz0VEwiGIALUaGGZmx5nZEcBUYGV8ATNrAH4GXOzuGwM4Zlr5tIKUfi4iEg55p6S5e5uZXQ6sAg4HVrh7k5nNjD1/D/Bt4GhgqZkBtCXL2AhKfCsol5WS0618LuGUbrkjSL43lJZCEgm3il2LL4iFYLVGWPlIt2AsJB8U1WKyIuEQucVitVJytChAiZSvyAWoDmoFRYMClEj5ShWgKr4p0TEJV0REyktkF4sVEZFwU4CSilBXU4NBl1tdTU3a50QkvCq+i0+iQeni5UljxJKOWlAiUnS5rpUp0aIAJSJFN+equax69jlu/skvtWOApKQAJSJF09bWxhVXzmbZvyxj144dLPjvk7jvlhuprT9aOwZIFwpQIlI08xcs5NcvreXuJ17g+088zx2/eIY317/Kj5bcrB0DpAsFKBEpio5dBi5fdEenXQauWHQHT//0Aba98bp2DJBOFKBEpCjS7TLQq/Yoll6b2Y4BmeySLZVBAUpEiiLdXmu7d2znU2POSLtjgDL/okcBSkSKItVea3ctnM2MGTO463t3pl3IOdddsqV8VfxisSISHrnuMtDa2sqgwUNY8vNfddnjbd7kCTRveUsTfctYZBeLFZHw6NGjB0tuu5Ubb7g+qxUkMtklW4tCV55Id/FpsFWkNDp2Gci01ZNu/EqZf5UrkgFKg60i5SXV+NXSazLL/JPyFEiAMrPzzGyDmW02sy4jltbuztjzL5vZqCCOmysNtoqUn1sWL2Lc6EbmTZ7AnPPOYt7kCYwb3Zg28086q6+txcy63Oprawv6u7nKO0nCzA4HNgLnAs3AauBL7v5qXJlJwBXAJOB04A53P7271y5EkkS6wdarPjeejevX079//0CPKeFTX1vL7r17uzxeV1OjldFDTiug5y7dztPdxYJ8fjeDeiVNkgiiBTUa2Ozur7n7fuABYEpCmSnAD73dC0BvMxsQwLGzlm6wtccRVXxk2HB191WgxG9/u/fuxaHLLVnQknDJdvxKylcQAWoQsCXufnPssWzLFEW6wdYD+//G4oceU3dfBUoMSCISfkEEKEvyWOJnQCZl2guazTCzNWa2pqWlJe/KJUo12Pq9hbOZ8I9TGTD0OK2qLCISAkEEqGZgSNz9wUDicsSZlAHA3Ze5e6O7N/bt2zeA6nXVMdh61efG89WzT2PO5z7FsR89iX+aew2AVlUWEQmBIALUamCYmR1nZkcAU4GVCWVWApfEsvnOAN5z920BHDsnHZMFN65fz/6/fsCin/ySS+d/m8NjM9k1t0JEKlFdTQ0GXW49odsMvVS/W1dTU7D65h2g3L0NuBxYBawDHnT3JjObaWYzY8UeAV4DNgP/AszK97hB6N+/PzNmzGDFd67R3IqIqaPrG63QbzaRUtu1Zw/u3uV2gK4JQ4lJQ6l+t5BZr5Ffiy/XtcGkvCitXCS1QqaQZ3j8pGnmkQ9QHTS3QkSiKqwBSk2EmI65FSIiEg6RXItPKl8plmURkWApQElF0koRIpkrRYZeJtTFJyIScWFNFFILSkREQkkBSkREQkkBSiQDSroQKT4FKKlIQQ/6KulCpPiUJCEVKayDviKSObWgREQklBSgREQklBSgREQklBSgRDIQ1pn2IpVMAUokA6XYC0ckH5UwNUIBSqQblfBGl+iphKkRSjMX6UbHGz2RldEbXaQcqQUlIiKhpAAlkid1AYoURl4ByszqzewJM9sU+7cuSZkhZvYrM1tnZk1mNjufY4qETbq+fgUpKTdh+sKVbwtqIfCUuw8DnordT9QGzHP3E4EzgMvM7KQ8jytSFsppQFoqS65TI8KUXJFvgJoC3Bf7+T7g84kF3H2bu/8+9vNeYB0wKM/jihSN5kCJlEa+WXzHuPs2aA9EZtYvXWEzGwqcBvwuz+OKFI3mOkk5qoTs024DlJk9CfRP8tS12RzIzKqBnwJz3D3lO97MZgAzABoaGrI5hIiIVJBuA5S7fzrVc2a23cwGxFpPA4AdKcr1pD04/Zu7/6yb4y0DlgE0NjYm+wIgEip1NTVJv5XWAbuLXx2JqPra2i7jREb73+GuktQof/mOQa0ELo39fCnwcGIBMzNgObDO3ZfkeTyR0Nm1Z0/S8ajdaJxKiidlckOWrxOmMdd8A9Qi4Fwz2wScG7uPmQ00s0diZcYCFwMTzGxt7DYpz+OKhIrW6pNKEaa/5bySJNx9J3BOkse3ApNiP/+a9gAsEjnJul2g/duogpcUS/wHcDm16rUWn0gBVUImlZQ/9/IcztdSRyIiEkoKUCIiFSBMyQ1BURefiEgFqMQxTbWgREQklBSgRAqoErtdJBhHJFkx3Mw4wpT03EEBSqSAcplTEqbtDqRwDtB1Uq3HHpd2GoMSKaFU86SUmi6iACVSUsnmSamDR6SduvhERCSUFKBERCSU1MUnIlICPUnenduz2BUJMQUokZCpI/kHl7u4z3oAAAgdSURBVFLTK8v+Ml0fr5gUoERKKOVmh1rtXEQBSqSUFIREUlOShIiIhJIClIhICGgFka4UoEREslSIYNIxaTvxlmylkajQGJSISJa0U3Jx5NWCMrN6M3vCzDbF/q1LU/ZwM/uDmf0in2OKiEg05NvFtxB4yt2HAU/F7qcyG1iX5/FEJCbXbiaNdUi5yDdATQHui/18H/D5ZIXMbDDwWeDePI8nIjG5jllorEPKRb4B6hh33wYQ+7dfinK3A/OBg3keT0SkImlzy666TZIwsyeB/kmeujaTA5jZ3wM73P0lMxufQfkZwAyAhoaGTA4hIlJU6VYAyZUmbXfVbYBy90+nes7MtpvZAHffZmYDgB1Jio0FJpvZJKAKqDWzH7n7P6U43jJgGUBjY6MWq5JISrWRoZZASq9Y/2/ZvpauZ27y7eJbCVwa+/lS4OHEAu7+LXcf7O5DganA06mCk4i00zhRbsL6/xbWeoVdvgFqEXCumW0Czo3dx8wGmtkj+VZORFLLdcxCYx1SLvKaqOvuO4Fzkjy+FZiU5PFngGfyOaaItMu1a0hdSlIutNSRiIiEkgKUiIiEkgKUSAiVepyoXFebKPX/W7nVK+y0WKxICOU6ThRUOnOyxVDrY4+bdd6QPkyp0mGpR6Kw1ivs1IISqSCFTGfeneR1w5oqXa4tQOlMLSgRqTjaDqMyqAUlIiKhpAAlIiWl7jhJRQFKRLpIlnVWKMVeBig+8Ck4hpsClEgFCSqdedeePbh7p1vQqdIdwYGE16vP6dUyFx/4tEZeuClJQqSCFDKdOejXTpnIEMBrp9wOg/ZsRCkPakGJSEGUsvusI5gmtox2pasvnbsz1d1XempBiUhBlFuqd8c8r0RhrW8UqAUlIqGjZYAEFKBEJITcPe8xr0ySOuLLSPioi09ESiJlIkNALadMAlx8mcQ1BqX0FKBEpCS0gKp0R118IlIQ5bbFRLnVNwrUghKRgii3FlK51TcK8mpBmVm9mT1hZpti/9alKNfbzB4ys/Vmts7MxuRzXBERqXz5dvEtBJ5y92HAU7H7ydwBPObuHwVOBdbleVwREalw+QaoKcB9sZ/vAz6fWMDMaoFPAssB3H2/u7+b53FFRKTC5RugjnH3bQCxf/slKXM80AL8wMz+YGb3mlmvPI8rIiIVrtsAZWZPmtkrSW5TMjxGD2AUcLe7nwa8T+quQMxshpmtMbM1LS0tGR5CREQqTbdZfO7+6VTPmdl2Mxvg7tvMbACwI0mxZqDZ3X8Xu/8QaQKUuy8DlgE0NjYmWxpLREQiIN8uvpXApbGfLwUeTizg7m8DW8zshNhD5wCv5nlcERGpcPkGqEXAuWa2CTg3dh8zG2hmj8SVuwL4NzN7GRgJ3JzncUVEpMLlNVHX3XfS3iJKfHwrMCnu/lqgMZ9jiYhItJh7eId5zKwFeLNAL98HeKdArx0mUTjPKJwjROM8o3COEI3zzOYcj3X3vokPhjpAFZKZrXH3im/VReE8o3COEI3zjMI5QjTOM4hz1GKxIiISSgpQIiISSlEOUMtKXYEiicJ5RuEcIRrnGYVzhGicZ97nGNkxKBERCbcot6BERCTEFKBERCSUIhOgzOyLZtZkZgfNLGXqo5mdZ2YbzGyzmaVcMzCssthE8g0z+5OZrTWzNcWuZy66uzbW7s7Y8y+b2ahS1DMfGZzjeDN7L3bd1prZt0tRz3yY2Qoz22Fmr6R4vuyvI2R0npVwLYeY2a9iG9E2mdnsJGVyv57uHokbcCJwAvAM0JiizOHAn2nfIuQI4I/ASaWue5bneQuwMPbzQmBxinJvAH1KXd8szqvba0P76iWPAgacAfyu1PUuwDmOB35R6rrmeZ6fpH2Hg1dSPF/W1zGL86yEazkAGBX7uQbYGOT7MjItKHdf5+4buik2Gtjs7q+5+37gAdo3ZSwn3W4iWaYyuTZTgB96uxeA3rFV9stFJfz9dcvdnwN2pSlS7tcRyOg8y567b3P338d+3kv7bumDEorlfD0jE6AyNAjYEne/ma7/2WGXySaSAA48bmYvmdmMotUud5lcm3K/fpnWf4yZ/dHMHjWzk4tTtaIq9+uYjYq5lmY2FDgN+F3CUzlfz7wWiw0bM3sS6J/kqWvdvctWIMleIsljocvDT3eeWbzMWHffamb9gCfMbH3sG19YZXJtyuL6pZFJ/X9P+7plrWY2Cfh/wLCC16y4yv06ZqpirqWZVQM/Bea4+57Ep5P8SkbXs6IClKfZXDFDzcCQuPuDga15vmbg0p1nhptI4u0rzuPuO8zsP2jvXgpzgMrk2pTF9Uuj2/rHv/nd/REzW2pmfdy9khYeLffrmJFKuZZm1pP24PRv7v6zJEVyvp7q4utsNTDMzI4zsyOAqbRvylhOut1E0sx6mVlNx8/AZ4CkmUYhksm1WQlcEssaOgN4r6O7s0x0e45m1t/MLPbzaNrfwzuLXtPCKvfrmJFKuJax+i8H1rn7khTFcr6eFdWCSsfM/gH4HtAX+KWZrXX3iWY2ELjX3Se5e5uZXQ6soj2jaoW7N5Ww2rlYBDxoZtOBt4AvQvsmksTOEzgG+I/Ye6MH8O/u/liJ6puRVNfGzGbGnr8HeIT2jKHNwAfAV0pV31xkeI4XAN8wszbgr8BUj6VKlQsz+zHtGWx9zKwZuA7oCZVxHTtkcJ5lfy2BscDFwJ/MbG3ssWuABsj/emqpIxERCSV18YmISCgpQImISCgpQImISCgpQImISCgpQImISCgpQImISCgpQImISCj9f424ovFy/xZgAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[y_db_gpu == 0, 0], X[y_db_gpu == 0, 1],\n",
    "            c='lightblue', marker='o', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 1')\n",
    "plt.scatter(X[y_db_gpu == 1, 0], X[y_db_gpu == 1, 1],\n",
    "            c='red', marker='s', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 2')\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarking\"></a>\n",
    "## Benchmarking: Comparing GPU and CPU\n",
    "\n",
    "RAPIDS uses GPUs to parallelize operations and accelerate computations. We saw porting an example from the traditional scikit-learn interface to cuML was trivial. So how much speedup do we get from using RAPIDS? \n",
    "\n",
    "The answer to this question varies depending on the size and shape of the data. As a good rule of thumb, larger datasets will benefit from RAPIDS. There is overhead associated with using a GPU; data has to be transferred from the CPU to the GPU, computations have to take place on the GPU, and the results need to be transferred back from the GPU to the CPU. However, the transactional overhead of moving data back and forth from the CPU to the GPU can quickly become negligible due to the performance speedup from computing on a GPU instead of a CPU.\n",
    "\n",
    "Feel free to change the number of rows and columns to see how this speedup might change depending on the size and shape of the data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "NumPy Version: 1.17.5\n",
      "(100000, 128)\n"
     ]
    }
   ],
   "source": [
    "import numpy as np; print('NumPy Version:', np.__version__)\n",
    "\n",
    "n_rows, n_cols = 100000, 128\n",
    "X = np.random.rand(n_rows, n_cols)\n",
    "print(X.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_df = pd.DataFrame({'fea%d'%i: X[:, i] for i in range(X.shape[1])})\n",
    "X_gdf = cudf.DataFrame.from_pandas(X_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingkmeans\"></a>\n",
    "## K-Means"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingkmeansgpu\"></a>\n",
    "### GPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "km_gpu = KMeans_GPU(n_clusters=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 361 ms, sys: 33 ms, total: 394 ms\n",
      "Wall time: 394 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "y_km_gpu = km_gpu.fit_predict(X_gdf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingkmeanscpu\"></a>\n",
    "### CPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "km = KMeans(n_clusters=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 28.4 s, sys: 7.32 s, total: 35.7 s\n",
      "Wall time: 26.8 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "y_km = km.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingdbscan\"></a>\n",
    "## DBSCAN"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingdbscangpu\"></a>\n",
    "### GPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "db_gpu = DBSCAN_GPU(eps=3, min_samples=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 825 ms, sys: 708 ms, total: 1.53 s\n",
      "Wall time: 1.52 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "y_db_gpu = db_gpu.fit_predict(X_gdf)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"benchmarkingdbscancpu\"></a>\n",
    "### CPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "db = DBSCAN(eps=3, min_samples=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "This takes a very long time!\n",
      "CPU times: user 115 µs, sys: 37 µs, total: 152 µs\n",
      "Wall time: 108 µs\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "print('This takes a very long time!')\n",
    "# y_db = db.fit_predict(X)  # uncomment this line to execute - run at your own risk of boredom!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"conclusion\"></a>\n",
    "## Conclusion\n",
    "\n",
    "In this notebook, we described K-Means and DBSCAN - several popular clustering algorithms - and showed how to use the GPU accelerated implementation of these algorithms in RAPIDS.\n",
    "\n",
    "To learn more about RAPIDS, be sure to check out: \n",
    "\n",
    "* [Open Source Website](http://rapids.ai)\n",
    "* [GitHub](https://github.com/rapidsai/)\n",
    "* [Press Release](https://nvidianews.nvidia.com/news/nvidia-introduces-rapids-open-source-gpu-acceleration-platform-for-large-scale-data-analytics-and-machine-learning)\n",
    "* [NVIDIA Blog](https://blogs.nvidia.com/blog/2018/10/10/rapids-data-science-open-source-community/)\n",
    "* [Developer Blog](https://devblogs.nvidia.com/gpu-accelerated-analytics-rapids/)\n",
    "* [NVIDIA Data Science Webpage](https://www.nvidia.com/en-us/deep-learning-ai/solutions/data-science/)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}