Crop Management in Controlled Environment Agriculture (CEA) Systems Using Predictive Mathematical Models †

Proximal sensors in controlled environment agriculture (CEA) are used to monitor plant growth, yield, and water consumption with non-destructive technologies. Rapid and continuous monitoring of environmental and crop parameters may be used to develop mathematical models to predict crop response to microclimatic changes. Here, we applied the energy cascade model (MEC) on green- and red-leaf butterhead lettuce (Lactuca sativa L. var. capitata). We tooled up the model to describe the changing leaf functional efficiency during the growing period. We validated the model on an independent dataset with two different vapor pressure deficit (VPD) levels, corresponding to nominal (low VPD) and off-nominal (high VPD) conditions. Under low VPD, the modified model accurately predicted the transpiration rate (RMSE = 0.10 Lm−2), edible biomass (RMSE = 6.87 g m−2), net-photosynthesis (rBIAS = 34%), and stomatal conductance (rBIAS = 39%). Under high VPD, the model overestimated photosynthesis and stomatal conductance (rBIAS = 76–68%). This inconsistency is likely due to the empirical nature of the original model, which was designed for nominal conditions. Here, applications of the modified model are discussed, and possible improvements are suggested based on plant morpho-physiological changes occurring in sub-optimal scenarios.


Introduction
Environmental control is a key factor to increase plant productivity in controlled environment agriculture (CEA) [1]. Recently, increased interest has been directed towards plant production in closed facilities (e.g., plant factories, vertical farms, indoor-growing modules) [2][3][4]. However, with the introduction of advanced monitoring and control technologies, it becomes necessary to properly discern plant/microclimate interaction, modulate environmental parameters, and manage cultivation factors. Indeed, in protected agriculture, environmental factors and plant responses are strictly interconnected: alteration of the microclimate can induce modifications in plants (both at morphological and physiological levels), affecting plant behavior especially in terms of transpiration and CO 2 /O 2 exchanges, which in turn re-modify the surrounding environment [5]. energy balance in closed systems. Moreover, such a model could easily be applied by stakeholders operating in controlled agriculture facilities for food production. Indeed, the model, requiring just a few environmental and cultivation parameters as inputs, can help forecasting changes happening during the cultivation after the modification of environmental factors, thus being possibly implemented as decision support system.
In the present study, we report an application of the MEC model on butterhead lettuce (Lactuca sativa L. var. capitate) cultivated under controlled conditions, also presenting a modification to the original model. More specifically, we introduced additional components to include the variation of canopy quantum yield of PSII (CQY) and carbon use efficiency (CUE) during leaf development. Such parameters, which represent respectively the moles of carbon fixed per mole of photons absorbed and the ratio of net carbon gain to gross carbon assimilation during growth, are critical to the model for the calculation of the net photosynthesis and the biomass production. Furthermore, we applied the modified version of the model on green-and red-leaf lettuce grown in a climatic chamber under two VPD scenarios, namely low and high VPD, corresponding to nominal and off-nominal conditions, respectively. This latter trial allowed the portrayal of the model parameters for nominal and off-nominal scenarios for green-and red-leaf plants, since non-identical behavior often occurs for different cultivars/varieties even under the same growth conditions.

The Original MEC Model
The original energy cascade model was developed for wheat by Volk et al. [25] and then calibrated for other crops like lettuce, rice, soybean, sweet potato and tomato [26]. It is an explanatory model, composed by multivariate equations whose coefficients have been determined through curve fitting of experimental data [27]. The input variables of the first version of the model were the light intensity (photosynthetic photon flux density; PPFD) and the photoperiod. All the model outputs, which mostly concerned the biomass and the growth rate, were function of these two parameters.
The original MEC model was based on an "energy cascade" with three fundamental steps: (1) the absorption of PPFD by the canopy; (2) the absorbed energy (A) used in the photosynthetic process to convert carbon into sucrose; (3) the conversion of sucrose into biomass. Such a simple model required only three crop parameters, namely (i) the time of canopy closure (t A ), (ii) the time of senescence onset (t Q ), and (iii) the time of harvesting (t M ). Eventually, the model was modified to add the following climatic parameters: air temperature, relative humidity, carbon dioxide concentration, dark period, and plant density, in order to improve the accuracy and robustness of the model [16]. In 2012, Boscheri et al. [12] implemented a modified version of the MEC model, for a multi-crop Lunar greenhouse prototype. The version modified by Boscheri et al. [12] also included a crop transpiration component, used to predict water and plant nutrient consumption.
The main model algorithm components were arranged according to twelve equations, sequentially computed at each time step to calculate the key variables as listed below.
The canopy quantum yield (CQY, mol −1 ) is defined by an empirical equation as function of the time t, expressed in days after emergence: where CQY MAX and CQY MIN are crop-specific parameters, while t Q and t M are time of the onset of senescence and time of harvesting, respectively. The carbon use efficiency (CUE), is expressed similarly to CQY, according to the following relation: where CUE MAX and CUE MIN are crop specific parameters, although CUE has been often assumed to be constant for most crops (i.e., CUE MAX = CUE MIN ), t Q and t M are time of the onset of senescence and time of harvesting, respectively. The parameter A is the fraction of PPFD absorbed by the top of the canopy and is assumed to increase with time t according to a power law equation, up to a maximum value A MAX at time t = t A , when canopy closure is established: where t A is the time of canopy closure and n is a crop dependent exponent, which is considered to be equal to 2.5 for lettuce [27]. The daily carbon gain (DCG, mol C m −2 d −1 ) is computed as follows: where H is the photoperiod, α = CUE·A·CQY, and 0.0036 is a constant used to convert µmol to mol and hours to seconds. The daily oxygen production (DOP, mol O 2 m −2 d −1 ) is then given by a fraction (oxygen production fraction, OPF) of DGC: where OPF is expressed in mol O 2 mol −1 C and is a crop-specific parameter. The crop growth rate (CGR, g m −2 d −1 ) is given by: where MWc = 12 g mol −1 is the carbon molecular weight, while BCF is a crop specific parameter representing the biomass carbon fraction. Thus, the total edible biomass (TEB), expressed as specific dry weight (g m −2 ), was calculated by integrating the crop growth rate, multiplied by the fraction of CGR allocated to edible biomass (XFRT): where t M is the time of harvesting, t E is the time of the onset of edible biomass formation and XFRT represents a partitioning coefficient for the edible biomass, which combines the effects of determinacy and temperature on storage organ growth rates [26]. The gross photosynthesis (P G , µmol CO 2 m −2 s −1 ) is computed as follows: where β = A·CQY and PPFD is the photosynthetic photon flux density (µmol m −2 s −1 ). The net photosynthesis (P N , µmol CO 2 m −2 s −1 ) is computed by accounting for the carbon use efficiency in the photoperiod: Sensors 2020, 20, 3110 where H is the photoperiod, α = CUE·A·CQY, β = A·CQY and PPFD is the photosynthetic photon flux density (µmol m −2 s −1 ). The stomatal conductance (g S , mol m −2 s −1 ) for planophile-type canopies (such as lettuce) is calculated according to Monje et al. [28] where T ( • C) is the mean air temperature during light cycle and [CO 2 ] is the air carbon dioxide concentration expressed as µmol CO 2 mol −1 .
The canopy surface conductance for water vapor (g c , mol m −2 s −1 ) is defined as follows: where g A = 2.5 mol m −2 s −1 is the aerodynamic conductance and g s the stomatal conductance. The daily canopy transpiration (DTR, L m −2 d −1 ) is also calculated as follows: where 3600 is a conversion constant from second to hours, MW W = 18 g mol −1 is the molecular weight of water, ρ W = 100 g L −1 is the water density, and P ATM (kPa) is the total atmospheric pressure, which was used to convert vapor pressure to mole fraction.

Limitations of the Original MEC Model Formulation
According to Equations (1)- (12), the original energy cascade model, used for advanced life support systems (ALSs) studies, predicted the biomass production and the photosynthetic rate based on three parameters: (i) the canopy light absorption (A), (ii) the crop quantum yield of PSII (CQY), and (iii) the carbon use efficiency (CUE). The physical and biological trends of these parameters were: a linear increase in PPFD absorption till the canopy closure; a constant CQY until the onset of senescence, followed by a linear decrease till the end of the cycle, and a constant CUE throughout the life cycle ( Figure 1, orange lines). However, for lettuce and a few other crops like sweet potato, which are harvested before the occurrence of senescence (t Q ), the CQY, as the CUE, were assumed to be constant during the entire growth cycle prior harvesting ( Figure 1, blue lines).
where gA = 2.5 mol m −2 s −1 is the aerodynamic conductance and gs the stomatal conductance. The daily canopy transpiration (DTR, L m −2 d −1 ) is also calculated as follows: where 3600 is a conversion constant from second to hours, MWW = 18 g mol −1 is the molecular weight of water, ρW = 100 g L −1 is the water density, and PATM (kPa) is the total atmospheric pressure, which was used to convert vapor pressure to mole fraction.

Limitations of the Original MEC Model Formulation
According to Equations (1)- (12), the original energy cascade model, used for advanced life support systems (ALSs) studies, predicted the biomass production and the photosynthetic rate based on three parameters: (i) the canopy light absorption (A), (ii) the crop quantum yield of PSII (CQY), and (iii) the carbon use efficiency (CUE). The physical and biological trends of these parameters were: a linear increase in PPFD absorption till the canopy closure; a constant CQY until the onset of senescence, followed by a linear decrease till the end of the cycle, and a constant CUE throughout the life cycle (Figure1, orange lines). However, for lettuce and a few other crops like sweet potato, which are harvested before the occurrence of senescence (tQ), the CQY, as the CUE, were assumed to be constant during the entire growth cycle prior harvesting ( Figure 1, blue lines). for lettuce (blue lines), which was considered to be constant throughout the whole crop cycle; and for other crops (orange lines), in which CQY (a) was considered to be constant until the onset of senescence (tQ), to linearly decrease till the time of harvesting (tM), while CUE (b) was considered to be constant.

Experiments to Retrieve CQY Temporal Pattern
The canopy quantum yield (CQY) represents the moles of carbon fixed per mole of photons absorbed [29][30][31] and can be assessed in different ways such as: (1) Dividing the daily PG (mol C m 2 ·d -1 ) by the total absorbed photons (mol m 2 ·d -1 ) [32,33]; (2) From the initial slope of saturation-photosynthetic curves [29]; (3) By means of a fluorimeter, an instrument which measures the proportion of the light absorbed by the chlorophyll associated with the photosystem II (PSII), thus indicating the efficiency of the carbon fixation and of the overall photosynthesis [34].

Experiments to Retrieve CQY Temporal Pattern
The canopy quantum yield (CQY) represents the moles of carbon fixed per mole of photons absorbed [29][30][31] and can be assessed in different ways such as: (1) Dividing the daily P G (mol C m 2 ·d -1 ) by the total absorbed photons (mol m 2 ·d -1 ) [32,33]; (2) From the initial slope of saturation-photosynthetic curves [29]; (3) By means of a fluorimeter, an instrument which measures the proportion of the light absorbed by the chlorophyll associated with the photosystem II (PSII), thus indicating the efficiency of the carbon fixation and of the overall photosynthesis [34].
In this study, we designed a calibration experiment for assessing the CQY temporal pattern in butterhead 'Salanova' lettuce (Lactuca sativa L. var. capitata), by means of a fluorimeter. The experiment was performed from 10 July 2019 to 6 August 2019, at the Department of Agricultural Sciences (University of Naples Federico II). Two weeks after sowing, the plants of butterhead lettuce cultivars, green-and red-leaf Salanova ® (Rijk Zwaan, Der Lier, The Netherlands) were transplanted into 10 cm pots filled with a peat:soil mixture (1:1 volume ratio) and exposed to solar light. Air temperature (T) was kept at 23 • C and plants were irrigated at 2-day intervals in order to reach the container capacity (till the beginning of drainage). Chlorophyll "a" fluorescence analyses were performed by means of a portable fluorimeter (ADC Bioscientific) on 3 expanded leaves per 4 green-and 4 red-lettuce plants every day, in order to highlight the leaf-age-driven variations in CQY, according to Genty et al. (1989) [35]. Fluorescence analyses were conducted at steady-state photosynthesis under a light intensity of about 400 µmol m −2 s −1 , with a saturation pulse duration of 0.8 s, by keeping the orientation of the leaf relative to the actinic light source when taking CQY measurements. The leaves were chosen from three different positions within the lettuce head (top-t, medium-m, and bottom-b) in order to gain representative data for retrieving the general temporal pattern of CQY throughout the canopy.

Experiments in a Controlled Environment Growth Chamber to Validate the Model
Nine green-and nine red-leaf Salanova lettuce plants were grown in a controlled environment growth chamber (KBP-6395F, Termaks, Bergen, Norvegia) ( Figure 2), in two consecutive trials. In the first trial, 1-week old plants were transplanted into 10 cm pots filled with peat:perlite substrate (1:1 volume ratio) and incubated with an average VPD of 0.69 kPa, corresponding to nominal conditions of Low VPD. In the second trial, 9 green-and 9 red-leaf lettuce plants were incubated with an average VPD of 1.76 kPa, corresponding to off-nominal conditions of High VPD. The two different VPDs were achieved by keeping air temperature (T) constant at 24 • C, while changing the relative humidity (RH) accordingly. Temperature and RH were monitored and recorded every 10 min by means of mini sensors equipped with a data logger (Testo 174H). All other microclimate parameters and agricultural practices were the same in the two consecutive experiments. The lighting system was an RGB LED panel, with a light intensity of 315 PPFD µmol m −2 s −1 at the canopy level (12 h photoperiod; 13.6 daily light integral, DLI). Daily rotation of the trays was performed to ensure homogenous light and humidity across the shelf surface. Plants were daily weighted to assess the loss of water by transpiration (DTR) and were re-watered to field-capacity. Evaporation losses from the substrate were minimized by covering the substrate with a plastic film. Plant growth was assessed by imaging, measuring canopy total area every day, and counting the number of leaves. Furthermore, dry weight was recorded at the beginning and at the end of both trials. These measurements were used to reconstruct the daily total edible biomass (TEB). Changes in leaf temperature were monitored with an infrared thermometer on three leaves per plant (H-1020; Helect). These measurements were averaged and used instead of T in Equation (10) of the MEC model, to obtain more precise information of the leaf-to-air VPD which influences stomata conductance the most. After 23 days, on fully developed leaves, eco-physiological analyses in terms of gas exchanges (LCA 4; ADC BioScientific Ltd., Hoddesdon, UK) and chlorophyll "a" fluorescence (through the above-reported portable fluorimeter), were performed. Gas-exchange analyses were carried out on fully expanded leaves, using 9 replicates per condition to assess plants physiological behavior (P N and g s ) in response to different VPD conditions. During the measurements PAR, RH and carbon dioxide concentrations were set at ambient value and the flow rate of air was set to 400 mL s −1 . P N and g s values were averaged and also used to evaluate the corresponding model prediction performances.

Model Structure and Parameter Identification
The results of the experiments conducted to retrieve CQY temporal pattern suggested the opportunity to redesign the MEC model by modifying the temporal variability of some key parameters, without increasing the model complexity. Indeed, fundamental model outputs like TEB and PN, directly depend on these key parameters, which temporal patterns are influenced by the VPDs levels (nominal and off-nominal) and often result to be cultivar-specific [13]. Based on Equations (1)- (12), the three variables A, CQY, and CUE were reduced to the variable α = CUE·A·CQY and β = A·CQY. The developmental stages observed for CQY were then assumed to be valid also for α and β, i.e., under the assumption of functional similarity between the variables involved, as also suggested by other studies [36][37][38][39]. The key temporal parameters, explaining the different development for α and β, were set equal to those observed for CQY. The other parameters, identifying the minimum and maximum values for α and β, were calibrated by minimizing the root mean square error (RMSE) of the predicted DTR with respect to the corresponding observations over the entire simulation period, with the generalized reduced gradient optimization algorithm. The calibrated model was then validated against TEB measured over the entire simulation period. Simulated and modelled gs and PN were also compared on day 23 after transplanting (DAT), the day of the experiment when gas exchange analyses were performed to experimentally determine PN and gs. Calibration and validation were performed for each of the four examined scenarios: green lettuce under nominal VPD conditions (G-N); green lettuce under off-nominal conditions (G-ON); red lettuce under nominal VPD conditions (R-N); red lettuce under off-nominal conditions (R-ON).
The statistical performance indices for CQY, DTR and TEB were the linear correlation between predictions and observations (r), the average difference between prediction and observation (BIAS), the root mean square error (RMSE) and the ratio of performance to deviation (RPD). BIAS and RMSE were computed as follows: where N is the number of simulation days considered for the computation of the performance index, X , is the prediction on the j-th day, X o,j is the observation day on the j-th day, averaged among the 9 sample plants for each examined scenario. Meanwhile, the RPD was calculated as the ratio of the

Model Structure and Parameter Identification
The results of the experiments conducted to retrieve CQY temporal pattern suggested the opportunity to redesign the MEC model by modifying the temporal variability of some key parameters, without increasing the model complexity. Indeed, fundamental model outputs like TEB and P N , directly depend on these key parameters, which temporal patterns are influenced by the VPDs levels (nominal and off-nominal) and often result to be cultivar-specific [13]. Based on Equations (1)-(12), the three variables A, CQY, and CUE were reduced to the variable α = CUE·A·CQY and β = A·CQY. The developmental stages observed for CQY were then assumed to be valid also for α and β, i.e., under the assumption of functional similarity between the variables involved, as also suggested by other studies [36][37][38][39]. The key temporal parameters, explaining the different development for α and β, were set equal to those observed for CQY. The other parameters, identifying the minimum and maximum values for α and β, were calibrated by minimizing the root mean square error (RMSE) of the predicted DTR with respect to the corresponding observations over the entire simulation period, with the generalized reduced gradient optimization algorithm. The calibrated model was then validated against TEB measured over the entire simulation period. Simulated and modelled g s and P N were also compared on day 23 after transplanting (DAT), the day of the experiment when gas exchange analyses were performed to experimentally determine P N and g s . The statistical performance indices for CQY, DTR and TEB were the linear correlation between predictions and observations (r), the average difference between prediction and observation (BIAS), the root mean square error (RMSE) and the ratio of performance to deviation (RPD). BIAS and RMSE were computed as follows: where N is the number of simulation days considered for the computation of the performance index, X p,j is the prediction on the j-th day, X o,j is the observation day on the j-th day, averaged among the 9 sample plants for each examined scenario. Meanwhile, the RPD was calculated as the ratio of the standard deviation of the experimental data to the standard error of the model predictions. Although considered redundant by some authors [40], it has been suggested that RPD > 3 are good for screening purpose, RPD > 5 are suitable for quality control and RPD values > 8 are considered optimal for every kind of analytical application [41]. The model performance for P N and g s were assessed by the relative BIAS, i.e., the BIAS normalized by the corresponding average measured value on the 23rd DAT.

Model Equations and Parameters
Results from fluorescence analyses on green-and red-leaf lettuce plants are reported in Figure 3. In both cultivars, top leaves always presented the highest values followed by medium leaves, while the lowest values were recorded in the bottom leaves. Based on the analyzed experimental data, in green-leaf plants, we distinguished three stages in CQY temporal patterns: (1) A period of CQY monotonically increasing, starting from the initial leaf lamina development till the beginning of the maturity stage (t Mi ). As shown in Table 1, predicted CQY BIAS was close to zero, while RMSE varied from 0.028 to 0.072 for green-and 0.020-0.021 for red-leaf lettuce. However, the linear correlation between observed and predicted CQY was high in all conditions, always being larger than 0.92. Furthermore, Furthermore, for the red-pigmented lettuce, the first phase was preceded by a period of stationary CQY.
In the light of these considerations, instead of Equation (1), we designed a new mathematical relation to analytically describe the temporal evolution of CQY, consistently with the experimental data. The first stage was modelled by a linear equation, that fits the observed data between the initial time of development (t D ) and maturity (t Mi ). The second stage was modelled by constant value, with CQY = CQY MAX , from t Mi till end of maturity (t M ). The third stage was modelled by a linear decreasing equation, between t Mi and the time t S , corresponding to the end of the senescence period, with CQY = CQY S . The third stage is not relevant from a practical perspective, since it is beyond the time of harvesting, which coincides with t M . Moreover, for the red-leaf cultivar, the initial period of stationary CQY was modelled by a constant value, with CQY = CQY MIN , from the first day after transplanting till t D . This stage is only relevant for the red lettuce, with t D = 8 days, since green lettuce does not present this "adaptation" period (i.e., t D = 1). As illustrated in Figure 3, the other relevant times t Mi and t M resulted to be equal to 16 and 23 for all examined scenarios and leaves.
As shown in Table 1, predicted CQY BIAS was close to zero, while RMSE varied from 0.028 to 0.072 for green-and 0.020-0.021 for red-leaf lettuce. However, the linear correlation between observed and predicted CQY was high in all conditions, always being larger than 0.92. Furthermore, the values for RPD always showed values around 5, indicating the robustness and reliability of CQY prediction model. Table 1. BIAS, root mean square error (RMSE), linear correlation (r) and ratio of performance to deviation (RPD) for CQY of PSII of green-(G) and red-leaf (R) lettuce cultivars for top-t, medium-m, and bottom-b leaves. In accordance with the assumptions stated in Section 2.5, variables α and β were modelled to change in time according to the following relation:

BIAS RMSE r RPD
where χ denotes the generic variable (either α or β), while χ MIN and χ MAX the corresponding minimum and maximum values. Parameters α min , α max, β min and β max were derived from experimental gas-exchange and chlorophyll "a" fluorescence measurement performed in the growth chamber experiment. These parameters were therefore differentiated for the nominal and off-nominal scenarios and for green and red lettuce cultivars, since these cultivars showed different behaviors under the same growth conditions. Afterwards, these parameters were calibrated by minimizing the RMSE with respect to the measured DTR values.

Model Performance
From our results, it is evident how DTR, TEB, g S , and P N followed similar trends in green-and red-leaf lettuce cultivars grown under both nominal and off-nominal conditions, but with different absolute values and magnitude of changes during the cultivation cycles. Figure 4 shows the observed and predicted DTR values after model calibrations. The irregular pattern of DTR with time, is due to high sensitivity of DTR to slight perturbations in VPD levels during the diurnal hours of experiment. This sensitivity was higher under off-nominal conditions. Indeed, the model was able to reproduce the observed DTR under nominal conditions better than under off-nominal conditions. As reported in Table 3, under nominal conditions, predicted DTR BIAS was almost null, while RMSE was 0.09 L m −2 and 0.10 L m −2 for green-and red-leaf lettuce, respectively, i.e., less than 20% of the average DTR observed during the entire experiment. Under off-nominal condition, DTR BIAS was still low (max 0.04 L m −2 for red-leaf lettuce) but the RMSE increased to 0.21 L m −2 for the green lettuce and to 0.35 for the red lettuce. These RMSE values were still acceptable, since they are below the 30% of the observed DTR. The linear correlation between observed and predicted DTR was almost always high (larger than 0.70), except for the R-ON scenario which exhibited a linear correlation equal to 0.56. However, the RPD values varied from 5.10 to 6.40 for nominal conditions and 5.25 to 5.34 for off-nominal conditions. Being always higher than 5, the RPD overall indicates a good quality of the model predictions.
The calibrated model was validated with the observed total edible biomass (TEB), representative of the lettuce daily growth. Total edible biomass was also influenced by VPD conditions: plants under nominal condition developed more biomass than those grown under off-nominal scenarios, both in green-and red-leaf lettuce cultivars ( Figure 5). The temporal evolution of TEB was more regular and less sensitive to the perturbations of the experimental settings. The TEB predictions curves fully reflected what was expected by lettuce grown under those different VPD conditions, showing an almost linear increment in biomass till the time of harvesting ( Figure 5). Furthermore, the predicted growth curves accurately simulated the lettuce biomass accumulation. In Table 4 Figure 4 shows the observed and predicted DTR values after model calibrations. The irregular pattern of DTR with time, is due to high sensitivity of DTR to slight perturbations in VPD levels during the diurnal hours of experiment. This sensitivity was higher under off-nominal conditions. Indeed, the model was able to reproduce the observed DTR under nominal conditions better than under off-nominal conditions. As reported in Table 3, under nominal conditions, predicted DTR BIAS was almost null, while RMSE was 0.09 L m −2 and 0.10 L m −2 for green-and red-leaf lettuce, respectively, i.e., less than 20% of the average DTR observed during the entire experiment. Under offnominal condition, DTR BIAS was still low (max 0.04 L m −2 for red-leaf lettuce) but the RMSE increased to 0.21 L m −2 for the green lettuce and to 0.35 for the red lettuce. These RMSE values were still acceptable, since they are below the 30% of the observed DTR. The linear correlation between observed and predicted DTR was almost always high (larger than 0.70), except for the R-ON scenario which exhibited a linear correlation equal to 0.56. However, the RPD values varied from 5.10 to 6.40 for nominal conditions and 5.25 to 5.34 for off-nominal conditions. Being always higher than 5, the RPD overall indicates a good quality of the model predictions.     growth curves accurately simulated the lettuce biomass accumulation. In Table 4, BIAS, RMSE, r and RPD for the TEB are reported. The linear correlation coefficient (r) was always close to 1, BIAS varied in the range 0.19-1.11 under nominal conditions, and 0.12-0.40 under off-nominal conditions, whereas RMSE and RPD values varied in the range 4.46-6.87 and 4.87-5.00 under nominal conditions and 2.98-3.60 and 4.88-4.96 under off-nominal conditions, indicating a good reliability of model predictions. Overall, the results showed that the TEB prediction errors were always below 10%.   Figure 6 shows the box-plot distribution of the stomatal conductance (gS) as well as of the net photosynthesis (PN). Stomatal conductance and PN were significantly higher under nominal conditions compared to off-nominal, in both butterhead cultivars. Unfortunately, these data are only available for the 23rd DAT, when leaf gas-exchange analyses were performed. As illustrated in Table  5, gS relative BIAS (rBIAS) was equal to 39% and −0.1% in green-and red-leaf lettuce, respectively; PN relative BIAS (rBIAS) was instead equal to 34.1% and −10.7%. Under off nominal conditions, the predicted gS and PN were much less accurate: rBIAS for gS was 68.2% and 48.6% for green-and red-  Figure 6 shows the box-plot distribution of the stomatal conductance (g S ) as well as of the net photosynthesis (P N ). Stomatal conductance and P N were significantly higher under nominal conditions compared to off-nominal, in both butterhead cultivars. Unfortunately, these data are only available for the 23rd DAT, when leaf gas-exchange analyses were performed. As illustrated in Table 5, g S relative BIAS (rBIAS) was equal to 39% and −0.1% in green-and red-leaf lettuce, respectively; P N relative BIAS (rBIAS) was instead equal to 34.1% and −10.7%. Under off nominal conditions, the predicted g S and P N were much less accurate: rBIAS for g S was 68.2% and 48.6% for green-and red-leaf lettuce respectively, while rBIAS for P N was 75.9% and 70.9%. These larger overestimations of g S and P N testifies the limit of some empirical formulations (e.g., Equation (10)) adopted by the MEC model to reproduce the impact of off-nominal conditions on the stomatal conductance. leaf lettuce respectively, while rBIAS for PN was 75.9% and 70.9%. These larger overestimations of gS and PN testifies the limit of some empirical formulations (e.g., Equation (10)) adopted by the MEC model to reproduce the impact of off-nominal conditions on the stomatal conductance.

Discussion and Conclusions
The proposed version of the MEC model, validated against experimental data on green-and red-leaf lettuce cultivars, grown under both nominal (low VPD) and off-nominal (high VPD) scenarios, proved to be reliable in predicting crop growth and transpiration rate. All previous versions of the MEC model considered CQY and CUE to be constant during plant growth [12,24,26,27]. This assumption for CQY and CUE constant behavior is not consistent with recent literature in which contrasting results on the topic are reported. For instance, in rice, Xu et al. (2019) [39] observed a decline in photosynthetic rate, dark respiration and quantum yield according to leaf aging. In the latter study, both parameters rapidly increased to a maximum around 15 days, to linearly decline as a response to plant aging. Similar findings were also reported in a study on Rhododendron maximum L., were the decline of CQY during leaf aging was exacerbated by the exposure to high light intensity [42].
The issue is even more complex for CUE, since the carbon use efficiency has been less characterized for horticultural plants and little information exists for lettuce under different environmental conditions [31]. Although many models still rely on a fixed value of CUE set around 0.5 [43,44], this topic has been questioned and more studies contrasting this theory have been reported. For example, Winzeler et al. (1976) [45] showed that CUE of barley increased during the early phases of the growing cycle, while a decrease was reported during the second half of the cultivation period. In forest species, Amthor (2000) [46] showed that CUE is reported to vary sharply with aging, within and among different species and environmental conditions, due to different respiratory needs for growth and maintenance [37,46]. Indeed, it should be noted that CUE represents how efficiently a plant incorporates carbon into biomass and can be defined as follows: where R D and R N are the daylight and night respiration, respectively. Thus, a constant CUE would indicate that plants always present a constant positive respiration rate and that changes in photosynthetic activity would determine limited variations in growth and respiration, both these scenarios being quite unlikely [36]. Many studies have found that changes in the respiration rate during the plant growth cycle are species/cultivar-specific and maintain similar trends as net photosynthesis [38,39,47]. Furthermore, the situation can be different for the same plant species under different environmental conditions. Plant growth under near-optimal conditions have been reported to have smaller changes in CUE than plant grown under limited conditions, because the relative growth rate, higher under optimal conditions, would minimize the effect of maintenance respiration coefficients on the carbon use efficiency [36].
Given these uncertainties in the determination of CUE and the observed variability for CQY, in the present study, we suggested a modified version of the MEC model structure, by aggregating the variables A, CUE, CQY into two variables (α and β) and by assuming for these two variables the same temporal patterns observed for CQY till maturity, under the assumption of physiological similarity. Thus, the number of model parameters to be calibrated was reduced to four. The calibration was performed against the DTR data, while the remaining TEB g s and P N experimental data were used for validation.
In the present study, we proved that in lettuce a temporal pattern exists in CQY changing during plant aging. Furthermore, our findings highlighted differences between green-and red-leaf lettuce plants. More specifically, a first stage of stationary CQY was observed just for the red 'Salanova' lettuce. This period could be attributed to the time required by this red-pigmented cultivar to adapt to the new environmental condition after transplanting. Indeed, green and red Salanova lettuce, although belonging to the same species, have been proven to have a different behavior even under the same environmental conditions. For instance, under both nominal and off-nominal conditions, red-leaf lettuce exhibited highest values of net photosynthesis, stomatal conductance, as well as a higher value of edible biomass [13].
Therefore, it was interesting to observe that all leaves of both green and red lettuce plants, notwithstanding the position inside the canopy level (top, medium, bottom), exhibited the same timing for CQY variation, except for the initial time of development (t D ). The time t D in red lettuce corresponds with the time of canopy closure (t A ), which was also observed to be equal to eight days for both green and red Salanova. Thus, it is feasible that the red cultivar requires an initial time interval for adapting after the transplant. However, red plants completely recover from their initial "delay" by reaching CQY MAX at time t Mi = 16 days, as it occurs for green-leaf lettuces, also showing the highest values of net photosynthesis, stomatal conductance and edible biomass, overall suggesting a better physiological performance [48].
In our study, the experimental data and particularly TEB, g s and P N assumed the highest values under nominal condition, as a result of the lower evaporative demand. Under a high VPD (dry air), the evaporative demand increases, and plants try to counteract dehydration by closing their stomata, thus decreasing photosynthetic rates and stomatal conductance [5,49,50]. Indeed, under high VPD conditions, transpiration rates were enhanced, and a plant might lose water from tissues with negative consequences on the whole plant-hydraulic system. Thus, these plants require more water to reach the field-capacity, compared to those grown under nominal-conditions. Generally, the cultivation of crops under high VPD results in yield drop-off [51] and often in quality loss [52][53][54][55], which are considered major problems for crop production.
This calibrated model was able to reproduce the observed transpiration and biomass growth, under both nominal and off-nominal conditions. This capability of prediction could represent an added value for the cultivation management in CEA because it may allow the prediction of any yield loss, consequent to sudden changes in the microenvironment. Furthermore, the reliability in the model prediction concerning the daily transpiration rates could allow the set-up of a precise irrigation schedule, according to changes in the environmental condition, similarly to what is done in other agricultural sectors, by developing decision support systems (DSS) based on the optimal combination of sensors and prediction models [56]. However, the model tended to overestimate stomatal conductance and photosynthesis under off-nominal conditions. A feasible technical explanation to these overestimations is that the empirical model was initially calibrated only for nominal conditions and that a "big leaf" approach is used to calibrate the model equations. A plant, especially when grown in a sub-optimal environment, triggers a cascade of biological processes leading to the development of leaves with different anatomical traits (especially those linked with conductance and hydraulics), thus influencing plant photosynthetic performance and the whole physiological behavior [57][58][59][60]. Therefore, in this specific case, overestimation of stomatal conductance and photosynthesis, while maintaining comparable values of transpiration under off-nominal conditions compared to nominal ones, can be explained by the lack of consideration of structural plasticity (e.g., mesophyll density and vein distribution) which can differentially establish the limits of different physiological processes [5]. In light of the above results, by applying this implemented version of the model to cultivation trials, it was possible to simulate variations in environmental parameters which can be due to sensor failure, power loss, and other problems related to environmental control. The present modified version of the MEC model can simulate crop growth, photosynthesis, and transpiration over a different range of environments and is therefore suitable to be implemented in decision support systems (DSS) for forecasting variations triggered by anomalies in the environmental control. However, the model still has a "big-leaf" approach and can therefore overestimate some processes happening at the crop morpho-physiological level. To increase the functionality of the model, a further step could be to modify the relation used to calculate g s and P N by considering morpho-physiological modifications that would affect plant gas exchanges under off-nominal conditions. Advanced