{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Localization of Findings\n\nIn this competition, the goal is to examine a chest x-ray image and plot a bounding box around findings of significance, along with a confidence rating. There are 13 classes that are of clinical significance with bounding boxes that span more than a single pixel in the chest x-ray. There is one class (`No finding`) that means that no signficiant finding was observed that has a bounding box that spans a single pixel. \n\nThere are many approaches to solving this problem, however, I believe the first step is to determine what the bounding box information actually looks like for each class. In particular, the purpose of this notebook is to see if bounding boxes for various target classes can be partially localized to one or more non-random areas on the chest x-ray. In this notebook we will:\n\n* Load the training data and augment it with information available from the dicom files.\n* Extract bounding box information for each class and overlay each on a heatmap that approximates the dimensions of a chest x-ray.\n* Examine patterns for each class of clinical finding.\n\nThe hope with this exercise is that we can reveal patterns as to where various clinical findings are likely to be found on a chest x-ray. We may be able to use these localizations to build classifiers that target those areas only, instead of concentrating on the x-ray image as a whole."},{"metadata":{},"cell_type":"markdown","source":"# 1. Reading and Augmenting Training Data\n\nIn this section, we will read in the training data, and then augment it with image meta-data from the associated dicom files. First, let's load up the training CSV file as we normally would."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import pandas as pd\nfrom sklearn.cluster import KMeans\n\ntrain_aug = pd.read_csv('../input/vinbigdata-chest-xray-abnormalities-detection/train.csv')\ntrain_aug","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Our next step is to pull in some basic information about each dicom file. In particular, what we're interested in is the height and width of the image. This gives us context for our bounding boxes. We'll use `pydicom` to read the image file and add some useful information to our training data. Let's leave out anything related to `class_id` 14, which is our `No finding` class. This will make loading meta-data faster since we don't have to load files for classes with no bounding boxes."},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"from pathlib import Path\nfrom pydicom import dcmread\n\ndef add_image_dimensions_gender(df):\n    path_spec = \"../input/vinbigdata-chest-xray-abnormalities-detection/train/{}.dicom\"\n    height = []\n    width = []\n    gender = []\n    age = []\n    for _, row in df.iterrows():\n        dcm = dcmread(Path(path_spec.format(row[\"image_id\"])), stop_before_pixels=True)\n        height.append(dcm.Rows)\n        width.append(dcm.Columns)\n        gender.append(dcm[0x10, 0x40].value)\n        patient_age = dcm[0x10, 0x1010].value if [0x10, 0x1010] in dcm else \"\"\n        age.append(patient_age)\n    df[\"image_height\"] = height\n    df[\"image_width\"] = width\n    df[\"gender\"] = gender\n    df[\"age\"] = age\n\nadd_image_dimensions_gender(train_aug)\ntrain_aug","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 2. Normalize Bounding Box Values\n\nNow that we have width and height information, we have context for the bounding box information we find in the rest of the training dataset. Our next step is to normalize the location of the bounding boxes. Why? Because each image may be a different size. We want to know relatively speaking where each bounding box appears in an idealized chest x-ray image."},{"metadata":{"trusted":true},"cell_type":"code","source":"def scale_bounding_boxes(df):\n    df[\"x_min_norm\"] = df[\"x_min\"] / df[\"image_width\"]\n    df[\"y_min_norm\"] = df[\"y_min\"] / df[\"image_height\"]\n    df[\"x_max_norm\"] = df[\"x_max\"] / df[\"image_width\"]\n    df[\"y_max_norm\"] = df[\"y_max\"] / df[\"image_height\"]\n    \nscale_bounding_boxes(train_aug)\ntrain_aug","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 3. Density Heatmaps\n\nNow we can create a density heatmap for each class of findings. We'll make our heatmap size have an aspect ratio of 4:5, which matches what we would expect on a chest x-ray image. This will help us localize where we expect to see each class of finding."},{"metadata":{"trusted":true},"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\n\ndef draw_box_on_array(row, np_array):\n    xy = [\n        int(row[\"x_min_norm\"] * 400), \n        int(row[\"y_min_norm\"] * 500),\n        int(row[\"x_max_norm\"] * 400),\n        int(row[\"y_max_norm\"] * 500),\n    ]\n    np_array[xy[1]:xy[3], xy[0]:xy[2]] += 1\n    \ndef get_bbox(df, class_id):\n    np_array = np.zeros(shape=(500, 400))\n    if class_id == 14:\n        return np_array\n    for _, row in df[df[\"class_id\"] == class_id].iterrows():\n        draw_box_on_array(row, np_array)\n    return np_array","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now let's loop through each class of finding and generate a heatmap."},{"metadata":{"trusted":true},"cell_type":"code","source":"classes = [\n    \"0 - Aortic enlargement\",\n    \"1 - Atelectasis\",\n    \"2 - Calcification\",\n    \"3 - Cardiomegaly\", \n    \"4 - Consolidation\",\n    \"5 - ILD\",\n    \"6 - Infiltration\",\n    \"7 - Lung Opacity\",\n    \"8 - Nodule/Mass\",\n    \"9 - Other lesion\",\n    \"10 - Pleural effusion\",\n    \"11 - Pleural thickening\",\n    \"12 - Pneumothorax\",\n    \"13 - Pulmonary fibrosis\",\n    \"14 - No finding\",\n]\n\nfig, axs = plt.subplots(nrows=5, ncols=3, figsize=(20, 30))\nfor i, ax in enumerate(axs.flatten()):\n    ax.imshow(get_bbox(train_aug, i), cmap='hot', interpolation='nearest')\n    _ = ax.set_title(classes[i], fontweight=\"bold\", size=15)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 3.1 Localization by Gender"},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=14, ncols=3, figsize=(20, 60))\n\nbboxes = [\n    get_bbox(train_aug[(train_aug[\"gender\"] == gender)], i)\n    for i in range(14)\n    for gender in [\"M\", \"F\", \"O\"]\n]\n\nnames = [\n    \"{} ({})\".format(classes[i], gender) \n    for i in range(14) \n    for gender in [\"M\", \"F\", \"O\"]\n]\n\nfor i, (name, bbox, ax) in enumerate(zip(names, bboxes, axs.flatten())):\n    ax.imshow(bbox, cmap='hot', interpolation='nearest')\n    _ = ax.set_title(name, fontweight=\"bold\", size=15)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 3.3 Per Class Exploration\n\nWe should take a closer look at some of the classes to see where their bounding boxes originate, and see if there are any outliers. A simple way of looking for outliers is to take the center point of the bounding box and see how many standard deviations away from the mean it falls. This will give us a sense of how many exceptional points we have to deal with when trying to fit our model."},{"metadata":{"trusted":true},"cell_type":"code","source":"import seaborn as sns\n\nfrom scipy.spatial import distance\nfrom scipy import stats\n\ndef calculate_centers(df):\n    x_centers = []\n    y_centers = []\n    for _, row in df.iterrows():\n        x_center = (row[\"x_min_norm\"] + ((row[\"x_max_norm\"] - row[\"x_min_norm\"]) / 2)) * 400\n        y_center = (row[\"y_min_norm\"] + ((row[\"y_max_norm\"] - row[\"y_min_norm\"]) / 2)) * 500\n        x_centers.append(x_center)\n        y_centers.append(y_center)\n    df[\"x_center\"] = x_centers\n    df[\"y_center\"] = y_centers\n\ndef get_probability_descriptor(x):\n    if x < 2.0:\n        return \"95% probability\"\n    if x >= 2.0 and x < 3.0:\n        return \"5% probability\"\n    if x >= 3.0 and x < 4.0:\n        return \"1% probability\"\n    return \"< 1% probability\"\n\ndef calculate_distances_to_center(df):\n    x_mean = df[\"x_center\"].mean()\n    y_mean = df[\"y_center\"].mean()\n    distances = []\n    for _, row in df.iterrows():\n        distances.append(distance.euclidean((row[\"x_center\"], row[\"y_center\"]), (x_mean, y_mean)))\n    df[\"distances\"] = distances\n    df[\"std_dev_distances\"] = np.abs(stats.zscore(df[\"distances\"]))\n    df[\"Probability Grouping\"] = df[\"std_dev_distances\"].apply(lambda x: get_probability_descriptor(x))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.0 Aortic Enlargement\n\nBelow we're plotting the density heatmap on the left with the bounding box centers on the right. The Probability Grouping describes how many standard deviations from the mean each bounding box center occurs. As we can see, outliers colored in both green and red represent bounding box centers that are extremely rare events. These bounding box centers (particularly the red ones) will be difficult to generalize in any model, as they exist very far away from the majority of our training examples."},{"metadata":{"trusted":true},"cell_type":"code","source":"n_clusters = 1\nclass_id = 0\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == class_id)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\n\nkmeans = KMeans(n_clusters=n_clusters)\nkmeans.fit(coi[[\"x_center\", \"y_center\"]])\ncentroids = kmeans.predict(coi[[\"x_center\", \"y_center\"]])\ncenters = kmeans.cluster_centers_\n\nfig, axs = plt.subplots(nrows=1, ncols=3, figsize=(15, 10), sharey=True)\n\naxs[0].imshow(get_bbox(train_aug, 0), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=centroids, \n    ax=axs[2],\n)\nsns.scatterplot(\n    centers[:, 0], centers[:, 1], c=['black'], s=200, alpha=0.5, ax=axs[2],\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.1 Atelectasis\n\nHere we should consider both the left and right sides separately. This is because we may see variations between sides."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 1) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 1), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 1) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.2 Calcification\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 2) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 2), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 2) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.3 Cardiomegaly\n\nSimilar to Aortic Enlargment, there are certain bounding box centers that exist far outside where we would expect to see them. In this case, groupings that exist in red represent less than a 1% probability of occurring (based on the distribution of bounding box centers in our test data). Again, red and green points are likely going to impact our ability to fit a model accordingly. "},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=2, figsize=(15, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 3), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.4 Consolidation\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 4) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 4), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 4) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.5 ILD\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 5) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 5), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 5) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.6 Infiltration\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 6) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 6), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 6) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.7 Lung Opacity\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 7) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 7), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 7) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.8 Nodule / Mass\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 8) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 8), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 8) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.9 Other Lesion\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=2, figsize=(15, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 9)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 9), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.10 Pleural Effusion\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 10) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 10), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 10) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.11 Pleural Thickening\n\nAgain, we have both the right and left side of the body to consider as separate entities. We'll also look to see if we can group them into specific segments using KMeans to search for centroids."},{"metadata":{"trusted":true},"cell_type":"code","source":"coi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 11)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\n\nfig, axs = plt.subplots(nrows=1, ncols=4, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 11) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 11), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 11) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 11)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nkmeans = KMeans(n_clusters=4)\nkmeans.fit(coi[[\"x_center\", \"y_center\"]])\ncentroids = kmeans.predict(coi[[\"x_center\", \"y_center\"]])\ncenters = kmeans.cluster_centers_\n\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=centroids, \n    ax=axs[3],\n)\nsns.scatterplot(\n    centers[:, 0], \n    centers[:, 1], \n    c=['black'], \n    s=200, \n    alpha=0.5, \n    ax=axs[3],\n)\n_ = axs[3].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[3].set(xlim=(0, 400), ylim=(500, 0))\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.12 Pnemothorax\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=3, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 12) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 12), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 12) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 3.3.13 Pulmonary Fibrosis\n\nAgain, we have both the right and left side of the body to consider as separate entities."},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axs = plt.subplots(nrows=1, ncols=4, figsize=(20, 10), sharey=True)\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 13) & (train_aug[\"x_min_norm\"] < 0.3)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\naxs[0].imshow(get_bbox(train_aug, 13), cmap='hot', interpolation='nearest', extent=[0,400,500,0], aspect=\"auto\")\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[1]\n)\n_ = axs[1].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[1].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 13) & (train_aug[\"x_min_norm\"] > 0.5)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=coi[\"Probability Grouping\"], \n    hue_order=[\"95% probability\", \"5% probability\", \"1% probability\", \"< 1% probability\"], \n    ax=axs[2]\n)\n_ = axs[2].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[2].set(xlim=(0, 400), ylim=(500, 0))\n\ncoi = pd.DataFrame(train_aug[(train_aug[\"class_id\"] == 13)])\ncalculate_centers(coi)\ncalculate_distances_to_center(coi)\nkmeans = KMeans(n_clusters=4)\nkmeans.fit(coi[[\"x_center\", \"y_center\"]])\ncentroids = kmeans.predict(coi[[\"x_center\", \"y_center\"]])\ncenters = kmeans.cluster_centers_\n\nsns.scatterplot(\n    x=\"x_center\", \n    y=\"y_center\", \n    data=coi[[\"x_center\", \"y_center\"]], \n    hue=centroids, \n    ax=axs[3],\n)\nsns.scatterplot(\n    centers[:, 0], \n    centers[:, 1], \n    c=['black'], \n    s=200, \n    alpha=0.5, \n    ax=axs[3],\n)\n_ = axs[3].set(xlabel=\"X position\", ylabel=\"Y position\")\n_ = axs[3].set(xlim=(0, 400), ylim=(500, 0))\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 3.4 Discussion\n\nOur heatmaps provide us with some interesting results:\n\n1. The heatmaps can be used to provide a confidence rating for any bounding box that we place. If the box is in a hot zone for a particular clinical finding, we can rate it with a higher degree of confidence. In essence, we would use the heatmap as a confidence prior.\n2. The heatmaps are nicely showing us where various clinical findings tend to be localized. We can see for class `0 - Aortic Enlargement` that we expect to see our bounding boxes to appear near the center of the x-ray towards the top of the image. This can help us out greatly. Instead of building a single monolithic classifier, we can build a highly targetted classifier that examines that single area. In essence, we can _zoom in_ on an area as indicated by the heatmap, and focus learning techniques on that specific area instead of the entire image.\n3. Localizations tell us that there are generally more findings in the right lung for issues such as `1 - Atelectasis`, `4 - Consolidation`, `7 - Lung Opacity` and `12 - Pneumothorax`. This may be significant possibly due to the differences in lung sizes - the right lung is larger than the left.\n4. Based on gender, `12 - Pneumothorax` appears to occur equally on the right and left lung for females. Similarly, `2 - Calcification` appears to occur more frequently on the left lung for females when compared to males. \n5. Certain heatmaps can give us clues how to build additional classifiers. For example, for heatmaps that have symmetrical localizations such as `5 - ILD`, we can very likely build a single classifier that recognizes a single side of the chest (say the right side). Then, we can simply flip the image and see if there is a finding using the same classifier. In this way, we can build single, highly targetted classifiers, but cut down training time (and increase training samples) by using both the left and right hand sides of the images as input sources."},{"metadata":{},"cell_type":"markdown","source":"# 4 Aspect Ratios\n\nAnother important aspect to consider for the bounding boxes and their localization, is the approximate aspect ratio that each bounding box may take. While we were able to see for some classes already that the aspect ratio reflects only a handful of options (generally 1:1 or 2:1), it would be a good idea to get an idea of how varied aspect ratios actually are on a per-class basis.\n\nThe code below will generate an approximate aspect ratio for each of our bounding boxes, and add it to the dataframe."},{"metadata":{"trusted":true},"cell_type":"code","source":"import math\n\ndef get_aspect_ratio(width, height):\n    gcd = math.gcd(width, height)\n    w_ratio = int(width / gcd)\n    h_ratio = int(height / gcd)\n    return w_ratio, h_ratio\n\ndef get_approximate_aspect_ratio(width, height):\n    w_ratio, h_ratio = get_aspect_ratio(width, height)\n    while w_ratio > 10:\n        w_ratio = round(w_ratio / 10)\n        h_ratio = round(h_ratio / 10)\n    return get_aspect_ratio(w_ratio, h_ratio)\n    \ndef calculate_aspect_ratio(df):\n    ar_width = []\n    ar_height = []\n    for _, row in df.iterrows():\n        if row[\"class_id\"] != 14:\n            height = int(row[\"y_max\"] - row[\"y_min\"])\n            width = int(row[\"x_max\"] - row[\"x_min\"])\n            width_ratio, height_ratio = get_approximate_aspect_ratio(width, height)\n            ar_width.append(width_ratio)\n            ar_height.append(height_ratio)\n        else:\n            ar_width.append(1)\n            ar_height.append(1)\n    df[\"width_aspect_ratio\"] = ar_width\n    df[\"height_aspect_ratio\"] = ar_height\n    \ncalculate_aspect_ratio(train_aug)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"With the dataframe updated to contain aspect information, we can now add that and plot via a heatmap. The resulting heatmaps will be a plot of _width_ by _height_ that denotes the specific aspect ratio. For example, a width of 4 and a height of 3 indicates a _4:3_ aspect ratio. Intensity will indicate how many findings of that class have a bounding box with that aspect ratio."},{"metadata":{"trusted":true},"cell_type":"code","source":"import seaborn as sns\n\ndef generate_heatmap_per_class(class_id, df):\n    dataframe = df[(df[\"class_id\"] == class_id)]\n    heatmap_array = np.zeros(shape=(dataframe[\"height_aspect_ratio\"].max() + 1, dataframe[\"width_aspect_ratio\"].max() + 1))\n    for _, row in dataframe.iterrows():\n        heatmap_array[row[\"height_aspect_ratio\"], row[\"width_aspect_ratio\"]] += 1\n    return heatmap_array\n\nnames = [name for name in classes]\nheatmaps = [generate_heatmap_per_class(index, train_aug) for index, _ in enumerate(classes)]\n\nfig, axs = plt.subplots(nrows=5, ncols=3, figsize=(20, 30))\nfig.suptitle(\"Approximate Bounding Box Aspect Ratios\", fontweight=\"bold\", size=15)\nfor i, (name, heatmap, ax) in enumerate(zip(names, heatmaps, axs.flatten())):\n    hm = sns.heatmap(heatmap, annot=False, ax=ax)\n    _ = ax.set_xlabel(\"Width\")\n    _ = ax.set_ylabel(\"Height\")\n    _ = ax.set_title(name, size=15)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# 4.1 Discussion\n\nFrom our aspect ratios, we can start to see some simple trends:\n\n1. Most classes have a strong affinity towards a _1:1_ aspect ratio. For example `0 - Aortic Enlargement` and `8 - Nodule / Mass` contain bounding boxes with an aspect ratio that are primarily _1:1_ with few other variants occurring.\n2. `3 - Cardiomegaly` tend to have bounding boxes that are wider than they are tall.\n3. `12 - Pneumothorax` have a wide variety of aspect ratios.\n4. `9 - Other lesion` can have any number of conceiveable aspect ratios.\n\nAgain, we can use these aspect ratio heatmaps to help us establish a confidence prior to our bounding box regressor. In general, we can look to see if bounding boxes we generate are within an established norm, and help us determine if they need refinement."},{"metadata":{},"cell_type":"markdown","source":"# 5. Conclusions\n\nBy plotting heatmaps of our bounding boxes, we can easily localize where various findings are expected. We can use heatmap information to establish confidence priors to help with our predictions. Aspect ratio information across classes can help us determine if our predicted bounding boxes make sense, and again can help establish confidence priors.\n\nIf you find this visualization useful, please comment or upvote the notebook! Thanks!"}],"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":4,"nbformat_minor":4}