George P Tiley
email: g.tiley@kew.org

Tobias Van Elst

Helena Teixeira

Dominik Schüßler

Jordi Salmona

Marina B Blanco

José M Ralison

Blanchard Randrianambinina

Rodin M Rasoloarison

Amanda R Stahlke

| Tobias Van Elst

| Dominik Schüßler

| Blanchard Randrianambinina

| Paul

A Hohenlohe

| Lounès Chikhi

Edward E Louis

Ute Radespiel

Anne D Yoder

Paul A Hohenlohe

Population genomic structure in Goodman's mouse lemur reveals long-standing separation of Madagascar's Central Highlands and eastern rainforests

établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

| INTRODUC TI ON

Madagascar is an iconic biodiversity hotspot with high levels of species diversity and endemism [START_REF] Myers | Biodiversity hotspots for conservation priorities[END_REF]. The majority of its biodiversity is at extreme risk of extinction due to deforestation and forest fragmentation [START_REF] Green | Deforestation history of the eastern rain forests of Madagascar from satellite images[END_REF][START_REF] Harper | Fifty years of deforestation and forest fragmentation in Madagascar[END_REF][START_REF] Vieilledent | Combining global tree cover loss data with historical national forest cover maps to look at six decades of deforestation and forest fragmentation in Madagascar[END_REF]. It is widely accepted that the rapid rate of biodiversity loss is due largely to human activity and agricultural practices [START_REF] Farris | Threats to a rainforest carnivore community: A multi-year assessment of occupancy and co-occurrence in Madagascar[END_REF][START_REF] Styger | Influence of slash-and-burn farming practices on fallow succession and land degradation in the rainforest region of Madagascar[END_REF]. Even so, the timing of and degree to which humans transformed and fragmented Madagascar's forests have been points of contention.

Much of central Madagascar is now covered with grasslands, traditionally thought to be strictly the product of fires from slash-andburn agriculture and cattle farming [START_REF] Gade | Deforestation and its effects in highland Madagascar[END_REF][START_REF] Humbert | La destruction d'une flore insulaire par le feu. Principaux aspects de la vegetation à Madagascar[END_REF] Perrier de la [START_REF] Perrier De La Bâthie | La vegetation Malgache[END_REF]. The best available archaeological evidence places human arrival in Madagascar between 2000 years ago (ka), consistent with human population genomics [START_REF] Pierron | Genomic landscape of human diversity across Madagascar[END_REF], and 10 ka [START_REF] Dewar | Stone tools and foraging in northern Madagascar challenge Holocene extinction models[END_REF][START_REF] Dewar | The culture history of Madagascar[END_REF][START_REF] Gommery | Les plus anciennes traces d'activités anthropiques de Madagascar sur des ossements d'hippopotames subfossiles d'Anjohibe (Province de Mahajanga)[END_REF][START_REF] Hansford | Early Holocene human presence in Madagascar evidenced by exploitation of avian megafauna[END_REF][START_REF] Mitchell | Settling Madagascar: when did people first colonize the world's largest Island?[END_REF]. This timeline coincides with evidence of increased frequency of fire based on isotope records from stalagmites [START_REF] Railsback | Relationships between climate change, human environmental impact, and megafaunal extinction inferred from a 4000-year multi-proxy record from a stalagmite from northwestern Madagascar[END_REF] and charcoal records from sediment cores covering the late to mid Holocene (Burney, 1987a;[START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF][START_REF] Virah-Sawmy | Evidence for drought and forest declines during the recent megafaunal extinctions in Madagascar[END_REF]. Pollen abundances from these same sediment cores also imply that the Central Highland vegetation transitioned from woody taxa to grasses during this timeframe, supporting the hypothesis that the grassy dominance of the Central Highlands is the consequence of human activities (reviewed in [START_REF] Burney | A chronology for late prehistoric Madagascar[END_REF][START_REF] Joseph | How expansive were Malagasy Central Highland forests, ericoids, woodlands and grasslands? A multidisciplinary approach to a conservation conundrum[END_REF]. However, it becomes difficult to disentangle anthropogenic and paleoclimatic effects given that the transition to grasslands probably began during the climatic fluctuations of the Pleistocene [START_REF] Burney | Late Quaternary stratigraphic charcoal records from Madagascar[END_REF][START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF], long before the arrival of humans. Changes in plant composition from woody taxa to grasses have also been reported outside of the Central Highlands in southwestern [START_REF] Burney | Late Holocene environmental changes in arid Southwestern Madagasar[END_REF] and southeastern [START_REF] Virah-Sawmy | Threshold response of Madagascar's littoral forest to sea-level rise[END_REF] Madagascar. Taken together, the paleoecological data imply that a grassy Central Highland vegetation may not strictly be the product of human agricultural practices, but instead, the consequence of global climate fluctuations, perhaps facilitated by the natural processes of fire ecology [START_REF] Bowman | Fire ecology[END_REF].

In support of this second scenario, ecological and phylogenetic studies of Madagascar's grasses indicate that Central Highland grasslands originated well before the Holocene, with levels of endemism comparable to the grasslands of South Africa [START_REF] Bond | The antiquity of Madagascar's grasslands and the rise of C 4 grassy biomes[END_REF] and tropical east Africa [START_REF] Vorontsova | Madagascar's grasses and grasslands: anthropogenic or natural[END_REF] where ancient origins are well-supported. And though few C 4 grass lineages such as Aristida and Loudetia dominate the Central Highlands [START_REF] Koechlin | Natural grasslands: eastern hemisphere and resume[END_REF][START_REF] Kull | Fire and the management of highland vegetation[END_REF][START_REF] Solofondranohatra | Grass functional traits differentiate forest and savanna in the Madagascar Central Highlands[END_REF], diversity among endemic C 4 grass clades would be more consistent with the global diversification of C 4 grasses during the Miocene (e.g., [START_REF] Edwards | Phylogenetic analyses reveal the shady history of C 4 grasses[END_REF] rather than as a consequence of many independent migration events [START_REF] Hackel | Grass diversification in Madagascar: In situ radiation of two large C 3 shade clades and support for a Miocene to Pliocene origin of C 4 grassy biomes[END_REF]. We are therefore left with a complex model of Pleistocene climate change, perhaps driving grassland expansions, but overlaid with deforestation and fragmentation, probably associated with and intensified by early human colonization. We hypothesize that present-day distributions of forest-dwelling species, their genetic diversity, and connectivity among populations will reflect both ancient climatic change as well as recent human-mediated processes. Here, we perform a population genomic investigation of Goodman's mouse lemur (Microcebus lehilahytsara) to test whether signatures of pre-and post-human vegetative change in Madagascar are detectable with extant molecular data and demographic models.

The geographic distribution of M. lehilahytsara is uniquely suited for this purpose. It is the only mouse lemur species found in isolated Central Highland forest patches [START_REF] Andriambeloson | Living in tiny fragments: a glimpse at the ecology of Goodman's mouse lemurs (Microcebus lehilahytsara) in the relic forest of Ankafobe, Central Highlands, Madagascar[END_REF][START_REF] Yoder | Geogenetic patterns in mouse lemurs (genus Microcebus) reveal the ghosts of Madagascar;s forests past[END_REF] as well as in the humid rainforests along the escarpment of Madagascar's east coast [START_REF] Blanco | Evidence of prolonged torpor in Goodman's mouse lemurs at Ankafobe forest, central Madagascar. Primates[END_REF][START_REF] Radespiel | First indications of a highland specialist among mouse lemurs (Microcebus spp.) and evidence for a new mouse lemur species from eastern Madagascar[END_REF].

It is also found in the rainforests of northeastern Madagascar, as the formerly recognized M. mittermeieri was recently recommended as a junior synonym of M. lehilahytsara based on both morphological [START_REF] Schüßler | Ecology and morphology of mouse lemurs (Microcebus spp.) in a hotspot of microendemism in northeastern Madagascar, with the description of a new species[END_REF] and molecular evidence (Poelstra Though our findings support the hypothesis that Central Highland grasslands predate human arrival, we also find evidence for human-mediated population declines. This highlights the extent to which species imminently threatened by human-mediated deforestation may already be vulnerable from paleoclimatic conditions.

K E Y W O R D S

Central Highland Savanna, coalescent, fragmentation, historical demography, Microcebus, RAD sequencing et al., 2021). Of additional paleoecological interest, M. lehilahytsara is found at higher altitudes than other mouse lemur species, perhaps due to its ability to enter a state of prolonged torpor [START_REF] Blanco | Evidence of prolonged torpor in Goodman's mouse lemurs at Ankafobe forest, central Madagascar. Primates[END_REF], thus reducing its metabolic demands during periods of extended low temperatures during the austral winter. To summarize, the species distribution and life history of M. lehilahytsara suggest that it has long been present in the Central Highlands and that demographic studies of these forest-dwelling primates may illuminate prehuman conditions in Madagascar. Population divergences can reflect forest fragmentation, and population size fluctuations through time will reflect pulses in forest distributions. Although loss of genetic diversity in endangered lemurs species can be attributed to human activities [START_REF] Olivieri | From genetic diversity and structure to conservation: Genetic signature of recent population declines in three mouse lemur species (Microcebus spp.)[END_REF][START_REF] Salmona | Climate change and human colonization triggered habitat loss and fragmentation in Madagascar[END_REF][START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF], a better informed reconstruction of prehuman conditions can reveal if some populations are especially vulnerable due to the effects of natural climate change [START_REF] Quéméré | Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region[END_REF][START_REF] Salmona | Climate change and human colonization triggered habitat loss and fragmentation in Madagascar[END_REF][START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF].

| MATERIAL S AND ME THODS

| Study sites and sampling

Microcebus lehilahytsara individuals were captured at five field sites (Figure 1). Two locations, Ambohitantely and Ankafobe, are located within the Central Highlands. Their vegetational components are characteristic of moist evergreen forests, although more open and drier. Both Ambohitantely and Ankafobe are protected forested areas well-characterized by the extensive grasslands surrounding them.

The three other locations, Anjozorobe, Ambatovy, and Tsinjoarivo, are mid-elevation rainforests at the southern end of the species range situated along the eastern escarpment (hereafter referred as southern [rainforest] populations; Figure S1). We complemented our sampling with an existing data set of five M. lehilahytsara populations at humid rainforest sites in northeastern Madagascar (Marojejy, Anjanaharibe-Sud, Anjiahely, Ambavala, Riamalandy) reported in Poelstra et al. (2021) (hereafter referred as northern [rainforest] populations), where it was recommended that M. mittermeieri should be synonymized as M. lehilahytsara. Consequently, our sample set included 54 M. lehilahytsara individuals covering the entire known geographic range of the species. A detailed description of climate and topography of the sampling sites can be found in the Supporting Information Methods "Study Sites". In addition, we included seven M. tanosi individuals from southeastern Madagascar and three M. murinus individuals from western Madagascar as outgroups. For a description of live-capture methods see [START_REF] Blanco | Evidence of prolonged torpor in Goodman's mouse lemurs at Ankafobe forest, central Madagascar. Primates[END_REF]. Collection information for all specimens is given in Table S1.

| RADseq library preparation and sequencing

We generated new RADseq data for 38 individuals while previously published sequences using the same methods were available for 26 individuals [START_REF] Poelstra | Cryptic patterns of speciation in cryptic primates: Microendemic mouse lemurs and the multispecies coalescent[END_REF]. RADseq library preparations generally followed [START_REF] Ali | Rad capture (Rapture): Flexible and efficient sequence-based genotyping[END_REF], using a single-digest protocol S2.

| Data processing and genotyping

Raw RAD reads were demultiplexed with the process_radtags function in STACKS version 2.0b [START_REF] Rochette | Stacks 2: Analytical methods for paired-end sequencing improve RADseqbased population genomics[END_REF], trimmed using Trimmomatic with the parameters Leading: 3, Trailing: 3, Slidingwindow: 4:15, and Minlen: 60 [START_REF] Bolger | Trimmomatic: a flexible trimmer for Illumina sequence data[END_REF], and aligned against the Microcebus murinus reference genome (Mmur 3.0; [START_REF] Larsen | Hybrid de novo genome assembly and centromere characterization of the gray mouse lemur (Microcebus murinus)[END_REF] using BWA version 7.0.17 with the BWA-MEM algorithm and default settings [START_REF] Li | Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM[END_REF]. Trimmed reads of technical replicates, that were available for 16 samples (Tables S2), were combined prior to reference alignment to increase coverage. Subsequently, BAM files were filtered with SAMtools version 1.11 [START_REF] Li | The Sequence Alignment/Map format and SAMtools[END_REF] to retain only reads mapping to autosomal scaffolds with a minimum mapping quality of 20. Paired-end sequenced data was additionally filtered for proper pairing and deduplicated. SAMtools was also used to count the number of RAD loci sequenced for each individual and estimate locus coverage as the forward read depth at the respective SbfI cutting site.

Two genotyping approaches were employed to evaluate robustness of our results to genotyping errors. First, we used the referencebased approach of STACKS version 2.53 [START_REF] Rochette | Stacks 2: Analytical methods for paired-end sequencing improve RADseqbased population genomics[END_REF], including only loci present in at least 75% of all individuals. Genotype calls were filtered using VCFtools version 0.1.17 [START_REF] Danecek | The variant call format and VCFtools[END_REF] following "FS6" recommendations of O' [START_REF] Oksanen | These aren't the loci you'e looking for: Principles of effective SNP filtering for molecular ecologists[END_REF] with modified thresholds (Methods "Genotype filtering"). A minor allele count filter of three was applied for analyses of population structure.

Second, genotype likelihoods were estimated with the SAMtools model in ANGSD version 0.934 [START_REF] Korneliussen | ANGSD: Analysis of next generation sequencing data[END_REF], using filtering thresholds applied in closely-related mouse lemurs [START_REF] Poelstra | Cryptic patterns of speciation in cryptic primates: Microendemic mouse lemurs and the multispecies coalescent[END_REF]) and excluding M. lehilahytsara individuals that did not pass genotype filtering and outgroups (Supporting Information Methods "Genotype likelihood estimation with ANGSD"). Genotype likelihoods were used for analyses of population structure and to infer minor allele frequency (MAF) spectra for all populations and population pairs (i.e., joint MAF spectra) for demographic analyses, using realSFS implemented in the ANGSD framework.

| Phylogenetic inference

We estimated phylogenies using three different approaches with genotype calls, which allowed us to explore uncertainty in phylogeographic history due to model assumptions for each method:

1. Maximum likelihood (ML; [START_REF] Felsenstein | Evolutionary trees from DNA sequences: A maximum likelihood approach[END_REF] analysis with RAxML-NG [START_REF] Kozlov | RAxML-NG: A fast, scalable and user-friendly tool for maximum likelihood phylogenetic inference[END_REF] was carried out on the unpartitioned concatenated alignment of SNPs using the GTR + Γ model of molecular evolution [START_REF] Tavaré | Some probabilistic and statistical problems in the analysis of DNA sequences[END_REF][START_REF] Yang | Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: Approximate methods[END_REF] while correcting for ascertainment bias [START_REF] Lewis | A likelihood approach to estimating phylogeny from discrete morphological character data[END_REF]. We conducted 20 unconstrained ML searches and performed 100 bootstrap replicates.

2. Topological estimation from invariants was performed with SVDquartets [START_REF] Chifman | Quartet Inference from SNP data under the coalescent model[END_REF] implemented through PAUP* version 4.0a (build 168;[START_REF] Swofford | PAUP*. Phylogenetic Analysis Using Parsimony (*and Other Methods)[END_REF], using a thinned data set with a minimum distance of 10,000 bp between SNPs to potentially avoid treating SNPs in linkage disequilibrium as independent. We evaluated all quartets with 100 standard bootstrap replicates and estimated phylogenies where each individual was a tip and where individuals were assigned to populations a priori based on sampling location.

3. SNAPP [START_REF] Bryant | Inferring species trees directly from biallelic genetic markers: Bypassing gene trees in a full coalescent analysis[END_REF], a Bayesian implementation of the multispecies coalescent (MSC; [START_REF] Rannala | Bayes estimation of species divergence times and ancestral population sizes using DNA sequences from multiple loci[END_REF] model for SNP data, was run using BEAST2 version 2.4.8 [START_REF] En Bouckaert | BEAST 2: A software platform for Bayesian evolutionary analysis[END_REF]. Similar to SVDquartets, we used a thinned data set to avoid site independence violations. To reduce the computational burden, we selected only two individuals per sampling location based on mean sequencing depth (Tables S1 andS2). Four replicates were run to evaluate convergence. Mixing and effective sample sizes were analysed with Tracer version 1.7 [START_REF] Rambaut | Posterior summarization in Bayesian pylogenetics using Tracer 1.7[END_REF]. Details regarding priors and MCMC options are in the Supporting Information Methods "SNAPP priors and settings".

To evaluate alternative topological hypotheses for ML and MSC methods, we used probabilistic topology hypothesis testing approaches. Specifically, we tested if the earliest split among M. lehilahytsara populations was between northern rainforest populations and southern rainforest plus Central Highland populations or between the Central Highlands and all rainforest populations. For ML, constrained optimizations of topologies and site likelihoods were performed with RAxML-NG, and support between the two competing topologies was evaluated with an approximately unbiased (AU; [START_REF] Shimodaira | Approximately unbiased tests of regions using multistep-multiscale bootstrap resampling[END_REF] test in IQ-TREE version 2.1.3 [START_REF] Minh | IQ-TREE 2: New models and efficient methods for phylogenetic inference in the genomic era[END_REF], using the GTR + Γ model and 10,000 RELL [START_REF] Kishino | Maximum likelihood inference of protein phylogeny and the origin of chloroplasts[END_REF] replicates. For SNAPP, we used stepping-stone integration [START_REF] Xie | Improving marginal likelihood estimation for Bayesian phylogenetic model selection[END_REF] to obtain marginal likelihoods with 24 steps. Because differences in marginal likelihoods were large, we presented natural log Bayes factors [START_REF] Kass | Bayes factors[END_REF] as opposed to computing model probabilities.

| Population structure

Population structure was evaluated through clustering analyses, principal component analysis (PCA), and fixation index (F ST ) statistics, leveraging both genotype calls and genotype likelihoods. Individual ancestries were inferred from genotype calls with STRUCTURE version 2.3.4 [START_REF] Pritchard | Inference of population structure using multilocus genotype data[END_REF] and from genotype likelihoods with ngsAdmix version 32 [START_REF] Skotte | Estimating individual admixture proportions from next generation sequencing data[END_REF]. We conducted 20 independent runs for genotype calls and 10 runs for genotype likelihoods with the number of a priori clusters (K) ranging from 1-10 and calculated ΔK for estimating the most likely number of Hardy-Weinberg groups following [START_REF] Evanno | Detecting the number of clusters of individuals using the software structure: A simulation study[END_REF]. PCA was conducted on genotype calls using the R package adegenet version 2.1.3 [START_REF] Jombart | adegenet 1.3-1: New tools for the analysis of genome-wide SNP data[END_REF] and on genotype likelihoods using PCAngsd version 1.01 [START_REF] Meisner | Inferring population structure and admixture proportions in low-depth NGS data[END_REF]. Pairwise F ST between populations were calculated from genotype calls as described by [START_REF] Weir | Estimating F-statistics for the analysis of population structure[END_REF] using the R package hierfstat version 0.5-7 [START_REF] Goudet | hierfstat, a package for r to compute and test hierarchical F-statistics[END_REF] and from genotype likelihood-based joint MAF spectra using realSFS with weighted F ST calculations [START_REF] Reynolds | Estimation of the coancestry coefficient: basis for a short-term genetic distance[END_REF]. The populations Riamalandy and Anjozorobe were excluded for F ST estimation because only two and one individuals were available for them, respectively. Mantel tests were conducted at the population level using the R package vegan version 2.5-7 [START_REF] Oksanen | These aren't the loci you'e looking for: Principles of effective SNP filtering for molecular ecologists[END_REF] to test the presence of isolationby-distance using 9,999 permutations. We tested three alternative sets of populations: (1) northern rainforests only, (2) southern rainforests plus Central Highlands, and (3) all populations. Physical distances between populations were calculated with the Geographic Distance Matrix Generator version 1.2.3 [START_REF] Ersts | Geographic distance matrix generator[END_REF]. Pairwise differences between individuals were used to perform an AMOVA analysis [START_REF] Excoffier | Analysis of molecular variance inferred from metric distances among DNA haplotypes: Application to human mitochondrial DNA restriction data[END_REF] with the R package poppr version 2.8.7 [START_REF] Kamvar | Poppr: An R package for genetic analysis of populations with clonal, partially clonal, and/or sexal reproduction[END_REF], grouping populations into the northern rainforests, southern rainforests, and Central Highlands to evaluate sources of variation.

| Demographic modelling

Genotype likelihoods can alleviate known biases of RADseq data and improve demographic inference from allele frequencies [START_REF] Heller | A reference-free approach to analyze non-model RADseq data using standard Next Generation Sequencing toolkits[END_REF][START_REF] Pedersen | A southern African origin and cryptic structure in the highly mobile plains zebra[END_REF][START_REF] Warmuth | Genotype-free estimation of allele frequencies reduces bias and improves demographic inference from RADSeq data[END_REF] because they include information about the uncertainty associated with genotype calls. Especially for low coverage data it may be important to account for this uncertainty, as stochastic under-and overrepresentation of alleles can impair genotype calling [START_REF] Nielsen | Genotype and ANP calling from next-generation sequencing data[END_REF]. Therefore, joint MAF spectra estimated from genotype likelihoods were used to test models of population splitting, population size change, and gene flow with fastsimcoal2 version 2.6.0.3 [START_REF] Excoffier | Robust demographic inference from genomic and SNP data[END_REF]. For each model, we used 100 independent approximations with 200,000 simulations and 40 cycles of the expectation maximization algorithm. The best likelihood for each model was used to calculate its Akaike information criterion (AIC; [START_REF] Akaike | A new look at the statistical model identification[END_REF] score. We used the AIC scores to calculate model probabilities such that model w has probability exp{AIC w -AIC min } ∑ n i exp{AIC i -AIC min } . Joint MAF spectra models were calibrated to absolute time using a per-generation mutation rate of 1.52 × 10 -8 based on a recent de novo mutation study with a M. murinus pedigree [START_REF] Campbell | Substitution rate analysis validates surprising results from a pedigree-based measurement of the de novo mutation rate in gray mouse lemurs[END_REF] and a generation time of 3.5 years based on parentage and survival estimations from wild populations of M. murinus [START_REF] Radespiel | An empirical estimate of the generation time of mouse lemurs[END_REF] and M. rufus [START_REF] Zohdy | Teeth, sex, and testosterone: Aging in the world's smallest primate[END_REF]. We consider both studies here by averaging the two (i.e., taking the midpoint) since generation time differences between species may reflect present-day resource availability rather than inherited life history traits.

Because joint modelling of all populations at once would lead to an untenable computational burden for model selection, populations were split into three sets of models that explored demographic histories of (1) Ambavala from the northern rainforests versus previously recognized M. mittermeieri (Anjanaharibe and Anjiahely), (2) Ambavala from the northern rainforests versus Ambatovy and Tsinjoarivo from the southern rainforests, and (3) Ambatovy and Tsinjoarivo from the southern rainforests versus Ankafobe and Ambohitantely in the Central Highlands. Complexity was added to models in a hierarchical fashion to constrain the total number of model comparisons (Supporting Information Methods "Demographic model selection"). A final set of models including all seven populations evaluated among subsets was run for parameter estimation. We used 200 nonparametric block-bootstraps (Supporting Information Methods "Block-bootstrapping MAF spectra") to estimate percentile-based 95% confidence intervals for the maximum pseudo-likelihood estimates.

| Characterizing population size changes through time

Population size change through time was evaluated with stairway plots [START_REF] Liu | Exploring population size changes using SNP frequency spectra[END_REF]. Similar to our other demographic analyses, we used genotype likelihood-based MAF spectra since there was no suitable outgroup to determine the ancestral state of each allele [START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF]. Because not all populations had appropriate sample sizes for stairway plot analyses, we focused on four representative populations: Ankafobe for the Central Highlands, Ambatovy for the southern rainforests, Ambavala for the northern rainforests, and a metapopulation of all formerly named M. mittermeieri individuals from the northern limit of the species distribution. For each population, we also generated 199 bootstrap replicates to obtain 95% confidence intervals, using the script provided by stairway plot version 2.0 [START_REF] Liu | Stairway Plot 2: demographic history inference with folded SNP frequency spectra[END_REF]. Analyses were rescaled from generations to absolute time given a per-generation mutation rate of 1.52 × 10 -8 and generation time of 3.5 years.

Because stairway plots are sensitive to the number of coalescent events, which is dependent on the number of individuals and SNPs in the data set [START_REF] Foote | Genome-culture coevolution promotes rapid divergence of killer whale ecotypes[END_REF][START_REF] Lapierre | Accuracy of demographic inferences from the site frequency spectrum: the case of the Yoruba population[END_REF][START_REF] Patton | Contemporary demographic reconstruction methods are robust to genome assembly quality: A case study in Tasmanian devils[END_REF]Teixeira, Salmona, et al., 2021), analyses were initially performed considering the entire population data sets (n ranging between 7 and 12) and repeated with an equal sample size (n = 7) to avoid differences in N e among populations related to sample size.

For the subsampled analyses, individuals were randomly selected (Table S1).

| RE SULTS

| RADseq data statistics

On average per individual, 8,181,656 raw reads were sequenced, of which 35.6% passed filtering and aligned to autosomal scaffolds (Table S2). Genotyping with STACKS and subsequent filtering resulted in a data set with 217,273 sites across 59 individuals (including 9 outgroup individuals) with an average 15× coverage per individual.

After the minor allele count filter, 111,104 sites were retained for downstream analyses. Genotype likelihoods estimated with ANGSD were available for 514,309 sites after filtering. Genotyping statistics can be found in Table S3 and percentage of missing data per individual is given in Table S4.

| Phylogenetic inference

All phylogenetic analyses recovered a well-supported split between northern rainforest and southern rainforest plus Central Highland M. lehilahytsara populations (Figure 2; Figures S2 andS3). There were minor incongruences between methods relating to the topology of the latter clade: MSC analyses (SNAPP and SVDquartets) found a well-supported split between Central Highland (Ambohitantely and Ankafobe) and southern rainforest (Ambatovy, Anjozorobe, and Tsinjoarivo) localities (Figure 2c; Figures S2 andS3) while the ML analysis (RaxML) found a well-supported clade comprising Anjozorobe and Central Highland populations (Figure 2a). ML analyses recovered sampling locations as clades, except for Ambatovy, which was paraphyletic with weak support, such that Tsinjoarivo was nested within (Figure 2a). In contrast, SVDquartets analyses of all individuals found Ambatovy and Tsinjoarivo to form two wellsupported clades (Figure S2). Furthermore, MSC-based analyses were not in complete agreement: SNAPP recovered Riamalandy and Ambavala as sister (Figure 2c) while SVDquartets recovered Riamalandy and Ambavala as a grade, albeit, with poor (Figure S2) or weak support (Figure S3).

| Delimiting population structure

Population structure mirrored findings from phylogenetic analyses and revealed admixture between M. lehilahytsara populations in the Central Highlands and southern rainforests, and between the southern and northern rainforest populations (Figure 2b). An optimal K of four was selected by STRUCTURE analyses of hard-filtered genotypes (Table S5), which yielded very similar results to the optimal K of three for genotype likelihoods from ngsAdmix (Figure S4; Table S6).

Notably, Ambohitantely had near 50:50 assignment to clusters from the Central Highlands and rainforest populations, represented by Ankafobe and Ambatovy, respectively (Figure 2b and Figure S5).

Populations previously identified as M. mittermeieri largely formed their own cluster, although some admixture between Anjiahely (south of the Antainambalana river) and Ambavala and between Anjiahely and Anjanaharibe (north of the Antainambalana river) was prevalent at higher K. Increasing K did not always recover the expected geographic groupings as separate clusters (Figures S5-S7).

PCA supported the STRUCTURE and ngsAdmix clustering results. PC1 generally separated northern rainforest, southern rainforest, and Central Highland populations. PC2 further separated southern rainforest populations from the northern rainforests and Central Highlands (Figure 3a). Ambohitantely took an intermediate position between Ankafobe and Ambatovy while Ambavala and Riamalandy fell between other northern rainforest populations and Ambatovy. PC3 further separated Tsinjoarivo from Ambatovy, and showed clustering among individual northern rainforest populations (Figure 3b). Analyses with genotype likelihoods found similar patterns (Figure S8). Significant geographic structure was detected between populations with AMOVA analyses (Table S7) although geographic sampling location explained only about 8% of genetic variation.

Pairwise F ST values (Table S8) regressed against Euclidean distances (Table S9) indicated the presence of isolation-by-distance when comparing all populations with a Mantel test (Pearson's r = .78, p = 7 × 10 -4 ; Table S10) (Figure 3c). This relationship was weaker when only considering F ST estimates between northern rainforests and southern/Central Highland populations (Pearson's r = .49, p = .15; Table S9). Several data points stood out as inconsistent with a simple IBD pattern. For example, genetic distance between Ankafobe and Ambohitantely was higher than between Ambatovy and Ambohitantely despite their phylogenetic placement as sisters and being only 12 km apart as opposed to 138 km. Similarly, the genetic distance between Ambatovy and Ambavala was much lower than expected from the regression model, which was congruent with evidence of admixture from STRUCTURE, ngsAdmix, and PCA results. Previously mentioned patterns were not always consistent with pairwise F ST estimates from genotype likelihoods (Figure S9), but results of Mantel tests (Table S9) similarly imply distance is not the sole explanation of observed genetic differentiation.

| Demographic modelling

Plausible models based on phylogenetic and population structure analyses were explored to understand historical and contemporary gene flow among northern populations (Ambavala, Anjiahely, and Anjanaharibe; Figure S10), between north and south (Ambavala vs. Ambatovy, and Tsinjoarivo; Figure S11), and between the Central Highlands and southern rainforests (Ankafobe and Ambohitantely vs. Ambatovy and Tsinjoarivo; Figures S12 andS13). The best models for each data set (Tables S11-S13) revealed that divergences between sister populations probably occurred sometime before human colonization of the Central Highlands and eastern Madagascar but after the LGM, even though confidence intervals were too wide to reject post-human arrival divergences (Figure 4). The models for northern populations (Figure 4a) and for the comparison of northern and southern populations (Figure 4b) detected a two-order of magnitude decline in N e , except for Ambatovy, which probably occurred after human colonization, but again, confidence intervals did not rule out early Holocene or late Pleistocene declines. No such crash was detected for the analysis of Central Highland and southern populations (Figure 4c). Rather, there was evidence for both historical and recent gene flow between the Central Highlands and Ambatovy. However, confidence intervals on migration rate parameters were typically wide, and there was no discernable pattern among migration rates through time to determine whether gene flow has actually decreased towards the present between the Central Highlands and Ambatovy (Figure 4c). Model selection results suggested that only historical gene flow was detectable between the north and south, which stopped shortly after the divergence of Ambatovy and Tsinjoarivo (Figure 4b), and that historical gene flow occurred between Ambavala and the common ancestor of Anjanaharibe and Anjiahely (Figure 4a).

Because confidence intervals were typically wide causing ambiguous interpretations and root divergence times were ultimately in conflict with the population tree (Figure 2c) split order (Tables S14-S16), we aimed to improve parameter estimation by jointly modelling all seven populations at once. This required some additional model testing to determine the population splitting order with gene flow (Figures S14 andS15; Table S17), and gene flow patterns observed in smaller models were assumed correct. Ancient population structure (i.e., the northern, southern, and Central Highland populations) formed at least before the LGM and probably prior to LIG (Figure 5). Subsequent population divergences occurred near the LGM, or at least the LGM was contained within confidence intervals while plausible human effects were excluded. The end of gene flow between northern and southern populations (i.e., Ambavala and Ambatovy) as well as an N e crash experienced by almost all populations occurred closer to the timing of human colonization, but potentially in the late Pleistocene after the LGM. Nearly a 10-fold N e decrease was estimated for all populations from the modelled recent N e change with the exception of Ambatovy (Figure 5; Table S18). There was no remarkable pattern among the migration rate parameters except that confidence intervals imply a high degree of uncertainty in migration rate point estimates. For a few parameters, the maximum pseudolikelihood estimates fall outside of their confidence intervals, but not for divergence times or more recent N e estimates.

| Characterizing population size changes through time

Similar demographic trends were recovered for the representative populations when using the same number of individuals per population (Figure 6). Assuming that undetected population structure is negligible, stairway plot inferred a period of constant N e that predated the LIG and lasted until the LGM, with N e lowest for the northern metapopulation formerly described as M. mittermeieri and highest for Ambatovy. After the LGM, there was a continuous decline in N e until present for all forest sites except for the former M. mittermeieri metapopulation, where declines were delayed until the mid-Holocene at approximately 7.5 ka. All populations show a second decline during or just before a hypothesized period of human colonization and anthropogenic landscape changes (c. 2 ka; Figure 6). Despite the large confidence intervals observed for the analyses considering seven individuals, our results support the presence of both pre-and post-human colonization population declines (Figure S16). The stairway plot analyses including more than seven individuals where possible confirmed sensitivity of N e estimates to sample size. Nevertheless, the same general patterns were observed compared to analyses with seven individuals (Figure S17), while recovering more fine-scale changes. Because of site filters for missing data, analyses with more than seven individuals ultimately had less data available for MAF spectra estimation than analyses subsampled to seven individuals (Table S19).

| Detecting potential migration corridors between the Central Highlands and southern rainforest populations

The Central Highlands of Madagascar are traversed by westward flowing rivers (Figure 7k), which were characterized by riparian forest or shrub vegetation before widespread transformation into agricultural fields [START_REF] Kull | Air photo evidence of historical land cover change in the highlands: Wetlands and grasslands give way to crops and woodlots[END_REF]. The Betsiboka river with its tributaries connects the humid rainforest population of M. lehilahytsara at Anjozorobe with the now isolated highland populations in Ambohitantely and Ankafobe. Remnants of riparian forest are still present in a few places, but most forests are converted into rice paddies today [START_REF] Kull | Air photo evidence of historical land cover change in the highlands: Wetlands and grasslands give way to crops and woodlots[END_REF]. The exchange of individuals using historically forested corridors along the Betsiboka river would also explain observed gene flow between southern rainforest populations such as Ambatovy and the Central Highlands.

However, the southern populations of M. lehilahytsara are separated by one large and several smaller rivers (Figure 7k). While it is likely that populations at Ambatovy and Anjozorobe (c. 60 km apart) could have been connected across the headwaters of the Mangoro river (Figure 7j

m b 0 2 m b 2 0 migration rate (N m /N e ) m a 1 2 m a 2 1 m b 0 2 m b 2 0 m a 1 2 m a 2 1 m a 1 3 m a 3 1 m a 2 3 m a 3 2 m b 1 2 m b 2 1 m b 1 3 m b 3 1 m b 2 3 m b 3 2 m c 1 3 m c 3 1 (a)
A n ja n a h a r ib e A n ji a h e ly A m b a v a la 5). Charcoal records (Burney, 1987a[START_REF] Burney | Late Quaternary stratigraphic charcoal records from Madagascar[END_REF] and stratigraphic pollen occurrence data [START_REF] Burns | Rapid humaninduced landscape transformation in Madagascar at the end of the first millennium of the common era[END_REF][START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF][START_REF] Virah-Sawmy | Evidence for drought and forest declines during the recent megafaunal extinctions in Madagascar[END_REF] have shown a period of increased drought and fire from about 6 ka to the present.

The LGM also induced cool and dry conditions in Madagascar [START_REF] Clark | The last glacial maximum[END_REF][START_REF] Scholz | Scientific drilling in the great rift valley: The 2005 Lake Malawi scientific drilling project -an overview of the past 145,000 years of climate variability in southern hemisphere East Africa[END_REF][START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF], which were characterized by more open habitats in the Central Highlands [START_REF] Burney | A chronology for late prehistoric Madagascar[END_REF][START_REF] Samonds | A new late Pleistocene subfossil site (Tsaramody, Sambaina basin, central Madagascar) with implications for the chronology of habitat and megafaunal community change on Madagascar's Central Highlands[END_REF], even though it is debatable what the vegetational components of that open habitat were (e.g., [START_REF] Joseph | How expansive were Malagasy Central Highland forests, ericoids, woodlands and grasslands? A multidisciplinary approach to a conservation conundrum[END_REF]. During this dry period, M. lehilahytsara, which is currently associated with humid environments, likely retreated to watersheds buffered from aridification [START_REF] Muldoon | Early Holocene fauna from a new subfossil site: A first assessment from Christmas River, south central Madagascar[END_REF][START_REF] Wilmé | Biogeographic evolution of Madagascar's microendemic biota[END_REF]. This would explain the near simultaneous population separations in the north, south, and Central Highland lineages. Although the ML estimate and confidence interval show uncertainty in the divergence time between Ankafobe and Ambohitantely, this may reflect the complex patterns of postdivergence gene flow between the Central Highlands and southern rainforest populations.

Rivers are known to function as barriers to gene flow in mouse lemurs [START_REF] Olivieri | The ever-increasing diversity in mouse lemurs: three new species in north and northwestern Madagascar[END_REF] and other lemurs [START_REF] Quéméré | Landscape genetics of an endangered lemur (Propithecus tattersalli) within its entire fragmented range[END_REF], and our data suggest that they are more significant in et al., 2016), but this would presumably not be possible after humanmediated fire-driven expansion of grasslands [START_REF] Burney | A chronology for late prehistoric Madagascar[END_REF][START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF]. Our gene flow model is supported by the identification of riparian forest corridors (Figure 7) that would have enabled movements between the Central Highlands and the eastern escarpment rainforests. Such corridors may also have been key to previous mouse lemur dispersal events across a grassland and S3) while filling a vital sampling gap along the proposed river corridors (Figure 7).

| Making conservation biology predictions from models of molecular evolution

Here, we provided a range of present-day N e estimates that are sensitive to model choice, though the direction of change in N e through time is not. Every analysis shows that M. lehilahytsara populations sizes have declined. Our SFS model selection and parameter estimation have been rigorous, and population substructure was accounted for in stairway plot analyses, since historical admixture or structure will manifest as N e changes in single-population methods (e.g., [START_REF] Chikhi | The confounding effects of population structure, genetic diversity and the sampling scheme on the detection and quantification of population size change[END_REF][START_REF] Mazet | On the importance of being structured: instantaneous coalescence rates and human evolution -lessons for ancestral population size inference[END_REF]. All evidence points to drastic population declines following the LGM as well as more recent declines coincident with human activities. These findings are consistent with previous studies of an endemic olive tree [START_REF] Salmona | Late Miocene origin and recent population collapse of the Malagasy savanna olive tree[END_REF], a mouse lemur [START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF], and two sifaka species [START_REF] Quéméré | Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region[END_REF][START_REF] Salmona | Climate change and human colonization triggered habitat loss and fragmentation in Madagascar[END_REF], where population structure, as a proxy for forest fragmentation, and population declines were attributed to both prehuman climatic and later human effects in northern Madagascar. The most recent IUCN evaluation of M. lehilahytsara considers the species to be vulnerable [START_REF] Dolch | Microcebus lehilahytsara. The IUCN Red List of Threatened Species 2020[END_REF] while other mouse lemur species in the eastern rainforest are considered endangered or critically endangered (Andriaholinirina et al., 2014a[START_REF] Andriaholinirina | Microcebus marohita. The IUCN Red List of Threatened Species[END_REF][START_REF] Borgerson | Microcebus jollyae (amended version of 2020 assessment). The IUCN Red List of Threatened Species 2020[END_REF]. We surmise that though M. lehilahytsara appears to be more robust to climate-driven and human-mediated disturbances than other species, this is more likely a consequence of its large ancestral geographic range combined with its physiological ability to inhabit both the lowlands and higher altitudes. It was previously regarded as a highland specialist [START_REF] Radespiel | First indications of a highland specialist among mouse lemurs (Microcebus spp.) and evidence for a new mouse lemur species from eastern Madagascar[END_REF] but is now also recognized in lowland areas in the northeast after synonymizing M. mittermeieri (Schüßler, Mantilla-Contreras, et al., 2020). M. lehilahytsara has been documented to enter a state of prolonged torpor [START_REF] Blanco | Evidence of prolonged torpor in Goodman's mouse lemurs at Ankafobe forest, central Madagascar. Primates[END_REF] which may either be an adaptive trait that has facilitated its increased geographic range compared to most other mouse lemur species, or a plastic one that is only evident at higher altitudes. The importance of torpor for sustaining population sizes despite limited resources remains unclear.

Regardless, the downward trends in N e across the Central Highlands and the studied rainforest populations imply a negative consequence of continued habitat loss and fragmentation.

| The importance of interpreting biodiversity through evolution

The processes that shaped current M. lehilahytsara population structure must have occurred before human arrival and mid-Holocene droughts. Although charcoal records (Burney, 1987a[START_REF] Burney | Late Quaternary stratigraphic charcoal records from Madagascar[END_REF], stratigraphic pollen occurrence [START_REF] Burns | Rapid humaninduced landscape transformation in Madagascar at the end of the first millennium of the common era[END_REF][START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF], and megafauna subfossil extinctions [START_REF] Hansford | Early Holocene human presence in Madagascar evidenced by exploitation of avian megafauna[END_REF] all point to human-mediated large-scale ecological changes in the Central Highlands [START_REF] Virah-Sawmy | Evidence for drought and forest declines during the recent megafaunal extinctions in Madagascar[END_REF], our analyses suggest that M. lehilahytsara population divergences coincided with Pleistocene climate fluctuations both before and near the

LGM. Our results, especially when combined with the narrative from other analyses of lemur populations [START_REF] Quéméré | Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region[END_REF][START_REF] Salmona | Climate change and human colonization triggered habitat loss and fragmentation in Madagascar[END_REF][START_REF] Teixeira | Past environmental changes affected lemur population dynamics prior to human impact in Madagascar[END_REF] and endemic olive trees [START_REF] Salmona | Late Miocene origin and recent population collapse of the Malagasy savanna olive tree[END_REF], reveal a timeline for vegetational change in Madagascar from a population genetic perspective. Analyses of multiple species and localities agree with ecological [START_REF] Bond | The antiquity of Madagascar's grasslands and the rise of C 4 grassy biomes[END_REF] and evolutionary [START_REF] Hackel | Grass diversification in Madagascar: In situ radiation of two large C 3 shade clades and support for a Miocene to Pliocene origin of C 4 grassy biomes[END_REF][START_REF] Vorontsova | Madagascar's grasses and grasslands: anthropogenic or natural[END_REF] studies of grasses that demonstrate Madagascar's Central Highland grasslands are of ancient origin, and further, that the Central Highlands have long persisted as a mix of grassland and woodland.

Although fires due to human activity are prominent across the Central Highlands (Burney, 1987a;[START_REF] Burns | Rapid humaninduced landscape transformation in Madagascar at the end of the first millennium of the common era[END_REF][START_REF] Gasse | A 40,000-yr pollen and diatom record from Lake Tritrivakely, Madagascar, in the southern tropics[END_REF][START_REF] Samonds | A new late Pleistocene subfossil site (Tsaramody, Sambaina basin, central Madagascar) with implications for the chronology of habitat and megafaunal community change on Madagascar's Central Highlands[END_REF][START_REF] Virah-Sawmy | Evidence for drought and forest declines during the recent megafaunal extinctions in Madagascar[END_REF], leading to devastating forest loss in Madagascar [START_REF] Vieilledent | Combining global tree cover loss data with historical national forest cover maps to look at six decades of deforestation and forest fragmentation in Madagascar[END_REF],

the

  Memorial Foundation; Laboratoire d'Excellence TULIP, Grant/Award Number: ANR-10-LABX-41 and ANR-11-IDEX-0002-02; North-German Supercomputing Alliance; Primate Action Fund Handling Editor: Kimberly Gilbert structure, demographic models of post-divergence gene flow, and population size change through time. The results offer evidence that open habitats are an ancient component of the Central Highlands, and that widespread forest fragmentation occurred naturally during a period of decreased precipitation near the last glacial maximum. Models of gene flow suggest that migration across the Central Highlands has been possible from the Pleistocene through the recent Holocene via riparian corridors.

  Here, we sampled M. lehilahytsara populations from multiple sites in the Central Highlands as well as from eastern rainforests, which include formerly recognized M. mittermeieri populations in the northeast and sites that are further south but along the escarpment where rainforests transition to the Central Highlands and more open vegetation. Using phylogeographic approaches and coalescent methods to analyse a large restriction site associated DNA sequencing (RADseq) data set, we explored patterns of genetic connectivity among populations. We also estimated divergence times between populations and changes in effective population size (N e ) through time. Our findings show that population divergences occurred well-before humans arrived in Madagascar. The initial isolation of the Central Highlands from eastern rainforests occurred even before the Last Interglacial period (LIG; c. 132-112 ka). Post-divergence gene flow among populations in the Central Highlands and eastern rainforests suggests that migration was possible despite hypothesized large-scale vegetative turnover from more closed to open ecosystems. We find that there have been two notable episodes of population size contraction; the first coinciding with the last glacial kaximum (LGM; c. 26.5-19 ka) and the second coinciding with recent human colonization (c. 2 ka). These patterns suggest that a combination of paleoclimatic and anthropogenic effects have synergistically increased extinction risk in M. lehilahytsara -a phenomenon that undoubtedly applies to other forest-dwelling species.

F

  I G U R E 1 Sampling sites for Microcebus lehilahytsara populations. Forest cover in green is taken from Schüßler, Mantilla-Contreras, et al. (2020), which reflects modern-day presence rather than historical connectivity. Water ways are shown in blue with relevant rivers named in italics. The x-axis is longitude, and the y-axis is latitude. [Colour figure can be viewed at wileyonlinelibrary. com] with the Sbfl restriction enzyme. A median insert size of 380 bp was selected for paired-end sequencing at the Duke Sequencing Core on an Illumina Hi-Seq 4000. Two libraries were prepared with the single-SbfI protocol described in Etter et al. (2011) and Genomic Resources Development Consortium et al. (2015). These two libraries were singleend sequenced on an Illumina HiSeq 2000 at the University of Oregon Core facility. Information on library preparation and sequencing for all new and previously published samples are given in Table

2. 8 |

 8 Identifying potential migration corridors along rivers We used ArcMap of ArcGis Desktop version 10.6.1 (ESRI Inc.) to visualize potential migration corridors along rivers by creating 10 m resolution false-colour composites (infrared, red green, bands 8-4-3) from tiles T38KQE, T38KQF, T38KRD and T38KRE of Sentinel-2 satellite images (Multi-Spectral Instrument, Bottomof-Atmosphere reflectance [S2MSI-2A]) taken on 23 August 2021 and acquired from the European Space Agency sci-hub (https:// scihub.coper nicus.eu/).

  Because the split order between populations is of consequence for downstream demographic modelling, we also explored the strength of evidence in favour of a north-versus-south/highlands hypothesis relative to a highlands-versus-humid rainforests hypothesis. ML-based AU tests rejected the alternative scenario of a highlands-versus-humid rainforests basal split (p-value = .0179).Similarly, marginal likelihoods from SNAPP decisively favoured the north-versus-south/highlands model (marginal ML = -82860.30) over the highlands-versus-humid rainforests model (marginal ML = -82889.96) with a Bayes factor of -59.32.

F

  I G U R E 2 Well-supported phylogenetic lineages of Microcebus lehilahytsara reflect genetic structure. (a) The ML phylogeny of all individuals indicates well-supported nodes (bootstrap support ≥70%) with an asterisk (*), while poorly supported nodes (bootstrap support <70%) are blank. Branch colours correspond to population colours in Figure 1. Branch lengths are expected substitutions per site. (b) A STRUCTURE analysis of all individuals reflects supported clades that correspond to geography at an optimal K of four. (c) A full likelihood MSC analysis with SNAPP further supports the ML topology and major splits with strong posterior probabilities, with the exception of Anjozorobe, which is represented by a single individual. Branch lengths are expected substitutions per site. [Colour figure can be viewed at wileyonlinelibrary.com] Principal component analysis (PCA) of genotype data and pairwise F ST . (a) The first and second axes generally separate northern, southern, and Central Highland populations. (b) The third axis further separates northern populations and Tsinjoarivo from Ambatovy. (c) Pairwise F ST values between populations reveal evidence for IBD, but several points deviate from the trend that may be explained by gene flow. The dashed line represents the fitted linear regression. [Colour figure can be viewed at wileyonlinelibrary.com]

  ), Anjozorobe and Tsinjoarivo are already separated by c. 140 km of rainforest traversed by several smaller eastward flowing rivers with their sources in the Central Highlands. Finally, Ambatovy and Tsinjoarivo are separated by the very wide Mangoro river which is, apart from its source region, probably too large to be crossed by mouse lemurs (Figure 7h-j).

F

  Demographic models reveal the history of divergence, gene flow, and population size among studied populations. Selected SFS models for evaluating population divergence, population size change, and gene flow (a) among northern, (b) between northern and southern, and (c) between the Central Highland and southern rainforest populations. All parameter estimates are shown on a log 10 scale to reflect changes in orders of magnitude. Points are maximum pseudo-likelihood estimates and error bars are percentile-based 95% confidence intervals. Colours for divergence times are the LGM (blue) and the time since human colonization (yellow). [Colour figure can be viewed at wileyonlinelibrary.com]

  in jo a r iv o x x x θ A m b a to v y θ A m b a v a la θ T s in jo a r iv o x x x θ A m b a to v y θ A m b a v a la θ r θ s θ t θ A n k a fo b e θ A m b o h it a n te ly θ A m b a to v y θ T s in jo a

|

  Population structure of Microcebus lehilahytsara predates human colonization of MadagascarOur results indicated that M. lehilahytsara Central Highland populations separated from rainforest populations long before human colonization of these regions (c. 2 ka; e.g.,[START_REF] Douglass | A critical review of radiocarbon dates clarifies the human settlement of Madagascar[END_REF] and predate even the oldest potential human presence in Madagascar (c. 10 ka;[START_REF] Hansford | Early Holocene human presence in Madagascar evidenced by exploitation of avian megafauna[END_REF]. Our demographic models indicated that this divergence and the distinct genetic clusters (Figure2b; Figure3a) could have originated as far back as 500 ka, which may be associated with older Pleistocene interglacial cycles and the fragmentation of a previously more continuous habitat (Figure5). Divergences among populations within the northern, southern, and Central Highland clusters likely originated sometime between the LGM and mid-Holocene (Figures 4 and

F I G U R E 5

 5 Parameter estimates for the seven-opulation model inferred with fastsimcoal2. Each node is labelled r through w, which has a divergence time and N e parameter associated with it. N e parameters are indicated by to avoid the extra subscript. Points are maximum pseudo-likelihood estimates and error bars are 95% confidence intervals. For the divergence time plot, the yellow colour indicates posthuman colonization of Madagascar, and the blue colour indicates the timing of the LGM. All y-axes were transformed to a log 10 scale to better visualize order-of-magnitude changes. [Colour figure can be viewed at wileyonlinelibrary.com]

  r θ s θ t θ u θ v θ w θ A n ja n a h a r ib e θ A n ji a h e ly θ A m b a v a la θ A n k a fo b e θ A m b o h it a n te ly θ A m b a to v y θ T s in jo a r iv o θ A n ja n a h a r ib e θ A n ji a h e ly θ A m b a v a la θ A n k a fo b e θ A m b o h it a n te ly θ A m b a to v y θ T s in jo a r iv o i a h e l y A m b a v a l a A n k a f o b e A m b o h i t a n t e l y A m b a t o v y T s i n j o a r i v o

  Figure S8). Similarly prominent was the divergence of Anjiahely from Anjanaharibe and Marojejy in the north. Anjiahely lies south of the Antainambalana river system while Anjanaharibe and Marojejy are to the north (Figure1). There was no evidence of post-divergence gene flow between Anjiahely and Anjanaharibe (Figure4), and although this region is characterized by heavy deforestation since the 1990s(Schüßler, Mantilla-Contreras, et al., 2020), their divergence time was

F I G U R E 7

 7 Potential riparian corridors along the Betsiboka and Mangoro Rivers. Riparian corridors potentially link rainforest and Central Highland populations of Microcebus lehilahytsara (k). The inlays (a-j) illustrate the presence of riparian forests (rf), riparian areas that were converted into rice paddies (rp) but were forested before[START_REF] Kull | Air photo evidence of historical land cover change in the highlands: Wetlands and grasslands give way to crops and woodlots[END_REF], transitional forests (tf) between rainforest and riparian highland vegetation and degraded vegetation (dv) along the two rivers. Inlays are produced from Sentinel-2 false colour composites. Rivers were derived fromSchüßler, Mantilla-Contreras, et al. (2020) and forest cover was taken from[START_REF] Vieilledent | Combining global tree cover loss data with historical national forest cover maps to look at six decades of deforestation and forest fragmentation in Madagascar[END_REF]. The images show in blueturquoise colour the river courses, in light green-beige grassland vegetation and bare soil, in pale-slightly red colour with regular shapes rice paddies, in light red shrubland, thin riparian vegetation narrower than the resolution of 10 m, and in rich red transitional forest vegetation or riparian forests. [Colour figure can be viewed at wileyonlinelibrary.com] woodland mosaic. It is unlikely that substantial migration between the studied Central Highland and southern rainforest populations persists today as the rapid deforestation of Madagascar's eastern rainforests in recent years has resulted and will further result in the fragmentation of forest corridors(Schüßler, Mantilla-Contreras, et al., 2020;[START_REF] Vieilledent | Combining global tree cover loss data with historical national forest cover maps to look at six decades of deforestation and forest fragmentation in Madagascar[END_REF]. Some caution is needed in the interpretation of our demographic models, as there is little discernable direction in gene flow, and confidence intervals for migration rate parameters typically span several orders of magnitude. This shows some limitations of our analyses, as we are obtaining estimates under a parameter-rich model, where there can be identifiability problems given the many combinations of divergence times, N e , and migration rates that can explain the observed MAF spectra. However, several other observations are consistent with our model such as the higher-than-expected genetic distance between Ankafobe and Ambohitantely (Figure3; Figure S9), model-based clustering analyses (Figure 2; Figures S5 and S6), and the intermediate position of Ambohitantely between Ankafobe and Ambatovy in PCAs (Figure 3; Figure S8). Additional sampling of individuals may be necessary to potentially reduce confidence interval widths. In addition, having population-level representation of Anjozorobe would resolve the site's current phylogenetic ambiguity (Figure 2; Figures S2 and

  effects of ongoing habitat fragmentation are occurring on a background of habitat mosaicism that is probably associated with Pleistocene climate change. AUTH O R CO NTR I B UTI O N S George P. Tiley and Anne D. Yoder conceived the study. Marina B. Blanco, Edward E. Louis, José M. Ralison, Blanchard Randrianambinina, and Rodin M. Rasoloarison performed the fieldwork. Amanda R. Stahlke and Paul A. Hohenlohe constructed RADseq libraries. George P. Tiley, Tobias van Elst, Helena Texeira, and Dominik Schüßler performed analyses. George P. Tiley, Tobias van Elst, and Anne D. Yoder wrote the manuscript. Jordi Salmona, Lounès Chikhi, and Ute Radespiel contributed to the intellectual and

Table S12

 S12 Population size change of Microcebus lehilahytsara through time estimated with stairway plots. The same sample size (n = 7) was used for all populations. The metapopulation of Anjiahely, Anjanaharibe and Marojejy northern rainforest sites is labelled M. mittermeieri to ease legend formatting although these individuals are now recognized as M. lehilahytsara. Absolute time calibrations were a generation time of 3.5 years and a pergeneration mutation rate of 1.52 × 10 -8 . The vertical bars identify the LGM between 19 and 26.5 ka, and the start of a significant human impact in Madagascar between 1 and 2 ka. [Colour figure can be viewed at wileyonlinelibrary.com]

) suggests movement across the Central Highland habitat through time, even compared to models that allowed for some post-divergence gene flow

F I G U R E 6

but not until the present. Previous research has also suggested that mouse lemurs were capable of moving across the Central Highlands during the Pleistocene based on phylogeographic inference

(Yoder

* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * -9 e -0 4 -6

Funding information

Alexander von Humboldt-Stiftung; Deutsche Forschungsgemeinschaft, Grant/Award Number: Ra 502/23-1; Duke Tropical Conservation Initiative; Fundação para a Ciência e Tecnologia, Grant/Award Number: PTDC-BIA-EVL/30815/2017; International Research Project-Bioinformatics, Ecology, Evolution, Genomics and Behaviour; Investissement d'Avenir, Grant/Award Number: ANR-10-LABX-25-01; John Simon Guggenheim of grasslands, interspersed with patches of forest. The historical perspective was that Madagascar's grasslands

DATA AVA I L A B I L I T Y S TAT E M E N T

All new sequencing data have been made available through NCBI BioProject PRJNA560399. Individual BioSample accessions are given in Table S1. VCFs, alignment files, phylogenies, and control files for demographic analyses will be available from the Dryad Digital Repository (https://doi.org/10.5061/dryad.dncjs xkvz).

Additional scripts for data processing and analysis are available from TvE's GitHub page (https://github.com/t-vane/Resea rchSu pplem ents).

creative direction of the study. All authors read, revised, and approved the final manuscript.

ACK N OWLED G EM ENTS

This is a Duke Lemur Centre publication #1510. We thank two anonymous reviewers and the editor for helpful comments that improved the manuscript.

CO N FLI C T O F I NTE R E S T

The authors declare no competing interests.

S U PP O RTI N G I N FO R M ATI O N

Additional supporting information can be found online in the Supporting Information section at the end of this article.