With the introduction of polytomous trees, the computation of the tree topology prior becomes more interesting but also much more complicated. We now need to not only know the number of binary tree topologies (the highest resolution class), we need to know the number of tree topologies in all the other resolution classes.
There are also at least two interesting variations on the concept of a flat prior with respect to tree topology. The prior could be flat across all tree topologies, or it could be flat across resolution classes. In the first case, the star tree and any particular fully-resolved tree topology have the same prior probability. In the second case, the star tree (which is the only representative of its resolution class), gets much greater weight than any particular fully-resolved tree topology. For example, with just 10 taxa, the star tree would have (in the flat resolution class prior) a prior probabilility 2,027,025 times that of any resolved tree!
To help us compute the topology prior, create the following class TopoPriorCalculator
in a header file named topo_prior_calculator.hpp as follows:
#pragma once
#include "tree.hpp"
#include "lot.hpp"
namespace strom {
class PolytomyTopoPriorCalculator {
public:
typedef std::shared_ptr<PolytomyTopoPriorCalculator> SharedPtr;
PolytomyTopoPriorCalculator();
virtual ~PolytomyTopoPriorCalculator();
bool isResolutionClassPrior() const;
bool isPolytomyPrior() const;
bool isRooted() const;
bool isUnrooted() const;
void setNTax(unsigned n);
unsigned getNTax() const;
void chooseRooted();
void chooseUnrooted();
double getLogCount(unsigned n, unsigned m);
double getLogSaturatedCount(unsigned n);
double getLogTotalCount(unsigned n);
std::vector<double> getLogCounts();
std::vector<double> getCountsVect();
std::vector<int> getNFactorsVect();
void chooseResolutionClassPrior();
void choosePolytomyPrior();
void setC(double c);
double getC() const;
void setLogScalingFactor(double lnf);
double getLogScalingFactor() const;
virtual double getLogTopoProb(Tree::SharedPtr t);
double getLogTopologyPrior(unsigned m);
double getLogNormalizedTopologyPrior(unsigned m);
double getLogNormConstant();
std::vector<double> getTopoPriorVect();
std::vector<double> getRealizedResClassPriorsVect();
unsigned sample(Lot::SharedPtr rng);
void clear();
void reset();
private:
void recalcCountsAndPriorsImpl(unsigned n);
void recalcPriorsImpl();
unsigned _ntax;
bool _is_rooted;
bool _is_resolution_class_prior;
double _C;
bool _topo_priors_dirty;
bool _counts_dirty;
double _log_scaling_factor;
std::vector<int> _nfactors;
std::vector<double> _counts;
double _log_total_count;
std::vector<double> _topology_prior;
};
// Member function bodies go here
}
The constructor calls the clear
member function to perform initialization, and the destructor does nothing, as usual.
PolytomyTopoPriorCalculator::PolytomyTopoPriorCalculator() {
//std::cout << "PolytomyTopoPriorCalculator being created" << std::endl;
clear();
}
PolytomyTopoPriorCalculator::~PolytomyTopoPriorCalculator() {
//std::cout << "PolytomyTopoPriorCalculator being destroyed" << std::endl;
}
inline void PolytomyTopoPriorCalculator::clear() {
_topo_priors_dirty = true;
_is_rooted = false;
_is_resolution_class_prior = true;
_C = 1.0;
_ntax = 4;
_counts_dirty = true;
_log_scaling_factor = 10.0;
}
This is the function that delivers what you need: the log of the prior probability of a tree topology having m
internal nodes. This represents the resolution class prior if
_is_resolution_class_prior
is true, otherwise it represents the polytomy prior.
A flat resolution class prior places equal prior probability mass on every possible value of m
. For example, for N=7
taxa, the resolution class variable m
can range from 1 (star tree) to N-2 = 5
(fully-resolved tree) for the unrooted tree case. The prior probability of each resolution class is thus 0.2 and a particular tree in each class equals 0.2 divided by the number of distinct labeled tree topologies in that class:
m | ntopo | ptree | ntopo*ptree |
---|---|---|---|
1 | 1 | 0.2/1 = 0.200000000 | 0.2 |
2 | 56 | 0.2/56 = 0.003571429 | 0.2 |
3 | 490 | 0.2/490 = 0.000408163 | 0.2 |
4 | 1260 | 0.2/1260 = 0.00015873 | 0.2 |
5 | 945 | 0.2/945 = 0.00021164 | 0.2 |
total | 2752 | 1.0 |
where m
is the resolution class (number of internal nodes), ntopo
is the number of tree topologies in resolution class m
(computed in the recalcCountsAndPriorsImpl
function), and ptree
is the probability of any given labeled tree topology in resolution class m
.
The flat polytomy prior is simpler: it assumes the same prior probability for each tree, regardless of the resolution class. Because the total number of tree topologies over all resolution class is 1+56+490+1260+945 = 2752
in the 7-taxon example above, the prior for each tree would be 1/2752 = 0.000363372
under the flat polytomy prior.
The priors described above were flat priors. It is possible to make some resolution classes more probable in a resolution class prior, and trees with more or fewer internal nodes more probable in a polytomy prior. For example, suppose you wished resolution class m = 1
(the star tree) to be twice as probable as the class of trees having 2 internal nodes (m = 2
), class m = 2
to be twice as probable as class m = 3
, and so on. The data member _C
determines the relative probability of resolution class m
over resolution class m+1
. Setting _C = 2
yields this resolution class prior distribution:
m | ntopo | ptree | ntopo*ptree |
---|---|---|---|
1 | 1 | (16/31)/1 = 0.516129032 | 0.516129032 |
2 | 56 | (8/31)/56 = 0.004608295 | 0.258064516 |
3 | 490 | (4/31)/490 = 0.000263331 | 0.129032258 |
4 | 1260 | (2/31)/1260 = 0.000051203 | 0.064516129 |
5 | 945 | (1/31)/945 = 0.000034136 | 0.032258065 |
total | 2752 | 1.0 |
The data member _C
in the case of a polytomy prior ensures that the probability of any particular tree in resolution class m
is _C
times more probable than any particular tree in resolution class m+1
. To see how this prior is calculated, note that have these 5 constraints (again for the 7-taxon unrooted case):
n_1 p_1 + n_2 p_2 + n_3 p_3 + n_4 p_4 + n_5 p_5 = 1.0
p_1 = C p_2
p_2 = C p_3
p_3 = C p_4
p_4 = C p_5
These constraints leave only p_5 to be determined:
p_1 = C^4 p_5
p_2 = C^3 p_5
p_3 = C^2 p_5
p_4 = C^1 p_5
(n_1 C^4 + n_2 C^3 + n_3 C^2 + n_4 C^1 + n_5 C^0) p_5 = 1.0
p_5 = 1/[n_1 C^4 + n_2 C^3 + n_3 C^2 + n_4 C^1 + n_5 C^0]
For _C = 2
and a polytomy prior, we have
m | ntopo | qtree | ptree | ntopo*ptree |
---|---|---|---|---|
1 | 1 | 1*2^4 = 16 | 2^4/5889 = 0.00271693 | 1*2^4/5889 = 0.00271693 |
2 | 56 | 56*2^3 = 448 | 2^3/5889 = 0.001358465 | 56*2^3/5889 = 0.076074036 |
3 | 490 | 490*2^2 = 1960 | 2^2/5889 = 0.000679232 | 490*2^2/5889 = 0.332823909 |
4 | 1260 | 1260*2^1 = 2520 | 2^1/5889 = 0.000339616 | 1260*2^1/5889 = 0.427916454 |
5 | 945 | 945*2^0 = 945 | 2^0/5889 = 0.000169808 | 945*2^0/5889 = 0.16046867 |
total | 2752 | 16+448+1960+2520+945 = 5889 | 1.0 |
where qtree
is the unnormalized probability of a particular tree topology (i.e. ptree is qtree divided by the sum of qtree values).
With this introduction, here is the deceptively simple getLogNormalizedTopologyPrior
function. The hard part is computing the _topology_prior
vector, which is described below in the explanation of the recalcPriorsImpl
function. The log of the normalized topology prior is obtained as _topology_prior[m]
minus _topology_prior[0]
(the 0th element of _topology_prior
holds the log of the normalization constant).
double PolytomyTopoPriorCalculator::getLogNormalizedTopologyPrior(unsigned m) {
if (_topo_priors_dirty)
reset();
assert(m < _topology_prior.size());
return (_topology_prior[m] - _topology_prior[0]);
}
This function is responsible for recomputing the _topology_prior
vector. If _is_resolution_class_prior
is true, this function requires knowledge of the number of tree topologies in each resolution class, so it is called inside the recalcCountsAndPriorsImpl
after the counts are calculated.
void PolytomyTopoPriorCalculator::recalcPriorsImpl() {
_topology_prior.clear();
_topology_prior.push_back(0.0); // This will hold the normalizing constant in the end
// Figure out the maximum possible value for m, the number of internal nodes
unsigned maxm = _ntax - (_is_rooted ? 1 : 2);
if (_is_resolution_class_prior) {
// _counts vector should have length equal to maxm if everything is ok
assert(maxm == (unsigned)_counts.size());
double logC = std::log(_C);
double log_normalization_constant = 0.0;
if (_C == 1.0)
log_normalization_constant = std::log(maxm);
else if (_C > 1.0) {
// factor out largest value to avoid overflow
double term1 = logC*maxm;
double term2a = std::exp(-logC*maxm);
double term2b = std::log(1.0 - term2a);
double term3 = std::log(_C - 1.0);
log_normalization_constant = term1 + term2b + term3;
}
else {
log_normalization_constant = std::log(1.0 - std::exp(logC*maxm)) - std::log(1.0 - _C);
}
_topology_prior[0] = log_normalization_constant;
for (unsigned m = 1; m <= maxm; ++m) {
double logCterm = (double)(maxm - m)*logC;
double log_count_m = std::log(_counts[m - 1]) + _log_scaling_factor*(double)_nfactors[m - 1];
double log_v = logCterm - log_count_m;
_topology_prior.push_back(log_v);
}
}
else {
double total = 0.0;
double logC = std::log(_C);
for (unsigned m = 1; m <= maxm; ++m) {
double logCterm = (double)(maxm - m)*logC;
total += std::exp(logCterm);
_topology_prior.push_back(logCterm);
}
_topology_prior[0] = std::log(total);
}
}
This function recomputes the _counts
vector for the supplied number of internal nodes (n
) using the method outlined by Joe Felsenstein in his 2004 book and also in Felsenstein (1978). Below I’ve reproduced the table from Felsenstein (2004) showing the number of possible rooted tree topologies for any number of internal nodes (rows) and up to 8 taxa (columns):
nodes | 2 | 3 | 4 | 5 | 6 | 7 | 8 | nfactors |
---|---|---|---|---|---|---|---|---|
1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 0 |
2 | 3 | 10 | 25 | 56 | 119 | 246 | 0 | |
3 | 15 | 105 | 490 | 1918 | 6825 | 0 | ||
4 | 0.0105 | 0.1260 | 0.945 | 5.6980 | 1 | |||
5 | 0.0945 | 1.7325 | 19.0575 | 1 | ||||
6 | 1.0395 | 27.0270 | 1 | |||||
7 | 13.5135 | 1 |
The recalcCountsAndPriorsImpl
function works from left to right, calculating each column in turn. Within a column, it works down. Each cell except the first in a given column requires knowledge of two cells from the column to its left: the cell to its immediate left as well as the cell above the cell to its immediate left. This is because in order to know how many tree topologies there are for N
taxa, one needs to know how many places a taxon could be added to a tree with N-1
taxa. The N-1
taxon tree could have the same number of internal nodes as the new tree (the new taxon was added to an existing node, either creating a new polytomy or enlarging an existing one), or the N-1
taxon tree could have one fewer internal nodes (in which case the new taxon inserts a new node).
Note that this function calls recalcPriorsImpl
before returning.
void PolytomyTopoPriorCalculator::recalcCountsAndPriorsImpl(unsigned n) {
if (_is_resolution_class_prior)
_counts_dirty = true;
if (_counts_dirty) {
double scaling_factor = exp(_log_scaling_factor);
_counts.clear();
_counts.push_back(1.0); // _counts are always 1 for m = 1
int last_factor = 0;
_nfactors.clear();
_nfactors.push_back(last_factor); // never need to scale this one
// temporary variables
double a, b, c;
double max_log_double = log(DBL_MAX);
// The value of epsilon is arbitrary, but must be larger than
// zero and less than scaling_factor
double epsilon = scaling_factor/10.0;
// Compute the vector of _counts for the number of internal nodes specified
// This is the main loop over columns. z is the number of taxa minus 1
// for rooted trees, and is the number of taxa minus 2 for unrooted trees
for (unsigned z = 2; z <= n; ++z) {
// This column is one element longer than the column to its left
_counts.push_back(0.0);
_nfactors.push_back(last_factor);
// _counts[0] is always 1.0 because there is only one star tree topology
b = _counts[0];
// This is the loop over rows within the current column.
// m + 1 is the number of internal nodes.
for (unsigned m = 1; m < z; ++m) {
unsigned num_internal_nodes = m + 1;
double diff = (double)(_nfactors[m - 1] - _nfactors[m]);
double log_factor = diff*_log_scaling_factor;
assert(log_factor < max_log_double);
a = b*exp(log_factor);
b = _counts[m];
c = a*((double)(z + num_internal_nodes - 1));
if (num_internal_nodes < z) {
c += b*(double)num_internal_nodes;
}
if (c > scaling_factor) {
double incr = floor(log(c - epsilon)/_log_scaling_factor);
_nfactors[m] += (unsigned)incr;
last_factor = _nfactors[m];
_counts[m] = exp(log(c) - incr*_log_scaling_factor);
b = exp(log(b) - incr*_log_scaling_factor);
}
else
_counts[m] = c;
}
}
// Now compute the log of the total number of tree topologies
// over all possible resolution classes (i.e. number of internal nodes)
// Begin by creating a vector of log _counts and finding the
// largest value (this will be factored out to avoid overflow)
std::vector<double> v;
unsigned sz = (unsigned)_nfactors.size();
assert(sz == _counts.size());
double max_log_count = 0.0;
for (unsigned i = 0; i < sz; ++i) {
double num_factors = (double)_nfactors[i];
double log_count = num_factors*_log_scaling_factor + log(_counts[i]);
if (log_count > max_log_count)
max_log_count = log_count;
v.push_back(log_count);
}
// Compute log sum of _counts by factoring out the largest count.
// Underflow will occur, but only for _counts that are so much
// smaller than the dominant _counts that the underflow can be
// ignored for all practical purposes
double sum = 0.0;
for (std::vector<double>::const_iterator it = v.begin(); it != v.end(); ++it) {
double diff = (*it) - max_log_count;
sum += exp(diff);
}
assert(sum > 0.0);
_log_total_count = log(sum) + max_log_count;
_counts_dirty = false;
}
else {
_nfactors.clear();
_counts.clear();
// _counts_dirty ensures that _counts will be calculated if requested,
// for example by calling getCount
_counts_dirty = true;
}
// Recalculate the _topology_prior vector too
recalcPriorsImpl();
_topo_priors_dirty = false;
}
This function forces recalculation of polytomy_prior
if _is_resolution_class_prior
is false, and both _counts
and polytomy_prior
if _is_resolution_class_prior
is true (or if _counts_dirty
is true).
void PolytomyTopoPriorCalculator::reset() {
unsigned num_internal_nodes = (_is_rooted ? (_ntax - 1) : (_ntax - 2));
recalcCountsAndPriorsImpl(num_internal_nodes);
}
This function returns the count of the number of distinct labeled tree topologies for n
taxa and m
internal nodes. It depends on the _counts vector being accurate, so it calls recalcCountsAndPriors
function if n
is not equal to _ntax
. If _is_rooted
is true, assumes m
is less than _ntax
. If _is_rooted
is false, assumes m
less than _ntax
- 1.
double PolytomyTopoPriorCalculator::getLogCount(unsigned n, unsigned m) {
assert((_is_rooted && (m < n)) || (!_is_rooted && (m < n - 1)));
if (n != _ntax)
setNTax(n);
if (_counts_dirty)
reset();
double nf = (double)(_nfactors[m - 1]);
double log_count = nf*_log_scaling_factor + log(_counts[m - 1]);
return log_count;
}
Returns the number of saturated (i.e. fully-resolved and thus having as many internal nodes as possible) trees of n
taxa. Calls recalcCountsAndPriors
function if n
is not equal to _ntax
because that means that the _counts
vector is no longer accurate and must be recalculated.
double PolytomyTopoPriorCalculator::getLogSaturatedCount(unsigned n) {
if (n != _ntax)
setNTax(n);
if (_counts_dirty)
reset();
unsigned last = (unsigned)(_counts.size() - 1);
double nf = (double)(_nfactors[last]);
double log_count = nf*_log_scaling_factor + log(_counts[last]);
return log_count;
}
This function constructs a vector in which the element having index m
(i = 0, 1, …, maximum number of internal nodes) represents the natural logarithm of the number of tree topologies having m
internal nodes. If _counts_dirty
is true, it recomputes the _counts
vectors first. The 0
th element of the returned vector holds the natural log of the total number of tree topologies (log of the sum of all other elements).
std::vector<double> PolytomyTopoPriorCalculator::getLogCounts() {
if (_is_resolution_class_prior)
_counts_dirty = true;
if (_counts_dirty)
reset();
std::vector<double> v;
unsigned sz = _ntax - (_is_rooted ? 0 : 1);
v.reserve(sz);
v.push_back(_log_total_count);
for (unsigned i = 1; i < sz; ++i) {
double log_Tnm = log(_counts[i - 1]) + _log_scaling_factor*(double)(_nfactors[i - 1]);
v.push_back(log_Tnm);
}
return v;
}
This function returns the natural log of the total number of trees for n
taxa, including all resolution classes from the star tree to fully resolved (saturated) trees. Calls recalcCountsAndPriors
function if n
is not equal to _ntax
or if not using the resolution class prior (in which case _counts
have not been calculated).
double PolytomyTopoPriorCalculator::getLogTotalCount(unsigned n) {
if (n != _ntax)
setNTax(n);
if (_counts_dirty)
reset();
return _log_total_count;
}
Constructs a vector of realized resolution class priors from the values in the _topology_prior
vector. If _topo_priors_dirty
is true, it recomputes the _topology_prior
vectors first. The m
th element of the returned vector is set to T_{n,m}*_topology_prior[m]
for m > 0
. The 0
th element of the returned vector holds the normalization constant (sum of all other elements). This function is not efficient because it is intended only to be used for providing information to the user on request. Table 2, p. 248, in the Lewis, Holder, and Holsinger (2005) presented (normalized) values from this vector, as do the columns labeled ntopo*ptree
in the tables presented above in the documentation for the getLogNormalizedTopologyPrior
member function.
std::vector<double> PolytomyTopoPriorCalculator::getRealizedResClassPriorsVect() {
if (!_is_resolution_class_prior)
_counts_dirty = true;
if (_topo_priors_dirty || _counts_dirty)
reset();
std::vector<double> v;
v.reserve(_topology_prior.size());
v.push_back(0.0);
unsigned sz = (unsigned)_topology_prior.size();
// First loop will be to determine largest value, which will be factored out
// the second time through so that the total does not overflow
double log_factored_out = 0.0;
for (unsigned i = 1; i < sz; ++i) {
double c = _counts[i - 1];
double nf = (double)_nfactors[i - 1];
double log_Tnm = log(c) + _log_scaling_factor*nf;
double log_prior = log_Tnm + _topology_prior[i];
v.push_back(log_prior);
if (log_prior > log_factored_out)
log_factored_out = log_prior;
}
// Now we can compute the total
double total = 0.0;
std::vector<double>::const_iterator it = v.begin();
for (++it; it != v.end(); ++it) {
total += exp((*it) - log_factored_out);
}
v[0] = log(total) + log_factored_out;
return v;
}
Samples a resolution class (i.e. number of internal nodes) from the realized resolution class distribution. This function is not particularly efficient because it calls PolytomyTopoPriorCalculator::getRealizedResClassPriorsVect
, resulting in an unnecessary vector copy operation.
unsigned PolytomyTopoPriorCalculator::sample(Lot::SharedPtr rng) {
std::vector<double> v = getRealizedResClassPriorsVect();
double u = rng->uniform();
double z = v[0];
double cum = 0.0;
for (unsigned i = 1; i < v.size(); ++i) {
cum += exp(v[i] - z);
if (u <= cum)
return i;
}
assert(0);
return (unsigned)(v.size() - 1);
}
The setNTax
function returns immediately if _ntax
equals new_ntax
; however, if new_ntax
differs from _ntax
, setNTax
sets _ntax
to new_ntax
and sets _topo_priors_dirty
to true so that subsequent requests for prior probabilities or counts will trigger recalculation. The getNTax
function returns the value of the data member _ntax
.
void PolytomyTopoPriorCalculator::setNTax(unsigned new_ntax) {
if (_ntax != new_ntax) {
// Set _ntax to the new value
assert(new_ntax > (unsigned)(_is_rooted ? 1 : 2));
_ntax = new_ntax;
_counts_dirty = true;
_topo_priors_dirty = true;
}
}
unsigned PolytomyTopoPriorCalculator::getNTax() const {
return _ntax;
}
The chooseRooted
function sets the _is_rooted
data member to true
. There are more rooted than unrooted trees for the same value of _ntax
, so this setting is important when asking questions that require knowledge of the numbers of possible trees. The chooseUnrooted
function sets the _is_rooted
data member to false
.
void PolytomyTopoPriorCalculator::chooseRooted() {
if (!_is_rooted) {
_is_rooted = true;
_topo_priors_dirty = true;
}
}
void PolytomyTopoPriorCalculator::chooseUnrooted() {
if (_is_rooted) {
_is_rooted = false;
_topo_priors_dirty = true;
}
}
These two functions can be used to choose between a resolution class prior or the alternative (the polytomy prior).
void PolytomyTopoPriorCalculator::chooseResolutionClassPrior() {
if (!_is_resolution_class_prior) {
_is_resolution_class_prior = true;
_topo_priors_dirty = true;
}
}
void PolytomyTopoPriorCalculator::choosePolytomyPrior() {
if (_is_resolution_class_prior) {
_is_resolution_class_prior = false;
_topo_priors_dirty = true;
}
}
These two functions either set or return the value of the data member _C
, which determines the nature of the resolution class or polytomy prior. For the resolution class prior, each resolution class m
has probability _C
times greater than resolution class m+1
. For the polytomy prior, a particular tree with m
internal nodes has probability _C
times greater than a particular tree with m+1
internal nodes.
void PolytomyTopoPriorCalculator::setC(double c) {
assert(c > 0.0);
if (c != _C) {
_C = c;
_topo_priors_dirty = true;
}
}
double PolytomyTopoPriorCalculator::getC() const {
return _C;
}
These two functions set or return the value of the data member _log_scaling_factor
data member.
void PolytomyTopoPriorCalculator::setLogScalingFactor(double lnf) {
assert(lnf > 0.0);
if (lnf != _log_scaling_factor) {
_log_scaling_factor = lnf;
_counts_dirty = true;
_topo_priors_dirty = true;
}
}
double PolytomyTopoPriorCalculator::getLogScalingFactor() const {
return _log_scaling_factor;
}
Returns copy of the _counts
vector, which contains in its (m-1)
th element the number of tree topologies having exactly m
internal nodes. Note that you will need to also call getNFactorsVect
if there is any chance that some of the _counts
are larger than exp(_log_scaling_factor)
. In such cases, the actual log count equals _nfactors[m]*_log_scaling_factor + log(count[m - 1]
.
std::vector<double> PolytomyTopoPriorCalculator::getCountsVect() {
if (_counts_dirty)
reset();
return _counts;
}
Returns copy of the _nfactors
vector, which contains in its (m-1)
th element the number of times _counts[m]
has been rescaled by dividing by the scaling factor (the log of which is _log_scaling_factor
).
std::vector<int> PolytomyTopoPriorCalculator::getNFactorsVect() {
if (_counts_dirty)
reset();
return _nfactors;
}
Returns copy of the _topology_prior
vector, which contains in its m
th element the unnormalized prior for tree topologies having exactly m
internal nodes. The 0
th element of _topology_prior
holds the normalizing constant.
std::vector<double> PolytomyTopoPriorCalculator::getTopoPriorVect() {
if (_topo_priors_dirty)
reset();
return _topology_prior;
}
Returns result of call to getLogTopologyPrior(m)
, where m
is the number of internal nodes in the supplied tree t
.
double PolytomyTopoPriorCalculator::getLogTopoProb(Tree::SharedPtr t) {
unsigned n = t->numLeaves();
assert(n > 3);
if (n != getNTax())
setNTax(n);
unsigned m = t->numInternals();
double topo_prior = getLogTopologyPrior(m);
return topo_prior;
}
Returns the natural logarithm of the unnormalized topology prior. This represents the resolution class prior if _is_resolution_class_prior
is true, otherwise it represents the polytomy prior.
double PolytomyTopoPriorCalculator::getLogTopologyPrior(unsigned m) {
if (_topo_priors_dirty)
reset();
assert(m < _topology_prior.size());
return _topology_prior[m];
}
Returns the natural logarithm of the normalizing constant for the topology prior. This value is stored in _topology_prior[0]
.
double PolytomyTopoPriorCalculator::getLogNormConstant() {
if (_topo_priors_dirty)
reset();
return _topology_prior[0];
}
The isResolutionClassPrior
function returns the current value of the data member _is_resolution_class_prior
and isPolytomyPrior
returns !isResolutionClassPrior()
.
bool PolytomyTopoPriorCalculator::isResolutionClassPrior() const {
return _is_resolution_class_prior;
}
bool PolytomyTopoPriorCalculator::isPolytomyPrior() const {
return !_is_resolution_class_prior;
}
The isRooted
function returns the current value of the data member _is_rooted
, and isUnrooted
returns !isRooted()
.
bool PolytomyTopoPriorCalculator::isRooted() const {
return _is_rooted;
}
bool PolytomyTopoPriorCalculator::isUnrooted() const {
return !_is_rooted;
}
Felsenstein, J. 1978. The number of evolutionary trees. Systematic Zoology 27:27-33 (see also the 1981 correction Systematic Zoology 30:122). https://doi.org/10.2307/2412810/
Felsenstein, J. 2004. Inferring Phylogeny. Sinauer, Sunderland, Massachusetts. ISBN: 9780878931774
Lewis, P. O., Holder, M. T., and Holsinger, K. E. 2005. Polytomies and bayesian phylogenetic inference. Systematic Biology 54(2):241–253. https://doi.org/10.1080/10635150590924208