custom inception cnn.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This code is meant to train and evaluate a custom Inception CNN model toward detecting pneumonia and differentiating bacterial and viral pneumonia in pediatric chest radiographs. You can use your own custom data to train and validate these models."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load libraries\n",
    "from __future__ import print_function\n",
    "\n",
    "import numpy as np\n",
    "from keras.utils import plot_model\n",
    "from keras.models import model_from_json\n",
    "from keras.models import Sequential\n",
    "from keras.models import Model\n",
    "from keras import layers\n",
    "from keras.layers import Conv2D, Dense, MaxPooling2D, GlobalAveragePooling2D, AveragePooling2D, Concatenate\n",
    "from keras.layers import Input\n",
    "import time\n",
    "from keras.optimizers import SGD\n",
    "from sklearn.metrics import roc_curve, auc\n",
    "from sklearn.preprocessing import label_binarize\n",
    "from scipy import interp\n",
    "from itertools import cycle\n",
    "from sklearn.metrics import classification_report,confusion_matrix,accuracy_score\n",
    "from keras.preprocessing.image import ImageDataGenerator\n",
    "import matplotlib.pyplot as plt\n",
    "from evaluation import plot_confusion_matrix\n",
    "from keras.callbacks import ModelCheckpoint, TensorBoard\n",
    "from keras.utils import plot_model\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Loading the data\n",
    "train_data_dir = 'cxr_normal_pneumonia_1024/train'\n",
    "test_data_dir = 'cxr_normal_pneumonia_1024/test'\n",
    "img_width = 1024\n",
    "img_height = 1024\n",
    "channel = 3\n",
    "input_img = Input(shape = (img_width, img_height, channel))\n",
    "epochs = 100\n",
    "batch_size = 8 #vary this parameter depending on your GPU capacity\n",
    "num_classes= 2 #[pneumonia, normal] [bacterial, viral]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Load the InceptionV3 based model layers.\n",
    "tower_1a = Conv2D(4, (1,1), padding='same', activation='relu', strides=2)(input_img)\n",
    "tower_2a = Conv2D(4, (1,1), padding='same', activation='relu')(input_img)\n",
    "tower_2a = Conv2D(4, (3,3), padding='same', activation='relu', strides=2)(tower_2a)\n",
    "tower_3a = AveragePooling2D((3,3), strides=(2,2), padding='same')(input_img)\n",
    "tower_3a = Conv2D(4, (3,3), padding='same', activation='relu')(tower_3a)\n",
    "tower_4a = Conv2D(4, (1,1), padding='same', activation='relu')(input_img)\n",
    "tower_4a = Conv2D(4, (3,3), padding='same', activation='relu')(tower_4a)\n",
    "tower_4a = Conv2D(4, (3,3), padding='same', activation='relu', strides=2)(tower_4a)\n",
    "output1 = layers.concatenate([tower_1a, tower_2a, tower_3a, tower_4a], axis = -1)\n",
    "\n",
    "#module2\n",
    "tower_1b = Conv2D(8, (1,1), padding='same', activation='relu', strides=2)(output1)\n",
    "tower_2b = Conv2D(8, (1,1), padding='same', activation='relu')(output1)\n",
    "tower_2b = Conv2D(8, (3,3), padding='same', activation='relu', strides=2)(tower_2b)\n",
    "tower_3b = AveragePooling2D((3,3), strides=(2,2), padding='same')(output1)\n",
    "tower_3b = Conv2D(8, (3,3), padding='same', activation='relu')(tower_3b)\n",
    "tower_4b = Conv2D(8, (1,1), padding='same', activation='relu')(output1)\n",
    "tower_4b = Conv2D(8, (3,3), padding='same', activation='relu')(tower_4b)\n",
    "tower_4b = Conv2D(8, (3,3), padding='same', activation='relu', strides=2)(tower_4b)\n",
    "output2 = layers.concatenate([tower_1b, tower_2b, tower_3b, tower_4b], axis = -1)\n",
    "\n",
    "#module3\n",
    "tower_1c = Conv2D(16, (1,1), padding='same', activation='relu', strides=2)(output2)\n",
    "tower_2c = Conv2D(16, (1,1), padding='same', activation='relu')(output2)\n",
    "tower_2c = Conv2D(16, (3,3), padding='same', activation='relu', strides=2)(tower_2c)\n",
    "tower_3c = AveragePooling2D((3,3), strides=(2,2), padding='same')(output2)\n",
    "tower_3c = Conv2D(16, (3,3), padding='same', activation='relu')(tower_3c)\n",
    "tower_4c = Conv2D(16, (1,1), padding='same', activation='relu')(output2)\n",
    "tower_4c = Conv2D(16, (3,3), padding='same', activation='relu')(tower_4c)\n",
    "tower_4c = Conv2D(16, (3,3), padding='same', activation='relu', strides=2)(tower_4c)\n",
    "output3 = layers.concatenate([tower_1c, tower_2c, tower_3c, tower_4c], axis = -1)\n",
    "\n",
    "#module4\n",
    "tower_1d = Conv2D(32, (1,1), padding='same', activation='relu', strides=2)(output3)\n",
    "tower_2d = Conv2D(32, (1,1), padding='same', activation='relu')(output3)\n",
    "tower_2d = Conv2D(32, (3,3), padding='same', activation='relu', strides=2)(tower_2d)\n",
    "tower_3d = AveragePooling2D((3,3), strides=(2,2), padding='same')(output3)\n",
    "tower_3d = Conv2D(32, (3,3), padding='same', activation='relu')(tower_3d)\n",
    "tower_4d = Conv2D(32, (1,1), padding='same', activation='relu')(output3)\n",
    "tower_4d = Conv2D(32, (3,3), padding='same', activation='relu')(tower_4d)\n",
    "tower_4d = Conv2D(32, (3,3), padding='same', activation='relu', strides=2)(tower_4d)\n",
    "output4 = layers.concatenate([tower_1d, tower_2d, tower_3d, tower_4d], axis = -1)\n",
    "\n",
    "#module5\n",
    "tower_1e = Conv2D(64, (1,1), padding='same', activation='relu', strides=2)(output4)\n",
    "tower_2e = Conv2D(64, (1,1), padding='same', activation='relu')(output4)\n",
    "tower_2e = Conv2D(64, (3,3), padding='same', activation='relu', strides=2)(tower_2e)\n",
    "tower_3e = AveragePooling2D((3,3), strides=(2,2), padding='same')(output4)\n",
    "tower_3e = Conv2D(64, (3,3), padding='same', activation='relu')(tower_3e)\n",
    "tower_4e = Conv2D(64, (1,1), padding='same', activation='relu')(output4)\n",
    "tower_4e = Conv2D(64, (3,3), padding='same', activation='relu')(tower_4e)\n",
    "tower_4e = Conv2D(64, (3,3), padding='same', activation='relu', strides=2)(tower_4e)\n",
    "output5 = layers.concatenate([tower_1e, tower_2e, tower_3e, tower_4e], axis = -1)\n",
    "\n",
    "#module6\n",
    "tower_1f = Conv2D(128, (1,1), padding='same', activation='relu', strides=2)(output5)\n",
    "tower_2f = Conv2D(128, (1,1), padding='same', activation='relu')(output5)\n",
    "tower_2f = Conv2D(128, (3,3), padding='same', activation='relu', strides=2)(tower_2f)\n",
    "tower_3f = AveragePooling2D((3,3), strides=(2,2), padding='same')(output5)\n",
    "tower_3f = Conv2D(128, (3,3), padding='same', activation='relu')(tower_3f)\n",
    "tower_4f = Conv2D(128, (1,1), padding='same', activation='relu')(output5)\n",
    "tower_4f = Conv2D(128, (3,3), padding='same', activation='relu')(tower_4f)\n",
    "tower_4f = Conv2D(128, (3,3), padding='same', activation='relu', strides=2)(tower_4f)\n",
    "\n",
    "#concatenating the outputs\n",
    "output6 = layers.concatenate([tower_1f, tower_2f, tower_3f, tower_4f], axis = -1)\n",
    "\n",
    "#adding global average pooling layer\n",
    "output = GlobalAveragePooling2D()(output6)\n",
    "\n",
    "#adding the final dense layer\n",
    "out = Dense(num_classes, activation='softmax', name='predictions')(output)\n",
    "\n",
    "#creating the model\n",
    "model = Model(inputs = input_img, outputs = out)\n",
    "\n",
    "#generating model summary\n",
    "model.summary()\n",
    "\n",
    "#plot the model\n",
    "plot_model(model, to_file='custom_inceptionv3_model.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We performed five-fold cross validation in this study. We have however shown running the script with a sample data fold. We performed no augmentation other than rescaling.We allocated 20% of the training data at random for validation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#declaring image data generators\n",
    "train_datagen = ImageDataGenerator(\n",
    "        rescale=1./255,\n",
    "        validation_split=0.2) #taking 20% of training for validation \n",
    "\n",
    "test_datagen = ImageDataGenerator(rescale=1./255)\n",
    "\n",
    "train_generator = train_datagen.flow_from_directory(\n",
    "        train_data_dir,\n",
    "        target_size=(1024, 1024),\n",
    "        batch_size=batch_size,\n",
    "        class_mode='categorical',\n",
    "        subset='training')\n",
    "\n",
    "validation_generator = train_datagen.flow_from_directory(\n",
    "        train_data_dir,\n",
    "        target_size=(1024, 1024),\n",
    "        batch_size=batch_size,\n",
    "        class_mode='categorical',\n",
    "        subset='validation')\n",
    "\n",
    "test_generator = test_datagen.flow_from_directory(\n",
    "        test_data_dir,\n",
    "        target_size=(1024, 1024),\n",
    "        batch_size=batch_size,\n",
    "        class_mode='categorical',shuffle=False)\n",
    "\n",
    "#count the number of samples\n",
    "nb_train_samples = len(train_generator.filenames)\n",
    "nb_validation_samples = len(validation_generator.filenames)\n",
    "nb_test_samples = len(test_generator.filenames)\n",
    "\n",
    "#check the class indices\n",
    "train_generator.class_indices\n",
    "validation_generator.class_indices\n",
    "test_generator.class_indices"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Allocate balanced weights to penalize over-represented classes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.utils import class_weight\n",
    "class_weights = class_weight.compute_class_weight(\n",
    "               'balanced',\n",
    "                np.unique(train_generator.classes), \n",
    "                train_generator.classes)\n",
    "print(class_weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#fix the optimizer\n",
    "sgd = SGD(lr=0.01, decay=1e-4, momentum=0.95, nesterov=True) #optimize to your requirements\n",
    "\n",
    "#compile the model\n",
    "model.compile(optimizer=sgd,\n",
    "              loss='categorical_crossentropy',\n",
    "              metrics=['accuracy'])\n",
    "\n",
    "#compute the training time\n",
    "start = time.time()\n",
    "\n",
    "#give the path to store the model weights\n",
    "filepath = 'weights/' + model.name + '.{epoch:02d}-{val_acc:.4f}.hdf5'\n",
    "\n",
    "#save only the best weights\n",
    "checkpoint = ModelCheckpoint(filepath, monitor='val_acc', verbose=1, \n",
    "                             save_weights_only=True, save_best_only=True, mode='max', period=1)\n",
    "\n",
    "#visualize model performance using TensorBoard\n",
    "tensor_board = TensorBoard(log_dir='logs/', histogram_freq=0, batch_size=batch_size)\n",
    "callbacks_list = [checkpoint, tensor_board]\n",
    "\n",
    "#train model\n",
    "custom_inception_history = model.fit_generator(\n",
    "      train_generator,\n",
    "      steps_per_epoch=nb_train_samples // batch_size, \n",
    "      epochs=epochs,\n",
    "      validation_data=validation_generator,\n",
    "      callbacks=callbacks_list,\n",
    "      class_weight = class_weights,\n",
    "      validation_steps=nb_validation_samples // batch_size, \n",
    "      verbose=1)\n",
    "\n",
    "#print the total time taken for training\n",
    "print(time.time()-start)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Testing the model's performance by loading the model weights\n",
    "\n",
    "model.load_weights('weights/model_1.01-0.7423.hdf5') #change this to your path and model weights\n",
    "custom_inception_y_pred = model.predict_generator(test_generator, nb_test_samples/batch_size, workers=1)\n",
    "\n",
    "#true labels\n",
    "Y_test=test_generator.classes\n",
    "\n",
    "#print the shape of y_pred and Y_test\n",
    "print(custom_inception_y_pred.shape)\n",
    "print(Y_test.shape)\n",
    "\n",
    "#measure accuracy\n",
    "custom_inception_model_accuracy=accuracy_score(Y_test,custom_inception_y_pred.argmax(axis=-1))\n",
    "print('The accuracy of custom Inception model is: ', custom_inception_model_accuracy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#print classification report\n",
    "target_names = ['class 0(normal)', 'class 1(pneumonia)'] #from the generator.class_indices\n",
    "print(classification_report(Y_test,custom_inception_y_pred.argmax(axis=-1),\n",
    "                                                              target_names=target_names, digits=4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute confusion matrix\n",
    "cnf_matrix = confusion_matrix(Y_test,custom_inception_y_pred.argmax(axis=-1))\n",
    "np.set_printoptions(precision=4)\n",
    "plt.figure(figsize=(20,10), dpi=100)\n",
    "plot_confusion_matrix(cnf_matrix, classes=target_names,\n",
    "                      title='Confusion matrix for custom Inception model without normalization')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#store the predictions in .csv files\n",
    "print(custom_inception_y_pred.argmax(axis=-1))\n",
    "print(Y_test)\n",
    "\n",
    "np.savetxt('custom_Inception_y_pred.csv',custom_inception_y_pred.argmax(axis=-1),fmt='%i',delimiter = \",\")\n",
    "np.savetxt('Y_test.csv',Y_test,fmt='%i',delimiter = \",\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "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, custom_inception_y_pred[:, i])\n",
    "    roc_auc[i] = auc(fpr[i], tpr[i])\n",
    "fig=plt.figure(figsize=(15,10), dpi=100)\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "major_ticks = np.arange(0.0, 1.0, 0.10)\n",
    "minor_ticks = np.arange(0.0, 1.0, 0.10)\n",
    "ax.set_xticks(major_ticks)\n",
    "ax.set_xticks(minor_ticks, minor=True)\n",
    "ax.set_yticks(major_ticks)\n",
    "ax.set_yticks(minor_ticks, minor=True)\n",
    "ax.grid(which='both')\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.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()"
   ]
  }
 ],
 "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
}