Journal of
Hydrology ELSEVIER
Journal of Hydrology 157 (1994) 177-195
[4]
Improvements to the modelling of British rainfall using a modified Random Parameter Bartlett-Lewis Rectangular Pulse Model Christian Onof*, Howard
S. W h e a t e r
Department of Civil Engineering, Imperial College of Science, Technology and Medicine, London SI,V7 2BU, UK (Received 26 August 1993; revision accepted 5 November 1993)
Abstract The Random Parameter Bartlett-Lewis Rectangular Pulse Model (RPBLRPM) was shown in a previous paper to provide generally good results when applied to a 38.5 year record of hourly rainfall from the Elmdon raingauge, Birmingham, UK: second-order properties of the depth and the interevent time distributions and the proportion of dry periods for different timescales were on the whole well reproduced by the optimized parameters. However, two main problems arose, namely an overestimation of the autocorrelations for daily data in particular, and relatively poor performance in reproducing the extreme values of the depth distribution. The first of these problems is here addressed by superposition of a jitter process on the rectangular pulse, thus providing a more realistic representation of the continuous time process. Applying the jitter to each rainfall cell appears more effective than applying it to the instantaneous rainfall depth and enables a more satisfactory reproduction of the autocorrelations. The use of a gamma variable to replace the exponential distribution of the cell intensity provides a good reproduction of the extreme rainfall depths at different time-scales. Moreover, the values of the parameters providing an optimal reproduction of the extreme values of the depth distribution are in most cases close to those obtained by optimizing on the interevent time features, which provides a method for estimation and ensures that the two objectives do not conflict.
1. Introduction S t o c h a s t i c m o d e l s o f p o i n t rainfall time-series present m a n y p o t e n t i a l a d v a n t a g e s over m o r e t r a d i t i o n a l a p p r o a c h e s to rainfall representation: in p a r t i c u l a r , they enable * Corresponding author. 0022-1694/94/$07.00 © 1994 - Elsevier Science B.V. All rights reserved SSD1 0022-1694(93)02429-2
178
C. OnoJ~ H.S. Wheater / Journal of Hydrology 157 (1994) 177 195
extensive time-series to be generated and first and higher order statistics to be calculated while requiring the input of a small number of parameters only. These may be used for such engineering applications as reservoir design (daily to monthly durations), flood studies (hourly to daily rainfalls) and design of sewerage systems (subhourly durations). Models which generate continuous series of rainfall are advantageous in that they readily allow for temporal disaggregation of rainfall data of short duration from longer duration statistics. Poisson-cluster models, which suppose that storms arrive according to a Poisson distribution within which cells arrive according to another Poisson process (BartlettLewis process) or have exponentially distributed arrivals (Neyman Scott process), have been examined by various authors (Waymire and Gupta, 1981) and in particular, models with rectangular pulses for each cell have been successfully applied by Rodriguez-Iturbe et al. (1987, 1988) and Entekhabi et al. (1989) to Denver Summer data, and by Cowpertwait (1991) to British data. In a previous paper (Onofand Wheater, 1993), the Random Parameter Bartlett-Lewis Rectangular Pulse Model (RPBLRPM) was applied to 38.5 years of hourly rainfall data from the Elmdon raingauge (Birmingham, UK). In this model, storms arrived according to a Poisson process of parameter A; for each storm, a parameter ~/was chosen from a gamma distribution of shape parameter a and scale parameter l/u; conditionally upon this choice, cells arrived according to a Poisson process of parameter ~r/and with each cell arrival a rectangular pulse was associated which had an exponentially distributed duration of parameter 7) and an exponentially distributed depth of mean #x; conditionally upon ~/, the storm had an exponentially distributed duration parameter ~brl. This was therefore a model with six parameters (A, ~, u, ~, ~, Ux). Apart from a very good reproduction of mean, variance and autocorrelations of the depth distribution at all time-scales, the model was particularly successful in reproducing the proportion of dry periods at different time-scales, thus providing an improvement on the BLRPM (non-randomized Bartlett-Lewis Rectangular Pulse Model). The optimization of the parameters provided a good reproduction of the interevent time distribution. However, two issues remain to be examined: first, an overestimation for some months of depth autocorrelations for time scales of 12 h and above, when the optimized parameters were used; and secondly an inadequate reproduction of extreme values for long return periods. These problems are addressed in this paper by examining two new types of model: the RPBLJM 1, RPBLJM2 models, in which a jitter process is added to the rainfall intensity, i.e. the rainfall profile acquires an irregular shape resulting from the addition of a high-frequency random perturbation and the RPBLGM model, in which a gamma model is used for rainfall depth. The notations used in this paper are the same as in Onof and Wheater (1993): yfhl represents the ith value in the time-series of h hourly rainfall depths and the characteristic variables are (for h = 1, 6, 12 and 24 h); the average rainfall depth at time-scale h h E[Y~h)];the standard deviation of the same - - S[Y~h)]; the autocorrelation of lag k of rainfall depths for time-scale h h
= c°v[r'!h)' ri+ ]
179
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177 195
the average duration o f an event in the hourly data - - mw, where an event is a wet period preceded and followed by at least one dry time-interval; the mean interevent time in the hourly data - - md; the standard deviation o f the same - - s0; the average number o f rainfall events per month - - mn; the proportion o f dry h hour time intervals - - p(h). The analytical estimates o f these variables are primed: m~i, s~, ! I mw, m, and p(h)'. The historical values of these variables are shown in Table 1.
2. Introduction of a jitter
The representation o f rainfall which is suggested in the types of model which have Table 1
Elmdon historical data synoptic table Hourly data
m',.
6 hourly data
(h)
J F M A M J J A S O N D
E[Y~ ')] (mm)
S(Y~ l)) (mm)
A(1)
A(2)
A(3)
S(Y~ 6)) (mm)
A(I)
A(2)
A(3)
0.0777 0.0701 0.0694 0.0634 0.0766 0.0757 0.0674 0.0963 0.0807 0.0730 0.0920 0.0811
0.3256 0.3070 0.2946 0.3095 0.3976 0.4297 0.5021 0.6126 0.4414 0.3798 0.4063 0.3453
0.6191 0.5800 0.5868 0.5355 0.4131 0.4643 0.3902 0.4324 0.5356 0.5689 0.6059 0.5849
0.4033 0.3741 0.3967 0.3224 0.2828 0.2847 0.1675 0.2206 0.3258 0.3512 0.3897 0.3650
0.2922 0.2859 0.3085 0.2669 0.2037 0.1960 0.1094 0.1518 0.224 0.2453 0.2540 0.2669
1.3858 1.2844 1.2479 1.2466 1.5071 1.6652 1.7876 2.2053 1.7878 1.5644 1.7027 1.4553
0.3679 0.3532 0.3856 0.3602 0.2913 0.2835 0.1845 0.2824 0.2515 0.3032 0.2782 0.3137
0.1044 0.1163 1.1551 0.0970 0.0738 0.0807 0.0549 0.1015 0.0561 0.0914 0.0838 0.1011
0.0363 0.0897 0.1076 0.0314 0.0514 0.0592 0.0382 0.1062 0.0360 0.0439 0.0756 0.0751
12 hourly
J F M A M J J A S O N D
24 hourly data
S(Y~ '2)) A(1) (mm)
S(Y~ 24)) A(1) (mm)
2.3257 2.1379 2.0387 2.0677 2.4935 2.6950 2.7240 3.4660 2.8077 2.6010 2.7364 2.3584
3.6246 3.3837 3.3271 3.2137 3.8620 4.0733 4.0853 5.3548 4.3279 4.0050 4.3548 3.6320
0.1947 0.2372 0.3152 0.1969 0.1521 0.1817 0.1638 0.2311 0.1787 0.1544 0.2091 0.2458
0.1015 0.2546 0.2388 0.1778 0.1075 0.2219 0.2311 0.1352 0.1507 0.1181 0.1210 0.2012
mid
s'~
(h)
(h)
20.12 21.69 22.43 23.97 22.76 26.30 29.80 24.97 27.01 27.26 20.50 20.80
34.61 42.71 50.06 48.74 46.68 52.11 61.42 49.04 55.89 53.82 36.85 39.15
m'n
31.16 25.76 29.50 25.97 28.21 24.26 23.21 26.51 24.13 24.79 29.85 30.77
3.05 3.05 2.97 2.83 2.68 2.62 2.42 2.66 2.65 2.85 3.01 2.89
p(h)' h=l
h=6
h=12
h=24
0.87 0.89 0.88 0.90 0.90 0.91 0.92 0.91 0.91 0.90 0.88 0.88
0.71 0.75 0.74 0.76 0.76 0.78 0.81 0.77 0.78 0.78 0.72 0.71
0.59 0.64 0.62 0.66 0.65 0.69 0.72 0.67 0.69 0.68 0.59 0.59
0.43 0.50 0.48 0.52 0.51 0.55 0.59 0.51 0.53 0.53 0.43 0.42
180
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
been examined so far is obviously highly simplified since it is clear that the precipitation generating phenomenon is complex at all time-scales. The hypothesis of a rectangular pulse is crude in its simplification of the continuous-time phenomenon. This may be improved by representing our absence of knowledge of the variations at very small time-scales (a few minutes) by a jitter process. The introduction of a jitter was suggested for the BLRPM model by RodriguezIturbe et al. (1987) referred to hereafter as RCI. Here, it is applied to the RPBLRPM model. Two types of jitter are considered depending on whether the jitter is added to the continuous-time process or to each cell and the models referred to as Random Parameter Bartlett-Lewis Jitter Models 1 and 2 (RPBLJM 1 and RPBLJM2).
2.1. First type of jitter We therefore first consider the RPBLJM1 for which the continuous-time rainfall depth becomes Y(t) = exp[Z(t)] Y(t) where Z(t) is a stationary Gaussian process of mean #z, variance a~ and autocovariance function cz(t). For a given t, the normal distribution Z(t) represents the sum of different influences which, according to Central Limit theorems (e.g. theorem of Liapunov) converges to a normal distribution. RCI suggest that we choose E{exp[Z(t)]} = 1, which means that 1
2
=0
(l)
With this condition E[Y(t)] = E[Y(t)] Since for a bivariate Gaussian distribution
E{exp[Z(u) + Z(u + T)]} = exp[cz(T)] for az 2 small the following are obtained var[]~(t)] : %2 + O.z2(//,y2 + O.y2) q.. O(o.z 4)
Cov[Y(/), ¢(t + T)] = ey(T) -+ Cz(T)[~y 2 + £y(T)] -~- O(6rz 4) where ~ry2 and Cy(T) are, respectively, the variance and autocovariance of lag-~- of the process { Y(t)}. Therefore, a decrease of the autocorrelations is obtained if and only if
Cy(T)
var[IY(/)]
cov[ (t),
+
Cz(r )
< °r--z-T Ozy <
Cy(T)
(2)
The effect of the jitter on the aggregated process at time-scale h may therefore be calculated E l L (h)] = h#y
(3a)
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
181
var[Yff )] = 2 [h(h - u)cov[l~(t), l?(t + u)]du j0 or
var[])i (h)] = var[Yi (h)] + 2 Ii(h - bl)Cz(U){Cy(U) + py2}du
(3b)
and
J
cov[:~(ht, /-/+ (~/] = -~ +h(h
-
Ivl)cov[Y(t), l)(t + v)]dv
or
cov[?(h/, ~/+~(h/] = cov[Y?/, y,+ (h/] +
I+2
_ (h - Ivl)c~(kh + v ) { c y ( k h + v) +
.yN}dv
(3c)
If it is assumed that cz(r) decreases much more rapidly than Cy(T), and thereby that Eq. (2) is fulfilled, then an approximation of the second-order properties will be (RCI) var[l?i(h)] ~ var[Yi(h)] + 2h{ay2 + #y2} [ ? Cz(U)du dO and
COV[]~i(h), ~..i+k(h)]
I
oc,
Cov[Yi(h), yi+k(h)] + h{ayZ + p)2} 0 c~(kh + v)dv
where the second term is negligible for k > 1. So the jitter reduces all autocorrelations by the same factor, except that of lag 1. To pursue this, a special form for the autocorrelation cz(r) has to be adopted
cz(r) = 10 -2 exp(-~r) where ~ is sufficiently large so that
Cz(t) << cy(t) This implies that crz2 = cz(0 ) = 10 2 is small. Thus, the second-order properties may be written E[ ri (h)] : hlAy
(4a)
var[~i(h/] ~ var[Yi(h)] + 2h{o.y2 + py2} 10-2
(4b)
cov[L (h), L+l (h)] ~ cov[Yi (h), Yi+l (h)] + h{a~,2 + py2} exp(-~h)10 -2 . ~-
(4c)
182
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
and cov[~ .(h), I~i+k(h)] ~ cov[Yi (h), Yi+k(h)] for k > 1
(4d)
and the condition is 10 -2 e x p ( - Q ) << cy(t) or
e x p ( - Q ) << 10 a -
~#x
j-1
1
2#x2 + ~
0_
#x2
/( ~
'~ '
u
(4e)
The other properties relative to the time-distribution of rainfall events are not affected by the introduction of the jitter.
2.2. Parameters and results for the first type of jitter As with the B L R P M and R P B L R P M models, parameter estimation is carried out by identifying characteristic features of the aggregated process and solving the equations thus obtained with a modification of the Powell hybrid method. Since a reduction in the autocorrelations for time-scales smaller than 12 h is to be avoided, the variances and lag-1 autocovariances for smaller time-scales are used among the estimation variables. Thus, the following set of variables is chosen (E[ ri 0)], var[ Yi (1)],coy[ yi(1)
Yi+l(1)],p(1), var[ yi(6)], p(24), cov[ yi(6), Yi+l(6)]}
Solving this set of equations is a difficult task, and in most cases, only pseudosolutions are found, e.g. such that the characteristic variables are reproduced with a relative error smaller than 5%. However, a reduction in the autocorrelations for time-scales greater than 12 is obtained for this method (Figs. l(a) and (b)). Because of the difficulties encountered in solving the system, parameter ( is now fixed, thus reducing the number of equations to six. These are now
{E[Yi(I)],var[Yi(l)],cov[Yi(1), Y (l)l"rla kP',J ),t'~ n124~l, c°v[Yi(6)'Yi+l(6)] var[Yi(6)] An optimal choice of the fixed value ff = fiT is then required, so as to obtain an improved reproduction of the 12 and 24 hourly A(1): this is achieved by varying the value of ff and taking that for which the sum of the squares of the deviations between analytical and historical values for the two autocorrelations is minimal. The values of these autocorrelations are shown in Figs. l(a) and (b). The only months where these autocorrelations are not satisfactorily reproduced are months such as February and June for which the 24 hourly lag-1 autocorrelation is larger than the 12 hourly one. This shows an improvement on R P B L R P M in the reproduction of
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177 195
"'~]
.... RPB~M~BUM~ ~]~
i
^
e.~q [
/ /
"~
k
s
//
.~
I'\
'<
---
l//k\ I, ~ ~,
~
/ ....
",
,'
1' /',.\
,'- ....
',
,'
/',, ",,:
k
(e equ,tio-- - op",,,,,ed)
Ktstorieal data RPBLJM1 a n a l y t i c a l (7 e q u a t i o n s )
- -
//~
.'
',,
\
,'
',
,,'
2
3
6
,' ,,'
',:,;;it,/
\
4
,'
.', ,,,',,
;
',,
,
/
k
e.ls 1
183
6
7
8
9
16
I 11
I 12
Months 0.:]0--
RPBIJM1 RPBLRPM Historical -- - P.,PBLJM1
a n a l y t i c a l (6 e q u a t i o n s - o p t i m ~ e d ) analytical data a n a l y t i c a l (7 e q u a t i o n s )
....
/
O
I
/
\\--/
~"
',./ .,,' /
,'/
,',v'/
0.16
0.10 1
I 2
I 3
I 4
J
K
I
6
J 7
J
R
J
9
J
10
11
I
12
Months
Fig. 1. (a) Analytical autocorrelations of lag-I of 12 hourly rainfall for RPBLJM 1. (b) Analytical autocorrelations of lag-I of 24 hourly rainfall for RPBLJMI.
184
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
the autocorrelations, although for some months (February, June and October) the deviations are still large. For these months, it appears that there is a tendency for storms to occur in the afternoon, thus reducing the 12 hourly A(1). This bias is not represented by the RPBLJM1 model. Finally, it is possible to attempt, by fixing both u and (, to provide an optimal representation of the variables " mw, ' md' and mn' (Onof and Wheater, 1993) as well as an improved reproduction of the 12 and 24 hourly A(1). This is particularly important in view of the fact that the optimization of fff leads to a poorer representation of these variables. The system now has five equations
Ely.(')1 c°v[ Yi(l)' Yi+l (l)] .I t , J,
var[Yi(i)]
c°v[ Yi(6)' Yi+l(6)] )
,t,~l),p(24),
var[Yi(6)]
It appears, however, that an optimization of u is not compatible with an improved reproduction of the autocorrelations. So this RPBLJM1 model is best used for the improved reproduction of 12 and 24 hourly A(1) alone. The performance of the RPBLJM1 model in reproducing the other variables of Table 1 is similar to that of the RPBLRPM model in its non-optimized form (Onof and Wheater, 1993).
2.3. Second type of jitter The conflicting aims of reducing the 12 and 24 hourly autocorrelations suggest that another approach to the introduction of the jitter, i.e. by adding it to each rainfall cell, might be more germane to this type of model. Thus, the RPBLJM2 model is one in which for each cell, the depth Xt(u) at a time u after the cell arrival at time t becomes
Xt(u) = X exp[Zt(u)] where {Zt(u)} is a Gaussian process with the properties of {Z(u)} used above. Again, it is assumed that E(exp[Zt(u)]} = 1, which means that #z + l/2a~ 2 = 0. This yields (for ~z~ small) a-
1 (1 + t~/~b)
var[l?(t)] = ~AlZcU E(X2)(1 + a.2). + ~T-~#x
cov[l?(t), Y(t + 7-)] -- A#cU
+ O(crz4)
E(X2)[1 + c.(7-)] +
#x 2
+
which means that the variance of ]~(t) is deduced from that of Y(t) by adding A#cu Crz:E(X 2)
C. Onof, H.S. Wheater/ Journalof Hydrology157 (1994) 177-195
185
and the covariance of lag 7- of l>(t) from that of Y(t) by adding the term
~--i cz(7-)\ - - ~ ,
e(x 2)
The second-order properties of the aggregated process may be deduced from the knowledge of cov[I?(t), I¢(t + 7-)] (with the assumption that a~2 is small) ElY/(h)] = E[ri (h)] = hA#x
u
o~--1
(1 + ~/0)
var[L(h/] ~ var[L.(h/l + 2 Amy E(X 2) - 1 COV[ ])i (h) ,
(5a)
Cz(U) \TT-@
(h - u)du
(5b)
Yi+k(h)] ~ cov[Yi (h) , Yi+k (h)]
i
Am~ E(X 2) cz(kh + u) ~ + ~ + u+c~- 1 -h
(h-
lul)du
(5c)
To pursue this, as above a particular form for the covariance function cz(T) must be chosen. The aim of the choice must be to represent the fact that the jitter process has an autocorrelation function which decreases very rapidly with the lag, so that the autocorrelation of {/~i(h)} is smaller than that of { Yi(h)}. Choosing Cz(7-) = 10-2 exp(-~PT-) ( u - - ~ r ) '-~ which has the advantage of providing a simplification of the integrals above - - an increase is obtained in the autocorrelations. Therefore, the following choice is made
Cz(7-) = v exp(-~bT-) where v is small.
2.4. Parameters and results for the second type of.jitter Parameters are estimated with the same method as previously, using the set of equations
{
E[YiO)],var[Yi(')],cov[Yi (1), Yi+l
(')1 ~ll' c°v[Yi(6)'yi+'(6)] J,et J,
var[~.t0~'-~'] ,p(24)
}
with ~ v ~f'10-2, For for January, the following parameters are obtained: A = 0.0297; #x = 0.7306; o~ = 3.0748; u = 1.1822; ~ = 0.6633; 0 = 0.1254; ~ = 1.0000. These yield 12 hourly A(1)=0.2771 and 24 hourly A(1)=0.1946, which represents no improvement on the results obtained with the RPBLRPM parameter sets.
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
186
A larger jitter covariance function is therefore required. Bearing in mind the condition that az 2 should be small, v = 10 -l may be chosen. Parameters are estimated for values of~pf ranging from 0.01 to 20 s 1. Optimal values of this parameter are chosen to yield the best reproduction of the 12 and 24 hourly A(1). This optimal reproduction is on the whole superior to that of the RPBLJM1 model. The analytical values of the characteristic variables are shown in Table 2. These parameters are also more effective than those of the RPBLJM1 model in reproducing the features of the time distribution of rainfall events: m~w,m/n and m~ because the optimization of ~f does not imply a poorer performance of those variables. Generally, this model performs as the non-optimized RPBLRPM model in reproducing the other characteristic variables. On the whole, therefore, the RPBLJM2 model is superior to the RPBLJM1 model in that it provides an improved reproduction of the 12 and 24 hourly A(1) (Figs. 2(a) and (b)). Table 2 RPBLJM2 analytical synoptic table Hourly data
6 hourly data
m',,
(h)
J F M A M J J A S O N D
E[Yi (1)]
S(Y~ '))
(mm)
(mm)
0.0777 0.0701 0.0694 0.0634 0.0766 0.0757 0.0609 0.0963 0.0807 0.0730 0.0920 0.0811
0.3256 0.3070 0.2946 0.3095 0.3976 0.4297 0.4807 0.6126 0.4414 0.3798 0.4063 0.3453
12 hourly
J F M A M J J A S O N D
A(l)
A(2)
A(3)
S(Y~ 6)) (ram)
A(1)
A(2)
A(3)
0.6192 0.5823 0.5860 0.5555 0.4132 0.4704 0.4133 0.4321 0.5356 0.5684 0.6050 0.5859
0.4039 0.3691 0.3758 0.3457 0.2455 0.2307 0.1805 0.2292 0.2582 0.3211 0.3348 0.3460
0.2971 0.2772 0.2881 0.2632 0.1916 0.1634 0.1163 0.1617 0.1636 0.2136 0.2205 0.2460
1.3928 1.2823 1.2443 1.2594 1.4683 1.6129 1.7018 2.2465 1.7081 1.5261 1.6683 1.4258
0.3677 0.3513 0.3849 0.3607 0.2903 0.2819 0.2161 0.2815 0.2503 0.2938 0.2782 0.3127
0.1379 0.1218 0.1637 0.1392 0.0870 0.1169 0.0990 0.1413 0.0839 0.1201 0.0768 0.0980
0.0801 0.0578 0.0925 0.0793 0.0407 0.0652 0.0801 0.1118 0.0404 0.0750 0.0318 0.0453
3.9691 4.3564 4.4332 4.0263 3.3299 3.0339 4.0009 2.6919 2.9608 3.3422 3.3649 3.7626
h=24
24 hourly
S(Y~ 12))
A(I)
S(Y~ 24)) A(I)
2.2747 2.1034 2.0692 2.0549 2.3309 2.5798 2.6497 3.5279 2.6985 2.4209 2.6645 2.3088
0.2494 0.2386 0.2896 0.2516 0.1802 0.2249 0.2006 0.2389 0.1819 0.2151 0.1800 0.2100
3.5705 3.3061 3.3216 3.2329 3.5562 4.0344 4.1004 5.4858 4.1456 3.7340 4.0908 3.5898
0.1615 0.1421 0.1985 0.1642 0.0979 0.1553 0.2227 0.2176 0.1117 0.1402 0.1009 0.1262
m'd
S'd
(h)
(h)
24.55 33.94 31.62 34.60 29.66 30.46 55.03 27.26 29.20 29.06 23.11 26.19
31.88 39.63 37.79 41.54 39.22 44.44 263.87 43.64 41.52 42.56 31.02 30.80
m',
26.09 19.43 20.64 19.26 22.55 22.21 12.60 24.84 23.13 22.96 28.10 24.84
p(h)' h= 1
h=6
h
0.86 0.89 0.88 0.90 0.90 0.91 0.93 0.91 0.91 0.90 0.87 0.87
0.73 0.78 0.76 0.79 0.79 0.80 0.77 0.75 0.79 0.77 0.74 0.74
0.61 0.67 0.65 0.69 0.68 0.70 0.66 0.65 0.69 0.68 0.61 0.61
12
0.43 0.50 0.48 0.52 0.51 0.55 0.57 0.51 0.53 0.53 0.43 0.42
C. Onof, H.S. Wheater I Journal of Hydrology 157 (1994) 177-195
"=1 I
^
"-:1 / .,I
t/%,
b, /",
e.~=-I
."
L
RPBIJi~ ~ e a l RPBLRPM a n a l y t i c a l
....
/,,:~ 1' 6/
"-'~
___-::"<"°""~_
\\ \';\
I '7 ,'
•
,' ,'
-
,,,
,,
",
i
',
/-
\ \
I
,,
', \
I
,,
',
k x
..:l-~
.... •
\'~
I/
\ ,'
..-]
2
3
//'\x
/,
~
I\
I
6
G
7
8
9
18
I'
I
I
11
12
Mouths
RPBI~I~ ~ e a l
....
,v,
I "..=,-I t .,
4
,' ,'
"
I
\1 I
,,
e.16
1
187
RPBLRPM l a l ~ c a l Hlatorioel d a t a
"~
,",\
i ,',\
/ ,'l
I1\',\
I,l
,',
I,
\\\
II
\
,, I
'
\, \
#
i7
",
",,,,'1
~
\/I
4/
1.1E:
/\ \ \1
t.
e'~e/
I
1
2
I 3
I 4
E
I
6
I
7
I
El
I
9
I
le
I
11
I
12
Mont.l~
Fig. 2. (a) Analytical autocorrelations of lag-1 of 12 hourly rainfall for RPBLJM2. (b) Analytical autocorrelations of lag-1 of 24 hourly rainfall for RPBLJM2.
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
188
3. Use of a Gamma depth distribution
The poor performance of the RPBLRPM model in the reproduction of extreme values suggests that another distribution be chosen for the cell depth. A twoparameter distribution would have the advantage of introducing a greater flexibility in the reproduction of extreme values. The Gamma distribution which was first used by Thom (1958) for rainfall depths, has been successfully applied by Acreman (1990) for Farnborough hourly precipitation and by Eagleson (1978) for Boston hourly rainfall. This latter author underlines the advantage of this choice in that the Gamma distribution: (1) is 'self-preservative' since the sum of n Gamma variables I'(a, bi) (a, scale parameter; b, shape parameter) is a i=n
F(a, Z b i ) i=1
distribution; (2) is simple enough for analytical tractability; (3) provides a good fit for rainfall depth observations at different time-scales and for different climatic types. We therefore choose as density function for X
fx(t)
6PtP-! exp(-6t), = F(p)
for t > 0
= 0 , for t = 0 This model is referred to as the RPBLGM model.
3.1. Model equations The continuous-time model equations therefore become
E[Y(t)]
//
= A/zx(1 + ~/~b) a---Z-~
var[Y(t)]
~52(-a-- 1)
cy(~-) = 62(~ _ 1)
(6a) t~
2
(6b)
42
\~--U~} (6c)
Similarly, the first- and second-order properties of the aggregated process are easily derived v ~[r~ (h/] = ~h~x(1 + ,~/~) ,T~- I
(7a)
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
189
and with the expressions al =
1)--~----(OeA"cV~ 2) - 3) [ p0, + 1) + - ~p"
62(a
and A2 =
A#cn#x2U'~ ~2(~b2 - 1)(ct - 1)(c~ - 2)(c~ - 3)
the following are obtained var[Yi/h}] = 2Al[(c~ - 3)hu 2-'~ - u 3-`~ + (u + h) 3 ~] - 2A2[qS(c~ - 3)hu 2-a - u 3-'~ + (u + ~bh)3-a]
(7b)
and cov[Yi (h), Yi+k (h)] = Al{[u + (k + 1)hi 3-a - 2(u + kh) 3-~ + [u + (k - 1)hi 3-a }
- A2{[u + (k + 1)q~h]3-'~ - 2(u + k(gh) 3-a + [u + (k - 1)@h]3-a}
(7c)
The properties of the time-distribution of rainfall events are not affected by the change of rainfall cell depth distribution. 3.2. P a r a m e t e r s and analytical results
The R P B L G M model is a seven-parameter model and estimating all these parameters is a difficult task, as was observed for the R P B L J M models. The parameter 6 must therefore be fixed. This is found to be preferable to the fixing of parameter p because it provides a greater range of parameter values with few iterations (Onof, 1992). The parameters are estimated using the same methods as previously for the identification of the following set of characteristic variables
var [yi(1)], cov[ yi(l), Yi+l (l)],P(1), var[ Yi(6)],p(24)}
{E[Yi(1)],
If the month of January is examined, the parameter sets obtained for different values of t5 = 6f are such that, as 6f increases, so does
E(X)
=P_ 6/
Consequently, the cell arrival rate E(3) =
OOg -
-
V
decreases and so does the average storm duration
L/
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177 195
190
The average number of cells per storm, lzc, reaches a minimum for a finite value of 6f as Fig. 3 shows. The analytical results enable the estimation of the parameter set which performs optimally for the objective function
2/[md(H) - m~(i)] 2 d;= v -~d~it-~
[mn(H) - m~n(i)]2
mn(n)2
+
where rnd(t ' "), mn(, ' ') are the analytical (iteration i) values of the two variables rn~, mn, and md(H) and mn(H) their historical values. The optimal set of parameters is obtained for #c approximately minimal: 6f = 1.5 and p = 1.77 (Fig. 3). This is practically equivalent to the optimal parameter set obtained for the Bartlett-Lewis RPBLRPM optimization (Onof and Wheater, 1993) with #x(RPBLRPM) ~ P (RPBLGM) To obtain this value, a small mesh (10 -2) is required for the incrementation of 6f. The optimization, when carried out for other months of the year, yields similar results. The analytical values of the characteristic variables are in Table 3. This optimization (which we shall refer to as O1) may of course be compared with that of u = uf for the (/z.,d'l)
40~
\
=
f(6)
--1.0
./UZXUL"7 R P B X . I ~ M (L fixed) \
36--
32--
28--
\ \
tO.8
\
\ \
24--
\
-0.6
\ \
20--
I \
i --e.4
16--
\
\
II
12--
--8.2
8-
4
8
1.0
I
I
I
1.2
1.4
1.6
Q.O 1.8
L (""-t.h) Fig. 3. N u m b e r of cells per storm and objective function di as a function of 6f for R P B L G M (January).
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
191
Table 3 B L R G analytical synoptic table
Hourly data
mwt
6 hourly data
(h)
J F M A M J J A S O N D
E[Y~'1]
S(Y~ '))
(mm)
(mm)
0.0777 0.0701 0.0766 0.0634 0.0766 0.0757 0.0674 0.0963 0.0807 0.0730 0.0920 0.0811
0.3256 0.3070 0.3976 0.3095 0.3976 0.4297 0.5021 0.6126 0.4414 0.3798 0.4063 0.3453
A(2)
A(3)
S(Y{ 6)) (ram)
A(1)
A(2)
A(3)
0.6197 0.5804 0.4131 0.5355 0.4131 0.4641 0.3903 0.4325 0.5353 0.5692 0.6058 0.5846
0.3759 0.3529 0.2525 0.3058 0.2526 0.2681 0.1802 0.1812 0.3148 0.3333 0.3555 0.3639
0.2757 0.2677 0.2026 0.2268 0.2026 0.2061 0.1262 0.1153 0.2368 0.2457 0.2546 0.2792
.3858 .2844 .5071 .2466 .5071 .6652 .7876 2.2053 1.7878 1.5644 1.7027 1.4553
0.3802 0.4032 0.3850 0.3654 0.3844 0.3717 0.2603 0.2274 0.3750 0.3701 0.3531 0.4140
0.1907 0.2281 0.2272 0.1988 0.2259 0.2194 0.1292 0.0997 0.2020 0.1936 0.1635 0.2355
0.1317 0.1685 0.1667 0.1439 0.1651 0.1669 0.0898 0.0654 0.1451 0.1377 0.1075 0.1751
2.7701 2.4736 3.0999 2.3604 2.5544 2.8875 2.2310 2.1389 2.3384 2.5404 2.7014 2.3050
h= 1
h=6
h = 12
h=24
0.87 0.89 0.88 0.90 0.90 0.91 0.92 0.91 0.91 0.90 0.88 0.88
0.70 0.74 0.73 0.76 0.77 0.77 0.80 0.76 0.78 0.76 0.72 0.70
0.59 0.64 0.62 0.66 0.67 0.68 0.72 0.66 0.68 0.67 0.60 0.58
0.43 0.50 0.48 0.52 0.51 0.55 0.59 0.51 0.53 0.53 0.43 0.42
24 hourly data
12 hourly
J F M A M J J A S O N D
A(I)
S(Y~12)) A(I)
S(Y~24)) A(1)
(mm)
(ram)
2.3024 2.1516 2.1105 2.0600 2.5078 2.7581 2.8381 3.4552 2.9619 2.5896 2.8010 2.4473
0.3236 0.3662 0.3947 0.3321 0.3616 0.3563 0.2414 0.2005 0.3285 0.3266 0.2910 0.3749
3.7460 3.5567 3.5250 3.3624 4.1384 4.5424 4.4720 5.5359 4.8279 4.2182 4.5009 4.0583
0.2708 0.3267 0.3564 0.2941 0.3231 0.3351 0.2145 0.1669 0.2827 0.2840 0.2336 0.3374
m~
s~
(h)
(h)
22.19 26.23 24.03 27.52 28.42 29.23 33.14 26.66 29.53 28.83 22.53 21.60
33.61 42.64 39.83 44.22 41.04 49.32 53.20 40.77 43.21 45.54 32.47 33.13
m~n
35.04 33.08 27.24 31.53 29.11 23.17 26.65 31.32 26.99 29.35 33.31 38.68
p(h)'
R P B L R P M model (Onof and Wheater, 1993) which also led us to consider the values of the fixed parameter for which #c is minimal. 3.3. Simulated results." the extreme values
The main purpose of the use of a gamma distribution is the improved reproduction of extreme rainfall depths. Therefore, 240 years of data were simulated for the month of January for different values of 6f. We may seek the value of 6 = ~Sfwhich yields the best reproduction of the extreme hourly and daily depths. This will be referred to as optimization 02. The choice of ~Sf is made by considering the best approximation of the mean values 17"max(I,1) and Ymax(4, 1) of the simulated distribution of annual extremes of hourly and daily data, respectively, for the month of January. We obtain
192
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
//
2
I
-1.8
I
e.e
I
l,e
I
2.0
1
3.0
I
4.0
I
s.e
6.e
Gumbcl r e d u c e d variate Return
i period (yrs) S
i
i
i
5
10
BO
i
I
100
200
IIS4SS~J~SIsIS
40--
EJIA~Jmd14sta 18--
e -1 . e
I
I
e.O
1 ,e
I
I
2,0
3,0
I 4.0
I
I
E.e
6.0
Gumbel r e d u c e d variate Return
i period (yre) 2
i 5
i 10
i 20
i 100
i 200
Fig. 4. (a) January hourly extreme values for RPBLGM. (b) January daily extreme values for RPBLGM. (c) June hourly extreme values for RPBLGM. (d) June daily extreme values for RPBLGM.
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
193
10-
14--
~
s s s j s S~SpS rS
!2-10-
o
js
•
~ ~ d a t a
•
I 0.0
-1.0
I 1 .O
I 2.0
I 3.0
I 4.0
q
I 5.0
6.0
Gu.mbel r e d u c e d varlate Return
period
I (yrs) Z
1 5
i 10
i 20
i 100
i 200
pJ~
s SS
i e -~.e
~s ~.s
S
.
.
.
.
sJPs~"
.
6 opt~.T
I
I
I
I
I
I
e.o
l.e
2,o
3.0
4.e
s.o
I
6.0
(;umbel reduced vartate Return
period
(yrs)
! ' Z
i 5
i tO
Fig. 4. (continued).
i 20
i tO0
data
i 200
194
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
a value of 6f approximately equal to 1.35. The optimal distributions of maxima obtained by simulation also reproduce accurately the skewness of the historical maxima (1.292 for the simulated hourly and 1.311 for the historical data). These results obtained for January may be generalized to the other months of the year. Figs. 4 (a)-(d) show the reproduction of hourly and daily extreme values for the months of January and July. For all months of the year, it is possible to obtain a very good reproduction of the extreme rainfall depths at these two timescales. The only remaining problem is as to the possibility of determining analytically the values of 6f which provide these very good approximations. So far they have been determined by simulation. Fig. 5, however, shows a comparison of the values of 6 = 6f which are optimal for O1 and for 02. For most months of the year, these optima are very close so that in most cases, the approximation
, r(ol )
6/(02)
(8)
may be used to estimate the parameters which will provide a good reproduction of the extreme values at both the hourly and the daily time-scales. These models provide a similar reproduction of the other characteristic variables to that of the optimized version of the RPBLRPM model (Onof and Wheater, 1993). 3,e-
/•
2, S
/
Optimalfor ~ propertle. / Optimalfor extreme-vaAue 1
L\~A
2,0-
~ 1,E;-
1.(~-
0.5-
0,0
1
2
I
3
1
I
I
I
I
I
4
6
6
7
B
9
Mot't~.s
Fig.5.Optimalvaluesof6t forRPBLGM.
I
10
I
11
I
12
C. Onof, H.S. Wheater / Journal of Hydrology 157 (1994) 177-195
195
4. Conclusions The refinements of the RPBLRPM model which have been presented above allow for improvements in the reproduction of the autocorrelations at time-scales larger than 12 h and of the extreme values of the rainfall depths. In this way, models are obtained which satisfactorily reproduce all the characteristic features of the rainfall process at different time-scales. The number of parameters has, in both the RPBLJM and RPBLGM models, been increased to seven. From the difficulties encountered in estimating parameter sets to fit these models, it appears that it is not advisable to seek models with more parameters than these for each month of the year. Rather, the success of these models suggests that further research should focus on more widespread application, including, for example, a method of regionalization of the parameters to enable application to areas where short-duration data are sparse or non-existent.
5. References Acreman, M.C., 1990. A simple stochastic model of hourly rainfall for Farnborough. Hydrol. Sci. J., 35(2): 119-148.
Cowpertwait, P.S.P., 1991. Further developments of the Neyman Scott clustered point process for modelling rainfall. Water Resour. Res., 27(7): 1431-1438. Eagleson, P.S., 1978. Climate, soil and vegetation, 2. The distribution of annual precipitation derived from observed storm sequences. Water Resour. Res., 14(5): 713--721. Entekhabi, D., Rodriguez-Iturbe, I. and Eagleson, P.S., 1989. Probabilistic representation of the temporal rainfall process by a modified Neyman Scott rectangular pulses model: parameter estimation and validation. Water Resour. Res., 25(2): 295 302. Onof, C., 1992. Stochastic modelling of British rainfall using Poisson processes. Ph.D. thesis, Imperial College, London. Onof, C. and Wheater, H.S., 1993. Modelling of British rainfall using a Random Parameter Bartlett-Lewis Rectangular Pulse Model. J. Hydrol., 149: 67-95. Rodriguez-Iturbe, I., Cox, D.R. and Isham, V., 1987. Some models for rainfall based on stochastic point processes. Proc. R. Soc. London, Ser. A, 410: 269-288. Rodriguez-Iturbe, I., Cox, D.R. and Isham, V., 1988. A point process model for rainfall: further developments. Proc. R. Soc. London, Ser. A, 417: 283-298. Thom, H.G., 1958. A note on the Gamma distribution. Month. Weather Rev., 86(4): 117 122. Waymire, E. and Gupta, V.K., 1981. The mathematical structure of rainfall representation, 3. Some applications of the point process theory to rainfall processes. Water Resour. Res., 17(5): 1287-1294.