An IPSO-BP neural network for estimating wheat yield using two remotely sensed variables in the Guanzhong Plain, PR China

Computers and Electronics in Agriculture 169 (2020) 105180

Contents lists available at ScienceDirect

Computers and Electronics in Agriculture journal homepage: www.elsevier.com/locate/compag

An IPSO-BP neural network for estimating wheat yield using two remotely sensed variables in the Guanzhong Plain, PR China

T

⁎

Huiren Tiana,b, Pengxin Wanga,b, , Kevin Tanseyc, Shuyu Zhangd, Jingqi Zhanga,b, Hongmei Lid a

College of Information and Electrical Engineering, China Agricultural University, East Campus, Beijing, PR China Key Laboratory of Remote Sensing for Agri-Hazards, Ministry of Agriculture and Rural Affairs, Beijing, PR China c School of Geography, Geology and the Environment, Centre for Landscape and Climate Research, University of Leicester, Leicester LE1 7RH, UK d Remote Sensing Information Center for Agriculture of Shaanxi Province, Shaanxi Meteorological Bureau, Xi'an, PR China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Vegetation temperature condition index Leaf area index Integrated index BP neural network Winter wheat yield estimation

Early and accurate information of crop growth condition is vital for agricultural industry and food security, which gives rise to a strong demand for timely monitoring crop growth condition and estimating crop yields. This study selected the remotely sensed leaf area index (LAI) and vegetation temperature condition index (VTCI) which closely relate to crop growth and crop water stress as two key variables for indicating crop growth condition and estimating crop yields in the Guanzhong Plain, PR China. The single VTCI, the single LAI and the combination of VTCI and LAI at four growth stages of winter wheat (the turning green, jointing, heading-filling, and dough stages) were used as three input variable schemes of the back propagation (BP) neural network and the improved particle swarm optimization algorithm (IPSO)-BP neural network using a nonlinear decreasing inertia weight, respectively. The relative importance of the input variables to the output variable, yield of winter wheat, was used to determine the weight values of input variables at each growth stage. Based on the weights, the integrated index (I) was established, and then three linear regression models (weighted VTCIs, weighted LAIs, and I values) were established with yield data to estimate winter wheat yields. By calculating several statistical functions, i.e., coefficient of determination (R2) and probability value (P), the model between the I values and wheat yield performed better than those between the weighted VTCIs or weighted LAIs and wheat yields. The yield estimation model of I values by using the IPSO-BP neural network (R2 = 0.342) was found to be better than that using the BP neural network (R2 = 0.310). Therefore, we applied the model with better performance (R2 = 0.342) to map the regional winter wheat yields pixel by pixel in the Guanzhong Plain during 2011–2018, and analyzed the spatial and temporal characteristics of the estimated yields. Regarding the spatial distribution, the yields in the west part of the plain are the highest, followed by the central part, and the yields in the east part are lowest, consistent with previous studies. The estimated yields showed inter-annual fluctuations along with an increasing trend on the whole. Winter wheat yields were most depleted in 2013 and most abundant in 2015. These results were consistent with the actual situation of winter wheat production in the plain, which indicated that I can be used to provide a better quantification for monitoring regional winter wheat growth conditions and estimating crop yield. Thus, the approach of this study can provide significant benefit for regional crop production monitoring.

1. Introduction Agricultural production is the foundation of human society. In recent years, the concerns about national food security and sustainable agricultural development have increased (Gaso et al., 2019). Improving crop yield production of supply and demand while reducing operating costs and environmental pollution is a key component in precision

agriculture. Therefore, a fast and accurate yield estimation is essential for crop management, food security, food trade, and policy-making, in particular for major crops such as wheat (Yang et al., 2019). Existing approaches rely on ground measured data and variables (such as remotely sensed data, weather and soil properties) related to crop growth to model crop yield (You et al., 2017). These approaches had been widely applied for developing regression models between

⁎

Corresponding author at: P. O. Box 116, China Agricultural University, East Campus, Qinghua East Road No. 17, Haidian, Beijing 100083, PR China. E-mail addresses: [email protected] (H. Tian), [email protected] (P. Wang), [email protected] (K. Tansey), [email protected] (S. Zhang), [email protected] (J. Zhang). https://doi.org/10.1016/j.compag.2019.105180 Received 9 September 2019; Received in revised form 6 December 2019; Accepted 25 December 2019 0168-1699/ © 2019 Published by Elsevier B.V.

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

thereby modeling the complex relationship between the input and the output. It has the unique advantages of distributed parallel processing, nonlinear mapping and adaptive learning, which provides an effective approach for the simulation and evaluation of nonlinear complex system processes (Chlingaryan et al., 2018). In recent years, the research of neural networks has been applied to agriculture (Lokers et al., 2016; Kamilaris et al., 2017; Kamilaris et al., 2018). Fortin et al. (2011) compared the accuracy of multiple linear regression and ANN models in potato yield forecasting in eastern Canada, and the results showed that predicted yields based on ANN were more consistent with field measurements. Pantazi et al. (2016) introduced three self-organizing map models of counter-propagation artificial neural networks (CP-ANNs), XY-fused networks (XY-Fs) and supervised Kohonen networks (SKNs) for predicting wheat yield in a 22 ha field in Bedfordshire, UK. Results showed that the SKN model had the best overall performance (81.65%). The particle swarm optimization (PSO) algorithm has attracted attention from the research community since it tends to find the global solution and has shown substantial promise as a robust method for solving continuous nonlinear optimization problems resulting in a more stable convergence. However, this approach has been mainly applied in soil nutrient management (Fu et al., 2010), forest fire susceptibility (Tien Bui et al., 2017), plant growth (Qi et al., 2010), rice farm quality inspection (Camci et al., 2018) and agriculture water resources allocation (Habibi Davijani et al., 2016) to date. These works considered the impact of the factors related to the target objects. None of these approaches had been utilized in crop yield modeling or used to quantify the impact of the factors on the target objects. The PSO algorithm is easily trapped into local optima with low accuracy in convergence (Chen et al., 2014). In this study, we propose an IPSO using a nonlinear decreasing inertia weight to train BP neural networks effectively for optimizing the performance of models in estimating yield of winter wheat. In addressing this issue, we will overcome the limitations described above and obtain comprehensive growth information of the main growth stage of winter wheat. The objectives of this study are to: (1) calculate the weight coefficients and thresholds of the VTCI and LAI at the four growth stages of winter wheat by using the BP neural network and the IPSO-BP neural network, respectively, and establish an integrated index (I) during the main growth period for improving the estimation accuracy of regional winter wheat yields, (2) construct winter wheat yield estimation models using the weighted LAIs, the weighted VTCIs and I respectively, and select the optimal estimation model to verify its applicability, and (3) present the results of estimation yields in the Guanzhong Plain, PR China from 2011 to 2018.

crop yield and Vegetation Index (VI) directly, which is a basic modeling approach for estimating crop yields at regional scales (Becker-Reshef et al., 2010; Bolton et al., 2013). Becker-Reshef et al., (2010) developed a winter wheat yield model based on generalized regression, which used the seasonal maximum normalized difference vegetation index (NDVI) as the main predictor and achieved accuracy between 7% and 10% of regional wheat estimated yields in Ukraine six weeks before harvest. Prasad et al. (2006) developed a linear model that used AVHRR-NDVI, soil moisture, land surface temperature and rainfall to predict Iowa maize and soybean yields. The results indicated that predicted values very close to observed ones for maize (R2 = 0.78) and for soybean crop (R2 = 0.86) in Iowa, USA. The National Agricultural Statistics Service (NASS) monitored crop growth based on the NDVI obtained by NOAA-AVHRR during the crop growth period, and provided independent national scale crop growth information for decision makers in Department of Agriculture, USA. However, in addition to vegetation indices, crop yields were also affected by temperature, soil moisture and crop growth status (Wang et al., 2018; Johnson, 2016). Many studies have shown that the degree of water stress on crops is closely related to the increase or decrease of grain yield. In order to obtain a more accurate yield estimation model, various indices such as water stress index and growth state index during crop growth and development period can be established. Among them, the vegetation temperature condition index (VTCI) is based on whether the two–dimensional scatterplot of the land surface temperature (LST) and the NDVI is triangular at the regional scale. The VTCI is a near real-time drought monitoring index and has been widely used in regional drought monitoring and yield estimation (Wang et al., 2001; Sun et al., 2008; Tian et al., 2016). Leaf Area Index (LAI) is an important agronomic indicator of crop productivity, photosynthetic capacity, and degree of stress, and has a stable correlation with crop yield (Liang et al., 2015; Darvishzadeha et al., 2008). Wang et al. (2018) constructed an integrated maize growth monitoring index (G) based on the weight values of VTCI and LAI, and then linear regression model between maize yields and G values was established in the central plain of Hebei Province from 2010 to 2015. The results showed that the correlation between G values and maize yields was higher than those between maize yields and VTCI or LAI alone. Wang et al. (2016) found that under the conditions of rainfed farming, the yield estimation model based on the two variables of LAI and VTCI assimilated in the main growth period of winter wheat was with the highest precision (R2 = 0.531). The LAI has the ability to drive CO2 assimilation and dry matter accumulation, and is a key indicator of potential grain yield. The VTCI is an effective index for monitoring crop water stress during the crop growing season, and water stress is one of the important factors affecting crop yield. Combining LAI and VTCI provides complementary information for monitoring crop growth and estimating crop yield. Therefore, the monitoring of crop growth with two or more variables can reflect the agricultural situation in a timely and comprehensive manner. Farmland ecosystems are complex. Crop yield estimation is typically complex and is a non-linear process that is difficult to assess. The main disadvantages of these regression models is that they fail to consider nonlinear relationships of factors which influence crop growth and are not so good for modeling complex data in comparison to artificial intelligence based modeling (Kitchen et al., 2003; Miao et al., 2006; Papageorgious et al., 2013; Singh et al., 2018). Mkhabela et al. (2011) found that there was possibly a nonlinear relationship between MODISNDVI and barley, canola, pea as well as spring wheat yields. Therefore, current models to estimate crop yields are not always satisfactory. To establish a more practical model of crop yield estimation, deep learning methods have been explored and developed, such as artificial neural networks (ANN), convolutional neural networks (CNN), deep neural networks (DNN), and deep belief networks (DBN). In the ANN, the observed data sets of selected variables are fitted, and the connection weights of input and output variables are continuously adjusted,

2. Materials and methods 2.1. Study area The Guanzhong Plain is located in the Wei River Basin in the central part of Shaanxi Province, PR China, stretching from Baoji city in the west to Tongguan county in the east, the Qinling Mountains to the south and Loess Plateau to the north of Shaanxi Province, with coordinates of 106°22′E–110°24′E and 33°57′N–35°39′N, and it includes five prefecture-level cities, which are Xi'an city, Tongchuan city, Baoji city, Xianyang city and Weinan city (Fig. 1). The Yangling Agricultural Hitech Industries Demonstration Zone is also located in the study area. The study area has relatively flat terrain and has fertile soil. The climate is continental monsoon. The annual average temperature is between 6 and 13 ℃, and the annual average rainfall ranges from 550 to 700 mm. The plain is in a warm temperate transitional zone between semi-humid and semi-arid climates, and has good natural conditions for planting crops such as wheat and maize. The main prevailing planting pattern in irrigated areas of the plain is winter wheat-summer maize rotation, while the prevailing planting pattern in rainfed areas is winter wheat. The plain is the most important farming area in Shaanxi Province and 2

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

Fig. 1. The location and county boundaries of the Guanzhong Plain in Shaanxi Province, PR China.

the initial LST and surface reflectance images were spliced and resampled, and their projection was converted from sinusoidal to a Lambert azimuthal projection using the MODIS reprojection tool (MRT) developed by NASA. The reflectances in bands 1 and 2 from MYD09GA were selected to calculate the daily NDVI products, and then applying the maximum value composite method to generate the NDVI and LST maximum synthetic images at 10-day intervals. Thus, VTCI time series data generated at 10-day intervals were generated according to Eqs. (1)–(3). The value of VTCI ranges from 0 to 1, and VTCI has capability of drought stress classification which therein lower VTCI is for drought and higher one for wet conditions. According to the crop phenological information obtained from the China Agricultural Information Network (http://www.agri.cn/), the planting and harvest dates for winter wheat growth and development stages in different parts of the study area are similar. Winter wheat is sowed in middle October and is harvested in early June of the following year. There are four main growth stages: the turning green stage (from 1st March to 20th March), the jointing stage (from 21st March to 20th April), the heading-filling stage (from 21st April to 10th May) and the dough stage (from 11st May to 31th May). The 10-day VTCI images were converted to VTCI images at the four growth stages by calculating the average values of the 10-day VTCIs for intervals belonging to each growth stage pixel by pixel. Then, the VTCI values of counties were generated by calculating the average VTCI values of pixels located within each county.

one of the main grain production areas in China (Li et al., 2014). Spring drought and early summer drought occur frequently in this region because of the low rainfall and its uneven distribution, which makes it difficult to meet the water demand of crops throughout the growing season. Effective water resource and crop management in this study area are important practices that influence yield. 2.2. Data description and preprocessing 2.2.1. Remotely sensed VTCI The VTCI is based on the shape of the scatter plots of LST and NDVI being triangular at a regional level (Gillies et al., 1995, 1997). It is defined as follows (Wang et al., 2001; Sun et al., 2008; Wan et al., 2004):

VTCI =

LSTmax (NDVIi ) − LST (NDVIi ) LSTmax (NDVIi ) − LSTmin (NDVIi )

(1)

LSTmax (NDVIi ) = a + bNDVIi

(2)

LSTmin (NDVIi ) = a′ + b′NDVIi

(3)

where LSTmax (NDVIi ) and LSTmin (NDVIi ) are the warm and cold edges, i.e., the maximum and minimum LST values of the pixels that have the same NDVI value in a study region, respectively. a, b, a′ and b′are coefficients to be determined according to the scatter plots of the LST and NDVI. An important issue in using the VTCI drought monitoring is to determine the warm and cold edges. The warm edge was determined using the multiyear maximum value composite (MVC) LST and NDVI products at 10-day intervals and the cold edge was determined using the multiyear maximum-minimum value composite LST products and the MVC NDVI products at 10-day intervals (Sun et al., 2008). The remotely sensed data products used in this study are MODIS LST data product with a spatial resolution of 1 km (MYD11A1, MODIS/Terra Land Surface Temperature/Emissivity Daily L3 Global 1 km SIN Grid V006) for tiles h26v05 and h27v05 (h for horizontal and v for vertical) and MODIS surface reflectance data product with a spatial resolution of 1 km (MYD09GA, MODIS/Terra Surface Reflectance Daily L2G Global 1 km SIN Grid V006) for tiles h26v05 and h27v05 during the main winter wheat growth period from March to May of 2011–2018. Then,

2.2.2. Remotely sensed LAI The annual time series LAI of the years 2011–2018 were generated using 4-day MODIS LAI data product with a spatial resolution of 500 m (MCD15A3H, MODIS/Terra + Aqua Leaf Area Index/FPAR 4-Day L4 Global 500 m SIN Grid V006) for tiles h26v05 and h27v05. Compared with MOD15A2 (MODIS/Terra Leaf Area Index/FPAR 8-Day L4 Global 1000 m SIN Grid V006) and MYD15A2 (MODIS/Aqua Leaf Area Index/ FPAR 8-Day L4 Global 1000 m SIN Grid V006) products, the product of MCD15A3H has higher temporal resolution, which is more advantageous to monitor the growth and phenology of crops. However, due to presence of clouds and seasonal snow, there are gaps and noise in the MODIS LAI time series data that may result in discontinuities in the ability to monitor the spatial and temporal status of vegetation and 3

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

growing areas of counties (districts) in the Guanzhong Plain were obtained. The unit area winter wheat yield data of the counties in the Guanzhong Plain from 2011 to 2017 were recorded in the Shaanxi Rural Yearbooks. As Tongchuan City is located in the transition zone between the Guanzhong Plain and the Loess Plateau in northern Shaanxi, the area of winter wheat is relatively small and mainly distributed in the Weibei Tableland in the north. Therefore, the winter wheat yield of the counties, including the other (four) remaining cities in the Guanzhong Plain from 2011 to 2017, were selected. 2.3. Methods 2.3.1. BP neural network The BP (back propagation) neural network was proposed by Rumelhart et al. (1986), and is the most widely used algorithm for supervised learning with multi-layered feed-forward networks. Its basic idea is to revise the weights and thresholds of the network by back propagation to minimize the error between the actual output value and the expected output value. Neural networks with at least one hidden layer are necessary and sufficient for arbitrary nonlinear function approximation. In practice, neural networks with one or two hidden layers, that is, three-layer or four-layer perceptrons (including input and output layers) are commonly used. The topology of BP neural network with one hidden layer is shown in Fig. 3. There are two processes in the BP neural network learning: the forward propagation of the input signal and the reverse propagation of the error. In the forward propagation, the input signal acts on the output node through the hidden layer to generate an output signal. The state of the neuron in each layer only affects the state of the neuron in the next layer. If the actual output does not match the expected output, the error is reversed. The error back propagation is to pass the output error back to the input layer through the hidden layer, and to minimize the error signal by modifying the weights of each layer of neurons. The processes can be described as the following steps:

Fig. 2. The Savitzky-Golay filtered time series LAI of winter wheat in a survey site of Fufeng County.

thereby reduce accuracy (Yang et al., 2007; Yuan et al., 2011). To resolve this, the upper envelope Savitzky-Golay filter was applied to smooth out in time series MODIS LAI, which was proposed by Savitzky and Golay (1964) and applied in the reconstitution of the original time series MODIS LAI pixel by pixel. The maximum value of the LAI over the 10-day period was taken as the 10-day time series MODIS LAI. Then, the 10-day time series MODIS LAI was generated to composite the LAI images at the four winter wheat growth and development stages using averaging methods. An example of a survey site from day of year 61 to day of year 161 in 2017 illustrated the effect of filtering, as shown in Fig. 2. The Savitzky-Golay filtered time series LAI removed the sudden drops, and reflected the growth status and temporal variation of the crops, which was consistent with the characteristics of crop growth conditions and improved the quality of the original LAI images significantly (Xun et al., 2018). In order to make LAI and VTCI with the same spatial resolution, the SG filtered time series LAI is doubled, pixel by pixel to generate LAI with 1 km resolution. 2.2.3. Extraction of winter wheat planting area and wheat yield data In this paper, the MODIS land cover type product (MCD12Q1, MODIS/Terra + Aqua Land Cover Type Yearly L3 Global 500 m SIN Grid V051) was used for the extraction of winter wheat planting area. According to the classification scheme of International GeosphereBiosphere Program (IGBP), which were superimposed with the administrative boundary vector map of the study area, and the wheat

(1) Network initialization. There are m input neurons, n hidden neurons, and one output neuron. The first step in training is to initialize the weight parameters w, and small random values are usually suggested. wjk (j = 1, 2, …, m; k = 1, 2, …, n) in Fig. 3 represents the connection weight between the jth node in the input layer and

Fig. 3. Structure of the back propagation neural network. 4

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

network training error. The formula for calculating the error is the same as (7) and can effectively reduce verification times and thus find the optimal number of nodes in the hidden layer as quickly as possible. For example, when the number of nodes in the input layer is 4, the optimal number of nodes in the hidden layer is 9 and when the number of nodes in the input layer is 8, the optimal number of nodes in the hidden layer is 11.

the kth node in the hidden layer, hk is the output of the kth node in the hidden layer, wk represents the weighted value between the kth node in the hidden layer and the output layer, and y is the neuron’s output in the output layer. The computation is given by: m

hk = f

⎛ ⎞ w a − θk ⎜ ∑ jk j ⎟ ⎝ j=1 ⎠

(4)

m

⎞ ⎛ y = f ⎜ ∑ wn hn − θ⎟ ⎠ ⎝n=1

2.3.2. IPSO-BP neural network Particle swarm optimization (PSO) is a population-based stochastic optimization technique inspired by swarm intelligence theory proposed by Kennedy and Eberhart (1995) to achieve optimal process through collective collaboration between bird flocking. Compared to other intelligent algorithms, PSO algorithms have simple structure and less parameters and are easy to describe and implement (Neri et al., 2013; Tan et al., 2002; Zitzler et al., 2003). However, the standard PSO algorithm also has shortcomings such as premature convergence and bad local searching ability similar to other intelligent algorithms (Del Valle et al., 2008; Tsekouras and Tsimikas, 2013; Nickabadi et al., 2011). For example, in the optimization of complex problems in high-dimension, the population may have accumulated to a certain point of stagnation without finding the global optimization point and forming premature convergence. Meanwhile, in the search process of the PSO algorithm, the convergence speed is obviously slow when the particle is approaching or entering the most advantageous region. That is, in the later period of particle optimization, the search ability is poor. In order to prevent standard PSO from falling in a local optimum, improved PSO algorithms have been subsequently developed, among which improvements have been made to the inertia weight in PSO algorithms (Shi and Eberhart, 1998). The inertia weight w is a parameter that controls the impact of previous velocity on the current one, its size represents how fast the particle inherits its parent particle so that controls the algorithm convergence speed. Thus proper control of the inertia weight is very important to find the optimum solution efficiently. The concept of a nonlinear decreasing inertia weight was developed by Chatterjee et al. (2006) to better control local exploitation and global exploration, which can be described as:

(5)

where θk is the bias of the kth node in the hidden layer, and θ is the bias of the neuron in the output layer. Biases are assigned with random values between 0 and 1 before the forward propagation of the operating signal. The most commonly–used hidden neuron activation function is the sigmoid given by:

f (x ) =

1 1 + e−x

(6)

(2) Error calculation. The neural network performance is evaluated by computing the difference between actual neural network outputs and excepted outputs for all the training samples. The difference, also known as the error, is quantified by:

E = (y − y')2 /2

(7)

where y' presents the excepted output of the network. (3) Weights and bias updating. According to the definition of E in (7), wk, wjk, θk and θ are updated along the negative direction of the gradient E, until E becomes small enough.

wk (t − 1) = wk (t ) + ηδyt + α [wk (t ) − wk (t − 1) ]

(8)

wjk (t − 1) = wjk (t ) + ηδk yj + α [wjk (t ) − wjk (t − 1) ]

(9)

θ(t − 1) = θ(t ) + ηδ + α [θ(t ) − θ(t − 1) ]

(10)

θk (t − 1) = θk (t ) + ηδk + α [θk (t ) − θk (t − 1) ]

(11)

t π w = wmax − (wmax − wmin) × tan ⎛ × ⎞ t 4⎠ ⎝ max

Here, the parameter η is called the learning rate, α is called the momentum constant, and t is called the iterations. In this study, three schemes of input variables, that is, single-index VTCI or single-index LAI or double-index (VTCI and LAI) at the four growth stages in 25 counties from 2012 to 2016 were selected as inputs of the input layer. Among them, the VTCIs and LAIs from 2012 to 2015 were used as the training samples and the VTCIs and LAIs in 2016 were used as the testing samples. The yield data were converted to 0–1 by using the Min-Max normalization transformer, and then the normalized yield data during 2012 to 2016 of each county were used as the expected output vector. Similarly, the yield data from 2012 to 2015 were used as training data, and the yield data in 2016 were used as testing data. Therefore, when single-index VTCI or single-index LAI were selected as inputs of the input layer, the number of nodes in the input layer is 4; when VTCI and LAI were selected as inputs of the input layer, the number of nodes in the input layer is 8. No matter how many input nodes, the number of nodes in the output layer is 1. At present, the number of hidden nodes is comparatively more difficult to determine. Thus, to compare the number of different neurons, to determine an appropriate number of neurons in the hidden layer an empirical formula is commonly used:

n=

p+m +b

⎜

⎟

(13)

where wmax and wmin respectively represent the maximum and minimum values of inertia weight w, t is the number of iteration times, and tmax is the maximum number of iteration steps. Under the circumstances, the inertia coefficient determines the search step length. When t is smaller, it is helpful for global search. When t is larger, it is good for local search. Therefore, the algorithm can flexibly adjust the global search and local search ability. In addition, in order to enable the PSO algorithm to expand the shrinking population search during the iterative process, the particles can jump out of the previously searched optimal position, and then search in a larger space, which could maintain the diversity of the population. The mutation operation is introduced into the PSO algorithm, that is, the particles are initialized with a certain probability after each update of the particles. The IPSO-BP algorithm optimizes the initial connection weights and thresholds of the BP neural network by using the IPSO algorithm. When the algorithm terminates, the points around the global optimal position can be found. On this basis, a local search is carried out from the point by using the BP algorithm with strong local optimization ability, which is performed to achieve the training goal of the network. The IPSO algorithm combines the advantages of two algorithms (PSO and BP) to improve the deficiencies of conventional neural networks, which can be described as follows:

(12)

where b is a constant, and the value ranges from 1 to 10. Therefore, the number of nodes in the hidden layer distribution interval of n is [ 5 + 1, 5 + 10 ] or [4, 13]. By setting the loop of the number of hidden layer nodes when the training samples are the same, the optimal number of nodes in the hidden layer were selected according to the

(1) Determine the dimension D of the particle swarm.

D = nh + no + ni × nh + nh × no 5

(14)

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

stage, which require a lot of water and nutrients, and are key stages for winter wheat yield formation. Therefore, if drought occurs during these two growth stages, it would have a greater impact on winter wheat yield. The turning green stage would have a stronger tolerance to water deficit, and the effect of drought on yield is relatively small. This is consistent with previous studies, such as Li et al. (2010) analyzed the correlation between soil moisture and the field-measured yields at the different growth stages of winter wheat and concluded that irrigation was more efficient for producing high grain yield at the jointing stage, followed by the heading-filling, dough and turning green stages. The effect of LAI on the yield of winter wheat is larger at the heading-filling stage and the dough stage. Because the heading-filling stage is the reproductive growth stage of winter wheat, which mainly determines the grain weight of winter wheat. The growth of winter wheat at the flowering and fertilization stage has a great influence on the seed setting rate. Moreover, the yield formation of winter wheat is closely related to the accumulation and distribution of dry matter. In addition, the normal growth of winter wheat green leaves is one of the important sources of assimilates especially from the flowering stage to the dough stage. Therefore, as one of the important population parameters, the LAI of winter wheat from flowering to dough stage, is important for the accumulation of dry matter and the yield ahead of the harvest. Xie et al. (2017) and Huang et al. (2015) demonstrated that LAI at the heading-filling stage played the most critical role in wheat yield estimation, which is similar to our assessment of the weights presented above. As shown in Table 2, the VTCI at the jointing stage and the LAI at the heading-filling stage and dough stage were more closely related to the yield of winter wheat, which were the same as the contribution rates of the single index at each growth stage. This suggested that the I may combine the useful information of VTCI and LAI, which probably reflect the growth of winter wheat comprehensively. The weight values of VTCI and LAI are different at different growth stages, and the weight values of VTCI are higher than LAI at the jointing stage. This is because the jointing stage of winter wheat is a key stage of spikelet (defined as a small or secondary spike where each spikelet has one or more florets) differentiation. The drought will cause nutrient deficiency and at winter wheat jointing stage. The normal development of spikelet is affected, resulting in smaller spikes at the end of growth period, fewer grains per panicle and serious effects on yield. Therefore, the growth of winter wheat at the jointing stage is more sensitive to water stress. The weight values of LAI are higher than VTCI at the dough stage, which suggested that LAI was a key parameter contributing to winter wheat yield at the dough stage. Regardless of whether the BP neural network or the IPSOBP neural network were used, the proportions of VTCI and LAI were between 40% and 60%, indicating that VTCI and LAI both play an important role in winter wheat growth period and are closely related to winter wheat yield. The constructed index (I) takes into account the effects of VTCI and LAI on the growth and yield of winter wheat, and is more sensitive to moisture content and temperature relative to the VTCI and LAI, which are able to monitor winter wheat growth condition more accurate. In summary, these results indicated that I is more suitable index for monitoring winter wheat growth at the main growth stage than based on the VTCI or LAI alone.

where ni is the number of neurons in the input layer, nh is the number of neurons in the hidden layer, and no is the number of neurons in the output layer. (2) Set the fitness function of the particle swarm. In this paper, the mean square error (MSE) of BP neural network was used as the fitness function.

E=

1 m (Σk = 1 (y0 − y¯0 )2) m

(15)

where yo is the excepted output of the oth node in the output layer, y¯o is the actual output, and m is the number of samples. (3) Update the weights and thresholds. The IPSO algorithm was used to optimize the weights and thresholds of the BP neural network between adjacent layers. Then, the optimized weights and thresholds are used as the initial connection weights and thresholds of the BP algorithm to be trained in the network. Finally, the weights and thresholds are adjusted according to the training of the BP algorithm until the mean square error is less than e, which e is a predetermined expected accuracy. 2.3.3. Calculation of weights and establishment of integrated index A technique proposed by Garson (1991) is used to calculate the relative importance of the input variables by examining the connection weights of the trained BP neural network and IPSO-BP neural network respectively. The calculation formula is:

|wij ∗ wjk |

n

Si =

∑j = 1

Ri =

Si n ∑i = 0 Si

m

∑i = 1 |wij ∗ wjk |

(16)

(17)

where Si represents the absolute contribution margin of the ith input to the kth output and Ri represents the relative contribution margin, which the ratio of the absolute contribution margin of the ith input to the absolute contribution margin of all inputs. This study calculated the absolute contribution of the single index VTCI, the single index LAI, and the double index VTCI and LAI at each growth stage to the yield of winter wheat, and then the weight values of single index VTCI, single index LAI and double index VTCI and LAI at each growth stage were determined by formula (17). On this basis, the winter wheat of the weighted VTCIs (V), weighted LAIs (L) and integrated index (I) were established:

V = w11 × V1 + w21 × V2 + w31 × V3 + w41 × V4

(18)

L = w12 × L1 + w22 × L2 + w32 × L3 + w42 × L4

(19)

I = w11 × V1 + w21 × V2 + w31 × V3 + w41 × V4 + w12 × L1 + w22 × L2 + w32 × L3 + w42 × L4

(20)

where V1-V4 represent the VTCIs of counties at the turning green, jointing, heading-filling, and dough stages; L1 − L4 represent the LAIs of counties at the four growth stages; w11 − w41 represent the weight coefficients of the VTCIs at the four growth stages; and w12 − w42 represent the weight coefficients of the LAIs at the four growth stages.

3.2. Establishment of yield estimation models 3. Results and analysis The linear regression models between different variables (L, V and I) and the winter wheat yield at the county level from 2012 to 2016 were established. The results presented in Table 3 clearly suggested that the regression models based on the IPSO-BP neural network were better than the regression models based on the BP neural network. Because the factors affecting winter wheat yields are complex and multiple, and that using the IPSO-BP neural network can better deal with nonlinear problems and therefore assign the importance of different input variables to output variables more reasonably. There were good correlations

3.1. Weights of each growth stage of winter wheat Table 1 presented the weight values of VTCI at the jointing stage and the heading-filling stage are higher than those of VTCI at the turning green stage and the dough stage based on the BP neural network and the IPSO-BP neural network. The reason is that the jointing stage and the heading-filling stage of winter wheat are in the transition stage from the vegetative growth stage to the reproductive growth 6

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

Table 1 The weights of single index VTCI and single index LAI at the four growth stages of winter wheat in the Guanzhong Plain. Method

BP neural network

IPSO-BP neural network

Growth stage Index

Turning green

Jointing

Heading-filling

Dough

Turning green

Jointing

Heading-filling

Dough

VTCI LAI

0.16 0.21

0.35 0.20

0.25 0.29

0.24 0.30

0.15 0.25

0.34 0.21

0.28 0.26

0.23 0.28

between I and winter wheat yields in the Guanzhong Plain, where the values of the significance test (P) were less than 0.001, indicating that the correlation between I and winter wheat yields reach the very significant level. As a result, the regression models between I and winter wheat yields are better than V or L, which indicated that the single index does not reflect the crop growth information comprehensively. The I combined the useful information of VTCI and LAI in estimating crop yields effectively, at the same time making up for the shortcomings of different variables in estimating yields. Neural network parameters were adjusted with different combination of inputs, leading to different training results. We selected input data for different combinations but the same amount from 2011 to 2015, from 2012 to 2016 and from 2013 to 2017, which further verify the applicability of I in winter wheat yield estimation and performance of the IPSO-BP neural network with different training samples. In addition, due to the management practice and scientific and technological improvement, models for monitoring crop growth condition and estimating crop yields were developed using 5 years of data. The estimated yield models of the whole study area from 2011 to 2015, from 2012 to 2016 and from 2013 to 2017 were established by using the BP neural network and IPSO-BP neural network (Table 4). The results showed that the correlations between I and winter wheat yield when using IPSO-BP neural network were higher than using BP neural network, and these models passed the significance test. In conclusion, the I constructed by using the IPSO-BP neural network was a more viable index to comprehensively estimate the winter wheat yields with more accuracy.

Table 3 Linear regression models of the relation between the different indices and yields of winter wheat in the Guanzhong Plain from 2012 to 2016 based on the BP neural network and the IPSO-BP neural network. Methods

Index

Linear regression model

R2

P

BP

V L I

Y = 4813.2 V + 1537.0 Y = 6529.2L + 3007.3 Y = 6495.5I + 1570.6

0.240 0.307 0.310

P < 0.001 P < 0.001 P < 0.001

IPSO-BP

V L I

Y = 4898.8 V + 1477.8 Y = 6654.3L + 3011.6 Y = 7130.7I + 1818.9

0.244 0.311 0.342

P < 0.001 P < 0.001 P < 0.001

Table 4 Linear regression models of the relation between the integrated index and yields of winter wheat on the whole study area in different years. Year

Method

Linear regression model

R2

P

2011–2015

BP IPSO-BP

Y = 7440.0I + 1474.5 Y = 7638.3I + 1421.3

0.390 0.400

P < 0.001 P < 0.001

2012–2016

BP IPSO-BP

Y = 6495.5I + 1570.6 Y = 7130.7I + 1818.9

0.310 0.342

P < 0.001 P < 0.001

2013–2017

BP IPSO-BP

Y = 7147.3I + 1793.4 Y = 7473.0I + 1936.4

0.349 0.360

P < 0.001 P < 0.001

County, Fufeng County, Mei County and Pucheng County, the relative errors are greater than 15%. These results suggest that the validated models are indicative of the models capability to estimate yields, since the winter wheat yield data of 2017 are independent of the data used for model development. To validate the proposed method for monitoring regional winter wheat growth and the yield results, linear regression analysis was carried out based on the actual yields and estimated yields of winter wheat in 25 counties from 2011 to 2017 using the optimal estimation model between I and winter wheat yield of 2012–2016. The results showed that there was a significant positive correlation between the estimated yields and the actual yields (R2 = 0.435, P < 0.001, RMSE = 673.41 kg/ha). The models reported here are appropriate for estimating winter wheat yields in the Guanzhong Plain. In summary, the yield estimation models based on the I were more in accordance with the actual growth conditions, which can provide comprehensive winter wheat growth-related information with good monitoring accuracy and early prediction.

3.3. Validation of the winter wheat yield The estimation model between I and winter wheat yield of 2012–2016 based on IPSO-BP neural network (R2 = 0.342, P < 0.001) was used to estimate the winter wheat yield of each county in the Guanzhong Plain in 2017, and the relative errors between the estimated yields and actual yields of 25 counties were analyzed. The results showed that the relative errors in 17 counties are less than 10%, of which the actual yield of Chunhua County is 4297.0 kg/ha and the estimated yield is 4299.7 kg/ha, the relative error is the smallest (0.57%) among all counties. Most of these 17 counties are concentrated in the western and central parts of the Guanzhong plain, which suggests that the established model is more reliable to moderate and high-yield areas. The relative errors of the 3 counties, Changan County, Zhouzhi County and Liquan County, were lower than 15%, and the values were 11.9%, 12.0% and 11.2% respectively. In Yanliang County, Gaoling

Table 2 The weights of double index VTCI and LAI at the four growth stages of winter wheat in the Guanzhong Plain by using the BP neural network and the IPSO-BP neural network respectively. VTCI

LAI

Growth stage Method

Turning green

Jointing

Heading-filling

Dough

Turning green

Jointing

Heading-filling

Dough

BP IPSO-BP

0.10 0.06

0.20 0.13

0.18 0.11

0.10 0.10

0.09 0.16

0.08 0.12

0.10 0.16

0.15 0.16

7

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

Fig. 4. Estimates of winter wheat yields (kg/ha) in the Guanzhong Plain from 2011 to 2018.

winter wheat yields estimation results, the distribution characteristics of winter wheat yields in the Guanzhong Plain were analyzed. For the spatial distribution of winter wheat yields, the yields are highest in the west part of the plain, followed by the central part, and in the east part the yields are lowest. These results were consistent with the fact that

3.4. Regional winter wheat yield estimation Based on the estimation models between I and winter wheat yields of 2012–2016, the yields of winter wheat in the study area from 2011 to 2018 were estimated pixel by pixel (Fig. 4). Through the regional 8

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

5. Conclusions

central part and west part are the main grain producing regions in the Guanzhong Plain. In the central part of Guanzhong Plain, the average winter wheat yields ranged from 4000 kg/ha to 5400 kg/ha, of which the highest yields in Fuping County in 2015 were 5353.2 kg/ha and the lowest yields in Gaoling county in 2013 were 4055.9 kg/ha. The average yields in the central part were 4625.8 kg/ha. In the west part, the average yields were distributed at 3800 kg/ha-6000 kg/ha, among which Fengxiang County in 2015 had the highest yield per unit, which was 5943.9 kg/ha, and in 2013, Qianxian County had the lowest yield (3860.7 kg/ha). The average yields in the western part were 4814.6 kg/ ha. In the east part the average winter wheat yields ranged from 3200 kg/ha to 4600 kg/ha, among which, the highest yields in Pucheng County in 2015 were 4583.6 kg/ha, and the lowest yields in Chengcheng County in 2013 were 3291.2 kg/ha. The average yields in the east were 3984.9 kg/ha. For the temporal characteristics of yields from 2011 to 2018, it can be seen from inter-annual variation that the estimated yields showed inter-annual fluctuations along with an increasing trend, of yield, on the whole. In 2013, due to the reduction in precipitation during the winter wheat growth season, the yields were the lowest. In 2015, due to abundant precipitation the yields of winter wheat were the highest, which was consistent with the field investigation results of actual yield in the Guanzhong Plain.

In this study, the single VTCI, the single LAI and double index VTCI and LAI at the turning green, jointing, heading-filling, and dough stages of winter wheat were derived from remotely sensed data and selected as model inputs. Then, the BP neural network and the IPSO-BP neural network were used to calculate the weight coefficients and thresholds of the VTCI and LAI at the four growth stages and to establish an integrated index, I, during the main growth period. The VTCI, LAI and I were selected as variables, which established estimation models between different variables and the winter wheat yield. The results showed that correlations between the I values and winter wheat yields for models based on the IPSO-BP neural network in the years from 2011 to 2017 at the county scale were higher, which can accurately reflect the winter wheat growth condition and yield information over the estimated values of winter wheat yield with the weighted VTCI or weighted LAI alone. Crop yield models were used to assess how well remote-sensing based neural networks are able to capture county-level geographic and inter-annual variability in crop yields in the Guanzhong Plain. The spatial distribution of winter wheat yields in the plain from 2011 to 2018 indicated that the yields from west part to east part are gradually decreasing. The inter-annual variation trend of the estimated yields was consistent with the records of rural yearbooks and the impacts of historical natural disasters, such as severe drought in 2013 resulted in low yields. Therefore, the estimation model was effective for estimating winter wheat yields.

4. Discussion The VTCI and LAI are two important indicators that can reflect the near real time information related crop water stress, crop growth condition and yield of winter wheat. In this study, the VTCIs and LAIs at each growth stage of winter wheat from 2012 to 2016 were used as variables to simulate the values of the I by calculating the weighting values of VTCI and LAI at each growth stage of winter wheat. The correlation between the I and the winter wheat yield was better than those use the VTCI or LAI alone, which greatly improved the accuracy of winter wheat yield estimation and can be used to estimate winter wheat yields in the Guanzhong Plain. We assumed that the start and end dates of winter wheat growth stages were the same from year to year and also from region to region. In reality, this might not be the case and may lead to a deviation between the winter wheat growth stage in non-normal years and the actual situation. This will, in turn, affect the winter wheat growth monitoring and the accuracy of the yield estimation, especially in 2013 where drought was severe. Future research will focus on investigating time series data products generated by fusing MODIS data and Sentinel series satellite data, and combine with the field survey data to achieve accurate identification of winter wheat phenophase. Our results showed that as a result of likely improvements in crop cultivation practices and associated technological advances, crop production per unit area (yield) has been on the upward trend. For the time series estimation of yields, influence of such improvements also should be considered (Kuwata and Shibasaki, 2015). We show that the yield estimation model based on the IPSO-BP neural network is better than that based on the BP neural network. IPSO-BP neural network can better reflect the nonlinear search process, but it also has the disadvantage of falling into the local minimum value in the late search. How to combine the BP algorithm with other evolutionary algorithms and optimize the network structure to derive a more characteristic and practical value of the hybrid algorithm remain to be further studied. Due to the small number of samples, the training neural network model has been affected to some extent, resulting in lower precision of yield prediction in some counties. In the future, the neural network can be retrained by increasing the number of training samples to develop a more accurate estimation model. Furthermore, additional prior knowledge of other factors affecting crop yield need to be considered. These factors can be investigated to train algorithms using reinforcement learning methods to extract important features for estimating crop yields.

CRediT authorship contribution statement Huiren Tian: Conceptualization, Methodology, Data curation, Writing original draft. Pengxin Wang: Supervision, Conceptualization, Methodology, Data curation, Resources, Writing review & editing, Funding acquisition. Kevin Tansey: Resources, Writing - review & editing, Funding acquisition. Shuyu Zhang: Resources, Data curation. Jingqi Zhang: Data curation. Hongmei Li: Resources, Data curation. Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgments This work was supported by the National Natural Science Foundation of China under Grants 41871336 and 41811530303. This work was supported by a UK Science & Technology Facilities Council (STFC) Agri-Tech in China Newton Network+ (ATCNN) grant administered through Rothamsted Research. The work was further supported by a Royal Society-Newton Mobility Grant (UK). References Becker-Reshef, I., Vermote, E., Lindeman, M., Justice, C., 2010. A generalized regressionbased model for forecasting winter wheat yields in Kansas and Ukraine using MODIS data. Remote Sens. Environ. 114 (6), 1312–1323. https://doi.org/10.1016/j.rse. 2010.01.010. Bolton, D.K., Friedl, M.A., 2013. Forecasting crop yield using remotely sensed vegetation indices and crop phenology metrics. Agric. Forest Meteorol. 173, 74–84. https://doi. org/10.1016/j.agrformet.2013.01.007. Camci, E., Kripalani, D.R., Ma, L., Kayacan, E., Khanesar, M.H., 2018. An aerial robot for rice farm quality inspection with type-2 fuzzy neural networks tuned by particle swarm optimization-sliding mode control hybrid algorithm. Swarm Evol. Comput. 41, 1–8. https://doi.org/10.1016/j.swevo.2017.10.003. Chatterjee, A., Siarry, P., 2006. Nonlinear inertia weight variation for dynamic adaptation in particle swarm optimization. Comput. Oper. Res. 33 (3), 859–871. https://doi.org/ 10.1016/j.cor.2004.08.012.

9

Computers and Electronics in Agriculture 169 (2020) 105180

H. Tian, et al.

Papageorgiou, E.I., Aggelopoulou, K.D., Gemtos, T.A., Nanos, G.D., 2013. Yield prediction in apples using fuzzy cognitive map learning approach. Comput. Electron. Agric. 91, 19–29. https://doi.org/10.1016/j.compag.2012.11.008. Pantazi, X.E., Moshou, D., Alexandridis, T., Whetton, R.L., Mouazen, A.M., 2016. Wheat yield prediction using machine learning and advanced sensing techniques. Comput. Electron. Agric. 121, 57–65. https://doi.org/10.1016/j.compag.2015.11.018. Qi, R., Ma, Y., Hu, B., Reffye, P.D., Cournède, Paul-Henry, 2010. Optimization of source–sink dynamics in plant growth for ideotype breeding: a case study on maize. Comput. Electron. Agric. 71 (1), 96–105. https://doi.org/10.1016/j.compag.2009. 12.008. Prasad, A.K., Chai, L., Singh, R.P., Kafatos, M., 2006. Crop yield estimation model for Iowa using remote sensing and surface parameters. Int. J. Appl. Earth Obs. 8 (1), 26–33. https://doi.org/10.1016/j.jag.2005.06.002. Rumelhart, D.E., Hinton, G.E., Williams, R.J., 1986. Learning representations by backpropagating errors. Nature 323 (6088), 533–536. https://doi.org/10.1038/ 323533a0. Savitzky, A., Golay, M.J.E., 1964. Smoothing and differentiation of data by simplified least squares procedures. Anal. Chem. 36 (8), 1627–1639. https://doi.org/10.1021/ ac60214a047. Shi, Y., Eberhart, R.C., 1998. A modified particle swarm optimization. In: Proc. IEEE Congr. Evol. Comput. pp. 69–73. https://doi.org/10.1109/ICEC.1998.699146. Singh, V.K., Singh, B.P., Kisi, O., Kushwaha, D.P., 2018. Spatial and multi-depth temporal soil temperature assessment by assimilating satellite imagery, artificial intelligence and regression based models in arid area. Comput. Electron. Agric. 150, 205–219. https://doi.org/10.1016/j.compag.2018.04.019. Sun, W., Wang, P., Zhang, S., Zhu, D., Liu, J., Chen, J., Yang, H., 2008. Using the vegetation temperature condition index for time series drought occurrence monitoring in the Guanzhong plain, PR China. Int. J. Remote Sens. 29 (17–18), 5133–5144. https:// doi.org/10.1080/01431160802036557. Tan, K.C., Lee, T.H., Khor, E.F., 2002. Evolutionary algorithms for multi-objective optimization: performance assessments and comparisons. Artif. Intell. Rev. 17, 251–290. https://doi.org/10.1023/A:1015516501242. Tian, M., Wang, P., Khan, J., 2016. Drought forecasting with vegetation temperature condition index using ARIMA models in the Guanzhong Plain. Remote Sens. 8, 690. https://doi.org/10.3390/rs8090690. Tien Bui, D., Bui, Q.T., Nguyen, Q.P., Pradhan, B., Nampak, H., Trinh, P.T., 2017. A hybrid artificial intelligence approach using GIS-based neural-fuzzy inference system and particle swarm optimization for forest fire susceptibility modeling at a tropical area. Agric. Forest Meteorol. 233, 32–44. https://doi.org/10.1016/j.agrformet.2016. 11.002. Tsekouras, G.E., Tsimikas, J., 2013. On training RBF neural networks using input–output fuzzy clustering and particle swarm optimization. Fuzzy Set. Syst. 221, 65–89. https://doi.org/10.1016/j.fss.2012.10.004. USDA NASS, 2009. Limited use for crop condition and crop yield. http://www.nass.usda. gov/surveys/remotely sensed data crop yield/index.asp. Wan, Z., Wang, P., Li, X., 2004. Using MODIS Land surface temperature and normalized difference vegetation index products for monitoring drought in the southern Great Plains, USA. Int. J. Remote Sens. 25, 61–72. https://doi.org/10.1080/ 0143116031000115328. Wang, L., Wang, P., Li, L., Xun, L., Kong, Q., Liang, S., 2018. Developing an integrated indicator for monitoring maize growth condition using remotely sensed vegetation temperature condition index and leaf area index. Comput. Electron. Agric. 152, 340–349. https://doi.org/10.1016/j.compag.2018.07.026. Wang, P., Gong, J., Li, X., 2001. Vegetation temperature condition index and its application for drought monitoring. Geomat. Inform. Sci. Wuhan Univ. 26 (5), 412–418. https://doi.org/10.13203/j.whugis2001.05.007. (in Chinese with English abstract). Wang, P., Sun, H., Xie, Y., Wang, L., Zhang, S., Li, L., 2016. Winter wheat yield estimation based on particle filter assimilation algorithm and remotely sensed LAI and VTCI. Trans. Chin. Soc. Agric. Mach. 47 (4), 248–256. https://doi.org/10.6041/j.issn.10001298.2016.04.033. (in Chinese with English abstract). Xie, Y., Wang, P., Bai, X., Khan, J., Zhang, S., Li, L., Wang, L., 2017. Assimilation of the leaf area index and vegetation temperature condition index for winter wheat yield estimation using Landsat imagery and the CERES-Wheat model. Agric. Forest Meteorol. 246, 194–206. https://doi.org/10.1016/j.agrformet.2017.06.015. Xun, L., Wang, P., Li, L., Wang, L., Kong, Q., 2018. Identifying crop planting areas using Fourier-transformed feature of time series MODIS leaf area index and sparse-representation-based classification in the North China Plain. Int. J. Remote Sens. 1–19. https://doi.org/10.1080/01431161.2018.1492181. Yang, P., Shibasaki, R., Wu, W., Zhou, Q., Chen, Z., Zha, Y., Shi, Y., Tang, H., 2007. Evaluation of MODIS land cover and LAI products in cropland of North China plain using in situ measurements and Landsat TM images. IEEE Trans. Geosci. Remote. 45 (10), 3087–3097. https://doi.org/10.1109/TGRS.2007.902426. Yang, Q., Shi, L., Zha, Y., Zhu, P., 2019. Deep convolutional neural networks for rice grain yield estimation at the ripening stage using UAV-based remotely sensed images. Field Crop. Res. 235, 142–153. https://doi.org/10.1016/j.fcr.2019.02.022. You, J., Li, X., Low, M., Lobell, D., Ermon, S., 2017. Deep Gaussian process for crop yield prediction based on remote sensing data. Proc. Thirty-First AAAI Conf. Artif. Intel. . https://aaai.org/ocs/index.php/AAAI/AAAI17/paper/view/14435. Yuan, H., Dai, Y., Xiao, Z., Ji, D., Wei, S., 2011. Reprocessing the MODIS leaf area index products for land surface and climate modelling. Remote Sens. Environ. 115 (5), 1171–1187. https://doi.org/10.1016/j.rse.2011.01.001. Zitzler, E., Thiele, L., Laumanns, M., Fonseca, C.M., TanFonseca, V.G.D., 2003. Performance assessment of multiobjective optimizers: an analysis and review. IEEE Trans. Evolut. Comput. 7 (2), 117–132. https://doi.org/10.1109/TEVC.2003. 810758.

Chen, S., Xu, Z., Tang, Y., Liu, S., 2014. An improved particle swarm optimization algorithm based on centroid and exponential inertial weight. J. Netw. Comput. Appl. 2014 (1), 1–14. https://doi.org/10.1155/2014/976486. Chlingaryan, A., Sukkarieh, S., Whelan, B., 2018. Machine learning approaches for crop yield prediction and nitrogen status estimation in precision agriculture: a review. Comput. Electron. Agric. 151, 61–69. https://doi.org/10.1016/j.compag.2018.05. 012. Darvishzadeha, R., Schlerf, M., Atzberger, C., Corsi, F., Cho, M., 2008. LAI and chlorophyll estimation for a heterogeneous grassland using hyperspectral measurements. ISPRS J. Photogramm. 63 (4), 409–426. https://doi.org/10.1016/j.isprsjprs.2008.01. 001. Del Valle, Y., Venayagamoorthy, G.K., Mohagheghi, S., Hernandez, J.C., Harley, R.G., 2008. Particle swarm optimization: basic concepts, variants and applications in power systems. IEEE Trans. Evolut. Comput. 12, 171–195. https://doi.org/10.1109/ TEVC.2007.896686. Fortin, J.G., Anctil, F., Parent, L., Bolinder, M.A., 2011. Site-specific early season potato yield forecast by neural network in Eastern Canada. Precis. Agric. 12 (6), 905–923. https://doi.org/10.1007/s11119-011-9233-6. Fu, Q., Wang, Z., Jiang, Q., 2010. Delineating soil nutrient management zones based on fuzzy clustering optimized by PSO. Math. Comput. Model. 51 (11–12), 1299–1305. https://doi.org/10.1016/j.mcm.2009.10.034. Garson, G.D., 1991. Interpreting neural-network connection weights. AI Expert 6, 47–51. https://doi.org/10.1207/s15327752jpa8502. Gaso, D.V., Berger, A.G., Ciganda, V.S., 2019. Predicting wheat grain yield and spatial variability at field scale using a simple regression or a crop model in conjunction with Landsat images. Comput. Electron. Agric. 159, 75–83. https://doi.org/10.1016/j. compag.2019.02.026. Gillies, R.R., Carlson, T.N., 1995. Thermal remote sensing of surface soil water content with partial vegetation cover for incorporation into climate models. J. Appl. Meteorol. 34 (4), 745–756. https://doi.org/10.1175/1520-0450(1995) 0342.0.CO;2. Gillies, R.R., Kustas, W.P., Humes, K.S., 1997. A verification of the 'triangle' method for obtaining surface soil water content and energy fluxes from remote measurements of the normalized difference vegetation index (NDVI) and surface radiant temperature. Int. J. Remote Sens. 18 (15), 3145–3166. https://doi.org/10.1080/ 014311697217026. Habibi Davijani, M., Banihabib, M.E., Nadjafzadeh, Anvar A., Hashemi, S.R., 2016. Optimization model for the allocation of water resources based on the maximization of employment in the agriculture and industry sectors. J. Hydrol. 533, 430–438. https://doi.org/10.1016/j.jhydrol.2015.12.025. Huang, J., Tian, L., Liang, S., Ma, H., Becker-Reshef, I., Huang, Y., Su, W., Zhang, X., Zhu, D., Wu, W., 2015. Improving winter wheat yield estimation by assimilation of the leaf area index from Landsat TM and MODIS data into the WOFOST model. Agric. Forest Meteorol. 204, 106–121. https://doi.org/10.1016/j.agrformet. Johnson, D.M., 2016. A comprehensive assessment of the correlations between field crop yields and commonly used MODIS products. Int. J. Appl. Earth Obs. 52, 65–81. https://doi.org/10.1016/j.jag.2016.05.010. Kamilaris, A., Kartakoullis, A., Prenafeta-Boldú, Francesc, X., 2017. A review on the practice of big data analysis in agriculture. Comput. Electron. Agric. 143, 23–37. https://doi.org/10.1016/j.compag.2017.09.037. Kamilaris, A., Prenafeta-Boldu, F.X., 2018. Deep learning in agriculture: a survey. Comput. Electron. Agric. 147, 70–90. https://doi.org/10.1016/j.compag.2018.02. 016. Kennedy, J., Eberhart, R., 1995. Particle swarm optimization. In: Proc. IEEE Int. Conf. Neur. Netw. pp. 1942–1948. https://doi.org/10.1109/ICNN.1995.488968. Kitchen, N.R., Drummond, S.T., Lund, E.D., Sudduth, K.A., Buchleiter, G.W., 2003. Soil electrical conductivity and topography related to yield for three contrasting soil-crop systems. Agron. J. 95, 483–495. https://doi.org/10.2134/agronj2003.0483. Kuwata, K., Shibasaki, R., 2015. Estimating crop yields with deep learning and remotely sensed data. IEEE Int. Geosci. Remote Sens. Symp. 858–861. https://doi.org/10. 1109/IGARSS.2015.7325900. Li, Q., Dong, B., Qiao, Y., Liu, M., Zhang, J., 2010. Root growth, available soil water, and water-use efficiency of winter wheat under different irrigation regimes applied at different growth stages in North China. Agric. Water Manag. 97 (10), 1676–1682. https://doi.org/10.1016/j.agwat.2010.05.025. Li, Y., Wang, P., Liu, J., Zhang, S., Li, L., 2014. Evaluation of drought monitoring effects in the main growth and development stages of winter wheat using vegetation temperature condition index III-Impact evaluation of drought on wheat yield. Agric. Res. Arid Areas 69 (4), 697–702 (in Chinese with English abstract). Liang, L., Di, L., Zhang, L., Deng, M., Qin, Z., Zhao, S., Lin, H., 2015. Estimation of crop LAI using hyperspectral vegetation indices and a hybrid inversion method. Remote Sens. Environ. 165, 123–134. https://doi.org/10.1016/j.rse.2015.04.032. Lokers, R., Knapen, R., Janssen, S., Randen, Y.V., Jansen, J., 2016. Analysis of big data technologies for use in agro-environmental science. Environ. Modell. Softw. 84, 494–504. https://doi.org/10.1016/j.envsoft.2016.07.017. Miao, Y., Mulla, D.J., Robert, P.C., 2006. Identifying important factors influencing corn yield and grain quality variability using artificial neural networks. Precis. Agric. 7 (2), 117–135. https://doi.org/10.1007/s11119-006-9004-y. Mkhabela, M.S., Bullock, P., Raj, S., Wang, S., Yang, Y., 2011. Crop yield forecasting on the Canadian Prairies using MODIS NDVI data. Agric. Forest Meteorol. 151 (3), 385–393. https://doi.org/10.1016/j.agrformet.2010.11.012. Neri, F., Mininno, E., Lacca, G., 2013. Compact particle swarm optimization. Inf. Sci. 239, 96–121. https://doi.org/10.1016/j.ins.2013.03.026. Nickabadi, A., Ebadzadeh, M.M., Safabakhsh, R., 2011. A novel particle swarm optimization algorithm with adaptive inertia weight. Appl. Soft Comput. 11 (4), 3658–3670. https://doi.org/10.1016/j.asoc.2011.01.037.

10