{"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":"markdown","source":"## Use stacked images (3D) and MedicalNet model + volumentations\n\n\nIn this notebook you will find:\n\n* Use of volumentations with PyTorch\n* Training with pretrained MedicalNet model\n* k-fold cross validation\n\nAcknowledgements:\n\n- https://www.kaggle.com/ihelon/brain-tumor-eda-with-animations-and-modeling\n- https://www.kaggle.com/furcifer/torch-efficientnet3d-for-mri-no-train\n- https://github.com/shijianjian/EfficientNet-PyTorch-3D\n- https://www.kaggle.com/rluethy/efficientnet3d-with-one-mri-type\n    \n    \nUse models with only one MRI type, then ensemble the 4 models \n\nV14: add image rotation augmentation","metadata":{}},{"cell_type":"markdown","source":"TODO: \n\n* measure time volumentations -> DONE\n* use efficietnet3d with in_channels=4 -> UNDONE\n* natural sorting, sigmoid -> DONE\n* Review configuration: batch_norm, architecture, etc.\n* Add augmentation + external data\n* Do EDA to see if possible improve data quality.\n* Unfreeze layers in medical net\n\n\nResults\n\n","metadata":{}},{"cell_type":"code","source":"# thanks to https://www.kaggle.com/ipythonx\n!pip install ../input/keras-3d-model-and-3d-augmentation/volumentations_3D-1.0.3-py3-none-any.whl -q","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:26:47.622893Z","iopub.execute_input":"2021-10-05T10:26:47.623308Z","iopub.status.idle":"2021-10-05T10:27:15.041710Z","shell.execute_reply.started":"2021-10-05T10:26:47.623226Z","shell.execute_reply":"2021-10-05T10:27:15.040564Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nimport sys \nimport json\nimport glob\nimport random\nimport collections\nimport time\nimport re\nimport torchvision\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\n\nimport torch\nfrom torch import nn\nfrom torch.utils import data as torch_data\nfrom sklearn import model_selection as sk_model_selection\nfrom torch.nn import functional as torch_functional\nimport torch.nn.functional as F\nimport glob\n\nfrom sklearn.model_selection import StratifiedKFold\nfrom sklearn.metrics import roc_auc_score","metadata":{"papermill":{"duration":1.048295,"end_time":"2021-07-14T20:26:46.309722","exception":false,"start_time":"2021-07-14T20:26:45.261427","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:15.045583Z","iopub.execute_input":"2021-10-05T10:27:15.045864Z","iopub.status.idle":"2021-10-05T10:27:17.904944Z","shell.execute_reply.started":"2021-10-05T10:27:15.045832Z","shell.execute_reply":"2021-10-05T10:27:17.904129Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data_directory = '/kaggle/input/rsna-miccai-brain-tumor-radiogenomic-classification'\ninput_monaipath = \"/kaggle/input/monai-v060-deep-learning-in-healthcare-imaging/\"\nmonaipath = \"/kaggle/tmp/monai/\"\n\n!mkdir -p {monaipath}\n!cp -r {input_monaipath}/* {monaipath}\n\nLEARNING_RATE=0.0005\n\nsys.path.append(monaipath)\n\nfrom monai.networks.nets.densenet import DenseNet121","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:17.908352Z","iopub.execute_input":"2021-10-05T10:27:17.908605Z","iopub.status.idle":"2021-10-05T10:27:24.571006Z","shell.execute_reply.started":"2021-10-05T10:27:17.908578Z","shell.execute_reply":"2021-10-05T10:27:24.570015Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"if os.path.exists(\"../input/rsna-miccai-brain-tumor-radiogenomic-classification\"):\n    data_directory = '../input/rsna-miccai-brain-tumor-radiogenomic-classification'\n    pytorch3dpath = \"../input/efficientnetpyttorch3d/EfficientNet-PyTorch-3D\"\n    medicalpath = \"../input/medicalnet\" \n    sys.path.append(medicalpath)\nelse:\n    data_directory = '/media/roland/data/kaggle/rsna-miccai-brain-tumor-radiogenomic-classification'\n    pytorch3dpath = \"EfficientNet-PyTorch-3D\"\n    \nmri_types = ['FLAIR','T1w','T1wCE','T2w']\nPRETRAINED_WEIGHTS='../input/medicalnet-pretrained-weights/resnet_34.pth'\nSIZE = 256\nNUM_IMAGES = 64\nVOLUMENTATIONS = True # quite slow\nEPOCHS = 6\nFOLDS = 3\nPATCH_SIZE = (SIZE, SIZE, NUM_IMAGES)\nRESIZING_VOLUMENTATIONS = True\nDATA_PATH = '/kaggle/input/rsna-processed-voxels-64x256x256-clahe'\nsys.path.append(pytorch3dpath)\n\n\nfrom efficientnet_pytorch_3d import EfficientNet3D\nfrom models import resnet","metadata":{"lines_to_end_of_cell_marker":2,"lines_to_next_cell":2,"papermill":{"duration":0.05565,"end_time":"2021-07-14T20:26:46.486521","exception":false,"start_time":"2021-07-14T20:26:46.430871","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:24.572754Z","iopub.execute_input":"2021-10-05T10:27:24.573131Z","iopub.status.idle":"2021-10-05T10:27:24.619013Z","shell.execute_reply.started":"2021-10-05T10:27:24.573093Z","shell.execute_reply":"2021-10-05T10:27:24.618261Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Functions to load images","metadata":{}},{"cell_type":"code","source":"split = 'train'\ntrain_voxels = sorted(glob.glob(f\"{DATA_PATH}/voxels/*/*.npy\"))\n\ndf_train = pd.DataFrame(train_voxels, columns=['voxel_paths'])\ndf_train['BraTS21ID'] = df_train.voxel_paths.map(lambda path:path.split('/')[-1].strip('.npy'))\ndf_train['MRI_Type'] = df_train.voxel_paths.map(lambda path:path.split('/')[-2])\ndf_train_labels = pd.read_csv('../input/rsna-miccai-brain-tumor-radiogenomic-classification/train_labels.csv',\n                             dtype={'BraTS21ID':np.object,\n                                   'MGMT_value':np.int32})\ndf_train = df_train.set_index('BraTS21ID').join(df_train_labels.set_index('BraTS21ID'), on='BraTS21ID', how='left')\ndf_train = df_train.reset_index()\ndf_train.to_csv(\"/kaggle/working/df_train_meta.csv\")\ndf_train\n\n\n","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:24.620415Z","iopub.execute_input":"2021-10-05T10:27:24.620798Z","iopub.status.idle":"2021-10-05T10:27:24.809482Z","shell.execute_reply.started":"2021-10-05T10:27:24.620760Z","shell.execute_reply":"2021-10-05T10:27:24.808737Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Augmentations","metadata":{}},{"cell_type":"code","source":"from functools import partial\nfrom volumentations import *\n\ndef get_augmentation(patch_size):\n    return Compose([\n        #Rotate((-5, 5), (0, 0), (0, 0), p=0.5), # slow\n        #GridDropout(ratio=0., holes_number_y=10, p=1.0),\n        RandomCropFromBorders(crop_value=0.25, p=0.3), # p=0.3\n        #ElasticTransform((0, 0.15), interpolation=2, p=0.5), # slow\n        #Resize(patch_size, interpolation=1, always_apply=True, p=1.0),\n        Flip(0, p=0.5),\n        Flip(1, p=0.5),\n        RandomRotate90((0, 1), p=0.6),\n        #GaussianNoise(var_limit=(0, 5), p=1.0), # slow\n        RandomGamma(gamma_limit=(0.5, 1.5), p=0.7),\n    ], p=1.0)\n\nvolume3D = get_augmentation((SIZE, SIZE, NUM_IMAGES))\n\ndef load_dicom_image(path, img_size=SIZE, voi_lut=True, rotate=0):\n    dicom = pydicom.read_file(path)\n    data = dicom.pixel_array\n    if voi_lut:\n        data = apply_voi_lut(dicom.pixel_array, dicom)\n    else:\n        data = dicom.pixel_array\n        \n    if rotate > 0:\n        rot_choices = [0, cv2.ROTATE_90_CLOCKWISE, cv2.ROTATE_90_COUNTERCLOCKWISE, cv2.ROTATE_180]\n        data = cv2.rotate(data, rot_choices[rotate])\n        \n    data = cv2.resize(data, (img_size, img_size))\n    return data\n\n\ndef natural_sort(l): \n    convert = lambda text: int(text) if text.isdigit() else text.lower()\n    alphanum_key = lambda key: [convert(c) for c in re.split('([0-9]+)', key)]\n    return sorted(l, key=alphanum_key)\n\n\ndef load_normalized_images_3d(path, target_size=(64, 256, 256)):\n    \n    voxel = np.load(path).astype(np.float32) / 255.0\n    voxel = np.reshape(voxel, target_size)\n    voxel = np.transpose(voxel, (1, 2, 0))\n    voxel = np.expand_dims(voxel,0)\n    return voxel\n\ndef load_dicom_images_3d(scan_id, num_imgs=NUM_IMAGES, img_size=SIZE, \n                         mri_type=\"FLAIR\", split=\"train\", \n                         rotate=0, volumentations=True):\n    \n    files = natural_sort(glob.glob(f\"{data_directory}/{split}/{scan_id}/{mri_type}/*.dcm\"))\n\n    middle = len(files)//2\n    num_imgs2 = num_imgs//2\n    p1 = max(0, middle - num_imgs2)\n    p2 = min(len(files), middle + num_imgs2)\n    img3d = np.stack([load_dicom_image(f, rotate=rotate) for f in files[p1:p2]]).T \n    \n    if VOLUMENTATIONS and (rotate != 0): # rotate = 0 means test set\n        data = volume3D(**{\"image\":img3d})\n        data = data['image']\n\n        if RESIZING_VOLUMENTATIONS:\n            new_data = []\n            for z in range(data.shape[2]):\n                new_data.append(cv2.resize(data[:, :, z], (SIZE, SIZE)))\n            img3d = np.transpose(np.array(new_data), (1, 2, 0))\n        \n    if img3d.shape[-1] < num_imgs: # Fill gaps\n        n_zero = np.zeros((img_size, img_size, num_imgs - img3d.shape[-1]))\n        img3d = np.concatenate((img3d,  n_zero), axis = -1)\n        \n    if np.min(img3d) < np.max(img3d): # normalize\n        img3d = img3d - np.min(img3d)\n        img3d = img3d / np.max(img3d)\n    \n    img3d = np.expand_dims(img3d,0)\n        \n    return img3d\n\ninitial = time.time()\na = load_normalized_images_3d(df_train.loc[0].voxel_paths)\n#a = load_dicom_images_3d(\"00000\", mri_type=\"T1wCE\", rotate=1)\n\nprint(f'Time for reading a voxel: {time.time() - initial}')\n\nprint(a.shape)\nprint(np.min(a), np.max(a), np.mean(a), np.median(a))","metadata":{"papermill":{"duration":0.035761,"end_time":"2021-07-14T20:26:46.726756","exception":false,"start_time":"2021-07-14T20:26:46.690995","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:24.810748Z","iopub.execute_input":"2021-10-05T10:27:24.811097Z","iopub.status.idle":"2021-10-05T10:27:24.994283Z","shell.execute_reply.started":"2021-10-05T10:27:24.811048Z","shell.execute_reply":"2021-10-05T10:27:24.993278Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(16, 5))\ndataplot = a[0][:, :, :4]\nfor i in range(dataplot.shape[2]):\n    plt.subplot(1, 4, i + 1)\n    plt.imshow(dataplot[:, :, i], cmap=\"gray\")\n    plt.title(\"Sample\", fontsize=16)\n    plt.axis(\"off\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:24.995728Z","iopub.execute_input":"2021-10-05T10:27:24.996351Z","iopub.status.idle":"2021-10-05T10:27:25.258013Z","shell.execute_reply.started":"2021-10-05T10:27:24.996290Z","shell.execute_reply":"2021-10-05T10:27:25.257015Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_slices(num_rows, num_columns, width, height, data):\n    \"\"\"Plot a montage of 20 CT slices\"\"\"\n    data = np.rot90(np.array(data))  \n    data = np.transpose(data)\n    data = np.reshape(data, (num_rows, num_columns, width, height))\n    rows_data, columns_data = data.shape[0], data.shape[1]\n    heights = [slc[0].shape[0] for slc in data]\n    widths = [slc.shape[1] for slc in data[0]]\n    fig_width = 12.0\n    fig_height = fig_width * sum(heights) / sum(widths)\n    f, axarr = plt.subplots(\n        rows_data,\n        columns_data,\n        figsize=(fig_width, fig_height),\n        gridspec_kw={\"height_ratios\": heights},\n    )\n    for i in range(rows_data):\n        for j in range(columns_data):\n            axarr[i, j].imshow(data[i][j], cmap=\"gray\")\n            axarr[i, j].axis(\"off\")\n    plt.subplots_adjust(wspace=0, hspace=0, left=0, right=1, bottom=0, top=1)\n    plt.show()\n# Visualize montage of slices.\n# 5 rows and 10 columns for 100 slices of the CT scan.\nplot_slices(5, 10, SIZE, SIZE, a[0][:, :, :50])","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:25.260759Z","iopub.execute_input":"2021-10-05T10:27:25.261079Z","iopub.status.idle":"2021-10-05T10:27:27.387346Z","shell.execute_reply.started":"2021-10-05T10:27:25.261040Z","shell.execute_reply":"2021-10-05T10:27:27.386572Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def set_seed(seed):\n    random.seed(seed)\n    os.environ[\"PYTHONHASHSEED\"] = str(seed)\n    np.random.seed(seed)\n    torch.manual_seed(seed)\n    if torch.cuda.is_available():\n        torch.cuda.manual_seed_all(seed)\n        torch.backends.cudnn.deterministic = True\n\nSEED = 12\nset_seed(SEED)","metadata":{"papermill":{"duration":0.668331,"end_time":"2021-07-14T20:27:48.114522","exception":false,"start_time":"2021-07-14T20:27:47.446191","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:27.389173Z","iopub.execute_input":"2021-10-05T10:27:27.389638Z","iopub.status.idle":"2021-10-05T10:27:27.436978Z","shell.execute_reply.started":"2021-10-05T10:27:27.389602Z","shell.execute_reply":"2021-10-05T10:27:27.436166Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## train / test splits","metadata":{}},{"cell_type":"code","source":"display(df_train)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.438267Z","iopub.execute_input":"2021-10-05T10:27:27.438750Z","iopub.status.idle":"2021-10-05T10:27:27.453091Z","shell.execute_reply.started":"2021-10-05T10:27:27.438713Z","shell.execute_reply":"2021-10-05T10:27:27.452152Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train.groupby('BraTS21ID').count().describe()","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.454540Z","iopub.execute_input":"2021-10-05T10:27:27.455041Z","iopub.status.idle":"2021-10-05T10:27:27.482716Z","shell.execute_reply.started":"2021-10-05T10:27:27.455005Z","shell.execute_reply":"2021-10-05T10:27:27.481948Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"skf = StratifiedKFold(n_splits=FOLDS, random_state=SEED, shuffle=True)\n\npatient_df = df_train.groupby('BraTS21ID').MGMT_value.max().reset_index()\n\nprint('Class Ratio:',sum(patient_df['MGMT_value'])/len(patient_df['MGMT_value']))\n\ntarget = patient_df.loc[:,'MGMT_value']\n\nfold_no = 1\ntrain_fold_dict = {}\nval_fold_dict = {}\ntrain_indices = []\nval_indices = [] \nfor train_index, val_index in skf.split(patient_df, target):\n    train = patient_df.loc[train_index,:]\n    val = patient_df.loc[val_index,:]\n    train_indices.append(df_train.merge(train, on='BraTS21ID').index)\n    val_indices.append(df_train.merge(val, on='BraTS21ID').index) # add indices from the general dataframe with all MRI types\n    print('Fold',str(fold_no),'Class Ratio:',sum(patient_df.iloc[val_index]['MGMT_value'])/len(patient_df.iloc[val_index]),\n          ',\\t len train, val, sum:',len(train), len(val), len(train)+len(val))\n    fold_no += 1\n\n","metadata":{"papermill":{"duration":0.633753,"end_time":"2021-07-14T20:27:49.350524","exception":false,"start_time":"2021-07-14T20:27:48.716771","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:27.483944Z","iopub.execute_input":"2021-10-05T10:27:27.484289Z","iopub.status.idle":"2021-10-05T10:27:27.532490Z","shell.execute_reply.started":"2021-10-05T10:27:27.484257Z","shell.execute_reply":"2021-10-05T10:27:27.531738Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"assert sum([len(x) for x in val_indices]) == df_train.shape[0]","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.533646Z","iopub.execute_input":"2021-10-05T10:27:27.533960Z","iopub.status.idle":"2021-10-05T10:27:27.540339Z","shell.execute_reply.started":"2021-10-05T10:27:27.533927Z","shell.execute_reply":"2021-10-05T10:27:27.539409Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_indices","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.543520Z","iopub.execute_input":"2021-10-05T10:27:27.543788Z","iopub.status.idle":"2021-10-05T10:27:27.553539Z","shell.execute_reply.started":"2021-10-05T10:27:27.543746Z","shell.execute_reply":"2021-10-05T10:27:27.552580Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.554967Z","iopub.execute_input":"2021-10-05T10:27:27.555306Z","iopub.status.idle":"2021-10-05T10:27:27.569984Z","shell.execute_reply.started":"2021-10-05T10:27:27.555274Z","shell.execute_reply":"2021-10-05T10:27:27.569035Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Model and training classes","metadata":{}},{"cell_type":"code","source":"DEVICE = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(DEVICE)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.571492Z","iopub.execute_input":"2021-10-05T10:27:27.571924Z","iopub.status.idle":"2021-10-05T10:27:27.578413Z","shell.execute_reply.started":"2021-10-05T10:27:27.571891Z","shell.execute_reply":"2021-10-05T10:27:27.577231Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Dataset(torch_data.Dataset):\n    def __init__(self, paths, targets=None, mri_type=None, label_smoothing=0.01, split=\"train\", augment=False):\n        self.paths = paths\n        self.targets = targets\n        self.mri_type = mri_type\n        self.label_smoothing = label_smoothing\n        self.split = split\n        self.augment = augment\n          \n    def __len__(self):\n        return len(self.paths)\n    \n    def __getitem__(self, index):\n        path = self.paths[index]\n        data = load_normalized_images_3d(path)        \n        if self.targets is None:\n            return {\"X\": torch.tensor(data).float(), \"id\": scan_id}\n        else:\n            y = torch.tensor(abs(self.targets[index]-self.label_smoothing), dtype=torch.float)\n            return {\"X\": torch.tensor(data).float(), \"y\": y}\n","metadata":{"papermill":{"duration":0.634322,"end_time":"2021-07-14T20:27:50.594701","exception":false,"start_time":"2021-07-14T20:27:49.960379","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:27.581472Z","iopub.execute_input":"2021-10-05T10:27:27.581972Z","iopub.status.idle":"2021-10-05T10:27:27.593592Z","shell.execute_reply.started":"2021-10-05T10:27:27.581942Z","shell.execute_reply":"2021-10-05T10:27:27.592686Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from models import resnet\nfrom collections import OrderedDict\n\nclass MedicalNet(nn.Module):\n\n    def __init__(self, path_to_weights=None, device='cuda'):\n        super(MedicalNet, self).__init__()\n        self.model = resnet.resnet34(sample_input_D=1, sample_input_H=256, sample_input_W=256, num_seg_classes=2)        \n        if path_to_weights:\n            net_dict = self.model.state_dict()\n            pretrained_weights = torch.load(path_to_weights, map_location=torch.device(device))\n            pretrain_dict = {\n                k.replace(\"module.\", \"\"): v for k, v in pretrained_weights['state_dict'].items() if k.replace(\"module.\", \"\") in net_dict.keys()\n              }\n            net_dict.update(pretrain_dict)\n            self.model.load_state_dict(net_dict)\n        self.model.conv_seg = nn.Sequential(OrderedDict([ \n            ('adapt1', nn.AdaptiveMaxPool3d(output_size=(1, 1, 1))),\n            ('flatten1', nn.Flatten(start_dim=1)),\n            ('dropout1', nn.Dropout(0.1))\n        ]))\n        self.model = self.model.to(device)\n        self.fc = nn.Linear(512, 1).to(device)\n\n    def forward(self, x):\n        features = self.model(x)\n        return self.fc(features)\n\nclass EfficientNetModel(nn.Module):\n    def __init__(self, device='cuda'):\n        super().__init__()\n        self.net = EfficientNet3D.from_name(\"efficientnet-b0\", override_params={'num_classes': 2}, in_channels=1)\n        n_features = self.net._fc.in_features\n        self.net._fc = nn.Linear(in_features=n_features, out_features=1, bias=True)\n        self.net = self.net.to(device)\n    \n    def forward(self, x):\n        out = self.net(x)\n        return out\n    \n\n    ","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.595111Z","iopub.execute_input":"2021-10-05T10:27:27.595751Z","iopub.status.idle":"2021-10-05T10:27:27.608986Z","shell.execute_reply.started":"2021-10-05T10:27:27.595711Z","shell.execute_reply":"2021-10-05T10:27:27.608037Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Trainer:\n    def __init__(\n        self, \n        model, \n        device, \n        optimizer, \n        criterion\n    ):\n        self.model = model\n        self.device = device\n        self.optimizer = optimizer\n        self.criterion = criterion\n        self.scheduler = torch.optim.lr_scheduler.ReduceLROnPlateau(optimizer, 'min', patience=5)\n        self.best_valid_score = np.inf\n        self.n_patience = 0\n        self.lastmodel = None\n        \n    def fit(self, epochs, train_loader, valid_loader, save_path, patience, kfold):        \n        for n_epoch in range(1, epochs + 1):\n            self.info_message(\"EPOCH: {}\", n_epoch)\n            \n            train_loss, train_time = self.train_epoch(train_loader)\n            valid_loss, valid_auc, valid_time = self.valid_epoch(valid_loader)\n            \n            self.info_message(\n                \"[Epoch Train: {}] loss: {:.4f}, time: {:.2f} s            \",\n                n_epoch, train_loss, train_time\n            )\n            \n            self.info_message(\n                \"[Epoch Valid: {}] loss: {:.4f}, auc: {:.4f}, time: {:.2f} s\",\n                n_epoch, valid_loss, valid_auc, valid_time\n            )\n\n            self.scheduler.step(valid_loss)\n            # if True:\n            # if self.best_valid_score < valid_auc: \n            if self.best_valid_score > valid_loss: \n                self.save_model(n_epoch, save_path, valid_loss, valid_auc, kfold)\n                self.info_message(\n                     \"auc improved from {:.4f} to {:.4f}. Saved model to '{}'\", \n                    self.best_valid_score, valid_loss, self.lastmodel\n                )\n                self.best_valid_score = valid_loss\n                self.n_patience = 0\n            else:\n                self.n_patience += 1\n            \n            \n    def train_epoch(self, train_loader):\n        self.model.train()\n        t = time.time()\n        sum_loss = 0\n\n        for step, batch in enumerate(train_loader, 1):\n            X = batch[\"X\"].to(self.device)\n            targets = batch[\"y\"].to(self.device)\n            self.optimizer.zero_grad()\n            outputs = self.model(X).squeeze(1)\n            \n            loss = self.criterion(outputs, targets)\n            loss.backward()\n\n            sum_loss += loss.detach().item()\n\n            self.optimizer.step()\n            message = 'Train Step {}/{}, train_loss: {:.4f}'\n            self.info_message(message, step, len(train_loader), sum_loss/step, end=\"\\r\")\n        \n        return sum_loss/len(train_loader), int(time.time() - t)\n    \n    def valid_epoch(self, valid_loader):\n        self.model.eval()\n        t = time.time()\n        sum_loss = 0\n        y_all = []\n        outputs_all = []\n\n        for step, batch in enumerate(valid_loader, 1):\n            with torch.no_grad():\n                X = batch[\"X\"].to(self.device)\n                targets = batch[\"y\"].to(self.device)\n\n                outputs = self.model(X).squeeze(1)\n                loss = self.criterion(outputs, targets)\n\n                sum_loss += loss.detach().item()\n                y_all.extend(batch[\"y\"].tolist())\n                outputs_all.extend(outputs.tolist())\n\n            message = 'Valid Step {}/{}, valid_loss: {:.4f}'\n            self.info_message(message, step, len(valid_loader), sum_loss/step, end=\"\\r\")\n            \n        y_all = [1 if x > 0.5 else 0 for x in y_all]\n        auc = roc_auc_score(y_all, outputs_all)\n        \n        return sum_loss/len(valid_loader), auc, int(time.time() - t)\n    \n    def save_model(self, n_epoch, save_path, loss, auc, kfold):\n        self.lastmodel = f\"{save_path}-k{kfold}-e{n_epoch}-loss{loss:.3f}-auc{auc:.3f}.pth\"\n        torch.save(\n            {\n                \"model_state_dict\": self.model.state_dict(),\n                \"optimizer_state_dict\": self.optimizer.state_dict(),\n                \"best_valid_score\": self.best_valid_score,\n                \"n_epoch\": n_epoch,\n            },\n            self.lastmodel,\n        )\n    \n    @staticmethod\n    def info_message(message, *args, end=\"\\n\"):\n        print(message.format(*args), end=end)","metadata":{"papermill":{"duration":0.637077,"end_time":"2021-07-14T20:27:58.09407","exception":false,"start_time":"2021-07-14T20:27:57.456993","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:27.611033Z","iopub.execute_input":"2021-10-05T10:27:27.611683Z","iopub.status.idle":"2021-10-05T10:27:27.632891Z","shell.execute_reply.started":"2021-10-05T10:27:27.611621Z","shell.execute_reply":"2021-10-05T10:27:27.632094Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## train models","metadata":{}},{"cell_type":"code","source":"# Select model\n\nMODEL_TYPE = 'medicalnet' \nFEATURES_EXTRACTION = False\n\ndef create_model(device):\n    if MODEL_TYPE == 'medicalnet':\n        model = MedicalNet(path_to_weights=PRETRAINED_WEIGHTS, device=device)\n        for name,param in model.named_parameters():\n            if name.startswith(\"fc\") or not FEATURES_EXTRACTION:\n                param.requires_grad = True\n            else:\n                param.requires_grad = False\n\n        print(\"Params to learn:\")\n        params_to_update = []\n        for name,param in model.named_parameters():\n            if param.requires_grad == True:\n                params_to_update.append(param)\n\n    elif MODEL_TYPE == 'efficientnet':\n        model = EfficientNetModel(device=device)\n        params_to_update = model.parameters()\n    else:\n        model = DenseNet121(spatial_dims=3, in_channels=1, out_channels=1)\n        model = model.to(device)\n        params_to_update = model.parameters()\n    return model, params_to_update","metadata":{"execution":{"iopub.status.busy":"2021-10-05T10:27:27.634462Z","iopub.execute_input":"2021-10-05T10:27:27.634836Z","iopub.status.idle":"2021-10-05T10:27:27.645593Z","shell.execute_reply.started":"2021-10-05T10:27:27.634801Z","shell.execute_reply":"2021-10-05T10:27:27.644726Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def train_mri_type(df_train, df_valid, mri_type, kfold, device):\n    print(f'kfold: {kfold}')    \n    train_data_retriever = Dataset(\n        paths=df_train.voxel_paths.values, \n        targets=df_train[\"MGMT_value\"].values, \n        mri_type=df_train[\"MRI_Type\"].values,\n        augment=True\n    )\n\n    valid_data_retriever = Dataset(\n        paths=df_valid.voxel_paths.values, \n        targets=df_valid[\"MGMT_value\"].values, \n        mri_type=df_valid[\"MRI_Type\"].values,\n    )\n\n    train_loader = torch_data.DataLoader(\n        train_data_retriever,\n        batch_size=4,\n        shuffle=True,\n        num_workers=8,\n    )\n\n    valid_loader = torch_data.DataLoader(\n        valid_data_retriever, \n        batch_size=4,\n        shuffle=False,\n        num_workers=8,\n    )\n\n# Observe that all parameters are being optimized\n    model, params_to_update = create_model(device=device)\n    optimizer = torch.optim.Adam(params_to_update, lr=LEARNING_RATE)\n    #optimizer = torch.optim.SGD(model.parameters(), lr=0.001, momentum=0.9)\n\n    criterion = torch_functional.binary_cross_entropy_with_logits\n\n    trainer = Trainer(\n        model, \n        device, \n        optimizer, \n        criterion\n    )\n\n    print(f'Fitting model for mri_type: {mri_type}')\n    history = trainer.fit(\n        EPOCHS, \n        train_loader, \n        valid_loader, \n        f\"{mri_type}\", \n        8,\n        kfold=kfold\n    )\n    \n    return trainer.lastmodel\n\nmodelfiles = None\n\n# all mri types\nif not modelfiles:\n    modelfiles = [train_mri_type(df_train.loc[train, :], df_train.loc[val, :], mri_type='all', kfold=i, device=DEVICE) \n                  for i, (train, val) in enumerate(zip(train_indices, val_indices))]\n    print(modelfiles)","metadata":{"lines_to_next_cell":2,"papermill":{"duration":447.387602,"end_time":"2021-07-14T20:35:26.110421","exception":false,"start_time":"2021-07-14T20:27:58.722819","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T10:27:27.647023Z","iopub.execute_input":"2021-10-05T10:27:27.647421Z","iopub.status.idle":"2021-10-05T17:01:59.180591Z","shell.execute_reply.started":"2021-10-05T10:27:27.647385Z","shell.execute_reply":"2021-10-05T17:01:59.178801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Predict function","metadata":{}},{"cell_type":"code","source":"def predict(modelfile, df, mri_type, split, device):\n    print(\"Predict:\", modelfile, mri_type, df.shape)        \n    data_retriever = Dataset(\n        df.voxel_paths.values, \n        mri_type=df[\"MRI_Type\"].values,\n        split=split\n    )\n        \n    data_loader = torch_data.DataLoader(\n        data_retriever,\n        batch_size=4,\n        shuffle=False,\n        num_workers=8,\n    )\n    \n    model, _ = create_model(device=device)\n    checkpoint = torch.load(modelfile, map_location=device)\n    model.load_state_dict(checkpoint[\"model_state_dict\"])\n    model.eval()\n    \n    y_pred = []\n    ids = []\n\n    for e, batch in enumerate(data_loader,1):\n        print(f\"{e}/{len(data_loader)}\", end=\"\\r\")\n        with torch.no_grad():\n            tmp_pred = torch.sigmoid(model(batch[\"X\"].float().to(device))).cpu().numpy().squeeze()\n            if tmp_pred.size == 1:\n                y_pred.append(tmp_pred)\n            else:\n                y_pred.extend(tmp_pred.tolist())\n            ids.extend(batch[\"id\"].numpy().tolist())\n            \n    preddf = pd.DataFrame({\"BraTS21ID\": ids, \"MGMT_value\": y_pred}) \n    preddf = preddf.set_index(\"BraTS21ID\")\n    return preddf","metadata":{"papermill":{"duration":0.990911,"end_time":"2021-07-14T20:35:30.482254","exception":false,"start_time":"2021-07-14T20:35:29.491343","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T17:01:59.182409Z","iopub.execute_input":"2021-10-05T17:01:59.182778Z","iopub.status.idle":"2021-10-05T17:01:59.194411Z","shell.execute_reply.started":"2021-10-05T17:01:59.182737Z","shell.execute_reply":"2021-10-05T17:01:59.193587Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Ensemble for validation","metadata":{}},{"cell_type":"code","source":"#modelfiles = ['../input/fake-model-flair/FLAIR-e1-loss0.707-auc0.570.pth'] * FOLDS # TODO: remove it!","metadata":{"execution":{"iopub.status.busy":"2021-10-05T17:01:59.198204Z","iopub.execute_input":"2021-10-05T17:01:59.198466Z","iopub.status.idle":"2021-10-05T17:01:59.231822Z","shell.execute_reply.started":"2021-10-05T17:01:59.198440Z","shell.execute_reply":"2021-10-05T17:01:59.230829Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# K-FOLD CV\npredictions_by_fold = []\nfor m, indices in zip(modelfiles, val_indices):\n    df_valid = df_train.loc[indices, :].set_index(\"BraTS21ID\")\n    df_valid[\"MGMT_pred\"] = 0 # add new colum\n    pred = predict(m, df_valid, mri_type='all', split=\"train\", device=DEVICE)\n    pred = pred.groupby(pred.index).mean() # give a single vote for the same patient study using all mri types\n    df_valid[\"MGMT_pred\"] = pred[\"MGMT_value\"]\n    predictions_by_fold.append(df_valid.copy())","metadata":{"execution":{"iopub.status.busy":"2021-10-05T17:01:59.233014Z","iopub.execute_input":"2021-10-05T17:01:59.234218Z","iopub.status.idle":"2021-10-05T17:02:04.575419Z","shell.execute_reply.started":"2021-10-05T17:01:59.234191Z","shell.execute_reply":"2021-10-05T17:02:04.572656Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_valid = pd.concat(predictions_by_fold)\nauc = roc_auc_score(df_valid[\"MGMT_value\"], df_valid[\"MGMT_pred\"])\nprint(f\"Validation ensemble AUC: {auc:.4f}\")\nsns.displot(df_valid[\"MGMT_pred\"])","metadata":{"execution":{"iopub.status.busy":"2021-10-05T17:02:04.576871Z","iopub.status.idle":"2021-10-05T17:02:04.577288Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Ensemble for submission","metadata":{}},{"cell_type":"code","source":"# TODO: train the model in the whole dataset\n#modelfiles = train_mri_type(df_train, df_train, mri_type='all') \n#                  for (train, val) in zip(train_indices, val_indices)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T17:02:04.578305Z","iopub.status.idle":"2021-10-05T17:02:04.578874Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission = pd.read_csv(f\"{data_directory}/sample_submission.csv\", index_col=\"BraTS21ID\")\nsubmission[\"MGMT_value\"] = 0\n\nfor m in modelfiles:\n    pred = predict(m, submission, 'all', split=\"test\", device=DEVICE)\n    pred = pred.groupby(pred.index).mean()\n    submission[\"MGMT_value\"] += pred[\"MGMT_value\"]\n    \nsubmission['MGMT_value'] /= len(modelfiles)\nsubmission[\"MGMT_value\"].to_csv(\"submission.csv\")","metadata":{"papermill":{"duration":0.990911,"end_time":"2021-07-14T20:35:30.482254","exception":false,"start_time":"2021-07-14T20:35:29.491343","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T17:02:04.580046Z","iopub.status.idle":"2021-10-05T17:02:04.580718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission","metadata":{"papermill":{"duration":0.972383,"end_time":"2021-07-14T20:35:41.939464","exception":false,"start_time":"2021-07-14T20:35:40.967081","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2021-10-05T17:02:04.581866Z","iopub.status.idle":"2021-10-05T17:02:04.582683Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.displot(submission[\"MGMT_value\"])","metadata":{"execution":{"iopub.status.busy":"2021-10-05T17:02:04.583962Z","iopub.status.idle":"2021-10-05T17:02:04.584764Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}