Introduction
During the twentieth century, due to the high demand for skin and meat in national and international trade, Mexican crocodile species became almost extinct because of the intensive hunting. At the same time, ecosystems began to transform, and crocodile habitats became fragmented, contaminated, and seriously limited (Casas-Andreu, 1995; Ross, 1998; SEMARNAP, 2000). In 1970, the Mexican government declared a total and permanent ban on the three crocodilian species distributed in Mexico; the American crocodile (Crocodylus acutus), Morelet’s crocodile (C. moreletii) and the caiman (Caiman crocodilus) (Casas-Andreu, 1995), which are currently protected by the Norma Oficial Mexicana Nom-059- SEMARNAT-2010 (SEMARNAT, 2010). Recently, C. moreletii was transferred from Appendix II to I of the Convention on International Trade in Endangered Species of Wild Fauna and Flora (CITES, 2021).
The American crocodile (C. acutus), that is distributed mainly along the Mexican coast of the Pacific Ocean, is considered a vulnerable species by the International Union for Conservation of Nature (IUCN) Red List. On the other hand, the Morelet’s crocodile is distributed in the marshy areas of the coastal region of the Caribbean and Gulf of Mexico and is considered by the IUCN as a least concern species (IUCN, 2021).
In the Yucatan Peninsula, both species are sympatric and eventually interbreed originating hybrids (Cedeño-Vázquez et al., 2008; González-Trujillo et al., 2012; Hekkala et al. 2015; Pacheco-Sierra et al., 2016; Ray et al., 2004; Rodríguez et al. 2008). Due to the need for generating genetic information on crocodile populations, molecular markers have been used as tools to understand the effect of various genetic factors causing the decline of population sizes. It is known that population fragmentation has an impact on the effective size of the population and therefore on the loss of genetic diversity and the ability of individuals to adapt to new environmental changes (Rocha & Gasca, 2007). Genetic analysis allows the study of genetic diversity and gene flow among populations, as well as to determine their structure and measure the effect of inbreeding crosses within a population (Amos & Balmford, 2001; Rocha & Gasca, 2007). Existing legislation and treaties governing the trade in wildlife, such as the Convention on the International Trade of Endangered Species (CITES) and the United States Endangered Species Act (ESA), are based on the recognition of distinct population or taxonomic units.
So far, genetic studies using microsatellites and mitochondrial deoxyribonucleic acid (mtDNA) have been conducted in the region of the caribbean and Gulf of Mexico by Cedeño-Vázquez et al. (2008), Dever et al. (2002), Machkour-M’Rabet et al. (2009), Ray et al. (2004) and Rodríguez et al. (2008). While Pacheco-Sierra et al. (2016) found evidence that these two species have hybridization areas and maintain genetic flow among their populations along distribution range.
A database of single gene ‘‘barcodes’’ has been proposed to classify the complete diversity of life (Hebert et al., 2003; Ratnasingham & Hebert, 2007). The region of the mitochondrial cytochrome c oxidase I (COX1) gene has been recommended as a standard for DNA barcoding as well as for the evaluation of genetic diversity and monitoring of legal and illegal trade of species (Eaton et al., 2010; Hebert et al., 2003; Ratnasingham & Hebert, 2007). Therefore, the objective of this study was to evaluate genetic diversity and some phylogenetic relationships in wild and captive populations of C. acutus and C. moreletii located in Southeastern Mexico using the Barcode of Life Data System (COX1, cytochrome C oxidase subunit 1 gene).
Materials and methods
Collection of samples: Species were identified according to the key characters as the scales embedded in the base of the tail, which correspond to fewer pronounced osteoderms and the wide snout in C. moreletii compared to C. acutus as Platt and Rainwater (2005) suggest. In 2013, blood and/or tissue samples were collected from five animals at the UMA-CICEA in Tabasco and nine samples were collected from a private farm in Campeche (Fig. 1). Blood was drawn from the supravertebral venous sinus and the blood samples were stored in Vacutainer® tubes with 4 ml of EDTA K3; tissue samples were preserved in ethanol. In addition, el Colegio de la Frontera Sur provided seven samples from Quintana Roo; and the samples from Oaxaca (N = 12) and Guerrero (N = 12) were provided by Serrano-Gómez (2010); The total number of samples was 45 individuals, 28 from C. acutus and 17 samples from C. moreletii (Fig. 1).
DNA Extraction and Sequencing: DNA was purified from blood or scales following the protocol from the DNeasy Blood and Tissue kit, (QIAGEN, Austin, Texas). Using the polymerase chain reaction (PCR), a 641 base pair (bp) fragment of the corresponding COX1 gene was amplified; this fragment corresponds to position 52-693 of the sequence of a C. acutus (GenBank of GQ144571, (Eaton et al., 2010)). The modified primers of Meusnier et al. (2008) were used: forward 5’-TCCACTAATCACAARGATATTGGTAC-3’ and reverse primer 5’-CCTCCGGGTGCCCGAAGAATCAG-3’.
PCR reactions was carried out in a final volume of 20 μl containing 100 ng of total DNA, 0.6 pmol of each primer, 0.2 mM of dNTP’s, 1.5 mM of MgCl2, 10 mM of TrisHCl pH 8.4, 50 mM of KCl, 10 µg/ml of gelatin; 150 µ/ml of BSA, Triton X100 0.1 % and 1 U of Taq polymerase (BioTecMol, México).The amplification conditions included an initial stage at 94 °C for 5 minutes followed by 30 cycles of 94 °C for 30 seconds, 56 °C for 30 seconds and 72 °C for 1 min and a final stage of 72 °C for 5 min. The amplified COX1fragments were purified with the Potassium iodide technique after electrophoresis in 3 % agarose gel (Vogelstein & Gillespie, 1979). The DNA fragments were purified, and quantified. Sanger-type sequencing of both DNA strands was performed using a commercial kit following the manufacturer’s protocol (BigDye Terminator v3.1 Cycle Sequencing Kit PE Applied Biosystems, Foster, CA). The primers used for PCR were also used for sequencing. Subsequently, the fragments obtained were read by a matrix of 16 ABI 3100 Genetic Analyzer capillaries (Applied Biosystems Inc.). Forward and reverse sequences were verified and edited in the trace editor implemented in the MEGA X (Kumar et al., 2018). The sequences obtained in the present study were submitted to GenBank (accession numbers: C. acutus: KY994087, KY994088, KY994090, KY994093, KY994094 and C moreletii: KY994097, KY994089, KY994091, KY994092, KY994095, KY994096, KY994098, KY994099).
Data analysis: A multiple alignment was performed through a Clustal W implemented in MEGA X (Kumar et al., 2018) using the C. acutus sequence in GenBank as reference (GQ144571, (Eaton et al., 2010)). To evaluate the relationship between the sequences in this study with the ones previously reported from different species, seven GenBank sequences were used (KF273840, KF273836, KF273834, KF273841, KF273838 (Bloor et al., 2015)), belonging to C. acutus from Colombia; HM636894, (Man et al., 2011), belonging to C. acutus of unreported origin; HQ585889 (Meganathan et al., 2011), belonging to C. moreletti from Belize.
The nucleotidic (π) and haplotypic (Ĥ) diversity (Nei & Kumar, 2000) by species was estimated using the software DnaSP 5.10 (Librado & Rozas, 2009). The neutrality test D (Tajima, 1989) and F test (Fu & Li, 1993) were used to measure the effect of the demographic changes of the populations on the DNA sequences using the software DnaSP 5.10 (Librado & Rozas, 2009).
Phylogenetic analysis: The best-fitting evolutionary model for the data and the gamma distribution parameters were estimated using the AICc value (Akaike Information Criterion, corrected). Non-uniformity of evolutionary rates among sites was modelled using a discrete gamma distribution with five rate categories assuming that a certain fraction of sites is evolutionarily invariable for which the routines implemented in MEGA X were used (Kumar et al., 2018).
The number of base substitutions per site was estimated as the average of all sequence pairs within groups using the Kimura’s two parameter model (Kimura, 1980; Nei & Kumar, 2000). The site rate variation was modelled using a gamma distribution (shape parameter 0.17017). The differences in the composition bias among sequences were considered in evolutionary comparisons (Tamura & Kumar, 2002). All positions containing gaps and missing data were eliminated. Standard errors were computed with the bootstrap method using 1 000 replicates (Felsenstein, 1985).
Phylogenetic trees were constructed with neighbour-joining methodology (Saitou & Nei, 1987) incorporated in the software MEGA X (Kumar et al., 2018). Additionally, a Bayesian analysis (Drummond & Rambaut, 2007) using the BEAST 1.8 program (Drummond et al., 2012). To choose the optimal model for nucleotide substitution, Modeltest 3.7 (Posada & Crandall, 1998) was used, using the value of the Bayesian Information Criterion (BIC) as a selection criterion. The “Yule Speciation” model was used to estimate the evolutionary history of the sequences. The MCMC chain ran for 150 x 106 generations, sampling every 100 generations, using a priori the normal distribution and a 95 % probability for this range (Oaks, 2011). To evaluate the convergence values, the effective sample size, and the “burnout” estimates, the Tracer 1.7 program was used. (Rambaut et al., 2018). The final tree was viewed and edited with the FigTree v.1.2.2 program.
In order to evaluate the relationship between the sequences in this study with the four New World Crocodylus species previously reported, nine GenBank sequences were used, where six belonged to C. acutus (KF273834, KF273838, KF273836, KF273840, KF273841 (Bloor et al., 2015); HM636894 (Man et al., 2011), one to C. moreletii (HQ585889; (Meganathan et al., 2011)), one to C. intermedius (JF502242, (Meredith et al., 2011)) and one to C. rhombifer (JF502247.1, (Meredith et al., 2011)); in order to obtaining an adequate tree topology, we include C. niloticus (JF502246, (Meredith et al., 2011)), C. porosus (DQ273698, (Li et al., 2007)) and C. johnsoni (HM488008, (Meganathan et al., 2011)) as external groups, considering the genus Crocodylus as a monophyletic group (Brochu, 2003).
Network analysis: Phylogenetic methods may not lead to the desired resolution at the intraspecific level due to the lower genetic diversity; complementary network approaches might be valuable alternatives for studying phylogenetic structures and haplogroups at the population level; therefore, the program PopART (Leigh & Bryant, 2015) was used to reconstruct median-joining network (Bandelt et al., 1999).
Results
The final sequences obtained began at position 52 of the COX1 gene and corresponded to the 5 370 position of the complete mitogenome of C. acutus (HM636894). A total of 28 sequences of 641 bp for C. acutus COX1 gene were obtained with 24 variable sites, 22 of them were informative and two were single mutations (singletons). For C. moreletii, 17 sequences were obtained with 24 variable sites, where 20 were informative sites and 4 were singletons. The polymorphic sites of the obtained sequences are shown in Table 1.
1 | 1 | 1 | 2 | 2 | 2 | 3 | 3 | 4 | 4 | 4 | 4 | 4 | 4 | 5 | 5 | 5 | 5 | 5 | 6 | 6 | 6 | 6 | 6 | 6 | 6 | 6 | 6 | ||||
8 | 8 | 9 | 1 | 6 | 6 | 2 | 4 | 6 | 4 | 8 | 0 | 2 | 6 | 6 | 7 | 8 | 0 | 5 | 7 | 8 | 9 | 0 | 1 | 1 | 2 | 2 | 6 | 6 | 7 | 7 | |
7 | 8 | 1 | 2 | 2 | 5 | 8 | 0 | 1 | 5 | 2 | 8 | 3 | 2 | 5 | 4 | 6 | 1 | 8 | 9 | 6 | 7 | 0 | 2 | 8 | 1 | 5 | 6 | 8 | 2 | 8 | |
GQ144571 | A | G | A | C | C | C | A | C | T | A | G | A | A | A | G | C | T | G | A | C | T | G | A | G | C | T | C | A | C | A | T |
Ca I | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . |
Ca II | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | A |
Ca III | . | . | T | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . |
Ca IV | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | G | . |
Ca V | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | C | . | . | . | . | A |
Cac08 | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . | T | . | . | . | . | . | . | . | . | . | . | . | . | . | . | . |
Cac03 | . | . | . | . | . | . | . | . | C | . | . | G | . | . | . | . | . | . | G | . | . | . | . | A | . | . | . | . | . | . | . |
Cac01 | . | . | . | . | . | . | . | . | . | . | . | G | . | . | . | . | . | . | G | . | . | A | . | A | . | . | . | . | . | . | . |
Cac05 | . | . | . | . | . | . | . | . | C | . | . | . | . | . | . | . | . | . | G | . | . | . | . | A | . | . | . | . | . | . | . |
Cm I | . | . | . | . | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm II | . | . | . | . | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | G | . | A |
Cm III | . | A | . | . | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm IV | G | . | . | . | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm V | . | . | . | G | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm VI | . | . | . | . | . | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm VII | . | . | . | . | . | T | G | T | . | C | A | . | G | G | A | . | C | A | G | T | C | A | C | A | T | C | T | G | . | . | A |
Cm VIII | . | . | . | . | T | T | G | T | . | C | . | . | G | G | A | . | C | A | G | T | . | A | . | A | . | C | T | G | . | . | A |
The upper numbers correspond to the position (pb) of the sites. The dots indicate equal nucleotides. Reference sequence of C. acutus GenBank GQ144571. Cac are haplotypes of C. acutus reported by Bloor et al. (2015).
Five haplotypes for C. acutus (CaI-CaV) and eight for C. moreletii (CmI-CmVIII) were identified in crocodile populations studied in Mexico. CaI haplotype was the most frequent and was found in both the Pacific and the Caribbean. The CaII haplotype was the second most frequent haplotype and was found in an organism in the Ventanilla lagoon (Oax.) in the Pacific, but an individual from the farm in Tabasco and another individual in Río Hondo (Quintana Roo), who were phenotypically classified as C. moreletii, they presented this haplotype. While the CaIII haplotype was found in the Chacahua lagoon (Oax), and the CaIV and CaV haplotypes appear in Fortuna (Gro).
For C. moreletii, the most frequent haplotype was CmI (Fig. 1), and it is is equal to a sequence of an animal sampled in Belize (HQ585889, (Meganathan et al., 2011)). Also with this haplotype, a phenotypically identified organism was found as C. acutus, which was sampled in Boca Paila, Quintana Roo. The CmVI haplotype was presented in Quintana Roo and Chacahua Lagoon (Oax) and CmVII and CmVIII haplotypes were from organisms identified as C. acutus, two individuals from Chacahua, Oaxaca.
Table 2 shows the results of the haplotype and nucleotide diversity of the two species. A first analysis was made including all the sequences; but since hybridization has an influence on the estimates, a second analysis was carried out, discarding the sequences of organisms with a phenotype different from their genotype. For C. acutus, the haplotype diversity was 0.455 (± 0.123) and nucleotide diversity was 0.0009 (± 0.0003); while for C. moreletii the haplotype diversity was 0.505 (± 0.158) and the nucleotide diversity was 0.0009 (± 0.0003).
Species | site | n | S | h | Ĥ ± SE | π ± SE |
C. acutus | ||||||
Ventanilla | 4 | 1 | 2 | 0.667 ± 0.2040 | 0.0010 ± 0.0003 | |
Chacahua* | 8 | 23 | 5 | 0.786 ± 0.1510 | 0.0177 ± 0.0039 | |
Chacahua** | 5 | 1 | 2 | 0.400 ± 0.2370 | 0.0006 ± 0.0003 | |
La Fortuna | 7 | 3 | 3 | 0.667 ± 0.1600 | 0.0016 ± 0.0006 | |
El Medano | 5 | 0- | 1 | 0 | 0 | |
Q Roo* | 4 | 21 | 3 | 0.833 ± 0.2220 | 0.0221 ± 0.0064 | |
Q Roo ** | 2 | 0 | 1 | 0 | 0 | |
Total* | 28 | 24 | 9 | 0.632 ± 0.1030 | 0.0103 ± 0.0029 | |
Total** | 23 | 4 | 5 | 0.455 ± 0.1230 | 0.0009 ± 0.0003 | |
C. moreletii | ||||||
L. Ilusiones | 5 | 1 | 2 | 0.400 ± 0.2370 | 0.0006 ± 0.0003 | |
Buenavista* | 3 | 20 | 2 | 0.667 ± 0.3140 | 0.0213 ± 0.0098 | |
Buenavista** | 1 | - | - | - | - | |
Campeche | 6 | 2 | 3 | 0.600 ± 0.2150 | 0.0011 ± 0.0004 | |
Q. Roo* | 3 | 21 | 3 | 1.000 ± 0.2720 | 0.0223 ± 0.0090 | |
Q. Roo** | 2 | 2 | 2 | 1.000 ± 0.5000 | 0.0016 ± 0.0007 | |
Total* | 17 | 24 | 6 | 0.675 ± 0.1170 | 0.0111 ± 0.0038 | |
Total** | 14 | 4 | 5 | 0.505 ± 0.1580 | 0.0008 ± 0.0003 |
*Localities with atypical haplotypes. **Analysis without atypical haplotypes, n: Number of sampled individuals, S: Polymorphic sites, h: Number of haplotypes, Ĥ: Haplotype diversity, π: Nucleotide diversity, SE: standard error. Q. Roo: Quintana Roo.
No significant differences were found between the values of nucleotide diversity among both species (H = 3.30, P = 0.770). Furthermore, no differences were found between nucleotide diversity between C. acutus and C. moreletii populations. The same statistical value was obtained for both tests (H = 2, P = 0.368). With the values of haplotype diversity, no significant differences were found between the C. acutus populations (H = 2.7, P = 0.440), nor between C. moreletii populations (H = 2, P = 0.368) and neither among the total of the two species (H = 1, P = 0.317). For C. acutus, the neutrality tests showed negative values for D and F, but not different from zero (D = -1.2951, P > 0.05, F = -1.1139, P > 0.05) suggesting that the analyzed sequences belonging to the COX1 gene for this species are not subject to selection. On the other hand, both indices were negative and significant for C. moreletii (D = -1.7975, P < 0.05, F = -2.4488, P < 0.05); this may be indicative that populations of this species are expanding or emerging from a recent bottleneck.
A phylogenetic analysis was carried out that includes 13 haplotypes of the partial COX1 fragment found in Mexico plus 10 haplotypes retrieved from the data bank. Fig. 2A shows the C. acutus and C. moreletii sequences in two well-defined clades. The Colombian C. acutus sequences obtained from GenBank (KF273834, KF273838, KF273836) together with HM636894 seem to form a subgroup within the C. acutus clade apart from the sequences found in Mexico, while the sequence of a Colombian C. acutus (KF273840, Cac07) is equal to the CaI and the reference sequence GQ14457; also, the sequence KF273841 (Cac08) showed that Colombian specimens are close to the Guerrero and Oaxaca populations. For its part, the CmI haplotype found in crocodiles of Tabasco, Campeche and Quintana Roo is equal to a sequence of an animal sampled in Belize (HQ585889, (Meganathan et al., 2011)). Fig. 2B was reconstructed using Bayesian analysis and the sequence of C. intermedius (HM636895) is within a specific clade of C. acutus like the neighbor phylogenetic tree. This analysis shows too the C. acutus and C. moreletii clades similar to Fig. 2A.
The median junction network (MJ) shows the two haplogroups formed by C. moreletii and C. acutus (Fig. 3), the most frequent haplotypes were CmI and CaI, respectively. The Colombian Cac07 haplotype and CaI haplotype are the same, and from this haplotype the Mexican haplotypes and the Colombian Cac05 and Cac08 haplotypes are derived. The connection between C. acutus and C. moreletii haplogroups appears in two lines that join the CmVIII haplotype. On the other hand, C. intermedius derives from CaI presenting four mutations, while C. rhombifer derives from CmVI with 20 mutations difference.
Discussion
Mitochondrial DNA is by far the most widely used population genetic marker in animals. The reasons for its use depend on its high level of variability, maternal (clonal) inheritance and mtDNA diversity that can be correlated with demographic effects (variations in population size between species or populations), which makes it reliable in the context of biological conservation (Harrison, 1989; Nabholz et al., 2008; Roman & Palumbi, 2003). In the phylogenetic tree obtained during this study using a partial sequence of the COX1 gene, the relationships between the New World crocodiles: C. acutus (American crocodile), C. moreletii (Morelet’s crocodile), C. intermedius (Orinoco crocodile) and C. rhombifer (Cuban crocodile) are consistent with other studies using the mitogenome sequence or partial ND6-tRNAglu-cytb, COX1-Cytb (Meganathan et al., 2011; Meredith et al., 2011).
The CaI haplotype in this COX1 region, is equal to the colombian haplotype KF273840 and to GQ144571 and, according to Eaton et al. (2010), tissue from an organism collected in Oaxaca, Mexico (possibly in the Chacahua Lagoon) was used. This haplotype occurs on the coasts of the Pacific and in the Mexican Caribbean; another Colombian haplotype (KF273841, Cac08) also appears in the Mexican clade; in addition, the haplotypes Cac01, Cac3 and Cac05 that correspond to the sequences KF273834, KF273836, KF273838 reported for Colombia by Bloor et al. (2015) form another branch in the clade of C. acutus, as well as the sequence HM636894 reported by Man et al. (2011), which does not indicate the locality of sampling, form another branch in the clade of C. acutus and it is very likely that the crocodile that was sampled comes from Colombia. This supports what was mentioned by Bloor et al. (2015), who consider the possibility of two lineages: The North American haplogroup, where CaI is located and the Central American haplotypes (C. acutus) and Southern american haplogroup. The North American haplogroup, where CaI is located and the Central American haplotypes (C. acutus) and Southern american haplogroup. Unfortunately we were unable to include C. acutus sequences from United State of America, Central America and the Caribbean reported in other studies (Milián-García et al., 2011; Milián-García et al., 2015; Milián-García et al., 2018) since they used another region of the COX1 gene sequence. Recently, studies were published that used fragments of control mtDNA sequences and microsatellites (Milián-García et al., 2020) and SNPs (Rossi et al., 2020), where there is a differentiation between the populations of the Caribbean (Mexico and Belize) and the populations of South America that coincide with those observed in this study using only information from the COX1 gene.
On the other hand, it was observed that the values of genetic diversity are like those obtained by Cedeño-Vázquez et al. (2008), Glenn et al. (2002), Ray et al. (2004) and Rodríguez (2007) using the control region of the mitochondria as well as what was observed in the COX1 sequence studied by Milián-García et al. (2018); these values have been characterized as low by these authors. The haplotype diversity was also like the one reported by Ray et al. (2004); these authors found values of 0.251 for C. moreletii, whereas Rodríguez (2007) found values of 0.518 for C. acutus. Milián-García et al. (2018) found, in C. acutus from captives and wild in Cuba, four haplotypes of the COX1 fragment with a haplotype diversity per population that varies from 0.286-0.558.
The low genetic diversity found with mitochondrial markers may be due to the low mutation rate that they present for being conserved genes or be interpreted as consequence of bottleneck when populations remain reduced in a demographic context. Vasconcelos et al. (2008), suggest that severe declines in crocodiles can be masked in mitochondrial DNA information, first, because previous population expansions may have left signature on mitochondrial DNA; second, signal cannot be detected on mitochondrial because of long generation periods in crocodiles, when compared to 100 years of continuous exploitation; and third, because crocodile counts can still be large, even when they represent only a fraction of ancient abundance. In this sense, population recovery may still be an artifact. Our neutrality test values did not detect significant deviations indicative of non-neutrality; also, we did not find any positive D values that could indicate a positive selection, balancing or a reduction of the population size (Perfectti et al., 2009). We think that a possible hypothesis in this regard is that crocodiles have remained largely unchanged throughout their evolutionary history, causing a low mutation rate in conserved genes as it is mentioned by Field (1988) and Li and Graur (2000). In this regard, Laird et al. (1969), suggested that the rate of nucleotide substitution may be negatively correlated with generational time, that is, species with short generational periods have high substitution rates, this is supported by Kimura (1983).
Particularly in poikilotherm organisms, Mooers and Harvey (1994) found low nucleotide substitution rates in mitochondrial DNA compared to homeoterm organisms of the same size, indicating that the metabolic rate influences the rate of molecular evolution (Allen et al., 2006; Martin & Palumbi, 1993), which could also respond to the low rate of nucleotide substitution found in the crocodiles of this study. Ray et al. (2004) using the mitochondrial DNA control region calculated a nucleotide substitution rate of 1.66 x 10-8, which is 3.6 times lower than the standard rate of 5.9 x 10-8 of mitochondrial genes for hominids calculated by Brown et al. (1982). Recently the mutation rate in the complete genome of crocodiles was estimated at 7.9 x 10-9 (Green et al., 2014) which is very low compared to the rate of nuclear genetic mutation in humans that are 1-10 x 10-5.
The variation in this gene is lower compared to nuclear markers, such as microsatellites; an example of this is the study of Dever et al. (2002), who found heterozygosity values of 0.49 in C. moreletii populations from Belize. Therefore, the highly polymorphic pattern of nuclear markers can detect higher polymorphism in crocodile populations. Flint et al. (2000) and Lawson et al. (1989) mentioned that patterns of low genetic variation in mitochondrial genes found in crocodile populations have been attempted to explain because of directional selection, which favors a particular phenotype in response to long periods of environmental stability due to its adaptation to relatively unchanging aquatic environments.
The analysis of the obtained sequences for both species evidenced individuals phenotypically characterized as C. acutus with C. moreletii haplotypes and viceversa. Specimens from the UMA in Tabasco, Quintana Roo and Oaxaca showed a different phenotype. Regarding the captive individuals from Buenavista, Tabasco, that had haplotypes of another species, at the time of collection, they were part of a breeding female group from a leather and meat marketing crocodile farm in Campeche. It is well known that these females have had offspring previously; therefore, this supports the analyzes of Rodríguez et al. (2008) in which they show that the hybrid individuals grouped in intermediate positions in their phylogenetic reconstructions are backcrosses towards one or the other species, so interspecific hybrid individuals are reproductively viable.
Of the crocodiles that were sampled in Quintana Roo, an individual phenotypically identified as C. moreletii corresponded to the CaII haplotype (characteristic of C. acutus) and two phenotypically C. acutus individuals had haplotype of C. moreletii had CmI y CmVI haplotype (characteristic of C. moreletii); this is very likely since they could correspond to hybrid organisms of C. acutus x C. moreletii. Hybridization within the Crocodylus genus has been reported in several species, mainly in captive organisms (FitzSimmons et al., 2002; Milián-García, et al., 2015; Milián-García, et al., 2018; Rodríguez et al., 2011; Tabora et al., 2012; Weaver et al., 2008). The presence of hybrids in Belize is the result of crosses between C. moreletii males and C. acutus females, which have a fertile offspring, as reported by Hekkala (2004) and Ray et al. (2004). However, in this study, we found that the interbreed between these species has been bidirectional, as also reported by Cedeño-Vázquez et al. (2008); since C. acutus haplotypes were grouped into organisms phenotypically determined as C. moreletii and viceversa.
The Pacific coast of Mexico is the ancestral habitat of C. acutus. It is known that in the last century (1970 approximately) in the lagoon of Chacahua (Oaxaca), a farm of C. morelleti was established started its operations with a batch of 40 animals of C. moreletii from the Atlanta Zoo, but by the scourge of a hurricane the crocodiles escaped towards the lagoon and eventually crossed with the C. acutus (Muñíz et al., 1997).
The results obtained in this study showed that there are reproductive individuals with different haplotypes to those of the species. It is important to take this into account for future reintroductions, since having the genetically identified individuals is paramount for the conservation of these species. The translocation of specimens into non-original distribution zones has also resulted in hybridization processes, such as in the case of the Chacahua Lagoons (Cedeño-Vázquez et al., 2008; Muñíz et al., 1997). Another case reported by Sánchez-Vilchis (2007) corresponded to the Alcazahue lagoon in Colima, where C. moreletii individuals were introduced in 1985 to perform intensive breeding. Currently, it is believed that C. moreletii has displaced C. acutus in this locality, even some C. acutus individuals with atypical nuchal osteoderms patterns have been observed, reason why it is believed that hybridization took place in this area; however, until now the individuals of these populations have not been studied with nuclear markers. Therefore, it is important to determine the actual conservation status of both species and the effect of anthropogenic hybridization on populations before developing exploitation strategies.
Management decisions should be based on the strongest possible evidence to allow a reliable estimate of the genetic processes that are occurring in crocodile populations. This study provides a small but significant advance in the genetic knowledge of both crocodile species and the use of mitochondrial markers, which in this case, the COX1 gene allowed the detection of hybrid organisms in wild and captive populations. Conservation efforts for both crocodile species should prevent interbreeding of both threatened species and would require the genetic identification of pure populations, to design effective conservation strategies considering the possibility of natural hybridization in sympatry areas.
Ethical statement: the authors declare that they all agree with this publication and made significant contributions; that there is no conflict of interest of any kind; and that we followed all pertinent ethical and legal procedures and requirements. All financial sources are fully and clearly stated in the acknowledgements section. A signed document has been filed in the journal archives.