Understanding species limits through the formation of phylogeographic lineages

Abstract

The outcomes of speciation across organismal dimensions (e.g., ecological, genetic, phenotypic) are often assessed using phylogeographic methods. At one extreme, reproductively isolated lineages represent easily delimitable species differing in many or all dimensions, and at the other, geographically distinct genetic segments introgress across broad environmental gradients with limited phenotypic disparity. In the ambiguous gray zone of speciation, where lineages are genetically delimitable but still interacting ecologically, it is expected that these lineages represent species in the context of ontology and the evolutionary species concept when they are maintained over time with geographically well-defined hybrid zones, particularly at the intersection of distinct environments. As a result, genetic structure is correlated with environmental differences and not space alone, and a subset of genes fail to introgress across these zones as underlying genomic differences accumulate. We present a set of tests that synthesize species delimitation with the speciation process. We can thereby assess historical demographics and diversification processes while understanding how lineages are maintained through space and time by exploring spatial and genome clines, genotype-environment interactions, and genome scans for selected loci. Employing these tests in eight lineage-pairs of snakes in North America, we show that six pairs represent 12 "good" species and that two pairs represent local adaptation and regional population structure. The distinct species pairs all have the signature of divergence before or near the mid-Pleistocene, often with low migration, stable hybrid zones of varying size, and a subset of loci showing selection on alleles at the hybrid zone corresponding to transitions between distinct ecoregions. Locally adapted populations are younger, exhibit higher migration, and less ecological differentiation. Our results demonstrate that interacting lineages can be delimited using phylogeographic and population genetic methods that properly integrate spatial, temporal, and environmental data.

Keywords: genome cline; genome selection; hybrid zone; spatial cline; speciation.

FIGURE 1

(a) Genomic admixture maps showing the locations of lineages within species pairs given clusters and admixture estimates using the SuperSom (self‐organizing maps) method with interpolated contour lines to reveal the center of the hybrid zone (purple) at the intersections of ecoregions and cline estimates using admixture frequency estimates for Pantherophis alleghaniensis/P. quadrivittatus (Pa), (b) Pantherophis slowinskii/P. emoryi (Ps), (c) P. meahllmorum/P. emoryi (Pe), (d) Lampropeltis triangulum/L. gentilis (Lt), (e), Crotalus atrox E/C. atrox W (Ca), (f) Lampropeltis splendida/ L. californiae (Ls). (g) Map showing the intersection between the following ecoregions: 1—southeastern coastal plains (blue) and the southeastern plains (aqua), 2—temperate (olive) and semiarid prairies (tan), 3—southeastern plains/Texas–Louisiana plains (dark green) and the semiarid prairies (tan), 4—semiarid prairies (tan) and Tamailipas–Texas semiarid plain (yellow), and 5—Chihuahuan Desert (orange) and Sonoran Desert (Red) at the Western Continental Divide/Cochise Filter Barrier (brown).

FIGURE 2

Best‐fit isolation and migration models using GADMA2 for all species pairs showing distinct clusters from the SuperSOM analysis. Time of divergence and changes in population sizes (width of blue lines) and migration rates (black arrows) are shown for all species pairs.

FIGURE 3

Histograms showing variable importance from general dissimilarity models (GDMs) for each species pair. Asterisks above variables indicate significance.

FIGURE 4

(a) Genome clines for all loci for each species pair. Gray lines indicate loci not showing a significant slope (v) > 1. Colored solid lines indicate loci with significant slope and colored dashed lines indicate the location along the genome‐wide hybrid index. (b) Proportion of significant genome cline loci as a percentage of all loci sampled for each species pair. (c) Correlation between the total number of loci sampled across all pairs and the number of significant loci discovered using the genome cline method.

FIGURE 5

Venn diagrams showing the number (and percentage) of loci that overlap among the following categories: RDA—environmental‐locus correlations, significant genome clines, significant genome scans, and 1% of loci showing contributions to the principal components for identifying clusters in DAPC >99%.