{"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":"**Key ideas** in this notebook:\n* create RGB (3D) images from spatially adjacent slices using a stride\n* for slice s, the corresponding RGB image is cretead by stacking [s-stride, s, s+stride] in the color channel\n* other winning solutions from previous RSNA Kaggle competitions used this technique to preprocess images and capture information from neighboring slices (e.g. [RSNA Intracranial Hemorrhage Detection 5th place solution](https://www.kaggle.com/c/rsna-intracranial-hemorrhage-detection/discussion/117232))\n* cervical vertebrae and fractions **might be more visible** in RGB images that display information from neighboring slices\n* the notebook displays the same slices for a given patient: 1) original slices from the training images, 2) original slices with the segmentation masks, 3) RGB images corresponding to these slices","metadata":{}},{"cell_type":"markdown","source":"## Libraries","metadata":{}},{"cell_type":"code","source":"# online\n! pip install python-gdcm\n! pip install pylibjpeg pylibjpeg-libjpeg pydicom","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:48:56.561125Z","iopub.execute_input":"2022-08-29T16:48:56.561687Z","iopub.status.idle":"2022-08-29T16:49:23.970384Z","shell.execute_reply.started":"2022-08-29T16:48:56.561574Z","shell.execute_reply":"2022-08-29T16:49:23.969132Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os\nimport glob\nimport pydicom\nimport nibabel as nib\nimport pandas as pd\nimport numpy as np\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\nimport matplotlib.pyplot as plt\nimport matplotlib.patches as patches\nimport seaborn as sns\nfrom joblib import Parallel, delayed\nfrom tqdm.notebook import tqdm\ntqdm.pandas()","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:23.973222Z","iopub.execute_input":"2022-08-29T16:49:23.973757Z","iopub.status.idle":"2022-08-29T16:49:24.964318Z","shell.execute_reply.started":"2022-08-29T16:49:23.973696Z","shell.execute_reply":"2022-08-29T16:49:24.963164Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# constants\nSTRIDE = 2\nDATA_DIR = \"../input/rsna-2022-cervical-spine-fracture-detection\"\nTRAIN_DIR = os.path.join(DATA_DIR, \"train_images\")\nSEGM_DIR = os.path.join(DATA_DIR, \"segmentations\")\n\ntrain_df = pd.read_csv(os.path.join(DATA_DIR, \"train.csv\"))\nprint(train_df.shape)\ntrain_df.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:24.965748Z","iopub.execute_input":"2022-08-29T16:49:24.966194Z","iopub.status.idle":"2022-08-29T16:49:25.005281Z","shell.execute_reply.started":"2022-08-29T16:49:24.966152Z","shell.execute_reply":"2022-08-29T16:49:25.004033Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Utilities","metadata":{}},{"cell_type":"code","source":"def load_dicom(dcm_path):\n    ds = pydicom.dcmread(dcm_path)\n    img = apply_voi_lut(ds.pixel_array, ds)\n    img -= img.min()\n    img /= img.max()\n    if not (ds.Rows == 512 and ds.Columns == 512):\n        img = cv2.resize(img, dsize=(512, 512))\n    return img\n\n\ndef show_segm_for_patient(patient_id, nrows, ncols, is_rgb, is_segm_on=True, step=1):\n    \"\"\"Shows slices with segmentation masks for the patient, skips slices that have empty segmentation masks.\"\"\"\n    if is_rgb:\n        slices_dir = os.path.join(\"./rgb_train_images\", patient_id)\n    else:\n        slices_dir = os.path.join(TRAIN_DIR, patient_id)\n    \n    num_slices = len(glob.glob(f\"{slices_dir}/*\"))\n    print(f\"Number of slices for patient: {num_slices}\")\n    \n    print(f\"Fracture information for patient: \")\n    print(f\"{train_df[train_df['StudyInstanceUID'] == patient_id].iloc[:, 1:-1]}\")\n    \n    segm_path = os.path.join(SEGM_DIR, f\"{patient_id}.nii\")\n    # source: https://www.kaggle.com/code/kretes/segmentations-fracture-zoom-in\n    segm_mask = nib.load(segm_path).get_fdata()\n        \n    # flip in z axis\n    segm_mask = np.flip(segm_mask, axis=-1)\n    # rotate 90 degrees in xy\n    segm_mask = np.rot90(segm_mask, axes=(0, 1))\n        \n    fig, axs = plt.subplots(nrows, ncols, figsize=(30, 20))\n    axs = axs.flatten()\n    \n    count = 0\n    idx = 0\n    while count < nrows*ncols:\n        if segm_mask[:, :, idx].max() > 0:\n            # filenames start from 1\n            if is_rgb:\n                img = np.load(os.path.join(slices_dir, f\"{idx+1}.npy\"))\n            else:\n                img = load_dicom(os.path.join(slices_dir, f\"{idx+1}.dcm\"))\n            axs[count].imshow(img, cmap=\"gray\")\n            axs[count].axis(\"off\")\n            \n            if is_segm_on:\n                segm_im = axs[count].imshow(segm_mask[:, :, idx], alpha=0.4)\n                segm_max = segm_mask[:, :, idx].max()\n\n                segm_labels = np.unique(segm_mask[:, :, idx][segm_mask[:, :, idx].nonzero()])\n                segm_colors = [segm_im.cmap(label/segm_max) for label in segm_labels]\n                ptch = [patches.Patch(color=segm_colors[i], label=f\"C{int(segm_labels[i])}\") for i in range(len(segm_labels))]\n                axs[count].legend(handles=ptch)\n            count += 1\n        idx += step","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:25.008539Z","iopub.execute_input":"2022-08-29T16:49:25.009307Z","iopub.status.idle":"2022-08-29T16:49:25.026361Z","shell.execute_reply.started":"2022-08-29T16:49:25.009259Z","shell.execute_reply":"2022-08-29T16:49:25.025392Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Original slices","metadata":{}},{"cell_type":"code","source":"show_segm_for_patient(\"1.2.826.0.1.3680043.10921\", nrows=4, ncols=5, is_rgb=False, is_segm_on=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:25.02775Z","iopub.execute_input":"2022-08-29T16:49:25.028134Z","iopub.status.idle":"2022-08-29T16:49:28.809137Z","shell.execute_reply.started":"2022-08-29T16:49:25.028097Z","shell.execute_reply":"2022-08-29T16:49:28.807787Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Original slices with segmentation masks","metadata":{}},{"cell_type":"code","source":"show_segm_for_patient(\"1.2.826.0.1.3680043.10921\", nrows=4, ncols=5, is_rgb=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:28.810988Z","iopub.execute_input":"2022-08-29T16:49:28.811403Z","iopub.status.idle":"2022-08-29T16:49:31.81974Z","shell.execute_reply.started":"2022-08-29T16:49:28.811365Z","shell.execute_reply":"2022-08-29T16:49:31.818615Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## RGB images corresponding to the slices","metadata":{}},{"cell_type":"code","source":"def save_rgb_imgs_for_patient(patient_id):\n    def create_rgb_with_adjacent_slices(slice_idx):\n        img_rgb = np.zeros((512, 512, 3))\n        prev_slice_idx = max(1, slice_idx - STRIDE)\n        next_slice_idx = min(num_slices, slice_idx + STRIDE)\n        img_rgb[:, :, 0] = load_dicom(os.path.join(slices_dir, f\"{prev_slice_idx}.dcm\"))\n        img_rgb[:, :, 1] = load_dicom(os.path.join(slices_dir, f\"{slice_idx}.dcm\"))\n        img_rgb[:, :, 2] = load_dicom(os.path.join(slices_dir, f\"{next_slice_idx}.dcm\"))\n        return img_rgb\n       \n    save_dir = os.path.join(\"./rgb_train_images\", patient_id)\n    os.makedirs(save_dir, exist_ok=True)\n    slices_dir = os.path.join(TRAIN_DIR, patient_id)\n    num_slices = len(glob.glob(f\"{TRAIN_DIR}/{patient_id}/*\"))\n    for slice_idx in np.arange(start=1, stop=num_slices+1):\n        img_rgb = create_rgb_with_adjacent_slices(slice_idx)\n        np.save(f\"{save_dir}/{slice_idx}.npy\", img_rgb)","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:31.82127Z","iopub.execute_input":"2022-08-29T16:49:31.821654Z","iopub.status.idle":"2022-08-29T16:49:31.83134Z","shell.execute_reply.started":"2022-08-29T16:49:31.821606Z","shell.execute_reply":"2022-08-29T16:49:31.830369Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"patient_ids = [\"1.2.826.0.1.3680043.10921\"]\n_ = Parallel(n_jobs=-1, backend=\"threading\")(delayed(save_rgb_imgs_for_patient)(patient_id) \\\n                                         for patient_id in tqdm(patient_ids, total=len(patient_ids)))","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:31.832522Z","iopub.execute_input":"2022-08-29T16:49:31.83291Z","iopub.status.idle":"2022-08-29T16:49:44.521102Z","shell.execute_reply.started":"2022-08-29T16:49:31.832861Z","shell.execute_reply":"2022-08-29T16:49:44.520166Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"show_segm_for_patient(\"1.2.826.0.1.3680043.10921\", nrows=4, ncols=5, is_rgb=True, is_segm_on=False)","metadata":{"execution":{"iopub.status.busy":"2022-08-29T16:49:44.522059Z","iopub.execute_input":"2022-08-29T16:49:44.522349Z","iopub.status.idle":"2022-08-29T16:49:47.434734Z","shell.execute_reply.started":"2022-08-29T16:49:44.522322Z","shell.execute_reply":"2022-08-29T16:49:47.433673Z"},"trusted":true},"execution_count":null,"outputs":[]}]}