Estimation of Actual Crop Coefficients Using Remotely Sensed Vegetation Indices and Soil Water Balance Modelled Data

Abstract

A new procedure is proposed for estimating actual basal crop coefficients from vegetation indices (K_{cb VI}) considering a density coefficient (K_d) and a crop coefficient for bare soil. K_d is computed using the fraction of ground cover by vegetation (f_{c VI}), which is also estimated from vegetation indices derived from remote sensing. A combined approach for estimating actual crop coefficients from vegetation indices (K_{c VI}) is also proposed by integrating the K_{cb VI} with the soil evaporation coefficient (K_e) derived from the soil water balance model SIMDualKc. Results for maize, barley and an olive orchard have shown that the approaches for estimating both f_{c VI} and K_{cb VI} compared well with results obtained using the SIMDualKc model after calibration with ground observation data. For the crops studied, the correlation coefficients relative to comparing the actual K_{cb VI} and K_{c VI} with actual K_cb and K_c obtained with SIMDualKc were larger than 0.73 and 0.71, respectively. The corresponding regression coefficients were close to 1.0. The methodology herein presented and discussed allowed for obtaining information for the whole crop season, including periods when vegetation cover is incomplete, as the initial and development stages. Results show that the proposed methods are adequate for supporting irrigation management.

1. Introduction

The accurate estimation of crop water requirements plays an important role in the improvement of crops water use and irrigation performance. This issue is particularly relevant considering the need for intensification of irrigated agriculture and increased water scarcity in several regions of the world.

A common approach for the estimation of crop water requirements is the K_c-ET_o approach adopted by FAO56 [1] where a reference crop evapotranspiration (ET_o) is multiplied by a crop coefficient (K_c) for estimating crop evapotranspiration, ET_c. The ET_o represents the climatic demand of the atmosphere while the K_c represents the differences distinguishing the reference crop and the considered crop in terms of ground cover, canopy properties and aerodynamic resistance, thus in terms of crop ET. Alternatively, crop evapotranspiration may be computed directly from ground observations using a combination equation such as the Penman-Monteith equation [2,3]. However, this approach is more demanding than the K_c-ET_o approach and is not used for operational purposes but limited to research.

Another approach for estimating crop ET is based on surface energy balance models that use remote sensing thermal infrared data (e.g., [4,5,6]). These models estimate crop ET by subtracting the soil heat flux (G) and sensible heat flux (H) from the net radiation (R_n) at the surface. Such models allow to directly integrating the effects related with soil water deficit or water vapor pressure in the crop ET estimation [7]. Crop coefficients can be further derived considering the ET_o [8,9,10]. However, these models present greater complexity and larger amount of input data than the K_c-ET_o approach [7]. The K_c-ET_o approach is commonly accepted for operational and research objectives [11]. A single or dual K_c approach may be used [1]. In the single approach both crop transpiration and soil evaporation are timely averaged into a single coefficient (K_c), whereas in the dual approach a daily basal crop coefficient (K_cb), representing primarily the plant transpiration, and a daily soil evaporation coefficient (K_e) are considered separately, i.e., K_c = K_cb + K_e.

For transferability purposes, K_c, K_cb and ET_c in FAO56 represent ET rates under optimal, well-watered conditions [1]. However, in the field and in common practice, crop conditions are often not optimal due to insufficient or non-uniform irrigation, crop density, soil salinity and/or agronomic management. Potential ET_c must then be replaced by the actual ET_c (ET_{c act}), and the resulting K_c is renamed K_{c act} and K_cb is also renamed K_{cb act} [11]. The term actual is adopted in this study. As commented by Pereira et al. [11], this concept better supports estimation and transferability avoiding the need to define multiple K_c values for the same crop as has occurred in the past for several crops, e.g., vines and orchards. Moreover, adopting K_{c act} and K_{cb act} rather than using the potential crop coefficients allows to perform irrigation management closer to the actual conditions, which is a main motivation of the present study.

Standard, potential K_c and K_cb values are defined and tabulated for a wide range of agricultural crops [1,12], but appropriate corrections may be required adopting a stress coefficient (K_s) to obtain the actual K_c (K_{c act} = K_s K_c or K_{cb act} = K_s K_cb). In addition, the transferability of K_c and K_cb values requires appropriate adjustment to local climate [1]. The adjustment of tabulated K_c and K_cb values for local conditions is also necessary when differences occur in planting density and geometry, vegetation height, and canopy architecture particularly in case of tree and vine crops [1,12].

In the last two decades, remote sensing (RS) data has increasingly been used for monitoring and mapping the spatial and temporal variation of ET and thus for computing crop water requirements (e.g., [13,14,15]). One of the most widely used approaches considers the relationship between vegetation indices (VI) derived from RS reflectance data and (actual) crop coefficients, either K_c or K_cb. The basis of this approach relies on the close correlation of several VI and various biophysical characteristics of the plants, e.g., leaf area index (LAI), ground cover fraction (f_c), biomass, and physiological processes depending on light absorption by the canopy, including ET (e.g., [16,17,18,19]). With the VI approach, because these indices reflect the actual vegetation cover conditions, the estimated K_c or K_cb represent actual rather than potential K_c or K_cb. As proposed by Pereira et al. [11], these coefficients should then be referred to as K_{c act} and K_{cb act} despite many authors have not adopted this conceptual difference and that most VI-based methods are unable to accurately observe reductions in K_c and K_cb caused by acute water or salinity stress.

One advantage of using VI-based actual crop coefficients is the ability to account for variations in plant growth due to abnormal weather conditions, e.g., the impact of a frost occurrence [20]. Using VI-based crop coefficients also allows obtaining the spatial variation of actual K_c (or K_cb) within fields. In addition they provide for a field-to-field description of the variation of actual K_c (or K_cb) due to variations in planting dates, plant spacing and cultivars [11]. A review of the advantages and disadvantages of VI-based crop coefficients was presented by Allen et al. [21].

The Normalized Difference Vegetation Index (NDVI; [22]) and the Soil Adjusted Vegetation Index (SAVI; [23]) are the most commonly VI used to estimate actual K_c and K_cb. The formulation of both VI combines the reflected light in the red and near infrared (NIR) bands, thus providing an indirect measure of the absorption of red light by chlorophylls (a and b) and reflectance of NIR by the mesophyll structure in leaves [17,24]. Choudhury et al. [25] refer that using SAVI instead of NDVI allows extending the range over which the VI respond to the increase in vegetation amount/density beyond a LAI around 3, which is related with NDVI saturation problems for high LAI values. In addition, the NDVI is considered more sensitive than SAVI to soil background reflectance changes due to the moisture of soil surface [7,13,26].

Several relationships between K_cb (or K_c) and VI are reported in the literature for different vegetation types, with the most recent studies focusing the K_cb estimation because plant transpiration is more directly related with the VI. Some examples of this approach consider a linear relationship between K_cb and NDVI [13,20,25,27,28], and between K_cb and SAVI [13,25,26] but an exponential relationship between NDVI and K_c was also considered [29]. A few studies report a K_cb estimation based on f_c or LAI obtained through VI, including for field and tree crops [7,30,31]. However, in most studies, the relationship between K_cb and VI was established for non-stressed conditions and for conditions of dry soil surface, which often do not represent the actual conditions of crop management. VI derived crop coefficients often do not account for the reduction in K_cb due to water or salinity stress [7,11]. Then the above referred stress coefficient K_s should be considered to obtain an actual K_c or K_cb, i.e., K_{c act} and K_{cb act} [1]. Most of those studies were established for a single crop and its application for other crops, particularly for discontinuous tree crops, is still limited. Differently, the approach proposed by Mateos et al. [7], which approximates K_cb using f_c and VI, refers to several crops and shows good results. However, this application does not always allow estimating K_cb for periods when the soil is not yet fully covered by the crop.

Although several equations were proposed for approaching K_cb via VI, a consensual equation that can be consistently used for different vegetation types and conditions is still lacking. In such context, a density coefficient (K_d) [12] may be used in the formulation of K_cb through VI to help incorporating the impact of vegetation density and height. Such approach would allow K_cb to be computed for a range of row crops, orchards and vines, thus including conditions when soil is not fully covered by the crop. For further estimating K_c, a water balance modeling approach adopting the dual crop coefficient approach [1] is required to estimate soil evaporation and thus to estimate the evaporation coefficient K_e that can be combined with the VI based K_cb, i.e., K_c = K_cb + K_e. The SIMDualKc model [32], which performs a daily water balance of the surface soil layer, provides for a daily K_e. Thus, by applying the dual coefficient approach and performing the soil water balance, the SIMDualKc model provides information on the stress coefficient K_s, on potential and actual K_cb, on K_e and ultimately on potential and actual K_c.

Considering the advances and limitations discussed above, particularly the need to represent actual rather than potential crop and evapotranspiration conditions for irrigation management purposes, the main objectives of this study are: (i) developing and testing a new equation for estimating actual K_cb based on VI and considering the stress coefficient K_s and the density coefficient K_d computed with the fraction of ground cover f_c estimated with VI; (ii) developing and testing a combined approach for estimation of the actual K_c using the VI-based K_cb and the evaporation coefficient K_e obtained with a daily soil water balance model; and (iii) testing the adequacy of the approach to different crop types, namely maize, barley, and an olive orchard.

2. Material and Methods

2.1. Study Areas

Three crops were considered to test the approach of actual K_c estimate based on reflectance VI in combination with the SIMDualKc soil water balance model: a super high density olive orchard, maize, and barley.

The super high density olive orchard test site is located in Viana do Alentejo, in South of Portugal (38°24′46′′ N, 7°43′38′′ W, 143 m a.s.l.; Figure 1). This olive orchard occupies a total area of 78 ha in an undulating terrain. The olive orchard was planted in 2006 in a hedgerow system with 1.35 m × 3.75 m spacing (1975 trees ha⁻¹), and an orientation South-North. The olive trees are of cultivar Arbequina. The fraction of ground covered by the vegetation (f_c) was approximately 0.35 and tree height was around 3.5 m. The procedure to obtain f_c was based upon the measurements of the projection of crown diameters in row direction and perpendicular to it in 51 trees. Soils are sandy loam, with soil water content averaging 0.24 cm³·cm⁻³ at field capacity and 0.12 cm³·cm⁻³ at the wilting point. Ground data from two consecutive years (2011 and 2012) were collected and used to validate information of several agronomic and biophysical parameters as detailed by Paço et al. [8] and Pôças et al. [33]. In 2012, between February 20th and 25th, an extremely heavy frost affected the orchard causing a strong leaf fall. A severe pruning of the trees was therefore applied following the frost occurrence. Some areas of the olive orchard were more affected by the frost, and consequently the pruning was not uniform along the orchard, which increased the variability of the vegetation conditions within the orchard in 2012.

Figure 1. Study areas location in Viana do Alentejo and Alpiarça, Portugal.

Figure 1. Study areas location in Viana do Alentejo and Alpiarça, Portugal.

Three farmer’s fields of Quinta da Lagoalva de Cima located in Alpiarça, Ribatejo region, were surveyed for the maize crop [34]: (i) field 1 (39°17′33′′ N, 8°34′19′′ W, 14 m a.s.l.); (ii) field 2 (39°17′45′′ N, 8°34′01′′ W, 14 m a.s.l.); and (iii) field 3 (39°18′04′′ N, 8°32′22′′ W, 16 m a.s.l.), (Figure 1). Fields have a cropped area of approximately 30 ha. In 2010, the maize was cropped in both fields 1 and 2, while in 2011 only field 1 was sowed; in 2012, fields 2 and 3 were cropped with maize. The PR33Y74 maize hybrid (FAO 600) was cultivated with a density of approximately 82,000 plants ha⁻¹. The length of each crop growth stages for all maize fields and seasons is described in Table 1. The crop growth stages were defined according to Allen et al. [1]: (i) initial stage, from the sowing/planting until 10% of ground cover (Ini); (ii) development stage, from 10% of ground cover until maximum ground cover by vegetation (Dev); (iii) mid-season stage, from full cover until maturity (Mid), and (iv) late-season stage, from the beginning of senescence and leaves yellowing until harvesting (Late). The observations of f_c were visually performed and adjusted with the help of photographs of the ground shadow by the crop near solar noon [34].

Barley was sown in field 1; two seasons were studied: (i) 2012, a dry year, and (ii) 2012/2013, a wet year [35]. The barley field was cropped with 200 kg·ha⁻¹ of malting barley (Cv. Publican) seeds using an inter-row spacing of 0.15 m. The plants density was measured before tillering, when the crop attained tree leaves, averaging 342 and 319 plants per m² respectively in 2012 and 2013. The lengths of the crop growth stages are given in Table 1. The barley f_c observations were performed similarly to the maize crop; furthermore, LAI measurements were also used to estimate f_c and therefore to confirm related observations as referred by Pereira et al. [35].

Table 1. Crop growth stages of maize and barley crops in each field and campaign.

Crop Growth Stages	Maize Study Areas	Barley Study Areas
	Field 1—Year 2010	2012 season
Initial (Ini)	25/05–25/06	16/01–06/02
Development (Dev)	26/06–22/07	07/02–03/04
Mid-season (Mid)	23/07–04/09	04/04–20/05
Late-season (Late)	05/09–13/10	21/05–26/06
	Field 2—Year 2010	2012/2013 season
Initial (Ini)	25/05–25/06	06/12–12/01
Development (Dev)	26/06–17/07	13/01–10/03
Mid-season (Mid)	18/07–02/09	11/03–05/05
Late-season (Late)	03/09–13/10	06/05–06/06
	Field 1—Year 2011
Initial (Ini)	20/04–17/05
Development (Dev)	18/05–28/06
Mid-season (Mid)	29/06–17/08
Late-season (Late)	18/08–20/09
	Field 2—Year 2012
Initial (Ini)	16/04–08/05
Development (Dev)	09/05 24/06
Mid-season (Mid)	25/06–20/08
Late-season (Late)	21/08–20/09
	Field 3—Year 2012
Initial (Ini)	30/05–12/06
Development (Dev)	13/06–15/07
Mid-season (Mid)	16/07–13/09
Late-season (Late)	14/09–12/10

Table 1. Crop growth stages of maize and barley crops in each field and campaign.

Crop Growth Stages	Maize Study Areas	Barley Study Areas
	Field 1—Year 2010	2012 season
Initial (Ini)	25/05–25/06	16/01–06/02
Development (Dev)	26/06–22/07	07/02–03/04
Mid-season (Mid)	23/07–04/09	04/04–20/05
Late-season (Late)	05/09–13/10	21/05–26/06
	Field 2—Year 2010	2012/2013 season
Initial (Ini)	25/05–25/06	06/12–12/01
Development (Dev)	26/06–17/07	13/01–10/03
Mid-season (Mid)	18/07–02/09	11/03–05/05
Late-season (Late)	03/09–13/10	06/05–06/06
	Field 1—Year 2011
Initial (Ini)	20/04–17/05
Development (Dev)	18/05–28/06
Mid-season (Mid)	29/06–17/08
Late-season (Late)	18/08–20/09
	Field 2—Year 2012
Initial (Ini)	16/04–08/05
Development (Dev)	09/05 24/06
Mid-season (Mid)	25/06–20/08
Late-season (Late)	21/08–20/09
	Field 3—Year 2012
Initial (Ini)	30/05–12/06
Development (Dev)	13/06–15/07
Mid-season (Mid)	16/07–13/09
Late-season (Late)	14/09–12/10

Soils in fields 1 and 2 are loamy sand soils, with total available water TAW = 171 and 149 mm·m⁻¹ respectively; field 3 is a silty-loam soil, with TAW = 209 mm·m⁻¹. Ground data collected in maize and barley study areas included LAI and fc, whose measurements were used to calibrate the SIMDualKc model as detailed by Paredes et al. [34] and Pereira et al. [35].

2.2. Satellite Imagery

Vegetation indices used for estimating the basal crop coefficients were obtained from Landsat 5 TM and Landsat 7 ETM+ satellite images. Landsat images of Path 203/Row 033 were considered for the olive orchard study area, while for the other study areas the Path 204/Row 033 was used. Table 2 presents the summary of the cloud-free satellite image dates considered for each study area and crop. The satellite image of 03 January 2013 was eliminated from the set of considered images due to the unfavorable light conditions during this period of the year, which affected the VI values.

Table 2. Satellite image dates considered for each study area and crop.

Olive Study Area (203/033)	Maize Study Areas (204/033)			Barley Study Areas (204/033)
	Field 1	Field 2	Field 3	2012	2012–2013
31/01/2011 (1)	06/07/2010 (2)	06/07/2010 (2),*	11/07/2012 (2)	02/02/2012 (1)	03/01/2013(1)
20/03/2011 (2)	22/07/2010 (3),*	22/07/2010 (3)	13/09/2012 (3)	18/02/2012 (2)	25/04/2013(3)
05/04/2011 (2)	30/07/2010 (3),*	30/07/2010 (3),*	29/09/2012 (4)	05/03/2012 (2)	11/05/2013(4)
23/05/2011 (3)	15/06/2011 (2)	11/07/2012 (3)		21/03/2012 (2)
24/06/2011 (3)	25/07/2011 (3)	13/09/2012 (4)		24/05/2012 (4),*
26/07/2011 (3)	18/08/2011 (3)
27/08/2011 (3),*	19/09/2011 (4)
12/09/2011 (3)
06/10/2011 (4),*
11/02/2012 (1)
15/04/2012 (2),*
20/07/2012 (3),*
21/08/2012 (3),*
06/09/2012 (3)
08/10/2012 (4)

⁽¹⁾ Ini; ⁽²⁾ Dev; ⁽³⁾ Mid; ⁽⁴⁾ Late (crop development stages). * Dates when water stress occurrence was detected using the soil water balance model.

Table 2. Satellite image dates considered for each study area and crop.

Olive Study Area (203/033)	Maize Study Areas (204/033)			Barley Study Areas (204/033)
	Field 1	Field 2	Field 3	2012	2012–2013
31/01/2011 (1)	06/07/2010 (2)	06/07/2010 (2),*	11/07/2012 (2)	02/02/2012 (1)	03/01/2013(1)
20/03/2011 (2)	22/07/2010 (3),*	22/07/2010 (3)	13/09/2012 (3)	18/02/2012 (2)	25/04/2013(3)
05/04/2011 (2)	30/07/2010 (3),*	30/07/2010 (3),*	29/09/2012 (4)	05/03/2012 (2)	11/05/2013(4)
23/05/2011 (3)	15/06/2011 (2)	11/07/2012 (3)		21/03/2012 (2)
24/06/2011 (3)	25/07/2011 (3)	13/09/2012 (4)		24/05/2012 (4),*
26/07/2011 (3)	18/08/2011 (3)
27/08/2011 (3),*	19/09/2011 (4)
12/09/2011 (3)
06/10/2011 (4),*
11/02/2012 (1)
15/04/2012 (2),*
20/07/2012 (3),*
21/08/2012 (3),*
06/09/2012 (3)
08/10/2012 (4)

⁽¹⁾ Ini; ⁽²⁾ Dev; ⁽³⁾ Mid; ⁽⁴⁾ Late (crop development stages). * Dates when water stress occurrence was detected using the soil water balance model.

The images were geometrically and atmospherically corrected using the procedure described by Allen et al. [15] and Tasumi et al. [36]. The procedure included image calibration using the coefficients proposed in the literature for Landsat 5 TM and Landsat 7 ETM+ [37,38]. The VI were calculated on a pixel-by-pixel basis and averaged for the entire study plot in each study area. Pixels from the edges of the fields and from roads within the olive orchard were not considered. The VI computed for the study were the Normalized Difference Vegetation Index (NDVI; [22]) and the Soil Adjusted Vegetation Index (SAVI; [23]).

2.3. SIMDualKc

The SIMDualKc is a soil water balance model that applies the dual crop coefficient approach, thus separately computing the daily soil evaporation and crop transpiration [32,39]. This model has been successfully applied to estimate ET and assess K_c of a wide range of field crops [34,35,39], including tree crops [8,40]. Data from SIMDualKc were used as benchmark, i.e., a standard against which comparisons can be done, mainly for the comparison of K_c and K_cb results. SIMDualKc was previously calibrated using ground data for the three crops, as described later in this section.

In SIMDualKc, K_cb is computed with the equation below [12,32] where impacts of plant density and/or leaf area are taken into consideration by a density coefficient:

$$K_{cb}= K_{c min}+ K_d\left(K_{cbfull}- K_{c min}\right)$$

(1)

where K_d is the crop density coefficient, K_{cb full} is the estimated basal K_cb for peak plant growth conditions having nearly full ground cover (or LAI > 3), and K_{c min} is the minimum K_c for bare soil (in the absence of vegetation). The K_{c min} value is about 0.15 under typical agricultural conditions and ranges 0.0–0.15 for native vegetation depending on rainfall frequency. K_cb is adjusted by the model for local climatic conditions where the average minimum relative humidity differs from 45% and/or the average wind speed is different from 2 m·s⁻¹ [1,12,32]. K_d is computed with the equation proposed by Allen and Pereira [12]:

$$K_d=min\left(1,M_L{f}_{c eff},f_{c eff}^{ \left(\frac{1}{1+h}\right)}\right)$$

(2)

where f_{c eff} is the effective fraction of ground covered or shaded by vegetation near solar noon [ ], M_L is a multiplier [ ] on f_{c eff} describing the effect of canopy density on shading and on maximum relative ET per fraction of ground shaded (to simulate the physical limits imposed on water flux through the plant root, stem and leaf systems), and h is the mean height of the vegetation [m]. Therefore, the use of K_d allows incorporating the impact of vegetation density and height in the estimation of K_cb. K_d has also been used in other studies to adjust K_cb to actual density conditions for full cover crops [41] and in partial cover crops using both ground [42] or remote sensing data [43].

When soil water deficit occurs, the stress coefficient K_s [ ] is computed by the model using a daily soil water balance for the entire root zone. K_s is expressed as a linear function of the root zone depletion D_r [1,32]:

$$K_s= \frac{TAW- D_r}{TAW-RAW}= \frac{TAW- D_r}{\left(1-p\right)TAW}, \text{for} D_r > \text{RAW}$$

(3)

$$K_s=1, \text{for} D_r \le \text{RAW}$$

(4)

where TAW and RAW are, respectively, the total and readily available soil water [mm], D_r is the root zone depletion [mm], and p is the depletion fraction for no stress [ ]. K_cb is multiplied by K_s to account for the effects of soil water stress and thus obtaining K_{cb act}.

The calculation of the soil evaporation coefficient K_e [ ] follows the procedures proposed by [1,44] which include the daily water balance of the evaporation soil surface layer that allows to compute the evaporation reduction coefficient K_r [ ] from the cumulative depth of water depleted (evaporated) from the topsoil [1]. It results that K_e is then computed daily with:

$$K_e= K_r\left(K_{c max}-K_{cb}\right)\le f_{ew}{K}_{c max}$$

(5)

as a function of K_{c max}, the maximum value of K_c following rain or an irrigation event and K_cb in that day (Equation (1)) after adjustment with K_r. K_e is limited to f_ew K_{c max} that represents the fraction of K_{c max} that refers to the fraction f_ew of the soil that is both exposed and wetted, i.e., the fraction of soil surface from which most evaporation occurs. The background and computation of K_cb, K_s and K_e in SIMDualKc is described in detail by Rosa et al. [32].

The calibration of SIMDualKc aims at optimizing the crop parameters K_cb and p relative to the various crop growth stages, as well as the soil evaporation parameters, the deep percolation parameters, and the runoff curve number (CN) using trial and error procedures until small errors were achieved. Further information is provided by [8,34,35]. For the olive orchard, the calibration was performed comparing daily transpiration data simulated with those obtained with sapflow measurements (and also eddy covariance and energy balance measurements as auxiliary data for calibration and data quality screening), as described in Paço et al. [8]. For the maize and barley crops, the calibration was performed by minimizing the differences between observed and simulated soil water content, as described in Paredes et al. [34] and Pereira et al. [35]. In all cases, the results of goodness of fit indicators showed good model performance and small errors of estimates [8,34,35]. Furthermore, ET partitioning was performed in other studies, e.g., testing the soil evaporation algorithm through comparing simulated with microlysimeter observations [45,46], crop transpiration against sap flow measurements [8,40] and crop ET compared with eddy covariance ET data [47]. An example of the results of the model calibration for barley is presented in Figure 2, which shows that the simulated soil water content accurately follows the dynamics of observations throughout the barley crop season [35].

Figure 2. Observed ($•$) and simulated (—–) daily soil water content (SWC) after the calibration of SIMDualKc for barley in the dry year (error bars correspond to the standard deviation of SWC observations); θ_Sat, θ_FC, θ_p and θ_WP are respectively the SWC at saturation, field capacity, wilting point and at non-stressed depletion fraction p). Precipitation is represented by dark bars and irrigation by the light ones (source: [35]).

Figure 2. Observed ($•$) and simulated (—–) daily soil water content (SWC) after the calibration of SIMDualKc for barley in the dry year (error bars correspond to the standard deviation of SWC observations); θ_Sat, θ_FC, θ_p and θ_WP are respectively the SWC at saturation, field capacity, wilting point and at non-stressed depletion fraction p). Precipitation is represented by dark bars and irrigation by the light ones (source: [35]).

2.4. Basal Crop Coefficients Derived from Reflectance Vegetation Indices

The proposed methodology is based on the dual crop coefficient approach, thus integrating daily K_cb and K_e, and considering K_s when water stress occurs [1,12]. Considering this approach, a new equation was developed to estimate K_cb from the density coefficient K_d and K_{c min} as proposed by Allen and Pereira [12] and from vegetation indices derived from satellite images, K_{cb VI}:

$$K_{cb VI}=K_{c min}+ K_d\left(\frac{V{I}_i-V{I}_{min}}{V{I}_{max}-V{I}_{min}}\right)$$

(6)

where VI_i corresponds to the VI for a specific date and pixel, VI_max is the VI for maximum vegetation cover and VI_min is VI for minimum vegetation cover (bare soil). This equation adapts to crops with continuous ground cover and to tree crops with bare soil by including K_d. For tree crops having an active ground cover, K_{c min} is replaced by a K_{cb cover}, corresponding to the K_cb of the ground cover in the absence of tree foliage [12]. K_d in Equation (6) is computed with Equation (2) using the effective fraction of ground covered estimated with the Equation (7) presented below. Table 3 presents K_d ranges for the three crops during the mid-season stage.

Table 3. Range of values of density coefficient computed for the satellite image dates of the mid-season stage in the three crops.

Kd	Maize Study Areas	Barley Study Areas (3)	Olive Study Areas
Kd (fc field) (1)	[0.93–0.99]	0.90	[0.62–0.78]
Kd (fc VI) (2)	[0.97–1.00]	0.90	[0.66–0.70]

⁽¹⁾ K_d computed using effective fraction of ground covered or shaded by vegetation derived from ground data; ⁽²⁾ K_d computed using effective fraction of ground covered or shaded by vegetation estimated with vegetation indices using Equation (7); ⁽³⁾ A single value is presented because a single satellite image was available for the mid-season stage.

Table 3. Range of values of density coefficient computed for the satellite image dates of the mid-season stage in the three crops.

Kd	Maize Study Areas	Barley Study Areas (3)	Olive Study Areas
Kd (fc field) (1)	[0.93–0.99]	0.90	[0.62–0.78]
Kd (fc VI) (2)	[0.97–1.00]	0.90	[0.66–0.70]

⁽¹⁾ K_d computed using effective fraction of ground covered or shaded by vegetation derived from ground data; ⁽²⁾ K_d computed using effective fraction of ground covered or shaded by vegetation estimated with vegetation indices using Equation (7); ⁽³⁾ A single value is presented because a single satellite image was available for the mid-season stage.

Both NDVI and SAVI vegetation indices may be used. The SAVI index is considered less sensitive to soil background reflectance variability than other VI, thus with higher potential for application in tree crops with discontinuous ground cover. NDVI_max was set to 0.75–0.85 for the three crops, and the value 0.10 was set for NDVImin, both based on reference values for, respectively, an irrigated agricultural field with full cover and a bare soil [48,49,50]. When SAVI was considered, the SAVI_min was set to 0.09 while the SAVI_max was adjusted to 0.75 according to Mateos et al. [7].

The K_{c min} was set to a value between 0.10 and 0.15 (0.13 was considered in this study) according to values proposed by Allen and Pereira [12]. For the olive orchard, K_{cb cover} was considered instead of K_{c min} and set to a value between 0.15 and 0.18 (0.17 was considered in this study) due to the presence of sparse vegetation cover between tree rows. Reference values for several parameters considered in the computation of Equation (6) can be found in Allen and Pereira [12]. K_d was computed with Equation (2) where, in a first approach, the f_{c eff} was the same considered in SIMDualKc, which was based in field observations (f_{c field}). In a second approach, f_{c eff} was computed using an equation based on vegetation indices derived from remote sensing reflectance data (f_{c VI}):

$$f_{c VI}= {\beta }_1\left(\frac{V{I}_i-V{I}_{min}}{V{I}_{max}-V{I}_{min}}\right)+{\beta }_2$$

(7)

where β₁ is an empirical coefficient (ranging between 0 and 1) depending upon the maximum VI value in each crop stage. VI_i is the average VI of the study area pixels for each date, VI_max and VI_min correspond to the VI respectively for maximum and minimum vegetation cover. Both NDVI and SAVI may be used, with the latter producing better results for incomplete cover crops, e.g., olive orchards. NDVI_max and NDVI_min were set according to references in the literature [48,49,50] (Table 4). For the super high density olive orchard, SAVI_max and SAVI_min were considered following values proposed by Mateos et al. [7] (Table 4). The parameter β₂ corresponds to an adjustment coefficient associated with crop senescence and leaves yellowing. Both NDVI and SAVI are sensitive to leaf senescence thus yielding smaller values at the end of the growing cycle than before peak development was reached [51]. Therefore, the β₂ adjustment is considered to compensate the decrease in VI due to the senescence and/or yellowing of plants, which is independent from f_c. Values of β₂ are presented in Table 4.

Table 4. Values of parameters considered in the estimation of fraction of ground cover from vegetation indices.

Parameters	Crop Growth Stage	Value
		Maize	Barley	Olive
NDVImax		0.75–0.85	0.75–0.85	0.75–0.85
NDVImin		0.1	0.1	0.1
SAVImax		0.75	0.75	0.75
SAVImin		0.09	0.09	0.09
β2		0–0.5 (2)	0–0.5 (2)	0–0.5 (2)
β1(1)	Ini	0.3	0.3	0.4–0.7
	Dev	0.3–1	0.3–1	1
	Dev (1st sub-stage)	0.6	0.5	1
	Dev (2nd sub-stage)	0.9	0.7	-
	Dev (3rd sub-stage)	-	0.9	-
	Mid	1	1	-
	End	1	1	1

⁽¹⁾ The values were defined according to stages of the crop growth. For the development stage were considered sub-stages by dividing the growth stage in two (spring-summer crops; maize) or three (winter-spring crops; barley) periods; ⁽²⁾ 0 is considered when the crop is in its maximum development and 0.5 close to the harvest (a value of 0.2 was set in the first half of the late stage).

Table 4. Values of parameters considered in the estimation of fraction of ground cover from vegetation indices.

Parameters	Crop Growth Stage	Value
		Maize	Barley	Olive
NDVImax		0.75–0.85	0.75–0.85	0.75–0.85
NDVImin		0.1	0.1	0.1
SAVImax		0.75	0.75	0.75
SAVImin		0.09	0.09	0.09
β2		0–0.5 (2)	0–0.5 (2)	0–0.5 (2)
β1(1)	Ini	0.3	0.3	0.4–0.7
	Dev	0.3–1	0.3–1	1
	Dev (1st sub-stage)	0.6	0.5	1
	Dev (2nd sub-stage)	0.9	0.7	-
	Dev (3rd sub-stage)	-	0.9	-
	Mid	1	1	-
	End	1	1	1

⁽¹⁾ The values were defined according to stages of the crop growth. For the development stage were considered sub-stages by dividing the growth stage in two (spring-summer crops; maize) or three (winter-spring crops; barley) periods; ⁽²⁾ 0 is considered when the crop is in its maximum development and 0.5 close to the harvest (a value of 0.2 was set in the first half of the late stage).

The four main stages of the crop growth, described by Allen et al. [1], were considered in the definition of β₁: (i) Ini, (ii) Dev, (iii) Mid, and (iv) Late. Due to the great variability of the ground cover by the vegetation throughout Dev, a differentiation of β₁ values according to sub-stages was considered: (i) two sub-stages in spring-summer crops (maize); (ii) three sub-stages in winter-spring crops (barley), where the Dev is longer. For the olive orchard, a crop with persistent leaves and where the soil cover by the vegetation is incomplete, the β₁ value was set to 1 throughout the year, except for the period of pruning and in the subsequent weeks, when it was slightly decreased to account for the decrease in canopies vegetation. Values of β1 are presented in Table 4.

Most of the optical VI track the effects of long-term water stress on plants but do not allow the early detection of water stress [52]. Therefore, Equation (6) expresses the potential K_cb, thus transpiration occurring at maximum rate. However, the occurrence of water stress cause a decrease of the potential value and K_s must be applied to obtain actual K_{cb VI}. In the current approach, K_s computed with SIMDualKc was used. The dates when water stress was detected through the daily soil water balance are given in Table 2. Therefore, hereafter, the actual K_{cb VI} refers to the K_{cb VI} computed by Equation (6) multiplied by the K_s computed with SIMDualKc.

For K_c estimation from VI, hereafter named K_{c VI}, the K_e computed with SIMDualKc was combined with the K_{cb VI}, i.e., K_{c VI} = K_e + K_{cb VI}. A similar combined approach was reported by Mateos et al. [7]. The results of K_{cb VI} and K_{c VI} were compared with the actual K_cb and K_c computed with SIMDualKc based on field data, later referred as K_{cb SIMDualKc} and K_{c SIMDualKc}, respectively.

Because K_s and K_e used in the estimation of actual K_{c VI} are computed with SIMDualKc, the results of comparing actual K_{c VI} with K_{c SIMDualKc} are somehow influenced by the fact that values of K_s and K_e used in both cases are the same. Nevertheless, this fact does not impede the assessment of the computation of the actual K_{cb VI} because this is the procedure that gives the possibility to compare the actual crop coefficients from both approaches and retrieves useful information from satellite imagery aiming at irrigation management and farmer advising.

A linear regression forced through the origin was used for the comparisons and the corresponding regression and determination coefficients were consequently used as indicators. Additionally, the root-mean-square deviation (RMSD; Equation (8)) and the relative mean deviation (RMD; Equation (9)) were used to compare two sets of data, with Xi relative to estimations with SIMDualKc, having mean $\overline{X}$ , and Y_i relative to the VI estimates (i = 1,2, ..., n).

$$RMSD= \sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(X_i-Y_i\right)}^2}{n}}$$

(8)

$$RMD=\frac{1}{n}{\displaystyle \sum }_{i=1}^n\frac{\left|X_i-Y_i\right|}{\overline{X}} 100$$

(9)

Since the database available for this study only comprises 15 VI observations for maize and olive and seven for barley (Table 2), the validation of K_cb estimated based on VI was performed using the “leave-one-out” (LOO) cross-validation [53]. Applications of the LOO cross-validation in remote sensing studies are provided by [54,55]. The LOO cross-validation evaluates the model performance for observations not considered in the estimation step, thus providing independent estimates of the predictive capability of the selected models. This technique consists in the removal of one observation from the dataset used and the estimation of a new regression model with the remaining observations. This new regression model is used to estimate the K_{cb VI} of the observation withdrawn (K_{cb VI LOO}). Equations (8) and (9) were also used to estimate the goodness of cross-validation results, with X_i and Y_i representing respectively K_{cb SimDualKc} and K_{cb VI LOO} data. A similar approach was considered for the K_c and f_c estimated based on VI.

3. Results

3.1. Estimation of the Fraction of Ground Cover from VI

The results of f_c estimated from NDVI (f_{c VI}; Equation (7)) were compared with f_c based on field observations (f_{c field}) for the maize and barley crops (Figure 3). The comparison of f_{c VI} and f_{c field} have shown a good agreement for both crops, with R² ≥ 0.81 and the regression coefficient b close to 1.0 (Figure 3). For maize, RMSD was 0.10 and RMD 10.6%, while for barley RMSD was 0.06 and RMD was 8.1%.

Figure 3. Comparison of the fraction of ground cover by vegetation derived from NDVI (f_{c VI}) with the f_c obtained from field observations (f_{c field}) for (a) maize and (b) barley.

Figure 3. Comparison of the fraction of ground cover by vegetation derived from NDVI (f_{c VI}) with the f_c obtained from field observations (f_{c field}) for (a) maize and (b) barley.

For the olive crop, due to the low variability of f_{c field} values, a regression for comparing f_{c VI} and f_{c field} was not considered. The fraction of ground cover for this tree crop with persistent leaves is nearly constant throughout the year, except for the pruning period. As such, for the analysis of the f_{c VI} results, only the indicators of residual estimation deviations were used: RMSD = 0.22 and RMD = 69.2%. When SAVI was used for the estimation of f_{c VI} results of all indicators improved much, with RMSD= 0.02 and RMD = 5.5%. The absolute mean differences between the f_{c VI} and the f_{c field} for the olive orchard were small, ranging from 0.003 to 0.046. Contrarily, for the annual crops (maize and barley), the use of SAVI in the estimation of f_{c VI} did not improve the residual estimation deviations.

Results of cross validation show a good performance of the f_{c VI} estimation, with RMSD of 0.06. 0.04, and 0.01 respectively for maize, barley and olive. RMD results for the cross validation varied between 2.6% for olive and 7.5% for maize.

3.2. Estimation of K_cb from Reflectance-Based Vegetation Indices (K_{cb VI})

The K_{cb VI} (Equation (6)) was adjusted using the K_s derived from SIMDualKc computations for the dates when stress occurred (Table 2) to obtain the actual K_{cb VI}. The values of K_{cb SIMDualKc} used in comparisons with K_{cb VI} were also adjusted to water stress using K_s. Both NDVI and SAVI were tested for the estimation of K_{cb VI} and the results obtained with the two VI were compared.

The comparison of actual K_{cb VI} with actual K_{cb SIMDualKc} for the studied crops is presented in Figure 4. The K_{cb VI} values in Figure 4 were estimated using f_{c NDVI} (Equation (7)) with the objective of fully estimating actual K_cb from reflectance-based vegetation indices. Figure 4a shows K_{cb VI} for maize when using NDVI, which has a determination coefficient (R²) lower than using SAVI (Figure 4b) but has a higher regression coefficient (b = 0.98) relative to SAVI (Table 5). Thus, K_{cb SAVI} underestimates K_cb values. For barley (Figure 4c,d) both b and R² are better when using NDVI (Table 5). For the olive orchard estimating K_{cb VI} with SAVI leads to b equal to 1.22 and a R² equal to 0.83; when using NDVI there is an overestimation of the K_cb values, b = 1.8 (Figure 4e,f).

Figure 4. Comparison of actual K_cb derived from NDVI (K_{cb NDVI}), on left, or from SAVI (K_{cb SAVI}), on right, with the actual K_cb obtained with SIMDualKc (K_{cb SIMDualKc}) for maize (a,b), barley (c,d) and the olive crop (e,f).

Figure 4. Comparison of actual K_cb derived from NDVI (K_{cb NDVI}), on left, or from SAVI (K_{cb SAVI}), on right, with the actual K_cb obtained with SIMDualKc (K_{cb SIMDualKc}) for maize (a,b), barley (c,d) and the olive crop (e,f).

The statistical indicators relative to the comparison between K_{cb VI} and K_{cb SIMDualKc} for all crops are presented in Table 5. The RMSD between K_{cb VI} and K_{cb SIMDualKc} for the set of image dates considered for maize was equal to 0.16 when using NDVI for the computation of K_{cb VI} and 0.22 when SAVI was adopted. RMD was also higher when adopting SAVI. The regression coefficient b = 0.98 when NDVI is used indicates a statistical similarity between K_{cb VI} and K_{cb field}. For barley, b and R² are higher for K_{cb NDVI} and RMSD and RMD are smaller than for K_{cb SAVI}. For the olive crop there is a clear overestimation of K_{cb NDVI} and K_{cb SAVI} when the f_{c VI} is computed with NDVI (1.81 and 1.19, respectively; Table 5).

Table 5. Statistical indicators relative to the comparison between K_{cb VI} and K_{cb SIMDualKc} for maize and barley when f_{c VI} computed with NDVI and for olive crops when f_{c VI} is computed with NDVI and SAVI *.

		b	R2	n	p	Estimation		Cross Validation
						RMSD ( )	RMD (%)	RMSD ( )	RMD (%)
Maize	Kcb NDVI	0.98	0.74	15	0.0000389	0.16	14.8	0.08	8.4
	Kcb SAVI	0.79	0.79	15	0.0000094	0.22	21.1	0.18	20.9
Barley	Kcb NDVI	0.99	0.95	7	0.0001935	0.07	9.6	0.04	4.7
	Kcb SAVI	0.67	0.89	7	0.0014203	0.23	32.0	0.22	31.3
Olive	Kcb NDVI	1.81	0.56	15	0.0013317	0.31	81.1	0.30	82.5
	Kcb SAVI	1.19	0.72	15	0.0000023	0.08	21.7	0.07	19.0
	Kcb SAVI *	1.03	0.88	15	0.0000002	0.02	5.6	0.01	3.3

* f_c estimated from SAVI.

Table 5. Statistical indicators relative to the comparison between K_{cb VI} and K_{cb SIMDualKc} for maize and barley when f_{c VI} computed with NDVI and for olive crops when f_{c VI} is computed with NDVI and SAVI *.

		b	R2	n	p	Estimation		Cross Validation
						RMSD ( )	RMD (%)	RMSD ( )	RMD (%)
Maize	Kcb NDVI	0.98	0.74	15	0.0000389	0.16	14.8	0.08	8.4
	Kcb SAVI	0.79	0.79	15	0.0000094	0.22	21.1	0.18	20.9
Barley	Kcb NDVI	0.99	0.95	7	0.0001935	0.07	9.6	0.04	4.7
	Kcb SAVI	0.67	0.89	7	0.0014203	0.23	32.0	0.22	31.3
Olive	Kcb NDVI	1.81	0.56	15	0.0013317	0.31	81.1	0.30	82.5
	Kcb SAVI	1.19	0.72	15	0.0000023	0.08	21.7	0.07	19.0
	Kcb SAVI *	1.03	0.88	15	0.0000002	0.02	5.6	0.01	3.3

* f_c estimated from SAVI.

For the olive crop, however, estimating K_{cb VI} using f_{c SAVI} (see Equations (6) and (7)) provided much better results (Figure 5 and Table 5; K_{cb SAVI *}) than using K_{cb VI} with f_{c NDVI} (Table 5). Comparing K_{cb SAVI} with K_{cb SIMDualKc} it results a higher R² = 0.88 with b of 1.03 (Table 5). In addition, RMSD decreased from 0.08 to 0.02 and RMD from 21.7 to 5.6% (Table 5). RMSD and RMD obtained in the cross validation procedure for the three crops are also presented in Table 5. The results obtained for cross validation are consistent with those obtained for the estimation procedure.

Figure 5. Comparison of actual K_cb obtained with SIMDualKc (K_{cb SIMDualKc}) with K_cb derived from SAVI (K_{cb SAVI}) for olive when f_c is estimated from SAVI.

Figure 5. Comparison of actual K_cb obtained with SIMDualKc (K_{cb SIMDualKc}) with K_cb derived from SAVI (K_{cb SAVI}) for olive when f_c is estimated from SAVI.

Table 6 presents the average and standard error of the mean or the range of the K_{cb VI} obtained for all fields and all crops. Results refer to the use of both NDVI and SAVI with Equation (6). As for the previous analysis, the best results refer to K_{cb NDVI} in case of maize and barley and K_{cb SAVI} in case of the olive orchard. Comparing results for both VIs in this Table 6 it becomes apparent that estimations using NDVI lead to higher K_{cb VI} results, which relates to the nature of the VI index.

Table 6. Average and standard error of the mean or the range of the actual daily basal crop coefficients estimated with vegetation indices for the various crops and crop growth stages.

Crop Growth Stages	Maize	Barley	Olive
Kcb NDVI
Initial	(1)	0.18 ( ) (2)	0.78 (±0.13)
Development (4)	[0.44–1.00]	[0.25–0.91]	[0.53–0.77] (3)
Mid-season	0.92 (±0.11) (5)	0.87 ( ) (2)	0.61 (±0.05) (6)
Late-season (4)	[0.84–0.46]	[0.86–0.69]	[0.75–0.55]
Kcb SAVI
Initial	(1)	0.15 ( ) (2)	0.33 (±0.01)
Development (4)	[0.33–0.76]	[0.20–0.63]	[0.30–0.44] (3)
Mid-season	0.76 (±0.10) (5)	0.54 ( ) (2)	0.37 (±0.03) (6)
Late-season (4)	[0.52–0.33]	[0. 56–0.52]	[0.40–0.33]

Note: f_{c VI} were obtained from NDVI for maize and barley and from SAVI for olive. ⁽¹⁾ No Landsat images available for this crop growth stage; ⁽²⁾ Only one Landsat image available for this crop growth stage; ⁽³⁾ Period after pruning; ⁽⁴⁾ Range of values for this crop growth stage; ⁽⁵⁾ Stress occurred in three out of the eight image dates considered; ⁽⁶⁾ Stress occurred in four out of the eight image dates considered.

Table 6. Average and standard error of the mean or the range of the actual daily basal crop coefficients estimated with vegetation indices for the various crops and crop growth stages.

Crop Growth Stages	Maize	Barley	Olive
Kcb NDVI
Initial	(1)	0.18 ( ) (2)	0.78 (±0.13)
Development (4)	[0.44–1.00]	[0.25–0.91]	[0.53–0.77] (3)
Mid-season	0.92 (±0.11) (5)	0.87 ( ) (2)	0.61 (±0.05) (6)
Late-season (4)	[0.84–0.46]	[0.86–0.69]	[0.75–0.55]
Kcb SAVI
Initial	(1)	0.15 ( ) (2)	0.33 (±0.01)
Development (4)	[0.33–0.76]	[0.20–0.63]	[0.30–0.44] (3)
Mid-season	0.76 (±0.10) (5)	0.54 ( ) (2)	0.37 (±0.03) (6)
Late-season (4)	[0.52–0.33]	[0. 56–0.52]	[0.40–0.33]

Note: f_{c VI} were obtained from NDVI for maize and barley and from SAVI for olive. ⁽¹⁾ No Landsat images available for this crop growth stage; ⁽²⁾ Only one Landsat image available for this crop growth stage; ⁽³⁾ Period after pruning; ⁽⁴⁾ Range of values for this crop growth stage; ⁽⁵⁾ Stress occurred in three out of the eight image dates considered; ⁽⁶⁾ Stress occurred in four out of the eight image dates considered.

Figure 6 shows selected examples of the seasonal variation of the daily K_cb, K_{cb act}, K_e and K_{c act} obtained with SIMDualKc for maize [34], barley [35] and olive [8] and the fitting of actual K_{cb VI} obtained for those sets of data. K_{cb VI} used in Figure 6 refer to K_{cb NDVI} in case of maize and barley and to K_{cb SAVI} in case of olive. For the olive orchard, the f_{c SAVI} (Equation (7)) were used in the computation of K_d while f_{c NDVI} were used for maize and barley. Unfortunately, there were few images for maize crop in the year when water stress occurred. Results show that actual K_{cb VI} match quite well the actual K_cb curve for the three crops. The variability of weather conditions in those years are well apparent through the numerous peaks of the evaporation coefficient K_e due to both rainfall and irrigation events. The effect of water stress on K_cb is well apparent by the deviations of K_{cb act} from the potential K_cb, thus when the K_{cb act} curve lays below the K_cb curve. Summarizing, Figure 6 demonstrates the goodness of the approach used to estimate K_{cb act} from a vegetation index.

Figure 6. Seasonal variation of the daily coefficients K_cb, actual K_cb (K_{cb act}), K_e and actual K_c (K_{c act}) obtained with SIMDualKc, and actual K_{cb VI} and actual K_c estimated by the combined approach (K_{c VI}) for: (a) maize, field 1—Year 2010; (b) barley, year 2012; and (c) olive orchard, year 2012.

Figure 6. Seasonal variation of the daily coefficients K_cb, actual K_cb (K_{cb act}), K_e and actual K_c (K_{c act}) obtained with SIMDualKc, and actual K_{cb VI} and actual K_c estimated by the combined approach (K_{c VI}) for: (a) maize, field 1—Year 2010; (b) barley, year 2012; and (c) olive orchard, year 2012.

3.3. Estimation of Actual K_c from Reflectance-Based Vegetation Indices (K_{c VI})

The actual K_{c VI} values were obtained by combining the actual K_{cb VI} with K_e derived from SIMDualKc. The results for K_{c VI} presented in Figure 7 and Table 7 were obtained using the K_{cb VI} that provided the best estimation for each crop as presented in Section 3.2. Thus, K_{cb NDVI} was considered for the computation of K_{c VI} for maize and barley crops, while K_{cb SAVI} was adopted for olive. Figure 6 shows the K_{c VI} corresponding to the same dates when the K_{cb VI} were obtained for selected data sets referred in Section 3.2. It may be observed that despite the enormous variability of K_e due to rainfall and irrigation wetting events the computed K_{c VI} match well the actual K_c curves for the three crops analyzed.

Table 7. Statistical indicators relative to the comparison between actual K_{c VI} and K_{c SIMDualKc}.

	b	R2	n	p	Estimation		Cross Validation
					RMSD ( )	RMD (%)	RMSD ( )	RMD (%)
Maize	0.99	0.72	15	0.0000637	0.16	12.7	0.08	6.9
Barley	0.99	0.83	7	0.0043194	0.07	6.4	0.03	2.6
Olive	1.01	0.99	15	0.0000000	0.02	3.3	0.01	2.1

Table 7. Statistical indicators relative to the comparison between actual K_{c VI} and K_{c SIMDualKc}.

	b	R2	n	p	Estimation		Cross Validation
					RMSD ( )	RMD (%)	RMSD ( )	RMD (%)
Maize	0.99	0.72	15	0.0000637	0.16	12.7	0.08	6.9
Barley	0.99	0.83	7	0.0043194	0.07	6.4	0.03	2.6
Olive	1.01	0.99	15	0.0000000	0.02	3.3	0.01	2.1

Figure 7. Comparison of actual K_c derived from SIMDualKc (K_{c SIMDualKc}) and actual K_c derived from the combined approach of K_{cb VI} and K_e from SIMDualKc (K_{c VI}) for (a) maize; (b) barley; and (c) olive.

Figure 7. Comparison of actual K_c derived from SIMDualKc (K_{c SIMDualKc}) and actual K_c derived from the combined approach of K_{cb VI} and K_e from SIMDualKc (K_{c VI}) for (a) maize; (b) barley; and (c) olive.

When comparing K_{c VI} and K_{c SIMDualKc} a high determination coefficient was obtained (R² ≥ 0.72) together with a regression coefficient very close to 1.0 for the three crops studied (Figure 7 and Table 7). The best statistical indicators were obtained for the olive orchard (Table 7), particularly with quite small RMSD and RMS. RMSD and RMD obtained in the cross validation procedure for the three crops show results very similar to those obtained for the estimation procedure (Table 7).

The average and range of the actual crop coefficients estimated with vegetation indices combined with the computed K_e (K_{c VI}) are given in Table 8 for the various crops.

Table 8. Average and standard error of the mean or the range of actual crop coefficients estimated with vegetation indices approach (K_{c VI}) for the various crop growth stages of each crop.

Crop Growth Stages	Maize	Barley	Olive
Initial	(1)	0.78 ( ) (2)	0.76 (±0.41) (3)
Development (4)	[0.90–1.30] (5)	[0.75–1.18]	[0.44–0.82] (6)
Mid-season	0.98 (±0.11) (7)	0.88 ( ) (2)	0.62 (±0.08) (8)
Late-season (4)	[0.84–0.55]	[1.03–0.73]	[0.53–0.40]

Note: f_{c VI} were obtained from NDVI for maize and barley and from SAVI for olive. ⁽¹⁾ No Landsat images available for this crop growth stage; ⁽²⁾ Only one Landsat image available for this crop growth stage; ⁽³⁾ High evaporation occurrence in one of the image dates considered in the initial stage due to precipitation occurrence. ⁽⁴⁾ Range of daily Kc values (not time-averaged); ⁽⁵⁾ High evaporation in two of the dates considered in this crop growth stage; ⁽⁶⁾ Period after pruning; ⁽⁷⁾ Stress occurred in three out of the eight image dates considered; ⁽⁸⁾ Stress occurred in four out of the eight image dates considered.

Table 8. Average and standard error of the mean or the range of actual crop coefficients estimated with vegetation indices approach (K_{c VI}) for the various crop growth stages of each crop.

Crop Growth Stages	Maize	Barley	Olive
Initial	(1)	0.78 ( ) (2)	0.76 (±0.41) (3)
Development (4)	[0.90–1.30] (5)	[0.75–1.18]	[0.44–0.82] (6)
Mid-season	0.98 (±0.11) (7)	0.88 ( ) (2)	0.62 (±0.08) (8)
Late-season (4)	[0.84–0.55]	[1.03–0.73]	[0.53–0.40]

Note: f_{c VI} were obtained from NDVI for maize and barley and from SAVI for olive. ⁽¹⁾ No Landsat images available for this crop growth stage; ⁽²⁾ Only one Landsat image available for this crop growth stage; ⁽³⁾ High evaporation occurrence in one of the image dates considered in the initial stage due to precipitation occurrence. ⁽⁴⁾ Range of daily Kc values (not time-averaged); ⁽⁵⁾ High evaporation in two of the dates considered in this crop growth stage; ⁽⁶⁾ Period after pruning; ⁽⁷⁾ Stress occurred in three out of the eight image dates considered; ⁽⁸⁾ Stress occurred in four out of the eight image dates considered.

4. Discussion

4.1. Estimation of the Fraction of Ground Cover from Vegetation Indices

The statistical results relative to the estimation of f_{c VI} indicate a good performance of Equation (7) for the annual crops studied (maize and barley). The R² values observed when comparing f_{c VI} with f_{c field} are highly significant for both crops (Figure 3). These R² values are within the range of the values reported by Trout et al. [56] for the estimation of canopy cover using NDVI derived from a multispectral camera for a set of 11 annual and perennial horticultural crops. Purevdorj et al. [57] reported similar R² results when estimating the percentage of vegetation cover in grasslands and lawn grass using NDVI derived from AVHRR satellite images, thus with a coarser spatial resolution. Estimated deviations are small, for both maize and barley, with RMSD of 0.10 and 0.06, respectively.

The estimation of f_c in the olive orchard also showed good estimation, with low RMD < 10%. When the two years studied (2011 and 2012) were analyzed separately, the average absolute deviations between f_{c VI} and f_{c field} were slightly higher in 2012 than in 2011. These results are likely due to the larger variability in the vegetation conditions within the orchard in 2012 because a severe frost occurred in February which was followed by an heavy and non-uniform pruning that increased the variability of f_c within the orchard as discussed by Paço et al. [8].

The better performance of SAVI for the estimation of f_{c VI} in the high density olive orchard, having a large fraction of soil exposed, is likely due to the lower sensitivity of this VI to soil background influence [23]. Similar conclusions were reported by Purevdorj et al. [57] relative to the better performance of SAVI to estimate the vegetation cover for low crop densities. Differently, in a study with savannah-type open woodlands dominated by evergreen oak species (montados or dehesa), the NDVI performed slightly better than SAVI for the estimation of tree canopy cover [58] but the study areas were characterized by an extremely variable understory including bare soil, dry grass and a few evergreen shrubs that impacted the VIs values as discussed by the authors.

4.2. Estimation of Actual K_cb from Reflectance-Based Vegetation Indices (K_{cb VI})

The R² results relative to the regression comparing the actual K_{cb VI} with K_{cb SIMDualKc} for maize were highly significant using NDVI (Table 5). Although using NDVI led to a lower R² than SAVI, the regression coefficient comparing K_{cb SAVI} with K_{cb SIMDualKc} (b = 0.79; Table 5) indicates an underestimation of K_{cb SAVI} values. The regression coefficient obtained when comparing K_{cb NDVI} with K_{cb SIMDualKc} is much better (b = 0.98; Table 5), which may be related with the tendency of NDVI reaching a maximum value at about the same time as K_cb, i.e., at a LAI around 3 [14,59]. Contrarily, SAVI tends to continue increasing for LAI above 3 [21,25,60]. Thus, according to Allen et al. [21], the NDVI supports better than SAVI the estimation of K_c and K_cb in crops reaching high LAI values as for the full cover crops considered in the current study.

The base data on maize refers to three different seasons and three different fields, which results in a large variability of soil, climate and crop management conditions. A large variability of K_{cb act} was therefore observed despite the potential/standard K_cb was the same after appropriate calibration of SIMDualKc [34]. This variability contributed for a larger variance of the actual K_{cb VI} and K_{cb SIMDualKc} and made it more difficult to parameterize a single value for K_{c min} and β₁ (Equations (6) and (7), respectively) when considering all three maize seasons and three cropped fields. The variability conditions observed in the three different fields and seasons during the initial and development stages may justify why comparing K_{cb NDVI} with K_{cb SIMDualKc} leads to a smaller R² (Table 5).

The actual K_{cb NDVI} obtained for the maize mid-season stage averaged 0.92 (Table 6). This value relates with the water stress that occurred during the mid-season stage of 2010 in both surveyed fields which was detected in two of the images, particularly that of 30 July 2010. The SIMDualKc computed K_s values for this date were 0.55 and 0.29 respectively for field 1 and field 2. Such low K_s values highly impacted actual K_cb for this stage. Figure 6a shows well the impact of water stress during that period, with the K_{cb act} curve largely below the K_cb curve. The actual K_{cb NDVI} are within the range of values presented by Padilla et al. [31]. Relative to the late season, the range of K_{cb NDVI} values in Table 6 refer to three images available for this period. One of these images was captured the day prior to harvest (19 September 2011 in field 1) and the K_{cb NDVI} was then equal to 0.48, therefore similar to the K_{cb end} proposed by Allen et al. [1], K_{cb end} = 0.5. However, values for K_{cb end} relate to crop management and therefore are difficult to be compared with other values in literature.

For the regression comparing the actual K_{cb VI} and K_{cb SIMDualKc} for barley the statistical indicators show that NDVI performed better, both in terms of the regression and determination coefficients and RMSD and RMD (Table 5 and Figure 4c,d). The actual K_{cb NDVI} obtained for the initial crop growth stage (actual K_{cb ini NDVI} = 0.18; Table 6) is just slightly higher than the value presented by Allen et al. [1], K_{cb ini} = 0.15, which relates with the fact that no water stress occurred in that stage (Figure 6b). K_{cb NDVI} values for the development growth stage ranged between 0.25 and 0.91 (Table 6), which refer to values commonly accepted [1]. The actual K_{cb mid NDVI} = 0.87 (Table 6) is lower than the standard K_cb = 1.10 proposed by Allen et al. [1] and values obtained for crops with similar canopy [30,31,61]. This is likely due to different vegetation conditions in the 2012/2013 crop season and to the adjustment to climate [35]. The K_{cb NDVI} presented in Table 6 for the late season ranged between 0.86 and 0.69 and refer to two images available for that stage. These images provided for likely appropriate actual K_cb values and compared well with values obtained with SIMDualKc (Figure 6b, [35]).

For the olive orchard, results obtained with SAVI are good (Figure 5) and in agreement with the above mentioned lower sensitivity of this VI to changes in soil background [23], which is particularly relevant for tree crops with discontinuous ground cover. Furthermore, SAVI improves the linearity between the VI and biophysical parameters as observed by other researchers [25,56]. Mateos et al. [7] also obtained good results when using SAVI for estimating K_cb in peach and mandarin orchards. Contrarily, in a study with grapevine, Campos et al. [13] did not find significant improvements in the accuracy when using SAVI but the application excluded days with wet soil. Accurate results for olive in the current study are likely due to the use of K_d in Equation (6).

The average K_{cb SAVI} value obtained for the mid-season stage in olive (0.37; Table 6 and Figure 5c) was affected by the occurrence of water stress in four out of the eight dates considered in this stage, with K_s < 0.77 in two of these dates. The impacts of water stress on K_cb is clearly shown in the example presented in Figure 6c, with the K_{cb act} curve laying bellow the K_cb curve during the mid-season. The actual K_{cb SAVI} for this stage would increase to 0.42 when only dates without water stress or with very mild stress (K_s > 0.95) are considered. The actual K_{cb SAVI} values for olive orchards for the mid- and late-seasons, respectively 0.37 and ranging 0.40 to 0.33 (Table 6), are within the range of values for olive orchards with an effective fraction of ground cover (fc eff) between 0.25 and 0.5, with a K_{cb mid} ranging between 0.35 and 0.55 and K_{cb end} ranging from 0.30 to 0.50 [12]. During the initial stage, some sparse vegetation was growing in the inter-rows, thus K_cb was only slightly influenced by this sparse vegetation. Consequently, K_{cb ini SAVI} = 0.33 was an appropriate value for this orchard.

Overall results obtained have shown a good performance of the K_{cb VI} Equation (6) for all three crops, particularly integrating results of f_{c VI} when computing the density coefficient K_d used in Equation (6). It must be stressed that the use of K_d, which integrates the information on f_c and crop height, and of K_{c min} provided good results in modeling the actual K_{cb VI} (Equation (6)) for crops with very different canopy structures and ground cover conditions as those used in this study (continuous versus discontinuous ground cover) for the whole crop season including when water stress was observed.

4.3. Estimation of Actual K_c by Combining K_{cb VI} with K_e from SIMDualKc (K_{c VI})

The statistical indicators from comparing the actual K_{c VI} with K_{c SIMDualKc} show a good performance of the model for the three crops (Table 7 and Figure 7). The deviations of estimate are very low for olive and barley and slightly higher for maize. Results are in agreement with those discussed for K_{cb VI}, which is a component of K_{c VI}.

For the maize crop, the actual K_{c NDVI} values obtained for the development stage ranged between 0.90 and 1.30, which relate with a large variability of K_e values, as for the example in Figure 6a. The average K_{c NDVI} = 0.98 ± 0.11 for the mid-season (Table 8) reflects also a large variability as exemplified in Figure 6a, where one value exceeds 1.20 and the other, due to water stress, is as low as 0.40. When only cases without water stress are considered, K_{c VI} increases to 1.14, a value close to the K_c reported by Gonzalez-Piqueras et al. [26]. Regarding the K_{c NDVI} values for the late-season, values obtained are likely to occur because they are within the interval limited by the K_c values at mid-season and those expected for end season.

The actual K_{c NDVI} estimated for the initial and development stages for barley (0.78 and 0.75–1.18; Table 8) reflects the abundance of rainfall during these stages, with various K_e peaks occurring then (Figure 6b). For the mid-season stage the value obtained is lower than the values estimated with SIMDualKc (K_{c SIMDualKc} = 0.99), thus indicating an underestimation of K_c. The average value obtained for the late-season stage (0.88) is close to that obtained by Liu et al. [62] for the period between wax ripeness and harvesting (equal to 0.83).

The results obtained for the olive orchard show a very good performance of the model for estimating K_{c SAVI} (Figure 7c). The average results of K_{c SAVI} (Table 8) are within the range of those proposed by Allen and Pereira [12] for an f_{c eff} between 0.25 and 0.50, which vary between 0.80 for the K_{c ini} (considering an active ground cover), 0.40–0.60 for the K_{c mid} (no active ground cover), and 0.35–0.75 for the K_{c end} (depending on the ground cover). The results obtained throughout the year follow a pattern similar to results reported by other authors [43,63], as discussed by Paço et al. [8].

It is important to highlight that results presented herein correspond to actual K_c, thus reflecting the management options observed in the study areas.

5. Conclusions

The overall results for the three crops—Maize, barley and olive—indicate that the newly proposed Equation (7) adequately models the fraction of soil covered by vegetation throughout the crop season. The results obtained indicate a better performance of SAVI for the computation of the fraction of ground cover for the incomplete cover crop, olive, and of the NDVI for the annual crops, maize and barley.

Using the novel Equation (6), which applies a density coefficient, the estimation of actual basal crop coefficients through reflectance-based vegetation indices was successful, with good accuracy. The use of NDVI was preferred for maize and barley and SAVI was more accurate for olive. These results are consistent with the lower sensitivity of SAVI to changes of reflectance in soil background due to soil moisture. This approach allowed considering the whole crop season, including periods when the soil cover by vegetation is low, as in the initial and development growth stages, which was not considered in several studies previously published in the literature. The use of the density coefficient allowed considering crops with different canopy structures and incorporating the impact of vegetation density and height. The use of the stress coefficient computed with the soil water balance SIMDualKc made it possible to apply the methodology not only in days when soil water was adequate to sustain full plant transpiration but also in days when heavy water stress conditions occurred, contrarily to most of the approaches previously published that do not account for the crop water stress in their formulation.

The proposed methodology for estimating actual crop coefficients follows a dual crop coefficient approach by combining a remotely sensed basal crop coefficient and an evaporation coefficient computed with the soil water balance SIMDualKc. Therefore, the approach allows integrating the effect of plant transpiration with that of soil evaporation.

The approaches proposed in this study, combining VI data with soil water balance data, provide for an operational methodology to support irrigation management considering the actual cropping conditions and associating real time satellite information with irrigation scheduling models potentially applicable to diverse crop types. Nevertheless, further studies are desirable to better test and implement the methodology, particularly for other crops and using a larger number of image dates.