A toolkit for working with phylogenetic data.
v0.24.0
tree/function/functions.cpp
Go to the documentation of this file.
1 /*
2  Genesis - A toolkit for working with phylogenetic data.
3  Copyright (C) 2014-2020 Lucas Czech and HITS gGmbH
4 
5  This program is free software: you can redistribute it and/or modify
6  it under the terms of the GNU General Public License as published by
7  the Free Software Foundation, either version 3 of the License, or
8  (at your option) any later version.
9 
10  This program is distributed in the hope that it will be useful,
11  but WITHOUT ANY WARRANTY; without even the implied warranty of
12  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  GNU General Public License for more details.
14 
15  You should have received a copy of the GNU General Public License
16  along with this program. If not, see <http://www.gnu.org/licenses/>.
17 
18  Contact:
19  Lucas Czech <lucas.czech@h-its.org>
20  Exelixis Lab, Heidelberg Institute for Theoretical Studies
21  Schloss-Wolfsbrunnenweg 35, D-69118 Heidelberg, Germany
22 */
23 
32 
37 #include "genesis/tree/tree.hpp"
39 
41 
42 #include <algorithm>
43 #include <cassert>
44 #include <functional>
45 #include <unordered_set>
46 #include <vector>
47 
48 #ifdef GENESIS_OPENMP
49 # include <omp.h>
50 #endif
51 
52 namespace genesis {
53 namespace tree {
54 
55 // =================================================================================================
56 // Node Properties
57 // =================================================================================================
58 
59 bool is_leaf( TreeLink const& link )
60 {
61  return &( link.next() ) == &link;
62 }
63 
64 bool is_leaf( TreeNode const& node )
65 {
66  return is_leaf( node.link() );
67 }
68 
69 bool is_leaf( TreeEdge const& edge )
70 {
71  return is_leaf( edge.secondary_link() );
72 }
73 
74 bool is_inner( TreeLink const& link )
75 {
76  return &( link.next() ) != &link;
77 }
78 
79 bool is_inner( TreeNode const& node )
80 {
81  return is_inner( node.link() );
82 }
83 
84 bool is_inner( TreeEdge const& edge )
85 {
86  return is_inner( edge.secondary_link() );
87 }
88 
89 bool is_root( TreeLink const& link )
90 {
91  return is_root( link.node() );
92 }
93 
94 bool is_root( TreeNode const& node )
95 {
96  // The link_ is always the one pointing towards the root. Also, the edge of that link always has
97  // the primary link set to that it points towards the root.
98  // At the root itself, however, this means we are pointing to ourselves. Use this to check
99  // if this node is the root.
100  return &( node.link().edge().primary_link() ) == &( node.link() );
101 }
102 
103 size_t degree( TreeLink const& link )
104 {
105  return degree( link.node() );
106 }
107 
108 size_t degree( TreeNode const& node )
109 {
110  size_t dgr = 0;
111  TreeLink const* lnk = &node.link();
112 
113  do {
114  ++dgr;
115  lnk = &lnk->next();
116  } while( lnk != &node.link() );
117 
118  return dgr;
119 }
120 
121 // =================================================================================================
122 // Node Count Properties
123 // =================================================================================================
124 
125 size_t max_degree( Tree const& tree )
126 {
127  size_t max = 0;
128  for( size_t i = 0; i < tree.node_count(); ++i ) {
129  max = std::max( max, degree( tree.node_at(i) ) );
130  }
131  return max;
132 }
133 
134 bool is_bifurcating( Tree const& tree, bool loose )
135 {
136  // Iterate all nodes and verify their degree.
137  for( size_t i = 0; i < tree.node_count(); ++i ) {
138  auto const deg = degree( tree.node_at(i) );
139 
140  // Any degree > 3 is multifurcating.
141  if( deg > 3 ) {
142  return false;
143  }
144 
145  // A degree of 2 is okay of we are loose, and always okay for the root.
146  if( deg == 2 ) {
147  auto const isroot = ( tree.node_at(i).index() == tree.root_node().index() );
148  assert( isroot == is_root( tree.node_at(i) ));
149 
150  if( ! isroot && ! loose ) {
151  return false;
152  }
153  }
154  }
155  return true;
156 }
157 
158 bool is_binary( Tree const& tree, bool loose )
159 {
160  return is_bifurcating( tree, loose );
161 }
162 
163 bool is_rooted( Tree const& tree )
164 {
165  return degree( tree.root_node() ) == 2;
166 }
167 
168 size_t leaf_node_count( Tree const& tree )
169 {
170  size_t sum = 0;
171  for (size_t i = 0; i < tree.node_count(); ++i) {
172  auto const& n = tree.node_at(i);
173  if( is_leaf(n) ) {
174  ++sum;
175  }
176  }
177  return sum;
178 }
179 
180 size_t inner_node_count( Tree const& tree )
181 {
182  return tree.node_count() - leaf_node_count( tree );
183 }
184 
185 size_t node_count( Tree const& tree )
186 {
187  return tree.node_count();
188 }
189 
190 size_t leaf_edge_count( Tree const& tree )
191 {
192  size_t sum = 0;
193  for( auto const& edge : tree.edges() ) {
194  if( is_leaf( edge.primary_node() ) || is_leaf( edge.secondary_node() ) ) {
195  ++sum;
196  }
197  }
198  return sum;
199 }
200 
201 size_t inner_edge_count( Tree const& tree )
202 {
203  size_t sum = 0;
204  for( auto const& edge : tree.edges() ) {
205  if( is_inner( edge.primary_node() ) && is_inner( edge.secondary_node() ) ) {
206  ++sum;
207  }
208  }
209  return sum;
210 }
211 
212 size_t edge_count( Tree const& tree )
213 {
214  return tree.edge_count();
215 }
216 
217 std::vector<size_t> inner_edge_indices( Tree const& tree )
218 {
219  std::vector<size_t> result;
220  for( auto const& edge : tree.edges() ) {
221  if( is_inner( edge.secondary_node() ) ) {
222  result.push_back( edge.index() );
223  }
224  }
225  return result;
226 }
227 
228 std::vector<size_t> leaf_edge_indices( Tree const& tree )
229 {
230  std::vector<size_t> result;
231  for( auto const& edge : tree.edges() ) {
232  if( is_leaf( edge.secondary_node() ) ) {
233  result.push_back( edge.index() );
234  }
235  }
236  return result;
237 }
238 
239 std::vector<size_t> inner_node_indices( Tree const& tree )
240 {
241  std::vector<size_t> result;
242  for( auto const& node : tree.nodes() ) {
243  if( is_inner( node )) {
244  result.push_back( node.index() );
245  }
246  }
247  return result;
248 }
249 
250 std::vector<size_t> leaf_node_indices( Tree const& tree )
251 {
252  std::vector<size_t> result;
253  for( auto const& node : tree.nodes() ) {
254  if( is_leaf( node )) {
255  result.push_back( node.index() );
256  }
257  }
258  return result;
259 }
260 
261 // =================================================================================================
262 // Tree Sides
263 // =================================================================================================
264 
266 {
267  // Get a quadratic matrix: for each edge, it gives a value whether each other edge is
268  // proximal (1) or distal (-1) relative to itself (0).
269  auto result = utils::Matrix<signed char>( tree.edge_count(), tree.edge_count(), 0 );
270 
271  // Helper function that traverses the subtree starting at a link,
272  // and for each edge in the subtree, sets its entry in the matrix to the given sign.
273  auto traverse = [&]( TreeLink const& start_link, size_t i, signed char sign ){
274  auto link = &( start_link.next() );
275  while( link != &start_link ) {
276  result( i, link->edge().index() ) = sign;
277  link = &( link->outer().next() );
278  }
279  };
280 
281  // For each edge, do the traversal in both directions and set the signs.
282  // This can probably be done more efficiently with only one smart traversal of the whole tree,
283  // but this function is not needed often enough right now.
284  for( size_t i = 0; i < tree.edge_count(); ++i ) {
285  traverse( tree.edge_at(i).primary_link(), i, -1 );
286  traverse( tree.edge_at(i).secondary_link(), i, 1 );
287  }
288 
289  return result;
290 }
291 
293 {
294  auto mat = utils::Matrix<signed char>( tree.node_count(), tree.node_count(), 0 );
295 
296  // Fill every row of the matrix.
297  #pragma omp parallel for
298  for( size_t i = 0; i < tree.node_count(); ++i ) {
299  auto const& row_node = tree.node_at(i);
300  auto const row_index = row_node.index();
301  auto const primary_link = &row_node.primary_link();
302 
303  // Fill root side subtree with 1s.
304  // We do set the value of inner nodes multiple times, but that's no problem.
305  // Also, we need to do an extra check for the root here, in order to set all
306  // subtrees of the root to -1.
307  signed char const value = is_root( row_node ) ? -1 : 1;
308  auto current_link = &( primary_link->outer() );
309  while( current_link != primary_link ) {
310  mat( row_index, current_link->node().index() ) = value;
311  current_link = &( current_link->next().outer() );
312  }
313 
314  // Fill all non-root side subtrees with -1s.
315  // We explicitly go through all non-root links of the node, in order to be really clear
316  // about our intentions. It would also work to simply follow the link chain until
317  // we reach the primary link again. However, this would also set -1 to our row node,
318  // thus we'd have to reset it, making the algorithm a bit messy.
319  // So, to be clear and clean, we avoid this.
320  auto sub_link = &( primary_link->next() );
321  while( sub_link != primary_link ) {
322 
323  // Now, for a given non-root subtree, set everything to -1.
324  current_link = &( sub_link->outer() );
325  while( current_link != sub_link ) {
326  mat( row_index, current_link->node().index() ) = -1;
327  current_link = &( current_link->next().outer() );
328  }
329 
330  // Go to next subtree.
331  sub_link = &( sub_link->next() );
332  }
333 
334  // Check that the diagonal element is untouched.
335  assert( mat( row_index, row_index ) == 0 );
336  }
337 
338  return mat;
339 }
340 
341 // =================================================================================================
342 // Subtrees
343 // =================================================================================================
344 
345 size_t subtree_size( Tree const& tree, TreeLink const& link )
346 {
347  if( ! belongs_to( tree, link )) {
348  throw std::runtime_error(
349  "Cannot caluclate subtree_size, as the given Link does not belong to the Tree."
350  );
351  }
352 
353  // TODO This is a quick and dirty solution. Traverse the whole subtree, add all nodes to a set
354  // and simply return the size of that set.
355 
356  std::unordered_set< TreeNode const* > visited_nodes;
357 
358  auto cur_link = &link.outer();
359  while( cur_link != &link ) {
360  visited_nodes.insert( &cur_link->node() );
361  cur_link = &cur_link->next().outer();
362  }
363 
364  return visited_nodes.size();
365 }
366 
367 std::vector<size_t> subtree_sizes( Tree const& tree, TreeNode const& node )
368 {
369  if( ! belongs_to( tree, node )) {
370  throw std::runtime_error(
371  "Cannot calculate subtree_sizes(), as the given Node does not belong to the Tree."
372  );
373  }
374 
375  // TODO this is an overly complex and thus error prone solution, maybe there is a better way?!
376 
377  // Prepare result vector.
378  std::vector<size_t> result;
379  result.resize( tree.node_count(), 0 );
380 
381  // We use a stack to track the subtree sizes.
382  // We store the entry link of the preorder traversal of the nodes. The entry link is the one
383  // that is given when visiting the node first while doing a eulertour traversal of the tree.
384  // This is always the next() link after the towards-the-starting-node/root link.
385  std::vector< TreeLink const* > stack;
386  stack.push_back( &node.link() );
387 
388  // Traverse the tree.
389  for( auto it : eulertour( node )) {
390 
391  // If this is the last time we visit that node on our way back up the tree.
392  // (The second part of the condition checks whether it is the starting node, because
393  // in this case, we do not want to remove it.)
394  if( &it.link().next() == stack.back() && stack.back() != &node.link() ) {
395 
396  // We finished with a subtree. Add the cummulative number of children of that subtree
397  // to the parent node, and remove the parent from the stack (as we are done with it).
398  auto st_size = result[ stack.back()->node().index() ];
399  stack.pop_back();
400  result[ stack.back()->node().index() ] += st_size;
401 
402  // If this node is already the current top stack element.
403  } else if( &it.node() == &stack.back()->node() ) {
404 
405  // Do nothing.
406 
407  // If it is a leaf.
408  } else if( is_leaf( it.link() )) {
409 
410  // Simply increment its parent's counter.
411  ++result[ stack.back()->node().index() ];
412 
413  // If we will visit that node in the future again.
414  } else {
415 
416  // Add a count for the immediate child (i.e., the current node) to the current stack
417  // end (i.e., increment the counter of children of that node),
418  // then add the current node itself to the stack, so that in the next iteration,
419  // we will increase its counts.
420  ++result[ stack.back()->node().index() ];
421  stack.push_back( &it.link() );
422  }
423  }
424 
425  // The stack now should contain only a single node, which is the starting node itself.
426  assert( stack.size() == 1 && stack.back() == &node.link() );
427 
428  // The size of the subtree of the starting node is always the number of nodes in the tree
429  // minus one for that node itself (as it is not counted as part of its subtree).
430  assert( result[ node.index() ] == tree.node_count() - 1 );
431 
432  return result;
433 }
434 
435 std::vector<size_t> subtree_sizes( Tree const& tree )
436 {
437  return subtree_sizes( tree, tree.root_node() );
438 }
439 
440 size_t subtree_max_path_height( Tree const& tree, TreeLink const& link )
441 {
442  if( ! belongs_to( tree, link )) {
443  throw std::runtime_error(
444  "Cannot calculate subtree_max_path_height(), "
445  "as the given Link does not belong to the Tree."
446  );
447  }
448 
449  // TODO make more efficient. no need for full dist vector.
450  auto dists = node_path_length_vector( tree, link.outer().node() );
451  size_t max = 0;
452 
453  auto cur_link = &link.outer();
454  while( cur_link != &link ) {
455  max = std::max( max, dists[ cur_link->node().index() ] );
456  cur_link = &cur_link->next().outer();
457  }
458  return max;
459 }
460 
461 std::vector<size_t> subtree_max_path_heights( Tree const& tree, TreeNode const& node )
462 {
463  if( ! belongs_to( tree, node )) {
464  throw std::runtime_error(
465  "Cannot calculate subtree_max_path_heights(), "
466  "as the given Node does not belong to the Tree."
467  );
468  }
469 
470  auto result = std::vector<size_t>( tree.node_count(), 0 );
471 
472  // Recursive helper function that evaluates the wanted size for a given subtree,
473  // stores the result in the vector and returns it for recursive usage.
474  std::function< size_t( TreeLink const* )> rec_subtree_height = [&]( TreeLink const* l ) {
475  size_t link_max = 0;
476  TreeLink const* cl = &l->next();
477  while( cl != l ) {
478  link_max = std::max( link_max, 1 + rec_subtree_height( &cl->outer() ));
479  cl = &cl->next();
480  }
481 
482  result[ l->node().index() ] = link_max;
483  return link_max;
484  };
485 
486  // Loop all subtrees of the given node and find the highest.
487  // This loop is a bit different from the one in the recursive function, as we need to evaluate
488  // all links of the given starting node, instead of just the ones away from the start node.
489  size_t node_max = 0;
490  TreeLink const* cur_l = &node.link();
491  do {
492  node_max = std::max( node_max, 1 + rec_subtree_height( &cur_l->outer() ));
493  cur_l = &cur_l->next();
494  } while( cur_l != &node.link() );
495  result[ node.index() ] = node_max;
496 
497  return result;
498 }
499 
500 std::vector<size_t> subtree_max_path_heights( Tree const& tree )
501 {
502  return subtree_max_path_heights( tree, tree.root_node() );
503 }
504 
505 utils::Matrix<signed char> sign_matrix( Tree const& tree, bool compressed )
506 {
507  // Edge cases and input checks.
508  if( tree.empty() ) {
510  }
511  if( ! is_rooted( tree )) {
512  throw std::invalid_argument( "Tree is not rooted. Cannot calculate its sign matrix." );
513  }
514  if( ! is_bifurcating( tree )) {
515  throw std::invalid_argument( "Tree is not bifurcating. Cannot calculate its sign matrix." );
516  }
517 
518  // Prepare a result matrix of the full size. For the compressed version,
519  // we later replate it again.
520  auto result = utils::Matrix<signed char>( tree.node_count(), tree.node_count(), 0 );
521 
522  // Helper function that fills all columns of a subtree with a given sign.
523  auto fill_subtree_indices = [&]( size_t row_idx, Subtree const& st, signed char sign ){
524  for( auto const& it : preorder(st) ) {
525  result( row_idx, it.node().index() ) = sign;
526  }
527  };
528 
529  // Fill every row of the matrix.
530  #pragma omp parallel for
531  for( size_t i = 0; i < tree.node_count(); ++i ) {
532  auto const& row_node = tree.node_at(i);
533  auto const row_idx = row_node.index();
534 
535  if( row_idx == tree.root_node().index() ) {
536 
537  // The root node is special: we use its two subtrees directly.
538  assert( &row_node.link().next().next() == &row_node.link() );
539  fill_subtree_indices( row_idx, Subtree{ row_node.link().outer() }, +1 );
540  fill_subtree_indices( row_idx, Subtree{ row_node.link().next().outer() }, -1 );
541 
542  } else if( is_inner( row_node )) {
543 
544  // All other inner nodes are filled using their subtrees.
545  assert( &row_node.link().next().next().next() == &row_node.link() );
546  fill_subtree_indices( row_idx, Subtree{ row_node.link().next().outer() }, +1 );
547  fill_subtree_indices( row_idx, Subtree{ row_node.link().next().next().outer() }, -1 );
548  }
549  }
550 
551  // For the compressed version, we re-use the previous result matrix,
552  // and simply fill a new one with the needed rows and columns.
553  // The data is not too big, and this is way easier and cleaner to implement.
554  if( compressed ) {
555  // Create a matrix with rows for each inner node and columns for each tip node.
556  auto const in_node_idcs = inner_node_indices( tree );
557  auto const lf_node_idcs = leaf_node_indices( tree );
558  auto result_cmpr = utils::Matrix<signed char>( in_node_idcs.size(), lf_node_idcs.size(), 0 );
559 
560  // Fill the matrix at the indices that belong to inner nodes (for rows) and
561  // leaf nodes (for columns).
562  for( size_t r = 0; r < in_node_idcs.size(); ++r ) {
563  for( size_t c = 0; c < lf_node_idcs.size(); ++c ) {
564  result_cmpr( r, c ) = result( in_node_idcs[r], lf_node_idcs[c] );
565  }
566  }
567 
568  // Replace the result matrix efficiently.
569  using std::swap;
570  swap( result, result_cmpr );
571  }
572 
573  return result;
574 }
575 
576 // =================================================================================================
577 // Misc
578 // =================================================================================================
579 
580 std::vector< TreeLink const* > path_to_root( TreeNode const& node )
581 {
582  std::vector< TreeLink const* > path;
583 
584  // Move towards the root and record all links in between.
585  TreeLink const* cur_link = &node.primary_link();
586  while( &cur_link->edge().secondary_link() == cur_link ) {
587 
588  // The above while condition means: is it the root?! Assert, that the default way of
589  // checking for the root by using the node gives the same result.
590  assert( ! is_root( cur_link->node() ));
591 
592  // Assert that the primary direction is correct.
593  assert( cur_link == &cur_link->edge().secondary_link() );
594 
595  // Add the primary link of the current node to the list.
596  path.push_back( cur_link );
597 
598  // Move one node towards the root.
599  // Assert that the default way of finding the next node towards the root (by using
600  // the edge) gives the same result as simply using the link's outer node.
601  // This is the case because the cur link is the one that points towards the root
602  // (which was asserted above).
603  assert( &cur_link->edge().primary_link() == &cur_link->outer() );
604  cur_link = &cur_link->outer().node().primary_link();
605  }
606 
607  // Now finally add the root itself and return the list.
608  assert( is_root( cur_link->node() ));
609  path.push_back( cur_link );
610  return path;
611 }
612 
613 TreeNode const& lowest_common_ancestor( TreeNode const& node_a, TreeNode const& node_b )
614 {
615  // Speedup and simplification.
616  if( &node_a == &node_b ) {
617  return node_a;
618  }
619 
620  auto path_a = path_to_root( node_a );
621  auto path_b = path_to_root( node_b );
622 
623  // We must have at least the two original links in the front and the root in the back.
624  assert( path_a.size() > 0 && path_b.size() > 0 );
625  assert( path_a.front() == &node_a.link() );
626  assert( path_b.front() == &node_b.link() );
627  assert( path_a.back() == path_b.back() );
628 
629  // Remove from back as long as the last two elements are the same.
630  // At the end of this, the remaining links are the ones on the path between
631  // the two original links.
632  while(
633  path_a.size() > 1 &&
634  path_b.size() > 1 &&
635  path_a.at( path_a.size() - 1 ) == path_b.at( path_b.size() - 1 ) &&
636  path_a.at( path_a.size() - 2 ) == path_b.at( path_b.size() - 2 )
637  ) {
638  path_a.pop_back();
639  path_b.pop_back();
640  }
641 
642  // Now, the last elements need to be the same (the LCA of the start and finish node).
643  assert( path_a.size() > 0 && path_b.size() > 0 );
644  assert( path_a.back() == path_b.back() );
645 
646  return path_a.back()->node();
647 }
648 
650 {
651  auto const& c_node_a = static_cast< TreeNode const& >( node_a );
652  auto const& c_node_b = static_cast< TreeNode const& >( node_b );
653  return const_cast< TreeNode& >( lowest_common_ancestor( c_node_a, c_node_b ));
654 }
655 
657 {
658  auto res = utils::Matrix<size_t>( tree.node_count(), tree.node_count() );
659 
660  // This is not the best way to calculate all pairwise LCAs.
661  // In the Quartet Scores code, we use range minimum queries and eulertours to achive the
662  // same result in less time. But for now, this code is good enough.
663 
664  // Parallel specialized code.
665  #ifdef GENESIS_OPENMP
666 
667  // We only need to calculate the upper triangle. Get the number of indices needed
668  // to describe this triangle.
669  size_t const max_k = utils::triangular_size( tree.node_count() );
670 
671  #pragma omp parallel for
672  for( size_t k = 0; k < max_k; ++k ) {
673 
674  // For the given linear index, get the actual position in the Matrix.
675  auto const rc = utils::triangular_indices( k, tree.node_count() );
676  auto const r = rc.first;
677  auto const c = rc.second;
678 
679  auto const& lca = lowest_common_ancestor( tree.node_at(r), tree.node_at(c) );
680  res( r, c ) = lca.index();
681  res( c, r ) = lca.index();
682  }
683 
684  // Lastly, because the trinangular indices exluce the diagonale, we need to fill this
685  // by hand. Luckily, those are always the indices themselves, as the LCA of a node
686  // and itself is again itself.
687  #pragma omp parallel for
688  for( size_t d = 0; d < tree.node_count(); ++d ) {
689  res( d, d ) = d;
690  }
691 
692  // If no threads are available at all, use serial version.
693  #else
694 
695  for( size_t r = 0; r < tree.node_count(); ++r ) {
696 
697  // The result is symmetric - we only calculate the upper triangle.
698  for( size_t c = r; c < tree.node_count(); ++c ) {
699 
700  auto const& lca = lowest_common_ancestor( tree.node_at(r), tree.node_at(c) );
701  res( r, c ) = lca.index();
702  res( c, r ) = lca.index();
703  }
704 
705  // See above: the diagonale contains its indices.
706  assert( res( r, r ) == r );
707  }
708 
709  #endif
710 
711  return res;
712 }
713 
714 } // namespace tree
715 } // namespace genesis
std::vector< size_t > leaf_node_indices(Tree const &tree)
Get a list of the node indices of all leaf TreeNodes.
utils::Matrix< signed char > sign_matrix(Tree const &tree, bool compressed)
Compute the sign matrix or Sequential Binary Partition (SBP) of a Tree.
TreeLink & secondary_link()
Return the TreeLink of this TreeEdge that points away from the root.
Definition: edge.hpp:130
bool is_rooted(Tree const &tree)
Return whether the Tree is rooted, that is, whether the root node has two neighbors.
bool is_root(TreeLink const &link)
Return whether the link belongs to the root node of its Tree.
bool is_inner(TreeLink const &link)
Return true iff the node of the given link is an inner node.
void swap(SequenceSet &lhs, SequenceSet &rhs)
std::vector< size_t > node_path_length_vector(Tree const &tree, TreeNode const &node)
Return a vector containing the depth of all nodes with respect to the given start node...
Tree operator functions.
bool empty() const
Return whether the Tree is empty (i.e., has no nodes, edges and links).
Definition: tree/tree.hpp:194
TreeLink & primary_link()
Return the TreeLink that points towards the root.
Definition: node.hpp:110
size_t inner_edge_count(Tree const &tree)
Return the number of Edges of a Tree that do not lead to a leaf Node.
size_t inner_node_count(Tree const &tree)
Count the number of inner Nodes.
utils::Matrix< signed char > node_root_direction_matrix(Tree const &tree)
Calculate a Matrix that indicates the nodes on the root side of a given node.
utils::Matrix< signed char > edge_sides(Tree const &tree)
Create a Matrix that indiciaces the relative position of the Edges of a Tree, i.e., whether they are on the root side or non-root side.
size_t edge_count() const
Return the number of TreeEdges of the Tree.
Definition: tree/tree.hpp:272
size_t node_count(Tree const &tree)
Return the number of Nodes of a Tree. Same as Tree::node_count().
bool is_bifurcating(Tree const &tree, bool loose)
Return whether the Tree is bifurcating.
double sum(const Histogram &h)
Container namespace for all symbols of genesis in order to keep them separate when used as a library...
size_t subtree_max_path_height(Tree const &tree, TreeLink const &link)
Calculate the height of a subtree, that is, the maximum path length to a leaf of that subtree...
Reference to a subtree of a Tree.
Definition: subtree.hpp:69
utils::Range< IteratorEulertour< true > > eulertour(ElementType const &element)
Definition: eulertour.hpp:206
size_t subtree_size(Tree const &tree, TreeLink const &link)
Return the size of the subtree defined by the given TreeLink, measured in number of nodes...
Matrix operators.
std::vector< TreeLink const *> path_to_root(TreeNode const &node)
Helper function that finds all TreeLinks between a given TreeNode and the root of the Tree...
size_t edge_count(Tree const &tree)
Return the number of Edges of a Tree. Same as Tree::edge_count().
utils::Range< IteratorPath< true > > path(ElementType const &start, ElementType const &finish)
Definition: path.hpp:337
size_t leaf_node_count(Tree const &tree)
Count the number of leaf Nodes of a Tree.
utils::Range< IteratorNodes > nodes()
Definition: tree/tree.hpp:418
Class for representing phylogenetic trees.
Definition: tree/tree.hpp:97
TreeLink & link()
Return the TreeLink that points towards the root.
Definition: node.hpp:129
size_t triangular_size(size_t n)
Calculate the number of linear indices needed for a triangular Matrix of size n.
utils::Range< IteratorEdges > edges()
Definition: tree/tree.hpp:452
TreeNode & node_at(size_t index)
Return the TreeNode at a certain index.
Definition: tree/tree.hpp:220
TreeNode const & lowest_common_ancestor(TreeNode const &node_a, TreeNode const &node_b)
Return the lowest common ancestor of two TreeNodes.
std::vector< size_t > subtree_max_path_heights(Tree const &tree, TreeNode const &node)
std::vector< size_t > leaf_edge_indices(Tree const &tree)
Get a list of the edge indices of all leaf edges, that is, all TreeEdges that lead to a leaf TreeNode...
TreeEdge & edge_at(size_t index)
Return the TreeEdge at a certain index.
Definition: tree/tree.hpp:238
bool belongs_to(Tree const &tree, TreeNode const &node)
Return whether the TreeNode belongs to the Tree, i.e., whether it is owned by the Tree...
TreeLink const & link() const
Get the TreeLink that separates the subtree from the rest of the tree.
Definition: subtree.hpp:125
utils::Matrix< size_t > lowest_common_ancestors(Tree const &tree)
Return the lowest common ancestor of each pair of TreeNodes for a given tree, in form of a Matrix of ...
size_t node_count() const
Return the number of TreeNodes of the Tree.
Definition: tree/tree.hpp:264
utils::Range< IteratorPreorder< true > > preorder(ElementType const &element)
std::pair< size_t, size_t > triangular_indices(size_t k, size_t n)
Given a linear index in a upper triangular Matrix, find the corresponding Matrix indices.
TreeNode & root_node()
Return the TreeNode at the current root of the Tree.
Definition: tree/tree.hpp:300
Header of Tree class.
bool is_binary(Tree const &tree, bool loose)
Alias for is_bifurcating().
Header of Tree distance methods.
size_t index() const
Return the index of this Node.
Definition: node.hpp:102
TreeLink & primary_link()
Return the TreeLink of this TreeEdge that points towards the root.
Definition: edge.hpp:114
size_t degree(TreeLink const &link)
Return the degree of the node for a given TreeLink, i.e. how many neighbouring nodes it has...
std::vector< size_t > inner_node_indices(Tree const &tree)
Get a list of the node indices of all inner TreeNodes.
size_t max_degree(Tree const &tree)
Return the highest degree of the Nodes of a Tree.
std::vector< size_t > inner_edge_indices(Tree const &tree)
Get a list of the edge indices of all inner edges, that is, all TreeEdges that do not lead to a leaf ...
bool is_leaf(TreeLink const &link)
Return true iff the node of the given link is a leaf node.
std::vector< size_t > subtree_sizes(Tree const &tree, TreeNode const &node)
Calculate the sizes of all subtrees as seen from the given TreeNode.
size_t leaf_edge_count(Tree const &tree)
Return the number of Edges of a Tree that lead to a leaf Node.