{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import cv2\nfrom tqdm import tqdm_notebook as tqdm\nimport fastai\nfrom fastai.vision import *\nimport os\nfrom mish_activation import *\nimport warnings\nwarnings.filterwarnings(\"ignore\")\nimport skimage.io\nimport numpy as np\nimport pandas as pd\nsys.path.insert(0, '../input/semisupervised-imagenet-models/semi-supervised-ImageNet1K-models-master/')\n#hubconf = '../input/semisupervised-imagenet-models/semi-supervised-ImageNet1K-models-master/hubconf.py'\nfrom hubconf import *\nfrom tqdm.notebook import tqdm\nfrom sklearn.utils import shuffle\nimport matplotlib.pyplot as plt\nfrom skimage.transform import resize","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install tifffile","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install ../input/installed-packages/imagecodecs-2020.5.30-cp37-cp37m-manylinux2014_x86_64.whl","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install ../input/installed-packages/spams-2.6.1-cp37-cp37m-linux_x86_64.whl","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import tifffile\nimport imagecodecs","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"TRAIN_DATA = '../input/prostate-cancer-grade-assessment/train_images/'\nTEST_DATA = '../input/prostate-cancer-grade-assessment/test_images/'\nTRAIN = '../input/prostate-cancer-grade-assessment/train.csv'\nTEST = '../input/prostate-cancer-grade-assessment/test.csv'\nSAMPLE = '../input/prostate-cancer-grade-assessment/sample_submission.csv'\nMODELS = [f'../input/panda-stain-norm-downsampling-12x64x64-models/RNXT50_{i}.pth' for i in range(4)]\n\nsz = 128\nbs = 1\nN = 12\nnworkers = 2\n\ndownsize = (N,64,64,3)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Model","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"def _resnext(url, block, layers, pretrained, progress, **kwargs):\n    model = ResNet(block, layers, **kwargs)\n    #state_dict = load_state_dict_from_url(url, progress=progress)\n    #model.load_state_dict(state_dict)\n    return model\n\nclass Model(nn.Module):\n    def __init__(self, arch='resnext50_32x4d', n=6, pre=True):\n        super().__init__()\n        #m = torch.hub.load('facebookresearch/semi-supervised-ImageNet1K-models', arch)\n        m = _resnext(semi_supervised_model_urls[arch], Bottleneck, [3, 4, 6, 3], False, \n                progress=False,groups=32,width_per_group=4)\n        self.enc = nn.Sequential(*list(m.children())[:-2])       \n        nc = list(m.children())[-1].in_features\n        self.head = nn.Sequential(AdaptiveConcatPool2d(),Flatten(),nn.Linear(2*nc,512),\n                Mish(),nn.BatchNorm1d(512),nn.Dropout(0.5),nn.Linear(512,n))\n        \n    def forward(self, x):\n        shape = x.shape\n        n = shape[1]\n        x = x.view(-1,shape[2],shape[3],shape[4])\n        x = self.enc(x)\n        shape = x.shape\n        x = x.view(-1,n,shape[1],shape[2],shape[3]).permute(0,2,1,3,4).contiguous()\\\n          .view(-1,shape[1],shape[2]*n,shape[3])\n        x = self.head(x)\n        x = F.softmax(x,dim=1)\n        return x","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"models = []\nfor path in MODELS:\n    state_dict = torch.load(path,map_location=torch.device('cpu'))\n    model = Model()\n    model.load_state_dict(state_dict)\n    model.float()\n    model.eval()\n    model.cuda()\n    models.append(model)\n\ndel state_dict","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Stain Normalization","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"# STAIN NORMALIZATION FUNCTIONS\nimport spams\nclass TissueMaskException(Exception):\n    pass\n\n######################################################################################################\n\ndef is_uint8_image(I):\n\n    if not is_image(I):\n        return False\n    if I.dtype != np.uint8:\n        return False\n    return True\n######################################################################################################\n\ndef is_image(I):\n\n    if not isinstance(I, np.ndarray):\n        return False\n    if not I.ndim == 3:\n        return False\n    return True\n######################################################################################################\n\ndef get_tissue_mask(I, luminosity_threshold=0.8):\n\n    I_LAB = cv2.cvtColor(I, cv2.COLOR_RGB2LAB)\n    L = I_LAB[:, :, 0] / 255.0  # Convert to range [0,1].\n    mask = L < luminosity_threshold\n\n    # Check it's not empty\n    if mask.sum() == 0:\n        raise TissueMaskException(\"Empty tissue mask computed\")\n\n    return mask\n\n######################################################################################################\n\ndef convert_RGB_to_OD(I):\n\n    mask = (I == 0)\n    I[mask] = 1\n    \n\n    #return np.maximum(-1 * np.log(I / 255), 1e-6)\n    return np.maximum(-1 * np.log(I / 255), np.zeros(I.shape) + 0.1)\n\n######################################################################################################\n\ndef convert_OD_to_RGB(OD):\n\n    assert OD.min() >= 0, \"Negative optical density.\"\n    \n    OD = np.maximum(OD, 1e-6)\n    \n    return (255 * np.exp(-1 * OD)).astype(np.uint8)\n\n######################################################################################################\n\ndef normalize_matrix_rows(A):\n\n    return A / np.linalg.norm(A, axis=1)[:, None]\n\n######################################################################################################\n\n\ndef get_concentrations(I, stain_matrix, regularizer=0.01):\n\n    OD = convert_RGB_to_OD(I).reshape((-1, 3))\n    return spams.lasso(X=OD.T, D=stain_matrix.T, mode=2, lambda1=regularizer, pos=True).toarray().T\n\n######################################################################################################\n\ndef get_stain_matrix(I, luminosity_threshold=0.8, angular_percentile=99):\n    \n    #assert is_uint8_image(I), \"Image should be RGB uint8.\"\n    # Convert to OD and ignore background\n    tissue_mask = get_tissue_mask(I, luminosity_threshold=luminosity_threshold).reshape((-1,))\n    OD = convert_RGB_to_OD(I).reshape((-1, 3))\n    \n    OD = OD[tissue_mask]\n\n    # Eigenvectors of cov in OD space (orthogonal as cov symmetric)\n    _, V = np.linalg.eigh(np.cov(OD, rowvar=False))\n\n    # The two principle eigenvectors\n    V = V[:, [2, 1]]\n\n    # Make sure vectors are pointing the right way\n    if V[0, 0] < 0: V[:, 0] *= -1\n    if V[0, 1] < 0: V[:, 1] *= -1\n\n    # Project on this basis.\n    That = np.dot(OD, V)\n\n    # Angular coordinates with repect to the prinicple, orthogonal eigenvectors\n    phi = np.arctan2(That[:, 1], That[:, 0])\n\n    # Min and max angles\n    minPhi = np.percentile(phi, 100 - angular_percentile)\n    maxPhi = np.percentile(phi, angular_percentile)\n\n    # the two principle colors\n    v1 = np.dot(V, np.array([np.cos(minPhi), np.sin(minPhi)]))\n    v2 = np.dot(V, np.array([np.cos(maxPhi), np.sin(maxPhi)]))\n\n    # Order of H and E.\n    # H first row.\n    if v1[0] > v2[0]:\n        HE = np.array([v1, v2])\n    else:\n        HE = np.array([v2, v1])\n\n    return normalize_matrix_rows(HE)\n\n######################################################################################################\n\ndef mapping(target,source):\n    \n    stain_matrix_target = get_stain_matrix(target)\n    target_concentrations = get_concentrations(target,stain_matrix_target)\n    maxC_target = np.percentile(target_concentrations, 99, axis=0).reshape((1, 2))\n    stain_matrix_target_RGB = convert_OD_to_RGB(stain_matrix_target) \n    \n    stain_matrix_source = get_stain_matrix(source)\n    source_concentrations = get_concentrations(source, stain_matrix_source)\n    maxC_source = np.percentile(source_concentrations, 99, axis=0).reshape((1, 2))\n    source_concentrations *= (maxC_target / maxC_source)\n    tmp = 255 * np.exp(-1 * np.dot(source_concentrations, stain_matrix_target))\n    return tmp.reshape(source.shape).astype(np.uint8)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def trans(img):\n    ### Stain Transformation\n    img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n    I_LAB = cv2.cvtColor(img, cv2.COLOR_RGB2LAB)\n    L = I_LAB[:, :, 0] / 255.0  # Convert to range [0,1].\n    mask = L < 0.8\n    if mask.sum() == 0:\n        trans = img\n        #empty_img.append(name)\n        #print('empty img')\n    elif mask.sum() == 1:\n        trans = img\n        #almost_empty_img.append(name)\n        print('almost empty img')\n    else:\n        trans = mapping(target,img)    \n    \n    return trans","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Data","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"target = cv2.imread('../input/panda-tiles-16x128x128/train/002a4db09dad406c85505a00fb6f6144_0.png')\ntarget = cv2.cvtColor(target, cv2.COLOR_BGR2RGB)\nplt.imshow(target)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def tile(img):\n    shape = img.shape\n    pad0,pad1 = (sz - shape[0]%sz)%sz, (sz - shape[1]%sz)%sz\n    img = np.pad(img,[[pad0//2,pad0-pad0//2],[pad1//2,pad1-pad1//2],[0,0]],\n                 constant_values=255)\n    img = img.reshape(img.shape[0]//sz,sz,img.shape[1]//sz,sz,3)\n    img = img.transpose(0,2,1,3,4).reshape(-1,sz,sz,3)\n    if len(img) < N:\n        img = np.pad(img,[[0,N-len(img)],[0,0],[0,0],[0,0]],constant_values=255)\n    idxs = np.argsort(img.reshape(img.shape[0],-1).sum(-1))[:N]\n    img = img[idxs]\n   \n    return img\n\n\n\n#mean = torch.tensor([1.0-0.90949707, 1.0-0.8188697, 1.0-0.87795304])\n#std = torch.tensor([0.36357649, 0.49984502, 0.40477625])\n\n### Image Stats for stain normalized tiles 12x64x64\nmean = torch.tensor([1.0-0.90008685, 1.0-0.80557228, 1.0-0.88988811])\nstd = torch.tensor([0.3247047, 0.41417728, 0.33721329])\n\n\n\nclass PandaDataset(Dataset):\n    def __init__(self, path, test):\n        self.path = path\n        self.names = list(pd.read_csv(test).image_id)\n\n    def __len__(self):\n        return len(self.names)\n\n    def __getitem__(self, idx):\n        name = self.names[idx]\n        #img = skimage.io.MultiImage(os.path.join(TEST_DATA,name+'.tiff'))[-1]\n        img = skimage.io.MultiImage(os.path.join(self.path,name+'.tiff'))[-1]\n        \n        tiles = tile(img)\n        \n        ### Stain Normalization ###\n        TMP = np.zeros(((N,sz,sz,3)))\n        #print(TMP.shape)\n        for i in range(N):\n            tran = trans(tiles[i])\n            #print(tran.shape)\n            TMP[i] = tran\n        ############################\n        \n        img_resized = resize(TMP, downsize)\n        img_resized = (255*img_resized).astype(np.uint8)\n        \n        \n        #tiles = torch.Tensor(1.0 - TMP/255.0)\n        tiles = torch.Tensor(1.0 - img_resized/255.0)\n        tiles = (tiles - mean)/std\n        return tiles.permute(0,3,1,2), name","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Show stain normalization examples","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"fig = plt.figure()\nfig, ax = plt.subplots(4,5, figsize=(20,20))\ni = 0\nnames = os.listdir(TRAIN_DATA)[:4]\n#img = skimage.io.MultiImage(os.path.join(TRAIN_DATA,names))[-1]\n#mask = skimage.io.MultiImage(os.path.join(MASKS,name+'_mask.tiff'))[-1]\n#tiles = tile(img,mask)\n\nfor name in tqdm(names):\n    img = skimage.io.MultiImage(os.path.join(TRAIN_DATA,name))[-1]\n    tiles = tile(img)\n    for j in range(2):\n        #j = 2j + 1\n        img = tiles[j]\n        tran = trans(img)\n\n        ax[i][0].imshow(target)\n        ax[i][0].set_title(\"Target Image\",fontsize=10)\n        ax[i][2*j+1].imshow(img)\n        ax[i][2*j+1].set_title(\"Source Image\",fontsize=10)\n        ax[i][2*j+2].imshow(tran)\n        ax[i][2*j+2].set_title(\"Transformed Image\",fontsize=10)\n    i = i + 1\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Prediction","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df = pd.read_csv(SAMPLE)\nif os.path.exists(TEST_DATA):\n    print('Starting predictions')\n    ds = PandaDataset(TEST_DATA,TEST)\n    dl = DataLoader(ds, batch_size=bs, num_workers=nworkers, shuffle=False)\n    names,probs,preds = [],[],[]\n    \n    with torch.no_grad():\n        for x,y in tqdm(dl):\n            x = x.cuda()\n            \n            #dihedral TTA\n            x = torch.stack([x,x.flip(-1),x.flip(-2),x.flip(-1,-2),\n              x.transpose(-1,-2),x.transpose(-1,-2).flip(-1),\n              x.transpose(-1,-2).flip(-2),x.transpose(-1,-2).flip(-1,-2)],1)\n            sz = 64 #after resize\n            x = x.view(-1,N,3,sz,sz)\n            p = [model(x) for model in models]\n            p = torch.stack(p,1)\n            #p = p.view(bs,8*len(models),-1).mean(1).argmax(-1).cpu()\n            prob = p.view(bs,8*len(models),-1).mean(1).cpu()\n            pred = prob.argmax(-1)\n            names.append(y)\n            probs.append(prob)\n            preds.append(pred)\n    \n\n    \n    names = np.concatenate(names)\n    probs = torch.cat(probs).numpy()\n    preds = torch.cat(preds).numpy()\n    sub_df = pd.DataFrame({'image_id': names, 'isup_grade': preds})\n    sub_df.to_csv('submission.csv', index=False)\n    sub_df.head()\n    \n    \nelse:\n    print('Found No test data')\n\n    ds = PandaDataset(TRAIN_DATA,TRAIN)\n    dl = DataLoader(ds, batch_size=bs, num_workers=nworkers, shuffle=False)\n    names,probs,preds = [],[],[]\n\n    #with torch.no_grad():\n    #    for x,y in tqdm(dl):\n    #        x = x.cuda()\n    #        #dihedral TTA\n    #        x = torch.stack([x,x.flip(-1),x.flip(-2),x.flip(-1,-2),\n    #          x.transpose(-1,-2),x.transpose(-1,-2).flip(-1),\n    #          x.transpose(-1,-2).flip(-2),x.transpose(-1,-2).flip(-1,-2)],1)\n    #        sz = 64 #after resize\n    #        x = x.view(-1,N,3,sz,sz)\n    #        p = [model(x) for model in models]\n    #        p = torch.stack(p,1)\n    #        #p = p.view(bs,8*len(models),-1).mean(1).argmax(-1).cpu()\n    #        prob = p.view(bs,8*len(models),-1).mean(1).cpu()\n    #        pred = prob.argmax(-1)\n    #        names.append(y)\n    #        probs.append(prob)\n    #        preds.append(pred)\n        \n    #names = np.concatenate(names)\n    #probs = torch.cat(probs).numpy()\n    #preds = torch.cat(preds).numpy()\n    #sub_df = pd.DataFrame({'image_id': names, 'isup_grade': preds})\n    #sub_df.to_csv('submission.csv', index=False)\n    #sub_df.head()\n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df.to_csv(\"submission.csv\", index=False)\nsub_df.head()","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":4}