{"cells":[{"metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5"},"cell_type":"markdown","source":"# Intro\nIn this notebook I'd like to check what we can get directly from a classifier. It can be considered the most simple baseline for weakly supervised segmentation based on image-level labels.\n\nI'd use the classifier trained in [this notebook](https://www.kaggle.com/samusram/cloud-classifier-for-post-processing). I'd use explainability technique [Gradient-weighted Class Activation Mapping](http://gradcam.cloudcv.org/) to extract masks from the classifier."},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","collapsed":true},"cell_type":"markdown","source":"# Plan\n1. [Libraries](#Libraries)\n2. [Preparing data](#Preparing-data)\n  * [One-hot encoding classes](#One-hot-encoding-classes)\n  * [Stratified split into train/val](#Stratified-split-into-train/val)\n  * [Data preprocessing](#Data-preprocessing)\n      * [Validation images](#Validation-images)\n      * [Masks](#Masks)\n3. [Grad-CAM routines](#Grad-CAM-routines)\n4. [Precision to threshold mapping per class](#Precision-to-threshold-mapping-per-class)\n5. [Vizually comparing GradCAM's with masks.](#Vizually-comparing-GradCAM's-with-masks)\n6. [Estimating performance](#Estimating-performance)\n7. [Predicting test masks](#Predicting-test-masks)\n8. [Conclusion](#Conclusion)\n9. [Credits](#Credits)"},{"metadata":{},"cell_type":"markdown","source":"# Libraries"},{"metadata":{"_kg_hide-output":true,"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import os, glob\nimport random\nfrom sklearn.model_selection import train_test_split\nimport cv2\nimport numpy as np\nimport pandas as pd\nimport multiprocessing\nfrom copy import deepcopy\nfrom sklearn.metrics import precision_recall_curve, auc\nimport tensorflow.keras as keras\nfrom tensorflow.keras.models import load_model\nfrom tensorflow.keras.utils import Sequence\nfrom tensorflow.keras import backend as K\nimport matplotlib.pyplot as plt\nfrom IPython.display import Image\nfrom tqdm import tqdm_notebook as tqdm\nfrom tensorflow.python.framework import ops\nfrom numpy.random import seed\nseed(10)\nfrom tensorflow import set_random_seed\nimport tensorflow as tf\nset_random_seed(10)\nimport gc\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"test_imgs_folder = '../input/understanding_cloud_organization/test_images/'\ntrain_imgs_folder = '../input/understanding_cloud_organization/train_images/'\nimage_width, image_height = 224, 224\nnum_cores = multiprocessing.cpu_count()\nbatch_size = 64\nmodel_class_names =  ['Fish', 'Flower', 'Sugar', 'Gravel']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Preparing data"},{"metadata":{},"cell_type":"markdown","source":"## One-hot encoding classes"},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"train_df_orig = pd.read_csv('../input/understanding_cloud_organization/train.csv')\ntrain_df = pd.read_csv('../input/understanding_cloud_organization/train.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"train_df = train_df[~train_df['EncodedPixels'].isnull()]\ntrain_df['Image'] = train_df['Image_Label'].map(lambda x: x.split('_')[0])\ntrain_df['Class'] = train_df['Image_Label'].map(lambda x: x.split('_')[1])\nclasses = train_df['Class'].unique()\ntrain_df = train_df.groupby('Image')['Class'].agg(set).reset_index()\nfor class_name in classes:\n    train_df[class_name] = train_df['Class'].map(lambda x: 1 if class_name in x else 0)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"# dictionary for fast access to ohe vectors\nimg_2_ohe_vector = {img:vec for img, vec in zip(train_df['Image'], train_df.iloc[:, 2:].values)}","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Stratified split into train/val"},{"metadata":{},"cell_type":"markdown","source":"The same split was used to train the classifier."},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"train_imgs, val_imgs = train_test_split(train_df['Image'].values, \n                                        test_size=0.2, \n                                        stratify=train_df['Class'].map(lambda x: str(sorted(list(x)))), # sorting present classes in lexicographical order, just to be sure\n                                        random_state=10)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Data preprocessing"},{"metadata":{},"cell_type":"markdown","source":"### Validation images"},{"metadata":{},"cell_type":"markdown","source":"Reading val images, resizing and storing in memory to speed up experiments."},{"metadata":{"trusted":true},"cell_type":"code","source":"val_imgs_np = np.empty((len(val_imgs), image_height, image_width, 3))\nfor img_i, img_name in enumerate(tqdm(val_imgs)):\n    img_path = os.path.join(train_imgs_folder, img_name)\n    val_imgs_np[img_i, :, :, :] = cv2.resize(cv2.imread(img_path), (image_height, image_width)).astype(np.float32)/255.0","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Masks"},{"metadata":{},"cell_type":"markdown","source":"Analogously, storing in memory validation masks."},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"# helper functions\n# credits: https://www.kaggle.com/artgor/segmentation-in-pytorch-using-convenient-tools?scriptVersionId=20202006\ndef rle_decode(mask_rle: str = '', shape: tuple = (1400, 2100)):\n    '''\n    Decode rle encoded mask.\n    \n    :param mask_rle: run-length as string formatted (start length)\n    :param shape: (height, width) of array to return \n    Returns numpy array, 1 - mask, 0 - background\n    '''\n    s = mask_rle.split()\n    starts, lengths = [np.asarray(x, dtype=int) for x in (s[0:][::2], s[1:][::2])]\n    starts -= 1\n    ends = starts + lengths\n    img = np.zeros(shape[0] * shape[1], dtype=np.uint8)\n    for lo, hi in zip(starts, ends):\n        img[lo:hi] = 1\n    return img.reshape(shape, order='F')\n\ndef mask2rle(img):\n    '''\n    Convert mask to rle.\n    img: numpy array, 1 - mask, 0 - background\n    Returns run length as string formated\n    '''\n    pixels= img.T.flatten()\n    pixels = np.concatenate([[0], pixels, [0]])\n    runs = np.where(pixels[1:] != pixels[:-1])[0] + 1\n    runs[1::2] -= runs[::2]\n    return ' '.join(str(x) for x in runs)\n\n\ndef make_mask(df, image_label, shape: tuple = (1400, 2100)):\n    \"\"\"\n    Create mask based on df, image name and shape.\n    \"\"\"\n    df = df.set_index('Image_Label')\n    encoded_mask = df.loc[image_label, 'EncodedPixels']\n    mask = np.zeros((shape[0], shape[1]), dtype=np.float32)\n    if encoded_mask is not np.nan:\n        mask = rle_decode(encoded_mask)\n            \n    return cv2.resize(mask, (image_height, image_width))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"val_masks_np = np.empty((len(model_class_names), len(val_imgs), image_height, image_width))\nfor class_i, class_name in enumerate(tqdm(model_class_names)):\n    for img_i, img_name in enumerate(val_imgs):\n        mask = make_mask(train_df_orig, img_name + '_' + class_name)\n        val_masks_np[class_i][img_i] = mask","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Grad-CAM routines"},{"metadata":{"trusted":true,"_kg_hide-output":true,"_kg_hide-input":true},"cell_type":"code","source":"model = load_model('../input/clouds-classifier-files/classifier_densenet169_epoch_21_val_pr_auc_0.8365921057512743.h5')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# gradcam functions source/inspiration: https://github.com/eclique/keras-gradcam/blob/master/grad_cam.py\nlayer_name='conv5_block16_concat'\n\n# for batch gradcam\ngradient_fns = []\nfor class_i in range(len(model_class_names)):\n    class_predictions = tf.slice(model.output, [0, class_i], [-1, 1])\n    conv_layer_output = model.get_layer(layer_name).output\n    grads = K.gradients(class_predictions, conv_layer_output)[0]\n    gradient_fns.append(K.function([model.input, K.learning_phase()], [conv_layer_output, grads]))\n\ndef grad_cam_batch(images, class_i):\n    \"\"\"GradCAM method for visualizing input saliency.\n    Same as grad_cam but processes multiple images in one run.\"\"\"\n    conv_output, grads_val = gradient_fns[class_i]([images, 0])    \n    weights = np.mean(grads_val, axis=(1, 2))\n    cams = np.einsum('ijkl,il->ijk', conv_output, weights)\n    \n    # Process CAMs\n    new_cams = np.empty((images.shape[0], image_height, image_width))\n    for i in range(images.shape[0]):\n        cam_i = cams[i] - cams[i].mean()\n        cam_i = (cam_i + 1e-10) / (np.linalg.norm(cam_i, 2) + 1e-10)\n        new_cams[i] = cv2.resize(cam_i, (image_height, image_width), cv2.INTER_LINEAR)\n        new_cams[i] = np.maximum(new_cams[i], 0)\n        new_cams[i] = new_cams[i] / new_cams[i].max()    \n    return new_cams","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# for a given class probability finds corresponding precision\ntabular_precisions =  np.concatenate((np.arange(0, 0.35, 0.1), np.arange(0.35, 0.6, 0.05), np.arange(0.6, 1, 0.025)))\ndef get_tabular_precisions(class_i, class_probs, class_probability_2_precision):\n    prob_2_precision = class_probability_2_precision[class_i]\n    return np.array([prob_2_precision[min(prob_2_precision.keys(), \n                                          key=lambda x: abs(x - class_prob))] \n                     for class_prob in class_probs])\n\n# credits https://www.kaggle.com/artgor/segmentation-in-pytorch-using-convenient-tools?scriptVersionId=20202006\ndef remove_small_mask_components(masks, min_size):\n    cleaned_masks = np.zeros_like(masks, np.float32)\n    for mask_i in range(masks.shape[0]):\n        mask = masks[mask_i]\n        num_component, component = cv2.connectedComponents(mask.astype(np.uint8))\n        for c in range(1, num_component):\n            p = (component == c)\n            if p.sum() > min_size:\n                cleaned_masks[mask_i][p] = 1\n    return cleaned_masks\n    \ndef generate_gradcam_masks(precision_2_binarization_threshold,\n                        min_size,\n                           imgs_np,\n                        class_probability_2_precision,\n                        model=model, model_class_names=model_class_names\n                       ):\n    gradcam_masks_np = np.empty((len(model_class_names),) + imgs_np.shape[:3], np.float32)\n    zero_mask = np.zeros((image_height, image_width), np.float32)\n    print(f'Generating masks using Grad-CAM')\n    num_batches = imgs_np.shape[0]//batch_size + 1\n    for batch_i in tqdm(range(num_batches)):\n        imgs_batch = imgs_np[batch_i*batch_size: (batch_i + 1)*batch_size]\n        predictions = model.predict(imgs_batch)\n        for class_i, class_name in enumerate(model_class_names):\n            gradcams_batch = grad_cam_batch(imgs_batch, class_i)\n            class_probs = predictions[:, class_i]            \n            class_precisions = get_tabular_precisions(class_i, class_probs, class_probability_2_precision)\n            binarization_thresholds = np.array([precision_2_binarization_threshold[class_precision]\n                                                for class_precision in class_precisions])\n            # to perform elementwise binarization on all images at once\n            binarization_thresholds_batch = np.repeat(binarization_thresholds, \n                                                      image_width*image_height).reshape(imgs_batch.shape[0], image_height, image_width)\n            binarized_gradcams_batch = gradcams_batch.copy()\n            binarized_gradcams_batch = np.where(binarized_gradcams_batch > binarization_thresholds_batch, 1.0, 0.0)\n            binarized_gradcams_batch[imgs_batch[:,:,:,0]==0] = 0\n            binarized_gradcams_batch = remove_small_mask_components(binarized_gradcams_batch, min_size)\n            gradcam_masks_np[class_i][batch_i*batch_size: (batch_i + 1)*batch_size] = binarized_gradcams_batch\n    return gradcam_masks_np\n\ndef visualize_img_gradcam_mask(img_idx, imgs_np, gradcam_masks_np, masks_val_np, names=None):\n    img = imgs_np[img_idx]\n    fig, axes = plt.subplots(2, 2, figsize=(10,10))\n    for class_i, class_name in enumerate(model_class_names):\n        ax = axes[int(class_i>1), class_i%2]\n        ax.set_title(f'GradCAM: {class_name}')\n        ax.axis('off')\n        ax.imshow(img)\n        ax.imshow(gradcam_masks_np[class_i][img_idx], cmap='jet', alpha=0.7)\n        ax.imshow(masks_val_np[class_i][img_idx], alpha=0.35)\n    if not names is None:\n        plt.suptitle(f'Image {names[img_idx]}', fontsize=15)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Precision to threshold mapping per class"},{"metadata":{},"cell_type":"markdown","source":"Using validation data, I map confidence to val-precision. As we've found out [here](https://www.kaggle.com/samusram/cloud-classifier-for-post-processing), some classes are harder than others. I try to unify thresholding while preserving class-specific characteristics by thresholding on precision. "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"class DataGenenerator(Sequence):\n    def __init__(self, images_list=None, folder_imgs=train_imgs_folder, \n                 batch_size=32, shuffle=True, augmentation=None,\n                 resized_height=224, resized_width=224, num_channels=3):\n        self.batch_size = batch_size\n        self.shuffle = shuffle\n        self.augmentation = augmentation\n        if images_list is None:\n            self.images_list = os.listdir(folder_imgs)\n        else:\n            self.images_list = deepcopy(images_list)\n        self.folder_imgs = folder_imgs\n        self.len = len(self.images_list) // self.batch_size\n        self.resized_height = resized_height\n        self.resized_width = resized_width\n        self.num_channels = num_channels\n        self.num_classes = 4\n        self.is_test = not 'train' in folder_imgs\n        if not shuffle and not self.is_test:\n            self.labels = [img_2_ohe_vector[img] for img in self.images_list[:self.len*self.batch_size]]\n\n    def __len__(self):\n        return self.len\n    \n    def on_epoch_start(self):\n        if self.shuffle:\n            random.shuffle(self.images_list)\n\n    def __getitem__(self, idx):\n        current_batch = self.images_list[idx * self.batch_size: (idx + 1) * self.batch_size]\n        X = np.empty((self.batch_size, self.resized_height, self.resized_width, self.num_channels))\n        y = np.empty((self.batch_size, self.num_classes))\n\n        for i, image_name in enumerate(current_batch):\n            path = os.path.join(self.folder_imgs, image_name)\n            img = cv2.resize(cv2.imread(path), (self.resized_height, self.resized_width)).astype(np.float32)\n            if not self.augmentation is None:\n                augmented = self.augmentation(image=img)\n                img = augmented['image']\n            X[i, :, :, :] = img/255.0\n            if not self.is_test:\n                y[i, :] = img_2_ohe_vector[image_name]\n        return X, y\n\n    def get_labels(self):\n        if self.shuffle:\n            images_current = self.images_list[:self.len*self.batch_size]\n            labels = [img_2_ohe_vector[img] for img in images_current]\n        else:\n            labels = self.labels\n        return np.array(labels)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data_generator_val = DataGenenerator(val_imgs, shuffle=False)\ny_pred = model.predict_generator(data_generator_val, workers=num_cores)\ny_true = data_generator_val.get_labels()\n\ndef get_probability_for_precision_threshold(y_true, y_pred, class_i, precision_threshold):\n    precision, recall, thresholds = precision_recall_curve(y_true[:, class_i], y_pred[:, class_i])\n    # consice, even though unnecessary passing through all the values\n    probability_threshold = [thres for prec, thres in zip(precision, thresholds) if prec >= precision_threshold][0]\n    return probability_threshold\n\nclass_probability_2_precision = [dict() for _ in range(len(model_class_names))]\nfor class_i in tqdm(range(len(model_class_names))):\n    for precition_thres in tabular_precisions:\n        class_prob = get_probability_for_precision_threshold(y_true, y_pred, class_i, precition_thres)\n        class_probability_2_precision[class_i][class_prob] = precition_thres","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Vizually comparing GradCAM's with masks"},{"metadata":{},"cell_type":"markdown","source":"I've prepared codes to experiment with different binarization thresholds for different precisions. Yet, for illustration purposes let's produce no masks for precisions under 0.85 and otherwise I use binarization threshold of 0.2 if precision is under 0.9 or 0.01 otherwise."},{"metadata":{"trusted":true},"cell_type":"code","source":"precision_2_binarization_threshold = {prec: 1.0 if prec < 0.85 else 0.2 if prec < 0.9 else 0.01 for prec in tabular_precisions}\nmin_size = 2000\ngradcam_masks_np = generate_gradcam_masks(precision_2_binarization_threshold, min_size, val_imgs_np, class_probability_2_precision)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for img_i in random.sample(list(range(len(val_imgs))), 15):\n    visualize_img_gradcam_mask(img_i, val_imgs_np, gradcam_masks_np, val_masks_np, names=val_imgs)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Estimating performance"},{"metadata":{"trusted":true},"cell_type":"code","source":"# credits: https://www.kaggle.com/artgor/segmentation-in-pytorch-using-convenient-tools?scriptVersionId=20202006\ndef dice(img1, img2):\n    img1 = np.asarray(img1).astype(np.bool)\n    img2 = np.asarray(img2).astype(np.bool)\n\n    intersection = np.logical_and(img1, img2)\n\n    return 2. * intersection.sum() / (img1.sum() + img2.sum())","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def estimate_dice(gradcam_masks_np, val_masks_np=val_masks_np):\n    assert(gradcam_masks_np.shape==val_masks_np.shape)\n    scores_all = []\n    scores_per_class = [[] for _ in range(len(model_class_names))]\n    for class_i, class_name in enumerate(model_class_names):        \n        for img_idx in range(gradcam_masks_np.shape[1]):\n            if gradcam_masks_np[class_i][img_idx].sum() != 0:\n                score = dice(gradcam_masks_np[class_i][img_idx], val_masks_np[class_i][img_idx])\n                scores_all.append(score)\n                scores_per_class[class_i].append(score)\n            else:\n                scores_all.append(val_masks_np[class_i][img_idx].sum() == 0)\n                scores_per_class[class_i].append(val_masks_np[class_i][img_idx].sum() == 0)\n    return np.mean(scores_all), [np.mean(scores) for scores in scores_per_class]  \n\nmean_dice, mean_dice_per_cls = estimate_dice(gradcam_masks_np)\nprint(f\"\"\"Val Mean Dice overall: {mean_dice: .3f}\n\n{chr(10).join([chr(9) + 'Val Mean Dice for class ' + \nclass_name + ' is {:.3f}'.format(mean_dice_per_cls[class_i]) for class_i, class_name in enumerate(model_class_names)])}\"\"\")","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Predicting test masks"},{"metadata":{"trusted":true},"cell_type":"code","source":"del val_imgs_np\ngc.collect()\n\ntest_imgs = os.listdir(test_imgs_folder)\ntest_imgs_np = np.empty((len(test_imgs), image_height, image_width, 3))\nfor img_i, img_name in enumerate(tqdm(test_imgs)):\n    img_path = os.path.join(test_imgs_folder, img_name)\n    test_imgs_np[img_i, :, :, :] = cv2.resize(cv2.imread(img_path), (image_height, image_width)).astype(np.float32)/255.0","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"gradcam_masks_np = generate_gradcam_masks(precision_2_binarization_threshold, 2000, test_imgs_np, class_probability_2_precision)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img_label_list = []\nenc_pixels_list = []\nfor test_img_i, test_img in enumerate(tqdm(test_imgs)):\n    for class_i, class_name in enumerate(model_class_names):\n        img_label_list.append(f'{test_img}_{class_name}')\n        mask = gradcam_masks_np[class_i][test_img_i]\n        if mask.sum() == 0:\n            enc_pixels_list.append(np.nan)\n        else:\n            mask = cv2.resize(mask, dsize=(525, 350), interpolation=cv2.INTER_LINEAR)\n            mask = np.where(mask > 0.5, 1.0, 0.0)\n            enc_pixels_list.append(mask2rle(mask))\nsubmission_df = pd.DataFrame({'Image_Label': img_label_list, 'EncodedPixels': enc_pixels_list})\nsubmission_df.to_csv('sub_gradcam.csv', index=None)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Conclusion"},{"metadata":{},"cell_type":"markdown","source":"It's interesing to experiment with explainability techniques and it's cool to see that a classifier itself has some localization capabilities. "},{"metadata":{},"cell_type":"markdown","source":"# Credits\n* Grad-CAM function was taken from [here](https://github.com/eclique/keras-gradcam/blob/master/gradcam_vgg.ipynb) and was slightly reorganized, \n* mask routines were taken from [Andrew's notebook](https://www.kaggle.com/artgor/segmentation-in-pytorch-using-convenient-tools?scriptVersionId=20202006)."}],"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}