{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Overview & Remarks (Important Disclaimer!)\n\n### - There is nothing novel here as I forked the notebook from [here](https://www.kaggle.com/leoisleo1/efficientnets-quantile-regression-inference) where it was originally created by my good bro [Wei Hao](https://www.kaggle.com/khoongweihao)'s [kernel](https://www.kaggle.com/khoongweihao/efficientnets-quantile-regression-inference).\n### - However, here I want to highlight that a small perturbation of the `mloss` function will change your score quite drastically on the LB.\n### - There will be fluctuations in the score."},{"metadata":{},"cell_type":"markdown","source":"# Acknowledgements\n\n- Main inspiration has been taken from the notebook\" LB score by tuning mloss (around -6.811)\" done by Hongnan Gao and the team.\n   - Model for the ensembling of qunatile regression and effnet. Can be found at https://www.kaggle.com/reighns/higher-lb-score-by-tuning-mloss-around-6-811/output?select=submission_regression.csv\n\n- Ulrich GOUE's Osic-Multiple-Quantile-Regression-Starter\n    - Model that uses images can be found at: https://www.kaggle.com/miklgr500/linear-decay-based-on-resnet-cnn\n- Michael Kazachok's Linear Decay (based on ResNet CNN)\n    - Model that uses tabular data can be found at: https://www.kaggle.com/ulrich07/osic-multiple-quantile-regression-starter\n\n- OSIC feature extract from CT by hfut_ybx \n     - Model to extract features from the CT scan https://www.kaggle.com/hfutybx/osic-feature-extract-from-ct/notebook\n     \n- Inspired by metadata analysis done by anarthal in \"medalDICOM metadata EDA\"\n     - Referenced for the understanding of the Metadata. https://www.kaggle.com/anarthal/dicom-metadata-eda\n\n- I have only combined few ways to come up with some additional features\n"},{"metadata":{},"cell_type":"markdown","source":"# Imports"},{"metadata":{"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"!pip install ../input/kerasapplications/keras-team-keras-applications-3b180cb -f ./ --no-index\n!pip install ../input/efficientnet/efficientnet-1.1.0/ -f ./ --no-index","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"conda install -c conda-forge -y gdcm","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install segmentation-models-pytorch","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"import os\nimport cv2\nimport pydicom\nimport pandas as pd\nimport numpy as np \nimport tensorflow as tf \nimport matplotlib.pyplot as plt \nimport random\nfrom tqdm.notebook import tqdm \nfrom sklearn.model_selection import train_test_split, KFold\nfrom sklearn.metrics import mean_absolute_error\nfrom tensorflow_addons.optimizers import RectifiedAdam\nfrom tensorflow.keras import Model\nimport tensorflow.keras.backend as K\nimport tensorflow.keras.layers as L\nimport tensorflow.keras.models as M\nfrom tensorflow.keras.optimizers import Nadam\nimport seaborn as sns\nfrom PIL import Image\nimport scipy.ndimage as ndimage\nfrom scipy.ndimage import zoom\nfrom scipy.stats import kurtosis\nfrom scipy.stats import skew\n\nimport torch\nimport torch.nn as nn\nfrom torch.utils.data import DataLoader\nfrom torch.utils.data.dataset import Dataset\n\nfrom sklearn.model_selection import KFold\n\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\ndef seed_everything(seed=2020):\n    random.seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    np.random.seed(seed)\n    tf.random.set_seed(seed)\n    \nseed_everything(42)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"config = tf.compat.v1.ConfigProto()\nconfig.gpu_options.allow_growth = True\nsession = tf.compat.v1.Session(config=config)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import sys\nimport glob\nsys.path.append('../input/efficientnet-pytorch/EfficientNet-PyTorch-master/')\nsys.path.append('../input/pretrainedmodels/pretrainedmodels-0.7.4/')\nsys.path.append('../input/segmentation-models-pytorch/')\nimport segmentation_models_pytorch as smp","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/train.csv') ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Linear Decay (based on EfficientNets)"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_tab(df):\n    print(df)\n    vector = [(df.Age.values[0] - 30) / 30] \n    \n    if df.Sex.values[0] == 'male':\n       vector.append(0)\n    else:\n       vector.append(1)\n    \n    if df.SmokingStatus.values[0] == 'Never smoked':\n        vector.extend([0,0])\n    elif df.SmokingStatus.values[0] == 'Ex-smoker':\n        vector.extend([1,1])\n    elif df.SmokingStatus.values[0] == 'Currently smokes':\n        vector.extend([0,1])\n    else:\n        vector.extend([1,0])\n    return np.array(vector) ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"A = {} \nTAB = {} \nP = [] \nfor i, p in tqdm(enumerate(train.Patient.unique())):\n    sub = train.loc[train.Patient == p, :] \n    fvc = sub.FVC.values\n    weeks = sub.Weeks.values\n    c = np.vstack([weeks, np.ones(len(weeks))]).T\n    a, b = np.linalg.lstsq(c, fvc)[0]\n    \n    A[p] = a\n    TAB[p] = get_tab(sub)\n    P.append(p)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## CNN for coeff prediction"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_img(path):\n    d = pydicom.dcmread(path)\n    return cv2.resize(d.pixel_array / 2**11, (512, 512))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tensorflow.keras.utils import Sequence\n\nclass IGenerator(Sequence):\n    BAD_ID = ['ID00011637202177653955184', 'ID00052637202186188008618']\n    def __init__(self, keys, a, tab, batch_size=32):\n        self.keys = [k for k in keys if k not in self.BAD_ID]\n        self.a = a\n        self.tab = tab\n        self.batch_size = batch_size\n        \n        self.train_data = {}\n        for p in train.Patient.values:\n            self.train_data[p] = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/')\n    \n    def __len__(self):\n        return 1000\n    \n    def __getitem__(self, idx):\n        x = []\n        a, tab = [], [] \n        keys = np.random.choice(self.keys, size = self.batch_size)\n        for k in keys:\n            try:\n                i = np.random.choice(self.train_data[k], size=1)[0]\n                img = get_img(f'../input/osic-pulmonary-fibrosis-progression/train/{k}/{i}')\n                x.append(img)\n                a.append(self.a[k])\n                tab.append(self.tab[k])\n            except:\n                print(k, i)\n       \n        x,a,tab = np.array(x), np.array(a), np.array(tab)\n        x = np.expand_dims(x, axis=-1)\n        return [x, tab] , a","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tensorflow.keras.layers import (\n    Dense, Dropout, Activation, Flatten, Input, BatchNormalization, GlobalAveragePooling2D, Add, Conv2D, AveragePooling2D, \n    LeakyReLU, Concatenate \n)\nimport efficientnet.tfkeras as efn\n\ndef get_efficientnet(model, shape):\n    models_dict = {\n        'b0': efn.EfficientNetB0(input_shape=shape,weights=None,include_top=False),\n        'b1': efn.EfficientNetB1(input_shape=shape,weights=None,include_top=False),\n        'b2': efn.EfficientNetB2(input_shape=shape,weights=None,include_top=False),\n        'b3': efn.EfficientNetB3(input_shape=shape,weights=None,include_top=False),\n        'b4': efn.EfficientNetB4(input_shape=shape,weights=None,include_top=False),\n        'b5': efn.EfficientNetB5(input_shape=shape,weights=None,include_top=False),\n        'b6': efn.EfficientNetB6(input_shape=shape,weights=None,include_top=False),\n        'b7': efn.EfficientNetB7(input_shape=shape,weights=None,include_top=False)\n    }\n    return models_dict[model]\n\ndef build_model(shape=(512, 512, 1), model_class=None):\n    inp = Input(shape=shape)\n    base = get_efficientnet(model_class, shape)\n    x = base(inp)\n    x = GlobalAveragePooling2D()(x)\n    inp2 = Input(shape=(4,))\n    x2 = tf.keras.layers.GaussianNoise(0.2)(inp2)\n    x = Concatenate()([x, x2]) \n    x = Dropout(0.35)(x) \n    x = Dense(1)(x)\n    model = Model([inp, inp2] , x)\n    \n    weights = [w for w in os.listdir('../input/osic-model-weights') if model_class in w][0]\n    model.load_weights('../input/osic-model-weights/' + weights)\n    return model\n\nmodel_classes = ['b5'] #['b0','b1','b2','b3',b4','b5','b6','b7']\nmodels = [build_model(shape=(512, 512, 1), model_class=m) for m in model_classes]\nprint('Number of models: ' + str(len(models)))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from sklearn.model_selection import train_test_split \n\ntr_p, vl_p = train_test_split(P, \n                              shuffle=True, \n                              train_size= 0.8) ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(list(A.values()));","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def score(fvc_true, fvc_pred, sigma):\n    sigma_clip = np.maximum(sigma, 70) # changed from 70, trie 66.7 too\n    delta = np.abs(fvc_true - fvc_pred)\n    delta = np.minimum(delta, 1000)\n    sq2 = np.sqrt(2)\n    metric = (delta / sigma_clip)*sq2 + np.log(sigma_clip* sq2)\n    return np.mean(metric)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subs = []\n# mid_df = pd.DataFrame([['patient_id', 'actual_fvc', 'pred_fvc', 'score']])\nfor model in models:\n    metric = []\n    for layer in model.layers:\n        if len(layer.weights) > 0:\n            print(layer.name, layer.weights[0].shape, layer.weights[0])\n    for q in tqdm(range(1, 10)):\n        m = []\n        for p in vl_p:\n            x, y = [] , []\n            tab = [] \n\n            if p in ['ID00011637202177653955184', 'ID00052637202186188008618']:\n                continue\n\n            ldir = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/')\n            for i in ldir:\n#                 print(len(ldir), int(i[:-4]), i, p)\n                if int(i[:-4]) / len(ldir) < 0.8 and int(i[:-4]) / len(ldir) > 0.15:\n                    x.append(get_img(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/{i}')) \n                    y.append(i)\n                    tab.append(get_tab(train.loc[train.Patient == p, :])) \n            if len(x) < 1:\n                continue\n            tab = np.array(tab)\n#             print(tab)\n\n            x = np.expand_dims(x, axis=-1) \n            _a = model.predict([x, tab]) \n            a = np.quantile(_a, q / 10)\n\n            percent_true = train.Percent.values[train.Patient == p]\n            fvc_true = train.FVC.values[train.Patient == p]\n            weeks_true = train.Weeks.values[train.Patient == p]\n\n            fvc = a * (weeks_true - weeks_true[0]) + fvc_true[0]\n            percent = percent_true[0] - a * abs(weeks_true - weeks_true[0])\n            m.append(score(fvc_true, fvc, percent))\n#             df2 = pd.DataFrame({\"patient_id\": p,  \"image_num\": y , 'actual_fvc': fvc_true , 'pred_fvc': fvc ,'percent_true': percent_true, 'percent_pred': percent, 'weeks_true': weeks_true, 'score': score })\n#             df2 = pd.DataFrame({\"patient_id\": p, 'actual_fvc': fvc_true , 'pred_fvc': fvc, 'score': score(fvc_true, fvc, percent) })\n#             print(df2)\n#             mid_df.append(df2, ignore_index = True)\n        print(np.mean(m))\n        metric.append(np.mean(m))\n\n    q = (np.argmin(metric) + 1)/ 10\n\n    sub = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/sample_submission.csv') \n    test = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv') \n    A_test, B_test, P_test,W, FVC= {}, {}, {},{},{} \n    STD, WEEK = {}, {} \n    for p in test.Patient.unique():\n        x = [] \n        tab = [] \n        ldir = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/test/{p}/')\n        for i in ldir:\n            if int(i[:-4]) / len(ldir) < 0.8 and int(i[:-4]) / len(ldir) > 0.15:\n                x.append(get_img(f'../input/osic-pulmonary-fibrosis-progression/test/{p}/{i}')) \n                tab.append(get_tab(test.loc[test.Patient == p, :])) \n        if len(x) <= 1:\n            continue\n        tab = np.array(tab) \n\n        x = np.expand_dims(x, axis=-1) \n        _a = model.predict([x, tab]) \n        a = np.quantile(_a, q)\n        A_test[p] = a\n        B_test[p] = test.FVC.values[test.Patient == p] - a*test.Weeks.values[test.Patient == p]\n        P_test[p] = test.Percent.values[test.Patient == p] \n        WEEK[p] = test.Weeks.values[test.Patient == p]\n\n    for k in sub.Patient_Week.values:\n        p, w = k.split('_')\n        w = int(w) \n\n        fvc = A_test[p] * w + B_test[p]\n        sub.loc[sub.Patient_Week == k, 'FVC'] = fvc\n        sub.loc[sub.Patient_Week == k, 'Confidence'] = (\n            P_test[p] - A_test[p] * abs(WEEK[p] - w) \n    ) \n\n    _sub = sub[[\"Patient_Week\",\"FVC\",\"Confidence\"]].copy()\n    subs.append(_sub)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"mid_df","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Averaging Predictions"},{"metadata":{"trusted":true},"cell_type":"code","source":"N = len(subs)\nsub = subs[0].copy() # ref\nsub[\"FVC\"] = 0\nsub[\"Confidence\"] = 0\nfor i in range(N):\n    sub[\"FVC\"] += subs[0][\"FVC\"] * (1/N)\n    sub[\"Confidence\"] += subs[0][\"Confidence\"] * (1/N)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub[[\"Patient_Week\",\"FVC\",\"Confidence\"]].to_csv(\"submission_img.csv\", index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img_sub = sub[[\"Patient_Week\",\"FVC\",\"Confidence\"]].copy()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Osic-Multiple-Quantile-Regression"},{"metadata":{"trusted":true},"cell_type":"code","source":"ROOT = \"../input/osic-pulmonary-fibrosis-progression\"\nBATCH_SIZE=128\n\ntr = pd.read_csv(f\"{ROOT}/train.csv\")\ntr.drop_duplicates(keep=False, inplace=True, subset=['Patient','Weeks'])\nchunk = pd.read_csv(f\"{ROOT}/test.csv\")\n\nprint(\"add infos\")\nsub = pd.read_csv(f\"{ROOT}/sample_submission.csv\")\nsub['Patient'] = sub['Patient_Week'].apply(lambda x:x.split('_')[0])\nsub['Weeks'] = sub['Patient_Week'].apply(lambda x: int(x.split('_')[-1]))\nsub =  sub[['Patient','Weeks','Confidence','Patient_Week']]\nsub = sub.merge(chunk.drop('Weeks', axis=1), on=\"Patient\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_root_path = '../input/osic-pulmonary-fibrosis-progression/train/'\nPatients_id = os.listdir(dicom_root_path)\nn_dicom_dict = {\"Patient\":[],\"n_dicom\":[],\"list_dicom\":[]}\n\nfor Patient_id in Patients_id:\n    dicom_id_path = glob.glob(dicom_root_path + Patient_id + \"/*\")\n    n_dicom_dict[\"n_dicom\"].append(len(dicom_id_path))\n    n_dicom_dict[\"Patient\"].append(Patient_id)\n    list_dicom_id = sorted([int(i.split(\"/\")[-1][:-4]) for i in dicom_id_path])\n    n_dicom_dict[\"list_dicom\"].append(list_dicom_id)\n\ndicom_pd = pd.DataFrame(n_dicom_dict)\ndicom_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(f\"min dicom number is {min(dicom_pd['n_dicom'])}\\n\\\nmax dicom number is {max(dicom_pd['n_dicom'])}\")\n\nplt.hist(dicom_pd['n_dicom'], bins=20)\nplt.title('Number of dicom per patient');","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd['height'],dicom_pd['width'], dicom_pd['kvp'] = -1,-1, -1\nfor Patient_id in Patients_id:\n    dicom_id_path = glob.glob(dicom_root_path + Patient_id + \"/*\")\n    for patient_dicom_id_path in dicom_id_path:\n        dicom = pydicom.dcmread(patient_dicom_id_path)\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'height'] = dicom.Rows\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'width'] = dicom.Columns\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'kvp'] = dicom.KVP\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'PatientPosition'] = dicom.PatientPosition\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'Manufacture'] = dicom.Manufacturer\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImageType0'] = dicom.ImageType[0] if len(dicom.ImageType) >= 1 else np.nan\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImageType1'] = dicom.ImageType[1] if len(dicom.ImageType) >= 2 else np.nan\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImageType2'] = dicom.ImageType[2] if len(dicom.ImageType) >=3 else np.nan\n        dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImageType3'] = dicom.ImageType[3] if len(dicom.ImageType) >=4 else np.nan\n#         print(dicom.ImageType)\n#         dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImagePositionPatientX'] = dicom.PatientPosition[0] if type(dicom.PatientPosition) is tuple else np.nan\n#         dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImagePositionPatientY'] = dicom.PatientPosition[1] if type(dicom.PatientPosition) is tuple else np.nan\n#         dicom_pd.loc[dicom_pd.Patient==Patient_id,'ImagePositionPatientZ'] = dicom.PatientPosition[2] if type(dicom.PatientPosition) is tuple else np.nan\n        break","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.hist(dicom_pd['PatientPosition'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reshape_dicom_pd = dicom_pd.loc[(dicom_pd.height!=512) | (dicom_pd.width!=512),:]\nreshape_dicom_pd = reshape_dicom_pd.reset_index(drop=True)\nreshape_dicom_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, ax = plt.subplots(len(reshape_dicom_pd.head()),2, figsize=(15, 18))\nfor idx,patient_id in enumerate(reshape_dicom_pd.head()['Patient']):\n    paths = random.sample(glob.glob(dicom_root_path + patient_id + \"/*\"),2)\n    dicom1 = pydicom.dcmread(paths[0])\n    dicom2 = pydicom.dcmread(paths[1])\n    ax[idx,0].set_title(f\"{patient_id}-{paths[0].split('/')[-1][:-4]}-{reshape_dicom_pd.loc[idx,'height']}-{reshape_dicom_pd.loc[idx,'width']}\")\n    ax[idx,0].imshow(dicom1.pixel_array, cmap=plt.cm.bone)\n    ax[idx,1].set_title(f\"patient id is {patient_id}-{paths[1].split('/')[-1][:-4]}-{reshape_dicom_pd.loc[idx,'height']}-{reshape_dicom_pd.loc[idx,'width']}\")\n    ax[idx,1].imshow(dicom2.pixel_array, cmap=plt.cm.bone)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"crop_id = ['ID00240637202264138860065','ID00122637202216437668965','ID00086637202203494931510',\n            'ID00419637202311204720264','ID00014637202177757139317','ID00094637202205333947361',\n            'ID00067637202189903532242',]\nreshape_dicom_pd['resize_type'] = 'resize'\nreshape_dicom_pd.loc[reshape_dicom_pd.Patient.isin(crop_id),'resize_type'] = 'crop'","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd['resize_type'] = 'no'\nfor idx,i in enumerate(reshape_dicom_pd['Patient']):\n    dicom_pd.loc[dicom_pd.Patient==i,'resize_type'] = reshape_dicom_pd.loc[idx,'resize_type']\ndicom_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_pd = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/train.csv')\ntemp_pd = pd.DataFrame(columns=train_pd.columns)\nfor i in range(len(dicom_pd)):\n    patient_pd = train_pd[train_pd.Patient==dicom_pd.iloc[i].Patient]\n    zeroweek = patient_pd['Weeks'].min()\n    #if sum(patient_pd.Weeks==zeroweek)>1:\n    #    print(pd.unique(patient_pd.Patient))\n    temp_pd = temp_pd.append(patient_pd[patient_pd.Weeks==zeroweek].iloc[0])\ndicom_pd = pd.merge(dicom_pd, temp_pd, on=['Patient'])\ndicom_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd[dicom_pd.resize_type!='no'].head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd['ImageType0'] = dicom_pd['ImageType0'].apply(lambda x: 'ORIGINAL' if x not in ['ORIGINAL', 'DERIVED'] else x)\ndicom_pd['ImageType1'] = dicom_pd['ImageType1'].apply(lambda x: 'PRIMARY' if x not in ['PRIMARY', 'SECONDARY'] else x)\ndicom_pd['ImageType2'] = dicom_pd['ImageType2'].apply(lambda x: 'AXIAL' if x not in ['AXIAL', 'REFORMATTED', 'OTHER'] else x)\ndicom_pd[dicom_pd['Patient']=='ID00421637202311550012437']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import seaborn as sns\ndicom_pd.isna().sum()\n# dicom_pd.shape\nsns.countplot(x ='Manufacture', hue = 'ImageType3', data = dicom_pd)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"d = dicom_pd[dicom_pd['Manufacture'].isin(['TOSHIBA', 'GE MEDICAL SYSTEMS'])]\nsns.countplot(x ='Manufacture', hue = 'ImageType3', data = d)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd_mf = dicom_pd[dicom_pd['ImageType3'].isna()][['Manufacture', 'ImageType3']]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# dicom_pd.shape\nsns.countplot(x ='Manufacture',  data = dicom_pd_mf)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# dicom_pd_mf.loc[(dicom_pd_mf['Manufacture']=='SIEMENS') & (dicom_pd['ImageType3'].isna()==True), 'ImageType3'] = 'CT_SOM5_SPI'\n# dicom_pd_mf.loc[(dicom_pd_mf['Manufacture']=='Philips') & (dicom_pd['ImageType3'].isna()==True), 'ImageType3'] = 'HELIX'\n# dicom_pd_mf.loc[(dicom_pd_mf['Manufacture']=='GE MEDICAL SYSTEMS') & (dicom_pd['ImageType3'].isna()==True), 'ImageType3'] = 'AVERAGE'\n\ndicom_pd_mf.loc[dicom_pd_mf['Manufacture']=='GE MEDICAL SYSTEMS']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_pd","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def load_scan(path,resize_type='no'):\n    \"\"\"\n    Loads scans from a folder and into a list.\n    \n    Parameters: path (Folder path)\n    \n    Returns: slices (List of slices)\n    \"\"\"\n    slices = [pydicom.read_file(path + '/' + s) for s in os.listdir(path)]\n    slices.sort(key = lambda x: int(x.InstanceNumber))\n    \n    try:\n        slice_thickness = abs(slices[-1].ImagePositionPatient[2] - slices[0].ImagePositionPatient[2])/(len(slices))\n    except:\n        try:\n            slice_thickness = abs(slices[-1].SliceLocation - slices[0].SliceLocation)/(len(slices))\n        except:\n            slice_thickness = slices[0].SliceThickness\n        \n    for s in slices:\n        s.SliceThickness = slice_thickness\n        if resize_type == 'resize':\n            s.PixelSpacing = s.PixelSpacing*(s.Rows/512)  \n    return slices","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def transform_to_hu(slices):\n    \"\"\"\n    transform dicom.pixel_array to Hounsfield.\n    Parameters: list dicoms\n    Returns:numpy Hounsfield\n    \"\"\"\n    \n    images = np.stack([file.pixel_array for file in slices])\n    images = images.astype(np.int16)\n\n    # convert ouside pixel-values to air:\n    # I'm using <= -1000 to be sure that other defaults are captured as well\n    #images[images <= -1000] = 0\n    \n    # convert to HU\n    for n in range(len(slices)):\n        \n        intercept = slices[n].RescaleIntercept\n        slope = slices[n].RescaleSlope\n        \n        if slope != 1:\n            images[n] = slope * images[n].astype(np.float64)\n            images[n] = images[n].astype(np.int16)\n            \n        images[n] += np.int16(intercept)\n    \n    return np.array(images, dtype=np.int16)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def crop_image(img: np.ndarray):\n    edge_pixel_value = img[0, 0]\n    mask = img != edge_pixel_value\n    return img[np.ix_(mask.any(1),mask.any(0))]\n\ndef resize_image(img: np.ndarray,reshape=(512,512)):\n    img = cv2.resize(img,(512,512))\n    return img\n\ndef preprocess_img(img,resize_type):\n    if resize_type == 'resize':\n        img = [resize_image(im) for im in img]\n    if resize_type == 'crop':\n        img = [crop_image(im) for im in img]\n        \n    return np.array(img, dtype=np.int64)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class Test_Generate(Dataset):\n    def __init__(self,imgs_dicom,resize_type='no'):\n        self.imgs_dicom = imgs_dicom\n        self.resize_type = resize_type\n    def __getitem__(self,index):\n        \n        slice_img = self.imgs_dicom[index].pixel_array\n        slice_img = (slice_img-slice_img.min())/(slice_img.max()-slice_img.min())\n        slice_img = (slice_img*255).astype(np.uint8)\n        if self.resize_type == 'crop':\n            slice_img = crop_image(slice_img)\n        elif self.resize_type == 'resize':\n            slice_img = cv2.resize(slice_img,(512,512))\n            \n        slice_img = slice_img[None,:,:]\n        slice_img = (slice_img/255).astype(np.float32)\n        return slice_img\n        \n    def __len__(self):\n        return len(self.imgs_dicom)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"device =  torch.device('cuda:0')\nmodel = smp.Unet('densenet121', classes=1, in_channels=1,activation='sigmoid',encoder_weights=None).to(device)\nmodel.load_state_dict(torch.load('../input/2020osic/best_lung_Unet_densenet121.pth'))\nbatch = 8\n\ndef Unet_mask(model: nn.Module,input_data: DataLoader):\n    model.eval()\n    outs = []\n    for idx, sample in enumerate(test_loader):\n        image = sample\n        image = image.to(device)\n        with torch.no_grad():\n            out = model(image)\n        out = out.cpu().data.numpy()\n        out = np.where(out>0.5,1,0)\n        out = np.squeeze(out,axis=1)\n        outs.append(out)\n\n    outs = np.concatenate(outs)\n    return outs","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, ax = plt.subplots(4,2, figsize=(14, 14))\n\nfor i in range(4):\n    path = os.path.join(dicom_root_path,dicom_pd.iloc[i].Patient)\n    patient_scans = load_scan(path)\n    \n    test_db = Test_Generate(patient_scans,dicom_pd.iloc[i].resize_type)\n    test_loader = DataLoader(test_db, batch_size=batch, shuffle=False, num_workers=4)\n    \n    masks = Unet_mask(model,test_loader)\n    \n    #patient_images = transform_to_hu(patient_scans)\n    \n    #patient_images = preprocess_img(patient_images,dicom_pd.iloc[i])\n    \n    num_slices = len(masks)\n    patient_image = test_db[num_slices//2][0]\n    patient_mask = masks[num_slices//2]\n    \n    #Mask = generate_internal_mask(patient_image)\n    \n    ax[i,0].set_title(f\"{dicom_pd.iloc[i].Patient}-{dicom_pd.iloc[i].FVC}\")\n    ax[i,0].imshow(patient_image,cmap='gray')\n    ax[i,1].imshow(patient_mask)\n    \nplt.show()\nplt.close()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"thresh = [-1000,0]\nf, ax = plt.subplots(2,2, figsize=(18, 18))\nsampler = random.sample(range(len(dicom_pd)),4)\nfor i in range(4):\n    path = os.path.join(dicom_root_path,dicom_pd.iloc[sampler[i]].Patient)\n    patient_scans = load_scan(path)\n    \n    test_db = Test_Generate(patient_scans,dicom_pd.iloc[sampler[i]].resize_type)\n    test_loader = DataLoader(test_db, batch_size=batch, shuffle=False, num_workers=4)\n    \n    masks = Unet_mask(model,test_loader)\n    \n    patient_images = transform_to_hu(patient_scans)\n    patient_images = preprocess_img(patient_images,dicom_pd.loc[sampler[i],'resize_type'])\n    \n    num_slices = len(patient_images)\n    #patient_images = patient_images[int(num_slices*0.1):int(num_slices*0.9)]\n    #patient_masks = masks[int(num_slices*0.1):int(num_slices*0.9)]\n    \n    #patient_images = patient_images[num_slices//2]\n    #patient_masks = pool.map(generate_internal_mask,patient_images)\n    #patient_masks = patient_masks[int(num_slices*0.1):int(num_slices*0.9)]\n    patient_images = masks*patient_images\n    patient_images_nonzero = patient_images[np.nonzero(patient_images)]\n    #patient_images_mean = np.mean(patient_images,0)\n    \n    s_pixel = patient_images_nonzero.flatten()\n    s_pixel = s_pixel[np.where((s_pixel>thresh[0])&(s_pixel<thresh[1]))]\n    \n    ax[i//2,i%2].set_title(f\"{dicom_pd.iloc[i].Patient}-{dicom_pd.iloc[i].FVC}\")\n    ax[i//2,i%2].hist(s_pixel, bins=20)\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#def func_volume(patient_scan,patient_mask):\n    \n\ndef caculate_lung_volume(patient_scans,patient_masks):\n    \"\"\"\n    caculate volume of lung from mask\n    Parameters: list dicom scans,list patient CT Mask\n    Returns: volume cm³　(float)\n    \"\"\"\n    lung_volume = 0\n    for i in range(len(patient_masks)):\n        \n        pixel_spacing = patient_scans[i].PixelSpacing\n        slice_thickness = patient_scans[i].SliceThickness\n        lung_volume += np.count_nonzero(patient_masks[i])*pixel_spacing[0]*pixel_spacing[1]*slice_thickness\n        \n    return lung_volume*0.001","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def caculate_histgram_statistical(patient_images,patient_masks,thresh = [-600,0]):\n    \"\"\"\n    caculate hisgram kurthosis of lung hounsfield\n    Parameters: list patient CT image 512*512,thresh divide lung\n    Returns: histgram statistical characteristic(Mean,Skew,Kurthosis)\n    \"\"\"\n    statistical_characteristic = dict(Mean=0,Skew=0,Kurthosis=0)\n    num_slices = len(patient_images)\n    \n    #patient_images = patient_images[int(num_slices*0.1):int(num_slices*0.9)]\n    #patient_masks = patient_masks[int(num_slices*0.1):int(num_slices*0.9)]\n    patient_images = patient_masks*patient_images\n    patient_images_nonzero = patient_images[np.nonzero(patient_images)]\n    \n    s_pixel = patient_images_nonzero.flatten()\n    s_pixel = s_pixel[np.where((s_pixel>thresh[0])&(s_pixel<thresh[1]))]\n    \n    statistical_characteristic['Mean'] = np.mean(s_pixel)\n    statistical_characteristic['Skew'] = skew(s_pixel)\n    statistical_characteristic['Kurthosis'] = kurtosis(s_pixel)\n    \n    return statistical_characteristic","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import gdcm\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"lung_stat_pd = pd.DataFrame(columns=['Patient','Volume','Mean','Skew','Kurthosis'])\n\nfor i in tqdm(range(len(dicom_pd))):\n    path = os.path.join(dicom_root_path,dicom_pd.iloc[i].Patient)\n    lung_stat_pd.loc[i,'Patient'] = dicom_pd.iloc[i].Patient\n    patient_scans = load_scan(path)\n    \n    test_db = Test_Generate(patient_scans,dicom_pd.iloc[i].resize_type)\n    test_loader = DataLoader(test_db, batch_size=batch, shuffle=False, num_workers=4)\n    masks = Unet_mask(model,test_loader)\n    \n    \n    patient_images = transform_to_hu(patient_scans)\n    patient_images = preprocess_img(patient_images,dicom_pd.loc[i,'resize_type'])\n    \n    lung_stat_pd.loc[i,'Volume'] = caculate_lung_volume(patient_scans,masks)                           \n    #patient_images = resize_image(patient_images) if dicom_pd.iloc[i].resize_type=='resize' else patient_images\n    #patient_images = resize_image(patient_masks) if dicom_pd.iloc[i].resize_type=='resize' else patient_images\n    \n    statistical_characteristic = caculate_histgram_statistical(patient_images,masks,thresh)\n    lung_stat_pd.loc[i,'Mean'] = statistical_characteristic['Mean']\n    lung_stat_pd.loc[i,'Skew'] = statistical_characteristic['Skew']\n    lung_stat_pd.loc[i,'Kurthosis'] = statistical_characteristic['Kurthosis']\n    \nlung_stat_pd.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_feature = pd.merge(dicom_pd, lung_stat_pd, on=['Patient'])\ndicom_feature.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_feature = dicom_feature.drop(['list_dicom', 'height','width','resize_type','n_dicom'], axis=1)\ndicom_feature = dicom_feature.drop(['ImageType3'], axis=1)\ndicom_feature.head(20)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_feature","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dicom_feature.to_csv('./CT_feature.csv',index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"ROOT = \"../input/osic-pulmonary-fibrosis-progression\"\ndevice = torch.device('cuda')\n\ntest_df = pd.read_csv(f\"{ROOT}/test.csv\")\nsub = pd.read_csv(f\"{ROOT}/sample_submission.csv\")\nsub['Patient'] = sub['Patient_Week'].apply(lambda x:x.split('_')[0])\nsub['Weeks'] = sub['Patient_Week'].apply(lambda x: int(x.split('_')[-1]))\nsub =  sub[['Patient','Weeks','Confidence','Patient_Week']]\nsub = sub.merge(test_df.drop('Weeks', axis=1), on=\"Patient\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class Inference_Generate(Dataset):\n    def __init__(self,imgs_dicom):\n        self.imgs_dicom = imgs_dicom\n        \n    def __getitem__(self,index):\n        metainf = self.imgs_dicom[index]\n        slice_img = metainf.pixel_array\n        slice_img = (slice_img-slice_img.min())/(slice_img.max()-slice_img.min())\n        slice_img = (slice_img*255).astype(np.uint8)\n        if metainf.Rows!=512 or metainf.Columns!=512:\n            slice_img = cv2.resize(slice_img,(512,512))\n            \n        slice_img = slice_img[None,:,:]\n        slice_img = (slice_img/255).astype(np.float32)\n        return slice_img\n        \n    def __len__(self):\n        return len(self.imgs_dicom)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"thresh = [-1000,0]\nct_root_path = '../input/osic-pulmonary-fibrosis-progression/test/'\nlung_test_feature = pd.DataFrame(columns=['Patient','Volume','Mean','Skew','Kurthosis', 'kvp', 'PatientPosition','Manufacture', 'ImageType0',\n                                          'ImageType1', 'ImageType2'])\nfor idx,i in enumerate(pd.unique(test_df['Patient'])):\n    lung_test_feature.loc[idx,'Patient'] = i\n    patient_scans = load_scan(ct_root_path + i)\n    dicom_id_path = glob.glob(ct_root_path + i + \"/*\")\n    for patient_dicom_id_path in dicom_id_path:\n        sc = pydicom.dcmread(patient_dicom_id_path)\n        lung_test_feature.loc[lung_test_feature.Patient==i,'kvp'] = sc.KVP\n        lung_test_feature.loc[lung_test_feature.Patient==i,'PatientPosition'] = sc.PatientPosition\n        lung_test_feature.loc[lung_test_feature.Patient==i,'Manufacture'] = sc.Manufacturer\n        lung_test_feature.loc[lung_test_feature.Patient==i,'ImageType0'] = sc.ImageType[0] if len(sc.ImageType) >= 1 else np.nan\n        lung_test_feature.loc[lung_test_feature.Patient==i,'ImageType1'] = sc.ImageType[1] if len(sc.ImageType) >= 2 else np.nan\n        lung_test_feature.loc[lung_test_feature.Patient==i,'ImageType2'] = sc.ImageType[2] if len(sc.ImageType) >=3 else np.nan\n    test_db = Inference_Generate(patient_scans)\n    test_loader = DataLoader(test_db, batch_size=batch, shuffle=False, num_workers=4)\n    masks = Unet_mask(model,test_loader)\n    \n    patient_images = transform_to_hu(patient_scans)\n  \n    if patient_images[0].shape!=(512,512):\n        patient_images = preprocess_img(patient_images,'resize')\n  \n    lung_test_feature.loc[idx,'Volume'] = caculate_lung_volume(patient_scans,masks)     \n    statistical_characteristic = caculate_histgram_statistical(patient_images,masks,thresh)\n    lung_test_feature.loc[idx,'Mean'] = statistical_characteristic['Mean']\n    lung_test_feature.loc[idx,'Skew'] = statistical_characteristic['Skew']\n    lung_test_feature.loc[idx,'Kurthosis'] = statistical_characteristic['Kurthosis']\nlung_test_feature.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"lung_test_feature['ImageType0'] = lung_test_feature['ImageType0'].apply(lambda x: 'ORIGINAL' if x not in ['ORIGINAL', 'DERIVED'] else x)\nlung_test_feature['ImageType1'] = lung_test_feature['ImageType1'].apply(lambda x: 'PRIMARY' if x not in ['PRIMARY', 'SECONDARY'] else x)\nlung_test_feature['ImageType2'] = lung_test_feature['ImageType2'].apply(lambda x: 'AXIAL' if x not in ['AXIAL', 'REFORMATTED', 'OTHER'] else x)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_df","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_df = test_df.merge(lung_test_feature, on='Patient')\nsub = sub.merge(lung_test_feature, on='Patient')\n\ntrain_df = pd.read_csv(f\"{ROOT}/train.csv\")\n#train_df = train_df.merge(lung_stat_pd, on='Patient')\nfeature_ct = pd.read_csv('CT_feature.csv',usecols=['Patient', 'kvp', 'PatientPosition', 'Manufacture', 'ImageType0',\n       'ImageType1', 'ImageType2', 'Volume', 'Mean', 'Skew', 'Kurthosis'])\ntrain_df = train_df.merge(feature_ct, on='Patient')\n\ntrain_df['WHERE'] = 'train'\ntest_df['WHERE'] = 'val'\nsub['WHERE'] = 'test'\ndata = train_df.append([test_df, sub])\n\ndata.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data = pd.read_csv('../input/final-data/final_data.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"list(data.columns)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data['min_week'] = data['Weeks']\ndata.loc[data.WHERE=='test','min_week'] = np.nan\ndata['min_week'] = data.groupby('Patient')['min_week'].transform('min')\n\nbase = (data\n    .loc[data.Weeks == data.min_week][['Patient','FVC']]\n    .rename({'FVC': 'min_FVC'}, axis=1)\n    .groupby('Patient')\n    .first()\n    .reset_index())\n\ndata = data.merge(base, on='Patient', how='left')\ndata['base_week'] = data['Weeks'] - data['min_week']\ndel base\n\nFE = list(data.Sex.unique()) + list(data.SmokingStatus.unique()) + list(data.PatientPosition.unique()) + list(data.Manufacture.unique()) + list(data.ImageType0.unique()) + list(data.ImageType1.unique()) + list(data.ImageType2.unique())\ndata = pd.concat([\n    data,\n    pd.get_dummies(data.Sex),\n    pd.get_dummies(data.SmokingStatus),\n    pd.get_dummies(data.PatientPosition),\n    pd.get_dummies(data.Manufacture),\n    pd.get_dummies(data.ImageType0),\n    pd.get_dummies(data.ImageType1),\n    pd.get_dummies(data.ImageType2),\n    pd.get_dummies(data.kvp),\n], axis=1)\n\ndata['age'] = (data['Age'] - data['Age'].min() ) / ( data['Age'].max() - data['Age'].min() )\ndata['BASE'] = (data['min_FVC'] - data['min_FVC'].min() ) / ( data['min_FVC'].max() - data['min_FVC'].min())\ndata['week'] = (data['base_week'] - data['base_week'].min() ) / ( data['base_week'].max() - data['base_week'].min() )\ndata['percent'] = (data['Percent'] - data['Percent'].min() ) / ( data['Percent'].max() - data['Percent'].min())\n\ndata['volume'] = (data['Volume'] - data['Volume'].min() ) / ( data['Volume'].max() - data['Volume'].min())\ndata['mean'] = (data['Mean'] - data['Mean'].min()) / (data['Mean'].max() - data['Mean'].min())\ndata['skew'] = (data['Skew'] - data['Skew'].min())/(data['Skew'].max() - data['Skew'].min())\ndata['kurthosis'] = (data['Kurthosis'] - data['Kurthosis'].min())/(data['Kurthosis'].max() - data['Kurthosis'].min())\n\nFE += ['age','percent','week', 'BASE','volume','mean','skew','kurthosis']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data.rename(columns={100: 'kvp_100', 110: 'kvp_110', 120: 'kvp_120', 130: 'kvp_130', 135: 'kvp_135', 140: 'kvp_140'}, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df = data.loc[data.WHERE=='train']\ntest_df = data.loc[data.WHERE=='val']\nsub = data.loc[data.WHERE=='test']\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"FE = ['Male',\n 'Female',\n 'Ex-smoker',\n 'Never smoked',\n 'Currently smokes',\n#  'FFS',\n#  'HFS',\n#  'FFP',\n#  'HFP',\n#  'GE MEDICAL SYSTEMS',\n#  'SIEMENS',\n#  'TOSHIBA',\n#  'Philips',\n#  'PACSMATT',\n#  'Hitachi Medical Corporation',\n#  'PACSGEAR',\n#  'ORIGINAL',\n#  'DERIVED',\n#  'PRIMARY',\n#  'SECONDARY',\n#  'AXIAL',\n#  'REFORMATTED',\n#  'OTHER',\n 'age',\n 'percent',\n 'week',\n 'BASE',\n 'volume',\n 'mean',\n 'skew',\n 'kurthosis',\n#  'kvp_100',\n#  'kvp_110',\n#  'kvp_120',\n#  'kvp_130',\n#  'kvp_135',\n#  'kvp_140'\n     ]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data.to_csv('./final_data.csv',index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# data = pd.read_csv('../input/final-data/final_data.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"tr = data.loc[data.WHERE=='train']\nchunk = data.loc[data.WHERE=='val']\nsub = data.loc[data.WHERE=='test']\n# del data","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"tr.shape, chunk.shape, sub.shape","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# The change of mloss\n\nHere is where I tuned the `mloss` from 0.8 to 0.65. You can try a grid-search to maybe find an optimal value - however, I have only tried a very few choices like 0.65, 0.7 and 0.75."},{"metadata":{"trusted":true},"cell_type":"code","source":"C1, C2 = tf.constant(70, dtype='float32'), tf.constant(1000, dtype=\"float32\")\n\ndef score(y_true, y_pred):\n    tf.dtypes.cast(y_true, tf.float32)\n    tf.dtypes.cast(y_pred, tf.float32)\n    sigma = y_pred[:, 2] - y_pred[:, 0]\n    fvc_pred = y_pred[:, 1]\n#     print(\"fvc_pred\", fvc_pred, 'fvc_true',y_true[:, 0], y_true, y_pred)\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 qloss(y_true, y_pred):\n    # Pinball loss for multiple quantiles\n    qs = [0.2, 0.50, 0.8]\n    q = tf.constant(np.array([qs]), dtype=tf.float32)\n    e = y_true - y_pred\n    v = tf.maximum(q*e, (q-1)*e)\n    return K.mean(v)\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\n# def make_model(nh):\n#     z = L.Input((nh,), name=\"Patient\")\n#     x = L.Dense(100, activation=\"relu\", name=\"d1\")(z)\n#     x = L.Dense(100, activation=\"relu\", name=\"d2\")(x)\n#     p1 = L.Dense(3, activation=\"linear\", name=\"p1\")(x)\n#     p2 = L.Dense(3, activation=\"relu\", name=\"p2\")(x)\n#     preds = L.Lambda(lambda x: x[0] + tf.cumsum(x[1], axis=1), \n#                      name=\"preds\")([p1, p2])\n    \n#     model = M.Model(z, preds, name=\"CNN\")\n#     model.compile(loss=mloss(0.65), optimizer=tf.keras.optimizers.Adam(lr=0.1, beta_1=0.9, beta_2=0.999, epsilon=None, decay=0.01, amsgrad=False), metrics=[score])\n#     return model\n\ndef make_model(nh):\n    z = L.Input((nh,), name=\"Patient\")\n    y = L.LSTM(25, return_sequences=True)(tf.expand_dims(z, axis=-1))\n    y = L.LSTM(25, return_sequences=True)(y)\n    y = L.Reshape([25*nh])(y)\n    x = L.Dense(100, activation=\"relu\", name=\"d1\")(y)\n    x = L.Dense(100, activation=\"relu\", name=\"d2\")(x)\n    x = L.Dense(100, activation=\"relu\", name=\"d3\")(x)\n    p1 = L.Dense(3, activation=\"linear\", name=\"p1\")(x)\n    p2 = L.Dense(3, activation=\"relu\", name=\"p2\")(x)\n    preds = L.Lambda(lambda x: x[0] + tf.cumsum(x[1], axis=1), \n                     name=\"preds\")([p1, p2])\n    \n    model = M.Model(z, preds, name=\"CNN\")\n    #model.compile(loss=qloss, optimizer=\"adam\", metrics=[score])\n    model.compile(loss=mloss(0.65), optimizer=tf.keras.optimizers.Adam(lr=0.1, beta_1=0.9, beta_2=0.999, epsilon=None, decay=0.01, amsgrad=False), metrics=[score])\n    return model","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y = tr['FVC'].values\nz=tr[FE].values\nze=sub[FE].values\nnh = z.shape[1]\npe = np.zeros((ze.shape[0], 3))\npred = np.zeros((z.shape[0], 3))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"net = make_model(nh)\nprint(net.summary())\nprint(net.count_params())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"NFOLD = 15 # originally 5\nkf = KFold(n_splits=NFOLD)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\ncnt = 0\nEPOCHS = 800\nBATCH_SIZE = 128\nfor tr_idx, val_idx in kf.split(z):\n    cnt += 1\n    print(f\"FOLD {cnt}\")\n    net = make_model(nh)\n    print(\"working\")\n    net.fit(z[tr_idx], y[tr_idx], batch_size=BATCH_SIZE, epochs=EPOCHS, \n            validation_data=(z[val_idx], y[val_idx]), verbose=0) #\n    print(\"train\", net.evaluate(z[tr_idx], y[tr_idx], verbose=0, batch_size=BATCH_SIZE))\n    print(\"val\", net.evaluate(z[val_idx], y[val_idx], verbose=0, batch_size=BATCH_SIZE))\n    print(\"predict val...\")\n    pred[val_idx] = net.predict(z[val_idx], batch_size=BATCH_SIZE, verbose=0)\n    print(y[val_idx], pred[val_idx])\n    print(\"predict test...\")\n    pe += net.predict(ze, batch_size=BATCH_SIZE, verbose=0) / NFOLD","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sigma_opt = mean_absolute_error(y, pred[:, 1])\nunc = pred[:,2] - pred[:, 0]\nsigma_mean = np.mean(unc)\nprint(sigma_opt, sigma_mean)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"idxs = np.random.randint(0, y.shape[0], 100)\nplt.plot(y[idxs], label=\"ground truth\")\nplt.plot(pred[idxs, 0], label=\"q25\")\nplt.plot(pred[idxs, 1], label=\"q50\")\nplt.plot(pred[idxs, 2], label=\"q75\")\nplt.legend(loc=\"best\")\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(unc.min(), unc.mean(), unc.max(), (unc>=0).mean())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.hist(unc)\nplt.title(\"uncertainty in prediction\")\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# PREDICTION\nsub['FVC1'] = 1.*pe[:, 1]\nsub['Confidence1'] = pe[:, 2] - pe[:, 0]\nsubm = sub[['Patient_Week','FVC','Confidence','FVC1','Confidence1']].copy()\nsubm.loc[~subm.FVC1.isnull()].head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm.loc[~subm.FVC1.isnull(),'FVC'] = subm.loc[~subm.FVC1.isnull(),'FVC1']\nif sigma_mean<70:\n    subm['Confidence'] = sigma_opt\nelse:\n    subm.loc[~subm.FVC1.isnull(),'Confidence'] = subm.loc[~subm.FVC1.isnull(),'Confidence1']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm.describe().T","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"otest = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv')\nfor i in range(len(otest)):\n    subm.loc[subm['Patient_Week']==otest.Patient[i]+'_'+str(otest.Weeks[i]), 'FVC'] = otest.FVC[i]\n    subm.loc[subm['Patient_Week']==otest.Patient[i]+'_'+str(otest.Weeks[i]), 'Confidence'] = 0.1","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm[[\"Patient_Week\",\"FVC\",\"Confidence\"]].to_csv(\"submission_regression_0.65_only_pos.csv\", index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reg_sub = subm[[\"Patient_Week\",\"FVC\",\"Confidence\"]].copy()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Ensemble (Simple Blend)"},{"metadata":{"trusted":true},"cell_type":"code","source":"img_sub = pd.read_csv('../input/img-sub/Submitted_op.csv')\ndf1 = img_sub.sort_values(by=['Patient_Week'], ascending=True).reset_index(drop=True)\ndf2 = reg_sub.sort_values(by=['Patient_Week'], ascending=True).reset_index(drop=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df = df1[['Patient_Week']].copy()\ndf['FVC'] = 0.25*df1['FVC'] + 0.75*df2['FVC']\ndf['Confidence'] = 0.26*df1['Confidence'] + 0.74*df2['Confidence']\ndf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df = pd.read_csv('../input/output/submission.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df.to_csv('./sample_submission.csv', 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}