{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np\nimport random\nimport pandas as pd\nimport pydicom\nimport os\nimport matplotlib.pyplot as plt\nfrom tqdm import tqdm\nfrom PIL import Image\nfrom sklearn.metrics import mean_absolute_error\nfrom sklearn.model_selection import KFold","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import tensorflow as tf\nimport tensorflow.keras.backend as K\nimport tensorflow.keras.layers as L\nimport tensorflow.keras.models as M","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def 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":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"ROOT = \"../input/osic-pulmonary-fibrosis-progression\"\n#DESIRED_SIZE = 256 # Memory issue\nDESIRED_SIZE = 128","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"tr = 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":"tr['WHERE'] = 'train'\nchunk['WHERE'] = 'val'\nsub['WHERE'] = 'test'\ndata = tr.append([chunk, sub])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(tr.shape, chunk.shape, sub.shape, data.shape)\nprint(tr.Patient.nunique(), chunk.Patient.nunique(), sub.Patient.nunique(), \n      data.Patient.nunique())\n#","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')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"base = data.loc[data.Weeks == data.min_week]\nbase = base[['Patient','FVC']].copy()\nbase.columns = ['Patient','min_FVC']\nbase['nb'] = 1\nbase['nb'] = base.groupby('Patient')['nb'].transform('cumsum')\nbase = base[base.nb==1]\nbase.drop('nb', axis=1, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data = data.merge(base, on='Patient', how='left')\ndata['base_week'] = data['Weeks'] - data['min_week']\ndel base","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"COLS = ['Sex','SmokingStatus']\nFE = []\nfor col in COLS:\n    for mod in data[col].unique():\n        FE.append(mod)\n        data[mod] = (data[col] == mod).astype(int)\n#=================","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#\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() )\nFE += ['age','percent','week','BASE']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","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']\ndel 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":"### Quick Image processing","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_images(df, how=\"train\"):\n    xo = []\n    p = []\n    w  = []\n    for i in tqdm(range(df.shape[0])):\n        patient = df.iloc[i,0]\n        week = df.iloc[i,1]\n        try:\n            img_path = f\"{ROOT}/{how}/{patient}/{week}.dcm\"\n            ds = pydicom.dcmread(img_path)\n            im = Image.fromarray(ds.pixel_array)\n            im = im.resize((DESIRED_SIZE,DESIRED_SIZE)) \n            im = np.array(im)\n            xo.append(im[np.newaxis,:,:])\n            p.append(patient)\n            w.append(week)\n        except:\n            pass\n    data = pd.DataFrame({\"Patient\":p,\"Weeks\":w})\n    return np.concatenate(xo, axis=0), data","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### BASELINE CNN ","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"from sklearn.linear_model import Ridge, ElasticNet","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def metric( trueFVC, predFVC, predSTD ):\n    \n    clipSTD = np.clip( predSTD, 70 , 9e9 )  \n    \n    deltaFVC = np.clip( np.abs(trueFVC-predFVC), 0 , 1000 )  \n\n    return np.mean( -1*(np.sqrt(2)*deltaFVC/clipSTD) - np.log( np.sqrt(2)*clipSTD ) )\n#","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"C1, C2 = tf.constant(70, dtype='float32'), tf.constant(1000, dtype=\"float32\")\n#=============================#\ndef kloss(y_true, y_pred):\n    tf.dtypes.cast(y_true, tf.float32)\n    tf.dtypes.cast(y_pred, tf.float32)\n    sigma = y_pred[:, 1]\n    fvc_pred = y_pred[:, 0]\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 kmae(y_true, y_pred):\n    tf.dtypes.cast(y_true, tf.float32)\n    tf.dtypes.cast(y_pred, tf.float32)\n    spread = tf.abs( (y_true[:, 0] -  y_pred[:, 0])  / (y_pred[:, 0] + 1.) )\n    #spread = tf.abs( (y_true[:, 0] -  y_pred[:, 0])  / y_true[:, 0] )\n    return K.mean(spread)\n#=============================#\n\ndef mloss(_lambda):\n    def loss(y_true, y_pred):\n        return _lambda * kloss(y_true, y_pred) + (1 - _lambda)*kmae(y_true, y_pred)\n    return loss\n#=================\ndef make_model():\n    z = L.Input((9,), name=\"Patient\")\n    x = L.Dense(100, activation=\"relu\", name=\"d1\")(z)\n    x = L.Dense(100, activation=\"relu\", name=\"d2\")(x)\n    #x = L.Dense(100, activation=\"relu\", name=\"d3\")(x)\n    preds = L.Dense(2, activation=\"relu\", name=\"preds\")(x)\n    \n    model = M.Model(z, preds, name=\"CNN\")\n    model.compile(loss=mloss(0.99), \n                  optimizer=tf.keras.optimizers.Adam(lr=0.1, beta_1=0.9, beta_2=0.999, \n                                                     epsilon=None, decay=0.01, amsgrad=False), \n                  metrics=[kloss]) #, kmae\n    #model.compile(loss=kmae, optimizer=\"adam\", metrics=[kloss])\n    #model.compile(loss=kloss, optimizer=\"adam\", metrics=[kmae])#\n    return model","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#net = make_model()\n#print(net.summary())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"tr.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y = tr['FVC'].values\nz = tr[FE].values\nze = sub[FE].values","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"NFOLD = 5\nkf = KFold(n_splits=NFOLD)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pe0 = np.zeros((ze.shape[0], 2))\npred0 = np.zeros((z.shape[0], 2))\n\n\npe1 = np.zeros((ze.shape[0], 2))\npred1 = np.zeros((z.shape[0], 2))\n\ncnt = 0\nfor tr_idx, val_idx in kf.split(z):\n    cnt += 1\n    print(f\"FOLD {cnt}\")\n    print(\"=====================  NEURAL NET =============================\")\n    net = make_model()\n    net.fit(z[tr_idx], y[tr_idx], batch_size=200, epochs=500, \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=500))\n    print(\"val\", net.evaluate(z[val_idx], y[val_idx], verbose=0, batch_size=500))\n    print(\"predict val...\")\n    pred0[val_idx] = net.predict(z[val_idx], batch_size=500, verbose=0)\n    print(\"predict test...\")\n    pe0 += net.predict(ze, batch_size=500, verbose=0) / NFOLD\n    print(\"=====================  RIDGE REG =============================\")\n    clf = Ridge(alpha=0.05)\n    clf.fit(z[tr_idx], y[tr_idx]) #\n    #print(\"predict val...\")\n    pred1[val_idx, 0] = clf.predict(z[val_idx])\n    pred_std = np.mean(np.abs(y[val_idx] - pred1[val_idx, 0])) * np.sqrt(2)\n    pred1[val_idx, 1] = pred_std\n    print(\"val\", metric(y[val_idx], pred1[val_idx, 0], pred1[val_idx, 1]))\n    #print(\"predict test...\")\n    pe1[:, 0] += clf.predict(ze) / NFOLD\n    pe1[:, 1] += pred_std / NFOLD    \n#==============\npred0[:, 1] = pred0[:, 1] + 70.","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"w = 0.5\npred = (1-w) * pred0 + w * pred1\npe = (1-w) * pe0 + w * pe1\npe[:, 1] = pe1[:, 1]\npred[:, 1] = pred1[:, 1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(\"oof neural net\", metric(y, pred0[:, 0], pred0[:, 1]))\nprint(\"oof ridge\", metric(y, pred1[:, 0], pred1[:, 1]))\nprint(\"oof ensemble\", metric(y, pred[:, 0], pred[:, 1]))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sigma_opt = mean_absolute_error(y, pred[:, 0])\nsigma_mean = np.mean(pred[:, 1])\nprint(sigma_opt, sigma_mean)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.plot(y)\nplt.plot(pred[:, 0])\n#plt.plot(pred[:, 1])\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pred[:, 1].min(), pred[:, 1].max()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.hist(pred[:, 1])\nplt.title(\"uncertainty in prediction\")\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### PREDICTION","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub['FVC1'] = pe[:, 0]\nsub['Confidence1'] = pe[:, 1]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm = sub[['Patient_Week','FVC','Confidence','FVC1','Confidence1']].copy()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm.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":"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\n#","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":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm[[\"Patient_Week\",\"FVC\",\"Confidence\"]].to_csv(\"submission.csv\", index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","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}