view DESeq_diff_ann_wrapper.py @ 0:8c5de60b3c04 draft

Uploaded

#!/usr/bin/env python
# -*- coding: utf-8 -*-

'''
Created on 10 sep. 2012

@author: Vladimir Daric, Rachel Legendre, Coline Billerey, Alban Ott, Claire Wallon
@copyright: eBio [e-mail]
@license: GPL v3
'''

from optparse import OptionParser
from os import path, mkdir, getcwd, walk
from shutil import copyfile
import subprocess
import sys

def stop_err( msg ):
    sys.stderr.write( "%s\n" % msg )
    sys.exit()    
    
def __main__():
    #Parse Command Line
    parser = OptionParser()
    
    # DESeq-count options.
    parser.add_option( "-a", "--alpha", type="float", dest="alpha", 
                       default = 0.05 )
    parser.add_option( "-f", "--foldchange", type="float", dest="foldchange", 
                       default = 0.05 )
    parser.add_option( "-b", "--header", action="store_true", dest="header", 
                       default = False )
    parser.add_option( "-r", "--replicates", action="store_true", dest="replicates", 
                       default = False )
    parser.add_option( "-l", "--logbase", type="int", dest="logbase", 
                       default = 10 )
    parser.add_option( "-n", "--geneNumber", type="int", dest="geneNumber", 
                       default = 30 )
    parser.add_option( "-d", "--subfolder", type="string", dest="subfolder")
    
    parser.add_option( "-c", "--cond1", type="string", dest="cond1")
    
    parser.add_option( "-C", "--cond2", type="string", dest="cond2")
    
    parser.add_option( "-s", "--sharingMode", type="string", dest="mode", default = "fit-only"  )
    parser.add_option( "-t", "--fitType", type="string", dest="fitType", default = "local"  )                       
    parser.add_option( "-m", "--method", type="string", dest="method", default = "pooled" )
    # Wrapper / Galaxy options.
    parser.add_option( '-o', '--out', dest='outfile', help="file containing results for all counts" )
    parser.add_option( '-u', '--up', dest='up_file', help="file containing results for up regulates counts" )
    parser.add_option( '-w', '--down', dest='down_file', help="file containing results for down regulated counts" )
    
    (options, args) = parser.parse_args()
    #make a folder for results
    outfile = options.outfile
    workfolder, filename = path.split(outfile)
    #mksubdir
    newdirname = options.subfolder
    if path.exists(newdirname):
        raise
    try:
        mkdir(newdirname)
    except:
        raise

    #launch script_deseq in subfolder
    exebinpath = path.join(path.dirname(__file__), 'script_deseq.R')

    cmd = "%s --OutputDir=%s/DESeq_out --alpha=%f --foldchange=%f --logBase=%d --geneNumber=%d --sharingMode=%s --method=%s --fitType==%s" % (exebinpath, newdirname, options.alpha, options.foldchange, options.logbase, options.geneNumber,options.mode, options.method, options.fitType)
    if options.replicates:
        cmd += " --IsReplicate=TRUE"
    if options.header:
        cmd += " --IsHeader=TRUE"
    for myarg in args:
                
        dataset, cond, rep = myarg.split(',')
        ## replace group with correct name
        if (cond == "1") :
                cond = options.cond1
        else :
                cond = options.cond2
        cmd += ' --condition=%s --replicate=%s --HTSeqcount=%s' % (cond, rep, dataset)

    try:
        proc = subprocess.Popen( args=cmd, shell=True, stdout=None, stderr=subprocess.PIPE )
        returncode = proc.wait()
        
        stderr = ''
        buffsize = 1048576
        try:
            while True:
                stderr += proc.stderr.read( buffsize )
                if not stderr or len( stderr ) % buffsize != 0:
                    break
        except OverflowError:
            pass
        if returncode != 0:
            raise Exception, stderr
    except Exception, e:
        stop_err( 'Error running DESeq analyse: %s'% e )
        
    ##copy up and down files in galaxy dataset :
    copyfile(path.join(newdirname, 'DESeq_out_up_genes.csv'), options.up_file)
    copyfile(path.join(newdirname, 'DESeq_out_down_genes.csv'), options.down_file)
    copyfile(path.join(path.dirname(__file__), 'deseq_result_html.html'), outfile)
    
if __name__=="__main__": 
    __main__()