Frakcionálne integrované procesy

1. Frakcionálne integrované procesy ACF stacionárnych procesov ARMA s rastúcim oneskorením exponenciál-ne klesá  korelácia náhodných premenných, ktoré sú od seba časovo vzdialené je štatisticky nevýznamná. V hydrologickej praxi sa však často stretávame s časovými radmi tvorenými stacionárnymi procesmi, v ktorých sú aj časovo veľmi vzdialené náhodné premenné pomerne silno korelované. Na tento jav v hydrologických časových radoch prvý raz upozornil už v roku 1951 Hurst. Časové rady s touto vlastnosťou nazývame rady s dlhou pamäťou (generujú ich stochastické procesy, ktoré nazývame procesy s dlhou pamäťou). Ich charakteristickou vlastnosťou je, že hodnoty ACFneklesajú s rastúcim oneskorením exponenciálne ale hyperbolicky. Uvažujme proces d Xt =t. Ak je d = 0, je to proces bieleho šumu. Ak je d celé kladné číslo, ide o integrovaný proces I(d), t. j. nestacionárny proces. Ak je d celé záporné číslo, je to neinvertibilný proces.

2. Ak d nie je celé číslo, potom sa operácia d = (1 - B)d nazýva frakcionálne diferencovanie. Modely s neceločíselným d sa nazývajú frakcionálne integrované procesy rádu d: (1 – B)d Xt = t kde t je proces bieleho šumu. Binomický rozvoj (1 – B)d Gamma funkcia (v systéme Mathematica Gamma[z]): Pre gamma funkciu plati: (k – d) = ((k – d – 1) (k – d – 2) …, (2 – d) (1 – d) (-d)) (– d) Autoregresná reprezentácia frakcionálne integrovaného procesu rádu d:

3. Binomický rozvoj (1 – B)-d Pretože môžeme písať: Frakcionálne integrovaný proces rádu d môžeme vyjadriť aj v tvare reprezentácie kĺzavých priemerov: Reprezentácia frakcionálne integrovaného procesu rádu d v tvare kĺzavých priemerov:

4. Stirlingova aproximácia pre veľké k: Môžeme teda písať: Autokorelačná funkcia frakcionálne integrovaného procesu rádu d: Parciálna autokorelačná funkcia frakcionálne integrovaného procesu rádu d: Pre d < 0.5 je proces stacionárny, pretože: Pre d > -0.5 je proces invertibilný, pretože:

5. Pre 0 < d < 0.5 má proces všetky autokorelácie kladné a platí: Pre -0.5 < d < 0 platí: Pre d  0.5 je proces nestacionárny a pre d  -0.5 je neinvertibilný. V tomto prípade nazývame proces perzistentný alebo proces s dlhou pamäťou (long memory process). V tomto prípade nazývame proces antiperzistentný alebo proces so strednou pamäťou.

6. 2. ACF má tvar t. j. pre k   , takže pre 0.5 > d > 0 je (k) > 0 a s rastúcim posunutím k klesá ACF monotónne a hyperbolicky k 0 (tento pokles je omnoho pomalší než pre d = 0). 3. PACF je daná vzťahom , k = 1, 2, .... Jej pokles je určený číslom k-1, ktoré nezávisí od d. Uvažujme proces ARIMA(0, d, 0). Pre | d | < 1/2 platí: 1. {Xt} je proces stacionárny a invertibilný

7. Proces ARIMA(0, d, 0) je pre: 0 < d < 0.5 proces s dlhou pamäťou (  |(k)|  ) -0.5 < d < 0 proces so strednou pamäťou (  |(k)|  konšt.) d  0.5 nestacionárny proces d  -0.5 neinvertibilný proces d < 0.5 stacionárny proces d > -0.5 invertibilný proces Kombináciou frakcionálneho diferencovania a procesu ARMA(p, q) sa získa trieda procesov ARFIMA(p, d, q). Dôležitou vlastnosťou týchto procesov je, že vplyv parametra d na vzdialené premenné s rastúcim oneskorením klesá hyperbolicky, ale vplyv AR a MA para-metrov klesá exponenciálne; d teda charakterizuje korelačnú štruk-túru pre veľké oneskorenia, ostatné parametre charakterizujú kore-lačnú štruktúru v malých oneskoreniach. Aj keď je ARFIMA(p, d, q) pre 0.5 < d < 1 nestacionárny proces, vygenerované časové rady sú pritahované k strednej hodnote procesu.

8. Korelačná funkcia Parciálna korelačná funkcia Príklad 1: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;-0,4;0) v tvare (1 – B)-0,4 Xt = t. Časový rad nemá žiadnu vývojovú tendenciu z hľadiska strednej hodnoty ani rozptylu, môžeme ho teda považovať za stacionárny. Odhad autokorelačnej funkcie je štatisticky nevýznamne rôzny od 0 už pre posunutie k = 2. Súčet absolútnych hodnôt autokorelácií generujúceho procesu je konečné číslo. Ide teda o proces so strednou pamäťou.

9. Korelačná funkcia Parciálna korelačná funkcia Príklad 2: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;0,45;0) v tvare (1 – B)0,45 Xt = t. Ani tento časový rad nemá žiadnu vývojovú tendenciu z hľadiska strednej hodnoty ani rozptylu, môžeme ho teda považovať za stacionárny. Hodnoty autokorelačnej funkcie klesajú k 0 veľmi pomaly a hyperbolicky, dá sa teda predpokladať, že súčet absolútnych hodnôt autokorelácií generujúceho procesu konverguje do nekonečna. Ide teda o proces s dlhou pamäťou.

10. Korelačná funkcia Parciálna korelačná funkcia Príklad 3: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;0,75;0) v tvare (1 – B)0,75 Xt = t. Tento časový rad má už jasnú vývojovú tendenciu, nemôžeme ho teda považovať za stacionárny. Hodnoty autokorelačnej funkcie klesajú k 0 veľmi pomaly, ale už nie hyperbolicky (skôr lineárne).

11. Testovanie hypotézy pre modely s dlhou pamäťou Uvažujme proces Xt =  + t, t = 1, …, n kde  je ľubovoľný parameter t je stochastický proces s nulovou strednou hodnotou Nulová hypotéza H0: Xt je proces s krátkou pamäťou Ak proces t spĺňa nasledujúce 4 podmienky: 1. E(t) = 0 2. 0 < D(t) <  3. E(t)k<  pre k > 2 4. Korelácia medzi vzdialenými premennými je štatisticky nevýznam-ná potom H0 pre daný proces nemôžeme zamietnuť. Zamietnutie H0 ešte neznamená, že proces je s dlhou pamäťou, ale len to, že všetky 4 predchádzajúce podmienky nie sú splnené súčasne.

12. kde a S(m) je aritmetický priemer a smerodajná odchýlka časového radu Xt, t = 1, …, m, 1  m  n (odporúčaná hodnota je m = n/4, resp. m = n/2). Používajú sa dve testovacie štatistiky. Pri prvej sa najprv vypočíta výraz: Potom sa vypočíta testovacia štatistika Vm: Ak vypočítaná hodnota testovacej štatistiky patrí do intervalu (0.809; 1.862), na hladine významnosti  = 0.05 nemôžeme zamietnuť nulovú hypotézu, že proces má len krátku pamäť. Ak je však hodnota mimo tohto intervalu, neznamená to ešte jednoznačne, že proces má dlhú pamäť. Potvrdíme (alebo vyvrátime) tento predpoklad použitím ďalšej testovacej štatistiky.

13. kde je odhad výberového rozptylu a autokovariancie s posunutím j časového radu Xt, t = 1, …, m, 1  m  n. Opäť vypočítame hodnotu normovanej testovacej štatistiky Vm(q) = . Ak je vypočítaná hodnota testovacej štatistiky mimo intervalu (0.809; 1.862), na hladine významnosti  = 0.05 zamietame nulovú hypotézu, že proces má len krátku pamäť. Najprv vypočítame výraz (odporúčané hodnoty q pre n > 125 sú n/10 a n/5):

14. Príklad 1: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;-0,4;0) v tvare (1 – B)-0,4 Xt = t. Prvá testovacia štatistika: Hodnota testovacej štatistiky je mimo (0.809; 1.862). Vypočítame preto aj druhú testovaciu štatistiku. Druhá testovacia štatistika: Vypočítané hodnoty testovacej štatistiky sú (až na jednu) mimo (0.809; 1.862). Na hladine významnosti  = 0.05 zamietame nulovú hypotézu, že proces má len krátku pamäť.

15. Príklad 2: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;0,45;0) v tvare (1 – B)0,45 Xt = t. Prvá testovacia štatistika: Hodnota testovacej štatistiky je mimo (0.809; 1.862). Vypočítame preto aj druhú testovaciu štatistiku. Druhá testovacia štatistika: Vypočítané hodnoty testovacej štatistiky sú mimo (0.809; 1.862). Na hladine významnosti  = 0.05 zamietame nulovú hypotézu, že proces má len krátku pamäť.

16. Príklad 3: Časový rad dĺžky 3000 simulovaný na základe procesu ARFIMA(0;0,75;0) v tvare (1 – B)0,75 Xt = t. Prvá testovacia štatistika: Vypočítaná hodnota testovacej štatistiky je mimo (0.809; 1.862). Neznamená to jednoznačne, že proces má dlhú pamäť. Potvrdíme (alebo vyvrátime) tento predpoklad použitím druhej testovacej štatistiky. Druhá testovacia štatistika: Dve hodnoty testovacej štatistiky sú mimo a dve v (0.809; 1.862). Na hladine významnosti  = 0.05 nevieme rozhodnúť, či má proces len krátku pamäť.

17. 2. Vypočítame priemernú hodnotu pre každý mesiac: 3. Vypočítame smerodajnú odchýlku pre jednotlivé mesiace: Postup pri testovaní mesačného hydrologického časového radu: Tieto časové rady obsahujú výraznú sezónnu zložku (s periódou 12 mesiacov), ktorú odstránime nasledovne: 1. Upravíme časový rad do tvaru matice 12 x N (N - počet rokov) s prvkami Xij, i = 1, …, 12, j = 1, …, N.

18. 4. Vypočítame „normované“ hodnoty: 5. Označme Yt = Mij, i = 1, …, 12 pre každé j = 1, …, N. Tento časový rad sa nazýva štandardizovaný časový rad. Testovacie štatistiky Vm a Vm(q) sa počítajú pre štandardizovaný časový rad. 6. Najprv sa vypočíta testovacia štatistika Vm (odporúčané hodnoty pre m sú n/4 a n/2). Ak je hodnota Vm vo vnútri intervalu [0.809, 1.862], nulovú hypotézu na hladine významnosti  = 0.05 nezamie-tame  časový rad je generovaný procesom s krátkou pamäťou. V tomto prípade už nie je nutné počítať aj druhu testovaciu štatistiku Vm(q). 7. Ak je hodnota Vm mimo intervalu [0.809, 1.862], nulovú hypotézu ešte nemôžeme zamietnuť. To, že časový rad nie je generovaný procesom len s krátkou pamäťou treba overiť vypočítaním testovacej štatistiky Vm(q) (odporúčané hodnoty q pre n > 125 je n/10 a n/5). Ak je aj Vm(q) mimo intervalu [0.809, 1.862], na hladine významnosti  = 0.05 môžeme zamietnuť nulovú hypotézu  časový rad nie je generovaný procesom len s krátkou pamäťou.

19. Odhad parametrov modelov ARFIMA Odhad parametrov pri modeloch ARFIMA je podstatne komplikova-nejší ako pri modeloch ARMA. Neceločíselný parameter d má v porovnaní s ostatnými parametrami špecifický význam, pretože podľa jeho hodnoty je možné určiť, či ide o proces s dlhou alebo krátkou pamäťou, stacionárny alebo nestacionárny. Existujú dve skupiny metód na odhad parametrov modelov ARFIMA. Prvú skupinu tvoria dvojkrokové metódy: najprv sa odhadne parame-ter d a potom sa odhadnú parametre modelu ARMA pre transformova-ný časový rad (filtrovaný frakcionálnym diferencovaním pre vypočíta-ný parameter d). Tieto metódy sú použiteľné len pre dostatočne dlhé časové rady. Druhú skupinu tvoria jednokrokové metódy, v ktorých sa súčasne odhadnú všetky parametre modelu ARFIMA. Väčšina z nich je založená na rôznych variantoch metódy maximálnej vierohodnosti či už v časovej alebo vo frekvenčnej oblasti. Hlavným nedostatkom týchto metód je ich výpočtová zložitosť, nutnosť poznať počiatočné hodnoty pre jednotlivé parametre a možná existencia lokálneho extrému funkcie vierohodnosti.

20. 1. Dvojkrokové metódy a. Hurstov koeficient Najprv musíme rozdeliť časový rad do m spojitých úsekov dĺžky N, pričom m x N = n (celková dĺžka časového radu). Pre každý tento j-ty úsek (j = 1, …, m): 1. Vypočítame aritmetický priemer (označíme ho Ej) a smerodajnú odchýlku Sj 2. Centrujeme dáta odčítaním aritmetického priemeru Yi, j = Xi, j - Ej pre i = 1, …, N 3. Vytvoríme kumulatívny časový rad postupným sčítaním po sebe idúcich centrovaných premenných:

21. 9. Hurstov koeficient H, pre ktorý platí H = d + 0.5 odhadneme ako regresný parameter v lineárnej regresii, kde závisle premenná je a nezávisle premenná je ln(N): 4. Určíme rozsah kumulatívneho časového radu odčítaním minimálnej hodnoty od maximálnej, čím získame štatistiku Rj: Rj = max (Z1, j, …, ZN, j) - min (Z1, j, …, ZN, j) 5. Vydelíme rozsah Rj smerodajnou odchýlkou Sj (reškálovaný rozsah) 7. Vypočítame priemer reškálovaných rozsahov pre všetky úseky dĺžky N: 8. Dĺžku úseku N zväčšíme na najbližšie číslo, pre ktoré existuje také celé číslo m, že m x N = n. Kroky 1 - 7 opakujeme, pokiaľ N  n/2. Aby sme dostali korektný odhad Hurstovho koeficientu H, musí byť v regresii minimálne 10 bodov.

22. b. Semiparametrický odhad d Táto metóda je založená na odhade parametra d vo frekvenčnej oblasti. Vychádza z teórie lineárnych filtrov, ktorá umožní vyjadriť proces (1 - B)d Xt = t v tvare fX() = | 1 - e-i | -2d f(), kde fX() a f() sú spektrálne hustoty procesov Xt a t. Po úprave môžeme písať (v spojitom tvare): v diskrétnom tvare: Odhad parametra d sa získa z lineárnej regresie vychádzajúcej z predchádzajúcej rovnice pre rôzne frekvencie 1, 2, …, N:

23. Predpokladá sa, že uj sú nezávislé rovnako rozdelené náhodné premenné s nulovou strednou hodnotou a rozptylom 2/6. Ak je {uj} proces bieleho šumu, regresiou sa získa dobrý odhad parametra d. Ak sú náhodné premenné autokorelované, platí vzťah regresie len pre frekvencie blízke 0  len v tomto prípade je odhad d konzistentný. Základnou otázkou v prípade korelovaných uj je teda určenie čísla N (napr. N , resp. sa volí také N, aby bol reziduálny rozptyl regresie rovný približne 2/6). Na testovanie parametra d sa používa štandardný t-test (jeden z výstupov príkazu Regress v systéme Mathematica). Napriek tomu, že semiparametrický odhad parametra d je veľmi jednoduchý, nejasnosti spojené s určením čísla N pri silnej autoko-relácii procesu {uj} znižujú jeho atraktívnosť. Pri nevhodnej voľbe N môžeme totiž získať veľmi vychýlený odhad. Preto sa stále častejšie v praxi používajú jednokrokové metódy.

24. 2. Jednokrokové metódy Metóda maximálnej vierohodnosti Táto metóda je určená k odhadu parametrov modelu ARFIMA v tvare p(B) (1 - B)d Xt = q(B) t za predpokladu, že t sú nezávislé náhodné premenné s rovnakým normálnym rozdelením pravdepodobnosti. Potom logaritmus presnej vierohodnostnej funkcie má tvar: kde  je vektor parametrov modelu AR,  je vektor parametrov modelu MA, n je autokovariančná matica procesu {Xt} typu n x n a X je stĺpcový vektor náhodných premenných typu 1 x n. Funkciu vierohodnosti môžeme prepísať pomocou autokorelačnej matice Pn typu n x n do tvaru:

25. Deriváciou podľa 2 získame odhad tohoto parametra Potom už môžeme uvažovať funkciu vierohodnosti bez parametra 2 Maximalizáciou tejto funkcie získame maximálne vierohodné odhady parametrov d, , . Hlavným problémom pri tejto metóde je výpočet autokorelačnej matice Pn a jej invertovanie. Boli vyvinuté mnohé iteračné metódy, ktoré okrem odhadov počítajú aj ich smerodajné chyby, potrebné k t-testom.