京都大学防災研究所年報 第 48 号 B 平成 17 年 4 月
Annuals of Disas. Prev. Res. Inst., Kyoto Univ., No. 48 B, 2005

Uncertainty Recognition and Quantification of Hydrologic Prediction

Shen CHIANG*, Yasuto TACHIKAWA, Kaoru TAKARA*Graduate School of Engineering, Kyoto University

Synopsis
We propose a methodology to identify prediction uncertainty through recognizing and quantifying
the different uncertainty sources in a hydrologic model. Statistical second moment is used as ameasure of uncertainty; also an index which originated from Nash coefficient of efficiency named
Model Structure Indicating Index (MSII) is proposed to quantify model structure uncertainty. The
results show that MSII can well reflect the goodness of model structure, while a larger value of MSII
indicating a poorer structure of hydrologic model. The index can be used as a tool for implementing
model quantitative comparison (selection).

Keywords: prediction uncertainty, Monte Carlo simulation, hydrological modeling

1. Introduction
Hydrological modeling is the discovery of generallaws and principles that govern the natural
phenomenon under observation.
A hydrologic
model is an integration of mathematical descriptionsof conceptualized hydrologic processes, which serves
for a specific purpose. Consequently, the spatial
scale, temporal scale, structure, architecture, andapplicability of a model are restricted a lot by the
hypothesis of the hydrologic model in most of thecases. As a result, there are numerous hydrologic
models developed for various aspects, and thedevelopment of new hydrologic models or
improvement of previously developed models
continues in Japan and elsewhere. With rapid
advances in computing technology, remote sensing,
GIS and DBMS, the role of hydrologic models isenhanced as a tool in planning, decision making andtends to incorporating with other process models such
like economic, social, political, administrative, andjudicial models. Thus, the watershed hydrologic
models will become a component in a larger
management strategy. Furthermore, these models
will become more global, not only in the sense of
spatial scale but also in the sense of hydrologic details(Singh, 2002). Therefore, a methodology to assess

the error, uncertainty and adequacy of adoptinghydrologic models in specific purpose is needed.
The present study provides a methodology for modelcomparison and selection through model uncertainty
recognition and quantification.
For many years hydrologists have been interested
in the effects of various uncertainties on the accuracy
and reliability of the estimation of catchmenthydrological variables such as peak flow and flood
volume (Høybye and Rosbjerg, 1999).
Amongearly contributions, some of them focus on the rainfall
uncertainty and its influence to the runoff (e.g., Stormet al., 1989); some of them focus on the sensitivity ofmodel structure due to the input error (e.g., Singh and
Woolhiser, 1976). Recent researches relating tohydrologic model uncertainty most refer to parameter
uncertainty identification, the procedure of parameter
calibration, and their impact to simulation result (e.g.,Freer et al., 1989; Kuczera and Parent, 1998).
Among them, Generalized Likelihood Uncertainty
Estimation (GLUE) Methodology (Beven and Binley,1992) offers a path of identifying parameter
uncertainty. Nevertheless, parameter equifinality
became the conclusion of GLUE; uncertainty related
to input data and other factors are excluded. Even it
is said that they could also be included in GLUE but
this has not normally been done (Beven, 2001).

The performance of hydrologic models is
profoundly affected by the sources of uncertainty,
briefly they are:(a) Observed data,(b) Data for model calibration,(c) Parameter space, and(d) Model structure.
Among those, data uncertainty occupies the mostand contaminates other sources of uncertainty.
Underestimating or misunderstanding of these
uncertainty sources and the interrelation among them
may cause tremendous misleading on interpreting theresult of hydrologic models. However, recognizing
and quantifying these different sources of uncertainty
in hydrologic models has received little attention in
the research literature.
In this study, a methodology is proposed torecognize and quantify the different uncertainty
sources. Firstly, Monte Carlo simulation method isapplied to add bias item in model input data series(rainfall), then rainfall realizations, parameter space,and model outcomes (outflow discharge) under
different bias level are acquired. Secondly, by
examining the counter relationship between model
simulation outcomes, calibration outcomes andobserved watershed response series (discharge), an
uncertainty structure is recognized.
Finally,
parameter uncertainty, calibration uncertainty, andmodel structure uncertainty caused by input data
uncertainty are recognized, separated, and quantified
through the methodology.
Statistical second moment is used as a measureof uncertainty, also an index which originated fromNash coefficient named Model Structure Indicating
Index (MSII) is proposed to quantify model structure
uncertainty which can be used as a tool for
implementing model quantitative comparison(selection). For the demostration of the proposed
method, several hydrologic models are applied toperform model comparison. They are: Storage
Function Method (SFM, Kimura, 1961), Topmodel(Kirkby and Beven, 1986) and KW-GIUH (Lee andYen, 1996). Through fixing the value of one
parameter of SFM, a poorer structure model isformed as an extream contradistinction.
The results show that a larger value of MSII
indicating a poorer structure of hydrologic model in a
dynamic manner, that is, incorporating the MSII to
the input uncertainty.
2. Uncertainty
modeling

recognitioninhydrologicalIn this study, prediction uncertainty which camefrom the four kinds of sources mentioned previously isclassified into four categories: system uncertainty,

entire uncertainty, inherent uncertainty, and structure
uncertainty. The definition and the procedure torecognize them are described below.
2.1 System uncertainty
Even it is well known that a hydrologic model is anapproximation of the real phenomena based on thehydrologic cycle, still it is needed to stress that there
always exiting discrepancy between model outcome
and observed data, no matter how precise the model isand how perfect the model calibrated. That is themodel predicting limitation underneath the hypothesis
and architecture of the model, which is defined as thesystem uncertainty in this study.
The system uncertainty can be recognized byevaluating the discrepancy between observed
watershed response series and the model outcome
during the process of model parameter calibration.
This can be clarified by the statement written below:~
y = f (~x ,θ )(1)x denote outcome and input series of awhere ~y and ~hydrologic model f(・) respectively, and θ denotesa parameter set. Denoting X and Y are observed
watershed input and response series data accordingly~
and assuming x = X, then θ can be determined byadjusting or tuning the value of parameters to makethe difference between model outcome ~y andobserved watershed response data Y in an acceptable
range according to some objective function. The
procedure of adjusting the value of parameter set iswell known as “model calibration”. In this study,
least sum of square error (LSE) method is selected as
the objective function:N

S LS (θ m ) = ∑ (Yn − f ( X n ,θ m )) 2

(2)

n =1

where n denotes time step of the simulation timeseries; m denotes the parameter set to be examined.
The discrepancy between observed watershed
response and model outcome during the process ofcalibration is supposed to be the minimum value bycomparison to the possible coming events. This hasoften been observed that the goodness-of-fit between
observed data and estimated data during calibration isbetter than what in validation, not to mention in
implementation. Under this fact, we can say that theuncertainty occurring here denotes the predicting
limitation of the model. Since this is the best that amodel can achieve.
2.2 Entire uncertainty
The discrepancy of the model outcome to theobserved data can be shown as below:

Y = f ( X ,θ ) + ε( 3)

where ε denotes system uncertainty mentioned
above, that is, the predicting limitation of the model
under perfect calibration. In agreement with thisdefinition, it is clear that the uncertainty came fromthe process of calibration. But the real uncertainty
source is supposed to be the data which is used forcalibration.
After calibrating the model parameter, thecalibrated parameter space will reflects its uncertainty
through the model structure and propagates to themodel outcome. This uncertainty can be recognizedby examining the discrepancy between observed
watershed response data and model outcome by using
input data and parameter sets. Also this is the way ofmost of the researches dealing with parameter
sensitivity. This uncertainty is defined as the entire
uncertainty in this study, since actually this is theutmost uncertainty a model could has under existing
input uncertainty.
2.3 Inherent uncertainty
Another categorized uncertainty what we called“inherent uncertainty” is very difficult to berecognized since it always contaminated by theuncertainty of model structure and input data.
Inherent uncertainty is different to what so calledparameter uncertainty since the propagation effect of

Fig. 1 Schematic diagram of the uncertainty structure

Fig. 2 Schematic diagram of the uncertainty
acquisition processmodel structure is excluded here intentionally.
Inherent uncertainty represents the sensitivity ofparameter space which determined according to the
input uncertainty and reflects to model outcomes.
This can be examined by the discrepancy amongmodel outcomes derived from different best fit
parameter sets. Watershed response data is not used
here, which indicates the model structure uncertainty
is eliminated as much as possible.
2.4 Structure uncertainty
The distance between entire uncertainty andinherently uncertainty is structure uncertainty. Since
the parameter set used for extracting both of theuncertainty is the same. The different is the data
used for comparison: observed data and outcome
during calibration. The propagation effect plays animportant role here, which is why the difference
between entire uncertainty and inherent uncertainty iscalled structure uncertainty in this study.
The schematic diagram of uncertainty structure
proposed here is depicted in Fig.1, which reveals astatic view of model uncertainty structure, that is, theinput uncertainty (data uncertainty) is neglected or
fixed to a certain magnitude. In order to know the
behavior of model structure uncertainty under
different input uncertainty level, a dynamic view ofmodel uncertainty structure is introduced byincreasing the input uncertainty during the processesof uncertainty recognition. The behavior of model
structure uncertainty under different input uncertainty
can be examined by observing the discrepancy among
different model outcomes accordingly. Fig.2 shows
the schematic diagram of the uncertainty acquisition
process. By increasing the magnitude of the bias
item (input uncertainty) which located in the centerof the diagram, the subsequent change of different
categorized uncertainty can be observed andevaluated. It is expected that within certain level of
input uncertainty, model with a better structure still
capable to reflect the true catchment hydrologic
behavior in certain level.
3. Uncertainty quantification
Several indexes are selected to perform
uncertainty quantification in this study; the basic formof the index is statistical second moment and Nashcoefficient.
The categorized prediction uncertainties definedin previous section are mathematically described as
below:

3.1 Index for quantifying system uncertainty∑ (Q − Q )Nash = 1 −∑ (Q − Q )
o

2

b

s

o

(4)

2

o

where Qo and Qb indicate the observed watershed
response series and model outcome by using the best
fit parameter set, n is the time step of the time series.
SDs quantifies the error variance, the residuals,
between observed watershed response and the best
fitted simulated outcomes; while Nashs elucidates theperformance of the best fitted simulated outcomes.
SD always larger than 0, while the range of Nashcoefficient is: − ∞ < Nash ≤ 1 .
3.2 Index for quantifying entire uncertainty

Nashe = 1 −∑ (Q − Q )∑ (Q − Q )
o

2

e

o

2

(5)

o

where Qe is the model outcomes acquired by using aparameter set within the whole parameter space andentire rainfall realizations.
Consequently, SDe
explains the total prediction error.
3.3 Index for quantifying inherent uncertainty

Nashi = 1 −is used as the numerator in the equation, while system
uncertainty is used as a denominator to form acriterion to evaluate model structure. The index
expresses the total variance unexplained by thestandard deviation in a dimensionless form. The
range of MSII is: 0 ≤ MSII < ∞ .

∑ (Q − Q )∑ (Q − Q )
b

b

2

e

2

(6)

b

Instead of sampling the parameter space directly likewhat GLUE did, the study here generates theparameter set space by introducing noise item into
input data with specified probability distribution.
Here Normal distribution with mean equals to zeroand standard deviation from 1.0 to 9.0 (mm/hr) is used
to acquire model parameter space and outcomes under
different input uncertainty. For each iteration, 10000
model outcomes for each specified input uncertainty
were derived from the combination of rainfall series
and parameter set generate output series through themodel. The system uncertainty and the prediction
ability were identified and recognized bycorresponding parameter set. Detail algorithm for
uncertainty recognition and quantification is listed in
Fig.3.
5. Model description and process of observed data
To demonstrate uncertainty identification, Storage
Function Model (SFM) proposed by Kimura (1975),
Topmodel (Kirkby and Beven) and KW-GIUH are
applied at the Yasu River basin (355km2) in this study.
Below are the brief description of the models.
5.1 Storage Function Model:
The form of SFM is as:

SDi and Nashi explain the inherent uncertainty because
Qe is calculated with a possible parameter set
identified by other model outcomes, thus, thediscrepancy between Qb and Qe indicates the error
inevitable in a model formulation.
3.4 Index for structure uncertainty
In order to implement dynamic view of therelationship among system uncertainty, structure
uncertainty and parameter uncertainty caused by input
data uncertainty, the Nash coefficient is used toformulate a Model Structure Indicating Index (MSII)
and is defined as:

Nashi − Nashe
MSII =Nashs

4. Algorithm for uncertainty recognition &quantificationdS= re (t − Tl ) − q, S = kq p
dt
re =f × r,r,if
if∑r ≤ R∑r > R

SA

(8)

SA

where S = water storage height; r = rainfall intensity;
q= runoff height; t = time step; and Tl = the lag-time.
This model is often used for the flood runoff
calculation in a basin with an area of less than five
hundred square kilometers in Japan (Takara and
For input uncertainty σ x = 1.0, 2.0,…, 8.0, 9.0 (unit:(7)

The difference between entire and inherent
uncertainty denotes the structure uncertainty. This

mm/hr)
1. Use Monte Carlo simulation to sample 100
rainfall realizations according to a real recorded
event by adding noise item to it.
2. Use LSE to determine 100 parameter setsregarding to its corresponding rainfall
realization.
N

S LS (θ j ) = min(∑ (Yn − f n ( x + ξ j , θ j )) 2 ), j = 1,...,100
n =1

3. System uncertainty ε s is acquired by examining
the discrepancy between 100 best fit modeloutcomes during calibration and observed
watershed response. Mean value of SDs andNashs for 100 cases are used for uncertainty
quantification.
Y = f ( x + ξ j , θ j ) + ε s , j = 1,...,100
4. Entire uncertainty ε e is acquired by examining
the discrepancy between observed watershed
response and 10000 model outcomes. Meanvalue of SDe and Nashe for 10000 cases are used
for uncertainty quantification.
Y = f ( x + ξ j ,θ j ) + ε e ;j = 1,...,100; k = 1,...,100

5. Inherent uncertainty ε i is acquired by examining
the discrepancy between 100 best fit modeloutcome during calibration and the rest 9900
model outcome.
Mean value of SDi andNashi for 10000 cases are used for uncertainty
quantification.
y l = f ( x + ξ l ,θ k ) + ε i , l ≠ k ;l = 1,...100; k = 1,...,100
yl = f l ( x + ξ l ,θ l ), l = 1,...,100
6. Apply Eq. (10) to acquire MSII of the Model.
End
where ξ 〜N (0, σ x 2 ) , n denotes time step of thesimulation time series, x denotes the observed
input data series.
Fig. 3 Algorithm for uncertainty recognition andquantification
Takasao, 1985).
Parameter p is a constant,commonly the value is 0.6, f is the ratio of contribution
area of the watershed which generates outflow, RSA isaccumulating saturated rainfall. Parameter k issolved by the equation proposed by Nagai et al.
(1982).
A fully functional SFM with parameter Tl, f, and
RSA is used.
For contradistinction, apoorer-structured model with comparison to theoriginal SFM by fixing the value of RSA to 0.0 ismanipulated in this study.
5.2 Topmodel
Topmodel is a set of programs for rainfall-runoff
modeling in single or multiple subcatchment in asemi-distributed way and using gridded elevation
data for the catchment area. It is considered aphysically based model as its parameters can be,theoretically measured in situ (Beven and Kirkby,

1979, Beven et al., 1984).
In this study,
TOPMODEL 95.02, written in FORTRAN 77 isused.
Topmodel uses the distribution of the topographic
index λof hydrological similarity:

λi = ln(ai)
tan β i(9)where ai (unit: m) is the specific contribution area,
which denotes the area draining though a grid cell perunit length of contour. tan β i denotes the local
surface slope. Flow is separated into surface runoff
generated by rainfall on saturated contributing areaand subsurface downhill flow. There are four basic
assumptions to relate local downslope flow from apoint to discharge at the catchment outlet (Campling
et al., 2002):
A1. The dynamics of the saturated zone areapproximatedbysuccessive
steadystate
representation.
A2. The recharge rate r (m/h) entering the watertable is spatially homogeneous.
A3. The effective hydraulic gradient of thesaturated zone is approximated by the local surface
topographic gradient tan β i .
A4. The distribution of downslope transmissivity
To (m2/h) with depth is a function of storage deficit.
From A1 and A2, the recharge rate of the flux
drains into the water table has the form like:qi = rai = To e − Si / m tan β i(10)where qi is the flux drains into the water table for
topographic index class i. r denotes recharge rate
the water drains to the water table, and Si is a local
storage deficit, and m is a model parameter
controlling the rate of decline of transmissivity withincreasing storage deficit.
By arranging Eq.(10), the local storage deficit
Si can be calculated:

S i = − m ln(rai)
To tan β i(11)The mean catchment storage deficit S can beobtained by integrating Eq.(11) over the entire
catchment area A:

1
rai
)
Ai (−m ln∑
To tan β i
A i

S=

(12)

where Ai is the fractional area of the topographic
index class i.
Then storage deficit of each
topographic index i can be expressed like below:ai)
S i = S + m(λ − lntan β i
where

(13)

λ is the areal average of the topographic

index lnλ =aof whole catchment.
tan β

a
1 A
ln( i ) Ai dA∫
0Atan β i

(14)

The vertical drainage qv from unsaturated zone
storage at any point I is calculated by Eq.(15) below:

qv =

Suz
Sitd

(15)

where Suz is the storage in the unsaturated zone, and td
is a time delay constant. The total vertical drainage
then can be calculated by summing up the flux alongcatchment area.

Qv = ∑ q v Ai

(16)

Output discharge from saturated zone cancalculatedusingasubsurface
storage
deficit-discharge function of the form:

Qb = Q0 e − S / m

(17)

−λ

is the discharge when S is zero.

where Q0 = Ae

At each topographic index class i, unsaturated
and saturated zone fluxes are simulated. Initial baseflow Q0 and the initial root zone storage deficit Sr0
are specified at the start of simulation. For thefollowing time step, the average storage deficit iscalculated according to Eq. (18). By which thecatchment average storage deficit S is updated bysubtracting the unsaturated zone recharge and adding
the baseflow from the previous time step:S t = S t −1 + [Qbt −1 − Qvt −1 ](18)

Topographic index derivation was obtained byusing DEM algorithm.
Infiltration excess
mechanism is not included in the study.
Subcatchment discharge is routed to the catchment
outlet by using time-area curve which generated fromdem-based algorithm with constant velocity all overthe catchment area.
5.3 KW-GIUH
If one unit rainfall excess fall on the watershedinstantaneously at time t=0, by assuming that theraindrops are independent and isolated from eachother and ignoring the raindrops fall on the river, thenthe distribution of the number of raindrop appears at
the outlet to time is the instantaneous unit hydrographof the watershed.
For a watershed of Ω order which categorizedunder the Horton-Strahler ordering scheme, the pathof a raindrop drains to the outlet of the watershed will
follow the sequence from a lower order overland flow
phase to a stream flow phase then a higher orderstream finally outlet of the watershed. If w denotes a
specified path xoi→xi→xj→・ ・ ・ →xΩ, the probabilityof a drop of rainfall excess adopting this path can be
expressed as

P ( w) = POAi ⋅ Pxoi xi ⋅ ⋅ ⋅ Pxi x j ⋅ ⋅ ⋅ Pxk xΩ(19)where POAi is the ratio of the ith-order overland flow
area to the total watershed area.

Pxoi xi is theprobability that rainwater drains from ith-order
overland flow to ith-order stream, the value supposed
to be 1 according to the assumption mentioned above;
Pxi x j is the transitional probability of a raindrop froman ith-order stream to a jth-order stream, which isequal to the ratio of the number of ith-order stream
drains to jth-order stream divided by the number of
ith-order stream.
Since the raindrops are assumed to behave
independently and isolated to each other, the storage
inside the watershed at time t equals to the number ofraindrops that the travel time T to the outlet larger thant.
Let f x j (t ) represents the distribution of travel
time that a raindrop moving in x j phase, and theaverage of f x j (t ) is Tx j , then the storage can beshown as (Feller, 1968):t

P (T ≤ t ) = ∑ [ ∫ f xoi (t ' ) * f xi (t ' ) * f x j (t ' ) * ⋅ ⋅ ⋅ * f xΩ (t ' )dt ' ] ⋅ P ( w)
w∈W

0

(20)

inflow rate joining the main stream flow qoi .

α o can be expressed as:
α o =

S oi
no

(24)

where S oi is the mean slope of the ith-order overland

Fig. 4 V-shape subbasin
among it, * denotes convolution integral. Since theraindrops are assumed to behave independently andisolated to each other, the storage inside the watershedat time t equals to the number of raindrops that thetravel time T to the outlet larger than t. For one unitof rainfall excess, I (t ) = 1.0 when t=0; andI (t ) = 0 when t ≠ 0 , then the outflow hydrographu(t) at the outlet of the watershed is the instantaneous
unit hydrograph of the watershed. From Eq. (20) we
can get (Rodriguez-Iturbe and Valdes, 1979; Gupta etal., 1980).

flow, and no is an effective roughness coefficient forthe overland flow plane, which is a parameter of themodel need to be calibrated. The flow rate at the endof a plane increases with time until it reaches
equilibrium. The longest time for a raindrop to travel
through the ith-order overland plane Txoi , is the time
for the flow to reach equilibrium in the plane. Thus,
the discharge for the ith-order overland subbasin at
any time t < Txoi is:

qoi = qL Loi =u (t ) = ∑ [ f xoi (t ) * f xi (t ) * f x j (t ) * ⋅ ⋅ ⋅ * f xΩ (t )dt ] ⋅ P ( w)

Soi m
hoi
no

(25)

w∈W

(21)
The following procedure is the determination of the
distribution function for the travel time.
KW-GIUH (Lee and Yen, 1997) is a refinement of
the GIUH method by using kinematic wave
approximation to calculate the travel time of overland
flow inside the catchment.
In KW-GIUH, an ith-order subbasins of thewatershed conceptually simplified as consisting of
two identical rectangular overland-flow planes. Eachplane contributes a lateral discharge into a channel ofconstant cross section and slope as shown as in Fig.4.
Assuming rainfall excess with intensity q L fall on
the plane homogeneously, the continuity equationand simplified momentum equation can be expressedas:then travel time Txoi can be calculated by(Wooding,1997):

Txoi =hoinL= ( 1 /o2 oim−1 )1 / m
qL
S oi qL(26)

Since the overland flow will concentrate to thechannel on the center of the V-shape plain eventually,
for a channel with channel width Bi , the lateraldischarge per unit width when it is equilibrium status
is 2q L Loi / Bi , then the travel time of ith-orderchannel flow can be expressed as(Lee and Yen,
1997):
1⎤⎡
2q L nc Loi Lci ⎞ m
Bi ⎢⎛⎜ m⎥⎟
T xi =
hh− coi ⎥co +2q L Loi ⎢⎜⎝ j
S c1i / 2 Bi ⎟⎠
⎥
⎢
⎦
⎣

(27)

∂hoi ∂qoi
+
= qL
∂t
∂x

(22)

qoi = α h(23)

the channel flow, which is a parameter need to becalibrated, hcoi is the water depth that drains from

where x is the distance along flow direction; hoi is the

upstream channel to the ith-order channel. For an Ωorder watershed, when i=1, hcoi =0; when 1< I <Ω,mo oi

where S ci is the mean slope of the ith-order channel
flow, and nc is an effective roughness coefficient for

ith-order overland flow depth; qoi is ith-order
overland flow discharge per unit width; m is constant,
α o is the parameter which reflects the hydrauliccharacteristics of the overland flow, and qL is lateralhcoi can be expressed as(Lee and Yen, 1997):

1⎡ q L nc ( N i Ai − APOAi ) ⎤ mhcoi = ⎢⎥N i Bi S c1i / 2
⎦⎥
⎣⎢

(28)

among it, Ni is the number of the ith-order channels;

Ai denotes the average area of ith-order subbasin, thearea is not ith-order subbasin but also its upstream
area is included. A = total area of the watershed;
POAi is the ratio of the ith-order overland flow area tothe total watershed area..
The mean length of the ith-order V-shape overland
flow planes is:

APOAi
2 N i Lci

Loi =

(29)

Lci is the mean channel length of the ith-ordersubbasins. The channel width is derived by linear
decreasing from the watershed outlet(Lee and Yen,
1997):

6. Results
Fig.5 through Fig.8 are the plot of TOPMODEL,
SFM, parameter-constrained SFM and KW-GIUH.
It can be seen that the entire uncertainty becoming
larger while input uncertainty increasing. Also themagnitude of standard deviation of the categorized
behaves the same. It is always expected that entire
uncertainty increases while input uncertainty
increasing; where the discrepancy between entire andinherent uncertainty is in the same tendency. System
uncertainty is supposed located in the middle between
entire and inherent uncertainty. The range between
entire and inherent uncertainty indicates the goodnessof a model structure; that is, the ability of a model
simulates the behavior of the watershed throughobserved data. The broader the distance is the worse
structure the model has.
The distance between entire uncertainty andsystem uncertainty indicates the divergence between
calibration outcome and simulation outcome. Fig.9
SFM

2

i

1

l =1

(30)

Ω

∑L

ci

l =1

Nash coefficientBi =BΩ ∑ Lci0
Entire Uncertainty-1

Inherent Uncerainty-2
System Uncertainty
-3
-4
-5
1

2

TOPMODEL

1

4

5

6

7

8

9

Fig. 6 Uncertainty quantification result of SFM

0 .8

Nash coefficien t

3

input uncertainty (mm/hr)

Entire Unc ertainty

0 .6

Inh ere nt Unc erainty
0 .4

System Unc ertainty

parameter-constrained SFM

0 .2
0
1

2

3

4

5

6

7

8

1.2

9

1

in put u nc ertainty (mm/ hr)

0.8

quantification

resultofBΩ is the channel width on the outlet of thewatershed, which can be acquired by field survey.
The rest of all the geomorphic factors can be acquired
manually from topographic map or derived
automatically through raster elevation data and DEMs
by some algorithms.
Rainfall data was collected from four rainfall
gauging stations inside the Yasu River basin(355km2); they are Yasu, Minakuchi, Kouka andOogawara. The average precipitation was obtainedby using Thiessen polygon method.

0.6

Nash coefficientFig. 5 Uncertainty
TOPMODEL

0.4
Entire Uncertainty

0.2

Inherent Uncerainty0

System Uncertainty
-0.2
-0.4
-0.6
-0.8
-1
1

2

3

4

5

6

7

8

9

input uncertainty (mm/hr)

Fig. 7 Uncertainty quantification
parameter-constrained SFM

resultofshows the MSII of each model.
Exceptparameter-constrained SFM, there is no apparent
distinction between the four candidate models withinsmall input uncertainty. However, along with theincreasing of the input uncertainty, KW-GIUH is
structurally more stable than the TOPMODEL and
SFM in this study.

7. Conclusions
For specific catchment, a methodology for uncertainty
recognition and quantification is proposed which canbe used as a tool for hydrologic models quantitative
KW-GIUH1
0.9

Nash coefficient

0.8
0.7
0.6

Entire Uncertainty

0.5

Inherent Uncerainty

0.4

System Uncertainty

0.3

structure. While a good model structure is defined
as: Within certain level of input uncertainty, themodel still capable to reflect the true catchment
behavior. The results show that a larger value of
MSII indicates a poorer structure of hydrologic model,
within increasing input uncertainty the tendency
becomes more apparently. The goodness of theprocess of model calibration is also included in theindex explicitly.
The index can be used as a tool for
implementing model quantitative comparison(selection).
Acknowledgments
A part of this study was supported by
Grant-in-Aid for
Scientific Research (A)(1)
16206050 (PI: Prof. K. Takeuchi, Yamanashi Univ.)
and CREST program of JST (PI: Prof. K. Takara,
Kyoto Univ.).

0.2
0.1

References

0
1

2

3

4

5

6

7

8

9

input uncertainty (mm/hr)

Fig. 8 Uncertainty quantification result of KW-GIUH

MSII

9
8

MSII

7
6

TOPMODEL

5

SFM

4

SFM constrained

3

KW-GIUH

2
1
0
1

2

3

4

5

6

7

8

9

Input uncertainty (mm/hr)

Fig. 9 MSII of TOPMODEL,
parameter-constrained SFM and KW-GIUH

SFM,

comparison or evaluation. Instead of sampling theparameter space directly like what GLUE did, thestudy here generates the parameter set space byintroducing noise item into input data with specified
probability distribution. This reflects the truth thatparameter uncertainty came from uncertainty of datato hand and the way the model structure responses it.
The results show that within increasing input
uncertainty, the distance between entire and inherent
uncertainty is also increased. A smaller magnitudeof the ratio of inherent uncertainty to entire
uncertainty indicates the system uncertainty (i.e. theprediction limitation) is larger, which means lackingof the ability to simulate the true watershed behavior.
MSII is proposed to evaluate the goodness of model

Beven, K. (2001): How far can we go in distributedhydrological modeling? Hydrology and Earth
System Science, Vol. 5(1), pp. 1-12.
Beven, K. J. and Binley, A. M. (1992): The future of
distributed models: Model calibration anduncertainty prediction, Hydrol. Proces., Vol. 6,
pp. 279-298.
Campling, P., Gobin, A., Beven, K. and Feyen, J.
(2002): Rainfall-runoff modeling of a humid
tropical catchment: the TOPMODEL approach,Hydrological Process, Vol. 16, pp231-253.
Freer, J., Beven, K. and Ambroise, B. (1996):
Bayesian estimation of uncertainty in runoff
prediction and the value of data: An application of
the GLUE approach, Water Resources Research,
Vol. 32(7), pp. 2161-2173.
Gupta, V. K., Waymire, E., and Wang, C. T. (1980): A
representation of an instantaneous unit hydrographfrom geomorphology, Water Resour. Res., Vol.
16(5), pp.855-862.
Høybye, J. and Rosbjerg, D. (1999): Effect of input
and parameter uncertainties in rainfall-runoff
simulation, Journal of Hydrologic Engineering,Vol. 4(3), pp.214-224.
Kimura, T. (1975): Storage Function Method.
Kawanabe Books Co., Ltd., Tokyo (in Japanese).
Kuczera, G. and Parent, E. (1998): Monte Carloassessment of parameter uncertainty in conceptual
catchment models: the Metropolis algorithm,
Journal of Hydrology, Vol. 211, pp. 69-85.
Lee, K. T. and Ben, C. Y. (1997): Geomorphology
and kinematic-wave-based hydrograph derivation,
Journal of Hydraulic Engineering, Vol. 123(1),pp73-80.
Nagai, A., Kadoya M., Sugiyama H. and Suzuki K.
(1982): Synthesizing storage function model forflood runoff analysis. Disaster Prevention

Research Institute Annuals, No. 25 B-2, pp.
207-220.
Rodriguez-Iturbe, I. and Valdes, J. B. (1979): The
geomorphologic structure of hydrologic response,
Water Resour. Res., Vol. 15(6), pp. 1409-1420.
Sherman, L. K. (1932): Streamflow from rainfall by
the unit-graph method, Eng. New-Rec., 108,
pp.501-505.
Singh, V. P. and Woolhiser, D. A. (1976): Sensitivityof linear and nonlinear surface runoff models to
input errors, Journal of Hydrology, Vol. 29, pp.
243-249.

Singh, V. P. and Woolhiser, D. A. (2002):
Mathematical modeling of watershed hydrology:
Journal of Hydrologic Engineering, Vol. 7, pp.
270-292.
Storm, B., Høgh Jensen, K. and Refsgaard, J. C.
(1989): Estimation of catchment rainfall
uncertainty and its influence on runoff prediction,
Nordic Hydro., Vol. 19, pp. 77-88.
Takara, K. and Takasao, T. (1985): A new attempt toevaluate rainfall-runoff models from theviewpoint of stochastic transformation, Journal of
Hydroscience and Hydraulic Engineering, Vol.
3(2), pp. 61-72.

水文モデル予測の不確かさの認出と定量化
江申*, 立川

康人, 宝

馨

*

京都大学大学院工学研究科土木工学専攻

要旨
本研究では，水文予測の不確かさを生み出すもととなる異なる要因を認識し定量化することにより，水文予測の不確かさを同定する手法を提案する。予測の不確かさを計量する手法として統計的な二次モーメ
ントを利用する。またモデル構造の不確かさを定量化する指標として，Nash 指標を基にしたモデル評価指標 MSII (Model Structure Indicating Index) を提案する。MSII の値は，モデル構造が不適当なモデルほど大きな値を示す。この指標は水文モデルを予測の不確かさの観点から定量的に評価するために有効に利用することができる。
キーワード: 予測の不確かさ，モンテカルロシミュレーション，水文モデリング，PUB