Multispecies Coalescent Process is a stochastic process model that describes the genealogical relationships for a sample of DNA sequences taken from several species. It represents the application of coalescent theory to the case of multiple species. The multispecies coalescent results in cases where the relationships among species for an individual gene (the gene tree) can differ from the broader history of the species (the species tree). It has important implications for the theory and practice of phylogenetics and for understanding genome evolution.
A gene tree is a binary graph that describes the evolutionary relationships between a sample of sequences for a non-recombining locus. A species tree describes the evolutionary relationships between a set of species, assuming tree-like evolution. However, several processes can lead to discordance between gene trees and species trees. The Multispecies Coalescent model provides a framework for inferring species phylogenies while accounting for ancestral polymorphism and gene tree-species tree conflict. The process is also called the Censored Coalescent.
Besides species tree estimation, the multispecies coalescent model also provides a framework for using genomic data to address a number of biological problems, such as estimation of species divergence times, population sizes of ancestral species, species delimitation, and inference of cross-species gene flow.Bruce Rannala, Scott V. Edwards, Adam Leaché, and Ziheng Yang (2020). The Multispecies Coalescent Model and Species Tree Inference. In Scornavacca, C., Delsuc, F., and Galtier, N., editors, Phylogenetics in the Genomic Era, chapter No. 3.3, pp. 3.3:1–3.3:21. No commercial publisher | Authors open access book.
Gene tree-species tree congruence
If we consider a rooted three-taxon tree, the simplest non-trivial phylogenetic tree, there are three different tree topologies but four possible gene trees. The existence of four distinct gene trees despite the smaller number of topologies reflects the fact that there are topologically identical gene tree that differ in their coalescent times. In the type 1 tree the alleles in species A and B coalesce after the speciation event that separated the A-B lineage from the C lineage. In the type 2 tree the alleles in species A and B coalesce before the speciation event that separated the A-B lineage from the C lineage (in other words, the type 2 tree is a deep coalescence tree). The type 1 and type 2 gene trees are both congruent with the species tree. The other two gene trees differ from the species tree; the two discordant gene trees are also deep coalescence trees.
The distribution of times to coalescence is actually continuous for all of these trees. In other words, the exact coalescent time for any two loci with the same gene tree may differ. However, it is convenient to break up the trees based on whether the coalescence occurred before or after the earliest speciation event.
Given the internal branch length in coalescent units it is straightforward to calculate the probability of each gene tree. For diploid organisms the branch length in coalescent units is the number of generations between the speciation events divided by twice the effective population size. Since all three of the deep coalescence tree are equiprobable and two of those deep coalescence tree are discordant it is easy to see that the probability that a rooted three-taxon gene tree will be congruent with the species tree is:
Where the branch length in coalescent units (T) is also written in an alternative form: the number of generations (t) divided by twice the effective population size (Ne). Pamilo and Nei also derived the probability of congruence for rooted trees of four and five taxa as well as a general upper bound on the probability of congruence for larger trees. Rosenberg followed up with equations used for the complete set of topologies (although the large number of distinct phylogenetic trees that becomes possible as the number of taxa increases makes these equations impractical unless the number of taxa is very limited).
The phenomenon of hemiplasy is a natural extension of the basic idea underlying gene tree-species tree discordance. If we consider the distribution of some character that disagrees with the species tree it might reflect homoplasy (multiple independent origins of the character or a single origin followed by multiple losses) or it could reflect hemiplasy (a single origin of the trait that is associated with a gene tree that disagrees with the species tree).
The phenomenon called incomplete lineage sorting (often abbreviated ILS in the scientific literatures) is linked to the phenomenon. If we examine the illustration of hemiplasy with using a rooted four-taxon tree (see image to the right) the lineage between the common ancestor of taxa A, B, and C and the common ancestor of taxa A and B must be polymorphic for the allele with the derived trait (e.g., a transposable element insertion) and the allele with the ancestral trait. The concept of incomplete lineage sorting ultimately reflects on persistence of polymorphisms across one or more speciation events.
Mathematical description of the multispecies coalescent The probability density of the gene trees under the multispecies coalescent model is discussed along with its use for parameter estimation using multi-locus sequence data.
Assumptions In the basic multispecies coalescent model, the species phylogeny is assumed to be known. Also, complete isolation after species divergence, with no migration, hybridization, or introgression is also assumed. Additionally, we assume no recombination so that all the sites within the locus share the same gene tree (topology and coalescent times). However, the basic model can be extended in different ways to accommodate migration or introgression, population size changes, recombination. Authors open access book.
Data and model parameters The model and implementation of this method can be applied to any species tree. As an example, the species tree of the great apes: humans (H), chimpanzees (C), gorillas (G) and orangutans (O) is considered. The topology of the species tree, , is assumed known and fixed in the analysis (Figure 1). Let be the entire data set, where represent the sequence alignment at locus , with for a total of loci.
The population size of a current species is considered only if more than one individual is sampled from that species at some loci.
The parameters in the model for the example of Figure 1 include the three divergence times , and and population size parameters for humans; for chimpanzees; and , and for the three ancestral species.
The divergence times ('s) are measured by the expected number of mutations per site from the ancestral node in the species tree to the present time (Figure 1 of Rannala and Yang, 2003).
Therefore, the parameters are .
Distribution of gene genealogies The joint distribution of is derived directly in this section. Two sequences from different species can coalesce only in one populations that are ancestral to the two species. For example, sequences H and G can coalesce in populations HCG or HCGO, but not in populations H or HC. The coalescent processes in different populations are different.
For each population, the genealogy is traced backward in time, until the end of the population at time , and the number of lineages entering the population and the number of lineages leaving it are recorded. For example, and , for population H (Table 1). This process is called a censored coalescent process because the coalescent process for one population may be terminated before all lineages that entered the population have coalesced. If the population consists of disconnected subtrees or lineages.
With one time unit defined as the time taken to accumulate one mutation per site, any two lineages coalesce at the rate . The waiting time until the next coalescent event, which reduces the number of lineages from to has exponential density
If , the probability that no coalescent event occurs between the last one and the end of the population at time ; i.e. during the time interval . This probability is and is 1 if .
(Note: One should recall that the probability of no events over time interval for a Poisson process with rate is . Here the coalescent rate when there are lineages is .)
In addition, to derive the probability of a particular gene tree topology in the population, if a coalescent event occurs in a sample of lineages, the probability that a particular pair of lineages coalesce is .
Multiplying these probabilities together, the joint probability distribution of the gene tree topology in the population and its coalescent times as
.
The probability of the gene tree and coalescent times for the locus is the product of such probabilities across all the populations. Therefore, the gene genealogy of Figure 1, we have
Likelihood-based inference The gene genealogy at each locus is represented by the tree topology and the coalescent times . Given the species tree and the parameters on it, the probability distribution of is specified by the coalescent process as
,
where is the probability density for the gene tree at locus locus , and the product is because we assume that the gene trees are independent given the parameters.
The probability of data given the gene tree and coalescent times (and thus branch lengths) at the locus, , is Felsenstein's phylogenetic likelihood. Due to the assumption of independent evolution across the loci,
The likelihood function or the probability of the sequence data given the parameters is then an average over the unobserved gene trees
where the integration represents summation over all possible gene tree topologies () and, for each possible topology at each locus, integration over the coalescent times . This is in general intractable except for very small species trees.
In Bayesian inference, we assign a prior on the parameters, , and then the posterior is given as
where again the integration represents summation over all possible gene tree topologies () and integration over the coalescent times . In practice this integration over the gene trees is achieved through a Markov chain Monte Carlo algorithm, which samples from the joint conditional distribution of the parameters and the gene trees
The above assumes that the species tree is fixed. In species-tree estimation, the species tree () changes as well, so that the joint conditional distribution (from which the MCMC samples) is
where is the prior on species trees.
As a major departure from two-step summary methods, full-likelihood methods average over the gene trees. This means that they make use of information in the branch lengths (coalescent times) on the gene trees and accommodate their uncertainties (due to limited sequence length in the alignments) at the same time. It also explains why full-likelihood methods are computationally much more demanding than two-step summary methods.
Markov chain Monte Carlo under the multispecies coalescent