{"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":"# Simple EDA","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\nimport matplotlib.pyplot as plt\nimport matplotlib\nmatplotlib.rcParams['animation.html'] = 'jshtml'\nimport seaborn as sns\n\nimport pydicom\nfrom skimage.transform import resize\n\nimport os\nfrom sklearn.decomposition import IncrementalPCA\nfrom sklearn.ensemble import GradientBoostingClassifier\nfrom sklearn.metrics import auc, roc_curve","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:11.68152Z","iopub.execute_input":"2021-09-21T04:22:11.681936Z","iopub.status.idle":"2021-09-21T04:22:11.688162Z","shell.execute_reply.started":"2021-09-21T04:22:11.681902Z","shell.execute_reply":"2021-09-21T04:22:11.687221Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## data loading stuff","metadata":{}},{"cell_type":"code","source":"root_dir = '../input/rsna-miccai-brain-tumor-radiogenomic-classification/'\ndf = pd.read_csv(root_dir+'train_labels.csv')","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:12.273834Z","iopub.execute_input":"2021-09-21T04:22:12.274174Z","iopub.status.idle":"2021-09-21T04:22:12.293065Z","shell.execute_reply.started":"2021-09-21T04:22:12.274146Z","shell.execute_reply":"2021-09-21T04:22:12.292323Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Add the full paths for each id for different types of sequences to the csv \ndef full_ids(data):\n    zeros = 5 - len(str(data))\n    if zeros > 0:\n        prefix = ''.join(['0' for i in range(zeros)])\n    \n    return prefix+str(data)\n        \n\ndf['BraTS21ID_full'] = df['BraTS21ID'].apply(full_ids)\n\n# Add all the paths to the df for easy access\ndf['flair'] = df['BraTS21ID_full'].apply(lambda file_id : root_dir+'train/'+file_id+'/FLAIR/')\ndf['t1w'] = df['BraTS21ID_full'].apply(lambda file_id : root_dir+'train/'+file_id+'/T1w/')\ndf['t1wce'] = df['BraTS21ID_full'].apply(lambda file_id : root_dir+'train/'+file_id+'/T1wCE/')\ndf['t2w'] = df['BraTS21ID_full'].apply(lambda file_id : root_dir+'train/'+file_id+'/T2w/')\ndf\n\n# first column is ID, second is the label, third is another ID, the rest are paths to dicoms","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:12.601354Z","iopub.execute_input":"2021-09-21T04:22:12.601698Z","iopub.status.idle":"2021-09-21T04:22:12.635559Z","shell.execute_reply.started":"2021-09-21T04:22:12.601669Z","shell.execute_reply":"2021-09-21T04:22:12.634578Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df['MGMT_value']","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:12.805332Z","iopub.execute_input":"2021-09-21T04:22:12.805699Z","iopub.status.idle":"2021-09-21T04:22:12.813264Z","shell.execute_reply.started":"2021-09-21T04:22:12.805663Z","shell.execute_reply":"2021-09-21T04:22:12.8126Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## function to load an image","metadata":{}},{"cell_type":"code","source":"def get_image(img_path):\n    '''\n    Returns the image data as a numpy array.\n    '''\n    data = pydicom.dcmread(img_path)\n    img = data.pixel_array\n#     if np.max(data.pixel_array)==0:\n#         img = data.pixel_array\n#     else:\n#         img = data.pixel_array/np.max(data.pixel_array)\n    \n    img = resize(img, (128, 128)) # resize it to a more compact size\n        \n    return img\n\nimg = get_image('../input/rsna-miccai-brain-tumor-radiogenomic-classification/train/00000/T1w/Image-24.dcm')\nplt.imshow(img, cmap='gray')","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:13.327127Z","iopub.execute_input":"2021-09-21T04:22:13.327721Z","iopub.status.idle":"2021-09-21T04:22:13.556244Z","shell.execute_reply.started":"2021-09-21T04:22:13.327686Z","shell.execute_reply":"2021-09-21T04:22:13.555552Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Make a dataframe for images around the midpoint for an MRI type","metadata":{}},{"cell_type":"code","source":"def get_mri_img_df(mri_type):\n\n    all_img_files = []\n    all_img_labels = []\n    all_img_patient_ids = []\n    \n    for row in df.iterrows():\n        if row[1]['BraTS21ID_full'] == '00109' and mri_type == 'flair':\n            continue\n        if row[1]['BraTS21ID_full'] == '00123' and mri_type == 't1w':\n            continue\n        if row[1]['BraTS21ID_full'] == '00709' and mri_type == 'flair':\n            continue\n        img_dir = row[1][mri_type]\n        img_files = os.listdir(img_dir)\n        img_nums = sorted([int(ele.replace('Image-', '').replace('.dcm', '')) for ele in img_files])\n        mid_point = int(len(img_nums)/2)\n        start_point = max(mid_point-5, 0)\n        end_point = min(mid_point+5, len(img_nums))\n        img_names = [f'Image-{img_nums[i]}.dcm' for i in range(start_point, end_point)]\n        img_paths = [img_dir+ele for ele in img_names]\n        img_labels = [row[1]['MGMT_value']]*len(img_paths)\n        img_patient_ids = [row[1]['BraTS21ID']]*len(img_paths)\n        all_img_files.extend(img_paths)\n        all_img_labels.extend(img_labels)\n        all_img_patient_ids.extend(img_patient_ids)\n\n    mri_img_df = pd.DataFrame({'patient_ids': all_img_patient_ids,\n                  'labels': all_img_labels,\n                  'file_paths': all_img_files})\n    \n    return mri_img_df","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:14.373033Z","iopub.execute_input":"2021-09-21T04:22:14.373653Z","iopub.status.idle":"2021-09-21T04:22:14.384824Z","shell.execute_reply.started":"2021-09-21T04:22:14.373602Z","shell.execute_reply":"2021-09-21T04:22:14.384109Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## First let's take a look at the MRI images","metadata":{}},{"cell_type":"code","source":"img_df = get_mri_img_df('t2w')","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:15.421355Z","iopub.execute_input":"2021-09-21T04:22:15.421815Z","iopub.status.idle":"2021-09-21T04:22:16.575546Z","shell.execute_reply.started":"2021-09-21T04:22:15.421785Z","shell.execute_reply":"2021-09-21T04:22:16.57457Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"len(img_df)","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:22:17.177814Z","iopub.execute_input":"2021-09-21T04:22:17.178288Z","iopub.status.idle":"2021-09-21T04:22:17.183572Z","shell.execute_reply.started":"2021-09-21T04:22:17.178258Z","shell.execute_reply":"2021-09-21T04:22:17.182533Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# lets load in all images\nall_img_array = np.zeros((len(img_df), 128*128))\nfor i, ele in enumerate(img_df['file_paths']):\n    if i%100 == 0:\n        print(i)\n    a = get_image(ele)\n    a = a.flatten()\n    all_img_array[i] = a","metadata":{"scrolled":true,"execution":{"iopub.status.busy":"2021-09-21T04:22:17.449968Z","iopub.execute_input":"2021-09-21T04:22:17.450332Z","iopub.status.idle":"2021-09-21T04:23:49.645611Z","shell.execute_reply.started":"2021-09-21T04:22:17.45029Z","shell.execute_reply":"2021-09-21T04:23:49.644467Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(df['MGMT_value'].values[:20])","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:23:49.651435Z","iopub.execute_input":"2021-09-21T04:23:49.654106Z","iopub.status.idle":"2021-09-21T04:23:49.662561Z","shell.execute_reply.started":"2021-09-21T04:23:49.654053Z","shell.execute_reply":"2021-09-21T04:23:49.661475Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"img_df","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:23:49.664843Z","iopub.execute_input":"2021-09-21T04:23:49.665265Z","iopub.status.idle":"2021-09-21T04:23:49.688981Z","shell.execute_reply.started":"2021-09-21T04:23:49.665214Z","shell.execute_reply":"2021-09-21T04:23:49.688002Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"imgs = all_img_array[img_df['labels'] == 1][:100].reshape((100, 128, 128))\n\n_, axs = plt.subplots(10, 10, figsize=(30, 30))\naxs = axs.flatten()\n\nfor img, ax in zip(imgs, axs):\n    ax.imshow(img)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-09-21T03:30:19.596743Z","iopub.execute_input":"2021-09-21T03:30:19.597129Z","iopub.status.idle":"2021-09-21T03:30:30.976957Z","shell.execute_reply.started":"2021-09-21T03:30:19.597096Z","shell.execute_reply":"2021-09-21T03:30:30.975684Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"imgs = all_img_array[img_df['labels'] == 0][:100].reshape((100, 128, 128))\n\n_, axs = plt.subplots(10, 10, figsize=(30, 30))\naxs = axs.flatten()\n\nfor img, ax in zip(imgs, axs):\n    ax.imshow(img)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-09-21T03:31:00.290429Z","iopub.execute_input":"2021-09-21T03:31:00.290798Z","iopub.status.idle":"2021-09-21T03:31:11.222345Z","shell.execute_reply.started":"2021-09-21T03:31:00.290764Z","shell.execute_reply":"2021-09-21T03:31:11.221227Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Run PCA in batches\nbatch_size = 500\nnc = 200\ntransformer = IncrementalPCA(n_components=nc)\n\nfor i in range(0, len(all_img_arrays), batch_size):\n    print(i)\n    batch = np.array(all_img_arrays[i:i+batch_size])\n    if len(batch) > nc:\n        transformer.partial_fit(batch)","metadata":{"scrolled":true,"execution":{"iopub.status.busy":"2021-09-20T18:38:19.920544Z","iopub.execute_input":"2021-09-20T18:38:19.921033Z","iopub.status.idle":"2021-09-20T18:40:18.678618Z","shell.execute_reply.started":"2021-09-20T18:38:19.920997Z","shell.execute_reply":"2021-09-20T18:40:18.677521Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sum(transformer.explained_variance_ratio_)","metadata":{"execution":{"iopub.status.busy":"2021-09-20T18:43:43.985651Z","iopub.execute_input":"2021-09-20T18:43:43.986047Z","iopub.status.idle":"2021-09-20T18:43:43.992733Z","shell.execute_reply.started":"2021-09-20T18:43:43.986011Z","shell.execute_reply":"2021-09-20T18:43:43.991407Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pcs = []\nfor i in range(0, len(all_img_arrays), batch_size):\n    print(i)\n    batch = np.array(all_img_arrays[i:i+batch_size])\n    batch_pcs = transformer.transform(batch)\n    pcs.append(batch_pcs)\npcs = np.concatenate(pcs, axis=0)","metadata":{"scrolled":true,"execution":{"iopub.status.busy":"2021-09-20T18:43:46.12959Z","iopub.execute_input":"2021-09-20T18:43:46.129931Z","iopub.status.idle":"2021-09-20T18:43:53.271128Z","shell.execute_reply.started":"2021-09-20T18:43:46.129897Z","shell.execute_reply":"2021-09-20T18:43:53.270086Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"del all_img_arrays","metadata":{"execution":{"iopub.status.busy":"2021-09-20T18:02:08.821532Z","iopub.execute_input":"2021-09-20T18:02:08.821895Z","iopub.status.idle":"2021-09-20T18:02:08.896367Z","shell.execute_reply.started":"2021-09-20T18:02:08.821863Z","shell.execute_reply":"2021-09-20T18:02:08.895476Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pixel_cols = [f'pixel-{i}' for i in range(128*128)]\nimg_df[pixel_cols] = all_img_array\ncollapsed_imgs = img_df.groupby('patient_ids').apply(lambda x: x[pixel_cols].mean())","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:23:49.690961Z","iopub.execute_input":"2021-09-21T04:23:49.691356Z","iopub.status.idle":"2021-09-21T04:25:23.041191Z","shell.execute_reply.started":"2021-09-21T04:23:49.691318Z","shell.execute_reply":"2021-09-21T04:25:23.04002Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"collapsed_imgs['patient_ids'] = img_df.groupby('patient_ids').apply(lambda x: x.iloc[0]['patient_ids'])\ncollapsed_imgs['labels'] = img_df.groupby('patient_ids').apply(lambda x: x.iloc[0]['labels'])","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:25:23.042714Z","iopub.execute_input":"2021-09-21T04:25:23.043194Z","iopub.status.idle":"2021-09-21T04:25:38.925395Z","shell.execute_reply.started":"2021-09-21T04:25:23.043131Z","shell.execute_reply":"2021-09-21T04:25:38.924262Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"collapsed_imgs","metadata":{"execution":{"iopub.status.busy":"2021-09-20T23:14:34.368391Z","iopub.execute_input":"2021-09-20T23:14:34.368781Z","iopub.status.idle":"2021-09-20T23:14:34.421069Z","shell.execute_reply.started":"2021-09-20T23:14:34.36874Z","shell.execute_reply":"2021-09-20T23:14:34.420117Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import random\nrandom.seed(369)\n\nunique_pids = collapsed_imgs['patient_ids'].unique()\ntrain_pids = random.sample(list(unique_pids), int(0.8*len(unique_pids)))\ntrain_df = collapsed_imgs[collapsed_imgs['patient_ids'].isin(train_pids)]\nval_df = collapsed_imgs[~collapsed_imgs['patient_ids'].isin(train_pids)]","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:25:38.926929Z","iopub.execute_input":"2021-09-21T04:25:38.927392Z","iopub.status.idle":"2021-09-21T04:25:38.981112Z","shell.execute_reply.started":"2021-09-21T04:25:38.927345Z","shell.execute_reply":"2021-09-21T04:25:38.980325Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from catboost import CatBoostClassifier","metadata":{"execution":{"iopub.status.busy":"2021-09-20T22:13:29.184737Z","iopub.execute_input":"2021-09-20T22:13:29.18521Z","iopub.status.idle":"2021-09-20T22:13:29.190154Z","shell.execute_reply.started":"2021-09-20T22:13:29.185173Z","shell.execute_reply":"2021-09-20T22:13:29.189181Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"clf = GradientBoostingClassifier(verbose=1, learning_rate=0.1, n_estimators=500)\nclf.fit(train_df[pixel_cols], train_df['labels'])","metadata":{"scrolled":true,"execution":{"iopub.status.busy":"2021-09-21T04:25:38.982102Z","iopub.execute_input":"2021-09-21T04:25:38.98245Z","iopub.status.idle":"2021-09-21T04:32:28.821783Z","shell.execute_reply.started":"2021-09-21T04:25:38.982424Z","shell.execute_reply":"2021-09-21T04:32:28.820709Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"clf.score(val_df[pixel_cols], val_df['labels'])","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:32:28.824579Z","iopub.execute_input":"2021-09-21T04:32:28.824938Z","iopub.status.idle":"2021-09-21T04:32:28.995834Z","shell.execute_reply.started":"2021-09-21T04:32:28.824909Z","shell.execute_reply":"2021-09-21T04:32:28.994779Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"y_pred = clf.predict_proba(val_df[pixel_cols])\ny_pred = y_pred[:, 1]","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:32:28.997407Z","iopub.execute_input":"2021-09-21T04:32:28.997844Z","iopub.status.idle":"2021-09-21T04:32:29.158068Z","shell.execute_reply.started":"2021-09-21T04:32:28.997796Z","shell.execute_reply":"2021-09-21T04:32:29.157206Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(y_pred)","metadata":{"execution":{"iopub.status.busy":"2021-09-21T04:32:29.15942Z","iopub.execute_input":"2021-09-21T04:32:29.159838Z","iopub.status.idle":"2021-09-21T04:32:29.323439Z","shell.execute_reply.started":"2021-09-21T04:32:29.159798Z","shell.execute_reply":"2021-09-21T04:32:29.322197Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"val_df['y_pred'] = y_pred\n\nmean_pred = val_df.mean()\ny_pred_agg = val_df.groupby('patient_ids').apply(lambda x: x['y_pred'].mean())\ny_true = val_df.groupby('patient_ids').apply(lambda x: x['labels'].iloc[0])","metadata":{"execution":{"iopub.status.busy":"2021-09-20T21:45:02.372034Z","iopub.execute_input":"2021-09-20T21:45:02.372478Z","iopub.status.idle":"2021-09-20T21:45:02.565755Z","shell.execute_reply.started":"2021-09-20T21:45:02.372438Z","shell.execute_reply":"2021-09-20T21:45:02.564108Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(y_pred_agg)","metadata":{"execution":{"iopub.status.busy":"2021-09-20T21:45:04.145245Z","iopub.execute_input":"2021-09-20T21:45:04.145697Z","iopub.status.idle":"2021-09-20T21:45:04.318146Z","shell.execute_reply.started":"2021-09-20T21:45:04.145659Z","shell.execute_reply":"2021-09-20T21:45:04.317415Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fpr, tpr, thresholds = roc_curve(y_true, y_pred_agg, pos_label=1)\nauc(fpr, tpr)","metadata":{"execution":{"iopub.status.busy":"2021-09-20T21:45:06.198964Z","iopub.execute_input":"2021-09-20T21:45:06.199629Z","iopub.status.idle":"2021-09-20T21:45:06.213574Z","shell.execute_reply.started":"2021-09-20T21:45:06.199556Z","shell.execute_reply":"2021-09-20T21:45:06.212178Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}