{"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_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"from IPython.core.magic import register_cell_magic\nimport os\nfrom pathlib import Path\n\n\n@register_cell_magic\ndef write_and_run(line, cell):\n    argz = line.split()\n    file = argz[-1]\n    mode = 'w'\n    if len(argz) == 2 and argz[0] == '-a':\n        mode = 'a'\n    with open(file, mode) as f:\n        f.write(cell)\n    get_ipython().run_cell(cell)\n    \nPath('/kaggle/working/scripts').mkdir(exist_ok=True)\nmodels_dir = Path('/kaggle/working/models')\nmodels_dir.mkdir(exist_ok=True)","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:24.881947Z","iopub.execute_input":"2021-08-31T13:28:24.882370Z","iopub.status.idle":"2021-08-31T13:28:24.893327Z","shell.execute_reply.started":"2021-08-31T13:28:24.882333Z","shell.execute_reply":"2021-08-31T13:28:24.892199Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\nimport os\nimport json\nimport glob\nimport random\nimport collections\n\nimport numpy as np\nimport pandas as pd\nimport pydicom\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\nimport cv2\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nplt.style.use('seaborn-talk')\n\nimport torch\n\nfrom sklearn.metrics import roc_auc_score, roc_curve, auc\n\nfrom torchvision import transforms\n\n\nimport time\n\nimport torch\nfrom torch import nn\nfrom torch.utils.data import DataLoader, Dataset\nfrom sklearn import model_selection as sk_model_selection\nfrom torch.nn import functional as F\n\nfrom sklearn.model_selection import StratifiedKFold\nimport pytorch_lightning as pl\nfrom transformers import DeiTFeatureExtractor, DeiTForImageClassification, AutoConfig\nfrom pytorch_lightning.core.memory import ModelSummary\nimport sys\nsys.path.append('../input/efficientnetpyttorch3d/EfficientNet-PyTorch-3D')\nfrom efficientnet_pytorch_3d import EfficientNet3D\nimport joblib\nimport nibabel as nb\nimport tarfile\nfrom pathlib import Path\n!pip install --upgrade https://github.com/VincentStimper/mclahe/archive/numpy.zip\n\nfrom mclahe import mclahe","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:24.898551Z","iopub.execute_input":"2021-08-31T13:28:24.898985Z","iopub.status.idle":"2021-08-31T13:28:26.859575Z","shell.execute_reply.started":"2021-08-31T13:28:24.898946Z","shell.execute_reply":"2021-08-31T13:28:26.858392Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Config:\n    seed = 42\n    lr = 2e-8\n    data_dir = Path('/tmp/input/nifti')\n    mri_types = ['FL', 'T1', 'T1CE', 'T2']\n\npl.utilities.seed.seed_everything(Config.seed, workers=True)","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:26.861639Z","iopub.execute_input":"2021-08-31T13:28:26.861989Z","iopub.status.idle":"2021-08-31T13:28:26.878276Z","shell.execute_reply.started":"2021-08-31T13:28:26.861953Z","shell.execute_reply":"2021-08-31T13:28:26.877304Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# tar = tarfile.open('../input/tumorbrainnifty/input.tar')\n# tar.extractall('/tmp/input')\n# tar.close()","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:26.880295Z","iopub.execute_input":"2021-08-31T13:28:26.880936Z","iopub.status.idle":"2021-08-31T13:28:26.886167Z","shell.execute_reply.started":"2021-08-31T13:28:26.880899Z","shell.execute_reply":"2021-08-31T13:28:26.885234Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tqdm import tqdm\ntqdm.pandas()\n\nkfold_df = pd.read_csv(\"../input/braintumor-sampling/brain_tumor_kfold.csv\", dtype={'BraTS21ID': str})\nprint(kfold_df.shape)\n\ndef filt_func(name):\n    files = glob.glob(f'/tmp/input/nifti/train/{name}/*')\n    if len(files) == 0:\n        return False\n    \n    for file in files:\n        data = nb.load(file).get_fdata()\n        if data.max() == 0:\n            return False\n    return True\n\n\nkfold_df['present'] = kfold_df.BraTS21ID.apply(func=filt_func)\nkfold_df = kfold_df[kfold_df.present]\nprint(kfold_df.shape)\n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T14:03:57.743898Z","iopub.execute_input":"2021-08-31T14:03:57.744230Z","iopub.status.idle":"2021-08-31T14:04:11.866322Z","shell.execute_reply.started":"2021-08-31T14:03:57.744202Z","shell.execute_reply":"2021-08-31T14:04:11.863380Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\ndef load_volume(dp_id, split=\"train\"):\n    mri_volumes = []\n    for mri_type in Config.mri_types:\n        volume = nb.load(Config.data_dir / split / dp_id / f'{mri_type}_to_SRI.nii.gz').get_fdata()\n#         print(volume.max())\n\n        volume = volume - volume.mean()\n        volume = volume / volume.std()\n#         print(volume.max())\n        mri_volumes.append(volume)\n        \n    mri_volumes = np.stack(mri_volumes)\n#     print(mri_volumes.shape, mri_volumes.dtype, mri_volumes.max())\n    \n\n#     print(mri_volumes.shape, mri_volumes.dtype, mri_volumes.mean(), mri_volumes.min(), mri_volumes.max())\n#     print(mri_volumes.shape)\n    return mri_volumes","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:26.989049Z","iopub.execute_input":"2021-08-31T13:28:26.989655Z","iopub.status.idle":"2021-08-31T13:28:26.997246Z","shell.execute_reply.started":"2021-08-31T13:28:26.989617Z","shell.execute_reply":"2021-08-31T13:28:26.996304Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%write_and_run scripts/dataset.py\n\nfrom torchvision.transforms import ToTensor\n\nclass DataRetriever(Dataset):\n    def __init__(self, paths, files_dir, targets=None):\n        \n        self.paths = paths\n        self.targets = targets\n        self.files_dir = files_dir\n        if self.targets is None:\n            self.split = 'val'\n        else:\n            self.split = 'train'\n        \n    def __len__(self):\n        return len(self.paths)\n    \n    def __getitem__(self, index):\n        dp_id = self.paths[index]\n        print(dp_id)\n        \n        volume = load_volume(dp_id)\n        if not self.targets is None:\n            y = torch.tensor(abs(self.targets[index]), dtype=torch.float)\n        else:\n            y = []\n            \n        return {\"X\": torch.tensor(volume).float(), 'y': y}\n    \n    \n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:33:20.857605Z","iopub.execute_input":"2021-08-31T13:33:20.857967Z","iopub.status.idle":"2021-08-31T13:33:20.867002Z","shell.execute_reply.started":"2021-08-31T13:33:20.857935Z","shell.execute_reply":"2021-08-31T13:33:20.866184Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df = kfold_df[kfold_df.fold!=0]\ntrain_ds = DataRetriever(\n    train_df[\"BraTS21ID\"].values, \n    '../input/rsna-miccai-brain-tumor-radiogenomic-classification/train',\n    targets=train_df[\"MGMT_value\"].values\n)\n\nX = train_ds[48]['X']\n\n\n# for i in range(40, 100):\n#     X = train_ds[i]['X']\n#     print(i, X.mean().item(), X.sum().item())\n\n# train_dl = DataLoader(\n#     train_ds,\n#     batch_size=1,\n#     shuffle=True,\n#     num_workers=1,\n# )\n\n# next(iter(train_dl))['X'].shape","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:33:22.526621Z","iopub.execute_input":"2021-08-31T13:33:22.526999Z","iopub.status.idle":"2021-08-31T13:33:23.980258Z","shell.execute_reply.started":"2021-08-31T13:33:22.526966Z","shell.execute_reply":"2021-08-31T13:33:23.979097Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for file in glob.glob('/tmp/input/nifti/train/00142/*'):\n    img = nb.load(file).get_fdata()\n    print(img.shape)\n#     print(np.nonzero(img))\n    plt.imshow(img[:,:,77])\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:56:38.988102Z","iopub.execute_input":"2021-08-31T13:56:38.988440Z","iopub.status.idle":"2021-08-31T13:56:40.814197Z","shell.execute_reply.started":"2021-08-31T13:56:38.988411Z","shell.execute_reply":"2021-08-31T13:56:40.813159Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for mri_type in glob.glob('../input/rsna-miccai-brain-tumor-radiogenomic-classification/train/00142/*'):\n    files = glob.glob(f'{mri_type}/*')\n    print(len(files))\n    mid_file = sorted(files, key=lambda n: int(n.split('-')[-1].split('.')[0]))[len(files) // 2+10]\n    print(mid_file)\n    img = pydicom.read_file(mid_file).pixel_array\n#     img = nb.load(file).get_fdata()\n    print(img.shape)\n    plt.imshow(img)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:56:21.565293Z","iopub.execute_input":"2021-08-31T13:56:21.565656Z","iopub.status.idle":"2021-08-31T13:56:22.303395Z","shell.execute_reply.started":"2021-08-31T13:56:21.565626Z","shell.execute_reply":"2021-08-31T13:56:22.302576Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def calc_roc_auc(y_true, y_pred):\n#     print(y_true, y_pred)\n    \n    return roc_auc_score(y_true, y_pred) \n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.221163Z","iopub.status.idle":"2021-08-31T13:28:27.221750Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%write_and_run scripts/model.py\n\n\nclass Model(pl.LightningModule):\n    def __init__(self, net, lr):\n        super().__init__()\n        self.net = net\n        self.lr = lr\n    \n    def forward(self, x):\n        out = self.net(x)\n        return out\n    \n    def training_step(self, batch, batch_idx):\n        X, y = batch['X'], batch['y']\n        y_hat = self(X).squeeze(1)\n        loss = F.binary_cross_entropy_with_logits(y_hat, y)  \n        \n        print('train:')\n        print('X:', X.mean(), X.sum())\n\n        print('y:', y)\n        print('y_hat:', y_hat)\n        print('loss:', loss.item())\n\n        return {\n            'loss': loss,\n            'y': y,\n            'y_hat': y_hat\n        }\n    \n    \n    def validation_step(self, batch, batch_idx):\n        self.train()\n        X, y = batch['X'], batch['y']\n        y_hat = self(X).squeeze(1)\n        loss = F.binary_cross_entropy_with_logits(y_hat, y)  \n\n        \n        print('val:')\n        print('X:', X.mean(), X.sum())\n        print('y:', y)\n        print('y_hat:', y_hat)\n        print('loss:', loss.item())\n        print('\\n\\n')\n        \n        return {\n            'loss': loss,\n            'y': y,\n            'y_hat': y_hat\n        }\n    \n    def predict_step(self, batch, batch_idx: int, dataloader_idx: int = None):\n        self.train()\n        X = batch['X']\n        return self(X)    \n    \n    def configure_optimizers(self):\n        optimizer = torch.optim.Adam(self.parameters(), lr=self.lr )\n        return optimizer\n    ","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.227017Z","iopub.status.idle":"2021-08-31T13:28:27.227622Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from time import time\nclass Metrics:\n    def __init__(self):\n        self.losses = []\n        self.reduced_losses = []\n        self.y_batches = []\n        self.y_hat_batches = []\n        self.roc_auc_list = []\n        \n        self.train_epoch_start_time = None\n        self.validation_epoch_start_time = None\n\n\nclass MetricsCallback(pl.callbacks.Callback):\n    def __init__(self):\n        self.train_metrics = Metrics()\n        self.validation_metrics = Metrics()\n        self.best_validation_roc_auc = float('-inf')\n        \n    def on_train_epoch_start(self, trainer, pl_module):\n        self.train_epoch_start_time = time()\n        \n    def on_validation_epoch_start(self, trainer, pl_module):\n        self.validation_epoch_start_time = time()\n    \n        \n    def on_train_batch_end(self, trainer, pl_module, outputs, batch, batch_idx, dataloader_idx):\n        self.after_batch(self.train_metrics, outputs)\n        \n    def on_validation_batch_end(self, trainer, pl_module, outputs, batch, batch_idx, dataloader_idx):\n        self.after_batch(self.validation_metrics, outputs)\n        \n        \n    def on_train_epoch_end(self, trainer, pl_module):\n        print('train epoch took', time() - self.train_epoch_start_time)\n        self.train_epoch_start_time = None\n        \n        \n            \n    def on_validation_epoch_end(self, trainer, pl_module):\n        train_loss = self.get_avg_loss(self.train_metrics)\n        validation_loss = self.get_avg_loss(self.validation_metrics)\n\n        \n        train_roc_auc = self.get_roc_auc(self.train_metrics)\n        validation_roc_auc = self.get_roc_auc(self.validation_metrics)\n        \n        self.log('roc_auc', validation_roc_auc)\n        if validation_roc_auc > self.best_validation_roc_auc:\n            self.best_validation_roc_auc = validation_roc_auc\n        \n        print('validation epoch took', time() - self.validation_epoch_start_time)\n        self.validation_epoch_start_time = None\n        print('train loss:', train_loss)\n        print('validation loss:', validation_loss)\n        print('train roc_auc:', train_roc_auc)\n        print('validation roc_auc:', validation_roc_auc)\n        print()\n        \n        \n    def get_avg_loss(self, metrics):\n        avg_loss = np.array(metrics.losses).mean()\n        metrics.reduced_losses.append(avg_loss)\n        metrics.losses = []\n        return avg_loss\n        \n    def after_batch(self, metrics, outputs):\n        metrics.losses.append( outputs['loss'].item())\n        metrics.y_batches.append(outputs['y'])\n        metrics.y_hat_batches.append(outputs['y_hat'])\n        \n    def get_roc_auc(self, metrics):   \n        if not metrics.y_batches:\n            return None\n        y_np = torch.hstack(metrics.y_batches).detach().cpu().numpy()\n        y_hat_np = torch.hstack(metrics.y_hat_batches).detach().cpu().numpy()\n        roc_auc = calc_roc_auc(y_np, y_hat_np)\n        metrics.roc_auc_list.append(roc_auc)\n        \n        metrics.y_batches = []\n        metrics.y_hat_batches = []\n        return roc_auc","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.232527Z","iopub.status.idle":"2021-08-31T13:28:27.233403Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_metrics(metrics_callback):\n    train_losses = metrics_callback.train_metrics.reduced_losses\n    validation_losses = metrics_callback.validation_metrics.reduced_losses\n\n    train_roc_aucs = metrics_callback.train_metrics.roc_auc_list\n    validation_roc_aucs = metrics_callback.validation_metrics.roc_auc_list\n\n    plt.plot(train_losses, label='train_loss')\n    plt.plot(validation_losses, label='val_loss')\n    plt.legend()\n    plt.show()\n\n    plt.plot(train_roc_aucs, label='train_roc_auc')\n    plt.plot(validation_roc_aucs, label='val_roc_auc')\n    plt.legend()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.238466Z","iopub.status.idle":"2021-08-31T13:28:27.239072Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"roc_aucs = []\nfor fold_n in range(5):\n    print('Fold:', fold_n+1)\n    train_df = kfold_df[kfold_df.fold!=fold_n]\n#     train_df = kfold_df.iloc[:10]\n    train_ds = DataRetriever(\n    train_df[\"BraTS21ID\"].values, \n    '../input/rsna-miccai-brain-tumor-radiogenomic-classification/train',\n    targets=train_df[\"MGMT_value\"].values\n)\n\n    val_df = kfold_df[kfold_df.fold==fold_n]\n#     val_df = kfold_df.iloc[:10]\n    val_ds = DataRetriever(\n    val_df[\"BraTS21ID\"].values, \n    '../input/rsna-miccai-brain-tumor-radiogenomic-classification/train',\n    targets=val_df[\"MGMT_value\"].values\n)\n\n    train_dl = DataLoader(\n        train_ds,\n        batch_size=2,\n        shuffle=True,\n        num_workers=8,\n        pin_memory=True\n    )\n\n    val_dl = DataLoader(\n        val_ds, \n        batch_size=2,\n        shuffle=False,\n        num_workers=2,\n        pin_memory=True\n    )\n\n    net = EfficientNet3D.from_name(\"efficientnet-b3\", override_params={'num_classes': 1}, in_channels=4)\n    model = Model(net, Config.lr)\n    checkpoint_callback = pl.callbacks.ModelCheckpoint(dirpath='models', filename=f'model_{fold_n}_' + '{epoch}_{roc_auc:.3}', monitor='roc_auc', mode='max', save_weights_only=True)\n    metrics_callback = MetricsCallback()\n    print(ModelSummary(model))\n    trainer = pl.Trainer(fast_dev_run=False, max_epochs=10, gpus=1,\n                         auto_lr_find=True, precision=16, limit_train_batches=1.0, limit_val_batches=1.0, \n                         num_sanity_val_steps=0, val_check_interval=0.333, \n                         accumulate_grad_batches=4,\n                         callbacks=[metrics_callback, checkpoint_callback])\n\n\n\n#     lr_find = trainer.tune(model, train_dl, val_dl)['lr_find']\n#     fig = lr_find.plot(suggest=True)\n#     plt.show()\n#     new_lr = lr_find.suggestion()\n#     print('new_lr', new_lr)\n\n    trainer.fit(model, train_dl, val_dl)\n    plot_metrics(metrics_callback)\n    roc_aucs.append(metrics_callback.best_validation_roc_auc)\n\n\nprint('roc_auc_s:', roc_aucs)\nprint('roc_auc mean:', np.array(roc_aucs).mean())\nprint('roc_auc std:', np.array(roc_aucs).std())\n\n\n\n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.241635Z","iopub.status.idle":"2021-08-31T13:28:27.242218Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!ls /tmp/input/nifti/train/0008*","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.244587Z","iopub.status.idle":"2021-08-31T13:28:27.245173Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import gc\n# del model, trainer, train_dl, val_dl, train_ds, val_ds\ngc.collect()\ntorch.cuda.empty_cache()","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.246531Z","iopub.status.idle":"2021-08-31T13:28:27.247100Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# !rm -rf models/*","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.250650Z","iopub.status.idle":"2021-08-31T13:28:27.251263Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!ls --block-size=M -l models","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.254581Z","iopub.status.idle":"2021-08-31T13:28:27.255169Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# base_model = Model(deit_model)\n\n# for ckpt_name in glob.glob('models/*.ckpt'):\n#     state_dict = torch.load(ckpt_name)['state_dict']\n#     base_model.load_state_dict(state_dict)\n#     fold_n = int(ckpt_name.split('_')[0][-1])\n    \n#     val_df = kfold_df[kfold_df.fold==fold_n]\n#     val_ds = DataRetriever(\n#         val_df[\"BraTS21ID\"].values, \n#         '../input/rsna-miccai-brain-tumor-radiogenomic-classification/train',\n#         targets=val_df[\"MGMT_value\"].values\n#     )\n\n#     val_dl = DataLoader(\n#         val_ds, \n#         batch_size=8,\n#         shuffle=False,\n#         num_workers=16,\n#         collate_fn=FeatureExtractorCollate(feats_extractor)\n#     )\n#     metrics_callback = MetricsCallback()\n\n#     trainer = pl.Trainer(gpus=1, num_sanity_val_steps=0, callbacks=[metrics_callback])\n#     trainer.validate(base_model, val_dl)\n\n\n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.257611Z","iopub.status.idle":"2021-08-31T13:28:27.258244Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# base_model = Model(deit_model)\n# test_df = pd.read_csv(\"../input/rsna-miccai-brain-tumor-radiogenomic-classification/sample_submission.csv\")\n\n# for ckpt_name in glob.glob('models/*.ckpt'):\n#     state_dict = torch.load(ckpt_name)['state_dict']\n#     base_model.load_state_dict(state_dict)\n#     fold_n = int(ckpt_name.split('_')[0][-1])\n    \n#     test_ds = DataRetriever(\n#         test_df[\"BraTS21ID\"].values, \n#         '../input/rsna-miccai-brain-tumor-radiogenomic-classification/test',\n#         targets=test_df[\"MGMT_value\"].values\n#     )\n\n#     test_dl = DataLoader(\n#         test_ds, \n#         batch_size=8,\n#         shuffle=False,\n#         num_workers=16,\n#         collate_fn=FeatureExtractorCollate(feats_extractor)\n#     )\n\n#     trainer = pl.Trainer(gpus=1, num_sanity_val_steps=0)\n#     preds = trainer.predict(base_model, test_dl)\n#     print([p.flatten().tolist() for p in preds])\n#     print()\n","metadata":{"execution":{"iopub.status.busy":"2021-08-31T13:28:27.260594Z","iopub.status.idle":"2021-08-31T13:28:27.261181Z"},"trusted":true},"execution_count":null,"outputs":[]}]}