[3] | 1 | #!/usr/bin/env python2.4 |
---|
| 2 | """ |
---|
| 3 | Returns all positions of a maf with any pwm score > threshold |
---|
| 4 | The positions are projected onto human coordinates |
---|
| 5 | """ |
---|
| 6 | |
---|
| 7 | import psyco_full |
---|
| 8 | from bx.align import maf as align_maf |
---|
| 9 | import bx.pwm.position_weight_matrix as pwmx |
---|
| 10 | from bx.pwm.pwm_score_maf import MafMotifSelect |
---|
| 11 | import sys |
---|
| 12 | from bx import intervals |
---|
| 13 | |
---|
| 14 | def isnan(x): |
---|
| 15 | return not x==x |
---|
| 16 | |
---|
| 17 | def main(): |
---|
| 18 | |
---|
| 19 | if len(sys.argv) < 5: |
---|
| 20 | print >>sys.stderr, "%s transfac|basic pwmfile inmaf threshold [motif]" % sys.argv[0] |
---|
| 21 | sys.exit(2) |
---|
| 22 | |
---|
| 23 | r = pwmx.Reader(open(sys.argv[2]),format=sys.argv[1]) |
---|
| 24 | pwm = iter(r).next() |
---|
| 25 | inmaf = open(sys.argv[3]) |
---|
| 26 | threshold = float(sys.argv[4]) |
---|
| 27 | if len(sys.argv) > 5: motif = sys.argv[5] |
---|
| 28 | else: motif = None |
---|
| 29 | |
---|
| 30 | for maf in align_maf.Reader(inmaf): |
---|
| 31 | for mafmotif,pwm_score,motif_score in MafMotifSelect(maf, pwm, motif, threshold): |
---|
| 32 | #mafwrite( mafmotif,pwm_score,motif_score) |
---|
| 33 | print mafmotif, pwm_score, motif_score |
---|
| 34 | print 'zzzzzzzzzzzzzzzzzzzzzzzzzzzzz' |
---|
| 35 | |
---|
| 36 | def mafwrite(alignment,kvec=None,jvec=None,file=sys.stdout): |
---|
| 37 | file.write( "a score=" + str( alignment.score ) ) |
---|
| 38 | for key in alignment.attributes: |
---|
| 39 | file.write( " %s=%s" % ( key, alignment.attributes[key] ) ) |
---|
| 40 | file.write( "\n" ) |
---|
| 41 | rows = [] |
---|
| 42 | if not kvec: kvec = [0 for c in alignment.components] |
---|
| 43 | if not jvec: jvec = [0 for c in alignment.components] |
---|
| 44 | for c,x,y in zip(alignment.components,kvec,jvec): |
---|
| 45 | #for c in alignment.components: |
---|
| 46 | rows.append( ( "s", c.src, str( c.start ), str( c.size ), c.strand, str( c.src_size ), c.text, "%.2f" % x, str(y) ) ) |
---|
| 47 | #rows.append( ( "s", c.src, str( c.start ), str( c.size ), c.strand, str( c.src_size ), c.text ) ) |
---|
| 48 | file.write( format_tabular( rows, "llrrrrrrr" ) ) |
---|
| 49 | |
---|
| 50 | def format_tabular( rows, align=None ): |
---|
| 51 | if len( rows ) == 0: return "" |
---|
| 52 | lengths = [ len( col ) for col in rows[ 0 ] ] |
---|
| 53 | for row in rows[1:]: |
---|
| 54 | for i in range( 0, len( row ) ): |
---|
| 55 | lengths[ i ] = max( lengths[ i ], len( row[ i ] ) ) |
---|
| 56 | rval = "" |
---|
| 57 | for row in rows: |
---|
| 58 | for i in range( 0, len( row ) ): |
---|
| 59 | if align and align[ i ] == "l": |
---|
| 60 | rval += row[ i ].ljust( lengths[ i ] ) |
---|
| 61 | else: |
---|
| 62 | rval += row[ i ].rjust( lengths[ i ] ) |
---|
| 63 | rval += " " |
---|
| 64 | rval += "\n" |
---|
| 65 | return rval |
---|
| 66 | |
---|
| 67 | |
---|
| 68 | |
---|
| 69 | if __name__ == '__main__': main() |
---|