Bt cotton area contraction drives regional pest resurgence, crop loss, and pesticide use.

Genetically-modified crops expressing Bacillus thuringiensis (Bt) proteins have been widely cultivated, permitting an effective non-chemical control of major agricultural pests. While their establishment can enable an area-wide suppression of polyphagous herbivores, no information is available on the impact of Bt crop abandonment in entire landscape matrices. Here, we detail a resurgence of the cosmopolitan bollworm Helicoverpa armigera following a contraction of Bt cotton area in dynamic agro-landscapes over 2007-2019 in North China Plain. An 80% reduction in Bt cotton was mirrored in a 1.9-fold increase of ambient H. armigera population levels, culminating in 1.5-2.1-fold higher yield loss and a 2.0-4.4-fold increase in pesticide use frequency in non-Bt crops (i.e. maize, peanut, soybean). Our work unveils the fate of herbivorous insect populations following a progressive dis-use of insecticidal crop cultivars, and hints at how tactically deployed Bt crops could be paired with agro-ecological measures to mitigate the environmental footprint of crop production.

Introduction

Since its first commercial use in 1996, genetically‐modified (transgenic) crops expressing insecticidal Bacillus thuringiensis (Bt) proteins (hereafter Bt crops) have been adopted worldwide and are presently grown on more than 100 million hectares (ISAAA, 2018). Bt crops have provided a lasting suppression of several (target) herbivorous pests, thereby protecting crop yields, curbing synthetic pesticide use, and enhancing farming profitability (Brookes, 2020, 2020; Klümper and Qaim, 2014). In addition to the net positive impacts on plant protection and societal well‐being (Brookes and Barfoot, 2020; Li et al., 2020; Smyth and Kerr, 2015), Bt crops safeguard biodiversity and bolster ecological resilience of agro‐ecosystems by alleviating the agrochemical burden on cropland (Lu et al., 2012; Romeis and Meissle, 2006).

Scientific research has documented a myriad of ecological and socio‐economic benefits of Bt crops at field, farm and agro‐landscape levels (National Academies of Sciences and Medicine, 2016). As Bt crops establish prevalence in agricultural landscapes, they can lower the abundance of target herbivores (Romeis and Shelton, 2008), attaining parallel declines in specialist feeder populations in both Bt and non‐Bt cultivars of a given crop (Carrière et al., 2003; Tabashnik et al., 2021; Wan et al., 2012). For generalist pests such as the corn borer Ostrinia nubilalis and the earworm Helicoverpa zea, a widespread adoption of Bt maize reduces infestation pressure in several non‐maize crops (e.g. high‐value horticultural commodities) in the mid‐Atlantic United States (Dively et al., 2018; Hutchison et al., 2010). Similarly, in the North China Plain, cotton is the preferred host crop of the polyphagous bollworm Helicoverpa armigera, and Bt cotton cultivation has achieved an area‐wide suppression of H. armigera in multiple crops (Wu et al., 2008). The extent of such pest suppression is directly tied to the spatial cover of Bt crops, with their increased prevalence in farming landscapes resulting in lowered pest pressure (Carrière et al., 2003; Gao and Feng, 2010). In principle, a year‐by‐year reduction in the spatial cover of insecticidal cultivars is expected to lead to a resurgence of target pest populations and a rebound of farmers’ pesticide use. Yet, these important social‐ecological features of Bt crop cultivation so far have not been examined.

In North China Plain, surging labour costs, lower farmer profits, and less favourable cotton prices recently led to an 80% contraction of the cotton planting area as compared to the early 2000s (Qiao and Yao, 2015). This provides an ideal case to investigate an eventual pest population resurgence due to Bt crop dis‐adoption.

Results

First, we examined the extent to which Bt cotton area affected H. armigera adult abundance (through pheromone trapping) in 29 different 2‐km radius landscape sectors over 2016–2017 (Figure 1; Tables S1 and S2). Generalized linear mixed model (GLMM) analysis showed that a declining share of Bt cotton among all host crops is the main driver for a heightened H. armigera abundance (effect size = −2.44, P < 0.001, Figure 2; Tables S3 and S4), whereas other land‐use variables, for instance, the percentage of maize (P = 0.929), woodlots (P = 0.146), and water cover (P = 0.909), had no significant effect on adult abundance (Tables S3 and S4) at the landscape level.

Figure 1

Geographical location of sampling sites in North China Plain. Helicoverpa armigera populations were surveyed through light‐trapping and field sampling in 50 different counties. Statistics of pesticide use and yield loss due to H. armigera during 2007–2019 were obtained from Henan, Hebei, and Shandong provinces. The pink triangle indicates the location of the 2016–2017 landscape‐level study.

Figure 2

Adult abundance of Helicoverpa armigera is associated with the relative spatial cover of Bt cotton within the agro‐landscape. Ambient population levels increased with a declining share of Bt cotton area. Adult abundance is log10‐transformed. Regression line with 95% confidence intervals (shaded areas) is shown. Data points represent the mean observed values for each trapping round in 29 sites.

Second, we assessed the variation in H. armigera adult abundance (through light‐trapping) in 29 cotton‐growing counties over time and recorded the respective larval infestation pressure in fields of three non‐Bt crops: maize, peanut, and soybean. Infestation sampling in the above non‐Bt crops was done in a respective 39, 30 and 19 counties of three provinces (Hebei, Henan and Shandong) over 2007–2019 (Figure 1; Tables S1 and S5). While the share of Bt cultivars within the total cotton area attained ~100% by 2011, its relative spatial cover among all H. armigera host crops decreased from 10% (2007) to less than 2% (2019) (Figure 3a). Meanwhile, the relative share of maize crops in agro‐landscape sectors increased from 41% (2007) to 53% (2019), whereas other non‐cotton host crops experienced no notable variation in cultivated area over time (Figure S1). This gradual decline in the Bt cotton spatial cover was associated with a 1.9‐fold increase in H. armigera adult abundance by 2019 (P = 0.001, Figure 3b and Table S6). Further, the relative cover of Bt cotton among all host crops negatively affected H. armigera abundance (Figure 3c and Table S6). Drawing upon a fixed‐effect panel data model (Zhang et al., 2018) (Table S6), we compellingly attribute increases in H. armigera abundance to the relative extent of Bt cotton cropping within local agro‐landscapes.

Figure 3

Temporal patterns in Bt cotton adoption, Helicoverpa armigera abundance, and infestation pressure in three non‐cotton host crops. Over time, the relative share of Bt cotton increased among all cultivated cotton cultivars and decreased among all host crops (a); H. armigera adult abundance increased from 2007 to 2019 (b); and H. armigera abundance related to a gradual declining share of Bt cotton area over time (c). Adult abundance records are log10‐transformed. A respective decrease or increase in H. armigera infestation level in three non‐cotton host crops (d); and larval infestation levels in the three non‐cotton crops increased with a declining share of Bt cotton area (e). Regression lines with 95% confidence intervals (shaded areas) are shown.

When assessed through a (1–5) qualitative score with 5 being extreme infestation, average infestation level of H. armigera larvae increased in maize (2.2–2.9), peanut (1.8–2.7) and soybean (1.6–2.4) from 2007 to 2019 (Figure 3d). Fixed‐effect panel data models further accentuate how, with a declining share of Bt cotton among all host crops, infestation levels increased significantly in all three crops (Figure 3e; Table S7). Further, larval infestation pressure on different crops was positively related to ambient adult abundance (maize: P = 0.034; peanut: P = 0.032; soybean: P = 0.021) (Figure S2; Table S8). Overall, regional H. armigera adult abundance and larval infestation pressure increased during 2007–2019 while the share of Bt cotton among host crops declined.

Third, we investigated chemical pesticide use against H. armigera and pest‐induced yield loss in non‐cotton crops in three provinces over 2007–2019 (Table S1). A decreasing share of Bt cotton among all host crops was mirrored by a 2.0‐, 2.0‐, and 4.4‐fold increase in insecticide application frequency for H. armigera in maize, peanut, and soybean, respectively (Figure 4a). The extent of Bt cotton cultivation (among all crops) was negatively related to the number of synthetic insecticide applications against H. armigera for maize (P < 0.001), peanut (P < 0.001), and soybean (P = 0.003) (Figure 4b; Table S9). Similarly, a gradual decline in Bt cotton cultivation may have aggravated H. armigera‐induced yield loss. From 2007 to 2019, absolute yield loss increased 1.5‐ to 2.l‐fold in maize (from 35.4 to 74.4 kg/ha), peanut (43.0–73.1 kg/ha), and soybean (38.5–57.2 kg/ha) (Figure 4c). Yield loss was negatively related to the share of Bt cotton among all host crops for maize (P < 0.001), peanut (P = 0.008), and soybean (P = 0.003) (Figure 4d; Table S10). The yield loss ratio (i.e. yield loss divided by average yield per ha) increased for maize (from 0.7% to 1.2%), peanut (1.1%–1.7%) and soybean (2.1%–2.2%) (Figure S3a). Regional cover of Bt cotton negatively related to this yield loss ratio for maize (P < 0.001), peanut (P = 0.020), and soybean (P = 0.012) (Figure S3b; Table S10). Hence, over a 13‐year time period, a contraction of Bt cotton area exacerbated pest‐induced yield loss and insecticide use in three prime food crops, thus entailing important socio‐economic and ecological impacts.

Figure 4

Chemical pest control and pest‐induced yield loss in maize, peanut, and soybean are associated with Bt cotton adoption. Temporal patterns in insecticide use over 2007–2019 (a); the number of insecticide applications in maize, peanut, or soybean is negatively related to the share of Bt cotton among all host crops (b); temporal patterns in yield loss (kg/ha) over 2007–2019 (c); declining yield loss (kg/ha) in maize, peanut and soybean with an increased share of Bt cotton among all host crops (d). Regression lines with 95% confidence intervals (shaded areas) are shown.

Discussion

In North China Plain, H. armigera completes four generations—with its 1st generation larvae feeding primarily on wheat in the spring. During summer (2nd, 3rd, and 4th generation), H. armigera preferably ovipositors on cotton plants and its larvae attain high fitness on non‐Bt cotton. While maize is a high‐quality host only for 4th generation bollworms, other crops (e.g. peanut, soybean, vegetables) are deemed to be optimum hosts for its 3rd generation. Thus, second generation H. armigera feed mainly on cotton while third and fourth generation larvae feed on a variety of crop species including maize, peanut, soybean, vegetables, and cotton (Guo, 1998). Since its commercial introduction in 1997, Bt cotton hereby effectively acted as a dead‐end trap crop for H. armigera and suppressed its population levels in multiple crops (Wu et al., 2008), thereby reducing pesticide use and bolstering biological control (Lu et al., 2010, 2012; Zhang et al., 2018). In this study, we demonstrate how a gradual dis‐adoption of Bt cotton cultivars can release this major agricultural pest from Bt‐mediated control and inflate pest‐induced crop losses.

An elevated H. armigera infestation pressure can inflict both direct and indirect negative effects on environmental health, food security, and societal well‐being. First and foremost, the increase in synthetic pesticide use can negatively impact upon biodiversity, beneficial arthropods, associated ecosystem services and the broader environment (Dainese et al., 2019; Desneux and Decourtye, 2007). The associated human health impacts are also substantial, involving direct exposure of farm workers and dietary intake of pesticide‐tainted produce (Boedeker et al., 2020; Wyckhuys et al., 2020). Feeding damage by lepidopteran pests such as H. armigera equally can increase mycotoxin load in harvested maize kernels, posing hazards for food safety and animal health (Li et al., 2015).

However, a sole reliance upon Bt technologies might carry potential risks; safeguards need to be put in place to circumvent resistance development and avert secondary pest outbreaks (Tabashnik, 2016; Zeilinger and Olson, 2016). So far, the generalist feeding habits of H. armigera and ample presence of non‐transgenic crops (as natural refuges) have counter‐acted Bt resistance development in North China Plain. While Bt cotton still effectively suppresses H. armigera (Zhang et al., 2019), Bt resistant genotypes have made their local appearance (Jin et al., 2015) and sap‐feeding pests have undergone population resurgence in Bt cotton landscapes (Lu et al., 2010; Zhang et al., 2018). Those potential risks can be resolved through a tactical (i.e. scientifically‐guided) deployment of Bt cotton cultivars, coupled with a full repertoire of agro‐ecological practices and a drastic reduction in pesticide use (e.g. Bommarco and Kleijn, 2013). By, thus, combining several (non‐chemical) management strategies, pest prevention can be maximized, agrochemical pollution can be curbed, and non‐negligible social‐ecological benefits can be reaped (Gurr et al., 2017; Lewis et al., 1997; Veres et al., 2020).

For polyphagous herbivores, the presence of one or few suitable host crops within a landscape mosaic can have—often dramatic—effects on their overall occurrence and population dynamics (Meisner and Zaviezo, 2017; Rand et al., 2014). Specifically for H. armigera, the impacts of (Bt) cotton are magnified due to its season‐long suitability and superior host‐plant quality as compared to, for example maize, peanut, or soybean (Guo, 1998). These attributes make (Bt) cotton a pivotal element in H. armigera area‐wide management; among all land‐cover categories (including maize and non‐crop habitat), Bt cotton exhibits the most pronounced effects on abundance levels of the polyphagous H. armigera (Table S4). These results continue to hold when extending our assessment to a broader spatio‐temporal scale. The latter analyses (i.e. Bt cotton contraction relates to H. armigera build‐up) equally illuminate how relative Bt cropping area explains pest‐induced crop losses and associated pesticide usage levels. By 2001, the percentual share of Bt cotton varieties (total cotton area) had surpassed 70% and Bt cotton occupied 7%–11% of all host crops in North China Plain, thus resulting in an effective area‐wide control of H. armigera (Wu et al., 2008). Though most cotton fields continued to be planted with Bt varieties over the past decade, Bt cotton spatial cover dropped from 10% to 2% (Figure 3a); this being the main driver of pest population resurgence in non‐cotton crops such as maize, peanut, and soybean. Bollworm populations rebounded in these crops, though remain below those of the early 1990s (Guo, 1998). Long‐term H. armigera suppression potentially can be attained once Bt crop cultivation surpasses a certain threshold with non‐Bt cotton serving as an effective refuge, as evidenced in Arizona (USA; Carrière et al., 2003). Area‐wide pest suppression is further shaped by H. armigera dietary breadth, seasonal host use, (long‐range) dispersal and the extent of cultivation, management intensity and nutritional quality of alternative (non‐Bt) host crops. Hence, a full integration of ecological knowledge, landscape‐level dynamics, and farming‐system perspectives is critical to design robust and environmentally sustainable plant health strategies (Lewis et al., 1997).

The establishment of Bt crops can generate clear benefits beyond pest suppression (Dively et al., 2018; Hutchison et al., 2010), though their landscape‐level dis‐adoption can compromise agricultural sustainability and exacerbate ecological and health risks (e.g. as tied to pesticide use). Given these profound social‐ecological implications of technological change or shifting land‐use patterns, we conclude that a judicious deployment of Bt transgenics within (a fraction of) the agro‐landscape can have long‐lasting and desirable effects on sustainable food production and societal well‐being.

Experimental procedures

Our analysis consists of three components (Table S1). First, we monitored adult H. armigera populations in agricultural landscapes with varying levels of Bt cotton cultivation, in order to relate (ambient) H. armigera abundance and regional Bt cotton adoption. Second, we used area‐wide monitoring to examine long‐term trends in H. armigera adult abundance and (in‐field) larval infestation on three non‐cotton host crops, and assessed the eventual effects of temporal shifts in Bt cotton adoption. Third, we used province‐level production statistics to gauge the socio‐economic and ecological repercussions (i.e. insecticide application frequency, yield loss) of increasing H. armigera infestation levels on the above crops. In earlier work, we reported how Bt cotton cultivation achieved an area‐wide suppression of H. armigera in multiple crop systems of North China Plain over 1997–2006 (Wu et al., 2008). Here, we present findings from the same geographical localities (i.e. North China Plain) and agro‐ecosystems over 2007–2019.

Data gathering & observational assays

Population monitoring at landscape level

Insect population surveys were conducted in a total of 29 different commercial cotton (100% Bt cotton) fields (Figure 1; Table S2), covering 15 and 14 fields in 2016 and 2017, respectively. Fields were selected along a gradient of increasing cotton area within the surrounding (2000 m radius) landscape sector—that is, a suitable spatial scale to study the response of pestiferous insects to landscape features (Karp et al., 2018; Tscharntke et al., 2012). Geo‐referenced individual fields were 7.02 ± 6.44 ha (mean ± SD) in size and spaced at least 4000 m apart within a given year. Land uses were mapped and assessed for each landscape sector using open‐access satellite imagery through Google Earth. Land‐use maps were validated through extensive ground‐truthing in each study year. Within each landscape sector, land‐use cover was classified into following categories: Bt cotton (i.e. focal crop), maize, other small crops, woodlots (i.e. forest, shrubland, and ruderal vegetation), water‐bodies, roads, and urban dwellings. For each land‐use category, the percent spatial cover was determined using ArcGIS 10.5 (ESRI, Redlands, CA).

Pheromone trapping was used to monitor the abundance of H. armigera adults. In each field, we deployed a total of five bucket traps (Cuerrero et al., 2014) baited with a commercial sex pheromone lure (Pherobio Technology Co., Ltd, Beijing, China). Trapping was conducted during July–August in each year, thus covering the 2nd and 3rd generations of local H. armigera populations. Traps were randomly positioned in each study field, installed at approximately 15 cm above the foliage and spaced at a minimum distance of 20 m. In each year, six successive week‐long trapping rounds were carried out from July 6–8 to August 24. At the end of each trapping round, lures in each bucket trap were replaced, the number of field‐caught H. armigera adults was recorded and traps were emptied. For a given field and trapping round, the mean number of adults captured across traps represented the in‐field H. armigera population abundance.

Multi‐year population monitoring at regional level

From 2007 to 2019, H. armigera populations were monitored through light‐trapping at 29 different sites within the main cotton‐growing counties of Henan, Hebei and Shandong provinces (Figure 1; Table S5). At each site, a 20 W ultraviolet light trap was installed at a height of 1.5 m above ground, and light‐trapping was performed every night under the National Standards of the People’s Republic of China (GB/T15800‐1995). Trapping was conducted from sunset through sunrise, except during rain, extreme weather or power outage, and carried out from June 1 to August 31 each year. Each day, the total number of field‐caught H. armigera adults was recorded. Light‐trapping was carried out over the span of 13 years, thus yielding unique insights into H. armigera population dynamics during a time period in which relative cover of Bt cotton in the agro‐landscapes gradually declined.

Concurrently, field‐level infestation levels of H. armigera larvae were assessed in maize, peanut and soybean in a respective 39, 30 and 19 different counties. Field surveys were conducted through the national pest monitoring network of the Ministry of Agriculture and Rural Affairs (MARA) (Figure 1; Table S5). Field surveys were performed in one or more crops per county, based upon their respective cropping area and economic importance. For each county, 10–20 fields were also sampled for each crop based on the standardized protocols of the Agricultural Technology Extension Service Center, China (Liu and Wang, 1996; Zhang et al., 2018). Infestation levels were recorded using a categorical score ranging from 1 (no infestation) to 5 (extreme infestation) (Qu, 1992; Zhang et al., 2018).

Annual cropping areas of the different H. armigera host crops and Bt cotton were either extracted from the China Agriculture Yearbook (MoA, 2008–2019) or provided by the Agricultural Technology Extension Service Center, MARA. On an annual basis, the relative share of Bt transgenic cultivars in the local cotton crop was calculated by dividing the area of Bt cotton by the total acreage of cotton crop. The share of Bt cotton among all key host crops was calculated by dividing the acreage of Bt cotton by the total acreage of H. armigera host crops.

Chemical control and yield loss at regional level

Aside from examining the effect of Bt cotton adoption on H. armigera ambient population levels and in‐field pest pressure, we assessed its effects on chemical pesticide use and (pest‐induced) yield loss.

Chemical control (i.e. cumulative area receiving synthetic pesticide applications in all cultivated areas) for H. armigera on maize, peanut, and soybean were assessed by the standardized proposals (Liu et al., 1996) and provided for Hebei, Henan and Shandong province (2007–2019) by the Agricultural Technology Extension Service Center, MARA (Figure 1). For each crop, the ‘cumulative area under chemical control’ represents the total (i.e. summed) area receiving single pesticide applications targeting H. armigera throughout the growing season. In a given crop, pesticide use intensity (i.e. number of pesticide applications) against H. armigera can be deduced by dividing the cumulative area under chemical control by the total cultivated area per crop.

Yield loss (i.e. total yield loss in all cultivated area) records for H. armigera on maize, peanut and soybean in Hebei, Henan and Shandong province (2007–2019) were provided by the Agricultural Technology Extension Service Center, MARA (Figure 1). For a given crop, H. armigera‐inflicted yield loss was annually estimated by comparing the yield levels of fields where H. armigera infestation was left uncontrolled versus pesticide‐treated fields (both of them had the similar irrigation and fertilization conditions) in a paired design according to standardized protocols of the Agricultural Technology Extension Service Center, China (Liu et al., 1996). For each crop, two indices were calculated: yield loss (kg/ha) and yield loss ratio (%). More specifically, Yield_loss (kg/ha) = total yield loss for all cultivated area/total cultivated area, and Yield_loss_ratio (%) = Yield_loss/average crop yield per ha.

Data analysis

Generalized linear mixed models were used to assess the effects of landscape variables (explanatory factors) on H. armigera adult abundance. The percent spatial cover of Bt cotton among the main H. armigera host crops (i.e. Bt cotton%/[Bt cotton% + other host crops %]) was the focal fixed effect, while individual land‐cover categories, that is, maize, woodlot, and water were also included as fixed effects. The proportion of road and urban dwellings were omitted in order to lower the degree of collinearity. These variables also provide minor or no nutritional resources to H. armigera and thus are of limited biological relevance. Study site was considered as a random effect. As a response variable, we used H. armigera adult abundance—that is, the average adult abundance across all sampling rounds for a given focal field in each site. GLMM analysis was performed with the negative binormal error distribution as this provided the best fit to the abundance data. Prior to analysis, adult abundance was log10‐transformed, we used the R package ‘spdep’ to conduct Moran’s I statistics to detect eventual spatial autocorrelation among adult abundance from different sites (Bivand and Wong, 2018), and no spatial autocorrelation were found here.

A multi‐model inference procedure was used to reduce model selection bias and obtain robust results. Model selection was done using bias‐corrected Akaike’s information criterion (AICc) to account for small sample size (Burnham and Anderson, 2004), and all models were selected with ΔAIC < 4 (Grueber et al., 2011). The model averaging procedure revealed the relative importance of explanatory variables, as quantified by the sum of their respective Akaike weights among the selected set of models (Grueber et al., 2011). Furthermore, marginal R 2 (i.e. variation explained by fixed factors) and conditional R 2 (i.e. variation explained by both fixed and random terms) were computed for each model (Nakagawa and Schielzeth, 2013).

We used the 'lme4' package for GLMM analysis (Bates et al., 2015), ‘r.squaredGLMM’ and ‘model.avg’ function in 'MuMIn' package for calculating R 2 and multi‐model inference (Kamil, 2020). All analyses were executed by using R version 4.0.2 software (R Development Core Team, 2020).

An econometric method of fixed‐effect panel data modelling (Larsen, 2013; Wooldridge, 2010; Zhang et al., 2018) was used to assess whether the share of Bt cotton among all host crops affects adult H. armigera abundance (i.e. as obtained through light‐trapping) over the 13‐year period (Table S6). More specifically, models were composed of dummy variables for each site (cross‐sectional unit) or period, with inclusion of Bt cotton adoption among all host crops (called Bt_ratio in the model) as an explanatory variable. We excluded other host crop variables (such as maize) from the fixed effects due to its strong collinearity and correlation with the focal crop‐ Bt cotton (i.e. variance inflation factor VIF = 10.85, r = −0.95). Dummy variables were included for year and data collection site. H. armigera abundance was expressed as the number of adults caught per month (i.e. 30 nights) for each site and year.

The fixed‐effect panel data approach helps capture the effects of explanatory variables using year‐to‐year variation within sites, controlling for unobserved heterogeneity between sites and years (Zhang et al., 2018). As such, the potential for confounding due to variable omission is minimized because time‐invariant characteristics drop out of the model (Larsen, 2013). This approach accounts for site‐specific heterogeneity that is difficult to observe, with the ensuing findings being more robust. Furthermore, the inclusion of a ‘year’ dummy variable accounts for time‐dependent effects shared by all sites such as national pesticide regulation or weather anomalies (Larsen, 2013). To test the robustness of model estimations, we built the fixed‐effect panel data model using three estimation methods: (i) linear regression using ordinary least squares (OLS) (column 1 in Table S6), (ii) tobit regression (column 2 in Table S6), and (iii) generalized linear model (GLM) regression with gamma error distribution and log link, which is suited to dealing with heteroskedasticity in non‐negative data (column 3 in Table S6), all conducted in Stata 15 (StataCorp).

Following Zhang et al. (2018), we note that repeated observations over time may have correlated disturbance terms, which could result in artificially small standard errors (SEs) (MacArthur, 1965). Hence, when estimating the OLS and GLM fixed‐effect regression models, we used cluster‐robust SEs clustered at the site to allow for arbitrary autocorrelation between observations within the same site. For the nonlinear tobit regression, nonparametric bootstrap estimation was used to obtain the correct SEs, taking into account clustering.

This same approach was also used to assess whether the relative area under Bt cotton affects (field‐level) H. armigera infestation pressure in maize, peanut, and soybean across multiple counties of the three provinces. Since the response variable (i.e. larval infestation level) is measured with ordered categorical data, nonlinear ordered probit regression is appropriate, with an additional OLS regression as a robustness check (Table S7). The analyses were performed for each of three non‐cotton host crops (maize, peanut, and soybean), respectively.

For each dependent variable (i.e. the number of pesticide application, absolute yield loss and yield loss ratio), data were shown per province in each year. Over 2007–2010, no data were available for Shandong Province on chemical control and yield loss in maize or soybean. Yet, as residuals for the other two provinces met normality distributions, absence of these records likely did not impact the robustness of statistical results. Linear regression analysis was performed to relate the different response variables to Bt cotton adoption among all host crops (Tables S8–S10). Yield loss data were log10‐transformed prior to analysis. Linear regression analyses were conducted using the ‘car’ package (John and Sanford, 2011) within R 4.0.2 software (R Development Core Team, 2020).

Conflict of interest

The authors declare no conflict interests.

Author contributions

Y.‐H.L. and K.‐M.W.: conceived the idea. Y.‐H.L., L.Y., J.Z., and Y.‐Y.J.: performed the research. Y.‐H.L., K.A.G.W., L.Y., B.L., and W.Z.: analysed the data. Y.‐H.L., K.A.G.W., L.Y., B.L., N.D., and W.Z.: wrote the draft. All authors contributed to the revision of the manuscript.

Figure S1. Yearly changes in the relative shares of cotton, maize, peanut, soybean, and others crops among all Helicoverpa armigera host crops in North China Plain.

Figure S2. Larval infestation pressure of Helicoverpa armigera in three non‐cotton host crops relates to its ambient adult abundance.

Figure S3. Pest‐induced yield loss in maize, peanut, and soybean is associated with Bt cotton adoption.

Table S1. Summary table of all data in this study.

Table S2. Geographical location records for all sampling sites included in the landscape‐level study in 2016 (sites 1–15) and 2017 (sites 16–29).

Table S3. Summary of candidate models (ΔAICc < 4) for multi‐model inference to assess the relative effects of Bt cotton and other land‐use variables on Helicoverpa armigera abundance in Bt cotton fields.

Table S4. Model‐averaged inference of the GLMM (based on ΔAICc <4, Table S3) to assess the relative effects of Bt cotton ratio among all host crops and other landscape variables (i.e. proportion cover of maize, woodlots, and water‐bodies within a 2‐km radius) on Helicoverpa armigera abundance in Bt cotton fields.

Table S5. Geographical records for all sampling sites included in the 2007–2019 region‐level monitoring of Helicoverpa armigera adult abundance and larval infestation pressure.

Table S6. Fixed‐effect (FE) panel data model analysis of the effect of Bt cotton cover (i.e. share of Bt cotton among all host crops) on Helicoverpa armigera adult abundance over a 13‐year period.

Table S7. Fixed‐effect (FE) panel data model analysis of the effect of Bt cotton cover (i.e. relative area under Bt cotton) on Helicoverpa armigera larval infestation level in three main non‐cotton host crops across multiple sites over a 13‐year period.

Table S8. Linear regression analysis to assess the relationship between Helicoverpa armigera infestation level in three non‐cotton crops and adult abundance (as obtained through light‐trapping) over a 13‐year period.

Table S9. Linear regression analysis to assess whether chemical control frequency on three non‐cotton crops is affected by Bt cotton cover (i.e. proportion of Bt cotton adoption among all host crops) over a 13‐year period.

Table S10. Linear regression analysis to assess whether Helicoverpa armigera‐inflicted yield loss in three non‐cotton host crops is affected by Bt cotton cover (i.e. proportion of Bt cotton area among all host crops) over a 13‐year period.

Click here for additional data file.