{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Import Packages"},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install ../input/kerasapplications/keras-team-keras-applications-3b180cb -f ./ --no-index\n!pip install ../input/osic-manual-bbox/efficientnet-1.1.0/efficientnet-1.1.0/ -f ./ --no-index","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import os\nimport gc\nimport glob\nimport math\nimport random\nfrom functools import partial, reduce\nfrom tqdm.notebook import tqdm\nimport warnings\n\nimport pydicom\n\nimport scipy\nimport numpy as np\nimport pandas as pd\n\nimport lightgbm as lgb\nimport tensorflow as tf\nprint(f\"TF version: {tf.__version__}\")\nimport efficientnet.tfkeras as efn\n\nfrom sklearn.preprocessing import OneHotEncoder, OrdinalEncoder\nfrom sklearn.model_selection import GroupKFold, KFold, train_test_split\nfrom sklearn.metrics import mean_squared_error\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"import re\n\n\ndef tryint(s):\n    try:\n        return int(s)\n    except ValueError:\n        return s\n    \ndef alphanum_key(s):\n    \"\"\" Turn a string into a list of string and number chunks.\n        \"z23a\" -> [\"z\", 23, \"a\"]\n    \"\"\"\n    return [ tryint(c) for c in re.split('([0-9]+)', s) ]\n\ndef sort_nicely(l):\n    \"\"\" Sort the given list in the way that humans expect.\n    \"\"\"\n    l.sort(key=alphanum_key)\n\nclass OSICTrainDataset:\n    def __init__(self, df):\n        self.df = df\n        self.create_base_df()\n        self._clean_dataset()\n        self._add_base_features()\n        self._add_col_id()\n        self.__sort_by_id()\n        \n    def create_base_df(self):\n        temp_dff = self.df.copy()\n        temp_dff['rank'] = temp_dff.groupby('Patient')['Weeks'].rank(method='min')\n        temp_dff = temp_dff[temp_dff['rank'] == 1]\n        temp_dff = temp_dff.drop_duplicates(subset='Patient')\n        self.base_df = temp_dff\n\n    def _clean_dataset(self):\n        \"\"\"\n        Preprocessing steps:\n            1. Drop duplicate Patient-Weeks combination\n        \"\"\"\n        self.__drop_duplicates()\n        \n    def __drop_duplicates(self):\n        before = self.df.shape[0]\n        self.df = self.df.drop_duplicates(subset=['Patient', 'Weeks'], keep='first').reset_index(drop=True)\n        after = self.df.shape[0]\n        print(f\"Dropped {before-after} rows of duplicate 'Patient-Weeks' values.\")\n        \n    def _add_base_features(self):\n        before = self.df.shape\n        temp_dff = self.df.copy()\n        temp_dff['rank'] = temp_dff.groupby('Patient')['Weeks'].rank(method='min')\n\n        all_dfs = []\n        for rank in sorted(temp_dff['rank'].unique()):\n            all_dfs.append(self.__get_ranked_base_features(temp_dff, rank))\n        self.df = pd.concat(all_dfs)\n        self.df = self.df.reset_index(drop=True)\n        after = self.df.shape\n        print(f\"Before-After shape adding base features: {before} {after}\")\n        \n    def __get_ranked_base_features(self, temp_dff, rank):\n        temp_df = temp_dff[temp_dff['rank'] == rank].reset_index(drop=True)\n        temp_df = temp_df.drop(['Sex', 'SmokingStatus', 'rank'], axis=1)\n        temp_df = temp_df.rename(columns={\n            \"FVC\": \"FVC_base\",\n            \"Percent\": \"Percent_base\",\n            \"Age\": \"Age_base\",\n            \"Weeks\": \"Weeks_base\"\n        })\n        temp_df = self.df[['Patient', 'Weeks', 'FVC', 'Sex', 'SmokingStatus']].merge(\n            temp_df,\n            how='inner',\n            on='Patient',\n        )\n        temp_df['Weeks_passed'] = temp_df['Weeks'] - temp_df['Weeks_base']\n        temp_df = temp_df[temp_df['Weeks_passed'] != 0]\n        return temp_df\n    \n    def _add_col_id(self):\n        col_id = 'Patient_Week'\n        self.df[col_id] = self.df['Patient'] + '_' + self.df['Weeks'].astype(str)\n        print(f\"ID column '{col_id}' added. After adding shape: {self.df.shape}\")\n        \n    def __sort_by_id(self):\n        self.df = self.df.sort_values(by='Patient').reset_index(drop=True)\n        \n        \nclass OSICTestDataset:\n    def __init__(self, test_df, submission_df):\n        self.df = test_df\n        self.submission_df = submission_df\n        self._prepare_test_df()\n        self.__sort_by_id()\n    \n    def _prepare_test_df(self):\n        before = self.df.shape\n        self.submission_df[['Patient', 'Weeks']] = self.submission_df['Patient_Week'].str.split('_', expand=True)\n        self.df = self.submission_df.drop(['FVC', 'Confidence'], axis=1).merge(\n            self.df.rename(columns={\n                \"FVC\": \"FVC_base\",\n                \"Percent\": \"Percent_base\",\n                \"Age\": \"Age_base\",\n                \"Weeks\": \"Weeks_base\"\n            }),\n            how='left',\n            on='Patient'\n            )\n        self.df['Weeks'] = self.df['Weeks'].astype(int)\n        self.df['Weeks_passed'] = self.df['Weeks'] - self.df['Weeks_base']\n        self.df = self.df.reset_index(drop=True)\n        after = self.df.shape\n        print(f\"Before-After shape adding base features: {before} {after}\")\n        \n    def __sort_by_id(self):\n        self.df = self.df.sort_values(by='Patient').reset_index(drop=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"basepath = \"../input/osic-pulmonary-fibrosis-progression/\"\ntrain_df = pd.read_csv(f\"{basepath}train.csv\")\ntest_df = pd.read_csv(f\"{basepath}test.csv\")\nsubmission_df = pd.read_csv(f\"{basepath}sample_submission.csv\")\nprint(train_df.shape, test_df.shape, submission_df.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_dataset = OSICTrainDataset(train_df)\ntest_dataset = OSICTestDataset(test_df, submission_df)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"bbox_map = pd.read_csv('../input/osic-manual-bbox/threshold_all.csv')\nprint(bbox_map.shape)\ncode_filter = (bbox_map['x'] == 0) & (bbox_map['y'] == 0) & (bbox_map['width'] == 1) & (bbox_map['height'] == 1)\nbbox_map = bbox_map[~code_filter]\nprint(bbox_map.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for x in test_df.Patient:\n    print(x in bbox_map.patient.tolist())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import collections\n\ndef create_scaler(min_, max_):\n    def scalar_scaler(val):\n        result = (val - min_) / (max_ - min_)\n        result = max(0.0, result)\n        result = min(1.0, result)\n        return result\n    return scalar_scaler\n\ndef osic_cat_encoder(cat):\n    mapper = {\n        'Male': 0,\n        'Female': 1,\n        'Never smoked': [0, 0],\n        'Ex-smoker': [0, 1],\n        'Currently smokes': [1, 0],\n    }\n    return mapper.get(cat, [1, 1])\n\ndef flatten(l):\n    for el in l:\n        if isinstance(el, collections.Iterable) and not isinstance(el, (str, bytes)):\n            yield from flatten(el)\n        else:\n            yield el","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## CNN-MLP Dataset Generator"},{"metadata":{"trusted":true},"cell_type":"code","source":"base_df = train_dataset.df.copy().reset_index(drop=True)\nbase_test_df = test_dataset.df.copy().reset_index(drop=True)\nprint(base_df.shape, base_test_df.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"BAD_ID = ['ID00011637202177653955184', 'ID00052637202186188008618']\nbase_df = base_df[~base_df.Patient.isin(BAD_ID)]\nbase_test_df = base_test_df[~base_test_df.Patient.isin(BAD_ID)]\nprint(base_df.shape, base_test_df.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cols_pred_num = ['FVC_base', 'Percent_base', 'Age_base', 'Weeks_passed']\ncols_pred_cat = ['Sex', 'SmokingStatus']\ncols_pred = cols_pred_num + cols_pred_cat\n\nfor col in cols_pred_num:\n    print(f\"Predictor {col}\")\n    min_ = min(base_df[col].min(), test_dataset.df[col].min())\n    max_ = max(base_df[col].max(), test_dataset.df[col].max())\n    print(min_, max_)\n    scaler = create_scaler(min_, max_)\n    base_df[col] = base_df[col].apply(scaler)\n    base_test_df[col] = base_test_df[col].apply(scaler)\nfor col in cols_pred_cat:\n    base_df[col] = base_df[col].apply(osic_cat_encoder)\n    base_test_df[col] = base_test_df[col].apply(osic_cat_encoder)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Modelling"},{"metadata":{"trusted":true},"cell_type":"code","source":"import tensorflow.keras.backend as K\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.layers import Dense, Input, Flatten, Conv2D, MaxPooling2D, GlobalAveragePooling2D, Concatenate\nfrom tensorflow.keras.optimizers import Adam, RMSprop, SGD\nfrom tensorflow.keras.callbacks import EarlyStopping, ModelCheckpoint, LearningRateScheduler\n\n\nquantiles = [0.2, 0.5, 0.8]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def score(y_true, y_pred):\n    y_true = tf.dtypes.cast(y_true, tf.float32)\n    y_pred = tf.dtypes.cast(y_pred, tf.float32)\n    C1, C2 = tf.constant(70, dtype='float32'), tf.constant(1000, dtype=\"float32\")\n    sigma = y_pred[:, 2] - y_pred[:, 0]\n    fvc_pred = y_pred[:, 1]\n    \n    #sigma_clip = sigma + C1\n    sigma_clip = tf.maximum(sigma, C1)\n    delta = tf.abs(y_true[:, 0] - fvc_pred)\n    delta = tf.minimum(delta, C2)\n    sq2 = tf.sqrt( tf.dtypes.cast(2, dtype=tf.float32) )\n    metric = (delta / sigma_clip)*sq2 + tf.math.log(sigma_clip* sq2)\n    return -K.mean(metric)\n\ndef mloss(_lambda):\n    def loss(y_true, y_pred):\n        return _lambda * qloss(y_true, y_pred) + (1 - _lambda)*score(y_true, y_pred)\n    return loss\n\ndef qloss(y_true, y_pred):\n    q = tf.constant(quantiles, dtype=tf.float32)\n    e = y_true - y_pred\n    v = tf.maximum(q*e, (q-1)*e)\n    return K.mean(v)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def effnet_mlp_model(input_shape, output_shape):\n    inp = Input(shape=input_shape)\n    base_model = efn.EfficientNetB5(\n        input_shape=input_shape,\n        weights=None,\n        include_top=False,\n    )\n    x1 = base_model(inp)\n    x1 = GlobalAveragePooling2D()(x1)\n    \n    inp2 = Input(shape=(7,))\n    x2 = Dense(100, activation='relu')(inp2)\n    x2 = Dense(100, activation='relu')(x2)\n    \n    x = Concatenate()([x1, x2])\n    output = Dense(output_shape, activation='linear', name='output')(x)\n    \n    model = Model(inputs=[inp, inp2], outputs=output, name='cnn_mlp_only_mid')\n    return model","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model = effnet_mlp_model(input_shape=(512, 512, 1), output_shape=len(quantiles))\nmodel.compile(loss=mloss(0.65),\n              optimizer=Adam(learning_rate=0.02),\n              metrics=[score])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model.load_weights(\"../input/osic-pretrained/cnn_mlp_only_mid.h5\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Make Prediction"},{"metadata":{"trusted":true},"cell_type":"code","source":"from skimage.filters import threshold_otsu, median\nfrom skimage.segmentation import clear_border\nfrom skimage.transform import resize\nfrom skimage import morphology\nfrom scipy.ndimage import binary_fill_holes\n\n\ndef preprocess_slice(slice, bbox_map, image_type, patient_id, img_width=512, img_height=512):\n    slice_data = to_HU(slice)\n\n    slice_data = segment_lung(slice_data, image_type, morphological_segmentation)\n    if patient_id not in bbox_map.patient.values:\n        bbox = infer_bounding_box(slice_data)\n    else:\n        x, y, width, height = bbox_map.loc[bbox_map.patient == patient_id, ['x', 'y', 'width', 'height']].values[0]\n        bbox = BoundingBox((x, y), width, height)\n\n    slice_data = crop_recenter(slice_data, bbox)\n    slice_data = rescale(slice_data)\n    resized_slice_data = resize(slice_data, (img_height, img_width), anti_aliasing=True)\n    resized_slice_data = np.expand_dims(resized_slice_data, axis=-1)\n    resized_slice_data = resized_slice_data.astype(np.float32)\n    return resized_slice_data\n\n\ndef segment_lung(slice_data, image_type, segment_func):\n    if image_type == 'zero':\n        slice_data[slice_data == 0] = -1000\n    segmented_image = segment_func(threshold_slices_data(slice_data, low=-1000, high=-400))\n    return segmented_image\n\ndef morphological_segmentation(img):\n    segmented_img, _ = lung_segment(img)\n    return segmented_img\n\ndef lung_segment(img):\n    try:\n        thresh = threshold_otsu(img)\n        binary = img <= thresh\n\n        lungs = median(clear_border(binary))\n        lungs = morphology.binary_closing(lungs, selem=morphology.disk(7))\n        lungs = binary_fill_holes(lungs)\n\n        final = lungs*img\n        final[final == 0] = np.min(img)\n        return final, lungs\n    except:\n        print(\"Segmentation failed. Returning original image.\")\n        return img, img\n\n\ndef infer_bounding_box(segmented_image):\n    try:\n        y_match, x_match = np.where(segmented_image != -1000)\n        y_min, x_min = y_match.min(), x_match.min()\n        y_max, x_max = y_match.max(), x_match.max()\n        width = abs(x_max - x_min)\n        height = abs(y_max - y_min)\n    except:\n        print(\"Inferring boundigng box failed, returning whole image\")\n        x_min = 0\n        y_min = 0\n        height, width = segmented_image.shape\n        height = height - 1\n        width = width - 1\n    return BoundingBox((x_min, y_min), width, height)\n\nclass BoundingBox:\n    \"\"\"Initiation of bbox follows matplotlib Rectangle patch\"\"\"\n    def __init__(self, xy, width, height):\n        self.x, self.y = xy\n        self.width = width\n        self.height = height\n        \n    @property\n    def attribute_list(self):\n        return [(self.x, self.y), self.width, self.height]\n    \n    def __repr__(self):\n        return f\"Bbox (bottom left width height): {self.x} {self.y} {self.width} {self.height}\"\n\n    \ndef crop_recenter(image, bbox, pad_value=-1000):\n    x, y, width, height = bbox.x, bbox.y, bbox.width, bbox.height\n    cropped_image = image[ y:y+height, x:x+width ]\n    out_height, out_width = image.shape\n    \n    padded_image = np.ones(image.shape, dtype=np.int16) * pad_value\n    x_start = (out_width - width) // 2\n    y_start = (out_height - height) // 2\n    padded_image[ y_start:y_start+height, x_start:x_start+width ] = cropped_image\n    return padded_image\n\n\ndef threshold_slices_data(slices_data, low=-1000, high=-400):\n    copy = slices_data.copy()\n    copy[copy < low] = low\n    copy[copy > high] = high\n    return copy\n\n\ndef rescale(slice_data):\n    min_, max_ = slice_data.min(), slice_data.max()\n    rescaled = (slice_data-min_) / (max_-min_)\n\n    if np.isfinite(rescaled).all():\n        return rescaled\n    else:\n        print(\"Rescaling failed, returning np.zeros() with original shape.\")\n        return np.zeros(slice_data.shape)\n\n\ndef to_HU(slice):\n    intercept, slope = slice.RescaleIntercept, slice.RescaleSlope\n\n    slice_data = slice.pixel_array.astype(np.int16)\n    slice_data[slice_data <= -1000] = 0\n\n    if slope != 1:\n        slice_data = slope * slice_data.astype(np.float64)\n        slice_data = slice_data.astype(np.int16)\n\n    slice_data += np.int16(intercept)\n    return slice_data","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"DOUBLE_IDS = ['ID00078637202199415319443']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"image_mapper = {}\nfor patient_id in base_test_df.Patient.unique():\n    slices = []\n    filepaths = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/test/{patient_id}/')\n    sort_nicely(filepaths)\n    for filepath in filepaths:\n        slices.append(pydicom.dcmread(f'../input/osic-pulmonary-fibrosis-progression/test/{patient_id}/{filepath}'))\n    if patient_id in DOUBLE_IDS:\n        x = slices[(len(slices)-1) //4]\n    else:\n        x = slices[(len(slices)-1) //2]\n    \n    intercept = x.RescaleIntercept\n    if intercept == 0:\n        image_type = 'zero'\n    else:\n        image_type = 'not-zero'\n    x = preprocess_slice(x, bbox_map, image_type, patient_id, 512, 512)\n    \n    image_mapper.update({\n        f'{patient_id}': x,\n    })","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"base_test_df['FVC'] = 2000\nbase_test_df['Confidence'] = 300\nfor patient_week in base_test_df['Patient_Week'].unique():\n    selector = base_test_df['Patient_Week'] == patient_week\n    dff = base_test_df[selector]\n    patient_id = dff['Patient'].values[0]\n    x = image_mapper[patient_id]\n    vector = np.array(list(flatten(dff[cols_pred].values[0])))\n    x = np.expand_dims(x, axis=0)\n    vector = np.expand_dims(vector, axis=0)\n    \n    pred = model.predict([x, vector])\n    FVC = pred[0][1]\n    confidence = pred[0][2] - pred[0][0]\n    base_test_df.loc[selector, 'FVC'] = FVC\n    base_test_df.loc[selector, 'Confidence'] = confidence","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission_df = submission_df[['Patient_Week']].merge(\n    base_test_df[['Patient_Week', 'FVC', 'Confidence']],\n    how='inner',\n    on='Patient_Week'\n)\nsubmission_df.to_csv(\"submission.csv\", header=True, index=False)","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}