|
9
|
1 #purpose: python wrapper to run spp
|
|
|
2 #author: Ziru Zhou
|
|
|
3 #Date: November 2012
|
|
|
4
|
|
|
5 import sys, subprocess, tempfile, shutil, glob, os, os.path, gzip
|
|
|
6 from galaxy import eggs
|
|
|
7 import pkg_resources
|
|
|
8 pkg_resources.require( "simplejson" )
|
|
|
9 import simplejson
|
|
|
10
|
|
|
11 CHUNK_SIZE = 1024
|
|
|
12
|
|
|
13 def main():
|
|
|
14 options = simplejson.load( open( sys.argv[1] ) )
|
|
|
15 output_narrow_peak = sys.argv[2]
|
|
|
16 output_region_peak = sys.argv[3]
|
|
|
17 output_peakshift_file = sys.argv[4]
|
|
|
18 output_rdata_file = sys.argv[5]
|
|
|
19 output_plot_file = sys.argv[6]
|
|
|
20 output_default_file = sys.argv[7]
|
|
|
21 script_path = sys.argv[8]
|
|
|
22
|
|
|
23 cmd = "Rscript %s/deploy_spp.R" % script_path
|
|
|
24 os.system(cmd)
|
|
|
25
|
|
|
26 #set file extensions and set mandatory options
|
|
|
27 #========================================================================================
|
|
|
28 experiment_name = '_'.join( options['experiment_name'].split() ) #save experiment name
|
|
|
29
|
|
|
30 chip_file = "%s.bam" % (options['chip_file'])
|
|
|
31 subprocess.call(["cp", options['chip_file'], chip_file])
|
|
|
32
|
|
|
33 cmdline = "Rscript %s/run_spp.R -c=%s" % (script_path, chip_file )
|
|
|
34 if 'input_file' in options:
|
|
|
35 input_file = "%s.bam" % (options['input_file'])
|
|
|
36 subprocess.call(["cp", options['input_file'], input_file])
|
|
|
37 cmdline = "%s -i=%s" % ( cmdline, input_file )
|
|
|
38
|
|
|
39 #test = "%s_VS_%s.narrowPeak.gz" %(chip_name, input_name)
|
|
|
40 #print test
|
|
|
41 #set additional options
|
|
|
42 #========================================================================================
|
|
|
43 if (options['action'] == "cross_correlation"):
|
|
|
44 cmdline = "%s %s %s %s > default_output.txt" % ( cmdline, options['savp'], options['out'], options['rf'] )
|
|
|
45 elif (options['action'] == "peak_calling"):
|
|
|
46 cmdline = "%s -fdr=%s -npeak=%s %s %s %s %s %s > default_output.txt" % ( cmdline, options['fdr'], options['npeak'], options['savr'], options['savd'], options['savn'], options['savp'], options['rf'] )
|
|
|
47 elif (options['action'] == "idr"):
|
|
|
48 cmdline = "%s -npeak=%s %s %s %s %s > default_output.txt" % ( cmdline, options['npeak'], options['savr'], options['savp'], options['out'], options['rf'] )
|
|
|
49 elif (options['action'] == "custom"):
|
|
|
50 cmdline = "%s -s=%s %s -x=%s -fdr=%s -npeak=%s %s %s" % ( cmdline, options['s'], options['speak'], options['x'], options['fdr'], options['npeak'], options['filtchr'], options['rf'] )
|
|
|
51 cmdline = "%s %s %s %s %s %s > default_output.txt" % ( cmdline, options['out'], options['savn'], options['savr'], options['savp'], options['savd'] )
|
|
|
52
|
|
|
53 #run cmdline
|
|
|
54 #========================================================================================
|
|
|
55 tmp_dir = tempfile.mkdtemp()
|
|
|
56 stderr_name = tempfile.NamedTemporaryFile().name
|
|
|
57 proc = subprocess.Popen( args=cmdline, shell=True, cwd=tmp_dir, stderr=open( stderr_name, 'wb' ) )
|
|
|
58 proc.wait()
|
|
|
59
|
|
|
60 #Do not terminate if error code, allow dataset (e.g. log) creation and cleanup
|
|
|
61 #========================================================================================
|
|
|
62 if proc.returncode:
|
|
|
63 stderr_f = open( stderr_name )
|
|
|
64 while True:
|
|
|
65 chunk = stderr_f.read( CHUNK_SIZE )
|
|
|
66 if not chunk:
|
|
|
67 stderr_f.close()
|
|
|
68 break
|
|
|
69 sys.stderr.write( chunk )
|
|
|
70
|
|
|
71
|
|
|
72 #determine if the outputs are there, copy them to the appropriate dir and filename
|
|
|
73 #========================================================================================
|
|
|
74 created_default_file = os.path.join( tmp_dir, "default_output.txt" )
|
|
|
75 if os.path.exists( created_default_file ):
|
|
|
76 shutil.move( created_default_file, output_default_file )
|
|
|
77
|
|
|
78 chip_name = os.path.basename(options['chip_file'])
|
|
|
79 input_name = os.path.basename(options['input_file'])
|
|
|
80 created_narrow_peak = os.path.join( "/mnt/galaxyData/files/000", "%s_VS_%s.narrowPeak.gz" % (chip_name, input_name) )
|
|
|
81 if os.path.exists( created_narrow_peak ):
|
|
|
82 shutil.move( created_narrow_peak, output_narrow_peak )
|
|
|
83
|
|
|
84 created_region_peak = os.path.join( "/mnt/galaxyData/files/000", "%s_VS_%s.regionPeak.gz" % (chip_name, input_name) )
|
|
|
85 if os.path.exists( created_region_peak ):
|
|
|
86 shutil.move( created_region_peak, output_region_peak )
|
|
|
87
|
|
|
88 created_peakshift_file = os.path.join( tmp_dir, "peakshift.txt" )
|
|
|
89 if os.path.exists( created_peakshift_file ):
|
|
|
90 shutil.move( created_peakshift_file, output_peakshift_file )
|
|
|
91
|
|
|
92 created_rdata_file = os.path.join( tmp_dir, "%s.Rdata" % options['chip_file'] )
|
|
|
93 if os.path.exists( created_rdata_file ):
|
|
|
94 shutil.move( created_rdata_file, output_rdata_file )
|
|
|
95
|
|
|
96 created_plot_file = os.path.join( tmp_dir, "%s.pdf" % options['chip_file'] )
|
|
|
97 if os.path.exists( created_plot_file ):
|
|
|
98 shutil.move( created_plot_file, output_plot_file )
|
|
|
99
|
|
|
100
|
|
|
101 os.unlink( stderr_name )
|
|
|
102 #os.rmdir( tmp_dir )
|
|
|
103
|
|
|
104 if __name__ == "__main__": main()
|
