Previous changeset 0:0fa83c466e9d (2013-08-29) Next changeset 2:e69151df910c (2013-09-09) |
Commit message:
Populate naive_variant_caller repository. |
added:
dependency_configs/tool_dependencies.xml test-data/fake_phiX174_reads_1.bam test-data/fake_phiX174_reads_1_test_out_1.vcf test-data/phiX174.fasta tool-data/sam_fa_indices.loc.sample tool-data/tool_data_table_conf.xml.sample tools/naive_variant_caller.py tools/naive_variant_caller.xml |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba dependency_configs/tool_dependencies.xml --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/dependency_configs/tool_dependencies.xml Thu Aug 29 10:54:14 2013 -0400 |
b |
@@ -0,0 +1,11 @@ +<tool_dependency> + <package name="numpy" version="1.7.1"> + <repository toolshed="http://toolshed.g2.bx.psu.edu" name="package_numpy_1_7" owner="iuc" changeset_revision="74c21f9bdc39" /> + </package> + <package name="pyBamParser" version="0.0.1"> + <repository toolshed="http://toolshed.g2.bx.psu.edu" name="package_pybamparser_0_0_1" owner="blankenberg" changeset_revision="ae1828dd8bbc" /> + </package> + <package name="pyBamTools" version="0.0.1"> + <repository toolshed="http://toolshed.g2.bx.psu.edu" name="package_pybamtools_0_0_1" owner="blankenberg" changeset_revision="524e3de06009" /> + </package> +</tool_dependency> |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba test-data/fake_phiX174_reads_1.bam |
b |
Binary file test-data/fake_phiX174_reads_1.bam has changed |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba test-data/fake_phiX174_reads_1_test_out_1.vcf --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/fake_phiX174_reads_1_test_out_1.vcf Thu Aug 29 10:54:14 2013 -0400 |
b |
@@ -0,0 +1,51 @@ +##fileformat=VCFv4.1 +##INFO=<ID=AC,Number=A,Type=Integer,Description="Allele count in genotypes, for each ALT allele, in the same order as listed"> +##INFO=<ID=AF,Number=A,Type=Float,Description="Allele Frequency, for each ALT allele, in the same order as listed"> +##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype"> +##FORMAT=<ID=AC,Number=.,Type=Integer,Description="Allele count in genotypes, for each ALT allele, in the same order as listed"> +##FORMAT=<ID=AF,Number=.,Type=Float,Description="Allele Frequency, for each ALT allele, in the same order as listed"> +##FORMAT=<ID=NC,Number=.,Type=String,Description="Nucleotide and indel counts"> +#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT A Fake phiX Sample +phiX174 1411 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=1, +phiX174 1412 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=3, +phiX174 1413 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=5, +phiX174 1414 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=6, +phiX174 1415 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=7, +phiX174 1416 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=8, +phiX174 1417 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=9, +phiX174 1418 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=10, +phiX174 1419 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1420 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1421 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1422 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=10, +phiX174 1423 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1424 . C A . . AC=7;AF=0.7 GT:AC:AF:NC 1/0:7:0.7:A=7,C=3, +phiX174 1425 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1426 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=10, +phiX174 1427 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1428 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1429 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1430 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1431 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1432 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=10, +phiX174 1433 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1434 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1435 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1436 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1437 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1438 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1439 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=10, +phiX174 1440 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1441 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=10, +phiX174 1442 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1443 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=10, +phiX174 1444 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1445 . C . . . AC=;AF= GT:AC:AF:NC 0/0:::C=10, +phiX174 1446 . C T . . AC=3;AF=0.3 GT:AC:AF:NC 0/1:3:0.3:C=7,T=3, +phiX174 1447 . T A . . AC=2;AF=0.222222222222 GT:AC:AF:NC 0/0:2:0.222222222222:A=2,T=7, +phiX174 1448 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=7, +phiX174 1449 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=5, +phiX174 1450 . T . . . AC=;AF= GT:AC:AF:NC 0/0:::T=4, +phiX174 1451 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=3, +phiX174 1452 . A . . . AC=;AF= GT:AC:AF:NC 0/0:::A=2, +phiX174 1453 . G . . . AC=;AF= GT:AC:AF:NC 0/0:::G=1, |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba test-data/phiX174.fasta --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/phiX174.fasta Thu Aug 29 10:54:14 2013 -0400 |
b |
@@ -0,0 +1,79 @@ +>phiX174 +GAGTTTTATCGCTTCCATGACGCAGAAGTTAACACTTTCGGATATTTCTGATGAGTCGAAAAATTATCTT +GATAAAGCAGGAATTACTACTGCTTGTTTACGAATTAAATCGAAGTGGACTGCTGGCGGAAAATGAGAAA +ATTCGACCTATCCTTGCGCAGCTCGAGAAGCTCTTACTTTGCGACCTTTCGCCATCAACTAACGATTCTG +TCAAAAACTGACGCGTTGGATGAGGAGAAGTGGCTTAATATGCTTGGCACGTTCGTCAAGGACTGGTTTA +GATATGAGTCACATTTTGTTCATGGTAGAGATTCTCTTGTTGACATTTTAAAAGAGCGTGGATTACTATC +TGAGTCCGATGCTGTTCAACCACTAATAGGTAAGAAATCATGAGTCAAGTTACTGAACAATCCGTACGTT +TCCAGACCGCTTTGGCCTCTATTAAGCTCATTCAGGCTTCTGCCGTTTTGGATTTAACCGAAGATGATTT +CGATTTTCTGACGAGTAACAAAGTTTGGATTGCTACTGACCGCTCTCGTGCTCGTCGCTGCGTTGAGGCT +TGCGTTTATGGTACGCTGGACTTTGTGGGATACCCTCGCTTTCCTGCTCCTGTTGAGTTTATTGCTGCCG +TCATTGCTTATTATGTTCATCCCGTCAACATTCAAACGGCCTGTCTCATCATGGAAGGCGCTGAATTTAC +GGAAAACATTATTAATGGCGTCGAGCGTCCGGTTAAAGCCGCTGAATTGTTCGCGTTTACCTTGCGTGTA +CGCGCAGGAAACACTGACGTTCTTACTGACGCAGAAGAAAACGTGCGTCAAAAATTACGTGCAGAAGGAG +TGATGTAATGTCTAAAGGTAAAAAACGTTCTGGCGCTCGCCCTGGTCGTCCGCAGCCGTTGCGAGGTACT +AAAGGCAAGCGTAAAGGCGCTCGTCTTTGGTATGTAGGTGGTCAACAATTTTAATTGCAGGGGCTTCGGC +CCCTTACTTGAGGATAAATTATGTCTAATATTCAAACTGGCGCCGAGCGTATGCCGCATGACCTTTCCCA +TCTTGGCTTCCTTGCTGGTCAGATTGGTCGTCTTATTACCATTTCAACTACTCCGGTTATCGCTGGCGAC +TCCTTCGAGATGGACGCCGTTGGCGCTCTCCGTCTTTCTCCATTGCGTCGTGGCCTTGCTATTGACTCTA +CTGTAGACATTTTTACTTTTTATGTCCCTCATCGTCACGTTTATGGTGAACAGTGGATTAAGTTCATGAA +GGATGGTGTTAATGCCACTCCTCTCCCGACTGTTAACACTACTGGTTATATTGACCATGCCGCTTTTCTT +GGCACGATTAACCCTGATACCAATAAAATCCCTAAGCATTTGTTTCAGGGTTATTTGAATATCTATAACA +ACTATTTTAAAGCGCCGTGGATGCCTGACCGTACCGAGGCTAACCCTAATGAGCTTAATCAAGATGATGC +TCGTTATGGTTTCCGTTGCTGCCATCTCAAAAACATTTGGACTGCTCCGCTTCCTCCTGAGACTGAGCTT +TCTCGCCAAATGACGACTTCTACCACATCTATTGACATTATGGGTCTGCAAGCTGCTTATGCTAATTTGC +ATACTGACCAAGAACGTGATTACTTCATGCAGCGTTACCGTGATGTTATTTCTTCATTTGGAGGTAAAAC +CTCTTATGACGCTGACAACCGTCCTTTACTTGTCATGCGCTCTAATCTCTGGGCATCTGGCTATGATGTT +GATGGAACTGACCAAACGTCGTTAGGCCAGTTTTCTGGTCGTGTTCAACAGACCTATAAACATTCTGTGC +CGCGTTTCTTTGTTCCTGAGCATGGCACTATGTTTACTCTTGCGCTTGTTCGTTTTCCGCCTACTGCGAC +TAAAGAGATTCAGTACCTTAACGCTAAAGGTGCTTTGACTTATACCGATATTGCTGGCGACCCTGTTTTG +TATGGCAACTTGCCGCCGCGTGAAATTTCTATGAAGGATGTTTTCCGTTCTGGTGATTCGTCTAAGAAGT +TTAAGATTGCTGAGGGTCAGTGGTATCGTTATGCGCCTTCGTATGTTTCTCCTGCTTATCACCTTCTTGA +AGGCTTCCCATTCATTCAGGAACCGCCTTCTGGTGATTTGCAAGAACGCGTACTTATTCGCCACCATGAT +TATGACCAGTGTTTCCAGTCCGTTCAGTTGTTGCAGTGGAATAGTCAGGTTAAATTTAATGTGACCGTTT +ATCGCAATCTGCCGACCACTCGCGATTCAATCATGACTTCGTGATAAAAGATTGAGTGTGAGGTTATAAC +GCCGAAGCGGTAAAAATTTTAATTTTTGCCGCTGAGGGGTTGACCAAGCGAAGCGCGGTAGGTTTTCTGC +TTAGGAGTTTAATCATGTTTCAGACTTTTATTTCTCGCCATAATTCAAACTTTTTTTCTGATAAGCTGGT +TCTCACTTCTGTTACTCCAGCTTCTTCGGCACCTGTTTTACAGACACCTAAAGCTACATCGTCAACGTTA +TATTTTGATAGTTTGACGGTTAATGCTGGTAATGGTGGTTTTCTTCATTGCATTCAGATGGATACATCTG +TCAACGCCGCTAATCAGGTTGTTTCTGTTGGTGCTGATATTGCTTTTGATGCCGACCCTAAATTTTTTGC +CTGTTTGGTTCGCTTTGAGTCTTCTTCGGTTCCGACTACCCTCCCGACTGCCTATGATGTTTATCCTTTG +AATGGTCGCCATGATGGTGGTTATTATACCGTCAAGGACTGTGTGACTATTGACGTCCTTCCCCGTACGC +CGGGCAATAATGTTTATGTTGGTTTCATGGTTTGGTCTAACTTTACCGCTACTAAATGCCGCGGATTGGT +TTCGCTGAATCAGGTTATTAAAGAGATTATTTGTCTCCAGCCACTTAAGTGAGGTGATTTATGTTTGGTG +CTATTGCTGGCGGTATTGCTTCTGCTCTTGCTGGTGGCGCCATGTCTAAATTGTTTGGAGGCGGTCAAAA +AGCCGCCTCCGGTGGCATTCAAGGTGATGTGCTTGCTACCGATAACAATACTGTAGGCATGGGTGATGCT +GGTATTAAATCTGCCATTCAAGGCTCTAATGTTCCTAACCCTGATGAGGCCGCCCCTAGTTTTGTTTCTG +GTGCTATGGCTAAAGCTGGTAAAGGACTTCTTGAAGGTACGTTGCAGGCTGGCACTTCTGCCGTTTCTGA +TAAGTTGCTTGATTTGGTTGGACTTGGTGGCAAGTCTGCCGCTGATAAAGGAAAGGATACTCGTGATTAT +CTTGCTGCTGCATTTCCTGAGCTTAATGCTTGGGAGCGTGCTGGTGCTGATGCTTCCTCTGCTGGTATGG +TTGACGCCGGATTTGAGAATCAAAAAGAGCTTACTAAAATGCAACTGGACAATCAGAAAGAGATTGCCGA +GATGCAAAATGAGACTCAAAAAGAGATTGCTGGCATTCAGTCGGCGACTTCACGCCAGAATACGAAAGAC +CAGGTATATGCACAAAATGAGATGCTTGCTTATCAACAGAAGGAGTCTACTGCTCGCGTTGCGTCTATTA +TGGAAAACACCAATCTTTCCAAGCAACAGCAGGTTTCCGAGATTATGCGCCAAATGCTTACTCAAGCTCA +AACGGCTGGTCAGTATTTTACCAATGACCAAATCAAAGAAATGACTCGCAAGGTTAGTGCTGAGGTTGAC +TTAGTTCATCAGCAAACGCAGAATCAGCGGTATGGCTCTTCTCATATTGGCGCTACTGCAAAGGATATTT +CTAATGTCGTCACTGATGCTGCTTCTGGTGTGGTTGATATTTTTCATGGTATTGATAAAGCTGTTGCCGA +TACTTGGAACAATTTCTGGAAAGACGGTAAAGCTGATGGTATTGGCTCTAATTTGTCTAGGAAATAACCG +TCAGGATTGACACCCTCCCAATTGTATGTTTTCATGCCTCCAAATCTTGGAGGCTTTTTTATGGTTCGTT +CTTATTACCCTTCTGAATGTCACGCTGATTATTTTGACTTTGAGCGTATCGAGGCTCTTAAACCTGCTAT +TGAGGCTTGTGGCATTTCTACTCTTTCTCAATCCCCAATGCTTGGCTTCCATAAGCAGATGGATAACCGC +ATCAAGCTCTTGGAAGAGATTCTGTCTTTTCGTATGCAGGGCGTTGAGTTCGATAATGGTGATATGTATG +TTGACGGCCATAAGGCTGCTTCTGACGTTCGTGATGAGTTTGTATCTGTTACTGAGAAGTTAATGGATGA +ATTGGCACAATGCTACAATGTGCTCCCCCAACTTGATATTAATAACACTATAGACCACCGCCCCGAAGGG +GACGAAAAATGGTTTTTAGAGAACGAGAAGACGGTTACGCAGTTTTGCCGCAAGCTGGCTGCTGAACGCC +CTCTTAAGGATATTCGCGATGAGTATAATTACCCCAAAAAGAAAGGTATTAAGGATGAGTGTTCAAGATT +GCTGGAGGCCTCCACTATGAAATCGCGTAGAGGCTTTACTATTCAGCGTTTGATGAATGCAATGCGACAG +GCTCATGCTGATGGTTGGTTTATCGTTTTTGACACTCTCACGTTGGCTGACGACCGATTAGAGGCGTTTT +ATGATAATCCCAATGCTTTGCGTGACTATTTTCGTGATATTGGTCGTATGGTTCTTGCTGCCGAGGGTCG +CAAGGCTAATGATTCACACGCCGACTGCTATCAGTATTTTTGTGTGCCTGAGTATGGTACAGCTAATGGC +CGTCTTCATTTCCATGCGGTGCATTTTATGCGGACACTTCCTACAGGTAGCGTTGACCCTAATTTTGGTC +GTCGGGTACGCAATCGCCGCCAGTTAAATAGCTTGCAAAATACGTGGCCTTATGGTTACAGTATGCCCAT +CGCAGTTCGCTACACGCAGGACGCTTTTTCACGTTCTGGTTGGTTGTGGCCTGTTGATGCTAAAGGTGAG +CCGCTTAAAGCTACCAGTTATATGGCTGTTGGTTTCTATGTGGCTAAATACGTTAACAAAAAGTCAGATA +TGGACCTTGCTGCTAAAGGTCTAGGAGCTAAAGAATGGAACAACTCACTAAAAACCAAGCTGTCGCTACT +TCCCAAGAAGCTGTTCAGAATCAGAATGAGCCGCAACTTCGGGATGAAAATGCTCACAATGACAAATCTG +TCCACGGAGTGCTTAATCCAACTTACCAAGCTGGGTTACGACGCGACGCCGTTCAACCAGATATTGAAGC +AGAACGCAAAAAGAGAGATGAGATTGAGGCTGGGAAAAGTTACTGTAGCCGACGTTTTGGCGGCGCAACC +TGTGACGACAAATCTGCTCAAATTTATGCGCGCTTCGATAAAAATGATTGGCGTATCCAACCTGCA + |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba tool-data/sam_fa_indices.loc.sample --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tool-data/sam_fa_indices.loc.sample Thu Aug 29 10:54:14 2013 -0400 |
b |
@@ -0,0 +1,28 @@ +#This is a sample file distributed with Galaxy that enables tools +#to use a directory of Samtools indexed sequences data files. You will need +#to create these data files and then create a sam_fa_indices.loc file +#similar to this one (store it in this directory) that points to +#the directories in which those files are stored. The sam_fa_indices.loc +#file has this format (white space characters are TAB characters): +# +#index <seq> <location> +# +#So, for example, if you had hg18 indexed stored in +#/depot/data2/galaxy/sam/, +#then the sam_fa_indices.loc entry would look like this: +# +#index hg18 /depot/data2/galaxy/sam/hg18.fa +# +#and your /depot/data2/galaxy/sam/ directory +#would contain hg18.fa and hg18.fa.fai files: +# +#-rw-r--r-- 1 james universe 830134 2005-09-13 10:12 hg18.fa +#-rw-r--r-- 1 james universe 527388 2005-09-13 10:12 hg18.fa.fai +# +#Your sam_fa_indices.loc file should include an entry per line for +#each index set you have stored. The file in the path does actually +#exist, but it should never be directly used. Instead, the name serves +#as a prefix for the index file. For example: +# +#index hg18 /depot/data2/galaxy/sam/hg18.fa +#index hg19 /depot/data2/galaxy/sam/hg19.fa |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba tool-data/tool_data_table_conf.xml.sample --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tool-data/tool_data_table_conf.xml.sample Thu Aug 29 10:54:14 2013 -0400 |
b |
@@ -0,0 +1,7 @@ +<tables> + <!-- Location of SAMTools indexes and other files --> + <table name="sam_fa_indexes" comment_char="#"> + <columns>line_type, value, path</columns> + <file path="tool-data/sam_fa_indices.loc" /> + </table> +</tables> |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba tools/naive_variant_caller.py --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tools/naive_variant_caller.py Thu Aug 29 10:54:14 2013 -0400 |
[ |
@@ -0,0 +1,64 @@ +#Dan Blankenberg +import sys +import optparse + +from pyBamParser.bam import Reader +from pyBamTools.genotyping.naive import VCFReadGroupGenotyper + +def main(): + #Parse Command Line + parser = optparse.OptionParser() + parser.add_option( '-b', '--bam', dest='bam_file', action='append', type="string", default=[], help='BAM filename, optionally index filename. Multiple allowed.' ) + parser.add_option( '-i', '--index', dest='index_file', action='append', type="string", default=[], help='optionally index filename. Multiple allowed.' ) + parser.add_option( '-o', '--output_vcf_filename', dest='output_vcf_filename', action='store', default = None, type="string", help='Output VCF filename' ) + parser.add_option( '-r', '--reference_genome_filename', dest='reference_genome_filename', action='store', default = None, type="string", help='Input reference file' ) + parser.add_option( '-v', '--variants_only', dest='variants_only', action='store_true', default = False, help='Report only sites with a possible variant allele.' ) + parser.add_option( '-s', '--use_strand', dest='use_strand', action='store_true', default = False, help='Report counts by strand' ) + parser.add_option( '-p', '--ploidy', dest='ploidy', action='store', type="int", default=2, help='Ploidy. Default=2.' ) + parser.add_option( '-d', '--min_support_depth', dest='min_support_depth', action='store', type="int", default=0, help='Minimum number of reads needed to consider a REF/ALT. Default=0.' ) + parser.add_option( '-q', '--min_base_quality', dest='min_base_quality', action='store', type="int", default=None, help='Minimum base quality.' ) + parser.add_option( '-m', '--min_mapping_quality', dest='min_mapping_quality', action='store', type="int", default=None, help='Minimum mapping.' ) + parser.add_option( '-t', '--coverage_dtype', dest='coverage_dtype', action='store', type="string", default='uint8', help='dtype to use for coverage array' ) + parser.add_option( '--allow_out_of_bounds_positions', dest='allow_out_of_bounds_positions', action='store_true', default = False, help='Allows out of bounds positions to not throw fatal errors' ) + parser.add_option( '--region', dest='region', action='append', type="string", default=[], help='region' ) + (options, args) = parser.parse_args() + + if len( options.bam_file ) == 0: + print >>sys.stderr, 'You must provide at least one bam (-b) file.' + parser.print_help( sys.stderr ) + sys.exit( 1 ) + if options.index_file: + assert len( options.index_file ) == len( options.bam_file ), "If you provide a name for an index file, you must provide the index name for all bam files." + bam_files = zip( options.bam_file, options.index_file ) + else: + bam_files = [ ( x, ) for x in options.bam_file ] + if not options.reference_genome_filename: + print >> sys.stderr, "Warning: Reference file has not been specified. Providing a reference genome is highly recommended." + if options.output_vcf_filename: + out = open( options.output_vcf_filename, 'wb' ) + else: + out = sys.stdout + + regions = [] + if options.region: + for region in options.region: + region_split = region.split( ":" ) + region = region_split.pop( 0 ) + if region_split: + region_split = filter( bool, region_split[0].split( '-' ) ) + if region_split: + if len( region_split ) != 2: + print >> sys.stderr, "You must specify both a start and an end, or only a chromosome when specifying regions." + cleanup_before_exit( tmp_dir ) + sys.exit( 1 ) + region = tuple( [ region ] + map( int, region_split ) ) + regions.append( region ) + + coverage = VCFReadGroupGenotyper( map( lambda x: Reader( *x ), bam_files ), options.reference_genome_filename, dtype=options.coverage_dtype, + min_support_depth=options.min_support_depth, min_base_quality=options.min_base_quality, min_mapping_quality=options.min_mapping_quality, + restrict_regions=regions, use_strand=options.use_strand, allow_out_of_bounds_positions=options.allow_out_of_bounds_positions ) + for line in coverage.iter_vcf( ploidy=options.ploidy, variants_only=options.variants_only ): + out.write( "%s\n" % line ) + out.close() + +if __name__ == "__main__": main() |
b |
diff -r 0fa83c466e9d -r ae6edc0012ba tools/naive_variant_caller.xml --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tools/naive_variant_caller.xml Thu Aug 29 10:54:14 2013 -0400 |
b |
b'@@ -0,0 +1,208 @@\n+<tool id="naive_variant_caller" name="Naive Variant Caller" version="0.0.1">\n+ <description> - tabulate variable sites from BAM datasets</description>\n+ <requirements>\n+ <requirement type="package" version="1.7.1">numpy</requirement>\n+ <requirement type="package" version="0.0.1">pyBamParser</requirement>\n+ <requirement type="package" version="0.0.1">pyBamTools</requirement>\n+ </requirements>\n+ <stdio>\n+ <exit_code range="1:" err_level="fatal" />\n+ <exit_code range=":-1" err_level="fatal" />\n+ </stdio>\n+ <command interpreter="python">naive_variant_caller.py\n+ -o "${output_vcf}"\n+ \n+ #for $input_bam in $reference_source.input_bams:\n+ -b "${input_bam.input_bam}"\n+ -i "${input_bam.input_bam.metadata.bam_index}"\n+ #end for\n+ \n+ #if $reference_source.reference_source_selector != "history":\n+ -r "${reference_source.ref_file.fields.path}"\n+ #elif $reference_source.ref_file:\n+ -r "${reference_source.ref_file}"\n+ #end if\n+ \n+ #for $region in $regions:\n+ --region "${region.chromosome}:${region.start}-${region.end}"\n+ #end for\n+ \n+ ${variants_only}\n+ \n+ ${use_strand}\n+ \n+ --ploidy "${$ploidy}"\n+ \n+ --min_support_depth "${min_support_depth}"\n+ \n+ #if str($min_base_quality):\n+ --min_base_quality "${min_base_quality}"\n+ #end if\n+ \n+ #if str($min_mapping_quality):\n+ --min_mapping_quality "${min_mapping_quality}"\n+ #end if\n+ \n+ --coverage_dtype "${coverage_dtype}"\n+ \n+ --allow_out_of_bounds_positions\n+ \n+ </command>\n+ <inputs>\n+ <conditional name="reference_source">\n+ <param name="reference_source_selector" type="select" label="Choose the source for the reference list">\n+ <option value="cached">Locally cached</option>\n+ <option value="history">History</option>\n+ </param>\n+ <when value="cached">\n+ <repeat name="input_bams" title="BAM file" min="1" >\n+ <param name="input_bam" type="data" format="bam" label="BAM file">\n+ <validator type="unspecified_build" />\n+ <validator type="dataset_metadata_in_data_table" table_name="sam_fa_indexes" metadata_name="dbkey" metadata_column="value" message="Sequences are not currently available for the specified build." /> <!-- fixme!!! this needs to be a select -->\n+ </param>\n+ </repeat>\n+ <param name="ref_file" type="select" label="Using reference genome" >\n+ <options from_data_table="sam_fa_indexes">\n+ <!-- <filter type="data_meta" key="dbkey" ref="input_bam" column="dbkey"/> does not yet work in a repeat...--> \n+ </options>\n+ <validator type="no_options" message="A built-in reference genome is not available for the build associated with the selected input file"/>\n+ </param>\n+ </when>\n+ <when value="history"> <!-- FIX ME!!!! -->\n+ <repeat name="input_bams" title="BAM file" min="1" >\n+ <param name="input_bam" type="data" format="bam" label="BAM file" >\n+ </param>\n+ </repeat>\n+ <param name="ref_file" type="data" format="fasta" label="Using reference file" optional="True" />\n+ </when>\n+ </conditional>\n+\n+ <repeat name="regions" title="Restrict to regions" min="0" >\n+ <param name="chromosome" type="text" value="" optional="False" label="Chromosome" />\n+ <param name="start" type="integer" value="" optional="True" label="Start" />\n+ <param name="end" type="integer" value="" optional="True" label="End" />\n+ </repeat>\n+\n+ <!-- TODO: enhance filtering -->\n+ <param name="min_support_depth" type="integer" value="0" min="0" label="Minimum number of reads needed to consider a REF/ALT" />\n+ <param name="min_base_quality" type="integer" value="" label="Minimum base quality" optional="True" />\n+ <param name="min_mapping_quality" type="integer" value="" label="Minimum m'..b'strand, if the strandedness option was specified.\n+ \n+\n+------\n+\n+**Inputs**\n+\n+Accepts one or more BAM input files and a reference genome from the built-in list or from a FASTA file in your history.\n+\n+\n+**Outputs**\n+\n+The output is in VCF format.\n+\n+Example VCF output line, without reporting by strand:\n+ ``chrM\t16029\t.\tT\tG,A,C\t.\t.\tAC=15,9,5;AF=0.00155311658729,0.000931869952371,0.000517705529095\tGT:AC:AF:NC\t0/0:15,9,5:0.00155311658729,0.000931869952371,0.000517705529095:A=9,C=5,T=9629,G=15,``\n+\n+Example VCF output line, when reporting by strand:\n+ ``chrM\t16029\t.\tT\tG,A,C\t.\t.\tAC=15,9,5;AF=0.00155311658729,0.000931869952371,0.000517705529095\tGT:AC:AF:NC\t0/0:15,9,5:0.00155311658729,0.000931869952371,0.000517705529095:+T=3972,-A=9,-C=5,-T=5657,-G=15,``\n+\n+**Options**\n+\n+Reference Genome:\n+\n+ Ensure that you have selected the correct reference genome, either from the list of built-in genomes or by selecting the corresponding FASTA file from your history.\n+\n+Restrict to regions:\n+\n+ You can specify any number of regions on which you would like to receive results. You can specify just a chromosome name, or a chromosome name and start postion, or a chromosome name and start and end position for the set of desired regions. \n+\n+Minimum number of reads needed to consider a REF/ALT:\n+\n+ This value declares the minimum number of reads containing a particular base at each position in order to list and use said allele in genotyping calls. Default is 0.\n+\n+Minimum base quality:\n+\n+ The minimum base quality score needed for the position in a read to be used for nucleotide counts and genotyping. Default is no filter.\n+\n+Minimum mapping quality:\n+\n+ The minimum mapping quality score needed to consider a read for nucleotide counts and genotyping. Default is no filter.\n+\n+Ploidy:\n+\n+ The number of genotype calls to make at each reported position.\n+\n+Only write out positions with with possible alternate alleles:\n+\n+ When set, only positions which have at least one non-reference nucleotide which passes declare filters will be present in the output.\n+\n+Report counts by strand:\n+\n+ When set, nucleotide counts (NC) will be reported in reference to the aligned read\'s source strand. Reported as: <strand><BASE>=<COUNT>.\n+\n+Choose the dtype to use for storing coverage information:\n+\n+ This controls the maximum depth value for each nucleotide/position/strand (when specified). Smaller values require the least amount of memory, but have smaller maximal limits.\n+\n+ +--------+----------------------------+\n+ | name | maximum coverage value |\n+ +========+============================+\n+ | uint8 | 255 |\n+ +--------+----------------------------+\n+ | uint16 | 65,535 |\n+ +--------+----------------------------+\n+ | uint32 | 4,294,967,295 |\n+ +--------+----------------------------+\n+ | uint64 | 18,446,744,073,709,551,615 |\n+ +--------+----------------------------+\n+\n+------\n+\n+**Citation**\n+\n+If you use this tool, please cite Blankenberg D, et al. *In preparation.*\n+\n+ </help>\n+ <tests>\n+ <test>\n+ <param name="reference_source_selector" value="history" />\n+ <param name="input_bam" value="fake_phiX174_reads_1.bam" ftype="bam" /> \n+ <param name="ref_file" value="phiX174.fasta" ftype="fasta" />\n+ <param name="regions" value="0" />\n+ <param name="min_support_depth" value="0" />\n+ <param name="min_base_quality" value="" />\n+ <param name="min_mapping_quality" value="" />\n+ <param name="ploidy" value="2" />\n+ <param name="variants_only" value="False" />\n+ <param name="use_strand" value="False" />\n+ <param name="coverage_dtype" value="uint8" />\n+ <output name="output_vcf" file="fake_phiX174_reads_1_test_out_1.vcf" compare="contains" />\n+ </test>\n+ </tests>\n+ \n+</tool>\n' |