Mercurial > repos > rlegendre > ribo_tools

comparison metagene_readthrough.py @ 15:702e60e819c2 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author rlegendre
date Mon, 11 May 2015 09:53:08 -0400
parents 344bacf6acb8
children c87c40e642af
comparison
equal deleted inserted replaced
14:344bacf6acb8 15:702e60e819c2
36 def compute_rpkm(length,count_gene,count_tot) : 36 def compute_rpkm(length,count_gene,count_tot) :
37 37
38 try : 38 try :
39 rpkm = "{0:.4f}".format(count_gene*1000000.0/(count_tot*length)) 39 rpkm = "{0:.4f}".format(count_gene*1000000.0/(count_tot*length))
40 except ArithmeticError : 40 except ArithmeticError :
41 stop_err( 'Illegal division by zero') 41 stop_err( 'Illegal division by zero during computing RPKM')
42 return float(rpkm) 42 return float(rpkm)
43 43
44 44
45 def find_stop(seq) : 45 def find_stop(seq) :
46 ''' 46 '''
680 680
681 with open(os.path.join(subfolder,'readthrough_result.csv'), 'rb') as csvfile: 681 with open(os.path.join(subfolder,'readthrough_result.csv'), 'rb') as csvfile:
682 spamreader = csv.reader(csvfile, delimiter='\t') 682 spamreader = csv.reader(csvfile, delimiter='\t')
683 ## skip the header 683 ## skip the header
684 next(spamreader, None) 684 next(spamreader, None)
685 for row in spamreader: 685 ##test if file is empty or not
686 if row[2] == '-' : 686 if next(spamreader, None):
687 gene_table += '<tr><td>%s</td><td><a href="%s.png" data-lightbox="%s"><img src="%s_thumbnail.png" /></a></td><td><a title="%s">%s</a></td><td>%s</td><td>%s:%s-%s</td><td>%s</td><td>%s</td><td>%s</td><td>%s</td></tr>' %(row[0], row[0], row[0], row[0], row[11], row[1], row[3], row[4], row[5], row[6], row[8], row[9], row[10], row[12]) 687 for row in spamreader:
688 688 if row[2] == '-' :
689 gene_table += '</tbody></table>' 689 gene_table += '<tr><td>%s</td><td><a href="%s.png" data-lightbox="%s"><img src="%s_thumbnail.png" /></a></td><td><a title="%s">%s</a></td><td>%s</td><td>%s:%s-%s</td><td>%s</td><td>%s</td><td>%s</td><td>%s</td></tr>' %(row[0], row[0], row[0], row[0], row[11], row[1], row[3], row[4], row[5], row[6], row[8], row[9], row[10], row[12])
690 690
691 gene_table += '</tbody></table>'
692 else :
693 gene_table = 'Sorry, there are no stop codon readthrough genes in your data\n'
691 694
692 695
693 html_str = """ 696 html_str = """
694 <!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" 697 <!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
695 "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd"> 698 "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
928 help="Bam Ribo-Seq alignments ", metavar="FILE") 931 help="Bam Ribo-Seq alignments ", metavar="FILE")
929 932
930 parser.add_option("-d", "--dirout", dest="dirout", type= "string", 933 parser.add_option("-d", "--dirout", dest="dirout", type= "string",
931 help="write report in this html file and in associated directory", metavar="STR,STR") 934 help="write report in this html file and in associated directory", metavar="STR,STR")
932 935
933 parser.add_option("-e", "--extend", dest="extend", type= "int", 936 parser.add_option("-e", "--extend", dest="extend", type= "int",default = 300 ,
934 help="Lenght of extension after stop in number of base pairs (depends on your organisme)", metavar="INT") 937 help="Lenght of extension after stop in number of base pairs (depends on your organisme)", metavar="INT")
935 938
936 parser.add_option("-q", "--quiet", 939 parser.add_option("-q", "--quiet",
937 action="store_false", dest="verbose", default=True, 940 action="store_false", dest="verbose", default=True,
938 help="don't print status messages to stdout") 941 help="don't print status messages to stdout")