{"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":"<div style=\"height:200px;width:100%;margin: 0;\">\n    <img src=\"https://storage.googleapis.com/kaggle-competitions/kaggle/52279/logos/header.png?t=2023-05-04-23-04-20\" style=\"width:100%;\" />\n</div>","metadata":{}},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"hubmap\"><center>HuBMAP - Vasculature</center></h3>","metadata":{}},{"cell_type":"markdown","source":" The proper functioning of your body's organs and tissues depends on the interaction, spatial organization, and specialization of your cells—all 37 trillion of them. With so many cells, determining their functions and relationships is a monumental undertaking.\n\nCurrent efforts to map cells involve the Vasculature Common Coordinate Framework (VCCF), which uses the blood vasculature in the human body as the primary navigation system. The VCCF crosses all scale levels--from the whole body to the single cell level--and provides a unique way to identify cellular locations using capillary structures as an address. However, the gaps in what researchers know about microvasculature lead to gaps in the VCCF. If we could automatically segment microvasculature arrangements, researchers could use the real-world tissue data to begin to fill in those gaps and map out the vasculature.\n\nCompetition host Human BioMolecular Atlas Program (HuBMAP) hopes to develop an open and global platform to map healthy cells in the human body. Using the latest molecular and cellular biology technologies, HuBMAP researchers are studying the connections that cells have with each other throughout the body.\n\n____\n\n\n##### HuBMAP Base\n[iafoss - HuBMAP Starter](https://www.kaggle.com/code/iafoss/hubmap-pytorch-fast-ai-starter) \n<br>\n##### Dataset\n[Dataset Segmentations 768x768](https://www.kaggle.com/alincijov/dataset-hubmap-vasc-768x768-segments)\n<br>\n##### Training\n[Training - FastAI Baseline](https://www.kaggle.com/alincijov/training-hubmap-vasc-fastai-baseline)\n<br>\n##### Blood and sweat of the Devastator for HuBMAP Body Competition\n[The Devastator](https://www.kaggle.com/thedevastator)\n____","metadata":{}},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"pycocotools\"><center>Pycocotools</center></h3>","metadata":{}},{"cell_type":"code","source":"# Install pycocotools package\nimport os\n!mkdir /kaggle/working/packages\n!cp -r /kaggle/input/pycocotools/* /kaggle/working/packages\nos.chdir(\"/kaggle/working/packages/pycocotools-2.0.6/\")\n!python setup.py install\n!pip install . --no-index --find-links /kaggle/working/packages/\nos.chdir(\"/kaggle/working\")","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:00:15.588694Z","iopub.execute_input":"2023-05-31T22:00:15.589078Z","iopub.status.idle":"2023-05-31T22:00:59.111466Z","shell.execute_reply.started":"2023-05-31T22:00:15.589047Z","shell.execute_reply":"2023-05-31T22:00:59.110185Z"},"_kg_hide-output":true,"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"config\"><center>Configuration</center></h3>","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport cv2\n\nimport base64\nimport numpy as np\nfrom pycocotools import _mask as coco_mask\nimport typing as t\nimport zlib\n\nimport os\nimport gc\nimport cv2\nimport rasterio\nimport numpy as np\nimport pandas as pd\nfrom PIL import Image\nfrom tqdm import tqdm\nimport tifffile as tiff\nimport matplotlib.pyplot as plt\nfrom fastai.vision.all import *\nfrom rasterio.windows import Window\nfrom torch.utils.data import Dataset, DataLoader\nimport warnings; warnings.filterwarnings(\"ignore\")","metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2023-05-31T22:01:36.478037Z","iopub.execute_input":"2023-05-31T22:01:36.478491Z","iopub.status.idle":"2023-05-31T22:01:44.061724Z","shell.execute_reply.started":"2023-05-31T22:01:36.478453Z","shell.execute_reply":"2023-05-31T22:01:44.060568Z"},"jupyter":{"source_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"bs = 32\nsz = 768    # the size of tiles\nreduce = 4  # reduce the original images by 4 times\nTH = 0.0  # threshold for positive predictions\nDATA = '/kaggle/input/hubmap-hacking-the-human-vasculature/test/'\nMODELS = [f'/kaggle/input/hubmap-vasc-fastai/model_{i}.pth' for i in range(3)]\ndf_sample = pd.read_csv('/kaggle/input/hubmap-hacking-the-human-vasculature/sample_submission.csv')\ndevice = torch.device('cuda' if torch.cuda.is_available() else 'cpu')","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:02:35.853276Z","iopub.execute_input":"2023-05-31T22:02:35.853760Z","iopub.status.idle":"2023-05-31T22:02:35.902698Z","shell.execute_reply.started":"2023-05-31T22:02:35.853718Z","shell.execute_reply":"2023-05-31T22:02:35.901681Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"additionals\"><center>Additionals</center></h3>","metadata":{}},{"cell_type":"code","source":"\"\"\"\nLovasz-Softmax and Jaccard hinge loss in PyTorch\nMaxim Berman 2018 ESAT-PSI KU Leuven (MIT License)\n\"\"\"\n\nimport torch\nfrom torch.autograd import Variable\nimport torch.nn.functional as F\nimport numpy as np\ntry:\n    from itertools import  ifilterfalse\nexcept ImportError: # py3k\n    from itertools import  filterfalse\n\n\ndef lovasz_grad(gt_sorted):\n    \"\"\"\n    Computes gradient of the Lovasz extension w.r.t sorted errors\n    See Alg. 1 in paper\n    \"\"\"\n    p = len(gt_sorted)\n    gts = gt_sorted.sum()\n    intersection = gts - gt_sorted.float().cumsum(0)\n    union = gts + (1 - gt_sorted).float().cumsum(0)\n    jaccard = 1. - intersection / union\n    if p > 1: # cover 1-pixel case\n        jaccard[1:p] = jaccard[1:p] - jaccard[0:-1]\n    return jaccard\n\n\ndef iou_binary(preds, labels, EMPTY=1., ignore=None, per_image=True):\n    \"\"\"\n    IoU for foreground class\n    binary: 1 foreground, 0 background\n    \"\"\"\n    if not per_image:\n        preds, labels = (preds,), (labels,)\n    ious = []\n    for pred, label in zip(preds, labels):\n        intersection = ((label == 1) & (pred == 1)).sum()\n        union = ((label == 1) | ((pred == 1) & (label != ignore))).sum()\n        if not union:\n            iou = EMPTY\n        else:\n            iou = float(intersection) / union\n        ious.append(iou)\n    iou = f_mean(ious)    # mean accross images if per_image\n    return 100 * iou\n\n\ndef iou(preds, labels, C, EMPTY=1., ignore=None, per_image=False):\n    \"\"\"\n    Array of IoU for each (non ignored) class\n    \"\"\"\n    if not per_image:\n        preds, labels = (preds,), (labels,)\n    ious = []\n    for pred, label in zip(preds, labels):\n        iou = []    \n        for i in range(C):\n            if i != ignore: # The ignored label is sometimes among predicted classes (ENet - CityScapes)\n                intersection = ((label == i) & (pred == i)).sum()\n                union = ((label == i) | ((pred == i) & (label != ignore))).sum()\n                if not union:\n                    iou.append(EMPTY)\n                else:\n                    iou.append(float(intersection) / union)\n        ious.append(iou)\n    ious = map(f_mean, zip(*ious)) # mean accross images if per_image\n    return 100 * np.array(ious)\n\n\n# --------------------------- BINARY LOSSES ---------------------------\n\n\ndef lovasz_hinge(logits, labels, per_image=True, ignore=None):\n    \"\"\"\n    Binary Lovasz hinge loss\n      logits: [B, H, W] Variable, logits at each pixel (between -\\infty and +\\infty)\n      labels: [B, H, W] Tensor, binary ground truth masks (0 or 1)\n      per_image: compute the loss per image instead of per batch\n      ignore: void class id\n    \"\"\"\n    if per_image:\n        loss = f_mean(lovasz_hinge_flat(*flatten_binary_scores(log.unsqueeze(0), lab.unsqueeze(0), ignore))\n                          for log, lab in zip(logits, labels))\n    else:\n        loss = lovasz_hinge_flat(*flatten_binary_scores(logits, labels, ignore))\n    return loss\n\n\ndef lovasz_hinge_flat(logits, labels):\n    \"\"\"\n    Binary Lovasz hinge loss\n      logits: [P] Variable, logits at each prediction (between -\\infty and +\\infty)\n      labels: [P] Tensor, binary ground truth labels (0 or 1)\n      ignore: label to ignore\n    \"\"\"\n    if len(labels) == 0:\n        # only void pixels, the gradients should be 0\n        return logits.sum() * 0.\n    signs = 2. * labels.float() - 1.\n    errors = (1. - logits * Variable(signs))\n    errors_sorted, perm = torch.sort(errors, dim=0, descending=True)\n    perm = perm.data\n    gt_sorted = labels[perm]\n    grad = lovasz_grad(gt_sorted)\n    #loss = torch.dot(F.relu(errors_sorted), Variable(grad))\n    loss = torch.dot(F.elu(errors_sorted)+1, Variable(grad))\n    return loss\n\n\ndef flatten_binary_scores(scores, labels, ignore=None):\n    \"\"\"\n    Flattens predictions in the batch (binary case)\n    Remove labels equal to 'ignore'\n    \"\"\"\n    scores = scores.view(-1)\n    labels = labels.view(-1)\n    if ignore is None:\n        return scores, labels\n    valid = (labels != ignore)\n    vscores = scores[valid]\n    vlabels = labels[valid]\n    return vscores, vlabels\n\n\nclass StableBCELoss(torch.nn.modules.Module):\n    def __init__(self):\n         super(StableBCELoss, self).__init__()\n    def forward(self, input, target):\n         neg_abs = - input.abs()\n         loss = input.clamp(min=0) - input * target + (1 + neg_abs.exp()).log()\n         return loss.mean()\n\n\ndef binary_xloss(logits, labels, ignore=None):\n    \"\"\"\n    Binary Cross entropy loss\n      logits: [B, H, W] Variable, logits at each pixel (between -\\infty and +\\infty)\n      labels: [B, H, W] Tensor, binary ground truth masks (0 or 1)\n      ignore: void class id\n    \"\"\"\n    logits, labels = flatten_binary_scores(logits, labels, ignore)\n    loss = StableBCELoss()(logits, Variable(labels.float()))\n    return loss\n\n\n# --------------------------- MULTICLASS LOSSES ---------------------------\n\n\ndef lovasz_softmax(probas, labels, only_present=False, per_image=False, ignore=None):\n    \"\"\"\n    Multi-class Lovasz-Softmax loss\n      probas: [B, C, H, W] Variable, class probabilities at each prediction (between 0 and 1)\n      labels: [B, H, W] Tensor, ground truth labels (between 0 and C - 1)\n      only_present: average only on classes present in ground truth\n      per_image: compute the loss per image instead of per batch\n      ignore: void class labels\n    \"\"\"\n    if per_image:\n        loss = f_mean(lovasz_softmax_flat(*flatten_probas(prob.unsqueeze(0), lab.unsqueeze(0), ignore), only_present=only_present)\n                          for prob, lab in zip(probas, labels))\n    else:\n        loss = lovasz_softmax_flat(*flatten_probas(probas, labels, ignore), only_present=only_present)\n    return loss\n\n\ndef lovasz_softmax_flat(probas, labels, only_present=False):\n    \"\"\"\n    Multi-class Lovasz-Softmax loss\n      probas: [P, C] Variable, class probabilities at each prediction (between 0 and 1)\n      labels: [P] Tensor, ground truth labels (between 0 and C - 1)\n      only_present: average only on classes present in ground truth\n    \"\"\"\n    C = probas.size(1)\n    losses = []\n    for c in range(C):\n        fg = (labels == c).float() # foreground for class c\n        if only_present and fg.sum() == 0:\n            continue\n        errors = (Variable(fg) - probas[:, c]).abs()\n        errors_sorted, perm = torch.sort(errors, 0, descending=True)\n        perm = perm.data\n        fg_sorted = fg[perm]\n        losses.append(torch.dot(errors_sorted, Variable(lovasz_grad(fg_sorted))))\n    return f_mean(losses)\n\n\ndef flatten_probas(probas, labels, ignore=None):\n    \"\"\"\n    Flattens predictions in the batch\n    \"\"\"\n    B, C, H, W = probas.size()\n    probas = probas.permute(0, 2, 3, 1).contiguous().view(-1, C)  # B * H * W, C = P, C\n    labels = labels.view(-1)\n    if ignore is None:\n        return probas, labels\n    valid = (labels != ignore)\n    vprobas = probas[valid.nonzero().squeeze()]\n    vlabels = labels[valid]\n    return vprobas, vlabels\n\ndef xloss(logits, labels, ignore=None):\n    \"\"\"\n    Cross entropy loss\n    \"\"\"\n    return F.cross_entropy(logits, Variable(labels), ignore_index=255)\n\n\n# --------------------------- HELPER FUNCTIONS ---------------------------\n\ndef f_mean(l, ignore_nan=False, empty=0):\n    \"\"\"\n    nanmean compatible with generators.\n    \"\"\"\n    l = iter(l)\n    if ignore_nan:\n        l = ifilterfalse(np.isnan, l)\n    try:\n        n = 1\n        acc = next(l)\n    except StopIteration:\n        if empty == 'raise':\n            raise ValueError('Empty mean')\n        return empty\n    for n, v in enumerate(l, 2):\n        acc += v\n    if n == 1:\n        return acc\n    return acc / n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:02:39.248697Z","iopub.execute_input":"2023-05-31T22:02:39.249130Z","iopub.status.idle":"2023-05-31T22:02:39.305204Z","shell.execute_reply.started":"2023-05-31T22:02:39.249094Z","shell.execute_reply":"2023-05-31T22:02:39.303467Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"dataset\"><center>Dataset</center></h3>","metadata":{}},{"cell_type":"code","source":"# functions to convert encoding to mask and mask to encoding\ndef enc2mask(encs, shape):\n    img = np.zeros(shape[0]*shape[1], dtype=np.uint8)\n    for m,enc in enumerate(encs):\n        if isinstance(enc,np.float) and np.isnan(enc): continue\n        s = enc.split()\n        for i in range(len(s)//2):\n            start = int(s[2*i]) - 1\n            length = int(s[2*i+1])\n            img[start:start+length] = 1 + m\n    return img.reshape(shape).T\n\ndef mask2enc(mask, n=1):\n    pixels = mask.T.flatten()\n    encs = []\n    for i in range(1,n+1):\n        p = (pixels == i).astype(np.int8)\n        if p.sum() == 0: encs.append(np.nan)\n        else:\n            p = np.concatenate([[0], p, [0]])\n            runs = np.where(p[1:] != p[:-1])[0] + 1\n            runs[1::2] -= runs[::2]\n            encs.append(' '.join(str(x) for x in runs))\n    return encs\n\n#https://www.kaggle.com/bguberfain/memory-aware-rle-encoding\n#with transposed mask\ndef rle_encode_less_memory(img):\n    #the image should be transposed\n    pixels = img.T.flatten()\n    \n    # This simplified method requires first and last pixel to be zero\n    pixels[0] = 0\n    pixels[-1] = 0\n    runs = np.where(pixels[1:] != pixels[:-1])[0] + 2\n    runs[1::2] -= runs[::2]\n    \n    return ' '.join(str(x) for x in runs)","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:02:39.851635Z","iopub.execute_input":"2023-05-31T22:02:39.854800Z","iopub.status.idle":"2023-05-31T22:02:39.873810Z","shell.execute_reply.started":"2023-05-31T22:02:39.854751Z","shell.execute_reply":"2023-05-31T22:02:39.872670Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# https://www.kaggle.com/datasets/thedevastator/hubmap-2022-256x256\nmean = np.array([0.62405882, 0.40748174, 0.68217106])\nstd = np.array([0.15112928, 0.20741797, 0.13163469])\n\ns_th = 40  #saturation blancking threshold\np_th = 1000*(sz//256)**2 #threshold for the minimum number of pixels\nidentity = rasterio.Affine(1, 0, 0, 0, 1, 0)\n\ndef img2tensor(img,dtype:np.dtype=np.float32):\n    if img.ndim==2 : img = np.expand_dims(img,2)\n    img = np.transpose(img,(2,0,1))\n    return torch.from_numpy(img.astype(dtype, copy=False))\n\nclass HuBMAPDataset(Dataset):\n    def __init__(self, idx, sz=sz, reduce=reduce):\n        self.data = rasterio.open(os.path.join(DATA,idx+'.tif'), transform = identity,\n                                 num_threads='all_cpus')\n        # some images have issues with their format \n        # and must be saved correctly before reading with rasterio\n        if self.data.count != 3:\n            subdatasets = self.data.subdatasets\n            self.layers = []\n            if len(subdatasets) > 0:\n                for i, subdataset in enumerate(subdatasets, 0):\n                    self.layers.append(rasterio.open(subdataset))\n        self.shape = self.data.shape\n        self.reduce = reduce\n        self.sz = reduce*sz\n        self.pad0 = (self.sz - self.shape[0]%self.sz)%self.sz\n        self.pad1 = (self.sz - self.shape[1]%self.sz)%self.sz\n        self.n0max = (self.shape[0] + self.pad0)//self.sz\n        self.n1max = (self.shape[1] + self.pad1)//self.sz\n        \n    def __len__(self):\n        return self.n0max*self.n1max\n    \n    def __getitem__(self, idx):\n        # the code below may be a little bit difficult to understand,\n        # but the thing it does is mapping the original image to\n        # tiles created with adding padding, as done in\n        # https://www.kaggle.com/iafoss/256x256-images ,\n        # and then the tiles are loaded with rasterio\n        # n0,n1 - are the x and y index of the tile (idx = n0*self.n1max + n1)\n        n0,n1 = idx//self.n1max, idx%self.n1max\n        # x0,y0 - are the coordinates of the lower left corner of the tile in the image\n        # negative numbers correspond to padding (which must not be loaded)\n        x0,y0 = -self.pad0//2 + n0*self.sz, -self.pad1//2 + n1*self.sz\n        # make sure that the region to read is within the image\n        p00,p01 = max(0,x0), min(x0+self.sz,self.shape[0])\n        p10,p11 = max(0,y0), min(y0+self.sz,self.shape[1])\n        img = np.zeros((self.sz,self.sz,3),np.uint8)\n        # mapping the loade region to the tile\n        if self.data.count == 3:\n            img[(p00-x0):(p01-x0),(p10-y0):(p11-y0)] = np.moveaxis(self.data.read([1,2,3],\n                window=Window.from_slices((p00,p01),(p10,p11))), 0, -1)\n        else:\n            for i,layer in enumerate(self.layers):\n                img[(p00-x0):(p01-x0),(p10-y0):(p11-y0),i] =\\\n                  layer.read(1,window=Window.from_slices((p00,p01),(p10,p11)))\n        \n        if self.reduce != 1:\n            img = cv2.resize(img,(self.sz//reduce,self.sz//reduce),\n                             interpolation = cv2.INTER_AREA)\n        #check for empty imges\n        hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)\n        h,s,v = cv2.split(hsv)\n        if (s>s_th).sum() <= p_th or img.sum() <= p_th:\n            #images with -1 will be skipped\n            return img2tensor((img/255.0 - mean)/std), -1\n        else: return img2tensor((img/255.0 - mean)/std), idx","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:04:13.326797Z","iopub.execute_input":"2023-05-31T22:04:13.327254Z","iopub.status.idle":"2023-05-31T22:04:13.361145Z","shell.execute_reply.started":"2023-05-31T22:04:13.327215Z","shell.execute_reply":"2023-05-31T22:04:13.360214Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"preds\"><center>Predictions</center></h3>","metadata":{}},{"cell_type":"code","source":"#iterator like wrapper that returns predicted masks\nclass Model_pred:\n    def __init__(self, models, dl, tta:bool=True, half:bool=False):\n        self.models = models\n        self.dl = dl\n        self.tta = tta\n        self.half = half\n        \n    def __iter__(self):\n        count=0\n        with torch.no_grad():\n            for x,y in iter(self.dl):\n                if ((y>=0).sum() > 0): #exclude empty images\n                    x = x[y>=0].to(device)\n                    y = y[y>=0]\n                    if self.half: x = x.half()\n                    py = None\n                    for model in self.models:\n                        p = model(x)\n                        p = torch.sigmoid(p).detach()\n                        if py is None: py = p\n                        else: py += p\n                    if self.tta:\n                        #x,y,xy flips as TTA\n                        flips = [[-1],[-2],[-2,-1]]\n                        for f in flips:\n                            xf = torch.flip(x,f)\n                            for model in self.models:\n                                p = model(xf)\n                                p = torch.flip(p,f)\n                                py += torch.sigmoid(p).detach()\n                        py /= (1+len(flips))        \n                    py /= len(self.models)\n\n                    py = F.upsample(py, scale_factor=reduce, mode=\"bilinear\")\n                    py = py.permute(0,2,3,1).float().cpu()\n                    \n                    batch_size = len(py)\n                    for i in range(batch_size):\n                        yield py[i],y[i]\n                        count += 1\n                    \n    def __len__(self):\n        return len(self.dl.dataset)","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:04:14.239544Z","iopub.execute_input":"2023-05-31T22:04:14.240190Z","iopub.status.idle":"2023-05-31T22:04:14.261171Z","shell.execute_reply.started":"2023-05-31T22:04:14.240128Z","shell.execute_reply":"2023-05-31T22:04:14.259789Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"model\"><center>Model</center></h3>","metadata":{}},{"cell_type":"code","source":"class FPN(nn.Module):\n    def __init__(self, input_channels:list, output_channels:list):\n        super().__init__()\n        self.convs = nn.ModuleList(\n            [nn.Sequential(nn.Conv2d(in_ch, out_ch*2, kernel_size=3, padding=1),\n             nn.ReLU(inplace=True), nn.BatchNorm2d(out_ch*2),\n             nn.Conv2d(out_ch*2, out_ch, kernel_size=3, padding=1))\n            for in_ch, out_ch in zip(input_channels, output_channels)])\n        \n    def forward(self, xs:list, last_layer):\n        hcs = [F.interpolate(c(x),scale_factor=2**(len(self.convs)-i),mode='bilinear') \n               for i,(c,x) in enumerate(zip(self.convs, xs))]\n        hcs.append(last_layer)\n        return torch.cat(hcs, dim=1)\n\nclass UnetBlock(Module):\n    def __init__(self, up_in_c:int, x_in_c:int, nf:int=None, blur:bool=False,\n                 self_attention:bool=False, **kwargs):\n        super().__init__()\n        self.shuf = PixelShuffle_ICNR(up_in_c, up_in_c//2, blur=blur, **kwargs)\n        self.bn = nn.BatchNorm2d(x_in_c)\n        ni = up_in_c//2 + x_in_c\n        nf = nf if nf is not None else max(up_in_c//2,32)\n        self.conv1 = ConvLayer(ni, nf, norm_type=None, **kwargs)\n        self.conv2 = ConvLayer(nf, nf, norm_type=None,\n            xtra=SelfAttention(nf) if self_attention else None, **kwargs)\n        self.relu = nn.ReLU(inplace=True)\n\n    def forward(self, up_in:Tensor, left_in:Tensor) -> Tensor:\n        s = left_in\n        up_out = self.shuf(up_in)\n        cat_x = self.relu(torch.cat([up_out, self.bn(s)], dim=1))\n        return self.conv2(self.conv1(cat_x))\n        \nclass _ASPPModule(nn.Module):\n    def __init__(self, inplanes, planes, kernel_size, padding, dilation, groups=1):\n        super().__init__()\n        self.atrous_conv = nn.Conv2d(inplanes, planes, kernel_size=kernel_size,\n                stride=1, padding=padding, dilation=dilation, bias=False, groups=groups)\n        self.bn = nn.BatchNorm2d(planes)\n        self.relu = nn.ReLU()\n\n        self._init_weight()\n\n    def forward(self, x):\n        x = self.atrous_conv(x)\n        x = self.bn(x)\n\n        return self.relu(x)\n\n    def _init_weight(self):\n        for m in self.modules():\n            if isinstance(m, nn.Conv2d):\n                torch.nn.init.kaiming_normal_(m.weight)\n            elif isinstance(m, nn.BatchNorm2d):\n                m.weight.data.fill_(1)\n                m.bias.data.zero_()\n\nclass ASPP(nn.Module):\n    def __init__(self, inplanes=512, mid_c=256, dilations=[6, 12, 18, 24], out_c=None):\n        super().__init__()\n        self.aspps = [_ASPPModule(inplanes, mid_c, 1, padding=0, dilation=1)] + \\\n            [_ASPPModule(inplanes, mid_c, 3, padding=d, dilation=d,groups=4) for d in dilations]\n        self.aspps = nn.ModuleList(self.aspps)\n        self.global_pool = nn.Sequential(nn.AdaptiveMaxPool2d((1, 1)),\n                        nn.Conv2d(inplanes, mid_c, 1, stride=1, bias=False),\n                        nn.BatchNorm2d(mid_c), nn.ReLU())\n        out_c = out_c if out_c is not None else mid_c\n        self.out_conv = nn.Sequential(nn.Conv2d(mid_c*(2+len(dilations)), out_c, 1, bias=False),\n                                    nn.BatchNorm2d(out_c), nn.ReLU(inplace=True))\n        self.conv1 = nn.Conv2d(mid_c*(2+len(dilations)), out_c, 1, bias=False)\n        self._init_weight()\n\n    def forward(self, x):\n        x0 = self.global_pool(x)\n        xs = [aspp(x) for aspp in self.aspps]\n        x0 = F.interpolate(x0, size=xs[0].size()[2:], mode='bilinear', align_corners=True)\n        x = torch.cat([x0] + xs, dim=1)\n        return self.out_conv(x)\n    \n    def _init_weight(self):\n        for m in self.modules():\n            if isinstance(m, nn.Conv2d):\n                torch.nn.init.kaiming_normal_(m.weight)\n            elif isinstance(m, nn.BatchNorm2d):\n                m.weight.data.fill_(1)\n                m.bias.data.zero_()","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:04:16.097596Z","iopub.execute_input":"2023-05-31T22:04:16.098014Z","iopub.status.idle":"2023-05-31T22:04:16.138108Z","shell.execute_reply.started":"2023-05-31T22:04:16.097980Z","shell.execute_reply":"2023-05-31T22:04:16.136894Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from torchvision.models.resnet import ResNet, Bottleneck\nclass UneXt50(nn.Module):\n    def __init__(self, stride=1, **kwargs):\n        super().__init__()\n        #encoder\n        m = ResNet(Bottleneck, [3, 4, 6, 3], groups=32, width_per_group=4)\n        #m = torch.hub.load('facebookresearch/semi-supervised-ImageNet1K-models',\n        #                   'resnext50_32x4d_ssl')\n        self.enc0 = nn.Sequential(m.conv1, m.bn1, nn.ReLU(inplace=True))\n        self.enc1 = nn.Sequential(nn.MaxPool2d(kernel_size=3, stride=2, padding=1, dilation=1),\n                            m.layer1) #256\n        self.enc2 = m.layer2 #512\n        self.enc3 = m.layer3 #1024\n        self.enc4 = m.layer4 #2048\n        #aspp with customized dilatations\n        self.aspp = ASPP(2048,256,out_c=512,dilations=[stride*1,stride*2,stride*3,stride*4])\n        self.drop_aspp = nn.Dropout2d(0.5)\n        #decoder\n        self.dec4 = UnetBlock(512,1024,256)\n        self.dec3 = UnetBlock(256,512,128)\n        self.dec2 = UnetBlock(128,256,64)\n        self.dec1 = UnetBlock(64,64,32)\n        self.fpn = FPN([512,256,128,64],[16]*4)\n        self.drop = nn.Dropout2d(0.1)\n        self.final_conv = ConvLayer(32+16*4, 1, ks=1, norm_type=None, act_cls=None)\n        \n    def forward(self, x):\n        enc0 = self.enc0(x)\n        enc1 = self.enc1(enc0)\n        enc2 = self.enc2(enc1)\n        enc3 = self.enc3(enc2)\n        enc4 = self.enc4(enc3)\n        enc5 = self.aspp(enc4)\n        dec3 = self.dec4(self.drop_aspp(enc5),enc3)\n        dec2 = self.dec3(dec3,enc2)\n        dec1 = self.dec2(dec2,enc1)\n        dec0 = self.dec1(dec1,enc0)\n        x = self.fpn([enc5, dec3, dec2, dec1], dec0)\n        x = self.final_conv(self.drop(x))\n        x = F.interpolate(x,scale_factor=2,mode='bilinear')\n        return x","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:04:16.981438Z","iopub.execute_input":"2023-05-31T22:04:16.984670Z","iopub.status.idle":"2023-05-31T22:04:17.005408Z","shell.execute_reply.started":"2023-05-31T22:04:16.984623Z","shell.execute_reply":"2023-05-31T22:04:17.004456Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"models = []\nfor path in MODELS:\n    state_dict = torch.load(path,map_location=torch.device('cpu'))\n    model = UneXt50()\n    model.load_state_dict(state_dict)\n    model.float()\n    model.eval()\n    model.to(device)\n    models.append(model)\n\ndel state_dict\ngc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:03:02.404767Z","iopub.execute_input":"2023-05-31T22:03:02.405218Z","iopub.status.idle":"2023-05-31T22:03:15.242779Z","shell.execute_reply.started":"2023-05-31T22:03:02.405181Z","shell.execute_reply":"2023-05-31T22:03:15.241448Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"<h3 class=\"list-group-item list-group-item-action active\" data-toggle=\"list\" style=\"color:#c3448b; background:#efe9e9; border:1px dashed #efe50b;\" role=\"tab\" aria-controls=\"submission\"><center>Submission</center></h3>","metadata":{}},{"cell_type":"markdown","source":"## Encode Binary Mask","metadata":{}},{"cell_type":"code","source":"def encode_binary_mask(mask: np.ndarray) -> t.Text:\n    \"\"\"Converts a binary mask into OID challenge encoding ascii text.\"\"\"\n\n    # check input mask --\n    if mask.dtype != np.bool:\n        raise ValueError(\"encode_binary_mask expects a binary mask, received dtype == %s\" % mask.dtype)\n\n    mask = np.squeeze(mask)\n    if len(mask.shape) != 2:\n        raise ValueError(\"encode_binary_mask expects a 2d mask, received shape == %s\" % mask.shape)\n\n    # convert input mask to expected COCO API input --\n    mask_to_encode = mask.reshape(mask.shape[0], mask.shape[1], 1)\n    mask_to_encode = mask_to_encode.astype(np.uint8)\n    mask_to_encode = np.asfortranarray(mask_to_encode)\n\n    # RLE encode mask --\n    encoded_mask = coco_mask.encode(mask_to_encode)[0][\"counts\"]\n    \n    binary_str = zlib.compress(encoded_mask, zlib.Z_BEST_COMPRESSION)\n    base64_str = base64.b64encode(binary_str)\n    return base64_str","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:04:20.683710Z","iopub.execute_input":"2023-05-31T22:04:20.684178Z","iopub.status.idle":"2023-05-31T22:04:20.700813Z","shell.execute_reply.started":"2023-05-31T22:04:20.684122Z","shell.execute_reply":"2023-05-31T22:04:20.699255Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Predictions","metadata":{}},{"cell_type":"code","source":"test_path = \"/kaggle/input/hubmap-hacking-the-human-vasculature/test\"\nsample_submission = pd.read_csv('/kaggle/input/hubmap-hacking-the-human-vasculature/sample_submission.csv')\n\nids = []\nheights = []\nwidths = []\nprediction_strings = []\nfor img_name in os.listdir(test_path):\n    img = cv2.imread(f\"{test_path}/{img_name}\", 0)\n    idx = img_name.split('.')[0]\n    \n    ds = HuBMAPDataset(idx)\n    dl = DataLoader(ds,bs,num_workers=0,shuffle=False,pin_memory=True)\n    \n    h, w = ds.sz, ds.sz\n    \n    mp = Model_pred(models,dl)\n    \n    #generate masks\n    mask = torch.zeros(len(ds),ds.sz,ds.sz,dtype=torch.int8)\n    for p,i in iter(mp): mask[i.item()] = p.squeeze(-1) > TH\n    \n    #reshape tiled masks into a single mask and crop padding\n    mask = mask.view(ds.n0max,ds.n1max,ds.sz,ds.sz).\\\n        permute(0,2,1,3).reshape(ds.n0max*ds.sz,ds.n1max*ds.sz)\n    mask = mask[ds.pad0//2:-(ds.pad0-ds.pad0//2) if ds.pad0 > 0 else ds.n0max*ds.sz,\n        ds.pad1//2:-(ds.pad1-ds.pad1//2) if ds.pad1 > 0 else ds.n1max*ds.sz]\n    mask = mask.numpy().astype(np.bool)\n    \n    encoded = encode_binary_mask(mask)\n    \n    if(idx not in ids):\n        ids.append(idx)\n        heights.append(h)\n        widths.append(w)\n        prediction_strings.append(f\"0 1.0 {encoded.decode('utf-8')}\")\n    else:\n        new_prediction_string = f\" 0 1.0 {encoded.decode('utf-8')}\"\n        prediction_strings[ids.index(idx)] += prediction_strings[ids.index(idx)] + new_prediction_string","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:14:28.638808Z","iopub.execute_input":"2023-05-31T22:14:28.639276Z","iopub.status.idle":"2023-05-31T22:14:29.599617Z","shell.execute_reply.started":"2023-05-31T22:14:28.639237Z","shell.execute_reply":"2023-05-31T22:14:29.598583Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Submission","metadata":{}},{"cell_type":"code","source":"submission = pd.DataFrame()\nsubmission['id'] = ids\nsubmission['height'] = heights\nsubmission['width'] = widths\nsubmission['prediction_string'] = prediction_strings\n\n\nsubmission = submission.set_index('id')\nprint(submission)\nsubmission.to_csv(\"submission.csv\")","metadata":{"execution":{"iopub.status.busy":"2023-05-31T22:14:32.445758Z","iopub.execute_input":"2023-05-31T22:14:32.446235Z","iopub.status.idle":"2023-05-31T22:14:32.471116Z","shell.execute_reply.started":"2023-05-31T22:14:32.446195Z","shell.execute_reply":"2023-05-31T22:14:32.469960Z"},"trusted":true},"execution_count":null,"outputs":[]}]}