Abstract

This paper describes a fully automated scheme that has provided calibrated 1–10-day ensemble river discharge forecasts and predictions of severe flooding of the Brahmaputra and Ganges Rivers as they flow into Bangladesh; it has been operational since 2003. The Bangladesh forecasting problem poses unique challenges because of the frequent life-threatening flooding of the country and because of the absence of upstream flow data from India means that the Ganges and Brahmaputra basins must be treated as if they are ungauged. The meteorological–hydrological forecast model is a hydrologic multimodel initialized by NASA and NOAA precipitation products, whose states and fluxes are forecasted forward using calibrated European Centre for Medium-Range Weather Forecasts ensemble prediction system products, and conditionally postprocessed to produce calibrated probabilistic forecasts of river discharge at the entrance points of the Ganges and Brahmaputra into Bangladesh. Forecasts with 1–10-day horizons are presented for the summers of 2003–07. Objective verification shows that the forecast system significantly outperforms both a climatological and persistence forecast at all lead times. All severe flooding events were operationally forecast with significant probability at the 10-day horizon, including the extensive flooding of the Brahmaputra in 2004 and 2007, with the latter providing advanced lead-time warnings for the evacuation of vulnerable residents.

1. Introduction

Accurate and timely forecasts of river flow have the potential of providing critical information for water resource management, agriculture practice optimization, and disaster mitigation. Nowhere is the need for such forecasts more urgent than in the Bangladesh delta that lies at the confluence of two of the largest river systems in the world: the Brahmaputra and the Ganges (Fig. 1). Combined, these rivers can possess discharges of well over 100 000 m3 s−1. Once these flood-stage waters reach Bangladesh—a small country the size of England (140 000 km2) but with a population of 155 million—they spread out onto a flat and wide floodplain occupying roughly 80% of the country, making Bangladesh part of the world’s largest river delta. The combined catchment area of the Brahmaputra–Ganges system (approximately 1.5 × 106 km2) ranks tenth in the world, with only the Amazon and the Congo surpassing the combined discharge of the two rivers.

Fig. 1.

(a) The Brahmaputra and Ganges River catchments with the blue and red lines denoting the major channels of each river, respectively. Bangladesh (yellow) sits at the confluence of the two rivers. (b) Detail of Bangladesh showing the paths of the two rivers. Major river gauging stations are located at Bahadurabad on the Brahmaputra and Hardinge Bridge on the Ganges. Daily streamflow measurements have been made at these two stations at least since the early 1950s. They also represent the only in situ data available to CFAB from which to develop extended predictions. Also shown are the five unions (green: Rajpur, etc.) selected for dissemination of the CFAB forecasts in 2006, which received advanced flood danger warnings during the two severe 2007 Brahmaputra flooding events.

Fig. 1.

(a) The Brahmaputra and Ganges River catchments with the blue and red lines denoting the major channels of each river, respectively. Bangladesh (yellow) sits at the confluence of the two rivers. (b) Detail of Bangladesh showing the paths of the two rivers. Major river gauging stations are located at Bahadurabad on the Brahmaputra and Hardinge Bridge on the Ganges. Daily streamflow measurements have been made at these two stations at least since the early 1950s. They also represent the only in situ data available to CFAB from which to develop extended predictions. Also shown are the five unions (green: Rajpur, etc.) selected for dissemination of the CFAB forecasts in 2006, which received advanced flood danger warnings during the two severe 2007 Brahmaputra flooding events.

Each year, short-lived flooding occurs throughout the summer and early autumn with sufficient irregularity to produce adverse social and agricultural effects. In these normal years, about one-fifth of the country can be inundated for days. The timing of the flooding is critical for the harvesting of the winter rice crop, the planting of summer rice, and the viability of other crops planted during the monsoons. Major flooding events occur in Bangladesh with a return period of 4–5 yr. For example, during the summer of 1998, more than 60% of Bangladesh was inundated for nearly three months (Ahmad and Mirza 2000). Again in July 2004, widespread flooding occurred as the Brahmaputra exceeded dangerous flood levels twice over a two-week period. In July and September 2007, the Brahmaputra severely flooded twice with inundation lasting weeks.

There is considerable utility in producing timely and accurate river flood forecasts for Bangladesh (RTI and EGIS 2000). Such forecasts can impact short-term agricultural and societal decision making, specifically when to evacuate flood-prone regions and ensure protection of food and water supplies and livestock safety. Preemptive actions taken can significantly lower relocation and reconstruction costs following a period of flooding. RTI and EGIS estimates that a 7-day forecast has the potential of reducing postflood costs by as much as 20% in household and agricultural costs compared to 3% with a 2-day forecast.

Predicting river discharge in Bangladesh is little different from that existing in many developing nations but with one very particular problem: the Brahmaputra and Ganges Rivers exit India to flow into Bangladesh, and there has been little consistent river discharge data sharing between these two countries. The only reliable river streamflow data come from what Bangladesh measures once the rivers cross the India–Bangladesh border (see Fig. 1). As a result, Bangladesh forecast lead times have been limited to 2–3 days for the interior of the country and have essentially no lead time in areas close to the border.

In the absence of upstream river discharge data from India, a different modeling approach would need to be taken to extend warning times of severe flooding. The Ganges and Brahmaputra catchments upstream of the India–Bangladesh border would have to be treated as if they were ungauged and effectively naturalized rivers. Can the latter be an accurate approximation for these highly utilized river basins?

The Brahmaputra has yet to have a major hydraulic structure built along its course (Singh et al. 2004) and as such, it can be modeled as a naturalized river. Although there are several irrigation projects along its course, their design discharges are typically small (well below 100 m3 s−1). However, the management of the Ganges River presents a further complication, primarily as a result of at least 34 barrages and other structures functional in India and Nepal (Mirza 2004). The actual amount of the water diverted by these structures is not known. However, these structures were designed primarily for use in the dry season and not during the monsoonal flood season, which is the period of interest of this study. For example, the Farakka Barrage, a major diversion located 10 km from the India–Bangladesh border, is diverted primarily to reduce silt in the Calcutta harbor during the dry season. During the monsoon season (June–October), there appears to be little difference between discharge–precipitation relationships before or after the 1974/75 construction of the barrage (Jian et al. 2009). On the basis of the Jian et al. (2009) study and the Flood Forecasting and Warning Centre (FFWC; 2000, personal communication), it is assumed Farakka diversions do not affect discharge into Bangladesh beyond 15 June.

The catastrophic Bangladesh flooding in 1998 encouraged the U.S. Agency for International Development (USAID) to fund an exploratory project called the Climate Forecast Applications in Bangladesh (CFAB) in 2000. The primary goal of CFAB was to provide advanced warning of Brahmaputra and Ganges flooding and drought in Bangladesh on time scales ranging from daily to seasonal, in partnership with the FFWC and the European Centre for Medium-Range Weather Forecasts (ECMWF). By utilizing ECMWF ensemble weather forecasts and conditional estimates of hydrologic modeling error, the CFAB discharge forecasts were designed to also provide forecasts of their own uncertainty. Operational (probabilistic) forecast dissemination to Bangladesh began during the monsoon season of 2003, based on ECMWF’s medium- and seasonal-range ensemble forecasts. Table 1 describes the CFAB three-tier system and illustrates the concept of overlapping forecasts to allow a nesting of strategic decisions and subsequent tactical adjustments. In this paper, we focus on the medium-term forecasting tier. (The project is described at http://cfab.eas.gatech.edu, and the three-tier forecasts generated as per Table 1 presented at http://cfab2.eas.gatech.edu.)

Table 1.

The three-tier overlapping CFAB forecasting scheme developed for Bangladesh. The table lists the forecast period of each tier, when the forecast is updated, and the type of economical/agricultural decision that can be made at each tier.

The three-tier overlapping CFAB forecasting scheme developed for Bangladesh. The table lists the forecast period of each tier, when the forecast is updated, and the type of economical/agricultural decision that can be made at each tier.
The three-tier overlapping CFAB forecasting scheme developed for Bangladesh. The table lists the forecast period of each tier, when the forecast is updated, and the type of economical/agricultural decision that can be made at each tier.

Section 2 provides a description of the medium-range forecast scheme. Five components are introduced and described in detail. Examples of the forecasts for 2003, 2004, and 20071 for the Brahmaputra and the Ganges are presented in section 3 along with forecast verifications. Section 4 summarizes the utility of the model and why it is necessary to include meteorological forecasts in a hydrological scheme. Possible extensions of the scheme are also described in section 4.

2. The hydrological prediction scheme

Figure 2 provides a flowchart of the Bangladesh multimodel discharge forecasting procedure. Step I describes the initial data inputs from ECMWF, National Aeronautics and Space Administration (NASA) and National Oceanic and Atmospheric Administration (NOAA) satellite and rain gauge products, river discharge measured at the India–Bangladesh border by the FFWC, and updated parameters for the hydrological models. In step II, the ECMWF ensemble predictions of precipitation are statistically modified using the satellite and rain gauge data to remove systematic errors. Together with the discharge data, these data are ingested into a set of hydrological models in step III. Step IV takes into account forecast uncertainties while preserving the useful information of the ensemble weather forecast dispersion, and it makes a final error correction to produce the final probabilistic forecast of river discharge. Step V is the dissemination of the forecasts, which is not discussed in detail in this paper. We note that the complete forecasting procedure shown in Fig. 2 is fully automated.

Fig. 2.

Flowchart of the CFAB short-term forecasting scheme. Five steps (I–V) are outlined schematically and within each is a brief description of the data used, tasks performed, and so on. The arrows show the path of the data through the system. Step I, clear box: these components represent daily inputs into the scheme. Step II, dotted box: the ECMWF ensemble forecasts are corrected statistically to reduce systematic error. Step III, hatched boxes: these represent the hydrologic modeling process, where all inputs are integrated to produce Brahmaputra and Ganges river flow. Step IV, checkered box: this represents the process by which all hydrologic uncertainties are accounted for in the probabilistic forecasting process. Step V: the ensemble forecasts are tailored to produce probabilistic information necessary for user needs. Note that all processes below the horizontal dashed line are done independently for each day and each 1–10-day forecast lead time.

Fig. 2.

Flowchart of the CFAB short-term forecasting scheme. Five steps (I–V) are outlined schematically and within each is a brief description of the data used, tasks performed, and so on. The arrows show the path of the data through the system. Step I, clear box: these components represent daily inputs into the scheme. Step II, dotted box: the ECMWF ensemble forecasts are corrected statistically to reduce systematic error. Step III, hatched boxes: these represent the hydrologic modeling process, where all inputs are integrated to produce Brahmaputra and Ganges river flow. Step IV, checkered box: this represents the process by which all hydrologic uncertainties are accounted for in the probabilistic forecasting process. Step V: the ensemble forecasts are tailored to produce probabilistic information necessary for user needs. Note that all processes below the horizontal dashed line are done independently for each day and each 1–10-day forecast lead time.

a. Step I: Initial data

There are four major sets of initial data in Fig. 2: the ensemble forecasts from ECMWF (box E), satellite and rain gauge data (box S), river discharge data (box Q), and hydrological model parameters (box MP). The arrow at the base of the boxes shows where the product or derivative of the product is used.

1) ECMWF (EPS) forecasts (box E)

A catchment integrates spatially and temporally all weather factors (precipitation, evaporation, etc.), together with hydrological effects (such as aquifer flow and soil moisture changes), into river flow and water storage and eventually river discharge out of the basin. Because of the spatial scale of the Ganges and Brahmaputra basins, the relatively large-grid global circulation weather model output from the ECMWF ensemble prediction scheme (EPS) can be used to generate useful discharge forecasts for these essentially ungauged basins. The ECMWF weather variables used in this study are surface fields of wind, humidity, and precipitation. The weather forecasts are provided as 51 ensemble members for each variable and for each forecast lead time. The model resolution was initially at approximately 80 km × 80 km grid (TL255L62) going out to 240 hours at 6-h increments; however, the current operational EPS is at approximately 50 km × 50 km grid (TL399L62 from 0 to 10 days), with the forecasts horizon also extending to 15 days (days 9–15 at TL255L62; Buizza et al. 2007). All forecast fields are interpolated to the same ½° × ½° grid used by the semidistributed hydrological model (SDM).

The joint CFAB–ECMWF effort was one of the first efforts to utilize ensemble medium-range numerical weather forecast input in the generation of operational short-to-medium-range ensemble discharge forecasts. The ECMWF system provides not only weather variable forecast information but also estimates of the likely forecast error through both 51-member medium-range ensembles (Molteni et al. 1996). The ECMWF ensemble predictions are passed in step II for statistical evaluation and correction.

2) Satellite and rain gauge precipitation estimates (box S)

For forecasting discharge of large spatial-scale catchments, and therefore long hydrological temporal scales, it is essential that the initial state of the hydrologic system (soil moisture content, groundwater storage, in-stream flows, etc.) be assessed as accurately as possible. The importance of setting these initial conditions correctly can be understood through the following simplified argument: i-day lead-time discharge forecasts require i-day weather forecasts and (ni) days of weather observations, for watersheds of time scale n, where n can be appreciable for large river basins.

Three precipitation observations were integrated in the operational forecast system in 2004. One source is ground-based ½° × ½° gridded daily averaged rainfall estimates from the NOAA Climate Prediction Center (CPC) utilizing World Meteorological Organization (WMO) global telecommunication system (GTS) reporting rain gauges and using Shepard’s inverse-distance algorithm for interpolation (Xie et al. 1996). Note, however, that the reporting gauge density is sparse over the Brahmaputra and Ganges catchments, with roughly 40 gauges reporting each day. The other two estimates are near-real-time ¼° × ¼° 3-hourly satellite-derived precipitation accumulations from the NASA Tropical Rainfall Measuring Mission (TRMM; Huffman et al. 2005, 2007) and the NOAA CPC morphing technique (CMORPH; Joyce et al. 2004) programs, both of which utilize infrared- and microwave-frequency satellite data but distinctly different algorithms.

These datasets are continuously updated as they become available and are interpolated to the same hydrologic computational temporal and spatial grid, with a simple grid-by-grid arithmetic average of the available data taken as the best estimate of the observed precipitation field at the time of daily forecast initialization. This observational estimate is used in both the catchment lumped model (CLM) and in the initialization (soil moisture states and in-stream flows) of the SDM. All other observed weather variables used in SDM evapotranspiration calculations are derived from the ECMWF analysis fields.

3) River discharge data (box Q)

An essential data source for skillful forecast generation are discharge and stage measurements at the near-border gauging stations of Bahadurabad for the Brahmaputra and Hardinge Bridge for the Ganges (Fig. 1b) provided by the FFWC. These stations have observed discharge records derived from standard river cross-section survey techniques completed approximately once per week for more than 13 years and daily river stage records of more than 52 years. These historical discharge records are essential for model calibration, and the near-real-time observations are utilized for hydrologic model data assimilation and error corrections considered later in step III of Fig. 2 and for the final discharge forecast probability distribution function (PDF) generation considered in step IV of Fig. 2.

Daily discharge estimates were derived from fitting rating curves to the FFWC observational data using rating curve functions of form

 
formula

where Qt is the estimated discharge for day t, Ht is the observed river stage (relative to a fixed datum), and i represents the stage interval for which the equation applies. The Ganges used a two-part rating curve i ∈ [1, 2], and the Brahmaputra, a three-part form. Because sediment transport and bed mobilities of both the Brahmaputra and Ganges are substantial, the parameter sets of these equations were fitted independently for each year over the historical record.2 This was done by weighting the residual error inversely proportional to the year separation and minimizing it using a nested version of the downhill simplex method in multidimensions [Nelder and Mead 1965; A Multidirectional Optimum Ecotope-Based Algorithm (AMOEBA) procedure; Press et al. 1992].

Absolute rating curve error ɛ is estimated to be |ɛ| = 0.10Q for Bahadurabad and |ɛ| = 0.07Q for Hardinge Bridge. By accounting for additional systematic measurement errors in the FFWC method to measure the observed river discharges of order |ɛ| − 0.1Q (rough) and by assuming independence of these sources of error to combine them in quadrature (Taylor 1982), we estimate total error bounds on the observed discharges (shown in Figs. 6, 8, 9, and 10) of |ɛ| ≈ ±0.15Q.

Fig. 6.

Multimodel forecasts of the (a) Brahmaputra and (b) Ganges Rivers. To obtain these forecasts, each of the 51-member ECMWF ensemble precipitation forecasts are run through the hydrological model to produce the ensemble (colored traces) with 1-day precipitation forecasts producing 1-day discharge forecasts, 2-day precipitation forecasts for the 2-day discharge forecasts and so on. The horizontal dashed line is the “critical discharge” flood level for the Brahmaputra at Bahadurabad and for the Ganges at Hardinge Bridge. Forecasts for 1, 4, 7, and 10 days are shown. The abscissa shows the dates for which the forecast is made. These results have not been “error corrected.” Notice that it is only at about the day-4 forecasts that the ensemble spread of the precipitation forecasts becomes noticeable for the Brahmaputra. For the Ganges, the spread does not become noticeable until after day 7.

Fig. 6.

Multimodel forecasts of the (a) Brahmaputra and (b) Ganges Rivers. To obtain these forecasts, each of the 51-member ECMWF ensemble precipitation forecasts are run through the hydrological model to produce the ensemble (colored traces) with 1-day precipitation forecasts producing 1-day discharge forecasts, 2-day precipitation forecasts for the 2-day discharge forecasts and so on. The horizontal dashed line is the “critical discharge” flood level for the Brahmaputra at Bahadurabad and for the Ganges at Hardinge Bridge. Forecasts for 1, 4, 7, and 10 days are shown. The abscissa shows the dates for which the forecast is made. These results have not been “error corrected.” Notice that it is only at about the day-4 forecasts that the ensemble spread of the precipitation forecasts becomes noticeable for the Brahmaputra. For the Ganges, the spread does not become noticeable until after day 7.

Fig. 8.

Brahmaputra discharges and forecasts for the most recent year with critical stage flooding, 2007. Shown are forecast horizons of (column I) 5 and (column II) 10 days. (row i) The complete ensemble distributions. Abscissa dates denote the forecast target dates. For example, the first date (20 Jun) in column I refers to a forecast made five days earlier ; in column II, the first date refers to a 10-day forecast made on . The bold black line denotes the observed flow matching the forecast dates, the horizontal dashed line denotes the river’s critical discharge level, and the distance between the plots’ right edge and the vertical line designates the forecasting interval. (row ii) From the ensemble set shown in row (i) 95% (red) and 50% (gold) confidence bounds are calculated and plotted. (row iii) The cumulative probability (red line) the river will exceed its danger level (dashed line), with the black line again denoting observed discharge.

Fig. 8.

Brahmaputra discharges and forecasts for the most recent year with critical stage flooding, 2007. Shown are forecast horizons of (column I) 5 and (column II) 10 days. (row i) The complete ensemble distributions. Abscissa dates denote the forecast target dates. For example, the first date (20 Jun) in column I refers to a forecast made five days earlier ; in column II, the first date refers to a 10-day forecast made on . The bold black line denotes the observed flow matching the forecast dates, the horizontal dashed line denotes the river’s critical discharge level, and the distance between the plots’ right edge and the vertical line designates the forecasting interval. (row ii) From the ensemble set shown in row (i) 95% (red) and 50% (gold) confidence bounds are calculated and plotted. (row iii) The cumulative probability (red line) the river will exceed its danger level (dashed line), with the black line again denoting observed discharge.

Fig. 9.

Same as Fig. 8 but for the Ganges for its most recent critical flood stage year, 2003.

Fig. 9.

Same as Fig. 8 but for the Ganges for its most recent critical flood stage year, 2003.

Fig. 10.

Brahmaputra discharges, forecasts, and above-danger-level probabilities for flood year 2004, presented as in Fig. 8. (i),(iii) The season’s 10-day lead forecasts and severe flooding probabilities, respectively, with (ii) and (iv) zooming in on the July severe flooding. (i)–(iv) The distance between the plots’ right edge and the vertical black line designates the forecasting interval, showing severe floods were forecast as the hydrograph was descending . (v),(vi) The 1–10-day ensemble plume and associated above-danger-level probabilities, highlighting both the skill and relative increase in uncertainty of the forecasts with lead time; forecasts issued 5 and 16 Jul, respectively, and both issued while the observed hydrograph was decreasing (left axis). The discharge first exceeded critical level on 11 Jul, and the two peak discharges occurred on 13–14 Jul (estimated 86 600 m3 s−1) and 23 Jul (81 600 m3 s−1). Black vertical lines bound the flooding period. Solid red lines give the probability forecasts (derived from the ensembles) the discharge will exceed its danger level (right axis).

Fig. 10.

Brahmaputra discharges, forecasts, and above-danger-level probabilities for flood year 2004, presented as in Fig. 8. (i),(iii) The season’s 10-day lead forecasts and severe flooding probabilities, respectively, with (ii) and (iv) zooming in on the July severe flooding. (i)–(iv) The distance between the plots’ right edge and the vertical black line designates the forecasting interval, showing severe floods were forecast as the hydrograph was descending . (v),(vi) The 1–10-day ensemble plume and associated above-danger-level probabilities, highlighting both the skill and relative increase in uncertainty of the forecasts with lead time; forecasts issued 5 and 16 Jul, respectively, and both issued while the observed hydrograph was decreasing (left axis). The discharge first exceeded critical level on 11 Jul, and the two peak discharges occurred on 13–14 Jul (estimated 86 600 m3 s−1) and 23 Jul (81 600 m3 s−1). Black vertical lines bound the flooding period. Solid red lines give the probability forecasts (derived from the ensembles) the discharge will exceed its danger level (right axis).

We note that both of these rivers have low gradients around the forecasting locations, so it was anticipated that pressure gradient effects due to variable surface water slopes, along with effects due to the filling and draining of overbank storage, may limit the utility of fixed stage–discharge relationships via this potential hysteresis. However, an attempt to minimize the 10% rating curve error just discussed for the Brahmaputra by including an additional hysteresis correction term in Eq. (1) never affected the estimated discharge more than 0.2%.

4) Updated SDM parameters (box MP)

The SDM model parameters are continuously updated through an automated background calibration process that uses a modified version of the root-finding AMOEBA procedure (Nelder and Mead 1965; Press et al. 1992) to minimize the square error between observed and modeled discharge. The search is performed over a reduced nine-parameter set to maintain a small ratio of adjustable parameters to daily discharge calibration points, and thus reduce the risk of overcalibration. These nine parameters are discussed in the SDM model, section 2c. If, at the time of daily discharge forecast initialization, a new optimal root has been found, soil moisture states and in-stream flows are automatically updated before the forecast proceeds.

b. Step II: Statistical rendering and downscaling of forecasts (dotted box, Fig. 2)

Generally, all practical hydrologic models require some form of calibration arising from reasons such as incomplete knowledge of watershed properties and the necessary parameterization of transport at unresolved scales. Through the hydrologic model parameter calibration process, some of the additional errors in the watershed inputs (rainfall, in particular) can also be implicitly reduced (e.g., runoff errors arising from an overbias in rainfall forcing can be implicitly minimized by increasing the parameterization of evapotranspiration). However, such bias reduction through calibration does not occur if there are relative errors between more than one input of the same physical quantity. In particular, using both weather forecasts and observationally based estimates of precipitation concurrently to generate unbiased discharge forecasts requires that these two data sources maintain their statistical similarity, such that there are no relative biases between their statistical moments (mean, variance, and skewness, in particular).

To minimize systematic differences in hydrologic model inputs leading to systematic hydrologic forecasting errors, a quantile-to-quantile (q-to-q) mapping technique was implemented. A similar approach has been applied to adjust nonsimultaneous radar data to rain gauge observations (Calheiros and Zawadzki 1987). As applied here, the q-to-q technique forces each ECMWF ensemble precipitation forecast at each ½° × ½° grid point (input E, Fig. 2) to be statistically sampled from the cumulative distribution function (CDF) of the associated observational historical record for the same grid point (input S, Fig. 2). In essence, this converts the numerical weather prediction product from a physical variable to a probability forecast by performing a bias correction at each quantile of the forecast CDF. This bias correction approach was implemented for the 2004 monsoon season for catchment lumped values of the ECMWF 40-member seasonal forecasts (Webster et al. 2006) and for each forecast grid of the ECMWF medium-range forecasts (Hopson 2005) the same year. The manner in which the technique is applied to the medium-range forecasts is described below.

Calculate the monsoon season “climatological” empirical CDFs for both the observations and the forecast model at each ½° × ½° grid location i. The observational (subscript o) climatological CDF for i and precipitation amount y, Coi(y), is derived from the previously discussed rain gauge (available from 1978 onward) and merged gauge–satellite data (available from 2002 onward). The forecast (subscript M) CDF’s at f and i, CMi,f(y) are derived from stored forecasts (available from 2003 onward), which are first downscaled through simple two-dimensional linear interpolation to the ½° × ½° computational grid of the rainfall observations and hydrologic models, and the CDF computation is performed independently for each forecast f (24, 48 h, etc.). Formally, these CDF’s are given by

 
formula
 
formula

where P{α} denotes the probability of the event α and x (y) is a specific value of the generic random variable X (Y) representing the observed (forecast) weather variable of interest. Using these CDFs, a direct quantile-by-quantile bias correction is made to each grid point i, f, and precipitation ensemble member fcst of the forecast, , of the forecast by simply setting (1) equal to (2) and determining the value of x = padj that satisfies the equation. In this manner, the “observation space” quantile is matched with the “forecast model space” quantile as shown in Fig. 3, and padj is used in replace of in the hydrometeorological application.

Fig. 3.

The q-to-q correction system. Both modeled and observed precipitation are binned into quantiles. The model precipitation is mapped onto the observed precipitation fields by setting respective modeled quantiles to observed quantiles. In the figure, pfcst of the 80th quantile is set to padj 80th quantile. The method ensures that the forecasts produce the same number of no-rain events as the observations.

Fig. 3.

The q-to-q correction system. Both modeled and observed precipitation are binned into quantiles. The model precipitation is mapped onto the observed precipitation fields by setting respective modeled quantiles to observed quantiles. In the figure, pfcst of the 80th quantile is set to padj 80th quantile. The method ensures that the forecasts produce the same number of no-rain events as the observations.

Note that this technique ensures that the forecasts produce the same climatological number of “no rain” events as observed. It also preserves the spatial and temporal (ranked) covariances of the weather variable forecast fields (as generated by the numerical weather prediction model). Figure 4 shows an example of both the spatial- and temporal-preserving ability of this approach. The top panel is a 144-h accumulated field of the original ECMWF 24-h forecasts from 20 to 25 July 2007 preceding the severe Brahmaputra flooding event from 26 July to 6 August. The lower panel shows the same accumulated field after being “q-to-q mapped” to the observational PDF. Note that structures in the original forecast precipitation field have been preserved.

Fig. 4.

Accumulated precipitation for 0–144-h 1200–1200 UTC 20–26 Jul 2007 preceding the 26 Jul–6 Aug severe Brahmaputra flooding event (see Fig. 8). The Brahmaputra catchment is shown in outline directly above the Bay of Bengal. Note that the spatial scale of the precipitation event is approximately the same as that of the Brahmaputra catchment itself, arguing why the CLM shows good forecasting skill.

Fig. 4.

Accumulated precipitation for 0–144-h 1200–1200 UTC 20–26 Jul 2007 preceding the 26 Jul–6 Aug severe Brahmaputra flooding event (see Fig. 8). The Brahmaputra catchment is shown in outline directly above the Bay of Bengal. Note that the spatial scale of the precipitation event is approximately the same as that of the Brahmaputra catchment itself, arguing why the CLM shows good forecasting skill.

If there is no contemporaneous forecast and observational data, we suggest that the q-to-q technique may be the best one for bringing forecasts into statistical similarity with the observations. However, although the mapping improves the forecast bias at each quantile by steering the forecasts into having the same frequency of events as the observation climatological record, this process still does not necessarily ensure that the forecasts’ skill will be increased (which we discuss elsewhere in detail). Such skill degradation can (arbitrarily) arise because the technique does not directly improve the forecast model’s reliability by neglecting the conditional relationship between observations and forecasts. Although not discussed here, to improve the forecast’s conditional bias and ensemble dispersion, following the gridded q-to-q bias correction, we employ quantile regression to the catchment-averaged values of the precipitation forecasts.

c. Step III: Hydrological modeling (hashed boxes, Fig. 2)

Next we discuss the two hydrological model approaches we employ along with their relative attributes and drawbacks: the CLM and the SDM. These are combined into a multimodel following the philosophy of Krishnamurti et al. (1999), Krishnamurti et al. (2000), and Georgakakos et al. (2004), who noted collectively that some models are better equipped in certain forecast circumstances than others and that the combined performance of the multimodel generally outperforms any one model. The aim of the multimodel approach applied here is to minimize hydrological modeling error, with total hydrologic forecast uncertainty dealt with separately in step IV. The availability of near-real-time discharge observations allows for autoregressive (AR) error corrections to be applied to the SDM forecasts before the two hydrological models are combined into the multimodel prediction, which is also discussed. A flowchart of each of the hydrological models and the manner in which they are combined is shown in Fig. 2.

1) CLM

The lumped modeling technique employed here is based on the work of Young and Beven (1994) and Beven (2000), whose methodology has been applied successfully in a flood forecasting scheme for the CI6 catchment at Llyn Briane in Wales and for the catchment of the River Nith upstream of the town of Dumfries in Scotland (Lees et al. 1994). This approach formed the basis of the 2003 CFAB operational ensemble forecasts for the Ganges and Brahmaputra.

The CLM employs a multiple linear regression between observed discharge and observed and forecasted precipitation, together with a nonlinear “effective rainfall” filter. The model is based on the idea of linear storage reservoirs, and the model structure is similar to the unit hydrograph approach (Sherman 1932), the transfer function plus noise (TFN) models of Box and Jenkins (1970), and to an autoregressive moving average model with exogenous variables (ARMAX) time series forecasting model (Haltiner and Salas 1988; Marco et al. 1993), where here the “exogenous variable” is a filtered (simultaneous and forecasted) catchment-averaged rainfall time series but the “moving average” process is not being explicitly modeled. The benefits of the CLM are that it has a flexible structure, it can be easily recalibrated daily for different forecast time horizons and for different in situ conditions, and it easily assimilates near-real-time discharge measurements. Drawbacks of the CLM include limited long-time-scale base-flow modeling and heterogeneities in soil moisture and the distribution of rainfall–runoff across the watershed for which they are not accounted. However, the CLM approach can still be quite effective under certain settings, such as when precipitation events are of a similar spatial scale to the catchment itself, or less strictly, when there is a sufficiently large number of precipitation events compared to the catchment’s inherent time scale, and their spatial distribution tends to maintain statistical similarity over time.

The CLM approach, employed operationally, uses the time series of both observed satellite–rain gauge combined product, as discussed earlier, and the daily ECMWF ensemble forecasts of catchment-averaged rainfall. These are adjusted to an estimated fraction of effective rainfall contributing to storm runoff ui using the nonlinear filter

 
formula

where Qi is the observed discharge at time i (or taken as the most recently observed discharge in the case of a precipitation forecast), Qmax is an arbitrary normalization constant, c is a power to be calibrated, and Ri is the catchment-averaged rainfall (observed and forecast). This filter provides a simple way to estimate effective rainfall dependent on watershed discharge used as a surrogate for antecedent soil moisture conditions.

Using the derived effective rainfall time series of (4) and the time series of daily observed catchment discharge, the linear forecast equation for the future discharge Qt+f is given by

 
formula

where t is the current day, f is the forecast lead time in days, the sets {ai} and {bj} are coefficients, and ɛt+f is the model error. For model calibration, all ui are derived from observed precipitation. However, in forecasting, the set {ut+f , ut+f−1, … , ut+1} is replaced by individual members of the 51-member ECMWF ensemble precipitation forecasts, which preserves the temporal covariability of the forecasted precipitation. This then produces a 51-member ensemble of discharge forecasts whose dispersion represents the effects on discharge uncertainty of input forecast precipitation uncertainty.

The number and value of coefficients {ai} and {bj}, the power c in (4), along with the length of a boxcar smoothing window performed on the time series of ui, are recalibrated before each daily forecast initialization using updated observed discharge data. This calibration is done independently for each forecast lead time, with the best model determined by having the minimum Bayesian information criterion (BIC; Schwarz 1978; Katz 1981). The daily calibrated models are then used to generate forecasts, as discussed earlier, and also to calculate hindcasts (box CLM, Fig. 2), which are used in the multimodel selection process (box M, Fig. 2). Note that observed precipitation is used in generating the hindcasts to isolate model uncertainty from weather variable forecast uncertainty in the multimodel selection, discussed further in step IV. Both the forecasts and hindcasts are then passed to the multimodel step.

2) SDM

In 2004, an SDM approach was included into the operational system (box SDM, Fig. 2) under the multimodel framework with the CLM. In contrast to the CLM, the strengths of the SDM approach are that water-balance physics (evapotranspiration, soil storage, etc.) are more explicitly modeled locally and the time delays of spatially distributed precipitation events are taken into account. The drawbacks are that model recalibration is more costly and therefore cannot be easily tuned to the current conditions of the watershed. Similarly, data assimilation of observed discharge into the forecast structure is not as directly applied as with the CLM.

The SDM employed here is a simple nine-parameter two-layer representative area soil model parameterizing most of the physical processes of water storage and transport with a slow and a fast time-scale response. Water balance calculations are carried out on a 0.5° × 0.5° computational grid, with each grid treated as a subcatchment. At this scale, heterogeneous water storage and transport processes are homogenized, with the model providing only conceptual representations of these processes. This model was designed to be of intermediate complexity but with a limited number of parameters, because both the large spatial scale and the lack of available data in the Ganges and Brahmaputra catchments for both discharge and soil and vegetation data limits the benefits of employing more complex models. Furthermore, the use of more complicated models increases the inherent risk of model overcalibration, as was discussed previously. The coupled nonlinear differential equations for the two-layer model are

 
formula
 
formula

where S1 and S2 are the amounts (length) of water in the upper and lower soil layers, respectively; rs1 and rs2 are the sizes (length) of the upper and lower soil layer reservoirs, respectively; ts1 and ts2 are the time constants for lateral outflow from the upper and lower soil layers, respectively; tp is the time constant for percolation from the upper to the lower soil layer; P is the precipitation rate (length/time); and Ep is the potential evaporation rate (length/time), which is described below. However, if S2 > r2, then the two-layer model simplifies to

 
formula
 
formula

The outflow of water (length/time) is given by

 
formula

The discharge calculations for each grid are calculated as a daily rate using Eq. (10) along with integrations of the two rate equations, Eqs. (6) and (7). These integrations are performed on 6-h intervals using a fourth-order Runge–Kutta method to generate the final daily discharge value.

The potential evaporation Ep is calculated using the aerodynamic method (Thornthwaite and Holzman 1939) using the ECMWF analysis and control (deterministic member) forecast fields and with an additional calibration coefficient. Although the aerodynamic method is derived for conditions where moisture supply is not limited, the evaporation losses here are made to depend linearly on the moisture availability of the upper soil surface through the factor (S1/rs1) in (6).

To route the discharge out of each subcatchment and down to the forecasting location and to simulate a typical watershed discharge response, the subcatchment is treated as a series of identical linear reservoirs (e.g., Nash 1958; Chow et al. 1988). For our application, this is done using an n-reservoir instantaneous unit hydrograph with its own decay constant tc (both n and tc are calibrated). For the jth subcatchment, this is given by

 
formula

where Γ(n) is the gamma function, m is the number of days prior to t presumed to generate outflow (operationally taken as 40), and Tlagdj represents the subcatchment travel-time delay from the subcatchment outlet to our forecasting location, where dj is the (normalized) distance from the subcatchment outlet to our forecasting location, and Tlag is a calibration coefficient (prescribed over the whole watershed). Summing over all subcatchments, j then gives the total discharge at the forecasting location.

Taken all together, there are five parameters of the two-layer model (rs1, rs2, ts1, ts2, tp), one evapotranspiration parameter, and the three routing parameters (n, tc, and Tlag), comprising nine calibration parameters for the whole watershed. The calibration of these parameters is discussed in section 2a(4). For each subcatchment, there are the two state variables S1 and S2 and four forcing variables (precipitation, wind, humidity, and temperature).

Note that snowmelt is not explicitly modeled in the SDM as presented here. Although this is a model deficiency, monsoon rainfall is by far the dominant source of river discharge during the summer months. For the Ganges, snow meltwater provides much less than 1% of the annual river flow volume (Central Water Commission 1988). However, for the Brahmaputra, roughly 10% of the monsoon runoff is derived from snowmelt. If one estimates natural variability in snowpack on the order of 10%, then only a few percent of the seasonal discharge variability is neglected without including snow runoff processes.

For each forecast initialization, discharge hindcasts are generated using the updated model parameters of box MP in Fig. 2 and using the ECMWF weather variable analysis fields. Precipitation, however, comes from the hourly updated satellite and rain gauge observations discussed earlier (box S, Fig. 2). Forecasts are then calculated forward from the updated current hydrologic state using the adjusted ECMWF EPS 51-member precipitation ensemble (output from step II, Fig. 2) and the deterministic forecast for all other weather variables.

As with the CLM, ensemble precipitation forecasts are transformed into ensemble discharge forecasts via the hydrologic modeling process, and the effect of the temporal covariances of the precipitation fields is preserved in the resultant ensemble discharge forecasts by tracking single-ensemble-member precipitation forecast fields through the hydrologic model from 1 to 10 days. Both the discharge forecasts and hindcasts are then passed to the error correction step discussed in the next section, which assimilates the FFWC discharge observations at Bahadurabad and Hardinge Bridge (box Q, Fig. 2), available at the time of forecast initialization.

3) Autoregressive error corrections

Young (2002) noted that the accuracy of operational hydrological forecasts could be improved if real-time corrections are applied to improve model performance utilizing recent discharge observations and hydrologic model errors. For our SDM discharge forecasts, we employ an AR data assimilation technique to forecast combined hydrologic model and weather observation error. This model is refit daily for each forecast lead time using SDM hindcasts in which estimated satellite–rain gauge precipitation fields are used instead of weather forecasts (discussed later). Using estimated precipitation fields isolates hydrologic forecasting errors due to deficiencies in hydrologic modeling and estimated precipitation from errors due to weather forecast uncertainty (which will be combined in step IV). The technique is similar to the stepwise AR technique discussed by Granger and Newbold (1986, their section 5.4) and is given by

 
formula

where Rf(t + f ) = Qo(t + f ) − Qf(t) is the hydrologic model residual (error) at time (day) t + f for a modeled discharge forecast Qf(t) initialized at time t with a lead time of f days in advance, where Qo(t + f ) is the contemporaneous observed discharge at time t + f, μ is the residual time series mean, and ɛt+f is the leftover error from the AR fit.

Note that for the AR model (14) to provide the best linear unbiased estimator of Rf(t + f ), the errors ɛt+f need to be homoscedastic and (generally) independent (see, e.g., Helsel and Hirsch 1992). However, stationarity in the errors is generally not valid (not shown), because when the modeled (and observed) discharges are larger (smaller), the variance of ɛt+f is generally also larger (smaller). Although a Box–Cox transform could be applied in this case to remove such heteroscedastic behavior, the utility of this is not always effective (see, e.g., Chatfield 1996). And more importantly, because (14) is updated and automatically applied to the operational forecasts each day, our motivation is only to make an initial improvement to the SDM, and then in a later step, we utilize a nonparametric nearest-neighbor approach to account for any remaining nonoptimal, correlative, and heteroscedastic behavior in the final hydrologic (multi)model time series forecasts, as discussed in step IV (Fig. 2). This fitting of (14) is carried out using stepwise regression using forward selection based on the minimization of the BIC. In addition, the selected AR model is also compared with up to fifth-order no-fit Taylor series models.

The optimal error correction model for each forecast lead time then “error corrects” the SDM hindcasts exactly as if they were actual forecasts (but driven with observed precipitation and analysis weather variable fields), as well as forecasting the error of (and then correcting) the current ensemble hydrologic forecasts

 
formula

Both the error-corrected hindcasts and corrected ensemble forecasts Qfc( f ) are then sent to the multimodel forecasting step (box M, Fig. 2). It should be noted that because the error correction is applied to the hindcast time series driven by observed fields, it is correcting for the combined hydrologic model and weather observational error and not for the hydrologic error due to forecasted weather variable uncertainty.

4) Multimodel approach

A natural extension of the earlier incorporation of the CLM in 2003 and the SDM in 2004 was to combine the two models under a multimodel framework, capitalizing on the strengths and weaknesses of both approaches as a function of forecast lead time based on the two models’ past performance.

The multimodel approach is automatically recalibrated daily utilizing both the updated CLM and SDM (AR error corrected) hydrologic model hindcasts as follows:

  1. Rank models based on their total hindcast residual error (square difference between modeled and observed discharges).

  2. Perform multiple linear regression of the CLM and SDM hindcast time series against the observed discharge time series to minimize multimodel residual error ɛMM(t + h): 
    formula
    where h is the hindcast “lead time.”
  3. Evaluate resultant residuals of the highest-ranked model and the multimodel using the BIC to avoid overcalibration. If BIC of the multimodel is minimized, then use the resultant regression coefficients {a1, a2} to generate “multimodel” forecasts and hindcasts; otherwise, use the highest-ranked model.

It should be noted that by using hindcasts derived from observed weather fields, as opposed to the forecasted fields, we are trying to isolate and minimize hydrologic model error in the multimodel process. Multimodel discharge forecasts are generated by applying the regression equations to the same ensemble member of each hydrologic model achieved by tracking each ECMWF ensemble member through both CLM and SDM hydrologic models and then individually combining them into a multimodel ensemble member.

The multimodel regression coefficients are essentially weighting coefficients and, as such, can be used to evaluate the skill of each hydrologic forecast model relative to one another. Figure 5 plots the actual regression coefficients (or weights) applied to each of the two hydrologic models (CLM, red; SDM, blue) as a function of forecast lead time for the Ganges (upper panel) and the Brahmaputra (lower panel) for a single day (15 October 2003). Because the Ganges catchment has a longer rainfall-response time scale than the Brahmaputra,3 for the particular meteorological and hydrological state of 15 October 2003, the CLM outperforms the SDM for the Ganges for short lead times. However, the performance reverses at longer lead times as the importance of modeling the longer time-scale groundwater response captured in the physics of the SDM begins to dominate. For the Brahmaputra (lower panel), the CLM outperforms the SDM at all lead times, requiring further investigation. Figure 4 gives an example of one argument why the CLM has shown good skill for severe flooding events. As discussed previously, the CLM should become more effective as large-magnitude precipitation events’ spatial scale approaches that of the catchment itself, as seen preceding the Brahmaputra 26 July–6 August flooding event, especially evident in the ECMWF forecasts of the top panel.

Fig. 5.

Weighting (regression) coefficients for the lumped and distributed models applied in the multimodel scheme as a function of forecast lead time here for 15 Oct 2003. Ganges and Brahmaputra weights for the CLM (red) and the SDM (blue) are shown in the upper and lower panels, respectively. Weighting factors will change depending on the meteorological and hydrological situation.

Fig. 5.

Weighting (regression) coefficients for the lumped and distributed models applied in the multimodel scheme as a function of forecast lead time here for 15 Oct 2003. Ganges and Brahmaputra weights for the CLM (red) and the SDM (blue) are shown in the upper and lower panels, respectively. Weighting factors will change depending on the meteorological and hydrological situation.

d. Step IV: Generation of probabilistic forecasts (gridded box, Fig. 2)

We have argued that by driving the hydrologic multimodel with calibrated ECMWF ensemble weather forecasts, we can project the effect of weather forecast uncertainty onto hydrologic forecast uncertainty. How important are the remaining sources of uncertainty?

Figures 6a and 6b show the operational multimodel discharge forecasts using all 51 ECMWF EPS members (colored lines) compared to observations (solid black lines) for the Brahmaputra (year 2004) and Ganges catchments (year 2006), respectively. Notice how in the day-1 lead-time forecasts for both of these catchments [Figs. 6a(i) and 6b(i), respectively], it appears that there is a single forecast, even though all 51 ensemble members are plotted, showing that numerical weather prediction (NWP) uncertainty provides a negligible source of discharge forecast uncertainty. However, as the lead times increase to 4, 7, and 10 days, the ensemble bundle widens further, showing that NWP uncertainty contributes more and more uncertainty to the hydrologic forecasts for longer and longer lead times. Contrasting the same lead-time Brahmaputra and Ganges ensembles dispersions, notice how the ensemble bundle is much narrower for the Ganges as a result of its longer transportation time for water to migrate through the catchment.

If all sources of uncertainty were properly accounted for, then the spread of the discharge ensemble would provide a representation of forecast uncertainty, and the observed discharge would generally fall within the ensemble bundle.4 However, in Fig. 6, the observed discharge is noticeably outside the ensemble bundle for a large fraction of days for both river basins, in particular for the Ganges. This points out how NWP forecast ensemble spread is not sufficient to account for discharge forecast uncertainty; other hydrologic modeling errors are of significance (e.g., watershed initial state uncertainties, hydrologic model deficiencies) and clearly they need to be accounted for, and their relative significance compared to NWP uncertainty varies for different lead times and catchments.

Our approach to optimally account for all sources of forecast uncertainty is to separate the effect of weather forecast uncertainty on hydrologic forecast uncertainty from that of other hydrologic modeling and initial state errors. This retains the useful information in the dynamical weather forecast dispersion. Hindcast errors of the remaining sources of uncertainty are then generated and conditionally sampled based on hydrologic conditions at the time of forecast initialization, which also provide an error correction of any residual forecast bias. These two sets of uncertainties are then combined to create a weighted ensemble set of hydrologic forecasts that provide a representation of the expected total hydrologic forecast error. These steps are represented schematically in Fig. 7, which will be discussed in more detail later.

Fig. 7.

Combining weather forecast uncertainty with hydrologic uncertainty to create the final set of discharge ensembles accounting for all aspect for discharge forecast uncertainty. (step I) Shown are seven ensemble members (colored bars). Each has an equal probability of , as all ensemble members are equally likely to occur. (step II) The residual hindcast time series, with the dashed segments being selected analogs using the Mahalanobis distance method. These analogs then generate an ensemble of forecasted model errors to the right of the forecast horizon. This set is weighted to produce a (step III) hydrological error PDF that is mapped onto the meteorologically derived discharge set of step I to provide an adjusted set of step IV ensemble members and a broader PDF.

Fig. 7.

Combining weather forecast uncertainty with hydrologic uncertainty to create the final set of discharge ensembles accounting for all aspect for discharge forecast uncertainty. (step I) Shown are seven ensemble members (colored bars). Each has an equal probability of , as all ensemble members are equally likely to occur. (step II) The residual hindcast time series, with the dashed segments being selected analogs using the Mahalanobis distance method. These analogs then generate an ensemble of forecasted model errors to the right of the forecast horizon. This set is weighted to produce a (step III) hydrological error PDF that is mapped onto the meteorologically derived discharge set of step I to provide an adjusted set of step IV ensemble members and a broader PDF.

For a given forecast day with a given current observed discharge, state of the catchment, and weather forecast, symbolically what we seek is a final CDF of the future discharge Qf , at a given forecast lead time f, that will tell us what the probability is that the river’s discharge will exceed a severe flooding level:

 
formula

where FE(Qfqf) gives the (cumulative) probability that the future discharge Qf will be less than the specific value qf , F1(Qfqf|qm) represents the (cumulative) probability that Qfqf given the discharge forecast model forecasts the value qm with no weather forecast uncertainty, and F2(Qmqm|pf) gives the (cumulative) probability that the modeled discharge Qm will be less than the specific value of qm given the NWP model forecasts the specific weather (precipitation) fields pf . In this representation, we have then decomposed the future discharge total uncertainty into two parts: one part due to hydrologic uncertainty (F1) and second part of the total error due to weather (precipitation) forecast uncertainty (F2).

Although the uncertainty calculation given by Eq. (16) in theory represents integrations over continuous functions, in practice, the determinations of both F1 and F2 are performed here by first discretely representing their associated PDFs ( f1 and f2, respectively) by (weighted) ensembles and then discretely convolving the two PDF’s together:

 
formula

where π represents a parameter set for which f1 is dependent on (which includes qm itself), qmi in the last equality is the discharge derived by running each weather forecast ensemble pfi through the hydrologic model, and is a set of conditionally selected hindcast errors (ɛjj = qojqHj) between observed qo and model hindcast discharge qH (derived using observed weather fields). The derivation of the final CDF FE (represented by the integration ) is done by a simple linear interpolation of the cumulative probability between each ensemble member. The specific details of this process are described next with reference to Fig. 7, with further details in Hopson (2005).

  • I: The first step is to generate the PDF F2 of Eq. (17) that represents the total error due to forecast weather uncertainty, which is done by running the adjusted 51-member ECMWF ensemble weather forecasts through the discharge (multi)model. The transformed ensembles represent the discretized PDF version of F2 that we seek, separating hydrologic uncertainty from the uncertainty due to probabilistic weather forecasts. These ensembles are shown schematically in the first block of Fig. 7, where each discharge ensemble member Qp is equally weighted with probability (although only seven ensembles are shown here). This discretized PDF is saved and is combined with F1 in Fig. 7 (IV) described later.

  • IIa: The first step in producing F1 is to generate an historic simulated discharge time series using the current updated hydrologic multimodel for each forecast lead time in an identical fashion to how the forecasts are made, except using “observed” satellite–rain gauge precipitation in place of weather forecasts. By doing so, we separate weather forecast variable uncertainty from the hydrologic uncertainty. This multimodel time series is then differenced with the historic rating curve–derived discharge observation record to generate the desired multimodel residual time series shown. Under the assumption that these historic errors (a nonparametric error “prior”) are statistically representative of future errors,5 this set of residual errors empirically accounts for all forecasting errors except for errors in the weather variable forecasts.6

  • IIb: Because of this dependence of hydrologic uncertainty on a number of measurable factors, the next step is to make a refinement on the hydrologic error prior by making a conditional selection of the hindcasted error ensembles that occurred within a similar region of “state space” as that for the system at the time of forecast initiation.

    To make this conditional selection of past error, we first must specify the “state” of the system (i.e., hydrologic error). In this application, we determine this state by a simple set of selector variables: the last observed hydrologic error eH (observed minus hindcast discharge), along with the slope (with respect to time) s(eH), and curvature c(eH) of this error, at the time of forecast (hindcast) initiation (the slope and curvature are discretely calculated).7

    With this set of state variables, steps Figs. 7(IIc)–7(IIf) below describe the conditional selection procedure, in reference to Fig. 7(II), where the selection (feature vector) in this figure is derived from the hindcasted hydrologic error time series.

  • IIc: At the time of forecast initiation t, with forecast horizon f (1 day in advance, 2 days, … , 10 days), a “feature” vector of the variables above is compiled: 
    formula
  • IId: For all periods i in the historic record, feature vectors of the observed variables are created: 
    formula
    where here, i represents the hindcast initialization day.
  • IIe: The next step is to compare the current (and forecasted) values of the feature vector (VE) with the feature vectors in the historic record (Vi), and those vectors (i.e., days i) most similar to the current conditions for which they are selected. The measure used in the selection process is the Mahalanobis distance Mi (see Davis 1986; Yates et al. 2003 for other applications). Specifically, we calculate 
    formula
    for each day i, selecting out the subset of days for which Mi is lowest, where 𝗦−1 is the inverse of the covariance matrix of the selector variables themselves. Use of Mi diminishes the weighting of the selector variables based on the value of their correlation.
  • IIf: The subset of historic days for which Mi are lowest is then selected out. The number taken is n = N, where N is the total number of historic days available. Referring to the visual schematic of Fig. 7(II), here the most recent residual string (the piece of the residual time series that precedes the solid vertical line of the forecast horizon) is compared with the rest of the historical residual error time series via their value, slope, and curvature, and “n nearest neighbors” for which they are selected. The arrows highlight the residual strings (dashed) that follow favorable comparisons in the historical record. These dashed residual time series then represent forecasted error estimates of the hydrologic error of the current hydrologic forecast, and as such, they forecast the hydrologic error into the forecast horizon, as shown, where the dashed forecasted model residuals are appended to the most recent residual error at the point of the forecast horizon, providing an ensemble estimate of hydrologic model error.

  • IIIa: This subset n of days is then weighted based on their value of Mi, where their weighting Wi is given by 
    formula
  • IIIb: From this weighted set of hindcast initialization days i, the associated hindcasted hydrologic error that occurs at day i + f days is then assigned the same weighting. The result is that we now have constructed a forecast of the estimated PDF of the hydrologic error that will occur at day t + f.

  • IV: The final step is to convolve the forecasted PDF of hydrologic error (steps II and III) with the discharge forecasts derived from weather forecast uncertainty (step I). This is done by mapping the inverse of the hydrologic error PDF onto the PDF of forecast input uncertainty ( f2); this is done by subtracting the value of each hydrologic error ensemble from each weather forecast ensemble–derived discharge ensemble Qp, and we then end up with a new set of ensembles with size [51 × n], and each new ensemble has a probability given by [(1/51) × Wi]. This set of ensembles represents a combined PDF of both weather forecast uncertainty along with hydrologic error uncertainty, and it also has simultaneously provided an additional error correction to the forecast.

To see how this last step provides a correction, recall that the hydrologic error ensembles represent error forecasts. So, we want to remove this error forecast from the final discharge forecast, and in essence, the hydrologic error PDF then represents a final hydrologic error correction. From this PDF, the final discharge forecast CDF FE is calculated.

3. Results

Figures 8 and 9 and the online supplemental figures (available at http://dx.doi.org/10.1175/2009JHM1006.s1) display discharge forecasts for the Brahmaputra at Bahadurabad and the Ganges at Hardinge Bridge based on the forecast system described in section 2. All forecasts shown are operational forecasts disseminated within the respective years with the following caveats: 1) all rating curves have been objectively recomputed and the forecasts rerun using the new rating curve–derived discharges; and 2) the operational forecasts for 2003 and part of 2004 were based only on the CLM model, they did not utilize the TRMM and CMORPH precipitation estimates, and they did not include the final error accounting algorithm of step IV, but they have been rerun to include these features as shown in the figures.8

Of the years of operational forecast dissemination discussed in this paper (2003–079), the Brahmaputra reached flood stage five times, in July 2003 (3–6; 9–17), July 2004 (10–26), and two distinct periods in 2007 (26 July 26–6 August; 8–17 September), with the latter shown in Fig. 8. The Ganges, on the other hand, only briefly exceeded critical flood level once in 2003 (19–21 September), shown in Fig. 9, but extensive flooding never occurred. All other operational years are given in the figures in the online supplement. In all these figures, the left column (I) of each panel summarizes the 5-day discharge forecasts, and the right (II) shows the 10-day discharge forecasts (Note that 1- to 10-day forecasts were operationally produced and disseminated but only 5- and 10-day forecasts are shown here for brevity). For presentation purposes, the forecasts have been sorted, removing the useful temporal covariability of the ensemble traces. The dates along the abscissa denote the forecast day. The bold black line shows the observed discharge verification curve. The distance between each plot’s right edge and the vertical black line is a metric of the forecasting lead time.

In Figs. 8 and 9, and those in the electronic appendix, the discharge ensembles [row (i)] were used to calculate the forecast 95% (red) and 50% (gold) confidence intervals [row (ii)]. Perhaps the most important metric is the probability that the rivers’ flow into Bangladesh will exceed its respective danger level (corresponding to discharges of approximately 62 000 and 55 000 m3 s−1 for the Brahmaputra and Ganges, respectively). These 5- and 10-day lead-time forecasted danger-level exceedance probabilities are shown in row (iii) for each river. These curves correspond in essence to the last step of the technical forecasting process (step V of Fig. 2).

To highlight further the error growth with lead time, in Fig. 10 we zoom in on the 10-day lead-time forecasts during the 2004 severe flooding event and isolate two days for which operational forecasts were issued (although forecasts were issued every day) during periods of extreme flood onset and recession. Shown in Figs. 10(v) and 10(vi) are plume plots of the ensemble forecasts and derived danger-level exceedance probabilities from two separate issuing dates out to the end lead time of 10 days. Vertical black lines bound the severe flooding periods. The two July-issued forecasts shown were chosen specifically for days when the observed hydrograph was decreasing on the issuing date (for which a simple “trend” forecast would have no skill). The first forecast [Fig. 10(v)] was issued 6 days before the discharge crossed the critical threshold and was forecasting 31% likelihood of exceedance for that date, but it was forecasting 66% likelihood at 7 days, and 84%–87% at 8–10 days. The second plume forecast shown [Fig. 10(vi)] points outs the forecast system’s ability to estimate flood duration. Again, as the observed hydrograph was descending, the forecasts instead warned of a second flood peak (maintaining severe flood probabilities of 81%–100% during this above-critical-level period).

Overall, both rivers’ discharge was well forecast at the 5-day and 10-day limit for the Brahmaputra but appeared to be conservative in forecasting the precise date of crossing danger-level thresholds. For the Brahmaputra, the flood forecasts of 2003 (refer to the electronic appendix) showed only moderate skill. However, the major flooding in 2004 was well forecast at 10 days, including the double peak and the duration of the forecast (Fig. 10). In particular, using the 10-day forecast lead-time metric (the gap between the vertical line and the plots’ right edge) in Fig. 10, we note that the first July 2004 flood wave was being forecast as the hydrograph was rapidly descending (continuously for eight days), implying a strictly discharge time series forecast model would not have captured this event. No floods were forecast during 2006 and none occurred. For the first major flood wave of 2007 (Fig. 8), the forecast amplitude was lower than occurred, but the timing and forecast critical-level probabilities were skillful. The second 2007 flood wave was generally well forecast. Overall, except for moderate skill in 2003, all major flooding events were well forecast at 10-day lead times with high probability of occurrence and with well-captured peak amplitudes and dates of flood recession. Notably, during the 5-yr period, there were no false positives of major flooding for either basin.

For the Ganges, the forecasts were generally less skillful than the Brahmaputra, although the modeling system has not been tested yet for severe extended flooding events. Low probabilities of short-duration flooding occurred throughout the forecast record, but no high probability false-positive flooding events were forecast. Overall though, forecast skill begins to lose utility past five days and as discussed previously, one reason for this may be that the forecast system presented here does not explicitly account for the Ganges’s water management (unlike the generally “naturalized” Brahmaputra) despite the Jian et al. (2009) study. Current improvements to the forecast system are being implemented to account for this management, a topic of future work.

All sources of discharge forecasting uncertainty were taken into account in each of the forecasts shown in Figs. 8 –10. Modeling the effect of weather forecast uncertainty, each of the 51-member ECMWF ensemble precipitation forecasts was run through the multimodel (using 1-, 2-, … , 5-day precipitation forecasts for the 5-day discharge forecasts and using 1-, 2-, … , 10-day precipitation forecasts for the 10-day discharge forecasts) and then merged with hydrologic and weather observation uncertainty, as discussed in step IV. Although not shown, there is a great deal of uncertainty—as depicted by the ensemble spread—in even the 1-day precipitation forecasts. However, after the transformation through the modeled Brahmaputra and Ganges watersheds with their basin time scales of order 5 days and 13 days, respectively, 1-day forecasts of precipitation only weakly affect the 1-day discharge forecasts. However, this is not the case for the 10-day discharge forecasts, where the weather forecast’s uncertainty clearly affects the discharge uncertainty. The probabilistic forecasts, so produced, provide the necessary information to decide if there is “sufficient risk” upon which to take action.

Objectively, how successful was the forecast system at providing skillful, calibrated forecasts? The rank histograms shown in Fig. 11 (Brahmaputra only; Ganges not shown, but they give similar results) provide one way to assess the calibration of the ensemble forecasts’ own uncertainty estimates. These panels were generated for the final calibrated 51-member ensemble (with 52 intervals) 5- and 10-day forecasts. The red dashed lines are 95% confidence intervals (CIs) on the expected counts in each bin for a perfect model. For 52 intervals, we would expect roughly 2–3 interval counts to fall outside the CI; instead, roughly 5–6 intervals for both 5- and 10-day forecasts for both basins lie outside the CI (with similar results for other forecast lead times). Although the forecast calibration is not that of a perfect model (at 95% confidence), overall the rank histograms are close to being statistically uniform.

Fig. 11.

Rank histograms of the 51-member (52-interval) Brahmaputra forecasts. Red dotted lines show 95% confidence limits for a perfectly calibrated forecast.

Fig. 11.

Rank histograms of the 51-member (52-interval) Brahmaputra forecasts. Red dotted lines show 95% confidence limits for a perfectly calibrated forecast.

The top panel of Fig. 12 provides an objective Brier scores (BSs) assessment of forecasting skill for the event of Brahmaputra discharge exceeding critical level, corresponding to a discharge of 62 400 m3 s−1, or the 94th percentile of the climatological record (q0.94). (BSs of the Ganges are not provided because only three days in the forecasting period exceeded critical level.) The results are presented as a function of forecast lead time for the forecast BSs (solid black) and 95% CI (black dash), along with a persistence forecast BS CI (red dash) and a climatological forecast BS CI (green dash). A persistence forecast refers to using the river discharge in day t to forecast the river discharge in day t + f. For the BS, a climatological forecast is given by the probability to exceed the river stage danger level on any given day during the monsoon season (calculated to be 0.06 for the Brahmaputra). The CIs are generated by “bootstrapping” (Huber 1981) the days (and the associated forecasts and observations) in the forecast period (years 2003–07) but preserving the temporal covariance of the time series (which widens the CI). Forecast skill analysis is performed on each bootstrap sample (here, performed 2000 times) and the results are ranked to determine the values corresponding to the 2.5 and 97.5 percentiles (giving the 95% CI). The discharge forecasts are superior to both climatology and persistence at all lead times, with only minimal-to-moderate overlap of their CI. [If one calculates the conditional ranking of the BS for each bootstrap sample (not shown), then there is more than 99% certainty the forecasts are superior to both.]

Fig. 12.

(top) BS for the 1–10-day Brahmaputra forecasts as a function of lead time. The scores are derived for the above-critical-level threshold (climatological 94th percentile). Solid black line is the forecast BS, with black dashed lines giving the 95% CIs on the results (see text). Green dashed lines give 95% CI for climatology (q = 0.94), and red dashed lines for persistence. Notice the most likely forecast BS (solid black) is essentially superior to both climatology and persistence at all lead times with high confidence. (bottom) The same as (top) but for BSSs. The solid green line gives the BSS for a climatological (q = 0.94) reference forecast, and solid red line applies to a persistence reference forecast. Dashed lines give 95% confidence limits. Notice the most likely BSSs (solid lines) show the forecasts outcompete both persistence and climatology at all lead times at >60% improvement over the respective reference forecast.

Fig. 12.

(top) BS for the 1–10-day Brahmaputra forecasts as a function of lead time. The scores are derived for the above-critical-level threshold (climatological 94th percentile). Solid black line is the forecast BS, with black dashed lines giving the 95% CIs on the results (see text). Green dashed lines give 95% CI for climatology (q = 0.94), and red dashed lines for persistence. Notice the most likely forecast BS (solid black) is essentially superior to both climatology and persistence at all lead times with high confidence. (bottom) The same as (top) but for BSSs. The solid green line gives the BSS for a climatological (q = 0.94) reference forecast, and solid red line applies to a persistence reference forecast. Dashed lines give 95% confidence limits. Notice the most likely BSSs (solid lines) show the forecasts outcompete both persistence and climatology at all lead times at >60% improvement over the respective reference forecast.

We further examine these comparisons through Brier skill scores (BSS) in the bottom panel of Fig. 12, where the BSS is calculated as

 
formula

where BSf is the Brier score for the forecasts, BSref is for a “reference forecast” [here, either “climatology” (green) or “persistence” (red)], and BSperf for a “perfect” forecast. The most likely BSSs (solid lines) show the forecasts outcompete both persistence (red solid and dashed CI) and climatology (green solid and dashed CI) at all lead times at >60% improvement over these respective reference forecasts.

Further skill score assessments were done for both the Brahmaputra and Ganges using other measures. Figures 13 and 14 show skill scores for the root-mean-square error of the ensemble mean forecasts [where RMSE is substituted for BS in Eq. (23)] and continuous rank probability skill score of the full ensemble [CRPSS substituted for BS in Eq. (23)] of the Brahmaputra (top panels) and Ganges (bottom panels) river discharge forecasts, respectively. The persistence forecast referred to in these figures is as defined earlier. The RMSE climatological forecast is the mean of the monsoon season historical discharge values, whereas the CRPS climatological forecast is given by the whole PDF of the same (i.e., not just the mean value). Again, the Brahmaputra forecasts (top panels of both figures) significantly outperform both climatology (45%–98% improvements) and persistence (55%–80%). The Ganges forecasts (bottom panels, both figures) also show significant skill gains over climatology (60%–98% improvements) and persistence, except for the RMSE of the 24-h Ganges forecast where we see only slight improvements over persistence. [However, the conditional ranking of the RMSE for each bootstrap sample (not shown) gives >67% certainty of the 24-h forecast outperforming persistence.]

Fig. 13.

RMSE skill scores and 95% confidence intervals (dashed) of the (top) Brahmaputra and (bottom) Ganges ensemble mean forecasts compared to persistence (red) and the climatological mean (green). Results show that for all lead times, the Brahmaputra forecasts outperform both persistence and climatology at >60% improvement, but only after 1 day do the Ganges forecasts significantly outperform persistence (>40% improvement).

Fig. 13.

RMSE skill scores and 95% confidence intervals (dashed) of the (top) Brahmaputra and (bottom) Ganges ensemble mean forecasts compared to persistence (red) and the climatological mean (green). Results show that for all lead times, the Brahmaputra forecasts outperform both persistence and climatology at >60% improvement, but only after 1 day do the Ganges forecasts significantly outperform persistence (>40% improvement).

Fig. 14.

CRPSS and 95% CIs (dashed) of the (top) Brahmaputra and (bottom) Ganges ensemble forecasts compared to persistence (red) and climatological CDF (green). Results show that for all lead times, forecasts for both catchments outperform both reference forecasts at high confidence (persistence and climatology) at >40% improvement.

Fig. 14.

CRPSS and 95% CIs (dashed) of the (top) Brahmaputra and (bottom) Ganges ensemble forecasts compared to persistence (red) and climatological CDF (green). Results show that for all lead times, forecasts for both catchments outperform both reference forecasts at high confidence (persistence and climatology) at >40% improvement.

4. Conclusions

In this work, we have presented an approach to produce probabilistic river discharge forecasts in the medium range for two large-scale river basins. The ability to produce accurate and timely forecasts at these spatial scales depends primarily on three factors: the availability of quality ensemble weather predictions [produced at a few centers, such as ECMWF, NOAA/National Centers for Environmental Prediction (NCEP), etc.], the availability of a range of satellite precipitation products (TRMM, CMORPH, etc.) to initialize modeled watershed states and to reduce systematic error in the meteorological forecasts, and a method of incorporating these data sources into a hydrological modeling system. In addition, an infrastructure to disseminate and utilize the forecasts is essential. The CFAB methodology presents an operational system that bridges the gap between the production of meteorological observations and forecasts, and the delivery of a useable hydrological product.

The CFAB system was developed in its present configuration because of the absence of upstream streamflow data beyond the borders of Bangladesh. The system depends critically on current and historical discharge data at the India–Bangladesh border for calibration and forecast data assimilation. The quality of the forecasts shown in Figs. 8 –10 may speak to the utility of similar schemes applied to other regions where upstream data may also be inaccessible or just not measured. However, if river flow measurements higher up in the catchment were available and could be routed downriver to the forecast location, errors in rainfall–runoff modeling (due to imprecise satellite-derived precipitation estimates and hydrologic modeling) could be reduced, potentially leading to significant improvements in forecasts at shorter lead times. However, hydrologic runoff process modeling coupled to statistically rendered weather forecast information (as per step II in Fig. 2) would play an increasingly greater role (and in situ river data less so) as the forecast lead times increase, or as catchment flow-through time scales decrease.

CFAB began producing operational 1–10-day advance operational flood forecasts in 2003, with the system presented here completed in 2004. During July 2004, the system produced skillful real-time forecasts of the severe Brahmaputra flooding event, motivating the establishment of a pilot program in 2006 to disseminate the forecasts down to vulnerable citizens in five selected districts with a total population of approximately 110 000 (Fig. 1, right panel; CEGIS 2006; Webster et al. 2010). During 2007, the Brahmaputra flooded severely twice, and the dissemination network was activated to warn vulnerable citizens, leading to the early evacuation of many thousands of people and livestock. As tragic as these floods were to the country as a whole, having the forecast system described here be of service to the people of Bangladesh was tremendously gratifying.

Acknowledgments

ECMWF has enthusiastically supported CFAB since its inception. In particular, we thank Philippe Bougeault, Roberto Buizza, the late Tony Hollingsworth, Dominique Marbouty, Martin Miller, and Tim Palmer for specific assistance. The Asian Disaster Preparedness Centre has directed the in-country work, and we are especially indebted to Sri. A. R. Subbiah and Sri. S. Ramasamy for their inspirational enthusiasm and dedication in helping the poorest of the poor in South Asia. We are also appreciative of the partnership with the Flood Forecasting and Warning Centre and the Bangladesh Meteorological Department, along with the Center for Environmental and Geographic Information Services, the Institute for Water Modelling, and the many other stakeholder institutions within Bangladesh for supporting this project. We would also like to thank Keith Beven for his very valuable assistance in early 2003 as the project went operational, and thanks to the inspiration of Gilbert White, whose life’s work focused on scientific research meeting the needs of the societal greater good. Initial support for the project was given by USAID under Grant USAID/OFDA AOT-A-00-00-00262-00 (2000–03) and most recently by CARE Bangladesh (2005–07). Basic science has been supported by National Science Foundation Grant ATM-0531771. We also appreciate the support of Georgia Tech’s School of Earth and Atmospheric Sciences, which provided the means for the project to continue during 2005, and the Advanced Study Program and the Research Applications Laboratory of the National Center for Atmospheric Research.

REFERENCES

REFERENCES
Ahmad
,
A. U.
, and
M. M. Q.
Mirza
,
2000
:
Review of causes and dimensions of floods with particular reference to Flood ‘98: National perspectives.
Perspectives on Flood 1998, Q. K. Ahmad et al., Eds., University Press Limited, 67–84
.
Akaike
,
H.
,
1974
:
A new look at the statistical model identification.
IEEE Trans. Autom. Control
,
19
,
716
723
.
Beven
,
K. J.
,
2000
:
Rainfall-Runoff Modelling: The Primer.
John Wiley and Sons, Ltd., 360 pp
.
Box
,
G. E. P.
, and
G. M.
Jenkins
,
1970
:
Time Series Analysis: Forecasting and Control.
Holden-Day, Inc., 553 pp
.
Buizza
,
R.
,
J-R.
Bidlot
,
N.
Wedi
,
M.
Fuentes
,
M.
Hamrud
,
G.
Holt
, and
F.
Vitart
,
2007
:
The new ECMWF VAREPS (Variable Resolution Ensemble Prediction System).
Quart. J. Roy. Meteor. Soc.
,
133
,
681
695
.
Calheiros
,
R. V.
, and
I. I.
Zawadzki
,
1987
:
Reflectivity–rain rate relationships for radar hydrology in Brazil.
J. Climate Appl. Meteor.
,
26
,
118
132
.
CEGIS
,
2006
:
Sustainable end-to-end climate/flood forecast application through pilot projects showing measurable improvements.
CEGIS Base Line Rep., 78 pp
.
Central Water Commission
,
1988
:
Water resources of India.
CWC Publication 30/88
.
Chatfield
,
C.
,
1996
:
The Analysis of Time Series: An Introduction.
Chapman & Hall, 283 pp
.
Chow
,
V. T.
,
D. R.
Maidment
, and
L. W.
Mays
,
1988
:
Applied Hydrology.
McGraw-Hill, 572 pp
.
Davis
,
J. C.
,
1986
:
Statistics and Data Analysis in Geology.
2nd ed. John Wiley, 646 pp
.
Georgakakos
,
K. P.
,
D-J.
Seo
,
H.
Gupta
,
J.
Schaake
, and
M. B.
Butts
,
2004
:
Towards the characterization of streamflow simulation uncertainty through multimodel ensembles.
J. Hydrol.
,
298
,
222
241
.
Granger
,
C. W. J.
, and
P.
Newbold
,
1986
:
Forecasting Economic Time Series.
2nd ed. Academic Press, 338 pp
.
Haltiner
,
J. P.
, and
J. D.
Salas
,
1988
:
Short-term forecasting of snowmelt runoff using ARMAX models.
Water Resour. Bull.
,
24
,
1083
1089
.
Helsel
,
D. R.
, and
R. M.
Hirsch
,
1992
:
Statistical methods in water resources.
Statistical Methods in Water Resources, Studies in Environmental Science Series, Vol. 49, Elsevier, 522 pp
.
Hopson
,
T. M.
,
2005
:
Operational flood-forecasting for Bangladesh. Ph.D. thesis, University of Colorado, 225 pp
.
Huber
,
P. J.
,
1981
:
Robust Statistics.
Wiley, 308 pp
.
Huffman
,
G. J.
,
R. F.
Adler
,
S.
Curtis
,
D. T.
Bolvin
, and
E. J.
Nelkin
,
2005
:
Global rainfall analyses at monthly and 3-hr time scales.
Measuring Precipitation from Space: EURAINSAT and the Future, V. Levizzani, P. Bauer, and J. F. Turk, Eds., Springer, 722 pp
.
Huffman
,
G. J.
, and
Coauthors
,
2007
:
The TRMM Multisatellite Precipitation Analysis (TMPA): Quasi-global, multiyear, combined-sensor precipitation estimates at fine scales.
J. Hydrometeor.
,
8
,
38
55
.
Jian
,
J.
,
P. J.
Webster
, and
C. D.
Hoyos
,
2009
:
Large-scale controls on Ganges and Brahmaputra river discharge on intraseasonal and seasonal time-scales.
Quart. J. Roy. Meteor. Soc.
,
135
(
639
)
353
370
.
Joyce
,
R. J.
,
J. E.
Janowiak
,
P. A.
Arkin
, and
P.
Xie
,
2004
:
CMORPH: A method that produces global precipitation estimates from passive microwave and infrared data at high spatial and temporal resolution.
J. Hydrometeor.
,
5
,
487
503
.
Katz
,
R. W.
,
1981
:
On some criteria for estimating the order of a Markov chain.
Technometrics
,
23
,
243
249
.
Krishnamurti
,
T. N.
,
C. M.
Kishtawal
,
T. E.
LaRow
,
D. R.
Bachiochi
,
Z.
Zhang
,
C. E.
Williford
,
S.
Gadgil
, and
S.
Surendran
,
1999
:
Improved Weather and Seasonal Climate Forecasts from Multimodel Superensemble.
Science
,
285
,
1548
1550
.
Krishnamurti
,
T. N.
,
C. M.
Kishtawal
,
Z.
Zhang
,
T.
LaRow
,
D.
Bachiochi
, and
E.
Williford
,
2000
:
Multimodel ensemble forecasts for weather and seasonal climate.
J. Climate
,
13
,
4196
4216
.
Lees
,
M.
,
P.
Young
,
S.
Ferguson
,
K.
Beven
, and
J.
Burns
,
1994
:
An adaptive flood warning scheme for the River Nith at Dumfries.
Second International Conference on River Flood Hydraulics, W. R. White and J. Watts, Eds., Wiley and Sons, 65–75
.
Marco
,
J. B.
,
R.
Harboe
, and
J. D.
Salas
,
Eds.
1993
:
Stochastic Hydrology and Its Use in Water Resources Systems Simulation and Optimization.
NATO Advanced Study Institute Series, Vol. 237, Kluwer Academic Publishers, 483 pp
.
Mirza
,
M. M. Q.
,
2003
:
The Choice of Stage-Discharge Relationship for the Ganges and Brahmaputra Rivers in Bangladesh.
Nord. Hydrol.
,
34
,
321
342
.
Mirza
,
M. M. Q.
,
2004
:
The Ganges Water Diversion: Environmental Effects and Implications.
Water Science and Technology Library Series, Vol. 49, Kluwer Academic Publishers, 364 pp
.
Molteni
,
F.
,
R.
Buizza
,
T. N.
Palmer
, and
T.
Petroliagis
,
1996
:
The ECMWF Ensemble Prediction System: Methodology and validation.
Quart. J. Roy. Meteor. Soc.
,
122
,
73
119
.
Nash
,
J. E.
,
1958
:
The form of the instantaneous unit hydrograph.
Surface Water, Prevision, Evaporation, IASH Publication 45, Vol. 3–4, IAHS, 114–121
.
Nelder
,
J. A.
, and
R.
Mead
,
1965
:
A simplex method for function minimization.
Comput. J.
,
7
,
308
313
.
Press
,
W. H.
,
S. A.
Teukolsky
,
W. T.
Vetterling
, and
B. P.
Flannery
,
1992
:
Numerical Recipes in FORTRAN: The Art of Scientific Computing.
Cambridge University Press, 963 pp
.
Riverside Technology Inc. and EGIS
,
2000
:
Main report.
Vol. 1, Information for flood management in Bangladesh, Report for CARE-Bangladesh, 131 pp
.
Schwarz
,
G.
,
1978
:
Estimating the dimension of a model.
Ann. Stat.
,
6
,
461
464
.
Sherman
,
L. K.
,
1932
:
Streamflow from rainfall by unit-graph method.
Eng. News-Rec.
,
108
,
501
505
.
Singh
,
V. P.
,
N.
Sharma
, and
C. S. P.
Ojha
,
2004
:
Brahmaputra Basin Water Resources.
Water Science and Technology Library Series, Vol. 47, Kluwer Academic Publishers, 632 pp
.
Taylor
,
J. R.
,
1982
:
An Introduction to Error Analysis: The Study of Uncertainties in Physical Measurements.
University Science Books, 270 pp
.
Thornthwaite
,
C. W.
, and
B.
Holzman
,
1939
:
The determination of evaporation from land and water surfaces.
Mon. Wea. Rev.
,
67
,
4
11
.
Webster
,
P. J.
,
T.
Hopson
,
C.
Hoyos
,
A.
Subbiah
,
H-R.
Chang
, and
R.
Grossman
,
2006
:
A three-tier overlapping prediction scheme: Tools for strategic and tactical decisions in the developing world.
Predictability of Weather and Climate, T. N. Palmer and R. Hagedorn, Eds., Cambridge University Press, 645–673
.
Webster
,
P. J.
, and
Coauthors
,
2010
:
Extended-range probabilistic forecasts of Ganges and Brahmaputra floods in Bangladesh.
Bull. Amer. Meteor. Soc.
,
in press
.
Xie
,
P. P.
,
B.
Rudolf
,
U.
Schneider
, and
P. A.
Arkin
,
1996
:
Gauge-based monthly analysis of global land precipitation from 1971 to 1994.
J. Geophys. Res.
,
101
, (
D14
).
19023
19034
.
Yates
,
D.
,
S.
Gangopadhyay
,
B.
Rajagopalan
, and
K.
Strzepek
,
2003
:
A technique for generating regional climate scenarios using a nearest-neighbor algorithm.
Water Resour. Res.
,
39
,
1199
.
doi:10.1029/2002WR001769
.
Young
,
P. C.
,
2002
:
Advances in real-time flood forecasting.
Philos. Trans. Roy. Soc. London
,
A360
,
1433
1450
.
Young
,
P. C.
, and
K. J.
Beven
,
1994
:
Data-based mechanistic modelling and the rainfall-flow non-linearity.
Environmetrics
,
5
,
335
363
.

Footnotes

+ Current affiliation: Advanced Study Program and Research Applications Laboratory, National Center for Atmospheric Research, Boulder, Colorado

Corresponding author address: Thomas M. Hopson, NCAR, P.O. Box 3000, Boulder, CO 80307-3000. Email: hopson@ucar.edu

1

Years for which severe flooding occurred in one or both basins.

2

However, Mirza (2003) concluded for the same gauging locations presented here that year-to-year variation in river stage for any common threshold discharge was minimal.

3

The (calibrated) Ganges SDM upper-layer (ts1), lower-layer (ts2), and unit hydrograph (tc) time constants are 4.2, 20.7, and 1.9 days; for the Brahmaputra, they are 0.3, 8.5, and 0.3 days.

4

Specifically, the observed discharge should have equal probability (here, ) of falling between any two ensembles or above (below) the largest (smallest) ensemble.

5

Because such residual errors are often generated from a “tuned” discharge model, strictly speaking the error estimates should come from cross validation or be inflated by calculated prediction intervals (see, e.g., Helsel and Hirsch 1992).

6

Although not discussed in detail here, we argue this prior also includes an implicit (although under) estimate of discharge rating curve error.

7

Note that there are many other hydrologic states that could also be utilized here, such as forecasted or observed catchment-averaged precipitation, soil moisture, forecast or observed discharge, among others.

8

Also, the Akaike Information Criteria (Akaike 1974) has been replaced with the BIC.

9

Operational forecasts were not made for 2005 as a result of a lapse in funding.

* Supplemental information related to this paper is available at the Journals Online Web site: http://dx.doi.org/10.1175/2009JHM1006.s1.

Supplemental Material