{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!pip install keras-vis","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:55:30.473284Z","iopub.execute_input":"2021-07-17T17:55:30.473790Z","iopub.status.idle":"2021-07-17T17:55:39.388102Z","shell.execute_reply.started":"2021-07-17T17:55:30.473759Z","shell.execute_reply":"2021-07-17T17:55:39.387062Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\nimport pickle\n\nfrom tqdm import tqdm\nimport glob\nimport os\nimport matplotlib.pyplot as plt\nimport matplotlib.pylab as pylab\nimport seaborn as sns\nimport pprint\nimport pydicom as dicom\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\nimport albumentations as A \nimport cv2\n\nfrom sklearn.metrics import roc_curve,roc_auc_score, auc\nimport sklearn\n\nfrom PIL import Image\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nos.listdir('../input')\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:50:45.407694Z","iopub.execute_input":"2021-07-17T17:50:45.408298Z","iopub.status.idle":"2021-07-17T17:50:48.033100Z","shell.execute_reply.started":"2021-07-17T17:50:45.408209Z","shell.execute_reply":"2021-07-17T17:50:48.032174Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from kaggle_datasets import KaggleDatasets\nfrom sklearn.metrics import roc_curve\n\nimport tensorflow.keras as K\nimport tensorflow as tf\ntry:\n    import tensorflow_io as tfio\nexcept:\n    print('tensorflow_io not installed')\nimport re\nimport vis ## keras-vis\nfrom vis.utils import utils","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:55:39.389906Z","iopub.execute_input":"2021-07-17T17:55:39.390167Z","iopub.status.idle":"2021-07-17T17:55:39.526810Z","shell.execute_reply.started":"2021-07-17T17:55:39.390138Z","shell.execute_reply":"2021-07-17T17:55:39.525805Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## History analysis\nOra che ho fatto training per 100 epochs posso analizzare la history e il comportamento del modello","metadata":{}},{"cell_type":"code","source":"hf, wf = [], []\nfor root, dirs, files in os.walk('../input/siim-covid-pooled/'):\n    for file in files:\n        if 'history' in file:\n            hf.append(os.path.join(root,file))\n        elif 'weights' in file:\n            wf.append(os.path.join(root,file))\n            \nhf, wf","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:00.847832Z","iopub.execute_input":"2021-07-17T17:53:00.848185Z","iopub.status.idle":"2021-07-17T17:53:00.858292Z","shell.execute_reply.started":"2021-07-17T17:53:00.848152Z","shell.execute_reply":"2021-07-17T17:53:00.857165Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"hf = '../input/reti-siim-covid/resmodel_pooled_history_16_07_2021__20_40.pkl'\nwf = '../input/reti-siim-covid/resmodel_pooled_weights_16_07_2021__20_40.h5'","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:28.015781Z","iopub.execute_input":"2021-07-17T17:53:28.016183Z","iopub.status.idle":"2021-07-17T17:53:28.020839Z","shell.execute_reply.started":"2021-07-17T17:53:28.016147Z","shell.execute_reply":"2021-07-17T17:53:28.019819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"history = pickle.load(open(hf, 'rb'))\nmodel = K.models.load_model(wf)\nmodel.summary()","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:29.140602Z","iopub.execute_input":"2021-07-17T17:53:29.141109Z","iopub.status.idle":"2021-07-17T17:53:34.409496Z","shell.execute_reply.started":"2021-07-17T17:53:29.141076Z","shell.execute_reply":"2021-07-17T17:53:34.408510Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"history['val_accuracy'][-1]","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:38.177950Z","iopub.execute_input":"2021-07-17T17:53:38.178462Z","iopub.status.idle":"2021-07-17T17:53:38.184565Z","shell.execute_reply.started":"2021-07-17T17:53:38.178429Z","shell.execute_reply":"2021-07-17T17:53:38.183696Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(history.keys())\nfig, ax = plt.subplots(2, 2, figsize=[15, 8], gridspec_kw=dict(left=0.1, right=0.9, bottom=0, top=1))\n\n# summarize history for accuracy\nax[0,0].plot(history['accuracy'], label='train')\nax[0,0].plot(history['val_accuracy'], label='test')\nax[0,0].set_title('model accuracy')\nax[0,0].set_ylabel('accuracy')\nax[0,0].set_xlabel('epoch')\nax[0,0].legend(loc='upper left')\n# summarize history for loss\n\nax[0,1].plot(history['loss'], label='train')\nax[0,1].plot(history['val_loss'], label='test')\nax[0,1].set_title('model loss')\nax[0,1].set_ylabel('loss')\nax[0,1].set_xlabel('epoch')\nax[0,1].legend(loc='upper left')\n\nax[1,0].plot(history['auc'], label='train')\nax[1,0].plot(history['val_auc'], label='test')\nax[1,0].set_title('model auc')\nax[1,0].set_ylabel('auc')\nax[1,0].set_xlabel('epoch')\nax[1,0].legend(loc='upper left')\n\nax[1,1].plot(history['prc'], label='train')\nax[1,1].plot(history['val_prc'], label='test')\nax[1,1].set_title('model prc')\nax[1,1].set_ylabel('prc')\nax[1,1].set_xlabel('epoch')\nax[1,1].legend(loc='upper left')\n\nplt.savefig('stats.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:39.570617Z","iopub.execute_input":"2021-07-17T17:53:39.570978Z","iopub.status.idle":"2021-07-17T17:53:40.718464Z","shell.execute_reply.started":"2021-07-17T17:53:39.570947Z","shell.execute_reply":"2021-07-17T17:53:40.717478Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x_train = np.load('../input/reti-siim-covid/x_train.npy')\ny_train = np.load('../input/reti-siim-covid/y_train.npy')","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:40.720083Z","iopub.execute_input":"2021-07-17T17:53:40.720617Z","iopub.status.idle":"2021-07-17T17:53:48.417441Z","shell.execute_reply.started":"2021-07-17T17:53:40.720576Z","shell.execute_reply":"2021-07-17T17:53:48.416508Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dist = [len(y_train[y_train[:, i] == 1]) for i in range(4)]\nlabels = ['Negative for Pneumonia', 'Typical Appearance', 'Indeterminate Appearance', 'Atypical Appearance']\n\nfig, ax = plt.subplots(1,2, figsize=[15, 7])\n\nsns.barplot(x=[l.split()[0] for l in labels], y=dist, ax=ax[0])\nax[0].set_title('Category distribution in the training dataset')\nax[1].imshow(x_train[0])\nax[1].set_title('Example image')\nplt.savefig('cat_distribution.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:48.419272Z","iopub.execute_input":"2021-07-17T17:53:48.419675Z","iopub.status.idle":"2021-07-17T17:53:48.984392Z","shell.execute_reply.started":"2021-07-17T17:53:48.419643Z","shell.execute_reply":"2021-07-17T17:53:48.983239Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"VAL_SPLIT = 0.15\nBATCH_SIZE = 200\ndatagen = K.preprocessing.image.ImageDataGenerator(rotation_range=30,\n                                                   width_shift_range=0.1,\n                                                   height_shift_range=0.1,\n                                                   brightness_range=(0.8, 1.2),\n                                                   shear_range=15,\n                                                   horizontal_flip=True,\n                                                   vertical_flip=True,\n                                                   validation_split=VAL_SPLIT,\n                                                  )\n\n\ntrain_it = datagen.flow(x_train,y_train, batch_size=BATCH_SIZE,subset='training')\n\nvalidation_it = datagen.flow(x_train, y_train, batch_size=BATCH_SIZE, subset='validation')","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:48.986403Z","iopub.execute_input":"2021-07-17T17:53:48.986945Z","iopub.status.idle":"2021-07-17T17:53:50.216076Z","shell.execute_reply.started":"2021-07-17T17:53:48.986898Z","shell.execute_reply":"2021-07-17T17:53:50.215179Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x_pred, y_true = validation_it.next()\ny_true = np.array([np.argmax(x) for x in y_true])\n\npred = model.predict(x_pred)\ny_val = np.array([np.argmax(x) for x in pred])\ncm = sklearn.metrics.confusion_matrix(y_true, y_val, labels=[0,1,2,3], normalize='pred')","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:53:50.217385Z","iopub.execute_input":"2021-07-17T17:53:50.217667Z","iopub.status.idle":"2021-07-17T17:54:07.112188Z","shell.execute_reply.started":"2021-07-17T17:53:50.217641Z","shell.execute_reply":"2021-07-17T17:54:07.111022Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize = (10,7))\nsns.heatmap(cm, annot=True, cmap=\"YlGnBu\")\nplt.savefig('confusion_matrix.png')","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:54:07.116897Z","iopub.execute_input":"2021-07-17T17:54:07.117301Z","iopub.status.idle":"2021-07-17T17:54:07.507079Z","shell.execute_reply.started":"2021-07-17T17:54:07.117263Z","shell.execute_reply":"2021-07-17T17:54:07.505047Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"x_pred, y_true = validation_it.next()\n\ny_pred = model.predict(x_pred)","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:54:07.509475Z","iopub.execute_input":"2021-07-17T17:54:07.509885Z","iopub.status.idle":"2021-07-17T17:54:23.078452Z","shell.execute_reply.started":"2021-07-17T17:54:07.509842Z","shell.execute_reply":"2021-07-17T17:54:23.077308Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"n_classes = 4\n\n# Compute ROC curve and ROC area for each class\nfpr = dict()\ntpr = dict()\nroc_auc = dict()\nfor i in range(n_classes):\n    fpr[i], tpr[i], _ = roc_curve(y_true[:, i], y_pred[:, i])\n    roc_auc[i] = auc(fpr[i], tpr[i])\n\n# Plot of a ROC curve for a specific class\n\nfig, ax = plt.subplots(2, 2, figsize=[15, 8], gridspec_kw=dict(left=0.1, right=0.9, bottom=0, top=1))\nfor i in range(n_classes):\n    ax_t = ax[i//2, i%2]\n    ax_t.plot(fpr[i], tpr[i], label='ROC curve (area = %0.2f)' % roc_auc[i])\n    ax_t.plot([0, 1], [0, 1], 'k--')\n    ax_t.set_xlabel('False Positive Rate')\n    ax_t.set_ylabel('True Positive Rate')\n    ax_t.set_title(f'ROC for class {labels[i]}')\n    ax_t.legend(loc=\"lower right\")\n    \nplt.savefig('roc.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-07-17T17:54:23.080063Z","iopub.execute_input":"2021-07-17T17:54:23.080330Z","iopub.status.idle":"2021-07-17T17:54:23.929531Z","shell.execute_reply.started":"2021-07-17T17:54:23.080304Z","shell.execute_reply":"2021-07-17T17:54:23.928773Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"single_it = datagen.flow(x_train, y_train, batch_size=1, subset='validation')\nimg, y = single_it.next()\nclass_idxs_sorted = np.argsort(y.flatten())[::-1]\nclasslabel = labels","metadata":{"execution":{"iopub.status.busy":"2021-07-17T18:19:04.623988Z","iopub.execute_input":"2021-07-17T18:19:04.624385Z","iopub.status.idle":"2021-07-17T18:19:04.785936Z","shell.execute_reply.started":"2021-07-17T18:19:04.624348Z","shell.execute_reply":"2021-07-17T18:19:04.785127Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Utility to search for layer index by name. \n# Alternatively we can specify this as -1 since it corresponds to the last layer.\n\nlayer_idx = -1\nmodel.layers[layer_idx].activation = K.activations.linear\nlinear_model = utils.apply_modifications(model)\nclass_idx = class_idxs_sorted[0]\nclass_idx = class_idxs_sorted[0]\n\ndef get_saliency_map(model, img, class_idx):\n    image = tf.convert_to_tensor(img, np.float32)\n    with tf.GradientTape() as tape:\n        tape.watch(image)\n        predictions = model(image)\n        \n        loss = predictions[:, class_idx]\n    \n    # Get the gradients of the loss w.r.t to the input image.\n    gradient = tape.gradient(loss, image)\n    \n    # take maximum across channels\n    gradient = tf.reduce_max(gradient, axis=-1)\n    \n    # convert to numpy\n    gradient = gradient.numpy()\n    \n    # normaliz between 0 and 1\n    min_val, max_val = np.min(gradient), np.max(gradient)\n    smap = (gradient - min_val) / (max_val - min_val + K.backend.epsilon())\n    \n    return smap[0]\n\ndef plot_map(grads, image):\n    img = image[0, :, : , 0]\n    fig, axes = plt.subplots(1,2,figsize=(14,5))\n    axes[0].imshow(img, cmap=plt.cm.gray)\n    axes[1].imshow(img, cmap=plt.cm.gray)\n    i = axes[1].imshow(grads,cmap=\"jet\",alpha=0.5)\n    fig.colorbar(i)\n    plt.suptitle(\"Pr(class={}) = {:5.2f}\".format(\n                      classlabel[class_idx],\n                      y_pred[0,class_idx]))","metadata":{"execution":{"iopub.status.busy":"2021-07-17T18:19:22.344291Z","iopub.execute_input":"2021-07-17T18:19:22.344641Z","iopub.status.idle":"2021-07-17T18:19:25.428890Z","shell.execute_reply.started":"2021-07-17T18:19:22.344611Z","shell.execute_reply":"2021-07-17T18:19:25.427935Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for class_idx in class_idxs_sorted:\n    plot_map(get_saliency_map(model, img, class_idx), img)","metadata":{"execution":{"iopub.status.busy":"2021-07-17T18:19:09.433046Z","iopub.execute_input":"2021-07-17T18:19:09.433391Z","iopub.status.idle":"2021-07-17T18:19:14.001960Z","shell.execute_reply.started":"2021-07-17T18:19:09.433361Z","shell.execute_reply":"2021-07-17T18:19:14.000856Z"},"trusted":true},"execution_count":null,"outputs":[]}]}