Introduction

Wheat ranks second among the main cereals consumed in Mexico, accounting for about 10 % of caloric consumption (U.S. Department of Agriculture - Foreign Agricultural Service [USDA-FAS], 2019). However, production costs have increased in recent years, with fertilization being the most expensive cultural activity, accounting for approximately 30 % of the total production cost. The balance between wheat price, yields and production costs generates little profit margin for the producer (Retes-López et al., 2020).

To determine whether a crop is economically profitable, it is essential to estimate the expected yield in order to identify and establish farming practices that favor production. This requires rapid and precise methods, as well as non-destructive techniques to manage resources, such as the application of nitrogen fertilization, on which the increase in yield depends (Guo et al., 2020; Qiao et al., 2022; Shibayama et al., 2012).

Yield is associated with leaf chlorophyll content; therefore, the use of optical sensors allows determining different crop parameters, such as biomass, canopy fraction or reflectance indices (such as the greenness index) (Yue et al., 2021), at various scales and different parts of the crop. The SPAD502+ (Minolta Ltd, Japan) is a commercial chlorophyll meter, and its use has shown a significant positive correlation with grain yield (Monostori et al., 2016); however, its ability to detect nitrogen status varies with the stage of crop development (Rahman et al., 2020). On the other hand, the GreenSeeker (GS) (Trimble®, USA) is an optical sensor that measures crop canopy reflectance. The use of this sensor has generated good results in estimating grain yield at different crop growth stages (Ali et al., 2020; Kaur et al., 2018; Yegül et al., 2020).

Technological progress has enabled the development of new optical sensors to indirectly estimate crop yields. These include sensors mounted on aerial devices, such as satellites, manned aircraft and, recently, unmanned aerial vehicles (UAVs), also known as drones. UAVs are affordable aerial platforms with several maneuverability advantages. Optical sensors are installed on these devices, which allow generating spectral information without direct contact with objects in the area of interest (Du & Noguchi, 2016; Saravia et al., 2023; Tsouros et al., 2019).

Cameras are optical sensors that are classified as multispectral, hyperspectral and thermal, and allow determining vegetation indices (VI). These indices have been used to indirectly determine agronomic variables such as yield, leaf chlorophyll content, biomass, leaf area index, crop height and canopy cover, among others (dos Santos et al., 2021; Feng et al., 2020; Gilliot et al., 2021; Li et al., 2022).

Walsh et al. (2022) evaluated the accuracy of an optical sensor placed on a UAV to estimate the average yield of three wheat genotypes, and compared it with GS sensor readings. These authors report that GS explained up to 79 % of yield, and the normalized difference vegetative index (NDVI) derived from UAV images explained 67 % of this variable; they also reported a strong correlation between GS and UAV (R² = 85 %).

Technology continues to evolve; therefore, the use of optical sensors under different climatic, nutritional and operating conditions should continue to be explored. Considering the above, this research aimed to evaluate the performance of the GS and SPAD502+ sensors, and of a multispectral camera mounted on a UAV to estimate the yield of a wheat crop under different nitrogen levels and at two phenological stages close to flowering (thickened sheath and heading).

Materials and methods

Description of the experimental site

The experiment was established in the experimental field of the Colegio de Postgraduados Campus Montecillo (19° 27’ 36.0” N and 98° 54’ 00” W, at 2 250 m a. s. l.) (Figure 1), characterized by a temperate climate with rainfall in summer and the dry season in winter [(Cwo)(w)b(i’) (g)] (García, 2004), an average annual temperature of 15.2 °C and 650 mm of annual precipitation.

The Nana 2007 wheat variety was established in mid-February and harvested at the end of June 2019. This variety was released by the National Rainfed Wheat Program of the National Institute of Forestry, Agriculture and Livestock Research (INIFAP), and was obtained by hybridization and selection through the combination of mass methods and families (Villaseñor-Mir et al., 2014). Seed density was 100 kg∙ha^-1, and drip irrigation was applied. Plant care was carried out according to conventional practices in the area.

The experimental design was completely randomized with seven treatments (nitrogen dose [kg∙ha^-1]: T1 = 0, T2 = 40, T3 = 60, T4 = 80, T5 =100, T6 =140 and T7 = 180) and four replications (Figure 2). The experimental unit size was 4 × 10.5 m, which covered an area of 1 600 m². N was applied at two times: during sowing and at the final phenological stage of stem elongation.

Experimental data acquisition with ground-based optical sensors

At each phenological stage (thickened sheath and heading) a systematic sampling was carried out: three on one side, three in the center and three on the other side (Figure 3), resulting in nine readings for each experimental unit, both with GS and SPAD502+.

Measurements with the GS sensor were taken, on average, at 1 m above the crop canopy (Figure 4a), whereas the SPAD502+ readings were taken at a midpoint of the fully developed flag leaf (Figure 4b) as indicated by Yue et al. (2020). The recorded readings were averaged to determine the representative value of the experimental unit and, subsequently, the average value of each treatment was determined; this was done for each data set obtained with each sensor.

Multispectral image acquisition and processing

The flight schedule and dates were synchronized with the samplings carried out with the GS and SPAD502+ sensors. The UAV (multirotor 3DR X8+, 3D Robotics, USA) was fitted with a multispectral camera (Canon S110 NIR, bands: blue [400-495 nm], green [490-550 nm], near-infrared [680-760 nm]; Haghighattalab et al., 2016) to obtain images at 120 m above ground level, resulting in a spatial resolution of 3 cm∙pixel^-1. Flights were conducted at midday to avoid shadowing in the collected images. The sky was clear during all sessions.

The images were processed in Pix4D software (Pix4D SA, Switzerland), and orthomosaics of the surface reflectance were obtained for each channel of the multispectral camera, as well as a point cloud and a digital surface model. The processing of the imagery derived from UAVs is mainly based on the combination of photogrammetric and computer vision algorithms (Remondino et al., 2014). Berra and Peppa (2020) state that image processing can be summarized in three phases: 1) sparse point cloud reconstruction, 2) georeferencing and 3) dense point cloud reconstruction.

High-resolution images collected from UAVs contain too much surface information, which requires equipment with great ability to extract and analyze the information (Ye et al., 2023). The images may contain vegetation, shadows, soil and weeds, among other objects; therefore, it is necessary to segregate the information to analyze exclusively the object of interest. One method that allows this action is object-based image analysis (OBIA), which uses a set of pixels with similar characteristics, based on the texture, shape, spatial structure and other multidimensional characteristics of adjacent pixels (Filippi et al., 2022).

In the present study, the multiresolution segmentation algorithm included in eCognition software (Trimble Inc., USA) was used to classify vegetation pixels into objects. This tool is intended to reduce the bias contributed by other objects to the VI values for each experimental unit, especially in treatments with low nitrogen doses, which have a less developed canopy and do not have full cover. This pixel segregation has been shown to improve the accuracy of VI estimation by minimizing the influence of the soil (Duan et al., 2017), reaching accuracies of up to 90 % in the classification of herbaceous crop vegetation (Torres-Sánchez et al., 2015). At the end of the analysis, eCognition delivers the objects that it classified as vegetation in the orthoimage in vector format.

Vegetation indices (VI)

It is known that VIs are obtained from plant reflectance using the electromagnetic spectrum; however, in order to relate them to vigor and productivity, it is necessary to understand the optical properties of the leaves, mainly the role of chlorophyll, carotenoids and xanthophylls, as well as the mesophyll cells in light reflectance and absorbance (Taddeo et al., 2019).

The selection of the spectral vegetation index should be made from the temporal spectral variability approach, so that it can be correlated with the physiological variable of interest, such as grain yield. This variability occurs in wheat plants with discrete flowers, which turn green at the beginning of the reproductive stages, and yellow, or even brown, at the maturation stage (Sulik & Long, 2016).

This research used the green normalized difference vegetation index (GNDVI; Gitelson et al., 1996) and the blue normalized difference vegetation index (BNDVI; Wang et al., 2007), which were calculated as follows:

where NIR is the near-infrared wavelength, G is the green wavelength and B is the blue wavelength.

The selection of these indices was made mainly based on the spectral channels available in the multispectral camera (green, blue and near infrared). The GNDVI performs well in predicting yield for different crops (Wahab et al., 2018; Jewan et al., 2021; Yang et al., 2022); in addition, Gitelson et al. (1996) state that the green electromagnetic region has shown higher sensitivity to a wider range of chlorophyll concentration than the red region. As for the BNDVI, it requires the blue of the visible spectrum to monitor areas sensitive to chlorophyll content; likewise, good results in yield estimation have been obtained with this index (Lukas et al., 2022). However, Zeng et al. (2021) evaluated different indices at different phenological stages of wheat and observed that indices including the blue band of the spectrum have a low correlation with yield. In particular, at the flowering stage they obtained a value of R² = 0.7 (average of all indices with blue band).

With the information from each VI, a new image (raster) was generated for the two phenological stages using the map algebra technique, which allows combining spectral bands from mathematical operators, since it treats the spatial data layers as variables (Mali et al., 2005). Finally, the average value of the BNDVI and GNDVI was estimated for each experimental unit. To do this, the vector layer classified as vegetation by eCognition and the vector layer of the perimeter of the experimental units in the Qgis program (OSGeo Foundation, USA) were used. Subsequently, the average of each treatment was determined.

Grain yield measurement

At the end of the wheat growth cycle, grain yield was determined by sampling with a 1.0 × 1.0 m frame, which was randomly dropped within each experimental unit. The plants remaining within the frame were counted, and the grain was extracted from these plants and weighed on a precision balance to obtain the weight per m². The average for each treatment was obtained from the values of the corresponding replicates.

Data analysis and relationship of sensor values with yield

Data obtained with the multispectral camera (BNDVI and GNDVI), the GS and the SPAD502+, as well as the grain yield values, were subjected to a process of eliminating outliers using the interquartile range technique. Outliers have a significant influence on the arithmetic mean, and can generate averages that are not representative of the experimental units.

The average values of BNDVI, GNDVI, GS and SPAD502+ of each treatment were correlated with the grain yield data by means of a regression analysis performed with R studio 4.2.2 software and the ggplot2 tool for generating graphs. Two regression models (linear and polynomial) were applied to determine the goodness of fit of each sensor at the two phenological stages. The coefficient of determination (R²) was also estimated to determine the prediction quality of the models.

Results and discussion

Ground-based optical sensor measurements and vegetation indices

Table 1 shows the average yield per treatment, as well as the values obtained with the sensors and VIs at the thickened sheath phenological stage.

Table 1. Average yield per treatment, and values obtained with the sensors (GreenSeeker and SPAD502+) and vegetation indices (BNDVI and GNDVI) at the thickened sheath stage.

Table 2 shows the average yield of each treatment, as well as the values obtained with the sensors and VIs at the heading stage.

Table 2. Average yield per treatment, and values obtained with the sensors (GreenSeeker and SPAD502+) and vegetation indices (BNDVI and GNDVI) at the heading stage.

Table 3 shows the coefficients of determination obtained with the models (linear and polynomial) when comparing grain yield with the values obtained with the sensors and VIs, and it can be observed that the second-degree polynomial model better explains the yield variable, by presenting higher R² values.

Table 3. Coefficients of determination of the comparison of grain yield with sensor measurements and vegetation indices at two phenological stages.

Second-degree polynomial models have been used to evaluate the potential of optical sensors in grain yield estimation at different stages of wheat crop development (Holzman et al., 2014; Zhang et al., 2019). On the other hand, Hassan et al. (2019), Zeng et al. (2021) and Walsh et al. (2022) obtained a good correlation when using a linear model to estimate yield at various stages of crop development (from the tillering stage to harvest). The differences between the present work and those mentioned above could be due to different factors, such as crop variety, plant density, irrigation, nutrition, climate and even the type of sensor (Zhang et al., 2019).

Figure 5 shows a positive trend when comparing the values obtained with the sensors and the yield for the different treatments at the thickened sheath phenological stage; that is, the higher the value recorded by the sensors, the higher the yield. A similar trend was observed between the levels of nitrogen applied and yield.

The comparison of the sensor values and yield for the different nitrogen treatments during the heading stage is shown in Figure 6.

GS showed the best performance (R² = 0.8678) at the thickened sheath stage, while at the heading stage the best yield prediction was obtained with GNDVI, followed by GS (R² of 0.8713 and 0.8932, respectively). The SPAD502+ sensor had the lowest performance at both phenological stages, with an R² around.

The results obtained coincide with those reported by Walsh et al. (2022), who observed that the GS sensor and the NDVI (derived from images from a multispectral camera mounted onto a UAV) perform better in yield estimation than the SPAD502+, which measures point readings on the leaf. The differences observed between the sensors may be due to the N concentration in the leaves (Cartelat et al., 2005) and the crop canopy, as it is not homogeneous (Eichelmann et al., 2005). Monostori et al. (2016) recommend calibrating the SPAD502+ sensor values for each variety and crop type, which could improve the accuracy of yield estimates.

The GS performance results agree with those reported by Zhang et al. (2019), who used this sensor to analyze the wheat crop canopy during the stages after full cover (stages 8 to 10 on the Feekes scale), and observed high correlations with yield (R² = 0.9). Similarly, Zeng et al. (2021) obtained reliable accuracy when using the GNDVI to estimate yield at various stages of crop development, from thickened sheath to grain filling, which coincides with the findings in the present work.

Conclusions

Due to the great diversity of factors that influence the signal received by the sensors, more variables must be considered to develop the models, so that they can be representative of different crop conditions.

The use of indices, such as the GNDVI, has the potential to develop reliable models for yield prediction, since they detect the variation in the optical characteristics of the wheat crop. The use of spectral images obtained from cameras mounted onto UAVs allows obtaining coefficients of determination comparable to those achieved with the GS and SPAD502+. The advantage of using UAVs is that a larger area can be covered in less time, which is why they can be recommended as a viable alternative to replace the GS and SPAD502+. In this sense, the methodological proposal used in this work is a reliable alternative for acquiring images using UAVs at field scale; however, its effectiveness should be evaluated on a larger scale and with other crops to strengthen the results.