Can a Combination of UAV-Derived Vegetation Indices with Biophysical Variables Improve Yield Variability Assessment in Smallholder Farms?

The rapid assessment of maize yields in a smallholder farming system is important for understanding its spatial and temporal variability and for timely agronomic decision-support. We assessed the predictability of maize grain yield using unmanned aerial/air vehicle (UAV)-derived vegetation indices (VI) with (out) biophysical variables on smallholder farms. High-resolution imageries were acquired with UAV-borne multispectral sensor at four and eight weeks after sowing (WAS) on 31 farmer managed fields (FMFs) and 12 nearby nutrient omission trials (NOTs) sown with two genotypes (hybrid and open-pollinated maize) across five locations within the core maize region of Nigeria. Acquired multispectral imageries were post-processed into three VIs, normalized difference VI (NDVI), normalized difference red-edge (NDRE), and green-normalized difference VI (GNDVI) while plant height (Ht) and percent canopy cover (CC) were measured within georeferenced plot locations. Result shows that the nutrient status had a significant effect on the grain yield (and variability) in NOTs, with a maximum grain yield of 9.3 t/ha, compared to 5.4 t/ha in FMFs. Generally, there was no relationship between UAV-derived VIs and grain yield at 4WAS (r < 0.02, p > 0.1), but significant correlations were observed at 8WAS (r ≤ 0.3; p < 0.001). Ht was positively correlated with grain yield at 4WAS (r = 0.5, R2 = 0.25, p < 0.001) and more strongly at 8WAS (r = 0.7, R2 = 0.55, p < 0.001), while the relationship between CC and yield was only significant at 8WAS. By accounting for withinand between-field variations in NOTs and FMFs (separately), predictability of grain yield from UAV-derived VIs was generally low (R2 ≤ 0.24); however, the inclusion of ground-measured biophysical variable (mainly Ht) improved the explained yield variability (R2 ≥ 0.62, Root Mean Square Error of Prediction, RMSEP ≤ 0.35) in NOTs but not in FMFs. We conclude that yield prediction with UAV-acquired imageries (before harvest) is more reliable under controlled experimental conditions (NOTs), compared to actual farmer managed fields where various confounding agronomic factors can amplify noise-signal ratio.


Introduction
Rapid assessment of crop yield is critical in order to monitor and address yield gaps [1,2], yet robust sampling is often cost-intensive across agricultural fields or landscapes. Hence, technologies and methods that can provide leverage for quick and non-destructive data collection are critical for monitoring changes both spatially and temporally. The evolution of remote-sensing technologies for the acquisition and processing of remotely-sensed imageries is considered valuable for the assessment of yield or other agronomic variables at various scales [2,3]. However, access to high-quality imageries is generally constrained by associated costs and time-lag in processing, posing limitations to the prospects of utilizing imagery-derived data for various agronomic use-cases, both in-season and out-of-season. Yet, there is a critical need to evolve reliable and rapid approaches for timely decision support within smallholder farming systems. Smallholder farming systems of sub-Saharan Africa (SSA) are often characterized by fragmented farmlands and differentiated management practices [3][4][5][6]. Most landscapes are complex mosaics with diffuse field boundaries and trees. Therefore, imageries of the landscape need to have sufficient spatial and temporal resolution to mask out areas of vegetation or undesired features. Spatially-explicit data acquired over farming landscapes can improve the understanding of the variability and dynamics of agronomic processes and variables, especially in smallholder cropping systems where changes may be more frequent at smaller scales, including within farm fields. Quite often, these changes are influenced by the management preferences of the farmers whose decisions are mostly driven by various external factors, including accessibility and affordability of inputs [7,8]. Most smallholder farming landscapes are defined by different varieties that are sown on different dates, with different soil nutrient application or status. Consequently, they are often characterized as mosaic(s) of individual fields that have contrasting vegetation structures or types within a very small area (often, within tens of meters).
The recent advances in satellite-based remote-sensing of global land-cover coincides with the emergence of unmanned aerial/air vehicles (UAVs) for crop monitoring and yield assessment, with increasing application for agronomic decision-support across large-scale mono-cropped and smallholder multi-cropped farming systems [9][10][11]. UAVs were initially developed for military use but have become recognized as a tool to acquire high-resolution imageries that can be post-processed and analyzed to understand spatially varying agronomic factors at the field scale. Within the past five years, several researchers have reported the applicability of UAV for monitoring agronomic variables, e.g., [11][12][13], in different cropping systems and across diverse geographies [10,[14][15][16]. Many of these variables are considered as potential proxies for yield estimation [10,14,15], especially at the plot and field level (typically up to 1000 ha). There are several existing methods for estimating crop yields with remote-sensing. A popular approach is to relate measured location-specific yields to vegetation indices (VI) derived from RGB, multi-spectral, or hyper spectral camera sensors. VIs often provide strong expression of the ground cover and the chlorophyll content of green material [17,18]. Many VIs have been developed, with the most common being the normalized difference vegetation index (NDVI) [10,11,[19][20][21]. Based on varied relationships between different reflectance spectral bands, other VIs have also evolved with a strong relevance for agricultural systems. These include the normalized difference red-edge (NDRE), green normalized difference vegetation index (GNDVI), green canopy vegetation index (GCVI), red vegetation index (RVI), red-edge canopy index (RECI), and many others [22,23]. Because these VIs represent spectral (and to a lesser extent, structural) characteristics of the vegetation, they are potential proxies for rapid assessment of yield and yield variability. Further, when these VIs are derived from spatially-explicit remotely-sensed imageries, they can provide a very useful understanding of yield variability at varying spatial scales. For instance, the use of individual VIs for the diagnosis of nutrient constraints and crop yield differences between fields/plots have been promising, e.g., [12,14]. Further, researchers have proposed that a combination of VIs can provide additive sensitivity effect and improve the detectability of nuanced vegetational characteristics to improve the assessment of agronomic parameters, including yield, e.g., [23,24]. This is because single VIs can be constrained by vegetation structure and composition, which may be undetectable at specific wavelengths of the electromagnetic spectrum. Researchers [21,25,26] have suggested that the greenness of plants grown under adequate nutrient conditions can compromise the accuracy of the remotely-sensed NDVI (from the multispectral sensor) due to saturation within the green spectral band. This type of limitation can be avoided by using other VIs that rely on spectral information from other reflectance bands within the electromagnetic spectrum.
While the use of UAVs for agronomic decision-support is fast evolving [11], there are limitations. For instance, Watanabe et al. [27] indicated that sorghum plant height was overestimated by UAV and that high fertilization affected the relationship between UAV and ground-based measurements. Schut et al. [28] reported that vegetation indices did not capture all management and biophysical factors that can aid the accurate assessment of yield within fields. Yet, new generation UAV-borne sensors may offer improved assessment accuracies for crop monitoring, especially in combination with (few) ground level data [29,30]. According to Sibley et al. [31] and Schut et al. [28], repeated in-season measurements and good field-level accuracy are important criteria that can be used to derive useful information from remotely-sensed imageries for rapid yield (and other agronomic) assessments. Given the complexity of smallholder farming systems [2,6], it is important to assess the field applicability of UAV beyond experimental plot conditions and within complex farming landscapes, where they can be used for rapid farm-level decision support. Schut et al. [28] indicated yield variability explained by selected VIs within specific crop farms reduced greatly across an array of fields compared to within fields where there is a higher homogeneity, and they noted that accurate ground reference data may improve the assessment of in-season yield variability. These ground data can include biophysical variables, such as plant height (Ht) and canopy cover (CC), which have been reported as valuable for non-destructive yield(-variability) assessment in smallholder farmers' fields. These biophysical variables represent morphological characteristics and are useful for understanding the allometric characteristics of plants [2]. Therefore, this study was conducted within the maize-producing savanna region of Nigeria to assess in-season predictability of grain yield in multilocational smallholder maize farmers' fields using UAV-derived VI with and without ancillary observations of biophysical variables. Our objective is to advance knowledge on the application of UAVs for in-season monitoring of maize yield in farmers' fields by combining ground-measured biophysical variables with UAV-derived VIs to assess yield variability.

Study Area
This study was carried out at multiple locations within the core maize production region of Nigeria; namely, Bunkure, Doguwa, Funtua, Ikara, and Soba ( Figure 1). The locations are within the Sudan and Northern Guinea savanna agroecological zones of the Country, which is the major cropping regions for grains and legumes such as maize (Zea mays L.), cowpea (Vigna unguiculata L.), peanut (Arachis hypogaea L.), and soybeans (Glycine max L.). The target areas for UAV-based data collection were selected based on: (i) the location of nutrient omission trials (NOTs), which were established under a different research activity to identify and understand nutrient constraints that are limiting maize yield among smallholder maize farmers [32]; (ii) the willingness of proximal farmers to grant access for ground-truthing and yield assessment in their farms, and (iii) the advisory guidance of the National Space Research and Development Agency (NASRDA) in Nigeria, which is critical for compliance with regulatory requirements on UAV use in the Country. The study area is characterized by a growing period of 151-160 days, mean annual rainfall of 800-1250 mm, and daily mean minimum and maximum temperature of 12 • C and 32 • C, respectively. Rainfall is unimodal and lasts from June to October, is considerably optimal for grain crops, and is well distributed throughout the growing season, while the dry season starts from late October to May.

Nutrient Omission Trial Field Establishment
At the onset of the planting season for the region, from mid to late June 2016, NOTs (n = 100) were established to assess the impact of varying soil nutrient limitations on maize yield within the maize-based system of Nigeria under the Taking Maize Agronomy to Scale in Africa (TAMASA) project (www.tamasa.cimmyt.org; [32]). A subset of the NOTs (n = 12) were covered within the target locations for UAV flight. Each experimental unit was comprised of 12 contiguous plots (5.2 m × 4 m) planted with maize in two blocks of six plots, with one block sown with open-pollinated genotype (OPV) and the other sown with hybrid genotype (HV). Within each of the genotype blocks, the nutrient treatments comprised a combination of major nutrients required for maize production, including nitrogen (N), phosphorus (P), potassium (K), and micronutrients (+). Hence, each genotype block received six fertilizer treatments (i.e., Control, PK, NP, NK, NPK, and NPK+) applied on six plots across twelve (12) multilocational NOTs, which were applied at recommended optimal dosage [32] based on previous soil tests in the region. Therefore, a total of 144 plots were covered during flight missions at 4 and 8 weeks after sowing (WAS) across the entire study location. All nutrients were applied at the establishment stage of NOT, except for N, which was applied in 3 splits (at establishment, 3WAS, and 6WAS). Other details on rates and management of NOTs are presented by Shehu et al. [32].

Farmers' Fields
Because farmers made their farm-level decisions independent of our research interests, we screened prospective volunteer farmers to select only farmlands that were sown with maize within about 3 days of NOT establishment. The selected farmers' fields (n = 32) differed in size and management; this is typical within smallholder maize-based systems. The specific genotype choice of the farmers is generally unknown; however, within the maize-based area, farmers sow both the open-pollinated (OPV) and hybrid (HV) genotypes.

Nutrient Omission Trial Field Establishment
At the onset of the planting season for the region, from mid to late June 2016, NOTs (n = 100) were established to assess the impact of varying soil nutrient limitations on maize yield within the maize-based system of Nigeria under the Taking Maize Agronomy to Scale in Africa (TAMASA) project (www.tamasa.cimmyt.org; [32]). A subset of the NOTs (n = 12) were covered within the target locations for UAV flight. Each experimental unit was comprised of 12 contiguous plots (5.2 m × 4 m) planted with maize in two blocks of six plots, with one block sown with open-pollinated genotype (OPV) and the other sown with hybrid genotype (HV). Within each of the genotype blocks, the nutrient treatments comprised a combination of major nutrients required for maize production, including nitrogen (N), phosphorus (P), potassium (K), and micronutrients (+). Hence, each genotype block received six fertilizer treatments (i.e., Control, PK, NP, NK, NPK, and NPK+) applied on six plots across twelve (12) multilocational NOTs, which were applied at recommended optimal dosage [32] based on previous soil tests in the region. Therefore, a total of 144 plots were covered during flight missions at 4 and 8 weeks after sowing (WAS) across the entire study location. All nutrients were applied at the establishment stage of NOT, except for N, which was applied in 3 splits (at establishment, 3WAS, and 6WAS). Other details on rates and management of NOTs are presented by Shehu et al. [32].

Farmers' Fields
Because farmers made their farm-level decisions independent of our research interests, we screened prospective volunteer farmers to select only farmlands that were sown with maize within about 3 days of NOT establishment. The selected farmers' fields (n = 32) differed in size and management; this is typical within smallholder maize-based systems. The specific genotype choice of the farmers is generally unknown; however, within the maize-based area, farmers sow both the open-pollinated (OPV) and hybrid (HV) genotypes.

UAV-Based Acquisition of Imageries and Post-Processing
We used an eBee UAV (SenseFly Inc., Switzerland, www.sensefly.com/drone/ebee.html) mounted with multispectral 4C sensor (Airinov, France, www.airinov.fr) to acquire fine resolution imageries at each target location. The e-Bee is a light-weight fixed-wing UAV that can cover up to 600 ha in a single flight at 1000 m altitude and is equipped with an onboard global positioning system (GPS), solar irradiance sensor, and a ground calibration target. The multi-spec 4C camera has four passive sensors that record reflectance in four spectral bands: red (R), green (G), red-edge (RE), and near infra-red (NIR) at 1.2 megapixels per sensor. It has a global shutter, instant field of view (IFOV) of 0.9 mrad, and low luminosity (>3000 lux). The UAV was flown to acquire fine-resolution imageries that covered the area of interest (including NOT and farmers' fields) at 4 and 8WAS, coinciding with the onset of vegetative (V7) and tasseling (VT) growth stages, respectively [33], and within the mid-season growth period where best indication of post-harvest grain yield is obtainable [34].

Post-Processing
After the completion of each flight mission, the imageries were exported from the UAV along with the flight log files for post-processing with Pix4D software (Pix4D v.3.1., Switzerland, www.pix4d.com). The software provided a platform for end-to-end post-processing of the acquired imageries and offered the needed flexibility to configure processing parameters based on desired output quality and target end-use of the products. Overall, for each successful flight mission,~400 imageries were geotagged and processed through several stages to generate final outputs in geotiff formats, including reflectance bands corresponding to the four spectral reflectance domains of the multi-spec sensor, digital surface model, digital orthomosaics, and VIs. The VI imageries were computed from corresponding spectral bands, based on Equations (1)-(3). The consideration of VIs to be computed was limited to those that had been reported as promising for agronomic application at the canopy level and in relation to vegetation status in croplands, especially maize, with a focus on selecting VIs that were based on red, red-edge, and green bands [23,[35][36][37][38].
where ρnir, ρred, ρred.edge, and ρgreen are the spectral reflectance of the near infrared band, red band, red edge band, and green band, respectively.

Ground-Truth Data Collection
We conducted in-situ measurement of NDVI with Greenseeker Handheld Crop Sensor HCS 100 (Trimble Ltd., Sunnyvale, CA; https://agriculture.trimble.com/precision-ag/products/greenseeker). The Greenseeker was held above the canopy (0.6 m) while walking for 30-60 s through each NOT plot, or through marked quadrats (4 m 2 , n = 5) within each farmers' field. The device proximally scans leaf greenness (within a swath of~0.25 m) through its infrared sensors and displays an NDVI value averaged over the duration of the scan. We did not acquire ground-truth measurement for other VIs due to cost (money and labor) and because our goal was not to recalibrate the sensor, but rather to test the application of the UAV-borne multispectral sensor for assessing yield variability based on the indices derived.
We adopted the crop-cut method to quantify grain yield, as recommended by the FAO and generally regarded as the most objective method for yield estimation [14,39]. Harvest was conducted within 9 m 2 quadrats in each NOT treatment plot, based on standard NOT protocols, as presented by Nziguheba et al. [40]. In farmers' fields, maize cobs were harvested and grain yield was quantified in five (2 m × 2 m) quadrats positioned along a diagonal transect within the field. In both farmers' fields and NOTs, the harvested cobs were shelled, and the grain was oven dried to determine moisture content. The weight of grain yield per sampling quadrat were converted to yields in metric tons per hectare (t/ha) at 12% moisture content. The geographic coordinates were recorded as degrees latitude and longitude (position precision of~3 m) at the center of each plot/quadrat in NOT and farmers' field using a Garmin eTrex 20 GPS device (https://buy.garmin.com/en-US/US/p/87771#overview). Using a height ruler, we measured the height of three randomly selected maize stands within each quadrat, and the recorded values were later averaged per quadrat at 4WAS and 8WAS. Similarly, at both growth stages, we used the smartphone-based Canopeo app (http://canopeoapp.com/) to measure percent canopy cover in each quadrat, based on standard protocols [41].

Georeferenced Locations and Data Extraction
The recorded geographic coordinates (X, Y) were processed and exported into point shapefiles in WGS1984 datum in ArcMap. Using the UAV-derived imageries as reference, some misaligned coordinates were noted and corrected by editing the point shapefiles in ArcGIS. These misalignments were likely due to the precision level of the recreational GPS unit (~4 m); hence, the editing process was critical to ensure that each point shapefile location rested within the designated field and referenced the appropriate NOT plot or farmers' field. Using the geoprocessing "Create Regular Polygon" add-on tool in ArcGIS, 4-sided polygons (2 m 2 for farmers' field and 3 m 2 for NOTs) were generated for each ground-truthing point, using the point datasets as the centroid. These polygons represent the sampling support unit for the ground-truthing process and were subsequently used in R (Raster package [42]), to compute the zonal average statistics of UAV-derived VI cell values for each location at the matching support scale of the yield measurements.

Statistical Analyses of Data
Due to the slight skewness of the yield data, we applied log-linear transformation to normalize the dataset. The output from the spatial zonal statistics, which was computed from the UAV-derived imageries, was processed as a table and imported into R-Studio to assess the correlation of VIs with measured grain yield across the study locations, and separately within NOTs and farmer managed fields (FMFs). In the first level of the analyses, we used a linear multivariate regression approach (Equation (4)) to assess the variability of yield that is explained by ground-measured variables in NOTs.
where, β 0 is the intercept and β 1,2,i are slopes of the observed variables X 1 , X 2 , X i , respectively. In the second level, we independently and randomly split the NOT and FMF data into 70% calibration data (93 datapoints for NOT and 72 datapoints for FMF) and 30% validation data (37 datapoints for NOT and 27 datapoints for FMF). Some data records were excluded from further analyses in NOTs and FMFs due to irresolvable missing datapoints or incomplete records. The random split of each dataset was implemented at the Farm (i.e., field) level for independence of locations in the model calibration and validation stages. Notwithstanding the assumption that both NOT and FMF data belong to the same population (i.e., smallholder farmers), separate models were fitted for NOTs vs. FMFs. The calibration datasets were used to fit linear mixed effect models (Equation (5)) for yield prediction, using UAV-derived VIs at 4WAS and 8WAS as input variables, with and without the significant biophysical variables (as determined from the first step). The validation datasets were used to assess the accuracy of the prediction. Similar to Burke and Lobell [1], model fit was also applied to time series combination of 4WAS and 8WAS data to assess potential applicability for an improved prediction outcome (Equation (6)). The linear mixed effect modeling framework included a random parameter that accounts for the potential indeterminable effect of differing management practices among farmers. Statistical parameters that were evaluated to make inference include descriptive statistics (mean, range, coefficient of variation), correlation coefficient (r), coefficient of determination (R 2 ), root means square error of prediction (RMSEP), and significance, P (tested at 95% confidence level). The best relationship between measured yield and UAV-derived VIs is indicated by r closer to 1, R 2 closer to 1, RMSEP closer to 0, and an acceptable significance level (p < = 0.05).
where, β 0 is the intercept, β 1 is slope of UAV-derived Vegetation Index, VI 1 , β i is the slope of measured biophysical variable X i (when included in the model), and t is the growth stage at which the imageries were acquired (in weeks after sowing (WAS)). The additional term, 1|Farm, denotes assignment of Farm as random variable where several factors can influence yield outcome.

Estimated Grain Yield and Ground-Truth Biophysical Variables (gNDVI, Ht, and CC)
The overall average maize grain yield was 3.12 t/ha, with a lower average yield estimated in FMFs (2.75 t/ha) than the average yield across the NOT plots (3.6 t/ha; Table 1). The maximum grain yield (9.3 t/ha) was attained by optimizing nutrient and genotype combination in NOTs, exceeding the maximum estimated yield in farmers' fields (5.4 t/ha; Figure 2). In NOT plots, average ground-measured gNDVI, Ht, and CC increased from 4WAS to 8WAS across all treatments (Figure 3). Comparing treatments, the ground-measured variables were similar at 4WAS, but at 8WAS the NOT plots that received either NP, NPK, or NPK+ treatments showed significantly higher gNDVI, Ht, and CC (p > 0.001) compared to the Control, NK, and PK treatments.
The initial analyses of univariate relationships suggest that no correlation exists between grain yield and ground-measured variables (R 2 < 0.02; r < 0.14, Figure 4). However, by including the known explanatory variables (treatment, location, and genotype) into the linear multivariate model, the explained variation in yield greatly improved (R 2 = 0.45 for all data, R 2 = 0.67 for NOT, and R 2 = 0.14 for FMF). Additionally, the multivariate analysis of the NOT data showed that treatment and location, but not genotype, significantly explained the observed yield variability (R 2 = 0.50; p < 0.001). The average estimated yields for HV (3.41 t/ha) was comparable to OPV (3.68 t/ha), considering all locations and treatments.    The initial analyses of univariate relationships suggest that no correlation exists between grain yield and ground-measured variables (R 2 < 0.02; r < 0.14, Figure 4). However, by including the known explanatory variables (treatment, location, and genotype) into the linear multivariate model, the explained variation in yield greatly improved (R 2 = 0.45 for all data, R 2 = 0.67 for NOT, and R 2 = 0.14 for FMF). Additionally, the multivariate analysis of the NOT data showed that treatment and location, but not genotype, significantly explained the observed yield variability (R 2 = 0.50; p < 0.001). The average estimated yields for HV (3.41 t/ha) was comparable to OPV (3.68 t/ha), considering all locations and treatments.

UAV-Derived VIs and Their Correlation with Yield in NOT and FMF
The imagery and extracted data from the UAV-derived VIs showed nutrient-induced differences in vegetational characteristics in the NOTs (e.g., Figures 5 and 6) and general field-level variations across FMFs, with an observable transition between 4WAS and 8WAS growth stages (e.g., Table 2). The NDRE exhibited highest variation across farms, with the highest coefficient of variation (CV) of 1.13 in NOT plots at 4WAS, compared to GNDVI, which consistently had the lowest CV at each growth stage in both NOTs and FMFs. Pooling NOT and FMF data, the correlation analyses indicate that none of the VIs had a significant relationship with yield at 4WAS (r < 0.02, p > 0.1); however, a very weak correlation emerged at 8WAS for NDVI and GNDVI (r ≤ 0.3; p < 0.001). In the NOTs, the VIs and grain yield were weakly correlated at 4WAS (r = 0.23 and 0.33 for GNDVI and NDVI, respectively, p < 0.001), but this improved at the later growth stage at 8WAS (r = 0.40 and 0.47 for GNDVI and NDVI, respectively, p < 0.001). Contrastingly, in FMF, no meaningful relationship could be established between VIs and grain yield at 4WAS, while the significant relationships at 8WAS were weak (r ≤ 0.2, p < 0.001).   Agronomy 2020, 10, x FOR PEER REVIEW 14 of 24 Figure 6. Nutrient-induced variation of UAV-derived vegetation indices (VIs) within multilocational nutrient omission trials (NOTs) in smallholder maize farming systems at 8 weeks after sowing. NDVI denotes Normalized difference VI, NDRE is Normalized difference red-edge, and GNDVI is green normalized difference VI. * and ** denotes outlier data points in the boxplot.

Predictability of Grain Yield Variability with(Out) Biophysical Variables
Relating grain yield to each VI separately, the highest explained yield variability at both growth stages were observed in NOTs (49% at 4WAS and 54% at 8WAS), compared to FMFs (42% at 4WAS) in the calibration dataset. The validation metrics (Table 3a,b) showed low yield predictability with UAV-derived VIs alone in both NOTs and FMFs. The maximum prediction of yield variability (r = Figure 6. Nutrient-induced variation of UAV-derived vegetation indices (VIs) within multilocational nutrient omission trials (NOTs) in smallholder maize farming systems at 8 weeks after sowing. NDVI denotes Normalized difference VI, NDRE is Normalized difference red-edge, and GNDVI is green normalized difference VI. * and ** denotes outlier data points in the boxplot.

Predictability of Grain Yield Variability with(Out) Biophysical Variables
Relating grain yield to each VI separately, the highest explained yield variability at both growth stages were observed in NOTs (49% at 4WAS and 54% at 8WAS), compared to FMFs (42% at 4WAS) in the calibration dataset. The validation metrics (Table 3a,b) showed low yield predictability with UAV-derived VIs alone in both NOTs and FMFs. The maximum prediction of yield variability (r = 0.56; R 2 = 0.29) was achieved by fitting the mixed model with GNDVI values at 8WAS in NOTs (Table 3a). There was no meaningful improvement in assessed yield variability after merging VI data from both phenological stages (i.e., 4WAS + 8WAS) and analyzing them as time-series data for either NOT or FMF.
In NOTs, the results from the initial null model testing shows that the UAV-derived VIs are viable explanatory variables, considering that R 2 increases from 0.03 (null model) to the highest value of 0.29 (for GNDVI). In contrast, there was no improvement in the assessed yield variability across the FMFs, suggesting that the VIs are poor explanatory variables for yield variability. Further, the combination of measured Ht with each VI (Figure 7) showed that the predictability of grain yield at 4WAS and 8WAS increased significantly in NOTs (Table 3a) but not in FMFs (Table 3b). For instance, in NOTs, assessed variability of grain yield improved when NDRE was combined with Ht, with R 2 increasing from 0.03 to 0.63 (p < 0.001). Similarly, other VIs became better predictors of grain yield variability when combined with the Ht variable, with the R 2 value peaking at 0.64. The overall prediction error, RMSEP, mainly decreased (up to~79%) after the inclusion of the height parameter, except when NDRE data for both growth stages were combined (Table 3a). In FMFs, the predictability of grain yield with VIs did not improve meaningfully (at 4WAS) after the inclusion of the Ht variable in each model run, but rather declined (~20%) at 8WAS, while the RMSEP values were stable, hovering around 0.03-0.07 at both of the growth stages that were assessed (Table 3b).   ns denotes non-significant at significance level, α, = 0.05; WAS denotes weeks after sowing; a is the model intercept, b is the slope, r is the correlation, Ht is plant height. Coefficient, R 2 is the adjusted coefficient of determination and RMSEP is the root mean square error of prediction.   Table 3a,b.  Table 3a,b.

Grain Yield Relationship with Measured Biophysical Variables
In this study, the set-up of NOTs next to farmers fields provided important context for the understanding of grain yield gap and predictability in smallholder farming systems. Our findings also have important ramification for agronomic decision support and aspirations for rapid (in-season) assessment of grain yield with UAV-derived data in farmers' fields, especially in Africa. Across the multilocational NOT plots, the attainment of high yield (up to 9.3 t/ha) in plots where nutrient limitations are fully addressed, in contrast to control plots (where grain yield was as low as 0.49 t/ha), supports the notion that proper soil nutrient management can reduce the existing yield gap within the smallholder farming systems [4]. The availability of required nutrients for uptake at various growth stages of the maize crop affects overall plant development, including height, greenness, and canopy formation, which are represented by the measured biophysical variables (height, gNDVI, and percent canopy, respectively). The significant relationships assessed between grain yield and measured variables (mainly gNDVI and Ht) aligns with previous findings that support the selection of relevant proxies for rapid in-season assessment of yield variability [2,43]. However, the assessed relationships were weakly expressed (R 2 ≤ 0.17 and r ≤ 0.42), deviating from reported findings where strong relationships were established between maize grain yield and proximally-sensed biophysical variables in farmers' fields, with R 2 -values typically greater than 0.5 (e.g., [43]). It should be noted that the relationship between yield and proxy variables can be influenced by artefacts of location-dependent soil nutrient conditions and environmental factors, such as short-term drought conditions. Such artefacts can negatively affect the final yield outcome by compromising plant vigor during the reproductive/grain-filling stages (before-and after−8WAS). In the study area, soil nutrient limitations and poor soil management practices are common [8,32], and we observed that the few farmers who applied nutrients predominantly adopt early growth stage application practice (mostly by surface dressing or broadcasting) without considering the high potential for rapid nutrient losses due to the sandy textural characteristics of the soil. This prevalent practice may explain the poor crop performance (lower yield) in farmers' fields compared to NOTs where the nutrient management schedule included proper timings relative to growth stages and site conditions. More importantly, in FMFs, the expected waning-off of nutrient availability with plant growth bears potential implication for the usability of the UAV-acquired spectral information at two growth stages in this study. While earlier growth stages may have progressed well, with UAV-sensed spectral signals indicating the potential for realizable grain yield, nutrient deficiencies at later stages can alter reproductive processes and result in poor grain yields, thereby translating into a poor relationship between estimated grain yield and remotely-sensed spectral data. Further research exploration may be needed to unravel the effects of farmers' nutrient management practices at various growth stages on yield predictability; however, this is beyond the scope of this article.

Nutrients, Not Genotype, May Influence UAV-Derived Vis-Insights from NOT
Within the multilocational NOTs, nutrient variations were more strongly expressed in the UAV-sensed VIs than genotype variations. Similar to other crops, maize genotypes are often released with the aim of improving resilience or tolerance to stress (such as drought and weed), and they consequently improve quality and quantity of yield outcomes. In some instances, such improvement may include modification of traits related to leaf morphology and canopy characteristics, to achieve desired responses to target stress conditions. Despite the expected difference in the phenotypic traits of different genotypes [44], the implementation of this study during the wet season (i.e., no drought stress), with full control of weeds in the experimental plots (i.e., no weed stress), may have addressed the major stresses that are likely to compromise the growth and development of the OPV genotype in this study. Therefore, in the absence of other stressors, soil nutrient status exerted a dominant effect on the overall yield outcomes for both OPV and Hybrid genotypes, without noticeable contrast. This is further supported by our observed data on percent canopy coverage (not reported), which showed that canopy coverage was similar within the OPV and Hybrid genotype plots at both of the growth stages evaluated. Because N and P are the most prevalent soil nutrient deficiencies within the savannah maize-based system of Nigeria [45,46], the application of appropriate fertilizer likely supported similar leaf formation and canopy morphology of maize stands across the plots, with no significant difference between genotypes. This suggests that the inherent genotype diversity in smallholder farming systems may not be a major concern for canopy-level light interception or reflectance that underlies the quality (including signal-noise ratio) of UAV-derived indices across (the complex mosaics of) farmers' fields.

In-Season Grain Yield Variability Assessment with UAV-Derived VI: A Nuanced Outcome
The comparison of assessed grain yield variability between NOTs and FMFs unravels an often-neglected limitation to the applicability/transferability of agronomic predictive tools from experimental "controlled" conditions to usually "complex" farming systems. Despite the set-up of multilocational NOTs next to fields that were managed by smallholder farmers, we observed a strongly contrasting difference in the potential to assess grain yield variability based on UAV-derived VIs under both conditions. The observation that a combination of Ht with VIs improved the assessed grain yield variability in NOTs but not in FMFs contradicts our expectations that the biophysical variable will greatly improve outcomes for yield variability assessment, irrespective of the management condition. On the contrary, ground measured Ht (on a few stands within each farm) may not be universally useful to implement spatially explicit yield prediction within farmers' fields, and seemingly requires further in-depth enquiry. However, UAV-derived canopy height data may be very relevant for spatially-explicit prediction of grain yield if the mounted sensor is properly calibrated (with ground control references) to provide reliable elevational/surface height (z) data. Field data collection in this study did not include instrumentation to calibrate the UAV-sensed surface model of the farms, thereby precluding the generation of gridded high-resolution height or elevation data. In agreement with Burke and Lobell [1], broadening the use of remotely-sensed imageries (including UAV-derived) for rapid yield assessment, especially across complex smallholder farmers' fields, requires rigorous ground data for proper calibration of sensors and models. Arguably, an interplay of factors can also confound the usefulness of high-resolution UAV-derived VIs, including the aforementioned (amplified) noise that may be generated at the canopy level, especially in smallholder farming systems where canopy formation and structure are usually indeterminate across farms [2,34].
The grain yield variability explained by the UAV-derived VIs in FMFs (maximum of 42% in calibration and 26% in validation) is lower than previously reported in other studies that are based on ground-level proximate-sensing of crop fields over several time points throughout the growth and reproductive stages of the crop [12,21,43,47]. Recently, similar findings from complex farming systems in sub-Saharan Africa have been reported (e.g., [14]), where the explained maize yield variability hovered around 40%. It is probable that single (or sparse) time-step imageries are insufficient to fully elucidate in-season yield variability because growth processes that modulates final yield outcomes occur and can change at more frequent temporal scales between the time-points of imagery acquisition. Generally, VI imageries that are acquired during or close to the reproductive stages have been reported to be more suitable for yield assessment [28]. For instance, around 8WAS, the maize plants are closer to anthesis, and plant health at this stage is more likely to determine grain yield outcomes [33]. However, unanticipated environmental stress after anthesis can negatively impact the grain-filling process, with consequent implication for the final grain yield. The potential to achieve better predictability of grain yield by combining time series VIs from UAV-acquired imageries (e.g., at 4WAS + 8WAS) deserves further consideration. Our comparison of results between combined time-series data and single timestamp (especially at 8WAS) showed negligible or no improvement in the explained yield variability across farmers' fields. Although this may obviate the suggestion that time-series imageries may improve yield forecasting in cropping systems [28], it aligns with idea that vegetation sensing should be implemented close to the most critical reproductive stage of the target crop [21,29,48,49]. In reality, justifying the additional costs associated with the acquisition of multiple/time-series imageries (time and resources for flights and analytics) may be contingent on gains in accuracy and a broader demand for agronomic decision support.
Notwithstanding the low predictability of grain yield in farmers' fields, high resolution UAV-imageries are useful for the generation of agronomically relevant information about crop health and nutrient status in experimental plots and in farmers' fields. For instance, the UAV-derived VIs distinctly showed major differences between plots based on nutrients and genotypes in NOTs (Figure 5), and this lends credence to the overall usefulness of UAV-derived VIs' for the assessment of relative differences in nutrient status and overall crop health status in non-uniformly managed smallholder farming systems. We were unable to fully elucidate the other factors that may account for yield variability across the farmer managed fields compared to NOTs due to incomplete information from the volunteer farmers. Although the farmers agreed to grant access for data collection and provided basic information on approximate sowing date, they were unavailable and unwilling to document actual tending operations (such as weeding, type and quantity of fertilizer, organic manure application, genotype, etc.). Future study on complex farming systems should include careful assessment of these potential factors.

Conclusions
Successful acquisition of quality high-resolution imageries and processing of the agronomically-relevant vegetation indices is an important step towards understanding within-and between-farm variability of yield and related indicators at the field-scale. This study provides further insight into the potential use of UAV-derived vegetational indices to assess yield variability for rapid agronomic monitoring and robust decision-support in smallholder farming systems. By setting up multilocational nutrient omission trials close to several farmers' fields, our findings showed that nutrients, not genotype, significantly explained the observed yield variability. The weak predictability of maize grain yield variability based on selected indices (NDVI, NDRE, and GNDVI) indicates a lingering gap in UAV application for rapid yield assessment within complex smallholder farming systems beyond experimental conditions and large-scale (mono-cropped) field conditions. Because the combination of UAV-derived VIs with the ground-measured biophysical variable (mainly height) improved assessed yield variability in NOTs but not in FMFs, we believe that VIs acquired at only two timesteps do not provide a sufficient basis to reliably predict grain yield within farmers' fields. However, it is advantageous that the UAV-based sensing of croplands can generate continuous high-resolution canopy height data that may be useful for spatially-explicit relative yield prediction, contingent on accurate calibration. While the demand for in-season prediction of crop yield in smallholder farming systems remains topical, especially between-and within-fields, further explorations of UAV application should consider and account for potential confounding factors in farmers' fields (such as nutrient application regimes, soil characteristics, and planting density), which can vary between farms and influence canopy formation, light interception, and spectral reflectance at various stages of growth.