BMI260 2019 Problem Set 3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# BIOMEDIN 260/RAD260: Problem Set 3 - Mammogram Project\n",
    "\n",
    "## Spring 2020"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Name: Yiheng Li"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Deadline = May 29 11:00 PM**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "Breast cancer has the highest incidence and second highest mortality rate for women in the US.\n",
    "\n",
    "Your task is to utilize machine learning to study mammograms in any way you want (e.g. classification, segmentation) as long as you justify why it is useful to do whatever it is you want to do. Turning in a deep dream assignment using mammograms might be amusing, for example, but not so useful to patients. That being said, choose something that interests you. As the adage goes, \"do what you love, and you’ll never have to work another day in your life, at least in BMI 260.\"\n",
    "\n",
    "Treat this as a mini-project. We highly encourage working with 1 other person, possibly someone in your main project team. \n",
    "\n",
    "In addition to the mammograms themselves, the dataset includes \"ground-truth\" segmentations and `mass_case_description_train_set.csv`, which contains metadata information about mass shapes, mass margins, assessment numbers, pathology diagnoses, and subtlety in the data. Take some time to research what all of these different fields mean and how you might utilize them in your work. You dont need to use all of what is provided to you.\n",
    "\n",
    "Some ideas:\n",
    "\n",
    "1. Use the ROI’s or segmentations to extract features, and then train a classifier based on those features using the algorithms presented to you in the machine learning lectures (doesn't need to use deep learning).\n",
    "\n",
    "2. Use convolutional neural networks. Feel free to use any of the code we went over in class or use your own (custom code, sklearn, keras, Tensorflow etc.). If you dont want to place helper functions and classes into this notebook, place them in a `.py` file in the same folder called `helperfunctions.py` and import them into this notebook.\n",
    "\n",
    "## Data\n",
    "\n",
    "The data is here:\n",
    "\n",
    "https://wiki.cancerimagingarchive.net/display/Public/CBIS-DDSM\n",
    "\n",
    "## Grading and Submission\n",
    "\n",
    "This assignment has 3 components: code, figures (outputs/analyses of your code), and a write-up detailing your mini-project. You will be graded on these categories.\n",
    "\n",
    "If you're OK with Python or R, please place all three parts into this notebook/.Rmd file that we have provided where indicated. We have written template sections for you to follow for simplicity/completeness. When you're done, save as a `.pdf` (please knit to `.pdf` if you are using `.Rmd`, or knit to `.html` and use a browser's \"Print\" function to convert to `.pdf`).\n",
    "\n",
    "If you don't like Python OR R, we will allow you to use a different language, but please turn your assignment in with: 1) a folder with all your code, 2) a folder with all your figures, and 3) a `.tex`/`.doc`/`.pdf` file with a write-up."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Using TensorFlow backend.\n"
     ]
    }
   ],
   "source": [
    "# Loading packages\n",
    "from tqdm import tqdm\n",
    "from keras_tqdm import TQDMCallback\n",
    "import h5py\n",
    "import os\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib.image import imsave\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from sklearn.preprocessing import MinMaxScaler\n",
    "from sklearn.metrics import confusion_matrix, accuracy_score, roc_auc_score\n",
    "from tensorflow.keras.models import load_model\n",
    "from skimage.transform import resize\n",
    "from sklearn.calibration import calibration_curve"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mammogram Classification and Uncertainty Analysis to Optimize Reveiw Workflow"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Util functions\n",
    "def getImageList(rel_path, verbose=False):\n",
    "    \"\"\"\n",
    "    Input relative path, get all the images from 'rel_path\\\\patient ID\\\\image names'\n",
    "    input: relative path\n",
    "    output: df(patient number, image file name list, image file dir list)\n",
    "    \"\"\"\n",
    "    patient_num_list = []\n",
    "    img_name_list = []\n",
    "    img_path_list = []\n",
    "    # Get all the image paths and names\n",
    "    cur_path = os.path.abspath('')\n",
    "    data_path = os.path.join(cur_path, rel_path)\n",
    "    p_folder_list = os.listdir(data_path)\n",
    "    for p_folder in p_folder_list:\n",
    "        img_dir = os.path.join(data_path, p_folder)\n",
    "        img_list = os.listdir(img_dir)\n",
    "        for img_name in img_list:\n",
    "            img_path = os.path.join(img_dir, img_name)\n",
    "            patient_num_list += [p_folder]\n",
    "            img_name_list += [img_name]\n",
    "            img_path_list += [img_path]\n",
    "    assert len(patient_num_list) == len(img_path_list), 'Lengths differ.'\n",
    "    df = pd.DataFrame(list(zip(patient_num_list, img_name_list, img_path_list)),\n",
    "                columns = ['patient number', 'img name', 'img path'])\n",
    "    return df\n",
    "\n",
    "def loadOneImage(img_path):\n",
    "    hf = h5py.File(img_path, 'r')\n",
    "    img = hf['data'][()]\n",
    "    return img\n",
    "\n",
    "# intergrated function to load images and do \n",
    "# preprocessing to the images loaded\n",
    "def loadImages(img_path_list, size = None):\n",
    "    assert len(img_path_list)>0, 'No files to load'\n",
    "    # image example\n",
    "    example = loadOneImage(img_path_list[0])\n",
    "    img_shape = list(example.shape)\n",
    "    if size is not None: \n",
    "        ar_shape = tuple([len(img_path_list)]+list(size)+[3])\n",
    "        imgs = np.zeros(ar_shape)\n",
    "        for i in tqdm(range(len(img_path_list))):\n",
    "            imgs[i] = makeRgb(resize(image = loadOneImage(img_path_list[i])[:,:,0], output_shape = size))\n",
    "    else: \n",
    "        ar_shape = tuple([len(img_path_list)] + img_shape[0:2]+[3])\n",
    "        imgs = np.zeros(ar_shape)\n",
    "        for i in tqdm(range(len(img_path_list))):\n",
    "            imgs[i] = makeRgb(loadOneImage(img_path_list[i])[:,:,0])\n",
    "    return imgs\n",
    "\n",
    "# Preprocess Image\n",
    "# Scale to 0 - 255\n",
    "def imgScale(img):\n",
    "    img -= np.min(img)\n",
    "    img *= 255/np.max(img)\n",
    "    return img.astype('int')\n",
    "\n",
    "# \"Make\" RGB\n",
    "def makeRgb(img):\n",
    "    \"\"\" input: 2D img; output: img w/ 3 channels identical \"\"\"\n",
    "    new_img = np.repeat(np.expand_dims(img, axis=-1), repeats=3, axis=-1)\n",
    "    return new_img\n",
    "\n",
    "def plotReliability(probs, y_test, n_bins = 10):\n",
    "    # reliability diagram\n",
    "    fop, mpv = calibration_curve(y_test, probs, n_bins=n_bins)\n",
    "    # plot perfectly calibrated\n",
    "    plt.plot([0, 1], [0, 1], linestyle='--')\n",
    "    # plot model reliability\n",
    "    plt.plot(mpv, fop, marker='.')\n",
    "    plt.show()\n",
    "    \n",
    "def getFeatrues(scaled_imgs):\n",
    "    features = np.zeros((len(scaled_imgs), 1024))\n",
    "    for i in tqdm(range(len(scaled_imgs))):\n",
    "        feature = getOutput(scaled_imgs[i:i+1])\n",
    "        features[i] = feature\n",
    "    return features"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**1. Describe what you are doing and why it matters to patients using at least one citation.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Traditional mammogram classification requires extra annotations such as bounding box of detection or like in our dataset, mask ground truch for segmentation, which requires a lot of work by the radiologists. And plus, classification of the lesions after having ROI or segmentation is also difficult which requires high expertise and time[<sup>1</sup>](#s1). Nowadays, machine learning appoaches, particularly deep learning, are widly utilized to improve computer-aided diagnosis (CAD) and decision report[<sup>2</sup>](#s2). Most popular methods seek to improve detection performances, but they didn't necessarily provide reliability of their prediction, which stil requires the interaction by an expert raiologist[<sup>3</sup>](#s3).\n",
    "\n",
    "This study aims to provide a potential feasible direction to ease the burden of radiologists by providing not the exact label classifications from machine learning models. But utilize models to excluded consirable proportion of mammograms which were believed to be very confidently classified. And left the rest to be further reviewed by radiologists."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**2. Describe the relevant statistics of the data. How were the images taken? How were they labeled? What is the class balance and majority classifier accuracy? How will you divide the data into testing, training and validation sets?**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Image statistics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1. Perspectives of images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>patient number</th>\n",
       "      <th>img name</th>\n",
       "      <th>img path</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>00001</td>\n",
       "      <td>LEFT_CC.h5</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>00001</td>\n",
       "      <td>LEFT_MLO.h5</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>00004</td>\n",
       "      <td>LEFT_CC.h5</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>00004</td>\n",
       "      <td>LEFT_MLO.h5</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>00004</td>\n",
       "      <td>RIGHT_CC.h5</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  patient number     img name  \\\n",
       "0          00001   LEFT_CC.h5   \n",
       "1          00001  LEFT_MLO.h5   \n",
       "2          00004   LEFT_CC.h5   \n",
       "3          00004  LEFT_MLO.h5   \n",
       "4          00004  RIGHT_CC.h5   \n",
       "\n",
       "                                            img path  \n",
       "0  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  \n",
       "1  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  \n",
       "2  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  \n",
       "3  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  \n",
       "4  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  "
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "img_info = getImageList('data_fixed_crop_w_mask')\n",
    "img_info.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Actually we can see that for each patient, there are multiple images related to them, and each of them are named based on the perspective. In order to have the complete list of all the perspectives of the images, we take the unique values of image names."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['LEFT_CC.h5' 'LEFT_MLO.h5' 'RIGHT_CC.h5' 'RIGHT_MLO.h5']\n"
     ]
    }
   ],
   "source": [
    "print(img_info['img name'].unique())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So all the images are taking from left or the right with bilateral craniocaudal (CC) or mediolateral oblique (MLO) views."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(2085, 4)\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>patient_id</th>\n",
       "      <th>patient number</th>\n",
       "      <th>img path</th>\n",
       "      <th>pathology</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>P_00001</td>\n",
       "      <td>00001</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "      <td>MALIGNANT</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>P_00001</td>\n",
       "      <td>00001</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "      <td>MALIGNANT</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>P_00001</td>\n",
       "      <td>00001</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "      <td>MALIGNANT</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>P_00001</td>\n",
       "      <td>00001</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "      <td>MALIGNANT</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>P_00004</td>\n",
       "      <td>00004</td>\n",
       "      <td>G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...</td>\n",
       "      <td>BENIGN</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  patient_id patient number  \\\n",
       "0    P_00001          00001   \n",
       "1    P_00001          00001   \n",
       "2    P_00001          00001   \n",
       "3    P_00001          00001   \n",
       "4    P_00004          00004   \n",
       "\n",
       "                                            img path  pathology  \n",
       "0  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  MALIGNANT  \n",
       "1  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  MALIGNANT  \n",
       "2  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  MALIGNANT  \n",
       "3  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...  MALIGNANT  \n",
       "4  G:\\Code\\BIOMEDIN260\\Problem Set 3\\data_fixed_c...     BENIGN  "
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Build a DataFrame\n",
    "description = pd.read_csv('mass_case_description_train_set.csv')\n",
    "fn = lambda x: 'P_' + x\n",
    "patient_id = img_info['patient number'].apply(fn)\n",
    "patient_id.name = 'patient_id'\n",
    "pred_info = pd.concat([patient_id, img_info['patient number'], img_info['img path']], axis = 1)\n",
    "pred_info = pd.merge(pred_info, description[['patient_id', 'pathology']], on='patient_id', how='inner')\n",
    "print(pred_info.shape)\n",
    "pred_info.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see after matching the labels with their images, we 2085 images with labels."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2. Labels of images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['MALIGNANT' 'BENIGN' 'BENIGN_WITHOUT_CALLBACK']\n"
     ]
    }
   ],
   "source": [
    "print(pred_info['pathology'].unique())\n",
    "labels = pred_info['pathology'] == 'MALIGNANT'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "All the images contain three kinds of labels. Benign without callback indicates that a physician marked some aspect of the case as worth tracking, but did not require further followup[<sup>4</sup>](#s4). We set class `'MALIGNANT'` as 1, while the other two classes are set to 0."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 3. Majority classifier accuracy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The majority classifier accuracy is given by the proportion of the majority cases in the dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "proportion of positive cases =  0.47817745803357314\n",
      "So the accuracy =  0.5218225419664269\n"
     ]
    }
   ],
   "source": [
    "# proportion of positive cases (malignant)\n",
    "pos_ratio = np.mean(pred_info['pathology'] == 'MALIGNANT')\n",
    "print('proportion of positive cases = ', pos_ratio)\n",
    "# So the accuracy would be\n",
    "print('So the accuracy = ', 1 - pos_ratio)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4. Data division"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `train_test_split` function is used to do the feature selection. I chose test set to have 15% of the data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**3. Describe your data pipeline (how is the data scrubbed, normalized, stored, and fed to the model for training?).**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. Image resize: In `loadImages` function, if size parameter is set, then the image will be loaded as the specified size. Multiple kinds of sizes are considered and the original size (256, 256) works fine.\n",
    "2. Mask remove: Masks of the image are automatically removed by `loadOneImage` function. To just keep the greyscale image for further analysis.\n",
    "3. RGB conversion: Greyscale images are converted to 3 channel-RGB images by duplicating the grey channel three times.\n",
    "4. Image rescale: All images are rescale to 0-255 integer data type.\n",
    "5. Image Storage: All images are loaded as an `ndarray` by `numpy`. There are four dimensions in total, which are respectively (n_examples, image_width, image_height, channels)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████████████████████████████████████████████████████████████████████████| 2085/2085 [00:22<00:00, 92.71it/s]\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(2085, 256, 256, 3)\n"
     ]
    }
   ],
   "source": [
    "# laod all images in a ndarry and rescale all the images to be `scaled_imgs`.\n",
    "imgs = loadImages(pred_info['img path'])\n",
    "scaled_imgs = imgScale(imgs)\n",
    "# The shape of X\n",
    "print(scaled_imgs.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**4. Explain how the model you chose works alongside the code for it. Add at least one technical citation to give credit where credit is due.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This study started with a model mensioned in an published article, which is an \"end-to-end\" classification model that requires only image level label to train and can achieve high performance[<sup>5</sup>](#s5). A pre-trained model is taken and used to give predictions and image features."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:Error in loading the saved optimizer state. As a result, your model is starting with a freshly initialized optimizer.\n"
     ]
    }
   ],
   "source": [
    "# load the pre-trained model\n",
    "resmodel = load_model('ddsm_resnet50_s10_[512-512-1024]x2.h5')\n",
    "# And try to make prediction using this model\n",
    "preds_matrix = resmodel.predict(scaled_imgs)\n",
    "preds = preds_matrix[:,1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The confusion matrix of the model predictions on our dataset.\n",
      "[[860 228]\n",
      " [376 621]]\n",
      "The accuracy of the model prediction.\n",
      "0.7103\n",
      "The AUCROC of the model prediction.\n",
      "0.7067\n"
     ]
    }
   ],
   "source": [
    "# Brief summary of the prediction results\n",
    "from sklearn.metrics import roc_auc_score, accuracy_score, confusion_matrix, plot_confusion_matrix\n",
    "print('The confusion matrix of the model predictions on our dataset.')\n",
    "print(confusion_matrix(labels, preds > 0.5))\n",
    "print('The accuracy of the model prediction.')\n",
    "print(round(accuracy_score(labels, preds > 0.5), 4))\n",
    "print('The AUCROC of the model prediction.')\n",
    "print(round(roc_auc_score(labels, preds > 0.5), 4))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that this model is doing fairly good in three results. However, when we look at the calibration reliability curve, the result is not that satisfying."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO2dd3xV9fnH398kBEggIYuRkBAgYW/CEGSj4N4W3K3WWrW22lpXt9X6s3W1tSruUVfBgYpgtbLDli0jBEjYkECA7OR+f398c+ESMu44dz/v1yuvO8655zwnN/nc5z7fZyitNYIgCELwE+FvAwRBEARrEEEXBEEIEUTQBUEQQgQRdEEQhBBBBF0QBCFEiPLXiZOTk3VmZqa/Ti8IghCUrF69+ojWOqWhbX4T9MzMTFatWuWv0wuCIAQlSqndjW2TkIsgCEKIIIIuCIIQIoigC4IghAgi6IIgCCGCCLogCEKI0KygK6VeU0odUkptbGS7Ukr9XSmVp5Rar5QaYr2ZgiAIQnM446G/AUxtYvsFQHbdz+3AC56bJQhOUrgCFj1lbgUhzGk2D11rvVApldnELpcBb2nTh3eZUqqdUqqT1nq/RTYKQsMUroA3LgZbNUS2hJtnQ/pwf1slCI1SVlVD0ckq0hNjvHJ8K2LoaUChw+M9dc+dhVLqdqXUKqXUqsOHD1twaiGs2bUIaitB26C2yjwWhABlad4Rpj67iDveWY3N5p05FFYIumrguQat1VrP0FrnaK1zUlIarFwVBOfp7OCNR0ZD5hj/2SIIjVBSXs2Ds9Zz3SvLiVDw24v7EBHRkGx6jhWl/3uAdIfHnYF9FhxXEJpGOfgj09+VcIsQcNTaNFe9sJT8wyf5ybhu3Du5B61aRHrtfFYI+mzgbqXU+8AIoETi54JPKFx2+n7bTv6zQxDqcbS0inYxLYiMUPzq/J6ktmvFgM7tvH5eZ9IW3wNygZ5KqT1KqVuVUncope6o22UOkA/kAS8Dd3rNWkFwpGDZaS/9+F7/2iIIgNaaj7/bw4Sn5vP+SrO0OLVfR5+IOTiX5TK9me0auMsyiwTBGWw2KFwOXcdB/rdwXKJ8gn/Zd6ycRz7ewLdbDzM4ox05XRJ8boPf2ucKgkcc3gIVJdD3csifL4Iu+JVP1+7lkY83UmvT/O7iPtw8KpNILy18NoUIuhCc2OPnmWOgTXso2eNfe4SwJr51Cwalt+MvV/b3Wo65M4igC8FJwXKITYHEbhCXJh664FNqam28ungn1bU27p6Yzfie7RnXIwWlfO+VOyKCLgQnBbmQMRKUgrhUKNrhb4uEMGHzvuM8MGs9G/aWcNGATmitUUr5XcxBBF0IRo7vh2O7Yfjt5nFcGuyUKlHBu1TW1PLP/+XxwvwdtItpwb+uH8IF/ToGhJDbEUEXgg97/DxjpLmNS4XKEqg8AS3b+s8uIaTZdaSMFxfs4NJBqfz2oj4kxEb726SzEEEXgo+C5RDVGjoOMI/j6loHHd8HKT39Z5cQcpRW1vDfzQe5fHAaPTu25Zv7xpOR5L9Fz+YQQReCj4JcSBsKUXUeUrxd0PeKoAuWsWj7YR76aAN7j5XTLy2OrPZtA1rMQSYWCcFG5Uk4sOF0uAVMyAUk00WwhJKyan49cx03vrqC6MgIPrj9HLLaB0coTzx0IbjYuxp07ZmCbu/jIoIueEitTXPVi0vZeaSUO8d3555J2V5tpmU1IuhCcFGwDFDQedjp56Jampx06eciuElxaRXtWptmWvdP6Ulau9b0S4v3t1kuIyEXIbgoXAbt+0Dres2O4lKhRARdcA2tNbNW72HC304305rSt2NQijmIhy4EE7ZaKFwJA645e1tcGhzd7XubhKBlz9EyHv54Iwu3HWZolwSGd030t0keI4IuBA8HN0HVCcg45+xtcWmwe6nvbRKCko+/28NvPt6IBv54aV9uHNnFa1OEfIkIuhA8FC43t+kjzt4WlwoVx6CqFKJjfWuXEHQkxrZkaGYij1/Rj84JgZ2K6Aoi6ELwUJALbVOhXcbZ204VF+2H5Czf2iUEPNW1Nl5elE9NreaeSdmM65HC2OzkgCrbtwIRdCF4KFgOGSNMQ676nMpF3yOCLpzBxr0lPDBrPZv2HeeSgakB1UzLakTQheDgWKER6/SfNbxdiouEelRU1/L3b7bz0sJ8EmKiefGGIUztF9qzZ0XQheDAHj93LChyJM6h/F8QgN1FZby8KJ8rB6fxm4v6EB/Twt8meR0RdCE4KMiFFrHQoV/D21u0gpgk8dDDnNLKGuZtOsCVQzrTs2Nb/vfL8X6dIORrRNCF4KBgOaQPg8gm/mTjUkXQw5gF2w7z8Ecb2FdSzoDO8WS1bxtWYg5SKSoEAxUlcGgTpDcSbrETlybVomHI0dIq7vtwLTe/toJWLSL4z0+Cp5mW1YiHLgQ+e1aCtjUeP7cTlwqFK3xjkxAQ2Jtp7S4q4+4JWdw9MSuommlZjQi6EPgULAcVAZ1zmt4vLg3Ki6G6HFq09o1tgvfYvQx2L4GuYyB9+Bmbik5WkhATTWSE4sGpvUhLaE3f1ODsv2IlEnIRAp+CXLMY2tx4OcfJRUJws3spvD4F/vcovHnpqW9eWms+XFXIhL/N572VBQCc37ejiHkdIuhCYFNbbXqgN9S/pT6Six46bPmi7o6G2irYtYjC4jJuem0Fv565nl4d4zinW5JfTQxEJOQiBDYH1kN1makQbQ7JRQ8d2nQ8fT8ymvmVPbnz2YUo4NHL+3H98IyQaKZlNSLoQmBTYG/I1cyCKDh46CLoQU9EXfAgqhXcPJuI8q4M77qTx67oT1o7WR9pDBF0IbApXAbxGacHQTdFdAy0TpCQSwhQe2Q7kQA1FZCUxdiYRMb2SPG3WQGPxNCFwEVrM3KuuXRFR+LSRNCDnI17S9iwbjU2bUIq+sg2P1sUPIigC4HL0V1w8qBz8XM7cakScglSKqpreeLLLVz2/BI61eylJHkQAKpoh58tCx6cEnSl1FSl1FalVJ5S6sEGtmcopb5VSn2nlFqvlLrQelOFsKNgmbl1Jn5uR2aLBi0FxWW8ujif6YOS6UARCf2mQEQLKNrub9OChmYFXSkVCTwPXAD0AaYrpfrU2+03wIda68HANOBfVhsqhCGFy6BlPLTv7fxr4tKg7AhUV3jPLsEyTlRU859VZjhzjw5t+fZX4/nz2Lr+Kyk9IbErFOX50cLgwhkPfTiQp7XO11pXAe8Dl9XbRwNxdffjAQliCp5jb8gV4UIptz118cR+79gkWMa3Ww4x5ZmFPDBrPXmHTgCYcXB2AU/KMj9HRNCdxRlBTwMKHR7vqXvOkT8ANyil9gBzgAanECilbldKrVJKrTp8+LAb5gphQ1kxHP7etQVRkOKiIKC4tIp7P1jLD99YSWzLKGb+dNSZzbTsgp7YzQh6cT7Yan1nYOEKWPRUUPYFciZtsaHsfV3v8XTgDa31U0qpc4C3lVL9tNa2M16k9QxgBkBOTk79YwjCafasNLeuxM9Byv8DnFqb5uoXllJQXMY9k7K5a0J3WkbV+wZWtMPMjm3Zxgh6bSWUFEJCpvcNLFwBb15i0iUjIuGiZ2DITQ2PPQxAnBH0PUC6w+POnB1SuRWYCqC1zlVKtQKSgUNWGCmEIQW5EBEFaUNde11c3Yix43ust0lwm8MnKkmKNc20Hr6wN2kJrendKa7hnYvyIKm7uZ+cffo5Xwj6rkVQU2nu22rhs3tg8TMw4Acw4NrTdgUozoRcVgLZSqmuSqlozKLn7Hr7FACTAJRSvYFWgMRUBPcpWA6dBppiIVdo2dYspIqHHhBorflgZQETn5rPuytMM63JfTo0LuZQJ+h1g77tt76Ko2eOOb1mE9UKxvwK2mXAgv+DfwyBlyfB8hlQesQ39rhIsx661rpGKXU3MA+IBF7TWm9SSv0JWKW1ng38EnhZKXUvJhxzi9ZaQiqCe9RUmoZcw25z7/XxXi4uKlxhPLnMs9u6CqcpKCrjwY/Ws3RHESO6JnJuVnLzLyorNi2Q7UIem2I+oH2V6ZI+3Lyv+76D6/9z+v09vg82zIT1H8KX98PcByFrsvHae17ouuPhJZwq/ddaz8Esdjo+9zuH+5uB0daaJoQt+9eZuKmrC6J2vFlc5BhjjWoJN38uot4AM1fv4befbCQyQvHYFf2YPszJZlr2IiK7oCtlwhy+zEWvKIHUwWe+r3GpMPoe83NwkxH2Df+BWfMgug30vsSIe9dxrmVlWYz0chECD3tBkSeCfmCDdfY44hhjrTFtXUXQz6ZDXEtGdU/iz1f0o1O8C820HFMW7SRnm/7ovqI4H/pf3fj2Dn3hvD/CpN+bARzrP4DNn8K690yXyP5XG3HvOMDni6ki6ELgUbDMpKy1ae/e6+PS4OQhI7hR0dbaljnGTE/StYCGjFHWHj9Iqaqx8cL8Hdi05t7zejAmO4Ux2W400yrKAxUJCV1OP5eUZUSzqsz7oY2yYqg4Zv7+miMiwkxT6joGLvwbbJtrvPblL0HuPyGllxH2/teYOLwPEEEXAgutTYVo9hT3jxGXCmhTXOQoDFbQeRi0bmfG3FWXmZ8wZ13hMX49cz1bD57gysFpaK1R7nqm9myWyBann7N768X50LGfx/Y2SXG+uXVG0B1p0Qr6Xm5+yoph8ycmLPPNn8xPxigj7u3STUjRS+svIuhCYFGUB2VF7odb4MziIqsF/dhuY9+Ux03mw7r3IGuStecIEsqrann6v1t5dfFO2rdtxSs35TC5TwfPDlq048xwC5x+XLQ9cAXdkZhEyPmR+Tm6y3jt6z6Az39xep+o1nDzbMtFXbotCoGFp/FzgLjO5tYbC6O7c81tt/HQ7yr4/nOoOG79eYKAwqNlvLl0N9OGZ/DVfWM9F3ObDYobEvS63G9fZLoU5wMK2lnkCCRkwtj74e6VMOzHnKrTrBurZzUi6EJgUbjMDKlIynb/GN4s/y/INWl0Kb1h4HVQU24WxMKE4xXVfOjQTGv+/eN5/Ir+xLVq0cwrneDEfhPCql+8Ex1r1kV8kYtenA/x6SaEYiVKmZBLVCuzRhAZbcIuFiMhFyGwKFhmyv0jPPA1WsVBdFvveOgFy0x/9ogI6JxjvMl178GQG60/V4Dxvy0HefijjRw6UcGQjASy2rch1cpxcA1luNhJyvKNh160w3R49Abpw02YxYs1DOKhC4FD6RHzT+vKQIvG8EYuemkRHNl6OhykFAycZlLXju6y9lwBRNHJSn7+/nf86I1VxLduwUd3jiarfRvrT1Rsz0FvoLw+KcvE0L1dr1ic71n8vDnSh8OYX3ot1VUEXQgcTsXPz/H8WHGp1odcCu32OaQqDpgGKLPoFYLU2jTXvJjLnA37uXdyDz772bkMSm/nnZMV7TCLhW1Tz96WnG0KfsqKvHNugPKjpkrVm4LuZSTkIgQOhctMbLHTIM+PFZ8GeVs8P44jBbnGvtTBp59rl27ykNe9B+N+HTRd+Zrj0IkKkmNbEhmheOSi3nROiKFnx7bNv9AT7E25Ggq3nerpsh1inWgh4A7FO81tEAu6eOhC4FCwDFKHWLMgFZcGJw5AbbXnx7LTmH0Dp8PRnVC43Lpz+QmbTfPv5buZ+LcF/Luumdak3h28L+ZwZpfF+pxKXfRiHN2eshjgHRWbQgRdCAyqy2HfWmvi5+BQXHTAmuNVlZmGTQ2lU/a+FFrEwNp3rTmXn9h1pJTrXlnGIx9vZEDneMa5U+npLrXVZh2ioQVRMJWWkdHe7eliF3RftOn1EiLoQmCw7zuwVbs+0KIxrB50sXc12GqgSwOl/i3bGFHf9In5YApCPlxVyJRnF7Jp73GeuLI//75tBBlJPuwgeKzA/H4bE/SISBMK8WbqYnG++btpYWHmjo8RQRcCg4K6gp10Kz10rMt0sS/YNpadMGg6VJbA1jkNbw9w0tq1ZmyPFP573zimDc9wv3TfXZpKWbTj7dRFb2e4+AARdCEwKFgOyT0gNsma41ntoRfkQvs+puipITLHmgrVde9bcz4vU1lTyzP/3cbTX20FYHRWMi/flEPHeIsLapzFWUEvzofaGu/YUJzvvRx0HyGCLvgfm81kuHhS7l+fVvHQItYaQbfVmj7oTaVTRkSYSsC8b+DEQc/P6UW+KzjKJf9YzHPfbGfvsQoCYhZNUZ75sIxJbHyfpCwTlispsP78Fceh9LB46ILgMUe2mhxjq+LnYNIH41KtmS16cCNUnWg+P37gdNNWd8OHnp/TC5RV1fDo55u58oWlnKio4bVbcnjq2oG+D680hOPYucawzxf1Rhz9VFOu4M1wARF0IRCwx8+t9NDBuuIiZxuGpfQwQ60DNOyy92g5by/bzfUjMvjq3rFM7OVhMy0raajLYn2SHAZGW40VXRYDABF0wYQTFj1lbv3B91+Y8IjVg3fjLJotWpBr4uPt0pvfd+B049HvX+/5eS2gpLya9+vyybM7tGXB/eP58+X9aWtFMy2rqCo1i9fN5X/HJEKrdt5JXTwl6BJDF4KZwhXwxsXwzaPw5qW+F/Ud82HHN1BdCm9dZu354+3FRR4somltWuZ2cbIdQb+rIKKFqRz1M19tOsB5Ty/gkU82knfoJIBr4+B8xamCnmY8dKVM2MUrHvpOMz4uOtb6Y/sQEfRwJ3+BGciMNoOPvdCjufFzz4cPbzLnBut7RMelmph26SH3j3F0F5w84Hw4KCYRek41Qw2srFJ1gSMnK7n73TXc/vZqEmOj+fjOUd5ppmUVzmS42EnK8l4MPcjDLSCCLtQ4FsJo01HQ21SehC9+aTzy6DYQ2dI7PaLtqYslHuSiu9MwbOB1JmMi7xv3z+smtTbN1S8s5atNB/nV+aaZ1oDOXmqmZRV2QXdGUJOy4MQ+8zdkJcX5kBT8gi7NucKZmirjSSb3gP7Xwo7/wbJ/QeehJnTgDXYthk/uNJWBI++Cib8xMWdv9Ig+o7homHvHKMg1KZApvZ1/TdZkiEkyYZeeU907r4scPF5BShvTTOv3l/Slc0Jrsjv4oP+KFRTVVWg6E+44NV90B3QaaM35K0+ab2HioQtBzdp3jLBOeRzG3Q83fmQ80Y9+Aju+tfZcVaXw5QPwxkWgIuCHc2Dq42aKu7d6RFtRXFSQ6/rAjahoM+l96xzTktVTmli0ttk0by/bzaSnFvDv5bsBmNCrvbVi7u1F86I858U02QuZLkeDv8uiHRH0cKW6Ahb+DToPNx4lmB4W098zHvsHN8DeNdaca3cuvDAalr8II+6Any5puCeK1bROMP213S3/Lz0CR7a5l045cLpZE9j4kXvntrPtK3h9qpkc/+YlZ4hq/uGTTHt5Gb/9ZCOD0tsxvmd7z87VEIUr4PUL687vpUVzZ3LQ7SR2A5S1cfQQSVkEEfTwZc2bRugmPnJmD+/W7eCGWdA6Ef59jckPdpfqcpj7MLx+AWgb3Pw5XPB/vsskOFVc5KaHbm+H687AjU4DTZjG3Zz0qjLjFX9wg6lUBaipPLVo/MHKAi54bhFb9h/nyasH8Patw0lP9EIzrbxvTHUmmMVzqxfNy4rNUAlnBb1FazPz00oPXQRdCGrsYtFlNHQdd/b2uE5w48eAhrcvd68FbeEKePFcWPY8DLsVfrrUDILwNZ6Motu91CzYpg1x/bVKmYZde1a45k3aauG7d+AfQ41XnDYEolrWbdTmGxXQOSGG8T1T+Pq+cVybk+69as+KYw4PlPWDje0Og7OCDiZf3cpc9OJ8iG0PLYNkzaEJRNDDkVWvwsmDMOGRxifsJGfB9f8xWS/vXAXlxxrerz7VFfDVb+G1KcajvOlTuOgp02LWH3hSXFSwrJ6gukj/a816wXonvHStYfvX8OIY+PQu86F6yxz40Vy4+XNqBt0EwNcLFwKmmdZLN+bQPs6LzbRqq2HLF9Chv/kgUcp6L9aVlEU7ydnmg8CqHjRFoZGyCCLo4UflSVj8LHQbD5mjm943bSj84G04vBXev86IdVPsWQ0vjYWlf4chNxmvvNt4a+x2l7hUOLH/dNjCWarKYP9az9oRxHWCbhNM2MVma3y//etMCue/rzIFVle/Drd9c+r9WW3LYkrelayy9WDInrfQNVXu2+QKmz6BkkITlrv0H6Zf+XfvWHuOojyTsprQxfnXJGVB5XE46UF9gSMhkoMOIujhx4oZUHYEJvzGuf2zJsEVL5rJ9rNubVgYayrh6z/Cq5Oh6iTc8BFc8hy0irPWdneITzNCVHrYtdfZB1p4OrB60HVGFHcvPnvbsQKTUfTSODiwAab+H9y1EvpdCUpRWlnDH2Zv4uoXc6mo0bSccD+J1QdQG2d5ZpMzaA1LnoPknpA9Bdr3gi7nwurXm/5wcpWiPDMhKNKFVgRWjqOrKjN57eEk6EqpqUqprUqpPKXUg43sc61SarNSapNSKrhncYUqFceN95x9PqS7kJfd/2ojNls+h8/vPfOr7r7vYMZ4WPy0Ea87c82HQKBwKnXRxTh6QS6gPE+l7HkhRLc9c3G0/JgJS/0jBzZ/Auf+An6+FkbeYVIe69h3rJx3VxRw08guzLt3LP3HXwMd+sHiZ6wV1YbI/xYOboBRPzudsjnsR6Zydsf/rDuPM0256nNK0C2Iox/dZW6DvIeLnWYLi5RSkcDzwHnAHmClUmq21nqzwz7ZwEPAaK31UaWUF/KnBI9Z9oLJi57wsOuvHXmHibsvftp4rgmZ5p9h3fvQpj1c9x/ocb7VFnuOvbioZK8JITlLcwMtnCU6BvpeDhtmGlsqT8D6D4yoD7rOvBfxnU/tXlJWzRcb9nPdiAyyO7Rl0a8n0MExTn7uveab0tYvoPclntnWFEueM71NBlx7+rlel5jFw1WvQvZkz89hs5kCoa5jXXtdfLpZrLbCQw+BwdCOOFMpOhzI01rnAyil3gcuAzY77PNj4Hmt9VEArbVFwS3BMsqPQu7z0OtiSB3s3jEm/Q4ObYa1/z79XNZ5cNXLnguft3CnuKi2xmTpDPiBNTakDobv3oaFfzWP03LgkmehY/8zdpu78QC//XQjxaVVjOiWSPeUNmeKOUDfK+Dbx0yWUq+LG1/U9oT960yfncl/OHNBOCoahtxoviEcK3Su+2RTnNgP1WWui2lEhHmNFbnoxXVZNgmh4aE7E3JJAwodHu+pe86RHkAPpdQSpdQypVSD9c5KqduVUquUUqsOH3Yxpil4xtJ/mpmX4x9y/xhKGTE69TjCdCEMVDEHU4IfGe1ayOXQJrMW4Gn83I5jtaiKgF4XniHmh05UcOe/V3PHO6tJadOST+8aTfeURrKCIiJh9C9MqCvf4mpeO0v+bsJEQ3949raht5iQ25o3PT+POxkudqyaL1qcb/5GWgd4vxsncUbQG3IB6ucLRQHZwHhgOvCKUuqs35DWeobWOkdrnZOSkuKqrYK7lBaZKs2+V0DHfp4dq9s4U32pIs3XXqvzkq3GneKi3XUDN5xtmdscXcc2+jurtWmufTGXr78/xP1TevLp3aPplxbf9PEGToO2qbDoaWvsc+Tobtj0MQy9uWGRa5cBPabAmrc87ybpqaAf3em5DSGU4QLOhVz2AI7frToD9f879gDLtNbVwE6l1FaMwK+0xErBM5Y8a77aeuKd20kfDjfP9k4zLW8R19k1QS/INXFah9i2RzTwO9tfUk6Htq1MM61L+5KeEON8i9uoljDqbpj3sBmunTHCGjvBNGdTCkbe2fg+OT+Cd681i+R9r3D/XEU7oEUMtO3k+muTs81azrECz+LfxTt904bCRzjjoa8EspVSXZVS0cA0YHa9fT4BJgAopZIxIZh8Kw0V3OTEQVjxsmkWldLTmmN6q5mWt3BltqjWpqDI6nF4db8zW9ow3liyk0lPLeAdezOtnu1d71c+5GbTnmGxhV56WbHxvPtfY9I9GyNrMsRnwMpXPTtfUZ6Z4elK4zM7dq/+iAeZLtUVULIn6OeIOtLsb1JrXQPcDcwDvgc+1FpvUkr9SSl1ad1u84AipdRm4Fvgfq21DxprC82y+BnTJGrcA/62xH/EpcLx/c6l+p0aaGFRuMWBvEMnufalXP7w2WZyMhOZ2MuDZLCWbWDkT2HbXDiw0RoDV75qvsmN+lnT+0VEQs4t5hvH4W3un68oz33v2opc9KO7AB1SIRenPhq11nO01j201t211o/VPfc7rfXsuvtaa32f1rqP1rq/1jowp+SGGyV7YdVrpqdIiKRluUVcmmkwVebEzNJTA6utFfT3VxRw4d8Xsf3QSZ66ZiBv/nAYnRM8bKY1/MdmQMjiZzw3sLocVrxkspY69G1+/8E3mVF7q19373y11XBst3vxczCToVonepaLHkJNuexIpWgos+gp0+Vw7K/9bYl/OWPQRTOcGmjRy1ITMpJimNy7PV/fN46rhna2pplW6wTT+GzTR551xQQzjKP0MIz+uXP7t0mBPpeaFNaqMtfPd6zAxMA9cTTsPV3cJUQGQzsigh6qHCsw8dAhN7rWJyMUiXchF71gmesDLRqgorqWJ+du4cm5WwAY1T2Zf10/lJS2bjb6aoyRdxlPeclz7h/DVmvSWlMHQ+a5zr8u51aoKDEfKK7iSYaLnaQsz2LoxfnQqp3x9kMEEfRQZcGTJud5zK/8bYn/cXa2qH2ghYfpiqt2FXPh3xfxr/k7KC6tQlvVFbAh2naAwTfA2nfd7yq55QtTYDP6564VKnUZZb7JuLM4apWgnzxgqm/dIcRSFkEEPTQp2mH+wXN+2HS2QrgQk2y82OZCLu4MhHbgZGUNv/90I9e8lEtVjY23fjScJ64a4L1e5XZG32NCa0v/6fpr7U24EjKh96XN7n4GSpkUxn1rTKGTKxTlmZCRJ96xpwujxTtCbm1JBD0UWfCkqY489z5/WxIYRESYVrbNebAFuabwx83WCAdKynl/ZSE3n5PJvF+MZWwPHxXPJWSaVMPVr5siMlcoyIW9q+Ccu032iqsMnGZyyV310l0ZO9cYp+aLuhFHr6msS1kUD12wGiuH8B7eBhs+hOG3ma/jgsGZQRduDLQ4WlrF28tMPnlWe9NM6w+X9iW2pTM1exZy7r0m5XDFS669bsnfTen7oOvdO2+reNONc753nv4AABy3SURBVMNM54eggHtdFuuT0BUzX9SNOPqxAvOtRgRdsJTCFWb47//+bM0Q3vl/MWXmo39hjX2hQlxa0yGXqtK6gRbOhVu01szZsJ/znlnAH2dvYsfhkwDenSDUFO17mWZdy190PqZ8aAts+xKG3266QrpLzq1QU+78/NSqUvNeeBruaNHKtCJwJ3UxBFMWQQTd/+xcCDUVxluoqTCtVd3l4CaTcTDyDohNts7GUMDez6WxBUoXBlocOl7BHe+s5s5/r6FTfGtm331u4820fMmY+0zWyarXnNs/9x/mw3/Yjz07b+og05p41WvOjYU71bLWQw/dfgx3Yugi6IJXiHWMs2pY+Qq8MhnW/wdcHTX27ePQMs7EQ4UziUszU+vLGokxFyzDmYEWtTbNNS/lMn/rYR66oBcf3zmKPqkBMJkJjKh2G28WR5sbF3h8P6z7wGTIxCZ5fu6cW+HIVtjVwGSm+liR4WLH3fmiRTugZbwJN4UQIuj+pqQQUDD2frjxYzMZqKwYProNnu0H858w/ViaY99a0yxp5J0hlVdrGc0VFxXkmgrJRtqo7jtWjs2miYxQ/Omyfnz58zH8ZFx3oiID7F9ozC+h9BCsbWb25/IXQdfCOXdZc96+V5h4ujPfDuyCboV3nJRlWh2fOODa64rzTUGRtzOQfEyA/TWGIdvmmq/5E38D3SeacMndq+D6WdBxgImJP9MXZt0Ge1Y1fpxvHzdFEuc00SUvnGlq0IV9oEUDDblqbZrX6zXTGtcjhW6BEGJpiMwx0HmYSUVsrLVsxXFY9Tr0ucy6KsnoGLOw+v1nzQ9vLtph3o/oWM/P6+44uhDMQQcRdP9SstcMB+4x5cznIyLMiK8bZsLP1sCw22DbPHhlEsyomyJfU3l6/8KVsH2eaarUqple2uFKfBOzRQ9ubHCgRd6hE1zz4lL++NlmRnRLZFLvIMgaUsp46ccKoLFh0mveNMNORt1j7blzfmR65qx5q+n9PGnKVR93ctFrq83vRwRdsJTt88xtjwYHPBmSusMFT8B9m+HCvxnh+fgnxmv/32OwZQ58dLuJnY+4wzd2ByOxKRAR1XC16KmCotMe+rvLC7jwucXsPFLKMz8YyOu3DCOtXWsfGesh2VOgfV8zAKN+h8maKjNbNnOMSdG0kuRsM8xj9RumnUBjWJGDbicuzSzsujKO7liBCTeJoAuWsu0raNfFuT7lLdua7np3rTCx9rQcWPgkvD8djuabHORDm5s/TrgSEWkGKTQUcinINf29HQZaZCbHcH7fDvz3vnFcMdiiZlq+IiLC5KUf2Qpb55y5beMs8y3F2SZcrpJzq1kX2v7fhreXFZuRfFYJun2+qCseeogNhnZEBN1fVJebQbw9pri2MKOUibVf9z6c8zNOTQjU2vSnFhonLvXskIvWUJBLbfoI/vLl9zzx5elmWv+8bgjJbSxupuUr+l5hKkgXPXU6A0RrWPp3aN/HDKnwBr0ugjYdYVUjlaNWZrjYScpyLYYeoimLIILuP3YuMsUY9ePnrtDnUohqVTerMjrw53v6m4Zmix7dCScP8tz2JF5akM+JimrvNtPyFZFRdcOk1xjHASDva/MtbtQ93svuiGwBQ24yHvrRXWdvt5fpWy3oR3c7n+ZbnG/6yMeG3lxjEXR/sW0utIiFLi60K62PfVblxEfMbbCMhPMX9vL/OsE+UVHNrE9mArBO9ebd20bw2BX9gyu80hSDrjPe8qKnzOMlz5nfQb+rvHveoTebD4zVb5y9rSjPOCDtMqw7X3K2iYk39AHSECGasggi6P5Ba5O10n2CKV/2hGCb7+lP4tLMt6LyowAcPF6JbddSyiPb8sK91zEqK8Sqa6NamsynXYtg+QxzO/KnEBXt3fPGdzYL/WvePjMbC4ygJ2QaT94qXM10CdGURRBB9w+HNpuhxZ6EWwTXqSsu+nzRSgCy2rfhiqRCWncfTUxLL4ucvxh6i2lT++X9ppNkh36+OW/OrWbk3/efnfm8FU256mNf3HQmjl5bY8IzITQY2hERdH+wba65zT7fv3aEEVprFh0yov3popXkHz4JpUeIOprXYEFRyNCyjWnaBWZY+HvTrenq2RzdJxpP3LFy1Gar60FusaC3TjA9753x0EsKTa68eOiCZWybB50GQduO/rYkLDh4vIIfv7Wa++eZIdGPTkw0lZ6n8s9H+dE6HxDXCZMNpY2o+yIbKiIChv4Qdi+BQ9+b507sN+m13kgXTM52Lhc9hDNcQATd95QWGQ+pqWIiwTJqbZprX8pl0fbD3DZ1JFpF0JFis/HUQItB/jXS22Sd559sqME3mPPZvXRvpCzacTYXPcQF3cdd+AXy/gtoiZ97mT1Hy+gU35rICMWjl/UjIzGGzORYWNPpdLVoQa7pUOjCQIugxJ4NtWuREXNfLaDHJkOfy02rikm/97KgZ0PpO6Z9cFPtL4p3mglLIfrtWDx0X7NtHrTpYEIuguXU2jSvLMpn8tMLeKduktDYHilGzOF0cVFVKexfF9rxc0f8lQ2V8yOoPA4bZ5oF0RYxpmLXapzNdLFnuIRgyiKIh+5baqsh7xvoc4mJMQqWsvXACX49az3rCo8xqVd7zu/bQDOtuFQ4uPn0QIsuIR4/9zcZI01l6spXjVec2N07f/v2+aJH8sy3rsYo3uFcq40gRVTFlxQsM13uJH5uOe8s283F/1hEYXEZz00bxCs359ApvoFmWvbiot25gDKtZgXvoZTx0g+sN9XR3uqfkpAJKqJpD91WV3wUovFzEEH3LdvmmkWibuP9bUnIYC/Tz2rfhgv7d+K/947lskFpjVd7xqVCdal5L5oYaCFYyIAfmKromnIzMcobaZNRLU2ju6Zy0Y/vNVk+IuiCJWybB5nnms6JgkeUV9Xy2BebeWKuaaY1slsSz00bTFJzzbTsgy72rQmf+Lm/aRV32onZtdiaYegN0dx80RDPcAERdN9RtMN4DxJu8ZjcHUVMfW4hLy/aSVllrWvNtOyCDk4NhBYsIqFL3R0v5sLb54vW7wFvRwRdsIztX5lbqQ51m+MV1Tz00Qamv2wKgt798Qgevbyfa8207LNFQQTdl/S9wgyi8GYufFKWKVw6sb/h7UU7TD5+29SGt4cAkuXiK7bNheSe1s1wDEMOHa/kk+/2cvvYbtw7uQetoyNdP4g9/7hlvImpxqc1vb9gDb7IhXecL9rQ+1q8ExK6hnSGmVNXppSaqpTaqpTKU0o92MR+VyultFIqxzoTQ4CK47BriRQTuUHRyUreWLITMAufix+YwMMX9nZPzAH2fQcok23krViu0DDezoW3py42FkcP4S6LdpoVdKVUJPA8cAHQB5iulOrTwH5tgXuA5VYbGfTkf2saAkn83Gm01ny6di+Tn17AY3O+N820oPlFz+bYteh0UYmv+poIvqFtJ5NN01BPF5vNDDMJ8W/Iznjow4E8rXW+1roKeB+4rIH9HgWeBCostC802DbPlCOnj/C3JUHBvmPl3PrmKn7+/lq6JMXyxT1jTDMtK8gcY/q3yJSn0EOpxnu6nNgPNRUh76E7E0NPAwodHu8BzlAmpdRgIF1r/blS6leNHUgpdTtwO0BGhoUTSwIZm80siGZNNmPBhCapqbUxbcYyDp+o5LcX9+GWUZlERlhYpu2vviaCb0jKMimp9Sm2j74LzT7odpxRmIb+m07liSmlIoBngFuaO5DWegYwAyAnJycEBjc6wb7voPSwhFuaobC4jNR2rYmKjODxK/qTkRhDRlKMd06WPlyEPFRJzobNn5hJSY5N18IgZRGcC7nsAdIdHncGHCfttgX6AfOVUruAkcBsWRitY9tcU5LsrSnrQU5NrY0ZC3cw+ekFvJ27C4Bzs5O9J+ZCaJOUBdpmMlocKc43Iba40M5qcsZDXwlkK6W6AnuBacB19o1a6xLg1DBGpdR84Fda61XWmhqkbJtrYucxif62JOD4fv9xHpi1nvV7SjivTwcu6O+FLnxCeOHYdbF9r9PPF+ebfi8RbmZHBQnNCrrWukYpdTcwD4gEXtNab1JK/QlYpbWe7W0jg5bj+0xTosl/8LclAcfbubv442ebiW/dgn9eN5iL+ndyrUBIEBrCMRfdkeKdIR9uAScLi7TWc4A59Z77XSP7jvfcrBDhVHWo5J/b0VqjlKJHh7ZcMjCV317ch8TYEB3QLPieVnFm3oBj6qLWxkPvOs5/dvkISbvwJtvmQXwGtO/tb0v8TllVDX+bt42oSMXDF/ZmRLckRnRL8rdZQihSv0nXiQOmJUCI56CD9HLxHtXlkD/fVIeGeShhSd4Rpjy7kNeW7KSqxuZaMy1BcJWkrDNDLmGS4QLioXuPXYuNVxDG6Yol5dU8/sX3fLCqkK7JsXz4k3MY3lUWhwUvk5Rl+q6XFZtkBBF0wWO2zTXzEzPP9bclfuPIyUo+W7+PO8Z15xeTs2nVIrQzDIQA4VRPlx2nBT0iCuLTm35dCCAhF2+gNWz7yjT1b9HK39b4lMMnKnltsckB7p7ShsUPTOTBC3qJmAu+o/7A6OIdJmUxDCq1Q/8K/cGh76GkAMb+0t+W+AytNZ+s3csfP9tMWWUtE3q1p2tyrGSwCL4nIdP06rHH0cOgy6IdEXRvsG2uuQ2TYRZ7j5XzyMcbmL/1MEMy2vHk1QPomhzrb7OEcCWyhRH1ory6lMWd0GW0v63yCSLo3mDbPOg08MzpOCGKaaaVS9HJKv5wSR9uPMfiZlqC4A7J2SYXvfQwVJ0UD11wk7Ji2LMCxt7vb0u8SkFRGWkJppnWE1cOICMxhvRE6b8iBAhJWSZt2B5HDxNBl0VRq8n72jQHCtHq0JpaGy/M38HkZxbwVu4uAEZnJYuYC4FFUpbpf76zboBJmAi6eOhWs20uxKZA6mB/W2I5m/aV8MCs9Wzce5wpfTtwkTTTEgIVe6bL9nlmgbRdeMxfEEG3ktpq46H3uiTkBtG+uXQXj36+mXYx0bxw/RDpjCgENvZc9L1r6lIWW/jVHF8hgm4lhcuhoiSkhkHbm2n16tiWywal8duLe9MuRlIRhQCnTQeIbhNWC6Iggm7QGr7/HAqXmfSmzNGmytPVT/XVr5thFi3jvGOnDymtrOGv87bSIlLxyEV9pJmWEFwoZcIu+9eKoIc0lSdN4c/BjXU/m2D/OtN3BSD3n6f3jYgywt6idd2t/X5riI51eL618cw3zAI0vDfNzK0M0jFnC7cd5qGPNrCvpJybz8k85aULQlDROsHchvhQC0dCV9BtNji2ywj2wU1wYIO5Peowmiq6LXToC+37wN7VmFGpdePi0ocbka8ur7t1vF9uhldUl59+rvI4p0at1laZIcRBJuglZdU8+sVmZq7eQ7cU00xrWKY00xKCkMIV5n8QYNVr0O+qoPt/dIfgFXT7G5Y5BlJ6wsHNpz3ugxvN4+rSup2VmfbdaSAMut6IeIe+ZuVbKXOsNy81QhwZDePud/3Nr3+MzDGWX7K3OVJayZcb9nPn+O7cM0maaQlBzK5FxqkDsNUGpYPlDsEp6IUr4I2LobYSUJzyjAFatYOO/WHIjaeFO6U3RDeRJ50+3IRI7B8Q7rzxVhzDDxw6UcHstfu4bUy3U820EqT/ihDsZI6BqJZB7WC5Q3AK+q5FYKuue6Ch+0QY8VMj3nGp7g2USB/uuQhbcQwfobVm1pq9PPr5Zsqra5nUuwNdk2NFzIXQIEgdLE8JTkHPHAORDp++4x8KmzfMCgqLy3j44w0s2n6EnC4JPHGVNNMSQpAgcrCsIjgFPUw/fa2gptbG9JeXcbS0ikcv68v1I7oQIc20BCEkCE5Bh7D89PWEXUdKSU+MISoygievNs20OidI/xVBCCVCqz5dOIvqWhvPf5vH+c8sPNVMa1T3ZBFzQQhBgtdDF5pl494Sfj1zPZv3H+ei/p24eEDo92cXhHBGBD1EeX3JTv78xfckxkbz4g1Dmdqvo79NEgTBy4ighxj2Mv2+qfFcOTiN31zUh/iY8Og0Jwjhjgh6iHCysoYn524hOjKC31zch+FdExneVcr2BSGckEXREGD+1kNMeWYhby/bjcZ46YIghB/ioQcxR0urePSLzXy0Zi9Z7dsw845RDO2S4G+zBEHwEyLoQczRsiq+2nSQeyZmcdfELFpGSTMtQQhnnAq5KKWmKqW2KqXylFIPNrD9PqXUZqXUeqXUN0qpLtabKgAcOl7BjIU70FrTLaUNSx6YyH3n9xQxFwSheUFXSkUCzwMXAH2A6UqpPvV2+w7I0VoPAGYCT1ptaLijtebDlYVMenoBT321jV1FZiCHZLAIgmDHmZDLcCBPa50PoJR6H7gM2GzfQWv9rcP+y4AbrDQy3CksLuOhjzawOO8Iw7sm8sSV/aWZliAIZ+GMoKcBhQ6P9wAjmtj/VuDLhjYopW4HbgfIyMhw0sTwxt5M61hZNX++vB/XDc+QZlqCIDSIM4LekHo0mBenlLoByAHGNbRdaz0DmAGQk5MjuXVNsPNIKRl1zbT+evVAuiTFkNqutb/NEgQhgHFmUXQPkO7wuDOwr/5OSqnJwCPApVrrSmvMCz+qa23845vtTHlmIW8u3QXAOd2TRMwFQWgWZzz0lUC2UqorsBeYBlznuINSajDwEjBVa33IcivDhPV7jvHrmevZcuAElwxM5dJB0kxLEATnaVbQtdY1Sqm7gXlAJPCa1nqTUupPwCqt9Wzgr0Ab4D/KjH8r0Fpf6kW7Q47XFu/kz19sJqVtS16+KYfz+nTwt0mCIAQZThUWaa3nAHPqPfc7h/uTLbYrbLA30xrQOZ4fDEvnwQt6E99aUhEFQXAdqRT1Eycqqnniyy20jIrkd5f0ISczkZxMaaYlCIL7SHMuP/DtlkOc/8xC3ltRQFSkkmZagiBYgnjoPqS4tIo/fbaJT9buo0eHNvzr+lEMzpBmWoIgWIMIug8pKa/mm+8P8fNJ2dw1IYvoKPmCJAiCdYige5kDJRV8snYvPxnbja7JsSx+cKIsegqC4BVE0L2E1pr3Vxby+BffU22zMbVvRzKTY0XMBUHwGiLoXmB3USkPztpAbn4RI7sl8sSVA8iUZlqCIHgZEXSLqam1cd3Lyykpr+bxK/ozbVi6NNMSBMEniKBbxI7DJ+lS10zrqWtNM61O8dJ/RRAE3yFpFh5SVWPj2a+3MfXZhbyVuxuAkd2SRMwFQfA54qF7wNrCYzwwcz1bD57gskGpXD44zd8mCYIQxoigu8mri3fy2Bebad+2Fa/enMOk3tJMSxAE/yKC7iL2ZlqD0uOZNjyDBy/oRVwrSUUUBMH/iKA7yfGKav4yZwutWkTw+0v6MrRLIkO7SDMtQRACB1kUdYKvNx/kvKcX8MHKAqKjIqSZliAIAYl46E1QdLKSP362mdnr9tGrY1tm3JjDwPR2/jZLEAShQUTQm+BERQ3fbj3EvZN78NPx3aWZliAIAY0Iej32HSvn4+/2cuf47mQmx7LkwYmy6CkIQlAggl6HzaZ5d0UBT3y5hVqb5qL+nchMjhUxFwQhaBBBB3YeKeXBWetZvrOY0VlJ/OWKAWQkxfjbLEEQBJcIe0GvqbVxwyvLOV5RzZNXDeCanM4oJc20BEEIPsJW0PMOnSAzKZaoyAie+cEguiTF0CGulb/NEgRBcJuwS9uorKnl6f9uY+qzi3izrpnW8K6JIuaCIAQ9YeWhryk4ygMz17P90EmuHJzGldJMSxCEECJsBP3lhfk8/uX3dIprxes/HMaEnu39bZIgCIKlhLyg22yaiAjFkC7tuH5EBg9M7UVbSUUUBCEECVlBLymv5rEvNtO6RSR/vKyfNNMSBCHkCclF0XmbDnDe0wuYtWYvsS2jpJmWIAhhQUh56EdOVvL7TzfxxYb99OkUx2u3DKNfWry/zRIEQfAJISXoJytqWLT9MPdP6cntY7vRIjIkv4AIgiA0SNAL+t5j5Xy8Zg93TcgiMzmWpQ9Nok3LoL8sQRAEl3HKhVVKTVVKbVVK5SmlHmxge0ul1Ad125crpTKtNrQ+Npvm7dxdnP/0Ap7/dge7i8oARMwFQQhbmlU/pVQk8DxwHrAHWKmUmq213uyw263AUa11llJqGvB/wA+8YTDAjsMneWjWBlbsKmZMdjKPX9Gf9ERppiUIQnjjjDs7HMjTWucDKKXeBy4DHAX9MuAPdfdnAv9USinthfSSmlobN726ghMV1fz16gFcPVSaaQmCIIBzgp4GFDo83gOMaGwfrXWNUqoESAKOOO6klLoduB0gIyPDPYMjI3h22iC6JMbQXvqvCIIgnMKZGHpD7m99z9uZfdBaz9Ba52itc1JSUpyxr0GGZSaKmAuCINTDGUHfA6Q7PO4M7GtsH6VUFBAPFFthoCAIguAczgj6SiBbKdVVKRUNTANm19tnNnBz3f2rgf95I34uCIIgNE6zMfS6mPjdwDwgEnhNa71JKfUnYJXWejbwKvC2UioP45lP86bRgiAIwtk4lbSttZ4DzKn33O8c7lcA11hrmiAIguAKUhsvCIIQIoigC4IghAgi6IIgCCGCCLogCEKIoPyVXaiUOgzsdvPlydSrQg0D5JrDA7nm8MCTa+6itW6wMtNvgu4JSqlVWuscf9vhS+SawwO55vDAW9csIRdBEIQQQQRdEAQhRAhWQZ/hbwP8gFxzeCDXHB545ZqDMoYuCIIgnE2weuiCIAhCPUTQBUEQQoSAFvRAHE7tbZy45vuUUpuVUuuVUt8opbr4w04rae6aHfa7WimllVJBn+LmzDUrpa6te683KaXe9bWNVuPE33aGUupbpdR3dX/fF/rDTqtQSr2mlDqklNrYyHallPp73e9jvVJqiMcn1VoH5A+mVe8OoBsQDawD+tTb507gxbr704AP/G23D655AhBTd/+n4XDNdfu1BRYCy4Acf9vtg/c5G/gOSKh73N7fdvvgmmcAP6273wfY5W+7PbzmscAQYGMj2y8EvsRMfBsJLPf0nIHsoZ8aTq21rgLsw6kduQx4s+7+TGCSCu6J0c1es9b6W611Wd3DZZgJUsGMM+8zwKPAk0CFL43zEs5c84+B57XWRwG01od8bKPVOHPNGoirux/P2ZPRggqt9UKantx2GfCWNiwD2imlOnlyzkAW9IaGU6c1to/WugawD6cOVpy5ZkduxXzCBzPNXrNSajCQrrX+3JeGeRFn3uceQA+l1BKl1DKl1FSfWecdnLnmPwA3KKX2YOYv/Mw3pvkNV//fm8WpARd+wrLh1EGE09ejlLoByAHGedUi79PkNSulIoBngFt8ZZAPcOZ9jsKEXcZjvoUtUkr101of87Jt3sKZa54OvKG1fkopdQ5mClo/rbXN++b5Bcv1K5A99HAcTu3MNaOUmgw8Alyqta70kW3eorlrbgv0A+YrpXZhYo2zg3xh1Nm/7U+11tVa653AVozAByvOXPOtwIcAWutcoBWmiVWo4tT/uysEsqCH43DqZq+5LvzwEkbMgz2uCs1cs9a6RGudrLXO1FpnYtYNLtVar/KPuZbgzN/2J5gFcJRSyZgQTL5PrbQWZ665AJgEoJTqjRH0wz610rfMBm6qy3YZCZRorfd7dER/rwQ3s0p8IbANszr+SN1zf8L8Q4N5w/8D5AErgG7+ttkH1/w1cBBYW/cz2982e/ua6+07nyDPcnHyfVbA08BmYAMwzd82++Ca+wBLMBkwa4Hz/W2zh9f7HrAfqMZ447cCdwB3OLzHz9f9PjZY8Xctpf+CIAghQiCHXARBEAQXEEEXBEEIEUTQBUEQQgQRdEEQhBBBBF0QBCFEEEEXBEEIEUTQBUEQQoT/BzAFfbO8VVTlAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plotReliability(preds, labels, n_bins=20)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**5. There are many ways to do training. Take us through how you do it (e.g. \"We used early stopping and stopped when validation loss increased twice in a row.\").**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The output layer was removed from the original model. The output of the last hidden layer which consist of 1024 values were taken as codings of images. Using the coding information, a support vector machine was trained on training set and evaluated on test set. The goal is that we want this SVM to be well calibratred, using `CalibratedClassifierCV` function by `sklearn`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████████████████████████████████████████████████████████████████████████| 2085/2085 [01:21<00:00, 25.65it/s]"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(2085, 1024)\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "# get output\n",
    "from tensorflow.keras.models import clone_model, Sequential\n",
    "from keras import backend as K\n",
    "getOutput = K.function(resmodel.layers[0].input, resmodel.layers[-2].output)\n",
    "features = getFeatrues(scaled_imgs)\n",
    "# Look at the shape of features extracted\n",
    "print(features.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Train a classifier and calibrate it\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.calibration import CalibratedClassifierCV\n",
    "clf = SVC()\n",
    "x_train, x_test, y_train, y_test = train_test_split(features, labels, test_size=0.2, random_state=9001)\n",
    "calibrator = CalibratedClassifierCV(clf, method = 'sigmoid', cv=5)\n",
    "_ = calibrator.fit(x_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get predict probability from SVM\n",
    "probs = calibrator.predict_proba(x_test)[:, 1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.7170263788968825"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accuracy_score(y_test, probs > 0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**6. Make a figure displaying your results.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd3hUVf7H8fc3DQg9BJCSEHoRETUCgtJFUBRBVMC662/RVVdX14J1Lavr2lfXRVFZu6g0gyBWmgJSFANEUHoo0mt6Ob8/JmjEQCZkJlPyeT1PHuZm7sz9XpJ8cnLuueeYcw4REQl9EYEuQEREfEOBLiISJhToIiJhQoEuIhImFOgiImEiKlAHjo+Pd0lJSYE6vIhISFq6dOku51z9kp4LWKAnJSWxZMmSQB1eRCQkmdnGoz2nLhcRkTChQBcRCRMKdBGRMKFAFxEJEwp0EZEwUWqgm9l4M9thZiuO8ryZ2XNmtsbMUs3sVN+XKSIipfGmhf4aMPAYzw8CWhd9jAbGlr8sEREpq1ID3Tk3F9hzjF2GAG84j4VAHTNr5KsCRUTCRWZuPul7Mv32/r7oQ28CpBfb3lz0ud8xs9FmtsTMluzcudMHhxYRCQ3z1+xi4LPzuO6tpRQW+mcdCl8EupXwuRKrdc6Nc84lO+eS69cv8c5VEZGwsj8rjzGTUhn1yjdEGNw3uAMRESXFZvn54tb/zUBCse2mwFYfvK+ISEgrKHRcNHY+63Ye4tpeLbilfxuqRkf67Xi+CPQU4EYzmwB0BfY757b54H1FRELS3oxc6sRGExlh3DagLY3rVKVT0zp+P26pgW5m7wK9gXgz2wz8HYgGcM69CMwAzgXWAJnAH/xVrIhIMHPOMXXZFh6clsadA9sxsksiAzueUGHHLzXQnXMjS3neATf4rCIRkRC0dV8W90xZzqzVOzklsQ7JzepWeA0Bmz5XRCRcfLhsC/dMWUFBoeP+wR24qnsSkX668HksCnQRkXKqXS2azgl1+Oewk0iIiw1YHQp0EZEyyi8o5NWv1pNXUMiNfVvTu20DerWpj1nFt8qLU6CLiJRB2tYD3DkpleVb9nNep0Y45zCzgIc5KNBFRLySk1/Af75cw9jZa6kTG81/LzuVQR1PCIogP0yBLiICkL4INsyDpLMgocvvnt6wK5MX56zlgs6Nue+8DtStHhOAIo9NgS4ikr4IXj8f8nMgqipclQIJXcjIyeeztO1ceEoT2p5Qky9u7U1ivcBd9CyNAl1EKrfCApj/PORne7YLcmDDPOZlN+euycvZsi+Ljk1q0apBzaAOc9CKRSJSma39El7qCT+kgBlYBC4yhv+sP4ErXl1ETGQE740+g1YNaga6Uq+ohS4ilc/2lfDpfbD2C6jTDIaPh1pNKdzwFTcvrMGMVXW5vncLburX2q+TafmaAl1EKo8D22DWI7DsbahSCwY8Al3+xJ4co061aCISuzK43s9cW6caHZvUDnS1ZaZAF5Hwl3MI5j/n6SsvyINu18NZf8NVq8vkb7fw0EeeybRGdU3knBMrbjItX1Ogi0j4KsiH796EWY9Cxg44cRj0ux/imrN5byZ3T1jM3B93clqzunRpHhfoastNgS4i4cc5+OlT+Ox+2LkKEs+Ake9C02QApny3mXunrMABD15wIld0a+a3VYQqkgJdRMJH+iJIfR82L4ZtyyCuBVz6FrQb7BnFUiSuehVOS4rj0aEdaVo3uIciloUCXUTCw7o58OZQcAWe7W43QP8HICqGvIJCXp63lvwCx039WtOrTX16to4Pqtv2fUGBLiLh4bP7fg1zi4Tq9SAqhhVb9nPnpFRWbj3A+Sc3DqrJtHxNgS4ioe+7t2Hb9xAR5ek/j4whp2l3/j1zFS/NXUfd2BhevPxUBnZsFOhK/UqBLiKhbccqmHGbZ1KtPvfApvmQdBYbotvz8rx5DDulCfee14HasdGBrtTvFOgiErpyM+CDqyCmOlz0Chkx8XyyO4FhCU1pC3z5t94BXUGooinQRSR0zbgDdq6GKyYzZ1skd0+ey9b9WXRqWptWDWpWqjAHTc4lIqFq2Tuw7C2yzriVW5fGcdX4RVSNjuCDa0NnMi1fUwtdRELPjlUw/W+4Zj04f3kP1u/Zyo19WnFj31YhNZmWrynQRSS05GaS/95VREbHYhe9yh3pRpO61TixcehNpuVrCnQRCRnOOda/eT1Ju1fzZfJY+tVqxIATA11V8FAfuoiEhPQ9mbz0/CO0SJ/ClOojaN71/ECXFHTUQheRoDf5282MnzqT9+05tsedxtAbniciKvzHlZeVAl1Egl6DqoW8WOV5qkTVIPYPb4PCvEQKdBEJOnkFhbw0Zy0FhXBz/9ac+dPjkLcRLp0EtcL79v3yUKCLSFBZsWU/t09M5YdtBxjSuTFu2bvYd2/CWbdBq36BLi+oKdBFJChk5xXw7Oc/8fK8dcRVj+GlK07jnAb7Ydyt0KwH9L4r0CUGPa9GuZjZQDNbbWZrzGxMCc8nmtksM/vOzFLN7Fzflyoi4WzTnkxe/Wodw09tyue39OKc1rXgg6shuhpc9ApEqv1ZmlL/h8wsEngBOBvYDCw2sxTnXFqx3e4F3nfOjTWzDsAMIMkP9YpIGDmYncfMFT9zcXICbRrWZNZtvX9dQSjlVtiRBpdPglqNA1toiPDmV14XYI1zbh2AmU0AhgDFA90BtYoe1wa2+rJIEQk/s1bt4J4py/n5QDanJNahVYOav4b59+/Bt2/AWX+DVv0DW2gI8SbQmwDpxbY3A12P2OcB4FMz+wtQHSjxK2Bmo4HRAImJiWWtVUTCwJ6MXB7+KI0p322hdYMaTPxz999OprXzR/joFkjsDr3vDlyhIcibQC9pnSZ3xPZI4DXn3FNmdgbwppl1dM4V/uZFzo0DxgEkJycf+R4iEuYKCh3Dx85n055MburXmhv6tKRKVLHJtHIzPfObR1eF4a+q37yMvPnf2gwkFNtuyu+7VK4BBgI45xaYWVUgHtjhiyJFJLTtPJhDveoxREYYd5/bniZ1q9G+Ua3f7zjzTk+/+WXqNz8e3oxyWQy0NrPmZhYDjABSjthnE9APwMzaA1WBnb4sVERCj3OO9xZvou9Ts3ln0SYA+ndoWHKYp77v6Tc/81ZorX7z41FqC905l29mNwKfAJHAeOfcSjN7CFjinEsB/ga8bGa34OmOudo5py4VkUps0+5MxkxOZf7a3XRtHseZreKPvvOun2DaXyHxDM+6oHJcvOqgcs7NwDMUsfjn7i/2OA3o4dvSRCRUTVy6mfumriAywnhkaEdGnp5IRERJl+OAvCx4v6jf/CL1m5eH/udExOca1qpC95b1+MfQjjSqXe3YO398J+xY6ek3r92kYgoMUwp0ESm33PxCxs5eS6Fz3HJ2G85qXZ+zWtcv/YWpH8C3r8OZt6jf3AcU6CJSLt+n7+OOiams3n6QYac0wTmH2VG6V4rb9RN89FdI6AZ97vV/oZWAAl1EjktWbgFPf7aaV79aT4OaVXnlymT6d2jo3YvzsjzztETGwPDx6jf3Ef0vishxSd+byevzNzKiSyJjBrWjVtUyLDoxcwxsXwGXTVS/uQ8p0EXEaweKJtO6pGgyrdm396ZxnVIueh5p+URY+hr0+Cu0PtsvdVZWCnQR8cqXq7Zz9+QV7DiYzamJdWnVoEbZw3zXGph2s6ffvK/6zX1NgS4ix7T7UA4PfZTGh8u20rZhTV684jRaNahR9jf6Tb/5qxCpdUF9TYEuIkdVUOi4+MUFpO/N5Jb+bfhz75bERHm1Ls7vzbwLti+HUR9A7aa+LVQABbqIlGDHwWziq1chMsK457z2NK0bS9sTapb+wqNZPhGW/g963AxtBviuUPmN4/xVKyLhqLDQ8fY3G+n75BzeLppMq1/7huUL81/6zbtC3/t8VKmURC10EQFgw64MxkxOZeG6PXRvWY9e3tzpWZq87KJ+8+ii8ebqN/cnBbqI8P6SdO6buoKYyAgeG3YSl56e4N3dnqX55HC/+fvqN68ACnQRoUmdavRsU5+Hh3TkhNpVy/+G6Ytg4VhYORm63wRtzin/e0qpFOgilVBOfgH/nbUW5xy3DmhLj1bx9DjWfOVlkb4IXh8M+TlgEdB2kG/eV0qli6Iilcx3m/Zy/vNf8e8vfmLLvmx8vhbN+nmeMAfAYNMC376/HJVa6CKVRGZuPk99+iPjv17PCbWqMv7qZPq283IyrbIoyC16YJ6biJLO8v0xpEQKdJFKYsveLN5cuJHLuiZy58B21CzLZFreytwDi1+B+HbQ6WJo3hMSuvj+OFIiBbpIGNuflcfHy7cxoksirRvWZM7tvUtfQag8Pv87ZO2FK6fCCSf57zhSIgW6SJj6dOXP3Dt1BbszcklOiqNVgxr+DfON8+HbNzyjWhTmAaFAFwkzuw7l8EDKSj5K3Ua7E2ryylXJxzeZVlnk58K0v0LtROg9xr/HkqNSoIuEkYJCx/Cx89m6L5vbBrTh2l4tiY6sgMFs8/8Nu1bDyPcgprr/jyclUqCLhIHtB7KpX8Mzmdbfzz+RpnWr0bphOeZfKYvda2HOE9BhCLQdWDHHlBJpHLpICCssdLy5cCP9nprD299sBKBPuwYVF+bOwfRbPcMTB/6rYo4pR6UWukiIWrfzEGMmL2fR+j2c2Sqe3m0bVHwRyz+AdbPh3CehVqOKP778hgJdJAS9t3gT93+4kipRETw+vBMXn9bUN5NplUXmHs+iFU1Og+Q/VuyxpUQKdJEQ1LRuLL3beibTalDLB5NpHY/PH/CMOb9iCkREBqYG+Q0FukgIyMkv4Pkv1gBw2zk+nkzreGxcAN++Dt3/Ao06Ba4O+Q0FukiQW7pxD3dMTGXtzgwuSW6Kc67iu1eKy8+Fj/4KtROg912Bq0N+R4EuEqQycvJ54pPVvL5gA41rV+P1P3ahVxsfrCJUXvOfg52rNOY8CHk1bNHMBprZajNbY2Yl3gZmZpeYWZqZrTSzd3xbpkjls3VfFu8s2sSV3ZrxyS09gyPMd6+FuRpzHqxKbaGbWSTwAnA2sBlYbGYpzrm0Yvu0Bu4Cejjn9ppZAMZPiYS+/Zl5TF++jVFdPZNpzbujDw0DddHzSM7B9L9BRLTGnAcpb7pcugBrnHPrAMxsAjAESCu2z5+AF5xzewGcczt8XahIuJu54mfu+3AFezJy6doijpb1awRPmAMsnwjrZmnMeRDzpsulCZBebHtz0eeKawO0MbOvzWyhmZX4t5iZjTazJWa2ZOfOncdXsUiY2XEwm+vfXsp1by2lfo0qfHhDD1rW9/NkWmWVucez4LPGnAc1b1roJV1OP3LNqiigNdAbaArMM7OOzrl9v3mRc+OAcQDJyck+XvdKJPQUFDoueXEBW/dnc/s5bRnds0XFTKZVVp8/4An1yydrzHkQ8ybQNwMJxbabAltL2Gehcy4PWG9mq/EE/GKfVCkSZrbtz6JhzaqeybQuOJGEurH+n+L2eGnMecjwpimwGGhtZs3NLAYYAaQcsc9UoA+AmcXj6YJZ58tCRcJBYaHjta/X0++pObx1eDKttg2CN8w15jyklNpCd87lm9mNwCdAJDDeObfSzB4CljjnUoqeG2BmaUABcLtzbrc/CxcJNWt2HGLMpFSWbNxLzzb16dsuBAaDacx5SDHnAtOVnZyc7JYsWRKQY0uYSl8EG+Z5VpkPsoWJJyzaxP0pK6kWHcn9gzsw7NQmgb3b0xu718LY7tDmHLjkjUBXI0XMbKlzLrmk53SnqISHHz+Fdy/1jJWOqgpXpQRVqCfWi6V/+wY8eEFH6tesEuhySqcx5yFJgS6hLz8XZtwGrrBoOwtWTA5ooGfnFfDcFz8BcMfAdnRvGU/3lgGcTKusDo85H/SExpyHkCAcHyVSRjPvhH0bITKaX76lF78Mc5+EgvwKL2fJhj2c+9w8/jt7LXsycglUt+ZxOzzmvPGpcPo1ga5GykAtdAlti1+BJeOhx83QbrCnD73hSfD9O/Dlw/DDNLjwv9DwRL+XcignnydmruKNhRtpUqcab/yxCz2DYf6VstKY85ClQJfQteEr+PhOaD0A+v3dEz6Hu1naDIAOF3r6gV/qBb3ugDNvKWrF+8fP+7OYsDidq85I4vZz2lK9Sgj+eB0ec37GjRpzHoI0ykVC096N8HIfqBYHf/oCqtYueb+M3fDxHbBiIpxwEgz5r0+Dam9GLh8t38YV3ZoBsONAduBWECqv/Fx46SzIzYDrF0KVIB0bX8kda5SL+tAl9OQcggmjPP3jIyccPcwBqteD4a/CpW/Dwe2eXwKz/ukJr3JwzjFj+TbOfmYOD6asZO3OQwChG+bw65jzc59UmIcoBbqElsJCmHod7EiDi8dDfCvvXtd+MNzwDXS8COY85gn2rcuOq4QdB7K57q2lXP/2tzSqXY2UG88Mvsm0ymrPOs885+0v0DznIUyBLqFl7uOeC51nPwSt+pfttbFxMGycp1WfsQte7gtfPAz5OV6/RUGh4+KXFjB79U7uGtSOKdd3p0PjWmU8iSBTfMz5II05D2UheNVGKq20FJj9Tzh5pOei3fFqOwgSu8En98C8J2HVdM9ImCanHvUlW/dlcUItz2RaDw3pSELdarQI9Vb5YSsmwdovi8acNw50NVIOaqFLaPh5BUy5Fpokw+Bnoby3zVer6wnxUR9A9n54pb9nuF5e9m92Kyh0/O+IybR6takfPmGetRdmjtGY8zChQJfgl7EbJoyEKrXg0rcg2ocXHtsMgBsWQudR8NUz8FJPSPfM+rxmx0EufnE+D05Lo2uLOPq1b+i74waLw2POz/+3xpyHAXW5SHAryIMPrvKMUPnDx/65Db1qbRjyHzjxQki5GcYPIC3pSi79sS9RVarxzKUnc2HnEJhMq6w2LYSlr2nMeRhRoEtwmznGc/fn0HHQ9DT/HqtVf7h+AXx2Hx2WvsbnsV9QZfhY6rRt6t/jBkJ+LkzTPOfhRl0uEryWjPfc2t/9L3DypX49VHZeAf/8+Acem7XV0/1wxVQaxhp13j0fZt4NuZl+PX6FW/A87PxBY87DjAJdgtOGr2HG7Z5Wc/8H/Xqob9btZtC/5/HSnHUczM7zTKbVsg9cP99zoXDhC/BiD9g43691VJg962DO4xpzHoYU6BJ89m2C96+Aus3holf9drHuYHYe905dzqXjFlJQ6Hjn/7ryyNCTfu0rr1ITznsKrpoGhQXwv3M9c8fkZvilngqhMedhTX3oElxyM+DdkUW39b8L1er47VDbD+Qwcelm/u/M5tw6oA2xMUf5cWjeE/48H754CL55EX6cCd1ugNyDQbk60jFpzHlY0+RcEjyc84xo+WEajHofWp/t80PsychleupWrjgjCYCdB3PKtoLQhq9h0jVwcJtnOyIazrwVGnf2tOir1vL8W6Xo36ggWp0oay/853TPhdD/+1zDFEOUlqCT0DD3CUj7EM5+2Odh7pzjo9RtPJCykgPZefRoFU+L+jXKvhxcUg847WqY/RjgoDAP5h6j6yIypijgD4d8rWLbxT6q1i7h88V+OURVKf/NVJ8/AJm74fJJCvMwpUCX4PDDNJj1CHS61DOqxYe2H8jmnikr+PyH7XRqWpu3h3ct352eLfvCV89CQa5nfvULx0Jcc8g+ADkHi30cKPnxgc2/Ps4+4PmlUJqI6N8HfdWSfjkc8RfC4cc/zvSMOe94ETQ6+fjPXYKaAl0Cb/tKmHwtNDkNzn+u/C3RYgoKHZe8tICf92dzz7nt+UOPJKIiyzkWIKGLZxHqDfN804eel13CL4DivwhK+vxBOLD1t9sFXkwytmo6pC8KrX5/8ZoCXQIrY7fnImiVmp45y310W//mvZk0ql2NyAjj4SEdSYyLJSm+uk/eG/AEoq9CMbqq56NGOZery8/xzBV/5C+B1Pdg5VTAee683TBPgR6mFOgSOL/c1v8z/GGGT27rPzyZ1pOfruauQe25qntSaK7reTyiqng+qtf77edj68HqmUVdRDGevyokLCnQJXBm3uVpLV74IjQt8aJ9maz++SB3TErl+/R99GvXgAEnhuFkWsfD111EErQU6BIYS/4Hi1/2TAzVeWS53+6thRt5cNpKalaN5t8jOnPByY3DbzKt8vBlF5EELQW6VLyN82HGbdCyn2floXJwzmFmtGpQg3NPasT9gztQr0YQjf0WqUAKdKlY+zbBe1dA3SQYPv64x0Nn5Rbw9GeriYgw7hrUnm4t6tGtRb3SXygSxjSXi1Sc3AyYMMpzcW7khOO+rX/B2t0M/PdcXp63nsycAgJ1t7NIsFELXSqGczD1es9ScqPeh/jWZX6LA9l5/HPGKt5dtIlm9WJ5509d6d4y3g/FioQmBbpUjLlPQtpUT595mwHH9RY7DuQw9bstjO7Zglv6t6FajG5fFynOqy4XMxtoZqvNbI2ZjTnGfsPNzJlZ+cegSfhYNR1m/aPotv6byvTS3YdyeO3r9QC0alCDr+7sw93ntleYi5Sg1Ba6mUUCLwBnA5uBxWaW4pxLO2K/msBNwDf+KFRC1PY0mDwaGp/iWQnIy6GEzjlSvt/KAykrOZSTT8829WlRv4ZGsIgcgzct9C7AGufcOudcLjABGFLCfg8DjwPZPqxPQlnmHnh3BMRUhxHvQHQ1r162dV8W17y+hJsnLKNZvepMv+ms8k2mJVJJeNOH3gRIL7a9GehafAczOwVIcM59ZGa3He2NzGw0MBogMTGx7NVK6Pjltv5tcPUMrxdTyC8oZMS4hew8mMN9gztwdfckIiN0g5CIN7wJ9JJ+mn4ZJ2ZmEcAzwNWlvZFzbhwwDjwLXHhXooSkT+6B9XM9U8smnF7q7ul7MmlcpxpRkRE8OvQkEuNiSawXWwGFioQPb7pcNgMJxbabAluLbdcEOgKzzWwD0A1I0YXRSmzp67DopaLb+kcdc9f8gkLGzV1L/6fn8OaCDQCc2TpeYS5yHLxpoS8GWptZc2ALMAL45afUObcf+GUwsJnNBm5zzml9ucpo4wLPIsQt+0L/B4+56w/bDnDnpFRSN+/n7A4NGXRS+WdbFKnMSg1051y+md0IfAJEAuOdcyvN7CFgiXMuxd9FSojYlw7vXwF1Ej239Uce/dvrzQUbeHBaGrWrRfOfUadw3kmNNJmWSDl5dWORc24GMOOIz91/lH17l78sCTnr5sLkP3lu7796OlSrW+JuhyfTatOwJuef3Jj7BncgrnpMBRcrEp50p6iU38aF8OYQcIWeBRSy9/9ul8zcfJ785EeiIo27z21P1xb16KrJtER8SpNzSfl9+bAnzAEKCzwLKRTz9ZpdnPPsXMZ/vZ7c/EJNpiXiJ2qhS/l89xZs/Aqs6Fb8Ykuc7c/K49HpP/DeknSax1fn/WvPoEvzuAAWKxLeFOhy/NbPhWk3Q4s+0OsO2LTgN0uc7TqUw7TUrVzXqyV/7d+aqtGaf0XEnyxQf/4mJye7JUs0sjFk7foJXukHNRvBHz/5ZW7znQdzmPb9Vv54ZnMA9mTk6qKniA+Z2VLnXIn3+aiFLmWXsRvevhgiomHUe1CtDs45pi7bwoPT0sjMKaBPuwY0j6+uMBepQAp0KZv8HHjvMjiw1TM8sW4SW/Zlcc+U5cxevZNTE+vw+PBONI+vHuhKRSodBbp4zzlI+Yunr3z4eEg4vWgyrQXsPpTLA+d34IozNJmWSKAo0MV7cx6H1Peg771sajSIJoWOqMgIHhvWicS4WBLiNP+KSCBpHLp4J/UDmP0ohZ1GMrZgKP2fmcMbCzYA0KNVvMJcJAiohS6l27QQPryejEbduGzzJSxbtJpzTmzIeZpMSySoKNDl2Pasgwmj2F/lBPpt/CPEFjL2slM1M6JIEFKgy9Fl7cW9cynmCll/zv/o9WM09w1uT51YDUUUCUYKdClRRmYW28deRLND64i8KoXOScl07hzoqkTkWHRRVH5n7uodfPnU5bQ4uJSUxLtwzboHuiQR8YJa6PKL/Zl5PDw9jXrLxnJX9Ods7XQjQ4cddc1vEQkyaqHLL3Zl5JC3fCp3Rb9LQYehNL7w4UCXJCJloECv5HYczOaVeesAaJm7mmdjxkLT04kcOhYi9O0hEkrU5VJJOeeY9O0WHv4ojay8AgY0zSNx0kisRn0Y8S5EVwt0iSJSRgr0Sih9TyZ3T1nOvJ92kdysLv86vzmJKRdBXhZcmQI16ge6RBE5Dgr0Sia/oJCRLy9kb0YuDw85kctOb0LEhJGwcxVcPhEatAt0iSJynBTolcSGXRkkxMUSFRnB48M9k2k1rRsLM26HNZ/B4GehZd9Alyki5aCrXmEur6CQF2atYcAzc3+ZTKt7y3hPmH/zEiwaB93/Asl/CGidIlJ+aqGHsRVb9nPHxFTSth3gvJMaMbhT41+f/Pp5+Ow+aNYD+j8YuCJFxGcU6GHqf1+v5x/TfyCuegwvXn4aAzue8OuTK6fAZ/d6Hm9Z6vkoWthZREKXAj3MOOcwM05sXJthpzTh3vM6UDs2+rc77V4LGOCgIA82zFOgi4QBBXqYOJSTz+MzVxETGcG9gzvQpXkcXZrHlbxz854QVRUKciEyBpLOqthiRcQvFOhhYPbqHdwzZQVb92fxxx7Nf2mlH1VCF7gqxdMyTzpLrXORMKFAD2F7M3J5eHoak7/dQqsGNZh4XXdOa1bXuxcndFGQi4QZBXoI25uZy6crt3NT31bc0LcVVaIiA12SiASQV+PQzWygma02szVmNqaE5281szQzSzWzL8ysme9LFYAdB7IZN3ctzjla1K/B13f25dYBbRXmIlJ6oJtZJPACMAjoAIw0sw5H7PYdkOyc6wRMBB73daGVnXOO9xen0+/pOTz16Y9s2J0J8PsRLCJSaXnT5dIFWOOcWwdgZhOAIUDa4R2cc7OK7b8QuNyXRVZ26XsyuWvycr5as4suzeN4bNhJNI+vHuiyRCTIeBPoTYD0Ytubga7H2P8a4OOSnjCz0cBogMTERC9LrNwOT6a1LzOPf1zYkVFdEomIOMYIFhGptLwJ9JLSw5W4o9nlQDLQq6TnnXPjgHEAycnJJb6HeKzflUFi0WRaTww/mWb1YmlcR3OUi8jReXNRdKBmKcIAAAllSURBVDOQUGy7KbD1yJ3MrD9wD3CBcy7HN+VVPnkFhTz/xU+c88xcXp+/AYAzWtZTmItIqbxpoS8GWptZc2ALMAIYVXwHMzsFeAkY6Jzb4fMqK4nUzfu4Y2Iqq34+yPknN+aCzo1Lf5GISJFSA905l29mNwKfAJHAeOfcSjN7CFjinEsBngBqAB8U3aG4yTl3gR/rDjvjv1rPP6anUb9mFV6+MpmzOzQMdEkiEmK8urHIOTcDmHHE5+4v9ri/j+uqNA7fpt+paW0uPT2BMYPaU7uahiKKSNnpTtEAOZidx2Mfr6JKVCT3n9+B5KQ4kpOOMpmWiIgXtGJRAMxatYMBz8zl3UWbiIo0nNOAHxEpP7XQK9CejFwemraSqcu20qZhDf57WXdOSfRyMi0RkVIo0CvQ/qw8vvhhBzf3a80NfVoRE6U/kETEdxTofvbz/mymLtvCtT1b0Dy+Ol+N6auLniLiFwp0P3HOMWFxOo9O/4G8wkIGnngCSfHVFeYi4jcKdD/YuDuDMZOWs2Ddbrq1iOOxYZ1I0mRaIuJnCnQfyy8oZNTL37A/K49Hh57EiNMTNJmWiFQIBbqPrN15iGZFk2k9dYlnMq1GtTX/iohUHA2zKKfc/EKe/fxHBj47lzcWbASgW4t6CnMRqXBqoZfDsvR93DkxldXbDzKkc2MuPKVJoEsSkUpMgX6cXv1qPY9MT6NBzaq8elUy/dprMi0RCSwFehkdnkyrc0JtRnRJZMygdtSqqqGIIhJ4CnQvHcjO458zVlE1OoK/n38ipzWL47RmmkxLRIKHLop64fO07Zz99BzeW7yJmKgITaYlIkFJLfRj2H0ohwenpZHy/VbanVCTcVckc3JCnUCXJSJSIgX6MRzMzmfW6h3c0r8Nf+7dUpNpiUhQU6AfYeu+LKZ8t4Xre7ckKb46X4/pq4ueIhISFOhFCgsd7yzaxGMfr6Kg0HHeSY1Iiq+uMBeRkKFAB9bvymDMpFS+Wb+HHq3q8c+hnUisFxvoskREyqTSB3p+QSGXv/INB7LzePyiTlyc3BQzTaYlIqGn0gb6mh0HSapXnajICJ65tDPN6sXSsFbVQJclInLcKt2wjZz8Ap7+7EcGPjuP14sm0+rSPE5hLiIhr1K10L/dtJc7J6by045DDDulCcM0mZaIhJFKE+gvz13Hox//QKNaVfnfH06nT9sGgS5JRMSnwj7QCwsdERHGqc3qcFnXRO4c2I6aGoooImEobAN9f1Yej0xPo1p0JA8O6ajJtEQk7IXlRdFPVv7M2U/PYdK3W6heJUqTaYlIpRBWLfRdh3L4+4crmb58Gx0a1WL81afTsUntQJclIlIhwirQD2XnM++nndx+TltG92xBdGRY/gEiIlKikA/0LfuymPLtZm7o04qk+OrMv6sfNaqE/GmJiJSZV01YMxtoZqvNbI2ZjSnh+Spm9l7R89+YWZKvCz1SYaHjzQUbGPD0HF6YtZaNuzMBFOYiUmmVmn5mFgm8AJwNbAYWm1mKcy6t2G7XAHudc63MbATwL+BSfxQMsHbnIe6atJxFG/ZwVut4Hh16EglxmkxLRCo3b5qzXYA1zrl1AGY2ARgCFA/0IcADRY8nAv8xM3N+GF6SX1DIla8u4mB2Hk8M78Tw0zSZlogIeBfoTYD0Ytubga5H28c5l29m+4F6wK7iO5nZaGA0QGJi4vEVHBnBsyM60ywulgaaf0VE5Bfe9KGX1Pw9suXtzT4458Y555Kdc8n169f3pr4SnZ4UpzAXETmCN4G+GUgott0U2Hq0fcwsCqgN7PFFgSIi4h1vAn0x0NrMmptZDDACSDlinxTgqqLHw4Ev/dF/LiIiR1dqH3pRn/iNwCdAJDDeObfSzB4CljjnUoBXgTfNbA2elvkIfxYtIiK/59WgbefcDGDGEZ+7v9jjbOBi35YmIiJloXvjRUTChAJdRCRMKNBFRMKEAl1EJExYoEYXmtlOYONxvjyeI+5CrQR0zpWDzrlyKM85N3POlXhnZsACvTzMbIlzLjnQdVQknXPloHOuHPx1zupyEREJEwp0EZEwEaqBPi7QBQSAzrly0DlXDn4555DsQxcRkd8L1Ra6iIgcQYEuIhImgjrQg3Fxan/z4pxvNbM0M0s1sy/MrFkg6vSl0s652H7DzcyZWcgPcfPmnM3skqKv9Uoze6eia/Q1L763E81slpl9V/T9fW4g6vQVMxtvZjvMbMVRnjcze67o/yPVzE4t90Gdc0H5gWeq3rVACyAG+B7ocMQ+1wMvFj0eAbwX6Lor4Jz7ALFFj/9cGc65aL+awFxgIZAc6Lor4OvcGvgOqFu03SDQdVfAOY8D/lz0uAOwIdB1l/OcewKnAiuO8vy5wMd4VnzrBnxT3mMGcwv9l8WpnXO5wOHFqYsbArxe9Hgi0M9Ce8XoUs/ZOTfLOZdZtLkQzwpSocybrzPAw8DjQHZFFucn3pzzn4AXnHN7AZxzOyq4Rl/z5pwdUKvocW1+vzJaSHHOzeXYK7cNAd5wHguBOmbWqDzHDOZAL2lx6iZH28c5lw8cXpw6VHlzzsVdg+c3fCgr9ZzN7BQgwTn3UUUW5kfefJ3bAG3M7GszW2hmAyusOv/w5pwfAC43s8141l/4S8WUFjBl/XkvlVcLXASIzxanDiFen4+ZXQ4kA738WpH/HfOczSwCeAa4uqIKqgDefJ2j8HS79MbzV9g8M+vonNvn59r8xZtzHgm85px7yszOwLMKWkfnXKH/ywsIn+dXMLfQK+Pi1N6cM2bWH7gHuMA5l1NBtflLaedcE+gIzDazDXj6GlNC/MKot9/bHzrn8pxz64HVeAI+VHlzztcA7wM45xYAVfFMYhWuvPp5L4tgDvTKuDh1qedc1P3wEp4wD/V+VSjlnJ1z+51z8c65JOdcEp7rBhc455YEplyf8OZ7eyqeC+CYWTyeLph1FVqlb3lzzpuAfgBm1h5PoO+s0CorVgpwZdFol27AfufctnK9Y6CvBJdylfhc4Ec8V8fvKfrcQ3h+oMHzBf8AWAMsAloEuuYKOOfPge3AsqKPlEDX7O9zPmLf2YT4KBcvv84GPA2kAcuBEYGuuQLOuQPwNZ4RMMuAAYGuuZzn+y6wDcjD0xq/BrgOuK7Y1/iFov+P5b74vtat/yIiYSKYu1xERKQMFOgiImFCgS4iEiYU6CIiYUKBLiISJhToIiJhQoEuIhIm/h8h6KDw8JLmXgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plots of calibrated SVM\n",
    "plotReliability(probs, y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What if we want the algorithm to exclude all the mammograms that score higher than 0.8. Let's compare the results before and after calibration."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There are 417 patients in 'y_test', and 96 patients have > 0.8 score.\n",
      "The excluded mammogram will have 0.8958 accuracy and reduce the workload of expert radiologists by 23.02%\n"
     ]
    }
   ],
   "source": [
    "print('There are {} patients in \\'y_test\\', and {} patients have > 0.8 score.'.format(len(y_test), np.sum(probs > 0.8)))\n",
    "print('The excluded mammogram will have {} accuracy and reduce the workload of expert radiologists by {}%'.format(round(np.mean(y_test[probs > 0.8]), 4), round(np.sum(probs > 0.8) / len(y_test) * 100, 2)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**7. Discuss pros and cons of your method and what you might have done differently now that you've tried or would try if you had more time.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This study did provide a potential feasible direction to explore which could help reduce the workload of expert radiologists so that they can spend more energy and time to focus on the more difficult mammograms to reduce misclassification finally.\n",
    "\n",
    "There are multiple things to try as next steps. Here I want to provide two wide fields which each of them contains a lot of potentials to explore.\n",
    "1. Use a better classifier at the first place or fine tune the model to yield better results. If this method is built on better classifiers, which could have more accurate predictions, then imaginably, the corresponding calibrated model would have more ability to exclude more images with high confidence.\n",
    "2. Incorporating `features_matrix`. In `features_matrix` there are scores which indicate the difficulty levels and confidence levels that radiologists can classify the mammogram. Build another CNN to predict those scores, and includes the predicted scores into the finall model of calibration, the results would be better. This is actually trying to predict how difficult human would think that this mammogram is to be classified."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**You will not be graded on the performance of your model. You'll only be graded on the scientific soundness of your claims, methodology, evaluation (i.e. fair but insightful statistics), and discussion of the strengths and shortcomings of what you tried. Feel free to reuse some of the code you are/will be using for your projects. The write-up doesn't need to be long (~1 page will suffice), but please cite at least one clinical paper and one technical paper (1 each in questions 1 and 4 at least, and more if needed).**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<span id='s1'> 1. Oeffinger, K. C., Fontham, E. T., Etzioni, R., Herzig, A., Michaelson, J. S., Shih, Y. C. T., ... & Wolf, A. M. (2015). Breast cancer screening for women at average risk: 2015 guideline update from the American Cancer Society. Jama, 314(15), 1599-1614.</span> \n",
    "\n",
    "<span id='s2'> 2. Wernick, M. N., Yang, Y., Brankov, J. G., Yourganov, G., & Strother, S. C. (2010). Machine learning in medical imaging. IEEE signal processing magazine, 27(4), 25-38.</span>\n",
    "\n",
    "<span id='s3'> 3. Rodriguez-Ruiz, A., Mordang, J. J., Karssemeijer, N., Sechopoulos, I., & Mann, R. M. (2018, July). Can radiologists improve their breast cancer detection in mammography when using a deep learning based computer system as decision support?. In 14th International Workshop on Breast Imaging (IWBI 2018) (Vol. 10718, p. 1071803). International Society for Optics and Photonics. </span>\n",
    "\n",
    "<span id='s4'> 4. Yi, D., Sawyer, R. L., Cohn III, D., Dunnmon, J., Lam, C., Xiao, X., & Rubin, D. (2017). Optimizing and visualizing deep learning for benign/malignant classification in breast tumors. arXiv preprint arXiv:1705.06362.</span>\n",
    "\n",
    "<span id='s5'> 5. Shen, L., Margolies, L. R., Rothstein, J. H., Fluder, E., McBride, R., & Sieh, W. (2019). Deep learning to improve breast cancer detection on screening mammography. Scientific reports, 9(1), 1-12. </span>"
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}