bx-python / scripts / maf_covered_ranges.py

#!/usr/bin/env python2.3

usage: %prog species1,species2,... nrequired < maf 

import psyco_full

import bx.align.maf
import copy
import sys

from itertools import *
from bx.cookbook import doc_optparse

SPAN = 100
MIN = 100

def main():

    options, args = doc_optparse.parse( __doc__ )

        species = args[0].split(',')
        nrequired = int( args[1] )

    maf_reader = bx.align.maf.Reader( sys.stdin )
    interval_start = None
    interval_end = None   

    for m in maf_reader:            
        ref = m.components[0]
        # Does this alignment have enough of the required species        
        if nrequired <= len( [ comp for comp in m.components if comp.src.split('.')[0] in species ] ):
            if interval_start is None:
                interval_start = ref.start
                interval_end = ref.end
                if ref.start - interval_end < SPAN:
                    interval_end = ref.end
                    if interval_end - interval_start >= MIN: print ref.src.split('.')[1], interval_start, interval_end
                    interval_start = ref.start
                    interval_end = ref.end    
            if interval_start != None and interval_end - interval_start >= MIN: print ref.src.split('.')[1],interval_start,interval_end
            interval_start = None
            interval_end = None

if __name__ == "__main__": main()
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for ProjectModifiedEvent.java.
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.