6
|
1 #!/usr/bin/env python
|
10
|
2 #Gregoire Seguin-Henry (Engineer IT)
|
|
3 #Amine Sbitti (Data Scientist)
|
|
4 #Ludovic Marie-Sainte (Project Manager)
|
|
5 #For Geviteam 2014
|
6
|
6
|
|
7 """
|
|
8 A wrapper script for running the GenomeAnalysisTK.jar commands.
|
|
9 """
|
|
10 from __future__ import print_function
|
|
11 import sys, argparse, os, tempfile, subprocess, shutil
|
|
12 from binascii import unhexlify
|
|
13 from string import Template
|
|
14 from galaxy import eggs
|
|
15
|
|
16 def cleanup_before_exit( tmp_dir ):
|
|
17 if tmp_dir and os.path.exists( tmp_dir ):
|
|
18 shutil.rmtree( tmp_dir )
|
|
19
|
|
20 def _create_config(args, config_path):
|
|
21 conf_file = open(config_path, "w")
|
|
22 conf_file.write("[user]\n")
|
|
23 for option in args:
|
7
|
24 if not option in ["tumorBam", "normalBam", "refFile", "configFile", "scriptPath"] and args[option]!=None:
|
6
|
25 conf_file.write("%s=%s\n" % (option, args[option]))
|
|
26 conf_file.close()
|
|
27
|
|
28 def my_Popen(cmd, prefix_for_stderr_name, tmp_dir, msg_error):
|
|
29 stderr_name = tempfile.NamedTemporaryFile( prefix = prefix_for_stderr_name ).name
|
|
30 proc = subprocess.Popen( args=cmd, shell=True, stderr=open( stderr_name, 'wb' ) )
|
|
31 return_code = proc.wait()
|
|
32 if return_code:
|
|
33 for line in open( stderr_name ):
|
|
34 print(line, file=sys.stderr)
|
|
35 os.unlink( stderr_name ) #clean up
|
|
36 cleanup_before_exit( tmp_dir )
|
|
37 raise Exception( msg_error )
|
|
38 else:
|
|
39 os.unlink( stderr_name )
|
|
40
|
|
41 def index_bam_files( bam_filenames, tmp_dir ):
|
|
42 for bam_filename in bam_filenames:
|
|
43 bam_index_filename = "%s.bai" % bam_filename
|
|
44 print("bam_filename is: " + bam_filename + " bam_index_filename is: " + bam_index_filename + " test is: %s" % os.path.exists(bam_index_filename))
|
|
45 if not os.path.exists( bam_index_filename ):
|
|
46 #need to index this bam file
|
|
47 command = 'samtools index %s %s' % ( bam_filename, bam_index_filename )
|
|
48 my_Popen( command, "bam_index_stderr", tmp_dir, "Error during indexation of fasta file :" + bam_filename)
|
|
49
|
|
50 def index_fasta_files( fasta_filenames, tmp_dir ):
|
|
51 for fasta_filename in fasta_filenames:
|
|
52 fasta_index_filename = "%s.fai" % fasta_filename
|
|
53 print("fasta_filename is: " + fasta_filename + " fasta_index_filename is: " + fasta_index_filename + " test is: %s" % os.path.exists(fasta_index_filename))
|
|
54 if not os.path.exists( fasta_index_filename ):
|
|
55 #need to index this bam file
|
|
56 command = 'samtools faidx %s %s' % ( fasta_filename, fasta_index_filename )
|
|
57 my_Popen( command, "fasta_index_stderr", tmp_dir, "Error during indexation of fasta file :" + fasta_filename)
|
|
58
|
|
59 def __main__():
|
10
|
60 #Manage options
|
17
|
61 print(os.environ['PATH'])
|
6
|
62 parser = argparse.ArgumentParser()
|
|
63 parser.add_argument( '-t', '--tumorBam', help='path to tumor bam file', required = False )
|
18
|
64 parser.add_argument( '-n', '--normalBam', help='', required = False )
|
|
65 parser.add_argument( '-r', '--refFile', help='', required = False )
|
|
66 parser.add_argument( '-c', '--configFile', help='', required = False )
|
|
67 parser.add_argument( '--depthFilterMultiple', help='', required = False )
|
|
68 parser.add_argument( '--snvMaxFilteredBasecallFrac', help='', required = False )
|
|
69 parser.add_argument( '--snvMaxSpanningDeletionFrac', help='', required = False )
|
|
70 parser.add_argument( '--indelMaxRefRepeat', help='', required = False )
|
|
71 parser.add_argument( '--indelMaxWindowFilteredBasecallFrac', help='', required = False )
|
|
72 parser.add_argument( '--indelMaxIntHpolLength', help='', required = False )
|
|
73 parser.add_argument( '--ssnvPrior', help='', required = False )
|
|
74 parser.add_argument( '--sindelPrior', help='', required = False )
|
|
75 parser.add_argument( '--ssnvNoise', help='', required = False )
|
|
76 parser.add_argument( '--sindelNoise', help='', required = False )
|
|
77 parser.add_argument( '--ssnvNoiseStrandBiasFrac', help='', required = False )
|
|
78 parser.add_argument( '--minTier1Mapq', help='', required = False )
|
|
79 parser.add_argument( '--minTier2Mapq', help='', required = False )
|
|
80 parser.add_argument( '--ssnvQuality_LowerBound', help='', required = False )
|
|
81 parser.add_argument( '--sindelQuality_LowerBound', help='', required = False )
|
|
82 parser.add_argument( '--isWriteRealignedBam', help='', required = False )
|
6
|
83 parser.add_argument( '--binSize', help='path to tumor bam file', required = False )
|
18
|
84 parser.add_argument( '--extraStrelkaArguments', help='', required = False )
|
|
85 parser.add_argument( '--isSkipDepthFilters', help='', required = False )
|
|
86 parser.add_argument( '--maxInputDepth', help='', required = False )
|
|
87 parser.add_argument( '--scriptPath', help='', required = False )
|
6
|
88 args = parser.parse_args()
|
|
89
|
8
|
90 root_dir= args.scriptPath
|
7
|
91 expected_dir="for_tests"
|
|
92 job_dir=os.getcwd()
|
|
93 analysis_dir=job_dir + "/StrelkaAnalysis"
|
|
94 config_script=root_dir + "/configureStrelkaWorkflow.pl"
|
9
|
95 tmp_dir = tempfile.mkdtemp( prefix='tmp-strelkaAnalysis-' )
|
7
|
96 config_ini = "%s/config.ini" % (tmp_dir)
|
|
97
|
14
|
98 print("root_dir: " + root_dir + "\njob_dir :" + job_dir + "\nanalysis_dir :" + analysis_dir + "\nconfig_script :" + config_script + "\ntmp_dir :" + tmp_dir + "\nconfig_ini :" + config_ini)
|
7
|
99
|
|
100
|
6
|
101 #verifying eveything's ok
|
|
102 if not os.path.isfile(config_script):
|
|
103 sys.exit("ERROR: The strelka workflow must be built prior to running. See installation instructions in '$root_dir/README'")
|
|
104 print("configuring...", file=sys.stdout)
|
|
105 if os.path.exists(analysis_dir):
|
|
106 sys.exit("'" + analysis_dir + "' already exist, if you are executing this tool from galaxy it should not happen")
|
|
107
|
|
108
|
|
109 # creating index if needed
|
|
110 bam_filenames = [ args.tumorBam, args.normalBam ]
|
|
111 index_bam_files( bam_filenames, tmp_dir )
|
|
112 fasta_files = [ args.refFile ]
|
|
113 index_fasta_files( fasta_files, tmp_dir )
|
|
114
|
|
115 #creating config file if needed
|
|
116 if args.configFile == "Custom":
|
|
117 _create_config(vars(args), config_ini)
|
18
|
118 elif args.configFile in ["strelka_config_bwa_default.ini", "strelka_config_isaac_default.ini", "strelka_config_eland_default.ini"]:
|
|
119 cmdbash="cp %s %s" % (root_dir + "/lib/" + args.configFile, config_ini)
|
6
|
120 my_Popen(cmdbash, "copy_default_file_err", tmp_dir, "Error during the copy of default config file, maybe it was removed")
|
|
121 else:
|
|
122 if not os.path.exists(args.configFile):
|
|
123 print( "The path to your configuration File seems to be wrong, use another one or custom option", file=sys.stderr)
|
|
124 cmdbash="cp %s %s" % (args.configFile, config_ini)
|
9
|
125 my_Popen(cmdbash, "copy_default_file_err", tmp_dir, "Error during the copy of the selected config file")
|
6
|
126
|
|
127
|
|
128
|
|
129
|
|
130 #configuration of workflow
|
|
131 cmd="%s --tumor=%s --normal=%s --ref=%s --config=%s --output-dir=%s" % (config_script, args.tumorBam, args.normalBam, args.refFile, config_ini, analysis_dir)
|
|
132 print( "**** Starting configuration.")
|
|
133 print( "**** Configuration cmd: '" + cmd + "'")
|
|
134 my_Popen( cmd, "cinfugation_stderr", tmp_dir, "Error during configuration !")
|
|
135 print("completed configuration")
|
|
136
|
|
137 #run the workflow !
|
|
138 cmd="make -C " + analysis_dir
|
|
139 print("**** starting workflow.")
|
|
140 print("**** workflow cmd: '" + cmd + "'")
|
|
141 my_Popen( cmd, "workflow_stderr", tmp_dir, "Error during workflow execution !")
|
|
142 print("**** completed workflow execution")
|
17
|
143
|
|
144 cmdbash="cp %s %s" % (config_ini, analysis_dir + "/config.ini")
|
|
145 my_Popen(cmdbash, "copy_final_conf_file_err", tmp_dir, "Error during the copy of conf file after job is done, quite strange...")
|
|
146
|
6
|
147
|
14
|
148 if __name__=='__main__':
|
|
149 __main__()
|