{"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":"This notebook is to train <a href=https://www.kaggle.com/c/sartorius-cell-instance-segmentation>Sartorius Competition</a> as Semantic Segmentation (<span style=\"color: red; \">Not</span> Instance Segmentation)\n- Model : Unet (backboned: EfficientNet, package: segmetation_models-1.0.1)  \n- Image : png images -> ndarray -> tfrecord (split into tiles)\n- Annotation(target) : tran.csv(run-lenght str) -> adarray -> tfrecord (split into tiles)\n\n\nRefs.  \nhttps://www.kaggle.com/wrrosa/hubmap-tf-with-tpu-efficientunet-512x512-train  \nhttps://www.kaggle.com/ammarnassanalhajali/sartorius-segmentation-keras-u-net-training","metadata":{}},{"cell_type":"markdown","source":"<a class='anchor' id='TOC'></a>\n# Table of Contents\n\n1. [Packages](#1)\n1. [Accelerator](#2)\n1. [Parameters](#3)\n1. [Data](#4)\n1. [Modeling](#5)\n1. [Train](#6)\n1. [Evaluation](#7)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19"}},{"cell_type":"markdown","source":"<a class='anchor' id='1'></a>\n# 1. Packages\n[Back to Table of Contents](#TOC)","metadata":{}},{"cell_type":"code","source":"!pip install segmentation_models==1.0.1 -q","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:18.509076Z","iopub.execute_input":"2021-11-09T04:24:18.50965Z","iopub.status.idle":"2021-11-09T04:24:28.146818Z","shell.execute_reply.started":"2021-11-09T04:24:18.509608Z","shell.execute_reply":"2021-11-09T04:24:28.14583Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os, glob, gc, re, yaml, json\nfrom tqdm.notebook import tqdm\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport matplotlib.gridspec as gridspec\nimport cv2\nimport seaborn as sns\nimport tensorflow as tf\nfrom datetime import datetime\nfrom pprint import pprint\nfrom PIL import Image, ImageEnhance\n\nimport warnings\nwarnings.filterwarnings('ignore')","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:30.144215Z","iopub.execute_input":"2021-11-09T04:24:30.144517Z","iopub.status.idle":"2021-11-09T04:24:35.048557Z","shell.execute_reply.started":"2021-11-09T04:24:30.144474Z","shell.execute_reply":"2021-11-09T04:24:35.047801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<a class='anchor' id='2'></a>\n# 2. Accelarator\n[Back to Table of Contents](#TOC)</br>\nSelect processor (priority: TPU>GPU>CPU)","metadata":{}},{"cell_type":"code","source":"def set_strategy():\n    try:\n        tpu = tf.distribute.cluster_resolver.TPUClusterResolver.connect() # TPU detection\n    except ValueError:\n        tpu = None\n        gpus = tf.config.experimental.list_logical_devices(\"GPU\")\n\n    if tpu:\n        strategy = tf.distribute.TPUStrategy(tpu)\n        print('Running on TPU ', tpu.cluster_spec().as_dict()['worker'])\n    elif len(gpus) > 1:\n        strategy = tf.distribute.MirroredStrategy([gpu.name for gpu in gpus])\n        print('Running on multiple GPUs ', [gpu.name for gpu in gpus])\n    elif len(gpus) == 1:\n        strategy = tf.distribute.get_strategy() # default strategy that works on CPU and single GPU\n        print('Running on single GPU ', gpus[0].name)\n    else:\n        strategy = tf.distribute.get_strategy() # default strategy that works on CPU and single GPU\n        print('Running on CPU')\n\n    print(\"Number of accelerators: \", strategy.num_replicas_in_sync)\n\n    return strategy\n\nstrategy = set_strategy()","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:36.525534Z","iopub.execute_input":"2021-11-09T04:24:36.526325Z","iopub.status.idle":"2021-11-09T04:24:38.636663Z","shell.execute_reply.started":"2021-11-09T04:24:36.52626Z","shell.execute_reply":"2021-11-09T04:24:38.633466Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<a class='anchor' id='3'></a>\n# 3. Parameters\n[Back to Table of Contents](#TOC)","metadata":{}},{"cell_type":"code","source":"# Path\nINPUT_PATH = '../input/sartorius-cell-instance-segmentation/'\n\n# Paramaters\nP = {}\nP['DEBUG'] = False # If true, number of epochs and fold calculation are minimized.\nP['MODEL'] = 'Unet' # UNet, FPN, Linknet, PNPNet\nP['BASE_TILE'] = [128, 128] # Tile sile to split original image\nP['RESIZED_TILE'] = [128, 128] # Image compressison for quick calculation (If no compression, P['RESIZED_TILE'] = P['BASE_SIZE'])\nP['MIN_OVERLAP'] = 32 # Overlap width of each tile (Note: Edge image may overlap more than MIN_OVERLAP)\nP['BACKBONE'] = 'efficientnetb1' \nP['WEIGHT'] = '../input/efficientnetb0b7-keras-weights/efficientnet-b1_weights_tf_dim_ordering_tf_kernels_autoaugment_notop.h5' # Downloaded from https://github.com/Callidior/keras-applications/releases\nP['SEED'] = 2021 # Random seed\nP['VERBOSE'] = 1 # Show the detail of train history, or not\nP['BATCH_COEF'] = 64\nP['BATCH_SIZE'] = P['BATCH_COEF'] * strategy.num_replicas_in_sync\nP['LR'] = 5e-4 # Learing rate\nP['STEPS_COEF'] = 3 # step_per_epoch = P['STEPS_COE'] * #TILES // BATCH_SIZE (Nominal: 1)\nP['NFOLDS'] =6 # Number of folds\nif P['DEBUG']==True:\n    P['EPOCHS'] = 1 # Number of epochs\n    P['CALC_FOLDS'] = 1 # One fold calcuration is only performed\nelif P['DEBUG']==False:\n    P['EPOCHS'] = 5 # Number of epochs\n    P['CALC_FOLDS'] = P['NFOLDS'] # Full fold calcuration is performed\n\n\nAUTO = tf.data.experimental.AUTOTUNE","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:38.638365Z","iopub.execute_input":"2021-11-09T04:24:38.638636Z","iopub.status.idle":"2021-11-09T04:24:38.647102Z","shell.execute_reply.started":"2021-11-09T04:24:38.638601Z","shell.execute_reply":"2021-11-09T04:24:38.646225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4. Data\n[Back to Table of Contents](#TOC)","metadata":{}},{"cell_type":"markdown","source":"## Tabular data","metadata":{}},{"cell_type":"code","source":"df_train = pd.read_csv(INPUT_PATH + 'train.csv')\ndisplay(df_train)","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:39.39714Z","iopub.execute_input":"2021-11-09T04:24:39.397774Z","iopub.status.idle":"2021-11-09T04:24:39.933848Z","shell.execute_reply.started":"2021-11-09T04:24:39.397734Z","shell.execute_reply":"2021-11-09T04:24:39.933137Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Anotations are grouped by id\ndf_tmp = df_train.drop_duplicates('id').reset_index(drop=True).sort_values('id')\ndf_tmp[\"annotation\"] = df_train.groupby('id')['annotation'].agg(list).reset_index(drop=True)\ndf_train = df_tmp.copy()\ndisplay(df_train)","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:39.935425Z","iopub.execute_input":"2021-11-09T04:24:39.935835Z","iopub.status.idle":"2021-11-09T04:24:40.004979Z","shell.execute_reply.started":"2021-11-09T04:24:39.935797Z","shell.execute_reply":"2021-11-09T04:24:40.004194Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Feature distributions","metadata":{}},{"cell_type":"code","source":"def fig_layout(seaborn_plot):\n    plt.xlabel('')\n    plt.ylabel('')\n    plt.yticks([])\n    # Hide frame\n    for l in ['right', 'top', 'left']:\n        seaborn_plot.spines[l].set_visible(False)\n    # Count record length\n    record_length = 0\n    for rectangle in seaborn_plot.patches:\n        record_length += rectangle.get_height()\n    # Add annotation of ratio \n    for rectangle in seaborn_plot.patches:\n        height = rectangle.get_height()\n        width = rectangle.get_width()\n        ratio = round(height/record_length*100,1)\n        # Ratio\n        sns_plot.annotate(f'{ratio}%',\n                          xy=(rectangle.get_x()+width/2, height),\n                          ha='center', va='center', size=8,\n                          xytext=(0, 10), textcoords='offset points')\n\n\ndf_tmp = df_train.copy()\n# dtype change\n# # plate_time: 11h30m00s -> 11.5\nhour   = pd.to_datetime(df_train['plate_time'], format='%Hh%Mm%Ss').dt.hour\nminute = pd.to_datetime(df_train['plate_time'], format='%Hh%Mm%Ss').dt.minute\ndf_tmp['plate_time'] = round(hour + minute/60, 2)\n# # sample_date: str -> datetime\ndf_tmp['sample_date'] = pd.to_datetime(df_train['sample_date'], format='%Y-%m-%d').dt.date\n\n# Plot\nfig = plt.figure(figsize=(22, 5))\nfeatures = ['cell_type', 'plate_time', 'sample_date']\nfor i, feature in enumerate(features):\n    plt.subplot(1,3,i+1)\n    plt.title(feature, size=20)\n    sns_plot = sns.countplot(x=feature, data=df_tmp)\n    fig_layout(sns_plot)\n    fig.autofmt_xdate(rotation=90)\n    \nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:40.971544Z","iopub.execute_input":"2021-11-09T04:24:40.971798Z","iopub.status.idle":"2021-11-09T04:24:41.465544Z","shell.execute_reply.started":"2021-11-09T04:24:40.971771Z","shell.execute_reply":"2021-11-09T04:24:41.464851Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Image data","metadata":{}},{"cell_type":"code","source":"train_imgs = INPUT_PATH + 'train/' + df_train['id'] + '.png'\nprint(f'train_images: {len(train_imgs)} files')\ndisplay(train_imgs)","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:41.516778Z","iopub.execute_input":"2021-11-09T04:24:41.517086Z","iopub.status.idle":"2021-11-09T04:24:41.525975Z","shell.execute_reply.started":"2021-11-09T04:24:41.517055Z","shell.execute_reply":"2021-11-09T04:24:41.525161Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# 3 cell type images are displayed.\ncell_types = df_train['cell_type'].unique()\ngs = gridspec.GridSpec(1, 3)\nplt.figure(figsize = (25, 20))\nfor i, cell_type in enumerate(cell_types):\n    idx = df_train[df_train['cell_type'] == cell_type].index[0]\n    img_id = df_train['id'][idx]\n    img = cv2.imread(INPUT_PATH + 'train/' + img_id + '.png')\n    ax = plt.subplot(gs[i])\n    ax.set_title(f'id: {img_id},  cell_type: {cell_type}')\n    ax.imshow(img)\n    ax.set_aspect('equal')\n    plt.axis('on')   \n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:41.767207Z","iopub.execute_input":"2021-11-09T04:24:41.767486Z","iopub.status.idle":"2021-11-09T04:24:42.497962Z","shell.execute_reply.started":"2021-11-09T04:24:41.767457Z","shell.execute_reply":"2021-11-09T04:24:42.49726Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Decode annotations (string -> ndarray)","metadata":{}},{"cell_type":"code","source":"# Decode run-length string -> ndarry\ndef rle_decode(annotation, shape):\n    '''\n    annotation: string\n    shape: (height, width)\n    return: ndarray, mask: 1, background: 0\n    '''\n    rle = annotation.split() # Even elements are starts, odd elements are the lengths.\n    starts  = np.asarray(rle[0:][::2], dtype=int)\n    lengths = np.asarray(rle[1:][::2], dtype=int)\n    starts -= 1 # Run-length start is numbered from one, on the other hand, list is numbered from zero.\n    ends = starts + lengths\n    \n    mask = np.zeros(shape[0]*shape[1], dtype=np.uint8)\n    for start, end in zip(starts, ends):\n        mask[start:end] = 1\n\n    return mask.reshape(shape)\n\n# Build mask image from all annotations with same id\ndef build_masks(annotations, shape, distinguish_objects=False):\n    '''\n    annotation_list: List[string]\n    shape: (height, width)\n    return: ndarray, mask: integer 1,2,3,..., background:0\n    '''\n    masks = np.zeros(shape, dtype=np.uint8)\n    for i, annotation in enumerate(annotations):\n        mask = rle_decode(annotation, shape)\n        if distinguish_objects:\n            masks = np.where(mask==0, masks, i+1)\n        else:\n            masks = np.where(mask==0, masks, 1)\n    \n    return masks\n\ndef plot_image_and_mask(img, mask, title=None):\n    fig, ax = plt.subplots(1, 4, figsize=(20,4))\n    \n    ax[0].set_title('Original image')\n    ax[0].imshow(img)\n    \n    ax[1].set_title('High contrasted image')\n    img_hc = img.max() - img\n    img_hc = np.asarray(ImageEnhance.Contrast(Image.fromarray(img_hc)).enhance(24))\n    ax[1].imshow(img_hc)\n       \n    ax[2].set_title('Mask')\n    ax[2].imshow(mask, cmap='inferno')\n    \n    ax[3].set_title('Image + Mask')\n    mask_ = np.tile(np.expand_dims(mask, 2), 3) # shape: (height, width) -> (height, width, 3) \n    mask_ = np.clip(mask_,0,1)*255 # mask: (255,255,255), background: (0,0,0)\n    mask_[:,:,2] = 0 # mask: (255,255,0): yellow\n    mask_ = mask_.astype(np.uint8) # type: np.uint16 -> np.unit8\n    merge_img_mask = cv2.addWeighted(img_hc, 0.80, mask_, 0.20, gamma=0.0)\n    ax[3].imshow(merge_img_mask)\n    \n    fig.suptitle(title, fontsize=14)\n    ","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:42.499702Z","iopub.execute_input":"2021-11-09T04:24:42.499968Z","iopub.status.idle":"2021-11-09T04:24:42.668056Z","shell.execute_reply.started":"2021-11-09T04:24:42.499929Z","shell.execute_reply":"2021-11-09T04:24:42.667118Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for cell_type in cell_types:\n    idx = df_train[df_train['cell_type'] == cell_type].index[0]\n    img_id = df_train['id'][idx]\n    img = cv2.imread(INPUT_PATH + 'train/' + img_id + '.png')\n    mask = build_masks(annotations=df_train.annotation[idx],\n                       shape=(df_train.height[idx], df_train.width[idx]),\n                       distinguish_objects=True)\n    plot_image_and_mask(img, mask, title=f'id: {img_id},  cell_type: {cell_type}')","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:42.669951Z","iopub.execute_input":"2021-11-09T04:24:42.670507Z","iopub.status.idle":"2021-11-09T04:24:45.634972Z","shell.execute_reply.started":"2021-11-09T04:24:42.670455Z","shell.execute_reply":"2021-11-09T04:24:45.634199Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Generate TFRecord","metadata":{}},{"cell_type":"code","source":"# Utility functions\n\n# Cast datatypes into 1 of the type lists (integer,float and bytes)\ndef _bytes_feature(value):\n    \"\"\"Returns a bytes_list from a string / byte.\"\"\"\n    if isinstance(value, type(tf.constant(0))):\n        value = value.numpy() # BytesList won't unpack a string from an EagerTensor.\n    return tf.train.Feature(bytes_list=tf.train.BytesList(value=[value]))\n\ndef _float_feature(value):\n    \"\"\"Returns a float_list from a float / double.\"\"\"\n    return tf.train.Feature(float_list=tf.train.FloatList(value=[value]))\n\ndef _int64_feature(value):\n    \"\"\"Returns an int64_list from a bool / enum / int / uint.\"\"\"\n    return tf.train.Feature(int64_list=tf.train.Int64List(value=[value]))\n\n# Serialization\ndef serialize_example(feature0, feature1, feature2):\n    # Create a feature dictionary which will be the contents of message\n    feature = {'image': _bytes_feature(feature0),\n               'mask' : _bytes_feature(feature1),\n               'cell_type' : _bytes_feature(feature2)}\n    # Serialization: convert the features into to bytes\n    example_proto = tf.train.Example(features=tf.train.Features(feature=feature))\n    return example_proto.SerializeToString()\n\n# Count images in a tfrecord file\ndef count_data_items(filenames):\n    n = [int(re.compile(r\"-([0-9]*)\\.\").search(f).group(1)) for f in filenames]\n    return np.sum(n)","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:45.636708Z","iopub.execute_input":"2021-11-09T04:24:45.637358Z","iopub.status.idle":"2021-11-09T04:24:45.64847Z","shell.execute_reply.started":"2021-11-09T04:24:45.637317Z","shell.execute_reply":"2021-11-09T04:24:45.64784Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Config of tiles\nBASE = P['BASE_TILE'][0] # Base tile size (Not original image size)\nRESIZE = P['RESIZED_TILE'][0] # Re-sized tile size\nreduce = BASE//RESIZE # Reduce base image size\n\n# Path to save tfrecords\nTFREC_PATH = f'./tfrec-{len(df_train)}-data_{RESIZE}x{RESIZE}-tile/'\nP['DATASET'] = TFREC_PATH\nif not os.path.exists(TFREC_PATH):\n    os.mkdir(TFREC_PATH)\n    os.mkdir(TFREC_PATH + 'train')\n    os.mkdir(TFREC_PATH + 'test')\n\n# For statistics\nx_tot, x2_tot  = [], [] \n\nprint('Generating tfrcords...')\nfor idx in range(len(df_train)):\n    \n    image_id = df_train.loc[idx, 'id']\n    # Load image and decode mask\n    img = cv2.imread(INPUT_PATH + 'train/' + image_id + '.png')\n    img = cv2.cvtColor(img, cv2.COLOR_RGB2BGR)\n    mask = build_masks(df_train.loc[idx,'annotation'], (img.shape[0], img.shape[1])) \n    \n    # Padding to make the image dividable by tile size \n    pad0 = (reduce*RESIZE - img.shape[0]%(reduce*RESIZE))%(reduce*RESIZE)\n    pad1 = (reduce*RESIZE - img.shape[1]%(reduce*RESIZE))%(reduce*RESIZE)\n    img  = np.pad(img,\n                  [[pad0//2, pad0-pad0//2], [pad1//2, pad1-pad1//2], [0, 0]],\n                  constant_values=0)\n    mask = np.pad(mask,\n                  [[pad0//2, pad0-pad0//2], [pad1//2,pad1-pad1//2]],\n                  constant_values=0)\n\n    # Tiling image and mask using the reshape + transpose trick\n    # image\n    img = cv2.resize(img, (img.shape[1]//reduce, img.shape[0]//reduce),\n                     interpolation=cv2.INTER_AREA)\n    img = img.reshape(img.shape[0]//RESIZE, RESIZE, img.shape[1]//RESIZE, RESIZE, 3)\n    img = img.transpose(0, 2, 1, 3, 4).reshape(-1, RESIZE, RESIZE, 3)\n    # mask\n    mask = cv2.resize(mask, (mask.shape[1]//reduce, mask.shape[0]//reduce),\n                      interpolation=cv2.INTER_NEAREST)\n    mask = mask.reshape(mask.shape[0]//RESIZE, RESIZE, mask.shape[1]//RESIZE, RESIZE)\n    mask = mask.transpose(0, 2, 1, 3).reshape(-1, RESIZE, RESIZE)\n    # cell_type\n    cell_type = df_train.loc[idx, 'cell_type']\n    cell_type = [cell_type.encode('utf-8') for _ in range(len(mask))]  \n    \n    # Generate TFRecord\n    num_tiles = 0\n    filename = TFREC_PATH + f'train/{image_id}.tfrec'\n    with tf.io.TFRecordWriter(filename) as writer:\n\n        for i, (img_, mask_, cell_type_) in enumerate(zip(img, mask, cell_type)):\n                       \n            x_tot.append((img_/255.0).reshape(-1, 3).mean(0))\n            x2_tot.append(((img_/255.0)**2).reshape(-1, 3).mean(0))\n            img_ = cv2.cvtColor(img_, cv2.COLOR_RGB2BGR)\n            \n            # Create tf.example\n            example = serialize_example(img_.tobytes(), mask_.tobytes(), cell_type_)\n            writer.write(example)\n            num_tiles +=1\n    \n    os.rename(filename, TFREC_PATH + 'train/'+ image_id + '-' + str(num_tiles) + '.tfrec')\n    if (idx+1)%100==0:\n        print(f'{idx+1} tfrecords generated')\n\nprint(f'Total {idx+1} tfrecords generated')\nprint(f'Saved in {TFREC_PATH}\\n')\n\ntrain_tfrecs = glob.glob(TFREC_PATH + 'train/*.tfrec')\nprint(f'Number of TFRecord files: {len(train_tfrecs)}')\nprint(f'Number of total tiles: {count_data_items(train_tfrecs)}')\nprint('Statistics of pixel values:')\nprint('Mean:', np.array(x_tot).mean(0))\nprint('STD :', np.sqrt(np.array(x2_tot).mean(0) - np.array(x_tot).mean(0)**2))    ","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:24:45.649752Z","iopub.execute_input":"2021-11-09T04:24:45.650499Z","iopub.status.idle":"2021-11-09T04:25:44.168754Z","shell.execute_reply.started":"2021-11-09T04:24:45.650453Z","shell.execute_reply":"2021-11-09T04:25:44.167966Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Confirm the generated records\nHere, confirm 3 tfrecord samples for 3 cell types","metadata":{}},{"cell_type":"code","source":"from skimage.segmentation import mark_boundaries\n\nDIM = RESIZE\n# Dataloader\ndef _parse_image_function(example_proto):\n    image_feature_description = {\n        'image': tf.io.FixedLenFeature([], tf.string),\n        'mask': tf.io.FixedLenFeature([], tf.string),\n        'cell_type': tf.io.FixedLenFeature([], tf.string)\n    }\n    single_example = tf.io.parse_single_example(example_proto, image_feature_description)\n    # image\n    image = tf.reshape(tf.io.decode_raw(single_example['image'], out_type=np.dtype('uint8')), (DIM,DIM,3))\n    image = tf.image.resize(image, (DIM, DIM))/255.0\n    # mask\n    mask =  tf.reshape(tf.io.decode_raw(single_example['mask'], out_type='bool'),(DIM,DIM,1))\n    mask = tf.image.resize(tf.cast(mask,'uint8'),(DIM, DIM))\n    # cell_type\n    cell_type = single_example['cell_type']\n    return image, mask, cell_type\n\ndef load_dataset(filenames):\n    dataset = tf.data.TFRecordDataset(filenames)\n    dataset = dataset.map(lambda ex: _parse_image_function(ex))\n    return dataset\n\ndef get_dataset(filename, n):\n    dataset = load_dataset(filename)\n    dataset = dataset.batch(n)\n    return dataset\n\n# Utility functions\ndef get_cell_type_from_tfrec(tfrec):\n    num_tiles = count_data_items([tfrec])\n    for i, m, cell_types in get_dataset(tfrec, num_tiles).take(1):\n        break\n    return cell_types[0].numpy().decode() # string\n\ndef get_unique_cell_samples(tfrecs):\n    sample_idx = []\n    sample_cell_type = []\n    for idx in range(100):\n        cell_type = get_cell_type_from_tfrec(tfrecs[idx])\n        if cell_type not in sample_cell_type:\n            sample_idx.append(idx)\n            sample_cell_type.append(cell_type)\n        if len(sample_idx) == 3:\n            break\n    return sample_idx, sample_cell_type\n\n# Plot 3 samples\nsample_idx, sample_cell_type = get_unique_cell_samples(train_tfrecs)\nprint(f'Sample index: {sample_idx}')\nprint(f'Sample cell type: {sample_cell_type}\\n')\nfor idx in sample_idx:\n    num_tiles = count_data_items([train_tfrecs[idx]])\n    for imgs, masks, cell_types in get_dataset(train_tfrecs[idx], num_tiles).take(1):\n        break\n    print(f'Sample image: {train_tfrecs[idx].split(\"/\")[-1]}')\n    print(f'image shape: {imgs.shape}, mask shape: {masks.shape},  cell type: {cell_types[0].numpy().decode()}')\n    \n    gs = gridspec.GridSpec(5, num_tiles//5)\n    plt.figure(figsize = (12, 10))\n    for i in range(num_tiles):\n        ax1 = plt.subplot(gs[i])\n        ax1.set_xticklabels([])\n        ax1.set_yticklabels([])\n        ax1.set_aspect('equal')\n        ax1.set_axis_off()\n        ax1.imshow(mark_boundaries(imgs[i], masks[i].numpy().squeeze().astype('bool')))\n    plt.subplots_adjust(wspace=0.01, hspace=0.01)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:25:44.171117Z","iopub.execute_input":"2021-11-09T04:25:44.171386Z","iopub.status.idle":"2021-11-09T04:25:50.174405Z","shell.execute_reply.started":"2021-11-09T04:25:44.17135Z","shell.execute_reply":"2021-11-09T04:25:50.173619Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# ","metadata":{}},{"cell_type":"markdown","source":"<a class='anchor' id='5'></a>\n# 5. Modeling\n[Back to Table of Contents](#TOC)","metadata":{}},{"cell_type":"markdown","source":"### Build model","metadata":{}},{"cell_type":"code","source":"# Metric\ndef iou_coef(y_true, y_pred, smooth=1):\n    intersection = K.sum(K.abs(y_true * y_pred), axis=[1,2,3])\n    union = K.sum(y_true,[1,2,3]) + K.sum(y_pred,[1,2,3]) - intersection\n    iou = K.mean((intersection + smooth) / (union + smooth), axis=0)\n    return iou\n\n# Loss\ndef bce_dice_loss(y_true, y_pred):\n    \n    def dice_loss(y_true, y_pred):\n        smooth = 1.\n        y_true_f = K.flatten(y_true)\n        y_pred_f = K.flatten(y_pred)\n        intersection = y_true_f * y_pred_f\n        score = (2. * K.sum(intersection) + smooth) / (K.sum(y_true_f) + K.sum(y_pred_f) + smooth)\n        return 1. - score\n    \n    bce_loss_ = tf.keras.losses.binary_crossentropy(tf.cast(y_true, tf.float32), y_pred)\n    dice_loss_ = dice_loss(tf.cast(y_true, tf.float32), y_pred)\n    return bce_loss_ * 0.5 + dice_loss_*0.5\n \n\n# Build model\nimport segmentation_models as sm\nmodel = sm.Unet(P['BACKBONE'], encoder_weights = P['WEIGHT'])\nmodel.compile(optimizer = tf.keras.optimizers.Adam(lr = P['LR']),\n              loss = bce_dice_loss,#'focal_tversky',\n              metrics=[iou_coef])\n\nprint('Model:', P['MODEL'])\nprint('Model backbone', P['BACKBONE'])\nprint('Total params: ', model.count_params())\ntrainable_params = sum(np.prod(w.shape) for w in model.trainable_weights)\nprint('Trainable params:', trainable_params)\nnon_trainable_params = sum(np.prod(w.shape) for w in model.non_trainable_weights)\nprint('Non-trainable params:', non_trainable_params)","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:25:50.175867Z","iopub.execute_input":"2021-11-09T04:25:50.176137Z","iopub.status.idle":"2021-11-09T04:25:54.566954Z","shell.execute_reply.started":"2021-11-09T04:25:50.1761Z","shell.execute_reply":"2021-11-09T04:25:54.566196Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Data parser","metadata":{}},{"cell_type":"code","source":"DIM = RESIZE\ndef _parse_image_function_argument(example_proto, augment=True):\n    image_feature_description = {\n        'image': tf.io.FixedLenFeature([], tf.string),\n        'mask': tf.io.FixedLenFeature([], tf.string),\n    }\n    single_example = tf.io.parse_single_example(example_proto, image_feature_description)\n    image = tf.reshape(tf.io.decode_raw(single_example['image'], out_type=np.dtype('uint8')), (DIM,DIM,3))\n    mask =  tf.reshape(tf.io.decode_raw(single_example['mask'], out_type='bool'), (DIM,DIM,1))\n    \n    if augment: # https://www.kaggle.com/kool777/training-hubmap-eda-tf-keras-tpu\n\n        if tf.random.uniform(()) > 0.5:\n            image = tf.image.flip_left_right(image)\n            mask  = tf.image.flip_left_right(mask)\n\n        if tf.random.uniform(()) > 0.4:\n            image = tf.image.flip_up_down(image)\n            mask  = tf.image.flip_up_down(mask)\n\n        if tf.random.uniform(()) > 0.5:\n            image = tf.image.rot90(image, k=1)\n            mask  = tf.image.rot90(mask, k=1)\n\n        if tf.random.uniform(()) > 0.45:\n            image = tf.image.random_saturation(image, 0.7, 1.3)\n\n        if tf.random.uniform(()) > 0.45:\n            image = tf.image.random_contrast(image, 0.8, 1.2)\n    \n    return tf.cast(image, tf.float32), tf.cast(mask, tf.float32)\n\ndef load_train_dataset(filenames, ordered=False, augment=True):\n    ignore_order = tf.data.Options()\n    if not ordered:\n        ignore_order.experimental_deterministic = False\n    dataset = tf.data.TFRecordDataset(filenames, num_parallel_reads=AUTO)\n    dataset = dataset.with_options(ignore_order)\n    dataset = dataset.map(lambda ex: _parse_image_function_argument(ex, augment=augment), num_parallel_calls=AUTO)\n    return dataset\n\ndef get_train_dataset(train_filenames, ordered=True, augment=True, batch_size=P['BATCH_SIZE']):\n    dataset = load_train_dataset(train_filenames, ordered=ordered, augment=augment)\n    dataset = dataset.repeat()\n    dataset = dataset.shuffle(128, seed=P['SEED'])\n    dataset = dataset.batch(batch_size, drop_remainder=True)\n    dataset = dataset.prefetch(AUTO)\n    return dataset\n\ndef get_valid_dataset(valid_filenames, ordered=True, augment=False, batch_size=P['BATCH_SIZE']):\n    dataset = load_train_dataset(valid_filenames, ordered=ordered, augment=augment)\n    dataset = dataset.batch(batch_size, drop_remainder=True)\n    #dataset = dataset.cache()\n    dataset = dataset.prefetch(AUTO)\n    return dataset\n\n","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:25:54.568085Z","iopub.execute_input":"2021-11-09T04:25:54.568337Z","iopub.status.idle":"2021-11-09T04:25:54.586724Z","shell.execute_reply.started":"2021-11-09T04:25:54.568302Z","shell.execute_reply":"2021-11-09T04:25:54.585627Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<a class='anchor' id='6'></a>\n# 6. Train\n[Back to Table of Contents](#TOC)","metadata":{"execution":{"iopub.status.busy":"2021-11-04T00:54:41.011555Z","iopub.execute_input":"2021-11-04T00:54:41.013514Z","iopub.status.idle":"2021-11-04T00:54:41.0247Z","shell.execute_reply.started":"2021-11-04T00:54:41.013432Z","shell.execute_reply":"2021-11-04T00:54:41.023082Z"}}},{"cell_type":"code","source":"from sklearn.model_selection import KFold\nfrom tensorflow.keras import backend as K\n\n# Saved model path\nMODEL_PATH = './model/'\nif not os.path.exists(MODEL_PATH):\n    os.mkdir(MODEL_PATH)\n\n# Saved metrics\nmetrics = ['loss', 'iou_coef','accuracy']\nM = {}\nfor m in metrics:\n    M['train_'+ m] = []\n    M['valid_'+ m] = []\n\n# Train\nkfold = KFold(n_splits=P['NFOLDS'], shuffle=True, random_state=P['SEED'])\nfor fold, (train_idx, valid_idx) in enumerate(kfold.split(df_train)):\n    \n    print('#'*35); print('############ FOLD ',fold+1,' #############'); print('#'*35);\n    print(f'Tile Size: {DIM}, Batch Size: {P[\"BATCH_SIZE\"]}')\n    \n    # Split into train and validation\n    train_split = [train_tfrecs[i] for i in train_idx]\n    valid_split = [train_tfrecs[i] for i in valid_idx]\n    STEPS_PER_EPOCH = P['STEPS_COEF'] * count_data_items(train_split) // P['BATCH_SIZE']\n    \n    # Build model\n    K.clear_session()\n    with strategy.scope():   \n        model = sm.Unet(P['BACKBONE'], encoder_weights = P['WEIGHT'])\n        model.compile(optimizer = tf.keras.optimizers.Adam(lr = P['LR']),\n                      loss = bce_dice_loss,\n                      metrics = [iou_coef,'accuracy'])\n    \n    # Callbacks\n    checkpoint = tf.keras.callbacks.ModelCheckpoint(\n        filepath = MODEL_PATH + f'/model-fold{fold}',\n        verbose = P['VERBOSE'],\n        monitor ='val_loss',\n        patience = 10,\n        mode='max',\n        save_best_only=True\n    )\n    early_stop = tf.keras.callbacks.EarlyStopping(monitor='val_iou_coef',\n                                                  mode='max',\n                                                  patience=10,\n                                                  restore_best_weights=True)\n    reduce = tf.keras.callbacks.ReduceLROnPlateau(monitor='val_loss',\n                                                  factor=0.1,\n                                                  patience=8,\n                                                  min_lr=0.00001)\n    # Fit\n    print(f'Training Model Fold {fold+1}...')\n    history = model.fit(\n        get_train_dataset(train_split),\n        validation_data = get_valid_dataset(valid_split),\n        epochs = P['EPOCHS'],\n        steps_per_epoch = STEPS_PER_EPOCH,\n        callbacks = [checkpoint, reduce, early_stop],\n        verbose = P['VERBOSE'],\n    )   \n    \n    # Load best model\n    with strategy.scope():\n        model = tf.keras.models.load_model(MODEL_PATH + f'model-fold{fold}',\n                                           custom_objects = {'iou_coef'     : iou_coef,\n                                                             'bce_dice_loss': bce_dice_loss})\n\n    # Save metrics (loss, iou, accuaracy)\n    train_metric = model.evaluate(get_valid_dataset(train_split), return_dict=True)\n    valid_metric = model.evaluate(get_valid_dataset(valid_split), return_dict=True)\n    for m in metrics:\n        M['train_'+m].append(train_metric[m])\n        M['valid_'+m].append(valid_metric[m])\n    \n    \n    # Plot train result\n    plt.figure(figsize=(15,5))\n    plt.xlabel('Epoch', size=14)\n    plt.ylabel('IoU', size=14)\n    ## Epoch-IoU\n    epochs = np.arange(len(history.history['iou_coef']))\n    plt.plot(epochs, history.history['iou_coef'],    '-o', label='Train IoU', color='black')\n    plt.plot(epochs, history.history['val_iou_coef'],'-o', label='Valid IoU', color='red')\n    ## Max IoU point\n    iou_max = np.max(history.history['val_iou_coef'] )\n    x_iou_max = np.argmax(history.history['val_iou_coef'] )\n    plt.scatter(x_iou_max, iou_max, s=200, color='red')\n    ## Max IoU text\n    xdist = plt.xlim()[1] - plt.xlim()[0]\n    ydist = plt.ylim()[1] - plt.ylim()[0]\n    plt.text(x_iou_max-0.03*xdist, iou_max-0.13*ydist, 'Max IoU\\n%.2f'%iou_max, size=12)\n    ## IoU Legend\n    plt.legend(loc=2)\n    ## Epoch-Loss\n    plt2 = plt.gca().twinx() # Secondary axis\n    plt.ylabel('Loss', size=14)\n    plt2.plot(epochs, history.history['loss'],     '-o', label='Train Loss', color='black', linestyle=\"dashed\")\n    plt2.plot(epochs, history.history['val_loss'], '-o', label='Valid Loss', color='red',   linestyle=\"dashed\")\n    # Loss legend\n    plt.legend(loc=3)\n    plt.show()\n    \n    # DEBUG\n    if (fold+1) == P['CALC_FOLDS']:\n        break\n\nprint('Metrics mean over folds')\nfor m in metrics:\n    M['mean_train_'+m] = np.mean(M['train_'+m])\n    M['mean_valid_'+m] = np.mean(M['valid_'+m])\n    print('Train '+ m + ': '+ str(round(M['mean_train_'+m], 4)))\n    print('Valid '+ m + ': '+ str(round(M['mean_valid_'+m], 4)))","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:25:54.589462Z","iopub.execute_input":"2021-11-09T04:25:54.589813Z","iopub.status.idle":"2021-11-09T04:27:24.160345Z","shell.execute_reply.started":"2021-11-09T04:25:54.589776Z","shell.execute_reply":"2021-11-09T04:27:24.159152Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Save parameters\nwith open(MODEL_PATH + 'params.yaml', 'w') as file:\n    yaml.dump(P, file)\n\n# Save metrics\nwith open(MODEL_PATH + 'metrics.json', 'w') as outfile:\n    json.dump(M, outfile)\n\nprint(f'Model paramaters and metrics saved in {MODEL_PATH}/params.yaml, metrics.json')","metadata":{"execution":{"iopub.status.busy":"2021-11-09T04:27:24.161484Z","iopub.status.idle":"2021-11-09T04:27:24.162207Z","shell.execute_reply.started":"2021-11-09T04:27:24.161957Z","shell.execute_reply":"2021-11-09T04:27:24.161988Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<a class='anchor' id='7'></a>\n# 7. Evaluation\n[Back to Table of Contents](#TOC)\nHere, last fold model and validataion set are used.","metadata":{}},{"cell_type":"code","source":"# Load model\nwith strategy.scope():\n    model = tf.keras.models.load_model(MODEL_PATH + f'model-fold{fold}', # Use last fold model\n                                       custom_objects = {'iou_coef'     : iou_coef,\n                                                         'bce_dice_loss': bce_dice_loss})\n    print(f'model-fold{fold} is loaded')","metadata":{"execution":{"iopub.status.busy":"2021-11-09T02:35:43.409019Z","iopub.execute_input":"2021-11-09T02:35:43.409577Z","iopub.status.idle":"2021-11-09T02:35:59.789821Z","shell.execute_reply.started":"2021-11-09T02:35:43.409539Z","shell.execute_reply":"2021-11-09T02:35:59.788937Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def predict(tfrecord, model, threshold=0.5):\n    \n    # Read tfrecord\n    num_tiles = count_data_items([tfrecord])\n    dataset = get_valid_dataset(tfrecord, batch_size=num_tiles)\n    for imgs, masks in dataset.take(1):\n        break\n    \n    imgs  = imgs/255\n    masks = masks.numpy().squeeze().astype('bool')\n    preds = model.predict_generator(dataset, verbose=1) # 0 to 1 values\n    preds = np.where(preds<threshold, False, True).squeeze() # binarize at threshold\n    \n    return imgs, masks, preds\n\n# Plot true mask and prediction\nsample_idx, sample_cell_type = get_unique_cell_samples(valid_split) # Use last fold validation set\n\nfor idx in sample_idx:\n    valid_imgs, valid_masks, valid_preds = predict(valid_split[idx], model) \n    image_id = valid_split[idx].split('/')[-1].split('-')[0]\n    cell_type = df_train[df_train['id']==image_id]['cell_type'].values[0]\n    num_of_objects = len(df_train[df_train['id']==image_id]['annotation'].values[0])\n    print(f'ID: {image_id}, Cell type: {cell_type}, Number of objects: {num_of_objects}')\n    \n    # Show true mask and prediction\n    gs = gridspec.GridSpec(5, 2*num_tiles//5) # 5row x 12col\n    fig = plt.figure(figsize = (24, 10.5))\n    for i in range(num_tiles):\n        row = i//6\n        col = i%6\n        \n        # true mask (left 5row x 6col)\n        ax1 = plt.subplot(gs[row, col])\n        ax1.set_xticklabels([])\n        ax1.set_yticklabels([])\n        ax1.set_aspect('equal')\n        ax1.set_axis_off()\n        ax1.imshow(mark_boundaries(valid_imgs[i], valid_masks[i]))\n        # prediction (right 5row x 6col)\n        ax2 = plt.subplot(gs[row, col+6])\n        ax2.set_xticklabels([])\n        ax2.set_yticklabels([])\n        ax2.set_aspect('equal')\n        ax2.set_axis_off()\n        ax2.imshow(mark_boundaries(valid_imgs[i], valid_preds[i]))\n\n    fig.suptitle('Ground Truth' + ' '*90 + 'Prediction',  fontsize=25)\n    plt.subplots_adjust(wspace=0.01, hspace=0.01)\n    plt.tight_layout()\n    plt.show()\n","metadata":{"execution":{"iopub.status.busy":"2021-11-09T02:38:23.231829Z","iopub.execute_input":"2021-11-09T02:38:23.232098Z","iopub.status.idle":"2021-11-09T02:38:36.48589Z","shell.execute_reply.started":"2021-11-09T02:38:23.232068Z","shell.execute_reply":"2021-11-09T02:38:36.48516Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Future work\nFor instance segmentation, cell types of \"shsy5y\" and \"cort\" are likely to separate into object, for example using cv2.connectedComponentsWithStats, but \"astro\" may not be able to. \n*Please upvoke, if useful for you.*","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}