{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":18613,"sourceType":"datasetVersion","datasetId":5839}],"dockerImageVersionId":30588,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Classification de radiographies pulmonaires à l'aide de réseaux de neurone convolutifs.\n**Gilbert Jabbour et Fayçal Zine-Eddine**\n\nLe but de ce projet est de classifier des radiographies pulmonaires à l'aide d'apprentissage profond.\n\nNous nous sommes intéressés plus particulièrement à la détection de pneumonie et d'oedème pulmonaire, car c'est des conditions vues fréquemment en clinique, qu'on veut généralement traiter rapidement et qui sont régulièrement confondues. \n\nCe *notebook* détaille la visualisation des données, la transformation et l'extraction des caractéristiques pertinentes et la production des jeux de données qui seront utilisées par notre algorithme.\n\nNous avons choisi de présenter ce notebook directement sur Kaggle, pour éviter de devoir télécharger la large base de données du NIH pour laquelle nous allons seulement utiliser une fraction des images.\n\nLe reste du code est accessible sur Github : https://github.com/Fay-Zine/NIH-Chest-XRay","metadata":{"execution":{"iopub.status.busy":"2023-12-01T11:12:02.886825Z","iopub.execute_input":"2023-12-01T11:12:02.887239Z","iopub.status.idle":"2023-12-01T11:12:03.151972Z","shell.execute_reply.started":"2023-12-01T11:12:02.88721Z","shell.execute_reply":"2023-12-01T11:12:03.151002Z"}}},{"cell_type":"markdown","source":"# Import and data loading","metadata":{}},{"cell_type":"code","source":"import os\nfrom glob import glob\nimport io\nfrom itertools import chain\nfrom collections import Counter\nimport pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport torch\nfrom torch.utils.data import Dataset, DataLoader\nfrom torchvision import transforms\nimport torch.nn as nn\nimport torchvision.models as models\nimport torch.optim as optim\nfrom PIL import Image\nfrom sklearn.model_selection import GroupShuffleSplit\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.metrics import classification_report\nimport pyarrow as pa\nimport pyarrow.parquet as pq\nimport warnings\nimport copy\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\n\n### Load dataset and visualize classes\nDATA_FILE_PATH = \"../input/data/Data_Entry_2017.csv\"\nBASE_IMAGE_DIR = \"/kaggle/input/data\"\n\nfull_df = pd.read_csv(DATA_FILE_PATH)\nfull_df = full_df[full_df['Patient Age'] <= 100]  #Présence d'outliers lors de la visualisation du df\nfull_df.describe()","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:17.471132Z","iopub.execute_input":"2023-12-06T13:19:17.471481Z","iopub.status.idle":"2023-12-06T13:19:20.423958Z","shell.execute_reply.started":"2023-12-06T13:19:17.471452Z","shell.execute_reply":"2023-12-06T13:19:20.422977Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Data visualization","metadata":{}},{"cell_type":"markdown","source":"## Classes distribution\nOn observe qu'il y a un débalancement majeur entre les classes. \nPour nos maladies d'intérêt, il y a environ 4000  images sur les 120 000 disponibles dans le dataset.","metadata":{}},{"cell_type":"code","source":"all_labels = np.unique(list(chain(*full_df['Finding Labels'].map(lambda x: x.split('|')).tolist())))\nall_labels = [x for x in all_labels if len(x) > 0]\n\nfor c_label in all_labels:\n    if len(c_label) > 0:\n        full_df[c_label] = full_df['Finding Labels'].map(lambda finding: 1.0 if c_label in finding else 0)\n\nlabel_counts = (\n    full_df[\"Finding Labels\"].str.split(\"|\", expand=True).stack().value_counts()\n)\nprint(label_counts)\nlabel_columns = [label for label in all_labels if label in full_df.columns]\nlabel_df = full_df[label_columns]\n\ncorrelation_matrix = label_df.corr()\n\n# Overall Label Counts\nplt.figure(figsize=(12, 8))\nsns.barplot(y=label_counts.index, x=label_counts.values)\nplt.title(\"Frequency of Each Label in the Dataset\")\nplt.xlabel(\"Counts\")\nplt.ylabel(\"Labels\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:20.425754Z","iopub.execute_input":"2023-12-06T13:19:20.426096Z","iopub.status.idle":"2023-12-06T13:19:22.660476Z","shell.execute_reply.started":"2023-12-06T13:19:20.426064Z","shell.execute_reply":"2023-12-06T13:19:22.659578Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Classes correlations\nÉtant donné que plusieurs \"label\" comme la cardiomégalie ou l'épanchement peuvent être associé à de l'insuffisance cardiaque, nous avons exploré l'interdépendance des maladies pour voir si des classes étaient fortement liées et pouvaient être utilisées pour améliorer notre modèle.\n\nAucun des coefficients n'était particulièrement haut, le plus grand étant la combinaison des deux maladies d'intérêt. Les probabilités conditionnelles variaient, mais les inclure aurait grandement augmenté la complexité du modèle avec un gain qu'on suspecte limité sur les performances.","metadata":{}},{"cell_type":"code","source":"# Calculate correlations for Pneumonia and Edema\ncorrelations_pneumonia = {label: correlation_matrix.at[label, 'Pneumonia'] for label in label_columns if 'Pneumonia' in label_columns}\ncorrelations_edema = {label: correlation_matrix.at[label, 'Edema'] for label in label_columns if 'Edema' in label_columns}\n\n# Visualize correlations with Pneumonia\ncorrelation_df_pneumonia = pd.DataFrame.from_dict(correlations_pneumonia, orient='index', columns=['Correlation']).drop(['Pneumonia']).drop(['No Finding'])\ncorrelation_df_pneumonia = correlation_df_pneumonia.reset_index().rename(columns={'index': 'Condition'})\ncorrelation_df_pneumonia = correlation_df_pneumonia.sort_values(by='Correlation', ascending=False)\n\nplt.figure(figsize=(12, 8))\nsns.barplot(x=correlation_df_pneumonia['Correlation'], y=correlation_df_pneumonia['Condition'], palette=\"RdYlGn_r\")\nplt.title(\"Correlations with Pneumonia\")\nplt.xlabel(\"Correlation\")\nplt.ylabel(\"Condition\")\nplt.show()\n\n# Visualize correlations with Edema\ncorrelation_df_edema = pd.DataFrame.from_dict(correlations_edema, orient='index', columns=['Correlation']).drop(['Edema']).drop(['No Finding'])\ncorrelation_df_edema = correlation_df_edema.reset_index().rename(columns={'index': 'Condition'})\ncorrelation_df_edema = correlation_df_edema.sort_values(by='Correlation', ascending=False)\n\nplt.figure(figsize=(12, 8))\nsns.barplot(x=correlation_df_edema['Correlation'], y=correlation_df_edema['Condition'], palette=\"RdYlGn_r\")\nplt.title(\"Correlations with Pulmonary Edema\")\nplt.xlabel(\"Correlation\")\nplt.ylabel(\"Condition\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:22.661714Z","iopub.execute_input":"2023-12-06T13:19:22.662092Z","iopub.status.idle":"2023-12-06T13:19:23.43493Z","shell.execute_reply.started":"2023-12-06T13:19:22.662054Z","shell.execute_reply":"2023-12-06T13:19:23.434042Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Target labelling\nPour notre problème de différentiation entre l'oedème pulmonaire et les pneumonies, plusieurs approches peuvent être utilisés. \n\nNous avons réfléchi notamment à utiliser une approche *un contre tous* où on entraîne 2 modèles à détecter chaque maladie, mais cela entraîne beaucoup de overload et de complexité dans l'application pratique.\n\nLa classification multiclasse est plus directe et nous semblait plus intéressante. \n\nNous avons initialement séparé notre problème de classification en 5 classes :  \"Pneumonia\", \"Edema\", \"Both\", \"Normal\" and \"Other\", mais il y avait un trop gros débalancement pour que le modèle performe bien sur notre tâche initiale. Nous avons testé une fonction de perte *weighted*, un *undersampling* massif des catégories dominantes avec un oversampling des autres, mais même à cela des tests préliminaires donnait une sensibilité <50% pour toutes nos classes d'intérêts.\n\nNous avons donc opté pour une classification plus simple à 3 classes : \"Pneumonia\", \"Edema\" et \"Both\".","metadata":{}},{"cell_type":"code","source":"def processing(df, label_col, conditions, undersample_size=2000):\n    \"\"\"\n    Categorizes labels, undersamples 'Normal' and 'Other' classes.\n    \"\"\"\n    def categorize_label(label):\n        if conditions[0] in label and conditions[1] in label:\n            return \"Both\"\n        elif conditions[0] in label:\n            return \"Pneumonia\"\n        elif conditions[1] in label:\n            return \"Edema\"\n        elif \"No Finding\" in label:\n            return \"Normal\"\n        else:\n            return \"Other\"\n\n    df[\"Target\"] = df[label_col].apply(categorize_label)  \n    return df\n\nprocessed_df = processing(full_df, \"Finding Labels\", [\"Pneumonia\", \"Edema\"])\nclasses = processed_df[\"Target\"].unique().tolist()\n\nclass_counts = processed_df[\"Target\"].value_counts()\n\n# Classes in final dataset\nplt.figure(figsize=(8, 6))\nplt.bar(class_counts.index, class_counts.values)\nplt.xlabel(\"Classes\")\nplt.ylabel(\"Count\")\nplt.title(\"Class Distribution in processed_df\")\nplt.xticks(rotation=45)  \nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:23.437953Z","iopub.execute_input":"2023-12-06T13:19:23.438358Z","iopub.status.idle":"2023-12-06T13:19:23.712897Z","shell.execute_reply.started":"2023-12-06T13:19:23.438322Z","shell.execute_reply":"2023-12-06T13:19:23.711964Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Séparation des groupes ","metadata":{}},{"cell_type":"code","source":"# Séparation des groupes basé sur l'identifiant patient pour éviter le \"data leakage\".\npatient_col = \"Patient ID\"\n\ntest_size = 0.2\nvalid_size = 0.1\n\ngss = GroupShuffleSplit(n_splits=1, test_size=test_size + valid_size, random_state=50)\n\n# Séparation en donnée d'entraînement, de validation et de test\nfor train_idx, test_valid_idx in gss.split(processed_df, groups=processed_df[patient_col]):\n    train_df = processed_df.iloc[train_idx]\n    test_valid_df = processed_df.iloc[test_valid_idx]\n\ngss_valid = GroupShuffleSplit(\n    n_splits=1, test_size=test_size / (test_size + valid_size), random_state=42\n)\n\nfor valid_idx, test_idx in gss_valid.split(\n    test_valid_df, groups=test_valid_df[patient_col]\n):\n    valid_df = test_valid_df.iloc[valid_idx]\n    test_df = test_valid_df.iloc[test_idx]\n    \ndef plot_class_distributions(dataframes, titles):\n    \"\"\"\n    Plots the class distributions of multiple DataFrames side by side.\n    \"\"\"\n    plt.figure(figsize=(15, 6))\n\n    num_dfs = len(dataframes)\n\n    for i, df in enumerate(dataframes):\n        plt.subplot(1, num_dfs, i+1)\n        class_counts = df['Target'].value_counts()\n        sns.barplot(x=class_counts.values, y=class_counts.index)\n        plt.title(titles[i])\n        plt.xlabel('Counts')\n        if i == 0:\n            plt.ylabel('Target Classes')\n        else:\n            plt.ylabel('')\n\n    plt.tight_layout()\n    plt.show()\n\nplot_class_distributions(\n    [train_df, valid_df, test_df], \n    ['Training Set Class Distribution', 'Validation Set Class Distribution', 'Test Set Class Distribution']\n)\n","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:23.71397Z","iopub.execute_input":"2023-12-06T13:19:23.714249Z","iopub.status.idle":"2023-12-06T13:19:24.37673Z","shell.execute_reply.started":"2023-12-06T13:19:23.714224Z","shell.execute_reply":"2023-12-06T13:19:24.375799Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Vu le débalancement extrême, on garde seulement un problème à 3 classes\ntrain_df_disease = train_df[train_df['Target'].isin([\"Edema\", \"Pneumonia\", \"Both\"])]\nval_df_disease = valid_df[valid_df['Target'].isin([\"Edema\", \"Pneumonia\", \"Both\"])]\ntest_df_disease = test_df[test_df['Target'].isin([\"Edema\", \"Pneumonia\", \"Both\"])]","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:24.377939Z","iopub.execute_input":"2023-12-06T13:19:24.378268Z","iopub.status.idle":"2023-12-06T13:19:24.397307Z","shell.execute_reply.started":"2023-12-06T13:19:24.378239Z","shell.execute_reply":"2023-12-06T13:19:24.396399Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Visualisation des caractéristiques démographiques et cliniques","metadata":{}},{"cell_type":"code","source":"def plot_distributions(dataframes, titles, columns):\n    \"\"\"\n    Plots the distributions of specified columns for multiple DataFrames.\n    \"\"\"\n    num_dfs = len(dataframes)\n    num_cols = len(columns)\n\n    plt.figure(figsize=(15, 5 * num_cols))\n\n    for i, df in enumerate(dataframes):\n        for j, col in enumerate(columns):\n            plt.subplot(num_cols, num_dfs, num_dfs * j + i + 1)\n            if df[col].dtype == 'object':\n                sns.countplot(y=df[col], order = df[col].value_counts().index)\n            else:\n                sns.histplot(df[col], bins=20, kde=False)\n                mean_val = df[col].mean()\n                std_val = df[col].std()\n                plt.axvline(mean_val, color='k', linestyle='dashed', linewidth=1)\n                plt.axvline(mean_val - std_val, color='r', linestyle='dashed', linewidth=1)\n                plt.axvline(mean_val + std_val, color='r', linestyle='dashed', linewidth=1)\n                min_ylim, max_ylim = plt.ylim()\n                min_xlim, max_xlim = plt.xlim()\n                plt.text(max_xlim-20, max_ylim*0.9, 'Mean: {:.2f}'.format(mean_val))\n                plt.text(max_xlim-20, max_ylim*0.85, 'Std: {:.2f}'.format(std_val))\n            \n            if j == 0:\n                plt.title(titles[i])\n            \n            if i == 0:\n                plt.ylabel(col)\n            else:\n                plt.ylabel('')\n\n    plt.tight_layout()\n    plt.show()\n\ncolumns_to_plot = ['Patient Age', 'Patient Gender', 'View Position']\n\nplot_distributions(\n    [train_df_disease, val_df_disease, test_df_disease], \n    ['Training Set', 'Validation Set', 'Test Set'],\n    columns_to_plot\n)","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:24.398494Z","iopub.execute_input":"2023-12-06T13:19:24.39888Z","iopub.status.idle":"2023-12-06T13:19:26.210939Z","shell.execute_reply.started":"2023-12-06T13:19:24.398849Z","shell.execute_reply":"2023-12-06T13:19:26.20998Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Rebalancement des classes","metadata":{}},{"cell_type":"code","source":"def create_inverse_weighted_sample(df, target_column, sample_size):\n    # Calculate class weights (inverse of class frequencies)\n    class_counts = df[target_column].value_counts()\n    total_samples = df.shape[0]\n    class_weights = {cls: total_samples/count for cls, count in class_counts.items()}\n\n    # Assign a weight to each row in the DataFrame\n    weights = df[target_column].map(class_weights)\n\n    # Perform weighted random sampling\n    return df.sample(n=sample_size, weights=weights, replace=True)\n\nequalized_train_df_disease = create_inverse_weighted_sample(train_df_disease, 'Target', sample_size=2500)\nequalized_val_df_disease = create_inverse_weighted_sample(train_df_disease, 'Target', sample_size=500)\n\nplot_class_distributions(\n    [equalized_train_df_disease, equalized_val_df_disease, test_df_disease], \n    ['Training Set Class Distribution', 'Validation Set Class Distribution', 'Test Set Class Distribution'])","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:26.212134Z","iopub.execute_input":"2023-12-06T13:19:26.212444Z","iopub.status.idle":"2023-12-06T13:19:26.840841Z","shell.execute_reply.started":"2023-12-06T13:19:26.212415Z","shell.execute_reply":"2023-12-06T13:19:26.839855Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"base_image_dir = \"/kaggle/input/data\"\n\n# Build the dictionary\nimage_folder_path = {\n    os.path.basename(x): x\n    for x in glob(os.path.join(base_image_dir, \"images*\", \"images\", \"*.png\"))\n}\n\n# Define a function to load and transform images\ndef load_transform_image(image_path, transform):\n    image = Image.open(image_path)\n    image = transform(image)\n    return image\n\n# Define a function to save images and labels to a Parquet file\ndef save_images_to_parquet(dataset, filename):\n    table = pa.Table.from_pandas(dataset)\n    pq.write_table(table, filename)\n    \n# Define the transforms for each dataset\ntrain_transform = transforms.Compose([\n    transforms.Resize((256, 256)), \n    transforms.RandomResizedCrop((224, 224)),  \n    transforms.RandomRotation(5),  \n    transforms.RandomAffine(degrees=0, translate=(0.05, 0.05), scale=(0.95, 1.05)),\n    transforms.ToTensor()\n])    \n    \nval_test_transform = transforms.Compose([\n    transforms.Resize((224, 224)),  \n    transforms.ToTensor(),\n])\n\n# Function to process and save dataset\ndef process_and_save_dataset(df, image_paths_dict, transform, filename1, filename2):\n    images = []\n    labels = []\n    \n    for _, row in df.iterrows():\n        img_path = image_paths_dict.get(row['Image Index'])\n        if img_path:\n            image = load_transform_image(img_path, transform)\n            image_np = np.array(image).flatten()/255 #Scale 0-1\n            images.append(image_np)\n            labels.append(row['Target'])\n           \n    images = pd.DataFrame({\n        'images': images })\n    labels = pd.DataFrame({\n        'Target': labels })\n        \n    save_images_to_parquet(images, filename1)\n    save_images_to_parquet(labels, filename2)","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:22:02.958542Z","iopub.execute_input":"2023-12-06T13:22:02.958903Z","iopub.status.idle":"2023-12-06T13:22:03.500114Z","shell.execute_reply.started":"2023-12-06T13:22:02.95887Z","shell.execute_reply":"2023-12-06T13:22:03.49932Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Visualisation de radiographies pulmonaires après transformation","metadata":{}},{"cell_type":"code","source":"def plot_images(df, image_paths_dict, figsize=(20, 20)):\n\n    grouped = df.groupby('Target')\n    fig, axs = plt.subplots(len(grouped), 5, figsize=figsize) \n\n    for i, (target, group) in enumerate(grouped):\n        random_samples = group.sample(n=5, replace=False)  \n        for j in range(5):\n            img_name = random_samples.iloc[j]['Image Index']\n            img_path = image_paths_dict.get(img_name)\n            if img_path:\n                image = Image.open(img_path).convert(\"RGB\")\n                image = np.asarray(image)\n                axs[i, j].imshow(image)\n                axs[i, j].set_title(f\"Target: {target}\")\n                axs[i, j].axis('off')\n\n    plt.tight_layout()\n    plt.show()\n\nplot_images(processed_df, image_folder_path)","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:19:27.398488Z","iopub.execute_input":"2023-12-06T13:19:27.398819Z","iopub.status.idle":"2023-12-06T13:19:34.636624Z","shell.execute_reply.started":"2023-12-06T13:19:27.398789Z","shell.execute_reply":"2023-12-06T13:19:34.635478Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Show characteristics of images before saving, used to ensure we have the correct shape as input\nimage_np = load_transform_image(\"/kaggle/input/data/images_001/images/00000001_000.png\", val_test_transform)  \n\nimage_shape = image_np.shape  \nimage_dtype = image_np.dtype  \n\nprint(\"Image Shape:\", image_shape)\nprint(\"Data Type:\", image_dtype)","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:22:07.482938Z","iopub.execute_input":"2023-12-06T13:22:07.483302Z","iopub.status.idle":"2023-12-06T13:22:07.505196Z","shell.execute_reply.started":"2023-12-06T13:22:07.483273Z","shell.execute_reply":"2023-12-06T13:22:07.504289Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Extraction des données sous format tabulaire en parquet.\nLe jeu de données pré-séparé et transformé est enregistré sous 6 fichiers parquet.\nX_train, y_train, X_val, y_val, X_test et y_test.\n\nChaque rangée des jeux de données X correspond à une image qui a été transformé, réduite en 224x224 (format communément utilisé en processing d'images) et aplatie en une seule rangée de 50 176 valeurs donnée sous forme de liste.","metadata":{}},{"cell_type":"code","source":"process_and_save_dataset(equalized_train_df_disease, image_folder_path, train_transform, \"X_train.parquet\", \"y_train.parquet\")\nprocess_and_save_dataset(equalized_val_df_disease, image_folder_path, val_test_transform, \"X_val.parquet\", \"y_val.parquet\")\nprocess_and_save_dataset(test_df_disease, image_folder_path, val_test_transform, \"X_test.parquet\", \"y_test.parquet\")","metadata":{"execution":{"iopub.status.busy":"2023-12-06T13:22:11.25299Z","iopub.execute_input":"2023-12-06T13:22:11.253937Z","iopub.status.idle":"2023-12-06T13:23:43.546265Z","shell.execute_reply.started":"2023-12-06T13:22:11.253885Z","shell.execute_reply":"2023-12-06T13:23:43.545196Z"},"trusted":true},"execution_count":null,"outputs":[]}]}