{"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":"<h1 style=\"font-family: Verdana; font-size: 28px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\"><center><br>[WSI Preprocessing] 👀: Fast Tiling + Tissue Segmentation </center></h1>\n                                                      \n<center><img src = \"https://drive.google.com/uc?id=1pbIvjTlhGywfhiMTqcsdOB5LSHlklM90\"/></center>   \n\n<h5 style=\"text-align: center; font-family: Verdana; font-size: 12px; font-style: normal; font-weight: bold; text-decoration: None; text-transform: none; letter-spacing: 1px; color: black; background-color: #ffffff;\">CREATED BY: NGHI HUYNH</h5>","metadata":{}},{"cell_type":"markdown","source":"<p id=\"toc\"></p>\n<h2 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" role=\"tab\" aria-controls=\"home\"><center><br>CONTENTS</center></h2>\n\n<h3 style=\"text-indent: 10vw; font-family: Verdana; font-size: 16px; font-style: normal; font-weight: normal; text-decoration: none; text-transform: none; letter-spacing: 2px; color: black; background-color: #ffffff;\"><a href=\"#preprocessing\">0&nbsp;&nbsp;&nbsp;&nbsp;WHY DATA PREPROCESSING?</a></h3>\n\n---\n\n<h3 style=\"text-indent: 10vw; font-family: Verdana; font-size: 16px; font-style: normal; font-weight: normal; text-decoration: none; text-transform: none; letter-spacing: 2px; color: black; background-color: #ffffff;\"><a href=\"#scaling\">1&nbsp;&nbsp;&nbsp;&nbsp;IMAGE SCALING</a></h3>\n\n---\n\n<h3 style=\"text-indent: 10vw; font-family: Verdana; font-size: 16px; font-style: normal; font-weight: normal; text-decoration: none; text-transform: none; letter-spacing: 2px; color: black; background-color: #ffffff;\"><a href=\"#segmentation\">2&nbsp;&nbsp;&nbsp;&nbsp;TISSUE SEGMENTATION</a></h3>\n\n---\n\n<h3 style=\"text-indent: 10vw; font-family: Verdana; font-size: 16px; font-style: normal; font-weight: normal; text-decoration: none; text-transform: none; letter-spacing: 2px; color: black; background-color: #ffffff;\"><a href=\"#retrieval\">3&nbsp;&nbsp;&nbsp;&nbsp;FAST TILING + TILE SELECTION + RETRIEVAL</a></h3>\n\n---\n\n<h3 style=\"text-indent: 10vw; font-family: Verdana; font-size: 16px; font-style: normal; font-weight: normal; text-decoration: none; text-transform: none; letter-spacing: 2px; color: black; background-color: #ffffff;\"><a href=\"#conclusion\">4&nbsp;&nbsp;&nbsp;&nbsp;CONCLUSION</a></h3>\n\n---\n\n","metadata":{}},{"cell_type":"markdown","source":"<div class=\"list-group\" id=\"list-tab\" role=\"tablist\">\n<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style='background:#F08080; border:0; color:black' role=\"tab\" aria-controls=\"home\"><center><br>If you find this notebook useful, do give me an upvote, it motivates me a lot.<br><br> This notebook is still a work in progress. Keep checking for further developments!😊</center></h3>","metadata":{}},{"cell_type":"markdown","source":"<a id=\"preprocessing\"></a>\n\n<h2 style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" id=\"preprocessing\"><left><br>&nbsp0. WHY DATA PREPROCESSING? <a href=\"#toc\">&#10514;</a><br></left> </h2>\n\n## Problems:\n* **WSIs with high resolution**: 3000x3000 for training data (HPA images) -> unable to fit into any SOTA/contemporary deep learning models in one go\n* **Small training data**: 351 images -> easily overfit\n* **Large portion of background, shadows, artefacts, etc**: -> add little to no information to the learning process, waste computation resources, and can cause severe detrimental effects on deep learning models\n    \n=> **Data preprocessing** is a crucial process in addressing these problems. Data preprocessing can increase the quality of the image data and facilitate our model training.  In brief, we will scale down whole-slide images (WSIs), apply filters to these scaled-down images for tissue segmentation, break the slides into tiles, select tiles containing tissue regions, and then retrieve and save those tiles into folders based on organs.\n\n## Data Preprocessing:\n### Steps:\n\n**1. Image Scaling**: scale down whole-slide images (WSIs) in the training data from 3000x3000 to 1024x1024 resolution\n\n**2. Tissue Segmentation**: identify tissue regions in scaled-down images using the following foreground/background filtering approaches:\n* [Simple Thresholding](https://docs.opencv.org/4.x/d7/d4d/tutorial_py_thresholding.html)\n* [Otsu's Binarization](https://docs.opencv.org/4.x/d7/d4d/tutorial_py_thresholding.html)\n* [Triangle Binarization](https://subscription.packtpub.com/book/data/9781789344912/9/ch09lvl1sec80/the-triangle-binarization-algorithm)\n    \n**3. Fast Tiling + Tile Selection + Retrieval**: select tiles based on thresholds calculated from ***step 2***. Note that tile size should be large enough that feature relevant to the task are visible to the model to learn\n\n[**References**: Apply filters for tissue segmentation](https://developer.ibm.com/articles/an-automatic-method-to-identify-tissues-from-big-whole-slide-images-pt2/)","metadata":{}},{"cell_type":"markdown","source":"<a id=\"scaling\"></a>\n\n<h2 style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" id=\"scaling\"><left><br>&nbsp1. IMAGE SCALING <a href=\"#toc\">&#10514;</a><br></left> </h2>","metadata":{}},{"cell_type":"markdown","source":"## Imports","metadata":{}},{"cell_type":"code","source":"import gc\nimport os\nimport cv2\nimport zipfile\nimport rasterio\nimport numpy as np\nimport math\nimport pandas as pd\nfrom PIL import Image\nimport tifffile as tiff\nimport seaborn as sns\nfrom tqdm.notebook import tqdm\nimport matplotlib.pyplot as plt\nfrom rasterio.windows import Window\nfrom torch.utils.data import Dataset","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:08.322255Z","iopub.execute_input":"2022-08-08T20:46:08.322872Z","iopub.status.idle":"2022-08-08T20:46:10.418386Z","shell.execute_reply.started":"2022-08-08T20:46:08.322753Z","shell.execute_reply":"2022-08-08T20:46:10.416797Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Config","metadata":{}},{"cell_type":"code","source":"TRAIN_PATH = '../input/hubmap-organ-segmentation/train_images/'\ntrain_df   = pd.read_csv('../input/hubmap-organ-segmentation/train.csv')\nfull_img_size = 3000","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:10.420820Z","iopub.execute_input":"2022-08-08T20:46:10.421480Z","iopub.status.idle":"2022-08-08T20:46:10.620499Z","shell.execute_reply.started":"2022-08-08T20:46:10.421439Z","shell.execute_reply":"2022-08-08T20:46:10.618685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Helper functions","metadata":{}},{"cell_type":"code","source":"# functions to convert encoding to mask and mask to encoding\n# https://www.kaggle.com/paulorzp/rle-functions-run-length-encode-decode\ndef mask2rle(img): # encoder\n    '''\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 \ndef rle2mask(mask_rle, shape=(1600,256)): # decoder\n    '''\n    mask_rle: run-length as string formated (start length)\n    shape: (width,height) of array to return \n    Returns numpy array, 1 - mask, 0 - background\n\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).T","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:10.622660Z","iopub.execute_input":"2022-08-08T20:46:10.624372Z","iopub.status.idle":"2022-08-08T20:46:10.638698Z","shell.execute_reply.started":"2022-08-08T20:46:10.624307Z","shell.execute_reply":"2022-08-08T20:46:10.637508Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# rescale to a desire img size\n# keep in mide that the size needs to be large enough\n# to keep all important features for model training\ndef rescale(img, mask, size=(1024,1024)):\n    scaled_img = cv2.resize(img, size)\n    scaled_mask = cv2.resize(mask, size, interpolation=cv2.INTER_NEAREST)\n    return scaled_img, scaled_mask","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:10.641036Z","iopub.execute_input":"2022-08-08T20:46:10.641563Z","iopub.status.idle":"2022-08-08T20:46:10.650307Z","shell.execute_reply.started":"2022-08-08T20:46:10.641526Z","shell.execute_reply":"2022-08-08T20:46:10.649184Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# visualize scaled-down images + masks based on organ\ndef visualize(df, organ):\n    tmp_df = df.loc[train_df['organ'] == organ].reset_index(drop=True)\n    plt.figure(figsize=(16,4))\n    for i in range(5):\n        img  = tiff.imread(TRAIN_PATH + str(tmp_df['id'][i]) +'.tiff')\n        mask = rle2mask(tmp_df['rle'][i], (full_img_size, full_img_size))\n        scaled_img, scaled_mask = rescale(img, mask)\n        plt.subplot(1,5,i+1)\n        plt.imshow(scaled_img)\n        plt.imshow(scaled_mask, cmap='seismic', alpha=0.4)\n        plt.title('ID:'+ str(tmp_df['id'][i]))\n    plt.suptitle(organ, fontsize=20)\n    plt.tight_layout()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:10.651523Z","iopub.execute_input":"2022-08-08T20:46:10.652534Z","iopub.status.idle":"2022-08-08T20:46:10.664585Z","shell.execute_reply.started":"2022-08-08T20:46:10.652479Z","shell.execute_reply":"2022-08-08T20:46:10.663487Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Scaling + Visualization","metadata":{}},{"cell_type":"code","source":"organs = train_df['organ'].unique()\nfor organ in organs:\n    visualize(train_df, organ)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:10.666216Z","iopub.execute_input":"2022-08-08T20:46:10.666917Z","iopub.status.idle":"2022-08-08T20:46:22.714031Z","shell.execute_reply.started":"2022-08-08T20:46:10.666881Z","shell.execute_reply":"2022-08-08T20:46:22.712545Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### *ID 12784 from kidney has no FTUs -> we need to thoroughly scan through all images to remove those faulted images if they exist*","metadata":{}},{"cell_type":"markdown","source":"<a id=\"segmentation\"></a>\n\n<h2 style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" id=\"segmentation\"><left><br>&nbsp2. TISSUE SEGMENTATION <a href=\"#toc\">&#10514;</a><br></left> </h2>\n\n## Goal: \nMask out non-tissue by setting non-tissue pixels to 0 for their red, green, and blue channels.\n\nConceptually and mathematically, it is often useful to have background values close to or equal to 0 (complement image)\n\nComplement image: simply subtract each pixel value from the maximum pixel value supported by the class (class uint8, the max value of pixel can be 255) and store in the output image array. In the output image, dark areas become lighter and light areas become darker\n\n## Foreground/Background Thresholding Techniques:\n* [Otsu](https://docs.opencv.org/4.x/d7/d4d/tutorial_py_thresholding.html): global image thresholding algorithm, applied to images which are bimodal, which means the image is having 2 peaks in the histogram\n* [Triangle binarization](https://subscription.packtpub.com/book/data/9781789344912/9/ch09lvl1sec80/the-triangle-binarization-algorithm): a line is drawn from the highest bar to the end of the histogram. Then for choosing the optimal threshold, distance of each bar is calculated from the line, whichever is least becomes the threshold value\n\n\n**Side Note**: [Understanding Histograms in Image Processing](https://towardsdatascience.com/histograms-in-image-processing-with-skimage-python-be5938962935)","metadata":{}},{"cell_type":"markdown","source":"## Helper Functions","metadata":{}},{"cell_type":"code","source":"def thresholding(img, method='otsu'):\n    # convert to grayscale complement image\n    grayscale_img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)\n    img_c = 255 - grayscale_img\n    thres, thres_img = 0, img_c.copy()\n    if method == 'otsu':\n        thres, thres_img = cv2.threshold(img_c, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)\n    elif method == 'triangle':\n        thres, thres_img = cv2.threshold(img_c, 0, 255, cv2.THRESH_BINARY+cv2.THRESH_TRIANGLE)\n    return thres, thres_img, img_c","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:22.715859Z","iopub.execute_input":"2022-08-08T20:46:22.716294Z","iopub.status.idle":"2022-08-08T20:46:22.725226Z","shell.execute_reply.started":"2022-08-08T20:46:22.716255Z","shell.execute_reply":"2022-08-08T20:46:22.723934Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def histogram(img, thres_img, img_c, thres):\n    \"\"\"\n    style: ['color', 'grayscale']\n    \"\"\" \n    plt.figure(figsize=(15,15))\n    \n    plt.subplot(3,2,1)\n    plt.imshow(img)\n    plt.title('Original Image')\n    \n    plt.subplot(3,2,2)\n    sns.histplot(img.ravel(), bins=np.arange(0,256), color='orange', alpha=0.5)\n    sns.histplot(img[:,:,0].ravel(), bins=np.arange(0,256), color='red', alpha=0.5)\n    sns.histplot(img[:,:,1].ravel(), bins=np.arange(0,256), color='Green', alpha=0.5)\n    sns.histplot(img[:,:,2].ravel(), bins=np.arange(0,256), color='Blue', alpha=0.5)\n    plt.legend(['Total', 'Red_Channel', 'Green_Channel', 'Blue_Channel'])\n    plt.ylim(0,0.3e6)\n    plt.xlabel('Intensity value')\n    plt.title('Color Histogram')\n    \n    plt.subplot(3,2,3)\n    plt.imshow(img_c, cmap='gist_gray')\n    plt.title('Complement Grayscale Image')\n    \n    plt.subplot(3,2,4)\n    sns.histplot(img_c.ravel(), bins=np.arange(0,256))\n    plt.axvline(thres, c='red', linestyle=\"--\")\n    plt.ylim(0,0.3e6)\n    plt.xlabel('Intensity value')\n    plt.title('Grayscale Complement Histogram')\n    \n    plt.subplot(3,2,5)\n    plt.imshow(thres_img, cmap='gist_gray')\n    plt.title('Thresholded Image')\n    \n    plt.subplot(3,2,6)\n    sns.histplot(thres_img.ravel(), bins=np.arange(0,256))\n    plt.axvline(thres, c='red', linestyle=\"--\")\n    plt.ylim(0,0.3e6)\n    plt.xlabel('Intensity value')\n    plt.title('Thresholded Histogram')\n    \n    plt.tight_layout()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:22.727417Z","iopub.execute_input":"2022-08-08T20:46:22.727969Z","iopub.status.idle":"2022-08-08T20:46:22.749007Z","shell.execute_reply.started":"2022-08-08T20:46:22.727929Z","shell.execute_reply":"2022-08-08T20:46:22.746874Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Otsu's Binarization","metadata":{}},{"cell_type":"code","source":"img_test2 = tiff.imread('../input/hubmap-organ-segmentation/train_images/10666.tiff')","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:22.751036Z","iopub.execute_input":"2022-08-08T20:46:22.751894Z","iopub.status.idle":"2022-08-08T20:46:22.801410Z","shell.execute_reply.started":"2022-08-08T20:46:22.751847Z","shell.execute_reply":"2022-08-08T20:46:22.799961Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"thres_otsu, thres_img, img_c = thresholding(img_test2, method='otsu')","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:22.807665Z","iopub.execute_input":"2022-08-08T20:46:22.809903Z","iopub.status.idle":"2022-08-08T20:46:22.880402Z","shell.execute_reply.started":"2022-08-08T20:46:22.809850Z","shell.execute_reply":"2022-08-08T20:46:22.878368Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"histogram(img_test2, thres_img, img_c, thres_otsu)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:46:22.882150Z","iopub.execute_input":"2022-08-08T20:46:22.882665Z","iopub.status.idle":"2022-08-08T20:47:24.776920Z","shell.execute_reply.started":"2022-08-08T20:46:22.882613Z","shell.execute_reply":"2022-08-08T20:47:24.775440Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Triangle Binarization","metadata":{}},{"cell_type":"code","source":"thres, thres_img, img_c = thresholding(img_test2, method='triangle')","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:47:24.778771Z","iopub.execute_input":"2022-08-08T20:47:24.779832Z","iopub.status.idle":"2022-08-08T20:47:24.809129Z","shell.execute_reply.started":"2022-08-08T20:47:24.779786Z","shell.execute_reply":"2022-08-08T20:47:24.807841Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"histogram(img_test2, thres_img, img_c, thres)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:47:24.811313Z","iopub.execute_input":"2022-08-08T20:47:24.812080Z","iopub.status.idle":"2022-08-08T20:48:22.727311Z","shell.execute_reply.started":"2022-08-08T20:47:24.812038Z","shell.execute_reply":"2022-08-08T20:48:22.724420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### *Otsu vs. Triangle binarization can only separate background and foreground but cannot detect artefacts. For example, Otsu treats lighter-colored tissue regions as background, whereas triangle binarization treats artefact as tissue regions. Thus, we need to dive deeper into these two thresholding techniques to have better tissue segmentation*","metadata":{}},{"cell_type":"markdown","source":"<a id=\"retrieval\"></a>\n\n<h2 style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" id=\"retrieval\"><left><br>&nbsp3. FAST TILING + TILE SELECTION + RETRIEVAL <a href=\"#toc\">&#10514;</a><br></left> </h2>","metadata":{}},{"cell_type":"markdown","source":"## Helper Functions","metadata":{}},{"cell_type":"code","source":"# adapted from: https://www.kaggle.com/code/analokamus/a-fast-tile-generation\ndef make_tiles(img, mask, tile_size=256):\n    '''\n    img: np.ndarray with dtype np.uint8 and shape (width, height, channel)\n    mask: np.ndarray with dtype np.uint9 and shape (width, height)\n    '''\n    w_i, h_i, ch = img.shape\n    w_m, h_m     = mask.shape\n    \n    pad0, pad1 = (tile_size - w_i%tile_size) % tile_size, (tile_size - h_i%tile_size) % tile_size\n    \n    padding_i = [[pad0//2, pad0-pad0//2], [pad1//2, pad1-pad1//2], [0, 0]]\n    padding_m = [[pad0//2, pad0-pad0//2], [pad1//2, pad1-pad1//2]]\n    \n    img = np.pad(img, padding_i, mode='constant', constant_values=255)\n    img = img.reshape(img.shape[0]//tile_size, tile_size, img.shape[1]//tile_size, tile_size, ch)\n    img = img.transpose(0, 2, 1, 3, 4).reshape(-1, tile_size, tile_size, ch)\n    \n    mask = np.pad(mask, padding_m, mode='constant', constant_values=255)\n    mask = mask.reshape(mask.shape[0]//tile_size, tile_size, mask.shape[1]//tile_size, tile_size)\n    mask = mask.transpose(0, 2, 1, 3).reshape(-1, tile_size, tile_size)\n    \n    num_tiles = len(mask)\n    #     if len(img) < num_tiles: # pad images so that the output shape be the same\n    #         padding = [[0, num_tiles-len(img)], [0, 0], [0, 0], [0, 0]]\n    #         img = np.pad(img, padding, mode='constant', constant_values=255)\n    #idxs = np.argsort(img.reshape(img.shape[0], -1).sum(-1))[:num_tiles] # pick up Top N dark tiles\n    #img = img[idxs]\n    return img, mask","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.730952Z","iopub.execute_input":"2022-08-08T20:48:22.732562Z","iopub.status.idle":"2022-08-08T20:48:22.753693Z","shell.execute_reply.started":"2022-08-08T20:48:22.732430Z","shell.execute_reply":"2022-08-08T20:48:22.751677Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def visualize_tiles(img, mask, img_tiles, mask_tiles):\n    plt.figure(figsize=(9.5,9.5))\n    plt.imshow(img)\n    plt.imshow(mask, cmap='seismic', alpha=0.4)\n    plt.title(f'Scaled Image + Mask\\nImage Size: {img.shape}')\n    plt.xticks([])\n    plt.yticks([])\n    plt.show()\n    \n    for i, imgs in enumerate(img_tiles):\n        plt.figure(figsize=(7.65,7.70))\n        num_tiles, size, _, _ = imgs.shape\n        rows = cols = int(math.sqrt(num_tiles))\n        for j, img_crop in tqdm(enumerate(imgs)):\n            plt.subplot(rows,cols,j+1)\n            plt.imshow(img_crop)\n            plt.imshow(mask_tiles[i][j], cmap='seismic',alpha=0.4)\n            plt.xticks([])\n            plt.yticks([])\n            plt.tight_layout()\n        plt.suptitle(f'Image + Mask \\nNum Tiles: {len(imgs)}\\nTile Size: {img_crop.shape}')\n        plt.tight_layout()\n        plt.show()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.755960Z","iopub.execute_input":"2022-08-08T20:48:22.756452Z","iopub.status.idle":"2022-08-08T20:48:22.774815Z","shell.execute_reply.started":"2022-08-08T20:48:22.756406Z","shell.execute_reply":"2022-08-08T20:48:22.773485Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# visualize selected tiles\ndef tile_selection(img, mask, thres, tile_size=256):\n    img_tiles, mask_tiles = make_tiles(img, mask, tile_size)\n    num_tiles, size, _, _ = img_tiles.shape\n    rows = cols = int(math.sqrt(num_tiles))\n    plt.figure(figsize=(7.65,7.70))\n    \n    for i, img_crop in tqdm(enumerate(img_tiles)):\n            plt.subplot(rows,cols,i+1)\n            plt.imshow(img_crop)\n            plt.imshow(mask_tiles[i], cmap='seismic',alpha=0.4)\n            plt.xticks([])\n            plt.yticks([])\n            plt.tight_layout()\n    plt.suptitle(f'Image + Mask \\nNum Tiles: {len(img_tiles)}\\nTile Size: {img_tiles[0].shape}')\n    plt.tight_layout()\n    plt.show()\n    \n    plt.figure(figsize=(7.65,7.70))\n    selected_tiles = 0\n    for i, img_crop in tqdm(enumerate(img_tiles)):\n        img_c = 255-img_crop # image complement\n        plt.subplot(rows, cols, i+1)\n        plt.xticks([])\n        plt.yticks([])\n        if img_c.mean() > thres: # pixel value mean greater threshold -> select tile \n            selected_tiles += 1\n            plt.imshow(img_crop)\n            plt.imshow(mask_tiles[i], cmap='seismic', alpha=0.4)\n            plt.tight_layout()\n    plt.suptitle(f'Num Tiles: {len(img_tiles)}\\nNum Selected Tiles: {selected_tiles}\\nThreshold: {thres}')\n    plt.tight_layout()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.776648Z","iopub.execute_input":"2022-08-08T20:48:22.778832Z","iopub.status.idle":"2022-08-08T20:48:22.800248Z","shell.execute_reply.started":"2022-08-08T20:48:22.778724Z","shell.execute_reply":"2022-08-08T20:48:22.798135Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def tile_saving(ID, df, thres, img_folder, mask_folder, rescale_size=(1024,1024), tile_size=256):\n    # load original img and mask\n    img  = tiff.imread(TRAIN_PATH + str(ID) +'.tiff')\n    mask = rle2mask(df.loc[df['id'] == ID].rle.values[0], (full_img_size, full_img_size))\n    \n    # rescale to 1024x1024\n    scaled_img, scaled_mask = rescale(img, mask, size=rescale_size)\n    \n    # make tiles + select based on given threshold\n    img_tiles, mask_tiles = make_tiles(scaled_img, scaled_mask, tile_size)\n    for i, img_crop in tqdm(enumerate(img_tiles)):\n        img_c = 255-img_crop # image complement\n        if img_c.mean() > thres: # pixel value mean greater threshold -> select tile \n            cv2.imwrite(os.path.join(img_folder,f'{ID}_{i}.png'), img_crop)\n            cv2.imwrite(os.path.join(mask_folder, f'{ID}_{i}.png'), mask_tiles[i])\n    print('-------------Done------------')","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.802637Z","iopub.execute_input":"2022-08-08T20:48:22.803845Z","iopub.status.idle":"2022-08-08T20:48:22.822286Z","shell.execute_reply.started":"2022-08-08T20:48:22.803778Z","shell.execute_reply":"2022-08-08T20:48:22.820379Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Tiling","metadata":{}},{"cell_type":"code","source":"# resize img, mask to 1024x1024 before tiling\nID = 10666\nimg = tiff.imread(f'../input/hubmap-organ-segmentation/train_images/{ID}.tiff')\nmask = rle2mask(train_df.loc[train_df['id'] == ID].rle.values[0], (full_img_size, full_img_size))\nscaled_img, scaled_mask = rescale(img, mask, size=(1024,1024))","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.825614Z","iopub.execute_input":"2022-08-08T20:48:22.826753Z","iopub.status.idle":"2022-08-08T20:48:22.893350Z","shell.execute_reply.started":"2022-08-08T20:48:22.826682Z","shell.execute_reply":"2022-08-08T20:48:22.891564Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimg_tiles, mask_tiles = make_tiles(scaled_img, scaled_mask, tile_size=128)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.896038Z","iopub.execute_input":"2022-08-08T20:48:22.897119Z","iopub.status.idle":"2022-08-08T20:48:22.909850Z","shell.execute_reply.started":"2022-08-08T20:48:22.897070Z","shell.execute_reply":"2022-08-08T20:48:22.907287Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"img_tiles_256, mask_tiles_256 = make_tiles(scaled_img, scaled_mask, tile_size=256)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.911970Z","iopub.execute_input":"2022-08-08T20:48:22.912795Z","iopub.status.idle":"2022-08-08T20:48:22.924187Z","shell.execute_reply.started":"2022-08-08T20:48:22.912749Z","shell.execute_reply":"2022-08-08T20:48:22.922325Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tiles = [img_tiles_256, img_tiles ]\nmasks = [mask_tiles_256, mask_tiles]","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.926622Z","iopub.execute_input":"2022-08-08T20:48:22.928351Z","iopub.status.idle":"2022-08-08T20:48:22.935277Z","shell.execute_reply.started":"2022-08-08T20:48:22.928283Z","shell.execute_reply":"2022-08-08T20:48:22.933576Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"visualize_tiles(scaled_img, scaled_mask, tiles, masks)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:22.937705Z","iopub.execute_input":"2022-08-08T20:48:22.939466Z","iopub.status.idle":"2022-08-08T20:48:39.851409Z","shell.execute_reply.started":"2022-08-08T20:48:22.939413Z","shell.execute_reply":"2022-08-08T20:48:39.850022Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### *There is a trade-off between different tile sizes. We can notice that a larger tile size contains more information about FTUs. In contrast, a smaller tile size leads to a larger number of tiles. However, many of these tiles only contain background or tissue regions with no FTUs. Thus, we need to be cautious when choosing tile size to yield a better result for training our model.*","metadata":{}},{"cell_type":"markdown","source":"## Tile Selection\n\nBased on the thresholds from **section 2** (Otsu, or Triangle), we will select tiles that contain tissue regions and discard background tiles. Since image pixels are now separated into 2 clusters with intensity values 0 and 255 (black, white), tiles consist pixel values' mean greater than a given threshold are selected for further analysis. A tissue region is determined by thresholding the grayscale image complement in **section 2**.","metadata":{}},{"cell_type":"code","source":"tile_selection(scaled_img, scaled_mask, thres_otsu, tile_size=128)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:48:39.853529Z","iopub.execute_input":"2022-08-08T20:48:39.854754Z","iopub.status.idle":"2022-08-08T20:49:02.909142Z","shell.execute_reply.started":"2022-08-08T20:48:39.854692Z","shell.execute_reply":"2022-08-08T20:49:02.907993Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### *We can discard many unnecessary tiles such as background and artefacts by a simple thresholding technique. However, there are many tiles containing tissue regions without FTUs. We can work on this a bit more and obtain a suitable threshold to increase the quality of our image training data*","metadata":{}},{"cell_type":"markdown","source":"## Tile Saving\n\nLet's save our selected image and mask tiles into the following folder structure:\n\n```\ntrain\n├── image\n└── mask  \n```","metadata":{}},{"cell_type":"code","source":"!mkdir train train/image train/mask\nimg_folder  = './train/image'\nmask_folder = './train/mask'","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:49:02.911246Z","iopub.execute_input":"2022-08-08T20:49:02.912088Z","iopub.status.idle":"2022-08-08T20:49:04.123208Z","shell.execute_reply.started":"2022-08-08T20:49:02.912042Z","shell.execute_reply":"2022-08-08T20:49:04.121503Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ID = 10666\ntile_saving(ID, train_df, thres_otsu, img_folder, mask_folder, tile_size=128)","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:49:04.125725Z","iopub.execute_input":"2022-08-08T20:49:04.126166Z","iopub.status.idle":"2022-08-08T20:49:04.311023Z","shell.execute_reply.started":"2022-08-08T20:49:04.126126Z","shell.execute_reply":"2022-08-08T20:49:04.309956Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!zip -r train.zip train","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-08-08T20:49:04.312466Z","iopub.execute_input":"2022-08-08T20:49:04.314088Z","iopub.status.idle":"2022-08-08T20:49:05.593439Z","shell.execute_reply.started":"2022-08-08T20:49:04.314032Z","shell.execute_reply":"2022-08-08T20:49:05.591046Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Sanity Check","metadata":{}},{"cell_type":"code","source":"img_path = './train/image'\nmask_path  = './train/mask'\nfor _, _, files in os.walk(img_path):\n    plt.figure(figsize=(10,6))\n    plt.subplots_adjust(left=0.5,\n                        right=0.6,\n                        top=0.6,\n                        bottom=0.5\n                       )\n    for i, name in enumerate(files):\n        selected_img = cv2.imread(os.path.join(img_path, name), cv2.IMREAD_UNCHANGED)\n        selected_msk = cv2.imread(os.path.join(mask_path, name), cv2.IMREAD_UNCHANGED)\n        plt.subplot(4,8,i+1)\n        plt.imshow(selected_img)\n        plt.imshow(selected_msk, cmap='seismic', alpha=0.4)\n        plt.xticks([])\n        plt.yticks([])\n    plt.tight_layout()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-08-08T20:49:05.595856Z","iopub.execute_input":"2022-08-08T20:49:05.596277Z","iopub.status.idle":"2022-08-08T20:49:07.757750Z","shell.execute_reply.started":"2022-08-08T20:49:05.596241Z","shell.execute_reply":"2022-08-08T20:49:07.756750Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<a id=\"conclusion\"></a>\n\n<h2 style=\"font-family: Verdana; font-size: 24px; font-style: normal; font-weight: bold; text-decoration: none; text-transform: none; letter-spacing: 3px; background-color: #CCCCFF; color: black;\" id=\"conclusion\"><left><br>&nbsp4. CONCLUSION <a href=\"#toc\">&#10514;</a><br></left> </h2>\n\n## WSI preprocessing framework for this challenge:\n\n### Steps:\n**1.** **Scale** down images\n\n**2.** **Filter** foreground/background (tissue/non-tissue)\n\n**3.**  **Break** down images into **tiles**\n\n**4.** **Select** tiles containing large portion of **tissue regions**\n\n**5.** **Retrieve** and **save** those tiles for further analysis\n\n### Remarks:\n* ID 12784 from kidney has no FTUs -> we need to thoroughly **scan** through all images to **remove** those faulted images if they exist\n* **Otsu** vs. **Triangle** binarization can only separate background and foreground but cannot detect artefacts. For example, Otsu treats lighter-colored tissue regions as background, whereas triangle binarization treats artefact as tissue regions. Thus, we need to dive deeper into these two thresholding techniques to have better tissue segmentation.\n* There is a **trade-off** between different tile sizes. We can notice that a larger tile size contains more information about FTUs. In contrast, a smaller tile size leads to a larger number of tiles. However, many of these tiles only contain background or tissue regions with no FTUs. Thus, we need to be cautious when choosing tile size to yield a better result for training our model.\n* We can **discard** many **unnecessary tiles** such as background and artefacts by a simple thresholding technique. However, there are many tiles containing tissue regions without FTUs. We can work on this a bit more and obtain a suitable threshold to increase the quality of our image training data\n\n","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}