{"cells":[{"metadata":{"_uuid":"d93a4b4073344f91c117d105278396c1ccf9718c"},"cell_type":"markdown","source":"I read the stuff that was posted under Acknowledgements namely this paper: Veeling et al \"Rotation Equivariant CNNs for Digital Pathology\"\nI tried to implement it (mainly made it useable with the fastai library since it was already implemented and available).\n <img src=\"https://raw.githubusercontent.com/basveeling/keras-gcnn/master/model.png\" >\nMicroscopy pictures are usually rotation and flip invariant. So a CNN that uses the image and all rotations by 90° as well as flipped left right and upside down, should be better suited than standard CNN.\nI largely used the code by [Adam Bielsky](https://github.com/adambielski) and [Taco Cohen's groupy](https://github.com/tscohen/GrouPy) but changed the batchnorm from 2D to 3D (this might be wrong but otherwise I couldnt get it to work). \n\nI originally intended this for the [Human Protein Atlas Image Classification](https://www.kaggle.com/c/human-protein-atlas-image-classification) but since it uses 8 Layers to include the rotations the memory, and training time is just too much for that data set.\nTheres also a keras implementation of a  [gDensenet ](https://github.com/basveeling/keras-gcnn/blob/master/keras_gcnn/applications/densenetnew.py) by Veeling online if you dont like the fastai stuff\n"},{"metadata":{"_uuid":"7e928da10247b9de66f9ee1114acedff2fe6f1f6"},"cell_type":"markdown","source":"# Imports (hidden)"},{"metadata":{"trusted":true,"_uuid":"a503853adef770ddb1513e85ad2d377f8965c1db","_kg_hide-input":true,"_kg_hide-output":false},"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\n\nimport numpy as np\nimport os\nfrom sklearn.metrics import f1_score\n\nfrom fastai import *\nfrom fastai.vision import *\n\nimport torch\nimport torch.nn as nn\nimport torchvision\nimport cv2\n\nfrom tqdm import tqdm\nfrom skmultilearn.model_selection import iterative_train_test_split\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.preprocessing import MultiLabelBinarizer\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\n%load_ext autoreload\n%autoreload","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"bdce7eee3f29d6489ac6487fbbd4e065b11038fe","_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"!git clone https://github.com/adambielski/GrouPy.git\n\nimport sys\nsys.path.insert(0,'./GrouPy/')\n\nimport groupy\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\n\nfrom torch.autograd import Variable\nfrom groupy.gconv.pytorch_gconv import P4ConvZ2, P4ConvP4\nfrom groupy.gconv.pytorch_gconv.pooling import plane_group_spatial_max_pooling\n\n\nfrom torch.autograd import Variable\nfrom groupy.gconv.pytorch_gconv import P4MConvZ2, P4MConvP4M","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"_uuid":"a33c8d5ad8998d9ee2a38d6547c383e833b90308"},"cell_type":"markdown","source":"# Definition of the group invariant Resnet (hidden)"},{"metadata":{"trusted":true,"_uuid":"ae007e0ae7658d364a10542a17c4c91c37f746b1","_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"class BasicBlock(nn.Module):\n    expansion = 1\n\n    def __init__(self, in_planes, planes, stride=1):\n        super(BasicBlock, self).__init__()\n        self.conv1 = P4MConvP4M(in_planes, planes, kernel_size=3, stride=stride, padding=1, bias=False)\n        self.bn1 = nn.BatchNorm3d(planes)\n        self.conv2 = P4MConvP4M(planes, planes, kernel_size=3, stride=1, padding=1, bias=False)\n        self.bn2 = nn.BatchNorm3d(planes)\n\n        self.shortcut = nn.Sequential()\n        if stride != 1 or in_planes != self.expansion*planes:\n            self.shortcut = nn.Sequential(\n                P4MConvP4M(in_planes, self.expansion*planes, kernel_size=1, stride=stride, bias=False),\n                nn.BatchNorm3d(self.expansion*planes)\n            )\n\n    def forward(self, x):\n        out = F.relu(self.bn1(self.conv1(x)))\n        out = self.bn2(self.conv2(out))\n        out += self.shortcut(x)\n        out = F.relu(out)\n        return out\n\n\nclass Bottleneck(nn.Module):\n    expansion = 4\n\n    def __init__(self, in_planes, planes, stride=1):\n        super(Bottleneck, self).__init__()\n        self.conv1 = P4MConvP4M(in_planes, planes, kernel_size=1, bias=False)\n        self.bn1 = nn.BatchNorm3d(planes)\n        self.conv2 = P4MConvP4M(planes, planes, kernel_size=3, stride=stride, padding=1, bias=False)\n        self.bn2 = nn.BatchNorm3d(planes)\n        self.conv3 = P4MConvP4M(planes, self.expansion*planes, kernel_size=1, bias=False)\n        self.bn3 = nn.BatchNorm3d(self.expansion*planes)\n\n        self.shortcut = nn.Sequential()\n        if stride != 1 or in_planes != self.expansion*planes:\n            self.shortcut = nn.Sequential(\n                P4MConvP4M(in_planes, self.expansion*planes, kernel_size=1, stride=stride, bias=False),\n                nn.BatchNorm3d(self.expansion*planes)\n            )\n\n    def forward(self, x):\n        out = F.relu(self.bn1(self.conv1(x)))\n        out = F.relu(self.bn2(self.conv2(out)))\n        out = self.bn3(self.conv3(out))\n        out += self.shortcut(x)\n        out = F.relu(out)\n        return out\n\n\nclass ResNet(nn.Module):\n    def __init__(self, block, num_blocks, num_classes=10):\n        super(ResNet, self).__init__()\n        self.in_planes = 23\n\n        self.conv1 = P4MConvZ2(3, 23, kernel_size=3, stride=1, padding=1, bias=False)\n        self.bn1 = nn.BatchNorm3d(23)\n        self.layer1 = self._make_layer(block, 23, num_blocks[0], stride=1)\n        self.layer2 = self._make_layer(block, 45, num_blocks[1], stride=2)\n        self.layer3 = self._make_layer(block, 91, num_blocks[2], stride=2)\n        self.layer4 = self._make_layer(block, 181, num_blocks[3], stride=2)\n        self.linear = nn.Linear(181*8*block.expansion, num_classes)\n\n    def _make_layer(self, block, planes, num_blocks, stride):\n        strides = [stride] + [1]*(num_blocks-1)\n        layers = []\n        for stride in strides:\n            layers.append(block(self.in_planes, planes, stride))\n            self.in_planes = planes * block.expansion\n        return nn.Sequential(*layers)\n\n    def forward(self, x):\n        out = F.relu(self.bn1(self.conv1(x)))\n        out = self.layer1(out)\n        \n        out = self.layer2(out)\n        \n        out = self.layer3(out)\n        \n        out = self.layer4(out)\n        \n        outs = out.size()\n        print(out.size(),out.size(0))\n        \n        out = out.view(outs[0], outs[1]*outs[2], outs[3], outs[4])\n        print(out.size(),out.size(0))\n\n        out = F.avg_pool2d(out, 4)\n        print(out.size(),out.size(0))\n        out = out.view(out.size(0), -1)\n        print(out.size())\n        out = self.linear(out)\n        print(out.size())\n        return out\n\n\ndef ResNet18(pretrained=False):\n    return ResNet(BasicBlock, [2,2,2,2])\n\ndef ResNet34(pretrained=False):\n    return ResNet(BasicBlock, [3,4,6,3])\n\ndef ResNet50(pretrained=False):\n    return ResNet(Bottleneck, [3,4,6,3])\n\ndef ResNet101(pretrained=False):\n    return ResNet(Bottleneck, [3,4,23,3])\n\ndef ResNet152(pretrained=False):\n    return ResNet(Bottleneck, [3,8,36,3])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e3a5c8e5cc0d06b10ef26c12e05f70340f1312ff"},"cell_type":"markdown","source":"### Defining necessary variables and load the data"},{"metadata":{"trusted":true,"_uuid":"5caf6cf236cf9d098331763b434806cd73d13823"},"cell_type":"code","source":"model_path='.'\npath='../input/'\ntrain_folder=f'{path}train'\ntest_folder=f'{path}test'\ntrain_lbl=f'{path}train_labels.csv'\nORG_SIZE=96\n\nbs=32\nnum_workers=None # Apprently 2 cpus per kaggle node, so 4 threads I think\nsz=96","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"ecb840384eb1e0aa90fe881f8c212beb135ad6c4"},"cell_type":"code","source":"df_trn=pd.read_csv(train_lbl)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7982793320648eb1ecf956b46da3f2da0d81aaef"},"cell_type":"code","source":"tfms = get_transforms(do_flip=False, flip_vert=False, max_rotate=.0, max_zoom=.3,\n                      max_lighting=0.15, max_warp=0.10)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"b1c089e47c2325cf783092cc235b9b5e9831022a"},"cell_type":"markdown","source":"\nIn Case I want to run quick tests use a subsample:\nI check if the last batch contains one element because otherwise batchnorm crashes."},{"metadata":{"trusted":true,"_uuid":"74e11bfba18cd641ada39ff8a78443d85492e434","scrolled":true},"cell_type":"markdown","source":"p=.1\nn=int(p*df_trn.shape[0])\nwhile True:\n    if p<1.0:\n        idx=np.random.choice(len(df_trn.index),n-np.random.randint(0,100))\n    else:\n        idx=range(len(df_trn.index))\n    X_train, y_train, X_val, y_val = iterative_train_test_split(df_trn.index[idx,None], df_trn.label[idx,None], test_size = 0.2)\n    print(X_train.shape,(X_train.shape[0]%bs),X_val.shape,(X_val.shape[0]%bs))\n    if  ((X_train.shape[0]%bs!=1) and (X_val.shape[0]%bs!=1)):\n        break\nX_train.shape,(X_train.shape[0]%bs),X_val.shape,(X_val.shape[0]%bs)\n\n\nsrc = (ImageItemList.from_csv(path,train_lbl, folder='train', suffix='.tif')\n        .split_by_idxs(X_train[:,0],X_val[:,0])\n        .label_from_lists(y_train.tolist(),y_val.tolist()))\n\ndata = (src.transform(tfms, size=sz)\n        .databunch(num_workers=0,bs=bs))\nstats=data.batch_stats()        \ndata.normalize(stats)"},{"metadata":{"_uuid":"6c091170357f0037a32675f8335f75f851766380"},"cell_type":"markdown","source":"### Otherwise the fulll dataset including the test set"},{"metadata":{"trusted":true,"_uuid":"546c687ba7770602c688c01a58427c47f48c80ac"},"cell_type":"code","source":"data = ImageDataBunch.from_csv(path,csv_labels=train_lbl,folder='train',bs=bs, ds_tfms=tfms, size=sz, suffix='.tif',test=test_folder);\nstats=data.batch_stats()        \ndata.normalize(stats)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"scrolled":true,"_uuid":"e8fec9df84b343e53c9e4379881967f713f17873"},"cell_type":"code","source":"data.show_batch(rows=3, figsize=(12,9))","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"25394f60865c798e0d91837ab9d48fb2ffd04512"},"cell_type":"markdown","source":"Defining the metric that is used in the competition as well as the loss function"},{"metadata":{"trusted":true,"_uuid":"b83ed47a0109f7d51f1c475958a5191f87c63b01"},"cell_type":"code","source":"from sklearn.metrics import roc_auc_score\ndef auc_score(y_pred,y_true):  \n    return roc_auc_score(to_np(y_true),to_np(F.sigmoid(y_pred))[:,1])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9098feefc03dad43016c7871df7dc086d3b13757"},"cell_type":"code","source":"#f1=Fbeta_binary(beta2=1)\nf1=partial(fbeta,beta=1)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"5414edd4805da95791fabed8ec7042213d71e170"},"cell_type":"markdown","source":"## Define a custom head for the gCNN\nThe original resnet has al ot more categories than 2 so a new head is needed.\nThe fastai package already delivers an awesome function for this. I had to define 2 little functions to adapt it to the group equivariant, namely make the pooling 3D not 2D"},{"metadata":{"trusted":true,"_uuid":"f77af9103304cd864a462159d3454a033b46a698","_kg_hide-output":true},"cell_type":"code","source":"class AdaptiveConcatPool3d(nn.Module):\n    \"Layer that concats `AdaptiveAvgPool3d` and `AdaptiveMaxPool3d`.\"\n    def __init__(self, sz:Optional[int]=None):\n        \"Output will be 2*sz or 2 if sz is None\"\n        super().__init__()\n        sz = sz or 1\n        self.ap,self.mp = nn.AdaptiveAvgPool3d(sz), nn.AdaptiveMaxPool3d(sz)\n    def forward(self, x): return torch.cat([self.mp(x), self.ap(x)], 1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"99ae81d4b5eb761e8a2d7e1856ae420c5f89acae","_kg_hide-output":true},"cell_type":"code","source":"def create_head_custom(nf:int, nc:int, lin_ftrs:Optional[Collection[int]]=None, ps:Floats=0.5, bn_final:bool=False):\n    \"Model head that takes `nf` features, runs through `lin_ftrs`, and about `nc` classes.\"\n    lin_ftrs = [nf, 512, nc] if lin_ftrs is None else [nf] + lin_ftrs + [nc]\n    ps = listify(ps)\n    if len(ps)==1: ps = [ps[0]/2] * (len(lin_ftrs)-2) + ps\n    actns = [nn.ReLU(inplace=True)] * (len(lin_ftrs)-2) + [None]\n    layers = [AdaptiveConcatPool3d(), Flatten()]\n    for ni,no,p,actn in zip(lin_ftrs[:-1],lin_ftrs[1:],ps,actns):\n        layers += bn_drop_lin(ni,no,True,p,actn)\n    if bn_final: layers.append(nn.BatchNorm1d(lin_ftrs[-1], momentum=0.01))\n    return nn.Sequential(*layers)\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a3da74e49fe40d07946dc00c612b711f654987fe"},"cell_type":"markdown","source":"## Now the data and model are merged in the fastai learner object\nSince pretrained is set to False the network is unfrozen so we dont need a warm up period where we only train the head of the network"},{"metadata":{"trusted":true,"_uuid":"e87abf1fbf9b061cb26fb05afc6932e253d4135c"},"cell_type":"code","source":"learn = create_cnn(\n    data,\n    ResNet18,\n    path='.',    \n    metrics=[fbeta], \n    ps=0.5,\n    pretrained=False\n)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"ba719e2b6518b7e9eab4757395aba7b8c6634af9"},"cell_type":"markdown","source":"## Now add the custom head"},{"metadata":{"trusted":true,"_uuid":"04ce6c3c0361d71c40fd232521c1e60fe4eb589e","_kg_hide-output":true},"cell_type":"code","source":"learn.model[1]=create_head_custom(362,2)\nlearn.model.to(learn.data.device)\nprint(learn.summary())\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"5082e4f5d5c9d00b663a0054962fd887dccb004c"},"cell_type":"markdown","source":"## And check for an appropriate learning rate\nThe maximum learning rate should be somewhere before the divergence/minimum and the minimum of the learning rate about a 10th of the maximum"},{"metadata":{"trusted":true,"_uuid":"d8f689bfd03e79e5f8f130032468a85788b7acce"},"cell_type":"code","source":"learn.lr_find(num_it=500,end_lr=1000)\nlearn.recorder.plot()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"fcfa743233605eafacfcbc383da90c7bc5fbd473","scrolled":true},"cell_type":"code","source":"learn.fit_one_cycle(2,slice(5e-5,1e-4))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1aeb8db4461f7110665b52869527604c04295fdf"},"cell_type":"code","source":"learn.save('gCNN_1cycle')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"61562941a9c242376136f4774e0117f1cf0b470b"},"cell_type":"markdown","source":"Theres an error I don't understand yet in plot_loss() when only doing one epoch (I think)"},{"metadata":{"trusted":true,"_uuid":"fcfa743233605eafacfcbc383da90c7bc5fbd473"},"cell_type":"code","source":"learn.recorder.plot()\ntry:\n    learn.recorder.plot_losses()\n    learn.recorder.plot_lr()\nexcept:\n    pass","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"5c4e009e09fb77b08997ac8e9ee69660c28676be"},"cell_type":"markdown","source":"## Check how well gCNN are doing\nNow that we learned a lot we check how well we are doing by predicting on the validation set.\nFor just one epoch this is pretty good!\nI achieved a LB score >0.96 with this, without TTA."},{"metadata":{"trusted":true,"_uuid":"b6917d2f83141d472e04073b59a9b5e1e6b26c57"},"cell_type":"code","source":"preds,y=learn.get_preds()\npred_score=auc_score(preds,y)\npred_score","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"bdfd3b80d33057d155e7da32391817a7e88ff80f"},"cell_type":"markdown","source":"### Now predict on test set"},{"metadata":{"trusted":true,"_uuid":"6c1890c9509ca13a01b6bba3bba296b5a9b34fa7"},"cell_type":"code","source":"preds_test,y_test=learn.get_preds(ds_type=DatasetType.Test)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6a85d720a61b996e3ee308c4a012345ef69fc461","scrolled":true},"cell_type":"code","source":"preds_test_tta,y_test_tta=learn.TTA(ds_type=DatasetType.Test)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"88a566de210643baa4a44519e67b6e377c9d26d7"},"cell_type":"markdown","source":"### prepare submission\nI now load in the sample submission and put my predictions in the label column and save to a new file."},{"metadata":{"_uuid":"c7f2e327a7943ed6e8d7f8bfc7bd3f9682186b7c"},"cell_type":"markdown","source":"Sometimes its important in which order the ids in the submissions are so to make sure I don't mess up I put them in the same order. My first submission had a 50% score so I somewhere messed up the order oder the matching of id to label.\nsince fname_clean is the id we can just use that as index when adding the correct label in our dataframe. "},{"metadata":{"trusted":true,"_uuid":"7c64426066b5d158024e4bbf05a6a05c10777ce7"},"cell_type":"code","source":"sub=pd.read_csv(f'{path}/sample_submission.csv').set_index('id')\nsub.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"20c89b9c62d4d98fe4cae6195e80338abfb3a553"},"cell_type":"code","source":"clean_fname=np.vectorize(lambda fname: str(fname).split('/')[-1].split('.')[0])\nfname_cleaned=clean_fname(data.test_ds.items)\nfname_cleaned=fname_cleaned.astype(str)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"25f186dca18d7855cf8b3910a45acd3de9f1f34e"},"cell_type":"markdown","source":"## I add the score to the name of the file so I can later plot the leaderboard score versus my validation score\nIn the fastai course Jeremy mentions that if you have a monotonic relation between validation and LB score the way you set up your validation set matches what the test set consists of."},{"metadata":{"trusted":true,"_uuid":"1e9f92573996df4ff7a5f047e7ffae6ef6e53b73"},"cell_type":"code","source":"sub.loc[fname_cleaned,'label']=to_np(preds_test[:,1])\nsub.to_csv(f'submission_{pred_score}.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6a934db6adbf52935aab26a9c6400bf0530d3bb9"},"cell_type":"code","source":"sub.loc[fname_cleaned,'label']=to_np(preds_test_tta[:,1])\nsub.to_csv(f'submission_TTA_{pred_score}.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"15af9ff7c16ef0e5bda33d4595b482c91d21f305"},"cell_type":"code","source":"!rm -r GrouPy/","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"6d57af97ba7c0404c55d7b05926d317bf364f5ce"},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}