{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Intro\n## Problem being solved\n  * when displaying pixel data, CT scans don't look like one would expect."},{"metadata":{"_kg_hide-input":true,"trusted":true,"_kg_hide-output":false},"cell_type":"code","source":"from IPython.display import Image\nImage('../input/illustrations/preproc.jpg')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Plan\n1. [Libraries](#Libraries)\n2. [Data](#Data)\n3. [Raw images visualization](#Raw-images-visualization)\n  * [Comparison with expected CT images](#Comparison-with-expected-CT-images)\n4. [Investigating lack of brain tissue details](#Investigating-lack-of-brain-tissue-details)\n  * [Intensities range](#Intensities-range)\n    * [Challenge with fine-grained intensities](#Challenge-with-fine-grained-intensities)\n    * [Manual range reduction](#Manual-range-reduction)\n5. [Dataset diversity](#Dataset-diversity)\n6. [Clustering images](#Clustering-images)\n  * [Generating intensity histograms](#Generating-intensity-histograms)\n  * [HDBSCAN](#HDBSCAN)\n  * [Visual check of the clusters](#Visual-check-of-the-clusters)\n7. [Range selection for each cluster](#Range-selection-for-each-cluster)\n  * [Cluster-0 image type](#Cluster-0-image-type)\n  * [Cluster-1 image type](#Cluster-1-image-type)\n  * [Cluster-2 image type](#Cluster-2-image-type)\n  * [Cluster-3 image type](#Cluster-3-image-type)\n  * [Cluster-4 image type](#Cluster-4-image-type)\n8. [Cluster-type Classifier](#Cluster-type-Classifier)\n  * [Data preprocessing](#Data-preprocessing)\n  * [Stratified train/test split](#Stratified-train/test-split)\n  * [SVM](#SVM)\n  * [Performance check](#Performance-check)\n9. [Check: preprocessing of new images](#Check:-preprocessing-of-new-images)\n10. [Distribution of hemorrhages per image-type](#Distribution-of-hemorrhages-per-image-type)\n11. [Next steps](#Next-steps)"},{"metadata":{},"cell_type":"markdown","source":"P.S. Metadata in the DICOM files cannot be used per competition rules. That's why I rushed into the detective work you can find below :)\n\nOtherwise, with metadata it'd be possible to extract the most relevant info (so that the scans would look like ones in the right picture). More details can be found [here](https://www.kaggle.com/c/rsna-intracranial-hemorrhage-detection/discussion/109328). It must be said that in the end metadata might be allowed, I'd recommend to check [the *Welcome!* thread](https://www.kaggle.com/c/rsna-intracranial-hemorrhage-detection/discussion/109258).\n\nP.P.S. In theory deep nets are able to learn preprocessing, however in practice lack of preprocessing unnecessarily complicates training\n\nP.P.P.S. :D \nI've checked 100.000 randomly sampled training images, and `RescaleSlope` was always equal to 1. I believe, it's reasonable to assume that pixel data are in Hounsfield units, and there's no need of intensity rescaling. "},{"metadata":{},"cell_type":"markdown","source":"# Libraries"},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"!pip install hdbscan","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"import os\nimport numpy as np\nimport pydicom\nimport pandas as pd\nfrom random import sample\nfrom tqdm import tqdm_notebook as tqdm\nimport hdbscan\nfrom scipy.spatial.distance import jensenshannon\nfrom sklearn.svm import LinearSVC\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.metrics import f1_score\nimport pickle\nimport random\nrandom.seed(1)\nfrom numpy.random import seed\nseed(1)\nimport matplotlib.pyplot as plt\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"train_folder = '../input/rsna-intracranial-hemorrhage-detection/stage_1_train_images/'","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Data"},{"metadata":{},"cell_type":"markdown","source":"## Organizing labels"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_df = pd.read_csv('../input/rsna-intracranial-hemorrhage-detection/stage_1_train.csv')\ntrain_df.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"img_labels_df = train_df[train_df['Label'] == 1].copy()\nimg_labels_df['Label_String'] = img_labels_df['ID'].map(lambda x: x.split('_')[-1])\nimg_labels_df['Image'] = img_labels_df['ID'].map(lambda x: 'ID_' + x.split('_')[1])\nimg_2_labels = img_labels_df.groupby('Image')['Label_String'].agg(list)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plt.hist(img_2_labels.map(len))\nplt.xlabel('Number of conditions per positive image')\nplt.ylabel('Count')\nplt.title('Counts of labels per positive image')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"No positive image has a single label. Let's check if the label ```any``` is redundant."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"print(f\"Label any is present in {100*sum(img_2_labels.map(lambda x: 'any' in x))/len(img_2_labels):.2f}% of positive images.\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's remove the redundant label."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"img_2_labels = img_2_labels.map(lambda x: [i for i in x if i != 'any'])\nimg_2_labels = img_2_labels.map(lambda x: ', '.join(x))\nimg_2_labels.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Images for analysis"},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_img_names = os.listdir('../input/rsna-intracranial-hemorrhage-detection/stage_1_train_images/')\nprint(f'There are {len(train_img_names)} train images.')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_img_names_subsample = sample(train_img_names, 1000)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Raw images visualization"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def get_img(img_name, folder=train_folder):\n    pydicom_filedataset = pydicom.read_file(os.path.join(folder, img_name))\n    return pydicom_filedataset.pixel_array\n\ndef vizualize_tuple(imgs, img_names, grid_size=3):\n    fig, axes = plt.subplots(grid_size, grid_size, figsize=(grid_size+4, grid_size+4))\n    fig.tight_layout()\n    for img_i, img in enumerate(imgs):\n        ax = axes[img_i//grid_size, img_i%grid_size]\n        ax.imshow(img, cmap='bone')\n        ax.axis('off')\n        ax.set_title(f\"{img_names[img_i]}:\\n{img_2_labels.get(img_names[img_i].replace('.dcm', ''), '')}\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Random subsample"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"img_names_subsample = sample(train_img_names_subsample, 9)\nvizualize_tuple([get_img(img_name) for img_name in img_names_subsample], img_names_subsample)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Comparison with expected CT images"},{"metadata":{},"cell_type":"markdown","source":"We can see the visualizations we get don't look like ones on the explanation page:\n![](https://cdn.discordapp.com/attachments/507208726864855060/623862702955036702/RSNA_Intracranial_Hemorrhage_Detection_Kaggle.png)"},{"metadata":{},"cell_type":"markdown","source":"## Summary on raw images"},{"metadata":{},"cell_type":"markdown","source":"To sum up,\nraw images don't have detalization of example images (e.g. in the figure with the explanations). Due to the lack of details, pictures of different hemorrhages then look alike. Even though the required info likely is somewhere in the image, pulling it to the surface before feeding to the network would be desirable. "},{"metadata":{},"cell_type":"markdown","source":"# Investigating lack of brain tissue details"},{"metadata":{},"cell_type":"markdown","source":"I'll focus on a pair of images which are look-alike to me, even though one has ```Intraventricular``` and another one has ```Intraparenchymal``` in it."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"imgs_to_compare = ['ID_df5ae8f49.dcm', 'ID_8cb74b318.dcm']\n\ndef plot_img_pair(imgs, suptitle='', img_names=imgs_to_compare):\n    fig, axes = plt.subplots(1, 2, figsize=(7, 6))\n    fig.tight_layout()\n    for img_i, img in enumerate(imgs):\n        ax = axes[img_i]\n        ax.imshow(img, cmap='bone')\n        ax.axis('off')\n        ax.set_title(f\"{img_names[img_i]}:\\n{img_2_labels.get(img_names[img_i].replace('.dcm', ''), '')}\")\n    plt.suptitle(suptitle, fontsize=15)\n        \nplot_img_pair([get_img(img_name) for img_name in imgs_to_compare])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Intensities range"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def plot_intensity_hists(imgs, bins=None, img_names=imgs_to_compare):\n    min_intesity = min([np.min(img) for img in imgs])\n    max_intesity = max([np.max(img) for img in imgs])\n    for img_i, img in enumerate(imgs):\n        plt.figure(figsize=(4, 4))\n        plt.hist(img.flatten(), bins=bins)\n        plt.xlabel('Image intensities', fontsize=11)\n        plt.ylabel('Pixels count', fontsize=11)\n        plt.xlim((min_intesity, max_intesity))\n        plt.title(f\"{img_names[img_i]}:\\n{img_2_labels.get(img_names[img_i].replace('.dcm', ''), '')}\", fontsize=14)\n        \nplot_intensity_hists([get_img(img_name) for img_name in imgs_to_compare])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Challenge with fine-grained intensities"},{"metadata":{},"cell_type":"markdown","source":"We can see, the ranges of intensities are very wide: 5000 instead of standard 255. As a result, approximately 20 different values of raw intensities become effectively indistinguishable from each other to a human eye. And for a deep net the relevant information is hidden in a pile of less relevant data.\n\nLet's manually select intensities range with the most important information for the competition task."},{"metadata":{},"cell_type":"markdown","source":"### Manual range reduction"},{"metadata":{},"cell_type":"markdown","source":"Let's first clip away negative intensities."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plot_img_pair([get_img(img_name) for img_name in imgs_to_compare], 'Initial Images')\n\ndef clip_negatives(img, threshold=0):\n    img[img < threshold] = threshold\n    return img\n\nplot_img_pair([clip_negatives(get_img(img_name)) for img_name in imgs_to_compare], 'Without Negative Intensities')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We got rid of the redundant background and we can already start seeing some brain-tissue details. Let's see new intensities ranges."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plot_intensity_hists([clip_negatives(get_img(img_name)) for img_name in imgs_to_compare], bins=50)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's leave values between 700 and 1500."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plot_img_pair([get_img(img_name) for img_name in imgs_to_compare], 'Initial Images')\ndef clip_positives(img, threshold=1500):\n    img[img > threshold] = threshold\n    return img\n\nlb = 700\nub = 1500\nclipped_imgs = [clip_positives(clip_negatives(get_img(img_name), lb), ub) for img_name in imgs_to_compare]\nplot_img_pair(clipped_imgs, f'Intensities between {lb} and {ub}')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It's exciting to start seeing brain tissue."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plot_intensity_hists(clipped_imgs, bins=50)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's tighten the intensities range."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"plot_img_pair([get_img(img_name) for img_name in imgs_to_compare], 'Initial Images')\n\nlb = 900\nub = 1200\nclipped_imgs = [clip_positives(clip_negatives(get_img(img_name), lb), ub) for img_name in imgs_to_compare]\nplot_img_pair(clipped_imgs, f'Intensities between {lb} and {ub}')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Dataset diversity"},{"metadata":{},"cell_type":"markdown","source":"Based on the [initial sample of images](#Random-subsample), we saw that there are quite diverse images in the dataset:"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"img_names_diverse = ['ID_c96b7ba2a.dcm', 'ID_8cb74b318.dcm']\nimgs_diverse = [get_img(img_name) for img_name in img_names_diverse]\nplot_img_pair(imgs_diverse, img_names=img_names_diverse)\nplot_intensity_hists(imgs_diverse, img_names=img_names_diverse)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Hence, naive intensity range selection is not an option."},{"metadata":{},"cell_type":"markdown","source":"# Clustering images"},{"metadata":{},"cell_type":"markdown","source":"Let's cluster images based on histograms of intensities. We'll use [HDBSCAN](https://hdbscan.readthedocs.io) for clustering."},{"metadata":{},"cell_type":"markdown","source":"## Generating intensity histograms"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def generate_hist_vec(img_name, bins=list(range(-50, 1500, 25))):\n    img = get_img(img_name)\n    return np.histogram(img[(img > -50) & (img < 1500)].flatten(), bins=bins)[0]\n\nhistogram_vectors = []\nfor img_name in tqdm(train_img_names_subsample, desc='Generation histogram vectors..'):\n    histogram_vectors.append(generate_hist_vec(img_name))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"histograms_df = pd.DataFrame(np.array(histogram_vectors), columns=[f'bin_{i}' for i in range(len(histogram_vectors[0]))])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## HDBSCAN"},{"metadata":{},"cell_type":"markdown","source":"As we are clustering histograms, we'd like to use metric capturing similarity of two distributions. Let's use [Jensen–Shannon divergence](https://en.wikipedia.org/wiki/Jensen%E2%80%93Shannon_divergence). One of HDBSCAN advantages is its ability to work with arbitrary similarity measures."},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"# clusterer = DBSCAN(eps=0.1, min_samples=3, metric=jensenshannon)\nclusterer = hdbscan.HDBSCAN(metric=jensenshannon)\nclusterer.fit(histograms_df)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"print(f'There are {max(clusterer.labels_) + 1} clusters.')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"print(f\"\"\"Sizes of the clusters:\n{pd.Series(clusterer.labels_).value_counts()}.\"\"\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Visual check of the clusters"},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_cluster_sample(cluster_i, sample_size=9):\n    cluster_img_names = [img for img, is_in in zip(train_img_names_subsample, clusterer.labels_==cluster_i) \n                               if is_in]\n    cluster_sample = sample(cluster_img_names, min(sample_size, len(cluster_img_names)))\n    return cluster_sample\n\n                \ndef check_cluster_vizually(cluster_sample, lb=-4000, ub=5000, hist=True): \n    clipped_imgs = [clip_positives(clip_negatives(get_img(img_name), lb), ub) for img_name in cluster_sample]\n                    \n    def vizualize_nine_hists(imgs):\n        min_intesity = min([np.min(img) for img in imgs])\n        max_intesity = max([np.max(img) for img in imgs])\n        fig, axes = plt.subplots(3, 3, figsize=(7, 7))\n        for img_i, img in enumerate(imgs):\n            ax = axes[img_i//3, img_i%3]\n            ax.hist(img.flatten(), bins=20)\n            ax.set_xlim((min_intesity, max_intesity))\n            ax.set_yticks([], [])\n        fig.suptitle('Histograms of pixel intensities')\n        \n    vizualize_tuple(clipped_imgs, cluster_sample)\n    if hist:\n        vizualize_nine_hists(clipped_imgs)\n    \n    \ncluster_samples = [get_cluster_sample(cluster_i) for cluster_i in range(max(clusterer.labels_) + 1)]\nfor cluster_i, cluster_sample in enumerate(cluster_samples):\n    plt.figure(figsize=(9, 1))\n    plt.plot(np.arange(20), np.ones(20))\n    plt.title(f'Cluster {cluster_i}')\n    check_cluster_vizually(cluster_sample, hist=False)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Range selection for each cluster\n\nLet's visualization of pixel data together with pixel intensity histograms.\nApproach: gradually clip values in a search for main class mode of the histograms, then try to center the mode (between range bounds). Then check pixels around boarders (which are clipped)."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"cluster_2_intensity_limits = []","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Cluster-0 image type"},{"metadata":{},"cell_type":"markdown","source":"### With chosen range clipping"},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[0], lb=930, ub=1150)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Checking that nothing very informative was thrown away\nThe lost information would be gray/brain-tissue-like (between zero and one)"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[0], lb=830, ub=930, hist=False)\ncheck_cluster_vizually(cluster_samples[0], lb=1150, ub=1250, hist=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"cluster_2_intensity_limits.append((930, 1150))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Cluster-1 image type"},{"metadata":{},"cell_type":"markdown","source":"### With range clipping"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[1], lb=920, ub=1200)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Checking that nothing very informative was thrown away\nThe lost information would be gray/brain-tissue-like (between zero and one)"},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[1], lb=820, ub=920, hist=False)\ncheck_cluster_vizually(cluster_samples[1], lb=1200, ub=1300, hist=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"cluster_2_intensity_limits.append((920, 1200))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Cluster-2 image type"},{"metadata":{},"cell_type":"markdown","source":"### With range clipping"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[2], lb=5, ub=300)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Checking that nothing very informative was thrown away\nThe lost information would be gray, brain-tissue-like (between zero and one)"},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[2], lb=-95, ub=5, hist=False)\ncheck_cluster_vizually(cluster_samples[2], lb=300, ub=400, hist=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"cluster_2_intensity_limits.append((5, 300))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Cluster-3 image type"},{"metadata":{},"cell_type":"markdown","source":"### Bones? "},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[3], lb=1500, ub=2000)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We observe higher HU values.\n\n> Dense materials such as bone have density values approaching +1000 HU.\n([source](https://www.sciencedirect.com/topics/medicine-and-dentistry/hounsfield-scale))\n\nTherefore, I'd consider this cluster to correspond to low brain visibility and strong skull pattern. However, the skull itself is of low interest (hopefully I'm not mistaken). So I'd clip its intensities."},{"metadata":{},"cell_type":"markdown","source":"### With range clipping"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[3], lb=970, ub=1300)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Checking that nothing very informative was thrown away\nThe lost information would be gray, brain-tissue-like (between zero and one)"},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[3], lb=870, ub=970, hist=False)\ncheck_cluster_vizually(cluster_samples[3], lb=1300, ub=1400, hist=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"cluster_2_intensity_limits.append((970, 1300))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Cluster-4 image type"},{"metadata":{},"cell_type":"markdown","source":"### Again Skull?"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[4], lb=1400, ub=2100)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### With range clipping"},{"metadata":{"trusted":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[4], lb=900, ub=1300)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Checking that nothing very informative was thrown away\nThe lost information would be gray, brain-tissue-like (between zero and one)"},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"check_cluster_vizually(cluster_samples[4], lb=800, ub=900, hist=False)\ncheck_cluster_vizually(cluster_samples[4], lb=1300, ub=1400, hist=False)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It seems this cluster is very similar to the previous one. Perhaps, a tiny a \"halo\" around the skull might be the main difference."},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"cluster_2_intensity_limits.append((900, 1300))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Cluster-type Classifier"},{"metadata":{},"cell_type":"markdown","source":"## Data preprocessing\n### Removing the noise cluster"},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"non_noise_bool_index = clusterer.labels_ != -1\nhistograms_df = histograms_df[non_noise_bool_index]\nlabels = clusterer.labels_[non_noise_bool_index]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Scaling data"},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"scaler = StandardScaler()\nhistograms_df_scaled = scaler.fit_transform(histograms_df)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Stratified train/test split"},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df, test_df, train_labels, test_labels = train_test_split(histograms_df_scaled, labels, test_size=0.2, stratify=labels)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## SVM\nLinear SVM for it to be fast."},{"metadata":{"trusted":true},"cell_type":"code","source":"svm = LinearSVC(multi_class='crammer_singer')\nsvm.fit(train_df, train_labels)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Performance check"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"test_predictions = svm.predict(test_df)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"print(f\"Average weighted F1 score on test data is: {f1_score(test_labels, test_predictions, average='weighted')}.\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"with open('scaler.pkl', 'wb') as f:\n    pickle.dump(scaler, f)\n    \nwith open('svm_cluster_type.pkl', 'wb') as f:\n    pickle.dump(svm, f)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Check: preprocessing of new images"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"used_set = set(train_img_names_subsample)\nunseen_img_names = [img for img in train_img_names if img not in used_set]\nunseen_img_sample = sample(unseen_img_names, 9)\n\nhistogram_vectors = [generate_hist_vec(img) for img in unseen_img_sample]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"used_set = set(train_img_names_subsample)\nunseen_img_names = [img for img in train_img_names if img not in used_set]\nunseen_img_sample = sample(unseen_img_names, 16)\n\nhistogram_vectors = scaler.transform([generate_hist_vec(img) for img in unseen_img_sample])\nimg_cluster_classes = svm.predict(histogram_vectors)\nclipping_limits = [cluster_2_intensity_limits[class_i] for class_i in img_cluster_classes]\nclipped_imgs = [clip_positives(clip_negatives(get_img(img_name), lb), ub) for img_name, (lb, ub) in zip(unseen_img_sample,\n                                                                                                       clipping_limits)]\nvizualize_tuple(clipped_imgs, unseen_img_sample, 4)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Distribution of hemorrhages per image-type"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def plot_pie_per_cluster_type(cluster_class_i):\n    all_known_cluster_imgs = get_cluster_sample(cluster_class_i, sample_size=float('inf'))\n    all_deceaseas = []\n    for img in all_known_cluster_imgs:\n        all_deceaseas.extend(img_2_labels.get(img.replace('.dcm', ''), 'nothing').split(', '))\n    all_deceaseas_counts = pd.Series(all_deceaseas).value_counts()\n    \n    fig, ax = plt.subplots(figsize=(7, 7))\n    # credits: https://stackoverflow.com/questions/6170246/how-do-i-use-matplotlib-autopct\n    def make_autopct(values):\n        def my_autopct(pct):\n            total = sum(values)\n            val = int(round(pct*total/100.0))\n            return '{p:.0f}%  ({v:d})'.format(p=pct,v=val)\n        return my_autopct\n    ax.pie(all_deceaseas_counts, labels=all_deceaseas_counts.index, autopct=make_autopct(all_deceaseas_counts), shadow=True, startangle=90)\n    ax.axis('equal')\n    plt.title(f'Distribution of deceases per class {cluster_class_i}', fontsize=15)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"for cluster_class_i in range(max(clusterer.labels_)+1):\n    plot_pie_per_cluster_type(cluster_class_i)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"To sum up, based on the allowed-to-use intensities we've arrived at cluster classes which in addition to the rules-complient preprocessing provide us with some info about final labels."},{"metadata":{},"cell_type":"markdown","source":"# Next steps\n  1. Generation of the pre-processed dataset,\n  2. Comparison of deep net training convergence with and without the preprocessing,\n  3. Improvements of the pre-processing if needed."}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.6.7"}},"nbformat":4,"nbformat_minor":1}