Taxonomy of Horned Lizards, Genus Phrynosoma (Squamata, Phrynosomatidae)

: In this article, I revise the taxonomy of the species and subspecies of the genus Phrynosoma through phylogenetic and species delimitation approaches based on four mtDNA markers (ND1, ND2, ND4, and 12S). The resulting taxonomy recognizes 12 species ( P. asio , P. bracconieri , P. cornutum , P. coronatum , P. douglasii , P. hernandesi , P. mcallii , P. modestum , P. orbiculare , P. platyrhinos , P. solare , and P. taurus ). Several of these species are divided into subspecies as follows: P. coronatum ( P. c. coronatum , P. c. blainvillii , P. c. cerroense , and P. c. frontale ), P. cornutum ( P. c. cornutum and P. c. bufonium ), P. hernandesi ( P. h. hernandesi , P. h. ditmarsi , and P. h. ornatum ), P. orbiculare ( P. o. orbiculare , P. o. bradti , P. o. boucardii , P. o. cortezii , P. o. dugesii , and P. o. durangoensis ), P. platyrhinos ( P. p. platyrhinos and P. p. goodei ), P. taurus ( P. t. taurus and P. t. sherbrookei ). In this coherent and objective approach, those taxa treated here as subspecies have diverged to a much lesser degree than those that are herein recognized as separate species. Typically, those taxa recognized as subspecies are one another’s closest relatives (i.e., they together form a monophyletic group that represents the species) and are distributed allopatrically. In this approach, all separate evolutionarily signiﬁcant units are recognized as named taxa—either species or subspecies—thereby reﬂecting the importance of identifying and naming such units for conservation. I provide a checklist of the recognized species and subspecies of Phrynosoma along with synonymies and distribution maps.


Introduction
The lizards of the genus Phrynosoma form a conspicuous component of the herpetofauna of parts of Canada, the USA, and Mexico. Due to their unique morphology, especially the more or less well developed spines on the head and the stout body shape, they are easy to recognize at a generic level. However, on the species level, the picture emerging from the available literature is somewhat more complicated.
In the course of preparing a new and updated edition of the book Krötenechsen (the German word for "horned lizards"), originally published by Betrand Baur and Richard Montanucci in 1998 [1], I undertook a literature review in order to bring the book up to date with respect to the current taxonomy of this genus. A large body of publications has accrued since the first edition, including several in-depth revisions of certain species complexes (e.g., [2,3]) and studies on the phylogenetic relationships among species based not only upon various mtDNA and nuclear DNA markers (e.g., [4][5][6][7][8]) but also on genomic sequence data (e.g., [9,10]). However, substantial uncertainty exists as to the recognition of various taxa that are apparently closely related. This is particularly true for the taxa in the P. orbiculare, P. douglasii, P. coronatum, and P. platyrhinos species complexes [2][3][4]11,12]. A wealth of morphological and genetic information has been produced by the combined efforts of various herpetologists [13][14][15][16][17][18]. A summary of recent changes in the taxonomy of horned lizards was published by Sherbrooke [19] and a list of the scientific and standard English names of Phrynosoma species occurring north of Mexico can be found in Crother [20]. However, evaluation of the literature revealed that the current taxonomy of the genus Phrynosoma has several weaknesses. The most obvious problem is the lack of a coherent analysis of evolutionary divergence among all taxa in the genus. There are species that are well-differentiated both morphologically and genetically whereas other taxa (e.g., some of those included in the species complexes listed above) seem to have much lower levels of differentiation.
Frost and Hillis [21] (p. 92) proposed that the subspecies category could be used for intraspecific historical sublineages that are "not incontrovertibly removed from the possibility of interaction with other sublineages" if the "historical relationships of the subpopulations" can be recovered by phylogenetic analyses. Implicitly, they acknowledged that such subspecies would not occur sympatrically. They argued against the use of the subspecies category for populations that represent arbitrary slices of a reproductively continuous geographic cline or for those that occur in sympatry. Recently, Hillis [22] has presented some stimulating thoughts regarding the use of the species versus subspecies categories. He explained that the geographic variation found in many wide-ranged species is caused by the presence of historical sublineages that had historically diverged in isolation but, after coming back into contact, had begun to merge again. Thus, in these species, different geographically segregated phenotypes and genetic entities can be identified but, along the contact zones, these sublineages often show evidence of hybridization. He concluded that "subspecies designations are ideal for this application" [22] (p. 54). In the following, I will refer to this as the Frost and Hillis subspecies concept. I propose the following definitions: a species level unit is discovered as a monophyletic group through genetic analyses with a degree of divergence typical for most of the taxa in the genus that are unambiguously recognized as species level units. A subspecies level unit is discovered as a monophyletic group through genetic analyses with a degree of divergence significantly lower than is typical for most of the taxa in the genus that are unambiguously recognized as species level units. It is allopatrically distributed relative to the other subspecies level taxa within the taxon recognized as the species to which they all belong. Together these subspecies level taxa form a monophyletic group which represents the species level unit. In practice, an evaluation of the differentiation pattern that is evident among the studied organisms will yield the needed information for the calibration of the taxonomic system. In cases in which all studied taxa are at same level of genetic differentiation, they should all be given the same taxonomic rank. If, however, there are two or more categories of significantly different genetic differentiation, they should not all be assigned the same taxonomic rank, but rather should be assessed differently from a taxonomic perspective. For the decision as to whether two different categories of genetic differentiation should then be classified as species and subspecies or as separate genera and species, the outgroup comparison can be helpful, i.e., a look at the current taxonomic system of closely related taxa.
I suggest the use of genetic data for the evaluation of differentiation because morphological data are usually much more subjective. This is not to say that morphological data are useless for phylogenetic analyses. On the contrary, I strongly believe that morphological data are essential for the definition and diagnosis of taxa as well as for identifying morphological synapomorphies that can be used in phylogenetic analyses. However, the history of taxonomy-in the genus Phrynosoma and elsewhere-shows that judgments of the usefulness of certain morphological characters for diagnosing certain taxa often vary drastically depending upon the individual researcher. The usefulness of mtDNA data for the detection of phylogenetic lineages within a genus of amphibians or reptiles has been demonstrated for many groups of organisms (e.g., [23][24][25]). Since most of the commonly used mtDNA markers seem to be relatively fast evolving, they are well-suited to the detection of genetic divergences among recently split taxa. We even tend to find more clades in an mtDNA tree than are supported by other evidence lines, such as external morphology, hemipenial morphology in squamate reptiles, or bioacoustics in frogs [26]. On the other hand, if the mtDNA data do not indicate any genetic structuring among the studied populations, then it can safely be concluded that there is no evidence for the recognition of more than one taxon for this set of populations. I know of no example where nuclear DNA data yielded more resolution for recent splits than mtDNA data. In this sense, the analysis of mtDNA data provides important information that can be used to evaluate the taxonomy of a given group of organisms.
In the genus Phrynosoma, we deal with several taxa that show a high degree of geographic variation. As is typical for these cases, a great deal of controversy exists regarding the recognition of taxa and regarding their respective taxonomic level-species versus subspecies. For example, in the P. coronatum complex, the number of recognized taxa vary from one to seven depending on the author (e.g., [3,11,27]). Similar problems are evident for taxa in the P. orbiculare, P. douglasii, and P. platyrhinos species complexes, respectively. There is no modern and coherent analysis and taxonomic assessment for the genus Phrynosoma that can use the Frost and Hillis subspecies concept to categorize the taxa within the genus into species and subspecies. The published taxonomic revisions are restricted to certain species groups, making it impossible to compare the obtained results on a larger scale, i.e., across the whole genus.
With this study I attempt to provide a coherent and objective approach, applying the Frost and Hillis subspecies concept to categorize the taxa within the genus Phrynosoma into species and subspecies.

Materials and Methods
I used phylogenetic and species delimitation approaches based on four mtDNA markers (12S, ND1, ND2, and ND4) to generate information on the degree of divergence of the studied taxa. A survey on Genbank showed that Phrynosoma sequences of the following markers were available for download (number of sequences for each marker in parentheses): ND1 (89), ND2 (86), ND4 (374), 12S (137), 16S (9), COI (34), and cytB (118), plus a few more markers with very small sample sizes. The DNA sequences used in this study were mostly downloaded from GenBank, except for a set of ND1 and ND2 sequences used in Leaché et al. (2009) that I obtained from Adam D. Leaché. I aligned the sequences with MUSCLE [28] using the default settings in Geneious 6.1.2. [29]. For software applications, sequence data formatting was converted using the online server Alter [30]. The best substitution model for each gene was identified using PartitionFinder2 [31], with linked branch lengths via PhyML 3.0 [32]. Model selection used the corrected (for finite sample size) Akaike information criterion (AICc) [33]. Bayesian inference (BI) analyses used MrBayes 3.2 [34,35] with five runs with eight chains. The first 25% of trees were discarded as burn-in. MCMC runs used an initial set of 1,000,000 generations, with sampling every 500 generations and addition of 500,000 generations until chains reached convergence. Convergence was considered achieved when the standard deviation of split frequencies was 0.015 or less. Additionally, convergence was diagnosed with the potential scale reduction factor (PRSF), which should approach 1.0 as runs converge [36]. I used MEGA 7.0.26 [37] to run maximum likelihood (ML) analyses using a bootstrap estimation method of 10,000 replications. I used FigTree 1.3.1 (http://tree.bio.ed.ac.uk/software/figtree/; accessed on 4 February 2021) to view the trees. I estimated evolutionary genetic divergence, computing uncorrected pairwise distances with MEGA 7.0.26 [37] to assess the degree of intra-and interspecific differences, using a bootstrap estimation method of 10,000 replications. I built a species tree using BEAST 2.4.7 [38] with 1,000,000 generations for the MCMC model. The resulting tree was visualized using DensiTree 2.2.6 [39]. Automatic barcode gap discovery (ABGD) [40] was run through its webserver (http://wwwabi.snv.jussieu.fr/public/abgd/abgdweb.html; accessed on 16 February 2021), using simple distance and default values for the prior intraspecific divergence, except for the relative gap width which I set to an intermediate value of 1.0. All trees were rooted on the Phrynosoma asio branch since this species was found to be the sister taxon to the rest of the species in this genus [9,41].
The genetic markers were analyzed separately in order to obtain a maximum coverage of populations to be treated equally in a given analysis. In concatenated data sets, there are always voucher specimens that are not represented by the full set of genetic markers.
To create distribution maps for the recognized taxa, I downloaded the available Phrynosoma records from VertNet (www.vertnet.org; accessed on 18 February 2021) and GBIF (GBIF.org; GBIF Occurrence Download: https://doi.org/10.15468/dl.mpyes6; accessed on 17 February 2021). Also, I received detailed locality information for the Phrynosoma specimens housed in the following institutions: AMNH, ANSP, ROM, SMF, UMMZ, USNM, and UTA. Finally, I included the locality data of Phrynosoma species obtained from the Geographic Distribution section of Herpetological Review. The maps were created using ArcMap 10.4. Records that were obviously out of the distributional range of the respective species were suspected to be based either on misidentified specimens or on horned lizards that had been released in areas away from their natural range of occurrence, and thus these records were omitted from the maps.

Results
The final alignments were 780 (12S), 975 (ND1), 1033 (ND2), and 839 (ND4) nucleotide positions long and contained 137 (12S), 89 (ND1), 86 (ND2), and 373 (ND4) samples. Partition schemes were as follows: 12S (GTR+I+G), ND1 (first codon position TVM+I+G; second HKY+I+G; and third GTR+G), ND2 (first codon position GTR+I+G; second HKY+I+G; and third GTR+I+G), and ND4 (GTR+I+G). The trees obtained by BI, ML, and BEAST analyses show a high degree of congruence at well-supported nodes, with some differences in branch arrangement at poorly supported nodes (Figures 1-5). The analyses based on ND4 covered most taxa currently recognized in the genus Phrynosoma except for those of the coronatum complex. The latter were well-represented by the other three markers. Taxa that were represented by more than one marker were recovered, indicating the same levels of genetic differentiation and phylogenetic position among the various markers.      The analyses recover all currently recognized taxa as monophyletic groups. For the assignment of clades to taxa in the P. orbiculare complex I follow the authors of [7], except that I combined their X and XI clades, since these had low bootstrap values in the ML analyses and were not recovered as separate clades in the BI analyses ( Figure 1). For the species related to douglasii I follow the authors of [8], except for naming their "GB/CP" clade "ornatum", since the type locality of the latter (see checklist below) coincides very closely with one of the genetically tested samples of this clade. The P. hernandesi clade includes specimens from the vicinities of the type localities of the following nominal taxa, which are therefore considered synonyms of Phrynosoma hernandesi Girard 1858: Phrynosoma ornatissima Girard 1858, Phrynosoma brevirostris Girard 1858, Phrynosoma douglassii brachycercum Smith 1942, Phrynosoma bauri Montanucci 2015, and Phrynosoma diminutum Montanucci 2015. Among the samples of cornutum, the two clades (eastern and western) found by the authors of [6] were also recovered by my analyses. I apply the name P. c. bufonium Wiegmann, 1828 to the western clade since this is the oldest available name for this clade. I recovered the same four clades in the coronatum complex as did the authors of [5]. However, because of its distinctness and high statistical support, I kept the "Northern Baja California" (terminology from [5]) clade separate from the combined "Northern California" and "Southern California" clade. Due to their respective type localities, I apply the names P. c. frontale Van Denburgh, 1894 to the California clade, P. c. blainvillii Gray, 1839 to the "Northern Baja California" clade, P. c. cerroense Stejneger, 1893 to the "Central Baja California" clade, and P. c. coronatum (Blainville, 1835) to the "Southern Baja California" clade. Finally, in my analyses I referred to the platyrhinos-like samples from Yuma Proving Grounds as "P. cf. platyrhinos". However, since no name is available for this population and its taxonomic status needs further evaluation, I do not consider it here any further.
Once combined to form respective single units in the analyses, these species complexes and sister taxa pairs have mean genetic distances from any other taxon in the ranges of 5.3-10.7% (12S), 9.7-15.6% (ND1), 11.2-15.9% (ND2), and 10.0-20.9% (ND4), and thus would fall into my categories A and B.
I provide the following checklist of the species and subspecies that I recognize. Geographic distribution: Northern Nayarit along the Pacific versant of Mexico to at least the Isthmus of Tehuantepec, at elevations from near sea level to 2500 m asl; there seem to be no recent verified records of this species east of the Isthmus of Tehuantepec ( Figure 6). Historically, the species was possibly present in the Rio Grijalva Valley of Chiapas and adjacent Guatemala [44] (J. A. Campbell, personal communication, 23 February 2021). According to A. Ramirez (personal communication, 19 March 2021), "there is no evidence for the ancient or current existence of P. asio in Chiapas. Dr. Álvarez del Toro referred to it in the state of Chiapas, but the data were not consistent. I personally visited the locality for more than ten years without corroborating the presence of the species. However, I was able to verify that the Phrynosoma specimens that were donated to Dr. Álvarez del Toro were purchased at the Isthmus of Tehuantepec, in the neighboring state of Oaxaca".

Checklist of the Species of
Content: No subspecies are recognized.
Geographic distribution: From southeastern Arizona and northeastern Sonora across New Mexico to eastern Texas and from Kansas southward into Mexico as far as Durango, Zacatecas, Nuevo León, and extreme northern Veracruz, from near sea level to 1865 m asl ( Figure 8). There are scattered records from northern Colorado, southern Nebraska, and southwestern Arkansas.  Geographic distribution: Definitely known (based on genetically confirmed records) from most of Texas, except for the extreme western portion, as well as from southeastern Colorado, from near sea level to 1485 m asl ( Figure 8); based on this distributional pattern, it is likely that the populations in the states of Oklahoma, Kansas, Nebraska, and Arkansas belong to this subspecies. The subspecific identity of the populations of the Mexican states of Coahuila, Nuevo León, and Tamaulipas needs to be clarified through genetic analyses.
Geographic distribution: Definitely known (based on genetically confirmed records) from southeastern Arizona and southwestern New Mexico to extreme western Texas, as well as in the Mexican states of Chihuahua and Sonora, from 390 to 1865 m asl (Figure 8).
Geographic distribution: Northern California to the Cape region of Baja California, Mexico, at elevations from near sea level to 1655 m asl (Figure 9). frontale; orange circles = genetically confirmed records of P. c. blainvillii; green squares = genetically confirmed records of P. c. cerroense; red circles = genetically confirmed records of P. c. coronatum; gray circles = records that have not been genetically confirmed; blue star = type locality of P. frontale; orange star = restricted type locality of P. blainvillii; gray star = type locality of P. ochoterrenai; green star = type locality of P. cerroense and P. schmidti; purple star = type locality of P. jamesi; black star = type locality of P. wigginsi; red star = restricted type locality of Agama coronata.
Geographic distribution: San Diego region and the northern portion of Baja California, Mexico; from about 31.8 • to 33.1 • latitude, at elevations from near sea level to 1135 m asl ( Figure 9).  (Figure 9).
Geographic distribution: Northern California southward to about 33.1 • latitude, about 30 km north of San Diego, at elevations from near sea level to 1655 m asl ( Figure 9).
Comments: In the original description of this taxon, Bell [62] used two different spellings of the specific epithet: in the main text and formal introduction of the species, the spelling is with a double "s" (i.e., Agama Douglassii) whereas he spelled the name with a single "s" in the labeling of his figure (i.e., Agama Douglasii). Hammerson and Smith [65], acting as First Revisor according to Art. 24.2. ICZN, fixed the name with a single "s". As the correct original spelling ends with "-ii", the spelling with a single "i" is considered to be an incorrect subsequent spelling (Art. 33.4 ICZN). Therefore, the correct scientific name of this species is Phrynosoma douglasii. See further comments in support of this spelling by Montanucci [2] (p. 45). The variant of the name spelled with a single "i", as proposed by Hammerson   Blue squares = genetically confirmed records of P. douglasii; green circles = genetically confirmed records of P. hernandesi; orange circles = genetically confirmed records of P. h. ornatum; red squares = genetically confirmed records of P. h. ditmarsi; gray circles = records that have not been genetically confirmed; blue star = type locality of P. pygmaea; gray star = restricted type locality of P. exilis; orange star = type locality of P. ornatum; green star = restricted type locality of P. hernandesi; purple star = restricted type locality of P. brevirostris; black star = restricted type locality of Phrynosoma douglassii brachycercum; red star = restricted type locality of P. ornatissima; yellow star = type locality of P. bauri; white star = type locality of P. diminutum. The type localities of P. douglasii and P. ditmarsi are too imprecise to plot. Geographic distribution: Definitely known (based on genetically confirmed records) from southern Canada (southeastern Alberta, southwestern Saskatchewan) to southeastern Arizona, New Mexico, and western Texas, at elevations of 600 to 2350 m asl (Figures 10 and 11); the taxonomic identity of the Mexican populations needs to be clarified. Figure 11. Map showing distributional records of Phrynosoma h. ditmarsi and P. h. hernandesi in the area of contact. Green circles = genetically confirmed records of P. hernandesi; red squares = genetically confirmed records of P. h. ditmarsi; gray circles = records that have not been genetically confirmed; green star = restricted type locality of P. hernandesi. The type locality of P. ditmarsi is too imprecise to plot.
Geographic distribution: Definitely known (based on genetically confirmed records) from the northern portion of the state of Sonora, Mexico, at elevations from 1050 to 1545 m asl (Figures 10 and 11). The true geographic distribution of P. h. ditmarsi is still poorly known.

Phrynosoma mcallii (Hallowell 1852) [53]
Geographic distribution: Highlands of Mexico from Chihuahua and eastern Sonora to central Veracruz and northern Oaxaca, at elevations from 800 to 3190 m asl (Figures 14 and 15). o. orbiculare; symbols with black edges (except for stars) are localities from where specimens have been included in the analyses; gray circles = records that have not been genetically confirmed; white star = restricted type locality of Lacerta orbicularis and P. wiegmanni; green star = type locality of P. o. bradti; blue star = type locality of P. o. durangoensis; orange star = type locality of Tapaya dugesii; purple star = type locality of P. o. orientale; yellow star = type locality of P. o. alticola; black star = type locality of Tapaya cortezii; red star = restricted type locality of Tapaya boucardii. The type locality of Tapaya orbicularis longicaudatus is too imprecise to plot. Phrynosoma orbiculare orbiculare (Linnaeus 1758) [71] Phrynosoma wiegmanni Gray 1839 [55] (p. 96). Type locality: "Mexico"; restricted to "México City, Distrito Federal" by Smith and Taylor [49]  Geographic distribution: Definitely known (based on genetically confirmed records) from Distrito Federal and the central and southern portions of the state of México, at elevations from 2020 to 3190 m asl (Figures 13 and 14); the taxonomic identities of the populations in the northern and eastern portions of the state of México, as well as of those in Tlaxcala and western Puebla, need to be clarified.
Geographic distribution: Populations in the south-central portion of the Mexican state of Hidalgo are currently assigned to this subspecies, at elevations from 800 to 2390 m asl ( Figures 14 and 15), although no genetically verified samples have been studied from near the restricted type locality (in northwestern Hidalgo) of this taxon.
Comments: For Bocourt as the single author of this taxon, see comments in the account for Phrynosoma braconnieri.
Geographic distribution: Highlands of Chihuahua and eastern Sonora, Mexico, at elevations from 1165 to 3165 m asl (Figures 14 and 15). Geographic distribution: Definitely known (based on genetically confirmed records) from a broad belt across north-central Mexico from Nuevo León, extreme southwestern Tamaulipas, and southeastern Coahuila across northern Hidalgo, Querétaro, San Luis Potosí, Aguascalientes, and Zacatecas to central Jalisco and northern Nayarit, at elevations from 1045 to 2600 m asl (Figures 14 and 15). The taxon's type locality, "Colima", should probably be understood as the state of Colima, as the city of the same name is situated at around 500 m elevation, too low for any member of the Phrynosoma orbiculare complex.
Comments: For Bocourt as the single author of this taxon, see comments in the account for Phrynosoma braconnieri.
Geographic distribution: According to the genetically confirmed records, this subspecies is widely distributed in the Mexican state of Durango and adjacent Nayarit and Zacatecas, at elevations from 800 to 2390 m asl (Figures 14 and 15). The taxonomic identities of the populations in northern Jalisco need to be clarified.
Geographic distribution: From southern Oregon and Idaho southward along eastern California, Nevada, Utah, and Arizona to northern Baja California and Sonora, Mexico ( Figure 16). platyrhinos; green circles = genetically confirmed records of P. p. goodei; yellow square = genetically confirmed records of P. cf. platyrhinos (Yuma Proving Grounds, California); gray circles = records that have not been genetically confirmed; black star = type locality of P. platyrhinos; blue star = type locality of Anota calidiarum; green star = restricted type locality of P. p. goodei.
Geographic distribution: Central and southern portions of Arizona and southward across Sonora to central Sinaloa, Mexico, from near sea level to about 1700 m asl ( Figure 17).  Comments: For Bocourt as the single author of this taxon, see comments in the account for Phrynosoma braconnieri.
Content: Two subspecies are recognized, P. t. taurus and P. t. sherbrookei.

Phrynosoma taurus taurus Bocourt 1870 [45]
Geographic distribution: Highlands of southern Mexico in the states of Morelos, Puebla, and Oaxaca and the Distrito Federal, as well as disjunct in east-central Guerrero, at elevations from 630 to 2440 m asl ( Figure 18). Geographic distribution: Definitely known (based on genetically confirmed records) from central northeastern Guerrero, Mexico, at elevations from 980 to 1795 m asl ( Figure 18).

Discussion
There is nothing wrong with being a subspecies! Except for parthenogenetic species, speciation is a gradual process and thus, along the continuum from the beginning of the process (i.e., 0% speciation) to fully fledged species (i.e., 100% speciation), the diverging populations are in an intermediate state of the speciation process; thus, in a kind of gray zone for taxonomists. Of course, the process of speciation is not inexorable and such diverging populations can come back into contact and merge again. Typically, the populations in this gray zone have a turbulent taxonomic history, as exemplified by the Phrynosoma coronatum complex. The taxonomy of this group of horned lizards has been revised at least 20 times in the past 180 years, and the various authors have recognized, for instance, between one and seven taxa among the coast horned lizards (see [10] and references therein). Thus, the persistence of this incongruence among some authors over many decades is an indication that the populations involved have not yet reached the 100% speciation level, and thus most probably qualify for an infraspecific classification. In the past, various terms have been proposed for these intraspecific lineages, such as evolutionarily significant units (ESUs; [82,83]). However, I agree with Hillis [22] (p. 54), who suggested that "subspecies designations are ideal for this application". The subspecies is a real evolutionary unit-a population or set of populations that has diverged to a certain degree from other such populations but has not yet reached a level of divergence typical for most of the taxa in the genus that are unambiguously recognized as species level units.
One might argue that "devaluing" a species level taxon to subspecies level would result in an underestimation of diversity. On the contrary: neglecting infraspecific recognition of historical phylogenetic lineages incurs the risk that these effectively disappear from the taxonomic system. One example is the recognition of taxa in the Phrynosoma coronatum complex. Obviously, there is a distinct genetic structure in the populations that compose this species complex, with four monophyletic genetic groups, although the mean genetic distances are relatively small; much smaller than among fully differentiated species such as P. asio, P. modestum, P. solare, P. cornutum, P. mcallii, and P. braconnieri. Thus, a perfect case for the use of the subspecies category proposed by Frost and Hillis [21]! As subspecies, these intraspecific phylogenetic lineages are recognized as a Linnaean category and thus made "visible" for political and conservation issues. Not using the subspecies category for such phylogenetic lineages can lead to them being ignored, as done by, e.g., Brattstrom [27,84], who recognized a single monotypic species, P. coronatum, in this species complex, thereby indeed underestimating diversity! If one were to treat all these "species in the making" (and yes, they could just as easily merge in the future) as full species level units, then this would result in an unbalanced system of full species and those populations that have diverged to a much lesser degree, thereby combining subjectively units of two different levels in a single category. By recognizing those populations that have diverged to a much lesser degree as subspecies, one has a taxonomic system with more phylogenetic information because then what belongs together-as a species level unit-is grouped together as subspecies.
One might ask, "why only use mtDNA data?" Well, most of all, this is because these data are available for more than 500 specimens of Phrynosoma, covering a huge portion of the geographic range of the genus. Also, mtDNA data seem to suit this purpose (i.e., evaluating the levels of genetic divergence among populations) very well because they evolve fast enough to identify younger splits. Nuclear DNA data and genomic datasets are great for phylogenetic reconstructions. However, the commonly used nuclear markers mostly evolve much more slowly and thus do not yield the resolution needed to identify historical phylogenetic lineages below the species level. Furthermore, genomic data are available only for a dozen or so individuals, i.e., they do not represent the geographic variation of the populations over the entire distribution area. The latter is necessary for an assessment of the possible existence of historical phylogenetic lines within a "species". The mtDNA data depict the genetic differentiation of the Phrynosoma populations very well, have a high resolution, and practically cover the entire spectrum of published names-with very few exceptions, but perfection does not exist.
Morphology is important-no doubt about it. However, morphological features are assessed much more subjectively than molecular genetic data. For example, look again at the Phrynosoma coronatum complex: Reeve [11] recognized six taxa in this complex, two of which (P. coronatum and P. cerroense) he recognized as species level units and the rest of which he recognized as subspecies of P. coronatum. Brattstrom [84] examined the "diagnostic" characteristics of the taxa and came to the conclusion that only a single monotypic species should be recognized in this complex. Montanucci [3] examined 24 "potentially informative characters" and analyzed them using principal component analysis. Taken individually, none of the traits he examined are diagnostic. Rather, all characters show broad overlap between the four species he distinguished in the P. coronatum complex. The only taxon that shows partially non-overlapping character states compared to the other "species" is P. wigginsi. However, the sample size for this taxon was very small, with only six or seven specimens (for the relevant characters). Additional specimens would most likely increase the documented range of variation and thus result in an overlap of the range with the other species. Brattstrom [27] was apparently not impressed by Montanucci's [3] taxonomic conclusions either, as he only recognized a single monotypic species in this complex: Phrynosoma coronatum.
After 25 years of experience working with a multitude of taxonomic issues involving salamanders, frogs, lizards, and snakes, and applying an integrative taxonomic approach including external morphological, genital morphological, and molecular genetic (both mtDNA and nuclear DNA, and recently genomic, data) features, as well as bioacoustics for frogs, I am convinced that, for the "first sorting" of the operational taxonomic units (OTUs) above all, mtDNA data are the most useful (but also bioacoustical and genital morphological data, if they are available), and only subsequently does the external morphology contribute to the understanding of the differentiation of the a priori sorted groups (i.e., the OTUs).
Incidentally, the phylogentic trees based on genomic data (e.g., [9]) are largely congruent with the trees based on mtDNA data presented here: the P. orbiculare and P. douglasii species complexes form a clade (Tapaya subgenus) in the analyses of both datasets, as do P. p. platyrhinos and P. p. goodei (Doliosaurus subgenus); P. t. taurus and P. t. sherbrookei (Brevicauda subgenus); and P. c. coronatum, P. c. cerroense, and P. c. blainvillii (Anota subgenus). However, there is some incongruence regarding the deeper splits, as can be expected. Finally, the genomic data show the same gradation in genetic differentiation as the mtDNA data, i.e., those taxa treated here as subspecies have diverged significantly less from one another than, e.g., P. asio, P. modestum, P. solare, P. cornutum, P. mcallii, and P. braconnieri. However, combined to represent a unit at species level, they show a similar differentiation as these well-separated species.
Typically, those taxa recognized as subspecies together form a monophyletic group (that represents the species) and are distributed allopatrically. In my approach, all separate evolutionary significant units are recognized as named taxa-either species or subspeciesthereby reflecting the importance of identifying and naming such important units for conservation.
Funding: This research received no external funding.
Institutional Review Board Statement: Ethical review and approval were waived for this study, due to the exclusive use of DNA sequence data downloaded from Genbank.

Data Availability Statement:
No new data were created or analyzed in this study. All analyses are based on DNA sequence data downloaded from Genbank. Data sharing is not applicable to this article.
Ontario Museum (ROM), Toronto; Greg Schneider, Museum of Zoology, University of Michigan (UMMZ), Michigan; Esther M. Langan, Smithsonian Institution, National Museum of Natural History (USNM), Washington; and Gregory G. Pandelis and Jonathan A. Campbell, Amphibian and Reptile Diversity Research Center, University of Texas at Arlington (UTA), Arlington. I thank my colleague at Senckenberg, Michael Oertel, for support at various levels during the preparation of the maps for this work. Jonathan A. Campbell (Arlington, TX, USA) and Antonio Ramirez (Tuxtla Gutiérrez, Chiapas, Mexico) shared information on the possible occurrence of Phrynosoma asio in Chiapas with me. I am grateful to Aaron Bauer (Villanova, PA, USA) for helpful comments on various nomenclatural issues. Finally, I am grateful to Helmut Diethert and three anonymous reviewers for constructive criticism and helpful comments on an early draft of this article.

Conflicts of Interest:
The author declares no conflict of interest.