Species delimitation and introgression analyses
Of the five species models tested in BFD*, the “4 species” model (S. zerene, S. atlantis, northern S. hesperis , and southern S. hesperis ) had the highest marginal likelihood estimate (MLE) and the most strongly negative Bayes Factor (BF), indicating that this model was the best supported by the data (Table 1). The maximum clade credibility (MCC) tree indicated strong support for the northern S. hesperis and S. atlantis clade (posterior probability of 1), but the clade containing S. zerene and southern S. hesperis was less supported with a posterior probability of 0.87 (Fig. 3a, left panel). The MCC tree additionally indicated that the southern S. hesperis/S. zerene clade diverged slightly earlier than the clade containing S. atlantis and northern S. hesperis .
DensiTree visualization of the BFD* results indicated discordance in the relationship between S. zerene and the remaining three species on the tree. The major recovered topology, which accounted for 87.4% of the sampled trees (shown in blue in Fig. 3a, right panel), was the same as the topology presented in the MCC tree, however the second most common topology (7.1% of the sampled trees, shown in purple) depictedS. zerene as basal to the other three species. A third topology accounting for the remaining 5.4% of the sampled trees (shown in orange) depicted S. zerene as the sister taxon to the S. atlantis/ northern S. hesperis clade. The S. atlantis /northern S. hesperis sister relationship was consistently recovered in each sampled topology.
Given the discordant relationship between S. zerene and S. atlantis/S. hesperis in our phylogenetic and BFD* analyses, we used TreeMix to test for introgression between these species. The unrooted maximum likelihood phylogeny produced in TreeMix (Fig. 3b, Fig. S2) was largely consistent with the SNP phylogeny in Fig. 1. Heatmaps containing the pairwise population residuals for each migration model generally indicated improved fit with the addition of migration events (Fig. S2), but this seemed to plateau for models m 3-5; the residual plots for three, four, and five migration events were highly similar, suggesting that continued addition of migration events past this point did not greatly improve fit. The m 3 model indicated a substantial migration event from the interior of the branch containing S. zerene to the central S. hesperis population (Fig. 3b; migration weight = 0.41, p = 0); a migration edge plotted along the interior of the branch rather than the terminus may indicate historical admixture or admixture from unsampled populations (Pickrell & Pritchard 2012). This model also indicated two less substantial, but still statistically significant, migration events between S. zerene and the clade containing the central and New Mexico/Arizona S. hesperispopulations (migration weight = 0.07, p = 0.037), and betweenS. zerene and S. atlantis (migration weight = 0.07,p = 0.017). Though models m 4 and m 5 had similar residuals (Fig. S2), they indicated non-significant migration weights between some populations: in m 4 the migration edge between the New Mexico/Arizona S. hesperis population and S. atlantiswas non-significant (migration weight: 0.02, p = 0.169), and inm 5 the edges between S. zerene and the clade containing the New Mexico/Arizona and central S. hesperis populations (migration weight = 0.09, p = 0.058), between the New Mexico/Arizona S. hesperis populations and S. atlantis(migration weight = 0.02, p = 0.092), and between the central and northern S. hesperis populations (migration weight = 0.01,p = 0.345) were non-significant.
Of the 61 f 3 tests for introgression that we computed, only three were significant (Table S2); one indicated admixture from S. zerene and the southern Utah population of S. hesperis into the central S. hesperis population (p = 0), and another fromS. zerene and the New Mexico/Arizona population of S. hesperis into the central S. hesperis population (p = 0.003). The third test suggested gene flow from the New Mexico/Arizona and northern S. hesperis populations into the central S. hesperis population (p = 0.035).