Identification of the taxonomic status of Scinax nebulosus and Scinax constrictus (Scinaxinae, Anura) based on molecular markers

30, 2021 (With 4 figures) Abstract The validation of many anuran species is based on a strictly descriptive, morphological analysis of a small number of specimens with a limited geographic distribution. The Scinax Wagler, 1830 genus is a controversial group with many doubtful taxa and taxonomic uncertainties, due a high number of cryptic species. One example is the pair of species Scinax constrictus and Scinax nebulosus , which share a similar morphology. Scinax constrictus is restricted to the Brazilian Cerrado savanna, while S. nebulosus is widely distributed throughout northern South America. Despite the validation of many anuran species, discriminations based only on morphological traits is quite difficult due to the high conservative morphology of some groups. In this context, the present study uses mitochondrial and nuclear genes to provide a more consistent diagnosis and test the validity of S. constrictus as a distinct species from S. nebulosus , as well as evaluate the position of these taxa within the Scinax genus. The topologies obtained herein uphold the monophyletic status of Scinax based on all molecular markers assessed in this study, in all analytical approaches, with


Introduction
Phylogenetic systematics seeks evidence on the evolutionary relationships among species in order to comprehend diversity patterns (Lourenço et al., 2015). Neotropical amphibians are characterized by considerable phylogeographic structuring and profound levels of genetic divergence, which require extended periods of time for the evolution of phenotypical differences (Sá et al., 2014;Suk et al., 2020).
As a group, Neotropical amphibians present relatively conserved morphological characteristics and, in many cases, reliable diagnostic traits are difficult to identify (Nogueira et al., 2016), which often implies the existence of cryptic species, which have been increasingly identified in recent years (Guarnizo et al., 2015;Pinto et al., 2019). This means that the species-level diversity of many genera has long been underestimated. In this context, molecular analyses can be extremely useful for the identification of cryptic variations in morphologically identical populations (Lima, 2011;Takazone, 2015;Bessa-Silva et al., 2016).
The Hylidae family is one of the most diverse anuran groups, found in Europe, northern Africa and Asia, South and Central America, western India, Australia, and New Guinea (Salducci et al., 2002). The Hylidae group is currently one of the most diverse anuran families, comprising about 701 species, around 14.5% of all anurans. Almost half of these species (338 species in 27 genera) are found in Brazil (Frost, 2017).
In the recently proposed arrangement by Duellman et al. (2016), the Scinax (Wagler, 1830) genus was excluded from the Hylinae subfamily and assigned to a new subfamily, Scinaxinae, which contains 70 described species (Frost, 2017). This analysis also supported a new arrangement for species belonging to the Scinax genus, which were originally categorized into two groups, the Scinax ruber and Scinax catharinae clades. In this new classification, the entire catharinae clade was reclassified into the Ololygon genus (Duellman et al., 2016).
The diagnostic characters traditionally used to identify S. nebulosus and S. constrictus are inadequate for the reliable differentiation of these species, leading to the uncertain taxonomic status of some specimens from populations historically identified as belonging to one of the two species. The increasing use of molecular data has reinforced the often enigmatic nature of amphibian morphology evolution (Menezes et al., 2016;Ferrão et al., 2016), and has revitalized the study of amphibian taxonomy. Many groups of amphibians are morphologically conservative and this, alongside their high convergence degree, has resulted in numerous erroneous interpretations of anuran phylogeny when based solely on morphological characteristics (Bossuyt and Milinkovitch, 2000;Guarnizo et al., 2015).
As morphological characters are not adequate to reliably distinguish these species, the genetic validity of populations diagnosed as S. constrictus or S. nebulosus was tested in the present study using both nuclear and mitochondrial data, as the lack of genetic data for these species reinforces the need for molecular analyses, allowing for assessments concerning the species arrangement within the Scinax genus.

Sample collection and molecular procedures
S. nebulosus and S. constrictus specimens were obtained through donations and exchanges with scientific institutions, and new specimens were also sampled in the field ( Figure 1, Table 1). Their collection was authorized by the Chico Mendes Institute for Biodiversity Conservation (ICMBio), through permanent license number 52593-3. The sampling covers a wide range of localities within the geographic distribution of the species (Figure 1). The samples were complemented by Scinax sequences available in GenBank. It is important to note that this study does not aim to generate threats of extinction.
In the case of the newly collected specimens, muscle tissue samples were extracted from the posterior members and conserved in absolute alcohol in 1.5 mL microtubes. The specimens themselves were fixed in 10% formalin and stored in 70% ethanol and subsequently catalogued and deposited at the Federal University of Maranhão (UFMA) scientific collection.
Pseudis paradoxa, representing a closely-related hylid genus, was used as the outgroup, based on the recent phylogenies published for the genus (Faivovich et al., 2005;Pyron and Wiens, 2011;Duellman et al., 2016). Genomic DNA was isolated from muscle tissue using a phenol-chloroform protocol modified from Sambrook and Russel (2001). As negative amplification of some samples was noted, a Wizard® Genomic DNA Purification kit was applied, using the mouse tail protocol. DNA extraction confirmation and the quantification and quality determination of the extracted material was based on electrophoresis runs in 1% agarose gels.
The sequencing reactions were prepared using the ABI Dye Terminator kit (Applied Biosystems), following the manufacturer's recommendations. Once amplified, the samples were injected into an ABI 3500XL automatic sequencer, resulting in a file with base sequences and the chromatogram.
It is important to note that, while an attempt was made to extract genomic DNA from all the samples and amplify all three markers, in some cases specimen conditions or tissue storage altered the quality of the extracted DNA, which was reflected in the amplification resolution.

Phylogenetic analyses
The DNA sequences were visualized using the Codoncode 3.0 (Codoncode Corporation), and aligned in MAFFT (Katoh and Standley, 2013), followed by visual  inspection to detect possible base identification errors. Variable sites were checked using the MEGA 7.0 software (Kumar et al., 2016), and genetic distance values were also calculated using this program Two different phylogenetic inferences were applied herein, namely Bayesian Inferences (BI), performed using the Mr. Bayes v.3.2 software (Ronquist et al., 2012) and a Maximum Likelihood (ML) analysis, using the PhyML v.3.0 program (Guindon et al., 2010). The most appropriate evolutionary model for each database was identified using jModeltest 2 (Darriba et al., 2012). For the ML analysis, the reliability of each node was estimated using a nonparametric bootstrap procedure comprising 1000 replicates (Felsenstein, 1985). The Bayesian Inferences were based on MCMC sampling in four simultaneous runs, each with four chains (one hot and three heated) with a total run of 10 million generations. The Bayesian probabilities were selected based on a consensus of 60%. As a conservative measure, 25% of the initial trees were discarded as burn-in. The log-likelihood scores were plotted in the Tracer v.1.5 program (Rambaut and Drummond, 2009) to confirm the validity of the initial burn-in period used in Mr. Bayes. The post burn-in samples were used to build the trees used in the present study.
As the number of sequences obtained for the different markers varied, the results were presented in two different ways: (i) mitochondrial gene trees and (ii) a tree combining the 16S and Rhodopsin genes. In the latter case, only S. constrictus and S. nebulosus specimens in which both the 16S and Rhodopsin genes were amplified were included in the analyses.
Inconsistencies among gene trees and in species analyses based on multiple loci are relatively common in evolutionary studies (Nichols, 2001). In the present study, the phylogenetic reconstructions were compared using the BEAST v. 1.6.2 software (Rambaut and Drummond, 2007), in order to evaluate inter-specimen relationships (*BEAST) (Heled and Drummond, 2010). Runs comprising 80 million generations were conducted, with samples being retrieved every 10,000 generations, of which 10% were discarded as burn-in. The runs were analyzed for possible cuts, using Tracer v.1.6, available in the BEAST program (Rambaut and Drummond, 2007), to verify whether the Effective Sample Size (ESS) was adequate, in all cases. All the trees generated in the present study were edited in Figtree v.1.3.1.

Results
The database included 477 base pairs (bps) for the 16S gene, 615 bps for COI, and 295 bps for Rhodopsin. The concatenated 16S+Rhod database contained 772 bps.
The greatest genetic distance was recorded for the COI marker (Table 2), with a divergence of 12.7% between S. nebulosus and S. constrictus. In the case of the 16S gene (Table 3), the divergence was of 4.7%, whereas for the Rhodopsin gene, the genetic distance between the two species was of only 1% (Table 4).
The same evolutionary model was selected by jModelTest 2 for the COI and Rhodopsin genes, with a different model selected for the 16S gene ( Table 5). As the topologies retrieved by both analytical approaches (BI and ML) were highly consistent and very similar to one another, only the BI trees are presented herein. The ML trees can be provided upon request.
Scinax was monophyletic in all topologies derived from all different genetic markers, and in all analytical approaches. Regarding the BI phylogenetic tree obtained for the 16S marker, a number of the Scinax species also formed well-supported monophyletic groups. In this topology, S. nebulosus and S. constrictus form a well-supported clade in which both species tend to form natural groups ( Figure 2).   Poorly-identified specimens (labeled as S. gr. nebulosus) were grouped with either S. nebulosus or S. constrictus. Some structuring was observed for S. nebulosus samples. Specimens from northeastern (Alagoas and Maranhão) and central (Goiás) Brazil were separated with high support from both the BI (0.97) and ML (98%) analyses, as were those from two central Brazilian Cerrado portions, namely Mato Grosso and Tocantins (BI = 1.00, ML = 99%) and Mato Grosso alone (BI = 1.00, ML = 99%). Scinax constrictus appears as a sister species for S. nebulosus in both the BI and ML analyses. Similar findings were obtained in the COI trees (Figure 3), in which S. constrictus and S. nebulosus appear as genetically distinct species, with high support values for the phylogenetic analyses (BI = 0.98, ML = 100%). As in the 16S trees, the COI tree recovered S. nebulosus as a polyphyletic species, formed by at least three groups, with S. constrictus as a sister group. One group contained the specimen from Roraima, the second group encompassed the specimens from Mato Grosso, while the third group included the sequences from Alagoas, Pará, and Maranhão.
The species tree obtained from BEAST also indicates that S. constrictus and S. nebulosus are separate genetic lineages, with high support values (1). The analysis nevertheless confirmed the close phylogenetic relationship between both species (Figure 4).

Discussion
The genetic distances obtained from all three markers indicate that S. constrictus and S. nebulosus are distinct species. The divergence value recorded for the COI gene was of over 10%, a comparatively high value for   differences between amphibian species (Vences et al., 2005a;Guarnizo et al., 2015). This differentiation was also supported by more conservative markers used in other studies (Nair et al., 2012;Nogueira et al., 2016) and in the present study, i.e., 16S and Rhodopsin. Based on an experimental analysis, Vences et al. (2005b) concluded that the mitochondrial 16S rRNA gene satisfies the prerequisites for a universal DNA barcode marker for amphibians and may be superior to COI sequences for the identification of vertebrate clades. This study observed that, while COI varied considerably in amphibians, 16S was highly conserved in these vertebrates. The authors thus recommend the use of 16S as an additional marker for vertebrate DNA barcoding. Xia et al. (2012) also support the use of the 16S rRNA gene for amphibian DNA barcoding, given that the authors indicated that genetic distances based on the Kimura 2-parameter (K2P) model revealed that the mean intraspecific distance was 1.4% in the COI gene, but only 0.3% in the 16S rRNA. The level of 16S rRNA variation obtained in the present study was higher than this aforementioned value, indicating a high degree of genetic variation among the study species, reinforcing their taxonomic validity.
Since the description of S. constrictus , many specimens of both species have been difficult to classify reliably, given the lack of precision in the traditional morphological approach for taxa differentiation. In general, taxa with uncertain morphological variations are considered to be species complexes, with cryptic forms. This is probably the case for both S. nebulosus and S. constrictus, which present high levels of genetic differentiation in morphologically similar forms. In such   cases, the molecular approach is essential for reliable species differentiation, as in the case of the Rhinella marina species complex, comprising distinct lineages composed of cryptic species (Bessa-Silva et al., 2016).
In addition to confirming Scinax as a monophyletic group, the topologies obtained herein support the validation of S. constrictus and S. nebulosus as distinct species. Some aspects of the arrangement within some of the clades were unexpected, however, which could invalidate the species as a monophyletic group. Most of the incongruent samples were specimens recently assigned to the S. nebulosus group, probably as a consequence of the unreliable morphological definition of the species. Scinax constrictus can be distinguished from S. nebulosus by its smaller size, the presence of a discrete anal flap in males (absent in S. nebulosus), and the presence of a row of small tubercles on the lower jaw margin, which is absent or indistinct in S. nebulosus . The advertisement call of S. constrictus is also quite distinct concerning dominant frequency and number of notes , even though these traits vary considerably in S. nebulosus. Despite these diagnoses, the analysis of specimens diagnosed previously as S. nebulosus or S. constrictus indicate that these traits are not adequate for reliable species distinction (Weber, L.N. personal communication).
It is, therefore, clear, that morphological analyses are extremely vulnerable to error and taxa misidentification, made apparent by molecular tools. The geographic distributions of both taxa appear to support the molecular diagnosis, given that S. constrictus appears to be restricted to the Cerrado savanna of central Brazil, whereas S. nebulosus is more widely distributed, and may even be sympatric with S. constrictus.
Genetic analyses have been widely applied to validate anuran species described on the basis of morphological characteristics. Bruschi et al. (2013), for example, evaluated the taxonomic status of Phyllomedusa hypochondrialis (Anura, Hylidae, Phyllomedusinae) applying cytogenetic and molecular data to differentiate this species from Phyllomedusa azurea. Their results contributed to the extension of the geographic range of P. hypochondrialis, and confirmed its sympatry with P. azurea, given that incorrect specimen identification appeared to have resulted in inconsistencies in the original P. azurea geographic distribution definition. In the present study, the analysis of the 16S gene ( Figure 2) indicated the grouping of S. nebulosus and S. constrictus specimens with S. constrictus specimens originally assigned to the nebulosus group, which appears to have been the result of species misidentification. In the case of both Phyllomedusa Bruschi et al. (2013) and Scinax (present study), molecular analyses indicated clearly a monophyletic position within the Hylidae family. Brusquetti et al. (2014) reviewed the taxonomy of Scinax fuscomarginatus (Lutz, 1925), and evaluated morphologically similar species exhibiting controversial taxonomy, such as S. fuscumarginatus, S. parkeri, S. trilineatus, S. lutzorum and S. pusillus, by molecular analyses, with paraphyletic clades being recuperated for S. fuscomarginatus in comparison to the other species, similar to the paraphyletism between S. nebulosus and S. constrictus reported in the present study (Figures 2 and 3). Duellman et al. (2016) subsequently redistributed Scinax species, and among the species analyzed by Brusquetti et al. (2014), only S. fuscomarginatus was maintained in the genus.
In their phylogenetic and biogeographic analyses of the Hylidae family, Duellman et al. (2016) proposed major changes in the Scinax genus arrangement, ranging from a new subfamily to species redistribution. Scinax nebulosus and S. constrictus remained in the genus, and the results of the present study confirm their status as valid Scinax genus lineages.
A number of the species described in the Scinax Wagler (Anura: Hylidae) genus require genetic validation, and the application of molecular markers has not only contributed to the identification of new species, but also in the differentiation of cryptic lineages (Nunes et al., 2012;Lourenço et al., 2013;Sturaro and Peloso, 2014;Nogueira et al., 2016).
The molecular phylogeny of Fouquet et al. (2007) indicates that Scinax ruber, found in the French Guiana region, is represented by six divergent lineages. The molecular phylogeny obtained by Lourenço et al. (2013) confirmed that the Scinax specimens collected at the Chapada dos Guimarães National Park, in central Brazil, which are morphologically similar to species from the S. catharinae group (sensu Faivovich et al., 2005), do, in fact, represent a new species. Nogueira et al. (2016) emphasize the need for further comparative analyses between these species, integrating morphological and bioacoustics data.
In this context, many researchers have highlighted the importance of a multiple approach, i.e., integrating molecular, morphological, bioacoustic, and ecological data, in taxonomic anuran assessments, in facilitating species identification (Padial et al., 2014;Silva, 2017;Pinto et al., 2019). Overall, molecular analyses aided in differentiating populations diagnosed as S. constrictus and S. nebulosus, despite the morphological similarities of both species.