{"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":"!pip install -qU python-gdcm pydicom pylibjpeg","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:07:58.813519Z","iopub.execute_input":"2023-01-12T17:07:58.8142Z","iopub.status.idle":"2023-01-12T17:08:10.005417Z","shell.execute_reply.started":"2023-01-12T17:07:58.814122Z","shell.execute_reply":"2023-01-12T17:08:10.003874Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\nimport pydicom\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\nimport gdcm\nimport pylibjpeg\nimport matplotlib.pyplot as plt\nimport matplotlib.cm as cm\nimport seaborn as sns\n\nimport os\nfrom os import listdir\n\nimport glob\nfrom tqdm.notebook import tqdm\nfrom joblib import Parallel, delayed\nimport cv2","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:23:03.366369Z","iopub.execute_input":"2023-01-12T17:23:03.366931Z","iopub.status.idle":"2023-01-12T17:23:03.374805Z","shell.execute_reply.started":"2023-01-12T17:23:03.366894Z","shell.execute_reply":"2023-01-12T17:23:03.373304Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Custom colors for charts\n\nc_0 = np.array([2, 48, 71,256])/256\nc_1 = np.array([251, 133, 0,256])/256\nc_2 = np.array([255, 183, 3,256])/256\nc_3 = np.array([33, 158, 188,256])/256\nc_4 = np.array([142, 202, 230,256])/256\n\n# Creating a custom colormap\n# The code comes from https://matplotlib.org/3.1.1/tutorials/colors/colormap-manipulation.html\nfrom matplotlib.colors import ListedColormap\n\nN = 256\nvals = np.ones((N, 4))\nvals[:, 0] = np.linspace(c_0[0], c_1[0], N)\nvals[:, 1] = np.linspace(c_0[1], c_1[1], N)\nvals[:, 2] = np.linspace(c_0[2], c_1[2], N)\ncustom_cmp1 = ListedColormap(vals)\n\n# A colormap with white in the middle\nN = 256\nvals = np.ones((N, 4))\nvals[:, 0] = np.concatenate((np.linspace(c_0[0], 1, 128), np.linspace(1, c_1[0], 128)), axis = None)\nvals[:, 1] = np.concatenate((np.linspace(c_0[1], 1, 128), np.linspace(1, c_1[1], 128)), axis = None)\nvals[:, 2] = np.concatenate((np.linspace(c_0[2], 1, 128), np.linspace(1, c_1[2], 128)), axis = None)\ncustom_cmp2 = ListedColormap(vals)\n\n# Allow value-dependant colorbar\n\ndef color_bar(y, cmap=custom_cmp1):\n    \n    # Create the matrix for colors\n    color_list = []\n    # Normalise the values\n    y = (y-min(y))/(max(y)-min(y))\n    \n    for i in range(len(y)):\n        color_list.append(cmap(int(N*y[i])))\n    \n    return color_list","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.776263Z","iopub.execute_input":"2023-01-12T17:08:10.776605Z","iopub.status.idle":"2023-01-12T17:08:10.79298Z","shell.execute_reply.started":"2023-01-12T17:08:10.776572Z","shell.execute_reply":"2023-01-12T17:08:10.791947Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# RSNA Breast Cancer Detection Project (Part 1)\n\nThe description of the competition project is the following :\n\n*According to the WHO, breast cancer is the most commonly occurring cancer worldwide. In 2020 alone, there were 2.3 million new breast cancer diagnoses and 685,000 deaths. Yet breast cancer mortality in high-income countries has dropped by 40% since the 1980s when health authorities implemented regular mammography screening in age groups considered at risk. Early detection and treatment are critical to reducing cancer fatalities, and your machine learning skills could help streamline the process radiologists use to evaluate screening mammograms.*\n\nThus, for obvious reasons the detection of a probable cancer in mammography scans could be a precious tool for fighting breast cancer. The data provided by the RSNA are distributed between a csv file containing **metadata** as the age of the patient or the presence or not of an implant and scans of mammography.","metadata":{}},{"cell_type":"markdown","source":"# 1. Data import and EDA\n\nAs explainer earlier, the dataset considered in this project is constituted of two parts : a folder with images and a csv file containing meta data. First, the metadata are expored. The procedure for importing the data have been inspired from the EDA performed by Laura Fink (<a href =https://www.kaggle.com/code/allunia/rsna-breast-cancer-eda>see</a>).\n\n## 1.1 Metadata\n\nThe metadata are imported from the project.","metadata":{}},{"cell_type":"code","source":"df_train = pd.read_csv('/kaggle/input/rsna-breast-cancer-detection/train.csv')\ndf_test = pd.read_csv('/kaggle/input/rsna-breast-cancer-detection/test.csv')\ndf_train.head()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.797111Z","iopub.execute_input":"2023-01-12T17:08:10.797622Z","iopub.status.idle":"2023-01-12T17:08:10.905686Z","shell.execute_reply.started":"2023-01-12T17:08:10.797586Z","shell.execute_reply":"2023-01-12T17:08:10.904493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train.shape","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.90768Z","iopub.execute_input":"2023-01-12T17:08:10.908294Z","iopub.status.idle":"2023-01-12T17:08:10.91549Z","shell.execute_reply.started":"2023-01-12T17:08:10.908249Z","shell.execute_reply":"2023-01-12T17:08:10.914353Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test.head()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.917465Z","iopub.execute_input":"2023-01-12T17:08:10.917919Z","iopub.status.idle":"2023-01-12T17:08:10.942207Z","shell.execute_reply.started":"2023-01-12T17:08:10.917877Z","shell.execute_reply":"2023-01-12T17:08:10.940927Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test.shape","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.943618Z","iopub.execute_input":"2023-01-12T17:08:10.943935Z","iopub.status.idle":"2023-01-12T17:08:10.955376Z","shell.execute_reply.started":"2023-01-12T17:08:10.943907Z","shell.execute_reply":"2023-01-12T17:08:10.95437Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The train dataset contains **54706 lines** for **14 columns**. However multiple rows corresponds to a single patient id.\n\nThe test set contains 4 rows for only 9 columns. The patient ID corresponds to a single patient. The columns *cancer, biopsy, invasive, BIRADS* and *difficult_negative_case* are not present in the test set while *prediction_id* have been added.","metadata":{}},{"cell_type":"markdown","source":"### Data type","metadata":{}},{"cell_type":"markdown","source":"The features types and statistics are looked closer thanks to the function info and describe.","metadata":{}},{"cell_type":"code","source":"df_train.info()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.956716Z","iopub.execute_input":"2023-01-12T17:08:10.957076Z","iopub.status.idle":"2023-01-12T17:08:10.981615Z","shell.execute_reply.started":"2023-01-12T17:08:10.957039Z","shell.execute_reply":"2023-01-12T17:08:10.980537Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train.describe()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:10.982898Z","iopub.execute_input":"2023-01-12T17:08:10.983202Z","iopub.status.idle":"2023-01-12T17:08:11.039323Z","shell.execute_reply.started":"2023-01-12T17:08:10.983175Z","shell.execute_reply":"2023-01-12T17:08:11.037956Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_train.describe(include='O')","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:11.040996Z","iopub.execute_input":"2023-01-12T17:08:11.04136Z","iopub.status.idle":"2023-01-12T17:08:11.077084Z","shell.execute_reply.started":"2023-01-12T17:08:11.041325Z","shell.execute_reply":"2023-01-12T17:08:11.075854Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The *site_id* indicates the number of hospitals from which the data have been taken. \nThe ID indicates that the data are coming from two hospitals with a majority of data coming from the hospital 1.\n\nThe patients *age* is between 26 and 89 years.\n\nThe target *cancer* in a int value of 0 or 1 with a mean at 0.021 meaning a low proportion of cancerous analysis.\n\nThe *laterality*, *view* and *density* are caterogrical values indicated as strings with a limited number of unique values.","metadata":{}},{"cell_type":"markdown","source":"As observed earlier, multiple picture can be attributed to a single patient. Thus, the number of unique patients is sutided.","metadata":{}},{"cell_type":"code","source":"df_train['patient_id'].nunique()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:11.078801Z","iopub.execute_input":"2023-01-12T17:08:11.080008Z","iopub.status.idle":"2023-01-12T17:08:11.087226Z","shell.execute_reply.started":"2023-01-12T17:08:11.079967Z","shell.execute_reply":"2023-01-12T17:08:11.086123Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**11913** unique patients are present in the train dataset.","metadata":{}},{"cell_type":"markdown","source":"### Missing values\n\nThe percentage of missing values is studied.","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(figsize = (4,4))\n\nax.barh(y = df_train.columns,\n        width = df_train.isna().sum(axis = 0),\n        log = True,\n        color = color_bar(df_train.isna().sum(axis = 0)))\n\nax.set_title('Missing values in the train dataset')\nax.set_xlabel('Number of missing values (log)')\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:11.088868Z","iopub.execute_input":"2023-01-12T17:08:11.089211Z","iopub.status.idle":"2023-01-12T17:08:11.735786Z","shell.execute_reply.started":"2023-01-12T17:08:11.089181Z","shell.execute_reply":"2023-01-12T17:08:11.734783Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Only three fields contains missing values : *age*, *BIRADS* and *density*.\n\nThe *BIRADS* and *density* are not present in the test set and thus can not be used for prediction. However, the *age* could be an imortant paramter and withdrawing or imputing the missing values could be necessary.","metadata":{}},{"cell_type":"markdown","source":"### Distribution, noisiness and outliers.\n\nThe distribution of the data as well as the presence of outliers is studied.","metadata":{}},{"cell_type":"code","source":"# Numerical and categorical data are studied separately\n\ndf_train_num = df_train.select_dtypes(exclude=['O', bool])\ndf_train_cat = df_train.select_dtypes(include='O')","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:11.740148Z","iopub.execute_input":"2023-01-12T17:08:11.741081Z","iopub.status.idle":"2023-01-12T17:08:11.751678Z","shell.execute_reply.started":"2023-01-12T17:08:11.74104Z","shell.execute_reply":"2023-01-12T17:08:11.7501Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, axes = plt.subplots(df_train_num.shape[1],\n                         2,\n                         figsize = (10, 16))\n\nfor i, col in enumerate(list(df_train_num.columns)):\n    sns.histplot(data = df_train[col],\n                 ax = axes[i, 0],\n                 color = c_1,\n                 alpha = 1,)\n    axes[i, 0].set_title(f'Histogram for {col}')\n    \n    sns.boxplot(data = df_train[col],\n                ax = axes[i, 1],\n                color = c_0,\n                orient = 'h')\n    axes[i, 1].set_title(f'Boxplot for {col}')\n    \nplt.tight_layout()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:11.753317Z","iopub.execute_input":"2023-01-12T17:08:11.753705Z","iopub.status.idle":"2023-01-12T17:08:15.978252Z","shell.execute_reply.started":"2023-01-12T17:08:11.753674Z","shell.execute_reply":"2023-01-12T17:08:15.976993Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As descibed earlier, only two sie are indicated in the *site_id*.\n\nThe *patient_id* and *image_id* do not show extreme deviations.\n\nThe *age* repartition of the patients could form a normal distribution centered on 60. It is worth noticing that a high frequency is observed at 50 (the age at which cancer screening is recommanded).\nMoreover, the frequency flattens rapidly bellow 40. Fiew outliers are observed bellow 30 and around 89.\n\nAs observed earlier, the proportion of cancerous analysis is way bellow the proprotion of begnin ones. The same trend can be observed for the other analysis.\n\nFinaly, the *machine_id* indicated that a large proprotion of analysis have been performed on the same machine.\n\n\nCategorical features are also studied :","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(df_train_cat.shape[1], 1, figsize = (5, 10))\n\nfor i, col in enumerate(list(df_train_cat.columns)):\n    axes[i].bar(x = df_train_cat[col].value_counts().index,\n                height = df_train_cat[col].value_counts(),\n                color = color_bar(df_train_cat[col].value_counts()))\n    # The view is better represented with a log scale\n    if col == 'view':\n        axes[i].set_yscale('log')\n        \n    axes[i].set_title(col)\nplt.show()\n    ","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:15.979384Z","iopub.execute_input":"2023-01-12T17:08:15.980286Z","iopub.status.idle":"2023-01-12T17:08:16.834316Z","shell.execute_reply.started":"2023-01-12T17:08:15.980249Z","shell.execute_reply":"2023-01-12T17:08:16.833155Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"First, the laterality shows that a balanced proportion of L and R scans are present in the dataset. Then, the majority of the scans view correspond to MLO and CC. This value corresponds to the orientation of the image ([source](https://two-views.com/mammograms/angles-views.html)):\n- MLO : Mediolateral oblique view (from side to side)\n- CC : Cranialcaudal (from above, the most comon view)\n- AT : *Information on AT view are yet to be found.*\n- LM : Lateromedial View (from the sides like the ML, only in reverse)\n- ML : Mediolateral View (from the center of the chest between the breasts, outward)\n- LMO : Lateral-medial oblique","metadata":{}},{"cell_type":"markdown","source":"The number of images by patient is looked more in details.","metadata":{}},{"cell_type":"code","source":"# The frequency of a given number of picture per patient\npic_num_freq = df_train['patient_id'].value_counts().value_counts()\n\nfig, ax = plt.subplots(figsize = (7,5))\n\nax.bar(x = pic_num_freq.index,\n       tick_label = pic_num_freq.index,\n       height= pic_num_freq,\n       color = color_bar(pic_num_freq.reset_index().patient_id),\n       log = True)\n\nax.set_title('Frequency of the number of picture per patient')\nax.set_xlabel('Number of picture per patient')\nax.set_ylabel('Frequency (log)')\nplt.show()\n\nprint('Values in a table :')\npic_num_freq","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:12:51.067832Z","iopub.execute_input":"2023-01-12T17:12:51.068963Z","iopub.status.idle":"2023-01-12T17:12:51.554025Z","shell.execute_reply.started":"2023-01-12T17:12:51.068921Z","shell.execute_reply":"2023-01-12T17:12:51.552848Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"For the vast majority of patients, 4 pictures where taken. According to the protocol descibred <a href= https://www.kaggle.com/competitions/rsna-breast-cancer-detection/discussion/369262>here</a>, the pictures could correspond to left and right CC and MLO.\n\nFor a lower proportion of patient a larger proportion of pictures have been taken. The number of patients seems to decrease exponetialy for an increase of picture number. A larger number of pictures could be related to a suspicion of cancer or to a difficulty of analysis (given in *difficult_negative_case*).","metadata":{"execution":{"iopub.status.busy":"2023-01-02T15:03:03.302552Z","iopub.execute_input":"2023-01-02T15:03:03.304167Z","iopub.status.idle":"2023-01-02T15:03:03.315507Z","shell.execute_reply.started":"2023-01-02T15:03:03.30408Z","shell.execute_reply":"2023-01-02T15:03:03.314461Z"}}},{"cell_type":"code","source":"df_patient = df_train[['patient_id',\n                       'site_id',\n                       'age', 'cancer',\n                       'implant',\n                       'difficult_negative_case']].groupby('patient_id').aggregate({'patient_id':'count',\n                                                                                    'site_id':'mean',\n                                                                                    'age' : 'mean',\n                                                                                    'cancer' : 'mean',\n                                                                                                                                 'implant' : 'mean',\n                                                                                                                                        'difficult_negative_case': 'mean'})\ndf_patient = df_patient.rename(columns = {'patient_id':'pic_num'})","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:17.306926Z","iopub.execute_input":"2023-01-12T17:08:17.307374Z","iopub.status.idle":"2023-01-12T17:08:17.32848Z","shell.execute_reply.started":"2023-01-12T17:08:17.307331Z","shell.execute_reply":"2023-01-12T17:08:17.327575Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_patient","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:17.329705Z","iopub.execute_input":"2023-01-12T17:08:17.330402Z","iopub.status.idle":"2023-01-12T17:08:17.350621Z","shell.execute_reply.started":"2023-01-12T17:08:17.330365Z","shell.execute_reply":"2023-01-12T17:08:17.349349Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def stacked_kde_plot(df, col1, col2):\n    '''\n    Create a stacked KDE plot where each plot correspond to\n    the frequency of col2 for each unique values of col1.\n    \n    If the number of unique values in col1 is above 6, the data are\n    bined.\n    '''\n    if df[col1].unique().shape[0] > 3:\n        bins = np.linspace(df_patient[col1].min(),\n                           df_patient[col1].max(),\n                           4)\n        y_data = pd.cut(df_patient[col1],\n                        bins,\n                        include_lowest = True)\n\n    else:\n        y_data = df[col1]\n        \n    fig, axes = plt.subplots(y_data.unique().shape[0], 1)\n    color = [c_0, c_1, c_2, c_3]\n    for index, unique in enumerate(y_data.unique()):\n        sns.kdeplot(df.loc[y_data == unique, col2],\n                    ax = axes[index],\n                    fill = True,\n                    color = color[index],\n                    alpha = 1)\n        sns.kdeplot(df.loc[y_data == unique, col2],\n                    ax = axes[index],\n                    color = 'w',\n                    lw = 3)\n\n        axes[index].get_xaxis().set_visible(False)\n        axes[index].get_yaxis().set_visible(False)\n        axes[index].spines[['right', 'top', 'left', 'bottom']].set_visible(False)\n        axes[index].set_facecolor((1,1,1,0))\n        axes[index].set_title(unique,\n                              loc = 'right',\n                              y = 0)\n        axes[index].set_xlim((df[col2].min()-1,\n                              df[col2].max()))\n\n    fig.subplots_adjust(hspace=-.25)\n    axes[y_data.unique().shape[0]-1].get_xaxis().set_visible(True)\n    fig.suptitle(f'{col2} kde plot depending on the {col1}')\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:14:38.875195Z","iopub.execute_input":"2023-01-12T17:14:38.875629Z","iopub.status.idle":"2023-01-12T17:14:38.889184Z","shell.execute_reply.started":"2023-01-12T17:14:38.875596Z","shell.execute_reply":"2023-01-12T17:14:38.888004Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for col in ['site_id', 'age', 'cancer', 'implant', 'difficult_negative_case']:\n    stacked_kde_plot(df_patient, col, 'pic_num')","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:14:40.238483Z","iopub.execute_input":"2023-01-12T17:14:40.239038Z","iopub.status.idle":"2023-01-12T17:14:42.145286Z","shell.execute_reply.started":"2023-01-12T17:14:40.239001Z","shell.execute_reply":"2023-01-12T17:14:42.143831Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The only parmeter that seems to highly affect the number of picture taken is the *implant*. Indeed, it seems logical that implants leads to further analysis. However, it is interessting that the number of pictures for difficult to analyse samples is only slightly higher than for other patients. The number of pictures is in average around 4 and thus, taking 14 (the maximum) pictures is exceptional.","metadata":{"execution":{"iopub.status.busy":"2023-01-03T09:46:43.992601Z","iopub.execute_input":"2023-01-03T09:46:43.993033Z","iopub.status.idle":"2023-01-03T09:46:44.001651Z","shell.execute_reply.started":"2023-01-03T09:46:43.992999Z","shell.execute_reply":"2023-01-03T09:46:44.000356Z"}}},{"cell_type":"markdown","source":"## 1.2 Pictures\n\nPictures present in the dataset are DICOM (Digital Imaging and Communications in Medicine) files. The DICOM files are specific to medical imaging and contains within the same file images as well as metadata about the patient and analysis. Again, the import of the data was inspired by the EDA performed by Laura Fink ([see](https://www.kaggle.com/code/allunia/rsna-breast-cancer-eda)).\n\nThe image processing have also been inspired of the tutorial of Laura Fink on Dicom images ([see](https://www.kaggle.com/code/allunia/pulmonary-dicom-preprocessing)) and on video explaining CT imaging [here](https://www.youtube.com/embed/KZld-5W99cI).","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-01-12T17:08:19.399469Z","iopub.execute_input":"2023-01-12T17:08:19.400076Z","iopub.status.idle":"2023-01-12T17:08:19.404906Z","shell.execute_reply.started":"2023-01-12T17:08:19.400041Z","shell.execute_reply":"2023-01-12T17:08:19.403711Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As a test the slices from a single patient are loaded.","metadata":{}},{"cell_type":"code","source":"test_slices = [pydicom.dcmread(train_path + '/'+ '10006'+'/'+filepath) for filepath in listdir(train_path + '/'+ '10006')]","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:19.406355Z","iopub.execute_input":"2023-01-12T17:08:19.407574Z","iopub.status.idle":"2023-01-12T17:08:19.43844Z","shell.execute_reply.started":"2023-01-12T17:08:19.407528Z","shell.execute_reply":"2023-01-12T17:08:19.437424Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The metadata are displayed","metadata":{}},{"cell_type":"code","source":"test_slices[0]","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:19.439515Z","iopub.execute_input":"2023-01-12T17:08:19.440485Z","iopub.status.idle":"2023-01-12T17:08:19.449331Z","shell.execute_reply.started":"2023-01-12T17:08:19.44045Z","shell.execute_reply":"2023-01-12T17:08:19.448284Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The images contained in the slices are plotted","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(2,2, figsize = (10,10))\nfor index, file in enumerate(test_slices):\n    img = file.pixel_array\n    axes[index//2,index%2].imshow(img, cmap=\"gray\")","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:19.450878Z","iopub.execute_input":"2023-01-12T17:08:19.451705Z","iopub.status.idle":"2023-01-12T17:08:33.307118Z","shell.execute_reply.started":"2023-01-12T17:08:19.451668Z","shell.execute_reply":"2023-01-12T17:08:33.305799Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"As explained by Laura Fink and in [this comment](https://www.kaggle.com/code/radek1/how-to-process-dicom-images-to-pngs/comments#2056978), Dicom scans contains raw data (as arrays) that should be converted to the Hounsfield scale (HU). The conversion is performed thanks to the parameters given in the metadata as:\n\n\n    HU = Pixel_value * Rescale Slope - Rescale Intercept\n    \nAfter a rapid study of the scans metadata, it seems that most of the scans have a Slope equal to one and an Intercept of 0.\n\nHowerver, to be sure that no sample is on another scale, a function to convert the images is created (the function was originaly created by Radek Osmulski [see](https://www.kaggle.com/code/radek1/eda-training-a-fast-ai-model-submission).","metadata":{}},{"cell_type":"code","source":"def rescale_img_to_hu(img, dicom):\n    \"\"\"Rescales the image to Hounsfield unit.\"\"\"\n    return img * dicom.RescaleSlope + dicom.RescaleIntercept","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:33.308601Z","iopub.execute_input":"2023-01-12T17:08:33.309071Z","iopub.status.idle":"2023-01-12T17:08:33.31401Z","shell.execute_reply.started":"2023-01-12T17:08:33.309039Z","shell.execute_reply":"2023-01-12T17:08:33.31312Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The samples presents different min and max values. The background value is given by the *Pixel Padding Value*.\n\nAccording to the study of Laura Fink, the background value depends on the machine on which the scan have been done. Moreover, on a single machine, the background can be different from a patient to the other.\n\nThis difference of background level can lead to difficulties in visualising the samples. Indeed dicom raw data have value varying with a wider range (about 1500) than the range of RGB picture (256). For a better interpretation, a specific window can be selected for increasing details in the scan. As an examples values bellow 500 can be set to 0, values above 2000 to 1 and the values between 500 and 2000 can be scaled between 0 and 1. This way, objects that are absorbing between 500 and 2000 present better readability.","metadata":{}},{"cell_type":"markdown","source":"As an example, the pixel values of the first sample have been investigated.","metadata":{}},{"cell_type":"code","source":"print(f'max value :, {np.max(test_slices[0].pixel_array)}')\nprint(f'min value :, {np.min(test_slices[0].pixel_array)}')","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:08:33.315383Z","iopub.execute_input":"2023-01-12T17:08:33.315881Z","iopub.status.idle":"2023-01-12T17:08:33.380716Z","shell.execute_reply.started":"2023-01-12T17:08:33.315849Z","shell.execute_reply":"2023-01-12T17:08:33.379433Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.histplot(test_slices[0].pixel_array.flatten(),log_scale = (False, True))\nplt.title('Histogram of the intensity of pixels in the grayscale scans')","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:16:07.234732Z","iopub.execute_input":"2023-01-12T17:16:07.235149Z","iopub.status.idle":"2023-01-12T17:16:26.039395Z","shell.execute_reply.started":"2023-01-12T17:16:07.235119Z","shell.execute_reply":"2023-01-12T17:16:26.03842Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"It seems that most of the image information is contained between 1500 and 3100. The background value seems to be at 3204.\n\nDifferent window value are displayed for this sample.","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(2,2, figsize = (10,10))\n\n\nfor index, windows in enumerate([(1504, 3204), (1504, 2250), (2000, 2500), (2500, 3204)]):\n    img = test_slices[0].pixel_array.copy()\n    img = (img - windows[0])/(windows[1]-windows[0])\n    img[img>1] = 1\n    img[img<0] = 0\n    img = 1 - img\n    mappable = axes[index//2,index%2].imshow(img, cmap=\"jet\")\n\nplt.colorbar(mappable)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:25:16.553158Z","iopub.execute_input":"2023-01-12T17:25:16.553564Z","iopub.status.idle":"2023-01-12T17:25:23.655832Z","shell.execute_reply.started":"2023-01-12T17:25:16.553532Z","shell.execute_reply":"2023-01-12T17:25:23.654893Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"According to the [documentation on Dicom files](https://dicom.innolitics.com/ciods/vl-endoscopic-image/vl-image/00281051), the center and window are mentioned in the file metadata as *Window Center* and *Window Width*. The optimal value can then be calculated as follow with center = c and width = w.\n  \n    if (x <= c - 0.5 - (w-1) /2), then y = ymin\n\n    else if (x > c - 0.5 + (w-1) /2), then y = ymax\n\n    else y = ((x - (c - 0.5)) / (w-1) + 0.5) * (ymax- ymin) + ymin\n\nAnother method for windowing the dicom files described [here](https://www.kaggle.com/code/radek1/how-to-process-dicom-images-to-pngs/comments#2061939) is to use the VOI LUT conversion. For this purpose, a function is directly available in the pydicom library :\n\n    from pydicom.pixel_data_handlers.util import apply_voi_lut\n\nFirst, a convertion using the dicom documentation and the different windows values is performed in order to compare the results.\n","metadata":{}},{"cell_type":"code","source":"fig, axes = plt.subplots(2,2, figsize = (10,10))\nfor index in range(len(test_slices[0].WindowCenter)):\n    img = test_slices[0].pixel_array.copy()\n    img = ((img - (test_slices[0].WindowCenter[index] - 0.5)) / (test_slices[0].WindowWidth[index]-1) + 0.5)\n    img[img>1] = 1\n    img[img<0] = 0\n    img = 1-img\n    mappable = axes[index//2,index%2].imshow(img, cmap=\"jet\")\n\nplt.colorbar(mappable)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:25:55.265557Z","iopub.execute_input":"2023-01-12T17:25:55.265977Z","iopub.status.idle":"2023-01-12T17:26:02.435912Z","shell.execute_reply.started":"2023-01-12T17:25:55.265944Z","shell.execute_reply":"2023-01-12T17:26:02.434687Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"When selecting the first window and width a better view of the internal tissues can be obtained.\nRandom pictures from the dataset are displayed with and without the first window parameters.","metadata":{}},{"cell_type":"code","source":"# To ease import process, a function for importing image for a given patient is created\n\ndef load_scans(path, patient_id):\n    '''\n    Function for loading a given patien scans in slices\n    '''\n    dcm_path = path + str(patient_id)\n    slices = [pydicom.dcmread(dcm_path + '/' + file) for file in listdir(dcm_path)]\n    return slices","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:09:11.998933Z","iopub.execute_input":"2023-01-12T17:09:11.999413Z","iopub.status.idle":"2023-01-12T17:09:12.006682Z","shell.execute_reply.started":"2023-01-12T17:09:11.999366Z","shell.execute_reply":"2023-01-12T17:09:12.00538Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, axes = plt.subplots(6,2, figsize = (15, 20))\npatient_id_sample = df_train['patient_id'].sample(6, random_state= 0)\n\nfor index, patient_id in enumerate(patient_id_sample):\n    slices = load_scans(train_path, patient_id)\n    img = slices[0].pixel_array.copy()\n    img_win = img.copy()\n    \n    # Window center and width can contain multiple values\n    img_win = apply_voi_lut(img_win, slices[0])\n    \n    if slices[0].PhotometricInterpretation == \"MONOCHROME1\":\n        img = np.max(img)-img\n        img_win = 1-img_win\n        \n    axes[index,0].imshow(img, cmap=\"jet\")\n    axes[index,0].set_title(f'Patient {patient_id} without window correction')\n    \n    axes[index,1].imshow(img_win, cmap=\"jet\")\n    axes[index,1].set_title(f'Patient {patient_id} with window correction')\n\nplt.tight_layout()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:09:12.008299Z","iopub.execute_input":"2023-01-12T17:09:12.009303Z","iopub.status.idle":"2023-01-12T17:09:33.305243Z","shell.execute_reply.started":"2023-01-12T17:09:12.009243Z","shell.execute_reply":"2023-01-12T17:09:33.304056Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"In general the window correction seems to allow a better visualisation of the breasts internal tissues.","metadata":{}},{"cell_type":"markdown","source":"# 2. Data preprocessing\n\nAs described in the image-focused part of the EDA, the data should be preprocessed through the following steps :\n- Convertion of the dicom file to an array\n- Windowing of the values with the Voi lut convertion\n- Scaling of the raw values between 0 and 1 or 0 and 255\n- Re-invert inverted scans\n\nAfter these steps, the convertion of the images to a lighter format as PNG could greatly improve the ML models performances. As a lot of scans present empty regions, the centering of the pictures on the objet of interest could also be highly beneficial.","metadata":{}},{"cell_type":"markdown","source":"## 2.1 PNG conversion\n\nStrong of all the previous observations, a function is created for converting the dicom files to lighter PNG files. This procedure is inspired from the notebooks of Theo Viel ([see](https://www.kaggle.com/code/theoviel/dicom-resized-png-jpg)) and Radek Osmulski ([see](https://www.kaggle.com/code/radek1/how-to-process-dicom-images-to-pngs/notebook)).","metadata":{}},{"cell_type":"code","source":"def process(path, size=512, save_folder=\"output/\", extension=\"png\"):\n    # Get patient_id and image id\n    patient = path.split('/')[-2]\n    image = path.split('/')[-1][:-4]\n    \n    # The scan is read and windiwing is applied\n    dicom = pydicom.dcmread(path)\n    img = apply_voi_lut(dicom.pixel_array, dicom)\n    \n    # The image is rescaled between 0 and 1\n    img = (img - img.min()) / (img.max() - img.min())\n    \n    # If the scan is negative, the value are inverted\n    if dicom.PhotometricInterpretation == \"MONOCHROME1\":\n        img = 1 - img\n    \n    # Resize the scan to a png\n    img = cv2.resize(img, (size, size))\n    \n    img = (img * 255).astype(np.uint8)\n    \n    # Save the image (255 based)\n    # cv2.imwrite(save_folder + f\"{patient}_{image}.{extension}\", img)\n    \n    return img","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:09:33.307076Z","iopub.execute_input":"2023-01-12T17:09:33.307424Z","iopub.status.idle":"2023-01-12T17:09:33.315934Z","shell.execute_reply.started":"2023-01-12T17:09:33.307393Z","shell.execute_reply":"2023-01-12T17:09:33.315022Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Convertion of the whole dataset requiers a long porcessing time due to the time requiered to open the heavy dicom files and to the high number of files (above 54000). Moreover, this export have already been performed by Radek Osmulski ([see](https://www.kaggle.com/datasets/radek1/rsna-mammography-images-as-pngs)) and thus his dataset will directly be used for the next steps of this study.\n\nIn accordance with the observation present in the EDA, the dataset selected is the one using pydicom, with VOI lUT conversion and a size of 512x512 pixels.","metadata":{}},{"cell_type":"markdown","source":"In order to ease the procedure of importing the png files, a function is created.","metadata":{}},{"cell_type":"code","source":"def import_png(patient_id):\n    train_png_path = '/kaggle/input/rsna-mammography-images-as-pngs/images_as_pngs_cv2_vl_512/train_images_processed_cv2_vl_512/'\n    path = train_png_path + str(patient_id)\n    slice_path_list = [path + '/' + img_fn for img_fn in os.listdir(path)]\n    slices = [cv2.imread(slice_path) for slice_path in slice_path_list]\n    \n    return slices","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:09:33.317371Z","iopub.execute_input":"2023-01-12T17:09:33.317695Z","iopub.status.idle":"2023-01-12T17:09:33.328528Z","shell.execute_reply.started":"2023-01-12T17:09:33.317666Z","shell.execute_reply":"2023-01-12T17:09:33.327567Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"patient_id_sample = df_train['patient_id'].sample(3, random_state = 42).to_numpy()\n\nfig, axes = plt.subplots(3,2, figsize = (10, 10))\n\nfor i in range(3):\n    train_png_path = '/kaggle/input/rsna-mammography-images-as-pngs/images_as_pngs_cv2_vl_512/train_images_processed_cv2_vl_512/'\n    img_path = train_png_path + str(patient_id_sample[i])\n    img_fn = os.listdir(img_path)[0].replace('png', 'dcm')\n    \n    \n    patient_path = train_path + str(patient_id_sample[i])\n    img_path = patient_path + '/' + img_fn\n    \n    axes[i, 0].imshow(process(img_path), cmap = 'gray')\n    axes[i, 0].set_title('Raw image processed')\n    \n    mappable = axes[i, 1].imshow(import_png(patient_id_sample[i])[0], cmap = 'Greys')\n    axes[i, 1].set_title('Imported PNG')\n    \nfig.colorbar(mappable)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-12T17:28:59.723098Z","iopub.execute_input":"2023-01-12T17:28:59.723524Z","iopub.status.idle":"2023-01-12T17:29:04.091352Z","shell.execute_reply.started":"2023-01-12T17:28:59.723487Z","shell.execute_reply":"2023-01-12T17:29:04.090349Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}