(D) ahpC and ahpC promoter analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# This notebook was created to 1. collect (SNP & INDEL) variants in *ahpC* and the *ahpC* promoter and 2. detect instances of co-occurence of LoF variants in *ahpC* and the *ahpC* promoter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style>.container { width:100% !important; }</style>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from IPython.core.display import display, HTML\n",
    "display(HTML(\"<style>.container { width:100% !important; }</style>\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import vcf\n",
    "\n",
    "%matplotlib inline\n",
    "import os\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "import matplotlib.ticker as ticker\n",
    "from itertools import compress\n",
    "from pylab import MaxNLocator\n",
    "import seaborn as sns; sns.set()\n",
    "from matplotlib.colors import LogNorm\n",
    "from matplotlib import gridspec\n",
    "import ast\n",
    "import itertools\n",
    "import seaborn as sns\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "import fastcluster\n",
    "from sklearn import cluster, datasets\n",
    "import scipy.cluster.hierarchy as hier\n",
    "from sklearn.cluster import KMeans\n",
    "import time\n",
    "import sys\n",
    "import pickle\n",
    "\n",
    "import Bio\n",
    "from Bio.Alphabet import IUPAC\n",
    "from Bio.Blast.Applications import NcbiblastnCommandline\n",
    "from Bio.Blast import NCBIXML\n",
    "from Bio.Seq import Seq\n",
    "from Bio.SeqRecord import SeqRecord\n",
    "from Bio.SeqFeature import SeqFeature, FeatureLocation\n",
    "from Bio import pairwise2\n",
    "from Bio import SeqIO\n",
    "from Bio.Graphics import GenomeDiagram\n",
    "from Bio.SeqUtils import GC\n",
    "from Bio import Phylo\n",
    "\n",
    "from Bio.Align.Applications import MuscleCommandline\n",
    "from StringIO import StringIO\n",
    "from Bio import AlignIO\n",
    "from Bio.Align import AlignInfo\n",
    "from Bio.Seq import MutableSeq\n",
    "import itertools\n",
    "import gzip\n",
    "\n",
    "import networkx as nx\n",
    "import scipy\n",
    "from collections import Counter"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# *ahpC*\n",
    "\n",
    "promoter: 2,726,088 - 2,726,192\n",
    "\n",
    "coding: 2,726,193 - 2,726,780\n",
    "\n",
    "H37Rv region: 2,726,088 - 2,726,780"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "H37Rv_region_start = 2726088\n",
    "H37Rv_region_end = 2726780\n",
    "regulator_id = 'Rv2427Ac_Rv2428'\n",
    "regulon_id = 'Rv2428'\n",
    "gene_symbol = 'ahpC'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### *Function* to convert lineage-calls to new/simpler lineage call scheme"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def convert_lineage_calls(isolate_annotation_DF):\n",
    "\n",
    "    lineage_hierarchincal_to_lineage_map_df = pd.read_csv('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/CSV files/lineage_hierachical_to_lineage_map.csv')\n",
    "    lineage_hierarchincal_to_lineage_map_series = pd.Series(list(lineage_hierarchincal_to_lineage_map_df.loc[:, 'lineage']), index = lineage_hierarchincal_to_lineage_map_df.lineage_hierarchical)\n",
    "\n",
    "    lineage_call_update_list = []\n",
    "\n",
    "    for hierarchical_lineage_call in isolate_annotation_DF.lineage_call:\n",
    "\n",
    "        # if lineage call has an updated shorter call\n",
    "        if hierarchical_lineage_call in lineage_hierarchincal_to_lineage_map_series.index:\n",
    "            lineage_call_update_list.append(lineage_hierarchincal_to_lineage_map_series[hierarchical_lineage_call])\n",
    "\n",
    "        # if lineage call doesn't have an updated shorter call\n",
    "        else:\n",
    "            lineage_call_update_list.append(hierarchical_lineage_call)\n",
    "\n",
    "    # replace the hierarchical calls with the shorter ones\n",
    "    isolate_annotation_DF.loc[:,'lineage_call'] = lineage_call_update_list\n",
    "    \n",
    "    return isolate_annotation_DF"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [1] INDELs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [1.1] Load INDEL genotype matrix and Annotation Files"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 31,428 isolates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load isolate annotation file (columns of Genotype Matrix)\n",
    "INDEL_isolate_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape_indels/Genotypes_Filtered_2/genotypes_isolate_annotation.pkl')\n",
    "INDEL_isolate_annotation_DF = INDEL_isolate_annotation_DF.loc[:, ['isolate_ID','lineage_call']] #drop columns\n",
    "INDEL_isolate_annotation_DF = convert_lineage_calls(INDEL_isolate_annotation_DF)\n",
    "\n",
    "#load INDEL annotation file (rows of Genotype Matrix) with gene annotation information\n",
    "INDEL_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape_indels/Genotypes_Filtered_2/genotypes_INDEL_functional_annotation.pkl')\n",
    "INDEL_annotation_DF.reset_index(inplace = True, drop = False)\n",
    "\n",
    "#load Genotypes Matrix\n",
    "INDEL_genotypes_array =  np.load('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape_indels/Genotypes_Filtered_2/genotypes_matrix.npy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>SAMN13051687</td>\n",
       "      <td>2.2.1.1.1.i3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>SAMN09100245</td>\n",
       "      <td>4.3.i4.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>SAMN08732238</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>SAMN07658260</td>\n",
       "      <td>3.1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>SAMN03648003</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     isolate_ID  lineage_call\n",
       "0  SAMN13051687  2.2.1.1.1.i3\n",
       "1  SAMN09100245      4.3.i4.1\n",
       "2  SAMN08732238     2.2.1.1.1\n",
       "3  SAMN07658260         3.1.1\n",
       "4  SAMN03648003     2.2.1.1.1"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_isolate_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(31428, 2)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_isolate_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>ACCGACGAAG_313_A</td>\n",
       "      <td>313</td>\n",
       "      <td>ACCGACGAAG</td>\n",
       "      <td>A</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>313</td>\n",
       "      <td>del</td>\n",
       "      <td>inframe</td>\n",
       "      <td>105</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>TC_1549_T</td>\n",
       "      <td>1549</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv0001_Rv0002</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>T_1552_TAA</td>\n",
       "      <td>1552</td>\n",
       "      <td>T</td>\n",
       "      <td>TAA</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv0001_Rv0002</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>TAA_1552_T</td>\n",
       "      <td>1552</td>\n",
       "      <td>TAA</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv0001_Rv0002</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>T_1552_TA</td>\n",
       "      <td>1552</td>\n",
       "      <td>T</td>\n",
       "      <td>TA</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv0001_Rv0002</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                key   pos         ref  alt gene_category gene_name  \\\n",
       "0  ACCGACGAAG_313_A   313  ACCGACGAAG    A     Essential      dnaA   \n",
       "1         TC_1549_T  1549          TC    T          None      None   \n",
       "2        T_1552_TAA  1552           T  TAA          None      None   \n",
       "3        TAA_1552_T  1552         TAA    T          None      None   \n",
       "4         T_1552_TA  1552           T   TA          None      None   \n",
       "\n",
       "         gene_id gene_pos ins_del  INDEL_type codon_pos  \n",
       "0         Rv0001      313     del     inframe       105  \n",
       "1  Rv0001_Rv0002     None     del  frameshift      None  \n",
       "2  Rv0001_Rv0002     None     ins  frameshift      None  \n",
       "3  Rv0001_Rv0002     None     del  frameshift      None  \n",
       "4  Rv0001_Rv0002     None     ins  frameshift      None  "
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50260, 11)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       ...,\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0]], dtype=int8)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_genotypes_array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50260, 31428)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### additional 12 eis C-14T mutants with AG MICs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load isolate annotation file (columns of Genotype Matrix)\n",
    "isolate_annotation_DF_extra_strains = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/eis_promoter_mutant_strains_from_Milan/Genotypes/Genotypes_Filtered/genotypes_isolate_annotation.pkl')\n",
    "isolate_annotation_DF_extra_strains = isolate_annotation_DF_extra_strains.loc[:, ['isolate_ID','lineage_call']]\n",
    "isolate_annotation_DF_extra_strains = convert_lineage_calls(isolate_annotation_DF_extra_strains)\n",
    "\n",
    "#load Genotypes Matrix\n",
    "genotypes_array_extra_strains =  np.load('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/eis_promoter_mutant_strains_from_Milan/Genotypes_indels/Genotypes_Filtered/genotypes_matrix.npy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>168-19</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>622-19</td>\n",
       "      <td>2.2.1.1.1.i3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>655-19</td>\n",
       "      <td>4.3.i3.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>IT1070</td>\n",
       "      <td>2.2.1.1.1.i3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>IT123</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  isolate_ID  lineage_call\n",
       "0     168-19     2.2.1.1.1\n",
       "1     622-19  2.2.1.1.1.i3\n",
       "2     655-19      4.3.i3.1\n",
       "3     IT1070  2.2.1.1.1.i3\n",
       "4      IT123     2.2.1.1.1"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isolate_annotation_DF_extra_strains.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       ...,\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 9, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0]], dtype=int8)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "genotypes_array_extra_strains"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50260, 12)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(genotypes_array_extra_strains)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Merge *isolate annotation files* and *genotypes matrices* together"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merge isolate annotation dataframes\n",
    "INDEL_isolate_annotation_DF = INDEL_isolate_annotation_DF.append(isolate_annotation_DF_extra_strains)\n",
    "INDEL_isolate_annotation_DF.reset_index(inplace = True, drop = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(31440, 2)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_isolate_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merge genotypes matrices\n",
    "INDEL_genotypes_array = np.hstack((INDEL_genotypes_array, genotypes_array_extra_strains))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50260, 31440)"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [1.2] Subset to targeted chomosomal region"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "# keep only indels that occur within chromosomal region\n",
    "indels_in_targeted_region = list(INDEL_annotation_DF[(INDEL_annotation_DF.pos >= H37Rv_region_start) & (INDEL_annotation_DF.pos <= H37Rv_region_end)].key)\n",
    "\n",
    "# construct boolean filter\n",
    "indels_to_keep_filter = [indel_i in indels_in_targeted_region for indel_i in INDEL_annotation_DF.key]\n",
    "\n",
    "# filter Genotypes Matrix\n",
    "INDEL_genotypes_array = INDEL_genotypes_array[indels_to_keep_filter , :]\n",
    "\n",
    "# filter INDEL annotation file\n",
    "INDEL_annotation_DF = INDEL_annotation_DF[indels_to_keep_filter]\n",
    "INDEL_annotation_DF.reset_index(drop = True , inplace = True) #re-index new filtered INDEL annotation DF (so new index matches indexing of genotypes matrix rows)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A_2726096_ATATG</td>\n",
       "      <td>2726096</td>\n",
       "      <td>A</td>\n",
       "      <td>ATATG</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>G_2726105_GAT</td>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>GAT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>CT_2726115_C</td>\n",
       "      <td>2726115</td>\n",
       "      <td>CT</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>G_2726146_GT</td>\n",
       "      <td>2726146</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>AT_2726148_A</td>\n",
       "      <td>2726148</td>\n",
       "      <td>AT</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               key      pos ref    alt gene_category gene_name  \\\n",
       "0  A_2726096_ATATG  2726096   A  ATATG          None      None   \n",
       "1    G_2726105_GAT  2726105   G    GAT          None      None   \n",
       "2     CT_2726115_C  2726115  CT      C          None      None   \n",
       "3     G_2726146_GT  2726146   G     GT          None      None   \n",
       "4     AT_2726148_A  2726148  AT      A          None      None   \n",
       "\n",
       "           gene_id gene_pos ins_del  INDEL_type codon_pos  \n",
       "0  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "1  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "2  Rv2427Ac_Rv2428     None     del  frameshift      None  \n",
       "3  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "4  Rv2427Ac_Rv2428     None     del  frameshift      None  "
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(14, 11)"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       ...,\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0],\n",
       "       [0, 0, 0, ..., 0, 0, 0]], dtype=int8)"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_genotypes_array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(14, 31440)"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [2] SNPs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 31,428 isolates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [2.1] Load SNP genotype matrix and Annotation Files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load isolate annotation file (columns of Genotype Matrix)\n",
    "SNP_isolate_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape/Genotypes_Filtered_2/genotypes_isolate_annotation.pkl')\n",
    "SNP_isolate_annotation_DF = SNP_isolate_annotation_DF.loc[:, ['isolate_ID','lineage_call']] #drop columns\n",
    "SNP_isolate_annotation_DF = convert_lineage_calls(SNP_isolate_annotation_DF)\n",
    "\n",
    "#load INDEL annotation file (rows of Genotype Matrix) with gene annotation information\n",
    "SNP_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape/Genotypes_Filtered_2/genotypes_SNP_annotation.pkl')\n",
    "SNP_annotation_DF = SNP_annotation_DF.loc[:, ['pos','ref','alt']] # drop columns\n",
    "SNP_annotation_DF.reset_index(inplace = True, drop = False)\n",
    "\n",
    "#load Genotypes Matrix\n",
    "SNP_genotypes_array =  np.load('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape/Genotypes_Filtered_2/genotypes_matrix.npy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>SAMEA3558733</td>\n",
       "      <td>4.3.i4.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>SAMN03648641</td>\n",
       "      <td>4.4.1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>SAMN03647419</td>\n",
       "      <td>3.1.1.i1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>SAMEA3671418</td>\n",
       "      <td>4.3.i2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>SAMN07659096</td>\n",
       "      <td>1.1.3</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     isolate_ID lineage_call\n",
       "0  SAMEA3558733     4.3.i4.2\n",
       "1  SAMN03648641      4.4.1.1\n",
       "2  SAMN03647419     3.1.1.i1\n",
       "3  SAMEA3671418       4.3.i2\n",
       "4  SAMN07659096        1.1.3"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_isolate_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(31428, 2)"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_isolate_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>index</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>48</td>\n",
       "      <td>C</td>\n",
       "      <td>[T]</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>64</td>\n",
       "      <td>G</td>\n",
       "      <td>[C]</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>67</td>\n",
       "      <td>G</td>\n",
       "      <td>[A, T]</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>69</td>\n",
       "      <td>C</td>\n",
       "      <td>[T]</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>71</td>\n",
       "      <td>C</td>\n",
       "      <td>[T]</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   index  pos ref     alt\n",
       "0      0   48   C     [T]\n",
       "1      1   64   G     [C]\n",
       "2      2   67   G  [A, T]\n",
       "3      3   69   C     [T]\n",
       "4      4   71   C     [T]"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(782565, 4)"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 1, 1, ..., 1, 1, 1],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       ...,\n",
       "       [1, 9, 9, ..., 1, 1, 1],\n",
       "       [1, 9, 1, ..., 1, 1, 1],\n",
       "       [1, 9, 1, ..., 1, 1, 1]], dtype=int8)"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_genotypes_array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(782565, 31428)"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### additional 12 eis C-14T mutants with AG MICs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load isolate annotation file (columns of Genotype Matrix)\n",
    "isolate_annotation_DF_extra_strains = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/eis_promoter_mutant_strains_from_Milan/Genotypes/Genotypes_Filtered/genotypes_isolate_annotation.pkl')\n",
    "isolate_annotation_DF_extra_strains = isolate_annotation_DF_extra_strains.loc[:, ['isolate_ID','lineage_call']]\n",
    "isolate_annotation_DF_extra_strains = convert_lineage_calls(isolate_annotation_DF_extra_strains)\n",
    "\n",
    "#load Genotypes Matrix\n",
    "genotypes_array_extra_strains =  np.load('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/eis_promoter_mutant_strains_from_Milan/Genotypes/Genotypes_Filtered/genotypes_matrix.npy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>168-19</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>622-19</td>\n",
       "      <td>2.2.1.1.1.i3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>655-19</td>\n",
       "      <td>4.3.i3.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>IT1070</td>\n",
       "      <td>2.2.1.1.1.i3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>IT123</td>\n",
       "      <td>2.2.1.1.1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  isolate_ID  lineage_call\n",
       "0     168-19     2.2.1.1.1\n",
       "1     622-19  2.2.1.1.1.i3\n",
       "2     655-19      4.3.i3.1\n",
       "3     IT1070  2.2.1.1.1.i3\n",
       "4      IT123     2.2.1.1.1"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isolate_annotation_DF_extra_strains.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 1, 1, ..., 1, 1, 1],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       ...,\n",
       "       [1, 1, 1, ..., 1, 1, 1],\n",
       "       [1, 1, 1, ..., 1, 1, 1],\n",
       "       [1, 1, 1, ..., 1, 1, 1]], dtype=int8)"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "genotypes_array_extra_strains"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(782565, 12)"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(genotypes_array_extra_strains)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Merge *isolate annotation files* and *genotypes matrices* together"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merge isolate annotation dataframes\n",
    "SNP_isolate_annotation_DF = SNP_isolate_annotation_DF.append(isolate_annotation_DF_extra_strains)\n",
    "SNP_isolate_annotation_DF.reset_index(inplace = True, drop = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(31440, 2)"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_isolate_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merge genotypes matrices\n",
    "SNP_genotypes_array = np.hstack((SNP_genotypes_array, genotypes_array_extra_strains))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(782565, 31440)"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [2.2] Subset to targeted chomosomal region"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "# keep only SNPs that occur within chromosomal region\n",
    "SNPs_in_targeted_region = list(SNP_annotation_DF[(SNP_annotation_DF.pos >= H37Rv_region_start) & (SNP_annotation_DF.pos <= H37Rv_region_end)].pos)\n",
    "\n",
    "# construct boolean filter\n",
    "SNPs_to_keep_filter = [SNP_i in SNPs_in_targeted_region for SNP_i in SNP_annotation_DF.pos]\n",
    "\n",
    "# filter Genotypes Matrix\n",
    "SNP_genotypes_array = SNP_genotypes_array[SNPs_to_keep_filter , :]\n",
    "\n",
    "# filter SNP annotation file\n",
    "SNP_annotation_DF = SNP_annotation_DF[SNPs_to_keep_filter]\n",
    "SNP_annotation_DF.reset_index(drop = True , inplace = True) #re-index new filtered SNP annotation DF (so new index matches indexing of genotypes matrix rows)\n",
    "SNP_annotation_DF.loc[:,'array_index'] = SNP_annotation_DF.index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>index</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>array_index</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>491723</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>[A]</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>491724</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>[T]</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>491725</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>[A]</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>491726</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>[A]</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>491727</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>[C]</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    index      pos ref  alt  array_index\n",
       "0  491723  2726088   G  [A]            0\n",
       "1  491724  2726089   C  [T]            1\n",
       "2  491725  2726091   G  [A]            2\n",
       "3  491726  2726093   T  [A]            3\n",
       "4  491727  2726095   A  [C]            4"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(133, 5)"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Re-format SNP annotation DataFrame so that each alternate allele is its own row"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [],
   "source": [
    "SNP_annotation_DF_mod = pd.DataFrame(columns = SNP_annotation_DF.columns)\n",
    "new_SNP_index = 0\n",
    "for SNP_i_row in SNP_annotation_DF.index:\n",
    "    \n",
    "    SNP_i_info = SNP_annotation_DF.loc[SNP_i_row, :]\n",
    "    \n",
    "    #iterate through alternate alleles, store a row for each one\n",
    "    for SNP_i_alt_allele in SNP_i_info.alt:\n",
    "        \n",
    "        SNP_annotation_DF_mod.loc[new_SNP_index, :] = list(SNP_i_info.loc[['index','pos','ref']]) + [SNP_i_alt_allele] + list(SNP_i_info.loc[['array_index']])\n",
    "        new_SNP_index += 1\n",
    "\n",
    "SNP_annotation_DF = SNP_annotation_DF_mod\n",
    "#create new column to merge func annotation on\n",
    "SNP_annotation_DF.loc[:, 'key'] = [str(SNP_i_pos) + '_' + SNP_i_alt for SNP_i_pos, SNP_i_alt in zip(SNP_annotation_DF.pos, SNP_annotation_DF.alt)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>index</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>array_index</th>\n",
       "      <th>key</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>491723</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>0</td>\n",
       "      <td>2726088_A</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>491724</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>1</td>\n",
       "      <td>2726089_T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>491725</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>2</td>\n",
       "      <td>2726091_A</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>491726</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>3</td>\n",
       "      <td>2726093_A</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>491727</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>4</td>\n",
       "      <td>2726095_C</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    index      pos ref alt array_index        key\n",
       "0  491723  2726088   G   A           0  2726088_A\n",
       "1  491724  2726089   C   T           1  2726089_T\n",
       "2  491725  2726091   G   A           2  2726091_A\n",
       "3  491726  2726093   T   A           3  2726093_A\n",
       "4  491727  2726095   A   C           4  2726095_C"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(153, 6)"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Merge functional annotation data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "SNP_func_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape/Genotypes_Filtered_2/genotypes_SNP_functional_annotation.pkl')\n",
    "SNP_func_annotation_DF.loc[:, 'key'] = SNP_func_annotation_DF.index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>SNP_type</th>\n",
       "      <th>AA_change</th>\n",
       "      <th>key</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>48_T</th>\n",
       "      <td>48</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>48</td>\n",
       "      <td>S</td>\n",
       "      <td>V16V</td>\n",
       "      <td>48_T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>64_C</th>\n",
       "      <td>64</td>\n",
       "      <td>G</td>\n",
       "      <td>C</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>64</td>\n",
       "      <td>N</td>\n",
       "      <td>G22R</td>\n",
       "      <td>64_C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>67_A</th>\n",
       "      <td>67</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>67</td>\n",
       "      <td>N</td>\n",
       "      <td>D23N</td>\n",
       "      <td>67_A</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>67_T</th>\n",
       "      <td>67</td>\n",
       "      <td>G</td>\n",
       "      <td>T</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>67</td>\n",
       "      <td>N</td>\n",
       "      <td>D23Y</td>\n",
       "      <td>67_T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>69_T</th>\n",
       "      <td>69</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>Essential</td>\n",
       "      <td>dnaA</td>\n",
       "      <td>Rv0001</td>\n",
       "      <td>69</td>\n",
       "      <td>S</td>\n",
       "      <td>D23D</td>\n",
       "      <td>69_T</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      pos ref alt gene_category gene_name gene_id gene_pos SNP_type AA_change  \\\n",
       "48_T   48   C   T     Essential      dnaA  Rv0001       48        S      V16V   \n",
       "64_C   64   G   C     Essential      dnaA  Rv0001       64        N      G22R   \n",
       "67_A   67   G   A     Essential      dnaA  Rv0001       67        N      D23N   \n",
       "67_T   67   G   T     Essential      dnaA  Rv0001       67        N      D23Y   \n",
       "69_T   69   C   T     Essential      dnaA  Rv0001       69        S      D23D   \n",
       "\n",
       "       key  \n",
       "48_T  48_T  \n",
       "64_C  64_C  \n",
       "67_A  67_A  \n",
       "67_T  67_T  \n",
       "69_T  69_T  "
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_func_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [],
   "source": [
    "SNP_annotation_DF = SNP_annotation_DF.merge(SNP_func_annotation_DF.loc[: , ['gene_category','gene_name','gene_id','gene_pos','SNP_type','AA_change','key']], how = 'left', on = 'key')\n",
    "SNP_annotation_DF = SNP_annotation_DF.loc[:, ['key','pos','ref','alt','gene_category','gene_name','gene_id','gene_pos','SNP_type','AA_change','array_index']] #re-arrange columns\n",
    "SNP_annotation_DF.rename(columns={'SNP_type':'mut_type', 'AA_change':'codon_pos'}, inplace = True) #change column names"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>array_index</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726088_A</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726089_T</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726091_A</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726093_A</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2726095_C</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         key      pos ref alt gene_category gene_name          gene_id  \\\n",
       "0  2726088_A  2726088   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "1  2726089_T  2726089   C   T          None      None  Rv2427Ac_Rv2428   \n",
       "2  2726091_A  2726091   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "3  2726093_A  2726093   T   A          None      None  Rv2427Ac_Rv2428   \n",
       "4  2726095_C  2726095   A   C          None      None  Rv2427Ac_Rv2428   \n",
       "\n",
       "  gene_pos mut_type codon_pos array_index  \n",
       "0     None        I      None           0  \n",
       "1     None        I      None           1  \n",
       "2     None        I      None           2  \n",
       "3     None        I      None           3  \n",
       "4     None        I      None           4  "
      ]
     },
     "execution_count": 51,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(153, 11)"
      ]
     },
     "execution_count": 52,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[2, 2, 2, ..., 2, 2, 2],\n",
       "       [1, 1, 1, ..., 1, 1, 1],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       ...,\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       [2, 2, 2, ..., 2, 2, 2],\n",
       "       [1, 1, 1, ..., 1, 1, 1]], dtype=int8)"
      ]
     },
     "execution_count": 53,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_genotypes_array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(133, 31440)"
      ]
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_genotypes_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [3] Pull info for isolates with mutant alleles"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [3.1] Get genotypes for SNPs of interest"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### *Function* to get the genotypes for a specifc SNP and get attributes for isolates that support SNP"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_genotypes_for_SNP(SNP_i_pos, SNP_i_alt_allele):\n",
    "    \n",
    "    #get the number code that corresponds to the alt allele, 9 indicated bad quality call\n",
    "    base_code_dict = {'A':0, 'C':1, 'G':2, 'T':3}\n",
    "    SNP_i_alt_allele_code = base_code_dict[SNP_i_alt_allele]\n",
    "    SNP_i_genotypes = SNP_genotypes_array[SNP_annotation_DF[SNP_annotation_DF.pos == SNP_i_pos].array_index.values[0] , :] \n",
    "\n",
    "    #count the number of isolates that support this SNP call\n",
    "    num_isolates_SNP_i = np.sum(SNP_i_genotypes == SNP_i_alt_allele_code)\n",
    "\n",
    "    #create a boolean filter and extract information for the isolates that support this call\n",
    "    SNP_i_isolate_filter = SNP_i_genotypes == SNP_i_alt_allele_code\n",
    "    isolate_annotation_SNP_i = SNP_isolate_annotation_DF[SNP_i_isolate_filter]\n",
    "\n",
    "    #find number of sub-lineages that have SNP support in at least 1 isolate\n",
    "    num_sublineages_with_SNP_i = len(set(list(isolate_annotation_SNP_i.lineage_call.values)))\n",
    "    \n",
    "    #get list of sublineages w/ at least 1 isolates that supports SNP call\n",
    "    sublineages_with_SNP_count_dict = Counter(isolate_annotation_SNP_i.lineage_call.values)\n",
    "    sublineages_with_SNP_list = []\n",
    "    for sublineage_i in sublineages_with_SNP_count_dict.keys():\n",
    "\n",
    "        sublineage_i_with_SNP = sublineage_i + '({0})'.format(str(sublineages_with_SNP_count_dict[sublineage_i]))\n",
    "        sublineages_with_SNP_list.append(sublineage_i_with_SNP)\n",
    "\n",
    "    sublineages_with_SNP = ' - '.join(sublineages_with_SNP_list)\n",
    "\n",
    "    return [num_isolates_SNP_i , num_sublineages_with_SNP_i , sublineages_with_SNP]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>array_index</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726088_A</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726089_T</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726091_A</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726093_A</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2726095_C</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         key      pos ref alt gene_category gene_name          gene_id  \\\n",
       "0  2726088_A  2726088   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "1  2726089_T  2726089   C   T          None      None  Rv2427Ac_Rv2428   \n",
       "2  2726091_A  2726091   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "3  2726093_A  2726093   T   A          None      None  Rv2427Ac_Rv2428   \n",
       "4  2726095_C  2726095   A   C          None      None  Rv2427Ac_Rv2428   \n",
       "\n",
       "  gene_pos mut_type codon_pos array_index  \n",
       "0     None        I      None           0  \n",
       "1     None        I      None           1  \n",
       "2     None        I      None           2  \n",
       "3     None        I      None           3  \n",
       "4     None        I      None           4  "
      ]
     },
     "execution_count": 56,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(153, 11)"
      ]
     },
     "execution_count": 57,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_isolates_with_SNP = []\n",
    "num_sublineages_with_SNP = []\n",
    "sublineages_with_SNP = []\n",
    "\n",
    "for SNP_i_pos, SNP_i_alt_allele in zip(SNP_annotation_DF.pos, SNP_annotation_DF.alt):\n",
    "    \n",
    "    num_isolates_with_SNP_i , num_sublineages_with_SNP_i , sublineages_with_SNP_i = get_genotypes_for_SNP(SNP_i_pos, SNP_i_alt_allele)\n",
    "    \n",
    "    num_isolates_with_SNP.append(num_isolates_with_SNP_i)\n",
    "    num_sublineages_with_SNP.append(num_sublineages_with_SNP_i)\n",
    "    sublineages_with_SNP.append(sublineages_with_SNP_i)\n",
    "    \n",
    "SNP_annotation_DF.loc[: , 'num_isolates'] = num_isolates_with_SNP\n",
    "SNP_annotation_DF.loc[: , 'num_sublineages'] = num_sublineages_with_SNP\n",
    "SNP_annotation_DF.loc[: , 'sublineages'] = sublineages_with_SNP\n",
    "\n",
    "#drop SNPs present in 0 isolates\n",
    "SNP_annotation_DF = SNP_annotation_DF[SNP_annotation_DF.num_isolates > 0]\n",
    "SNP_annotation_DF.reset_index(inplace = True, drop = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>array_index</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726088_A</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.1.1.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726089_T</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726091_A</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1.1.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726093_A</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>6(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2726095_C</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(2)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         key      pos ref alt gene_category gene_name          gene_id  \\\n",
       "0  2726088_A  2726088   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "1  2726089_T  2726089   C   T          None      None  Rv2427Ac_Rv2428   \n",
       "2  2726091_A  2726091   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "3  2726093_A  2726093   T   A          None      None  Rv2427Ac_Rv2428   \n",
       "4  2726095_C  2726095   A   C          None      None  Rv2427Ac_Rv2428   \n",
       "\n",
       "  gene_pos mut_type codon_pos array_index  num_isolates  num_sublineages  \\\n",
       "0     None        I      None           0             2                1   \n",
       "1     None        I      None           1             1                1   \n",
       "2     None        I      None           2             1                1   \n",
       "3     None        I      None           3             2                1   \n",
       "4     None        I      None           4             2                1   \n",
       "\n",
       "       sublineages  \n",
       "0     2.2.1.1.1(2)  \n",
       "1  4.1.i1.1.1.1(1)  \n",
       "2         1.1.2(1)  \n",
       "3             6(2)  \n",
       "4       4.10.i1(2)  "
      ]
     },
     "execution_count": 59,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SNP_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(151, 14)"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(SNP_annotation_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [3.2] Get genotypes for indels of interest"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### *Function* to get the genotypes for a specifc indel and get attributes for isolates that support indel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_genotypes_for_indel(indel_i_key):\n",
    "    \n",
    "    #0 supports Ref, 1 supports Alt, 9 indicated bad quality call\n",
    "    indel_i_genotypes = INDEL_genotypes_array[INDEL_annotation_DF[INDEL_annotation_DF.key == indel_i_key].index[0] , :] \n",
    "\n",
    "    #count the number of isolates that support this indel call\n",
    "    num_isolates_indel_i = np.sum(indel_i_genotypes == 1)\n",
    "\n",
    "    #create a boolean filter and extract information for the isolates taht support this call\n",
    "    indel_i_isolate_filter = indel_i_genotypes == 1\n",
    "    isolate_annotation_indel_i = INDEL_isolate_annotation_DF[indel_i_isolate_filter]\n",
    "\n",
    "    #find number of sub-lineages that have indel support in at least 1 isolate\n",
    "    num_sublineages_with_indel_i = len(set(list(isolate_annotation_indel_i.lineage_call.values)))\n",
    "    \n",
    "    #get list of sublineages w/ at least 1 isolates that supports indel call\n",
    "    sublineage_with_indel_count_dict = Counter(isolate_annotation_indel_i.lineage_call.values)\n",
    "    sublineages_with_indel_list = []\n",
    "    for sublineage_i in sublineage_with_indel_count_dict.keys():\n",
    "\n",
    "        sublineage_i_with_indel = sublineage_i + '({0})'.format(str(sublineage_with_indel_count_dict[sublineage_i]))\n",
    "        sublineages_with_indel_list.append(sublineage_i_with_indel)\n",
    "\n",
    "    sublineages_with_indel = ' - '.join(sublineages_with_indel_list)\n",
    "    \n",
    "    return [num_isolates_indel_i , num_sublineages_with_indel_i , sublineages_with_indel]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A_2726096_ATATG</td>\n",
       "      <td>2726096</td>\n",
       "      <td>A</td>\n",
       "      <td>ATATG</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>G_2726105_GAT</td>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>GAT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>CT_2726115_C</td>\n",
       "      <td>2726115</td>\n",
       "      <td>CT</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>G_2726146_GT</td>\n",
       "      <td>2726146</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>AT_2726148_A</td>\n",
       "      <td>2726148</td>\n",
       "      <td>AT</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               key      pos ref    alt gene_category gene_name  \\\n",
       "0  A_2726096_ATATG  2726096   A  ATATG          None      None   \n",
       "1    G_2726105_GAT  2726105   G    GAT          None      None   \n",
       "2     CT_2726115_C  2726115  CT      C          None      None   \n",
       "3     G_2726146_GT  2726146   G     GT          None      None   \n",
       "4     AT_2726148_A  2726148  AT      A          None      None   \n",
       "\n",
       "           gene_id gene_pos ins_del  INDEL_type codon_pos  \n",
       "0  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "1  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "2  Rv2427Ac_Rv2428     None     del  frameshift      None  \n",
       "3  Rv2427Ac_Rv2428     None     ins  frameshift      None  \n",
       "4  Rv2427Ac_Rv2428     None     del  frameshift      None  "
      ]
     },
     "execution_count": 62,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(14, 11)"
      ]
     },
     "execution_count": 63,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_annotation_DF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_isolates_with_indel = []\n",
    "num_sublineages_with_indel = []\n",
    "sublineages_with_indel = []\n",
    "\n",
    "for indel_i_key in INDEL_annotation_DF.key:\n",
    "    \n",
    "    num_isolates_with_indel_i , num_sublineages_with_indel_i , sublineages_with_indel_i = get_genotypes_for_indel(indel_i_key)\n",
    "    \n",
    "    num_isolates_with_indel.append(num_isolates_with_indel_i)\n",
    "    num_sublineages_with_indel.append(num_sublineages_with_indel_i)\n",
    "    sublineages_with_indel.append(sublineages_with_indel_i)\n",
    "    \n",
    "INDEL_annotation_DF.loc[: , 'num_isolates'] = num_isolates_with_indel\n",
    "INDEL_annotation_DF.loc[: , 'num_sublineages'] = num_sublineages_with_indel\n",
    "INDEL_annotation_DF.loc[: , 'sublineages'] = sublineages_with_indel\n",
    "\n",
    "#drop indels present in 0 isolates\n",
    "INDEL_annotation_DF = INDEL_annotation_DF[INDEL_annotation_DF.num_isolates > 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A_2726096_ATATG</td>\n",
       "      <td>2726096</td>\n",
       "      <td>A</td>\n",
       "      <td>ATATG</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>G_2726105_GAT</td>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>GAT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>4</td>\n",
       "      <td>4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>CT_2726115_C</td>\n",
       "      <td>2726115</td>\n",
       "      <td>CT</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>10</td>\n",
       "      <td>3</td>\n",
       "      <td>3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>G_2726146_GT</td>\n",
       "      <td>2726146</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>37</td>\n",
       "      <td>4</td>\n",
       "      <td>4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>AT_2726148_A</td>\n",
       "      <td>2726148</td>\n",
       "      <td>AT</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.2(2)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               key      pos ref    alt gene_category gene_name  \\\n",
       "0  A_2726096_ATATG  2726096   A  ATATG          None      None   \n",
       "1    G_2726105_GAT  2726105   G    GAT          None      None   \n",
       "2     CT_2726115_C  2726115  CT      C          None      None   \n",
       "3     G_2726146_GT  2726146   G     GT          None      None   \n",
       "4     AT_2726148_A  2726148  AT      A          None      None   \n",
       "\n",
       "           gene_id gene_pos ins_del  INDEL_type codon_pos  num_isolates  \\\n",
       "0  Rv2427Ac_Rv2428     None     ins  frameshift      None             1   \n",
       "1  Rv2427Ac_Rv2428     None     ins  frameshift      None             4   \n",
       "2  Rv2427Ac_Rv2428     None     del  frameshift      None            10   \n",
       "3  Rv2427Ac_Rv2428     None     ins  frameshift      None            37   \n",
       "4  Rv2427Ac_Rv2428     None     del  frameshift      None             2   \n",
       "\n",
       "   num_sublineages                                        sublineages  \n",
       "0                1                                      4.1.i1.2.1(1)  \n",
       "1                4  4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...  \n",
       "2                3               3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)  \n",
       "3                4       4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)  \n",
       "4                1                                           2.2.2(2)  "
      ]
     },
     "execution_count": 65,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(14, 14)"
      ]
     },
     "execution_count": 66,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(INDEL_annotation_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [4] Check co-occurrence between *regulator* SNPs/INDELs & *regulon* SNPs/INDELs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### SNPs in regulator"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*NOTE:* **2726105_A** is the ahpC -88 mutation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_SNPs_DF = SNP_annotation_DF[SNP_annotation_DF.gene_id == regulator_id]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>array_index</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726088_A</td>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.1.1.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726089_T</td>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726091_A</td>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1.1.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726093_A</td>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>6(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2726095_C</td>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(2)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         key      pos ref alt gene_category gene_name          gene_id  \\\n",
       "0  2726088_A  2726088   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "1  2726089_T  2726089   C   T          None      None  Rv2427Ac_Rv2428   \n",
       "2  2726091_A  2726091   G   A          None      None  Rv2427Ac_Rv2428   \n",
       "3  2726093_A  2726093   T   A          None      None  Rv2427Ac_Rv2428   \n",
       "4  2726095_C  2726095   A   C          None      None  Rv2427Ac_Rv2428   \n",
       "\n",
       "  gene_pos mut_type codon_pos array_index  num_isolates  num_sublineages  \\\n",
       "0     None        I      None           0             2                1   \n",
       "1     None        I      None           1             1                1   \n",
       "2     None        I      None           2             1                1   \n",
       "3     None        I      None           3             2                1   \n",
       "4     None        I      None           4             2                1   \n",
       "\n",
       "       sublineages  \n",
       "0     2.2.1.1.1(2)  \n",
       "1  4.1.i1.1.1.1(1)  \n",
       "2         1.1.2(1)  \n",
       "3             6(2)  \n",
       "4       4.10.i1(2)  "
      ]
     },
     "execution_count": 68,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulator_SNPs_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(49, 14)"
      ]
     },
     "execution_count": 69,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulator_SNPs_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### INDELs in regulator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_INDELs_DF = INDEL_annotation_DF[INDEL_annotation_DF.gene_id == regulator_id]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A_2726096_ATATG</td>\n",
       "      <td>2726096</td>\n",
       "      <td>A</td>\n",
       "      <td>ATATG</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>G_2726105_GAT</td>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>GAT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>4</td>\n",
       "      <td>4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>CT_2726115_C</td>\n",
       "      <td>2726115</td>\n",
       "      <td>CT</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>10</td>\n",
       "      <td>3</td>\n",
       "      <td>3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>G_2726146_GT</td>\n",
       "      <td>2726146</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>37</td>\n",
       "      <td>4</td>\n",
       "      <td>4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>AT_2726148_A</td>\n",
       "      <td>2726148</td>\n",
       "      <td>AT</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.2(2)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "               key      pos ref    alt gene_category gene_name  \\\n",
       "0  A_2726096_ATATG  2726096   A  ATATG          None      None   \n",
       "1    G_2726105_GAT  2726105   G    GAT          None      None   \n",
       "2     CT_2726115_C  2726115  CT      C          None      None   \n",
       "3     G_2726146_GT  2726146   G     GT          None      None   \n",
       "4     AT_2726148_A  2726148  AT      A          None      None   \n",
       "\n",
       "           gene_id gene_pos ins_del  INDEL_type codon_pos  num_isolates  \\\n",
       "0  Rv2427Ac_Rv2428     None     ins  frameshift      None             1   \n",
       "1  Rv2427Ac_Rv2428     None     ins  frameshift      None             4   \n",
       "2  Rv2427Ac_Rv2428     None     del  frameshift      None            10   \n",
       "3  Rv2427Ac_Rv2428     None     ins  frameshift      None            37   \n",
       "4  Rv2427Ac_Rv2428     None     del  frameshift      None             2   \n",
       "\n",
       "   num_sublineages                                        sublineages  \n",
       "0                1                                      4.1.i1.2.1(1)  \n",
       "1                4  4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...  \n",
       "2                3               3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)  \n",
       "3                4       4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)  \n",
       "4                1                                           2.2.2(2)  "
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulator_INDELs_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(8, 14)"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulator_INDELs_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop **INDEL_type** column and rename **ins_del** column to fit with SNPs DataFrames"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/rv76/anaconda2/lib/python2.7/site-packages/ipykernel_launcher.py:1: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n",
      "/home/rv76/anaconda2/lib/python2.7/site-packages/pandas/core/frame.py:2746: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  **kwargs)\n"
     ]
    }
   ],
   "source": [
    "regulator_INDELs_DF.drop(['INDEL_type'], axis=1, inplace = True)\n",
    "regulator_INDELs_DF.rename(columns={'ins_del':'mut_type'}, inplace = True) #change column names"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### SNPs in regulon"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We're going to include SNPs in a gene that either \n",
    "\n",
    "- cause a **premature STOP codon** (nonsense mutation)\n",
    "- **synonymous mutations at the first codon that abolish the start codon** (i.e. eis, whiB7 and mmpR which start with a **valine** as this is not possible with ahpC, mmpS5 or mmpL5, which start with a **methionine**)\n",
    "- **non-synonymous mutations at the first codon that abolish the start codon** (i.e. that change the start to an AA other than **valine** and **methionine**)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_SNPs_DF = SNP_annotation_DF[SNP_annotation_DF.gene_id == regulon_id]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>array_index</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>49</th>\n",
       "      <td>2726194_C</td>\n",
       "      <td>2726194</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>2</td>\n",
       "      <td>N</td>\n",
       "      <td>M1T</td>\n",
       "      <td>39</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50</th>\n",
       "      <td>2726195_A</td>\n",
       "      <td>2726195</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>3</td>\n",
       "      <td>N</td>\n",
       "      <td>M1I</td>\n",
       "      <td>40</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.5(1)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          key      pos ref alt          gene_category gene_name gene_id  \\\n",
       "49  2726194_C  2726194   T   C  Antibiotic Resistance      ahpC  Rv2428   \n",
       "50  2726195_A  2726195   G   A  Antibiotic Resistance      ahpC  Rv2428   \n",
       "\n",
       "   gene_pos mut_type codon_pos array_index  num_isolates  num_sublineages  \\\n",
       "49        2        N       M1T          39             1                1   \n",
       "50        3        N       M1I          40             1                1   \n",
       "\n",
       "    sublineages  \n",
       "49  4.3.i3.1(1)  \n",
       "50       4.5(1)  "
      ]
     },
     "execution_count": 75,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulon_SNPs_DF[[(AA_change[-1]=='*') | ((len(AA_change)==3) and (AA_change[1]=='1')) for AA_change in regulon_SNPs_DF.codon_pos]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_SNPs_DF = regulon_SNPs_DF[[(AA_change[-1]=='*') | ((len(AA_change)==3) and (AA_change[1]=='1')) for AA_change in regulon_SNPs_DF.codon_pos]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(2, 14)"
      ]
     },
     "execution_count": 77,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulon_SNPs_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Save SNPs that occur in the first codon seperately"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_SNPs_DF[[((len(AA_change)==3) and (AA_change[1]=='1')) for AA_change in regulon_SNPs_DF.codon_pos]].to_excel('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/CSV files/START_codon_SNVs/ahpC_START_codon_SNVs.xlsx')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "non-synonymous SNV in *ahpC* that changes the *methionine i* (ATG) to a *isoleucine i* (ATA) does not abolish the START codon so we will exclude this SNV from analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_SNPs_DF.drop(labels = [50], axis = 0, inplace = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1, 14)"
      ]
     },
     "execution_count": 79,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulon_SNPs_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### INDELs in regulon"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### keep only *frameshift* indels for gene body"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>TC_2726221_T</td>\n",
       "      <td>2726221</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>29</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>10</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>AGATCG_2726383_A</td>\n",
       "      <td>2726383</td>\n",
       "      <td>AGATCG</td>\n",
       "      <td>A</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>191</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>64</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>TG_2726440_T</td>\n",
       "      <td>2726440</td>\n",
       "      <td>TG</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>248</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>83</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1.2.2.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>TA_2726509_T</td>\n",
       "      <td>2726509</td>\n",
       "      <td>TA</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>317</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>106</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>TC_2726515_T</td>\n",
       "      <td>2726515</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>323</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>108</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 key      pos     ref alt          gene_category gene_name  \\\n",
       "8       TC_2726221_T  2726221      TC   T  Antibiotic Resistance      ahpC   \n",
       "9   AGATCG_2726383_A  2726383  AGATCG   A  Antibiotic Resistance      ahpC   \n",
       "10      TG_2726440_T  2726440      TG   T  Antibiotic Resistance      ahpC   \n",
       "11      TA_2726509_T  2726509      TA   T  Antibiotic Resistance      ahpC   \n",
       "12      TC_2726515_T  2726515      TC   T  Antibiotic Resistance      ahpC   \n",
       "\n",
       "   gene_id gene_pos ins_del  INDEL_type codon_pos  num_isolates  \\\n",
       "8   Rv2428       29     del  frameshift        10             1   \n",
       "9   Rv2428      191     del  frameshift        64             1   \n",
       "10  Rv2428      248     del  frameshift        83             2   \n",
       "11  Rv2428      317     del  frameshift       106             1   \n",
       "12  Rv2428      323     del  frameshift       108             1   \n",
       "\n",
       "    num_sublineages    sublineages  \n",
       "8                 1    4.3.i3.1(1)  \n",
       "9                 1     4.10.i1(1)  \n",
       "10                1     1.2.2.1(2)  \n",
       "11                1  4.1.i1.2.1(1)  \n",
       "12                1       3.1.1(1)  "
      ]
     },
     "execution_count": 80,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_annotation_DF[(INDEL_annotation_DF.gene_id == regulon_id) & (INDEL_annotation_DF.INDEL_type == 'frameshift')].head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_INDELs_DF = INDEL_annotation_DF[(INDEL_annotation_DF.gene_id == regulon_id) & (INDEL_annotation_DF.INDEL_type == 'frameshift')]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>key</th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>ins_del</th>\n",
       "      <th>INDEL_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>TC_2726221_T</td>\n",
       "      <td>2726221</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>29</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>10</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>AGATCG_2726383_A</td>\n",
       "      <td>2726383</td>\n",
       "      <td>AGATCG</td>\n",
       "      <td>A</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>191</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>64</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>TG_2726440_T</td>\n",
       "      <td>2726440</td>\n",
       "      <td>TG</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>248</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>83</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1.2.2.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>TA_2726509_T</td>\n",
       "      <td>2726509</td>\n",
       "      <td>TA</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>317</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>106</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>TC_2726515_T</td>\n",
       "      <td>2726515</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>323</td>\n",
       "      <td>del</td>\n",
       "      <td>frameshift</td>\n",
       "      <td>108</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 key      pos     ref alt          gene_category gene_name  \\\n",
       "8       TC_2726221_T  2726221      TC   T  Antibiotic Resistance      ahpC   \n",
       "9   AGATCG_2726383_A  2726383  AGATCG   A  Antibiotic Resistance      ahpC   \n",
       "10      TG_2726440_T  2726440      TG   T  Antibiotic Resistance      ahpC   \n",
       "11      TA_2726509_T  2726509      TA   T  Antibiotic Resistance      ahpC   \n",
       "12      TC_2726515_T  2726515      TC   T  Antibiotic Resistance      ahpC   \n",
       "\n",
       "   gene_id gene_pos ins_del  INDEL_type codon_pos  num_isolates  \\\n",
       "8   Rv2428       29     del  frameshift        10             1   \n",
       "9   Rv2428      191     del  frameshift        64             1   \n",
       "10  Rv2428      248     del  frameshift        83             2   \n",
       "11  Rv2428      317     del  frameshift       106             1   \n",
       "12  Rv2428      323     del  frameshift       108             1   \n",
       "\n",
       "    num_sublineages    sublineages  \n",
       "8                 1    4.3.i3.1(1)  \n",
       "9                 1     4.10.i1(1)  \n",
       "10                1     1.2.2.1(2)  \n",
       "11                1  4.1.i1.2.1(1)  \n",
       "12                1       3.1.1(1)  "
      ]
     },
     "execution_count": 82,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulon_INDELs_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 83,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(6, 14)"
      ]
     },
     "execution_count": 83,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulon_INDELs_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop **INDEL_type** column (since all indels are frameshifts) and rename **ins_del** column to fit with SNPs DataFrames"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 84,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/rv76/anaconda2/lib/python2.7/site-packages/ipykernel_launcher.py:1: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n"
     ]
    }
   ],
   "source": [
    "regulon_INDELs_DF.drop(['INDEL_type'], axis=1, inplace = True)\n",
    "regulon_INDELs_DF.rename(columns={'ins_del':'mut_type'}, inplace = True) #change column names"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Save DataFrame for *regulator* & *regulon* mutations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 85,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_SNPs_DF_for_export = regulator_SNPs_DF.drop(['array_index'], axis=1)\n",
    "regulon_SNPs_DF_for_export = regulon_SNPs_DF.drop(['array_index'], axis=1)\n",
    "regulator_regulon_mutations_df = regulator_SNPs_DF_for_export.append(regulon_SNPs_DF_for_export.append(regulator_INDELs_DF.append(regulon_INDELs_DF)))\n",
    "\n",
    "regulator_regulon_mutations_df.drop(['key'], axis=1, inplace = True) #drop col\n",
    "regulator_regulon_mutations_df.sort_values(by = 'pos', inplace = True) #sort by reference position\n",
    "regulator_regulon_mutations_df.reset_index(inplace = True, drop = True) #reset index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 86,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(64, 12)"
      ]
     },
     "execution_count": 86,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(regulator_regulon_mutations_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 87,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726088</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.1.1.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726089</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726091</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1.1.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726093</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>6(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2726095</td>\n",
       "      <td>A</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>2726096</td>\n",
       "      <td>A</td>\n",
       "      <td>ATATG</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>2726100</td>\n",
       "      <td>G</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>4.10.i1(2) - 3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>2726100</td>\n",
       "      <td>G</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>5(1) - 4.2.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2726100</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>3.1.1(2) - 1.1.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>2726103</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "      <td>4.10.i1(1) - 4.2.2(1) - 3.1.1(1) - 2.2.1(1) - ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>2726104</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>3.1.1(3) - 1.2.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3350</td>\n",
       "      <td>12</td>\n",
       "      <td>3.1.1(2240) - 4.2.1.1.1.2.2.1(1) - 4.6.1.2(1) ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>GAT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>4</td>\n",
       "      <td>4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>2726106</td>\n",
       "      <td>A</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>2726107</td>\n",
       "      <td>T</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>2726109</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>2726110</td>\n",
       "      <td>A</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>2726112</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>48</td>\n",
       "      <td>21</td>\n",
       "      <td>4.5(3) - 2.2.1.1.1.i3(2) - 3.1.1(2) - 2.2.1.1....</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>2726113</td>\n",
       "      <td>A</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>6</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(6)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>2726114</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>2726114</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>8</td>\n",
       "      <td>1</td>\n",
       "      <td>1.2.2.2(8)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>2726115</td>\n",
       "      <td>CT</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>None</td>\n",
       "      <td>10</td>\n",
       "      <td>3</td>\n",
       "      <td>3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>2726115</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(4)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>2726116</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>7</td>\n",
       "      <td>2</td>\n",
       "      <td>3.1.1.i1(1) - 3.1.1(6)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>2726116</td>\n",
       "      <td>T</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>10</td>\n",
       "      <td>4</td>\n",
       "      <td>3.1.1(6) - 3.1.1.i2(1) - 3.2.2(1) - 3.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>2726117</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>13</td>\n",
       "      <td>3</td>\n",
       "      <td>3.1.1(11) - 3.1.1.i2(1) - 3.2.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>2726117</td>\n",
       "      <td>T</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(4)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27</th>\n",
       "      <td>2726119</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>11</td>\n",
       "      <td>7</td>\n",
       "      <td>4.3.i3.1(1) - 4.1.i1.1.1.1(3) - 4.3.i4.1(2) - ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>2726119</td>\n",
       "      <td>G</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.2(3)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29</th>\n",
       "      <td>2726121</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>32</td>\n",
       "      <td>15</td>\n",
       "      <td>4.5(3) - 4.3.i3.2(1) - 4.3.i3.1(1) - 4.1.i1.1....</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34</th>\n",
       "      <td>2726139</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>35</th>\n",
       "      <td>2726139</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>28</td>\n",
       "      <td>13</td>\n",
       "      <td>4.5(3) - 4.11(2) - 4.1.i1.1.1.1(3) - 2.2.1.1.1...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>36</th>\n",
       "      <td>2726141</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>92</td>\n",
       "      <td>25</td>\n",
       "      <td>4.11(1) - 2.2.1.1.1.i3(2) - 2.2.1.1.1.i1(1) - ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>37</th>\n",
       "      <td>2726141</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>91</td>\n",
       "      <td>10</td>\n",
       "      <td>4.11(74) - 2.2.1.1.1.i3(1) - 4.5(1) - 4.1.i1.1...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>38</th>\n",
       "      <td>2726142</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>49</td>\n",
       "      <td>16</td>\n",
       "      <td>4.5(2) - 2.2.1.1.1.i3(2) - 4.3.i3.1(4) - 4.1.i...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>39</th>\n",
       "      <td>2726144</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>40</th>\n",
       "      <td>2726145</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>85</td>\n",
       "      <td>20</td>\n",
       "      <td>4.11(1) - 2.2.1.1.1.i3(3) - 4.3.i3.1(3) - 4.5(...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>41</th>\n",
       "      <td>2726146</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>None</td>\n",
       "      <td>37</td>\n",
       "      <td>4</td>\n",
       "      <td>4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>42</th>\n",
       "      <td>2726147</td>\n",
       "      <td>A</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>43</th>\n",
       "      <td>2726148</td>\n",
       "      <td>AT</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.2(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>44</th>\n",
       "      <td>2726149</td>\n",
       "      <td>T</td>\n",
       "      <td>TC</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.5(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>45</th>\n",
       "      <td>2726151</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>12</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(12)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46</th>\n",
       "      <td>2726151</td>\n",
       "      <td>T</td>\n",
       "      <td>G</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>47</th>\n",
       "      <td>2726152</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>48</th>\n",
       "      <td>2726153</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>4.10.i1(1) - 2.2.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>49</th>\n",
       "      <td>2726154</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>4.1.i1.2.1(2) - 4.1.i1.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50</th>\n",
       "      <td>2726157</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2.2.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>51</th>\n",
       "      <td>2726158</td>\n",
       "      <td>CA</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>del</td>\n",
       "      <td>None</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>1.2.2.1(1) - 1.2.2.2(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>52</th>\n",
       "      <td>2726174</td>\n",
       "      <td>T</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53</th>\n",
       "      <td>2726182</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.5(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>54</th>\n",
       "      <td>2726190</td>\n",
       "      <td>G</td>\n",
       "      <td>GT</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>ins</td>\n",
       "      <td>None</td>\n",
       "      <td>7</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i4.2(7)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>55</th>\n",
       "      <td>2726190</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.1.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>56</th>\n",
       "      <td>2726190</td>\n",
       "      <td>G</td>\n",
       "      <td>C</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>57</th>\n",
       "      <td>2726194</td>\n",
       "      <td>T</td>\n",
       "      <td>C</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>2</td>\n",
       "      <td>N</td>\n",
       "      <td>M1T</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>58</th>\n",
       "      <td>2726221</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>29</td>\n",
       "      <td>del</td>\n",
       "      <td>10</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.3.i3.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>59</th>\n",
       "      <td>2726383</td>\n",
       "      <td>AGATCG</td>\n",
       "      <td>A</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>191</td>\n",
       "      <td>del</td>\n",
       "      <td>64</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.10.i1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>60</th>\n",
       "      <td>2726440</td>\n",
       "      <td>TG</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>248</td>\n",
       "      <td>del</td>\n",
       "      <td>83</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1.2.2.1(2)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>61</th>\n",
       "      <td>2726509</td>\n",
       "      <td>TA</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>317</td>\n",
       "      <td>del</td>\n",
       "      <td>106</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.1.i1.2.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>62</th>\n",
       "      <td>2726515</td>\n",
       "      <td>TC</td>\n",
       "      <td>T</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>323</td>\n",
       "      <td>del</td>\n",
       "      <td>108</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>63</th>\n",
       "      <td>2726725</td>\n",
       "      <td>GGCGC</td>\n",
       "      <td>G</td>\n",
       "      <td>Antibiotic Resistance</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>533</td>\n",
       "      <td>del</td>\n",
       "      <td>178</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>4.5(1)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>64 rows × 12 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "        pos     ref    alt          gene_category gene_name          gene_id  \\\n",
       "0   2726088       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "1   2726089       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "2   2726091       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "3   2726093       T      A                   None      None  Rv2427Ac_Rv2428   \n",
       "4   2726095       A      C                   None      None  Rv2427Ac_Rv2428   \n",
       "5   2726096       A  ATATG                   None      None  Rv2427Ac_Rv2428   \n",
       "6   2726100       G      C                   None      None  Rv2427Ac_Rv2428   \n",
       "7   2726100       G      T                   None      None  Rv2427Ac_Rv2428   \n",
       "8   2726100       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "9   2726103       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "10  2726104       T      A                   None      None  Rv2427Ac_Rv2428   \n",
       "11  2726105       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "12  2726105       G    GAT                   None      None  Rv2427Ac_Rv2428   \n",
       "13  2726106       A      G                   None      None  Rv2427Ac_Rv2428   \n",
       "14  2726107       T      G                   None      None  Rv2427Ac_Rv2428   \n",
       "15  2726109       T      C                   None      None  Rv2427Ac_Rv2428   \n",
       "16  2726110       A      G                   None      None  Rv2427Ac_Rv2428   \n",
       "17  2726112       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "18  2726113       A      G                   None      None  Rv2427Ac_Rv2428   \n",
       "19  2726114       C      A                   None      None  Rv2427Ac_Rv2428   \n",
       "20  2726114       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "21  2726115      CT      C                   None      None  Rv2427Ac_Rv2428   \n",
       "22  2726115       C      A                   None      None  Rv2427Ac_Rv2428   \n",
       "23  2726116       T      A                   None      None  Rv2427Ac_Rv2428   \n",
       "24  2726116       T      G                   None      None  Rv2427Ac_Rv2428   \n",
       "25  2726117       T      A                   None      None  Rv2427Ac_Rv2428   \n",
       "26  2726117       T      G                   None      None  Rv2427Ac_Rv2428   \n",
       "27  2726119       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "28  2726119       G      T                   None      None  Rv2427Ac_Rv2428   \n",
       "29  2726121       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "..      ...     ...    ...                    ...       ...              ...   \n",
       "34  2726139       C      A                   None      None  Rv2427Ac_Rv2428   \n",
       "35  2726139       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "36  2726141       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "37  2726141       C      A                   None      None  Rv2427Ac_Rv2428   \n",
       "38  2726142       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "39  2726144       T      C                   None      None  Rv2427Ac_Rv2428   \n",
       "40  2726145       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "41  2726146       G     GT                   None      None  Rv2427Ac_Rv2428   \n",
       "42  2726147       A      G                   None      None  Rv2427Ac_Rv2428   \n",
       "43  2726148      AT      A                   None      None  Rv2427Ac_Rv2428   \n",
       "44  2726149       T     TC                   None      None  Rv2427Ac_Rv2428   \n",
       "45  2726151       T      C                   None      None  Rv2427Ac_Rv2428   \n",
       "46  2726151       T      G                   None      None  Rv2427Ac_Rv2428   \n",
       "47  2726152       C      A                   None      None  Rv2427Ac_Rv2428   \n",
       "48  2726153       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "49  2726154       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "50  2726157       C      T                   None      None  Rv2427Ac_Rv2428   \n",
       "51  2726158      CA      C                   None      None  Rv2427Ac_Rv2428   \n",
       "52  2726174       T      A                   None      None  Rv2427Ac_Rv2428   \n",
       "53  2726182       T      C                   None      None  Rv2427Ac_Rv2428   \n",
       "54  2726190       G     GT                   None      None  Rv2427Ac_Rv2428   \n",
       "55  2726190       G      A                   None      None  Rv2427Ac_Rv2428   \n",
       "56  2726190       G      C                   None      None  Rv2427Ac_Rv2428   \n",
       "57  2726194       T      C  Antibiotic Resistance      ahpC           Rv2428   \n",
       "58  2726221      TC      T  Antibiotic Resistance      ahpC           Rv2428   \n",
       "59  2726383  AGATCG      A  Antibiotic Resistance      ahpC           Rv2428   \n",
       "60  2726440      TG      T  Antibiotic Resistance      ahpC           Rv2428   \n",
       "61  2726509      TA      T  Antibiotic Resistance      ahpC           Rv2428   \n",
       "62  2726515      TC      T  Antibiotic Resistance      ahpC           Rv2428   \n",
       "63  2726725   GGCGC      G  Antibiotic Resistance      ahpC           Rv2428   \n",
       "\n",
       "   gene_pos mut_type codon_pos  num_isolates  num_sublineages  \\\n",
       "0      None        I      None             2                1   \n",
       "1      None        I      None             1                1   \n",
       "2      None        I      None             1                1   \n",
       "3      None        I      None             2                1   \n",
       "4      None        I      None             2                1   \n",
       "5      None      ins      None             1                1   \n",
       "6      None        I      None             3                2   \n",
       "7      None        I      None             2                2   \n",
       "8      None        I      None             3                2   \n",
       "9      None        I      None             5                5   \n",
       "10     None        I      None             4                2   \n",
       "11     None        I      None          3350               12   \n",
       "12     None      ins      None             4                4   \n",
       "13     None        I      None             2                1   \n",
       "14     None        I      None             2                1   \n",
       "15     None        I      None             1                1   \n",
       "16     None        I      None             1                1   \n",
       "17     None        I      None            48               21   \n",
       "18     None        I      None             6                1   \n",
       "19     None        I      None             1                1   \n",
       "20     None        I      None             8                1   \n",
       "21     None      del      None            10                3   \n",
       "22     None        I      None             4                1   \n",
       "23     None        I      None             7                2   \n",
       "24     None        I      None            10                4   \n",
       "25     None        I      None            13                3   \n",
       "26     None        I      None             4                1   \n",
       "27     None        I      None            11                7   \n",
       "28     None        I      None             3                1   \n",
       "29     None        I      None            32               15   \n",
       "..      ...      ...       ...           ...              ...   \n",
       "34     None        I      None             1                1   \n",
       "35     None        I      None            28               13   \n",
       "36     None        I      None            92               25   \n",
       "37     None        I      None            91               10   \n",
       "38     None        I      None            49               16   \n",
       "39     None        I      None             1                1   \n",
       "40     None        I      None            85               20   \n",
       "41     None      ins      None            37                4   \n",
       "42     None        I      None             1                1   \n",
       "43     None      del      None             2                1   \n",
       "44     None      ins      None             1                1   \n",
       "45     None        I      None            12                1   \n",
       "46     None        I      None             1                1   \n",
       "47     None        I      None             1                1   \n",
       "48     None        I      None             2                2   \n",
       "49     None        I      None             3                2   \n",
       "50     None        I      None             1                1   \n",
       "51     None      del      None             2                2   \n",
       "52     None        I      None             1                1   \n",
       "53     None        I      None             1                1   \n",
       "54     None      ins      None             7                1   \n",
       "55     None        I      None             1                1   \n",
       "56     None        I      None             1                1   \n",
       "57        2        N       M1T             1                1   \n",
       "58       29      del        10             1                1   \n",
       "59      191      del        64             1                1   \n",
       "60      248      del        83             2                1   \n",
       "61      317      del       106             1                1   \n",
       "62      323      del       108             1                1   \n",
       "63      533      del       178             1                1   \n",
       "\n",
       "                                          sublineages  \n",
       "0                                        2.2.1.1.1(2)  \n",
       "1                                     4.1.i1.1.1.1(1)  \n",
       "2                                            1.1.2(1)  \n",
       "3                                                6(2)  \n",
       "4                                          4.10.i1(2)  \n",
       "5                                       4.1.i1.2.1(1)  \n",
       "6                               4.10.i1(2) - 3.1.1(1)  \n",
       "7                                   5(1) - 4.2.1.1(1)  \n",
       "8                                 3.1.1(2) - 1.1.2(1)  \n",
       "9   4.10.i1(1) - 4.2.2(1) - 3.1.1(1) - 2.2.1(1) - ...  \n",
       "10                            3.1.1(3) - 1.2.1.1.1(1)  \n",
       "11  3.1.1(2240) - 4.2.1.1.1.2.2.1(1) - 4.6.1.2(1) ...  \n",
       "12  4.10.i1(1) - 4.10.i2(1) - 4.6.1.1.1.1.2(1) - 2...  \n",
       "13                                           3.1.1(2)  \n",
       "14                                           3.1.1(2)  \n",
       "15                                           3.1.1(1)  \n",
       "16                                           3.1.1(1)  \n",
       "17  4.5(3) - 2.2.1.1.1.i3(2) - 3.1.1(2) - 2.2.1.1....  \n",
       "18                                           3.1.1(6)  \n",
       "19                                             3.1(1)  \n",
       "20                                         1.2.2.2(8)  \n",
       "21               3.1.1.i1(7) - 3.1.1(2) - 3.1.1.i2(1)  \n",
       "22                                           3.1.1(4)  \n",
       "23                             3.1.1.i1(1) - 3.1.1(6)  \n",
       "24         3.1.1(6) - 3.1.1.i2(1) - 3.2.2(1) - 3.1(2)  \n",
       "25                 3.1.1(11) - 3.1.1.i2(1) - 3.2.2(1)  \n",
       "26                                           3.1.1(4)  \n",
       "27  4.3.i3.1(1) - 4.1.i1.1.1.1(3) - 4.3.i4.1(2) - ...  \n",
       "28                                           2.2.2(3)  \n",
       "29  4.5(3) - 4.3.i3.2(1) - 4.3.i3.1(1) - 4.1.i1.1....  \n",
       "..                                                ...  \n",
       "34                                       2.2.1.1.1(1)  \n",
       "35  4.5(3) - 4.11(2) - 4.1.i1.1.1.1(3) - 2.2.1.1.1...  \n",
       "36  4.11(1) - 2.2.1.1.1.i3(2) - 2.2.1.1.1.i1(1) - ...  \n",
       "37  4.11(74) - 2.2.1.1.1.i3(1) - 4.5(1) - 4.1.i1.1...  \n",
       "38  4.5(2) - 2.2.1.1.1.i3(2) - 4.3.i3.1(4) - 4.1.i...  \n",
       "39                                        4.3.i3.1(1)  \n",
       "40  4.11(1) - 2.2.1.1.1.i3(3) - 4.3.i3.1(3) - 4.5(...  \n",
       "41       4.5(1) - 4.2.2(33) - 2.2.1.1.1(2) - 4.4.2(1)  \n",
       "42                                         4.10.i2(1)  \n",
       "43                                           2.2.2(2)  \n",
       "44                                             4.5(1)  \n",
       "45                                        4.10.i1(12)  \n",
       "46                                         4.10.i2(1)  \n",
       "47                                        4.3.i3.1(1)  \n",
       "48                          4.10.i1(1) - 2.2.1.1.1(1)  \n",
       "49                    4.1.i1.2.1(2) - 4.1.i1.1.1.1(1)  \n",
       "50                                       2.2.1.1.1(1)  \n",
       "51                            1.2.2.1(1) - 1.2.2.2(1)  \n",
       "52                                           3.1.1(1)  \n",
       "53                                             4.5(1)  \n",
       "54                                        4.3.i4.2(7)  \n",
       "55                                    4.1.i1.1.1.1(1)  \n",
       "56                                           3.1.1(1)  \n",
       "57                                        4.3.i3.1(1)  \n",
       "58                                        4.3.i3.1(1)  \n",
       "59                                         4.10.i1(1)  \n",
       "60                                         1.2.2.1(2)  \n",
       "61                                      4.1.i1.2.1(1)  \n",
       "62                                           3.1.1(1)  \n",
       "63                                             4.5(1)  \n",
       "\n",
       "[64 rows x 12 columns]"
      ]
     },
     "execution_count": 87,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulator_regulon_mutations_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_regulon_mutations_df.to_excel('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/CSV files/mutations_detected_within_isolates/promoter_ahpC_gene_body_mutations_in_31440_isolates.xlsx')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Save DataFrame for *regulator* & *regulon* mutations for *Table 1*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**REGULATOR** Subset to the 4 most commonly occurring mutations in the regulator regions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 89,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_mutations_df_subset = regulator_regulon_mutations_df[regulator_regulon_mutations_df.gene_id == regulator_id].sort_values(by = 'num_isolates', ascending = False).head(n=4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**REGULON** Subset to mutations that occured in at least 2 Mtb sub-lineages in the regulon regions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulon_mutations_df_subset = regulator_regulon_mutations_df[regulator_regulon_mutations_df.gene_id == regulon_id]\n",
    "regulon_mutations_df_subset = regulon_mutations_df_subset[regulon_mutations_df_subset.num_sublineages >= 2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Export as excel file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 91,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_regulon_subset_mutations_df = regulator_mutations_df_subset.append(regulon_mutations_df_subset)\n",
    "regulator_regulon_subset_mutations_df.sort_values(by = 'pos', ascending = True, inplace = True)\n",
    "regulator_regulon_subset_mutations_df.reset_index(inplace = True, drop = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 92,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>pos</th>\n",
       "      <th>ref</th>\n",
       "      <th>alt</th>\n",
       "      <th>gene_category</th>\n",
       "      <th>gene_name</th>\n",
       "      <th>gene_id</th>\n",
       "      <th>gene_pos</th>\n",
       "      <th>mut_type</th>\n",
       "      <th>codon_pos</th>\n",
       "      <th>num_isolates</th>\n",
       "      <th>num_sublineages</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>3350</td>\n",
       "      <td>12</td>\n",
       "      <td>3.1.1(2240) - 4.2.1.1.1.2.2.1(1) - 4.6.1.2(1) ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2726141</td>\n",
       "      <td>C</td>\n",
       "      <td>T</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>92</td>\n",
       "      <td>25</td>\n",
       "      <td>4.11(1) - 2.2.1.1.1.i3(2) - 2.2.1.1.1.i1(1) - ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2726141</td>\n",
       "      <td>C</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>91</td>\n",
       "      <td>10</td>\n",
       "      <td>4.11(74) - 2.2.1.1.1.i3(1) - 4.5(1) - 4.1.i1.1...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2726145</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>85</td>\n",
       "      <td>20</td>\n",
       "      <td>4.11(1) - 2.2.1.1.1.i3(3) - 4.3.i3.1(3) - 4.5(...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       pos ref alt gene_category gene_name          gene_id gene_pos mut_type  \\\n",
       "0  2726105   G   A          None      None  Rv2427Ac_Rv2428     None        I   \n",
       "1  2726141   C   T          None      None  Rv2427Ac_Rv2428     None        I   \n",
       "2  2726141   C   A          None      None  Rv2427Ac_Rv2428     None        I   \n",
       "3  2726145   G   A          None      None  Rv2427Ac_Rv2428     None        I   \n",
       "\n",
       "  codon_pos  num_isolates  num_sublineages  \\\n",
       "0      None          3350               12   \n",
       "1      None            92               25   \n",
       "2      None            91               10   \n",
       "3      None            85               20   \n",
       "\n",
       "                                         sublineages  \n",
       "0  3.1.1(2240) - 4.2.1.1.1.2.2.1(1) - 4.6.1.2(1) ...  \n",
       "1  4.11(1) - 2.2.1.1.1.i3(2) - 2.2.1.1.1.i1(1) - ...  \n",
       "2  4.11(74) - 2.2.1.1.1.i3(1) - 4.5(1) - 4.1.i1.1...  \n",
       "3  4.11(1) - 2.2.1.1.1.i3(3) - 4.3.i3.1(3) - 4.5(...  "
      ]
     },
     "execution_count": 92,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "regulator_regulon_subset_mutations_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "metadata": {},
   "outputs": [],
   "source": [
    "regulator_regulon_subset_mutations_df.to_excel('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/CSV files/mutations_detected_within_isolates/promoter_ahpC_gene_body_mutations_in_31440_isolates_for_Table_1.xlsx')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Construct DataFrame (which displays each pair of co-occuring mutations as a seperate row with info for each pair of mutations)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 94,
   "metadata": {},
   "outputs": [],
   "source": [
    "variant_mutation_pair_df = pd.DataFrame(columns=['A_type', 'A_pos', 'A_ref', 'A_alt', 'A_mut_type', 'A_codon_pos', 'A_gene_id', 'A_gene_name', 'A_gene_pos' ,'A_num_isolates', 'B_type', 'B_pos', 'B_ref', 'B_alt', 'B_mut_type', 'B_codon_pos', 'B_gene_id', 'B_gene_name', 'B_gene_pos', 'B_num_isolates', 'num_isolates_co_occur', 'isolate_IDs', 'sublineages'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get the ordering of the INDELs genotypes columns in terms of the SNPs genotypes columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "metadata": {},
   "outputs": [],
   "source": [
    "SNP_isolate_annotation_DF['isolate_order'] = SNP_isolate_annotation_DF.index\n",
    "SNP_isolate_annotation_DF.set_index('isolate_ID', drop = True, inplace = True)\n",
    "INDEL_isolate_ordering_for_SNP_genotypes_array = SNP_isolate_annotation_DF.loc[INDEL_isolate_annotation_DF.isolate_ID, :].isolate_order.values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  188,  2733, 24409, ..., 31437, 31438, 31439])"
      ]
     },
     "execution_count": 96,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "INDEL_isolate_ordering_for_SNP_genotypes_array"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fill in dataframe by quering the genotypes matrices"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {},
   "outputs": [],
   "source": [
    "#SNPs or INDELs in regulator x SNPs or INDELs in regulon\n",
    "regulator_mut_keys = list(regulator_SNPs_DF.key) + list(regulator_INDELs_DF.key)\n",
    "regulon_mut_keys = list(regulon_SNPs_DF.key) + list(regulon_INDELs_DF.key)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 98,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "#index for DataFrame\n",
    "variant_mutation_pair_df_index = 0\n",
    "\n",
    "for variant_i in regulator_mut_keys: #iterate through regulator SNPs & INDELs\n",
    "    \n",
    "    if variant_i.count('_') == 1: # if variant in regulator is a SNP\n",
    "\n",
    "        SNP_i_pos = int(variant_i.split('_')[0])\n",
    "        SNP_i_alt_allele = variant_i.split('_')[1]\n",
    "\n",
    "        #get the number code that corresponds to the alt allele, 9 indicated bad quality call\n",
    "        base_code_dict = {'A':0, 'C':1, 'G':2, 'T':3}\n",
    "        SNP_i_alt_allele_code = base_code_dict[SNP_i_alt_allele]\n",
    "\n",
    "        #genotypes for each SNP site\n",
    "        genotypes_for_SNP_i = SNP_genotypes_array[SNP_annotation_DF[SNP_annotation_DF.pos == SNP_i_pos].array_index.values[0] , :]\n",
    "\n",
    "        #re-order isolates according to the order of the columns in the INDEL genotypes matrix\n",
    "        genotypes_for_SNP_i = genotypes_for_SNP_i[INDEL_isolate_ordering_for_SNP_genotypes_array]\n",
    "\n",
    "        for variant_j in regulon_mut_keys: #iterate through mutations in regulon\n",
    "\n",
    "            if variant_j.count('_') == 1: #if variant in regulon is a SNP\n",
    "\n",
    "                SNP_j_pos = int(variant_j.split('_')[0])\n",
    "                SNP_j_alt_allele = variant_j.split('_')[1]\n",
    "\n",
    "                #get the number code that corresponds to the alt allele, 9 indicated bad quality call\n",
    "                SNP_j_alt_allele_code = base_code_dict[SNP_j_alt_allele]\n",
    "\n",
    "                #genotypes for each SNP site\n",
    "                genotypes_for_SNP_j = SNP_genotypes_array[SNP_annotation_DF[SNP_annotation_DF.pos == SNP_j_pos].array_index.values[0] , :]\n",
    "\n",
    "                #re-order isolates according to the order of the columns in the INDEL genotypes matrix\n",
    "                genotypes_for_SNP_j = genotypes_for_SNP_j[INDEL_isolate_ordering_for_SNP_genotypes_array]\n",
    "\n",
    "                #boolean filter of isolates with both variants\n",
    "                isolates_with_both_variants_filter = [(SNP_i_called and SNP_j_called) for SNP_i_called , SNP_j_called in zip(list(genotypes_for_SNP_i == SNP_i_alt_allele_code) , list(genotypes_for_SNP_j == SNP_j_alt_allele_code))]\n",
    "\n",
    "                #number of isolates that have both variants\n",
    "                num_isolates_with_both_variants = np.sum(isolates_with_both_variants_filter)\n",
    "\n",
    "                #check to see if the number of isolates carrying both variants was >= 1\n",
    "                if num_isolates_with_both_variants >= 1:\n",
    "\n",
    "                    # get info for mutation in regulator\n",
    "                    regulator_SNP_info = list(regulator_SNPs_DF[(regulator_SNPs_DF.pos == SNP_i_pos) & (regulator_SNPs_DF.alt == SNP_i_alt_allele)].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    # get info for mutation in regulon\n",
    "                    regulon_SNP_info = list(regulon_SNPs_DF[(regulon_SNPs_DF.pos == SNP_j_pos) & (regulon_SNPs_DF.alt == SNP_j_alt_allele)].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    strain_names = [' - '.join(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].isolate_ID))]\n",
    "\n",
    "                    # sublineage calls\n",
    "                    strain_lineage_calls_count_dict = Counter(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].lineage_call))\n",
    "                    strain_lineage_calls_list = []\n",
    "                    for sublineage_i in strain_lineage_calls_count_dict.keys():\n",
    "                        sublineage_i_with_mutations = sublineage_i + '({0})'.format(str(strain_lineage_calls_count_dict[sublineage_i]))\n",
    "                        strain_lineage_calls_list.append(sublineage_i_with_mutations)\n",
    "                    strain_lineage_calls = ' - '.join(strain_lineage_calls_list)\n",
    "\n",
    "                    co_occurring_SNP_i_and_SNP_or_indel_j_info = ['SNP'] + regulator_SNP_info + ['SNP'] + regulon_SNP_info + [num_isolates_with_both_variants] + strain_names + [strain_lineage_calls]\n",
    "                    variant_mutation_pair_df.loc[variant_mutation_pair_df_index, :] = co_occurring_SNP_i_and_SNP_or_indel_j_info\n",
    "                    variant_mutation_pair_df_index += 1\n",
    "\n",
    "            if variant_j.count('_') == 2: #if variant in regulon is an INDEL\n",
    "\n",
    "                #genotypes for each INDEL\n",
    "                genotypes_for_indel_j = INDEL_genotypes_array[INDEL_annotation_DF[INDEL_annotation_DF.key == variant_j].index[0] , :]\n",
    "\n",
    "                #boolean filter of isolates with both variants\n",
    "                isolates_with_both_variants_filter = [(SNP_i_called and indel_j_called) for SNP_i_called , indel_j_called in zip(list(genotypes_for_SNP_i == SNP_i_alt_allele_code) , list(genotypes_for_indel_j == 1))]\n",
    "\n",
    "                #number of isolates that have both variants\n",
    "                num_isolates_with_both_variants = np.sum(isolates_with_both_variants_filter)\n",
    "\n",
    "                #check to see if the number of isolates carrying both variants was >= 1\n",
    "                if num_isolates_with_both_variants >= 1:\n",
    "\n",
    "                    # get info for mutation in regulator\n",
    "                    regulator_SNP_info = list(regulator_SNPs_DF[(regulator_SNPs_DF.pos == SNP_i_pos) & (regulator_SNPs_DF.alt == SNP_i_alt_allele)].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    # get info for mutation in regulon\n",
    "                    regulon_INDEL_info = list(regulon_INDELs_DF[regulon_INDELs_DF.key == variant_j].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    strain_names = [' - '.join(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].isolate_ID))]\n",
    "\n",
    "                    # sublineage calls\n",
    "                    strain_lineage_calls_count_dict = Counter(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].lineage_call))\n",
    "                    strain_lineage_calls_list = []\n",
    "                    for sublineage_i in strain_lineage_calls_count_dict.keys():\n",
    "                        sublineage_i_with_mutations = sublineage_i + '({0})'.format(str(strain_lineage_calls_count_dict[sublineage_i]))\n",
    "                        strain_lineage_calls_list.append(sublineage_i_with_mutations)\n",
    "                    strain_lineage_calls = ' - '.join(strain_lineage_calls_list)\n",
    "\n",
    "                    co_occurring_SNP_i_and_SNP_or_indel_j_info = ['SNP'] + regulator_SNP_info + ['INDEL'] + regulon_INDEL_info + [num_isolates_with_both_variants] + strain_names + [strain_lineage_calls]\n",
    "                    variant_mutation_pair_df.loc[variant_mutation_pair_df_index, :] = co_occurring_SNP_i_and_SNP_or_indel_j_info\n",
    "                    variant_mutation_pair_df_index += 1\n",
    "    \n",
    "    elif variant_i.count('_') == 2: # if variant in regulator is an INDEL\n",
    "        \n",
    "        #genotypes for each INDEL\n",
    "        genotypes_for_indel_i = INDEL_genotypes_array[INDEL_annotation_DF[INDEL_annotation_DF.key == variant_i].index[0] , :]\n",
    "\n",
    "        for variant_j in regulon_mut_keys: #iterate through mutations in regulon\n",
    "\n",
    "            if variant_j.count('_') == 1: #if variant in regulon is a SNP\n",
    "\n",
    "                SNP_j_pos = int(variant_j.split('_')[0])\n",
    "                SNP_j_alt_allele = variant_j.split('_')[1]\n",
    "\n",
    "                #get the number code that corresponds to the alt allele, 9 indicated bad quality call\n",
    "                SNP_j_alt_allele_code = base_code_dict[SNP_j_alt_allele]\n",
    "\n",
    "                #genotypes for each SNP site\n",
    "                genotypes_for_SNP_j = SNP_genotypes_array[SNP_annotation_DF[SNP_annotation_DF.pos == SNP_j_pos].array_index.values[0] , :]\n",
    "\n",
    "                #re-order isolates according to the order of the columns in the INDEL genotypes matrix\n",
    "                genotypes_for_SNP_j = genotypes_for_SNP_j[INDEL_isolate_ordering_for_SNP_genotypes_array]\n",
    "\n",
    "                #boolean filter of isolates with both variants\n",
    "                isolates_with_both_variants_filter = [(indel_i_called and SNP_j_called) for indel_i_called , SNP_j_called in zip(list(genotypes_for_indel_i == 1) , list(genotypes_for_SNP_j == SNP_j_alt_allele_code))]\n",
    "\n",
    "                #number of isolates that have both variants\n",
    "                num_isolates_with_both_variants = np.sum(isolates_with_both_variants_filter)\n",
    "\n",
    "                #check to see if the number of isolates carrying both variants was >= 1\n",
    "                if num_isolates_with_both_variants >= 1:\n",
    "\n",
    "                    # get info for mutation in regulator\n",
    "                    regulator_INDEL_info = list(regulator_INDELs_DF[regulator_INDELs_DF.key == variant_i].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    # get info for mutation in regulon\n",
    "                    regulon_SNP_info = list(regulon_SNPs_DF[(regulon_SNPs_DF.pos == SNP_j_pos) & (regulon_SNPs_DF.alt == SNP_j_alt_allele)].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    strain_names = [' - '.join(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].isolate_ID))]\n",
    "\n",
    "                    # sublineage calls\n",
    "                    strain_lineage_calls_count_dict = Counter(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].lineage_call))\n",
    "                    strain_lineage_calls_list = []\n",
    "                    for sublineage_i in strain_lineage_calls_count_dict.keys():\n",
    "                        sublineage_i_with_mutations = sublineage_i + '({0})'.format(str(strain_lineage_calls_count_dict[sublineage_i]))\n",
    "                        strain_lineage_calls_list.append(sublineage_i_with_mutations)\n",
    "                    strain_lineage_calls = ' - '.join(strain_lineage_calls_list)\n",
    "\n",
    "                    co_occurring_SNP_i_and_SNP_or_indel_j_info = ['INDEL'] + regulator_INDEL_info + ['SNP'] + regulon_SNP_info + [num_isolates_with_both_variants] + strain_names + [strain_lineage_calls]\n",
    "                    variant_mutation_pair_df.loc[variant_mutation_pair_df_index, :] = co_occurring_SNP_i_and_SNP_or_indel_j_info\n",
    "                    variant_mutation_pair_df_index += 1\n",
    "\n",
    "            if variant_j.count('_') == 2: # if variant in regulon is an INDEL\n",
    "\n",
    "                #genotypes for each INDEL\n",
    "                genotypes_for_indel_j = INDEL_genotypes_array[INDEL_annotation_DF[INDEL_annotation_DF.key == variant_j].index[0] , :]\n",
    "\n",
    "                #boolean filter of isolates with both variants\n",
    "                isolates_with_both_variants_filter = [(indel_i_called and indel_j_called) for indel_i_called , indel_j_called in zip(list(genotypes_for_indel_i == 1) , list(genotypes_for_indel_j == 1))]\n",
    "\n",
    "                #number of isolates that have both variants\n",
    "                num_isolates_with_both_variants = np.sum(isolates_with_both_variants_filter)\n",
    "\n",
    "                #check to see if the number of isolates carrying both variants was >= 1\n",
    "                if num_isolates_with_both_variants >= 1:\n",
    "\n",
    "                    # get info for mutation in regulator\n",
    "                    regulator_INDEL_info = list(regulator_INDELs_DF[regulator_INDELs_DF.key == variant_i].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    # get info for mutation in regulon\n",
    "                    regulon_INDEL_info = list(regulon_INDELs_DF[regulon_INDELs_DF.key == variant_j].loc[:, ['pos','ref','alt','mut_type','codon_pos','gene_id','gene_name','gene_pos','num_isolates']].values[0])\n",
    "\n",
    "                    strain_names = [' - '.join(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].isolate_ID))]\n",
    "\n",
    "                    # sublineage calls\n",
    "                    strain_lineage_calls_count_dict = Counter(list(INDEL_isolate_annotation_DF[isolates_with_both_variants_filter].lineage_call))\n",
    "                    strain_lineage_calls_list = []\n",
    "                    for sublineage_i in strain_lineage_calls_count_dict.keys():\n",
    "                        sublineage_i_with_mutations = sublineage_i + '({0})'.format(str(strain_lineage_calls_count_dict[sublineage_i]))\n",
    "                        strain_lineage_calls_list.append(sublineage_i_with_mutations)\n",
    "                    strain_lineage_calls = ' - '.join(strain_lineage_calls_list)\n",
    "\n",
    "                    co_occurring_SNP_i_and_SNP_or_indel_j_info = ['INDEL'] + regulator_INDEL_info + ['INDEL'] + regulon_INDEL_info + [num_isolates_with_both_variants] + strain_names + [strain_lineage_calls]\n",
    "                    variant_mutation_pair_df.loc[variant_mutation_pair_df_index, :] = co_occurring_SNP_i_and_SNP_or_indel_j_info\n",
    "                    variant_mutation_pair_df_index += 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>A_type</th>\n",
       "      <th>A_pos</th>\n",
       "      <th>A_ref</th>\n",
       "      <th>A_alt</th>\n",
       "      <th>A_mut_type</th>\n",
       "      <th>A_codon_pos</th>\n",
       "      <th>A_gene_id</th>\n",
       "      <th>A_gene_name</th>\n",
       "      <th>A_gene_pos</th>\n",
       "      <th>A_num_isolates</th>\n",
       "      <th>...</th>\n",
       "      <th>B_alt</th>\n",
       "      <th>B_mut_type</th>\n",
       "      <th>B_codon_pos</th>\n",
       "      <th>B_gene_id</th>\n",
       "      <th>B_gene_name</th>\n",
       "      <th>B_gene_pos</th>\n",
       "      <th>B_num_isolates</th>\n",
       "      <th>num_isolates_co_occur</th>\n",
       "      <th>isolate_IDs</th>\n",
       "      <th>sublineages</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>SNP</td>\n",
       "      <td>2726105</td>\n",
       "      <td>G</td>\n",
       "      <td>A</td>\n",
       "      <td>I</td>\n",
       "      <td>None</td>\n",
       "      <td>Rv2427Ac_Rv2428</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>3350</td>\n",
       "      <td>...</td>\n",
       "      <td>T</td>\n",
       "      <td>del</td>\n",
       "      <td>108</td>\n",
       "      <td>Rv2428</td>\n",
       "      <td>ahpC</td>\n",
       "      <td>323</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>SAMEA1118691</td>\n",
       "      <td>3.1.1(1)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>1 rows × 23 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "  A_type    A_pos A_ref A_alt A_mut_type A_codon_pos        A_gene_id  \\\n",
       "0    SNP  2726105     G     A          I        None  Rv2427Ac_Rv2428   \n",
       "\n",
       "  A_gene_name A_gene_pos A_num_isolates     ...     B_alt B_mut_type  \\\n",
       "0        None       None           3350     ...         T        del   \n",
       "\n",
       "  B_codon_pos B_gene_id B_gene_name B_gene_pos B_num_isolates  \\\n",
       "0         108    Rv2428        ahpC        323              1   \n",
       "\n",
       "  num_isolates_co_occur   isolate_IDs sublineages  \n",
       "0                     1  SAMEA1118691    3.1.1(1)  \n",
       "\n",
       "[1 rows x 23 columns]"
      ]
     },
     "execution_count": 99,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "variant_mutation_pair_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 100,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1, 23)"
      ]
     },
     "execution_count": 100,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(variant_mutation_pair_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Save dataframe as CSV file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {},
   "outputs": [],
   "source": [
    "variant_mutation_pair_df.to_excel('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/CSV files/co_occurring_mutations_within_isolates/ahpC_promoter_gene_body_co_occurence_mutations_in_31440_isolates.xlsx')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}