AIforHealth_DL1.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "2JWhrS6VjEXz"
      },
      "source": [
        "# **AI for Health Deep Learning workshop**\n",
        "\n",
        "Welcome to the AI for Health Deep Learning workshop! In the first part of this workshop, you will be training a convolutional neural network (CNN) to classify handwritten digits. \n",
        "\n",
        "In the second part, you will employ a CNN to determine the malignancy of cells in pathology images. Let's see if the neural network performs better than you!\n",
        "\n",
        "# How to use this notebook\n",
        "This notebook contains cells with snippets of code that help you to use to load, process and visualize data, to train CNNs and to visualize results. The code snippets are accompanied by snippets of text that explain what the code does. \n",
        "\n",
        "Some text cells also contain questions for you to answer and discuss. \n",
        "\n",
        "To execute a cell of code, move the cursor into the cell by clicking and press ctrl+Enter. Some code snippets build on results of previous cells, so make sure that you run all code cells preceding the current cell. \n",
        "\n",
        "Feel free to experiment by altering the code and to ask any questions you have!"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Z3DLTdf1h5V0"
      },
      "source": [
        "# MNIST: training a network to 'read' handwritten digits"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "etAMMkNYOzNs"
      },
      "source": [
        "## **Understanding the data** "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "NLsSgQN1dXp1"
      },
      "source": [
        "First we import the libraries necessary for the code in this notebook to run. There are different libraries available that make it easier to create and train deep learning networks without having to code things like back propagation or convolutions yourself. Here we use PyTorch. Another much used library is Tensorflow."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {
        "cellView": "both",
        "id": "eoC_ZUZ2FXOJ"
      },
      "outputs": [],
      "source": [
        "from __future__ import print_function\n",
        "\n",
        "import os\n",
        "os.environ['TF_CPP_MIN_LOG_LEVEL'] = '3'\n",
        "\n",
        "import torch\n",
        "import torch.optim as optim\n",
        "import torch.nn as nn\n",
        "import torch.nn.functional as F\n",
        "\n",
        "import torchvision\n",
        "import torchvision.transforms as transforms\n",
        "\n",
        "from matplotlib import pyplot as plt\n",
        "from matplotlib.colors import LogNorm\n",
        "\n",
        "import numpy as np\n",
        "from sklearn.metrics import confusion_matrix, ConfusionMatrixDisplay\n",
        "from math import sqrt, ceil\n",
        "import pandas as pd\n",
        "import seaborn as sn"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "62xfpmbzdW5H"
      },
      "source": [
        "MNIST is a famous dataset of images of handwritten digits. We use it as toy example. And the good thing is: Pytorch has this dataset available via torchvision!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {
        "id": "plGqX5_tFauC"
      },
      "outputs": [
        {
          "ename": "OSError",
          "evalue": "[Errno 30] Read-only file system: '/data'",
          "output_type": "error",
          "traceback": [
            "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
            "\u001b[0;31mOSError\u001b[0m                                   Traceback (most recent call last)",
            "\u001b[1;32m/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb Cell 7\u001b[0m line \u001b[0;36m7\n\u001b[1;32m      <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=0'>1</a>\u001b[0m transform \u001b[39m=\u001b[39m transforms\u001b[39m.\u001b[39mCompose(\n\u001b[1;32m      <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=1'>2</a>\u001b[0m     [transforms\u001b[39m.\u001b[39mToTensor()]\n\u001b[1;32m      <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=2'>3</a>\u001b[0m     )\n\u001b[1;32m      <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=4'>5</a>\u001b[0m batch_size \u001b[39m=\u001b[39m \u001b[39m128\u001b[39m\n\u001b[0;32m----> <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=6'>7</a>\u001b[0m mnist_train \u001b[39m=\u001b[39m torchvision\u001b[39m.\u001b[39mdatasets\u001b[39m.\u001b[39mMNIST(\u001b[39m'\u001b[39m\u001b[39m/data/mnist\u001b[39m\u001b[39m'\u001b[39m, download\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m, train\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m, transform \u001b[39m=\u001b[39m transform)\n\u001b[1;32m      <a href='vscode-notebook-cell:/Users/junda_huang/Work/radboudumc_projects/ai4h_course/AIforHealth_DL1.ipynb#W6sZmlsZQ%3D%3D?line=8'>9</a>\u001b[0m mnist_test \u001b[39m=\u001b[39m torchvision\u001b[39m.\u001b[39mdatasets\u001b[39m.\u001b[39mMNIST(\u001b[39m'\u001b[39m\u001b[39m/data/mnist\u001b[39m\u001b[39m'\u001b[39m, download\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m, train\u001b[39m=\u001b[39m\u001b[39mFalse\u001b[39;00m, transform \u001b[39m=\u001b[39m transform)\n",
            "File \u001b[0;32m~/anaconda3/envs/ai4h_dl/lib/python3.11/site-packages/torchvision/datasets/mnist.py:99\u001b[0m, in \u001b[0;36mMNIST.__init__\u001b[0;34m(self, root, train, transform, target_transform, download)\u001b[0m\n\u001b[1;32m     96\u001b[0m     \u001b[39mreturn\u001b[39;00m\n\u001b[1;32m     98\u001b[0m \u001b[39mif\u001b[39;00m download:\n\u001b[0;32m---> 99\u001b[0m     \u001b[39mself\u001b[39m\u001b[39m.\u001b[39mdownload()\n\u001b[1;32m    101\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39mnot\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_check_exists():\n\u001b[1;32m    102\u001b[0m     \u001b[39mraise\u001b[39;00m \u001b[39mRuntimeError\u001b[39;00m(\u001b[39m\"\u001b[39m\u001b[39mDataset not found. You can use download=True to download it\u001b[39m\u001b[39m\"\u001b[39m)\n",
            "File \u001b[0;32m~/anaconda3/envs/ai4h_dl/lib/python3.11/site-packages/torchvision/datasets/mnist.py:179\u001b[0m, in \u001b[0;36mMNIST.download\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    176\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_check_exists():\n\u001b[1;32m    177\u001b[0m     \u001b[39mreturn\u001b[39;00m\n\u001b[0;32m--> 179\u001b[0m os\u001b[39m.\u001b[39mmakedirs(\u001b[39mself\u001b[39m\u001b[39m.\u001b[39mraw_folder, exist_ok\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m)\n\u001b[1;32m    181\u001b[0m \u001b[39m# download files\u001b[39;00m\n\u001b[1;32m    182\u001b[0m \u001b[39mfor\u001b[39;00m filename, md5 \u001b[39min\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39mresources:\n",
            "File \u001b[0;32m<frozen os>:215\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n",
            "File \u001b[0;32m<frozen os>:215\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n",
            "File \u001b[0;32m<frozen os>:215\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n",
            "File \u001b[0;32m<frozen os>:225\u001b[0m, in \u001b[0;36mmakedirs\u001b[0;34m(name, mode, exist_ok)\u001b[0m\n",
            "\u001b[0;31mOSError\u001b[0m: [Errno 30] Read-only file system: '/data'"
          ]
        }
      ],
      "source": [
        "transform = transforms.Compose(\n",
        "    [transforms.ToTensor()]\n",
        "    )\n",
        "\n",
        "batch_size = 128\n",
        "\n",
        "mnist_train = torchvision.datasets.MNIST('/data/mnist', download=True, train=True, transform = transform)\n",
        "\n",
        "mnist_test = torchvision.datasets.MNIST('/data/mnist', download=True, train=False, transform = transform)\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "sUYcIuC2mWOM"
      },
      "outputs": [],
      "source": [
        "print('Training set input size:', mnist_train.data.shape)\n",
        "print('Training set output size:', mnist_train.targets.shape)\n",
        "print('Testing set input size:', mnist_test.data.shape)\n",
        "print('Testing set output size:', mnist_test.targets.shape)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "z9_7v160d7bq"
      },
      "source": [
        "*Exercise:* \n",
        " - *What do you think these numbers mean?* \n",
        " - *How many images does each set contain?* \n",
        " - *How large are these images (in pixels)?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dY1EEUAFqlka"
      },
      "source": [
        "We are going to split the train set into a smaller training set and a validation set. The validation set will be used during training to see intermediate performance. Also, we can change the model architecture or hyperparameters and check back with the validation performance what worked best. In this way the test set will be left untouched during the training and development process and will be a real final test to see how our model works."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "JcbRmWzmqk7S"
      },
      "outputs": [],
      "source": [
        "#from sklearn.model_selection import train_test_split\n",
        "#x_train, x_val, y_train, y_val = train_test_split(x_train, y_train, test_size=1/6, random_state=44, shuffle=True)\n",
        "\n",
        "train_set, val_set = torch.utils.data.random_split(mnist_train, [50000,10000])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "bPUwKcWtrqOU"
      },
      "outputs": [],
      "source": [
        "print('Training set input length:', len(train_set))\n",
        "print('Training set output length:', len(train_set))\n",
        "print('Validation set input length:', len(val_set))\n",
        "print('Validation set output length:', len(val_set))\n",
        "print('Testing set input length:', len(mnist_test.data))\n",
        "print('Testing set output length:', len(mnist_test.targets))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "G60SAVpuoL-l"
      },
      "source": [
        "*Exercise:*\n",
        "\n",
        "- *The following snippet of code prints the `i`th image in the train set. Change the value of index `i` below (between 0 and the size of the train set) to see different samples.*"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "3DcnimSdmvlE"
      },
      "outputs": [],
      "source": [
        "i=5\n",
        "plt.figure()\n",
        "plt.imshow(train_set[i][0][0])\n",
        "plt.title('label: {}'.format(train_set[i][1]));"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "oMJg-4DIPgop"
      },
      "source": [
        "## **Preparing the data**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "mmwyixIlsRcj"
      },
      "source": [
        "An image dataset that can be fed to a neural network usually has 4 dimensions.\n",
        "- The number of images\n",
        "- The number of color channels\n",
        "  - When using color images there are three channels for red, green and blue (RGB).\n",
        "  - When using black and white there is only 1 channel.\n",
        "- The number of pixel rows\n",
        "- The number of pixel columns\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "RKbbiVUjLyQL"
      },
      "source": [
        "*Exercise:*\n",
        "- *What should the dimensions of the train set be?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "mRKSgueeYkP8"
      },
      "source": [
        "Pytorch has very useful DataLoaders that can be used to feed data to the network when training, validating and testing. As we will see later you can loop over these DataLoader objects to feed the images and labels to your model during train, validation and testing. When creating these DataLoaders we already define the size of the batches that we feed into the model. A batch contains multiple images (e.g. 128) that are fed to the model at the same time. This helps us making more accurate updates to the weights of the model."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "-_zpkwhZ-o37"
      },
      "outputs": [],
      "source": [
        "# Create DataLoaders\n",
        "\n",
        "from torch.utils.data import DataLoader\n",
        "\n",
        "batch_size = 128\n",
        "\n",
        "trainloader = DataLoader(train_set, batch_size=batch_size,\n",
        "                                          shuffle=True, num_workers=2)\n",
        "\n",
        "valloader = DataLoader(val_set, batch_size=batch_size,\n",
        "                                          shuffle=False, num_workers=2)\n",
        "\n",
        "testloader = DataLoader(mnist_test, batch_size=batch_size,\n",
        "                                         shuffle=False, num_workers=2)\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MnboYbqbPqWP"
      },
      "source": [
        "## **Building a convolutional neural network**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Zp87Mn6yhAWX"
      },
      "source": [
        "Now that the data is ready, let's build a convolutional neural network (CNN)! In Pytorch, neural networks are class objects. In the '__init__' function you define all the layers that your network will have. In the 'forward' function the different calculations that happen when you feed an image in the network are defined. \n",
        "\n",
        "In the network below we use 2D convolution layers (Conv2d), max pooling layers (MaxPool2d) and linear or dense layers (Linear). Generally, every layer also has an activation function associatied with it. Here we use the Rectified Linear Unit (relu)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "XMLy_GwmLDMR"
      },
      "outputs": [],
      "source": [
        "import torch.nn as nn\n",
        "import torch.nn.functional as F\n",
        "from torchsummary import summary\n",
        "\n",
        "# Define model\n",
        "class Net(nn.Module):\n",
        "    def __init__(self):\n",
        "        super().__init__()\n",
        "        self.conv1 = nn.Conv2d(1, 6, 5)\n",
        "        self.pool = nn.MaxPool2d(2, 2)\n",
        "        self.conv2 = nn.Conv2d(6, 16, 5)\n",
        "        self.fc1 = nn.Linear(256, 120)\n",
        "        self.fc2 = nn.Linear(120, 84)\n",
        "        self.fc3 = nn.Linear(84, 10)\n",
        "\n",
        "    def forward(self, x):\n",
        "        x = self.pool(F.relu(self.conv1(x)))\n",
        "        x = self.pool(F.relu(self.conv2(x)))\n",
        "        x = torch.flatten(x, 1) # flatten all dimensions except batch\n",
        "        x = F.relu(self.fc1(x))\n",
        "        x = F.relu(self.fc2(x))\n",
        "        x = self.fc3(x)\n",
        "        return x\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "XS1PfM0NmGtM"
      },
      "source": [
        "*Exercise:*\n",
        "\n",
        "- Check out the above neural network. Do you understand the steps that happen in the 'forward' function?\n",
        "- Why is torch.flatten() needed, do you think?"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HqYN5oBzrjbr"
      },
      "source": [
        "Now, we can compile our network which makes it ready for training. While doing so, we state what loss metric and optimizer we want to use during training.\n",
        "- The optimizer is the algorithm we use to train the neural network. It determines what weigths have to be increased and decreased. We use stochastic gradient descent (SGD) as optimization algorithm.   \n",
        "- The loss is the metric that the optimizer uses to determine the performance of the model (we use categorical cross-entropy). The larger the loss, the worse the performance. Zero indicates perfect performance.\n",
        "\n",
        "To speed up training of neural networks, GPUs are very important. Where CPUs perform all operations in sequence, GPUs are able to do many similar computations in parallel. The computations that happen in neural networks are very well suited for parallel computations, especially the convolution computations used in many deep learning architectures.\n",
        "\n",
        "To make use of a GPU we have to explicitly send our model and our data to the GPU. Below we first detect whether there is a GPU available, and then initiate the model and send it to the GPU."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "9EvnadnnIeqC"
      },
      "outputs": [],
      "source": [
        "\n",
        "device = torch.device('cuda:0' if torch.cuda.is_available() else 'cpu')\n",
        "\n",
        "# Assuming that we are on a CUDA machine, this should print a CUDA device:\n",
        "\n",
        "print(device)\n",
        "\n",
        "# Initiate model and send model to device\n",
        "model = Net()\n",
        "model.to(device)\n",
        "\n",
        "# Define loss function and optimizer\n",
        "criterion = nn.CrossEntropyLoss()\n",
        "optimizer = optim.SGD(model.parameters(), lr=0.01, momentum=0.9)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DW7nTbS3P0DT"
      },
      "source": [
        "## **Training and testing a convolutional neural network**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TWhu7X64h8XM"
      },
      "source": [
        "Before we start training the network, let's see what happens if we let the current, randomly initialized model predict the digits for the validation set.\n",
        "Run the next 2 cells to see the accuracy. The first cell defines a useful 'predict_output' function that loops over the testset and calculates the predictions. In the second cell we calculate the accuracy between these predictions and the actual labels."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ZJhu8CqbHBRS"
      },
      "outputs": [],
      "source": [
        "# Define function to predict the digits for the test set\n",
        "def predict_output(testset):\n",
        "\n",
        "    x_test = []\n",
        "    y_test_pred = []\n",
        "\n",
        "    # since we're not training, we don't need to calculate the gradients for our outputs\n",
        "    with torch.no_grad():\n",
        "        for i, data in enumerate(testset, 0):\n",
        "            # get the inputs; data is a list of [inputs, labels]\n",
        "            inputs = data[0].to(device)\n",
        "            # calculate outputs by running images through the network\n",
        "            outputs = model(inputs)\n",
        "            # the class with the highest energy is what we choose as prediction\n",
        "            _, predicted = torch.max(outputs.data, 1)\n",
        "\n",
        "            x_test.extend(inputs.data.cpu().numpy())\n",
        "            y_test_pred.extend(predicted.data.cpu().numpy())\n",
        "\n",
        "    return np.array(x_test), np.array(y_test_pred)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "OT1B5ICpKMJX"
      },
      "outputs": [],
      "source": [
        "# Get preditions\n",
        "x_test, y_test_pred = predict_output(testloader)\n",
        "\n",
        "# Get labels\n",
        "y_test = []\n",
        "for i, data in enumerate(testloader):\n",
        "    y_test.extend(data[1].data.numpy())\n",
        "y_test = np.array(y_test)\n",
        "\n",
        "# Calculate accuracy\n",
        "correct = (y_test_pred == y_test)\n",
        "print(correct)\n",
        "print(f'Accuracy of the network on the 10000 test images: {100 * np.sum(correct)/len(correct)} %')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "STg_--nqd9QR"
      },
      "source": [
        "We can also visualize the predictions that are made. Below is a useful funtion 'visualize_predictions'. You do not have to understand what it does exactly. Run the next two cells to get a nice visualization of the predictions."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ZgXVjKW2qb5U"
      },
      "outputs": [],
      "source": [
        "# Define function to visualize predictions\n",
        "\n",
        "def visualize_predictions(x, y_true, y_pred, correct, label_to_name={}):\n",
        "    per_row = int(0.5 + sqrt(len(x)))\n",
        "    nr_rows = int(ceil(len(x) / float(per_row)))\n",
        "    fig, grid = plt.subplots(nr_rows, per_row, figsize=(10,10))  \n",
        "    cmap = None\n",
        "    if x.shape[-1] == 1:\n",
        "        cmap=\"gray\"    \n",
        "    for r in range(nr_rows):\n",
        "        for c in range(per_row):\n",
        "            if (r*per_row + c < len(x)):\n",
        "                grid[r,c].imshow(np.squeeze(np.moveaxis(x[r*per_row + c],0,2)), cmap=cmap)\n",
        "                grid[r,c].set_xticklabels([])\n",
        "                grid[r,c].set_xticklabels([])\n",
        "                grid[r,c].set_yticklabels([])\n",
        "                if label_to_name:\n",
        "                    grid[r,c].set_xlabel(\"Label = \" + str(label_to_name[y_true[r*per_row + c]]) + \"\\nPredicted = \" + str(label_to_name[y_pred[r*per_row + c]]))\n",
        "                else:\n",
        "                    grid[r,c].set_xlabel(\"Label = \" + str(y_true[r*per_row + c]) + \"\\nPredicted = \" + str(y_pred[r*per_row + c]))\n",
        "                for tic in grid[r][c].xaxis.get_major_ticks():\n",
        "                    tic.tick1line.set_visible(False)\n",
        "                    tic.tick2line.set_visible(False)\n",
        "                for tic in grid[r][c].yaxis.get_major_ticks():\n",
        "                    tic.tick1line.set_visible(False)\n",
        "                    tic.tick2line.set_visible(False)\n",
        "                color = \"green\"\n",
        "                if correct[r*per_row + c] == False:\n",
        "                    color = \"red\"\n",
        "                grid[r,c].spines['bottom'].set_color(color)\n",
        "                grid[r,c].spines['bottom'].set_linewidth(2)\n",
        "                grid[r,c].spines['top'].set_color(color)\n",
        "                grid[r,c].spines['top'].set_linewidth(2)\n",
        "                grid[r,c].spines['left'].set_color(color)\n",
        "                grid[r,c].spines['left'].set_linewidth(2)\n",
        "                grid[r,c].spines['right'].set_color(color)\n",
        "                grid[r,c].spines['right'].set_linewidth(2)\n",
        "            else :\n",
        "                grid[r,c].axis('off')\n",
        "    plt.tight_layout()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "nvL9QTQ5qdb7"
      },
      "outputs": [],
      "source": [
        "# Select some random images\n",
        "random_indices = np.random.randint(0,len(x_test), 16)\n",
        "\n",
        "# Visualize predictions\n",
        "visualize_predictions(x_test[random_indices], y_test[random_indices], y_test_pred[random_indices], correct[random_indices])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "JcmgpWKhOkeB"
      },
      "source": [
        "*What do you think, how does the model do?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dwKwdPCre3gW"
      },
      "source": [
        "To train our model we use the 'train_model' function below. It contains multiple loops. On the highest level you loop over epochs. During every epoch the model first loops over the training batches and subsequently over the validation batches. The loss and accuracy are printed for both the train and test set, after every epoch.\n",
        "\n",
        "To train the network, the loss is calculated for every train batch based on the loss function, i.e. 'criterion(outputs, labels)'. Subsequently the gradients are calculated using backpropagation with 'loss.backward'. Finally the weights are updated based on the gradients during 'optimizer.step()'. These steps are the core of training every deep learning model."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "EXs9dPG9F824"
      },
      "outputs": [],
      "source": [
        "def train_model(model,criterion,optimizer,trainloader,valloader,epochs):\n",
        "    train_loss = []\n",
        "    val_loss = []\n",
        "    train_acc = []\n",
        "    val_acc = []\n",
        "\n",
        "    for epoch in range(epochs):  # loop over the dataset multiple times\n",
        "\n",
        "        # train step\n",
        "        model.train(True)\n",
        "        running_loss = 0.0\n",
        "        running_acc = 0.0\n",
        "        total = 0.0\n",
        "\n",
        "        for i, data in enumerate(trainloader, 0):\n",
        "            # get the inputs; data is a list of [inputs, labels]\n",
        "            inputs, labels = data[0].to(device), data[1].to(device)\n",
        "\n",
        "            # zero the parameter gradients\n",
        "            optimizer.zero_grad()\n",
        "\n",
        "            # forward + backward + optimize\n",
        "            outputs = model(inputs)\n",
        "            loss = criterion(outputs, labels)\n",
        "            loss.backward()\n",
        "            optimizer.step()\n",
        "\n",
        "\n",
        "            # print statistics\n",
        "            running_loss += loss.item()\n",
        "            _, predicted = torch.max(outputs.data, 1)\n",
        "            running_acc += (predicted == labels).sum().item()\n",
        "            total += outputs.shape[0]\n",
        "\n",
        "        print(f'epoch: {epoch + 1}, loss: {running_loss / i:.3f}')\n",
        "        train_loss.append(running_loss / i)\n",
        "        train_acc.append(running_acc / total)\n",
        "\n",
        "        # validation step\n",
        "        model.train(False)\n",
        "        running_vloss = 0.0\n",
        "        running_vacc = 0.0\n",
        "        vtotal = 0.0\n",
        "\n",
        "        for i, data in enumerate(valloader, 0):\n",
        "            # get the inputs; data is a list of [inputs, labels]\n",
        "            inputs, labels = data[0].to(device), data[1].to(device)\n",
        "\n",
        "            # forward + backward + optimize\n",
        "            outputs = model(inputs)\n",
        "            loss = criterion(outputs, labels)\n",
        "\n",
        "            running_vloss += loss.item()\n",
        "            _, predicted = torch.max(outputs.data, 1)\n",
        "            running_vacc += (predicted == labels).sum().item()\n",
        "            vtotal += outputs.shape[0]\n",
        "\n",
        "        print(f'epoch: {epoch + 1}, vloss: {running_vloss / i:.3f}')\n",
        "        val_loss.append(running_vloss / i)\n",
        "        val_acc.append(running_vacc / vtotal)\n",
        "\n",
        "    print('Finished Training')\n",
        "\n",
        "    return train_loss, train_acc, val_loss, val_acc"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Hqcu2-HcqrLR"
      },
      "source": [
        "Now we are ready to train our model! We use 5 epochs here to get some quick results. Training for more epochs probably gives you better results, but at some point you will start to overfit."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "M1uFFQcpYRfc"
      },
      "outputs": [],
      "source": [
        "epochs = 5\n",
        "\n",
        "train_loss, train_acc, val_loss, val_acc = train_model(model,criterion,optimizer,trainloader,valloader,epochs)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "U17dGVh_isA1"
      },
      "source": [
        "We can visualize what happened during the training. Here we plot the accuracy and loss after each epoch."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "O_rn8Kj1W9yJ"
      },
      "outputs": [],
      "source": [
        "# Plot training & validation accuracy values\n",
        "plt.plot(train_acc)\n",
        "plt.plot(val_acc)\n",
        "plt.title('Model accuracy')\n",
        "plt.ylabel('Accuracy')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()\n",
        "\n",
        "# Plot training & validation loss values\n",
        "plt.plot(train_loss)\n",
        "plt.plot(val_loss)\n",
        "plt.title('Model loss')\n",
        "plt.ylabel('Loss')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "i7dUTJDa6v8N"
      },
      "source": [
        "*Exercise:*\n",
        "- *What happens to the model loss and model accuracy during training. Why?*\n",
        "- *How can differences between the train and validation performance be explained?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "es90Klo4i2Pc"
      },
      "source": [
        "Let's check out that accuracy and some of the results again after having trained the model!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ETdmoaEoUDDq"
      },
      "outputs": [],
      "source": [
        "# Let the model predict the test set labels\n",
        "x_test,y_test_pred = predict_output(testloader)\n",
        "\n",
        "y_test = []\n",
        "for i, data in enumerate(testloader):\n",
        "    y_test.extend(data[1].data.numpy())\n",
        "y_test = np.array(y_test)\n",
        "\n",
        "correct = (y_test_pred == y_test)\n",
        "\n",
        "print(f'Accuracy of the network on the 10000 test images: {100 * np.sum(correct)/len(correct)} %')\n",
        "\n",
        "random_indices = np.random.randint(0,x_test.shape[0], 16)\n",
        "visualize_predictions(x_test[random_indices], y_test[random_indices], y_test_pred[random_indices], correct[random_indices])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "3XFK2KHvp7bj"
      },
      "source": [
        "\n",
        "A nice way to visualize your results on a classification task is to create a confusion matrix. This shows for every class how often it is correctly classified or confused with another class."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "uArEqKo_0wdP"
      },
      "outputs": [],
      "source": [
        "def plt_confusion_matrix(y_test, y_test_pred, labels, normalize='true', log_cmap=False):\n",
        "  cm = confusion_matrix(y_test, y_test_pred, normalize = None)\n",
        "  if log_cmap:\n",
        "    cm2 = cm.copy().astype(float)\n",
        "    cm2[cm==cm.min()] = 0.1\n",
        "    norm  = LogNorm(vmin=cm2.min(), vmax=cm2.max())\n",
        "  else:\n",
        "    cm2 = cm\n",
        "    norm = None\n",
        "  df_cm = pd.DataFrame(cm2, index=labels, columns=labels)\n",
        "  plt.figure(figsize = (len(labels),len(labels)))\n",
        "  sn.heatmap(df_cm, annot=cm, norm=norm, cbar=False, fmt='.4g')\n",
        "  plt.xlabel(\"Predicted label\"); \n",
        "  plt.ylabel(\"True label\")\n",
        "  plt.title(\"Confusion matrix\")\n",
        "\n",
        "  \n",
        "plt_confusion_matrix(y_test, y_test_pred, labels=range(10), log_cmap=True)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-ipZtQjdQm-b"
      },
      "source": [
        "*Exercise:*\n",
        "\n",
        "*It is always good to have a look at the outcomes and errors as a sanity check if nothing went wrong.*\n",
        "- *Which digits got confused by the neural network most often? Does that make sense to you?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MzOrr6pZ3mrR"
      },
      "source": [
        "*Bonus exercise:*\n",
        "\n",
        "*The architecture of the CNN we trained is not the only one that could work. Variants on it could also work well for the task of classifying digits.* \n",
        "- *What variables could you change in the training process and what effects do you expect these changes to have?*\n",
        " - *E.g.: What are the trade-offs to keep in mind when increasing/decreasing the size or complexity of the model?*\n",
        "- *Try changing the model by altering the number of parameters (filters) of the convolutional (Conv2D) layers or adding or removing layers in the code that builds the model and see what happens!* \n",
        "- *Play around with other variables such as number of epochs, learning rate (lr) and others!* \n",
        "\n",
        "*(Note: If you change the code in a cell, the code cells below it have to be re-executed to work with the changes. So for example, don't forget to re-compile and re-train your model after rebuilding it!*)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "zPYedqB0jAKz"
      },
      "source": [
        "# PatchCamelyon: Classifying pathology patches as benign or malignant"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MFCWFqeZ7nmD"
      },
      "source": [
        "Now that we know how to train and test a convolutional neural network, let's use the same methods to determine malignancy of pathology patches! Will the neural network beat you at it?"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "W-6rzNs5QNCU"
      },
      "source": [
        "## **Understanding the data**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4OQ9VHJO8TEv"
      },
      "source": [
        "First, we download the pathology data we will use. Again we have a train, validation and test set. Because of time and memory constraints we use only 5% of the original dataset. If you are into this you can download the complete open source dataset here: https://patchcamelyon.grand-challenge.org/. Also we train for a limited amount of epochs to not keep you waiting for too long.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "DGaaDr-Yy9xW"
      },
      "outputs": [],
      "source": [
        "# load dataset\n",
        "!pip install --upgrade --no-cache-dir gdown\n",
        "import h5py    \n",
        "\n",
        "if not os.path.exists('./data'):\n",
        "  os.mkdir('./data')\n",
        "\n",
        "def download_pathology_data(identifier):\n",
        "    #gdd.download_file_from_google_drive(file_id=identifier, dest_path='./data/temp.h5.gz', overwrite=True, showsize=True)\n",
        "    !gdown --id $identifier -O './data/temp.h5.gz'\n",
        "    !gunzip -f ./data/temp.h5\n",
        "    x = h5py.File('./data/temp.h5','r+')['x'][()]\n",
        "    y = h5py.File('./data/temp.h5','r+')['y'][()].flatten()\n",
        "    os.remove('./data/temp.h5')\n",
        "\n",
        "    return x,y\n",
        "\n",
        "images_train, labels_train = download_pathology_data('1tLw33To0BplzqTx8OVG8okjmXaJjVE_K')\n",
        "images_val, labels_val = download_pathology_data('1iEouljzAbk8AtwoufRt8WK0aqxPO5MYd')\n",
        "images_test,  labels_test  = download_pathology_data('14eQ32RHVbj112zf7ud01hndyqVrAITw1')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ttdLrV0fDFv0"
      },
      "source": [
        "We use the same patches for the test set as were used in the reader study earlier today. This way, in the end you can compare your performance with that of the neural network."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "6pEsJFk2zHmh"
      },
      "outputs": [],
      "source": [
        "# print data shape\n",
        "print('Training set input size:', images_train.shape)\n",
        "print('Training set output size:', labels_train.shape)\n",
        "print('Validation set input size:', images_val.shape)\n",
        "print('Validation set output size:', labels_val.shape)\n",
        "print('Testing set input size:', images_test.shape)\n",
        "print('Testing set output size:', labels_test.shape) "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "o0dbV-2TFS_5"
      },
      "source": [
        "*Exercise:* \n",
        "- *What differences do you spot between the dimensions of the MNIST data and the dimensions of the pathology data?*\n",
        "- *What changes would we need to make to the neural network we used earlier to account for these differences?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GzppkY2SGSoa"
      },
      "source": [
        "*Exercise:*\n",
        "- *Again, we should take a look at the input images to make sure we have done everything right so far. Change the index `i` to view different patches.*"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "TvW-3pnG9HbG"
      },
      "outputs": [],
      "source": [
        "i=0\n",
        "plt.figure()\n",
        "plt.imshow(images_train[i])\n",
        "plt.title('label: {}'.format(labels_train.flatten()[i]));"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7kcfc4C-HdVU"
      },
      "source": [
        "*Exercise:* \n",
        "- *What does the image label (title above the plot) mean?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eGNx5IgpI88o"
      },
      "source": [
        "Below we create the dataset and dataloader object. A TensorDataset is a useful way to organize your data. It allows you to automatically keep track of images and their labels and you can perform transformations on the data if you want, such as normalization or augmentations. \n",
        "\n",
        "You can also create your own custom Dataset object that allows you to further specify how your data is modified before it is used to train the network."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "u1jXhRbqv3Cg"
      },
      "outputs": [],
      "source": [
        "from torch.utils.data import TensorDataset\n",
        "from torch.utils.data import DataLoader\n",
        "\n",
        "batch_size = 128  \n",
        "\n",
        "patchcam_train = TensorDataset(torch.tensor(np.moveaxis(images_train,3,1), dtype = torch.float32),torch.tensor(labels_train, dtype = torch.uint8))\n",
        "patchcam_val = TensorDataset(torch.tensor(np.moveaxis(images_val,3,1), dtype = torch.float32),torch.tensor(labels_val, dtype = torch.uint8))\n",
        "patchcam_test = TensorDataset(torch.tensor(np.moveaxis(images_test,3,1), dtype = torch.float32),torch.tensor(labels_test, dtype = torch.uint8))\n",
        "\n",
        "trainloader = DataLoader(patchcam_val, batch_size = batch_size, shuffle = True)\n",
        "valloader = DataLoader(patchcam_val, batch_size = batch_size, shuffle = False)\n",
        "testloader = DataLoader(patchcam_test, batch_size = batch_size, shuffle = False)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "1_-nfMhtQXQX"
      },
      "source": [
        "## **Training and testing a convolutional neural network**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TgjEQP1vQ8ae"
      },
      "source": [
        "### **Building the model** "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "O21hi326J4GC"
      },
      "source": [
        "Execute the following three cells to build and compile the same model architecture that we trained for the MNIST."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "c-NIPutIoEGa"
      },
      "outputs": [],
      "source": [
        "# define model\n",
        "\n",
        "class Net(nn.Module):\n",
        "    def __init__(self):\n",
        "        super().__init__()\n",
        "        self.conv1 = nn.Conv2d(3, 6, 5)\n",
        "        self.pool = nn.MaxPool2d(2, 2)\n",
        "        self.conv2 = nn.Conv2d(6, 16, 5)\n",
        "        self.fc1 = nn.Linear(7056, 120)\n",
        "        self.fc2 = nn.Linear(120, 84)\n",
        "        self.fc3 = nn.Linear(84, 2)\n",
        "\n",
        "    def forward(self, x):\n",
        "        x = self.pool(F.relu(self.conv1(x)))\n",
        "        x = self.pool(F.relu(self.conv2(x)))\n",
        "        x = torch.flatten(x, 1) # flatten all dimensions except batch\n",
        "        x = F.relu(self.fc1(x))\n",
        "        x = F.relu(self.fc2(x))\n",
        "        x = self.fc3(x)\n",
        "        return x\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "CejHsJKkoYQ4"
      },
      "outputs": [],
      "source": [
        "device = torch.device('cuda:0' if torch.cuda.is_available() else 'cpu')\n",
        "\n",
        "# assuming that we are on a CUDA machine, this should print a CUDA device:\n",
        "\n",
        "print(device)\n",
        "\n",
        "model = Net()\n",
        "model.to(device)\n",
        "criterion = nn.CrossEntropyLoss()\n",
        "optimizer = optim.SGD(model.parameters(), lr=0.01, momentum=0.9)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GaVPNaO4Qz2s"
      },
      "source": [
        "### **Training and testing**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "wSa9lvQrSQAf"
      },
      "source": [
        "Execute the following two cells to see that the model performance is worthless again before training."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "kSTlvg1boYQ5"
      },
      "outputs": [],
      "source": [
        "# Let the model predict the test set labels\n",
        "x_test,y_test_pred = predict_output(testloader)\n",
        "\n",
        "y_test = labels_test\n",
        "correct = (y_test_pred == y_test)\n",
        "\n",
        "print(f'Accuracy of the network on the 10000 test images: {100 * np.sum(correct)/len(correct)} %')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "L4h5Ih07aGSf"
      },
      "outputs": [],
      "source": [
        "cm = confusion_matrix(y_test, y_test_pred, normalize = None)\n",
        "disp = ConfusionMatrixDisplay(confusion_matrix=cm,\n",
        "                              display_labels=['benign' , 'malignant'])\n",
        "disp.plot(values_format='d'); "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "QeF3ljIVScSe"
      },
      "source": [
        "Now let's train our malignancy classification model!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "y_T6n_XKrTOP"
      },
      "outputs": [],
      "source": [
        "epochs = 5\n",
        "\n",
        "train_loss, train_acc, val_loss, val_acc = train_model(model,criterion,optimizer,trainloader,valloader,epochs)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "lRaQ-dT9S88S"
      },
      "source": [
        "Just like with training the MINIST model, we visualize the training process here."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "sjRTvjhimsAP"
      },
      "outputs": [],
      "source": [
        "# Plot training & validation accuracy values\n",
        "plt.plot(train_acc)\n",
        "plt.plot(val_acc)\n",
        "plt.title('Model accuracy')\n",
        "plt.ylabel('Accuracy')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()\n",
        "\n",
        "# Plot training & validation loss values\n",
        "plt.plot(train_loss)\n",
        "plt.plot(val_loss)\n",
        "plt.title('Model loss')\n",
        "plt.ylabel('Loss')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "8jgqYmXKmxfW"
      },
      "source": [
        "*Exercise:*\n",
        "- *What is happening here? Is the model learning/getting better just as with the MNIST dataset?*\n",
        "- *What could be reasons that the model is performing this way?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "JHVtwI1pN_Co"
      },
      "source": [
        "When you understand what is happening, continue to the next section."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "BQd_nhOJRIaA"
      },
      "source": [
        "## **A different neural network**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OkBVzidVWBxR"
      },
      "source": [
        "So our first try did not work. Run the next two blocks to train a different kind of network architecture. Run them first, then read the code and text around it while you wait for it to learn ;)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "HshOZrHPE1P0"
      },
      "outputs": [],
      "source": [
        "model = torch.hub.load('pytorch/vision:v0.10.0', 'mobilenet_v2', pretrained=True)\n",
        "\n",
        "model.to(device)\n",
        "criterion = nn.CrossEntropyLoss()\n",
        "optimizer = optim.SGD(model.parameters(), lr=0.01, momentum=0.9)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "gtuO9DU1E-du"
      },
      "outputs": [],
      "source": [
        "epochs = 20\n",
        "\n",
        "train_loss, train_acc, val_loss, val_acc = train_model(model,criterion,optimizer,trainloader,valloader,epochs)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "WpFPR0EVOiux"
      },
      "source": [
        "The previous cell of code trains the MobileNetV2 model: a standardized network architecture. As you can see in the summary, this model has more and also different layers than our previous model.\n",
        "\n",
        "Instead of training it from scratch with random initial weights, we train it using transfer learning. As our starting point for training, we use the architecture and weights of the MobileNetV2 model that was trained on a completely different image dataset (ImageNet). This way, we only need to fine-tune the weights to make the model work for our malignancy classification problem. \n",
        "\n",
        "When the model is done training, you can run the next code cell to visualize the training process."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "2xTcw9JwrTOU"
      },
      "outputs": [],
      "source": [
        "# Plot training & validation accuracy values\n",
        "plt.plot(train_acc)\n",
        "plt.plot(val_acc)\n",
        "plt.title('Model accuracy')\n",
        "plt.ylabel('Accuracy')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()\n",
        "\n",
        "# Plot training & validation loss values\n",
        "plt.plot(train_loss)\n",
        "plt.plot(val_loss)\n",
        "plt.title('Model loss')\n",
        "plt.ylabel('Loss')\n",
        "plt.xlabel('Epoch')\n",
        "plt.legend(['Train', 'Val'], loc='upper left')\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "JdrV2TPkVGLA"
      },
      "source": [
        "That's looking better! "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "45HRELLHeQTn"
      },
      "source": [
        "*Exercise:*\n",
        "\n",
        "- *Is there a difference between the performance on the train and validation set? If so, how can you explain this difference?* "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "NqSnon5pR9Xi"
      },
      "source": [
        "## **The neural network's performance on the reader study**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Vxz5iZG2YmEP"
      },
      "source": [
        "So far, the model was only trained and tested with data from the training and validation data sets. \n",
        "\n",
        "Now let's see how well the model on the data you scored during the reader study earlier today!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "B-Jsx6x1qcaw"
      },
      "outputs": [],
      "source": [
        "# Let the model predict the test set labels\n",
        "x_test,y_test_pred = predict_output(testloader)\n",
        "\n",
        "y_test = labels_test\n",
        "correct = (y_test_pred == y_test)\n",
        "\n",
        "print(f'Accuracy of the network on the 10000 test images: {100 * np.sum(correct)/len(correct)} %')\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "caZbvGJ4xQlZ"
      },
      "outputs": [],
      "source": [
        "# Plot confusion matrix\n",
        "cm = confusion_matrix(y_test, y_test_pred, normalize = None)\n",
        "disp = ConfusionMatrixDisplay(confusion_matrix=cm,\n",
        "                              display_labels=['benign', 'malignant'])\n",
        "disp.plot(values_format='d');"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "e3XhdlnqX-JS"
      },
      "source": [
        "*Exercise:*\n",
        "- *So who would be better suited as a pathologist? You or the trained neural network? (Did you forget your own score, go to https://grand-challenge.org/reader-studies/)*\n",
        "- *What kind of mistakes were made? If you would not want to miss any malginancies with the network, which number should be the lowest of them all?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AxXwl8mhXxcp"
      },
      "source": [
        "The network outputs a probability value between 0 and 1. Rounding this number gives us a binary classification of benign (0) and malignant (1).\n",
        "Check out some examples of your model's predictions."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "hdQjzn7CoYQ7"
      },
      "outputs": [],
      "source": [
        "# Run again if you want to see a different set of example images\n",
        "random_indices = np.random.randint(0, len(x_test), 16)\n",
        "visualize_predictions(np.array(x_test)[random_indices]/255, np.array(y_test)[random_indices], np.array(y_test_pred)[random_indices], np.array(correct)[random_indices])"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "A-7A1Q4uaEtw"
      },
      "source": [
        "To show the distribution of the predicted probability scores, we plot histograms."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "BXWSoRpFlQom"
      },
      "outputs": [],
      "source": [
        "plt.figure()\n",
        "plt.hist(np.squeeze(y_test_pred), bins=10)\n",
        "plt.title('Histogram of all predictions')\n",
        "plt.ylabel('Number of images')\n",
        "plt.xlabel('Probability')\n",
        "plt.figure()\n",
        "plt.hist(np.squeeze(y_test_pred)[correct], bins=10, color='green')\n",
        "plt.title('Histogram of correct predictions')\n",
        "plt.ylabel('Number of images')\n",
        "plt.xlabel('Probability')\n",
        "plt.figure()\n",
        "plt.hist(np.squeeze(y_test_pred)[~correct], bins=10, color='red')\n",
        "plt.title('Histogram of incorrect predictions')\n",
        "plt.ylabel('Number of images')\n",
        "plt.xlabel('Probability');"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "MJBXMgZLaVeY"
      },
      "source": [
        "*Exercise:*\n",
        "- *In what probability range were most mistakes made?*\n",
        "- *Think of measures you could take as an AI-software developer or clinician to make less mistakes based on the probability values. How can clinicians and machines work together?*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "IIXgBb-z5OM6"
      },
      "source": [
        "*Bonus exercise:*\n",
        "- *Within the Patch Camelyon section, go back to the first model architecture block (the one that is the same as we used for the MNIST). Can you create an architecture yourself that performs better? Think of changing layers, using drop-out, etc. You could also play around with parameters as batch size and learning rate. NB: Make sure it does not only perform well on the trainset but also on the validation set!*"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "WIBR_hgizLTY"
      },
      "source": [
        "## Camelyon heatmap code"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hYU__VoywWwc"
      },
      "source": [
        "The network we trained predicts for a patch of 96 by 96 whether or not there is malignant tissue present. However, in real life, pathologists look at way larger images. We load a larger tile of an orginal image here and use our network to locate the malignant tissue."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "y0Ca0DuD_yMP"
      },
      "outputs": [],
      "source": [
        "import zipfile\n",
        "import numpy as np\n",
        "from matplotlib import pyplot as plt\n",
        "from PIL import Image\n",
        "from google_drive_downloader import GoogleDriveDownloader as gdd\n",
        "\n",
        "gdd.download_file_from_google_drive(file_id='1G7rBoPoma5cmsQSPGF05NoS4x4SQ4Qa1', dest_path='./data/camelyon-tiles.zip', overwrite=True)\n",
        "with zipfile.ZipFile('./data/camelyon-tiles.zip', \"r\") as zip_ref:\n",
        "  zip_ref.extractall(\"./data\")\n",
        "os.remove('./data/camelyon-tiles.zip')\n",
        "\n",
        "npz = np.load(\"./data/camelyon-test-tile-1.npz\")\n",
        "\n",
        "tile = npz['tile']\n",
        "annotations = np.asarray(Image.open(\"./data/camelyon-test-tile-1.png\"))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "R57rzbJdsh-1"
      },
      "outputs": [],
      "source": [
        "plt.figure(figsize=(30,30))\n",
        "plt.subplot(1, 3, 1); plt.imshow(tile); plt.axis('off');\n",
        "plt.title('1. Larger pathology image')\n",
        "plt.subplot(1, 3, 2); plt.imshow(annotations); plt.axis('off');\n",
        "plt.title('2. Manual annotations for malignancy')\n",
        "plt.subplot(1, 3, 3); plt.imshow(tile); plt.ylim([tile.shape[0], 0]), plt.xlim([0,tile.shape[0]])\n",
        "plt.vlines(range(96,tile.shape[0], 96), ymin=0, ymax=tile.shape[0], linestyle='-', linewidth=1)\n",
        "plt.hlines(range(96,tile.shape[0], 96), xmin=0, xmax=tile.shape[0], linestyle='-', linewidth=1)\n",
        "plt.title('3. Grid showing the patch size we trained the network on')\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "C2rnsBUbzFit"
      },
      "source": [
        "Above we see the original image (1) and the manual annotation of the malignancy (2). To get a feel for the scale we are looking at, the most right image (3) has a grid overlay where each square is 96 by 96 such as the patches we trained the network on."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "8O-BuCzA0a8U"
      },
      "source": [
        "To get predictions for the larger tile we slide over the image and make a prediction for each 96 by 96 patch. We do this with overlap to create a more detailed mapping (`stride=16`).\n",
        "\n",
        "![alt text](https://raw.githubusercontent.com/iamaaditya/iamaaditya.github.io/master/images/conv_arithmetic/full_padding_no_strides_transposed.gif)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "rn22GEPzC98l"
      },
      "outputs": [],
      "source": [
        "tile_size = tile.shape[0]\n",
        "patch_size = 96\n",
        "stride = 16\n",
        "\n",
        "# Slide over the large tile\n",
        "patches = []\n",
        "for hor in range(0, tile_size-patch_size, stride):\n",
        "  for vert in range(0, tile_size-patch_size, stride): \n",
        "    patch = tile[hor:hor+patch_size, vert:vert+patch_size]\n",
        "    patches.append(patch)\n",
        "\n",
        "patches_test = TensorDataset(torch.tensor(np.moveaxis(patches,3,1), dtype = torch.float32))\n",
        "\n",
        "patchloader = DataLoader(patches_test, batch_size = batch_size, shuffle = False)\n",
        "\n",
        "# Predict the outcome of all patches in the slide at once.\n",
        "_, patch_pred = predict_output(patchloader)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "NaMGZDzeJVB1"
      },
      "outputs": [],
      "source": [
        "# Some resizing to make the heatmap overlay the original tile\n",
        "import scipy.ndimage\n",
        "heatmap_size=np.int((tile_size-patch_size)/stride)  \n",
        "heat_map = np.array(patch_pred).reshape([heatmap_size,heatmap_size])\n",
        "heat_map = scipy.ndimage.zoom(heat_map, stride, order=1)\n",
        "heat_map_pad = np.pad(heat_map, np.int(patch_size/2), mode='constant')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "8KjR1YLeHNPO"
      },
      "outputs": [],
      "source": [
        "# Plotting the results\n",
        "plt.figure(figsize=(30,30))\n",
        "plt.subplot(1, 2, 1); plt.imshow(tile); plt.axis('off');\n",
        "plt.imshow(heat_map_pad, 'jet', alpha=0.5); \n",
        "plt.title('Heatmap of probability values for malignancy likelihood predicted by the network')\n",
        "plt.subplot(1, 2, 2); plt.imshow(annotations); plt.axis('off');\n",
        "plt.title('Manual annotations for malignancy')\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "r7RHCQyB11Vi"
      },
      "source": [
        "This shows how a neural network trained on small images (a simpler task) can aid in the classification or even segmentation of larger images."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "F0fS0r0f9b5o"
      },
      "source": [
        "# Wrap up\n",
        "\n",
        "Great! You did it!\n",
        "You can download the notebook to your own pc or save it to your drive to continue playing with it at home. Try your own data or different networks and see if you can get better at it!"
      ]
    }
  ],
  "metadata": {
    "accelerator": "GPU",
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3",
      "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.11.4"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}