[2] | 1 | #!/usr/bin/env python |
---|
| 2 | |
---|
| 3 | """ |
---|
| 4 | Read a maf file and write out a new maf with only blocks having all of |
---|
| 5 | the passed in species, after dropping any other species and removing columns |
---|
| 6 | containing only gaps. This will attempt to fuse together any blocks |
---|
| 7 | which are adjacent after the unwanted species have been dropped. |
---|
| 8 | |
---|
| 9 | usage: %prog input_maf output_maf species1,species2 |
---|
| 10 | """ |
---|
| 11 | #Dan Blankenberg |
---|
| 12 | import sys |
---|
| 13 | from galaxy import eggs |
---|
| 14 | import pkg_resources; pkg_resources.require( "bx-python" ) |
---|
| 15 | import bx.align.maf |
---|
| 16 | |
---|
| 17 | from bx.align.tools.thread import * |
---|
| 18 | from bx.align.tools.fuse import * |
---|
| 19 | |
---|
| 20 | def main(): |
---|
| 21 | input_file = sys.argv.pop( 1 ) |
---|
| 22 | output_file = sys.argv.pop( 1 ) |
---|
| 23 | species = sys.argv.pop( 1 ).split( ',' ) |
---|
| 24 | |
---|
| 25 | try: |
---|
| 26 | maf_reader = bx.align.maf.Reader( open( input_file ) ) |
---|
| 27 | except: |
---|
| 28 | print >> sys.stderr, "Unable to open source MAF file" |
---|
| 29 | sys.exit() |
---|
| 30 | try: |
---|
| 31 | maf_writer = FusingAlignmentWriter( bx.align.maf.Writer( open( output_file, 'w' ) ) ) |
---|
| 32 | except: |
---|
| 33 | print >> sys.stderr, "Unable to open output file" |
---|
| 34 | sys.exit() |
---|
| 35 | try: |
---|
| 36 | for m in maf_reader: |
---|
| 37 | new_components = m.components |
---|
| 38 | if species != ['None']: |
---|
| 39 | new_components = get_components_for_species( m, species ) |
---|
| 40 | if new_components: |
---|
| 41 | remove_all_gap_columns( new_components ) |
---|
| 42 | m.components = new_components |
---|
| 43 | m.score = 0.0 |
---|
| 44 | maf_writer.write( m ) |
---|
| 45 | except Exception, e: |
---|
| 46 | print >> sys.stderr, "Error steping through MAF File: %s" % e |
---|
| 47 | sys.exit() |
---|
| 48 | maf_reader.close() |
---|
| 49 | maf_writer.close() |
---|
| 50 | |
---|
| 51 | print "Restricted to species: %s." % ", ".join( species ) |
---|
| 52 | |
---|
| 53 | if __name__ == "__main__": main() |
---|