- Research article
- Open Access
Multilocus coalescent analyses reveal the demographic history and speciation patterns of mouse lemur sister species
BMC Evolutionary Biology volume 14, Article number: 57 (2014)
Debate continues as to whether allopatric speciation or peripatric speciation through a founder effect is the predominant force driving evolution in vertebrates. The mouse lemurs of Madagascar are a system in which evolution has generated a large number of species over a relatively recent time frame. Here, we examine speciation patterns in a pair of sister species of mouse lemur, Microcebus murinus and M. griseorufus. These two species have ranges that are disparately proportioned in size, with M. murinus showing a much more extensive range that marginally overlaps that of M. griseorufus. Given that these two species are sister taxa, the asymmetric but overlapping geographic ranges are consistent with a model of peripatric speciation. To test this hypothesis, we analyze DNA sequence data from four molecular markers using coalescent methods. If the peripatric speciation model is supported, we predict substantially greater genetic diversity in M. murinus, relative to M. griseorufus. Further, we expect a larger effective population size in M. murinus and in the common ancestor of the two species than in M. griseorufus, with a concomitant decrease in gene tree/species tree incongruence in the latter and weak signs of demographic expansion in M. murinus.
Our results reject a model of peripatric divergence. Coalescent effective population size estimates were similar for both extant species and larger than that estimated for their most recent common ancestor. Gene tree results show similar levels of incomplete lineage sorting within species with respect to the species tree, and locus-specific estimates of genetic diversity are concordant for both species. Multilocus demographic analyses suggest range expansions for M. murinus, with this species also experiencing more recent population declines over the past 160 thousand years.
Results suggest that speciation occurred in allopatry from a common ancestor narrowly distributed throughout southwest Madagascar, with subsequent range expansion for M. murinus. Population decline in M. murinus is likely related to patterns of climate change in Madagascar throughout the Pleistocene, potentially exacerbated by continual anthropogenic perturbation. Genome-level data are needed to quantify the role of niche specialization and adaptation in shaping the current ranges of these species.
In recent years, there has been growing interest in using genealogical tree structure to reconstruct the demographic and temporal context of diverging populations and species [1–4]. This tree-based approach permits the examination of lineages as they have diverged in the past using principles derived from coalescent theory, since ancestral polymorphisms that are shared among lineages are sorted during population segregation and speciation. A comparative demographic approach utilizing tree-based methods has been implemented in a variety of biological sub-fields, including systematics, phylogeography, conservation, and life history. Further, estimates of divergence times, migration rates, and effective population sizes among groups of organisms are now inferred with consideration of the genealogical structure of multilocus data sets [5–9]. These estimates may be synthesized and compared among sister taxa and used to infer the impact of the biogeographic or climatic context or other historical processes on patterns of genetic differences among groups [10–12].
In general, it is assumed that a complex suite of geologic, climatic, and population genetic forces have led to divergence and speciation in allopatry, where a reproductive barrier effectively divides an ancestral species into two populations of roughly equal size. With the continual progress of next-generation DNA sequencing and the increasing abundance of sequence data from across the genome, researchers are beginning to test alternate models of speciation including parapatric divergence along ecological clines [13–15]. Few studies, however, have used recently developed multilocus coalescent methods to fully understand the geography of speciation and the frequency of peripatric events in nature. Peripatric speciation can be seen as a subset of allopatric speciation, with a founder effect leading to the formation of a new species as a small population becomes physically separated from a broadly distributed ancestor. By definition then, the effective population size (Ne) of the diverging population will be smaller than that of the common ancestor in a peripatric scenario. Although evidence of peripatric divergence is commonly inferred throughout archipelagos using traditional phylogenetic approaches (e.g. [16, 17]), examples from strictly mainland taxa are relatively rare.
New algorithms can effectively estimate demographic parameters on a species tree to help disentangle the historical context of divergence and speciation and to test the plausibility of a peripatric scenario in comparison with the generally accepted null model of allopatric divergence [2, 3, 9]. Mitochondrial DNA (mtDNA) has been commonly used in studies of speciation given its small effective population size relative to nuclear DNA (nDNA) loci. Accordingly, mitochondrial gene trees undergo lineage sorting much more rapidly, and thus have a higher probability of yielding gene trees that are reciprocally monophyletic for the hypothesized species. However, the stochasticity inherent in single-locus gene trees makes this a controversial approach for estimating species histories, as previous authors have discussed (refer to [18–20]). Therefore, the implementation of a comparative demographic approach incorporating multiple unlinked loci is now regarded as a more powerful approach.
Mouse lemurs (Microcebus spp.) are small, nocturnal primates that inhabit a wide range of habitats throughout Madagascar. As many as 21 species have been described within the genus [21–25], thus a remarkable increase from the two species taxonomy that was recognized as recently as 1994 [26, 27]. The geography of the island has historically been considered to play a significant role in explaining divergence patterns among lemurs and other vertebrate species (e.g. [26, 28–32]) and specifically Microcebus species [33–36]. Within the genus, however, M. murinus is unique in having a relatively vast geographic distribution. Whereas all other Microcebus species have geographically limited ranges --- some restricted to isolated forest fragments --- M. murinus is distributed along most of the western half and along the south of the island (Figure 1A). Over the extensive range of M. murinus, the species has been found to overlap with multiple congeners [37, 38], many of which exhibit a high degree of endemicity. One example of this pattern of sympatry couples M. murinus with M. griseorufus, with multiple lines of evidence for distribution overlap and potential hybridization between the two species in southern Madagascar [39–41]. Within their respective ranges, M. murinus and M. griseorufus also exhibit marked differences in abundance. Microcebus murinus is a common, widespread species throughout the west (Figure 1A) that has been studied extensively since its description [26, 27, 42]. Conversely, M. griseorufus has been recorded from relatively few locations within a narrow range in the southwest of the island in the recent past [21, 24].
Despite their partially overlapping distributions and large genetic distances , M. murinus and M. griseorufus are sister taxa [34, 43, 44]. Given the current ranges of M. murinus and M. griseorufus and their status as sister species, we consider two possible scenarios for the speciation event within this small clade. One scenario is that M. griseorufus was derived from the most recent common ancestor (MRCA) as a small founding population, with divergence the result of random genetic drift (Figure 2A). Alternatively, we consider that the MRCA was a small, narrowly distributed population located in southwestern Madagascar. Under this alternative hypothesis, the two modern species resulted from an event that partitioned that ancestral species into comparably sized populations prior to a substantial range expansion in M. murinus (Figure 2B). These two hypotheses differ most critically in their assumptions regarding the ancestral range and population size: was the MRCA more like modern M. murinus (hypothesis 1: peripatric speciation) or more like modern M. griseorufus (hypothesis 2: allopatric speciation)?
To test these alternative scenarios, we analyze multilocus DNA sequence data with coalescent-based methods to estimate relevant demographic parameters. If M. griseorufus diverged as a small founding population from a large ancestral population (Figure 2A), then we expect 1) that M. murinus has a greater effective population size and higher levels of genetic diversity than M. griseorufus, 2) that estimates of effective population size for the MRCA are more similar to those for M. murinus than to those for M. griseorufus, 3) that historical changes in population size in M. murinus were minimal, and finally, 4) that M. griseorufus is more likely to segregate as a clade in individual gene trees. Under the alternative hypothesis of recent allopatric speciation (Figure 2B), we expect that 1) that M. murinus and M. griseorufus have similar effective population sizes and levels of genetic diversity, 2) that Ne estimates for the MRCA are smaller than those for contemporary species or similar to M. griseorufus, 3) that significant population expansion has occurred in M. murinus, and finally, 4) that individual gene trees are more likely to show incomplete lineage sorting between the two species.
To estimate relevant demographic parameters for both species, we mined GenBank for all available high-quality sequences as both M. murinus and M. griseorufus have been used in a number of recent molecular studies (e.g. [41, 44]). We sought to utilize genes with sufficient sequence variation to reliably estimate coalescent parameters. We also chose data sets with samples covering a large portion of the known geographic range of both species to assemble the most geographically complete data set possible (Figure 1B; Additional file 1: Table S1). Although recent studies have suggested that M. murinus may possibly contain at least three undocumented species, [43–45], we treat these populations as M. murinus sensu lato for the purpose of this study. The final data set consisted of the following nuclear loci: alpha enolase intron (ENOL: 916 bp), alpha fibrinogen intron (FIB: 608 bp), von Willebrand factor intron (VWF: 795 bp). Although the Adora 3 gene has been sequenced numerous times within Microcebus, levels of variation were too low to provide any meaningful information for this study. We also assembled and concatenated information from two mitochondrial loci including cytochrome b (cytb) and cytochrome c oxidase subunit II (COII) for a total of 1,141 bp. The total number of sequences per locus was as follows: ENOL = 110 sequences; FIB = 124 sequences; VWF = 86 sequences; mtDNA = 55 concatenated sequences. All of the nuclear sequences represented phased haplotypes, which were included for all subsequent analyses. Cheirogaleus major was used as an outgroup for all phylogenetic and demographic analyses (see Additional file 1: Table S1). Although Cheirogaleus can be considered a relatively distant outgroup, we were interested in adding a temporal component to our analyses (see below) and sought to maintain consistency with taxon sampling across analyses. All sequence manipulation was performed in Geneious v. 6.1.5 . Multiple sequence alignments were conducted using MAFFT v. 7.017  within Geneious.
Sequence diversity and effective population sizes
We first tested each marker for signs of recombination using the program RDP v. 4.0 . Each marker was tested for recombination events using the GENECONV , MaxChi , and RDP methods . Default settings were used for all analyses. Because no signs of recombination were detected, all subsequent analyses utilized the entire read of each locus. General sequence diversity statistics for each locus and species including number of haplotypes, haplotype diversity, number of segregating sites, and nucleotide diversity were calculated using DnaSP v. 5.10.1 .
Using the species tree topology ((M. murinus, M. griseorufus), Cheirogaleus major) prescribed by both nuclear and mitochondrial phylogenies [43, 44, 53], we used the program Bayesian Phylogenetics and Phylogeography (BPP v. 2.2; ) to estimate coalescent-scaled population sizes (θ = 4Neμ) and time of divergence (τ = μt). This method accommodates the species tree divergence patterns as well as gene tree lineage sorting processes. A gamma prior G(2,1000), with mean 2/1000 = 0.002 was used for the population size parameters (θ values). The age of the root in the species tree (τ0) was also assigned a gamma prior G(2,1000), while the other divergence time parameters were assigned the Dirichlet prior (; Equation 2). Locus-specific mutation rates were allowed to vary, and we specified a heredity multiplier value G(4,5) to account for the combined mtDNA and nDNA data. We also performed analyses using larger priors for both τ0 and θ (G(2,100)) to evaluate the sensitivity of our results to the choice of prior. For this study, the parameters of interest were θM, θG, θMG, and τMG, where M = M. murinus, G = M. griseorufus, and MG = the MRCA of M. murinus and M. griseorufus. The default numbers of generations and sampling intervals were used for all analyses. Due to computational issues with the full data set, all BPP analyses were implemented on a reduced data set consisting of approximately 30 sequences per locus per Microcebus species and one to six Cheirogaleus sequences per locus. Sequences were randomly sampled to encompass a broad geographic range for each species. Previous simulation-based studies have shown that similar sample sizes can be sufficient to infer speciation processes . Additional analyses on further reduced data sets (e.g. 10 sequences per locus per species) yielded similar results. All analyses were run at least twice to check for consistency among runs.
To further test between the alternative speciation hypotheses, and to determine if the relatively large geographic range of M. murinus was due to population expansion following the split with M. griseorufus, we tested for signs of demographic change through time for both species. Under a peripatric speciation model, an ancestral species with a large range similar to that of modern M. murinus would have diverged into M. murinus and M. griseorufus populations, with the former experiencing little to no demographic change associated with the speciation event and the latter experiencing a significant population bottleneck event. Conversely, under the allopatric speciation model, a narrowly-distributed ancestral species would have diverged into M. murinus and M. griseorufus populations, with the former subsequently experiencing a significant population growth event and the latter experiencing a relatively smaller but possibly detectable demographic change after the speciation event.
To test for signs of demographic expansion or contraction we implemented multilocus Bayesian methods. Specifically, we constructed extended Bayesian skyline plots (EBSPs; ) for both species using BEAST v. 1.7.5 . Studies have shown that doubling the number of independent loci can reduce error and 95% credible intervals in demographic reconstruction by √2 . Thus, compared to single locus estimates, our multilocus data provided a powerful approach for estimating demographic trends in mouse lemurs. First, we used jModeltest v. 2.1.4 [57, 58] to calculate model likelihood scores for each locus and to estimate optimal models using BIC (Table 1). Because mtDNA is inherited as a single linked unit, and to minimize the computational burden for BEAST, the concatenated mtDNA data were treated as a single partition. The three nuclear loci were each specified as a separate partition for model fitting. We tested the likelihood of 24 commonly used models in BEAST.
We were also interested in adding a temporal component to the demographic analyses. Although there is no fossil record for lemurs, recent studies have utilized fossil information from more distantly related groups to date divergence times within the lemurs (e.g. [53, 59]). Using multiple calibration points outside the clade and relaxed clock methods, these authors estimated the split between Cheirogaleus and Microcebus to be approximately 25 million years ago (Ma). We first used this information to estimate the substitution rate for each of the test loci in BEAST using a reduced data set of both species and Cheirogaleus (see Additional file 1: Table S1). Best-fitting models were calculated and used for all divergence dating and rate estimation. Because alignments contained representative alleles from both multiple species and multiple individuals within species, we used *BEAST  to estimate the posterior distribution of substitution rates. We defined three species (Cheirogaleus major, M. griseorufus, M. murinus) and grouped alleles accordingly. For each analysis, the root node of the species tree was calibrated with a normal distribution around a mean of 25 Ma and standard deviation of 5 Ma, which encompassed the 95% HPD estimates from previous studies . To increase computational efficiency, we ran four independent *BEAST analyses (by locus) to estimate the posterior distribution of rates. Both strict and relaxed clock (lognormal; ) models were tested. All analyses were run between 10–50 million generations; sampling was chosen at intervals to utilize 10,000 draws from the posterior. Following analyses, the program Tracer v1.4  was used to examine effective sample size (ESS) values (target > 200) and examine the posterior distribution of relevant parameters. We used the 95% HPD of the substitution rate for each locus as a uniform prior to add a temporal component to all EBSP analyses. All EBSP analyses used a strict clock. Operators were modified according to author recommendations and analyses were run for 50 million generations (M. griseorufus) or 200 million generations (M. murinus) to obtain adequate ESS values. All EBSP and *BEAST analyses were implemented via the Duke Shared Cluster Resource (DSCR).
Under a scenario of peripatric speciation, we expected to detect differing signals in the degree of incomplete lineage sorting in gene trees for M. griseorufus and M. murinus. For example, some coalescence times within M. murinus would predate speciation, whereas M. griseorufus would show a higher degree of reciprocal monophyly among different loci due to its smaller effective population size, particularly during the bottleneck that was hypothesized to be associated with the peripatric speciation event. Conversely, under an allopatric speciation model we would expect the degree of incomplete lineage sorting among gene trees to be similar for both species. Thus, we performed maximum likelihood (ML) phylogenetic analysis of each locus using RAxML v. 7.6.0 . Because all mitochondrial genes are linked, we performed a single ML analysis for the concatenated mtDNA loci. For each gene we ran a full ML analysis followed by rapid bootstrapping  using the autoMRE bootstopping criterion.
Sequence diversity and effective population sizes
In general, sequence diversity characteristics showed moderate values for each gene for both species (Table 2). Although average values were slightly higher, M. murinus did not exhibit consistently higher diversity values than M. griseorufus based on haplotype diversity, nucleotide diversity, or average number of nucleotide differences. Multiple runs of BPP gave similar results indicating adequate sampling of the posterior. ESS values were also high for all parameters (Additional file 2: Table S2, Additional file 3: Table S3). Using the prior G(2,1000) for θ and τ0, mean effective population size was slightly greater for M. murinus (θM = 0.0099) versus M. griseorufus (θG = 0.0060). However, 95% HPDs for the two species almost completely overlapped (Figure 3). The MRCA had a significantly smaller mean population size than either species (θMG = 0.0032; 95% HPD 0.0011–0.0054). BPP results using a gamma prior of G(2,100) for θ and τ0 resulted in slightly larger parameter estimates (Additional file 3: Table S3). With these priors, the mean effective population size of the MRCA was intermediate between M. murinus and M. griseorufus. (θMG = 0.0208; 95% HPD 0.0087–0.0351; θM = 0.0316; 95% HPD 0.0183–0.0470; θG = 0.0185; 95% HPD 0.01026–0.02789). However, confidence intervals again overlapped substantially. In all analyses, effective population size estimates of the MRCA were more similar to M. griseorufus than they were to M. murinus.
For all *BEAST analyses, the 95% HPD for the coefficient of variation parameter included zero for all relaxed clock analyses, indicating that a strict clock was sufficient to explain the data. Mean estimated rates of nucleotide substitution (substitutions per site per million years) and 95% HPDs for each locus were as follows: mtDNA = 0.0132 [0.0060862–0.0222]; ENOL = 0.0025731 [0.00037366–0.0046715]; FIB = 0.0010495 [0.00027744–0.0019127]; VWF = 0.0023935 [0.00037916–0.0045323]). EBSP results for M. murinus indicated an increase in effective population size at about 1 Ma, with a subsequent rapid decline in size starting approximately 160 Ka and continuing to the present (Figure 4A,B). The mean number of population size changes throughout the history of M. murinus was estimated as 2.4 (95% HPD 2–4). Conversely, results for M. griseorufus showed signs of relatively constant population size through time (Figure 4C,D) with an estimated mean of 0.59 size changes (95% HPD 0–2).
Maximum likelihood phylogenetic analyses of the mtDNA data revealed reciprocal monophyly for both species with strong bootstrap support (Additional file 4: Figure S1). Gene copies within Microcebus griseorufus were also monophyletic based on ML analysis of the ENOL locus (bootstrap support = 81; Additional file 5: Figure S2). Conversely, none of the nuclear markers showed reciprocal monophyly for M. murinus. Moderate geographic structure among populations was detected from the mtDNA analysis, particularly for M. murinus (Additional file 4: Figure S1). However, nuclear gene tree analyses suggested little to no signal of geographic population structure within species as haplotypes were shared among many localities (Additional file 5: Figure S2, Additional file 6: Figure S3, Additional file 7: Figure S4). Species tree analyses suggested that divergence of M. griseorufus and M. murinus occurred approximately 3–6 Ma.
Our motivation for this study was to differentiate between competing models of speciation that can potentially explain the historical divergence between a sister species pair of mouse lemurs. These primates are of particular interest for such a study given their cryptic species diversity, highly threatened status, and their endemic distribution in Madagascar, one of Earth’s hottest biodiversity hotspots. Given the difference in the sizes of their geographic distributions (Figure 1A), we initially hypothesized that M. griseorufus diverged peripatrically from a geographically widespread common ancestor (Figure 2A). The observation that two sister species occur with different yet overlapping ranges in Madagascar raises obvious questions regarding the driving mechanisms behind their divergence. Given that one species, M. griseorufus, shows a more limited though contiguous range with M. murinus is on its face entirely congruent with a peripatric model of speciation. To test this hypothesis, we formed a set of predictions that rely on a series of historical demographic variables including comparisons of effective population size in both the modern species and their common ancestor, as well as patterns of population size change in the history of the focal species. The majority of our results reject a model of peripatric speciation in favor of a model of allopatric divergence with subsequent range expansion for M. murinus.
Contemporary and ancestral Ne
Our BPP results using both large and small priors for divergence times and population sizes strongly suggest that contemporary Ne values are similar for both species. For example, although mean estimates of Ne for M. murinus are slightly larger than those for M. griseorufus, 95% HPDs overlapped significantly. In analyses using the priors G(2,1000), Ne estimates of the MRCA are substantially smaller than for either contemporary species, although using larger priors results in Ne estimates of the MRCA that are intermediate to the two contemporary species. However, regardless of the prior used, Ne estimates for the MRCA were more similar to M. griseorufus than to the more widespread M. murinus. These results support a model in which the MRCA was a species with a relatively small geographic distribution likely situated within southern Madagascar. These conclusions are congruent with our phylogenetic results and with estimates of genetic diversity for both species. Although genetic diversity is slightly higher in M. murinus, this may be an artifact of both our sampling regime and/or recent demographic trends for this species (see below). Under a peripatric scenario we would expect to find substantially larger diversity values in M. murinus, a pattern that was not recovered with any marker. The reciprocal monophyly and population structure we observed with the mtDNA and incomplete lineage sorting in nuclear markers is also congruent with other recent studies of these species (e.g. [43, 44]).
These results strongly support the allopatric rather than the peripatric speciation model. Of particular importance for this conclusion is the estimation of Ne for the MRCA of M. murinus and M. griseorufus. The model employed in BPP allows for the combined analysis of multiple genetic markers in a coalescent framework, a necessary approach as individual loci may suffer from rate heterogeneities and idiosyncratic gene genealogies [65, 66]. Because a single genetic locus provides a limited and highly stochastic perspective, we find that examining four independent loci is effective for estimates when combined, but may remain insufficient for reconciling historical demographic processes when analyzed individually given the limited power of single locus analyses [5, 53]. A rejection of a peripatric model for mouse lemurs is also similar to a recent study on mantellid frogs based on range overlap analysis, where the authors found that range size differences among sister species increased with evolutionary age . For example, many sister species of frogs were composed of microendemics encompassing similarly-sized geographic ranges. Under a peripatric scenario, range asymmetry would be high in younger species.
Our results regarding population size changes provide further support for the allopatric model. Although only four independent loci are used to infer demographic changes, four loci are predicted to reduce the error by one-half as compared to single-locus estimates . EBSP results suggest demographic expansion for M. murinus in the Quaternary around 1 Ma. Demographic expansion during the Quaternary has also been documented previously for populations of M. murinus in northwest Madagascar  as well as for Malagasy rodents . However, our results differ from previous single-locus studies of M. murinus that have utilized different analytical methods and found evidence for more recent episodes of expansion during the Pleistocene and Holocene . Evidence from palynological records indicates that the Pleistocene climate and vegetation of Madagascar, like most of the world, was quite different than the climate of today . For example, during the Last Glacial Maximum (LGM) ~40–20,000 years ago, humid forest was likely restricted to isolated refugia scattered throughout the island, whereas dry, xeric vegetation was allowed to expand. Indeed, a wide body of evidence is available that suggests that a large portion of Madagascar experienced substantially drier conditions during the LGM than the present (see ). Both M. murinus and M. griseorufus are common in dry environments and both were likely affected by climate conditions associated with Quaternary Madagascar.
Following the demographic expansion of M. murinus during the Quaternary, our results suggest recent and substantial population decline of this species beginning approximately 160 Ka and continuing to the present. It has been proposed that vegetation shifts associated with Pleistocene climate change were more substantial in western dry forest versus the arid spiny forest to the south , which may partly explain why no evidence of recent population decline is indicated for M. griseorufus. Furthermore, evidence suggests that humans first colonized Madagascar ~2,000 years ago and subsequently had a rapid and profound impact on the native biota and their habitats [70, 71]. A variety of hypotheses have been put forth to explain the decline of Malagasy flora and fauna subsequent to human colonization, including increased frequency of fire , drought , hunting , invasive species , disease , and synergistic anthropogenic influences . Regardless of the exact mechanism(s), it is highly probable that the recent and rapid population decline inferred from our data for M. murinus has been exacerbated by subsequent anthropogenic influences beginning around 2,000 years ago. The relatively constant population size of M. griseorufus suggests that human impacts and habitat fragmentation throughout the southern spiny forests may have been less severe than impacts throughout western dry forests.
Alternative speciation hypotheses
Our results indicate that the large geographic range of M. murinus seems to be a uniquely derived feature of this species. This phenomenon begs investigation. Even so, the evidence presented in this study fails to specifically explain the mechanism either promoting the geographic range expansion of M. murinus or limiting the range of M. griseorufus. The initial divergence between the two species may have resulted from a geographically-based vicariant event or an ecological niche separation. There is no discernible extrinsic barrier separating the two species, but there is evidence that suggests ecological segregation [39–41, 77, 78]. Although some studies suggest that M. murinus preferentially inhabits dry forest habitat in northwestern Madagascar , it could be argued that M. murinus is a generalist species as it is often found in both dry deciduous forest and wet, gallery forest habitats in the southeast [38, 39]. These habitat associations contrast with M. griseorufus, which is more common in xeric, spiny forest [21, 39], though recently M. griseorufus has been shown to inhabit both spiny and gallery forests at the Beza Mahafaly Private Reserve . Therefore, while these habitat preferences may have been important in the original subdivision of the ancestral species, there is no direct evidence supporting the hypothesis that this drove their initial divergence. The use of next-generation DNA sequencing methods, in particular RAD-Seq or whole transcriptomes, may be useful for investigating and quantifying genomic islands of divergence and adaptation in this system. Additionally, multilocus demographic methods such as those used here should be combined with future projections of species distribution models to better inform conservation practices. This is of utmost importance for highly threatened taxa inhabiting areas that are experiencing high rates of habitat loss, such as the case with Madagascar’s lemur fauna [80, 81].
We estimated historical demographic parameters in a multilocus coalescent framework to test the predictions associated with two models of speciation that may have driven the divergence of M. griseorufus and M. murinus. The majority of our results reject the hypothesis of peripatric speciation. Our results instead favor a model of allopatric divergence from a range-restricted common ancestor in southwestern Madagascar, with subsequent range expansions for M. murinus. Whether due to ecological constraint or interspecific competition, M. griseorufus is presently restricted to the arid spiny forest in the south, whereas M. murinus has successfully expanded throughout much of western Madagascar and limited areas to the southeast. The methods used here can be easily applied to address similar evolutionary questions in other systems to help elucidate the geographic context of divergence and speciation. In turn, these approaches can help guide conservation priorities when synthesized with complex geospatial methods and species distribution models.
CB is presently a postdoctoral associate in the laboratory of Anne Yoder at Duke University. His research interests include combining next-generation DNA sequencing data with spatially explicit modeling to understand the contemporary and historical processes responsible for shaping patterns of genomic variation in natural populations. He is also interested in using next-generation sequencing of environmental samples to quantify patterns of biodiversity.
Effective population size
Most recent common ancestor
Alpha enolase intron
Alpha fibrinogen intron
Von Willebrand factor intron
Cytochrome c oxidase subunit II
Bayesian Phylogenetics and Phylogeography
Million years ago
Effective sample size
Duke Shared Cluster Resource
Highest posterior density
Effective population size of Microcebus murinus
Effective population size of Microcebus griseorufus
Effective population size of the most recent common ancestor
Divergence time of the most recent common ancestor
Thousand years ago
Last Glacial Maximum.
Hey J, Nielsen R: Multilocus methods for estimating population sizes, migration rates and divergence time, with applications to the divergence of Drosophila pseudoobscura and D. persimilis. Genetics. 2004, 167 (2): 747-760. 10.1534/genetics.103.024182.
Hey J: Isolation with migration models for more than two populations. Mol Biol Evol. 2010, 27 (4): 905-920. 10.1093/molbev/msp296.
Hey J: The divergence of chimpanzee species and subspecies as revealed in multipopulation isolation-with-migration analyses. Mol Biol Evol. 2010, 27 (4): 921-933. 10.1093/molbev/msp298.
Edwards S, Beerli P: Perspective: gene divergence, population divergence, and the variance in coalescence time in phylogeographic studies. Evolution. 2000, 54 (6): 1839-1854.
Felsenstein J: Estimating effective population size from samples of sequences: inefficiency of pairwise and segregating sites as compared to phylogenetic estimates. Genet Res. 1992, 59 (02): 139-147. 10.1017/S0016672300030354.
Kuhner MK, Yamato J, Felsenstein J: Estimating effective population size and mutation rate from sequence data using metropolis-hastings sampling. Genet. 1995, 140 (4): 1421-1430.
Beerli P, Felsenstein J: Maximum-likelihood estimation of migration rates and effective population numbers in two populations using a coalescent approach. Genet. 1999, 152 (2): 763-773.
Rosenberg NA, Feldman MW: The relationship between coalescence times and population divergence times. Modern Developments in Theoretical Population Genetics. Edited by: Slatkin M, Veuille M. 2002, Oxford: Oxford University Press, 130-164.
Yang Z, Rannala B: Bayesian species delimitation using multilocus sequence data. Proc Natl Acad Sci. 2010, 107 (20): 9264-9269. 10.1073/pnas.0913022107.
Knowles LL, Maddison WP: Statistical phylogeography. Mol Ecol. 2002, 11 (12): 2623-2635. 10.1046/j.1365-294X.2002.01637.x.
Jennings WB, Edwards SV: Speciational history of Australian grass finches (Poephila) inferred from thirty gene trees. Evol. 2005, 59 (9): 2033-2047.
Russell A, Ranivo J, Palkovacs E, Goodman S, Yoder AD: Working at the interface of phylogenetics and population genetics: a biogeographical analysis of Triaenops spp.(Chiroptera: Hipposideridae). Mol Ecol. 2007, 16 (4): 839-851. 10.1111/j.1365-294X.2007.03192.x.
Rundle HD, Nosil P: Ecological speciation. Ecol Lett. 2005, 8 (3): 336-352. 10.1111/j.1461-0248.2004.00715.x.
Nosil P, Gompert Z, Farkas TE, Comeault AA, Feder JL, Buerkle CA, Parchman TL: Genomic consequences of multiple speciation processes in a stick insect. Proc Biol Sci. 2012, 279 (1749): 5058-5065. 10.1098/rspb.2012.0813.
Nosil P, Funk DJ, Ortiz-Barrientos D: Divergent selection and heterogeneous genomic divergence. Mol Ecol. 2009, 18 (3): 375-402. 10.1111/j.1365-294X.2008.03946.x.
Malay MCMD, Paulay G: Peripatric speciation drives diversification and distributional pattern of reef hermit crabs (Decapoda: Diogenidae: Calcinus). Evol. 2010, 64 (3): 634-662. 10.1111/j.1558-5646.2009.00848.x.
Bernardi G, Bucciarelli G, Costagliola D, Robertson DR, Heiser J: Evolution of coral reef fish Thalassoma spp.(Labridae). 1. Molecular phylogeny and biogeography. Mar Biol. 2004, 144 (2): 369-375. 10.1007/s00227-003-1199-0.
Maddison WP: Gene trees in species trees. Syst Biol. 1997, 46 (3): 523-536. 10.1093/sysbio/46.3.523.
Moritz C, Cicero C: DNA barcoding: promise and pitfalls. PLoS Biol. 2004, 2: 1529-1531.
Blair C, Murphy RW: Recent trends in molecular phylogenetic analysis: where to next?. J Hered. 2011, 102 (1): 130-138. 10.1093/jhered/esq092.
Rasoloarison RM, Goodman SM, Ganzhorn JU: Taxonomic revision of mouse lemurs (Microcebus) in the western portions of Madagascar. Int J Primatol. 2000, 21 (6): 963-1019. 10.1023/A:1005511129475.
Rasoloarison RM, Weisrock DW, Yoder AD, Rakotondravony D, Kappeler PM: Two new species of mouse lemurs (Cheirogaleidae: Microcebus) from eastern Madagascar. Int J Primatol. 2013, 34 (3): 455-469. 10.1007/s10764-013-9672-1.
Kappeler P, Rasoloarison R, Razafimanantsoa L, Walter L, Roos C: Morphology, behaviour and molecular evolution of giant mouse lemurs (Mirza spp.) Gray, 1870, with description of a new species. Primate Report. 2005, 71: 3-
Louis EE, Coles MS, Andriantompohavana R, Sommer JA, Engberg SE, Zaonarivelo JR, Mayor MI, Brenneman RA: Revision of the mouse lemurs (Microcebus) of eastern Madagascar. Int J Primatol. 2006, 27 (2): 347-389. 10.1007/s10764-006-9036-1.
Radespiel U, Ratsimbazafy JH, Rasoloharijaona S, Raveloson H, Andriaholinirina N, Rakotondravony R, Randrianarison RM, Randrianambinina B: First indications of a highland specialist among mouse lemurs (Microcebus spp.) and evidence for a new mouse lemur species from eastern Madagascar. Primates. 2012, 53 (2): 157-170. 10.1007/s10329-011-0290-2.
Martin RD: Review lecture: Adaptive radiation and behaviour of the Malagasy lemurs. Phil Trans Roy Soc Lond B Biol Sci. 1972, 264: 295-352. 10.1098/rstb.1972.0013.
Tattersall I: Primates of Madagascar. 1982, New York: Columbia University Press, 52
Pearson RG, Raxworthy CJ, Nakamura M, Townsend Peterson A: Predicting species distributions from small numbers of occurrence records: a test case using cryptic geckos in Madagascar. J Biogeogr. 2007, 34 (1): 102-117.
Russell AL, Goodman SM, Fiorentino I, Yoder AD: Population genetic analysis of Myzopoda (Chiroptera: Myzopodidae) in Madagascar. J Mammal. 2008, 89 (1): 209-221. 10.1644/07-MAMM-A-044.1.
Wilmé L, Goodman SM, Ganzhorn JU: Biogeographic evolution of Madagascar’s microendemic biota. Sci. 2006, 312 (5776): 1063-1065. 10.1126/science.1122806.
Townsend TM, Vieites DR, Glaw F, Vences M: Testing species-level diversification hypotheses in Madagascar: the case of microendemic Brookesia leaf chameleons. Syst Biol. 2009, 58 (6): 641-656. 10.1093/sysbio/syp073.
Vences M, Wollenberg KC, Vieites DR, Lees DC: Madagascar as a model region of species diversification. Trends Ecol Evol. 2009, 24 (8): 456-465. 10.1016/j.tree.2009.03.011.
Mittermeier RA, Tattersall L, Konstant WR, Meyers DM, Mast RB: Lemurs of Madagascar. 1994, Washington D.C: Conservation International, 1994
Yoder AD, Rasoloarison RM, Goodman SM, Irwin JA, Atsalis S, Ravosa MJ, Ganzhorn JU: Remarkable species diversity in Malagasy mouse lemurs (primates, Microcebus). Proc Natl Acad Sci. 2000, 97 (21): 11325-11330. 10.1073/pnas.200121897.
Yoder AD, Heckman KL: Mouse lemur phylogeography revises a model of ecogeographic constraint in Madagascar. Primate Biogeography. 2006, New York: Springer, 255-268.
Olivieri G, Zimmermann E, Randrianambinina B, Rasoloharijaona S, Rakotondravony D, Guschanski K, Radespiel U: The ever-increasing diversity in mouse lemurs: three new species in north and northwestern Madagascar. Mol Phylogenet Evol. 2007, 43 (1): 309-327. 10.1016/j.ympev.2006.10.026.
Schmid J, Kappeler P: Sympatric mouse lemurs (Microcebus spp.) in western Madagascar. Folia Primatol. 1994, 63 (3): 162-170. 10.1159/000156812.
Zimmermann E, Cepok S, Rakotoarison N, Radespiel U: Sympatric mouse lemurs in north-west Madagascar: a new rufous mouse lemur species (Microcebus ravelobensis). Folia Primatol. 1998, 69 (2): 106-114. 10.1159/000021571.
Yoder AD, Burns MM, Génin F: Molecular evidence of reproductive isolation in sympatric sibling species of mouse lemurs. Int J Primatol. 2002, 23 (6): 1335-1343. 10.1023/A:1021187106641.
Gligor M, Ganzhorn J, Rakotondravony D, Ramilijaona O, Razafimahatratra E, Zischler H, Hapke A: Hybridization between mouse lemurs in an ecological transition zone in southern Madagascar. Mol Ecol. 2009, 18 (3): 520-533. 10.1111/j.1365-294X.2008.04040.x.
Hapke A, Gligor M, Rakotondranary SJ, Rosenkranz D, Zupke O: Hybridization of mouse lemurs: different patterns under different ecological conditions. BMC Evol Biol. 2011, 11 (1): 297-10.1186/1471-2148-11-297.
Martin RD: Origins, diversity and relationships of lemurs. Int J Primatol. 2000, 21 (6): 1021-1049. 10.1023/A:1005563113546.
Heckman KL, Mariani CL, Rasoloarison R, Yoder AD: Multiple nuclear loci reveal patterns of incomplete lineage sorting and complex species history within western mouse lemurs (Microcebus). Mol Phylogenet Evol. 2007, 43 (2): 353-367. 10.1016/j.ympev.2007.03.005.
Weisrock DW, Rasoloarison RM, Fiorentino I, Ralison JM, Goodman SM, Kappeler PM, Yoder AD: Delimiting species without nuclear monophyly in Madagascar’s mouse lemurs. PLoS One. 2010, 5 (3): e9883-10.1371/journal.pone.0009883.
Heckman KL, Rasoazanabary E, Machlin E, Godfrey LR, Yoder AD: Incongruence between genetic and morphological diversity in Microcebus griseorufus of Beza Mahafaly. BMC Evol Biol. 2006, 6 (1): 98-10.1186/1471-2148-6-98.
Kearse M, Moir R, Wilson A, Stones-Havas S, Cheung M, Sturrock S, Buxton S, Cooper A, Markowitz S, Duran C: Geneious basic: an integrated and extendable desktop software platform for the organization and analysis of sequence data. Bioinformatics. 2012, 28 (12): 1647-1649. 10.1093/bioinformatics/bts199.
Katoh K, Misawa K, Kuma K, Miyata T: MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Res. 2002, 30 (14): 3059-3066. 10.1093/nar/gkf436.
Martin DP, Lemey P, Lott M, Moulton V, Posada D, Lefeuvre P: RDP3: a flexible and fast computer program for analyzing recombination. Bioinformatics. 2010, 26 (19): 2462-2463. 10.1093/bioinformatics/btq467.
Padidam M, Sawyer S, Fauquet CM: Possible emergence of new geminiviruses by frequent recombination. Virol. 1999, 265 (2): 218-225. 10.1006/viro.1999.0056.
Smith JM: Analyzing the mosaic structure of genes. J Mol Evol. 1992, 34 (2): 126-129.
Martin D, Rybicki E: RDP: detection of recombination amongst aligned sequences. Bioinformatics. 2000, 16 (6): 562-563. 10.1093/bioinformatics/16.6.562.
Librado P, Rozas J: DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009, 25 (11): 1451-1452. 10.1093/bioinformatics/btp187.
Horvath JE, Weisrock DW, Embry SL, Fiorentino I, Balhoff JP, Kappeler P, Wray GA, Willard HF, Yoder AD: Development and application of a phylogenomic toolkit: resolving the evolutionary history of Madagascar’s lemurs. Genome Res. 2008, 18 (3): 489-499. 10.1101/gr.7265208.
Knowles LL, Carstens BC: Delimiting species without monophyletic gene trees. Syst Biol. 2007, 56 (6): 887-895. 10.1080/10635150701701091.
Heled J, Drummond A: Bayesian inference of population size history from multiple loci. BMC Evol Biol. 2008, 8 (1): 289-10.1186/1471-2148-8-289.
Drummond AJ, Suchard MA, Xie D, Rambaut A: Bayesian phylogenetics with BEAUti and the BEAST 1.7. Mol Biol Evol. 2012, 29 (8): 1969-1973. 10.1093/molbev/mss075.
Guindon S, Gascuel O: A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood. Syst Biol. 2003, 52 (5): 696-704. 10.1080/10635150390235520.
Darriba D, Taboada GL, Doallo R, Posada D: jModelTest 2: more models, new heuristics and parallel computing. Nature Methods. 2012, 9 (8): 772-
Yoder AD, Yang Z: Divergence dates for Malagasy lemurs estimated from multiple gene loci: geological and evolutionary context. Mol Ecol. 2004, 13 (4): 757-773. 10.1046/j.1365-294X.2004.02106.x.
Heled J, Drummond AJ: Bayesian inference of species trees from multilocus data. Mol Biol Evol. 2010, 27 (3): 570-580. 10.1093/molbev/msp274.
Drummond AJ, Ho SY, Phillips MJ, Rambaut A: Relaxed phylogenetics and dating with confidence. PLoS Biol. 2006, 4 (5): e88-10.1371/journal.pbio.0040088.
Rambaut A, Drummond A: Tracer v1. 4. 2007, http://beast.bio.ed.ac.uk/Tracer,
Stamatakis A: RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics. 2006, 22 (21): 2688-2690. 10.1093/bioinformatics/btl446.
Stamatakis A, Hoover P, Rougemont J: A rapid bootstrap algorithm for the RAxML web servers. Syst Biol. 2008, 57 (5): 758-771. 10.1080/10635150802429642.
Yang Z: On the estimation of ancestral population sizes of modern humans. Genet Res. 1997, 69 (2): 111-116. 10.1017/S001667239700270X.
Chen F-C, Li W-H: Genomic divergences between humans and other hominoids and the effective population size of the common ancestor of humans and chimpanzees. Am J Hum Genet. 2001, 68 (2): 444-456. 10.1086/318206.
Wollenberg KC, Vieites DR, Glaw F, Vences M: Speciation in little: the role of range and body size in the diversification of Malagasy mantellid frogs. BMC Evol Biol. 2011, 11 (1): 217-10.1186/1471-2148-11-217.
Schneider N, Chikhi L, Currat M, Radespiel U: Signals of recent spatial expansions in the grey mouse lemur (Microcebus murinus). BMC Evol Biol. 2010, 10 (1): 105-10.1186/1471-2148-10-105.
Shi JJ, Chan LM, Rakotomalala Z, Heilman AM, Goodman SM, Yoder AD: Latitude drives diversification in Madagascar’s endemic dry forest rodent Eliurus myoxinus (subfamily Nesomyinae). Biol J Linn Soc. 2013, 110 (3): 500-517. 10.1111/bij.12143.
Burney DA, Burney LP, Godfrey LR, Jungers WL, Goodman SM, Wright HT, Jull A: A chronology for late prehistoric Madagascar. J Hum Evol. 2004, 47 (1): 25-63.
Perez V, Burney D, Godfrey L, Nowak-Kemp M: Box 4. Butchered sloth lemurs. Evol Anthropol. 2003, 12 (6): 260-
Humbert H: La destruction d’une flore insulaire par le feu: principaux aspects de la végétation à Madagascar, documents photographiques et notices, vol. 5. 1927, Tananarive, Madagascar: Impr. Moderne de l’Emyrne, G. Pitot et Cie
Mahé J, Sourdat M: Sur l’extinction des vertébrés subfossiles et l’aridification du climat dans le Sud-ouest de Madagascar. Bulletin de la Société de Géologie de France. 1972, 14: 295-309.
Martin PS: Prehistoric overkill: the global model. Quaternary Extinctions: a Prehistoric Revolution. Edited by: Martin PS, Klein RG. 1984, Tucson: University of Arizona Press, 354-403.
Dewar RE: Extinctions in Madagascar: the loss of the subfossil fauna. Quaternary Extinctions: a Prehistoric Revolution. Edited by: Martin PS, Klein RG. 1984, Tucson: University of Arizona Press, 574-593.
Burney DA: Rates, patterns, and processes of landscape transformation and extinction in Madagascar. Extinctions in Near Time: Causes, Contexts, and Consequences. Edited by: MacPhee RDE. 1999, New York: Springer Science+Business Media, 145-164.
Rakotondranary SJ, Hapke A, Ganzhorn JU: Distribution and morphological variation of Microcebus spp. along an environmental gradient in southeastern Madagascar. Int J Primatol. 2011, 32 (5): 1037-1057. 10.1007/s10764-011-9521-z.
Rakotondranary SJ, Ganzhorn JU: Habitat separation of sympatric Microcebus spp. in the dry spiny forest of south-eastern Madagascar. Folia Primatol. 2012, 82 (4–5): 212-223.
Rakotondravony R, Radespiel U: Varying patterns of coexistence of two mouse lemur species (Microcebus ravelobensis and M. murinus) in a heterogeneous landscape. Am J Primatol. 2009, 71 (11): 928-938. 10.1002/ajp.20732.
Harper GJ, Steininger MK, Tucker CJ, Juhn D, Hawkins F: Fifty years of deforestation and forest fragmentation in Madagascar. Environ Conservat. 2007, 34 (4): 325-
Irwin MT, Wright PC, Birkinshaw C, Fisher BL, Gardner CJ, Glos J, Goodman SM, Loiselle P, Rabeson P, Raharison J-L: Patterns of species change in anthropogenically disturbed forests of Madagascar. Biol Conservat. 2010, 143 (10): 2351-2362. 10.1016/j.biocon.2010.01.023.
This research was supported by NSF award DEB-0516276 and a Biodiversity Leadership Award from the Bay & Paul Foundations to A.D.Y. The manuscript was significantly improved thanks to comments from M. Noor and discussion with J. Wall.
The data sets supporting the results of this article are available in the Dryad repository, [Provisional DOI: doi:10.5061/dryad.9v43h].
The authors declare that they have no competing interests.
CB, KH, AR, ADY designed the study. CB collected and analyzed the data. CB, AR, ADY wrote the paper. All authors read an approved the final manuscript.
Electronic supplementary material
Additional file 1: Table S1: Database of all individuals and haplotypes used for this study. The last four columns indicate which haplotypes were used for each analysis. (XLSX 91 KB)
Additional file 2: Table S2: Results from BPP analysis using priors G(2,1000) for both τ0 and θ. Parameters were estimated on a fixed species tree (Cheirogaleus major (M. murinus, M. griseorufus)). (XLSX 39 KB)
Additional file 3: Table S3: Results from BPP analysis using priors G(2,100) for both τ0 and θ. Parameters were estimated on a fixed species tree (Cheirogaleus major (M. murinus, M. griseorufus)). (XLSX 26 KB)
Additional file 4: Figure S1: Maximum likelihood mtDNA gene tree (concatenated cytochrome b and cytochrome c oxidase II) for all Microcebus griseorufus (purple) and M. murinus (blue) sequences used for this study. Values at nodes represent bootstrap support values >50 calculated using the autoMRE function in RAxML. For ease of visualization the outgroup taxon (Cheirogaleus major) is removed. (PDF 274 KB)
Additional file 5: Figure S2: Maximum likelihood alpha enolase gene tree for all Microcebus griseorufus (purple) and M. murinus (blue) sequences used for this study. Values at nodes represent bootstrap support values >50 calculated using the autoMRE function in RAxML. For ease of visualization the outgroup taxon (Cheirogaleus major) is removed. (PDF 346 KB)
Additional file 6: Figure S3: Maximum likelihood alpha fibrinogen gene tree for all Microcebus griseorufus (purple) and M. murinus (blue) sequences used for this study. Values at nodes represent bootstrap support values >50 calculated using the autoMRE function in RAxML. For ease of visualization the outgroup taxon (Cheirogaleus major) is removed. (PDF 457 KB)
Additional file 7: Figure S4: Maximum likelihood von Willebrand factor gene tree for all Microcebus griseorufus (purple) and M. murinus (blue) sequences used for this study. Values at nodes represent bootstrap support values >50 calculated using the autoMRE function in RAxML. For ease of visualization the outgroup taxon (Cheirogaleus major) is removed. (PDF 327 KB)
About this article
Cite this article
Blair, C., Heckman, K.L., Russell, A.L. et al. Multilocus coalescent analyses reveal the demographic history and speciation patterns of mouse lemur sister species. BMC Evol Biol 14, 57 (2014). https://doi.org/10.1186/1471-2148-14-57
- Coalescent methods
- Modes of speciation
- Historical demography
- Lemur evolution
- Peripatric speciation