Journal of Japan Society of Civil Engineers, Ser.
B1 (HydraulicEngineering),
Vol.
68, No.
4, I_25-I_30,
2012.
Annual Journalof HydraulicEngineering,
JSCE,
Vol.56,2012, February

SHORT TERM PREDICTION OF WATER LEVELAND DISCHARGE USING A 2D DYNAMIC WAVE
MODEL WITH PARTICLE FILTERS
Yeonsu KIM1, Yasuto TACHIKAWA2, Sunmin KIM2,
Michiharu SHIIBA2, Kazuaki YOROZU2, Seong Jin NOH1
1Student Member of JSCE, Graduate School of Eng., Civil and Earth Resources Eng., Kyoto Univ.
2Member of JSCE, M. Eng., Professor, Dept. of Civil and Earth Resources Eng., Kyoto Univ.

A short-term flood prediction algorithm with a 2D dynamic wave model and a particle filter is
proposed to consider the uncertainties of hydrologic input data and channel roughness. The particle filter
makes it possible to utilize a non-linear and non-Gaussian model for estimating time variant channelroughness and inlet flow uncertainty by considering sequentially updated water stages. The proposed
method was applied to the Katsura River located in Kyoto, Japan, and it was verified first through asynthetic experiment. The experiment result shows that the algorithm successfully traces the hidden true
values, which are the correct inlet discharges and Manning’s roughness coefficient, on a real-time basis.
The prediction results were also compared with observed water stages, and they showed good agreementswith the observed water stage.
Key Words: River discharge, short-term prediction, particle filter, 2D dynamic wave model

1. INTRODUCTION
Flood forecasting is crucial information for
mitigating and/or protecting the damages induced byflood to our properties and human lives. In general,a hydrological model is utilized to predict thevolume and peak time of flood, based on a given
rainfall information. Then, a hydraulic model
calculates water stage profile for a specific river
channel to provide more detailed flood information.
The most successful method of accurate flood
forecasting may need precise rainfall information
and a well-organized hydrologic model, as well as
hydraulic models. However, it is not avoidable tosuffer much erroneous forecasting information
through a series of conversion, which is fromrainfall to discharge and from discharge to waterstage. Thus, we are still lacking in our forecasting
methods, not only in accurate rainfall forecasting
but also in proper hydrologic and hydraulicmodeling.
Among many error sources of forecasting in theabove mentioned processes, this study focuses ondata assimilation into a hydraulic model to properly
convert uncertain river discharge information into
correct water stage information on a specific river

channel. Generally, input data into a hydraulicmodel, such as inflow from the upper boundary,already include a certain amount of error, andadditional system error that comes from a hydraulicmodel will be added during simulations. To avoid
these conventional error sources, this study proposesa new method of hydraulic model utilization, which
considers the input error, system error, andobservation error using a recursive Monte Carlosimulation algorithm.
Real-time updating of model state variables has
already been adopted in many researches withstochastic approaches like the Kalman filter andparticle filters. Shiiba et al.1) introduced the Kalmanfilter into the 1D dynamic wave model to improveforecasting accuracy considering observed waterstage and discharges. Hsu et al. 2) showed that
updating of the channel roughness coefficient during
a simulation considering observed water stage
improves prediction capability, and a fixed channelroughness coefficient results in an inaccurate
prediction. Discharge also includes uncertainties, as
pointed out by Dottori et al.3) and Baldasarre andMontari4) even observed that discharge data include

I_25

Fig. 1 Real-time water stage prediction algorithm proposed in this study

many uncertainties. In addition, discharge data
converted by a rating curve from observed waterstages includes many uncertainties.
Arico et al.5) presented a simultaneous
estimation method for discharge and channelroughness, and proved that it is essential to consider
heterogeneous channel roughness. It is obvious thatchannel roughness, inlet flow, and the interaction ofthese factors are critical items for accurate waterstage estimation. A stochastic model or a real-time
based calibration seems necessary for considering
many errors during hydrologic and hydraulic model
simulations.
Recent researches have introduced a new type ofrecursive updating scheme, called particle filters(PFs), into 1D hydraulic models to consider the
non-linearity of system models (e.g., Montanari etal.6); Matgen et al.7); Giustarini et al.8)). Tachikawa
et al.9) also introduced PFs into the 1D hydraulicmodel, and they improved the predictability of waterstages by using the sequentially updated waterstages. However, the 1D model was not able to fully
consider the geomorphologic characteristics ofchannel, and the tracking ability of water stage intheir study was limited. In terms of the method
incorporating noise, the conventional Kalman filter
algorithm has limitations in its application to
unknown non-normal variances in the state orobservational equation (Bradley et al.10)). However,PFs are applicable in the non-linear system without
any Gaussian assumption.
In this study, we introduce an improved
prediction algorithm that is based on a 2D dynamic
wave model and particle filters considering input,
system, and observation errors. To verify the
proposed algorithm, we tested our method on a shortreach on the Katsura River in Kyoto, Japan. A 2D
dynamic wave model is adopted in this study to
reproduce the relation of water stage, discharge, andparameters on a complex river bed more precisely.
The paper is composed of the following
sections. Section 2 describes the proposed

methodology. Section 3 is devoted to particle filtersand noise generation methods. The presented
algorithm is verified with a synthetic experiment insection 4. Then, in section 5, we implement theperformance analysis of the proposed method.
Finally, section 6 summarizes the study.

2. METHODOLOGY(1) Prediction algorithm
The prediction algorithm largely consists of aparticle filtering system (estimation process) and a
prediction system (Fig.1). The prediction algorithm
runs on results through an estimation process. The
estimation process is composed of a perturbation
step, an update simulation step, and a resamplingstep. First, the boundary condition, such as the
upstream discharge and the downstream water stage,
and model parameter values (channel roughness) are
disturbed to consider their uncertainties at the
perturbation step. Then, the state variables, which
are calculated by the 2D dynamic wave model, andparameter values are sequentially updated in theresampling step according to the weight against thesequentially updated water stage. Next, the statevariables and parameter values are transferred to theprediction process to reflect the current state moreaccurately after updating. The prediction process iscalculated up to 6 hours, with the updated statevariables and parameter values every hour.
The proposed method is applied to the Katsura
River located in Kyoto, Japan. The reach length isabout 2km and the study reach is covered with 500
structure grids in the calculation domain (Fig. 2).
There are two water stage stations at both ends
without tributaries, so the characteristics of thereach are good for application of the method since
we neglect the lateral flow in consideration ofdischarge uncertainties. In addition, the flood plain
exists on both sides of the main channel. The Kamo
River joins at the upstream of the reach and the Uji
River joins at the downstream of the reach.

I_26

h M N0t x yZ M uM vM  gh s  bx  (u  2 h)  (u v h)
t
x
y
x
 x
y
Z  by N uN vN gh s  (u vh)  (v 2 h)
t
x
y
y
 x
y

(1)
(2)

(3)

where h is water depth and Z s is water stage; tis time; u and v are depth averaged velocity in x
and y direction; M and N denote discharge flux inx and y direction;  is density of water; and  bx
and  by are the bed shear stress in x and y

Fig. 2 Location of the study area and the 2D model structure

Although they do not join inside the study reach, theeffect of the tributaries and flood plain makes flow
conditions very complex.
We will firstly verify our proposed method witha synthetic experiment, which is based on artificial
true values. The method is going to be tested on anatural river channel with the observed data to
confirm its predictability. The synthetic experimentconsists of two steps. The first step is to generate―True Values‖ with the existing upstream discharge,
downstream water level, and arbitrarily generated
Manning roughness coefficient. The water stage anddischarge of each cross section are simulated
through this step, and these simulation results are
regarded as ―True Values.‖ The second step is totest the particle filter algorithm—namely, whether
or not the algorithm successfully traces the hidden―True Values.‖ In this testing step, upstream
discharges and downstream water stage, including acertain level of error, are used as the boundary
condition. Then some errors were considered in
Manning roughness coefficient as well. Based on
this condition, the particle filter algorithm with thecorrect water stage information was utilized to trace
the correct discharges and proper Manningroughness coefficient on a real-time basis.
It is believed that a large number of particles
improves the accuracy of the estimation andprediction. The algorithm proposed in this study has
been tested with different particle numbers, such as100, 300, and 500 particles. We have determinedthat 100 particles in our algorithm are sufficient in
the sense of calculation time and simulation
accuracy.

direction.

3. PARTICLE FILTERING SYSTEM
PFs perform the sequential Monte Carlo (SMC)
estimation based on particle representations ofprobability densities within Bayesian theorem(Ristic et al.12)). The sequential process of particles
in time using the nonlinear model is preceded up to
the next available measurement (Salmon andFeyen13)). Among the various PFs, the sequential
importance resampling (SIR) method is introducedin this study since the SIR can reduce themeaningless calculations and estimate current state
more exactly.
(1) Perturbation process
From the assumption that all errors come fromuncertainties of channel roughness and inlet inflow,
we incorporate some errors to Manning’s roughness
coefficient, inlet discharge, and downstream waterstage before the simulation step.
The general method, which adds white noise to
the chosen particles like Eq. (6), shows thelimitations in tracking rapidly changing parametervalues and state variables. Therefore, the modified
extrapolation method, which showed appropriate
ability when tracking abruptly changed channelroughness with Kalman filter in the study of
Crissman et al.14), is introduced to generate the noiseof Manning roughness coefficient and inlet flow, asin Eq. (5) and Eq. (4), respectively.
The polynomial extrapolation method reduces
the variation of particles, but it shows the weakness
in inflection point, such as peak time of dischargehydrograph. As an alternative to these problems, we
introduced the two types of perturbation equations
to generate noise of Manning’s roughness
coefficient. In order to choose the perturbationobs
equation, the discharge ratio ( (Qtobs  Qtobs),
1 ) / Qt(2) 2D Dynamic wave model
The 2D dynamic wave model is composed of thecontinuity equation (Eq. (1)) and the momentum
equation (Eq. (2) and (3)). Then, the equations are
numerically solved with Finite Volume Method(FVM) on a structure grid. In dealing with theconvection term, the simple first order upwind
scheme is utilized. The Adams-Bashforth method is
introduced for time integration (Nagata11)).

where Qtobs is the discharge at the current time and
Qtobs
1

3

I_27

is the discharge at the previous time, isintroduced. Several cases were tested to determine
the criteria to choose the perturbation equation, andit has been confirmed that a 20% discharge ratio isthe proper value in our subject river channel. In thecase when the ratio is less than 20%, Eq. (5) isapplied, and Eq. (6) is applied when the ratio is
more than 20%. It should be noted that this value is
not a general index for other river channels.
obs
i
Qti1  Qti  (Qtobs
1  Qt )   t 1

 ti1

N (0,0.1Qtobs
1 ) (4)ti1

N (0, 0.005) (5)nti1  nti  ti1

ti1

N (0, 0.01) (6)hti1  htobs  ti1

ti1

N (0, 0.1) (7) nˆt 1 )  ti1

nti1  nti  (nˆtFig. 3 The comparison of averaged water stage in the synthetic
experiment

where Qti , nti and hti indicate the discharge,
Manning roughness n, and downstream water stageof i-th particle at time t, respectively.  ti , ti ,  ti
and  ti indicate the errors drawn from the normal
distribution. Qtobs and htobs indicate the observed
discharge and the observed downstream water stage.
nˆt indicates the average roughness coefficient ineach time step according to the weight in Eq. (11).

Fig. 4 The comparison of inlet discharge in the synthetic
experiment(2) Resampling process
The key idea in the particle filters is to represent
the posterior pdf p( xt | yt ) with a set of randomdraws, called particles (Salmon and Feyen13)). The
posterior density at time t is approximated as Eq. (8)and the updated equation can be shown as Eq. (9).
Fig. 5 The comparison of Manning roughness coefficient in thesynthetic experimentNp( xt | Yt )   wti ( xt  xti )i 1

wti  wti1p( zt | xti ) p( xti | xti1 )q( xti | xti1 , yt )

(8)

wi
t
Wi 
t
N
i
 w
t
j 1

(9)

where xti , w ti denote i th particle and its weight,
respectively.  indicates the Direc delta function
and xt , yt are the state vector and the measurement
vector at time t.
At each updating step, the weight, which is
calculated by Eq. (9), is calculated against theobserved water stage of 2 points (upstream anddownstream ends), respectively. Then, we
calculated joint probability by multiplying the
likelihood of 2 points to represent the weight ofeach particle (Eq. (10)) (Matgen et al.7)). The joint
probability of each particle was normalized asfollows (Eq. (11)):Npi, jwi   w
t
t
j 1

(11)

where wi,j indicates the weight of i th particle atttime t calculated against

j th observation.

Witdenotes the normalized weight of each particle attime t , N p is the number of updating points.
According to the normalized weight, theparticles are removed or multiplied in theresampling step. The systematic resampling method
presented in Kitagawa15) is introduced among
various resampling methods.
(3) Constraints on noise generation
This method considers various cases with the
disturbed channel roughness and inlet discharge, butparticles are controlled by only the normalized

(10)

4

I_28

weight calculated against the updated water stage.
Therefore, it is not enough to control the two
randomly generated variables—the Manningroughness and inlet discharge—with only the weightcalculated against updated water stage. In order toreduce the meaningless calculation, we assume theerrors of observed discharge are fewer than 35%.

4. EVALUATION WITH SYNTHETIC DATA
(1) Synthetic experimental design
We verify the estimation method with thesynthetic data. The synthetic data is based on thereal event, which occurs from 6:00 on October 20,
2004, to 15:00 on October 21, 2004, at the naturalriver reach from Hazukashi station to Nosou station(See Fig. 2). In general, Manning’s roughness isdependent on water depth and ranges from 0.02 to
0.07. Thus, we present a simple relationshipequation of Manning’s roughness to verify the
proposed algorithm. The equation is based on theassumption that Manning’s roughness coefficient( nt ) is linearly inverse proportional to water stage( ht ) at time t as follows:nt  (18  ht ) / 250  0.02

Fig. 6 The comparison of prediction water stage at Hazukashi

Fig. 7 The comparison of prediction discharge at Hazukashi(12)

The water stage of upstream is generated fromthe first step of the synthetic experiment. It isconsidered to be the synthetic truth instead of themeasured upstream water stage. Then predeterminedsynthetic truths, which are water stage, input
discharge, and Manning n, are utilized asverification data and marked as ―TRUE‖ in the Figs.
3, 4, and 5. In addition, only the water stage ofsynthetic truth is utilized as updating data for thesecond step of the synthetic experiment.

Fig. 8 The estimated tendency of Manning roughness
Table 1 Water stage comparison of the estimation and prediction at Hazukashi
RMSE of water stage Maximum absolute
Classification(m)error (m)(2) Evaluation results
Averaged values of water stage, inlet discharge,
and Manning’s roughness coefficient are marked asa black line. The values of each particle are marked
as red points in Figs. 3, 4, and 5. According to thegraphs, the averaged values are similar to ―True,‖
which is synthetic truth in each graph. Root meansquare error (RMSE) is utilized to compare theaveraged values with synthetic truth. The RMSE forwater stage and discharge are 0.11m and 31.75 m3/s,
respectively. Then, the range of particles covers the―True.‖
Particularly, the initial errors of Manning’s
roughness were also corrected during the tracking
procedure by the learning process of particle filters.
Such a simultaneous estimation of the Manning’s n
and inlet discharge could present approximately
averaged values and reduce the errors.

Estimation

0.024

0.075

1hr prediction

0.097

0.258

3hr prediction

0.291

0.758

6hr prediction

0.565

1.376

5. PREDICTION RESULTS AND ANALYSISUSING OBSERVED DATAIn this section, the predictability of the algorithmwith observed data is presented. Most of the applied
conditions for this simulation are the same withsynthetic experiment, but the simulation is different
in that this test utilizes real observed water stage for
updating and true values do not exist.
In the prediction process, Manning’s roughness
coefficient is fixed as the value chosen in theresampling step at the current time step, and the

5

I_29

predicted inlet flows up to 6 hours are determinedby the proportions of particle’s inlet discharge
against the given discharge at the current time step.
The prediction process performs up to 6 hours.
Then the results at 1 hour (1HR), 3 hours (3HR),
and 6 hours (6HR) are compared with the observed
water stage and discharge (OBS) at Hazukashi
station (Figs. 6 and 7). EST indicates the averagevalues with the weight, while the prediction results
are averaged with even weight. Table 1 shows thecomparison results of the estimation and predictionof the water stage.
The estimated water stage and 1 hour aheadwater stage shows good agreement with theobserved one. With increased lead time, theaccuracy of the prediction result is less than the 1hour ahead prediction result due to the fixedroughness coefficient chosen in the resampling step,
while the real channel roughness varies according totime.
In addition, the estimated Manning’s roughness
coefficient is plotted in Fig. 8 with the values ofeach particle. Manning’s roughness coefficient is sovaried due to the flood plain and backwater
increasing the channel roughness when flooding
occurs. Finally, Fig. 7 also shows that the estimateddischarges are similar to the observed discharges,
even though there is some discrepancy at the peak
point and the recession limb of the hydrograph.
However, the discrepancy in Fig. 7 may come fromthe rating curve conversion because the discharge isconverted from the observed water stage.

2) Hsu, M.H., Fu, J.C., and Liu, W.C.: Dynamic routing modelwith real time roughness updating for flood forecasting,
Journal of Hydraulic Engineering, Vol.132, pp.605-619,
2006.
3) Dottori, F., Martina, M. L. V., and Todini, E.: A dynamic
rating curve approach to indirect discharge measurement,
Hydrol. Earth Syst. Sci., Vol.13, pp.847-863, 2009.
4) Di Baldassarre, G. and Montanari, A.: Uncertainty in riverdischarge observations: a quantitative analysis, Hydrol.
Earth Syst. Sci., Vol.13, pp.913-921, 2009.
5) Aricò, C., Nasello, C., and Tucciarelli, T.: Using
unsteady-state water level data to estimate channelroughness and discharge hydrograph. Adv. Wat. Res.,
Vol.32, pp.1223-1240, 2009.
6) Montanari, M., Hostache, R., Matgen, P., Schumann, G.,
Pfister, L., and Hoffmann, L.: Calibration and sequential
updating of a coupled hydrologic model using remotesensing-derived water stage, Hydrol. Earth Syst. Sci.,
Vol.13, pp.367-380, 2009.
7) Matgen, P., Montanari, M., Hostache, R., Pfister, L.,
Hoffmann, L., Plaza, D., Pauwels, V.R.N., De Lannoy,
G.J.M., De Keyser, R., and Savenije, H.H.G.: Towards thesequential assimilation of SAR-derived water stages into
hydraulic models using the particle filter: proof of concept,
Hydrol. Earth Syst. Sci., Vol.14, pp.1773-1785, 2010.
8) Giustarini, L., Matgen, P., Hostache, R., Montanari, M.,
Plaza, D., Pauwels, V.R.N., de Lannoy, G.J.M., de Keyser,
R., Pfister, L., Hoffmann, L., and Savenije, H.H.G.:
Assimilating SAR-derived water level data into a hydraulicmodel: a case study, Hydrol. Earth Syst. Sci., Vol.15,
pp.2349-2365, 2011.
9) Tachikawa, Y., Sudo, J., Shiiba, M., Yorouzu, K., and Kim,
S.: Development of a real-time rivers stage forecasting
method using particle filter, Annual Journal of HydraulicEngineering, JSCE, Vol.55, pp. S_511-S_516, 2011.
10) Bradely, P. C., Nicholas, G., Polson, N.G., and David, S.S.:
A Monte Carlo Approach to Nonnormal and Nonlinear
State-Space Modeling, Journal of the American Statistical
Association, Vol.87, pp.493-500, 1992.
11) Nagata, T.: Hydraulics formulae: Hydraulics worked
examples with CD-ROM, Hydraulic committee, JSCE,
pp.16-19, 2002.
12) Ristic, B., Arulampalam, S., and Gordon, N.: Beyond theKalman Filter, Artech House, pp. 35-65, 2004.
13) Salamon, P. and Feyen L.: Disentangling uncertainties in
distributed hydrological modelling using multiplicative
error models and sequential data assimilation, Water
Resour. Res., Vol.46, 2010.
14) Crissman, R.D., Chiu, C-L., Yu, W., Mizumura, K., andCorbu, I.: Uncertainties in Flow Modeling and Forecastingfor Niagara River, American Society of Civil Engineers
Journal of Hydraulic Engineering, ASCE, Vol.119,
pp.1231–1250, 1994.
15) Kitagawa, G.: Non-Gaussian State-Space Modeling of
Nonstationary Time Series, Journal of the American
Statistical Association, Vol.82, pp.1032-1063, 1996.

6. CONCLUSION
Stochastic approaches are introduced to the 2D
dynamic wave model to incorporate errors toconsider the uncertainties of inlet discharge andchannel roughness. The estimation method isverified with the synthetic experiment, andprediction results were compared with the observed
water stage. The results of tests implemented in thisstudy show that the proposed algorithm is not only
an efficient prediction tool for water stage, but alsoa simultaneous estimation tool of inlet discharge andchannel roughness.
In further research, the lateral flow will beconsidered and the method will be applied to more
events and other reaches for more general approachto natural channels.
REFERENCES
1) Shiiba, M., Laurenson, X., and Tachkawa, Y.: Real-timestage and discharge estimation by a stochastic-dynamic
flood routing model, Hydrological Process, Vol.14,
pp.481-495, 2000.

(Received September 30, 2011)6

I_30