{"cells":[{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport pydicom as dcm\nimport cv2\nimport os\nimport tensorflow as tf\nfrom tensorflow.keras.applications.vgg16 import preprocess_input\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Mirror Statergy Setup"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_available_gpus():\n    local_device_protos = tf.python.client.device_lib.list_local_devices()\n    return [x.name for x in local_device_protos]\n\nget_available_gpus()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"strategy = tf.distribute.MirroredStrategy()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from vtk import vtkDICOMImageReader\nfrom vtk import vtkImageShiftScale\nimport vtk\nfrom vtk.util import numpy_support\nfrom vtk import vtkPNGWriter\n\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\n\ndef window_image(img, window_center,window_width, intercept, slope, rescale=True):\n    img = (img*slope +intercept) #for translation adjustments given in the dicom file. \n    img_min = window_center - window_width//2 #minimum HU level\n    img_max = window_center + window_width//2 #maximum HU level\n    img[img<img_min] = img_min #set img_min for all HU levels less than minimum HU level\n    img[img>img_max] = img_max #set img_max for all HU levels higher than maximum HU level\n    if rescale: \n        img = (img - img_min) / (img_max - img_min)*255.0 \n    return img\n    \ndef get_first_of_dicom_field_as_int(x):\n    #get x[0] as in int is x is a 'pydicom.multival.MultiValue', otherwise get int(x)\n    if type(x) == dcm.multival.MultiValue: return int(x[0])\n    else: return int(x)\n    \ndef get_windowing(data):\n    dicom_fields = [data[('0028','1050')].value, #window center\n                    data[('0028','1051')].value, #window width\n                    data[('0028','1052')].value, #intercept\n                    data[('0028','1053')].value] #slope\n    return [get_first_of_dicom_field_as_int(x) for x in dicom_fields]\n\n\ndef decode_dcm(file, dim = 1,  rescale = True):\n    data = dcm.dcmread(file)\n    try:\n        image = data.pixel_array\n    except:\n        image = get_img(file)\n    window_center , window_width, intercept, slope = get_windowing(data)\n    output = window_image(image, window_center, window_width, intercept, slope, rescale = rescale)\n    if dim != 1:\n        output = np.repeat(output[..., np.newaxis], dim, -1)\n    return output","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def tags(data):\n    tags = tio.image.dicom_tags()\n    window_center = tio.image.decode_dicom_data(data, tags = tags.WindowCenter)\n    WindowWidth = tio.image.decode_dicom_data(data, tags = tags.WindowWidth)\n    RescaleIntercept = tio.image.decode_dicom_data(data, tags= tags.RescaleIntercept)\n    RescaleSlope = tio.image.decode_dicom_data(data, tags = tags.RescaleSlope)\n    return tf.strings.to_number(window_center), tf.strings.to_number(WindowWidth), tf.strings.to_number(RescaleIntercept), tf.strings.to_number(RescaleSlope)\n\ndef dcm_to_image(filename):\n    return tio.image.decode_dicom_image(filename, on_error = 'lossy', dtype = tf.float32, )\n\ndef window_image_tf(img, window_center,window_width, intercept, slope, rescale=True):\n    img = (img*slope +intercept)\n    img_min = window_center - window_width//2\n    img_max = window_center + window_width//2\n    img = tf.clip_by_value(img, img_min, img_max)\n    if rescale: \n        img = (img - img_min) / (img_max - img_min)*255.0 \n    return img\n\ndef string_to_image(filename):\n    encodedata = tf.io.read_file(filename)\n    image = dcm_to_image(encodedata)\n    window_center,window_width, intercept, slope = tags(encodedata)\n    return window_image_tf(image, window_center,window_width, intercept, slope)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def preprocess_image(x, y):\n    x = x.numpy()\n    if isinstance(x, np.ndarray):\n        x = list(map(lambda a : decode_dcm(a.decode('utf-8'), dim = 3), x))\n    else:\n        x = decode_dcm(x.decode('utf-8'), dim = 3)\n    return (np.array(x), y)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def process_string(x, y):\n    return tf.strings.unicode_decode(x,input_encoding='UTF-8'), y","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df = pd.read_csv(\"../input/rsna-str-pulmonary-embolism-detection/train.csv\")\nx_columns = train_df.columns[:3]\noutput = ['pe_present_on_image','negative_exam_for_pe','indeterminate', 'rv_lv_ratio_gte_1', 'rv_lv_ratio_lt_1',\n            'chronic_pe','acute_and_chronic_pe', 'indeterminate',\n            'leftsided_pe', 'rightsided_pe', 'central_pe']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"root_path = \"../input/rsna-str-pulmonary-embolism-detection/train/\"\ntrain_df_x = list(map(lambda x : root_path + \"/\".join(x) +\".dcm\", train_df[x_columns].values.tolist()))\ntrain_df_y = train_df[output].values.astype('float32')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df_x = tf.data.Dataset.from_tensor_slices(train_df_x)\n# train_df_x = train_df_x.map(lambda x: tf.py_function(preprocess_image, [x], [tf.float32]))\ntrain_df_y = tf.data.Dataset.from_tensor_slices(train_df_y)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"ds = tf.data.Dataset.zip((train_df_x, train_df_y))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"BUFFER_SIZE=1300\nPER_REPLICA_BATCH_SIZE=64\ntry:\n    NUM_REPLICAS=strategy.num_replicas_in_sync\nexcept:\n    NUM_REPLICAS = 1\nGLOBAL_BATCH_SIZE = PER_REPLICA_BATCH_SIZE * NUM_REPLICAS\ntrain_set = int((np.round(train_df.shape[0] * 0.70)) // GLOBAL_BATCH_SIZE)\n\nds = ds.shuffle(BUFFER_SIZE)\nds = ds.batch(GLOBAL_BATCH_SIZE)\ntrain_ds = ds.take(train_set)\ntest_ds = ds.skip(train_set)\n#train_ds = train_ds.interleave(train_ds, num_parallel_calls = tf.data.experimental.AUTOTUNE)\ntrain_ds = train_ds.map(lambda x, y : tf.py_function(preprocess_image, [x, y], [tf.float32, tf.float32]))\ntest_ds = test_ds.map(lambda x, y : tf.py_function(preprocess_image, [x, y], [tf.float32, tf.float32]))\ntest_ds = test_ds.prefetch(tf.data.experimental.AUTOTUNE)\ntrain_ds = train_ds.prefetch(tf.data.experimental.AUTOTUNE)\n\ntrain_example_ds = train_ds.take(100)\ntest_example_ds = test_ds.take(100)\ntrain_dist_dataset = strategy.experimental_distribute_dataset(train_example_ds)\ntest_dist_dataset = strategy.experimental_distribute_dataset(test_example_ds)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"it = iter(train_example_ds)\nnext(it)\nnext(it)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"del train_df_x, train_df_y\n\nimport gc \ngc.collect()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"x, y = next(iter(train_dist_dataset))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"print(x.numpy().shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Create Model"},{"metadata":{"trusted":true},"cell_type":"code","source":"feature_extraction = tf.keras.applications.ResNet50(include_top = False, input_shape= (512, 512, 3))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for layer in feature_extraction.layers:\n    layer.trainable = False","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"kernel_init = tf.keras.initializers.glorot_uniform()\nbias_init = tf.keras.initializers.Constant(value=0.2)\ndef inception_module(x,\n                     filters_1x1,\n                     filters_3x3_reduce,\n                     filters_3x3,\n                     filters_5x5_reduce,\n                     filters_5x5,\n                     filters_pool_proj,\n                     name=None):\n    \n    if filters_1x1 == 0:\n        conv_1x1 = tf.keras.layers.Conv2D(filters_1x1, (1, 1), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(x)\n    \n    conv_3x3 = tf.keras.layers.Conv2D(filters_3x3_reduce, (1, 1), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(x)\n    conv_3x3 = tf.keras.layers.Conv2D(filters_3x3, (3, 3), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(conv_3x3)\n\n    conv_5x5 = tf.keras.layers.Conv2D(filters_5x5_reduce, (1, 1), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(x)\n    conv_5x5 = tf.keras.layers.Conv2D(filters_5x5, (5, 5), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(conv_5x5)\n\n    pool_proj = tf.keras.layers.MaxPool2D((3, 3), strides=(1, 1), padding='same')(x)\n    pool_proj = tf.keras.layers.Conv2D(filters_pool_proj, (1, 1), padding='same', activation='relu', kernel_initializer=kernel_init, bias_initializer=bias_init)(pool_proj)\n\n    if filters_1x1 == 0:\n        output = tf.keras.layers.concatenate([conv_1x1, conv_3x3, conv_5x5, pool_proj], axis=3, name=name)\n    else:\n        output = tf.keras.layers.concatenate([conv_3x3, conv_5x5, pool_proj], axis=3, name=name)\n    return output","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def create_model():\n    inputs = tf.keras.Input(shape=(512, 512, 3))\n    features = feature_extraction(inputs)\n    inception_downsample = tf.keras.layers.MaxPool2D((2, 2), strides=(2, 2), name='max_pool_1')(features)\n    inception_downsample = inception_module(inception_downsample,filters_1x1=32,filters_3x3_reduce=32,filters_3x3=32,filters_5x5_reduce=32,\n                     filters_5x5=32,filters_pool_proj=10,name='inception_1')\n    inception_downsample = tf.keras.layers.BatchNormalization(name = 'BatchNormalization_1')(inception_downsample)\n    inception_downsample = tf.keras.layers.MaxPool2D((2, 2), name='max_pool_2')(inception_downsample)\n    # Classification Layer\n    x = tf.keras.layers.Flatten()(inception_downsample)\n    x = tf.keras.layers.Dense(230)(x)\n    x = tf.keras.layers.LeakyReLU()(x)\n    x = tf.keras.layers.Dense(128)(x)\n    x = tf.keras.layers.LeakyReLU()(x)\n    y1 = tf.keras.layers.Dense(1,'sigmoid',name = 'pe_present_on_image')(x)\n    y2 = tf.keras.layers.Dense(1,'sigmoid',name = 'negative_exam_for_pe')(x)\n    y3 = tf.keras.layers.Dense(1,'sigmoid',name = 'indeterminate')(x)\n\n    model = tf.keras.Model(inputs = inputs, outputs = [y1, y2,y3])\n    print(model.summary())\n    return model","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Loss Function"},{"metadata":{"trusted":true},"cell_type":"code","source":"with strategy.scope():\n    loss_object = tf.keras.losses.BinaryCrossentropy(\n      reduction=tf.keras.losses.Reduction.NONE)\n\n    def compute_loss(labels, predictions):\n        labels1 = labels[:, 0:1]\n        labels2 = labels[:, 1:2]\n        labels3 = labels[:, 2:3]\n        predictions1 = predictions[0]\n        predictions2 = predictions[1]\n        predictions3 = predictions[2]\n        per_example_loss1 = loss_object(labels1, predictions1)\n        per_example_loss2 = loss_object(labels2, predictions2)\n        per_example_loss3 = loss_object(labels3, predictions3)\n        return tf.nn.compute_average_loss(per_example_loss1, global_batch_size=GLOBAL_BATCH_SIZE), tf.nn.compute_average_loss(per_example_loss2, global_batch_size=GLOBAL_BATCH_SIZE), tf.nn.compute_average_loss(per_example_loss3, global_batch_size=GLOBAL_BATCH_SIZE)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Metrics"},{"metadata":{"trusted":true},"cell_type":"code","source":"with strategy.scope():\n    test_loss1 = tf.keras.metrics.Mean(name='test_loss1')\n    test_loss2 = tf.keras.metrics.Mean(name='test_loss2')\n    test_loss3 = tf.keras.metrics.Mean(name='test_loss3')\n    train_accuracy1 = tf.keras.metrics.BinaryAccuracy(name='pe_present_on_image_accuracy')\n    train_accuracy2 = tf.keras.metrics.BinaryAccuracy(name='negative_exam_for_pe')\n    train_accuracy3 = tf.keras.metrics.BinaryAccuracy(name='indeterminate')\n    test_accuracy1 = tf.keras.metrics.BinaryAccuracy(name='test_accuracy1')\n    test_accuracy2 = tf.keras.metrics.BinaryAccuracy(name='test_accuracy2')\n    test_accuracy3 = tf.keras.metrics.BinaryAccuracy(name='test_accuracy3')\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Training Loop"},{"metadata":{"trusted":true},"cell_type":"code","source":"checkpoint_dir = './training_checkpoints'\ncheckpoint_prefix = os.path.join(checkpoint_dir, \"ckpt\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# model, optimizer, and checkpoint must be created under `strategy.scope`.\nwith strategy.scope():\n    model = create_model()\n    optimizer = tf.keras.optimizers.Adam()\n    checkpoint = tf.train.Checkpoint(optimizer=optimizer, model=model)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def train_step(inputs):\n    images, labels = inputs\n    with tf.GradientTape(persistent=True) as tape:\n        predictions = model(images, training=True)\n        loss1, loss2, loss3 = compute_loss(labels, predictions)\n    gradients1 = tape.gradient(loss1, model.trainable_variables)\n    gradients2 = tape.gradient(loss2, model.trainable_variables)\n    gradients3 = tape.gradient(loss3, model.trainable_variables)\n    optimizer.apply_gradients(zip(gradients1, model.trainable_variables))\n    optimizer.apply_gradients(zip(gradients2, model.trainable_variables))\n    optimizer.apply_gradients(zip(gradients3, model.trainable_variables))\n    labels1 = labels[:, 0:1]\n    labels2 = labels[:, 1:2]\n    labels3 = labels[:, 2:3]\n    predictions1 = predictions[0]\n    predictions2 = predictions[1]\n    predictions3 = predictions[2]\n    train_accuracy1.update_state(labels1, predictions1)\n    train_accuracy2.update_state(labels2, predictions2)\n    train_accuracy3.update_state(labels3, predictions3)\n    return loss1, loss2, loss3\n\ndef test_step(inputs):\n    images, labels = inputs\n    predictions = model(images, training=False)\n    t_loss1, t_loss2, t_loss3 = compute_loss(labels, predictions)\n    labels1 = labels[:, 0:1]\n    labels2 = labels[:, 1:2]\n    labels3 = labels[:, 2:3]\n    predictions1 = predictions[0]\n    predictions2 = predictions[1]\n    predictions3 = predictions[2]\n    test_loss1.update_state(t_loss1)\n    test_loss2.update_state(t_loss2)\n    test_loss3.update_state(t_loss3)\n    test_accuracy1.update_state(labels1, predictions1)\n    test_accuracy2.update_state(labels2, predictions2)\n    test_accuracy3.update_state(labels3, predictions3)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"@tf.function\ndef distributed_train_step(dataset_inputs):\n    per_replica_losses1, per_replica_losses2, per_replica_losses3 = strategy.run(train_step, args=(dataset_inputs,))\n    return strategy.reduce(tf.distribute.ReduceOp.SUM, per_replica_losses1, axis=None),strategy.reduce(tf.distribute.ReduceOp.SUM, per_replica_losses2,\n                         axis=None), strategy.reduce(tf.distribute.ReduceOp.SUM, per_replica_losses3,\n                         axis=None)\n\n@tf.function\ndef distributed_test_step(dataset_inputs):\n    return strategy.run(test_step, args=(dataset_inputs,))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"EPOCHS = 5\nfor epoch in range(EPOCHS):\n    # TRAIN LOOP\n    total_loss1 = 0.0\n    total_loss2 = 0.0\n    total_loss3 = 0.0\n    num_batches = 0\n    for x in train_dist_dataset:\n        loss1, loss2, loss3 = distributed_train_step(x)\n        total_loss1 += loss1\n        total_loss1 += loss2\n        total_loss1 += loss3\n        num_batches += 1\n        train_loss1 = total_loss1 / num_batches\n        train_loss2 = total_loss2 / num_batches\n        train_loss3 = total_loss3 / num_batches\n\n    # TEST LOOP\n    for x in test_dist_dataset:\n        distributed_test_step(x)\n\n    if epoch % 10 == 0:\n        checkpoint.save(checkpoint_prefix)\n\n    template = (\"Epoch {},\\n\\t Loss: \\n\\t\\t pe_present_on_image {} \\n\\t\\t negative_exam_for_pe {} \\n\\t\\t indeterminate {}, \\n\\tAccuracy:  \\n\\t\\t pe_present_on_image {} \\n\\t\\t negative_exam_for_p {} \\n\\t\\t indeterminate {}, \\n\\t Test Loss: \\n\\t\\t pe_present_on_image {} \\n\\t\\t negative_exam_for_pe {} \\n\\t\\t indeterminate {}, \\n\\t Test Accuracy: \\n\\t\\t pe_present_on_image {} \\n\\t\\t negative_exam_for_pe {} \\n\\t\\t indeterminate {}\")\n    print (template.format(epoch+1, train_loss1, train_loss2, train_loss3,\n                         train_accuracy1.result()*100, train_accuracy2.result()*100, train_accuracy3.result()*100\n                        ,test_loss1.result(),test_loss2.result(),test_loss3.result(),\n                         test_accuracy1.result()*100, test_accuracy2.result()*100, test_accuracy3.result()*100))\n\n    test_loss1.reset_states()\n    train_accuracy1.reset_states()\n    test_accuracy1.reset_states()\n    test_loss2.reset_states()\n    train_accuracy2.reset_states()\n    test_accuracy2.reset_states()\n    test_loss3.reset_states()\n    train_accuracy3.reset_states()\n    test_accuracy3.reset_states()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"x, y = next(iter(test_example_ds))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"x = x[0]\ny = y[0][:3]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"tf.round(model(x[np.newaxis, ...]))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model.save(\"model.h5\")","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}