Multidecadal, county-level analysis of the effects of land use, Bt cotton, and weather on cotton pests in China.

Long-term changes in land use, climate, and agricultural technologies may affect pest severity and management. The influences of these major drivers can only be identified by analyzing long-term data. This study examines panel data on land use, adoption of genetically modified Bacillus thuringiensis (Bt) insect-resistant cotton, weather, pest severity, and insecticide use on three major cotton pests for 51 counties in China during 1991-2015. Bt cotton had pervasive effects on the whole pest complex in cotton and its management. Adoption resulted in major reductions in insecticide use for bollworm control. The resulting restoration of aphid biological control decreased aphid severity. However, mirid bugs, which have few effective natural enemies in cotton, increased in severity with warming May and reduced insecticide spraying against bollworm. The effects of landscape on pest severity were pest specific. The severity of cotton aphid and mirid bugs decreased with higher land use diversity, but the severity of highly polyphagous cotton bollworm was unrelated to land use diversity. Shares of forest, water body, and unused land area were negatively associated with the severity of mirid bugs, whereas cotton bollworm responded positively to the shares of water body and unused land area. Farmers sprayed insecticides at mild infestation levels and responded aggressively to severe bollworm outbreaks. Findings support the usefulness of Bt-based plant resistance as a component of integrated pest management (IPM) but highlight the potential for unexpected outcomes resulting from agro-ecosystem feedback loops as well as the importance of climate.

Despite tremendous improvements in breeding and other technologies for robust yield increase, crop pests remain an important cause of considerable yield loss, triggering the use of insecticides that affect farm profit and the health of humans and their environment. Cotton is the most heavily treated agricultural commodity, accounting for one-third of all pesticide use globally, and the widespread use of highly toxic insecticides in cotton poses significant threats to human health and the environment (1).

Complex landscapes generally provide improved pest suppression ecosystem services that support crop production and reduce the need for external inputs compared with simplified landscapes. Simplified landscapes tend to have fewer plant and invertebrate species, larger fields, and less noncrop habitat where natural enemies of crop pests may find important resources such as nectar, pollen, other alternative foods, and shelter (2–4). While natural habitat in agricultural landscapes has been shown to increase pest control in many systems, the net effect of landscape complexity on pest severity has been mixed, context dependent, and system specific (2–4). The relative importance of natural habitat for biological pest control can vary dramatically, depending on type of crop, pest, predator, land management, and landscape structure (5). This needs to be considered when designing measures aimed at enhancing biocontrol services through restoring or maintaining natural habitat (5).

Landscape-based ecosystem services interact with crop management (6). This interaction is characterized by feedback: If the landscape-based ecosystem service of biocontrol is effective at crop level, a farmer may refrain from using pesticides, and vice versa. However, if pesticide use in crops is indiscriminate, it should be expected that natural enemies may not be effective, even if they are present in a landscape, while the biocontrol service may be destroyed in the long term due to the interruption of natural enemies’ life cycles. Because of this linkage, farmers can develop a “lock-in” syndrome where continued heavy spraying is necessary to compensate for the lack of enemies that this same spraying has caused, a syndrome described as a “pesticide treadmill” (7). Restoration of natural enemy communities by broadscale adoption of judicious pesticide use can revert this syndrome (7, 8). The introduction of genetically engineered insect resistance in crop opens an opportunity to escape from this treadmill and integrate agroecology with biotechnology (9). The impacts of crop management versus landscape-based ecosystem services on pest severity have been little explored in empirical studies because they require large-scale intensive data gathering on land use and crop management.

Weather is an important factor affecting insect population and movement. Positive physiological responses to increasing temperatures allow for faster insect population growth and facilitate movement, and milder winters allow for earlier commencement of colonization of crops and a reduction in winter mortality (10–14). Most analyses show that, in a warmer climate, pests may become more abundant and may expand their geographical range (15). Precipitation also affects crop–pest interactions (15). Both direct and indirect effects of moisture stress on crops make them more vulnerable to damage by pests, especially in the early stages of growth (15). Therefore, improved understanding of climate-induced risks on pest severity over extended periods is essential for the discussion of climate adaptation policies related to agriculture. China’s climate has warmed since the 1960s, with stronger warming in the north and increased rainfall contrast between northeastern and southern China (16–18). Weather effects on agricultural pests need to be assessed in the light of climate change, which is projected to bring significant warming to large parts of China over the coming decades (19, 20).

China is the largest cotton producer in the world (21). It is also one of the largest pesticide consumers worldwide. For instance, China used 4.4 times more pesticides (in tons of active ingredients) than the United States in 2012 (22). An estimated 30–40% of all pesticides applied in China are used on cotton, making it the most heavily treated agricultural crop (1). Nearly 40% of the pesticides used by Chinese cotton farmers contain active ingredients that are classified as extremely or highly hazardous by the World Health Organization (23), contributing to around 400–500 cotton farmer deaths every year from pesticide poisoning (24). Genetically modified Bacillus thuringiensis (Bt) insect-resistant cotton was introduced in China’s agriculture in 1997, and adoption progressed at different rates in different regions depending on the timing of Bt cotton varieties approved for commercialization and the availability of Bt cotton seed in local markets (25). The practice of applying excessive amounts of highly toxic pesticides has continued even after the adoption of Bt cotton (26–28). Research has shown that microlevel behavioral factors such as risk aversion (23) and lack of knowledge by Chinese cotton farmers (29) are important factors driving pesticide use, while others have suggested that Bt cotton seed quality (23, 30) and secondary pests (31) are also at play.

The effects of insecticide use on both human health and ecological systems need to be mitigated. Different agricultural, land use, conservation, and climate adaptation policy measures may address different drivers of pest severity and insecticide use (32). For effective policy development, it is necessary to understand how abiotic or biotic factors jointly and quantitatively drive pest severity and pesticide use at broad geographic scales at which experimentation is not feasible. Using a unique long-term panel (i.e., longitudinal) dataset for 51 counties in eight provinces (Fig. 1), this study applies rigorous econometric methods to examine the main drivers of cotton pest severity and insecticide use at the county scale in China over a 25-y period, from 1991 to 2015. The eight provinces included in the study (Anhui, Hebei, Henan, Hubei, Jiangsu, Shaanxi, Shandong, and Shanxi) cover two of China’s three major cotton production regions, the Yangtze River valley and Yellow River valley regions, and accounted for one-half of the national cotton production in 2010 (33).

Fig. 1.

Locations of counties included in the study and land use diversity (Shannon diversity) for 1990. Data used to make the map are reported in .

We consider three major cotton pests: cotton aphid (Aphis gossypii), mirid bugs (mainly including Apolygus lucorum, Adelphocoris suturalis, Adelphocoris lineolatus, and Adelphocoris fasciaticollis), and cotton bollworm (Helicoverpa armigera). Together, these three pests accounted for 81% and 75% of cotton crop loss caused by arthropod pests in the presence of insect control in 1994 and 2001, respectively (34). We conclude with a discussion of policy implications on sustainable solutions to pest management.

Results

Descriptive Analysis.

Cotton acreage in the sampled counties varied over time, with a steady decline from 1991 to 2001, a restoration from 2001 until 2004, and a slowly decreasing trend afterward (Fig. 2). The adoption rate of Bt cotton varied considerably across provinces, with Anhui, Hubei, and Jiangsu in the more southern Yangtze River valley region introducing Bt cotton at a slower pace than the provinces in the Yellow River valley region (Hebei, Henan, Shaanxi, Shandong, and Shanxi). Six of the eight provinces had reached almost full adoption by 2008, while all eight provinces had full adoption by 2015 (see Bt cotton adoption rate by province in ).

Fig. 2.

Proportion of cotton area in total cultivated area and proportion of Bt cotton area in the total cotton area across counties in the sample, 1991–2015.

Pest severity was measured using a five-point scale: (i) level I, no pest or very low level of infestation; (ii) level II, slight; (iii) level III, moderate; (iv) level IV, severe; and (v) level V, extremely severe infestation. Across the sampled counties, cotton aphid severity showed considerable decline over the early 1990s and remained largely constant from1995 until 2010, and declined slowly thereafter (Fig. 3). Cotton bollworm severity declined remarkably from an average rating of 4.4 during 1991–1996 to 2.2 during 2005–2010, and continued to decrease afterward (Fig. 3). In contrast, mirid bugs severity was low until around 1997 before it started to increase drastically, from an average level of 1.1 during 1991–1996 to 3.1 during 2005–2010, after which it slightly decreased (Fig. 3).

Fig. 3.

Cotton pest severity (A) and number of insecticide applications targeted at cotton aphids, mirid bugs, and cotton bollworm (B) across all counties in the sample, 1991–2015.

Trends in insecticide use largely corresponded to those of infestation levels of the respective pests at which the insecticides were targeted (Fig. 3). Severe bollworm outbreaks were met with frequent sprays in the early 1990s, but overall treatment declined from an average of 9.3 sprays during 1991–1996 to less than 3 sprays per year after 2006. Largely untreated before 1997, mirid bugs were on average sprayed as many as 3.7 times during 2005–2015. Overall, the severity of pests and frequency of insecticide sprays showed a slightly decreasing trend since 2008 for all three pests.

Cultivated land was the dominant land use in the surveyed counties throughout the period (1990–2015), accounting for over 55% of the total land area, followed by built-up land (26%), water body (10%), forest (5%), grassland (4%), and unused land (1%). Cultivated land, grassland, and unused land declined in area proportions between 1990 and 2015, whereas the area shares of built-up land, water, and forest increased (). Land use diversity changed marginally during the period (), but there was a large variation across counties (Fig. 1).

Monthly average temperatures in May, June, and August increased during 1991–2015 (Fig. 4). Monthly total precipitation decreased in June but increased in August during the same period (Fig. 4). The months of May and June became significantly warmer over time, while June also became dryer (see estimated time trends from regression analysis in ), and these changes are likely to affect cotton production adversely, by stimulating pests (15).

Fig. 4.

Monthly average temperature (in degrees Celsius) (A) and monthly total precipitation (in millimeters) (B) for the months of May, June, July, and August 1991–2015. Time trends were smoothed using local polynomial smoothers.

Econometric Analysis Results.

Drivers of pest severity.

Land use diversity, as measured by Shannon diversity (, ), negatively affected cotton aphid and mirid bug severity (Table 1; full regression results are reported in ). For cotton aphid, pest severity tended to be higher when there was a greater share of crop land devoted to cotton, ceteris paribus. In addition, the proportion of built-up area was negatively correlated with cotton aphid severity, relative to the proportion of cultivated land. For mirid bugs, the proportions of forest, water body, and unused land area were negatively correlated with pest severity.

Table 1.

Estimated coefficients for the severity of cotton aphid, mirid bugs, and cotton bollworm

	Cotton aphid		Mirid bugs		Cotton bollworm
Explanatory variables	(1)†	(2)†,‡	(1)†	(2)†,‡	(1)†	(2)†,‡
Land use diversity (Shannon diversity)	−2.486**		−2.071*		−0.705
	(1.129)		(1.213)		(1.208)
Share of cotton area in total cultivated area, %	0.008*	0.011**	0.001	0.002	−0.008	−0.006
	(0.005)	(0.005)	(0.005)	(0.005)	(0.006)	(0.006)
Share of Bt cotton area in total cotton area, %	0.003	0.003	0.001	0.001	−0.012***	−0.012***
	(0.002)	(0.002)	(0.002)	(0.003)	(0.002)	(0.003)
No. of insecticide sprays against cotton bollworm	0.068***	0.063***	−0.044*	−0.043*
	(0.021)	(0.017)	(0.025)	(0.024)
Percentage of forest area, %		0.018		−0.763*		0.035
		(0.380)		(0.461)		(0.323)
Percentage of grassland area, %		0.030		−0.152		0.097
		(0.125)		(0.152)		(0.148)
Percentage of built-up area, %		−0.071*		0.003		−0.022
		(0.039)		(0.030)		(0.026)
Percentage of water body area, %		−0.066		−0.435***		0.312***
		(0.141)		(0.159)		(0.113)
Percentage of unused land area, %		−0.029		−0.274**		0.135*
		(0.098)		(0.109)		(0.073)
Monthly mean temperature: May	0.058	0.073	0.154***	0.159**	−0.014	−0.007
	(0.046)	(0.048)	(0.052)	(0.066)	(0.055)	(0.053)
Monthly mean temperature: June	0.005	0.008	0.049	0.052	−0.212***	−0.201***
	(0.048)	(0.054)	(0.066)	(0.068)	(0.062)	(0.069)
Monthly mean temperature: July	0.129*	0.135**	−0.088	−0.087	0.027	0.033
	(0.066)	(0.067)	(0.066)	(0.061)	(0.069)	(0.073)
Monthly mean temperature: August	−0.145**	−0.156**	−0.066	−0.067	−0.012	−0.023
	(0.067)	(0.075)	(0.079)	(0.063)	(0.069)	(0.059)
Monthly total precipitation: May	0.002*	0.002*	−0.000	−0.000	0.001	0.001
	(0.001)	(0.001)	(0.001)	(0.001)	(0.001)	(0.001)
Monthly total precipitation: June	0.000	0.000	0.000	0.000	−0.001*	−0.001
	(0.001)	(0.001)	(0.001)	(0.001)	(0.001)	(0.001)
Monthly total precipitation: July	0.001	0.001	0.000	0.000	0.000	0.000
	(0.000)	(0.000)	(0.001)	(0.000)	(0.000)	(0.000)
Monthly total precipitation: August	−0.000	−0.000	−0.001**	−0.001*	−0.000	−0.000
	(0.000)	(0.000)	(0.000)	(0.001)	(0.000)	(0.001)
Observations	1,228	1,228	1,228	1,228	1,228	1,228
Log likelihood	−1,328	−1,324	−1,003	−998.2	−1,202	−1,197

Notes: (i) SEs are in parentheses: ***P < 0.01, **P < 0.05, and *P < 0.1. (ii) Estimated coefficients for dummy variables for counties and years are reported in . (iii) Monthly mean temperature has units of degrees Celsius. Monthly total precipitation has units of millimeters.

Variables representing land use are different between models 1 and 2. Model 1 uses Shannon diversity to characterize land use diversity, whereas model 2 uses data on proportion of land use for six land use classes (one omitted from the model because the fractions sum to 1).

The reference land use class is cultivated land.

Land use diversity had no significant effect on the severity of cotton bollworm (Table 1). The proportion of cotton area that was genetically resistant to bollworm had a significant negative impact on bollworm severity, indicating that this pest control strategy was effective. The proportions of water body and unused land area, compared with cultivated land, had significant positive effects on bollworm.

The number of insecticide applications against bollworm had a significant positive impact on cotton aphid severity but a negative impact on mirid bug severity. This result indicates that the severity of cotton aphid decreased as treatments against bollworm declined over time, likely due to the restoration of aphid suppression by natural enemies. In contrast, reduced insecticide use against cotton bollworm benefited mirid bugs.

The marginal effect of a given explanatory variable calculated from an ordered probit regression captures the change in the probability of observing different classes of pest severity (levels I–V) due to a unit change in the explanatory variable (). Table 2 reports the marginal effects of explanatory variables that were found significant in the regressions reported in Table 1, evaluated at the sample means of covariates. For instance, the probability of observing the most severe class of cotton aphid severity (level V) increased by 5.5% (marginal effect of 0.0553) for a unit of increase in land use diversity (i.e., if there is one more type of land use in the county). Considering the marginal effect of land area share, for instance, a 1% increase of forest area in the county made it 28.3% more likely to have a nondamaging level of mirid bug severity (level I) while a 1% increase in water body area made it 16.1% more likely to have a nondamaging mirid bug level (level I), but the latter change made a severe bollworm outbreak (level IV) 8.4% more likely. These differences in sign of significant coefficients point to diverging effects of land use on mirid bugs and cotton bollworm, even though both are polyphagous and able to exploit multiple hosts across the landscape.

Table 2.

Estimated marginal effects† of selected explanatory variables in models for pest severity

	Model 1, with Shannon diversity					Model 2, with area proportions of different land uses
Variables	I	II	III	IV	V	I	II	III	IV	V
Cotton aphid
Land use diversity (Shannon diversity)	0.2084*	0.7449**	−0.3523**	−0.5458**	−0.0553**
Share of cotton area in total cultivated area, %	−0.0007	−0.0024*	0.0011	0.0017*	0.0002	−0.0009**	−0.0032**	0.0015**	0.0023**	0.0002**
No. of insecticide sprays against cotton bollworm	−0.0057***	−0.0204***	0.0096***	0.0149***	0.0015***	−0.0053***	−0.0189***	0.0090***	0.0138***	0.0014***
Percentage of built-up area, %						0.0060***	0.0214**	−0.0102**	−0.0157***	−0.0015**
Monthly mean temperature: July	−0.0108**	−0.0388*	0.0183*	0.0284*	0.0029*	−0.0113**	−0.0405***	0.0193**	0.0296***	0.0029**
Monthly mean temperature: August	0.0122**	0.0436**	−0.0206**	−0.0319**	−0.0032**	0.0131**	0.0467**	−0.0222**	−0.0342**	−0.0034**
Monthly total precipitation: May	−0.0001*	−0.0005*	0.0002*	0.0003*	0.0000	−0.0001**	−0.0005**	0.0002*	0.0003**	0.0000*
Mirid bugs
Land use diversity (Shannon diversity)	0.7713*	−0.2792*	−0.4219*	−0.0690*	−0.0012
No. of insecticide sprays against cotton bollworm	0.0162**	−0.0059*	−0.0088**	−0.0014**	−0.0000	0.0159*	−0.0060*	−0.0086**	−0.0014*	−0.0000
Percentage of forest area, %						0.2831**	−0.1061**	−0.1525**	−0.0241**	−0.0004*
Percentage of water body area, %						0.1606***	−0.0602***	−0.0865***	−0.0137***	−0.0002**
Percentage of unused land area, %						0.1028***	−0.0385***	−0.0554***	−0.0087***	−0.0001**
Monthly mean temperature: May	−0.0565***	0.0205***	0.0309***	0.0051***	0.0001*	−0.0587***	0.0220***	0.0316***	0.0050***	0.0001*
Monthly total precipitation: August	0.0004**	−0.0001**	−0.0002**	−0.0000**	−0.0000	0.0004**	−0.0001**	−0.0002**	−0.0000**	−0.0000
Cotton bollworm
Share of Bt cotton area in total cotton area, %	0.0004**	0.0030***	0.0006*	−0.0033***	−0.0007***	0.0004***	0.0030***	0.0006**	−0.0033***	−0.0007***
Percentage of water body area, %						−0.0102***	−0.0770***	−0.0142**	0.0838***	0.0175***
Percentage of unused land area, %						−0.0044**	−0.0334**	−0.0061*	0.0363**	0.0076**
Monthly mean temperature: June	0.0071**	0.0524***	0.0095*	−0.0567***	−0.0123***	0.0065***	0.0496***	0.0091*	−0.0539***	−0.0113***
Monthly total precipitation: June	0.0000*	0.0002*	0.0000	−0.0003*	−0.0001*

Notes: (i) †The marginal effect of a given explanatory variable calculated from an ordered probit regression captures the change in the probability of observing different classes of pest severity (levels I–V) due to a unit change in the explanatory variable. (ii) I, II, III, IV, and V in the first row refer to pest infestation levels I–V. (iii) ***P < 0.01, **P < 0.05, and *P < 0.1.

Weather, particularly temperature, had substantial effects on the severity of cotton aphid and mirid bugs during the period of 1991–2015, compared with the other significant factors (Table 2). The likelihood of a severe cotton aphid infestation (level IV) occurring was increased by 2.8% for 1 °C increase in the average temperature in July and was decreased by 3.2% for 1 °C increase in the average temperature in August. The likelihood of a low aphid infestation (level II) was decreased by 3.9% and increased by 4.4% for these same temperature changes, respectively. Mirid bug severity was generally higher in years with higher May temperature, with the chance of a no-outbreak year (level I) lowered by 5.7% for each degree Celsius increase in May temperature. Temperature in June was an important driver of cotton bollworm severity, with the chance of having a severe infestation (level IV) lowered by 5.7% for each degree Celsius increase.

The effect of pest severity on insecticide use intensity.

Cotton farmers applied on average 0.44, 1.10, 2.06, and 2.60 more insecticide sprays for cotton aphid control, and 0.93, 1.64, 2.68, and 3.73 more sprays for mirid bugs control at infestation levels II, III, IV, and V, respectively, than at level I (Table 3; full regression results are reported in ). For cotton bollworm, farmers reacted vigorously in the event of an extremely severe infestation (level V), applying on average 4.69 more sprays than at severity level I. During the study period, farmers on average sprayed 1.65 times (±1.23), 0.38 times (±0.70), and 1.01 times (±1.07) against cotton aphid, mirid bugs, and cotton bollworm at the nondamaging infestation levels (level 1), respectively. The results indicate that farmers use insecticides prophylactically for cotton aphid and bollworm, but not for mirid bugs.

Table 3.

Estimated coefficients for the frequency of insecticide sprays against cotton aphid, mirid bugs, and cotton bollworm

	No. of insecticide applications against each pest
Explanatory variables	Cotton aphid	Mirid bugs	Cotton bollworm
Infestation level: category II	0.438***	0.929***	0.540*
	(0.159)	(0.173)	(0.281)
Infestation level: category III	1.095***	1.641***	1.440***
	(0.248)	(0.185)	(0.331)
Infestation level: category IV	2.059***	2.679***	2.343***
	(0.324)	(0.282)	(0.404)
Infestation level: category V	2.598***	3.731***	4.685***
	(0.373)	(0.585)	(0.451)
Share of Bt cotton area in total cotton area, %	0.004	0.008***	−0.017***
	(0.003)	(0.003)	(0.005)
Constant	1.943***	−0.223*	5.314***
	(0.380)	(0.126)	(0.737)
Observations	1,228	1,228	1,228
R2	0.704	0.771	0.779
Adjusted R2	0.684	0.755	0.764

Notes: (i) SEs are in parentheses. ***P < 0.01, **P < 0.05, and *P < 0.1. (ii) Infestation level I is the reference category for the infestation level explanatory variables. (iii) Estimated coefficients for dummy variables for counties and years are reported in .

Discussion

Cotton aphids exploit a variety of host plants, but cotton is the single most important host. Our results indicate that counties with a greater proportion of cotton in total cultivated land had higher severities of cotton aphid. Increases in land use diversity in a county decreased cotton aphid severity, and increases in built-up area (urbanization) were associated with lower cotton aphid severity. The role of urban areas in supporting a diverse predator community in the cotton ecosystem in northern China was also found (35).

Mirids exploit a wide variety of crop hosts including cotton (31), fruit trees, and over 200 other host species (36). Therefore, it is anticipated that mirid bugs will profit rather than suffer from diverse land uses (31). However, land use diversity negatively affected mirid bug severity, while greater area proportions of forest, water body, and unused land were associated with lower mirid bug severity. These findings contradict results of earlier analyses indicating that diverse land use, through the provision of multiple host habitats, stimulates mirid bugs in the Chinese context (31). Results of the current county-level analysis, based on the largest database collected thus far, indicate that in the Chinese context advantages of diversifying the cropping landscape are more important for mirid bug control than disadvantages.

Land use diversity had no effect on cotton bollworm severity. However, the area proportions of water body and unused land were positively associated with cotton bollworm severity, which may signal the importance of a broad range of host plants during different seasons for the survival and population dynamics at the landscape level of this generalist herbivore. Another key result of this analysis is that cotton bollworm depends negatively on Bt cotton proportion. However, the estimated average (negative) effect of Bt cotton area proportion for the period of 1991–2010 (n = 1,020; mean, −0.013; SD, 0.003) was significantly lower than that for the period of 1991–2015 (n = 1,028; mean, −0.012; SD, 0.003) (t test, P value of 0.0000). This signals a slight decline in the effectiveness of Bt cotton in suppressing cotton bollworm in recent years, likely due to resistance development. Zhang et al. (37) detected cotton bollworm resistance in laboratory bioassays but noted that control failures of Bt cotton have not been reported in China. These early warnings may spur proactive countermeasures, including a switch to transgenic cotton producing two or more toxins distinct from Cry1A toxins (37).

Our results show that the adoption of Bt cotton has not only effectively brought down severity levels of the cotton bollworm but also led to significant spillover effects on cotton aphid and mirid bugs because of reduced insecticide use against bollworm. Specifically, reducing insecticide use against bollworm had helped to conserve aphid natural enemies, which subsequently regulated cotton aphid populations. This finding is consistent with that of Lu et al. (38), who found that there was a marked increase in the abundance of three types of generalist arthropod predators (ladybirds, lacewings, and spiders) and a decreased abundance of aphid pests associated with widespread adoption of Bt cotton and lower insecticide use in this crop. Previous studies identified insecticide use in cotton as an important factor affecting the population levels of cotton aphid during the 1990s and 2000s in China (34, 39, 40). Before the 1970s, cotton aphids were kept well under control by natural enemies and treating seeds with insecticides (34). During the 1970s and 1980s, aphid damage to cotton became more serious and frequent because insecticide use against cotton bollworm was wiping out most natural enemies of cotton aphid (32, 39, 41–43). Natural enemies play an essential role in population suppression of cotton aphid (44–47), and measures that allow conservation of predators in cotton fields can greatly help to control the pest.

In contrast, reduced insecticide use against bollworm was responsible for increasing mirid bug severity in cotton. This is plausible, because common insecticides that are used to control bollworm can effectively suppress mirid bug populations, while at the same time the effect of natural enemies on mirid bugs is very limited (48). Our analysis of long-term panel data confirms the report of Lu et al. (31), based on a much smaller dataset, that mirid bugs increased in importance as a pest in cotton because of increased adoption of Bt cotton and associated reduction in insecticide use.

Findings of this study have the following policy implications. First, landscape structure plays an important role in providing biological pest control ecosystem services and our analysis at the county level confirms this. However, the role of land use diversity in providing pest control services is likely to be idiosyncratic, depending on the biological features of and interactions among species, and should not be considered a “silver bullet” for managing all insect pest problems. Policies that aim at enhancing natural habitat and biodiversity conservation should anticipate the possibility of varying impacts on different insect species. Second, the large-scale expansion of Bt cotton in China has played a critical role in regulating cotton bollworm. However, the associated reduction in insecticide use against cotton bollworm has also led to effects on other cotton pests such as cotton aphid and mirid bugs. Such linkage among different pest species and unintended effects of management practices need to be accounted for when designing agricultural policies and insect management methods. Third, weather is an important driver of the system and the projected climate warming is likely to heighten pest outbreak risks. Compared with Bt cotton whose adoption nearly saturated about a decade ago, weather is gaining its momentum as an important threat. Policy makers and agricultural practitioners need to consider the impact of climate change on pest management and implications for chemical insecticide use. Fourth, cotton farmers sprayed against cotton aphid and bollworm even at nondamaging levels of infestation. This prophylactic spraying has been attributed to farmers’ risk-averse attitude and lack of knowledge (23, 29).

Alarmingly, Huang et al. (49) showed that one additional kilogram of insecticide per hectare from the current usage level by cotton farmers in North China Plain would reduce farmers’ income by 63.14 CNY/ha, indicating that the “true” marginal value of insecticides is negative. Furthermore, there is a clear lack of awareness among farmers on the value of natural enemies for crop protection, both in China (49, 50) and in other developing countries (e.g., refs. 51 and 52). Better transfer of information to farmers can support the health of crops, the environment, and farmers themselves. Policies that encourage farmers to account for the human health and environmental costs of insecticides would help to incentivize the use of natural enemies to regulate pest problems. Further research is needed to identify how policies may be designed to encourage the use of ecosystem services as a public good in farming communities and reduce the reliance on chemical insecticides (52). To understand the true (net) costs and benefits of alternative pest management approaches, the overall use of insecticides in cotton and associated economic and environmental costs to farmers and society need to be better disentangled.

We note some limitations of our study. Our land use data are based on six main categories of land use/land cover. A more disaggregated land use classification that enables the identification of main crop types would allow for deeper analysis. Additionally, as data become more available, it will be useful to explicitly address the effects of varying cropping systems and the use of other inputs in cotton on pest severity.

Materials and Methods

Data.

We compiled a unique long-term panel dataset that consists of data on (i) pest severity and insecticide applications per annum per county by pest species, (ii) land cover/use, and (iii) monthly average temperature and total precipitation data. The counties in the database are the 51 most important cotton growing counties, by production, in the Yangtze River valley and Yellow River valley cotton production regions, while the data covered the years 1991–2015, with complete coverage of counties in all years when cotton was cultivated. Between 2011 and 2013, eight counties in our sample stopped cultivating cotton. The number increased to 11 and 12 counties in 2014 and 2015, respectively, resulting in 47 missing records (8 × 3 + 11 + 12). Despite the decline in cotton area in the eight provinces since 2011 due to increasing labor cost and decreasing cotton price (53), the two production regions still accounted for 37% of the national cotton production or one-half of national cotton acreage in 2016 (54).

The national cotton pests monitoring network, maintained by the Ministry of Agriculture, mandates the main cotton-producing counties to collect yearly data on pest infestation levels and insecticide applications for key cotton pests following national standardized monitoring and categorization methods (50). Tailored scouting methods were used for different pests. In each county, 10–20 fields were selected for pest monitoring in each year. Insect populations were recorded every 3–10 d from early June to late August each year (31, 38), and the seasonal average abundance across the surveyed fields was used for scoring using a five-point scale of levels I–V (55). Data on the number of insecticide applications targeted at specific pests were collected by interviewing farmers at each scouting to estimate yearly pest-specific total number of sprays for each county. While the detailed data collection methods and protocols should inspire confidence in the data, the reliability of our pest level data depends on the accuracy, knowledge, and honesty of the respondents, as is the case with any non–first-hand data. The results of analysis are robust because of (i) the large number of counties and years in the dataset, and (ii) the use of a fixed-effect panel data model (), which controls for unobserved county-level heterogeneity.

County-level land use data were drawn from a national land cover/use database developed by the Chinese Academy of Sciences, using satellite remote-sensing data from the Landsat Thematic Mapper/Enhanced Thematic Mapper images (56, 57). The database offers the most comprehensive coverage of China’s land use/cover and has been used in a number of published studies (e.g., refs. 58–60). Six main land use classes are used in this study: cultivated land, forest, grassland, water body, built-up land, and unused land (description of land use classes is provided in ). We extracted land use data for 6 y (1990–2015 at 5-y intervals) and computed the proportional area for each land use as well as the Shannon index for land use diversity for each county in each year (data for 1990 and 2015 are provided in ). Land use proportions in intermediate years (e.g., 1991, 1992, 1993, 1994, 1996, etc.) were calculated by linear interpolation between the data.

Daily temperature and precipitation data from meteorological stations (China Meteorological Data Service Center, data.cma.cn/) were used to derive monthly data for May, June, July, and August for each year. For counties without a meteorological station, we used data collected at the nearest station.

Econometric Analysis.

We adopt a two-step approach to analyze the main drivers of cotton pest severity and insecticide use for 51 counties in China over a 25-y period. First, we identify how land use, crop management, and weather affect pest severity using a fixed-effect ordered probit regression model (61). We then estimate the relationship between pest infestation levels and the frequency of insecticide applications with a fixed-effect linear regression model. As Larsen (32) showed, landscape simplification has inconsistent effects on insecticide use over years. Therefore, multiyear studies are key to unlocking the true drivers of variation in insecticide application.

Like Larsen (32), we use a fixed-effect panel data approach to control for unobserved heterogeneity between counties and periods. Throughout, we follow the econometric use of the term fixed effects to describe models with dummy variables for each county (cross-sectional unit) or period. The fixed-effect panel data approach allows our multiyear analysis to identify the effects of explanatory variables using year-to-year variation within counties. In comparing a unit with itself, the potential for confounding due to omitted variables is minimized because time-invariant characteristics, such as slope, soil quality, location, and cultural norms or agricultural traditions, drop out of the model (32). The approach addresses county-specific heterogeneity that is difficult to observe. As a result, the estimation results are robust against bias caused by these unobserved factors. Furthermore, the inclusion of dummy variables for years account for time effects shared by all counties in a given year such as national pesticide regulation or weather anomalies (32).

Our two-step analysis approach differs from two previous studies examining the relationship between landscape simplification and insecticide use in the United States—those by Meehan et al. (62) and Larsen (32)—in three main respects. First, the use of the proportion of cropland in a county treated with insecticides as an indicator for insecticide use is reasonable for the US context but does not apply to the Chinese context because virtually all agricultural land is treated with pesticides, although at varying degree. We therefore use as a variable for the intensity of pesticide use the average number of sprays per year. Second, while Meehan et al. (60) and Larsen (32) focused on the proportion of county in cropland as a proxy for landscape simplification, we consider both land use diversity and proportions for different land uses to examine the effect of these land uses on pest severity. Third, we consider the effects of weather and the interactions among the three key pests of cotton, via Bt cotton-induced change in the frequency of insecticide use against cotton bollworm.

The nonlinear ordered probit regression is appropriate for examining the determinants of pest severity because severity is measured by five infestation levels, that is, ordered categorical outcomes. The ordered probit regression for pest severity of each pest is estimated with two specifications, which differ in how the effects of land use are captured. In the first specification, Shannon diversity index, H, as defined in Eq. , is used to quantify land use diversity, denoted by L (63, 64):

where P is the area proportion of land use i, and m is the total number of land uses (65). To capture the true diversity represented by the value of the index, we translate it back to the scale of “effective number of land uses” (65) using , which is used in the subsequent analysis. The resulting regression model is as follows:

where i and t are indices denoting county and year, respectively; S denotes pest severity measured with an ordinal variable of five infestation levels; is the land use diversity in county i in year t; is the proportion of cotton in the total agricultural area in county i in year t; is the proportion of the cotton area grown with Bt cotton in county i in year t; is the average number of insecticide applications against cotton bollworm per field in county i in year t (in the cotton aphid and mirid bug models only); is a set of four explanatory variables representing May, June, July, and August average temperature (in degrees Celsius) in county i in year t; is a set of four explanatory variables representing monthly total precipitation (in millimeters) in May, June, July, and August in county i in year t; is a set of 50 county dummy variables (one county dropped from the regression as a reference) to represent county specific effects; and is a set of 24-y dummy variables (with one year dropped from the regression as a reference) to represent year-specific effects. ε is a random error term. Table 4 reports the descriptive statistics of all of the variables.

Table 4.

Sample means of variables

Variable	No. of observations	Mean	SD	Min	Max
Infestation level: cotton aphid	1,228	2.77	1.14	1	5
Infestation level: mirid bugs	1,228	2.14	1.19	1	5
Infestation level: cotton bollworm	1,228	3.05	1.31	1	5
No. of insecticide applications: cotton aphid	1,228	3.25	2.16	0	17
No. of insecticide applications: mirid bugs	1,228	2.17	2.17	0	14
No. of insecticide applications: cotton bollworm	1,228	5.07	4.11	0	33
Share of cotton area in total cultivated area, %	1,228	19.82	13.25	0.05	67.29
Share of Bt cotton area in total cotton area, %	1,228	57.55	45.50	0	100
Percentage of cultivated area, %	1,326	54.64	10.07	18.72	72.66
Percentage of forest area, %	1,326	5.02	9.41	0.03	60.52
Percentage of grassland area, %	1,326	4.27	6.69	0	33.18
Percentage of built-up area, %	1,326	25.61	7.52	7.79	56.58
Percentage of water body area, %	1,326	9.53	8.44	0.88	37.45
Percentage of unused land area, %	1,326	0.92	3.06	0	25.79
Land use diversity (Shannon diversity)	1,326	2.98	0.59	2.04	4.18
Monthly mean temperature, °C
May	1,275	20.46	1.82	14.31	24.75
June	1,275	24.79	1.90	18.1	28.6
July	1,275	26.93	1.79	20.28	31.28
August	1,275	25.67	1.87	20	30.78
Monthly total precipitation, mm
May	1,275	76.68	67.50	0	390.5
June	1,275	108.16	95.67	1.30	645.7
July	1,275	169.27	113.60	5.7	806.3
August	1,275	146.63	105.90	0.70	639.2

In the second specification, instead of a single explanatory variable for land use diversity, , we include percentages of area for five different land uses (the sixth, cultivated land area, being one minus the other five). The share of cultivated land area in total land area is excluded from regressions as a reference, and regression coefficients are thus estimated for the effects of increasing the shares of forest, grassland, water body, built-up area, and unused land, relative to the share of cultivated area. The empirical model for pest severity (Eq. ) is estimated separately for each of the three pests: cotton aphid, mirid bugs, and cotton bollworm.

To examine the effects of pest severity on the frequency of sprays targeted at each pest, , we estimate a fixed-effect linear regression model using the ordinary least squares estimation method. We include share of Bt cotton in total cotton area (R) and dummy variables for counties (D) and years (E). Since here pest severity is a categorical explanatory variable, we use four dummy variables for infestation levels II, III, IV, and V (Slevel II, , Slevel III, , Slevel IV, , and Slevel V, ) in the regression, with level I being the reference. The model equation is as follows:

Repeated observations over time may have correlated disturbance terms. Although estimates from the fixed-effects model would remain unbiased and consistent (59), such autocorrelation in the error terms could result in artificially small SEs (66). Here, we used cluster-robust SEs clustered at the county to allow for arbitrary autocorrelation between observations within the same county when estimating Eqs. and . For nonlinear ordered probit regressions on pest severity, nonparametric bootstrap estimation was used to obtain the correct SEs, taking into account clustering. All data analyses were conducted in Stata 15 (StataCorp).