tcr_cnn_batch.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import random\n",
    "import data_processing\n",
    "import models\n",
    "import plotting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "chain = 'chain1'\n",
    "select_epitopes = ['NLVPMVATV', 'GILGFVFTL', 'GLCTLVAML', 'LLWNGPMAV']\n",
    "\n",
    "sub, other_ep_tcrs = data_processing.load_and_prepare_data(iedb_csv='./data/iedb_HLA-A2_1.csv', chain=chain,\n",
    "                                           select_epitopes=select_epitopes, show_summary=False, V_and_J_only=True)\n",
    "if chain == 'full_chain':\n",
    "    sub = sub[[chain, 'chain1', 'chain2', 'epitope', 'bind', 'set']]\n",
    "    other_ep_tcrs = other_ep_tcrs[[chain, 'chain1', 'chain2', 'epitope', 'bind', 'set']]\n",
    "else:\n",
    "    sub = sub[[chain, 'epitope', 'bind', 'set']]\n",
    "    other_ep_tcrs = other_ep_tcrs[[chain, 'epitope', 'bind', 'set']]    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Decide on whether or not to subset the chain1 or chain2 analysis based on what's available for the full_chain analysis.\n",
    "Looking at the full_chain requires samples that have both alpha and beta sequences. This significantly reduces the number of overall samples."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "use_same_tcrs_as_full_chain = False\n",
    "\n",
    "if use_same_tcrs_as_full_chain:\n",
    "    ### First run the full_chain analysis and save the alpha and beta chains\n",
    "    if chain == 'full_chain':\n",
    "        alpha = sub['chain1']\n",
    "        beta = sub['chain2']\n",
    "        beta_others = other_ep_tcrs['chain2']\n",
    "        alpha_others = other_ep_tcrs['chain1']\n",
    "\n",
    "    ### For chain1 analysis, subset on alpha chains from the full_chain analysis\n",
    "    elif chain == 'chain1':\n",
    "        sub = sub[sub[chain].isin(alpha.values)]\n",
    "        other_ep_tcrs = other_ep_tcrs[other_ep_tcrs[chain].isin(alpha_others.values)]\n",
    "\n",
    "    ### For chain2 analysis, subset on beta chains from the full_chain analysis\n",
    "    elif chain == 'chain2':\n",
    "        sub = sub[sub[chain].isin(beta.values)]\n",
    "        other_ep_tcrs = other_ep_tcrs[other_ep_tcrs[chain].isin(beta_others.values)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train / test split"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Save model results\n",
    "dict_of_models = {}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model_results = models.BatchModelResults()\n",
    "dict_of_models[chain] = model_results\n",
    "\n",
    "tcr_max_len = sub[chain].str.len().max()\n",
    "epitopes = sub['epitope']\n",
    "encoded_epitopes = data_processing.encode_epitopes(sub['epitope'].unique())\n",
    "\n",
    "for seed in range(20,27):\n",
    "    train, test = models.train_test_split(sub, classes=sub['epitope'].unique(), chain=chain, class_col='epitope', seed=seed, train_size=0.7,\n",
    "                                          put_dups_in_train=True)\n",
    "    other_ep_train, other_ep_test = models.train_test_split(other_ep_tcrs, classes=other_ep_tcrs['epitope'].unique(), chain=chain,\n",
    "                                                            class_col='epitope', seed=seed, train_size=0.7 , put_dups_in_train=True)\n",
    "    \n",
    "    ### Make synthetic train data, assign false epitopes to 'other_epitopes', combine all train data\n",
    "    ### Synthetic types: 'swapped_epitopes', 'scrambled_seqs_and_random_eps', 'scrambled_seqs_and_same_eps'\n",
    "    synthetics_to_make = ['swapped_epitopes']\n",
    "    synth_train = data_processing.make_synthetic_data(train, chain, epitopes, synthetics_to_make, tile_sequences=False)\n",
    "    other_ep_train['epitope'] = other_ep_train['epitope'].apply(lambda current_ep: data_processing.swap_epitopes(current_ep, epitopes))\n",
    "    all_train = pd.concat([train, synth_train, other_ep_train], sort=True)\n",
    "    all_train = all_train.sample(frac=1.0).reset_index(drop=True)\n",
    "\n",
    "    \n",
    "    ### Make synthetic test data, assign false epitopes to 'other_epitopes', combine all test data\n",
    "    synthetics_to_make = ['swapped_epitopes']\n",
    "    synth_test = data_processing.make_synthetic_data(test, chain, epitopes, synthetics_to_make)\n",
    "    other_ep_test['epitope'] = other_ep_train['epitope'].apply(lambda current_ep: data_processing.swap_epitopes(current_ep, epitopes))\n",
    "    all_test = pd.concat([test, synth_test, other_ep_test], sort=True)\n",
    "    all_test = all_test.sample(frac=1.0).reset_index(drop=True)\n",
    "\n",
    "    print('Training samples: {} || Testing samples: {}'.format(len(all_train), len(all_test)))\n",
    "\n",
    "    ### Encode sequences and epitopes\n",
    "    print('Encoding sequences...')\n",
    "    scrambled_train_tcrs = all_train[chain].sample(frac=1)\n",
    "    scrambled_test_tcrs = all_test[chain].sample(frac=1)\n",
    "    randomized_test_tcrs = all_test[chain].apply(lambda x: data_processing.randomize_seq(x))\n",
    "\n",
    "    \n",
    "    ### One hot encode sequences\n",
    "    train_tcrs = data_processing.encode_seqs(all_train[chain], target_length=tcr_max_len)\n",
    "    test_tcrs = data_processing.encode_seqs(all_test[chain], target_length=tcr_max_len)\n",
    "    scrambled_train_tcrs = data_processing.encode_seqs(scrambled_train_tcrs, target_length=tcr_max_len)\n",
    "    scrambled_test_tcrs = data_processing.encode_seqs(scrambled_test_tcrs, target_length=tcr_max_len)\n",
    "    randomized_test_tcrs = data_processing.encode_seqs(randomized_test_tcrs, target_length=tcr_max_len)\n",
    "    \n",
    "    \n",
    "    ### One hot encode epitopes\n",
    "    all_train['oh_ep'] = all_train['epitope'].apply(lambda x: encoded_epitopes[x])\n",
    "    all_test['oh_ep'] = all_test['epitope'].apply(lambda x: encoded_epitopes[x])\n",
    "    train_eps = np.vstack(all_train['oh_ep'])\n",
    "    test_eps = np.vstack(all_test['oh_ep'])\n",
    "\n",
    "    \n",
    "    ### Set up class weights - optional\n",
    "    all_train['weight'] = all_train['bind'].apply(lambda x: 1 if x == 0 else 1)\n",
    "    class_weights = dict(zip(all_train.index, all_train['weight']))\n",
    "    \n",
    "    \n",
    "    ### Modeling\n",
    "    print('Training...')\n",
    "    model, callbacks = models.build_model(train_tcrs, train_eps)\n",
    "    model.fit([train_tcrs, train_eps], all_train['bind'].values, epochs=20, class_weight=class_weights, validation_split=0.1, callbacks=callbacks, verbose=0)\n",
    "\n",
    "    \n",
    "    ### Model evaluation\n",
    "    all_test['prob'] = model.predict([test_tcrs, test_eps])\n",
    "    all_test['pred'] = all_test['prob'].apply(lambda x: 1 if x >= 0.5 else 0)\n",
    "    model_results.add_test_output(all_test)\n",
    "\n",
    "model_results.plot_all_auroc()\n",
    "model_results.plot_set_accuracy()\n",
    "model_results.plot_metrics_boxplots()\n",
    "model_results.plot_epitope_performance()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Alpha chain performance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dict_of_models['chain1'].plot_all_auroc()\n",
    "dict_of_models['chain1'].plot_set_accuracy()\n",
    "dict_of_models['chain1'].plot_metrics_boxplots()\n",
    "dict_of_models['chain1'].plot_epitope_performance()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Beta chain performance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dict_of_models['chain2'].plot_all_auroc()\n",
    "dict_of_models['chain2'].plot_set_accuracy()\n",
    "dict_of_models['chain2'].plot_metrics_boxplots()\n",
    "dict_of_models['chain2'].plot_epitope_performance()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Full chain performance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dict_of_models['full_chain'].plot_all_auroc()\n",
    "dict_of_models['full_chain'].plot_set_accuracy()\n",
    "dict_of_models['full_chain'].plot_metrics_boxplots()\n",
    "dict_of_models['full_chain'].plot_epitope_performance()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "____________________________"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot model architecture"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "models.plot_model(model, to_file='./images/model_architecture.png', show_layer_names=False, rankdir='LR')"
   ]
  }
 ],
 "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.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}