Phylogenetic Trees

Introduction

Phylogenetic trees represent the evolutionary history of some species of consideration. Here we consider phylogenetic trees with branch lengths as defined in [SS03].

Construction

phylogenetic_treeFunction
phylogenetic_tree(T::Type{<:Union{Float64, QQFieldElem}}, newick::String)

Constructs a rooted phylogenetic tree with Newick representation newick. T indicates the numerical type of the edge lengths.

Examples

Make a phylogenetic tree with 4 leaves from its Newick representation and print its taxa and cophenetic matrix.

julia> phylo_t = phylogenetic_tree(Float64, "((H:3,(C:1,B:1):2):1,G:4);");

julia> taxa(phylo_t)
4-element Vector{String}:
 "B"
 "C"
 "G"
 "H"

julia> cophenetic_matrix(phylo_t)
4×4 Matrix{Float64}:
 0.0  2.0  8.0  6.0
 2.0  0.0  8.0  6.0
 8.0  8.0  0.0  8.0
 6.0  6.0  8.0  0.0
source
phylogenetic_tree(M::Matrix{T}, taxa::Vector{String}) where T <: Union{Float64, QQFieldElem}

Constructs a phylogenetic tree with cophenetic matrix M and taxa taxa. The matrix M must be ultrametric, otherwise an error will be thrown.

Examples

Make a phylogenetic tree on 4 taxa with given cophenetic matrix and print one Newick representation.

julia> mat = [0. 2 8 6; 2 0 8 6; 8 8 0 8; 6 6 8 0]
4×4 Matrix{Float64}:
 0.0  2.0  8.0  6.0
 2.0  0.0  8.0  6.0
 8.0  8.0  0.0  8.0
 6.0  6.0  8.0  0.0

julia> tax = ["Bonobo", "Chimpanzee", "Gorilla", "Human"]
4-element Vector{String}:
 "Bonobo"
 "Chimpanzee"
 "Gorilla"
 "Human"

julia> tree_mat = phylogenetic_tree(mat, tax);

julia> newick(tree_mat)
"Gorilla:4,(Human:3,(Bonobo:1,Chimpanzee:1):2):1;"
source

Some Helpful Functions

adjacency_treeFunction
adjacency_tree(ptree::PhylogeneticTree)

Returns the underlying graph of the phylogenetic tree ptree.

Examples

Make a phylogenetic tree with given Newick format and print its underlying graph.

julia> ptree = phylogenetic_tree(Float64, "((H:3,(C:1,B:1):2):1,G:4);");

julia> adjacency_tree(ptree)
Directed graph with 7 nodes and the following edges:
(1, 2)(1, 7)(2, 3)(2, 4)(4, 5)(4, 6)
source
is_equidistantFunction
is_equidistant(ptree::PhylogeneticTree)

Checks if the phylogenetic tree ptree is equidistant.

Examples

Make a phylogenetic tree with given Newick format and check if it is equidistant.

julia> ptree = phylogenetic_tree(Float64, "((H:3,(C:1,B:1):2):1,G:4);");

julia> is_equidistant(ptree)
true
source
cophenetic_matrixFunction
cophenetic_matrix(ptree::PhylogeneticTree)

Returns the cophenetic matrix of the phylogenetic tree ptree.

Examples

Make a phylogenetic tree with given Newick format and print its cophenetic matrix.

julia> ptree = phylogenetic_tree(Float64, "((H:3,(C:1,B:1):2):1,G:4);");

julia> cophenetic_matrix(ptree)
4×4 Matrix{Float64}:
 0.0  2.0  8.0  6.0
 2.0  0.0  8.0  6.0
 8.0  8.0  0.0  8.0
 6.0  6.0  8.0  0.0
source
taxaFunction
taxa(ptree::PhylogeneticTree)

Returns the taxa of the phylogenetic tree ptree.

Examples

Make a phylogenetic tree with given Newick format and print its taxa.

julia> ptree = phylogenetic_tree(Float64, "((H:3,(C:1,B:1):2):1,G:4);");

julia> taxa(ptree)
4-element Vector{String}:
 "B"
 "C"
 "G"
 "H"
source
newickFunction
newick(ptree::PhylogeneticTree)

Returns a Newick representation of the phylogenetic tree ptree.

Examples

Make a phylogenetic tree from a matrix and print a Newick representation of it.

julia> mat = [0. 2 8 6; 2 0 8 6; 8 8 0 8; 6 6 8 0]
4×4 Matrix{Float64}:
 0.0  2.0  8.0  6.0
 2.0  0.0  8.0  6.0
 8.0  8.0  0.0  8.0
 6.0  6.0  8.0  0.0

julia> tax = ["Bonobo", "Chimpanzee", "Gorilla", "Human"]
4-element Vector{String}:
 "Bonobo"
 "Chimpanzee"
 "Gorilla"
 "Human"

julia> tree_mat = phylogenetic_tree(mat, tax);

julia> newick(tree_mat)
"Gorilla:4,(Human:3,(Bonobo:1,Chimpanzee:1):2):1;"
source
tropical_median_consensusFunction
tropical_median_consensus(arr::Vector{PhylogeneticTree{T}})

Computes the tropical median consensus tree of the phylogenetic trees from the vector arr.

Examples

Compute the tropical median consensus of three trees and print one of its Newick representations.

julia> t1 = phylogenetic_tree(Float64, "((H:30,(C:10,B:10):20):10,G:40);");

julia> t2 = phylogenetic_tree(Float64, "(((H:10,C:10):20,B:30):10,G:40);");

julia> t3 = phylogenetic_tree(Float64, "((H:25,C:25):15,(B:15,G:15):25);");

julia> arr = [t1, t2, t3];

julia> tc = tropical_median_consensus(arr);

julia> newick(tc)
"G:40,(B:35,(C:30,H:30):5):5;"
source
tropical_median_consensus(trees::Vararg{PhylogeneticTree, N}) where {N}

Computes the tropical median consensus tree of any number of phylogenetic trees given as parameters.

Examples

Compute the tropical median consensus of three trees and print one of its Newick representations.

julia> t1 = phylogenetic_tree(Float64, "((H:30,(C:10,B:10):20):10,G:40);");

julia> t2 = phylogenetic_tree(Float64, "(((H:10,C:10):20,B:30):10,G:40);");

julia> t3 = phylogenetic_tree(Float64, "((H:25,C:25):15,(B:15,G:15):25);");

julia> tc = tropical_median_consensus(t1, t2, t3);

julia> newick(tc)
"G:40,(B:35,(C:30,H:30):5):5;"
source