Water Use and Yield of Soybean under Various Irrigation Regimes and Severe Water Stress. Application of AquaCrop and SIMDualKc Models

: Data relative to two soybean seasons, several irrigation scheduling treatments, including moderate and severe deﬁcit irrigation, and rain-fed cropping were used to parameterize and assess the performance of models AquaCrop and SIMDualKc, the latter combined with the Stewart’s yield model. SIMDualKc applies the FAO56 dual crop coefﬁcient approach for computing and partitioning evapotranspiration (ET) into actual crop transpiration (T c act ) and soil evaporation (E s ), while AquaCrop uses an approach that depends on the canopy cover curve. The calibration-validations of models were performed by comparing observed and predicted soil water content (SWC) and grain yield. SIMDualKc showed good accuracy for SWC estimations, with normalized root mean square error (NRMSE) ≤ 7.6%. AquaCrop was less accurate, with NRMSE ≤ 9.2%. Differences between models regarding the water balance terms were notable, and the ET partition revealed a trend for under-estimation of T c act by AquaCrop, mainly under severe water stress. Yield predictions with SIMDualKc-Stewart models produced NRMSE < 15% while predictions with AquaCrop resulted in NRMSE ≤ 23% due to under-estimation of T c act , particularly for water stressed treatments. Results show the appropriateness of SIMDualKc to support irrigation scheduling and assessing impacts on yield when combined with Stewart’s model.


Introduction
Uruguay is characterized by a warm temperate and humid climate, where summer crops are commonly rain-fed. Due to rainfall uncertainty, supplemental irrigation is often required for achieving high yields [1]. Thus, adequate irrigation scheduling has to be considered for soybean production. Predicting soybean yield response to water is required for assessing irrigation management strategies to be adopted by farmers. Attention should be paid to the crop stages where water stress is most critical, with several studies having identified the period from flowering to grain filling as the most sensitive to water stress [2][3][4].
Crop growth and yield models are often used. The CROPGRO-Soybean model is probably the most used to simulate soybean growth and yield. It is one of the Decision Support System for Agrotechnology Transfer-Cropping System Models (DSSAT-CSM) whose features are discussed in detail by Jones et al. [5]. Because DSSAT-CSM are oriented to represent the growth and yield processes considering a variety of constraints and stresses, they are rarely used for assessing water use or for Southern Oscillation [28] and the Pacific Decadal Oscillation [29]. Local climate is warm temperate, with humid and hot summers: Cfa according to the Köppen-Geiger classification [30]. Weather daily data including maximum and minimum air temperature ( • C), solar radiation (MJ·m −2 ·d −1 ), air relative humidity (%), wind speed (m·s −1 ), and precipitation (mm) were collected by an automatic meteorological station (Vantage Pro 2TM, Davis Instruments, Hayward, CA, USA) located near the experimental fields. These data were used to compute daily reference evapotranspiration (ET o ) with the FAO Penman Monteith (FAO-PM) equation [31]. The variability of daily rainfall and ET o during the soybean crop seasons is given in Figure A1.
The soil in the experimental fields is a Eutric Cambisol, loamy in the top layer and clay loamy underneath. The total available water (TAW), which represents the difference between the water storage in the root zone at field capacity (33 kPa) and permanent wilting point (1500 kPa), is 176 mm and 144 mm for soils 1 and 2, respectively. The respective main soil hydraulic properties are presented in Table A1. The soil water content (SWC) was measured with a calibrated neutron probe (503DR HYDROPROBE, InstroTek Inc., Martinez, CA, USA). Measurements were performed every 0.10 m until a maximum depth of 1.00 m. Soil sampling was used for the upper 0.10 m layer. Plots were cropped with the soybean variety "Don Mario 5.1i RR" (maturity group V) that is of indeterminate growth and has high yield potential. Each plot was 5 m × 2 m, with five rows spaced 0.4 m. The plant density was 30 plants m −2 . Cropping practices were those recommended locally by the extension services. The irrigation system consisted of pressure compensating in-line drippers spaced 0.20 m and discharging 1.5 L·h −1 . Irrigations were scheduled by performing a simple daily soil water balance applied to a depth of 1.0 m using the computed ET o and the measured SWC data. The irrigation trigger was a depletion of 60% of TAW during periods when water stress was induced, and a depletion of 40% of TAW otherwise. Irrigation depths were set to refill SWC up to 90% of θ FC in the periods when water stress was not allowed and up to 60% of θ FC otherwise.
The following treatments were adopted: (a) FI, full irrigation, aimed at fully satisfying crop water requirements, thus to avoid water stress in all crop growth stages; (b) DI GFill , deficit irrigation during the flowering to grain filling periods; (c) DI Veg , deficit irrigation during the vegetative period; (d) DI Veg-GFill , deficit irrigation during the vegetative to the grain filling periods; (e) Rain-fed.
Water deficits were induced by withholding irrigation or precipitation using rain shelters to allow for water deficits to be induced at desired timings in the crop season. Three replications of the referred five irrigation treatments were adopted. Completely random blocks were used. To assure good crop establishment, no stress was allowed during emergence. The irrigation depths applied during both crop seasons and all irrigation treatments are presented in Table A2.
The dates of each crop growth stage as defined in FAO56 [31] and the respective cumulated growing degree days (CGDD, • C) are presented in Table A3. Measurements of the photosynthetically active radiation (PAR) were performed in the treatment FI using a ceptometer (Decagon AccuPar LP 80). Following Farahani et al. [32], these measurements were converted into canopy cover (CC) and fraction of ground cover (f c ) for use with AquaCrop and SIMDualKc, respectively. The crop height (h, m) and rooting depths (Z r , m) were randomly measured, and the maximum root depth observed was 1 m. The final above ground biomass and soybean grain yield were obtained from harvesting all experimental plots, thus, three samples per irrigation treatment were used; samples were oven dried to a constant weight at 65 ± 5 • C.

Modelling
Two modelling approaches were used: (a) the SIMDualKc [33] soil water balance model that uses the FAO56 dual crop coefficient approach for partitioning crop ET and was combined with the modified Stewart's global water-yield model [17] for yield predictions; and (b) the crop growth and yield model AquaCrop, that partitions ET based upon the canopy cover (CC).
As revised previously [34,35], the FAO56 dual crop coefficient approach (dual-K c , [31,36]) accurately models and partitions ET as described in several studies (e.g., [37,38]) and when compared with the dual-source Shuttleworth-Wallace model [39]. The SIMDualKc model has been positively tested for actual transpiration using sap-flow measurements [40,41] and for soil evaporation using micro-lysimeters [42,43] including soybeans [26]. The SIMDualKc model computes crop evapotranspiration (ET c ) under standard/potential, non-limiting conditions as where ET o (mm) is the reference evapotranspiration, K cb (dimensionless) is the potential basal crop coefficient that describes transpiration (T c ), and K e (dimensionless) is the soil water evaporation coefficient that describes soil evaporation (E s ). The model provides for separately computing potential transpiration T c = K cb ET o (mm) and soil evaporation E s = K e ET o (mm). The actual crop ET (ET c act , mm) is computed by the model as a function of the available soil water in the root zone (ASW): when soil water extraction is smaller than the depletion fraction for no stress (p) then ET c act = ET c , otherwise ET c act < ET c and decreases with decreasing ASW. The ET c act and the T c act are, therefore, defined as where K s (dimensionless) is the water stress coefficient (0-1). K s is computed through a soil water balance applied to the entire root zone (SWB). Soil evaporation is given as with K e depending on the fraction of ground cover by vegetation (f c ) and the SWC in the soil layer with depth Z e of 10-15 cm. K e is computed daily through an SWB of the evaporation layer, which is characterized by the readily and total evaporable water (REW, TEW, mm); REW and TEW may be computed from the soil textural and water holding characteristics of the top-layer [31,36]. K e is adjusted for mulches and for the fraction of soil wetted by irrigation and exposed to radiation. The SWB of the root zone is performed by computing the soil water depletion D r,i at the end of every day i: where the depletion D r,i−1 of the precedent day is i − 1 and the precipitation P, runoff RO, net irrigation depth I, capillary rise CR, deep percolation DP, and crop ET c act are in mm and refer to day i. CR was not considered because the water table was deep. RO was computed using the curve number (CN) approach [44]. DP was computed with a parametric equation [45] requiring two parameters, a D characterizing storage and b D referring to the velocity of vertical drainage, both estimated from soil physical characteristics [45]. The SIMDualKc model calibration consists of searching the model crop parameters-basal crop coefficients K cb and depletion fraction for no stress p, soil evaporation parameters Z e , TEW and REW, runoff curve number CN, and DP parameters a D and b D -that minimize the deviations between the simulated and observed SWC values. The calibration is performed through an iterative procedure of searching the best parameter values within a reasonable range until SWC errors stabilize, as discussed by Pereira et al. [21]. This procedure is first applied to the crop parameters and, after, to the remaining parameters and, finally, to all parameters together. Validation consists of testing the model using the calibrated set of parameters with one or more sets of independent field data collected in the same or different years. However, if validation is performed in a soil with different characteristics, then parameters Z e , TEW, REW, a D, and b D have to be adjusted as described by Giménez et al. [46] for maize in Paysandú. Model calibration was performed using SWC values observed in the FI treatment in 2011-2012. The validation used all other datasets of 2011-2012 and 2012-2013.
As stated above, the SIMDualKc model was combined with a modified version of the water-yield model proposed by Stewart et al. [17] to assess the impacts of water deficits on yields. The version used in the present study assumes a linear variation of the relative yield loss with the relative crop transpiration deficit [19]: where Y a and Y m are the actual and maximum yields (kg·ha −1 ) corresponding, respectively, to the seasonal T c act and T c (mm), and K y is the water-yield response factor. The Y a values consist of observed dry grain. Values for Y m were obtained from maximum yields observed, further using the Wageningen method [18] and checking results against maximal yields achieved by best farmers. The resulting Y m are 6.15 and 5.22 t·ha −1 , respectively, for 2011-2012 and 2012-2013. The value K y = 1.25 was adopted from solving Equation (6) relative to K y using all experimental data available. After knowing K y and Y m , yield predictions were performed by solving Equation (6) in relation to Y a for all T c act results of SIMDualKc. The AquaCrop model is a crop growth and yield model used for a variety of field crops, including soybean, mainly aiming at yield prediction. The model is described by Raes et al. [14] and Vanuytrecht et al. [47], and its open source is described by Foster et al. [48], as well as in various papers quoted there. T c is computed as where CC* is the crop canopy cover (%) adjusted for micro-advective effects, and K cTr,x is the maximum standard crop transpiration coefficient (dimensionless) that corresponds to the K cb mid parameter in FAO56 [31]. T c act is obtained by adjusting T c using the water stress coefficient K s (0-1), as T c act = K s T c (8) K s in AquaCrop is, however, more complex than in FAO56 because it describes the effects of the soil water stress on various processes and the depletion fractions p are inputs of the model that, contrary to SIMDualKc, do not require calibration [14].
Soil evaporation is also obtained from CC* as where K ex is the maximum soil evaporation coefficient (non-dimensional) and K r is the evaporation reduction coefficient (0-1), with K r < 1 when insufficient water is available in the top soil to respond to the evaporative demand of the atmosphere [14]. The product K r (1 − CC*) K ex corresponds to K e as defined in FAO56 as described above. The canopy cover (CC) is similar to f c in FAO56 but while SIMDualKc uses observed f c for adjusting K e , in AquaCrop the CC observations are used to parameterize a CC* curve which is performed in three phases and focuses on four parameters that describe the curve: canopy cover at 90% emergence (cc o ), maximum canopy cover (CC x ), canopy growth coefficient (CGC), and canopy decline coefficient (CDC) [14]. The above ground dry biomass (B, t·ha −1 ) is estimated by the model using the water transpired by the crop throughout the season and the normalized biomass water productivity (BWP*, g·m −2 ). BWP* represents B produced per unit of area considering the cumulative transpiration and ET o [14]. The crop yield (Y, t·ha −1 ) is computed from B as Y = f HI HI o B (10) where HI o is the reference harvest index, describing the harvestable proportion of B, and f HI is an adjustment factor integrating five water stress factors [14]. The model parameterization was initialized using the parameter values proposed by Raes et al. [14]. The parameterizations of the CC curves were first performed using a trial and error procedure. Once these curves were properly parameterized, the trial and error procedure was applied to search the K cTr,x value that leads to a better fit of SWC. In this search, the CN and REW values found for SIMDualKc were used. Growth and yield parameters of AquaCrop were obtained using the above-ground biomass observations. The parameters retained after parameterization using FI data of 2011-2012 were used for model testing using all data sets.
"Goodness-of-fit" indicators were used to assess the accuracy of model simulations at calibration and validation of SIMDualKc and parameterization and testing of AquaCrop. Indicators, following Legates and McCabe Jr. [49], Moriasi et al. [50], and described by Pereira et al. [21], were computed from the pairs of observed and predicted values, respectively, O i and P i (i = 1, 2, ..., n) with means O and P. The regression coefficient b 0 of a regression forced to the origin relating O i and P i was used to verify the similarity between the simulated and observed values. The determination coefficient R 2 of the ordinary least squares regression of the same variables was used to assess the dispersion of pairs of O i and P i values along the regression line, with large R 2 indicating that a large fraction of the variance of observations was explained by the model. The root mean square error (RMSE) and the normalized root mean square error relative to the mean of observations (NRMSE) were adopted to assess modelling errors. In addition, the Nash and Sutcliff [51] modelling efficiency (EF) was adopted to express the relative magnitude of the mean square error (MSE = RMSE 2 ) in relation to the variance of the observed data [49].

Soil Water Simulation and Models Calibration and Parameterization
Simulations with both models are presented in Figure 1: Figure 1a,b refer to DI Veg-GFill in 2011-2012, when a severe water deficit was applied from the vegetative growth to grain filling, a sensitive period to water stress; Figure 1c,d refer to FI in 2012-2013, where water stress was avoided; and Figure 1e,f are relative to rain-fed cropping in 2012-2013, where only a limited stress occurred during pod formation. Water stress for the rain-fed crop is smaller than for that of the DI Veg-GFill treatment because, contrarily to the latter, rainfall was not avoided during any period. Results show that both models behaved well and in a similar way, which is due to their careful calibration/parameterization.
The "goodness-of-fit" indicators relative to all SWC simulations with SIMDualKc and AquaCrop (Table 1) show a better model performance when SIMDualKc is used. Regression coefficients (b 0 ) ranged from 0.95 to 1.01 and R 2 varied from 0.65 to 0.94 for SIMDualKc applications indicating that the predicted and observed values were statistically similar and a large fraction of the total variance of the observed SWC values was explained by the model. Wider but acceptable values were obtained for AquaCrop, with b 0 ranging from 0.92 to 1.06 and R 2 varying from 0.61 to 0.92. The estimation errors were small for SIMDualKc (RMSE < 0.025 cm 3 ·cm −3 and NRMSE < 7.6%) and slightly larger for AquaCrop (RMSE < 0.029 cm 3 ·cm −3 and NRMSE < 9.2%). Model efficiency was high for SIMDualKc, with EF ranging from 0.61 to 0.91, indicating that simulation errors MSE were much smaller than the variance of SWC observations. In contrast, EF values obtained for AquaCrop showed a wider range of variation, 0.16 to 0.93, indicating that MSE varied widely relative to the variance of observations. Overall, results indicate that though both models are appropriate for simulating daily SWC, SIMDualKc performed better.
The SIMDualKc calibrated parameters-K cb , p, TEW, REW, Z e , CN, a D, and b D -are presented in Table 2. CN, Z e , REW, TEW, a D, and b D are the same as those previously obtained by Giménez et al. [46] for the same experimental area because they essentially depend upon the soil characteristics rather than the crop. The K cb and p values are equal to those proposed by Allen et al. [31], K cb ini = 0.15, K cb mid = 1.10 and K cb end = 0.33. Slightly lower K cb mid values were obtained by Odhiambo and Irmak [52] and Wei et al. [26]. K cb ini and K cb end reported by those authors are about the same as for the current study. Calibrated K cb mid values are also coherent to the single crop coefficients K c mid reported by Karam et al. [3], Tabrizi et al. [53] and Payero and Irmak [54]. Thus, results relative to potential K cb and p values confirm those proposed in FAO56 [31].
Relative to AquaCrop, the "goodness-of-fit" of CC curves for FI in both seasons have shown a slight under-estimation trend, with b 0 = 0.93, but other goodness-of-fit indicators were generally high, with an R 2 of 0.99 and RMSE of 6.8% and 6.4% for 2011-2012 and 2012-2013 seasons, respectively. These values are within the range of other AquaCrop applications to soybean [15,16,55]; thus, one may consider the parameterization of the CC curves in the current study adequate. and (e,f) rain-fed in 2012-2013 (error bars indicate the standard deviation of SWC observations; θSat, θFC, θWP, and θp are, respectively, the SWC at saturation, field capacity, wilting point, and at the threshold depletion for no stress).
The conservative and non-conservative parameters used in AquaCrop are also presented in Table 2. The value for KcTr,x = 1.10 equals the Kcb mid calibrated with SIMDualKc (Table 2). Similar values were reported by Abi Saab et al. [15] and Paredes et al. [16]. BWP* = 14 g·m −2 equals the one reported by Abi Saab et al. [15] and Khoshravesh et al. [55]; a higher value was reported by Paredes et al. [16]. For no-stress conditions, HIo observed in both seasons averaged 0.38. That HIo value equals that reported by Paredes et al. [16]; slightly smaller values were reported by Andrade [2] and larger values by Abi Saab et al. [15] and Khoshravesh et al. [55]. Differences in BWP* and HIo values may relate to soybean varieties. Results analyzed show that parameters in Table 2 ) daily average soil water content (SWC) in the soil root zone using the models SIMDualKc (left) and AquaCrop (right) for: (a,b) deficit irrigation during the vegetative to the grain filling periods (DI Veg-GFill ) in 2011-2012; (c,d) full irrigation (FI) in 2012-2013; and (e,f) rain-fed in 2012-2013 (error bars indicate the standard deviation of SWC observations; θ Sat , θ FC , θ WP , and θ p are, respectively, the SWC at saturation, field capacity, wilting point, and at the threshold depletion for no stress).
The conservative and non-conservative parameters used in AquaCrop are also presented in Table 2. The value for K cTr,x = 1.10 equals the K cb mid calibrated with SIMDualKc (Table 2). Similar values were reported by Abi Saab et al. [15] and Paredes et al. [16]. BWP* = 14 g·m −2 equals the one reported by Abi Saab et al. [15] and Khoshravesh et al. [55]; a higher value was reported by Paredes et al. [16]. For no-stress conditions, HI o observed in both seasons averaged 0.38. That HI o value equals that reported by Paredes et al. [16]; slightly smaller values were reported by Andrade [2] and larger values by Abi Saab et al. [15] and Khoshravesh et al. [55]. Differences in BWP* and HI o values may relate to soybean varieties. Results analyzed show that parameters in Table 2 are appropriate for use in Uruguay.  Notes: K cb ini , K cb mid and K cb end are respectively the basal crop coefficients for the initial, mid and end-season stages; p ini, p dev, p mid, and p end are the depletion fractions for no stress for the initial, crop development, mid and end-seasons stages; REW and TEW are the readily and total evaporable water; Z e is the depth of the soil evaporation layer; CN is the curve number; a D and b D are the parameters of the deep percolation equation [46]. * different values were obtained due to the spatial heterogeneity of the soil among experimental plots.

Water Balance and Water Use Components
The actual ET c computed with SIMDualKc and AquaCrop were quite similar (Table 3), which agree with the results of SWC simulation discussed above. However, its partition on T c act and E s produced different values, with AquaCrop generally giving a smaller T c act and a larger E s . Comparing Equations (3)- (7) and Equations (4)- (9), it is apparent that differences mainly stem from procedures used to compute the actual K cb and K e . In fact, the daily K cb and K cb act curves obtained with SIMDualKc and AquaCrop are quite different (Figure 2), particularly under severe water stress (Figure 2a,b). Differences largely stem from the form of the potential K cb curve, with SIMDualKc using the typical linear variation of K cb for the four crop growth stages adopted in FAO56 [31], i.e., a K cb curve defined with only three values-K cb ini , K cb mid , and K cb end - (Figure 2a,c,e), while a curvilinear variation of K cb dictated by the parameterized CC curve is adopted in AquaCrop (Figure 2b,d,f). Without a very severe stress, the variation of K cb are somewhat similar for both models (Figure 2c,d, and Figure 2e,f) but when a severe water stress occurs, e.g., DI Veg-GFill in 2011-2012 (Figure 2b), the K cb curve of AquaCrop is far from representing the potential K cb defined in FAO56 [31,35] because this model does not use K cb mid but just the maximum K cTr,x . When water stress occurs but it does not affects crop development noticeably, as is the case of the rain-fed treatment, both Figure 2e,f show a similar behavior of K cb act until the end of February, but not afterwards, likely due to the model approach used to compute the stress coefficient K s in AquaCrop which includes various stresses in addition to soil water deficits.

Water Balance and Water Use Components
The actual ETc computed with SIMDualKc and AquaCrop were quite similar (Table 3), which agree with the results of SWC simulation discussed above. However, its partition on Tc act and Es produced different values, with AquaCrop generally giving a smaller Tc act and a larger Es. Comparing Equations (3)- (7) and Equations (4)-(9), it is apparent that differences mainly stem from procedures used to compute the actual Kcb and Ke. In fact, the daily Kcb and Kcb act curves obtained with SIMDualKc and AquaCrop are quite different (Figure 2), particularly under severe water stress (Figure 2a,b). Differences largely stem from the form of the potential Kcb curve, with SIMDualKc using the typical linear variation of Kcb for the four crop growth stages adopted in FAO56 [31], i.e., a Kcb curve defined with only three values-Kcb ini, Kcb mid, and Kcb end- (Figure 2a,c,e), while a curvilinear variation of Kcb dictated by the parameterized CC curve is adopted in AquaCrop (Figure 2b,d,f). Without a very severe stress, the variation of Kcb are somewhat similar for both models (Figure 2c    ) and irrigation ( Water 2017, 9, 393 9 of 18

Water Balance and Water Use Components
The actual ETc computed with SIMDualKc and AquaCrop were quite similar (Table 3), which agree with the results of SWC simulation discussed above. However, its partition on Tc act and Es produced different values, with AquaCrop generally giving a smaller Tc act and a larger Es. Comparing Equations (3)-(7) and Equations (4)-(9), it is apparent that differences mainly stem from procedures used to compute the actual Kcb and Ke. In fact, the daily Kcb and Kcb act curves obtained with SIMDualKc and AquaCrop are quite different (Figure 2), particularly under severe water stress (Figure 2a,b). Differences largely stem from the form of the potential Kcb curve, with SIMDualKc using the typical linear variation of Kcb for the four crop growth stages adopted in FAO56 [31], i.e., a Kcb curve defined with only three values-Kcb ini, Kcb mid, and Kcb end- (Figure 2a,c,e), while a curvilinear variation of Kcb dictated by the parameterized CC curve is adopted in AquaCrop (Figure 2b,d,f). Without a very severe stress, the variation of Kcb are somewhat similar for both models (Figure 2c    ) are also depicted. Figure 2 shows a similar behavior during the initial and early vegetative crop stages, though AquaCrop shows a tendency to estimate a larger K e . In contrast, K e tends to be larger afterwards when water stress occurs particularly during mid-season. Differences in K e computed by both models increase when water deficits are larger (Figure 2a,b). Differences between models are due to the fact that K e in SIMDualKc varies with the observed f c and the daily computed depletion of the soil evaporation layer [33], while K e in AquaCrop depends upon the fitted CC curve. Therefore, K e tends to be higher with AquaCrop when water deficits occur. Under no stress conditions, differences are negligible (Figure 2c,d) as observed by Paredes et al. [16].

The daily variation of the K e in all examples of
Analyzing the ET estimates and partition into E s and T c act during the initial period (Table 3) it was observed that E s simulated by SIMDualKc represented 81% to 85% of ET c act while AquaCrop simulated a larger E s corresponding to 92% to 97% of ET c act , thus resulting in a small T c act during this period.  Throughout the crop development stage, E s progressively decreased, as shown in Figure 2, due to the progressive decrease of the soil surface fraction exposed to solar radiation. During this period, E s /ET c act falls, in average, to 29% and 33% when computed with SIMDualKc and AquaCrop, respectively. During the mid-season, the soil is nearly fully shadowed by the crop and the energy available for evaporation drops to minimum values. Thus, E s /ET c act falls to 2% and 6% on average when computed with SIMDualKc and AquaCrop, respectively (Table 3). However, estimated E s /ET c act using AquaCrop had a very wide range, from 1% to 26%, likely due to the heavy dependency of E s on CC (Equation (7)), i.e., whenever the model simulated high impacts of water stress and reduced CC, as for DI Veg-GFill and the rain-fed treatments during 2011-2012, higher E s /ET c act values were estimated. Thus, differences between models relative to T c act amounted to up to 40% when water stress occurred, with higher T c act values being estimated by SIMDualKc (Table 3). During the late season, despite lower coverage of the soil due to leaf senescence, because watering events were small and infrequent, E s /ET c act increased slightly with SIMDualKc but to a higher value averaging 15% when using AquaCrop. Farahani et al. [32] also reported high E s /ET c act ratios with AquaCrop under water stress. Consequently, it could be concluded that AquaCrop tends to underestimate T c act throughout the crop season, mainly under water deficit conditions. Differences relative to the non-consumptive use terms, runoff, and deep percolation are notable, particularly for the 2012-2013 season ( Table 4). RO and DP, whose sum equals the difference between the water input and ET c act , differ between models, with differences stemming from computational approaches as also observed by Pereira et al. [21]. The CN value used for RO computations was the same with both models but related computational processes are different [14,33], thus RO values are also different.  SIMDualKc  821  354  16  266  207  718  606  112  16  AquaCrop  821  354  23  245  266  688  574  114  17  DI GFill  SIMDualKc  676  288  37  190  132  679  577  102  15  AquaCrop  676  288  −8  109  159  688  571  117  17  DI Veg  SIMDualKc  773  162  15  211  204  535  450  85  16  AquaCrop  773  162  19  248  259  448  275  173  39  DI Veg-GFill  SIMDualKc  628  90  22  100  129  511  429  82  16  AquaCrop  628  90  6  80  138  507  370  137  27  Rain-fed  SIMDualKc  821  0  17  98  205  535  457  78  15  AquaCrop  821  0  5  89  217  520  384  136  26   2012-2013  FI  SIMDualKc  786  342  10  408  164  566  447  119  21  AquaCrop  786  342  20  306  273  569  454  115  20  DI GFill  SIMDualKc  666  306  39  304  158  549  434  115  21  AquaCrop  666  306  56  216  274  539  424  115  21  DI Veg  SIMDualKc  746  216  45  330  164  512  407  105  21  AquaCrop  746  216  −2  159  274  527  403  124  24  DI Veg-GFill  SIMDualKc  668  306  39  318  149  546  437  109  20  AquaCrop  668  306  30  180  257  567  451  116  20  Rain-fed  SIMDualKc  786  0  61  183  164  501  400  101  20  AquaCrop  786  0  75  139  219  504  383  121  24 Notes: P is precipitation, I is net irrigation depths, ∆SWC is variation in stored soil water, DP is deep percolation, RO is runoff, ET c act is actual crop evapotranspiration, T c act is the actual crop transpiration, E s is the soil evaporation.
DP values computed with SIMDualKc were generally higher, up to 171 mm, than those estimated with AquaCrop (Table 4) due to differences in the computation of DP: SIMDualKc uses a parametric function (Liu et al. 2006) whose parameters a D and b D are calibrated, as per this application. In contrast, in AquaCrop, DP is estimated using a quasi-deterministic redistribution and drainage module based on the hydraulic characteristics of the soil [14] but does not use calibrated parameters. Possible deficiencies in that DP module were referred by Pereira et al. [21] and Iqbal et al. [56], as well as Farahani et al. [32] who compared computed with field observed DP. As analyzed by several authors (e.g., [57]), AquaCrop had not been tested for severe water stress conditions yet. Results herein relative to the soil water balance components and the insufficiencies in partition of ET c act support the need for improving that model.

Yield Predictions
The available data on water use and transpiration, biomass, and yield covering four seasons, 2009-2010 to 2012-2013, were used to test the biomass and yield predictions by AquaCrop and the Stewart's model combined with SIMDualKc (Stew-SIM). Yields of all treatments were significantly different as per an application of ANOVA (data not shown).
Yield predictions often show better results with the Stew-SIM combined approach relative to AquaCrop ( Table 5). The Stew-SIM approach shows a tendency for slightly over-predicting yields (b 0 = 1.04), with relative deviations between predicted and simulated yields ranging from 1% to 66% (Table 5). AquaCrop results do not show any tendency for under-or over-estimation (b 0 = 0.99) but deviations vary in a wider range, from 1% to 103%. Deviations between observed and simulated yields using the Stew-SIM approach are in the range of those reported by Ma et al. [6] using the CROPGRO-Soybean and the hybrid RZWQM-CROPGRO-Soybean model, and by Banterng et al. [58] when using the CROPGRO-Soybean model. Note: * dried at 65 ± 5 • C; the standard deviation is presented between brackets.
The "goodness-of-fit" indicators relative to all yield predictions with AquaCrop were poor, with RMSE = 1.01 t ha −1 , NRMSE = 22.8%, and EF = −0.41. The negative EF indicates that the MSE is larger than the variance of observations, thus, modelling predictions are poor and there is no effective advantage in using this model. Nevertheless, results in the current study relative to AquaCrop applications are in the range of those reported by Mercau et al. [7] using CROPGRO-Soybean and Cera et al. [12] using SoySim. However, better results using AquaCrop for soybean were reported by Abi Saab et al. [15], Paredes et al. [16], and Battisti et al. [27] whose studies only considered small water stress levels. The above referred results are likely due to the previously discussed poor estimation of actual transpiration when water stress occurs. Katerji et al. [57] and Pereira et al. [21] also reported that AquaCrop biomass and yield predictions were poor under severe water stress because they were hampered by poor estimations of T c act . Good predictions were, however, obtained with AquaCrop for vining pea [59], which was cultivated without water stress, thus confirming that the use of AquaCrop predictions is only recommended when severe water stress is not considered.
The "goodness-of-fit" indicators relative to yield predictions with the Stew-SIM combined approach were RMSE = 0.65 t·ha −1 , NRMSE = 14.5%, and EF = 0.43, which are much better than the indicators relative to the AquaCrop predictions. These RMSE and NRMSE values are in the range of those reported for other model applications, e.g., with the CROPGRO-Soybean model [8]. However, much lower NRMSE were reported when using DSSAT CSM CROPGRO-Soybean [5] and with the hybrid RZWQM-CROPGRO model for soybean [6]. Lower RMSE values were also reported by Setiyono et al. [9] when using the SoySim model in a comparative study using the models CROPGRO-Soybean, Sinclair-Soybean, and WOFOST. Results for these models [9] resulted in a much higher RMSE than the one obtained with the combined Stew-SIM approach. Therefore, the latter is adequate to predict yields aimed at assessing impacts of alternative irrigation scheduling strategies even when a severe stress is considered.

Conclusions
Experimental results relative to various deficit irrigation scheduling treatments confirm that the crop growth stage from flowering to grain filling is the most sensitive to water stress. However, the highest impacts of water stress were observed when deficits were imposed from the vegetative to the grain filling period.
Both SIMDualKc and AquaCrop models were successfully calibrated and validated for soybean using SWC data relative to all treatments and two soybean seasons. The accuracy for simulating the SWC dynamics along the crop seasons was better for SIMDualKc and lower for AquaCrop mainly for the treatments subjected to severe water stress. The water balance terms resulting from the application of both models were quite different, mainly due to different procedures for computing the daily actual basal crop coefficient and the evaporation coefficient, resulting in different values of T c act and E s . Computations of potential and actual K cb in SIMDualKc follow the well-established FAO56 dual-K c methodology while maximum and actual K cb values in AquaCrop depend heavily on the fitted CC curve which only works well for non-stressed crops. Relative to E s , there are large computational differences, also due to the strong dependency of K e on the fitted CC curve in AquaCrop, while K e in SIMDualKc is obtained after calibration of the parameters characterizing the evaporative top soil layer and considering the observed f c fraction.
Differences between models are quite evident in terms of non-consumptive water use, RO and DP. Differences in RO, computed with the same CN, resulted from differences in the algorithms used for the calculations by the models. Relative to DP, the computation modules are very different: in AquaCrop a quasi-deterministic module is used but without calibration; on the contrary, a parametric function is used in SIMDualKc but after calibration of its parameters.
It can be concluded that the calibrated parameters of both SIMDualKc and AquaCrop may be further used for soybean in this region and that SIMDualKc performed more accurately in computing the soil water balance, mainly in estimating T c act , thus proving to be more appropriate to support advising farmers on supplemental irrigation scheduling. Both the AquaCrop model and the SIMDualKc-Stewart's combined approach may be used for yield predictions. However, AquaCrop responded poorly when severe water stress was imposed, which relates to the above referred poor estimation of T c act under those conditions. Thus, whenever the model fitting of CC is worse, the model poorly estimates T c act and, as a consequence, biomass and yield are under-estimated. Results herein clearly identified the main weaknesses of AquaCrop, thus, the need for its further improvement for high water deficit situations. Contrarily, yield predictions with the Stew-SIM approach were good because T c act was predicted accurately and the empirical yield response factor K y was calibrated. Thus, that simple approach can be further used when devising supplemental irrigation strategies for soybean in Uruguay.
Based on the current study, the next step is to design supplemental irrigation strategies to cope with climate variability in line with previous studies [13,60] and to consider water productivity and economic farmers' returns. Further research should also assess the usability of weather forecasts for supporting real time irrigation scheduling. thanks the postdoctoral fellowship (SFRH/BPD/102478/2014) provided by the Fundação para a Ciência e a Tecnologia (FCT). The support of FCT through the research unit LEAF-Linking Landscape, Environment, Agriculture and Food (UID/AGR/04129/2013) is also acknowledged.
Author Contributions: Luis Giménez contributed for the current study by developing the field experiments and related data analysis. Paula Paredes contributed by supporting the testing of both model applications and manuscript writing. Luis S. Pereira was responsible for the theory regarding water balance, water use, and evapotranspiration and for directing manuscript writing.

Conflicts of Interest:
The authors declare no conflict of interest.  Notes: θsat, θFC,, and θWP are, respectively, the soil water content at saturation, field capacity, and wilting point; Ksat is the saturated hydraulic conductivity.   Notes: θ sat , θ FC, , and θ WP are, respectively, the soil water content at saturation, field capacity, and wilting point; K sat is the saturated hydraulic conductivity.