[Bioperl-l] Check the hierarchy (level) of nodes in a Newick tree
Ross KK Leung
ross at cuhk.edu.hk
Thu Oct 6 07:07:11 UTC 2011
Hi Tom,
With Dave's and your help, I'm moving forward to the core part.
The previous problem was indeed due to old versions of Perl/Bioperl. Sorry
for being unaware that perl interpreter can be a problem... -_-
Rutger tipped previously:
" You should probably pre-compute the distances from each node to all its
ancestors (post-order traversal) and attach those values as hashes to the
internal nodes. Then, for any assigned node, you just traverse from it to
the root, look up those hashes along the way and compute the distance from
the assigned node to all others that way."
Yet, we still do not have any information about the level of hierarchy. For
the following tree, we say nodes E and F
((A:0.1,B:0.2)E:0.3,(C:0.3,D:0.4)F:0.5)G;
are at the same level. get_lca, ancestor, get_descendant etc cannot help
solve the problem easily because they either return LCA regardless of level
(get_lca(A,C) will return G) or throw exception (by the following codes and
the above tree). I notice that splice() may actually be dealing one level
(ancestor/descent) by one level but it's a bit complicated for me to hack
it...
use Bio::TreeIO;
use Bio::Tree::TreeFunctionsI;
my ($testtree) = @ARGV;
my $treein = Bio::TreeIO->new( -file => $testtree, -format => 'newick' );
while( my $tree = $treein->next_tree ){
%dist_matrix = ();
my @leaves = $tree->get_leaf_nodes;
foreach my $leaf1( @leaves ){
my $id1 = $leaf1->id;
foreach my $leaf2( @leaves ){
my $id2 = $leaf2->id;
if ($id1 eq $id2) {
$dist_matrix{$id1}->{$id2} = 0;
next;
}
my $distance = $tree->distance( -nodes => [$leaf1,
$leaf2] );
print "L1: ";
print $id1; <STDIN>;
print "L2: ";
print $id2; <STDIN>;
my $lca = $tree->get_lca(-nodes => [$leaf1,
$leaf2] );
my $anc = $lca->ancestor;
print "LCA: "; print $lca->id; <STDIN>;
my $distance2 = $tree->distance( -nodes => [$lca,
$anc] );
print "LCA-L2: ";
print $distance2; <STDIN>;
$dist_matrix{$id1}->{$id2} = $distance;
}
}
}
More information about the Bioperl-l
mailing list