{"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_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"<img src=\"https://markdown.data-ensta.fr/uploads/upload_62aae64efd54f4a28389617d2ff98804.png\" alt=\"drawing\" width=\"1400\"/>\n<center><h1 style=\"color:Blue;font-size:40px;\">HuBMAP 2022</h1> </center>\n<center><h1 style=\"color:Blue;font-size:27px;\">Vanilla Submission with Keras</h1> </center>","metadata":{}},{"cell_type":"markdown","source":"## Gist & acknowledgements","metadata":{}},{"cell_type":"markdown","source":"+ This notebook proposes a simple Keras pipeline to make a quick submission\n+ It relies heavily on @quvbel's [segmentation library](https://github.com/qubvel/segmentation_models)\n+ Thanks for @thedevastator's heavy lifting, this year's first [tiled dataset](https://www.kaggle.com/datasets/thedevastator/hubmap-2022-256x256) can be used early in the competition\n+ This notebook is an adaptation of team @kishan98joshi's awesome [notebook](https://www.kaggle.com/code/joshi98kishan/hubmap-keras-pipeline-training-inference).\n+ Cover picture from [NIH's website](https://commonfund.nih.gov/hubmap/image-of-the-week).","metadata":{}},{"cell_type":"code","source":"# Keras utilities (don't forget to \"Add data\" in the side menu if content not there, or else failure)\n!pip install -U ../input/kerasapplications/Keras_Applications-1.0.8-py3-none-any.whl\n!pip install ../input/qubvel/efficientnet-1.0.0-py3-none-any.whl\n!pip install ../input/qubvel/image_classifiers-1.0.0-py3-none-any.whl\n!pip install ../input/qubvel-segmentation-model-keras-v101/segmentation_models-master\n\n%env SM_FRAMEWORK=tf.keras","metadata":{"_kg_hide-input":false,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-06-23T05:33:01.622875Z","iopub.execute_input":"2022-06-23T05:33:01.623215Z","iopub.status.idle":"2022-06-23T05:34:50.818236Z","shell.execute_reply.started":"2022-06-23T05:33:01.623155Z","shell.execute_reply":"2022-06-23T05:34:50.817032Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Imports","metadata":{}},{"cell_type":"code","source":"import os, gc\nimport numpy as np \nimport pandas as pd \nimport cv2\nimport glob\nimport numba\nimport pathlib\nfrom tqdm.notebook import tqdm\nimport tifffile as tiff\nimport tensorflow as tf\nimport tensorflow.keras.backend as K\nimport matplotlib.pyplot as plt\nfrom sklearn.model_selection import GroupKFold\n\nfrom albumentations import *\nimport segmentation_models as sm\nimport rasterio\nfrom rasterio.windows import Window\n\nimport warnings\nwarnings.filterwarnings('ignore')","metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-23T05:34:50.820906Z","iopub.execute_input":"2022-06-23T05:34:50.821474Z","iopub.status.idle":"2022-06-23T05:34:53.322744Z","shell.execute_reply.started":"2022-06-23T05:34:50.821426Z","shell.execute_reply":"2022-06-23T05:34:53.321851Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Loading and formatting utilities ","metadata":{}},{"cell_type":"code","source":"def get_settings(batch_size = 32, \n                 encoder = 'resnet34',\n                 epochs = 5):\n    \n    return {'BATCH_SIZE' : batch_size,\n           'ENCODER' : encoder,\n           'EPOCHS' : epochs}\n\n\ndef get_path(pathof):\n    \"\"\"\n    Returns mapping function based on the input argument\n    \n    Arguments:\n            pathof : value can be either 'image' or 'mask'\n    \"\"\"\n\n    if pathof=='image':\n        return lambda fname: tf.strings.join([path_train, fname], \n                                             separator = '/')\n    \n    elif pathof=='mask':\n        return lambda fname: tf.strings.join([path_masks, fname], \n                                              separator = '/') \n\ndef get_normalised_tensor(path):\n    \"\"\"\n    Reads the image and scale it to [0, 1]\n    \n    Arguments:\n            path : path of the image to be read\n            \n    Returns:\n            Returns normalized image tensor\n    \"\"\"\n    \n    im = tf.io.read_file(path)\n    im = tf.io.decode_png(im)\n    im = tf.image.convert_image_dtype(im, \n                                      tf.float32)\n    return im\n\ndef get_tensor(path): \n    \"\"\"\n    Reads the image or mask\n    \n    Arguments:\n            path : path of the image or mask to be read\n            \n    Returns:\n            Returns image or mask tensor\n    \"\"\"\n\n    mask = tf.io.read_file(path)\n    mask = tf.io.decode_png(mask)\n    \n    return mask\n\ndef augment_data(image, mask):\n    \"\"\"\n    This is a mapping function on a zipped dataset.\n    \n    It takes tf tensors, do numpy based albumentation transformations \n    and returns augmented tf tensors.\n    \n    Reference: https://albumentations.ai/docs/examples/tensorflow-example/\n    \"\"\"\n    \n    def _fn(image, mask):\n        sample = transforms(image = image, \n                            mask = mask)\n        aug_img = sample['image']\n        aug_msk = sample['mask']\n        \n        aug_img = tf.cast(aug_img/255.0, tf.float32)\n        \n        return aug_img, aug_msk\n\n    aug_img, aug_msk = tf.numpy_function(func = _fn, \n                                         inp = [image, mask],\n                                         Tout = [tf.float32, tf.uint8])\n    \n    return aug_img, aug_msk\n\n\ndef dice_coeff(y_true, y_pred, epsilon=1.):\n    \n    \"\"\"\n    Calculates dice coefficient\n\n    Arguments: \n            y_true : tensor of ground truth values.\n            y_pred : tensor of predicted values.\n            epsilon : constant to avoid divide by 0 errors.\n    \n    Returns:\n            dice_coefficient\n    \"\"\"\n    \n    y_true_f = K.flatten(y_true)\n    y_pred_f = K.flatten(y_pred)\n    intersection = K.sum(y_true_f * y_pred_f)\n    score = (2. * intersection + epsilon) / (K.sum(y_true_f) + K.sum(y_pred_f) + epsilon)\n    return score\n\n\n# https://www.kaggle.com/leighplt/pytorch-fcn-resnet50\ndef make_grid(shape, window=256, min_overlap=32):\n    \"\"\"\n        Return Array of size (N,4), where N - number of tiles,\n        2nd axis represente slices: x1,x2,y1,y2 \n    \"\"\"\n    x, y = shape\n    nx = x // (window - min_overlap) + 1\n    x1 = np.linspace(0, x, num=nx, endpoint=False, dtype=np.int64)\n    x1[-1] = x - window\n    x2 = (x1 + window).clip(0, x)\n    ny = y // (window - min_overlap) + 1\n    y1 = np.linspace(0, y, num=ny, endpoint=False, dtype=np.int64)\n    y1[-1] = y - window\n    y2 = (y1 + window).clip(0, y)\n    slices = np.zeros((nx,ny, 4), dtype=np.int64)\n    \n    for i in range(nx):\n        for j in range(ny):\n            slices[i,j] = x1[i], x2[i], y1[j], y2[j]    \n    return slices.reshape(nx*ny,4)\n\n@numba.njit()\ndef rle_numba(pixels):\n    size = len(pixels)\n    points = []\n    if pixels[0] == 1: points.append(0)\n    flag = True\n    for i in range(1, size):\n        if pixels[i] != pixels[i-1]:\n            if flag:\n                points.append(i+1)\n                flag = False\n            else:\n                points.append(i+1 - points[-1])\n                flag = True\n    if pixels[-1] == 1: points.append(size-points[-1]+1)    \n    return points\n\ndef rle_numba_encode(image):\n    pixels = image.flatten(order = 'F')\n    points = rle_numba(pixels)\n    return ' '.join(str(x) for x in points)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-06-23T05:34:53.324367Z","iopub.execute_input":"2022-06-23T05:34:53.324740Z","iopub.status.idle":"2022-06-23T05:34:53.356472Z","shell.execute_reply.started":"2022-06-23T05:34:53.324700Z","shell.execute_reply":"2022-06-23T05:34:53.355535Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Training set up","metadata":{}},{"cell_type":"code","source":"# TRAINING SETTINGS\nsettings = get_settings(batch_size = 32,\n                        epochs = 30,\n                        encoder = 'seresnext50')\n\nDATA_ORIG_PATH = '../input/hubmap-organ-segmentation'\nDATA_PATH = '../input/hubmap-2022-256x256'\n\npath_train = os.path.join(DATA_PATH, 'train')\npath_masks = os.path.join(DATA_PATH, 'masks')\n\npath_submission_file = os.path.join(DATA_ORIG_PATH, 'sample_submission.csv')\n\nprint(f'No. of training images: {len(os.listdir(path_train))}')\nprint(f'No. of masks: {len(os.listdir(path_masks))} \\n')","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:34:53.358536Z","iopub.execute_input":"2022-06-23T05:34:53.359383Z","iopub.status.idle":"2022-06-23T05:34:53.379473Z","shell.execute_reply.started":"2022-06-23T05:34:53.359301Z","shell.execute_reply":"2022-06-23T05:34:53.378803Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fnames = np.array(os.listdir(path_train))\ngroups = [fname[:9] for fname in fnames]\ngroup_kfold = GroupKFold(n_splits = 4)\n\n# https://www.tensorflow.org/api_docs/python/tf/data/Dataset#shuffle\nBUFFER_SIZE = 1000\nBATCH_SIZE = settings['BATCH_SIZE']","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:34:53.382574Z","iopub.execute_input":"2022-06-23T05:34:53.382899Z","iopub.status.idle":"2022-06-23T05:34:53.394450Z","shell.execute_reply.started":"2022-06-23T05:34:53.382865Z","shell.execute_reply":"2022-06-23T05:34:53.393719Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Model","metadata":{}},{"cell_type":"code","source":"# We will be using pretrained resnet34 as our encoder.\n\n!mkdir -p ~/.keras/models\n\nif settings['ENCODER']=='resnet34':\n    !cp ../input/keras-pretrained-imagenet-weights/resnet34_imagenet_1000_no_top.h5 ~/.keras/models/resnet34_imagenet_1000_no_top.h5\nelif settings['ENCODER']=='resnet50':\n    !cp ../input/keras-pretrained-imagenet-weights/resnet50_imagenet_1000_no_top.h5 ~/.keras/models/resnet50_imagenet_1000_no_top.h5\nelif settings['ENCODER']=='seresnext50':\n    !cp ../input/keras-pretrained-imagenet-weights/seresnext50_imagenet_1000_no_top.h5 ~/.keras/models/seresnext50_imagenet_1000_no_top.h5","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:34:53.397547Z","iopub.execute_input":"2022-06-23T05:34:53.397813Z","iopub.status.idle":"2022-06-23T05:34:55.223822Z","shell.execute_reply.started":"2022-06-23T05:34:53.397786Z","shell.execute_reply":"2022-06-23T05:34:55.222590Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Data Augmentation","metadata":{}},{"cell_type":"code","source":"# https://www.kaggle.com/iafoss/hubmap-pytorch-fast-ai-starter\ntransforms = Compose([\n            HorizontalFlip(),\n            VerticalFlip(),\n            RandomRotate90(),\n            ShiftScaleRotate(shift_limit=0.0625, scale_limit=0.2, rotate_limit=15, p=0.9, \n                             border_mode=cv2.BORDER_REFLECT),\n            OneOf([\n                OpticalDistortion(p=0.3),\n                GridDistortion(p=.1),\n                IAAPiecewiseAffine(p=0.3),\n            ], p=0.3),\n            OneOf([\n                HueSaturationValue(10,15,10),\n                CLAHE(clip_limit=2),\n                RandomBrightnessContrast(),            \n            ], p=0.3),\n        ], p = 1.0)","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:34:55.226470Z","iopub.execute_input":"2022-06-23T05:34:55.227058Z","iopub.status.idle":"2022-06-23T05:34:55.235782Z","shell.execute_reply.started":"2022-06-23T05:34:55.227009Z","shell.execute_reply":"2022-06-23T05:34:55.234674Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Training","metadata":{}},{"cell_type":"code","source":"for fold, (t_idx, v_idx) in enumerate(group_kfold.split(fnames, \n                                                        groups = groups)):\n    \n    print(f'Fold: {fold+1}')\n    \n    t_fnames_ds = tf.data.Dataset.from_tensor_slices(fnames[t_idx])\n    v_fnames_ds = tf.data.Dataset.from_tensor_slices(fnames[v_idx])\n\n    t_img_ds = t_fnames_ds.map(get_path('image')).map(get_tensor)\n    t_msk_ds = t_fnames_ds.map(get_path('mask')).map(get_tensor)\n    \n    v_img_ds = v_fnames_ds.map(get_path('image')).map(get_normalised_tensor)\n    v_msk_ds = v_fnames_ds.map(get_path('mask')).map(get_tensor)\n    \n    train_ds = tf.data.Dataset.zip((t_img_ds, t_msk_ds))\n    train_ds = train_ds.map(augment_data)  \n    val_ds = tf.data.Dataset.zip((v_img_ds, v_msk_ds))\n    del t_fnames_ds, v_fnames_ds, t_img_ds, t_msk_ds, v_img_ds, v_msk_ds\n    gc.collect()\n    train_ds = train_ds.shuffle(BUFFER_SIZE)\\\n                       .batch(BATCH_SIZE)\\\n                       .repeat()\n    val_ds = val_ds.batch(BATCH_SIZE)\\\n                   .repeat()\n    model = sm.Unet(settings['ENCODER'], \n                encoder_weights='imagenet')\n    model.compile(optimizer = 'adam',\n                  loss = tf.keras.losses.BinaryCrossentropy(),\n                  metrics = [dice_coeff])\n    \n    checkpoint_callback = tf.keras.callbacks.ModelCheckpoint(f'./saved_models/fold_model_{fold+1}.pb',\n                                                            save_best_only = True)\n    EPOCHS = settings['EPOCHS']\n\n    TOTAL_TRAIN_SAMPLES = len(t_idx)\n    TOTAL_VAL_SAMPLES = len(t_idx)\n    STEPS_PER_EPOCH = TOTAL_TRAIN_SAMPLES//BATCH_SIZE\n    \n    VALIDATION_STEPS = TOTAL_VAL_SAMPLES//BATCH_SIZE\n\n    history = model.fit(train_ds, \n                        epochs = EPOCHS,\n                        steps_per_epoch = STEPS_PER_EPOCH,\n                        validation_data = val_ds,\n                        validation_steps = VALIDATION_STEPS,\n                        callbacks = [checkpoint_callback])\n    del train_ds, val_ds, model\n    gc.collect()\n    \n    # Not running the last iteration\n    if fold==0:\n        break","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-06-23T05:34:55.239064Z","iopub.execute_input":"2022-06-23T05:34:55.239443Z","iopub.status.idle":"2022-06-23T05:40:11.756229Z","shell.execute_reply.started":"2022-06-23T05:34:55.239409Z","shell.execute_reply":"2022-06-23T05:40:11.755354Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Inference","metadata":{}},{"cell_type":"code","source":"identity = rasterio.Affine(1, 0, 0, 0, 1, 0)\n\n# Size of the tile to be read by rasterio\nWINDOW = 1024\n\nMIN_OVERLAP = 32\n\n# Tiles will be resized to NEW_SIZE, which is the size of the image\n# on which the model was trained\nNEW_SIZE = 256","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:40:11.760694Z","iopub.execute_input":"2022-06-23T05:40:11.761666Z","iopub.status.idle":"2022-06-23T05:40:11.769755Z","shell.execute_reply.started":"2022-06-23T05:40:11.761632Z","shell.execute_reply":"2022-06-23T05:40:11.769038Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!ls ./saved_models","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:40:11.771167Z","iopub.execute_input":"2022-06-23T05:40:11.771692Z","iopub.status.idle":"2022-06-23T05:40:13.571904Z","shell.execute_reply.started":"2022-06-23T05:40:11.771647Z","shell.execute_reply":"2022-06-23T05:40:13.571017Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fold_models = []\n\nfor fold_model_path in glob.glob(os.path.join('./saved_models/*.pb')):\n    fold_models.append(tf.keras.models.load_model(fold_model_path, custom_objects={'dice_coeff': dice_coeff}))","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:40:13.573918Z","iopub.execute_input":"2022-06-23T05:40:13.574249Z","iopub.status.idle":"2022-06-23T05:45:40.355874Z","shell.execute_reply.started":"2022-06-23T05:40:13.574218Z","shell.execute_reply":"2022-06-23T05:45:40.354905Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"p = pathlib.Path(DATA_ORIG_PATH)\nsubm = {}\n\nfor i, filename in tqdm(enumerate(p.glob('test_images/*.tiff')), \n                        total = len(list(p.glob('test_images/*.tiff')))):\n    \n    print(f'{i+1} Predicting {filename.stem}')\n    \n    dataset = rasterio.open(filename.as_posix(), transform = identity)\n    slices = make_grid(dataset.shape, window=WINDOW, min_overlap=MIN_OVERLAP)\n    preds = np.zeros(dataset.shape, dtype=np.uint8)\n    \n    for (x1,x2,y1,y2) in slices:\n        image = dataset.read([1,2,3],\n                    window=Window.from_slices((x1,x2),(y1,y2)))\n        image = np.moveaxis(image, 0, -1)\n        \n        image = tf.image.convert_image_dtype(image, \n                                 tf.float32)\n        image = cv2.resize(image.numpy(), (NEW_SIZE, NEW_SIZE))\n        image = np.expand_dims(image, 0)\n        \n        pred = None\n        \n        for fold_model in fold_models:\n            if pred is None:\n                pred = np.squeeze(fold_model.predict(image))\n            else:\n                pred += np.squeeze(fold_model.predict(image))\n        \n        pred = pred/len(fold_models)\n        \n        pred = cv2.resize(pred, (WINDOW, WINDOW))\n        preds[x1:x2,y1:y2] = (pred > 0.5).astype(np.uint8)\n            \n    subm[i] = {'id':filename.stem, 'rle': rle_numba_encode(preds)}\n    del preds\n    gc.collect();","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:45:40.360788Z","iopub.execute_input":"2022-06-23T05:45:40.362960Z","iopub.status.idle":"2022-06-23T05:45:49.891019Z","shell.execute_reply.started":"2022-06-23T05:45:40.362920Z","shell.execute_reply":"2022-06-23T05:45:49.890138Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission = pd.DataFrame.from_dict(subm, orient='index')\nsubmission.to_csv('../working/submission.csv', index=False)\n\nsubmission.head()","metadata":{"execution":{"iopub.status.busy":"2022-06-23T05:45:49.892386Z","iopub.execute_input":"2022-06-23T05:45:49.892931Z","iopub.status.idle":"2022-06-23T05:45:50.023203Z","shell.execute_reply.started":"2022-06-23T05:45:49.892888Z","shell.execute_reply":"2022-06-23T05:45:50.022445Z"},"trusted":true},"execution_count":null,"outputs":[]}]}