{"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":"# 1. Introduction\n\nThis is a notebook for a past kaggle competition [**HuBMAP - Hacking the Human Vasculature**.](https://www.kaggle.com/competitions/hubmap-hacking-the-human-vasculature) The goal of this competition is to detect Blood Vessels from images of kidney tissue taken by microscope, and detecton mask shall have [IoU (Intersection over Union)](https://learnopencv.com/intersection-over-union-iou-in-object-detection-and-segmentation/) greater than 0.6. The kaggle score is calculated by Average Precision Over Confidence which is the same as [Open Images 2019 - Instance Segmentation](https://www.kaggle.com/c/open-images-2019-instance-segmentation/overview/evaluation). \n\nThis HuBMAP competition seems to be difficult for few reasons. First, it looks almost impossible to correctly identify blood vessels by non-experts(refer to EDA). It would be hard for deep learning too. Second, only 1633 label data are given for images which is insufficient for such a complex segemtation. Another reason is imbalanced label: only 3.3% of image pixels are positive. Therefore, I defined following targets for this project. \n\n**Target**\n\n* To develop an effective training loop which includes **data augmentation** and **custom loss function**\n* To compare score and calculation time between two well-known model (**FCN and U-NET**)\n* **Average IoU > 0.3** in validation data when **20% of images are selected for validation**.\n* To submit prediction to kaggle","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport os\nimport cv2\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nimport gc\nimport time\nimport math\nimport json\n\nimport tensorflow as tf\nfrom keras import backend as K\nfrom tensorflow.keras import layers\nfrom tensorflow.keras import Model","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-10-02T16:33:08.363362Z","iopub.execute_input":"2023-10-02T16:33:08.363952Z","iopub.status.idle":"2023-10-02T16:33:18.888157Z","shell.execute_reply.started":"2023-10-02T16:33:08.363911Z","shell.execute_reply":"2023-10-02T16:33:18.887138Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 2. Data\n\nThere are 1633 training images with their label (`polygon`), tile information (`tile_df`) which indicate datasource and source wsi. Source wsi is profiles of human subjects described in`wsi_df`.\n\n**Data Source**\n\nhttps://www.kaggle.com/competitions/hubmap-hacking-the-human-vasculature/data","metadata":{}},{"cell_type":"code","source":"train_folder = \"/kaggle/input/hubmap-hacking-the-human-vasculature/train/\"\ntest_folder = \"/kaggle/input/hubmap-hacking-the-human-vasculature/test/\"\nwsi_fpath = \"/kaggle/input/hubmap-hacking-the-human-vasculature/wsi_meta.csv\"\ntile_fpath = \"/kaggle/input/hubmap-hacking-the-human-vasculature/tile_meta.csv\"\nsample_fpath = \"/kaggle/input/hubmap-hacking-the-human-vasculature/sample_submission.csv\"\npolygon_fpath = \"/kaggle/input/hubmap-hacking-the-human-vasculature/polygons.jsonl\"","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:18.889679Z","iopub.execute_input":"2023-10-02T16:33:18.890299Z","iopub.status.idle":"2023-10-02T16:33:18.895564Z","shell.execute_reply.started":"2023-10-02T16:33:18.890273Z","shell.execute_reply":"2023-10-02T16:33:18.894294Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"with open(polygon_fpath) as f:\n    polygon = [json.loads(line) for line in f]\n    \nNP = len(polygon)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:18.897362Z","iopub.execute_input":"2023-10-02T16:33:18.898171Z","iopub.status.idle":"2023-10-02T16:33:23.012058Z","shell.execute_reply.started":"2023-10-02T16:33:18.898135Z","shell.execute_reply":"2023-10-02T16:33:23.011082Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"N of train data = \", NP)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.014518Z","iopub.execute_input":"2023-10-02T16:33:23.015511Z","iopub.status.idle":"2023-10-02T16:33:23.021185Z","shell.execute_reply.started":"2023-10-02T16:33:23.015476Z","shell.execute_reply":"2023-10-02T16:33:23.019995Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_files = os.listdir(train_folder)\ntile_df = pd.read_csv(tile_fpath)\nwsi_df = pd.read_csv(wsi_fpath)\nsample_df = pd.read_csv(sample_fpath)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.022711Z","iopub.execute_input":"2023-10-02T16:33:23.023475Z","iopub.status.idle":"2023-10-02T16:33:23.702726Z","shell.execute_reply.started":"2023-10-02T16:33:23.023423Z","shell.execute_reply":"2023-10-02T16:33:23.701791Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tile_df","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.704283Z","iopub.execute_input":"2023-10-02T16:33:23.704958Z","iopub.status.idle":"2023-10-02T16:33:23.723368Z","shell.execute_reply.started":"2023-10-02T16:33:23.704922Z","shell.execute_reply":"2023-10-02T16:33:23.722211Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"wsi_df","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.725294Z","iopub.execute_input":"2023-10-02T16:33:23.725648Z","iopub.status.idle":"2023-10-02T16:33:23.738172Z","shell.execute_reply.started":"2023-10-02T16:33:23.725616Z","shell.execute_reply":"2023-10-02T16:33:23.736981Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_fnames = []\nfor i in range(NP):\n    train_fnames.append(polygon[i][\"id\"])\n    \ntrain_name_df = pd.DataFrame({\"id\":train_fnames})\ntrain_name_df[\"count\"] = 1\n\ntile_df = pd.read_csv(tile_fpath)\ntile_df = tile_df.merge(train_name_df, on = \"id\", how = \"left\")\ntile_df = tile_df.query(\"count > 0\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.739703Z","iopub.execute_input":"2023-10-02T16:33:23.740702Z","iopub.status.idle":"2023-10-02T16:33:23.779478Z","shell.execute_reply.started":"2023-10-02T16:33:23.740664Z","shell.execute_reply":"2023-10-02T16:33:23.778397Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Reading Images**\n\nThis code read images from folder. ","metadata":{}},{"cell_type":"code","source":"time1 = time.time()\ni = 0\nL = 512\nX_train = np.zeros((NP, L, L, 3), dtype =  np.uint8)\n\nfor i in range(NP):\n    #cv2 read images as BGR that should be converted into RGB\n    img = cv2.imread(train_folder + train_fnames[i] + \".tif\")[:,:,::-1]    \n    X_train[i,:,:,:] = img\n    \ntime2 = time.time()\ntime3 = np.round(time2 - time1)\nprint(time3, \"sec\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:33:23.781401Z","iopub.execute_input":"2023-10-02T16:33:23.781796Z","iopub.status.idle":"2023-10-02T16:34:37.089781Z","shell.execute_reply.started":"2023-10-02T16:33:23.781759Z","shell.execute_reply":"2023-10-02T16:34:37.088805Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 3. Data Cleaning\n\nQuality fo image is good and unnecessary to clean. The problem is that the label is given as polygon (geometry) data. For image segmentation, it shall be converted to pixelwise label. This operation is accomplished by [other notebook](https://www.kaggle.com/code/hidetaketakahashi/hubmap-create-mask). Put simply, it checks whether pixels are in polygon or not by breath-first-seach. BFS (starts from one of polygon element) significantly reduced calculation time compared with checking all the elements. ","metadata":{}},{"cell_type":"code","source":"# This data is pixel wise label creatated by other notebook\nlabel_mat = np.load(\"/kaggle/input/hubmap-label/mask_mat.npy\")\nlabel_mat[label_mat > 1] = 0\nlabel_mat = label_mat.reshape(NP, 512, 512, 1)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:37.094812Z","iopub.execute_input":"2023-10-02T16:34:37.095421Z","iopub.status.idle":"2023-10-02T16:34:39.640009Z","shell.execute_reply.started":"2023-10-02T16:34:37.095387Z","shell.execute_reply":"2023-10-02T16:34:39.638665Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#This is code to copy polygon data into matrix \n\nnv = len(polygon[0][\"annotations\"])\n\nvmat =np.zeros((512, 512), dtype =  np.uint8)\n\nfor i in range(nv):\n    type1 = polygon[0][\"annotations\"][i][\"type\"]\n    if type1 == \"blood_vessel\":\n        crd = polygon[0][\"annotations\"][i][\"coordinates\"][0]\n        \n        for x, y in crd:\n            #pixels in polygon line is 1, otherwize 0\n            vmat[y, x] = 1","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:39.641663Z","iopub.execute_input":"2023-10-02T16:34:39.642339Z","iopub.status.idle":"2023-10-02T16:34:39.659145Z","shell.execute_reply.started":"2023-10-02T16:34:39.642303Z","shell.execute_reply":"2023-10-02T16:34:39.651974Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Following figure shows one of given image, its corresponding label, and pixel wise label.","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(1,3, figsize = (10, 3))\n\nax[0].imshow(X_train[0])\nax[1].imshow(vmat)\nax[2].imshow(label_mat[0])\n\nax[0].set_title(\"Image\")\nax[1].set_title(\"given label polygon\")\nax[2].set_title(\"pixel wise label (preprocessed)\")\n","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:39.661384Z","iopub.execute_input":"2023-10-02T16:34:39.662212Z","iopub.status.idle":"2023-10-02T16:34:40.415228Z","shell.execute_reply.started":"2023-10-02T16:34:39.662179Z","shell.execute_reply":"2023-10-02T16:34:40.414260Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 4. EDA","metadata":{}},{"cell_type":"markdown","source":"## 4.1 Basic Information\n\nImages and labels are provided from two datasets and four sources (persons). Number of data is not equal between them. Kaggle test data is taken from only dataset1, but they are invisible for participants. ","metadata":{}},{"cell_type":"code","source":"select = [\"dataset\", \"source_wsi\", \"count\"]\n\nsum_data_df = tile_df[select].groupby(select[0:-1]).sum().astype(int)\nsum_data_df.plot(kind = \"bar\", title = \"count of (datast, source wsi)\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:40.416304Z","iopub.execute_input":"2023-10-02T16:34:40.416643Z","iopub.status.idle":"2023-10-02T16:34:40.713239Z","shell.execute_reply.started":"2023-10-02T16:34:40.416607Z","shell.execute_reply":"2023-10-02T16:34:40.712285Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_photo(img, title = \"\"):\n    \n    N = img.shape[0]\n    \n    NC = 5\n    NR =  math.ceil(N/NC)\n    fig, ax = plt.subplots(NR, NC, figsize = (12, NR*2.3))\n    \n    for k in range(N):\n        i =  int(k/NC)\n        j = k % NC\n        \n        if N <= 5:\n            ax[j].imshow(img[k])\n            ax[j].tick_params(left = False, right = False , labelleft = False ,\n                        labelbottom = False, bottom = False)\n        else:\n            ax[i,j].imshow(img[k])\n            ax[i,j].tick_params(left = False, right = False , labelleft = False ,\n                        labelbottom = False, bottom = False)\n    \n    if title != \"\":\n        if N <= 5:\n            ax[2].set_title(title)\n        else:\n            ax[0,2].set_title(title)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:40.714578Z","iopub.execute_input":"2023-10-02T16:34:40.715526Z","iopub.status.idle":"2023-10-02T16:34:40.725093Z","shell.execute_reply.started":"2023-10-02T16:34:40.715493Z","shell.execute_reply":"2023-10-02T16:34:40.724105Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 4.2 Images and Labels\n\nNext, images and labels are visualized by each dataset and source (person). Size and quantity of blood vessels are different between sources. ","metadata":{}},{"cell_type":"code","source":"#This code select some samples from each dataset and source\nnp.random.seed(1)\n\nsamples_ds = []\nfor ds in [1,2]:\n    filter1 = tile_df[\"dataset\"]==ds\n    samples = []\n    for wsi in [1,2,3,4]:\n        \n        if ds ==1:\n            if wsi in [1,2]:\n                filter2 = filter1 & (tile_df[\"source_wsi\"] == wsi)\n                s_idx = np.random.choice(np.where(filter2)[0], 50)\n                \n        else:\n            filter2 = filter1 & (tile_df[\"source_wsi\"] == wsi)\n            s_idx = np.random.choice(np.where(filter2)[0], 50)\n            \n        samples.append(s_idx)\n        \n    samples_ds.append(samples)\n","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:40.726478Z","iopub.execute_input":"2023-10-02T16:34:40.727127Z","iopub.status.idle":"2023-10-02T16:34:40.744984Z","shell.execute_reply.started":"2023-10-02T16:34:40.727092Z","shell.execute_reply":"2023-10-02T16:34:40.743754Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Dataset1**\n\nIn the dataset 1, source 1 contains larger size but small number of blood vessels, whereas source 2 has many smaller size blood vessels. ","metadata":{}},{"cell_type":"code","source":"plot_photo(X_train[samples_ds[0][0][0:5]], \"dataset 1, source 1\")\nplot_photo(label_mat[samples_ds[0][0][0:5]], \"dataset 1, source 1, Label\")\nplot_photo(X_train[samples_ds[0][1][0:5]], \"dataset 1, source 2\")\nplot_photo(label_mat[samples_ds[0][1][0:5]], \"dataset 1, source 2, Label\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:40.746908Z","iopub.execute_input":"2023-10-02T16:34:40.747288Z","iopub.status.idle":"2023-10-02T16:34:44.517055Z","shell.execute_reply.started":"2023-10-02T16:34:40.747253Z","shell.execute_reply":"2023-10-02T16:34:44.516119Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Dataset 2**\n\nImage and labels of dataset 2 looks similar to dataset 1, but their color is slightly different. ","metadata":{}},{"cell_type":"code","source":"plot_photo(X_train[samples_ds[1][0][0:5]], \"dataset 2, source 1\")\nplot_photo(label_mat[samples_ds[1][0][0:5]], \"dataset 2, source 1, Label\")\n\nplot_photo(X_train[samples_ds[1][1][0:5]], \"dataset 2, source 2\")\nplot_photo(label_mat[samples_ds[1][1][0:5]], \"dataset 2, source 2, Label\")\n\nplot_photo(X_train[samples_ds[1][2][0:5]], \"dataset 2, source 3\")\nplot_photo(label_mat[samples_ds[1][2][0:5]], \"dataset 2, source 3, Label\")\n\nplot_photo(X_train[samples_ds[1][3][0:5]], \"dataset 2, source 4\")\nplot_photo(label_mat[samples_ds[1][3][0:5]], \"dataset 2, source 4, Label\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:44.518641Z","iopub.execute_input":"2023-10-02T16:34:44.519227Z","iopub.status.idle":"2023-10-02T16:34:51.887159Z","shell.execute_reply.started":"2023-10-02T16:34:44.519192Z","shell.execute_reply":"2023-10-02T16:34:51.886246Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 4.3 Data Imbalance\n\nSince it was hard to clearly explain difference of distribution, following items are calculated by each image. Then distribution are plotted. \n\n* ratio of positive pixel (blood vessel)\n* Mean size and quantity of blood vessel\n\nThe function `label_stat` calculates number of blood vessel and mean size of them. OpenCV's function `cv2.connectedComponents` takes boolean image as input, and count number of connected area, then it assigns label into each of them in 2D matrix. `densty_plot` and `scatter_plot` visualize ratio of positive pixels, mean size and quantity of blood vessels respectively.","metadata":{}},{"cell_type":"code","source":"def label_stat(Y_mat):\n    \n    #n_labels = number of labels\n    # labels_mat= 2D matrix in which label of connected areas are assigned. \n    \n    n_labels, labels_mat = cv2.connectedComponents(Y_mat)\n    size_list = []\n    for i in range(1, n_labels):\n        #Calculate number of pixels for each label\n        filter1 = labels_mat == i\n        size = filter1.sum()\n        size_list.append(size)\n\n    return n_labels - 1, np.mean(size_list)\n\n\ndef density_plot(dataset, n_source):\n    \n    fig, ax = plt.subplots(1,n_source, figsize = (3.5*n_source,2.7))\n    \n    for i in range(n_source):\n        filter1 = (tile_df[\"dataset\"] == dataset) & (tile_df[\"source_wsi\"] == i + 1)\n        sns.kdeplot(tile_df[\"positive_ratio\"].loc[filter1], ax = ax[i])\n        ax[i].set_xlim(0, 0.5)\n        ax[i].set_title(\"(dataset, source) = (\" +str(dataset) + \", \" + str(i+1)  + \")\" )\n        if i > 0:\n            ax[i].set_ylabel(\"\")\n            \n    plt.show()\n    \n\ndef scatter_plot(dataset, n_source):\n    fig, ax = plt.subplots(1,n_source, figsize = (3.5*n_source,2.7))\n\n    for i in range(n_source):\n        filter1 = (tile_df[\"dataset\"] == dataset) & (tile_df[\"source_wsi\"] == i + 1)\n        sns.scatterplot(data = tile_df.loc[filter1], x =  \"n_BloodVessel\", y = \"meanSize_BloodVessel\", ax = ax[i], alpha = 0.4)\n        ax[i].set_xlim(0, 35)\n        ax[i].set_ylim(0, 15)\n        if i > 0:\n            ax[i].set_ylabel(\"\")\n        ax[i].set_title(\"(dataset, source) = (\" +str(dataset) + \", \" + str(i+1)  + \")\" )\n    \n    ax[0].set_ylabel(\"mean size (1000 pixel)\")\n    plt.show()\n    \n","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:51.888595Z","iopub.execute_input":"2023-10-02T16:34:51.889170Z","iopub.status.idle":"2023-10-02T16:34:51.899872Z","shell.execute_reply.started":"2023-10-02T16:34:51.889135Z","shell.execute_reply":"2023-10-02T16:34:51.898964Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"n_label_list = []\nmean_size_list = []\nratio_list = []\nfor i in range(NP):\n    n_label, mean_size = label_stat(label_mat[i])\n    n_label_list.append(n_label)\n    mean_size_list.append(mean_size)\n    ratio_list.append(np.mean(label_mat[i]))","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:51.901232Z","iopub.execute_input":"2023-10-02T16:34:51.901890Z","iopub.status.idle":"2023-10-02T16:34:56.931453Z","shell.execute_reply.started":"2023-10-02T16:34:51.901858Z","shell.execute_reply":"2023-10-02T16:34:56.930401Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tile_df[\"n_BloodVessel\"] = n_label_list\ntile_df[\"meanSize_BloodVessel\"] = np.array(mean_size_list)/1000\ntile_df[\"positive_ratio\"] = ratio_list","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:56.932819Z","iopub.execute_input":"2023-10-02T16:34:56.933183Z","iopub.status.idle":"2023-10-02T16:34:56.942946Z","shell.execute_reply.started":"2023-10-02T16:34:56.933149Z","shell.execute_reply":"2023-10-02T16:34:56.942123Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Ratio of positive pixel\n\nOverall, only 3.3% of pixels are positive (blood vessel). When it is calculated by each dataset and source, dataset 1 has more positive pixels than dataset 2.","metadata":{}},{"cell_type":"code","source":"pct = np.round(tile_df[\"positive_ratio\"].mean()*100, 1)\nprint(pct, \"% of pixels are positive in all data\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:56.944654Z","iopub.execute_input":"2023-10-02T16:34:56.945518Z","iopub.status.idle":"2023-10-02T16:34:56.965249Z","shell.execute_reply.started":"2023-10-02T16:34:56.945477Z","shell.execute_reply":"2023-10-02T16:34:56.964038Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tile_df[[\"dataset\", \"source_wsi\", \"positive_ratio\"]].groupby([\"dataset\", \"source_wsi\"]).mean()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:56.966890Z","iopub.execute_input":"2023-10-02T16:34:56.968129Z","iopub.status.idle":"2023-10-02T16:34:56.994200Z","shell.execute_reply.started":"2023-10-02T16:34:56.968087Z","shell.execute_reply":"2023-10-02T16:34:56.992910Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The ratio positive pixel is calculated by each image, and the distribution is plotted in following density plots.","metadata":{}},{"cell_type":"markdown","source":"\n\n**Dataset 1**","metadata":{}},{"cell_type":"code","source":"density_plot(1, 2)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:56.996396Z","iopub.execute_input":"2023-10-02T16:34:56.997243Z","iopub.status.idle":"2023-10-02T16:34:57.427107Z","shell.execute_reply.started":"2023-10-02T16:34:56.997202Z","shell.execute_reply":"2023-10-02T16:34:57.426153Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Dataset 2**","metadata":{}},{"cell_type":"code","source":"density_plot(2, 4)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:57.428601Z","iopub.execute_input":"2023-10-02T16:34:57.429265Z","iopub.status.idle":"2023-10-02T16:34:58.150146Z","shell.execute_reply.started":"2023-10-02T16:34:57.429211Z","shell.execute_reply":"2023-10-02T16:34:58.149238Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Size and quantity of Blood Vessels\n\nDataset 1 source 1 has different distribution from others. It has less number of blood vessels, but their size are larger. Please note that **one dot corresponds to one image** in scatterplots.","metadata":{}},{"cell_type":"markdown","source":"**dataset1**","metadata":{}},{"cell_type":"code","source":"scatter_plot(1, 2)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:58.151645Z","iopub.execute_input":"2023-10-02T16:34:58.152707Z","iopub.status.idle":"2023-10-02T16:34:58.577385Z","shell.execute_reply.started":"2023-10-02T16:34:58.152670Z","shell.execute_reply":"2023-10-02T16:34:58.576494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**dataset2**","metadata":{}},{"cell_type":"code","source":"scatter_plot(2, 4)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:58.578981Z","iopub.execute_input":"2023-10-02T16:34:58.579626Z","iopub.status.idle":"2023-10-02T16:34:59.360010Z","shell.execute_reply.started":"2023-10-02T16:34:58.579588Z","shell.execute_reply":"2023-10-02T16:34:59.359086Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 5. Models\n\nFor image segmentation, Fully Convolutional Network and U-net are developed. Both of them are trained, then validation results are compared. ","metadata":{}},{"cell_type":"markdown","source":"## 5.1 FCN (Fully Convolutional Network)\n\nFully Convolutional Network (FCN) has two type of elements. First,`ConvBlock` which consists of `Conv2D` layer followed by `BatchNormalization` and `Relu`. It downsamples input by `MaxPooling2D` upon request. `ConvBlock` are allocated sequantially just like VGG16 or AlexNet architecture. \n\nUnlike architecture for image classificaiton, however, FCN does not have fully connected layer. Instead, it has `UpsampleBlock` which takes input of size 64x64x512, then generates 512x512x1 by `Conv2DTranspose`.","metadata":{}},{"cell_type":"code","source":"def ConvBlock(channel, X, ksize = 3, downsample = True):\n    \n    if downsample:\n        X =  layers.MaxPooling2D(pool_size=(2, 2), strides = (2,2))(X)        \n    \n    X = layers.Conv2D(channel, kernel_size = ksize, strides = 1, padding = \"same\")(X)        \n    X = layers.BatchNormalization()(X)\n    X = layers.ReLU()(X)\n    \n    return X\n\ndef UpsampleBlock(channel, X, ksize, stride):\n    \n    X = layers.Conv2DTranspose(channel, kernel_size=ksize, strides=stride, padding =\"same\")(X)\n        \n    return X","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:59.361316Z","iopub.execute_input":"2023-10-02T16:34:59.362203Z","iopub.status.idle":"2023-10-02T16:34:59.369533Z","shell.execute_reply.started":"2023-10-02T16:34:59.362162Z","shell.execute_reply":"2023-10-02T16:34:59.368651Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def create_FCN():\n    \n    L = 512\n    Input =  layers.Input(shape=(L, L, 3))\n    \n    X = layers.Rescaling(scale = 1./127.5, offset= -1, )(Input)\n    \n    KS = 3\n    channel = 64\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#512\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    \n    channel = channel*2\n    X = ConvBlock(channel, X, ksize = KS, downsample = True)#256\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    channel = channel*2\n    X = ConvBlock(channel, X, ksize = KS, downsample = True)#128\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    channel = channel*2\n    X = ConvBlock(channel, X, ksize = KS, downsample = True)#64\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    X = ConvBlock(channel, X, ksize = KS, downsample = False)#\n    \n    \n    X = UpsampleBlock(1, X, ksize = 8, stride = 8)\n\n    model = Model(inputs = Input, outputs = X)\n    \n    return model","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:59.377577Z","iopub.execute_input":"2023-10-02T16:34:59.378051Z","iopub.status.idle":"2023-10-02T16:34:59.387679Z","shell.execute_reply.started":"2023-10-02T16:34:59.378024Z","shell.execute_reply":"2023-10-02T16:34:59.386686Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model_FCN = create_FCN()\nmodel_FCN.summary()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:34:59.389297Z","iopub.execute_input":"2023-10-02T16:34:59.389725Z","iopub.status.idle":"2023-10-02T16:35:03.130326Z","shell.execute_reply.started":"2023-10-02T16:34:59.389679Z","shell.execute_reply":"2023-10-02T16:35:03.129492Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tf.keras.utils.plot_model(model_FCN)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:03.131392Z","iopub.execute_input":"2023-10-02T16:35:03.131711Z","iopub.status.idle":"2023-10-02T16:35:03.516977Z","shell.execute_reply.started":"2023-10-02T16:35:03.131679Z","shell.execute_reply":"2023-10-02T16:35:03.516020Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 5.2 U-NET\n\nU-net is basically similar to convolutional autoencoders that have downsampling decoder and upsampling encoder. What makes U-net unique is **skip connection from decoder to encoder at each level**. This structure looks similar to [V-model](https://en.wikipedia.org/wiki/V-model). \n\n`LeftBlock` written in the below code has `Conv2D` followed by `BatchNormalization` and `Relu`. It downsamples input by `MaxPooling2D` upon request. `RightBlock` upsamples input by `Conv2DTranspose` and it is concatenated with skip connection from `LeftBlock`. It is actually simple when plotted in the figure. ","metadata":{}},{"cell_type":"code","source":"def LeftBlock(channel, X, ksize = 3, downsample = True):\n    \n    if downsample:\n        X =  layers.MaxPooling2D(pool_size=(2, 2), strides = (2,2))(X)\n            \n    X = layers.Conv2D(channel, kernel_size = ksize, strides = 1, padding = \"same\")(X)\n    X = layers.BatchNormalization()(X)\n    X = layers.ReLU()(X)\n    \n    return X\n\n\ndef RightBlock(channel, X, ksize = 3, X_skip = None, upsample = True):\n    \n    if upsample:\n        X = layers.Conv2DTranspose(channel, kernel_size=4, strides=2, padding=\"SAME\")(X)\n    \n    if X_skip is not None:\n        X = layers.Concatenate()([X, X_skip])\n    \n    X =  layers.Conv2D(channel, kernel_size = ksize, strides = 1, padding = \"same\")(X)\n    X = layers.BatchNormalization()(X)\n    X = layers.ReLU()(X)\n\n    return X","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:03.518398Z","iopub.execute_input":"2023-10-02T16:35:03.519412Z","iopub.status.idle":"2023-10-02T16:35:03.529364Z","shell.execute_reply.started":"2023-10-02T16:35:03.519377Z","shell.execute_reply":"2023-10-02T16:35:03.528257Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def create_UNET():\n    \n    L = 512\n    Input =  layers.Input(shape=(L, L, 3))\n\n    X0 = layers.Rescaling(scale = 1./127.5, offset= -1, )(Input)\n    \n    KS = 5\n    \n    channel = 48\n    X1 = LeftBlock(channel, X0, ksize = KS, downsample = False)#512\n\n    channel = channel*2 #128\n    X2 = LeftBlock(channel, X1, ksize = KS, downsample = True)#256\n\n    channel = channel*2 #256\n    X3 = LeftBlock(channel, X2, ksize = KS, downsample = True)#128\n    \n    channel = channel*2 #512\n    X4 = LeftBlock(channel, X3, ksize = KS, downsample = True)#64\n\n    channel = channel #512\n    X5 = LeftBlock(channel, X4, ksize = KS, downsample = True)#32\n    \n    XR = RightBlock(channel, X5, ksize = KS, X_skip = X4) #64\n\n    channel = int(channel/2) #256\n    XR = RightBlock(channel, XR, ksize = KS, X_skip = X3) #128\n\n    channel = int(channel/2) #128\n    XR = RightBlock(channel, XR, ksize = KS, X_skip = X2) #256\n\n    channel = int(channel/2) #64\n    XR = RightBlock(channel, XR, ksize = KS, X_skip = X1) #512\n\n    channel = 1\n    XR = layers.Conv2D(channel, kernel_size = 1, strides = 1, padding = \"same\")(XR)\n\n    model = Model(inputs = Input, outputs = XR)\n    \n    return model","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:03.531035Z","iopub.execute_input":"2023-10-02T16:35:03.531797Z","iopub.status.idle":"2023-10-02T16:35:03.543026Z","shell.execute_reply.started":"2023-10-02T16:35:03.531744Z","shell.execute_reply":"2023-10-02T16:35:03.541972Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model_UNET = create_UNET()\nmodel_UNET.summary()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:03.544436Z","iopub.execute_input":"2023-10-02T16:35:03.545198Z","iopub.status.idle":"2023-10-02T16:35:03.969921Z","shell.execute_reply.started":"2023-10-02T16:35:03.545162Z","shell.execute_reply":"2023-10-02T16:35:03.969146Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Figure of U-net**","metadata":{}},{"cell_type":"code","source":"tf.keras.utils.plot_model(model_UNET)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:03.970955Z","iopub.execute_input":"2023-10-02T16:35:03.971315Z","iopub.status.idle":"2023-10-02T16:35:04.339338Z","shell.execute_reply.started":"2023-10-02T16:35:03.971281Z","shell.execute_reply":"2023-10-02T16:35:04.338469Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 6. Training and Hyperparameter Tuning","metadata":{}},{"cell_type":"markdown","source":"## 6.1 Split Dataset\n\n80% of image is randomly selected for training. ","metadata":{}},{"cell_type":"code","source":"batch_size = 4","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:04.340514Z","iopub.execute_input":"2023-10-02T16:35:04.341421Z","iopub.status.idle":"2023-10-02T16:35:04.345704Z","shell.execute_reply.started":"2023-10-02T16:35:04.341384Z","shell.execute_reply":"2023-10-02T16:35:04.344812Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.random.seed(1)\n\nn_data = NP\nn_train = int(n_data*0.8)\n\nall_idx = np.arange(n_data)\nnp.random.shuffle(all_idx)\n\ntrain_idx = all_idx[:n_train]\nval_idx = all_idx[n_train:]\n\nprint(\"N of data for training = \", train_idx.shape[0]) \nprint(\"N of data for validation = \", val_idx.shape[0])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:04.346967Z","iopub.execute_input":"2023-10-02T16:35:04.347704Z","iopub.status.idle":"2023-10-02T16:35:04.361370Z","shell.execute_reply.started":"2023-10-02T16:35:04.347668Z","shell.execute_reply":"2023-10-02T16:35:04.360221Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_ds = tf.data.Dataset.from_tensor_slices((X_train[train_idx], label_mat[train_idx])).shuffle(1000).batch(batch_size)\n\nX_val = X_train[val_idx]\nY_val = label_mat[val_idx]\nn_val = X_val.shape[0]","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:04.362903Z","iopub.execute_input":"2023-10-02T16:35:04.363439Z","iopub.status.idle":"2023-10-02T16:35:08.213105Z","shell.execute_reply.started":"2023-10-02T16:35:04.363407Z","shell.execute_reply":"2023-10-02T16:35:08.212067Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"`val_tile_df` is later used to analyze validation result.","metadata":{}},{"cell_type":"code","source":"val_tile_df = tile_df.iloc[val_idx].reset_index()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.214654Z","iopub.execute_input":"2023-10-02T16:35:08.215005Z","iopub.status.idle":"2023-10-02T16:35:08.220740Z","shell.execute_reply.started":"2023-10-02T16:35:08.214972Z","shell.execute_reply":"2023-10-02T16:35:08.219828Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.2 Custom Loss Function\n\nSince the label is imbalanced (refer to chapter 4 EDA), simple binary cross entropy would not work well. Instead, this `custom_loss` function calculate binary cross entropy for positive and negative pixels separately. Because tesorflow's `BCE` contains reduce_mean function, positive and negative losses can be balanced in this way. It is same as **weighted log loss**\nThen, positve loss and negative loss are summed with **weight of negative loss**. When the weight of negative loss is 2 or 3 , the result was good. ","metadata":{}},{"cell_type":"code","source":"# Custom Loss Function\n\nBCE = tf.keras.losses.BinaryCrossentropy(from_logits=True)\n\ndef custom_loss(y_true, y_pred):\n\n    filter1 = y_true == 1\n\n    p_loss = BCE(y_true[filter1], y_pred[filter1])\n\n    filter1 = y_true == 0\n    n_loss = BCE(y_true[filter1], y_pred[filter1])\n\n    loss =  p_loss + n_loss*3 #weighted log loss\n\n    return loss","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.222132Z","iopub.execute_input":"2023-10-02T16:35:08.222694Z","iopub.status.idle":"2023-10-02T16:35:08.236561Z","shell.execute_reply.started":"2023-10-02T16:35:08.222661Z","shell.execute_reply":"2023-10-02T16:35:08.235636Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.3 Image Augmentation\n\nThe funciton `augment` is applied to both image and label, because when image is flipped, label shall follow it. First part of augmentation is **horizontal and vertical flip** at probability 0.5 respectively. Next is image rotation. **The rotation angle is selected from 0, 90, 180 and 270**. Finally, **the image and the label is zoomed** at 0.3 probability. It select uppler left pixel randomly. Then the length is selected. It is clipped and resized to 512x512. ","metadata":{}},{"cell_type":"code","source":"#Image Augment Function\ndef augment(X, Y):\n    \n    #1. random flip--------------\n    #1.1 horizontal\n    if tf.random.uniform(shape=[1]) > 0.5:\n        X = X[:,:,::-1]\n        Y = Y[:,:,::-1]\n        \n    #1.2 vertical\n    if tf.random.uniform(shape=[1]) > 0.5:\n        X = X[:,::-1]\n        Y = Y[:,::-1]\n    \n    #2. rotation------------------\n    #2.1 set angle \n    angle = np.random.randint(4) \n    X = tf.image.rot90(X, k=angle)\n    Y = tf.image.rot90(Y, k=angle)\n    \n    #3. zoom -----------------------------\n    L = 512\n    if tf.random.uniform(shape=[1]) < 0.3:\n        \n        #left upper\n        x1 = np.random.choice(np.arange(0, int(L*0.2)), 1)[0]\n        y1 = np.random.choice(np.arange(0, int(L*0.2)), 1)[0]\n        \n        #Length\n        L2 = min(L-x1, L-y1)\n        L3 = np.random.choice(np.arange(int(L2*0.8), L2), 1)[0]\n        \n        X = X[:,y1:(y1+L3), x1:(x1+L3),:] \n        Y = Y[:,y1:(y1+L3), x1:(x1+L3),:] \n        \n        X = tf.cast(tf.image.resize(X,[L, L]), dtype = tf.uint8)\n        Y = tf.cast(tf.image.resize(Y,[L, L]), dtype = tf.uint8)\n    \n    return X, Y","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.237845Z","iopub.execute_input":"2023-10-02T16:35:08.238841Z","iopub.status.idle":"2023-10-02T16:35:08.249232Z","shell.execute_reply.started":"2023-10-02T16:35:08.238790Z","shell.execute_reply":"2023-10-02T16:35:08.248203Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.4 Custom train loop\n\nThe custom train loop `train_loop` consists of `augment` fuction, `train_step` which calculates loss and apply gradients, and some functions to calculate validation results (`cal_IoU`, `predict_probability`, `val_score`). The statement `@tf.function` before the `train_step` significantly accelerates training by creating computation graph. Furthermore, to keep the best model, `train_loop` saves trained model when validation loss is the lowest. ","metadata":{}},{"cell_type":"code","source":"@tf.function\ndef train_step(X, Y, model):\n        \n    with tf.GradientTape() as tape:\n        Y_pred = model(X)\n        loss = custom_loss(Y, Y_pred)\n        \n        \n    grads = tape.gradient(loss, model.trainable_weights)\n    optimizer.apply_gradients(zip(grads, model.trainable_weights))\n    \n    return loss","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.250704Z","iopub.execute_input":"2023-10-02T16:35:08.251404Z","iopub.status.idle":"2023-10-02T16:35:08.266830Z","shell.execute_reply.started":"2023-10-02T16:35:08.251361Z","shell.execute_reply":"2023-10-02T16:35:08.265858Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#calculate IoU from predicted prob\ndef cal_IoU(Y_true, prob, cutoff = 0):\n    \n    Y_pred = (prob > cutoff).astype(int)\n    \n    score_list = []\n    for k in range(prob.shape[0]):\n    \n        and_score =  np.sum(Y_pred[k][Y_true[k] == 1])\n        or_score = np.sum(Y_true[k]) + np.sum(Y_pred[k]) - and_score\n        \n        \n        if or_score == 0:\n            score = 1\n        else:\n            score = and_score/or_score\n        score_list.append(score)\n    \n    return score_list #np.round(np.mean(np.array(score_list)), 3) \n\n\n#calculation of probability from logit\ndef predict_probability(X, model):\n    pred = model.predict(X, verbose = 0)\n    odd = np.exp(pred)\n    prob = odd/(1+odd)\n    \n    return prob\n    \n# validation loss \ndef val_score(X, Y, model):\n    \n    Y_pred = model.predict(X)\n    loss =  custom_loss(Y, Y_pred)\n    \n    return loss","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.268334Z","iopub.execute_input":"2023-10-02T16:35:08.268681Z","iopub.status.idle":"2023-10-02T16:35:08.281246Z","shell.execute_reply.started":"2023-10-02T16:35:08.268649Z","shell.execute_reply":"2023-10-02T16:35:08.280322Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def train_loop(n_epoch, history, model, model_name, print_result = True):\n    \n    best_val_loss = 100000.\n    \n    for k in range(n_epoch):\n        \n        time1 = time.time()\n        \n        loss_list = []\n        for _, ds in enumerate(train_ds):\n            X, Y = ds #to extract X and Y from dataset\n\n            X, Y = augment(X, Y) #image augmentation\n            \n            loss = train_step(X, Y, model) #calculate loss and apply gradients\n            loss_list.append(loss)\n        \n        #recording loss\n        train_loss = np.mean(loss_list)    \n        val_loss = val_score(X_val, Y_val, model)\n        \n        #calculate IoU\n        prob = predict_probability(X_val, model)\n        val_IoU = cal_IoU(Y_val, prob, 0.8)\n        \n        #records loss and IoU in history\n        history[\"train_loss\"].append(train_loss)\n        history[\"val_loss\"].append(val_loss)\n        history[\"val_IoU\"].append(val_IoU)\n        \n        time2 = time.time()\n        time3 = np.round(time2- time1)\n        \n        if k >= 0:\n            #keep saving the best model\n            if val_loss < best_val_loss:\n                best_val_loss = val_loss\n                model.save_weights(model_name + \"/ckpt1\")\n                print(\"write model at epoch \", k)\n        \n        if print_result:\n            print(k, \"train loss\", train_loss, \", val loss \", val_loss, \", val_IoU\", val_IoU, \" time[s] = \", time3)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.282665Z","iopub.execute_input":"2023-10-02T16:35:08.283445Z","iopub.status.idle":"2023-10-02T16:35:08.298365Z","shell.execute_reply.started":"2023-10-02T16:35:08.283414Z","shell.execute_reply":"2023-10-02T16:35:08.297337Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#This is the function to save train and validation history.\ndef save_history(history, model_name):\n    \n    np.save(\"train_loss_\" + model_name, np.array(history[\"train_loss\"]))\n    np.save(\"val_loss_\" + model_name, np.array(history[\"val_loss\"]))\n    np.save(\"val_IoU_\" + model_name, np.array(history[\"val_IoU\"]))","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.299877Z","iopub.execute_input":"2023-10-02T16:35:08.300569Z","iopub.status.idle":"2023-10-02T16:35:08.316482Z","shell.execute_reply.started":"2023-10-02T16:35:08.300535Z","shell.execute_reply":"2023-10-02T16:35:08.315599Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.5 Hyperparameter Tuning\n\nAfter running notebooks many times, following hyperparameters were tuned. Only the results are shown here, because demonstrating all the process is not possible within limited calculation time. \n\n**General**\n\n* adam optimizer with learning rate = 0.0002. (much better than default lr 0.001) with batch size = 4\n* weight of negative loss (false positive) in `custom_loss` = 3\n* best epoch: 20 ~ 30 \n\n**FCN**\n* Number of conv block at each level = 3\n* Convolution kernel size = 3\n\n**U-NET**\n* Number of conv block at each level = 1\n* Convolution kernel size = 5","metadata":{}},{"cell_type":"markdown","source":"## 6.6 Training FCN","metadata":{}},{"cell_type":"markdown","source":"**This model is previously trained in the notebook to avoid OutOfMemoryError**, because training two models in one notebook was not possible. It loads trained results.","metadata":{}},{"cell_type":"code","source":"history_FCN = {}\nhistory_FCN[\"train_loss\"] = []\nhistory_FCN[\"val_loss\"] = []\nhistory_FCN[\"val_IoU\"] = []","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.317696Z","iopub.execute_input":"2023-10-02T16:35:08.318423Z","iopub.status.idle":"2023-10-02T16:35:08.332598Z","shell.execute_reply.started":"2023-10-02T16:35:08.318391Z","shell.execute_reply":"2023-10-02T16:35:08.331571Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"optimizer =  tf.keras.optimizers.Adam(learning_rate=0.0002)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.333980Z","iopub.execute_input":"2023-10-02T16:35:08.334474Z","iopub.status.idle":"2023-10-02T16:35:08.350201Z","shell.execute_reply.started":"2023-10-02T16:35:08.334444Z","shell.execute_reply":"2023-10-02T16:35:08.349224Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#train_loop(20, history_FCN, model_FCN, \"HuBMAP-FCN-model\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.351755Z","iopub.execute_input":"2023-10-02T16:35:08.352301Z","iopub.status.idle":"2023-10-02T16:35:08.357861Z","shell.execute_reply.started":"2023-10-02T16:35:08.352268Z","shell.execute_reply":"2023-10-02T16:35:08.356819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#save_history(history_FCN, \"FCN\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.359423Z","iopub.execute_input":"2023-10-02T16:35:08.360171Z","iopub.status.idle":"2023-10-02T16:35:08.369490Z","shell.execute_reply.started":"2023-10-02T16:35:08.360137Z","shell.execute_reply":"2023-10-02T16:35:08.368540Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**import trained model and result**","metadata":{}},{"cell_type":"code","source":"model_FCN.load_weights(\"/kaggle/input/hubmap-reportfiles/fcn_model/ckpt1\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:08.371077Z","iopub.execute_input":"2023-10-02T16:35:08.371702Z","iopub.status.idle":"2023-10-02T16:35:09.178778Z","shell.execute_reply.started":"2023-10-02T16:35:08.371668Z","shell.execute_reply":"2023-10-02T16:35:09.177765Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"history_FCN[\"train_loss\"] = np.load(\"/kaggle/input/hubmap-reportfiles/fcn_loss/train_loss_FCN.npy\")\nhistory_FCN[\"val_loss\"] = np.load(\"/kaggle/input/hubmap-reportfiles/fcn_loss/val_loss_FCN.npy\")\nhistory_FCN[\"val_IoU\"] = np.load(\"/kaggle/input/hubmap-reportfiles/fcn_loss/val_IoU_FCN.npy\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.180323Z","iopub.execute_input":"2023-10-02T16:35:09.180689Z","iopub.status.idle":"2023-10-02T16:35:09.217168Z","shell.execute_reply.started":"2023-10-02T16:35:09.180654Z","shell.execute_reply":"2023-10-02T16:35:09.216327Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_history(history):\n    \n    fig, ax = plt.subplots(figsize = (5,4))\n    ax.plot(history[\"train_loss\"], label = \"train loss\")    \n    ax.plot(history[\"val_loss\"], label = \"val loss\")\n    ax.set_title(\"loss\")\n    ax.legend()\n    ax.grid()\n    #fig, ax = plt.subplots(1,2, figsize = (10,4))\n    #ax[0].plot(history[\"train_loss\"], label = \"train loss\")    \n    #ax[0].plot(history[\"val_loss\"], label = \"val loss\")\n    \n    #ax[1].plot(history[\"val_IoU\"])\n    #ax[0].set_title(\"loss\")\n    #ax[0].legend()\n    #ax[1].set_title(\"val IoU (in all images)\")\n    #ax[0].grid()\n    #ax[1].grid()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.218644Z","iopub.execute_input":"2023-10-02T16:35:09.219246Z","iopub.status.idle":"2023-10-02T16:35:09.225823Z","shell.execute_reply.started":"2023-10-02T16:35:09.219212Z","shell.execute_reply":"2023-10-02T16:35:09.224296Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_history(history_FCN)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.230479Z","iopub.execute_input":"2023-10-02T16:35:09.231877Z","iopub.status.idle":"2023-10-02T16:35:09.558142Z","shell.execute_reply.started":"2023-10-02T16:35:09.231839Z","shell.execute_reply":"2023-10-02T16:35:09.557092Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 6.7 Training U-NET","metadata":{}},{"cell_type":"markdown","source":"**This model is previously trained in the notebook to avoid OutOfMemoryError**, because training two models in one notebook was not possible. It loads trained results. ","metadata":{}},{"cell_type":"code","source":"history_UNET = {}\nhistory_UNET[\"train_loss\"] = []\nhistory_UNET[\"val_loss\"] = []\nhistory_UNET[\"val_IoU\"] = []","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.559465Z","iopub.execute_input":"2023-10-02T16:35:09.560405Z","iopub.status.idle":"2023-10-02T16:35:09.566431Z","shell.execute_reply.started":"2023-10-02T16:35:09.560364Z","shell.execute_reply":"2023-10-02T16:35:09.565068Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#train_loop(20, history_UNET, model_UNET, \"HuBMAP-UNET-model\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.567963Z","iopub.execute_input":"2023-10-02T16:35:09.568645Z","iopub.status.idle":"2023-10-02T16:35:09.577957Z","shell.execute_reply.started":"2023-10-02T16:35:09.568609Z","shell.execute_reply":"2023-10-02T16:35:09.576934Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#save_history(history_UNET, \"UNET\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.579294Z","iopub.execute_input":"2023-10-02T16:35:09.580260Z","iopub.status.idle":"2023-10-02T16:35:09.590905Z","shell.execute_reply.started":"2023-10-02T16:35:09.580223Z","shell.execute_reply":"2023-10-02T16:35:09.589803Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**import trained result**","metadata":{}},{"cell_type":"code","source":"model_UNET.load_weights(\"/kaggle/input/hubmap-reportfiles/unet_model/ckpt1\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:09.592297Z","iopub.execute_input":"2023-10-02T16:35:09.593097Z","iopub.status.idle":"2023-10-02T16:35:11.109627Z","shell.execute_reply.started":"2023-10-02T16:35:09.593064Z","shell.execute_reply":"2023-10-02T16:35:11.108606Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"history_UNET[\"train_loss\"] = np.load(\"/kaggle/input/hubmap-reportfiles/unet_loss/train_loss_UNET.npy\")\nhistory_UNET[\"val_loss\"] = np.load(\"/kaggle/input/hubmap-reportfiles/unet_loss/val_loss_UNET.npy\")\nhistory_UNET[\"val_IoU\"] = np.load(\"/kaggle/input/hubmap-reportfiles/unet_loss/val_IoU_UNET.npy\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.111009Z","iopub.execute_input":"2023-10-02T16:35:11.111576Z","iopub.status.idle":"2023-10-02T16:35:11.150118Z","shell.execute_reply.started":"2023-10-02T16:35:11.111541Z","shell.execute_reply":"2023-10-02T16:35:11.149104Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_history(history_UNET)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.151430Z","iopub.execute_input":"2023-10-02T16:35:11.152308Z","iopub.status.idle":"2023-10-02T16:35:11.485573Z","shell.execute_reply.started":"2023-10-02T16:35:11.152268Z","shell.execute_reply":"2023-10-02T16:35:11.484667Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 7. Analysis of Validation Result\n\n## 7.1 Validation Loss\n\nThere are lots of ways to analyse results. First, validation loss of both models are compared. According to the plots, both model learned well but U-NET has slightly smaller loss than FCN.","metadata":{}},{"cell_type":"code","source":"def compare_val(history_FCN, history_UNET):\n    \n    fig, ax = plt.subplots(figsize = (6, 5))\n    ax.plot(history_FCN[\"val_loss\"], label = \"FCN val loss\")\n    ax.plot(history_UNET[\"val_loss\"], label = \"U-NET val loss\")\n    ax.set_title(\"val loss\")\n    ax.grid()\n    ax.legend()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.486845Z","iopub.execute_input":"2023-10-02T16:35:11.487732Z","iopub.status.idle":"2023-10-02T16:35:11.494505Z","shell.execute_reply.started":"2023-10-02T16:35:11.487697Z","shell.execute_reply":"2023-10-02T16:35:11.493493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"compare_val(history_FCN, history_UNET)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.495794Z","iopub.execute_input":"2023-10-02T16:35:11.496653Z","iopub.status.idle":"2023-10-02T16:35:11.794411Z","shell.execute_reply.started":"2023-10-02T16:35:11.496617Z","shell.execute_reply":"2023-10-02T16:35:11.793493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 7.2 Predicted labels\n\nSince both of FCN and U-net trained well, they must be able to predict label with good accuracy. The funciton `plot_result` visualize an image, its label, and prediction results by FCN and U-net.","metadata":{}},{"cell_type":"code","source":"def plot_result(X, Y_true, Y_pred1, Y_pred2, cutoff = 0.8):\n    \n    N = X.shape[0]\n    \n    fig, ax = plt.subplots(N,6, figsize = (13,2.2*N))\n    \n    for k in range(N):\n        \n        cutoff_img1 = (Y_pred1[k,:,:,0] > cutoff).astype(int)\n        cutoff_img2 = (Y_pred2[k,:,:,0] > cutoff).astype(int)\n\n        true_img = np.zeros((512, 512, 3), dtype = np.uint8)\n        true_img[:,:,1] = Y_true[k,:,:,0]*200\n        \n        cutoff1 = np.zeros((512, 512, 3), dtype = np.uint8)\n        cutoff2 = np.zeros((512, 512, 3), dtype = np.uint8)\n        \n        cutoff1[:,:,0] = cutoff_img1*230\n        cutoff2[:,:,0] = cutoff_img2*230\n        \n        cutoff1[:,:,1] = cutoff_img1*50\n        cutoff2[:,:,1] = cutoff_img2*50\n        cutoff1[:,:,2] = cutoff_img1*50\n        cutoff2[:,:,2] = cutoff_img2*50\n        \n        diff_photo1 = cutoff1.copy()\n        diff_photo2 = cutoff2.copy()\n        diff_photo1[:,:,1] += Y_true[k,:,:,0]*200\n        diff_photo2[:,:,1] += Y_true[k,:,:,0]*200\n        \n        ax[k, 0].imshow(X[k])\n        ax[k, 1].imshow(true_img, cmap = \"gray\")\n        ax[k, 2].imshow(cutoff1, cmap = \"gray\")\n        ax[k, 3].imshow(diff_photo1)\n        ax[k, 4].imshow(cutoff2, cmap = \"gray\")\n        ax[k, 5].imshow(diff_photo2)\n        \n        for j in range(6):\n            ax[k,j].set_xticks([])\n            ax[k,j].set_yticks([])\n    \n        if k == 0:\n            ax[k, 0].set_title(\"val img\")\n            ax[k, 1].set_title(\"true label\")\n            ax[k, 2].set_title(\"FCN (cutoff at 0.8)\")\n            ax[k, 3].set_title(\"Compare (Y:tp)\")\n            ax[k, 4].set_title(\"UNET (cutoff at 0.8)\")\n            ax[k, 5].set_title(\"Compare (Y:tp)\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.795856Z","iopub.execute_input":"2023-10-02T16:35:11.796430Z","iopub.status.idle":"2023-10-02T16:35:11.809215Z","shell.execute_reply.started":"2023-10-02T16:35:11.796396Z","shell.execute_reply":"2023-10-02T16:35:11.808206Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.random.seed(5)\nval_sample = np.random.choice(n_val, 10)\nval_sample","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.810933Z","iopub.execute_input":"2023-10-02T16:35:11.811529Z","iopub.status.idle":"2023-10-02T16:35:11.831115Z","shell.execute_reply.started":"2023-10-02T16:35:11.811490Z","shell.execute_reply":"2023-10-02T16:35:11.829971Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Prediction by FCN**","metadata":{}},{"cell_type":"code","source":"time1 = time.time()\n\n#prob_val_FCN = predict_probability(X_val[val_sample], model_FCN)\nprob_val_FCN = predict_probability(X_val, model_FCN)\ntime2 = time.time()\n\nval_time_FCN = np.round((time2 - time1), 2)\nprint(val_time_FCN, \"sec\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:11.832233Z","iopub.execute_input":"2023-10-02T16:35:11.833024Z","iopub.status.idle":"2023-10-02T16:35:40.130501Z","shell.execute_reply.started":"2023-10-02T16:35:11.832998Z","shell.execute_reply":"2023-10-02T16:35:40.129539Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Prediction by U-NET**","metadata":{}},{"cell_type":"code","source":"time1 = time.time()\n#prob_val_UNET = predict_probability(X_val[val_sample], model_UNET)\nprob_val_UNET = predict_probability(X_val, model_UNET)\ntime2 = time.time()\n\nval_time_UNET = np.round((time2 - time1), 2)\nprint(val_time_UNET, \"sec\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:35:40.132172Z","iopub.execute_input":"2023-10-02T16:35:40.132823Z","iopub.status.idle":"2023-10-02T16:36:16.501446Z","shell.execute_reply.started":"2023-10-02T16:35:40.132789Z","shell.execute_reply":"2023-10-02T16:36:16.500406Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Results**\n\nTen of prediction resutls are shown below. Apparently, both models predict somewhat well, although it is not simple task. The true label is **green** and the predicted label is **red**. Predicted label is created by cutting off predicted probability at 0.8. When both are plotted in the same image (refer to **Compare (Y:tp)**), true positive becomes **yellow**. Precisions of prediction are varied, some are very well predicted but some are not. The next question would be \"Does the result depend on dataset?\" ","metadata":{}},{"cell_type":"code","source":"#plot_result(X_val[val_sample], Y_val[val_sample], prob_val_FCN, prob_val_UNET)\nplot_result(X_val[val_sample], Y_val[val_sample], prob_val_FCN[val_sample], prob_val_UNET[val_sample])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:16.507852Z","iopub.execute_input":"2023-10-02T16:36:16.508149Z","iopub.status.idle":"2023-10-02T16:36:21.197084Z","shell.execute_reply.started":"2023-10-02T16:36:16.508126Z","shell.execute_reply":"2023-10-02T16:36:21.196261Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 7.3 Average IoU\n\nFinally, the average IoU is calculated: predict label, calculate IoU of an image, then take mean of them. Overall, both **FCN and U-NET achieved the target: Average IoU > 0.3 in validation**. When it is calculated by datasets/sources, however, both models did not achieve 0.5 in dataset1 source2. The reason could be lack of quantity of dataset1 source2 (refer to the chapter 4 EDA, section 4.1 Basic Information). To figure out the relationship between number of images and validation scores, scatterplot `Average IoU vs train image quantity` are created. Those are kind of correlated.  ","metadata":{}},{"cell_type":"code","source":"#IoU_FCN_list = []\n#IoU_UNET_list = []\n#for i in range(n_val):\n#    val1 = cal_IoU(Y_val[i], prob_val_FCN[i], cutoff = 0.8)\n#    val2 = cal_IoU(Y_val[i], prob_val_UNET[i], cutoff = 0.8)\n#    IoU_FCN_list.append(val1)\n#    IoU_UNET_list.append(val2)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:21.198235Z","iopub.execute_input":"2023-10-02T16:36:21.199070Z","iopub.status.idle":"2023-10-02T16:36:21.203529Z","shell.execute_reply.started":"2023-10-02T16:36:21.199034Z","shell.execute_reply":"2023-10-02T16:36:21.202778Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"IoU_FCN_list = cal_IoU(Y_val, prob_val_FCN, cutoff = 0.8)\nIoU_UNET_list = cal_IoU(Y_val, prob_val_UNET, cutoff = 0.8)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:21.204574Z","iopub.execute_input":"2023-10-02T16:36:21.205531Z","iopub.status.idle":"2023-10-02T16:36:22.009378Z","shell.execute_reply.started":"2023-10-02T16:36:21.205501Z","shell.execute_reply":"2023-10-02T16:36:22.008320Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"val_tile_df[\"IoU_FCN\"] = IoU_FCN_list\nval_tile_df[\"IoU_UNET\"] = IoU_UNET_list","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.010767Z","iopub.execute_input":"2023-10-02T16:36:22.011769Z","iopub.status.idle":"2023-10-02T16:36:22.018581Z","shell.execute_reply.started":"2023-10-02T16:36:22.011701Z","shell.execute_reply":"2023-10-02T16:36:22.017584Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Average IoU**","metadata":{}},{"cell_type":"code","source":"val_tile_df[[\"IoU_FCN\", \"IoU_UNET\"]].mean()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.020204Z","iopub.execute_input":"2023-10-02T16:36:22.020962Z","iopub.status.idle":"2023-10-02T16:36:22.041473Z","shell.execute_reply.started":"2023-10-02T16:36:22.020926Z","shell.execute_reply":"2023-10-02T16:36:22.040241Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Average IoU by dataset/source**","metadata":{}},{"cell_type":"code","source":"mean_IoU_df = val_tile_df[[\"dataset\", \"source_wsi\", \"IoU_FCN\", \"IoU_UNET\"]].groupby([\"dataset\", \"source_wsi\"]).mean()\nmean_IoU_df","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.043075Z","iopub.execute_input":"2023-10-02T16:36:22.043793Z","iopub.status.idle":"2023-10-02T16:36:22.063013Z","shell.execute_reply.started":"2023-10-02T16:36:22.043755Z","shell.execute_reply":"2023-10-02T16:36:22.061822Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"select = [\"dataset\", \"source_wsi\", \"count\"]\nsum_train_df = tile_df[select].iloc[train_idx].groupby(select[0:-1]).sum().astype(int)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.064590Z","iopub.execute_input":"2023-10-02T16:36:22.065635Z","iopub.status.idle":"2023-10-02T16:36:22.080833Z","shell.execute_reply.started":"2023-10-02T16:36:22.065589Z","shell.execute_reply":"2023-10-02T16:36:22.079666Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sum_train_df = sum_train_df.reset_index()\nsum_train_df[\"text\"] = sum_train_df.apply(lambda x: \"d\" +str(x[\"dataset\"]) +\", s\" + str(x[\"source_wsi\"]) , axis = 1) \nsum_train_df = sum_train_df.merge(mean_IoU_df, on = [\"dataset\", \"source_wsi\"])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.082502Z","iopub.execute_input":"2023-10-02T16:36:22.083173Z","iopub.status.idle":"2023-10-02T16:36:22.099402Z","shell.execute_reply.started":"2023-10-02T16:36:22.083132Z","shell.execute_reply":"2023-10-02T16:36:22.098169Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig, ax = plt.subplots(figsize = (5, 4))\nsns.scatterplot(data =sum_train_df, x = \"count\", y = \"IoU_FCN\", ax =ax, label = \"FCN\")\nsns.scatterplot(data =sum_train_df, x = \"count\", y = \"IoU_UNET\", ax = ax, label = \"U-NET\")\n\nfor i in range(sum_train_df.shape[0]):\n    ax.text(x = sum_train_df[\"count\"].iloc[i] + 12, y = sum_train_df[\"IoU_FCN\"].iloc[i] + 0.04, s= sum_train_df[\"text\"].iloc[i])\nax.set_xlabel(\"number of train images in dataset/source\")\nax.set_ylabel(\"Average IoU\")\nax.grid()\nax.set_xlim(0, 400)\nax.set_ylim(0, 0.7)\nax.set_title(\"Average IoU vs number of train image\")\nax.legend(loc = 4)\nfig.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.101005Z","iopub.execute_input":"2023-10-02T16:36:22.101930Z","iopub.status.idle":"2023-10-02T16:36:22.535151Z","shell.execute_reply.started":"2023-10-02T16:36:22.101838Z","shell.execute_reply":"2023-10-02T16:36:22.534098Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 7.4 IoU Distribution\n\nThose density plots shows IoU Distribution. First one is overall distribution, and latter ones are by dataset/source. For both models, They are normally distributed.","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(figsize =(5, 4))\nsns.kdeplot(val_tile_df[\"IoU_FCN\"], ax = ax, label = \"FCN\")\nsns.kdeplot(val_tile_df[\"IoU_UNET\"], ax = ax, label = \"U-NET\")\nax.grid()\nax.legend()\nax.axvline(x = 0.3, ls = \"-.\", color = \"black\", alpha = 0.5)\nax.set_title(\"Validation IoU distribution\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.536517Z","iopub.execute_input":"2023-10-02T16:36:22.537473Z","iopub.status.idle":"2023-10-02T16:36:22.964607Z","shell.execute_reply.started":"2023-10-02T16:36:22.537435Z","shell.execute_reply":"2023-10-02T16:36:22.963764Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def val_density_plot(dataset, n_source):\n    \n    fig, ax = plt.subplots(1,n_source, figsize = (3.5*n_source,2.7))\n    \n    for i in range(n_source):\n        filter1 = (val_tile_df[\"dataset\"] == dataset) & (val_tile_df[\"source_wsi\"] == i + 1)\n        sns.kdeplot(val_tile_df[\"IoU_FCN\"].loc[filter1], label = \"FCN IoU\", ax = ax[i])\n        sns.kdeplot(val_tile_df[\"IoU_UNET\"].loc[filter1], label = \"UNET IoU\", ax = ax[i])\n        ax[i].set_xlim(0, 1)\n        ax[i].set_title(\"(dataset, source) = (\" +str(dataset) + \", \" + str(i+1)  + \")\" )\n        ax[i].legend()\n        ax[i].set_xlabel(\"Val IoU\")\n        ax[i].axvline(x = 0.3, ls = \"-.\", color = \"black\", alpha = 0.5)\n        if i > 0:\n            ax[i].set_ylabel(\"\")\n            \n            \n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.965824Z","iopub.execute_input":"2023-10-02T16:36:22.967061Z","iopub.status.idle":"2023-10-02T16:36:22.975320Z","shell.execute_reply.started":"2023-10-02T16:36:22.966993Z","shell.execute_reply":"2023-10-02T16:36:22.974317Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"for dataset1/source1, U-NET has peak greater than 0.6.","metadata":{}},{"cell_type":"code","source":"val_density_plot(1, 2)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:22.976621Z","iopub.execute_input":"2023-10-02T16:36:22.977581Z","iopub.status.idle":"2023-10-02T16:36:23.553272Z","shell.execute_reply.started":"2023-10-02T16:36:22.977547Z","shell.execute_reply":"2023-10-02T16:36:23.552289Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"For dataset2/source1, both models have peak greather than 0.6.","metadata":{}},{"cell_type":"code","source":"val_density_plot(2, 4)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T16:36:23.554802Z","iopub.execute_input":"2023-10-02T16:36:23.555853Z","iopub.status.idle":"2023-10-02T16:36:24.590809Z","shell.execute_reply.started":"2023-10-02T16:36:23.555817Z","shell.execute_reply":"2023-10-02T16:36:24.589821Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 8. Kaggle Score\n\n## 8.1 Evaluation Metric\n\nThis competition adopted different metric for scoring. \n\n*Submissions are evaluated by computing the **Average Precision over confidence** scores. ... Segmentation is calculated using **IoU with a threshold of 0.6**.*\n\nhttps://www.kaggle.com/competitions/hubmap-hacking-the-human-vasculature/overview/evaluation\n\nUnlike Average IoU, this metric calculate **true positive and false positive for each blood vessel**. **True positive of a blood vessel means IoU > 0.6**. For example, when there are multiple blood vessels closely exist and predicted connected labels are connected, then it counted as false positive. ","metadata":{}},{"cell_type":"markdown","source":"## 8.2 How I submitted\n\nUsing U-NET or FCN accompanies following additional procedures (**step1 - step3**). Those who applied Mask-RCNN or YOLO does not need step1 - step3. \n\nStep1: To cutoff low probability in predicted mask at threshold 0.8~0.84 (tuned before every submission by IoU). \n\nStep2: To identify all predicted blood vessels indivisually by `cv2.connectedComponents`\n\nStep3: To calculate **mean probability of each component**. It is confidence of prediction required by the competition. \n\nStep4: To encode every predicted label using the provided function [`encode_binary_mask`](https://www.kaggle.com/competitions/hubmap-hacking-the-human-vasculature/overview/evaluation)\n\nStep5: To submit the encoded labels with confidences.","metadata":{}},{"cell_type":"markdown","source":"## 8.3 Kaggle Score\n\nBoth models are trained by 95% of given image, and max_epoch is 35, and the models with lowest validation loss are selected.\nThen kaggle private scores are: \n\n* FCN: 0.264\n* U-NET: 0.255","metadata":{}},{"cell_type":"markdown","source":"# 9. Conclusion\n\nThis project achieved avrage IoU in validation > 0.3 by both FCN and U-NET. First, **EDA discovered imbalance in label** and difference of distribution between datasets/source. Based on EDA, **custom loss function is developed to increase importance of positive labels** which are small in images. **Image augmentation consists horizontal&vertical flip, rotation, and random zoom**, because those operations do not collapse the context (there must be no definitive direction of tissues). It was effective to privent overfitting when there are no plenty of images (only 1600 images are given). \n\nBoth FCN and U-NET took nearly similar computation time in validation: 25~45sec. However, training of FCN was much faster than U-NET. FCN took one epoch for 100sec, whereas U-NET took 210sec, because architectre of FCN is much simpler and number of parameters are be smaller. On the otherhand, U-NET got smaller validation loss at every epoch. Thus, the best architecture depends on circumstances.\n\n**What did not work well**\n\nFollowing attemps did not improve results:\n\n* random brightness change \n* transfer learning from VGG16 trained with imagenet dataset\n\nRandom brightness change is one of typical auguentation, but merely worsened both validation and kaggle score. The reason should be small variation of brightness in dataset.\n\n**Ideas for improvement**\n\nChapter 7.3 revealed that validation score was worst in dataset1/source2, likely due to its small number of images. Since four images were randomly selected for a batch of training input (batch size = 4), dataset1/source2 had relatively less chance to be selected. In this case, **selecting images by equal probability from datasets/source** rectifies this imbalance of training that could improve precision. \n\n","metadata":{}},{"cell_type":"markdown","source":"# References\n\n[1] Fully Convolutional Networks for Semantic Segmentation, Jonathan Long, Evan Shelhamer, Trevor Darrell, UC Berkeley: https://arxiv.org/abs/1411.4038\n\n[2] U-Net: Convolutional Networks for Biomedical Image Segmentation, Olaf Ronneberger, Philipp Fischer, Thomas Brox, Computer Science Department and BIOSS Centre for Biological Signalling Studies,\nUniversity of Freiburg, Germany: https://arxiv.org/abs/1505.04597","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}