Heatmap_CAM_Visualization.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This code helps to visualize the salient regions in the image used by the model to arrive at the predictions. We perform grad-CAM visualizations to identify the salient image regions. The process involves computing the output feature map of the deepst convolutional layer in the trained model and compute the Gradient of the expected class (abnormal/normal) with regard to the output feature map of this layer. We then multiply each channel in the feature-map array by “how important this channel is” with regard to the “abnormal/normal” class and obtain the heatmap. For visualization purposes, we normalize the heatmap between 0 and 1. An image is then generated by superimposing the original image on the heatmap obtained. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To begin with, let us define a few functions to load the data and convert them to Keras compatible targets. We will load the libraries to begin with."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# load libraries\n",
    "import cv2\n",
    "import numpy as np\n",
    "import os\n",
    "from keras.utils import np_utils\n",
    "import matplotlib.pyplot as plt\n",
    "import itertools\n",
    "import time\n",
    "from keras.models import Sequential\n",
    "from keras.models import Model\n",
    "from keras import applications\n",
    "from keras.callbacks import ModelCheckpoint, TensorBoard\n",
    "from keras.layers import Conv2D, Activation, Dense, MaxPooling2D, Flatten, Dropout, GlobalAveragePooling2D\n",
    "from sklearn.metrics import log_loss\n",
    "from sklearn.utils import class_weight\n",
    "from keras.optimizers import SGD\n",
    "from sklearn.metrics import roc_curve, auc\n",
    "from sklearn.metrics import classification_report,confusion_matrix, accuracy_score\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.metrics import average_precision_score\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We performed 5-fold cross validation at the patient level. we had train and test splits for each fold to ensure that none of the patienet information in the training data leaks into the test data. We randomly split 10% of the training data for validation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#define data directories\n",
    "train_data_dir = 'f1_mal/train' #path to your data\n",
    "valid_data_dir = 'f1_mal/valid'\n",
    "test_data_dir = 'f1_mal/test'\n",
    "\n",
    "# declare the number of samples in each category\n",
    "nb_train_samples = 19808 #  modify for your dataset\n",
    "nb_valid_samples = 4952 #  modify for your dataset\n",
    "nb_test_samples = 2730 # modify for your dataset\n",
    "num_classes = 2 # binary classification \n",
    "img_rows_orig = 100 # modify these values depending on your requirements\n",
    "img_cols_orig = 100"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lets define functions to load and resize the training, validation and test data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def load_training_data():\n",
    "    labels = os.listdir(train_data_dir)\n",
    "    total = len(labels)\n",
    "    X_train = np.ndarray((nb_train_samples, img_rows_orig, img_cols_orig, 3), dtype=np.uint8)\n",
    "    Y_train = np.zeros((nb_train_samples,), dtype='uint8')\n",
    "    i = 0\n",
    "    print('-'*30)\n",
    "    print('Creating training images...')\n",
    "    print('-'*30)\n",
    "    j = 0\n",
    "    for label in labels:\n",
    "        image_names_train = os.listdir(os.path.join(train_data_dir, label))\n",
    "        total = len(image_names_train)\n",
    "        print(label, total)\n",
    "        for image_name in image_names_train:\n",
    "            img = cv2.imread(os.path.join(train_data_dir, label, image_name), cv2.IMREAD_COLOR)\n",
    "            img = np.array([img])\n",
    "            X_train[i] = img\n",
    "            Y_train[i] = j\n",
    "            if i % 100 == 0:\n",
    "                print('Done: {0}/{1} images'.format(i, total))\n",
    "            i += 1\n",
    "        j += 1    \n",
    "    print(i)                \n",
    "    print('Loading done.')\n",
    "    print('Transform targets to keras compatible format.')\n",
    "    Y_train = np_utils.to_categorical(Y_train[:nb_train_samples], num_classes)\n",
    "    np.save('imgs_train.npy', X_train, Y_train) #save as numpy files\n",
    "    return X_train, Y_train\n",
    "    \n",
    "def load_validation_data():\n",
    "    # Load validation images\n",
    "    labels = os.listdir(valid_data_dir)\n",
    "    X_valid = np.ndarray((nb_valid_samples, img_rows_orig, img_cols_orig, 3), dtype=np.uint8)\n",
    "    Y_valid = np.zeros((nb_valid_samples,), dtype='uint8')\n",
    "    i = 0\n",
    "    print('-'*30)\n",
    "    print('Creating validation images...')\n",
    "    print('-'*30)\n",
    "    j = 0\n",
    "    for label in labels:\n",
    "        image_names_valid = os.listdir(os.path.join(valid_data_dir, label))\n",
    "        total = len(image_names_valid)\n",
    "        print(label, total)\n",
    "        for image_name in image_names_valid:\n",
    "            img = cv2.imread(os.path.join(valid_data_dir, label, image_name), cv2.IMREAD_COLOR)\n",
    "            img = np.array([img])\n",
    "            X_valid[i] = img\n",
    "            Y_valid[i] = j\n",
    "            if i % 100 == 0:\n",
    "                print('Done: {0}/{1} images'.format(i, total))\n",
    "            i += 1\n",
    "        j += 1\n",
    "    print(i)            \n",
    "    print('Loading done.')\n",
    "    print('Transform targets to keras compatible format.');\n",
    "    Y_valid = np_utils.to_categorical(Y_valid[:nb_valid_samples], num_classes)\n",
    "    np.save('imgs_valid.npy', X_valid, Y_valid) #save as numpy files\n",
    "    return X_valid, Y_valid\n",
    "\n",
    "def load_test_data():\n",
    "    labels = os.listdir(test_data_dir)\n",
    "    X_test = np.ndarray((nb_test_samples, img_rows_orig, img_cols_orig, 3), dtype=np.uint8)\n",
    "    Y_test = np.zeros((nb_test_samples,), dtype='uint8')\n",
    "    i = 0\n",
    "    print('-'*30)\n",
    "    print('Creating test images...')\n",
    "    print('-'*30)\n",
    "    j = 0\n",
    "    for label in labels:\n",
    "        image_names_test = os.listdir(os.path.join(test_data_dir, label))\n",
    "        total = len(image_names_test)\n",
    "        print(label, total)\n",
    "        for image_name in image_names_test:\n",
    "            img = cv2.imread(os.path.join(test_data_dir, label, image_name), cv2.IMREAD_COLOR)\n",
    "            img = np.array([img])\n",
    "            X_test[i] = img\n",
    "            Y_test[i] = j\n",
    "            if i % 100 == 0:\n",
    "                print('Done: {0}/{1} images'.format(i, total))\n",
    "            i += 1\n",
    "        j += 1\n",
    "    print(i)            \n",
    "    print('Loading done.')\n",
    "    print('Transform targets to keras compatible format.');\n",
    "    Y_test = np_utils.to_categorical(Y_test[:nb_test_samples], num_classes)\n",
    "    np.save('imgs_test.npy', X_test, Y_test) #save as numpy files\n",
    "    return X_test, Y_test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will define functions to resize the original images to that dimensions required for the pretrained models using the functions defined below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def load_resized_training_data(img_rows, img_cols):\n",
    "\n",
    "    X_train, Y_train = load_training_data()\n",
    "    X_train = np.array([cv2.resize(img, (img_rows,img_cols)) for img in X_train[:nb_train_samples,:,:,:]])\n",
    "    \n",
    "    return X_train, Y_train\n",
    "    \n",
    "def load_resized_validation_data(img_rows, img_cols):\n",
    "\n",
    "    X_valid, Y_valid = load_validation_data()\n",
    "    X_valid = np.array([cv2.resize(img, (img_rows,img_cols)) for img in X_valid[:nb_valid_samples,:,:,:]])\n",
    "        \n",
    "    return X_valid, Y_valid   \n",
    "\n",
    "def load_resized_test_data(img_rows, img_cols):\n",
    "\n",
    "    X_test, Y_test = load_test_data()\n",
    "    X_test = np.array([cv2.resize(img, (img_rows,img_cols)) for img in X_test[:nb_test_samples,:,:,:]])\n",
    "    \n",
    "    return X_test, Y_test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An evaluation script has been written to compute the confusion matrix for the performance of the trained model. This function prints and plots the confusion matrix. Normalization can be applied by setting 'normalize=True'."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def plot_confusion_matrix(cm, classes,\n",
    "                          normalize=False, #if true all values in confusion matrix is between 0 and 1\n",
    "                          title='Confusion matrix',\n",
    "                          cmap=plt.cm.Blues):\n",
    "    plt.imshow(cm, interpolation='nearest', cmap=cmap)\n",
    "    plt.title(title)\n",
    "    plt.colorbar()\n",
    "    tick_marks = np.arange(len(classes))\n",
    "    plt.xticks(tick_marks, classes, rotation=45)\n",
    "    plt.yticks(tick_marks, classes)\n",
    "\n",
    "    if normalize:\n",
    "        cm = cm.astype('float') / cm.sum(axis=1)[:, np.newaxis]\n",
    "        print(\"Normalized confusion matrix\")\n",
    "    else:\n",
    "        print('Confusion matrix, without normalization')\n",
    "\n",
    "    print(cm)\n",
    "    thresh = cm.max() / 2.\n",
    "    for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n",
    "        plt.text(j, i, cm[i, j],\n",
    "                 horizontalalignment=\"center\",\n",
    "                 color=\"white\" if cm[i, j] > thresh else \"black\")\n",
    "    plt.tight_layout()\n",
    "    plt.ylabel('True label')\n",
    "    plt.xlabel('Predicted label')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will now proceed to extract the features from our dataset using the pretrained models."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "img_rows=224 \n",
    "img_cols=224\n",
    "channel = 3 #RGB\n",
    "num_classes = 2 #binary classification\n",
    "batch_size = 32 # modify based on the GPUs in your system\n",
    "num_epoch = 100 # modify depending on the model's convergence with your data\n",
    "\n",
    "#load data\n",
    "X_train, Y_train = load_resized_training_data(img_rows, img_cols)\n",
    "X_valid, Y_valid = load_resized_validation_data(img_rows, img_cols)\n",
    "X_test, Y_test = load_resized_test_data(img_rows, img_cols)\n",
    "\n",
    "\n",
    "#print the shape of the data\n",
    "print(X_train.shape, Y_train.shape)\n",
    "print(X_valid.shape, Y_valid.shape)\n",
    "print(X_test.shape, Y_test.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us now configure our pretrained model. This code uses VGG16 as a feature extractor.\n",
    "\n",
    "you can use the rest of the models like:\n",
    "\n",
    "ResNet50:\n",
    "\n",
    "feature_model = applications.ResNet50((weights='imagenet', include_top=False, input_shape=(img_rows, img_cols, 3)) \n",
    "feature_model = Model(input=feature_model.input, output=feature_model.get_layer('res5c_branch2c').output) \n",
    "\n",
    "Xception:\n",
    "\n",
    "feature_model = applications.Xception((weights='imagenet', include_top=False, input_shape=(img_rows, img_cols, 3))\n",
    "feature_model = Model(input=feature_model.input, output=feature_model.get_layer('block14_sepconv1').output) \n",
    "\n",
    "DenseNet121:\n",
    "\n",
    "For DenseNet, the main file densenet121_model is included to this repository. The model can be used as :\n",
    "feature_model = densenet121_model(img_rows=img_rows, img_cols=img_cols, color_type=channel, num_classes=num_classes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "base_model = applications.VGG16(weights='imagenet', include_top=False, input_shape=(img_rows, img_cols, channel))\n",
    "\n",
    "#extract feature from the optimal layer for your data\n",
    "base_model = Model(input=base_model.input, output=base_model.get_layer('block5_conv2').output) \n",
    "\n",
    "#get the model summary\n",
    "base_model.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lets modify the architecture by adding a global spatial average pooling layer and a fully-connected layer with a dropout ratio of 0.5 to prevent overfitting and help model generalization. We will train only the top layers which are randomly initialized, freeze all the convolutional layers to prevent large gradient updates wrecking the learned weights. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x = base_model.output\n",
    "x = GlobalAveragePooling2D()(x)\n",
    "x = Dense(1024, activation='relu')(x)\n",
    "x = Dropout(0.5)(x)\n",
    "predictions = Dense(num_classes, activation='softmax', name='predictions')(x)\n",
    "\n",
    "# this is the model we will train\n",
    "model = Model(inputs=base_model.input, outputs=predictions)\n",
    "\n",
    "# Train only the top layers (which were randomly initialized)\n",
    "# i.e. freeze all convolutional layers to prevent large gradient updates wrecking the learned weights\n",
    "for layer in base_model.layers:\n",
    "    layer.trainable = False\n",
    "\n",
    "#fix the optimizer\n",
    "sgd = SGD(lr=0.00001, decay=1e-6, momentum=0.9, nesterov=True) \n",
    "\n",
    "#compile the model\n",
    "model.compile(optimizer=sgd,\n",
    "              loss='categorical_crossentropy',\n",
    "              metrics=['accuracy'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Its time to train the model. We will store only the best model weights by initializing callbacks. Also we can view the performance of our model during run-time by visualizing the performance graphs with Tensorboard. Create a log directory named 'logs' to store the training logs and a separate folder named 'weights' to store the model weights. You can visualize tensorboard graphs simply by navigating to your working directory and do:\n",
    "\n",
    "$tensorboard --logdir=path/to/log-directory/ --port 6006\n",
    "\n",
    "Then open localhost:6006 in your browser to view the performance graphs, model architecture and other parameters of your interest."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "filepath = 'weights/' + model.name + '.{epoch:02d}-{val_acc:.4f}.hdf5'\n",
    "checkpoint = ModelCheckpoint(filepath, monitor='val_acc', verbose=1, \n",
    "                             save_weights_only=True, save_best_only=True, mode='max', period=1)\n",
    "tensor_board = TensorBoard(log_dir='logs/', histogram_freq=0, batch_size=batch_size)\n",
    "callbacks_list = [checkpoint, tensor_board]\n",
    "\n",
    "#compute training time\n",
    "t=time.time()\n",
    "hist = model.fit(X_train, Y_train, batch_size=batch_size, \n",
    "                 callbacks=callbacks_list,\n",
    "                 epochs=num_epoch, verbose=1, \n",
    "                 shuffle=True, validation_data=[X_valid, Y_valid])\n",
    "\n",
    "#compute the training time\n",
    "print('Training time: %s' % (time.time()-t))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to visualize the performance of the model in the console other than with Tensorboard, you can use the following script."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "train_loss=hist.history['loss']\n",
    "val_loss=hist.history['val_loss']\n",
    "train_acc=hist.history['acc']\n",
    "val_acc=hist.history['val_acc']\n",
    "xc=range(num_epoch)\n",
    "\n",
    "plt.figure(1,figsize=(20,10), dpi=100)\n",
    "plt.plot(xc,train_loss)\n",
    "plt.plot(xc,val_loss)\n",
    "plt.xlabel('num of Epochs')\n",
    "plt.ylabel('loss')\n",
    "plt.title('train_loss vs val_loss')\n",
    "plt.grid(True)\n",
    "plt.legend(['train','val'])\n",
    "plt.style.use(['classic'])\n",
    "\n",
    "plt.figure(2,figsize=(20,10), dpi=100)\n",
    "plt.plot(xc,train_acc)\n",
    "plt.plot(xc,val_acc)\n",
    "plt.xlabel('num of Epochs')\n",
    "plt.ylabel('accuracy')\n",
    "plt.title('train_acc vs val_acc')\n",
    "plt.grid(True)\n",
    "plt.legend(['train','val'])\n",
    "plt.style.use(['classic'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the model is trained, load the best model weights to predict on the test data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.load_weights('weights/model_2.01-0.8546.hdf5') #modify for your own model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#predict on the test data\n",
    "X_test, Y_test = load_resized_test_data(img_rows, img_cols)\n",
    "print(X_test.shape, Y_test.shape)\n",
    "print('-'*30)\n",
    "print('Predicting on the test data...')\n",
    "print('-'*30)\n",
    "y_pred = model.predict(X_test, batch_size=batch_size, verbose=1)\n",
    "\n",
    "# compute the accuracy\n",
    "Test_accuracy = accuracy_score(Y_test.argmax(axis=-1),y_pred.argmax(axis=-1))\n",
    "print(\"Test_Accuracy = \",Test_accuracy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us now compute the performance metrics for the pretrained VGG16 model with the test data. The performance metrics involve computing the ROC-AUC values, cross-entropy loss score, average precision score, prediction probabilities and storing these values and plotting the ROC curves."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#compute the ROC-AUC values\n",
    "fpr = dict()\n",
    "tpr = dict()\n",
    "roc_auc = dict()\n",
    "for i in range(num_classes):\n",
    "    fpr[i], tpr[i], _ = roc_curve(Y_test[:, i], y_pred[:, i])\n",
    "    roc_auc[i] = auc(fpr[i], tpr[i])\n",
    "    \n",
    "# Compute micro-average ROC curve and ROC area\n",
    "fpr[\"micro\"], tpr[\"micro\"], _ = roc_curve(Y_test.ravel(), y_pred.ravel())\n",
    "roc_auc[\"micro\"] = auc(fpr[\"micro\"], tpr[\"micro\"])\n",
    "\n",
    "#Plot ROC curves\n",
    "plt.figure(figsize=(20,10), dpi=100)\n",
    "lw = 1\n",
    "plt.plot(fpr[1], tpr[1], color='red',\n",
    "         lw=lw, label='ROC curve (area = %0.4f)' % roc_auc[1])\n",
    "plt.plot([0, 1], [0, 1], color='black', lw=lw, linestyle='--')\n",
    "plt.xlim([0.0, 1.0])\n",
    "plt.ylim([0.0, 1.05])\n",
    "plt.xlabel('False Positive Rate')\n",
    "plt.ylabel('True Positive Rate')\n",
    "plt.title('Receiver operating characteristics')\n",
    "plt.legend(loc=\"lower right\")\n",
    "plt.show()\n",
    "\n",
    "# computhe the cross-entropy loss score\n",
    "score = log_loss(Y_test,y_pred)\n",
    "print(score)\n",
    "\n",
    "# compute the average precision score\n",
    "prec_score = average_precision_score(Y_test,y_pred)  \n",
    "print(prec_score)\n",
    "\n",
    "# transfer it back\n",
    "y_pred = np.argmax(y_pred, axis=1)\n",
    "Y_test = np.argmax(Y_test, axis=1)\n",
    "print(y_pred)\n",
    "print(Y_test)\n",
    "\n",
    "#save the predictions as a CSV file for further analysis\n",
    "np.savetxt('vgg16_model_y_pred.csv',y_pred,fmt='%i',delimiter = \",\")\n",
    "np.savetxt('vgg16_model_Y_test.csv',Y_test,fmt='%i',delimiter = \",\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, let us plot the confusion matrix of the model's performance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "target_names = ['class 0(abnormal)', 'class 1(normal)'] #decide the labels for your own data\n",
    "print(classification_report(Y_test,y_pred,target_names=target_names))\n",
    "print(confusion_matrix(Y_test,y_pred))\n",
    "cnf_matrix = (confusion_matrix(Y_test,y_pred))\n",
    "np.set_printoptions(precision=4)\n",
    "\n",
    "# Plot non-normalized confusion matrix\n",
    "plt.figure(figsize=(20,10), dpi=100)\n",
    "plot_confusion_matrix(cnf_matrix, classes=target_names,\n",
    "                  title='Confusion matrix')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load image to be predicted\n",
    "img_path = 'C182P143NThinF_IMG_20151201_172216_cell_151.png' #path to the image\n",
    "img = image.load_img(img_path)\n",
    "\n",
    "#preprocess the image\n",
    "x = image.img_to_array(img)\n",
    "x = np.expand_dims(x, axis=0)\n",
    "\n",
    "#predict on the image\n",
    "preds = model.predict(x)[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#grad-CAM visualization and heat-map\n",
    "malaria_output = model.output[:, 0] #index of the positive (abnormal) class\n",
    "\n",
    "#deepest convolutional layer in the trained model\n",
    "model.summary()\n",
    "last_conv_layer = model.get_layer('block5_conv2')\n",
    "\n",
    "#compute the Gradient of the expected class (abnormal) with regard to the \n",
    "#output feature map of the deepst convolutional layer\n",
    "grads = K.gradients(malaria_output, last_conv_layer.output)[0]\n",
    "\n",
    "#Vector of shape (512,), where each entry is the mean intensity of the gradient over a specific feature-map channel\n",
    "pooled_grads = K.mean(grads, axis=(0, 1, 2))\n",
    "\n",
    "#Lets you access the values of the quantities you just defined: pooled_grads and the\n",
    "#output feature map of block5_conv2, given the input image\n",
    "iterate = K.function([model.input],[pooled_grads, last_conv_layer.output[0]])\n",
    "\n",
    "#Values of these two quantities, as Numpy arrays, given the input abnormal image\n",
    "pooled_grads_value, conv_layer_output_value = iterate([x])\n",
    "\n",
    "#Multiplies each channel in the feature-map array by “how important this channel is” with regard to the “abnormal” class\n",
    "for i in range(512):\n",
    "    conv_layer_output_value[:, :, i] *= pooled_grads_value[i]\n",
    "    \n",
    "#For visualization purposes, we normalize the heatmap between 0 and 1.\n",
    "heatmap = np.mean(conv_layer_output_value, axis=-1)\n",
    "heatmap = np.maximum(heatmap, 0)\n",
    "heatmap /= np.max(heatmap)\n",
    "plt.matshow(heatmap)\n",
    "\n",
    "# we OpenCV to generate an image that superimposes the original image\n",
    "#on the heatmap obtained\n",
    "img = cv2.imread(img_path)\n",
    "\n",
    "#Resizes the heatmap to be the same size as the original image\n",
    "heatmap = cv2.resize(heatmap, (img.shape[1], img.shape[0]))\n",
    "\n",
    "#Converts the heatmap to RGB \n",
    "heatmap = np.uint8(255 * heatmap)\n",
    "\n",
    "#Applies the heatmap to the original image\n",
    "heatmap = cv2.applyColorMap(heatmap, cv2.COLORMAP_JET)\n",
    "superimposed_img = heatmap * 0.4 + img #0.4 here is a heatmap intensity factor\n",
    "\n",
    "#Saves the image to disk\n",
    "cv2.imwrite('CAM_parasitic_cell.jpg', superimposed_img)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These visualization techniques answer two important questions: Why did the network think this is a parasitic cell? Where is the parasite seen in the cell image? In particular, it’s interesting to note that the areas showing the parasites are strongly activated: this is probably how the network can tell the difference between normal and abnormal (parasitic) cells."
   ]
  }
 ],
 "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.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}