Estimation of Cry3Bb1 resistance allele frequency in field populations of western corn rootworm using a genetic marker.

Following the discovery of western corn rootworm (WCR; Diabrotica virgifera virgifera) populations resistant to the Bacillus thuringiensis (Bt) protein Cry3Bb1, resistance was genetically mapped to a single locus on WCR chromosome 8 and linked SNP markers were shown to correlate with the frequency of resistance among field-collected populations from the US Corn Belt. The purpose of this paper is to further investigate the relationship between one of these resistance-linked markers and the causal resistance locus. Using data from laboratory bioassays and field experiments, we show that one allele of the resistance-linked marker increased in frequency in response to selection, but was not perfectly linked to the causal resistance allele. By coupling the response to selection data with a genetic model of the linkage between the marker and the causal allele, we developed a model that allowed marker allele frequencies to be mapped to causal allele frequencies. We then used this model to estimate the resistance allele frequency distribution in the US Corn Belt based on collections from 40 populations. These estimates suggest that chromosome 8 Cry3Bb1 resistance allele frequency was generally low (<10%) for 65% of the landscape, though an estimated 13% of landscape has relatively high (>25%) resistance allele frequency.

Introduction

Western corn rootworm (Diabrotica virgifera virgifera; WCR) is a significant maize pest in the United States. In its larval phase, WCR feeds on maize roots and can cause substantial yield loss during heavy infestations (Gray ; Gassmann 2016). Various control options are available for WCR, including crop rotation, chemical insecticides, and transgenically expressed Cry toxins from Bacillus thuringiensis (Bt). However, WCR has a history of evolving resistance to control measures, and field-evolved resistance has been detected to each of these control options (Levine ; Gray ; Gassmann ). As a recent example, field-evolved resistance to the transgenically expressed Bt toxin Cry3Bb1 was reported by Gassmann . This was a striking discovery, as plant-expressed Cry3Bb1 is widely used in WCR control products (for example, VT Triple PRO® and SmartStax® corn) and resistance apparently evolved within 6 years of the first commercial use of this protein.

The discovery of Cry3Bb1 resistance in WCR prompted several research activities. Follow-up work showed that Cry3Bb1-resistant WCR were cross-resistant to mCry3A and eCry3.1Ab, other commercially available Bt toxins that belong to the same Cry3 family as Cry3Bb1 (Gassmann ). Moreover, phenotypic assays demonstrated that insects harboring Cry3Bb1 resistance suffer little fitness cost (Hoffmann ; Paolino and Gassmann 2017), suggesting that resistance allele frequencies might not be held in check by natural selection. Despite this work, uncertainty remains about the distribution and trajectory of resistance alleles across the US Corn Belt. Maintaining WCR susceptibility to Cry3Bb1 remains a major priority among all stakeholders (Tabashnik and Gould 2012; Cullen ; Head ; Gassmann 2016). To this end, insect resistance management (IRM) practices have been initiated.

A key feature of IRM is the accurate assessment of resistance allele frequency (Roush and Daly 1990; Head and Greenplate 2012). It is common to use bioassays to monitor resistance allele frequencies, and there has been success using this approach in WCR (Gassmann ; Gassmann ; Shrestha ). However, bioassays can also be costly, time-consuming, and sensitive to subtle variation in environmental factors (Meihls ). A promising alternative approach is to use genetic markers that tag resistance alleles. Genetic markers have several advantages over bioassays. First, because of advances in genotyping technologies, genetic markers are affordable and can be rapidly deployed. Genetic markers can also scale to a large number of populations more readily than bioassays. Finally, markers can be used identify historical resistance allele frequencies in preserved specimens. This type of analysis can yield helpful insights about resistance allele frequencies in materials collected before resistance was first detected (e.g., Banerjee ). Despite these benefits, genetic markers have a higher upfront cost than bioassays because they require specific genetic knowledge of the resistance locus.

Fortunately, much of this upfront work has been done for Cry3Bb1 resistance in WCR. Using a Cry3Bb1 resistant colony collected near Hopkinton, Iowa (Gassmann ), resistance was mapped to a single locus on WCR chromosome 8 (Flagel ). Resistance at this locus is nearly recessive, and markers tightly linked to the trait were described and shown to correlate with Cry3Bb1 resistance in other field populations. Moreover, the identified chromosome region contains a candidate gene known as ABCb1 (Flagel ). This gene has been implicated in Cry3Aa resistance in Chrysomela tremula (Pauchet ), a coleopteran species like WCR. Though this gene is a promising candidate, to date, the causal Cry3Bb1 resistance allele in WCR has not been identified. Thus, the link between the identified markers and the true resistance locus remains correlative, which needs to be taken into account to utilize genetic markers for estimating Cry3Bb1 resistance allele frequencies in WCR.

The purpose of this paper is to further investigate the relationship between a linked resistance marker identified by Flagel and the causal resistance locus itself by coupling response to selection data from laboratory and field experiments with genetic models. Our objectives were to (1) confirm the association between the marker and the causal locus, (2) test whether the marker allele associated with resistance is perfectly linked with the causal resistance allele, and (3) describe and fit a theoretical model accounting for imperfect linkage that allows mapping from marker allele frequency to causal allele frequency.

Materials and methods

Marker development and genotyping assays

We used the CRW1862 TaqMan marker described in Flagel . This marker tags a C/T SNP polymorphism on WCR chromosome 8 at approximately 126.7 cM on the genetic map, which is very near the Cry3Bb1 resistance QTL peak described in Flagel and approximately 7 cM from the ABCb1 candidate gene. To assay this marker, we followed the DNA extraction and genotyping protocols described in Flagel .

Plant bioassays

Adult WCR were collected from 16 grower fields (Figure 1) with confirmed greater-than-expected damage (GTED) to Cry3Bb1-expressing maize. For each population, the collected adults were mated and their offspring tested in a plant assay as previously described (Flagel ). These assays included two treatments—Cry3Bb1-expressing maize and near-isoline non-Cry3Bb1-expressing maize (control)—using 10 replicates of 30 neonate larvae for each treatment for each population. WCR shows low soil movement (<1 m; Strnad and Bergman 1987), and plants from both treatments were grown in tents with isolation distances of several meters. Individuals that survived to adulthood were genotyped for the CRW1862 marker (Table 1).

Figure 1

Locations of WCR populations used in this study.

Table 1

CRW1862 genotype (CC, CT, and T) and allele (C and T) frequencies for individuals recovered from Cry3Bb1-expressing plants and control plants in plant assay.

	Control-exposed insects							Cry3Bb1-exposed insects
Population	N	CC	CT	TT	C	T	TAF	N	CC	CT	TT	C	T	TAF
14016	44	13	18	13	44	44	0.5	44	2	12	30	16	72	0.82
14016r	13	5	3	5	13	13	0.5	44	6	22	16	34	54	0.61
150014651	44	22	15	7	59	29	0.33	19	7	8	4	22	16	0.42
15002	30	1	13	16	15	45	0.75	44	1	18	25	20	68	0.77
150023420	44	8	20	16	36	52	0.59	42	5	24	13	34	50	0.6
150027567	44	26	15	3	67	21	0.24	36	12	17	7	41	31	0.43
15003	22	8	8	6	24	20	0.45	19	1	6	12	8	30	0.79
15005	30	4	14	12	22	38	0.63	41	1	12	28	14	68	0.83
15008	43	12	24	7	48	38	0.44	44	3	15	26	21	67	0.76
15010	26	12	12	2	36	16	0.31	37	7	19	11	33	41	0.55
15011	44	17	16	11	50	38	0.43	44	6	18	20	30	58	0.66
15014	24	9	12	3	30	18	0.38	40	4	25	11	33	47	0.59
15015	31	4	14	13	22	40	0.65	44	1	14	29	16	72	0.82
15016	32	5	18	9	28	36	0.56	21	1	10	10	12	30	0.71
15020	43	18	19	6	55	31	0.36	29	4	19	6	27	31	0.53
15021	44	5	21	18	31	57	0.65	44	1	18	25	20	68	0.77

N is the number of genotyped individuals. TAF is the CRW1862 T-allele frequency. In all populations, TAF in the Cry3Bb1-exposed insects was equal to or greater than that in the control-exposed insects.

If CRW1862 is linked to a causal gene for Cry3Bb1 resistance, there should be a difference in allele frequencies between the two treatments. Specifically, because the T allele of CRW1862 was previously found to be associated with Cry3Bb1 resistance (Flagel ), we would expect the frequency of the T allele to be greater in individuals recovered from Cry3Bb1-expressing plants than in individuals recovered from control plants.

Logistic regression was used to compare the relative frequencies of T alleles between individuals recovered from Cry3Bb1 with those in the control group. Let the random variable if an individual was recovered from a Cry3Bb1 plant and if individual was recovered from a control plant, and let be the probability that an allele recovered from treatment is a T allele. Allele shifts within populations were assessed using the model: for the population and the jth individual within that population, where would suggest that the marker is linked to the causal gene.

To test for an average allele frequency shift across populations, the following model with constant slope parameter was used:

Again, would suggest that on average, across populations, the marker is linked to the causal gene.

GTED field experiment

The plant assay was used to statistically test for a shift in marker allele frequency following Cry3Bb1 selection. If there is a shift, it remains to be determined whether the marker allele is perfectly linked to the causal allele (i.e., marker allele frequency equals resistance allele frequency), or whether the linkage is imperfect. If the marker allele is perfectly linked with the causal allele, it should be possible to predict changes in genotype frequencies given knowledge of selection parameters. If are genotype frequencies for the susceptible homozygote, heterozygote, and resistant homozygote, respectively, before Cry3Bb1 selection, and are their corresponding relative fitness in the presence of Cry3Bb1 selection, then the expected genotype frequencies after selection are: with The corresponding resistance allele frequency before selection is and after selection is .

To investigate the relationship between marker genotype frequencies before and after selection, we used data from field experiments where the selection parameters are well characterized. Trials were conducted at seven grower fields (Figure 1) with previously observed GTED fields. In each trial, non-Cry3Bb1-expressing corn and Cry3Bb1-expressing corn (either MON88017 or MON87411) were grown in separate plots, with cages placed over each plot to capture emerging adults. Surviving adult beetles were genotyped for the CRW1862 marker.

A likelihood approach was used to formally test the hypothesis that the CRW1862 T allele is perfectly linked with resistance. Let and be the number of genotyped individuals emerging from control and Cry3Bb1 plots, respectively; and let and be vectors of corresponding observed genotype counts. If are genotype frequencies before selection and are selection parameters, then follows the multinomial distribution: and, if the marker allele is perfectly linked with the causal allele, follows the multinomial distribution: where .

Thus, under the assumption that the marker allele is the causal allele or is perfectly linked to the causal allele, the likelihood function is:

To test the consistency of GTED field trial data with the hypothesis that the marker allele is (or is perfectly linked with) the causal allele, we built the following formal statistical hypothesis test: vs . Under field conditions, the efficacy of Cry3Bb1 in commercial hybrid maize is well studied, and it has been consistently found that (Head ). Thus, rejecting in favor of strongly implies that the marker allele is not perfectly linked with the causal Cry3Bb1 resistance allele. To test the statistical hypothesis, the five parameters need to be estimated to maximize the likelihood , though under the parameter space of is restricted. By taking the logarithm of and equating the first derivative with 0, we obtain five nonlinear equations of the five unknowns, from which estimates are obtained by numerically solving these equations. To test vs , the maximized likelihood is compared with the maximized likelihood .

Causal linkage model

If the suspected resistance allele is not perfectly linked to the causal allele, then it is necessary to find a model that maps the marker allele frequency to the causal allele frequency. We developed a three-haplotype model to search for this mapping. For example, if (T allele of CRW1862) and (C allele of the CRW1862) are the marker alleles associated with Cry3Bb1 resistance and susceptibility, respectively, and and are the corresponding causal Cry3Bb1 resistance and susceptibility alleles, respectively, then one possibility is:

In this model, is the probability of seeing the resistance-linked T allele at CRW1862 when the Cry3Bb1 susceptibility allele is present at the causal locus. If then the marker allele frequency is greater than the causal allele frequency. Marker genotype probabilities, as a function of causal genotype probabilities, are given by:

We use a maximum likelihood approach to estimate parameters. Using the same notation as in previous section, with , then follows the multinomial distribution: and follows the multinomial distribution:

It follows that the likelihood function is:

Note that if , then marker genotype frequencies are same as causal genotype frequencies, distributions (7) and (8) are equivalent to distributions (4) and (5) from the previous section, and likelihood (9) is equivalent to likelihood (6). We used this model to estimate with the constraint .

This model assumes the existence of three haplotypes: , , and . The haplotype represents the CRW1862 T allele () linked with the causal allele () for Cry3Bb1 resistance. Similarly, is the CRW1862 T allele linked to the susceptibility allele () for Cry3B1 resistance, and is the CRW1862 C allele linked to the susceptibility allele (Figure 2). The model assumes the remaining haplotype () is absent or so rare that it can be ignored. Assuming the haplotype (C allele at CRW1862 linked to casual allele for Cry3Bb1 resistance) does not exist before Cry3Bb1 selection began, then the only way it can emerge it is from crossover in an X heterozygote. Flagel estimate the average fitness of the heterozygote at the causal locus (i.e. ) is 12% of the homozygote (i.e. ) under Cry3Bb1 selection. They also estimate that CRW1862 is within a few centimorgans of the resistance allele (conservatively we will use 5 cM). If we assume Hardy–Weinberg genotype proportions, and that only and haplotypes exist (neglecting , which cannot produce haplotypes, and thereby making our work below conservative), then the frequency of the X heterozygote individual in the population can be modeled as 2*(1−)*. Accounting for linkage and the fitness of this heterozygote, we would expect the influx of new gametes under Cry3Bb1 selection to be 2*(1−)**0.12*0.05 each generation. An interesting value to extract from this function is max{/}. This occurs as approaches zero and / = 1.2%. So the production of haplotypes under Hardy–Weinberg is limited by the frequency of in the population and cannot exceed 1.2% of the frequency of in a Cry3Bb1 environment. In a non-Cry3Bb1 environment, we take away the relative fitness term, and max{/} = 0.1, again at the point where approaches zero. The landscape is a mix of these two selective environments. This model suggests that the production of haplotypes is too rare in the short timeframe between the emergence of Cry3Bb1 resistance and our population sampling to impact our estimates of frequency in a meaningful way.

Figure 2

Temporal representation of hypothesized three-haplotype model. The x-axis represents time from past to present. The genealogical tree tracks the birth of new haplotypes near the CRW1862 marker and the closely linked Cry3Bb1 resistance locus. Each present-day haplotype is listed at the tips of the tree. According to this model, first, the CRW1862 SNP occurred and created two haplotypes: the C (green) and T (purple) alleles at this marker. These alleles spread, and in the recent past, another mutation occurred causing Cry3Bb1 resistance in the T-allele background of CRW1862, creating two new haplotypes (for a total of three). S, susceptibility-linked marker allele, susceptibility to Cry3Bb1; R, resistance-linked marker allele, resistance to Cry3Bb1; R, resistance-linked marker, susceptibility to Cry3Bb1.

Field collected populations

To estimate Cry3Bb1 resistance allele frequencies, adult insects from 40 populations of WCR were field collected across the US Corn Belt (Figure 1). All insects we collected at adulthood and some likely came from Cry3Bb1 expressing maize as that is prevalent in the landscape. These populations are a random sample representative of the mating population for the next generation of WCR.

Data availability

All data generated for this project (genotype and allele frequencies) are available in the tables or Supplementary Table S1.

Supplementary material is available at G3 online.

Results

Plant bioassay

To estimate Cry3Bb1 resistance allele frequencies in WCR in the US Corn Belt, we used a Cry3Bb1 resistance-linked TaqMan marker called CRW1862. This marker is tightly linked with a Cry3Bb1 resistance locus on chromosome 8, as first described by Flagel . It tags a C/T SNP where the T allele was shown to be positively correlated with resistance in the field. Our first goal was to determine if the CRW1862 marker is linked with Cry3Bb1 resistance in these populations.

CRW1862 genotypes were obtained for individuals recovered in plant bioassays from Cry3Bb1 plants and near-isoline control plants for each of the 16 populations (Figure 1). Table 1 shows the genotype and allele frequencies for each treatment, and CRW1862 T-allele frequencies are shown in Figure 3. The T-allele frequency was greater in the Cry3Bb1 exposure group for all populations, suggesting that it had responded to selection in these populations. In an across-population logistic regression analysis (equation 2), was estimated to be 0.813 with standard error 0.091 , indicating a statistically significant increase in T-allele frequency across populations upon exposure to Cry3Bb1. In within-population logistic regression analyses (equation 1), was estimated to be significantly greater than 0 (indicating a statistically significant increase in T-allele frequency) for 10 of the 16 populations (Table 2). Consistent with the across-population regression analysis, the six nonsignificant populations nevertheless showed an increase in the T allele in the Cry3Bb1-exposed group. These data confirm that the CRW1862 T allele is positively correlated with Cry3Bb1 resistance in the field.

Figure 3

CRW1862 T-allele frequencies for larvae recovered from Cry3Bb1-expressing plants (red triangles) and control plants (blue triangles) in plant assay. Population IDs are listed on the y-axis, with populations sorted by T-allele frequencies in control plants.

Table 2

Estimated logistic regression parameters (Model 1) for plant bioassay.

Population	α^	SE (α^)	β^	SE(β^)	P-value
14016	0	0.213	1.504	0.349	0
14016r	0	0.392	0.463	0.449	0.303
150014651	−0.71	0.227	0.392	0.399	0.326
15002	1.099	0.298	0.125	0.392	0.749
150023420	0.368	0.217	0.018	0.311	0.954
150027567	−1.16	0.25	0.881	0.345	0.011
15003	−0.182	0.303	1.504	0.5	0.003
15005	0.547	0.268	1.034	0.397	0.009
15008	−0.234	0.217	1.394	0.331	0
15010	−0.811	0.3	1.028	0.381	0.007
15011	−0.274	0.215	0.934	0.311	0.003
15014	−0.511	0.298	0.864	0.375	0.021
15015	0.598	0.265	0.906	0.383	0.018
15016	0.251	0.252	0.665	0.424	0.117
15020	−0.573	0.225	0.711	0.346	0.04
15021	0.609	0.223	0.615	0.338	0.069

Reported P-values are for the test Small P-values suggest that the CRW1862 marker is linked to the causal gene for Cry3Bb1 resistance. Across all populations, .

Flagel posited that CRW1862 is itself likely not the causal Cry3Bb1 resistance causing mutation nor that it is perfectly linked with the causal allele. To test this, CRW1862 genotypes were obtained for adult beetles recovered from Cry3Bb1-expressing plants and near-isoline control plants for seven GTED field populations. Table 3 shows the genotype and allele frequencies for each treatment. To informally investigate the relationship between marker genotype frequencies before and after selection, we compared the observed T-allele frequencies after selection with predicted T-allele frequencies under the assumption that the marker allele is perfectly linked to the causal allele and that (where is the relative fitness of the susceptible homozygote at the causal locus) and the dominance of resistance is or (Table 4). In all cases but one (Field D with ), the observed frequency was lower than the predicted frequency, suggesting either that is substantially larger than 0.05 or that the marker allele is not perfectly linked to the causal resistance allele.

Table 3

CRW1862 genotype (CC, CT, and T) and allele (C and T) frequencies for adult beetles emerging from Cry3Bb1 expressing plants and control plants in GTED field trials.

	Control-exposed insects							Cry3Bb1-exposed insects
Field	N	CC	CT	TT	C	T	TAF	N	CC	CT	TT	C	T	TAF
A	47	14	25	8	53	41	0.44	54	13	31	10	57	51	0.47
B	48	9	15	24	33	63	0.66	48	2	15	31	19	77	0.8
C	44	24	16	4	64	24	0.27	45	15	16	14	46	44	0.49
D	32	9	23	0	41	23	0.36	46	13	28	5	54	38	0.41
E	44	16	23	5	55	33	0.38	44	11	20	13	42	46	0.52
F	41	8	15	18	31	51	0.62	43	0	19	24	19	67	0.78
G	44	1	13	30	15	73	0.83	43	1	12	30	14	72	0.84

N is the number of genotyped individuals. TAF is the CRW1862 T-allele frequency. In all populations, TAF in Cry3Bb1-exposed insects was equal to or greater than that in the control-exposed insects.

Table 4

Observed T-allele frequencies for Control and Cry3Bb1-exposed populations, compared with predicted assuming complete linkage T-allele frequencies for Cry3Bb1-exposed populations assuming and dominance or .

			Predicted Cry3Bb1
Field	Control	Cry3Bb1	h = 0.05	h = 0.50
A	0.44	0.47	0.87	0.68
B	0.66	0.8	0.97	0.87
C	0.27	0.49	0.73	0.61
D	0.36	0.41	0.35	0.48
E	0.38	0.52	0.8	0.63
F	0.62	0.78	0.96	0.84
G	0.83	0.84	0.99	0.91

To more formally test this relationship, we constructed a likelihood ratio test. This test relies on computing two likelihoods given the data, where and where (see equation 6 in Materials and Methods). Noting that is approximately distributed as Table 5 shows the result of likelihood ratio tests conducted separately for fields A–G. All but one of the seven GTED fields had a small P-value (P < 0.1), indicating that is preferred over A composite likelihood ratio test across all fields also suggested preference for (P < 1e−7). This in turn implies that if the CRW1862 T allele is perfectly linked to the causal allele, then (i.e., Cry3Bb1 mortality is <95% on homozygous susceptible insects). Again, this result is incongruous with multiple independent field experiments showing that is 0.05 (Head ). An alternative and simpler explanation is that the CRW1862 T allele is not in perfect linkage with the causal locus for Cry3Bb1 resistance. This explanation is also consistent with the findings of Flagel , where they demonstrated a moderate correlation between resistance frequency and an the CRW1862 T-allele frequency among a different set of field populations. Together, these results argue the need for a model that can map the marker allele frequency to the causal allele frequency.

Table 5

Likelihood ratio test statistics and P-values for comparing hypotheses vs .

Field	-2(L^_0-L^_A)	P-value
A	35.621	0
B	2.966	0.085
C	8.498	0.004
D	30.768	0
E	16.568	0
F	0	1
G	3.441	0.064

Small P-value suggests that .

Results of the causal linkage model

Based on the evidence above, we concluded that the CRW1862 T allele is not the causal Cry3Bb1 resistance allele nor is it in perfect linkage with it. To estimate the Cry3Bb1 resistance allele frequency on WCR chromosome 8 given the CRW1862 marker allele observations, we developed a three-haplotype model that accounts for incomplete linkage (see Materials and Methods). Using this three-haplotype model, our goal was to estimate the proportion of CRW1862 T alleles that are linked to the causal resistance allele. We used a likelihood approach to fit the model (equation 9). The parameter we estimated, θ, is the probability of seeing the T allele at CRW1862 when the susceptibility allele is present at the causal locus. Our estimate of θ (i.e. ) was made separately for each GTED field (Table 6 and Supplementary Table S1). Also, because is a function of the estimated T-allele frequency at the CRW1862 marker (), we rescaled to for each population. The estimate of is greater than 0.5 for all but one population, suggesting that in general, the causal resistance allele frequency is substantially lower than the observed T-allele frequency for marker CRW1862.

Table 6

Maximum likelihood estimates of θ in Causal Linkage Model () under the assumption that , and ratio of estimated θ () to estimated CRW1862 T-allele frequency ().

Field	p^	θ^	θ^/p^
A	0.436	0.433	0.993
B	0.656	0.565	0.861
C	0.273	0.151	0.553
D	0.359	0.354	0.986
E	0.375	0.3	0.8
F	0.622	0	0
G	0.83	0.829	0.999

The ratio , which takes values between 0 and 1, is a relative measure of deviation between T-allele frequency and causal allele frequency. A high ratio suggests a large deviation, i.e., that the causal allele frequency is substantially lower than the T-allele frequency.

Estimates of were used to estimate the distribution of resistance allele frequencies for 40 populations collected between 2013 and 2015 in fields with no reported Cry3Bb1 performance issues (Figure 1). The purpose of this analysis was to estimate the levels of resistance allele frequencies across the US Corn Belt. For each population, each of the seven estimated values (corresponding to the seven GTED populations described above) was used to compute the adjusted causal resistance allele frequency using the equation:

Thus, for each population, there are seven adjusted causal allele frequency estimates (Supplementary Table S1). The distribution of all 40 • 7 = 280 adjusted allele frequencies is shown in Figure 4 (solid line), which has mean value of 0.098 and median 0.042. For comparison, the distribution of the observed T-allele distribution for marker CRW1862 is shown as a dashed line, which has mean of 0.32. The estimated causal allele frequency distribution is substantially lower than the CRW1862 T-allele frequency with 65% of its density <10%.

Figure 4

Estimated distribution of resistance allele frequency (RAF) in 2013–2015 (dashed line) based on 40 fields without product performance inquiries (non-PPI fields) compared to distribution of observed T-allele frequency distribution (solid line). Both distributions are drawn using kernel density estimation from 280 allele frequency estimates.

Discussion

Genetic markers can be a valuable tool for IRM. Once developed, they can offer a cost-effective and scalable platform for diagnosing resistance allele frequencies in natural populations. Moreover, genetic markers can be run on DNA extracted from dead insects, such as historical specimens or preserved insects collected from the field. This flexibility opens access to new sources of information and can simplify field collection strategies. Despite their advantages, however, genetic markers require more upfront work and an understanding of the genetics of resistance.

Here we have shown that the T allele of the CRW1862 marker identified in a previous study (Flagel ) responds to Cry3Bb1 selection across 16 WCR populations, confirming that CRW1862 is linked to Cry3Bb1 resistance. However, we also found compelling evidence that the T allele for CRW1862 is not itself the causal allele or in perfect linkage with the causal allele, as also predicted by Flagel . Specifically, we showed that if the T allele was itself the causal allele or in perfect linkage with it, then the field efficacy of Cry3Bb1 would have to be substantially less than 95% (i.e., > 0.05). Were this true, it would contradict considerable evidence about Cry3Bb1 performance in the field (Head ). Thus, the most parsimonious explanation of the data is to accept that the T allele of CRW1862 is only partially linked with the causal resistance allele on WCR chromosome 8. Accepting this hypothesis makes the observed allele frequency shifts under Cry3Bb1 selection fit better with established work on Cry3Bb1, but also complicates estimation of the true underlying resistance allele frequency with the CRW1862 marker.

To address the lack of perfect linkage between the CRW1862 marker and the causal Cry3Bb1 resistance allele, we developed a three-haplotype model. This model was fitted to our observed allele frequency responses to Cry3Bb1 selection among the 7 GTED populations using a maximum likelihood approach. From this model, we estimated a parameter (θ) for each population, which serves as a map between marker allele frequencies and underlying causal Cry3Bb1 resistance allele frequencies. Finally, with our estimates of θ, we mapped CRW1862 marker frequencies to causal resistance allele frequencies in 40 populations collected between 2013 and 2015 in fields with no reported Cry3Bb1 performance issues (Figure 4). From this mapping, we estimated the Cry3Bb1 resistance allele frequency between 0% and 20% in most populations, with a mean of approximately 10% and a substantial fraction of the estimates at <10%.

If—as these results suggest—a substantial portion of WCR populations still have relatively low Cry3Bb1 resistance allele frequencies, then plant-expressed Cry3Bb1 should still provide some level of WCR control in most maize fields where it is used and, where combined with other modes of action such as the Cry34/35Ab1 proteins (as in the product SmartStax®), WCR control generally should be sufficient to keep WCR populations below economic thresholds (Head ). This is consistent with the low percentage of SmartStax® maize fields reported to have product performance issues (<0.1%, Head ), and low overall root damage ratings in SmartStax® fields (Head ). For this reason, one focus of IRM strategies for WCR has been to transition from single mode of action products containing only Cry3Bb1 or a similar Bt protein to pyramided products like SmartStax® maize with two or more Bt proteins (EPA 2017). A second IRM strategy has been crop rotation, which has been shown to substantially decrease problem fields associated with Cry3Bb1 resistance (Carrière ). Our observation of relatively low Cry3Bb1 resistance allele frequencies is consistent with the joint beneficial effects of Bt trait pyramiding and crop rotation, which have both been implemented aggressively since Cry3Bb1 resistance was discovered in 2011.

The variability in Cry3Bb1 resistance allele frequency observed among WCR populations (Figure 4 and Supplementary Table S1) is consistent with reported low levels of WCR adult movement among maize fields (Spencer ). Resistance in a given field appears to be primarily a function of the products and practices used in the field over time and less a function of immigration of resistant beetles (Gassmann ). This observation supports a field-by-field IRM strategy for WCR in which additional control measures such as crop rotation are differentially applied to fields with the highest resistance allele frequencies and product performance issues (Head ; EPA 2017).

For the time being, the genetic basis of Cry3Bb1 resistance in WCR is not understood. Thus, to estimate Cry3Bb1 resistance allele frequency with genetic markers, we had to make several assumptions and we must acknowledge the limitations of this approach. First, our estimates apply only to the resistance locus on WCR chromosome 8 that was previously mapped by Flagel . It is possible that there are other resistance loci or modifiers of resistance in the WCR genome, and our approach cannot estimate the impact of these should they exist. Flagel attempted to map modifiers in a Cry3Bb1-resistant WCR colony from Hopkinton, Iowa, and were unable to detect any. On the other hand, Ingber and Gassmann (2015) found that Cry3Bb1 resistance in populations from Cresco and Hopkinton, Iowa, was considerably less recessive than that reported by Flagel . One explanation for this could be that the colonies studied by Ingber and Gassmann (2015) have different genetics than those studied in Flagel ; e.g., they may contain novel resistance alleles or genetic modifiers. Were this true, it could mean that our analyses only capture some of the total Cry3Bb1 resistance segregating in the 40 field sampled populations with no reported Cry3Bb1 performance issues.

Our second assumption was that there are only three appreciable haplotypes containing the CRW1862 marker and the Cry3Bb1 resistance allele on chromosome 8. Excluding the haplotype (i.e. the C allele of CRW1862 linked with the causal Cry3Bb1 resistance allele) made it possible to develop likelihood equations to estimate θ, which in turn allowed us to map marker allele frequencies to causal resistance allele frequencies. However, it is possible that this assumption does not hold for all populations studied, and if so could cause us to underestimate Cry3Bb1 resistance allele frequencies. One reason may be a recombination event that creates the haplotype, though as discussed in the Materials and Methods, we expect this to be a rare event. That said, another mechanism that could create the haplotype is either of two varieties of a “soft sweep.” The first variety is often referred to as a multiple origin soft sweep, which is when multiple independent resistance-producing mutations arise in the same gene (Pennings and Hermisson 2006). This scenario could give rise to some populations in which the T allele of CRW1862 is linked to a particular Cry3Bb1 resistance allele at the causal locus, and others in which the C allele of CRW1862 is linked to a different Cry3Bb1 resistance allele at the same locus but with an independent mutational origin. In lepidopteran insects, there are several examples of multiple origin soft sweeps in response to various insecticides (Flagel ; Walsh ; Boaventura ). The second variety of soft sweep could occur if the Cry3Bb1 resistance allele arose prior to the introduction of Cry3Bb1 expressing maize, which could give the resistance allele sufficient time to recombine on to diverse haplotypic backgrounds (Hermisson and Pennings 2005; Pennings and Hermisson 2006). The net result of either of these soft sweep scenarios is that different populations could harbor different resistance haplotypes or even mixtures of resistance haplotypes. Among the GTED populations in Table 3, population D showed a slight increase in CC genotype frequency in response to Cry3Bb1 selection. This could be simply be a sampling anomaly, but it could also be a population that harbors the haplotype and may be a good candidate for future follow-up.

Finally, the third assumption of our methods was to assume similar patterns of linkage between the CRW1862 marker and the Cry3Bb1 resistance allele across a broad geography (Figure 1). We used GTED populations to estimate θ, and these populations came from Iowa where WCR resistance was first reported (Gassmann ). We then extrapolate from these Iowa-based θ estimates to adjacent states. Evidence from our Plant Bioassay populations suggests that the CRW1862 T allele increases in response to Cry3Bb1 selection in many of the adjacent states, which would suggest similar patterns of linkage (Table 1vsTable 3). We also used estimates from seven GTED populations to sample any potential variation in linkage patterns (Figure 4 and Supplementary Table S1). However, we may under- or overestimate Cry3Bb1 resistance allele frequencies in any location where the linkage patterns strongly differ from those found in the seven GTED populations.

Monitoring Cry3Bb1 resistance in WCR remains critical to IRM strategies for this species. Here we show that the CRW1862 marker can be used to estimate resistance allele frequencies, with the caveats provided in this report. While we can continue to use CRW1862 for this purpose, future research should investigate potential improvements. For example, one approach may be to discover additional linked markers for the purposes of defining a multimarker haplotype or haplotypes linked to the causal resistance allele. This could increase precision to detect true resistance haplotypes as compared to the CRW1862 marker alone, particularly if soft sweeps are present. Another step that should be taken is to further investigate the WCR ABCb1 gene in resistant and susceptible colonies. There is compelling evidence that this gene may be involved in Cry3Bb1 resistance in WCR (Flagel ; Pauchet ), and if this were to be confirmed, it could allow the direct detection of the causal Cry3Bb1 resistance allele(s).

Click here for additional data file.