root/galaxy-central/eggs/bx_python-0.5.0_dev_f74aec067563-py2.6-macosx-10.6-universal-ucs2.egg/EGG-INFO/scripts/maf_to_axt.py

リビジョン 3, 2.4 KB (コミッタ: kohda, 14 年 前)

Install Unix tools  http://hannonlab.cshl.edu/galaxy_unix_tools/galaxy.html

行番号 
1#!/usr/bin/python2.6
2"""
3Application to convert MAF file to AXT file, projecting to any two species.
4Reads a MAF file from standard input and writes an AXT file to standard out; 
5some statistics are written to standard error.  The user must specify the
6two species of interest.
7
8usage: %prog primary_species secondary_species < maf_file > axt_file
9"""
10
11__author__ = "Bob Harris (rsharris@bx.psu.edu)"
12
13import sys
14import copy
15import bx.align.maf
16import bx.align.axt
17
18def usage(s=None):
19        message = """
20maf_to_axt primary_species secondary_species < maf_file > axt_file
21"""
22        if (s == None): sys.exit (message)
23        else:           sys.exit ("%s\n%s" % (s,message))
24
25
26def main():
27
28        # parse the command line
29
30        primary   = None
31        secondary = None
32
33        args = sys.argv[1:]
34        while (len(args) > 0):
35                arg = args.pop(0)
36                val = None
37                fields = arg.split("=",1)
38                if (len(fields) == 2):
39                        arg = fields[0]
40                        val = fields[1]
41                        if (val == ""):
42                                usage("missing a value in %s=" % arg)
43
44                if (primary == None) and (val == None):
45                        primary = arg
46                elif (secondary == None) and (val == None):
47                        secondary = arg
48                else:
49                        usage("unknown argument: %s" % arg)
50
51        if (primary == None):
52                usage("missing primary species")
53
54        if (secondary == None):
55                usage("missing secondary species")
56
57        # read the alignments and other info
58
59        out = bx.align.axt.Writer(sys.stdout)
60
61        axtsRead = 0
62        mafsWritten = 0
63        for mafBlock in bx.align.maf.Reader(sys.stdin):
64                axtsRead += 1
65
66                p = mafBlock.get_component_by_src_start(primary)
67                if (p == None): continue
68                s = mafBlock.get_component_by_src_start(secondary)
69                if (s == None): continue
70
71                axtBlock = bx.align.Alignment (mafBlock.score, mafBlock.attributes)
72                axtBlock.add_component (clone_component(p))
73                axtBlock.add_component (clone_component(s))
74
75                remove_mutual_gaps (axtBlock)
76                if (axtBlock.text_size == 0):
77                        continue
78
79                out.write (axtBlock)
80                mafsWritten += 1
81
82        sys.stderr.write ("%d blocks read, %d written\n" % (axtsRead,mafsWritten))
83
84
85def clone_component(c):
86        return bx.align.Component (c.src, c.start, c.size, c.strand, c.src_size, \
87                                   copy.copy(c.text))
88
89
90def remove_mutual_gaps (block):
91
92        if (len(block.components) == 0): return
93
94        nonGaps = []
95
96        for c in block.components:
97                for ix in range(0,block.text_size):
98                        if (ix not in nonGaps) and (c.text[ix] != "-"):
99                                nonGaps.append(ix)
100
101        nonGaps.sort()
102
103        for c in block.components:
104                c.text = "".join([c.text[ix] for ix in nonGaps])
105
106        block.text_size = len(nonGaps)
107
108
109if __name__ == "__main__": main()
110
Note: リポジトリブラウザについてのヘルプは TracBrowser を参照してください。