Calibration of Hybrid-Maize Model for Simulation of Soil Moisture and Yield in Production Corn Fields

Abstract

Model calibration is essential for acceptable model performance and applications. The Hybrid-Maize model, developed at the University of Nebraska-Lincoln, is a process-based crop simulation model that simulates maize growth as a function of crop and field management and environmental conditions. In this study, we calibrated and validated the Hybrid-Maize model using soil moisture and yield data from eight commercial production fields in two years. We used a new method for the calibration and multi-parameter optimization (MPO) based on kriging with modified criteria for selecting the parameter combinations. The soil moisture-related parameter combination (SM-PC3) improved simulations of soil water dynamics, but improvement in model performance is still required. The grain yield-related parameter combination significantly improved the yield simulation. We concluded that the calibrated model is good enough for irrigation water management at the field scale. Future studies should focus on improving the model performance in simulating total soil water (TSW) dynamics at different soil depths by including more soil water processes in a more dynamic manner.

1. Introduction

Crop models are essential tools for understanding and predicting crop growth and yield in response to weather variation and management [1,2,3,4,5,6,7]. With climate change being one of the significant uncertainties to crop productivity [8], simulation models have been used extensively to investigate the potential impacts of climate change on crop growth and yield [9,10,11]. However, successful applications are still a challenge [12,13], because the mathematical functions that describe crop growth and its responses to the environment are still relatively insufficient in capturing all interactions of the processes [14].

The Hybrid-Maize model is a process-based crop model [15] developed using experimental data largely from the US Corn Belt. It has been used to simulate maize yield responses to changes in climate and crop management at plot and field scales in different geographic regions [15,16,17,18,19,20,21,22,23,24]. Although the Hybrid-Maize model captures corn growth processes and yield satisfactorily in most reported studies, many of the crop physiology parameters are calibrated based on data that are a few decades old. With farmers’ constant adoption of newer crop cultivars [25,26,27,28,29,30], it becomes important to know how well the Hybrid-Maize model simulates corn yield and other processes when using the default model parameters, in production fields under current cultivar and management conditions, and how to calibrate model parameters effectively if needed.

Several calibration methods have been used in modeling, including least squares [13,31,32], Bayesian parameter estimation [33,34], Markov chain Monte Carlo parameter estimation [35,36], generalized likelihood uncertainty estimation (GLUE) [37,38,39], and manual trial-and-error method [3,40]. While most studies calibrated their models for yield only [40,41,42,43], some also calibrated for variables that are important to yield, including total biomass, LAI, and soil moisture. [44,45,46,47]. Since crop yield and many model outputs often depend on multiple model parameters [34], there is a need to calibrate all sensitive parameters, either separately or simultaneously. A good calibration method should not only reduce the uncertainty of parameters [14,48] but also mitigate the computational burden of the calibration process.

The Hybrid-Maize model does not have a built-in calibration routine. As a result, it requires a separate external process of model calibration. Because values of some model parameters are inter-correlated, it requires efficient and effective methods for this type of model calibration. A kriging-based multi-parameter optimization method for calibrating the Hybrid-Maize model was recently developed [20]. Although the stepwise approach was robust enough to cater to uncertainties in parameter estimation, the model was only calibrated for crop yield on plot-scale experiments. In addition, the Hybrid-Maize model was not calibrated for soil moisture, which is one of the major variables that affect crop yield in the western US Corn Belt. The objectives of this study were the following: (1) demonstrate a stepwise kriging approach for calibration and multi-parameter optimization (MPO) for daily soil moisture and final grain yield in production fields, and (2) validate the model outputs using independent datasets from the same farm.

2. Materials and Methods

2.1. The Study Fields

The field data were collected in two years from 2019 to 2020 from eight privately owned maize production fields in Elgin, Nebraska (Figure 1). The region has a humid temperate climate with a maize growing season from late April to late September. The long-term average annual precipitation is 686 mm. Figure 2 shows the monthly and long-term precipitation and grass-reference evapotranspiration (ETo) of crop growing seasons of 2019 and 2020. The predominant soil textures across all the fields are loamy sand and sandy loam (https://websoilsurvey.nrcs.usda.gov/app/, accessed on 10 June 2019). The soil is well-drained, with a slope of the fields ranging from 0 to 6%.

2.2. Crop Management

The eight fields are referred to as Home Field (HF), East Field (EF), Links Field (LF), Kelly Field (KF), North Koinzan Field (NKF), South Koinzan Field (SKF), Johnson Field (JF), and North Field (NF). Due to crop rotation between maize and soybean, only six fields were used in each of the growing seasons. The size of the fields ranges from 49 to 55 ha. Prior to planting, all fields were tilled using a disc harrow. Two row spacings, 51 and 76 cm, were used in all fields with a seeding depth of 6 cm. The farmer also used variable seeding rate technology for achieving optimal profits. The seeding rates were based on normalized yield return of the past 10 years.

Table 1. (a) Crop and management practices for the four continuous maize fields in 2019 and 2020. (b) Crop and management practices of the two maize after soybean rotation fields in 2019 and 2020.

(a)
		Field
	KF	HF	NKF	SKF
2019
Planting date	25 April 2019	4 May 2019	26 April 2019	26 April 2019
Harvest date	8 October 2019	15 October 2019	10 October 2019	11 October 2019
Hybrid	Channel 213-19stxrib; Pioneer P1197AM and P1197AMT	Channel 210-79stxrib 212-90stxrib, and 213-19vt2prib	Channel 213-19vt2prib; Pioneer P1197AM and P1370Q	Channel 213-19vt2prib; Pioneer P1197AM and P1370Q
Weighted average seeding rate (seed/ha)	81,500	80,220	90,100	87,290
Nitrogen fertilizer (kg N/ha)	281	269	278	280
Irrigation amount (mm)	227	214	193	193
Rainfall amount (mm)	518	503	514	514
2020
Planting date	1 May 2020	1 May 2020	25 April 2020	25 April 2020
Harvest date	20 October 2020	19 October 2020	12 October 2020	13 October 2019
Seed brand Cultivar	Channel 211-66stx and 213-93stxrib; Pioneer P1108Q	Channel 213-19stxrib; Pioneer P1108Q	Channel 213-19vt2 and 216-36stxrib; Pioneer P1415Q	Channel 213-19vt2 and 216-36stxrib; Pioneer P1415Q
Weighted average seeding rate (seed/ha)	81,390	78,340	86,800	87,070
Nitrogen fertilizer (Kg N/ha)	280	280	280	270
Irrigation amount (mm)	278	278	298	302
Rainfall amount (mm)	360	360	361	361
(b)
Field	EF	LF	JF	NF
	2019		2020
Planting date	2 May 2019	24 April 2019	26 April 2020	26 April 2020
Harvest date	16 October 2019	8 October 2019	14 October 2020	15 October 2020
Seed brand Cultivar	Channel 213-19stxrib and 213-19vt2prib	Channel 213-19stxrib; and Pioneer P1197AM	Golden Harvest 13H15 and Pioneer P1415Q	Golden Harvest 13H15 Channel 209-51VT2PRIB, 211-66stx, and 213-19VT2PRIB.
Weighted average seeding rate (seed/ha)	79,680	80,900	78,980	82,490
Fertilizer (Kg/ha)	268	280	283	268
Irrigation amount (mm)	202	218	298	284
Rainfall amount (mm)	503	518	360	334

a summarizes the crop management for the same four fields out of the six fields in each growing season in both years, while Table 1b contains similar information but for other four fields, two from each year due to crop rotation.

2.3. Measurement of Soil Moisture and Maize Yield

Acclima Time Domain Reflectometer (TDR-315 and TDR-315L, Acclima, Inc., Meridian, ID, USA) soil moisture sensors were installed at depths of 0.25, 0.46, 0.66, and 0.86 m at the dominant soil texture in each field. The sensors were connected to data loggers and soil moisture data were recorded at 30 s intervals. For data analysis, however, the data were aggregated to a daily interval. Corn grain yield for each field was determined by the combine yield monitor and adjusted to a 15.5% grain moisture content.

2.4. Model Simulation Setup

Input data required by the Hybrid-Maize model include daily weather data (including solar radiation, maximum and minimum temperature, precipitation, humidity, and reference evapotranspiration), planting date, hybrid brand and maturity, and plant population. The daily weather data were obtained from the Elgin weather station, which is less than 1 km away from the farthest field used in our study. This weather station (Figure 1) is part of the Automated Weather Data Network (AWDN) of the High Plains Regional Climate Center (HPRCC). For user-set irrigation, the model also requires soil texture, maximum rooting depth, bulk density, surface residual coverage, and soil moisture at planting, and record of irrigation date and amount were also required. The model assumes optimal management for nutrients, and control of insects and diseases.

The input settings were the following: loam sand for soil texture, maximum soil rooting depth of 150 cm, soil bulk density of 1.3 g/cm³, available water of 75% at planting for topsoil and at 100% for subsoil, soil surface residue coverage of 50%. A generic seed brand with the estimate of growing degree days (GDD) to silking and total GDD to physiological maturity were used as inputs as recommended in the model user manual [49] for a situation where the seed brand cultivated in the study is unavailable in the Hybrid-Maize database. The weighted average plant population for each field was used in the model simulation. Other agricultural management data, such as planting dates and irrigation amounts and dates (Table 1a,b), were obtained from the farmer.

2.5. Sensitivity Analysis of Model Parameters

Seven soil moisture-related parameters and twelve growth-related parameters were selected for both one-at-a-time sensitivity (OAT) and global sensitivity analysis (GSA) based on expert knowledge of the model structure (

Table 2. Hybrid-Maize model (version 2018) parameters related to soil water content (SWC) and crop yield [49].

Parameter Abbreviation	Parameter Description	Unit	Default Value
SWC parameters
POR	Porosity	%	0.4400
GAM	Texture-specific constant	cm−2	0.0330
PSImax	Texture-specific suction boundary	cm	200
Ksat	Saturated hydraulic conductivity	cm/d	26.50
Alfa	Texture-specific geometry constant	cm−1	0.0398
Ak	Texture-specific empirical constant	cm−2.4 d−1	16.40
BD	Bulk density	g/cm3	1.3
Yield parameters
G5	Potential kernel filling rate	mg kernel−1 day−1	8.70
G2	Potential number of kernels per ear	kernel ear−1	675
ILUE	Initial light use efficiency	g CO2 MJ−1 PAR	12.5
GRG	Growth respiration coefficient of grain	g CH2O g−1 dry matter	0.490
GRL	Growth respiration coefficient of leaf	g CH2O g−1 dry matter	0.470
MPR	Maximum photosynthetic rate	g CO2 m−2 leaf h−1	7.0
K	Light extinction coefficient	-	0.55
ECT	Efficiency of carbohydrate translocation from stem or leaf to grain	-	0.260
MRG	Maintenance respiration coefficient for grain	g CH2O g−1 dry matter d−1	0.0050
GRS	Growth respiration coefficient of stem	g CH2O g−1 dry matter	0.520
MFB	Maximum fraction of leaf biomass at silking that can be translocated as carbohydrate from leaf to grain	-	0.15
SLW	Empirical parameter determining the relative contribution of a soil layer to water uptake	-	3.0

).

For soil moisture OAT analysis, selected parameters were allowed to vary between the minimum and maximum values of each parameter and across the soil texture spectrum, as documented in the Hybrid-Maize model. For crop growth OAT, changes prescribed at ±5%, ±10%, ±15%, and ±20% of the default values were used. Since soil moisture outputs from Hybrid-Maize are a time series, an objective function [deviation sensitivity (DS) in our case] comparing simulated daily soil moisture time series (based on relative changes in one parameter at a time) to that using default parameters was used for the OAT analysis. As an objective function is used to maximize or minimize variables when dealing with mathematical optimization problems, the objective function is to minimize DS. If a change in a soil moisture-related parameter does not change the simulated soil moisture values from the default simulation values, then the DS would be zero. For yield-related parameters, the OAT sensitivity was determined by the slope of the relative percentage change of simulated yield versus the percentage change of that parameter.

For both soil moisture GSA and yield GSA, forward stepwise selection (regression-based method) and Sobol (variance-based method) similar to the approach of [20] were used but with soil moisture time series converted to DS values. Our forward stepwise selection method used the Akaike Information Criterion (AIC) values to rank and select the model parameters based on their order of importance. It started with the most influential parameter and ended with the least influential one. The Monte Carlo estimation of Sobol indices as reported by [50,51] was used since Sobol GSA could handle nonlinear responses and evaluate the influence of interactions in non-additive systems. The sobol2007 function in the “sensitivity” package was developed in R with 10,000 runs to determine the total Sobol GSA indices for Hybrid-Maize model parameters. It is important to point out that although OAT and the more robust (forward stepwise and Sobol) SA methods were compared, only the more robust SA methods were eventually used to select the soil moisture and yield-related parameters.

2.6. Multi-Parameter Calibration and Validation

The model was first run using default parameter values. It produced poor results with large discrepancies between observed and simulated yield values (Figure 3). For model calibration, combinations of six fields per year for the two years were mixed and split into two datasets, the calibration and validation datasets, to accommodate the climatic variability between the wet and dry years. This prevented biased model calibration for a specific year. The four fields present in both years (HF, KF, NKF, and SKF) were first split into two, so that both calibration and validation datasets had two wet and two dry fields. That is, two fields (HF and KF) in 2019 were combined with two fields (NKF and SKF) in 2020, making a total of four field-years used for the calibration. The remaining four fields (two per year) that were not present in both years (i.e., EF, LK, JF, and NK) were added to the validation dataset, making a total of eight field-years. With this, the calibration and validation datasets were drawn from the same underlying “target population” (i.e., same environmental conditions and management practices), which could be referred to as “interpolation studies”.

Table 3. Calibration and validation dataset for wet and dry years.

	2019 Fields	2020 Fields
Calibration dataset	KF, HF	NKF, SKF
Validation dataset	NKF, SKF, EF, LF	KF, HF, NF, JF

shows the calibration and validation datasets used in this study.

Multi-Parameter Optimization for Soil Moisture and Yield Simulations

As soil water directly affects crop yield, soil water simulation was first calibrated, followed by yield. We adopted the MPO approach used by [20], but with modifications such as using a standard deviation of the spatial yield variation within each field as a criterion for selecting the best parameter combinations instead of using ±8% of the observed yield as a constraint when calibrating yield parameters. The pairwise calibration and MPO approach for each response (soil moisture and yield) involved selecting the top two most sensitive parameters, calibrating, and fixing them, and then calibrating the next two sensitive parameters until all the selected parameters had been calibrated in pairs. Each pairwise calibration and MPO approach was able to produce ten thousand yield simulation outputs using the ordinary kriging interpolation type in Surfer 16 software (Golden Software, LLC, Golden, CO, USA) with the best variogram model. For better understanding, the steps for the pairwise calibration and MPO are:

(a)

The top two parameters that were most sensitive to soil moisture were selected based on the results of OAT and GSA sensitivity analyses.

(b)

A grid made of 64 two-parameter combinations (8 × 8 grid nodes) was generated based on the parameter range of each parameter in the pair.

(c)

For the 64 two-parameter combinations in (b), the Hybrid-Maize model was run while keeping other sensitive parameters at default values. Soil moisture response (i.e., RMSE between each simulated and observed time series) was recorded for every parameter pair.

(d)

The grid generated in (b) and the corresponding soil moisture response in (c) were used to create a 3-dimensional response surface using ordinary kriging interpolation in Surfer 20.1 software. Gaussian and wave components in Surfer are the two variogram models that were used for kriging to give the best output grids with the lowest errors and best cross-validation results. In order to deal with potential trends in the model parameters, AutoFit tool in Surfer 20.1 software was used. This tool takes a user-specified variogram model and an initial set of parameters and attempts to find a better set of parameter values. The response surface generated was made up of 10,000 parameter pairs (100 × 100 nodes) and their corresponding DS values. These are equivalent to 10,000 Hybrid-Maize model simulations.

(e)

Cross-validation was carried out to determine the accuracy of the soil moisture DS response surface by randomly selecting 30 nodes (n = 30) from the 10,000 kriged nodes with the exclusion of the 64 nodes in (b). Each of the pair parameters corresponding to the 30 selected output grid was then run using the Hybrid-Maize model and the simulated output was compared with that of the kriged node to ascertain the accuracy of the kriged response surface.

(f)

The MPO process was carried out by selecting the parameter pair with the lowest objective function (DS) on the response surface. The lower the DS, the more the simulated daily soil moisture matched the observed daily soil moisture time series.

(g)

The averages of the top two most sensitive soil moisture parameters from (f) were fixed and steps (b–f) were repeated for the next two sensitive parameters until all the sensitive parameters had been calibrated. If there was an odd number of sensitive parameters, the least sensitive parameter was calibrated based on OAT approach.

(h)

After calibrating for soil moisture, steps (a–g) were repeated for the parameters that were most sensitive to crop yield based on the results of OAT and GSA sensitivity analyses. However, the MPO process for crop yield was carried out by averaging all the parameter combinations that met three constraints similar to [22], but with simulated yield varying between one standard deviation (±SD) of the observed yield instead of ±8% of the observed yield. This approach was deemed appropriate, as SD was indicative of the uncertainty from the natural/inherent yield variation across each field.

The crop yield constraints mentioned above were used to prevent overfitting the model parameters for each calibration field, and to also reduce the uncertainty due to equifinality of models since many simulations using multiple combinations of parameters (i.e., equal pathways) could give the exact observed crop yield (i.e., finality) for each calibration field. Mathematically, the modified yield MPO is written as:

$$\begin{array}{c}\underset{A^{*},B^{*}}{\mathrm{arg\; min}}\left|Y_{sim}-Y_{obs}\right|, \mathrm{subject}\mathrm{to}:\\ -SD\left(Y_{obs}\right)<Y_{sim}<SD\left(Y_{obs}\right),\\ A^{*}\in \left[0.85\left(A_{def}\right),1.15\left(A_{def}\right)\right],\\ B^{*}\in \left[0.85\left(B_{def}\right),1.15\left(B_{def}\right)\right]\end{array}$$

where the optimal input parameter argument represents the parameter pair {$A^{*}, B^{*}$} that minimizes the value of the objective function $\left|Y_{sim}-Y_{obs}\right|$ with the added constraints; $Y_{sim}$ is the simulated yield in Mg ha⁻¹; $Y_{obs}$ is the observed yield in kg ha⁻¹; $A^{*}$ and $B^{*}$ are the optimized model parameters; $A_{def}$ and $B_{def}$ are the default model parameters; and SD is the standard deviation. It is important to note that the second and third yield-related constraints were based on the premise that the re-calibrated yield-related parameters for each field would not be too far (±15% of default values) from the default values since the current Hybrid-Maize model had already been calibrated for the US Corn Belt region. The stepwise MPO approach is also illustrated in Figure 4.

2.7. Model Evaluation

We used the following indices to evaluate model performance [3,21,52], including percent deviation (PD), root-mean-square error (RMSE), normalized root-mean-square error (NRMSE), mean absolute error (MAE), Nash–Sutcliffe efficiency (NSE), and index of agreement (d).

$$PD=\frac{\left(S_i-M_i\right)}{M_i}\times 100$$

(1)

$$RMSE=\sqrt{\frac{{\sum }_{i=1}^n{\left(S_i-M_i\right)}^2}{n}}$$

(2)

$$NRMSE=\frac{RMSE}{\overline{M}}$$

(3)

$$MAE=\frac{1}{n}{\sum }_{i=1}^n\left|S_i-M_i\right|$$

(4)

$$NSE=1-\frac{\sum _{i=1}^n{\left(M_i-S_i\right)}^2}{\sum _{i=1}^n{\left(S_i-\overline{S}\right)}^2}$$

(5)

$$d=1-\frac{{\sum }_{i=n}^n{\left(M_i-S_i\right)}^2}{{\sum }_{i=1}^n\left(|S_i\right.-\overline{M}|+|M_I-{\left.\overline{M}|\right)}^2}$$

(6)

where n is the number of paired samples, M_i and S_i are the measured and simulated values of the ith observation, respectively, and $\overline{M}$ and $\overline{S}$ are the average measured values, respectively.

PD measures the percent deviation between observed and simulated values. RMSE ranges from zero to positive infinity with the same unit as both observed and simulated variables [53]. RMSE can be standardized using the mean of the measured values and regarded as NRMSE, which is expressed as a fraction. As proposed by [54], model simulations were considered excellent, good, fair, and poor based on the NRMSE values of <10%, 10–20%, 20–30%, and >30%, respectively. MAE measures the average magnitude of errors in a set of predictions, without considering their direction. Both MAE and RMSE express average model prediction errors in the units of the variables of interest and range from zero to positive infinity. Neither is affected by the direction of errors. NSE describes the relative magnitude of the residual variance in comparison with the measured data variance. Its values range from negative infinity to 1, where NSE = 1 represents perfect correlation between simulations and observations or perfect model fit, and NSE < 0 implies that the model is a weaker predictor than the mean of the observations [55]. The index of agreement (d) represents the ratio of the mean square error and the potential error because it determines the additive and proportional differences in the measured and simulated means and variance. However, d is overly sensitive to extreme values owing to the squared differences. d varies between 0 and 1, where d = 1 indicates a perfect match, and d = 0 implies no agreement [56].

3. Results and Discussion

3.1. Weather Conditions of the Experimental Years

The 2019 growing season (May to September) was wetter and cooler than 2020 and the long-term average (Figure 2a,b). The 2019 growth season had a total of 466 mm rainfall, which was 63 and 106 mm more than the long-term average and 2020, respectively (Figure 2a,b). The cumulative ET_o of the 2020 growing season was 209 mm higher than 2019 but 67 mm lower than the long-term average. The average daily maximum temperature of the 2020 growing season was higher than 2019, while the daily minimum temperature for 2020 was lower than 2019.

3.2. Sensitivity Analysis

The results of the sensitivity analysis are shown in

Table 4. One-at-a-time (OAT) and global sensitivity analysis (GSA) of parameters related to soil moisture and yield.

Response	Parameter	OAT Slope	Stepwise GSA AIC	Sobol GSA Total Index
Soil moisture	GAM	15.00	−438	0.90
	PSImax	28.00	−479	0.07
	BD	0.27	−507	0.03
	Ksat	0.00	-	0.00
	AK	0.00	-	0.00
	Alfa	0.00	-	0.00
	Porosity	0.00	-	0.00
Yield	G5	0.69	94	0.45
	G2	0.66	57	0.41
	ILUE	0.65	9	0.12
	MPR	0.13	2	0.01
	GRG	0.22	-	0.01
	K	0.09	-	-
	ECT	0.09	-	-
	MRG	0.04	-	-
	GRS	0.03	-	-
	GRL	0.02	-	-
	MFB	0.00	-	-
	SLW	0.00	-	-

Notes: Model performance increases as the AIC value decreases. Parameters whose AIC values were not in the table were not chosen in the stepwise regression.

. PSImax and G5 were the most sensitive parameters for soil water content and grain yield, respectively. Of all the soil moisture-related parameters considered in the OAT sensitivity approach, PSImax, GAM, and BD were the only parameters sensitive to simulated soil moisture with PSImax, GAM, and BD having sensitive indexes of 15.00, 28.00, and 0.27, respectively. Except for MFB and SLW, the yield was sensitive to changes to all other parameters considered in this study based on OAT alone. G5 had the highest OAT slope of 0.69, while GRL had the lowest OAT slope of 0.02. Previous OAT sensitivity studies also observed that G5, G2, and ILUE are the most sensitive parameters related to yield. Refs. [17,42] reported that, in a decreasing order, G5, G2, and ILUE were the three most sensitive yield-related parameters. Ref. [22], on the other hand, reported a different order of importance with ILUE being the most sensitive parameter, then followed by both G5 and G2 with an equal sensitivity index.

For the stepwise GSA, soil moisture was again sensitive to changes in only three parameters (GAM, PSImax, and BD), with GAM also being the most sensitive one. This is contrary to what was observed for the OAT approach, where PSImax was the most sensitive parameter and ranked above GAM. The stepwise GSA picked just three soil moisture parameters before terminating, starting with GAM and ending with BD, with each consecutive parameter selection improving the model as shown by the reducing AIC values in Table 4. Adding extra parameters beyond that did not improve the model’s performance. However, a similar order of importance when compared to OAT was observed for the yield-related parameters. The stepwise GSA picked just four model parameters before terminating, starting with G5 and ending with MPR, with each subsequent parameter selection improving the model as shown by the reducing AIC values in Table 4. Since the lower the AIC value, the better the model performance, the stepwise GSA for yield showed that a combination of the top four most sensitive yield parameters resulted in a lower AIC value of 2.64 when compared to that with three sensitive parameters (AIC = 9.04). However, adding extra parameters after MPR did not improve the model’s performance.

Results from Sobol GSA for soil moisture-related parameters depicted a similar order of parameter importance when compared to OAT, while yield-related parameters depicted a similar order of parameter importance to both OAT and stepwise GSA. The only difference was that, whilst the Sobol GSA method selected five sensitive yield parameters, stepwise GSA and OAT selected four and ten parameters, respectively. For soil moisture parameters, GAM, PSImax, and BD had total indices of 0.90, 0.07, and 0.03, respectively, while for yield parameters, G5, G2, ILUE, MPR, and GRG had total indices of 0.45, 0.41, 0.12, 0.01, and 0.01, respectively. It is important to note that a parameter may be highly sensitive during OAT but may not be significant when in combination with other parameters in a GSA. This is due to possible interactions among parameters in a GSA. Based on the result of the sensitivity analysis, GAM and PSImax were selected for the calibration of soil water content while G5, G2, ILUE, and GRG were chosen for the calibration of grain yield.

3.3. Soil Water Content Calibration and MPO

The kriged DS surfaces (with 10,000 nodes or simulations) created with the top parameter pair {GAM, PSImax} were used for calibrating the model for soil moisture at the calibration fields based on the multi-parameter optimization steps discussed earlier. As previously stated, the MPO process for the {GAM, PSImax} pair was carried out by minimizing the DS between the observed and simulated daily soil moisture, and selecting the parameter pair with the lowest value on the 3-D DS surfaces. Only three fields were selected from the soil water content calibration and multi-parameter optimization. This is because these are the only fields with measured soil water content across the two growing seasons. The optimized {GAM, PSImax} pair produced the lowest soil moisture DS nodes with values of 0.02 mm, 0.03 mm, and 0.04 mm for the FK, NKF, and SKF fields, respectively. Since the {GAM, PSImax} pair accounted for about 97% of the Sobol GSA total indices (Table 4), a further calibration for soil moisture using the third most sensitive parameter (BD) resulted in insignificant changes in the daily soil water content dynamics.

On each DS response surface in Figure 5, the red dot shows the DS value using the default {GAM, PSImax} pair, while the white square shows the lowest DS value corresponding to the best {GAM, PSImax} pair. Although each response surface has distinctive attributes peculiar to the calibration field under consideration, there seem to be similarities in the shapes of the response surfaces for KF 2019 (wet year) and for SFK 2020 and NFK 2020 (dry year). According to each of these three response surfaces in Figure 5, there was a shift between the red and white dots, such that GAM reduced while PSImax increased for the MPO-calibrated pair.

GAM is an empirical soil texture-specific constant, while PSImax is soil texture-specific maximum water potential at the field capacity. The default values of GAM and PSImax for loamy sand in the Hybrid-Maize model are 0.0330 cm⁻² and 200 cm, respectively. The GAM is derived from an empirical function with little to no very good physical meaning and is impossible to measure. PSImax, on the other hand, can be measured independently. Regardless of the field and across the two growing seasons, the default GAM shifted to a new optimized value lower than the default value (Figure 5). The optimized GAM values for the three calibration fields have a close range of values (0.0193–0.0280 cm⁻²) as shown in

Table 5. Default and calibrated parameter scenarios.

Variable	Parameter	Field Calibrated Parameter				SM/GY-PC1	SM/GY-PC2				SM/GY-PC3
		2019 Wet Year		2020 Dry Year			HF	KF	NKF	SKF
		HF	KF	NKF	SKF
Soil	GAM		0.0280	0.0263	0.0193	0.0330		0.0280	0.0263	0.0193	0.0245
	PSImax		300	300	300	200		300	300	300	300
Yield	G5	7.8	7.7	8.1	8.0	8.7	7.8	7.7	8.1	8.0	7.9
	G2	603	592	625	621	675	603	592	625	621	611
	ILUE	12.2	12.0	12.2	12.4	12.5	12.2	12.0	12.2	12.4	12.2
	GRG	0.495	0.494	0.490	0.495	0.490	0.495	0.494	0.490	0.495	0.494

Notes: SM/GY-PC1: Soil moisture/grain yield-related parameter combination based on default value; SM/GY-PC2: Soil moisture/grain yield-related parameter combination for each field-year calibrated values; SM/GY-PC3: Soil moisture/grain yield-related parameter combination based on overall pooled average calibrated values; GAM: Texture-specific constant (cm⁻²); PSImax: Maximum water potential at field capacity (cm); G5: Potential kernel filling rate (mg kernel⁻¹ day⁻¹); G2: Potential number of kernels per ear (kernel ear⁻¹); ILUE: Initial light use efficiency (g CO₂ MJ⁻¹ PAR); GRG: Growth respiration coefficient of grain (g CH₂O g⁻¹ dry matter).

. For the same soil texture across different fields, the GAM value should be constant. However, there is a possibility of having a narrow GAM value as shown in Table 5 for the sample soil textural class, considering that there is a range of particle size distributions of silt, sand, and clay for the same soil texture class as well as different levels of compaction across the various fields for the same particle size distribution of the same soil texture. In addition, GAM may be closely related to pore space tortuosity, which is a geometric parameter that describes interconnected pore spaces. A lower GAM value implies a decrease in the ratio of the actual flow path length to the straight distance between the ends of the flow path. This ratio could be reduced by a reduced total pore fraction, and the total pore fraction for the topsoil (loamy sand) found in our study area could have been reduced by tillage, the use of heavy equipment, and heavy irrigation [57,58].

In addition, regardless of the field and across the two growing seasons, the default PSImax value shifted to a new optimized value of 300 cm, which is the upper limit for PSImax in the calibration space. This upper limit for PSImax was based on the maximum value for all soil textures in the Hybrid-Maize model. The higher calibrated PSImax value (300 cm) relative to the default value (200 cm) indicates an increase in the suction boundary, which could consequently reduce the total water readily available in the topsoil. Similar to a lower GAM value, a higher PSImax could also be due to the effect of topsoil compaction processes related to tillage, heavy irrigation, or compaction due to the traffic of heavy equipment during field operations [59].

It is important to note that the default values of soil parameters in Hybrid-Maize are merely approximations or representative values. Firstly, the soil texture triangle represents a range of particle size distributions of sand, silt, and clay. Since a range of combinations of sand, silt, and clay could end up in the same texture class, soils within the same texture class can have different particle combinations and, hence, different parameter values.

Secondly, soil organic matter (SOM) content has a strong influence on the water-holding capacity of any soil texture, especially for coarse textures such as sandy soils [60,61]. The Hybrid-Maize model, however, does not consider the effect of SOM on soil moisture. This implies that two soils with the same soil texture but different SOM content could have different water-holding properties in reality, whereas the model would not make a distinction. In addition, water retention, porosity, and saturated hydraulic conductivity are very sensitive to compaction and thus are spatially sensitive. Therefore, two fields with the same soil texture but different compactions may have different soil properties. Lastly, as stated earlier, different soil layers or horizons could affect the net or equivalent parameter values in the root zone on the field when compared to values from the laboratory analysis of soil samples.

3.4. Yield Calibration and Multi-Parameter Optimization

The MPO results based on the kriged yield surfaces for both the {G5, G2} and {ILUE, GRG} pairs are shown in Figure 6a,b. The most sensitive yield parameter pairs were calibrated in a sequential order after calibrating {GAM, PSImax} for soil water, as discussed in the previous section. The red contour line in each yield surface in Figure 6 represents the average crop yield observed in each calibration field, while the two black lines depict one standard deviation (±SD) about the average observed crop yield based on the spatial variation of crop yield measured by the combine harvester. The red dot shows the model-simulated yield using the default parameter values, while the white square shows the model-simulated yield using the calibrated parameters via the MPO approach.

As shown on the yield surfaces in Figure 6a, multiple parameter combinations {G5, G2} could result in yield simulations within ±SD boundaries (black contour lines) of the average crop yield in each field. In fact, there are multiple {G5, G2} combinations that could result in yield simulations with the exact average crop yield (red contour lines) in each calibration field. This is also true for the yield surfaces in Figure 6b with respect to {ILUE, GRG} pairs. Owing to this issue of equifinality (i.e., equal pathways), the average of all the parameter pair combinations (i.e., nodes on the kriged yield surfaces) that met the three yield MPO constraints were used to determine the calibrated parameter pair values for each calibration field as shown by the white squares. It is important to note that aggregating only the parameter pairs that are on the red contour line is the same as calibrating the model parameters such that the simulated yield matches the observed average yield exactly. This would amount to overfitting the model while not capturing the spatial variation of crop yield within the ±SD boundaries (black contour lines) of the average crop yield in each calibrating field. With every pairwise step of the yield MPO (i.e., {G5, G2} first and then {ILUE, GRG} later), the simulated yield value, based on the calibrated parameter pair (white square), moves closer to the average observed yield (Figure 6a,b). This is based on the contribution of each parameter pair to the overall yield simulation.

Although each yield surface in Figure 6a or Figure 6b is unique to the calibration field under consideration, the shape of all the yield surfaces for each parameter pair ({G5, G2} or {ILUE, GRG}) appears to be comparable across the different fields and years. As observed for each yield surface in Figure 5a, the position of simulated yield using the default parameter values (red dot) fell outside the black contour lines (±SD yield constraint) on the yield surfaces for all the calibration fields, while the white square falls within the two black contour lines. On the other hand, for each yield surface in Figure 6b, the positions of simulated yields using both the default (red dot) and the calibrated (white dot) parameter values fell within the black contour lines (±SD yield constraint) on the yield surfaces for all the calibration fields. When compared to the shift from the red dot to the white square in each yield surface in Figure 6a, the shift in Figure 6b is relatively smaller, because each red dot in Figure 6b is a simulation based on the combination of the default {ILUE, GRG} values and the previously calibrated {G5, G2} values from the preceding step. In other words, the white squares in the {G5, G2} parameter space shown in Figure 5a are the same as the red dots in the {ILUE, GRG} parameter space shown in Figure 6b. This is expected since {G5, G2} is a more sensitive parameter pair than {ILUE, GRG} as shown in Table 4.

As indicated in Table 5, G2, with a default of 675 kernels per plant, indicates the potential number of kernels per plant; G5, with a default value of 8.7 mg d⁻¹, represents a potential daily grain filling rate (mg d⁻¹ kernel⁻¹); ILUE, a conservative parameter with a default value of 12.5 g CO₂ MJ⁻¹ PAR (photosynthetically active radiation), denotes the photosynthetic rate at the very low level of radiation; and GRG, which is set at a default value of 0.49, signifies the efficiency when converting carbohydrate to grain biomass. Although these default values are used for high-yielding maize hybrids at plant densities ranging from 60,000 to 85,000 plants/ha prevalent in the US Corn Belt, they were mostly derived or calibrated using data from two decades ago. Advancements in crop breeding and crop management would call for new examinations of those parameters in the model.

Lower calibrated values of G2 across the fields and years than the default values were observed in Table 5. This could be related to the fact that the currently higher-yielding maize hybrids cultivated at currently higher plant densities have a tendency of bearing smaller cobs with potentially fewer kernels per plant. Recent plant breeding is more tailored towards improving yield on a per ground area basis instead of individual plants, because the best individual plant performance may not result in the best yield based on per ground area [62]. Again, as opposed to the G5 default value of 8.7 mg d⁻¹ kernel⁻¹ in the current Hybrid-Maize model, lower calibrated values of G5 were observed across the different fields and years (Table 5). This could be a result of higher plant density and thus stronger competition for resources among individual plants.

ILUE describes the CO₂ assimilation (carbon fixation) at a very low light intensity and is regarded as a conservative parameter [63]. Table 5 shows that a narrow range of lower calibrated values of ILUE was observed across the different fields and years. Note that ILUE is different from radiation use efficiency (RUE) used in some crop models. RUE is the average net biomass production per unit of light intercepted over the entire season, with a value ranging from 3.3 to 3.8 g CH² MJ⁻¹ PAR [64], while ILUE is the maximum CO2 assimilation efficiency without accounting for respiration losses. Like ILUE, the default value of 0.49 of GRG was adapted from [65]. This value is related to the composition of grains, especially the protein content. From the calibration results, the values of GRG obtained for the different fields across the growing seasons were closer to the default value of 0.490.

Since the farmer cultivated more than one hybrid brand but with comparable genetic properties in the same fields, it is worth noting that the new yield-related parameters are not associated in absolute terms with any of the hybrids. However, it may reflect a close sense of the biological attributes of most recent hybrids because hybrids within the same seed company (i.e., the parent company owes several subsidiaries) may share a lot of the same genetic materials and have comparable biological characteristics. This is important because sensitivity analysis or model calibration does not have to be very specific to cultivars in order not to limit the potential applications of the model. In addition, emphasis should not be placed on the true meaning of the calibrated parameters but on the need for the model to first reproduce as much of the field observations as possible in order to use the model for crop management, such as irrigation scheduling. This is important because no matter how good a model theory is, it has to be able to simulate close to reality before it can be applied effectively.

3.5. Simulation of Soil Water Content

Typical pivot irrigation, as used in this study, takes three to four days to cover a regular field. The model, however, assumes that the irrigation is completed in one day on the specified date of irrigation. This creates a potential mismatch between model-simulated and measured soil moisture. In order to harmonize this inherent difference in response between the simulated and measured soil moisture, we used a three-day moving average of measured soil moisture. To ensure that the Hybrid-Maize model was well calibrated for irrigation management in the study area, it is important to evaluate the effectiveness of the calibrated model in simulating soil water dynamics in the major rooting depth [0–0.3 m (SD1) and 0.3–0.6 m (SD2)] as well as in the entire 1 m [i.e., 0–1.0 m (SD3)].

Figure 7, Figure 8, Figure 9 and Figure 10 show the comparisons between the measured and simulated three-day moving average 1 m depth total soil water (TSW) using calibrated soil moisture parameter combination (SM-PC2) and overall pooled soil moisture calibrated parameter combination values (SM-PC3). The red lines depict simulated daily TSW using SM-PC2, while the purple lines show the simulated daily TSW using SM-PC3 for the calibrating fields (Figure 7).

Table 6. (a) Goodness-of-fit indicators for total soil water (TSW) for different soil depths using the calibrated Hybrid-Maize parameters for the calibrated fields. (b) Goodness-of-fit indicators for total soil water (TSW) for different soil depths using the calibrated Hybrid-Maize parameters in the validation fields.

(a)
	Year	Field	Goodness-of-Fit Indicators
			0–0.3 m Soil Depth (SD1)					0.3–0.6 m Soil Depth (SD2)					0–1.0 m Soil Depth (SD3)
			RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d	RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d	RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d
SM-PC1	2019	KF	6.76	11.80	5.12	−19.04	0.11	12.68	21.00	12.10	−10.65	0.13	30.64	15.60	27.98	−7.23	0.13
	2020	NKF	12.20	21.80	9.48	−13.15	0.22	21.26	34.30	19.11	−38.27	0.09	47.50	24.10	42.24	−13.69	0.08
		SKF	23.96	33.90	21.11	−25.52	0.02	27.57	39.80	25.53	−31.41	0.01	74.11	32.80	69.75	−21.55	0.00
SM-PC2	2019	KF	6.28	11.00	4.81	−7.85	0.17	11.93	19.70	11.30	−9.35	0.15	27.96	14.20	25.12	−5.85	0.16
	2020	NKF	10.21	18.20	7.91	−8.91	0.26	12.87	20.70	8.30	−13.39	0.19	25.15	12.80	19.15	−3.12	0.16
		SKF	11.36	16.10	9.52	−4.95	0.07	11.02	15.90	6.67	−4.17	0.08	23.68	10.50	19.17	−1.30	0.17
SM-PC3	2019	KF	8.53	14.90	7.44	−14.99	0.08	6.22	10.30	5.47	−1.80	0.39	12.99	6.60	10.77	−0.48	0.45
	2020	NKF	9.96	17.80	8.12	−8.43	0.29	12.68	20.40	8.15	−12.96	0.19	24.12	12.20	18.76	−2.79	0.18
		SKF	13.27	18.80	9.33	−7.13	0.02	11.20	16.10	6.71	−4.35	−0.07	27.62	12.20	20.13	−2.13	0.17
(b)
	Year	Field	Goodness-of-Fit Indicators
			0–0.3 m Soil Depth (SD1)					0.3–0.6 m Soil Depth (SD2)					0–1.0 m Soil Depth (SD3)
			RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d	RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d	RMSE (mm)	NRMSE (%)	MAE (mm)	NSE	d
SM-PC3	2019	NFK	8.59	15.30	7.53	−10.62	0.31	3.81	6.70	2.92	−1.49	0.60	20.81	11.90	18.42	−8.75	0.24
		SKF	10.26	18.90	8.87	−3.50	0.40	6.7	11.00	5.65	−0.76	0.49	12.62	6.80	9.46	−0.41	0.59
		EF	16.38	34.30	15.45	−20.78	0.10	16.6	42.60	16.36	−39.86	0.07	53.63	3.80	52.78	−44.97	0.05
		LF	8.98	15.40	7.63	−5.73	0.22	6.36	11.10	5.36	−0.39	0.21	29.05	13.80	23.22	−2.56	0.20
	2020	KF	10.31	17.00	8.26	−5.75	0.13	17.05	26.00	14.29	−8.34	0.07	36.20	17.00	29.72	−4.94	0.06
		JF	10.76	16.60	7.82	−8.08	0.35	13.78	22.10	10.58	−19.77	0.24	23.89	11.80	18.34	−4.84	0.35
		NF	12.19	23.40	10.33	−3.04	0.19	10.15	17.30	7.40	−2.25	0.36	19.84	11.00	17.20	−0.10	0.57

Notes: SM-PC1: Soil moisture-related parameter combination based on default value; SM-PC2: Soil moisture-related parameter combination for each field-year calibrated values; SM-PC3: Soil moisture-related parameter combination based on overall pooled average calibrated values.

a,b shows the goodness-of-fit statistics of the model calibration and validation.

Prior to calibration, the simulations were run using default model parameters (SM-PC1) for the calibration fields. For SD1, results of using SM-PC1 showed poor simulations with high RMSEs ranging from 6.7 mm to 23.9 mm, NRMSEs in the range of 11.8% to 33.9%, and MAEs in the range of 5.1 mm to 21.0 mm as well as low NSEs ranging from −25.3 to −9.1, and d values varying from 0.02 to 0.22. For SD2, high RMSEs ranged from 12.7 mm to 27.6 mm, NRMSEs in the range of 21.0% to 39.8%, and MAEs in the range of 12.1 mm to 25.5 mm, as well as low NSEs ranging from −38.3 to −10.7, and d values varying from 0.01 to 0.13. In addition, for SD3, the model-simulated TSW had very low accuracy as indicated by the high RMSEs ranging from 30.6 mm to 74.1 mm, NRMSEs in the range of 15.6% to 32.8%, and MAEs in the range of 28.0 mm to 69.7 mm, as well as low NSEs ranging from −21.5 to −7.2, and d varying from 0.00 to 0.08. The large errors observed in the 1 m rooting depth were due to the additive nature of errors from the different soil depths of SD1 and SD2, as well as the remaining soil depths between 0.6 m and 1 m across the various calibrating fields.

Results of TSW simulations slightly improved when SM-PC2 and SM-PC3 were used with both simulated daily TSW having similar daily fluctuations across the two years compared to measured TSW for the calibrating fields throughout the growing seasons. When the SM-PC2 was used, the performance metrics for the SD1 resulted in moderate values ranging from 6.3 mm to 11.4 mm, 11.0% to 18.2%, 4.8 mm to 9.5 mm, −8.9 to −5.0, and 0.07 to 0.26 for RMSEs, NRMSEs, MAEs, NSEs, and d, respectively, while for SD2, the values ranged 11.0 mm to 12.8 mm, 15.9% to 20.7%, 6.7 mm to 11.3 mm, −13.4 to −4.2, and 0.08 to 0.19, correspondingly (Table 6a). For SD3, the performance metrics improved with a moderate accuracy of RMSEs ranging from 23.7 mm to 28.0 mm, NRMSEs in the range of 10.5% to 14.2%, MAEs in the of 19.2 mm to 25.1 mm, NSEs ranging from −5.8 to −1.3, and d varying from 0.16 to 0.17. In addition, when the SM-PC3 was used, the model simulation of TSW slightly improved for SD1, resulting in values varying from 8.5 mm to 13.3 mm, 14.9% to 18.8%,7.4 mm to 9.3 mm, −15.0 to −7.1, and 0.02 to 0.29 for RMSEs, NRMSEs, MAEs, NSEs, and d, in that order, while for SD2, the values ranged 6.2 mm to 12.7 mm, 10.3% to 20.4%, 5.5 mm to 8.2 mm, −13.0 to −1.8, and 0.07 to 0.39 for RMSEs, NRMSEs, MAEs, NSEs, and d, respectively. For SD3 also, performance metrics improved with RMSEs ranging from 13.0 mm to 27.6 mm, NRMSEs in the range of 7.0% to 12.2%, MAEs in the range of 10.8 mm to 20.1 mm, NSEs ranging from −2.8 to −0.5, and d varying from 0.17 to 0.46.

Except for SKF 2020, in most cases for the various soil depths across the calibrating fields, the simulated TSW using SM-PC3 performed better than using SM-PC2. This could be the result of a lower GAM value of 0.0245 for SM-PC3 as opposed to the values of 0.0280 and 0.0263 for SM-PC2 for KF and NKF, respectively, as the same soil texture, loamy sand, was used for the model input. It is worth noting that the irrigation amount and distribution in terms of their order of magnitude across the two calibrating fields (NKF and SKF) in the 2020 growing season were almost comparable but the sensors at the same soil depths across the two fields were responding in a different order of magnitude, either partly because of differences in their calibration and sensitivity of the sensors or variation in soil texture at the locations where the sensors were installed. However, considering that we cannot transfer each unique field calibrated parameter combination values (SM-PC2) to another field during validation and in line with our intention to have a single set of calibrated parameter combinations to be used by the farmer across the various climatic conditions in the study area, only the SM-PC3 was used for comparing the field-measured TSW during validation.

During validation, however, the statistical indicators also depict a slightly improved model performance, as indicated by relatively low RMSEs ranging from 8.5 mm to 16.4 mm and NRMSEs in the range of 15.3% to 34.3%, MAEs ranging from 7.53 mm to 15.45 mm, as well as relatively high NSEs varying from −20.8 to −3.1, and d in the range of 0.10 to 0.40 for SD1 (Table 6b). In addition, for SD2, the statistical indicators resulted in values ranging from 3.8 mm to 17.1 mm, 6.7% to 42.7%, 2.92 mm to 16.36 mm, −39.9 to −0.39, and 0.07 to 0.60 for RMSEs, NRMSEs, MAEs, and d, respectively. On the other hand, for SD3, the statistical indicators show a poor agreement between measured and simulated TSW with RMSEs ranging from 12.6 mm to 53.6 mm, NRMSEs in the range of 6.8% to 38.0%, MAEs ranging from 3.5 mm to 52.8 mm, NSE varying from −44.9 to −0.1, and d in the range of 0.05 to 0.59. The model did not exhibit any consistent and systematic trend of over- or under-estimation of the TSW in all fields for SD1 and SD2, as well as SD3 in both 2019 (wet year) and 2020 (dry year) growing seasons (Figure 8, Figure 9 and Figure 10). In general, the performance of the model in simulating TSW at the different layers (SD1 and SD2) as well as the entire rooting depth (SD3) was very poor during validation. In part, this could be a result of the model’s poor performance in simulating evapotranspiration.

Similar findings have been reported. A slightly better agreement between the model-simulated and field-measured data without calibration for both irrigated and rainfed treatments with NRMSEs varying from 11.0 to 19.5% and absolute mean prediction error (MPE) in the range of 0.5% to 3.7% was reported [21]. However, it should be noted that the relatively low values of the above statistical indicators resulting in high model performance were due to the small number of samples used for comparison, and the experiment was conducted in carefully monitored experimental plots where the potency of data collection should be relatively high and accurate as opposed the large volume of continuous SWC data collection mechanism used by the farmer in this study. Ref. [66] reported that the smaller the number of samples used for comparing the model-simulated and field-measured data, the tendency for high model performance. Ref. [52] reported highly unsatisfactory performance when the Hybrid-Maize model was used in simulating water balance components for rainfed, limited, and full irrigation conditions. The authors observed negative values of NSE coupled with high NRMSEs in the range of 47–62% across the two soil layers of SD1 and SD2 and noted that rainfed treatments had the greatest disparities between simulated and measured data, indicating that the model is unsuitable for rainfed/dry locations. Although each model structure is unique and parameterized differently and, in most cases, should not be compared with other crop models, there is reported evidence of the poor performance of other crop models in simulating TSW in standard-size fields solely managed by farmers, where the fields are not carefully monitored, and the potency of data collection seems to be relatively poor. Ref. [65] conducted research in farmer’s fields with an average size of about 35 ha to evaluate the performance of the AquaCrop model in simulating field observations relating to leaf area index (LAI), crop evapotranspiration, soil water content, biomass, and final yield.

Despite the careful calibration of each field, the authors observed poor fitting with an RMSE ranging from 8.4% to 11.7%, and NSE in the range of 0.03 to 0.72. They attributed less goodness-of-fit indicators or the bias of estimates of available soil water (ASW) to the overestimation of actual transpiration and underestimation of soil evaporation by the model. The authors finally suggested that the maximum standard crop transpiration coefficient (K_cTr) and actual crop canopy cover (CC) curve proportionality should be improved, as the latter is sensitive primarily to water stress during the vegetative stage and not to daily water stress throughout the season.

The Hybrid-Maize model assumes that soil water is not retained above field capacity (FC) and drains to the subsequent soil layers underneath at the end of a day. This may not be true with the continuously field-measured soil moisture because there may be readings when the soil water content is above field capacity or even saturated, especially during major rainfall or irrigation events. That may be responsible for some of the underestimations of the TSW at the top-soil depths. Therefore, if the model is revised in relation to the above problem to accommodate periods when saturation may occur, this may improve the model simulation by reducing the peak gap between simulated and measured soil moisture values as well as the error in computing or estimating soil water balance.

3.6. Simulation of Yield

The performance of the grain yield simulation using default parameter combinations (GY-PC1) and calibrated parameters is presented in

Table 7. (a) Measured and default-simulated yield of all the fields in 2019 and 2020 growing seasons. (b) Goodness-of-fit for yield using different parameter combinations.

(a)
Yield (Mg ha−1)
Year		Field	Meas	Std. Dev (±)		GY-PC1		Difference		PD (%)
2019		HF	13.9	1.97		18.4		4.5		32.4
		KF	13.4	1.68		18.6		5.2		38.8
		NKF	15.4	2.06		18.8		4.6		31.5
		SKF	15.6	1.28		18.8		4.3		28.9
		EF	13.4	1.97		18.5		3.4		22.1
		LK	12.7	1.96		15.9		3.2		20.5
2020		NKF	14.6	3.19		19.2		5.1		38.1
		SKF	14.9	2.62		19.2		3.2		25.2
		HF	13.6	2.29		18.6		5.0		36.8
		KF	13.2	2.6		15.6		2.4		18.2
		JF	13.3	2.74		17.3		4		30.1
		NF	15.0	2.55		18.5		3.5		23.3
MAE (Mg/ha)						4.03
RMSE (Mg/ha)						4.12
NRMSE (%)						29.00
(b)
Year	Calibration Yield (Mg ha−1)					Validation (Mg ha−1)
	Field	Meas.	GY-PC1	GY-PC2	GY-PC3	Year	Field	Meas.	GY-PC1	GY-PC3
2019	HF	13.9	18.4	15.0	15.4	2019	NKF	15.4	18.8	15.8
							SKF	15.6	18.8	15.8
	KF	13.4	18.6	14.6	15.6		EF	13.4	18.5	15.6
							LK	12.7	15.9	12.9
2020	NKF	14.6	19.2	16.3	15.6	2020	HF	13.6	18.6	15.1
							KF	13.2	15.6	12.7
	SKF	14.9	19.2	16.1	15.6		JF	13.3	17.3	14.0
							NF	15.0	18.5	15.0
	MAE (Mg/ha)		4.65	1.30	1.35		MAE		3.73	0.71
	RMSE (Mg/ha)		4.66	1.32	1.47		RMSE		3.83	1.00
	nRMSE (%)		33.00	9.30	10.30		nRMSE		27.30	7.20

Notes: PD, percent deviation of default-simulated yield from the field-measured yield; Meas., field-measured yield; Std. dev., standard deviation of field-measured yield; GY-PC1, grain yield-related parameter combination based on the default value; GY-PC2, grain yield-related parameter combination for each field-year calibrated values; GY-PC3, grain yield-related parameter combination based on overall pooled average calibrated values.

a,b. The default model predicted grain yields across all the fields and years with low accuracy and with percentage deviations (PD) ranging from 18.2–38.8%. The PD range corresponds to an overprediction error range of 2.4–5.2 Mg/ha and the pooled data MAE of 4.0 Mg/ha, RMSE of 4.1 Mg/ha, and NRMSE of 29.0% before calibration (Table 7a). Ref. [64] reported a moderate to poor simulation accuracy with an underestimation of 1.2 Mg/ha for fully irrigated treatments, an overestimation of 2.0 Mg/ha for limited irrigation treatments, and up to 5.7 Mg/ha underestimation for rainfed treatments, with R² of 0.69, mean bias error (MBE) of 0.0–0.44 mm, NRMSE of 14%, and NSE of 0.57 using pooled data without calibrating the model. These values are better than those obtained in this current research prior to calibration. This could be because of the small size of the fields used in the referenced study, with relatively fewer yield variations coupled with a carefully monitored experiment. Also, the GY-PC1 in the Hybrid-Maize model may be close to the true parameters (if they were measured) in the referenced study. In another research conducted by [21], the authors reported PDs of yield ranging from 2.3–11.7% under both irrigation and rainfed treatments in a semiarid environment. Contrary to the results in the present study using the GY-PC1, results from the previous investigation showed that the model tends to underestimate grain yield prior to calibration for irrigated treatments. The poor default simulation results (GY-PC1) in our study (Table 7b) forced the need to calibrate the model. The simulation results of each field show a good match between the field-measured and model-simulated grain yield using each field’s calibrated parameter combination (GY-PC2). However, the goal of this study was to have a good combination of calibrated parameters for the study area. Therefore, we evaluated two different grain yield-related parameter combinations (GY-PC2 and GY-PC3) from all the calibrated parameter values of the calibrated fields across the two growing seasons (Table 7b). The results showed improved yield simulations with GY-PC2 significantly reducing the GY-PC1 MAE, RMSE, and NRMSE values by 72%, 72%, and 72%, respectively, while GY-PC3 reduced the ME, RMSE, and NRMSE using GY-PC1 values considerably by 71%, 69%, and 69%, respectively. Besides using the GY-PC1 values, the largest discrepancies between measured and simulated grain yield were observed using GY-PC3, wherein the MAE = 1.4 Mg/ha, RMSE = 1.5 Mg/ha, and NRMSE = 10.3%, while the lowest discrepancies between measured and simulated grain yield were detected using GY-PC2 with MAE = 1.3 Mg/ha, RMSE = 1.3 Mg/ha, and NRMSE = 9.3%. With the above evaluation results, GY-PC2 was the best parameter combination considering that its corresponding reduction estimation errors were the highest. However, even though the simulation results using GY-PC2 were considered the best parameter combination during calibration, it is important to note that GY-PC3 still fell within the acceptance range. Additionally, considering the need to have a single set of calibrated yield-related parameters for all fields across years, GY-PC3 was the only calibrated parameter combination evaluated during validation. This decision supports the need to have a single parameter combination that would accommodate the climatic variations between years and across fields for regional decision-making.

Validation was carried out using GY-PC3 across the two growing seasons for the validation fields. Again, some correlational regression-based and deviation-based statistical methods were also used to evaluate the model performance since each statistical method would provide unique information relating to the accuracy of the simulation. The validation results in Table 7b depict that GY-PC3 significantly reduced the MAE, RMSE, and NRMSE of using GY-PC1 by 84%, 74%, and 74%, respectively. Figure 10 shows the results of the Hybrid-Maize model grain yield simulation using GY-PC3 for the pooled data for the validation fields across the two growing seasons of 2019 and 2020. R² of 0.55 as shown in Figure 11 indicated that the model simulated grain yield across the various fields well. Although there was no consistent or systematic trend of under- or over-estimation using GY-PC3 for all the fields across the two growing seasons, only one was under-estimated and three fields closely match the 1:1 line. This shows an improved simulation of the grain yield across the fields as opposed to using the default parameter (GY-PC1). In a study where four different parameter combinations were compared in simulating grain yield using the Hybrid-Maize model under different agronomic management practices from different farmers in pilot-scale fields, Ref. [20] reported an improved Hybrid-Maize model simulation of grain yield with an MAE of 1.1 Mg/ha, RMSE of 1.4 Mg/ha, and NRMSE of 7% using the best parameter combination. The authors suggested that using in situ sensors combined with the calibrated model’s precision/adequacy would be more beneficial for improving in-season crop management and yield predictions for field-scale study.

3.7. Limitations, Practical Considerations, and Recommendations for Improvement

Our study attempted to use data collected from commercial-scale production fields for the calibration and validation of the Hybrid-Maize model, instead of first using data from small research plots before moving to large fields. Although large fields may produce uncertainty in calibration and validation, as in our study, they provide a direct idea of the magnitude of errors that may be encountered in real-field situations when using the model.

One of the limitations of this study is that only eight fields were used for the new method of model calibration, which may not be enough to test the validity of the innovation in the research. With this limited field dataset, there is a high risk of overfitting a crop simulation model with the new method of calibration and the overfitted model will not be able to generalize well to new datasets. Another limitation is that the current Hybrid-Maize model does not consider the upward capillary flux of soil water when simulating soil water dynamics, especially for locations with possibilities of high-water tables, as in the case of our study area with considerable upward flux. It is also worth noting that the assumption of instantaneous drainage of soil water above field capacity used in the Hybrid-Maize model is another serious drawback. This tends to support why the model showed declining patterns in soil water, whereas the soil water sensors showed almost no decline, because it may take some time for the soil water above field capacity to drain down the soil profile. Another possible reason for the no significant decline in the soil water measured by the soil water sensors could be that the electronic soil moisture sensors do not work very well.

Although the re-calibrated model produced a relatively moderate overall average RMSE of less than 20 mm in relation to the total average applied water (precipitation + irrigation) of 651 mm, there is still a need to reprove the part of the model responsible for simulating soil water balance.

4. Summary and Conclusions

In this study, we demonstrated the use of calibration and MPO approach for enhancing the Hybrid-Maize simulation of v and grain yields in production maize fields. The sensitivity analyses revealed that the texture-specific constant (GAM) and texture-specific suction boundary (PSImax) had the greatest impact on the simulation of SWC, while grain yield simulations were most sensitive to potential kernel filling rate (G5), the potential number of kernels per ear (G2), initial light use efficiency (ILUE), and growth respiration coefficient of grain (GRG) out of all the model parameters taken into consideration. We concluded the GY-PC3 (G5 = 7.9 mg kernel⁻¹ day⁻¹, G2 = 611 kernel ear⁻¹, ILUE = 12.2 g CO₂ MJ⁻¹ PAR, and GRG = 0.494 g CH₂O g⁻¹ dry matter) was a good, calibrated parameter combination for simulating grain yield, while SM-PC3 (GAM = 0.0245, PSImax = 300) is a good, calibrated parameter combination for soil moisture. Future studies should also consider applying and evaluating the MPO approach to other distinct groups of hybrids with unique characteristics and farmers’ fields across different agroclimatic zones.