Estimating the impact of divergent mating phenology between residents and migrants on the potential for gene flow

Abstract Gene flow between populations can allow the spread of beneficial alleles and genetic diversity between populations, with importance to conservation, invasion biology, and agriculture. Levels of gene flow between populations vary not only with distance, but also with divergence in reproductive phenology. Since phenology is often locally adapted, arriving migrants may be reproductively out of synch with residents, which can depress realized gene flow. In flowering plants, the potential impact of phenological divergence on hybridization between populations can be predicted from overlap in flowering schedules—the daily count of flowers capable of pollen exchange—between a resident and migrant population. The accuracy of this prospective hybridization estimate, based on parental phenotypes, rests upon the assumptions of unbiased pollen transfer between resident and migrant active flowers. We tested the impact of phenological divergence on resident–migrant mating frequencies in experiments that mimicked a single large gene flow event. We first prospectively estimated mating frequencies two lines of Brassica rapaselected or early and late flowering. We then estimated realized mating frequencies retrospectively through progeny testing. The two estimates strongly agreed in a greenhouse experiment, where procedures ensured saturating, unbiased pollination. Under natural pollination in the field, the rate of resident–migrant mating, was lower than estimated by phenological divergence alone, although prospective and retrospective estimates were correlated. In both experiments, differences between residents and migrants in flowering schedule shape led to asymmetric hybridization. Results suggest that a prospective estimate of hybridization based on mating schedules can be a useful, although imperfect, tool for evaluating potential gene flow. They also illustrate the impact of mating phenology on the magnitude and symmetry of reproductive isolation.


| INTRODUC TI ON
Gene flow between populations has important consequences for the evolution, survival, and spread of species. It can introduce adaptive alleles into a declining resident populations, potentially saving it from extinction through genetic rescue (Frankham, 2015;Tallmon, Luikart, & Waples, 2004). One conservation measure-Assisted Gene Flow (AGF), has been proposed to affect genetic rescue of populations threatened by global change by introducing adaptive alleles from related population in situ (Aitken & Whitlock, 2013). On the other hand, a small locally adapted population with distinct genetic combinations can be "swamped out" by the inflow of genetic material from a larger population (Ellestrand, Prentice, & Hancock, 1999;Lenormand, 2002). Finally, gene flow into a species invading a new community can increase the species' invasiveness by introducing alleles that increase their survival, reproductive rate, or dispersal ability in their newly invaded environment (Ellestrand et al., 1999, Ellstrand & Schierenbeck 2000. Crucially, consequences of gene flow dependent not only on migration rates, but also on successful interbreeding between residents and migrants. Successful interbreeding, in turn, depends on synchrony of resident and migrant reproductive phenologies. Many species show geographic variation in reproductive phenology (e.g., Conklin, Battley, Potter, & Fox, 2010;Dambroski & Feder, 2007;Eckhart, Geber, & McGuire, 2004;Guttman, Bramble, & Sexton, 1991, Hall & Willis, 2006 These populations are thus not only isolated by distance, but also by time (Wright, 1943). In seasonally breeding species, even adjacent populations can be in reproductive isolation if their mating periods are asynchronous. For instance, flowering periods do not overlap between cultivated sunflowers and the weedy sunflowers growing in the vicinity, impeding gene flow between the two (Burke, Gardner, & Rieseberg, 2002).The apparent widespread nature of phenological divergence makes it likely to impact the dynamics of invasiveness, the efficacy of conservation practices, and the evolution of new species by reducing the potential of populations to interbreed.
Two plants can exchange pollen only if they are in flower at the same time. Most plants produce multiple flowers over a period of days, and so the mating probability between any two individuals can be estimated by the overlap in their flowering schedules. Wadgymar and Weis (2017) explored this by mimicking the first stage of an AGF program. Wadgymar and Weis grew northern "resident" Minnesota population of the prairie annual Chamecrista fasiculata alongside "migrants" from several potential source populations and estimated the potential for pollen exchange form their flowering schedules. The southern-most migrant population, sourced from North Carolina, came into flower ~24 days behind the residents when grown at the northern latitude. By the time the migrants reached peak flower production, the residents had already advanced to fruit maturation. Even though the residents and migrants were planted in equal numbers, fewer than 3% of the resident flowers could have received North Carolina migrant pollen, and reciprocally, less than 1% of the migrant flowers could have received resident pollen. Repeating the same procedure with a migrant population sourced from Pennsylvania, which flowered ~10 days after the residents, they found higher expected hybridization rates of 31% and 38% for residents and migrants, respectively, but these were still well short of the 50% expected under phenological matching.
The flowering schedule method devised by Weis and Wadgymar (2017) is simple to employ when designing AGF programs, but its predictive value rests on the assumption that on each day of the flowering season, every open flower has an equal probability of exchanging pollen with every other open flower. If pollinators transfer pollen between residents and migrants non-randomly, realized hybridization can differ from that predicted by phenological overlap alone. Gametic incompatibilities between residents and migrants could further cause realized hybridization rates to diverge from phenological expectations. The degree of bias can be evaluated by comparing the resident-migrant hybridization rate estimated prospectively form flowering schedules to a retrospective estimate based on progeny testing. Quite simply, with unbiased pollination (and complete gametic compatibility) resident plants should produce resident-migrant offspring it proportion to the number of contemporaneous mating opportunities they shared with migrants.
We performed a pair of two-generation experiments to mimic the translocation of phenologically mismatched migrants into a resident population. These experiments were designed to mimic an AGF program, whereby a single massive gene flow event is imposed upon a resident population in need of genetic rescue. In both experiments, the potential for hybridization was estimated in the first generation from overlap in flowering schedules of the phenologically divergent lines-the prospective estimate. The realized hybridization rate was revealed in the progeny generation through quantitative genetic methods. Specifically, the genotypes of the progeny generation were expected to include purebred residents, purebred migrants, and F 1 hybrids between the two. The proportion of successful hybrid matings by residents was determined retrospectively by comparing the flowering time distribution of their progeny to the distributions of known residents and hybrids, as detailed below. We found that realized hybridization rates were accurately predicted form phenological overlap in a greenhouse experiment, where procedures ensured random and saturating pollen loads to all flowers. Although prospective and retrospective estimates were correlated under natural pollination in the field, the rate of resident-migrant mating was lower than predicted by phenological divergence alone.

| Resident and migrant strains
We used two lines of Brassica rapa selected for divergent bolting time ) as a model system for phenological mismatch between residents and migrants in an AGF program. This species is a self-incompatible annual. The founding stock for the lines was collected from a large natural Quebec population (see . After two generations of selection the flowering times of the two lines differed by about two weeks under field conditions, with the earlier line flowering ~30 d after sowing and the later ~43 d. We will refer to the earlier flowering line as the "resident" population in need of evolutionary rescue through AGF, with the "migrant" population flowering later.

| Experiment 1: Hybridization under controlled pollination
The first experiment was designed to test if prospective estimates of hybridization, derived solely from flowering schedules, agreed with retrospective estimates, derived through progeny testing. It was performed in the greenhouse, and employed daily, random, hand pollen transfer. These are the conditions most likely to produce agreement between the two estimates. Failure to find agreement in this experiment would conclusively show that the Wadgymar and Weis (2017)

| Experiment 2: Hybridization under natural pollination
The second experiment, performed in the field, tested if flowering schedules adequately predict hybridization under natural pollination. Deviance of realized from potential hybrid matings could arise in several ways. For instance, pollinator preference for large inflorescences could cause plants at peak flower to contribute disproportionately to the successful pollen pool, unless there is a compensatory decline in the per-flower visitation rate with display size (e.g., Brys & Jacquemyn, 2010;Dudash, Hassler, Stevens, & Fenster, 2011;Karron & Mitchell, 2011). At low plant density, pollinators typically move between adjacent plants (Fenster, 1991;Thomson & Thomson, 1989), but at high plant density, skip over near neighbors (Cresswell, 1997). Finally, flowers produced late in the flowering schedule are less likely to set seed (Austen, Forrest, & Weis, 2015;Ison & Wagenius, 2014), so that pollen donated late in the season is less likely to sire offspring. Thus, on any one day, the pool of successful gametes may be a biased subset of the overall gamete pool, decaying the relationship between potential and realized hybridization.
This experiment also tested the notion that hybridization declines as phenological mismatch increases. We manipulated mismatch between populations by staggering the planting date for the In the summer of 2015, parental plants from the resident and migrant lines were exposed to natural pollinators at the Koffler Scientific Reserve at Joker's Hill (KSR), 58 km north of the University of Toronto. As of the summer of 2015, no wild population of Brassica rapa has been observed on the property. We applied three levels of mismatch. For the "Asynchronous Flowering" (mismatched) treatment, the two strains were planted on the same day (May 19th), so that the resident line would start to flower two weeks ahead of the migrant, as they naturally would. In the "Synchronous Flowering" treatment, migrants were planted 2 weeks ahead (May 6th) of the residents so that the two would come into flower on the same calendar date. For the "Intermediate" treatment, the migrant strain was planted only one week ahead of the resident. In addition to the three synchrony treatments, there were two control treatments, resident only and migrant only (expected hybridization = 0) planted on May 19th. These five treatments were replicated three times, for a total of 15 plots.
The experimental plants were started from seed in conetainer pots in the KSR greenhouse, using the same soil as above. The seeds were derived from the purebred controlled crosses performed as part of Experiment 1. Once all plants were established across all treatments, they were moved, in the containers, and sunk into 1.3 × 1.3 m experimental plots filled with local sand (see  on May 23rd. Each plot in each treatment held 21 plants per line, randomly dispersed in a hexagonal grid, with one plant from each purebred full-sibling family, so that experimental plots had a total of 42 plants and control plots had 21 plants. Thus, each plant in a plot had a full sibling in every other plot but no siblings in its own plot. Plots were separated by at least 200 meters to prevent pollen exchange between them, as in previous pollination experiments (Kunin, 1993;Mustajärvi, Siikamäki, Rytkönen, & Lammi, 2001). families from the controlled crosses in experiment 1 were planted: 12 resident and 12 migrant purebred sibships, and F 1 hybrid sibships from 12 resident and 12 migrant mothers, with 12 plants per sibship, for a total of 576 plants of known parentage. We used the same conetainer pots and soil as for the field generation.

| Prospective estimate of hybridization
Our notation for the several hybridization estimates are as follows: h iR is the proportion of flowers on individual i of the resident population expected to be pollinated by migrants; H R is the mean of h iR , and constitutes the prospective estimate of the frequency of hybrids among the seed progeny produced by the resident population; H' R is the retrospective estimate for the frequency of hybrids among the seed progeny of residents, derived from progeny testing. The corresponding hybridization rates for migrant seed parents use the subscript M.
We calculated a prospective estimate of hybridization rate for each plant (Wadgymar & Weis 2017;Weis, 2015)

| Retrospective estimate of hybridization
Our strategy to retrospectively estimate hybridization rates, H' R and H' M , used flowering time itself as the marker for paternity. The two and 97.5th percentile values to obtain confidence intervals. Known hybrids and purebreds were run through the LDA as a "training set" to determine how accurately it assigned genotype. The LDA was run using the MASS package in R (Venables & Ripley, 2002).
To determine the relationship of fruiting success to plant age (days

| Hybridization under controlled pollination
In the parental generation of the controlled pollination experi-   We note that when the minimum χ 2 goodness-of-fit method was applied to the control plots, neither H' R nor H' M were different from the true value of zero (Table 2).

| Hybridization under field conditions
Overall, the prospective estimate of hybridization was a reasonable, but biased predictor of the retrospective estimate ( Figure 10). The Spearman's rank correlation between the two was 0.84 and 0.68 (one-tailed p < 0.0005 and p < 0.01) for residents and F I G U R E 9 χ 2 goodness-of-fit distributions for migrants from each experimental plot with the hybridization value at the minimum χ 2 , confidence intervals, and prospective estimate plotted. Plot number is in the top left migrants respectively, when the control plots are included. If only experimental plots are considered, the rank correlations are 0.64 and 0.27 (one-tailed p < 0.05 and p < 0.24) for resident and migrant plants, respectively. However, H R > H' R in 8 of the 9 experimental plots, and H M > H' M in 6 of the 9.
We evaluated age-dependent fruit set as a source of bias.
Although the flowers that plants produce late in life are less likely to set fruit, the effect was trivial.
The Linear Discriminant Analysis gave the same qualitative results for retrospective estimates of H' R and H' M as the minimum χ 2 goodness-of-fit method (Table 3) and correlated with the prospective estimates. These appear biased in the complementary way, overestimating hybridization for migrant plants, but not for residents. However, the LDA method estimated hybridization rates in the control plots to be ~0.2, far greater than the true value of zero.
Therefore, this second estimation method appears to be fundamentally biased.  In the field experiment, where plants were exposed to natural pollinators only, realized hybridization rates were correlated to those predicted from flowering schedules, but were lower.

| D ISCUSS I ON
Beyond divergence in the onset of flowerings, schedules for residents and migrants also differed in size and shape. These differences caused asymmetry in hybridization. In the greenhouse, the early-blooming residents produced more flowers across their schedule than the later migrants, which increased opportunities for the former to transfer pollen to the latter. In the field, the later-blooming residents flowered over more days than residents, which increased opportunities for migrant-migrant pollen exchange.
Together, these results suggest that while the opportunity for mating can be a useful indicator of potential gene flow, it may underestimate actual rates of pollen exchange in some environments.
Below we explore reasons for quantitative differences between prospective and retrospective hybridization estimates, and more broadly, discuss implications of phenological mismatch.

| Lost mating opportunities
The quantitative discrepancy between prospective and retrospective estimates for residents (H R vs. H' R ) may reflect the loss of mating opportunities under natural pollination. As noted above, non-random pollen exchange could arise though pollinator preferences for large inflorescences, decline in the per-flower visitation rate with display size, or foraging responses to flower density. These could cause the proportion of migrant pollen delivered to receptive stigmas to differ from its proportion in the pollen pool on a given day (q d , Equation (1) and (2). Day to day changes in weather could have caused pollination intensity to fluctuate, both by directly effecting pollinator foraging, and indirectly by altering floral rewards (Pleasants, 1983;Vicens & Bosch, 2000). Inclement weather during the narrow window of overlap in resident and migrant flowering schedules could thus reduce hybridization. Conversely, unfavorable conditions before and after the overlap window could decrease purebred mating frequency.
Herbivory in the field experiment was very low, but in other cases it could affect mating opportunities, especially if the strength of herbivory varies over time. Many studies have found that herbivory often peaks around peak-flowering (Elzinga et al., 2007). This could reduce hybridization opportunities when plants flower in synchrony, or increase it when they flower asynchronously, by reducing seed-set in a biased way. Disease could have similar on-peak effects, with transmission rates of floral diseases being higher when more flowers are open, or diseases in general which are carried by pollinators (Elzinga et al., 2007). The anther smut Microbotryum violaceum, attacks more often earlier in the season, countering the reproductive advantage early-flowerers have in male function (Biere & Antonovics, 1996). Finally, pollinators may bypass inflorescences with several damaged flowers, leaving the neighboring intact flowers unmated (Krupnick, Weis, & Campbell, 1999).
We note that hybridization rates were asymmetrical in the field experiment; that is, migrant plants produced about as many hybrid offspring as predicted by flowering schedules, while resident plants produced fewer. This asymmetry may have emerged from the selection regime that created the two phenological lines. The linear discriminant analysis showed that plants that flower early have larger flowers (see Appendix S1B). In Brassica rapa flower size has been correlated with pollen size (Sarkissian & Harder, 2001), which in turn is correlated with increased siring ability (McCallum & Chang, 2016;Sarkissian & Harder, 2001) and with increased competitiveness of sired seeds within the ovary (Cruzan, 1990). Earlier experiments with the base population for the selection lines also suggested a siring advantage for early flowering (Austen & Weis, 2016). Asymmetric hybridization could alter the evolutionary trajectory if cytoplasmic genes, which are maternally inherited, are under local selection (e.g., Galloway 2005). Additionally, migrants were expected to produce more hybrid offspring in the greenhouse, whereas the opposite was true in the field, owing to them producing more late-season flowers.
This asymmetry in hybrid siring rate could vary across years and locales due to changes in growing season length.
Reduced initial hybridization, as we detected, would not necessarily prevent large-scale introgression of migrant genes. If the purebred progeny of migrants have sufficient absolute fitness (net reproductive rate), adaptive loci can continue to introgress into the resident genetic background through both migrant × resident and migrant × hybrid matings at a rate dependent on their phenological mismatch. However, local selection may act against migrants, and introgression may be slowed considerably if it acts against migrant phenology itself. The residents that mate with migrants will be those with the most migrant-like phenology (Weis, 2015). If selection acts against migrant phenology, it will be acting against migrant-like residents as well, meaning that hybrids will tend to inherit maladaptive phenology alleles from both parents. This "narrows the bridge" over which adaptive migrant alleles pass into the resident population. Any adaptive alleles the migrant may have will likely be in linkage disequilibrium with their maladaptive phenology loci in early generations, further slowing introgression (Soularue & Kremer, 2014).

| Implications
The field experiment showed initial hybridization rates can decline with increasing phenological divergence between residents and migrants. Being able to predict potential gene flow between populations based on phenology has implications particularly for agriculture and conservation, but also for invasion biology. Understanding what allows and inhibits hybridization can do much to elucidate how new species arise and diverge. Introgression can also break-down barriers between species, resulting in the birth of a new species from the collapse of previous ones. Many such species complexes are breaking down on observable time-scales. The three-spined stickleback in some lakes has collapsed from two biological species into one due to hybridization (Taylor et al., 2006). On the other hand, weedy sunflowers in highly-infested fields have evolved to be more divergent phenologically from their croppy cousins (Roumet, Noilhan, Latreille, David, & Muller, 2013). Flowering phenology can pose an effective barrier to hybridization, allowing early-and late-flowering varieties to exist in sympatry with minimal gene flow. (Soliva & Widmer 1999).
These divergences in reproductive timing could eventually lead to the development of new species. Crop-weed gene flow is an important problem in modern agriculture. The transfer of genes into weedy populations is well documented (Roumet et al. 2013, Chèvre et al. 2000, Langevin, Clay, & Grace, 1990, Ellstrand et al., 1999. In particular, gene flow can lead to the spread of herbicide resistance to weedy populations, reducing the effectiveness of herbicides in weed control (Kreiner, Stinchcombe, & Wright, 2017). Multiple herbicide resistance alleles have spread into escaped populations of canola, even outside of crop fields (Knispel, McLachlan, Van Acker, & Friesen, 2008). Of the world's 13 most important food crops, 12 have been found to hybridize with wild populations somewhere in their range (Ellestrand et al., 1999). If selection for crop traits is advantageous, it can lead to genetic swamping in wild populations (Ellestrand et al., 1999;Haygood, Ives, & Andow, 2003). Transgene escape in sunflowers has been shown to lead to increased fecundity and decreased herbivory in wild populations (Snow et al. 2003).
Demographic swamping may also occur in instance where large populations of crops successfully breed with wild populations, flooding the ecosystem with unfit hybrids (Haygood et al., 2003).
Crops planted in large monocultures are likely to spread their genes to receptive wild populations. Choosing cultivars that flower asynchronously from their wild relatives (Jenczewski, Ronfort, & Chèvre, 2003) can be employed to slow crop-weed gene flow through phenological divergence (Roumet et al. 2012).
Gene flow is an important factor in the survival of colonizing species. Gene flow to a founder population will prevent founder effects and improve genetic diversity. Understanding how species colonizing novel ecosystems become invasive is critical to preserving global biodiversity (Vitousek, D'antonio, Loope, Rejmanek, & Westbrooks, 1997). Many invasive species do not become invasive immediately after colonization, suggesting that many species are not simply preadapted to invasion (Ewel et al. 1999). Hybridization between native species and an introduced species may cause a native species to become invasive through the input of both additional heterozygosity and novel adaptive variation (Ellstrand & Schierenbeck 2000).
Multiple introductions of a species may augment the amount of genetic diversity an invading species has, provided these introductions readily interbreed (Dlugosch & Parker, 2008).
Rising global temperatures may send some narrowly-adapted populations into decline. A new, unfavorable climate can reduce mean fitness of a local population, reducing net reproductive rate below replacement level, eventually leading to local extinction. If rapid enough, adaptation can restore mean fitness. This process, called "evolutionary rescue", depends critically on standing genetic variation in the traits under selection (Bell & Gonzalez, 2009;Gomulkiewicz & Shaw, 2013). Evolutionary rescue is at the core of Assisted Gene Flow (AGF), which depends on interbreeding between populations. Hybridization between species and populations can lead to increased diversity and introduce novel traits, and potentially promote evolutionary rescue (Aitken & Whitlock, 2013, Janes & Hamilton, 2017. Gene flow can be inhibited by factors other than phenological divergence, however many of these factors are unlikely in most scenarios where AGF is relevant (Aitken & Whitlock, 2013). Phenological divergence poses a much more significant problem and is rooted in the biology of the conservation targets.
Migrant phenology in its habitat may not reliably indicate phenology after translocation, however. Flowering in many species of plants occurs after accumulating a number of degree-days above a basal temperature (Forrest & Thomson, 2011). Other plant species have a vernalization requirement, that is, the heat accumulation mechanism is activated in spring after accumulating a particular number of winter chilling degree-days (e.g., Mimura & Aitken, 2010;Stinchcombe et al., 2004). Fruiting trees, such as almonds, have different chilling requirements, causing them to be very phenologically divergent (Egea, Ortega, & Martıńez-Gómez, P., & Dicenta, F., 2003). In many species, heat accumulation mech- Our experiment demonstrates that reproductive timing has a large, yet predictable, effect on hybridization. While field results are likely to differ from estimates of mating opportunities, our prospective estimate will still give insight into potential for populations to hybridize.

ACK N OWLED G M ENTS
We thank the following individuals for assistance with plant care and data collection in field and greenhouse: ChuChu Lu, Susan Gordon, Sara Bowman, Madeline Peters, Julia Holden, Khadija Jafarova, and Hyerin Jeong. Stephan Schnieder and John Jensen assisted in field setup. The manuscript benefited from the helpful advice from Benjamin Gilbert, Spencer Barrett, James Thomson, and Stephen Wright.

AUTH O R CO NTR I B UTI O N S
Colin Bonner and Arthur E. Weis conceived of the presented idea and wrote the manuscript. Sally Erin Westover bred the phenological lines of Brassica rapa. Michelle Ho and Nina A. Sokolov contributed to data collection and analysis and assisted in editing the manuscript.