Trees and Phylogenetics (skbio.tree
)#
This module provides functionality for working with trees, including phylogenetic trees and various types of hierarchies with relevance in biology or beyond. It supports trees that are multifurcating and nodes that have single descendants. Functionality is provided for constructing trees, for traversing in multiple ways, comparisons, fetching subtrees, and more.
See the Tutorial section for working with trees using scikit-bio.
Tree structure and operations#
TreeNode
is the data structure for tree representation. A tree consists of an
arbitrary number of interconnected TreeNode
objects representing individual nodes.
Each node has pointers to its parent and children. The class provides a large number of
methods to achieve various tasks for tree analysis and manipulation.
|
Represent a node within a tree. |
Tree Construction#
Algorithms for phylogenetic reconstruction based on distance matrices.
|
Perform unweighted pair group method with arithmetic mean (UPGMA) or its weighted variant (WPGMA) for phylogenetic reconstruction. |
|
Perform neighbor joining (NJ) for phylogenetic reconstruction. |
|
Perform greedy minimum evolution (GME) for phylogenetic reconstruction. |
|
Perform balanced minimum evolution (BME) for phylogenetic reconstruction. |
Tree Rearrangement#
Algorithms for improving an existing phylogenetic tree via rearrangement.
|
Perform nearest neighbor interchange (NNI) to improve a phylogenetic tree. |
Tree Comparison#
Metrics for assessing the topological or length dissimilarity among trees.
|
Calculate Robinson-Foulds (RF) distances among trees. |
|
Calculate weighted Robinson-Foulds (wRF) distances or variants among trees. |
|
Calculate path-length distances or variants among trees. |
Tree utilities#
|
Determine consensus trees from a list of rooted trees. |
Exceptions#
|
General tree error. |
|
Missing length when expected. |
|
Duplicate nodes with identical names. |
|
Expecting a node. |
|
Missing a parent. |
Tutorial#
In scikit-bio, the TreeNode
class is the main and only data structure for tree
represention. It is simple, flexible and powerful.
>>> from skbio import TreeNode
Loading a tree#
A tree can be constructed from a Newick string. Newick is a common file format used to represent tree structure based on nesting with parentheses. For instance, the following string describes a 3-taxon tree, with one internal node:
((A,B)C,D)root;
Where A, B, and D are tips of the tree, and C is an internal node that unites tips A and B. C and D are in turn united by the root of the tree.
If this string is stored in a file input.nwk
, we can read it into scikit-bio with:
tree = TreeNode.read("input.nwk")
Alternatively, we can directly convert a Newick string into a tree with:
>>> tree = TreeNode.read(["((A,B)C,D)root;"])
>>> tree
<TreeNode, name: root, internal node count: 1, tips count: 3>
Let’s display an ASCII representation of the tree:
>>> print(tree.ascii_art())
/-A
/C-------|
-root----| \-B
|
\-D
Optionally, we can attach branch lengths to the nodes in a Newick string. In phylogenetics, the branch length usually represents the amount of evolution accumulated from one node to another. For example:
>>> tree = TreeNode.read(["((A:2,B:3)C:1,D:4)root;"])
To save a tree to a Newick file:
tree.write("output.nwk")
Traversing a tree#
The process of visiting each and every node within a tree is called tree
traversal. There are a few common ways to traverse a tree, and
depending on your use, some methods are more appropriate than others. The wiki page
linked above will go into further depth on this topic. Here, we’re only going to cover
two of the commonly used traversal methods: preorder
and
postorder
.
A preorder traversal moves from root to tips, looking at the left most child first. For instance:
>>> for node in tree.preorder():
... print(node.name)
root
C
A
B
D
A postorder traversal moves from tips to root, accessing the left subtree tips first before walking back up the tree:
>>> for node in tree.postorder():
... print(node.name)
A
B
C
D
root
If you are not certain (or don’t mind) what method should be used to traverse the tree,
simply use traverse
.
Note, by default, those traversal methods will include self, i.e., the root node from
which traversal starts. One can pass the include_self=False
flag to the method call
if you wish to exclude it.
TreeNode
provides multiple methods for iterating over certain nodes. For examples,
tips
and non_tips
iterate over tips and internal
nodes, respectively:
>>> for node in tree.tips():
... print(node.name)
A
B
D
>>> for node in tree.non_tips():
... print(node.name)
C
If we just want to know the total number of nodes (or just tips) in a tree, we can do:
>>> tree.count()
5
>>> tree.count(tips=True)
3
Instead, if we want to obtain a set of taxa (tip names) from a tree, we can use the
subset
method (it is called the “subset” because one can call it at
any node and get a subset of taxa descending from it).
>>> taxa = tree.subset()
>>> sorted(taxa)
['A', 'B', 'D']
Editing a tree#
As discussed above, a tree consists of multiple TreeNode
objects connected to each
other. To manipulate a tree, one just needs to change the connections between nodes.
TreeNode
provides several standard methods that let you work with nodes like Python
lists.
Create a new tip “E” and attach it to the internal node “C” as a child of it:
>>> node = TreeNode("E")
>>> tree.find("C").append(node)
>>> print(tree.ascii_art())
/-A
|
/C-------|--B
| |
-root----| \-E
|
\-D
Create two tips “F” and “G” and attach both of them to tip “D”. Now “D” becomes an internal node because it has children:
>>> D = tree.find("D")
>>> D.extend([TreeNode("F"), TreeNode("G")])
>>> print(tree.ascii_art())
/-A
|
/C-------|--B
| |
-root----| \-E
|
| /-F
\D-------|
\-G
Prune clade “D” and graft it to node “B”. Once a node is attached to a new parent, the connection to its old parent is automatically removed.
>>> B = tree.find("B")
>>> B.append(D)
>>> print(tree.ascii_art())
/-A
|
| /-F
-root---- /C-------|-B------- /D-------|
| \-G
|
\-E
Remove clade “D” from node “B”. This operation detaches the two nodes, but the removed node still exists in the memory and once can still use it.
>>> B.remove(D)
True
>>> print(tree.ascii_art())
/-A
|
-root---- /C-------|--B
|
\-E
>>> print(D.ascii_art())
/-F
-D-------|
\-G
These methods automatically clear lookup tables and other caches because they may be
obsolete for the modified tree. This safety measure could impact performance when
working with very large trees. If caches are not present or relevant, we can add
uncache=False
to these methods to skip this check. See has_caches
for details.
TreeNode
provides several advanced methods for editing a tree as a batch. For
examples, bifurcate
converts an arbitrary tree into a bifurcating
tree (binary tree), shear
refines a tree to a given set of taxa,
unpack_by_func
unpacks (contracts) branches that suffice certain
criteria (e.g., bootstrap support below a threshold).
Building a tree#
scikit-bio provides several common algorithms for building a tree structure based on
a distance matrix (DistanceMatrix
) that stores pairwise
distances between taxa.
In the field of phylogenetics, the evolutionary relationships
among organisms are inferred using genetic data.
Distance-based methods are a common category
of methods for such purpose. The distances are usually computed from sequence alignment
(see skbio.alignment
), but they can also be inferred using alternative methods
such as _k_-mer frequency (see kmer_frequencies
) and
community diversity (see skbio.diversity.beta
), and the applications of
distance-based tree building are not limited to evolutionary biology, but can be
extended to other fields (e.g., hierarchical clustering in data science).
To begin with, let’s create a distance matrix of seven taxa.
>>> from skbio import DistanceMatrix
>>> dm = DistanceMatrix([
... [0. , 0.005, 0.033, 0.387, 0.626, 0.635, 0.665],
... [0.005, 0. , 0.028, 0.387, 0.626, 0.635, 0.665],
... [0.033, 0.028, 0. , 0.396, 0.635, 0.64 , 0.679],
... [0.387, 0.387, 0.396, 0. , 0.611, 0.611, 0.66 ],
... [0.626, 0.626, 0.635, 0.611, 0. , 0.147, 0.758],
... [0.635, 0.635, 0.64 , 0.611, 0.147, 0. , 0.754],
... [0.665, 0.665, 0.679, 0.66 , 0.758, 0.754, 0. ],
... ], ids=[
... 'human', 'chimp', 'monkey', 'pig', 'mouse', 'rat', 'chicken'
... ])
The following code applies neighbor joining (nj
), a
classic distance-based tree-building algorithm:
>>> from skbio.tree import nj
>>> tree = nj(dm)
Display the reconstructed phylogenetic tree:
>>> print(tree.ascii_art())
/-human
|
|--chimp
---------|
| /-monkey
| |
\--------| /-pig
| |
\--------| /-mouse
| /--------|
\--------| \-rat
|
\-chicken
Other methods include upgma
(simplest), gme
(highly scalable) and
bme
. All of them can be applied in the same way.
These four methods create a tree from scratch through a fixed set of procedures.
Additionally, starting from an existing tree, tree rearrangement methods will explore the tree space further in order to
improve it. Nearest neighbor interchange (nni
) is one such method.
>>> from skbio.tree import nni
>>> tree = nni(tree, dm)
Rooting a tree#
The phylogenetic tree created above seems to properly reflect the relatedness among the seven organisms. However, the “root” of the tree doesn’t appear to be aligned with our knowledge of evolution. This is because the tree is in fact unrooted. That is, the placement of branches does not inform the direction of evolution.
In scikit-bio, the distinguishment between rooted and unrooted trees is implicit: Every tree has a root node. A tree is considered as “rooted” if its root node has exactly two children. In contrast, an “unrooted” tree may have three (the most common case), one, or more than three children attached to its root node.
There are several approaches to root an unrooted tree to reflect the correct direction
of evolution. Midpoint rooting (root_at_midpoint
) is a simple method
that places the root at the midpoint between the two most distant tips in the tree.
>>> rooted = tree.root_at_midpoint()
>>> print(rooted.ascii_art())
/-chicken
|
-root----| /-mouse
| /--------|
| | \-rat
\--------|
| /-pig
\--------|
| /-monkey
\--------|
| /-chimp
\--------|
\-human
If we already know which taxon or taxa (called an “outgroup”) were separated from all
other taxa in evolution, we can perform outgroup rooting
(root_by_outgroup
):
>>> rooted = tree.root_by_outgroup(["chicken"])
Or we can directly root a tree at a given taxon, node or branch:
>>> rooted = tree.root_at("chicken")
To convert a rooted tree into unrooted, we can apply the unroot
method.
Analyzing topology#
The topology of a tree describes the grouping pattern among taxa. For example, in the rooted tree displayed above, human, chimp and monkey form a clade (a monophyletic group), while mouse and rat form another clade. A clade can be referred to by its root node, which is the lowest common ancestor (LCA) of all descending taxa:
>>> primate = rooted.lca(["human", "chimp", "monkey"])
>>> rodent = rooted.lca(["mouse", "rat"])
Note that one doesn’t need to list all descending taxa in order to find the LCA. Two or more taxa that join at the LCA node on their paths to the root are sufficient for locating the LCA.
>>> primate = rooted.lca(["human", "monkey"])
A clade regardless of its internal structure can be simply defined by the
subset
of taxa descending from it (see above). All
subsets
within a rooted tree describes all possible taxon groups.
Assume we have another rooted tree tree2
of the same taxa. To assess whether
that tree also support human, chimp and monkey grouped together, we can simply do:
primate.subset() in tree2.subsets()
When working with unrooted trees, we should replace the notion of subset with
bipartition (bipart
), which defines the separation of all taxa into
two parts by a given branch, regardless of the direction of evolution. For example,
both the rooted and unrooted trees support the separation of primates from other taxa:
>>> primate.bipart() in tree.biparts()
True
The topology of two trees with shared taxa can be compared using the Robinson-Foulds (RF) distance, which is the number of bipartitions (unrooted trees) or subsets (rooted trees) that are distinct between the two trees. For example, let’s define three trees, the first two of which have the same topology despite different orders of clades.
>>> tree1 = TreeNode.read(["((A,B),(C,D),(E,F));"])
>>> tree2 = TreeNode.read(["((F,E),(C,D),(B,A));"])
>>> tree3 = TreeNode.read(["((C,B),(A,D),(E,F));"])
Calculate RF distances using (compare_rfd
):
>>> tree1.compare_rfd(tree1) # identity case
0.0
>>> tree1.compare_rfd(tree2) # same tree but different clade order
0.0
>>> tree1.compare_rfd(tree3) # only 1 of 3 common clade
4.0
We can calculate pairwise RF distances between any number of trees using a single call
of rf_dists
, which returns a distance matrix:
>>> from skbio.tree import rf_dists
>>> dm = rf_dists([tree1, tree2, tree3])
>>> print(dm)
3x3 distance matrix
IDs:
'0', '1', '2'
Data:
[[ 0. 0. 4.]
[ 0. 0. 4.]
[ 4. 4. 0.]]
Analyzing distances#
With branch lengths annotated, a phylogenetic tree describes the evolutionary distances
among taxa. This important property is the basis for multiple downstream analyses, such
as community diversity (see skbio.diversity
).
The distance
between two nodes in a tree (a.k.a., patristic distance)
is the total length of branches constituting the path connecting them. For example:
>>> human, mouse = tree.find("human"), tree.find("mouse")
>>> human.distance(mouse).round(3)
0.628
The maximum distance between any pair of taxa (maxdist
) and the total
branch length connecting all taxa (total_length
) can be used to
measure the biodiversity of groups of organisms.
The distances between all pairs of tips (taxa) in a tree can be efficiently calculated
using cophenet
, which returns a distance matrix.
dm = tree.cophenet()
The evolutionary distances reflected by two trees can be compared. First, let’s create two trees with the same taxa but different branch lengths.
>>> tree1 = TreeNode.read(["((A:0.1,B:0.2):0.3,C:0.4,D:0.5);"])
>>> tree2 = TreeNode.read(["((A:0.4,B:0.8):0.3,C:0.1,D:0.5);"])
Calculate the distance between the cophenetic distance matrics using
(compare_cophenet
):
>>> tree1.compare_cophenet(tree2).round(3)
0.12
Likewise, we can calculate this metric between multiple trees using path_dists
.