A Comparative Analysis of Genetic Diversity and Structure in Jaguars (Panthera onca), Pumas (Puma concolor), and Ocelots (Leopardus pardalis) in Fragmented Landscapes of a Critical Mesoamerican Linkage Zone

With increasing anthropogenic impact and landscape change, terrestrial carnivore populations are becoming more fragmented. Thus, it is crucial to genetically monitor wild carnivores and quantify changes in genetic diversity and gene flow in response to these threats. This study combined the use of scat detector dogs and molecular scatology to conduct the first genetic study on wild populations of multiple Neotropical felids coexisting across a fragmented landscape in Belize, Central America. We analyzed data from 14 polymorphic microsatellite loci in 1053 scat samples collected from wild jaguars (Panthera onca), pumas (Puma concolor), and ocelots (Leopardus pardalis). We assessed levels of genetic diversity, defined potential genetic clusters, and examined gene flow for the three target species on a countrywide scale using a combination of individual- and population-based analyses. Wild felids in Belize showed moderate levels of genetic variation, with jaguars having the lowest diversity estimates (HE = 0.57 ± 0.02; AR = 3.36 ± 0.09), followed by pumas (HE = 0.57 ± 0.08; AR = 4.20 ± 0.16), and ocelots (HE = 0.63 ± 0.03; AR = 4.16 ± 0.08). We observed low to moderate levels of genetic differentiation for all three target species, with jaguars showing the lowest degree of genetic subdivision across the country, followed by ocelots and pumas. Although levels of genetic diversity and gene flow were still fairly high, we detected evidence of fine-scale genetic subdivision, indicating that levels of genetic connectivity for wild felids in Belize are likely to decrease if habitat loss and fragmentation continue at the current rate. Our study demonstrates the value of understanding fine-scale patterns of gene flow in multiple co-occurring felid species of conservation concern, which is vital for wildlife movement corridor planning and prioritizing future conservation and management efforts within human-impacted landscapes.

With increasing anthropogenic impact and landscape change, terrestrial carnivore populations are becoming more fragmented. Thus, it is crucial to genetically monitor wild carnivores and quantify changes in genetic diversity and gene flow in response to these threats. This study combined the use of scat detector dogs and molecular scatology to conduct the first genetic study on wild populations of multiple Neotropical felids coexisting across a fragmented landscape in Belize, Central America. We analyzed data from 14 polymorphic microsatellite loci in 1053 scat samples collected from wild jaguars (Panthera onca), pumas (Puma concolor), and ocelots (Leopardus pardalis). We assessed levels of genetic diversity, defined potential genetic clusters, and examined gene flow for the three target species on a countrywide scale using a combination of individual-and population-based analyses. Wild felids in Belize showed moderate levels of genetic variation, with jaguars having the lowest diversity estimates (H E = 0.57 ± 0.02; A R = 3.36 ± 0.09), followed by pumas (H E = 0.57 ± 0.08; A R = 4.20 ± 0.16), and ocelots (H E = 0.63 ± 0.03; A R = 4.16 ± 0.08). We observed low to moderate levels of genetic differentiation for all three target species, with jaguars showing the lowest degree of genetic subdivision across the country, followed by ocelots and pumas. Although levels of genetic diversity and gene flow were still fairly high, we detected evidence of fine-scale genetic subdivision, indicating that levels of genetic connectivity for wild felids in Belize are likely to decrease if habitat loss and fragmentation continue at the current rate. Our study demonstrates the value of understanding fine-scale patterns of gene flow in multiple co-occurring felid species of conservation concern, which Introduction Carnivores are particularly vulnerable to loss and fragmentation of natural habitats because of their space and prey requirements, occurrence at low densities, and dependence on forest and/ or habitats that provide sufficient cover [1,2]. Landscape changes due to forest loss, degradation, and anthropogenic development can severely impact animal movement and the degree of genetic connectivity, consequently decreasing reproductive fitness and adaptive potential (e.g., [3,4,5]). Assessing genetic diversity and connectivity is important for understanding responses of wild populations to fragmented landscapes and anthropogenic disturbance. Sensitivity to these threats is species-specific, depending on dispersal characteristics, habitat specialization, trophic level, and other ecological traits (e.g., [2,6]). Thus, comparative analysis of genetic connectivity among multiple species in the same landscape is particularly valuable as it may reveal common factors driving gene flow.
Neotropical carnivore guilds, which are among the most threatened worldwide [7] are understudied. The jaguar (Panthera onca), the largest Neotropical felid, which co-occurs with several other felids [e.g., puma (Puma concolor), ocelot (Leopardus pardalis)], has been extirpated from more than half of its historic range during the last 100 years, and its distribution continues to contract, mainly due to severe deforestation and habitat fragmentation, direct persecution, and loss of main prey species (e.g., [1,[8][9][10]). Mesoamerican jaguar populations are reduced to one third of their historic range, and 75% of the remaining jaguar populations have declined in numbers and are potentially highly fragmented [10,11]. To prevent further loss of these predators, it is crucial to conduct efficient monitoring of their wild populations, and gather valuable scientific data to inform conservation and management of this internationally protected species (listed as 'Near Threatened' under the International Union for Conservation of Nature [IUCN], as an Appendix I species under the Convention on International Trade of Endangered Species of Wild Fauna and Flora [CITES]). Co-occurring pumas and ocelots face similar challenges to their survival. Pumas have also disappeared from large portions of their historic range and are listed as a 'Least Concern' species under IUCN and as an Appendix I (eastern and Central American subspecies) and an Appendix II (remaining subspecies) species under CITES [1,12]. Ocelots, the only medium-sized felid found in the Neotropics, are listed as a 'Least Concern' species under IUCN and are included as an Appendix I species under CITES [1,8,13]. Wild populations of all three Neotropical felid species are considered understudied across Mesoamerica, a region that faces one of the highest deforestation rates worldwide [14,15].
Molecular genetic approaches provide powerful tools to assess the conservation status of multiple species of concern by monitoring genetic diversity and connectivity. Advancements in noninvasive genetic monitoring techniques, including the use of molecular scatology and scat detector dogs [16][17][18][19][20][21], have made it more powerful and feasible to genetically study multiple elusive forest carnivores simultaneously in tropical environments, without having to physically capture and handle animals. The use of professionally trained scat detector dogs significantly increased scat-collection rates for several recent molecular scatology studies conducted in tropical environments (e.g., [17,20,21]). Implementing conservation genetic studies for multiple species simultaneously is more efficient and enables better identification of landscapes of high conservation concern to all species, thus aiding in wildlife movement corridor planning and conservation management efforts (e.g., [5,22]). Nonetheless, molecular population genetics studies of jaguars and other Neotropical felids are still relatively rare (e.g., [23,[24][25][26][27][28]). Only a handful of current noninvasive genetic studies sample wild populations of multiple felid species and most of them are based in South America (e.g., [20,21,28,29,30]).
To address this knowledge gap, we conducted a 4-year (2007-2010) noninvasive genetic study on jaguars and two co-occurring felids, pumas and ocelots, using fecal DNA samples collected with the aid of a scat detector dog across several study sites in Belize, Central America. Belize and its neighboring countries, Guatemala, and Mexico, are part of La Selva Maya (the Maya Forest), the largest remaining tropical forest in Mesoamerica, and the largest intact forest north of the Amazon, representing a critical link in the Mesoamerican Biological Corridor. Like other Mesoamerican countries, Belize has been experiencing widespread land conversion due to agricultural and urban development, logging, cattle-ranching, and natural disasters (e.g., hurricanes, forest fires), that have increased deforestation rates (2.3% yearly) above the Central American average (1.2% yearly) over the last two decades [31,32]. We hypothesized that genetic diversity and connectivity would be reduced for felids in the most northern and southern protected areas, which are more isolated and face one of the highest deforestation rates within the country [31]. We also expected to see lowered levels of genetic connectivity between protected areas separated by human-constructed barriers (e.g., intensive agriculture, roads, and urban development) such as in central and northern Belize. Furthermore, we hypothesized that jaguars would exhibit higher levels of genetic differentiation than pumas, since pumas, especially males reportedly are long-distance dispersers and are more likely to move through disturbed and fragmented areas than jaguars (e.g., [33,34]). To test these hypotheses and assess the conservation status of wild felids in Belize, we specifically aimed to: (a) estimate levels of genetic diversity within different regions of the country, (b) examine patterns and spatial scale of genetic structure using multiple individual-and population-based analyses, and (c) assess contemporary gene flow and dispersal movements for all three target species. We discuss the implications of our results to support current and future conservation and management efforts, including wildlife movement corridors for wild felids in Belize.

Ethics Statement
Permission to undertake fieldwork and sample collection was obtained from the Belize Forest Department. No review from the ethics committee was required, as our research work applied noninvasive genetic sampling methods and did not involve any direct manipulation or disturbance of animals. barriers to gene flow, including natural (e.g., mountain ranges) and anthropogenic landscape (e.g., urban and agricultural areas, roads) features. For jaguars and pumas, we evaluated five geographical regions in Belize, which included the following study sites: north (FB, SNP), north-central (RBCMA, BF, MFR, TSB), central (MPR, CFRNP, BFR, HVPR), south-central (CBWS), and south (BNR, BC, GS, MH, STNP) (Fig 1). For ocelots, due to low sample size, we evaluated felids detected in the north (FB, SNP, RBCMA, BF, MFR, TSB), and south (MPR, CFRNP, BFR, HVPR, BNR, BC, GS, MH, STNP) (Fig 1). Across study sites, elevation ranges from 0 to 1120 m, and mean annual rainfall varies from 1524 mm in the north to 4064 mm in the south with a pronounced wet season from June to December. Average annual temperatures fluctuate between 17.7 and 31.3°C. A high diversity of native habitat types occurs within the study sites, including lowland and submontane broadleaf moist and wet forests (both primary and secondary growth), lowland and submontane pine forests, mangrove and littoral forests, lowland savannah, shrub land, and wetland swamps called bajos.

Fecal DNA sampling and genotyping
Fecal samples were detected opportunistically in the field by using a professionally trained scat detector dog (PackLeader LLC, Gig Harbor, WA, US) following the study design and sampling protocols described by Wultsch, Waits [35]. DNA extraction was carried out for all fecal samples detected by the scat detector dog following protocols in Wultsch, Waits [20]. Fecal samples were assigned to feline species using species-specific microsatellite alleles identified in 14 polymorphic microsatellite loci, confirmed by mitochondrial DNA sequencing [20]. The microsatellite genotypes were also used to identify individuals, and by combining and arranging these loci with two additional DNA markers associated with the Y sex chromosome carried by males (Zn, Zn-finger; Amel, Amelogenin), but not by females [36], into three multiplexes (S1 Materials and Methods), we were also able to determine the sex of individuals. All extractions and polymerase chain reactions (PCR) included negative controls. To finalize consensus genotypes, we conducted an average of 4.29 (SD ± 1.84) PCR replicates per locus and sample. In addition, a multi-tube approach was used, where at least three identical PCR results were required for homozygote genotypes, and each allele had to be observed at least twice to confirm heterozygote genotypes. The genetic analysis was conducted at a laboratory specialized in the analysis of noninvasive genetic samples (Laboratory for Ecological, Evolutionary, and Conservation Genetics, University of Idaho, Moscow, ID, US).

Genetic diversity
We used GenAlEx, version 6.41 [37] to assess genetic variation per species and geographical region (north, north-central, central, south-central, and south for jaguars and pumas; north and south for ocelots) across all loci by estimating the number of alleles (N A ), the number of private alleles (A P ), observed (H O ) and expected heterozygosities (H E ), and inbreeding coefficients (F IS ). Additionally, we determined allelic richness (A R ) using the rarefaction method with HP-RARE, version 1.0 [38]. Statistical differences between groups were evaluated using non-parametric Kruskal-Wallis rank-sum tests in R, version 3.1.3 [39]. We tested for linkage disequilibrium and departures from Hardy-Weinberg proportions using exact tests in GENE-POP, version 4.1 [40] with default settings for Markov chain parameters. Results for multiple significance tests were adjusted by applying a sequential Bonferroni approach [41].

Fine-scale genetic structure and contemporary gene flow
To determine whether felid species in Belize represent single panmictic populations or several subdivided populations, we used multiple approaches to assess genetic structure and levels of gene flow indirectly. We used F-statistics to examine the degree of genetic differentiation by calculating pairwise F ST values [42] between groups of felids detected across the five geographical regions. To account for high mutation rates of microsatellites, which may underestimate genetic differentiation using F ST values [43], we also calculated pairwise standardized G' ST estimates [44] using GenAlEx, version 6.41 [37]. Genetic differentiation was further examined via analysis of molecular variation (AMOVA, [45]) in Arlequin, version 3.5.1.3 [46]. Different hypotheses about hierarchical genetic differentiation for each species were tested, including genetic subdivision by sampling geographical To complement these analyses, we also conducted exploratory multivariate analyses, including discriminant analysis of principal components (DAPC), using adegenet, version 1.4.2 [47] in R, version 3.1.3 [39].
To assess gene flow more directly and gain insight into contemporary rates of genetic change, we conducted several different types of individual-based assignment tests using aspatial [48,49] and spatial Bayesian clustering techniques [50][51][52]. We applied aspatial Bayesian clustering in STRUCTURE, version 2.3.4 [48,49] to determine the optimal number of genetic groups/clusters (K) by running 2,000,000 Markov chain Monte Carlo (MCMC) iterations after a burn-in of 200,000 replicates with 10 independent runs per K (ranging from 1 to 10 for each species), using the admixture model with correlated allele frequencies, which estimates fractions of individual genomes that belong to different ancestry groups. We also applied the same model adding sampling locations as prior (LOCPRIOR, [49]). The optimal K value was chosen by calculating the mean posterior probability for each K value, LnP(D), which is based on estimated maximum log-likelihood values [48]. We also calculated ΔK values (the rate of change in the log probability of data between successive K values) as suggested by Evanno et al. [53] using STRUCTURE-HARVESTER, version 0.6.94 [54]. Individual membership assignments estimated in STRUCTURE were averaged using CLUMPP, version 1.2.2 [55] with a Greedy algorithm and 10,000 random permutations. Furthermore, we implemented spatially-explicit Bayesian clustering using GENELAND, version 4.0.3 [52] in R, version 3.1.3 [39]. We used the spatially explicit model (1,000,000 MCMC iterations, thinning = 100), and set the number of potential populations to K equals 1 to 10. We applied the uncorrelated frequency model, which uses the Dichirilet distribution [48] to model allele frequencies as recommended by Guillot et al. [52]. The uncertainty associated with spatial coordinates was set to 20 m. Next, we performed ten additional runs for the selected K value, using the same parameters. To verify the consistency of this analysis, we repeated the analysis three times total and compared the results with each other. After finding the optimal K value with STRUCTURE and GENELAND, individuals were assigned to distinct genetic clusters using the percentage of genotypes' ancestry (Q; scores 70% indicate assignment to a genetic cluster and scores < 70% admixture).
To characterize countrywide spatial genetic structure, we also examined isolation by distance (IBD) and spatial autocorrelation patterns. First, we determined whether a significant correlation existed between pairwise codominant genotypic and geographical distances by applying simple Mantel tests for each species and sex with 10,000 permutations using ecodist, version 1.2.9 [56] in R, version 3.1.3 [39]. Second, spatial autocorrelation analysis was conducted in GenAlEx, version 6.41 [37] to examine the spatial extent of genetic structure and to determine if dispersal patterns were sex-biased. We correlated pairwise geographical and genetic distance matrices for each sex and species and generated autocorrelation coefficients (r) similar to Moran's I coefficient for each spatial distance class [57], which were visualized as correlograms. The distance classes ranging from 5 to 250 km were chosen based on the distribution of geographic distances among individuals. Spatial genetic structure is tested against the null hypothesis of no autocorrelation (r = 0) by generating 95% confidence intervals (CI) for each distance class via permutation (9,999 simulations) and bootstrapping (999 repeats). Within a given correlogram, significant spatial autocorrelation was confirmed only when a positive r-value fell outside the 95% CI (derived from the permutation test), and when the 95% CI about r (derived from bootstrapping) did not intercept the axis of r = 0 as described by Peakall, Ruibal [57]. The extent of non-random spatial genetic structure was based on the location of the first x-intercept (e.g., [57]).
Furthermore, to indirectly infer patterns of contemporary gene flow, we calculated maximum-likelihood estimates of pairwise relatedness coefficients (r) ranging from 0 (unrelated) to 0.5 (parent-offspring, full-siblings) among individual felids for all three target species within and among different sites in ML-RELATE [58]. Statistical differences between groups were evaluated using non-parametric Wilcoxon and Kruskal-Wallis rank-sum tests in R, version 3.1.3 [39].

Detection of migrants and contemporary migration rates
To assess levels of contemporary migration (i.e. within the last few generations) per species and among the five different geographical regions within Belize, we used two different approaches. First, we assessed contemporary dispersal by determining first-generation migrants in GENE-CLASS, version 2.0, which applies MCMC resampling algorithms to compute individual probabilities of being a resident (i.e., not a first-generation migrant) to each reference population [59,60]. The analysis is based on maximum likelihood estimation (Lhome:Lmax, the ratio of L_home to the highest likelihood value among all available population samples including the population where the individual was sampled [L_max]) [60]. We applied the Bayesian criterion of Rannala and Mountain [61] in combination with the MCMC resampling method of Paetkau, Slade [60]. We simulated 10,000 individuals and selected an alpha level of 0.01. Second, we estimated contemporary directional migration rates (i.e. migration events that occurred within the last few generations) using Bayesian inference framework implemented in BayesAss +, version 3.0 [62]. We conducted the analysis using 3,000,000 iterations, 1,000,000 iterations burn-in, and a sampling frequency of 2000. We used the default delta value of 0.15 for allele frequencies, migration rate, and inbreeding. We conducted ten independent runs of the analysis to confirm the consistency of migration rate estimates. Bidirectional migration rates were visualized as circos plots using circlize, version 0.3.4 [63] in R, version 3.1.3 [39]. Pumas had significantly lower allelic richness estimates at the northern site compared with most other regions (Kruskal-Wallis rank-sum test, A R , H = 11.10, P = 0.03), but did not differ for other diversity estimates (Kruskal-Wallis rank-sum test, H e , H = 8.10, P = 0.09; A P, H = 4.58, P = 0.33). F IS values for all target species and geographical regions were positive and at low to moderate levels, ranging from 0.00 to 0.08 for jaguars, 0.01 to 0.12 for pumas, and 0.04 to 0.14 for ocelots (Table 2). After sequential Bonferroni corrections, loci FCA043 (P < 0.000) and F98 (P = 0.011) for pumas, and loci FCA391 (P = 0.001), FCA275 (P < 0.000), and FCA741 (P = 0.011) for ocelots significantly deviated from Hardy-Weinberg expectations. Significant linkage disequilibrium after sequential Bonferroni correction (P 5.00E-04) was only detected among one pair of loci (FCA096 and FCA441) in jaguars.

Fine-scale genetic structure and contemporary gene flow
Pairwise F ST values between most geographical regions were significantly different from zero for all three species (Table 3). Our results suggested low to moderate genetic differentiation for North-Central 16 13  Table 2. Summary statistics of genetic diversity for Panthera onca, Puma concolor, and Leopardus pardalis in Belize.

Region
Panthera onca (n = 65) Puma concolor (n = 54) Leopardus pardalis (n = 30)  (Table 3). The AMOVA results for jaguars, pumas, and ocelots between genetic clusters defined by geographical regions within Belize are given in Table 4. All three comparisons indicated low but significant genetic differentiation. The AMOVA test showed highest levels of genetic variation within sampling sites and regions (> 85%), indicating that genetic differentiation among sampling sites and regions is low for all three target species. DAPC analysis suggested that jaguars within Belize form two groups with low levels of genetic differentiation as north and north-central jaguars were slightly separated from the remaining three regions (Fig 2A). For pumas, DAPC grouped individuals in three groups corresponding to the north, north-central, and central/central-south/south geographical regions (Fig 2B). Non-spatial Bayesian clustering analysis in STRUCTURE indicated no genetic structure for jaguars with K = 1 as the most probable number of genetic clusters using the admixture model with correlated allele frequencies (Fig 3A, S1A Fig). However, when prior location information of the origin of the scat samples was included (LOCPRIOR model), the number of genetically distinguishable groups detected by STRUCTURE was K = 2, roughly grouping sites into a northern and southern cluster consistent with DAPC analyses (Fig 3B, S1B Fig). Jaguars detected within the north/north-central regions had a high membership assignment to the northern genetic cluster (mean Q = 0.98, range 0.92-0.99). More than half of the jaguars (56%)  For pumas, we identified one genetic cluster using the admixture model with correlated allele frequencies (Fig 4A, S2A Fig) and two genetic clusters using the LOCPRIOR model in STRUCTURE, separating most individuals (77%) detected within the northern sampling sites  Using spatiallyexplicit Bayesian clustering in GENELAND, we identified four genetic clusters and one ghost population for pumas, corresponding to our main a priori geographical regions within Belize, including the north (100% residents with mean Q = 1), north-central (79% residents with mean Q = 1; four migrants with mean Q = 0.02), central and central-south (71% residents with mean Q = 1; two migrants with mean Q = 0.07; three admixed genotypes) and south (80% residents with mean Q = 0.98; one migrant with Q = 0.24) (Fig 4C). For ocelots, Bayesian clustering analysis in STRUCTURE (admixture and LOCPRIOR models) did not reveal any population subdivision (K = 1) (Fig 5A and 5B, S3 Fig). Spatially-explicit Bayesian clustering in GENELAND revealed seven genetic clusters, through which individuals were roughly grouped into two main genetic clusters (northern and southern Belize) (Fig 5C). While the northern cluster showed a high degree of admixture, the southern cluster consisted of several individuals detected at Cockscomb Basin Wildlife Sanctuary that were strongly assigned (Q > 0.95) to this site.
Examination of IBD patterns in jaguars showed no significant correlation (all jaguars, r = 0.010, P = 0.510, Fig 6A; males, r = 0.099, P = 0.100; females, r = -0.201, P = 0.200) between-individual genetic and geographic distances, verifying that geographic distance is not driving genetic differentiation at this spatial scale. Significantly positive relationships between genetic and geographic distance were detected for pumas (Fig 6B), which were most pronounced in females (all pumas, r = 0.266, P = 0.010; males, r = 0.116, P = 0.050; females, r = 0.379, P = 0.010). This suggested that gene flow in female pumas is spatially limited and may be one of the factors driving genetic differentiation in pumas in Belize. In ocelots, we found a significant, but weak signature of IBD (all ocelots, r = 0.134, P = 0.020, Fig 6C; males, r = -0.029, P = 0.590; females, r = 0.131, P = 0.070), indicating that geographic distance potentially has a small effect on gene flow. However, sample sizes in female jaguars and male ocelots were low (> 10 individuals), thus careful interpretation of results is essential.
To examine the spatial extent of genetic structure, we conducted spatial autocorrelation analysis in all three target species and by sex. Female jaguars and male ocelots were not analyzed as individual groups due to low sample sizes. The autocorrelogram for all jaguars showed significantly positive autocorrelation in the first three distance classes (5 km, r = 0.050, P = 0.002; 10 km, r = 0.053, P = 0.001; 15 km, r = 0.042, P = 0.012) and an x-intercept of r at 20 km (S4A Fig), empirically confirming the presence of nonrandom spatial structuring and genetic association among individuals at distances < 20 km. Individuals below this threshold, share a higher proportion of genes, than spatially distant individuals. Patterns of spatial  km for all samples and at~83 km for females, indicating an overall larger spatial extent of positive autocorrelation in ocelots than in the two larger felids.
Relatedness among sites was significantly different (Kruskal-Wallis rank-sum test, H = 38.17, P < 0.000) among the three target species, with highest mean relatedness coefficients for jaguars (0.07 ± 0.11, SD), followed by pumas (0.05 ± 0.09, SD), and ocelots (0.03 ± 0.07, SD). Mean relatedness coefficients within sites were slightly higher for jaguars (0.11 ± 0.15, SD), followed by pumas (0.10 ± 0.14, SD), and ocelots (0.08 ± 0.13, SD). However, mean relatedness coefficients within sites were not significantly different among species (Kruskal-Wallis rank-sum test, H = 3.38, P = 0.184). In all three species, within-site relatedness was significantly higher than among-site relatedness (jaguars [Wilcoxon rank-sum test, W = 24726, P = 0.024], pumas [Wilcoxon rank-sum test, W = 8127, P < 0.000], ocelots [Wilcoxon ranksum test, W = 2180, P < 0.000]). However, a skewed sex-ratio towards males in jaguars and females in ocelots may introduce bias into the analysis. Using GENECLASS, we identified three first-generation migrants (i.e. dispersers) for jaguars. The three animals were located, and originally came from, respectively: (1) the central from the south, (2) the central-south from the south, and (3) the south from the north sites. For pumas, we detected a total of four first-generation migrants. The animals were located, and originally came from respectively: (1) north from the north-central, (2) north-central from the north, (3) central from the north, and (4) south from the south-central. For ocelots, we identified two first-generation migrants. One ocelot was located in the northern cluster (Fireburn Nature Reserve) and genetically assigned to the southern cluster. The other ocelot originally came from the southern cluster (Cockscomb Basin Wildlife Sanctuary) and was signed to the northern cluster (Table 5). All first-generation migrants detected in this study were males.
Contemporary and bidirectional migration rates estimated between all sampling sites in BayesAss+ suggested an average migration rate of 0.07 (range 0.02-0.18) for jaguars across Belize ( Fig 7A, Table 6). Migration rates, which represent the proportion of individuals that move from corresponding source population each generation, were highest between the central and south-central sites (0.18). Migration rates for pumas were similar (mean = 0.06, range 0.01-0.12; Fig 7B) to jaguars across Belize, with the highest estimates between central and north, south-central and central, and south and north-central sites (0.12). Mean migration rates for ocelots were relatively high (0.15) and two-to-three times higher from the southern to  [37] versus pairwise geographic distances (km). Statistical significance was assessed using simple Mantel tests in ecodist, version 1.2.9 [56] in R, version 3.1.3 [39]. Each point represents a pairwise comparison among individual felids.

Discussion
Preserving and restoring genetic connectivity in anthropogenically altered and fragmented landscapes is a key aspect of conservation (e.g., [5,64]). Comparative conservation genetics studies increase our understanding of species-specific responses to habitat loss and fragmentation, and have been invaluable for developing conservation and management strategies for multiple species of conservation concern (e.g., [65,[66][67][68][69][70]). The main objective of our study was to use noninvasive genetic sampling to examine levels of genetic diversity, population structure, and contemporary gene flow of multiple elusive felid species living in fragmented tropical forest habitats. Ultimately, we aim to support current and future conservation and management efforts for wild felids within Belize and throughout their range. Our results provide the first comparative genetic assessment of wild populations of jaguars, pumas, and ocelots within Mesoamerica. The use of scat detector dogs to increase detection rates of fecal samples and the application of standardized methods for sample collection and storage of fecal DNA [35] have been instrumental in successfully applying these methods on multiple felid species in tropical environments, where relatively little is known about genetic diversity and structure of these species.

Genetic diversity of three Neotropical felids
Genetic diversity across Belize was moderate for all species and lowest for jaguars, followed by pumas and ocelots (Table 2), which likely reflects differences in effective population sizes for In comparison with other genetic studies, using Fca microsatellite markers [71], levels of genetic diversity for Belizean jaguars were relatively low. For example, Brazilian (H e = 0.73) [25] and Colombian jaguars (H e = 0.85) [23] were genetically more diverse. But genetic diversity estimates for Central American jaguars obtained through a range-wide genetic study [24] were concordant with our findings. For pumas, diversity estimates in Belize (H e = 0.65) were similar [72] or higher [73,74] compared to North American pumas and lower than South American pumas (21,27) using Fca microsatellite markers [71]. Culver, Johnson [75] observed similar levels of genetic diversity for pumas in Central America. Belizean ocelots had lower levels of genetic diversity compared to South American populations (e.g., [21,76]), but showed higher estimates than ocelots studied in northeastern Mexico and southern parts of the United States [77,78]. The differences in genetic diversity estimates for Central and South American felids can be explained by a combination of historical and contemporary factors, including range contraction/expansion, and changes in population size and gene flow (e.g., [24,75,77,79]). Nonetheless, a strict direct comparison of diversity measures between these studies should be viewed cautiously since different sets of Fca microsatellite markers were used by most studies.

Fine-scale genetic structure and contemporary gene flow
Using a combination of individual-and population-based genetic structure analyses revealed low to moderate levels of genetic differentiation for Neotropical felids within Belize, with jaguars showing the lowest levels of genetic subdivision (i.e. highest genetic connectivity) across the country, followed by ocelots and pumas. Considering the dispersal capabilities, particularly of the two larger felids (e.g., [9,33]), the presence of large tracts of forests, and the relatively small total area of Belize (~22,966 km 2 ; Fig 1), we were not surprised to see relatively high levels of gene flow and low genetic differentiation. However, our fine scale genetic analyses revealed that jaguars in Belize, may not exist in total panmixia, but rather exhibit low levels of genetic differentiation when genetic divergence was estimated indirectly (F statistics, AMOVA, DAPC) or through non-spatial Bayesian clustering methods. Spatial autocorrelation analysis detected genetic structuring in jaguars at fine geographic scale (< 20 km) across all samples and in males (females were not analyzed separately due to small sample size). This indicated high genetic associations among individuals sampled in close proximity, which could be caused by several factors including natal philopatry, habitat fragmentation, and/or clustered and highdensity sampling efforts (e.g., [74,80]). Since male philopatry in jaguars is unlikely, spatial autocorrelation detected at close distance classes was either a sampling artifact or caused by limited gene flow resulting from restricted movement of individuals across the landscape as described in other studies (e.g., [81,82]). When using spatially informed Bayesian clustering, we detected moderate levels of genetic differentiation. Spatial Bayesian clustering methods, which integrate geographically explicit prior distributions into the inference process, have been shown to be more powerful in detecting weak genetic discontinuities in small data sets (e.g., [66,72,83]). However, careful consideration is necessary when interpreting the Bayesian clustering results, since opportunistic sampling, as well as sampling density, may cause inconsistencies in the clustering analysis [77,78]. Multiple methods roughly grouped jaguars into a northern (north and north-central sites) and a southern (central, central-south, south sites) cluster, which are separated by a~50 km wide stretch of open and human-dominated landscape bordering north-central and central Belize. Our results suggested that the genetic subdivision detected corresponded largely to patterns of habitat fragmentation and human disturbance. This finding was concordant with another genetic study on jaguars showing that human-dominated landscapes have the potential to limit gene flow in jaguars on a fine geographic scale [25]. Although jaguars are considered generalist species and utilize a variety of habitat types (e.g., [1,9,84,85]), they prefer areas with tree cover [86], close to water bodies [11], and are less likely to use human-modified landscapes such as agricultural [87] and more populated areas [88]. Former studies also showed that physical barriers such as roads limit jaguar movement, with females exhibiting a high degree of road avoidance [88,89]. However, spatial Bayesian clustering analysis also identified several admixed individuals (Fig 3C), which indicates that jaguars did not only move across the region, but were also able to reproduce. In addition, our contemporary gene flow analysis found evidence for recent dispersal events by classifying three male jaguars as first-generation migrants. Two of the jaguars were detected in central (Chiquibul Forest Reserve and National Park) and south-central (Cockscomb Basin Wildlife Sanctuary) Belize, and were genetically assigned to southern Belize. One male jaguar that was originally sampled in the Golden Stream Corridor Preserve in southern Belize, was genetically assigned to the most northern site, which is located approximately 250 km away. These findings were supported by contemporary migration rates, which suggest that movement of jaguars between some sites within Belize is relatively common, although asymmetric migration rates were detected in some cases (Fig 7A, Table 6). Interestingly, we detected asymmetric migration rates within central Belize (Chiquibul Forest Reserve and National Park, Mountain Pine Ridge Forest Reserve) to Cockscomb Basin Wildlife Sanctuary, suggesting that about 18% of the jaguars moved out from these areas during this generation, whereas immigration rates into these areas were not higher than 7% (Fig 7A). Central Belize, including the Chiquibul Forest Reserve and other areas close to the Guatemalan border (e.g., Columbia River Forest Reserve), have experienced high rates of incursions by illegal Guatemalan immigrants during the last forty years (e.g., [31,90,91]). Illegal activities within these protected areas include farming, hunting, harvesting of timber and non-timber products, and drug trafficking. In addition, industrial activities such as agroindustry, hydroelectricity, mining, logging, and road expansions into more remote areas have increased the rates of forest loss and human impact [90]. The elevated levels of human disturbance and land use change within these protected forest habitats could potentially cause emigration rates to increase for jaguars. In contrast, outward migration rates from Cockscomb Basin Wildlife Sanctuary, which represents undisturbed forest habitat and is also part of the Maya Mountain Massif, a large block of tropical forest in Belize, are relatively low (5% to south, 3% to north, 7% to central, and 7% to northcentral Belize). Besides human-induced impacts (e.g., human developments, hunting, humanwildlife conflict), migration rates may be also driven by ecological (e.g., population density, competition, social system, prey density), and/or behavioral (e.g., dispersal, habitat use) factors as observed in several other carnivore studies (e.g., [67,72,75,[92][93][94][95]). Jaguars sampled at the Rio Bravo Conservation Management Area, showed relatively high outward migration rates towards northern and central Belize (Fig 7A). This could be explained by Rio Bravo Conservation Management Area being part of La Selva Maya, which likely represents one of the most important source populations for jaguars within the region, potentially exhibiting dispersal induced by positive density-dependence. In conclusion, our results suggested that jaguars within Belize are genetically still relatively well connected; nonetheless, there is fine-scale genetic differentiation occurring and if habitat loss and fragmentation continue at the current rate, jaguar populations may consequently exhibit much stronger genetic structure, as shown in a previous fine-scale study [25]. For pumas, which are found in a wide range of natural environments (e.g., [1,9,96]), we detected moderate levels of genetic differentiation within Belize with sub-structuring between most sites surveyed. Spatial, individual-based Bayesian clustering analysis showed that genetic clusters are more distinct and less admixed in areas experiencing higher levels of human disturbance. Genetic subdivision was most pronounced between Fireburn Nature Reserve, one of the geographically most isolated protected areas in northern Belize, and all other areas surveyed. These findings matched with previous fine-scale genetic studies of pumas (e.g., [1,9]), which concluded that despite the pumas' ability to disperse over long distances [72,[97][98][99][100], genetic substructure increases in areas where habitat is less contiguous and is altered due to anthropogenic impacts (e.g., [33,101]), or is negatively influenced by habitat barriers such as rivers [97], open deserts and grasslands [98]. However, the detection of admixed individuals and dispersers suggested that pumas were able to move across Belize and reproduce to some extent. Protected areas in central (Chiquibul Forest Reserve and National Park, Mountain Pine Ridge Forest Reserve) and south-central (Cockscomb Basin Wildlife Sanctuary) Belize formed one genetic cluster with several admixed individuals and dispersers from northern Belize (Fig 4C). Furthermore, our study also showed evidence for fine-scale genetic structure in female pumas, inferring that females living in close geographic proximity are on average genetically closely related, exhibiting short dispersal distances and female philopatry as described by other puma studies [102,103] and in mammals in general (e.g., [104,105]). In contrast, and different from jaguars, male pumas did not show evidence of positive spatial autocorrelation at any geographic distance tested, suggesting that dispersal in pumas across this landscape is male-biased, which is consistent with previous studies [106,107]. However, the spatial extent of positive autocorrelation described by our study was smaller in comparison to North American pumas (e.g., [74,106,108]), suggesting that dispersal movements in Neotropical pumas, especially in females may be more restricted. Furthermore, contemporary gene flow analysis found evidence for four male first-generation migrants (i.e. dispersers) in Belizean pumas. We detected two dispersers between north-central and the most northern site. One puma was detected in the Mountain Pine Ridge Forest Reserve in central Belize and genetically assigned to the most northern sites within Belize. The fourth migrant was detected in the Golden Stream Corridor Preserve in southern Belize, and originated from the neighboring Cockscomb Wildlife Basin Sanctuary in south-central Belize. Our results suggested that most dispersal events were directional and that dispersers moved away from either the most northern or the most southern sites, which both overlap with areas experiencing high levels of deforestation and land conversion (e.g., [31,109]). Contemporary migration rate analysis, which suggested that on average about 6% of Belizean pumas move among different geographical areas, complemented these findings. An interesting observation was that outward migration rates for pumas from the Cockscomb Basin Wildlife Sanctuary moving into central Belize (Chiquibul Forest Reserve and National Park, Mountain Pine Ridge Forest Reserve) were more than double in comparison to their inward migration rates, which stands in contrast to our findings for jaguars. Despite the fact that both areas are separated by the Maya Mountain range, we believe that high population densities for jaguars and corresponding levels of interspecific competition at Cockscomb Basin Wildlife Sanctuary may cause pumas to have a two to three times higher outward migration rate in comparison to levels of incoming pumas. Furthermore, Rio Bravo Conservation Management Area showed moderate puma outward migration rates (9%) into central Belize, whereas inward migration was only one third of it. Similar to jaguars, Rio Bravo Conservation Management Area may also represent an important source population for pumas within the region. Migration rates into the most southern sites were generally limited ( 3%), which most likely was caused by higher levels of land conversion and human disturbance in this region. Despite the detection of several first-generation migrants within the country, we conclude that pumas exhibited moderate levels of genetic differentiation within Belize, which may be primarily driven by spatial autocorrelation patterns in female pumas.
For ocelots, non-spatial Bayesian clustering analysis rendered support for no or low levels of genetic differentiation. Furthermore, autocorrelation analysis suggested a large spatial extent of genetic structure (~84 km) in a mesocarnivore, indicating that gene flow for ocelots across this landscape may be relatively high. Given ocelots' smaller body size and their tolerance to human activity to some extent (e.g., [13]) they may be more successful moving through human-dominated landscapes than the two larger felids. However, spatial autocorrelation analysis also revealed fine-scale genetic structure over short distances (< 5 km) in female ocelots, suggesting female philopatry. In addition, spatial Bayesian clustering resulted in moderate genetic subdivision, roughly grouping Belizean ocelots into a northern and southern cluster, indicating that human-dominated areas bordering north-central and central Belize presumably restrict gene flow between these sites, but still allow for occasional dispersal events. From previous studies it is known that ocelots are habitat specialists preferring closed forest types and dense thorn scrubs (e.g., [9,96,110,111]), thus movement across central Belize may be limited by open areas (e.g., savannah, agricultural areas). Nonetheless, among the northern and southern cluster, we documented movement into and out of both areas. First, we observed two male first-generation migrants in ocelots going from the north to the south and one from south to north (Table 5). Second, migration rates between northern and southern sites were also found to be directional, with higher levels of gene flow going from the north to the south (21% of individuals for this generation), than from southern to northern regions (9%) (Fig 7; Table 6). Although ocelots are still considered widely abundant across Central America, they are understudied, and we found that the effects of habitat loss and fragmentation may have the potential to severely reduce genetic connectivity, which has been documented in ocelots at the northern extent of their range (e.g., [78]) and other small felid populations (e.g., [112,113]).

Comparative genetic structure of two sympatric felids
This study directly compared patterns of fine-scale genetic connectivity in multiple Neotropical felids. For jaguars and pumas, which differ in size across regions of coexistence, jaguars are the larger, dominating species [114]. We observed regional differences in genetic connectivity levels within Belize, with a tendency for pumas to be less connected than jaguars. We did not include ocelots in this section due to the smaller sample size obtained for this species. Our hypothesis that pumas would experience higher levels of genetic connectivity relative to sympatric jaguars was not supported. Varying levels of genetic differentiation for jaguars and pumas were especially evident between the most northern and most southern sites within Belize, where felids are likely confined to small and fragmented forest patches surrounded by less suitable habitat and human-modified areas of higher disturbance (e.g., human-wildlife conflict, farming). This pattern of genetic divergence was more pronounced in pumas. However, we also detected evidence for first-generation migrants in both species (3 dispersers in jaguars, 4 dispersers in pumas). Considering that this analysis depends on sample coverage and several other factors, the results only represent an approximation of the actual dispersal events. We found a general trend that pumas were more likely to disperse to geographically close or neighboring areas in comparison to jaguars. Comparative relatedness analysis supported the findings related to genetic divergence by revealing that jaguars had significantly higher pairwise relatedness values among sites than pumas, and showed a generally higher tendency to disperse from their natal areas. Autocorrelation analysis revealed significant genetic structure in jaguars and female pumas at close distance classes with positive correlations almost three times as high in female pumas. Consequently, we reject the null hypothesis of a random distribution of genotypes at this spatial scale. The spatial structuring we detected can be associated with a multitude of species-specific processes, including social and mating systems, sex-biased dispersal, and/or restricted movement and consequently gene flow (e.g., [57,67]). Although it is difficult to differentiate between the effects of social organization and reduced gene flow, differences in autocorrelation patterns between sexes in pumas indicated that female philopatry and male-biased dispersal may be primarily driving the spatial structure detected, which is common in felids and mammals in general (e.g., [9,104]). Alternatively, for jaguars, positive spatial structure at close geographic distances, especially for males, may be the result of restricted gene flow across the landscape, which stands in contrast to our findings for male pumas. In addition, the spatial extent of genetic structure detected in jaguars and pumas (with the exception of male pumas) also suggested that effective dispersal distances across this heterogeneous landscape are potentially limited.
Gene flow may be also driven by various other ecological factors, including the degree to which habitat use and human disturbance limit movement in wide-ranging species such as jaguars and pumas. Former studies described pumas as opportunistic in their habitat use since they are known to use a wider variety of habitats in comparison to jaguars (e.g., [9,96]). In contrast to jaguars, which have a tendency to avoid open areas and prefer dense forest habitats, previous studies also reported movement of pumas through disturbed and human-developed areas (e.g., [25,72,97,98,115]). Opposite to these findings and in agreement with our study, research on habitat use of jaguars and pumas in Belize concluded that pumas were less likely to be found outside of protected areas than jaguars due to a differential tolerance to human disturbance and resource limitation [116]. Davis, Kelly [117] also reported that pumas, and jaguars to a lesser extent, were sensitive to human disturbance, even within protected areas of Belize. In addition, interspecific spatial interactions between sympatric jaguars and pumas resulting in spatial avoidance of jaguars by pumas as described by a few studies (e.g., [118,119]), may also affect movement and consequently gene flow patterns of these two sympatric species across the landscape.

Conservation and management implications
With increasing anthropogenic impact and landscape change, it is crucial to genetically monitor wild felid populations. Understanding levels of fine-scale gene flow and genetic structure patterns of multiple co-occurring felid species reveal behavioral and other factors that influence gene flow across the landscape, which is vital for planning and prioritizing future conservation and management efforts on a countrywide scale. This study provides comprehensive baseline conservation genetic data and demonstrates that noninvasive genetic sampling is an efficient research approach to simultaneously assess levels of genetic diversity and differentiation of multiple Neotropical felid species in the wild. We found that genetic diversity for wild felids in Belize is moderate and that levels of genetic connectivity within the country are moderate to high. Although Belize has a high proportion of forest cover (~62.7%) and protected areas (3 6%) compared to other Mesoamerican countries, we believe that levels of genetic connectivity are likely to decrease if habitat loss and fragmentation continue at the current rate. Despite the dispersal capabilities of the two large felids, our study detected evidence for fine-scale genetic differentiation particularly between northern and southern sites, indicating that the moredeveloped and human-dominated areas adjacent to north-central and central Belize may be subtlety-affecting movements of felids. Consequently, we recommend prioritizing countrywide conservation and management efforts for wild felids, such as maintaining and enhancing biological connectivity among the national protected area system within Belize through strengthening formerly identified wildlife movement corridors (Northern Belize Corridor, Central Belize Corridor, Southern Belize Corridor) [120,121]. The fledgling Central Belize Corridor, which connects Belize's two largest forest blocks (north: RBCMA, Yalbac, Laguna Seca, Gallon Jug; south: Maya Mountain Massif), represents a critical link in the Mesoamerican Biological Corridor within this region, and deserves special attention. In addition, we also recommend conducting continued genetic monitoring efforts, and assessing functional connectivity of movement corridors through a landscape genetics approach. Furthermore, we also encourage large-scale research efforts focusing on multiple, sympatric species to increase understanding of species interactions and responses to fragmented landscapes and to develop regional multispecies conservation and management strategies. With Mesoamerica bearing one of the highest deforestation rates in the world, and lacking in genetic studies of wild felids in general, geneticbased monitoring focusing on wildlife of conservation concern is needed to assess connectivity and genetic health across the entire region.