(C) Phylogeny Construction for sub-lineage 4.11.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# This notebook was created to prepare the data and construct a phylogeny for the sub-lineage 4.11 (also known as 4.2.1.1.1.1.2) isolates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "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": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import os\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import itertools\n",
    "import gzip\n",
    "import vcf\n",
    "from slurmpy import Slurm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##################################################################################################################################################################################################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [1] Interactive Session to prepare data for IQTree"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##################################################################################################################################################################################################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Different isolates groups: '1', '2', '3', '4A', '4B', '4C', '5', '6'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Inputs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "sub_lineage = '4.2.1.1.1.1.2'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Load in the isolate annotation DF for Genotypes Matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "#load isolate annotation file (columns of Genotype Matrix)\n",
    "isolate_annotation_DF = pd.read_pickle('/n/data1/hms/dbmi/farhat/Roger/homoplasy_project/rolling_DB_scrape/Genotypes_Filtered_2/genotypes_isolate_annotation.pkl')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "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>lineage_1</th>\n",
       "      <th>lineage_2</th>\n",
       "      <th>lineage_3</th>\n",
       "      <th>lineage_4</th>\n",
       "      <th>lineage_5</th>\n",
       "      <th>lineage_6</th>\n",
       "      <th>lineage_7</th>\n",
       "      <th>lineage_8</th>\n",
       "      <th>lineage_9</th>\n",
       "      <th>lineage_10</th>\n",
       "      <th>lineage_11</th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "      <th>group</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>i3</td>\n",
       "      <td>1</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>SAMEA3558733</td>\n",
       "      <td>4.2.1.2.1.1.i3.1</td>\n",
       "      <td>4B</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>SAMN03648641</td>\n",
       "      <td>4.2.1.2.2.1.1</td>\n",
       "      <td>4B</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>i1</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>SAMN03647419</td>\n",
       "      <td>3.1.1.i1</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>i1</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>SAMEA3671418</td>\n",
       "      <td>4.2.1.2.1.1.i1</td>\n",
       "      <td>4B</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>SAMN07659096</td>\n",
       "      <td>1.1.1.2</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  lineage_1 lineage_2 lineage_3 lineage_4 lineage_5 lineage_6 lineage_7  \\\n",
       "0         4         2         1         2         1         1        i3   \n",
       "1         4         2         1         2         2         1         1   \n",
       "2         3         1         1        i1       NaN       NaN       NaN   \n",
       "3         4         2         1         2         1         1        i1   \n",
       "4         1         1         1         2       NaN       NaN       NaN   \n",
       "\n",
       "  lineage_8 lineage_9 lineage_10 lineage_11    isolate_ID      lineage_call  \\\n",
       "0         1       NaN        NaN        NaN  SAMEA3558733  4.2.1.2.1.1.i3.1   \n",
       "1       NaN       NaN        NaN        NaN  SAMN03648641     4.2.1.2.2.1.1   \n",
       "2       NaN       NaN        NaN        NaN  SAMN03647419          3.1.1.i1   \n",
       "3       NaN       NaN        NaN        NaN  SAMEA3671418    4.2.1.2.1.1.i1   \n",
       "4       NaN       NaN        NaN        NaN  SAMN07659096           1.1.1.2   \n",
       "\n",
       "  group  \n",
       "0    4B  \n",
       "1    4B  \n",
       "2     3  \n",
       "3    4B  \n",
       "4     1  "
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isolate_annotation_DF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(31428, 14)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(isolate_annotation_DF)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Create a text file with a list of all of the isolates that belong to a particular lineage"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "isolates_belonging_to_sub_lineage = isolate_annotation_DF[isolate_annotation_DF.lineage_call == sub_lineage]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "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>lineage_1</th>\n",
       "      <th>lineage_2</th>\n",
       "      <th>lineage_3</th>\n",
       "      <th>lineage_4</th>\n",
       "      <th>lineage_5</th>\n",
       "      <th>lineage_6</th>\n",
       "      <th>lineage_7</th>\n",
       "      <th>lineage_8</th>\n",
       "      <th>lineage_9</th>\n",
       "      <th>lineage_10</th>\n",
       "      <th>lineage_11</th>\n",
       "      <th>isolate_ID</th>\n",
       "      <th>lineage_call</th>\n",
       "      <th>group</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>49</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>01R0685</td>\n",
       "      <td>4.2.1.1.1.1.2</td>\n",
       "      <td>4C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>280</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>00R1405</td>\n",
       "      <td>4.2.1.1.1.1.2</td>\n",
       "      <td>4C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>348</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>02R1942</td>\n",
       "      <td>4.2.1.1.1.1.2</td>\n",
       "      <td>4C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>621</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>Peru4582</td>\n",
       "      <td>4.2.1.1.1.1.2</td>\n",
       "      <td>4C</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>770</th>\n",
       "      <td>4</td>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>Peru5115</td>\n",
       "      <td>4.2.1.1.1.1.2</td>\n",
       "      <td>4C</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    lineage_1 lineage_2 lineage_3 lineage_4 lineage_5 lineage_6 lineage_7  \\\n",
       "49          4         2         1         1         1         1         2   \n",
       "280         4         2         1         1         1         1         2   \n",
       "348         4         2         1         1         1         1         2   \n",
       "621         4         2         1         1         1         1         2   \n",
       "770         4         2         1         1         1         1         2   \n",
       "\n",
       "    lineage_8 lineage_9 lineage_10 lineage_11 isolate_ID   lineage_call group  \n",
       "49        NaN       NaN        NaN        NaN    01R0685  4.2.1.1.1.1.2    4C  \n",
       "280       NaN       NaN        NaN        NaN    00R1405  4.2.1.1.1.1.2    4C  \n",
       "348       NaN       NaN        NaN        NaN    02R1942  4.2.1.1.1.1.2    4C  \n",
       "621       NaN       NaN        NaN        NaN   Peru4582  4.2.1.1.1.1.2    4C  \n",
       "770       NaN       NaN        NaN        NaN   Peru5115  4.2.1.1.1.1.2    4C  "
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isolates_belonging_to_sub_lineage.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(188, 14)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(isolates_belonging_to_sub_lineage)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "isolate_IDs_for_sub_lineage = [isolate_ID for isolate_ID in list(isolates_belonging_to_sub_lineage.isolate_ID) ]\n",
    "\n",
    "#save each element as a new row in a text file\n",
    "with open('/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/sublineage_isolate_lists/isolate_list_lineage_' + sub_lineage + '.txt', 'w') as f:\n",
    "    for item in isolate_IDs_for_sub_lineage:\n",
    "        f.write(\"%s\\n\" % item)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Number of isolates that have a lineage call at the sub-lineage level (**global lineage.sub lineage**)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "188"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.shape(isolates_belonging_to_global_lineage)[0] - sum( np.isnan(isolates_belonging_to_global_lineage.lineage_2.values.astype(float)) )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### [interactive session] BASH function to prepare data for iqtree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "prep_data_for_iqtree(){\n",
    "    \n",
    "# Load BCFtools\n",
    "module load bcftools/1.9\n",
    "    \n",
    "## The default ulimit is 1024. I should increase it. I think on o2 the hard limit is\n",
    "ulimit -n 50000  #'ulimit -Hn' to check the hard limit, do this for lineage 4 and find login node with high ulimit capacity\n",
    "    \n",
    "# I initialize the variables\n",
    "local OUTDIR=$1\n",
    "local TAG=$2\n",
    "local LIST_ISOL=$3\n",
    "\n",
    "echo \"OUTDIR: ${1}\"\n",
    "echo \"TAG: ${2}\"\n",
    "echo \"LIST_ISOL: ${3}\"\n",
    "\n",
    "# I create the output directory\n",
    "rm -rf ${OUTDIR}\n",
    "mkdir -p ${OUTDIR}\n",
    "    \n",
    "# Change the current directory to the output directory\n",
    "cd ${OUTDIR}\n",
    "\n",
    "# I clean the directory with the bcf data\n",
    "rm -rf bcf/\n",
    "mkdir -p bcf/\n",
    "\n",
    "# I generate the bcf files, change the headers and index them\n",
    "for i in `cat ${LIST_ISOL}`;do\n",
    "bcftools view /n/data1/hms/dbmi/farhat/rollingDB/genomic_data/${i}/pilon/${i}.vcf --types snps -O b -o ${i}_snps.bcf\n",
    "bcftools index ${i}_snps.bcf\n",
    "bcftools view ${i}_snps.bcf --regions-file /home/lf61/lf61/mic_assemblies/40-full-analysis/lin-sp-var-10k/results/data_cleaning/regions_to_include_final.bed -o ${i}.bcf -O b;\n",
    "bcftools reheader -s <(echo \"${i}\") ${i}.bcf -o bcf/${i}_renamed.bcf;\n",
    "bcftools index bcf/${i}_renamed.bcf;\n",
    "done\n",
    "\n",
    "# I add canettii\n",
    "bcftools view /n/data1/hms/dbmi/farhat/lfreschi/mic_assemblies/40-full-analysis/lin-sp-var-10k/results/vcf_canettii/canettii.vcf --types snps -O b -o canettii_snps.bcf\n",
    "bcftools index canettii_snps.bcf\n",
    "bcftools view canettii_snps.bcf --regions-file /home/lf61/lf61/mic_assemblies/40-full-analysis/lin-sp-var-10k/results/data_cleaning/regions_to_include_final.bed -o canettii.bcf -O b;\n",
    "bcftools reheader -s <(echo \"canettii\") canettii.bcf -o bcf/canettii_renamed.bcf;\n",
    "bcftools index bcf/canettii_renamed.bcf;\n",
    "\n",
    "# cleaning\n",
    "rm -rf *.bcf\n",
    "rm -rf *.bcf.csi\n",
    "\n",
    "# I merge the bcf files\n",
    "ls bcf/*|grep \"bcf$\" > list_bcfs.txt\n",
    "bcftools merge -l list_bcfs.txt -o ${OUTDIR}/${TAG}.vcf -O v -0\n",
    "\n",
    "## I use vcf2phylip (from luca's home directory)\n",
    "/home/lf61/sw/vcf2phylip/1.5/vcf2phylip.py -i ${OUTDIR}/${TAG}.vcf -f\n",
    "\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Get variables for data preperation function & for iqtree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "OUTDIR = '/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/tree_output_files/phylogeny_lineage_' + sub_lineage +  '/'\n",
    "TAG = 'lineage_' + sub_lineage\n",
    "LIST_ISOL = '/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/sublineage_isolate_lists/isolate_list_lineage_' + sub_lineage + '.txt'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/tree_output_files/phylogeny_lineage_4.2.1.1.1.1.2/'"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "OUTDIR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'lineage_4.2.1.1.1.1.2'"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "TAG"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'/n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/sublineage_isolate_lists/isolate_list_lineage_4.2.1.1.1.1.2.txt'"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "LIST_ISOL"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### [interactive session] Call function above in an interactive session to prepare the data for iqtree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "prep_data_for_iqtree /n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/tree_output_files/phylogeny_lineage_4.2.1.1.1.1.2/ lineage_4.2.1.1.1.1.2 /n/data1/hms/dbmi/farhat/Roger/mmpR_BDQ_mutant_project/phylogenies/sublineage_isolate_lists/isolate_list_lineage_4.2.1.1.1.1.2.txt\n"
     ]
    }
   ],
   "source": [
    "print 'prep_data_for_iqtree' + ' ' + OUTDIR + ' ' + TAG + ' ' + LIST_ISOL"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##################################################################################################################################################################################################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [2] Submit jobs to O2 to create trees"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##################################################################################################################################################################################################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Use code below to submit a iqtree job (uses ModelFinder)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iq_L4.2.1.1.1.1.2\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "submitted: Submitted batch job 16619247\n"
     ]
    }
   ],
   "source": [
    "#iqtree_job_p1 = 'VIRTUAL_ENV_DISABLE_PROMPT=true source activate iqtree_virtualenv'\n",
    "iqtree_job_p1 = 'set +eu\\nVIRTUAL_ENV_DISABLE_PROMPT=true source activate iqtree_virtualenv\\nset -eu' # to fix \"/home/rv76/anaconda2/etc/profile.d/conda.sh: line 55: PS1: unbound variable\" ERROR\n",
    "iqtree_job_p2 = 'time iqtree -s {0}{1}.min4.fasta -mem 7G -pre {0}tree_{1}_iqtree_FINAL -m MFP -mset GTR -bb 1000 -alrt 1000 -nt AUTO -ntmax 8 -redo'.format(OUTDIR, TAG)\n",
    "\n",
    "iqtree_job = iqtree_job_p1 + '\\n' + iqtree_job_p2\n",
    "\n",
    "#directory where you want output + error files\n",
    "os.chdir(OUTDIR)\n",
    "\n",
    "job_name = 'iq_L' + TAG.split('_')[1]  \n",
    "\n",
    "s = Slurm(job_name , {'partition':'priority' , 'N':'1' , 'c':'8' , 't':'0-12:00:00' , 'mem-per-cpu':'8G' , 'mail-type':'ALL' , 'mail-user':'roger_vargas@g.harvard.edu' , 'o':'out_run_iqtree_{}.txt'.format(TAG)})\n",
    "\n",
    "#submits the job\n",
    "job_id = s.run(iqtree_job)\n",
    "\n",
    "print job_name"
   ]
  }
 ],
 "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
}