Modelling water requirements of greenhouse spinach for irrigation management purposes

The spread of greenhouse crops has significantly increased in recent decades in many countries around the world, concerning different latitudes and climatic zones (Pardossi et al. 2004; FAO 2013). In the coastal areas of the Mediterranean basin, mild temperatures and high levels of solar radiation even during winters have supported the recent expansion of protected crops (Pardossi et al. 2004; Bonachela et al. 2006). The development of this type of farming system has been promoted by the growing demand for high-quality fresh food throughout the year, the introduction of very efficient post-harvest storage systems, the development of transportation networks that quickly connect production and consumption locations, and the adoption of new manufacturing technologies and automated systems for monitoring and controlling most of the environmental and agronomic parameters (Fitz-Rodrìguez et al. 2009; Pardossi & Incrocci 2009). Although no reliable data on this issue exist, the Food and Agriculture Organization (FAO) estimates that approximately 405,000 hectares are covered by greenhouses worldwide, with different levels of sophistication and automation. This figure almost rises to 2 million hectares if it is expanded to include crops under tunnels or temporary protections. Greenhouse typologies span from simple plastic greenhouses with no real climate control systems to the most sophisticated systems with metal and glass structures equipped with modern technologies for the precise adjustment of the indoor environmental parameters to achieve the optimal microclimate and agronomic conditions for each crop and cultivar (Pardossi & Incrocci 2009). Low technology greenhouse systems and very simple plastic shelters represent the predominant typology in the Mediterranean basin (FAO 2013).

Greenhouse cropping systems are destined to play an increasingly important role in the future due to the possibility of offering sustainable intensive agricultural productions, even in marginal or urban areas or wherever there is a shortage of agricultural land (FAO 2013). Greenhouse cultivations, which mainly focus on high-income crops such as vegetables and freshly cut flowers, often permit a more efficient use of production factors compared to full field crops, and also offer a significant economic return, albeit in the context of a high initial investment (Pardossi & Incrocci 2009).

The sustainability of these agricultural systems, both in environmental and economic terms, is attained as long as each agronomic input (water, fertilizers, plant protection products and energy) is used as efficiently as possible (FAO 2013; Masseroni et al. 2016). Therefore, the development and the application of methods and tools able to provide accurate crop irrigation requirements and irrigation scheduling for each crop, are crucial to increase water use efficiency and to minimize crop water stress in greenhouse environments (FAO 2013). Furthermore, the need to plan and design greenhouse irrigation systems in advance calls for methods able to quantify the irrigation requirements of the different crops.

Crop water requirements correspond to the crop evapotranspiration computed in well-watered conditions and in the absence of all other kinds of stress (ETc). To calculate how much of the crop water requirement shall be provided by irrigation (i.e. crop irrigation requirement) and, conversely, how much water rather comes from rain, soil moisture and capillary rise and is lost due to actual evapotranspiration, percolation and runoff, a water balance equation needs to be developed for the root zone volume (Allen et al. 1998; Orgaz et al. 2005; Facchi et al. 2013). Inside the greenhouse, ETc is the cornerstone of the hydrological balance, due to the fact that it is usually the only outgoing water flux from the crop root zone, since runoff and percolation can be very often considered null or negligible (FAO 2013). Therefore, the amount of water that a hypothetical irrigation system must be able to replenish coincides almost exactly with the amount of water evapotranspirated in the period of time elapsed since the last irrigation, since the efficiency of greenhouse irrigation systems could be often considered close to one (Pardossi et al. 2004; FAO 2013). A precise estimate of ETc is therefore the starting point for a rational irrigation scheduling and a proper irrigation system sizing. The most commonly adopted method for the estimation of ETc is the FAO-56 ‘single crop coefficient’ approach proposed by Allen et al. (1998), where ETc is obtained by multiplying the reference evapotranspiration (ETo) and the crop coefficient (K_c), and a water balance equation is used to describe the variation of soil water content in the root zone with the aim of computing irrigation water requirements. In particular, an irrigation is due when the soil water content in the root zone drops below a critical threshold depending on the type of soil and crop, and the irrigation replenishes soil water content to field capacity.

The literature offers several models to support the schedule of irrigation, developed and applied by different authors for open-field crops; some of these approaches were later adapted to protected crops (Farahani et al. 2007; Fitz-Rodrìguez et al. 2009; Fernández et al. 2010). However, many of these approaches were calibrated and validated over short periods (sometimes less than a week), and may require many parameters for describing the evapotranspiration processes (e.g., stomatal response to inlet solar radiation); moreover, they are typically applied for polytunnel greenhouses with crops cultivated directly in the ground (Fernández et al. 2008).

In this work, a simplified version of the FAO-56 method designed for greenhouse crops (FAO-56-GH) is presented. The model was implemented using data coming from two different experiments carried out in 2015 on two spinach varieties (Verdi F1 and SV2157VB) grown in pots in a glass-greenhouse located in Milan (northern Italy), respectively January 28th to March 20th and from April 18th to May 24th. ETo was calculated by applying the FAO Penman–Monteith equation to the climatic data registered in the greenhouse during the two periods. Irrigation was applied daily to the pots to replenish evapotranspiration losses. Crop coefficient values used in the model were calibrated using the irrigation data collected during the first experiment. Data from the second experiment, carried out considering another spinach variety grown in a different period, were used for the model validation. Three statistical indices (R², root mean square error (RMSE) and Nash–Sutcliffe efficiency (NSE)) were computed to evaluate the model performance, while the analysis of variance (ANOVA) test was carried out to exclude that the nitrogen management could influence significantly crop development and irrigation requirements.

Two separate experiments (EX1 and EX2) were carried out in two different periods of the year. Two spinach cultivars were grown in a total of 18 pots (six for EX1 and 12 for EX2). In particular, for EX1 the spinach cultivar Verdi F1 (ISI Sementi), suitable for an early-spring greenhouse growing, was selected. During EX2, the SV2157VB cultivar (Seminis Seeds Vegetables) was grown due to its suitability to spring–summer greenhouse production. The two cultivars have a similar phenotypic development and canopy structure. The sowing dates were January 28th and April 18th, while the germination dates were February 4th and April 24th for the first and second cultivars, respectively.

The substrate of each pot was prepared by mixing 50% in volume of silica sand with neutral reaction (with granulometry between 0.4 and 0.8 mm) with 50% of soil coming from an agricultural field. The final substrate had the following characteristics: sand 52%, silt 39%, clay 9%, pH_water 7.4, pH_KCl 6.6, total nitrogen content 0.13%, total carbon content 1.4%, nitrogen nitric content 11.28 ppm, ammonia nitrogen content 3.63 ppm and bulk density 1.28 t m⁻³. The substrate volumetric water content at the field capacity (−10 kPa) and at the wilting point (−1,500 kPa), as obtained through conducting measurements by a Richards apparatus over disturbed soil samples collected in three different pots, was found to be respectively 0.14 m³ m⁻³ and 0.05 m³ m⁻³. The substrate was used to fill 18-litre volume pots. A layer of approximately 5 cm of expanded clay was placed on the bottom of the pots, and a thin filter fabric with fine mesh was placed between the substrate and the expanded clay layer to limit the drainage of the fine soil particles on the bottom of the pots, and to separate the expanded clay and substrate layers.

In both experiments, in order to evaluate whether crop evapotranspiration (and consequently irrigation water management) could be significantly influenced by a different nitrogen fertilization supply, pots were treated with four nitrogen levels: not fertilized (N0), limited concentration (N1), optimal concentration (N2) and redundant concentration (N3) (Table 1).

The optimal level of nitrogen fertilization (N2) was obtained for EX1 by a generalized equation for the C3 species, while the N1 and N3 levels were obtained by decreasing and increasing the optimal treatments by 50%, respectively. In EX1, the nitrogen application rate for N2 was 66 kg ha⁻¹. For EX2, the amount of nitrogen applied for N2 was increased to 132 kg ha⁻¹, with the aim to differentiate more significantly the four nitrogen treatments. Other macro- and micronutrients were added to the substrate during its preparation, to discount these nutrients as limiting factors for the crop growth (28 kg ha⁻¹ of P₂O₅, 118 kg ha⁻¹ of K₂O, 9 kg ha⁻¹ of Mg and 1 kg ha⁻¹ of Fe). The sowing was performed by manually placing the seeds in parallel furrows, covering the seeds with a thin layer of soil and applying a gentle pressure on the substrate. In EX1, the seeding density was approximately 100 plants m⁻² (2 rows per pot), while in EX2, the seeding density was increased to 200 plants m⁻² (3 rows per pot) to eventually accelerate the symptoms of water and nitrogen scarcity in pots under stress. The position of each pot belonging to the four different nitrogen treatments on the benches was randomly assigned. Further details on the experimental setup can be found in Corti et al. (2016).

Pots were maintained in well-watered conditions, replacing water lost by evapotranspiration daily. The water replacement was performed by weighing each pot daily during EX1 (using a weighing balance with an accuracy of ±1 g), and replenishing the soil water content to the field capacity value. The registration of the daily weight variation of each pot provided an affordable measure of the daily ETc, successively used for the calibration of the crop coefficient values considered in the model. Due to the higher number of pots in EX2, weighing was performed once every 2 days. During the day following the one in which the weighing was conducted, in order to still return a volume of water that would approximately replenish evapotranspiration losses without originating percolation, the minimum volume of water registered the previous day over all the pots was provided. In the successive day, the conduction of the weighing operations could ensure that volumes of water added to the pots were enough to replenish soil to the field capacity value. The irrigation procedure adopted during the two experiments allowed the spinach crop to avoid water stress conditions. In fact, soil water content of each pot never dropped below the critical soil water content threshold, corresponding to the depletion of the readily available water (RAW). RAW for spinach, according to Allen et al. (1998), was assumed equal to 20% of the total available water (TAW; i.e. soil water content between field capacity and wilting point).

where θ_t and θ_t−1 are the soil water content values at the time t and t−1, respectively; P_t is the precipitation; I_t is the irrigation amount; RO_t is the net surface runoff; G_t is the net subsurface water flow; DP_t is the deep percolation; CI_t is the rainfall (or irrigation) canopy interception; CR_t is the capillary rise; and ETc is the crop evapotranspiration in well-watered conditions. All the quantities in the equation can be expressed as volumetric or height units. In particular, I_t is the amount of water that shall be provided by irrigation to avoid ET_t to drop below ETc_t (i.e. to ensure the absence of crop water stress).

Since the unit volume for the application of the soil water balance equation in the case of greenhouse crops in pots is the pot volume, Equation (1) can be simplified. In particular, in the case of closed pots not irrigated by sprinklers, deep percolation, capillary rise, crop interception and surface/subsurface runoff can be neglected. Thus, I_t and ETc_t are, respectively, the only water input and water output in the soil water balance.

where e_a is the air vapour pressure (kPa) and T is the air temperature (K). The term ɛ_c was considered constant throughout the experiments, being calculated as the average of the surface emissivities of three main components of the greenhouse crop system, namely the spinach crop (ɛ = 0.80, Pieters & Deltour 1997), the plastic mulch on the bench (ɛ = 0.53, Zhu et al. 1998) and the soil in the pots (ɛ = 0.95, Valdés-Gómez et al. 2009).

To define the K_c curve along the crop cycle, seven parameters are necessary. Four parameters are related to the length of the crop phenological stages (Δt_ini is the initial stage, Δt_dev is the development stage, Δt_mid is the middle-season stage and Δt_end is the late-season stage), while the other three parameters correspond to the values that K_c assumes in the same stages (K_{c_ini}, K_{c_mid}, K_{c_end}). Values of these parameters are reported in Allen et al. (1998) for several crops, but site-specific experimental values are highly preferred when they are available. This is particularly true in the case of greenhouse crops, where these values can be markedly altered by the specific climatic conditions of greenhouses (Orgaz et al. 2005; Farahani et al. 2007).

Microclimatic parameters measured inside the greenhouse were used as input variables for the estimation of daily reference evapotranspiration (ETo). Horizontal wind velocity inside the structure proved to be the most uncertain variable, since it was not measured by the instrumentation installed in the greenhouse. To set a value for this variable, a measuring campaign was conducted on January 28th at different locations inside the structure by positioning a Young 81000 sonic anemometer at 2 m above the ground, and considering only the two horizontal components of the wind speed (R.M. Young Company, PA, USA). The measured average value of horizontal wind velocity was equal to 0.4 m s⁻¹, which was considered constant over the two experiments. Albedo (α) was calculated as the ratio between the inlet and outlet shortwave radiation in the PAR range. The two radiation terms were measured at a distance of about 30 cm from the top of the canopy by a PAR sensor (Apogee, USA) in five different points over the bench on March 8th, and the resulting average albedo value was 0.2.

For each phenological stage of the EX1 period, the three K_c values (K_c-ini, K_c-mid, K_c-end) were obtained averaging the daily ratio between observed evapotranspiration (ETc_obs) and the corresponding reference evapotranspiration (ETo). In particular, ETc_obs corresponded to the water loss measured by the daily weighing. K_c values obtained in the calibration phase EX1 (K_c_cal) were then applied in EX2 (validation phase) considering the specific phenological stage lengths determined for the second experimental period.

When the estimation is perfect (i.e. ETc_mod = ETc_obs) RMSE is equal to zero; the higher the RMSE, the higher the deviation between the simulated and the measured values. With respect to NSE, simulation is perfect if NRMSE is zero; predictions are worse than using the mean of observed values if NRMSE is greater than one.

Finally, the one-way ANOVA test was applied to check: (a) whether the difference in fraction cover values reached during the middle-season stage for the experiments EX1 and EX2 could be considered statistically not significant over the four nitrogen treatments; (b) whether the difference in fraction cover values reached during the middle-season stage for the experiments EX1 and EX2 could be considered not significant despite of the seeding density; (c) whether the difference in cumulative irrigation amount values at the end of the crop cycle measured for EX1 and EX2 in the case of the four nitrogen treatments could be considered statistically not significant. For objective (a), the ANOVA test was applied separately for EX1 and EX2 considering the fraction covers obtained through the image analysis for the four nitrogen treatments respectively on March 19th and May 23rd for the two experiments. In the case of objective (b), the ANOVA test was applied considering the two fraction covers obtained for the experiments EX1 and EX2 by averaging, for each experiment, the cover fractions for the four nitrogen treatments retrieved on March 19th and May 23rd. Finally, for objective (c), the ANOVA test was applied, for each experiment (EX1 and EX2), over the four cumulative irrigation values obtained at the end of each experiment from the corresponding four nitrogen treatments.

The average daily relative air humidity (not shown in the graph) in EX2 ranged from 60% to 80%, while in EX1 was found to be between 45% and 60%. The higher relative humidity in EX2 may be partially justified by the transpiration flux coming from the larger number of plants in the greenhouse (belonging to different experiments), and by the use of the cooling system.

The internal global radiation (Figure 3(c) and 3(d)) showed a typical sinusoidal trend with highest and lowest peaks around midday and midnight. The radiation inside the greenhouse reached a maximum of approximately 200 W m⁻² in both campaigns, while the external radiation reached maximum values of 600 W m⁻² in EX1 and 900–1,000 W m⁻² in EX2. The difference between the internal and external radiation values was likely due to the cover glass decreasing the incoming radiation. The effect of the lamps is clearly visible in Figure 3(c) and 3(d) during early morning and evening hours, i.e. when the external radiation was around zero while the internal one was 50 W m⁻².

The daily average ETo for the entire EX1 and EX2 periods is shown in Figure 3(e) and 3(f). ETo pattern is obviously determined by the pattern of the other climatic variables, according to Equation (1). For instance, March 5th, which showed the highest ETo value over the entire EX1 period (3.5 mm day⁻¹), was also the day characterized by the maximum values of air temperature and global radiation (35 °C and approximately 650 W m⁻²), as well as by the lowest value of relative humidity (below 20%). Overall, ETo values ranged between 1.5 and 3.5 mm day⁻¹ in EX1 and between 2 and 3.5 mm day⁻¹ in EX2.

For both EX1 and EX2, the result of the ANOVA test confirmed that a not significant difference (p-value >0.05) was found between the fraction covers reached in the mid-season stage in the case of the four nitrogen treatments. Thus, for each experimental campaign, a single fraction cover curve was calculated averaging all the RGB images acquired at each single date; the two curves were then used to derive the growth stage lengths reported above.

A further ANOVA test was conducted considering the two average fraction covers (one for EX1 and one for EX2) thereby obtained. The result showed that no significant differences can be evidentiated between the two values (p-value >0.05), so that the sowing density proved to be a factor not affecting cover fraction (and thus K_c values) in the mid-season stage, at least in this study.

Figure 5 highlights the high values of K_{c_ini} and K_{c_mid} obtained experimentally from the calibration procedure. In particular, the value of K_{c_ini} = 0.9 is high with respect to that tabulated in FAO-56 (K_c__ini_FAO-56 = 0.7), also after applying the corrections for the wetting frequency recommended by Allen et al. (1998) (K_c__ini_FAO-56_adj = 0.8). In addition, also the K_{c_mid}= 1.5 is significantly higher than the values obtained applying what was reported in the FAO-56 handbook, both before (K_c__mid_FAO-56 = 1.0) and after the suggested corrections (K_c__mid_FAO-56_adj = 1.05). At first, the K_{c_mid} value seems to contradict the FAO-56 prescriptions for small herbaceous crops, that suggest that this value should not exceed 1.3. However, the literature illustrates many cases in which strong differences are registered between calibrated and FAO-56 K_c values, especially in the case of greenhouse environments (Cereković et al. 2010,). In a greenhouse experiment on cucumbers, Blanco & Folegatti (2003) determined a K_{c_mid} value of approximately 1.6, in contrast to the value of 1.0 reported in the FAO-56 handbook. Moreover, Orgaz et al. (2005) showed examples of K_{c_mid} values for beans and peppers respectively of 1.4 and 1.3, which are higher than those reported in the FAO-56 (K_c__mid_FAO-56 = 1.15). In the study K_{c_end} was assumed to be equal to that tabulated in FAO-56 (K_{c_end_}adj = 0.90), as a consequence of the absence of irrigation interventions in the senescence phase before the harvesting.

Part of the discrepancy between experimental and FAO-56 K_c values may be surely explained by the particular agro-climatic conditions of greenhouses in contrast to open fields. Another factor to be taken into account in this study is the fact that spinach was cultivated in pots. Soil inside the black plastic pots was probably more subjected to heating due to the higher solar radiation absorption with respect to a soil extending with continuity under vegetation. This circumstance could have resulted in an increased soil evaporation, and thus in higher K_c values.

During the model validation phase (EX2), K_c_cal values obtained from EX1 were applied considering the phenological stage lengths measured for the second spinach variety. In EX2, from emergence to April 29th K_{c_ini} was 0.9, from May 11th to May 19th K_{c_mid} was 1.5, and finally, on May 24th K_{c_end} was 0.95.

Figure 6(b) compares single and cumulative irrigation amounts provided to the pots during EX1. Cumulative irrigations, both from the model and from observational data were approximately 130 mm, denoting a good ability of the model to predict the total amount of water needed by the crop. When looking at the daily irrigation patterns, it can be noted that the model suggested irrigation, on average, every 2 days instead of daily, as performed during the experimental activity to avoid any kind of plant water stress. This is explained by the fact that the model scheduled an irrigation event only when RAW (20% of TAW for spinach) was completely depleted. The experimental irrigation heights were in the range of 2–7 mm, while the model suggested an average irrigation height of 5–7 mm with peaks of 8–9 mm.

The ANOVA test performed considering the cumulative irrigation amount over the entire crop cycle showed no significant differences between irrigation values supplied to pots with different nitrogen treatments (p > 0.05), further confirming what had already been found for fraction cover; namely, that from the point of view of crop development and irrigation consumption, no statistically significant differences between the nitrogen treatments can be noted in EX1.

The lower model performance with respect to EX1 could be partially explained by the fact that a different weighing procedure was adopted in EX2 compared to EX1. In EX2, for practical reasons, the weighing was performed every 2 days. On the weighing day, soil water content of each pot was reported to field capacity, while in the day between 2 weighing days, the irrigation amount provided was equal to the minimum volume registered the day before. The value obtained averaging the irrigation amounts for the 2 days was than attributed to each of the two. This explains the pairs of equal evapotranspiration values for 2 consecutive days in Figure 6(a), and surely introduced an uncertainty in the daily values of ETc_obs for EX2. Additionally, the modelled evapotranspiration is notably different from the observed evapotranspiration on May 12th and May 13th (circled in Figure 7(a)). On these days, the observed evapotranspiration is approximately 7 mm, which is much higher than that expected by applying the model. This difference is very probably due to an error made in recording data of water amount provided to the pots (values are very high with respect to other measured ETc, also considering the climatic conditions in the greenhouse which were not dissimilar to those of the period). For that reasons, this pair of ETc_obs value was discarded in the calculation of performance indices for EX2. At the harvesting time, the cumulative irrigation was 115 mm with respect to 108 mm simulated by the model. Irrigation frequencies and heights were similar to those simulated by the model in EX1.

As for EX1, the ANOVA test showed a not significant difference between irrigation provided to pots with different nitrogen contents (p > 0.05), further confirming the possibility to develop and apply a single model to support the irrigation scheduling independently from the crop nitrogen management, at least in the case study.

In this work, we demonstrated the reliability of a simple model based on the FAO-56 method (FAO-56-GH) in the prediction of daily evapotranspiration fluxes, daily irrigation scheduling and total irrigation requirements of a greenhouse spinach crop cultivated in pots. The K_c pattern used in the model was calibrated and validated using two independent datasets, obtained by cultivating, in two separate experiments, two varieties of spinach (Verdi F1 and SV2157VB) in a greenhouse in different periods of the year. Despite the differences in crop variety, seeding density and growth period, the performance of the model was high even in the validation period (R² = 0.80, RMSE = 0.41 mm day⁻¹, NSE = 0.78 for the validation period; R² = 0.84, RMSE = 0.21 mm day⁻¹, NSE = 0.83, for the calibration period). Lengths of crop development stages were obtained in the two experiments by post-elaborating RGB images acquired at different times for the crop. Values of K_c calibrated in the first experiment and then validated in the second were the following: K_{c_ini} = 0.9; K_{c_mid} = 1.5; K_{c_end} = 0.95. The high K_c-mid value is confirmed by other studies carried out in greenhouses and may be explained by the particular agro-climatic conditions of these controlled environments, as well as, in the case of this specific study, by the plastic pots in which the crop grew which may likely have increased the soil evaporation. These considerations reinforce the belief that, even in greenhouse environments, effort must be placed in obtaining experimental K_c values, since they could be very different from the FAO tabulated ones. During the two experiments, statistically significant differences in crop development and irrigation requirements for pots subjected to different nitrogen treatments were excluded by applying the ANOVA tests. Consequently, crop data used to implement the model and model results in terms of irrigation scheduling were found to be essentially independent from the nitrogen management adopted.

The model designed and applied in this study may be useful for a rational irrigation scheduling in greenhouse environments, representing a valuable tool for decreasing water waste and increasing water use efficiency. This can become particularly important in productive conditions where a specific crop is cultivated in the same environment in a nearly continuous process, and thus the effort placed for the site-specific determination of K_c could be more exploited. Moreover, the prediction of crop water requirements could be crucial in supporting greenhouse irrigation systems designing to supply water to specific crops.

The authors wish to thank the DiSAA (Department of Agricultural and Environmental Sciences of the University of Milan) for founding the STeP project (Sensing Technologies for Precision Agriculture), in which this research was conducted. Moreover, special thanks go to Dr Pietro Marino, Dr Martina Corti, Antonia Moreno and all the other colleagues who gave invaluable help in projecting and managing the greenhouse experiment.