{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"markdown","source":"![](http://)# **RSNA Pneumonia Detection Challenge**"},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","collapsed":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":false},"cell_type":"markdown","source":"## **Overview**\n\n\n![alt text](https://i.pinimg.com/564x/ff/cd/c4/ffcdc4d74eed036d029a84c381604a10.jpg)\n\n## **Symptoms to detect Pneumonia**\n\nThe Diagnosis of pneumonia on CXR( is complicated because of a number of other nditions in the lungsuch as uid overload (pulmonary edema), bleeding, volume loss (atelectasis or collapse), lung cancer, or post-radiation or surgical changes."},{"metadata":{"_uuid":"c63d5668351eec6135d443fde229c44944461927"},"cell_type":"markdown","source":"## **1-Exploring the data**"},{"metadata":{"trusted":true,"_uuid":"a8702c6c864c8f4639e5b7897c1a9705318bc7bb","collapsed":true},"cell_type":"code","source":"import glob, pylab, pandas as pd\nimport pydicom, numpy as np\n\nimport os\nimport csv\nimport random\nfrom skimage import measure\nfrom skimage.transform import resize\n\nimport tensorflow as tf\nfrom tensorflow import keras\n\nfrom matplotlib import pyplot as plt\n\n!ls ../input","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a23bf9a1f1aae2bcec212483d5e315b12beff5a8","collapsed":true},"cell_type":"code","source":"df = pd.read_csv('../input/stage_1_train_labels.csv')\nprint(df.iloc[0])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"bee4507e172fa8955effc6314a9abf2baaa830b4","collapsed":true},"cell_type":"code","source":"print(df.iloc[8])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"81f51bb0577a92160d901d7fb5de092536a9a9f4"},"cell_type":"markdown","source":"## **Data Summary**\n\n**Stage 1 Images** - stage_1_train_images.zip and stage_1_test_images.zip\nimages for the current stage. Filenames are also patient names.\n\n**Stage 1 Labels** - stage_1_train_labels.csv and Stage 1 Sample Submission stage_1_sample_submission.csv\nWhich provides the IDs for the test set, as well as a sample of what your submission should look like\n\n**Stage 1 Detailed Info** - stage_1_detailed_class_info.csv\ncontains detailed information about the positive and negative classes in the training set, and may be used to build more nuanced models."},{"metadata":{"_uuid":"7fe06ecd7e5d97bb59131ccb0317d869a0f1fea2"},"cell_type":"markdown","source":"## **File descriptions**\n\n**stage_1_train.csv** - the training set. Contains patientIds and bounding box / target information.\n\n**stage_1_sample_submission.csv** - a sample submission file in the correct format.\n\n**stage_1_detailed_class_info.csv** - provides detailed information about the type of positive or negative class for each image."},{"metadata":{"_uuid":"25a969732444b85f3fc3c73c16f595d79bfeadf4"},"cell_type":"markdown","source":"## **Data fields**\n**patientId _**- A patientId. Each patientId corresponds to a unique image.\n\n**x_ **- the upper-left x coordinate of the bounding box.\n\n**y_ **- the upper-left y coordinate of the bounding box.\n\n**width_** - the width of the bounding box.\n\n**height_** - the height of the bounding box.\n\n**Target_** - the binary Target, indicating whether this sample has evidence of pneumonia."},{"metadata":{"trusted":true,"_uuid":"cd930fe22947651cae0329ae8ffb55c71551d7a0","collapsed":true},"cell_type":"code","source":"from os import listdir\nfrom os.path import isfile, join\n\n\ndet_class_path = '../input/stage_1_detailed_class_info.csv'\nbbox_path = '../input/stage_1_train_labels.csv'\ndicom_dir = '../input/stage_1_train_images/'\n\ntrain_images_dir = '../input/stage_1_train_images/'\ntrain_images = [f for f in listdir(train_images_dir) if isfile(join(train_images_dir, f))]\ntest_images_dir = '../input/stage_1_test_images/'\ntest_images = [f for f in listdir(test_images_dir) if isfile(join(test_images_dir, f))]\n\nprint('Number of train images:', len(train_images))\nprint('Number of test images:', len(test_images))","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"aa0437c54b260410bebabee2d7454ea3dfee2ab9"},"cell_type":"markdown","source":"## **DICOM files:**"},{"metadata":{"trusted":true,"scrolled":true,"_uuid":"d51193846e4c17ee960deabdb0d9040f55fcedc9","collapsed":true},"cell_type":"code","source":"patientId = df['patientId'][4]\ndcm_file = '../input/stage_1_train_images/%s.dcm' % patientId\ndcm_data = pydicom.read_file(dcm_file)\n\nprint(dcm_data)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"2a16d5665eb6e2ad843a2622be221dccaf6343d9","collapsed":true},"cell_type":"code","source":"patientId2 = df['patientId'][55]\ndcm_file2 = '../input/stage_1_train_images/%s.dcm' % patientId2\ndcm_data2 = pydicom.read_file(dcm_file2)\n\nprint(dcm_data2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e7c450e6975a4a6694c6e91ffbdcd982b0600556","collapsed":true},"cell_type":"code","source":"im = dcm_data.pixel_array\nprint(type(im))\nprint(im.dtype)\nprint(im.shape)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"df23716ba2c54d59faae613443617f4f710612ca"},"cell_type":"markdown","source":"* high bit-depth original images have been rescaled to 8-bit encoding (256 grayscales)\n* Thel image matrices (typically acquired at >2000 x 2000) have been resized to the shape of 1024 x 1024"},{"metadata":{"trusted":true,"_uuid":"53fdb2d65d9706e47438e9c7805b799e199d6cf1","collapsed":true},"cell_type":"code","source":"pylab.imshow(im, cmap=pylab.cm.gist_gray)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"153a108ce223f2d2d07cbeac2b3dd82c11a9ce9c","collapsed":true},"cell_type":"code","source":"im2 = dcm_data2.pixel_array\npylab.imshow(im2, cmap=pylab.cm.gist_gray)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"362609c777e6587da10a48907257ec00a1c53cf3"},"cell_type":"markdown","source":"## **CSV into a data structure with unique entries (the patient ID):**\n\nAny patient may  have many boxes if there are several different suspicious areas of pneumonia. "},{"metadata":{"trusted":true,"_uuid":"3a1c7afeaa3887a738ae950f6962cec77468a92d","collapsed":true},"cell_type":"code","source":"def parse_data(df):\n    \"\"\"\n      parsed = {\n        \n        'patientId-00': {\n            'dicom': path/to/dicom/file,\n            'label': either 0 or 1 for normal or pnuemonia, \n            'boxes': list of box(es)\n        },\n        'patientId-01': {\n            'dicom': path/to/dicom/file,\n            'label': either 0 or 1 for normal or pnuemonia, \n            'boxes': list of box(es)\n        }, ...\n\n      }\n\n    \"\"\"\n    # Define lambda to extract coords in list [y, x, height, width]\n    extract_box = lambda row: [row['y'], row['x'], row['height'], row['width']]\n\n    parsed = {}\n    for n, row in df.iterrows():\n        pid = row['patientId']\n        if pid not in parsed:\n            parsed[pid] = {\n                'dicom': '../input/stage_1_train_images/%s.dcm' % pid,\n                'label': row['Target'],\n                'boxes': []}\n\n        # Add box if opacity is present\n        if parsed[pid]['label'] == 1:\n            parsed[pid]['boxes'].append(extract_box(row))\n\n    return parsed\n\n\n\nparsed = parse_data(df)\n\n\n\nprint(parsed['00436515-870c-4b36-a041-de91049b9ab4'])\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"3ce3a306f843df6c9b51a65d26946e5954a0f2f4"},"cell_type":"markdown","source":"## **Overlay color boxes on the original grayscale DICOM files:****"},{"metadata":{"trusted":true,"_uuid":"9064982148b3bb222cd49e2ba755cc59d41e9c40","collapsed":true},"cell_type":"code","source":"def draw(data):\n    #Open DICOM file\n    d = pydicom.read_file(data['dicom'])\n    im = d.pixel_array\n\n    #Convert from single-channel grayscale to 3-channel RGB\n    im = np.stack([im] * 3, axis=2)\n\n    #Add boxes with random color if present\n    for box in data['boxes']:\n        rgb = np.floor(np.random.rand(3) * 256).astype('int')\n        im = overlay_box(im=im, box=box, rgb=rgb, stroke=6)\n\n    pylab.imshow(im, cmap=pylab.cm.gist_gray)\n    \n    \n\ndef overlay_box(im, box, rgb, stroke=1):\n    #Convert coordinates to integers\n    box = [int(b) for b in box]\n    \n    #Extract coordinates\n    y1, x1, height, width = box\n    y2 = y1 + height\n    x2 = x1 + width\n\n    im[y1:y1 + stroke, x1:x2] = rgb\n    im[y2:y2 + stroke, x1:x2] = rgb\n    im[y1:y2, x1:x1 + stroke] = rgb\n    im[y1:y2, x2:x2 + stroke] = rgb\n\n    return im\n\n\n\ndraw(parsed['00436515-870c-4b36-a041-de91049b9ab4'])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"bc01ee9bed986c97bb458e17dd5fc971cfe5f743"},"cell_type":"markdown","source":"## **Exploring Detailed Labels**\n\nIn addition to the binary classification (the presence or absence of pneumonia), each bounding box without pneumonia is further categorized into normal or no lung opacity / not normal (abnormality on the image)"},{"metadata":{"trusted":true,"_uuid":"a333f128f9ac221c796bfdcbd311b3268020df20","collapsed":true},"cell_type":"code","source":"df_detailed = pd.read_csv('../input/stage_1_detailed_class_info.csv')\nprint(df_detailed.iloc[6])\nprint(df_detailed.iloc[80])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"85c55ee17a6d87e1fba5943e07e544da20b746da"},"cell_type":"markdown","source":"**Model**"},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"72b218c91de23b6fb51aa3acf2030f92cab201c0"},"cell_type":"code","source":"# empty dictionary\nnodule_locations = {}\n# load table\nwith open(os.path.join('../input/stage_1_train_labels.csv'), mode='r') as infile:\n    reader = csv.reader(infile)\n    # skip header\n    next(reader, None)\n\n    for rows in reader:\n        filename = rows[0]\n        location = rows[1:5]\n        nodule = rows[5]\n        # if row contains a nodule add label to dictionary\n        # which contains a list of nodule locations per filename\n        if nodule == '1':\n            # convert string to float to int\n            location = [int(float(i)) for i in location]\n            # save nodule location in dictionary\n            if filename in nodule_locations:\n                nodule_locations[filename].append(location)\n            else:\n                nodule_locations[filename] = [location]\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"aa5006c4161e86854aa3adaa11d500357a263fba","collapsed":true},"cell_type":"code","source":"folder = '../input/stage_1_train_images'\nfilenames = os.listdir(folder)\nrandom.shuffle(filenames)\n# split into train and validation filenames\nn_valid_samples = 2000\ntrain_filenames = filenames[n_valid_samples:]\nvalid_filenames = filenames[:n_valid_samples]\nprint('n train samples', len(train_filenames))\nprint('n valid samples', len(valid_filenames))\nn_train_samples = len(filenames) - n_valid_samples","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"f13efc45f0e0416e71283d0e129cc3957129def4"},"cell_type":"code","source":"class generator(keras.utils.Sequence):\n    \n    def __init__(self, folder, filenames, nodule_locations=None, batch_size=32, image_size=256, shuffle=True, predict=False):\n        self.folder = folder\n        self.filenames = filenames\n        self.nodule_locations = nodule_locations\n        self.batch_size = batch_size\n        self.image_size = image_size\n        self.shuffle = shuffle\n        self.predict = predict\n        self.on_epoch_end()\n        \n    def __load__(self, filename):\n        # load dicom file as numpy array\n        img = pydicom.dcmread(os.path.join(self.folder, filename)).pixel_array\n        # create empty mask\n        msk = np.zeros(img.shape)\n        # get filename without extension\n        filename = filename.split('.')[0]\n        # if image contains nodules\n        if filename in nodule_locations:\n            # loop through nodules\n            for location in nodule_locations[filename]:\n                # add 1's at the location of the nodule\n                x, y, w, h = location\n                msk[y:y+h, x:x+w] = 1\n        # resize both image and mask\n        img = resize(img, (self.image_size, self.image_size), mode='reflect')\n        msk = resize(msk, (self.image_size, self.image_size), mode='reflect') > 0.5\n        # add trailing channel dimension\n        img = np.expand_dims(img, -1)\n        msk = np.expand_dims(msk, -1)\n        return img, msk\n    \n    def __loadpredict__(self, filename):\n        # load dicom file as numpy array\n        img = pydicom.dcmread(os.path.join(self.folder, filename)).pixel_array\n        # resize image\n        img = resize(img, (self.image_size, self.image_size), mode='reflect')\n        # add trailing channel dimension\n        img = np.expand_dims(img, -1)\n        return img\n        \n    def __getitem__(self, index):\n        # select batch\n        filenames = self.filenames[index*self.batch_size:(index+1)*self.batch_size]\n        # predict mode: return images and filenames\n        if self.predict:\n            # load files\n            imgs = [self.__loadpredict__(filename) for filename in filenames]\n            # create numpy batch\n            imgs = np.array(imgs)\n            return imgs, filenames\n        # train mode: return images and masks\n        else:\n            # load files\n            items = [self.__load__(filename) for filename in filenames]\n            # unzip images and masks\n            imgs, msks = zip(*items)\n            # create numpy batch\n            imgs = np.array(imgs)\n            msks = np.array(msks)\n            return imgs, msks\n        \n    def on_epoch_end(self):\n        if self.shuffle:\n            random.shuffle(self.filenames)\n        \n    def __len__(self):\n        if self.predict:\n            # return everything\n            return int(np.ceil(len(self.filenames) / self.batch_size))\n        else:\n            # return full batches only\n            return int(len(self.filenames) / self.batch_size)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"de6c835be6485e25403f752959dfbeeb254ef18f"},"cell_type":"code","source":"def create_downsample(channels, inputs):\n    x = keras.layers.BatchNormalization(momentum=0.9)(inputs)\n    x = keras.layers.LeakyReLU(0)(x)\n    x = keras.layers.Conv2D(channels, 1, padding='same', use_bias=False)(x)\n    x = keras.layers.MaxPool2D(2)(x)\n    return x\n\ndef create_resblock(channels, inputs):\n    x = keras.layers.BatchNormalization(momentum=0.9)(inputs)\n    x = keras.layers.LeakyReLU(0)(x)\n    x = keras.layers.Conv2D(channels, 3, padding='same', use_bias=False)(x)\n    x = keras.layers.BatchNormalization(momentum=0.9)(x)\n    x = keras.layers.LeakyReLU(0)(x)\n    x = keras.layers.Conv2D(channels, 3, padding='same', use_bias=False)(x)\n    return keras.layers.add([x, inputs])\n\ndef create_network(input_size, channels, n_blocks=2, depth=4):\n    # input\n    inputs = keras.Input(shape=(input_size, input_size, 1))\n    x = keras.layers.Conv2D(channels, 3, padding='same', use_bias=False)(inputs)\n    # residual blocks\n    for d in range(depth):\n        channels = channels * 2\n        x = create_downsample(channels, x)\n        for b in range(n_blocks):\n            x = create_resblock(channels, x)\n    # output\n    x = keras.layers.BatchNormalization(momentum=0.9)(x)\n    x = keras.layers.LeakyReLU(0)(x)\n    x = keras.layers.Conv2D(1, 1, activation='sigmoid')(x)\n    outputs = keras.layers.UpSampling2D(2**depth)(x)\n    model = keras.Model(inputs=inputs, outputs=outputs)\n    return model","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"18ff18e1b02deafebcf23d0a4a8a0c74e8b2e2f3","collapsed":true},"cell_type":"code","source":"# define iou or jaccard loss function\ndef iou_loss(y_true, y_pred):\n    y_true = tf.reshape(y_true, [-1])\n    y_pred = tf.reshape(y_pred, [-1])\n    intersection = tf.reduce_sum(y_true * y_pred)\n    score = (intersection + 1.) / (tf.reduce_sum(y_true) + tf.reduce_sum(y_pred) - intersection + 1.)\n    return 1 - score\n\n# combine bce loss and iou loss\ndef iou_bce_loss(y_true, y_pred):\n    return 0.5 * keras.losses.binary_crossentropy(y_true, y_pred) + 0.5 * iou_loss(y_true, y_pred)\n\n# mean iou as a metric\ndef mean_iou(y_true, y_pred):\n    y_pred = tf.round(y_pred)\n    intersect = tf.reduce_sum(y_true * y_pred, axis=[1, 2, 3])\n    union = tf.reduce_sum(y_true, axis=[1, 2, 3]) + tf.reduce_sum(y_pred, axis=[1, 2, 3])\n    smooth = tf.ones(tf.shape(intersect))\n    return tf.reduce_mean((intersect + smooth) / (union - intersect + smooth))\n\n# create network and compiler\nmodel = create_network(input_size=256, channels=32, n_blocks=2, depth=4)\nmodel.compile(optimizer='adam',\n              loss=iou_bce_loss,\n              metrics=['accuracy', mean_iou])\n\n# cosine learning rate annealing\ndef cosine_annealing(x):\n    lr = 0.001\n    epochs = 25\n    return lr*(np.cos(np.pi*x/epochs)+1.)/2\nlearning_rate = tf.keras.callbacks.LearningRateScheduler(cosine_annealing)\n\n# create train and validation generators\nfolder = '../input/stage_1_train_images'\ntrain_gen = generator(folder, train_filenames, nodule_locations, batch_size=32, image_size=256, shuffle=True, predict=False)\nvalid_gen = generator(folder, valid_filenames, nodule_locations, batch_size=32, image_size=256, shuffle=False, predict=False)\n\nhistory = model.fit_generator(train_gen, validation_data=valid_gen, callbacks=[learning_rate], epochs=2, shuffle=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f2ee52cc1785521c8d206b2f8b32a8481d8de37b","collapsed":true},"cell_type":"code","source":"from sklearn import feature_selection, linear_model, metrics, preprocessing\n\n\nfolder = '../input/stage_1_test_images/'\nfilenames = os.listdir(folder)\n#random.shuffle(filenames)\n\n\n\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"56d58de85b91f74f28bf6bc97ab2658d3d05cab2","collapsed":true},"cell_type":"code","source":"n_y_samples = 1 #must be at least 500 (it is currently used to test 1 image only)\ntest_X = filenames[999:]  #must be: test_X = filenames[n_y_samples:]\ntest_Y = filenames[:n_y_samples]\nprint('n test_x', len(test_X))\nprint('n test_y', len(test_Y))\nn_x_samples = len(filenames) - n_y_samples","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"20076a22bdbf47fe541f013c0317bbc1029d8059","collapsed":true},"cell_type":"code","source":"#should return a shape of (1, 32, 32, 3)  ???\n\ntest_x_gen = generator(folder, test_X, nodule_locations, batch_size=32, image_size=256, shuffle=True, predict=False)\ntest_y_gen = generator(folder, test_Y, nodule_locations, batch_size=32, image_size=256, shuffle=True, predict=False)\n\n#print (test_x_gen)\n#print (test_y_gen)\n    \n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"aaa85a8843852d943821499679fd420347a4a604","collapsed":true},"cell_type":"code","source":"#score = model.evaluate(np.expand_dims(test_x_gen, axis=3), test_y_gen, batch_size=32)\n#print (score)\n\n#model.predict(test_x_gen)\n              \n#metrics.accuracy_score(np.array(test_y_gen), model.predict(np.array(test_x_gen)))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"c06abd384d144d43c991c2494c3c86821b9baf4e"},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.4","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}