Estimate edge lengths, in coalescent units, on an unrooted species tree from a table of resolved quartet counts from a collection of gene trees.
estimateEdgeLengths(tree, rqt, terminal = 1, method = "simpleML", lambda = 1/2)
a phylo object, giving a resolved tree on which to estimate edge lengths
a resolved quartet table, as from
an edge length to assign to terminal edges, whose lengths cannot be estimated
a positive parameter for the
While the argument
tree may be supplied as rooted or unrooted, metric or topological,
only its unrooted topology will be used for determining new metric estimates.
Counts of quartets for all those quartets which define a single edge
on the tree (i.e., whose internal edge is the
single edge on the unrooted input tree) are summed, and from this an
estimate of the branch length is computed. If
method= "simpleML" this is the maximum likelihood estimate.
method="simpleBayes" this is the Bayesian estimate of Theorem 2
of \insertCiteSayMir16;textualMSCquartets, using parameter
gives a flat prior on [1/3,1] for the probability of the quartet displayed on the species tree.
These methods are referred to as ‘simple’ since they use only the quartets defining a single edge of the species tree. Quartets with central edges composed of several edges in the species tree are ignored.
Note that branch length estimates may be 0 (if the count for the quartet
displayed on the input tree is not dominant),
the counts for quartet topologies not displayed on the input tree are all 0, and
an unrooted metric tree with the same topology as
tree, of type phylo
1 2 3 4 5 6 7 8 9 10 11 12 13 14
gtrees=read.tree(file=system.file("extdata","dataGeneTreeSample",package="MSCquartets")) taxanames=taxonNames(gtrees) QT=quartetTable(gtrees,taxanames[1:6]) RQT=quartetTableResolved(QT) DQT=quartetTableDominant(RQT) tree=QDS(DQT) write.tree(tree) plot(tree) metricMTree=estimateEdgeLengths(tree,RQT,method="simpleML") write.tree(metricMTree) plot(metricMTree) metricBTree=estimateEdgeLengths(tree,RQT,method="simpleBayes") write.tree(metricBTree) plot(metricBTree)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.