{"cells":[{"metadata":{},"cell_type":"markdown","source":"# OSIC Pulmonary Fibrosis Progression\n## Predict lung function decline\nThe aim of this competition is to predict a patient’s severity of decline in lung function based on a CT scan of their lungs. Lung function is assessed based on output from a spirometer, which measures the forced vital capacity (FVC), i.e. the volume of air exhaled.\n\nIn the dataset, you are provided with a baseline chest CT scan and associated clinical information for a set of patients. A patient has an image acquired at time Week = 0 and has numerous follow up visits over the course of approximately 1-2 years, at which time their FVC is measured.\n\nIn the training set, you are provided with an anonymized, baseline CT scan and the entire history of FVC measurements.\nIn the test set, you are provided with a baseline CT scan and only the initial FVC measurement. You are asked to predict the final three FVC measurements for each patient, as well as a confidence value in your prediction.\nThere are around 200 cases in the public & private test sets, combined. This is split roughly 15-85 between public-private.\n\nSince this is real medical data, you will notice the relative timing of FVC measurements varies widely. The timing of the initial measurement relative to the CT scan and the duration to the forecasted time points may be different for each patient. This is considered part of the challenge of the competition. To avoid potential leakage in the timing of follow up visits, you are asked to predict every patient's FVC measurement for every possible week. Those weeks which are not in the final three visits are ignored in scoring.\n\n#### Files\nThis is a synchronous rerun code competition. The provided test set is a small representative set of files (copied from the training set) to demonstrate the format of the private test set. When you submit your notebook, Kaggle will rerun your code on the test set, which contains unseen images.\n\n* train.csv - the training set, contains full history of clinical information\n* test.csv - the test set, contains only the baseline measurement\n* train/ - contains the training patients' baseline CT scan in DICOM format\n* test/ - contains the test patients' baseline CT scan in DICOM format\n* sample_submission.csv - demonstrates the submission format\n\n#### Columns\ntrain.csv and test.csv\n* Patient- a unique Id for each patient (also the name of the patient's DICOM folder)\n* Weeks- the relative number of weeks pre/post the baseline CT (may be negative)\n* FVC - the recorded lung capacity in ml\n* Percent- a computed field which approximates the patient's FVC as a percent of the typical FVC for a person of similar characteristics\n* Age\n* Sex\n* SmokingStatus\n\n#### sample submission.csv\n* Patient_Week - a unique Id formed by concatenating the Patient and Weeks columns (i.e. ABC_22 is a prediction for patient ABC at week 22)\n* FVC - the predicted FVC in ml\n* Confidence - a confidence value of your prediction (also has units of ml)"},{"metadata":{},"cell_type":"markdown","source":"## Data Importing"},{"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 typing as tp\nimport pydicom\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.model_selection import GridSearchCV\nfrom sklearn.metrics import make_scorer\n\nfrom xgboost import XGBRegressor\nfrom lightgbm import LGBMRegressor\nfrom sklearn.linear_model import HuberRegressor\n\nfrom statsmodels.formula.api import quantreg","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/train.csv')\ntest_df = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv')\ndf = pd.concat([train_df, test_df], ignore_index=True)\ndf['Patient_Week'] = df['Patient'].astype(str) + '_ '+ df['Weeks'].astype(str)\ndf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Shape of Training data: ', train_df.shape)\nprint('Shape of Test data: ', test_df.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Data Preprocessing"},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_height(data)->'dataframe':\n    data['Height'] = 0\n    data['Height'] = data.apply(lambda x: (x.FVC / (21.78 - (0.101 * x.Age))) if x.Sex == 1 else (x.FVC / (27.63 - (0.112 * x.Age))), axis=1)\n    \ndef add_norm(data)->'dataframe':\n    return (data - data.mean()) / data.std()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df['Sex'] = df['Sex'].map({'Female': 0, 'Male': 1})\ndf['SmokingStatus'] = df['SmokingStatus'].map({'Currently smokes': 0, 'Never smoked': 1, 'Ex-smoker': 2})\ndf = df.drop('Patient_Week', axis=1)\ndf = df.set_index('Patient')\nadd_height(df)\ndf[df.columns[~df.columns.isin(['FVC', 'Sex', 'Weeks', 'SmokingStatus'])]] = add_norm(df[df.columns[~df.columns.isin(['FVC', 'Sex', 'Weeks', 'SmokingStatus'])]])\ndf.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Submission part"},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/sample_submission.csv')\nsub_df.drop(['FVC'], axis=1, inplace=True)\nsub_df['Patient'] = sub_df['Patient_Week'].apply(lambda x: x.split('_')[0])\nsub_df['pred_Weeks'] = sub_df['Patient_Week'].apply(lambda x: x.split('_')[1]).astype(int)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_FVC = pd.merge(sub_df, df, how='left', on=['Patient'])\ntest_FVC = test_FVC[['Patient_Week', 'pred_Weeks', 'Percent', 'Age', 'Sex', 'SmokingStatus', 'Height']]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_FVC.rename(columns={'pred_Weeks':'Weeks'}, inplace=True)\ntest_FVC = test_FVC.groupby('Patient_Week').mean()\ntest_FVC.sort_values(['Weeks'], inplace=True)\ntest_FVC","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(test_FVC.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Train&Validation part"},{"metadata":{"trusted":true},"cell_type":"code","source":"X = df.iloc[:, df.columns != \"FVC\"]\ny = df[\"FVC\"]\n\nX_train = X[:-5]\ny_train = y[:-5]\nX_val = X[-5:]\ny_val = y[-5:]\n\nprint(X.shape)\nprint(y.shape)\nprint(X_train.shape)\nprint(y_train.shape)\nprint(X_val.shape)\nprint(y_val.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Visualization & EDA"},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(2, 2, figsize=(15, 10))\ntrain_df.boxplot('FVC', by='SmokingStatus', ax = axs[0, 0])\ntrain_df.boxplot('Percent', by='SmokingStatus', ax = axs[0, 1])\ntrain_df.boxplot('FVC', by='Sex', ax = axs[1, 0])\ntrain_df.boxplot('Percent', by='Sex', ax = axs[1, 1])\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img = \"../input/osic-pulmonary-fibrosis-progression/train/ID00009637202177434476278/100.dcm\"\nds = pydicom.dcmread(img)\nplt.figure(figsize = (7,7))\nplt.imshow(ds.pixel_array, cmap=plt.cm.bone)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img_1 = \"../input/osic-pulmonary-fibrosis-progression/train/ID00009637202177434476278/100.dcm\"\nimg_2 = \"../input/osic-pulmonary-fibrosis-progression/train/ID00012637202177665765362/10.dcm\"\n\nfig, ax = plt.subplots(1, 2, figsize=(10, 10))\nds = pydicom.dcmread(img_1)\nax[0].set_title('Patient 1: Ex-Smoker')\nax[0].imshow(ds.pixel_array, cmap=plt.cm.bone)\n\nds = pydicom.dcmread(img_2)\nax[1].set_title('Patient 2: Never smoked')\nax[1].imshow(ds.pixel_array, cmap=plt.cm.bone)\n\nplt.show","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Modeling"},{"metadata":{},"cell_type":"markdown","source":"## compute loss (and also applied to the Evaluation metric)\nRef: https://www.kaggle.com/ttahara/osic-baseline-lgbm-with-custom-metric"},{"metadata":{"trusted":true},"cell_type":"code","source":"def baseline_loss_metric(trueFVC, predFVC, predSTD=100):\n    clipSTD = np.clip(predSTD, 70 , 9e9)  \n    deltaFVC = np.clip(np.abs(trueFVC - predFVC), 0 , 1000)  \n#     error = -1 * (np.sqrt(2) * deltaFVC / clipSTD) - np.log(np.sqrt(2) * clipSTD)\n    error = np.mean(-1 * (np.sqrt(2) * deltaFVC / clipSTD) - np.log(np.sqrt(2) * clipSTD))\n    return error","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class model_selection(): \n    \n    def __init__(self): \n        self.y_pred_FVC = pd.DataFrame()\n#         self.my_scorer = make_scorer(baseline_loss_metric, greater_is_better=False)\n        self.best_param = None\n        self.scoring_train = 0\n        self.scoring_val = 0\n    \n    \n    def xgboost(self, X, y, X_val, y_val,test): \n        parameters = {'learning_rate': [0.1], 'n_estimators':[10000], 'max_depth':[3], 'reg_alpha':[0.005]}\n        clf = GridSearchCV(XGBRegressor(min_child_weight=0, gamma=0, \n                                        colsample_bytree=0.7, objective='reg:squarederror', nthread=-1,\n                                        scale_pos_weight=1, subsample=.7, seed=27), \n                           param_grid=parameters, \n                           scoring=make_scorer(baseline_loss_metric, greater_is_better=True))\n        clf.fit(X, y)\n        self.best_param = clf.best_params_\n        self.scoring_train = clf.score(X, y)\n        self.scoring_val = clf.score(X_val, y_val)\n        y_pred_xgb_FVC = clf.predict(test)\n        \n        self.y_pred_FVC = pd.concat([pd.Series(test.index), pd.Series(y_pred_xgb_FVC)], axis=1)\n#         self.y_pred_FVC = y_pred_FVC.groupby('Patient_Week').mean()\n        return self.y_pred_FVC, self.best_param, self.scoring_train, self.scoring_val\n        \n    \n    def lightgbm(self, X, y, X_val, y_val,test): \n        parameters = {'learning_rate': [0.05], 'n_estimators':[3000], 'num_leaves':[45]}\n        clf = GridSearchCV(LGBMRegressor(boosting_type='gbdt', \n                                         objective='regression',\n                                         bagging_fraction=0.8,\n                                         bagging_freq=1, \n                                         verbose=-1,), \n                           param_grid=parameters, \n                           scoring=make_scorer(baseline_loss_metric, greater_is_better=True))\n        clf.fit(X, y)\n        self.best_param = clf.best_params_\n        self.scoring_train = clf.score(X, y)\n        self.scoring_val = clf.score(X_val, y_val)\n        y_pred_lgb_FVC = clf.predict(test)\n        \n        self.y_pred_FVC = pd.concat([pd.Series(test.index), pd.Series(y_pred_lgb_FVC)], axis=1)\n#         self.y_pred_FVC = self.y_pred_FVC.groupby('Patient_Week').mean()\n        return self.y_pred_FVC, self.best_param, self.scoring_train, self.scoring_val\n    \n    \n    def HuberRegressor(self, X, y, test): \n        hbr = HuberRegressor(max_iter=200)\n        hbr.fit(X, y)\n        y_pred_hbr_FVC = hbr.predict(test)\n        \n        self.y_pred_FVC = pd.concat([pd.Series(test.index), pd.Series(y_pred_hbr_FVC)], axis=1)\n#         self.y_pred_FVC = self.y_pred_FVC.groupby('Patient_Week').mean()\n        return self.y_pred_FVC","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### FVC Prediction"},{"metadata":{},"cell_type":"markdown","source":"#### XGBoost"},{"metadata":{"trusted":true},"cell_type":"code","source":"model = model_selection()\noutput = model.xgboost(X_train, y_train, X_val, y_val,test_FVC)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(output[1]) # best params on previous model {'learning_rate': 0.1, 'max_depth': 3, 'n_estimators': 10000, 'reg_alpha': 0.005}\nprint(output[2]) # error on training -4.954089770686458\nprint(output[3]) # error on validation -4.954233496482487","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y_pred_xgb_FVC = output[0]\ny_pred_xgb_FVC.rename(columns={0:'FVC'}, inplace=True)\ny_pred_xgb_FVC","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### LightGBM"},{"metadata":{"trusted":true},"cell_type":"code","source":"model = model_selection()\noutput = model.lightgbm(X_train, y_train, X_val, y_val,test_FVC)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(output[1]) # best params on previous model {'learning_rate': 0.05, 'n_estimators': 3000, 'num_leaves': 45}\nprint(output[2]) # error on training -4.97098194152317\nprint(output[3]) # error on validation -4.958517234575529","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y_pred_lgb_FVC = output[0]\ny_pred_lgb_FVC.rename(columns={0:'FVC'}, inplace=True)\ny_pred_lgb_FVC","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Confidence Prediction (on process)"},{"metadata":{"trusted":true},"cell_type":"code","source":"formula = 'FVC ~ Weeks+Percent+Age+Sex+SmokingStatus'\ndef QuantileRegression(train, test, printsumL=False,quanL=0.15, quanH=0.85): \n        modelL = quantreg(formula, train).fit(q=quanL)\n        modelH = quantreg(formula, train).fit(q=quanH)\n        test['FVC_L'] = modelL.predict(test).values\n        test['FVC_H'] = modelH.predict(test).values\n        if printsumL: \n            print(modelL.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_conf = pd.merge(sub_df, df, how='left', on=['Patient'])\ntest_conf = test_conf.rename(columns={'Patient_Week_x':'Patient_Week'})\ntest_conf = test_conf.drop(['Confidence', 'Weeks'], axis=1)\ntest_conf = test_conf.groupby('Patient_Week').mean()\ntest_conf.sort_values(['Patient_Week'], inplace=True)\n\ntest_conf = pd.merge(test_conf, y_pred_xgb_FVC, how='left', on=['Patient_Week']).rename(columns={'FVC_x':'FVC', 'FVC_y':'pred_FVC', 'pred_Weeks':'Weeks'})\ntest_conf.set_index('Patient_Week', inplace=True)\ntest_conf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"QuantileRegression(df, test_conf,True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# test_conf = test_conf[['Weeks', 'FVC', 'Age', 'Sex', 'SmokingStatus', 'Height']]\ntest_conf['Confidence'] = 0.55*np.abs(test_conf['FVC_H'] - test_conf['pred_FVC'])+0.45*np.abs(test_conf['pred_FVC'] - test_conf['FVC_L'])\ntest_conf.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# y_pred_conf = pd.DataFrame(competition_metric(test_conf['FVC'], test_conf['pred_FVC'], 100))\ny_pred_conf = test_conf[['Confidence']]\ny_pred_conf.rename(columns={0:'Confidence'}, inplace=True)\ny_pred_conf","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Evaluation"},{"metadata":{"trusted":true},"cell_type":"code","source":"def competition_metric(trueFVC, predFVC, predSTD):\n    clipSTD = np.clip(predSTD, 70 , 9e9)  \n    deltaFVC = np.clip(np.abs(trueFVC - predFVC), 0 , 1000)  \n#     error = -1 * (np.sqrt(2) * deltaFVC / clipSTD) - np.log(np.sqrt(2) * clipSTD)\n    error = np.mean(-1 * (np.sqrt(2) * deltaFVC / clipSTD) - np.log(np.sqrt(2) * clipSTD))\n    return error","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"competition_metric(test_conf['FVC'], test_conf['pred_FVC'], test_conf['Confidence'])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Merge pred_FVC, pred_Confidence into submission"},{"metadata":{"trusted":true},"cell_type":"code","source":"submission = pd.merge(sub_df,y_pred_xgb_FVC, how='left', on=['Patient_Week'])\nsubmission = pd.merge(submission,y_pred_conf, how='left', on=['Patient_Week']) # Predicted Confidence\nsubmission = submission[['Patient_Week', 'FVC', 'Confidence_y']].rename(columns={'Confidence_y':'Confidence'})\n\n# #baseline Confidence\n# submission = pd.merge(submission, test_FVC, how='left', on=['Patient_Week']) # Raw Confidence (Percent)\n# submission = submission[['Patient_Week', 'FVC', 'Confidence']].rename(columns={'Percent':'Confidence'}) # for baseline\n\nsubmission.set_index('Patient_Week', inplace=True)\nsubmission","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"submission.to_csv('./submission.csv')","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}