Mercurial > repos > xuebing > sharplabtool
comparison tools/maf/maf_limit_to_species.py @ 0:9071e359b9a3
Uploaded
author | xuebing |
---|---|
date | Fri, 09 Mar 2012 19:37:19 -0500 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:9071e359b9a3 |
---|---|
1 #!/usr/bin/env python | |
2 | |
3 """ | |
4 Read a maf file and write out a new maf with only blocks having the | |
5 required species, after dropping any other species and removing | |
6 columns containing only gaps. | |
7 | |
8 usage: %prog species,species2,... input_maf output_maf allow_partial min_species_per_block | |
9 """ | |
10 #Dan Blankenberg | |
11 from galaxy import eggs | |
12 import pkg_resources; pkg_resources.require( "bx-python" ) | |
13 import bx.align.maf | |
14 from galaxy.tools.util import maf_utilities | |
15 import sys | |
16 | |
17 assert sys.version_info[:2] >= ( 2, 4 ) | |
18 | |
19 def main(): | |
20 | |
21 species = maf_utilities.parse_species_option( sys.argv[1] ) | |
22 if species: | |
23 spec_len = len( species ) | |
24 else: | |
25 spec_len = 0 | |
26 try: | |
27 maf_reader = bx.align.maf.Reader( open( sys.argv[2],'r' ) ) | |
28 maf_writer = bx.align.maf.Writer( open( sys.argv[3],'w' ) ) | |
29 except: | |
30 print >>sys.stderr, "Your MAF file appears to be malformed." | |
31 sys.exit() | |
32 allow_partial = False | |
33 if int( sys.argv[4] ): allow_partial = True | |
34 min_species_per_block = int( sys.argv[5] ) | |
35 | |
36 maf_blocks_kept = 0 | |
37 for m in maf_reader: | |
38 if species: | |
39 m = m.limit_to_species( species ) | |
40 m.remove_all_gap_columns() | |
41 spec_in_block_len = len( maf_utilities.get_species_in_block( m ) ) | |
42 if ( not species or allow_partial or spec_in_block_len == spec_len ) and spec_in_block_len > min_species_per_block: | |
43 maf_writer.write( m ) | |
44 maf_blocks_kept += 1 | |
45 | |
46 maf_reader.close() | |
47 maf_writer.close() | |
48 | |
49 print "Restricted to species: %s." % ", ".join( species ) | |
50 print "%i MAF blocks have been kept." % maf_blocks_kept | |
51 | |
52 if __name__ == "__main__": | |
53 main() |