t/Tree.t

   1 # -*-Perl-*- Test Harness script for Bioperl
   2 # $Id$
   3
   4 use strict;
   5
   6 BEGIN {
   7     use lib 't/lib';
   8     use BioperlTest;
   9
  10     test_begin(-tests => 60);
  11
  12     use_ok('Bio::TreeIO');
  13 }
  14
  15 my $verbose = test_debug();
  16
  17 my $treeio = Bio::TreeIO->new(-verbose => $verbose,
  18                              -format => 'nhx',
  19                              -file   => test_input_file('test.nhx'));
  20 my $tree = $treeio->next_tree;
  21
  22 # tests for tags
  23 ok ! $tree->has_tag('test');
  24 is $tree->add_tag_value('test','a'), 1;
  25 ok $tree->has_tag('test');
  26 is $tree->add_tag_value('test','b'), 2;
  27 my @tags = $tree->get_tag_values('test');
  28 is scalar @tags, 2;
  29 is scalar $tree->get_tag_values('test'), 'a', 'retrieve the first value';
  30 is $tree->remove_tag('test2'), 0;
  31 is $tree->remove_tag('test'), 1;
  32 ok ! $tree->has_tag('test');
  33 is $tree->set_tag_value('test',('a','b','c')), 3;
  34 is $tree->remove_all_tags(), undef;
  35 ok ! $tree->has_tag('test');
  36
  37
  38 my @nodes = $tree->find_node('ADH2');
  39 is(@nodes, 2,'Number of nodes that have ADH2 as name');
  40
  41 if( $verbose ) {
  42     $treeio = Bio::TreeIO->new(-verbose => $verbose,
  43                               -format => 'nhx',
  44                               );
  45     $treeio->write_tree($tree);
  46     print "nodes are: \n",
  47     join(", ", map {  $_->id . ":". (defined $_->branch_length ?
  48                                      $_->branch_length : '' ) } @nodes), "\n";
  49 }
  50
  51 $treeio = Bio::TreeIO->new(-format => 'newick',
  52                           -file   => test_input_file('test.nh'));
  53 $tree = $treeio->next_tree;
  54
  55
  56 if( $verbose ) {
  57     my $out = Bio::TreeIO->new(-format => 'tabtree');
  58
  59     $out->write_tree($tree);
  60 }
  61
  62 my @hADH = ( $tree->find_node('hADH1'),
  63              $tree->find_node('hADH2') );
  64 my ($n4) = $tree->find_node('yADH4');
  65
  66 is($tree->is_monophyletic(-nodes    => \@hADH,
  67                           -outgroup => $n4),1,'Test Monophyly');
  68
  69 my @mixgroup = ( $tree->find_node('hADH1'),
  70                  $tree->find_node('yADH2'),
  71                  $tree->find_node('yADH3'),
  72                  );
  73
  74 my ($iADHX) = $tree->find_node('iADHX');
  75
  76 # test height
  77 is($iADHX->height, 0,'Height');
  78 is($iADHX->depth,0.22,'Depth');
  79 isnt( $tree->is_monophyletic(-nodes   => \@mixgroup,
  80                             -outgroup=> $iADHX),1, 'non-monophyletic group');
  81
  82 # binary tree?
  83 is $tree->is_binary, 0, 'not a binary tree';
  84 is scalar $tree->get_nodes, 12, '12 nodes';
  85 $tree->verbose(-1);
  86 $tree->force_binary;
  87 is $tree->is_binary, 1, 'after force_binary() it is';
  88 is scalar $tree->get_nodes, 17, 'and there are more nodes (17)';
  89
  90 my $in = Bio::TreeIO->new(-format => 'newick',
  91                          -fh     => \*DATA);
  92 $tree = $in->next_tree;
  93 my ($a,$b,$c,$d) = ( $tree->find_node('A'),
  94                      $tree->find_node('B'),
  95                      $tree->find_node('C'),
  96                      $tree->find_node('D'));
  97
  98 is($tree->is_monophyletic(-nodes => [$b,$c],
  99                           -outgroup => $d),1, 'B,C are Monophyletic');
 100
 101 is($tree->is_monophyletic(-nodes => [$b,$a],
 102                           -outgroup => $d),1,'A,B are Monophyletic');
 103
 104 $tree = $in->next_tree;
 105 my ($e,$f,$i);
 106 ($a,$b,$c,$d,$e,$f,$i) = ( $tree->find_node('A'),
 107                            $tree->find_node('B'),
 108                            $tree->find_node('C'),
 109                            $tree->find_node('D'),
 110                            $tree->find_node('E'),
 111                            $tree->find_node('F'),
 112                            $tree->find_node('I'),
 113                            );
 114 isnt( $tree->is_monophyletic(-nodes => [$b,$f],
 115                             -outgroup => $d),1,'B,F are not Monophyletic' );
 116
 117 is($tree->is_monophyletic(-nodes => [$b,$a],
 118                           -outgroup => $f),1, 'A,B are Monophyletic');
 119
 120 # test for paraphyly
 121
 122 isnt(  $tree->is_paraphyletic(-nodes => [$a,$b,$c],
 123                            -outgroup => $d), 1,'A,B,C are not Monophyletic w D as outgroup');
 124
 125 is(  $tree->is_paraphyletic(-nodes => [$a,$f,$e],
 126                            -outgroup => $i), 1, 'A,F,E are monophyletic with I as outgroup');
 127
 128
 129 # test for rerooting the tree
 130 my $out = Bio::TreeIO->new(-format => 'newick',
 131                            -fh => \*STDERR,
 132                            -noclose => 1);
 133 $tree = $in->next_tree;
 134 $tree->verbose( -1 ) unless $verbose;
 135 my $node_cnt_orig = scalar($tree->get_nodes);
 136 # reroot on an internal node: should work fine
 137 $a = $tree->find_node('A');
 138 # removing node_count checks because re-rooting can change the
 139 # number of internal nodes (if it is done correctly)
 140 my $total_length_orig = $tree->total_branch_length;
 141 is $tree->total_branch_length, $tree->subtree_length,
 142     "subtree_length() without attributes is an alias to total_branch_lenght()";
 143 cmp_ok($total_length_orig, '>',$tree->subtree_length($a->ancestor),
 144        'Length of the tree is larger that lenght of a subtree');
 145 $out->write_tree($tree) if $verbose;
 146 is($tree->reroot($a),1, 'Can re-root with A as outgroup');
 147 $out->write_tree($tree) if $verbose;
 148 is($node_cnt_orig, scalar($tree->get_nodes), 'Count the number of nodes');
 149 my $total_length_new = $tree->total_branch_length;
 150 my $eps = 0.001 * $total_length_new;    # tolerance for checking length
 151 warn("orig total len ", $total_length_orig, "\n") if $verbose;
 152 warn("new  total len ", $tree->total_branch_length,"\n") if $verbose;
 153 # according to retree in phylip these branch lengths actually get larger
 154 # go figure...
 155 #ok(($total_length_orig >= $tree->total_branch_length - $eps)
 156 #   and ($total_length_orig <= $tree->total_branch_length + $eps));
 157 is($tree->get_root_node, $a->ancestor, "Root node is A's ancestor");
 158
 159 # try to reroot on an internal, will result in there being 1 less node
 160 $a = $tree->find_node('C')->ancestor;
 161 $out->write_tree($tree) if $verbose;
 162 is($tree->reroot($a),1, "Can reroot with C's ancsestor");
 163 $out->write_tree($tree) if $verbose;
 164 is($node_cnt_orig, scalar($tree->get_nodes), 'Check to see that node count is correct after an internal node was removed after this re-rooting');
 165 warn("orig total len ", $total_length_orig, "\n") if $verbose;
 166 warn("new  total len ", $tree->total_branch_length,"\n") if $verbose;
 167 cmp_ok($total_length_orig, '>=', $tree->total_branch_length - $eps,
 168        'Total original branch length is what it is supposed to be');
 169 cmp_ok($total_length_orig, '<=',$tree->total_branch_length + $eps,
 170        'Updated total branch length after the reroot');
 171 is($tree->get_root_node, $a->ancestor, 'Make sure root is really what we asked for');
 172
 173 # try to reroot on existing root: should fail
 174 $a = $tree->get_root_node;
 175 isnt( $tree->reroot($a),1, 'Testing for failed re-rerooting');
 176
 177 # try a more realistic tree
 178 $tree = $in->next_tree;
 179 $a = $tree->find_node('VV');
 180 $node_cnt_orig = scalar($tree->get_nodes);
 181 $total_length_orig = $tree->total_branch_length;
 182 $out->write_tree($tree) if $verbose;
 183 is($tree->reroot($a->ancestor),1, 'Test that rooting succeeded');
 184 $out->write_tree($tree) if $verbose;
 185 is($node_cnt_orig+1, scalar($tree->get_nodes), 'Test that re-rooted tree has proper number of nodes after re-rooting');
 186 $total_length_new = $tree->total_branch_length;
 187 $eps = 0.001 * $total_length_new;    # tolerance for checking length
 188 cmp_ok($total_length_orig, '>=', $tree->total_branch_length - $eps, 'Branch length before rerooting');
 189 cmp_ok($total_length_orig, '<=', $tree->total_branch_length + $eps,
 190        'Branch length after rerooting');
 191 is($tree->get_root_node, $a->ancestor->ancestor,'Root is really the ancestor we asked for');
 192
 193 # BFS and DFS search testing
 194 $treeio = Bio::TreeIO->new(-verbose => $verbose,
 195                              -format => 'newick',
 196                              -file   => test_input_file('test.nh'));
 197 $tree = $treeio->next_tree;
 198 my ($ct,$n) = (0);
 199 my $let = ord('A');
 200 for $n (  $tree->get_leaf_nodes ) {
 201     $n->id(chr($let++));
 202 }
 203
 204 for $n ( grep {! $_->is_Leaf } $tree->get_nodes ) {
 205     $n->id($ct++);
 206 }
 207 # enable for debugging
 208 Bio::TreeIO->new(-format => 'newick')->write_tree($tree) if( $verbose );
 209
 210 my $BFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'b')));
 211 is($BFSorder, '0,1,3,2,C,D,E,F,G,H,A,B', 'BFS traversal order');
 212 my $DFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'd')));
 213 is($DFSorder, '0,1,2,A,B,C,D,3,E,F,G,H', 'DFS travfersal order');
 214
 215
 216 # test some Bio::Tree::TreeFunctionI methods
 217 #find_node tested extensively already
 218 $tree->remove_Node('H');
 219 $DFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'd')));
 220 is($DFSorder, '0,1,2,A,B,C,D,3,E,F,G', 'DFS traversal after removing H');
 221 #get_lineage_nodes tested during get_lca
 222 $tree->splice(-remove_id => 'G');
 223 $DFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'd')));
 224 is($DFSorder, '0,1,2,A,B,C,D,3,E,F', 'DFS traversal after removing G');
 225 $tree->splice(-remove_id => [('E', 'F')], -keep_id => 'F');
 226 $DFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'd')));
 227 is($DFSorder, '0,1,2,A,B,C,D,F', 'DFS traversal after removing F');
 228 $tree->splice(-keep_id => [qw(0 1 2 A B C D)]);
 229 $DFSorder = join(",", map { $_->id } ( $tree->get_nodes(-order => 'd')));
 230 is($DFSorder, '0,1,2,A,B,C,D', 'DFS after removing all but 0,1,2,A,B,C,D');
 231 #get_lca, merge_lineage, contract_linear_paths tested in in Taxonomy.t
 232
 233
 234 # try out the id to bootstrap copy method
 235 $treeio = Bio::TreeIO->new(-format => 'newick',
 236                            -file   => test_input_file('bootstrap.tre'));
 237 $tree = $treeio->next_tree;
 238 my ($test_node) = $tree->find_node(-id => 'A');
 239 is($test_node->ancestor->id, 90,'Testing bootstrap copy');
 240 is($test_node->ancestor->ancestor->id, '25','Testing bootstrap copy');
 241 $tree->move_id_to_bootstrap;
 242 is($test_node->ancestor->id, '','Testing bootstrap copy');
 243 is($test_node->ancestor->bootstrap, '90', 'Testing bootstrap copy');
 244 is($test_node->ancestor->ancestor->id, '', 'Testing bootstrap copy');
 245 is($test_node->ancestor->ancestor->bootstrap, '25', 'Testing bootstrap copy');
 246
 247 # change TreeIO to parse
 248 $treeio = Bio::TreeIO->new(-format => 'newick',
 249                            -file   => test_input_file('bootstrap.tre'),
 250                            -internal_node_id => 'bootstrap');
 251 $tree = $treeio->next_tree;
 252 ($test_node) = $tree->find_node(-id => 'A');
 253 is($test_node->ancestor->id, '','Testing auto-boostrap copy during parse');
 254 is($test_node->ancestor->ancestor->id, '',
 255    'Testing auto-boostrap copy during parse');
 256 is($test_node->ancestor->bootstrap, '90',
 257    'Testing auto-boostrap copy during parse');
 258 is($test_node->ancestor->ancestor->bootstrap, '25',
 259    'Testing auto-boostrap copy during parse');
 260
 261
 262
 263 __DATA__
 264 (D,(C,(A,B)));
 265 (I,((D,(C,(A,B))),(E,(F,G))));
 266 (((A:0.3,B:2.1):0.45,C:0.7),D:4);
 267 (A:0.031162,((((((B:0.022910,C:0.002796):0.010713,(D:0.015277,E:0.020484):0.005336):0.005588,((F:0.013293,(G:0.018374,H:0.003108):0.005318):0.006047,I:0.014607):0.001677):0.004196,(((((J:0.003307,K:0.001523):0.011884,L:0.006960):0.006514,((M:0.001683,N:0.000100):0.002226,O:0.007085):0.014649):0.008004,P:0.037422):0.005201,(Q:0.000805,R:0.000100):0.015280):0.005736):0.004612,S:0.042283):0.017979,(T:0.006883,U:0.016655):0.040226):0.014239,((((((V:0.000726,W:0.000100):0.028490,((((X:0.011182,Y:0.001407):0.005293,Z:0.011175):0.004701,AA:0.007825):0.016256,BB:0.029618):0.008146):0.004279,CC:0.035012):0.060215,((((((DD:0.014933,(EE:0.008148,FF:0.000100):0.015458):0.003891,GG:0.010996):0.001489,(HH:0.000100,II:0.000100):0.054265):0.003253,JJ:0.019722):0.013796,((KK:0.001960,LL:0.004924):0.013034,MM:0.010071):0.043273):0.011912,(NN:0.031543,OO:0.018307):0.059182):0.026517):0.011087,((PP:0.000100,QQ:0.002916):0.067214,(RR:0.064486,SS:0.013444):0.011613):0.050846):0.015644,((TT:0.000100,UU:0.009287):0.072710,(VV:0.009242,WW:0.009690):0.035346):0.042993):0.060365);