Skin_Lesion_Analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Risk-Aware Machine Learning Classifier for Skin Lesion Diagnosis\n",
    "\n",
    "### Aryan Mobiny, Aditi Singh, Hien Van Nguyen, Journal of Clinical Medicine, 2019"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"imgs/fig1.png\" width=\"5000\" height=\"100\" >"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 0. Import the required packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Using TensorFlow backend.\n"
     ]
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "import os\n",
    "import h5py\n",
    "import numpy as np\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "from keras import Model\n",
    "from models import backbone\n",
    "from config import options\n",
    "from utils.py_utils import one_hot, softmax, load_data\n",
    "from utils.eval_utils import uncertainty_fraction_removal, normalized_uncertainty_toleration_removal\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Load the data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"imgs/fig2.png\" width=\"500\" height=\"100\" >"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## skin lesion categories"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Melanoma (MEL) \n",
    "- Melanocytic Nevi (NV) \n",
    "- Basal Cell Carcinoma (BCC) \n",
    "- Actinic Keratoses and Intraepithelial Carcinoma (AKIEC) \n",
    "- Benign Keratosis (BKL) \n",
    "- Dermatofibroma (DF) \n",
    "- Vascular (VASC) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "BASE_DIR = os.getcwd()\n",
    "data_dir = os.path.join(BASE_DIR, 'preprocessed_data.h5')\n",
    "x_train, y_train, x_valid, y_valid = load_data(data_dir)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Training images: (8010, 224, 224, 3)\n",
      "Training labels: (8010, 7)\n",
      "\n",
      "Validation images: (2005, 224, 224, 3)\n",
      "Validation labels: (2005, 7)\n"
     ]
    }
   ],
   "source": [
    "print('Training images: {}'.format(x_train.shape))\n",
    "print('Training labels: {}\\n'.format(y_train.shape))\n",
    "print('Validation images: {}'.format(x_valid.shape))\n",
    "print('Validation labels: {}'.format(y_valid.shape))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2. Create the Model\n",
    "\n",
    "### DenseNet-169\n",
    "\n",
    "\n",
    "\n",
    "<img src=\"imgs/fig3.png\" width=\"800\" height=\"100\" >\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "### In Keras:\n",
    "\n",
    "<img src=\"imgs/fig4.PNG\" width=\"800\" height=\"100\" >\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.1. Model parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# These can also be set up and imported from config.py\n",
    "backbone_name = 'densenet169'\n",
    "num_cls = 7\n",
    "dropout_rate = 0.5\n",
    "num_dense_layers = 1\n",
    "num_dense_units = 128\n",
    "pooling = 'avg'\n",
    "dense_layer_regularizer = 'L2'\n",
    "model_name = 'lesion_densenet169'\n",
    "lr = 1e-4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "save_dir = os.path.join(base_dir, 'default_save_dir')\n",
    "if not os.path.exists(save_dir):\n",
    "    os.mkdir(save_dir)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from models import backbone\n",
    "from utils.eval_utils import compute_class_weights\n",
    "from utils.callback_utils import config_cls_callbacks\n",
    "\n",
    "callbacks = config_cls_callbacks()\n",
    "\n",
    "model = backbone(backbone_name).classification_model(input_shape=x_train.shape[1:],\n",
    "                                                     num_classes=options.num_classes,\n",
    "                                                     num_dense_layers=options.num_dense_layers,\n",
    "                                                     num_dense_units=options.num_dense_units,\n",
    "                                                     pooling=options.pooling,\n",
    "                                                     dropout_rate=options.dropout_rate,\n",
    "                                                     kernel_regularizer='L2',\n",
    "                                                     save_to=save_dir,\n",
    "                                                     print_model_summary=False,\n",
    "                                                     plot_model_summary=False,\n",
    "                                                     lr=options.lr)\n",
    "\n",
    "class_weights = compute_class_weights(y_train, wt_type='balanced')\n",
    "\n",
    "batch_size = 2\n",
    "\n",
    "model.fit(x=x_train,\n",
    "          y=y_train,\n",
    "          batch_size=batch_size,\n",
    "          epochs=1000,\n",
    "          verbose=1,\n",
    "          validation_data=(x_valid, y_valid),\n",
    "          class_weight=class_weights,\n",
    "          shuffle=True,\n",
    "          callbacks=callbacks)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.2. Start the Evaluation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2.1. Normal Evaluation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Don't run this; run the following cell to load the results\n",
    "model = backbone(backbone_name).classification_model(load_weights_from=model_name,\n",
    "                                                     num_classes=num_cls,\n",
    "                                                     num_dense_layers=num_dense_layers,\n",
    "                                                     num_dense_units=num_dense_units,\n",
    "                                                     dropout_rate=0.,\n",
    "                                                     pooling=pooling,\n",
    "                                                     kernel_regularizer=dense_layer_regularizer)\n",
    "\n",
    "y_pred = model.predict(x_valid)\n",
    "y_prob = softmax(y_pred)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load from the saved file\n",
    "file_name = os.path.join(base_dir, 'predictions.h5')\n",
    "h5f = h5py.File(file_name, 'r')\n",
    "y_prob = h5f['y_prob'][:]\n",
    "h5f.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2.2. MC-Evaluation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "monte_carlo_simulations = 10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Don't run this; run the following cell to load the results\n",
    "model = backbone(backbone_name).classification_model(load_weights_from=model_name,\n",
    "                                                     num_classes=num_cls,\n",
    "                                                     num_dense_layers=num_dense_layers,\n",
    "                                                     num_dense_units=num_dense_units,\n",
    "                                                     pooling=pooling,\n",
    "                                                     dropout_rate=dropout_rate,\n",
    "                                                     kernel_regularizer=dense_layer_regularizer)\n",
    "\n",
    "y_prob_mc = np.zeros((monte_carlo_simulations, y_valid.shape[0], num_cls))\n",
    "\n",
    "for mc_iter in range(monte_carlo_simulations):\n",
    "    print('MC iteration #{}'.format(mc_iter))\n",
    "    y_prob_mc[mc_iter] = model.predict(x_valid)\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load from the saved array (instead of running the above loop)\n",
    "file_name = os.path.join(base_dir, 'MC_predictions.h5')\n",
    "h5f = h5py.File(file_name, 'r')\n",
    "y_prob_mc = h5f['y_prob'][:]\n",
    "h5f.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### - Compute the predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\mu_{pred} \\approx \\frac{1}{T}\\sum_{t=1}^T p(y^*|\\mathbf{x}^*, \\hat{\\mathbf{w}}_t)$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.display import Math\n",
    "Math(r'\\mu_{pred} \\approx \\frac{1}{T}\\sum_{t=1}^T p(y^*|\\mathbf{x}^*, \\hat{\\mathbf{w}}_t)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "mean_prob = y_prob_mc.mean(axis=0)        # prediction probabilities of shape [2005, 7]\n",
    "y_pred = mean_prob.argmax(axis=1)         # predicted class labels of shape [2005,]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### - Compute the uncertainty estimates (using predictive entropy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle H(\\mathbf{y^*\\mathbf{| \\mathbf{x}^*, \\mathcal{D}}}) = - \\sum_c p(y^*=c|\\mathbf{x}^*, \\mathcal{D})\\,\\text{log} \\, p(y^*=c|\\mathbf{x}^*, \\mathcal{D})$"
      ],
      "text/plain": [
       "<IPython.core.display.Math object>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Math(r'H(\\mathbf{y^*\\mathbf{| \\mathbf{x}^*, \\mathcal{D}}}) = - \\sum_c p(y^*=c|\\mathbf{x}^*, \\mathcal{D})\\,\\text{log} \\, p(y^*=c|\\mathbf{x}^*, \\mathcal{D})')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predictive_entropy(prob):\n",
    "    \"\"\"\n",
    "    Entropy of the probabilities (to measure the epistemic uncertainty)\n",
    "    :param prob: probabilities of shape [batch_size, C]\n",
    "    :return: Entropy of shape [batch_size]\n",
    "    \"\"\"\n",
    "    eps = 1e-5\n",
    "    return -1 * np.sum(np.log(prob+eps) * prob, axis=1)\n",
    "\n",
    "\n",
    "H = predictive_entropy(mean_prob)\n",
    "\n",
    "# Let's normalize the uncertainty values (in range [0, 1])\n",
    "H_norm = (H - H.min()) / (H.max() - H.min())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. Analyze the Results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.1. Compare the prediction performance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy of Standart DenseNet-169 = 88.23%\n",
      "Accuracy of Bayesian DenseNet-169 = 88.23%\n"
     ]
    }
   ],
   "source": [
    "acc_reg = 1 - np.count_nonzero(np.not_equal(y_prob.argmax(axis=1), y_valid.argmax(axis=1))) / mean_prob.shape[0]\n",
    "acc_mc = 1 - np.count_nonzero(np.not_equal(mean_prob.argmax(axis=1), y_valid.argmax(axis=1))) / mean_prob.shape[0]\n",
    "print('Accuracy of Standart DenseNet-169 = {:.02%}'.format(acc_reg))\n",
    "print('Accuracy of Bayesian DenseNet-169 = {:.02%}'.format(acc_mc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.2. Uncertainty of correct/incorrect predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x7f5aaa586cc0>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAD8CAYAAABXe05zAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xd4m/W5//H3Lcm2vHcSJ7Gzd8h0EjaUHaBQyl5ltZSu00IXtNc5Pe35lZ4u2p6uc3IKtAXaQ0vZoVAyIIEMkkAGWWSTxDO2bMlLsq3v74/Hck3iIduSHsm+X9flKx6PpRthf/LN/XyHGGNQSimVOBx2F6CUUqp/NLiVUirBaHArpVSC0eBWSqkEo8GtlFIJRoNbKaUSjAa3UkolGA1upZRKMBrcSimVYFzReNCCggIzfvz4aDy0UkoNSVu2bDlhjCkM59qoBPf48ePZvHlzNB5aKaWGJBE5Eu612ipRSqkEo8GtlFIJRoNbKaUSjAa3UkolGA1upZRKMBrcSimVYMIKbhHJEZFnRGSPiOwWkTOiXZhSSqnuhTuP+xfAq8aY60QkGUiLYk0JyRiDiNhdhlJqGOhzxC0i2cC5wKMAxpiAMaYu2oUlkse3PM6Fj15IY6DR7lKUUsNAOK2SCUA18LiIvCcivxOR9CjXlVCe2/kcR+qO8PiWx+0uRSk1DIQT3C5gAfBbY8x8oBF48OSLROReEdksIpurq6sjXGb8qvBVsLNqJ8nOZJa9s4zaplq7S1JKDXHhBPcx4JgxZmPHx89gBflHGGOWGWNKjTGlhYVh7ZMyJLxx6A0APrv4szS3NfPbjb+1tyCl1JDXZ3AbYyqAoyIyreNTFwK7olpVAll9YDWF6YUsGbuEhaMX8sreV+wuSSk1xIU7j/tLwFMish2YBzwcvZISR0trC28deYt5RfMQESbkTqCioYKGQIPdpSmlhrCwpgMaY7YCpVGuJeFsPLaRlrYW5hfNB2BM9hgADtYcZE7RHDtLU0oNYbpychD2Vu8FYEr+FABGZ44GYF/NPttqUkoNfRrcg1DuKyctKY20ZGs90siMkbgcLg7UHrC5MqXUUKbBPQhlvjLy0/I7P3Y6nBRlFrG/Zr+NVSmlhjoN7kEo95Z/JLjBapdoq0QpFU0a3INQ5i0jLzXvI58bkzWGY/XH8Lf5bapKKTXUaXAPUEtrC54WDwVpBR/5/Ois0QRNkEOeQzZVppQa6jS4B6jcVw5AXtqpI27QmSVKqejR4B6gMl8ZwCkj7qLMIhzi4ECNzixRSkWHBvcAlXutEffJNyeTnEmMSB+hwa2UihoN7gEq95UjCLmpuad8rSC9gOO+4zZUpZQaDjS4B6jMV0a2O5skZ9IpX8tPze8ckSulVKRpcA9Que/UOdwheWl5VDdW09reGuOqlFLDgQb3AB33Hj9lRklIXmoeBkN14/A5UEIpFTsa3ANgjKHcV37KjJKQ0Eg8NGVQKaUiSYN7AOpb6mlubSY/tedWCWhwK6WiQ4N7AEKBnJ/efXCHAr3CVxGzmpRSw4cG9wBUNVYBkOPO6fbrqUmpuF1uHXErpaJCg3sAPM0eADJTMrv9uoiQn5avI26lVFRocA9AZ3Andx/cALmpuZ3L4pVSKpI0uAfA0+zBIY7Ok2+6k5+Wr60SpVRUaHAPgKfZQ0ZyBg7p+eXLS83jROMJXYSjlIo4De4BqGuu67G/HZKflo/BUNVQFaOqlFLDhQb3ANQ215KRnNHrNaGTccobtF2ilIosDe4BCLVKeqOrJ5VS0aLBPQCeZk+frZLQiFunBCqlIs0VzkUichjwAe1AmzGmNJpFxTNjTFgj7rTkNNwutwa3UiriwgruDh8zxpyIWiUJoqm1idZga58jbrDmcodWWSqlVKRoq6Sf+lo12VWOO0dnlSilIi7c4DbAP0Rki4jcG82C4l0ouPtqlYAV3JUNldEuSSk1zIQb3GcbYxYAS4EviMi5J18gIveKyGYR2VxdPXQPEOjXiDs1h+rGaowx0S5LKTWMhBXcxpjjHX9WAc8Bi7u5ZpkxptQYU1pYWBjZKuNIf0bcue5cWtpaaAg0RLsspdQw0mdwi0i6iGSG3gcuAd6PdmHxytPSvxE3oH1upVREhTPiHgm8JSLbgHeA5caYV6NbVvzyNFkbTKUnpfd5bWi/bp1ZopSKpD6nAxpjDgJzY1BLQvA0e0hPTsfh6PvvvNCIWw8NVkpFkk4H7CdPi6fXfbi7ynXnAtoqUUpFlgZ3P4WzajIkNSmVZGeyBrdSKqI0uPvJ0+whIyW84BYRXT2plIo4De5+CmeDqa6y3dka3EqpiNLg7odwN5jqKtedq60SpVREaXD3Q3NrM4H2QL9G3KHVk0opFSka3P1Q11IHhLdqMiTHnUNDoIGmQFO0ylJKDTMa3P3g9XsBSE/ue/FNSOfqSe1zK6UiRIO7H7wtVnCnJaWF/T2hudzaLlFKRYoGdz90jrjDWO4eovuVKKUiTYO7H0LBnZYc/ohb9ytRSkWaBnc/dAZ3P1olGckZuBwuHXErpSJGg7sffC0+oH/BrasnlVKRpsHdD16/F7fLjdPh7Nf35bh1LrdSKnI0uPvB6/f2aypgSLY7m0qfnj2plIoMDe5+8Pl9/WqThOS4c7RVopSKGA3ufvD6vQML7tQcvH4v/jZ/FKpSSg03Gtz94G0ZWHDrIhylVCRpcPeD1+/t1xzukNAinMoG7XMrpQZPg7sfvH5vv1ZNhoQW4eiIWykVCRrcYTLGDPjmZG6qnj2plIocDe4wNbY2EjTBAbVKMlMycYhDZ5YopSJCgztMoZ0BB9IqcYjDmhKoI26lVARocIfJ57eWu6cmpQ7o+3Uut1IqUjS4wzSQQxS6yknVEbdSKjI0uMM0kJ0Bu9L9SpRSkRJ2cIuIU0TeE5GXo1lQvBrI6Tdd5aTmUNtcS6A9EMmylFLDUH9G3F8GdkerkHg32FZJaPXkicYTEatJKTU8hRXcIjIWuAL4XXTLiV+h4B7wzclUXYSjlIqMcEfcPwe+AQR7ukBE7hWRzSKyubp66IWTr8WH2+XG5XAN6PtDi3AqGioiWZZSahjqM7hF5EqgyhizpbfrjDHLjDGlxpjSwsLCiBUYLwa6M2BIXmoeABU+DW6l1OCEM+I+C7hKRA4D/wdcICJPRrWqODTQDaZCMlMycTlcutGUUmrQ+gxuY8xDxpixxpjxwE3AKmPMbVGvLM4MdsTtEAd5qXk64lZKDZrO4w5TfUv9oIIbrD53ua88QhUppYarfgW3MeYNY8yV0Somnvn8vgHtU9JVXmqeBrdSatB0xB2mwfa4AfLS8qhsqMQYE6GqlFLDkQZ3GIwxNPgbBt0qyUvNI9AewNPsiVBlSqnhSIM7DE2tTbSb9sG3StKsKYE6s0QpNRga3GHo3GBqsK0SncutlIoADe4wDHaDqZDQ6snyBr1BqZQaOA3uMIQOURjoBlMhOe4cHOLQEbdSalA0uMMw2L24Q5wOJznuHA1updSgaHCHIVLBDVafWzeaUkoNhgZ3GDoPCh5kqwQgN01XTyqlBkeDOwyD3Yu7q7zUPCp9Oh1QKTVwGtxh8Pq9g9qLu6u81DwaWxs7b3gqpVR/aXCHwef3RWS0DZCflg9AmbcsIo+nlBp+NLjD4G3xDnrVZEhhunXIxDHvsYg8nlJq+NHgDoPX743YiLszuOs1uJVSA6PBHYZ6f31EZpQAZKVkkexM5nj98Yg8nlJq+NHgDoO3ZXCn33QlIhSkF2irRCk1YBrcYfD5fRELboDCtEJtlSilBkyDuw/GGOv0mwi1SsDqc+uIWyk1UBrcfQjtxR3JEXdBegH1LfU6l1spNSAa3H2I5D4lIaGZJce9eoNSKdV/Gtx9iOQ+JSGFaTolUCk1cBrcfQi1MyLdKgFdhKOUGhgN7j5E6tiyrrJSskhxpuhcbqXUgGhw9yEU3JFa8g7WXO7C9EKO1h+N2GMqpYYPDe4+RKNVAla7RHvcSqmB6DO4RcQtIu+IyDYR2Ski341FYfEiGrNKAArSrOA2xkT0cZVSQ184I24/cIExZi4wD7hMRE6Pblnxw9viJcWZgss5+L24uxqVMQpfwIen2RPRx1VKDX19BrexNHR8mNTxNmyGiV6/N6JTAUNGZY4C4JDnUMQfWyk1tIXV4xYRp4hsBaqA140xG6NbVvyI9D4lIUWZRYAGt1Kq/8IKbmNMuzFmHjAWWCwis0++RkTuFZHNIrK5uro60nXaxuuP3M6AXRWmF+IUJ4c9hyP+2Eqpoa1fs0qMMXXAauCybr62zBhTaowpLSwsjFR9tqtvqY/YIQpdOR1ORmaM5GDtwYg/tlJqaAtnVkmhiOR0vJ8KXAzsiXZh8aK+JXKHKJxsVOYoDtVqq0Qp1T/hjLiLgNUish3YhNXjfjm6ZcWPaN2cBGtmyeG6wwRNMCqPr5Qamvqc42aM2Q7Mj0EtcSdoglG7OQnWDcpAe4BybzljssdE5TmUUkOPrpzsRWOgkaAJRm3EHZpZctCjfW6lVPg0uHtR31IPRHafkq4653Jrn1sp1Q8a3L3oDO4ojbhz3Dm4XW6dy62U6hcN7l7U+6Mb3CLC6MzRGtxKqX7R4O6Fr8XaGTBawQ1Wu2Rfzb6oPb5SaujR4O5F54g7Sj1ugJLsEip8FZ1HpCmlVF80uHsR7R43QHFOMQB7T+yN2nMopYYWDe5e1LfU4xAHbpc7as9Rkl0CwJ7qYbMYVSk1SBrcvfC2WKsmRSRqz5GbmktGcoYGt1IqbBrcvaj310e1vw3WzJKSnBINbqVU2DS4exHNDaa6Ks4uZu+JvbpniVIqLBrcvahvqY/aPiVdlWSX0NzarKe+K6XCosHdi5iNuDtmlmi7RCkVDg3uXvj8vpgE99issQjC3mqdEqiU6psGdw+MMdZe3FG+OQmQ4kphVOYodlXtivpzKaUSnwZ3D5pam2gLtsVkxA0wIXcC28q3xeS5lFKJTYO7B6FVk7G4OQkwKW8SVY1VVPgqYvJ8SqnEpcHdA6/f2jskViPuyfmTAXTUrZTqkwZ3D2KxT0lXJTkluBwutpZvjcnzKaUSlwZ3D2Id3MnOZMbljNPgVkr1SYO7B6FtVmMxqyRkUt4kdlTuoD3YHrPnVEolHg3uHkT79JvuTMqfRHNrsx6soJTqlQZ3D7wtXgQhNSk1Zs85Oc+6QantEqVUbzS4e1Dvt5a7OyR2L9HIjJFkpmSy5fiWmD2nUirxaHD3wNvijdkc7hARYXrBdNZ/uB5jTEyfWymVOPoMbhEpFpHVIrJLRHaKyJdjUZjd6lrqyEjOiPnzzhwxk3Jfue4UqJTqUTgj7jbgq8aYmcDpwBdEZGZ0y7Kfp9lDRoo9wQ2w4eiGmD+3Uiox9BncxphyY8y7He/7gN3AmGgXZrfaployUzJj/rxjssaQ7c5mw4ca3Eqp7vWrxy0i44H5wMZoFBNPPM0eMpNjH9wiwozCGdrnVkr1yBXuhSKSAfwN+IoxxtvN1+8F7gUoKSmJWIF28Lf5aWxttGXEDVa7ZMPRDRz2HGZC3gRbaoip+lp45004sBuOHoDaauvzTicUjobiCTB5Fiw4GzLs+X+iVDwJK7hFJAkrtJ8yxjzb3TXGmGXAMoDS0tKEHip6mj0AtgY3wLoP1w3t4D6wG5b/Gd7fDMEgZOdBYREUTwQEgu1QWwXbNsK6FfDUr2F2KVx6HUydbXf1Stmmz+AWEQEeBXYbYx6Jfkn2szu4R2WMYmTGSFYdWMWt8261pYaoqqmEZx6DTW9CWgaUnguzFkL+iO6vNwYqj8OebbDrXdi2AaaeBtfeBZOG/H1ypU4Rzoj7LOB2YIeIhJb0fcsY80r0yrJXbXMtgC09brD63POL5rPq4CqaAk2kJcd2PnlUbVgFT/4S2trgjAth0XmQnNL794jAqLHW21kXw/Z3rND/wQNw9qVw7d2QmR2b+pWKA30GtzHmLUBiUEvcsHvEDbBg9AJe3fcqbx15i0umXGJbHRET8MMTv4T1K2DsBFh6I2Tn9v9xkpJh4dlw2iJYvxLWvQ5b18Od98O8MyJft1JxSFdOdiMegnta4TTSktJYdWCVbTVEjLcOfvwN2LASzrgIbvjMwEK7q+QUOO9yuP1frHbLr74Lf/wF+FsiU7NScUyDuxuh4LZj5WSIy+Fi7qi5rDq4iqAJ2lbHoFUcg4e/DEcPwlW3Wa0OhzNyj19YBLd8wWq5rH0V/vMBOKHHv6mhTYO7G7XNtWQkZ+CMZMAMwPzR86lpquG9svdsrWPAjh+GH30Nmhrhxs/ClCjNBHG5rNH3NXdCVTn8x5dgt+6wqIYuDe5ueJo9trZJQuaPnk+yM5mX9rxkdyn9d/Sg1R4JBq3QLiqO/nNOnA63fRHcafCzb8Hrz1kzUpQaYjS4u+Fp9tjaJglJS0pjwegFvLznZVrbW+0uJ3xlR+An3wSHwwrtnqb5RUNuAdz6BZg0A57+H3jsJ9AaiN3zKxUDGtzdsGufku6cNe4sPM0e1h5ea3cp4TlRAY98y5rCd8O9VpDGWnLKP/vp61fCj78J9Z7Y16FUlGhwd6O2uda2OdwnmzNqDpkpmTy/63m7S+lbvQd++hC0NMO190BOvn21iMOawXLV7fDhfvj+l62eu1JDgAb3SYwxtm3p2h2Xw8XpxaezYv+KzgOM41JTg9VXrquBT94JhaPsrsgydTbcdB8EWuAH91uLd5RKcBrcJ2lqbSLQHoibVgnA+RPOx9/u56/v/9XuUrrnb4Ff/BuUfQhX3w6jx9ld0UeNGmv1vbPz4JffgRXP601LldA0uE/SufgmTlolAONzxzO9cDq/f/f3tAXb7C7no9pa4bf/Dw7uhitugvFT7a6oe5k5cON91t4m//ff8OSvrGX3SiUgDe6TdO5TEkcjboDLplxGmbeMlftX2l3KPwXb4dGfWLv7XfxJmDbH7op6l5wMV98Gi8+HN5fDL/7VavEolWA0uE8SD8vdu7Nw9EIK0wt5bMtjdpdiMQae+o212dO5S2HOYrsrCo84rHovvQ4+2A4PfwUqy+yuSql+0eA+SbwGt8Ph4LIpl7H5+GbePvK2vcUYA8/8zhq1LjrPGsEmmtMWwXWftg5xePjL1paxSiUIDe6ThII7KyXL5kpOdcGkCyhIK+BHa35k7/4lLz4Jr/3N2o3v3KX21TFYxROtfU5S3PDTB+HlP1ntH6XinAb3SWqaanCIg9SkVLtLOUWyM5nrZl/H+5Xv88pem7ZD//tf4KWnrJNoLrzKWmiTyHIL4NYvwvS58Pwf4WfftkbhSsUxDe6TVDdWk+POwSHx+dKcVXIWJdkl/HDND2kIxPjG2orn4W+PWSF3ybVWv3goSHHD5TfBpdfCvp3w75+Dne/aXZVSPRoiv3mRU9lQSY47x+4yeuRwOLh74d1U+Cp4+I2HY/fEa/5uTaObPMs6BMExxH50ROC0xdboO9kNP/82PL1M9/dWcWmI/fYNXmVDJTmp8RvcAFMKpnDFtCt4evvTrD64OvpPuOpF65CCCdPgylus09eHqsJR1g6Dc5bA68/Cd+6D3Qm6ra4asjS4T1LVUEVeap7dZfTp2lnXUpJdwgPLH+Bg7cHoPdHf/wJ/+g1MnmmtinSFc0xpgktKhouvsXY2bG+z9l/5zX9AtR7QoOKDBncX/jY/dS11cT/iBkhyJnH/WfcjItz9t7upaaqJ7BMEg/CX/+3oac+Dj98GrqTIPke8K54Id3zFOpB4xyb4189YW8V66+yuTA1zGtxdVDVUAZDrHuR5iDEyImMED5z1AJUNldz69K1U+CI0ImwNwLIfwD/+BvPPhMtvHNrtkd64kuD0C+Dur1krQ1c8Dw/eYfW/9Yg0ZRMN7i4qGysByE1NjOAGmJw/ma+d8zWOeY9x/Z+vZ9+JfYN7wLoaa//qzWvhvCvggquG3o3IgcjMhqU3wF1ftfY7WfE8PHSX1UL54H3dtErFlP5GdlHpS7zgBpg1YhbfPu/bNAWa+MSTn+CvO/6KGUiQHNhlndd49IB1EMGicxN/nnak5RVam2l95pvWqtFd71nnan7vi9ZhxS3NdleohgEN7i6qGjtaJQkW3AAT8ibw/Yu/z6S8STz42oN89ZWvhj/P2xj4x7Pwo29YQX3L52HqadEtONFl5VirRj/7kLXBVlMD/OHn8LVb4In/gsOD/JePUr0YBlMEwlfZUInL4YqL8yYHIic1hwfPfZAX9rzAszufZWv5Vh654hHmFc3r+Zt8dfDYI7DjHWvmyKXXQ2pa7IpOdEnJMHeJtclW2RHYthHefh3efAVKJlmnzy8+H1LT7a5UDSEa3F1UNVSRm5qLJHB7wOFwcM3Ma5heOJ3fbvwt1//peu5bch9fOuNLJDuTP3rxlrfgiV9CcyNccDXMP0NbIwMlAmPGW28XXGXN/d7+jvX6Pr0MFp8H515uzYXX11gNkvTVCxWRx4ArgSpjzOxwHrS0tNRs3rw5AuXF1i1P34Kn2cN3LviO3aVERFOgiSe2PsGaw2uYUTiDn17+U6YVTgPPCStMNq+BkWPhsuvj56ixocQYqDgG2zdauw+2BmDsBGuWTuk54BimM3VUt0RkizGmNKxrwwjuc4EG4I9DPbgvfPRCRqaP5F/O/Be7S4moLce38OiWRwn4G/ld5pUs2XoAaW+zprktPn/4TvWLJX+LFd7vvg01lTBitHWTc8kFw2NRk+pTf4K7z58YY8waERk/2KISQVVDFdMLpttdRsQtHLOQJb5kil5+jjFNO9mcC6OW3snY0TPsLm34SHF39MIXwb5dsGEVPP4IvPAkXHkznH2JjsBV2CL2V72I3AvcC1BSUhKph40Zn99HU2tTQs4o6U1qVRXFq1eRt3cP/sxMlpcW8K2WtbRv+x4PBj7FreMuS+iefsIRh3Xy/JRZcGgvbFhp7QOz6kW4+XPxf/ybigsRC25jzDJgGVitkkg9bqx0rpocIsGdUudh7BtvULB9G8GkJKrmzadm1iwmOJ38MjCXXx16nu+8/7+8XvEOP5j7BUanFthd8vAiAhOnWzcrP9hhzUL58Tdg4Tlww6chf6TdFao4ps21DqFVk4mwT0lvknxexrz1FiO2bAYRambOomb2bNrd7s5rCpKz+c7UT/Fq9SYe+/BVLn/zK/zb7Hu4Zsz5OvqONRFrlD1xhnV+5ztvWDczr7gJLrtB+9+qW/pT0aFz1WSC7FNystSqKorWr6Ngx3bEGOomT6Z6zlza0rufPywiLB2xmHlZk/j5wWf5+tZf8uThV3lg2s2cVTBXAzzWkpLgzIusk4XeXG6dxrN5Ldz1AIybYnd1Ks70Gdwi8mfgfKBARI4B3zHGPBrtwmLtaP1RBKEgLYFaBsaQeeQIo9e/Te6+fQSdLuqmTKFm5ixaM8M77LjInc/DM+5hRfW7PF22mjs2fo9pmSXcPv5yrh5zLmkud98PoiInKwc+fivMmGfth/L9L1sj74/fYi32UYowpgMORCJOB7x/+f2sO7KOX1z5C7tL6VswSN6e3Yxe9zYZZWW0ud3UTpuOZ9q0j7RE+isQbOWNmm0sr9zIwaZyMl1pXFd8ATeVXMLkzLER/A9QYWlpgjeWw/uboajY2uBq4tCb9aQsEZ3HPRCJGNzXPHkNQRPkofMesruUHjlaAxRu3UrRhvW4PR78mVnUzpxB3aTJmAj2Qo0x7G74kOWVG1jn2UWbaWdR3kxuLrmYy4rOJMU5zPblttuhvdZpPL56a4fCj9+qo+8hSIO7n4wxzP/VfJaMXcJdC++yu5xTuBobGbXpHUZu3kRSUxNNhYXUzJyFr7g46luu1rU2sPLEe7xWtYlyfy0jU/K4b/InubHkIlJOXkKvosffDKtftkbfY8bDPV+39kJRQ4YGdz/VNNWw+DeLuW3ebSydutTucjql1NZStGE9hVvfw9nWhndsMbWzZtE0YkTM97sImiBb6w/wdPkb7PIdYURKLvdN/iQ3j7uEZIeOwGPmwG5rJ8fmRqvvvfRGnXkyRER05eRwcNhzGIBRGfGxX0f68WOMXreOvD27MSLUT5xIzcxZBHLsm6roEAcLcqYwP3sy230H+fPx1Xxv56M8deRVvjf7s5xeENZuCGqwJs2AO79iLdh54QnYugHu+RqMHmd3ZSqGNLj5Z3AXZRbZV4QxZB84wJi315J15AjtycnUzJpF7fQZtKXFzzarIsLcrEnMyZzIprq9LPtwObdu+Dc+MeY8Hpp5BwUpiT0PPiGkpsMVN8PkWbDyeesQh098Ci75pC6bHyY0uIFDnkM4xEFBug1TAU2Q3L17GbN2DRnl5bSmpVFRWkrdlKkEk+K3BSEiLM6dzpysify1/E2eLXuLVZWb+dbMO7mu+AKdBx4L0+ZYuw2ueA6eeRTeW2+djTlytN2VqSjT4AaOeI4wMmMkLkcMXw4TJH/nTsasXUNadTX+zEzKzjiDuomTEmq3PrczmdvHXsz5+fP49eEXeHD7r3mxbC3fP+0+StLjo/U0pKVnwlW3W/t/r3wRvvs5uO4eOP9KPSt0CNObk8AVf7gCt8vN18/5evSfrKMlUrLyddIrK2nJyeHE7NPwjh+f8L9oQRPkterN/P7oawSN4YHpt3DnhCtwSuL8RZTQfPXwj79Z0wenz4U77td91hOIzirpB2MMp/3XaZw3/jxun397VJ8rvayMkhX/IPvwYQIZGVTNn493/IQhdyLKiUA9vzn8Epvq9jAnezI/mPt5pmeNt7us4cEY6xi6N5Zb7195C1x6Lbjit+2mLDqrpB+qGqtobm1mVGb0RiauxgbG/eMfFO7YTltKChWLFuGZOg2TQC2R/ihIzuZfp9zK2todLDuynKvXfp3PTr6Gz0++Frczxe7yhjYRmLPE2nVw9cvw3O9h/Qq47UvWKFwNCcM+uPfX7AeiNKPEGAq2b2P8a6/iCAQ4Mfs0TsyeTTB56C9cERHOzZ/DvKzJPPrhK/x63zM88+Eqvjj1eq4d+zFdvBNtmTlw1W1wcA+sfAF+8k3rxKPr7oGcfLurU4M07IN7a/lWACbkTojo46Z4PEx4+SVyDh2kqbCQsjPOtHUetl2yktK4f9J1XFS4kD8ee51/3fE/PLLnT1xfciHHCxVaAAAOmklEQVSXjFrCnJzJffbAgyaIJ+DDE/ARCLbSGmyj1bThEAdZrnRGuHPJStJT1Ls1cToUT4SNq+GdN+HddbD0erjkWutUHpWQhn2P+56/3cP+2v388NIfRuYBg+2M2riR4tWrAUPVgoV4punJ3mDdT9juPcjLVRt4x7OHIIYMVypTMksYk1pImtNNu2nH29qIt7WRutYGagL1eAJe2k2w18fOSkpnRtYE5uVM4azCuSzKm6ErOk9WVwNr/m4d3JBbAJ+8C5Z8LOFvig8VenMyTEETZOGvFrJwzEI+XfrpQT9eWkUFE196gYzycnxjx1K+5PQe98Me7nxtTbxbv4+dviMcba6mttVLS3sApzhJd7lJd7rJcLrJSkonNymDHFcGWUlpJEsSLoeTJHHSboI0trdQE/BS7q9lf2MZh5rKaTPtpDtTuWL0mVxXfCELcqfpvPKujh2y+t+Vx6y9vj95J8xcoIMLm2lwh2l/zX4uffxS7l10L+dNOG/AjyOtrYxd8yaj171Ne0oKFYsWW9P79Bch5lraA2z3HmSdZydv1+6kJRhgQvpori++kOuKLyA/JdvuEuODCcKu9+Dt18HrgSmz4Zo7YOppdlc2bGlwh+kvO/7CQ689xE+W/mTANyezDh9iwssvkVpbi2fSZCpLSwmm6MyJeNDU7uft2vd5vXoLuxs+JElcLC06g1vHX8bC3Ok6Cgdoa4Mdm2DjKmjwWiPvT9wBE6fZXdmwo9MBw/Tu8XfJTMkc0OZSzuZmSla8zsj33iWQmcmRiy6mcbQuNY4nac4ULi5cyMWFCznSVMmrVZtYUbmJF8vWMjWzhFvGXcrVY84d3jc2XS6YfwbMXmhtWLXpDXj4yzDvdLj6U9aNTRV3hvWI++LHLibHncNXz/5qv74vb/cuxr+ynKSmJmpmzqR67ryIHmSgoqelPcCbNdt5teod9jeVkSQuzh0xjytHn80FIxeR4Uq1u0R7Bfzw7tuwaY21B/jc0+HyG61dCVVU6Yg7DCcaT3Cw9iA3nHZD2N+TXF/P+Ff/Tt7ePTTn5XHs/I/Rkq9zYhOJ25nMpSNKuaRwIfsby1hTu523at9nZeVmUhxJLM6fxTmF8zincB5TMoqHXzslOcWa7z3vdNjyNry3DrZtgGlzrQCfOV/v3cSBYRvcL+x+AYCFoxf2ea20tVG0YT1j1q5BgkEqFyykZuZMnUaVwESEKRljmJIxhruKL2V3w4esq93Je/X7WVttze0fkZLLovyZLMqbyaK8GUzNLMEhw+T/uTsNzroYFp0L2zbClrXws2/BuMlw0TVQeo4en2ajYdkqMcaw9PdLMcbw3Yu+2+u12fv3Mf7Vv5NaW4u3pITK0kW0ZmTEqFJlhyp/He/V72e79wC7Gj7kRKAegExXGgvzprMobyaleTM4LXvy8Dl/s60Ndm6xAry22jqN/rwr4PwrIDvP7uqGBJ1V0oft5du55qlruHvh3Vw46cJur0mtrqJ45UryPtiLPyuLikWLaRwzJsaVKrsZY6gK1LHTd5hdviPs9B3hWEs1ACmOJObkTOkckc/PnUZmUvwcehEVJgiH91ktlIN7wOmCBWdZo/MZ8xNqS+J4oz3uPjzz/jMkO5M5o/iMU77mrjnBmLVrKdi+jWBSEpXzF1A7c+aQ3RBK9U5EGJmSy8iUXC4omA9AfWtjR4gfZmfDEf57/7P8hiAOhBnZE1iUN4PS3Bksyp859E4EEoe1gdWEaeA5YQX4jk2w6U3IyrX642deZB3woKJm2I24K3wVXPLYJcwrmsfnT/+89UljyDh2lKL1661zHp1OaqdNp2b2bNrdup+D6l1zu5+9DUfZ6TvCroYj7G04ij/YCsD49CIW5c1kasey/tGphRSl5pObnDl09ilva7NG37vehYO7IRiEkWNh7hKYs9g6Yk1nXfVJWyU9MMZw5zN3sunYJh6+5GGKJZ38nTsZ8e4W0qqraU9OpnbaNGqnz6A9dZhPC1MD1hps40BTObt8h9npO8Luhg/xtTWdcl2WK43spAyykzPJScog56Q/C1JyKUrNZ5Q7nxHu3MTYe6WpAfZuh/27rKX17W3Wjc5ZC6xVmROmQ/EEvbHZjYgHt4hcBvwCcAK/M8b8Z2/Xx2tw/37L4/zp1R/y1dyPcXqtIevIEcQYmvILqJsyhfoJEzBxfM6jSkzGGHxtzVQH6qgK1HHCX09DezPetiYa2prxtTXT0NaEr72ZhjbrLcipv5f5ydmMcOeS6Uoj3ZXa+ZbhSiXTlWb9mZRGpiudQncuRR2Bb9vIPuCHI/utUfjhD6wTesAafRdPstotRSUwaow1Qs/JH9YztSIa3CLiBD4ALgaOAZuAm40xu3r6nrgIbmOgvhbKjtB2+AM+2PwSuWUVFLVZwdySnYNv3Di848bhz821t1alugiaIE3tfmpbfZwI1FMT8Hb+Wdvqo7k9QEvQT3N7gKZ2P83tflqCgW4fyykORqTkMTZtBCVpIylOG8nYtBEUp41kpDuPnKRMMlyp0Z+vbowV3OUfQvlRqDgKlcehtUvdSclQWAS5+ZBTADl5Vphn5UJahnW+ZugtxT3k5pNH+ubkYmC/MeZgx4P/H3A10GNwD4ox1j+v2tt7/rOtFZqboLnR+qdZc6O1z0JtNXiqrT9rq8Hf0vkfmZHUSnluOm0T59Eyplin9Km45RAHGR0j6ZLUEWF9T7tpp6ndT1O7n4a25s7Qr/bXWyN9fx1vNGyhJuA95Xtd4iQ7KYNkhwunOHGKA6c4cTkcOMSBA+tPpzhwINbnOj5Od7nJcKV1jPjTyEhK/ef7rlRSO3Z6THNZb6njx+GaMBGHCBhwNjXg8NSCp5q2mkraPdW0Vxwk7dAe3M0t1iyWbl8kJ6RnWIGekgruVCvMU9zWxynujs+d9PnkFGvE70rq8ub66OecLmvkLw5wCCDWXxKdb46TPo79XyDhBPcY4GiXj48BS6JSzReusZbZDlRGljW/NCsXxk609hwuLGKj1LAvUMmS/NNo6bhUGyJqKEkC3EBfM6oD7a1U+mupbK6hNuDF17H3ua+tkbZgO+0mSJAg7SZIu7E+NsYQNEGCWH+2myAB00a7CVLtr6OpvYXGtmYa21poN+2D+w9Js95OL5zLU+f8DBp91sCsubFjsNYELR2DtpYmaGm2Ru3+Fmioh0DAatGE3oKDrCdcoQDPzoMfPxn1p4vYrV4RuRe4t+PDBhHZO4CHKQBORKqmGEikehOpVkisehOpVkiAeg9ykD/xHCRArV1Ytf7kqYF+/7hwLwwnuI8DxV0+HtvxuY8wxiwDloX7xN0Rkc3h9njiQSLVm0i1QmLVm0i1QmLVq7V2L5xbuJuAKSIyQUSSgZuAF6NbllJKqZ70OeI2xrSJyBeB17CmAz5mjNkZ9cqUUkp1K6wetzHmFeCVKNcCg2y12CCR6k2kWiGx6k2kWiGx6tVauxGVlZNKKaWiZ/guU1JKqQRla3CLSJ6IvC4i+zr+PGUJo4jME5H1IrJTRLaLyI0xrvEyEdkrIvtF5MFuvp4iIk93fH2jiIyPZX3d1NNXvQ+IyK6O13KliIQ9BSnS+qq1y3XXiogREVtnF4RTr4jc0PH67hSRP8W6xi519PVzUCIiq0XkvY6fhcvtqLOjlsdEpEpE3u/h6yIi/9Xx37JdRBbEusaT6umr3ls76twhIutEZG7EizDG2PYG/Ah4sOP9B4EfdnPNVGBKx/ujgXIgJ0b1OYEDwEQgGdgGzDzpms8D/93x/k3A0za+nuHU+zEgreP9z9lVbzi1dlyXCawBNgClcf7aTgHeA3I7Ph4Rx7UuAz7X8f5M4LCNr+25wALg/R6+fjnwd0CA04GNdtUaZr1ndvkZWBqNeu1ulVwN/KHj/T8Anzj5AmPMB8aYfR3vlwFVQGGM6utc7m+MCQCh5f5ddf1veAa4UOw7qLDPeo0xq40xoa3qNmDNy7dDOK8twH8AP4TORa92CafezwC/NsZ4AIwxVTGuMSScWg2Q1fF+NlAWw/o+Wogxa4DaXi65GvijsWwAckSkKDbVnaqveo0x60I/A0Tpd8zu4B5pjCnveL8CGNnbxSKyGGsEcSDahXXobrn/ycfgdF5jjGkD6gG7ThAOp96u7sEaydihz1o7/klcbIxZHsvCehDOazsVmCoib4vIho5dNe0QTq3/DtwmIsewZox9KTalDUh/f67jSVR+x6K+u7mIrABGdfOlb3f9wBhjRKTHKS4df8M+AdxhTE87z6hwichtQClwnt21dEdEHMAjwJ02l9IfLqx2yflYo6w1InKaMabO1qq6dzPwe2PMT0XkDOAJEZmtv1uRIyIfwwrusyP92FEPbmPMRT19TUQqRaTIGFPeEczd/tNSRLKA5cC3O/6pFCvhLPcPXXNMRFxY/+ysiU15pwhrewIRuQjrL87zjDH+GNV2sr5qzQRmA290dJ5GAS+KyFXGGDv2DA7ntT2G1c9sBQ6JyAdYQb4pNiV2CqfWe4DLAIwx60XEjbXXhl3tnd6E9XMdT0RkDvA7YKkxJuJ5YHer5EXgjo737wBeOPmCjmX2z2H1uJ6JYW0Q3nL/rv8N1wGrTMddCRv0Wa+IzAf+B7jKxh4s9FGrMabeGFNgjBlvjBmP1Su0K7QhvJ+F57FG24hIAVbr5GAsi+wQTq0fAhcCiMgMrM0Fq2NaZfheBD7VMbvkdKC+S4s17ohICfAscLsx5oOoPInNd2fzgZXAPmAFkNfx+VKsk3YAbgNaga1d3ubFsMbLsQ6SOIA14gf4HlaIgPUD/1dgP/AOMNHm17SvelcAlV1eyxfjtdaTrn0DG2eVhPnaClZ7ZxewA7gpjmudCbyNNeNkK3CJjbX+GWu2WCvWv1ruAe4D7uvyuv66479lRxz8HPRV7+8AT5ffsc2RrkFXTiqlVIKxu1WilFKqnzS4lVIqwWhwK6VUgtHgVkqpBKPBrZRSCUaDWymlEowGt1JKJRgNbqWUSjD/HxaEGLrL2SRRAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.kdeplot(H_norm[y_valid.argmax(axis=1) == y_pred], shade=True, color='forestgreen')\n",
    "sns.kdeplot(H_norm[y_valid.argmax(axis=1) != y_pred], shade=True, color='tomato')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### This shows that the model is generally less confident when predictions are incorrect. \n",
    "\n",
    "### Now let's look at some of the predictions and the estimated uncertainty values:\n",
    "(figure from the original paper)\n",
    "\n",
    "<img src=\"imgs/fig5.PNG\" width=\"800\" height=\"100\" >\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3.3. Uncertainty-Aware Skin Lesion Referral\n",
    "\n",
    "\n",
    "### 3.3.1. Accuracy as a function of tolerated uncertainty"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xl8FPX9x/HXJzchECAXN3jggbdGvI9az1pR0bbgUbX+PPCqWm2xeBVLrdajWhFvKtaj2lrFaqWKB2rVEupBAaF4cUMgBAh3ks/vj5nAEkJ2A9lssvt+Ph77yM7szOxnNsl89jvfy9wdERGRxqQlOgAREWn9lCxERCQqJQsREYlKyUJERKJSshARkaiULEREJColCxERiUrJQkREolKyEBGRqDISHUBzKSws9L59+yY6DBGRNmXy5MlL3L0o2nZJkyz69u1LWVlZosMQEWlTzOzbWLbTbSgREYlKyUJERKJSshARkaiULEREJColCxERiUrJQkREokqaprPbav6SaoaPXsycRdX0Kslg5NBiuhem/MciIrKZlC9ZDH9wMd8uqKa2FmYvDBKHiIhsLuWTxZxF1Rufu8OchdWNbC0ikppSPln06pqB2ablWocHnq9g2cqaxAUlItLKpHyyGDm0mN5dM0hLg57FGXzngHa8NLGKc26ezx//XsmqNbWJDlFEJOHM3RMdQ7MoLS315hobavbCDTzxSiUTP1lDfl4aZ5/YkYFHdCAr06LvLCLShpjZZHcvjbqdksXWffHNOh57uZL/zFhHcZd0Lvh+PscOaE96mpKGiCQHJYtmNPmLtTz6UiUzZ6+nb7dMLjw1n0P3aoeZkoaItG1KFs3M3Zn4yRoeH1fJ3MXV9N8hi4tO68Q+/XLi9p4iIvGmZBEnNTXOPz5cxZOvLmfp8hoG7JHDaUfm8fBLlerYJyJtjpJFnK1bX8vf3q3i2fErWLl6U4spM+jdNYMxN3VvsVhERLZVrMki5ZvObqvsrDQGH9eRp0d036yfhnvQE7y6JjmSsIgIKFlst7zcNHrX69jnDmfdNJ+n/rGcihXq3CcibZ+SRTOI7NjXp2sGPzu7Czt0z2TMK8sZPHwevxmzhGlfryNZbvmJSOpRnUUczVm0gZffXcnrH61i9Vpn195ZnHZ0Ht85oL06+IlIq6AK7lZk9dpa3vj3Kl56t4pvF2ygY/s0Tj4sj4FH5lHSRa2mRCRxlCxaIXfnk5nreOmdlfzr8zUAHLp3O47YL5dnXl/OnMVqeisiLUvJopVbVFHNuIlVvPpBFStWqemtiCSGms62ciVdMrjotE48/5seW7SkipxjQ0SkNYhrsjCzE81shpnNMrNhDbzex8wmmNnnZvaOmfWMeO1OM5tqZtPN7H5L0oGYsjJti6a3melG1WoNjS4irUfckoWZpQOjgJOA/sAQM+tfb7O7gLHuvjcwArg93PdQ4DBgb2BP4EDgqHjFmmiRTW8LO6VRXeNcfc8illSqhCEirUM8a1EHALPc/SsAM3sOOBWYFrFNf+Da8PnbwEvhcwdygCzAgExgURxjTajuhZvXUUz+Yi03P1zOVXcv4o4riulVkpnA6ERE4nsbqgcwJ2J5brgu0mfAoPD56UAHMytw9w8JkseC8DHe3afXfwMzu9jMysysrLy8vNlPIFEO2C2He64uZu0656q7F/HFN+sSHZKIpLhEV3BfBxxlZp8Q3GaaB9SY2c7A7kBPggRzjJkdUX9nd3/E3UvdvbSoqKgl4467Xftkc/91JbTLNq69bzFl09ckOiQRSWHxTBbzgF4Ryz3DdRu5+3x3H+Tu+wHDw3WVBKWMj9y9yt2rgH8Ah8Qx1lapZ3Emf7iuK90LM/jlg+VMmLQq0SGJSIqKZ7KYBPQzsx3MLAsYDIyL3MDMCs2sLoYbgCfC57MJShwZZpZJUOrY4jZUKijIT+f315awx47ZjByzlL++tSLRIYlICopbsnD3auAKYDzBhf55d59qZiPMbGC42dHADDObCZQAI8P1fwG+BKYQ1Gt85u6vxCvW1i6vXRp3XFHM4fu0Y9RfKnns5UoNSigiLUo9uNuQmlrnvueW8ff3qzjpkPZce1YX0tOTsvuJiLSQWHtwawCiNiQ9zbhmSGc6d0jjqX+soLKqlpsuLCAnK9HtFEQk2ekq08aYGRec0omrftSZj/67hp//oXyzaV1FROJByaKNOu2oDtx0YSEzvl3HT+9eRLl6e4tIHClZtGFH75/L7ZcXs3hZNVfetYjZCzckOiQRSVJKFm3c/rvmcM/VJaxdV8tPblvAsZfP5oLb5jN/iUoaItJ8lCySwC69s+iQm06tQ63D7IXVDB+9ONFhiUgSUbJIEguWbipJeJgwRESai5JFkuhVkrHFJErP/nOFOu+JSLNQskgSkXNi9O6awcF75fDoS5WMemEZtbVKGCKyfdQpL0nUnxOjttYZ/WIlf31rJRUrahl2XgFZmertLSLbRskiSaWlGZed0YmC/HQe+VsllVU1jLikiLx2KkyKSNPpypHEzIzBx3Vk2HkFTJm1jmvuXcTS5TWJDktE2iAlixRw/EHtGTm0iHnl1Vx510LmLFLnPRFpGiWLFDFgj3bc89NNU7VO/1pTtYpI7JQsUshufYOpWnOzjZ/dt5iPp2qqVhGJjZJFiqmbqrVnSQbDR5cz/qOqRIckIm2AkkUK6pKfzr1Xl7Bvv2zuGFvBM+OXq/OeiDRKySJFtW+Xxu2XF3NMaS6PvbxcnfdEpFHqZ5HCMjOMX55fQOeO6eq8JyKNUrJIcXWd9wrz03n4b5V8+N81bKh2epVkMHJoMd0L9SciIroNJQSd9350XEcK8tNYt96prdUw5yKyOSUL2WjZyk1zedcNc16jegwRQclCIjQ0zPm19y5mXrl6fIukOiUL2ShymPM+XTO49PR8vpq/notGLuTliSvVvFYkhVmyXABKS0u9rKws0WEkncUV1fzuTxVM/mItB+yWw/XndKG4iyq9RZKFmU1299Jo20UtWZjZi2Z2spmpFJKCirtkcOeVRVw9uDNTv17Hhb9ewOsfVqmUIZJiYkkADwJnAf8zs9+a2a5xjklaGTNj4JEdePSXXdmxZxZ3PlXBjQ8toULDnYukjKjJwt3fdPezgf2Bb4A3zexfZnaBmWU2tq+ZnWhmM8xslpkNa+D1PmY2wcw+N7N3zKxnxGu9zeyfZjbdzKaZWd+mnpw0rx5FmdxzdTFDz+hE2fQ1/OTXC3hn8qpEhyUiLSCmW0tmVgCcD/wf8AlwH0HyeKORfdKBUcBJQH9giJn1r7fZXcBYd98bGAHcHvHaWOB37r47MABQo/9WID3N+MF3O/LIDd3oVpjBiMeXctvjS1hepVKGSDKLpc7ib8B7QC5wirsPdPc/u/uVQF4juw4AZrn7V+6+HngOOLXeNv2Bt8Lnb9e9HiaVDHd/A8Ddq9x9dRPOS+KsT7dMHriuhAtOyWfiJ6u58NcL+Nfn+hWJJKtYmrXc7+5vN/RClBr0HsCciOW5wEH1tvkMGERQUjkd6BCWYnYBKs3sRWAH4E1gmLvr62srkp5unHtSPofs2Y7fPrmUGx9awpH7tePr+RuYV16tIUNEkkgst6H6m1mnugUz62xmlzXT+18HHGVmnwBHAfOAGoIkdkT4+oHAjgS3wTZjZhebWZmZlZWXlzdTSNJUO/fK4sFfdOXsEzoy8ZM1zFlUrSFDRJJMLMniInevrFtw92XARTHsNw/oFbHcM1y3kbvPd/dB7r4fMDxcV0lQCvk0vIVVDbxEUEdCvf0fcfdSdy8tKiqKISSJl6xM48JTO5FWrwf4nEXViQtKRJpNLMki3WzTIBBhxXVWDPtNAvqZ2Q5mlgUMBsZFbmBmhRH9N24AnojYt5OZ1WWAY4BpMbynJFivrpsPGdKrWLegRJJBLMnideDPZvZdM/su8Gy4rlFhieAKYDwwHXje3aea2QgzGxhudjQww8xmAiXAyHDfGoJbUBPMbApgwKNNOjNJiLohQ+oSxvEHt09sQCLSLKIO9xF+878E+G646g3gsdZW2azhPlqX2lrn6nsXM3vhBv54czc6dUhPdEgi0oBmG+7D3WvdfbS7nxk+Hm5tiUJan7Q049qzurB6bS2j/7os0eGIyHaKpZ9FPzP7S9iL+qu6R0sEJ21b326ZDD6+I2/8ezWTv1ib6HBEZDvEUmcxBhgNVAPfIehZ/ad4BiXJ4+wTOtKjKIPfP1vBuvW10XcQkVYplmTRzt0nENRvfOvutwInxzcsSRbZWWlcM6QL88qr+dPrKxIdjohso1iSxbqwkvt/ZnaFmZ1O48N8iGxm/91yOG5ALn9+YwVfz1+f6HBEZBvEkix+SjAu1FXAAcA5wHnxDEqSz9AzOpObk8a9zy6jVvN6i7Q5jSaLsAPej8KB/Oa6+wXufoa7f9RC8UmS6NQhnUsGdeK/X67jtX9pWHORtqbRZBE2kT28hWKRJHfiwe3Zp182j/xtmSZOEmljYrkN9YmZjTOzc81sUN0j7pFJ0jEzrjmrC+s2OKPU90KkTYklWeQASwnGZzolfHw/nkFJ8updkslZJ+Tzdtlq/j11TaLDEZEYRR3lzd0vaIlAJHUMOb4jb5Wt4r7nKnj8pm7kZMU0YaOIJFDUZGFmY4Atmq+4+0/iEpEkvaxM49ohXbjm94sZ+9oKLj6tU/SdRCShYvlK93fg1fAxAegIVMUzKEl+++ySw0mHtOf5N1fw5Vz1vRBp7WIZSPCvEY+ngR8CUUcoFInm4tM70SE3jXueqVDfC5FWbltuFvcDips7EEk9+XnpXHZGZ6Z/s55X3lNhVaQ1i2XU2ZVmtqLuAbwC/CL+oUkqOHZALvvvms1jL1eypFJTsIq0VrHchurg7h0jHru4+19bIjhJfmbGNUO6sL7aeeAF9b0Qaa1iKVmcbmb5EcudzOy0+IYlqaRHcSbnnpTPxE/W8OEU9b0QaY1iqbO4xd2X1y24eyVwS/xCklT0o+M60qdbJvf9uYI1azXvhUhrE0uyaGibqP0zRJoiM8O4dkhnFlfU8MdXl0ffQURaVCzJoszM7jGzncLHPcDkeAcmqWevnXP4/uF5/PWtlfxvjvpeiLQmsSSLK4H1wJ+B54C1wOXxDEpS10WndSKvvXHl7xZy7BWzueC2+cxfolZSIokWS2uoVe4+zN1L3f1Ad/+lu2tCAomLDrlpZGUY66uhthZmL6xm+OjFiQ5LJOXF0hrqDTPrFLHc2czGxzcsSWUVKzZVcLvDnEUqWYgkWiy3oQrDFlAAuPsy1INb4qhXSQZmm5Z7Fqs9hUiixZIsas2sd92CmfWhgVFoRZrLyKHF9O66KWHssWN2YgMSkZiawA4H3jezdwEDjgAuiWtUktK6F2Yw5qbuANz3XAUvT6ziqP1yGbBHuwRHJpK6Yqngfh3Yn02toQ4I14nE3aWDOtG3WyZ3jF1KxQrN2y2SKDGNOuvuS9z978B0YKiZTY1lPzM70cxmmNksMxvWwOt9zGyCmX1uZu+YWc96r3c0s7lm9kAs7yfJJzsrjRt/UkDVmlrufGqphjIXSZBYWkN1N7NrzGwSMDXcZ3AM+6UDo4CTgP7AEDPrX2+zu4Cx7r43MAK4vd7rtwETo56FJLUde2Qx9IzO/HvqWl58Z2WiwxFJSVtNFmZ2sZm9DbwDFAAXAgvc/VfuPiWGYw8AZrn7V+6+nuAW1qn1tukPvBU+fzvydTM7ACgB/hnjuUgSO/XIPA7dux2PvlSp3t0iCdBYyeKB8PWz3P1Gd/+cprWC6gHMiVieG66L9BkwKHx+OtDBzArMLA24G7iusTcIE1qZmZWVl5c3ITRpa8yM68/pQsf26Yx8Yglr1mmwQZGW1Fiy6AY8C9wd1jvcBmQ28/tfBxxlZp8ARwHzgBrgMuA1d5/b2M7u/kjYs7y0qKiomUOT1iY/L50bzitgzuJqHvyL5r4QaUlbTRbuvtTdH3L3o4DvApXAIjObbma/ieHY84BeEcs9w3WR7zHf3Qe5+34ETXTrhkA/BLjCzL4hqNf4sZn9tgnnJUlq/91y+NFxHXn1g1W8+5/ViQ5HJGXE2hpqrrvf7e6lBPUKa2PYbRLQz8x2MLMsgkrxcZEbmFlheMsJ4AbgifD9znb33u7el6D0Mdbdt2hNJanpJ6fks2ufLO5+eimLKjQUiEhLiClZRHL3me4+IobtqoErgPEETW6fd/epZjbCzAaGmx0NzDCzmQSV2SObGo+knox048afFFBTC7/541Jq1JxWJO7MPTn+0UpLS72srCzRYUgL+ufHq/jtk0u54Pv5nPu9/Og7iMgWzGxyeNeoUU0uWYi0Fscf1J5jD8zlydeW898v1yU6HJGkFlOyMLMeZnaomR1Z94h3YCKx+OngLpR0TmfkmCVUrVZzWpF4iaUH9x3AB8CNwPXho9H+DyItpX27NIb/pJDyyhrufbaCZLmtKtLaxDLq7GnAru6ucr60Sv13yOaC7+fz+LjlHNg/hxMPyUt0SCJJJ5bbUF/R/J3xRJrV4OM7sm+/bO5/fhlzFm1IdDgiSSeWZLEa+NTMHjaz++se8Q5MpCnS04wbLiggK8MYOWYpG6p1O0qkOcWSLMYRjP76L2ByxEOkVSnqlMF1Z3dh5uz1PD6uMvoOIhKzqHUW7v5k2AN7l3DVDHdXOV9apcP3zeWUI/J4/s2VlO6eQ+numl1PpDlETRZmdjTwJPANwbSqvczsPHfXPBPSKg09oxOTv1jLsFHl4NCrawYjhxbTvTCW9hwi0pBYbkPdDRzv7ke5+5HACcC98Q1LZNvlZKXh7tTWQq3D7IXVDB+9ONFhibRpsSSLTHefUbfg7jNR6yhp5RZVbJqv2x3mLNKAgyLbI5ZkUWZmj5nZ0eHjUUCDMEmr1qskA7NNy106aGQbke0Ry3/QUGAacFX4mBauE2m1Rg4tpnfXDNLSIDvTqKyq5Ytv1K9UZFtp1FlJesurahh6x0I2VMNDw7pSkJ+e6JBEWo3tHnXWzJ4Pf04xs8/rP5ozWJF4ys9L57ZLili1tpabHy5n/Ybk+IIk0pIauw310/Dn94FTGniItBk79czihvMKmP7Neg04KLINGpuDe0H49DJ3/zbyAVzWMuGJNJ8j9s3lvJPzGf/RKv7y1spEhyPSpsRSwX1cA+tOau5ARFrCuSd15Mj92vHwi5VMmrYm0eGItBmN1VkMNbMpwG716iu+Bqa0XIgizSctzfjFuQX07Z7JbY8v0Qi1IjFqrGTxDEHdxMtsXldxgLuf3QKxicRFu5w0fn1pEenpxo0PlVO1RjPsiUTTWJ3Fcnf/BrgPqIior6g2s4NaKkCReOhakMGtFxUyv7yakU8soaZWFd4ijYmlzmI0UBWxXBWuE2nT9umXw5U/7MzHU9fy+Msa0lykMbEMw2ke0c7Q3WvNTMN3SlIYeGQHvpq3gefeWMmOPbI4dkD7RIck0irFNK2qmV1lZpnh46cEU62KJIUrftiZffpl87s/LdWQICJbEUuyuBQ4FJgHzAUOAi6OZ1AiLSkj3bjl/wopyE/npoeXsHR5TfSdRFJM1GTh7ovdfbC7F7t7ibuf5e6aHECSSqcOGhJEpDGN9bP4efjzD2Z2f/1Hy4Uo0jJ26pnFsB9rSBCRhjRWspge/iwDJjfwiMrMTjSzGWY2y8yGNfB6HzObEHb2e8fMeobr9zWzD81savjaj5p0ViLb6Mj9NCSISEO22qrJ3V8Jfz65LQc2s3RgFMFwIXOBSWY2zt2nRWx2FzDW3Z80s2OA24FzgdXAj939f2bWHZhsZuPdXe0bJe7OPakjX89fz8MvVtK3WyYH9m+X6JBEEm6rycLMXgG2Wg5394FRjj0AmOXuX4XHew44lWDypDr9gWvD528DL4XHnhnxPvPNbDFQBChZSNzVDQny9fyFDBtVjgG9umYwcmgx3QvValxSU2O3oe4C7ga+BtYAj4aPKuDLGI7dA5gTsTw3XBfpM2BQ+Px0oIOZFURuYGYDgKwY31OkWbTLSaO2Npi/u9Zh9sJqho9Wuw5JXY3dhnoXwMzurjeL0itm1lxT0l0HPGBm5wMTCZrnbmy3aGbdgKeA89x9iwF8zOxiwma8vXv3bqaQRAILllZvfO5hwhBJVbH0s2hvZjvWLZjZDkAs3VznAb0ilnuG6zZy9/nuPsjd9wOGh+sqw/fpCLwKDHf3jxp6A3d/xN1L3b20qKgohpBEYterJAOzTcvu8PybK9RKSlJSLMniGuCdsLXSuwR1C1fHsN8koJ+Z7WBmWcBgYFzkBmZWaGZ1MdwAPBGuzwL+RlD5/ZfYTkWkeY0cWkzvrhmkpUHvkgwO7J/NQy9Wcs8zFVTXKGFIaolaW+fur5tZP2C3cNUX7h51TAR3rzazK4DxQDrwhLtPNbMRQJm7jwOOBm43Mye4DXV5uPsPgSOBgvAWFcD57v5p7Kcmsn26F2Yw5qbuG5dra50xryzn6fErWLCkmlsuKqJDbizft0TaPotWpDazXIIWS33c/aIwcezq7n9viQBjVVpa6mVlzVWVIrJ1r39YxT3PVNC9MIORlxXRoygz0SGJbDMzm1yvXrpBsXwtGgOsBw4Jl+cBv96O2ETatBMPyeN3VxZTWVXL5XcuYsqstYkOSSTuYkkWO7n7ncAGAHdfDVjju4gkt312yeGB60vo2D6N6+5fzBsfr0p0SCJxFUuyWG9m7Qg76JnZToDGcZaU17M4kweuL2GPHbK5/cmljHmlUi2lJGnFkixuAV4HepnZ08AE4OdxjUqkjejYPp07rizmxEPa89Q/VvDrMUs1Yq0kpUZbQ5mZAV8Q9LI+mOD200/dfUkLxCbSJmRmGNef04VeJZk8+lIli5ZWM+KSIrp0TE90aCLNptGSRTid6mvuvtTdX3X3vytRiGzJzBhyfEduvaiQL+du4PI7F/L1/PWJDkuk2cRyG+o/ZnZg3CMRSQJH7pfLvdcUs6HauequRUyatibRIYk0i1j6WXwB9AO+AVYR3Ipyd9877tE1gfpZSGuyuKKaX44u5+v5G+jcIY3Kqlp6lWjkWml9mrOfxQnAjsAxwCnA98OfIrIVxV0yuP9nJeRkGRUraqmt1ci10rY1Np9FDnApsDMwBXjc3TXspkiMcnPSWBfRMsod5mjkWmmjGitZPAmUEiSKkwjmthCRJthi5Frg/U9XJywekW3VWLLo7+7nuPvDwJnAES0Uk0jSiBy5tkdRBn26ZXLzI0sY/ddlGrlW2pTGato21D0JR5BtgXBEkkv9kWvXb3AeenEZL0xYydSv1nHThYWUdFGFt7R+jZUs9jGzFeFjJbB33XMzW9FSAYokk6xM46ofdeHmCwv4ZsEGLrl9IR/9V81rpfXbarJw93R37xg+Orh7RsTzji0ZpEiyOfqA9oz+RVcKO6XzywfLeezlSmp0W0paMc3cIpIgvUoyGXV9CScf1p5nxq/gZ/ctZkmlWktJ66RkIZJA2Vlp/OzsAoadV8DM2eu5+DcLmfyF5seQ1kfJQqQVOP6g9jz4i67k56Xx8z8s5slXl1NTq9tS0nooWYi0En27ZfLgL7py7IG5PPnqcoY9UM6ylTWJDksEULIQaVXaZacx7LwCrju7C1O+XMfFv1nIZ//TbSlJPDXwFmllzIzvHZbHrn2y+NVjS7j294vJz0tjxSoNRiiJo5KFSCu1U88sRv+iK+2yjcqV4WCEC6q5YZQGI5SWp2Qh0oq1b5fG2vURgxECcxZVc9fTS5k5W5MrSctRWVakletVksHshdW4gxm0zzEm/Hs1r32wit37ZjHwyDyO3j+X7Cx995P40V+XSCsXORhh764ZPHRDN164vQeXn9mJqjW13DG2gh8Nn8/ovy5j3uIN0Q8osg2izpTXVmimPElF7s6nM9fx8sSVfPDZGmpq4YDdcjj1yDwO2asd6ekaAFQaF+tMeboNJdKGmRn77ZrDfrvmsHR5Da99UMUr71dx8yNLKOqUzsmH53HyYXkU5KcnOlRp41SyEEkyNTXOh1PWMO69KsqmryU9DfbfLZvZC6spr6xR81vZTHPOwb09QZxoZjPMbJaZDWvg9T5mNsHMPjezd8ysZ8Rr55nZ/8LHefGMUySZpKcbh++by51XFjP21m4M+k4HJk9fx6KKGs0FLtssbsnCzNKBUQRTsvYHhphZ/3qb3QWMdfe9gRHA7eG+XYBbgIOAAcAtZtY5XrGKJKuexZkMPaMzRE7tqrnAZRvEs2QxAJjl7l+5+3rgOeDUetv0B94Kn78d8foJwBvuXuHuy4A3gBPjGKtIUqs/F7gZfD1f/TQkdvFMFj2AORHLc8N1kT4DBoXPTwc6mFlBjPtiZhebWZmZlZWXlzdb4CLJJrL5bbeCdPJyjavuWsSkaZqlT2KT6Bqu64AHzOx8YCIwD4h5mE13fwR4BIIK7ngEKJIM6s8FvriimuEPlXPDg+Vc+YPOnHpUhwRGJ21BPEsW84BeEcs9w3Ubuft8dx/k7vsBw8N1lbHsKyLbrrhLBvddW8KA/jnc9+dlPPB8hebPkEbFM1lMAvqZ2Q5mlgUMBsZFbmBmhWZWF8MNwBPh8/HA8WbWOazYPj5cJyLNJDcnjdsuLeLMYzrw4jtV3Di6nNVraxMdlrRScUsW7l4NXEFwkZ8OPO/uU81shJkNDDc7GphhZjOBEmBkuG8FcBtBwpkEjAjXiUgzSk8zLjuzM9cM6cyk6Wu56q5FLKpQSynZkjrliQgAZdPX8KvHlpCVafz60iJ275ud6JCkBbSKTnki0naU7t6OB67rSk6mcc29i3nnP6sTHZK0IkoWIrJRn26ZjPp5V3bpncWIx5bwp38sJ1nuPsj2UbIQkc106pDOXVcVc+yAXJ54ZTl3jK1g/QYljFSX6H4WItIKZWUaN5xXQK+STMa8spwFS6oZcUkh+XkavTZVKVmISIPMjHNPyqdXcQa/HVvBxbcvJDPdWFhRrZFrU5BuQ4lIo44+oD33Xl3M0uU1zF9SrZFrU5SShYhEtfsO2RBRbeEeJAx14ksdShYiEpNeXTcfudYdzrppPk+9tpyq1UoayU7JQkRiEjlybZ9uGfzq4gL23CmbMX9fzpAb5/H4y5VUrozYReJIAAAPaElEQVR5HFBpY9SDW0S2y5dz1/P06yt495PVZGcapxyRxw+P7ah5v9uIWHtwK1mISLP4dsEGnhm/nAllq0lPg+8dlsfg4zpS0kUtplozJQsRSYh55Rt47p8rGP/RKtzhhIPbM+T4jvQozkx0aNIAJQsRSahFFdX8+Y0VvPpBFTU1cPBeOXwzv1r9NFoZJQsRaRWWLq/hhQkreGHCSuouNwZ0L8pg7K3dsMgmVtLiNOqsiLQKBfnpXDqoM5EpwYF55dUMHj6fkWOW8Mp7K/lmwQYNWtiKqQwoIi2iV9cMZi+sxh3MoHOHNPbcKZtPZ65jwqRgOPT8vDT22imbvftls/fOOezUM5P0NJU8WgMlCxFpESOHFjN89GLmLNq8zsLdmV9ezeez1m18vP/ZGgByc4w9d8xm752z2btfDh1y07j1sfItjiHxpzoLEWl1ypdtnjy+XbBhi23MoHdJBmNu7p6ACJOHKrhFJGksr6phyqx13PLIEupfsQ7bpx0H7p5Daf92KmVsg1iThT5ZEWn18vPSOXzfXHp3i6j3ANq3M/43Zz0ffLYGWEaPogxK++dQunsO++2SQ26O2vA0FyULEWkzGqr36FaQzpzF1ZRNW0vZ9DWM/3AVL79bRXoa7LFjNgf2z+HA/u3YuWcmaaos32a6DSUiSWX9BmfqV+som76WSdPXMGtOUN+Rn5fGAbvlsHPPTF771yrmL1ElOajOQkQEgGUra5g8fS2Tpq2h7Iu1LFux+XDqnTqk8fNzC9ixeyZFndNTrpOgkoWISD3uznFXzKF2K5e9vHbGDj2y2LFHJjt2z2THHlns0D0zqes+VMEtIlKPmW3RObBncQbXnd2Fr+Zv4Kt5G/hq3nre+HgVq9duyijdCtKDxNEjSCB57YxRLyxjzuJtv5U1f0l1g/1OWvoYsVLJQkRSSiwXWHdnUUUNX85bz9fzNiWRuYurGyyVZKRDcZcM3IN9ax1wqHU2rvPwed3+K1fXEnn5NYPc7KbdAlu9zjeNt2XQu2sGY25qWr8TlSxERBrQvTD6BdXM6FqQQdeCDA7be9P6detr+XZhNUPvWLjZhb66Bnbvm4URXLTrWl2lGVha0Mw3zQwzNj5efrdq8zd1OPHQvCady4tvrdy0u8OcRdVN2r8p4poszOxE4D4gHXjM3X9b7/XewJNAp3CbYe7+mpllAo8B+4cxjnX32+MZq4hINNlZaezSO4ve9W5l9e6awfALCpt0rE9mrt3iGJef2blJxyibvmazY/Qqid8lPW61NmaWDowCTgL6A0PMrH+9zW4Ennf3/YDBwIPh+h8A2e6+F3AAcImZ9Y1XrCIiTRE5H3nvrsGtrLZ6jFjFs2QxAJjl7l8BmNlzwKnAtIhtHOgYPs8H5kesb29mGUA7YD2wIo6xiojELJZbWW3lGLGKZ3uwHsCciOW54bpItwLnmNlc4DXgynD9X4BVwAJgNnCXu1fEMVYREWlEohsPDwH+6O49ge8BT5lZGkGppAboDuwA/MzMdqy/s5ldbGZlZlZWXl7eknGLiKSUeCaLeUCviOWe4bpIFwLPA7j7h0AOUAicBbzu7hvcfTHwAbBF0y53f8TdS929tKioKA6nICIiEN9kMQnoZ2Y7mFkWQQX2uHrbzAa+C2BmuxMki/Jw/THh+vbAwcAXcYxVREQaEbdk4e7VwBXAeGA6QaunqWY2wswGhpv9DLjIzD4DngXO96CX4Cggz8ymEiSdMe7+ebxiFRGRxqkHt4hICku5gQTNrBz4NtFxJEghsCTRQSRQqp8/6DPQ+W/7+fdx96iVvkmTLFKZmZXF8s0gWaX6+YM+A51//M8/0U1nRUSkDVCyEBGRqJQsksMjiQ4gwVL9/EGfgc4/zlRnISIiUalkISIiUSlZtCFmdqKZzTCzWWY2rIHXrzWzaWb2uZlNMLM+iYgzXqKdf8R2Z5iZm1lStY6J5fzN7Ifh38BUM3umpWOMtxj+B3qb2dtm9kn4f/C9RMQZD2b2hJktNrP/buV1M7P7w8/mczPbv1kDCKb706O1Pwgmh/oS2BHIAj4D+tfb5jtAbvh8KPDnRMfdkucfbtcBmAh8BJQmOu4W/v33Az4BOofLxYmOOwGfwSPA0PB5f+CbRMfdjOd/JMGEcP/dyuvfA/5BMDHfwcDHzfn+Klm0HRvnB3H39UDd/CAbufvb7r46XPyIYPDGZBH1/EO3AXcAa1syuBYQy/lfBIxy92UAHgzCmUxi+Qy2NkdOm+fuE4HGpmo4lWBWUXf3j4BOZtatud5fyaLtiGV+kEgXEnzLSBZRzz8sdvdy91dbMrAWEsvvfxdgFzP7wMw+Cqc1TibbM0dOKmjqNaJJ4joHtySGmZ1DMKT7UYmOpaWE86DcA5yf4FASKYPgVtTRBKXKiWa2l7tXJjSqllU3R87dZnYIwRw5e7p7baIDa+tUsmg7YpkfBDM7FhgODHT3dS0UW0uIdv4dgD2Bd8zsG4J7tuOSqJI7lt//XGCcB/PAfA3MJEgeyWJ75shJBTFdI7aVkkXbEXV+EDPbD3iYIFEk2/3qRs/f3Ze7e6G793X3vgR1NgPdPVmGIo5lfpiXCEoVmFkhwW2pr1oyyDjbnjlyUsE44Mdhq6iDgeXuvqC5Dq7bUG2Eu1ebWd38IOnAEx7ODwKUufs44HdAHvCCmQHMdveBWz1oGxLj+SetGM9/PHC8mU0jmJb4endfmriom1eMn8HPgEfN7BqCyu7zPWwq1NaZ2bMEXwYKwzqZW4BMAHd/iKCO5nvALGA1cEGzvn+SfI4iIhJHug0lIiJRKVmIiEhUShYiIhKVkoWIiESlZCEiIlEpWSSpcNTVuyOWrzOzW1s4hj+a2Znh88fMrP92Hq9v/RE3zWwvM/s0fFSY2dfh8zcbOc7OZvbp9sQScay0xkbAbWS//zOz329lvZvZ0RHrzgzXndaE4x9rZi9tyzYNxWZm75vZvrG+/7YIR03OiWG7MWa2a5RtBpnZbs0XnShZJK91wKCwc1aTmVmz9sFx9/9z92nNeczwuFPcfV9335egU9L14fKxzfUeUT6LNKDJySKKKQQdzuoMIRhhNWmZWTpwLUEnuka5+wXuPiPKZoMAJYtmpGSRvKoJhmu+pv4L4Tf0t2zTvBe9w/V/NLOHzOxj4E4zu9XMnjSz98zs2/Db2p1mNsXMXjezzHC/m81skpn918wesbBHYL33fMfMSs1sYERJYIaZfR2+foCZvWtmk81svIWjZYbrPzOzz4DLm/IBhN/67wnjmlJXyqm3TUa4zb/Dz+P/wvXHhjH/neDijZm9EsY3tW474LdAh/B8xobbnRce71Mze9CCcavqvrHPNLN/EwxHsjXvAIeGsXUEegMbS1Rmdnx47Clm9qgFvZkxs5PDz/Q/RIzGamZ54e/23xbM83BKUz7HBj6vSjP7bfh7+dDMisPXuprZy+Hn+JmZHbS1zyPiOL83s8+BnwPFwHsWlgrDv6Wy8PO+OSKG981s363FYmZHEHROuzd8z13NbFLE/ruHvwNpAiWL5DYKONvM8uut/wPwpLvvDTwN3B/xWk/gUHe/NlzeCTgGGAj8CXjb3fcC1gAnh9s84O4HuvueQDvg+1sLyN3HRZQEPgPuCpPOH4Az3f0A4AlgZLjLGOBKd99nG87/B8DuwD7AcQQXj+J621wMLHb3AcCBwOV1yZNgMMbL3H33cPm8ML4DgWvNrDNBqWJleE4/NrM9gdMJPsN9CUZJGGxmPYGbgEOAwwnGsdqaWoKEcWx4rI23iswsl+DzOSP8PeQCF4frHya4SB4AdI843s3A6+E5HgPcbTHc7mlEPvBu+Dv5EPhJuH4U8Eb4d3UAMH1rn0fEcSa6+97ufjuwGDgiolQ4zN1LCX9/1vBtzC1icff3CHozXxP+XmYAa8JYIOjZPGY7zj8lKVkkMXdfAYwFrqr30iFA3SxqTxFcvOq84O41Ecv/cPcNBN+u04HXw/VTgL7h8++Y2cdmNoXgYrRHtNjM7OfAGncfBexKcPF8w4K6hBuBnmbWCegUjuNfF2tTHA486+417r4QeJ8gAUQ6HrggfN+PgU5sGnzvQ3efHbHtNWEJ50OCpLpTA+95LEEyKQuPeVS43cHABHdfGs7F8HyU2J8juKgODp/X2R2Y6e5fhstjCSbF6V+3Phze4ul65zg8jOdtgls9vdm6rQ3rULd+jbvXDX8/mU1/B0cTJCzcvTr8+9va5wGwHvhbI3EMCUtJ/wnPu6FksbVY6nuc4PecQfAl4tlG3lcaoLGhkt/vCf7ZYv0mtare8joAd681sw0R4+zUAhnhN9QHCWalm2NBJXqj31otGBn3BwQXOQhm9prq7ofU265TjDFvDyMoPUxoIMZV9ZaPBA529zVm9j4Nn6cRjFl0U73jbXELLIoPCS68K9z9S9vyzl5TGHBaRIKpi2lrCWMp0Lneui7AkvD5+oj1NWx+HamfaLb2eWQQXOgbTExm1g/4KTDA3SvN7E80/Hk3FkukF4BfAh8QfAlIpWHbm4VKFknO3SsIvsVeGLH6X2y6FXA28N52vEXdP/ASM8sDGr0oWjAv+CjgB+6+Jlw9AyiyYP4BzCzTzPYI/6Erzayu5HN2E2N7j+AWUJqZlQCHAfVHoR0PXBZevAjvb7dr4Fj5QEWYKPYg+LaMu1eH+9VdpN4EfmhhwwIzKwgvyh8Bx5hZl7COodHPKbyIDiO4wEWaTjDy6o7h8jnAu8A0No3IagSV4pHnuHESIAtGJ27Mx8CREXURBxFc9KPNOvc2cGm4T3pY37K1z6MhKwmGmodgtruVwAoL6q9OiPLejR2LcAbJt4AH0C2obaJkkRruZvMx/a8kKJJ/DpxL8A1um4QX9EcJKmDHEwwj3ZjzgQLgpbDy8bXwtsyZwB3hbZ5PgUPD7S8ARoW3MJr69fovwBfA5wQXrWsbGLr9YeB/wKcWNMsdTcPfTl8Fci0Y0fXXBBfUOo8Dn5vZWHefAvwKeDP8fP8JlLj73HC/jwiSWNSWYe7+qru/W2/daoLE/2J4228d8Gi4/lKC2RHLgMihqX8FtLegQnwqwWxyjb3vAoLRW8eHn/tdwJAYRm+9AjghjKsM2G1rn8dW9n8k3O5NgtLwNILf31iCEkFTPAv8Mvwb6xuuexrYAEzY2k6ydRp1VkRSggX9YbLd/VeJjqUtUp2FiCQ9M3uFYBa5YxIdS1ulkoWIiESlOgsREYlKyUJERKJSshARkaiULEREJColCxERiUrJQkREovp/noQKukls7EMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "num_intervals = 20\n",
    "normalized_uncertainty_toleration_removal(y_valid.argmax(axis=1), y_pred, H_norm, num_intervals)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.3.2. Accuracy as a function of fraction of referrals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzs3XeYVOX1wPHvmbZ9l4WlN1l6ExQCUVAQG6JiSTQKdn/B3qImttgNlmjUaIy9BzVqDCKiqGCLhBIpAoKIIlVYYFmWLdPO74+5Mw4ryw6ws7PlfJ7nPnPnzr1zzwzLPfPet4mqYowxxuyOK9UBGGOMqf8sWRhjjKmRJQtjjDE1smRhjDGmRpYsjDHG1MiShTHGmBpZsjDGGFMjSxbGGGNqZMnCGGNMjTypDqC2FBQU6H777ZfqMIwxpkGZN29ekaq2rGm/RpMs9ttvP+bOnZvqMIwxpkERkVWJ7Ge3oYwxxtTIkoUxxpgaWbIwxhhTI0sWxhhjamTJwhhjTI2SlixE5BkR2SgiX1XzuojIwyKyQkQWisiBca+dLSLfOMvZyYrRGGNMYpJZsngOGL2b148BujvLBOAxABFpDtwCDAWGALeISH4S42Txykr+MW0bi1dWNsjjjTEm2ZLWz0JVPxGR/XazywnACxqZ13WWiDQTkbbASGC6qm4BEJHpRJLOpGTE+b9l5fzhr5sIhcHt2saRQzJpmZ/417Jpa5Dps8v2+fhwGLyebfz58lb065a+Nx/FGGOSJpWd8toDq+Oer3G2Vbf9Z0RkApFSCZ06ddqrIL5cVkkoHFkPhWHarDJEEj8+fgrzfT3eH4SrH95I3y5p9Ojso1dnHz06p9G2hRvZkzetY4tXVrJgeQUDeqTTtzAt1eEYY5KgQffgVtUngCcABg8erDXsvku/7JfB6x9uJxBSvG7hz1e02qML3uKVlVzz0MZ9Pz6ouFyReDYVh/jXzO0EgpF9crNc9Ojko2csgfho2az2/ul2dbEPh5VAUAkEwR+MrsdtC0TWV6zx8+RbxYTC4HVv466LWzKoV0atxWaMqR9SmSzWAh3jnndwtq0lcisqfvvMZAXRtzCNP1/Raq9/GSfr+EBQ+W5dgGWrKlm2ys+yH/xMer+EsFMKapHnpkenSPJI8wk/bAjQsbWX1s3dlPuVSr9SUalU+MNURJ9Hl8pw5HlAKd4eYsPmkBPNNnxeCIWIlbb2hD8I1z68ifQ0oUWum+Z57p8e89y0yHXF1pvnusnNciEiVjIxpgEQ1b36QZ7Ym0fqLKaoar9dvHYscCkwhkhl9sOqOsSp4J4HRFtH/Q8YFK3DqM7gwYO1sY8NVeEP8+2aQCR5OEnkhx+DNR6X5hXSfEK6s6T5hPQ0F+k+YdPWIN+v/+k9+nX1sX+3dLwewesBr0fwecR5Lni9ke3RbWs2Bvjra1sJhsDtgmOHZeP1CJtLQmzZFmKzs5RX/vzvzOuB7AwXxaVhVMHngfuvbG0Jw5g6JCLzVHVwTfslrWQhIpOIlBAKRGQNkRZOXgBV/TswlUiiWAGUAec6r20RkTuAOc5b3V5Tomgq0n0u+hamORfTHACem1LMi++WoAougRNHZvPrUbmx5JDmFVyu6us7qt5Gu+Ck/D26WA/skU6Xdr4aSwblFWE2l0QSx5ZtoVgymbu0gq3bI8UYfxDufKaIs47NY+QBmWSkWzcgY+qLpJYs6lJTKFnsyr7WmUTfI1W3geLrbEQit9c2bg2RniaMPDCTYw7Kol/XtHpdwW9MQ5ZoycKSRSPQ0O/5x8ffp4uPxSv9TPuilBnzyiivVDq08jD6l1kc+cusWq3YN8ZYsjCNQHllmE++LOPd/+xg4YpKXAK/6JPOMQdnc1D/DLweK20Ys68sWZhGZe3GANNm7eC9WTsoKg6Rm+XiiCFZHHNQFhV+bdAlK2NSyZKFaZRCYWXe0gqmfbGDzxeWEQiCSKRzo9cNt/62gIP2z0x1mMY0GJYsTKO3rTTEPS9sZtZXFTttL2jmpmdnHz2djow9OvnIy3anKEpj6reUN501Jtnyst2MH53Hl8sqCYQUtwuOPySbbaVhlq3y8/mC8ti+bQs89Ozkiw2j0r2jj6yMSNPcfW0g0NAbGBiTCEsWpkHbXQ/60vIwy3/ws3yVn69XVfL1qkpm/q8MiNy66tjKQ6vmbuYvr4wNBDl+dC4dWnkTPv+ajQFenlZCKAw+T8leNV02piGw21CmSSneHmL5D5HhU5at8jN/WQVlu+hdvrda5bsZcWBmrPNkizy7/WXqN6uzMCYBi1dWcvVDGwkGFY8bfn9WC7p19CV8/IrVfu59YTOBUKQHfee2Xlb/GIgNAtm6uZu+hWn06ZJGv65pFLb34nFbk19Tf1idhTEJ6FuYxv37MBBkp9ZeWjf37HS8P6B8u8bPVysrWfKdn4UrKvlobuT2V5pX6LWfjz5d0uhT6KNvlzTWbgpanYep96xkYUwd2LglyJLvKvlqpZ8l31WyYrWfoDPYrwBKZGDF+y5vxf42+ZWpQ3Ybyph6rNIfqXx/aVoJc5b81PQ3zSuMHJTJoQdkMqhXOj6v3bIyyWW3oYypx9J8Lvp3S+esMcLCbyJNf10C/bv5+GxBGe/N2kFmujC0XwaHDsxkSJ90G4XXpJQlC2NSaFdNfwNB5ctlFXw6v4zPF5QzY24ZPq/wiz7pHDIwk4P7Z5CdaYnD1C27DWVMPRYKKYu+reTT+WV8Or+couIQbhcc2CuSOIYPyLAKcrNPrM7CmEYmHFaWrfLzyZdlfDK/nPVFQQRiNeQeD9x/RSv6dbUKcpM4SxbGNGKqysq1AR59fSvzl1fGtmemCYcckMng3ukM6pVOsxzrFGh2r15UcIvIaOAhwA08pap3V3m9M/AM0BLYApyhqmuc1+4FjgVcwHTgCm0smc2YfSQidO3g4/yxzWIzDbpc0LOzjy8WlfPerB0AdO/oZXDvDAb3Tqdf1zSbA8TstaSVLETEDSwHjgTWEJlT+3RVXRK3zz+BKar6vIiMAs5V1TNF5GDgPuBQZ9fPgOtVdWZ157OShWmqqg5kGAor3/zgZ+7SCuYsrWDJysjYV+lpwsDuabHk0bG1x6arNfWiZDEEWKGqK52AXgFOAJbE7dMH+J2zPgN4y1lXIB3wEbkj6wV+TGKsxjRY0XGootwuodd+afTaL40zjsljR3mY+d9UMHdpBfOWVjDrq61AZByrwb3Tad3CQyik/KJPhlWQm2olM1m0B1bHPV8DDK2yzwLgZCK3qk4CckSkhap+ISIzgPVEksUjqro0ibEa02hlZbgYtn8mw5xJodYXBZm7tJy5Syv4aO4OKvyR/V5+r4T7Lm3JwJ4ZKYzW1Fepbqx9DTBCRL4ERgBrgZCIdAN6Ax2IJJ1RInJI1YNFZIKIzBWRuZs2barLuI1psNoWeDj+kBxum9CScUfnEb0TFQrBTY8XMeWzUkIhqx40O0tmslgLdIx73sHZFqOq61T1ZFU9ALjR2VZMpJQxS1VLVbUUeBc4qOoJVPUJVR2sqoNbtmyZrM9hTKN1QM90fB7B5YqMTdUq380D/9jCuXesZ+a8HYTDljRMRDJvQ80BuotIFyJJ4jRgXPwOIlIAbFHVMHA9kZZRAD8AvxWRiURuQ40AHkxirMY0SVV7kPfp4uM/C8t5evI2bn96M906lvB/Y5vxiz7pVhnexCW1n4WIjCFykXcDz6jqXSJyOzBXVSeLyK+BiUQqtD8BLlHVSqcl1d+ItIZSYJqq/m7XZ4mw1lDG1J5QWPlwThnPTSlmw+YQA7qn8X8nNLMK8EbIOuUZY/ZZIKhM+ayUl6ZtY2tJmIP6Z3D+2DwK2yc+QZSp3yxZGGNqTXllmDdnbOeV6SWUVSiH/yKTc45rRrsCG4u0obNkYYypdSU7QrzyfglvziwlHFaOHZ7NkN7pfLcuYAMZNlCWLIwxSVNUHOTFd0t457NSwhppheLzCn++opUljAYm0WSR6n4WxpgGqKCZh6tOb86vRuUAkVYolQHlPwvLUhuYSRpLFsaYvXboAZmkeSXWsW/yp6XMWlSe2qBMUthtKGPMPokOZNimhYdJ00v4dk2AX43K4bcnNLM5xBuA+jCQoDGmCYgfyHD4wEwe/9dW3vhoOwu/qeCP5xfQoZU3xRGa2mC3oYwxtcbnFS47tTl3XFDAj1tCXDBxA9P/uyPVYZlaYMnCGFPrhg3I5Mkb2tCto4+Jz2/m7uc3U14RTnVYZh9YsjDGJEXLfA8PXNGKs8bk8sHsHVxw9waW/+BPdVhmL1myMMYkjdstnHNcM+6/ohUVfuWyP2/gjY9KaCwNa5oSSxbGmKQb0COdJ29ow+DeGTz6ejE3PraJbaWhVIdl9oAlC2NMncjLdnPnhQVceko+876u4P/u2sD85RWpDsskyJrOGmPqjIhw8mE59O+Wxh1PF3H1Qxs5emgm7Vp6OaCnjS1Vn1nJwhhT57p39PH4dW0Y0iedabPKeObtbVz1lx9ZtMJKGvWVJQtjTEpkpLvo3zUtNlRIMAQ3P1HE9P/usDnA6yFLFsaYlBnQ46c5wD1uyEoXJj6/mbNvX8/Uz0sJBC1p1BfJnlZ1NPAQkWlVn1LVu6u83pnIvNstgS3AGaq6xnmtE/AU0JHIoJZjVPX76s5lY0MZ0zBFx5Ya0COd3vv5+GJROS++W8LyH/y0au7m9CNzOebgbBtnKklSPp+FM4/2cuBIYA0wBzhdVZfE7fNPYIqqPi8io4BzVfVM57WZwF2qOl1EsoGwqlY7/rElC2MaD1VlzpIKXnx3G4tX+mmR5+bUI3I4bng2GWl2Q6Q21dp8FiLypogcKyJ7+i80BFihqitV1Q+8ApxQZZ8+wEfO+ozo6yLSB/Co6nQAVS3dXaIwxjQuIsKQvhk8fHVr7r+iFR1be3jsjWLG/3Ed/3hvGzvKbeiQupZIAvgbMA74RkTuFpGeCb53e2B13PM1zrZ4C4CTnfWTgBwRaQH0AIqdRPWliNznlFSMMU2IiHBAz3QeuLI1D1/dmh6dfDz1722M++M6nn9nG9vLLGnUlRr7WajqB8AHIpIHnO6srwaeBF5S1cA+nP8a4BEROQf4BFgLhJy4DgEOAH4AXgXOAZ6OP1hEJgATADp16rQPYRhj6rt+XdO4+9JWLFtVyUvvlvD8O9v454clDBuQQevmHob2zbB+GkmUUJ2F82v/DOBMYB3wMjAc6K+qI6s55iDgVlU92nl+PYCqTqxm/2zga1XtICK/BO5R1RHOa2cCv1TVS6qL0eosjGlavl3j52+vb+XL5ZUApNkc4HulNuss/gV8CmQCx6vqWFV9VVUvA7J3c+gcoLuIdBERH3AaMLnKexfE1YVcT6RlVPTYZiLS0nk+CliCMcY4unbwMahXOtE2Uv6gssCGD0maROosHlbVPqo6UVXXx7+wu2ykqkHgUuA9YCnwmqouFpHbRWSss9tIYJmILAdaA3c5x4aI3KL6UEQWAULktpcxxsQM6JEea1KrCr27+FIcUeNV420oEbkEeFlVi53n+USawP6tDuJLmN2GMqZpWryykmlflPLO5zs47/g8zjgmL9UhNSi1dhsK+G00UQCo6lbgt/sSnDHG1Ja+hWlcPb4FhwzM4B/vl7Blmw19ngyJJAu3iMS6TjpNWK2sZ4ypVyac2IxgUHlmSnHNO5s9lkiymAa8KiKHi8jhwCRnmzHG1BvtW3k5cUQO0/6zg2/X2PSttS2RZPEHIr2rL3KWD4HfJzMoY4zZG2eOySM708Xf3yy2qVtrWY3JQlXDqvqYqv7aWR53WisZY0y9kpPp4qwxecz7uoL/LrZmtLUpkX4W3UXkdRFZIiIro0tdBGeMMXtq7KHZdGjl4e9vbCVo82LUmkRuQz0LPAYEgcOAF4CXkhmUMcbsLY9buODkZvzwY5B3PitNdTiNRiLJIkNVPyTSJ2OVqt4KHJvcsIwxZu8d3D+Dgd3TeO6dbZTaYIO1IpFkUekMyfGNiFwqIiex+2E+jDEmpUSEC3+VT8mOMC9P25bqcBqFRJLFFUTGhbocGERkQMGzkxmUMcbsqx6dfBw1NIs3Z25nfVEw1eE0eLtNFk4HvN84kw+tUdVzVfVXqjqrjuIzxpi9dt7YPNwu4cm3rKPevtptsnCayA6vo1iMMaZWtWzm4TdH5jLzf2V89W1lqsNp0BK5DfWliEwWkTNF5OTokvTIjDGmFpx6RA4t8tz87Y2thMPWlHZvJZIs0oHNROaUON5ZjktmUMYYU1sy0lycPzaPr7/3M2NeWarDabASmVb13LoIxBhjkiVa0f3kv4sZPiCDNF8iv5NNvBqThYg8C/ys7Kaq5yUlImOMqWUul3DRyflc/dBG3vhoO+NG25wXeyqR9DoFeMdZPgRyAesWaYxpUA7omc6w/Z05L0pseLs9lchAgm/ELS8DpwI1zqpkjDH1zYSTmlHpV56fYh319tTe3LjrDrRKZEcRGS0iy0RkhYhct4vXO4vIhyKyUERmikiHKq/nisgaEXlkL+I0xpiddGzt5YQRObzzeSnfrbM5L/ZEIqPObheRkugCvE1kjouajnMDjwLHAH2A00WkT5Xd/gy8oKr7A7cDE6u8fgfwSc0fwxhjEnPWmFwy04W/v2kd9fZEIrehclQ1N27poapvJPDeQ4AVqrpSVf3AK8AJVfbpA3zkrM+If11EBgGtgfcT+SDGGJOI3Cw3Z47JY86SCmYvLk91OA1GIiWLk0QkL+55MxE5MYH3bg+sjnu+xtkWbwEQ7eB3EpAjIi2cgQvvB66pIbYJIjJXROZu2rQpgZCMMQZOHJFD+5YeHnuzmJDNeZGQROosblHVWG2QqhYDt9TS+a8BRojIl8AIYC0QAi4Gpqrqmt0drKpPqOpgVR3csmXLWgrJGNPYeT3ChJOasWp9gHc+t8adiaixnwW7TiiJHLcW6Bj3vIOzLUZV1+GULEQkG/iVqhaLyEHAISJyMZHh0H0iUqqqP6skN8aYvTF8QAb7d0vjqX8Xs6UkxC/6ZNC3MC3VYdVbiZQs5orIAyLS1VkeAOYlcNwcoLuIdBERH3AaMDl+BxEpcG45AVwPPAOgquNVtZOq7kek9PGCJQpjTG0SEUYflEVpufLC1BKueWgji1faYIPVSSRZXAb4gVeJVFJXAJfUdJCqBoFLgfeApcBrqrpYRG4XkbHObiOBZSKynEhl9l17/AmMMWYvbdn2U+e8QFBZsLwihdHUb4mMDbUD2Ktf9ao6FZhaZdvNceuvA6/X8B7PAc/tzfmNMWZ3BvRIx+ctwR9QwhqZMMnsWiKtoaaLSLO45/ki8l5ywzLGmOTrW5jG/Ve0Yuwh2YjAtFk7ULXWUbuSSEV1gdMCCgBV3SoiCfXgNsaY+q5vYRp9C9No0czNs29vY0jfDI4ampXqsOqdROoswiLSKfpERDqzi1FojTGmIRt3dC79u6bx8KtbWGdzdv9MIsniRuAzEXlRRF4iMvzGDckNyxhj6pbbJVx/TgtE4E/PFllnvSoSGe5jGnAgP7WGGuRsM8aYRqVNCw9Xnd6cJd/5efFdG5k2XkKjzqpqkapOIdIE9iIRWZzcsIwxJjVGDc7iyCGZvPRuCV99a/0uohJpDdVORK4SkTnAYueY05IemTHGpMjlv2lO6+Zu7nq2iNLycKrDqReqTRbOIH0zgJlAC+B8YL2q3qaqi+ooPmOMqXNZGS5uOLeATcUhHn5lS6rDqRd2V7J4xHl9nKrepKoLsVZQxpgmom9hGmeNyeODOWV8MHtHqsNJud0li7bAJOB+Z7a7OwBv3YRljDGpN/7oXPp1TePBV7awvok3p602WajqZlX9u6qOAA4HioEfRWSpiPypziI0xpgUcbuFG85pgQB/eq5pN6dNtDXUGlW9X1UHE5nNzkbbMsY0CW1aeLjy9OYsXunnpWklqQ4nZRJKFvFUdbmq3p6MYIwxpj46/BdZHDEkkxenbmuyzWn3OFkYY0xTdPlvmtOquZuJzxWxowk2p7VkYYwxCcjOcHHDOQX8uDXEQ682vea0CSULEWkvIgeLyKHRJdmBGWNMfdOvaxpnHpPHB7PL+HBO02pOW+MQ5SJyD/AbYAkQnVZKiQwoaIwxTcoZo3OZu7ScBydtoW9hGm1aJDLTQ8OXSMniRKCnqo5R1eOdZWyNRwEiMtrpo7FCRH42256IdBaRD0VkoYjMFJEOzvaBIvKFiCx2XvvNnn0sY4xJjkhz2gIU+NNzm5tMc9pEksVK9qIznoi4gUeBY4A+wOki0qfKbn8GXlDV/YHbgYnO9jLgLFXtC4wGHoyfrc8YY1KpbYGHK37TnK++reTBV7bwj2nbWLyycbeSSqT8VAbMF5EPgdi3oaqX13DcEGCFqq4EEJFXiPTRWBK3Tx/gd876DOAt572Xx51nnYhsBFoS6RhojDEpd+TQLKbPLuWdz3cgAj5PCX++ohV9C9NSHVpSJFKymAzcAfwHmBe31KQ9sDru+RpnW7wFwMnO+klAjoi0iN9BRIYAPuDbBM5pjDF1pneXSGJQhUBQWbC88fZXrrFkoarPi4gP6OFsWqaqgVo6/zXAIyJyDpEK87X8VImOiLQFXgTOVtWfNWwWkQnABIBOnTpVfdkYY5JqSJ8MXnm/hEAw0uqnTyMtVUBi81mMBL4hUv/wN2B5gk1n1wId4553cLbFqOo6VT1ZVQ8gMn0rqlrsnDcXeAe4UVVn7eoEqvqEqg5W1cEtW7ZMICRjjKk9fQvTeODK1owanIkqvP/fHag2zgrvROos7geOUtVlACLSg8hotINqOG4O0F1EuhBJEqcB4+J3EJECYItTargeeMbZ7gP+RaTy+/XEP44xxtStvoVp9C1Mo0OrYl6YWkKXdl5OOTw31WHVukTqLLzRRAGxyucaW0epahC4FHiPyHSsr6nqYhG5XUSiTW9HAstEZDnQGrjL2X4qcChwjojMd5aBiX4oY4ypa2eNyePQAzJ4/M1i/ru4PNXh1DqpqcgkIs8AYeAlZ9N4wK2q5yU5tj0yePBgnTt3bqrDMMY0YeWVYS6//0c2FAV59Pdt6NSm/k8BJCLznBHFdyuRksVFRJq7Xu4sS5xtxhhj4mSkubjjgpb4vMKNj22iZEeo5oMaiBqThapWquoDTkX0yar6F1Vt3L1PjDFmL7Vp4eG2CS35cUuQO55uPD28q00WIvKa87jIGXJjp6XuQjTGmIalX9c0rhrXnHlfV/C3N7amOpxasbvWUFc4j8fVRSDGGNOYHHNQNt+vC/DPD7fTpZ2P44ZnpzqkfbK7ObjXO6sXq+qq+AW4uG7CM8aYhmvCSc0Y0iedh17Z0uB7dydSwX3kLrYdU9uBGGNMY+N2CTedV0C7lh5uebKI9UXBVIe013ZXZ3GRiCwCelWpr/gOWFR3ISZfOBymoqJhZ31jTP2Unenirotaogo3/X0TZRUNc0rW3ZUs/gEcD/zbeYwug1R1fB3EVmc2btzIW2+9xbZt21IdijGmEerQysvN/1fAqg0B/vTcZsLhhtdCand1FttU9XvgISJDckTrK4IiMrSuAqwL4XCYtWvXMnnyZLZv357qcIwxjdCgXulc/Kt8/rOwnGfebng/TBOps3gMKI17Xupsa1Sys7MJh8O8/fbb7NjRtObWNcbUjZNGZnPssCz+8V5Jg5vDO5FkIRo3Jogz6F+jnHS2efPmVFRUMGXKFMrKylIdjjGmkRERLv9Nc/bvlsa9L25m6fcNp39zQtOqisjlIuJ1liuITLXaKLVs2ZLS0lKmTp1qld7GmFrn9Qi3TSigRZ6bmx8vYlNxw2ghlUiyuBA4mMgw42uAoTgTDjVWrVq1YuvWrUybNo3KyoaT+Y0xDUNetps7L2xJWUWYax/ayAvvFNf7ObwTGRtqo6qepqqtVLW1qo5T1Y11EVwqtWnTho0bN/L+++8TCNTWxIDGGBNR2N7HGcfk8sOPQZ57p4RrHtpYrxNGtXUPIvJ7Vb1XRP5KZMbAnajq5UmNrB5o27Yt69atY/r06Rx55JF4vfV/uGFjTMMRP1l0dA7vvvV0atbdVVQvdR6b9CQR7dq1Y9WqVcyYMYPDDz8ct9ud6pCMMY3EgB7p+Lwl+AMae15fVZssVPVt5/H5ugunfmrfvj3ffvstbrebww47DJcrkaoeY4zZvb6Fadx/RSuef6eYuUvr7y0o2P1wH2+LyOTqlkTeXERGi8gyEVkhItft4vXOIvKhM4zITBHpEPfa2SLyjbOcvXcfr3aICO3bt+frr7/ms88+IxxumN31jTH1T9/CNG77bUvyc1w8/e9iapq9NFV29xP5z8D9wHdAOfCks5QC39b0xiLiBh4lMuhgH+B0Eemzi3O8oKr7A7cDE51jmwO3EGl5NQS4RUTyE/9YtU9E6NChA4sWLeI///lPvf0HNcY0PBnpLsaPzmP+N5X8b1n9LGHsbriPj1X1Y2CYqv5GVd92lnHAIQm89xBghaquVFU/8ApwQpV9+gAfOesz4l4/GpiuqltUdSswHRid+MdKDpfLRYcOHZg/fz6zZ8+2hGGMqTXHDc+mVXM3T9XT0kUiN9+zRKQw+kREugBZCRzXHlgd93yNsy3eAuBkZ/0kIEdEWiR4bEpEE8acOXP43//+l+pwjDGNhM8rnH1sHstW+flsQXmqw/mZRJLFVcBMp07hYyIlgCtr6fzXACNE5EtgBJGOfwnPcC4iE0RkrojM3bRpUy2FVDO3202HDh344osvmD9/fp2d1xjTuB01JIuOrT088/Y2QvVsZNpEOuVNA7oTmWb1cqCnqr6XwHuvBTrGPe/gbIt/73WqerKqHgDc6GwrTuRYZ98nVHWwqg5u2bJlAiHVHrfbTfv27fnss8/46quv6vTcxpjGye0Wzj2+GavWB/hwTv0an67GZCEimcC1wKWqugDoJCKJzMs9B+guIl1ExAecBuzUikpECkQkGsP1wDPO+nvAUSKS71RsH+Vsq1c8Hg/t2rVj5syZzJ49m5UrV7Ju3Tq2bNlCWVmZtZoyxuyxQwdm0L2jl+enFBMI1p/SRSKjxz4LzAMOcp6vBf4JTNndQaoaFJFLiVzk3cAzqrpYRG4H5qrqZGAkMFFEFPgEuMQ5dosVVKN6AAAgAElEQVSI3EEk4QDcrqpb9uiT1RGv10u7du2YP39+LDmISKyCKisri5ycHHJzc8nLyyM3N5f09HTS09PJyMggPT0dj6dRDuJrjNkLLpdw/thmXPfoJqZ+XsoJI3JSHRKQWLLoqqq/EZHTAVS1TEQkkTdX1anA1Crbbo5bfx14vZpjn+Gnkka95vV6adu27c+2qyqBQIDy8nJKSkpYuXIlwWBkhMn4hOLz+cjMzNzpuOpaQ1S3XURIS0sjPT2dtLQ0MjIyYo8+nw+Px4PX693lo8fjsY6GxtQjv+iTTv9uabz47jaOPiiLdF/q/38mkiz8IpKBMz6UiHQF6mdD4HpGRPD5fPh8vt3uFwwGCQaDJJiDd3r/KFUlGAyybds2gsEg4XCYYDBIKBTaKcFUPSb6mJOTQ5cuXejYsSMtW7YkIyNjj2IxxtQeEeH/xuZxxQMbeWtmKacdlZvqkBJKFrcA04COIvIyMAw4J5lBNTXRX/epVFlZybJly1i0aBGqSosWLejWrRtt27aloKDABlE0po7175bO0L7pTHq/hOOGZ5OdmdrSxW6vUM7tpq+J9IX4JSDAFapaVAexmTqUlpZGtEWZqlJeXs68efMIBoO4XC7atWtH165dadWqFc2bN7fbVsbUgfPGNuOCiRt47cMSzju+WUpj2W2yUFUVkamq2h94p45iMikmImRmZsbqUcLhMCUlJXzyySeoKj6fj86dO7PffvvRqlUrcnJy9vgWmjGmZt07+hh5YCavf7Sdk0bmkJ+TulGvE7n38T8R+YWqzql5V9MYuVwucnNzyc2N3DcNBoOsW7eOFStWoKpkZWXRrl074KfK+Wh9SDgcjj2PthaLrsfvJyJkZWWRlZVFZmYm2dnZsfqetLS02LrX67XEZJqUc47P45P5ZfzjvRIu+XXqhshLJFkMBc4Qke+BHURuRakz+J9pgjweD/n5+eTnR/5wKysrWbduHRC56Ecv5vEX9eq2RZ+rKjt27CAQCMQq5uP3jyYVl8tFRkZGLKFkZmaSlZVFRkYG4XCYUCgUq9SPrkcTU7TiP7oeTVqhUAiPx0OrVq1o0aJFLDHa3CWmPujU2svoX2Yx+ZPt/HpUDq2bp6Z+M5GzHp30KEyDlpaWRlpa3czuFW3lFQwGKSoqIhAIEAgEdurjUt0CxOpaqr4WDodZvXr1TkkqPz+f1q1b06ZNm1gfmfgmzsbUlTPH5DF99g5enLqNa85okZIYdjetajpwIdANWAQ8rarBugrMmF1xuVwJNUfeG3l5ebH1aCX/d999x9KlS2PbfT4frVu3pm3btrFSSE5OTspbs5nGrXVzD8cfksNbH2/n1CNz6dS67lsn7u4v/HkgAHzKT3NSXFEXQRmTalUr+aMCgQDFxcWsW7cuVgqBSKLJzs4mOzubnJwcsrOzYx0j4xe7tWX21vjRuUz9TynPT9nGH88vqPPz7y5Z9HFaQSEiTwOz6yYkY+ovr9eL1+uNVfZDpBRSWVnJ9u3b2bJlS6zeRVVjt73ie+tnZ2fHhoGJJpb4hJKeno7P57PmyWYn+Tlufj0qh5feLeH0o/x061j7pevd2V2yCERXnHGe6iAcYxoeEYmN91WTYDAYK51s2rSJQCCwy8p8VY1V4scvWVlZsSFdoovP57MWYk3EqYfn8u+PS3l6cjETL2lVp+feXbIYICIlzroAGc7zaGuo1Pc/N6aBifbWr2k4lWiLrUAgwKZNm1i3bt3PKvLjmx1nZmbGEkfVUkq0pBI9d3SJHxss/vZY/PAw8UPCVI0v/jE6tI3dZkuu7EwXpx2Vy5NvFbNoRQX9u9X8A6W2VJssVNX+1Y1JERGJ3fKqSbSFWDgcpry8nB07dsSaDcc3JVbVnVqGRdXmFJ5erzfWXyb+dlvVuptoIwUrEe25k0Zm88ZHJTw1eRsPXpVWZ9+hNeEwfP311yxatIj+/fvTq1evVIdT5xr654+2EKsPQqEQgUCA7du3s3Xr1p3qb6oOYhnfiCDaCTM6UnJ6enosWVZXGmqqoyWn+1yceUweD726lblLK/hFn7oZ9NOSRT2wrxerPT2+rKwMn8/H6tWrKSkp4bbbbiMcDuN2u7nlllvo27dvnd5OqOvPX/XYm266iWAwiMfj4c4779yr92jIyaY2ud1u3G53QvU38bfaioqKCIVCsQ6Z0c6TUdWNluzxeGK327xe7896/EdLMlVvuVVdd7vdeL3eBnMbbcywbF79oISn/l3MoF7puFzJL11YskiB6H+GhQsXEggEuPvuu2MtZy6//HJKS0s54IADWLZsGX379uX777+nW7durF27lk6dOrF582batGnD9u3b2bhxI3fccQfBYBC3282NN95IeXk5vXv3ZubMmRx88MH861//4thjj+XZZ59l3LhxPPfcc5x33nm88847ZGZmxubYCIfDzJ49m6effprf/e53TJo0ibPPPpsZM2Zw9NFHs3TpUgYOHMjmzZtp27Zt7MJQ9WIZCARwu91s2LCB1q1b89VXX9GnTx8++eQThg8fzr/+9S9OOOEEnnrqKYYPH86tt94auxjcfPPNzJ49mzPPPJPJkydz4oknMmPGDA4//HDmzJnDkCFDWL58OT179mT9+vVs376dW265Jfb5b7/9dsrKyjjwwAN5//33OfLII5k0aRKnnXYaf/3rX7n00ku59dZbueWWW7j00ksZNWoUfr8fAL/fz/z583nyySeZOHEiDzzwAFdffTUvvPACZ511FlOnTuXYY49l9uzZDB06lG+//ZZgMMjNN98cSzZ33HEHPXr0wOVyJXx7oKkmmz251Vad6P+lUChEZWVlbIbK+O27Szyw8204l8uF1+vd6XZZfN1PdD1ayokmmejz6JLspOP1COccm8fdL2zh0/nljDgw+Z1FLVkkWTgc5rvvviM7O5tZs2bRtWtX3n33XUaMGMGGDRvYsGHDTmMnrV69mh49euz0y2zjxo106NCBJUuWkJ+fz8cff8xhhx3GlClT8Pl8sYtdOBxm7ty5tG7dmoqKCtq2bYvP5+PII4+kRYsWXHXVVWRlZXHXXXcBcPnll/P1118zbdq02MVu+PDhnHfeeQCMGzeO5s2b0759e9LT0/nuu+/o168fL774IhdeeCE33HAD55xzDvfcc0+swvWGG27gySef5J577uHhhx/mD3/4Ax988AFdunRh9erVBAIBcnJyUFWGDBnCsmXLYrGLCN988w29e/fG5/ORlZWFy+Viy5bIJInz589n0KBBTJ48mSuvvJJHHnmEAQMG7PT558+fz7p16+jbty8bNmwgGAzSrFkzVJWDDjoIEeHcc8/F4/Fw7733sn79et54443Y5x8wYABDhw7F6/UyatQoPB4PrVu3xuPxsGXLFlwuF59//jlDhw7lqaeeYvDgwTslm0WLFnHDDTfw+uuvc+qpp/Laa69xxhln8NJLL3HeeefxzDPPcOGFF/L3v/+dyy+/nIsvvpgbbrgh9v396U9/4uOPP+aCCy7g9ddf55RTTuGjjz5i1KhRLFiwgIEDB7J69Wo6duxIWVkZmZmZsaHlU1EyS/Xx0R8stTWEfvwwMcFgkMrKSoqLi3dKPvH9a0SElStXsnz5cnr06EFhYSGqGitdxc+KGb9ExzmLNjyILxEl+iPj8CFZTJq+nTsf+oBHZR4nHDeK8b8eUSvfw65IbVZupdLgwYN17ty5e3XsunXrePvtt2OD4e2p6B97hw4dKCws5JtvviEzM5MFCxbQunVrwuEwHTp0IDMzk/z8fFq0aBG717qvt0FSeRslFArx5ptv8tJLL8VKRuPGjeOUU07Zo1/VDfXz13T+QCCAx+OhvLycjIwMiouLyc/PZ8OGDbRp04YffviB//73v7HvT0QYP348ubm5HHXUUbzyyiucdtppPProo1x88cXcfvvt3HzzzVxxxRU89NBDnHHGGdx0001cf/31uFwuwuEwEydO5C9/+QuPP/44l156KY888ghXXnklDz74INdeey333Xcf119/PRMnTuTmm29m3LhxP0tWzz//PBMnToztHz3+kksu4ZFHHuGCCy7g8ccf59xzz+W666772fn/9re/8fDDD3Pddddx9913x+K+//77ufrqq3nssce46KKLeOGFFxgyZAg33XRT7PPfdddd/Pjjj4wYMYJ58+YxaNAgli5dSp8+fVi1ahWdO3emqKiIgoKC2Pf69ddf89VXXyUlWcU3DAgEAni93th5S0tLWbNmDX/84x93Kll27NiRjIwMduzYgcfjobKycqckVF2y6d69O926dYt17MzJyYl19owv3URb01VUVHDP3z/ljmtPIBwO4HJ7eWHStD1OGCIyT1UH17RfUmuHRGS0iCwTkRUict0uXu8kIjNE5EsRWSgiY5ztXhF5XkQWichSEbk+mXEmKhgMUlZWxpo1a/j222+ZN28er732Gtdffz0vvvgi9913H9OnT6egoID8/HzGjRvHMcccw3HHHcfAgQPp0aMHLVu23KlSrlevXtx5552MHz9+ry50+3p89D1OOeWUvfpV179/f7xeLy6XC4/Hw/77779HrTMa8uev6fzREXIzMzNjY00BtGnTBoBOnTrt9P15vV72339/Ro8eHUu8LpeLyy67DLfbzW233Ybb7eaRRx7B7XYzadIklixZErtQu1wulixZwsSJE3G5XPzxj3/E5XJx7bXX4nK5uOSSSxARzj//fESE008/ncWLFwPEmuQuWbKE8847L7Z//PG33HILbrebe++9F7fbzaOPPvqz8y9evDiWfCZMmADAySefDMDw4cNj3xlA8+bN+eqrr3b61b5o0SIWLlxIKBRi+vTphEIhJk2aRDAY5KGHHiIQCHDjjTcSCASYMGECixYt4g9/+AMvv/wyv//971m4cCGnn346lZWVjB8/fqfHM844g8rKSs4880z8fj9nnXUWixYt4ve///0ujz/llFPw+/2ceuqp+P1+zj77bPx+PxMmTCAQCHDZZZexYMEC/H4/4XAYv9/PggULuOiiiwgGg1x00UW43W4uvPBCsrKyuOiii8jPz+d3v/sdBQUFXHvttRQVFXHvvfcyefJk7rvvPjZu3MiECRMoKipizJgxfP755/Tr14+33nqLDh068Pzzz9O8eXMeeeQRCgoKmD71FcKhCtAQ4VCAye98tMd/w4lKWslCRNzAcuBIYA0wBzhdVZfE7fME8KWqPiYifYCpqrqfiIwDxqrqaSKSCSwBRqrq99Wdb19KFpMnT+bZZ5/l4IMPpmvXrpSUlFBSUoLf72fTpk2ICKtXryY3N5eioiJ69+5NSUkJhYWFTJs2jffeey/2n2X8+PGccsopexVHQ9VU77nXllRW0Df04//5z3/y8ssv7/T/b+zYsfh8PiorK0lLS4uVBKKP0dt3O3bs4N133632+Gijj9qKP3qtDYVCuN1u/H4/kydP/tn5jz32WIBY/UtlZWWsDix+JID09HTmffUjD953A+FQ8ksWyayzGAKsUNWVTkCvACcQufBHKRDt3JcHrIvbniUiHiAD8AMlJMEXX3wR+wUxZcoUBgwYwNixY/n+++8ZOHAgGRkZdOvWjZ49e1JQ8PPxWA477DA+/PDD2B9L//79kxFmvdarVy9LEvtgX76/aMlmb5NNQz++f//+eDyenf7/RUdAjtb7Rcf3ij5mZWUBkJ2dvdvjE6mkjo+/X79+dO/ePVbSqFrJXvVRRGjTpk3sPG63mzZt2lBSUkJeXh7t2rUjLy+PZs2akZmZuVOdR/SYccCBBx7I5Hc+YuyxDbTOQkR+DYxW1f9znp8JDFXVS+P2aQu8D+QDWcARqjpPRLzAi8DhQCZwlao+sYtzTAAmAHTq1GnQqlWr9jjOiRMnctNNN+1TycB+WRuTOrVdwR6tX4gOhR8dkiX6vOrcLPH1GtGK67S0tFiz3uhzr9cba+Ibfd3tdrNgwQJmzZrFYYcdxogRI+q8z0x9KFkk4nTgOVW9X0QOAl4UkX5ESiUhoB2RRPKpiHwQLaVEOQnkCYjchtqbAEaOHBlrUbS3JQP7ZW1MzaJNWOMfdzVjYvzjrrZVfa2wsJDCwkLC4XDs1k38L/iq61WPz83NZdiwYQCsXbsWINY5MDpzY/x6/KyN8et7O0x9586dGTt27F4dW5eSmSzWAh3jnndwtsU7HxgNoKpfOHNoFBApXU1T1QCwUUQ+BwYDK6llBx10EK+++mqszsIu+qauxLeOib+IiUi1kzRV3VZV1Slrq7swxz+v+hhd311Dhd29tqu7FdHmpFU7wkUXj8cTq2SPn3636hIfY/zUvfDTiMDxF/D4X/XRbdWdP36fptgzvCbJTBZzgO4i0oVIkjiNSBKI9wORW03PiUhvIB3Y5GwfRaSkkQX8EngwWYEOHjyY9evX73XTWdP0xLfFrzoOU3ScJvj51LFR0d7H0Y5e0VsP0XGcou8R/VUcPwVs/EVzVxftXV0E4y/S0Qtm1Yt3fKu26FAa8Us0ie1qqfqa2+3+2T6mYUtasnCGNb8UeA9wA8+o6mIRuR2Yq6qTgauBJ0XkKiKV2ueoqorIo8CzIrKYyCi3z6rqwmTFapqmcDgcm5Y1fu7vqu3go+KHmYgfSiIzM3OnDljx7eHjh5mI7+lbG+MaVf217Xa7qy1xGLOvklpnoapTgalVtt0ct74EGLaL40qBptX+1NSa+DGH/H5/LBnE/xKPdiKMjoqanZ29U0/b+KEbqg7lUF8GsIsmhvoQi2n8Ul3BbZqQ6EW8srKSysrKWBPDff0lXDUJAGRkZJCTk0N+fj55eXnk5ubGRjONH+fHfoUbkxhLFqbWRHux+v3+WDKAnSfqycjIIC8vjzZt2pCfnx+7iEfbtsdXXMY/1rQtfra69PR0+7VtTC2zZGF+pmrlbdXH6MW/6v18l8tFbm4uLVq0IC8vj/z8/Go7ExljGhZLFk1UKBSirKyM0tLSWMudqPgRM6OT00TX40fMjF+iFb52W8eYxsmSRRNRUVHBjh07qKioACLj9rdp04bevXvTqlUrcnJydmqHbowx8SxZpJiqsmXLFioqKmJt76s2sdzT1jehUIjS0tJYb1ZVJTc3l8LCQtq3bx+r9LWkYIxJlCWLFAmHw2zevJmKigq6dOlCz549CQQCsVtDZWVl7Nixg7KysthFPypaZxCtJ4hO/FJZWYmI4PF4aNeuHQMGDIgNl57INJfGGFMdSxZ1LBQKsWnTJoLBID179mTAgAG0aNGixuOifQbiWxr5/X7Ky8vZvn07IkLbtm3Jz88nJyfHWgMZY2qVJYs6Ep2UXlXp168f/fr1Iy8vL+HjoxXJ0eGVjTGmLlmySLLKykqKiopwu90MGjSIXr162QXfGNPgWLJIkrKyMrZu3YrP5+Pggw+mR48eVm9gjGmwLFnUstLSUrZu3Upubi6HHXYYhYWFsQpoY4xpqCxZ1JJt27ZRUlJC8+bNGT16NJ07d7amqcaYRsOSRS3YsGEDzZs3Z+TIkbRr185aIhljGh1LFvto27ZtZGVlceyxx1qdhDGm0bKfwPvA7/ezfft2jj76aEsUxphGzZLFXgqHw2zYsIFRo0Yl1KnOGGMasqQmCxEZLSLLRGSFiFy3i9c7icgMEflSRBaKyJi41/YXkS9EZLGILBKRevXTfcOGDfTt25cePXqkOhRjjEm6pNVZiIgbeBQ4ElgDzBGRyc5UqlE3Aa+p6mMi0ofIFKz7iYgHeAk4U1UXiEgLIJCsWPfUli1byM/PZ9iwYTYktzGmSUhmyWIIsEJVV6qqH3gFOKHKPgrkOut5wDpn/ShgoaouAFDVzaoaSmKsCauoqKCyspKjjjrK+k8YY5qMZCaL9sDquOdrnG3xbgXOEJE1REoVlznbewAqIu+JyP9E5PdJjDNhoVCIH3/8kSOOOIJmzZqlOhxjjKkzqa7gPh14TlU7AGOAF0XEReT22HBgvPN4kogcXvVgEZkgInNFZO6mTZuSHuz69es58MADKSwsTPq5jDGmPklmslgLdIx73sHZFu984DUAVf0CSAcKiJRCPlHVIlUtI1LqOLDqCVT1CVUdrKqDW7ZsmYSP8JOioiLatGnDkCFDknoeY4ypj5KZLOYA3UWki4j4gNOAyVX2+QE4HEBEehNJFpuA94D+IpLpVHaPAJaQImVlZYRCIY444gg8HuvHaIxpepKWLFQ1CFxK5MK/lEirp8UicruIjHV2uxr4rYgsACYB52jEVuABIglnPvA/VX0nWbHuTigUoqioiKOPPpqcnJxUhGCMMSmX1J/JqjqVyC2k+G03x60vAYZVc+xLRJrPptTatWsZOnQoHTt2rHlnY4xppFJdwV2v/fjjj3Tq1IkDD/xZdYkxxjQpliyqUVpaisfjYdSoUTaKrDGmybOr4C4EAgG2bNnC6NGjbQpUY4zBksXPqCrr169n2LBhtGnTJtXhGGNMvWDJoooNGzbQvXt39t9//1SHYowx9YYlizglJSVkZmZy6KGHWj2FMcbEsSuiw+/3U1JSYhMZGWPMLliycASDQUaOHEmyhw0xxpiGyJIFkJGRwbBhw+jVq1eqQzHGmHrJBjoC8vPzbYBAY4zZDStZGGOMqZElC2OMMTWyZGGMMaZGliyMMcbUyJKFMcaYGlmyMMYYUyNLFsYYY2pkycIYY0yNRFVTHUOtEJFNwKpUx5FCBUBRqoNIIfv89vnt8++dzqpa4zhHjSZZNHUiMldVB6c6jlSxz2+f3z5/cj+/3YYyxhhTI0sWxhhjamTJovF4ItUBpJh9/qbNPn+SWZ2FMcaYGlnJwhhjTI0sWTQgIjJaRJaJyAoRuW4Xr/9ORJaIyEIR+VBEOqcizmSq6TuI2+9XIqIi0qhayCTy+UXkVOfvYLGI/KOuY0ymBP4PdBKRGSLypfP/YEwq4kwGEXlGRDaKyFfVvC4i8rDz3SwUkQNrNQBVtaUBLIAb+BYoBHzAAqBPlX0OAzKd9YuAV1Mdd11/B85+OcAnwCxgcKrjruO/ge7Al0C+87xVquOu48//BHCRs94H+D7Vcdfi5z8UOBD4qprXxwDvAgL8EvhvbZ7fShYNxxBghaquVFU/8ApwQvwOqjpDVcucp7OADnUcY7LV+B047gDuASrqMrg6kMjn/y3wqKpuBVDVjXUcYzIl8vkVyHXW84B1dRhfUqnqJ8CW3exyAvCCRswCmolI29o6vyWLhqM9sDru+RpnW3XOJ/IrozGp8Ttwit4dVfWdugysjiTyN9AD6CEin4vILBEZXWfRJV8in/9W4AwRWQNMBS6rm9DqhT29RuwRm4O7ERKRM4DBwIhUx1KXRMQFPACck+JQUslD5FbUSCIly09EpL+qFqc0qrpzOvCcqt4vIgcBL4pIP1UNpzqwhs5KFg3HWqBj3PMOzradiMgRwI3AWFWtrKPY6kpN30EO0A+YKSLfE7lvO7kRVXIn8jewBpisqgFV/Q5YTiR5NAaJfP7zgdcAVPULIJ3IuElNQULXiL1lyaLhmAN0F5EuIuIDTgMmx+8gIgcAjxNJFI3pXnXUbr8DVd2mqgWqup+q7kek3masqs5NTbi1rsa/AeAtIqUKRKSAyG2plXUZZBIl8vl/AA4HEJHeRJLFpjqNMnUmA2c5raJ+CWxT1fW19eZ2G6qBUNWgiFwKvEekVcgzqrpYRG4H5qrqZOA+IBv4p4gA/KCqY1MWdC1L8DtotBL8/O8BR4nIEiAEXKuqm1MXde1J8PNfDTwpIlcRqew+R52mQg2diEwi8kOgwKmTuQXwAqjq34nU0YwBVgBlwLm1ev5G8j0aY4xJIrsNZYwxpkaWLIwxxtTIkoUxxpgaWbIwxhhTI0sWxhhjamTJwtQZEQmJyPy4Zb9aeM8TRaRP3PPbnY6JSSMik5xRPa+qsv1WEVnrfLYlInJ6Au+1U/y72e9CETlrX+KOe6/nROTX1Wz/TkQWiMhyEXlBRGocX0xErhSRzNqIzdRf1s/C1KVyVR1Y3Ysi4lHV4B6+54nAFGAJgKrevA/x1UhE2gC/UNVu1ezyF1X9s4h0B+aJyOuqGtjNW+4Uf3WcdvR14VpVfV0iHXWuBD5yhsvw7+aYK4GXiLTtN42UlSxMSonIOSIyWUQ+Aj4UkWxnLo7/icgiETkhbt+znF/0C0TkRRE5GBgL3Of8mu8a/6tZRA535jVY5MwFkOZs/15Ebos7R69dxJUuIs86r38pIoc5L70PtHfOd0h1n0tVvyFy8cx33q+riEwTkXki8qmI9Kom/t+KyBznM74R/cXulFqucdZnisg9IjLbKQEc4mx3i8h9zvELReQCZ7uIyCMSmQfiA6BVTf8uzsilfwE2AMc47/OYiMyVyDwZtznbLgfaATNEZEZ1+5lGINVjtNvSdBYiPYrnO8u/nG3nEBnPqLnz3APkOusFRHqjCtCXyDhHBc5r0f2fA34dd47ngF8TGeZhNdDD2f4CcKWz/j1wmbN+MfDULmK9mkgPYYBeRIaRSAf2o/r5BG4FrnHWDwQ+jXvtQ6C7sz4U+Kia+FvErd8ZF2f8e88E7nfWxwAfOOsTgJuc9TRgLtAFOBmYTqTXczugOP6cVb+7KtseBP5Q5Tt3OzHsH/d9FsQds8v9bGnYi92GMnWputtQ01U1Ok6/AH8SkUOBMJEhllsDo4B/qmoRQNz+1ekJfKeqy53nzwOXELn4AbzpPM4jcjGtajjwV+dcX4vIKiLjLJXUcN6rRORcZ9/jAUQkGziYn4ZhgcjFfFf6icidQDMiQ7e8V81+8fHv56wfBewfVx+RR2QQwUOBSaoaAtY5pbhESdz6qSIygUhCb0tkcqGFuzgm0f1MA2LJwtQHO+LWxwMtgUGqGpDI6LHpSThndETeELX7/yBaZzEWeFpEuhK53VtcTaKs6jngRFVdICLn4AwKuAu7il+IlER2SjCyb1OLHkDk9mAX4Boi9b4tKhIAAAGUSURBVDVbReQ5dvHvkuh+puGxOgtT3+QBG51EcRgQnUf8I+AUEWkBICLNne3biQxNXtUyYD8RiVZEnwl8vAdxfEokcSEiPYBOznsmRCOD2s0FzlbVEuA7ETnFeT8RkQHVxJ8DrBcRb/T8e+A94CLnWESkh4hkEZli9jf/394dqkQURHEY/474DmL2FdYH8DkEi90gGLcIJrtBTIKCGKzbZIOIRVh0YdGnMAhi8BhmLoLIHW5z4fvVO3c4k/4MM8ypZxrrlPa7vWqNe5SdwYTSfe4deIuINeo5xh9r6BunJWZY6L+5AEYR8QTsAAuAzJwDR8A0ImaUJkdQWmse1EPojW6SzPygvLp5Xef6AobcKDoBVuq/V5TXS4f2BzkE9qM0ZdoGdmvtc37agf6ufww8AHfd2gc4o9yqeoyIZ8pz9avADfBav50D9z1zHNcaX4BNYCszPzNzRuntvQAua32dU2ASEbeNcVpivjorSWpyZyFJajIsJElNhoUkqcmwkCQ1GRaSpCbDQpLUZFhIkpoMC0lS0zc4ncdkSH8wUQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "reps_for_random = 40\n",
    "num_fractions = 20\n",
    "uncertainty_fraction_removal(y_valid.argmax(axis=1), y_pred, H_norm, num_fractions, reps_for_random)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}