Modeling All Exceedances Above a Threshold Using an Extremal Dependence Structure: Inferences on Several Flood Characteristics

Mo deling All Exceedances Ab o v e a Threshold Using an 1 Extremal Dep endence Structure: 2 Inferences on Sev eral Flo o d Characteris tics 3 Mathieu Ribatet ∗ , † T aha B.M.J. Ouarda † Eric Sauquet ∗ 4 Jean-Mic hel Gr´ esillon ∗ 5 Submitted to: Water R esour c es R ese ar ch 6 ∗ Cemagref Lyon, Unit ´ e de Recherch e Hydro logie-Hydraulique, 3 bis quai Chauveau, CP2 20, 7 69336 Lyon cedex 0 9, F r ance 8 † INRS-ETE, Universit y of Qu´ ebe c, 490, de la Couro nne Qu´ eb ec, Qc , G1 K 9A9, CANADA . 9 Corresp o nding author: M. Ribatet; Email: ribatet@lyon.cemagref.fr 10 Phone: +33 4 72 20 87 64; F ax: +33 4 78 47 78 75 11 Abstract 12 Floo d quantile estimation is of gr eat imp ortance for many engineering studies and policy deci- 13 sions. Ho wev er, p ractitioners must often deal with small data av ailable. Thus, the i nformation 14 must be used op t imally . In the last decades, t o redu ce the waste of data, inferen tial metho d- 15 ology has evolv ed from annual maxima mo deling to p eaks ov er a threshold one. T o mitigate 16 the lack of data, p eaks ov er a threshold are sometimes combined with add itional information 17 - mostly regional and historical information. H o wev er, whateve r the extra information is, the 18 most precious information for th e practitioner is found at the t arget site. In th is study , a mo del 19 that allows inferences on the whole time series is introduced. In particular, the prop osed mo del 20 takes into account the d ep endence b etw een successive extreme observ ations u sing an appropri- 21 ate ex tremal d ep endence structure. R esults s how t h at this mo del leads to more accurate flo od 22 p eak qu antil e estimates th an conv entio nal estimators. In ad d ition, as t h e time dep enden ce is 23 taken into account, inferences on other floo d c haracteristics can b e p erformed. An illustration 24 is given on fl oo d duration. Our analysis show s that the accuracy of the prop osed mo dels to 25 estimate the flo o d duration is related to sp ecific catc hment characteristics. Some suggestions 26 to increase the floo d duration predictions are introduced. 27 1 1 In tro duction 28 Estimation o f extreme flo o d event s is an imp ortant stage for man y enginee r ing designs a nd ris k 29 management. This is a consider a ble task as the amount of data av ailable is often small. Thus, 30 to incr ease the pr ecision a nd the qualit y of the estimates, several authors use extra information 31 in addition to the targ et s ite one. F or ex a mple, Riba tet et al. [2007a], Kjeldsen and Jones [200 6, 32 2007] a nd Cunderlik and O uarda [2006 ] add information from other homogeneo us gaging stations. 33 W err itt y et al. [2 006] and Reis Jr. and Steding er [2005] us e his to rical information to improv e 34 inferences. Incorp ora tion of extra information in the estimation pro cedur e is a ttractive but it 35 should no t be more prominent tha n the original data [Ribatet et a l., 2 007b]. Before look ing at 36 other kinds o f information, it s eems reasonable to use efficien tly the o ne a v aila ble at the targe t site. 37 Most often, practitioners have initially the who le time series, not only the extreme observ ations. In 38 particular, it is a considerable waste of information to r educe a time s eries to a sample of Annu al 39 Maxima ( AM ). 40 In this p ersp ective, the Peaks Over Threshold ( POT ) appro ach is less wasteful as mo re than one 41 even t p er y ear co uld b e inferred. How ever, the declustering method used to identify indep endent 42 even ts is quite sub jective. F urthermore, even tho ug h a “qua si automa tic” pro ce dure was recently 43 int ro duced by F erro and Seger s [200 3], there is still a waste of informa tion as only cluster maxima 44 are used. 45 Coles et al. [1994] and Smith et al. [1997 ] pro p ose an approach using Mar ko v c hain models 46 that uses a ll exceeda nc e s a nd accounts for tempo r al dep endence betw een succes sive observ ations. 47 Finally , the ent ire information a v ailable within the time series is taken into accoun t. Mo re recen tly , 48 F aw cett and W alshaw [2006] give an illustrativ e applica tion of the Ma rko v c hain mo del to ex tr eme 49 wind s pe ed mo de ling . 50 In this study , extreme flo o d ev ents are of interest. The perfo rmance of the Marko v chain model 51 is compared to the conv entional POT appro ach. The data analyzed consist of a collection of 50 52 F rench gaging statio ns. The a rea under study r anges fr om 2 ◦ W to 7 ◦ E and from 45 ◦ N to 51 ◦ N. The 53 drainage areas v ar y from 72 to 3 8300 km 2 with a median v a lue of 792 km 2 . Daily observ ations were 54 2 recorded from 39 to 105 years, with a mean v alue of 60 y ears. F or the remainder of this article, the 55 quantile benchmark v a lues are der ived from the maximum lik eliho o d estimates on the whole times 56 series using a conv entional POT ana lysis. 57 The pap er is or ganized as follows. Section 2 intro duces the theoretical asp ects for the Ma r ko v 58 chain mo del, while Section 3 chec ks the relev ance of the Marko vian mo del hypothesis. Section 4 59 and 5 analyze the pe rformance of the Marko vian mo de l to es timate the flo o d p eaks and durations 60 resp ectively . Finally , some conclusio ns and p ersp ectives are dr awn in Sec tio n 6. 61 2 A Mark o v Chain Mo del for Cluster Exceedances 62 In this section, the extremal Mar ko v chain mo del is pre sented. In the remainder of this a rticle, it 63 is assumed that the flow Y t at time t dep ends on the v alue Y t − 1 at time t − 1. The dependence 64 betw een tw o consecutive obser v ations is modeled by a first order Mark ov chain. Before introducing 65 the theoretical asp e c ts o f the mo del, it is w orth justifying and descr ibing the main adv antages of 66 the pr op osed approach. 67 It is now well-known that the univ ariate Extr eme V alue Theory ( EVT ) is relev a nt when mo d- 68 eling either AM or POT. Nevertheless, its ex tension to the mult iv a riate case is s ur prisingly rare ly 69 applied in practice. This work aims to motiv a te the us e of the Multiv ariate EVT ( MEVT ). In 70 our applicatio n, the m ultiv ar iate results a re used to mo del the dep endence b etw een a set of lagged 71 v alues in a times s e ries. Consequently , compared to the AM or the POT a pproaches, the amount 72 of observ a tions used in the inference pro cedure is clearly larger. F or instance, while only cluster 73 maxima a re used in a POT analysis, a ll exceedances a r e inferred using a Marko vian mo del. 74 2.1 Lik eliho o d function 75 Let Y 1 , . . . , Y n be a s tationary first- o rder Marko v c hain with a joint distribution function of t wo 76 consecutive observ ations F ( y 1 , y 2 ), and F ( y ) its marginal distribution. Th us, the likeliho o d function 77 L ev aluated at po int s ( y 1 , . . . , y n ) is : 78 3 L ( y 1 , . . . , y n ) = f ( y 1 ) n Y i =2 f ( y i | y i − 1 ) = Q n i =2 f ( y i , y i − 1 ) Q n − 1 i =2 f ( y i ) (1) where f ( y i ) is the marginal density , f ( y i | y i − 1 ) is the conditional density , and f ( y i , y i − 1 ) is the join t 79 density of tw o consecutive observ a tio ns. 80 T o mo del all exceedances a b ov e a sufficiently lar ge threshold u , the joint and mar ginal densities 81 m ust b e known. Standard univ a riate EVT arg uments [Coles, 2001] justify the use of a Genera lized 82 Pareto Distribution ( GPD ) for f ( y i ) - e.g. a term of the denominator in equation (1). As a 83 consequence, the marginal distr ibution is defined by: 84 F ( y ) = 1 − λ  1 + ξ y − u σ  − 1 /ξ + , y ≥ u (2) where x + = max(0 , x ), λ = Pr[ Y ≥ u ], σ and ξ are the scale a nd sha pe parameter s resp ectively . 85 Similarly , MEVT ar g uments [Res nick, 1987] argue for a biv a riate extreme v a lue distr ibution for 86 f ( y i , y i − 1 ) - e.g. a ter m o f the numerator in equation (1). Thus, the joint distr ibution is defined 87 by: 88 F ( y 1 , y 2 ) = exp [ − V ( z 1 , z 2 )] , y 1 ≥ u , y 2 ≥ u (3) where V is a homog eneous function of order -1 , e.g. V ( nz 1 , nz 2 ) = n − 1 V ( z 1 , z 2 ), sa tisfying 89 V ( z 1 , ∞ ) = z − 1 1 and V ( ∞ , z 2 ) = z − 1 2 , and z i = − 1 / log F ( y i ), i = 1 , 2. 90 Contrary to the univ aria te case, ther e is no finite parametriza tion fo r the V functions. Thus, it is 91 common to use sp ecific parametric families for V such as the log istic [Gumbel, 1 960], the asymmetric 92 logistic [T a wn, 19 88], the negative logistic [Gala m b os, 19 75] or the asy mmetric neg ative logistic [Jo e, 93 1990] mo dels. Some details for these parametris ations ar e rep o rted in Annex A. These mo dels , as 94 all mo dels o f the form (3) a re asymptotically dep endent , that is [Coles et al., 199 9] 95 4 χ = lim ω → 1 χ ( ω ) = lim ω → 1 Pr [ F ( Y 2 ) > ω | F ( Y 1 ) > ω ] > 0 (4) χ = lim ω → 1 χ ( ω ) = lim ω → 1 2 log (1 − ω ) log Pr[ F ( Y 1 ) >ω ,F ( Y 2 ) >ω ] − 1 = 1 (5) Other para metric families exist to consider simultaneously asymptotically de p endent a nd inde- 96 pendent cases [Bo r tot and T awn, 1998 ]. How ever, apart from a few particular cases (see Sectio n 3), 97 the d ata analyzed here seem to b elong to the asymptotically dependent class. Consequently , in this 98 work, only asy mptotically dep endent mo dels a re consider ed - i.e. of the form (1)– (3 ). 99 2.2 Inference 100 The Ma rko v c hain mo del is fitted using maximum censored lik eliho o d estimation [Le dfo r d and 101 T awn, 1996]. The contribution L n ( y 1 , y 2 ) of a p oint ( y 1 , y 2 ) to the numerator of equation (1) is 102 given b y: 103 L n ( y 1 , y 2 ) =                        exp [ − V ( z 1 , z 2 )] [ V 1 ( z 1 , z 2 ) V 2 ( z 1 , z 2 ) − V 12 ( z 1 , z 2 )] K 1 K 2 , if y 1 > u , y 2 > u exp [ − V ( z 1 , z 2 )] V 1 ( z 1 , z 2 ) K 1 , if y 1 > u , y 2 ≤ u exp [ − V ( z 1 , z 2 )] V 2 ( z 1 , z 2 ) K 2 , if y 1 ≤ u , y 2 > u exp [ − V ( z 1 , z 2 )] , if y 1 ≤ u , y 2 ≤ u (6) where K j = − λ j σ − 1 t 1+ ξ j z 2 j exp(1 /z j ), t j = [1 + ξ ( y j − u ) / σ ] − 1 /ξ + and V j , V 12 are the partia l der iv ative 104 with resp e ct to the co mpo nent j and the mixed partial deriv ative resp ectively . The contribution 105 L d ( y j ) o f a p oint y j to the denominator of equatio n (1) is given b y: 106 L d ( y j ) =        σ − 1 λ [1 + ξ ( y j − u ) /σ ] − 1 /ξ − 1 + , if y j > u, 1 − λ, otherwise. (7) 5 Finally , the log-likeliho o d is given b y: 107 log L ( y 1 , . . . , y n ) = n X i =2 log L n ( y i − 1 , y i ) − n − 1 X i =2 log L d ( y i ) (8) 2.3 Return lev els 108 Most o ften, the ma jor issue of a n extre me v alue a nalysis is the q uantile es timation. As for the POT 109 approach, return le vel estimates can b e computed. Howev er, as a ll exc eedances are inferred, this 110 is do ne in a different wa y a s the dep endence betw een successive observ ations must be taken into 111 account. F or a sta tionary sequence Y 1 , Y 2 , . . . , Y n with a marg ina l distr ibution function F , Lindgren 112 and Ro o tze n [1 987] have sho wn that: 113 Pr [max { Y 1 , Y 2 , . . . , Y n } ≤ y ] ≈ F ( y ) nθ (9) where θ ∈ [0 , 1] is the extremal index a nd can b e interpreted as the rec ipr o cal of the mean cluster 114 size [Lea dbe tter, 19 83] - i.e. θ = 0 . 5 means that extreme (enough) even ts ar e exp ected to o ccur by 115 pair. θ = 1 (resp. θ → 0) corr esp onds to the indep endent (re sp. per fect dependent) ca se. 116 As a consequence, the qua ntile Q T corres p o nding to the T -year return p erio d is obta ined b y 117 equating equa tion (9) to 1 − 1 /T and solving for T . By definition, Q T is the observ ation that is 118 exp ected to b e exceeded onc e every T y ear s, i.e, 119 Q T = u − σ ξ − 1  1 − n λ − 1 h 1 − (1 − 1 /T ) 1 / ( nθ ) io − ξ  (10) It is worth emphasizing equation (9) a s it has a larg e impact on b oth theoretical and practical 120 asp ects. Indeed, for the AM approa ch, equation (9) is repla ced by 121 Pr [ max { Y 1 , Y 2 , . . . , Y n } ≤ y ] ≈ G ( y ) (11) where G is the distribution function o f the ra ndo m v ariable M n = max { Y 1 , Y 2 , . . . , Y n } , that is a 122 generalized extr eme v a lue distribution. In particula r, the e quations (9) and (11 ) differ a s the fir st 123 6 one is fitted to the whole obser v ations Y i , while the la tter is fitted to the AM ones. By definition, the 124 nu mber n Y of the Y i observ ations is mu ch larger than the size n M of the AM data set. Esp ecially , 125 for da ily data , n Y = 365 n M . 126 [Figure 1 ab out her e .] 127 F ro m equation (10), the extr emal index θ must be known to obtain quantile estimates. The 128 metho dology applied in this study is similar to the one sug gested by F awcett and W alshaw [2006]. 129 Once the Marko vian mo del is fitted, 100 Markov chains of length 200 0 were gener ated. F or each 130 chain, the extremal index is estimated using the estimator prop osed by F erro and Segers [2003 ] to 131 av oid iss ues related to the choice of declustering par a meter. In par ticular, the extr emal index θ is 132 estimated using the following equations: 133 ˆ θ ( u ) =        max  1 , 2 [ P N − 1 i =1 ( T i − 1) ] 2 ( N − 1) P N − 1 i =1 T 2 i  , if max { T i : 1 ≤ i ≤ N − 1 } ≤ 2 max  1 , 2 ( P N − 1 i =1 T i ) 2 ( N − 1) P N − 1 i =1 ( T i − 1)( T i − 2)  , otherwise (12) where N is the num b er o f observ ations exce eding the threshold u , T i is the inter-exceedance time, 134 e.g. T i = S i +1 − S i and the S i is the i -th exce e da nce time. 135 Lastly , the extremal index related to a fitted Marko v chain model is estimated us ing the sample 136 mean of the 1 00 extremal index es timations. Figure 1 represents the histo gram of these 10 0 extremal 137 index e s timations. In this s tudy , as lo ts o f time serie s are involv ed, the num b er a nd length of the 138 simulated Markov chains may b e too small to lea d to the most accura te extr emal index estimations; 139 but av oid intractable CPU times . If less sites a re c onsidered, it is preferable to increas e these tw o 140 v alues. 141 A preliminary s tudy (not shown here) demonstrates that, for quantile estimation, this pro cedur e 142 was more a ccurate than es timating θ using the estimator o f [Leadb etter, 1983]. This co nfirms the 143 conclusions dr awn b y F aw cett [200 5 ] for the extreme wind sp eed data. 144 7 3 Extreme V alue Dep endence Structure Assessmen t 145 Prior to per forming an y estimations , it is necessary to test whether: (a) the first o rder Markov 146 chain ass umption a nd (b) the extreme v alue dependence structure (equation (3)) are appropr iate 147 to mo del succes sive observ ations ab ov e the thresho ld u . 148 [Figure 2 ab out her e .] 149 [Figure 3 ab out her e .] 150 Figures 2 and 3 plot the a uto -corr e lation functions and the sca tter plots b etw e en tw o co nsecutive 151 observ ations for tw o different gaging sta tions. As the partial autoco rrelatio n co efficient at lag 1 152 is large, Figure 2 a nd 3 (left panels) corrob or ate the (a) hypothes is . How ever, as some partial 153 auto-cor relation co efficients are significant b eyond lag 1, it may suggest that a higher -order mo del 154 may b e more appropriate but do es not necessar ily mean that a first-orde r assumption is completely 155 flaw e d. Simplex plots [Cole s and T awn, 199 1] (not s hown) ca n b e used to assess the suitability of 156 a second-order as sumption over a first-order one. F or our a pplica tion, it seems that a first-order 157 mo del seems to b e v alid - except for the five slow est dynamic catchmen ts. 158 Though it is an imp ortant stage b ecause of its co nsequences on quantile estimates [Ledford 159 and T awn, 1996; Bortot and Cole s , 2 000], verifying the (b) hypo thesis is a co ns iderable task. An 160 ov erwhelming dep endence b etw een consecutive observ ations a t finite levels is not sufficient as it 161 do es not g ive any information ab out the dependence rela tion at asymptotic levels. F o r instance, 162 the overwhelming dependence a t lag 1 (Figure 2 and 3, right panels) do es cer tainly no t justify the 163 use of an asymptotic dep endent mo del. 164 [Figure 4 ab out her e .] 165 [Figure 5 ab out her e .] 166 Figures 4 and 5 plot the e volution of the χ ( ω ) and χ ( ω ) statistics as ω increases for t wo differe nt 167 sites. F or these fig ures, the confidence in terv a ls a re der ived by b o ots tr apping contiguous blo cks 168 to ta ke into ac count the success ive obser v ations dep endence [Ledfor d and T awn, 2003]. The χ ( ω ) 169 8 and χ ( ω ) sta tis tics s eem to depict tw o different asymptotic extremal dep endence. F ro m Figure 4, 170 it seems that lim χ ( ω ) ≫ 0 a nd lim χ ( ω ) = 1 for ω → 1 . On the contrary , Figur e 5 advo c ates 171 for lim χ ( ω ) = 0 and lim χ ( ω ) < 1 for ω → 1. Consequently , Figure 4 seems to conclude for an 172 asymptotic dependent ca se while Figur e 5 for an a s ymptotic indepe ndent ca se. 173 In theory , asymptotic (in)dependence sho uld not be assessed using scatterplots. How ever, these 174 t wo differ e nt features can be deduced from Figures 2 a nd 3. F or Figure 2, the scatterplot ( Y t − 1 , Y t ) 175 is incr easingly less spr ead as the observ ations b ecomes larger ; while increasingly mo re spr ead for 176 Figure 3. In other words, for the first cas e, the dependence seems to b ecome stronger at la rger 177 levels while this is the contrary for the second ca se. 178 [T able 1 ab out here.] 179 [Figure 6 ab out her e .] 180 Two sp ecific cases for different as y mptotic dep endence structures were illustra ted. T able 2 shows 181 the evolution of the χ ( ω ) statistics as ω incre a ses for all the sites und er study . Most of the stations 182 hav e s ig nificantly p ositive χ ( ω ) v a lues. In addition, only 1 3 sites have a 95% confidence interv al 183 that contains the 0 v a lue. F or 9 of these stations , the 95 % co nfidence interv als co r resp ond to the 184 theoretical low er and upp er b ounds; so that uncer tainties are to o larg e to determine the ex tremal 185 depe ndence c la ss. F or the χ statistic, results are less clear- cut. Figure 6 repr e sents the histogra ms 186 for χ ( ω ) for successive ω v alues. Despite only a few obser v ations being close to 1, most of the 187 stations hav e a χ ( ω ) v alue gr eater than 0.75 . These v alues can b e considered a s significantly high 188 as − 1 < χ ( ω ) ≤ 1 , for all ω . Consequently , mo dels of the form (1)–(3) may b e suited to mo del the 189 extremal dependenc e b etw een success ive observ ations. 190 Other m etho ds exist to test the extremal dep endence but w ere uncon vincing for our application 191 [Ledford and T a wn, 200 3; F alk and Michel, 2 006]. Indeed, the appro ach of F alk and Michel [2006] 192 do es not take into a ccount the dep endence b etw een Y t − 1 and Y t ; while the test of Ledford and 193 T awn [200 3] app ear s to b e p o orly discriminator y for o ur ca se study . 194 9 4 P erformance of the M ark o vian Mo dels on Quan tile Esti- 195 mation 196 4.1 Comparison b etw een Marko vian estimators 197 In this s e c tion, the p erforma nce of six different extremal dependence s tructures is analyzed on the 198 50 gaging s tations in tro duced in se c tion 1. These mo dels are: l og for the logistic, nl og for the 199 negative logistic, mix for the mixed mo dels and their relative asymmetric counterparts - e.g. al og , 200 anl og and amix . T o assess the impact o f the dep endence structure o n flo o d p eak estimation, the 201 efficiency of ea ch mo del to estimate qua nt iles with r eturn p erio ds 2, 10, 2 0, 50 and 1 00 years is 202 ev aluated. 203 As practitioners often hav e to dea l with small r ecord lengths in practice, the p er fo rmance of the 204 Marko vian mo dels is analyzed on all sub time series of length 5 , 1 0, 15 and 20 y ears . Finally , to 205 assess the efficiency for a ll the ga ging stations consider ed in this study , the no rmalized bias ( nbias ), 206 the v arianc e ( v ar ) and the no rmalized mean squa red err o r ( nmse ) are computed: 207 nbias = 1 N N X i =1 ˆ Q i,T − Q T Q T (13) v ar = 1 N − 1 N X i =1 ˆ Q i,T − Q T Q T − nbia s ! 2 (14) nmse = 1 N N X i =1 ˆ Q i,T − Q T Q T ! 2 (15) where Q T is the b enchmark T - year return level and ˆ Q i,T is the i -th estimate o f Q T . 208 [Figure 7 ab out her e .] 209 Figure 7 depicts the nbi as densities fo r Q 20 with a record le ng th of 5 years. It is ov erwhelming 210 that the extremal dep endence str ucture has a grea t impact on the estimation of Q 20 . Comparing 211 the t wo panels, it can b e noticed that t he symmetric dependence structures give spreader densities; 212 10 that is, mor e v a riable es timates. Independently of the symmetry , Figure 7 sho ws that the mixed 213 depe ndence fa mily is more accura te. 214 [T able 2 ab out here.] 215 T able 3 sho ws the nbias , var and nmse statistics for all the Marko v ian es timators as the rec o rd 216 length incr eases for quant ile Q 50 . This table co nfir ms results derived from Figure 7 . Indeed, the 217 asymmetric depe ndence structures give less v ar iable a nd biased estimates - as their n b ias and v ar 218 statistics are smaller. In addition, whatever the record leng th is, the Mar ko vian mo dels p erform 219 with the same hierarch y; that is the mix and amix models are by fa r the most accurate estimators 220 - i.e. with the s mallest nmse v a lues. The same res ults (not shown) ha ve b een fo und for other 221 quantiles. 222 F ro m an h ydro logical p o int of view, these t wo res ults are not surpris ing. T he symmetric models 223 suppo se that the v ariables Y t and Y t +1 are exchangeable. In our c ontext, e xchangeabilit y mea ns 224 that the time s eries are rev ersible - e.g. the time vector dir ection has no im p orta nce. When dea ling 225 with AM or POT and statio na ry time se r ies, it is a reas onable hypothesis . F or example, the MLE 226 remains the s ame with any permutations of the AM/POT sample. Ho wev er, when mo deling all 227 exceedances, the time dir ection can not be consider ed as reversible as flo o d h ydrogr aphs ar e clea rly 228 non symmetr ic . 229 [Figure 8 ab out her e .] 230 The Pick ands’ dep endence function A ( ω ) [Pick ands, 198 1] is another representation for the 231 extremal dep endence structure for any extreme v a lue dis tr ibution. A ( ω ) is related to the V function 232 in e q uation (3 ) a s follows: 233 A ( ω ) = V ( z 1 , z 2 ) z − 1 1 + z − 1 2 , ω = z 1 z 1 + z 2 (16) Figure 8 represents the Pick ands’ dep endence function for all the gag ing stations a nd the thr ee 234 asymmetric Mark ovian mo dels. One ma jor s p e c ificit y of the mixed models is that these models can 235 not a c count for perfect dependence cases. In par ticula r, the Pick ands’ dep endence functions for 236 11 the mixed mo dels sa tisfy A (0 . 5) ≥ 0 . 75 while A (0 . 5) ∈ [0 . 5 , 1] for the logis tic and negative log istic 237 mo dels. F rom Figure 8, it can b e seen that only few stations hav e a dependence function that could 238 not be modeled b y the amix model. Therefore, the d ep endence range limit ation of the amix model 239 do es no t s eem to o r e strictive. 240 In this section, the effect of the extremal dep endence structure has b ee n as sessed. It has b een 241 established tha t the symmetric mo dels a re hydrologically inconsistent as they could not r epro duce 242 the flo o d event asymmetry . In addition, for all the quantiles a nalyzed, th e asymmetric mixed model 243 is the most accurate for flo o d p eak estimations . Therefore, in the remainder of this section, only 244 the a m i x mo del will b e compared to c onv ent iona l POT estimator s. 245 4.2 Comparison b etw een amix and conv en tional POT estimators 246 In this s e c tion, the p er formance of the amix estimator is compa r ed to the estimators usually used 247 in floo d freq ue ncy analys is. F o r this pur p o se, the quantile es timates deriv ed from the Maximum 248 Likelihoo d Estima tor ( MLE ), the Unbiased and Biased Probability W eighted momen ts estimators 249 [Hosking and W allis, 1987 ] ( PWU a nd PWB resp ectively) ar e cons idered. 250 [Figure 9 ab out her e .] 251 Figure 9 depicts the nbias densities for the amix , M L E , P W U and P W B estimators r elated 252 to the Q 5 , Q 10 and Q 20 estimations with a recor d length of 5 years. It can be seen that amix is the 253 most ac c urate mo del for a ll quantiles. Indeed, the amix nbias densities a re the mos t sha rp with a 254 mo de close to 0. F o cusing o nly on “cla s sical” estimato r s (e.g. M LE , P W U and P W B ), there is 255 no estimator that p erform better than any o ther anytime. These tw o results a dvocate the use of 256 the a m i x mo del. 257 [T able 3 ab out here.] 258 T able 4 s hows the performance of each estimator to estimate Q 50 as the reco rd length increases . 259 It can b e se e n t hat the amix model p erforms better than the con ven tional estimators for the whole 260 range o f r ecord lengths a nalyzed. First, amix has the sa me bias than the conv ent ional estimators. 261 12 Thu s, the amix dep endence str uc tur e s eems to b e suited to estimate flo o d quantile estimates. 262 Second, b ecause of its smaller v ariance, amix is more a ccurate than M LE , P W U and P W B 263 estimators. This smaller v a riance is mainly a result of a ll of the exceedances (not o nly cluster 264 maxima) being used in the inference procedur e. Co nsequently , th e amix model has a smaller nmse 265 - ar ound half of the conv entional models ones. 266 [Figure 1 0 a b out here.] 267 Figure 1 0 shows the evolution of the nmse as the r eturn p erio d increas e s for the ami x , M LE , 268 P W U and P W B models. This figur e co rrob o rates the conclusio ns drawn from Figure 9 a nd T a ble 4 . 269 It can be seen that the amix mo de l has the smallest nmse , indep endently of the r eturn per io d and 270 the r e cord leng th. In addition, the amix b ecomes incre a singly more efficient as the retur n p erio d 271 increases - mos tly for retur n per io ds gr eater than 20 years. While the c o nv en tional estimators 272 present an erratic nmse behavior as the return p erio d incr eases, the amix mo del is the only one 273 that has a smo oth evolution. T o conclude, these res ults confirm that the amix mo del clear ly 274 improv es flo o d pe ak quantile estimates - esp ecially for la rge retur n p erio ds . 275 5 Inference on Other Flo o d Characteristics 276 As all exceedances are mo deled using a first o rder Mar ko v c hain, it is po ssible to infer other quan- 277 tities than flo o d p eaks - e.g. volume and dur ation. In this se ction, the ability of these Marko vian 278 mo dels to repro duce the flo o d duration is analyzed. F or this purp os e, the most severe flo o d h y- 279 drogra phs within each y ear are considered and normalized b y their p eak v alues. Consequently , 280 from this obs e r ved normalized hydrograph set, tw o flo o d characteristics derived from a da ta s et of 281 hydrographs [Robson and Reed, 1 9 99; Sauquet et a l., 2008 ] are considered: (a) the duration d mean 282 ab ov e 0.5 of the nor ma lized hydrograph set mean and (b) the media n d med of the durations ab ov e 283 0.5 o f each norma lized hydrograph. 284 13 5.1 Global Performa nce 285 [Figure 1 1 a b out here.] 286 Figure 11 plots the flo o d durations d mean and d med biases derived fro m the thr ee asy mmetr ic 287 Marko vian mo dels in function o f their empiric a l estimates. It can b e seen that no mo del lea ds to 288 accurate flo o d dur ation estimations. In addition, the extremal dep endence structure has a clear 289 impact on these estimations. In particular , the anlog and amix mo dels seem to underestimate 290 the floo d durations, while the alog mo del leads to o verestimations. Consequently , t wo different 291 conclusions can be drawn. Fir st, a s large durations are p o or ly estimated, higher order Markov 292 chains may b e o f interest. How ever, this is a considera ble task as higher dimensio nal multiv ariate 293 extreme v alue distr ibutions often lead to numerical pro blems. I ns tead of considering higher order , 294 another alterna tive ma y b e to change daily obse rv atio ns for d -day observ ations - where d is larger 295 than 1. Seco nd, it is overwhelming that the extremal dep endence str ucture affects the flo o d dur ation 296 estimations. As noticed in Section 2.1, there is no finite parametr ization for the extrema l dep endence 297 structure V - see E quation (3). Consequently , it seems reaso nable to suppose tha t one suited for 298 flo o d hydrograph estimation may exist. 299 [Figure 1 2 a b out here.] 300 Figure 12 depicts the observed norma lized mean hydrographs and the ones predicted by the three 301 asymmetric Mar ko vian mo dels. F or the J0621610 station (left panel), the normalized hydrograph is 302 well estimated by the th ree mo dels; wher eas for the L0400610 sta tion (right panel), the normalize d 303 hydrograph is p o or ly predicted. This result confirms the ina bilit y o f the three Markovian mo dels 304 to repr o duce lo ng flo o d even ts with daily da ta and a first order Ma rko v chain. 305 [Figure 1 3 a b out here.] 306 Figure 13 r epresents the biases related to ea ch v alue of the nor malized mean h ydrog raph. In 307 addition, to help estima to r comparison, the nmse is repo rted at the r ight side. It can be see n 308 that the al og mo del dr a matically ov er estimates the h ydrogr aph rising lim b while giving rea s onable 309 estimations for the falling phase. The anlog mo del slight ly ov eres timates the rising part while 310 14 strongly underestimates the falling o ne. The a m i x mo del always leads to underes timations - this 311 is mo r e pr onounced for the falling lim b. How ever, despite these different b ehaviors, these three 312 estimators se ems to have a similar per formance - in terms o f nmse . 313 [Figure 1 4 a b out here.] 314 Figure 1 4 repr esents the spa tia l dis tr ibution o f the nmse on the normalized mean hydrograph 315 estimation for each Markovian mo del. It seems that there is a specific spatial distr ibution. In 316 particular, the worst cases a re r elated to the middle par t of F rance. In addition, for differen t 317 extremal dep endence structures , the b est nm se v alues corr esp ond to differen t spatial lo ca tio ns. 318 The alog model is more accurate for the extreme north par t o f F rance; the anl og model is more 319 efficient for the east part of F rance; while the amix mo del pe rforms best in the middle part of 320 F ra nce. Consequen tly , a s at a global scale no mo de l is ac curate to estimate the normalized mean 321 hydrograph, it is w o rth trying to identify which catchmen t type s are rela ted to the b est es timations. 322 F or our data set, this is a considerable task. No standar d sta tistical technique lea d t o r easona ble 323 results. In par ticular, the principal comp o ne nt analysis, hierarchical cla ssification, sliced inv erse 324 regres s ion lead to no conclusio n ab o ut whic h catchmen t types ar e more suitable for our models. 325 Only a r egressio n approa ch gives some first g uidelines. F or this purp ose, a regr ession b etw een the 326 nbias on t he d mean estimation for ea ch asymmetric model a nd so me geomorphologic a nd hydrologic 327 indices ar e p er formed. The effect o f the dr ainage area, a n index of catchmen t slop e der ived from 328 the hypsometric curve [Roche, 1963], the Bas e Flow Index ( BFI ) [T allaksen and V an Lanen, 2 004, 329 Section 5.3 .3 ] and an index characterizing the ra infall pe r sistence [V asko v a and F ra nc` es, 19 98] are 330 considered. 331 nbias ( d mean ; anl og ) = 0 . 89 − 2 . 19 B F I , R 2 = 0 . 40 (17) nbias ( d mean ; amix ) = 0 . 49 − 1 . 74 B F I , R 2 = 0 . 43 (18) F ro m equations (17) and (18), the B F I v ariable expla ins a round 40% of the v a riance. Despite 332 the f act that a large v a riance pr op ortion is not tak en into acco unt , the B F I is clearly related to the 333 15 d mean estimation per formance. These equations indicate that the an l og (re s p. ami x ) mo del is mo re 334 accurate to repr o duce the d mean v aria ble for ga ging stations with a B F I a round 0 .4 (res p. 0.28 ). 335 These B F I v alues corr esp ond to catchmen ts with mo dera te up to flash flow r egimes r esp ectively . 336 These results corr ob ora te the ones derived fro m Figure 13: the firs t order Marko vian mo dels with a 337 1-day la g c o nditioning are not appropriate for long flo o d duration estimations. Consequently , while 338 no physiographic characteristic is rela ted to the al og p erforma nce; it is sugg ested, for such 1- day 339 lag c onditioning, to use the anl og and amix mo dels fo r quick basins. 340 6 Conclusion 341 Despite that univ a r iate EVT is widely applied in en vironmental s ciences, its m ultiv ariate extension 342 is rar ely c onsidered. This work tries to promo te the use of the MEVT in hydrology . In this work, the 343 biv aria te case was considered as the dep endence b etw een t wo successive observ ations was mo deled 344 using a fir st or der Ma rko v chain. This appro a ch has tw o main adv antages for practitioner s as: (a) 345 the num b er of data to be inferred increases considerably and (b) o ther featur es ca n be e s timated - 346 flo o d duration, volume. 347 In this study , a co mparison betw een six different extremal dependence structures (including 348 bo th symmetric a nd asymmetric for ms) has b een per formed. Results show that an a symmetric 349 depe ndence structure is more relev a nt. F ro m a hydrological point of view, this asymmetry is 350 rational as flo o d hydrographs are asymmetr ic. In par ticular, fo r our data, the asymmetric mixed 351 mo del gives the most acc ur ate flo o d pea k estima tions and clea rly improv es flo o d p eak estimatio ns 352 compared to co nv entional estimators indepe ndently of the r e turn p erio d co nsidered. 353 The ability of these Markovian mo de ls to es timate the floo d duration was carrie d out. It has 354 bee n s hown that, at first sig ht , no dep endence structure is able to repro duce the flo o d hydrogra ph 355 accurately . How ever, it s e ems that the anl og and a mix mo dels ma y b e more appropr iate when 356 dealing with mo derate up to flash flow regimes. These res ults depend strongly on the co nditio ning 357 term (i.e. Pr[ Y t ≤ y t | Y t − δ = y t − δ ]) of the first order Ma rko v chain and o n the auto-corre la tion 358 within the time series . In our application, δ = 1 and daily time step was considered. 359 16 More general conclusions can b e drawn. The weakness of the prop osed mo dels to derive consis- 360 ten t floo d hydrogr a phs may not be r e lated to t he daily time step but to the inadequa cy b etw een the 361 conditioning term and the floo d dyna mics. T o ensur e b etter results, higher order Marko v chains 362 may b e o f interest [F aw cett and W alshaw, 2 006]. Ho wev er, as numerical pro blems may ar ise, an- 363 other alternative may b e to still consider a first order c hain but to change the “conditioning lag 364 v alue” δ . In particular, for some basins, it may b e more relev ant to condition the Marko v chain 365 with a lar g er but more appropr ia te lag v alue. 366 Another o ption to improve the pro p o sed mo dels for flo o d hydrograph estimation is to use a mo re 367 suitable dep endence function V . As there is no finite parametriza tion for the extremal dep endence 368 structure, it seems reas o nable tha t one mo re appropr iate for flo o d hydrographs may exist. In this 369 work, r esults s how that the anlog model is more a ble to reproduce the hydrogr aph rising par t, while 370 the a l og is b etter for the falling phase. Define 371 V ( z 1 , z 2 ) = αV 1 ( z 1 , z 2 ) + β V 2 ( z 1 , z 2 ) where V 1 (resp. V 2 ) is the extremal de p endence f unction for the al og (resp. anl og ) mo del and α and 372 β ar e real constants such as α + β = 1. By definition, V is a new extr emal dep e ndence function. 373 In particular , V ma y co m bine the accuracy of the al og and anl og mo dels for b oth the rising 374 and falling part of the flo o d hydrograph. Ano ther a lternative may be to lo o k a t non-parametric 375 Pick ands’ dependence function estima tors [Cap´ era` a et al., 19 9 7] but that will require tec hniques to 376 simulate Marko v chains from these non- parametric e stimations. 377 All statistical analysis were perfor med within the R Developmen t Core T eam [2007] framework. 378 In particular, t he P O T pack age [Ribatet, 2 007] integrates the to ols that were developed to carry out 379 the modeling effort presented in this pape r . This pac k ag e is av ailable, free of c ha rge, at the website 380 http:/ /www. R- project.org , section CRAN, Pack a ges or at its own w ebpage htt p://po t.r- fo rge.r- project.org/ . 381 17 Ac kno wledgmen ts 382 The authors wish to thank the F rench H YDRO database for pro viding the data. Benjamin Rena rd 383 is acknowledged for criticizing thoro ughly the data analyzed in this s tudy . 384 A P arametrization for the Extremal Dep endence 385 This a nnex presents some useful results for the six extremal dep endence mo dels that have b een 386 considered in this w ork. As first order Marko v c hains w e r e used, only the biv aria te results are 387 describ ed. 388 [T able 4 ab out here.] 389 18 T able 1: Partial and mixed partia l deriv atives, definition domain, total independent and perfect dep endent cas es for ea ch extremal asymmetric depe ndenc e function V . Mo del Asymmetric Mo dels al og anl og am i x V ( x, y ) 1 − θ 1 x + 1 − θ 2 y +   x θ 1  − 1 /α +  y θ 2  − 1 /α  α 1 x + 1 y − h x θ 1  α +  y θ 2  α i − 1 /α 1 x + 1 y − (2 α + θ ) x + ( α + θ ) y ( x + y ) 2 V 1 ( x, y ) − 1 − θ 1 x 2 − θ 1 α 1 x − 1 α − 1   x θ 1  − 1 /α +  y θ 2  − 1 /α  α − 1 − 1 x 2 + θ − α 1 x α − 1 h x θ 1  α +  y θ 2  α i − 1 /α − 1 − 1 x 2 − 2 α + θ ( x + y ) 2 + 2 ( 2 α + θ ) x + ( α + θ ) y ( x + y ) 3 V 2 ( x, y ) − 1 − θ 2 y 2 − θ 1 α 2 y − 1 α − 1   x θ 1  − 1 /α +  y θ 2  − 1 /α  α − 1 − 1 y 2 + θ − α 2 y α − 1 h x θ 1  α +  y θ 2  α i − 1 /α − 1 − 1 y 2 − α + θ ( x + y ) 2 + 2 ( 2 α + θ ) x + ( α + θ ) y ( x + y ) 3 V 12 ( x, y ) α − 1 α ( θ 1 θ 2 ) 1 α ( xy ) − 1 α − 1   x θ 1  − 1 /α +  y θ 2  − 1 /α  α − 2 − ( α + 1)( θ 1 θ 2 ) − α ( xy ) α − 1 h x θ 1  α +  y θ 2  α i − 1 /α − 2 6 α +4 θ ( x + y ) 3 − 6 (2 α + θ ) x + ( α + θ ) y ( x + y ) 4 A ( w ) (1 − θ 1 ) (1 − w ) + (1 − θ 2 ) w + h (1 − w ) 1 α θ 1 α 1 + w 1 α θ 1 α 2 i α 1 −   1 − w θ 1  − α +  w θ 2  − α  − 1 α θw 3 + αw 2 − ( α + θ ) w + 1 Independenc e α = 1 or θ 1 = 0 or θ 2 = 0 α → 1 or θ 1 → 0 or θ 2 → = 0 α = θ = 0 T ota l dep endence α → 0 α → + ∞ Never r e a c h e d Constraint 0 < α ≤ 1, 0 ≤ θ 1 , θ 2 ≤ 1 α > 0, 0 < θ 1 , θ 2 ≤ 1 α ≥ 0, α + 2 θ ≤ 1 , α + 3 θ ≥ 0 19 References 390 P . Bo r tot and S. Coles. The multiv ariate Gaussian tail mo del: An applica tion to o ceanogr aphic 391 data. Journal of the Ro yal Statistic al S o ciety. Series C: Applie d Statistics , 49 (1):31–49 , 2000. 392 P . Bo rtot and J.A. T awn. Mo dels for the extremes of Markov ch ains. Biometrika , 85(4):851– 867, 393 1998. ISSN 00 0634 4 4. 394 P . Cap´ era ` a, A.-L. F oug` ere s , a nd C. Genest. A nonparametric estimation pro cedure for biv ariate 395 extreme v alue copulas . Biometrika , 84(3):567 –577, 1997 . ISSN 00 0634 44. 396 S. Co les. A n Intr o duction to Statistic al Mo del ling of Extre me V alues . Springer Series in Statistics. 397 Springers Ser ies in Statistics, Lo ndon, 20 01. 398 S. Coles and J.A. T awn. Modelling Extreme Multiv ariate Even ts. Journal of the R oyal Statistic al 399 So ciety. Series B (Metho dolo gic al) , 53 (2):377– 392, 1991. ISSN 00 3 5-92 46. 400 S. Coles, J.A. T awn, and R.L. Smith. A seasonal Marko v mo del for ex tremely lo w tempera ture. 401 Envir onmetrics , 5 :221–2 39, 1994. 402 S. Cole s, J. Hefferna n, and J.A. T awn. Dep endence Mea sures for Ex tr eme Value Analy s es. Extr emes , 403 2(4):339– 365, Decem b er 19 99. 404 J.M. Cunderlik and T.B.M.J. O uarda. Regional floo d-dura tion-frequency mo deling in the changing 405 environmen t. Journal of H ydr olo gy , 3 18(1-4 ):276–2 91, 2006. 406 M. F a lk and R. Michel. T esting for tail indepe ndence in extreme v alue mo dels. Annal. In st. Stat. 407 Math. , 5 8(2):261 –290, 2006. ISSN 00 2031 57. 408 L. F aw cett. Statistic al Metho dolo gy for the Estimation of Envir onmental Extr emes . PhD thesis, 409 Univ ers it y of Newcastle up on Tyne, 200 5 . 410 L. F aw cett and D. W a lshaw. Markov chain mo dels for e xtreme wind sp eeds. Envir onmetrics , 17 411 (8):795–8 09, 2006. ISSN 1 18040 09. 412 C.A.T. F erro and J. Segers. Inference for clusters of e xtreme v alues. Journal of the R oyal St atistic al 413 So ciety. Series B: Statistic al Metho dolo gy , 65 (2):545–5 56, 2003. ISSN 13 69741 2. 414 J. Gala mbos. Order statistics of samples from multiv ariate distributio ns . Journal of the Americ an 415 Statistic al Asso ciation , 9:674 – 680, 1975. 416 E.J. Gumbel. Biv aria te exp onential distributions. Journal of t he Americ an Statistic al Asso ciation , 417 55(292 ):698–7 07, 1960. 418 J.R.M. Ho s king and J.R. W allis. P ara meter and Q uantile Es timation for the Gener alized P a reto 419 Distribution. T e chnometrics , 29(3):3 39–34 9, 1987. 420 H. Jo e . F amilies o f min-stable multiv a riate exp onential and multiv ariate extreme v alue distributions. 421 Statist. Pr ob ab. Le tt. , 9:75– 82, 1990 . 422 T.R. Kjeldsen a nd D. Jones. Estimatio n of an index flo o d us ing data tra nsfer in the UK. Hydr ol. 423 Sci. J. , 52 (1 ):86–98 , 2007. ISSN 026 2 6667 . 424 20 T.R. Kjeldsen and D.A. Jones. Prediction uncertaint y in a median- ba sed index flo o d metho d using 425 L moment s. Water Resour c es R ese ar ch , 4 2(7):–, 20 06. ISSN 0043 1397 . 426 M.R. Le adb etter. Extremes and lo cal dep endence in stationary s equences. Pr ob ability The ory and 427 R elate d Fields (H ist oric al Ar chive) , 65 (2):291– 306, 1983. 428 A.W. Ledford and J.A. T a wn. Statistics for near indepe ndenc e in m ultiv ar iate extreme v alues. 429 Biometrika , 83:169 –187 , 19 96. 430 A.W. Ledford and J.A. T awn. Diag nostics for dep endence within time ser ies ex tremes. Journal of 431 the Ro yal S tatistic al So ciety. Series B: Statistic al Metho dolo gy , 6 5 (2):521– 543, 2003. 432 G. Lindgr en and H. Rootzen. E xtreme v alues: theory a nd technical applications. Sc andinavian 433 journal of statistics , 14(4 ):241–2 79, 1987. 434 J. Pick ands. Multiv ariate E xtreme Va lue Distr ibutio ns. In Pr o c e e dings 43r d Session Intern ational 435 Statistic al Institu te , 198 1 . 436 R Developmen t Cor e T eam. R: A L anguage and Envir onment for Statistic al Co mputing . R F ounda- 437 tion for Statistical Computing, Vienna, Austria, 2007. URL htt p://w ww.R- p roject.org . ISBN 438 3-900 051-0 7-0. 439 D.S. Reis Jr . and J.R. Stedinger. Bay esian MCMC floo d frequency ana ly sis with historica l infor - 440 mation. Journal of Hydr olo gy , 313 (1-2):97 –116, 2005. ISSN 0022 1694 . 441 S.I. Resnick. Extr eme Values, Re gular Variation and Point Pr o c esses . New–Y o r k: Springer –V erlag, 442 1987. 443 M. Riba tet. POT: Mo de lling Pea ks Over a Threshold. R N ews , 7 (1):34–3 6, April 200 7. 444 M. Ribatet, E. Sauque t, J.-M. Gr´ es illo n, and T.B.M.J. Ouarda. A reg ional Bay esian P OT mo del 445 for floo d frequency analysis. Sto chastic Envir onmental Rese ar ch and Ri sk Assessment (SERRA) , 446 21(4):327 –339 , 2007a. 447 M. Ribatet, E. Sa uquet, J.-M. Gr´ esillon, and T.B.M.J. Ouar da. Usefulness of the Reversible Jump 448 Marko v Chain Mon te Carlo Mo del in Regional Flo o d Frequency Analysis. Water R esour c es 449 R ese ar ch , 4 3 (8):W0840 3 , 2007b. doi: 1 0.102 9/200 6WR0055 25. 450 A.J. Robson and D.W. Reed. Flo o d Estimation Handb o ok , v o lume 3. Institute o f Hydrolog y , 451 W allingfo rd, 1 9 99. 452 M. Ro che. Hydr olo gie de surfac e . Ga uthier-Villars, Paris, 1963. 453 E. Sauquet, M.-H. Ramos, L. Cha pe l, and P . Ber nardar a. Stream flow scaling pr o p erties: inv estigat- 454 ing c hara cteristic scales fro m differen t statistical approaches. Ac c epte d in Hydr olo gic al Pr o c esses , 455 2008. doi: 10.100 2/hyp.6952. 456 R.L. Smith, J.A. T awn, and S.G. Co les. Ma rko v chain mo dels for thresho ld excee dances. Biometrika , 457 84(2):249 –268 , 1997. ISSN 000 6344 4. 458 L. T allakse n and H. V an Lane n. Hydr olo gic al Dr ought: Pr o c esses and Estimation Metho ds for 459 Str e amflow and Gr oundwater , volume 48. Elsevier, 2004. 460 21 J.A. T a wn. Biv ar iate extreme v alue theor y: Mo dels and estimation. Biometrika , 75(3):397 –415 , 461 1988. 462 I. V asko v a and F. F r anc` es. Rainfall analysis and r egionaliza tion computing in tensity-duration- 463 frequency c urves. In Flo o d Awar e Final R ep ort , pages 95– 108. Cemagr ef edition, 1998. 464 A. W erritty , J.L. Paine, N. Ma cdonald, J.S. Row an, and L.J. McEwen. Use of multi-proxy flo o d 465 records to improv e estimates o f flo o d risk: Lo wer Riv er Tay , Scotland. Catena , 66(1- 2):107– 119, 466 2006. ISSN 03 4181 6 2. 467 22 List of Figures 468 1 Histogram of the extremal index estimations from the 100 s imulated Mar ko v Cha ins 469 of length 2000. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 470 2 Autoco rrela tio n plot (left panel) and scatterplot of the time s e ries at lag 1 (right 471 panel) for the Somme river at Abbeville (E647 0910 ). . . . . . . . . . . . . . . . . . . 25 472 3 Autoco rrela tio n plot (left panel) and scatterplot of the time s e ries at lag 1 (right 473 panel) for the Mos e lle river at Noirg ue ux (A420063 0). . . . . . . . . . . . . . . . . . 26 474 4 Plot o f the χ and χ s ta tistics and t he related 9 5% confidence in ter v als f or the Somme 475 river at Abb eville (E6 4709 10). The solid blue lines a re the theor etical b o unds . . . . 27 476 5 Plot of the χ and χ s ta tistics and the related 95% interv als for the Mo selle river a t 477 Noirgueux (A42006 30). The solid blue lines are the theor etical b ounds. . . . . . . . . 28 478 6 Histograms of the χ ( ω ) statistics for differen t ω v alue s . Left pa nel: ω = 0 . 98, middle 479 panel: ω = 0 . 985 and r ight panel: ω = 0 . 99 . . . . . . . . . . . . . . . . . . . . . . . . 29 480 7 Densities of the normalized biases of Q 20 estimates for the symmetric Marko vian 481 mo dels (left panel) and the asymmetric o nes (righ t panel). T arg et site record leng th: 482 5 years. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 483 8 Representation of the Pick ands’ dep endence functions fo r the 5 0 ga ging stations. 484 Left panel : alog , middle panel: anlog a nd r ight panel: amix . “+ ” re pr esents the 485 theoretical dep e ndence bo und for the amix mo del. . . . . . . . . . . . . . . . . . . . 31 486 9 Densities of the normalize d bias es for the ami x model and the M LE , P W U , and 487 P W B estimator s for quantiles Q 5 (left panel), Q 10 (middle panel) and Q 20 (right 488 panel). Record length: 5 years. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 489 10 Evolution o f the nmse as the return p erio d increa ses for the amix , M LE , P W U and 490 P W B estimators . Record leng th: (a) 5 years, (b) 1 0 years, (c) 15 years a nd (d) 20 491 years. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 492 11 d mean and d med normalized bias es in function of the theoretical v alues for the three 493 asymmetric Marko vian mo dels. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 494 12 Observed and sim ulated normalized mean hydrographs for the J0621 610 (left panel) 495 and the L 0 4006 10 (right panel) stations . . . . . . . . . . . . . . . . . . . . . . . . . . 35 496 13 Evolution of the bia s es for the norma lized mean hydrogr a ph estimations in function 497 of the dis ta nce of the flo o d p e a k time. . . . . . . . . . . . . . . . . . . . . . . . . . . 36 498 14 nmse spa tial distribution acc o rding for the three Marko vian models. Left panel: 499 al og , middle pa ne l: anl og and right panel: a mix . The radius is pr o p ortional to the 500 nmse v alue. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 501 23 θ Frequency 0.2 0.3 0.4 0.5 0.6 0.7 0 5 10 15 20 Figure 1: Histogram of the extr emal index estimations from the 10 0 simulated Markov Chains of length 2000. 24 Figure 2: Auto co rrelatio n plot (left pa nel) a nd scatterplo t of the time ser ie s at lag 1 (right panel) for the Somme river at Abbeville (E647 0910 ). 25 Figure 3: Auto co rrelatio n plot (left pa nel) a nd scatterplo t of the time ser ie s at lag 1 (right panel) for the Mos e lle river at Noirg ue ux (A420063 0). 26 0.0 0.2 0.4 0.6 0.8 1.0 −1.0 −0.5 0.0 0.5 1.0 ω χ 0.0 0.2 0.4 0.6 0.8 1.0 −1.0 −0.5 0.0 0.5 1.0 ω χ Figure 4: Plot of the χ and χ statistics a nd the related 95% confidence interv als for the So mme river at Abb eville (E6 4709 10). The solid blue lines a re the theoretical b ounds. 27 0.0 0.2 0.4 0.6 0.8 1.0 −1.0 −0.5 0.0 0.5 1.0 ω χ 0.0 0.2 0.4 0.6 0.8 1.0 −1.0 −0.5 0.0 0.5 1.0 ω χ Figure 5: Plot of the χ and χ statistics and the r elated 95% interv als for the Mo s elle river at Noirgueux (A42006 30). The solid blue lines a re the theor etical b ounds. 28 χ ( ω ) Frequency 0.6 0.7 0.8 0.9 1.0 0 5 10 15 20 χ ( ω ) Frequency 0.6 0.7 0.8 0.9 1.0 0 5 10 15 χ ( ω ) Frequency 0.6 0.7 0.8 0.9 1.0 0 5 10 15 Figure 6: Histograms of the χ ( ω ) s tatistics for different ω v alues. Le ft panel: ω = 0 . 98 , middle panel: ω = 0 . 9 85 and right panel: ω = 0 . 9 9. 29 −150 −100 −50 0 50 100 0.000 0.005 0.010 0.015 0.020 0.025 NBIAS Density MIX NLOG LOG −150 −100 −50 0 50 100 0.000 0.005 0.010 0.015 0.020 0.025 NBIAS Density AMIX ANLOG ALOG Figure 7: Densities of the normalized bia ses of Q 20 estimates for the symmetric Mar ko vian models (left panel) and the asy mmetric ones (right panel). T arget site reco r d length: 5 years. 30 0.0 0.2 0.4 0.6 0.8 1.0 0.5 0.6 0.7 0.8 0.9 1.0 w A (w) 0.0 0.2 0.4 0.6 0.8 1.0 0.5 0.6 0.7 0.8 0.9 1.0 w A (w) 0.0 0.2 0.4 0.6 0.8 1.0 0.5 0.6 0.7 0.8 0.9 1.0 w A (w) Figure 8: Repr esentation of the Pick ands’ dependence functions for the 50 g aging stations. Left panel : al og , middle panel: anl og a nd rig ht panel: ami x . “+” r epresents the theor e tical dep endence bo und for the amix mo del. 31 −100 −50 0 50 100 150 0.000 0.010 0.020 0.030 NBIAS Density AMIX MLE PWB PWU −150 −50 0 50 100 150 200 0.000 0.005 0.010 0.015 0.020 0.025 NBIAS Density AMIX PWU MLE PWB −50 0 50 100 150 0.000 0.005 0.010 0.015 NBIAS Density AMIX PWB MLE PWU Figure 9: Densities of the normalized bia s es for the amix mo del and the M LE , P W U , and P W B estimators for quantiles Q 5 (left panel), Q 10 (middle panel) and Q 20 (right panel). Reco rd leng th: 5 years. 32 Figure 10: Evolution of the nms e as the return p erio d increase s fo r the amix , M LE , P W U and P W B es timators. Record length: (a) 5 years, (b) 1 0 years, (c) 15 years and (d) 20 years. 33 Figure 11: d mean and d med normalized bia ses in function of the theoretical v alues for the three asymmetric Marko vian mo dels. 34 −2 0 2 4 6 0.0 0.2 0.4 0.6 0.8 1.0 (days) Normalised Time Normalised Flows −10 0 10 20 0.0 0.2 0.4 0.6 0.8 1.0 (days) Normalised Time Normalised Flows obs alog anlog amix Figure 12: O bserved a nd sim ulated no rmalized mean h ydrogr aphs for the J062 1610 (left pa nel) and the L0 4006 10 (right panel) stations . 35 t 0 − 4 t 0 − 3 t 0 − 2 t 0 − 1 t 0 t 0 + 1 t 0 + 2 t 0 + 3 t 0 + 4 t 0 + 5 t 0 + 6 t 0 + 7 alog anlog amix NBIAS (%) −40 −20 0 20 40 nmse 0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14 Figure 13: Evolution o f the biases fo r the no rmalized mean hydrograph estimations in function of the dis ta nce of the flo o d p eak time. 36 Figure 1 4: nmse spatial distribution accor ding for the three Ma rko vian mo dels. Left panel: al og , middle pa nel: anl og and rig ht panel: amix . The radius is pr op ortiona l to the nmse v alue . 37 List of T a bles 502 1 Partial and mixed par tial deriv atives, definition domain, total indep endent and p er- 503 fect dep endent cases for each extrema l asymmetric dep endence function V . . . . . . 19 504 2 χ ( ω ) statistics for all stations. ω = 0 . 98 , 0 . 9 85 , 0 . 99 . . . . . . . . . . . . . . . . . . . . 39 505 3 Several characteristics of the Markovian estimators on Q 50 estimation a s the record 506 length increases. Standard erro r s are r ep orted in brackets. . . . . . . . . . . . . . . . 40 507 4 Several characteristics of the amix , M LE , P W U and P W B estimators for Q 50 508 estimation as the r ecord leng th incr eases. Standard er r ors are rep or ted in brackets. . 41 509 5 Partial and mixed par tial deriv atives, definition domain, total indep endent and p er- 510 fect dep endent cases for each extrema l symmetric dep endence function V . . . . . . . 42 511 38 T able 2 : χ ( ω ) sta tis tics for all stations. ω = 0 . 9 8 , 0 . 98 5 , 0 . 99. Stations ω = 0 . 9 8 ω = 0 . 985 ω = 0 . 99 χ ( ω ) 95% C.I. χ ( ω ) 9 5% C.I. χ ( ω ) 95% C.I. A34720 10 0.67 (-0.0 2, 1 .0 0) 0.60 (-0.02, 1.00 ) 0.57 (-0.01, 1.00) A42006 30 0.53 ( 0.21 , 0.81 ) 0.45 ( 0.07, 0 .77) 0.38 (-0.01, 0.7 6) A42506 40 0.55 ( 0.27 , 0.82 ) 0.49 ( 0.18, 0 .76) 0.41 ( 0.0 2, 0.71 ) A54310 10 0.44 (-0.0 2, 1 .0 0) 0.44 (-0.02, 1.00 ) 0.41 (-0.01, 1.00) A57306 10 0.59 ( 0.25 , 0.94 ) 0.56 ( 0.20, 0 .90) 0.50 ( 0.0 7, 0.97 ) A69410 10 0.62 ( 0.22 , 0.99 ) 0.60 ( 0.16, 1 .00) 0.56 ( 0.0 6, 1.00 ) A69410 15 0.63 ( 0.29 , 0.95 ) 0.60 ( 0.20, 0 .96) 0.58 ( 0.1 7, 0.98 ) D01370 1 0 0.39 ( 0.04, 0.69 ) 0.33 (-0.02, 0.67) 0.2 8 (-0 .01, 0.6 9) D01565 1 0 0.59 ( 0.25, 0.88 ) 0.55 ( 0.20, 0 .86) 0.53 ( 0.1 4, 0.92 ) E1727 510 0 .62 ( 0.18, 0.91) 0.59 ( 0.16, 0 .9 3) 0.47 (-0.01, 0.89 ) E1766 010 0 .63 ( 0.23, 0.98) 0.59 ( 0.17, 0 .9 6) 0.54 ( 0.0 9 , 0.96) E3511 220 0 .59 ( 0.10, 1.00) 0.53 (-0.02 , 1.00) 0.50 (-0.01, 0.9 9) E4035 710 0 .77 ( 0.02, 1.00) 0.68 (-0.02 , 1.00) 0.60 (-0.01, 1.0 0) E5400 310 0 .88 ( 0.30, 1.00) 0.89 ( 0.29, 1 .0 0) 0.83 ( 0.1 3 , 1.00) E5505 720 0 .91 ( 0.24, 1.00) 0.87 ( 0.09, 1 .0 0) 0.86 ( 0.0 2 , 1.00) E6470 910 0 .96 ( 0.40, 1.00) 0.94 ( 0.25, 1 .0 0) 0.98 ( 0.0 0 , 1.00) H04000 10 0.84 ( 0.12 , 1.00 ) 0.83 ( 0.02, 1 .00) 0.78 (-0.01, 1.0 0) H15010 10 0.82 ( 0.36 , 1.00 ) 0.90 ( 0.39, 1 .00) 0.84 ( 0.2 6, 1.00 ) H23420 10 0.68 ( 0.31 , 1.00 ) 0.67 ( 0.25, 1 .00) 0.60 ( 0.1 1, 1.00 ) H50710 10 0.75 ( 0.30 , 1.00 ) 0.76 ( 0.22, 1 .00) 0.75 ( 0.1 5, 1.00 ) H51720 10 0.80 ( 0.47 , 1.00 ) 0.77 ( 0.42, 1 .00) 0.73 ( 0.3 0, 1.00 ) H62010 10 0.69 ( 0.29 , 1.00 ) 0.69 ( 0.14, 1 .00) 0.69 ( 0.0 8, 1.00 ) H74010 10 0.85 ( 0.46 , 1.00 ) 0.85 ( 0.38, 1 .00) 0.81 ( 0.2 7, 1.00 ) I92210 10 0.67 ( 0.23 , 1.0 0) 0.66 ( 0.19 , 1.00) 0 .5 9 ( 0.0 4, 1.00 ) J0621 610 0 .61 ( 0.25, 0.92) 0.58 ( 0.20, 0.9 4) 0.51 ( 0.08 , 0.91) K0433 010 0 .59 ( 0.22, 0.91) 0.54 ( 0.15, 0.8 9) 0.45 ( 0.0 0 , 0.85) K0454 010 0 .71 ( 0.37, 1.00) 0.67 ( 0.24, 1.0 0) 0.65 ( 0.1 4 , 1.00) K0523 010 0 .62 (-0.02, 1.0 0 ) 0.58 (-0.02, 1.00) 0 .53 (- 0 .01, 1.0 0) K0550 010 0 .61 ( 0.22, 0.94) 0.57 ( 0.15, 0.9 4) 0.54 ( 0.0 7 , 1.00) K0673 310 0 .67 ( 0.24, 1.00) 0.65 ( 0.18, 1.0 0) 0.66 ( 0.0 7 , 1.00) K0910 010 0 .65 (-0.02, 1.0 0 ) 0.61 (-0.02, 1.00) 0 .58 (- 0 .01, 1.0 0) K1391 810 0 .68 ( 0.27, 1.00) 0.64 ( 0.16, 0.9 8) 0.60 ( 0.0 6 , 0.96) K1503 010 0 .69 ( 0.38, 0.98) 0.67 ( 0.30, 0.9 8) 0.64 ( 0.2 3 , 1.00) K2330 810 0 .68 ( 0.29, 1.00) 0.66 ( 0.22, 1.0 0) 0.62 ( 0.0 9 , 1.00) K2363 010 0 .65 ( 0.26, 0.98) 0.66 ( 0.16, 1.0 0) 0.61 ( 0.0 1 , 1.00) K2514 010 0 .61 ( 0.24, 1.00) 0.61 ( 0.21, 1.0 0) 0.58 ( 0.1 2 , 1.00) K2523 010 0 .53 (-0.02, 1.0 0 ) 0.53 (-0.02, 1.00) 0 .51 (- 0 .01, 1.0 0) K2654 010 0 .68 ( 0.37, 1.00) 0.68 ( 0.31, 1.0 0) 0.60 ( 0.1 0 , 1.00) K2674 010 0 .60 ( 0.25, 0.89) 0.58 ( 0.22, 0.9 4) 0.54 ( 0.0 8 , 0.95) K2871 910 0 .62 ( 0.26, 0.95) 0.57 ( 0.15, 0.9 4) 0.56 ( 0.1 0 , 0.97) K2884 010 0 .62 ( 0.25, 1.00) 0.57 ( 0.17, 0.9 7) 0.59 ( 0.1 6 , 1.00) K3222 010 0 .56 ( 0.21, 0.90) 0.53 ( 0.18, 0.9 3) 0.46 ( 0.1 1 , 0.89) K3292 020 0 .59 ( 0.27, 0.91) 0.57 ( 0.17, 0.9 1) 0.48 ( 0.0 7 , 0.90) K4470 010 0 .76 ( 0.39, 1.00) 0.77 ( 0.40, 1.0 0) 0.73 ( 0.2 7 , 1.00) K5090 910 0 .64 ( 0.27, 0.93) 0.64 ( 0.26, 0.9 6) 0.58 ( 0.1 2 , 0.98) K5183 010 0 .57 ( 0.14, 0.91) 0.56 ( 0.15, 0.9 6) 0.53 ( 0.0 6 , 0.97) K5200 910 0 .63 ( 0.24, 0.93) 0.62 ( 0.20, 0.9 5) 0.56 ( 0.1 1 , 0.97) L01406 10 0 .7 3 ( 0.2 3, 1 .00) 0.66 ( 0.1 5, 1.00 ) 0.58 (-0.01, 1.00) L02315 10 0 .5 9 ( 0.1 6, 0 .91) 0.55 ( 0.1 1, 0.92 ) 0.53 (-0.01, 0.92) L04006 10 0 .7 4 (-0 .02, 1.00) 0.65 (-0.02, 1.0 0) 0.61 (-0.01, 1.00 ) 39 T able 3: Several characteristics o f the Marko vian estimator s on Q 50 estimation as the record length increases. Standard er rors are rep or ted in brackets. Mo del 5 years 10 years 15 years 20 yea r s nbias v ar n mse nbias v ar n mse nbia s v ar nmse nbias v ar n m s e l og -0.35 0.54 0.66 -0.32 0.32 0.42 -0.3 0 0.23 0.3 2 -0.28 0.17 0 .2 5 (16e-3) (22e-3) (18e-3) (12e-3) (12 e-3) (14e-3) (11e-3) (9e-3) (12 e-3) (9e-3) (7e-3) ( 1 1 e - 3 ) nl og -0.21 0.20 0.24 -0.20 0.11 0.15 -0.18 0.0 8 0.12 -0.18 0.06 0 .0 9 (10e-3) (7e-3) (11e-3 ) (7e-3) (4 e-3) (9e -3) (6e-3) (3e-3) (8e-3) (5e-3) (2e-3 ) ( 7 e - 3 ) mix -0.08 0.14 0.14 -0.07 0.08 0.08 -0.0 6 0.05 0.0 6 -0.05 0.04 0 .0 4 (8e-3) (5e-3) (8 e-3) (6e-3) (2e-3) (6e-3 ) (5e-3) (2e - 3) (5e-3) (4e-3) (1e -3) ( 5 e - 3 ) al og -0.15 0.39 0.41 -0.13 0.22 0.24 -0.11 0.16 0 .17 -0.10 0.12 0 .1 3 (14e-3) (15e-3) (14e-3) (10e-3) (9e-3 ) (11e- 3) (9e-3) (6 e -3) (9e-3 ) (8e-3) (4 e -3) ( 8 e - 3 ) anl og - 0 .10 0.20 0.21 -0.09 0.11 0.12 -0.08 0.08 0.09 -0.08 0 .06 0 .0 6 (10e-3) (7e-3) (10e-3 ) (7e-3) (4 e-3) (8e -3) (6e-3) (2e-3) (6e-3) (5e-3) (2e-3 ) ( 6 e - 3 ) amix -0.0 6 0.11 0.12 -0.05 0.06 0.06 -0.04 0.0 4 0.05 -0.03 0.0 3 0 .0 3 (7e-3) (4e-3) (7 e-3) (6e-3) (2e-3) (6e-3 ) (5e-3) (1e - 3) (5e-3) (4e-3) (1e -3) ( 4 e - 3 ) 40 T able 4: Sev eral c ha racteristics of the amix , M LE , P W U and P W B estimato rs for Q 50 estimation as the re c ord length increa ses. Standard erro rs are rep orted in brack ets. Mo del 5 years 10 years 15 years 20 years nbias v ar n mse nbias v ar nmse nbias v ar nmse nbias v ar nm s e amix -0.0 6 0.11 0.12 -0.05 0.06 0.0 7 -0.04 0.04 0.0 5 -0.04 0.03 0.0 3 (8e-3) (4e-3) (8 e-3) (6e-3 ) (2e - 3) (6e-3) (5e-3) (1 e-3) (5e-3) (4e-3) (1e-3) (4e - 3 ) M LE -0.13 0.2 5 0.27 -0.14 0 .13 0.1 4 -0.13 0.08 0.1 0 -0.11 0.05 0.0 7 (12e-3) (15e-3) (12e-3) (8e-3) (6e-3) (9e- 3 ) (7e-3) (3e-3) (7e- 3) (5e-3) (2e-3) (6 e - 3 ) P W U 0.08 0.30 0.31 -0.01 0.15 0 .15 -0.03 0.10 0 .10 -0.03 0.06 0 .0 6 (13e-3) (13e-3) (13e-3) (9e-3) (6e-3) (9e- 3 ) (7e-3) (3e-3) (7e- 3) (6e-3) (2e-3) (6 e - 3 ) P W B -0.07 0.20 0.21 -0.10 0.11 0.1 2 -0.11 0.08 0.0 9 -0.10 0.05 0.0 6 (10e-3) (8e-3) (11e-3 ) (7e-3) (4e-3) (8e- 3) (6e-3) (2e-3) (7e - 3) (5e-3) (1e-3) (6 e - 3 ) 41 T able 5: Partial and mixed partial deriv atives, definition domain, total indepe ndent and p er fect depe ndent ca ses for ea ch extremal symmetric dep endence function V . Mo del Symmetric Mo dels l og nl og mix V ( x, y )  x − 1 /α + y − 1 /α  α 1 x + 1 y − ( x α + y α ) − 1 /α 1 x + 1 y − α x + y V 1 ( x, y ) − x − 1 α − 1 V ( x, y ) α − 1 α − 1 x 2 + x α − 1 ( x α + y α ) − 1 α − 1 − 1 x 2 + α ( x + y ) 2 V 2 ( x, y ) − y − 1 α − 1 V ( x, y ) α − 1 α − 1 y 2 + y α − 1 ( x α + y α ) − 1 α − 1 − 1 y 2 + α ( x + y ) 2 V 12 ( x, y ) − ( xy ) − 1 α − 1 1 − α α V ( x, y ) α − 2 α − ( α + 1)( xy ) α − 1 ( x α + y α ) − 1 α − 2 − 2 α ( x + y ) 3 A ( w ) h (1 − w ) 1 α + w 1 α i α 1 − [(1 − w ) − α + w − α ] − 1 α 1 − w (1 − w ) α Independenc e α = 1 α → 0 α = 0 T ota l dep endence α → 0 α → + ∞ Never reached Constraint 0 < α ≤ 1 α > 0 0 ≤ α ≤ 1 42