{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import tensorflow as tf\nimport pandas as pd\nimport numpy as np\nimport os\nimport random\nimport pathlib\nimport pickle\nimport time\n\nimport pydicom\nimport matplotlib.pyplot as plt\n\nimport tensorflow.keras.backend\nfrom tensorflow import keras as K\nfrom tensorflow.keras import layers as L\nfrom skimage.transform import resize\nfrom scipy.ndimage import zoom\nimport scipy.ndimage as ndimage\nfrom skimage import measure, morphology, segmentation\n\nfrom math import ceil","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"raw_test = pd.read_csv('/kaggle/input/osic-pulmonary-fibrosis-progression/test.csv')\nX_prediction = pd.read_csv('/kaggle/input/osic-pulmonary-fibrosis-progression/sample_submission.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#Constant\nTEST_PATH = '/kaggle/input/osic-pulmonary-fibrosis-progression/test'\n\n#For data prep\nDESIRED_SIZE = (30,256,256)\nBATCH_SIZE = 256\n\nMASK_ITERATION = 4\n\nclip_bounds = (-1000, 200)\npre_calculated_mean = 0.02865046213070556","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# def create_submission(value):\n#     if value==0:\n#         sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':[0 for i in X_prediction.index], 'Confidence': [10000 for i in X_prediction.index]})\n#     elif value==1:\n#         sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':[100 for i in X_prediction.index], 'Confidence': [5000 for i in X_prediction.index]})\n#     elif value==2:\n#         sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':[500 for i in X_prediction.index], 'Confidence': [1000 for i in X_prediction.index]})\n#     elif value==3:\n#         sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':[1000 for i in X_prediction.index], 'Confidence': [500 for i in X_prediction.index]})\n#     elif value==4:\n#         sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':[2000 for i in X_prediction.index], 'Confidence': [100 for i in X_prediction.index]})\n#     return sub\n        ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Tabular prep"},{"metadata":{"trusted":true},"cell_type":"code","source":"X_prediction['Patient'] = X_prediction['Patient_Week'].str.extract(r'(.*)_.*')\nX_prediction['Weeks'] = X_prediction['Patient_Week'].str.extract(r'.*_(.*)').astype(int)\nX_prediction = X_prediction[['Patient', 'Weeks', 'Patient_Week']]\nrename_cols = {'Weeks_y':'Min_week', 'Weeks_x': 'Weeks', 'FVC':'Base_FVC'}\nX_prediction = X_prediction.merge(raw_test, how='left', left_on='Patient', right_on='Patient').rename(columns=rename_cols)[['Patient', 'Min_week', 'Base_FVC', 'Percent', 'Age', 'Sex', 'SmokingStatus', 'Weeks', 'Patient_Week']].reset_index(drop=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"X_prediction['Base_week'] = X_prediction['Weeks'] - X_prediction['Min_week']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from sklearn.preprocessing import OneHotEncoder as SklearnOneHotEncoder\n\n#Override OneHotEncoder to have the column names created automatically (Ex-smoker, Never Smoked...) \nclass OneHotEncoder(SklearnOneHotEncoder):\n    def __init__(self, **kwargs):\n        super(OneHotEncoder, self).__init__(**kwargs)\n        self.fit_flag = False\n\n    def fit(self, X, **kwargs):\n        out = super().fit(X)\n        self.fit_flag = True\n        return out\n\n    def transform(self, X, categories, index='', name='', **kwargs):\n        sparse_matrix = super(OneHotEncoder, self).transform(X)\n        new_columns = self.get_new_columns(X=X, name=name, categories=categories)\n        d_out = pd.DataFrame(sparse_matrix.toarray(), columns=new_columns, index=index)\n        return d_out\n\n    def fit_transform(self, X, categories, index, name, **kwargs):\n        self.fit(X)\n        return self.transform(X, categories=categories, index=index, name=name)\n\n    def get_new_columns(self, X, name, categories):\n        new_columns = []\n        for j in range(len(categories)):\n            new_columns.append('{}_{}'.format(name, categories[j]))\n        return new_columns\n\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.exceptions import NotFittedError\n\ndef standardisation(x, u, s):\n    return (x-u)/s\n\ndef normalization(x, ma, mi):\n    return (x-mi)/(ma-mi)\n\nclass data_preparation():\n    def __init__(self, bool_normalization=True,bool_standard=False):\n        self.enc_sex = LabelEncoder()\n        self.enc_smok = LabelEncoder()\n        self.onehotenc_smok = OneHotEncoder()\n        self.standardisation = bool_standard\n        self.normalization = bool_normalization\n        \n        \n    def __call__(self, data_untransformed):\n        data = data_untransformed.copy(deep=True)\n        \n        #For the test set/Already fitted\n        try:\n            data['Sex'] = self.enc_sex.transform(data['Sex'].values)\n            data['SmokingStatus'] = self.enc_smok.transform(data['SmokingStatus'].values)\n            data = pd.concat([data.drop(columns=['SmokingStatus']), self.onehotenc_smok.transform(data['SmokingStatus'].values.reshape(-1,1), categories=self.enc_smok.classes_, name='', index=data.index).astype(int)], axis=1)\n            \n            #Standardisation\n            if self.standardisation:\n                data['Base_week'] = standardisation(data['Base_week'],self.base_week_mean,self.base_week_std)\n                data['Base_FVC'] = standardisation(data['Base_FVC'],self.base_fvc_mean,self.base_fvc_std)\n                data['Base_percent'] = standardisation(data['Base_percent'],self.base_percent_mean,self.base_percent_std)\n                data['Age'] = standardisation(data['Age'],self.age_mean,self.age_std)\n                data['Weeks'] = standardisation(data['Weeks'],self.weeks_mean,self.weeks_std)\n            \n            #Normalization\n            if self.normalization:\n                data['Base_week'] = normalization(data['Base_week'],self.base_week_max,self.base_week_min)\n                data['Base_FVC'] = normalization(data['Base_FVC'],self.base_fvc_max,self.base_fvc_min)\n                data['Percent'] = normalization(data['Percent'],self.base_percent_max,self.base_percent_min)\n                data['Age'] = normalization(data['Age'],self.age_max,self.age_min)\n                data['Weeks'] = normalization(data['Weeks'],self.weeks_max,self.weeks_min)\n                data['Min_week'] = normalization(data['Min_week'],self.base_week_max,self.base_week_min)\n\n        #For the train set/Not yet fitted    \n        except NotFittedError:\n            data['Sex'] = self.enc_sex.fit_transform(data['Sex'].values)\n            data['SmokingStatus'] = self.enc_smok.fit_transform(data['SmokingStatus'].values)\n            data = pd.concat([data.drop(columns=['SmokingStatus']), self.onehotenc_smok.fit_transform(data['SmokingStatus'].values.reshape(-1,1), categories=self.enc_smok.classes_, name='', index=data.index).astype(int)], axis=1)\n            \n            #Standardisation\n            if self.standardisation:\n                self.base_week_mean = data['Base_week'].mean()\n                self.base_week_std = data['Base_week'].std()\n                data['Base_week'] = standardisation(data['Base_week'],self.base_week_mean,self.base_week_std)\n\n                self.base_fvc_mean = data['Base_FVC'].mean()\n                self.base_fvc_std = data['Base_FVC'].std()\n                data['Base_FVC'] = standardisation(data['Base_FVC'],self.base_fvc_mean,self.base_fvc_std)\n\n                self.base_percent_mean = data['Base_percent'].mean()\n                self.base_percent_std = data['Base_percent'].std()\n                data['Base_percent'] = standardisation(data['Base_percent'],self.base_percent_mean,self.base_percent_std)\n\n                self.age_mean = data['Age'].mean()\n                self.age_std = data['Age'].std()\n                data['Age'] = standardisation(data['Age'],self.age_mean,self.age_std)\n\n                self.weeks_mean = data['Weeks'].mean()\n                self.weeks_std = data['Weeks'].std()\n                data['Weeks'] = standardisation(data['Weeks'],self.weeks_mean,self.weeks_std)\n\n                \n            #Normalization\n            if self.normalization:\n                self.base_week_min = data['Base_week'].min()\n                self.base_week_max = data['Base_week'].max()\n                data['Base_week'] = normalization(data['Base_week'],self.base_week_max,self.base_week_min)\n\n                self.base_fvc_min = data['Base_FVC'].min()\n                self.base_fvc_max = data['Base_FVC'].max()\n                data['Base_FVC'] = normalization(data['Base_FVC'],self.base_fvc_max,self.base_fvc_min)\n\n                self.base_percent_min = data['Percent'].min()\n                self.base_percent_max = data['Percent'].max()\n                data['Percent'] = normalization(data['Percent'],self.base_percent_max,self.base_percent_min)\n\n                self.age_min = data['Age'].min()\n                self.age_max = data['Age'].max()\n                data['Age'] = normalization(data['Age'],self.age_max,self.age_min)\n\n                self.weeks_min = data['Weeks'].min()\n                self.weeks_max = data['Weeks'].max()\n                data['Weeks'] = normalization(data['Weeks'],self.weeks_max,self.weeks_min)\n                \n                self.base_week_min = data['Min_week'].min()\n                self.base_week_max = data['Min_week'].max()\n                data['Min_week'] = normalization(data['Min_week'],self.base_week_max,self.base_week_min)\n\n            \n        return data","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pickefile = open('/kaggle/input/data-preparation-for-osic/data_prep', 'rb')\ndata_prep = pickle.load(pickefile)\npickefile.close()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"X_prediction = data_prep(X_prediction).sort_values('Patient')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Data Generator"},{"metadata":{"trusted":true},"cell_type":"code","source":"class ConvertToHU:\n    def __call__(self, imgs, dicom):\n\n#         img_type = dicom.ImageType\n#         is_hu = img_type[0] == 'ORIGINAL' and not (img_type[2] == 'LOCALIZER')\n        # if not is_hu:\n        #     warnings.warn(f'Patient {data.PatientID} CT Scan not cannot be'\n        #                   f'converted to Hounsfield Units (HU).')\n\n        intercept = dicom.RescaleIntercept\n        slope = dicom.RescaleSlope\n        imgs = (np.array(imgs.to_list()) * slope + intercept).astype(np.int16)\n        return imgs\nconvertohu = ConvertToHU()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class Clip:\n    def __init__(self, bounds=(-1000, 500)):\n        self.min = min(bounds)\n        self.max = max(bounds)\n\n    def __call__(self, image):\n        image[image < self.min] = self.min\n        image[image > self.max] = self.max\n        return image\nclip = Clip(clip_bounds)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class MaskWatershed:\n    def __init__(self, min_hu, iterations):\n        self.min_hu = min_hu\n        self.iterations = iterations\n\n    def __call__(self, image, dicom):\n        \n        # Structuring element used for the filter\n        blackhat_struct = [[0, 0, 1, 1, 1, 0, 0],\n                           [0, 1, 1, 1, 1, 1, 0],\n                           [1, 1, 1, 1, 1, 1, 1],\n                           [1, 1, 1, 1, 1, 1, 1],\n                           [1, 1, 1, 1, 1, 1, 1],\n                           [0, 1, 1, 1, 1, 1, 0],\n                           [0, 0, 1, 1, 1, 0, 0]]\n\n        blackhat_struct = ndimage.iterate_structure(blackhat_struct, self.iterations)\n        stack = []\n        for slice_idx in range(image.shape[0]):\n            sliced = image[slice_idx]\n            stack.append(self.seperate_lungs(sliced, blackhat_struct, self.min_hu,\n                                             self.iterations))\n        \n        \n        return np.stack(stack)\n\n    @staticmethod\n    def seperate_lungs(image, blackhat_struct, min_hu=min(clip_bounds), iterations=2):\n#         h, w = image.shape[0], image.shape[1]\n\n        marker_internal, marker_external, marker_watershed = MaskWatershed.generate_markers(image)\n        \n        # Sobel-Gradient\n        sobel_filtered_dx = ndimage.sobel(image, 1)\n        sobel_filtered_dy = ndimage.sobel(image, 0)\n        sobel_gradient = np.hypot(sobel_filtered_dx, sobel_filtered_dy)\n        sobel_gradient *= 255.0 / np.max(sobel_gradient)\n\n        watershed = morphology.watershed(sobel_gradient, marker_watershed)\n\n        outline = ndimage.morphological_gradient(watershed, size=(3,3)).astype(bool)\n#         outline = outline.astype(bool)\n\n#         # Structuring element used for the filter\n#         blackhat_struct = [[0, 0, 1, 1, 1, 0, 0],\n#                            [0, 1, 1, 1, 1, 1, 0],\n#                            [1, 1, 1, 1, 1, 1, 1],\n#                            [1, 1, 1, 1, 1, 1, 1],\n#                            [1, 1, 1, 1, 1, 1, 1],\n#                            [0, 1, 1, 1, 1, 1, 0],\n#                            [0, 0, 1, 1, 1, 0, 0]]\n\n#         blackhat_struct = ndimage.iterate_structure(blackhat_struct, iterations)\n\n        # Perform Black Top-hat filter\n        outline += ndimage.black_tophat(outline, structure=blackhat_struct)\n\n        lungfilter = np.bitwise_or(marker_internal, outline)\n        lungfilter = ndimage.morphology.binary_closing(lungfilter, structure=np.ones((5,5)), iterations=3)\n\n        segmented = np.where(lungfilter == 1, image, min_hu * np.ones((image.shape[0], image.shape[1])))\n\n        return segmented  #, lungfilter, outline, watershed, sobel_gradient\n\n    @staticmethod\n    def generate_markers(image, threshold=-400):\n#         h, w = image.shape[0], image.shape[1]\n\n#         marker_internal = image < threshold\n#         marker_internal = segmentation.clear_border(image < threshold)\n        marker_internal_labels = measure.label(segmentation.clear_border(image < threshold))\n        \n        areas = [r.area for r in measure.regionprops(marker_internal_labels)]\n        areas.sort()\n        \n        if len(areas) > 2:\n            for region in measure.regionprops(marker_internal_labels):\n                if region.area < areas[-2]:\n                    for coordinates in region.coords:\n                        marker_internal_labels[coordinates[0], coordinates[1]] = 0\n        \n        marker_internal = marker_internal_labels > 0\n\n        # Creation of the External Marker\n        external_a = ndimage.binary_dilation(marker_internal, iterations=10)\n        external_b = ndimage.binary_dilation(marker_internal, iterations=55)\n        marker_external = external_b ^ external_a\n\n        # Creation of the Watershed Marker\n        marker_watershed = np.zeros((image.shape[0], image.shape[1]), dtype=np.int)\n        marker_watershed += marker_internal * 255\n        marker_watershed += marker_external * 128\n\n        return marker_internal, marker_external, marker_watershed\nmaskwatershed = MaskWatershed(min_hu=min(clip_bounds), iterations=MASK_ITERATION)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class Normalize:\n    def __init__(self, bounds=(-1000, 500)):\n        self.min = min(bounds)\n        self.max = max(bounds)\n\n    def __call__(self, image):\n        image = image.astype(np.float)\n        image = (image - self.min) / (self.max - self.min)\n        return image\n     \nclass ZeroCenter:\n    def __init__(self, pre_calculated_mean):\n        self.pre_calculated_mean = pre_calculated_mean\n\n    def __call__(self, image):\n        return image - self.pre_calculated_mean\n\nnormalize = Normalize(bounds=clip_bounds)\nzerocenter = ZeroCenter(pre_calculated_mean=pre_calculated_mean)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def sort_function(x): #Get the files in the right order\n    return int(x.split('.')[0])\n\ndef _read(path, patients=[], desired_size=(60,512,512)):\n    X = np.empty(np.concatenate(([len(patients), 1], np.array(desired_size))))\n    i=0\n    for patient in patients:\n        dicom = pydicom.dcmread(path+'/'+patient+'/'+os.listdir(path+'/'+patient)[0])\n        list_patient_files = sorted([i for i in os.listdir(path+'/'+patient)], key=sort_function)\n#         list_patient_files = list_patient_files[::max(int(len(list_patient_files)/40),1)][:30]\n        df = pd.DataFrame(list_patient_files).apply(lambda x:path+'/'+patient+'/'+x)\n        df = convertohu(df.iloc[:,0].apply(lambda x: pydicom.dcmread(x).pixel_array), dicom)\n        df = zoom(df, np.array(DESIRED_SIZE)/np.array(df.shape), mode='nearest')\n#         df = np.array(df.iloc[:,0].apply(lambda x: resize(pydicom.dcmread(x).pixel_array, DESIRED_SIZE[1:3])).to_list())\n#         df = zoom(df, np.array(DESIRED_SIZE)/np.array(df.shape), mode='nearest')\n        X[i,0,:,:,:] = zerocenter(normalize(maskwatershed(clip(df), dicom)))\n        i+=1\n            \n    #apply on resize on a batch and not on one picture\n    return X","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# X = _read('/kaggle/input/osic-pulmonary-fibrosis-progression/train', patients=['ID00007637202177411956430'], desired_size=DESIRED_SIZE)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class DataGenerator(K.utils.Sequence):\n    \n    def on_epoch_end(self):#Indices=np.arange(size of dataset)\n        self.indices = np.arange(len(self.list_IDs))\n    \n    def __len__(self):\n        return int(ceil(len(self.indices) / self.batch_size))\n    \n    def __init__(self, train, list_IDs, batch_size=1, desired_size=(10,512,512), img_path=TEST_PATH, *args, **kwargs):\n        self.train = train\n        self.list_IDs = list_IDs\n        self.batch_size = batch_size\n        self.desired_size = desired_size\n        self.img_path = img_path\n        self.on_epoch_end()\n    \n    def __getitem__(self, index):\n        indices = self.indices[index*self.batch_size:(index+1)*self.batch_size]\n        list_IDs_temp = [self.list_IDs[k] for k in indices]\n        \n        patients = self.train.loc[list_IDs_temp, 'Patient'].unique()\n        imgs = _read(self.img_path, patients=patients, desired_size=self.desired_size)\n        return self.train.loc[list_IDs_temp, :].reset_index(drop=True), np.transpose(imgs, (0, 2, 3, 4, 1))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pred_generator = DataGenerator(X_prediction, X_prediction.index, batch_size=BATCH_SIZE, desired_size=DESIRED_SIZE, img_path=TEST_PATH)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Load model and predict"},{"metadata":{"trusted":true},"cell_type":"code","source":"# model = K.models.load_model('/kaggle/input/cnn-for-latent-features/model',custom_objects={'mloss':mloss(LAMBDA_LOSS)})\nmodel = K.models.load_model('/kaggle/input/3d-cnn-mlp/model_4',compile=False)\nCNN = K.models.load_model('/kaggle/input/3d-cnn-mlp/CNN_4',compile=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# start_time = time.perf_counter()\nSELECTED_COLUMNS = ['Weeks', 'Percent', 'Age', 'Sex', 'Min_week', 'Base_FVC','Base_week', '_Currently smokes', '_Ex-smoker', '_Never smoked']\ny_prediction = np.array([[0,0,0]])\nfor X1, X2 in pred_generator:\n    out_imgs = CNN(tf.convert_to_tensor(X2))\n    X_imgs = tf.concat([tf.repeat(tf.reshape(out_imgs[i], (1,-1)), X1['Patient'].value_counts()[i], axis=0) for i in range(len(X1['Patient'].value_counts()))], axis=0)\n    X_tabular = model.get_layer('dense_to_freeze2')(model.get_layer('dense_to_freeze1')(tf.convert_to_tensor(np.asarray(X1[SELECTED_COLUMNS]))))\n    inp_mlp = model.get_layer('concatenate')([X_imgs, X_tabular])\n    inp_mlp = model.layers[-1](model.layers[-2](model.layers[-3](inp_mlp)))\n    y_prediction = np.append(y_prediction, inp_mlp.numpy(), axis=0)\ny_prediction = y_prediction[1:]\n# tf.print(\"Execution time:\", time.perf_counter() - start_time)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub = pd.DataFrame(data={'Patient_Week':X_prediction['Patient_Week'], 'FVC':y_prediction[:,1], 'Confidence': (y_prediction[:,2]-y_prediction[:,0])})\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}