{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport pydicom\nimport os\nimport random\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\n\nfrom tqdm import tqdm\nimport torch\n\nimport warnings\n\nwarnings.simplefilter('ignore')\n\nfrom torch.distributions import Normal\nfrom torch.utils.data import Dataset, DataLoader\nimport torch\nimport torchvision\nimport torch.nn.functional as F\nimport torch.nn as nn\nfrom torch.optim.lr_scheduler import ReduceLROnPlateau, StepLR, CosineAnnealingLR\n\n# Competation metric for numpy array\ndef metric(outputs, std, target):\n    confidence = std\n    clip = np.where(confidence > 70, confidence, 70)\n    delta = np.abs(outputs - target)\n    delta = np.where(delta > 1000, 1000, delta)\n\n    metrics = (delta * np.sqrt(2) / clip) + np.log(clip * np.sqrt(2))\n\n    return np.mean(metrics)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class LungDataset(Dataset):\n    def __init__(self, df, train=True, meta_features=None):\n\n        self.df = df\n        self.train = train\n        self.meta_features = meta_features\n\n    def __getitem__(self, index):\n\n        # print(self.df.iloc[index]['Patient'])\n        meta = np.array(self.df.iloc[index][self.meta_features].values, dtype=np.float32)\n\n        if self.train:\n            y = self.df.iloc[index]['target']\n            return (meta), y\n        else:\n            return (meta)\n\n    def __len__(self):\n        return len(self.df)\n\n\nclass Bayesian_Regression(nn.Module):\n\n    def __init__(self, input_size):\n        super(Bayesian_Regression, self).__init__()\n\n        self.fc1 = nn.Linear(input_size, 128)\n        self.fc2 = nn.Linear(128, 64)\n        self.fc3 = nn.Linear(64, 32)\n        self.std = nn.Linear(32, 1)\n        self.mean = nn.Linear(32, 1)\n\n    def forward(self, X):\n        h = self.fc1(X)\n        h = self.fc2(h)\n        h = self.fc3(h)\n\n        std = F.softplus(self.std(h))\n        mean = self.mean(h)\n\n        return mean, std\n\n\nfrom torch.utils.data import Dataset, DataLoader\n\n\nclass EarlyStopping:\n    \"\"\"Early stops the training if validation loss doesn't improve after a given patience.\"\"\"\n\n    def __init__(self, patience=7, verbose=False, delta=0):\n\n        self.patience = patience\n        self.verbose = verbose\n        self.counter = 0\n        self.best_score = None\n        self.early_stop = False\n        self.val_loss_min = np.Inf\n        self.delta = delta\n\n    def __call__(self, val_loss, model, path):\n\n        score = -val_loss\n\n        if self.best_score is None:\n            self.best_score = score\n            self.save_checkpoint(val_loss, model, path)\n        elif score < self.best_score - self.delta:\n            self.counter += 1\n            print(f'EarlyStopping counter: {self.counter} out of {self.patience}')\n            if self.counter >= self.patience:\n                self.early_stop = True\n        else:\n            self.best_score = score\n            self.save_checkpoint(val_loss, model, path)\n            self.counter = 0\n\n    def save_checkpoint(self, val_loss, model, path):\n        '''Saves model when validation loss decrease.'''\n        if self.verbose:\n            print(f'Validation loss decreased ({self.val_loss_min:.6f} --> {val_loss:.6f}).  Saving model ...')\n        torch.save(model, path)\n        self.val_loss_min = val_loss\n\n\ndef training(model, train_data, val_data, epochs,\n             idx, batch_size=64, lr=0.001, patience=5):\n\n    trainloader = DataLoader(train_data, batch_size=batch_size, shuffle=True)\n    valloader = DataLoader(val_data, batch_size=batch_size, shuffle=True)\n\n    if os.path.isfile('checkpoint_%s.pt' % idx):\n\n        # load the last checkpoint with the best model\n        model = (torch.load('checkpoint_%s.pt' % idx))\n\n        return model, valloader\n\n    else:\n\n        optimizer = torch.optim.Adam(model.parameters(), lr=lr)\n        scheduler = ReduceLROnPlateau(optimizer=optimizer, mode='min', patience=3, verbose=True)\n        train_losses = []\n        valid_losses = []\n        avg_train_losses = []\n        avg_valid_losses = []\n        val_loss_torch = []\n        early_stopping = EarlyStopping(patience=patience, verbose=True)\n\n        for e in range(epochs):\n\n            model.train()\n            for x, y in tqdm(trainloader):\n\n                trainX = torch.tensor(x, device=device, dtype=torch.float32)\n                trainY = torch.tensor(y, device=device, dtype=torch.double)\n\n                # Forward pass\n                mean, std = model(trainX)\n\n                gaussian = torch.distributions.normal.Normal(mean.double(), std.double())\n                pred_y = gaussian.sample([1000])\n                likelihood = gaussian.log_prob(trainY)\n\n                loss = torch.mean(torch.log(std)/2 + torch.square(trainY-mean)/(2*std))-torch.mean(likelihood)\n\n                # Backward and optimize\n                optimizer.zero_grad()\n                loss.backward()\n\n                # nn.utils.clip_grad_norm_(meta_reg.parameters(), 5)\n                optimizer.step()\n                train_losses.append(loss.item())\n            \n            pred_y_mean = []\n            pred_y_std = []\n            true_y = []\n            model.eval()\n            for val_x, val_y in tqdm(valloader):\n                valX = torch.tensor(val_x, device=device, dtype=torch.float32)\n                valY = torch.tensor(val_y, device=device, dtype=torch.double)\n\n                mean, std = model(valX)\n\n                gaussian = torch.distributions.normal.Normal(mean.double(), std.double())\n                pred_y = gaussian.sample([1000])\n                likelihood = gaussian.log_prob(trainY)\n\n                val_loss = torch.mean(torch.log(std)/2 + torch.square(trainY-mean)/(2*std)) - torch.mean(likelihood)\n\n                valid_losses.append(val_loss.item())\n                val_loss_torch.append(val_loss.unsqueeze(0))\n                \n                pred_y_mean.append(mean)\n                pred_y_std.append(std)\n                true_y.append(valY)\n            \n            pred_y_mean = torch.cat(pred_y_mean, dim=0).detach().cpu().numpy()\n            pred_y_std = torch.cat(pred_y_std, dim=0).detach().cpu().numpy()\n            true_y = torch.cat(true_y, dim=0).detach().cpu().numpy()\n\n            train_loss = np.average(train_losses)\n            valid_loss = np.average(valid_losses)\n            avg_train_losses.append(train_loss)\n            avg_valid_losses.append(valid_loss)\n\n            scheduler.step(torch.mean(torch.cat(val_loss_torch, dim=0)))\n            epoch_len = len(str(epochs))\n\n            print_msg = (f'[{e:>{epoch_len}}/{epochs:>{epoch_len}}] ' +\n                         f'train_loss: {train_loss:.5f} ' +\n                         f'valid_loss: {valid_loss:.5f}' + \n                        f'metric: {metric(pred_y_mean, pred_y_std, true_y)}')\n\n            print(print_msg)\n\n            early_stopping(valid_loss, model, path='checkpoint_%s.pt' % idx)\n\n            if early_stopping.early_stop:\n                print(\"Early stopping\")\n                break\n\n        # load the last checkpoint with the best model\n        model = torch.load('checkpoint_%s.pt' % idx)\n\n        return model, valloader\n\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"ROOT = \"../input/osic-pulmonary-fibrosis-progression\"\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":"tr['WHERE'] = 'train'\nchunk['WHERE'] = 'val'\nsub['WHERE'] = 'test'\ndata = tr.append([chunk, sub])\n\ndata['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.loc[data.Weeks == data.min_week]\nbase = base[['Patient', 'FVC', 'Percent']].copy()\nbase.columns = ['Patient','min_FVC', 'min_percent']\nbase['nb'] = 1\nbase['nb'] = base.groupby('Patient')['nb'].transform('cumsum')\nbase = base[base.nb==1]\nbase.drop('nb', axis=1, inplace=True)\n\ndata = data.merge(base, on='Patient', how='left')\ndata['base_week'] = data['Weeks'] - data['min_week']\ndel base\n\nCOLS = ['Sex','SmokingStatus'] #,'Age'\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\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['min_percent_norm'] = (data['min_percent'] - data['min_percent'].min() ) / ( data['min_percent'].max() - data['min_percent'].min() )\ndata['target'] = data['FVC']#/data['FVC'].max()\nFE += ['age','week','BASE', 'min_percent_norm']\n\nprint(FE)\n\n# data.to_csv('df.csv', index=False)\n\ntr = data.loc[data.WHERE=='train']\nchunk = data.loc[data.WHERE=='val']\nsub = data.loc[data.WHERE=='test']\n\nkf = KFold(n_splits=5)\n\nidx=0\n\npred_test = np.zeros((sub.shape[0], 2))\npred_train = np.zeros((tr.shape[0], 3))\n\nmodels = []\n\nfor tr_idx, val_idx in kf.split(tr.index):\n\n    model = Bayesian_Regression(input_size=len(FE))\n    model = model.to(device)\n\n    train_data = LungDataset(tr.loc[tr_idx], meta_features=FE)\n    val_data = LungDataset(tr.loc[val_idx], meta_features=FE)\n\n    model_nn, val_loader = training(model, train_data, val_data, 1000,\n             idx, batch_size=64, lr=0.025, patience=5)\n    \n    models.append(model_nn)\n\n    pred_y_val = []\n    std_y_val = []\n    true_y_val = []\n    model.eval()\n    for val_x, val_y in tqdm(val_loader):\n        valX = torch.tensor(val_x, device=device, dtype=torch.float32)\n        mean, std = model_nn(valX)\n\n        #gaussian = torch.distributions.normal.Normal(mean.double(), std.double())\n        #pred_y = gaussian.sample([1000])\n\n        pred_y_val.append(mean)\n        std_y_val.append(std)\n        true_y_val.append(val_y)\n\n    pred_y_mean = torch.cat(pred_y_val, dim=0).detach().cpu().numpy()\n    pred_y_std = torch.cat(std_y_val, dim=0).detach().cpu().numpy()\n    true_y = torch.cat(true_y_val, dim=0).detach().cpu().numpy()\n    true_y = true_y #* data['FVC'].max()\n\n    plt.plot(true_y)\n    plt.plot(pred_y_mean)\n    plt.plot(pred_y_mean+pred_y_std)\n    plt.show()\n    \n    print(metric(pred_y_mean, pred_y_std, true_y))\n\n\n    test_data = LungDataset(sub, train=False, meta_features=FE)\n    testloader = DataLoader(test_data, batch_size=256, shuffle=False)\n\n    pred_y_test = []\n    std_y_test = []\n    model.eval()\n    for test_x in tqdm(testloader):\n        testX = torch.tensor(test_x, device=device, dtype=torch.float32)\n        mean, std = model_nn(testX)\n\n        pred_y_test.append(mean)\n        std_y_test.append(std)\n\n    pred_y_mean = torch.cat(pred_y_test, dim=0).detach().cpu().numpy()\n    pred_y_std = torch.cat(std_y_test, dim=0).detach().cpu().numpy()\n\n    pred_test[:, 0] = pred_test[:, 0] + pred_y_mean[:, 0]\n    pred_test[:, 1] = pred_test[:, 1] + pred_y_std[:, 0]\n\n    idx=idx+1\n\npred_test = pred_test/5\n\nsub = pd.read_csv(f\"{ROOT}/sample_submission.csv\")\nsub[['FVC','Confidence']] = pred_test\n\nsub.to_csv(\"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}