Feasibility analysis of using inverse modeling for estimating field-scale evapotranspiration in maize and soybean fields from soil water content monitoring networks

In this study, the feasibility of using inverse vadose zone modeling for estimating field-scale actual evapotranspiration (ETa) was explored at a long-term agricultural monitoring site in eastern Nebraska. Data from both point-scale soil water content (SWC) sensors and the areaaverage technique of cosmic-ray neutron probes were evaluated against independent ETa estimates from a co-located eddy covariance tower. While this methodology has been successfully used for estimates of groundwater recharge, it was essential to assess the performance of other components of the water balance such as ETa. In light of recent evaluations of land surface models (LSMs), independent estimates of hydrologic state variables and fluxes are critically needed benchmarks. The results here indicate reasonable estimates of daily and annual ETa from the point sensors, but with highly varied soil hydraulic function parameterizations due to local soil texture variability. The results of multiple soil hydraulic parameterizations leading to equally good ETa estimates is consistent with the hydrological principle of equifinality. While this study focused on one particular site, the framework can be easily applied to other SWC monitoring networks across the globe. The value-added products of groundwater recharge and ETa flux from the SWC monitoring networks will provide additional and more robust benchmarks for the validation of LSM that continues to improve their forecast skill. In addition, the value-added products of groundwater recharge and ETa often have more direct impacts on societal decision-making than SWC alone. Water flux impacts human decision-making from policies on the long-term management of groundwater resources (recharge), to yield forecasts (ETa), and to optimal irrigation scheduling (ETa). Illustrating the societal benefits of SWC monitoring is critical to insure the continued operation and expansion of these public datasets.

Abstract. In this study, the feasibility of using inverse vadose zone modeling for estimating field-scale actual evapotranspiration (ET a ) was explored at a long-term agricultural monitoring site in eastern Nebraska. Data from both point-scale soil water content (SWC) sensors and the areaaverage technique of cosmic-ray neutron probes were evaluated against independent ET a estimates from a co-located eddy covariance tower. While this methodology has been successfully used for estimates of groundwater recharge, it was essential to assess the performance of other components of the water balance such as ET a . In light of recent evaluations of land surface models (LSMs), independent estimates of hydrologic state variables and fluxes are critically needed benchmarks. The results here indicate reasonable estimates of daily and annual ET a from the point sensors, but with highly varied soil hydraulic function parameterizations due to local soil texture variability. The results of multiple soil hydraulic parameterizations leading to equally good ET a estimates is consistent with the hydrological principle of equifinality. While this study focused on one particular site, the framework can be easily applied to other SWC monitoring networks across the globe. The value-added products of groundwater recharge and ET a flux from the SWC monitoring networks will provide additional and more robust benchmarks for the validation of LSM that continues to improve their forecast skill. In addition, the value-added products of groundwater recharge and ET a often have more direct im-pacts on societal decision-making than SWC alone. Water flux impacts human decision-making from policies on the long-term management of groundwater resources (recharge), to yield forecasts (ET a ), and to optimal irrigation scheduling (ET a ). Illustrating the societal benefits of SWC monitoring is critical to insure the continued operation and expansion of these public datasets.

Introduction
Evapotranspiration (ET) is an important component in terrestrial water and surface energy balance. In the United States, ET comprises about 75 % of annual precipitation, while in arid and semiarid regions, ET comprises more than 90 % of annual precipitation (Zhang et al., 2001;Glenn et al., 2007;Wang et al., 2009a). As such, an accurate estimation of ET is critical in order to predict changes in hydrological cycles and improve water resource management (Suyker and Verma, 2008;Anayah and Kaluarachchi, 2014). Given the importance of ET, an array of measurement techniques at different temporal and spatial scales have been developed (see Maidment, 1992;Zhang et al., 2014), including lysimeter, Bowen ratio, eddy covariance (EC), and satellite-based surface energy balance approaches. However, simple, low-cost, and accurate field-scale measurements of actual ET (ET a ) still remain a challenge due to the uncertainties of available estima-tion techniques (Wolf et al., 2008;Li et al., 2009;Senay et al., 2011;Stoy, 2012). For instance, field techniques, such as EC and Bowen ratio, can provide relatively accurate estimation of local ET a , but are often cost prohibitive for widespread use beyond research applications (Baldocchi et al., 2001;Irmak, 2010). By comparison, satellite-based remote sensing techniques are far less costly for widespread spatial coverage (Allen et al., 2007), but are limited by their accuracy, temporal sampling frequency (e.g., Landsat 8 has a 16-day overpass), and technical issues that further limit temporal sampling periods (e.g., cloud coverage during overpass) (Chemin and Alexandridis, 2001;Xie et al., 2008;Li et al., 2009;Kjaersgaard et al., 2012).
As a complement to the above-mentioned techniques, recent studies have used process-based vadose zone models (VZMs) for estimating field-scale ET a with reasonable success, particularly in arid and semiarid areas (Twarakavi et al., 2008;Izadifar and Elshorbagy, 2010;Galleguillos et al., 2011;Wang et al., 2016). Although VZMs are time and cost effective for estimating field-scale ET a , they generally require complex model parameterizations and inputs, some of which are not readily available (e.g., soil hydraulic parameters and plant physiological parameters; see Wang et al., 2016). In order to address the issue of missing soil hydraulic parameters, a common approach is to use pedotransfer functions to convert readily available soil information (texture, bulk density, etc.) to soil hydraulic parameters (Wösten et al., 2001); however, significant uncertainties are usually associated with this method for estimating local-scale water fluxes . In fact, Nearing et al. (2016) identified soil hydraulic property estimation as the largest source of information lost when evaluating different land surface modeling schemes versus a soil moisture benchmark. Poor and uncertain parameterization of soil hydraulic properties is a clear weakness of land surface model (LSM) predictive skill in sensible and latent heat fluxes (Best et al., 2015). This problem will continue to compound with the continuing spatial refinement of hyper-resolution LSM grid cells to less than 1 km (Wood et al., 2011).
In order to address the challenge of field-scale estimation of soil hydraulic properties, here we utilize inverse modeling for estimating soil hydraulic parameters based on field measurements of soil water content (SWC) (see Hopmans and Šimunek, 1999;Ritter et al., 2003). While VZM-based inverse approaches have already been examined for estimating groundwater recharge (e.g., Jiménez-Martínez et al., 2009;Andreasen et al., 2013;Min et al., 2015;Ries et al., 2015;Turkeltaub et al., 2015;Wang et al., 2016), their application for ET a estimation has not been adequately tested. Moreover, we note that simultaneous estimation of SWC states and surface energy fluxes within LSMs is complicated by boundary conditions, model parameterization, and model structure (Nearing et al., 2016). With the incorporation of regional soil datasets in LSMs like Polaris (Chaney et al., 2016), effective strategies for estimating ground truth soil hydraulic proper-ties from existing SWC monitoring networks (e.g., SCAN, CRN, COSMOS, State/National Mesonets; see Xia et al., 2015) will become critical for continuing to improve the predictive skill of LSMs.
The aim of this study is to examine the feasibility of using inverse VZMs for estimating field-scale ET a based on long-term local meteorological and SWC observations for an AmeriFlux (Baldocchi et al., 2001) EC site in eastern Nebraska, USA. We note that while this study focused on one particular study site in eastern Nebraska, the methodology can be easily adapted to a variety of SWC monitoring networks across the globe (Xia et al., 2015), thus providing an extensive set of benchmark data for use in LSMs. The remainder of the paper is organized as follows. In the methods section, we will describe the widely used VZM, Hydrus-1D (Šimunek et al., 2013), used to obtain soil hydraulic parameters. We will assess the feasibility of using both profiles of in situ SWC probes as well as the area-average SWC technique from cosmic-ray neutron probes (CRNPs). In the results section, we will compare simulated ET a resulted from calibrated VZM with independent ET a estimates provided by EC observations. Finally, a sensitivity analysis of key soil and plant parameters will be presented.  (Baldocchi et al., 2001) and has been operating continually since 2001. The regional climate is of a continental semiarid type with a mean annual precipitation of 784 mm yr −1 (according to the AmeriFlux US-Ne3 website). According to the Web Soil Survey data (Soil Survey Staff, 2016, http: //websoilsurvey.nrcs.usda.gov/), the soils at the site are comprised mostly of silt loam and silty clay loam ( Fig. 1b and Table 1). Soybean and maize are rotationally grown at the site under rain-fed conditions, with the growing season beginning in early May and ending in October (Kalfas et al., 2011). Since 2001, crop management practices (i.e., planting density, cultivars, irrigation, and herbicide and pesticide applications) have been applied in accordance with standard best management practices prescribed for production-scale maize systems (Suyker and Verma, 2008). More detailed information about site conditions can be found in Suyker et al. (2004) and Verma et al. (2005).
An EC tower was constructed at the center of the field (Figs. 1 and 2a) and continuously measures water, energy, and CO 2 fluxes (e.g., Baldocchi et al., 1988). At this field, sensors are mounted at 3.0 m above the ground when the canopy is shorter than 1.0 m. At canopy heights greater than 1.0 m, the sensors are then moved to a height of 6.2 m until harvest in order to have sufficient upwind fetch (in all directions) representative of the cropping system being studied (Suyker et al., 2004). In this study, hourly latent heat flux measurements were integrated to daily values and then used for calculating daily EC ET a integrated over the field scale. Detailed information on the EC measurements and calculation procedures for ET a are given in Suyker and Verma (2009). Hourly air temperature, relative humidity, horizontal wind speed, net radiation, and precipitation were also measured at the site. Destructive measurements of leaf area index (LAI) were made every 10 to 14 days during the growing season at the study site . We note that the LAI data were linearly interpolated to provide daily estimates. ThetaProbes (TPs) (Delta-T Devices, Cambridge, UK) were installed at four locations in the study field with measurement depths of 10, 25, 50, and 100 cm at each location to monitor hourly SWC in the root zone (Suyker and Verma, 2008 (Zreda et al., 2008, which are converted into SWC following standard correction procedures and calibration methods (see Zreda et al., 2012). In addition, the changes in aboveground biomass were removed from the CRNP estimates of SWC following Franz et al. (2015). The CRNP measurement depth  at the site varies between 15 and 40 cm, depending on SWC. Note that, for simplicity in this analysis, we assume the CRNP has an effective depth of 20 cm (mean depth of 10 cm) for all observational peri- ods. The areal footprint of the CRNP is ∼ 250 ± 50 m radius circle (see Zreda, 2013 andKöhli et al., 2015 for details). Here, we assume for simplicity that the EC and CRNP footprints are both representative of the areal-average field conditions.

Vadose zone model
The Hydrus-1D model (Šimunek et al., 2013), which is based on the Richards equation, was used to calculate ET a . The setup of the Hydrus-1D model is explained in detail by Jiménez-Martínez et al. (2009), Min et al. (2015, and Wang et al. (2016), and only a brief description of the model setup is provided here. Given the measurement depths of the ThetaProbes, the simulated soil profile length was chosen to be 175 cm with 176 nodes at 1 cm intervals. An atmospheric boundary condition with surface runoff was selected as the upper boundary. This allowed the occurrence of surface runoff when precipitation rates were higher than soil infiltration capacity or if the soil became saturated. According to a nearby USGS monitoring well (Saunders County, NE, USGS 411005096281502, ∼ 2.7 km away), the depth to water tables was greater than 12 m during the study period. Therefore, free drainage was used as the lower boundary condition.
Based on the ASCE Penman-Monteith equation, ET r values can be computed for either grass or alfalfa and then, using crop-specific coefficients, daily potential evapotranspiration (ET p ) can be calculated. Here, daily ET r values were calculated for the tall (0.5 m) ASCE reference (ASCE-EWRI, 2005), and daily potential evapotranspiration (ET p ) was calculated according to FAO 56 (Allen et al., 1998): where K c is a crop-specific coefficient at time t. The estimates of growth stage lengths and K c values for maize and soybean suggested by Allen et al. (1998) and Min et al. (2015) were adopted in this study. In order to partition daily ET p into potential transpiration (T p ) and potential evaporation (E p ) as model inputs, Beer's law (Šimunek et al., 2013) was used as follows: where k[−] is an extinction coefficient with a value set to 0.5 (Wang et al., 2009b) and LAI [L 2 L −2 ] is leaf area index described in the previous section. The root water uptake, S(h), was simulated according to the model of Feddes et al. (1978): where α(h)[−] is the root water uptake water stress response function and varies between 0 and 1 depending on soil matric potentials, and S p is the potential water uptake rate and assumed to be equal to T p . The summation of actual soil evaporation and actual transpiration is ET a .
Since the study site has annual cultivation rotations between soybean and maize, the root growth model from the Hybrid-Maize model (Yang et al., 2004) was used to model the root growth during the growing season: where D (cm) is plant root depth for each growing season day, MRD is the maximum root depth (assumed equal to 150 cm for maize and 120 cm for soybean in this study following Yang et al., 2004), AGDD is the accumulated growing degree days, and GDD Silking is the accumulated GDD at the silking point (e.g., accumulated plant GDD approximately 60-70 days after crop emergence). GDD for each growing season day was calculated as where T max and T min are the maximum and minimum daily temperature ( • C), respectively, and T base is the base temperature set to be 10 • C following McMaster and Wilhelm (1997) and Yang et al. (1997). Finally, the Hoffman and van Genuchten (1983) model was used to calculate root distribution. Further details about the model can be found in Šimunek et al. (2013).

Inverse modeling to estimate soil hydraulic parameters
Inverse modeling was used to estimate soil hydraulic parameters for the van Genuchten-Mualem model (Mualem, 1976;van Genuchten, 1980): where θ [L 3 L −3 ] is volumetric SWC; θ r [L 3 L −3 ] and θ s [L 3 L −3 ] are residual and saturated water content, respec- are unsaturated and saturated hydraulic conductivity, respectively; and S e (=(θ − θ r )/(θ s − θ r )) [−] is saturation degree. With respect to the fitting factors, α [1 L −1 ] is inversely related to air entry pressure, n[−] measures the pore size distribution of a soil with m = 1 − 1/n, and l [−] is a parameter accounting for pore space tortuosity and connectivity. Daily SWC data from the four TP locations and CRNP location were used for the inverse modeling. Based on the measurement depths of the TPs, the simulated soil columns were divided into four layers for TP locations (i.e., 0-15, 15-35, 35-75, and 75-175 cm), which led to a total of 24 hydraulic parameters (θ r ; θ s ; α, n, K s , and l) to be optimized based on observed SWC values. In order to efficiently optimize the parameters, we used the method outlined in Turkeltaub et al. (2015). Since Hydrus-1D is limited to optimizing a maximum of 15 parameters at once and that the SWC of the lower layers changes more slowly and over a smaller range than the upper layers, the van Genuchten parameters of the upper two layers were first optimized, while the parameters of the lower two layers were fixed. Then, the optimized van Genuchten parameters of the upper two layers were kept constant, while the parameters of the lower two layers were optimized. The process was continued until there were no further improvements in the optimized hydraulic parameters or until the changes in the lowest sum of squares were less than 0.1 %. Given the sensitivity of the optimization results to the initial guesses of soil hydraulic parameters in the Hydrus model, soil hydraulic parameters from six soil textures were used as initial inputs for the optimizations at each location (Carsel and Parish, 1988), including sandy clay loam, silty clay loam, loam, silt loam, silt, and clay loam. Based on the length of available SWC data from the TP measurements, the periods of 2007, 2008-2010, and 2011-2012 were used as the spin-up, calibration, and validation periods, respectively. Moreover, to minimize the impacts of freezing conditions on the quality of SWC measurements, data from January to March of each calendar year were removed (based on available soil temperature data) from the optimizations.
In addition to the TP profile observations, we used the CRNP area-average SWC in the inverse procedure to develop an independent set of soil parameters. The CRNP was assumed to provide SWC data with an average effective measurement depth of 20 cm at this study site. The observation point was therefore set at 10 cm. As a first guess and in the absence of other information, soil properties were assumed to be homogeneous throughout the simulated soil column with a length of 175 cm. Because the CRNP was installed in 2011 at the study site, the periods of 2011, 2012-2013, and 2014 were used as spin-up, calibration, and validation periods, respectively, for the optimization procedure.
The lower and upper bounds of each van Genuchten parameter are provided in Table 2. With respect to the goodness-of-fit assessment, root mean square error (RMSE) between simulated and observed SWC was chosen as the objective function to minimize in order to estimate the soil hydraulic parameters. The built-in optimization procedure in Hydrus-1D was used to perform parameter estimation. A sensitivity analysis of the six soil model parameters was performed. In addition, three additional performance criteria, including coefficient of determination (R 2 ), mean average error (MAE), and the Nash-Sutcliffe Efficiency (NSE) were used to further evaluate and validate the selected model behavior: ) 2 (10) where n is the total number of SWC data points; O i and P i are, respectively, the observed and simulated daily SWC on day i; andŌ i is the observed mean value. Based on the best scores (i.e., lowest RMSE values), the best optimized set of soil hydraulic parameters at each location was selected. Using the selected parameters, the Hydrus model was then run in a forward mode in order to estimate ET a between 2007 and 2012. Finally, we note that the years 2004-2006 were used as a model spin-up period for the forward model and evaluation of ET a because of the longer climate record length.
3 Results and discussions

Vadose zone inverse modeling results
The time series of the average SWC from the four TP locations along with 1 standard deviation at each depth are plotted in Fig. 4. Based on the large spatial standard deviation values (Fig. 4), despite the relatively small spatial scale (∼ 65 ha) and uniform cropping at the study site, SWC varies considerably across the site, particularly during the growing season. The comparison between SWC data from the CRNP and spatial average of SWC data at the four TP locations in the study field (i.e., average of 10 and 25 cm depths at TP locations) is presented in Fig. 5. The daily RMSE between the spatial average of the TPs and CRNP data is 0.037 cm 3 cm −3 , which is consistent with other studies that reported similar values in semiarid shrublands , German forests (Bogena et al., 2013;Baatz et al., 2014), montane forests in Utah (Lv et al., 2014), sites across Australia (Hawdon et al., 2014), and a mixed land use agricultural site in Austria . We note that we would expect lower RMSE (∼ < 0.02 cm 3 cm −3 ) with additional point sensors located at shallower depths and in more locations distributed across the study site. Nevertheless, the consistent be-  havior between the spatial mean SWC of TPs and the CRNP allows us to explore spatial variability of soil hydraulic properties within the footprint using inverse modeling. This will be described in the next sections. The study period (2007-2012; Fig. 6) contained significant interannual variability in precipitation. During the spin-up period in 2007, the annual precipitation (942 mm) was higher than the mean annual precipitation (784 mm); 2008 was a wet year (997 mm); 2009-2011 were near-average years (715 mm); and 2012 was a record dry year (427 mm) with widespread drought across the region. Therefore, both wet and dry years were considered in the inverse modeling simulation period.
As an illustration, Fig. 7 shows the daily observed and simulated SWC during the calibration (2008-2010) and validation (2011-2012) periods at the TP 1 location (the simulation results of the other three sites can be found in the Supplement Figs. S1, S2, and S3). The results of objective function criterion (RMSE) and the other three performance criteria (e.g., R 2 , MAE, and NSE) between simulated and observed SWC values at TP locations are presented in Table 3.  In this research, we define RMSE values less than 0.03 cm 3 cm −3 between observed and simulated SWC values as well matched and RMSE between 0.03 and 0.06 cm 3 cm −3 as fairly well matched. We note the target error range of satellite SWC products (e.g., SMOS and SMAP) is less than 0.04 cm 3 cm −3 (Entekhabi et al., 2010). Similar to previous studies (e.g., Jiménez-Martínez et al., 2009;Andreasen et al., 2013;Min et al., 2015;Wang et al., 2016), the results of all the performance criteria at TP locations show the capability of inverse modeling in estimation of soil hydraulic parameters. The results of the calibration period (2008-2010) indicate that the simulated and observed SWC values are in good agreement (i.e., well matched as defined above) throughout the entire period at most locations and depths ( Fig. 7 and Table 3). In addition, the simulated and observed SWC data are fairly well matched at most locations and depths during the validation period (2011-2012), with notable differences during the second half of 2012 during the extreme drought conditions ( Fig. 7 and Table 3). Reasons for this disagreement in the observed and simulated SWC data will be discussed in the following sections.
The results of inverse modeling using the CRNP data also indicate the feasibility of using these data to estimate effective soil hydraulic parameters ( Fig. 8 and Table 4). Based Table 3. Goodness-of-fit measures for simulated and observed SWC data at different depths during the calibration period (2008-2010) and validation period (2011-2012) at TP locations. Note that we assume a good fit as an RMSE between 0 and 0.03 cm 3 cm −3 and fair as between 0.03 and 0.06 cm 3 cm −3 .

Location Depth
Calibration period ( on the performance criteria (Table 4), the simulated data are fairly well matched with the observed SWC data during both the calibration and validation periods. Additional information from deeper soil probes or more complex modeling approaches, such as data assimilation techniques Renzullo et al., 2014), may be needed to fully utilize the CRNP data for the entire growing season. However, this was beyond the scope of the current study and merits further investigation given the global network of CRNP  dating back to ∼ 2011. Table 5 summarizes the optimized van Genuchten parameters for the four different depths of the four TP locations and the single layer for the CRNP location. The optimized parameters were then used to estimate ET a for the entire study period as an independent comparison to the EC ET a data. The results of the ET a evaluation will be discussed in the next section. According to the simulation results (Table 5), in most of the soil layers, the TP 4 location results in lower n and K s values and higher θ r values than the other three locations (TPs 1-3), suggesting either underlying soil texture variability in the field or texture-dependent sensor sensitivity/calibration. As a validation for the simulation results, the publicly available Web Soil Survey data (http: //websoilsurvey.nrcs.usda.gov/) was used to explore whether the optimized van Genuchten parameters from the inverse modeling ( Fig. 1b and Table 2) agreed qualitatively with the survey data. Based on the Web Soil Survey data, the soil at the TP 4 location contains higher clay percentage than the other locations. Meanwhile, the optimized parameters reflect the spatial pattern of soil texture in the field as shown by the Web Soil Survey data (e.g., lower n and K s values and higher θ r values at the TP 4 location with finer soil texture). Physically, finer-textured soils generally have lower K s and higher θ r values (Carsel and Parrish, 1988). Moreover, the shape factor n is indicative of pore size distributions of soils. In general, finer soils with smaller pore sizes tend to have lower n values (Carsel and Parrish, 1988). The observed SWC at the TP 4 location is consistently higher than the average SWC of the other three locations (Fig. S4), which can be partly attributed to the higher θ r values at the TP 4 location . Overall, the obtained van Genuchten parameters from the inverse modeling are in qualitatively good agreement with the available spatial distribution of soil texture in the study field, indicating the capability of using inverse VZM to infer soil hydraulic properties. Further work on validating the Web Soil Survey data soil hydraulic property estimates is of general interest to the LSM community.

Comparison of modeled ET a with observed ET a
Because a longer set of climatic data was available at the study site (as compared to SWC data), we used 2004-2006 as a spin-up period. Using the best-fit soil hydraulic parameters for the four TP locations and the single CRNP location, the Hydrus-1D model was then run in a forward mode to calculate ET a over the entire study period (2007)(2008)(2009)(2010)(2011)(2012). The simulated daily ET a was then compared with the independent EC ET a measurements using RMSE (Eq. 9) as the  evaluation criterion. In order to upscale TP ET a estimation to the field/EC scale, we used the soil textural boundaries and areas defined by the Web Soil Survey data map to compute a weighted average ET a . In this research, we consider RMSE values less than 1 mm day −1 between observed and simulated ET a values as well matched and RMSE values between 1 and 1.2 as fairly well matched ( Fig. 9 and Table 6). The performance criterion results indicate that the simulated daily ET a is in a better agreement with EC ET a measurements at the TP 1-3 locations than at the TP 4 and CRNP locations (Table 6). However, based on the performance criteria from inverse modeling results and on the Web Soil Survey data, we conclude that spatial heterogeneity of soil texture in the study field results in significant spatial variation in ET a rates across the field (e.g., less ET a occurs at the TP 4 location than in the other parts of the field). Here, smaller ET a rates at the TP 4 location are likely due to finer soil texture at this location, which makes it more difficult for the plant/roots to overcome potentials to extract water from the soil, thus leading to a lower ET a rate and greater plant stress. In addition, higher surface runoff can be expected at the TP 4 location due to finer-textured soils (as we observed during our field campaigns). According to the simulation results, the average surface runoff at the TP 4 location was about 44.8 mm yr −1 from 2007 to 2012, while the average surface runoff at the other three locations (TPs 1-3) was around 10.6 mm yr −1 , which partially accounts for the lower ET a rates. We note that future work using historic yield maps may also be used to further elucidate the soil hydraulic property differences given the direct correlation between transpiration and yield. Given that CRNPs have a limited observational depth and that only a single soil layer was optimized in the inverse model for the CRNP, one could expect the simulated daily ET a from the CRNP to have larger uncertainty. Here, we found an RMSE of 1.14 mm day −1 using the CRNP versus 0.91 mm day −1 for the upscaled TP locations. However, when the optimized soil parameters obtained from the CRNP data were used to estimate ET a , the model did simulate daily ET a fairly well during both non-growing and growing seasons in comparison to the EC ET a measurements.
On the annual scale, ET a measured by the EC tower accounted for 87 % of annual P recorded at the site during the study period (Fig. 6). Overall, the simulated annual ET a at all the TP and CRNP locations is comparable to the annual ET a measured by the EC tower, except during 2012 (Table 7), in which a severe drought occurred in the region. One explana- Table 5. Optimized van Genuchten parameters in different locations at the study site. Note that 95 % confidence intervals are in parentheses.

Location
Depth  sensitivity analysis was performed on the maximum rooting depth and presented in the following section. However, we note that given the fact that EC ET a estimation can have up to 20 % uncertainty (Massman andLee, 2002, andHollineger andRichardson, 2005), and accounting for the natural spatial variability of ET a due to soil texture and root depth growth uncertainties, the various ET a estimation techniques performed fairly well. In fact, it is difficult to identify which ET a estimation method is the most accurate method. These results are consistent with the concept of equifinality in hydrologic modeling given the complexity of natural systems (Beven and Freer, 2001). Moreover, the findings here are consistent with Nearing et al. (2016) who showed that information lost in model parameters greatly affects the soil moisture comparisons against a benchmark. However, soil parameterization was less important in the loss of information for the comparisons of ET/latent energy against a benchmark. Fully resolving these issues remains a key challenge to the land surface modeling community and the model's ability to make accurate predictions (Best et al., 2015). The following section provides a detailed sensitivity analysis of the soil hydraulic parameters and root depth growth functions in order to begin to understand the sources of error in estimating ET a from SWC monitoring networks.

Sensitivity analysis of soil hydraulic parameters and rooting depth
In this research, we compared simulated ET a with the measured EC ET a . As expected, some discrepancies between simulated and measured ET a values existed. In order to begin to understand the key sources of error, we performed a set of sensitivity analysis experiments on the estimated soil hydraulic parameters. Building on Wang et al. (2009b), a sensitivity analysis for a single homogeneous soil layer (6 parameters) and a four-layer soil profile (24 parameters) was performed over the study period (2007)(2008)(2009)(2010)(2011)(2012). Here, we performed a preliminary sensitivity analysis by changing a single soil hydraulic parameter one at a time while keeping the other parameters constant (i.e., at the average value). Figure 10 illustrates the sensitivity results on simulated ET a , indicating that the soil hydraulic parameters have a range of sensitivities with tortuosity (l) being the least. We found that n and α were the most sensitive, particularly in the shallowest soil layer. This sensitivity to the shallowest soil layer provides an opportunity to use the CRNP observations, particularly in the early growing season (i.e., when evaporation dominates latent energy flux), to help constrain estimates of n and α. As the crop continues to develop (and transpiration contributes a relatively larger component of latent energy), additional information about deeper soil layers should be used to estimate soil hydraulic parameters or perform data assimilation. Moreover, the CRNP may be useful in helping constrain and parameterize soil hydraulic functions in simpler evaporation models widely used in remote sensing (see Allen et al., 2007) and crop modeling (see Allen et al., 1998).
Following the sensitivity analysis, we repeated the optimization experiment using only α, n, and K s , and used model default estimates for the other parameters in each layer. We found that the RMSE values were significantly higher (1.511 vs. 0.911 mm day −1 ) than when considering all 24 parameters. We suspect, given the high correlation between soil hydraulic parameters (Carsel and Parrish, 1988), that fixing certain parameters leads to a degradation in overall performance. We suggest that further sensitivity analyses, in particular changing multiple parameters simultaneously or using multiple objective functions, be used to fully understand model behavior (see Bastidas et al., 1999 andRosolem et al., 2012).
A sensitivity analysis of ET a by varying rooting depth is summarized in Fig. 11. As would be expected with increasing rooting depth, higher ET a occurred. In addition, Fig. 11 illustrates a decreasing RMSE against EC observations for increases up to 200 %. Again, it is unclear if the EC observations are biased high or if in fact rooting depths are much greater than typically considered in these models. The high observed EC values in the drought year of 2012 indicate that roots likely take up water from below the 1 m observations. Certainly the results shown here further indicate the importance of root water uptake parameters in VZMs and LSMs,  even in homogeneous annual cropping systems. While it is beyond the scope of this paper, we refer the reader to the growing literature on the importance of root water uptake parameters on hydrologic fluxes (see Schymanski et al. 2008 andGuswa, 2012).

Applications and limitations of the vadose zone modeling framework
Given their simplicity and the widespread availability of ground data, ET r and K c values are often used in a wide variety of applications to estimate ET p and thus approximate ET a . It is well known that SWC is a limiting factor affecting the assumption that ET p ∼ ET a . On the other hand, we know that SWC observations are local in nature and not necessarily representative of ET a footprint estimates. The key questions include what the value of SWC observations is, how many profiles should be installed in a footprint, and which Figure 11. Sensitivity analysis of the effect of root depth on ET a estimation for a single homogeneous soil layer profile. Note that root depth is shown in terms of percent depth, as it is dynamic over the growing period.
depths should be used to constrain estimates of fluxes. The well-instrumented and long-term study presented here allows us to start to answer these key questions. First, we find that ET p has an average annual value of 1064.9 mm as compared to EC at 612.5 mm (Table 7). By including individual SWC profiles (TP 1-4) and the CRNP in the VZM framework, we are able to constrain our estimate of ET a to between 525.3 and 643.1 mm and reduce ET a RMSE from 1.992 to around 1 mm day −1 (Table 6). In addition, a range of soil hydraulic parameters for each depth and spatially averaged top layer can be estimated to help better constrain recharge fluxes simultaneously. Given the principle of equifinality in hydrologic systems, the VZM framework may lead to equally reasonable estimates of parameters, which is a limitation of the method and LSMs in general. Based on our sensitivity analysis (Fig. 10), the key parameters of α and n may greatly affect ET a . Although sparsely distributed, widespread state, national, and global meteorological observations paired with SWC profiles (Xia et al., 2015) and the VZM framework provide an opportunity to better constrain ET a and local soil hydraulic functions. Moreover, where multiple SWC profile information is available, a range of ET a and soil hydraulic parameters can be estimated and thus considered in LSM data assimilation frameworks. The combination of basic metrological observations with a CRNP in the VZM framework further allows for estimates of upscaled soil hydraulic parameters with similar estimates of ET a as found with individual SWC profiles. Moving forward, combining CRNP with deeper SWC observations from point sensors seems to be a reasonable strategy in order to average the inherent SWC variability in the near surface yet provide SWC constraints at depth, particularly as annual crops develop over the growing season.

Conclusions
In this study, the feasibility of using inverse vadose zone modeling for field-scale ET a estimation was explored at an agricultural site in eastern Nebraska. Both point SWC sensors (TPs) and area-average techniques (CRNPs) were explored. This methodology has been successfully used for estimates of groundwater recharge but it was critical to assess the performance of other components of the water balance such as ET a . The results indicate reasonable estimates of daily and annual ET a but with varied soil hydraulic function parameterizations. The varied soil hydraulic parameters were expected given the heterogeneity of soil texture at the site and consistent with the principle of equifinality in hydrologic systems. We note that while this study focused on one particular site, the framework can be easily applied to other networks of SWC monitoring across the globe (Xia et al., 2015). The value-added products of groundwater recharge and ET a flux from the SWC monitoring networks will provide additional and more robust benchmarks for the validation of LSM that continues to improve their forecast skill.
Data availability. The climatic and EC data used in this research can be found at http://ameriflux.lbl.gov/. The TP SWC and LAI data in the study site are provided by Andrew Suyker and CRNP SWC are provided by Trenton E. Franz, and both sets of data can be requested directly from the authors. The US soil taxonomy information is provided by Soil Survey Staff and is available online at http://websoilsurvey.nrcs.usda.gov/. The remaining datasets are provided in the Supplement associated with this paper.
The Supplement related to this article is available online at doi:10.5194/hess-21-1263-2017-supplement.