In [2]:

    

baseDir = '/home/nick/notebook/SIPSim/dev/bac_genome1210/'
genomeDir = os.path.join(baseDir, 'genomes')
rnammerDir = os.path.join(baseDir, 'rnammer')


    

In [3]:

    

import glob
from IPython.display import Image


    

In [4]:

    

%load_ext rpy2.ipython


    

In [5]:

    

%%R
library(ggplot2)
library(dplyr)
library(tidyr)
library(gridExtra)


    

    






Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

Loading required package: grid

In [6]:

    

if not os.path.isdir(rnammerDir):
    os.makedirs(rnammerDir)


    

In [13]:

    

%%bash -s "$genomeDir" "$rnammerDir"

find $1 -name "*fasta" | \
    perl -pe 's/.+\/|\.fasta//g' | \
    xargs -n 1 -I % -P 30 bash -c \
    "rnammer -S bac -m ssu -gff $2/%_rrn.gff -f $2/%_rrn.fna -xml $2/%_rrn.xml < $1/%.fasta"


    

In [14]:

    

## Summarizing the results

!cd $rnammerDir; \
    egrep -v "^#" *.gff | \
    grep "16s_rRNA" | \
    perl -pe 's/:/\t/' > ssu_summary.txt


    

In [28]:

    

%%bash -s "$rnammerDir"
cd $1

printf "ssu gene length distribution:\n"
cut -d $'\t' -f 7 ssu_summary.txt | NY_misc_perl stats_descriptive


    

    




ssu gene length distribution:
1	min	589.20
1	Q1	1873.20
1	mean	1895.42
1	median	1938.90
1	Q3	1975.30
1	max	2090.20
1	stdev	146.53

In [29]:

    

! cd $rnammerDir; \
    cat *_rrn.fna > bac_genome1210_16S.fna
    
!printf "Number of sequences: "
! cd $rnammerDir; \
    grep -c ">" bac_genome1210_16S.fna


    

    




Number of sequences: 4557

In [ ]: