Genetic population structure of endangered ring‐tailed lemurs (Lemur catta) from nine sites in southern Madagascar

Abstract Madagascar's ring‐tailed lemurs (Lemur catta) are experiencing rapid population declines due to ongoing habitat loss and fragmentation, as well as increasing exploitation for bushmeat and the illegal pet trade. Despite being the focus of extensive and ongoing behavioral studies, there is comparatively little known about the genetic population structuring of the species. Here, we present the most comprehensive population genetic analysis of ring‐tailed lemurs to date from across their likely remaining geographic range. We assessed levels of genetic diversity and population genetic structure using multilocus genotypes for 106 adult individuals from nine geographically representative localities. Population structure and F ST analyses revealed moderate genetic differentiation with localities being geographically partitioned into northern, southern, western and also potentially central clusters. Overall genetic diversity, in terms of allelic richness and observed heterozygosity, was high in the species (AR = 4.74, H O = 0.811). In fact, it is the highest among all published lemur estimates to date. While these results are encouraging, ring‐tailed lemurs are currently affected by ongoing habitat fragmentation and occur at lower densities in poorer quality habitats. The effects of continued isolation and fragmentation, coupled with climate‐driven environmental instability, will therefore likely impede the long‐term viability of the species.


| INTRODUC TI ON
Sustaining natural levels of genetic diversity within wildlife populations is a key concern for conservation biologists (Frankham, 1995(Frankham, , 2003(Frankham, , 2005. Pressures from climate change, anthropogenic habitat modification, overexploitation, and the introduction of novel competitors and infectious diseases are producing rapidly and ever-changing environments, forcing species to adapt and evolve or go extinct (Di Marco, Venter, Possingham, & Watson, 2018;Frankham, Ballou, & Briscoe, 2010). Genetic diversity, the variation of alleles and genotypes present within a population, is the foundation on which natural selection acts and is therefore necessary for adaptive evolutionary change to occur (Frankham, 1995(Frankham, , 2003(Frankham, , 2005Frankham et al., 2010). Populations with low levels of genetic diversity struggle to evolve in modified environments. For instance, the Tasmanian devil (Sarcophilus harrisii) has become vulnerable to the spread of devil facial tumor disease due to the lack of diversity across immune genes following human-induced population crashes from introduced diseases (Guiler, 1970;Morris, Wright, Grueber, Hogg, & Belov, 2015).
Habitats modified by human activity hold less genetic diversity than pristine environments, thus putting their inhabitants at high risk (Miraldo et al., 2016). That is because deforestation, fragmentation, and habitat degradation-all critical threats to biodiversityinteract to restrict the amount of viable habitat available to species, reduce carrying capacity and consequently maximum population size, and create isolated patches separated by matrix (i.e., inhospitable habitat) that impedes gene flow among remaining species' populations (Baden et al., 2019;Holmes et al., 2013;Stangel, Lennartz, & Smith, 1992). Combined, these processes lead to greater inbreeding, reduced genetic diversity, and ultimately an increased extinction risk (Frankham et al., 2010;Lino, Fonseca, Rojas, Fischer, & Pereira, 2019;Struebig et al., 2011).
Though not alone in its vulnerability to habitat loss, Madagascar's biodiversity is considered to be a top concern, in part because of its incredible concentration of species endemism (Myers, Mittermeier, Mittermeier, da Fonseca, & Kent, 2000). Since its colonization by humans as recently as 4,000 years ago, the island has undergone extensive forest cover loss and with it more than 17 species of large-bodied lemurs (Dewar et al., 2013;Godfrey & Irwin, 2007;Kistler et al., 2015;Myers et al., 2000). Unfortunately, deforestation in Madagascar continues unabated (Vieilledent et al., 2018) and scientists anticipate that remaining rainforest habitat will be lost before the end of this century (Morelli et al., 2020). When considered alongside the impacts of climate change, this threat poses significant risk to the persistence of remaining extant lemur species (Brown & Yoder, 2015;Morelli et al., 2020). It is therefore an urgent conservation priority to quantify the genetic variability present within Madagascar's only endemic primate radiation to assess whether and to what extent lemur species can cope with intensifying environmental pressures.
Despite being one of the most-studied lemur species, there is relatively little known about the genetic diversity and population structure of remaining wild ring-tailed lemur individuals.
Moreover, low F ST values among sites in the north (Clarke et al., 2015) and among those in the west (Parga et al., 2015) indicate minimal genetic differentiation, suggesting the presence of historical gene flow. While encouraging, there is also evidence that genetic erosion within the species has already begun to negatively impact their health and fitness (e.g., Charpentier, Williams, & Drea, 2008;Grogan, Sauther, Cuozzo, & Drea, 2017). It is therefore likely that at least some of Madagascar's remaining ring-tailed lemur populations are already experiencing a time-delayed response (i.e., extinction debt), as extinctions do not typically occur until several generations after a fragmentation event (Jackson & Sax, 2010;Tilman, May, Lehman, & Nowak, 1994).
In an effort to characterize the remaining genetic diversity present within the species and identify how this diversity is apportioned among remnant populations, we provide a preliminary population genetic assessment of ring-tailed lemurs across their remaining geographic range. We evaluate within-and across-site levels of genetic diversity and infer population genetic structure to better understand this species' adaptive potential and diagnose possible conservation priorities.

| Sample collection
Our sample included 106 adult ring-tailed lemurs from nine geographically representative localities from across their existing range (Table 1, Figure 2). This dataset includes previously published genetic data from 75 adult ring-tailed lemurs from five localities (Clarke et al., 2015;Parga et al., 2012), as well as 31 new individuals from an additional four sites (Table 1). Published genetic data were collected in May through August 2006 (Parga et al., 2012) and August to October 2012 (Clarke et al., 2015). Data for new individuals were generated from fecal samples collected from Isalo in July 2016 and from Ambirary (AMB), Beoloke (BLK), and Berenty (BER) in June and July 2017. Multiple individuals and groups were sampled at each locality (Table 1). Fecal samples were immediately stored in RNAlater (Ambion) to prevent DNA degradation and were banked within 1 month of collection at −80°C for long-term storage. Sample collection and export/import protocols adhered to Malagasy and International laws and were approved by Malagasy wildlife authorities and the US Fish and Wildlife Service.

F I G U R E 1 Photograph of ring-tailed lemur (Lemur catta) by M LaFleur
TA B L E 1 Sampling localities, geographic coordinates, and sample sizes of ring-tailed lemur populations used in this study

| Microsatellite genotyping
Extraction products were amplified via PCR in a 13 μl reaction volume using 6.25 μl HotStarTaq DNA polymerase Master Mix, 20 mg/ ml BSA, 10 μM primer pairs, and 3 μl (0.25-1 ng) gDNA using annealing temperatures outlined in the Appendix S1. The 5′ end of the forward primer was fluorescently labeled. PCR products were separated by capillary electrophoresis (ABI 3730xl Genetic Analyzer), and alleles were sized to an internal size standard (Rox-500) using GeneMarker software v.2.6.7 (SoftGenetics). Genotype assignment was based on multiple independent reactions, where heterozygotes were confirmed with at least two independent reactions and homozygotes were confirmed with five independent reactions (Morin, Chambers, Boesch, &amp;Vigilant, 2001;Taberlet et al., 1996). Individuals from earlier studies (n = 45, Parga et al., 2012;n = 30, Clarke et al., 2015) were regenotyped and scaled to ensure datasets were comparable. CERVUS v.3.0 (Kalinowski, Taper, & Marshall, 2007) was used to calculate probability of identity (P ID ), that is, the probability that two randomly drawn individuals from a population will show identical multilocus genotypes. using GenoDive (Meirmans & Van Tienderen, 2004). We standardized allelic richness (AR) to the smallest sample size in the dataset to account for uneven sampling between populations using HP-Rare 1.1 (Kalinowski, 2005).

| Population genetic structure
To assess the genetic distances between sampling localities, Wright's F ST (Weir & Cockerham, 1984) was calculated for all pairs of populations using GenoDive (Meirmans & Van Tienderen, 2004). F ST is a measure of genetic differentiation among subpopulations and illustrates whether and to what extent populations are considered genetically distinct (Frankham et al., 2010). Significance was calculated using 10,000 permutations corrected for multiple comparisons (Bonferroni adjusted p = 0.001).
The presence of isolation-by-distance (IBD) was evaluated using the program GenAlEx v.6.5 (Peakall & Smouse, 2012) and significance estimated with Mantel's test using 10,000 permutations  Table 1. Gray shading illustrates historic ring-tailed lemur distribution across southern Madagascar retrieved from IUCN website (Mantel, 1967). Genetic distances between populations were esti- To identify genetic clusters, we used three different methods.
This method is used to estimate the number of genetically distinct clusters (K) with no a priori information regarding the individuals' geographic sampling locations provided, so clusters were formed solely on genetic information. We evaluated the hypothesis K = 1-12, three more than the number of wild populations (Evanno, Regnaut, & Goudet, 2005), using 100,000 iterations of MCMC following a burn-in of 50,000 iterations, as longer burn-in or MCMC did not significantly change our results. We implemented 20 runs for each value of K and assumed admixture and correlated allele frequencies.
The admixture model allows individuals to have mixed ancestry, assuming that a portion of an individual's genome, q, comes from a subpopulation, k (where ∑ k q k = 1; Falush, Stephens, & Pritchard, 2003). To account for unbalanced sampling, the ALPHA value was changed from the default value (1.0) to 0.5 (Wang, 2016). The most likely number of genetic clusters (K) was assessed using the highest value of ∆K (Evanno et al., 2005) using the program Structure Harvester v0.6.94 (Earl & vonHoldt, 2012). Structure Harvester calculates the second-order rate of change of the likelihood distribution (ΔK), which indicates the most pronounced subdivision within the data and the optimal number of genetic clusters. We implemented a twostep approach to evaluate further substructure in the dataset. We first identified the most likely number of clusters within the overall sample (n = 106) and then ran subsequent analyses within each of the K clusters (Evanno et al., 2005).
To corroborate the Structure analysis, a discriminant analysis of principal components (DAPC) was performed in R, using the adegenet package (Jombart, 2008;R Core Team, 2017). This multivariate method identifies clusters of genetically related individuals that maximize between-group variability and minimize within-group variability by using a set of retained principal components (determined by the user to optimize variance explained; Jombart, Devillard, & Balloux, 2010). The optimal number of clusters is determined by the number of clusters with the smallest Bayesian information criterion (BIC) value.

| Genetic diversity
The nonexclusion combined probability of identity (P ID ; Paetkau & Strobeck, 1994) of the six markers used in this study was 5.22 × 10 −10 and for P ID-sibs was 1.13 × 10 −3 . These values demonstrate a low probability that two individuals would share the same multilocus genotype. The six loci were highly polymorphic, with the number of alleles ranging between 12 and 15 alleles ( Table 2). Individuals were pooled across sampling localities, and there was no evidence of linkage disequilibrium. One locus was found to significantly deviate from Hardy-Weinberg equilibrium (L-2; Table 2); however, it did not deviate at any one site specifically and was therefore included in further analyses.
As a species, ring-tailed lemurs maintain high levels of genetic diversity despite severe habitat fragmentation across their range. Mean number of alleles (MNA) ranged from 4.33 to 8.67 (Table 3). The mean observed heterozygosity across sampling sites was 0.811 ± 0.044, while mean expected heterozygosity was 0.775 ± 0.054. Overall F IS was −0.052, and values ranged from −0.194 at ISALO to 0.042 at BER (Table 3).    to Cluster 2, and the central locality (ISALO) to Cluster 3. Individuals in the western localities (BEZA and TNP) showed higher membership probabilities to Cluster 4. From our initial analyses, we could not confidently group ISALO, as our Structure results indicated that ISALO clustered with southern localities (BER-AMB-BLK) despite the distance (301 km) and our DAPC results grouped the site into its own central cluster. Upon further investigation, we did find support in our Structure analysis for higher K-values (K = 3-4, Figure S2), supporting the DAPC results and grouping ISALO into a separate central cluster.

11
This study.
15 Salmona et al. (2015).  (Jolly, 1966;Sussman, 1972Sussman, , 1974, ring-tailed lemurs can disperse more easily across nonforested areas than forest-dependent arboreal species. Furthermore, the largest river drainage systems in southern Madagascar are seasonal and therefore do not pose permanent dispersal barriers to this species; in fact, they may actually be used as dispersal corridors (Goodman et al., 2006).
There was strong evidence for isolation-by-distance (IBD), meaning there was a positive correlation between genetic and geographic distances among populations. There are currently no records of subfossil ring-tailed lemurs outside of their current distribution (Godfrey, Jungers, Simons, Chatrath, & Rakotosamimanana, 1999), suggesting that this broad geographical range has been stable through geological time (Goodman et al., 2006). If localities had been isolated for a long time, we would expect genetic drift to erase any pattern of IBD (Baden et al., 2014;Olivieri et al., 2008). Therefore, there is potentially movement and relatively recent interconnectivity between ring-tailed populations via river basins (e.g., Mandrare River; Goodman et al., 2006). Genetic structure can result from limited gene flow or from historical events such as fragmentation; however, distinguishing between these processes can be challenging, especially when demography is unknown and forest fragmentation is recent.
Measures of genetic diversity, gene flow, and population structure are subjected to time-lag effects (Epps & Keyghobadi, 2015).
For instance, F ST values typically reflect historic rather than current population structure if populations have not yet reached migrationdrift equilibrium (Whitlock & McCauley, 1999). Moreover, heterozygosity is slow to decline in previously large populations following a genetic health bottleneck (Cornuet & Luikart, 1996), a pattern which has been documented in the western localities of this species range (BEZA and TNP; Parga et al., 2012). Because deforestation has occurred within the last few decades (Brinkmann et al., 2014;Clarke et al., 2015;Gardner & Davies, 2014), it may therefore be too recent to gauge whether habitat loss has negatively impacted genetic diversity and gene flow in this species (Keyghobadi, Roland, Matter, & Strobeck, 2005;Nei, Maruyama, & Chakraborty, 1975).

| Conservation implications
Genetic diversity is lost more rapidly within fragmented and isolated habitats, elevating a species' extinction risk (Frankham, 1995(Frankham, , 2003(Frankham, , 2005. Our results indicate that ring-tailed lemur populations have high levels of genetic diversity. While this is encouraging for the conservation of the species, this may reflect past, not current population processes. Though this species is considered the most ecologically flexible lemur, exploiting anthropogenic landscapes and persisting in small fragments (Cameron & Gould, 2013;Gabriel, 2013;LaFleur & Gould, 2009;Sauther et al., 2006), they are significantly affected by fragmentation and occur at lower densities in poorer habitats (Eppley, Santini, Tinsman, & Donati, 2020;Gabriel, 2013;Kelley, 2011;Sussman et al., 2003). In addition, continued fragmentation and further isolation, coupled with climate change, may prove too much for this historically abundant lemur species. Climatic cycles have been shown to strongly affect mortality rate within this species; a 2-year drought period in southwestern Madagascar resulted in a tenfold increase (3%-27%) in mortality among adult populations in this region (Gould et al., 2003).
Our future work aims to increase sampling efforts in underrepresented and unprotected regions, to reflect the full geographic range of this species and provide a species-wide genetic health assessment (e.g., Calkins & Baden, in review). Moreover, this dataset forms the basis for future landscape genetics analyses which will be used to infer migration and gene flow across the species' remaining range.
Because ring-tailed lemurs are highly terrestrial and are suspected to utilize river basins as dispersal corridors (Goodman et al., 2006), they may have been able to disperse across Madagascar more easily than the more restricted arboreal lemur species. We can use landscape genetic analyses to test this hypothesis to not only better understand what facilitates and impedes gene flow, but also to develop targeted management plans moving forward.

ACK N OWLED G M ENTS
Logistical support in Madagascar was provided by Madagascar Biodiversity Partnership, MICET, ANGAP/MNP. Laboratory support was provided by the Hunter Primate Molecular Ecology Lab. University of New York (ALB).

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