annotate VCFFiltering_wrapper.py @ 1:cfd4eaadad42 draft

Uploaded
author urgi-team
date Tue, 15 Dec 2015 05:36:12 -0500
parents 3552a8d9f51c
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
1 #!/usr/bin/env python
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
2
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
3
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
4 import subprocess, tempfile, sys, os, glob, shutil, time
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
5 from optparse import OptionParser
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
6 from optparse import Option, OptionValueError
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
7
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
8 class VCFFilteringWrapper(object):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
9
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
10 def __init__(self):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
11 self._options = None
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
12
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
13
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
14 def stop_err(self, msg ):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
15 sys.stderr.write( "%s\n" % msg )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
16 sys.exit()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
17
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
18
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
19 def setAttributesFromCmdLine(self):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
20 description = "VCFFiltering_wrapper"
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
21 description += "\nWrapper for VCFFiltering ;\n VCFFiltering filters SNP on a VCF depending on depth (DP) allele number (AN), allele frequency (AF) and SNP quality.\n"
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
22 description += "example 1 : VCFFiltering.py -f myVCF.vcf -o FilteredVCF.vcf\n"
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
23 description += "example 2 : VCFFiltering.py -f myVCF.vcf -N 2 -F 0.87 -b bed1.bed bed2.bed -o FilteredVCF.vcf\n"
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
24 parser = OptionParser(description = description, version = "0.1")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
25 parser.add_option("-f", "--vcf", dest = "VCFFile", action = "store", type = "string", help = "Input VCF File name [compulsory] [format: VCF]", default = "")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
26 parser.add_option("-o", "--output", dest = "outFile", action = "store", type = "string", help = "output VCF File name [compulsory] [format: VCF]", default = "")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
27 parser.add_option("-m", "--minDP", dest = "minDP", action = "store", type = "int", help = "minimum of depth ; if both minDP and maxDP are set, optimal DP will not be calculated ", default = 0)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
28 parser.add_option("-M", "--maxDP", dest = "maxDP", action = "store", type = "int", help = "maximum of depth ; if both minDP and maxDP are set, optimal DP will not be calculated ", default = 0)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
29 parser.add_option("-N", "--AN", dest = "AN", action = "store", type = "int", help = "maximum number of allele for a SNP; default = 2", default = 2)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
30 parser.add_option("-F", "--AF", dest = "AF", action = "store", type = "float", help = "minimum frequency for the alternative allele of a SNP; default = 0.9", default = 0.9)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
31 parser.add_option("-b", "--bed", dest = "bedFiles", action = "append", type = "string", help = "bed files: list of coordinates to filter, multiple arguments allowed '-b file1 file2' ", default = [])
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
32 parser.add_option("-G", "--graphHTML", dest = "graphHTML", action = "store", type = "string", help = "name of the HTML linking to graphs ", default = "")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
33 parser.add_option("-d", "--dirGraphs", dest = "dirGraphs", action = "store", type = "string", help = "name of the folder containing graphs ", default = "")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
34 options = parser.parse_args()[0]
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
35 self._setAttributesFromOptions(options)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
36
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
37
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
38 def _setAttributesFromOptions(self, options):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
39 self._options = options
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
40
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
41 def run(self):
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
42 if self._options.minDP and self._options.maxDP :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
43 if self._options.minDP > self._options.maxDP :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
44 self.stop_err( 'error in options : minDP > max DP (%s > %s)' %(self._options.minDP,self._options.maxDP))
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
45
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
46 prg = "VCFFiltering.py -g -G 'png' "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
47 args = ""
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
48 args += "-f %s" % self._options.VCFFile
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
49 args += " "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
50 args += "-o %s" % self._options.outFile
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
51 if self._options.AF :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
52 args += " "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
53 args += "-F %s" % self._options.AF
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
54 if self._options.AN :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
55 args += " "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
56 args += "-N %s" % self._options.AN
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
57 if self._options.minDP :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
58 args += " "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
59 args += "-m %s" % self._options.minDP
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
60 if self._options.maxDP :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
61 args += " "
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
62 args += "-M %s" % self._options.maxDP
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
63 if self._options.bedFiles :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
64 if self._options.bedFiles == "":
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
65 pass
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
66 else :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
67 self._lBedFiles = self._options.bedFiles
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
68 cmd = "%s %s" %(prg, args)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
69
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
70 print cmd
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
71
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
72 try:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
73 tmp_err = tempfile.NamedTemporaryFile().name
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
74 tmp_stderr = open( tmp_err, 'wb' )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
75 proc = subprocess.Popen( args=cmd, shell=True, cwd=".", stderr=tmp_stderr )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
76 returncode = proc.wait()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
77 tmp_stderr.close()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
78 # get stderr, allowing for case where it's very large
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
79 tmp_stderr = open( tmp_err, 'rb' )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
80 stderr = ''
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
81 buffsize = 1048576
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
82 try:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
83 while True:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
84 stderr += tmp_stderr.read( buffsize )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
85 if not stderr or len( stderr ) % buffsize != 0:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
86 break
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
87 except OverflowError:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
88 pass
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
89 tmp_stderr.close()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
90 if stderr:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
91 raise Exception, stderr
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
92 except Exception, e:
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
93 self.stop_err( 'Error in VCFFiltering:\n' + str( e ) )
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
94
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
95 if True :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
96 html = open(self._options.graphHTML, "w")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
97
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
98 os.mkdir(self._options.dirGraphs)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
99 lGraphsFiles = glob.glob("VCFFiltering_graphs/*")
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
100 for file in lGraphsFiles :
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
101 baseName = os.path.basename(file)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
102 shutil.move( file ,"%s/%s" %(self._options.dirGraphs, baseName))
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
103 line = "<img src=\"%s\" > \n" %(baseName)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
104 html.write(line)
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
105
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
106
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
107 if __name__ == "__main__":
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
108 iWrapper = VCFFilteringWrapper()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
109 iWrapper.setAttributesFromCmdLine()
3552a8d9f51c Uploaded
urgi-team
parents:
diff changeset
110 iWrapper.run()