Quantifying ecological, morphological, and genetic variation to delimit species in the coast horned lizard species complex (Phrynosoma)

Abstract

Lineage separation and divergence form a temporally extended process whereby populations may diverge genetically, morphologically, or ecologically, and these contingent properties of species provide the operational criteria necessary for species delimitation. We inferred the historical process of lineage formation in the coast horned lizard (Phrynosoma coronatum) species complex by evaluating a diversity of operational species criteria, including divergence in mtDNA (98 specimens; 2,781 bp) and nuclear loci (RAG-1, 1,054 bp; BDNF 529 bp), ecological niches (11 bioclimatic variables; 285 unique localities), and cranial horn shapes (493 specimens; 16 landmarks). A phylogenetic analysis of mtDNA recovers 5 phylogeographic groups arranged latitudinally along the Baja California Peninsula and in California. The 2 southern phylogeographic groups exhibit concordance between genetic, morphological, and ecological divergence; however, differentiation is weak or absent at more recent levels defined by phylogeographic breaks in California. Interpreting these operational species criteria together suggests that there are 3 ecologically divergent and morphologically diagnosable species within the P. coronatum complex. Our 3-species taxonomic hypothesis invokes a deep coalescence event when fitting the mtDNA genealogy into the species tree, which is not unexpected for populations that have diverged recently. Although the hypothesis that the 3 phylogeographic groups distributed across California each represent distinctive species is not supported by all of the operational species criteria evaluated in this study, the conservation status of the imperiled populations represented by these genealogical units remains critical.

Fig. 1.

Mitochondrial DNA genealogy (A), nuclear allele networks (B), and multilocus nuclear network showing the genetic distance among specimens (C) for the P. coronatum species complex. The mtDNA genealogy is based on a partitioned Bayesian analysis of 2,781 base pairs, and asterisks denote major clades with posterior probability values ≥0.99. Allele networks for the RAG−1 and BDNF nuclear loci are based on statistical parsimony with a 95% connection significance. Specific locality data are provided in Table S5, and a detailed mtDNA genealogy and map illustrating the collecting locality of each sample is provided in Fig. S1.

Fig. 2.

Multivariate plots of the climate variables (A) and geometric morphometrics (B) for the P. coronatum complex. The climate data plots include 1,000 background points (shown in gray), drawn at random from the distribution of the P. coronatum complex (A). Variable loadings for the first, second, and third components for the climate data, respectively, are: mean temperature of coldest quarter, mean temperature of warmest quarter; temperature annual range, mean diurnal range; precipitation of driest month, precipitation of warmest quarter.

Fig. 3.

Geographic variation in cranial horn shapes. (A) Cranial horn shape divergence and geometric morphometric landmark positions between P. coronatum complex specimens from Contra Costa County, California (Left; MVZ 33623) and Southern Baja California (Right; MVZ 100473). (B) The Southern Baja California phylogeographic group differs notably from all others in cranial horn shape. Mean landmark positions (average Bookstein coordinates; CS = 1) for each phylogeographic group are shown in relation to the parietal eye (x = 0, y = 0). Males and females produced similar results, and only results for females are shown here.

Fig. 4.

Proposed species tree for the P. coronatum species complex and the contained mtDNA genealogy composed of 5 phylogeographic groups. The failure of all mtDNA haplotype lineages within P. blainvillii to coalesce between T1 and the present results in an opportunity for deep coalescence between P. blainvillii and P. cerroense.