{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Unsupervised Masks - CV 0.60\nIn this kernel, we create segmentation masks **without** using the annotators' training masks! \n\nInstead we build a classifier to classify images into `Fish`, `Flower`, `Gravel`, and `Sugar`. By looking at hundreds of images of `Fish`, the CNN learns to recognize `Fish` in the images. Even though we don't tell the network what specifically is a `Fish` within the image, the CNN determines where and what `Fish` are by understanding the similarity between many images containing `Fish`. This feels like magic! There is a great blog why this works [here][1]. Some of this notebook's code was taken from the blogger's GitHub [here][2]\n\n[1]: https://alexisbcook.github.io/2017/global-average-pooling-layers-for-object-localization/\n[2]: https://github.com/alexisbcook/ResNetCAM-keras"},{"metadata":{},"cell_type":"markdown","source":"# Load Image Labels\nFrom the training data, we will not use the annotators' masks. Instead we will only use a label for each image indicating whether `Fish`, `Flower`, `Gravel`, and/or `Sugar` clouds are present. Note that we don't tell the network where the `Fish`, `Flowers`, `Gravel` nor `Sugar` are located within the images!"},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"!pip install tensorflow-gpu==1.14.0 --quiet\n!pip install keras==2.2.4 --quiet\n\nimport keras\nimport numpy as np, pandas as pd, os \nfrom keras import layers\nfrom keras.models import Model\nfrom PIL import Image\nfrom keras import optimizers\nimport scipy, cv2   \nimport matplotlib.pyplot as plt \nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import roc_auc_score, accuracy_score","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":false,"trusted":true},"cell_type":"code","source":"train = pd.read_csv('../input/understanding_cloud_organization/train.csv')\ntrain['Image'] = train['Image_Label'].map(lambda x: x.split('.')[0])\ntrain['Label'] = train['Image_Label'].map(lambda x: x.split('_')[1])\ntrain2 = pd.DataFrame({'Image':train['Image'][::4]})\ntrain2['e1'] = train['EncodedPixels'][::4].values\ntrain2['e2'] = train['EncodedPixels'][1::4].values\ntrain2['e3'] = train['EncodedPixels'][2::4].values\ntrain2['e4'] = train['EncodedPixels'][3::4].values\ntrain2.set_index('Image',inplace=True,drop=True)\ntrain2.fillna('',inplace=True); train2.head()\ntrain2[['d1','d2','d3','d4']] = (train2[['e1','e2','e3','e4']]!='').astype('int8')\ntrain2[['d1','d2','d3','d4']].head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"os.listdir('../input/understanding_cloud_organization/test_images')\n# train['Image'] = train['Image_Label'].map(lambda x: x.split('.')[0])\n# train['Label'] = train['Image_Label'].map(lambda x: x.split('_')[1])\n# train2 = pd.DataFrame({'Image':train['Image'][::4]})\n# train2['e1'] = train['EncodedPixels'][::4].values\n# train2['e2'] = train['EncodedPixels'][1::4].values\n# train2['e3'] = train['EncodedPixels'][2::4].values\n# train2['e4'] = train['EncodedPixels'][3::4].values\n# train2.set_index('Image',inplace=True,drop=True)\n# train2.fillna('',inplace=True); train2.head()\n# train2[['d1','d2','d3','d4']] = (train2[['e1','e2','e3','e4']]!='').astype('int8')\n# train2[['d1','d2','d3','d4']].head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Helper Functions\nFunctions to help manipulate masks and generate data."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"def rle2maskX(mask_rle, shape=(2100,1400), shrink=1):\n    # Converts rle to mask size shape then downsamples by shrink\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).T[::shrink,::shrink]\n\ndef rle2mask2X(mask_rle, shape=(2100,1400), shrink=(512,352)):\n    # Converts rle to mask size shape then downsamples by shrink\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    img = img.reshape(shape).T\n    img = Image.fromarray(img)\n    img = img.resize(shrink)\n    img = np.array(img)\n    return img\n\ndef mask2contour(mask, width=5):\n    w = mask.shape[1]\n    h = mask.shape[0]\n    mask2 = np.concatenate([mask[:,width:],np.zeros((h,width))],axis=1)\n    mask2 = np.logical_xor(mask,mask2)\n    mask3 = np.concatenate([mask[width:,:],np.zeros((width,w))],axis=0)\n    mask3 = np.logical_xor(mask,mask3)\n    return np.logical_or(mask2,mask3) \n\ndef mask2rle(img, shape=(525,350)):    \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\ndef dice_coef6(y_true_rle, y_pred_rle, y_pred_prob, th):\n    if y_pred_prob<th:\n        if y_true_rle=='': return 1\n        else: return 0\n    else:\n        y_true_f = rle2maskX(y_true_rle,shrink=4)\n        y_pred_f = rle2maskX(y_pred_rle,shape=(525,350))\n        union = np.sum(y_true_f) + np.sum(y_pred_f)\n        if union==0: return 1\n        intersection = np.sum(y_true_f * y_pred_f)\n        return 2. * intersection / union\n\ndef dice_coef8(y_true_f, y_pred_f):\n    union = np.sum(y_true_f) + np.sum(y_pred_f)\n    if union==0: return 1\n    intersection = np.sum(y_true_f * y_pred_f)\n    return 2. * intersection / union\n\nclass DataGenerator(keras.utils.Sequence):\n    # USES GLOBAL VARIABLE TRAIN2 COLUMNS E1, E2, E3, E4\n    'Generates data for Keras'\n    def __init__(self, list_IDs, batch_size=8, shuffle=False, width=512, height=352, scale=1/128., sub=1., mode='train',\n                 path='../input/understanding_cloud_organization/train_images/', flips=False):\n        'Initialization'\n        self.list_IDs = list_IDs\n        self.shuffle = shuffle\n        self.batch_size = batch_size\n        self.path = path\n        self.scale = scale\n        self.sub = sub\n        self.path = path\n        self.width = width\n        self.height = height\n        self.mode = mode\n        self.flips = flips\n        self.on_epoch_end()\n        \n    def __len__(self):\n        'Denotes the number of batches per epoch'\n        ct = int(np.floor( len(self.list_IDs) / self.batch_size))\n        if len(self.list_IDs)>ct*self.batch_size: ct += 1\n        return int(ct)\n\n    def __getitem__(self, index):\n        'Generate one batch of data'\n        indexes = self.indexes[index*self.batch_size:(index+1)*self.batch_size]\n        X, y = self.__data_generation(indexes)\n        if (self.mode=='train')|(self.mode=='validate'): return X, y\n        else: return X\n\n    def on_epoch_end(self):\n        'Updates indexes after each epoch'\n        self.indexes = np.arange(int( len(self.list_IDs) ))\n        if self.shuffle: np.random.shuffle(self.indexes)\n\n    def __data_generation(self, indexes):\n        'Generates data containing batch_size samples' \n        # Initialization\n        lnn = len(indexes)\n        X = np.empty((lnn,self.height,self.width,3),dtype=np.float32)\n        y = np.zeros((lnn,4),dtype=np.int8)\n        \n        # Generate data\n        for k in range(lnn):\n            img = cv2.imread(self.path + self.list_IDs[indexes[k]] + '.jpg')\n            img = cv2.resize(img,(self.width,self.height),interpolation = cv2.INTER_AREA)\n            # AUGMENTATION FLIPS\n            hflip = False; vflip = False\n            if (self.flips):\n                if np.random.uniform(0,1)>0.5: hflip=True\n                if np.random.uniform(0,1)>0.5: vflip=True\n            if vflip: img = cv2.flip(img,0) # vertical\n            if hflip: img = cv2.flip(img,1) # horizontal\n            # NORMALIZE IMAGES\n            X[k,] = img*self.scale - self.sub      \n            # LABELS\n            if (self.mode=='train')|(self.mode=='validate'):\n                y[k,] = train2.loc[self.list_IDs[indexes[k]],['d1','d2','d3','d4']].values\n            \n        return X, y","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Train Classifier Model\nWe will build an Xception model by removing it's ImageNet top and adding our own top consisting of one Global Average Pooling Layer and one Dense Layer with 4 sigmoid units. Our classifier is now a fully convolutional classifier. The base is pretrained on ImageNet data and we will train our new top on clouds.\n\nThe layer preceeding the Global Average Pooling Layer (i.e. the top layer of base model) will have dimensions of the input shape divided by 32 because Xception does five 2x downsamplings. Our input shape will be `352x512` with 3 maps and this reduces to `11x16` with 2048 maps. Each of these 2048 maps is like a segmentation mask which specializes in spacially locating a certain type of pattern in the original image. (We'll discuss this more later)."},{"metadata":{"trusted":true},"cell_type":"code","source":"# USE KERAS XCEPTION MODEL\nfrom keras.applications.xception import Xception\nbase_model = Xception(weights='imagenet',include_top=False,input_shape=(None,None,3))\n# FREEZE NON-BATCHNORM LAYERS IN BASE\nfor layer in base_model.layers:\n    if not isinstance(layer, layers.BatchNormalization): layer.trainable = False\n# BUILD MODEL NEW TOP\nx = base_model.output\nx = layers.GlobalAveragePooling2D()(x)\nx = layers.Dense(4,activation='sigmoid')(x)\nmodel = Model(inputs=base_model.input, outputs=x)\n# COMPILE MODEL\nmodel.compile(loss='binary_crossentropy', optimizer = optimizers.Adam(lr=0.001), metrics=['accuracy'])\n\n# SPLIT TRAIN AND VALIDATE\nidxT, idxV = train_test_split(train2.index, random_state=42, test_size=0.2)\ntrain_gen = DataGenerator(idxT, flips=True, shuffle=True)\nval_gen = DataGenerator(idxV, mode='validate')\n\n# TRAIN NEW MODEL TOP LR=0.001 (with bottom frozen)\nh = model.fit_generator(train_gen, epochs = 2, verbose=2, validation_data = val_gen)\n# TRAIN ENTIRE MODEL LR=0.0001 (with all unfrozen)\nfor layer in model.layers: layer.trainable = True\nmodel.compile(loss='binary_crossentropy', optimizer = optimizers.Adam(lr=0.0001), metrics=['accuracy'])\nh = model.fit_generator(train_gen, epochs = 2, verbose=2, validation_data = val_gen)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"base_model.output","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Evaluate Validation Accuracy\nWe see that our classifier has a high accuracy of 75%"},{"metadata":{"trusted":true},"cell_type":"code","source":"# PREDICT HOLDOUT SET\ntrain3 = train2.loc[train2.index.isin(idxV)].copy()\noof_gen = DataGenerator(train3.index.values, mode='predict')\noof = model.predict_generator(oof_gen, verbose=2)\nfor k in range(1,5): train3['o'+str(k)] = 0\ntrain3[['o1','o2','o3','o4']] = oof\n\n# COMPUTE ACCURACY AND ROC_AUC_SCORE\ntypes = ['Fish','Flower','Gravel','Sugar']\nfor k in range(1,5):\n    print(types[k-1],': ',end='')\n    auc = np.round( roc_auc_score(train3['d'+str(k)].values,train3['o'+str(k)].values  ),3 )\n    acc = np.round( accuracy_score(train3['d'+str(k)].values,(train3['o'+str(k)].values>0.5).astype(int) ),3 )\n    print('AUC =',auc,end='')\n    print(', ACC =',acc) \nprint('OVERALL: ',end='')\nauc = np.round( roc_auc_score(train3[['d1','d2','d3','d4']].values.reshape((-1)),train3[['o1','o2','o3','o4']].values.reshape((-1)) ),3 )\nacc = np.round( accuracy_score(train3[['d1','d2','d3','d4']].values.reshape((-1)),(train3[['o1','o2','o3','o4']].values>0.5).astype(int).reshape((-1)) ),3 )\nprint('AUC =',auc, end='')\nprint(', ACC =',acc) ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Display Class Activation Maps\nEarlier we discussed how the layer preceeding the Global Average Pooling Layer has 2048 maps of size `11x16`. Each of these maps is like a segmentation mask that specializes in detecting a certain pattern. Let's imagine that detecting `Sugar` requires the use of maps 1, 45, 256, and 1039 where 1 detects small white shapes, 45 detects circular objects, 256 detects a lattice arrangement, and 1039 detects blurred edges. Then the segmentation map for `Sugar` is the addition of these 4 maps.\n\nBelow are 25 rows of images. The images on the left are the class activation maps (formed by summing all relevent pattern detection maps from among the 2048 possible maps). These CAMs result from feeding the associated image into our Xception CNN. For each row, we display the cloud type that activated the strongest. The images of the right are the true mask in yellow and the activation map converted to a mask in blue. Remember our network has never seen the annotators' training masks! But none-the-less, our model has found masks! Truly magical"},{"metadata":{"trusted":true},"cell_type":"code","source":" np.dot(last_conv_output.reshape((16*11, 2048)), layer_weights).shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# NEW MODEL FROM OLD TO EXTRACT ACTIVATION MAPS\nall_layer_weights = model.layers[-1].get_weights()[0]\ncam_model = Model(inputs=model.input, \n        outputs=(model.layers[-3].output, model.layers[-1].output)) \n\n# DISPLAY 25 RANDOM IMAGES\nPATH = '../input/understanding_cloud_organization/train_images/'\nIMGS = os.listdir(PATH)\nfor k in np.random.randint(0,5000,25):\n    \n    # LOAD IMAGE AND PREDICT CLASS ACTIVATION MAP\n    img = cv2.resize( cv2.imread(PATH+IMGS[k]), (512, 352))\n    x = np.expand_dims(img, axis=0)/128. -1.\n    last_conv_output, pred_vec = cam_model.predict(x) \n    last_conv_output = np.squeeze(last_conv_output) \n    pred = np.argmax(pred_vec)\n    layer_weights = all_layer_weights[:, pred] \n    final_output = np.dot(last_conv_output.reshape((16*11, 2048)), layer_weights).reshape(11,16) \n    final_output = scipy.ndimage.zoom(final_output, (32, 32), order=1) \n\n    # DISPLAY IMAGE WITH CLASS ACTIVATION MAPS\n    plt.figure(figsize=(12,6))\n    plt.subplot(1,2,1)\n    mx = np.round( np.max(final_output),1 )\n    mn = np.round( np.min(final_output),1 )\n    final_output = (final_output-mn)/(mx-mn)\n    mask0 = (final_output>0.3).astype(int)\n    contour0 = mask2contour(mask0,5)\n    plt.imshow(img, alpha=0.5)\n    plt.imshow(final_output, cmap='jet', alpha=0.5)\n    plt.title('Found '+types[pred]+'  -  Pr = '+str(np.round(pred_vec[0,pred],3)) )\n    \n    # DISPLAY IMAGE WITH MASKS\n    plt.subplot(1,2,2)\n    img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n    rle = train2.loc[IMGS[k].split('.')[0],'e'+str(pred+1)]\n    mask = rle2mask2X(rle,shrink=(512,352))\n    contour = mask2contour(mask,5)\n    img[contour==1,:2] = 255\n    img[contour0==1,2] = 255\n    diff = np.ones((352,512,3),dtype=np.int)*255-img\n    img=img.astype(int); img[mask0==1,:] += diff[mask0==1,:]//4\n    plt.imshow( img )\n    dice = np.round( dice_coef8(mask,mask0),3 )\n    plt.title('Dice = '+str(dice)+'  -  '+IMGS[k]+'  -  '+types[pred])\n    \n    plt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img.shape\nx.shape\nlast_conv_output.shape\ncam_model.predict(x)[0].shape\nall_layer_weights.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"last_conv_output.shape","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Evaluate Validation Dice\nIf we use these class activation maps as a prediction for the annotators' masks in Kaggle's Cloud competition, the Dice score validates over 0.600 !! Not bad for unsupervised learning. One can even argue that these segmentation masks are more accurate than the annotators' masks."},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Computing',len(train3),'masks...')\nfor i in range(1,5): train3['p'+str(i)] = ''\nfor i in range(1,5): train3['pp'+str(i)] = 0\n\nfor i,f in enumerate(train3.index.values):\n    \n    # LOAD IMAGE AND PREDICT CLASS ACTIVATION MAPS\n    img = cv2.resize( cv2.imread(PATH+f+'.jpg'), (512, 352))\n    x = np.expand_dims(img, axis=0)/128. -1.\n    last_conv_output, pred_vec = cam_model.predict(x) \n    last_conv_output = np.squeeze(last_conv_output) \n    \n    for pred in [0,1,2,3]:\n        # CREATE FOUR MASKS FROM ACTIVATION MAPS\n        layer_weights = all_layer_weights[:, pred]  \n        final_output = np.dot(last_conv_output.reshape((16*11, 2048)), layer_weights).reshape(11,16) \n        final_output = scipy.ndimage.zoom(final_output, (32, 32), order=1)\n        mx = np.round( np.max(final_output),1 )\n        mn = np.round( np.min(final_output),1 )\n        final_output = (final_output-mn)/(mx-mn)\n        final_output = cv2.resize(final_output,(525,350))\n        train3.loc[f,'p'+str(pred+1)] = mask2rle( (final_output>0.3).astype(int) )\n        train3.loc[f,'pp'+str(pred+1)] = pred_vec[0,pred]\n    if i%25==0: print(i,', ',end='')\nprint(); print()\n        \n# COMPUTE KAGGLE DICE\nth = [0.8,0.5,0.7,0.7]\nfor k in range(1,5):\n    train3['ss'+str(k)] = train3.apply(lambda x:dice_coef6(x['e'+str(k)],x['p'+str(k)],x['pp'+str(k)],th[k-1]),axis=1)\n    dice = np.round( train3['ss'+str(k)].mean(),3 )\n    print(types[k-1],': Kaggle Dice =',dice)\ndice = np.round( np.mean( train3[['ss1','ss2','ss3','ss4']].values ),3 )\nprint('Overall : Kaggle Dice =',dice)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# right, one more session in which to make a submission. then, back to those kpis","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print('Computing',len(train3),'masks...')\nfor i in range(1,5): train3['p'+str(i)] = ''\nfor i in range(1,5): train3['pp'+str(i)] = 0\n\nsub = pd.DataFrame(columns=['EncodedPixels'])\n\npred_names = {0: 'Fish', 1: 'Flower', 2: 'Gravel', 3: 'Sugar'}\n    \nth = [0.8,0.5,0.7,0.7]\n    \nfor i,f in enumerate(os.listdir('../input/understanding_cloud_organization/test_images')):\n    \n    # LOAD IMAGE AND PREDICT CLASS ACTIVATION MAPS\n    #print(f'../input/understanding_cloud_organization/test_images/{f}')\n    img = cv2.resize( cv2.imread(f'../input/understanding_cloud_organization/test_images/{f}'), (512, 352))\n    x = np.expand_dims(img, axis=0)/128. -1.\n    last_conv_output, pred_vec = cam_model.predict(x) \n    last_conv_output = np.squeeze(last_conv_output) \n    \n    for pred in [0,1,2,3]:\n        # CREATE FOUR MASKS FROM ACTIVATION MAPS\n        layer_weights = all_layer_weights[:, pred]  \n        final_output = np.dot(last_conv_output.reshape((16*11, 2048)), layer_weights).reshape(11,16) \n        final_output = scipy.ndimage.zoom(final_output, (32, 32), order=1)\n        mx = np.round( np.max(final_output),1 )\n        mn = np.round( np.min(final_output),1 )\n        final_output = (final_output-mn)/(mx-mn)\n        final_output = cv2.resize(final_output,(525,350))\n        #sub.loc[f,'p'+str(pred+1)] = mask2rle( (final_output>0.3).astype(int) )\n        \n        if pred_vec[0,pred] > th[pred]:\n            sub.loc[f'{f}_{pred_names[pred]}', 'EncodedPixels'] = mask2rle( (final_output>0.3).astype(int) )\n        else:\n            sub.loc[f'{f}_{pred_names[pred]}', 'EncodedPixels'] = ''\n    if i%25==0: print(i,', ',end='')\nprint(); print()\n        \n# COMPUTE KAGGLE DICE\n\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.reset_index().rename({'index': 'Image_Label'}, axis=1) # this is good for learning. this is how you learn to do ML - by doing it. so...keep doing this without guilt. once it's over, you can go back to the KPIs. besides, Kaggle\n# is how you've been learning all this stuff, and if they can't see that, they don't know what's good for them","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.reset_index().rename({'index': 'Image_Label'}, axis=1).to_csv('submission.csv', index=False)","execution_count":null,"outputs":[]}],"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":1}