{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import os\nimport csv\nimport random\nimport pydicom\nimport numpy as np\nimport pandas as pd\nfrom skimage import io\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\nimport matplotlib.patches as patches","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# model configurations\nclass Configuration:\n    EPOCHS = 5\n    LEARNING_RATE = 0.001\n    BATCH = 32\n    IMAGE_SIZE = 64\n    CHANNLES = 32\n    N_BLOCKS = 2\n    NN_DEPTH = 4\n    WORKERS = 4\n    TRAIN_LABELS_PATH = '../input/rsna-pneumonia-detection-challenge/stage_2_train_labels.csv'\n    TRAIN_IMAGES_PATH = '../input/rsna-pneumonia-detection-challenge/stage_2_train_images'\n    TEST_IMAGES_PATH = '../input/rsna-pneumonia-detection-challenge/stage_2_test_images'\n    TRAIN_IMAGES_SIZE = 500 \n    VALID_IMAGES_SIZE = 200","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# create empty disctionary to store the pneumonia locations for all the unique files\npneumonia_locations = {}\n# Open and load stage_2_train_labels.csv file to read all the records and add locations to dictionary for all uniqueue files\n# Fields \n    # patientid, x, y, width, height and target\nwith open(os.path.join(Configuration.TRAIN_LABELS_PATH), mode='r') as infile:\n    # open reader\n    reader = csv.reader(infile)\n    # skip header\n    next(reader, None)\n    # loop through rows\n    for rows in reader:\n        filename = rows[0] # get file name\n        location = rows[1:5] # location as (x, y) point and width, height \n        target = rows[5] # it shows if the patient have pheumonia or not\n        # if row contains pneumonia add location to dictionary\n        # which contains a list of pneumonia locations per filename\n        if target == '1':\n            # convert string to float to int\n            location = [int(float(i)) for i in location]\n            # save pneumonia location in dictionary\n            if filename in pneumonia_locations:\n                pneumonia_locations[filename].append(location)\n            else:\n                pneumonia_locations[filename] = [location]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"file_name = '00436515-870c-4b36-a041-de91049b9ab4'\nimg = pydicom.dcmread(os.path.join(Configuration.TRAIN_IMAGES_PATH, file_name +'.dcm')).pixel_array\nmsk = np.zeros(img.shape)\nif filename in pneumonia_locations:\n    # loop through pneumonia\n    for location in pneumonia_locations[filename]:\n        # add 1's at the location of the pneumonia\n        x, y, w, h = location\n        msk[y:y+h, x:x+w] = 1\nimg = resize(img, (Configuration.IMAGE_SIZE, Configuration.IMAGE_SIZE), mode='reflect')\nmsk = resize(msk, (Configuration.IMAGE_SIZE, Configuration.IMAGE_SIZE), mode='reflect') > 0.5\nimg = np.expand_dims(img, -1) \nmsk = np.expand_dims(msk, -1) \n\naxidx = 0\nf, axarr = plt.subplots(1, 2, figsize=(20,15))\naxarr = axarr.ravel()\naxarr[axidx].imshow(img[:, :, 0])\n# threshold true mask\ncomp = msk[:, :, 0] > 0.5\n# apply connected components\ncomp = measure.label(comp)\n# apply bounding boxes\npredictionString = ''\nfor region in measure.regionprops(comp):\n    # retrieve x, y, height and width\n    y, x, y2, x2 = region.bbox\n    height = y2 - y\n    width = x2 - x\n    axarr[axidx].add_patch(patches.Rectangle((x,y),width,height,linewidth=2,edgecolor='b',facecolor='none'))\naxidx += 1","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# load and shuffle train images\nfolder = Configuration.TRAIN_IMAGES_PATH\nfilenames = os.listdir(folder)\nrandom.shuffle(filenames)\n# split into train and validation image filenames\ntrain_filenames = filenames[:Configuration.TRAIN_IMAGES_SIZE]\nvalid_filenames = filenames[Configuration.TRAIN_IMAGES_SIZE:Configuration.VALID_IMAGES_SIZE + Configuration.TRAIN_IMAGES_SIZE]\nprint('number of train samples', len(train_filenames))\nprint('number of valid samples', len(valid_filenames))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Total train images:',len(filenames))\nprint('Images with pneumonia:', len(pneumonia_locations))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class generator(keras.utils.Sequence):\n    \n    # Construcor to initialize the object's properties\n    def __init__(self, folder, filenames, pneumonia_locations=None, batch_size=32, image_size=256, shuffle=True, augment=False, predict=False):\n        self.folder = folder\n        self.filenames = filenames\n        self.pneumonia_locations = pneumonia_locations\n        self.batch_size = batch_size\n        self.image_size = image_size\n        self.shuffle = shuffle\n        self.augment = augment\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 pneumonia\n        if filename in self.pneumonia_locations:\n            # loop through pneumonia\n            for location in self.pneumonia_locations[filename]:\n                # add 1's at the location of the pneumonia\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        # if augment then horizontal flip half the time\n        if self.augment and random.random() > 0.5: \n            img = np.fliplr(img)  \n            msk = np.fliplr(msk) \n        # add trailing channel dimension\n        img = np.expand_dims(img, -1) \n        msk = np.expand_dims(msk, -1) \n        # print('Load method img size {0} and msk size is {1}'.format(img.shape, msk.shape))\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        #print('LoadPredict method img size {0}'.format(img.shape))\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            #print('GetItem method and predict is true img size {0} and total filenames are {1}'.format(img.shape, len(filenames)))\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            #print('GetItem method and predict is false img size {0} and mask size {1}'.format(imgs.shape, msks.shape))\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            #print('Len method and predict is true and values {0}'.format(int(np.ceil(len(self.filenames) / self.batch_size))))\n            return int(np.ceil(len(self.filenames) / self.batch_size))\n        else:\n            # return full batches only\n            #print('Len method and predict is false and values {0}'.format(int(len(self.filenames) / self.batch_size)))\n            return int(len(self.filenames) / self.batch_size)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"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},"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#     if type(y_true) == 'bool':\n#         print(' type is bool')\n#     print(y_true.dtype == 'bool')\n#     print(y_true)\n#     print(y_pred)\n    if y_true.dtype == 'bool':\n        y_true = tf.cast(y_true, tf.float32)\n#     if y_pred.dtype == 'bool':\n#         y_pred = tf.cast(y_pred, tf.float32)\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 = Configuration.IMAGE_SIZE, channels = Configuration.CHANNLES, n_blocks = Configuration.N_BLOCKS, depth = Configuration.NN_DEPTH)\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 = Configuration.LEARNING_RATE\n    epochs = Configuration.EPOCHS\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 = Configuration.TRAIN_IMAGES_PATH \ntrain_gen = generator(folder, train_filenames, pneumonia_locations, batch_size=Configuration.BATCH, image_size=Configuration.IMAGE_SIZE, shuffle=True, augment=True, predict=False)\nvalid_gen = generator(folder, valid_filenames, pneumonia_locations, batch_size=Configuration.BATCH, image_size=Configuration.IMAGE_SIZE, shuffle=False, predict=False)\n\nhistory = model.fit_generator(train_gen, validation_data=valid_gen, callbacks=[learning_rate], epochs=Configuration.EPOCHS, workers=Configuration.WORKERS, use_multiprocessing=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=(12,4))\nplt.subplot(131)\nplt.plot(history.epoch, history.history[\"loss\"], label=\"Train loss\")\nplt.plot(history.epoch, history.history[\"val_loss\"], label=\"Valid loss\")\nplt.legend()\nplt.subplot(133)\nplt.plot(history.epoch, history.history[\"mean_iou\"], label=\"Train iou\")\nplt.plot(history.epoch, history.history[\"val_mean_iou\"], label=\"Valid iou\")\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# load and shuffle filenames\ntest_filenames = os.listdir(Configuration.TEST_IMAGES_PATH)\nprint('number of test samples:', len(test_filenames))\n\n# create test generator with predict flag set to True\ntest_gen = generator(Configuration.TEST_IMAGES_PATH, test_filenames, None, batch_size=Configuration.BATCH, image_size=256, shuffle=False, predict=True)\n\n# create submission dictionary\nsubmission_dict = {}\n# loop through testset\nfor imgs, filenames in test_gen:\n    # predict batch of images\n    preds = model.predict(imgs)\n    # create figure\n    f, axarr = plt.subplots(4, 8, figsize=(50,30))\n    axarr = axarr.ravel()\n    axidx = 0\n    # loop through batch\n    for img, pred, filename in zip(imgs, preds, filenames):\n         # plot image\n        axarr[axidx].imshow(img[:, :, 0])\n        # resize predicted mask\n        pred = resize(pred, (1024, 1024), mode='reflect')\n        # threshold predicted mask\n        comp = pred[:, :, 0] > 0.5\n        # apply connected components\n        comp = measure.label(comp)\n        # apply bounding boxes\n        predictionString = ''\n        for region in measure.regionprops(comp):\n            # retrieve x, y, height and width\n            y, x, y2, x2 = region.bbox\n            height = y2 - y\n            width = x2 - x\n            axarr[axidx].add_patch(patches.Rectangle((x,y),width,height,linewidth=2,edgecolor='b',facecolor='none'))\n            # proxy for confidence score\n            conf = np.mean(pred[y:y+height, x:x+width])\n            # add to predictionString\n            predictionString += str(conf) + ' ' + str(x) + ' ' + str(y) + ' ' + str(width) + ' ' + str(height) + ' '\n        # add filename and predictionString to dictionary\n        axidx += 1\n        filename = filename.split('.')[0]\n        submission_dict[filename] = predictionString\n    plt.show()\n    # only process one batch\n    break\nprint(submission_dict)","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}