Analysis and evaluation of a dynamic model for greenhouse lettuce growth

INTRODUCTION

Crop models are logical or quantificational algorithms which capture quantitative information on how a crop interacts with its environment during growth (Asseng et al., 2014). They simulate the behavior of a real crop by predicting its growth and development as influenced by environmental conditions and crop management (Li et al., 2012). Many crop models have been developed to predict growth changes and to assist in controlling the crop growing environment, particularly for the crop in greenhouses. This is because it´s easy to improve the quality of crop growth and yield by controlling the greenhouse environment (van Henten, 1994).

Many published studies in crop modeling (Ehret et al., 2011; Ramírez-Pérez et al., 2018; Gong et al., 2019; Wang et al., 2021) have been carried out in the field of greenhouse crops, which have accelerated the development of crop production’s automation and management. There are several greenhouse crop growth models. Jones et al. (1991) developed a physiological model of tomato development and growth called TOMGRO (TOMato GROwth), which accurately describes the differences between growth and yield of tomatoes, using a source-sink approach for the distribution of carbohydrate to different organs. TOMGRO was evaluated using boundary data (Shamshiri et al., 2016), and applied to predict the tomato yield in greenhouses based on controllable greenhouse environmental parameters (Lin et al., 2019). HORTISIM (HORTIcultural SIMulation) is a combined model for a variety of crops in greenhouses, and it is designed to predict the yields and timing of production related to crops, greenhouse characteristics and climate control (Gijzen et al., 1998). Körner & Holst (2017) developed an open-source modelling platform based on HORTISIM model. Marcelis (1994) studied a dynamic growth and development model of cucumber in greenhouses, which presents the accumulation and distribution of dry matter in various organs during the period from flowering to fruiting. The NICOLET (NItrate COntrol in LETtuce) model is a well-known and well-researched crop growth model. It is similar to the TOMGRO model in that it represents quantitative relationships between major environmental variables and the growth and development of specific crops, namely lettuce (Lactuca sativa L.) and tomato (Solanum lycopersicum L.). The HORTISIM model is a combined model of seven sub-models including the growth model. The initial NICOLET model was originally designed to predict the nitrate levels of greenhouse lettuce, to combat the health hazard of high nitrates to the population (Seginer et al., 1998). Among these crop models, the NICOLET model, having fewer equations and parameters, is relatively easy to be applied and has been studied by several researchers (Stigter & van Straten, 2000; Mathieu et al., 2006; López-Cruz et al., 2012). Many investigations focused on lettuce models (van Henten, 1994; Escobar-Gutiérrez & Burns, 2002; Zhang et al., 2004, 2008; Shimizu et al., 2008), such as NICOLET, due to the fact that this type of crop has a fast development characteristic and is easy to be managed in greenhouses (Juárez-Maldonado et al., 2012). It can also be calibrated to be applied in other crops, such as cauliflower (Seginer & Stützel, 2006). The NICOLET model has been then calibrated and developed by several researchers since 1990s (Seginer et al., 1999; van Straten et al., 1999; Seginer, 2003, 2004). What’s more, in subsequent researches, it was found that the NICOLET model was also accurate as a tool in predicting the dry and fresh weight of lettuce and had a good application prospect.

In general, the objectives of agricultural crop production are to maximize crop yields and obtain higher quality in order to increase profitability. It requires an appropriate model as a function of the environmental conditions to describe the growth of crops over time, such as solar radiation, maximum and minimum temperatures, relative humidity, carbon dioxide, and cultivar characteristics (Asseng et al., 2014). The NICOLET model is a suitable crop model, whose formulation was intended to satisfy the main experimental observations of general vegetative plants and especially for lettuce (Seginer et al., 2004). The study of such models requires an in-depth understanding of the biological process, and appropriate evaluation is necessary before a model can be successfully applied. Compared with other versions (Seginer et al., 2004) of the NICOLET model which have been studied in different periods, the NICOLET B3 (López-Cruz et al., 2004) contains fewer parameters and formulas, being then simpler to be analyzed and applied. What’s more, most of the studies on the NICOLET have focused on sensitivity analysis and optimal control. However, there are few researches, such as analysis and evaluation of model, on the NICOLET B3 than on other previous versions.

The aim of present study was to evaluate and calibrate a simple dynamic lettuce model, i.e., NICOLET B3, which is capable of predicting the fresh and dry matter of lettuce during its growth period. The novelty of this research is to discuss and identify the feasibility of the studied model applied in predicting lettuce growth with multi-dimensional evaluation criteria. Moreover, the study of the model robustness based on the optimal values of main parameters has been performed with parameter estimation. The possibility of reducing the measuring cost for the application of the NICOLET B3 model was also discussed.

MATERIAL AND METHODS

Description of model

The original NICOLET model (Seginer et al., 1998) is schematically shown in Fig. 1, where the compartments and the carbon fluxes are indicated. The two compartments are designated vacuoles (v), where soluble, non-structural material is stored, and structure (S), which is composed of structural material, including proteins. The various processes are affected by different factors in the environment: Photosynthesis by light and CO₂ concentration, and growth and respiration by temperature. A central element of the model is that a negative linear correlation between primary carbon compounds and nitrate in vacuole which maintain a constant osmotic pressure (Behr & Wiebe, 1992). The NICOLET model has two state variables: M_Cv and M_Cs, which represent carbon content in the vacuoles and cell structure, respectively.

The NICOLET B3 model is based on first principles of plant physiology and has been described in detail by López-Cruz et al. (2004), and only a brief introduction is given here. The core of this model is carbon balance, which is illustrated below:

Carbon in vacuoles

Carbon in structure

The term M_Cv[mol(C)m^-2] originates from assimilation by photosynthesis (F_Cav [mol(C)m^-2s^-1]) which is driven by light (I[μmol/m²/s]) and CO₂ (CC[ppm]). The production of M_Cs[mol(C)m^-2] is driven by growth (F_Cvs[mol(C)m^-2s^-1]). The terms F_Cm[mol(C)m^-2s^-1] and F_Cg[mol(C)m^-2s^-1] denote maintenance respiration and growth respiration, respectively. F_Cm and F_Cvs are controlled by temperature (T[℃]). The term h_g[dimensionless] is a growth inhibition function. It approaches 0 when vacuolar carbon content depletes, and increases asymptotically to 1 with the increase of vacuolar carbon content. There is another function h_p[dimensionless] called photosynthesis inhibition, whose behavior is qualitatively a mirror image of the behavior of h_g. Detailed equations for the NICOLET B3 are summarized in Appendix A [suppl].

The difference, between the NICOLET B3 model and the original NICOLET model, is that there is an extra term (h_g×F_Cm) in balance equations for NICOLET B3 model. This term accounts for the structural carbon in accordance to maintenance respiration, when the level of carbon content in the vacuoles is very low (López-Cruz et al., 2003). Another important difference is that the inhibition functions (h_g and h_p) are not calculated by an exponential function, but are equal to 1 when the concentration of carbon levels is non-inhibiting, that is, the functions have no inhibitory effect. And the reason for the NICOLET B3 modified two original inhibition functions is that they did not completely avoid growth when the level of carbon in the vacuoles is very low (López-Cruz et al., 2003).

Model evaluation method

In this study, three types of evaluation criteria, adopted to evaluate the goodness of fit between real measurement data of lettuce and those predicted by the studied model, are i) coefficient of determination (R²), ii) root mean squared error (RMSE) and iii) relative error.

˗ R² measures how successful the fit is between the observations and predictions (Neter et al., 1996). It is a property of the fitted model, defined as:

where SSR is the sum of squares due to regressions, SST is the sum of squares of total, SSE is the sum of squares for error.

˗ RMSE is a widely used goodness-of-fit measure which can assess bias and estimation accuracy (Harwell, 2018), defined as:

where y_i is the measured value, ŷ_i is the corresponding simulated value, and n is the number of measurements.

˗ Relative error is a measure of the uncertainty of measurement compared to the size of the measurement (Helmenstine, 2020), defined as:

where E_R is relative error, V is the measured value, V' is the corresponding simulated value.

R² approaches 1 (Rios-Moreno et al., 2007) and RMSE approaches 0, both of which are indicating a perfect match between the values observed and predicted by the model.

Methodology for the sensitivity analysis of model inputs

The objective of the sensitivity analysis in this study was to evaluate the relative importance of the input variables of the NICOLET B3. One way to evaluate the studied model sensitivity is calculating the normalized deviation ratio for three input variables of the model, whose influence on the model output is reflected. The normalized deviation ratio (%D) is defined as

where y_k is the original output of the model; ŷ_k is the output of the model with a changed input variable; N is the total number of the data.

Another way to evaluate the model sensitivity to changes in the input variables is the sensitivity function according to France & Thornley (1984).

where S_i is the value of the sensitivity function for a variation of input variable i; Y_dw(T) is the simulated dry or fresh matter at the time T using original inputs; p_i is the original input value; ∂_pi is a small variation in input i while keeping the other inputs constant; p_i and ∂_pi are all normalized; and ∂_Ydw(T) is the difference between the simulated dry or fresh matter at the time T with and without the variation in input i. The larger the S_i of the parameter, the greater the effect of the parameter on Y_dw(T).

Simulation

The simulation was carried out using MATLAB/Simulink R²020b (MathWorks, Inc., 2020), with a variable step-size solver (function ode45). The NICOLET B3 model was programmed in Simulink using photosynthetically active radiation (PAR), air temperature (T) and CO₂ concentration (C_Ca) as model inputs. The experimental data used for this work were obtained from Juárez-Maldonado et al. (2010). The article explains the conditions under which such an experiment was carried out, in which climate and lettuce fresh and dry data were sampled, showing the concentration of CO₂, PAR, air temperature and fresh and dry mass of lettuce. The parameter values of the model were collected from the literature (van Straten et al., 1999; Juárez-Maldonado et al., 2010), and are presented in Appendix B [suppl].

Case study description

The NICOLET B3 model has been used for a case study to predict the lettuce growth trend and its yield. The experiment was carried out at Beijing International Urban Agricultural Science and Technology Park, Tongzhou, Beijing, China, coordinates 116° 48’ N, 39° 52’ W, in a greenhouse of 200 m². Lettuce was cultivated from July 13^th 2020 to August 20^th 2020. The crops were cultivated in Deep Flow Technique (DFT) hydroponic system. The cultivar of lettuce was cream green (´Boston´ lettuce), and each plant was 2 cm apart. There were four planting frames, and each frame could grow 81 lettuces. A water-soluble fertilizer was used to provide the following (in mmol L^-1): NO₃^-, 13.5; NH₄⁺, 1.8; P, 1.3; K, 8.0; Ca, 2.5; Mg, 2.0. The internal environment of the greenhouse is demonstrated in Fig. 2.

With a data logger system connected to a real-time measurement and control management platform in the greenhouse, indoor climate data were recorded every two hours. There were four CO₂ sensors VMS-3002-GZ (Fig. 3a), which were used to measure CO₂ concentration, a sensor LightScout 3415FQF (Spectrum Technologies, Inc., Aurora, USA) PAR (Fig. 3b) and a sensor MZIOT-WS01-W of air temperature with relative humidity (Fig. 3c). The abovementioned sensors were installed inside the greenhouse (Fig. 4). Three lettuce plants were randomly picked every three days to measure the total dry matter. The average values of these measurements were used in the numerical simulation.

Optimization of model

Lopez-Cruz et al. (2004) proposed that only three parameters, maintenance respiration coefficient (k), leaf area closure (a) and growth rate coefficient (v), have a major impact on both state variables and outputs of the NICOLET B3 model. This is the result of a sensitivity analysis which was used to evaluate the effects of the parameters of the model on outputs. It is essential to evaluate and demonstrate model robustness, that is, the sensitivity of empirical results of model parameters to some changes in the model (Young & Holsteen, 2017). And it should be pointed out that each new source of simulation data for the NICOLET model requires a re-adjustment of two or three parameters (Seginer et al., 2004). An accurate crop model could be built up via estimating or calibrating the main parameters of the model. What’s more, it can represent current state of the model better and make the model satisfy different operating conditions more suitable for simulation under diverse experimental environment. So, a calibration process for the main parameters is required.

Parameter estimation can be regarded as an optimization problem. The optimization problem solution is the estimated parameter values set. In this study, maintenance respiration coefficient (k), leaf area closure parameter (a) and growth rate coefficient (v) were the parameters to be estimated. The objective function F(x) of this optimization problem is defined as

where y_sim(i) is the simulated response obtained on the software by tuning the model parameters; y_ref(i) is measured response, namely the experimental measurement data; and n is the number of data.

The process was performed in Simulink by means of least squares identification, which is a recommended method for parameter estimation.

RESULTS AND DISCUSSION

Evaluation of model

The ability of this model to describe the change of lettuce dry and fresh matter was investigated. The outputs of the model were affected by CO₂ concentration, air temperature and PAR in the greenhouse. The input data have been extracted from Juárez-Maldonado et al. (2010). The dry and fresh matter of lettuce solved by Simulink (see Fig. 5) were used to evaluate the NICOLET B3 model compared with measurement data (Juárez-Maldonado et al., 2010). According to Fig. 5, the simulated dry and fresh matter followed quite accurately the general dynamic trend in the mean of the measurement data of crop growth (van Holsteijn, 1980). Simultaneously, the curve change trend of dry weight (DW) and fresh weight (FW) was basically the same, because the ratio of DW and FW is stable at 0.059. This ratio corresponds to the constant DFR in the above model formula in Eq (9),

where the value of DFR is 0.05 in van Straten et al. (1999).

The comparisons between the simulation results in this study, real measurement data of lettuce crop and its corresponding model results from the literature (Juárez-Maldonado et al., 2010) are demonstrated in Figs. 6a (dry matter) and 6b (fresh matter), respectively. It can be seen that the simulation results of this study were very close to both from the literature (Juárez-Maldonado et al., 2010). Therefore, both of them could ensure that the studied model in this work has the credibility for the following investigation.

Figs. 7a and 7b demonstrate the values predicted by the studied model fit the collected data of dry and fresh matter (Juárez-Maldonado et al., 2010) very well, respectively. For dry matter weight, the coefficient R² (0.9844) of the studied model in this work is better than the coefficient from the literature (Juárez-Maldonado et al., 2010), which was equal to 0.974. It indicates that the studied model has a better prediction accuracy for DW of lettuce growth compared with Juárez-Maldonado’s. Regarding fresh matter weight, the correlation coefficient R² (0.9940) shows almost perfect fitting relation, which is very close to the result of 0.9975 in the literature (Juárez-Maldonado et al., 2010).

For the average values of relative errors between simulation results from the studied model and real data from Juárez-Maldonado et al. (2010), the error for DW data was 14.93% and for FW data 11.84%. At the same time, for the simulated results in Juárez-Maldonado et al. (2010) and real data, the average values of relative errors were 15.14% and 11.02%, respectively. The simulated lettuce weight was similar to real measurement data before harvest, with a RMSE = 0.3698 kg/m² for FW and a RMSE = 0.0803 kg/m² for DW. The simulated FW and DW in 39 days (13.68 kg/m² and 0.82 kg/m², respectively) were about 1:37 and 1:10 of the values of RMSE, respectively. It’s enough to prove the good reliability of the model simulation.

The NICOLET B3 model has been evaluated in a multi-dimensional way using three criteria. All these statistics suggest a good fit for lettuce growth prediction, which could fully prove that the NICOLET B3 model is effective and accurate for predicting lettuce yield and taking decisions concerning the harvest date selection. The abovementioned results are summarized in Table 1.

Sensitivity analysis of model inputs

After the NICOLET B3 model was evaluated and its credibility ensured, it was then further evaluated. The variation characteristics in values of lettuce fresh and dry matter were tested when model inputs change (Figs. 8a-c), where the blue line corresponds to the data measured within the greenhouse, meanwhile, the red dotted line considers a constant value of the average of the input data.

Through the two methods mentioned above to measure the model’s input sensitivity, the calculated results are as follows.

In the case of taking the average of the input temperature, the deviation degree of the model output from the original fresh matter and dry matter reached 20.29% and 0.93%, respectively. And the deviation degree was only 0.01% and 0.02% for fresh matter, both 0.01% for dry matter, when PAR level and CO₂ concentration were averaged. This difference can also be seen from the ending points of the polylines in Figs. 8a, 8b and 8c. The sensitivity of the model to the three inputs was compared using %D as an evaluation criterion. The analysis found that the effect of temperature on the model was much more obvious than other two inputs.

The average results of S_i were calculated over the planting days for fresh and dry matter. The results (Table 2) indicate that temperature was the most significant factor determining the growth of lettuce among three inputs of the NICOLET B3. Taking the result of FW as an example, compared with PAR and CO₂, S of T increased by 1329.7 ((2.1291-0.0016)/0.0016=1329.7) times and 7096 ((2.1291-0.0003)/0.0003=7096) times, respectively, which shows that the model sensitivity for variations in the other two inputs seems to be much smaller. Additionally, under the same input, S for fresh matter was always lower than S for dry matter. The S of CO₂ was minimal regardless of the DW or FW. In other words, CO₂ concentration had the least impact on the model. The aforementioned results just conform to the analysis and conclusions of normalized deviation ratio.

Based on both previous evaluation results for the studied model sensitivity, it can be thus concluded that temperature as input variable has the most obvious influence on the output of the NICOLET B3 model, i.e., dry and fresh matter. In view of this conclusion, the following research was conducted to simplify the input variables of the NICOLET B3 model when those input variables maintain within the ranges of crop growth requirement. That is, one or two input values were retained, and the remaining input values were derived from the average of the daily measurements to observe the error of the model outputs.

The inputs of the studied model were further optimized by comparing the model performance of one or two excluded input variables, respectively (Table 3). The influence of reducing different input variables in the model was different. It can be seen from the δ_r-m of both retaining two inputs and of all the inputs, that temperature was the most important input on the output of the model. Although the model also performed well when temperature was the only input of the model, however, too few input variables may make the outputs of the model uncertain and undermine the principles of plant physiology.

A larger number of inputs requires that more couplings between inputs and outputs be described, which usually leads to more complex models. And it is more difficult for the model to be applied and modified for control purposes. So, the temperature and the PAR, or the temperature and CO₂ concentration can be used as the inputs of the NICOLET B3 model to simplify the input variables when simplified input maintain within the ranges of crop growth requirement. The simplified model input was only replaced by the average of the normal value of the greenhouse, eliminating the need for precise measurement. What’s more, this simplified model could eliminate the need to measure an environment variable in the application, which helps to facilitate the application of the NICOLET B3 to a certain extent and decrease the measuring cost without affecting significantly the capability of the model.

Results of the case study

The measurement results of average CO₂ concentration, PAR and air temperature inside the greenhouse during the experiment from 13 July to August 20, 2020, are shown in Figs. 9a, 9b and 9c, respectively.

The growth trend of lettuce under greenhouse conditions in Tongzhou was established by the fresh and dry matter predicted by NICOLET B3 model. Fig. 10 shows the growth curves of FW and DW. They are consistent with the general trend of fresh and DW in the previous lettuce experiment. The ratio of DW to FW is also a fixed value, stable at about 0.059. Using the NICOLET B3 model, the growth of lettuce in the greenhouse can be well predicted, which is of benefit to the efficiency of the management of the greenhouse climate and yield maximization.

The aforementioned three main parameters were estimated to fit the real data of the lettuce from the literature (Juárez-Maldonado et al., 2010). Table 4 shows the values of three parameters before and after calibration. The reason for the difference between the original and the calibrated values is that the original values were calibrated using data from the lettuce grown in Germany (van Straten et al., 1999) and the Netherlands (Linker et al., 2004). However, the greenhouse experimental conditions used in this case study were completely different. It can be seen from this result that the robustness of model parameters is poor. Some main parameters of a nonlinear dynamic model may not be accurate for all scenarios due to the influence of environment, resulting in a mismatch problem between the model and specific application environment. Therefore, when the experimental environment and crop characteristics change greatly, it is better to readjust the main parameters of the model, which will help to fit correctly the plant growth experimental data from very different treatments.

For the remaining parameters, the currently adopted parameter values are sufficient for the application of the model. From the experimental results of the NICOLET model applied in three different places, namely Germany (van Straten et al., 1999), the Netherlands (Linker et al., 2004) and Mexico (Juárez-Maldonado et al., 2010), with different varieties of lettuce, it can be seen that the model has good performance by using the existing model parameter values. Some parameters were corrected according to the season (van Straten et al., 1999) and the current experimental environment (Juárez-Maldonado et al., 2010). Then, when the NICOLET model is applied in different environments, some key parameters can be corrected on the basis of known model parameter values, which is helpful to improve the robustness of the model. The method of calibration of parameters described above is a simple and feasible way.

Final considerations

This study has presented that the NICOLET B3 model can predict the general trend in fresh and dry matter of the lettuce with satisfactory accuracy. However, there is a big discrepancy between the sim ulated results and the real data observed at the end of lettuce growth curve, especially on the 43^rd day of lettuce growth. The reason for the discrepancy may be partly attributed to the limitations of the model simulation. The output data on 43^rd day were estimated by interpolation through input data at other days. An absence of input data for the model after 43^rd day would affect the estimation of that day. And it was also found in the test that the output value of this dynamic model was affected by the input data before and after its corresponding time. For example, when the input data of the model only ended on the 39^th day, there was also a large discrepancy between the simulation results and the real data of the fresh and DW of lettuce on that day. Therefore, in the previous experiment, the inputs data after the 39^th day were not used for the model simulation, and there was also due to no comparison of the measured lettuce data after that day. Greenhouse data collection for more prolonged periods of time will contribute to configure more representative climate files of the study site (Bojacá et al., 2009).

The evaluation of a model’s ability to describe the dynamic evolution of crop is also important, and the crop dry and FW is determined by the model inputs (van Henten, 1994). Temperature as a model input produced equally positive results, that is, the relative error between the predicted results of the simplified model inputs and the real data is close to that of the original model. This is because the other two input variables have less influence on the model than temperature. This behavior was similar to that reported by Juárez-Maldonado et al. (2012), who also obtained that temperature have the greatest impact on NICOLET model.

A new calibration is required to obtain an accurate prediction of lettuce growth when the NICOLET B3 model is simulated and applied in different locations and scenarios. The variations in parameters are mainly due to factors such as the area in which the experiment was carried out and the characteristics of the variety, as well as the species (Quesada-Roldán & Bertsch-Hernández, 2013). The model may need to be applied to different locations in the future to validate the parameter estimation method. In addition, the data assimilation method, such as unscented Kalman filter (UKF), can also be used to estimate parameters and improve the prediction performance of the NICOLET model by incorporating information coming from samples of destructive measurements of actual crop (Ruíz-García et al., 2014).

Despite of the limitations already mentioned, the results of the present work provide insight into the effectiveness of the NICOLET model for predicting yield of lettuce. The dynamic model initially applied to the tomato crop was calibrated and validated to adequately simulate the growth and mineral uptake of cucumber (Ramírez-Pérez et al., 2018). In consequence, it is possible to study the NICOLET model is applied to other crops with a similar approach.

In this paper, the knowledge about the mechanism of the NICOLET B3 model and the data information has been collected to analyze this nonlinear dynamic crop model. The mechanism of the model accords with plant physiology. The model is validated with multi-dimensional evaluation criteria. The correlation between simulated results from the proposed model and real data from the literature was very well and close to 1, with R² = 0.9939 for fresh matter and R² = 0.9858 for dry matter. It demonstrated that the NICOLET B3 model can predict the behavior of lettuce growth with a convincing degree of accuracy.

A sensitivity analysis has been performed to analyze the influence that three inputs may have on the prediction of lettuce yield. Model simulations indicate that the growth of the lettuce is less sensitive to the changes of PAR and CO₂ as inputs compared to temperature. This might indicate that from the control point of view, the numbers of the input variables in the NICOLET B3 model can be reduced. Test results demonstrate that the model is still effective after reducing PAR or CO₂ from input variables. This reduced model will help to decrease the measuring cost without affecting significantly the capability of the model to predict the dynamic growth and yield of the crop.

Meanwhile, The NICOLET B3 model was applied to a greenhouse in Tongzhou, Beijing. A method, which provides a solution for the mismatch between the model and the specific application environment, to modify the model to apply in different scenarios is also proposed, which help to improve the reliability and applicability of the model.

Further research is needed to test the performance of the NICOLET B3 model under different climate conditions in the greenhouse. At the same time, some model parameters should be calibrated according to the aforementioned different situations. The NICOLET model was originally developed to predict nitrate content of a greenhouse lettuce crop. Nitrate as one of the outputs of the model should be brought into the study, which also helps control high concentrations of nitrate in leafy vegetables constituting a health hazard for people. At the same time, nitrogen is a key element in the aquaponics system, the investigation of the model is thus valuable in this field.