Lowland Leopard Frogs (Rana yavapaiensis) have experienced extensive population declines over the last century. In California, this species was historically known to occur in scattered localities in the extreme southeastern portion of the state, but it has not been positively documented since 1965. Subsequent to this decline in California, nonnative Rio Grande Leopard Frogs (R. berlandieri) have expanded into localities previously occupied by R. yavapaiensis. The lack of extensive formal surveys and the difficulty distinguishing between these species using morphological characters have caused uncertainty about whether Lowland Leopard Frogs persist within their historical range in California. Recently, leopard frogs that could not be confidently identified to species have been observed at historical localities of R. yavapaiensis. Thus, we undertook a formal study of these populations to characterize their morphological and genetic variation, and conclusively determine to which species they belong. Our genetic analyses demonstrate that these frogs are R. berlandieri, but the morphological characters typically used to diagnose these species are largely overlapping. Further complicating field identifications, for some morphological characters, the California R. berlandieri are more similar to R. yavapaiensis than to native-range R. berlandieri. Additionally, invasive R. berlandieri show greater variation in a key character—the condition of the inset dorsolateral folds—than that found across much of the species' native range. These results demonstrate the potential for morphological change during rapid population expansions to confound species identifications. Our findings have implications for future efforts to resolve the status of R. yavapaiensis in California and to identify other native leopard frogs found within the expanding range of R. berlandieri. Our results also highlight the utility of genetic approaches for reliably identifying morphologically similar leopard frogs.

OVER the past century, ranid frog populations in the American Southwest have experienced tremendous flux. The native leopard frog species have declined dramatically, with the extirpation of many populations and the near extinction of at least one species, the Relict Leopard Frog, Rana onca (Clarkson and Rorabaugh, 1989; Jennings and Hayes, 1994a; Jaeger et al., 2001; Bradford et al., 2005; Savage et al., 2011). Simultaneously, nonnative American Bullfrogs (Rana catesbeiana) and Rio Grande Leopard Frogs (Rana berlandieri) have been introduced and have rapidly expanded their range (Rorabaugh et al., 2002; Stebbins, 2003; Goodward and Wilcox, 2019). Tracking these rapid range shifts has been especially challenging because the native and nonnative leopard frog species (Rana pipiens complex) are so morphologically similar that positive field identifications are often not possible. This presents a serious challenge for survey and monitoring efforts, especially in areas where native species have declined and documenting their persistence could affect conservation management efforts. As a result, researchers have increasingly turned to genetic analyses to identify species and track distributional changes (Platz et al., 1990; Benedict and Quinn, 1999; Rorabaugh et al., 2002).

The Lowland Leopard Frog (Rana yavapaiensis) presents one such challenge. This species is native to the southwestern United States and northern Mexico, ranging from southwestern New Mexico through central and southern Arizona, northern Sonora, southwestern Utah, southern Nevada, and southeastern California. It has declined severely throughout its range, notably in California where it has not been documented since 1965 and may be extirpated (Clarkson and Rorabaugh, 1989; Jennings and Hayes, 1994a, 1994b; Sredl, 2005). As a result of these declines, R. yavapaiensis is designated as a Species of Special Concern in California (Thomson et al., 2016). Because this species may persist at low densities, and because much of its historical range is remote and difficult to survey completely, it remains possible, though extremely unlikely, that the species persists in areas where it has not been documented for long periods. If extant populations of R. yavapaiensis were to be detected in California, this would pose both a dramatic herpetological rediscovery, as well as a critical conservation opportunity. A needed step following any such claim would be a convincing demonstration that the newly found leopard frogs were in fact R. yavapaiensis.

Rana yavapaiensis is difficult to distinguish morphologically from R. berlandieri, which has invaded the historical range of R. yavapaiensis. This poses a serious difficulty in developing a convincing positive species identification if remaining populations of R. yavapaiensis are found. While the two species differ slightly for some morphological characteristics (e.g., R. berlandieri has a slightly larger body size and a different configuration of the dorsolateral folds), their identifying characteristics broadly overlap (Stebbins, 2003). Any populations of R. yavapaiensis that remain in California, or other areas where the species is declining, therefore run the risk of being overlooked due to the close morphological similarities of these two species. Similarly, populations of R. berlandieri, especially if found in locations previously occupied by R. yavapaienis, may be misidentified and incorrectly assumed to be a declining native species of conservation concern. Positive identification of an individual as R. yavapaiensis is therefore dependent on genetic data.

Three populations of leopard frogs that appear to be morphologically intermediate between R. yavapaiensis and R. berlandieri were recently reported in California. In November 2014, B. J. Stacey and J. Keller photographed three leopard frogs southeast of the Salton Sea near Calipatria, Imperial County, California and uploaded these observations to the Reptiles and Amphibians of Southern California (RASCals) project on the iNaturalist citizen science platform (iNaturalist 1066114 and 1066936). These observations were made in a highly-modified agricultural area but within the historical range of R. yavapaiensis and in close proximity (ca. 11.5 km) to the collection locality of one of the last specimens of R. yavapaiensis collected in the state (LACM 91310, collected March 1956). One photograph showed a frog's entire left dorsolateral fold (DLF). Dorsolateral folds are lines of raised glandular skin, and in many ranid frogs, the length and configuration of these folds is a useful character for distinguishing among species (e.g., McAlister, 1962; Pace, 1974). Rana berlandieri and R. yavapaiensis have folds that are discontinuous, with the posterior section inset medially. This short, inset DLF is usually a single fold (i.e., a single, long segment) in R. berlandieri but a series of small dots of raised tissue in R. yavapaiensis (Fig. 1; Pace, 1974; Miera and Sredl, 2000; Stebbins, 2003). Surprisingly, in iNaturalist 1066936, the inset DLF is a series of small dots, consistent with this frog being R. yavapaiensis, which, if confirmed, would make it the first Lowland Leopard Frog documented in the state in over 50 years. Later, in February 2016, one of us (SS) observed leopard frogs in the San Felipe Creek watershed, southwest of the Salton Sea, at some of the last localities where R. yavapaiensis were observed in California (KU 194236–194241 from 1956 and 1959; see also Ruibal, 1959). Because San Felipe Wash is the only historical locality of R. yavapaiensis in California that has not been heavily altered by agricultural development, reservoir construction, and/or irrigation activities, and because this site is isolated and receives only infrequent surveys, there was optimism that these morphologically ambiguous frogs might be R. yavapaiensis. Later, in April 2016, leopard frogs of questionable identity were also observed in Indio, Riverside County, California, north of the Salton Sea, and more than 70 km from the two aforementioned localities. The discovery of leopard frogs that could not be confidently identified to species by multiple biologists highlights the need for a better understanding of the morphological variation in these frogs and a genetically confirmed species identification.

Fig. 1.

Types of dorsolateral folds observed in Rana berlandieri and Rana yavapaiensis. These representations depict the right dorsolateral folds. The color of the medially inset, posterior folds (states 1–6) indicate similarity of the tissue to the primary dorsolateral folds (see text for more information). Black depicts posterior folds that are similar in appearance to the primary fold, and lighter shades of gray indicate skin features that are less pronounced than the raised, glandular tissue of the primary folds. State 1 is typical of R. berlandieri across its native range, and state 5 is typical of R. yavapaiensis. Modeled after Pace (1974: fig. 23).

Fig. 1.

Types of dorsolateral folds observed in Rana berlandieri and Rana yavapaiensis. These representations depict the right dorsolateral folds. The color of the medially inset, posterior folds (states 1–6) indicate similarity of the tissue to the primary dorsolateral folds (see text for more information). Black depicts posterior folds that are similar in appearance to the primary fold, and lighter shades of gray indicate skin features that are less pronounced than the raised, glandular tissue of the primary folds. State 1 is typical of R. berlandieri across its native range, and state 5 is typical of R. yavapaiensis. Modeled after Pace (1974: fig. 23).

Given the extreme conservation importance of any remaining populations of R. yavapaiensis in California, we conducted a morphological and genetic study to both positively identify the morphologically intermediate frogs and to more fully characterize the extent of morphological and genetic variation present within these two species. We positively identified these morphologically intermediate populations as nonnative R. berlandieri. Our morphological and genetic datasets provide a baseline to which future leopard frog populations in California can be compared, thus easing positive identification for any future leopard frog discoveries.

MATERIALS AND METHODS

Field surveys.—

In May 2015, nocturnal searches were conducted along agricultural ditches in the vicinity (<3.5 km) of the iNaturalist observations. We also conducted diurnal surveys of a lower reach of San Felipe Creek near Highway 78. San Felipe Creek was targeted because it was home to one of the last documented populations of R. yavapaiensis in California, although the surveys were conducted downstream of the historical breeding localities (Ruibal, 1959). Frogs encountered during these surveys showed a great deal of variation in the inset DLFs, including several individuals with folds more consistent with R. yavapaiensis than R. berlandieri, as was also seen in iNaturalist 1066936. However, the presence of prominent external vocal sacs in adult males and the occurrence of individuals with inset DLFs typical of R. berlandieri convinced us these were not R. yavapaiensis, and two to nine individuals were collected opportunistically from each of four localities (lower San Felipe Creek and three localities near the iNaturalist observations). Additional surveys were conducted in March 2016 more broadly around Calipatria and the iNaturalist observations. For the 2016 surveys, collections were made at five localities in an area approximately 14 km N-S and 13 km E-W around Calipatria. Surveys in 2016 also included the upper reaches of the San Felipe Creek watershed including localities with historical records of R. yavapaiensis along San Felipe Creek, Fish Creek, and Carrizo Wash. Using the same criteria as in the 2015 surveys, one to 14 individuals were collected opportunistically from each locality. Captured frogs were euthanized via anesthetic overdose (5% benzocaine), and then liver samples were taken prior to preserving the specimen in 10% formalin. Tissue samples and the specimens were deposited at the Natural History Museum of Los Angeles County (LACM). Fieldwork was conducted under a scientific collecting permit issued to GBP by the California Department of Fish and Wildlife (SCP-4307).

Morphological analyses.—

We examined seven external morphological characters that have previously been suggested to be useful in differentiating R. berlandieri and R. yavapaiensis. These included the presence of prominent external vocal sacs in males, five morphometric characters, and condition of the posterior dorsolateral folds. The five morphometric characters were (i) snout–urostyle length (SUL), distance from the tip of the snout to the base of the urostyle; (ii) tibiofibula length; (iii) head width at the posterior margin of the tympanum, measured below the labial fold in cases where a pronounced fold was present (this follows Platz et al., 1990; but see Platz, 1976 for a similar measurement taken at the jaw articulation); (iv) head length, distance from the jaw articulation to the tip of the snout; and (v) tympanum diameter. We selected these characters because they were examined in previous studies of R. yavapaiensis (Platz, 1976; Platz and Frost, 1984; Platz et al., 1990), with characters iii–v being the most informative characters for differentiating R. yavapaiensis and R. berlandieri in a previous statistical analysis (Platz et al., 1990: fig. 5, table 2). Characters ii, iv, and v were measured from the right side of the specimen except in cases of injury or damage to the specimen, in which case the character was measured from the left side. All measurements were taken to the nearest 0.01 mm with a digital caliper, and all measurements were taken by one person (GBP).

We also examined the dorsolateral folds. Rana berlandieri is typically described as having a discontinuous dorsolateral fold in which the posterior section is interrupted and inset medially near the groin (Fig. 1; Pace, 1974; Powell et al., 1998; Miera and Sredl, 2000; Stebbins, 2003). Rana yavapaiensis has a similar configuration except that the small inset fold near the groin tends to be broken into a series of short segments or small dots of raised glandular tissue (Fig. 1; Pace, 1974; Miera and Sredl, 2000; Stebbins, 2003). We scored the inset DLFs as one of the following character states, all of which are depicted in Figure 1:

  • 0.

    DLF is continuous.

  • 1.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF is a single fold (i.e., a single, long segment) similar in appearance to the primary DLF. This is the typical condition in R. berlandieri (see also Pace, 1974: fig. 23d).

  • 2.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF has 2 longer segments, each segment consisting of raised tissue similar in appearance to the primary fold.

  • 3.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF has 1 or 2 segments and ≥1 raised dots of tissue similar in appearance to the primary fold.

  • 4.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF has a single short segment or a series of short segments and dots with the tissue showing some similarity to the raised, glandular tissue in the primary fold.

  • 5.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF has a series of dots with minimal or no similarity to the raised, glandular tissue in the primary fold. This is the typical condition in R. yavapaiensis (see also Pace, 1974: fig. 23e).

  • 6.

    DLF is discontinuous, with the posterior section inset medially. Inset DLF has some dots but not in a distinct row.

  • 7.

    DLF is absent in the groin region (i.e., no posterior DLF).

The left and right inset DLFs were found to often vary in the same individual, so both sides were scored. Although this could introduce pseudoreplication into statistical analyses of the DLFs, our goal was to examine the extent to which characters frequently used for identification were diagnostic; thus, understanding the overall levels of variability was a more central concern. The character states were treated as unordered, though states 1 through 5 can be considered points along a continuum of variation from the typical condition in R. berlandieri through the typical condition in R. yavapaiensis. However, how these states are related developmentally and how they relate to states 0, 6, and 7 are unknown.

Sex was determined by examination of the size and shape of the first digit of the hand especially presence/absence of nuptial pads, examination of the gonads in those specimens previously cut open or for which permission to examine the gonads was secured, and presence/absence of a prominent external vocal sac. This third character is only reliable for differentiating sexes in R. berlandieri because R. yavapaiensis does not have a prominent external vocal sac, although a less-conspicuous vocal sac can be observed in some adult males during the breeding season (e.g., LACM 13846).

The morphological analyses were largely restricted to adult individuals. Platz (1988) lists minimum sizes at sexual maturity for R. yavapaiensis as 46 mm for males and 53 mm for females, but these cut-offs are for SVL, not SUL, which was measured here. Depending on measurement technique, SVL is either the same or slightly larger than SUL. Similar cut-off data do not exist for R. berlandieri so specimens between 45 and 55 mm SUL were examined for development of secondary sexual characteristics in males (vocal sac and nuptial pads) and ovaries in females. To add to the existing knowledge of minimum size at sexual maturity for R. yavapaiensis, we similarly examined males and females between 46 and 53 mm SUL.

Sex and all morphological characters were examined for 84 museum specimens of R. berlandieri from their native range in Texas, New Mexico, and Mexico; 41 museum specimens of R. yavapaiensis; and 38 recently collected leopard frogs of unknown identity from southeastern California. One R. yavapaiensis was not scored for inset DLFs because of inadequate preservation. After finding high levels of variation in the inset DLFs of some South and West Texas populations of R. berlandieri, we scored 6–14 additional individuals from each of four localities (total of 50 additional individuals). We also scored inset DLFs for an additional four R. berlandieri, nine unknown individuals from southeastern California, and four R. yavapaiensis. In total, we scored inset DLFs from 229 specimens. Locality information, sex, morphological measurements, and DLF condition are provided in the supplemental materials (see Data Accessibility).

To assess whether morphometric trait values differed among R. berlandieri, R. yavapaiensis, and the frogs of uncertain identity, we performed analyses of variance (ANOVA). To do so, we log transformed each trait and then divided the transformed tibiofibula length, head width, head length, and tympanum diameter by transformed SUL. We treated the two species and the uncertain individuals as separate groups and performed a one-way ANOVA on SUL and each of the four morphometric ratios to ask whether trait means differed among the three groups. When significant differences were found, we identified which comparisons were different using a Tukey multiple pairwise comparison. We then tested for normality of residuals for each ANOVA using a Shapiro-Wilk test and by visualizing the residuals directly using q-q plots. These analyses were carried out in the R programming environment (R Core Team, 2017) and made use of the dplyr v0.7.4 (Wickham et al., 2018) and car v2.1-5 packages (Fox and Weisberg, 2019).

Molecular identification.—

We used molecular tools to determine the identity of the unknown individuals. To do so, we assembled a comparative panel from individuals that spanned much of the native range of R. berlandieri in the US (n = 5), two localities in the nonnative range of R. berlandieri in Arizona (n = 10), and two localities of R. yavapaiensis in Arizona (n = 4). We included only individuals for which morphological identifications were unequivocal. We also included samples of the unidentified California frogs, including two from Indio, four from the San Felipe Creek watershed (three from the upper watershed at Fish Creek, and one from the lower watershed), and 11 collected across seven sampling sites in the vicinity of Calipatria (Fig. 2). Many of the individuals from San Felipe Creek and the vicinity of Calipatria were specifically selected because they had inset DLFs more representative of R. yavapaiensis than R. berlandieri (i.e., inset DLF states 4, 5, and 7). Locality information and GenBank numbers are provided in the supplemental materials (see Data Accessibility).

Fig. 2.

Map of southeastern California, USA depicting historical localities of Rana yavapaiensis (black dots), contemporary localities of Rana berlandieri (white dots), and the uncertain specimens that we examine here (red dots).

Fig. 2.

Map of southeastern California, USA depicting historical localities of Rana yavapaiensis (black dots), contemporary localities of Rana berlandieri (white dots), and the uncertain specimens that we examine here (red dots).

We extracted genomic DNA from toe or liver tissue using a salt extraction protocol (Sambrook and Russell, 2001), then amplified and sequenced two mitochondrial markers (NADH subunit two and 12S) and four nuclear markers (NTF3, RAG1, Rhodopsin, and Tyrosinase) for all individuals (Table 1). We amplified each marker using PCR in either a 11 μl MangoTaq (Bioline) or 25 μl GoTaqGreen (Promega) mediated reaction. Each reaction was carried out with an initial denaturation of 1 min at 95°C, 38 subsequent cycles with denaturation for 30 s at 94°C, locus-specific annealing temperatures (see Table 1) for 45 s, and extension at 72°C for 1 min, followed by a final extension of 10 min at 72°C. Products were then purified with ExoSAP-IT (Affymetrix) and sequenced in both directions on an Applied Biosystems 3730XL DNA Analyzer at the University of Hawai'i.

Table 1.

Primers used for amplification and sequencing and the associated model of sequence evolution selected via jModelTest2 and used in the Bayesian analyses.

Primers used for amplification and sequencing and the associated model of sequence evolution selected via jModelTest2 and used in the Bayesian analyses.
Primers used for amplification and sequencing and the associated model of sequence evolution selected via jModelTest2 and used in the Bayesian analyses.

We edited, assembled, and aligned sequences using MUSCLE in Geneious v7.1.5 (Kearse et al., 2012). We checked each nuclear locus for heterozygous positions and coded these using IUPAC ambiguity codes. For each alignment, we selected a model of molecular evolution using Akaike Information Criterion (AICc) implemented in jModelTest2 (Darriba et al., 2012). We then estimated phylogeny using MrBayes v3.2.5 (Ronquist et al., 2012) for each individual marker and the concatenated dataset as a whole. These analyses used four independent runs that each used one cold and three incrementally heated Markov chains (temperature parameter = 0.05). We ran each analysis for 10 million generations, recording the state of the cold chain every 10,000th generation. We examined the Markov chain Monte Carlo (MCMC) output in Tracer v1.6 to verify that all chains mixed well, appeared to converge to the same target distribution, and had obtained an adequate number of effectively independent samples from the posterior (ESS > 1000 for all parameters). Finally, we discarded the first 25% of samples as ‘burnin,' combined the remaining samples from the four runs, and summarized the results.

RESULTS

Morphological analyses.—

We found a large amount of variation in most morphological traits that we measured. No single measurement was adequate to identify the two species or the unknown individuals with certainty. We did find differences in means for several traits, although their distributions often overlapped widely (Table 2, Fig. 3).

Table 2.

Mean (and standard deviation) for morphometric traits measured in this study. SUL: snout–urostyle length; TFL: tibiofibula length; HL: head length; HW: head width; TYMP: diameter of tympanum.

Mean (and standard deviation) for morphometric traits measured in this study. SUL: snout–urostyle length; TFL: tibiofibula length; HL: head length; HW: head width; TYMP: diameter of tympanum.
Mean (and standard deviation) for morphometric traits measured in this study. SUL: snout–urostyle length; TFL: tibiofibula length; HL: head length; HW: head width; TYMP: diameter of tympanum.
Fig. 3.

Distribution of morphometric characters for Rana berlandieri (Rb), uncertain individuals, and Rana yavapaiensis (Ry). Trait values are log transformed. Different letters indicate significant differences between groups based on Tukey multiple comparison tests (ANOVA, P < 0.05). Abbreviations as in Table 2.

Fig. 3.

Distribution of morphometric characters for Rana berlandieri (Rb), uncertain individuals, and Rana yavapaiensis (Ry). Trait values are log transformed. Different letters indicate significant differences between groups based on Tukey multiple comparison tests (ANOVA, P < 0.05). Abbreviations as in Table 2.

The most conspicuous morphological trait that differentiates R. berlandieri and R. yavapaiensis is the relatively larger size (measured by SUL) of adult R. berlandieri (Table 2). In addition, adult male R. berlandieri have well developed vocal sacs, and these are absent or only weakly developed in R. yavapaiensis (Platz and Frost, 1984). Adult males from the Southern California populations had well developed vocal sacs, suggesting that they are R. berlandieri. For large and male frogs, these two traits may suffice to reliably differentiate the two species and identify uncertain individuals.

The remaining morphometric traits were far less useful for differentiating the two species, or for identification of the uncertain specimens, owing to a large amount of overlap in trait values (Table 2, Fig. 3). One-way analysis of variance confirms that R. berlandieri has a larger mean SUL and a relatively smaller head width and length than R. yavapaiensis. We found no significant difference between the two species in relative length of the tibiofibula or size of the tympanum. Finally, we observed a large amount of variation in the condition of the dorsolateral folds within both R. berlandieri and R. yavapaiensis, as well as between the left and right dorsolateral fold of individual frogs (Fig. 4). The most common DLF character state differed between the two species, although both showed a large and overlapping amount of variation, such that it is not possible to identify an individual frog with this trait alone. We also saw an effect of locality in the DLF condition for R. berlandieri (Table 3). Twenty of the 35 R. berlandieri with dorsolateral folds in states 4 through 7 came from a single locality (Miller Ranch, Presidio County, TX). Thus, across most of the native US range of this species, states 4 through 7 are uncommon. Yet, these same character states make up 27% of the observations for the specimens with uncertain identification in California.

Fig. 4.

Character states for the left and right dorsolateral folds for R. berlandieri, uncertain individuals, and R. yavapaiensis. Points along the dark line indicate individuals with the same character state for the left and right dorsolateral folds. See text for description of character states. We have applied a small amount of random noise to the points to enhance visibility.

Fig. 4.

Character states for the left and right dorsolateral folds for R. berlandieri, uncertain individuals, and R. yavapaiensis. Points along the dark line indicate individuals with the same character state for the left and right dorsolateral folds. See text for description of character states. We have applied a small amount of random noise to the points to enhance visibility.

Table 3.

Variation in dorsolateral folds. Numbers given are counts of dorsolateral fold character states for the given species or locality. We scored left and right dorsolateral folds for each individual. Thus, the sample size (n) corresponds to twice the number of specimens examined.

Variation in dorsolateral folds. Numbers given are counts of dorsolateral fold character states for the given species or locality. We scored left and right dorsolateral folds for each individual. Thus, the sample size (n) corresponds to twice the number of specimens examined.
Variation in dorsolateral folds. Numbers given are counts of dorsolateral fold character states for the given species or locality. We scored left and right dorsolateral folds for each individual. Thus, the sample size (n) corresponds to twice the number of specimens examined.

The specimens with uncertain identification did not show clear affinities with either R. berlandieri or R. yavapaiensis. They were intermediate between the two in size and head width, and significantly more similar to R. yavapaiensis for head length (Fig. 3). The dorsolateral folds were more similar to those of R. berlandieri, although again a large amount of variation within and across individuals was present, and some specimens had character states identical to R. yavapaiensis. Taken together, their larger size, head width, and presence of well-developed vocal sacs suggests that they are R. berlandieri, although the widely overlapping morphometric data and dorsolateral folds suggest that caution is warranted when making identifications based on morphology, especially when vocal sac condition cannot be assessed.

Molecular identification.—

We obtained a total of 4066 bp (1706 bp mitochondrial, 2360 bp nuclear) of sequence data for the 36 individuals. All sequences have been deposited in GenBank (MT114437–MT114472 for 12S; MT112715–MT112858 for NTF3, Rhodopsin, Tyrosinase, and ND2; and MT124627–MT124661 for RAG1; see Data Accessibility). The aligned matrix was nearly complete and contained 2.2% missing data that arose from a small amount of low quality sequence data that we trimmed, and a small number of insertion–deletion events. The phylogenetic analyses mixed well and appeared to converge within the first 25% of MCMC samples.

The phylogenetic analysis of the concatenated alignment recovers a well-supported (Posterior Probability ∼1.0) bipartition separating R. berlandieri from R. yavapaiensis (Fig. 5). All frogs of questionable identity are included within the R. berlandieri clade with high posterior probability, indicating that they belong to this species. The independent gene tree estimates were less well resolved but were compatible with the concatenated tree. None of the gene trees gave any indication that the uncertain specimens could potentially be R. yavapaiensis, or be hybrids of R. yavapaiensis and R. berlandieri. The tree contained no strong geographic signal within either the R. berlandieri or R. yavapaiensis clades.

Fig. 5.

Majority rule consensus tree of the posterior distribution of trees from the Bayesian analysis of the concatenated alignment (two mitochrondrial and four nuclear markers, with the model of evolution for each marker specified in Table 1). Individuals with uncertain species identification are shown in bold. Numbers below nodes are estimated posterior probabilities. See Data Accessibility for tree file.

Fig. 5.

Majority rule consensus tree of the posterior distribution of trees from the Bayesian analysis of the concatenated alignment (two mitochrondrial and four nuclear markers, with the model of evolution for each marker specified in Table 1). Individuals with uncertain species identification are shown in bold. Numbers below nodes are estimated posterior probabilities. See Data Accessibility for tree file.

DISCUSSION

Morphological variation.—

Our results have several implications for ongoing efforts to document, conserve, and manage declining leopard frog populations in the southwestern United States. The dataset reported here provides a quantitative assessment of the extent of morphological variation within wild populations of an invasive species with an expanding range. Morphological traits that have previously been suggested as being useful for differentiating the two species do vary between R. yavapaiensis and R. berlandieri, but there is a great deal of overlap. Many individual frogs cannot be positively identified on the basis of morphology alone, using any known morphological character. Thus, these traits offer little utility for field identifications.

Most importantly, our morphological analyses confirm that the invasive populations of R. berlandieri are morphologically intermediate, and for some characters, are actually more similar to R. yavapaiensis than they are to native-range R. berlandieri. The condition of the DLFs is an especially useful character for rapid identification of many leopard frog species (McAlister, 1962; Pace, 1974; Stebbins, 2003), but in southeastern California, some nonnative R. berlandieri have DLFs similar to those observed in R. yavapaiensis and rarely observed in native-range populations of R. berlandieri (Table 3; Fig. 4; supplemental material; see Data Accessibility). Further, two of five morphometric traits are statistically distinguishable between introduced and native populations of R. berlandieri (HL and HW; Fig. 3), with nonnative R. berlandieri being more similar to R. yavapaiensis than would be expected based on analyses of native-range R. berlandieri (see also SUL; Fig. 3A).

These results highlight the potential for morphometric trait means to vary in different parts of a species' range, particularly when parts of the range are the result of recent expansions. Thus, we advise caution when using morphological traits to identify leopard frogs, particularly when these traits are being applied to populations outside of those in which the traits were originally assessed.

The unexpected morphological variation observed in nonnative R. berlandieri could, in principle, simply result from differences in life stage or growth rates in different parts of the species range, plastic responses to varying environmental cues across the range, or due to evolution of the trait itself. Rapid evolution of morphometric traits has been documented in other invasive populations. This is especially common on invasion fronts, where dispersal-adapted morphological traits are favored, potentially resulting in range-edge phenotypes that are not observed in the range-core area of introduced populations nor in the native range (Chuang and Peterson, 2016). Such patterns have been especially well studied in invasive anuran populations (Phillips et al., 2006; Shine et al., 2011; Perkins et al., 2013) and are also found in other taxa (Laparie et al., 2013). Adaptation to novel food sources, habitats, and climates can also result in rapid morphological shifts in introduced populations (Losos et al., 1997; Herrel et al., 2008). While speculative, it is possible that such a process could be occurring here as well. The available habitat for leopard frogs in Southern California has changed dramatically over the last half century, and it is possible that this change in habitat has led to plastic and/or evolutionary changes in the morphometric characters that we use to recognize species. The fact that these traits are useful for species identification (at least for other leopard frog species pairs) means that they are variable and therefore relatively labile. A better understanding of this variation would be most easily generated by assembling larger morphometric datasets across species ranges and directly assessing variation (and co-variation) that is present. Such data are conceptually simple to collect and would greatly help our understanding of the morphological diversity present in this threatened group of frogs.

Molecular identification and the need for more surveys.—

Our phylogenetic analysis establishes the identity of morphologically intermediate leopard frog populations in Southern California as Rana berlandieri. We find no evidence for hybrid ancestry or other potential explanations for the morphological intermediacy of these populations, although it is possible that this could exist and go undetected in an analysis of a small number of loci, as we have done here. These molecular tools provide a simple means of species identification for future surveys if additional morphologically intermediate populations are discovered.

It deserves emphasis that these results do not change our understanding of whether R. yavapaiensis is extirpated in California. Rather, they should encourage additional survey efforts and the collection of genetic samples. Leopard frogs are known to persist at low densities and to undergo rather marked fluctuations in population size in other areas of the desert southwest (Rorabaugh, 2005). If the species persists in California, it probably does so at low densities. In addition, surveys of tadpoles in suitable waterways could be undertaken to increase detectability. Non-lethal tissue samples are simple to collect for leopard frog tadpoles, and this may serve as our best chance to rapidly inventory and identify leopard frog populations within the state.

Native populations of Rana berlandieri encompass an extensive geographic range, occurring from Texas and New Mexico south through Mexico and into Central America. More recent introductions of this species into Arizona, California, and Baja California have led to the establishment of populations that appear to be expanding into additional suitable habitat in these regions (Rorabaugh et al., 2002; Goodward and Wilcox, 2019). This may partially be facilitated by the fact that R. berlandieri is successful in a wide range of habitats, and in particular, human modified habitats such as agricultural landscapes. By contrast, R. yavapaiensis has experienced extensive declines throughout its range over the past century, and appears to be less tolerant of anthropogenic disturbances (Clarkson and Rorabaugh, 1989). Thus, there is general concern about the potential displacement of any existing native leopard frog populations by introduced R. berlandieri, especially in the highly modified remaining habitats in California. This should also motivate surveys of additional localities to more definitively determine if any remaining populations of R. yavapaiensis might persist in California as well as to monitor changing distributions of these two species and potential threats posed by R. berlandieri in Arizona.

DATA ACCESSIBILITY

Supplemental material is available at https://www.copeiajournal.org/ch-19-222.

ACKNOWLEDGMENTS

We thank J. Keller and B. J. Stacey for their incredible efforts to document biodiversity; together they have contributed over 115,000 observations to iNaturalist, including the observations that first motivated this research effort. M. Sredl kindly provided tissues of both species from Arizona and was a tremendous source of natural history information. We also thank B. Zitt for contributing information and tissue samples from the Indio population. For assistance in the field, we thank B. Hardy, A. Backlin, and C. Demetropoulos. For permission to examine specimens, and assistance during museum visits, we thank T. LaDuc and D. Cannatella at the Texas Natural History Collections and B. Hollingsworth at the San Diego Natural History Museum. T. LaDuc and D. Cannatella of the TNHC also kindly loaned tissue samples of R. berlandieri. M. Sredl, J. Rorabaugh, and R. Jennings provided helpful advice. Funding for fieldwork and museum visits was provided by the Urban Nature Research Center (UNRC) at the Natural History Museum of Los Angeles County, and UNRC staff kindly provided feedback on earlier versions of the manuscript. Grants from the National Science Foundation (DEB-1354506 and DBI-1356796) provided funding for molecular and computational resources. M. Shaulsky was supported by a National Science Foundation Research Experience for Undergraduates site award (DBI-1560491). A. Barley was supported by a postdoctoral fellowship award from the Arnold and Mabel Beckman Foundation. This work was approved by the Institutional Animal Care and Use Committee of the University of Southern California.

LITERATURE CITED

LITERATURE CITED
Benedict,
N.,
and
Quinn.
T. W.
1999
.
Identification of Rio Grande Leopard Frogs by mitochondrial DNA analyses: a tool for monitoring the spread of a non-native species
.
Report submitted September 30, 1999 to the Arizona Game and Fish Department.
Bossuyt,
F.,
and
Milinkovitch.
M. C.
2000
.
Convergent adaptive radiations in Madagascan and Asian ranid frogs reveal covariation between larval and adult traits
.
Proceedings of the National Academy of Sciences of the United States of America
97
:
6585
6590
.
Bradford,
D. F.,
Jennings,
R. D.
and
Jaeger.
J. R.
2005
.
Rana onca Cope, 1875(b) Relict Leopard Frog
,
p.
567
568
.
In:
Amphibian Declines: The Conservation Status of United States Species
.
Lannoo
M. J.
(ed.)
University of California Press
,
Berkeley, California
.
Chuang,
A.,
and
Peterson.
C. R.
2016
.
Expanding population edges: theories, traits, and trade-offs
.
Global Change Biology
22
:
494
512
.
Clarkson,
R. W.,
and
Rorabaugh.
J. C.
1989
.
Status of leopard frogs (Rana pipiens complex: Ranidae) in Arizona and southeastern California
.
Southwestern Naturalist
34
:
531
538
.
Darriba,
D.,
Taboada,
G. L.
Doallo,
R.
and
Posada.
D.
2012
.
jModelTest 2: more models, new heuristics and parallel computing
.
Nature Methods
9
:
772
.
Fox,
J.,
and
Weisberg.
S.
2019
.
An R Companion to Applied Regression. Third edition
.
Sage
,
Thousand Oaks, California
.
Goodward,
D. M.,
and
Wilcox.
M. D.
2019
.
The Rio Grande leopard frog (Lithobates berlandieri) and other introduced and native riparian herpetofauna of the Coachella Valley, Riverside County, California
.
California Fish and Game
105
:
48
71
.
Herrel,
A.,
Huyghe,
K.
Vanhooydonck,
B.
Backeljau,
T.
Breugelmans,
K.
Grbac,
I.
Van Damme,
R.
and
Irschick.
D. J.
2008
.
Rapid large-scale evolutionary divergence in morphology and performance associated with exploitation of a different dietary resource
.
Proceedings of the National Academy of Sciences of the United States of America
105
:
4792
4795
.
Jaeger,
J. R.,
Riddle,
B. R.
Jennings,
R. D.
and
Bradford.
D. F.
2001
.
Rediscovering Rana onca: evidence for phylogenetically distinct leopard frogs from the border region of Nevada, Utah, and Arizona
.
Copeia
2001
:
339
354
.
Jennings,
M. R.,
and
Hayes.
M. P.
1994
a.
Decline of native ranid frogs in the desert southwest
,
p.
185
213
.
In:
Herpetology of the North American Deserts. Proceedings of a Symposium
.
Brown
P. R.
and
Wright
J. W.
(eds.)
Southwestern Herpetologists Society
,
Special Publications No. 5
.
Jennings,
M. R.,
and
Hayes.
M. P.
1994
b.
Amphibian and Reptile Species of Special Concern in California
.
California Department of Fish and Game, Inland Fisheries Division
,
Rancho Cordova, California
.
Kearse,
M.,
Moir,
R.
Wilson,
A.
Stones-Havas,
S.
Cheung,
M.
Sturrock,
S.
Buxton,
S.
Cooper,
A.
Markowitz,
S.
Duran,
C.
Thierer,
T.
Ashton,
B.
Meintjes,
P.
and
Drummond.
A.
2012
.
Geneious Basic: an integrated and extendable desktop software platform for the organization and analysis of sequence data
.
Bioinformatics
28
:
1647
1649
.
Laparie,
M.,
Renault,
D.
Lebouvier,
M.
and
Delattre.
T.
2013
.
Is dispersal promoted at the invasion front? Morphological analysis of a ground beetle invading the Kerguelen Islands, Merizodus soledadinus (Coleoptera, Carabidae)
.
Biological Invasions
15
:
1641
1648
.
Losos,
J. B.,
Warheit,
K. I.
and
Schoener.
T. W.
1997
.
Adaptive differentiation following experimental island colonization in Anolis lizards
.
Nature
387
:
70
73
.
McAlister,
W. H.
1962
.
Variation in Rana pipiens Schreber in Texas
.
American Midland Naturalist
67
:
334
363
.
Miera,
V.,
and
Sredl.
M. J.
2000
.
Museum photo voucher series and photo identification cards: status of the Rio Grande Leopard Frog in Arizona
.
Unpublished report submitted to the Nongame and Endangered Wildlife Program
,
Arizona Game and Fish Department
.
Newman,
C. E.,
Feinberg,
J. A.
Rissler,
L. J.
Burger,
J.
and
Shaffer.
H. B.
2012
.
A new species of leopard frog (Anura: Ranidae) from the urban northeastern US
.
Molecular Phylogenetics and Evolution
63
:
445
455
.
Pace,
A. E.
1974
.
Systematic and biological studies of the leopard frogs (Rana pipiens complex) of the United States
.
Miscellaneous Publications Museum of Zoology, University of Michigan
,
No. 148
.
Perkins,
A. T.,
Phillips,
B. L.
Baskett,
M. L.
and
Hastings.
A.
2013
.
Evolution of dispersal and life history interact to drive accelerating spread of an invasive species
.
Ecology Letters
16
:
1079
1087
.
Phillips,
B. L.,
Brown,
G. P.
Webb,
J. K.
and
Shine.
R.
2006
.
Invasion and the evolution of speed in toads
.
Nature
439
:
803
.
Platz,
J. E.
1976
.
Biochemical and morphological variation of leopard frogs in Arizona
.
Copeia
1976
:
660
672
.
Platz,
J. E.
1988
.
Rana yavapaiensis
.
Catalogue of American Amphibians and Reptiles
418
:
1
2
.
Platz,
J. E.,
Clarkson,
R. W.
Rorabaugh,
J. C.
and
Hillis.
D. M.
1990
.
Rana berlandieri: recently introduced populations in Arizona and southeastern California
.
Copeia
1990
:
324
333
.
Platz,
J. E.,
and
Frost.
J. S.
1984
.
Rana yavapaiensis, a new species of leopard frog (Rana pipiens complex)
.
Copeia
1984
:
940
948
.
Powell,
R.,
Collins,
J. T.
and
Hooper,
E. D.
Jr.
1998
.
A Key to the Amphibians and Reptiles of the Continental United States and Canada
.
University Press of Kansas
,
Lawrence, Kansas
.
R Core Team
.
2017
.
R: a language and environment for statistical computing
.
R Foundation for Statistical Computing
,
Vienna, Austria
.
Ronquist,
F.,
Teslenko,
M.
van der Mark,
P.
Ayres,
D. L.
Darling,
A.
Hohna,
S.
Larget,
B.
Liu,
L.
Suchard,
M. A.
and
Huelsenbeck.
J. P.
2012
.
MrBayes 3.2: efficient Bayesian phylogenetic inference and model choice across a large model space
.
Systematic Biology
61
:
539
542
.
Rorabaugh,
J. C.
2005
.
Rana berlandieri Baird, 1854(a) Rio Grande Leopard Frog
,
p.
530
532
.
In:
Amphibian Declines: The Conservation Status of United States Species
.
Lannoo
M. J.
(ed.)
University of California Press
,
Berkeley, California
.
Rorabaugh,
J. C.,
Sredl,
M. J.
Miera,
V.
and
Drost.
C. A.
2002
.
Continued invasion by an introduced frog (Rana berlandieri): Southwestern Arizona, Southeastern California, and Río Colorado, México
.
Southwestern Naturalist
47
:
12
20
.
Ruibal,
R.
1959
.
The ecology of a brackish water population of Rana pipiens
.
Copeia
1959
:
315
322
.
Sambrook,
J.,
and
Russell.
D. W.
2001
.
Molecular Cloning: A Laboratory Manual. Third edition
.
Cold Spring Harbor Laboratory Press
,
New York
.
Savage,
A. E.,
Sredl,
M. J.
and
Zamudio.
K. R.
2011
.
Disease dynamics vary spatially and temporally in a North American amphibian
.
Biological Conservation
144
:
1910
1915
.
Shine,
R.,
Brown,
G. P.
and
Phillips.
B. L.
2011
.
An evolutionary process that assembles phenotypes through space rather than through time
.
Proceedings of the National Academy of Sciences of the United States of America
108
:
5708
5711
.
Sredl,
M. J.
2005
.
Rana yavapaiensis Platz and Frost, 1984 Lowland Leopard Frog
,
p.
596
599
.
In:
Amphibian Declines: The Conservation Status of United States Species
.
Lannoo
M. J.
(ed.)
University of California Press
,
Berkeley, California
.
Stebbins,
R. C.
2003
.
A Field Guide to Western Reptiles and Amphibians. Third edition
.
Houghton Mifflin Company
,
New York
.
Thomson,
R. C.,
Wright,
A. N.
and
Shaffer.
H. B.
2016
.
California Amphibian and Reptile Species of Special Concern
.
University of California Press
,
Oakland, California
.
Wickham,
H.,
François,
R.
Henry,
L.
and
Müller.
K.
2018
.
dplyr: a grammar of data manipulation. R package version 0.7.6
.
Zaldívar-Riverón,
A.,
León-Regagnon,
V.
and
Nieto-Montes de Oca.
A.
2004
.
Phylogeny of the Mexican coastal leopard frogs of the Rana berlandieri group based on mtDNA sequences
.
Molecular Phylogenetics and Evolution
30
:
38
49
.

Author notes

1

Department of Herpetology, Natural History Museum of Los Angeles County, 900 Exposition Boulevard, Los Angeles, California 90007; Email: (GBP) gpauly@nhm.org. Send reprint requests to this address.

2

Urban Nature Research Center, Natural History Museum of Los Angeles County, 900 Exposition Boulevard, Los Angeles, California 90007.

3

School of Life Sciences, University of Hawai'i, Honolulu, Hawai'i 96822; Email: (MCS) mayacs@hawaii.edu; (AJB) ajbarley@hawaii.edu; (SKG) skennedygold@gmail.com; and (RCT) thomsonr@hawaii.edu.

4

Southwest Aquatic & Terrestrial Biology, 3225 Maine Avenue, Long Beach, California 90806; Email: bullfrogslayer@gmail.com.

5

California Department of Fish and Wildlife, 78078 Country Club Drive, Suite 109, Bermuda Dunes, California 92203; Email: sharon.keeney@wildlife.ca.gov.

Associate Editor: B. Stuart.