HW3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import Packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "import numpy as np\n",
    "from itertools import cycle\n",
    "from collections import defaultdict\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mean-Shift Implementation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "class MeanShift:\n",
    "    def __init__(self, band_width=2.0, min_fre=3, epsilon=None, bin_seeding=False, bin_size=None):\n",
    "        self.epsilon = epsilon if epsilon else 1e-3 * band_width\n",
    "        self.bin_size = bin_size if bin_size else self.band_width\n",
    "        self.band_width = band_width\n",
    "        self.min_fre = min_fre  # 可以作为起始质心的球体内最少的样本数目\n",
    "        self.bin_seeding = bin_seeding\n",
    "        self.radius2 = self.band_width ** 2  # 高维球体半径的平方\n",
    "\n",
    "        self.N = None\n",
    "        self.labels = None\n",
    "        self.centers = []\n",
    "\n",
    "    def init_param(self, data):\n",
    "        # 初始化参数\n",
    "        self.N = data.shape[0]\n",
    "        self.labels = -1 * np.ones(self.N)\n",
    "        return\n",
    "\n",
    "    def get_seeds(self, data):\n",
    "        # 获取可以作为起始质心的点（seed）\n",
    "        if not self.bin_seeding:\n",
    "            return data\n",
    "        seed_list = []\n",
    "        seeds_fre = defaultdict(int)\n",
    "        for sample in data:\n",
    "            seed = tuple(np.round(sample / self.bin_size))  # 将数据粗粒化，以防止非常近的样本点都作为起始质心\n",
    "            seeds_fre[seed] += 1\n",
    "        for seed, fre in seeds_fre.items():\n",
    "            if fre >= self.min_fre:\n",
    "                seed_list.append(np.array(seed))\n",
    "        if not seed_list:\n",
    "            raise ValueError('the bin size and min_fre are not proper')\n",
    "        if len(seed_list) == data.shape[0]:\n",
    "            return data\n",
    "        return np.array(seed_list) * self.bin_size\n",
    "\n",
    "    def euclidean_dis2(self, center, sample):\n",
    "        # 计算均值点到每个样本点的欧式距离（平方）\n",
    "        delta = center - sample\n",
    "        return delta @ delta\n",
    "\n",
    "    def gaussian_kel(self, dis2):\n",
    "        # 计算高斯核\n",
    "        return 1.0 / self.band_width * (2 * math.pi) ** (-1.0 / 2) * math.exp(-dis2 / (2 * self.band_width ** 2))\n",
    "\n",
    "    def shift_center(self, current_center, data):\n",
    "        # 计算下一个漂移的坐标\n",
    "        denominator = 0  # 分母\n",
    "        numerator = np.zeros_like(current_center)  # 分子, 一维数组形式\n",
    "        for sample in data:\n",
    "            dis2 = self.euclidean_dis2(current_center, sample)\n",
    "            if dis2 <= self.radius2:\n",
    "                d = self.gaussian_kel(dis2)\n",
    "                denominator += d\n",
    "                numerator += d * sample\n",
    "        if denominator > 0:\n",
    "            return numerator / denominator\n",
    "        else:\n",
    "            return None\n",
    "\n",
    "    def classify(self, data):\n",
    "        # 根据最近邻将数据分类到最近的簇中\n",
    "        center_arr = np.array(self.centers)\n",
    "        for i in range(self.N):\n",
    "            delta = center_arr - data[i]\n",
    "            dis2 = np.sum(delta * delta, axis=1)\n",
    "            self.labels[i] = np.argmin(dis2)\n",
    "        return\n",
    "\n",
    "    def fit(self, data):\n",
    "        # 训练主函数\n",
    "        self.init_param(data)\n",
    "        seed_list = self.get_seeds(data)\n",
    "        for seed in seed_list:\n",
    "            bad_seed = False\n",
    "            current_center = seed\n",
    "            # 进行一次独立的均值漂移\n",
    "            while True:\n",
    "                next_center = self.shift_center(current_center, data)\n",
    "                if next_center is None:\n",
    "                    bad_seed = True\n",
    "                    break\n",
    "                delta_dis = np.linalg.norm(next_center - current_center, 2)\n",
    "                if delta_dis < self.epsilon:\n",
    "                    break\n",
    "                current_center = next_center\n",
    "            if not bad_seed:\n",
    "                # 若该次漂移结束后，最终的质心与已存在的质心距离小于带宽，则合并\n",
    "                for i in range(len(self.centers)):\n",
    "                    if np.linalg.norm(current_center - self.centers[i], 2) < self.band_width:\n",
    "                        break\n",
    "                else:\n",
    "                    self.centers.append(current_center)\n",
    "        self.classify(data)\n",
    "        return"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generate Samples"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_samples(mean, var, num):\n",
    "    samples = []\n",
    "    samples.append(mean[0] + np.random.randn(num) * var[0])\n",
    "    samples.append(mean[1] + np.random.randn(num) * var[1])\n",
    "    return samples\n",
    "\n",
    "samples1 = generate_samples([-2, -2], [1, 1], 200)\n",
    "samples2 = generate_samples([2, 2], [1, 1], 200)\n",
    "samples = np.vstack([np.hstack([samples1[0], samples2[0]]) , np.hstack([samples1[1], samples2[1]])])\n",
    "samples = samples.T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clustering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "MS = MeanShift(band_width=3, min_fre=3, bin_size=4, bin_seeding=True)\n",
    "MS.fit(samples)\n",
    "labels = MS.labels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualization Raw Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAD4CAYAAADxeG0DAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2df2xd53nfv++9unTJylLAKw8pYvMyf+iPNVO9zEKwwUNXVEsX026cGIPR7sog4gJEpBWggAVGEwJTHIArmgCzCWyyILQ2GPECQdA6zTLLaBunAbZiC0KlcdTEbW10JOOsRSWqU6WIDSny3R9Hhzr33Pd9z/ue855f934/wIXEw/PjPYfk933O8z4/hJQShBBC6kuj7AEQQgjJBoWcEEJqDoWcEEJqDoWcEEJqDoWcEEJqzoEyLnrkyBE5PT1dxqUJIaS2XL58+ZqU8oH49lKEfHp6Gqurq2VcmhBCaosQYl21na4VQgipORRyQgipORRyQgipORRyQgipORRyQgipORRyQkqgd6WH6Ren0Xi+gekXp9G70it7SKTGUMgJKZjelR7mvjaH9RvrkJBYv7GOua/NFS7mnEyGBwo5IQWz8MYCbu/c7tt2e+c2Ft5YKGwMVZlMiB8o5IQUzMaNDafteVCFyYT4g0JOSMFMHZ5y2p4HVZhMiD8o5IQUzOKJRUy0Jvq2TbQmsHhisbAxVGEyIf6gkBNSMN1jXVz45QvoHO5AQKBzuIMLv3wB3WPdwsZQhcmE+EOU0bPz+PHjkkWzCCmX3pUeFt5YwMaNDUwdnsLiicVCJpOyrjsMCCEuSymPD2ynkBNC8kAl2AAw97W5voXWidZE4W8kdYVCTggpjDC8MS7Y4wfGsbm1ObB/53AHa2fWChxhPdEJOX3khFSYuibt6MIbVSIOBNEydb3XKkAhJyQHTr92Ggc+dwDieYEDnzuA06+ddj5HnZN2XMMYJ8cna3uvVYBCTohnTr92Gi+tvoRduQsA2JW7eGn1JWcxr3PSji6MsT3eVkbLAKjtvVYBCjkhnrlw+YLTdh0uSTtVc0vowhuXHltShl5e37quPA8TlOzw1rNTCNEEsArgR1LKJ3ydl5C6EVrittt1TB2ewvqNwRaNcWs3vrAYuiUAlBYJEl5XF2YYH9fCGwtW90rU+LTI5wG85fF8hNSSpmg6bddhm7RTVRdM91gXa2fWsHd2D2tn1oyTChOUsuFFyIUQDwJ4HMBv+zgfIXVm7pE5p+06bDNAde6H9RvrhbhbfLh1qpDtWme8xJELIX4XwG8CuB/Ap1SuFSHEHIA5AJiamnpkfX3wNYqQYeH0a6dx4fIF7MpdNEUTc4/M4dzj53K51vSL00q3hICAxL2/7zwSb3Tx4hThfMgtIUgI8QSAGSnlaSHEL0Aj5FGYEESIP1RiGhfxEN+JN7pJhAk++ZBnQtCjAD4qhFgD8CUAvyiEWPFwXkKGhjyjSlRuCZWIA/6jQFgOtxpkFnIp5aellA9KKacB/AqAb0gpT2YeGSFDQhGJPfGFxc7hjnI/31EgLIdbDRhHTkjO+I4qsbHui4oC0V1n5uhMpeLahx1vceQAIKX8JoBv+jwnIXXHp/vBNmY8KY7bF6rrzBydwfKby5WKax92WP2QkJzxuSCY9lxZaoC7HssF0Pxg9UNCSsKnmyONdZ/FR5/m2KosgFatbEGeUMgJyRmfyS5pFhdtfPQ60Uvj36/CAmidK0emgUJOSAG4pKubSGPdJ1nIJtFLY11XId2+qmUL8oJCTkiNcLXue1d6aAj1n3loIZtEL411nfYNxKcrpCrunaLwGrVCCLEjy+Jj91jXat/Q0lZVXYxayCbRu/jURWUKfpJ1bTvG+Fh9RbrYVo4cFmiRE1IwRflvVZY2EFRhjFrIJqs7tK7b4+397eMHxvfvw5cF7dsVUgX3TpFQyAkpmKL8tzpLe0/u9Vm5M0dnlPtFt2/d2dr//+bWJj7x+5/As1991ttk5NsVMmrVFOlaIaRgVK/8gH//ra174dLbl5THh9tVE8/O3s7A/uFkVBVXiKt7p87QIiekQHpXehAQyu/59t/auheSrGGXCWb9xnoqV8uouUJ8QyEnpEAW3lhQViYUEN5Fy9a9kBSZ4jrBpHG1jJorxDdM0SekQBrPN7QlZuXZ4v8WgeTmEKrv28K0fL8wRZ+QCqCzbnVlZwF9dIivqJEkazj6fVfSFgYbldR6X9AiJ6RAXFuj6faffXi2r8Kg6TxZYtZtxzN+YBybW5sD+7ta5GwdZya3Vm9poJCTUcZFWHWVBJuiqUz0aY+3cXDsoLakLJBdGFXjB5BagKPna4iG8r7oogmgkBNSAD6tX8DsU7chbe/O8D7Wb6zvTxqdwx3j/aS5d1v/u4DA3tk94z6jgE7IGUdOiCd8p5kD+vhqnUUeJ03vztOvncb51fP7x4bXSbqfNHHbuuzTOMOaWu8LLnYS4on51+e9Z2zq4qvnHpkb2O6CThh7V3p9Ih7HdwaqzWIo48mToZAT4oHelZ5ysQ+wj9xQRWvoIkrOPX5uYHu0HkqUeAKSSRh1ce5p7scG3YTSFE3GkztAHzkhlph8wLpFScBuoc5HtIYpwuXS25esfNc2PnmfC4/DFKXie31EBX3khGQgyf+d1GhBdb7oH/2t7Vtat4xLeVsgW8NlnU8+yvqNdRz5/BEsPbbkFJGiGk9RTaLzJo/1ERdokRNiQVJDYd332+NtXHvuWt82l0zJvKI1dALrMrax5hhefvJlYxRL/FxjzTHcP3Y/rm9dr61oqyiq4TQzOwlxJPRZi+dFYsVC3aLk0mNLA75v1aKojjyiNUz10EOffFM0E8+zvbttXPhURaRs725jc2uztD6aeWWNlt2RiEJOiIKo2JkIhVa3KAlgQDR1i6JxxppjA24ZH0KUVA+9e6yLPWn3FhDt+xkfl42IFdlHM8+GHmU3nKZrhRAFpsXLEJtFOZvz6Gg1WnjlY6/sn9/XwqBpQVNA7PvsbSacMEnIJW1fdc0ikn3ydH8UtWhL1wohDthYkzZ/pFlerXf2djD7ldl9izFLZ6GoxWwitFRvbt9Eq9Ey7hu+MejGFe6TRGi1xq3606+dHrDys7yR5On+KLsMLy1yQhQkWdK2VlwWizwktOxOvnpS+f0kizZtGdpo3ZbJ8Un8w51/wI93frz/vTBqRWfhCwhMjk8arfLw3oDBWi1xWo0WhBDY3t0eON5GMItakMwTWuSEOKBavAxRJdToLEXVecaaY1aLiSG3d25j/vX51J2FbNPg41zfuo61M2vYO7uHa89dw63P3II8KyHPSlx77ppV8+brW9e1549arTZj3Nnb6RNxwM3HPsxdiCjkhCiI1+AOhTf+yty70sORzx/ByVdPGqNAoq/cv/bBX4Prm3AY6RFH1VkoHJN4XhgjbpKwXagzCaSp/vrambX955jFvRE/Vjep+nB/VLVWOl0rhCQQjbmeHJ8EEFirk+OTuLl9c8BKDNG9svtwt0SJdhbqXenh2a8+qx2TLfEOQUkJOy5x6apsU9vFVRXR56xzI0VdQWmpQhYqy9gSkoIsbc5C33Vc5GyiYdI2bvA1SYTCB6SvMx4Sv39VjfSx5hiklNjZ29Gep9VoQULizt6dvuOiSUmm+88qulXwsdNHTkgK0vqXgcA1oYpdNhG+7sdf/5/+wNPK/WeOzvR9ncVFEfXBb25tYu5rc9qKjidfPWntWuge6+772tfOrOHS25eUiUKH7ju078qKrwcICPx85+cHtscNUdP9R/3paVwkZSf9mKCQE6Ig/ENPa92GPmKXiSAM54sLX/dYF5fevqQ8Jrq9d6WHhtD/SbfH29oFXFUDits7t43ujrQJNTrhCxdXO4c7A2ORkPjm2jcHLPadvZ2+xc4kv/7GjY3UiUFlJ/2YoJCTkcDFArPN6tTRHm/vv8K7WGtRF0F8vEklAsIx65pNtBotLD22pF3ATduFKE1mZpIg6p6Z7t6i+5uijcJrpI3Hr3LUC4WcDD2uFpiLFd1qtNAeb++7QVaeWrEKzYvTOdzZXxxURcEkhR6axtweb+OVj72yv9/GjQ10Dnew/PFlyLNy3wrWHZuU1OPqWkgSRFONchXR/cPIFFVt9vAaaV0kZSf9mKCQk6HHZIG51ghpj7f7hPuVj72Ca89d63ODRFk8sagV4ZBQYMIJR+XOkJDGBhG6MQuI/eqLpslMJ65Pf+Bp3Nm9M3DeKA3RcHKvJAmiS1cklUXcPdbFteeuYeWpFeU1srhIVG6vKpA5akUI8RCALwJ4L4A9ABeklEumYxi1QorEVFskS4SILeJ5vZBHGxrb+OQ7hztOjS2aook9uaftTh/P3gTQV2I2bMCceI8Q+OTxT+Lc4+cS97XBFM6YtnZ5tKF0fE2gLs0scgs/FEL8DICfkVJ+RwhxP4DLAD4mpfyB7hgKOSkSk8jpxG3rzpa3eGHbsDWT4Kv2j5IlTDJK/D5tOgaFCAhcfOpiJcVQ9XxCMY9OplUnt/BDKeVfSym/c/f/NwG8BeB9Wc9LiC90r+q6xbPrW9e9+kJtFsl6V3pGF4xq/6hLCEDfQmZawtDCA587APG8MEbBxJGQhZWkdUXlXgtFvEoukrR49ZELIaYBfBDAtxTfmxNCrAohVq9everzsoQY0flkdaI3dXjKqy/UZpHM1PQ4GgUD6BdvgeSoDVvCSU432ekw1ScvE5cFzqqN3QZvPTuFEAcB/B6AM1LKv49/X0p5AcAFIHCt+LouITZ0j3WVYqzKWiwjnMy0wBpvFadL0gmtYZN7JfSZu6bFR33K9zXvw092f6LcL5oEldS/sohmxdFxqdxb8QXOsntvpsWLRS6EaCEQ8Z6U8lUf5yTEBybrqqhwMpvwR1Nxqfi5dOK7cWPDOCFMtCaw/PHl/beMpceWrK33MGqmc7iDg2MHlfuEBbxs4rRdQkJ9WMi2MeBZar6XiY/FTgFgGcB1KeUZm2O42EmKwKXIUZ7Woc1ip+1YTZEtTdHEe37qPUqhb4omlj++rL1vX0W85FlprE8e1k23XQBWPZdWo4VD9x1ybuBs8zO2GXuZ5Bm18i8B/A8AVxCEHwLAZ6SU6pxiUMhJMWQRC5/haLbioBIa4F4ST1KTBiB98wUfUS9h1IpuYog+d1M0TDTE0sb94/NnVYXCWCbyjFr5n1JKIaX8OSnlP7370Yo4IUVhu8CV9+u0zm3SEI0+d0F8gfVPNv4Ez7z6zL77wcafvbO3g/vH7rdyF0VdFgtvLGD24dmB9H0XwqgVGzeG7pkICOdG1T5/VlVOwzfBzE4ytNhm8OVd1U4XSbIrd7X+4d6VHs6vnk9VAyXa2UcXcaPyUS+/uYzFE4uQZyXu/EdzNqeOjRsbVmsPiycWlan/aWu++PpZVTkN34S3qBVCqoauu3vcurKNaEhLKAKhi0SVZRlaldF904qazbhNbyFhBqWqIqLttXVRQlF89kLwWYHQZuxVgxY5GVpsrasiXqejbpM9qV40i1qVaS3MpHEnlecNr5tmInF5ZgtvLBibSERpj7f3f4aqIl51cH3kDYWcDDU2iT1Fv07buHzSWpjRYmBxbMrzTo5POtVhD7NRXZ+Zy0T19Aee7msC/fKTL3v5WdUx8UcHW70RUjBJdT9mjs7gy9//cuoeloB76CIQLHA2G01jv89okS1VZI1tOKDLZJFHxIhrpJKP8FQf52DPTkIqhKkSn4qwvraLuMcFMKkAVkM0tG4fQB2/Deh7egJ6gXcJd8wjhtslzNBHeKqvEFcKOSEV5Mjnj1iJc1ih7xO//wlr33JcALO0rmuPt3Fz++ZAfLqu5K+qgmS81G3vSg8nXz2ZeO08LHKXxB8fseW+4tPZfJmQimFKt48ThvUduu+Q9fnjfvakglq62PEwLT/ucjH19Nzc2lRWGzy/en7fF9091k2s1pjXQqZLcwkf4al5h7hSyAkpifnX5633DQXm+tZ1q/11nXNMbdBUHXgEBGaOzngTnHipW1MHpaZo5rbo7BKp5KPpct6NmynkhJSES+XBpH6WUUwCaGqDdu7xc5h9eLZPWCUklt9c3u8eFKc93nYum7t+Y72vjrrOb78n93KLHHKJVPIRnpp3iCt95ISURFJHoCjybPB3arNImGVxUOfLNXVNAmDl61aRR2s93/Su9DD/+vz+GNvjbSw9tuRceC3PqBVa5ISUhMrFoaJzuLMf8/zMq89g/MC48Vid1W4TN61zoVzfuo7Zh2f3/ehN0cTsw7P7WZBpOxOFE0NV65uoGmJv3dnS7mcqy5tn42YKOSElsfTYEg40zFUyJloTmDk60ycSm1ub2LqzhRPvPzHgX9YJoG39b90kMDk+ieU3l/s6B/3On/4ODv6ngxDPC6UVP9GawKnjp/bdFzp8t9bziW1BtbLrmFPICSmJ7rEuDt93WPv9UNAuvX1JKRLf+D/fGPAv63ps2gqNzpcb7h9le3cbP975sfJ67fE2xg+M4/zqeQDAxacu5tpaL68sTdtok7yjUpKgkBNSIqYolI0bG8amD6pFwlvbt/DMq89APC/6BM1WaFSRLTofto7Qn765tdln/c8cncnFheLSbcgV22iTvKNSkqCQE1Iipj/0UJRMbgndcQD6BM1VaKJ+4M2tTacxqGLIb+/cxqW3L+XiQsnTrWEbbVJ2HXMKOakkw1TQyMTiiUW0Gi3jPmG/zCi2whoKmovQqIRRNQZXwqSmtC4U3e9Enm4N2zDFsuuYM/yQVI68W69VhXhYWxLRFmgzR2ew/OayU60S2/A3mzZsk+OT+Lutv8Me+sMcx5pjuH/sfu/hhLrfidmHZ3Hh8oWB+u5Zr1dVGH5IakMREQBlW/yqsDYToSiFluy5x89pszTjuPppdftHx3DtuWv44lNf7Lt+e7yNl598GUuPLXl3M+h+J86vnleKeFXCF4uCQu6ZXg+YngYajeDf3nB6BHIl7wiAPBfHbFEJkw6dKMWzNIFBl0t4rMs9zxydsQprDK8vz0rIsxLXnru2H1eexs1gmlx1P3vVm0Oeqf1Vha3ePNLrAXNzwO27f5/r68HXANAdnd+pzOTdei2pzVkR2E5KYdVDm3EJiP1U+mip2e6xLqZfnLa6596VHpbfXO4TSAGxn/xji2u7tLjrJJxownPpfidU5JnaX1VokXtkYeGeiIfcvh1sJ/bkHQFQdswvYDcpCYi+BUGVxRq3tMNkoYtPXew71vaedQudl96+lOY2rUlyp6l+J3SLr0WF/FUJCrlHNjQ6oNtO1OQdAVB2zC+QXFIW6B+PzjUy//q81XqC7T2XNcklXVf1O/HJ45+sbGp/0dC14pGpqcCdotpO3Mizk/niiUVlBESRAhDem65LUHw8OotV52ePC6PtPeft1tJhc13V78SjU49mLkQ1DNAi98jiIjARM7ImJoLtpDqUHfMbHcfamTXIs3I/hV03HleLOC68tvdcVmJL2uvmWYiqTtAi90i4oLmwELhTpqYCEedCZ/XI0+JPQhfPbRqPy2JfVADj15o5OoNLb1/aT/8H0Hfd6JtCkVZuWdcdFijknul2KdxET1J0hg6VaySOgOgTQNW1Xlp9aX9/3bXLnORIOuhaISQnVFEmaZOdQteIqa9m3L1gE6teZKlVEz5i+8tO8ioTCjkhOaATJp17xMYH3j3Wxdwjc8rvzRydSXVOm/2KEMis2bxVSPIqE7pWCMkBnTA1RVOZUm4bFaKL51Ztt/Wrx68d9atPjk/i5vZNbO9uA7B3BbmSNeyxCkleZUKLnJAc0AnQrtzNFBXiIng2serxa6sSjEIRD7m9cxvzr88bz+tqxWeN7a9CkleZUMhJbahTHRtT4aksoY8ugqcKOYy2XlNd27YGzObWplac07g5soY9ViHJq0xYxjYDvR5DDYsiXscGCGL0L1yo5jPPqxRv3iV+TSVs4+jKxE6/OK106SSVlc3SZX5USh+zjK1nQmFZXwekvFcgq6pWYp2sWRV1q2OTV9JRWeULVLi6M5LcHFmSe6qS5FUWtMhTMj2tTsfvdIC1tWzn9m3p182aVdFoBBNmHCGAvb3B7UWSxZKsGirLVodvi5wkQ4vcM3kVyMrD0q+bNatCV6+m7Do2wxb2pvOru/ivy+5fOZJIKTN/AHwEwF8AeAfAbyTt/8gjj8i60+lIGUht/6fTqd55hVCfU4hsY11ZCcYlRPDvykq28yVda2Kif/wTE/le04bOCx2Jz2Lg03mhU+7APLPyvRXZeaEjxWeF7LzQkSvfMz941/2JHQBWpUJTM7tWhBBNAH8J4MMA3gXwbQC/KqX8ge6YYXCt5OWuyMOFoHMDNZvA8nK68Zbhrqni4rJucTDsk0mIT/J0rXwIwDtSyr+SUm4D+BKAJz2ct9J0u4FodTqByHY6fkQsDxeCqiojAOzupnfblOGu6XaD9Ye9veDfskUcYNgbqQY+hPx9AH4Y+frdu9v6EELMCSFWhRCrV69e9XDZ8slDWPIohRtOOk1FmY4k8dVFu7CJRgD9waQK+BByVb+lgXdNKeUFKeVxKeXxBx54wMNlh5O8LP1uV++a0YmvaeG1qouPRTPqYW+kGvjwkf8LAJ+VUv6bu19/GgCklL+pO2YYfOR1xDVk0rT/4mL9QxoJqRt5+si/DeCoEOL9QogxAL8C4L95OO/IUFSyjqvbxuQ+yevNgRDiTmYhl1LeAfDrAP4AwFsAviyl/H7W844KRWaIRsUXCHzmoY9cdb0k90kVFx+LYJTrXpNqwszOkskzQ1SHbejgMGSE+mZUanqQaqJzrVDIS6aM1HOXyaOKsdtlwvRzUiY6IWdjiZKZmlKLap7RHy6hg+xBeo9eD1j/fxvKOK1RqXtNqglrrZRMHnHjIbpFVIYOuhO6mXCDCUCkelDIS8Y1+sM2wsW0iJrn5OGTKpXe3c9kfWMR2GYCEKkYqgIseX+GoWhWGrIWmXIpHJVUfKvIgldpqFqRrL7CY8dWJM50JM4KiTMsCEWKA5qiWRTylLgKoQ9h0olzszl4nrwqHhZFXtUlh2U8ZDTRCTldKylIE/vto8iUbpFSVfyq7n7wqtVyqYs7qkow3r44KOQpSCPKJmEy+YKj32sYflrx69ddeHQTTqNRjq+cmaxuDFvDjarDOPIUpIn91sVut9vA1pY66QYYTMgxEb9+lWLAXceiSkYKGfWkpDrAePt8YKs3j6RxW+gsZEBv3assf5dx+Uihj74RHDkSfFyjSEyuKN3bSJbSu6R80jZgJumgkKcgjdtC92p+/bp6/40Nd3/wzIzb/knEBXhzM/i41oTRuaLm581rDWlK75JqwIYbxUIhT0Faf6nKQjZZ964Lk5cu2e9rE6Od9EZgaxnrRHdzM3mtoe6LtqMKG24UjCqUJe/PMIQfuqAKVQy3hSGBqrBEVcii6WMbWmgbCqkLYXS9pi50z+acVYsnJ/awAbN/wDjyclAJ0diYlK3WoHiFcclRkYoKftLHNqbZNiba5rrtdnI8vU6M2227cVQ9eYmQoqCQF0RcdHRi5SrENtaxrcDZJgslvRG0WsGkZGMt695KaG0TYg+FvABcXSE6MXdJuY9+bHHJUowK8MGD9yaBZjP4OsubQfz8tLYJMaMTci52esQ1XFCFLhpEFSkTxzYk0CXqJlygvXgxWKSVMti+uwvcuqU+v0tEyah2GSLEJxRyj9gK2NgY0Grpv6+KBgkjZYSiFnaIbUhgUtSNKqLFZZJiRAkhBaMy0/P+VNW1kvU1X+eyUC0IJi1i6qJBVlYGF0qzuDZU51f5rW1dQ/RxE5IfoI/cjI+FtzTnSFNVL+0kYINrqGA4UcUnP/q+CfGPTsjpWrmLj+qEaRKF8ihuNTlpv2/cjaKqB5PEwYP9Pm7b6pDhtYUADhwI/i27gQQhtUSl7nl/qmiRl1m/28V6XVkZDPlTWci2142/QdiEOereIJJCLqNvGaYIH7pnCFEDjUVem+qHeVfyc+ksXyZHjgSp7SZMVRiBe88yjfWdFSGCn9+tW+b7qNpzJ6QK1Lr6YZpGDq7UpX53kogD5prd0WdZBuHPL+k+bCOAqtTXk5DSUJnpeX9cXStFtdmqwwJd1uiRNHVP0rhasn5sfrbMDCWjBuq82FlU26+iklNsrEjdPu223TV0C7Wuz0zKe4u3qtrgOtrt4Lg0CGFXktfHAjUhQ4FK3fP+VNUiN+HDWl9ZUS8Exq1Ik6W5shKkx9ta03Fcar9En7NLNcbo/ZjeAML4eptrq34GLvdNyDCAOseRl/0KnVeMuW5SMk1cNglBponOJOQmN0p0IjGJaLPZ/1xM5wz3S5pcJiakPHXKPsKGne3JsFJrIZeyXP+1jzeCJMszakWaQiF159HVNI+T5PM2jTN6v7bhmrrz/fRP3/t52kxKurcQ2/smZBjQCXktfORAucWVfPjok/aN1icxdcXRnUdKu0QkUx2UTsccpRO9tm3nHlU0UKsF7Ozci0KyYXdXvd32vgkZZmoj5GXio92Yad94mKMpFFJ3njDuOpzoAPVi6eJiULQrTqsVfM+0UBi9tmqMQgTiHM3SBAazXQ8dAra39ddRoVtojd83RZyMJCozPe9PFTM7TeTpI2+3A/+vrqiWqoaJaSy6BdVoB6JTp/r3abfvHW/jeokvPpr81WlbyNn4yOlGIaMG6u4jLxtfUSs+uuS4iLytuIbYxJnHj086ptnsH6tu/3C/dpuFuAhRQSGvKEkLqS7i5ZLsE4pl3Cq3nQyi9VxcLGwX6zoq+uFiJwWcjDIU8opiiv5wtdazZmG2WoNuHdP+NrHiukkqaYKyLapFK52MErkIOYAvAPhzAN8D8BUA77E5ropCnrcg6M6vi6E2Jcvowh7T1BJPOndSQo+UgYXtMonYJOzYJgrRb05GibyE/JcAHLj7/98C8Fs2x1VNyPMWBNP5TUJuErLwvOHk0G7rmyG7fOIiu7Ji3l/lJgk/jYbbRBQlaWIwxdQzIYgMKzoh91bGVgjxcQD/VkqZGACWpoxtnuRdwtZ0/o2NQH7iCKHeHrKyElQxzNrsOU67DVy71r/NVDq32VTHeLfbwM2bg2GGrRbwyivJYYJJTS6Snp2pjC8hdaWIMrbPAnjd4/kKQ3bw8T0AAAwjSURBVJdks77up0SqKaFIFxeeNL/aNEPWxV632+pYciAQ7NOn+7ctLemvoUvU2dxUx4ofOmQX662KU49y65a+E5JLfD/L4JKhQGWmRz8Avg7gzxSfJyP7LCDwkQvDeeYArAJYnZqaKu5dxALXtHdXf3pS7RSX5sbAvTA9m311Lp1Tp8zHxYt4mcaSxXWjQhWtovq0WoPdklwWQuljJ3UDeUWtAJgF8L8ATNgeUwcfuU4o2+10cd9JSTwuC5XxhB7TRxeTbbOYGGKa6E6cUE94Nu3ebJ9V0v2lSZwy3Rd97KSq5CLkAD4C4AcAHnA5rmpCLqV9iVQb4VOdT5W9GcfGyj5xItjXpRytaqKxOc4m21M1AZ46ld7aTfvs4+e1Eeky+7QSkoa8hPwdAD8E8N27n/M2x1VRyOOk6aQTkqeIhULkGjMen2hsXCLhmE2ZmKZrpQnpTBsLH3++NiJNi5zUjdxcK2k+dRBynRjbuAzSCoSNWyEUoiwTjZRuE4DuWdheK6/s1Kw/A/rISd2gkKcgbW2ULK/spthy4F4STpaJRkq3RUrds8hDLFdW0lvlad6KTp269yyazeBrQqoKhdwjSRZmllf2JKs8WuPEdqIZGxtc8LQVx2bTbaw+FhRVmaLR6o22kxWjVsiwQSEvkCwC4dJJyHT9aMZnvDVcq6XPutRZ5Lrr6MrhhqR9OzGJsC8Bpo+c1A0KecGkrd2S5FbQhdvpyFp/Jd6DM3p/NpZ/XmLpozYOo1ZI3aCQ1wST8LZaUh44MCieJhGz9Tc3GoPJNTprd2XFPupFV7I2DMcMJ4tQ3E2+c99FzWwnGVZYJFWBQl5BbH3coSWuK4oV9ZvHsbXIw7K5NiGFLgk7qpK1pkJburrkvppvuJ6XfnRSJSjkFcMkEDoRMgmm6TpxH7lJqJPcDVnj611i5UNc3TMu4pvnwjUhvqGQe8THq3YagUgr5DqXiUrkksaVJTQw7hYyfbKEc/oU36L96HTjEBMUck/4etU2CaKudohu/zSulXgfTdv7052v0bCz/G0/LpOL7bNNI75FWuR045AkKOSeSPrDtrWoTIWoVH/IJpdEmsVOk6iZwgptXUI+xDx8prpKjbrkHV8/o6T79Q3dOCQJCrknfPXYVO2rO3eSOJowdSAKxxGPHlHFnrv2yfQh5El++TQ+8rwWTn3AcEiSBIXcEyZRSSM4UYEwCVoaay0p3V0l2DYWchK2maNJbxu2C7A2zzatm6ZIqjw2Ug0o5J4wWXRZLSrTH3KamiWuDStsLeS09wEEIZSutcPrvHjpAn3kJAkKuUfysvaS/pBdXvGzZnRmuY80bqAiUvJNz6UqVi+jVogJCnkB+BAcX3/IvhYco5+wqUUSeYilbYKPjyQgQqoKhbwgVIuHroKcZ5x6lk+jYTcu3UKuKspEd6+uz8BnEhAhVYVCXiBZrD5fFqPqPKpmxVk+OnGWUl2K1jb9XVefJY0rqSouE0J8QCEvkCyi4lOQdLVc4nVPsoh7GHaZ5j5MyUquz6DKi5iE+EIn5CL4XrEcP35crq6uFn7domg0AhmJIwSwt5ffsWnp9YD5eWBzM/h6bAzY3rY/vtMB1tb6t9nch24fHaZnMD0NrK/bjY2QuiKEuCylPB7f3ihjMMPO1JTbdl/HpqXbBa5du2fH/uQnbsevrweiPD0dTAqA3X3o9mk29cf2esF14tdbXAQmJvr3FwKYmbG8CUJqDIU8B1SiMjERbE9zrBCBWEaFyzdxgWy33Y6XMhjj3FxwLptnoNtnbk69fWYm+N76+uD1ul1gdjZ4VtExLS/n98wIqQwqf0ven2HxkSfFPqeNjIhGviQtGPoYr2rRcWwsfRGs8NxJbeCk1Dc/Vo01ye/OBU8y7ICLnX4pIh7ZhzDZTAq664Rt5dKIuc2zcX2GSQuaXPAkww6F3DNFWH9JBbrSJL+oxptWIHUf26gT12dIi5yMOjohp488JRsbbtvToFsMnJzU+4qjzM8Dt2+br7Gxkbww6bLQOjEB7O7qr2X6Oml7kt89y9oEIXWGQp6SIqJLdMIEDAr07dvAwsK9r3u9e+GEJqam0gkkECyInjoVhPgJEfx74ULwr+5apq+Ttne7984fvV63a/d9QoYWlZme92cYXCtF1exQuVBsfME2vm2bgly2JQeixyfVMy/6GRIyLIA+cv+UVbPDxhec5Nf2Wb9EF/XSbtvVZWHdE0Ls0Ak5MztrSK8X+MSj7pWJiX43gi7Tsd0Okn9ssM2WZFYlIcXAzM4hwsYXrPN7Ly3ZX0e36BgX7TwWfnUZnISQQSjkNaXbDazdvb3g3/iCno+FP92ioxD3hLXXC8TW5fgkwjeOpKgcQkgAXStES68HPPNMIKZxOp3A6o+7eELirh4X6KohRI3OtUIhJ0aitUvi26em1ILbbAY1TtKG/ZVRAZKQOkAfOUmFKSZc5wPf28sWu11GBUhC6gyFfEjxtVhoShbKS3CZoUmIGxTyIcTnYqFp0TQvwWWGJiFuePGRCyE+BeALAB6QUiZGKdNHni9FLhb2ekFpgLBmy+IiBZeQvMjNRy6EeAjAhwF4LBdFsuArrtvGPWMKg2QsOCHF4MO18gKA5wAUH/5ClOh81FLaC2pW9wxjwQkpjkxCLoT4KIAfSSnftNh3TgixKoRYvXr1apbLVpoqWKG6aoWAWVCjY5+dTa6waGJhIdvxhBAHVAVYoh8AXwfwZ4rPkwC+BeDw3f3WABxJOp8coqJZcXxW88taTMrU+UfVbMGmCYVLtx126yHEP/Bd/RDAMQB/e1fA1wDcQeAnf2/SscMq5L461PicEGwF1balm+292D4LVj8kxB6dkKd2rUgpr0gp/5GUclpKOQ3gXQD/TEr5N2nPWXd8LTL6dEvYxnrbjNEltNAmNJF+dEL8wDhyj/hKkPFZTdA21ls3xmYzXSy3TSw4/eiE+MGbkN+1zC0rXQ8nvhJkfGZM2ibX6Ma+vKyvsGhzbVOFxiL6nhIyCtAi94ivjETfGZNJgupz7C6wpgohfmD1w4oyChmTNp2OCCH30GV2HihjMCSZbnf4xSy8v2GfsAjJGwo5KZVRmLAIyRv6yAkhpOZQyAkhpOZQyAkhpOZQyAkhpOZQyAkhpOaUEkcuhLgKQNHDppIcATBqGaujeM8A73vUqON9d6SUD8Q3liLkdUIIsaoKwB9mRvGeAd532eMommG6b7pWCCGk5lDICSGk5lDIk7lQ9gBKYBTvGeB9jxpDc9/0kRNCSM2hRU4IITWHQk4IITWHQm6JEOJTQggphDhS9liKQAjxBSHEnwshvieE+IoQ4j1ljylPhBAfEUL8hRDiHSHEb5Q9niIQQjwkhPhjIcRbQojvCyHmyx5TUQghmkKIPxVC/Peyx+IDCrkFQoiHAHwYwCg1IfsjAP9ESvlzAP4SwKdLHk9uCCGaAP4rgMcA/CyAXxVC/Gy5oyqEOwD+g5TyHwP45wD+/YjcNwDMA3ir7EH4gkJuxwsAngMwMivDUso/lFLeufvl/wbwYJnjyZkPAXhHSvlXUsptAF8C8GTJY8odKeVfSym/c/f/NxEI2/vKHVX+CCEeBPA4gN8ueyy+oJAnIIT4KIAfSSnfLHssJfIsgNfLHkSOvA/ADyNfv4sRELQoQohpAB8E8K1yR1IILyIwzPbKHogv2CEIgBDi6wDeq/jWAoDPAPilYkdUDKb7llJ+9e4+CwhewXtFjq1ghGLbyLx9CSEOAvg9AGeklH9f9njyRAjxBIC/lVJeFkL8Qtnj8QWFHICU8l+rtgshjgF4P4A3hRBA4F74jhDiQ1LKvylwiLmgu+8QIcQsgCcAnJDDnXDwLoCHIl8/COD/ljSWQhFCtBCIeE9K+WrZ4ymARwF8VAgxA+CnABwSQqxIKU+WPK5MMCHIASHEGoDjUsq6VUxzRgjxEQD/GcC/klJeLXs8eSKEOIBgQfcEgB8B+DaAfyel/H6pA8sZEVgnywCuSynPlD2eorlrkX9KSvlE2WPJCn3kRMd/AXA/gD8SQnxXCHG+7AHlxd1F3V8H8AcIFvy+POwifpdHATwD4Bfv/oy/e9dSJTWDFjkhhNQcWuSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJz/j9xu//snpe/YAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = plt.figure()\n",
    "for i in range(200):\n",
    "    plt.plot(samples1[0][i], samples1[1][i], 'o', c='b')\n",
    "for i in range(200):\n",
    "    plt.plot(samples2[0][i], samples2[1][i], 'o', c='g')\n",
    "plt.savefig(\"fig1.png\")\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualization Clustering Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAD4CAYAAADxeG0DAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2dX4xc133fv7+ZnZV2SYkGhypsWNpZP/ChShk31sJo4YcGZpNaK9ly/BAkHRGEBGMhsgFWQA2hzj7QMjBIIwORCTSUsFAlsNoBDKGR46aikMRq89CgMbJyLTO20toIuLTdBCGXKCWKG3O5e/owvLt37pxz7zn3nvtv5vsBBtTeuXPuubPa7/nd378jSikQQgipL42yJ0AIISQbFHJCCKk5FHJCCKk5FHJCCKk5FHJCCKk5U2Vc9MiRI2p+fr6MSxNCSG15++23ryql7oseL0XI5+fnsb6+XsalCSGktojIhu44XSuEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEEFJzKOSEFEz/Yh/zX59H49kG5r8+j/7FftlTIjWHQk5IgfQv9rH0R0vYuL4BBYWN6xtY+qOlwsWci8l4QSEnpEBW3lrBze2bQ8dubt/Eylsrhc2hKosJ8QeFnJACuXz9stPxPKjCYkL8QiEnpEDmDs05Hc+DKiwmxC8UckIKpHe8h9nW7NCx2dYsesd7hc2hCosJ8QuFnJAC6R7rYvWzq+gc6kAg6BzqYPWzq+ge6xY2hyosJsQvopQq/KILCwuKbWwJKY/+xT5W3lrB5euXMXdoDr3jvdwXkzKuOW6IyNtKqYWR4xRyQohvoqK9eHQR5985PxRknW3NFv40Unco5ISQQgjSG8OiLRAojGpN51AHl56+VODs6o1JyOkjJ6Si1LVoR5feqBNxYD9Tpq73WhUo5IR45vQbpzH11SnIs4Kpr07h9Bunnceoc9GOSxrj3KG5Wt9rVaCQE+KR02+cxgvrL2BH7QAAdtQOXlh/wVnMXYp2qmbNmtIYBTL0c5ApwwKl7FDICfHI6turTsdN2BbtVNGaNaU3PrXwlDbtkgVK2ZnyNZCINAGsA/iZUupRX+MSUicCS9z2uIm5Q3PYuD66YXrU2o2zZsvKBgmua5tqaHuvxIxPi3wZwLsexyOkdjSl6XTchG3RTlWt2e6xLi49fQm7Z3Zx6elLsYsKC5Sy40XIReR+AI8AeMnHeITUlaWHlpyOm7CtADVZrYdnDhfiN/fhn69CtWvd8ZJHLiL/GcDvALgHwJd0rhURWQKwBABzc3MPbWyMPkoRMg6cfuM0Vt9exY7aQVOaWHpoCeceOZfLtXQ529PNaSilsL27vXcsj+Ib3bVZ5JMvuRUEicijABaVUqdF5JdhEPIwLAgixB/RKsobt25gc2tz5DzfxTfzX5/X+rZZ5JMfeRYEfQrA50TkEoBvAPi0iKx5GJeQsSDv9MCoP/ra1jXteb795lX1z08imYVcKfVlpdT9Sql5AL8B4L8ppR7PPDNCxgDf6YE2i0JRbWrL9s+TfZhHTkiO+Cx2sV0UisoC0V2n1Wjh/VvvVyqvfRLwKuRKqT9jDjkh+/h0P9guCtEskPZMGzNTMzjx+glnCznuCUCXbXLvXffi1s6txDkSv7D7ISE54jMg2Hi2oW0+JRDsntnVfiZLZkmaz6aZI7GH3Q8JKQGfbo40vu8srp00n63CNnJV6z1TBBRyQnLEZ7FLmkXBxrVjEr40bqGyqzSr2HumCLz1WiGE6Oke63opkHHtYdK/2EdDGto+L4GFHHWfBMIXnOPaA8V1jsEcfG0BV8XeM0VAHzkhBVPE3pU6/3ZA2M8d58PvHe/lXrnpuzp03H309JETUgGKevTXWabAoHlXWCTj3CdxbiFffmjfvcir4KMvAwo5IQVS1CYKJoHeVbtDlu7hmcPa84Lj3WNd9I73MHdoDpevX8bKWys4/cZpb4uR7+rQsn30ZUEhJ6RAdG4MwH9Zuy/LVPcE8eL6i94WI98W9KR2UqSQE1IQ/Yv9ke3OAnw/+ttapqa+LMFxl42UN65vOLta8rCgXXqhjwsUckIKYuWtFWMgzvejf9Z+5sFx1ycFV1fLpFrQvmHWCiEFYcqoAAB1pvi/QyA5a8SU1WID29n6h1krhJSMyfrtHOoYP2PKDvGVNZJkEetcH7ak8ftPYlWmD2iRE1IQrjnTpvNPfvwkzr9z3mocHznrwRgmy7wpTW3RkatFzh2Hkslth6A0UMjJpOIirCa3hkk42zNtHJw+uDf24tFFa8G3nXvWhSU6Xhk7G9UZCjkhOeO7YjPOp26DQLSfzyKMpnt0vfe4ylPdfYxDVaYPKOSE5EgebgFXi9yWJGEMu1KCawUl+75cHC5BVFrk+zDYSUiOLL+57L1i05RjvfTQUuoAJGAOuvYv9nHkuSN4/PXH90Q2WDB8txKwDYROQlWmDyjkhGSkf7Gv9e0C9oKly9YwZZSce+TcyPH2TFs7brQAySSMwROF6T4Av60ETItJe6bNnPIU0LVCiAVxPuA4N4GNW8CHWyYuEHnhRxcSfde2rg5f/upxylApoptlgMm1wn7khCQQ17O7e6ybuNGCbrxotkbWHtpp+oCHsX1yaEhj72khC1nnWxWS/t8oClrkhCSQtO+m6f32TBtXn7k6dKyq2RouwUeB4KmFp3DukXPGc6KL1eLRRasng7rhc09WGxjsJMSRwG+d1LHQFJQ8+/DZEd+3LihqIq8e2jp/vEsFp4LCi+svGgOfuo6JL6y/UOr2a3lVjPpuw5sWCjkhGsJiZCIQWlNQEsCIoMUFE8NMN6dH3DI+xMi0sQWAvXuwQUFh5a0V7ZxMm1qEyaMHu4k8N/OoykYWdK0QoiHJ1SAQvPqFV2PdA1kaTrUaLbzy+Vf2xvcVHIzLTd9Vu7EVljpmW7Mjc7J94ijKdZSn+6PooC1dK4Q4kPRorKAS/1CzPF5v727j5DdP7lmNWXYWClvNpoVlR+3sWavv33ofrUYrcdymNLVzakoz8bPA8AbQUave9pgNebo/qtKGlxY5IRqSrGkbay6LRR4QWHcnXj+RalNhl+BqmHDfltnWLD7Y/mBkXnFj2rwfdj+Fz51uTkMphe3d7b1jrUYLIoJbO7dGxkgSzaIDknlCi5wQB+KCf7bVhroxTDsEmbi5fRPLby6jIfo/1SRfrI2/Wse1rWt7u+zc+O0bWPvC2ojVafKnh2MEJgIB1s3v1s6tIREHBk8oYREH7J9IJmEfTwo5IRrCj8wA9twF0UfnoKxdnhXIs4Ijzx3Ze+TXPXZ/+mOfdp7L5tamtreKjRildR9EFwjd9mlxAtk91o0V+uD7y+reCH/e5Hrx4f6oep90ulYIcSCcH3145jCu//w6bu/eHjonGqgM48PdAgwWlvO/dl7bf3z5zeW9YKWpA2IcYZdFUtVi3Ps21aYNaWRqABa4R3TXCu49a8OvKlWhsvshISkJdwO0FcawwISFLknEbbM+dL7x/sU+nvzWkyMuCFfaM22cffgsgFH/dZrWAdHCoGjvchsaaGAXw/frsiVdFuGtko+dQk5ICtIGC4P0RJfPBpZjWPg2tzZx49aNkXN1VaNZrP3oAjXbmsXM1Iw2DdH0NGBDUmte3UJ5oHVgxEcerS616d2eVnhNY5fRJ53BTkIcCHyij7/+eKpg4dyhOadAY1AAFPVF39W8y3q+Nta+Dp143ty+acwl31E7qQtqTD7xXbWLzqGOVjD/4fY/jDxlKChc+NGFvZ9tCnCCa7v6u6tS9BMHhZyMPa5/uDZVnXG0Gi30jvecAnkvP/bykG85mK9JTK9tXRuZbxzRTJNw8DbNLkRpKzPjRNH0fZl86OHzbVoMzB2aS1XlWYesFwo5GWvS/OG6WNKtRgsHpw/u/dyeae8FOm0ttiCLI7qxQ5zAhsdOmm+r0cLi0cU9l03nUAfnf+081BmFS09fMmaXtGfamG5OG8dNk3ESJ4qm78tUYBQ+P5plZOrDnqawqipFP3FQyMlYE/eHa7LUkwQqEInOoQ5e+fwreP/L70OdUVBnFK4+c3XvD7x3vJeYNx4IjM3GDtHPBMTNtz3Txhc/8UWcf+e8cTEzieuv/8Kv4/bObd2wAAbuDddUvDhRdNkRSWcRB24pdUbh1S+8qr1G2ipPXfpllcgc7BSRBwD8JwAfBrALYFUpdTbuMwx2kqKIC4Lp+oSsfnZ1L0NFh2sqmzxrFvLwWDaBSoFoUwDTttkNB/90aYRx30OY6eb0kGsoC742dzaNaUp3rEuVZ25ZKyLyEQAfUUp9V0TuAfA2gM8rpX5o+gyFnBSF6wbGgbj6yhu2TV2LE3zd+WF0mTWtRgv33nUvrm1di3XRmBYHwC4TJECXRVMFbLKO6rQzUW5ZK0qpv1VKfffOf78P4F0AH806LiE+MD2uxwXQfPpEbQJl/Yv9WBeM7vywSwgYbkHbkAa2d7exubWZKMSBq+XJbz25V6E69dUpyLNibAugw7ZbYtGY4gdNaVbW350Grz5yEZkH8EsAvqN5b0lE1kVk/cqVKz4vS4gRkyibAnzhHuM+fKI2i8LKWytGwW3PtEdaApj6iQeLxq5yz22+tXNrT4yDRS5txWWVytnj0h2jv9sqzdsVb3t2ishBAH8A4Gml1HvR95VSqwBWgYFrxdd1CUmie6yrFWKd+6SMlLK4QFvUXaHbYSicdZGUbZOmZN/mM+2ZNgC7PSyL3KzYVE0bzZCpyt6bafFikYtICwMR7yulXvcxJiFZibOwikops0l/NKXdRZ8a+hf7RhfG5euXEzMvOoc62D2za70LUICCgkBi+4wHJf1J6X0u6aA+LGTbHPAs/d6rQGYhFxEB8B8BvKuU+r3sUyIkOzaCEbhPXv3CqwCAE6+f8P5IbSMQLmJjoiENHJ45bHw/PJ7L/pwBCirW1ZLUzTA4biuYut/fE3/4BI48d8RJ2G0X7KrsvZkWHxb5pwCcAPBpEfnendeih3EJSU0WwfC5MbCNQMTt+RlYpEeeOxKbCrijdvDez9/TFvBE/ezR67Vn2lY7ApkQyN73lVTObiuYut9fOIDr8nuyiXfUoQw/Dh9ZK/9DKSVKqV9USv3TO68LyZ8kJD+yCIbPR2qTEDSkMWRZRsXmzy//OU68fmJvgbHJCtne3cY90/cMLQhrX1gbKlLScXD6IL74iS86u1wCgo2YgeSnC9P3cXjm8JAbxSZ/3efvqQ5l+HGwspOMJbYWVt6P1CY3RniPzKhl2b/Yx4vrL6bqgRLe2cdkfeqeQs6/cx694z2oM8p5FyNg//tKcmX0jvdGnhqmGlN47+fvDc3Hdg6+fk91KMOPw1vWCiFVwlTUE7WwbLMa0hIIQVxlYWBZhs9NI+KA3bzjnkK6x7o4PHPYOS882vckTgCjRYjRjTmA/QBr0vfg0/WRNO8qQ4ucjCW2FlYRj9Rht4kpxztsWaa1MpPmHWSBmNwWeV03zMpbKyP7cZoIdvcJ/PhRS75Oro+8oZCTscUmyFX0I7WNyyetlRluBhbFpjVv4Ke2scbDjcNcvi+XxUIg6B3vYffMLq4+cxUvP/ayt99TnYt/dHCHIEIKJGlvycWji3jtB69lKnnX9Q5JCiC2Gi2ISOw2cU1pYlft7hXxAHAu7HHdxSiPZlaue3BmLWDyWQDFrd4IqQiue4AGVZMu4h4VwLgGWE1p4kN3fyh2/OhGxoB5P0/ALPCuW+flsZ2ayx6cWTde9r1xM4WckIpx5LkjVuIciOcTf/iEtX85KoBxlnDgJjEJvct+nu2ZNrZub408cYT31wwvZEnkYZG77MGZdeNl3xs3c89OQipEXLl9lKAj47133Ws9ftTPHrfJxdyhudjdeVz289zc2hyxthUUXlx/cc8PHcQukvLW8wpmuhT/ZE1PLapilEJOSAksv7lsfW4gMOF9OuMw7Z7z1MJTxi3QdNk7AkndATFKuGgoIG5xaUozt6CzS6ZS1orPoipGKeSElICtNR5kbgB2f/xxAnjukXPGLdC6x7o4+fGTQ8Ia5HLraM+0nfu1bFzfGKlmNblzdtVubplDLplKWdNTi6oYpY+ckBJI2hEojDoz+Bu1CRRmCQ6a/Lk6H3k4qOmShRIdw/T5Kmy9FvblBztK6bb6S8pKKSJrhRY5ISUQZKIk0TnU2ct5PvH6CcxMzcR+Vme12+ZMm/y2CmqvhW1Tmjj58ZN7VryNr1tHkPNe1R4n0bz7HbWzN6+oSNt22cxz42YKOSElcPbhs5hqxHfImG3NYvHo4pBQbG5tYuv2Fo5/7LjR3x3GpbujyXUT9pXvqB2cf+c8Tr9xem9rOJ1FPd2cRnumHdszxce2enkV9tg2U6tKH3O6Vggpibj0w+AR3uR60OWfP3jkQXyw/cHQI7yL66J/se+U4mgi6obIy33iO0c7jG2Koksqow/oWiGkYsRloVy+fjnW/6wTjx9e/eGI5e3aV2WwT0x6wtZ7MIfFo4u5uE/ytIZts02q0secQk5IScT9sbu2c9Vxc/umcXs23bVX3lqJLdG3QZdzfuFHF3Jxn+SZo23ru6+Kj59CTirJuDU10tE73kvcmUeXAugi7kGQLoxJaGwF0HVxCXzhaQJ+Jh//6TdOoyF6+fJhDdv67qvSx5z9yEnlqPuO5jb0L/ax/OaylT866HES+L4Xjy5abzwR+KmX31ze88fPTM1ozzX1Zm/PtHFw+uDQ9V/67kvWvvQswmpyn5ju36c1bNufvAp9zGmRk8qRdyZA2dZ+sFDZFgUFQcHAmj33yDltlWaUsKht3d7aO765tanNXDG5Cc4+fHbk+q98/pWhNMj2TBunFk55dzPEpURGybMatOpQyD3T7wPz80CjMfi3P34egdzJ0/eZ92bLNugWKhMmIdRVaZ5aOKV9xLddGIPqTl3OeJTusS6uPnMV6oyCOqNw9ZmrOPfIuVRuhriF1cWaz7MatOow/dAj/T6wtATcDP3NzM4Cq6tAdzL//0qF745xRY1tS1xL2TC6KsIoNlWDtilyeabzxc0/7ppx/dujVKEaNG+YflgAKyvDIg4Mfl4ptjag9uSZCVBUN7o4bKzMQJTCfbyjVqvt04VtilwZxS1J19QFE59aeKoSmSJVgkLukcsGLTAdJ3ryzASoQt6vbqEKExUlk2Avv7lsJby2C2MZi5zNNaMZL2ldOOMMs1Y8MjcHbGjqL+aKrQ0YC/LKBOgd72kf5Yu05oL7Clwih2cOAxgUCOncIyar1eRnj4pj9HomF4wpayXPRS7tNauQKVIlKOQe6fX0PvLe5D7xVQ5bUStiHnlsWAzsi6Br170yFrkqLKzjAIOdnun3Bz7xy5cHlnivx0An2SdNS1OXDYvDLWajAtlqtHDvXfcaLf+088tKGdesK6ZgJy1yz3S7FG6iJ22hk85qjSKQIRGc//r8yPnbu9t7ueuma9fVZTHpiwEtckJyQCcsWboA9i/2cfKbJ7Vbr+k+75LiWGbKno+UxzLSJsuC6YeEFIQpy8S1E2GY7rEulh5a0r63eHRx5JhtgDJ87WiK4+k3TudeAesj5bEqPcHLhEJOiGdMwuLSiVDHhR9dsD6elOIYvbZu8Xlh/YWhnx9//XEcee6IV0H3kfJYhdqAsqGQE+IZk4C4dCJ0GVd3PJqL355pY7o5bby2bdsAU5+WMC69bHzk9VehNqBsKOSkNtSlj41JQILClbSFLK6CFS6kufrMVbz82MvGa7tYr3FuC9deNj6qeKvSE7xMGOzMAFMNi6NOfWzyCr7lGdRzSXEEzFuZpell4yPjZFKyVkzBTgp5SuokLED9F535eX3VbKcDXLpU9Gz2MQlIXsKS57hJKY5hTMJc9B6WkwaF3DN5CUseglu3RUdHowHo/lcVAXZL0odxS3uLLhKLRxfx2g9eG+mbHnePVeguOc7kKuQi8hkAZwE0AbyklPr3ceePg5DnISx5CW5VrVkXqngPkyJaLk8B47a4VY3c8shFpAng9wE8DOBBAL8pIg9mHbfqmBphZWmQlVcb3Dy6MhYdeOz1BotamLL72ExK2pvLfptV2cNy0vCRtfJJAD9WSv2NUuoWgG8AeMzDuJUmD2HJqw2uaXFpNNIJcPDksLExeCrZ2Bj8nKeYd7uDJ5NOZ/DU0+mU7xpi2puetBstk/T4EPKPAvhJ6Oef3jk2hIgsici6iKxfuXLFw2XLJQ9hycPKB/SLDgDs7KQT4LI20Oh2B26U3d3Bv2X795n2RqqCDyHX7QA74j1WSq0qpRaUUgv33Xefh8uWj29hyct9ECw6TU1hYZIA61wo3EBjAN0IpCpkDnaKyD8H8BWl1L+68/OXAUAp9Tumz4xDsDMv8kwTdA3QmoKvMzPApmYD+DoFTwmpI3k2zfpLAEdF5GMiMg3gNwD8Fw/jTgxhq3dlZSDeebgPXF03JhcKUL3AIyGTTGYhV0rdBvBbAP4YwLsAXlNK/SDruJNCkYFDV9eNyVVy7Vr1Ao9F4dJHhJCiYEFQyRSdH93vA8vL+66Rdhs4e1YvwlXM3S4T5kiTsmE/8opSRuBwa2v/vzc3zU8AVczdLhP2vSZVhUJeMnmlHJpwSR2sYu52WfT7wMb/m4wCIFI/KOQlU7TV6/oEULXc7TII4hi4zgIgUk0o5CWTp9WrywEv+gkgLVXqPb73FPNWD7jFAiBSQZRShb8eeughNYmsrSnV6SglMvh3bS2/MdbWlJqdVWqQCzN4zc4qdeqU/niaueSFae5lzVEkNJdjawpPdxTOiMLTHbX2/Qp9cWTsAbCuNJpKIU+Jqyj7ECeXMTqd4fOCVzDXrAtKnsTNnfMhkwyF3CNpRNmHGJjGaDZHrz1kRYZeIiluuGCqNveqPSGQycUk5PSRpyBN0ygfaYamc3XNr+riC9dRtbkze8cNFk0VD4U8BWlEOU6c4gJ74fcaMb+t6EJSpRxw18ClqVvjjRvlBT2ZvWOH6+bLxBM6Mz3vV91dK2ncJGmCj7rPxL2irgefvvC0Y8W5JeLGXFtTqt0evUe6NKpN5/mOwlcw8uo83yl7amMB6CP3R1qfqU644hYF03umVx7Bt6yCarqHdjv5O2SQsX7IV0Qr5PKVGgRnaoBJyOlaSUFan6nu8TzOTeNapr+4aH+ujbsjKITRtay13UjCdA+bm8lxBvY9rx/cNakcKOQpcfGZ6kQzOKaU/jNzc+7BvQsX7M6z7bioC+qGsRFU13sIj1m1oCdJhrsmlYTOTM/7VXfXigs6N8z0tFKtltlF4stHbsLWZWFKA3RxcZjcUDp3TXRMpv3Vk7Xvr6nO8x0lXxHVeZ5FUz4BfeTlkMbPHQ362Y5h6zu2zdOOu24QqLWtMo2eZyvSVS9eIqRIKOQFERUeFxHXCXlAknUM2IucrUVueiJot/2U+lOkCXGDQl4AOuGzEWCTayVM0qLQbmebp0mEw08Ezea+4Nu4RgghfjEJOYOdHtEFB5UaZLaEmZ4GWi3zOLqMkF5vdJwwm5v2XQJdsm663f0CnZ2dwbGNDX0mC8CMEkLKgELuEZOIKTUsmi+/DLzyyuC/bcfqdoGnnooXc5f9PuOybqJZNsvL8dkrYZhRQkjxUMhDZO2BbRKxYI/LsGgGQmoSc91Y584Br766vyg0m6Pn2OZ3m9ClJpqs7yjhFgBV6idOyNij87fk/aqij7zoNrM+rhvnM3eZczjgaPJ9x/nmo8FK14yUqP+dQU9C9IDBznh8lYOnycRI09s8LvjZbNrP1SVP3STkLm0HbK7NXHFC9JiEXAbvFcvCwoJaX193+ky/P3AZXL48cDv0en470DUaAxmJIjJwiVSFfh948kng1q348+J+rcF3ubHhd27AwL0S508XGfz+btyId9kE7ihCyD4i8rZSaiF6fKqMybgS+G0DgQiCeoA/MZ+b0wtb1YJ3y8vJIt5um9+Lfpe+uXlz4LsPMlyiKGW3gDD7hRB7ahHsTLORgytV6t8dh03g8b33zMHFpP4pPtjZ0fcTd8FmAWVAlZABtRDyIrrgFbkLjG3nwbQitb1tXuRcv7Nm05whYyL47uKeDOKwWUBtG38RMhHoHOd5v1yDnVXoS+2jnNy2t3dc1seBA3ZBSFMDLdeslOB71pXkJwUqbapR486J+56r8P8EIUWDOmetlN0FL6/URJMAxW3GENc10UbQ0gh5cL9Bk6xgodDNL/ydxLUnCHd4tOkEGXyHwWLquoARMg7UWsiVKrfBkg/rL8k6DQtQmv4stotMlrGD+7X9PkznNZvDi4LNdW1TJWmRk3Gm9kJeJrZtX9OMoROgtF0TbRa5uLEPHLCzdm2/D5d9SpMWOZvvhPnnZNwxCXktgp1l42Onmrhzo8E9UwaNKXgYbQEAmAOlvd6gaVeUVgu4++6BJCbdg+leGo1BYHRqavDvygpw8uRoAPnCBbfMmbm5+CBt3sFpQiqPTt3zftXNIs/TRx71K4fPd92MwRRMbbWGS+lPnRo+L5hDkjUeWP4uFrXue3Jx7wSfZ3CTELNFTiG3xFfWim4Ml7HjxsgirkrZu3SiuwMFfVJMr2bTroQ/eDUa+yLt2r+FkHGGQl5RfAmUq189EN9wsyqdpW2ynsMbWbha2La7C0WfMA4c0DfpImRSMAk5feQlE1e16lIU5FroE5TQhzeLOH9+1KetlP7zm5v783GJFdy8OfCRJxVf9fvAE08MV7J+8AHw/vuDVr5JsQBCJgqduuf9qqJFnnd6o2n8JOs1euz4cf34aTJdbHzOSQU7wb25ZqEkYVMoRFcLmTSQh2sFwNcA/DWA7wP4JoAP2XyuakKetyjEjZ/kX9a9Tp3aHzdYHFyKhVxEdm0t/lxdT/Gke7IJUCYFXhn8JJNIXkL+qwCm7vz37wL4XZvPVU3I8xaFuPHTiG2zmWwFpy380d2zqS3AgQPueeK2C2SSRe4jt5+QumES8kw+cqXUnyilbt/58S8A3J9lvLIw+Zc3Nvz4YOOafsXt22liZye5i+FgbdUT15nw6tXR+7z7bv25W1t6//7qqn5uzaZ9rnevZ96g+sYN4PBh/XtVaztMSBH4DHY+CeBN05sisiQi6yKyfuXKFY+XzY7pj19E313PtdPFZd4AAAu6SURBVDNhXEGRrvgniWYzfefHTmcQ0DTxwQejXQSvXdOfa9pww9SLfHfXvmCn2x1sUH3gwOh7m5uDVr3RwqZoYVXS74ltcMnYoDPTwy8A3wbwV5rXY6FzVjDwkUvSeKqCrhWdmyIu7c73vpw2udXhV7Sgx/aVVFxjcrGYrmX6jkw+cltXlamwSfe7MAWobb5zBktJ3UBeeeQATgL4nwBmbT9TNSFXajSrJItv2VSVmZQVY+PXDrJWbIVct6mx7T0FuCwarZZ9nrjp9zA9bX890/eZFPdgsJTUkVyEHMBnAPwQwH0un6uikEdxFfMgyJbF0nOxlNOUuQfYZMoEGSmu1woKhdKmc6YNAEfvMSkYymApqSN5CfmPAfwEwPfuvF60+VwdhNwkyCbr1IelZ5OPHQiNq+CFr+8qyi7X0gmhi6j7aLNr83ugRU7qSG6ulTSvOgi5UukaV2W19JL8w2GLVzcPG4F1yV033bPpPqNC6PqEkqWwKXyPNj7yaN59q0UfOak2FHKPxFmYWS29JKs83ONEN4+4zRzS+P/D1Zvha9n6wV2/jzgfedCl0XbMuN+T7jrT0xRyUm0o5AWRNRvCZSch2+tnecVd79Sp4WrOoOI0TJonlOhTSbTVr4+ME7pWSB2hkBdIlr4tST7iuJQ73fXTtAAwPQFEr2Fyt4TnlZdgZu2Nw2AnqSMU8poQZ5G3WkpNTQ0fS3IH2AYPDx7UuzR0fmObHjHhzZWjgj89vW9x69Ijo+TR0MyHe4aQoqGQVxDbYGpgGR88qBcfk9WslL1PPGi8lSRuLq4bnX89rrmXqSd5mgIsm6eWpHFZNESqBoW8YsSJRJqWt3HXsemKaNuIKksqoo0lb9tGNy5Yaiu+SYJPPzqpGhRyj/h43E4jEmmFPKlS0mVfTJc87/C5Bw7YV2y6FPb4+F5NFO1HpxuHJEEh94Svx+04QdT9Mcf1BU/jWonuo2l7b1nyvG1f4Wuarme6Z5/iW6RFTjcOsYFC7omkP25bq8o0TlSIbBpdpQl2moQtfC1dINIkOC6bMdu8wt+nzjVkCvL6+v3E3Wse4ko3DrGBQu6JOGF09c+6VEvGWfBxmAKY0Z4oYeHWBSR13RpNYphVxKMLjU0QNu57jcug8RE49QHTIYkNFHJPxFlOaQJzYZGIE7U0FtvaWvwC4Lo9nG2/GFs/uk3vGqXSP1VExbfKVm+V50aqA4XcE3FWXVarKu6P2dWa9F3haXsfcQtSszkQbpfeNUnfiwtVtnrpIyc2UMg9kpfFZ9PoyfYxP4+gpM19xFnjaQt+fIlc1a1eZq2QJCjkBeBDcHz9MWdpB2t6BZtaxFFWSb6vIiBCqgyFvCCSsj5cxsgjT92HeyWpnF4XxNU11Iq717yyS2j1kjpDIS+QLJafL6tRN06r5baNWpKgm8T51ClzGqXNvbpuFVd1lwkhvqCQF0gWYfEpSqZeLuFjp06lF/fwdnBp7iGuWMnlO6hyEJMQn5iEXAbvFcvCwoJaX18v/LpF0WgMpCSKCLC7m99n09LvA8vLwObm4OfpaeDWLbvPdjrApUvDx2zvwXSeCdN3MD8PbGzYzY2QOiMibyulFqLHG2VMZtyZm3M77uuzael2gbNnB8InAnzkI/af3dgYCGmjMfi337e/B9N5zab++OHDo9cCgF4PmJ0dPlcEWFy0uwdCao/OTM/7Ne6uFd8+ct2GDXnPN23Wi4uP28VHrvPvh8e09csTUmdAH7l/4jIgsmRHhDNffIqTa/57WjFvt+O3agvPx3RedK5JFaAMeJJJgELumSJykn2Ik82ikNSJMY2Yp7XGTd9fUkCTAU8yCVDIPVOEBZjUoCtNAYxuvkn34qO4yDZjxfT9JZ1Pi5xMAiYhZ7AzJZcvux1PgykYePgwsLQ0CDQqNfh3aWk/+BewvAzcvBl/jcuX9cHC2dnB8bh5AIMMFxui34vr95c0x6T3CRlrdOqe94sWuR0m94Pt3poulnKSv99Urenqy87y/fko0yekzoCuFb8U1bdDJ06+9ta0ma9tywHXAC37nhDiDoU8B8qyAG2s2SS/ts18bcU2bcokLWhC3DAJOSs7a0i/P/CJh/3fs7PA6uqguAcwVzu228DVq3bXsa2YZGUlIcXAys4xotsdiHZQidnpDIs4YA7+nT1rfx1T4DEq2nkEfvt9fRUnIWSUqbInQNLR7Q4Lt+59AFhZGQjq3NxA3OM+E2VuTm9pi+wL68rKwJli+nwaok8cQVYO4DZ/QiYFulaIkX4fOHFCL9TtNrC1ZU5vjLp6XKCrhhA9dK0QZ7pds7W9uWkWcZ2rx4UicvQJGSfoWiGxdDp669iESHar2eTSybMDJCF1hhb5GOIzUGgKmrbb+vN9iC2rNAlxg0I+ZgSBwqTyfVtMGTJnz+YntjZZOYSQfbwEO0XkSwC+BuA+pVRiljKDnflRZKCw38+WFUMIcSO3YKeIPADgVwAwFFUB4nK/XV0tSS6abnewOOzuDv4NizjzwAkpDh+ulecBPAOg+DxGMkKcj9rF1ZLFRePbvUMIiSeTkIvI5wD8TCn1jsW5SyKyLiLrV65cyXLZSuPLEk07ji5QGOXmzYFLJO6aJ0+OpheaPhdlZSX9ZwkhKdA1YAm/AHwbwF9pXo8B+A6AQ3fOuwTgSNJ4aoyaZkXx1dEv6zjhZlSmplnRnXNsNqHQfU4Hd+shJB/gu/shgGMA/v6OgF8CcBsDP/mHkz47rkLuq0e5z17ntmPZbulmMwfba7L7ISFumIQ8tWtFKXVRKfWPlFLzSql5AD8F8Aml1N+lHbPu+KpI9FnZaJuTbTO2bXqhzTXpRyfEH8wj94gp0OhaJONrHMA+J9s0drPpnsttc0360QnxiM5Mz/s1rq6VqvjI01D0NelHJ8QdcPPl/PFVkVhGZWPR1/T51EHIpMM2tqQUbHY5IoQMwza2pFKwnwoh/mAbW1IaSbscEULsoEVOCCE1h0JOCCE1h0JOCCE1h0JOCCE1h0JOCCE1p5Q8chG5AsBhS99SOQIgcdejMYP3PDlM4n3X+Z47Sqn7ogdLEfI6ISLrugT8cYb3PDlM4n2P4z3TtUIIITWHQk4IITWHQp7MatkTKAHe8+Qwifc9dvdMHzkhhNQcWuSEEFJzKOSEEFJzKOQOiMiXRESJyJGy55I3IvI1EflrEfm+iHxTRD5U9pzyQkQ+IyL/W0R+LCL/ruz55I2IPCAi/11E3hWRH4jIctlzKgoRaYrI/xKR/1r2XHxCIbdERB4A8CsAUmyBXEv+FMA/UUr9IoD/A+DLJc8nF0SkCeD3ATwM4EEAvykiD5Y7q9y5DeDfKqX+MYB/BuDfTMA9BywDeLfsSfiGQm7P8wCeATAR0WGl1J8opW7f+fEvANxf5nxy5JMAfqyU+hul1C0A3wDwWMlzyhWl1N8qpb5757/fx0DYPlrurPJHRO4H8AiAl8qei28o5BaIyOcA/Ewp9U7ZcymJJwG8WfYkcuKjAH4S+vmnmABRCxCReQC/BOA75c6kEL6OgTG2W/ZEfMMdgu4gIt8G8GHNWysAfhvArxY7o/yJu2el1LfunLOCwaN4v8i5FYhojk3EU5eIHATwBwCeVkq9V/Z88kREHgXw90qpt0Xkl8uej28o5HdQSv1L3XEROQbgYwDeERFg4GL4roh8Uin1dwVO0Tumew4QkZMAHgVwXI1vwcFPATwQ+vl+AP+3pLkUhoi0MBDxvlLq9bLnUwCfAvA5EVkEcDeAe0VkTSn1eMnz8gILghwRkUsAFpRSde2eZoWIfAbA7wH4F0qpK2XPJy9EZAqDYO5xAD8D8JcA/rVS6gelTixHZGCRnAdwTSn1dNnzKZo7FvmXlFKPlj0XX9BHTkz8BwD3APhTEfmeiLxY9oTy4E5A97cA/DEGQb/XxlnE7/ApACcAfPrO7/Z7dyxVUlNokRNCSM2hRU4IITWHQk4IITWHQk4IITWHQk4IITWHQk4IITWHQk4IITWHQk4IITXn/wM3Gx0ZYHCVpgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def visualize(data, labels):\n",
    "    color = 'bgrymk'\n",
    "    unique_label = np.unique(labels)\n",
    "    for col, label in zip(cycle(color), unique_label):\n",
    "        partial_data = data[np.where(labels == label)]\n",
    "        plt.scatter(partial_data[:, 0], partial_data[:, 1], color=col)\n",
    "    plt.show()\n",
    "    return\n",
    "\n",
    "fig = plt.figure()\n",
    "visualize(samples, labels)\n",
    "plt.savefig(\"fig2.png\")\n",
    "plt.close()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}