Genetic diversity and structure of Lolium perenne ssp. multiflorum in California vineyards and orchards indicate potential for spread of herbicide resistance via gene flow

Abstract Management of agroecosystems with herbicides imposes strong selection pressures on weedy plants leading to the evolution of resistance against those herbicides. Resistance to glyphosate in populations of Lolium perenne L. ssp. multiflorum is increasingly common in California, USA, causing economic losses and the loss of effective management tools. To gain insights into the recent evolution of glyphosate resistance in L. perenne in perennial cropping systems of northwest California and to inform management, we investigated the frequency of glyphosate resistance and the genetic diversity and structure of 14 populations. The sampled populations contained frequencies of resistant plants ranging from 10% to 89%. Analyses of neutral genetic variation using microsatellite markers indicated very high genetic diversity within all populations regardless of resistance frequency. Genetic variation was distributed predominantly among individuals within populations rather than among populations or sampled counties, as would be expected for a wide‐ranging outcrossing weed species. Bayesian clustering analysis provided evidence of population structuring with extensive admixture between two genetic clusters or gene pools. High genetic diversity and admixture, and low differentiation between populations, strongly suggest the potential for spread of resistance through gene flow and the need for management that limits seed and pollen dispersal in L. perenne.


| INTRODUCTION
Weedy plants pose a major problem to agricultural production causing significant crop losses worldwide and economic damages estimated to total $33 billion annually in the United States (Oerke, 2006;Pimentel, Zuniga, & Morrison, 2005). Weeds are an ongoing challenge for farmers as weed control practices exert strong selection for the evolution of weed adaptations that render the management practices less effective over time (Barrett, 1983;Owen, Michael, Renton, Steadman, & Powles, 2011;Powles & Yu, 2010). One of the best examples of this process is the evolution of resistance to herbicides. In weed populations containing phenotypic variation for susceptibility to an herbicide, those individuals with an inherited ability to survive and reproduce following an herbicide application are favored and resistance increases in the population over time (Delye, Jasieniuk, & Le Corre, 2013;Neve, Vila-Aiub, & Roux, 2009). To date, over 470 cases of resistance in 250 species have been documented to a wide variety of herbicides worldwide (Heap, 2016).
Whether or not a weed population is able to adapt in response to management practices depends on whether that population contains the necessary genetic variation (Jasieniuk, Brule-Babel, & Morrison, 1996;Sakai et al., 2001). Population size, standing genetic variation, selection, and gene flow with other populations all play a role in the spatial distribution of evolved adaptive traits (Delye, Jasieniuk, et al., 2013;Lawton-Rauh, 2008). For studies of adaptation in agricultural weeds, strong selection pressures on weed populations such as tillage or herbicide application are usually known and population sizes are often large . Population sizes of common weeds vary across an agricultural landscape with some areas containing heavy infestations, allowing for high genetic diversity within a species across a region through the accumulation of mutations over time. In self-pollinating weeds, populations may be genetically uniform as individuals within populations often share nearly identical highly homozygous genotypes because of repeated inbreeding, but populations are likely to differ genetically (Ward & Jasieniuk, 2009). In contrast, obligately outcrossing weeds are expected to contain high genetic diversity within populations but low genetic differentiation among populations. The amount and distribution of phenotypic and genetic variation within weed populations influence the potential for adaptation in agricultural landscapes, which are variable in both space and time as a result of habitat fragmentation due to diverse crops and associated crop and weed management practices. Ultimately, the adaptation of weed populations to a variable environment across an agricultural landscape may lead to population structuring in both selfing and outcrossing weeds.
Genetic diversity in weed populations is required for weed adaptation, but is also impacted by it as a result of strong positive selection and population bottlenecks . Successful herbicide applications kill 95%-99% of individuals in susceptible weed populations. This substantial reduction in population size may mean that the alleles of only a small fraction of individuals are passed on to the next generation, potentially causing some alleles to be lost by genetic drift.
Alternatively, strong selection will favor selectively advantageous alleles, if present in the population, and reduce population genetic diversity. For instance, individuals which survive herbicide treatment due to heritable mechanisms will pass on their resistance-conferring alleles to their progeny, and resistance will increase in frequency in the population over time. As the frequency of resistant individuals increases in a population, further herbicide applications will become less effective in reducing population size, leading to restoration of populations to their original size but with decreased genetic diversity.
Strong selection for resistance may also be associated with a selective sweep at causative loci which not only results in the loss of susceptible alleles at the adaptive locus but also any alleles at loci in gametic disequilibrium with it (Maynard-Smith & Haig, 1974;Menchari, Delye, & Le Corre, 2007). In summary, weed populations with a high frequency of resistant individuals are expected to contain lower genetic diversity than populations with a low frequency of resistant individuals both due to population bottlenecks while an herbicide is still effective in controlling the weed and due to selection as resistance to the herbicide evolves.
It has been hypothesized that gene flow may spread herbicide resistance among weed populations within an agricultural landscape to a greater degree than novel mutations as rates of gene flow are generally believed to be higher than rates of mutation (Jasieniuk et al., 1996).
Herbicide resistance alleles may be present in populations prior to the onset of selection pressure by an herbicide (Delye, Deulvot, et al., 2013), and may spread by gene flow even before the trait is selectively advantageous. Evidence for the spread of herbicide resistance among populations by seed dispersal has been shown in several highly selfpollinating weed species, based on patterns of molecular marker and phenotypic variation (Okada et al., 2013(Okada et al., , 2014Osuna, Okada, Ahmad, Fischer, & Jasieniuk, 2011). Interestingly, however, neutral genetic and phenotypic variation in Ipomoea purpurea, a weed species with a mixed mating system, provided support for independent origins of resistance in multiple geographic locations (Kuester, Chang, & Baucom, 2015). In outcrossing weeds, analyses of neutral genetic variation revealed low population differentiation, and possible spread of resistance through local gene flow (Delye, Clement, Pernin, Chauvel, & Le Corre, 2010) but independent origins through novel mutations (Menchari et al., 2007).
The goal of the current study was to characterize genetic variation of northwestern California L. perenne populations where herbicide resistance evolution is very recent and likely ongoing. We examined the frequency of glyphosate-resistant plants in populations across the landscape along with microsatellite marker variation to address the following questions: (i) do populations of outcrossing weeds contain high genetic diversity and is this diversity reduced in populations with a high frequency of glyphosate-resistant individuals, (ii) is there evidence of genetic structuring and differentiation among populations of this widespread weed across an agricultural landscape, and (iii) is there potential for spread of resistance alleles across the landscape through gene flow?

| Population sampling
To determine whether glyphosate-resistant individuals are present in L. perenne ssp. multiflorum populations in northwest California, we sampled 13 orchards and vineyards in 2013 from Sonoma County and Lake County (Table 1) in the general regions where growers had reported difficulty controlling plants with glyphosate to farm advisors and in surrounding areas where populations may be experiencing gene flow with resistant plants. One population identified as resistant to glyphosate from Butte County was also sampled to serve as a comparison with an area which had evolved resistance greater than 10 years ago (Jasieniuk et al., 2008;Simarmata et al., 2003). Within each population, young leaf tissue and panicles with mature seed were collected from each of 30-40 individuals at least one meter apart from one another while walking randomly selected tree or vine rows. Leaf tissue was transported to the laboratory for DNA extraction. Seed panicles were stored in paper envelopes for 3 months to allow seeds to after-ripen and overcome dormancy before planting and testing plants for resistance to glyphosate.

| Phenotyping plant response to glyphosate
Eight seeds from each sampled plant were germinated on moistened filter paper in petri dishes at 20°C and a 12-hr photoperiod.
Germinated seedlings were transplanted into 8 × 8 cm square pots filled with UC soil mix (sand, compost, and peat in 1:1:1 ratio with 1.8 kg/m 3 dolomite) with two seedlings per pot and grown in the glasshouse at 27/15°C with ambient light conditions. At the tillering stage, individual plants were divided into genetically identical clones following the method described by Boutsalis (2001) and grown in the glasshouse to the 2-3 leaf stage. One clone of a genotype was treated with water, which served as a control. The second clone was treated with glyphosate (Roundup PowerMax, Monsanto, St. Louis, MO) at the rate of acid equivalent 1,681 g/ha, which is twice the recommended  (Jasieniuk et al., 2008) were included during each herbicide application to confirm herbicide activity.
Two sets of alleles, with one or two alleles in each of the two re-  (Table 2). PCR products were multiplexed into six pairs of PCR product and separated using an ABI 3100 Genetic Analyzer (Applied Biosystems, Foster City, CA) with GENESCAN 400HD as an internal size standard (Applied Biosystems).
Genotypes were also inspected manually.

| Microsatellite marker error rate
To assess genotyping errors, we tested for null alleles or genotyp-

| Genetic diversity and structure
To estimate the allelic diversity of each locus, we calculated the total number of alleles detected at each locus (N A ) using F STAT software T A B L E 2 Characteristics and sources of 12 microsatellite loci used to genotype Lolium perenne ssp. multiflorum. Microsatellite markers were selected based on polymorphism and consistent amplification of alleles version 2.9.3.2 (Goudet, 1995 (Cornuet & Luikart, 1997;Piry, Luikart, & Cornuet, 1999). A sequential Bonferroni correction was applied to adjust significance levels for multiple comparisons (Rice, 1989).
To investigate the spatial structuring of genetic variation among populations and counties, we performed multiple distance-based analyses using GenAlEx 6.5 (Peakall & Smouse, 2006. First, we calculated a matrix of genetic distances, based on Nei's D (Nei, 1978), between each pair of populations with GenAlEx software interpolating for missing data, and then used the pairwise genetic distances in the following analyses. We performed a principal coordinate analysis (PCA) using the standardized covariance method in GenAlEx to assess whether populations located within the same county and/ or geographically near each other were also genetically similar and thus may share an evolutionary history. To test for isolation by distance (IBD) and determine whether genetic structuring correlates with geographic structuring between populations, we conducted a Mantel test (Mantel, 1967) in GenAlEx with 1,000 permutations between the matrix of pairwise genetic distances described above and a matrix of pairwise geographic distances calculated using the GPS coordinates (latitude and longitude) of each population. We also performed a hierarchical analysis of molecular variance (AMOVA; Excoffier, Smouse, & Quattro, 1992) in GenAlEx, which examined the distribution of genetic variation at five hierarchical levels: among counties within the total sample, among populations within those counties, among populations in the total sample, among individuals within populations, and among individuals within the total sample. The AMOVA was performed with 1,000 permutations.
To assess population structure and determine the degree of admixture among populations, we used a model-based Bayesian clustering algorithm in STRUCTURE software version 2.3.4 (Pritchard, Stephens, & Donnelly, 2000). STRUCTURE infers genetic clusters or populations based on the multilocus genotypes of all individuals, independent of sampled location, by probabilistically assigning individuals to a cluster or jointly to multiple clusters if their genotypes indicate admixture between clusters, while simultaneously maximizing Hardy-Weinberg equilibrium and minimizing linkage disequilibrium within those clusters. STRUCTURE analysis was performed with a burn-in period of 1,000,000 iterations followed by 1,000,000 iterations for the number of genetic clusters (K) ranging from K = 1 to K = 12. Five independent runs at each value of K were performed using the population admixture model for potentially interbreeding populations and correlated allele frequencies. Likelihood values of ln P(D) were assessed for each run. The most likely value of K was inferred using the ΔK method (Evanno, Regnaut, & Goudet, 2005) in STRUCTURE Harvester online (Earl & von Holdt, 2012). Each individual's probability of assignment to each cluster (q), also interpreted as the proportion of an individual's genome that originated in each cluster (Pritchard et al., 2000), was visualized for all individuals using Distruct software version 1.1 (Rosenberg, 2004). To examine substructuring within genetic clusters, the multilocus genotypes of individuals with q > 0.6 to a cluster were analyzed independently, as suggested by Evanno et al. (2005), using the same parameters as above.

| Plant response to glyphosate
Within sampled populations of L. perenne ssp. multiflorum, resistance to glyphosate, estimated as the percentage of individuals surviving glyphosate treatment per population, varied from 9.7% to 89.0% ( Figure 1, Table 1

| Genetic diversity of microsatellite loci
Information on null alleles for each locus can be found in the Data S1 section. We detected 259 distinct alleles in 412 individuals of L.

| Genetic diversity of populations
Within populations, the average number of alleles detected per locus (N A ) was high ranging from 8.7 to 12.1. Correspondingly, allelic richness (A R ) ranged from 7.7 to 9.2 with an average of 8.4 over all populations (  (Table 4). However, the population (population 13) containing the highest frequency of glyphosate-resistant individuals (89% R) (Table 1)    As population clustering by county was not explained by geographic distance, an AMOVA was conducted to determine how much of the genetic variance could be attributed to county or population differences. The AMOVA revealed low, but significant, genetic differ- STRUCTURE analysis (Pritchard et al., 2000) was used to further examine genetic structuring. STRUCTURE revealed increasing values of ln P(D) with increasing K values ranging from 1 to 12 with no clear maximum likelihood (Figure 3a). ΔK (Evanno et al., 2005) clearly showed the highest value at K = 2, indicating two genetic clusters ( Figure 3b). The proportion of the genome, as represented by the 12 microsatellite loci, that assigns to each cluster, q, was calculated for each individual. Individuals assigning to cluster 1 with q > 0.7 were comprised of some individuals from Sonoma County and Butte County and a few individuals from Lake County, while most individuals from Lake County and some from Sonoma County and Butte County assigned to cluster 2 (Figure 4a). While individuals from Sonoma County and Butte County assigned to both genetic clusters, individuals from Lake County assigned highly to cluster 2. All populations contained some individuals that assigned partially to each cluster (q < 0.6) indicating admixture between genetic clusters (Figure 4a).

| Population structure
There is no apparent pattern of genetic structuring based on whether the individuals originated from field populations that were predominantly resistant or susceptible to glyphosate (Table 1). The majority of individuals from populations 12, 13, and 14, where the frequency of glyphosate resistance was 85%-89%, assigned to cluster 2 with q > 0.7, but most individuals from population 15 where resistance frequency was 21% also assigned to the same cluster. Of the individuals from population 1 genotyped, 48% and 52% assigned to clusters 1 and 2, respectively, but 77% of individuals phenotyped were resistant to glyphosate, indicating that glyphosate-resistant individuals likely assign to both clusters.
To examine patterns of hierarchical population structure, individuals assigning to each cluster with q > 0.6 were separated and analyzed independently. STRUCTURE analysis of each cluster revealed evidence of subclustering, with K = 3 within cluster 1 and K = 4 within cluster 2 ( Figure 3c-f). Among individuals that assigned to cluster 1, most individuals assigned highly to one of the three subclusters ( Figure 4b).
However, there is little apparent geographic substructuring of genetic variation, with the exception that most individuals from Butte County assigned to subcluster 1. Among individuals that assigned to cluster 2, most individuals were admixed, assigning to multiple subclusters ( Figure 4c).

| High genetic variation observed in weedy L. perenne regardless of resistance frequency
Our analyses indicate very high genetic diversity within populations of L. perenne ssp. multiflorum as expected based on the biology of this widespread obligately outcrossing species. However, given the large number of detected alleles and the allele frequencies observed, a higher level of heterozygosity was expected than was observed across all populations (   (Luikart & Cornuet, 1998), although the bottleneck imposed by intense selection by glyphosate and the evolution of glyphosate resistance likely has been recent, probably occurring within the last 20 years, that is 20 generations, based on when resistance to glyphosate first was identified in Lolium in the United States and worldwide (Powles, Lorraine-Colwill, Dellow, & Preston, 1998;Simarmata et al., 2003). Although resistance alleles may have been present at low frequencies previous to 20 years ago, they did not rise to frequencies high enough to be problematic until glyphosate  (Delye et al., 2010;Menchari et al., 2007;respectively), and populations of a glyphosate-resistant weed with a mixed mating system, Ipomoea purpurea, also had lower total genetic diversity (H E = 0.304) based on microsatellite markers (Kuester et al., 2015). Populations of species with an outcrossing mating system would be expected to have higher genetic diversity than those with predominantly self-pollinating mating systems. Glyphosateresistant populations of the two closely related selfing species Conyza canadensis and Conyza bonariensis in California displayed a wide range of genetic diversities (H E = 0.0-0.45 and H E = 0.009-0.513, respectively) based on microsatellite markers (Okada et al., 2013(Okada et al., , 2014. in California populations was also similar or higher than that seen in studies of L. perenne cultivars using some of the same microsatellite markers (N A = 19.4, 13.3, and 9.9; Kubik et al., 2001;Wang et al., 2009;Brazauskas et al., 2011, respectively). Interestingly, the forage cultivars of L. perenne do not show lower genetic diversity than their wild and weedy relatives despite many generations of breeding.
Rather, L. perenne seems to display a very high level of genetic diversity, regardless of origin.

| Spatial and genetic structuring of populations
Glyphosate-resistant plants were detected in L. perenne populations from all sampled areas (Table 2) (Table S3) both reveal larger genetic distance between populations within Lake County than between populations within Sonoma County.
This indicates that despite the relatively larger geographic distances, Sonoma County populations are closely related to each other due to either higher genetic exchange between them or less differentiation over time, possibly related to the relative homogeneity of the vineyard cropping system across Sonoma County compared to the mix of perennial crops grown in Lake County. Differences in water availability in primarily drip-irrigated vineyards compared to sprinkler-irrigated orchards may lead to local adaptation or phenological differences due to water stress may contribute to population differentiation in Lake County populations. The lack of correlation between genetic distance and geographic distance found by a Mantel test indicates that genetic distance between populations is not due mainly to factors associated with isolation by distance, but to some other factor such as local adaptation or differences in the strength of selection pressures across the landscape. Extensive long-distance gene flow may also erode the relationship between genetic distance and geographic distance, especially if gene flow is not homogeneous across the entire range. An analysis of molecular variance indicates a significant amount of genetic variation is distributed at the county and population level (Table 5).
However, both county and population differences are outweighed by the high genetic variation among individuals, as might be expected in this highly diverse outcrossing species (Brazauskas et al., 2011).
Bayesian clustering STRUCTURE analysis indicates the presence of two distinct gene pools or genetic clusters (Figure 3a). While the individuals assigning to cluster 1 are primarily from Butte County and Sonoma County, individuals assigning to cluster 2 come from all three counties and include almost all individuals from Lake County ( Figure 4a). This indicates that there is little admixture between individuals in Lake County and individuals in the other sampled areas or that L. perenne has been introduced into the region too recently for substantial admixture to have occurred. In the subclustering analysis, there is little apparent spatial structure among individuals assigning to subclusters (Figure 4b,c). Many individuals had admixed genomes assigning partially to multiple subclusters, indicating little differentiation or high gene flow between individuals assigning to these subclusters.

| Evolution of resistance and potential for spread of resistance alleles
Together, data on spatial patterns of population structuring and frequencies of glyphosate-resistant phenotypes allow comparison of hypotheses regarding single or multiple evolutionary origins and subsequent spread of glyphosate resistance in L. perenne populations in northwest California. Populations with moderate frequencies of resistant individuals in Sonoma County (30% > R > 80%, populations 7, 9, 10, 11) contain mostly individuals that assign to cluster 1 with q > 0.7, while populations with high frequencies of resistant individuals (R > 80%, populations 12, 13, 14) contain a large proportion of individuals that assign to cluster 2 with q > 0.7 (Figure 4a, Table 2).
This suggests that glyphosate resistance has likely evolved independently in individuals that assign to each cluster. There is also evidence that some plants are resistant due to a mechanism other than target-site mutations, indicating that multiple mechanisms of resistance are present in the region resulting from at least one additional independent origin of resistance (Karn and Jasieniuk, in review). It is possible that in the future, additional novel mutations may result in yet more independent origins of resistance.
STRUCTURE analysis revealed potential for future spread of resistance alleles through gene flow. Many admixed individuals with genotypes that assigned partially to each cluster were identified (Figure 4).
While the majority of individuals assign highly to a single cluster, the considerable number of admixed individuals in many populations indicates that gene flow is common. Localized gene flow between populations located near each other may be pollen-mediated, and pollen movement over distances of 3 km has been documented in L. perenne ssp. rigidum (Busi, Yu, Barrett-Lennard, & Powles, 2008). Short-and long-distance gene flow may also be mediated by seed movement on agricultural machinery and vehicles, or over short distances by wind or animals. The higher levels of admixture detected in this study compared to studies of other herbicide-resistant weeds (Kuester et al., 2015;Okada et al., 2013Okada et al., , 2014 likely relate to the outcrossing nature of L. perenne. Spread of resistance alleles through gene flow could also result in populations and individuals containing multiple mechanisms of glyphosate resistance. Successful management of glyphosate-resistant L. perenne populations in perennial cropping systems will likely require implementation of integrated pest management programs that include chemical and nonchemical techniques to not only control currently resistant plants, but reduce the intensity of selection pressure for future independent origins of resistance and limit the spread of resistance through gene flow. Increasing tillage where possible, mowing to reduce seed set, applying herbicide alternatives to glyphosate, or applying glyphosate in mixtures with other herbicides may help prevent or delay future origins of resistance by reducing population sizes and reducing the selection pressure of glyphosate in these systems. To limit spread of resistance through gene flow, cleaning weed seed from equipment and shoes moved between infested fields may help reduce long-distance seed transfer. Mowing or tillage may also reduce shortdistance gene flow through pollen dispersal by reducing the number of flowers resistant plants produce, although these management techniques may not be able to entirely eliminate pollen production. It is not yet known whether a fitness cost is associated with glyphosate resistance in these L. perenne populations, and whether the frequency of resistance would be maintained in the absence of continued glyphosate use, as no fitness studies have yet been conducted in California populations of L. perenne and fitness costs associated with herbicide resistance vary depending on mechanism, the measure of fitness used, and genetic background of the population (Giacomini, Westra, & Ward, 2014;Preston et al., 2009;Vila-Aiub, Neve, & Powles, 2009). If any of the mechanisms of resistance present in the sampled area do confer a fitness cost, discontinuing use of glyphosate could result in a gradual decrease in the frequency of resistant individuals in populations.
However, glyphosate is still effective on many other weed species, and will likely continue to be used in weed management in perennial crops.
If gene flow between populations acts to produce populations con-