Package Bio :: Package Phylo :: Module NexusIO
[hide private]
[frames] | no frames]

Source Code for Module Bio.Phylo.NexusIO

 1  # Copyright (C) 2009 by Eric Talevich (eric.talevich@gmail.com) 
 2  # This code is part of the Biopython distribution and governed by its 
 3  # license. Please see the LICENSE file that should have been included 
 4  # as part of this package. 
 5   
 6  """I/O function wrappers for Bio.Nexus trees.""" 
 7  __docformat__ = "epytext en" 
 8   
 9  from itertools import chain 
10   
11  from Bio.Nexus import Nexus 
12  from Bio.Phylo import Newick 
13   
14  import NewickIO 
15   
16  # Structure of a Nexus tree-only file 
17  NEX_TEMPLATE = """\ 
18  #NEXUS 
19  Begin Taxa; 
20   Dimensions NTax=%(count)d; 
21   TaxLabels %(labels)s; 
22  End; 
23  Begin Trees; 
24   %(trees)s 
25  End; 
26  """ 
27   
28  # 'index' starts from 1; 'tree' is the Newick tree string 
29  TREE_TEMPLATE = "Tree tree%(index)d=[&U]%(tree)s;" 
30   
31   
32 -def parse(handle):
33 """Parse the trees in a Nexus file. 34 35 Uses the old Nexus.Trees parser to extract the trees, converts them back to 36 plain Newick trees, and feeds those strings through the new Newick parser. 37 This way we don't have to modify the Nexus module yet. When we're satisfied 38 with Bio.Phylo, we can change Nexus to use the new NewickIO parser directly. 39 """ 40 nex = Nexus.Nexus(handle) 41 # NB: Once Nexus.Trees is modified to use Tree.Newick objects, do this: 42 # return iter(nex.trees) 43 # Until then, convert the Nexus.Trees.Tree object hierarchy: 44 def node2clade(nxtree, node): 45 subclades = [node2clade(nxtree, nxtree.node(n)) for n in node.succ] 46 return Newick.Clade( 47 branch_length=node.data.branchlength, 48 name=node.data.taxon, 49 clades=subclades, 50 support=node.data.support, 51 comment=node.data.comment)
52 53 for nxtree in nex.trees: 54 newroot = node2clade(nxtree, nxtree.node(nxtree.root)) 55 yield Newick.Tree(root=newroot, rooted=nxtree.rooted, name=nxtree.name, 56 weight=nxtree.weight) 57
58 -def write(obj, handle, **kwargs):
59 """Write a new Nexus file containing the given trees. 60 61 Uses a simple Nexus template and the NewickIO writer to serialize just the 62 trees and minimal supporting info needed for a valid Nexus file. 63 """ 64 trees = list(obj) 65 writer = NewickIO.Writer(trees) 66 nexus_trees = [TREE_TEMPLATE % {'index': idx+1, 'tree': nwk} 67 for idx, nwk in enumerate( 68 writer.to_strings(plain=False, plain_newick=True, 69 **kwargs))] 70 tax_labels = map(str, chain(*(t.get_terminals() for t in trees))) 71 text = NEX_TEMPLATE % { 72 'count': len(tax_labels), 73 'labels': ' '.join(tax_labels), 74 'trees': '\n'.join(nexus_trees), 75 } 76 handle.write(text) 77 return len(nexus_trees)
78