Conservation implications of high gene flow and lack of pronounced spatial genetic structure in elephants supported by contiguous suitable habitat in north‐western India

The western Terai Arc Landscape (wTAL) in Uttarakhand, India, marks the range limit for the Asian elephant in north‐western India. This region has been impacted by land‐use changes and infrastructure expansion for the last seven decades. To evaluate the impact of habitat deterioration on the population structure of elephants in the region, we characterized their genetic diversity and local genetic structure using mitochondrial (D loop) and nuclear DNA (microsatellites; n = 15) markers. We used tissue samples of 114 elephants from five different sub‐populations, collected between 2005 and 2014. The genetic variation was moderate (HO = 0.49–0.55) compared with other Indian elephant populations. Two mtDNA haplotypes were identified without strong spatial patterns across wTAL. Bayesian individual‐based clustering algorithm identified two genetic clusters (K = 2) with high admixture (50% at Q < 0.7) and no spatial adherence. Though K = 1 was not supported by the Bayesian algorithm, multivariate analysis and sibship patterns did not indicate genetic differentiation. The lack of spatial genetic structuring suggests high levels of gene flow, indicating that this population is still panmictic. This suggests that the life history traits of elephants as well as the ecological features of this landscape influence genetic connectivity. However, ongoing land use changes necessitate regular genetic monitoring in wTAL to identify incipient structuring caused by anthropogenic barriers to movement.


| INTRODUCTION
The Asian elephant (Elephas maximus L.) is classified as "Endangered" by the International Union for Conservation of Nature (Williams et al., 2020) and is listed in "Schedule I" of the Wildlife Protection Act (1972) of India, providing robust legal protection.However, Asian elephants now occupy only $5% of their historical range (Olivier, 1978), while habitat fragmentation is a serious threat to populations across their distribution in the 13 range countries (Leimgruber et al., 2003;Williams et al., 2020).A recent study has indicated that elephant habitats in India have lost c. 24.7% and 60% of forest cover and core areas respectively between 1930 and 2013 (Padalia et al., 2019).The recent extinction of an isolated population in east-central India shows that habitat fragmentation is a major threat to elephant populations in the Indian subcontinent (Parida et al., 2022).Habitat loss and fragmentation have been associated with increased incidents of crop raiding and human casualties (Puyravaud et al., 2019;Vasudev et al., 2023), and with higher levels of elephant mortality (Aziz et al., 2016;Joshi & Puri, 2019;Menon & Tiwari, 2019).
The western segment of the Terai belt (Terai Arc Landscape; TAL) in the state of Uttarakhand, north-western India has been undergoing rapid land transformation recently due to expansion of agriculture and settlements, and implementation of infrastructure projects, namely, roads and railway lines (Dinerstein, 2003).Loss of forest cover from 1930 to 2013 has been estimated to be around 46% in this landscape (Padalia et al., 2019).This region supports $2000 elephants, geographically isolated from other elephant populations in India (PED-MoEFCC, 2017).Field studies have indicated that five sub-populations of elephants occur within Uttarakhand (Johnsingh et al., 2006).As this landscape is only 20-25 km wide along the northsouth axis in most places, it restricts the movement of elephants along their established routes (Johnsingh, 1993).The effect of this linear landscape configuration is exacerbated by anthropogenic activities potentially affecting $3 generations of elephants (Wikramanayake et al., 2010).Impacts of such relatively recent fragmentation could be masked by high gene flow (Armbruster et al., 1999) and longer lifespan.On the other hand, the effects of impedance to movement on genetic structure are expected to differ between the sexes as adolescent males tend to disperse from natal herds while the females are philopatric (Desai & Johnsingh, 1995).Therefore, an in-depth study of fine-scale population genetic structure could shed light on the effect of anthropogenic impacts and sex-biased dispersal on the elephant population in this region and provide guidelines for the conservation and management of this population to ensure long-term viability in the face of accelerating climate change and human impacts.
Anthropogenic factors leading to fragmentation of animal populations can lead to loss of genetic diversity and decrease in gene flow (Frankham et al., 2010;Haddad et al., 2015;Young et al., 2000).Habitat fragmentation can also disrupt social interactions due to the Allee effect (Angulo et al., 2018), which in turn can decrease fitness (Luque et al., 2016) and increase the threat of extinction for isolated populations (Wittmann et al., 2018).Genetic data provides insight into critical issues, namely, prioritization of populations for conservation action, assessing census and effective population sizes, estimating adaptive potentials, and detecting hybridizations (Hohenlohe et al., 2021).Though significant progress has been achieved in the field of conservation genetics during the last few decades, the application of genetic data for on-ground management of threatened species remains limited (Willi et al., 2022).
We used 114 elephant tissue samples collected by the State Forest Department of Uttarakhand, India during autopsy examinations of dead elephants from five subpopulations (Johnsingh et al., 2006) across the western Terai Arc region between March 2005 and August 2014.These sub-populations are circumscribed by six major river channels (Figure 1).We isolated nuclear and mitochondrial DNA from these samples to assess (i) patterns of genetic variation and extent of genetic differentiation between the sub-populations, (ii) whether the population genetic structure is sex-specific, (iii) the level of gene flow facilitated by males and females, and (iv) genetic signatures of historical changes in the demographic structure of elephants.Considering the field-based observation of recent fragmentation of this population (Johnsingh et al., 2006) and habitat quality, we tested the null hypothesis of a panmictic population genetic structure with high levels of genetic admixture between proposed sub-populations.

| Study area and sample collection
This study was conducted across the Shivalik Elephant Reserve (SER) in the western Terai Arc Landscape (wTAL) stretching from river Yamuna to river Sarda (west to east) in Uttarakhand, India (Figure 1).The majority of the north-west Indian elephant population occurs within SER comprised of three parallel tracts-the Shivalik hills, the gently sloped Bhabar region with a low water table and the Terai alluvial plains.The total area is approximately 5400 km 2 , which includes $2100 km 2 of protected areas in the form of tiger reserves.The vegetation is an assortment of moist and dry deciduous forests, scrub savannah and grasslands in the floodplains (Johnsingh et al., 2004) and hosts many charismatic endangered fauna, including the Bengal tiger.
We used 114 tissue samples of elephants collected by the Uttarakhand Forest Department, Government of Uttarakhand, India (UKFD) during autopsy examinations of dead elephants between March 2005 and August 2014.We received information on the date of mortality, sex and place of death up to forest beats (smallest administrative unit of Indian forests), and approximate age for the majority of samples.The UKFD veterinarians determined the age by body size criteria and the presence of molars (Roth & Shoshani, 1988;Todd, 2010).Considering sample sizes and evenness, we used the age estimates to classify the samples into two broad age categories (Arivazhagan & Sukumar, 2008), (i) combined calf, juvenile and sub-adults (<15 years) and (ii) adults (>15 years).All samples were preserved in a cryogenic freezer at À80 C.

| DNA extraction, mtDNA sequencing, and microsatellite genotyping
We used DNeasy Blood and Tissue Kit (QIAGEN GmbH) to extract DNA from the samples following the manufacturerspecified protocol with minor modifications, that is, incubation of the samples in a shaking water bath at 56 C and addition of 20 μL proteinase K (20 mg/mL) every 6-8 h till the tissue dissolved in the lysis buffer.After the digestion step, isolation and purification of DNA were carried out in a separate high-yield DNA extraction facility using a spin column-based method.The samples were stored at À20 C until further use.Blank negative controls were used in every batch of DNA extraction to monitor for contamination.
We amplified a 630 base-pair (bp) mtDNA fragment targeting the D-loop region from a subset of 79 individuals using the primer pairs MDL5 and MDL3 (Fernando & Lande, 2000).For each sample, the polymerase chain reaction (PCR) comprised of Thermo 2Â Maxima Hot Start Green PCR Master Mix, 1.5 μg Bovine Serum Albumin, 3 pmol each of the forward and reverse primers, 3 μL tissue DNA isolate and nuclease-free water to make the reaction volume up to 15 μL.The thermocycling profile included an initial denaturation at 95 C for 5 min, followed by a total of 40 cycles of denaturation (95 C for 30 s), annealing (63 C for 40 s) and extension (72 C for 40 s).The final extension was at 72 C for 10 min before holding the reaction at 4 C until removal from thermocycler.All PCR batches included negative and positive controls.We visualized the resulting amplicons using 2% w/v agarose gel electrophoresis with ethidium bromide stain in an ultraviolet light-based gel documentation system.A 100-bp DNA size standard was used to verify the amplicon lengths.We performed enzymatic hydrolyzation on the amplicons using Exonuclease I and Shrimp Alkaline Phosphatase to remove excess primers and nucleotides from the PCR products.The resultant products were subjected to Sanger sequencing from both directions using BigDye Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems) following manufacturer-specified protocols.We purified the sequencing products using the ethanol precipitation method and dissolved the precipitate in Hi-Di™ Formamide (10 μL) (Invitrogen) before denaturation at 95 C for 5 min followed by capillary injection in an ABI 3130 Genetic Analyzer with POP-7 polymer (Invitrogen).
We selected 15 dinucleotide microsatellite markers, namely, EMU03, EMU04, EMU06, EMU07, EMU09, EMU10, EMU11, EMU12, EMU13, EMU14, EMU15, EMU17, LafMS02, LafMS03, and LafMS05, designed for Asian or African elephants (Kongrit et al., 2008;Nyakaana & Arctander, 1998) for genotyping nuclear DNA from the tissue samples.On the basis of the reported allele size range and dye labels, we optimized these markers for co-amplification in four multiplex reactions (Table S2).Reaction compositions were as follows: 5 μL Qiagen Multiplex PCR Master-mix, 10 μg bovine serum albumin (BSA), 1.0 μL of multiplexed primer panel (equal ratio of 10 μM primers; fluorescence-labeled forward primers), 2 μL genomic DNA extract having varying DNA concentration, and nucleasefree water to a volume of 10 μL.The thermal cycling profile for all four panels were similar which included initial denaturation at 95 C for 15 min followed by 40 cycles of denaturation at 95 C for 30 s, annealing at 58 C for 1 min and extension at 72 C for 40 s followed by a final extension at 60 C for 30 min before hold at 4 C.All reaction batches included negative and positive controls to track reagent contamination and PCR failure.To screen for PCR success, a fraction of the amplicons ($10%) across panels were run in 2% w/v ethidium bromide-stained agarose gel alongside a 100-bp DNA ladder and then visualized in a geldocumentation system.Each amplicon (1.0 μL) was then added to Hi-Di™ Formamide (8.93 μL) and GeneScan™ 500 LIZ™ size marker (0.07 μL) (Invitrogen), mixed thoroughly and denatured at 95 C for 5 min before fragment analysis in an ABI Prism 3130 Genetic Analyzer using POP-7 polymer (Invitrogen).

| Data analyses
The sequence chromatograms were aligned using the software SEQUENCHER v5.4.6 (Gene Codes Corporations) while verifying each mutation manually.The ends of the raw sequences were trimmed to uniform length and discrepancies in base calls were resolved by manual scrutiny.We compared the clean data with the NCBI GenBank sequence repository for the identification of the haplotypes.
We used GENEMAPPER v3.7 (Applied Biosystems) to examine the fragment analysis electropherograms followed by automated allele scoring and manually verifying each allele call.The Excel macro AUTOBIN v0.9 Excel Macro (Salin, 2010) was used to bin the raw allele sizing data.
We computed the mean number of alleles per locus (MNA), observed heterozygosity (H O ), expected heterozygosity (H E ), unbiased expected heterozygosity (uH E ) using the R package diveRsity (Keenan et al., 2013;R Core Team, 2019).The 95% confidence interval of inbreeding coefficient (F IS ) values were estimated by bootstrap re-samplings (n = 999).We used rarefied measures of allelic richness (AR) that accounted for samplesize variability across EUs.We computed all population genetic parameters for each EU as well as on pooled data across all EUs (i.e., all individuals pooled into one group).We tested allele frequencies at each locus for departure from Hardy-Weinberg expectations across EUs.We assessed the presence of null alleles at each locus using FreeNA (Chapuis & Estoup, 2007), employing the maximum-likelihood algorithm (Dempster et al., 1977).We checked the microsatellite data for potential genotyping errors, namely, stutter bands, and large allele dropouts using MicroChecker (van Oosterhout et al., 2004).
For the microsatellite data, we executed a hierarchical analysis of molecular variance (AMOVA) (Excoffier et al., 1992;Michdaki & Excoffier, 1995) using GenAlEx (Peakall & Smouse, 2006, 2012) to determine the significance of genetic variability across the sampling sites and within geographic locations.As differential diversity levels strongly influence pairwise F ST (Charlesworth, 1998), we also estimated pairwise D Jost values which are independent of within-population diversity (Jost, 2008).These two metrics of population differentiation were estimated using the diveRsity R package (Keenan et al., 2013;R Core Team, 2019) along with their 95% confidence intervals by bootstrapping (n = 999).We computed pairwise F ST and D Jost for the entire dataset, as well as separately for males and females, to determine if there is sex-biased structure or gene flow patterns.Genetic differentiation was also examined across groups of two age classes, that is, subadults including calf and juvenile and adults.
We investigated population genetic structure using the individual-based clustering (IBC) algorithm as implemented in STRUCTURE v. 2.3.4 (Pritchard et al., 2000) with assumptions of correlated allele frequencies and admixture without a-priory information on sampling locations.We tested for the number of genetic clusters (K) in the study population (K = 1 to 10) with 20 iterations for each K, while a total of 5 Â 10 5 Markov Chain Monte Carlo (MCMC) simulations were computed after discarding the initial 5 Â 10 4 simulations as burn-in for each independent run.The optimal K was identified by computing the rate of log probability change (ΔK) of the Bayesian simulation fitting the observed data across consecutive values of K (Evanno et al., 2005).We used the model choice statistics L(K) (Pritchard et al., 2000), which estimates the posterior probability of the observed data for any given K value to support the inference from ΔK.A web-based tool, STRUCTURE HARVESTER (Earl & vonHoldt, 2012), was used to calculate model selection statistics while CLUMPAK (Kopelman et al., 2015) was used for model averaging for the K across the independent runs.Expecting a sex-biased gene-flow pattern, we ran similar STRUCTURE algorithms using male and female individuals as two different groups.We also assessed the population structure for elephants using individuals sampled at different time points (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) to alleviate temporal effects.
We applied multivariate discriminant analysis of principal components (DAPC) for identifying clusters of genetically related individuals.Unlike Bayesian clustering algorithms, DAPC is devoid of any underlying population genetic model, namely, Hardy-Weinberg Equilibrium (HWE) (Vergara et al., 2015).DAPC uses a k-means clustering algorithm-transforming allele frequencies into uncorrelated components using principal component analysis (PCA) and then, discriminant analysis (DA) to maximize the genetic differential between the groups keeping intragroup differentiation to the minimum (Jombart, 2008).DAPC demonstrates better resolution than STRUCTURE at characterizing population segmentation (Jombart et al., 2010).We implemented DAPC using the R package ADEGENET (Jombart, 2008; R Core Team, 2019).

| Isolation by distance
We used Mantel test in GenAlEx v6.5 (Peakall & Smouse, 2006, 2012) to test for the presence of isolation by distance (IBD) across the study area.The analysis examined congruence between a pairwise genetic distance matrix and a Euclidian geographic distance matrix of the samples with 999 permutations to assess significance.

| Sibship analysis
We used the program COLONY v.2.0.6.6 (Jones & Wang, 2010) to identify sibling pairs among the sampled individuals under a full-likelihood framework.The assumptions of polygamy and inbreeding were considered for the three independent runs of medium length having weak sibship priors.

| Test for changes in effective population size, and dating with MSVAR
We used MSVAR 0.4 to test whether a signature of population bottleneck was present in the wTAL elephant population and estimated the timing of such historic demographic changes using coalescent modeling of microsatellite genotyped data following Sharma et al. (2012).To test for a genetic bottleneck MSVAR 0.4 runs were carried out under different models of population size changes (Beaumont, 1999).MSVAR 1.3 (Storz & Beaumont, 2002) estimates the current (N 0 ) and the ancestral effective population sizes (N 1 ) and the time (T; in generations) since the effective population size changed.To control the confounding effects of population structure and gene flow on the false assignment of bottlenecks (Chikhi et al., 2010), individuals sampled at different time points were not separately analyzed.Instead, a pooled data set was constructed by resampling the individuals randomly, for both MSVAR version 0.4 and MSVAR version 1.3.As the program is very slow, only randomly selected individuals (n = 80) from each of the 10 time points (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) were used to create a pooled data set for which we inferred the population demographic history.We carried out three runs with wide uninformative priors and used a 25-year generation length (Williams et al., 2020).Different sets of priors were used to test their influence on the posteriors.Supplementary materials (Figure S4) contain additional details on MSVAR analyses.

| Information on elephant spatial demography
To corroborate elephant life history with our results, we used demographic data based on field survey from across wTAL (UKFD, 2016) to visualize forest division-wise differential spatial distribution patterns of adult male and female individuals as well as male to female ratios using ArcMap 10.2 (ESRI, Redlands, CA, USA).

| RESULTS
3.1 | Genetic diversity in the wTAL elephant population

| Mitochondrial DNA (mtDNA) analyses
We successfully amplified and sequenced 440 bp of mtDNA D-loop in 79 tissue samples.On comparison with the available sequences at GenBank, we confirmed the presence of two previously described haplotypes AC and AH within the sampled individuals while the AH haplotype had previously not been reported from northwestern India (Vidya et al., 2009).The two haplotypes differ from each other by only two base-pair mutational events.There was no spatial pattern of mtDNA haplotype distribution along the landscape although the frequency of AC (74%) was much higher than AH.The haplotype diversity (h) in TAL was 0.397 ± 0.002 while nucleotide diversity (π) was 0.002 ± 0.000.

| Microsatellite analyses
We successfully genotyped all 114 tissue samples of elephants at 15 microsatellite loci.Of the 114 samples analyzed across the landscape, 79 were male individuals, and 35 were females.Total gaps in the multi-locus genotype data were 4.5%.Although significant deviation from Hardy-Weinberg equilibrium (HWE) was present at most loci, the degree of deviation varied across EUs (Table S3).The results were similar when analyses were performed on the pooled group instead of EUs and 12 out of 15 loci deviated from HWE (Table S3).When tested separately, 10 loci deviated from HWE for females while males were out of HWE for 12 loci (data not shown).Linkage disequilibrium was not evident between any locus pairs.Null allele frequencies were intermediate at all loci and their frequency estimates ranged between 4% and 18% (Table S3).We did not find any evidence of large allele dropouts and stutter band-related scoring errors.
We found congruent genetic variability patterns within and between the five EUs.The number of alleles was between four and 14 per locus.The observed heterozygosity (H O ) values were similar for all the EUs and ranged from 0.49 to 0.55, whereas unbiased expected heterozygosity (uH E ) ranged from 0.62 to 0.70 (Table 1).
The mean number of alleles for each EU ranged between 5.00 and 6.20.Applying correction for variable sample sizes from each sampling unit using rarefaction, allelic richness (AR) was found to be nearly similar across the EUs (between 4.04 and 4.63).The F IS values were all positive in the range of 0.17-0.25 and were statistically significant in all EUs except for EU4.These positive values indicate a heterozygote deficiency at all 15 polymorphic loci across the population.

| F-statistics and analysis of molecular variance
Both D Jost and F ST showed low estimates and indicated non-significant genetic differentiation between all pairs of EUs for the microsatellite data across all loci (Tables 2  and S4).Both statistical measures remained low and nonsignificant when the sexes were tested separately (Table 2).We also found non-significant differentiation between the two age classes, that is, between the subadults including calf and juvenile and adults (F ST <0.05; Hartl & Clark, 1997).
There were no significant differences between F ST and D Jost across the five EUs for any of the time points.Similarly, AMOVA results showed that most of the mtDNA and nuclear DNA-based genetic variation could be attributed to differences between individuals within EUs.

| Bayesian clustering methods
STRUCTURE detected a weak genetic substructure, that is, the presence of two clusters in the study area at the uppermost level of the hierarchy with no apparent geographic pattern (Figure 2a, b).Although K = 2 (mean LnProb = À4456) was obtained using the ΔK criterion (Figure S1a, b), the individual assignment showed only 50% of the 114 individuals belonged to any of these two clusters with the probability of assignment Q > 0.7 while the other half of the individuals were not assigned to any genetic cluster (Figure 2a) hence categorized as the admixture of the two ancestries.This observation supports an absence of clear population structuring while indicating that the ΔK based estimate is unable to define the true K (K = 1).Contrary to our prior expectations, STRUCTURE analyses displayed a similar pattern of two clusters with no spatial adherence (Figure 2c) for male and female individuals.

| Multivariate analyses
The multivariate DAPC displayed overlaps among the EUs sampled, augmenting the finding of a panmictic population pattern (Figure 3).DAPC with prior information on sampled EUs produced highly overlapping clusters except for EU4 (Ramnagar, Terai Central and part of Terai West FD), the centroid of which did not lie inside the ellipses of any other EU clusters.The DAPC with either We did not compute measures of differentiation for females of EU5 due to low sample size (n = 1).
males or females did not reveal a sex-biased pattern, indicating a lack of genetic structure (data not shown).

| Isolation-by-distance
We did not find any evidence for isolation-by-distance across the study area, as the effect of geographic distance on genetic differentiation was not statistically significant (Rxy = 0.06, p = .11)(Figure S2).The IBD pattern for males was weak but significant (Rxy = 0.10, p = .03)while it was non-significant in females (Rxy = 0.11, p = .13).This difference in the level of significance between the sexes could be attributed to the unequal sample size for male and female individuals.

| Relatedness analyses
Sibship analyses with COLONY identified in total 12 fullsibling relationships, out of which, five pairs were found within the EUs, seven across the EUs and five across the adjacent EUs.Moreover, the proportions of full-and halfsiblings to all sampled individuals within each of the EUs were less than the proportions observed across the study area (Figure S3).

| Population size changes and dating
The demographic analysis using MSVAR consistently supported a clear population decline under the exponential model of population size change.The posterior distributions of the effective population size change for the elephant population are shown in Figure S4a.Across all runs, the median log 10 (r) was À0.91 in the dataset.The negative value indicates an order of magnitude decline in the effective size of the population and does not support an increase or a stable population size.MSVAR 1.3 simulations also showed evidence for population decline (Figure S4b) consistent with a small current effective population size.The estimates of current (median N 0 ) and ancestral population size (median N 1 ) were 2100 (90% highest posterior density; HPD = 425/11,708) and 13,000 (90% HPD = 2756/66,205), respectively.Estimates of the time of population contraction (Figure S4c) showed support for population decline occurring long before the last two centuries of habitat fragmentation in north-western India (median T = 30,000 years before present (YBP) or ≥1200 generations ago).Wide 90% HPD limits were observed, and the broad values ranged between 1778 and 465,411 years.

| Spatial demography of elephants
As per the 2015 census data (UKFD, 2016), the spatial demographic profile of adult elephants was highly skewed and localized (Figure S5a,b).CTR alone contained 51% of males and 58% of total female individuals, respectively.Rajaji National Park (this excludes buffer areas of RTR; Lansdowne and Shyampur of Haridwar) and Lansdowne FD were inhabited by 22% of all males in wTAL and 27% of females.Haridwar and Terai Central FD had the highest adult male to adult female ratio (1:2.08 and 1:1.25, respectively), followed by Kalsi Soil Conservation Division (1:1) and Terai West FD (1:1) (Figure S5c).A male to female ratio between 1:0.55 and 1:0.61 was recorded in Terai East, Dehradun and Ramnagar FD.Lansdowne FD, CTR, Haldwani FD and Rajaji NP have male to female proportion of <0.5 (Figure S5c).

| DISCUSSION
Our results show that elephant sub-populations (EUs) in wTAL which is isolated from other mainland populations in India, are panmictic and not genetically subdivided, supporting the observations made in an earlier study (De et al., 2021).African elephants show that habitat modification influences spatial genetic structuring through changes in the species' life-history traits and behavior (Archie & Chiyo, 2012).Genetic diversity and the lack of differentiation in wTAL indicate that recent land use changes have not yet affected gene flow among EUs.The observed level of genetic diversity is comparable to wellconnected populations inhabiting a large landscape, for example, Nilgiri Biosphere Reserve ($5500 km 2 ), southern India (De et al., 2021;Vidya et al., 2005).However, a recent study by Lorenzana et al. (2020) on jaguars (Panthera onca) demonstrated moderate genetic diversity in the fragmented (since early 1950s) and declining (since late 1980s) (Beisiegel et al., 2012) Atlantic Forest population despite evidence of strong genetic drift while retaining heterozygosity comparable to the Amazon population having high demographic connectivity.Therefore, the severe anthropogenic pressure across our study area can potentially alter and eventually disconnect these populations over the long term in spite of the current levels of genetic diversity.Hence, our results must be interpreted with caution because a number of studies in mammals reported declines in genetic diversity across landscapes due to habitat modifications (Curry et al., 2021;Fünfstück & Vigilant, 2015;Lino et al., 2019).

| Genetic diversity and population genetic structure
The results showed significant departures from Hardy-Weinberg expectations for the microsatellite loci examined.We attribute this deviation to non-random mating and null alleles in the population.Deviation from random mating is expected and is tightly linked to the interplay of the dispersal of males away from their natal group and the mate choice of the individuals (Chelliah & Sukumar, 2015;Fernando & Lande, 2000;Vidya & Sukumar, 2005).The presence of null alleles has been documented in population genetics studies of elephants (De et al., 2021;Goossens et al., 2016;Parida et al., 2022).Null alleles are most likely to be found in species with large effective population sizes due to high nucleotide diversity in flanking regions of target microsatellites (Chapuis & Estoup, 2007).
The levels of nuclear genetic diversity (H O = 0.50; H E = 0.65) observed were moderate, supporting the previous findings of De et al. (2021) for the same population (H O = 0.46, H E = 0.62).Also, a global comparison revealed either a higher or comparable level of genetic diversity in wTAL elephants, except for Nepal and Lao PDR (Ahlering et al., 2011;De et al., 2021;Flagstad et al., 2012;Goossens et al., 2016;Vidya et al., 2005;Zhang et al., 2015).However, this comparison should be treated with caution as these studies amplified different combinations and numbers of microsatellite loci (Lorenzana et al., 2020).We also found that the elephants retained genetic variations over a 10-year period which is not unexpected as these samples cover less than a generation time of the Asian elephant.
Haplotype and nucleotide diversity for mtDNA was low.In addition, the haplotype AH, found for the first time in wTAL population, is known to co-occur in the northeast Indian elephant population along with haplotype AC (Vidya et al., 2009).This finding confirms historical connectivity between the wTAL and northeastern elephant populations and can be indicative of allopatric divergence followed by admixture due to natural processes (Fleischer et al., 2001;Vidya et al., 2009).
Bayesian and multivariate clustering algorithms concordantly showed no population structuring and indicated the presence of admixture individuals all over the study area.This is further supported by analyses using multiple metrics, that is, pairwise F ST , and D Jost which indicated high to moderate gene flow across EUs including different age-sex classes.Similar findings were previously recorded for other elephant populations in India, for example, the Nilgiris population in southern India and the populations in the north bank and southwestsouthcentral bank population across river Brahmaputra in northeast India (Vidya et al., 2005).
Sex-biased dispersal in wTAL elephants was not supported by our results.The population genetic indices were also similar for both females and males.This result suggests that dispersal for either sex is not limited, and the genetic structure is weaker than many other elephant populations (Fernando & Lande, 2000;Vidya & Sukumar, 2005).Our results also corroborate an earlier radio-telemetry-based study in wTAL, where the home range size varied from 188 to 407 km 2 , and there was no significant difference in the average home range size of male and female elephants (Williams et al., 2008).However, the sample size of females in our study is relatively low, and the results might change upon further sampling.
In addition, the following two factors may have also influenced genetic diversity and gene flow among the EUs in wTAL: 1. Habitat quality and relatively recent origin of anthropogenic fragmentation Western TAL is a highly productive and contiguous wildlife habitat due to high plant productivity and water availability.Most of the anthropogenic changes in this area occurred since the late 1950s spanning $3 generations of elephants (Wikramanayake et al., 2010).Simulations using the generation time for African savanna elephants (25 years) suggest that at least 100 years (4 generations) following a steep decline in gene flow is required to detect an alteration in population genetic structure (Epps et al., 2013).Hence, potential barriers to gene flow between the EUs in wTAL may not have been in place long enough to detect a strong genetic signature.In this context, the observed patterns might reflect conditions prior to habitat fragmentation, and hence, should be interpreted with caution.

Landscape configuration, ranging behavior, and presence of congregation sites
Genetic structure is linked to life history aspects, such as foraging strategy, movement patterns and habitat association that influence dispersal and gene flow.Reduction in connectivity due to habitat fragmentation impedes gene flow and, in turn, affects allele frequencies.It can lead to a loss of genetic diversity and fixation of certain alleles, with potential deleterious effects on populations.However, in our case, the landscape contiguity within the EUs, extensive ranging behavior (Williams et al., 2008), and robust sex ratio between males and females (Figure S5) may have facilitated panmixia within wTAL.
Geographic distance is expected to have a weaker effect on gene flow in elephants due to their high dispersal potential.The distance between the farthest populations (250 km) in wTAL is within the dispersal capabilities of elephants (Williams et al., 2008).Hence, the observed genetic integrity in this landscape may have been mainly manifested by similar patterns of mating and dispersal as documented earlier (Vidya & Sukumar, 2005).
Our results also suggest that the presence of rivers and roads did not influence patterns of relatedness among individuals.Results from sib-ship analyses suggest that despite the total numbers of full-sibs being low, their high proportions across adjacent EUs indicate uninterrupted gene flow.The absence of IBD between individuals suggests that both males and females have similar movement patterns across the landscape.While rivers may slow down the movement of elephants during high water levels in the monsoon, they do not act as absolute barriers to broadscale movement as animals have been observed crossing the river during strong water flow in the Ganga (Nigam, personal communication).The space use of elephants is known to be shaped by the distribution of water sources (Eisenberg & Lockhart, 1972).A number of perennial water sources (Johnsingh et al., 2004) located inside wTAL support elephant movement.These rivers, pools and reservoirs serve as elephant congregation points during summer, facilitating gene flow across EUs by intermixing of individuals from one clan to another (Johnsingh et al., 2004).However, the recent update in the tourism policy of the Government of Uttarakhand (https:// uttarakhandtourism.gov.in/sites/default/files/2023-06/Tourism-policy-2023.pdf) and therefore, likely improvement of road infrastructure may affect genetic cohesiveness of elephants in the future as intense road traffic impact genetic structure of mammals (Habrich et al., 2021).

| Population size changes and dating
Our results suggested a historical population contraction while there was no support for a substantial recent bottleneck.Based on pooled data analysis across wTAL, a onefold decrease in population size (N 0 , $2100; N 1 , $13,000) was detected before the beginning of the Last Glacial Maxima, about 30,000 years ago.However, genetic bottleneck estimates were unreliable, as the log(t a ) point estimates were consistently higher and factors like mutation rate could influence and possibly bias the dating estimates (Girod et al., 2011).
Studies on Asian elephants have shown that Pleistocene climatic oscillations have been paramount in shaping Asian elephant biogeography (Fleischer et al., 2001;Sharma et al., 2018;Vidya et al., 2009).However, elephants were documented to be widely distributed across the Indian subcontinent and extended through Afghanistan and Syria to the West.These populations were extirpated through anthropogenic impact and possibly climate change (Trautmann, 2015).Hence, the genetic signature of the bottleneck identified in our data could be related to an ancient event.This is also supported by the presence of only two maternal lineages (mtDNA haplotypes) at different frequencies in wTAL elephant population.At glacial maxima, the major proportion of the northern and western Indian landform was dry, and central TAL was drier compared to the western and eastern extremities of this region (Kale et al., 2003).More recently (i.e., Quaternary period), probably forest cover shrank and expanded due to cycles of multiple dry and wet periods (Karanth, 2003;Morley & Morley, 2022).Hence, eastern and wTAL could have been a refugia for elephants from west of Yamuna to the Brahmaputra in northeast India.Similar biogeographic processes have also affected the demographic history of other sympatric species found in the same ecosystem, for example, Swamp deer and greater one-horned rhinoceros, which are only found in the eastern and western part of TAL (Ellis & Talukdar, 2019;Groves, 1982).Though the observed results across different elephant populations (Sharma et al., 2018) indicate robust genetic signature of demographic changes, future studies using larger, highresolution datasets (e.g., genome-wide SNPs) and methods based on summary statistics (e.g., approximate Bayesian computation) can potentially provide narrower confidence intervals for parameters of interest.

| Implications for conservation and conclusion
The wTAL elephant population is isolated from other mainland populations in India.This population has adapted to the unique bioclimatic region, differs in preferred forage species and, therefore, of high conservation priority (Johnsingh et al., 2004;Williams et al., 2005).The lack of genetic structure in the wTAL provides robust support for a panmictic elephant population in the study area.However, ongoing habitat destruction and barriers to movement in important dispersal corridors (Johnsingh et al., 2004), likely to be exacerbated by the proposed tourism policy of the state, may contribute to population isolation and loss of genetic diversity over a period.For the long-term viability of the elephant population and retention of genetic structure in wTAL, it is important to: 1. Conduct similar, harmonized studies in other elephant populations in India and Nepal to identify habitat linkages that could be restored (Davidar et al., 2023).2. Carry out systematic genetic monitoring of elephants of all age classes every 10 years to detect genetic changes over generations.3. Protect movement paths and congregation sites for maintaining gene flow.4. Monitor the vegetational changes at regular intervals using high-resolution geospatial data and ground truthing to detect and thereby avert habitat loss and fragmentation.

F
I G U R E 3 Clusters of different elephant population units (EUs) in the discriminant space using discriminant analysis of principal components (DAPCs) in the individuals across Shivalik Elephant Reserve, western Terai Arc Landscape, India.The left inset figure represents percentages of cumulated variance (PCA eigenvalues) explained by sequential principal components and the right inset shows the discriminant analysis eigenvalues.
Genetic diversity in the elephant population units (EUs) of Shivalik Elephant Reserve, western Terai Arc Landscape, India.
T A B L E 1 a Denotes statistically significant results.b Refer Figure 1 for geographic location.