Embed Size (px)
Citation preview
This article was downloaded by: [University of New Hampshire]On: 09 October 2014, At: 00:01Publisher: Taylor & FrancisInforma Ltd Registered in England and Wales Registered Number: 1072954 Registeredoffice: Mortimer House, 37-41 Mortimer Street, London W1T 3JH, UK
Hydrological Sciences BulletinPublication details, including instructions for authors andsubscription information:http://www.tandfonline.com/loi/thsj19
Adaptive hydrological forecasting—areview / Revue des méthodes deprévision hydrologique ajustablesP. E. O'CONNELL a & R. T. CLARKE aa Institute of Hydrology , Wallingford, Oxfordshire, OX10 8BB,UKPublished online: 21 Dec 2009.
To cite this article: P. E. O'CONNELL & R. T. CLARKE (1981) Adaptive hydrological forecasting—areview / Revue des méthodes de prévision hydrologique ajustables, Hydrological SciencesBulletin, 26:2, 179-205, DOI: 10.1080/02626668109490875
To link to this article: http://dx.doi.org/10.1080/02626668109490875
PLEASE SCROLL DOWN FOR ARTICLE
Taylor & Francis makes every effort to ensure the accuracy of all the information (the“Content”) contained in the publications on our platform. However, Taylor & Francis,our agents, and our licensors make no representations or warranties whatsoeveras to the accuracy, completeness, or suitability for any purpose of the Content. Anyopinions and views expressed in this publication are the opinions and views of theauthors, and are not the views of or endorsed by Taylor & Francis. The accuracy ofthe Content should not be relied upon and should be independently verified withprimary sources of information. Taylor and Francis shall not be liable for any losses,actions, claims, proceedings, demands, costs, expenses, damages, and other liabilitieswhatsoever or howsoever caused arising directly or indirectly in connection with, inrelation to or arising out of the use of the Content.
This article may be used for research, teaching, and private study purposes. Anysubstantial or systematic reproduction, redistribution, reselling, loan, sub-licensing,systematic supply, or distribution in any form to anyone is expressly forbidden. Terms& Conditions of access and use can be found at http://www.tandfonline.com/page/terms-and-conditions
Hydrological Sciences-Builetin-des Sciences Hydrologiques, 26, 2, 6/1981
Adaptive hydrological forecasting — a review
P, E, O'CONNELL & R, T. CLARKE Institute of Hydrology, Wallingford, Oxfordshire OX10 8BB, UK
ABSTRACT With the increasing use of telemetry in the control of water resource systems, a considerable amount of effort is being devoted to the development of models and parameter estimation techniques for on-line use. A variety of models and parameter estimation algorithms have been considered, ranging from complex conceptual models of the soil moisture accounting type, which are traditionally calibrated off-line, to state-space/Kalman filter models which, perhaps, have enjoyed undue popularity in the recent literature due to their mathematical elegance. The fundamental assumptions underlying the various approaches are reviewed, and the validity of these assumptions in the hydrological forecasting context is assessed. The paper draws on some results obtained during a recent workshop at the Institute of Hydrology in making assessments of the relative merits of different models and parameter estimation algorithms; these results have been derived from an intercomparison of a number of real time forecasting models.
Revue des méthodes de prévision hydrologique ajustables RESUME Avec le développement de 1'emploi de la télémesure dans le contrôle des systèmes de ressources en eau, des efforts considérables sont consacrés à la mise au point de modèles et de techniques d'estimation des paramètres pour leur emploi avec possibilité d'ajustement en cours de prévision. On a considéré divers modèles et divers algorithmes d'estimation des paramètres depuis les modèles conceptuels complexes du type qui tient compte de l'humidité du sol, lesquels sont habituellement calibrés "off-line" jusqu'aux modèles à filtre de Kalman/état-espace qui ont bénéficiés peut être d'une popularité imméritée dans la literature récente par suite de leur élégance du point de vue mathématique. Les hypothèses fondamentales qui sont sous-jacentes aux diverses approches sont passées en revue et on estime la validité de ces hypothèses dans le contexte de la prévision hydrologique. Cet article part de certains résultats obtenus lors de deux "ateliers" récents qui se sont tenus à l'Institut d'Hydrologie en vue de procède à une estimation des mérites relatifs de différents modèles et algorithmes d'estimation des paramètres; ces résultats ont été déduits d'études fondamentales d'intercomparaison d'un certain nombre de modèles de prévision en temps réel.
179
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
180 P.E. O'Connell & R.T. Clarke
DEFINITIONS
For the purpose of this paper, adaptive hydrological forecasting is defined as follows. We suppose (a) that it is necessary to provide forecasts of river discharge in the near future (say for the next 6 h) for a basin instrumented with telemetering rainfall recorders and a telemetering water level gauge, for which a reliable stage discharge curve is available; (b) that water level is recorded at intervals At, commonly but not necessarily equal, and that the rainfall accumulated by each recorder is available for the same time intervals; (c) that a model q-£ = f ({pt},(5_) + et has been identified which describes the transformation by which mean areal rainfall p-(-in the intervals (t-At,t) becomes runoff qt; (d) that, at time t, estimates 6_ of the model parameters have been calculated from past observations of pt and q.(- (since these estimates will depend upon the records up to time t, we include the suffix t, writing 6 = 6j-) ; (e) that the model is to be used to make forecasts qt+l, §t+2' •••» of runoff during future time intervals (t,t+At), (t+At,t+2At), ... .
If (a) the estimates 9 are corrected at the end of the time interval (t,t+At), when pt+1, <ït+l become available for use to give a new estimate 8-t+i = ©t + correction; (b) the forecasts of future runoff qt+2 ' t+3 ' •••' a r e a^-so corrected by the use of pt+1, qt+l' then we have adaptive estimation of the parameters Q_, and adaptive forecasting of future runoff. The definition can clearly be extended to the case where water quality variables as well as runoff occur on the left-hand side of the model equation, giving
it = i< {Pt }' ••• i } + £t
where y_t, f_(-) and e_t are now vector variables. The distinction between adaptive and non-adaptive estimation and
forecasting can be further illustrated by consideration of the standard off-line procedure for model fitting and use. If records of rainfall and runoff are available up to time Tj, these may be used as before to estimate the parameters ^(i); holding these estimates unchanged, non-adaptive forecasting proceeds to use the fitted model q- = f({p^}, ... 8 (i)), together with measured rainfall pt, to provide forecasts of runoff qt+l/ <3t+2' •••' until some future time T]_+T2, when the entire records, now of length T1+T2, are used to recalculate estimates of model parameters, giving estimates 6_(2)- These are again held constant until time T1+T2+T3 when the re-estimation, using all the record, is repeated. Thus, for non-adaptive forecasting, recalculation of model parameters and forecasts occurs at irregular intervals, or possibly not at all, using the entire records of rainfall and runoff in the recalculation; for adaptive forecasting, recalculation of model parameters and forecasts occurs whenever new observations become available, and the new estimates of parameter values are obtained by using the most recent data to correct the old estimates.
The above definitions have been given within the context of rainfall-runoff models; where runoff at a downstream site is to be forecast given observations of upstream runoff, extension to this case is straightforward, the general form for the model equation now being written qD t = f ({qjj ^ } , 6) + £t.
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 181
HISTORICAL PERSPECTIVE
It is clear that adaptive forecasts have been made for as long as discharges have been plotted on graph paper and the curves drawn through them extrapolated; such curves have been modified by eye as new measurements of water level (and hence of estimated discharge) were passed to the plotter. Whilst adaptive forecasting is not new, the use of numerical procedures (algorithms) for the modification of forecasts, given new observations, is a relatively recent development. In the UK work on adaptive forecasting derived considerable impetus from the Dee research programme, begun in 1967 to investigate how modern methods of management, mathematical modelling and remote data collection could be combined for the more efficient operation of a water resource system. In 1972, a telemetric network had been put into operation; by 1975, rainfall-runoff models for catchment response and hydraulic models for channel routing had been selected, calibrated and tested. This was followed by much research to render the models adaptive (Lambert, 1978; Moore & Weiss, 1980; Jones & Moore, 1980; O'Connell, 1980; Cameron, 1980). Research of a basic nature also flourished at the International Institute of Applied Systems Analysis and was later disseminated by the group working there (Szollosi-Nagy, 1976a, b; Szollosi-Nagy et al., 1977; Todini & Wallis, 1978), whilst work on adaptive forecasting was also active in the USA (e.g. Kitanidis & Bras, 1978).
Tf€ USE OF CONCEPTUAL MODELS FOR ADAPTIVE FORECASTING
"Conceptual" rainfall-runoff models (those in which river basin behaviour is represented as a network of storages, with rainfall "input" partitioned into components which are routed through the storages either to the basin outfall as runoff, or to the atmosphere as evaporation) have received relatively little attention in adaptive forecasting. O'Connell (1980) has set out the possible reasons for this; by comparison with the simple structure of many input-output or black box models, conceptual models do not readily lend themselves to the updating of forecasts in a computationally simple manner. Four possibilities for adaptive forecasting exist with conceptual models:
(a) the estimates 9t may be adjusted so that the required agreement between forecast qt and observed qt is obtained. However, for models with an appreciable number of parameters, manual adjustment of their parameters would not be straightforward, whilst automatic parameter optimization, as usually applied, would require the calculation of an objective function derived from all data up to the present;
(b) the contents of the model reservoirs may be adjusted in each time interval to give the desired agreement between forecast qt and observed qt. However, this would probably need to be done empirically, and would require specialist knowledge of the model which may not always be readily to hand;
(c) the conceptual model may be cast into state-space form, so that adaptive estimation algorithms of the kind discussed below may be applied. Some success has been obtained with this approach
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
182 P.E. O'Connell & R.T. Clarke
(Kitanidis & Bras, 1978), but the nonlinearities in the system equation can give rise to formidable computational difficulties;
(d) the conceptual model may be fitted off-line, and an auto-regressive moving average (ARMA) model may be fitted to the residuals qt-qt. Using the methods of Box & Jenkins (1970), the "first-stage" forecast qt+l provided by the conceptual model may be adjusted by using the ARMA model to provide a forecast of the residual qt+i-<ït+l- This procedure has been explored by Jamieson et al. (1972); it should be noted, however, that the forecasts so obtained are not strictly adaptive in the sense defined above, unless the parameters of the ARMA model are adaptively estimated after receipt of each new residual.
Recently, a new approach to the use of conceptual models for adaptive forecasting has been proposed (Tucci & Clarke, 1980), in which only the M most recent values of observed rainfall and runoff (pt,qt) are retained in computer store together with the contents of the model storages at times t-M; at time t, an objective (least squares) function is calculated, and minimized with respect to the parameters 9, using as initial values the optimized values 6t_]_ calculated in the preceding time interval. The updated parameter values 6^ are used to calculate the forecasts qt+l, *3t+2' •-•' contents of model storages at time t-M are overwritten in the computer by those at time t-M+1, and the values Pt-M' ^t-M a r e
overwritten by Pt-M+1' ^t-M+l• This procedure shows some promise in that the lengths of rainfall and runoff records retained in computer storage are not excessive, and the automatic optimization procedure (although undertaken in each time interval) begins on each occasion with good starting values. However, further research is required to establish what difficulties are likely to arise in practice, and to explore how short a length of record M may be used.
ISO (INPUT-STORAGE-OUTPUT) MODELS
Models of this type were developed originally by Lambert (1969, 1972) and extended by McKerchar (1975), Lowing et al. (1975) and Green (1979). They start from the continuity equation
dS/dt = p - e - q (1)
in which S = S(t) is the mean depth of water stored in a basin, p = p(t) is mean areal rainfall "input" (units LT~"X), and e, q are mean areal evaporation loss and runoff respectively in the same units; e is set to zero (not an unreasonable assumption for winter months in the Dee basin, for which Lambert's model was originally developed), and S is expressed as a non-decreasing function of q, commonly S = k log q. The above equation can then be easily solved to give
5t+i= ^t/^tP-1 + (i-qtp-1)exP(-PAt/k)) (P>0> <2a^
= qt/(l+qtAt/k) (p=0) (2b)
where p is the accumulated rainfall in an interval of length At. Equation (2a) assumes that rainfall p is of constant intensity during each time interval; furthermore, both equations require the
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 183
use of p-^-L' t n e rainfall input over some earlier time period, so that L, the lag, is the second of two parameters (k and L) to be estimated.
The equations also include the quantity qt, observed runoff rate at the end of the preceding time interval, which is used as the initial value for the solution of equation (1) taken with S = k log q. Using equations (2) therefore, forecasts may easily be obtained by substituting the observed rainfalls Pt-L+1' Pt-L+2' •••î if forecasts of future rainfall Pt+1/ -••» can also be assumed, runoff can be forecast for times greater than the time of concentration of the basin.
As described in the above paragraph, Lambert's model is not adaptive in the sense that model parameters are re-estimated following each new observation of p t and q- ; forecasts are updated, however, by the use of new initial values qt given by the observed runoff values at the end of each time interval. Variants of the basic model (Green, 1979) have included the following:
(a) Replacement of Pt-L ky a weighted mean of the form
^a-eipt-L-i+Qpt-i+^u-^Pt-iw-i to allow for variation in the lag time taken by rainfall on different parts of the basin to reach the basin outfall. Green concluded that the benefits of rainfall "smoothing" by use of a weighted mean were small; there was some evidence that the forecast hydrograph with 0 = 0.6 was a better fit than with 6 = 1 for the basins of the Alwen (137 km2) and Ceiriog (113.7 km2) within the Dee basin.
(b) Automatic optimization of the parameters k, L. Lambert had used a graphical procedure to estimate these parameters; McKerchar (1975) calculated them by automatic optimization of a least squares objective function computed using several months of rainfall-runoff record, but the estimates of k, L obtained proved unsatisfactory for real time use.
(c) Replacement of the storage relation S = k log q by a piece-wise continuous curve. McKerchar (1975) experimented with the use of two linear relations S = kj log q, S = k2 log q such that k^ log q(t) = k2 log q (t) for some value q*. Three parameters, k]_, k2 and L, required estimation, but this variation on the basic model proved to be of small advantage.
(d) Replacement of the storage parameter k by a functional relation k = f(q). By expressing equations (2) as a pair in which the parameter k is given in terms of p, t, q-j- and qt+1' G r e e n w a s
able to derive real time estimates of the parameter k associated with the initial discharge q0 at which it occurred. The rising and falling limbs of the hydrograph were treated separately and a separate relation between k and q calculated for each limb. Real time estimates of the parameter k gave a significant improvement in the accuracy of forecasts for the Hirnant, Ceiriog, Gelyn and Alwen sub-basins of the River Dee basin.
INPUT-OUTPUT OR TRANSFER FUNCTION MODELS
Model formulation The models discussed in the preceding sections have to a greater or
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
184 P.E. O'Connell & R.T. Clarke
lesser extent been internally descriptive of the processes whereby system input (rainfall) is transformed to system output (runoff); however, for some real time forecasting applications, empirical relations between rainfall and runoff may be adequate. This section describes the forms which such relationships between rainfall and runoff might take, and the problems which the estimation of the parameters occurring in them present.
(a) A linear relationship between runoff gt and rainfall in current and earlier time periods pt, pt_j, .... If qt and pt are discrete time measurements of runoff and rainfall respectively, then a linear relationship between them may be written as
qt = aoPt + «lPt-1 + ... + ampt_m (3)
where a,0, a^, ---, otm are parameters; they can be regarded as the ordinates of a finite-memory impulse response for the river basin. Frequently, a pure time-delay between rainfall and runoff may exist, in which case the parameters {a} would be applied to the variables
(Pt-b' •••' Pt-b-m> • (b) A linear relationship between runoff gt and runoff in earlier
time intervals Qt-j, ... • Again writing qt for basin runoff, the form of this relationship is
qt = 3-L qt_! + B2 St-2 + • • - + Bn 3t-n (4>
where 3i, ..., Bn are parameters. (c) A linear relationship between runoff q-^r runoff in earlier
time periods g^-I' • - - , a.nd rainfall depths Pt, Pt-1' ••• • This relationship combines (3) and (4) to give
qt = 61 qt-1 + &2 9t-2 + - - - + 6r qt_r + 0)o pt_b + WX Pt_b_i
+ ws-l Pt-b-s-1 (5)
where {6} is a set of r parameters applied to r previous runoffs, {to} is a set of s parameters applied to s previous rainfalls delayed by b time intervals.
So far, neither the suitability of the above relationships for adaptive hydrological forecasting, nor the estimation of the parameters from the available hydrological data have been discussed. A discussion of the latter (estimation) problem is deferred to the next section; we now proceed to discuss the former.
If qt were "storm runoff" and p t were "effective rainfall", then (3) becomes the classic linear discrete-time unit hydrograph (UH) model. The unit hydrograph has been used extensively for hydrological design studies; however its use for real time forecasting in the form given by equation (3) requires reappraisal. Firstly, effective rainfall (as defined in the classical UH approach: e.g. Natural Environment Research Council, 1975) cannot be defined from total rainfall without prior knowledge of total storm runoff. Secondly, assuming that forecasts of storm runoff could be obtained, some assumption is then required about the magnitude of the baseflow contribution before a forecast of total flow qt could be derived. Thirdly, it is also desirable that corrections be applied to future forecasts according as telemetered streamflow information shows them to be in error; this could be achieved by allowing the ordinates of the UH to be time variant, involving their adaptive estimation at
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 185
each time point. However, this may not prove to be a wholly satisfactory solution; in any case, the number of parameters to be estimated (equal to the number of ordinates of the UH) may be large, suggesting that a parsimonious description of the rainfall-runoff relationship has not been achieved.
One way of overcoming the first and second difficulties discussed above is to apply a relationship of the form of (3) to total rainfall and total streamflow; this is no longer equivalent to the UH representation, and the set of ordinates {a} would not be expected to sum to unity since allowance must be made for losses. However, the rainfall-runoff process is recognized to be highly nonlinear (the existence of nonlinearity is at least implicitly recognized by the UH approach); hence a linear description of the rainfall-runoff process would be expected to exhibit severe shortcomings. In addition, the existence of baseflow in the absence of rainfall implies that the "memory" in equation (3) would have to be unduly large to avoid streamflow estimates becoming abruptly zero during a long recession prior to a storm event for which real time forecasts are required.
The use of a relationship of the form (4) which does not use rainfall, for forecasting streamflow in real time would require that telemetered flow information be available; however, while a high degree of interdependence usually exists between discrete-time measurements of streamflow over short time intervals, thus enabling future values to be forecasted, there is no means whereby the natural lag between the occurrence of rainfall and the response of streamflow can be utilized in making forecasts, nor can the model anticipate the occurrence of a peak discharge. Nonetheless, the model does have a natural self-correcting ability when telemetered streamflow information is available, and if rainfall telemetry proved unreliable, a relationship such as (4) could be used for forecasting.
The use of the relationship (4) for the real time forecasting of flow from rainfall offers the advantages of both models (3) and (4), in that it possesses a natural self-correcting ability as well as the ability to represent explicitly the lagged response of streamflow to rainfall. Thus, for real time use, it does not present any of the initialization problems associated with the use of (3); in fact, it can be shown that the ordinates of the impulse response CSJ_ are directly related to the parameters of (5) through the following equations (Box & Jenkins, 1970):
aj_ = 0 i < b
ai = 6lai-l + 62 a i -2 + - - - + <$rai-r + w0 i = b
a i = 6 l a i - l + 62 a i -2 + • • • + <5rai-r - aji.b i = b+1, b+2, . . . ,
b+s-1 (6)
«i = ôlai-l + 62ai-2 + ••• + ôra±_r i > b+s-1 Thus, b denotes the number of zero ordinates prior to the rising limb; the next (s-r) ordinates follow no fixed pattern and are absent if r>s-l. The remaining ordinates follow an rth order difference equation, i.e. ou for j>b+s-r with r starting values ou, b+s-rjjb+s-l. An important additional advantage of (5) is that it invariably results in considerably fewer parameters than (3).
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
186 P.E. O'Conneli & R.T. Clarke
The linearity inherent in (5) as a description of the rainfall-runoff process remains a shortcoming; however, there are a number of ways in which this can be overcome while still retaining the basic form of (5). One approach is to incorporate an auxiliary variable m t into the model as an additional "input" variable, which reflects antecedent basin conditions; equation (5) then becomes
9t = Sl^t-l + <S2gt-2 + • - - + <5rqt_r + W o P t_ b + WiPt-b-i + . . . +
"s-lPt-b-s-1 + 8 mt C7>
where 6 is an additional parameter. The variable m- can be defined in a number of ways, depending on the data available. Given that only antecedent rainfall data are available, an antecedent precipitation index (API) could be computed at each time point as
APIt = K API^JL + pt_]_ (8)
where K is a decay factor. If soil moisture deficit (SMD) data are available, then a catchment wetness index (CWI) can be computed (in mm) as
CWIt = 125 - SMD + APl5t (9)
where API5-|- denotes a 5-day API computed prior to and during a flood event as described in the Flood Studies Report (Natural Environment Research Council, 1975) and the constant, 125, is introduced to ensure that CWI remains positive since, in the UK, SMD rarely exceeds 125 mm. Variables such as APIt or CWIt can then be incorporated directly into (7) as the variable mt. The resulting model is then nonlinear in the systems theory sense but linear in the parameters (Clarke, 1973).
Instead of using m t as an additive term in (7), m t can be used to scale the precipitation input p t to yield a transformed input (or, in some sense, an effective rainfall input) as
pt* = m tp t (10)
the term pt* then replaces p t in (5). This approach has been adopted by Whitehead & Young (1975). Procedures for determining the efficacy of different transformations will be discussed later in the section on parameter estimation.
A further means of coping with the nonlinearity in basin response is through the use of thresholds; Todini & Wallis (1977) used this approach as follows. If m t again denotes some measure of antecedent basin conditions, and T is a threshold value of mt, then the threshold principle is applied to the rainfall input p t to yield two separate input series p ^ and P2t:
mt > T : Pit = P t; p2t = °
m t < T: p l t = 0; P 2 t = p t
The resulting model is then a multiple input type, i.e.
<*t = Ej=i 6jqt-j + Ei =i ^ ô 1 "ikP^t-bi-k ( 1 1 )
This model structure assumes that the response of a basin to rainfall under "wet" and "dry" conditions can be represented by separate impulse responses; this approach can be regarded as a piece-wise
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 187
linear representation of a nonlinear relationship. The type of nonlinearity represented by the above approach can be
regarded as "entry to storage" nonlinearity affecting the hydrograph rise (Moore, 1980); however, "release from storage" nonlinearity is also demonstrated in hydrograph recessions, where a steep fall from the hydrograph peak is followed by a long drawn-out tail. While such behaviour could perhaps be represented by making r in equations (5) and (11), for example, very large, this would not lead to efficient parameterization. An alternative approach is to apply a threshold also to the hydrograph recession; this could be the same threshold mfc applied to the rainfall, or some other threshold criterion which could even be some particular runoff value (Moore, 1980).
Transfer function relationships can also be useful for real time flow routing. Taking the case of a reach with a single upstream input, equation (5) can be applied to describe the relationship between upstream flow and downstream flow, and offers the same advantages over (4) and (3) as for the case of flow forecasting from rainfall. Equation (5) can easily be generalized to the multiple input case of several upstream tributaries. However, the assumption of linearity brings restrictions in the case of flow routing also. While equation (5) can represent the effect of flood wave attenuation, the variation in travel time with discharge observed for many rivers cannot be represented by a simple linear model. There are ways in which this difficulty might be overcome; for example the parameters could be allowed to follow some time variant relationship with flow if this could be identified from the available data, or alternatively, thresholds could be introduced as demonstrated above. In this case, a threshold might be identified with the discharge at which out-of-bank flow occurred.
Parameter estimation Least squares When applied to observed realizations of rainfall
and runoff, the relationships presented in the previous section cannot be expected to represent river basin behaviour exactly; accordingly, an error term will be required to complete the specification of each relationship. Thus, for example, equation (3) can then be written as
qt = a 0p t + a l P t_! + ... + ampt_m + e t (12.
where e t is an error term representing the effects of model inadequacy, and possibly measurement noise in the variables qt and p t. As will be seen subsequently, statistical assumptions to be attributed to the noise term £t, and the presence or absence of measurement noise in the rainfalls p t, will have a critical bearing on the approach to be adopted to parameter estimation. If p t is assumed to be a deterministic error free variable, and £ t is assumed to be a serially uncorrelated random variable, then the estimation of the parameter set {a} becomes a straightforward problem of linear regression analysis to which ordinary least squares (OLS) can be applied; the resulting least squares estimates will be unbiassed and will have sampling variance at least as small as that of any other linear unbiassed estimator, regardless of the distribution of £-£• Estimation of the parameters can be performed
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
188 P.E. O'Conneil & R.T. Clarke
using the conventional "block data" processing procedure, or, alternatively, the recursive least squares algorithm (see, e.g. Young, 1974; or Weiss, 1980) can be used to calculate the regression coefficients recursively. If the latter approach is to be adopted, then time variance can, if necessary, be introduced into the parameters {a}. Writing (12) in vector notation as
qt = h^6t + e t (13)
with
i}l = fct pt_! . - - pt_j and
§t = K al • • • aml then time variance can be introduced into the parameters by writing
§t = et-i + ït ( 1 4 )
where w t is a vector of independently distributed Gaussian random variables with variances chosen to impart variation to each of the parameters in the vector 9. The parameters thus follow a random walk in which their values in any particular time interval are independent of their values in earlier time intervals. The resulting algorithm is known as dynamic least squares (e.g. Beck, 1979). If wt = O and the parameters are assumed time invariant, then as the number of data increase, the variance of the parameters decreases and they eventually assume steady state values. There is then little point in continuing to estimate the parameters recursively on-line, and the steady state values can be used. If, on the other hand, it is desired to impart an amount of variation which is intermediate to the above alternatives, then the recursive least squares algorithm can be modified to give more weight to recent data and less weight to past data in estimating the parameter vector 6 at each time point (Beck, 1979). Rather than minimizing the function
JH- = E t , (q, - hfe,)2 (15) Jt i =l <*i - *i§i>
a modified loss function defined as
a't = £^ = 1 at _ 1 (q± - hje^)
2 o<a<l (16)
is minimized. The weighting factor a is usually chosen to be just less than unity; smaller values of a give a more rapid forgetting of past data, resulting in more time variation in the parameters.
The introduction of time variance into model parameters at the estimation stage is an option which recursive estimation algorithms can provide and which "block data" processing procedures cannot; the merits of the particular types of time variance considered here will be discussed later when reviewing some of the applications of adaptive estimation techniques to hydrological forecasting.
Turning again to the problems of parameter estimation, we now set out the effects of alternative assumptions about the statistical properties of the variables p t and £•(- on the estimation of the parameter vector 6. If p t is still regarded as a deterministic
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 189
variable but {E-^} is an autocorrelated sequence of random variables, then the application of OLS to (12) will result in unbiassed parameter estimates, but their sampling variances may be inflated, and as a consequence, forecasts made using these estimates will be imprecise. If, however, the autocorrelation structure of the noise term e . can be specified, then generalized least squares could be used to provide statistically efficient estimates of the parameters in the vector 6 (Johnston, 1972). Specification of the autocorrelation structure of e t will present difficulties in practice; approximate methods of overcoming these are discussed in Johnston (1972) for specific cases.
So far it has been assumed that p t, Pt-i/ •••/ i n (13) are
measured without error; if however, the true rainfall p t is regarded as subject to measurement noise Ç t such that
Pt = Pt + 5t (17>
and if the presumed model is now
St = ao^t+al^t-l + -•• + amPt-m+et <18>
with e t assumed to be serially uncorrelated, then substitution of (17) into (18) results in
gt = ZT=o aiPt-i+et-Ei=oaiÇt-i = Zi=0aiPt-i+et <19>
where {et} is now an autocorrelated noise sequence which is also cross-correlated with the error-corrupted rainfall observations p-j-(the measurement errors {Ç^ï a r e common to both). This cross-correlation causes the estimates of {a} obtained by OLS to be inconsistent. The parameters aj_ will in general by underestimated; for the simple case of m = 1 in equation (17) , the limiting estimate of ct]_ will be
- D al
a Z J—- (20) 1 l+a|/a^
where at denotes the variance of the measurement error in p t, the estimated rainfall, and Op is the variance of the true rainfall values, assuming that the measurement errors and true values are uncorrelated (Johnston, 1972). If a| is 10% of Op, then a]_ will be underestimated by about 10%. Thus measurement errors in rainfall can, if a model of the form (18) is to be postulated, pose a serious problem for least squares estimation, and alternative estimators are required for its solution.
So far, attention has been confined to the estimation problems presented by the relationship in equation (3) in order to demonstrate as simply as possible the effects of different assumptions on parameter estimation by least squares. It was noted earlier that models of the form (5) have more desirable characteristics for real time forecasting than (3); however the inclusion of lagged values of runoff q- on the right-hand side of equation (5) introduces additional complications for parameter estimation. If the model is written as
St = ôiqt-l + (52qt-2 + ••- M t - r + u0Pt-b + ulPt-b-l + -•• +
"s-lPt-b-s-1 + Et ( 2 1 )
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
190 P.E. O'Connell & R.T.Clarke
where qt, qt-1/ ---/ are measured values incorporating the noise terms £t, £t-l •••/ respectively, then, provided that the £ t are serially uncorrelated, least squares estimates will be consistent and asymptotically efficient. If, however, the noise term £t is autocorrelated, then £ t will be cross-correlated with qt-1» qt_2,•••, as these variables contain the noise terms Et-1» £t-2>- - - » respectively, and least squares estimates will be inconsistent. The problem is further accentuated if the variables Pt-b' Pt-b-1' •••' a r e subject to measurement noise.
The foregoing has attempted to illustrate the limitations of ordinary least squares estimation under different assumptions for the explanatory variables and noise terms in the various models considered; estimation techniques have been developed which overcome some of these problems, some of which are described in brief in thé following sections.
Transfer function-noise model estimation Two approaches which overcome the difficulties associated with applying ordinary least squares to the estimation of the parameters in equation (5) will be considered here. Both approaches assume the same overall model structure in which observed discharge is written as the sum of a deterministic component and a stochastic component, i.e.
qt = qt + Tlt (22)
where the deterministic component is given by the transfer function
model
<?t = Vt-l + 629t-2 + • • • + M t - r + 0)opt_b + UiPt.b.! + . . .
us-lPt-b-s-l (23)
Introducing the backward shift operator B defined as Bgt = ?-(--l' t*ie
transfer function model (5) may be written as
6(B) gt = u(B)pt_b (24)
where
6(B) = (1 - 6LB - 62B2 - ... <5rB
r)
s-1,
and
U)(B) = (u)0 + tO-jB + ... Ws_1B'
The stochastic component or noise term rit in (22) is assumed to be related to white noise at through
n t = * i n t - i + M t - 2 + ••• V i t - p + a t ^ e i a t - i - ••• -8qat_g (25)
i.e. an ARMA (p,q) model (Box & Jenkins, 1970) which can be written in concise form as
(26) <I>(B) n t =
<f>(B) = (1
9(B) - (1
0(B) a t
- (f)1B - $2B2 -
- e1B - e2B2 -... v P )
. . . eqBq)
Substituting for qt and rit fr o m (24) and (26) into (22) then yields
the composite transfer function-noise model
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 191
Fig, 1 The transfer function-noise model.
q t " 6 (B) p t - i 9(B) (()(B)
(27)
which is depicted graphically in Fig. 1. The following approach to the estimation of the parameters in the
polynomials 6(B), 03 (B) , 6(B) and cj> (B) is advocated by Box & Jenkins (1970). Assuming vectors of starting values q0, p 0, aQ prior to the commencement of the observed data series were available, then, given the data, for any choice of the parameter vector B = {<$]_, ..., ôr; 0)o, ..., ws_j_, b} values of at = a-(-(B|q0, p Q, a0) can be calculated. Under the assumption of normality for the a's, close approximations to the maximum likelihood estimates of the parameters can be obtained by minimizing the conditional sum of squares
so<§> = E t = 1 at 3|q0,P0'ar (26
such estimates will be consistent and asymtotically efficient. Minimization of the function S0(3) can be achieved using a standard nonlinear least squares computer program, while initial estimates of the parameters can be derived using cross-correlation and autocorrelation analyses as described in Box & Jenkins (1970).
The above method of estimating the parameters in (27) is a "block-data" processing procedure which would be performed "off-line" prior to using the model for real time forecasting. A procedure in which the parameters of the model (27) can be estimated recursively has been proposed by Young et al. (1971); their approach is to decouple the estimation of the transfer function parameters and of the noise model parameters, giving two separate estimation problems. The transfer function parameters are estimated as follows. If the measurement and parameter vectors in the linear model (13) are defined as
hj: = [qt-i/ qt-2'
ft = D5!' «2- •••
9t-r<' Pt-b' Pt-b-1' -••' Pt-b-s-lJ
Js-U
(29)
, u r , uiQ, o.1, ..
then, as noted in the previous section, least squares estimates will be inconsistent. If, however, equation (22) is used as the basis for the estimation, then the resulting model is
«It = £t§t + Tit (30)
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
192 P.E. O'Connell & R.T. Clarke
where
hi = [gt-l< St-2' •••' ît-r' Pt-b' Pt-b-1' •••' Pt-b-s-l]
As the elements of h-t are deterministic variables uncorrelated with Ht' this model provides a basis for consistent least squares estimates of the elements of 6-f However, the elements gt-1' gt-2' ••- of hj. are unknown; to overcome this problem, an estimate of ht, denoted by h t and referred to as an instrumental variable (IV) vector, can be provided which is defined to be highly correlated with h t but uncorrelated with the noise r|t.
The instrumental variable vector is generated from a linear transfer function model as follows. Firstly, the least squares estimation of the parameter vector 0 is formulated in recursive form so that the estimate is updated at each time point. After all the available data have been processed, an estimate of 6 is available; this can then be used to generate a set of estimates of gt
recursively as
qt = h^B (31)
This set of estimates is then used in conjunction with the (assumed) deterministic rainfall input as the IV variable for the next "pass" through the data, i.e.
hl = fct-1' ât-2' ••- St-r' Pt-b' pt-b-l Pt-b-s-l] ( 3 2 )
and this procedure is repeated until stability is achieved in the estimate of 6. Details of the recursive estimation algorithm are given in Young et al. (1971).
Once convergent parameter estimates and instrumental variables have been obtained, an estimate of the noise series T). in (22) can be obtained as
fit = q t - gfc = q t - §£g (33) An ARMA (p,q) model is then applied to the noise series rit; the parameter estimates are obtained recursively and approximate to maximum likelihood estimates. The overall estimation procedure is referred to as the Instrumental Variable-Approximate Maximum Likelihood (IVAML) algorithm (Young et al., 1971).
While the IVAML algorithm employs recursive estimation procedures, it cannot be used on-line in the form described above as several passes or iterations through the data are required before final parameter estimates are obtained. An alternative version has been outlined (Young et al., 1971) which is suitable for on-line use and which does not require successive passes through the data. However, for off-line estimation, the former version has been used more extensively (e.g. Venn & Day, 1978; Moore & O'Connell, 1978; Moore, 1980) and it may be that this provides more efficient parameter estimates than the on-line version. The efficiency of the parameter estimates derived from the use of instrumental variables depends on the correlation between the latter and underlying noise-free values. If the correlation is high, then the estimates will be efficient, although they will have larger sampling variance than least squares estimates; if the correlation is weak, then estimates with poor efficiency will result. Procedures for improving the efficiency of parameter estimates derived using instrumental
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 193
variables have been developed by Young & Jakeman (1978). Both the Box-Jenkins and IVAML estimation procedures provide
consistent parameter estimates; it is not clear how the efficiencies of parameter estimates from the two procedures would compare, although under conditions of normality, the Box-Jenkins estimates approximate closely to maximum likelihood estimates. If a model of the form of (27) is to be used to forecast streamflow from rainfall, then, as previously noted, some means of accounting for the non-linearity in basin response is desirable. One of the advantages of the IVAML recursive estimation approach is that the efficacy of different transformations to the input pt (equation (lO)) can be judged on the basis that time invariant parameter behaviour should result from the transformation (Whitehead & Young, 1975). This is a useful attribute of recursive estimation and is quite separate from considerations about on-line use.
In considering ways in which nonlinearity might be introduced into transfer-function models in the section on model formulation, multiple-input transfer function models were employed (equations (7) and (11)). The IVAML estimation procedure can be extended in a straightforward manner to cope with the case of multiple inputs by extending the measurement and parameter vectors in equation (29) to include the additional inputs and corresponding parameters; however, this approach does not account for the cross-correlation structure of the multiple inputs. The multivariable IVAML algorithm described by Young & Whitehead (1977) does account for the full cross-correlation structure; however, this is achieved at the expense of considerable computational complexity.
Self-tuning algorithms Self-tuning algorithms have been developed in the control theory field over a number of years; their application to problems in real time hydrological forecasting and control has recently been explored by Ganendra (1979, 1980) amongst others. The starting point for the development of self-tuning algorithms is the nth order stochastic difference equation representation of a system (i.e. river basin):
A(Ç_1) qfc = B(Ç_1) p t + CiC1) e t (34)
where Ç- is the backward shift operator Ç q. = qt_j_ (the symbol Ç is used here rather than the Box-Jenkins symbol B in order to preserve the notation used in the literature on self-tuning algorithms) and A(Ç - 1), B(Ç-1) and CtÇ""1) are polynomials in Ç^1
defined as
A(Ç_1) = 1 + a;LÇ-1 + ... + anÇ~
n
BiC1) = 1 + b 1 Ç _ 1 + . . . + b n Ç ~ n (35)
etc""1) = l + c 1ç_ 1 + ... + cnç~
n
We now discuss the self-tuning predictor algorithm which seeks to provide minimum mean square error, k-step-ahead forecasts of runoff qt. However, rather than trying to estimate the parameters of the system model in (34), the problem posed is the estimation of the parameters of the minimum mean square error predictor. The derivation of the self-tuning predictor is given elsewhere (see, e.g., Wittenmark, 1974) and will not be repeated here; the predictor for qt is
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
194 P.E, O'Connell & R.T.Clarke
St+k|t = f Pt + }& Vt C36)
where q^+kIt denotes the k-step-ahead prediction of qt+k given rainfall and runoff up to time t, G and F are further polynomials in the backward shift operator Ç-1 and V t is the k-step-ahead prediction error defined as
Vt = qt - qt|t-k
Equation (34) may be rearranged to give
%+k|t = ( 1~ A F ) <ît+k|t + BFPt + G vt ( 3 7 )
The estimation of the parameters in the predictor is approached by
(3E
wr i t i ng
( 1 -
BF
G =
and de
£ f
-AF) =
= h +
= Yi +
^fining
Ct]_Ç
• 32?
Y2S~
-1 +
r1 + -1 +
a2r2
. . . +
• measurement
= R t+k- i | v t ,
= [«!' • • * /
t - i ' v t -
' a n
. . - ,
•£+l]
," h'
+ . . . +
+ emrm
• y ^ - z
anrn
and parameter vec
q t+k -n t - n ' p t ' •
. . . , e m ; Yi, •••
: tors as
• • ' Pt-:
' YjJ
(39)
so that the predictor equation can now be written as
qt+k|t = 5? § (4o)
The predictor parameters are estimated by applying recursive least squares to (40); it has been shown that if C(Ç_1) = 1 in (34) then consistent parameter estimates are obtained and the predictor parameters will converge to those of the optimal predictor. The predictor performance is judged according to whether or not the prediction errors satisfy
CVV(T) = E{vfc V t + T} = 0 T = k, ..., k+l-1
CV-(T) = E{vt+T qt+k|t} = 0 T = k k+n (41)
However, even if C(Ç_1) ^ 0, the predictor also gives good results. Ideally, the sequence of prediction errors {vt} should obey the above relationships for all lags, since all relevant information has then been abstracted from the data and knowledge of one error gives no information about the next.
Additional inputs can be readily incorporated into the predictor by extending the data and parameter vectors in (39); Ganendra (1979) has incorporated an auxiliary input called the PERLOG into the self-tuning predictor as a measure of antecedent basin conditions. The PERLOG is defined as
p£(t) = -^ (loga(t)) = ^-~ A (q(t)) (42) dt 3 e a q(t) dt
where q(t) is continuous streamflow; a lagged finite difference approximation to (42) can be calculated for prediction purposes. The rationale behind the use of the PERLOG is that it amplifies the rate of change of streamflow with time; thus, for example, the
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 195
transition from recession to rising limb of a hydrograph would be expected to reflect antecedent conditions, which in turn would be reflected in the behaviour of the PERLOG.
STATE-SPACE MODELS AND KALMAN FILTERING
In recent years, state-space modelling and Kalman filtering have attracted considerable attention in hydrology. These techniques find their origin principally in the control theory field; because of the increasing use of telemetry in the management of water resources systems, hydrologists have been exploring the potential of these techniques for real time forecasting and control. This has led to a proliferation of papers in this subject area; no attempt will be made to review this extensive literature here. Instead, the main assumptions underlying the state-space/Kalman filtering approach will be discussed together with some of the problems encountered in applying these methods to the representation of river flow behaviour.
The merit of the state-space modelling approach is that it allows many different models (physics-based, conceptual or black box) to be cast within the mathematical framework of two equations, a system equation and a measurement equation. Equally both deterministic and stochastic models can be represented within the same framework; the stochastic representation is employed when it is required to estimate the "state" in the presence of both system and measurement noise. The concept of the system "state" is introduced primarily for mathematical convenience; although the variables comprising the states need not be physically measurable quantities, they are related to observed outputs through the measurement equation. State-space models may be formulated for systems which are linear or nonlinear and observed in discrete or continuous time. For a linear dynamic system the discrete-time system and measurement equations are
ït = tt-l*t-l + A tu t + Ttwt (43)
3t = !?t*t + Yt <44>
where xt, qt are (nxl) and (mxl) state and measurement vectors, u t
is an (rxl) vector of deterministic inputs or control variables, w t and vt are (mxl) vectors of system and measurement noise, $t-l i
s
an (nxn) transition matrix, and At, T-ç- and H t are (nxr) , (nxm) and (mxn) weighting matrices. In the standard exposition of Kalman filter theory, the system and measurement noises are independently and identically distributed Gaussian random variables with the following properties:
E(wt) = O; E(vt) = 0; E(wtw£) = Q 6 t f k
E^YtYk^ = ? ^t k' E tYÏD = 2 v t » k
where 6 t j , is the Kronecker delta defined as
«t,k = 1 fc = k
«t,k = ° t * k
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
196 P.E. O'Connell & R.T. Clarke
Finally, the system noise w t is assumed to be independent of xs
for s < t, and the measurement noise v t is assumed independent of xt. One of the advantages of the state-space formulation is that the system and measurement noises are separated; the problem posed then is the estimation of the system state xt in the presence of the measurement noise v- - This problem is solved by applying the linear Kalman filter to derive filtered state estimates, x
t t' recursively using the measurements obtained up to time t.
A prerequisite for the use of the linear Kalman filter is that the model be expressed in the form of the equations (43) and (44); however, the flexibility of this framework is such that models can be so expressed in a number of different ways. This is well illustrated by the linear transfer function model (5); taking the specific case of r = 2, s = 3 and b = 1 as an example, then one possible formulation is
9t
L t-lJ
IS, <5. <?t-l
L9t-2J
ktl = [i o] .fft-lJ
0 u l u 2
0 0 0
P t - 1
P t - 2
J ? t - 3 -
+ Wj_
0 J t
+ v.
If p t is regarded as rainfall, and q. as streamflow, then the above formulation assumes rainfall to be a deterministic error free input variable. However, rainfall is subject to measurement error and is computed in many cases as a lumped areal average which introduces further errors. Also, if forecasts in advance of one-step-ahead are to be made from the above model, then some assumption has to be made about the future behaviour of rainfall. An alternative approach is to treat rainfall as a stochastic process; for illustrative purposes a third order autoregressive (AR(3)) process will be used to represent rainfall. The state-space model is then written as
vt g t - i
Pt
Pt-i
Pt-2_
=
' « i
o o 0
0
62 Yo Yi ^
1 0 O 0
o a^ a 2 a 3
0 1 0 0
0 0 1 0
? t - l "
«Tt-2
* t - l
Pt-2
Pt-3
+
r i o"
0 1
0 0
0 0
0 0
V w2
tJ
1 0 0 0 o"
0 0 1 0 0_
"?t
* t - l
Pt-1
Pt-2
+ ~V1
- v 2.
Thus, short term rainfall forecasts can be provided within the model using the statistical structure of rainfall, and rainfall measurement
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 197
noise is also accounted for. Models with multiple inputs such as those discussed in the section on model formulation can also be easily accommodated within the state-space framework.
A third possible formulation of (5) is
'h $2 u>0
0J-j_
_w2_
=
t
Ï&1
62
w0
w l w 2 .
% = ht-l^t-2Pt-lPt-2Pt-3]
where the model parameters are assumed to constitute the state vector, and in the above formulation are assumed time invariant; time variance can be introduced by adding noise terms to the parameters in the system equation. If the Kalman filter is applied to this formulation, this is exactly equivalent to recursive least squares,; the noise term v t includes both the effects of model error and measurement noise. Thus, if v, is not an independently distributed random variable, inconsistent parameter estimates will result. This formulation is frequently adapted when applying Kalman filtering to models of the form of (3) and (5), and it is not always appreciated that consistent parameter estimates can result.
For the formulations (43) and (44), nothing has been said so far about model parameters such as $t-i; i n some control theory applications these may be assumed known, but for application to hydrological systems the parameters must always be estimated. There are a number of ways in which both the states x^ and the model parameters may be jointly estimated; one is to augment the state vector with the parameters and solve the resulting estimation problem using the extended Kalman filter. A second approach is to decouple the state/parameter problem into two linear estimation problems which can be solved using two Kalman filters interactively (Todini, 1978; Todini et al., 1980). Other approaches are also possible.
If the dynamics of system behaviour are assumed to be nonlinear, then the system equation can be written in discrete time as
*t = f t-l^t-l* + wt <45>
while if the measurement system is nonlinear this can be written as
qt = h(xt) + v (46)
where f (.) and h(.) are nonlinear functions. The extended Kalman filter can be applied to (45) and (46) by deriving Taylor series expansions of f(.) and h(.) retaining only terms linear in xt, and then applying the linear Kalman filter to the linearized problem.
&1
W2
+ v,
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
198 P.E. O'Connell & R.T. Clarke
As an example, a simple nonlinear rainfall runoff model formulated by Mandeville (1975) and Moore & Weiss (1980) can be written in differential equation form as
dhk b
— = a{c p* exp(-d.S75) - hk}hk (tk<t<tk+l) (47)
where hk(tk) = qk, p* is a lagged and distributed rainfall term, S75 is a term representing soil moisture deficit and a, b, c and d are model parameters. This nonlinear model was written in state-space form by Moore & Weiss (1980) as
ïk+1 = *k
^k+l = ^ ^ k + l ' + vk+l
where xk = [abed],- i.e. the parameters constitute the state vector and the dynamics of system behaviour are represented by the measurement equation. The extended Kalman filter is then applied to the linearized measurement equation to obtain recursive estimates of the model parameters; flow forecasts are then obtained at time points tk+j_, tk+2 by solving the differential equation (45) subject to the initial condition qk = h k(t k).
Finally, there is the problem of estimating the system and measurement noises Q and R, as again these will be unknown for hydrological systems. This is a problem which sometimes appears to be glossed over when Kalman filtering is applied in hydrology; yet the performance of this estimation procedure depends critically on the proper specification of Q and R. Algorithms for estimating them in the case of a discrete linear dynamic system are presented by Todini et al. (1980); they are applicable to systems where the dimensionality of the system and measurement noise matrices are equal as is frequently the case with models of the transfer function type.
CURRENT EXPERIENCE AND FUTURE DEVELOPMENTS
This paper has so far presented a spectrum of modelling procedures for real time application, a discussion of assumptions underlying them, and some problems encountered in their application. The procedures range from conceptual hydrological models which use relatively unsophisticated parameter estimation procedures (such as the River Dee model), to simple input-output type models which use relatively sophisticated state filtering and parameter estimation procedures from the fields of statistics and control theory. Most of the research effort in recent years has been spent in exploring the use of sophisticated estimation procedures with relatively little effort devoted to the improvement of the underlying models. It is perhaps a pertinent juncture to assess whether this balance of effort is justified in terms of improved accuracy of forecasts. Some evidence for this assessment comes from a workshop held at the Institute of Hydrology in 1977 to establish the relative merits of different models and parameter estimation algorithms. The results of the comparison study carried out during this workshop are described in O'Connell (1980); extracts from these results will be
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 199
Table 1 Means, variances and values of Q and S-tests applied to 1-h-ahead forecast errors for the Hirnant for the period 1 -30 December 1972 (after O'Connell, 1980)
Statistic Model A Model B Model C Model D
No. of degrees of freedom (n.d.f) X0.05 Q/(n.d.f)
No. of degrees of freedom (n.d.f) v 2
Xo.os S/(n.d.f)
0.007 0.078
308.4
28
41.34 11.01
225.5
28
41.34 8.04
0.041 0.075
288
27
40.11 10.70
191
27
40.11 7.07
0.018 0.269
100
21
32.67 4.76
228
21
32.67 10.86
0.042 0.153
83
23
35.17 3.61
80
23
35.17 3.48
discussed here. Four models were applied to 2 months of half-hourly rainfall and runoff data for the Hirnant sub-basin of the River Dee; the first month of data was used for model identification and for obtaining preliminary parameter estimates; 1-h-ahead forecasts were then provided by each of the four models for the second month. The four models compared were
Model A: a version of the River Dee sub-basin model; Model B: the simple nonlinear conceptual model used by Moore &
Weiss (1980) and described briefly in the last section; Model C: a linear transfer function model with time-variant
parameters; Model D: the self-tuning predictor with PERLOG input described
in the section on self-tuning algorithms. Several properties of 1-h-ahead forecast errors were computed;
some of these properties are summarized in Tables 1 and 2, while plots of the errors for selected flood events are shown in Figs 2 and 3. Values of the Q and S statistics (Box & Jenkins, 1970) testing (a) that the errors are uncorrelated and (b) that they are not cross-correlated with the input, are presented in Table 1
Table 2 Statistics of 1 -h-ahead forecasts for Hirnant events on 4-6 December (event 1 ) and 11-12 December (event 2) 1972 (after O'Connell, 1980)
% error at peak Timing error in forecasting peak Maximum % error on rising limb (as % of observed peak)
Model A
Event 1
5.0 0
9.4
Event 2
-7.5 0
10.4
Model B
Event 1
10.1 0
10.1
Event 2
-2.1 1 h late
9.9
Model C
Event 1 Event 2
-15.3 -10.5 0 0
-8 .4 -9 .6
Model D
Event 1 Event 2
2.2 -1.9 0 0
12.0 11.5
Note: positive errors denote underestimation.
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
200 P.E, O'Connell & R.T. Clarke
3 -
2
1 -
Ë-lJ
- 2 -
- 3 -
3 -
2 -
1
v 0 CO
n E 1
- 2
- 3 J
3 -
2 -
1 -
o °
- 2 -
- 3 -
3
2 -
1 -
o °
- 2 -
- 3 -
Fig.2 Hirnant
MODEL A
I 10
I 10
' I ' 10;
I 20
MODEL B
30 I
40 I
SO I
6 0
I 20
MODEL C
I 30
I 4 0
I I SO 6 0
I 20
MODEL D
I 30
I 50
I 70
I 70
•60 70
I 10 3 0
1 4 0
I 5 0 6 0
I 7 0
I 80
I 80
I 80
I 80
One-hour-ahead forecast errors for the four models for event 1 (4-6 December 1972) on the (after O'Connell, 1980). (Note: positive errors denote underestimation.)
together with the mean and variance of the 1-h-ahead forecast errors. Since the different models had different numbers of parameters, the Q and S statistics were divided by the number of degrees of freedom to allow comparison of the results. Table 1 illustrates that
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 201
2 -
1 -
'(/) E_,_
oJ
MODEL A
I ! ! / | ; •-•J
' • • '
I 10
I 3 0
I 4 0
I 50
1
V ° 'to
m
E - 1
-
-
MODEL B
A ; i
I 10
I 20
I 30
I 4 0
I 50
!] T 0 -
n
E .
I 10
I 20
I 30
I 40 50
2 - 1
0) M
e 1 -
i 10
,. . .- • ' . . . . • ' - . > - . ' f -
i
20
I
3 0
I
4 0
I 50
Fig. 3 One-hour-ahead forecast errors for the four models for event 2 (11-12 December 1972) on the Hirnant (after O'Connell, 1980). {Note: positive errors denote underestimation.)
Models A and B have the smallest error variances (minimum variance is clearly a desirable property for forecast errors); however Model D performs best in terms of the standardized Q and S statistics. The statistics presented in Table 2 pertain to two particular flood events within the 1-month test period; the results illustrate that the magnitudes and signs of the forecast errors computed for these particular events are all rather similar, with
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
202 P,E. O'Connell & R.T.Clarke
Model D proving marginally best. The plots in Figs 2 and 3 illustrate the similarity in structure of the forecast errors for each of the four models, all of which underestimate on the rising limb and overestimate on the peak and early recession. These plots illustrate (a) that the variance of the forecast errors is non-stationary and (b) that the errors are autocorrelated in the neighbourhood of flood events. For Models B and C, these results indicate that the assumption of independently and identically distributed random variables for the system noise in the state-space formulation (last section) is not justified; this is also true of the forecast errors in the self-tuning predictor (Model C). The underestimation-overestimation phases in forecasting are not confined to these particular results and are also noticeable, for example, in Todini & Wallis (1978) and Andjelic & Szollosi-Nagy (1980).
The apparent nonstationary nature of the forecast errors may derive from the manner in which noise enters these systems. When streamflow is in a state of recession, then modelling errors derive only from inaccuracy in modelling the recession and there is no rainfall forcing function. However, during flood events, noise in discrete "lumps" enters the system through the inadequate definition of rainfall inputs and through errors in modelling the system's response to rainfall. If, then, the system noise is not Gaussian, linear predictors will not be optimal. Weiss (1981) has recently investigated the performance of linear Gaussian predictors when applied to continuous linear systems activated by discrete shot noise; the one-step-ahead forecast errors were found to display remarkably similar characteristics to the errors in Figs 2 and 3. This may suggest that linear predictors are not appropriate for hydrological systems; however, the complexity introduced by nonlinear predictors is formidable.
The above results suggest that the standard assumptions of independence and identical distribution for system noise are not justified for hydrological systems, and that failure to recognize this may mean that the attention paid, e.g. to the problems of inconsistent parameter estimation may not be justified. Given that the River Dee model, which did not employ recursive estimation but which has a natural self-correcting mechanism, performed as well in forecasting flood peaks as the other models, the computational effort required to apply, for example, the extended Kalman filter to Model B can hardly be justified. In addition the steady state estimate of the error forecast variance provided by the Kalman filter
p t + i | t = <f>tpt|t *l + rt Q ^
will not have much meaning when the error forecast variance is clearly time variant. This might suggest that time variant parameter models are to be preferred; even if this is so, the random-walk time variance imparted to the parameters in Model C tends to lead to erratic forecast behaviour over flood events. Furthermore, the use of exponential weighting of past data (see section on parameter estimation) included in Model B does not appear to have led to any great improvement.
The above discussion suggests that there are still considerable
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 203
unsolved estimation problems in real time forecasting, but it is not clear to what extent their solution would result in improved forecasts. It may be more beneficial to seek a better representation of the spatial variation in rainfall and its effect on streamflow response, and in improving the structure of real time forecasting models than to expend effort in solving estimation problems. Information on where efforts will be best rewarded can only be obtained by feedback from case studies.
REFERENCES
Andjelic, M. & Szollosi-Nagy, A. (1980) On the use of stochastic structural models for real time forecasting of river flow on the River Danube. In: Hydrological Forecasting (Proc. Oxford Symp., April 1980), 371-380. IAHS Publ. no. 129.
Beck, B. (1979) System identification estimation and forecasting of water quality. Part I: Theory. WP-79-31, International Institute for Applied Systems Analysis, Laxenburg, Austria.
Box, G. E. P. & Jenkins, G. M. (1970) Time Series Analysis Forecasting and Control. Holden-Day, San Francisco.
Cameron, R. J. (1980) An updating version of the Muskingum-Cunge flow routing technique. In: Hydrological Forecasting (Proc. Oxford Symp., April 1980), 381-387. IAHS Publ. no. 129.
Clarke, R. T. (1973) A review of some mathematical models used in hydrology, with observations on their calibration and use. J. Hydrol. 19, 1-20.
Ganendra, T. (1979) Real-time forecasting and control in the operation of water resource systems. PhD Thesis, University of London.
Ganendra, T. (1980) The self-tuning predictor. In: Real-Time Hydrological Forecasting and Control (ed. by P. E. O'Connell), (Proc. 1st International Workshop, July 1977). Institute of Hydrology, Wallingford.
Green, C. G. (1979) An improved subcatchment model for the River Dee. Report no. 58, Institute of Hydrology, Wallingford.
Jamieson, D. G., Wilkinson, J. C. & Ibbitt, R. P. (1972) Hydrologie forecasting with sequential deterministic and stochastic stages. Proc. Int. Symp. on Uncertainties in Hydrologie and Water Resources Systems (Tuscon), Vol. 1, 177-187.
Johnston, J. (1972) Econometric Methods, 2nd edn. McGraw Hill, New York.
Jones, D. A. & Moore, R. J. (1980) A simple channel flow routing model for real time use. In: Hydrological Forecasting (Proc. Oxford Symp., April 1980), 397-408. IAHS Publ. no. 129.
Kitanidis, P. K. & Bras, R. L. (1978) Real-time forecasting of river flows. Report no. 235, R. M. Parsons Laboratory for Water Resources and Hydrodynamics, MIT, Cambridge, Mass.
Lambert, A. O. (1969) A comprehensive rainfall/runoff model for an upland catchment area. J. Instn Wat. Engrs 23 (4), 231-238.
Lambert, A. O. (1972) Catchment models based on ISO functions. J. Instn Wat. Engrs 26 (8), 413-422.
Lambert, A. 0. (1978) The River Dee regulation scheme: operational experience of on-line hydrological simulation. International
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
204 P.E, O'Connell & R.T. Clarke
Symposium on Logistics and Benefits of using Mathematical Models of Hydrological and Water Resources Systems, IAHS/WMO/IBM, Pisa.
Lowing, M. J., Price, R. K. & Harvey, R. A. (1975) Real-time conversion of rainfall to runoff for flow forecasting in the River Dee. Proc. Dee Weather Radar and Water Management Conference (Chester). Water Research Centre, Medmenham.
McKerchar, A. I. (1975) Subcatchment modelling for River Dee forecasting. Report no. 29, Institute of Hydrology, Wallingford.
Mandeville, A. N. (1975) Non-linear conceptual catchment modelling of isolated storm events. PhD Thesis, Dept. Environmental Sciences, Univ. of Lancaster.
Moore, R. J. (1980) Real-time Forecasting of Flood Events Using Transfer-Function Noise Models: Part 2, Report to Water Research Centre. Institute of Hydrology, Wallingford
Moore, R. J. & O'Connell, P. E. (1978) Real-time Forecasting of Flood Events Using Transfer-Function Noise Models: Part 1, Report to the Water Research Centre. Institute of Hydrology, Wallingford.
Moore, R. J. & Weiss, G. (1980) A simple non-linear conceptual model. In: Real-Time Hydrological Forecasting and Control (ed. by P. E. O'Connell) (Proc. 1st International Workshop, July 1977). Institute of Hydrology, Wallingford.
Natural Environment Research Council (1975) Flood Studies Report, Vol. 1, chap. 7, 513-531. Natural Environment Research Council, London.
O'Connell, P. E. (1980) Real-Time Hydrological Forecasting and Control (Proc. 1st International Workshop, July 1977). Institute of Hydrology, Wallingford.
Szollosi-Nagy, A (1976a) An adaptive identification and prediction algorithm for the real-time forecasting of hydrologie time series. Hydrol. Sci. Bull. 21, 163-176.
Szollosi-Nagy, A. (1976b) Introductory remarks on the state-space modelling of water resource systems. RM-78-73, International Institute for Applied Systems Analysis, Laxenburg.
Szollosi-Nagy, A., Todini, E. & Wood, E. F. (1977) A state-space model for real-time forecasting of hydrological time series. J. Hydrol. Sci. 4 (1), 61-76.
Todini, E. (1978) Mutually interactive state-parameter (MISP) estimation. In: Proc. AGU Chapman Conference on Applications of Kalman Filtering Theory and Techniques to Hydrology, Hydraulics and Water Resources, 135-151. Dept of Civil Engineering, Univ. of Pittsburgh, Pittsburgh.
Todini, E. O'Connell, P. E. & Jones, D. A. (1980) Basic methodology: Kalman filter estimation problems. In: Real-Time Hydrological Forecasting and Control (ed. by P. E. O'Connell) ((Proc. 1st International Workshop, July 1977). Institute of Hydrology, Wallingford.
Todini, E. & Wallis, J. R. (1977) Using CLS for daily or longer period rainfall-runoff modelling. In: Mathematical Models in Surface Water Hydrology (ed. by T. A. Ciriani, U. Maione & J. R. Wallis), 149-168. Wiley, London.
Todini, E. S Wallis, J. R. (1978) A real-time rainfall-runoff model for an on-line flood warning system. In: Proc. AGU Conference on Applications of Kalman Filtering Theory to Hydrology,
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
Adaptive hydrological forecasting 205
Hydraulics and Water Resources (ed. by Chao-lin Chiu), 355-368. Dept of Civil Engineering, Univ. of Pittsburgh, Pittsburgh.
Tucci, C. E. M. & Clarke, R. T. (1980) Adapative forecasting with a conceptual rainfall runoff model. In: Hydrological Forecasting (Proc. Oxford Symp., April 1980), 445-454. IAHS Publ. no. 129.
Venn, M. W. & Day, B. (1977) Computer aided procedure for time series analysis and identification of noisy processes (CAPTAIN): user manual. Report no. 39, Institute of Hydrology, Wallingford.
Weiss, G. (1980) Basic methodology: the Kalman filter. In: Real-Time Hydrological Forecasting and Control (ed. by P. E. O'Connell) (Proc. 1st International Workshop, Wallingford, July 1977) . Institute of Hydrology, Wallingford.
Weiss, G. (1981) Prediction in the presence of non-Gaussian noise. In: Proc. 2nd International Workshop of Real-time Hydrological Forecasting and Control (ed. by P. E. O'Connell), (Wallingford, July 1979). Institute of Hydrology, Wallingford (in prep).
Whitehead, P. G. & Young, P. C. (1975) A dynamic stochastic model for water quality in part of the Bedford-Ouse River system. In: Computer Simulation of Water Resources Systems (ed. by G. C. Vansteenkiste), 417-438. North-Holland, Amsterdam.
Wittenmark, B. (1974) A self-tuning predictor. Inst. Elect. Electron. Engrs Trans. Automat. Control AC-19, 848-851.
Young, P. C. (1974) A recursive approach to time-series analysis. Bull. Inst. Maths Applic. 10, 209-224.
Young, P. C. & Jakeman, A. (1978) Refined instrumental variable methods of recursive time series analysis, Part 1: Single-input, single output systems. CRES, Report no. AS/R12, Australian National University, Canberra.
Young, P. C , Shellswell, S. H. & Neethling, C. G. (1971) A recursive approach to time-series analysis. Dept of Engng, Univ. of Cambridge, Report no. CUED/B-C0NTR0L/TRI6.
Young, P. C. & Whitehead, P. G. (1977) A recursive approach to time series analysis for multivariable systems. Int. J. Control 25 (3), 457-482.
Paper presented at the Symposium on Hydrological Forecasting, Oxford, April 1980; revised 2 October 1980.
Dow
nloa
ded
by [
Uni
vers
ity o
f N
ew H
amps
hire
] at
00:
01 0
9 O
ctob
er 2
014
