doc/source/examples/cd_clf_cifar10.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Learned drift detectors on CIFAR-10\n",
    "\n",
    "Under the hood drift detectors leverage a function (also known as a test-statistic) that is expected to take a large value if drift has occurred and a low value if not. The power of the detector is partly determined by how well the function satisfies this property. However, specifying such a function in advance can be very difficult. In this example notebook we consider two ways in which a portion of the available data may be used to **learn** such a function before then applying it on the held out portion of the data to test for drift.\n",
    "\n",
    "## Detecting drift with a learned classifier\n",
    "\n",
    "The classifier-based drift detector simply tries to correctly distinguish instances from the reference data vs. the test set. The classifier is trained to output the probability that a given instance belongs to the test set. If the probabilities it assigns to unseen tests instances are significantly higher (as determined by a Kolmogorov-Smirnov test) to those it assigns to unseen reference instances then the test set must differ from the reference set and drift is flagged. To leverage all the available reference and test data, stratified cross-validation can be applied and the out-of-fold predictions are used for the significance test. Note that a new classifier is trained for each test set or even each fold within the test set.\n",
    "\n",
    "### Backend\n",
    "\n",
    "The method works with both the **PyTorch** and **TensorFlow** frameworks. Alibi Detect does however not install PyTorch for you. \n",
    "Check the [PyTorch docs](https://pytorch.org/) how to do this.\n",
    "\n",
    "### Dataset\n",
    "\n",
    "[CIFAR10](https://www.cs.toronto.edu/~kriz/cifar.html) consists of 60,000 32 by 32 RGB images equally distributed over 10 classes. We evaluate the drift detector on the CIFAR-10-C dataset ([Hendrycks & Dietterich, 2019](https://arxiv.org/abs/1903.12261)). The instances in\n",
    "CIFAR-10-C have been corrupted and perturbed by various types of noise, blur, brightness etc. at different levels of severity, leading to a gradual decline in the classification model performance. We also check for drift against the original test set with class imbalances. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import tensorflow as tf\n",
    "\n",
    "from alibi_detect.cd import ClassifierDrift\n",
    "from alibi_detect.datasets import fetch_cifar10c, corruption_types_cifar10c"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load data\n",
    "\n",
    "Original CIFAR-10 data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "(X_train, y_train), (X_test, y_test) = tf.keras.datasets.cifar10.load_data()\n",
    "X_train = X_train.astype('float32') / 255\n",
    "X_test = X_test.astype('float32') / 255\n",
    "y_train = y_train.astype('int64').reshape(-1,)\n",
    "y_test = y_test.astype('int64').reshape(-1,)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For CIFAR-10-C, we can select from the following corruption types at 5 severity levels:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['brightness', 'contrast', 'defocus_blur', 'elastic_transform', 'fog', 'frost', 'gaussian_blur', 'gaussian_noise', 'glass_blur', 'impulse_noise', 'jpeg_compression', 'motion_blur', 'pixelate', 'saturate', 'shot_noise', 'snow', 'spatter', 'speckle_noise', 'zoom_blur']\n"
     ]
    }
   ],
   "source": [
    "corruptions = corruption_types_cifar10c()\n",
    "print(corruptions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's pick a subset of the corruptions at corruption level 5. Each corruption type consists of perturbations on all of the original test set images."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "corruption = ['gaussian_noise', 'motion_blur', 'brightness', 'pixelate']\n",
    "X_corr, y_corr = fetch_cifar10c(corruption=corruption, severity=5, return_X_y=True)\n",
    "X_corr = X_corr.astype('float32') / 255"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We split the original test set in a reference dataset and a dataset which should not be flagged as drift. We also split the corrupted data by corruption type:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(5000, 32, 32, 3) (5000, 32, 32, 3)\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(0)\n",
    "n_test = X_test.shape[0]\n",
    "idx = np.random.choice(n_test, size=n_test // 2, replace=False)\n",
    "idx_h0 = np.delete(np.arange(n_test), idx, axis=0)\n",
    "X_ref,y_ref = X_test[idx], y_test[idx]\n",
    "X_h0, y_h0 = X_test[idx_h0], y_test[idx_h0]\n",
    "print(X_ref.shape, X_h0.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_corr = len(corruption)\n",
    "X_c = [X_corr[i * n_test:(i + 1) * n_test] for i in range(n_corr)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can visualise the same instance for each corruption type:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAGbCAYAAAAr/4yjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAfAUlEQVR4nO3de5BeBZnn8d95731Pd7pzIZAb5IIhMJARWUvArO4oZUbLRdxCTBmKWqHc8rKOVFlYVqF/YBVaBVVWWYtu5LIIAg4wsogsuzXDxRkvBYOdSYAQkiYhJJ2k0+nLe7+c/YPxGVgYeJ5sMgn6/fxn58nTp897+b2nw/mZpGmaCgAASZkTfQAAgJMHoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKGAP1rXX3+9kiQ5qr972223KUkSjY2NHduDep2xsTElSaLbbrvtuH0PIIpQwElp69at+uxnP6tFixapWCzqlFNO0RVXXKGtW7ee6EMD/qgldB/hZHP//ffr8ssv19DQkK666iotW7ZMY2Nj2rx5syYmJvTTn/5Un/zkJ99xT6vVUqvVUqlUCh9Du91Ws9lUsVg86quNdzI2NqZly5bp1ltv1aZNm47L9wCicif6AIDXe+mll7Rx40YtX75cTzzxhEZGRuzPvvzlL+vCCy/Uxo0bNTo6quXLl7/ljnK5rJ6eHuVyOeVyR/cUz2azymazR/V3gXczfn2Ek8p3v/tdVSoV/fCHP3xDIEjS8PCwbrnlFpXLZd14442S/uXfDbZt26bPfOYzGhwc1Ac+8IE3/NnrVatVfelLX9Lw8LD6+vr08Y9/XHv37lWSJLr++utt7q3+TWHp0qXasGGDnnrqKZ1//vkqlUpavny57rjjjjd8j8OHD+trX/ua1q5dq97eXvX39+uSSy7R73//+2N4poDjgysFnFQeeughLV26VBdeeOFb/vlFF12kpUuX6uGHH37D1y+77DKtWLFCN9xwg97uN6KbNm3Svffeq40bN+qCCy7Q448/ro997GPu49uxY4c+9alP6aqrrtLnPvc5/fjHP9amTZu0bt06rVmzRpK0c+dOPfjgg7rsssu0bNkyjY+P65ZbbtHFF1+sbdu26ZRTTnF/P+DfGqGAk8bU1JReffVVfeITn3jbubPPPls///nPNTMzY18755xzdNddd73t33vmmWd077336itf+YpuuukmSdIXvvAFXXnlle5P8S+88IKeeOIJC61Pf/rTOu2003Trrbfqe9/7niRp7dq12r59uzKZf7kQ37hxo1avXq3Nmzfrm9/8put7AScCvz7CSeMPb/J9fX1vO/eHP5+enravXXPNNe+4/5e//KWk14Lg9b74xS+6j/E973nPG65iRkZGtGrVKu3cudO+ViwWLRDa7bYmJibU29urVatW6ZlnnnF/L+BEIBRw0vjDm/3rrwDeyluFx7Jly95x/8svv6xMJvOm2TPOOMN9jIsXL37T1wYHBzU5OWn/u9Pp6KabbtKKFStULBY1PDyskZERjY6Oampqyv29gBOBUMBJY2BgQAsXLtTo6Ojbzo2OjmrRokXq7++3r3V1dR3vw5Okf/W/SHr9v2PccMMN+upXv6qLLrpId955px599FE99thjWrNmjTqdzr/JcQJHi39TwEllw4YN+tGPfqSnnnrK/iui13vyySc1Njamq6++Orx7yZIl6nQ62rVrl1asWGFf37Fjx//XMf+/fvazn2n9+vXavHnzG75+5MgRDQ8PH9PvBRxrXCngpHLttdeqq6tLV199tSYmJt7wZ4cPH9Y111yj7u5uXXvtteHdH/nIRyRJP/jBD97w9e9///tHf8BvIZvNvum/gLrvvvu0d+/eY/p9gOOBKwWcVFasWKHbb79dV1xxhdauXfumO5oPHTqku+++W6effnp497p163TppZfq5ptv1sTEhP0nqdu3b5ekY3bn8oYNG/Ttb39bV155pd7//vdry5Yt+slPfvKv3mwHnEwIBZx0LrvsMq1evVrf+c53LAjmzp2r9evX67rrrtNZZ5111LvvuOMOLViwQHfffbceeOABffjDH9Y999yjVatWHVUdxlu57rrrVC6Xddddd+mee+7Reeedp4cfflhf//rXj8l+4Hii+wh/8p599lmde+65uvPOO3XFFVec6MMBTij+TQF/UqrV6pu+dvPNNyuTyeiiiy46AUcEnFz49RH+pNx44416+umntX79euVyOT3yyCN65JFH9PnPf16nnXbaiT484ITj10f4k/LYY4/pW9/6lrZt26bZ2VktXrxYGzdu1De+8Y2jblQF/pgQCgAAw78pAAAMoQAAMO5fot6+sv+dh15nYa+/42VBKfb/cFVKWu7ZvlLst2Nzev03MGUzzdDudtJ2z2bysRupmv5TIkmaqfofn2o9dizt1P9ZI5uJfS5pyv94TpZjj0+5Efs508CxtPtjN9u11q1zz04//neh3Qdy/p9zvFEI7R4qH3TP7prMh3a3emPvQertdY+OVyqh1QN1/3yxXA7trmT9r81sJ/b+dvuhxjvOcKUAADCEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAADj7j5avzi2uD/Q3ZPNx4p7Zqvv3N/xB5k01quUtvzH3ej4O0okqdbwz2cysW7/eit2LNN1/2y5GetXaQXOS/DHVDtQTzRTjXUZBauS1Er8+yvlidDunb/4P+7ZgXQ2tDsNPJ5J9Jxk/eekt3c4tHtHb19ofsuRcffsQDvYkRZ4ahVib0FqJYHuo0DPmBdXCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAACMu2RgUd5fLSFJ7VbbPVtT7BbzSsO/ux1r0FAj8GO2m7FqiWqgWiKNnRI12rF8nw0cetl/uiVJkdOSzcWOu534T8xsM7a7Fqg4kaR64FhazdgTMdPxP1mmi7EHqLfj710oBH5GSTqY+Hfv7S+Gdm+bLofmd01W3LPLgz9nruj/OUtp7H2i0wkcS/B9woMrBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGHf30cSRWPdRrePPm3YxH9pdTUv+4VysX2V6eto924526wS6j5pprIenpdh8NfUfeyXSxSKpFTj2fDN23M3E3/NTD3T8SFKtHTznib/TJm3H+m+6Ak+tmWA31VTLf14ySSG0u1rwv95eafi7iSSpfaQWmp/Xcb+9aTAXO4l9gadKPvbQqxjoMWsnx/5zPVcKAABDKAAADKEAADCEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAIz7PvBXk57Q4tms/3b33nx3aHej7r8lvVyJ3b5emfXnZNqK1SLUApUOtSRYLZHEjqURqKKoB2/TTwOVG4VgPUcr459vBGal2DmRpMhD1AxUYkhSLvUvz7X8dQ6SVB9e5J4tzvXPStLUvn3u2XRyPLR7QWhamsn4z/mS7lgdTj4T6Kzp6grtzsz6j7vdaYV2u77/Md8IAHjXIhQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGHdpyr5Mb2hxOS24Z9OJWmh3bbrinq0GepIkKRuZTWOZWmv7+2zqwa6cdjDe08Tfl9MJ/pxJoEOopVjHUztSTxTsg4rO5yLnMPLEkpRt+ztterKx12Zp7Z+5Z19KSqHdB+tN9+xg6p+VpJnpQ6H5ub159+zi/li3W2/gNZFmYt1UtXrZPZs06T4CABxHhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAAAMoQAAMIQCAMC4779+bt9saHGzEagvaMaqDjptf3VFJhOrLugk/pwsBOocJKme+qsrMortzgQqNCQpk/XvzwbqHCQpk/qPJZvEjjvw8CibC37mCTyvJElpoOYieCj5QM1FOjQY2r2r5X/sf7NzV2j39OEJ9+yquXNDu/vSRmh+WaBapCeJ1UVka4HnSqMe2p2m/hqf6Puba+cx3wgAeNciFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGEIBAGAIBQCAIRQAAIZQAAAYd3nL+GQ1tLgof/FINlZ/oyTQ91HMBApQJLXlP5hOMFPTQM9PpONHkrKBviFJygW+QaAm6Z/n/X+hlMmHdrfk749Kc7HHvp6PdR9l8/5jzxRiP2fSLrlnD/X3hnY/t2+/e3bni8+Hdufq/t6eUns4tHtF1v/YS1JP1X8sjST22Lfq/q6kfLCXLBt4jncC77NeXCkAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAAAMoQAAMO6ai7Z/VJKUJP7br3OB+gdJSgJVFNEKjch4PthFkQa2Z4K786Ejl3KZQM1FoFZEklrd/udKa+5gaHdX018BUCwVQrtn5a8ukKRc6p9vtmOPZzXQXjDTiu0+cHDCPZukzdDuvrz/ebiw5q+hkKR5aS0030791RWdYJdLPfCelQZrYjId/+sn2KDh+/7HfiUA4N2KUAAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBg3CUb2Uw+tDgT6OTIBXt7JH+nSUax4pGkE+gnCkZqpG8oWpiSSWLz+Zy/F6hnsD+0u9bnf660+7tDu9OJqn933d+TJEntdqz7qNz29wJ1cr2h3fVSj3v2SDNQlCSpt2uOe3bpkq7Q7q7WtHs21/a/jiVpqhF7jncq/m6lXCfY7xXoyWoHeuAkKUn97xORPjUvrhQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGEIBAGDc3UfFWDWIcqm/kyMb3N0J9AIFV4f+Qhr4GV+bD8xGyqMktRXrV2nliu7Z2UBPkiQdmKm5Z0s5f8ePJFXy/q6k0mBfaHf/4oWh+SXLlrhnF572ntDu7NBc92zlqV+FdtcP+R+f8T17Qrv3bnvGPbt//pzQ7ul8rIMrN37IPTtnZja0u536e5vSYC9ZJtDD1E7oPgIAHEeEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwLhrLkqJ/7bu1xb7b9XOJu7DkCRFbuxOs8cv9zqd4HzglvTozetpLvZzTjf9j2erGTuanhVr3bOr//1fhHbPXXSqezbTG6vQKA7EahRakdm2v55Dkiaa/iqK5edfENp94eIz3LNbf/2b0O7/9rtfu2f/fuzl0O6+voHQ/MXLznTPprt3hXa3J/b6Z4NVFJlAH047WLXj+v7HfCMA4F2LUAAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBg3KVDxWDHRprxdx8lgVlJSlJ/lmWCuRdpYWplg+ck8R9LJ40VK1WyhdB8z/KV7tnhs88O7S4uXe6ePZCL9dls2f6Kf/f4gdDu6uRUaH5m9oh79vBkJbT7SMU//+cX/Hlo9/v/6oPu2d4LY6+fpy/w9zDd//gvQ7sPTe8Lzc/rG3LPnh/og5KkyvSsezbT9M9KUi7QqtWi+wgAcDwRCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAONudcgpG1rciMRNbLWStr8WI5vGKjSygZw8kvhvR5ekvPzH0kxKod0Dq88KzTeXnO6e/e3BWP3DkbFfu2c7ha7Q7q07d7pnd+/cEdrdnbZD8yODfe7ZfROHQ7vrib+25MKLLw7tLpcb7tmunuHQ7ov+8lL37D9s2xbaPbbnpdD81lf2uGcLXf2h3UnRX8/SV6+Fdg8G3leouQAAHFeEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAADj7j5Kcu7R1+bTjns2m4mVH6VtfzdIJ4l1g3QSf05WG8HdGf85zC9fHtp9uM/fxSJJW7f8k3v2yORMaPfQ8Dz3bGswdtztjr+3J1uIfeapzMR+TnUNukfzAz2h1avX/Jl79n0finUf1Zr+jqfcbOw5fvZ5/849+8EPXRLafe/d/yM0nzb87xOjO54P7e7L5d2zI1n/rCS1O/7j7lIztNuDKwUAgCEUAACGUAAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAxt270Mkksc3+O7WVSWO7I4fSSPx1G5LU6O13zw7NXxXaXav5j+XIyILQ7qd37Q7NFwK1JUND/nMiScNz/fOvtGO36Tda/vne/thxZ3pKofnhxUvds+vXvTe0+0Mf3eCeHVm0JLS7Ufc/D3OlQmh3rV5xzxaC1R9r15wdmt+/Y5t7dqI6G9pdHhxyz5511rrQ7pGq/xxObvltaLcHVwoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADDuApw02E+kjn8+TWOrpax7slkqhjYfnjPonh0644zQ7krL/4PumJgJ7Z5/5trQ/J6Xt7tn27nYY58mDfdspRHrPlpz1lnu2Y9+9KOh3SuWLw3NL1p0qnt2aF6sy6oT+Lx26PBUaLfy/tdPq1ENrf7Jbbe6Z5964P7Q7rXzlobmay3/OZxsxjrSzjzT/zz8wIcvCe3OjY+7Z3+1dTS024MrBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGHf3UaJY/00m8edNM1h+1My6D1uTxd7Q7i1VfxdPeeuLod1dc+a4Z/vnxrpypsuV0PzL+/a7Z1P/6ZYklSaPuGfLk+XQ7r/62qXu2f90+eWh3Y1mrIcpDXRZVWb9fVCSVK/X3bO5YHdYLvH3/Pzirx8I7f71Xfe5Z7sOHQ7trs7GnogL5y/xzy46N7T7fReud8/Om7cwtLvQ0++eLQ6MhHZ7cKUAADCEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwPjvG2/H7qVP/XfSqxm47V6SOv1z3LML33tBaPeW8Qn37Mz+mdDuxlSg6qAwG9q988XnY8cyU3PPptl8aPfcgUH3bH6wO7R7YGDYPbtv/1Ro9+GZ2Hy16j+Hwae4Bgf8VQe9PaXY8kCtzIIFi0Kr1645xz1bmZwM7Z63bEVofnjlavds/8hQaHc28HF6Zjb2PjHY7X88O4P+15oXVwoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCB7qNWbHOr7R6tzJ0XWv2+T1/hni2te19o99/ed797dnbnodDuTst/DvNdhdDu2akjofnm7LR7ttjdF9rdXepxz86dH+vWyRb9u8cnjoR2z1bLofl2oM9osH8gtLseeK5Mj+8P7e7t8fdNnfvB9aHdhcDP+cq+vaHd+TnBc5gm7tlMsxna3alV/LPt2O49+19xz45XY71KHlwpAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADDumot2moYW15v+motTL/4Pod3v3XSNe/Z3u/23jEtS/8hC92y+Z0dod5r6b3dvNuqh3ZWZ4O3ukWOpz4ZWv7jrJffsaaefGdqdKfrrP2qtRmh3I1h10BWo8yjPxM7h//rFz92zo1ueCe0emT/inv3IX3wstPv0VWe5Z3PzTwntnjkyGZqv1P1VFPVAbYUkNQJPrcp07LH/1RN/5559ZV/s/c2DKwUAgCEUAACGUAAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABh391GllYQWd7oH3bNdS1aGdj/6G3/Xy/6pWF/KnMEh92yxVAztTtr+/qj9e3eHdtfq5dB8oeg/9kKpN7S7u3+OezZfiJ3DTDbrnm1ECmoktVqt0Hzefyh66G8eDO2+88f/3T2bJrHjTnL+z4LbRreEdv/n//Jf3bMrAz1JkpQocMIlHZ447J6tlqdDu5vlKffsk//70dDu0d/8vXt2KBfrpPPgSgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGEIBAGAIBQCAcddcNIIVAKWRAffsU//4bGj3Q5vvcs+efd45od1nnOOfLwaqIiSpVa27ZyvB2+5zuVgFQKbQ7Z4967zzQ7uXnLHaPdvV5T8OScoGai7CtRX5fGj+4IFX3bOP/s8HQ7tLef/ntaG580O7q42qe3bnjhdCu//mZz91z37iP14e2j0z4z9uSZo4csg/3K6Fdv/D3z7mnh39rb+2QpKKqb+epasnVkHjwZUCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAACMu/uorViPTK3TdM/ufmUstDuXSd2zMzMzod2FQsE9O2fOnNDuF1992T3bbPn7TySp2B3rQOkeHHHP9s2ZG9pdLlfcs0NDsd3z5s0LzUfkAr1KkrR96z+6Z6emDod2z+nrc89OTsZ2t9OOe7a/tye0e+uzz7hnV648M7R7wanLQ/OR1/LOF2IdT9uf2+qeLWb851uSRvr63bM9pVJotwdXCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAACMu+ZiNpgf9dmye7Yz4q+tkKRli09zz7aTJLQ7lf9Yurq6Qrvb7bZ7NlsohnYPDPlrKyRpcMES92yaxs5htex/7E899dTQ7kzG/zysVPx1G5KUBJ8r4+Pj7tlcPh/a3ROouejujVWczAYen+kjk6HdMzP+yo0dz/9TaPfCxUtD80nify3vGRsL7W5V/c+tOSV/3YYklbKB99pOrELDgysFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGEIBAGAIBQCAIRQAAIZQAAAYd/fRgXYsPxp1fydHpe7vBJKktOTvNOkEu0FqtZp/d6yySfVm0z2b6+oJ7R4YWRCaP3XJcvfs8ODc0O4k0h9VKoV279v3qns2DT5AxVKsb6odeG4luVj3UTbn78vpH5gT2t3qHPTPVhuh3ZWZaffsy7t2hHafse/l0Pxs2f9a3rtnT2h3o+7f3Uxjz8NKq+4fLsT61zy4UgAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgHF3H03ngvmR9c9na7F+lVavvyspTWKdM5VAX0rvUKwT6JRlK92z/cNDod0rVp8Zml+1co179tQF80O7I0+VYnesb6hY8D+eaScJ7VYam+/p6nPPZoLPw3bg89rCRYtCu0fm+3uynhsdDe2u1Gfds/vH/T1WkrR9a+xYypWqe/bggb2h3Y1Wy38c0c/e+ax/tnDsP9dzpQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAuGsuVPKPSpKS1D2ar/hvjZek/pK/MmAmcMe4JDWmJ92zhw9PxJan/nqO6mw5tHr7c8+F5vfv3uOe7e3qCu3O5/yPT76rENqdkb+KotPyn+/obkmaOnTAfyzNZmh3Ie9/vb24fXtodzbjf1EcODge2l1v1t2zMzNTod2/+9WTsWNp+Ctr6tVKaHcuUONTC9anpKl/dy7wWHpxpQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAOMuWMnkYh01pby/L6esTmj3/t0vuWerxWJo96t7XvAfx4FYL0x5asY9mwa6VSTJ3zT1mkj/TfiTQ+LfnWRjnVqZxH80SRo8K8H5jFru2WbD3wkkScsXn+KeTZLYOTx0yN/ZtWjhgtDu5573vyY6rUZo99RkrGssDbwqMoGuNklKI/PZWPdRJ+t/P0yT2G4PrhQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGPf98fnu3tDiQs5fL9HJxG7Tr9f9lQH7p/3VEpJUbvp3F4IVGiMLF/qPo1oN7W512qH5XKDmIi5QARA87sh8p+WvoTia+Vbqr2noBH/OrVtH3bOrVq4J7V44319dsXv3jtDuWq3ink2C5SyBhpN/3h+ogAi2RSSB6ookHzvwQnfJvzt37D/Xc6UAADCEAgDAEAoAAEMoAAAMoQAAMIQCAMAQCgAAQygAAAyhAAAwhAIAwBAKAADjLh2qB/uJ0kCZSCO4uzQ04J5d2N0f2p0UC+7Z3v6+0O605e+/Gdu1K7S7Wot1JZW6utyz2WysJymbdNyzSTN23O1m0z0bOd+S1Gr4d0tSu1NzzzYb/p4kSSrX/B1c217YHtqdCbzepqcOhnangQ6hfOC1JklpJ1ZQlCTHr/so8prIFfKh3YWiv/uok8b6ozy4UgAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYAgFAIAhFAAAhlAAABhCAQBg3Pe7t4rF0OJU/lu7c4Mjod3zT13inu2etyi0uxmIyXK1Etp95NCEe7bQ46/ykKTeoXmh+UjFQJqJdQDks/6TmOu0QrvTjr9CI237ZyWpWfPXVkhSo1Z2z9YCs5IU+DFVKPgrSyRJgfPSTmPVH/XGjHs2E/xImsnE6iIiDRCdTqwSJZ/1V4WUCrE6j2zir9DopLH6FA+uFAAAhlAAABhCAQBgCAUAgCEUAACGUAAAGEIBAGAIBQCAIRQAAIZQAAAYQgEAYJI0jTSEAAD+mHGlAAAwhAIAwBAKAABDKAAADKEAADCEAgDAEAoAAEMoAAAMoQAAMP8XXQidP2kqUAUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAGbCAYAAAAr/4yjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAmSklEQVR4nO3deXTV9bnv8Wdn2hnJnJAQSEiAhBDmQRFkCgiKOLYOx3qVllar9Zx22em0p4K2arVVr3f1qFWU9mi1onUAFRURREHm0TBkZApkIgOZh53f/eMev9cIyvP0uGpr36+1ulaJHz/ujB+2JI8+z/M8AQBARIK+7AcAAPj7wSgAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcBgFAIDDKOArZ8aMGTJjxowv+2F84ZYsWSI+n+/Lfhj4imMUAACOj9tH+Krp6uoSEZGwsLAv+ZF8sXp6eqSnp0fCw8O/7IeCrzBGAQDg8K+PcEbr1q2TCRMmSHh4uOTk5Mjvf//70/6d9rJly2TWrFmSkpIifr9f8vPz5dFHHz2ty+fzyZIlS057eVZWltx4443u193d3XLnnXfK0KFDJTw8XBITE2Xq1KmyevVql6mqqpKFCxdKRkaG+P1+SUtLk0svvVQOHTrkMp/+M4Wuri654447ZPz48RIbGytRUVFy/vnny9q1a/s8nkOHDonP55Pf/va38vjjj0tOTo74/X6ZOHGibN261fz28/l8snz5crn77rslIyNDwsPDpbCwUEpLS0/Lv/DCCzJ+/HiJiIiQpKQk+cY3viGVlZV9Mmf6M4XVq1fL1KlTJS4uTqKjoyU3N1d+9rOf9cl0dnbK4sWLZciQIeL3+2XgwIHy4x//WDo7O02vE/45hHzZDwB/f3bu3Cnz5s2TtLQ0ufPOOyUQCMhdd90lycnJfXKPPvqojBgxQi655BIJCQmRlStXyi233CK9vb1y6623mv+5S5YskXvvvVcWLVokkyZNklOnTsm2bdtkx44dMmfOHBERufLKK6WoqEhuu+02ycrKkpqaGlm9erUcOXJEsrKyzth76tQpWbp0qVx77bXy7W9/W5qbm+XJJ5+UuXPnypYtW2TMmDF98s8++6w0NzfLTTfdJD6fT+6//3654oorpLy8XEJDQ02v069//WsJCgqSH/7wh9LU1CT333+/XHfddbJ582aX+cMf/iALFy6UiRMnyr333ivV1dXy8MMPy4YNG2Tnzp0SFxd3xu6ioiK5+OKLZdSoUXLXXXeJ3++X0tJS2bBhg8v09vbKJZdcIh988IF85zvfkeHDh8vevXvloYcekuLiYnnllVdMrw/+CXjApyxYsMCLjIz0Kisr3ctKSkq8kJAQ75MfMm1tbaf9vXPnzvWys7P7vExEvMWLF5+WzczM9G644Qb369GjR3vz58//zMfV0NDgiYj3m9/85nMf//Tp073p06e7X/f09HidnZ2ndaWmpnrf/OY33csqKio8EfESExO9+vp69/JXX33VExFv5cqVn/vP/aS1a9d6IuINHz68zz/74Ycf9kTE27t3r+d5ntfV1eWlpKR4BQUFXnt7u8u99tprnoh4d9xxh3vZ4sWL+7z9H3roIU9EvNra2s98HE8//bQXFBTkvf/++31e/thjj3ki4m3YsEH9OuGfA//6CH0EAgF555135LLLLpP09HT38iFDhsiFF17YJxsREeH+f1NTk9TV1cn06dOlvLxcmpqazP/suLg4KSoqkpKSkjP+9YiICAkLC5N169ZJQ0ODujc4ONj9oXNvb6/U19dLT0+PTJgwQXbs2HFa/uqrr5b4+Hj36/PPP19ERMrLyy2vjoiILFy4sM8feH+6a9u2bVJTUyO33HJLnz9Anj9/vuTl5cnrr7/+md0fP4N49dVXpbe394yZF154QYYPHy55eXlSV1fn/jdr1iwRkdP+FRrAKKCPmpoaaW9vlyFDhpz21z79sg0bNsjs2bMlKipK4uLiJDk52f377L9mFO666y5pbGyUYcOGyciRI+VHP/qR7Nmzx/11v98v9913n6xatUpSU1Nl2rRpcv/990tVVdVZu//4xz/KqFGj3J9VJCcny+uvv37Gxzlo0KA+v/54ICxDpO06fPiwiIjk5uae9vfm5eW5v34mV199tUyZMkUWLVokqampcs0118jy5cv7DERJSYkUFRVJcnJyn/8NGzZMRP7f+xv4JEYBf5WysjIpLCyUuro6efDBB+X111+X1atXyw9+8AMRkc/8nesnBQKBPr+eNm2alJWVyVNPPSUFBQWydOlSGTdunCxdutRlvv/970txcbHce++9Eh4eLr/4xS9k+PDhsnPnzs/85zzzzDNy4403Sk5Ojjz55JPy5ptvyurVq2XWrFlnfJzBwcFn7PH+im/U+yK7Pi0iIkLWr18v77zzjlx//fWyZ88eufrqq2XOnDnubdvb2ysjR46U1atXn/F/t9xyy//4ceCrhT9oRh8pKSkSHh5+xu+Q+eTLVq5cKZ2dnbJixYo+vxs+07+OiI+Pl8bGxj4v6+rqkhMnTpyWTUhIkIULF8rChQulpaVFpk2bJkuWLJFFixa5TE5Ojtx+++1y++23S0lJiYwZM0YeeOABeeaZZ874Or344ouSnZ0tL730Up/v3lm8ePFnvyH+RjIzM0VE5ODBg+5f6Xzs4MGD7q9/lqCgICksLJTCwkJ58MEH5Z577pGf//znsnbtWpk9e7bk5OTI7t27pbCwkJ+GhgrPFNBHcHCwzJ49W1555RU5fvy4e3lpaamsWrWqT06k7+94m5qaZNmyZad15uTkyPr16/u87PHHHz/tmcLJkyf7/Do6OlqGDBnivnWyra1NOjo6TuuOiYn53G+vPNNj3bx5s3z44Yef+ff8rUyYMEFSUlLkscce6/M6rFq1Svbv3y/z58//zL+3vr7+tJd9/J1UH3ddddVVUllZKU888cRp2fb2dmltbf0fvgb4quGZAk6zZMkSefvtt2XKlCny3e9+VwKBgPzud7+TgoIC2bVrl4iIXHDBBRIWFiYLFiyQm266SVpaWuSJJ56QlJSU054BLFq0SG6++Wa58sorZc6cObJ792556623JCkpqU8uPz9fZsyYIePHj5eEhATZtm2bvPjii/K9731PRESKi4ulsLBQrrrqKsnPz5eQkBB5+eWXpbq6Wq655prPfH0uvvhieemll+Tyyy+X+fPnS0VFhTz22GOSn58vLS0tX+wbzyg0NFTuu+8+WbhwoUyfPl2uvfZa9y2pWVlZ7l/Hncldd90l69evl/nz50tmZqbU1NTII488IhkZGTJ16lQREbn++utl+fLlcvPNN8vatWtlypQpEggE5MCBA7J8+XJ56623ZMKECX+rVxf/CL7cb37C36s1a9Z4Y8eO9cLCwrycnBxv6dKl3u233+6Fh4e7zIoVK7xRo0Z54eHhXlZWlnffffd5Tz31lCciXkVFhcsFAgHvJz/5iZeUlORFRkZ6c+fO9UpLS0/7ltRf/epX3qRJk7y4uDgvIiLCy8vL8+6++26vq6vL8zzPq6ur82699VYvLy/Pi4qK8mJjY71zzjnHW758eZ/H/ulvSe3t7fXuueceLzMz0/P7/d7YsWO91157zbvhhhu8zMxMl/v4W1LP9C2v8hnfVvtZPv6W1BdeeKHPyz/+ZyxbtqzPy59//nlv7Nixnt/v9xISErzrrrvOO3bsWJ/Mp78ldc2aNd6ll17qpaene2FhYV56erp37bXXesXFxX3+vq6uLu++++7zRowY4fn9fi8+Pt4bP368d+edd3pNTU3q1wn/HDhzAbXLLrvsc79lFMA/Pv5MAWfU3t7e59clJSXyxhtvfCVPUgP4/3imgDNKS0uTG2+8UbKzs+Xw4cPy6KOPSmdnp+zcuVOGDh36ZT+8L01XV9cZ/4D3k2JjY/v8YB/wj4Q/aMYZzZs3T5577jmpqqoSv98vkydPlnvuueefehBERDZu3CgzZ8783MyyZcv6HPoD/pHwTAEwaGhokO3bt39uZsSIEZKWlvY3ekTAF4tRAAA4/EEzAMBR/5nCSyNsPyIfbHj23N5p+88mzj8v6eyh/+ZtOX720CecaNNnc0dnm7qlulsdrUs7aqo+dPqxz881IkafrQkbYOqOjqw8e+i/+Sps/32ChF7972NKTtj+IzIna2y/R9qcP0ydjc/Qv+9FRIKunqjO5jz1Z1P39q4r1dny8GhTd33tX9TZpOxYU7e3y/bT11Wd+h/KOxA4+1HFT5rQqP8CV5+4+uyhTxh9Uv/ndscKak3dj+09+1FHnikAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAR30ldfsC2+2j+voMdXZOYKCpuzlkszq7rdX2uGdeEFBn218zVYtvoP5WTllmsak7L8xwzEhEqjqb1Vl/ualamlqGqLO7cktN3aNfz1Rnt4w5bOpOtp35kbYV+mxnzGhT96sj9DdtUqJt971iW/X/AaC2ze1nD31CS0D/ed8QYrt5tvcS2/2o6jcb1dlJI/Wf9yIi40R/JK2pOM7U3d3RqM6WxSWaut8+XnfWDM8UAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwQrTB/jLBVDw+/oA6e/CA7RRFdUqvOjuzw/Zj4N0fnVRnI0JCTd2bG5vU2Zh1pmrZOVZ/tkJEpM5wuqK5Zqypuyp9nzp7Tpf+LIKIyDu+BnV2UlWkqXuN4QyJiEjmmF3qbOe4s58X+KQRRVHqbEWRqVr6n9KfrkgqsL1Ndrfrz7PkxM8ydW/Ztd6Ur2vrUWfX2C5uSO7ubHU2PrrR1F3drv6yLON8+pMyWjxTAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAI7P8zxPEwwMst0neneAPjuixtb9UY7+Lkz/+mBT96hjAXX2cFqqqTu6qkyd3Vhle9xNg4eb8mnHP1JnyyebqmVA8GB1dvv2WlP3eY0t6uzGtBRT94U1tt8jrQzo700FjW41dTe1Z6mzJ4/EmLpLQ+rV2YtiK03drZE56uyqyY2m7px9+ve9iMhftsWrs1MTq0zd2YZP/WhbtYQ2JKizYaH696WIyP2dZ/9yzzMFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAc9ZmLojzbKYqc8BnqbFi/rabuoHj9GYVHWmyPe1iD/sfGYxJDTN3v7GtUZ3MjO03draURpnzHpU3qbGBdtKnbi2pTZ8NCw0zd1Ycz1NljWcWm7uT+oaZ8ly9SnR0cMNx9EZEw/z59uM1v6t6SF6vOTqw4x9S9plf/NszyNpi6P+ypNuUTtus/VlrOs/3+ODf6mDp7oGysqTsneq86O6Bfl6n7P9Zz5gIAYMAoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgqI/3bIxNMxXv7VmnzuYF+pu6V+6rUGczSltN3SVD9dnInfobTCIisZn6e0NN41JM3SeP1pnyISvC9eEL403dPW8E1Nl9Ynvff03K1Nl14VGm7hNd/Uz5c1NPqLNlq0zVEpaXpc7mDNC/vUVEJqYvUGff89vuQe1r1t8E6unW3yYSEYl4xXYnq3OM/mvW9E07Td3tkfp7UxMG2O6Y7SzW3xoL8myfP6rOL7wRAPAPi1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAcRgEA4DAKAACHUQAAOD7P8zxNcNkE249TH2iqVmeHHzVVS1mnfsv6Des1dffUjVFnfd4uU3e9/ifjJUv/k+4iIhIcYcvvPq7P5p20nTqoT+9WZ0cFRZq6tx2LU2erJ580dV/UbjtHcKxUnz0elW7qHtxUpc6WjbOdW6nqn63OFm/tsnWH6U9RjBtte/98Ld72SbFtk/4Ly+ge2wkNOXhYHw0fYqoeHKM/41Pvt53x+enRs3+555kCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcEK0wZ2V+ltGIiIpsVnqbFvXIFP3VNHfTCmNLjJ1r6/dpc4O7jBVy7We/vVc2mK7w5MeZTt+lBiUoQ/HHDN1Zx4fps5WDd1l6g69Sn+0aeTbto+rw8OPmPKJC/T3jKI+MhybEpHOy/TZ+JIoU/dvt/nU2dRGw4EnERnYPECdDa80HAMTkf0dNaZ80AL9+7OoKM/UnR+r/+RvH95o6m4v0d8z2h9su6ukwTMFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAc9ZmLI1W24smGfHV6man73eP16uzo8nBTt9eiPxkwtHegqftI9i51Nroo29Sd3qE//SEiUpB8SJ090WM7R1AT0a7OZtWkmLpD2vWnK0JG206zDOwabMoHDmSps88HTNUypzxRnS3yYmzlw99WR3sibd31m1vU2aT0JlN37ME0Uz7xPf0b3X/Kdg7nuWR99oqaVFP3jmh9tjDTdoZEg2cKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwFHfPlooEabi9+IOqbPBeReZumcc36bOvhVku32U3lOrzsaP2WXqfrd2iDqbKbabJk2hwab82jj9TZvJo0aauqMGVqqzO9pjTd39K/T3jFoCNabuzcVTTfm88zvU2fSd6k81ERF5N6ifOjtwbaipO/eqW9XZ5s3bTd2pKfqP26P1trtXKfGRpvxWw6fEeads94nGnNJ/HK5LLzF1z6vTf6yUrL/Y1K3BMwUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgqI9sPBKhv5UjIjImWn9zqG7vflP37lz9TZshEmbqrjTcYyk5qr+TJCISl9SkzlbH2t7ek/XVIiKytXu0OnuyfKepu9Hwvk850W7q7gj3q7P1U2aauhN8A0x5KdC/j77/s5tN1R/UTFZnSxNfNHXnxHykzobljzN1d6YeVmcDTXmm7j9V2z6Xp0Trf89b3ZVj6m6q3KzODtndYOpONJyZW55fZurW4JkCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAACO+szF5T360xIiIvHR+mzpsVxTd6ChQp09Fhdr6m6TY+rsBH+8qXvbQX02L6vZ1F0fbIpLS0OROlsyebypOyqvWp1tW3ChqfvSGRers15lj6nbf8UBU75fdbI6W96sf9wiIqE9leps0uxBpu7JGf+mzva89l+m7ndfXanObvOtNXVfcMT2cSiZ/dVRf802U3VHrP4cTua/2D4O33hU/7mfta/F1K3BMwUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDg+DzP8zTBlxMHmIrfbTuuzl7YYaqWTYn67KmTkabuhpA4dTY6oH8dRUTCB4eqs6FB3abuwfW2O0zdw69UZwfNG2LqPhk0UJ0t63/K1J3XXKXOvtlrux9VcqTVlB+zu586+3ZDmak7UfSP/dJ/nWnqvnzBbepsaMBv6v769x5QZ/e++O+m7sjY4aZ82tdr1dlzS2z313oaN6izJ5ptb8OkaH0+v8z2+fPLlrN/ueeZAgDAYRQAAA6jAABwGAUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAATog22BJVbyoe5dOfOqiJ0J8LEBEZkFCkzg7uiTJ192uKUGdXRg0ydV9Uc0Sd3dw53dQdnW+7FVIxoledfb9yr6l7/LFV6uyRfueauleUbFFn65M+MHW3rhptyrdcfkKdjVlbaeoekJmnzqYd6TR1n0rSf04E3jpg6r5tSJY6u3jSjabuE1tsH4dlS/W/570odYype1eY/gzJuOZSU3dIj/593x36kalbg2cKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwFHfPmqKtd3WyTjqU2d94pm6g5v13YeiG03drQlN6mxOe5KpuyRSf9OkYPpRU3dR9SxT/tC7h/TdwY2m7tYo/V2Y+iO2m0BBXaH68NZuU3fjONsdmbaD+s+J+NxEU3fWlCvU2ZiptjtZcRs2qbMJ44ebutvGtauzo/22m0ARW94w5ZskoM4+3XPQ1N3m26PO1qeYqiXipP5G2oxG29dODZ4pAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgqM9c1HZGmYpb+jWqs8FDq0zdIcX6H+0eM2aUqXu3t0udbe2JM3VLif7n3T+onmCqfrvlWVN+aJj+sQyNsP2c/uR549TZ42W1pu73tleos/3PKzR19waKTfl/GZGrzjZnLjB1F/4yTp0dVmQ7c9Hdvk2dbY/VPw4RkY4PBquz2XHJpu5TlxSY8rWb9OdWhtfvMHVvKchRZ9tS0k3dEXVt+nC+/nNNi2cKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwFHfPjoV1G0qHpvlU2ePNNnuKsWEBNTZP9cfMXVfFae/Z3O07SJTd+i8o+rs5oP7TN1Nwy415XtPbVBna2tiTd0b/qtdnT2cW2bqTiocqc4u/l+XmbpDAnGmfPbUWepsqtg+xqViozq6c0SJqTpixyR1tlz/qSYiItv+/A119u3Xwk3dSRm2z4lZNR3qbFlOhKk7x1N/6ZTCb803dc+Yqf98e/OXvzJ1a/BMAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADjqAx5hYfqbQCIiA47r7988e3mpqXvWTn22PcZULUu79bd19kS+aOoeWnqFOls/eoWpu3/lFFP+jdBWdTa52/AGF5HYhEh1tnSd7W7Pip9/S52dMvl6U3fQrkOmfIN8pM7u6TzH1N26W/9xmB2pv2MlIhLXnqzOrvmlrfv+P9ers9cG9LeJRER8xweZ8vuy9I8lKm+cqXtu5nnq7HmxM03dp4qq1dmythZTtwbPFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcNRnLsr27jUVt+SMUWcLVqebusOj/Ors8Jq5pu63C7aqs/u6fKbujojj6uzgDwpN3WtiNpryskX/4/GjPds5gv4D56mzVRm2t2GJv78622K7niKSaTvp8PTJ89XZfiG1pu7roveos8eOjzB1h+Tr3595B2zv+9umjlFnW7oPmrrjY4eY8oPm6B97RHGaqTtjzER1tvek/muKiMiIaUfU2ZFHbSc0NHimAABwGAUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABz17aMF+Qmm4nf27Vdnh8tkU/ehH1+kzl4wbpypu/Y/T6iz6w6bqmVC0yvq7KrIGFN3WqftbVgXq89uirLdhYnur7+rdGX3aFP3sI8i9OEF9abuTcW2O1kzxujvTVUduMTUvbI5V52NSOln6m47pH8bFs6dberO9Z2nzr7bucvUXTSuzJTP2pGhzu6Z0mXq9tfpP/kPjAwzdcdt0b8/NzS8YurW4JkCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAACO+syFVNhOBlguQPi+mWrq/sVPb1Zndx/bYeo+PMBTZ6dG+kzdBysX6LvzPzJ1rzpw3JTvH1atzvr3B5u6+01Yp85+UDzL1P21+qPq7N6DttMssY0dpvyrR/WnK2amLzd1/5/ni9TZyNe3m7rj9iaps/semGPqnnJVoTrb//3ppu74rY2mvG9mnDo74y39572IyNHaVnV2UkObqXvz2tfV2SMVpmoVnikAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAR3376EiI31Q8X6LV2ZrBpmr51q4P1dkZh2yPe0jH1ersys4HTN39YvW3ko6u0N8mEhGRENsdpqpQQzix0/ZQ1gxVZyd/rdzUnZCuv61zsFF/P0hEJHOU/hSYiEhBW391tug/9LdyREQqNr2rzqYFuk3dDYP078+2f/8PU7fvNv0tq6HxAVO3xOlvNomIyHP6u2etweNM1R1x+ltjz61caureuL1WnU3tNn7xVOCZAgDAYRQAAA6jAABwGAUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAAjvrn+lMHn2sqbjl+UJ0tXTPI1P3Koz9VZ49d/G+m7iHj4tXZOS+fZ+p+I0XfLbmbTN1J1fofuxcRqa/K14fnzDd1/7AzU50NC200ddeN159GmNEywNT94gnbyY2gXR+os//ZudnU7WsOU2cPe1Wm7uCd+u7W6JOm7p51+s+30MnfNnUPLE825bcF9B8raeNtn2+b/vS0Onu8/oCpOylLnx1ZY/g8VuKZAgDAYRQAAA6jAABwGAUAgMMoAAAcRgEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHDUt48699hut6SHDlFna6pXmLoHJ6Sqs9Xyv03d88d8R51t/CDC1J3a+JY6W10fberurOox5XtF3/+vg203gVrq9Fl/fJ6pOzo8V53dmLjL1D3i0HhTfmPuy+psxSN7Td1jQ1rV2dJS2+/tJl/2vjq7vS3D1H3klW36cKqtO3PKZFM+pqhNnS17wlQt723W3zPyx9q6B0Tqs17BPlu5As8UAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABw1Gcuoqbpz1aIiKzdflKdrQgZauquLGhQZ1PXZZu6S6f0qrNdWd2mbv/JcHU2b/MhU3dicJIp7105QJ3dVB1j6s5M1v+c/pzU/abujM5gdbZso+39U3+iwpRvDzE89iMfmbpb2q9VZy9JsZ3QWDcooM7mPWw7n7JD9G/D1trXTN2lh6pN+aiYRers29G/MXWfVH/lFMltijJ1Z4alqLOHdneZujV4pgAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAc9QWPDTX6m0AiIqWtfnW2rrHD1N2vSX8b5LzkYlN3R3WGOusdKjV1h+zX3+050K/G1D3wvDRTPtMXr85elZ5r6g4K1Wd3J4wxdYd0bFRnw6KPmbob8g0PXERqPxqnD7cPN3XnZ2xRZ1MHzjd1Zx3dqs7uHWj7uGo/tUqd3bTC9vZu6We7vzav9R11tuuE7U6WNOjve7WFJJqqPzLcj/Iao03dGjxTAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAI769tF7B/aZigsMexPaPc/UnbfneXW2eOIsU/fAvbPV2WHDwk3dq4sb1dlLv2W7CzPl2LdN+f5zo9TZxpE9pu5pSePV2Snv6u+8iIjs7zhPnc3PbDB194RuM+Xza/W3r0ZnV5q6y7LGqrO1Kba7PbNj9Xd7oo7qb02JiHwYHlBnO3uPmLqrdti+Bu1P1z+Wyu5qU7ecq7+RFva+7f5a8MEEdbZXckzdGjxTAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDA8Xme52mCo3zDTMWjpF6drQxrMXUH3TpYnc37yyBTd9ttBers/tJeU3f1s+qrIuKL6G/qjiywnTpIG5ipziZ0vGHqLh8/Qp2d01Bl6q5M0P8+ZkZ5k6l7X67tY2X7m1vV2fIdtlMHA0YPUWdL/fqTCyIiec1H1dn1NWGm7qQRp9TZjhWHTN0yIMUU943zqbPNxZttj2WvIZseY6rOS2pWZ+sPxpu6qzvO/nWZZwoAAIdRAAA4jAIAwGEUAAAOowAAcBgFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAUd8+un6Q7RZP74hqdbZNBpq6D3dMVGcbKg6Yur2Uc/SP4+SHpu5saVdny+v0bz8REfHFmuLhTfqbNgPSVR8izuEY/d2rnoNtpm4R/c2mkIw6U3PYMdudn7bRieps5u5SU7cv/2J1NnWs7fd27X86os7mn5Nq6n6ntlydrWstMXVnR5niEqw/ISQVtZGm7h7JUWcTzikydYdt9quzQSH6rykiIpXdZ/9c5pkCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAABOiDZY3xFuKo6oH63OpsbHmbr3NZxQZw/lNpi6Q959Vh8OyzZ1J8WMU2dPnXrZ1N0sHaZ8R3KKOltz3PY2HBaUrM7uy80wdU85WazO+o71M3XvTuo05aVcf4rklD/aVN098AN19tyiNFP3oK/r3/c7Xt1s6o7PGqTO1iWYqqW+I8mUb6zVf2yNDjpm6g6P15/zOLC519SdlaL//Am1XedQ4ZkCAMBhFAAADqMAAHAYBQCAwygAABxGAQDgMAoAAIdRAAA4jAIAwGEUAAAOowAAcHye53ma4Mgo232ViWlT1Nk32/R3REREkgv1Bz/CKmtN3cWjstTZmLxmU/fcR8arsy/E7DR1N7eoz1iJiEjS1Ap1NnZ3lKlbevTvn7re46bq6NJ0dbYjot7U3VQ/0JQvyCtTZ0PqE03dnUeC1dn9GbYjQl5SmDob2BUwdYu8rU6eK/Gm5k2xthtcqc3617O6X6ypW6IMX1dSC0zVCWX6z7fOWNvXiZbDZ7/vxTMFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAc9W2EYdGhpuLImv3qbGr+GFN3flCmOjto6ilT99BKw+u53VQt7aO3qrMdvmxT97Qi1bUSZ8CWAepsUYb+5IKIyNHGwerspO5uU3fzzHXqbOW+GFN36iTb2/Do+4PU2WFh+s8HEZHWGeers3ldtnMrR+qK1dkesX3+tGWmqLMhUfpzKCIiEjXKFD+5dZs+bPswlNhavzqbmmA7E9OWf0idPbU3z9StwTMFAIDDKAAAHEYBAOAwCgAAh1EAADiMAgDAYRQAAA6jAABwGAUAgMMoAAAcRgEA4Pg8z7MdfAEAfGXxTAEA4DAKAACHUQAAOIwCAMBhFAAADqMAAHAYBQCAwygAABxGAQDg/F8zNzbKgpvAsgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAGbCAYAAAAr/4yjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAemElEQVR4nO3dbYydBbnu8etZLzNrpjOdmb7QUgotlJdjQ+RsMbhJYFf33qhRE0RN1EhCEwoxRIkKH9AoBTURY9AaE41EQxH9QGJi/CCRTY4BjsaziRI1yDaA0lIobaftTOdlva/1nA8e7lh5u6+dFpDz/32y03vuPut51lrXrJbnsijLshQAAJIqr/UBAABePwgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBbxh3HLLLSqK4rU+DBVFoU984hOvOLd7924VRaE9e/ac/IMCkggF/ENpNpu65ZZb9MADD7zWhwK8IREK+IfSbDZ16623vmgofP7zn1er1Xr1Dwp4AyEU8IZRq9XUaDRe68N4zSwvL7/Wh4A3AEIBJ8zzf6f/+OOP68orr9TU1JTWrl2rL3zhCyrLUvv27dPll1+ulStXav369br99tuP+/5Dhw7p6quv1rp169RoNHTBBRforrvuit/fs2eP1q5dK0m69dZbVRSFiqLQLbfcctyf/7f6/b6+9KUvacuWLRodHdXmzZv1uc99Tp1O57i5zZs3633ve59++ctf6qKLLlKj0dBZZ52lH/zgB//t8/GjH/1I5513nhqNhi688EI99NBDr/g9f/t4/v74tm/fHr9+/t8jHnzwQV133XU65ZRTtHHjxv/2sQLPIxRwwn34wx/WcDjUbbfdpre97W368pe/rF27dumyyy7Taaedpq9+9as6++yzdeONN8YbZavV0tvf/nbdfffd+tjHPqavfe1rmpqa0vbt2/XNb35TkrR27Vp95zvfkSRdccUVuvvuu3X33XfrAx/4wEsey44dO3TzzTfrLW95i77xjW9o27Zt+spXvqKPfOQjL5h98skn9aEPfUiXXXaZbr/9ds3MzGj79u364x//aJ+DBx98UJ/61Kd05ZVX6otf/KKOHDmid7/73Xr00UftXS/nuuuu02OPPaabb75ZN9100wndjf9PlcAJsnPnzlJSee2118bX+v1+uXHjxrIoivK2226Lr8/NzZVjY2PlVVddVZZlWe7atauUVP7whz+MmW63W1588cXlxMREubCwUJZlWc7OzpaSyp07d77kn/+83/3ud6WkcseOHcfN3XjjjaWk8he/+EV8bdOmTaWk8qGHHoqvHTp0qBwdHS1vuOEG6zxIKiWVv/nNb+Jre/fuLRuNRnnFFVfE1+68885SUvnUU08d970v9tg2bdoU5+pvv/eSSy4p+/2+dXzAy+GTAk64HTt2xP+uVqt661vfqrIsdfXVV8fXp6endd555+kvf/mLJOnee+/V+vXr9dGPfjRm6vW6rr/+ei0tLenBBx+0j+Pee++VJH3mM5857us33HCDJOlnP/vZcV/funWrLr300vj12rVrjztGx8UXX6wLL7wwfn3GGWfo8ssv13333afBYGDveynXXHONqtXqCdsHEAo44c4444zjfj01NaVGo6E1a9a84Otzc3OSpL179+qcc85RpXL8U/JNb3pT/L5r7969qlQqOvvss4/7+vr16zU9Pf2CnX9/3JI0MzMTx+g455xzXvC1c889V81mU7Ozs/a+l3LmmWeesF2ARCjgJHixn1xf6qfZ8lX4f4PN3tD2Wh7jK3mpTxdjY2Ov8pHgjY5QwOvCpk2b9MQTT2g4HB739T/96U/x+1L+Df757xkOh3riiSeO+/rBgwc1Pz8fO0+Gv/8zJenxxx/X+Ph4/BdUL2ZmZkbz8/PHfa3b7eq555470YcIvChCAa8L73nPe3TgwAHdc8898bV+v69vfetbmpiY0LZt2yRJ4+PjkvSCN86X2ilJu3btOu7rX//61yVJ733ve0/Akb+4X//613rkkUfi1/v27dNPf/pTvfOd73zZfwPYsmXLC/7T1TvuuOOE/jsE8HJqr/UBAJJ07bXX6rvf/a62b9+u3/72t9q8ebN+/OMf61e/+pV27dqlyclJSX/965KtW7fqnnvu0bnnnqtVq1bp/PPP1/nnn/+CnRdccIGuuuoq3XHHHZqfn9e2bdv08MMP66677tL73/9+veMd7zhpj+f888/Xu971Ll1//fUaHR3Vt7/9bUl/vb/i5ezYsUMf//jH9cEPflCXXXaZfv/73+u+++57wb/HACcLoYDXhbGxMT3wwAO66aabdNddd2lhYUHnnXee7rzzzuNu2pKk733ve/rkJz+pT3/60+p2u9q5c+eLhsLzs2eddZZ2796tn/zkJ1q/fr0++9nPaufOnSf18Wzbtk0XX3yxbr31Vj399NPaunWrdu/erTe/+c0v+33XXHONnnrqKX3/+9/Xz3/+c1166aW6//779W//9m8n9XiB5xXl6+Ff0QAArwv8mwIAIPDXR0DSgQMHXvb3x8bGNDU19SodDXBy8NdHQNIr/eewV111lXbv3v3qHAxwkvBJAUi6//77X/b3N2zY8CodCXDy8EkBABD4h2YAQEj/9dGDF05ai2dG8x9AVo54H1ZWGLsbI8NXHvobtZoxX/F2D4yH2XWGJbW71rg6vXxdRH+Yn5WkgfLzHfNxNvv53a2B9zNPr+49x3vTL11X8YJjMc63JC0O88+tuSNeYV+r188fR8u7Pt1B/nEuF97/S167Yv4MO8w/zrLt/d+4VvvGC26QPw5JGhb5c26MSpLuOvzKx80nBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAAhHT30aoxb/FKo5S7YXYf1ar5XphCZjmIMW7U00iSBkbPz8Do+JGkfs97nH2jjqVvnsK+UbzbG5i7jYMZmNenb/6I1DeOvdvpWbuPHZtPz3bMJ2LZNS6+eYFKjaZnq6N1a3fP7D5abLbTs+PmORw3ZiveS9lyMlbzSQEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBQBASJdRnLbCuw28KPO3xw8Lb/fAuPO+NKsO1M1/w2Do9T/0jHaBttkt0fJaFNQxzovZoCHjYTqnW5LUHuRv7G+XXglAx+mtkLTQOZqeXW51rN3dQf6C9s3nodNbUi3ytRWS1Knn+3D2V40uHEkHmk1rvm5Ui4yanSg14z2ralbtWJfTfI5n8EkBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAAAhXT5SHXStxS2jwKNbutlk9H2Y3SDdTr7/pu+U/Ejq9o3eHqPjR5LaboeQcV7Mh6mBcX3Miifrcba91Wr1vH6i5TL/mui7JVxOH5h5Druqp2eXjC4jSdpfze9+cnnZ2j1Y8q7oeuOcj1W8k9hwhs1LPzRem+alT+GTAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAICQrrk42qlai1tFfr6o5G+Nl6RBL3/feKfTs3Z3W/lbzHt9L1M7w/xur3BB6rh1Hsb98QPzXnrjYdq7ncfZMXf3Su8bhsb8oOJ1HTitGP38y1iS1BpbmZ6db0xau2cXWunZoVlbsXLgncOVRnXFyqr3+hkxDmVY9a7PMN+0Yz0Hs/ikAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAkC7lmC0mrMUD5fuMeq2+tbu9nO9X6bS8vpTBIN/Z1Df7hpwunq612Z8fWsfuPU6njaVfeN0tTpOVUSEjSRpaR24ynyvD6mh6tlix2trdn1mbnj263LR2a5ifX13z+tTGq95reaae37+i5l2fyjD/8/RgYPaSDfPvhxV55yS3EwCA/4dQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAAhHTNxb6md0t6r9lOz3ZbTnmB1O/lbwMvzHoBp7qib7Yi9Ir8LelDu1rCm1eR/3mgsI8l/zjLiln/YIxXzN2FWXPhnPNaJV/7IknV6XwVxdIpG63dB4f51/KxuSVrt8p8ucj4ihFr9erCew86pZ6/PmNmXUS/nd897Hnvb0XhFrScWHxSAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBQBASHcf7T3odaCUPaNLZOh1zhRO/43R8SNJ+VYlqXS7coxepYrzICUV5uOslEb3kVurZByL209Uc8bN3X3zcapinPMVk9bqpZmZ9OwzZr/XnqPH0rPNZtPaPWHUE622LqZ0ei39diVJWtHvpmeHfeeVL5Xd/G4NvS6jSpF/X7E7zzJ//gnfCAD4h0UoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAQvq+8YWWd6t2zbj9umZWNDgq5m3gTnWFWy1RGLvdeo6qme9V47RU3Z6LmlGh0ah7u40KgKJmdC5IqlbN62nsb45OWLtnjZfb03NHrd0HZ/PzE6VR5yBpjdFEsXHg1cRM9nrW/KCfnx8OjFoeSUOjusJsIXGe4uZwDp8UAACBUAAABEIBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQ8k0lpdkLY/R9FGY/kcPvPsrPu7udc+J0R0lStWLOG709tRVj1m6nz6gyYpTlSCqMjhqnx0qSzFOotvLncLH0HueR5Xzn0PzCkrV7pMyfw/V1r5tqQ6Wfnl0xyM9K0mDozfcGRoGUWSHkNCU57yl//Yb8vPscz+CTAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAQrqQpVp4+VE1OjnsZDLqPipO4ZCs2hG/+8iYrVS8s1IfGbHmK5NGn5EzK6lX5i9Q3egPkqSyzPfZDM1zOBwbt+YrK6bSs6Oj+VlJmi7z52XN+GFrd9ntpGdXt5et3cOFI+nZJbNXqTD714Zl/nFW+l6v0tA4lqHzpiKvhunENx/xSQEA8DcIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQEjXXNTNG6qd8oKKmU1l4RyLd4u5xV3tVG5U05dGkjQca1jzbWN/b+g90NHxfKXD+Oq13u7p/O6RVdPW7rHV6635ofE4n2v2rN3HZo+mZzetPsXaXW0307PDg/ut3cc6+VqM+V7b2j3a996DVo9NpmfLVv6cSFLZyVdoDMwuCqMlRuVJeH/jkwIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAEK6AMdr4pHV81M4nUCSyooz7+0ujF6lgVk7UqnmG6F6xqwktUov3yuTq9OzM+eca+2e2XxWerY2vcraPWzkO57Ket3avVx68wvtbnr22U6+y0iSysZYenbjxpXe7sVj6dl9ywvW7rnR8fRstze0dq+ueSVCI4NBenbSOG5J6nXyu8syP/vXb8ifF6cnKYtPCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAABCur2iWnj5MXSqKMyaC6vkwqyiKI1vGJiR6pyToXnb/cQWr4pi7Oz8/HDVKdbuZyoj6dnWXMvaffTYwfTssfl5a3er2bbm271+enZQ9Z6IZ205Oz27acsma/f8wQPp2WcO7Ld2lxP52pLWUsfavVzxOh0qw/wLtGdUYkhSpZKvRKn2vd1OwU15Enou+KQAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAICQ7j5S9eTlR8XsVZLR9+G1jkgDo6KmO/T6bGoT0+nZxtnnWbtbp55mze+dW07PLh36s7W7PczPzh47Zu0+eng2PdttN63d3Y7XfVQbbaRnN2w+y9q9dt2p6dl1p2+2djsv+6ln8l1TkrTiSP56Li/ln4OSdGRhwZrvGN1HtbbXwTU+yD/Jx8wCtrpRZ2RWaqXwSQEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBQBAyN/vXvHup64YVQduMllVFDLuGZfUK0bSs/U1p1i7i3X5KoqDjRXW7qee3GPND5U/iVWjzkGSlgf99OxzR45Yu4/NHU3P1grv2ldqVWt+Zs2q9Ow//fPbrN3nX3hhena0Pmbtnl7VS89uOMOr55g9ejg9+8z+Z63dvUr+7UqSjjTztRj1fv45K0k9o/anMTll7a4Yr5/uolcTk/rzT/hGAMA/LEIBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQEiXiRSl130ko3PIqEmSJPUr+Y6a4YoJb/nU6vxxrNtgrT5UzR/33kNeJ1B9ZNSaVzlIj/bKfFeOJC10WunZbr9r7Z6czvfIbNp0hrX7zHPOtuYvuDDfZ/SmN/+TtbtWy/dNtZba1u7JVTPp2cbUuLV7cXkxPdtdyM9KUt3oBJKk5UH+OV4Z9x7nqlNPTc+uHs33qUlSfy7fH9VaWrJ2Z/BJAQAQCAUAQCAUAACBUAAABEIBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAId99JK/7qDTme1Uvm9rj+T6j7rp8R4kkLU5Op2cP9fP9TpK02De6WFbmO34kqdVtWvPtdr6fqNP1Hmff6Mk6dePp1u5LLrkkPXvR2y6ydm/Y6HUlTc3ke7K6Pe8cNpv5PqOxRt3aPfvcbHr2D7/+lbX7yH89np6d6phdRi2v42lqzdr07Jat/8Pafc6G09Kz3YP7rd2LnXwfWFH3Opsy+KQAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIKRrLkrvLn31i3zetOsNa/fc+GR6dnnCq4t4pjdMz3ZVtXbXRkfTs83WkrX78OED1vxgmH+cpfmzw5Yt56ZnL/2Xf7F2/+u//nt69lSjikCSBvlTIklqNvN1BL1uz9pdLfJVIYcPedf+v/7z/6Rnl/bss3avGuafK/1G/nUsSWs2nWXNn37+1vTsmWd5u+u9/PU8eOyYtbtqvL9Vzfe3DD4pAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgpLuPBmb3Ua+W7wXqrznF2j1qdOvM18at3csHD6dnu22vz6a/dDQ9u7Dk9aUsLi5a87VqPT07s2qNtfsMo6Pm9DPOtnYX9fz1nJ3z+qPa7Y41XzH6vVaMef1e7dZyenZ+dtbaXTc6ns40+4Y66zamZ2vT09bu6dO9LquJqXwvUGU4sHYPlo3X27j3HtSu5V+b3caYtTuDTwoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAQrrmolfN11ZIUn1D/nb3tdv+3drdXH96enbfY49bu3vDI/njaLes3e12vnahvexVNPSaXkXDxJrJ9OzKmRlrd6Wer3RYavet3fOLzfRsYT5nKyqs+XGjumJhccHavX/vnvTs0cP556wk1Sbz9Q/rtq62do+tWpWera5YYe0eDox+DkmDTv712Wvna0UkqaN8LcYx4zgk6dm5/PVcNndn8EkBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAAAh3X1UXbveWnzqP29Lz674nxdZux/Z91x6drnTs3Y7fTmDYb7/RJJ6nXZ6dnnB68opVXrHYhx7u+/t7hf5c9gden027V7+eo6YXUYqvJ+RDj6Xfx7+9j9/be1+6s9PpGdXmd1Uq9asTc9OTKw8abtr9RFrt9sH5nSNNVv5Ti1JOnhgf3r20cce9Xbv25OerZtP8Qw+KQAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAI6ZqLxmmbrcW96TXp2Uef2mft/uMTf07PLi55t8ZXqvmcHPS9Co1Wczk92+13rN3V+qg1PyjSl14T06us3TOr89d+bHzc2j00ajFKs0Jj9sgha/7hX/3v9OyD/+s/rN3Var6/4JR1XgXN6YP8eVlf86ooRo4cSc+Oj09YuxeOzVnzSwtH07OHntlj7f7NL/PXfv+fH7d2jyp/fYp6vlImi08KAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAI6QKc+a7X8/PoI79Lz+4ZltbuYSPfl1PK212p5HOyHA6s3e12K7+78PJ6Ynq1Nb9x85b87Bmbrd2Nxlh6tmI+zpF6PT3bauW7piTpiT89Zs0/+odH0rPN5UVr94oV+ef4gQP7rd194/VWmNdneamZnl0xMWntLkvv9XZw/zPp2T88/LC1e8+T+T6jscLr4BobzT/HaxW6jwAAJxGhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACOmai70HDliLn9Lh9Gx39Vpr9/TUdHq21+tbuwfddn5336v+UFGkR8cmp6zVazdstOZXn7I+PVuv5W+7l6RaLf200sqpldZux5HDh6z5/c8+bc13WvnakqpRnyJJvV7+uVWYuxfm5tKzzzy919q9ak3+nMysXmPtlllZ8/RTf0nPHnh2n7W7blRXTIw1rN2NWr66omq8p2TxSQEAEAgFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBQBAIBQAACFdUnPU6HmRpGJyVf4g6iPe7iKfZUXF6wYZDAbOgVi7R8ZXpGdXmr0wE9Orrfl6fTQ9Ozkxae1evSp/7SvmOZyfz/f2HD1yxNrdbDat+aKafx5WjD4oSRpxXhPG60GSlhYX07Nu79XY2Fh6tlb3zkmr5V2fwwefS88OOt7726jRT+Q9SknDfK+SivxxZPFJAQAQCAUAQCAUAACBUAAABEIBABAIBQBAIBQAAIFQAAAEQgEAEAgFAEAgFAAAIV3L0Sy9jpqymu/kKAZG14ekodFPNCyt1eoZx1I3uowkaWbFRHp27YYzrN0bTt9kza9ftz49u+l071imp6fTs91e19rdd+a9p6yqZs/PyEgjPTvW8J4r1Wq+Mada8X62qxpdSZ221wm0uHgsPdvseF1Gi/P53ZI0d/hQerZnPs6q8u9BndJ8fzPeO436rTQ+KQAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACAQCgCAQCgAAAKhAAAI6Xvph2ZlgJzqin7fWl2U+e6KilkBMDo2np6dLLyTMrFyKj176mletcT69fnaCklaNZU/FqMVQZJ0bH4uPVvK6yHptvM1F0dnj1q7jx4+Ys23ltvp2bpRiSFJU1Mr07PNpSVrd7OZn2+2lq3dvTL/WjYbaDR3KF9bIUnzswfSs5W+V7dS1PJVFO4LqG+8rdTc9+UEPikAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACIQCACCku48qQ6PLSFLd6OQYVLwCj67R3eL0iEhej0xvMLB2t5byfSmzB57xdi/k+4Ykaf9IPT1bc3peJBVG18ug753DTjvfNzR/xOsyOvCsd86Hg156dnw836klSZ1W+qWpVqtp7V5aXEjPDo0uI0k68Gy+b8q5lpLUWvY6nspB/tjrZo9Zx3o79HaXxrEM3WKyBD4pAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAjpe+lHavnb7iVpYFQAdJbzt91L0lIzPz/XNCsAjPlKxat/GBkdSc8eeMa7Nb5S8fK9ahx7Ye4uyzI9O+h5NQqDXv551Wu3rN29bsean5ycSM+O1LzreWQ2XwFRrXrPw0ZjND174MCstXtxYT49OzBrYgrln1eSVLHqIrzrMzBeEsOqWXNhXM+hWRGUwScFAEAgFAAAgVAAAARCAQAQCAUAQCAUAACBUAAABEIBABAIBQBAIBQAAIFQAACEfKFR1cuPTr+bnj3Wyfe8SNLcYJg/DnndIPVGIz1bq9et3U7XS7vl9fY4fUOS35djGeYf57CXf55IXvdRUeafJ5JUFN457PXzXUnNpvf6GWusSM92O15n0+zhg+nZhWPz1u7hMN9lVZh9Qydzvqh516dSz791Vsz3icLomSvNPqgMPikAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgEAoAgEAoAAACoQAACOn7qXsVs+ZCRo3CeP6WfkmaWjGZnh2dyM9KUsW4xXx5ecnafezYsfRstfRqKNzKjapTW1LmayskqRjkqw7Kfn5WkobGbpk1F9ZuSUOjYaDnHYq6y830bLvtVaK0jWqRqlHnIEmV0qiWMCtoCvM9qKjk97u1L7WR0fys+dqUUc9Rms/xDD4pAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAAqEAAAiEAgAgFGVZGg0uAIA3Mj4pAAACoQAACIQCACAQCgCAQCgAAAKhAAAIhAIAIBAKAIBAKAAAwv8FF5Yf8QJAq/sAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAGbCAYAAAAr/4yjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAgoklEQVR4nO3da5CddbXn8fU8+9q9e/e9O+mEXDsBkuYQCRDEMIAaMYwnyEAGGMYqUcYqTaJyyqIsyxcEHdSScqwhqAVaThWD4rGM4jlT4EhQZobrcEwIJIQknSu5dCfdSd+79+XZ/3lBnTUTCbIWEi5nvp8qqkxn9eqn997P/u2nO/tnFEIIAgCAiMTv9gEAAN47CAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFnFHr16+XKIpkYGDgLe+IokjWr1//lj933bp1b/lrA/+/IRQAEXn66adl/fr1MjQ09G4fCvCuSr/bBwC8mcnJSUmnz+xD9emnn5Y777xTbrnlFmlubj6jXwt4L+NKAe9JtVpNpqamREQkn8+f8VAA8BpCAe+IgYEBueGGG6SxsVHa2trky1/+sj7pi/zfn/3/7Gc/k56eHsnlcvK73/1O/+7Pf6fwxBNPyEUXXST5fF66u7vlvvvu099fnM7DDz8s5513nuRyOenp6dHdIq/93uP2228XEZF58+ZJFEUSRZHs37//lGP7Szv+2eHDh+Wzn/2sTJs2Ted++tOfvm5uw4YN0tPTI/X19dLS0iIXXXSR/PznP9e/Hx0dldtuu03mzp0ruVxOOjs75WMf+5hs3rzZdoMDbxEvv/COuOGGG2Tu3Lny7W9/W5599lm555575OTJk/LAAw/ozB/+8Af55S9/KevWrZP29naZO3fuaXdt2bJFVq5cKV1dXXLnnXdKkiTyjW98Qzo6Ok47/+STT8qvf/1rWbNmjRSLRbnnnnvk+uuvl4MHD0pbW5tcd911smvXLnnooYfk+9//vrS3t4uInLLvzXaIiPT398sHP/hBDZGOjg559NFH5dZbb5WRkRG57bbbRETkxz/+sXzpS1+S1atXazi++OKL8txzz8nNN98sIiKf//zn5Ve/+pWsW7dOFi9eLIODg/Lkk0/Kjh07ZOnSpX/t3QG8sQCcQXfccUcQkXDNNdec8vE1a9YEEQlbt24NIYQgIiGO47B9+/bX7RCRcMcdd+ifV61aFerr68Phw4f1Y7t37w7pdDr8+UNaREI2mw29vb36sa1btwYRCRs2bNCP3X333UFEwr59+0779S07br311tDV1RUGBgZO+fybbropNDU1hYmJiRBCCJ/85CdDT0/P677O/6upqSmsXbv2L84AZwI/PsI7Yu3ataf8+Ytf/KKIiDzyyCP6sSuuuEIWL178F/ckSSKbNm2Sa6+9VmbMmKEfX7BggVx99dWn/ZwVK1ZId3e3/vn888+XxsZG2bt3r/n432xHCEE2btwoq1atkhCCDAwM6H8f//jHZXh4WH/009zcLIcOHZLnn3/+Db9ec3OzPPfcc3LkyBHzMQJvB0IB74iFCxee8ufu7m6J41h/bi/y2s/z38yxY8dkcnJSFixY8Lq/O93HRERmz579uo+1tLTIyZMn3/TrWXccP35choaG5P7775eOjo5T/vvMZz6jxy4i8tWvflUaGhpk2bJlsnDhQlm7dq089dRTp+z+7ne/K9u2bZNZs2bJsmXLZP369a4QA94qQgHvitP9Qriuru6MfK1UKnXajwfH/xPtm+2o1WoiIvKpT31KHnvssdP+t3z5chERWbRokezcuVN+8YtfyGWXXSYbN26Uyy67TO644w7de8MNN8jevXtlw4YNMmPGDLn77rulp6dHHn30UfMxA28Fv2jGO2L37t2nXAn09vZKrVZ7w18mv5HOzk7J5/PS29v7ur873ces3uhfLVl1dHRIsViUJElkxYoVbzpfKBTkxhtvlBtvvFHK5bJcd911ctddd8nXvvY1yefzIiLS1dUla9askTVr1sixY8dk6dKlctddd73hj8mAtwNXCnhH/OAHPzjlzxs2bBARcT/BpVIpWbFihTz88MOn/Ly9t7f3r3oVXSgURETe8juaU6mUXH/99bJx40bZtm3b6/7++PHj+r8HBwdP+btsNiuLFy+WEIJUKhVJkkSGh4dPmens7JQZM2ZIqVR6S8cHWHGlgHfEvn375JprrpGVK1fKM888Iw8++KDcfPPNsmTJEveu9evXy+9//3tZvny5fOELX5AkSeTee++V8847T1544YW3dHwXXnihiIh8/etfl5tuukkymYysWrVKw8LiO9/5jvzxj3+USy65RD73uc/J4sWL5cSJE7J582bZtGmTnDhxQkRErrrqKpk+fbosX75cpk2bJjt27JB7771XPvGJT0ixWJShoSE566yzZPXq1bJkyRJpaGiQTZs2yfPPPy/f+9733tL3B5i9q//2Cf/i/fM/SX355ZfD6tWrQ7FYDC0tLWHdunVhcnJS50TkDf8JpvzZP0kNIYTHH388XHDBBSGbzYbu7u7wk5/8JHzlK18J+Xz+dZ97ur1z5swJn/70p0/52De/+c0wc+bMEMfxKf881bOjv78/rF27NsyaNStkMpkwffr08NGPfjTcf//9OnPfffeFyy+/PLS1tYVcLhe6u7vD7bffHoaHh0MIIZRKpXD77beHJUuWhGKxGAqFQliyZEn44Q9/eNrbB3g7RSE4ftsGvIdde+21sn37dtm9e/e7fSjA+xa/U8D70uTk5Cl/3r17tzzyyCNy5ZVXvjsHBPwLwZUC3pe6urrklltukfnz58uBAwfkRz/6kZRKJdmyZcvr3hMBwI5fNON9aeXKlfLQQw9JX1+f5HI5ufTSS+Vb3/oWgQD8lbhSAAAofqcAAFCEAgBAmX+nMHhJu2txW7FqH46d2RTK9tlMzbc756g7iCq+3XFin804axec36aUHZ9Q9VZAOO7PyHnfx46fdk457x/v9xk5jqX9A77d/9r+Tu/Rv3/ItfpY/el7nE6nv9rg2t16fI95dl9/vWt3tbPLNS9v8P+vcTr9E+Ou1U0n+s2zueMDrt0TOfvzRKrqO/H/7Z6JN53hSgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAMrcfdTW4eyRiRydHBnn7rKj+6hm73kREZHE0X9TcxYOeXpKYuf/1UXiPBZHDZMkznb14DiWtLNvyHN3lny7q96qpNi+f+LAUdfuvf/lZ+bZpqqvWyccs9+f3nqvasp+mzS0zXXt7p013zX/UvqEebZpxHNCiDQ7HuLZjO/8qUb23rhU4nx+M+BKAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAICydykkJd/mxPO2cWeNQtWx2/fudRHP2/o9tRUiIo52Du9NIjVnvkeOCgjvsXjm08636afsy8fKvt1TFV8tRim2H0u15OuLiEcnzbMjBXstgohIQ81+2mcj351/PJM1zx6e3+Ha/XJd0TX/ywn7/Xlj317X7nRlyjybT3z3T83zvOKt8THgSgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAMpegjLuKe4RceVNZO9LeW21/bAlm/ftHhmyz5admeqqj/J2GfnmK2LvTJmo+fpvqsHeOZNJOY87tncIlWoZ1+6pqq/7qBrbi7WCp69LROrS9tt81FetI8Nl+/kWR3Wu3ZO5RvPsoayvTy3Z5usn+sSYfbal5Lt/io7eq4yzIi1XsZ+bSUz3EQDgDCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAIAiFAAAyt4XUVd0bq63z+YKvt2et6RXnO8xnxh17HbUbYjIlKMWY8rxNnoRkWrkq2goO14PlBLfsQSxH0vW+bqk6qjFKKd9u8vBNx95HoYpX41CumZ/3KZLvpqY0vzzzLO5mYtcu4cPvWqeDTv3uHZPH/Y9DkcdLRpzGn23YSaZsA8Xm1y740H7YyVJvPVDhq//tm8EALxvEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAlL28p6HNtXgqsncfheEx3+6BEfPs5GTVtTvlqFdJBWf3UdW+vBT7jjtJucYlxPaul1rwLY8cHUJV8XVTJZ6XMbGvD0qCbz4d5cyztYzvUFIVe3FPIdPq2p3/yErz7J722a7dx3ftMs+2bPM9xkdTu13zbek68+xsZ7dbQ8XeTxRS9seJiMjUxEnzbDRF9xEA4AwiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAMrc07Dr4LBrcWVyyD485as6qFUr5tk45cu9WmyvrsjWfLtLYn9bf+yoihARiau+2zBO2/enHJUY3mNJOasoIkfjRirj7P4Q++NKRESCoyrEeSgZR81FmDnDtXtfvsk8+7udR1y7J05MmWcvmrfItbtY89XhzCtNmmcLwT4rIpIaczxWJsddu0PN/lwbp5xVLpadb/tGAMD7FqEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQJmLfvr7RlyLcyFjnk35anskcvQZ5VL24xARScR+MDXxFdqEKJhno5Rvd6rmuxHTjhIhb71KKrLfP/nYd/9UQ2KeDRnf7lLO3k0lIpLK58yzcV3etTuqNJhnB7o6Xbt3jNjP5W07d7t2zznysnk2n7X3jImILKyMuuYLYyfMs+WU776vjtm7qTIV37mZiuzHUhPfY9yCKwUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAIAiFAAAyvw+80Tsb+kXEYki+9uv086ei8hRReGt0LAXUYhkHFURr+22b4+duzORs+bCUaORcvZcVJuy9tmZXa7ddZP2CoBcoc61eyyuuObTtbJ5tlL13Z+TaXvNxWjNfnuLiBx7dcg825Tz1dsUq/3m2a4h323SWR5yzSfBUReR9b0+LsX2cznEvvMnTuzPtYnnCcv69d/+lQCA9ytCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAIAydx+lUnnX4jjxHISvG0Qc3Uexc3fkOO7YGanptKPrJTj7UmJnV1LGfn8Wpre7dk912DuHko4m1+5wcNi+e8JxZ4pI4ugyEhEZr9rna7lW1+5So31+qFZw7c6n7LsXdPl2tzbYH1fp3gOu3cOjvqKf2tCQ/VgS38lcjUrm2SS298CJiERV89OyBMdzoRVXCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUOaSjZyznijt6ORIOXfXgj3LvK1Knk8Iwdc7Emr25SHl63lJ7HeliIhUc/ZOm7F8vWv3seFJ82y+vs21e6Jg7+3Jn+XrG2o8d75rfvqCbvtwcY5rt7TZ+6Ymjh93rT7eVzXPpg686to9enSfebZv4QLX7pH0LNd8evtu82zzkX7X7iTYb8PgfIKLHT1MSex7njB9/bd9IwDgfYtQAAAoQgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKHM3Qj6qOBfb39qdirKu3Z43dodMyrXbo1b1vX29FttrMbxvXg9Z3/c5Uk7Ms9Wq77VD4ZIV5tmF/+ZG125pbLHPFu1VHiIiUl/nm89kzKPVqbxr9aTj7pzZfY5rd2fFfn9O9PkqNH76m9+YZ5+eXXLtLjbNdM1fsewq82x47knX7uTANvus47wXEYkd9TlJ4ttt+vpv+0YAwPsWoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFCEAgBAmbuPcs6OjZCyl7dEsS+bomDfHTtmRUTSjjqjasbXUBRi+7HUgr2bSERkIlvvmi9ceKl5tv2KK127c+ctMc8OtMxz7T48NGye7d+937W7PGjfLSIyPjFknj3aN+Ha3Tdun19+5UWu3csuvdg82zCzw7X74mtWmmcf3PIPrt3bjx9zzXem7T1Zy5Ze5to9MXDSPBtP+I47LfZOqGrwtqS9Oa4UAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAIAiFAAAilAAAChzzUXaPioiIuWUoy/Ct1qiqn13KjiOQ0RSjlqMoXTFtTsj9mOpxEXX7qblH3HNV5Z+0Dz7vyd8FSdD/2OzebZav9e1+9ld9vl923pdu88aO+6a7+jKm2ePHhl07e6rm22eXfGxK1y7y2X7bMV318ulN3zKPPvEzn2u3VMn9rjmnx2yV1FkZy9y7Y6mnWueLe4bc+1uiR01F876IQuuFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoMytQ1E261ocBXsnRyrlKz8KFXt5Sy0Ort212N59NDll7yh5bXedeTZz8TLX7hMzZ7nmt2/Zap4d6h927W49a455tjLfed/H9vs+1+Z7zTNxYsA1L51nm0cz7QXX6ovPvdg8u2DpYtduT/dRdtS1WtK5BvPs5Vdd7drdu+G/uubzSdU8+/SeXa7dXS2N5tmOk+2u3cmo/Xmlrjzu2m3BlQIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAAZe4YqKWc+VGx11zEtci1OnYcSjlOXLvL7dPMs63di1y7p8btlRtD5/S4dv+p74RrPpvJmGdbpze7drd3Fs2zh1IV1+6pxD5fbLZXEYiINC27wDXfcu4C8+zyCy9x7V5ygX1+yvcQl5rj/Kn4WkiklrfPplt91R8fWn6+a37omT+aZwdLju4PETk4t9s8e97F/961u+PgfvPsyX98wLXbgisFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoc7NJCM78SOx9RsHZryJi7+2pFH39KidmzzHPtn7oQ67dE45vtHcy5drdfO55rvm+w6+YZ5PI1wsTGu3HPlouuXZfsMT+ff7t1Stduzua7Z1NIiJir7KSsr0KTEREqo754Nyddjy00lnf7hc2P2ee3fLAg67df9MwwzX/Sv+4efZkbJ8VEbng4/bH4SUfvdq1u27/fvPsU4/+1rXbgisFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAocxlPJPYuIxGROLIXrFRqjhIZEalk7IUsJ5s6XbtfStebZ8dfPenanW5vM8/mWnzHXRmbcM2/eLTPPDuV9vUTNYyMmmfHj1dcu9fffr15trnJ12VUcnYISeIYrfpWe86ItO/0kZxjdqT3Vdfuzf/xHvNs3Ys7XLsnO892zXed/QH77MULXLuXXvYR82xHZ5drd0Ndg3k2N2Oha7cFVwoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAlLnmQiq+99IHRwVAJfb1C9TaOsyzXauuc+3eKnnz7P6d9joHEZH6ctk82xyNuXbv2vmKaz4Zth9LSey1IiIi0xdMM89mi74KjcZiu3nW2fwhFV+Ti9QcD9vIWaGRcbxcy9jPYrdp02e65v/m0uXm2Ym53a7dnedf6JpvPP9882xxWqtrd9re4iOjY77niXyD/TmoNneea7cFVwoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFD21pSqr6NGylXz6MTcBa7VF//d1+3Dl3/EtfvX/22TefbYyQHX7o6K/TbJFH19Q2ODQ675xmN95tlUi73LSESkIVswz2Y7fN06ZcmYZ5OKa7UkzpdIiaPfy37Ur4kcx+KsbJLg+IRK1nejfGjtGvPs0OiIa7c0N7rGxxz3T81zZ4pIMmUv1qoF3wPxYP8h82x/1deRZsGVAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABlrrlIasG1uDRlr3Q46+b/4NotH77aPHok9tVF1Ld3mWfTLb2u3UnK/nb3asVXKzLeN+qar6ufNM8OZHxvpT+0a495dsGCRa7dVc/LGOdLHudD3FUvEWq+3Vtf2GGe3fzSZtfurhkd5tmrPnyVa7f9rBepFHy1FZGztiSdss+GrGNYRCqOY5ka9p0/z//PJ8yzh3a+4NptwZUCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAACUuftoouLLj1rzTPNs3eKLXbv7q/Y+o6PDJ127W1tazbP1hZxrdya2l+v0vXrQtbsUxl3zQ/lm8+xEaHDtbmqy787mfLehp57IWWXk5qnLeWn7dtfu//yffmKejes9jUMiUZ39XH7hTy+5dq/78t+ZZ2Pna9Ka79uUyXLZPDs1MeJcPmwefeqx/+5a/cpvHzbPtgbnjWLAlQIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAAZa65KDveMi4ikp83zTz75JYXXbu/v+G35tkrPrzEtXvBUvt8fZ2voiFVLplnyyXf2+6nTXN0LojI/sl68+wV/2qZa/fcBeeaZ+vr7MchIhI5ZsMZ7rmIHTf5P/7mYdfuzgb767WumfZzTURkPJk0z27ZsdO1e+Pf/8I8e/3qf+faPTRsP24RkYHhAfNsKky5dj/z+GPm2X2PP+3anRvcb56t62h07bbgSgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAMrcfZSIr/toqmbv+Tm4f5drdyFrPmwZGh517c5ls+bZ9tZm1+4dLx0wz1ZTvtt7MmlwzTd3dJhni01trt2joxPm2Y5W3+70GXwZ41194uSgfXb4hGv3jJaieXZg0Le7LDXz7JzWgmv3n57bbJ5ddM4i1+7OmfNd89mc/Vzet9PX8bT9he3m2ZbJPa7dHbO6zLOFjP17tOJKAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAIhQAAIpQAAAoQgEAoAgFAIAy90WMRSnX4tLJEfOs/U33r/nAObPMs+NR5NodomCeLdTXuXZXq4l5diqVc+1unmavrRARKbTPMc8mie82nBgbN8/Onn2Wa7fnSLyPK/u985qj/f3m2Uw249pdaLLXXNQ3+ipOhkbt98/w4EnX7pEBe+XGnle2uXZPnzPXNR9X7efy/j37XbvrHVUuxXrfc2fe88itVl27LbhSAAAoQgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKAIBQCAMncfHav6ulvKk/bZiZKvG6RSsXea1NK+BpzJySnzrKNaRUREhqcq5tm4seDa3dA+3TU/Z+588+z0ljbX7tjTH5XPu3a7ml6c5UeR8yVSkti/QJL2nT9RKmuebW5udu0eLx83zzaky67dU2LvPNu3t9e1u7v/gGt+aNR+Lu/d86prd3rQ/n1WqsdcuyeC41Fe1+rabcGVAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAlLn7aKTO3sUiIiK5nHk0NTDkWh03zTDPVmq+zpkxR19Koc3XCTT77LPNsy3TfZ0mi3sWueZ7FvaYZxuzvtcOGcd45Lt7JB3ZZxPfagnOTyjki/ZhZ3fYVM1+Iy6aPdO1u6PL3pO17Z9edO0er42ZZ/f3HXHt3rXddyzDY/YCtgNHDrt2tw7Z759s5OsxE/tTp0hdo2+3AVcKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAAJS55kKa877Ncc08mul71bV67pIl5tn9zhqF2uhJ8+zgiUHX7thRvDA+Mu7ave2lHa75I/vst3lzfZ1rdzZlv9EzDb76lFSw91zUKr7eipQ4OjREZPT4MftwqeLaXZezn5o7Xt7l2p2KUubZQ0f6XbtLUjLP7h0bdu2Wx/+Xa3wksVfWjCUTrt3tZ9lfT0+NnePaHeqq5tl0k2u1CVcKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAEARCgAARSgAABShAABQ5oKVON/gWpzPNZpnx2tl1+6+7X8yz052TnPt3tq/0zy746ivF+bArlHz7ETRl9f2xpnXNMX2/ptC8O2uJPbdkrF3/IiIZB2vYzKJ78Djmm++ud7eUVMa991DFy6aYZ6t1Xy34dE+e2fX2fOmu3Zv320/J1KZo67dO/fYb28RkbEG+/15bvuAa3dd7OiCa/R1h9Xy9sd4CL7dFlwpAAAUoQAAUIQCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFDm98dnWjtdi7PpevNsLVdw7S5V7LUYOw/63r5+oGqvIygWcq7dS5Z0mWf7Riddu8ckcc03iKOKwilK2esFMsF33HHFMV/y1SJEU775psYR82y54Kty2bL1RfNszzk9rt1zZtirKw692uvaXZEJ82y7+M7NzLwO13yzRObZ0aji2j2Usj8OM3nfudZc32yejaq+ihMLrhQAAIpQAAAoQgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAKEIBAKDMxRklZz9RiDLm2XK+6Nqdmr3QPNueb3Ttbq3LmmebW3zHHRy9PcM797l2T034upJyDXXm2Uza192ST9fMswXxHXdStnfURCVfr1Jl0td/kyo02GfLvu6jg/32Dq7Bf9rl2h3ZT3s5Pnjctbtwlv3+zOamuXZHNV/Pz1hk7z6qRPY+KBGRfMr+/BZlfa+9szn780ot7XuMW3ClAABQhAIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAAAUoQAAUIQCAECZ3zdebWlxLQ41e41Ceu45rt2F+faai/M7Z7p2TzkaHUYmJ1y7B/sHzbOZ+ibX7jmtna75nKPOI0nZ6wJEROoy9tca9aHq2h1q9gqNULHPiohUJqd88+Vx82zVMSsi0py3z+ZS9nNNRCSu2m+XJOOr/qhm7CdQOtvs2i1V3/dZchx6LfLVRcTBXrnREtnPNRGRdGK/DSux73FlwZUCAEARCgAARSgAABShAABQhAIAQBEKAABFKAAAFKEAAFCEAgBAEQoAAEUoAABUFEII7/ZBAADeG7hSAAAoQgEAoAgFAIAiFAAAilAAAChCAQCgCAUAgCIUAACKUAAAqP8D/rHUR904mSkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAGbCAYAAAAr/4yjAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAVNUlEQVR4nO3dfYyedZno8euZ55mZTqcvUKhtFShjW5YFlWW7K8JKWV8iqIemKnIIQeBIQMNLZDecnJhzEsBA9Q/TiAYJPSnV0HiSXaueXZG4WwMbDerBgwKCQBHKm5S29H3eOjPPff4gXnGWt/tm+dnu2c8nMbHTa3695mmf5zu3M3PbqqqqCgCIiJ6DvQAAhw5RACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRR4N+dY489Ni6++OJi51988cVx7LHHFjsfDmWiAG+i1atXx/e+972DvQa8YZ2DvQA09eijj0ZPz6H5+czq1avjnHPOiVWrVh3sVeANEQX+3env7z/YK8D/tw7NT7f4D+m6666LVqsVjzzySJx77rkxZ86cOOKII+Jzn/tcjI2N5dwffk2hqqp43/veF/Pnz49t27blzIEDB+Kd73xnLFmyJIaHh/PtGzZsiOXLl8fAwEDMmzcvzjvvvHjmmWded7cvf/nLcdppp8URRxwRAwMDsXz58vj2t789babVasXw8HB885vfjFarFa1Wa9rXPp577rn49Kc/HQsWLIj+/v448cQT47bbbnuDjxaUIQoccs4999wYGxuLL37xi/GRj3wkvvrVr8Zll132irOtVituu+22GBsbi89+9rP59muvvTYeeuihWL9+fQwODkZExI033hgXXnhhLFu2LNasWRNXX311/OhHP4oVK1bE7t27X3Onm266KU4++eT4whe+EKtXr45OpxOf/OQn44477siZ22+/Pfr7++P000+P22+/PW6//fb4zGc+ExERL7zwQrznPe+JTZs2xZVXXhk33XRTLF26NC655JL4yle+8m97wODNVMEh4tprr60iolq5cuW0t19++eVVRFT3339/VVVVtXjx4uqiiy6aNnPrrbdWEVFt2LCh+tnPfla12+3q6quvzt/fsmVL1W63qxtvvHHa+z344INVp9OZ9vaLLrqoWrx48bS5kZGRab8+cOBA9Y53vKN6//vfP+3tg4ODL9utqqrqkksuqRYtWlTt2LFj2tvPO++8au7cuS87Hw4WVwoccq644oppv77qqqsiIuIHP/jBq77PZZddFmeeeWZcddVV8alPfSqWLFkSq1evzt//zne+E91uN84999zYsWNH/mfhwoWxbNmyuOuuu15zp4GBgfzvu3btij179sTpp58e99133+t+PFVVxcaNG+Pss8+Oqqqm/flnnnlm7Nmzp9Y58MfgC80ccpYtWzbt10uWLImenp7YsmXLa77funXrYsmSJbF58+a45557pr2Qb968OaqqetnZv9fb2/uaZ3//+9+PG264IX71q1/F+Ph4vr3Var3ORxOxffv22L17d6xduzbWrl37ijN/+PUQOJhEgUNenRfeiIi77747X7AffPDBOPXUU/P3ut1utFqtuPPOO6Pdbr/sfWfNmvWq5/74xz+OlStXxooVK+LrX/96LFq0KHp7e2P9+vXxrW9963X36na7ERFxwQUXxEUXXfSKM+9617te9xz4YxAFDjmbN2+OoaGh/PXjjz8e3W73NX/K+Pnnn4+rrroqPvShD0VfX19cc801ceaZZ8bixYsj4qWrjaqqYmhoKI477rhG+2zcuDFmzJgRP/zhD6d9O+z69etfNvtKAZs/f37Mnj07pqam4oMf/GCjPxv+2HxNgUPOzTffPO3XX/va1yIi4sMf/vCrvs+ll14a3W431q1bF2vXro1OpxOXXHJJVFUVEREf//jHo91ux/XXX59v+72qquLFF1981bPb7Xa0Wq2YmprKt23ZsuUVf3J5cHDwZd/J1G634xOf+ERs3Lgxfv3rX7/sfbZv3/6qfzb8sblS4JDz5JNPxsqVK+Oss86Kn/70p7Fhw4Y4//zz46STTnrF+fXr18cdd9wR3/jGN+Koo46KiJdCcsEFF8Qtt9wSl19+eSxZsiRuuOGG+PznPx9btmyJVatWxezZs+PJJ5+M7373u3HZZZfFNddc84rnf/SjH401a9bEWWedFeeff35s27Ytbr755li6dGk88MAD02aXL18emzZtijVr1sRb3/rWGBoailNOOSW+9KUvxV133RWnnHJKXHrppXHCCSfEzp0747777otNmzbFzp0739wHEd6og/idTzDN778l9eGHH67OOeecavbs2dXhhx9eXXnlldXo6GjO/eG3pD7zzDPV3Llzq7PPPvtl533sYx+rBgcHqyeeeCLftnHjxuq9731vNTg4WA0ODlbHH398dcUVV1SPPvpozrzSt6SuW7euWrZsWdXf318df/zx1fr163PfP/TII49UK1asqAYGBqqImPbtqS+88EJ1xRVXVEcffXTV29tbLVy4sPrABz5QrV279t/wqMGbq1VV/+paGg6S6667Lq6//vrYvn17HHnkkQd7HfgPydcUAEiiAEASBQCSrykAkFwpAJBEAYBU+4fXfnHanGJLLJjZ7H/B6ulOFtokYt5gvfvsRER0qwPF9ij3Eb5k10i39uzwgYKfOzT8v9Ucnqy/d1M7hsudPda3qNH8i4MzCm0SMbzr1X96+1/bN17uMXl24vBiZ0dEbOnUfy4Pjo+9/tAbtGB4T6P5sSj3mK/fPv66M64UAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoApE7dwT+bPV5siZGeqUbz2/ZXhTaJ2LGnfifHx7rF9hg90C52dkTEnsn65++bLLhIu9VofLRqNt/EroKP+djuvY3mX+x5sdAmEYMNHsLRKPd4PzZjRrGzIyJ+8uzztWff11vu43x7T7PXiVa33OtbHa4UAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQBS7XsfPb634L14ZsxuND5c8P43+3fuqT17YLy32B77psp9jBERo1X9zwdGyt3iKTqtZh/neKvZfbKaGC74cU5WzfZuF7zn0O4Zh9WeHe+bVWyPkZ27ip0dETFU1X/NWljw3ke9Pc1eJ/rHC/5DrMGVAgBJFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgNSpO/hYe16xJcZ3TTSa379ruNAmEd3JVu3Zqtsutsf+mCx2dkTEVE/9zweqqtznDj1V1Wh+qrf+309jrXJ/n93az7SX9Fflduk95fTas4/tHi22x9yRXxQ7OyJidqdbe/aYOQPF9uhtzWw0P/rC3kKb1ONKAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgCk2j98/8BvdxZbojXZ7FYHPa1ytzqoeurfXqBdldujVfCWCxERvQ0+H+gU3KWv4d/lVMPbRTTRLXgLjcneZot3+8rdduEXu+vfJuaeB35TbI+/GhwsdnZExOKJfbVnO6MjxfboGT3QbD6a3fbnzeZKAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQOrUHRydKNePWVXVaL5TMGXdbv3Zvp7aD19jrcK5Hui0a8925s4stkd3/mGN5ucMT5RZJCJGqgPFzp6oehvNP987WGiTiM1bXqg9O9CaLLbH0Mh4sbMjIgbG99eeHW81ew1qYqzBa0pERKuq/9wswZUCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoApE7dwf5oFVuit+HR5TaJ6EZVf7hqMNtQVbjX3b6ZtWf39g0U22Mq+hvNH5g9p9AmEUf95UnFzl7yzlMbzT+yf6TQJhFP3PdY7dnH7/s/xfZ4une02NkREQPbt9eenbtvd7E9plpTjearbslXuNfnSgGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoApNq3uRhodcst0eptNN/tKXd7iZ5W/bOnWuV+HL3bLnybi9mH1Z6dfdLJxfaYd8oZjebHO7X/yTa2bXys2Nm/fWZvo/mRnnLPt7+55m9rz3771luK7bHum7cWOzsi4vTjTqg9O9iztdgeB3b+rtF8q1vu9a0OVwoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQBSq6qqqs7gmrfMKrZET6fT7B26rTKLRES7wey+9lSxPXoXLC52dkREz7v/qvbsY7v2FdtjJHobzd/zf39RaJOIvtG9xc7eM7y/0fw7Tjmj0CYRX77lf9aeff6pJ4rtceF5q4qdHRExFd3as3++8Jhie8x+7vFG83O7o4U2ibht5+uf7UoBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRRACDVvulQq93w/kRNluhpdvZENVlok4jR/oHas1P95e4HtW3ewmJnR0Q8+cCDtWcHOn3F9hg//PBG8z295e571ekbLHb2ny0/qdH8ORf+l0KbRIxM1H/+zJt/VLE9/vwvTit2dkTEr3/+k/qz254rtseJRzW7j1nfrq2FNqnHlQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYDUqT1ZletHt9tsfrR/ZplFIuI3846sPds++u3F9tjXN6vY2RERT+56qPbszPZUsT3GfjfcaP6//ff/UWiTiPf/9RnFzh6cdVij+dGxA2UWiYjJyfpPuL+/bW2xParfbC52dkTEcTPn154dOGGo2B4fXfXJRvPP/tMdhTapx5UCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGA1Kk92S23xGg11Wi+d8lxhTaJ2N6ZWXt227O7i+3R6t9T7OyIiG2/21F7dta8I4vt0dff32h+8bF/UmiTiN/tGil29uTWvY3mF8yfV2iTiIFOu/ZsZ2Si2B5Df3pCsbMjIha/+7Tas4e9pdzjPdX09W3RokKb1ONKAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAg1b73UXfiQLkl3v2eRvOnXvlfC20S8fN/vLP27Jbn7i22x9TwvmJnR0SM7t1de7Z3YLDYHgsWL200PzxR7iZcE/vHi509s39Go/l/3rSp0CYR7fp3PIu/WLmq2B5HnfreYmdHRAxP1n/Nmiz4fBvfU/8+YxER9z/2UKFN6nGlAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIHXqDg63ZhRb4u0rzmo0/8COfYU2iVjwtmNqz24efKDYHk9vfqzY2RERE1HVn22V+9xhsmp29mTVLbRJRF+UO/sfvvN3jeY3rF9baJNo8KyPWLnqPxdbY9UnLih2dkTE1t8+Unt2ePf2Ynvc/YP/3Wh+870/L7RJPa4UAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEiiAECq/QPvU33tYkvc8S93N5r/zb6JMotExF+f/Z9qz3Y6De4X0FC3O1ns7IiIgbnza8/+5WlnFNvjqGOObTQ/Y0a5263s37e72Nn/sPF/NZo/fPasQptETHTrP3/u/uc7i+0xNHR8sbMjIvaPjdaevf9n/1Jsj1/+9CeN5hfN6i+0ST2uFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKAKRO3cGnRseLLbFt585G83MOW1Bok4hu1a0922qVa+rcBW8tdnZExNCfnlx7dkb/zGJ7LFu6tNH8xMREoU0itm7dWuzsqapqNL9wwaJCm0Ts3be/9uz2558ttsdDD95b7OyIiPlvW1x79uH7f1lsjzl9rUbz/T3N5t9srhQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgdeoObivYj+7+sWbzcwotEhF79+yvPXvYwrcV2+PoE04sdnZExKmnnVF79phFC4rtMTDY32i+v7ev0CYRM/oGip3d0272cR52xFsKbRIxtPRPas/eOzpSbI8Hf3VvsbMjIuY+s6X27LZtzxfbY2BqotF8u/7LchGuFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAqv3z1O2BcrcX6B/Z12h+anRPoU0innrovtqzTz/9dLE9Jquq2NkREY/8sv4tBmb0NbtFQxOd/mZn91StQptEdMeb3W6lib0vbms0P7y33K1FdmzdWnt2bGy42B5bnnqs2NkREVXUfw5VU5PF9hjraXbbipGeg/u5uisFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIASRQASKIAQBIFAJIoAJBEAYAkCgAkUQAgiQIAqVN3cGDWnGJL9Hb6G81vfvapQptEbH3s4dqzg4ODxfYYHR8vdnZExO5Wq+j5tXUnG41Xk1OFFomYOjBa7OzJyQON5rdvPaLQJhFHHzVUe/bpLSPF9piaLPtvvBX1/4232gX36Gv2XOuZ0Vtok5p//kH90wE4pIgCAEkUAEiiAEASBQCSKACQRAGAJAoAJFEAIIkCAEkUAEi173003OkrtkRrcG6j+aHjTiq0ScTiWbNrz05NNbtvTxNPPv54sbMjInr76/99djrlPnfoTDW7J1B3YqLQJhFjw/uLnT06srfR/IvD5e7DtKfBv619o+X2mDE4q9jZL2lw76OC9wKb0d/s3m7tjnsfAXCIEAUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkEQBgCQKACRRACC1qqqqDvYSABwaXCkAkEQBgCQKACRRACCJAgBJFABIogBAEgUAkigAkP4fmfGDn+gMtkoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "i = 6\n",
    "\n",
    "n_test = X_test.shape[0]\n",
    "plt.title('Original')\n",
    "plt.axis('off')\n",
    "plt.imshow(X_test[i])\n",
    "plt.show()\n",
    "for _ in range(len(corruption)):\n",
    "    plt.title(corruption[_])\n",
    "    plt.axis('off')\n",
    "    plt.imshow(X_corr[n_test * _+ i])\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Detect drift with a TensorFlow classifier\n",
    "\n",
    "#### Single fold\n",
    "\n",
    "We use a simple classification model and try to distinguish between the reference data and the corrupted test sets. The detector defaults to `binarize=False` which means a Kolmogorov-Smirnov test will be used to test for significant disparity between continuous model predictions (e.g. probabilities or logits). Initially we'll test at a significance level of $p=0.05$, use $75$% of the shuffled reference and test data for training and evaluate the detector on the remaining $25$%. We only train for 1 epoch."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.keras.layers import Conv2D, Dense, Flatten, Input\n",
    "\n",
    "tf.random.set_seed(0)\n",
    "\n",
    "model = tf.keras.Sequential(\n",
    "  [\n",
    "      Input(shape=(32, 32, 3)),\n",
    "      Conv2D(8, 4, strides=2, padding='same', activation=tf.nn.relu),\n",
    "      Conv2D(16, 4, strides=2, padding='same', activation=tf.nn.relu),\n",
    "      Conv2D(32, 4, strides=2, padding='same', activation=tf.nn.relu),\n",
    "      Flatten(),\n",
    "      Dense(2, activation='softmax')\n",
    "  ]\n",
    ")\n",
    "\n",
    "cd = ClassifierDrift(X_ref, model, p_val=.05, train_size=.75, epochs=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If needed, the detector can be saved and loaded with `save_detector` and `load_detector`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:Compiled the loaded model, but the compiled metrics have yet to be built. `model.compile_metrics` will be empty until you train or evaluate the model.\n",
      "WARNING:tensorflow:No training configuration found in the save file, so the model was *not* compiled. Compile it manually.\n"
     ]
    }
   ],
   "source": [
    "from alibi_detect.saving import save_detector, load_detector\n",
    "# Save detector\n",
    "filepath = 'tf_detector'\n",
    "save_detector(cd, filepath)\n",
    "\n",
    "# Load detector\n",
    "cd = load_detector(filepath)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check whether the detector thinks drift occurred on the different test sets and time the prediction calls:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from timeit import default_timer as timer\n",
    "\n",
    "labels = ['No!', 'Yes!']\n",
    "\n",
    "def make_predictions(cd, x_h0, x_corr, corruption):\n",
    "    t = timer()\n",
    "    preds = cd.predict(x_h0)\n",
    "    dt = timer() - t\n",
    "    print('No corruption')\n",
    "    print('Drift? {}'.format(labels[preds['data']['is_drift']]))\n",
    "    print(f'p-value: {preds[\"data\"][\"p_val\"]:.3f}')\n",
    "    print(f'Time (s) {dt:.3f}')\n",
    "    \n",
    "    if isinstance(x_corr, list):\n",
    "        for x, c in zip(x_corr, corruption):\n",
    "            t = timer()\n",
    "            preds = cd.predict(x)\n",
    "            dt = timer() - t\n",
    "            print('')\n",
    "            print(f'Corruption type: {c}')\n",
    "            print('Drift? {}'.format(labels[preds['data']['is_drift']]))\n",
    "            print(f'p-value: {preds[\"data\"][\"p_val\"]:.3f}')\n",
    "            print(f'Time (s) {dt:.3f}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "No corruption\n",
      "Drift? No!\n",
      "p-value: 0.963\n",
      "Time (s) 3.412\n",
      "\n",
      "Corruption type: gaussian_noise\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 3.194\n",
      "\n",
      "Corruption type: motion_blur\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 3.233\n",
      "\n",
      "Corruption type: brightness\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 3.162\n",
      "\n",
      "Corruption type: pixelate\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 3.210\n"
     ]
    }
   ],
   "source": [
    "make_predictions(cd, X_h0, X_c, corruption)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As expected, drift was only detected on the corrupted datasets and the classifier could easily distinguish the corrupted from the reference data.\n",
    "\n",
    "#### Use all the available data via cross-validation\n",
    "\n",
    "So far we've only used $25$% of the data to detect the drift since $75$% is used for training purposes. At the cost of additional training time we can however leverage all the data via stratified cross-validation. We just need to set the number of folds and keep everything else the same. So for each test set `n_folds` models are trained, and the out-of-fold predictions combined for the significance test:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Both `n_folds` and `train_size` specified. By default `n_folds` is used.\n"
     ]
    }
   ],
   "source": [
    "cd = ClassifierDrift(X_ref, model, p_val=.05, n_folds=5, epochs=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "No corruption\n",
      "Drift? No!\n",
      "p-value: 0.738\n",
      "Time (s) 11.437\n",
      "\n",
      "Corruption type: gaussian_noise\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 16.781\n",
      "\n",
      "Corruption type: motion_blur\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 19.431\n",
      "\n",
      "Corruption type: brightness\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 17.203\n",
      "\n",
      "Corruption type: pixelate\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 16.747\n"
     ]
    }
   ],
   "source": [
    "make_predictions(cd, X_h0, X_c, corruption)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Detecting drift with a learned kernel\n",
    "\n",
    "An alternative to training a classifier to output high probabilities for instances from the test window and low probabilities for instances from the reference window is to learn a kernel that outputs high similarities between instances from the same window and low similarities between instances from different windows. The kernel may then be used within an [MMD-test](https://docs.seldon.io/projects/alibi-detect/en/stable/examples/cd_mmd_cifar10.html) for drift. [Liu et al. (2020)](https://arxiv.org/abs/2002.09116) propose this learned approach and note that it is in fact a generalisation of the above classifier-based method. However, in this case we can train the kernel to directly optimise an estimate of the detector's power, which can result in superior performance.\n",
    "\n",
    "### Detect drift with a learned PyTorch kernel\n",
    "\n",
    "Any differentiable *Pytorch* or *TensorFlow* module that takes as input two instances and outputs a scalar (representing similarity) can be used as the kernel for this drift detector. However, in order to ensure that MMD=0 implies no-drift the kernel should satify a *characteristic* property. This can be guarenteed by defining a kernel as $$k(x,y)=(1-\\epsilon)*k_a(\\Phi(x), \\Phi(y)) + \\epsilon*k_b(x,y),$$ where $\\Phi$ is a learnable projection, $k_a$ and $k_b$ are simple characteristic kernels (such as a [Gaussian RBF](https://en.wikipedia.org/wiki/Radial_basis_function_kernel), and $\\epsilon>0$ is a small constant. By letting $\\Phi$ be very flexible we can learn powerful kernels in this manner.\n",
    "\n",
    "This can be implemented as shown below. We use *Pytorch* instead of *TensorFlow* this time for the sake of variety. Because we are dealing with images we give our projection $\\Phi$ a convolutional architecture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "import torch.nn as nn\n",
    "\n",
    "# set random seed and device\n",
    "seed = 0\n",
    "torch.manual_seed(seed)\n",
    "torch.cuda.manual_seed(seed)\n",
    "device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\n",
    "\n",
    "# define the projection\n",
    "proj = nn.Sequential(\n",
    "    nn.Conv2d(3, 8, 4, stride=2, padding=0),\n",
    "    nn.ReLU(),\n",
    "    nn.Conv2d(8, 16, 4, stride=2, padding=0),\n",
    "    nn.ReLU(),\n",
    "    nn.Conv2d(16, 32, 4, stride=2, padding=0),\n",
    "    nn.ReLU(),\n",
    "    nn.Flatten(),\n",
    ").to(device)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We may then specify a `DeepKernel` in the following manner. By default `GaussianRBF` kernels are used for $k_a$ and $k_b$ and here we specify $\\epsilon=0.01$, but we could alternatively set `eps='trainable'`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "from alibi_detect.utils.pytorch.kernels import DeepKernel\n",
    "kernel = DeepKernel(proj, eps=0.01)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since our *PyTorch* encoder expects the images in a *(batch size, channels, height, width)* format, we transpose the data. Note that this step could also be passed to the drift detector via the `preprocess_fn` kwarg:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(5000, 3, 32, 32) (5000, 3, 32, 32) (10000, 3, 32, 32)\n"
     ]
    }
   ],
   "source": [
    "def permute_c(x):\n",
    "    return np.transpose(x.astype(np.float32), (0, 3, 1, 2))\n",
    "\n",
    "X_ref_pt = permute_c(X_ref)\n",
    "X_h0_pt = permute_c(X_h0)\n",
    "X_c_pt = [permute_c(xc) for xc in X_c]\n",
    "print(X_ref_pt.shape, X_h0_pt.shape, X_c_pt[0].shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We then pass the kernel to the `LearnedKernelDrift` detector. By default $75\\%$ of the data is used to train the kernel and the MMD-test is performed on the other $25\\%$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "from alibi_detect.cd import LearnedKernelDrift\n",
    "cd = LearnedKernelDrift(X_ref_pt, kernel, backend='pytorch', p_val=.05, epochs=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, the detector can be saved and loaded:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "from alibi_detect.saving import save_detector, load_detector\n",
    "# Save detector\n",
    "filepath = 'torch_detector'\n",
    "save_detector(cd, filepath)\n",
    "\n",
    "# Load detector\n",
    "cd = load_detector(filepath)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, lets make some predictions with the detector:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "No corruption\n",
      "Drift? No!\n",
      "p-value: 0.890\n",
      "Time (s) 1.450\n",
      "\n",
      "Corruption type: gaussian_noise\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 1.760\n",
      "\n",
      "Corruption type: motion_blur\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 1.765\n",
      "\n",
      "Corruption type: brightness\n",
      "Drift? Yes!\n",
      "p-value: 0.000\n",
      "Time (s) 1.724\n",
      "\n",
      "Corruption type: pixelate\n",
      "Drift? No!\n",
      "p-value: 0.060\n",
      "Time (s) 1.739\n"
     ]
    }
   ],
   "source": [
    "make_predictions(cd, X_h0_pt, X_c_pt, corruption)"
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}