Improving streamflow predictions at ungauged locations 1 with real-time updating : Application of an EnKF-based 2 state-parameter estimation strategy 3

The challenge of streamflow predictions at ungauged locations is primarily attributed to various uncertainties in hydrological modelling. Many studies have been devoted to addressing this issue. The similarity regionalization approach, a commonly used strategy, is usually limited by subjective selection of similarity measures. This paper presents an application of a partitioned update scheme based on the ensemble Kalman filter (EnKF) to reduce the prediction uncertainties. This scheme performs real-time updating for states and parameters of a distributed hydrological model by assimilating gauged streamflow. The streamflow predictions are constrained by the physical rainfall-runoff processes defined in the distributed hydrological model and by the correlation information transferred from gauged to ungauged basins. This scheme is successfully demonstrated in a nested basin with real-world hydrological data where the subbasins have immediate upstream and downstream neighbours. The results suggest that the assimilated observed data from downstream neighbours have more important roles in reducing the streamflow prediction errors at ungauged locations. The real-time updated model parameters remain stable with reasonable spreads after short-period assimilation, while their estimation trajectories have slow variations, which may be attributable to climate and land surface changes. Although this real-time updating scheme is intended for streamflow predictions in nested basins, it can be a valuable tool in separate basins to improve hydrological predictions by assimilating multi-source data sets, including ground-based and remote-sensing observations.


Introduction
The streamflow prediction plays a central role in hydrology because it is an important element for water resources management, the design of hydraulic infrastructures and flood risk mapping (Srinivasan et al., 2010).Because it is an important component in the terrestrial water budget, streamflow is also a direct diagnostic variable measuring the impact of climate changes and human activities that act on a given watershed.Streamflow prediction depends highly on reliable hydrological data and sophisticated hydrological models.However, hydrological data are often insufficient due to ungauged or poorly gauged basins in many parts of the world (Sivapalan, 2003).Because of the scarcity of data, hydrological modelling is also plagued by various sources of uncertainties.To reduce uncertainties from those hydrological data and hydrological modelling, the International Association of Hydrological Sciences (IAHS) launched an initiative on Predictions in Ungauged Basins (PUB) (Sivapalan, 2003;Sivapalan et al., 2003).
Through the past PUB decade, major advances have been achieved including data acquisition and exploitation, modelling strategies and uncertainty analysis, and catchment classification and new theory (Hrachowitz et al., 2013).There is a growing consensus that remote sensing techniques provide valuable data for understanding the land surface hydrological system (Yang et al., 2013).Moreover, considerable progress has been made on hydrological models (typically the distributed hydrological models) to capture the physical process associated with the basin rainfall-runoff and Published by Copernicus Publications on behalf of the European Geosciences Union.
X. Xie et al.: Improving streamflow predictions at ungauged locations with real-time updating snowmelt-runoff responses.This progress has fostered specific problem areas in the field: uncertainty quantification with respect to model input forcing, model structures and parameters (Ajami et al., 2007;Vrugt et al., 2008;Gupta et al., 2012).To reduce the uncertainty from model parameters, one common practice is the parameter calibration by adjusting model parameters to make the simulated water discharges correspond to the observations (typically the data from the outlet of a watershed) (Duan et al., 1992(Duan et al., , 1994)).However, a calibrated parameter set with acceptable streamflow simulation performance at the watershed outlet does not guarantee the performance at interior locations (Zhang et al., 2008).
The essence of PUB is to transfer information from neighbouring basins to the basins of interest (Sivapalan et al., 2003).Such process is generally referred to as hydrological regionalization, based on either regression methods or measurable distances (with respect to physical similarity or spatial proximity) between gauged and ungauged locations (Hrachowitz et al., 2013).Regionalization techniques regarding model parameters are popular for discharge prediction in ungauged basins.Merz and Blöschl (2004) evaluated the performance of various regionalization methods for parameters of a conceptual catchment model, determining that spatial proximity is able to represent the unknown controls on the runoff regime and the relationships of model parameters within neighbouring basins.Sellami et al. (2014) presented a model parameter regionalization approach based on physical similarity between gauged and ungauged catchments, indicating that similar hydrological behaviour may appear due to physically similar catchments in the same geographic and climatic region.Parajka et al. (2013) reported that the spatial proximity and geostatistics probably perform better than the regression or regionalization with a simple averaging of model parameters from gauged catchments.One drawback of the regionalization of model parameters is that it often confronts an arbitrary criterion for selecting the "behavioural" model parameter sets from the gauged catchment (Sellami et al., 2014).Hrachowitz et al. (2013) provides a comprehensive review of the parameter regionalization and catchment similarity.
In addition to those parameter regionalization approaches, newly developed data assimilation methods are also encouraging and are capable to address some issues associated with PUB.They are generally based on physical correlations between the neighbouring basins, and they can combine multisource observations to transfer information from gauged to ungauged basins (Sivapalan et al., 2003;Troch et al., 2003;Chen et al., 2011).As a typical sequential data assimilation approach, the ensemble Kalman filter (EnKF) is popular in hydrology (Reichle et al., 2002;Evensen, 2003Evensen, , 2009)).The EnKF is attractive in hydrology primarily because it can perform real-time updating with simple implementation and it considers various uncertainties in modelling and observations (Blöschl et al., 2008).The feature of real-time updating is very important for flood forecasting (Norbiato et al., 2008).
In some current applications, EnKF is mainly dedicated to dynamic state estimations in which the model parameters are defined with prior values or calibrated in advance (Vrugt et al., 2005;Clark et al., 2008).
The EnKF method also provides a general framework to perform state-parameter estimation which is the core of PUB issues.It is has been successfully used for parameter estimation of hydrological models.Moradkhani et al. (2005b) proposed a dual state-parameter estimation of hydrological models and made an acceptable application of this method for a lumped hydrological model.Wang et al. (2009) presented three constrained schemes with EnKF to prevent the violation of parameter physical constraints.Most of these studies performed parameter estimations for lumped hydrological models with a small number of parameters to be estimated.Xie and Zhang (2010) successfully demonstrated a joint state-parameter estimation based on the EnKF for a distributed hydrological model, i.e.Soil and Water Assessment Tool (SWAT), focusing on one dominant parameter in SWAT.For multiple types of parameter estimation, Xie and Zhang (2013) developed a partitioned update scheme and indicated the potential of this scheme for streamflow predictions in ungauged basins based on distributed hydrological models.
In this study, we present the application of the partitioned update scheme to improve streamflow predictions in ungauged locations by assimilating gauged streamflow.This data assimilation algorithm is fully coupled with the distributed hydrological model, i.e.SWAT.The state vector and parameters in ungauged subbasins are estimated when information is transferred from gauged subbasins.To our knowledge, this study is the first one which explicitly employs a data assimilation method with state-parameter estimation to improve streamflow predictions in ungauged locations.Although a few applications of data assimilation methods are dedicated to streamflow predictions based on distributed models (Clark et al., 2008;Chen et al., 2011;Lee et al., 2012;Rakovec et al., 2012;McMillan et al., 2013), the model parameter estimation, which is important for PUB, is not systematically considered.In addition to the EnKFbased scheme, note that the other data assimilation methods, e.g. the particle filter (Moradkhani et al., 2005a;DeChant and Moradkhani, 2012), the Particle-DREAM (Vrugt et al., 2013) and the Maximum Likelihood Ensemble Filter (Tran et al., 2014), may also be optional for state-parameter estimation.
In the following sections, we first introduce the EnKFbased data assimilation scheme and give a brief description of the SWAT model.We then present an application case concerning a real-world problem in the Zhanghe River basin in China in which river channels are connected and subbasins have nested upstream and downstream neighbours.Three scenarios regarding different combinations of observed streamflow are designed to discuss the impact of Hydrol.Earth Syst.Sci., 18, 3923-3936, 2014 www.hydrol-earth-syst-sci.net/18/3923/2014/ gauged locations on streamflow predictions.Finally, conclusions are given in the last section.

EnKF-based state and parameter estimation scheme
To describe the information transfer process from gauged to ungauged locations, we define a joint state vector X that contains gauged (x g ) and ungauged (x u ) states: Moreover, we consider the diagnostic variables, i.e. the water discharge and the evapotranspiration, as model states and include them in the vector X to perform streamflow updating in the data assimilation.The joint state vector X and the parameter vector θ estimation at time t are conditioned on measurements (y t ) from gauged basins.The information transfer process, i.e. the posterior probability density function (pdf) p( X t , θ t | y t ), can be expressed within Bayes' framework, where p( y t | X t , θ t ) is the likelihood function of measurements given model estimations at time t.Moreover, p( X t , θ t | X t−1 , θ t−1 ) is the prior pdf of X and θ at time t that represents model forecasting and parameter evolutions.
The updating framework defined in Eq. ( 1) is well included in and effectively solved by sequential data assimilation strategies -typically, the EnKF strategy (Evensen, 1994).The EnKF strategy operates sequentially with a forecast step and a filter update step.In the forecasting process, uncertainty propagation is characterized by an ensemble of model realizations: where "−" and "+" denote the forecast and analysis for the state vectors X and the parameter vector θ, t is the time step, u is the input forcing vector, and N is the ensemble size.The model error vector ω is assumed to follow a Gaussian distribution with zero mean and covariance W t .Equation ( 2) is a general expression with representative errors for all state variables.In implementation, one may define errors for only a few of the state variables (e.g.soil moisture) to reflect realistic modeling uncertainties.Detailed prescription of the errors will be given in Sect.3.2.Prior to model forecasting using Eq. ( 2), the model parameters can be perturbed, similar to the forecast of the state vector, to avoid the shrinkage of the parameter ensemble during the updating (Wang et al., 2009).However, the parameter perturbation is susceptible to over-dispersion in sampling (Moradkhani et al., 2005b).A kernel smoothing technique is effective to address the over-dispersion while maintaining a reasonable ensemble spread for the parameters (Liu, 2000;Moradkhani et al., 2005b;Xie and Zhang, 2013).This technique is briefly expressed as where α is the shrinkage factor typically within [0.95, 0.99], h is the smoothing factor, and T t is the covariance constrained by the ensemble variance var θ + t .The smoothing factor h is defined as √ 1 − α 2 to maintain equal variances of the parameter before and after the perturbation.This kernel smoothing technique has been discussed based on synthetic cases (Liu, 2000;Moradkhani et al., 2005b;Xie and Zhang, 2013), so we do not provide any more experiments to demonstrate the properties of the kernel smoothing.The prescription of the shrinkage factor α is subject to trial and error experimentation, but it has limited impact on the parameter estimation (An illustrative case was shown in the response to the reviewers' comments at version 4 of this paper).In this study, it is specified with 0.98 according to the suggestions by Moradkhani et al. (2005b) and Xie and Zhang (2013).
With the forecast of the states and parameters, the filter update step is performed when observations are available.This updating is actually the solving process for Eq. ( 1).Here we intentionally create an explicit expression of the updating for gauged and ungauged states and parameters: where y i t is the observation vector, which is appropriately perturbed using covariance of R to account for uncertainties in observations, and H is the observation operator and it is linear in this study.The Kalman gain matrix K t is expressed as where cov(•) is the covariance operator that is computed from the ensembles of states and parameters.Note that the size of the matrix K t is n × m, where n is the total number of state variables and parameters and m is the number of observations.
The above two equations rely on EnKF with a stateaugmentation technique.This technique is valid and able to retrieve correct parameter estimates in real time primarily

X. Xie et al.: Improving streamflow predictions at ungauged locations with real-time updating
because it allows for parameter dynamics and performs the parameter evolution.Specifically, model parameters are assumed as an extension of state variables and they can travel slowly with time, in response to changes in environmental forcing inputs (Liu and Gupta, 2007).Like the model state forecasting, the parameters are perturbed/evolved using the kernel smoothing technique.In this way, the evolution of model parameters is consistent with the forecasting of model state variables.Thus the model parameters can be appended to the state vector (Moradkhani et al., 2005b;Xie andZhang, 2010, 2013).When observations are available, the parameters are updated along with state variables by assimilating these observations.Therefore, their estimates are expected to converge to the "correct" posterior target distribution (Xie and Zhang, 2013).This technique has been successfully used in many cases for real-time state and parameter estimation (Moradkhani et al., 2005b;Wang et al., 2009;Xie andZhang, 2010, 2013).
We can see that EnKF provides a general framework to transfer information from gauged to ungauged basins.However, when used for parameter estimations in distributed hydrological models, it is vulnerable to corruption due to spurious covariance computation in Eq. ( 7), primarily resulting from a large degree of freedom for high-dimensional vectors of the augmented state.To relieve this problem, Xie and Zhang (2013) proposed a partitioned forecast-update scheme (PU_EnKF) that is inspired by the dual state-parameter estimation algorithm (Moradkhani et al., 2005b).In the partitioned forecast-update scheme, the parameter set of a hydrological model is partitioned into different types (N p types in total) based on their sensitivities.Each type is estimated in an individual loop by repeated forecasting and updating.Here, the parameter type maintains an aggregation connotation.A parameter type can contain only one parameter (e.g. for lumped hydrological models) or many parameters associated with the same number of computational units in distributed hydrological models.For example, the parameter CN 2 in SWAT (will be introduced in Sect.2.2) is considered as a parameter type.
At time t, the PU_EnKF is iteratively applied as follows for N p loops: I. Perform parameter evolution using Eq. ( 3) for the j th parameter type, producing a new ensemble of parameters.
II. Run the model N times following Eq.( 2) to obtain ensemble predictions for gauged and ungauged state variables.In the prediction, the j th parameter type is prescribed with a member of the ensemble produced in step I, while the others are set with the ensemble means that are estimated from previous loops at this time step and from the previous time step.
III. Compute the Kalman gain matrix using Eq. ( 7) based on the ensembles of states and parameters when observations become available at time t.
IV. Update the state vector and the j th parameter type using Eq. ( 6).
V. Compute the ensemble means of the j th parameter type.
The means are the estimates of the parameters and will be used in step II in the subsequent loops to estimate the other parameter types.
VI.Return to step I if j <N p .Otherwise, go to the next time step t + 1.The updated state vector from the loop j = N p is considered as estimates of gauged and ungauged state variables; and all estimates of parameters are also obtained.
We can see that the partitioned update scheme employs an iterative algorithm to update each parameter type at each time step -not only is one parameter considered at a time.At time t, the new estimated parameter values from previous loops are used for the model forecasting (Eq.2) in the current loop in which a target parameter type (the j th parameter type) is estimated.This iterative update is expected to push the estimates towards their optimal values.Therefore, this scheme is suitable for distributed hydrological models to estimate high-dimensional parameters.Its capability has been demonstrated using synthetic cases and it has been successfully used in a real watershed for state and parameter estimation (Xie and Zhang, 2013).In this study, we apply this scheme to improve the streamflow prediction in ungauged sites and to estimate model parameters.

Model description
The distributed hydrologic model, SWAT, is a basin-scale hydrological model developed by the USDA Agricultural Research Service (Arnold et al., 1998;Arnold and Fohrer, 2005).In implementation of SWAT, a basin is partitioned into multiple subbasins that are then divided into hydrologic response units (HRUs), which consist of unique land cover, management, and soil characteristics (Neitsch et al., 2001;Gassman et al., 2007).The HRUs are the basic computational units in which the overall hydrologic balance is simulated, including precipitation partitioning, surface runoff generation, evapotranspiration (ET), soil water and groundwater movement.
The surface runoff generation is commonly simulated using the Soil Conservation Service (SCS) model (Rallison and Miller, 1981;Ponce et al., 1996).This model has only one parameter, i.e. the curve number at moisture condition II (CN 2 ), which is also the dominant parameter in SWAT.Actual ET is formulated based on potential ET to account for evaporation from the plant canopy, transpiration, sublimation and evaporation from the soil.The soil water movement is characterized by a storage routing technique that uses the field capacity to dominate redistribution of water between layers.By infiltration or percolation, a fraction of water below the soil profile enters groundwater storage as recharge and is partitioned between shallow and deep aquifers.Base flow from the shallow aquifer is also routed to river channels.Details regarding these processes can be found in the SWAT user's manual (Neitsch et al., 2001).SWAT contains a large number of spatially varying parameter types to be prescribed before hydrologic simulation and prediction.These parameters consist of the surface roughness, soil properties, land-cover pattern and hydraulic conditions of the river channel.Although their default values can be prescribed according to lookup tables, the optimal values must be calibrated on the basis of modelling behaviour and observations.To reduce the number of calibrating parameters, a sensitivity analysis is usually required (van Griensven et al., 2006).Considerable effort has been devoted to sensitivity analysis for SWAT; several parameters are recognized as the most influential ones that dominate the model behaviour (Holvoet et al., 2005;Muleta and Nicklow, 2005;van Griensven et al., 2006).Based on these studies, seven parameters (also called parameter types) are selected and shown in Table 1.They underpin different hydrologic processes in a basin involving the surface runoff, soil water, baseflow, groundwater, evapotranspiration and channel water processes.Their ranges are determined in terms of the lookup tables (Neitsch et al., 2001) and the specific soil and land use properties of the Zhanghe River basin (Post and Jakeman, 1999).
In addition to these sensitive parameter types, ten hydrologic variables are selected to be updated in data assimilation (Table 2).They can be divided into three groups: (1) quick water storage (marked with QW in Table 2) regarding surface runoff, (2) slow water storage (marked with SW) associated with baseflow and groundwater flow and soil moisture, and (3) river channel storage (marked with CW) and flow.The first nine variables are the dynamic states that characterize water storage status in HRUs or subbasins and partially influence the diagnostic variables, i.e.ET and the water discharge (Q r ).Therefore, along with both outputs, these states should be updated to guarantee consistent model behaviour.In this study, ET is excluded from the state vector because there are no ET observations and its passive update in data assimilation does not impact other state estimations.
The SWAT model is used for this study for two main reasons.First, SWAT is a very popular distributed hydrological model to predict water, sediment, and agricultural chemical yields in large, complex watersheds (Gassman et al., 2007).An improved version of this model has been used to simulate the water movement in the Zhanghe River basin, an irrigation district with paddy rice planting (Xie and Cui, 2011).Second, we have coupled it with the EnKF-based algorithms with a few successful applications (Xie, 2013;Xie andZhang, 2010, 2013).Therefore, such a coupled SWAT-EnKF data assimilation platform is expected to be more powerful and widely used for real-time hydrological predictions.SWAT requires a significant amount of data including model input and system response data (e.g.streamflow, evapotranspiration), which seems not consistent with effort of predictions in ungauged basins.But this issue can be eased to some degree because streamflow data from just a few locations downstream (e.g. the outlet) can favour estimation for the entire basin by the data assimilation scheme used in this study.
3 Application to a real case

Study area and database
The data assimilation scheme is applied in the Zhanghe River basin in Hubei Province, China (Fig. 1).The Zhanghe drains an area of 1129 km 2 , and the elevation difference between the north and the south is more than 400 m.It has a typical subtropical climate with an annual mean temperature of 17 • C. The annual rainfall in the catchment is approximately 970 mm per year, although rainfall varies substantially from year to year depending upon the monsoon strength.This basin is actually an agricultural irrigation area and its cultivated area accounts for 59 %.Paddy rice is the primary cultivated plant, which, from May to August, requires irrigation water from the Zhanghe reservoir and thousands of local ponds.Owing to intense human activities, including cultivation, irrigation and drainage, streamflow prediction in this basin is a challenge with large uncertainties (Cai, 2007;Xie and Cui, 2011).We chose the Zhanghe River basin as a study area because there are sufficient data sets associated with weather conditions, land use and soil properties, and hydrological information.This area has been chosen for a few modelling studies  (Cai, 2007;Xie and Cui, 2011).The land use classification with resolution of 14.25 m was retrieved based on remotely sensed data (Landsat ETM+) for the years 2000 and 2001 (Fig. 1b).The land use pattern in this basin has exhibited only small changes since 2000.Therefore, we assume that the land use pattern in the period 2004-2006 is the same as in 2000-2001.The soil map with soil properties, which is used to derive model parameters, is obtained from the local agriculture department.The weather data set, including daily temperature, radiation, wind speed and relative humidity, from January 2000 to December 2006 is available from five stations distributed in and around this basin as shown in Fig. 1c.Moreover, four streamflow gauges were installed, marked as A, B, C and D for simple referencing.Gauge D is the outlet of the basin.Gauge A is located at the outlet of a small source subbasin.Because these four gauges observe the river stages and then transform the data into streamflow according to calibrated rating curves, daily streamflow data for the period 2003-2006 are available.
The Zhanghe River basin is divided into 20 subbasins based on a digital elevation model (DEM) with a resolution of 90 m (Fig. 1c).Thereafter, 98 HRUs are obtained according to land use and the soil map.With this delineation, Gauge A drains runoff from a source subbasin, Gauge B drains four, Gauge C drains ten, and Gauge D drains all the basins.

Error quantification
The success of ensemble-based data assimilation methods depends partly on ensemble generations to quantify errors from model input forcing, parameters and model structures.Moreover, quantifying observation errors is also critical to account for uncertainties from measurements and derivations.Due to the dynamics of the SWAT model, the errors/uncertainties from the input forcing, parameters and the model structure are transferred to the water storages (e.g.soil moisture and channel storages) and diagnostic variables (e.g.streamflow).Although 10 selected variables require updating in SWAT, two of them (i.e.soil moisture and streamflow) are perturbed in this study to represent the modelling uncertainties, because the other variables are internal and their uncertainties are transferred to the soil moisture and the simulated streamflow (Xie and Zhang, 2013).Precipitation as a major forcing input is also perturbed to represent the uncertainty probably derived from weather forecasting and other sources.
Perturbations to the above three variables are conducted based on zero-mean Gaussian distributions.The standard deviation (σ ) for SWAT-simulated soil moisture is set as 0.03 m 3 m −3 as suggested by Chen et al. (2011).The standard deviations for streamflow and precipitation are assumed to be proportional to their values (Clark et al., 2008), where η is the fractional factor of the standard deviation to the variable x.Thus, there are three fractional factors corresponding to the simulated streamflow (η Qm ), observed streamflow (η Qo ) and precipitation (η p ). Therefore the PU_EnKF scheme used in this study is also applicable to hydrological prediction when measured rainfall data are unavailable but could be derived from various sources (e.g.weather forecasting).With this error quantification, the three standard deviations vary with time, depending on the magnitudes of the four variables.These fractional factors should not only represent the related uncertainties in modelling and the observations but also produce ensemble streamflow predictions with reasonable ensemble spread (Clark et al., 2008).Based on the uncertainty analysis by Xie and Cui (2011), the prediction errors with the SWAT model are more than 10 % of the variables due to the irrigation and drainage practices in the Zhanghe River basin; the measurement of precipitation also has the Values of fractional factor 0.10 0.15 0.10 same level of uncertainty.Therefore, various combinations of factor values are evaluated by running the data assimilation procedure.Table 3 presents the final choice of the three fractional factors.
Note that the error quantification remains challenging for land surface data assimilation.A few newly developed approaches may be a good attempt, e.g.adaptive filtering (Crow and Reichle, 2008;Reichle et al., 2008).However, we quantify the model and observation uncertainties in terms of an experiential and practical perspective in which large storm events normally induce larger uncertainties in modelling and observations.Moreover, an overestimation of uncertainties is a better practice than underestimation to avoid the ensemble shrinkage (Crow and Van Loon, 2006;Clark et al., 2008).

Assimilation set-up and scenario design
The assimilation process is performed with three successive periods (Xie and Zhang, 2013).First, the model is prescribed with prior parameters and spun-up within the period 1 January-30 June 2003 to initialize the model states.At the end of this period, the seven parameters of the SWAT model are perturbed using the Latin hypercube method (Helton and Davis, 2003) with Gaussian distributions.The parameter means regarding the Gaussian distributions are set according the lookup table suggested in SWAT (Neitsch et al., 2001); the associated variances are constrained to ensure that random samples are within their respective physically or model-required ranges in Table 2.The uniform distribution is more intuitive than the Gaussian and often also used in sampling (Moradkhani et al., 2005b).In this study, we use the Gaussian because the lookup table provides prior estimates for the parameters.The number of parameter samples (i.e. the ensemble size) is 80.After the parameter perturbations, the second period begins (1 July-31 December 2003) to perturb the model input forcing, model states and diagnostic variables as described in Sect.3.2.The aim of this perturbation period is to quantify the uncertainties in prediction and to generate reasonable ensemble spread for subsequent data assimilation.The third period is the data assimilation period (1 January 2004-31 December 2005) in which the streamflow observations are assimilated when data are available.Given that streamflow originates primarily from either surface runoff or subsurface runoff in different periods, the variables of quick water storage (QW in Table 2) are updated only when precipitation occurs.The variables of slow water storage (SW) are updated during dry periods (no To demonstrate the improvement of streamflow prediction in ungauged locations, we only assimilate streamflow from one or two of the four gauges, and the remaining gauges, regarded as pseudo-ungauged locations, are used to validate the performance of data assimilation.Three scenarios with different combinations of data from the four gauges are designed: I. ASS_D: the observed data of streamflow from Gauge D are assimilated; Gauges A, B and C are assumed as pseudo-ungauged.This scenario is similar to a common calibration practice for which only the outlet (Gauge D) discharge data are employed to calibrate the parameters and to extrapolate streamflow of ungauged subbasins.
II. ASS_BD: the observed data of streamflow from Gauge B and D are assimilated; the other two are regarded as pseudo-ungauged subbasins.This scenario adds the data from Gauge B at the upstream in this basin based on scenario ASS_D.
III. ASS_AB: the observed data of streamflow only from Gauge A and B are assimilated.This scenario only uses the streamflow from the two gauges in the upstream subbasins.

Prediction in ungauged locations
Ensemble streamflow predictions along with parameter estimations are performed for the three scenarios.To distinguish the improvement of streamflow prediction, a controlrun scenario is conducted in which the model parameters are prescribed with the calibrated estimates from Xie and Cui (2011).The data assimilation performance is evaluated by comparing with the four series of observed streamflow.
Although the observed streamflow series still contain uncertainties, we consider them to be a benchmark because the observations are commonly assumed to be the best estimates of "real" streamflow processes.Therefore, the series of streamflow prediction errors are computed (predictions minus observations).The root mean square error (RMSE) and the mean absolute error (MAE) are used as comprehensive indices for evaluations.To quantify the ensemble spread of streamflow in data assimilation, we define a measure, the ensemble coverage index (EnCI), which is the percentage of discharge data contained in the 95 % ensemble simulation intervals.
Figure 2 shows the streamflow errors from the controlrun prediction and scenario ASS_D.The reason for the errors being presented instead of the streamflow observations is that some of the streamflow observations are so large that the difference between the cases is not notable.The controlrun simulation clearly overestimates the peak flow (in wet periods of rainfall occurrence) for the four gauges, while it underestimates the base flow in some dry periods (e.g.230th-300th time steps).This poor performance is significantly improved by assimilating the observed streamflow and by considering the uncertainties from the input forcing and model states.It may not be surprising that the Gauge D streamflow errors in ASS_D are less than those in the controlrun scenario because the observed streamflow from Gauge D is assimilated to update the prediction.For the (pseudo-)ungauged locations, the streamflow predictions of Gauges A, B and C are also more acceptable than from the controlrun scenario.At Gauge C, for example, the RMSE decreases from 3.539 m 3 s −1 to 1.912 m 3 s −1 .Moreover, there is no notable biased prediction due to the slight overestimations and underestimations for peak flow.
The EnCI for Gauge D is up to 95.72 % (see Fig. 2).This means that 95.72 % discharge data are contained in the 95 % ensemble intervals, except that some discharge data with considerable magnitudes of flood are outside of the intervals.The lowest EnCI for Gauge A (75.21 %) is partly due to the fact that Gauge A is the farthest gauge to the outlet (Gauge D, its data are assimilated).Nevertheless, all ensemble spreads for the four gauges are reasonable to trace and to contain the discharge data.
Figure 3 shows the results for Gauge C from scenarios ASS_BD and ASS_AB.Adding an observed gauge (Gauge B) at the upstream in the basin, i.e. the ASS_BD scenario, provides better streamflow predictions in the pseudoungauged subbasins than the ASS_D scenario; the RMSE drops to 1.669 m 3 s −1 and the EnCI is up to 90.28 %.If assimilating the data from the upstream locations, i.e. the ASS_AB scenario, the improvement is degraded and the predictions are only slightly better than the control-run scenario.The improvement of streamflow prediction using the PU_EnKF scheme depends on the correlation of physical processes between gauged and ungauged locations.If the two locations are very close (which means the correlation of flow processes will be strong), quite favorable data assimilation performance will be shown.In addition to Gauge C (for pseudo-ungauged locations), Gauges A, B and D have encouraging streamflow predictions due to the fact the data from these gauges are assimilated to update the predicted streamflow (not shown in Fig. 3).
Along with the updating of model states and diagnostic variables, the model parameters are also estimated.Figure 4 shows examples of real-time parameter updating from the ASS_D scenario.After about 130 time steps, the ensemble trajectories are nearly stable with slow variations which are probably induced by the changes of land surface and river channel conditions for runoff generation and routing (Liu et al., 2008;Troch et al., 2013).At every time step in data assimilation, the parameter samples can be approximated with Gaussian distributions and they are constrained within the prior ranges (Min-Max, see Table 1) as shown in the histograms in Fig. 4.This property is favourable for parameter estimation with ensemble-based data assimilation.The uncertainties of parameter estimates at every time step are represented using the ensemble spread (EnSp), which is computed based on sample variances (see the illustration in caption of Fig. 5).At the beginning of the data assimilation, the parameters have broad ensemble spreads.The spreads quickly shrink after 100 time steps with the evolution of the streamflow assimilation, and remain stable after 400 time steps.Therefore, the estimate uncertainties of the parameters decrease with the data assimilation and state updating.Moreover, the relative stabilities of ensemble trajectories (Fig. 4) and the ensemble spreads (Fig. 5) imply an attractive potential that it is possible to use short-term data to retrieve optimal estimates of parameters.
Even though the three scenarios provide different parameter estimates due to the assimilation of different observations, encouraging properties of parameter estimations are achieved in the three scenarios.It is not sure so far whether the parameter estimates converge to their appropriate values in this real-word application, so the parameter estimates require a further validation to evaluate the effectiveness of the PU_EnKF scheme.

Validation for parameter estimates
It is difficult to directly validate the parameter estimates using measurements because the SWAT model is a conceptual hydrological model and most parameters do not have physical meanings.Only a few parameters (e.g. the SOL_AWC in Table 1) can be measured at local sites; those parameters regarding HRUs, subbasins and river channels remain difficult to obtain by sampling experiments.We perform single-run predictions using the parameter estimates from the three scenarios and evaluate the predicted streamflow against observed streamflow.This is a commonly used strategy to   Figure 6 shows the streamflow prediction errors from the four scenarios.Only the results of Gauges C and D are shown because they are located at the downstream locations in the Zhanghe River basin.The three scenarios using prescribed parameters with estimates from data assimilation achieve better predictions for the two gauges than the control-run scenario.The RMSE of Gauge D from the ASS_D scenario decreases from 5.550 m 3 s −1 to 2.324 m 3 s −1 .Moreover, the ASS_BD scenario provides the best predictions among the four scenarios.All of these improvements are attributable to the appropriate parameter estimates from the data assimilation.The ASS_BD scenario renders the most reasonable parameter estimates.Comparably, the parameter estimates from ASS_D are also satisfactory for streamflow predictions, while the estimates from the ASS_AB scenario lead to slight improvements for streamflow predictions.Therefore, the parameter estimation performance of the three scenarios is consistent with the prediction of diagnostic variables (i.e. the water discharge) as illustrated in Sect.3.4.The assimilated observations from downstream, especially the outlet of the basin, have more important roles than those from upstream for parameter estimation and streamflow predictions in ungauged subbasins.

Conclusions
We present an application of PU_EnKF for improving streamflow predictions at ungauged locations.This scheme features real-time updating and simultaneous state-parameter estimation, considering modelling and observing uncertainties.Moreover, the scheme constrains the predictions by the physical rainfall-runoff processes that are defined in the distributed hydrological model (i.e. the SWAT model), and it accounts for the correlations of states and parameters between gauged and ungauged subbasins.The correlations are represented by the covariance matrix in the Kalman gain.With the constraint and the correlation representation, the observed information is successfully transferred to ungauged locations and thereby improves streamflow prediction.
The real-word application case suggests that the PU_EnKF scheme performs better than the control-run simulation (with calibrated parameters) for streamflow predictions at gauged and ungauged locations.Although only the outlet-gauged data are assimilated, the streamflow predictions at ungauged sites are still acceptable, since they contain convergent flow information from all subbasins due to runoff routing.Generally, the downstream data (especially the data from the outlet) have important roles to reflect the runoff generation for the entire basin.This data assimilation scheme provides reasonable estimates of model parameters for all computational units (i.e.subbasins and HRUs), including both gauged and ungauged sites, as validated by the conventional single-run simulation.Moreover, the parameter estimates approach nearly stable levels after a small number of time steps (130 steps in this study).The parameter estimates show slow variations to trace parameter travels, which indicates the PU_EnKF has a potential advantage of identifying the changes in underlying land surface (e.g., the land use and land cover changes).
Although favourable performance to improve streamflow predictions is obtained using the EnKF-based scheme, the runoff routing is neglected within the PU_EnKF assimilation set-up because the travel time of generated runoff is less than 1 day in the Zhanghe River watershed.In fact, the time lag of runoff routing is an important factor for short-time (e.g. the hourly step) flood forecasting (Li et al., 2013;Pan and Wood, 2013).Moreover, this scheme is intent on PUB for the nested basins in which the correlations of states and parameters between neighbouring subbasins can be constructed.For separate basins in the same climatic regions and land surface conditions, assimilating other sources of data (e.g. the remotely sensed soil moisture and bright temperature) is expected to improve the predictions of hydrological variables (Troch et al., 2003).Nevertheless, this study provides an encouraging application for PUB by assimilating streamflow, which is generally regarded as quality observations compared with the remotely sensed data.There are optional methods to address PUB, e.g. the Particle-DREAM by Vrugt et al. (2013).It will be an encouraging attempt to compare these methods with distributed hydrological models for hydrological diagnosis and predictions.

Figure 1 .
Figure 1.Zhanghe River basin in China (a), the land use (b) and subbasin distribution with DEM (c).

Figure 2 .
Figure 2. Streamflow prediction errors from the control-run simulation (left column) and the data assimilation of scenario ASS_D (right column), i.e.only the observed streamflow from Gauge D (outlet) is assimilated to update model states and parameters.

Figure 3 .
Figure 3. Streamflow prediction errors from scenarios ASS_BD and ASS_AB.Only the results for Gauge C are shown because Gauge C is at the outlet of a pseudo-ungauged subbasin in both scenarios.

Figure 4 .
Figure 4. Estimations of two typical parameters (CN 2 and CH_K) from the ASS_D scenario.The histograms in each plot, fitted with the Gaussian distribution function, represent the ensemble distribution at three time steps.

Figure 5 .
Figure 5. Ensemble spreads (EnSp) of the seven parameters listed in Table 1: EnSp= 1 Nu Nui=1 VAR En (i), where "Nu" is the number of HRUs or subbasins and VAR En (i) denotes the ensemble variance at each HRU or subbasin with respect to each parameter.

Figure 6 .
Figure 6.Streamflow predictions using four scenarios of different parameter sets.Only results of Gauges C and D are shown.

Table 1 .
Model parameters to be estimated in data assimilation.
*The hydrologic variables are with respect to the scales to reflect the related hydrologic processes.

Table 2 .
Dynamic hydrologic states and outputs to be updated in data assimilation.Q surfstor Amount of surface runoff stored or lagged (mm) HRU QW Q latstor Amount of lateral flow stored or lagged (mm) HRU QW Q pregw Amount of groundwater flow into the main channel (mm) HRU QW W sol Amount of water stored in the soil layer for each HRU (mm) HRU ×N lay SW This column indicates the scale at which each variable is simulated.N lay is the number of soil layers (N lay = 4 for this study), and HRU ×N lay means the soil profile of each HRU is partitioned into N lay layers. 2 This column denotes water storage condition for each variable: QW, quick water storage; SW, slow water storage; and CW, river channel storage.