analysis.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline \n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import editdistance as ed\n",
    "from tqdm import tqdm_notebook\n",
    "from matplotlib import pyplot as plt\n",
    "from collections import Counter\n",
    "import time\n",
    "\n",
    "import torch\n",
    "from torch import nn\n",
    "from torch import optim\n",
    "from torch.utils.data import TensorDataset, DataLoader\n",
    "\n",
    "from sklearn.metrics import accuracy_score, f1_score, classification_report\n",
    "\n",
    "import uuid\n",
    "\n",
    "import utils\n",
    "import importlib\n",
    "importlib.reload(utils)\n",
    "from utils import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " - Dropped 3106 duplicates\n",
      " - Dropped 62656 duplicates\n",
      " - Removed 15566 sequences\n",
      "80214\n",
      "495796\n"
     ]
    }
   ],
   "source": [
    "df_pos = load_dataset(\"data/GSM3155092_P01_CRVstim_CD8_beta.txt.gz\")\n",
    "df_neg = load_dataset(\"data/GSM3155090_P01_unstim_CD8_beta.txt.gz\")\n",
    "\n",
    "df_neg = filter_specific(df_pos, df_neg)\n",
    "\n",
    "print(len(df_pos))\n",
    "print(len(df_neg))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 131,
   "metadata": {},
   "outputs": [],
   "source": [
    "# model_1.pth\n",
    "# class CMVModel(nn.Module):\n",
    "    \n",
    "#     def __init__(self):\n",
    "#         super(CMVModel, self).__init__()\n",
    "        \n",
    "#         self._hidden_size = 64\n",
    "#         self.rnn = nn.GRU(input_size=len(oh_dict[\"A\"]), hidden_size=self._hidden_size, num_layers=2, bidirectional=True, dropout=0, batch_first=True)\n",
    "#         self.final = nn.Sequential(nn.LeakyReLU(),\n",
    "#                                    nn.BatchNorm1d(self._hidden_size),\n",
    "#                                    nn.Linear(self._hidden_size, self._hidden_size), \n",
    "#                                    nn.LeakyReLU(),\n",
    "#                                    nn.BatchNorm1d(self._hidden_size),\n",
    "#                                    nn.Linear(self._hidden_size, 1), \n",
    "#                                    nn.Sigmoid())\n",
    "#         torch.nn.init.kaiming_uniform_(self.final[2].weight, nonlinearity=\"leaky_relu\")\n",
    "#         torch.nn.init.kaiming_uniform_(self.final[5].weight, nonlinearity=\"leaky_relu\")\n",
    "        \n",
    "    \n",
    "#     def forward(self, batch):\n",
    "#         x = self.rnn(batch)[1]\n",
    "#         x = (x[-2, :, :] + x[-1, :, :]).mul(.5)\n",
    "#         x = self.final(x)\n",
    "#         return x\n",
    "\n",
    "\n",
    "# model_2.pth\n",
    "# class CMVModel(nn.Module):\n",
    "    \n",
    "#     def __init__(self):\n",
    "#         super(CMVModel, self).__init__()\n",
    "        \n",
    "#         self._hidden_size = 64\n",
    "#         self.rnn = nn.GRU(input_size=len(oh_dict[\"A\"]), hidden_size=self._hidden_size, num_layers=2, bidirectional=True, dropout=0.0, batch_first=True)\n",
    "# #         self.final = nn.Sequential(nn.LeakyReLU(),\n",
    "# #                                    nn.BatchNorm1d(self._hidden_size),\n",
    "# #                                    nn.Linear(self._hidden_size, 1), \n",
    "# #                                    nn.Sigmoid())\n",
    "#         self.final = nn.Sequential(nn.Linear(self._hidden_size, 1), \n",
    "#                                    nn.Sigmoid())\n",
    "        \n",
    "    \n",
    "#     def forward(self, batch):\n",
    "#         x = self.rnn(batch)[1]\n",
    "#         x = (x[-2, :, :] + x[-1, :, :]).mul(.5)\n",
    "#         x = self.final(x)\n",
    "#         return x\n",
    "\n",
    "\n",
    "class CMVModel(nn.Module):\n",
    "    \n",
    "    def __init__(self, hidden_size=32, num_layers=1, rnn_dropout=0):\n",
    "        super(CMVModel, self).__init__()\n",
    "        \n",
    "        self._hidden_size = hidden_size\n",
    "        self.rnn = nn.GRU(input_size=len(oh_dict[\"A\"]), hidden_size=self._hidden_size, num_layers=num_layers, bidirectional=True, dropout=rnn_dropout, batch_first=True)\n",
    "        self.final = nn.Sequential(nn.Linear(self._hidden_size, 1), \n",
    "                                   nn.Sigmoid())\n",
    "        \n",
    "    \n",
    "    def forward(self, batch):\n",
    "        x = self.rnn(batch)[1]\n",
    "        x = (x[-2, :, :] + x[-1, :, :]).mul(.5)\n",
    "        x = self.final(x)\n",
    "        return x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 130,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train DFs sizes: 70214 485796\n",
      "Test DFs sizes: 10000 10000\n",
      "Train Tensors sizes: 70214 485796\n",
      "Test Tensors sizes: 10000 10000\n"
     ]
    }
   ],
   "source": [
    "seed = 42\n",
    "oh_dict = load_dict()\n",
    "\n",
    "TEST_SIZE = 10000\n",
    "\n",
    "indices = np.random.choice(len(df_pos), TEST_SIZE, replace=False)\n",
    "df_pos_test = df_pos.iloc[indices, :]\n",
    "df_pos_train = df_pos.drop(indices)\n",
    "\n",
    "indices = np.random.choice(len(df_neg), TEST_SIZE, replace=False)\n",
    "df_neg_test = df_neg.iloc[indices, :]\n",
    "df_neg_train = df_neg.drop(indices)\n",
    "\n",
    "print(\"Train DFs sizes:\", len(df_pos_train), len(df_neg_train))\n",
    "print(\"Test DFs sizes:\", len(df_pos_test), len(df_neg_test))\n",
    "\n",
    "\n",
    "X_pos_train = seq2vec(df_pos_train[\"CDR3.sequence\"], oh_dict)\n",
    "X_pos_test = seq2vec(df_pos_test[\"CDR3.sequence\"], oh_dict)\n",
    "X_neg_train = seq2vec(df_neg_train[\"CDR3.sequence\"], oh_dict)\n",
    "X_neg_test = seq2vec(df_neg_test[\"CDR3.sequence\"], oh_dict)\n",
    "\n",
    "print(\"Train Tensors sizes:\", len(X_pos_train), len(X_neg_train))\n",
    "print(\"Test Tensors sizes:\", len(X_pos_test), len(X_neg_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 132,
   "metadata": {},
   "outputs": [],
   "source": [
    "#\n",
    "# Update the model\n",
    "#\n",
    "model = CMVModel().to(\"cuda\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 167,
   "metadata": {},
   "outputs": [],
   "source": [
    "BATCH_SIZE = 64\n",
    "# LR = 0.0003\n",
    "LR = 0.001\n",
    "\n",
    "indices = np.random.choice(len(X_neg_train), len(X_pos_train), replace=False)\n",
    "X_neg_train_ss = X_neg_train[indices].contiguous()\n",
    "\n",
    "dl_train = DataLoader(TensorDataset(X_neg_train_ss.transpose(2, 1), X_pos_train.transpose(2, 1)), batch_size=BATCH_SIZE, shuffle=True, num_workers=1, pin_memory=True, drop_last=False)\n",
    "dl_test = DataLoader(TensorDataset(X_neg_test.transpose(2, 1), X_pos_test.transpose(2, 1)), batch_size=2048, shuffle=False, num_workers=1, pin_memory=True, drop_last=False)\n",
    "\n",
    "optimiser = optim.Adam(model.parameters(), lr=LR)\n",
    "criterion = nn.BCELoss()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 137,
   "metadata": {},
   "outputs": [],
   "source": [
    "def run_epoch(model, criterion, optimiser, dl, train_mode=True):\n",
    "    if train_mode:\n",
    "        model.train()\n",
    "    else:\n",
    "        model.eval()\n",
    "        \n",
    "    loss_list = []\n",
    "    pred_list = []\n",
    "    true_list = []\n",
    "    \n",
    "    start = time.time()\n",
    "    for batch_neg, batch_pos in dl:\n",
    "        optimiser.zero_grad()\n",
    "        loss = 0\n",
    "        \n",
    "        for batch, y_true in create_batches(batch_neg, batch_pos):\n",
    "            y_pred = model(batch.to(\"cuda\")).reshape((-1,))\n",
    "            loss += criterion(y_pred, y_true.to(\"cuda\"))\n",
    "            \n",
    "            y_pred_round = y_pred.cpu().detach().round().numpy()\n",
    "            pred_list.append(y_pred_round)\n",
    "            true_list.append(y_true.cpu().detach().numpy())\n",
    "        \n",
    "        if train_mode:\n",
    "            loss.backward()\n",
    "            optimiser.step()\n",
    "        \n",
    "        loss_list.append(loss.cpu().item())\n",
    "    end = time.time()\n",
    "    \n",
    "    return loss_list, pred_list, true_list, end - start"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 168,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Epoch:   1\n",
      "  - Train (in 8.8s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.56      0.58      0.57     70214\n",
      "        CMV+       0.57      0.55      0.56     70214\n",
      "  - Test (in 0.73s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.54      0.56      0.55     10000\n",
      "        CMV+       0.54      0.52      0.53     10000\n",
      "\n",
      "Epoch:   2\n",
      "  - Train (in 8.6s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.56      0.58      0.57     70214\n",
      "        CMV+       0.57      0.55      0.56     70214\n",
      "  - Test (in 0.64s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.53      0.55      0.54     10000\n",
      "        CMV+       0.54      0.52      0.53     10000\n",
      "\n",
      "Epoch:   3\n",
      "  - Train (in 8.5s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.57      0.58      0.57     70214\n",
      "        CMV+       0.57      0.56      0.56     70214\n",
      "  - Test (in 0.52s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.53      0.57      0.55     10000\n",
      "        CMV+       0.54      0.50      0.52     10000\n",
      "\n",
      "Epoch:   4\n",
      "  - Train (in 8.5s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.57      0.58      0.58     70214\n",
      "        CMV+       0.57      0.56      0.56     70214\n",
      "  - Test (in 1.0s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.54      0.50      0.52     10000\n",
      "        CMV+       0.53      0.57      0.55     10000\n",
      "\n",
      "Epoch:   5\n",
      "  - Train (in 8.5s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.57      0.58      0.58     70214\n",
      "        CMV+       0.57      0.56      0.57     70214\n",
      "  - Test (in 0.5s): \n",
      "              precision    recall  f1-score   support\n",
      "\n",
      "        CMV-       0.54      0.49      0.51     10000\n",
      "        CMV+       0.53      0.58      0.56     10000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "EPOCHS = 5\n",
    "\n",
    "for i_epoch in range(1, EPOCHS+1):\n",
    "    loss_list, pred_list, true_list, exec_time = run_epoch(model, criterion, optimiser, dl_train, True)\n",
    "    pred_list = np.concatenate(pred_list).ravel()\n",
    "    true_list = np.concatenate(true_list).ravel()\n",
    "    cr = classification_report(true_list, pred_list, target_names=[\"CMV-\", \"CMV+\"])\n",
    "    \n",
    "    print(\"Epoch:{0:4}\".format(i_epoch))\n",
    "    print(\"  - Train (in {0:3.2}s): \".format(exec_time), )\n",
    "    print(cr[:162])\n",
    "    \n",
    "    loss_list, pred_list, true_list, exec_time = run_epoch(model, criterion, optimiser, dl_test, False)\n",
    "    pred_list = np.concatenate(pred_list).ravel()\n",
    "    true_list = np.concatenate(true_list).ravel()\n",
    "    cr = classification_report(true_list, pred_list, target_names=[\"CMV-\", \"CMV+\"])\n",
    "    print(\"  - Test (in {0:3.2}s): \".format(exec_time), )\n",
    "    print(cr[:162])\n",
    "    print()\n",
    "    torch.save(model.state_dict(), \"./model.pth\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 143,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0., 1., 0., ..., 1., 0., 1.], dtype=float32)"
      ]
     },
     "execution_count": 143,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}