ISSN 1581-0267 UDK/UDC: 519.21:532.433:556 Prej eto/Received: 04. 02. 2014 Sprejeto/Accepted: 14. 05. 2014 Izvirni znanstveni članek - Original scientific paper PRIMERJAVA MED KLASIČNIMI UNIVARIATNIMI VERJETNOSTNIMI ANALIZAMI IN BIVARIATNIMI Z UPORABO FUNKCIJE KOPULA COMPARISON BETWEEN CLASSICAL UNIVARIATE FREQUENCY ANALYSIS AND BIVARIATE ANALYSIS WITH COPULA 1 Fakulteta za gradbeništvo in geodezijo, Univerza v Ljubljani, Hajdrihova 28, 1000 Ljubljana Izvleček Verjetnostne analize so osnova za določanje projektnih pretokov. Običajno se pri analizah upošteva le ena spremenljivka, večinoma letna konica pretoka. Ker pa so hidrološki pojavi določeni z več medsebojno odvisnimi spremenljivkami, je pri analizah smiseln multivariaten pristop. Primer takega postopka je funkcija kopula. Klasične, univariatne verjetnostne analize so še vedno predpogoj za izvedbo analiz z uporabo funkcije kopula. Bivariatno verjetnostno analizo z uporabo funkcije kopula smo naredili za letne konice pretokov in pripadajoče volumne vodomerne postaje Litija na reki Savi. Uporabili smo tri funkcije kopula iz Arhimedove družine in parametre ocenili s pomočjo Kendallovega koeficienta korelacije (metoda momentov). Izračunali smo skupne povratne dobe in jih primerjali z univariatnimi povratnimi dobami. Ugotovili smo, da razlike med povratnimi dobami niso zanemarljivo majhne. Pri dogodku iz leta 1990, ki je bil največji v opazovanem obdobju, povratna doba TAND znaša 92 let, TR pa 17 let. Univariatni povratni dobi konic in volumnov pa ležita v območju med omenjenima povratnima dobama. S pomočjo statističnih in grafičnih kriterijev ustreznosti smo ugotovili, da kopula Gumbel-Hougaard za obravnavani primer izkazuje boljše rezultate kot kopula Clayton ali kopula Frank. Ključne besede: funkcija kopula, bivariatne analize, Kendallov koeficient korelacije, Arhimedova družina kopul, skupna povratna doba OR, skupna povratna doba AND. Frequency analyses are a basis for designing discharge estimations. Univariate flood frequency analyses are usually applied in hydrological practice. Hydrological processes are multivariate, however multivariate analyses are needed. Copula function can be used for multivariate modelling. Classical univariate flood frequency analyses are a precondition for the copula analyses. Flood frequency analyses were made on the annual maximum series data from gauging station Litija on the Sava River. Three copulas from the Archimedean family were used; parameters were estimated with method of moments (based on the Kendall correlation coefficient). Some joint return periods were calculated and compared with the univarite return periods. Differences between return periods were not negligible. In the case of a flood event in 1990, which was the largest in the observed period, TAND was 92 years and T°R was 17 years. Univariate return periods lay * Stik / Correspondence: mojca.sraj@fgg.uni-lj.si © Šraj M. et al.; Vsebina tega članka se sme uporabljati v skladu s pogoji licence Creative Commons Priznanje avtorstva -Nekomercialno - Deljenje pod enakimi pogoji 4.0. © Šraj M. et al.; This is an open-access article distributed under the terms of the Creative Commons Attribution - Non Commercial -Share Alike 4.0 Licence. v 1 * 1 1 Mojca Sraj1 , Nejc Bezak1, Mitja Brilly1 Abstract between these two values. Statistical and graphical performance measures were used to choose the best fit copula function. Gumbel-Hougaard copula gave better results than Clayton and Frank copulas. Keywords: copula function, bivariate analyses, Kendall's correlation coefficient, Archimedean copula family, joint return period OR, joint return period AND. 1. Uvod Verjetnostne analize določajo razmerje med velikostjo in verjetnostjo nastanka opazovanega pojava. So osnova za določanje projektnih pretokov pri dimenzioniranju hidrotehničnih objektov, določanju poplavnih območij, pri upravljanju z vodami, obvladovanju tveganja zaradi poplav ipd. Ker gre večinoma za ukrepe, ki neposredno vplivajo na ogroženost določenega območja in prebivalstva (Đurović & Mikoš, 2004; Brilly & Polič, 2005; Kobold et al., 2005; Mikoš et al., 2004;), je potrebno verjetnostne analize izvesti s čim večjo zanesljivostjo in s hkratnim upoštevanjem več spremenljivk, ki določajo visokovodni val. V svetovni praksi se največkrat uporabljajo univariatne verjetnostne analize, kjer se večinoma upoštevajo le konice pretokov. To velja tudi za Slovenijo. Če želimo verjetnostne analize izvesti zanesljivo, je smiselno poleg konic pretokov statistično obravnavati tudi volumne in čase trajanja visokovodnih valov (Karmakar & Simonovic, 2008). Z upoštevanjem le ene spremenljivke lahko podcenimo dejansko povratno dobo opazovanega ekstremnega dogodka (Favre et al., 2004; Renard & Lang, 2007). Za izvedbo verjetnostne analize lahko uporabimo klasične bivariatne porazdelitvene funkcije, kot so: normalna, logaritemsko normalna, gama in generalizirana porazdelitev ekstremnih vrednosti (Genest & Favre, 2007). Pogoj za uporabo klasične bivariatne analize je, da sta obe spremenljivki definirani z isto porazdelitveno funkcijo (Genest & Favre, 2007). Kot alternativa klasičnemu bivariatnemu pristopu se v zadnjih letih v svetovni hidrološki praksi uveljavlja funkcija kopula, ki so jo pred tem že uspešno uporabljali zlasti v ekonomiji. V zadnjem desetletju so jo že uporabili tudi v nekaterih hidroloških študijah, predvsem pri analizah visokovodnih valov (Genest & Favre, 2007; Karmakar & Simonovic, 2008; Karmakar & Simonovic, 2009; Zhang & Singh, 2006; Zhang & Singh, 2007), ter tudi pri analizah suše (Ma et al., 2011; Wong et al., 2010), analizah padavin (Balistrocchi & Bacchi, 2011; Grimaldi & Serinaldi, 2006), preverjanju ustreznosti prelivnega objekta (De Michele et al., 2005), analizah tveganja (angl. risk) (Chen et al., 2012; Ganguli & Reddy, 2013). Glavna prednost funkcije kopula pred klasičnimi bivariatnimi analizami je, da lahko spremenljivke opišemo z različnimi robnimi porazdelitvenimi funkcijami (angl. marginal distribution functions) (Genest & Favre, 2007). Poleg tega sta izbira ustreznih robnih porazdelitev in uporaba funkcije kopula dva ločena koraka. Je pa res, da take kompleksne, multivariatne verjetnostne analize z uporabo funkcije kopula zahtevajo predhodno natančno analizo lastnosti in odvisnosti posameznih spremenljivk ter verjetnostne analize vsake posamezne spremenljivke (Bezak et al., 2013; Šraj et al., 2012), kar vsekakor pomeni mnogo več vloženega časa in dela. Zastavlja se vprašanje, kolikšna je razlika razmerij med velikostjo in verjetnostjo pojava pri uporabi univariatnih verjetnostnih analiz in bivariatnih analiz z uporabo funkcije kopula. V prispevku so prikazani rezultati bivariatnih analiz z uporabo funkcije kopula za konice in volumne, možno pa je tudi hkratno upoštevanje vseh treh spremenljivk (konice, volumni in časi trajanja), ki določajo visokovodne valove (trivariatne analize). Cilji raziskave so bili: (i) narediti klasične, univariatne verjetnostne analize konic in volumnov visokovodnih valov, (ii) izvesti bivariatne verjetnostne analize konic in volumnov z uporabo funkcije kopula, (iii) primerjati univariatne povratne dobe z različnimi skupnimi povratnimi dobami (angl. joint return periods), določenimi z bivariatno kopulo. 2. Podatki Za potrebe raziskave smo uporabili podatke z vodomerne postaje Litija I na reki Savi (slika 1). Hidrološka postaja Litija je najstarejša postaja Agencije RS za okolje, kjer vodostaje opazujejo od leta 1893, podatki pa so na razpolago od leta 1895 (Ulaga, 2011). Leta 1953 je bila postaja zaradi zasipavanja dna prestavljena približno 500 m gorvodno in opremljena z limnigrafom, kar pomeni, da so od takrat dalje na voljo podatki o celotnih hidrogramih. Za analize smo torej uporabili 58 let podatkov (1953-2010) dnevnih vrednosti pretokov z vključenimi maksimalnimi konicami, ki so bili izmerjeni z limnigrafom Seba Omega. Slika 1: Hidrološka postaja Litija I na reki Savi (ARSO, 2013;. Figure 1: Hydrologie station Litija I on the Sava River (ARSO, 2013;. Za uporabo dnevnih vrednosti pretokov z vključenimi absolutnimi konicami smo se odločili po pregledu podatkov in analizi ter primerjavi vseh treh spremenljivk pri uporabi urnih in dnevnih vrednosti (Šraj & Bezak, 2013). Ker smo za analize potrebovali tudi volumne visokovodnih valov, smo najprej izločili bazni odtok z uporabo grafične tri-točkovne metode. Nato smo določili vzorec letnih maksimumov (Bezak, 2012), ki so ga sestavljale maksimalne vrednosti pretokov in pripadajočih volumnov posameznih let. 3. Metode 3.1 Univariatne verjetnostne analize Vzorec smo določili po metodi letnih maksimumov (Bezak et al., 2013). Za vsako spremenljivko posebej (konice pretokov in pripadajoči volumni) smo najprej naredili univariatne verjetnostne analize (Šraj et al., 2012). Uporabili smo različne porazdelitvene funkcije, ki se v svetovni hidrološki praksi najpogosteje uporabljajo (logaritemsko normalno porazdelitev, Pearsonovo III porazdelitev, logaritemsko Pearsonovo III porazdelitev, Gumbelovo porazdelitev, normalno porazdelitev, generalizirano porazdelitev ekstremnih vrednosti (GEV) in generalizirano logistično porazdelitev (GL)) in tri metode ocenj evanj a parametrov (metodo momentov, metodo momentov L in metodo največjega verjetja) (Bezak, 2012; Turk, 2012; Šraj et al., 2012). Z uporabo statističnih testov, kriterijev ujemanja in grafičnih prikazov ujemanja smo določili najprimernejšo porazdelitev za izvedbo verjetnostnih analiz konic in volumnov (Šraj et al., 2012). 3.2 Funkcija kopula Funkcij a kopula predstavlj a alternativo klasičnemu, univariatnemu pristopu k verj etno stnim analizam. Po stopek analize j e sestavljen iz naslednjih korakov, ki jih lahko združimo v dve skupini: • Oblikovanje vzorca letnih maksimumov ali metode vrednosti nad izbranim pragom, izločanje baznega odtoka, izvedba univariatnih verjetnostnih analiz in izbira najustreznejše porazdelitvene funkcije (Bezak, 2012; Šraj et al., 2012). • Določitev robnih porazdelitvenih funkcij, ocena odvisnosti obravnavanih spremenljivk, pregled nabora funkcij in izbira ustreznih kopul, določitev parametrov izbranih funkcij, uporaba statističnih testov, kriterijev ujemanja in grafičnih prikazov ujemanja za določitev najprimernejše kopule ter izračun različnih skupnih ter pogojnih (angl. conditional) povratnih dob. Prva skupina korakov je identična tistim pri univariatnih verj etnostnih analizah. Zapišemo lahko, da je izvedba klasičnih verjetnostnih analiz, kjer upoštevamo le eno spremenljivko, predpogoj za izvedbo verjetnostnih analiz z uporabo funkcije kopula. Robne porazdelitvene funkcije so tiste, ki najbolje opišejo spremenljivke pri univariatnih verjetnostnih analizah (Karmakar & Simonovic, 2009). Za prikaz odvisnosti spremenljivk, ki določajo visokovodne valove, lahko uporabimo grafične prikaze kot so: diagram K (Kendallov diagram), diagram Chi (Genest & Favre, 2007) ali razsevni diagram podatkov (Brilly & Šraj, 2005). Fisher in Swatzer (2001) sta podala več informacij o diagramu Chi, Genest in Boies (2003) pa o diagramu K. Za oceno parametrov kopul lahko podobno kot pri univariatnih verjetnostnih analizah uporabimo različne metode: metodo momentov (angl. method of moments) (Genest & Favre, 2007), metodo največjega verjetja (angl. maximum likelihood method) (Dupuis, 2007), pseudo metodo največjega verjetja (angl. maximum pseudo-likelihood method) (Genest & Favre, 2007). Pri metodi momentov parametre določimo s pomočjo Kendallovega ali Spearmanovega koeficienta korelacije (Genest & Favre, 2007). Ta način ocenjevanja parametrov je v praksi pogost (Genest & Favre, 2007; Karmakar & Simonovic, 2009; Zhang & Singh, 2006). Kendallov koeficient korelacije temelji na vrstnem redu podatkov in ga lahko izračunamo z naslednjim izrazom (Genest & Favre, 2007): kjer Pn označuje število parov spremenljivk (X,Y) za katere velja: (xi-x])(yi-y])>0, n pa je velikost vzorca. Na podlagi Kendallovega koeficienta korelacije lahko določimo kopule, ki so primerne za modeliranje odvisnosti med obravnavanimi spremenljivkami. Tako lahko z nekaterimi funkcijami opišemo le pozitivne vrednosti Kendallovega koeficienta korelacije, nekatere pa so primerne le za vrednosti koeficienta blizu vrednosti 0 (neodvisni spremenljivki). Razvitih je že veliko različnih tipov kopul, ki pripadajo posameznim družinam kopul, kot so: Arhimedova družina (angl. Archimedean copulas), eliptična družina (angl. Elliptical copulas), družina kopul ekstremnih vrednosti (angl. Extreme value copulas), Farlie-Gumbel-Morgensternova družina kopul (Genest & Favre, 2007) in še nekatere druge (Nelsen, 1999). V prispevku smo uporabili le kopule iz Arhimedove družine kopul (preglednica 1). Kopula Gumbel-Hougaard je primerna le za pozitivne vrednosti Kendallovega koeficienta korelacije, medtem ko lahko kopuli Frank in Clayton uporabimo tudi pri negativnih vrednostih korelacijskega koeficienta. Do negativnih vrednosti lahko pride predvsem pri hkratnem obravnavanju konic pretokov in časov trajanja visokovodnih valov. t = ■ n(n-l) rPr, - 1 (1) Preglednica 1: Osnovne lastnosti nekaterih kopul iz Arhimedove družine. Table 1: Some properties of copulas from Archimedean family . Ime / Name Cg(u, v) в e Gumbel-Hougaard Clayton (Cook-Johnson) exp {-((- ln u)9 + (- ln v)s)l/S] {u-9 + v-e - 1} -e_irl/e 1-i 1 + 2 [1,») [-1,«0\{0} Frank 1 ( (exp(-eu) - 1)(exp(-ev) - 1) -^lnU+-exp(-в) - 1- 1 --ri 1 - 9 if- в) ef -1 -dt (-от, »)\{0} Multivariatno porazdelitveno funkcijo H, ki jo določa n spremenljivk, lahko zapišemo z naslednjim izrazom (Zhang & Singh, 2006): .n \Cn-Cg\ AD = may —, l