In [1]:
workDir = '/home/nick/notebook/SIPSim/dev/Ecoli/'
genomeDir = '/home/nick/notebook/SIPSim/dev/Ecoli/genomes/'
In [79]:
import glob
import nestly
from IPython.display import Image, display, HTML
In [2]:
%load_ext rpy2.ipython
%matplotlib inline
In [3]:
%%R
library(ggplot2)
library(dplyr)
library(tidyr)
library(gridExtra)
In [64]:
# building tree structure
nest = nestly.Nest()
## varying params
nest.add('percIncorp', range(0,101,20))
## set params
nest.add('np', [24], create_dir=False)
nest.add('percTaxa', [100], create_dir=False)
## input/output files
nest.add('fileName', ['ampFrags'], create_dir=False)
nest.add('genome_index', [os.path.join(genomeDir, 'genome_index.txt')], create_dir=False)
nest.add('genome_dir', [genomeDir], create_dir=False)
nest.add('primers', [os.path.join(workDir, '../', '515F-806R.fna')], create_dir=False)
# building directory tree
buildDir = os.path.join(workDir, 'incorp')
nest.build(buildDir)
In [65]:
bashFile = os.path.join(workDir, 'SIPSimRun.sh')
In [66]:
%%writefile $bashFile
#!/bin/bash
# simulating fragments
SIPSim fragments \
{genome_index} \
--fp {genome_dir} \
--fr {primers} \
--fld skewed-normal,9000,2500,-5 \
--flr None,None \
--nf 10000 \
--np {np} \
2> {fileName}.log \
> {fileName}.pkl
# converting to kde object
SIPSim fragment_kde \
{fileName}.pkl \
> {fileName}_kde.pkl
# plotting
SIPSim KDE_plot \
{fileName}_kde.pkl \
--xMin 1.66 --xMax 1.8 \
-o {fileName}_kde.png
SIPSim KDE_plot \
{fileName}_kde.pkl \
--xMin 1.66 --xMax 1.8 --logY 10 \
-o {fileName}_kde_log10.png
# adding diffusion
SIPSim diffusion \
{fileName}_kde.pkl \
--np {np} \
> {fileName}_kde_dif.pkl
# plotting
SIPSim KDE_plot \
{fileName}_kde_dif.pkl \
--xMin 1.66 --xMax 1.8 \
-o {fileName}_kde_dif.png
SIPSim KDE_plot \
{fileName}_kde_dif.pkl \
--xMin 1.66 --xMax 1.8 --logY 10 \
-o {fileName}_kde_dif_log10.png
# making incorp file
SIPSim incorpConfigExample \
--percTaxa {percTaxa} \
--percIncorpUnif {percIncorp} \
> {percTaxa}_{percIncorp}.config
# adding isotope incorporation to BD distribution
SIPSim isoIncorp \
{fileName}_kde_dif.pkl \
{percTaxa}_{percIncorp}.config \
--np {np} \
> {fileName}_kde_dif_incorp.pkl
# plotting
SIPSim KDE_plot \
{fileName}_kde_dif_incorp.pkl \
--xMin 1.66 --xMax 1.8 \
-o {fileName}_kde_dif_incorp.png
SIPSim KDE_plot \
{fileName}_kde_dif_incorp.pkl \
--xMin 1.66 --xMax 1.8 --logY 10 \
-o {fileName}_kde_dif_incorp_log10.png
# BD shift
SIPSim BD_shift \
{fileName}_kde_dif.pkl \
{fileName}_kde_dif_incorp.pkl \
> {fileName}_BD-shift.txt
In [67]:
!chmod 775 $bashFile
In [68]:
!cd $workDir; \
nestrun -j 1 --template-file $bashFile -d incorp --log-file log.txt
In [112]:
tmpDir = os.path.join(workDir, 'incorp')
incorp_L1 = glob.glob(os.path.join(tmpDir, '*/*incorp_1.png'))
incorp_L2 = glob.glob(os.path.join(tmpDir, '*/*incorp_2.png'))
for L1,L2 in zip(sorted(incorp_L1), sorted(incorp_L2)):
print L1
img1 = Image(L1, retina=True)
img2 = Image(L2, retina=True)
display(img1, img2)
In [113]:
# log10
tmpDir = os.path.join(workDir, 'incorp')
incorp_L1 = glob.glob(os.path.join(tmpDir, '*/*incorp_log10_1.png'))
incorp_L2 = glob.glob(os.path.join(tmpDir, '*/*incorp_log10_2.png'))
for L1,L2 in zip(sorted(incorp_L1), sorted(incorp_L2)):
print L1
img1 = Image(L1, retina=True)
img2 = Image(L2, retina=True)
display(img1, img2)
In [102]:
tmpDir = os.path.join(workDir, 'incorp')
incorp_L1 = glob.glob(os.path.join(tmpDir, '*/*incorp_1.png'))
incorp_L2 = glob.glob(os.path.join(tmpDir, '*/*incorp_2.png'))
s = """<table>
<tr>
<th><img src="{}"/></th>
<th><img src="{}"/></th>
</tr></table>"""
for L1,L2 in zip(sorted(incorp_L1), sorted(incorp_L2)):
s = s.format(str(L1), str(L2))
print(s)
t = HTML(s)
display(t)
In [ ]: