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