{"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":"import os\nimport sys \nimport json\nimport glob\nimport random\nimport collections\nimport time\nimport re\nimport math\nimport numpy as np\nimport pandas as pd\nimport cv2\nfrom tqdm import tqdm\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport pydicom\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\n\nfrom random import shuffle\nfrom sklearn import model_selection as sk_model_selection\nfrom sklearn.metrics import roc_auc_score\n\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nimport torch.optim as optim\nfrom torch.optim import lr_scheduler\nfrom torch.utils.data import DataLoader","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2021-10-05T15:21:40.019239Z","iopub.execute_input":"2021-10-05T15:21:40.019774Z","iopub.status.idle":"2021-10-05T15:21:40.027059Z","shell.execute_reply.started":"2021-10-05T15:21:40.019738Z","shell.execute_reply":"2021-10-05T15:21:40.026285Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data_directory = '../input/rsna-miccai-brain-tumor-radiogenomic-classification'\n \nmri_types = ['FLAIR','T1w','T1wCE','T2w']\nmri_types_id=0 # 0,1,2,3\n\nIMAGE_SIZE = 256\nNUM_IMAGES = 100\nBATCH_SIZE= 16\nNUM_EPOCHS = 50\nDEVICE = \"cuda:0\"\n\ntrain_df = pd.read_csv(\"../input/rsna-miccai-brain-tumor-radiogenomic-classification/train_labels.csv\")\ntrain_df['BraTS21ID5'] = [format(x, '05d') for x in train_df.BraTS21ID]\ntrain_df.head(5)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.120529Z","iopub.execute_input":"2021-10-05T15:21:40.121145Z","iopub.status.idle":"2021-10-05T15:21:40.140480Z","shell.execute_reply.started":"2021-10-05T15:21:40.121101Z","shell.execute_reply":"2021-10-05T15:21:40.139854Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"indexs = [not val in [\"00109\", \"00123\", \"00709\"] for val in train_df[\"BraTS21ID5\"]]\ntrain_df = train_df[indexs]","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.143811Z","iopub.execute_input":"2021-10-05T15:21:40.144357Z","iopub.status.idle":"2021-10-05T15:21:40.149183Z","shell.execute_reply.started":"2021-10-05T15:21:40.144328Z","shell.execute_reply":"2021-10-05T15:21:40.148322Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def load_dicom_image(path, img_size=IMAGE_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 load_dicom_images_3d(scan_id, num_imgs=NUM_IMAGES, img_size=IMAGE_SIZE, mri_type=mri_types[mri_types_id], split=\"train\", rotate=0):\n#     rotate = np.random.randint(4)\n\n    files = sorted(glob.glob(f\"{data_directory}/{split}/{scan_id}/{mri_type}/*.dcm\"), \n               key=lambda var:[int(x) if x.isdigit() else x for x in re.findall(r'[^0-9]|[0-9]+', var)])\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 img3d.shape[-1] < num_imgs:\n        s1 = (num_imgs - img3d.shape[-1]) // 2\n        s2 = (num_imgs - img3d.shape[-1] + 1) // 2\n        n_zero1 = np.zeros((img_size, img_size, s1))\n        n_zero2 = np.zeros((img_size, img_size, s2))\n        img3d = np.concatenate((n_zero1, img3d,  n_zero2), axis = -1)\n        \n    if np.min(img3d) < np.max(img3d):\n        img3d = img3d - np.min(img3d)\n        img3d = img3d / np.max(img3d)\n            \n    return np.expand_dims(img3d,0)\n\n# a = load_dicom_images_3d(\"00003\")\n# print(a.shape)\n# print(np.min(a), np.max(a), np.mean(a), np.median(a))\n# image = a[0]\n# print(\"Dimension of the CT scan is:\", image.shape)\n# plt.imshow(np.squeeze(image[:, :, 100]), cmap=\"gray\")","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.246777Z","iopub.execute_input":"2021-10-05T15:21:40.247005Z","iopub.status.idle":"2021-10-05T15:21:40.260086Z","shell.execute_reply.started":"2021-10-05T15:21:40.246966Z","shell.execute_reply":"2021-10-05T15:21:40.259375Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class Dataset(object):\n    def __init__(self,df,is_train=True, device=\"cpu\"):\n#         self.idx = df[\"BraTS21ID\"].values\n        self.paths = df[\"BraTS21ID5\"].values\n        self.y =  df[\"MGMT_value\"].values\n        self.is_train = is_train\n        self.len = df.shape[0]\n        self.device = device\n    def __len__(self):\n        return self.len\n   \n    def __getitem__(self,ids):\n        \n        if self.is_train:\n            X = load_dicom_images_3d(self.paths[ids],split=\"train\")\n            Y = self.y[ids]\n            X = torch.tensor(X, dtype=torch.float32, device=self.device)\n            Y = torch.tensor(Y, dtype=torch.float32, device=self.device)\n            return X, Y\n        else:\n            X = load_dicom_images_3d(self.paths[ids],split=\"train\")\n            X = torch.tensor(X, dtype=torch.float32, device=self.device)\n            return X","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.362604Z","iopub.execute_input":"2021-10-05T15:21:40.363093Z","iopub.status.idle":"2021-10-05T15:21:40.370457Z","shell.execute_reply.started":"2021-10-05T15:21:40.363065Z","shell.execute_reply":"2021-10-05T15:21:40.369628Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class View(nn.Module):\n    def __init__(self, size):\n        super(View, self).__init__()\n        self.size = size\n\n    def forward(self, tensor):\n        return tensor.view(self.size)\n    \nclass Model(nn.Module):\n    def __init__(self):\n        super(Model, self).__init__()\n        h_dim = [32, 64, 128]\n        self.net = nn.Sequential(\n            nn.Conv3d(1, h_dim[0], 4, 2, 1), # 128, 128, 50\n            nn.BatchNorm3d(h_dim[0]),\n            nn.Dropout(0.05),\n            nn.MaxPool3d((2,2,2)),  # 64, 64,50\n            nn.Conv3d(h_dim[0], h_dim[1], 4, 2, 1), # 32, 32,12\n            nn.BatchNorm3d(h_dim[1]),\n            nn.Dropout(0.05),\n            nn.MaxPool3d((2,2,2)),   # 16, 16,12\n            nn.Conv3d(h_dim[1], h_dim[2], 4, 2, 1),  # 8, 8, 3\n            nn.BatchNorm3d(h_dim[2]),\n            nn.Dropout(0.05),\n            nn.MaxPool3d((4,4,1)), # 2, 2, 3\n            View((-1, h_dim[-1]*2*2*3)),\n            nn.Linear(h_dim[-1]*2*2*3, 128),\n            nn.ReLU(True),\n            nn.Linear(128, 1))\n    \n    def forward(self, x):\n        #x: [bz, 1, 256, 256, 400]\n        x = self.net(x)\n        x = torch.sigmoid(x)\n        eps = 1e-6\n        x = 0.5*eps + (1-eps)*x\n        return x","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.462571Z","iopub.execute_input":"2021-10-05T15:21:40.462763Z","iopub.status.idle":"2021-10-05T15:21:40.474236Z","shell.execute_reply.started":"2021-10-05T15:21:40.462741Z","shell.execute_reply":"2021-10-05T15:21:40.473429Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_loss(y_pred, y):\n    y = y.unsqueeze(-1)\n    loss = - torch.log(y_pred) * y - torch.log(1 - y_pred) * (1 - y)\n#     print(loss.shape)\n    acc = (y_pred > 0.5) == y\n    acc = acc.to(torch.float32)\n#     print(y_pred)\n    auc = roc_auc_score(y.detach().cpu().numpy(), y_pred.detach().cpu().numpy())\n#     print(acc.shape)\n    return loss.mean(), acc.mean(), auc","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.475549Z","iopub.execute_input":"2021-10-05T15:21:40.476071Z","iopub.status.idle":"2021-10-05T15:21:40.485534Z","shell.execute_reply.started":"2021-10-05T15:21:40.476034Z","shell.execute_reply":"2021-10-05T15:21:40.484845Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_dataset = Dataset(train_df[:int(train_df.shape[0]*0.85)], device=DEVICE)\nval_dataset = Dataset(train_df[int(train_df.shape[0]*0.85):], device=DEVICE)\n\ntrain_loader = DataLoader(train_dataset, BATCH_SIZE, shuffle=True)\nval_loader = DataLoader(val_dataset, BATCH_SIZE, shuffle=False)\n\nmodel = Model()\nmodel = model.to(DEVICE)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.563733Z","iopub.execute_input":"2021-10-05T15:21:40.563922Z","iopub.status.idle":"2021-10-05T15:21:40.582719Z","shell.execute_reply.started":"2021-10-05T15:21:40.563900Z","shell.execute_reply":"2021-10-05T15:21:40.582100Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def train(train_dataloader, val_dataloader, device, num_epochs):\n    optimizer = optim.Adam(model.parameters(), lr=1e-4)\n    scheduler = lr_scheduler.StepLR(optimizer, step_size=5, gamma=0.75)\n    # record training\n    best_model_path = \"best_model.pt\"\n    best_val_loss = np.inf\n    torch.save(model, best_model_path)\n    \n    for epoch in tqdm(range(num_epochs)):\n        model.train()\n        train_loss = torch.tensor(0, device=device).float()\n\n        for X, y in train_dataloader:\n            y_pred = model(X)\n            loss, acc, auc = get_loss(y_pred, y)\n#                 print(y.shape)\n\n            optimizer.zero_grad()\n            loss.backward()\n            optimizer.step()\n            train_loss += loss * X.shape[0]\n\n            print(\"loss: %.4f\" % loss.item(), \"acc: %.4f\" % acc.item(), \"auc : %.2f\" % auc)\n            \n\n        model.eval()\n        val_loss = torch.tensor(0, device=device).float()\n        val_acc = torch.tensor(0, device=device).float()\n        val_auc = torch.tensor(0, device=device).float()\n        \n        with torch.no_grad():\n            for X, y in val_dataloader:\n                y_pred = model(X)\n                loss, acc, auc = get_loss(y_pred, y)\n                val_loss += loss * X.shape[0]\n                val_acc += acc * X.shape[0]\n                val_auc += auc * X.shape[0]\n            \n            val_data_num = len(val_dataloader.dataset)\n            print(\"Current val loss: %.4f\" % (val_loss/val_data_num))\n            print(\"Current acc: %.4f\" % (val_acc/val_data_num))\n            print(\"Current auc: %.4f\" % (val_auc/val_data_num))\n\n            if best_val_loss > val_loss/val_data_num:\n                best_val_loss = val_loss/val_data_num\n                print(\"Impove! saving the model. Current best val loss: %.4f\" % best_val_loss)\n                torch.save(model, best_model_path)\n        scheduler.step()","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.653200Z","iopub.execute_input":"2021-10-05T15:21:40.653390Z","iopub.status.idle":"2021-10-05T15:21:40.664239Z","shell.execute_reply.started":"2021-10-05T15:21:40.653368Z","shell.execute_reply":"2021-10-05T15:21:40.663388Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 用以下命令跑后台将输出存log\n# import sys\n# stdout_backup = sys.stdout\n# log_file = open(\"message.log\", \"w\")\n# sys.stdout = log_file\ntrain(train_loader, val_loader, DEVICE, num_epochs=NUM_EPOCHS)","metadata":{"execution":{"iopub.status.busy":"2021-10-05T15:21:40.665870Z","iopub.execute_input":"2021-10-05T15:21:40.666425Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}