The Meso-American slider turtle (Trachemys venusta) is a freshwater turtle that is widely distributed from Mexico to Colombia. Due to the overexploitation of populations of this species in Mexico, it has been placed within the “subject to special protection” category formulated by the Official Mexican Standard NOM-059-ECOL-2010. In the state of Tabasco, Mexico, Management Units for the Conservation of Wildlife (UMA) were created to reduce the impact of overexploitation of freshwater turtles bred in captivity. However, no genetic management plan was considered. The present study was carried out in an UMA in the state of Tabasco. We obtained the level of genetic diversity of the founder individuals of the UMA in order to develop a management plan which will optimize reproduction in the UMA. Genetic diversity was compared between captive (n = 86) and wild (n = 45) individuals using 14 microsatellite molecular markers. The genetic diversity parameter determined in this study was slightly higher for captive than for wild population (He = 0.606 and He = 0.594 respectively), reflecting the mix of genetic sources in captive group (founding individuals from different localities) and demonstrating that the captive population contains a diverse subset of alleles from representative populations. The analysis of genetic structure revealed a relationship between captive and wild populations, indicating the influence of the two principal river basins in this region on the populations structure of freshwater turtles. Finally, according to the results obtained from the relationship analysis, we recommend the use of 19 females and 13 males to constitute the appropriate breeding group, generating a potential of 247 dyads with no relationship. However, in order to improve breeding program and the genetic diversity of captive population, we suggest to introduce wild-caught individuals. These results are the first regarding genetic management in a Mexican UMA and demonstrate the importance of molecular approaches in the management and conservation of captive species.
Currently, the world faces a rapid loss of its biodiversity (Valiente-Banuet et al., 2015), predomi nantly due to anthropogenic activities that have resulted in climate change, habitat loss, the introduction of invasive species, pollution, habitat degradation among others; consequently, captive breeding programs are a priority for species conservation (Miller et al., 2018). During the last few decades there has been an increase in global awareness on the need to conserve species through the management of captive populations and develop integrative conservation and management strategies (Ballou et al., 2010; Ortega-Argueta et al., 2016). Accordingly, the World Association for Zoos and Aquariums promotes an expansion of animal management programs at the population level to establish viable populations while recognizing that this represents a biological and organizational challenge (Ballou et al., 2010). In addition, other initiatives are emerging such as the Community-Based Natural Resources Management approach that began in Africa and promotes the integration of complex relationships among society (e.g., local livelihoods), economic and political authorities, the environment, and sustainability principles for the management of natural resources (Ortega-Argueta et al., 2016). In Mexico, this type of environmental policy initiative commenced in 1997 with the creation of the National Program for Wildlife Conservation and Productive Diversification of the Rural Sector. Part of this program included the development of the System of Management Units for Wildlife Conservation (SUMA, by its Spanish acronym) created by the Ministry of Environment and Natural Resources (SEMARNAT, by its Spanish acronym) (Ortega-Argueta et al., 2016; Weber et al., 2006). Through this program, many Wildlife Management Units (UMAs, by its Spanish acronym) have been created (Valdez et al., 2006; Weber et al., 2006). The goal of UMAs is the appropriate and economically viable management of wildlife resources (fauna and flora) for conservation, rescue, and preservation purposes (Ministry of Environment and Natural Resources [SEMARNAT], 1997). UMAs incorporate a wide range of activities such as research, recreation, environmental education, game farms, and commercialization of wildlife with products that are subject to regulated laws (Valdez et al., 2006). Furthermore, these units of management create biological corridors by interconnecting protected natural areas and zones that implement the sustainable exploitation of species, thus contributing to biodiversity conservation (de Benito, 2009). Unfortunately, genetic aspects linked to captive breeding are rarely considered in the management of an UMA and consequently may lead to the failure of such programs.
One of the principal objectives of captive breeding programs is to obtain and maintain a population with a high level of genetic variation (Ballou & Lacy, 1995). Generally, captive populations are small and established with a few founder individuals, making these populations prone to significant genetic changes such as a loss of genetic diversity, inbreeding depression, accumulation of novel deleterious mutations, and genetic adaptation to captivity (Williams & Hoffman, 2009). However, when considering the small amount of research conducted in UMAs, we have observed that genetic considerations (for example, the genetic diversity of founders and kinship relations) are one of the least considered aspects when establishing and managing an UMA. In Mexico, studies that have implemented genetic tools to improve the management of captive populations within UMAs are scarce. Zarza et al. (2016) identified the origin of 24 individuals of the threatened, endemic species Ctenosaura pectinata Wiegmann 1834 (Squamata, Iguanidae), the Mexican spiny-tailed iguana, to two Mexican UMAs. The study was conducted using mtDNA and microsatellite molecular markers in addition to genetic assignment methods, with the aim of comparing these individuals with a database of 341 individuals from 49 localities. This type of research demonstrates the importance of genetic tools in identifying the origin of individuals before a possible release in the wild.
Slider turtles of the genus Trachemys Agassiz, 1857 (Testudines, Emydidae) are freshwater turtles distributed throughout the New World with 26 recognized extant forms (Seidel, 2002). In Mexico a total of 13 subspecies are recognized, belonging to between one and nine species depending on the author concerned and the species concept used to identify lineages (review in Parham et al., 2015). One species from the Trachemys genus is T. venusta (Gray, 1856) that presents a range that extends from Mexico to Colombia (Ceballos & Brand, 2014); the most recent genetic analysis of this species identified populations present from the Atlantic river basins of Mexico (and an isolated population in Acapulco) to Central America (Parham et al., 2015). Its nomenclature is still uncertain; therefore, genetic studies have been conducted in order to obtain a clearer understanding of its classification (Fritz et al., 2012; Parham et al., 2013, 2015; Vargas-Ramírez et al., 2017). Three subspecies are described for T. venusta (Ernst & Seidel, 2006), all found in Mexico: T. venusta venusta (Gray, 1856), T. venusta cataspila (Günther, 1885), and T. venusta grayi. Based on genetic analysis, the latter subspecies was recently proposed as a full species T. grayi (Bocourt, 1868) by Parham et al. (2013). In the state of Tabasco, Mexico, only one subspecies was reported, T. venusta venusta (Ernst & Seidel, 2006; Parham et al., 2015), commonly called hicotea. Although Vargas-Ramírez et al. (2017) suggested two lineages for T. venusta venusta, the origin of specimens used is not evident and the authors did not explain each lineage; consequently, we decided to maintain the general and confirmed nomenclature of T. venusta venusta for the subspecies in Tabasco until more elements confirmed the existence of two different lineages for this subspecies.
The hicotea has played an important role in the culture and gastronomy of the Maya people since the prehispanic period (Guevara Chumacero et al., 2016). One of the few studies that deals with the exploitation of turtles in the state of Tabasco shows that the ingrained consumption of turtles persists despite the fact that the species is considered endangered (Guevara Chumacero et al., 2016). Turtles are still in demand for commercial and consumption purposes. In addition to this culinary practice, urban growth has increased markedly since the mid-20th century leading to high rates of deforestation (Palomeque de la Cruz et al., 2017; Tudela, 1992). The city of Villahermosa, the largest urban area in the state of Tabasco, is a perfect example of such environmental degradation; over a period of only 40 years, 4008 ha of forest vegetation and 289 ha of wetlands have been lost as the city grew rapidly (Palomeque de la Cruz et al., 2017). The uncontrolled exploitation of freshwater turtles in Mexico, particularly in Tabasco, together with habitat loss has resulted in many species of turtles being placed in the endangered category in the Official Mexican Standard (Norma Oficial Mexicana), the NOM-059-SEMARNAT-2010. This regulation establishes the protection of the environment and native species of wild flora and fauna of Mexico considering the following categories: endangered, threatened, under special protection, and probably extinct in the wild. Trachemys venusta (previously known as T. scripta venusta (Flores-Villela & Canseco-Márquez, 2004) is classified as “under special protection” (SEMARNAT, 2010).
Molecular tools allow us identify the origin of individuals when unknown (confiscated pets, illegal seizure) and to infer kin relationships among captive founder individuals (Miller et al., 2018). For the successful management of breeding programs, knowledge on the relatedness between individuals is fundamental for success, because it allows a minimization of mean kinship and retains the maximum of genetic variation (review in Miller et al., 2018). There are two different methods that could be used to evaluate kin relationships: (1) the relatedness (r) which is a measure of the fraction of identical alleles shared by offspring/descent among individuals, and (2) the relationship category which is a particular pedigree (genealogical) relationship such as full siblings or half siblings (Blouin, 2003).
In this study, genetic characteristics (e.g., genetic diversity, genetic structure, kinship analysis) of the founder individuals of T. venusta venusta from one UMA in the state of Tabasco were characterized and compared to wild individuals, with the purpose of establishing an efficient management plan for reproduction in the UMA. Considering that mating between unrelated individuals and the maintenance of sufficient genetic variability in individuals reared at the UMA could optimize the management program, our particular objectives were (1) to determine the genetic diversity level of founder individuals of Trachemys venusta venusta in a UMA and compare it with wild individuals, (2) to evaluate intra- and inter- genetic structure of the UMA and wild population, (3) to determine kin relationships (relatedness and pedigree) between males and females in the UMA, and finally (4) our results will be discussed from a management and conservation perspective.
Sampling Sites and Collection
The samples used for this study were obtained from one Wildlife Management Unit (UMA) and three different wild localities (Figure 1). All samples were collected with trammel nets, Fyke nets, or hand captures during 2017 and 2018 (Vogt, 1980, 2016). Adult (n = 32) and juvenile (n = 13) wild turtles were collected and classified as follows: juveniles = shell length <20 cm with bright green coloration; adults = shell length >20 cm with green, brown or olive coloration (Ernst & Seidel, 2006; Gibbons, 1990). The UMA is located in the state of Tabasco, Mexico (municipality of Nacajuca; 18°11’23’’N – 92°59’37’’W) and was created in 1978. It is dedicated to the reproduction of seven freshwater turtles Dermatemys mawii Gray, 1947 (Testudines, Dermatemydidae), Chelydra serpentine Linnaeus, 1758 (Testudines, Chelytridae), Staurotypus triporcatus Wiegmann, 1828 (Testudines, Kinosternidae), Trachemys venusta venusta, Rhinoclemys aerolata Duméril & Duméril, 1851 (Testudines, Geoemydidae), Claudius angustatus Cope, 1865 (Testudines, Kinosternidae), and Kinosternon leucostomum Duméril & Duméril, 1851 (Testudines, Kinosternidae) for conservation purposes and to exchange individuals with other turtle farms. All initial founder individuals were obtained from the UMA (n = 86; nfemale = 73 and nmale = 13) and were obtained when the UMA was created from Pomposu and El Espino (C. Zenteno, personal communication, May 4, 2020; Figure 1). In 2008, this UMA reported a total of 4125 individuals of the T. venusta species, including 1650 newborn individuals (SEMARNAT, 2009). In addition, we obtained a total of 45 wild individuals (nfemale = 32 and nmale = 13) from (1) La Venta Park Museum located in Tabasco state (municipality of Centro; 18°00’02’’N–92°56’08’’W, n = 29), (2) Miguel Hidalgo (municipality of Centla; 18°20’00’’N-92°30’00’’W, n = 4), and (3) Bosques de Saloya (municipality of Nacajuca,18°11’22.85’’ N–92°58’39.37’’ W, n = 12). These individuals were considered as wild as they were collected from a natural environment (eMiguel Hidalgo and Bosques de Saloya), or had been recently rescued from the wild and donated directly to a center/farm (e.g., La Venta Park). Fewer individuals than recommended (25–30 for analysis using microsatellites; Hale et al., 2012) were collected from the localities of Miguel Hidalgo and Bosques de Saloya; therefore, all 45 wild individuals were grouped and analyzed as a single wild population. Furthermore, although Miguel Hidalgo and Bosques de Saloya localities are distant from each other, they belong to a geographical landscape consisting mainly of lacustrine and palustrine wetlands (areas with lagoons of more than 100 ha) (Barba-Macías et al., 2006). This type of landscape is very important for freshwater turtles as it helps maintain the regional stability of turtle populations (Mali et al., 2016). The individuals from La Venta Park Museum originate from seizures carried out by the Federal Attorney of Environmental Protection of Mexico (PROFEPA, by its Spanish acronym) which confiscates illegally trafficked animals. These seizures have been carried out mainly in the municipality of Centla which includes the locality of Miguel Hidalgo (Federal Attorney of Environmental Protection of Mexico [PROFEPA], 2014). For the aforementioned reasons, we consider that the wild individuals in our study are representative of the genetic diversity of the species at a regional level in the state of Tabasco.
Tissue samples were collected from the interdigital webbing by making small skin cuts (1 cm2). Before sampling, we disinfected the skin with 70% ethanol and subsequently applied an antiseptic (methylthioninium chloride) on the wound to avoid infection. The collected samples were preserved in a tissue buffer solution (salt-saturated 20% DMSO; Proebstel et al., 1993), and stored at −80°C until DNA extraction.
DNA Extraction and Microsatellite Amplification
DNA was extracted using the DNeasy Blood and Tissue kit (Qiagen); quality and concentration were verified using a UV-Vis spectrophotometer (NanoDrop One, Thermo Scientific). We used a total of 14 microsatellites, described as polymorphic for the genus Trachemys (Simison et al., 2013), one (TSC-108) was designed for this study using Primer3 software (Rozen & Skaletsky, 2000). DNA amplification was performed at a final volume of 20 µl per reaction, which contained: 20 ng of DNA, 2 µl of ultrapure water (Invitrogen), 100 µM of forward primer, 100 µM of reverse primer, 15 µl of HotStarTaq ® Master Mix Kit (Qiagen) (5 units/µl of polymerase, 20 mM of TRIS-Cl, 100 mM of KCl, 0.1 mM EDTA and 400 µM of dNTPs). The amplifications were carried out in a Bio-Rad T100 thermocycler under the following conditions: 94 °C for 7 min, followed by 45 cycles consisting of denaturation at 92 °C for 1 min, primer annealing temperature of 54 °C or 55 °C (see Table 1) for 1 min, an extension at 72 °C for 1 min, and a final extension at 72 °C for 7 min. To visualize the PCR products, horizontal electrophoresis was conducted using high definition agarose gels (agarose 1000, Invitrogen) at a concentration of 3.4%. Gels were digitalized in a UV transilluminator Bio-print CX4 (Vilber Loumar) and the results were interpreted using the BioVision software (Vilber Loumar). Microsatellite reproducibility was tested considering the allelic dropout (when a heterozygote is typed as a homozygote; Valière, 2002) and false allele (when homozygote is typed as a heterozygote; Valière, 2002) rates using Gimlet 1.3.3 (Valière, 2002). The criteria used to test the genotyping error were: (1) analysis based on 5% of samples (eight individuals chosen at random) from the total sample size as suggested by Bonin et al. (2004), (2) the same DNA extract was used in the reproducibility test and experiments, and (3) each DNA sample was replicated three times (twice in the same PCR run and once in the other PCR run) (Machkour-M’Rabet et al., 2017). In addition, the resolution level of agarose gel electrophoresis was controlled using a random subsample (10%) which was then processed by high-resolution capillary electrophoresis (QIAxcel Advanced System, Qiagen). Results showed an identical number and pattern of alleles than when using high definition agarose gel electrophoresis (Gómez-Carrasco et al., 2018).
Characteristics of the 14 Microsatellites Used for Trachemys venusta venusta in Southern Mexico Over the Whole Dataset.
Genetic diversity was estimated for the UMA and wild population, by calculating the number of alleles (NA), effective number of alleles (Ne), observed heterozygosity (Ho), and expected heterozygosity (He) using Genalex 6.502 (Peakall & Smouse, 2012). To test for differences in these parameters between the UMA and wild population, a Mann-Whitney U test was performed using Past 3.21 (Hammer et al., 2001). Allelic richness (NA), as with other measures, does not consider sample size; therefore, we determined allelic richness (AR) following the rarefaction method implemented in HP-Rare 1.0 (Kalinowski, 2004, 2005) which takes into account differences in sample size. Additionally, the index of the polymorphic content (PIC) was determined using Cervus 3.0.7 (Kalinowski et al., 2007). Levels of inbreeding were estimated using three different parameters: (1) at the population level, we determined the inbreeding coefficient (F) using Genalex 6.502, (2) at the individual level, we determined the internal relatedness (IR) which evaluates the relatedness of the parents of an individual (Amos et al., 2001) (outbred individuals present values equal to or below zero, and positive values indicate some relatedness of parents with a maximum value of 1), and (3) the homozygosity by loci (HL) which is a homozygosity index that weights the contribution of each locus depending on their allelic variability (Aparicio et al., 2006), with values that range from 0 (all loci are heterozygous) to 1 (all loci are homozygous). The Hardy-Weinberg equilibrium (HWE) for heterozygote deficit was tested with Genepop 4.0.10 (Raymond & Rousset, 1995; Rousset, 2008), and the presence of null alleles at each loci was estimated with FreeNa (Chapuis & Estoup, 2007). Finally, we evaluated the effective population size (NE) using NeEstimator 2.1 (Do et al., 2014) considering the linkage disequilibrium method with the lowest allele frequency at a critical level of 0.05, and the jackknife method to evaluate the 95% confidence interval (CIs) as recommended by Waples and Do (2008).
To test a possible bottleneck event, we compared the levels of He excess related to the expected equilibrium heterozygosity (Heq) using Bottleneck 1.2.02 (Cornuet & Luikart, 1996; Piry et al., 1999). Two models of mutational equilibrium were considered: a stepwise mutation model (SMM) and a two-phase model (TPM), with 95% of the mutation assumed in a single step for the TPM model and a variance among multiple steps of 12 as recommended by Piry et al. (1999). Significance was evaluated from the one-tailed Wilcoxon rank test which is more appropriate and powerful for less than 20 microsatellites (Piry et al., 1999). Bottleneck event will be considered only if both models (SMM and TPM) are significant.
To determine the level of genetic differentiation between the UMA and the wild population, the value of FST was estimated with Genalex 6.502. To identify the genetic structure between the UMA and the wild population, the Bayesian analysis implemented in Structure 2.3.1 (Pritchard et al., 2000) was applied. This method determines the optimal number of groups/clusters (K) and assigns to each individual a membership probability (qi) to new clusters. The admixture and allele correlated frequency models were used. To determine the optimal number of clusters, the program was executed 10 times for different numbers of K (K from 1 to 5), and for each run the Markov Chain Monte Carlo (MCMC) algorithm was executed with a burn-in period of 100,000 steps followed by 100,000 steps. We used the Evanno method (ΔK method; Evanno et al., 2005), implemented in the Structure Harvester website (Earl & Von Holdt, 2012) to determine the best value of K that fit with our data. This method assigns individuals to the number of clusters considered as optimal; therefore, we decided to conduct a principal coordinate analysis (PCoA) to investigate genetic structure without a priori assignment. Finally, a hierarchical analysis of molecular variance (AMOVA; 9999 permutations) was processed considering region level (UMA vs wild populations) and population level (different wild localities) using an allelic distance matrix as input. These final analyses (PCoA and AMOVA) were conducted using Genalex 6.502.
To evaluate pedigree relationships in the UMA and wild localities, we used the ML-Relate program (Kalinowski et al., 2007). This calculates the likelihood of each type of relationship for each pair of individuals considering the following pedigree relationships: unrelated (U), half-siblings (HS), full-siblings (FS), and parents-offspring (PO). For each pair of individuals, the program provides by default the highest likelihood of relationship, then suggests the highest relationship. ML-Relate provides a list of several possible relationships (the putative corresponding to the highest likelihood, and the alternatives) for each individual pair based on a test (statistical test and simulations) that determines which are consistent with the data (confidence set option; 0.05 level of significance and 1000 simulations) (Kalinowski et al., 2007). Based on these relationships for each female-male pair from the UMA, we used a likelihood ratio test (1000 simulations; specific hypotheses test option) to obtain a P value: if P is small (P < 0.05) the alternative hypothesis is rejected, thus we maintained the highest relationship proposed by the highest likelihood. If the P value is large (P > 0.05), putative and alternative relationships are consistent with the data and the lowest relationship was selected (U < HS < FS < PO) (Kalinowski et al., 2007).
Furthermore, relatedness coefficient (r) among each female-male pair in the UMA was evaluated using Storm (Frasier, 2008) and ML-Relate programs.
A total of 86 individuals from the UMA and 45 from the wild population were successfully amplified for 14 microsatellite loci. However, four loci contained null alleles (Table 1) and were therefore eliminated; the following analysis only considered the 10 loci without null alleles. Genotyping error rates could be considered null (average allelic dropout rate of 0% and average false allele rate of 0%; as errors were not significant) when compared with other studies (Valière et al., 2006; Table 4).
Globally, genetic diversity parameters (Na, AR, Ne, Ho, He, PIC) were slightly higher for the UMA than for the wild population but never significant (Table 2) and both parameters used to evaluate the allelic richness (Na and AR) were similar (Table 2). The effective population size (NE) was low, with a value of 91 for the UMA and 110 for wild populations. The inbreeding coefficients (F) were positive for both groups but with a lower value for the UMA than for the wild population. The HWE test for deficiency of heterozygotes was only highly significant for the wild population (Table 2). At the individual level, internal relatedness (IR) was very low for the UMA and high for the wild population with high individual variability for both groups (Supplemental Figure S1). The level of homozygosity (HL) was high for the UMA (a mean of 33% of loci are homozygous per individual) and very high for the wild population (a mean of 55% of loci are homozygous per individual) (Table 2), with a higher number of loci presenting a high level of homozygosity (>0.5) for the wild population (Supplemental Figure S2). No recent bottleneck was detected for any population as suggested by the L-shape graph (Figure 2) that shows a characteristic L-shape distribution (alleles with low frequency are the most numerous), nor by both heterozygote excess tests (SMM and TPM) (Table 2).
Statistical Summary for Trachemys venusta venusta in Southern Mexico Using 10 Polymorphic Microsatellites for Founder Individuals of the UMA and for the Wild Population.
The FST value between the UMA and wild population was low (0.031) but significant (P = 0.0001); the results of AMOVA at the region level presented only 2% of variance (P = 0.0017) and 4% for the population (P = 0.001), while the rest of the variance corresponded to individual variance (see details of analysis in Supplemental Table S1). The Bayesian analysis performed with Structure identified an optimal number of two groups (K = 2) using the ΔK method (Supplemental Figure S3). All individuals were relatively well partitioned into the two new clusters according to their Q probabilities (Table 3), showing good separation between individuals from the UMA group (purple in Figure 3) and individuals from the wild population (blue in Figure 3). Nevertheless, some individuals of the UMA present a genetic profile characteristic of a wild population. A subsequent Structure analysis for each new cluster did not permit the identification of a sub-structure (not shown). As observed for Bayesian analysis, PCoA analysis showed some separation between the UMA (purple in Figure 4) and wild individuals (blue in Figure 4). Wild individuals have a tendency to be more on the negative side of the x-axis and positive side of y-axis while captive individuals present a wide range of distribution on the x-axis. Individuals from the La Venta locality (blue triangle) are more separated than others (blue circle and square).
Probability Q of Membership of Individuals in Each Initial Group to Each of the Inferred Clusters Determined by Structure 2.3.1.
Relatedness coefficient (r) determined using Storm and ML-Relate clearly shows some kinship relations in the initial founders of the UMA. Storm and ML-Relate analysis presents values from −4.78 to 0.67 and 0 to 0.63 respectively, considering all types of dyad (female x female, male x male, and female × male). However, considering that we are interested in genetic management to optimize reproduction, we focused our analysis on female x male pairs only, which represent a total of 949 dyads. Using the relatedness coefficient generated by Storm we determined the proportion of each pedigree relation (PO, FS, HS, and U) based on the following criteria: values ≤0 are assigned to unrelated (U), values ≤0.25 are assigned to half-sibling (HS), and values > 0.25 are assigned to FS or PO without distinction (Kamel et al., 2012; Queller & Goodnight, 1989). We obtained 45% of U (n = 427), 39% of HS (n = 367), and 16% of FS or PO (n = 155). Our pedigree predictions differed considerably from ML-Relate relationships based on the highest likelihood with the following results; 79.7% of U, 15.4% of HS, 3.8% of FS, and 1.2% of PO for the same dataset (a total of 193 pairs of individuals presented half-sibling or greater kinship relation). Relatedness coefficients corresponding to these pedigree predictions were highly variable (Figure 5, full color histogram) with a mean of 0.039 (±0.067SD, min = 0.000, max = 0.474) for U (blue in Figure 5), 0.230 (±0.068SD, min = 0.087, max = 0.388) for HS (yellow in Figure 5), 0.383 (±0.083SD, min = 0.230, max = 0.625) for FS (green in Figure 5), and 0.536 (±0.046SD, min = 0.500, max = 0.611) for PO (red in Figure 5). Applying the likelihood ratio test to compare the putative relationship (highest value of likelihood) with alternatives proposed by the confidence interval, we obtained the following results: 91.0% of U, 7.5% of HS, 1.4% of FS, and 0.1% of PO for the same dataset (a total of 85 pairs of individuals presented half-sibling or greater kinship relation). The relatedness coefficient (Figure 5, hatched color histogram) produced some highly variable results with a mean of 0.060 (±0.085SD, min = 0.000, max = 0.500) for U (blue in Figure 5), 0.342 (±0.067SD, min = 0.228, max = 0.511) for HS (yellow in Figure 5), 0.509 (±0.099SD, min = 0.340, max = 0.625) for FS (green in Figure 5), and only one data for PO (red in Figure 5).
To optimize reproduction by using optimum females and males in terms of kinship, we analyzed all female x male dyads (total of 949) and identified which pairs presented half-sibling or a greater kinship relation. In the first instance, we considered kinship based on the highest likelihood prediction (bold in Table 4) and subsequently, we looked at which of these relationships changed to unrelated considering the likelihood ratio test results (grey in Table 4). Of 73 female founders from the UMA, only 19 have no strict kin relationship with the 13 males (females highlighted in grey in Table 4). All other females present a minimum of one kin relationship with the male founders of the UMA.
Pedigree Relationships for All Female × Male Dyads (n = 949) of the UMA Determined by ML-Relate.
Results for wild localities showed a very low number of pedigree relationship: 1) Bosques de Saloya; out of a total of 66 dyads (n = 12 individual), only 8 dyads were identified that demonstrated some type of relationship (two FS and six HS), 2) La Venta Park Museum; out of a total of 406 dyads (n = 29), 54 dyads were identified with some type of relationship (one PO, 20 FS, and 33 HS), and 3) Miguel Hidalgo, from a total of 6 dyads (n = 4 individuals) no pedigree relationship among individuals was observed.
There are numerous studies on Trachemys spp. related to topics such as general ecology (Gibbons, 1990; Gradela et al., 2017; Works & Olson, 2018), biogeography (Ennen et al., 2017), systematics (Ernst & Seidel, 2006; Fritz et al., 2012; McCord et al., 2010; Parham et al., 2015; Vamberger et al., 2020; Vargas-Ramírez et al., 2017), and even as an invasive species (Cadi & Joly, 2003; Martins et al., 2018; Rodrigues et al., 2016; Standfuss et al., 2016). One reason that this particular species has been subject to a great deal of research could be that many species in this genus are included within the IUCN Red List (International Union for Conservation Nature) and TFTSG list (Tortoise and Freshwater Turtle Specialist Group) or are not yet evaluated (Rhodin et al., 2018). Genetic studies are also largely represented, generally focusing on phylogenic or biogeographic considerations (Fritz et al., 2012; Parham et al., 2015; Vamberger et al., 2020; Vargas-Ramírez et al., 2017), hybridization process (Parham et al., 2013,), or even identifying genetic damage under radioactive effects (Lamb et al., 1991); however, population genetic studies still comparatively rare (Martínez et al., 2007; McGaugh, 2012; Scribner et al., 1984; Smith & Scribner, 1990). Particularly, for Trachemys venusta, no population genetic study is reported to date. We present, the first results related to the population genetics of wild population of T. venusta venusta in Mexico compared to founder individuals of a captive group (Mexican UMA).
With the exception of some studies that present low genetic diversity (value of He < 0.4; Echelle et al., 2010; Vargas-Ramírez et al., 2012), freshwater turtles generally exhibit high values of genetic diversity (generally He > 0.6–0.7; Vargas-Ramírez et al., 2012, Table 4 for examples; Davy et al., 2014) when using methodology based on microsatellites. Our values could be considered lower (He < 0.6) than many other freshwater turtles or even threatened species (Davy et al., 2014; Petre et al., 2015; Vargas-Ramírez et al., 2012, Table 4 for examples), and correspond with values for wild populations of Elusor macrurus Cann & Legler 1994 (Testudines, Chelidae) an endangered Australian freshwater turtle (Schmidt et al., 2018). Very few studies report values of genetic diversity for Trachemys spp. based on microsatellites. The first study (McGaugh, 2012) demonstrated low values for T. taylori Legler, 1960 (He from 0.442 to 0.569; Ho from 0.458 to 0.590) and high values for T. scripta elegans Wied-Neuwied, 1839 (He = 0.809; Ho = 0.828). The recent study by Vamberger et al. (2020) showed relatively similar results for T. scripta elegans (He = 0.841; Ho = 0.690), T. scripta scripta (He = 0.753; Ho = 0.720), and T. scripta troostii (He = 0.757; Ho = 0.663). Other population genetic studies that report levels of genetic diversity for Trachemys spp used allozymes which makes comparison with microsatellites data difficult. In this respect, all studies showed very low values of genetic diversity (H < 0.15; Martínez et al., 2007; Scribner et al., 1984; Smith & Scribner, 1990). The relatively low value of genetic diversity observed for wild populations of T. venusta venusta compared to other turtle species, suggests that this species deserves more attention and probably a higher level of protection considering the threshold of 0.54 proposed by Willoughby et al. (2015) as a value to consider a species as Critically Endangered (Schmidt et al., 2018). We cannot discard the possibility that T. v. venusta individuals have been and are still translocated, as mentioned for Dermatemys mawii (González-Porter et al., 2011), which could affect the level of genetic diversity. To confirm this hypothesis more ecological studies are necessary on this highly preyed species (Reynoso et al., 2016). Furthermore, the low genetic diversity observed could be the result of a Wahlund effect considering that wild individuals originate from different places and could belong to different populations. Extensive sampling over the distribution range of the species will be necessary to understand the genetic structure of T. v. venusta.
Values of genetic diversity for captive populations were slightly higher than wild population as reported for Elusor macrurus (Schmidt et al., 2018) and Dermatemys mawii Gray, 1847 (Testudines, Dermatemydidae) (Gallardo-Alvárez et al., 2019). Because founder individuals of the UMAs (captive population) originate from two distinct localities: Pomposu that is located within the Grijalva River Sub-Basin and El Espino located in the Usumacinta River Sub-Basin (see Figure 1), the genetic pool of these founders could be higher due to the mix of different genetic sources that will reflect the genetic structure of the wild populations related with the river basins (see below in the discussion for more detail). This situation was reported for Lithobates sevosus Goin & Netting, 1940 (Anura, Ranudae), a critically endangered frog from the southeastern USA (Hinkson et al., 2016), and for Dermatemys mawii a critically endangered freshwater turtle (Gallardo-Alvárez et al., 2019), suggesting that the captive population is genetically representative of natural populations (Hinkson et al., 2016) and could be used to found new wild populations (Edwards et al., 2014). Also, founder individuals are approximately 40 years of age, and genetic diversity represents the situation from several decades ago. The genetic diversity of the wild population of T. venusta venusta may have declined since the UMA was founded, probably due to habitat fragmentation and water pollution (Tudela, 1992). Indeed, the state of Tabasco has experienced widespread environmental degradation, losing around 60% of its wetlands at the beginning of the 21st century mainly due to anthropogenic activities such as the oil industry, the establishment of new crop areas and grassland for livestock use, road construction, population growth, and increased pollution derived from these activities (Melchor et al., 2017; Palomeque de la Cruz et al., 2017). Additionally, in the state of Tabasco, freshwater turtles have experienced intense hunting and illegal trading and trafficking (PROFEPA, 2014, 2015; Reynoso et al., 2016) that is likely to have contributed to a decline in population size and number, leading to a loss of genetic diversity in local populations. It is reported that in Villahermosa and surrounding cities (e.g., Nacajuca, Comalcalco, Jalpa de Méndez, and Cunduacán) freshwater turtles were abundant, but today it is difficult to find turtles such as T. venusta (Reynoso et al., 2016). Many species of turtles are characterized by a long generation time (often > 25 years), making a loss of genetic diversity or genetic differentiation among populations due to recent and/or anthropogenic events (e.g., habitat fragmentation) very difficult or impossible to detect (examples in Davy et al., 2014). However, sexual maturity of T. venusta has been evaluated at approximately four to seven years (Torres et al., 2011) which results in five to 10 generations (40 years) between UMA founder individuals and the wild population analyzed in this study, thus allowing the observation of a loss of genetic diversity in the wild population due to recent anthropogenic pressure.
Although captive population (UMA) did not demonstrate a significant loss of heterozygotes, the wild population shows a significant loss of heterozygotes that could indicate a high level of inbreeding and small population size that are probably a consequence of recent anthropogenic pressures (see above) on freshwater turtles in Tabasco state. Nevertheless, we cannot discard the influence of some level of pedigree relationship among individuals in two localities or the Wahlund effect. This hypothesis of higher inbreeding in wild populations is reinforced by our results that showed a significant HWE deviation and higher values of internal relatedness and homozygosity by loci for the wild population. Indeed, a positive value (maximum of 1) of internal relatedness (higher value of IR for the wild group than captive group) indicates that individuals are a result of consanguineous mating (O’Leary et al., 2013). Similarly, the higher value of HL (homozygosity by loci) for the wild population than the UMA indicates higher homozygosity in the wild group which could be result of higher inbreeding (Frankham et al., 2002). Several species of freshwater turtles, such as Apalone spinifera emoryi Agassiz, 1857 (Testudines, Trionychidae), Mesoclemmys dahli Zangerl & Medem, 1958 (Testudines, Chelidae), Chrysemys p. picta Schneider, 1783 (Testudines, Emydidae), and Clemmys guttata Schneider, 1792 (Testudines, Emydidae) have shown evidence of genetic isolation, genetic differentiation, as well as modest to high inbreeding rates, but surprisingly the values of heterozygosity in these species are medium to high (0.6–0.7) despite experiencing anthropogenic pressures (Buchanan et al., 2019; Gallego-García et al., 2018; Mali et al., 2015a). As previously mentioned, the decrease in genetic diversity of many species of turtles, even after a prolonged decrease in population size, may not be observed due to their long generation times and late maturity associated with chelonian life history (Buchanan et al., 2019; Kuo & Janzen, 2004; Willoughby et al., 2013).
Bayesian, and in a lesser extent PCoA, analyses of genetic structure among wild and captive populations revealed some good separation between both, probably reflecting the low-connectivity among sites where founder individuals originated. Certainly, the Pomposu ecosystem is largely dependent on the Grijalva River Sub-Basin (Figure 1) as well as wild populations considered in this work (captive individuals showing a dominant wild genetic profile in analyses; blue on Figure 3). El Espino is dependent on the Usumacinta River Sub-Basin (Figure 1) and could be represented by captive individuals with the alternative genetic profile shown in Structure analysis (purple profile on Figure 3). Therefore, the population structure of this freshwater turtle could be conditioned by the river basin structure; however, more extensive sampling along both rivers would be necessary to confirm this. Consequently, both genetic profiles observed in our analyses (blue or purple in Structure and PCoA analyses) could reflect the separation between the two principal river basins in Tabasco (Grijalva River and Usumacinta River Sub-Basins) more than between the UMA and wild groups.
The relationship coefficient evaluated for the captive population (Tabasco State Government UMA) suggests a relatively low level of kinship (considering all pairs: female-male, female-female, and male-male) which probably reflects the origin from different populations of the UMA founder individuals. However, we highlighted the differences observed between programs (Storm and ML-Relate) and between methods used in ML-Relate (application or non-application of the likelihood ratio test to evaluate the more probable relationship) which can give values of 2-fold difference (e.g., 45% for unrelated using Storm and 91% using ML-Relate with ratio test).
For all founder females of the UMA (n = 73), the stricter analysis provided by ML-Related enabled the identification of 19 females that present no kin relationship with males (n = 13). This could result in a total of 247 potential dyads to establish an optimal breeding program in the UMA. In order to develop optimal reproduction in a captive population, the number and identity of the founders is a key factor for the genetic pool (Witzenberger & Hochkirch, 2011). According to these authors, who reviewed 188 genetic studies of breeding programs, to maintain 90% of the natural genetic diversity for 100 years (meaning a He = 0.54) in captive population, 15 founder individuals and a captive population size of at least 100 individuals are necessary to preserve this heterozygosity threshold. Considering that the UMA have a sufficient unrelated founders for reproduction (>15) and more than 100 individuals, the UMA management has exceeded the benchmarks recommended by Witzenberger and Hochkirch (2011). Further, in order to improve breeding program and the genetic diversity of captive population, it will be recommended to introduce wild-caught individuals. On the other hand, if local T. venusta venusta wild populations decline to critically low numbers or show evidence of inbreeding depression, the captive population in Tabasco could be used for reintroduction to the wild. To implement a reintroduction program, it would need to take into consideration the population demographics and genetic diversity of the species in the reintroduction sites. In this scenario, we recommend the use of the best females (no kin relationship with all males of the UMA), and the introduction of new founders to avoid genetic erosion (Witzenberger & Hochkirch, 2011). To date, very few studies have used genetic information of founders to optimize breeding programs and avoid inbreeding in turtle programs. For example, J. M. Miller et al. (2018) evaluated the pedigree on three seasons of a captive breeding program of the Galápagos giant tortoise, revealing that the genetic diversity of the progeny was reduced in a single generation accompanied by a tendency towards the reduction of their physical shape when more related individuals were bred. In addition, Spitzweg et al. (2018) analyzed the genetic diversity and kin relationship of a turtle (Batagur baska Gray, 1830; Testudines, Geoemydidae) in captivity, and showed that most of the founder individuals that came from the wild demonstrated a high kinship relationship.
Successful management and conservation programs for long-lived organisms are those that recognize the need for protection and biological knowledge of all life stages of the species they breed (Congdon et al., 1993). For this reason, our results have important implications for the conservation and management of T. venusta venusta, as they contribute to a better selection of pairs, decreasing the possibility of inbreeding in the generations born in captivity. Furthermore, through good genetic management, it could be considered to develop a reintroduction program from the UMA individuals considering previous genetic studies of wild populations to avoid jeopardizing the genetic health of local populations. The program conducted with the highly endangered species of freshwater turtle Batagur trivittata Duméril & Bibron, 1835 (Testudines, Geoemydidae) is an example of a successful genetic management program, where the breeding, reproduction, and release of captive individuals has resulted in an improvement in the wild population. A captive management program has been established for this species, and it was thought that only 12 breeding turtles existed in the wild. Since 2002, more than 700 individuals of this species have been obtained and after selecting individuals with a high genetic diversity and reintroducing them to their natural habitat, the program reported an increase in the fertility of the eggs (Çilingir et al., 2017). Furthermore, as has been suggested for the turtle Trachemys scripta elegans Wied-Neuwied, 1839 (Testudines, Emydidae) that is principally exploited as a pet (Mali et al., 2015b), our results could be used to create predictive models of exploitation since T. venusta venusta is a highly valued species for human consumption and also as a pet.
Finally, we believe that it is convenient to analyze the offspring of the founders already in the UMA in order to improve the genetic management of this captive population. Another recommendation to consider in the implementation of a genetic management program for T. venusta could be the storage of sperm and the study of the multiple paternity, since it has been reported in several species of turtles (Edwards et al., 2014) as a reproductive strategy that serves to maximize the genetic diversity of the offspring of long-lived organisms (Davy et al., 2011).
To conclude and based on this first experience to optimize management in captivity through genetic tools, we propose a concise protocol to establish optimal mates. Initially, all individuals can be genotyped using an adequate molecular marker (e.g., microsatellites, SNP), followed by determining the genetic diversity (homozygosity by loci) for each individual and the pedigree relationship between female and male individuals. Based on these results, the optimal breeding colony could be established considering the likelihood ratio test to verify and confirm the level of kinship. Once the breeding colony is created, it is important to monitor the genetic diversity of the offspring. Furthermore, the inclusion of new individuals from illegal trafficking seizures could enrich genetic diversity by introducing new alleles to the UMA. This would require the known origin (to the extent possible) and genotyping of new individuals before entering the breeding group. Subsequently, for the release of individuals into the wild, a translocation program following IUCN guidelines (International Union for Conservation of Nature/Species Survival Commission, 2013) would be implemented. The results obtained from this study could be applied following a basic reproduction protocol such as suggested by D. A. Williams and Osentoski (2007, Figure 1).
Implications for Conservation
Wild populations of freshwater turtles have declined worldwide mainly due to human activities such as overexploitation, water pollution, flow modification, destruction or degradation of habitat and invasion by exotic species (Dudgeon et al., 2006; Lovich et al., 2018); in addition, 12 new emerging threats have been detected that are affecting freshwater environments (changing climates; E-commerce and invasions; infectious diseases; harmful algal blooms; expanding hydropower; emerging contaminants; engineered nanomaterials; microplastic pollution; light and noise; freshwater salinization; declining calcium and cumulative stressors) (Reid et al., 2019) which increase the vulnerability of freshwater turtle populations. Another pressure to which this group of chelonians has been subjected is exploitation for human consumption; the main consumers being Asian countries, of which China is reported as the main consumer of freshwater turtles in the world (Gong et al., 2009; Shiping et al., 2006). The consumption of chelonians is also recorded both in Central America and in southeastern Mexico, being a traditional activity inherited from Mesoamerican cultures (Guevara Chumacero et al., 2017). For example, traces of Trachemys venusta shells have been found frequently in archeological zones of Mayan culture, indicating that this species has been linked since pre-Hispanic times to human gastronomy. Today the communities of southeastern Mexico and Central America prefer to consume T. venusta (Páez et al., 2012; Vargas-Ramírez et al., 2017), due to the reduction of other species of freshwater turtles such as Dermatemys mawii (Reynoso et al., 2016).
However, despite the pressures that T. venusta has historically undergone, few studies address the population and genetic situation in the species. The preponderance of research associated with T. venusta has focused on clarifying its phylogeny given that T. venusta is often considered a subspecies of T. scripta (Ceballos & Brand, 2014; Ernst & Seidel, 2006; Seidel, 2002). This situation creates a conflict in establishing an appropriate risk categorization; for example, in Mexican legislation, the species T. venusta does not appear as a species, instead, this law considers the name T. scripta as a species (SEMARNAT, 2010) and its risk categorization is threatened; similarly, the IUCN does not consider T. venusta as a species. Within its risk classification, it considers the name T. scripta as the species, its conservation status as a minor concern and its populations stable (van Dijk et al., 2011).
Therefore, our work is important because it addresses two points that have not been studied in this species: 1) The information is relevant because the results found may indicate the first detrimental effects on wild populations of this species caused by anthropogenic pressures. Consequently, due to the data we obtained and the lack of existing studies on this turtle, we suggest more studies that contribute to an increase in knowledge and improved understanding of the population genetics and dynamics of T. venusta, in order to make a satisfactory update to the conservation status of this species. 2) Genetic evaluation in captivity and a proposal of reproductive breeding; captive breeding has been a strategy for several decades for the recovery and exploitation of wildlife species (Farquharson et al., 2017; Mandimbihasina et al., 2020; D. A. Williams & Osentoski, 2007; Witzenberger & Hochkirch, 2011). The turtle T. venusta has been successfully bred in the UMAs of southeastern Mexico; however, the reproduction of this species has been mainly for commercial purposes, thus our work would help enrich future captive management programs that may be created. Furthermore, it helps to identify the genetic diversity of the reproductive individuals of the studied UMA, one of the most important sites for the reproduction and breeding of freshwater turtles in southeastern Mexico.
The data obtained showed that founding individuals in the UMA presented greater genetic diversity than in wild populations. This knowledge can be used to obtain individuals with a high genetic diversity by applying the proposed crosses; these individuals can then be released to enrich the genetic diversity of wild populations if for example, a demographic study suggests a decrease in population size or a phylogeographical research demonstrates a need to improve connectivity among populations with new intermediate populations. Our research could also be used as an example for other captive breeding sites to replicate if required. The use of genetic and genomic methods to help endangered species could permit a reasonable selection of founding individuals to develop captive breeding programs; the proper selection of founding individuals would ensure the protection of a greater genetic diversity in the wild and would eventually allow the reintroduction of a progeny with a high genetic diversity (W. Miller et al., 2010).
The capture of organisms was done with the permit number SPGA/DGVS/01085/16 issued by the Mexican Ministry of Environment and Natural Resources (SEMARNAT). We would also like to thank Dr. Humberto Hernández Trejo for his contribution of information to strengthen the discussion. Thanks to Dr. Leon David for helping to produce the Figure 1. Finally, we are particularly grateful to the anonymous referees for their time and their valuable comments.
Data Accessibility Statement
The morphometric and weight data taken from Trachemys venusta venusta can be found at Knowledge Network for Biocomplexity Digital Repository: https://doi.org/10.5063/F1C53J44.
Declaration of Conflicting Interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the Mexican Commission for the Knowledge and Use of Biodiversity (CONABIO) with the funding number NE017-project: “Establishment of a genetic management program for the species Dermatemys mawii (white turtle) and Trachemys venusta (hicotea) in Wildlife Management Units (UMA) to increase the genetic flow and connectivity of the Mesoamerican Biological Corridor in Tabasco.” Thanks to the Mexican National Council of Science and Technology (CONACyT) to provide financial support to Elsi Beatriz Recino Reyes (scholarship #485105), Salima Machkour-M’Rabet (fellowship #217950), Julia M. Lesher Gordillo (fellowship #240997).