{"cells":[{"metadata":{},"cell_type":"markdown","source":"# **Load data**","execution_count":null},{"metadata":{"_uuid":"451276ce-0746-48ab-acc9-f0b2a2ac0e4d","_cell_guid":"738374b8-6220-4024-8ee7-7bf6f9b20b3f","trusted":true},"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\n# Loss function and lots of inspiration from https://www.kaggle.com/ulrich07/osic-multiple-quantile-regression-starter#kln-163\n\nimport os\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport torch as pt\nfrom torch.utils.data import Dataset, DataLoader, RandomSampler\nimport torch.optim as optim\nimport matplotlib.pyplot as plt\nimport pydicom\nfrom sklearn.model_selection import KFold\nimport seaborn as sns\nimport pydicom\nfrom glob import glob\nimport scipy.ndimage\nfrom skimage import morphology\nfrom skimage import measure\nfrom skimage.filters import threshold_otsu, median\nfrom scipy.ndimage import binary_fill_holes\nfrom skimage.segmentation import clear_border\nfrom scipy.stats import describe\n    \ntrainImagesPath = \"/kaggle/input/osic-pulmonary-fibrosis-progression/train/\"\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\ndtype = pt.float\nuse_cuda = pt.cuda.is_available()\ndevice = pt.device(\"cuda:0\" if use_cuda else \"cpu\")\n\n\n\n# You can write up to 5GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session\n#\n#inputCols = ['PercentIn','AgeIn',  'Smoke', 'Gender', 'WeekIn']\ninputCols = ['PercentIn', 'AgeIn', 'WeekIn', 'base_FVC', 'min_FVC', 'SmokingStatus_Currently smokes', 'SmokingStatus_Ex-smoker', 'SmokingStatus_Never smoked', 'Sex_Male', 'Sex_Female']\n#inputCols=['PercentIn']\ntest = pd.read_csv(\"/kaggle/input/osic-pulmonary-fibrosis-progression/test.csv\")\ntrain = pd.read_csv(\"/kaggle/input/osic-pulmonary-fibrosis-progression/train.csv\")\nsubms = pd.read_csv(\"../input/osic-pulmonary-fibrosis-progression/sample_submission.csv\")\nsubms['Patient'] = subms.Patient_Week.apply(lambda x: x.split(\"_\")[0])\nsubms['Weeks'] = subms.Patient_Week.apply(lambda x: int(x.split(\"_\")[-1]))\ntrain['Split'] = \"train\"\ntest[\"Split\"] = \"test\"\ntest['PatientDir'] = test.Patient.apply(lambda x: \"../input/osic-pulmonary-fibrosis-progression/test/\" + x)\ntrain['PatientDir'] = train.Patient.apply(lambda x: \"../input/osic-pulmonary-fibrosis-progression/train/\" + x)\nsubms['PatientDir'] = subms.Patient.apply(lambda x: \"../input/osic-pulmonary-fibrosis-progression/train/\" + x)\n\n#Add submissions\nsubms = subms[['Patient', 'Weeks', 'Patient_Week']]\nsubms = subms.merge(test.drop(\"Weeks\", axis=1), on=\"Patient\")\nsubms[\"Split\"] = \"subm\"\ndata = test.append([subms, train])\n\n# Scale week input relative to first week\ndata['first_week'] = data.Weeks\ndata['first_week'] = data.groupby('Patient')['first_week'].transform('min') # This assumes all patient ids in submission set exist in test set later. Set missing to train set mean?\ndata['first_week'] = data.Weeks - data.first_week\n\n# Minimum FVC and FVC on first test\ndata['min_FVC'] = data.groupby('Patient')['FVC'].transform('min')\ndata = data.sort_values(['Patient', 'Weeks'], ascending=True)\nb = data[data['first_week'] == 0]\nb['base_FVC'] = b.FVC\nb = b[['Patient', 'base_FVC']]\nb = b.drop_duplicates('Patient')\ndata = data.merge(b, on='Patient', how='left')\ndata['WeekIn'] = (data.first_week - data.first_week.min()) / (data.first_week.max() - data.first_week.min())\ndata = pd.concat([data, pd.get_dummies(data.SmokingStatus, prefix=\"SmokingStatus\")], axis=1)\ndata = pd.concat([data, pd.get_dummies(data.Sex, prefix=\"Sex\")], axis=1)\n\n# Normalize\ndata['Smoke'] = data.SmokingStatus.replace({'Ex-smoker' : 0.5, 'Never smoked' : 0, 'Currently smokes' : 1})\ndata['Gender'] = data.Sex.replace({'Male' : 1, 'Female' : 0})\ndata['min_FVC'] = (data.min_FVC - data.min_FVC.min()) / (data.min_FVC.max() - data.min_FVC.min())\ndata['base_FVC'] = (data.base_FVC - data.base_FVC.min()) / (data.base_FVC.max() - data.base_FVC.min())\ndata['WeekIn'] = (data.first_week - data.first_week.min()) / (data.first_week.max() - data.first_week.min())\ndata['AgeIn'] = (data.Age - data.Age.min()) / (data.Age.max() - data.Age.min())\ndata['PercentIn'] = (data.Percent - data.Percent.min()) / (data.Percent.max() - data.Percent.min())\n\ntrain = data.loc[data.Split == 'train'].copy()\nsubms = data.loc[data.Split == \"subm\"].copy()\ntest = data.loc[data.Split == \"test\"].copy()\ndata.corr()\n\n# I/O for nn\ninputs = [data.columns.get_loc(c) for c in inputCols if c in data.columns]\noutputColI = data.columns.get_loc('FVC')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# **Image features**","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"print(data['PatientDir'])\nprint((data.nunique()))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# **Preprocess data and define model**","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"class OSICDataSet(Dataset):\n    def __init__(self, data, mode = \"train\"):\n        self.data = data\n        self.mode = mode\n    \n    def __len__(self):\n        return len(self.data)\n    \n    def __getitem__(self, idx):\n        #image = self.data.loc[idx, ['Patient']].map(lambda filename: pydicom.dcmread(trainImagesPath + filename + \"/1.dcm\")).item().pixel_array.astype(np.float64)\n        otherData = pt.from_numpy(self.data.iloc[idx, inputs].values.astype(np.float32))\n        if self.mode == \"train\":\n            targets = pt.from_numpy(self.data.iloc[idx, outputColI].values.astype(np.float32))\n            return otherData, targets\n            #return otherData, targets\n        return otherData #TODO: Add image data\n    \nclass Model(pt.nn.Module):\n    def __init__(self):\n        super(Model, self).__init__()\n        self.start = pt.nn.Sequential(\n            pt.nn.Linear(len(inputCols), 128),\n            pt.nn.ReLU(),\n            pt.nn.Linear(128, 256),\n            pt.nn.ReLU(),\n            pt.nn.Linear(256, 512),\n            pt.nn.ReLU(),\n            pt.nn.Linear(512, 1024),\n            pt.nn.ReLU(),\n            pt.nn.Linear(1024, 512),\n            pt.nn.ReLU(),\n            pt.nn.Linear(512, 256),\n            pt.nn.ReLU(),\n        )\n        self.left = pt.nn.Sequential(\n            pt.nn.Linear(256, 128)\n        )\n        self.sigmoid = pt.nn.Sigmoid()\n        self.right = pt.nn.Sequential(\n            pt.nn.Linear(256,128)\n        )\n        self.last = pt.nn.Sequential(\n            pt.nn.Linear(128, 3),\n        )\n        self.lastRelu = pt.nn.Sequential(\n            pt.nn.Linear(128, 3),\n            pt.nn.ReLU()\n        )\n        \n        \n    def forward(self, x):\n        h = self.start(x)\n        l = self.left(h)\n        r = self.right(h)\n        h = l * self.sigmoid(r)\n        p1 = self.last(h)\n        p2 = self.lastRelu(h)\n        out = p1 + pt.cumsum(p2, 1)\n        return out\n\nmodel = Model()\n\nC1, C2 = pt.FloatTensor([70]), pt.FloatTensor([1000])\n#=============================#\ndef score(y_pred, y_true):\n    sigma = y_pred[:, 2] - y_pred[:, 0]\n    fvc_pred = y_pred[:, 1]\n\n    #sigma_clip = sigma + C1\n    sigma_clip = pt.max(sigma, C1.expand_as(sigma))\n    delta = pt.abs(y_true - fvc_pred)\n    delta = pt.min(delta, C2)\n    sq2 = pt.sqrt(pt.FloatTensor([2]))\n    metric = (delta / sigma_clip)*sq2 + pt.log(sigma_clip* sq2)\n    return pt.mean(metric)\n\ndef pinballLoss(pred, label, quant):\n    err = label.unsqueeze(1) - pred\n    m = pt.mean(pt.max(quant * err, (quant-1)*err))\n    return m\n\ndef loss(pred, label):\n    quantiles = pt.FloatTensor([0.2, 0.5, 0.8])\n    return 0.8 * pinballLoss(pred, label, quantiles) + 0.2 * score(pred, label)\n\ndef earlyStop(useEarlyStopping, earlyStopping):\n    if useEarlyStopping:\n        return not earlyStopping\n    return True","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# **Train**","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"trainDataSet = OSICDataSet(train) \noptimizer = optim.Adam(model.parameters(), lr=0.075, weight_decay = 0.1, eps=0.001)\nlosses = []\nvalLosses = []\nnSplits = 5\nkf = KFold(n_splits=5)\nc = 0\nuseEarlyStopping = True\nearlyStopping = False\nsubmDataSet = OSICDataSet(subms, \"submission\")\nstopping = 0\ni = 0\nwhile i < 800 and earlyStop(useEarlyStopping, earlyStopping):\n    for train_index, val_index in kf.split(trainDataSet):\n        train_x, train_y = trainDataSet[train_index]\n        val_x, val_y = trainDataSet[val_index]\n        optimizer.zero_grad()\n        pred = model(train_x)\n        los = loss(pred, train_y)\n        los.backward()\n        if c % 10 == 0:\n            pred = model(val_x)\n            valLos = loss(pred, val_y)\n            losses.append(los)\n            valLosses.append(valLos)\n            if useEarlyStopping and (len(valLosses) > 1 and valLosses[-1] < valLosses[-2]):\n                stopping += 1\n                if stopping > nSplits:\n                    earlyStopping = True\n                break\n        c += 1\n        i += 1\n        stopping = 0\n        optimizer.step()\n        c+=1\n    ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# **Plots, submission**","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.set(style=\"white\", palette=\"muted\", color_codes=True)\n\nplt.plot(losses, label=\"Train loss\")\nplt.plot(valLosses, label=\"Val loss\")\nplt.legend()\nplt.show()\n\nasd = pt.from_numpy(trainDataSet.data[inputCols].values.astype(np.float32))\npred = model(asd).detach()\nidxs = np.random.randint(0, len(trainDataSet), 100)\nplt.plot(trainDataSet.data.iloc[idxs, outputColI].values.astype(np.float32), 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()\nc = pred[:,2] - pred[:,0]\nf = pred[:,1]\n\n\nfig, axes = plt.subplots(1, 2, figsize=(13, 6.5))\nsns.distplot(c, color=\"g\", kde=False, kde_kws={\"shade\": True}, ax=axes[0]).set_title(\"Predicted Confidence on train set\")\nsns.distplot(f, color=\"g\", kde=False, kde_kws={\"shade\": True}, ax=axes[1]).set_title(\"Predicted FVC on train set\")\nplt.show()\n\nsubmDataSet = OSICDataSet(subms, \"submission\")\ninp = pt.from_numpy(submDataSet.data[inputCols].values.astype(np.float32))\nout = model(inp).detach()\nconfidence = out[:,2] - out[:,0]\nfvc = out[:,1]\nfig, axes = plt.subplots(1, 2, figsize=(13, 6.5))\nsns.distplot(confidence, color=\"g\", kde=False, kde_kws={\"shade\": True}, ax=axes[0]).set_title(\"Confidence on submission set\")\nsns.distplot(fvc, color=\"g\", kde=False, kde_kws={\"shade\": True}, ax=axes[1]).set_title(\"FVC on submision set\")\nplt.show()\n\nsubmDataSet.data['pFVC'] = fvc\nsubmDataSet.data['cConf'] = confidence\n\nsubmission = pd.DataFrame({'Patient_Week' : submDataSet.data.Patient_Week, 'FVC' : fvc, 'Confidence' : confidence})\notest = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv')\nfor i in range(len(otest)):\n    submission.loc[submission['Patient_Week']==otest.Patient[i]+'_'+str(otest.Weeks[i]), 'FVC'] = otest.FVC[i]\n    submission.loc[submission['Patient_Week']==otest.Patient[i]+'_'+str(otest.Weeks[i]), 'Confidence'] = 0.1\n\nsubmission.to_csv(\"submission.csv\", index=False)\nnull_columns=submission.columns[submission.isnull().any()]\n#print(submission[submission.isnull().any(axis=1)][null_columns].head())\nprint(submission.head())\nprint(len(submission))","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}