Assimilating in situ and radar altimetry data into a large-scale hydrologic-hydrodynamic model for streamflow forecast in the Amazon

In this work we introduce and evaluate a data assimilation framework for gauged and radar altimetry-based discharge and water levels applied to a large scale hydrologic-hydrodynamic model for stream ﬂow forecasts over the Amazon River basin. We used the process-based hydrological model called MGB-IPH coupled with a river hydro-5 dynamic module using a storage model for ﬂoodplains. The Ensemble Kalman Filter technique was used to assimilate information from hundreds of gauging and altimetry stations based on ENVISAT satellite data. Model state variables errors were generated by corrupting precipitation forcing, considering log-normally distributed, time and spatially correlated errors. The EnKF performed well when assimilating in situ discharge, 10 by improving model estimates at the assimilation sites and also transferring information to ungauged rivers reaches. Altimetry data assimilation improves results at a daily basis in terms of water levels and discharges with minor degree, even though radar altimetry data has a low temporal resolution. Sensitivity tests highlighted the importance of the magnitude of the precipitation errors and that of their spatial correlation, 15 while temporal correlation showed to be dispensable. The deterioration of model performance at some unmonitored reaches indicates the need for proper characterization of model errors and spatial localization techniques for hydrological applications. Finally, we evaluated stream ﬂow forecasts for the Amazon basin based on initial conditions produced by the data assimilation scheme and using the ensemble stream ﬂow predic-20 tion approach where the model is forced by past meteorological forcings. The resulting forecasts agreed well with the observations and maintained meaningful skill at large rivers even for long lead times, e.g. > 90 days at the Solim˜oes/Amazon main stem. Results encourage the potential of hydrological forecasts at large rivers and/or poorly monitored regions by combining models and remote sensing information.


Introduction
Land surface waters play an important role in global water cycle and earth system, regulating freshwater discharge from land into oceans (Oki and Kanae, 2006) and also land-atmosphere exchanges of water, energy (Krinner, 2003;Decharme et al., 2012) and gases such as methane (Gedney et al., 2004).Moreover, it affects directly society that uses it for drinking water and also transportation of people and goods, agriculture and energy production from hydropower.More specific to the Amazon basin, important extreme hydrological events have occurred recently, for instance, the 2009 and 2012 floods and the 1996, 2005 and 2010 droughts (Chen et al., 2010;Tomasella et al. 2010;Marengo et al., 2008;Espinoza et al., 2011;Marengo et al., 2011).These events caused several impacts on local population that strongly depends on the rivers and is very vulnerable to floods since most settlements lie along the rivers.
In situ measurements of river stage and discharge at stream gauges are the most conventional alternative for monitoring surface waters, although observation networks are rather sparse at several regions such as the Amazon River basin.Alternatively, radar altimetry techniques have been developed in past years to monitor water levels (e.g.Santos da Silva et al., 2010;Alsdorf et al., 2007) or discharges using rating curves (e.g.Leon et al., 2006;Papa et al., 2010a;Getirana and Peters-Lidard, 2013).If compared to in situ gauges in remote regions, these satellite instruments can provide observations with much better spatial resolution, but with worse temporal sampling.
Moreover, the forthcoming Surface Water and Ocean Topography (SWOT) mission (Durand et al., 2010a) is designed to provide high resolution images of inland water surface elevation, including rivers, lakes, wetlands and reservoirs, using a swath mapping radar altimeter with high frequency repeat orbit.Additionally, it may also be possible to derive discharge estimates from SWOT data by using specially developed algorithms (e.g.Durand et al., 2010b).
In contrast, there are several efforts on hydrological modeling to simulate processes as river and floodplain dynamics in large river basins such as the Amazon (Paiva et al., Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | 2012, 2013;Yamazaki et al., 2011;Getirana et al., 2012;Decharme et al., 2012;Coe et al., 2008;Wilson et al., 2007;Trigg et al., 2009).These models can potentially provide detailed information on surface waters, both spatially and temporally, but such estimates are somehow imperfect due to uncertainty in model structure, parameters and forcing data (Liu and Gupta, 2007).
Data assimilation (DA) methods are an alternative to optimally merge uncertain model predictions with both in situ and the newly remote sensing observations of surface waters.The aim of DA techniques is to "produce physically consistent representations or estimates of the dynamical behaviour of a system by merging the information present in imperfect models and uncertain data in an optimal way to achieve uncertainty quantification and reduction" (Liu and Gupta, 2007).Such methods can also be used to estimate balanced initial states of hydrological models for forecasting the aforementioned extreme events.There are already some hydrological regional/global forecast systems founded on physically-based hydrological models (e.g.Wood et al., 2002;Thielen et al., 2009;Alfieri et al., 2012), and also several physical modeling experiments in the Amazon basin, as previously mentioned.However, current attempts for developing hydrological forecasts in this basin are mostly based on statistical methods (e.g.Uvo and Grahan, 1998;Uvo et al., 2000).Furthermore, Paiva et al. (2012b) showed that, for lead times up to 3 months, uncertainty of initial conditions plays a major role for discharge predictability on main Amazonian Rivers, if compared to the importance of precipitation forcing, suggesting the importance of DA techniques for streamflow forecasts in this region.
Research on data assimilation applied to hydrology has increased in past years with various applications utilizing Kalman filters (e.g. the Ensemble Kalman Filter -EnKF, developed by Evensen, 2003), particle filters or variational methods, as extensively reviewed in Liu and Gupta (2007), Reichle (2008) and Liu et al. (2012).These applications include a wide range of observations, both in situ and remotely sensed, data assimilation methods and models representing different hydrological processes, at different spatial scales and with several objectives, such as: the assimilation of snow  (Andreadis and Lettenmaier, 2006) and soil moisture (Reichle et al., 2002) data into land surface models using the EnKF; assimilation of in situ water level measurements into a small scale 1-D hydrodynamic model for flood forecast using Kalman filtering methods (Neal et al., 2007;Ricci et al., 2011); assimilation of synthetic SWOT data into hydrodynamic models at restricted areas using the EnKF and some variations (Biancamaria et al., 2011;Andreadis et al., 2007;Durand et al., 2008); assimilation of discharge data into distributed hydrological models (Clark et al., 2008;McMillan et al., 2013;Lee et al., 2012;Thirel et al., 2010;Rakovec et al., 2012) using the EnKF or variational methods; simultaneous assimilation of soil moisture and discharge data into a distributed hydrological model using variational DA (Lee et al., 2011); assimilation of radar altimetry data of reservoir water levels using the EnKF (Pereira-Cardenal et al., 2011); development of a modelling platform (Land Information System -LIS) to merge multiple in situ and remotely sensed observations with land surface models (Kumar et al., 2008); merging water levels information derived from a satellite Synthetic Aperture Radar (SAR) image and digital terrain model (DTM) with a 1-D hydrodynamic model for estimating river discharge (Neal et al., 2009); among others.Although there is an extensive bibliography on hydrological data assimilation, the current state of the art regional/global hydrological prediction systems (e.g.Thielen et al., 2009;Alfieri et al., 2012) still do not incorporate advanced data assimilation systems for updating model initial states.Also, the assimilation of discharge and water levels from in situ and remotely sensed observations into regional/global hydrologic-hydrodynamic models is still uncommon.
In this paper, we present the development and evaluation of a data assimilation framework for both gauged and radar altimetry-based discharge and water levels into a large scale hydrologic-hydrodynamic model of the Amazon River basin using the EnKF.We also explore the usefulness of such system to provide streamflow forecasts when forced by past remotely sensed precipitation data and based mostly on model initial conditions.This paper is in the context of recent developments of techniques for integrating information from hydrological models with newly remotely sensed data such as the forthcoming SWOT mission and also in the context of regional/global hydrological forecast systems including large, poorly gauged river basins.Through our experimental results, we explore questions such as: is an EnKF-based DA scheme feasible for assimilating discharge and water level data into large scale hydrologichydrodynamic models?Is it able to improve discharge and water level estimates at sites where data were assimilated and also at ungauged rivers?Does the assimilation of radar altimetry data also improve model estimates at large river basins, considering that it has lower temporal resolution and accuracy if compared to gauged in situ data?Would it be possible to provide accurate streamflow forecasts at large basins such as the Amazon using a large scale hydrologic-hydrodynamic model based mostly on the initial hydrological states gathered by the DA scheme?At which spatial and temporal scales?

The hydrologic-hydrodynamic model
We used the MGB-IPH model (Collischonn et al., 2007), which is a large scale, distributed and process-based hydrological model with a hydrodynamic module described in Paiva et al. (2011).It simulates surface energy and water balance and also discharge, water level and flood inundation on a complex river network.Vertical hydrological processes include soil water budget using a bucket model, energy budget and evapotranspiration using the Penman Monteith approach, and also surface, subsurface and groundwater flow generation, among others.The flow generated within each catchment is routed to the stream network using a linear reservoir type model.River flow routing is performed using a combination of either a Muskingum-Cunge (MC) method or a hydrodynamic model (HD).
The hydrodynamic model of MGB-IPH (Paiva et al., 2011)  model assuming that the floodplains act only as storage areas.River-floodplain parameters (river width, bottom levels, roughness coefficient, floodplain bathymetry) are estimated using GIS-based algorithms from the Shuttle Radar Topography Mission (SRTM) Digital Elevation Model (DEM) (Farr et al., 2007) and using geomorphological relations.

The Ensemble Kalman Filter
The goal of data assimilation is to combine the uncertain and complementary information from measurements and simulation models into an optimal estimate of the hydrological fields of interest, providing a general framework for dealing with uncertainty from measurements and also input, output and model structure (Reichle, 2008;Liu and Gupta, 2007;Liu et al. 2012;Vrugt et al., 2005).
A great part of the hydrological applications of data assimilation methods uses schemes based on the Kalman filter (Kalman, 1960), specially the Ensemble Kalman Filter (Evensen, 2003(Evensen, , 2009)), which is also used in this study and is briefly described below.The model representing the dynamics of the simulated system can be represented in a discrete form by the process equation: where x represents the state variables and, u and θ represent model forcings and parameters, respectively, M is the nonlinear model operator that relates model states from time interval t k to t k+1 = t k + ∆t, and q k represents errors in model structure, parameter, forcings and antecedent states.In this study, x is a vector composed by all MGB-IPH state variables, i.e. soil moisture, storage and discharge from surface, subsurface and groundwater reservoirs, soil temperature, canopy storage and river discharge and water level.The measurement equation is defined by: where y is observation vector, ε is the vector of observation errors and H is the observation operator which relates the model state variables vector x to the observations vector y.In our case, observations are river discharges or water levels at selected sites.
In the context of forecast systems and sequential data assimilation methods, at each time interval, the model is integrated in time using Eq. ( 1) to provide a short-term forecast (or background) x f k+1 and whenever an observation is available, the forecast errors are computed as (y k+1 − H(x f k+1 )).Therefore, the goal of data assimilation is to obtain an optimal estimate of model state variables x a (called analysis) given model and observation errors.In the case of the original Kalman Filter (Kalman, 1960), the DA problem is solved using a linear estimator assuming that (i) model and observation operators (M and H) are linear; (ii) observation errors are unbiased and both temporally and spatially uncorrelated; (iii) model errors are unbiased and temporally uncorrelated; and (iv) there is no correlation between model and observation errors.Consequently, an unbiased and minimum variance estimate of model states is obtained by: where K is the Kalman gain, P is the covariance of model errors q, and R is the covariance of measurement errors ε.The Kalman filter also provides the maximum likelihood solution of the DA problem when model and measurement errors are assumed to be also Gaussian.However, the applicability of the KF is limited since most hydrological systems exhibit nonlinear dynamics (Liu and Gupta, 2007) and the assumptions about model errors (e.g.Gaussian, unbiased, among others) are not always valid.Moreover, both the original KF and also its non-linear version, the Extended KF (EKF), have additional drawbacks when applied in large and complex systems with lots of state variables (e.g.distributed hydrological models) due to extra programming and heavy computational requirements associated with the storage and forward integration of the error covariance matrix P (Vrugt et al., 2005).Evensen (2003) presented the Ensemble Kalman Filter (EnKF), which is a stochastic or Monte Carlo alternative for the deterministic EKF (Evensen, 2009).In this method, observations and model states are perturbed using a priori-known errors and by means of the model operator M, the algorithm generates an ensemble of model trajectories from which the time evolution of model errors and error covariance matrix can be sampled: Each ensemble member is then updated using the same analysis equation from the original KF (Eqs. 3 and 4).Alternatively, efficient computational implementations of the EnKF are presented in Evensen (2003Evensen ( , 2004Evensen ( and 2009) ) (see http://enkf.nersc.no/for Fortran codes) where the explicit computation and storage of P are not required.We used the square root scheme presented in Evensen (2004Evensen ( , 2009) ) where the perturbation of measurements is not performed, in order to reduce sampling errors.
For linear systems and with large ensemble sizes, the EnKF provides the same solution as the KF method.However, it is noteworthy that it does not fully take into account non-Gaussian errors nor solve the Bayesian update equation for non-Gaussian probability distribution functions.Still, it is a computational efficient analysis scheme for nonlinear models that provides a satisfactory solution, although it is suboptimal, that somehow lies between a linear Gaussian update and a full Bayesian computation (Evensen, 2009).

Uncertainty in precipitation forcing
We perturbed model states variables by adding a noise in precipitation forcing, considering (i) that this is the most uncertain model input (Liu et al., 2012) and possibly the most important source of model uncertainty and (ii) that this is a proper method to generate physically coherent model errors.A similar approach performed satisfactorily in other hydrological applications of DA methods such as Andreadis and Lettenmaier (2006) and Biancamaria et al. (2011).Precipitation values were corrupted using a lognormally distribution as presented by Nijjsen and Lettenmaier (2004) and also applied by Andreadis and Lettenmaier (2006): where P c (mm ∆ t −1 ) is the perturbed daily precipitation, P (mm ∆ t −1 ) is the unperturbed daily precipitation, E is the relative error (%), β is the relative bias and s ∼ N(0, 1) is a normally distributed and spatially correlated random variable with zero mean and unit variance.Spatially correlated pseudo random fields w were generated by means of the algorithm based on the two dimensional Fourier transform presented in Evensen (2003) (see http://enkf.nersc.no/for Fortran codes), having zero mean, unit variance and isotropic covariance function decreasing to the e −1 value at the distance τ x called spatial decorrelation length.At each spatial location, temporal correlation was also considered using the following equation for simulating the time evolution of errors (Evensen, 2003): where k is the time interval, s k is a sequence of time correlated errors with zero mean and unit variance (input for Eq. 6) and α determines the time decorrelation of the stochastic forcing, e.g.α = 0 generates a sequence of white noise while α = 1 removes the stochastic forcing.The parameter α is determined by: where τ t is the temporal decorrelation length, that determines that s decreases by the ratio e −1 after a time period t = τ t if the stochastic term w is excluded.

Measurement errors
Water level (z) and discharge (Q) observations errors were modeled using the following relations: where z c (m) and Q c (m 3 s −1 ) are the corrupted values of z(m) and Q (m 3 s −1 ), ε z (m) and ε Q (m 3 s −1 ) are the normally distributed errors with parameters σ z (m) and σ Q (%) of z and Q, respectively.The formulation of discharge errors allows representing larger uncertainties for high stage levels than for low flows due to uncertainties in discharge rating curves, as pointed out by Clark et al. (2008).Alternatively, simulated and observed discharges were also transformed into the log space before the assimilation, following Clark et al. (2008).In this case, observation errors were modeled by: where now ε ′ Q is a log-normally distributed error with unit mean and standard deviation σ ′ Q (%), similar to Eq. ( 10).At log space, standard deviation is given by 3 Experimental design

The Amazon River basin
The study area is the Amazon River basin, known as the largest hydrological system of the world (∼ 6 million km 2 of surface area) and contributing with ∼ 15 % to the total fresh water released into the oceans.The Amazon basin is characterized by extensive seasonally flooded areas (Hess et al., 2003;Papa et al., 2010b;Melack and Hess, 2010), which store and release large amounts of water from the rivers and consequently attenuate and delay flood waves in several days or months (Paiva et al., 2012a(Paiva et al., , 2013;;Yamazaki et al., 2011).Also, complex river hydraulics are present, where the low river slopes cause backwater effects that control part of river dynamics (Meade, 1991;Trigg et al., 2009;Tomasella et al., 2010;Paiva et al., 2012aPaiva et al., , 2013)).Additionally, this region presents high precipitation rates (average ∼ 2200 mm yr −1 ) with high spatial variability and contrasting rainfall regimes in the northern (rainfall peak at JJA) and southern (rainfall peak at DJF) parts of the basin, with more defined wet and dry seasons occurring in southern and eastern regions (Espinoza et al., 2009).

Model implementation
We used a MGB-IPH implementation on the Amazon basin developed by Paiva et al. (2013), as briefly described below.The model was forced using meteorological data obtained from the CRU CL 2.0 dataset (New et al., 2002) and remotely sensed precipitation estimates from the TRMM 3B42 v6 product (Huffman et al., 2007), with spatial resolution of 0.25 • × 0.25 • and daily time step for a period spanning 12 yr (1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009).The model parameters related to soil water budget were calibrated using daily discharge data from stream gauges (see next section for description of gauged data).Then, the model was validated against daily discharge and water level data from stream gauge stations, water levels derived from ENVISAT satellite altimetry data (Santos da Silva et al., 2010) (212 sites with 35-day repeat orbit), monthly Terrestrial Water Storage from GRACE mission (Frappart et al. 2010(Frappart et al. , 2011b) )  Since this study aimed at applications of data assimilation to hydrological forecasting, we also used a real time precipitation product to force the MGB-IPH model.We choose to use the TRMM Merge product (Rozante et al., 2010), which is a near to real time precipitation estimate based on TRMM 3B42RT (Huffman et al., 2007) merged with data from in situ gauges and provided by the Brazilian center for weather forecasts and climate studies CPTEC (Centro de Previs ão do Tempo e Estudos Clim áticos), a division of the Brazilian National Institute for Space Research INPE (Instituto Nacional de Pesquisas Espaciais).

Discharge and water level observations
We evaluated the assimilation of three types of data: (1) in situ discharge observations; (2) remotely sensed water levels derived from the ENVISAT radar altimeter; and (3) remotely sensed discharge estimates derived from radar altimetry water levels and rating curves.
In situ daily discharge from 109 stream gauges were provided by the Brazilian agency for water resources ANA (Ag ência Nacional das Águas), the Peruvian and Bolivian national meteorology and hydrology services SENAMHI (Servicio Nacional de Meteorología e Hidrología) and the French ORE-HYBAM program (Hydrologie, Biogeochimie and Geodynamique du Bassin Amazonien, http://www.ore-hybam.org).We also used stage data from 66 ANA gauge stations, but only for validation purposes.
Remotely sensed water levels were obtained from the ENVISAT satellite altimeter.
The ENVISAT satellite has a 35-day repeat orbit and an  2010) result in ∼ 10 to 40 cm water level accuracy.Due to differences in water levels datum reference, the comparisons between simulated and observed water levels were performed in terms of anomalies, i.e. after removing the long-term average.
Altimetry-based discharge data was developed by Getirana and Peters-Lidard (2013) for the Amazon basin, following the methodology first presented by Leon et al. (2006) in the Negro River sub-basin.This dataset was constructed using a rating-curve-based methodology deriving water discharge from ENVISAT altimetry data at 475 altimetric stations (AS).The stage-discharge relations at each AS were built based on satellite altimetry and outputs from a global flow routing (GFR) scheme (Getirana et al., 2012).A second experiment was performed in this study using observed discharges at gauge stations to force the GFR scheme at downstream reaches.Validation of the methodology against observed discharges at 90 sites showed a mean relative error of 27 % for the experiment using in situ discharge within the GFR scheme.We assimilated data only from the 287 ASs located downstream of a gauging station where results were improved in the second experiment.

Parameters of the DA scheme
The first sensitivity experiments used the following standard parameters of the DA scheme.Ensemble size of the EnKF was set as N = 200.Precipitation fields were corrupted considering the following error parameters: precipitation relative error E = 25 %, and precipitation relative bias β = 1.0 following Andreadis and Lettenmaier (2006); temporal decorrelation length of precipitation errors τ t = 10 days; and spatial decorrelation length of precipitation errors τ x = 1.0 • , similarly to Andreadis and Lettenmaier (2006) and Clark et al. (2008).The parameter of water level measurements error was set as σ z = 0.20 m, based on the accuracy of ENVISAT estimates provided by Santos da Silva et al. (2010).We computed the mean relative error between in situ discharge measurements and values provided by rating curves at 87 gauging stations from the ANA database as a surrogate of the discharge error parameter σ Q .The median value of all stations was 13 %, while Clark et al. (2008) used 10 % in its DA experiments.
Therefore, we choose to also use σ Q = 10 % for simplicity.We used σ Q = 27 % for assimilation of satellite based discharge data, based on the error value found in Getirana and Peters-Lidard (2013).

Data assimilation experiments
We performed three data assimilation experiments, namely: (i) in situ discharge assimilation (Exp 1) (ii) radar altimetry assimilation (Exp 2) and (iii) assimilation of discharge series based on satellite altimetry (Exp 3).
In the first experiment, we tested: (Exp.1a) the assimilation of discharge from almost all gauge stations (80 %) using a few of them for validation (20 %); (Exp.1b) the assimilation of only 12 stations (∼ 10 %) located at some of the major tributaries to emulate the situation of using only telemetric stream gauges for real time applications; (Exp.1c) the assimilation of discharge from almost all gauge stations, similar to (Exp.1a), but without transforming discharge into the log space (Sect.2.4).Moreover, we explored the sensibility of the DA scheme to some of its parameters, namely the ensemble size N, precipitation relative error E and temporal and spatial decorrelation lengths of precipitation errors τ t , and τ x .
The second experiment (Exp.2) evaluated the assimilation of ENVISAT radar altimetry water level anomalies.Stage data from in situ gauges were used for model verification.Simulations were also compared in terms of discharge using in situ data to evaluate the impact of water level assimilation in discharge estimates.
In the third experiment (Exp.3), we assessed the assimilation of discharge derived from radar altimetry water level.Discharge data from stream gauges were used for verification.
In all cases, simulations started in 1998 and ran to 2002 for model spin-up.The year of 2003 was used for the spin-up of the DA scheme, where no update was performed in the first months allowing the system to develop a coherent correlation structure, following Andreadis and Lettenmaier (2006).Results were evaluated for the two year period 2004-2005 using the following performance statistics: (i) the Nash-Suttcliffe coefficient E NS ranging from −∞ to 1 (optimum) and (ii) changes in root mean square error ∆rms ranging from −100 % (optimum) to ∞.

Prospects of streamflow forecasting
Hindcast streamflow forecasts were generated using an ensemble streamflow prediction (ESP) approach (Day, 1985), as described below.The model uses an estimate of initial conditions derived from the DA scheme and runs forced by an ensemble of observed meteorological data from past years.An estimate of initial conditions is computed during the spin-up period using a hydrological model driven by observed meteorological forcings, updated using data assimilation of observations up to the time of forecast (e.g.forecast starts with model states from 1 June 2010).Then, an ensemble forecast is obtained using observed meteorological data resampled from past years (e.g.meteorological data from 1 June to 1 September of years 1998, 1999, ..., 2009).
Precipitation from TRMM Merge was used during spin-up period, while during forecast the model was forced with TRMM 3B42 data for the period spanning 12 yr (1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009) and, consequently, the forecast ensemble had 12 members.The DA scheme used the configuration from Exp. 1b where in situ discharge data were assimilated to update model states before starting a forecast.ESP runs generated decadal forecasts up to 90 days lead time and starting at every 1st, 10th and 20th day of the month for the two year period of 2004-2005.
For simplicity reasons, forecasts were evaluated only by deterministic means by averaging ensemble values into a single forecast.We used the skill score SS cli which compares the performance of the model forecasts with a control forecast based on climatology (Wilks, 2006): where t is the time interval, Q obs is daily discharge observed at stream gauge stations,

In situ discharge assimilation
We start our analysis evaluating the sensibility of the DA scheme performance in terms of ∆rms to some of its parameters, as presented in Fig. 1.The objective of such examination is to verify which parameters are the most important ones and if the DA performance is improved by using values of these parameters that are different from the first guess ones based on previous studies (see Sect. 3.4).The configuration of Exp.1a was used, where in situ discharge data were assimilated.Results were evaluated in terms of mean changes in root mean squared error (∆rms) between observed and simulated discharges, computed for two samples, the first including stream gauges used for data assimilation and the latter only the validation ones.Larger decreases in the rms error indicate better performance of the DA scheme.
According to the analysis, the DA scheme strongly depends on the ensemble size N. Small N values produce small improvements in discharge results and larger values enhance the DA performance (smaller ∆rms values), although the improvement rate is small for N values larger than 150 members.Such behaviour is possibly due to numerical reasons, since a larger N enable a better sampling of model covariance errors from the ensemble, as discussed by Evensen (2009).The DA scheme is also very sensitive to precipitation relative error E and increasing E values improves DA performance.However, if E is larger than 50 %, ∆rms increases in validation sites causing worse results (see Fig. 1).Possibly, larger precipitation errors cause larger model uncertainty and consequently the DA scheme gives more weight to observations, but it starts to degrade model results at different locations after some point.A moderate dependence to the τ x parameter was found and spatial correlation of precipitation errors showed to be of importance, since the performance degrades for smaller decorrelation lengths.The best results were obtained for 1.5 • for both the assimilation and validation samples.Finally, a weak sensibility to the τ t parameter was found, which indicates that considering temporal correlation in precipitation errors is not as important as spatial correlation.
Based on the sensitivity tests, we used the following new parameter values for the further experiments: N = 200 (unchanged), E = 50 %, τ x = 1.5 • and τ t = 10 days (unchanged).However, it is noteworthy that these parameter values related to precipitation errors, although providing better results for data assimilation, may not realistically represent errors in the TRMM Merge dataset or the spatially variable satellite precipitation errors presented in Tian and Peters-Lidard (2010).That is possibly because we considered that model uncertainty comes from precipitation errors and neglected other sources such as parameter and model structural errors (Liu and Gupta, 2007).Therefore, and since the first guess values were not fully justified in the previous studies (Andreadis and Lettenmaier, 2006;Clark et al., 2008), we preferred to use the parameter values where the DA scheme performs better.
We first evaluate results from the Exp.1a.The DA scheme improves results by decreasing model errors in almost all stream gauges (blue sites in Fig. 2a), including both assimilation and validation sites.On average, E NS values increase from 0.71 to 0.94 and the rms error decreases by 49 % (Table 1).For example, at an assimilation site located on the Negro River (Fig. 3a), when the EnKF is used, the discharge estimates are much closer to observations if compared with the open-loop simulation.The E NS index increases from 0.62 to 0.91 and the rms error decreases by 51 %.Similarly, results also improve at validation sites, although with a smaller degree, and the E NS index increases from 0.60 to 0.73, with a reduction in rms error of −16 % (Table 1), as illustrated at a validation site located at upper Juru á River basin (Fig. 3b).Such results demonstrate that the DA scheme improves model discharge estimates, not only at sites where data were assimilated but possibly at ungauged rivers reaches as well.
In Exp.1b, results improve at assimilation sites −E NS increases from 0.89 to 0.98 and the rms error decreases by 50 % (Table 1).However, since data from only a few gauges were assimilated, there is no important improvement (∆rms = −3 %) if all validation sites are examined together.As expected, according to Fig. 2b the DA scheme improves discharge estimates mostly at large rivers (e.g.Fig. 3c), where E NS increases from 0.79 to 0.87 while ∆rms equals −23 %.But at smaller rivers, in most cases the DA scheme has minor effect on simulated discharges (green squares at Fig. 2b) or in some cases it degrades results.
In previous studies conducted over smaller basins (e.g. in Clark et al., 2008; and in others summarised by Lee et al., 2012), the attempt to transfer information to neighbour or upstream ungauged river reaches was unsuccessful and corrupted model results, while in our case (Exp.1b) the DA scheme degraded model outputs mostly at smaller basins and improved results at larger rivers.Such behaviour possibly happens because the state estimation in distributed hydrological models is subject to overfitting due to the large dimensionality of the model state space, and consequently, when limited data is available, the data assimilation may update state variables at some lumped fashion such as the sub-basin scale, as explained by Lee et al. (2012).
Finally, we compare the use (Exp.1a) or not (Exp.1c) of the transformation of discharge values into the log space before data assimilation.The performance of the DA scheme degrades if the log transformation is not used, and in this case ∆rms increases to −29 % and −10 % for the assimilation and validation samples, respectively, instead of the −49 % and −16 % values obtained in the Exp.1a.Clark et al. (2008) argue that the EnKF with log transformation performs better because relationships between streamflow and model states are non-linear and state updates are exceptionally large when differences between model and observed values are high.However, the worst performance was observed mostly at smaller river reaches (see Fig. 2c) as illustrated in Fig. 3c.Also, DA performs better at gauging stations in large rivers and ∆rms increases from −34 % (Exp.1a) to −40 % (Exp.1c).Apparently, when the log transformation is not used, the DA scheme gives more weight to large discharge values (∼ 10 3 to ∼ 10 5 m 3 s −1 ) at large rivers while observations at the smaller ones (<∼ 10 3 m 3 s −1 ) are not fully taken into account.These results indicate the importance of using the log space transformation also to deal with very different discharge magnitudes, including the ones arising from different spatial scales but also concerning to flood and drought flows.

Radar altimetry data assimilation
In this section, the assimilation of water levels derived from ENVISAT altimetry is evaluated (Exp.2).Stage and discharge data observed at in situ gauging stations were used for validation purposes.
The DA scheme improves water level simulations at altimetry stations used for assimilation (Fig. 4a), as illustrated in Fig. 5a.On average, rms decreases by 56 % and E NS values increase from 0.66 to 0.96 (Table 2).Simulated water level accuracy also increases when compared to in situ stage data (∆rms = −13 %).However, the improvement is more evident if only gauging stations located at rivers where altimetry data were assimilated (Figs.4b and 5b).In this case, mean ∆rms equals −43 % and E NS changes from 0.75 to 0.94 (Table 2), similar to what was obtained at altimetry stations.At other sites, the DA scheme has a minor effect on simulated water levels or even degrades results in some cases.Similar results were found for the in situ discharge validation sample (Fig. 4c).Assimilating water level data improves discharge estimates (∆rms = −15 %) mostly at the same rivers in which altimetry data is available (e.g.Fig. 5c).But it also degrades results at some of the other river reaches.
Furthermore, the DA scheme can degrade results in some reaches where no data were assimilated.Such a problem is possibly caused by spurious correlations in the model covariance matrix from the EnKF due to a poor sampling from the ensemble.Aiming to avoid spurious correlations, methods such as covariance localization or local analysis (Sakov and Bertino, 2010) could be used to constrain the influence of observations based on distance criteria as already used in atmospheric or ocean applications.
However, in our view, a particular localization criteria should be developed for hydrological applications, since the correlation between the model states can be a function of an Euclidean distance in some cases (e.g.soil moisture) or in others, of a distance measured following the rivers' path (e.g.river discharge and water levels).
Results from this experiment demonstrate that the assimilation of radar altimetry data into large scale hydrologic models can improve simulations, mainly in terms of water levels but also in discharge.Even though ENVISAT data is provided at a 35day temporal resolution, its assimilation can improve model results at a daily basis as illustrated in Fig. 5b, c possibly due to the low temporal variability of Amazonian hydrographs and the fact that ENVISAT measurements are non-simultaneous.

Assimilation of discharge series based on satellite altimetry
In the last DA experiment (Exp.3), we evaluate the assimilation of discharge data derived from ENVISAT water level.Therefore, the data assimilated into the model has the same high spatial coverage and low temporal sampling as altimetry water levels have, but it also contains discharge information which is the most important hydrological variable of the model.In situ discharge data were also used for validation.
The DA scheme was able to assimilate altimetry-based discharges increasing the agreement between these observations and model results in most of the altimetry stations (Fig. 6a), as exemplified in Fig. 7a.In average, the rms error between altimetry-based discharges and model results decreased 23 % (Table 3), which represents a smaller improvement if compared to the assimilation of in situ discharge (∆rms = −49 % in Exp. 1a).Since observation errors are larger in the altimetry-based discharges, the DA scheme gives more weight to background model results and updated discharge values are not so close to measurements.The comparison of model results with in situ discharge data (Fig. 6b) shows that errors decrease mostly at gauging stations located at rivers where altimetry data were assimilated (e.g.Fig. 7b).At these sites, E NS changes from 0.76 to 0.80 and the mean ∆rms is −14 % (Table 3), which is comparable to the improvement obtained in the altimetry data assimilation (Exp.2, Table 2) over discharge results (∆rms = −15 %), but smaller than the enhancement of water level results (∆rms = −43 %).Moreover, similarly to Exp. 2, at gauging stations located outside the assimilation domain, the DA scheme has a minor effect on simulated discharge or degrades results in some cases.
Results from this section show the potential of assimilating discharge data derived from satellite altimetry into large scale hydrological models instead of in situ discharge or satellite water levels, even though such data has lower temporal resolution and accuracy if compared to data gathered at in situ gauging stations.

Prospects of streamflow forecasting
We now assess the potential of a large scale hydrologic-hydrodynamic model coupled with a DA scheme to provide streamflow forecasts in the Amazon basin.Since hydrological initial states governs discharge predictability at the large Amazonian rivers (discussed by Paiva et al., 2012b), we have chosen to generate forecasts starting with initial states gathered by the DA scheme and then using the ensemble streamflow prediction approach -ESP (Day, 1985), where the model is run forced by an ensemble of observed meteorological data from past years.Since this is a first attempt, we have chosen to evaluate forecasts using only the DA scheme configuration from Exp. 1b, where the model is updated using discharge data from 12 gauging stations located on the Amazon and its main tributaries.
We first evaluate hindcast forecasts at the Solim ões/Amazon mainstem, including upper Solim ões River at Tamishiyacu, Solim ões River at Manacapuru and Amazon River at Óbidos (Fig. 8).The model was able to forecast discharges with relatively high accuracy even for very large lead times (90 days).In all cases, forecasts are markedly better than simply using discharge climatology, as shown by positive values of SS cli skill score (Fig. 8).As expected, the agreement between model values and observations decreases as function of lead time, and, for example, SS cli decreases from 0.90 to 0.49 for forecasts 10 and 90 days ahead at Óbidos station.But it remains very high, showing that it would be possible to produce accurate forecasts at the Amazon main river for even larger lead times.Model performance also increases from the upper to the lower part of the Solim ões/Amazon River, and at the same time, the spread of the ESP 2901 ensemble at large lead times increases upstream and decreases downstream.Such behaviour is explained by the fact that in larger rivers, the hydrological predictability is much more influenced by the current volumes of water stored upstream than by future precipitation forcing, as discussed by Paiva et al. (2012b).Analysis from Fig. 8 also demonstrates that the model successfully predicted the severe 2005 drought at the Solim ões/Amazon main stem.At this year, discharges dropped ∼ 1 month earlier than normal (Fig. 8) and river levels fell to historically low levels causing navigation to be suspended (Marengo et al., 2008).Even so, the model was able to predict this low flows ∼ 90 days ahead.We now evaluate forecasts at gauging stations located all over the Amazon basin.
Forecasts for a smaller lead time of 5 days or even 15 days (Fig. 9a, b) were relatively accurate with positive SS cli values at several gauging stations located at both upstream and downstream rivers.However, the quality of the forecasts decreased as a function of lead time (30 and 90 days, Fig. 9c, d).It becomes very poor at smaller rivers and remains meaningful with positive SS cli values mainly at gauging stations with large draining areas.For 90 days lead time, SS cli index remains positive at almost all stations along the Solim ões/Amazon main stem and in some of the main tributaries.This behavior is also illustrated at Fig. 10.SS cli values are usually higher at gauging stations located in rivers draining large areas and decrease with lead time.For instance, if only stations gauging rivers with drainage area larger than 10 5 km 2 or 4 × 10 5 km 2 are considered, on average, forecasts remain better than climatology (SS cli > 0) up to ∼ 15 and ∼ 25 days lead time, respectively.On the other hand, if only the largest rivers are taken into account (> 10 6 km 2 ), SS cli values are high and always positive, which demonstrates the good performance of the model forecasts.SS cli is also high at stream gauges used for data assimilation where forecasts are usually better, SS cli is close to one for small lead times, as expected, and becomes negative after ∼ 55 days.

Summary and conclusions
We presented the development and evaluation of a data assimilation scheme for both gauged and satellite altimetry-based discharge and water levels into a large scale hydrologic-hydrodynamic model of the Amazon River basin using the Ensemble Kalman Filter -EnKF.We also evaluated hindcast forecasts based on this system using the ensemble streamflow prediction approach, where the model was forced by an ensemble of past precipitation forcing from TRMM mission.According to our results, the data assimilation scheme performed well in assimilating in situ and remotely sensed discharge and water levels into the large scale hydrologichydrodynamic model.The assimilation of in situ discharge showed that EnKF can improve discharge estimates at assimilation gauges but, differently from previous studies at smaller basins (e.g.Clark et al., 2008; and others summarised by Lee et al., 2012), also transfer information to ungauged rivers by improving results at validation sites, although with a smaller degree.The assimilation of discharge data at a reduced number of gauging stations located at larger rivers improves results mostly at the large reaches but it degrades results at some smaller basins.Also, the transformation of discharge measurements into the log space proved to be important to deal with very different discharge magnitudes arising from different spatial scales or from contrasting flood and recession flows.
The assimilation of satellite altimetry data improved model water levels, and also discharges, mostly at the same river reaches where altimetry stations are located.
Assimilating altimetry-based discharge also improved model estimates, although with minor degree if compared to the in situ discharge assimilation, probably due to the larger errors in remotely sensed observations.However, in both cases, even though radar altimetry data has low temporal resolution (35 days), its assimilation can improve model results at a daily basis, possibly due to its higher spatial resolution and the low temporal variability of Amazonian hydrographs.The sensitivity analysis of the parameter from the DA scheme highlighted the importance of the magnitude of precipitation errors and that of their spatial correlation, while temporal correlation showed to be dispensable.
The deterioration of model performance at some unmonitored reaches may be due to the large dimensionality of state space in distributed hydrological models compared to the available information.Consequently, data assimilation may update state variables at some lumped fashion such as the sub-basin scale, as explained by Lee et al. (2012).This problem can be also due to spurious correlations that can arise by numerical reasons and could be avoided by using proper spatial localization methods (e.g.Sakov and Bertino, 2010) developed for hydrological applications to constrain the influence of measurements.Additionally, the DA scheme could benefit from a better characterization of model errors, where not only precipitation but other sources of uncertainty, such as in model parameters and structure could be included, as suggested by Liu et al. (2012).
Although limitations still exist, results are encouraging.This kind of DA scheme could also be easily employed to other similar regional/global scale hydrological models (e.g.Yamazaki et al., 2011;Decharme et al., 2012;Alfieri et al., 2012).It has also the potential to improve by assimilating remotely sensed water levels gathered by other satellite missions as the existing ones, or the altimetry missions to be launched in the coming years by the European Spatial Agency ESA, namely the Sentinel-3 constellation and the forthcoming SWOT mission (Durand et al., 2010a).Moreover, the altimetry-based discharge assimilation can improve when better discharge estimates become available, such as the ones under development for the future SWOT mission (Durand et al., 2010b).Finally, the model was able to provide relatively accurate streamflow forecasts in the Amazon basin.For smaller lead times (∼ 5 to 15 days), forecasts agreed with observations in lots of gauging stations and for larger lead times (> 30 days) they remained meaningful mostly at larger rivers.Forecasts were usually better at stream gauges used for data assimilation, especially for smaller lead times.Along the Solim ões/Amazon

2883
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | solves the full 1-D Saint-Venant equations for the river network and flood inundation is simulated using a simple 2885 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper |

2887
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

2889
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | and monthly flood inundation extent from Papa et al. (2010b).Simulations agreed with observations, with relatively high Nash and Suttcliffe index (E NS ) values: E NS > 0.6 in ∼ 70 % of discharge gauges, E NS > 0.6 in ∼ 60 % of the water level stations derived from satellite altimetry, E NS = 0.71 for total flood extent and E NS = 0.93 for terrestrial water storage.2891 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | 2893 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | the climatological value of discharge on day t 2895 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | computed from observations.SS cli ranges from −∞ to 1 (optimum) and positive values show an improvement over a forecast based on climatology.

2897
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

2899
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Discussion
Paper | Discussion Paper | Discussion Paper | Discussion Paper |

2903
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Fig. 4 .
Fig. 4. Evaluation of ENVISAT radar altimetry data assimilation.Spatial distribution of change in root mean square error (∆rms) at (a) altimetry stations used for data assimilation and stream gauges with (b) stage and (c) discharge data used for verification.