Digital Commons @ Trinity Digital Commons @ Trinity

: Turtles are one of the most threatened groups of vertebrates, with about 60% of species classiﬁed at some level of extinction risk. Compounding this extinction crisis are cryptic species and species complexes that are evaluated under a single species epithet but harbor multiple species, each of which needs to be evaluated independently. The Phrynops geoffroanus species group is a classic example. Described ﬁrst in 1812, it is currently thought to harbor multiple species. To test this hypothesis, we collected mitochondrial and nuclear genomic data, morphometric data, and distribution and associated biome information. We applied statistically rigorous species delimitation analyses, taxonomic hypotheses tests, and fully coalescent phylogenetic reconstruction methods, concluding that the Phrynops geoffroanus species complex comprises four geographically structured species/lineages that diverged during the Pleistocene and are currently geographically structured along the main South American river basins and biomes. These species/lineages show subtle and largely non-signiﬁcant differences in shape but are characterized by differences in coloration and patterns of marks on the head and plastron. Our results contribute to the understanding of species diversity and diversiﬁcation of biodiversity in South America and provide an important basis for the conservation of freshwater turtles.


Introduction
It is becoming ever more apparent that many Neotropical taxonomic groups harbor cryptic diversity, i.e., evolutionarily divergent lineages that do not present obvious or any morphological differences, e.g., the works of [1][2][3][4][5][6]. Coupled with this issue is the conceptualization of species as evolving lineages [7], which often are easily identified through the analysis of molecular data, and current taxonomic practices that treat species as immutable categories defined by morphological diagnoses, the hallmark of Linnaean taxonomy. Yet it is only species as evolutionary lineages that fully capture the evolutionary history and species diversity of any given taxonomic group. In some taxonomic groups, the levels of cryptic diversity may be high enough to cause a "taxonomic crisis" since traditional taxonomy often does not reflect the diversity and distribution of species [5,8]. Consequently, the non-inclusion of cryptic diversity causes biases in biodiversity studies [9,10]. Considering that ecological and biogeographic models are fundamental bases for conservation actions and policies and that biodiversity is being lost faster than scientists can describe new species, cryptic diversity has also reduced the efficiency of conservation programs [11].
Molecular markers have been used efficiently to quantify lineage diversity and refine taxonomy, even in organisms showing relatively little morphological diversity, including certain South American and Australasian freshwater turtle groups [4][5][6]12,13]. The analyses of molecular data also permit powerful inferences about the evolutionary history of organisms [1,2,12,[14][15][16][17], adding an additional dimension to traditional biodiversity studies [18]. A plethora of methods for obtaining and analyzing molecular data have been developed in the last decades, and mitochondrial markers used for DNA barcoding have been shown to be cost-effective for phylogenetic analyses focused on cryptic diversity [19]. Mitochondrial markers are particularly efficient in identifying cryptic diversity when analyzed under movement-or coalescent-based species delimitation methods, which objectively delimit boundaries between evolutionarily distinct lineages in non-arbitrary ways, e.g., the works of [20,21].
Ecological and biogeographic models biased by cryptic diversity are particularly problematic when focused on threatened organisms because different lineages may have different environmental requirements and sources of threat [22]. Turtles, along with primates, are the two most threatened groups of vertebrates, with about 60% of species classified at some level of extinction risk [23,24]. Several factors have been described as general sources of population declines, including habitat loss and environmental degradation, pet and food trade, traditional medicines, and climate change [25][26][27][28].
The intensity of these threats is expected to vary over the distribution of a species as a consequence of spatial heterogeneity of environment and sources of threat [2,29]. Yet broadly distributed species often are complexes of cryptic species that occur across a range of environments and geomorphological landscapes [1,2,6,29]. Due to the large geographic extent of these "species", they tend to be evaluated as not threatened or in low-risk categories, yet clearly, conservation status needs to be assessed for each cryptic species [30][31][32].
Previous studies have suggested that Phrynops geoffroanus comprises a complex of cryptic species [5,13,24,30,[32][33][34][35]. What are these cryptic species, where do they occur, and what are their phylogenetic relationships are, however, uncertain given the lack of sampling. It is likely, however, that these cryptic species will be structured along major river basins since this group depends on the drainage network to disperse.
In this study, we analyzed a geographically comprehensive sample of the Geoffroy's side-necked turtle Phrynops geoffroanus (Chelidae) species complex spanning all known major South American river basins in which it is known to occur. For all collected specimens, we collected mitochondrial DNA, genomic ddRAD, and morphometric data, which we then used to test two principal hypotheses: (1) Phrynops geoffroanus is a species complex comprised of multiple evolutionary lineages. (2) The distribution of these cryptic species coincides with principal South American river basins. Additionally, we provide new insights into the phylogenetic relationships of Phrynops and other Chelidae. Molecular data. We extracted mitochondrial DNA from striated muscle tissue using the phenol-chloroform protocol [36]. DNA extracts were used to generate sequences of three mitochondrial genes and a genomic library from which haplotypes and single nucleotide polymorphisms were extracted. We amplified partial fragments of the genes 16S rDNA, Cytochrome b (Cytb), and Cytochrome Oxidase I (COI) using a respective combination of the primers 16Sar(L) and 16Sbr(H) [37], L14725 [38], H15573 [39] and LturtCOIa, HturtCOIa [40]. The final volume of the PCR mixture was 12 µL and contained 4.8 µL ddH 2  We purified the PCR products using ExoSap, following the manufacturer's instructions (ThermoFisher, Waltham, MA, USA). Sequencing reactions for the three genes were carried out using the Big Dye Terminator kit following the directions of the manufacturer (ThermoFisher, Waltham, MA, USA). We used the amplification primers for bidirectional sequencing at 55 • C, followed by precipitating the product of the sequencing reaction with EtOH and EDTA. The precipitated PCR products were resuspended in 10 µL of formamide and sequenced using ABI 3500 automatic sequencer (ThermoFisher, Waltham, MA, USA). The generated sequences were edited in Geneious 8.1.8 [41], aligned in Clustal W [42], and concatenated in CodonCode Aligner v.3.5.2 (http://www.codoncode.com/ aligner/download.htm, accessed on 29 February 2018), resulting in a concatenated matrix of 1843 bp. The sequence data were deposited in GenBank under the following accession numbers: ON063048-ON063118 (COI), ON061103-ON061173 (Cytb), ON063137-ON063207 (16S); SUB11224878 (ddRAD).
We prepared the genomic library using a modified ddRADseq protocol [43] adapted for IonTorrentPGM (https://github.com/legalLab/protocols-scripts, accessed on 10 January 2016) for one individual per sampled locality plus other species of the genus (P. tuberosus and P. williamsi) and outgroup taxa (Mesoclemmys gibba, Ranacephala hogei, and Platemys platycephala). Briefly, after DNA extraction, the 200 ng of whole genomic DNA was digested with SdaI and Csp6I restriction enzymes, and barcoded sequencing adapters were added with T4 ligase. Following PCR, purification with AMPure Beads (Beckman Coulter, Brea, CA, USA) and fragment size selection (320-400 bp) on the Pippin Prep electrophoresis platform (Sage Science, Beverly, MA, USA), all samples were pooled in equimolar quantities and sequenced in the Genetic Sequencer Generation Ion Torrent PGM, using reagents and sequencing kits 400 pb and chip 318.
After demultiplexing our run, we added three representatives each of Chelus fimbriata and C. orinocensis generated for a previous study [6]; these individuals were outgroups and were used for calibrating divergence times of our ingroup taxa. We processed our demultiplexed runs in PyRAD [44] to generate data matrices for phylogenetic analyses. We filtered our reads on quality (nucleotides with PHRED scores <30 were changed to ambiguities, reads with >3 ambiguities were discarded) and coverage (a minimum 5× coverage per read). In the final data matrix, we included only those loci that were present in at least 75% of the individuals, resulting in a final alignment of 244 loci representing 83,703 bp. We also generated a matrix of SNPs using DiscoSnp-RAD [45]. We extracted single nucleotide polymorphisms (SNPs) from our reads using a minimum read depth of 5. The highest quality SNPs were then sampled from each locus, and the SNPs were further filtered on quality. We retained only those SNPs with a rank >0.9, a statistic incorporating the discriminant power and read coverage of each SNP, and those that were present in at least 80% of the samples. This resulted in 7668 SNPs and a matrix with 2.01% missing data.
Phylogenetic reconstruction-ddRAD. To infer phylogenetic relationships among individuals of the P. geoffroanus species complex, we ran a coalescent-based phylogenetic reconstruction using SVDQuartets [46] as implemented in PAUP* [47] using the SNP matrix as input. The result of this analysis was used in species hypotheses testing.
Single-locus species discovery analysis. We reconstructed the Bayesian inference mtDNA phylogeny of P. geoffroanus in the software BEAST v2.6 [48]. Markov chain Monte Carlo (MCMC) searches were made for 10 million generations, sampling every 1000th generation for a total of 10 000 trees using the coalescent constant population size and coalescent tree prior to implementing the GTR + G + I model of molecular evolution; the most likely model indicated in jModeltest [49]. We verified convergence of the Markov chain convergence using Tracer v1.7 [50] and summarized the topologies in TreeAnnotator v1.6.2 [48] using the maximum clade credibility criterion after discarding the first 10% of the trees as burn-in.
We subjected the resulting topology to four single-locus species discovery (SLSD) methods to partition our data set into putative species-like clusters and to assess the support for these clusters following Machado et al. [3]: GMYC, the general mixed Yule coalescent model; bGMYC, a Bayesian implementation of the GMYC [51]; mPTP, the Poisson tree process method [52]; and local minima (locMin), a distance threshold optimizing and clustering approach from the spider_1.3-0 software package [53].
For GMYC and bGMYC, we used splits_1.0-19 [51] and bGMYC_1.0.2 [54] packages in R [55]. We summarized the bGMYC posterior samples into putative species with a conservative posterior probability of conspecificity at 0.05. The locMin analyses were conducted as a point estimate and on a set of 1000 bootstrapped data sets to generate a confidence interval. For mPTP, the BEAST chronograms (ultrametric trees with branch lengths scaled by time) were transformed into phylograms using maximum likelihood optimization in phangorn_2.2.0 [56] under the GTR + G + I substitution model. We retained only those lineages that were discovered by at least three of the four methods. Finally, we calculated uncorrected p-distance among lineages using the concatenated sequences in MEGA 7.0 [57].
Species tree analysis-ddRAD. We ran PartitionFinder2 [59] to select an optimal number of partitions and their respective models of molecular evolution in our genomic alignment, resulting in 86 partitions. We then generated an XML file for species tree analyses where species were clusters of individuals delimited as such in the species hypotheses tests. We analyzed the data in BEAST v2.6 [48], implementing the calibrated Yule tree prior, generating 5 × 10 8 topologies, sampling every 5000th topology, and discarding the first 15% topologies as burn-in after verifying stationarity in Tracer v1.7 [50]. Topologies were visualized in DensiTree [48], and a final consensus maximum clade credibility tree was generated using TreeAnnotator v1.6.2 [48] after discarding 15% burn-in. The maximum clade credibility tree and divergence times calibrated using the Chelidae divergence prior (mean = 36.0 My; 95% HPD 26.5-46.8 My; Thomson et al. [60]) were then visualized in FigTree v1.3 [61].
Morphometric analyses. We collected 34 linear measurements from 198 individuals (127 individuals in addition to those used in the molecular analyses) of the Phrynops geoffroanus complex (Supplemental Table S2). In addition to linear measurements, we also recorded the sex, basin, and biome in which the individuals occurred, and posterior to the SLSD analysis, the lineage to which each individual belonged. Due to the undesirable properties of proportions in statistical analyses [62], we performed linear regression against total carapace length and retained residuals of each of the remaining 33 linear measurements; the residuals largely encompass components of shape. The residuals did not deviate from normality, so we used parametric analyses. After filtering all juveniles or individuals of undetermined sex, we retained 151 individuals, 73 females and 78 males. We tested for differences in shape between males and females. Since there are significant differences in shape between males and females, we analyzed males and females separately in subsequent analyses. We reduced the dimensionality of our data by performing principal component analysis and retaining the first 15 principal components; the first 15 principal components explain >90% of the total variance. We used ANOVA to test for shape differences among lineages among basins and among biomes. We subsequently applied the discriminant analysis of principal components-DAPC [63] to each data set to test if individuals of the different lineages could be discriminated morphologically. In the DAPC analysis groups (the cryptic species identified in the single-locus delimitation analyses) were identified a priory, and the first 15 principal components were used in the discriminant analyses.
Distribution analyses. Finally, we tested the association between lineages and their distribution classified as major river basin and biome using Chi-square and Fisher's exact test.
Genomic species trees. We reconstructed two dated species trees: the four species hypothesis (Figure 3) and the five species hypothesis (Figure 4). The topology of both phylogenies was congruent, as were divergence times. Phrynops williamsi diverged from the P. geoffroanus complex~2.7 ma, and Phrynops tuberosus was nested within the P. geoffroanus complex. The crown age of the P. geoffroanus complex was 1.2 ma representing the divergence of the Atlantic Rainforest lineage of P. geoffroanus complex (L4) and all other lineages of the P. geoffroanus complex, including P. tuberosus, which occurs in the Amazonian savannas of Roraima. In the four species phylogeny, P. tuberosus diverged from all remaining lineages of P. geoffroanus at~710 ka, and the Paraná lineage of P. geoffroanus diverged from the Amazon + São Francisco lineage of P. geoffroanus at~640 ka; although the entire clade was highly supported (pp = 0.95), the sister taxon relationship of the Paraná and Amazon + São Francisco lineages was only weakly supported (pp = 0.57). In the five species phylogeny, P. tuberosus diverged from all remaining lineages of P. geoffroanus at 1 ma, the Paraná lineage of P. geoffroanus diverged from the Amazon + São Francisco lineage of P. geoffroanus at~790 ka, and the Amazon and São Francisco lineages diverged at~730 ka of P. geoffroanus. The entire clade was weakly supported (pp = 0.70), and the sister relationships of the Paraná, Amazon, and São Francisco lineages were also weakly supported (pp = 0.92 and 0.90). Table 2. Test of taxonomic hypotheses using path sampling in BEAST2. Differences in marginal likelihoods were evaluated using Bayes factors. The four species hypothesis (P. tuberosus, P. geoffroanus L4, P. geoffroanus L1a, P. geoffroanus L1b + L2 + L3) was supported as the best taxonomic hypothesis, and the five species hypothesis (P. tuberosus, P. geoffroanus L4, P. geoffroanus L1a, P. geoffroanus L3, P. geoffroanus L1a + L2) as the second-best taxonomic hypothesis.

Hypothesis
Marginal Genomic species trees. We reconstructed two dated species trees: the four species hypothesis (Figure 3) and the five species hypothesis (Figure 4). The topology of both phylogenies was congruent, as were divergence times. Phrynops williamsi diverged from the P. geoffroanus complex ~2.7 ma, and Phrynops tuberosus was nested within the P. geoffroanus complex. The crown age of the P. geoffroanus complex was 1.2 ma representing the divergence of the Atlantic Rainforest lineage of P. geoffroanus complex (L4) and all other lineages of the P. geoffroanus complex, including P. tuberosus, which occurs in the Amazonian savannas of Roraima. In the four species phylogeny, P. tuberosus diverged from all remaining lineages of P. geoffroanus at ~710 ka, and the Paraná lineage of P. geoffroanus diverged from the Amazon + São Francisco lineage of P. geoffroanus at ~640 ka; although the entire clade was highly supported (pp = 0.95), the sister taxon relationship of the Paraná and Amazon + São Francisco lineages was only weakly supported (pp = 0.57). In the five species phylogeny, P. tuberosus diverged from all remaining lineages of P. geoffroanus at ~1 ma, the Paraná lineage of P. geoffroanus diverged from the Amazon + São Francisco lineage of P. geoffroanus at ~790 ka, and the Amazon and São Francisco lineages diverged at ~730 ka of P. geoffroanus. The entire clade was weakly supported (pp = 0.70), and the sister relationships of the Paraná, Amazon, and São Francisco lineages were also weakly supported (pp = 0.92 and 0.90).   Morphometric analyses. The shape of males and females is significantly different (one-way MANOVA, F(1149) = 13.75, p = 0.00030). However, only 2 of the 33 linear measurements are significantly different between males and females, namely plastron width (t-test Holm's p = 0.00196) and plastron length (t-test Holm's p = 0.02425). The size of adult males and females is also significantly different, with females being approximately 7% larger than males (t-test t(148.75) = −2.6384, p = 0.00922). All subsequent analyses with morphometric data were therefore carried out separately for each sex.
We investigated which variables explain differences in shapes of each sex. Of the three variables analyzed (lineage, basin and biome) lineage and basin were significant for males (lineage: one-way MANOVA, F(74) = 5.911, p = 0.00113; basin: one-way MANOVA, F(71) = 3.332, p = 0.00600) and basin and biome were significant for females (basin: oneway MANOVA, F(66) = 2.604, p = 0.02510; biome: one-way MANOVA, F(67) = 3.733, p = 0.00485). Tukey's honest significant differences test indicated this was caused by differences in the shape of males from the southern Atlantic forest basin where lineage L4 occurs when compared to males from other basins. In the case of females, Tukey's honest significant differences test indicated differences in female shape were between the Caatinga and other biomes, but not the Amazon Forest biome; the significant difference by basin was accounted for by a single comparison between the Paraná and the northeastern Atlantic forest drainages. Significant results were observed in analyses of PCA data but not the residuals. Differences in the shape of both males and females are clearly very subtle to non-existent.
Discriminant analysis of principal components of males easily discriminated individuals of lineage L4 from the southern Atlantic forest, the first diverging lineage of the Phrynops geoffroanus species complex, from all others ( Figure 5A). The other lineages could Morphometric analyses. The shape of males and females is significantly different (one-way MANOVA, F(1149) = 13.75, p = 0.00030). However, only 2 of the 33 linear measurements are significantly different between males and females, namely plastron width (t-test Holm's p = 0.00196) and plastron length (t-test Holm's p = 0.02425). The size of adult males and females is also significantly different, with females being approximately 7% larger than males (t-test t(148.75) = −2.6384, p = 0.00922). All subsequent analyses with morphometric data were therefore carried out separately for each sex.
We investigated which variables explain differences in shapes of each sex. Of the three variables analyzed (lineage, basin and biome) lineage and basin were significant for males (lineage: one-way MANOVA, F(74) = 5.911, p = 0.00113; basin: one-way MANOVA, F(71) = 3.332, p = 0.00600) and basin and biome were significant for females (basin: one-way MANOVA, F(66) = 2.604, p = 0.02510; biome: one-way MANOVA, F(67) = 3.733, p = 0.00485). Tukey's honest significant differences test indicated this was caused by differences in the shape of males from the southern Atlantic forest basin where lineage L4 occurs when compared to males from other basins. In the case of females, Tukey's honest significant differences test indicated differences in female shape were between the Caatinga and other biomes, but not the Amazon Forest biome; the significant difference by basin was accounted for by a single comparison between the Paraná and the northeastern Atlantic forest drainages. Significant results were observed in analyses of PCA data but not the residuals. Differences in the shape of both males and females are clearly very subtle to non-existent.
Discriminant analysis of principal components of males easily discriminated individuals of lineage L4 from the southern Atlantic forest, the first diverging lineage of the Phrynops geoffroanus species complex, from all others ( Figure 5A). The other lineages could not be discriminated, however. For females, DAPC was able to largely unable to discriminate individuals ( Figure 5B), with individuals of three of the four lineages overlapping in morphospace, while individuals of lineage L4 from the southern Atlantic forest were somewhat divergent, although this inference was hampered by low sample size. ty 2022, 13, x FOR PEER REVIEW not be discriminated, however. For females, DAPC was able to largely unable to dis inate individuals ( Figure 5B), with individuals of three of the four lineages overlappi morphospace, while individuals of lineage L4 from the southern Atlantic forest somewhat divergent, although this inference was hampered by low sample size.

Discussion
Testing taxonomic hypotheses using genomic data, we were able to provide sup for four or potentially five species/lineages within the Phrynops geoffroanus species plex: Phrynops tuberosus from northern Amazonian savannas of Roraima, one specie tributed in the Paraná-Paraguai basin (L1a), one species from southern and southea savannas of the Amazon basin (L1 + L2), one in the northeastern Atlantic forest-São cisco basin (L3) and one in the southeastern Atlantic forest (L4). In the case of the species hypothesis, animals from southern and southeastern savannas of the Amazo sin (L1 + L2) and the northeastern Atlantic forest-São Francisco basin (L3) would com just one taxon. Our results were consistent with the observation that many widely di uted Amazonian "species" are multiple geographically structured lineages/cryptic sp that hide under the same epithet. This idea has been widely supported based on mole evidence from terrestrial animals, e.g., the works of [1,2,14,64,65], and aquatic fauna the works of [3,6,66].
Based on the dated species tree, these species/lineages diverged between 1.2 an ma. These species/lineages are restricted to specific river basins and biomes (the two b largely correlated in our analysis) but are morphologically cryptic. While there are s icant morphological differences between males and females, there were no statist Distribution analyses. Individuals of the four lineages of the Phrynops geoffroanus species complex are not randomly distributed among major river basins (Pearson's Xsquared = 440.67, df = 18, p < 2.2 × 10 −16 ; Fisher's exact test p < 2.2 × 10 −16 ) and biomes (Pearson's X-squared = 300.89, df = 15, p < 2.2 × 10 −16 ; Fisher's exact test p < 2.2 × 10 −16 ) with major river basins being restricted to biomes.

Discussion
Testing taxonomic hypotheses using genomic data, we were able to provide support for four or potentially five species/lineages within the Phrynops geoffroanus species complex: Phrynops tuberosus from northern Amazonian savannas of Roraima, one species distributed in the Paraná-Paraguai basin (L1a), one species from southern and southeastern savannas of the Amazon basin (L1 + L2), one in the northeastern Atlantic forest-São Francisco basin (L3) and one in the southeastern Atlantic forest (L4). In the case of the four species hypothesis, animals from southern and southeastern savannas of the Amazon basin (L1 + L2) and the northeastern Atlantic forest-São Francisco basin (L3) would comprise just one taxon. Our results were consistent with the observation that many widely distributed Amazonian "species" are multiple geographically structured lineages/cryptic species that hide under the same epithet. This idea has been widely supported based on molecular evidence from terrestrial animals, e.g., the works of [1,2,14,64,65], and aquatic fauna, e.g., the works of [3,6,66].
Based on the dated species tree, these species/lineages diverged between 1.2 and 0.7 ma. These species/lineages are restricted to specific river basins and biomes (the two being largely correlated in our analysis) but are morphologically cryptic. While there are significant morphological differences between males and females, there were no statistically significant differences among the morphology of females and males of the different species, except for males from the southern Atlantic forest. In this region occurs species/lineage L4, a species/lineage that is sister to all other species/lineages of P. geoffroanus and P. tuberosus, from which it has diverged approximately 1.2 ma. Phrynops tuberosus and all other species/lineages of P. geoffroanus diverged at approximately 0.7 ma. The estimated divergence times are compatible with those of Thomson et al. [60], who estimated the divergence of P. tuberosus and P. geoffroanus at 1 ma.
An apparently surprising result of our analysis is the phylogenetic position of P. tuberosus. Phrynops tuberosus is generally considered the sister taxon of P. geoffroanus [60]; however, in our genomic analyses, it is sister to the clade comprising all species/lineages of P. geoffroanus except the southern Atlantic forest species/lineage. Males of the southern Atlantic forest species/lineage of P. geoffroanus are morphologically differentiated in the form of their plastron from other species/lineages of P. geoffroanus in addition to being the first diverging species of the P. geoffroanus species complex, and therefore likely represent a true case of mistaken identity, i.e., these individuals should not have been attributed to P. geoffroanus. Additionally, in the five species taxonomic hypothesis, the sister taxon relationship of P. tuberosus and all species/lineages except the southern Atlantic forest species/lineage was weakly supported (pp = 0.70), while in the four species taxonomic hypothesis, the posterior probability for this sister taxon relationship was 0.95. Similarly, in the SVDQuartet phylogeny (Figure 2), relationships among the species/lineages are also not well supported. In the mtDNA phylogeny (data not shown), P. tuberosus is a sister taxon to all species/lineages of the P. geoffroanus species complex. In conclusion, although there is strong support for the existence of multiple species/lineages in the P. geoffroanus complex, including P. tuberosus, phylogenetic relationships among these species/lineages are not well supported likely due to their relatively rapid and recent divergence.
Our data clearly show relatively recent divergences, from at least 1.2 to 0.7 ma, and very little divergence in morphology. The diversification of the P. geoffroanus species group is largely consistent with the evolutionary history of Chelidae in South America [6,33,[67][68][69] and Neotropical Squamata [70][71][72], suggesting the same set of events that similarly affected the evolutionary history of multiple taxonomic groups.
Morphologically cryptic diversification is also evident in many taxa that have divergence in the Late Miocene and the Pliocene. Possible explanations are simply just too little time and/or selective pressure for taxa that have speciated allopatrically but have not diverged ecologically to have accumulated any relevant morphological differences suitable for species-level diagnoses. The Neotropics in general and the Amazon specifically harbor large numbers of cryptic species, e.g., the works of [2,6,29,65,66,69,[73][74][75]. Oftentimes, however, once discovered through phylogenetic and/or species delimitation methods, many of these taxa can be diagnosed once the search for diagnostic morphological characters is guided by molecular evidence. In the case of the P. geoffroanus species group, there really are minimal differences in the shape of individuals of the different cryptic species, even when a large number of individuals were analyzed and robust statistical methods were employed. However, there are differences in non-morphometric characters, including the color of the carapace, plastron, and lateral bands and stripes on the head ( Figure 6). These findings are relevant to resolving the taxonomic uncertainties associated with the epithet P. geoffroanus. Historically this epithet has been applied to almost any and all individuals of the genus Phrynops from almost anywhere in South America [32,76]. Phrynops geoffroanus is also replete with junior synonyms that, although currently considered invalid, may, in fact, be valid taxa.
Rhodin and Mittermeier [34] were the first to investigate the taxonomic relationships within P. geoffroanus and to propose that P. geoffroanus is, in fact, a species complex, which has been accepted by other studies since then [24,30,32,73]. Our study supports this hypothesis for the first time based on broad geographic sampling, which revealed distinct and spatially structured species/lineages within P. geoffroanus. We provide clear evidence supporting four species/lineages that occupy distinct major river basins and biomes in South America. In addition, while taxonomic hypothesis tests support the existence of only one species/lineage shared between the Amazon and São Francisco basins, the second-best taxonomic hypothesis supports two species/lineages, one occupying the Amazon and one in the São Francisco basin. In addition, although not statistically significant, there are differences in shape and in color patterns of individuals from these two river basins. Therefore, animals from these two river basins may be in the early stages of divergence. We, therefore, propose that all these species/lineages be considered for conservation and decision-making purposes, although formal taxonomic revision needs to await a future study. As mentioned above, there are differences in color and pattern among the different species/lineages, and we are conducting a taxonomic revision of this group. However, any diagnostic characteristics are likely to be quite subtle, and therefore in situ species/lineage conservation efforts should focus less on identifying individuals and more on conserving individuals within a basin/biome. Rhodin and Mittermeier [34] were the first to investigate the taxonomic relationships within P. geoffroanus and to propose that P. geoffroanus is, in fact, a species complex, which has been accepted by other studies since then [24,30,32,73]. Our study supports this hypothesis for the first time based on broad geographic sampling, which revealed distinct and spatially structured species/lineages within P. geoffroanus. We provide clear evidence supporting four species/lineages that occupy distinct major river basins and biomes in South America. In addition, while taxonomic hypothesis tests support the existence of only one species/lineage shared between the Amazon and São Francisco basins, the second-best taxonomic hypothesis supports two species/lineages, one occupying the Amazon and one in the São Francisco basin. In addition, although not statistically significant, there are differences in shape and in color patterns of individuals from these two river basins. Therefore, animals from these two river basins may be in the early stages of divergence. We, therefore, propose that all these species/lineages be considered for conservation and decision-making purposes, although formal taxonomic revision needs to await a future study. As mentioned above, there are differences in color and pattern among the different species/lineages, and we are conducting a taxonomic revision of this group. However, any diagnostic characteristics are likely to be quite subtle, and therefore in situ species/lineage conservation efforts should focus less on identifying individuals and more on conserving individuals within a basin/biome.
Our results are relevant to future conservation decisions because decision making and policies are largely dependent on species and ecological and biogeographic models using these species [6,29,65,69]. Because reptiles have declined globally [24,27] and turtles are among the most endangered vertebrates [23,24], resolving the taxonomy of "P. geoffroanus" and other widely distributed taxa is particularly urgent. The IUCN has not yet evaluated P. geoffroanus, largely due to its large geographic distribution, but it is clear that "P. geoffroanus" is an epithet that encompasses multiple apparently cryptic species. These cryptic species are allopatrically distributed and occupy different portions of Our results are relevant to future conservation decisions because decision making and policies are largely dependent on species and ecological and biogeographic models using these species [6,29,65,69]. Because reptiles have declined globally [24,27] and turtles are among the most endangered vertebrates [23,24], resolving the taxonomy of "P. geoffroanus" and other widely distributed taxa is particularly urgent. The IUCN has not yet evaluated P. geoffroanus, largely due to its large geographic distribution, but it is clear that "P. geoffroanus" is an epithet that encompasses multiple apparently cryptic species. These cryptic species are allopatrically distributed and occupy different portions of anthropogenic disturbance gradients (e.g., human density, pollution, deforestation). Therefore, we argue that a single conservation approach considering just one species of "P. geoffroanus" would be minimally ineffective if not detrimental.