Mapping Phylogenetic Trees to Reveal Distinct Patterns of Evolution

Evolutionary relationships are frequently described by phylogenetic trees, but a central barrier in many fields is the difficulty of interpreting data containing conflicting phylogenetic signals. We present a metric-based method for comparing trees which extracts distinct alternative evolutionary relationships embedded in data. We demonstrate detection and resolution of phylogenetic uncertainty in a recent study of anole lizards, leading to alternate hypotheses about their evolutionary relationships. We use our approach to compare trees derived from different genes of Ebolavirus and find that the VP30 gene has a distinct phylogenetic signature composed of three alternatives that differ in the deep branching structure. Key words: phylogenetics, evolution, tree metrics, genetics, sequencing.


Introduction
Phylogenetic trees are fundamental tools for understanding evolution.Improvements in sequencing technology have meant that phylogenetic analyses are growing in size and scope.However, when a tree is inferred from data there are multiple sources of uncertainty.Competing approaches to tree estimation can produce markedly different trees.Trees may conflict due to signals from selection (e.g.convergent evolution), and/or when derived from different data (e.g. the organisms' mitochondrial vs nuclear DNA, individual genes or other subsets of sequence data [16]).Evolution is not always tree-like: species trees differ from gene trees, and many organisms exchange genes through horizontal gene transfer.It is therefore crucial to be able to compare trees to identify these signals.
Trees can be compared by direct visualization, aided by methods such as tanglegrams and software such as DensiTree [4], but this does not lend itself to detailed comparison of large groups of trees.Current quantitative methods for tree comparison suffer from the challenges of visualizing non-Euclidean distances [12] and from counter-intuitive behavior.For example, the nearest-neighbor interchange (NNI) distance of Robinson and Foulds (RF) [25], which is the most widely used, is hampered by the fact that large NNI distances do not imply large changes among the shared ancestry of most tips [31,20,17].In fact, two trees differing in the placement of a single tip can be a maximal NNI distance apart.
We introduce a metric which flexibly captures both tree structure and branch lengths.It can be used as a quantitative tool for comparing phylogenetic trees.Each metric on trees defines a tree space; this tree space lends itself to clear visualizations in low dimensions, and captures and highlights differences in trees according to their biological significance.
In Section 2 we formally define our distance function, prove that it is a metric, and explain its capacity to capture tree structure and branch lengths.We also provide a brief survey, explaining how our metric relates to and differs from existing metrics (Section 2.3).In Section 3 we explain some of the applications of our metric.We show how our metric enables visualization of tree space (Section 3.1) and detection of islands (Section 3.2), which we demonstrate with a simple application to Dengue fever phylogenies.We also explain how our metric provides a new suite of methods for selecting summary trees in Section 3.3.We conclude with some ideas for extensions to our metric in Section 4. 1 arXiv:1507.05211v3[q-bio.PE] 11 Sep 2015 2 Metrics

Our metric: definition and proof
Let T k be the set of all rooted trees on k tips with labels 1, . . ., k.In common with previous literature [11,25] we say that trees T a , T b ∈ T k have the same labeled shape or topology if the set of all tip partitions admitted by internal edges of T a is identical to that of T b , and we write this as T a ∼ = T b .We say that T a = T b if they have the same topology and each corresponding branch has the same length.
For any tree T a ∈ T k let m i,j be the number of edges on the path from the root to the most recent common ancestor (MRCA) of tips i and j, let M i,j be the length of this path, and let p i be the length of the pendant edge to tip i.Then, including all pairs of tips, we have two vectors: which captures the tree topology, and which captures the topology and the branch lengths.The vector M (T ) is similar to the vector of cophenetic values [29,5] (Section 2.3).We form a convex combination of these vectors, parameterized with λ ∈ [0, 1], to give .
Figure 1 provides an example of this calculation for two small trees.A metric is a mathematical notion of distance; specifying a metric gives structure and shape to a set of objects, forming a space.

Topology and lengths Topology
Proof.Since the Euclidean distance between vectors satisfies the conditions (1), ( 3) and ( 4) for being a metric, it remains to prove that d 0 (T a , T b ) = 0 ⇔ T a ∼ = T b (i.e. the distance is 0 with λ = 0 if and only if the trees have the same topology) and d λ (T a , T b ) = 0 ⇔ T a = T b for all λ ∈ (0, 1] (i.e. the distance is 0 for 0 < λ 1 if and only if the trees are identical).We will address this in three stages, showing that (1) the tree topology vector, (2) the branch-length focused vector, and (3) their convex combination each uniquely define a tree.That is, we show that for For ease of notation we restrict our attention here to binary trees; it is straightforward to extend these arguments to trees that are not binary.
1. We show that m(T ) characterizes a tree topology.Suppose that for T a , T b ∈ T k we have d 0 (T a , T b ) = 0, so m i,j (a) = m i,j (b) for all pairs i, j ∈ 1, . . ., k.Consider the tip partition created by the root of T a .That is, if the root and its two descendant edges were removed, then T a would be split into two subtrees, whose tip sets we label L and R. For all leaf pairs (i, j) with i ∈ L and j ∈ R we have m i,j (a) = 0, and therefore m i,j (b) = 0. Thus the root of T b also admits the partition {L, R}.
Similarly, any internal node n in T a partitions its descendant tips into non-empty sets L n , R n .Let the number of edges on the path from the root to n be x n .For all leaf pairs (i, j) with i ∈ L n , j ∈ R n we have m i,j (a) = x n = m i,j (b), and so there must also be an internal node in T b which partitions the leaves into the sets L n , R n .Since this is true for all internal nodes, and hence all internal edges, we have T a ∼ = T b , and d 0 is a metric on tree topologies.Note that the final k fixed entries of m(T ) are redundant for unique characterization of the topology of the tree, but are included to allow the convex combination of the topological and branch-length focused vectors.

2.
We show that M (T ) characterizes a tree using a similar argument to that of part (1).Suppose that for T a , T b ∈ T k we have d 1 (T a , T b ) = 0, so M i,j (a) = M i,j (b) for all pairs i, j ∈ 1, . . ., k.Let the length of the path from the root to internal node n be X n .Then for all i ∈ L n , j ∈ R n we have M i,j (T a ) = X n = M i,j (T b ), which means that T b also contains an internal node at distance X n from the root which admits the partition {L n , R n }.Since this holds for all internal nodes including the root (where X n = 0), we have that T a and T b have the same topology and internal branch lengths.
The final k elements of M (T ) correspond to the pendant branch lengths.When M (T a ) = M (T b ) we have that for each i ∈ 1, . . ., k the pendant branch length to tip i has length p i in both T a and T b .Thus T a and T b have the same topology and branch lengths, hence T a = T b and d 1 is a metric.3. Finally, we need to show that v λ (T ) characterizes a tree for λ ∈ (0, 1).Suppose that for

2
. It is clear that for the final k entries, that is for We therefore restrict our attention to the first k for all i, j ∈ 1, . . ., k.We show that, for any λ ∈ (0, 1), although it is possible for Equation 1to hold for some i, j ∈ 1, . . ., k it will only hold for all i, j ∈ 1, . . ., k when T a = T b .Suppose for a contradiction that we have and so d λ (T a , T b ) = 0 implies that T a and T b must share the same root partition.Now fix λ ∈ (0, 1) and consider a pair of tips x, y ∈ 1, . . ., k with m x,y (T a ) = m x,y (T b ), m x,y (T a ), m x,y (T b ) = 0, which must exist since T a = T b , using part (1).Without loss of generality, suppose that m x,y (a) − m x,y (b) = n, where n ∈ N. Then there exist at least n tips z 1 , . . ., z n for which, because the trees have the same root partition, we have for each i ∈ 1, . . ., n (see Figure 2).Pick z j so that m Our metric is fundamentally for rooted trees.A single unrooted tree, when rooted in two different places, produces two distinct rooted trees, and our distance between these will be positive.It will be large if the two distinct places chosen for the roots are separated by a long path in the original unrooted tree.However, it would be straightforward to check if two trees have the same (unrooted) topology in our metric: root both trees on the edge to the same tip and find the distance.Re-rooting a tree will induce systematic changes in v(T ), with some entries increasing and others decreasing by the same amount.The metric d λ is invariant under permutation of labels.That is, for trees T a and T b and a label permutation ).We note that alternative, similar definitions for a metric on T k are possible.In particular, the metric defined by gives similar behavior to the metric we have used.The difference between the two is that in D, the Euclidean distances are taken between the m and M vectors before they are weighted by λ.Rather than a Euclidean distance between two vectors (v for each tree), D is a weighted sum of two different metrics: the distance between m(T a ) and m(T b ) (first term in the above), and between M (T a ) and M (T b ) (second term).A benefit of D λ is that it is linear in λ, so that the changes as λ moves from 0 to 1 are more intuitive.A disadvantage is that D λ itself is not Euclidean, leading to (typically only slightly) poorer-quality visualization in MDS plots (Section 3.1).

The role of λ
The parameter λ allows the user to choose to what extent the branch lengths of a tree, vs its topology alone, contribute to the tree distance.The distance between two trees may increase or decrease as λ increases from 0 to 1. Since the topology-based vector, m, contains the number of edges along paths in the tree, and M contains the path lengths, the branch lengths are implicitly compared to 1 in the convex combination v.In other words, if the branch lengths are much larger than 1, then the entries of M will be much larger than the corresponding entries of m, and M will dominate in the expression for v even when λ is relatively small.Conversely, if the branch lengths are much less than 1, the entries of M will be much less than those of m, and a value of λ near 1 will be required in order for lengths to substantially change v.In the case when all branch lengths are equal to 1, m = M and the distance is independent of λ.The example in Figure 3 may provide some intuition.
In order to capture length-sensitive distances between trees, we may wish to use a value of λ such that neither (1 − λ)m nor λM dominate excessively, but naturally this will depend on the analysis.For a more gradual change in d λ as λ tends to 1, and for comparison of this change across different data sets, it is possible to rescale the branch lengths, for example by dividing all branch lengths by the median, or by changing the units.However, this should be done with caution because information is inevitably lost through rescaling.For example, if a phylogenetic analysis of multiple genes from the same organism had produced trees with similar topologies but different clock rates (e.g.branches in trees from gene 1 were typically twice as long as branches in trees from gene 2), this information would be obscured by rescaling.

Other metrics on labeled phylogenetic trees
Various metrics have been defined on phylogenetic trees.For a recent comparative survey, see [17].
The vector M (T ) is similar to the cophenetic vector of Cardona et al. [5], following Sokal and Rohlf [29], where M i,j is called the cophenetic value of tips i and j.Parts (1) and (2) of our proof follow directly from results in [5].Instead of the pendant branch lengths p i , Cardona et al. use the depth of each taxon, which can be considered as M i,i .This involves a repetition of information between M i,i , M j,j and M i,j whenever M i,j > 0. However, their definition does allow for the presence of nested taxa (taxa which are internal nodes of the tree).Cardona et al. also note that tree vectors such as these can be compared by any norm L p , but that the Euclidean norm L 2 , which we also use, has the benefits of being more discriminative than larger values of p, and enabling many geometrical and clustering methods.
The most widely used metric is that of Robinson-Foulds (RF) [25].However, RF and its branch-length weighted version [24] are fundamentally very different from our metric because they are defined on unrooted trees, whereas our metric emphasizes the placement of the root and all the descendant MRCAs.Similarly, the path difference metrics of Williams and Clifford [33] and Steel and Penny [31] are for unrooted trees.They compare the distance between each pair of tips in a tree; in essence, they consider the distance between tips and their MRCA, whereas our metric considers the distance between the root and the MRCA.These metrics therefore capture different characteristics of trees and are only loosely correlated with our metric.
The metric introduced by Billera, Holmes and Vogtmann (BHV) captures branch lengths as well as tree structure [3] on rooted trees.The BHV tree space is formed by mathematically 'gluing' together orthants.Each orthant corresponds to a tree topology and moving within an orthant corresponds to changing the tree's branch lengths.Moving from one orthant to an adjacent one corresponds to a nearest-neighbor interchange move.The metric is convex: for any two distinct trees T 1 and T 2 , there is a tree T 3 'in between' them, i.e. such that d BHV (T 1 , T 3 ) + d BHV (T 3 , T 2 ) = d BHV (T 1 , T 2 ).This is a mathematically appealing and useful property, in part because it allows averaging of trees [1].However, it does not allow the user to choose a balance between the topology of the tree and the branch lengths.We provide further comparisons in Figure 4.
Our metric compares trees with the same set of taxa (i.e. the same tips).As a consequence, it is suited for studies in which there is one set of taxa, and trees can be compared from different genes, inference methods, and sources of data.Our metric does not capture distances between trees with different taxa; where the taxa overlap between two trees, our approach can compare the subtrees restricted to the taxa present in both trees.In contrast, comparisons between unlabeled trees take a different form (e.g.kernel methods [22]), suitable to comparing trees on different sets of taxa.
Many phylogenetic analyses are, implicitly or explicitly, conducted in the context of a rooted tree.In the context of macroevolution, examples include estimates of times to divergence, ancestral relationships and ancestral character reconstruction.In more recent literature, most methods to link pathogen phylogenies to epidemic dynamics (phylodynamics) [30,23,7] are based on rooted phylogenetic trees.For these reasons, the fact that the relationships to the root of the tree play a central role in our metric allows it to capture intuitive similarities in groups of trees in a way that other metrics do not.

Exploring tree space
Tree spaces are large and complex.It is important to understand the 'shape' of a tree space before attempting to summarize it.Our metric creates a space which can be effectively visualized (Section 3.1) and where islands (distinct clusters) of tree topologies can be detected.We demonstrate these techniques on a sample dataset of BEAST posterior trees for Dengue fever.Finally, in Section 3.3 we describe how our metric can be used to make a principled selection of summary trees.

Visualizing tree space
Visualization techniques like multidimensional scaling (MDS) [6] have been used to explore tree space previously, but are challenged by poor-quality projections [14,2].When a set of distances is projected into a low-dimensional picture, there is typically some loss of information, which may result in a poor-quality visualization.For example, if 10 points are all 3 units away from each other, this will not project well into two dimensions; some will appear more closely grouped than others.However, if there are only 3 such points they can be arranged on a triangle, capturing the distances in two dimensions.
One approach to checking the quality of a visualization is a Shepard plot [28], which is a scatter plot of the true distance vs the MDS distance (i.e. the distance in the projection).Figure 4 shows the MDS plot of the space of trees on 6 tips (with unit branch lengths) under our metric and two others: RF [25] and BHV [3].Shepard plots are included as an indication of the quality of each projection.
Each metric captures differences in both shape (shown by color) and labeling.Our approach produces a wide range of tree distances and captures intuitive similarities (e.g. the similar chimp-human pairing in the yellow and gray triangles in Figure 4a).All 945 possible tree shapes and permutations of their labels are present in the input set of trees, and consequently there is no asymmetry that should lead to one group being separated from the rest.Our metric captures the symmetry in the space and illustrates this in the MDS projection (Figure 4a), whereas in RF and BHV (Figures 4b and 4c), poor-quality projections lead to apparent distinct tree islands where none exist.This makes detecting genuine islands in posterior sets of trees difficult using RF or BHV.The Euclidean nature of our metric means that it is well-suited to visualizations that project distances into twoor three-dimensional Euclidean space.The Shepard plots illustrate that the correspondence between the projected distances and true distances is better in our metric than the others, though the projection distance can be much smaller than the true distance (but not the converse).MDS projections are of higher quality for trees from data than in the space of all trees on 6 tips (e.g. Figure 5).

Islands in tree space
Tree inference methods explore the set of possible trees given the data, but there are many alternative trees.Bayesian Markov Chain Monte Carlo (MCMC) methods as implemented in BEAST [9] and MrBayes [15] produce a posterior set of trees, each with associated likelihoods.Distinct islands of trees within small NNI distance can share a high parsimony or likelihood [21,26].Complicating matters further, not all taxa in a dataset will have complete data at all loci.In this case, there are 'terraces' of many equally likely trees, with trees in a terrace all supporting the same subtrees for the taxa with data at a given locus [27].These facts have deep implications for tree inference and analysis, but the difficulty of detecting and interpreting tree islands has meant that the majority of analyses, particularly on large datasets, remain based on a single summary tree method such as the maximum clade credibility (MCC) tree with posterior support values illustrating uncertainty, or on maximum likelihood or parsimony trees with bootstrap supports.Our metric can detect distinct clusters or islands of close tree topologies (λ = 0) within a collection of trees.Since distance is defined by the metric that is used, these are different from previously described tree islands [21,26].We demonstrate our approach using the examples from the original paper introducing BEAST [8], where Drummond and Rambaut demonstrated their Bayesian analysis on 17 dengue virus serotype 4 sequences from [18] under varying priors for model and clock rate.As a means of comparing posterior tree distributions under different BEAST settings, we ran the xml files provided in [8] through BEAST v1.8 and analyzed the resulting trees.In Figure 5 we demonstrate MDS plots of two of these analyses: Figure 5a is a sample of the posterior under the standard GTR + Γ + I substitution model with uncorrelated lognormal-distributed relaxed molecular clock; Figure 5b is a sample from the posterior under the codon-position specific substitution model GTR + CP, with a strict clock.These analyses demonstrate some of the different signals which can be detected by visualizing the metric's tree distances: distinct islands are visible in (a), whereas in (b) there are some tight bunches of points but the posterior is not as clearly separated into distinct islands.Additionally, trees in (b) are more tightly grouped together, indicating that is less conflict in the phylogenetic signals in (b).We ran BEAST twice with the settings from (a) (using different random starting seeds), and found that the space of trees explored and accepted in each run was similar, with the same islands.It is also encouraging that the MCC tree from the first BEAST run had the same topology as that from the second run, and that this topology sits in the largest island (yellow triangle in Figure 5a).Similarly, the MCC tree is in the largest cluster in (b).
Islands are of concern for tree inference and for outcomes that require the topology of tree, which will affect ancestral character reconstruction and consequently the interpretation of many phylogenetic datasets [32].However, other analyses, and tree estimation methods themselves, take trees' branch lengths as well as topology into account.We find that islands typically merge together in the metric as λ approaches 1; the posterior becomes unimodal.

Summary trees
Summarizing groups of phylogenetic trees is challenging, particularly when there are different alternative and inconsistent topologies [13].MCC trees can summarize posterior distributions; they rely on including the clades with the strongest posterior support but where these are not concordant the resulting MCC trees can have negative branch lengths.Furthermore, the MCC tree itself may never have been sampled by the MCMC chain, casting doubt on its ability to reflect the relationships in the data.
Our metric allows us to find 'central' trees within any group of trees: a posterior set of trees, or any island or cluster of trees.To do this, we exploit the fact that our metric is  simply the Euclidean distance between the two vectors v λ (T a ) and v λ (T b ).Among N trees T i (i = 1, . . ., N ) in a posterior sample, we can find the tree closest to the average vector v = 1 N N i=1 v λ (T i ).The average vector v may not in itself represent a tree, but we can then find the tree vectors from our sample which are closest to this average.These vectors correspond to trees, T c , (not necessarily unique) which minimize the distance between v and v λ (T c ).This minimal distance is a measure of the quality of the summary: if it is small, T c is close to 'average' in the posterior.T c is known as the geometric median tree [10].The geometric median is one of a range of barycentric methods which can be used with our metric to select a tree as a representative of a group.It is also straightforward to weight trees by likelihood or other characteristics when finding the geometric median.This provides a suite of tools for summarizing collections of trees.Geometric median trees will always have been sampled by the MCMC, and will not have negative branch lengths.We found that within islands, geometric median trees are very close to the MCC tree for the island.

Concluding remarks
The fact that our metric is a Euclidean distance between two vectors whose components have an intuitive description means that simple extensions are straightforward to imagine and to compute.For example, it may be the case that the placement of a particular tip is a key question.This could occur, for example, in a real-time analysis of an outbreak, where new cases need to be placed on an existing phylogeny to determine the likely source of infection.We could form a metric that emphasizes differences in the placement of a particular tip (say, A), by weighting A's entries of m and M highly compared to all other entries.In this new metric, trees would appear similar if their placement of A was similar; patterns of ancestry among the other tips would contribute less to the distance.Indeed, it is possible to design numerous metrics, extending this one and others, and using linear combinations of existing metrics [19].
Our metric enables quantitative comparison of trees.It is relevant to viral, bacterial and higher organisms and can help to reveal distinct, likely patterns of evolution.It allows quantitative comparison of tree estimation methods and can provide a heuristic for convergence of tree estimates.There are also many applications in comparing trees derived from different data.For example, the metric can be used to detect informative sites which, when removed from sequence alignments, change the phylogeny substantially.More generally, our metric can find distances between any rooted, labeled trees with the same set of tips.It can be used to compare tree structures from a variety of scientific disciplines, including decision trees, network spanning trees, hierarchical clustering trees and language trees.

Figure 1 :
Figure 1: A tree is characterized by the vectors m and M , which are calculated as shown.These are used to calculate the distance between the trees for any λ ∈ [0, 1].Here, d 0 (T 1 , T 2 ) = 2 and d 1 (T 1 , T 2 ) = 1.96.

Figure 2 :
Figure 2: If d(T a , T b ) = 0 then T a and T b must share the same root partition, hence S 2 is the same set of tips in both trees.If m x,y (T a ) = m x,y (T b ), m x,y (T a ) − m x,y (T b ) = n (here m x,y (T a ) − m x,y (T b ) = 5 − 2 = 3), then there exist at least 3 tips z 1 , z 2 , z 3 between the root and the MRCA of x and y in T a , but positioned further from the root than the MRCA of x and y in T b .

Figure 3 :
Figure 3: Example trees from T 3 to illustrate the effect of changing λ.The distance between T a and T c (d λ (T a , T c )) is fixed for λ ∈ [0, 1] because their unmatched edges have the same length.d λ (T b , T d ) < d λ (T b , T c ) for λ ∈ (0, 1] because the edge which T c and T d share and which is not found in T b is shorter in T d than in T c .Most entries increase with λ.The only distance to decrease as λ → 1 is d λ (T a , T d ), because the difference between the lengths of their unmatched branches is less than one.

Figure 4 :
Figure 4: MDS projections of the shape of T 6 according to metrics as shown, with corresponding Shepard plots.Colors correspond to tree shapes, of which examples are shown with triangles.Symmetries correspond to permutations of the labels.In order to include the BHV metric in this comparison we assigned all branch lengths to be 1, with the result that m = M and our metric is invariant to λ ∈ [0, 1].

Figure 5 :
Figure 5: MDS plots of dengue fever trees sampled from posteriors demonstrate differences in the space of trees explored by BEAST under different settings.MCC trees are marked by yellow triangles.(a) GTR + Γ + I substitution model with uncorrelated lognormal-distributed relaxed molecular clock (b) Codon-position specific substitution model GTR + CP, with a strict clock.