Standing geographic variation in eclosion time and the genomics of host race formation in Rhagoletis pomonella fruit flies

Abstract Taxa harboring high levels of standing variation may be more likely to adapt to rapid environmental shifts and experience ecological speciation. Here, we characterize geographic and host‐related differentiation for 10,241 single nucleotide polymorphisms in Rhagoletis pomonella fruit flies to infer whether standing genetic variation in adult eclosion time in the ancestral hawthorn (Crataegus spp.)‐infesting host race, as opposed to new mutations, contributed substantially to its recent shift to earlier fruiting apple (Malus domestica). Allele frequency differences associated with early vs. late eclosion time within each host race were significantly related to geographic genetic variation and host race differentiation across four sites, arrayed from north to south along a 430‐km transect, where the host races co‐occur in sympatry in the Midwest United States. Host fruiting phenology is clinal, with both apple and hawthorn trees fruiting earlier in the North and later in the South. Thus, we expected alleles associated with earlier eclosion to be at higher frequencies in northern populations. This pattern was observed in the hawthorn race across all four populations; however, allele frequency patterns in the apple race were more complex. Despite the generally earlier eclosion timing of apple flies and corresponding apple fruiting phenology, alleles on chromosomes 2 and 3 associated with earlier emergence were paradoxically at lower frequency in the apple than hawthorn host race across all four sympatric sites. However, loci on chromosome 1 did show higher frequencies of early eclosion‐associated alleles in the apple than hawthorn host race at the two southern sites, potentially accounting for their earlier eclosion phenotype. Thus, although extensive clinal genetic variation in the ancestral hawthorn race exists and contributed to the host shift to apple, further study is needed to resolve details of how this standing variation was selected to generate earlier eclosing apple fly populations in the North.


| INTRODUC TI ON
The raw material for novel adaptation can come from new mutations or standing variation (Barrett & Schluter, 2008). Adaptation is expected to occur more rapidly when based on standing variation, because the process is not limited by wait time for new beneficial mutations to arise and standing variants are likely to have been filtered for negative pleiotropic fitness effects or detrimental epistatic interactions (Barrett & Schluter, 2008;Yeaman, 2015). Therefore, taxa harboring high levels of standing genetic variation may be more likely to adapt to rapid environmental shifts and experience ecological speciation than those requiring new favorable mutations to respond. Indeed, myriad examples are accumulating of rapid adaptive evolution fueled by standing variation, including coat color in the old field mouse, Peromyscus polionotus (Steiner, Weber, & Hoekstra, 2007), reduced defensive armor in the threespine stickleback fish, Gasterosteus aculeatus (Colosimo et al., 2005), mimicry in Heliconius butterflies (The Heliconius Genome Consortium, 2012), and beak morphology in Darwin's Finches (Lamichhaney et al., 2015). Standing variation may be especially relevant in cases where adaptation involves polygenic traits, and the majority of favorable alleles have only minor effects on fitness and genomic incompatibility (Barrett & Schluter, 2008;Hermisson & Pennings, 2005;Schluter & Conte, 2009). Large stores of standing genetic variation observed in some populations may be a product of complex evolutionary histories, including past gene flow coupled with ecological selection to maintain high levels of genetic variation resulting from tracking local, regional, or temporal differences in environmental or biotic conditions (Bergland, Tobler, González, Schmidt, & Petrov, 2016;Berner & Salzburger, 2015;Brawand et al., 2014;Feder, Berlocher, et al., 2003;Gompert, Fordyce, Forister, Shapiro, & Nice, 2006;Jeffery et al., 2017;Jones et al., 2012;Lamichhaney et al., 2015;Loh et al., 2013;Pease, Haak, Hahn, & Moyle, 2016;Roesti, Gavrilets, Hendry, Salzburger, & Berner, 2014;The Heliconius Genome Consortium, 2012).
In addition, hybridization need not have an immediate effect on generating new taxa (Abbott, Barton, & Good, 2016) but can also create and maintain extensive standing variation, enabling divergence at a later time when populations experience new ecological opportunities (Berner & Salzburger, 2015).
Rhagoletis pomonella (Diptera: Tephritidae), a model for population divergence in sympatry, provides an avenue for investigating the role of standing clinal variation in rapid adaptation and speciation in response to novel ecological opportunity (Berlocher & Feder, 2002). Ancestral R. pomonella infested the fruits of native North American hawthorns (Crataegus spp.) and shifted (<170 years ago) to domesticated apple (Malus domestica) after the plant was introduced by European settlers ~400 years ago (Bush, 1966;Walsh, 1867). phenology, alleles on chromosomes 2 and 3 associated with earlier emergence were paradoxically at lower frequency in the apple than hawthorn host race across all four sympatric sites. However, loci on chromosome 1 did show higher frequencies of early eclosion-associated alleles in the apple than hawthorn host race at the two southern sites, potentially accounting for their earlier eclosion phenotype. Thus, although extensive clinal genetic variation in the ancestral hawthorn race exists and contributed to the host shift to apple, further study is needed to resolve details of how this standing variation was selected to generate earlier eclosing apple fly populations in the North.

K E Y W O R D S
clinal variation, eclosion time, ecological speciation, host races, standing variation Previous studies have implied that hawthorn-infesting populations of R. pomonella possess large stores of phenotypic and genetic variation, notably for life history timing, which may have facilitated the host shift to apple Feder, Berlocher, et al., 2003).
In particular, differences in the timing of adult eclosion have been shown to be an important host-related ecological adaptation contributing to partial allochronic premating isolation between apple and hawthorn flies (Feder et al., 1994;Feder, Hunt, & Bush, 1993).
Fruit on apple varieties favorable for larval survivorship ripens about 3-4 weeks earlier than those of downy hawthorn (C. mollis), the primary host of R. pomonella in the Midwestern and Northeastern United States (Feder et al., 1994). Flies must synchronize breakage of pupal diapause and adult eclosion with the availability of ripe host fruit for mating and oviposition. This is critical because Rhagoletis is univoltine, adults take a week to reach sexual maturity, and flies live for a maximum of 1 month in nature (Dean & Chapman, 1973).
As a result, apple flies have evolved to eclose an average of 10 days earlier than hawthorn flies, as measured in field capture studies, reducing host race overlap to ~80% (Feder et al., 1993(Feder et al., , 1994. This partial allochronic isolation, in combination with host fidelity (e.g., host fruit odor preference), reduces gene flow between the host races to ~4%-6% per generation (Feder et al., 1993(Feder et al., , 1994. Eclosion time appears to have a complex evolutionary history associated with standing clinal variation in the ancestral hawthorn race Feder, Berlocher, et al., 2003;Feder, Chilcote, & Bush, 1990;Michel et al., 2010;Michel, Rull, Aluja, & Feder, 2007). It has been hypothesized that part of the variation originated in the Eje Volcánico Trans Mexicano (EVTM) ~1.5 million years ago, during a period of isolation from other Mexican and North American populations (see Figure 1; Michel et al., 2007;Xie et al., 2007). It then spread through the Sierra Madre Oriental Mountains (SMO) of Mexico and into North America over the past ~1 million F I G U R E 1 Map of the four paired sympatric collection sites for apple and hawthorn flies in the Midwestern United States. Also given are the ranges of the apple race and native hawthorn-infesting populations of Rhagoletis pomonella in the United States and Mexico (see Supporting Information Table S2 for numerical designations of populations and site information). Note that the 430-km transect through the Midwestern United States encompasses much of the latitudinal range of overlap of the apple and hawthorn host races in the region. The primary hawthorn host of R. pomonella in the Northeastern and Midwestern United States is Crataegus mollis. However, moving south from Urbana, apple is not infested and C. mollis becomes rare, although a variety of the species, C. mollis texana, exists in the state of Texas, United States (site 5). Other hawthorn species, with varying fruiting times, are the primary hosts of the fly in the southern United States and Mexico (see Figure 2a and Supporting Information Table S2; Lyons-Sobaski & Berlocher, 2009;Rull et al., 2006). DNA sequencing data imply that hawthorn-infesting R. pomonella populations from the Eje Volcánico Trans Mexicano (EVTM) and those in the Sierra Madre Oriental Mountains of Mexico (SMO) and United States have undergone cycles of allopatry followed by secondary contact and gene flow over the past ~1.5 My Feder, Berlocher, et al., 2003;Michel et al., 2007;Xie et al., 2007), contributing to the creation and maintenance of geographic genetic variation in eclosion time in the fly 0 200 400 600 800 km N years via episodes of secondary contact and introgression between hawthorn fly populations in Mexico and the United States Feder, Berlocher, et al., 2003;Michel et al., 2007;Xie et al., 2007). Subsequently, selection related to latitudinal, altitudinal, and species differences in hawthorn fruiting time in the SMO and United States likely maintained this adaptive variation in eclosion timing.
Consistent with this scenario, many species of hawthorns ripen later in the year with decreasing latitude, reflected in later dates of fly collection (Figure 2a; Lyons-Sobaski & Berlocher, 2009;Rull, Aluja, Feder, & Berlocher, 2006;Xie et al., 2007). As a result, eclosion time varies geographically, with hawthorn flies from further south requiring more time to eclose than those from further north, both in nature and in controlled rearing experiments (Figures 2b and 3a;Dambroski & Feder, 2007;Hood et al., 2015;Lyons-Sobaski & Berlocher, 2009;Xie et al., 2007). This standing geographic variation in eclosion timing of hawthorn flies is thought to have contributed to the adaptive radiation of the R. pomonella sibling species group by allowing these short-lived, univoltine flies to attack novel host plants with differing fruiting times, including the formation of a number of races and potentially species on different hawthorns in southern latitudes (Powell, Cha, Linn, & Feder, 2012;Powell, Forbes, Hood, & Feder, 2014) and most recently the apple race in the Eastern United States (Feder, Chilcote, & Bush, 1988).
Previous studies from our group provided evidence for a link between standing clinal variation and host race formation but had insufficient genomic resolution to infer the genetic architecture of either the underlying eclosion time phenotypes or host race differentiation Feder et al., 1990Feder et al., , 2005Feder, Berlocher, et al., 2003;Michel et al., 2007Michel et al., , 2010Xie et al., 2008Xie et al., , 2007. Ragland et al. (2017) recently laid a genomic foundation for understanding the architecture of eclosion timing and its association with host race differentiation in R. pomonella (see Supporting Information Appendix S1 for additional details). In a genome-wide association study (GWAS) of adult eclosion time in hawthorn and apple flies collected from a field site in Fennville, Michigan (MI), United States, Ragland et al. (2017) compared single nucleotide polymorphism (SNP) allele frequency differences between the earliest (≤3%) and latest (≥97%) eclosing quantiles of flies. They found that SNPs displaying significant allele frequency differences between early-and late-eclosing flies in both host races were concentrated on three of the five major chromosomes (numbers 1-3) constituting the R. pomonella genome (Supporting Information Table S1). Within chromosomes 1-3, SNPs displaying high levels of linkage disequilibrium (LD) with one another, presumably due to chromosomal inversions (Feder, Roethele, Filchak, Niedbalski, & Romero-Severson, 2003b), showed the greatest frequency responses (Table S1; see Supporting Information Appendix S1 for further discussion of inversions and their origins). However, a proportion of low LD SNPs on chromosomes 1-3 also displayed significant responses above random expectation in the GWAS (Supporting Information Table   S1). These low LD SNPs presumably represent loci in more freely recombining, colinear regions of chromosomes also contributing to eclosion timing.
The highly polygenic nature of eclosion time demonstrated by Ragland et al. (2017) provides additional evidence supporting the hypothesis that standing genetic variation, rather than new mutations, likely fueled R. pomonella's shift to the earlier fruiting apple host. However, a stronger case would be made if allele frequencies for SNPs responding in the eclosion time GWAS also showed significant and predictable geographic and host-related variation in nature. Here, we investigate genome-wide patterns of differentiation in Rhagoletis pomonella (Diptera: Tephritidae) for 10,241 SNPs by comparing results from the adult eclosion GWAS (Ragland et al., 2017) with a population survey of four sites distributed from north to south along a 430-km transect across the Midwestern United States, where populations of hawthorn and apple flies co-occur in sympatry (Figure 1). We determine the degree to which the responses of eclosion-associated SNPs from the GWAS predict geographic allele frequency differences within the host races across the Midwest and local, host-related differences between apple and hawthorn fly populations at the four sympatric sites.
In this study, we test four specific predictions of the standing

| Genotyping by sequencing
Methods for genotyping by sequencing (GBS) of individually barcoded double-digest restriction site-associated DNA libraries of flies (ddRAD-seq), de novo assembly of contigs, SNP calling, and allele frequency estimation for field-collected samples at the four sites surveyed in the current study were performed as in Ragland et al. (2017). We used custom scripts and the Genome Analysis Toolkit (GATK version 2.5-2; DePristo et al., 2011) to identify 10,241 variable sites passing quality filters that were genotyped in the eclosion study and at all four paired field sites. Average SNP coverage per individual was 6.2X in the eclosion time GWAS and 3.3X in the geographic survey of sympatric sites.

| Measuring linkage disequilibrium
To assess genome structure, Burrow's composite measures of linkage disequilibrium (∆) standardized to r values between −1 and 1 were estimated following Weir (1979) between pairs of SNPs within hawthorn and apple fly samples at the four collecting sites. Analysis of LD was restricted to a subset of 4,244 of the 10,241 total SNPs genotyped that were previously assigned to one of the five major chromosomes in the genome Ragland et al., 2017).
To further take LD and genome structure into account, we also sepa- Appendix S1 for details). These high LD clusters of linked SNPs presumably represent inversions, with eight different groups identified on chromosome 2 and one each on the remaining four chromosomes.
F I G U R E 2 Geographic variation in host fruiting phenology and fly eclosion time along a latitudinal transect of 18 different Rhagoletis pomonella populations analyzed at 13 sites from Grant, MI, United States, to the highlands of Chiapas, Mexico (MX) (Data are from Dambroski & Feder, 2007;Hood et al., 2015;Lyons-Sobaski & Berlocher, 2009;Xie et al., 2007; see Figure 1 for a map of the 13 sampled sites and Supporting Information Table S2 for population designations and additional information. Sites with the same number but different letters represent flies sampled from different hosts at the site). (a) Host fruiting time, as indicated by the collecting date of fly populations at sites, plotted against latitude; (b) mean time to eclosion measured as the average number of days following postwinter warming of fly populations determined from controlled laboratory rearing studies plotted against latitude; and (c) mean eclosion time of fly populations plotted against host fruiting time at sites. Grant, MI (site 1), and Urbana, IL (site 4), represent the approximate northernmost and southernmost ends, respectively, where apple (green circles) and Crataegus mollis (red triangles) host races geographically overlap in the Midwestern United States. Note that flies from site 5 infest C. mollis v. texana in Texas which is also depicted with a red triangle, while other hawthorn-infesting fly populations are shown in black triangles. Geographic and host-related differentiation in host fruiting time and fly eclosion in the Midwestern United States is subsumed within a larger pattern of variation across North America. There is a general trend for both host fruit to ripen and flies to eclose later in the year with decreasing latitude, as evidenced by significant negative correlations (r values are given considering the hawthorn data alone = r, and for all hosts including apples = r all). The relationship is complicated, however, by R. pomonella attacking different hawthorn species moving southward from the Midwest (see Supporting Information Table S2) that fruit at varying times during the season (Lyons-Sobaski & Berlocher, 2009;Rull et al., 2006;Xie et al., 2007), as well as altitudinal effects in Mexico

| Tests for allele frequency differences
Probabilities of single locus genotypes and allele frequencies for Carlo distribution of absolute differences was taken as significantly different than expected. To determine whether overall percentages of SNPs associated with eclosion time and geographic divergence were significantly greater than expected by chance, the total percentage of SNPs showing significant differences in each simulation was calculated and the distribution for all 10,000 runs compared to the observed percentages to assess whether these were in the upper 5% of simulated percentages. By using whole fly genotype probabilities to generate the null distribution, LD relationships among loci and, thus, the effects of inversions and genetic correlations between SNPs were accounted for by mirroring their associations in the randomized samples. Ragland et al. (2017) found that allele frequency differences between early-and late-eclosing apple vs. hawthorn flies from Fennville were significantly correlated genome wide in the GWAS (r = 0.54, p < 0.0001 for all 10,241 SNPs genotyped). We therefore used the mean of the frequency difference averaged between the apple and hawthorn races for the allele associated with later eclosion time for all analyses involving eclosion time, except where noted.

| Tests for genetic associations of eclosion time with geographic and host-related differentiation
We performed linear regression analyses (R Core & Team, 2017) to assess the degree to which allele frequency differences for SNPs between the earliest and latest quantiles of flies in the eclosion time GWAS were related to (a) allele frequency differences within the apple and the hawthorn host races between the Grant and Urbana sites (i.e., a measure of the extent of geographic variation for flies between the two most distant sites we sampled, essentially encompassing the latitudinal range of overlap of the host races in the Midwest; see Figure 1) and (b) allele frequency differences between the apple and hawthorn races at the four sympatric sites (i.e., local host-related F I G U R E 3 (a) Mean time to eclosion measured as the average number of days following postwinter warming for apple and hawthorn flies in controlled rearing experiments (data from Hood et al., 2015; see also Supporting Information Table S2 for additional site information) and b-d) mean standardized allele frequencies for SNP variants significantly associated with later eclosion time in the GWAS of Ragland et al. (2017) at the four latitudinally arrayed study sites surveyed across the Midwestern United States. Populations having higher mean standardized allele frequencies are therefore predicted to be genetically inclined to eclose later than populations having lower mean standardized allele frequencies. Panel (b) chromosome 2 and 3 SNPs together; (c) chromosome 1 SNPs; (d) chromosome 1-3 SNPs together; see Materials and Methods for details concerning the calculation of mean standardized allele frequencies for SNPs. n = number of significant SNPs on chromosome(s) contributing to mean standardized allele frequencies divergence). Significance levels were determined by Monte Carlo simulations in which correlation coefficients were calculated between two random samples of whole fly genotype probabilities taken with replacement from the appropriate pooled data sets of flies, as above.
Absolute values of observed correlation coefficients that exceeded the upper 5% of absolute r values in the simulated distributions of 10,000 replicates were taken as significantly different than expected.
Again, by using whole fly genotype probabilities to generate null distributions, the nonindependence of SNPs displaying high LD with one another was factored into the analyses testing for significance. To visualize how overall patterns of SNP allele frequencies significantly associated with eclosion time changed geographically and between the races, we calculated mean standardized allele frequencies for apple and hawthorn fly populations at all sites. For eclosion time-related SNPs (Ragland et al., 2017), we used the allele at a higher frequency in late-eclosing flies as the reference, set the average allele frequency for each race across all sites to zero, and calculated the mean deviation in allele frequency for each population. also highly correlated between the host races at Grant (r = 0.90, p < 0.0001, 4,243 df) and Urbana (r = 0.79, p < 0.0001, 4,243 df).

| Linkage disequilibrium
Likely due in large part to shared inversions and ongoing gene flow (Feder et al., 1994;Feder, Roethele, et al., 2003b), these similar patterns of LD imply that genome structure is concordant within and between the host races across the Midwest and allow for our broader application of the high, intermediate, and low LD SNPs classes categorized at Grant.

| Patterns of geographic variation within and between the host races
Both hawthorn and apple host races displayed substantial geo- Two major trends in geographic variation were apparent. First, geographic divergence varied among chromosomes, being greater for SNPs on chromosomes 1-3 than on chromosomes 4 and 5 (Table 1). In an exception to this pattern, chromosome 1 showed less geographic variation in the apple vs. hawthorn race (Table 1).  (Table 1). Together, these two trends support the first prediction of the standing variation hypothesis: High LD loci on chromosomes 1-3 displayed both the strongest associations with eclosion time (Ragland et al., 2017) and the highest levels of geographic divergence (Table 1).

| Eclosion time and geographic variation
Allele frequency differences between the early-and late-eclosing SNPs, the genome-wide correlation between the eclosion time GWAS and geographic differentiation in the hawthorn race was r = 0.60 (p < 0.0001; Figure 5a). The relationship was not as pronounced, but still significant, for the apple race (r = 0.35, p < 0.0001; Figure 5b).
For the hawthorn race, the genetic relationship between eclosion time and geographic SNP variation was significant for chromosomes 1-3 (r = 0.78, p < 0.0001, n = 2,620 SNPs) but not for the remainder of the genome (chromosomes 4 and 5, r = 0.09, p = 0.396, n = 1,624 SNPs; Table 2). The genetic correlation between eclosion time and geographic variation was also high among the eight different groups of high LD SNPs on chromosome 2 that presumably represent smaller inversions (r = 0.88, p = 0.0039, 7 df; Figure 6a).

| Host race, eclosion, and geographic differentiation across the genome
Genomic variation associated with eclosion time was also significantly related to host race divergence at the Grant, Fennville, Dowagiac, and Urbana sites. However, this relationship varied across chromosomes and was dependent upon the pattern of geographic F I G U R E 6 Relationships of the mean genetic responses of SNPs belonging to the eight different high LD groups on chromosome 2 in the eclosion time GWAS (i.e., their average allele frequency differences between early-and late-eclosing quantiles of flies) vs. (a) their mean frequency differences in the hawthorn race between Grant and Urbana; (b) their mean frequency differences in the apple race between Grant and Urbana; and (c) their mean frequency differences between the apple and hawthorn races at Fennville, MI. Numbers (1-8) in figures designate the eight different groups of high LD loci, presumably representing eight different relatively small-sized inversions residing on chromosome 2 variation displayed by each LD class (Table 3). In particular, host race divergence at sympatric sites was significantly associated with eclosion time for chromosomes 1-3 (Table 3). In contrast, no significant relationship was detected for chromosomes 4 or 5 at any site (Table 3). We therefore focus on chromosomes 1-3 below.
Patterns of host-related allele frequency divergence were similar between chromosomes 2 and 3, but differed for chromosome 1. When considered together, chromosomes 2 and 3 displayed a consistent pattern of divergence in mean eclosion-standardized allele frequencies between apple and hawthorn host races across sympatric sites and host-related divergence were significant and positive in sign at Urbana and Dowagiac, in line with prediction four, while they were significant and negative in sign at Grant and Fennville (Table 3). When chromosomes 1-3 were considered collectively, the effects of chromosome 1 outweighed those of chromosomes 2 and 3, such that standardized allele frequencies correctly predicted that apple flies should eclose earlier than hawthorn flies at Urbana and Dowagiac ( Figure 3d). Moreover, mean standardized allele frequencies predicted that the apple fly population at Dowagiac would deviate from the other sites and eclose relatively late (Figure 3d), as observed in nature ( Figure 3a). Nevertheless, mean standardized allele frequencies for chromosomes 1-3 still predicted that apple flies at Grant and Fennville would eclose later than hawthorn flies (Figure 3d), deviating from the observed phenotypic pattern (Figure 3a) and running counter to prediction four of the standing variation hypothesis.

| Standing variation and rapid adaptive evolution
Our results support the importance of standing variation in the recent formation of the apple-infesting host race of R. pomonella, a prime example of ecological divergence with gene flow (Berlocher & Feder, 2002). The patterns of eclosion-associated genomic divergence ob- Future studies are needed to identify the specific genes affecting eclosion time in R. pomonella, test for signatures of selection, and discern their demographic histories. In this regard, Ragland, Egan, Feder, Berlocher, and Hahn (2011) and Meyers et al. (2016) showed marked differences between host races in gene expression patterns prior to diapause termination and suggested genes in the Wnt and TOR signaling pathways as candidate loci. Whole-genome DNA sequencing is also needed examine the apple race for soft vs. hard selective sweeps around candidate loci (Hermisson & Pennings, 2005), the former predicted under the standing variation hypothesis, while the latter is expected if new mutations enabled the shift to apple.
In addition, previous analysis of cDNA loci implied that the putative inversions on chromosomes 1-3 containing eclosion loci are old (>1.5 Mya;Feder, Berlocher, et al., 2003;Feder et al., 2005). Thus, as discussed earlier, if secondary contact and gene flow with Mexican flies from the EVTM generated the standing variation in the hawthorn race, then we predict that genome sequencing studies will reveal that many alleles affecting eclosion will date to ~1.5 Mya.
However, this does not mean that variants of large effect (Colosimo et al., 2005;Hopkins & Rausher, 2011Lamichhaney et al., 2015;Steiner et al., 2007;Yuan, Sagawa, Young, Christensen, & Bradshaw, 2013) and of recent mutational origin (e.g., industrial melanism in the peppered moth; van't Hof, 2016) do not also play important roles in ecological adaptation and speciation. This could even be the case for Rhagoletis with respect to loci affecting other selected traits such as host plant choice (Dambroski et al., 2005).
Future work is needed to establish the genomic architecture of host odor preference in R. pomonella and infer the role of new mutations vs. standing variation in its contribution to host race divergence.

| Eclosion time and genomic divergence between the host races
Notably, our results differed from prediction four of the standing variation hypothesis: On balance, apple flies did not have higher frequencies of alleles associated with earlier eclosion than sympatric hawthorn flies. Although we found substantial standing geographic genetic variation for eclosion time on chromosomes 1-3 in the hawthorn race, which was significantly correlated with host-related divergence at all sympatric sites, the direction of the allele frequency differences between the races generally differed from expectation. Finally, gene-by-environment or nonadditive epistatic interactions may also play a role in determining eclosion time (Filchak, Roethele, & Feder, 2000). In the former case, differences in the preor overwintering period experienced by the host races at the more northern Grant and Fennville sites may result in alleles associated with later diapause termination in the laboratory having different effects in nature, causing apple flies to eclose earlier (Filchak et al., 2000). In insects, most latitudinal clines associated with diapause vary in response to photoperiod (Hut, Paolucci, Dor, Kyriacou, & Daan, 2013), including voltinism clines in the European corn borer that contribute to temporal isolation between pheromone strains (Levy, Kozak, & Dopman, 2018;Levy, Kozak, Wadsworth, Coates, & Dopman, 2015). Photoperiod also affects diapause in Rhagoletis (Prokopy, 1968). Thus, it is possible that longer photoperiods experienced by larval apple flies or by their mothers may cause them to eclose earlier than hawthorn flies in a manner not predicted by their genotypes. However, temperature and the length of the prewinter period following pupariation have been shown to exert a larger influence than photoperiod on eclosion phenotype in R. pomonella (Filchak, Roethele, & Feder, 2001). Moreover, Dambroski and Feder (2007) controlled for photoperiod and verified that there is a genetic basis for apple flies eclosing earlier than hawthorn flies, including at the Grant site. With respect to epistasis, it is possible that interactions between SNPs on chromosomes 1-3 might cause apple flies to eclose earlier than hawthorn flies. Elevated LD between these loci might be expected if strong epistatic interactions between chromosomes influenced the trait. However, considering the host races together, we did not find strong LD between eclosion-associated SNPs for other traits, such as host odor preference, are needed to fully characterize differentiation between host races at this site.

| Geographic variation not explained by eclosion time
Our results also raise questions concerning why genomic regions, other than those harboring eclosion variation (chromosomes 1-3), display significant geographic variation. For example, the high and intermediate LD classes of SNPs on chromosome 4 show significant geographic variation in the apple race (Table 1) but are not associated with eclosion time ( Table 2). The same is true in the hawthorn race for the low LD class of SNPs on chromosomes 4 and 5 (Tables 1 and   2). These results suggest that other factors, in addition to eclosion time, impose differential selection across the Midwest. In this regard, previous studies imply that the length of the prewinter period also differentially selects on the intensity at which flies initially enter pupal diapause. Feder, Roethele, Wlazlo, & Berlocher, 1997). The prewinter period is longer for apple than hawthorn flies and also varies with latitude (lengthening from north to south), resulting in both geographic and host-related differential selection on initial diapause intensity (Dambroski & Feder, 2007).
Moreover, initial diapause intensity may be under separate genetic control from eclosion time, at least in hawthorn flies (Ragland et al., 2017). Thus, selection on other aspects of the diapause phenotype might account for the significant geographic variation not explained by eclosion time.

| CON CLUS ION
In conclusion, we have shown that selection on eclosion time variation significantly sculpts geographic and host-related genomic differentiation in the apple and hawthorn races of R. pomonella.
Diapause life history adaptation in the recently formed apple race appears to have been extracted, in significant part, from standing genetic variation in the ancestral hawthorn race. Thus, our results highlight the potential inherent in geographic clines to contribute to the origins of new biodiversity when new resource opportunities become available. However, details concerning how this process has resulted in apple flies eclosing earlier than hawthorn flies, particularly at the more northern Grant and Fennville sites, remain to be resolved. Regardless, in Rhagoletis, polygenic standing variation, particularly in putative inverted regions, appears to be due to a history of geographic isolation and secondary contact Feder, Berlocher, et al., 2003). How general this is for other model systems of ecological speciation remains to be determined. Similarly, whether selection has widespread genomic effects on divergence in other organisms, as it does in Rhagoletis, or has simpler genetic underpinnings and consequences for facilitating speciation, is an open question. Nevertheless, it is clear that intraspecific clinal variation in life history timing is extensive in other organisms, including plants (Blackman, Michaels, & Rieseberg, 2011;Chen et al., 2012;Kawkakami et al., 2011;Keller, Levsen, Olson, & Tiffin, 2012;Lowry & Willis, 2010), vertebrates (Johnsen et al., 2007;O'Malley, Ford, & Hard, 2010), and insects (Adrion, Hahn, & Cooper, 2015;Bradford & Roff, 1995;Bradshaw, Quebodeaux, & Holzapfel, 2003;Kyriacou, Peixoto, Sandrelli, Costa, & Tauber, 2008;Lankinen, Tyukmaeva, & Hoikkala, 2013;Lehmann, Lyytinen, Piiroinen, & Lindström, 2015;Levy et al., 2015;Masaki, 1978;Paolucci, Zande, & Beukeboom, 2013;Roff, 1980;Wang et al., 2012). However, whether and how this intraspecific clinal variation is transformed into interspecific divergence are less well documented. It would seem that divergent selection on standing variation in life history timing could be a prime axis for contributing to the adaptation of organisms to new habitats and environments (review by Taylor & Friesen, 2017), potentially generating reproductive isolation and a wealth of new biodiversity.  -1560363).

CO N FLI C T S O F I NTE R E S T
All authors confirm that they have no conflicts of interest.

DATA ACCE SS I B I LIT Y
Data from the eclosion GWAS were previously deposited in DRYAD (https://doi.org/10.5061/dryad.kn568). All new sequence data generated for this study are deposited in DRYAD (https://doi.org/dryad. k42t7g2).