{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"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\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\n# import os\n# for dirname, _, filenames in os.walk('/kaggle/input'):\n#     for filename in filenames:\n#         print(os.path.join(dirname, filename))\n\n# You can write up to 5GB 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","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Import \n\nThis notebook is taken from https://www.kaggle.com/seraphwedd18/pe-detection-with-keras-model-creation/output?select=submission.csv"},{"metadata":{"trusted":true},"cell_type":"code","source":"import os\nimport gc\nimport time\nfrom IPython.display import clear_output\nfrom tensorflow.keras.models import load_model\nfrom tensorflow.keras.callbacks import ModelCheckpoint as MC\nfrom tensorflow.keras import backend as K\n\nroot = \"/kaggle/input/rsna-str-pulmonary-embolism-detection/\"\n\nfor item in os.listdir(root):\n    path = os.path.join(root, item)\n    if os.path.isfile(path):\n        print(path)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Load Dataset"},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Reading train data...')\ntrain = pd.read_csv(\"../input/rsna-str-pulmonary-embolism-detection/train.csv\")\nprint(train.shape)\ntrain.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Reading test data...')\ntest = pd.read_csv(\"../input/rsna-str-pulmonary-embolism-detection/test.csv\")\nprint(test.shape)\ntest.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Reading sample data...')\nss = pd.read_csv(\"../input/rsna-str-pulmonary-embolism-detection/sample_submission.csv\")\nprint(ss.shape)\nss.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Check Targets and Input Image"},{"metadata":{"trusted":true},"cell_type":"code","source":"ids = ss.id\n\ncounter = [1 for _ in range(10)]\nmapper = []\n\nfor i in ids:\n    n = '_'.join(i.split('_')[1:])\n    if n not in mapper:\n        mapper.append(n)\n    else:\n        counter[mapper.index(n)] += 1\nprint(\"List of keys:\")\nprint(mapper, sep='\\n')\nprint()\nprint(\"Count of items per key:\")\nprint(counter)\n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import vtk\nfrom vtk.util import numpy_support\nimport cv2\n\nreader = vtk.vtkDICOMImageReader()\ndef get_img(path):\n    reader.SetFileName(path)\n    reader.Update()\n    _extent = reader.GetDataExtent()\n    ConstPixelDims = [_extent[1]-_extent[0]+1, _extent[3]-_extent[2]+1, _extent[5]-_extent[4]+1]\n\n    ConstPixelSpacing = reader.GetPixelSpacing()\n    imageData = reader.GetOutput()\n    pointData = imageData.GetPointData()\n    arrayData = pointData.GetArray(0)\n    ArrayDicom = numpy_support.vtk_to_numpy(arrayData)\n    ArrayDicom = ArrayDicom.reshape(ConstPixelDims, order='F')\n    ArrayDicom = cv2.resize(ArrayDicom,(512,512))\n    return ArrayDicom","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#test read a dcom file and view it\nfpath = \"../input/rsna-str-pulmonary-embolism-detection/train/0003b3d648eb/d2b2960c2bbf/00ac73cfc372.dcm\"\nds = get_img(fpath)\n\nimport matplotlib.pyplot as plt\n\n#Convert dcom file to 8bit color\nfunc = lambda x: int((2**15 + x)*(255/2**16))\nint16_to_uint8 = np.vectorize(func)\n\ndef show_dicom_images(dcom):\n    f, ax = plt.subplots(1,2, figsize=(16,20))\n    data_row_img = int16_to_uint8(ds)\n    ax[0].imshow(data_row_img, cmap=plt.cm.bone)\n    ax[1].imshow(ds, cmap=plt.cm.bone)\n    #print(data_row_img)\n    ax[0].axis('off')\n    ax[0].set_title('8-bit DICOM Image')\n    ax[1].axis('off')\n    ax[1].set_title('16-bit DICOM Image')\n    plt.show()\n    \nshow_dicom_images(ds)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Model Creation"},{"metadata":{"trusted":true},"cell_type":"code","source":"import tensorflow as tf\nfrom tensorflow import keras\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.layers import Input, Dense, Dropout, Conv2D\n\ninputs = Input((512, 512, 3))\n\n#x = Conv2D(3, (1, 1), activation='relu')(inputs)\n\nbase_model = keras.applications.Xception(\n    include_top=False,\n    weights=\"imagenet\"\n)\n\nbase_model.trainable = False\n\noutputs = base_model(inputs, training=False)\noutputs = keras.layers.GlobalAveragePooling2D()(outputs)\n# outputs = keras.layers.MaxPooling2D()(outputs)\n\noutputs = Dropout(0.4)(outputs)\n# outputs = keras.layers.MaxPooling2D((3))(outputs)\n# outputs = keras.layers.GlobalAveragePooling2D()(outputs)\n# outputs = Dropout(0.3)(outputs)\n# outputs = Dense(2048, activation='relu')(outputs)\noutputs = Dense(1024, activation='relu')(outputs)\noutputs = Dense(256, activation='relu')(outputs)\noutputs = Dense(64, activation='relu')(outputs)\noutputs = Dense(32, activation='relu')(outputs)\nnepe = Dense(1, activation='sigmoid', name='negative_exam_for_pe')(outputs)\nrlrg1 = Dense(1, activation='sigmoid', name='rv_lv_ratio_gte_1')(outputs)\nrlrl1 = Dense(1, activation='sigmoid', name='rv_lv_ratio_lt_1')(outputs) \nlspe = Dense(1, activation='sigmoid', name='leftsided_pe')(outputs)\ncpe = Dense(1, activation='sigmoid', name='chronic_pe')(outputs)\nrspe = Dense(1, activation='sigmoid', name='rightsided_pe')(outputs)\naacpe = Dense(1, activation='sigmoid', name='acute_and_chronic_pe')(outputs)\ncnpe = Dense(1, activation='sigmoid', name='central_pe')(outputs)\nindt = Dense(1, activation='sigmoid', name='indeterminate')(outputs)\n\nmodel = Model(inputs=inputs, outputs={'negative_exam_for_pe':nepe,\n                                      'rv_lv_ratio_gte_1':rlrg1,\n                                      'rv_lv_ratio_lt_1':rlrl1,\n                                      'leftsided_pe':lspe,\n                                      'chronic_pe':cpe,\n                                      'rightsided_pe':rspe,\n                                      'acute_and_chronic_pe':aacpe,\n                                      'central_pe':cnpe,\n                                      'indeterminate':indt})\n\nopt = keras.optimizers.Adam(lr=0.001)\n\nmodel.compile(optimizer=opt,\n              loss='binary_crossentropy',\n              metrics=['accuracy'])\n\nmodel.summary()\nmodel.save('pe_detection_model.h5')\ndel model\nK.clear_session()\ngc.collect()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Training Model"},{"metadata":{"trusted":true},"cell_type":"code","source":"def convert_to_rgb(array):\n    array = array.reshape((512, 512, 1))\n    return np.stack([array, array, array], axis=2).reshape((512, 512, 3))\n    \ndef custom_dcom_image_generator(batch_size, dataset, test=False, debug=False):\n    \n    fnames = dataset[['StudyInstanceUID', 'SeriesInstanceUID', 'SOPInstanceUID']]\n    \n    if not test:\n        Y = dataset[['negative_exam_for_pe', 'rv_lv_ratio_gte_1', 'rv_lv_ratio_lt_1',\n                     'leftsided_pe', 'chronic_pe', 'rightsided_pe',\n                     'acute_and_chronic_pe', 'central_pe', 'indeterminate']]\n        prefix = 'input/rsna-str-pulmonary-embolism-detection/train'\n        \n    else:\n        prefix = 'input/rsna-str-pulmonary-embolism-detection/test'\n    \n    X = []\n    batch = 0\n    for st, sr, so in fnames.values:\n        if debug:\n            print(f\"Current file: ../{prefix}/{st}/{sr}/{so}.dcm\")\n\n        dicom = get_img(f\"../{prefix}/{st}/{sr}/{so}.dcm\")\n        image = convert_to_rgb(dicom)\n        X.append(image)\n        \n        del st, sr, so\n        \n        if len(X) == batch_size:\n            if test:\n                yield np.array(X)\n                del X\n            else:\n                yield np.array(X), Y[batch*batch_size:(batch+1)*batch_size].values\n                del X\n                \n            gc.collect()\n            X = []\n            batch += 1\n        \n    if test:\n        yield np.array(X)\n    else:\n        yield np.array(X), Y[batch*batch_size:(batch+1)*batch_size].values\n        del Y\n    del X\n    gc.collect()\n    return","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"history = {}\nstart = time.time()\ndebug = 0\nbatch_size = 500\ntrain_size = int(batch_size*0.8)\n\nmax_train_time = 3600 * 2 #hours to seconds of training\n\ncheckpoint = MC(filepath='../working/pe_detection_model.h5', monitor='val_loss', save_best_only=True, verbose=1)\n#Train loop\nfor n, (x, y) in enumerate(custom_dcom_image_generator(batch_size, train.sample(frac=1), False, debug)):\n    \n    if len(x) < 10: #Tries to filter out empty or short data\n        break\n        \n    clear_output(wait=True)\n    print(\"Training batch: %i - %i\" %(batch_size*n, batch_size*(n+1)))\n    \n    model = load_model('../working/pe_detection_model.h5')\n    hist = model.fit(\n        x[:train_size], #Y values are in a dict as there's more than one target for training output\n        {'negative_exam_for_pe':y[:train_size, 0],\n         'rv_lv_ratio_gte_1':y[:train_size, 1],\n         'rv_lv_ratio_lt_1':y[:train_size, 2],\n         'leftsided_pe':y[:train_size, 3],\n         'chronic_pe':y[:train_size, 4],\n         'rightsided_pe':y[:train_size, 5],\n         'acute_and_chronic_pe':y[:train_size, 6],\n         'central_pe':y[:train_size, 7],\n         'indeterminate':y[:train_size, 8]},\n\n        callbacks = checkpoint,\n\n        validation_split=0.2,\n        epochs=5,\n        batch_size=32,\n        verbose=debug\n    )\n    \n    print(\"Metrics for batch validation:\")\n    model.evaluate(x[train_size:],\n                   {'negative_exam_for_pe':y[train_size:, 0],\n                    'rv_lv_ratio_gte_1':y[train_size:, 1],\n                    'rv_lv_ratio_lt_1':y[train_size:, 2],\n                    'leftsided_pe':y[train_size:, 3],\n                    'chronic_pe':y[train_size:, 4],\n                    'rightsided_pe':y[train_size:, 5],\n                    'acute_and_chronic_pe':y[train_size:, 6],\n                    'central_pe':y[train_size:, 7],\n                    'indeterminate':y[train_size:, 8]\n                   }\n                  )\n    \n    try:\n        for key in hist.history.keys():\n            history[key] = np.concatenate([history[key], hist.history[key]], axis=0)\n    except:\n        for key in hist.history.keys():\n            history[key] = hist.history[key]\n            \n    #To make sure that our model don't train overtime\n    if time.time() - start >= max_train_time:\n        print(\"Time's up!\")\n        break\n        \n    model.save('pe_detection_model.h5')\n    del model, x, y, hist\n    K.clear_session()\n    gc.collect()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for key in history.keys():\n    if key.startswith('val'):\n        continue\n    else:\n        epoch = range(len(history[key]))\n        plt.plot(epoch, history[key]) #X=epoch, Y=value\n        plt.plot(epoch, history['val_'+key])\n        plt.title(key)\n        if 'accuracy' in key:\n            plt.axis([0, len(history[key]), -0.1, 1.1]) #Xmin, Xmax, Ymin, Ymax\n        plt.legend(['train', 'validation'], loc='upper right')\n        plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Prediction"},{"metadata":{"trusted":true},"cell_type":"code","source":"from tensorflow.keras import backend as K\n\npredictions = {}\nstopper = 3600 * 2 #4 hours limit for prediction\npred_start_time = time.time()\n\np, c = time.time(), time.time()\nbatch_size = 500\n    \nl = 0\nn = test.shape[0]\n\nfor x in custom_dcom_image_generator(batch_size, test, True, False):\n    clear_output(wait=True)\n    model = load_model(\"../working/pe_detection_model.h5\")\n    preds = model.predict(x, batch_size=8, verbose=1)\n    \n    try:\n        for key in preds.keys():\n            predictions[key] += preds[key].flatten().tolist()\n            \n    except Exception as e:\n        print(e)\n        for key in preds.keys():\n            predictions[key] = preds[key].flatten().tolist()\n            \n    l = (l+batch_size)%n\n    print('Total predicted:', len(predictions['indeterminate']),'/', n)\n    p, c = c, time.time()\n    print(\"One batch time: %.2f seconds\" %(c-p))\n    print(\"ETA: %.2f\" %((n-l)*(c-p)/batch_size))\n    \n    if c - pred_start_time >= stopper:\n        print(\"Time's up!\")\n        break\n    \n    del model\n    K.clear_session()\n    \n    del x, preds\n    gc.collect()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for key in predictions.keys():\n    print(key, np.array(predictions[key]).shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_ids = []\nfor v in test.StudyInstanceUID:\n    if v not in test_ids:\n        test_ids.append(v)\n        \ntest_preds = test.copy()\ntest_preds = pd.concat([test_preds, pd.DataFrame(predictions)], axis=1)\ntest_preds","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"IDS = []\nlabels = []\n\nfor label in ['negative_exam_for_pe', 'rv_lv_ratio_gte_1', 'rv_lv_ratio_lt_1',\n                 'leftsided_pe', 'chronic_pe', 'rightsided_pe',\n                 'acute_and_chronic_pe', 'central_pe', 'indeterminate']:\n    for key in test_ids:\n        temp = test_preds.loc[test_preds.StudyInstanceUID==key]\n        \n        IDS.append('_'.join([key, label]))\n        labels.append(np.max(temp[label]))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"IDS += test_preds.SOPInstanceUID.tolist()\nlabels += test_preds['negative_exam_for_pe'].tolist()\n\nsub = pd.DataFrame({\"id\":IDS, 'label':labels})\nsub","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.fillna(0.2, inplace=True)\nsub.to_csv('submission.csv', index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":4}