Satellite-Based Evapotranspiration in Hydrological Model Calibration

Abstract

Hydrological models are usually calibrated against observed streamflow (Q_obs), which is not applicable for ungauged river basins. A few studies have exploited remotely sensed evapotranspiration (ET_RS) for model calibration but their effectiveness on streamflow simulation remains uncertain. This paper investigates the use of ET_RS in the hydrological calibration of a widely used land surface model coupled with a source–sink routing scheme and global optimization algorithm for 28 natural river basins. A baseline simulation is a setup based on the latest model developments and inputs. Sensitive parameters are determined for Q_obs and ET_RS-based model calibrations, respectively, through comprehensive sensitivity tests. The ET_RS-based model calibration results in a mean Kling–Gupta Efficiency (KGE) value of 0.54 for streamflow simulation; 61% of the river basins have KGE > 0.5 in the validation period, which is consistent with the calibration period and provides a significant improvement over the baseline. Compared to Q_obs, the ET_RS calibration produces better or similar streamflow simulations in 29% of the basins, while further significant improvements are achieved when either better ET or precipitation observations are used. Furthermore, the model results show better or similar performance in 68% of the basins and outperform the baseline simulations in 90% of the river basins using model parameters from the best ET_RS calibration runs. This study confirms that with reasonable precipitation input, the ET_RS-based spatially distributed calibration can efficiently tune parameters for better ET and streamflow simulations. The application of ET_RS for global scale hydrological model calibration promises even better streamflow accuracy as the satellite-based ET_RS observations continue to improve.

1. Introduction

Physically-based hydrological models are an essential tool for water resources planning, drought, and flood prediction, climate change assessment, and many other precipitation-runoff related applications. To achieve optimized performance, hydrological models are often calibrated against streamflow observations (hereafter termed Q_obs), which measure integrated runoff fluxes [1,2,3]. However, Q_obs may not always be the best applicable option for hydrological model calibration in various circumstances. Firstly, streamflow prediction in ungauged basins often uses the regionalization approach which transfers hydrological parameters from Q_obs-based calibrations developed from donor gauged basins. Such regionalization methods tend to lose efficiency during the transfer process because the assumptions of spatial proximity and physical similarity are not always valid [4,5,6] Secondly, Q_obs-based calibration has long been criticized for its equifinality problem, i.e., different parameter sets can satisfy the objective functions during the calibration period, while the model performance deteriorates significantly in the validation period [7,8]. Thirdly, Q_obs calibration models tend to derive spatially homogeneous parameters across the river basin which neglects characteristic landscape heterogeneity [9]. Furthermore, to accomplish the “goodness-of-fit” between simulated and observed streamflow, Q_obs calibration models are computationally time-consuming in routing spatially distributed runoff through complex drainage networks to produce estimated streamflow at the river basin outlet for each model time step [10].

Research has been conducted to calibrate hydrological models using observations of other water balance variables, such as soil moisture and evapotranspiration (ET) with or without observed streamflow [2]. Satellite-based estimations of these water balance variables have demonstrated value in providing spatially distributed observations with large geographical coverage at regular time intervals and with sufficient record length. For example, several remotely sensed ET records (hereafter referred to as ET_RS) have been developed in recent years which differ in their governing equations, forcing data, and spatial and temporal scales of operation [11,12,13,14,15]. Many studies on the use of ET_RS for hydrological model calibration have also been conducted ([16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34], see

Table 1. Literature reviews on integrating ET_RS in hydrological model calibration. In the “ET_RS*” column, either available ET_RS products or user calculated ET using a specific equation with remote sensing data input are shown. In the “Model Configuration” column, the number of parameters involved in the calibration process is also included in the bracket for each study. The study of Nijzink et al. [27] employed five different conceptual models with multiple parameters; readers are recommended to investigate the original literature if interested; we used (*) to represent this complex case in the “Model Configuration” column.

Reference	River Basin	ETRS*	Model Configuration	Precipitation Input
Immerzeel and Droogers [16]	One regulated basin, India (45,678 km2)	SEBAL(MODIS)	Distributed SWAT model (6 to 53)	Field measurements
Winsemius et al. [17]	One regulated basin, Zambia (15,000 km2)	SEBAL(MODIS)	Distributed HBV model (7)	Satellite estimations
Zhang et al. [18]	120 basins, Australia (50 to 2000 km2)	P–M (MODIS)	Lumped SIMHYD model (7)	Field measurements
Rientjes et al. [19]	Seven basins, Iran (1286 to 9873 km2)	SEBS(MODIS)	Lumped HBV model (7)	Field measurements
Vervoort et al. [20]	Four basins, Australia (147 to 2183 km2)	MOD16 ET	Lumped SMA model (8)	Field measurements
Kunnath-Poovakka et al. [21]	11 basins, Australia (62 to 1028 km2)	CMRSET	Distributed AWRA-L model (5)	Field measurements
López et al. [22]	One basin, Morocco (38,025 km2)	GLEAM ET	Distributed PCR-GLOBWB model (6)	Multi-Source data
Tobin and Bennett [23]	One basin, the United States (610 km2)	GLEAM ET	Distributed SWAT model (16)	Multi-Source data
Demirel et al. [24]	One basin, Denmark (2500 km2)	TSEB(MODIS)	Distributed mHM model (26)	Field measurements
Herman et al. [25]	One basin, the United States (1110 km2)	ALEXI (GOES), SSEBop(MODIS)	Distributed SWAT model (3 to 15)	Field measurements
Kunnath-Poovakka et al. [26]	28 basins, Australia (55 to 193 km2)	CMRSET	Distributed AWRA-L model (5)	Field measurements
Nijzink et al. [27]	27 basins, Europe (144 to 1587 km2)	LSA-SAF ET, MOD16 ET	Lumped (*) models (*)	Multi-Source data
Pan et al. [28]	One basin, China (5996 km2)	SEBAL(MODIS)	Distributed DHSVM model (8)	Field measurements
Poméon et al. [29]	36 basins, Africa (473,718 km2 in total )	GLEAM ET, MOD16 ET	Distributed mHM model (29)	Multi-Source data
Rajib et al. [30]	One basin, the United States (1670 km2)	MOD16 ET	Distributed SWAT model (21 to 31)	Field measurements
Wambura et al. [31]	One basin, Tanzania (41,170 km2)	MOD16 ET	Distributed SWAT model (19)	Field measurements
Becker et al. [32]	One basin, Pakistan (15,000 km2)	SEBAL(MODIS)	Distributed SWAT model (44)	Field measurements
Odusanya et al. [33]	One basin, Nigeria (20,292 km2)	GLEAM ET, MOD16 ET	Distributed SWAT model (11)	Field measurements
Tobin and Bennett [34]	One basin, the United States (3356 km2)	GLEAM ET	Distributed SWAT model (23)	Field measurements
This study	28 basins, the United States (326 to 14,332 km2)	MOD16 ET	Distributed VIC model (4 to 6)	Field measurements

). As one of the pioneering studies in this field, Immerzeel and Droogers [16] optimized the parameters of the Soil and Water Assessment Tool (SWAT) model by tuning the modeled ET against ET_RS for a human-regulated catchment in India. Their results showed close correspondence between the calibrated model ET and ET_RS, although they were not able to provide direct validation of the model streamflow simulations due to heavy human-regulation of the study area. Subsequent studies modified the calibration approach by taking Q_obs or ET_RS as the calibration target independently using a single objective function or calibrating against both “Q_obs + ET_RS” simultaneously using a multi-objective function. Most of these studies concluded that, compared to the Q_obs calibration scheme, the “Q_obs + ET_RS” calibration was able to improve the ET simulation while maintaining suitable model streamflow performance [18,19,20,21,24,25,26,28,29,30,31]. However, compared to the Q_obs calibration model, some recent ET_RS calibration models seem to result in a significant reduction in streamflow simulation accuracy [19,20,24,25], while the improvement of streamflow prediction skill is very low [21,22,26].

As a critical hydrological component, ET represents moisture lost to the atmosphere from plants and soil, and thus significantly impacts the precipitation partitioning into different water balance components, particularly for runoff estimation. Hence, hydrological model performance should benefit from better ET accuracy by calibrating model parameters against high-quality ET observations. Previous studies using the ET_RS calibration method have generally not demonstrated corresponding improvements in streamflow prediction performance. The lack of progress in this area reflects the need for a systematic investigation of the effectiveness of ET_RS calibration and greater attention to the uncertainties in model inputs and model structure. The uncertainties in the ET_RS estimations over space and time will directly influence the efficiency of the ET_RS calibration method. Error in the input data, limitation of the algorithm in representing biophysical processes, and the uncertainty of biophysical and environmental parameters will together contribute to the biases in the ET_RS estimations. When precipitation uncertainty is unknown or significant for hydrological model calibration [17,22,23,27,29,35], the optimization algorithm will be forced to adjust ET to compensate for errors in the precipitation inputs, rather than reducing parametric uncertainty. The lumped or oversimplified representation of hydrological processes and storages may be the main reason for poor streamflow performance [18,19,20,27]. The model structure limitation can cause serious difficulties in the attribution of the model performance to the ET_RS calibration method, particularly in study areas with significant anthropogenic activities [16,17]. Furthermore, most prior studies have focused on only one small catchment, with the representativeness and transferability of the ET_RS calibration results remaining unknown [23,24,25,30].

Given the significant improvement in ET_RS estimation accuracy [36,37] and foreseeably higher-quality global ET products to be available in the near future, the ET_RS calibration method may provide a viable approach for large scale or global hydrological modeling. For example, the need for better calibration is one of the major challenges of large scale hydrological modeling for global flood and drought prediction [38,39,40,41]. Therefore, the purpose of this study is to investigate the utility of remote sensing-based ET estimation in hydrological model calibration, emphasizing the evaluation of ET_RS calibration impacts on streamflow simulation.

The following sections describe the methodology for this study, including the detailed calibration experimental design and a brief introduction of the hydrological model (Section 2); a description of the study area and datasets used (Section 3); a presentation of the major results (Section 4); followed by the conclusions of the implications of this work (Section 5).

2. Methodology

Previous studies involving ET_RS calibration for hydrological modeling have largely adopted a lumped calibration method, whereby the modeled basin average ET is tuned toward spatially averaged ET_RS observations. Relatively few studies have used the spatially distributed ET_RS calibration method with all grid cells or sub-basins being treated individually [30,32,33]. The lumped ET_RS calibration method relies heavily on the seasonal pattern of ET_RS observations. The ET_RS based spatially distributed calibration takes a further step than the lumped calibration by deriving an optimized spatially distributed model parameter set that accounts for the ET_RS spatial patterns [9]. It is technically convenient to implement parallel computing for spatially distributed ET_RS calibration on a per grid cell basis to reduce the computational load. For global applications and ungauged basins, this study tests the efficiency and effectiveness of the spatial ET_RS calibration method based on the large-scale, semi-distributed, Variable Infiltration Capacity (VIC) model [42], coupled with a dominant river-tracing-based source-sink hydrologic routing scheme [43]. A total of 28 natural river basins (with negligible human interference) across a diversity of hydroclimate in the Pacific Northwest region are selected. The model simulations start on 1 January 1999, and end on 31 December 2013. The years 2004–2008 are used as the calibration period, while the 2009–2013 period is used for validation.

2.1. Hydrological Model

The VIC is a semi-distributed, physically-based hydrological model, and has been extensively used in studies involving water resources management, land-atmosphere interactions, and climate change [3,35,40,44,45,46,47,48,49,50]. A key challenge in the application of the VIC model in many large-scale applications has involved the need for automated calibration. Therefore, the VIC model has been tested, improved, and optimized for many regional river basins, resulting in relatively less uncertainty in various model inputs and predictions. Wu et al. [3] developed dominant river-tracing-based streamflow and temperature model, integrated with the VIC model, and successfully produced streamflow and temperature predictions at 1/16 degree resolution for the Pacific Northwest (PNW) domain. Here, we implement a similar modeling framework as Wu et al. [3], including the same VIC model components for ET and runoff simulations, the runoff-routing model, and calibration algorithms. The VIC model balances both the water and energy fluxes at the grid cell scale at each time step. While a sub-daily time step water-and-energy balance simulation is crucial in land-atmosphere studies, satisfying surface energy budgets requires substantial computational time, but has little influence on streamflow or ET simulations [51]. Therefore, we chose the VIC model running in the water balance mode at a daily time step for this study, similar to Wu et al. [3].

Based on Wu et al. [3], we set the VIC model run at 1/16-degree spatial resolution with subgrid spatial variability in vegetation and soil [48]. Vegetation parameters profoundly affect all components of the water cycle in the VIC model. The predefined subgrid vegetation parameters (e.g., land cover types, root depths) are set generally identical to those in previous studies [3,48,52], which were derived from the 1-km-resolution University of Maryland land cover database [53]. The subgrid soil parameters (e.g., bulk density, field capacity) are obtained from previous studies [3,48,52], which were derived from the 1-km-resolution Pennsylvania State University database [54]. The VIC model calibration processes for the streamflow simulation involves six parameters (

Table 2. The six soil parameters calibrated for the Variable Infiltration Capacity (VIC) model simulations in this study.

Parameters	Unit	Description	Range
B	N/A	The variable infiltration curve parameter	0.0–0.4
Dsmax	mm/day	maximum baseflow that can occur in the third soil layer	0.0–30.0
Ds	fraction	The fraction of Dsmax where non-linear baseflow begins	0.0–1.0
Ws	fraction	The fraction of maximum soil moisture where non-linear baseflow occurs	0.0–1.0
D2	M	Depth of second soil layer	0.1–1.5
D3	M	Depth of third soil layer	1.0–3.0
Note: All six parameters (B, Dsmax, Ds, Ws, D2, and D3) will be used in Qobs calibration model, and four parameters (Dsmax, Ds, Ws, and D3) will be used in ETRS calibration model.

). A primary effort of calibrating the VIC model for the major river basins over the Continental United States at the 1/8 degree resolution was made by Maurer et al. [48] by matching the simulated streamflow to available observations. The resulting parameters were used in the study of Wu et al. [3] for the Columbia River Basin. Livneh et al. [52] further refined the parameters B, Dsmax, and D3 at the 1/16-degree resolution using a similar method as Troy et al. [1] from which the resulting parameter set is taken as the baseline soil information input for this study.

2.2. Model Calibration Experiments

A four-step workflow (Figure 1) was developed to perform the ET_RS based spatially distributed calibration in this study and evaluate the impact of the ET_RS calibration on the VIC model streamflow simulations. Major uncertainties for the model streamflow simulations are contributed from the model inputs (e.g., vegetation, soil, and precipitation inputs) and the simplified model structure that generally underrepresents the complexity of human-regulated river basins. Only upstream or headwater river basins without dams are selected for this study to emphasize natural flow conditions while minimizing model uncertainty impacts from river regulation.

2.2.1. Updating Model Vegetation Inputs

We first evaluate and update the VIC model vegetation inputs using a satellite observation-based updated Leaf Area Index (LAI) dataset to define variations in canopy cover influencing the model ET calculations. Two widely used LAI climatology datasets are evaluated by comparing their corresponding streamflow prediction skills. The first LAI product is monthly averaged from a quarter-degree resolution monthly global LAI database [55] based on Advanced Very High-Resolution Radiometer (AVHRR) observations extending from 1981 to 1994. The AVHRR LAI record has been used for delineating vegetation cover in the VIC model in many prior studies [48,52] and is still widely used in VIC and other hydrological modeling applications [56,57,58,59,60,61,62]. The AVHRR LAI dataset is therefore referred to as the Default LAI record for this study. However, recent studies have found that temporally updated LAI inputs, rather than an outdated LAI monthly climatology, can improve VIC model performance [51,63,64,65]. Therefore, we adopted the 8-day, 1-km resolution GLASS/MODIS LAI record [66] as an additional model LAI input in this study, which extends over the 2004 to 2013 period and is monthly averaged as the Updated LAI record. The VIC model was then configured with the up-to-date vegetation information-based on the comparison of the streamflow simulations using the Default LAI and Updated LAI datasets.

2.2.2. Identifying Sensitive Parameters for Model Calibration

In the second step, we identified key model parameters influencing the resulting streamflow and spatial ET simulations. Given the up-to-date vegetation information used in the model, spatially distributed soil parameters are of critical importance. The model sensitivity tests are performed both at the river basin level and on a per grid cell basis, by separately varying each parameter equally over 100 intervals throughout the predefined parameter space (Table 2). This setting is useful for conducting a local sensitivity analysis, providing diagnostic insights on the model’s behavior around a point [67].

The Nash–Sutcliffe coefficient (NSC) [68] is used to evaluate the monthly streamflow simulation. The absolute Percent Bias (PBIAS) [69] is adopted to measure the 8-day spatial ET simulation and abbreviated as APBIAS hereafter.

$$NSC\left(streamflow\right)=1-\frac{{{\displaystyle \sum }}_{i=1}^n{\left(Q_{sim}-Q_{obs}\right)}^2}{{{\displaystyle \sum }}_{i=1}^n{\left(Q_{sim}-{\mu }_Q{}_{obs}\right)}^2},$$

(1)

$$APBIAS\left(spatialET\right)=\left|\frac{{{\displaystyle \sum }}_{i=1}^n\left(E{T}_{sim}-E{T}_{RS}\right)}{{{\displaystyle \sum }}_{i=1}^n\left(E{T}_{RS}\right)}\right|*100,$$

(2)

where $Q_{obs}$ and $Q_{sim}$ are the monthly time series of streamflow observations and simulations respectively, and ${\mu }_Q{}_{obs}$ is the mean streamflow observation. $E{T}_{RS}$ and $E{T}_{sim}$ are the 8-day time series for the spatial ET_RS observations and model ET simulations, respectively. The sensitivity index is then defined as the ratio between the difference and the sum of the maximum and minimum NSC (streamflow) per catchment or APBIAS (spatial ET) per grid cell among each parameter space setting (Equation (3)). A parameter is considered as more sensitive if the sensitivity index is closer to 1 than other parameters.

$$SensitivityIndex\left(streamfloworspatialET\right)=\left|\frac{Ma{x}_{NSCorAPBIAS}-Mi{n}_{NSCorAPBIAS}}{Ma{x}_{NSCorAPBIAS}+Mi{n}_{NSCorAPBIAS}}\right|.$$

(3)

2.2.3. Calibration Model

In the third step, the spatial ET_RS observations and Q_obs calibration model are executed respectively using the Shuffled Complex Evolution (SCE-UA) global optimization method [70], which follows the modeling framework developed by Wu et al. [3]. The SCE-UA parameters (e.g., the number of complexes, the number of points in each complex) are recommended by [70] and 500 minimum iterations are set unless the parameters or objective function reach convergence prior to this threshold, within a ±1% level of further performance improvement. The objective functions of the Q_obs and spatial ET_RS calibration schemes are NSC (streamflow) and APBIAS (spatial ET), respectively.

The uncertainties in the ET_RS and precipitation inputs were investigated; whereby, a synthetic ET calculation was derived from the simulation using the Q_obs calibration optimized model parameters. A similar approach was used by [21] with favorable results. The synthetic ET results were then used for the spatially distributed model calibration, providing a reference to evaluate uncertainty contributed by the ET_RS observations. A high-quality precipitation dataset is used in this study [52], and a simple correction approach is applied to further refine the precipitation product according to [35]. Specifically, annual reference precipitation at the basin scale is defined as the sum of the observed streamflow and ET_RS for each hydrological year, assuming negligible annual soil water storage change for a water-budget closed river basin between hydrological years. The raw daily precipitation used in this study is then bias corrected by the annual reference precipitation. Such annual bias correction for precipitation input is proven to have positive impacts on streamflow simulations [35]. The further bias-corrected precipitation is then used to drive the model simulations for the comparison, and help understand the uncertainty contributed from the precipitation inputs.

2.2.4. Model Validation and Evaluation

In the fourth step, the effectiveness of the ET_RS calibration on the streamflow simulation is evaluated by comparing the model results against independent observations, and through inter-comparison of outputs from different model experiments. Three model validation runs were performed for evaluating the Q_obs calibration model (Q_VAL); spatial ET_RS calibration model (ET_VAL); the baseline simulation model with the latest community efforts optimized soil information (BL_VAL).

The Kling–Gupta Efficiency (KGE) metric [71] was used to measure the agreement between observed and simulated monthly streamflow time series. The value of KGE ranges from negative infinity (no agreement) to one (perfect agreement).

$$\begin{gathered} KGE\left(streamflow\right)=1-\sqrt{(\alpha -1{)}^2+(\beta -1{)}^2+(\gamma -1{)}^2} \\ \alpha =\frac{{\delta }_Q{}_{sim}}{{\delta }_{Q_{obs}}},\beta =\frac{{\mu }_Q{}_{sim}}{{\mu }_Q{}_{obs}},\gamma =\frac{1}{n}{\displaystyle \sum }_{i=1}^n\frac{\left(Q_{obs}-{\mu }_Q{}_{Obs}\right)*\left(Q_{sim}-{\mu }_Q{}_{Sim}\right)}{{\delta }_Q{}_{Sim}*{\delta }_{Q_{Obs}}}, \end{gathered}$$

(4)

where ${\mu }_Q{}_{obs}$ and ${\delta }_Q{}_{sim}$ are the mean and standard deviation of the streamflow observations; ${\mu }_Q{}_{sim}$ and ${\delta }_Q{}_{sim}$ are the mean and standard deviation of the streamflow simulations. We defined three general performance categories for the ET_VAL results, as improved (KGE improvement for streamflow > +10%) comparable (KGE improvement for streamflow within ± 10%), or degraded (KGE improvement for streamflow < −10% ) of their counterparts, respectively.

3. The Test Area and Datasets

The Columbia River Basin (CRB) encompasses approximately 668,000 km² of the Pacific Northwest region of North America including portions of the Rocky Mountains in British Columbia, Canada. The CRB encompasses a wide range of physical and ecological environments ranging from semiarid and desert areas of the central plateaus to relatively wet forests in the Cascade Mountains [72]. The Columbia River drains the CRB and represents the fourth largest river in the continental United States. The water resources in the CRB have been extensively exploited in the past 70 years to meet various needs such as hydropower production, flood control, agricultural withdrawals, instream flow requirements for fish, and recreational needs [73].

River gauges in the CRB domain maintained by the U.S. Geological Survey (USGS) are selected for this study and used for model calibration and validation. Three general criteria were used for stream gauge selection, including (1) having continuous monthly discharge records from January 2004 to December 2013; (2) being identified as “reference gauges” (least-disturbed by human activities) in the Geospatial Attributes of Gages for Evaluating Streamflow (GAGES-II) database [74,75]; and (3) having dimensional consistency with the dominant river tracing (DRT) hydrography dataset defined at 1/16-degree spatial resolution [76]. The DRT river network upscaling algorithm, utilized 1-km-resolution HydroSHEDS (Hydrological data and maps based on SHuttle Elevation Derivatives at multiple Scales) baseline hydrography inputs, and successfully preserves critical sub-grid scale basin attributes affecting runoff and river flows [77,78]. To meet the requirements for criteria (3) above, the absolute relative difference between drainage areas derived from the DRT datasets and reported by the USGS has to be less than 20%, while visual inspection is also used to ensure geolocation accuracy. As a result, a total of 28 sub-basins in the CRB were selected for this study with sizes ranging from 325.70 km² to 14,331.70 km² (Figure 2). Detailed characteristics of the 28 basins can be found in

Table 3. Gauge information and basin characteristics of the 28 natural catchments examined in this study. Gauge information is retrieved from the USGS database, while basin characteristics are obtained from the GAGE-II database.

No.	Gauge Information					Basin Characteristics
	ID (USGS)	Latitude (degree)	Longitude (degree)	Drainage area (km²)	Streamflow (mm/year)		Precipitation (mm/year)	Temperature (°C)	Elevation (m)
1	12413875	47.06	−115.35	325.7	946.1		1383.8	5.0	1640.3
2	14222500	45.84	−122.47	333.1	1785.2		2455.4	9.4	573.3
3	14236200	46.60	−122.46	385.8	1791.7		1772.6	8.6	672.9
4	13334450	46.27	−117.29	395.0	115.1		561.8	7.8	1271.8
5	14159200	44.05	−122.22	435.7	1336.5		1782.8	7.6	1280.8
6	14185000	44.39	−122.50	449.4	1573.7		1998.6	9.0	889.5
7	12458000	47.54	−120.72	532.9	973.2		1655.0	4.7	1547.5
8	14154500	43.74	−122.87	539.7	927.0		1561.4	9.9	857.5
9	12452800	47.82	−120.42	553.0	559.1		962.4	5.7	1519.8
10	12451000	48.33	−120.69	878.2	1419.1		1494.8	4.9	1534.0
11	13046995	44.06	−111.15	939.2	806.7		1175.2	2.8	2349.5
12	13161500	41.93	−115.68	994.2	96.0		450.4	6.0	2045.0
13	12447383	48.57	−120.39	1065.9	401.4		984.8	3.7	1685.7
14	12411000	47.71	−115.98	1078.6	536.3		1254.2	6.4	1202.5
15	13235000	44.09	−115.62	1187.7	570.8		936.6	4.9	2079.0
16	13010065	44.10	−110.67	1304.6	571.9		981.4	2.0	2508.5
17	14231000	46.53	−121.96	1492.9	1715.3		1804.6	6.5	1128.6
18	12448000	48.48	−120.19	1784.1	206.1		617.4	4.4	1613.0
19	13185000	43.67	−115.73	2144.7	455.1		883.8	5.9	1955.1
20	13296500	44.27	−114.73	2168.3	358.0		753.8	3.1	2375.2
21	12414500	47.27	−116.19	2656.6	724.4		1231.0	6.3	1381.3
22	13309220	44.72	−115.01	2666.8	455.6		897.8	3.8	2192.4
23	12413000	47.57	−116.25	2712.7	587.1		1207.0	6.7	1167.8
24	13338500	46.09	−115.98	3077.4	227.4		863.8	6.1	1384.6
25	12358500	48.50	−114.01	3239.1	771.9		1166.6	3.8	1723.7
26	13340600	46.84	−115.62	3354.2	846.4		1426.0	5.9	1442.7
27	13337000	46.15	−115.59	3962.2	579.0		1245.0	5.3	1584.1
28	13340000	46.48	−116.26	14,331.7	493.5		1021.2	6.0	1443.6

.

A recently updated gridded meteorological dataset for the Continental United States by Livneh et al. [52] is adopted for this study, which is an extended work of Maurer et al. [48] but has longer temporal coverage (1950 to 2013) and finer spatial resolution (1/16-degree). The quality-controlled dataset contains gridded precipitation, daily maximum and minimum temperature aggregated from point observations at National Weather Service Cooperative Observer (NWS COOP) stations, and surface winds from the National Centers for Environmental Prediction (NCEP)/ National Centers for Atmospheric Research (NCAR) reanalysis products. The gridded precipitation data are interpolated using the synergraphic mapping system (SYMAP) [79] and are scaled to match the Parameter-Elevation Regressions on Independent Slopes Model (PRISM) 1981–2010 precipitation climatology for topographic correction.

The widely used 8-day, 1-km resolution MODIS (MODerate resolution imaging spectroradiometer) MOD16 ET_RS operational product [12] was used in this study. The MOD16 ET algorithm uses a modified Penman–Monteith algorithm driven by daily meteorological inputs from the NASA Global Modeling and Assimilation Office (GMAO) reanalysis, along with MODIS derived land cover, albedo, LAI, and enhanced vegetation index observations for global ET mapping and monitoring. The quality of the MOD16 ET record has been validated using tower eddy covariance-based daily ET measurements over North America [12,80,81,82]. Like other ET_RS estimations, validation results of MOD16 products showed positive and negative bias in different places. Several studies have also used the MOD16 ET record for hydrological model calibration [20,27,30,31] and data assimilation [83]. For this investigation, the MOD16 ET record was geo-referenced and spatially aggregated to the 1/16-degree resolution PNW model grid.

4. Results

4.1. The Baseline Simulation

The hydrological model produced generally improved streamflow simulations in 68% (19 out of 28) of the test river basins when using the Updated LAI (mean KGE: 0.52) rather than the Default LAI (mean KGE: 0.45) as model inputs (Figure 3a). The differences in numerical streamflow modeling between the Default LAI and Updated LAI derived results can be large for some river basins. For river basin No.28, for example, Evergreen Needleleaf forest presents the highest LAI (Figure 3b,c), followed by Mixed Cover and Woodland vegetation types from both Default LAI and Updated LAI records. Significant discrepancies exist between the two LAI products, especially for the Evergreen Needleleaf class. The Default LAI fails to capture monthly variations in Evergreen Needleleaf forest, but generally captures similar monthly variations as the Updated LAI record for other land cover types, but with significant value differences, particularly for Mixed Cover. The LAI values are smaller and more spatially variable in the Updated LAI record than in the Default LAI for all the three major land cover classes, which is consistent with previous studies [84,85]. Lower LAI reduces transpiration from vegetation and evaporative losses from canopy interception, which promotes wetter soil conditions and greater surface and sub-surface runoff in the model [85].

Using a more temporally updated LAI record and more heterogeneous LAI data results in a marked improvement in the model streamflow predictions, as the seasonal LAI record is significant in the role of partitioning rainfall to runoff and ET in the VIC module. All of the following sections are based on model simulations using the Updated LAI record as a key input, allowing this study to build on recent advances in landscape parameterization.

4.2. Sensitivity Analysis of Soil Parameters

Model sensitivity tests were performed to identify the important parameters affecting both streamflow and the spatial ET simulations. Figure 4a suggests that all six of the soil parameters examined (Table 2) are essential for accurate streamflow simulations across all 28 study catchments, though the relative importance of the parameters varies across basins. However, the Ws and Ds parameters show the most significant individual influence on ET at the river basin scale (Figure 4b). The Dsmax and D3 parameters tend to have lower model sensitivity but are more effective than parameter D2. The model ET seems to show little sensitivity to parameter B in most of the test river basins.

Plots showing the sensitivity of the six model parameters at each individual grid cell in the four landscape groups are presented in Figure 4c; these results show the relative model sensitivity to the six parameters, which is consistent with the river basin scale results (Figure 4b). Most of the grid cell points are at the top of the Ws and Ds plots, indicating that the model ET calculations are highly sensitive to both parameters for all land cover classes (Figure 4c). The model is also sensitive to the other four parameters, but the sensitivity is more geolocation dependent. However, there were no significant differences in the model ET sensitivities to the model parameters among the four main landcover classes. However, the D2 parameter showed relatively less impact on the model ET calculation in woodlands due to the relatively thinner soil layer represented by this land cover class. Therefore, we allowed all six model parameters to vary in the Q_obs calibration consistent with the majority of previous studies, while assuming B and D2 have no substantial impact on spatial ET_RS calibration. Further research is warranted to investigate the validity of this assumption for other areas.

4.3. Uncertainty in the Calibrated Model Parameters

The model parameter ensembles paired with the NSC (streamflow) values from Q_obs calibration are shown in Figure 5, six significant basins are selected to show the uncertainties in the calibrated model parameters against Q_obs. For each of the six selected river basins, grey points denote all 500 evolutions in the Q_obs calibration runtime; the 25 black points in each plot refer to equifinal simulations with NSC values very close to each other (within 1% difference), while the red point denotes the best streamflow prediction. The distribution of black points in the figure plots shows that Q_obs calibration can yield very good streamflow predictions in most cases, with generally narrow ranges in the model parameter values. However, a lack of convergence in the calibrated parameters can be found in some river basins, indicating that “seemingly right” answers could be produced for the wrong reasons [30]. For example, at river basin No.16, all six model parameters except Ds, show large variations, indicating the low quality of model inputs (precipitation in particular), which may give the optimization algorithm more freedom to fit the objective function.

Keeping water balance closure can help ensure that a model with parameters constrained by reliable ET_RS observations yields favorable runoff and streamflow predictions [17]. As shown in Figure 6, a number of river basins (e.g., No.2, 5, 7, 9, 10, 15, 16) show that improved ET prediction propagates to more accurate streamflow simulations. The model parameter set resulting in the best streamflow prediction (green star) is not (even close to) the one with the best ET simulation in many cases (except No. 15). This reflects both the limitations in the ET_RS observations used for model ET calibration and the precipitation inputs. Effective spatial calibration of the parameter fields will be constrained as long as the errors in ET_RS and the precipitation inputs remain larger than the variability in the calibrating parameters [32]. However, in all river basins, the model streamflow performance (i.e., KGE metrics) increases or decreases consistently with the gradient in ET bias. The calibration scheme steadily modulates the interplay of the errors in precipitation and ET_RS, resulting in progressive changes of KGE values for streamflow simulation. For instance, if precipitation is biased high and ET_RS biased low, a larger positive bias in ET would lead to a high KGE value for the streamflow simulation, while a successful ET calibration based on a relatively low biased ET_RS would not necessarily produce a favorable streamflow simulation KGE result.

As shown in Figure 6, when the ET and streamflow simulations simultaneously achieve good calibration performance (i.e., KGE for streamflow greater than 0.7), APBIAS for ET is less than 10% and the KGE value is very close to the best run (green star) result, which would also indicate good quality ET_RS and precipitation inputs (e.g., the No. 15 river basin). On the other hand, if both the best run (green star) and corresponding streamflow simulation using the optimized model ET have a negative KGE value, the precipitation input may be considered as very low quality (see No.9 and 13). All of the rest (90%) of the test river basins in Figure 6 achieve positive KGE values with a mean of 0.50 for the model streamflow simulations after the ET calibration, indicating that the ET_RS observations and precipitation inputs are of reasonable quality for most of the study basins. However, the distance between the actual (red point) and best run (green star) KGE results indicates potential room for further model improvements using better quality precipitation and/or ET_RS records.

It is worth noting that incorporating the spatial distribution of the estimated parameters in the model calibration is a distinct advantage of the spatial ET_RS calibration, which is considered a major challenge for large-domain hydrological modeling [9]. As shown in Figure 7, taking river basin No.28 as an example, the spatial ET_RS optimized parameters show large spatial variability, while the Dsmax and Ds spatial distributions are similar and generally congruent with the pattern of soil texture (bulk density) across the catchment.

4.4. Validation and Evaluation

As shown in

Table 4. The KGE metrics of streamflow predictability in BL_VAL, ET_VAL, and Q_VAL setups over the 28 study basins (involved in Step-Four). The basins No. 16, 20, 22, 25, 26, and 28 are illustrated in Figure 9 and highlighted below (in bold). The “MOD16 ET*” represents the optimum parameter set retrieved offline when observed and simulated streamflow fit best (green star in Figure 6).

	ID	Raw Precipitation					Adjusted Precipitation
NO.	(USGS)	BL_VAL	ET_VAL (MOD16 ET)	ET_VAL (MOD16 ET*)	ET_VAL (Synthetic ET)	Q_VAL	BL_VAL	ET_VAL (MOD16 ET)	Q_VAL
1	12413875	0.53	0.63	0.55	0.69	0.56	0.71	0.74	0.62
2	14222500	0.48	0.47	0.51	0.47	0.49	0.57	0.57	0.55
3	14236200	0.68	0.66	0.78	0.93	0.86	0.84	0.77	0.87
4	13334450	0.49	0.31	0.50	0.63	0.13	0.16	0.54	0.57
5	14159200	−0.12	0.23	0.35	0.20	0.52	0.06	0.36	0.48
6	14185000	0.52	0.56	0.66	0.56	0.69	0.48	0.54	0.66
7	12458000	0.73	0.85	0.83	0.90	0.92	0.78	0.89	0.91
8	14154500	0.56	0.60	0.76	0.57	0.72	0.57	0.62	0.58
9	12452800	0.04	−0.11	−0.02	0.03	0.29	0.56	0.43	0.66
10	12451000	0.74	0.80	0.86	0.85	0.90	0.58	0.66	0.86
11	13046995	0.36	0.43	0.90	0.73	0.91	0.17	0.25	0.95
12	13161500	0.37	0.26	0.57	0.34	0.74	0.26	0.81	0.72
13	12447383	0.42	0.26	0.37	0.43	0.45	0.49	0.59	0.83
14	12411000	0.56	0.60	0.83	0.62	0.82	0.38	0.54	0.89
15	13235000	0.69	0.89	0.82	0.73	0.85	0.54	0.54	0.90
16	13010065	0.85	0.91	0.91	0.55	0.86	0.88	0.90	0.90
17	14231000	0.63	0.70	0.76	0.74	0.78	0.51	0.68	0.68
18	12448000	0.28	−0.02	0.31	0.47	0.48	0.32	0.46	0.79
19	13185000	0.66	0.95	0.91	0.68	0.83	0.77	0.79	0.92
20	13296500	0.46	0.44	0.52	0.68	0.58	0.67	0.78	0.76
21	12414500	0.83	0.72	0.80	0.73	0.77	0.77	0.72	0.85
22	13309220	0.74	0.76	0.84	0.51	0.85	0.59	0.87	0.80
23	12413000	0.32	0.40	0.68	0.45	0.82	0.26	0.38	0.81
24	13338500	0.49	0.48	0.79	0.48	0.78	0.27	0.41	0.88
25	12358500	0.66	0.61	0.77	0.88	0.92	0.62	0.61	0.88
26	13340600	0.71	0.62	0.85	0.61	0.86	0.67	0.73	0.80
27	13337000	0.62	0.56	0.70	0.51	0.82	0.60	0.56	0.87
28	13340000	0.62	0.59	0.81	0.52	0.89	0.69	0.63	0.83
Average		0.53	0.54	0.68	0.59	0.72	0.53	0.62	0.78

, in the “Raw Precipitation” setup, KGE metrics from the validation simulations derived using model parameters from the ET_RS calibration range from −0.11 to 0.95 (mean: 0.54). Approximately 61% of the study basins have KGE > 0.5, which is consistent with the calibration period results. Compared to the baseline validation, BL_VAL, in which the KGE values are between −0.12 and 0.85 (mean: 0.53).

The ET_VAL results show significant improvements in streamflow simulation, with 50% and 25% of the river basins showing improved and comparable performance, respectively. When compared to the Q_VAL, the ET_VAL results show improved and comparable performance for 18% and 11% of river basins, respectively, while 71% of the river basins show degraded performance (Figure 8a).

Using synthetic ET data for the model calibration causes significant improvements in model performance, with improved or comparable performance to BL_VAL for 64% and 11% of the study basins, respectively (Figure 8b). Thus, relatively closer model performance to Q_VAL is achieved using synthetic ET for the calibration, with a total of 43% of the river basins showing improved (14%) or comparable (29%) performance to Q_VAL. Using simply corrected precipitation in the ET_RS calibration and validation, the streamflow simulations are also improved significantly over the ET_VAL results, with 39% of the study basins showing improved (21%) or comparable (18%) KGE metrics against Q_VAL (Figure 8c). Surprisingly, the simply corrected precipitation causes further marked improvement in the model performance, with 96% of the total river basins showing improved (71%) and comparable (25%) performance than the BL_VAL results. Further improvements in the model performance for both ET and streamflow simulations are also expected with better quality ET or precipitation inputs as shown in Figure 8b,c.

Furthermore, almost all the best model performances (in terms of streamflow) are achieved by using the parameter sets from the best KGE runs (i.e., green stars in Figure 6) using the ET_RS calibration without involving streamflow observations. The best model runs represent 89% and 4% of the total study basins with improved and comparable performance to BL_VAL, respectively, and only 7% of the river basins with degraded metrics (Figure 8d). Interestingly, the ET_VAL results appear to be even better than the Q_VAL results, with 29% and 39% of the study basins showing improved and comparable KGE values, respectively. These results promote additional questions regarding how the best model runs can perform better than Q_VAL given the same maximum time of optimizing evolutions for the calibrations, and why the ET_RS calibration approach can achieve very good streamflow simulations more efficiently than the Q_obs calibration. Such improvements are likely due to the spatial distribution of ET_RS, which may help the calibration scheme efficiently tune the model ET calculations at each grid cell toward optimal magnitudes, thereby improving the spatial distribution of runoff relative to the Q_obs based model calibration. Understanding the internal processing structure to identify the best model run and optimal outputs under a diversity of climate and landscape conditions remains challenging and is warranted in future studies.

As shown in Figure 9, six river basins were selected to illustrate the ET_VAL simulated monthly hydrograph comparisons with Q_VAL, including two improved (No.20 and No.22), two comparable (No.16 and No.26), and two degraded (No.25 and No.28) basins from the adjusted precipitation setup (see Figure 8c). In the CRB domain, the snow melts in spring and accumulates in winter; thus, seasonal variations of modeled hydrograph are all well captured without consideration of parameter optimization method. For Basin No.20 and No.22 with improved performances of streamflow prediction, the ET_VAL results generally capture more flows than the Q_VAL results. For Basin No.16 and No.26 with comparable performances of streamflow prediction, the ET_VAL results for Basin No.16 overestimate flows in general and capture high flows well but fail in simulating low flow for Basin No.26. For both Basin No.25 and No.28 with degraded performances of streamflow prediction, ET_VAL underestimates high flows but overestimates low flows. All streamflow simulations do not capture the peak flows in the wet seasons of 2010 of all six basins. It is probable that the adjusted precipitation input here is insufficient to reflect discharge magnitudes.

5. Conclusions

In this study, we hypothesized that hydrological models calibrated against improved ET_RS observations would lead to better streamflow predictions. A modeling framework developed by Wu et al. [3], including the widely used Variable Infiltration Capacity (VIC) model [42], the dominant river-tracing-based source-sink hydrologic routing model [43] and the SCE-UA global optimization algorithm [70] was adopted for model calibration, validation, and evaluation for 28 natural river basins in the Pacific Northwest region of the continental United States. The widely used MODIS MOD16 operational ET product [12] provided a spatially and temporally dynamic ET observational record (ET_RS) for this study. A four-step workflow was designed to evaluate the impacts of the ET_RS based spatially distributed model calibration for the streamflow simulations, and to clarify how the model parameters determine reasonable ET estimates coherent with other water budget components including precipitation, soil moisture, and runoff.

The baseline model simulations for this study used the latest VIC model inputs including precipitation, vegetation, and soil datasets. The resulting model simulations showed significant improvement in ET and streamflow simulations using a recently updated Leaf Area Index (LAI) climatology [66] as a critical input. Comprehensive model sensitivity tests were performed indicating a set of six critical model parameters for calibrating against Q_obs, and four parameters for calibrating against ET_RS (Table 2). Q_obs based lumped calibration and ET_RS based spatially distributed calibration were designed for this study, while uncertainties of precipitation and ET were considered to improve model performances (Figure 1). The validation results show that the model using the optimized parameters from the ET_RS calibration results in a mean KGE value of 0.54, while 61% of the 28 natural river basins tested had KGE values exceeding 0.5. Importantly, the ET_RS calibration scheme not only resulted in significant improvements over the baseline model, with 75% of the basins showing improved or comparable KGE metrics, but also resulted in consistent performance between model calibration and validation periods. Our results also indicate that either improved ET observations (indicated from a synthetic ET calibration experiment) or better precipitation inputs (through further bias correction) can provide further improvements to the effectiveness of the ET based calibration for better streamflow simulation. Compared to the Q_obs calibration-based results, the synthetic ET based calibration produces better or comparable performance in 43% of the study basins, while the ET_RS calibration using further annual bias-corrected precipitation produced better or comparable performance in 39% of the basins.

Interestingly, using the parameter sets from the best model calibration runs also produced the best overall model performance, with 93% and 68% of the study basins showing improved or close performance compared to the baseline simulations and the Q_obs calibration-based results, respectively. Given reasonably accurate precipitation and ET_RS records, the spatially distributed calibration scheme can be efficient in parameter optimization and effective in tuning the model for accurate streamflow calculations. However, if the precipitation inputs have a significant bias, the model calibration may still provide accurate streamflow simulations, but at the expense of degraded ET performance as the model ET parameters may be adversely adjusted to compensate for the bias in precipitation.

Our validation and evaluation results show that the ET_RS calibration can efficiently tune the relevant model parameters for better ET and streamflow simulations within their physically meaningful ranges and with reasonable spatial distribution (e.g., consistent with soil texture patterns), which is likely due to the positive contribution from the spatial distribution information of ET_RS. In gauged river basins with reliable Q_obs data, the utility of ET_RS in hydrological model calibration can further improve model performance by providing additional observational information to constrain the calibration. Moreover, our results demonstrate efficient and effective model calibrations using ET observations from global satellite records. These results are expected to be particularly beneficial for global-scale applications and natural basins lacking river gauge networks. The relative effectiveness of the model calibrations is expected to improve further as the quality of the ET_RS observations continue to improve. A similar model calibration framework may apply to other satellite-based hydrological observations including soil moisture, surface inundation, snow cover, and river discharge. Further studies are needed to test the value of these observations for enabling further model calibration and performance enhancements.