EPA/600/R-XX/XXX | October 2010 | www.epa.gov/nrmrl

Risk Assessment of the Evolution of 

Resistance in Western corn rootworm 

(Diabrotica virgifera virgifera LeConte) to 

Smartstax® corn 

Letter Report 

Sustainable Technology Division, Cincinnati, OH 45268 

National Risk Management Research Laboratory 

Office of Research and Development
EPA/600/R-XX/XXX 

October 2010

Risk Assessment of the Evolution of 

Resistance in Western corn rootworm 

(Diabrotica vigifera vigifera LeConte) to 

Smartstax® corn 

Letter Report 

October 29, 2010 

 

Michael Caprio1, John A. Glaser2 

1

Mississippi State University, Department of Entomology and Plant Pathology, Mississippi 

State, MS 39762-9775 

2

US Environmental Protection Agency, National Risk Management Research Laboratory, 

Sustainable Technology Division, Cincinnati, OH 45268 

Sustainable Environments Branch 

Sustainable technology Division 

National Risk Management Research Laboratory 

Office of Research and Development 

Cincinnati, Ohio 45258 

U.S. Environmental Protection Agency 

Office of Research and Development 

Washington, DC 20460 

2
Notice

The U.S. Environmental Protection Agency through its Office of Research and Development 

funded and managed the research described here. It has been subjected to the Agency=s peer 

and administrative review and has been approved for publication as an EPA document. Mention 

of trade names or commercial products does not constitute endorsement or recommendation for 

use 

. 

The information provided in this paper is not intended for determining compliance with state or 

federal regulations for Plant Iincorporated Pesticide resistance management. For these purposes, 

information should be obtained from the appropriate local, state and federal regulatory agencies, 

the Federal Insecticide Fungicide and Rodenticide Act, the United States Code, the Code of 

Federal Regulations, the Federal Register and other sources. 

3
Abstract 

Genetically-Modified(GM) crops containing pesticidal traits are valued as environmental 

assets since they directly replace broad spectrum chemical pesticides having detrimental 

environmental and human health effects. As part of the risk assessment for each new GM 

crop registration requires an evaluation of the sustainable lifetime for each new GM crop 

variety. Models employing principles of population genetics have been used to simulate field 

conditions where the development of insect pest resistance to the control traits of the crop is 

encountered. The Monsanto Company’s submission of the Superstax® corn event for control 

of the Western corn rootworm (Diabrotica vigifera vigifera LeConte) has been evaluated for 

the development of insect pest resistance using the EPA developed deterministic two locus 

simulation model POPGEN D2L. The design of research is presented with the interpreted 

results. 

4
Foreword 

The U. S. Environmental Protection Agency (USEPA) is charged by Congress with 

protecting the Nation's land, air, and water resources. Under a mandate of national 

environmental laws, the Agency strives to formulate and implement actions leading to a 

compatible balance between human activities and the ability of natural systems to support and 

nurture life. To meet this mandate, USEPA’s research program is providing data and technical 

support for solving environmental problems today and building a science knowledge base 

necessary to manage our ecological resources wisely, understand how pollutants affect our 

health, and prevent or reduce environmental risks in the future. 

The National Risk Management Research Laboratory (NRMRL) is the Agency's center for 

investigation of technological and management approaches for preventing and reducing risks 

from pollution that threaten human health and the environment. The focus of the Laboratory's 

research program is on methods and their cost-effectiveness for prevention and control of 

pollution to air, land, water, and subsurface resources; protection of water quality in public water 

systems; remediation of contaminated sites, sediments and ground water; prevention and control 

of indoor air pollution; and restoration of ecosystems. NRMRL collaborates with both public and 

private sector partners to foster technologies that reduce the cost of compliance and to anticipate 

emerging problems. NRMRL's research provides solutions to environmental problems by: 

developing and promoting technologies that protect and improve the environment; advancing 

scientific and engineering information to support regulatory and policy decisions; and providing 

the technical support and information transfer to ensure implementation of environmental 

regulations and strategies at the national, state, and community levels. 

This publication has been produced as part of the Laboratory's strategic long-term 

research plan. It is published and made available by USEPA's Office of Research and 

Development to assist the user community and to link researchers with their clients. 

Sally Gutierrez, Director 

National Risk Management Research Laboratory 

5
Contents 

Notice............................................................................................................................... 

Abstract............................................................................................................................ 

Acknowledgements 

Foreword.......................................................................................................................... 

Contents............................................................................................................................ 

List of Figures.............................................................................................................................. 

List of Tables.............................................................................................................................. 

Terminology.................................................................................................................... 

CHAPTER 1 Introduction 

1.1 Background 

1.2 Purpose 

CHAPTER 2 Study design 

 2.1 RM Models 

 2.2 Parameters 

 2.3 Methods and materials 

 2.4 Results 

CHAPTER 3 Discussion and Conclusions 

References....................................................................................................................... 

Appendix 1 Resistance Management Models Used for this Analysis................................ 

Appendix 2 Deterministic Parameter Description................................................... 

Appendix 2 Sensitivity Analysis....................................................................................... 

6
List of Figures 

Figure 1: The distribution of the time taken for resistance to evolve in 1000 simulations of a 

5% seed-mixture. One value of 1000 generations was excluded for the purposes of this 

graph. 

Figure 2: The distribution of the time taken for resistance to evolve in 1000 simulations of a 

5% structured refuge 

Figure 3: Direction comparison of the densities of the distributions of the time until resistance 

evolved for 5% seed-mixtures and 5% structured refuges for the critical first 200 generations. 

7
List of Tables 

Table 1: Fitness values calculated from the most likely values for individual toxins (Table 2) 

and a dominance of 0.05 

Table 2: The PERT-Beta distributions used in the corn rootworm risk assessment 

Table 3: Quantiles for the distributions of time until resistance evolved for 1000 simulations 

each of the seed-mixture and structured refuge models. 

8
Risk Assessment of the Evolution of 

Resistance in Western Corn 

Rootworm(Diabrotica vigifera vigifera LeConte) 

to Smartstax® corn 

Letter Report 

 

 

Michael Caprio1, John A. Glaser2 

1

Misssissippi State University, Department of Entomology and Plant Pathology, Mississippi 

State, MS 39762-9775 

2

US Environmental Protection Agency, National Risk Management Research Laboratory, 

Sustainable Technology Division, Cincinnati, OH 45268 

Keywords: PIP corn variety, PIP crop registration, resistance management, simulation modeling, POPGEN D2L 

CHAPTER 1 Introduction 

The USEPA, under its administration of the Federal Insecticide Fungicide and Rodenticide 

Act (FIFRA), requires the registration of all pesticides and pesticidal materials. The GM or 

Plant Incorporated Protectant (PIP) crops containing pesticidal traits are subject to FIFRA 

registration requirements. Concerns relating to future environmental effects and sustainability 

of these crops can be investigated only using simulation models. Model application in this 

instance offers an analysis of the useful lifetime for the pesticidal traits. USEPA has clearly 

supported the use of these traits as replacements of broad spectrum pesticides that have 

significant environmental footprints. 

1.1 Background

Resistance management (RM) simulation models are designed for use as analytical 

methodologies to provide realistic assessment of the risk of resistance evolution given the 

9
allele frequencies used as initial conditions. New concepts of GM crop deployment such as 

refuge-in-the-bag require a detailed evaluation to ensure the claims associated with each 

concept offer the sustainable protection of the crop. Simulation models offer a means to 

evaluate the threat of crop loss through the development of resistance in the near future. 

1.2 Purpose

This research effort was designed to evaluate the relative merit of resistance management 

concepts as submitted by Monsanto Company relative to the Superstax® corn registration. 

This research specifically analyzed the relative risks related to the control of western corn 

rootworm in seed mixture scenarios 

CHAPTER 2 Study Design 

To understand the relative risks related to the evolution of resistance in western corn 

rootworm for the Superstax® corn variety as submitted by Monsanto Company, a two locus 

frequency-based deterministic model was developed and evaluated. 

2.1 RM Models 

POPGEN-D2L 

The deterministic model POPGEN-D2L is a modification of a previously developed simple 

frequency-based model that incorporates noncompliance, non-random mating in compliant 

fields and larval movement (Caprio & Glaser, 2010a). This dual locus simulation model 

designed to evaluate the population dynamics of crop pest insects relating to the development 

of resistance to pesticidal controls installed in the crop. The model includes a switch to vary 

between a block-refuge design and a seed-mixture design. Biological parameters were 

10
primarily derived from that model with the exception of mortality rates, block size, dispersal 

rates and spatial aspects of refuge placement and movement.

2.2 Methods and Materials 

The risk assessments utilized a two locus deterministic model, an extension of the model 

used previously (Caprio & Glaser 2010 a&b). This was desig1ned as a two compartment 

model, one compartment for hard-core non-compliant fields (continuously non-compliant) and 

a second for compliant and unintentionally non-compliant fields. The latter compartment was 

further divided into two compartments, one for the Bt-plants and the second for refuge (nonBt) plants. This model extended the model used in an earlier risk assessment to two 

independent resistance loci (Caprio & Glaser 2010a&b). The model could be run as either a 

seed-mixture model or a structured refuge model. The risk assessment was conducted by 

running each version 1000 times, each time resampling the parameters listed in Table 2 using 

PERT-Beta distributions (Vose 2000). The value of each parameter and the time until 

resistance evolved was recorded for each simulation. Several options exist for sensitivity 

analysis of this data, including Bayesian model analysis, stepwise regression and 

standardized regression coefficients. Although the use of standardized regression coefficients 

is somewhat controversial, they are appropriate in this case because the variance in each 

variable is determined by the PERT-Beta parameters and is related to the uncertainty in the 

parameter. 

The calculation of the fitnesses in this non-additive interaction calculated as follows. The 

dominance of locus 1 was h1 and the dominance at locus 2 was h2. Note that when locus 1 is 

RR, then the fitness only depends on the locus 2 genotype and vice versa when the locus 2 

genotype is RR. The exterior row and column of the table can therefore be completed from 

data submitted as part of the registration materials because one of the loci is fully resistant 

11
and only one locus is responsible for the mortality. The SSSS genotype fitness can be 

estimated by the overall survivorship to the dual toxin product (also provided). If we look at 

each locus separately, we know the survivorship of an SS is 0.0025, and RR is 1.0. Assuming 

a dominance of 0.2, a heterozygote fitness of 0.202 is determined. We also know that the 

SSSS fitness is 0.0005, and RRRR is 1.0. With these assumptions, the following table can be 

generated: 

 Locus 1 

Locus 2 SS RS RR 

SS 0.0005 A 0.0025 

RS B C (RRRRRRSS)*h2+RRSS

RR 0.0025 (RRRRSSRR)*h1+SSRR 

1.0000 

The remaining difficulty is filling out the cells marked A & B (single heterozygotes) and C (the 

double heterozygote class). Making the assumption that the dominance properties would 

hold for these genotypes as well, the fitness of the SSRS genotype (B) would be calculated 

as h2*(SSRR-SSSS)+SSSS, while the fitness of the RSSS genotype (A) would be calculated 

as h1*(RRSS-SSSS) + SSSS. Finally, the double heterozygote fitness was calculated by 

using an average dominance: RSRS = (h1+h2)/2 * (RSRR – RSSS) + RSSS. 

 Locus 1 

Locus 2 SS RS RR 

SS 0.01800 0.01900 0.03800 

RS 0.02000 0.02330 0.08610 

RR 0.05800 0.10510 1.00000 

Table 1: Fitness values calculated from the most likely values for individual toxins (Table 2) 

and a dominance of 0.05 

12
Variable Min Mode Max Comments

L1-fit 0.008 0.058 0.1

L2-fit 0.0008 0.038 0.08

Overall Survivorship 0.0003 0.018 0.03This will always be as low as lowest of L1 and L2

Dominance 0.01 0.05 0.2Drawn independently for both loci

Initial frequency 0.001 0.005 0.02

These are variables specific to spatial refuges

Refuge size 0 0.03 0.06 Always set to 0.05 if seed-mix

Non-compliance 0 0.3 0.6Always set to 0 for seed-mix

Non-compliance mixing 0.8 0.9 1

Non-random mating 0.2 0.5 0.7Always set to 0 for seed-mix

These are variables specific to seed-mixes

Base larval movement 0.05 0.2 0.4Always set to 0 for spatial refuges

Days of larval selection 5 10 15

Asymmetrical larval movement 0 0.3 0.5 Always set to 0 for spatial refuges

Table 2: The PERT-Beta distributions used in the corn rootworm risk assessment 

13
2.4 Results 

The mean time to resistance for the seed-mixture model was 81.9 generations, while it was 

52.9 generations for the structured refuge simulations. Histograms of the results from the 

1000 simulations of each version of the model were similar in shape (Figures 1 & 2), though 

bootstrapping revealed significant differences (p <<0.01, R statistical package, sm.density 

function). A comparative density plot for the two distributions is shown in Fig. 3. In the case 

of seed-mixtures, 2.5% of the simulations resulted in times until resistance evolved of less 

than 30 generations, while 17.0% of the structured refuge simulations fell below this limit. To 

assist in assessing different levels of risk, different quantiles of the two distributions are show 

in Table 3. The quantiles can be interpreted as levels of acceptable risk. If one accepted a 

5% chance of not meeting a target as an acceptable level of risk, the seed-mixture would be 

acceptable if the target were 33 generations or less, while a structured refuge would only 

meet a target of 22 generations or less. The sensitivity analysis (Appendix 1) indicated that 

all variables except non-random mating and non-compliant mixing had statistically significant 

impacts on the time it took for resistance to evolve. The Bayesian model averaging 

suggested that non-compliant mixing should be included in the model. For seed-mixtures, all 

variable except base larval movement, asymmetrical movement and the time period over 

which larval movement occurred had statistically significant impacts on the time it took for 

resistance to evolve. The Bayesian model averaging suggested that the base larval 

movement rate should be included as a factor. The difference can be seen in the BMA output. 

The best fit model, the one returned by the stepwise regression, did not include base larval 

movement, but the second best fit model did include that variable. In both of the BMA runs 

only the two top models were used as other models had significantly higher absolute 

Bayesian information contents (i.e., they did not fit the data nearly as well). 

14
Quantiles

0 5.00% 10.00% 25.00% 50.00%

Seed-mixture 23 33 39 51 70

Structured refuge 12 22 26 34 45

Table 3: Quantiles for the distributions of time until resistance evolved for 

1000 simulations each of the seed-mixture and structured refuge models. 

15
CHAPTER 3 Discussion

These results provide support for the use of a seed-mixture refuge strategy in a two locus 

corn rootworm system. There are several factors included in the model that made the 

structured refuge less effective than one might have expected. Corn rootworm movement 

and refuge placement, combined with source-sink dynamics, can alter the effective refuge 

size. If growers always plant refuges in the same place year after year, it is likely that that 

refuge will have more than 5% of the total egg lay within the 5% of the total area it contains. 

Alternatively, if the refuge is moved every year, it is possible that very few eggs would be laid 

under the new refuge site. In addition, the proportion of the growers were assumed to be 

non-compliant year after year was increased to a maximum of 60% to account for local areas 

where non-compliance might be higher than average. Some non-random mating 

(independent of effective refuge size) was also included in structured refuge simulations. 

Finally, asymmetrical movement of CRW larvae off of transgenic plants was included in the 

seed-mixture model. This can, at least in some cases, lead to an increased density of larvae 

on the non-Bt plants in a seed-mixture, somewhat increasing the efficacy of that strategy. The 

latter two suggestions, however, are not supported by our sensitivity analysis. It should be 

noted that several of these factors are dependent of specifics of CRW biology (low vagility, 

asymmetrical movement), and caution should be used when extrapolating these results to 

other species. 

It is relevant to note that these data present more risk for dual gene products than might have 

been expected from earlier research (Caprio and Glaser 2010a & b). This is due primarily to 

a change in our perception of the initial resistance allele frequency for western corn rootworm 

16
that was changed to a minimum of 0.001, a most likely value of 0.005 and a maximum value 

of 0.02. 

References 

Caprio, M. A. and D. Sumerford. 2007. Evaluating transgenic plants for suitability in pest and 

resistance management programs, 2nd edition. L. Lacey & H. Kaya [eds]. Field Manual of 

Techniques for the Application and Evaluation of Entomopathogens. 

Caprio, M.A. and J.A. Glaser 2010a, The Use of Simulation Models to Evaluate Pest 

Resistance Development to Refuge in the Bag Concepts related to Pioneer Submission: 

Letter Report, January 2010, US EPA National Risk Management Research Laboratory, 

Cincinnati, OH. 

Caprio, M.A. and J.A. Glaser 2010b, The Use of Simulation Models to Evaluate Pest 

Resistance Development to Refuge in the Bag Concepts related to Pioneer Submission: 

Letter Report Addendum, April 2010, US EPA National Risk Management Research 

Laboratory, Cincinnati, OH. 

Hibbard, BE, MI Higdon, DP Duran, YM Schweikert, and MR Ellersieck. 2004. Role of egg 

density on establishment and plant-to-plant m7(3): 871-882 

Hibbard, BE, TT Vaughn, IO Oyediran, TL Clark, MR Ellersieck. 2005. Effect of Cry3Bb1 

expressing transgenic corn on plant-to-plant movement by western corn rootworm larvae 

(Coleoptera Chrysomeldiae). J. Econ. Entomol. 98(4): 1126-1138 

Vose, 2000. Risk Analysis: A quantitative guide. Wiley, New York. 

17
Figure 1: The distribution of the time taken for resistance 

to evolve in 1000 simulations of a 5% seed-mixture. One 

value of 1000 generations was excluded for the purposes 

of this graph. 

Figure 2: The distribution of the time taken for resistance 

to evolve in 1000 simulations of a 5% structured refuge 

18
Figure 3: Direction comparison of the densities of the distributions of 

the time until resistance evolved for 5% seed-mixtures and 5% 

structured refuges for the critical first 200 generations.

19
Appendix 1 - Resistance Management Model Used for this Analysis 

The deterministic model is a relatively simple frequency-based model that incorporates noncompliance, non-random mating in compliant fields and larval movement. It includes a switch 

to vary between a block-refuge design and a seed-mixture design. Appendix 2 includes 

descriptions of the different parameters. The primary difference between the Bt population in 

the non-compliant population and the Bt-population in the non-random mating population is 

that no refuge was planted for the non-compliant fields. Non-compliant populations reduce 

both non-random mating and refuge size. (Caprio & Glaser, 2010 a&b). 

While we believe that this may overestimate movement of refuges by growers, this was 

incorporated into the model by defining the refuge size not as the physical area devoted to the 

refuge, but rather the proportion of the eggs in a field that start out underneath the new 

refuge. This is the effective refuge size to distinguish it from physical refuge size (the 

proportion of an area planted to a refuge). These small refuge sizes are independent of noncompliance and result from growers moving their refuge to areas of the field were few eggs 

were deposited the previous year. To estimate this important parameter, we used the 

stochastic model to estimate egg densities at the start of a year across our four different field 

sizes. In each case, we estimated a best and worst scenario. For the best case scenario, the 

grower moved their refuge to the opposite side of the field, next to where the refuge would 

have been planted in the neighboring field the previous year (these simulations were run with 

all growers planting refuges to the same sides of their fields and swapping them every year). 

The worst case scenario was based on planting the refuge in the area of the field where the 

egg density was the lowest. Because the egg density distribution across fields has never 

been measured, it is unlikely that growers could accomplish this, but it does represent and 

20
theoretical lower bound for effective refuge size. The worst case effective refuge sizes for the 

four different field sizes were (from smallest to largest field size): 15%, 11%, 3% and 0.055%. 

Few insects were noted to be moving into the center of large fields as reflected by the 

effective refuge size of 0.055% for 144 ha fields. Other estimates for dispersal are higher and 

may depend on the region of the country. To test the effect of higher dispersal rates we 

doubled the vector of daily movement (to 29m/day) and repeated the large field size 

simulation In this case the density of eggs in the center of the field was fifty times greater 

under the worst case scenario and the effective refuge size was 1%. Most fields in the corn 

belt, however, are smaller, in the range of 30-50 ha, and we decided to use a lower limit of 1% 

for this parameter. 

Asymmetrical movement was included in the model, though the amount movement of 

susceptible larvae from transgenic corn is uncertain (Hibbard et al 2004, Hibbard et al 2005). 

Asymmetrical movement was defined as the probability that larvae that have survived 

exposure on a Bt-plant will move. When larvae move, they move onto plant genotypes 

according to the ratio of those genotypes in the field. In this case, a larva leaving a transgenic 

plant had a 95% probability of moving onto another transgenic plant. As with other rates in the 

model, this movement rate is converted into a daily rate (depending on the days of larval 

movement parameter). 

21
Appendix 2 Deterministic Parameter Description 

Wss – The fitness of homozygous susceptible genotypes on the transgenic plants. Without 

better knowledge of field tests it's difficult to develop a range, but a specific factor over 

and under was used to develop our estimate. One question to be asked might be what 

is the lowest/highest performance under field conditions. Instead of 98.75% mortality, 

could it be 95% (as an example)? 

Dominance – The relative fitness of heterozygotes. The traditional population genetics 

equation was employed: 

 Wrs = (Wrr – Wss) * h + Wss 

where h is the dominance and Wrr = 1.0. Our distribution is derived from an estimate of 

h=0.05, but that seems low compared to many estimates based on LC50s. LC50s, however, 

are based on a log scale, however this research has used a linear scale. In general, for 

products with lower effectiveness, dominance tends to be higher ( Caprio et al 2007). 

IRF – The initial resistance allele frequency. My impression from dealing with CRW folks is 

that this value has tended to be higher in empirical estimates than we have seen for other 

insects. 

Block refuge model parameterization: 

Values cited above are shared with the larval movement model and changes in would effect 

both the block and larval movement models to a similar degree. 

22
Refuge proportion – Initially set this to 20%. However, it can be viewed as the proportion of 

eggs in a compliant field that end up being in a refuge (perhaps after being moved), and the 

non-random oviposition patterns can be addressed using this parameter. If more individuals 

survive in refuges, and they tend to lay eggs near that area, then the number of eggs in the 

soil to start the new growing season will be greatest near the refuge from the previous year. If 

the refuge is moved, it stands to reason that the new refuge, even if it covers 20% of the area, 

will receive less than 20% of the eggs oviposited in the field. Based on this modeling 

research, this fact and the likelihood that refuges are moved are identified as key 

assumptions of previous modeling efforts. Unfortunately, there is little empirical data on the 

non-random oviposition patterns of CRW. 

Non-compliance – This is the proportion of fields that are planted without a refuge. 

Individuals emerging from this area mate as a separate population. The most likely value is a 

literature value (CSPI). This parameter should be based on a field by field basis and 

incorporate some estimate of amount of non-compliance. A person planting a 3% refuge in a 

field is not the same as someone not planting any refuge. Alternatively, it shouldn't measure if 

a field is compliant (as yes./no answer), but how compliant a non-compliant field is. General 

comments have been encountered which indicated that 30% non-compliance was too low. It 

is not obvious how non-compliance was estimated. This parameter may be important 

because it determines how much Bt-corn is planted without any refuge. 

Non-compliance mixing – This determines how much movement there is between compliant 

and non-compliant fields. If non-compliant growers are accidentally non-compliant (i.e. they 

never got around to planting the refuge because of weather), then they will generally be 

compliant in the next year, and this parameter should be close to 1.0. On the other hand, if 

23
non-compliant growers are growers who purposefully don't plant refuges and continue this 

practice year-by-year, then this parameter should be low. How low would depend on how 

much movement naturally occurs between non-compliant and compliant fields. Inevitably this 

leads to discussion of the clustering of non-compliant fields. 

Non-random mating – This parameter addresses the proportion of the population in the Btportion of compliant fields that is far enough away from the refuges that they do not mate at 

random. These differ from non-compliant fields in that refuge is planted, but it was too far 

away. The non-random oviposition that results from this non-random mating (and limited 

dispersal capability) is covered by the reductions in refuge size discussed above. 

Larval movement model parameters 

The refuge proportion is fixed at 0.05 in the larval movement model. 

Base larval movement – This is the amount of movement that might occur in a non-transgenic 

field with populations low enough so that significant density-dependent dispersal is not a 

factor. In the stochastic models, population decline is found to be a greater problem than 

density-dependent factors such as mortality or dispersal. The model does include 

asymmetrical movement, with the movement of individuals on Bt plants increased by an 

amount related to their daily mortality rates. Susceptibles therefore move more than 

heterozygotes, and resistant homozygotes would only move at the base larval movement 

rate. Of course, an individual leaving a transgenic plant only has a 5% change of moving 

onto a non-Bt plant. Various estimates of larval dispersal rates have been encountered but 

mostly on the low side. Generally, this research assumes that populations will be low in seed24
mix fields and there will be little to no density-dependent larval dispersal. This could be an 

area of concern since the selection of alternate values is relatively unimportant compared to 

other parameters and a precise estimate would be a misapplication of resources. 

Days of larval movement – Movement and larval development are assumed to occur over this 

period of days. Movement occurs once per day. Daily survivorship is constant and is 

determined as the rate that raised to the days power results in the overall fitness listed for the 

genotype in the original table. In a previous case, that was 0.0125 for susceptible larvae. If 

the days of larval movement is 10, then the daily survivorship rate is 0.645. This means that 

on day 1, 36.5% of the SS larvae on Bt die (assume conservatively that mortality occurs 

before dispersal) , and of the remaining 64.5%, 36.5% will disperse (see above). Of these 

dispersants, only 5% will move onto non-Bt plants, increasing the refuge size from 5% to ca. 

6%. If the base larval movement is low (remember it is spread over days as well), then this 

process can lead to an increase in the refuge population and actually decrease the rate of 

resistance evolution. This feature was installed in the model just to see what effect it might 

have. The value was varied between 5, 10 and 15 and it had very little effect on the results. 

25
Appendix 3 – Sensitivity Analysis 

Structured Refuges 

Call: 

lm(formula = time ~ dominance + IRF + l2dom + l2wss + Non.compliance + 

 non.compliance.mixing + non.random.mating + Overall + Refuge + 

 wss, data = DatasetStruct) 

Residuals: 

 Min 1Q Median 3Q Max 

-33.035 -8.763 -3.390 4.458 229.222 

Coefficients: 

 Estimate Std. Error t value Pr(>|t|) 

(Intercept) 141.953 16.615 8.544 < 2e-16 *** 

dominance -115.819 18.192 -6.366 2.96e-10 *** 

IRF -3417.503 217.135 -15.739 < 2e-16 *** 

l2dom -386.607 18.362 -21.054 < 2e-16 *** 

l2wss -564.326 44.426 -12.703 < 2e-16 *** 

Non.compliance -28.121 5.747 -4.893 1.16e-06 *** 

non.compliance.mixing -27.456 17.147 -1.601 0.110 

non.random.mating 6.533 7.237 0.903 0.367 

Overall 2504.917 119.961 20.881 < 2e-16 *** 

Refuge 825.109 60.190 13.708 < 2e-16 *** 

wss -905.127 38.116 -23.746 < 2e-16 *** 

--- 

Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 18.27 on 989 degrees of freedom 

Multiple R-squared: 0.6627, Adjusted R-squared: 0.6593 

F-statistic: 194.3 on 10 and 989 DF, p-value: < 2.2e-16 

Stepwise regression: 

Df Sum of Sq RSS AIC 

<none> 331361 5865.4 

+ non.compliance.mixing 1 885 330476 5869.6 

+ non.random.mating 1 301 331060 5871.4 

- Non.compliance 1 7979 339340 5882.3 

- dominance 1 13531 344892 5898.5 

- l2wss 1 53272 384633 6007.6 

- Refuge 1 62456 393817 6031.1 

- IRF 1 83047 414408 6082.1 

- Overall 1 145166 476527 6221.8 

- l2dom 1 147904 479265 6227.5 

- wss 1 188771 520132 6309.3 

Call: 

lm(formula = time ~ dominance + IRF + l2dom + l2wss + Non.compliance + 

Overall + Refuge + wss, data = DatasetStruct) 

26
Coefficients: 

 (Intercept) dominance IRF l2dom 

l2wss 

 120.38 -115.78 -3423.40 -386.44 -

559.39 

Non.compliance Overall Refuge wss 

 -28.09 2500.98 823.09 -905.43 

Standardized regression coefficients: 

> lm.beta(RegModel.2) 

 dominance IRF l2dom 

 -0.11778716 -0.29144635 -0.39016175 

 l2wss Non.compliance non.compliance.mixing 

 -0.23694429 -0.09057749 -0.02969303 

 non.random.mating Overall Refuge 

 0.01670632 0.38795030 0.25349491 

 wss 

 -0.44004402 

iterative Bayesian Model Averaging summary: 

Call: 

iBMA.bicreg.matrix(x = as.matrix(DatasetStructnrs[, 1:10]), Y = 

DatasetStructnrs$time, sorted = TRUE, nbest = 10) 

completed 1 iteration 

number of iterations in last call = 100 

variables selected = 1 2 3 4 5 6 7 8 9 

statistics for variables examined so far 

 variable status max Probne0 n times in model 

1 wss selected 100 1 

2 dominance selected 100 1 

3 l2wss selected 100 1 

4 l2dom selected 100 1 

5 Overall selected 100 1 

27
6 IRF selected 100 1 

7 Refuge selected 100 1 

8 Non.compliance selected 100 1 

9 Non.compliance.mixing selected 10.8 1 

10 Non.random.mating dropped 0 1 

 Results of running BMA on final set of selected variables: 

Call: 

bicreg(x = currentX, y = cs$Y, wt = cs$wt, maxCol = cs$maxNvar + 1, 

nbest = 10) 

 2 models were selected 

 Best 2 models (cumulative posterior probability = 1 ): 

 p!=0 EV SD model 1 model 2 

Intercept 100.0 123.102 10.44 120.38 145.52 

wss 100.0 -905.259 38.11 -905.43 -903.86 

dominance 100.0 -115.823 18.20 -115.78 -116.15 

l2wss 100.0 -559.995 44.36 -559.39 -565.01 

l2dom 100.0 -386.431 18.37 -386.44 -386.37 

Overall 100.0 2501.336 120.03 2500.98 2504.25 

IRF 100.0 -3423.297 217.21 -3423.40 -3422.49 

Refuge 100.0 823.258 60.22 823.09 824.67 

Non.compliance 100.0 -28.085 5.75 -28.09 -28.03 

28
Non.compliance.mixing 10.8 -3.023 10.35 . -27.91 

 

nVar 8 9 

r2 0.661 0.662 

BIC -1027.91 -1023.69 

post prob 0.892 0.108 

29
Seed-Mixtures 

Call: 

lm(formula = time ~ AsymmMove + BaseLarvalMovement + dominance + 

 IRF + Overall + SelectionDays + w2dominance + w2ss + wss, 

 data = DatasetLM) 

Residuals: 

 Min 1Q Median 3Q Max 

-63.358 -16.039 -6.685 7.804 810.420 

Coefficients: 

 Estimate Std. Error t value Pr(>|t|) 

(Intercept) 244.0960 12.3854 19.708 < 2e-16 *** 

AsymmMove 17.0510 14.1070 1.209 0.2271 

BaseLarvalMovement -46.9418 20.5009 -2.290 0.0222 * 

dominance -226.3551 36.2211 -6.249 6.12e-10 *** 

IRF -4976.0100 438.0047 -11.361 < 2e-16 *** 

Overall 3576.2111 246.6545 14.499 < 2e-16 *** 

SelectionDays -0.6829 0.7174 -0.952 0.3414 

w2dominance -605.4791 37.3586 -16.207 < 2e-16 *** 

w2ss -1081.5755 90.2482 -11.984 < 2e-16 *** 

wss -1504.3279 76.6823 -19.618 < 2e-16 *** 

--- 

Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 

Residual standard error: 37.16 on 990 degrees of freedom 

Multiple R-squared: 0.5517, Adjusted R-squared: 0.5476 

F-statistic: 135.3 on 9 and 990 DF, p-value: < 2.2e-16 

Stepwise regression: 

Df Sum of Sq RSS AIC 

<none> 1377719 7276.5 

+ BaseLarvalMovement 1 7262 1370457 7278.2 

+ AsymmMove 1 2337 1375381 7281.7 

+ SelectionDays 1 1148 1376570 7282.6 

- dominance 1 57428 1435147 7310.5 

- IRF 1 179238 1556957 7391.9 

- w2ss 1 199628 1577347 7404.9 

- Overall 1 286797 1664515 7458.7 

- w2dominance 1 369432 1747151 7507.2 

- wss 1 538532 1916251 7599.6 

Call: 

lm(formula = time ~ dominance + IRF + Overall + w2dominance + w2ss + 

wss, data = DatasetLM) 

Coefficients: 

(Intercept) dominance IRF Overall w2dominance 

w2ss 

 234.4 -232.7 -4988.5 3543.8 -608.7 -

1084.9 

30
 wss 

 -1511.4 

Standardized regression coefficients 

> lm.beta(RegModel.1) 

 AsymmMove BaseLarvalMovement dominance IRF 

 0.02590690 -0.04890295 -0.13375587 -0.24243465 

 

 Overall SelectionDays w2dominance w2ss 

 0.31162764 -0.02038417 -0.34873674 -0.25629715 

 

 wss 

 -0.42115295 

iterative Bayesian Model Averaging summary: 

Call: 

iBMA.bicreg.matrix(x = as.matrix(DatasetLMnrs[, 1:9]), Y = 

DatasetLMnrs$time, sorted = TRUE, nbest = 10) 

completed 1 iteration 

number of iterations in last call = 100 

variables selected = 1 2 3 4 5 6 7 

statistics for variables examined so far 

 variable status max Probne0 n times in model 

1 wss selected 100 1 

2 dominance selected 100 1 

3 w2ss selected 100 1 

4 w2dominance selected 100 1 

5 Overall selected 100 1 

6 IRF selected 100 1 

7 BaseLarvalMovement selected 29.2 1 

8 SelectionDays dropped 0 1 

9 AsymmMove dropped 4.9 1 

31
 Results of running BMA on final set of selected variables: 

Call: 

bicreg(x = currentX, y = cs$Y, wt = cs$wt, maxCol = cs$maxNvar + 1, 

nbest = 10) 

 2 models were selected 

 Best 2 models (cumulative posterior probability = 1 ): 

 p!=0 EV SD model 1 model 2 

Intercept 100.0 236.96 8.92 234.35 242.85 

wss 100.0 -1509.29 76.75 -1511.40 -1504.51 

dominance 100.0 -231.47 36.20 -232.67 -228.77 

w2ss 100.0 -1084.60 90.39 -1084.93 -1083.85 

w2dominance 100.0 -607.94 37.30 -608.69 -606.27 

Overall 100.0 3546.71 246.37 3543.80 3553.29 

IRF 100.0 -4989.24 438.61 -4988.47 -4990.98 

BaseLarvalMovement 30.7 -14.43 24.46 . -46.97 

 

nVar 6 7 

r2 0.548 0.551 

BIC -753.00 -751.38 

post prob 0.693 0.307 

32