Spatiotemporal Dynamics of Bemisia tabaci MEAM1 (Hemiptera: Aleyrodidae) in Commercial Soybean Crops

Abstract

Spatiotemporal dynamics studies of crop pests enable the determination of the colonization pattern and dispersion of these insects in the landscape. Geostatistics is an efficient tool for these studies: to determine the spatial distribution pattern of the pest in the crops and to make maps that represent this situation. The aim of this study was to determine the spatiotemporal distribution of B. tabaci adults and nymphs in commercial soybean crops from planting to harvest using geostatistics. Infestation by adults and nymphs of B. tabaci started between 30 and 50 days after the emergence of the plants. The maximum population density of ten adults per plant and two nymphs per leaf was registered between 90 and 101 days after plant emergence. The colonization of soybean plants by B. tabaci may be divided into three stages: beginning infestation (at the outermost parts of the crop), whole area colonization, and dispersion colonization (when the whole crop area is infested). The density of adult insects was positively correlated with rainfall and relative humidity. Wind speed positively affected the dispersion of adult whiteflies. The distribution pattern of B. tabaci in the soybean crop was aggregated.

1. Introduction

Bemisia tabaci (Gennadius) Middle East -Asia Minor 1 - MEAM1 (Hemiptera: Aleyrodidae), commonly known as whitefly, is one of the most invasive and destructive pests for global agriculture, causing significant yield losses in commercial crops [1, 2]. This insect has high fecundity and fertility rates. It is polyphagous and occurs in areas of temperate, subtropical, and tropical climates [3, 4, 5]. The whitefly has a wide host range; thus, crops such as cotton, beans, squash, melon, and tomato may serve as alternate host plants for this species infesting soybean crops.

The whitefly causes several problems in Brazilian crops, and the insect may cause direct damage by sucking phloem sap. Indirect damage is linked to the excretion of honeydew, which serves as a substrate for the growth of opportunistic fungi (Capnodium sp.) [6, 7, 8]. B. tabaci may inhibit the gathering of carotenoids and chlorophyll affecting the photosynthetic rate of plants [9]. In soybean [G. max (L.) Merr., Fabaceae], B. tabaci is the vector for the cowpea mild mottle virus (CpMMV) [10].

In recent years, soybean has been highly attacked by high whitefly populations, mainly in the reproductive stages (R1 stage) Until very recently, this species was only seen as an occasional pest, however it is now considered a key pest in soybean. Although the importance of the whitefly B. tabaci MEAM1 has increased its incidence in soybean cultivars [11, 12, 13], there is no sampling method currently defined for this pest in this crop. The recommendations do not specify the sample numbers, size of area to be sampled, location of the different stages, or the most suitable forms of evaluation. In general, management is based on information obtained from other countries and what has been observed in other crops, such as cotton [14, 15, 16] and melon [17], in the state of Arizona, USA. In Brazil [18].

Understanding the spatiotemporal dynamics of pest insects in crops provides important information that can be incorporated into integrated pest management programs [19, 20, 21, 22]. The study of the spatiotemporal dynamics of insects using geostatistics can provide important information about the pattern of colonization, aggregation, and dispersion of pests in the field [20, 22, 23, 24, 25]. This statistical tool uses a method that characterizes spatial variation by comparing similarities between distant and proximal points.

This technique provides results that allow colonization maps to be produced, zoning the different densities and determining the pattern of spatial distribution of the insects in the field [20, 25]. A sequence of these maps during crop development may indicate areas that demand greater attention for the observation of pest infestation over time. In these spatiotemporal dynamic’s studies, plants in the landscape that can serve as pest sources can be identified and located. This may be useful for farmers when planning the location of crops in the landscape or even for predicting where pest attacks will start within the crop area [24, 26, 27].

Considering the importance of studying the spatiotemporal dynamics of insect pests and to the best of our knowledge, there are few published studies on B. tabaci in soybean. Thus, the aim of this study was to determine the spatiotemporal distribution of adults and nymphs of B. tabaci in soybean crops from planting to harvest by using geostatistics and verifying abiotic factors (temperature, relative humidity, wind, and rainfall) that affect the dispersion of these insects.

2. Materials and methods

The experiment was carried out at Celeiro Seeds Farm, and the farm’s total area was estimated at 16,000 hectares, where grains and soybean seeds are produced. The farm is located in the municipality of Monte Alegre do Piauí and Serra do Quilombo in the Brazilian state of Piaui (09°21′12”S; 45°07′42” W, 642 m) (Figure 1A). The region’s soil is characterized as a yellow Latosol (Oxisol)28, and the climate is tropical with a dry season in winter, Aw, according to the Köppen-Geiger classification 29. Climatic data such as the mean air temperature (°C), relative humidity (%), wind speed (m/s) and wind direction were obtained from NIMET (National Institute of Meteorology), and rainfall data were obtained with a rainfall meter during the period of execution of the experiment (Figure 2). This work was carried out for 2 seasons (2017/2018 and 2018/2019). During the first season, the cultivars Brasmax BonoIPRO® and BRS 9180IPRO® were used; in the second season, the cultivars were Brasmax Bonus IPRO® and Brasmax Extreme IPRO®. The characteristics of the cultivars are shown in Table 1.

Figure 1.

(A) Brazil’s map - Monte Alegre do Piauí, Piauí. (B) sample grid of the evaluated points and geographic location. (C) sampling unit (3 plants equidistant from the central point, total of 12 plants per unit).

Figure 2.

Rainfall, temperature and average relative humidity of the region. Celeiro seed farm in Serra do Quilombo, Monte Alegre do Piauí. National Institute of meteorology (INMET) data. * sowing month of cultivars.

Season	I (November 2017–April 2018)		II (November 2018–April 2019)
Variety	Brasmax BonusIPRO ®	BRS 9180IPRO®	Brasmax BonusIPRO ®	Brasmax ExtremaIPRO ®
Field	Field I	Field II	Field III	Field IV
Planting date and characteristics	Maturity group 7.9 - average cycle of 108 days. Planting on November 11, 2017	Maturation group 9.1 - average cycle 119–139 days. Planting on November 18, 2017.	Planting on November 27, 2018.	Maturation group 8.1 - average cycle of 110 days. Planting on December 4, 2018
Stand and spacing	390 thousand plants/ha. Spacing between lines of 0.50 m	200 thousand plants/ha. Spacing between lines of 0.45 m	390 thousand plants/ha. Spacing between lines of 0.50 m	200 thousand plants/ha Spacing between lines of 0.45 m
Fertilization	620 kg/ha of Simple SuperPhosphate, 200 kg/ha of Potassium Chloride, 3 kg/ha of copper sulphate, 2.5 kg/ha of Manganese Monoxide, 3 kg/ha of Ulexite and 20 kg/ha of MIB Granary.	200 kg/ha Simple SuperPhosphate, 200 kg/ha of Potassium Chloride and 20 kg/ha of MIB Granary in the planting line.	630 kg/ha of Simple SuperPhosphate, 200 kg/ha of Potassium Chloride, 20 kg/ha of Manganese Monoxide, 6 kg/ha of Ulexite and 13 kg/ha of Zinc Sulphate, by haul.
Nearest native forest	1.79 km to the North	4.09 km to the West	4.22 km North	0.32 km Northwest
Assessments days	14, 28, 42, 56, 70 and 84 DAEa	17, 31, 45, 59, 73, 87 and 101 DAE.	13, 27, 41, 55, 69, 83 and 87 DAE	20, 34, 48, 62, 76, 90 and 104 DAE
Geographic coordinates	9°23′27.21”S; 45° 6′56.62”W	9°25′3.54”S; 45° 0′23.63”W	9°24′10.88”S; 45° 8′37.50”W	9°22′25.27”S; 45° 6′56.69”W

Table 1.

Phenological cycle, fertilization, plant population, planting date, distance from cultivars to the nearest native forest, and evaluation days by cultivar in the two seasons studied.

^a

= Days After Emergence.

The area of 1.6 ha (16,000 m²), with 100 m length and 160 m width, was divided into 160 sampling units, with 10 m². In the center of each sample unit, sets of four ordered points, equally spaced 2 m from a central point, were evaluated to improve the estimate of the nugget effect. Each sampling point was georeferenced with the aid of a GPS device model GPSMAP 60CS®x (Garmin). The distance from each field to the nearest native forest is described in Table 1.

2.1 Evaluation of the B. tabaci density

The density of insects was evaluated in the sampling grid of 160 points previously georeferenced (Figure 1B). At each georeferenced point, 12 plants were evaluated; four subsamples of three plants were taken (sampling unit). The subsamples were 2.0 m equidistant from the central point (Figure 1C). The evaluated plants were positioned along a regular grid pattern throughout the crop cycle to obtain systematic sampling points and avoid directional trends.

The third leaflet from top to bottom (apical third) of each plant was evaluated by direct counting. Leaves were handled with care, and the nymphs and adults of B. tabaci were counted. A magnifying glass was used to count the number of nymphs. We evaluated these leaves in particular using a direct counting method because these are the ideal sample type and technique for assessing the density of B. tabaci nymphs and adults in soybean crops 14. The sampling schedule was standardized; the evaluations started at 7 am and ended at 11 am. The sampling schedules for each cultivar are described in Table 1.

2.2 Analysis of the spatial distribution of B. tabaci in soybean

The data on adult and nymph densities of B. tabaci were submitted to statistical analysis. Subsequently, principal component analysis was performed between the range and density of whitefly adults with climatic parameters (temperature, relative humidity, wind speed, and rainfall). These analyses were performed using R® software. 30 All geostatistical procedures were performed using the Geostatistical Analyst Tool for ArcGIS 10.5 (ESRI). The procedure in this tool can be simplified as follows: map and examination of data; preprocess of data if necessary (transform, detrend, decluster); definition of spatial structure model; definition of search strategy; prediction of values at unsampled locations; quantification of uncertainty of the predictions; checking if the model produces reasonable results for predictions and uncertainties and, using the information in risk analysis and decision making.

Subsequently, geostatistical analysis was performed; the spatial patterns and interpolations were determined using the parameters of adjusted experimental semivariograms and ordinary kriging, respectively. The semivariograms were calculated from the primary collected data, allowing us to detect differences between pairs of sampled points in relation to the distances (this procedure was used to adjust the theoretical semivariogram). Once the semivariance increases, there is a spatial dependence relationship between the densities of B. tabaci and the sampled points. The nugget effect and the sill value were calculated for each of the adjusted models (spherical, exponential, and Gaussian).

The experimental semivariograms were adjusted to the theoretical models, and the selection was made based on the cross-validation parameters, in which the measured and estimated values were compared using the standard mean error (SME) and the root mean square standardized error (RMSSE).

The kriging indicator was applied to spatial distribution maps of the insect in the crops, with the objective of modeling the probability of unsampled locations exceeding the quality reference values (QRVs). The kriging indicator does not use normal distribution assumptions, because it transforms the original data into binary values and, consequently, into cumulative distribution functions from pre-established values, in this case the QRVs. 31, 32 The data were transformed into log-normal values to minimize distribution errors and meet the requirements of common kriging based on the rejection of the null hypothesis of the Kolmogorov-Smirnov test for normal distribution. Spatial variability was determined from isotropic and anisotropic semivariograms. Anisotropic calculations were performed in four directions (0, 45, 90 and 135°).

To obtain reliable estimates, the theoretical model needs to show SME values close to 0 and RMSSE close to 1.33 The Akaike information criterion (AIC) was used as the last selection criterion. The spatial dependence rate (SDR) was calculated according to the formula: [C0/(C0 + C1) × 100]. 34 A nugget effect less than or equal to 25% of the plateau was considered strong. The value was considered moderate when it was between 25 and 75% and weak when it was above 75%. Variables that showed SDR less than one unit were not considered.

Spatial distribution maps of B. tabaci adults and nymphs were prepared for each crop. In these maps, the predominant direction of the winds is shown.

3. Results

According to the statistical analysis, the values of variance were within the range 0 to 7.52 and were predominantly positive, where the peak of the normality curve of the results was higher than the standard value, the distribution curves were classified with positive asymmetries, and the coefficient of variation values were predominantly greater than 30% (Online Resource 1). In the two seasons evaluated, infestation by adults and nymphs of B. tabaci generally started between 30 and 50 days after the emergence of the plants. The maximum population density of 10 adults per plant and approximately two nymphs per leaflet were recorded 101 days after plant emergence (Online Resource 1 and Figure 3).

Figure 3.

Occurrence of Bemisia tabaci adults and nymphs in the Brasmax BonoIPRO®, BRS 9180IPRO® and Brasmax ExtremaIPRO® cultivars for seasons I, II, and l.

The results that showed a strong or moderate degree of spatial dependence (SDR) were subjected to geostatistical models and chosen based on cross-validation parameters using the standard mean error (SME) and the root mean square standardized error (RMSSE). Of the 126 models that were processed, 33 were selected, of which 13 were exponential, 10 Gaussian, 10 spherical, and 36 pure nugget effects (Table 1). All selected models were isotropic (i.e., the spatial autocorrelation was the same in all directions).

Differences in the variogram parameters (nugget, sill, and range) were observed in the three cultivars in the two sampled seasons. A nugget effect was observed in the variograms, always at a low density of B. tabaci adults or nymphs (Figures 4 and 5). The SDR of the models ranged from 0 to 100. This showed a significant effect of the nugget effect on the interpolation for some sampled days. These proportions displayed that the spatial component accounted for 66% of the total spatial variance. From the selected models, 15.15% demonstrated a strong SDR (<0.25), 51.52% moderate SDR (between 0.25 and 0.75), and 33.33 weak SDR (Table 2).

Figure 4.

Kriging maps for the spatial variability of Bemisia tabaci adults in the Brasmax BonoIPRO®, BRS 9180IPRO® and Brasmax ExtremaIPRO® cultivars for the 2017/2018 season I and 2018/2019 season II.

Figure 5.

Kriging maps referring to spatial variability of Bemisia tabaci nymphs in the Brasmax BonoIPRO®, BRS 9180IPRO® and Brasmax ExtremaIPRO® cultivars for the 2017/2018 season I and 2018/2019 season II.

Variety	Stages	SAD	Model	C0a	C0+C1b	RDS (%) c	RDS	Range (m)	RMSSEd	SME and
Season I (November 2017–April 2018)
Brasmax BonusIPRO ®	Adult	56	Exponential	0.06	0.08	72.02	Moderate	24.91	0.98	0.29
		70	Exponential	0.09	0.12	74.19	Moderate	24.91	0.98	0.35
		84	Exponential	0.14	0.24	60.06	Moderate	45.41	1.07	0.46
	Nynph	70	Exponential	0.00	0.00	97.83	Weak	34.40	1.00	0.05
		84	Exponential	0.00	0.02	0.00	Strong	27.01	0.96	0.11
BRS 9180IPRO®	Adult	45	Exponential	0.00	0.16	0.00	Strong	27.74	0.99	0.33
		59	Exponential	0.00	0.25	0.00	Strong	52.13	1.02	0.32
		73	Exponential	0.00	2.01	0.00	Strong	41.37	0.81	1.02
		87	Exponential	0.00	2.89	0.00	Strong	33.73	0.94	1.33
		101	Exponential	2.28	8.20	27.78	Moderate	52.59	0.98	2.29
	Nymph	73	Exponential	0.04	0.11	40.00	Moderate	20.50	0.96	0.33
		87	Exponential	0.12	0.42	27.87	Moderate	54.50	0.96	0.52
		101	Exponential	0.53	1.27	41.24	Moderate	132.53	0.99	0.86
Season II (November 2018–April 2019)
Brasmax BonusIPRO ®	Adult	41	Spherical	0.00	0.00	100.0	Weak	72.82	1.06	0.03
		55	Spherical	0.04	0.06	70.0	Moderate	128.48	1.00	0.22
		69	Spherical	0.51	0.92	55.4	Moderate	45.18	0.96	0.83
		83	Spherical	0.28	0.53	52.8	Moderate	42.39	0.94	0.62
		97	Spherical	1.38	2.51	55.2	Moderate	25.20	0.89	1.48
	Nymph	41	Spherical	0.00	0.00	98.2	Weak	25.34	1.10	0.06
		55	Spherical	0.03	0.03	95.1	Weak	24.94	1.04	0.18
		69	Spherical	0.12	0.16	74.6	Moderate	20.82	0.96	0.40
		83	Spherical	0.54	0.55	96.5	Weak	45.18	0.94	0.76
		97	Spherical	0.73	0.73	100.0	Weak	128.48	0.94	0.88
Brasmax ExtremaIPRO ®	Adult	48	Gaussian	0.064997639	0.08	80.48	Weak	83.89	1.01	0.26
		62	Gaussian	0.108601774	0.17	65.00	Moderate	21.62	1.01	0.39
		76	Gaussian	1.167979313	1.75	66.74	Moderate	41.45	0.96	1.17
		90	Gaussian	2.574639589	3.34	77.07	Weak	24.10	0.99	1.79
		104	Gaussian	1.444284816	5.23	27.64	Moderate	67.53	1.02	1.32
	Nymph	48	Gaussian	0.01394965	0.02	87.38	Weak	13.03	0.95	0.13
		62	Gaussian	0.033389883	0.03	95.60	Weak	29.94	0.98	0.19
		76	Gaussian	0.167840405	0.41	41.16	Moderate	19.93	0.98	0.65
		90	Gaussian	1.675549944	2.05	81.93	Weak	35.23	0.99	1.38
		104	Gaussian	0.823629539	1.81	45.55	Moderate	22.01	0.96	1.18

Table 2.

Models and parameters estimated by semivariogram for Bemisia tabaci adults and nymphs in the cultivars Brasmax BonoIPRO®, BRS 9180IPRO® and Brasmax ExtremaIPRO.

^a

Nugget Effect.

^b

Threshold.

^c

(C0/C0 + C1) * 100) = Spatial Dependence Rate.

^d

Root Mean Square Standardized Error; and: Standard Mean Error.

The ranges of the models varied from 13.03 to 132.53 m, and the maximum range obtained for adults was 67 m at 104 DAE in the Brasmax ExtremaIPRO® cultivar and approximately 130 m for nymphs 101 DAE in the BRS 9180IPRO® cultivar (Table 2).

We observed that colonization of soybean plants by B. tabaci may be divided into three stages: the start of infestation, colonization of the area, and dispersion in the area. The adults of B. tabaci began to infest the experimental area from the outermost area; in season I, the insects came from the west, and in season II, the insects originated from the north (Figure 4).

Colonization by B. tabaci nymphs occurred as adults disperse in the area. Near the end of the culture cycle, adults, and nymphs of B. tabaci had already colonized the whole experimental area. With the colonization of the total area, adults began the process of migration to nearby areas (Figures 4 and 5).

The soybean cultivars were divided into three groups according to the density of the pest. In the first group (Brasmax Bonus IPRO®), the lowest densities of the pest (2 adults and 1 nymph per leaf) were observed. In the second group (Brasmax ExtremaIPRO®), lands with a higher density of nymphs (up to 10 nymphs per leaf) were observed. In the third group (BRS 9180IPRO®), plants exhibited the highest densities of adults (up to 15 adults per leaf) (Figures 4 and 5).

We observed two patterns of variation in the densities of B. tabaci throughout the development of the plants. The first occurred in Brasmax Bonus IPRO®, where the pest population varied little over time and remained at low density. The second pattern occurred in Brasmax ExtremaIPRO® and BRS 9180IPRO®, where the pest population increased over time and was distributed throughout the crop area (Figures 4 and 5).

The mean air temperature was high (29–33°C) during the cropping period. On the other hand, rainfall ranged from 10 to 127 mm per month during the cropping period (Figure 2). The winds during the experimental periods occurred predominantly in the northeast direction in 70% of cases with a maximum speed of 2.1 m/s. Performing the main component analysis, it was apparent that the two axes explained 74.79% of the data variation. There were positive correlations between the range and wind speed, and this trend was observed also, in addition, between adult whitefly densities and rainfall and relative humidity (Figure 6).

Figure 6.

Principal component analysis of the effects of rainfall and climatic elements on the range and abundance of Bemisia tabaci in soybean crops. The first and second axes explained 74.79% of the variance.

4. Discussion

The data presented here provide information on the spatiotemporal dynamics of B. tabaci in commercial soybean crops. Under the conditions evaluated, the initial foci of colonization by adults occurred near the end of the vegetative stage of the crop, and it is common to observe the occurrence of nymphs and adults of B. tabaci at this stage of plant development . The peak infestation by adults and nymphs recorded in the present study was relatively low compared to other studies [22, 28].

The cultivar Brasmax BonoIPRO® had high attractiveness for adults and nymphs of B. tabaci. However, there was low infestation in the two evaluated seasons, which may indicate that cultivars of an early cycle with early planting may escape the whitefly attack [29], as mentioned above. B. tabaci infestation was more severe in late February and early March.

B. tabaci showed an aggregated distribution in the soybean crops, indicating that there were factors that influenced this outcome. The spatial distribution of pest insects in crops is a consequence of the colonization and dispersion of these organisms at these sites [24, 30]. Climate elements, terrain topography, characteristics of insect species, and landscapes are among the factors influencing the spatial distribution of pest insects in crops [11, 24, 31, 32, 33, 34].

Some results obtained for B. tabaci adults and nymphs showed a weak degree of spatial dependence. This finding is characterized as a pure nugget effect (total absence of spatial dependence), probably due to the low infestation by adults and whitefly nymphs in each evaluation, making the presence of insects in geostatistical analyses imperceptible [25, 35]. These results are a common finding in this type of study; other authors with experiments on the spatial distribution of Euschistus heros (Fabricius, 1798) (Hemiptera: Pentatomidae) in soybean have made similar observations [25, 36].

According to the semivariograms made, whitefly adults and nymphs had an isotropic distribution. This finding confirms the premise that B. tabaci disperses across the field in a diffusion pattern, allowing us to achieve an equal or similar number of adults throughout the sampling points regardless of direction, with the center of origin as the beginning of the evaluation [3]. Isotropic and aggregate pattern distributions were observed for adults and nymphs, and a reduction in range (dispersion radius) was observed as the density of adults per plant increased, while for nymphs, the opposite effect was verified. Other authors have observed the same pattern of whitefly distribution in soybean and other crops [3, 20, 22]. In experiments focused on the spatial distribution of B. tabaci in commercial watermelon crops, it was demonstrated that the dispersion of the insect was not influenced by a physical barrier or height gradient [20].

With the kriging maps, it was possible to verify that colonization of the area by B. tabaci proceeded in the western and northern regions. Near the experimental area, there was a permanent preservation area, native forest. Therefore, when favorable environmental conditions existed, the infestation of B. tabaci in soybean crops was to be expected. Migration from the native forest to the crop occurred when there was abundant food and shelter [26]. However, the northeast winds were predominant in both seasons; that is, the adults of B. tabaci flew against the wind searching for food. This fact is common since wind helps in the dispersion of volatiles in the area. Insects use volatiles as olfactory cues to find food, so insects can fly upwind, cross-wind and downwind when having perceived volatiles [37].

Kriging and semivariogram maps are fundamental tools that can assist in the identification and determination of the spatiotemporal dynamics of pests, in addition to providing information such as the dispersion pattern and range of insects. Based on this information, an adequate and efficient pest management monitoring program may be defined [20, 22, 25]. In addition, the elaboration of the maps allows the identification of the places with the highest incidence of whiteflies in crops, as well as the delimitation of the direction of colonization by adults. This allows the application of efficient control measures [20].

Geostatistics allowed us to verify the movements of B. tabaci adults and nymphs during the evaluations and soybean cycle. As previously explored, there are factors that can affect the distribution of pests in the field, features linked to insect characteristics such as population growth (reproduction, mortality) and dispersion (immigration, colonization, emigration) [3, 22, 25, 27, 36]. B. tabaci may fly to a height of up to 7 m and a distance of 7 km [38]. Climatic factors such as rainfall and the relative humidity of the air favored an increase in the density of whitefly adults, while wind speed favored the dispersion of adults, according to the literature [3, 20]. The development of monitoring and management techniques aimed at the precise control of B. tabaci depends on a better understanding of the flight behavior of the whitefly, especially short-range migration. The present study provides important information about the pattern of aggregation and distribution of B. tabaci, which may be used in future work to assist in the sampling method.

Acknowledgments

We thank the National Council for Scientific and Technological Development — CNPq — Brazil for financial support. We also thank the Coordination for the Improvement of Higher Education Personnel (CAPES), Brazil and the Piauí State Research Foundation (FAPEPI) for the scholarships and resources provided. To the Federal University of Piauí, Brazil, for providing logistical support. We also thank the Celeiro Seeds Farm for the supply of the area and inputs.