{"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":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport os\nfilename = None\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    \n    for filename in filenames:\n        if filename == 'train.pt':\n            dataset_dir = dirname\n        \n        print(os.path.join(dirname, filename))\nprint(\"\\nInput directory: \", dirname)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2021-10-24T02:14:37.290907Z","iopub.execute_input":"2021-10-24T02:14:37.291251Z","iopub.status.idle":"2021-10-24T02:14:37.306398Z","shell.execute_reply.started":"2021-10-24T02:14:37.291214Z","shell.execute_reply":"2021-10-24T02:14:37.305257Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Our approach uses Pytorch Geometric. It is not currently included in Kaggle's environment, so we must install it via  `pip`. \n\nUnfortunately, building these wheel files for CUDA=11.0 is *reaaalllllly slow* (>1 hour) for some reason, and that is the version available for the Kaggle environment.","metadata":{}},{"cell_type":"code","source":"!pip install torch-scatter torch-sparse torch-cluster torch-spline-conv torch-geometric -f https://data.pyg.org/whl/torch-1.9.0+cu110.html\nimport torch\nfrom torch_geometric.data import DataLoader,Dataset,Data","metadata":{"_kg_hide-output":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#### Now that this is out of the way, let us *read the Pytorch Geometric version of our dataset* (created on a previous kernel: https://www.kaggle.com/isonettv/converting-kddbr-2021-data-to-pytorch-files/notebook). \n\nThere are 2 files, **'train.pt'** and **'test.pt'**. Both of them are lists, with each element of them being an instance of the TSP problem, represented as a ``pytorch_geometric.data.Data`` object. The attributes of each element are:\n1. **x**: tensor with shape [200,2]. The coordinates for up to 200 vertices. Padded out with zeros in case the number is < 200.\n2. **m**: tensor with shape [200]. A mask such that ``x[m]`` returns the coordinates of the actual vertices.\n3. **y_index**: tensor with shape [?,2]. Each row corresponds to some edge $(i,j)$ that is part of the solution. Only present in `train.pt`.","metadata":{}},{"cell_type":"code","source":"import torch\nfrom torch_geometric.data import DataLoader,Dataset,Data\nimport pandas as pd\nimport numpy as np\nimport gc\nimport torch.nn.functional as F\nimport os.path as osp\n\ndevice = 'cuda'\n\nL = torch.load(osp.join(dataset_dir,'train.pt'),\n               map_location=device)\nL_train = [L[i] for i in range(len(L)) if i % 10 > 0]\nL_val = [L[i] for i in range(len(L)) if i % 10 == 0]\nimport random\nrandom.Random(4).shuffle(L_val)\n\nL_test = torch.load(osp.join(dataset_dir,'test.pt'),\n                    map_location=device)\n\nprint(\"Training example: \", L_train[0])\nprint(\"Test example: \", L_train[0])","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Next, we implement the ``Dataset`` interface from Pytorch Geometric. When we index the dataset (`dataset['train'][i]`, for example), we look up the corresponding TSP problem and apply some small transformations only at runtime, saving up on memory:\n\n1. We dynamically compute the KNN graph (k=20 is a hyperparameter)\n2. We transform the target edge list (**y_index**) into a dense adjacency matrix, which we use for the cost function.","metadata":{}},{"cell_type":"code","source":"import torch\nfrom torch_geometric.data import DataLoader,Dataset,Data\nimport pandas as pd\nimport numpy as np\nimport gc\nimport torch.nn.functional as F\nfrom torch_geometric.utils import to_dense_adj\nfrom torch_geometric.transforms import KNNGraph,Distance\n\nfrom torch_geometric.nn import knn_graph\n\n\n\nclass KDD_Dataset(Dataset):\n    \"\"\"\n    Base class representing a dataset for the competition\n    \"\"\"\n    \n    distance_transf = Distance()\n    knn_transf = KNNGraph(k=20)\n    def __init__(self, L,mode='train'):\n        \"\"\"\n        Args:\n            csv_file (string): Path to the csv file with annotations.\n            root_dir (string): Directory with all the .pt files.\n            transform (callable, optional): Optional transform to be applied\n                on a sample.\n        \"\"\"\n        self.L = L\n        \n        \n    def __len__(self):\n        return len(self.L)\n    def __getitem__(self, idx):\n        sample = self.L[idx]\n        num_nodes = torch.where(sample.m==True)[0][-1] + 1\n        \n        num_nodes_vector = num_nodes * torch.ones(200,device=sample.x.device)\n        num_nodes_vector[sample.m.view(-1)==False] = 1\n        \n        m = sample.m.view(-1,1).bool()\n            \n        \n        data = Data(m=m, pos=sample.x[0:num_nodes,0:num_nodes],\n                         num_nodes=num_nodes,\n                         num_nodes_vector=num_nodes_vector)\n        if hasattr(sample,'y_index'):\n            data.y = to_dense_adj(sample.y_index.transpose(0,1).long(),max_num_nodes=200)\n        \n        data = self.knn_transf(data)\n        data = self.distance_transf(data)\n        \n        data.num_nodes=200\n        data.pos = sample.x\n        data.x = sample.x\n        \n        data.n_edges = data.edge_index.shape[1]\n\n        return data   \n\ndatasets = {\"train\":KDD_Dataset(L_train,mode='train'),\n            \"val\":KDD_Dataset(L_val,mode='val'),\n            \"test\":KDD_Dataset(L_test,mode='test'),\n           }\n\nbatch_size = 32\ndataloaders = {\"train\":DataLoader(datasets['train'],batch_size=batch_size,shuffle=True),\n            \"val\":DataLoader(datasets['val'],batch_size=batch_size,shuffle=True),\n            \"test\":DataLoader(datasets['test'],batch_size=batch_size),\n           }\nfor k in dataloaders.keys():\n    print(f\"Length of {k} dataset : {len(datasets[k])}\")","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### The next step is to define our model. \n\n#### We feed our input data to the **edge convolutional operator** from the **“Dynamic Graph CNN for Learning on Point Clouds”** paper :https://arxiv.org/abs/1801.07829. After a few layers of EdgeConv, we obtain our final prediction using a **Multilayer Perceptron**. LeakyRELUs serve as the activation function.\n\n#### We have not performed any sort of extensive hyperparameter tuning, though this configuration seemed to perform best from our (limited) experiments","metadata":{}},{"cell_type":"code","source":"import torch\nimport pandas as pd\nfrom torch_geometric.nn import MessagePassing\nfrom torch_geometric.utils import add_self_loops, degree\n\nfrom torch.nn import Sequential, Linear, BatchNorm1d, ReLU, LeakyReLU, Sigmoid\nfrom torch_geometric.nn.conv import GATConv,TransformerConv,GCNConv,GINConv,EdgeConv\nimport torch.nn.functional as F\nfrom torch_scatter import scatter\nclass ChallengeModel(torch.nn.Module):\n\n    \n    def __init__(self,num_iterations=5,H=64,**kwargs):\n        super(ChallengeModel, self).__init__(**kwargs)\n        self.num_iterations = num_iterations\n        \n        \n        self.H = H\n        self.convs = []\n        \n        for i in range(self.num_iterations):\n            _in = 2 if i == 0 else H\n            _out = H\n            mlp = Sequential(Linear(2 * _in, 128),\n                       LeakyReLU(),\n                       Linear(128, _out))\n            self.convs.append(EdgeConv(mlp).cuda())\n        \n        setattr(self,'convs',torch.nn.ModuleList(getattr(self,'convs')))\n        \n        \n        self.edgeMLP = Sequential(Linear(2 * H, 512),\n                       LeakyReLU(),\n                       Linear(512, 256),\n                       LeakyReLU(),\n                       Linear(256, 1))\n    \n    def forward(self,data):\n        edge_index = data.edge_index.clone().long()\n        x = sample.pos.float().clone()\n        \n        for i in range(len(self.convs)):\n            x = F.leaky_relu(self.convs[i](x,edge_index))\n        \n        \"\"\"\n            To create the input for our MLP edge classifier, we simply concatenate the outputs for\n            the source and target vertices.\n        \"\"\"\n        X_source = x[edge_index[0,:],:]\n        X_target = x[edge_index[1,:],:]\n        X = torch.cat([X_source,X_target],axis=1).clone()\n        def act(x):\n            return F.leaky_relu(x)\n        \n        X = self.edgeMLP(X)\n\n        L = []\n        num_batches = sample.batch.max()+1\n\n        batch = data.batch.clone().to('cuda')\n        \n        \"\"\"\n            We use this helper function to convert the prediction ´´X´´ to an adjacency matrix format.\n            This yields B 200x200 matrices in total. \n            \n            Note that any edge that is not present to the KNN graph will have its entry corresponding \n        to 0 in the dense format.\n        \n        \"\"\"\n        mat = to_dense_adj(edge_index,edge_attr=X\n                               ,max_num_nodes=200,\n                         batch=batch).squeeze(-1)\n        del batch\n\n        return mat\n               ","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### With our model defined, we will need to add code for the training loop. \n\n\n### Our custom loss function uses the same kind of weight balancing as previous approaches :https://arxiv.org/abs/1906.01227https://arxiv.org/abs/1906.01227. \n$$w_0 = \\frac{n^2}{c(n^2-2n)} $$\n$$w_1 = \\frac{n^2}{c\\times 2n} $$\nIn our case, $c=2$ because we perform binary classification. \n\n---\n\nThe final classification probabilities have to be converted to binary labels, maximizing the F1 score. One possibility is to select the best threshold hyperparameter using the validation data. However, currently we simply set a fixed hyperparameter ``THRESH``.  ","metadata":{}},{"cell_type":"code","source":"from tqdm import tqdm,trange\nimport sklearn.metrics\nfrom sklearn.metrics import f1_score\n\ndef _xent(preds,actuals):\n    \"\"\" Here we assume both preds and actuals are square matrices with the relevant nodes only.\n        preds sums up to 1, actuals sums up to 2.\n    \"\"\"\n    num_nodes = preds.shape[0]\n    preds = preds.reshape(-1)\n    actuals = actuals.reshape(-1)\n    \n    w0 =  (num_nodes**2) / (2*(num_nodes**2 - 2*num_nodes))\n    w1 = (num_nodes**2) / (2*(2*num_nodes))\n    \n    weights = w0 * torch.ones(num_nodes**2,device=preds.device)\n    weights[actuals > 0] = w1\n\n    \n    return F.binary_cross_entropy_with_logits(preds,actuals,\n                                              weight=weights,\n                                              reduction='mean')\n\ndef custom_XENT(preds,actuals,edge_index,m,num_nodes_vector):\n    \"\"\" edge_index[1,:] now points from 0 to \"\"\"\n    nb = actuals.shape[0]\n    m = m.clone().reshape(nb,200)\n    loss = 0\n    for i in range(nb):\n        num_nodes = torch.where(m[i,:]==True)[0][-1].item() + 1\n        loss += _xent(preds[i,0:num_nodes,0:num_nodes],\n                      actuals[i,0:num_nodes,0:num_nodes])\n\n    \n    loss /= nb\n    return loss\n\nTHRESH = 0.56\ndef _metrics(preds,actuals):\n    preds = np.array(preds > THRESH,dtype=np.int32)\n    #print(preds,actuals)\n    return sklearn.metrics.f1_score(preds,actuals)\n\ndef custom_metrics(preds,actuals,edge_index,m,num_nodes_vector):\n    \"\"\" edge_index[1,:] now points from 0 to \"\"\"\n    nb = actuals.shape[0]\n    m = m.clone().reshape(nb,200)\n    loss = 0\n    preds = torch.sigmoid(preds).clone().detach().cpu().numpy()\n    actuals = actuals.clone().detach().cpu().numpy()\n    \n    for i in range(nb):\n        num_nodes = torch.where(m[i,:]==True)[0][-1].item() + 1\n\n        \n        loss += _metrics(preds[i,0:num_nodes,0:num_nodes].reshape((-1,)),\n                      actuals[i,0:num_nodes,0:num_nodes].reshape((-1,)))\n\n    \n    loss /= nb\n    return loss\n\nimport numpy as np\nbest_f1 = -1\nmodel = ChallengeModel().cuda()\n\nfrom torch.utils.tensorboard import SummaryWriter\nwriter = SummaryWriter()\nnum_epochs = 10\nopt = torch.optim.Adam(lr=1e-3,params=model.parameters())\nstep = 0\n\n#model.load_state_dict('./model/model_3650.pt')\ntorch.manual_seed(420)\n        \n        \nfor epoch in range(num_epochs):\n    \n        \n    for mode in ['train','val']:\n        if mode == 'train':\n                model.train()\n        else:\n                model.eval()\n        stats = {'loss':[],\n                 'f1_score':[],\n                }\n                \n        running_loss = []\n        out_mean = []\n        occup_mean = []\n        actual_out_mean = []\n        cnt = 0\n        total=len(dataloaders[mode])//10 if mode == 'train' else len(dataloaders[mode])//10\n        #mode_maybe_shuffle = 'train' if mode == 'train' else mode\n        for i,sample in tqdm(enumerate(dataloaders[mode]),\n                             total=total):\n                \n                \n                with torch.set_grad_enabled(mode == 'train'):\n                    if i == total:\n                        break\n                    if mode == 'train':\n                        opt.zero_grad()\n                    cnt += 1\n                    \n                    out = model(sample)\n                    \n                    \n                    #sklearn.metrics.f1_score()\n                    \n                    if mode == 'train':\n                        loss = custom_XENT(out,sample.y,sample.edge_index.clone().long(),\n                                          sample.m,sample.num_nodes_vector)\n                        loss_np = loss.cpu().item()\n                        loss.backward()\n                        opt.step()\n                    \n                    f1_score = custom_metrics(out,sample.y,sample.edge_index.clone().long(),\n                                          sample.m,sample.num_nodes_vector)\n                    _stats = {'loss':loss,\n                              'f1_score':torch.as_tensor(f1_score*100.0,device=loss.device)\n                             }\n                    for k in _stats.keys():\n                        stats[k].append(_stats[k].cpu().item())\n                    del _stats\n               \n\n        print(\"======================================\")\n        print(f\"Epoch {epoch} - Avg stats ({mode})\")\n        for k in stats.keys():\n            print(f'Mean {k}: {np.array(stats[k]).mean()}')\n            epoch_f1 = np.array(stats['f1_score']).mean()\n            if best_f1 < epoch_f1:\n                best_f1 = epoch_f1\n                torch.save(model.state_dict(),f'./best_model.pt')\n                print(\"New best\")\n    print(\"Flushed\")\n    step += 1\n    torch.cuda.empty_cache()\n\n    ","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Lastly, we put our model into evaluation mode, and  save the predictions to a CSV file.","metadata":{}},{"cell_type":"code","source":"model.eval()\nmodel.load_state_dict(model.state_dict(),f'best_model.pt')\n\nmode='test'\n\n\ndef get_pred(preds,actuals,edge_index,m,num_nodes_vector,upload_file,cnt):\n    \"\"\" edge_index[1,:] now points from 0 to \"\"\"\n    nb = preds.shape[0]\n    m = m.clone().reshape(nb,200)\n    loss = 0\n    preds = torch.sigmoid(preds).clone().detach().cpu().numpy()\n    for i in range(nb):\n        num_nodes = torch.where(m[i,:]==True)[0][-1].item() + 1\n        for i1 in range(num_nodes):\n            for j1 in range(num_nodes):\n                if i1 == j1:\n                    continue\n                else:\n                    upload_file.write(f'{cnt+i}:{i1}-{j1},{1 if preds[i,i1,j1]>THRESH else 0}\\n')\nimport os \ncnt = 0\nwith open('submission.csv', \"w\") as upload_file:\n    \n    upload_file.write('Id,Predicted\\n')\n    \n    for i,sample in tqdm(enumerate(dataloaders[mode]),\n                                 total=len(dataloaders[mode])):\n        out = model(sample)\n        with torch.set_grad_enabled(False):\n            pred = get_pred(out,sample.y,sample.edge_index.clone().long(),\n                                              sample.m,sample.num_nodes_vector,\n                              upload_file,cnt)\n            \n            \n            cnt += sample.batch.max().item() + 1","metadata":{},"execution_count":null,"outputs":[]}]}