Skip to main content
Download PDF

Impact of model assumptions on demographic inferences: the case study of two sympatric mouse lemurs in northwestern Madagascar

BMC Ecology and Evolution volume 21, Article number: 197 (2021) Cite this article

Abstract

Background

Quaternary climate fluctuations have been acknowledged as major drivers of the geographical distribution of the extraordinary biodiversity observed in tropical biomes, including Madagascar. The main existing framework for Pleistocene Malagasy diversification assumes that forest cover was strongly shaped by warmer Interglacials (leading to forest expansion) and by cooler and arid glacials (leading to forest contraction), but predictions derived from this scenario for forest-dwelling animals have rarely been tested with genomic datasets.

Results

We generated genomic data and applied three complementary demographic approaches (Stairway Plot, PSMC and IICR-simulations) to infer population size and connectivity changes for two forest-dependent primate species (Microcebus murinus and M. ravelobensis) in northwestern Madagascar. The analyses suggested major demographic changes in both species that could be interpreted in two ways, depending on underlying model assumptions (i.e., panmixia or population structure). Under panmixia, the two species exhibited larger population sizes across the Last Glacial Maximum (LGM) and towards the African Humid Period (AHP). This peak was followed by a population decline in M. ravelobensis until the present, while M. murinus may have experienced a second population expansion that was followed by a sharp decline starting 3000 years ago. In contrast, simulations under population structure suggested decreasing population connectivity between the Last Interglacial and the LGM for both species, but increased connectivity during the AHP exclusively for M. murinus.

Conclusion

Our study shows that closely related species may differ in their responses to climatic events. Assuming that Pleistocene climatic conditions in the lowlands were similar to those in the Malagasy highlands, some demographic dynamics would be better explained by changes in population connectivity than in population size. However, changes in connectivity alone cannot be easily reconciled with a founder effect that was shown for M. murinus during its colonization of the northwestern Madagascar in the late Pleistocene. To decide between the two alternative models, more knowledge about historic forest dynamics in lowland habitats is necessary. Altogether, our study stresses that demographic inferences strongly depend on the underlying model assumptions. Final conclusions should therefore be based on a comparative evaluation of multiple approaches.

Peer Review reports

Background

Marked Quaternary climatic oscillations have been largely acknowledged as a major driver of evolutionary and biogeographical patterns of species worldwide [1, 2]. The present-day distribution, genetic diversity patterns and demography of many temperate and tropical species have been shaped by historical warming–cooling cycles that forced species to retract and expand according to their ecological requirements [2,3,4,5,6,7,8]. Accordingly, there is an increasing interest in reconstructing the climate, biome and fire regimes of the late Quaternary. Most studies focus on the well-pronounced climate fluctuations that occurred during the last Interglacial–Glacial cycle, which included the Last Interglacial (LIG; ca. 132–112 kyr; kyr = thousand years) and the Last Glacial Maximum (LGM; ca. 26.5–19 kyr) [9].

There is a general consensus that the climate during the LGM was globally cooler than today, even if the magnitude of the cooling was not spatially uniform across the globe [1, 10, 11]. However, there is less general knowledge about the effects of the last glaciation across the tropics [1, 11,12,13], partly due to the scarcity of high-resolution paleoenvironmental records in the southern tropics. In contrast to these uncertainties about the last Interglacial–Glacial cycle, the so-called African Humid Period (AHP; ca. 15 to 5 kyr [14,15,16], but timing differed slightly across Africa) represents a well-established climatic event in various African regions. The AHP was characterized by a sudden increase in summer precipitation that was followed by an abrupt shift toward more arid conditions, strongly impacting the vegetation cover across continental Africa and Madagascar [14, 17, 18]. The demographic history of a species, e.g., the timing and extent of population expansions or bottleneck events, should indirectly mirror past environmental fluctuations and can provide information about species resilience to past [19] and possibly future climatic oscillations.

Madagascar is a natural evolutionary laboratory, allowing to investigate how past climatic changes shaped species demographic history. First, the island is characterized by exceptional levels of species richness and endemism [20]. Second, Madagascar has been isolated from other landmasses for over 80 million years (Mya) [21] and exhibits marked environmental gradients. These conditions promoted multiple adaptive radiations and resulted in many cases of micro-endemism that evolved in response to a particular set of local or regional environmental and climatic conditions [22]. Third, Madagascar was one of the last major landmasses on Earth settled by humans (e.g., [23,24,25,26,27] but see [24, 28, 29]), which enables us to control for the confounding effect of the anthropogenic impact on endemic species demography until recent times. It has long been assumed that the climate in Madagascar was generally cooler and more arid during glaciations, and that the extent of the forest cover dramatically contracted during these periods [30,31,32], likely resulting in high levels of specialization [33]. However, solid evidence for historical forest cover dynamics across different regions and habitat types on the island is still missing.

More than 90% of the Malagasy species and almost all lemurs live exclusively in forests and woodlands [34]. Among the lemuriforms, mouse lemurs (Microcebus spp.) provide a suitable model for demographic studies, because they have a very young age at first reproduction (ca. 8 months) [35] and therefore a relatively short generation time (ca. 2.5 years) [36]. Moreover, as forest-dwelling species they should be susceptible to vegetation shifts, and their Pleistocene demographic dynamics can therefore be expected to correspond to past environmental changes. While mouse lemurs are typically microendemic species, M. murinus is the only mouse lemur species with a large geographic distribution, inhabiting various forest habitats from southern to northwestern Madagascar [37]. Across its range, M. murinus co-occurs with five other locally restricted mouse lemur species, including M. ravelobensis in one northern part of its distribution [38, 39]. M. ravelobensis occurs exclusively in the so-called Inter-River-System Ia (IRS Ia; Fig. 1a) delimited by the Betsiboka and Mahajamba rivers [40, 41]. The two species have most likely undergone very different evolutionary trajectories. While M. ravelobensis is thought to have diverged allopatrically from its sister species M. bongolavensis also distributed in northwestern Madagascar [40], M. murinus likely diverged allopatrically from its sister species M. griseorufus at about 3–6 Mya in southwestern [38] and colonized northwestern Madagascar only during the Late Pleistocene [42].

Fig. 1
figure 1

Study area and sampling strategy. a Distribution range of M. murinus and M. ravelobensis across northwestern Madagascar and location of study sites. Individual capture locations of mouse lemurs in b Ravelobe and c Ankomakoma along four forest transects that were arranged in proximity to one lake each (in blue). The two forest sites were approximately 10 km apart. a The river shape files were provided by [141] and the outlines of the Ankarafantsika National Park were obtained from the Protected Planet Database [142]. b and c Individual coordinates can be found in Additional file 1: Table S8. M.mur = M. murinus; M.rav = M. ravelobensis; ANP = Ankarafantsika National Park

Full size image

The present study aims to investigate the impact of the Late Quaternary environmental changes on the demographic dynamics of these two mouse lemur species (M. murinus and M. ravelobensis) living in partial sympatry in the lowland forests of northwestern Madagascar in spite of their different phylogeographical background. Notwithstanding recent advances in analytical tools to reconstruct populations dynamics of non-model organisms (e.g., [43, 44]), previous studies on lemurs were mostly based on a reduced number of molecular markers (but see [45,46,47]) or a single demographic method (but see [19, 48]). Also, with the exception of very few studies (e.g., [19, 48, 49]), possible effects of population structure (i.e., non-random mating) on demographic inferences have been largely neglected and results were often interpreted exclusively as size changes of panmictic populations (e.g., [45, 46]). However, it has been shown that population structure can generate spurious signals of population size changes, even when the populations were stationary through time [50,51,52,53,54]. To overcome these limitations, we used two types of genome-wide data (Restriction site Associated DNA sequencing (RADseq) and whole-genome sequences) and three complementary demographic approaches (Stairway Plot, PSMC and IICR-simulations) to model both changes in effective population size (Ne) and connectivity over time in M. murinus and M. ravelobensis in northwestern Madagascar. The Stairway Plot [55] and PSMC [56] methods have been widely used to detect population size changes. The first typically uses genomic data from independent loci obtained from a population sample of several individuals, whereas the second method uses the whole-genome of a single diploid individual. The methods also differ in the fact that the Stairway Plot method has been shown to perform best towards the recent past whereas the PSMC is more informative about events occurring in a more distant past [55, 57,58,59,60]. Results from both methods have been interpreted by assuming that the genomic data used for the analyses stem from a panmictic population. However, theoretical work and simulations have shown that the PSMC dynamics observed for many species that is interpreted as reflecting population size changes might also be caused by population structure and changes in connectivity. More specifically, the PSMC curves are actually estimates of a complex temporal and sample-based function called the IICR (Inverse Instantaneous Coalescence Rate). Under panmixia, the IICR and the PSMC should be interpreted as changes in population size, but the same PSMC plot may also have been generated under vastly different modelling conditions, such as a stationary n-island model that undergoes a series of changes in connectivity [53]. For this reason, we also tested the impact of different modeling parameters on the estimated IICR in order to evaluate whether the trajectories revealed by coalescent-based methods such as the PSMC [56] may also be the result of potential changes in connectivity in a structured population of a constant size [53, 61].

Based on our knowledge about Quaternary climatic oscillations in Africa/Madagascar and in the ecology of our study species, we hypothesize that (I) M. murinus should have undergone a founder effect during its relatively recent colonization of northwestern Madagascar and should subsequently have expanded its population size in the region but not during the LGM (hypothesis I); (II) M. ravelobensis shows signals of a demographic decline and/or reduced levels of population connectivity during the LGM (hypothesis II); (III) both M. murinus and M. ravelobensis show signatures of population size increase and/or higher levels of population connectivity during the AHP (hypothesis III); and finally (IV) both mouse lemur species underwent a population decline and/or a decrease in population connectivity after the termination of the AHP (hypothesis IV) (see Table 1). The results of our study will provide a first step towards a better understanding of species responses to past climatic changes in Malagasy lowland forests.

Table 1 Summary of the expected demographic dynamics for M. ravelobensis and M. murinus under population panmixia and population structure
Full size table

Results

Genomic resources

Both species of mouse lemurs were trapped in two forest sites within the Ankaranfantsika National Park (ANP), Ravelobe and Ankomakoma. Similar to previous studies in the region, the capture data revealed that the two species were not evenly distributed within the two forest sites, even though the trapping effort was the same along all transects [41] (Fig. 1b,c and Additional file 1: Table S1 for details). Two genomic datasets were generated based on the RADseq data for complementary analyses. The Analysis of Next Generation Sequencing Data (ANGSD; dataset 1) [62] pipeline resulted in genotype likelihood information from a total of 324,608 variable sites for M. murinus (n = 22) and 601,571 variable sites for M. ravelobensis (n = 56). The number of SNPs in the called genotype dataset (dataset 2) yielded 122,053 SNPs for M. murinus (n = 22) and 242,121 for M. ravelobensis (n = 56). The mean depth per individual ranged between 13.82X and 42.18X for M. murinus, and between 13.96X and 47.88X for M. ravelobensis (Additional file 1: Tables S2 and S3). Additionally, we generated whole-genome sequences of three female individuals. Whole-genome sequencing and mapping of the single M. murinus and the two M. ravelobensis individuals resulted in 2,197,400,000 (M. murinus), 2,154,480,000 (M. ravelobensis, Ravelobe) and 2,159,200,000 (M. ravelobensis, Ankomakoma) sites for the PSMC analyses. The mean depth of coverage ranged between 16.02X (M. murinus, Ankomakoma), 17.04X (M. ravelobensis, Ravelobe) and 18.79X (M. ravelobensis, Ankomakoma).

Present day population structure and isolation-by-distance

A weak to moderate genomic signal of population structure was detected in the two mouse lemur species. First, pairwise FST estimates suggested low levels of genetic differentiation between the two forest sites, although the value was significant for M. ravelobensis (M. murinus: FST = 0.011, p = N.S.; M. ravelobensis: FST = 0.015, p < 0.0001). Second, individuals of both species were rather continuously spread along the first axis (> 22% variation explained) of a Principal Component Analyses (PCA), although no overlap existed between both sites (Additional file 1: Fig. S1a, b). Third, a rather fine-scale population-genomic structure was revealed among the two sites for M. murinus and M. ravelobensis by the NGSadmix analyses, as the admixture plots revealed no complete separation under K = 2 (Fig. 2, but see Additional file 1: Figs. S2 and S3). Instead, substantial levels of admixture were observed in both species and both locations, indicating the occurrence of gene flow between the two close-by sites. Finally, a Mantel test based on all dyadic comparisons revealed a positive and significant correlation between geographical and genetic distance in both mouse lemur species (Mantel statistic r = 0.2344 for M. murinus and r = 0.2173 for M. ravelobensis, p < 0.001), supporting a pattern of isolation-by-distance.

Fig. 2
figure 2

Population genomic structure of the two mouse lemur species. a Clustering assignment of 22 M. murinus individuals, and b Clustering assignment of 56 M. ravelobensis individuals to two genetic clusters (K = 2) using dataset 1, respectively. Each single vertical bar represents an individual and each color a distinct genetic cluster. Samples are sorted according to sampling site and respective latitude. Animal illustrations copyright 2013 Stephen D. Nash / IUCN SSC Primate Specialist Group. Used with permission

Full size image

Demographic modelling

The demographic history of M. murinus and M. ravelobensis was inferred using three complementary approaches (Stairway Plot, PSMC, IICR-simulations). The Stairway Plot [55] and the PSMC (Pairwise Sequentially Markovian Coalescent, [56]) were used to infer population size changes from the RADseq data and whole-genome sequences, respectively. The Stairway Plot analyses with the entire dataset for each forest site (Additional file 1: Fig. S4) and for the two sites together (Fig. 3b) suggested the same overall demographic trend for M. murinus and M. ravelobensis. Both species exhibited an increase of population size < 100 kyr and unexpectedly reached their largest size between the LGM and the African Humid Period (Ne ~ 210,000 for M. ravelobensis and ~ 60,000 for M. murinus). These maxima were followed by a subsequent continuous decline in population size that started at the onset of the AHP and lasted until the present. Stairway Plot results consistently suggested a higher Ne for M. ravelobensis than for M. murinus, even after standardizing the number of individuals in the analyses (n = 7 per forest site and species, Additional file 1: Fig. S5). Despite the sex-biased dispersal patterns suggested for our study species [63,64,65], the Stairway Plot analyses performed with the M. ravelobensis dataset (n = 22 males and 22 females) confirmed an identical demographic history when considering males and females separately (Additional file 1: Fig. S6).

Fig. 3
figure 3

Reconstruction of history of M. murinus and M. ravelobensis using three complementary methods. The grey vertical bars identify three well-pronounced climatic events in Africa: LIG (Last Interglacial), LGM (Last Glacial Maximum) and AHP (African Humid Period). All analyses were performed considering 2.5 years as generation time and 1.2 × 10–8 as mutation rate. a Demographic history inferred by the PSMC method using “4 + 25*2 + 4 + 6” free atomic time intervals. The thick lines represent the inferred mean trajectories for three populations, and each light line represents 100 subsampled bootstrap replicates for each individual. The humps observed in PSMC plots of the two mouse lemur species during the last 2–5 kyr seem to be a common artefact (e.g., see [6, 56]) and do not correspond to a real demographic event. Note that those humps are no longer present when we consider a different free atomic time interval (see Additional file 1: Fig. S6). b Reconstruction of the demographic history of M. murinus (N = 22) and M. ravelobensis (N = 55) using the Stairway Plot method, considering the two forest sites together. c and d PSMC plots estimated from the genomic data and IICR obtained from simulated data for both species assuming (i) a n-island model of migration; (ii) a constant population size over the time, (iii) the occurrence of five changes in population connectivity. The color code for the top horizontal bars summarizes the inferred changes in connectivity across time for each species. According to the simulation results, population connectivity changed in M. murinus at ~ 129.1 kyr (LIG), 42.7 kyr, 30.7 kyr, 13.7 kyr (onset of AHP) and 5.1 kyr (termination of AHP, Fig. 3c). Accordingly, population connectivity changed in M. ravelobensis at ~ 338.9 kyr, 135.6 kyr (LIG), 27.1 kyr, 20.1 kyr (LGM) and 7.9 kyr (Fig. 3d). RAV = Ravelobe; ANK = Ankomakoma

Full size image

The PSMC suggested a similar demographic trend for M. ravelobensis as the previous method, but different dynamics for M. murinus (Fig. 3a). Specifically, it inferred a rather stable population size of M. ravelobensis between the LIG period (~ 130 kyr) and about 30 kyr which was followed by a population increase that reached its maximum around the onset of the AHP in Madagascar (~ 15 kyr; [18]; Ne ~ 45,000 for Ravelobe and ~ 55,000 for Ankomakoma). Population sizes then decreased towards the present. The recent peak observed in both curves could be interpreted as an increase of population size. However, this peak was no longer present in the PSMC when using a different free atomic time interval (Additional file 1: Fig. S7). Therefore, the recent peaks most likely represent a common artefact (also present in other studies; e.g., see [6, 56]) rather than a real demographic event.

The PSMC for M. murinus suggested an ancient massive population decline which started before the LIG period and reached a minimum at ~ 70 kyr (Ne ~ 10,000). This bottleneck was followed by a population recovery that reached an unexpected moderate maximum during the LGM (19–26.5 kyr, Ne ~ 25,000). Afterwards, M. murinus experienced a second population decline that lasted until the AHP. This event was followed by a population size increase until ~ 3 kyr that may be partly, but not entirely, due to the above mentioned artefact, and a subsequent population decline towards the present (see Additional file 1: Fig. S7).

The previous interpretations of the PSMC rely on the assumption that population structure can be neglected and that individuals were part of panmictic populations. Alternatively, it is possible that mouse lemur populations are structured into sub-populations (= islands) and that changes in the PSMC are caused by changes in migration rate between them over time (= changes in connectivity). As a first approximation we estimated the effects of connectivity changes in a n-island model under constant population size using simulations of the IICR (Inverse Instantaneous Coalescence Rate) as in [61]. These simulations suggested that models with a large number of islands (n = 61 for M. ravelobensis and n = 84 for M. murinus) and five major historical changes in population connectivity per species may best explain the dynamics inferred for the IICR (Fig. 3c, d). The congruence between the two IICR curves inferred by the PSMC and by simulations is relatively high for the dynamics older than 7 kyr in both species. The timing of the changes in population connectivity overlapped only partly between the two mouse lemur species (horizontal bars above curves in Fig. 3c, d). In both species, population connectivity appeared to have been higher during and after the LIG than during the LGM. However, whereas population connectivity stayed relatively low for M. ravelobensis across the last 30 kyr, population connectivity was again higher during the AHP for M. murinus, the species thought to have undergone a recent expansion into northwestern Madagascar [42].

Discussion

The baseline: present-day population structure and connectivity among forest sites

The relatively low levels of genetic differentiation and the relatively high levels of genomic admixture between sampling sites suggest the existence of gene flow and therefore genetic connectivity among the two sampling sites in both mouse lemur species at present times. Ravelobe and Ankomakoma are separated only by about 10 km straight-line distance, but some larger patches of savannah would preclude straight-line dispersal between them (e.g., [42, 66]). Although the two sub-populations showed an isolation-by-distance effect, our results suggest that the two forest sites were still connected via some forest corridors that surrounded the savannah (see Additional file 1: Fig. S8) and may facilitate the dispersal of small organisms such as mouse lemurs [67]. Nevertheless, dispersal events are likely limited to small geographical distances, and should by all means be smaller than the geographic distance between both sites, since dispersal distances of mouse lemurs were previously shown to reach a maximum of up to 1 km (M. murinus in western Madagascar; [68] but see also [66, 69, 70]). Furthermore, multiple genetic studies showed that most dispersal events in mouse lemurs (e.g., M. murinus; M. ravelobensis; M. tavaratra) [64, 66, 68,69,70,71], but also in larger body sized lemur species with longer dispersal distances (e.g., Propithecus tattersalli and P. perrieri) [72, 73] occur among neighboring social groups. Altogether, our results suggest that mouse lemurs from both forests are part of a larger set of sub-populations that were connected, and maybe still are.

Taking stock: historical ecological changes inferred for northwestern Madagascar

The longest paleoenvironmental record available for Madagascar stems from a sediment core from Lake Tritrivakely in the central highlands (154 kyr) [32]. This record is the only one from the island reaching back to the Last Interglacial, and suggested that the vegetation of the highlands during this period was dominated by a grassland/woodland mosaic. This vegetation was replaced by Ericaceae (adapted to survive strong seasonal droughts) [31, 32] during the LGM as a result of the abrupt decrease in temperature (~ 4 ℃ cooler than today) [31, 74]. Multiple paleoenvironmental records from other parts of the island and corresponding to more recent periods confirmed substantially drier and cooler conditions during the LGM (e.g., [18, 30, 75, 76] but see [77]), suggesting that a general contraction of forest habitats across Madagascar at that time period is likely [30]. Such conditions during the LGM probably induced lowland forests to contract, possibly to riverine refugia, as those areas would have secured access to water and generated rather mesic local conditions [33]. As a consequence, many forests would have become partly isolated or only connected via mosaic forest corridors. This period was followed by an abrupt warming and increasing levels of moisture across the island, largely in congruence with the African Humid Period in continental Africa [17, 18, 31, 32]. A paleoenvironmental record from the Anjohibe Cave in northwestern Madagascar confirmed that the climatic conditions were wetter in the early-Holocene (9.1–4.9 kyr; AHP). It can therefore be expected that the wetter conditions favored the expansion of the woodland/grassland vegetation in the Malagasy lowland forests [25]. Finally, the period that followed the AHP was marked by a climatic warming and aridification in northwestern Madagascar, resulting in a vegetation shift towards a grassland-dominated landscape at ~ 4.8 kyr [17]. This ecosystem shift was likely intensified by the introduction of swidden agriculture and spread of pastoralism in the region during the last two millennia [17, 75, 78].

Demographic history of M. ravelobensis

Assuming that population structure is negligible, the PSMC and Stairway Plot inferred a population expansion of M. ravelobensis starting before the LGM, whereas the IICR-simulations under population structure suggested that the same PSMC dynamics could be the result of a decreased population connectivity under constant population size. Considering the presumably rather cold and arid environmental conditions in northwestern Madagascar across the LGM [30, 33], a population expansion was rather unlikely. However, a decrease in connectivity between sub-populations of M. ravelobensis would be concordant with the hypothesized contraction of forests during this period, since animal populations in the resulting mosaic landscapes would be less connected than before (hypothesis II; see introduction).

The PSMC and Stairway Plot suggested that M. ravelobensis subsequently underwent a demographic decline that started before the onset of the AHP. Such a scenario is also rather unlikely, given that lowland forests were presumably at their maximum extension at these times due to warmer temperatures and sufficient water supply [14, 15]. Our IICR-simulations suggested that the same PSMC curve could be the result of low levels of population connectivity for M. ravelobensis during the AHP. This scenario also contradicts our predictions for mouse lemur populations during the AHP (hypothesis III). Moreover, they are not in congruence with the ecological preference of M. ravelobensis for mesic microhabitats that should have been widespread during this time (e.g., [41, 79]). One possible explanation for these conflicting results could be the potential interspecific competition with M. murinus after its expansion into an inter-river system that was previously only inhabited by one mouse lemur species, M. ravelobensis [40]. Previous field studies confirmed that M. murinus has a higher competitive potential than M. ravelobensis in an experimental setting [80]. The spatial expansion of M. murinus may have required new patterns of habitat partitioning and resulted in new direct or indirect interspecific competition (e.g., for food resources, sleeping sites). These species interactions may have precluded establishing large population sizes and/or higher population connectivity for M. ravelobensis during the AHP.

The PSMC and Stairway Plot inferred a continuing population decline for M. ravelobensis after the termination of the AHP, which would be in concordance with hypothesis IV. The increasing aridification during the mid-Holocene and the ecosystem shift towards a grassland-dominated landscape [17, 75, 78] might have contributed to this development. Such a decline was already inferred in previous studies on this and other lemur species inhabiting the northwestern region (e.g., Microcebus bongolavensis, M. danfossi, and Lepilemur edwardsi) [81, 82] and in multiple Malagasy species distributed across the island (e.g., frogs, birds of prey and rodents) [83,84,85] and was mostly attributed to intensified human pressures on the island during the last thousand years. It should be noted that the IICR-simulations did generate the same dynamics under a model of population structure and constant population size. This approach required similar population connectivity levels during recent times than in the LGM. Such reduced levels of population connectivity would be congruent with the loss of suitable mouse lemur habitats across the late-Holocene [17, 75, 78]. In conclusion, it appears that both changes in population connectivity (LGM) and size (late-Holocene) may have happened and shaped the demographic dynamics of M. ravelobensis populations over time.

Demographic history of M. murinus

The late Pleistocene demographic dynamics of M. murinus were expected to be shaped by its rather late colonization of northwestern Madagascar [42]. Such a colonization would be a classical example of a founder effect with a relatively small number of individuals arriving in this region via a highland corridor [86, 87]. Such an event would very likely mirror a population bottleneck and would have been followed by a rapid spatial and demographic expansion into the lowland forests between the Betsiboka and Mahajamba rivers (hypothesis I). Our demographic analyses, if interpreted in terms of population panmixia, indeed suggested that M. murinus underwent a strong population bottleneck at around 70 kyr, and that the effective population size (Ne) subsequently increased and reached a maximum during the LGM. In fact, our time estimate for the founder effect (~ 70 kyr) does agree very well with that of the previous study (26.5–33.5 kyr) [42] if taking into account that the generation time used in both studies differed by the factor 2.5 (2.5 years [this study, 18, 36] versus 1 year [42]; see methods for details on generation time). However, the IICR-simulations revealed that the IICR dynamics may also have been the result of changes in population connectivity under stable population size. Indeed, the IICR-simulations suggested higher levels of population connectivity during and after the LIG than during the LGM, which would be in concordance with the predicted forest contractions during the LGM [30, 33]. These two competing interpretations cannot be easily reconciled, since both scenarios do partially fit hypothesis I for this species. Information about the definite time point of colonization by M. murinus or details of the paleoenvironmental dynamics in this inter-river-system back to the LIG are needed to evaluate the two alternative scenarios.

Towards more recent times and in contrast to our expectations, the PSMC and Stairway Plot suggested that M. murinus underwent a temporary demographic decline that started before the onset of the AHP. Conversely, the IICR-simulations point towards higher connectivity levels during the AHP under constant population size, which would be in line with hypothesis III. Assuming that increasing humidity after the LGM likely resulted in an expansion of dry deciduous forests [17], higher levels of gene flow among sub-populations would be rather expected, and were previously documented for M. arnholdi in northern Madagascar [18].

Finally, both PSMC and the Stairway Plot revealed a population decline for M. murinus during the late-Holocene (hypothesis IV), which is congruent with the climatic warming and aridification that followed after the AHP [17] and the presumably increasing degree of habitat fragmentation in the region [17, 75, 78]. The timing of the M. murinus population decline coincided well with previously documented bottlenecks for other Malagasy species and with the collapse of the Malagasy megafauna (< 3 kyr; e.g., [25, 30, 88]).The lower levels of population connectivity during recent times suggested by the IICR-simulations would also be in concordance with the hypothesis IV.

Conclusions

The present study suggests that climatic fluctuations have been important drivers of evolutionary trajectories for mouse lemurs in northwestern Madagascar. Interestingly, our demographic reconstructions also revealed distinct dynamics for M. murinus and M. ravelobensis, suggesting that even closely related species may differ in their responses to the same climatic events. Different demographic scenarios emerged for both mouse lemur and the decision for one of possible alternative explanations was not always straight forward. Population structure and changes in connectivity very likely impacted the demographic dynamics of mouse lemurs in this region of Madagascar (see Additional file 1: Table S4), while forest contractions and expansions may have extensively shaped the history of lowland forests during the Late Pleistocene (e.g., [33]) like in other regions of the world (e.g., [89,90,91,92,93]). For example, a decrease in population connectivity may explain the IICR dynamics during the LGM better than population size changes in both species. However, it also became clear that changes in connectivity alone may not explain well all findings. For instance, if the PSMC dynamics before the LGM would be the result of changes in connectivity alone, this would imply that M. murinus must have colonized northwestern Madagascar already earlier than the LIG, because a founder event and subsequent colonization of this IRS would necessarily include population size changes. Systematic simulations are needed to clarify the consequences of a founder event and a subsequent spatial expansion on the IICR dynamics when considering a n-island model of migration.

The two alternative models studied here represent endpoints on a scale of possible events, and mixed scenarios including both population size and connectivity changes may ultimately explain best the historic population dynamics of mouse lemurs in northwestern Madagascar (e.g., see [18, 19]). Therefore, the implementation of model-based approaches such as the Approximate Bayesian Computation [94] or composite-likelihood methods [95] is crucial to test and compare alternative demographic hypotheses including both population size and connectivity changes (e.g., [19, 96, 97]). In addition, further high-resolution paleoenvironmental reconstructions for northwestern Madagascar are urgently needed to better understand the impact of the last Interglacial–Glacial cycle on the lowland dry forest dynamics. Altogether, our study shows that it is essential to consider the impact of different model assumptions (e.g., panmixia and population structure) when exploring, weighting and inferring alternative species demographic scenarios. In particular, the application of the PSMC should in the future be always complemented with a simulation approach in order to avoid oversimplification.

Methods

Study species

The sympatric M. murinus and M. ravelobensis display marked differences in their distribution, evolutionary history and ecology. It was previously estimated that M. murinus and M. ravelobensis diverged by about 9.60 Ma [98], although diversification dates for mouse lemurs were recently challenged [47]. While M. murinus is regarded as a habitat generalist due to its wide distribution across various forest types in western and northwestern Madagascar (reviewed in [39]), M. ravelobensis is assumed to be more specialized, although it was shown to reach higher local abundancies and a broader distribution within the Ankarafantsika National Park than its congener [41]. Both species occur in partial sympathy in northwestern Madagascar but were shown to respond differently to habitat fragmentation [79]. While M. murinus can be often found in even small forest fragments and seems to have a higher overall vagility across mixed and partially open landscapes, M. ravelobensis is typically found in larger forests and may be less able to connect across open landscapes [79]. Both species were shown to have different microhabitat preferences [41, 99], and direct competition between the two species is therefore not regarded as main mechanism of abundance regulation [41, 79]. Both species forage solitarily during the night, but form sleeping groups during the day that vary in composition between the species. M. murinus males sleep alone and females form groups of related individuals in wooden tree holes, while M. ravelobensis forms mixed‐sex sleeping groups consisting of matrilinear relatives of varying degree of relatedness who use diverse substrates as sleeping site (e.g., lianas and leaves) [65, 100, 101]. Natal dispersal in M. murinus is male-biased [63, 64], whereas M. ravelobensis displays only moderate and delayed male-biased dispersal which was suggested to lead to a higher risk of inbreeding [65].

Study area and sample collection

The Ankarafantsika National Park (ANP) is located in northwestern Madagascar, in an area of about 135,000 ha of dry deciduous forest that is delimited by the Betsiboka (western limit) and the Mahajamba rivers (eastern limit; IRS Ia; Fig. 1a) [40, 41]. The ANP is one of the largest remaining forest blocks in western Madagascar [42], although it shows some degree of forest fragmentation towards the edges [41, 79]. The locally sympatric M. murinus and M. ravelobensis were sampled around two natural lakes on the western part of the ANP, Ravelobe (Rav, –16.302413°N, 46.821346°E) and Ankomakoma (Ank, -16.342752°N, 46.740293°E, Fig. 1b, c). The sites are approximately 10 km apart and characterized by a mosaic of savannah and forest corridors. Ravelobe (85–176 m above sea level, a.s.l.) is situated near Ampijoroa next to the National Road RN4 that connects Antananarivo to Mahajanga, while Ankomakoma (105–185 m a.s.l.) is part of a mosaic landscape consisting of dry deciduous forest with interspersed savannah patches (Additional file 1: Fig. S8).

Fieldwork took place during the dry season of 2017 (April–July). Four 1 km transects were installed for field work in each site in order to cover various forest parts around the lakes. All transects followed pre-existing dirt roads or foot paths. A total of 1200 Sherman Traps (Sherman Traps Inc, Tallahassee, FL, USA), baited with banana, were installed overnight in Ravelobe across 12 nights, and 1000 traps were installed in Ankomakoma across 10 nights (see [102] for details). M. murinus and M. ravelobensis were distinguished in the field based on their head coloration (greyish in M. murinus vs. brownish in M. ravelobensis) [103] and their distinctive tail length (130.81 ± 6.15 mm in M. murinus vs. 155.48 ± 7.57 mm in M. ravelobensis) [40]. All animals were released at dusk of the same day at their individual capture position. Small ear biopsies (approx. 2–3 mm2) were taken from all captured animals for genomic analyses. Tissue samples were stored in Queen’s lysis buffer [104] during the field season and subsequently at − 20 °C in the laboratory.

RADseq library & whole-genome sequencing

A total of 24 M. murinus (7 from Ravelobe and 17 from Ankomakoma) and 60 M. ravelobensis (33 from Ravelobe and 27 from Ankomakoma) individuals were available for RADseq, based on the spatial distribution and abundance of each mouse lemur species per transect and sampling site. Total genomic DNA was extracted from the ear biopsies using the DNeasy Blood & Tissue Kit (Qiagen) following the manufacturer’s protocol with few modifications (see [46] for details), and the DNA concentration was estimated with the Qubit® Fluorometer (Life Technologies). DNA samples (~ 200 ng of DNA) were then digested with the restriction enzyme SbfI at the GenoToul-GeT-PlaGE Core Facility (Toulouse, France). RAD Libraries were prepared in sets of 24 samples sorted by original DNA concentration, where each sample was assigned to one of 48 unique barcode sequences during adapter ligation. Sub-libraries were randomly sheared [105, 106] using a Covaris M220 ultrasonicator, resulting in fragments with an average size of 550 bp. Sheared DNA fragments were ligated to the second adapter and all fragments with both adapters were amplified in 10 Polymerase Chain Reaction (PCR) cycles. DNA concentration and fragment sizes of the amplified libraries were verified using qPCR and Fragment Analyzer. The sub-libraries were sequenced using 150 bp paired-end reads on an Illumina HiSeq3000 platform at the GenoToul-GeT-PlaGE Core Facility (Toulouse, France). Raw data was initially demultiplexed using the tool splitbc implemented on the FASTX-toolkit (http://hannonlab.cshl.edu/fastx_toolkit/), and the quality of the raw data was verified with FastQC v0.11.7 (http://www.bioinformatics.bbsrc.ac.uk/projects/fastqc). Trimmomatic v0.36 [107] was used to remove Illumina adapters from the reads, remove reads with low quality bases, and to trim the reads to a minimum 4-base sliding window with quality score below 15. Additionally, reads with less than 60 bp length after the filtering steps were removed from the analyses. After quality filtering, BWA-MEM (http://bio-bwa.sourceforge.net/) was used to align the paired-end reads to a high quality genome assembly of Microcebus murinus (GCA_000165445.3) [108]. Finally, the software SAMtools v1.8 [109] was used to remove PCR duplicates (i.e., sequence reads that result from sequencing two or more copies of the exact same DNA fragment; [110]) and to convert the Sequence Alignment Map (SAM) format to the corresponding binary version (BAM). To ensure that only autosomal data was used for the analyses, the M. murinus scaffold NC_033692.1 (designated as the X-chromosome) as well as the M. murinus mitochondrial genome (NCBI Accession Number: KR911908.1) were excluded from the aligned BAM files. The autosomal BAM files were then used as input files for all downstream genomic analyses.

In addition to the RADseq dataset, one M. murinus (Ankomakoma) and two M. ravelobensis samples (Ravelobe and Ankomakoma) were selected for whole-genome sequencing. All samples were females which are the philopatric sex in both species [65, 69]. Given that no female M. murinus was caught at Ravelobe, no whole-genome sequence was generated at this site for that species. Libraries were prepared at the Institute for Animal Breeding and Genetics of the University of Veterinary Medicine Hannover using the NEBNext Ultra DNA Library Prep Kit from Illumina (New England BioLabs, Ipswich, MA, USA). DNA samples (~ 200 ng of DNA) were first sheared with an ultrasonicator (Covaris M220, Woburn, Massachusetts, USA) and the respective sizes were selected according to the manufacturer’s recommendations. Whole-genome sequencing was performed on an Illumina NextSeq 500 (Illumina, San Diego, CA, USA) for 300 cycles in paired‐end mode. Visual quality control of whole-genome sequencing data was performed using FastQC version 0.11.7. Reads were trimmed using PRINSEQ version 0.20.4 [111] and mapped to the M. murinus reference genome (“Microcebus_murinus.Mmur_3.0.dna.toplevel.fa.gz”) using the BWA-MEM algorithm implemented in the BWA version 0.7.17 [112]. Similarly to the RADseq dataset, the reads that mapped against the M. murinus X-chromosome and mitochondrial genome were discarded from our analyses (see details above). For all analyses with the whole-genome sequences, the read depth and quality of the variant sites were controlled by applying the following quality filters: base quality above 20, mapping quality above 30, minimum read depth of 3 and a maximum read depth of 100. See Additional file 1: Text S1 and S2 for details about the sequencing libraries.

RADseq datasets: genotype likelihoods & SNP calling

Next-Generation Sequencing platforms can generate large amounts of sequencing data but are prone to sequencing errors [113,114,115]. It is therefore advisable to keep the data in form of genotype likelihoods during downstream analyses, because genotype likelihoods retain information about uncertainty in base calling, which enable to control for some problems commonly associated with RADseq datasets (e.g., unevenness in sequencing depth and allelic dropout) [114, 116, 117]. Consequently, most of our analyses were carried out on genotype likelihoods estimated with ANGSD [115], considering the following filters: a minimum base quality of 20, a minimum mapping quality of 30, minimum Minor Allele Frequency below 0.5, a minimum mean depth of coverage of 4X [115, 117, 118], and sites present in at least 75% of the individuals (dataset 1: genotype likelihoods). In addition to the genotype likelihoods, genotypes were called for each mouse lemur species for complementary analyses. Genotypes were called with SAMtools v1.8 [109] using the same quality filters previously used in ANGSD (see above), but considering a minimum depth of coverage of 10× to ensure high-confidence genotype calls (dataset 2: genotype calls) [119].

Knowing that population genetics analyses can be biased by social structure and relatedness between individuals (e.g., [120, 121]), a relatedness analyses was performed next. Relatedness between two individuals is usually described by the concept of identity-by-descent (IBD), where two alleles are considered identical by descent if they recently descended from a common ancestral allele [117, 122]. NGSrelate [117] was used to calculate the IBD coefficients (k0, k1 and k2; probability of two individuals sharing 0, 1 or 2 alleles from a single ancestor at any locus, respectively) between pairs of individuals using dataset 1 [122]. First-degree relatives were then inferred based on the comparison of the obtained IBD coefficients with the expected IBD probabilities (i.e., parent-offspring: k0 = 0, k1 = 1 and k2 = 0; full sibs: k0 = 0.25, k1 = 0.50 and k2 = 0.25) [122]. Only one individual of each dyad of closely related individuals was retained in our dataset for downstream analyses (see Additional file 1: Text S3).

Population-genomic structure & isolation-by-distance

Population-genomic structure patterns were investigated in M. murinus (n = 22) and M. ravelobensis (n = 56) using three distinct approaches. First, genetic differentiation between both sites was estimated per species using Wright’s F-statistics FST [123]. The analyses were performed with ARLEQUIN [124] and significance was tested using 10,000 permutations (dataset 2). Second, signals of population genetic structure were evaluated by Principal Component Analysis (PCA) computed using PCAngsd (dataset 1) [125]. The PCA Eigenvalues that explained most of the genetic variation between individuals (PC1 and PC2) were extracted and individuals were plotted using R (R CoreTeam 2014). Third, NGSadmix [114], a maximum-likelihood clustering method, was used to assign individuals to a specific number of clusters (K) using dataset 1. For both species, the number of explored clusters ranged between 1 and 3, and a total of 10 independent runs were performed for each value of K. The most suitable value of K was determined with Clumpak [126], following the Evanno method [127].

Geographically restricted gene flow results in a significant correlation between genetic and geographic distance, known as isolation-by-distance [128], where the genetic dissimilarity increases with the increase of geographical distance [67]. The effect of geographic distance on genomic differentiation over our small geographic scale was investigated using the individual as the unit of the analyses [71]. Genetic dissimilarity between any two individuals was measured by the Rousset’s genetic distance (â) [129], an estimator analogous to the FST/(1-FST) ratio using pairs of individuals instead of populations. Geographic distance was measured as the linear geographical distance in km separating each pair of individuals. Both Rousset’s genetic distance and geographic distance were computed with SPAGeDi [130] based on dataset 2. The occurrence of isolation-by-distance was finally investigated with a Mantel test [131] using the VEGAN package [132] available in R (R Development Core Team 2005). Significance was determined via 10,000 permutations.

In addition, genetic summary statistics were calculated for each species and forest site (see Additional file 1: Text S4 and Table S5 for details). Inbreeding coefficients per individual (F) were also estimated using dataset 2 to detect deviations from Hardy–Weinberg equilibrium (Additional file 1: Tables S6 and S7) [133]. One M. ravelobensis individual captured at Ravelobe showed a highly negative F value and was therefore excluded from the demographic analyses (Additional file 1: Table S6) [133]. For details about the set of individuals used in each step of this study see Additional file 1: Table S8.

Mutation rate and generation time

The demographic history of both species was investigated using three complementary modeling approaches: Stairway Plot [55], PSMC [56] and IICR-simulations [53]. All analyses were performed assuming a mutation rate value of 1.2 × 10–8 [47, 134]. This mutation rate is the most accurate estimate available for mouse lemurs and was calculated from average pedigree-based estimates of seven primates species [47, 135]. During the last decade, generation time (GT) values between 1 and 4.5 years have been used for genetic studies on mouse lemurs [36, 40, 134]. Since a recent study on free-living M. murinus from the ANP estimated a 2.5 year generation time as the average age of parents [36], we recently compared the performance of three GT values (1, 2.5 and 4.5 years) for another mouse lemur species (M. arnholdi) using alternative demographic methods. As the results under 2.5 years fitted best to on-site high-resolution paleoenvironmental reconstructions [18], we decided to also use this estimate for M. murinus and M. ravelobensis in this study.

Stairway plot

The Stairway Plot [55] method uses the Site Frequency Spectrum (SFS; i.e., the distribution of the allele frequencies of a given set of SNPs in a population) [136] from population genomic sequence data to estimate a series of population mutation rates (θ = 4Neµ) following a multi-epoch demographic model, where epochs coincide with coalescent events [55]. The realSFS tool implemented in ANGSD [115] was used to estimate a folded one-dimensional SFS for each mouse lemur species (i.e., considering the two forest sites together) based on dataset 1, because we lack a suitable outgroup to determine the ancestral state of each allele [136]. The folded 1d-SFS of each species was used as input data to generate 199 additional SFS by bootstrap, following the software guidelines [55]. Inferences were then made based on the 200 1d-SFS for each species with Stairway Plot v2.0 [55]. Given the observed isolation-by-distance pattern in both M. murinus and M. ravelobensis, the Stairway Plots were also repeated considering the two forest sites separately. Knowing that changes in Ne through time are dependent on the number of possible coalescent events, this method is sensitive to the number of individuals and SNPs in the dataset [59, 60, 137]. Therefore, analyses were first performed considering the entire dataset of each forest site and species (n ranging between 7 and 30), and also repeated with an equal sample size (n = 7) per site and species. For this analysis individuals of each site were randomly selected using R (R Development Core Team 2005). Finally, to evaluate the effect of sex-biased dispersal in the demographic inferences of mouse lemurs, we reran the Stairway Plot for the M. ravelobensis dataset considering males and females separately. The analyses were performed considering all M. ravelobensis females available in our dataset (n = 22) and an equal number of males to avoid differences in Ne related to sample size. The males were randomly selected using R (R Development Core Team 2005).

PSMC & IICR-simulations

The PSMC [56] method makes use of the distribution of heterozygous sites along the entire genome of a single diploid individual to estimate the population size change history that best explains the corresponding coalescence times. As shown by [53] the PSMC method is actually inferring the IICR of the sampled individual. This IICR can be interpreted as a history of population size change under panmixia, but under structured models it is more difficult to interpret. Here we focused on the interpretation of PSMC plots under either panmixia or piecewise-stationary n-island models. In the latter, we allowed gene flow to be constant for certain periods of time and to change between periods of time that we call components, always assuming that the population structure can be approximated by an n-island model. For instance, if we assumed that there were three components this means that there were three period of time with respective migration rates M1, M2 and M3, and two times at which the migration rates changes, namely t1 and t2 (t0 is assumed to be the most recent time of sampling, i.e., t0 = 0). The three whole-genome autosomal sequences (one M. murinus and two M. ravelobensis) generated as part of this study were submitted to PSMC analyses, considering the following parameters: minimum read depth per site of 3 (-d3), maximum read depth per site of 100 (-D100), upper limit of the TMRCA and initial θ/ρ value of 5 (-t5 and -r5). The Ne was inferred using 4 + 25*2 + 4 + 6 free atomic time intervals, with a total of 100 bootstrap replicates (command line: ‐N30 –t5 –r5 –p “4 + 25*2 + 4 + 6” ‐D100 –d3 –q30). The PSMC analyses were repeated considering the same parameters but using -p “64*1” free atomic time intervals. Since the mean genome coverage of two of our sequences was close to but below the recommended 18X (see [6]), simulations for variant coverage divergence were performed to evaluate the impact of lower genome-wide coverages in the PSMC inferences for mouse lemurs. The analyses were performed for the M. ravelobensis individual sampled in Ankomakoma by varying the minimum read depth option per site (-d) between one and twelve. Results are discussed in the Additional file 1: Text S5 and Fig. S9.

In order to find a demographic scenario that could explain the observed PSMC curves under a piecewise stationary n-island model as mentioned above, we used the SNIF (Structured Non-stationary Inference Framework) inferential method of [138]. This method uses the PSMC plot as a summary statistic of genomic information and as a target for an optimization algorithm under a piecewise stationary n-island models with unknown parameter values (N, n, ti, Mi), where N is the deme size, n is the number of demes, and the ti values correspond to the times at which the migration rates Mi change. SNIF uses a search algorithm that explores the parameter space and uses an optimality criterion to select the structured scenario that best explains a given target IICR (simulated) or PSMC (observed or simulated) curve. The method has been validated using target IICRs generated under piecewise stationary n-island models of increasing complexity (i.e., number of components or periods during which the Mi can change) by comparing inferred and simulated parameter values. Technical details of the search algorithm and distance computation can be found in [138]. Here we simply note that the method finds the parameter values that minimize a distance between the IICR curve generated under a very large number of piecewise stationary n-island models, and the observed PSMC. The algorithm stops its search when it reaches a minimum distance set by the user or a pre-set number of optimization steps. The final scenario is then validated by generating an IICR curve under these "best parameter values" and by running SNIF on this IICR to test whether the method would indeed re-infer the same values. If SNIF does not manage to infer the right scenario, this suggests that the scenario is not to be trusted as even if it were correct the method would not infer it properly. If we infer the right parameters, this suggests that the scenario is inferable by our method, but it cannot be seen as a proof that it is correct. This is the best scenario under the piecewise stationary n-island model that (i) explains the data and (ii) is of reasonable complexity and (iii) can be inferred and thus trusted to some extent. In order to produce a well estimated IICR curve we simulated 100,000 T2 values for the M. murinus and the M. ravelobensis individuals sampled in Ankomakoma using the ms software [139]. The IICR was plotted with the observed PSMC using a python script available at: https://github.com/willyrv/IICREstimator [61]. The SNIF program and its documentation can be found in github.com/arredondos/snif. The ms commands used to generate the best-fitting IICR plots can be found in the Additional file 1: Text S6.

Impact of repeat regions in the demographic modelling

It has been recently shown that genomic repeat regions may affect demographic inferences using diverse methods, including the Stairway Plot and PSMC [140]. To evaluate the impact of repeat regions in our demographic inferences, we reran the Stairway Plot for the entire M. ravelobensis dataset (n = 55) and the PSMC analyses for the three whole-genome sequences without the repeat regions. We firstly removed all repeat regions from our BAM files following the RepeatMasker output file kindly provided by J. Rogers for the Microcebus murinus reference genome (GenBank Assembly Accession Number: GCA_000165445.3 [108]). We then reran the Stairway plot and PSMC using the same command options aforementioned. Results are presented in the Additional file 1: Text S7 and Fig. S10.

Availability of data and materials

The RADseq sequences generated under this study are publicly available at Sequence Read Archive (NCBI) in the BioProject PRJNA560399 (Number Accession: SAMN16955525–SAMN16955602). Whole-genome sequences are available in the BioProject PRJNA632451. Accession Numbers: SAMN15237929 (M. ravelobensis, Ravelobe), SAMN15246264 (M. ravelobensis, Ankomakoma), SAMN15237931 (M. murinus, Ankomakoma). Scripts used for all analyses are available upon reasonable request.

Abbreviations

a.s.l.:

Above sea level

AHP:

African Humid Period

ANK:

Ankomakoma

ANP:

Ankarafantsika National Park

IBD:

Identity-by-descent

IICR:

Inverse Instantaneous Coalescence Rate

IRS:

Inter-River-System

kyr:

Thousand years

LGM:

Last Glacial Maximum

LIG:

Last Interglacial

MIG:

Migration rate

M.mur:

M. murinus

M.rav:

M. ravelobensis

Mya:

Million years ago

Ne :

Effective population size

PCA:

Principal component analyses

PSMC:

Pairwise Sequentially Markovian Coalescent

RADseq:

Restriction site Associated DNA Sequencing

RAV:

Ravelobe

SFS:

Site frequency spectrum

References

  1. Hewitt GM. Genetic consequences of climatic oscillations in the Quaternary. Philos Trans R Soc B Biol Sci. 2004;359:183–95.

    Article  CAS  Google Scholar 

  2. Martínez-Freiría F, Crochet PA, Fahd S, Geniez P, Brito JC, Velo-Antón G. Integrative phylogeographical and ecological analysis reveals multiple pleistocene refugia for Mediterranean Daboia vipers in north-west Africa. Biol J Linn Soc. 2017;122:366–84.

    Article  Google Scholar 

  3. Martínez-Freiría F, Velo-Antón G, Brito JC. Trapped by climate: interglacial refuge and recent population expansion in the endemic Iberian adder Vipera seoanei. Divers Distrib. 2015;21:331–44.

    Article  Google Scholar 

  4. Velo-Antõn G, Parra JL, Parra-Olea G, Zamudio KR. Tracking climate change in a dispersal-limited species: reduced spatial and genetic connectivity in a montane salamander. Mol Ecol. 2013;22:3261–78.

    Article  PubMed  Google Scholar 

  5. Veríssimo J, Znari M, Stuckas H, Fritz U, Pereira P, Teixeira J, et al. Pleistocene diversification in Morocco and recent demographic expansion in the Mediterranean pond turtle Mauremys leprosa. Biol J Linn Soc. 2016;119:943–59.

    Article  Google Scholar 

  6. Nadachowska-Brzyska K, Burri R, Smeds L, Ellegren H. PSMC analysis of effective population sizes in molecular ecology and its application to black-and-white Ficedula flycatchers. Mol Ecol. 2016;25:1058–72.

    Article  PubMed  PubMed Central  Google Scholar 

  7. Chattopadhyay B, Garg KM, Gwee CY, Edwards SV, Rheindt FE. Gene flow during glacial habitat shifts facilitates character displacement in a Neotropical flycatcher radiation. BMC Evol Biol. 2017;17:1–15.

    Article  Google Scholar 

  8. Sow AS, Martínez-Freiría F, Dieng H, Fahd S, Brito JC. Biogeographical analysis of the Atlantic Sahara reptiles: Environmental correlates of species distribution and vulnerability to climate change. J Arid Environ. 2014;109:65–73.

    Article  Google Scholar 

  9. Kukla GJ, Bender ML, de Beaulieu J-L, Bond G, Broecker WS, Cleveringa P, et al. Last interglacial climates. USGS Staff Res. 2002;174.

  10. Sylvestre F. Moisture pattern during the last glacial maximum in South America. In: Past climate variability in South America and surrounding regions. Springer; 2009. p. 3–27.

  11. Hewitt G. The genetic legacy of the quaternary ice ages. Nature. 2000;405:907–13.

    Article  CAS  PubMed  Google Scholar 

  12. Garg KM, Chattopadhyay B, Koane B, Sam K, Rheindt FE. Last Glacial Maximum led to community-wide population expansion in a montane songbird radiation in highland Papua New Guinea. BMC Evol Biol. 2020;20:82.

    Article  PubMed  PubMed Central  Google Scholar 

  13. Piñeiro R, Dauby G, Kaymak E, Hardy OJ. Pleistocene population expansions of shade-tolerant trees indicate fragmentation of the African rainforest during the ice ages. Proc R Soc B Biol Sci. 2017;284:20171800.

    Article  Google Scholar 

  14. Demenocal P, Ortiz J, Guilderson T, Adkins J, Sarnthein M, Baker L, et al. Abrupt onset and termination of the African Humid Period: rapid climate responses to gradual insolation forcing. Quat Sci Rev. 2000;19:347–61.

    Article  Google Scholar 

  15. Tierney JE, DeMenocal PB. Abrupt shifts in Horn of Africa hydroclimate since the last glacial maximum. Science. 2013;342:843–6.

    Article  CAS  PubMed  Google Scholar 

  16. Burrough SL, Thomas DSG. Central southern Africa at the time of the African Humid Period: a new analysis of Holocene palaeoenvironmental and palaeoclimate data. Quat Sci Rev. 2013;80:29–46.

    Article  Google Scholar 

  17. Wang L, Brook GA, Burney DA, Voarintsoa NRG, Liang F, Cheng H, et al. The African Humid Period, rapid climate change events, the timing of human colonization, and megafaunal extinctions in Madagascar during the Holocene: evidence from a 2m Anjohibe Cave stalagmite. Quat Sci Rev. 2019;210:136–53. https://doi.org/10.1016/j.quascirev.2019.02.004.

    Article  Google Scholar 

  18. Teixeira H, Montade V, Salmona J, Metzger J, Bremonde L, Kasper T, et al. Past environmental changes affected lemur population dynamics prior to human impact in Madagascar. Comun Biol. 2021. https://doi.org/10.1038/s42003-021-02620-1.

    Article  Google Scholar 

  19. Salmona J, Heller R, Quéméré E, Chikhi L. Climate change and human colonization triggered habitat loss and fragmentation in Madagascar. Mol Ecol. 2017;26:5203–22.

    Article  PubMed  Google Scholar 

  20. Myers N, Mittermeier RA, Mittermeier CG, da Fonseca GAB, Kent J. Biodiversity hotspots for conservation priorities. Nature. 2000;403:853–8. https://doi.org/10.1038/35002501.

    Article  CAS  PubMed  Google Scholar 

  21. Poux C, Madsen O, Marquard E, Vieites DR, de Jong WW, Vences M. Asynchronous colonization of Madagascar by the four endemic clades of primates, tenrecs, carnivores, and rodents as inferred from nuclear genes. Syst Biol. 2005;54:719–30.

    Article  PubMed  Google Scholar 

  22. Vences M, Wollenberg KC, Vieites DR, Lees DC. Madagascar as a model region of species diversification. Trends Ecol Evol. 2009;24:456–65.

    Article  PubMed  Google Scholar 

  23. MacPhee RDE, Burney DA. Dating of modified femora of extinct dwarf Hippopotamus from Southern Madagascar: implications for constraining human colonization and vertebrate extinction events. J Archaeol Sci. 1991;18:695–706.

    Article  Google Scholar 

  24. Dewar RE, Radimilahy C, Wright HT, Jacobs Z, Kelly GO, Berna F. Stone tools and foraging in northern Madagascar challenge Holocene extinction models. Proc Natl Acad Sci. 2013;110:12583–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Burney DA, Robinson GS, Burney LP. Sporormiella and the late holocene extinctions in Madagascar. Proc Natl Acad Sci U S A. 2003;100:10800–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Perez VR, Godfrey LR, Nowak-kemp M, Burney DA, Ratsimbazafy J, Vasey N. Evidence of early butchery of giant lemurs in Madagascar. J Hum. 2005;2005:722–42.

    Google Scholar 

  27. Douglass K, Hixon S, Wright HT, Godfrey LR, Crowley BE. A critical review of radiocarbon dates clari fi es the human settlement of Madagascar. Quat Sci Rev. 2019;221:105878.

    Article  Google Scholar 

  28. Gommery D, Ramanivosoa B, Faure M, Guérin C, Kerloc P, Sénégas F, et al. Oldest evidence of human activities in Madagascar on subfossil hippopotamus bones from Anjohibe (Mahajanga Province). CR Palevol. 2011;10:271–8. https://doi.org/10.1016/j.crpv.2011.01.006.

    Article  Google Scholar 

  29. Hansford J, Wright PC, Rasoamiaramanana A, Pérez VR, Godfrey LR, Errickson D, et al. Early Holocene human presence in Madagascar evidenced by exploitation of avian megafauna. Sci Adv. 2018;1–7.

  30. Burney DA, Pigott L, Godfrey LR, Jungers WL, Goodman SM, Wright HT, et al. A chronology for late prehistoric Madagascar. J Hum Evol. 2004;25–63.

    Article  PubMed  Google Scholar 

  31. Gasse F, Van Campo E. A 40 000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics. Quat Res. 1998;49:299–311.

    Article  Google Scholar 

  32. Gasse F, Van Campo E. Late Quaternary environmental changes from a pollen and diatom record in the southern tropics (Lake Tritrivakely, Madagascar). Palaeogeogr Palaeoclimatol Palaeoecol. 2001;167:287–308.

    Article  Google Scholar 

  33. Wilmé L, Goodman SM, Ganzhorn JU. Biogeographic evolution of Madagascar’s microendemic biota. Science. 2006;312:1063–6.

    Article  PubMed  CAS  Google Scholar 

  34. Goodman SM, Benstead JP. Updated estimates of biotic diversity and endemism for Madagascar. Oryx. 2005;39:73–7.

    Article  Google Scholar 

  35. Perret M. Influence of social factors on sex ratio at birth, maternal investment and young survival in a prosimian primate. Behav Ecol Sociobiol. 1990;27:447–54.

    Article  Google Scholar 

  36. Radespiel U, Lutermann H, Schmelting B, Zimmermann E. An empirical estimate of the generation time of mouse lemurs. Am J Primatol. 2019;81:1–8. https://doi.org/10.1002/ajp.23062.

    Article  Google Scholar 

  37. Mittermeier RA, Jr. EEL, Richardson M, Schwitzer C, Langrand O, Rylands AB. Lemurs of Madagascar. Conservation. USA; 2010

  38. Blair C, Heckman KL, Russell AL, Yoder AD. Multilocus coalescent analyses reveal the demographic history and speciation patterns of mouse lemur sister species. BMC Evol Biol. 2014;14:1–12.

    Article  Google Scholar 

  39. Radespiel U. Can behavioral ecology help to understand the divergent geographic range sizes of mouse lemurs. Gremlins night Biol Behav Conserv Biogeogr Cheirogaleidae (SM Lehman, U Radespiel, E Zimmermann, eds) Cambridge Univ Press Cambridge, United Kingdom. 2016;498–519.

  40. Olivieri G, Zimmermann E, Randrianambinina B, Rasoloharijaona S, Rakotondravony D, Guschanski K, et al. The ever-increasing diversity in mouse lemurs: three new species in north and northwestern Madagascar. Mol Phylogenet Evol. 2007;43:309–27.

    Article  CAS  PubMed  Google Scholar 

  41. Rakotondravony R, Radespiel UTE. Varying patterns of coexistence of two mouse lemur species (Microcebus ravelobensis and M. murinus) in a heterogeneous landscape. Am J Primatoll. 2009;938:928–38.

    Article  Google Scholar 

  42. Schneider N, Chikhi L, Currat M, Radespiel U. Signals of recent spatial expansions in the grey mouse lemur (Microcebus murinus). BMC Evol Biol. 2010. https://doi.org/10.1186/1471-2148-10-105.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Salmona J, Heller R, Lascoux M, Shafer A. Inferring demographic history using genomic data. In: Population Genomics. Springer; 2017. p. 511–37.

  44. Beichman AC, Huerta-Sanchez E, Lohmueller KE. Using genomic data to infer historic population dynamics of nonmodel organisms. Annu Rev Ecol Evol Syst. 2018;49:433–56.

    Article  Google Scholar 

  45. Meyer WK, Venkat A, Kermany AR, Geijn BD, Zhang S, Przeworski M. Evolutionary history inferred from the de novo assembly of a nonmodel organism, the blue-eyed black lemur. Mol Ecol. 2015;24:4392–405.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Hawkins MTR, Culligan RR, Frasier CL, Dikow RB, Hagenson R, Lei R, et al. Genome sequence and population declines in the critically endangered greater bamboo lemur (Prolemur simus) and implications for conservation. BMC Genomics. 2018;19:1–15.

    Article  CAS  Google Scholar 

  47. Poelstra J, Salmona J, Tiley G, Schüßler D, Blanco M, Andriambeloson J, et al. Cryptic patterns of speciation in cryptic primates: microendemic mouse lemurs and the multispecies coalescent. Syst Biol. 2020.

  48. Chikhi L, Zaonarivelo JR, Jan F, Zaranaina R, Rakotonanahary A, Kun-Rodrigues C, et al. Genetic differentiation and demographic history of the Northern Rufous Mouse Lemur (Microcebus tavaratra) across a fragmented landscape in Northern Madagascar. Int J Primatol. 2018;39:65–89.

    Article  Google Scholar 

  49. Quéméré E, Amelot X, Pierson J, Crouau-Roy B, Chikhi L. Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region. Proc Natl Acad Sci. 2012;109:13028–33.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Beaumont MA. Recent developments in genetic data analysis: what can they tell us about human demographic history? Heredity (Edinb). 2004;92:365–79.

    Article  CAS  Google Scholar 

  51. Chikhi L, Sousa VC, Luisi P, Goossens B, Beaumont MA. The confounding effects of population structure, genetic diversity and the sampling scheme on the detection and quantification of population size changes. Genetics. 2010;186:983–95.

    Article  PubMed  PubMed Central  Google Scholar 

  52. Mazet O, Rodríguez W, Chikhi L. Demographic inference using genetic data from a single individual: separating population size variation from population structure. Theor Popul Biol. 2015;104:46–58. https://doi.org/10.1016/j.tpb.2015.06.003.

    Article  PubMed  Google Scholar 

  53. Mazet O, Rodríguez W, Grusea S, Boitard S, Chikhi L. On the importance of being structured: instantaneous coalescence rates and human evolution-lessons for ancestral population size inference? Heredity (Edinb). 2016;116:362–71.

    Article  CAS  Google Scholar 

  54. Heller R, Chikhi L, Siegismund HR. The confounding effect of population structure on Bayesian skyline plot inferences of demographic history. PLoS ONE. 2013;8:e62992.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Liu X, Fu Y-X. Exploring population size changes using SNP frequency spectra. Nat Genet. 2015;5:555–9.

    Article  CAS  Google Scholar 

  56. Li H, Durbin R. Inference of human population history from individual whole-genome sequences. Nature. 2011;475:493–6. https://doi.org/10.1038/nature10231.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Terhorst J, Song YS. Fundamental limits on the accuracy of demographic inference based on the sample frequency spectrum. Proc Natl Acad Sci U S A. 2015;112:7677–82.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Terhorst J, Kamm JA, Song YS. Robust and scalable inference of population history from hundreds of unphased whole genomes. Nat Genet. 2017;49:303–9.

    Article  CAS  PubMed  Google Scholar 

  59. Lapierre M, Lambert A, Achaz G. Accuracy of demographic inferences from the site frequency spectrum: the case of the yoruba population. Genetics. 2017;206:139–449.

    Article  Google Scholar 

  60. Patton AH, Margres MJ, Stahlke AR, Hendricks S, Lewallen K, Hamede RK, et al. Contemporary demographic reconstruction methods are robust to genome assembly quality: a case study in tasmanian devils. Mol Biol Evol. 2019;36:2906–21.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Chikhi L, Rodríguez W, Grusea S, Santos P, Boitard S, Mazet O. The IICR (inverse instantaneous coalescence rate) as a summary of genomic diversity: insights into demographic inference and model choice. Heredity (Edinb). 2018;120:13–24. https://doi.org/10.1038/s41437-017-0005-6.

    Article  Google Scholar 

  62. Heller R, Nursyifa C, Garcia Erill G, Salmona J, Chikhi L, Meisner J, et al. A reference-free approach to analyze RADseq data using standard Next Generation Sequencing toolkits. Mol Ecol Resour. 2021;21:1085.

    Article  CAS  PubMed  Google Scholar 

  63. Lutermann H, Zimmermann E, Bruford MW, Radespiel U, Schmelting B. Patterns and dynamics of sex-biased dispersal in a nocturnal primate, the grey mouse lemur, Microcebus murinus. Anim Behav. 2003;65:709–19.

    Article  Google Scholar 

  64. Fredsted T, Pertoldi C, Schierup MH, Kappeler PM. Microsatellite analyses reveal fine-scale genetic structure in grey mouse lemurs (Microcebus murinus). Mol Ecol. 2005;14:2363–72.

    Article  CAS  PubMed  Google Scholar 

  65. Radespiel U, Jurić M, Zimmermann E. Sociogenetic structures, dispersal and the risk of inbreeding in a small nocturnal lemur, the golden-brown mouse lemur (Microcebus ravelobensis). Behaviour. 2009;146:607–28.

    Article  Google Scholar 

  66. Radespiel U, Schulte J, Burke RJ, Lehman SM. Molecular edge effects in the Endangered golden-brown mouse lemur Microcebus ravelobensis. Oryx. 2019;53:716–26.

    Article  Google Scholar 

  67. Manel S, Holderegger R. 10 years of landscape genetics. Trends Ecol Evol. 2013;28:614–21.

    Article  PubMed  Google Scholar 

  68. Schliehe-Diecks S, Eberle M, Kappeler PM. Walk the line-dispersal movements of gray mouse lemurs (Microcebus murinus). Behav Ecol Sociobiol. 2012;66:1175–85.

    Article  PubMed  PubMed Central  Google Scholar 

  69. Radespiel U, Lutermann H, Schmelting B, Bruford MW, Zimmermann E. Patterns and dynamics of sex-biased dispersal in a nocturnal primate, the grey mouse lemur, Microcebus murinus. Anim Behav. 2003;65:709–19.

    Article  Google Scholar 

  70. Radespiel U, Rakotondravony R, Chikhi L. Natural and anthropogenic determinants of genetic structure in the largest remaining population of the endangered golden-brown mouse lemur, Microcebus ravelobensis. Am J Primatol. 2008;70:860–70.

    Article  PubMed  Google Scholar 

  71. Aleixo-Pais I, Salmona J, Sgarlata GM, Rakotonanahary A, Sousa AP, Parreira B, et al. The genetic structure of a mouse lemur living in a fragmented habitat in Northern Madagascar. Conserv Genet. 2019;20:229–43. https://doi.org/10.1007/s10592-018-1126-z.

    Article  Google Scholar 

  72. Quéméré E, Crouau-Roy B, Rabarivola C, Louis EE, Chikhi L. Landscape genetics of an endangered lemur (Propithecus tattersalli) within its entire fragmented range. Mol Ecol. 2010;19:1606–21.

    Article  PubMed  Google Scholar 

  73. Salmona J, Teixeira H, Rasolondraibe E, Aleixo-Pais I, Kun-Rodrigues C, Rakotonanahary AN, et al. Genetic diversity, population size, and conservation of the critically endangered Perrier’s Sifaka (Propithecus perrieri). Int J Primatol. 2015;36:1132.

    Article  Google Scholar 

  74. Burney DA. Modern pollen spectra from Madagascar. Palaeogeogr Palaeoclimatol Palaeoecol. 1988;66:63–75.

    Article  Google Scholar 

  75. Railsback LB, Dupont LA, Liang F, Brook GA, Burney DA, Cheng H, et al. Relationships between climate change, human environmental impact, and megafaunal extinction inferred from a 4000-year multi-proxy record from a stalagmite from northwestern Madagascar. Quat Sci Rev. 2020;234: 106244. https://doi.org/10.1016/j.quascirev.2020.106244.

    Article  Google Scholar 

  76. Mietton M, Gunnell Y, Nicoud G, Ferry L, Razafimahefa R, Grandjean P. ‘Lake’ aaotra, Madagascar: a late quaternary wetland regulated by the tectonic regime. CATENA. 2017;2018(165):22–41.

    Google Scholar 

  77. Burney D, James H, Grady F, Rafamantanantsoa J-G, Wright H, Cowart J. Environmental change, extinction and human activity: evidence from caves in NW Madagascar. J Biogeogr. 1997;24:755–67.

    Article  Google Scholar 

  78. Voarintsoa NRG, Wang L, Railsback LB, Brook GA, Liang F, Cheng H, et al. Multiple proxy analyses of a U/Th-dated stalagmite to reconstruct paleoenvironmental changes in northwestern Madagascar between 370 CE and 1300 CE. Palaeogeogr Palaeoclimatol Palaeoecol. 2017;469:138–55. https://doi.org/10.1016/j.palaeo.2017.01.003.

    Article  Google Scholar 

  79. Andriatsitohaina B, Ramsay MS, Kiene F, Lehman SM, Rasoloharijaona S, Rakotondravony R, et al. Ecological fragmentation effects in mouse lemurs and small mammals in northwestern Madagascar. Am J Primatol. 2020;82:1–11.

    Article  Google Scholar 

  80. Thorén S, Linnenbrink M, Radespiel U. Different competitive potential in two coexisting mouse lemur species in northwestern Madagascar. Am J Phys Anthropol. 2011;145:156–62.

    Article  PubMed  Google Scholar 

  81. Olivieri GL, Sousa V, Chikhi L, Radespiel U. From genetic diversity and structure to conservation: genetic signature of recent population declines in three mouse lemur species (Microcebus spp.). Biol Conserv. 2008;141:1257–71.

    Article  Google Scholar 

  82. Craul M, Chikhi L, Sousa V, Olivieri GL, Rabesandratana A, Zimmermann E, et al. Influence of forest fragmentation on an endangered large-bodied lemur in northwestern Madagascar. Biol Conserv. 2009;142:2862–71. https://doi.org/10.1016/j.biocon.2009.05.026.

    Article  Google Scholar 

  83. Sommer S. Effects of habitat fragmentation and changes of dispersal behaviour after a recent population decline on the genetic variability of noncoding and coding DNA of a monogamous Malagasy rodent. Mol Ecol. 2003;12:2845–51.

    Article  CAS  PubMed  Google Scholar 

  84. Johnson JA, Tingay RE, Culver M, Hailer F, Clarke ML, Mindell DP. Long-term survival despite low genetic diversity in the critically endangered Madagascar fish-eagle. Mol Ecol. 2009;18:54–63.

    PubMed  Google Scholar 

  85. Orozco-Terwengel P, Andreone F, Louis E, Vences M. Mitochondrial introgressive hybridization following a demographic expansion in the tomato frogs of Madagascar, genus Dyscophus. Mol Ecol. 2013;22:6074–90.

    Article  CAS  PubMed  Google Scholar 

  86. Ramstad KM, Colbourne RM, Robertson HA, Allendorf FW, Daugherty CH. Genetic consequences of a century of protection: serial founder events and survival of the little spotted kiwi (Apteryx owenii). Proc R Soc B Biol Sci. 2013;280:20130576.

    Article  Google Scholar 

  87. Seixas FA, Juste J, Campos PF, Carneiro M, Ferrand N, Alves PC, et al. Colonization history of Mallorca Island by the European rabbit, Oryctolagus cuniculus, and the Iberian hare, Lepus granatensis (Lagomorpha: Leporidae). Biol J Linn Soc. 2014;111:748–60.

    Article  Google Scholar 

  88. Crowley BE. A refined chronology of prehistoric Madagascar and the demise of the megafauna. Quat Sci Rev. 2010;29:2591–603. https://doi.org/10.1016/j.quascirev.2010.06.030.

    Article  Google Scholar 

  89. Anshari G, Kershaw AP, van der Kaars S. A late Pleistocene and Holocene pollen and charcoal record from peat swamp forest, Lake Sentarum Wildlife Reserve, West Kalimantan, Indonesia. Palaeogeogr Palaeoclimatol Palaeoecol. 2001;171:213–28.

    Article  Google Scholar 

  90. Bush MB, Silman MR. Observations on Late Pleistocene cooling and precipitation in the lowland Neotropics. J Quat Sci. 2004;19:677–84.

    Article  Google Scholar 

  91. Wurster CM, Bird MI, Bull ID, Creed F, Bryant C, Dungait JAJ, et al. Forest contraction in north equatorial Southeast Asia during the Last Glacial Period. Proc Natl Acad Sci. 2010;107:15508–11.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  92. Bohoussou KH, Cornette R, Akpatou B, Colyn M, Kerbis Peterhans J, Kennis J, et al. The phylogeography of the rodent genus Malacomys suggests multiple Afrotropical Pleistocene lowland forest refugia. J Biogeogr. 2015;42:2049–61.

    Article  Google Scholar 

  93. Mizerovská D, Nicolas V, Demos TC, Akaibe D, Colyn M, Denys C, et al. Genetic variation of the most abundant forest-dwelling rodents in Central Africa (Praomys jacksoni complex): Evidence for Pleistocene refugia in both montane and lowland forests. J Biogeogr. 2019;46:1466–78.

    Google Scholar 

  94. Beaumont MA, Zhang W, Balding DJ. Approximate Bayesian computation in population genetics. Genetics. 2002;162:2025–35.

    Article  PubMed  PubMed Central  Google Scholar 

  95. Nielsen R, Williamson S, Kim Y, Hubisz MJ, Clark AG, Bustamante C. Genomic scans for selective sweeps using SNP data. Genome Res. 2005;15:1566–75.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  96. Nater A, Greminger MP, Arora N, Van Schaik CP, Goossens B, Singleton IAN, et al. Reconstructing the demographic history of orang-utans using Approximate Bayesian Computation. Mol Ecol. 2015;24:310–27.

    Article  PubMed  Google Scholar 

  97. Meier JI, Sousa VC, Marques DA, Selz OM, Wagner CE, Excoffier L, et al. Demographic modelling with whole-genome data reveals parallel origin of similar Pundamilia cichlid species after hybridization. Mol Ecol. 2017;26:123–41.

    Article  CAS  PubMed  Google Scholar 

  98. Kumar S, Stecher G, Suleski M, Hedges SB. TimeTree: a resource for timelines, timetrees, and divergence times. Mol Biol Evol. 2017;34:1812–9.

    Article  CAS  PubMed  Google Scholar 

  99. Sehen L, Goetze D, Rajeriarison C, Roger E, Thoren S, Radespiel U, et al. Structural and floristic traits of habitats with differing relative abundance of the lemurs Microcebus murinus and M. ravelobensis in northwestern Madagascar. Ecotropica. 2010;16:15–30.

    Google Scholar 

  100. Radespiel U, Sarikaya Z, Zimmermann E, Bruford MW. Sociogenetic structure in a free-living nocturnal primate population: sex-specific differences in the grey mouse lemur (Microcebus murinus). Behav Ecol Sociobiol. 2001;50:493–502.

    Article  Google Scholar 

  101. Radespiel U, Ehresmann P, Zimmermann E. Species-specific usage of sleeping sites in two sympatric mouse lemur species (Microcebus murinus and M. ravelobensis) in northwestern Madagascar. Am J Primatol. 2003;59:139–51.

    Article  PubMed  Google Scholar 

  102. Rina Evasoa M, Zimmermann E, Hasiniaina AF, Rasoloharijaona S, Randrianambinina B, Radespiel U. Sources of variation in social tolerance in mouse lemurs (Microcebus spp.). BMC Ecol. 2019;19:1–16. https://doi.org/10.1186/s12898-019-0236-x.

    Article  CAS  Google Scholar 

  103. Zimmermann E, Cepok S, Rakotoarison N, Zietemann V, Radespiel U. Sympatric Mouse Lemurs in North-West Madagascar: a new rufous Mouse Lemur species (Microcebus ravelobensis). Folia Primatol. 1998;69:106–14.

    Article  CAS  Google Scholar 

  104. Seutin G, White BN, Boag PT. Preservation of avian blood and tissue samples for DNA analyses. Can J Zool. 1991;69:82–90.

    Article  CAS  Google Scholar 

  105. Davey JL, Blaxter MW. RADseq: next-generation population genetics. Brief Funct Genomics. 2010;9:416–23.

    Article  CAS  PubMed  Google Scholar 

  106. Davey JW, Hohenlohe PA, Etter PD, Boone JQ, Catchen JM, Blaxter ML. Genome-wide genetic marker discovery and genotyping using next-generation sequencing. Nat Rev Genet. 2011;12:499–510. https://doi.org/10.1038/nrg3012.

    Article  CAS  PubMed  Google Scholar 

  107. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  108. Larsen PA, Harris RA, Liu Y, Murali SC, Campbell CR, Brown AD, et al. Hybrid de novo genome assembly and centromere characterization of the gray mouse lemur (Microcebus murinus). BMC Biol. 2017;15:1–17.

    Article  CAS  Google Scholar 

  109. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  110. Ebbert MTW, Wadsworth ME, Staley LA, Hoyt KL, Pickett B, Miller J, et al. Evaluating the necessity of PCR duplicate removal from next-generation sequencing data and a comparison of approaches. BMC Bioinformatics. 2016;17:491–500.

    Article  Google Scholar 

  111. Schmieder R, Edwards R. Quality control and preprocessing of metagenomic datasets. Bioinformatics. 2011;27:863–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  112. Li H, Durbin R. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics. 2009;25:1754–60. https://doi.org/10.1093/bioinformatics/btp324.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  113. Metzker ML. Sequencing technologies—the next generation. Nat Rev Genet. 2010;11:31.

    Article  CAS  PubMed  Google Scholar 

  114. Skotte L, Korneliussen TS, Albrechtsen A. Estimating individual admixture proportions from next generation sequencing data. Genetics. 2013;195:693–702.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  115. Korneliussen TS, Albrechtsen A, Nielsen R. ANGSD: analysis of next generation sequencing data. BMC Bioinform. 2014;15:1–13.

    Article  Google Scholar 

  116. Korneliussen TS, Albrechtsen A, Nielsen R. Open access ANGSD: analysis of next generation sequencing data. BMC Bioinform. 2014;15:1–13.

    Article  Google Scholar 

  117. Korneliussen TS, Moltke I. Sequence analysis NgsRelate: a software tool for estimating pairwise relatedness from next-generation sequencing data. Bioinformatics. 2015;2015(31):4009–11.

    Google Scholar 

  118. Soraggi S, Wiuf C, Albrechtsen A. Powerful inference with the D-statistic on low-coverage whole-genome data. G3 Genes, Genomes, Genet. 2017;8:551–66.

    Google Scholar 

  119. Bagley RK, Sousa VC, Niemiller ML, Linnen CR. History, geography and host use shape genomewide patterns of genetic variation in the redheaded pine sawfly (Neodiprion lecontei). Mol Ecol. 2017;26:1022–44.

    Article  CAS  PubMed  Google Scholar 

  120. Parreira BR, Chikhi L. On some genetic consequences of social structure, mating systems, dispersal, and sampling. Proc Natl Acad Sci U S A. 2015;112:E3318–26.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  121. Parreira B, Quéméré E, Vanpé C, Carvalho I, Chikhi L. Genetic consequences of social structure in the golden-crowned sifaka. Heredity (Edinb). 2020;125:1–12.

    Article  Google Scholar 

  122. Blouin MS. DNA-based methods for pedigree reconstruction and kinship analysis in natural populations. Trends Ecol Evol. 2003;18:503–11.

    Article  Google Scholar 

  123. Weir BS, Cockerham CC. Estim F-statistics. Anal Popul Struct Evol. 1984;38:1358–70.

    CAS  Google Scholar 

  124. Excoffier L, Lischer HEL. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour. 2010;10:564–7.

    Article  PubMed  Google Scholar 

  125. Meisner J, Albrechtsen A. Inferring population structure and admixture proportions in low-depth NGS data. Genetics. 2018;210:719–31.

    Article  PubMed  PubMed Central  Google Scholar 

  126. Kopelman NM, Mayzel J, Jakobsson M, Rosenberg NA, Mayrose I. Clumpak: a program for identifying clustering modes and packaging population structure inferences across K. Mol Ecol Resour. 2015;15:1179–91.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  127. Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol Ecol. 2005;14:2611–20.

    Article  CAS  PubMed  Google Scholar 

  128. Wright S. Evolution in Mendelian populations. Genetics. 1931;16:97.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  129. Rousset F. Genetic differentiation between individuals. J Evol Biol. 2000;13:58–62.

    Article  Google Scholar 

  130. Hardy OJ, Vekemans X. SPAGeDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. Mol Ecol Notes. 2002;2:618–20.

    Article  CAS  Google Scholar 

  131. Mantel N. The detection of disease clustering and a generalized regression approach. Cancer Res. 1967;27(2 Part 1):209–20.

    CAS  PubMed  Google Scholar 

  132. Dixon P. VEGAN, a package of R functions for community ecology. J Veg Sci. 2003;14:927–30.

    Article  Google Scholar 

  133. Vieira FG, Fumagalli M, Albrechtsen A, Nielsen R. Estimating inbreeding coefficients from NGS data: impact on genotype calling and allele frequency estimation. Genome Res. 2013;23:1852–61.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  134. Yoder AD, Campbell CR, Blanco MB, Ganzhorn JU, Goodman SM. Geogenetic patterns in mouse lemurs (genus Microcebus) reveal the ghosts of Madagascar ’ s forests past. PNAS. 2016;113:8049–56.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  135. Campbell CR, Tiley GP, Poelstra JW, Hunnicutt KE, Larsen PA, dos Reis M, et al. Pedigree-based measurement of the de novo mutation rate in the gray mouse lemur reveals a high mutation rate, few mutations in CpG sites, and a weak sex bias. bioRxiv. 2019;724880.

  136. Gutenkunst RN, Hernandez RD, Williamson SH, Bustamante CD. Inferring the joint demographic history of multiple populations from multidimensional SNP frequency data. PLoS Genet. 2009;5:e1000695.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  137. Foote AD, Vijay N, Ávila-Arcos MC, Baird RW, Durban JW, Fumagalli M, et al. Genome-culture coevolution promotes rapid divergence of killer whale ecotypes. Nat Commun. 2016. https://doi.org/10.1038/ncomms11693.

    Article  PubMed  PubMed Central  Google Scholar 

  138. Arredondo A, Mourato B, Nguyen K, Boitard S, Rodriguez W, Noûs C, et al. Inferring number of populations and changes in connectivity under the n-island model. Heredity (Edinb). 2021;126:896–912.

    Article  Google Scholar 

  139. Hudson RR. Generating samples under a Wright-Fisher neutral model of genetic variation. Bioinformatics. 2002;18:337–8.

    Article  CAS  PubMed  Google Scholar 

  140. Patil AB, Vijay N. Repetitive genomic regions and the inference of demographic history. Heredity (Edinb). 2021. https://doi.org/10.1038/s41437-021-00443-8.

    Article  Google Scholar 

  141. Schüßler D, Mantilla-Contreras J, Stadtmann R, Ratsimbazafy JH, Radespiel U. Identification of crucial stepping stone habitats for biodiversity conservation in northeastern Madagascar using remote sensing and comparative predictive modeling. Biodivers Conserv. 2020;29(7):2161–84.

    Article  Google Scholar 

  142. UNEP-WCMC and IUCN. Protected planet: the World Database on Protected Areas (WDPA); 2020. Cambridge, UK: UNEP-WCMC and IUCN; 2020. www.protectedplanet.net.

Download references

Acknowledgements

We thank Madagascar National Parks and the Direction du Système des Aires Protégées (DSAP), the Direction Générale du Ministère de l’Environnement et des Forêts de Madagascar, Madagascar’s Ad Hoc Committee for Fauna and Flora, and the Committee for Environmental Research (CAFF/CORE) for the permission to conduct the field work in the Ankarafantsika National Park. We thank Madame Jacqueline, Planet Madagascar and Durrell Wildlife Conservation Trust for their help concerning field work logistics. The field work was possible thanks to the help of Adolphe Rakotomanantena, Augustin Rakotonandrianina, Hauke Henkel, the guides from ANP and thanks to the welcoming local community. We thank Ariel Rodriguez, Jordan Silvester and Jörn Wrede for their help with computational issues and Felix Pellerin for his help and guidance with custom scripts. We are grateful to the Get-Plage platform for library preparation and sequencing, to the Genotoul bioinformatics (Bioinfo Genotoul), and to the University of Veterinary Medicine Hannover for providing computing resources. We are also grateful to Dominik Schüßler and Sönke von den Berg for their help with Fig. 1, and to Jeffrey Rogers and Alan R. Harris for providing the RepeatMasker output file for the M. murinus reference genome. Finally, we acknowledge the two anonymous reviewers whose comments and suggestions improved the quality of our manuscript substantially.

Funding

Open Access funding enabled and organized by Projekt DEAL. This study was supported by the Deutsche Forschungsgemeinschaft (DFG Ra 502/20-1, Ra 502/20-3, ME 4517/2-1) and H. Teixeira was additionally funded by a scholarship of the Hannover Graduate School for Veterinary Pathobiology, Neuroinfectiology and Translational Medicine (HGNI). The work also benefited from several meetings with collaborators organized under the ERA-NET BiodivERsA initiative (No. 2015-138), project: INFRAGECO (Inference, Fragmentation, Genomics, and Conservation), financed on the German side by the German Bundesministerium für Bildung und Forschung (#01LC1617A).

Author information

Authors and Affiliations

  1. Institute of Zoology, University of Veterinary Medicine Hannover, Foundation, Bünteweg 17, 30559, Hannover, Germany

    Helena Teixeira & Ute Radespiel

  2. Laboratoire Évolution and Diversité Biologique (EDB UMR 5174), Université de Toulouse Midi-Pyrénées, CNRS, IRD, UPS, 118 Route de Narbonne, Bât. 4R1, 31062, Toulouse cedex 9, France

    Jordi Salmona, Sophie Manzi & Lounès Chikhi

  3. Instituto Gulbenkian de Ciência, Rua da Quinta Grande, 6, 2780-156, Oeiras, Portugal

    Armando Arredondo, Beatriz Mourato & Lounès Chikhi

  4. Université de Toulouse, Institut National des Sciences Appliquées, Institut de Mathématiques de Toulouse, Toulouse, France

    Armando Arredondo & Olivier Mazet

  5. Ecole Doctorale Ecosystèmes Naturels (EDEN), University of Mahajanga, 5 Rue Georges V - Immeuble KAKAL, Mahajanga Be, B.P. 652, 401, Mahajanga, Madagascar

    Romule Rakotondravony

  6. Faculté des Sciences, de Technologies et de l’Environnement, University of Mahajanga, 5 Rue Georges V - Immeuble KAKAL, Mahajanga Be, B.P. 652, 401, Mahajanga, Madagascar

    Romule Rakotondravony

  7. Institute of Animal Breeding and Genetics, University of Veterinary Medicine Hannover, Foundation, Bünteweg 17p, 30559, Hannover, Germany

    Julia Metzger

  8. Veterinary Functional Genomics, Max Planck Institute for Molecular Genetics, Ihnestrasse 73, 14195, Berlin, Germany

    Julia Metzger

Authors
  1. Helena Teixeira
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Jordi Salmona
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Armando Arredondo
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Beatriz Mourato
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Sophie Manzi
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Romule Rakotondravony
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Olivier Mazet
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Lounès Chikhi
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Julia Metzger
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Ute Radespiel
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

UR and HT conceptualized and designed the study and planned the field expedition together with RR. HT collected the mouse lemur samples for the genomic analyses, carried out the initial laboratory work, processed the RADseq data with guidance of JS and UR, and wrote the first draft of the manuscript. SM prepared the RADseq libraries and was involved in all laboratorial procedures related to the RAD sequencing. JS developed the bioinformatics pipeline used to analyze the raw RADseq data and provided supervision during genomic analyses. JM performed the whole-genomic sequencing and the PSMC analyses. AA and BM conducted the IICR-simulations with the guidance of OM and LC. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Helena Teixeira or Ute Radespiel.

Ethics declarations

Ethics approval and consent to participate

The samples used in these study stem from wild mouse lemur populations. Therefore, sample collection took place in Madagascar. Our research proposal was submitted to the Malagasy research and conservation authorities (Madagascar National Parks, Ministère de l’Environnement et des Forêts de Madagascar and Committee for Environmental Research) before the field work, and all scientific methods and ethical aspects of trapping and sampling were approved by them (Permission Number: No 78/17/MEEF/SG/DGF/DSAP/SCB.Re). The field work took place as initially planned and it was in full agreement with the legal requirements of Madagascar and with the ethical guidelines of the International Council of Laboratory Animal Science (ICLAS), the International Union for Conservation of Nature (IUCN), Policy Statement on Research Involving Species at Risk of Extinction, and the Principles for the Ethical Treatment of Non-Human Primates of the American Society of Primatologists. All capture and handling procedures followed routine protocols and were approved by the Malagasy Authorities and by the Institute of Zoology, University of Veterinary Medicine Hannover.

Consent for publication

Mouse lemur illustrations copyright 2013 Stephen D. Nash / IUCN SSC Primate Specialist Group. Used with permission.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Supplementary text (Text S1–Text S7), figures (Figure S1–Figure S10) and tables (Table S1–Table S10).

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Teixeira, H., Salmona, J., Arredondo, A. et al. Impact of model assumptions on demographic inferences: the case study of two sympatric mouse lemurs in northwestern Madagascar. BMC Ecol Evo 21, 197 (2021). https://doi.org/10.1186/s12862-021-01929-z

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12862-021-01929-z

Keywords

BMC Ecology and Evolution

ISSN: 2730-7182

Contact us