{"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":"# 🎗️[EDA] RSNA Screening Mammography Breast Cancer Detection\n---","metadata":{}},{"cell_type":"markdown","source":"Welcome to this notebook on **exploratory data analysis (EDA)** for the competition [*RSNA Breast Cancer Detection Competition*](https://www.kaggle.com/competitions/rsna-breast-cancer-detection). In this notebook, we will be using a dataset containing images of different patients in dicom format, along with the corresponding meta data for each patient in a csv table. Our goal is to understand the characteristics of the data and identify any potential issues or inconsistencies that may need to be addressed before moving on to the modeling phase.\n\nTo begin, we will load the data and take a look at the structure of the images and meta data. We will also perform some basic statistical analysis to get a sense of the overall distribution and range of values in the data. Additionally, we will visualize the data in various ways to help us gain a better understanding of the relationships between different variables and the underlying patterns present in the data.\n\nBy the end of this notebook, we should have a good understanding of the data and be prepared to move on to the next steps in the machine learning process. Let's get started!\n\n[*Link to the competition dataset*](https://www.kaggle.com/competitions/rsna-breast-cancer-detection/data)\n\n\n*This notebook has been inspired by the [work of Laura Fink](https://www.kaggle.com/code/allunia/rsna-breast-cancer-eda)*.","metadata":{}},{"cell_type":"markdown","source":"# Import packages","metadata":{}},{"cell_type":"code","source":"!pip install -qU python-gdcm pydicom pylibjpeg\n\n# Data manipulation and visualization libraries\nimport numpy as np\nimport pandas as pd\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nimport plotly.express as px\nimport cv2\n\n# File and directory handling libraries\nimport pydicom\nfrom os import listdir\n\n# Statistical and data processing libraries\nfrom scipy.stats import mode, skew\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.mixture import GaussianMixture\n\n# Progress tracking and warning suppression libraries\nfrom tqdm.notebook import trange\nimport warnings\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\n\n# Set seaborn style\nsns.set_style('darkgrid')","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2023-05-13T16:01:54.268397Z","iopub.execute_input":"2023-05-13T16:01:54.269195Z","iopub.status.idle":"2023-05-13T16:02:13.105195Z","shell.execute_reply.started":"2023-05-13T16:01:54.269153Z","shell.execute_reply":"2023-05-13T16:02:13.104198Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Exploring the patient metadata","metadata":{}},{"cell_type":"markdown","source":"The test set only contains an example for a single patient. Let's only focus on the train set.","metadata":{}},{"cell_type":"code","source":"filepath = '/kaggle/input/rsna-breast-cancer-detection/train.csv'\ndata = pd.read_csv(filepath)\ndata.head()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:13.107268Z","iopub.execute_input":"2023-05-13T16:02:13.1077Z","iopub.status.idle":"2023-05-13T16:02:13.244377Z","shell.execute_reply.started":"2023-05-13T16:02:13.10766Z","shell.execute_reply":"2023-05-13T16:02:13.243128Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Missing values ","metadata":{}},{"cell_type":"code","source":"data.isna().sum()\n#used to check the number of missing values (NaN or null values) in each column of my data","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:34.90933Z","iopub.execute_input":"2023-05-13T16:02:34.9098Z","iopub.status.idle":"2023-05-13T16:02:34.928525Z","shell.execute_reply.started":"2023-05-13T16:02:34.909759Z","shell.execute_reply":"2023-05-13T16:02:34.927406Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"There are missing values only for age, BIRADS and density.","metadata":{}},{"cell_type":"markdown","source":"## Basic information","metadata":{}},{"cell_type":"code","source":"num_patients = data['patient_id'].nunique()\nmin_patient_age = int(data['age'].min())\nmax_patient_age = int(data['age'].max())\ngroupby_id = data.groupby('patient_id')['cancer'].max()\nn_negative = (groupby_id == 0).sum()\nn_positive = (groupby_id == 1).sum()\n\nprint(f\"There are {num_patients} different patients in the train set.\\n\")\nprint(f\"The younger patient is {min_patient_age} years old.\")\nprint(f\"The older patient is {max_patient_age} years old.\\n\")\nprint(f\"There are {n_negative} patients negative to breast cancer. Ratio = {n_negative / num_patients}\")\nprint(f\"There are {n_positive} patients positive to breast cancer. Ratio = {n_positive / num_patients}\")","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:38.059708Z","iopub.execute_input":"2023-05-13T16:02:38.060446Z","iopub.status.idle":"2023-05-13T16:02:38.086691Z","shell.execute_reply.started":"2023-05-13T16:02:38.060407Z","shell.execute_reply":"2023-05-13T16:02:38.085489Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The classification problem is very unbalanced. About 4% of the patients are in the positive class which means that we will have a minority of cancer examples.","metadata":{}},{"cell_type":"markdown","source":"## Age distribution","metadata":{}},{"cell_type":"code","source":"ages = data.groupby('patient_id')['age'].apply(lambda x: x.unique()[0])\ncancer_ages = data[data['cancer'] == 1].groupby('patient_id')['age'].apply(lambda x: x.unique()[0])\nno_cancer_ages = data[data['cancer'] == 0].groupby('patient_id')['age'].apply(lambda x: x.unique()[0])\n\nplt.figure(figsize=(16, 10))\n\nplt.subplot(1, 2, 1)\nsns.histplot(ages, bins=63, color='orange', kde=True) # bins to indicate the number of devisions of age axis\nplt.title(\"All the patient\")\nplt.xlim(33, 89)\n\nplt.subplot(2, 2, 2)\nsns.histplot(cancer_ages, bins=51, color='red', kde=True)\nplt.title(\"Patients with cancer\")\nplt.xlim(33, 89)\n\nplt.subplot(2, 2, 4)\nsns.histplot(no_cancer_ages, bins=63, color='green', kde=True)\nplt.title(\"Patients without cancer\")\nplt.xlim(33, 89)\n\nplt.suptitle(\"Age distribution of the patients\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:42.29998Z","iopub.execute_input":"2023-05-13T16:02:42.300747Z","iopub.status.idle":"2023-05-13T16:02:44.675697Z","shell.execute_reply.started":"2023-05-13T16:02:42.300692Z","shell.execute_reply":"2023-05-13T16:02:44.674431Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Notice that the age distribution of the patients with a breast cancer is **negatively skewed** compared to the healthy ones.","metadata":{}},{"cell_type":"code","source":"# Statistics\nprint(\"Mean:\", ages.mean())\nprint(\"Std:\", ages.std())\nprint(\"Q1:\", ages.quantile(0.25))\nprint(\"Median:\", ages.median())\nprint(\"Q3:\", ages.quantile(0.75))\nprint(\"Mode:\", ages.mode()[0])","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:49.367703Z","iopub.execute_input":"2023-05-13T16:02:49.368141Z","iopub.status.idle":"2023-05-13T16:02:49.383722Z","shell.execute_reply.started":"2023-05-13T16:02:49.368104Z","shell.execute_reply":"2023-05-13T16:02:49.382432Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n- Most of the patients are **older than 40 years old**\n- There is a a **peak at the age of 50**\n- Then, there is a **plateau until 70 years old** before the count of patients drops\n- **Patients with cancer** are usually **older** than the others\n- **Age will be an important feature** for a future model","metadata":{}},{"cell_type":"markdown","source":"## Number of images per patient","metadata":{}},{"cell_type":"markdown","source":"","metadata":{}},{"cell_type":"code","source":"n_images_per_patient = data['patient_id'].value_counts()\nplt.figure(figsize=(16, 6))\nsns.countplot(n_images_per_patient, palette='Reds_r')\nplt.title(\"Number of images taken per patients\")\nplt.xlabel('Number of images taken')\nplt.ylabel('Count of patients')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:53.467977Z","iopub.execute_input":"2023-05-13T16:02:53.468735Z","iopub.status.idle":"2023-05-13T16:02:53.759202Z","shell.execute_reply.started":"2023-05-13T16:02:53.46868Z","shell.execute_reply":"2023-05-13T16:02:53.758038Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* Most of the patients have 4 images (2 views per side). However there are sometimes more of them.","metadata":{}},{"cell_type":"markdown","source":"## Image features","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(2, 3, figsize=(16, 8))\nsns.countplot(data['laterality'], palette='Blues_r', ax=ax[0, 0])\nsns.countplot(data['implant'], palette='Greens_r', ax=ax[0, 1])\nsns.countplot(data['difficult_negative_case'], palette='Reds_r', ax=ax[0, 2])\nsns.countplot(data['view'], palette='Oranges_r', ax=ax[1, 0])\nsns.countplot(data['density'], palette='Purples_r', order=['A', 'B', 'C', 'D'], ax=ax[1, 1])\nsns.countplot(data['site_id'], palette='Greys_r', ax=ax[1, 2])\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:02:57.008298Z","iopub.execute_input":"2023-05-13T16:02:57.008722Z","iopub.status.idle":"2023-05-13T16:02:57.825811Z","shell.execute_reply.started":"2023-05-13T16:02:57.008687Z","shell.execute_reply":"2023-05-13T16:02:57.824618Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* The pictures are balanced in term of laterality\n* Only a few images have implents\n* Some images were difficult to diagnose\n* There are usually only two types of views which are CC and MLO\n* There is a minority of images where the breast tissus density is very large (D) or very low (A). Most of the time, the density is in the middle (B) and (C)\n* Images where taken in two different sites in a balanced way","metadata":{}},{"cell_type":"markdown","source":"## Target features","metadata":{}},{"cell_type":"code","source":"biopsy_counts = data.groupby('cancer')['biopsy'].value_counts().unstack().fillna(0)\nbiopsy_perc = biopsy_counts.transpose() / biopsy_counts.sum(axis=1)\n\nfig, ax = plt.subplots(1, 2, figsize=(10, 4))\nsns.countplot(data['cancer'], palette='Greens', ax=ax[0])\nsns.heatmap(biopsy_perc, square=True, annot=True, fmt='.1%', cmap='Blues', ax=ax[1])\nax[0].set_title(\"Number of images showing cancer\")\nax[1].set_title(\"Percentage of images\\nresulting in a biopsy\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:01.608051Z","iopub.execute_input":"2023-05-13T16:03:01.608865Z","iopub.status.idle":"2023-05-13T16:03:01.98346Z","shell.execute_reply.started":"2023-05-13T16:03:01.608825Z","shell.execute_reply":"2023-05-13T16:03:01.982235Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* The number of sample images with cancer is very low compared to the number of healthy breast samples. It will be more difficult to make the model learn the desired pattern.\n* All patients with a cancer had a biopsy\n* Only a small percentage of images without cancer resulted in a biopsy","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(1, 3, figsize=(16, 4))\nsns.countplot(data[data['cancer'] == True]['invasive'], ax=ax[0], palette='Reds')\nsns.countplot(data[data['cancer'] == False]['BIRADS'], order=[0, 1, 2], ax=ax[1], palette='Blues')\nsns.countplot(data[data['cancer'] == True]['BIRADS'], order=[0, 1, 2], ax=ax[2], palette='Blues')\nax[0].set_title(\"Count of invasive cancer images\")\nax[1].set_title(\"BIRADS for healthy images\")\nax[2].set_title(\"BIRADS for cancer images\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:05.900715Z","iopub.execute_input":"2023-05-13T16:03:05.90116Z","iopub.status.idle":"2023-05-13T16:03:06.33294Z","shell.execute_reply.started":"2023-05-13T16:03:05.901122Z","shell.execute_reply":"2023-05-13T16:03:06.331809Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"BIRADS value signification:\n* 0: required follow-up\n* 1: rated as negative for cancer\n* 2: rated as normal","metadata":{}},{"cell_type":"markdown","source":"### Insights\n* Most of the images with cancer are invasive\n* A few healthy images led to a follow-up\n* All cancer images needed a follow-up","metadata":{}},{"cell_type":"markdown","source":"## Machine IDs","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(14, 5))\nsns.countplot(data['machine_id'])\nplt.title(\"Count of images taken by machine ID\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:11.421106Z","iopub.execute_input":"2023-05-13T16:03:11.421569Z","iopub.status.idle":"2023-05-13T16:03:11.679628Z","shell.execute_reply.started":"2023-05-13T16:03:11.421529Z","shell.execute_reply":"2023-05-13T16:03:11.678409Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* The images were taken from 10 different machines\n* Most of the images are from machine 49, 21, 29 and 48\n* Machine difference can lead to image distribution difference","metadata":{}},{"cell_type":"markdown","source":"# Exploring the image data","metadata":{}},{"cell_type":"markdown","source":"### DICOM files\nThe image data is in DICOM format. **DICOM (Digital Imaging and Communications in Medicine)** is a standard for storing and transmitting medical images and related information.\n\nIt consists of a set of data elements that are organized into a file structure. These data elements contain information about the medical image, such as the patient's name and medical record number, the image modality (e.g., CT, MRI, X-ray), the date and time the image was taken, and the image itself. The image data can be stored in various formats, such as 8-bit or 16-bit grayscale, or 24-bit color.\n\nIn addition to the image data, the DICOM format also includes metadata that describes the characteristics of the image, such as the image resolution, the size of the image in pixels, and the orientation of the image. This metadata is important for accurately displaying and interpreting the image.\n\nThe DICOM format is widely used in the medical community for storing, sharing, and analyzing medical images. It is supported by a wide range of medical devices, such as scanners, modalities, and workstations, and is used in hospitals, clinics, and research facilities around the world.","metadata":{}},{"cell_type":"code","source":"train_path = '/kaggle/input/rsna-breast-cancer-detection/train_images'\ntest_path = '/kaggle/input/rsna-breast-cancer-detection/test_images'","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:16.888593Z","iopub.execute_input":"2023-05-13T16:03:16.889488Z","iopub.status.idle":"2023-05-13T16:03:16.895414Z","shell.execute_reply.started":"2023-05-13T16:03:16.889434Z","shell.execute_reply":"2023-05-13T16:03:16.893449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Getting all scans of a single patient","metadata":{}},{"cell_type":"code","source":"def load_patient_scans(path, patient_id):\n    patient_path = path + '/' + str(patient_id)\n    return [pydicom.dcmread(patient_path + '/' + file) for file in listdir(patient_path)]","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:20.124807Z","iopub.execute_input":"2023-05-13T16:03:20.125548Z","iopub.status.idle":"2023-05-13T16:03:20.132114Z","shell.execute_reply.started":"2023-05-13T16:03:20.125493Z","shell.execute_reply":"2023-05-13T16:03:20.131187Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can now load all the scans for the first patient that we will use for exploring.","metadata":{}},{"cell_type":"code","source":"# Load all scans of the twelfth patient\npatient_ids = data['patient_id'].unique()\nscans = load_patient_scans(train_path, patient_ids[11])","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:23.387947Z","iopub.execute_input":"2023-05-13T16:03:23.389094Z","iopub.status.idle":"2023-05-13T16:03:25.047214Z","shell.execute_reply.started":"2023-05-13T16:03:23.38904Z","shell.execute_reply":"2023-05-13T16:03:25.046196Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Raw pixelarrays\n\nLet's check the pixel distribution of an image.","metadata":{}},{"cell_type":"code","source":"# Look at raw pixelarrays\nfig, ax = plt.subplots(1, 2, figsize=(20, 5))\nim = ax[0].imshow(scans[0].pixel_array, cmap='jet')\nax[0].grid(False)\nfig.colorbar(im, ax=ax[0])\nsns.histplot(scans[0].pixel_array.flatten(), ax=ax[1], bins=50)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:28.223623Z","iopub.execute_input":"2023-05-13T16:03:28.224427Z","iopub.status.idle":"2023-05-13T16:03:40.322096Z","shell.execute_reply.started":"2023-05-13T16:03:28.224385Z","shell.execute_reply":"2023-05-13T16:03:40.320747Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The value of the background is zero. This means that the raw pixelarray is not in HU. Some values are very large compared to what we would expect for soft tissues. The images should be converted to HU in order to have something consistant.\n\nLet's dive deep into the different machine IDs.","metadata":{}},{"cell_type":"code","source":"def get_scan_info():\n    modes, rows, cols = [], [], []\n    machine_ids = data['machine_id'].unique()\n    for m_id in machine_ids:\n        m_id_modes, m_id_rows, m_id_cols = [], [], []\n        print(f\"Machine id {m_id} in progress\")\n        patient_ids = data[data['machine_id'] == m_id]['patient_id'].unique()\n        for n in range(50):\n            try:\n                scan = load_patient_scans(train_path, patient_ids[n])[0]\n                m_id_modes.append(mode(scan.pixel_array.flatten())[0][0])\n                m_id_rows.append(scan.Rows)\n                m_id_cols.append(scan.Columns)\n            except IndexError:\n                break\n        modes.append(m_id_modes)\n        rows.append(m_id_rows)\n        cols.append(m_id_cols)\n    return modes, rows, cols","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:46.987785Z","iopub.execute_input":"2023-05-13T16:03:46.988715Z","iopub.status.idle":"2023-05-13T16:03:46.997752Z","shell.execute_reply.started":"2023-05-13T16:03:46.988669Z","shell.execute_reply":"2023-05-13T16:03:46.996451Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"modes, rows, cols = get_scan_info()\n\nmachine_ids = data['machine_id'].unique()\nmedians = [np.median(x) for x in modes]\nstds = [np.std(x) for x in modes]\nrows = [np.mean(x) for x in rows]\ncols = [np.mean(x) for x in cols]\ndf = pd.DataFrame(data={'Machine ID': machine_ids, 'Mode (median)': medians, 'Mode (std)': stds, 'Rows (mean)': rows, 'Cols (mean)': cols})\ndf.astype(int).set_index('Machine ID').T","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:03:55.908625Z","iopub.execute_input":"2023-05-13T16:03:55.909808Z","iopub.status.idle":"2023-05-13T16:15:40.429729Z","shell.execute_reply.started":"2023-05-13T16:03:55.909759Z","shell.execute_reply":"2023-05-13T16:15:40.427709Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* It seems that the pixel distribution is different based on the machine ID. In fact, the mode corresponds to the background pixels. Most of the time, this value is calibrated to 0 but for some machines (e.g. IDs 29 and 210), it is above a 3000 and 1000 respectively.\n* The value of Photometric Interpretation must be observed. There are 2 types:\n - In a *MONOCHROME1* image, the pixel values represent the grayscale values of the image, with higher values corresponding to brighter pixels and lower values corresponding to darker pixels.\n - In a *MONOCHROME2* image, the pixel values are reversed, with higher values corresponding to darker pixels and lower values corresponding to brighter pixels.\n* The image size also depends on the machine. Again, machines 29 and 210 have a high resolution compared to others like machine IDs 216 and 197.\n* It will be necessary to normalize the image size and the pixel values in order to train a robust model.","metadata":{}},{"cell_type":"markdown","source":"## Display examples for different machines","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(22, 8))\nfor i, m_id in enumerate(machine_ids):\n    patient_ids = data[data['machine_id'] == m_id]['patient_id'].unique()\n    scan = load_patient_scans(train_path, patient_ids[0])[0] # Load first scan of first patient\n    plt.subplot(2, 5, i+1)\n    plt.imshow(scan.pixel_array, cmap='jet')\n    plt.title(f\"Machine {m_id}\")\n    plt.colorbar()\n    plt.grid(False)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:16:17.635886Z","iopub.execute_input":"2023-05-13T16:16:17.636507Z","iopub.status.idle":"2023-05-13T16:16:33.296712Z","shell.execute_reply.started":"2023-05-13T16:16:17.636394Z","shell.execute_reply":"2023-05-13T16:16:33.295708Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* For machines 29 and 210, the background values seem to be switched. In other words, it appears to be the highest value (white color). For the other images, it appears to be the lowest value (black color). This should be taken into account during processing.\n* These images are just one example for each machine. It could be also different for other examples.\n\n### Tips\n* It is possible to fix the photometric interpretation by applying these formulas to the pixelarrays:\n - If MONOCRHOME1:\n >array = array.max() - array\n - If MONOCHROME2:\n >array = array - array.min()\n \nThis makes the way all the images are interpreted the same. The background of the images will be set to zero.","metadata":{}},{"cell_type":"markdown","source":"## Images showing implants","metadata":{}},{"cell_type":"markdown","source":"In this dataset, there are scans of breasts with implants. It can be interesting to have a look at them and compare them with other scans.","metadata":{}},{"cell_type":"code","source":"m_id_implants = data[data['implant'] == 1]['machine_id'].unique()\nprint(\"Scans showing implents are from machines\", m_id_implants)","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:16:45.010522Z","iopub.execute_input":"2023-05-13T16:16:45.010956Z","iopub.status.idle":"2023-05-13T16:16:45.020368Z","shell.execute_reply.started":"2023-05-13T16:16:45.010921Z","shell.execute_reply":"2023-05-13T16:16:45.019193Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"patient_ids = data[data['implant'] == 1]['patient_id'].unique()\n\n# Display scans showing implants\nplt.figure(figsize=(22, 8))\nfor i in range(10):\n    scan = load_patient_scans(train_path, patient_ids[i])[0] # Load first scan of the patient\n    plt.subplot(2, 5, i+1)\n    plt.imshow(scan.pixel_array, cmap='jet')\n    plt.title(f\"Patient {patient_ids[i]}\")\n    plt.grid(False)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:16:55.854597Z","iopub.execute_input":"2023-05-13T16:16:55.855774Z","iopub.status.idle":"2023-05-13T16:17:12.113781Z","shell.execute_reply.started":"2023-05-13T16:16:55.855727Z","shell.execute_reply":"2023-05-13T16:17:12.112573Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The implant is distinguished by the white mass in the breast. However, it is sometimes not very easy to see it like with patients 10960 and 13095.\n\nLet's check all the scans of a single patient having implant.","metadata":{}},{"cell_type":"code","source":"# Display scans of patient 13095\nplt.figure(figsize=(22, 8))\nscans = load_patient_scans(train_path, 13095)\nfor i in range(10):\n    plt.subplot(2, 5, i+1)\n    plt.imshow(scans[i].pixel_array, cmap='jet')\n    plt.grid(False)\nplt.suptitle(\"All scans of patient 13095\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:17:24.778851Z","iopub.execute_input":"2023-05-13T16:17:24.780151Z","iopub.status.idle":"2023-05-13T16:17:35.519331Z","shell.execute_reply.started":"2023-05-13T16:17:24.780079Z","shell.execute_reply":"2023-05-13T16:17:35.518423Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* All the scans showing implants are from machines 49 and 170.\n* The implants are sometimes only present on one side and not in both breasts for a given patient.\n* The metadata does not specify the presence of implent for a scan but for a given patient. It is therefore likely that an image without an implant indicates its presence.","metadata":{}},{"cell_type":"markdown","source":"## Images with cancer","metadata":{}},{"cell_type":"code","source":"# Choose to display images with or without cancer\ndef display_cancer_or_not(cancer=True):\n    cancer_scans = data[data['cancer'] == int(cancer)].sample(frac=1, random_state=0)\n    plt.figure(figsize=(22, 10))\n    for i in trange(10):\n        patient = str(cancer_scans.iloc[i][['patient_id']][0])\n        file = str(cancer_scans.iloc[i][['image_id']][0]) + '.dcm'\n        scan = pydicom.dcmread(train_path + '/' + patient + '/' + file)\n        plt.subplot(2, 5, i+1)\n        plt.imshow(scan.pixel_array, cmap='jet')\n        plt.title(f\"Patient {patient}\\nScan {file}\")\n        plt.grid(False)\n    plt.suptitle(f\"Cancer = {cancer}\")\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:18:59.55587Z","iopub.execute_input":"2023-05-13T16:18:59.557145Z","iopub.status.idle":"2023-05-13T16:18:59.566862Z","shell.execute_reply.started":"2023-05-13T16:18:59.557084Z","shell.execute_reply":"2023-05-13T16:18:59.565678Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Images with cancer\ndisplay_cancer_or_not(cancer=True)","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:19:05.250332Z","iopub.execute_input":"2023-05-13T16:19:05.25124Z","iopub.status.idle":"2023-05-13T16:19:25.620475Z","shell.execute_reply.started":"2023-05-13T16:19:05.251194Z","shell.execute_reply":"2023-05-13T16:19:25.61925Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Images without cancer\ndisplay_cancer_or_not(cancer=False)","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:20:02.680134Z","iopub.execute_input":"2023-05-13T16:20:02.681054Z","iopub.status.idle":"2023-05-13T16:20:20.318464Z","shell.execute_reply.started":"2023-05-13T16:20:02.681007Z","shell.execute_reply":"2023-05-13T16:20:20.317225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n* It is impossible for me to distinguish by eye an image with cancer from a healthy image\n* Again, there are images from different machine IDs which results in different pixel distributions","metadata":{}},{"cell_type":"markdown","source":"## Outliers analysis\nTo improve the performance of our future model, we also could focus on identifying and removing outliers. One approach to accomplish this is to create a table containing various statistical data on the distribution of pixels in the images. By comparing these statistics, we can quickly identify images that deviate significantly from the norm and eliminate them as outliers.\n\nThis task can be very computational expensive and this is why I will do it on just a small part of the patient images. Let's do it for the first scan of all the patients of machine ID 21.\n\nThe outlier analysis would require its own notebook to go through the images of the full dataset. This is something that I plan to do later.\n\n### Machine ID 21\nI am using the [preprocessed 256x256 images](https://www.kaggle.com/code/theoviel/dicom-resized-png-jpg) of Theo Viel in order to save computation time. The images do not have to be loaded with `pydicom`, and they are already smaller.\n\nThank you again to *Laura Fink* and her EDA [notebook](https://www.kaggle.com/code/allunia/rsna-breast-cancer-eda/notebook) that inspired me. She did the same work on patients of machine ID 49.","metadata":{}},{"cell_type":"code","source":"source = '../input/rsna-breast-cancer-256-pngs/'\ndata['path'] = source + data['patient_id'].astype(str) + \"_\" + data['image_id'].astype(str) + \".png\"","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:20:41.797793Z","iopub.execute_input":"2023-05-13T16:20:41.798242Z","iopub.status.idle":"2023-05-13T16:20:41.901342Z","shell.execute_reply.started":"2023-05-13T16:20:41.798205Z","shell.execute_reply":"2023-05-13T16:20:41.900255Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# All patients of machine ID 21\npatients = data[data.machine_id == 21].patient_id.unique()\n\nmeans, medians, stds, skews, paths, labels = [], [], [], [], [], []\n\nfor i in trange(len(patients)):\n    path = data[data['patient_id'] == patients[i]]['path'].iloc[0]\n    label = data[data['patient_id'] == patients[i]]['cancer'].iloc[0]\n    img = cv2.imread(path, cv2.IMREAD_GRAYSCALE)\n    means.append(np.mean(img))\n    medians.append(np.median(img))\n    stds.append(np.std(img))\n    skews.append(skew(img, axis=None))\n    paths.append(path)\n    labels.append(label)\n    \nstats = pd.DataFrame()\nstats['mean'] = means\nstats['median'] = medians\nstats['std'] = stds\nstats['skew'] = skews\nstats['path'] = paths\nstats['cancer'] = labels\nstats.head()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:20:45.74788Z","iopub.execute_input":"2023-05-13T16:20:45.748297Z","iopub.status.idle":"2023-05-13T16:21:07.926699Z","shell.execute_reply.started":"2023-05-13T16:20:45.748262Z","shell.execute_reply":"2023-05-13T16:21:07.925251Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig = px.scatter_3d(stats, x='mean', y='std', z='skew', color='median')\nfig.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:21:17.317627Z","iopub.execute_input":"2023-05-13T16:21:17.318897Z","iopub.status.idle":"2023-05-13T16:21:18.493298Z","shell.execute_reply.started":"2023-05-13T16:21:17.318851Z","shell.execute_reply":"2023-05-13T16:21:18.492105Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n- Looking at the medians, there are different groups. The biggest one has a median around zero, which means that most of the pixels are the background. Orange dots represent images where the breast covers more than the half of the pixels.\n- The mean of those images is also greater which makes sense.\n- The skewness seems correlated with the mean and so with the breast size.","metadata":{}},{"cell_type":"code","source":"sns.heatmap(np.abs(stats.corr()), annot=True,\n            square=True, vmin=0, vmax=1, cmap='Blues')\nplt.title(\"Correlation matrix of the statistics\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:21:55.348986Z","iopub.execute_input":"2023-05-13T16:21:55.3503Z","iopub.status.idle":"2023-05-13T16:21:55.668435Z","shell.execute_reply.started":"2023-05-13T16:21:55.350248Z","shell.execute_reply":"2023-05-13T16:21:55.66717Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As previously mentioned, means and skewness are highly correlated, while standard deviations and medians are least correlated. There is also no way to predict if a patient has cancer based on such statistics.\n\nA Gaussian mixture model can be used to group similar images based on the statistics that have been calculated. By utilizing the tools from the `scikit-learn` library, we can implement this model. It is possible to plot again the different clusters for visualization purposes, providing a clearer understanding of how the images have been grouped.","metadata":{}},{"cell_type":"code","source":"scaler = StandardScaler()\nX = stats.drop(['path', 'cancer'], axis=1)\nX = scaler.fit_transform(X)\n\ngmm = GaussianMixture(n_components=4, random_state=0)\nstats['cluster_label'] = gmm.fit_predict(X)","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:00.228036Z","iopub.execute_input":"2023-05-13T16:22:00.22886Z","iopub.status.idle":"2023-05-13T16:22:00.276622Z","shell.execute_reply.started":"2023-05-13T16:22:00.228818Z","shell.execute_reply":"2023-05-13T16:22:00.275547Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig = px.scatter_3d(stats, x='mean', y='std', z='skew', color='cluster_label')\nfig.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:06.709118Z","iopub.execute_input":"2023-05-13T16:22:06.709591Z","iopub.status.idle":"2023-05-13T16:22:06.779303Z","shell.execute_reply.started":"2023-05-13T16:22:06.70955Z","shell.execute_reply":"2023-05-13T16:22:06.778116Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's look at a couple of images belonging to the different clusters.","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(10,4,figsize=(20,40))\nfor l in range(4):\n    paths = stats[stats['cluster_label'] == l]['path'].values[0:10]\n    cancers = stats[stats['cluster_label'] == l]['cancer'].values[0:10]\n    for n in range(10):\n        img = cv2.imread(paths[n], cv2.IMREAD_GRAYSCALE)\n        ax[n,l].imshow(img, cmap='jet')\n        ax[n,l].set_title(f\"label {l} / cancer {cancers[n]}\")\n        ax[n,l].axis('off')","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:14.389051Z","iopub.execute_input":"2023-05-13T16:22:14.389576Z","iopub.status.idle":"2023-05-13T16:22:19.592703Z","shell.execute_reply.started":"2023-05-13T16:22:14.389531Z","shell.execute_reply":"2023-05-13T16:22:19.591443Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n- The images for a given label (along a column) are quite similar in their pixel distribution.\n- The contrasts vary between clusters.\n\nNow let's check some outliers since it was the purpose of this section. For this, let's calculate the log-likelihood of each sample. Then we check the value of the first 0.5% quantile.","metadata":{}},{"cell_type":"code","source":"stats['logL'] = gmm.score_samples(X)\nstats['logL'].quantile(0.005)","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:29.858521Z","iopub.execute_input":"2023-05-13T16:22:29.858942Z","iopub.status.idle":"2023-05-13T16:22:29.871252Z","shell.execute_reply.started":"2023-05-13T16:22:29.858906Z","shell.execute_reply":"2023-05-13T16:22:29.870339Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's say that the outliers we are interested in have a logL value that is less that -10 and we can display them.","metadata":{}},{"cell_type":"code","source":"outliers = stats[stats['logL'] < -10].sort_values(by='logL')\nprint(f\"{len(outliers)} outliers have been found!\")","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:33.955925Z","iopub.execute_input":"2023-05-13T16:22:33.956433Z","iopub.status.idle":"2023-05-13T16:22:33.965809Z","shell.execute_reply.started":"2023-05-13T16:22:33.956393Z","shell.execute_reply":"2023-05-13T16:22:33.964474Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"outliers['path'].iloc[0]","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:22:39.737738Z","iopub.execute_input":"2023-05-13T16:22:39.738726Z","iopub.status.idle":"2023-05-13T16:22:39.748023Z","shell.execute_reply.started":"2023-05-13T16:22:39.738676Z","shell.execute_reply":"2023-05-13T16:22:39.74669Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(20, 8))\nfor i in range(10):\n    plt.subplot(2, 5, i+1)\n    img = cv2.imread(outliers['path'].iloc[i], cv2.IMREAD_GRAYSCALE)\n    plt.imshow(img, cmap='jet')\n    plt.axis('off')","metadata":{"execution":{"iopub.status.busy":"2023-05-13T16:23:00.436838Z","iopub.execute_input":"2023-05-13T16:23:00.437763Z","iopub.status.idle":"2023-05-13T16:23:01.769541Z","shell.execute_reply.started":"2023-05-13T16:23:00.437715Z","shell.execute_reply":"2023-05-13T16:23:01.768422Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Insights\n- Outliers are images where there is either a lot of background or just a little surface.\n- It could also be images with a special contrast in the pixels.\n\nI am sure that there many more outliers from other machine IDs and by taking the full set of scans for a given patient. This long work will be achieved in a future notebook.","metadata":{}},{"cell_type":"markdown","source":"# Conclusion\n* The dataset is heavily unbalanced between scans with and without cancer.\n* Most of the patients are over 40 years old.\n* Images are quite large and will need to be rescaled during preprocessing.\n* Pixel distributions vary significantly depending on the machine ID used.\n* The dataset is also unbalanced in terms of images showing implants.\n* It is very difficult for a novice to distinguish a scan with cancer from a healthy one.\n\nPerforming EDA was very informative. It helped me understand what needs to be taken into consideration during the preprocessing steps.\n\n**Thank you for reading. I welcome any feedback you may have. 👋**","metadata":{}}]}