UNIVERZITA KARLOVA V PRAZE Matematicko fyzikální fakulta

UNIVERZITA KARLOVA V PRAZE Matematicko–fyzik´ aln´ı fakulta

´ PRACE ´ DIPLOMOVA

Tom´ aˇ s Senft Predikce popt´ avky po obˇ eˇ zivu v ekonomice z hlediska centr´ aln´ı banky

Katedra pravdˇepodobnosti a matematické statistiky Vedouc´ı diplomové práce: RNDr. Michael Koˇ na´k, PhD. Studijn´ı program: Matematika, Finanˇcn´ı a pojistná matematika

Na tomto m´ıstˇe bych rád podˇekoval RNDr. Michaelu Koˇ na´kovi, PhD. za veden´ı diplomové práce, za cenné pˇripom´ınky a za ˇcas stráven´ y na konzultac´ıch.

Prohlaˇsuji, ˇze jsem svou diplomovou práci napsal samostatnˇe a v´ yhradnˇe s pouˇzit´ım citovan´ ych pramen˚ u. Souhlas´ım se zap˚ ujˇcován´ım práce. V Praze dne 9. dubna 2006 Tomáˇs Senft

Obsah ´ 1 Uvod

6

2 Modelov´ an´ı objemu obˇ eˇ ziva 8 2.1 Funkce penˇez . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2 Centráln´ı bankovnictv´ı . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.3 Modelován´ı objemu obˇeˇziva . . . . . . . . . . . . . . . . . . . . . . 10 3 ARIMA model 3.1 Základy BJ metodologie . . . . . . . 3.2 ARIMA model poptávky po obˇeˇzivu 3.2.1 Stabilizace ˇrady . . . . . . . . 3.2.2 Struktura modelu . . . . . . . 3.2.3 Odhad parametr˚ u modelu . . 3.3 Verifikace modelu . . . . . . . . . . . 3.4 Predikce . . . . . . . . . . . . . . . . 4 GARCH model 4.1 Základn´ı principy . . . . . . . . . . . 4.2 GARCH model poptávky po obˇeˇzivu 4.2.1 Stabilizace ˇrady . . . . . . . . 4.2.2 Struktura modelu . . . . . . . 4.2.3 Odhad parametr˚ u modelu . . 4.3 Verifikace modelu . . . . . . . . . . . 4.4 Predikce . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . .

5 STS model 5.1 Teorie stavov´ ych model˚ u ˇcasov´ ych ˇrad 5.1.1 Základn´ı model . . . . . . . . . 5.1.2 Kalmanovy rekurze . . . . . . . 5.1.3 Inicializace Kalmanova filtru . . 5.1.4 Odhad parametr˚ u modelu . . . 5.2 STS model poptávky po obˇeˇzivu . . . . 5.2.1 Modelované sezónnosti . . . . . 5.2.2 Konstrukce kubick´ ych splin˚ u. .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . .

. . . . . . .

12 12 13 13 14 15 17 20

. . . . . . .

23 23 25 25 26 27 27 29

. . . . . . . .

32 32 32 33 36 39 41 42 44

5.2.3 Alternativn´ı modelován´ı mˇes´ıˇcn´ı sezónnosti . . . . . . . . . . 47 5.3 Verifikace modelu . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 5.4 Predikce . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 6 Kombinace model˚ u

54

7 Porovn´ an´ı model˚ u

57

8 Z´ avˇ er

60

Literatura

61

4

N´ azev pr´ ace: Predikce poptávky po obˇeˇzivu v ekonomice z hlediska centráln´ı banky Autor: Tomáˇs Senft Katedra: Katedra pravdˇepodobnosti a matematické statistiky Vedouc´ı diplomov´ e pr´ ace: RNDr. Michael Koˇ na´k, PhD. e-mail vedouc´ıho: [email protected] Abstrakt: Tato diplomová práce se zab´ yvá modelován´ım a predikc´ı objemu obˇeˇziva, které patˇr´ı mezi hlavn´ı autonomn´ı veliˇciny ovlivˇ nuj´ıc´ı likviditu trhu. V práci jsou objasnˇeny potˇreby jeho modelován´ı a prezentovány tˇri zkonstruované stochastické modely. Jsou jimi ARIMA a GARCH model vycházej´ıc´ı z Box-Jenkinsovy metodologie a STS model. STS model je strukturovan´ y model ˇcasové ˇrady vyuˇz´ıvaj´ıc´ı Kalmanovy rekurze. Pˇredpovˇedi jednotliv´ ych model˚ u jsou dále kombinovány a statisticky porovnány. V´ ysledky ukazuj´ı, ˇze nejvhodnˇejˇs´ım modelem pro pˇredpov´ıdán´ı objemu obˇeˇziva je kombinace STS a ARIMA modelu, kter´ y dosáhl stejné ˇ kvality pˇredpovˇed´ı jako expertn´ı odhad pouˇz´ıvan´ y v CNB. Lze jej tedy pouˇz´ıt ˇ pˇrinejmenˇs´ım jako podp˚ urn´ y prostˇredek pro ˇr´ızen´ı likvidity v CNB. Kl´ıˇ cov´ a slova: obˇeˇzivo, ARIMA, GARCH, STS model, Kalman˚ uv filtr

Title: Forecasting of the Demand of Currency in Circulation from the View of Central Bank Author: Tomáˇs Senft Department: Department of Probability and Mathematical Statistics Supervisor: RNDr. Michael Koˇ na´k, PhD. Supervisor’s e-mail address: [email protected] Abstract: This diploma thesis deals with modeling and forecasting of the daily series of currency in circulation, which is one of the main autonomous factors influencing the liquidity of financial markets. Reasons for its modeling are explained and three constructed stochastic models are presented. There are ARIMA and GARCH models based on Box-Jenkins methodology and STS model. STS model is structured time series model using Kalman equations. Forecasts of models are combined together and statistically compared. The results show that the combination of STS and ARIMA models is the best model for forecasting of the daily series of currency in circulation and it has the same forecasting performance as the current model-judgement practice in the Czech National Bank. The model might be also applied at least as a supportive tool for the liquidity management. Key words: currency in circulation, ARIMA, GARCH, STS model, Kalman filter

Kapitola 1 ´ Uvod Jedn´ım z hlavn´ıch u ´kol˚ u centráln´ı banky v ekonomickém prostˇred´ı je dosaˇzen´ı a udrˇzen´ı cenové stability, tj. vytváˇren´ı n´ızkoinflaˇcn´ıho prostˇred´ı v ekonomice. Pˇri ˇ a národn´ı banka (CNB) ˇ plnˇen´ı svého c´ıle vyuˇz´ıvá Cesk´ nˇekolik mˇenovˇepolitick´ ych nástroj˚ u. Efektivita vyuˇzit´ı nástroj˚ u mˇenové politiky závis´ı na kvalitˇe odhadu v´ yvoje poptávky a nab´ıdky po penˇez´ıch (tzv. likvidita trhu). Celková likvidita trhu se skládá z nˇekolika komponent, které lze rozdˇelit na dvˇe skupiny podle moˇznost´ı jejich kontroly. Prvn´ı skupinu tvoˇr´ı veliˇciny, které jsou pˇresnˇe pˇredpovˇeditelné, respektive pˇredem známé, napˇr´ıklad vypl´ yvaj´ı z jiˇz sjednan´ ych obchod˚ u. Druhou skupinu tvoˇr´ı takzvané autonomn´ı veliˇciny, které nejsou pod pˇr´ımou kontrolou centráln´ı banky. Jednou z nejh˚ uˇre predikovateln´ ych autonomn´ıch veliˇcin je právˇe objem obˇeˇziva, tedy celková hodnota bankovek a minc´ı v drˇzen´ı komerˇcn´ıch bank a nebankovn´ıch subjekt˚ u. Poptávka po obˇeˇzivu se mˇen´ı v závislosti na celé ˇradˇe faktor˚ u. Zpravidla stoupá pˇred v´ yplatami a klesá koncem mˇes´ıce. Naopak v pr˚ ubˇehu t´ ydne poptávka vzr˚ ustá a nejvˇetˇs´ı b´ yvá pˇred v´ıkendem. Tento efekt b´ yvá umocnˇen pokud je v pátek nebo v pondˇel´ı svátek. Faktor˚ u ovlivˇ nuj´ıc´ıch poptávku po obˇeˇzivu je vˇsak mnohem v´ıce, napˇr´ıklad turistick´ y ruch, podm´ınky pro v´ ybˇer z bankomat˚ u atd. V souˇcasné dobˇe ˇ jsou denn´ı predikce objem˚ u obˇeˇziva v CNB realizovány pomoc´ı tzv. expertn´ıch odhad˚ u. Jedná se o expertn´ı sledován´ı pˇredchoz´ıch hodnot a zmˇen v mˇenové bázi a následné vyhodnocen´ı. C´ılem práce je vyˇreˇsen´ı tohoto u ´kolu aplikac´ı ekonometrick´ ych a statistick´ ych teori´ı. Jedná se pˇredevˇs´ım o podrobnou anal´ yzu ˇcasové ˇrady obˇeˇziva, nalezen´ı relevantn´ıch exogenn´ıch promˇenn´ ych a následnˇe vytvoˇren´ı vhodného stochastického modelu v´ yvoje objemu obˇeˇziva pˇredevˇs´ım s d˚ urazem na jeho pˇredpov´ıdac´ı schopnost. V kapitole 2 práce je struˇcnˇe shrnuta problematika penˇez, centráln´ıho bankovnictv´ı a d˚ uvody modelován´ı objemu obˇeˇziva. V kapitolách 3, 4 a 5 jsou prezentovány sestavené stochastické modely objemu

6

obˇeˇziva. U kaˇzdého modelu jsou vysvˇetleny teoretické základy a je popsáno jeho sestaven´ı, vˇcetnˇe odhadu parametr˚ u. Následnˇe je provedena verifikace a testován´ı model˚ u. Na konci kaˇzdé kapitoly je ˇca´st vˇenovaná predikci pomoc´ı vytvoˇreného modelu a jsou shrnuty obdrˇzené v´ ysledky. Kapitola 6 je vˇenována kombinovan´ ym model˚ um, které d´ıky tomu, ˇze vyuˇz´ıvaj´ı informace obsaˇzené v d´ılˇc´ıch modelech, dávaj´ı lepˇs´ı v´ ysledky pˇredpovˇed´ı. Porovnán´ı model˚ u, vˇcetnˇe k tomu pouˇzit´ ych kritéri´ı, je popsáno v kapitole 7, kde je zároveˇ n vybrán fináln´ı model. Pro pˇr´ıpravu dat byl pouˇzit program Microsoft Office Excel 2002. Model ARIMA byl zpracován v programu Matlab 6.5, pˇredevˇs´ım s pouˇzit´ım aplikaˇcn´ıho bal´ıˇcku System Identifacation Toolbox. Pro GARCH model byl pouˇzit bal´ıˇcek GARCH Toolbox ze stejného programu. Nejvˇetˇs´ı ˇca´st práce – STS model byl naprogramován v systému Mathematica 5.0 s pouˇzit´ım aplikaˇcn´ıho bal´ıˇcku Time Series. Kombinován´ı model˚ u a jejich porovnán´ı je zpracováno taktéˇz v systému Mathematica 5.0. Vˇsechna data, v´ ysledky a pouˇzité skripty jsou uloˇzené v CD pˇr´ıloze.

7

Kapitola 2 Modelov´ an´ı objemu obˇ eˇ ziva Obsah této kapitoly je vˇenován vysvˇetlen´ı problematiky penˇez, centráln´ıho bankovnictv´ı a modelován´ı poptávky po obˇeˇzivu. Podrobné a aktuáln´ı informace ˇ lze nalézt na internetov´ ych stránkách CNB.

2.1

Funkce penˇ ez

Pen´ıze pˇredstavuj´ı vˇseobecnˇe pˇrij´ıman´ y prostˇredek smˇeny. Za pen´ıze m˚ uˇzeme tedy povaˇzovat vˇse, co je v urˇcitém ˇcase a m´ıstˇe vyuˇz´ıváno ekonomick´ ymi subjekty pˇri smˇenˇe statk˚ u. Dále maj´ı funkci z´ uˇctovac´ı jednotky, tj. lze jimi vyjádˇrit cenu ostatn´ıch statk˚ u, a slouˇz´ı jako uchovatel hodnoty. Nejstarˇs´ı pen´ıze se objevily zhruba pˇred jedenácti tis´ıci lety, kdy funkci penˇez zastával dobytek a zemˇedˇelské plodiny. Forma penˇez se bˇehem historie postupnˇe vyv´ıjela. Dneˇsn´ı pen´ıze jsou penˇezi nekryt´ ymi, neboli penˇezi s nucen´ ym obˇehem, které jsou pˇrij´ımány proto, ˇze jsou zákonem pˇredepsan´ ym platidlem. Maj´ı dvˇe základn´ı formy - hotovostn´ı a pˇrevládaj´ıc´ı bezhotovostn´ı. Pˇri tomto dˇelen´ı je nutné si uvˇedomit jejich velmi u ´zkou vzájemnou vazbu. V podstatˇe lze hotovostn´ı pen´ıze uloˇzen´ım do banky pˇremˇenit na bezhotovostn´ı a naopak. Obˇe formy jsou také emitovány v zásadˇe stejnˇe, pˇredevˇs´ım u ´vˇerem. Jedin´ ym zásadn´ım rozd´ılem je, ˇze hotovostn´ı pen´ıze, aˇz na drobné vyj´ımky v nˇekter´ ych zem´ıch, m˚ uˇze emitovat pouze centráln´ı banka. Bezhotovostn´ı pen´ıze oproti tomu mohou emitovat i jiné instituce, napˇr´ıklad obchodn´ı banky.

8

2.2

Centr´ aln´ı bankovnictv´ı

Pod pojmem centráln´ı bankovnictv´ı rozum´ıme dvoustupˇ nov´ y bankovn´ı systém, kde prvn´ım stupnˇem je centráln´ı banka, tzv. banka bank, a druh´ ym stupnˇem obchodn´ı a dalˇs´ı banky, které slouˇz´ı ostatn´ım subjekt˚ um a podnikaj´ı na rozd´ıl ˇ e republice je podle Ustavy ´ od centráln´ı banky za u ´ˇcelem zisku. V Cesk´ centráln´ı ˇ bankou CNB, v Evropské mˇenové unii pln´ı tuto u ´lohu Evropská centráln´ı banka ˇ (ECB). Hlavn´ım c´ılem CNB je cenov´ a stabilita, která je d˚ uleˇzitá pro správn´ y chod ˇ ekonomiky. Podobnˇe jako ostatn´ı centráln´ı banky se i CNB zamˇeˇruje pˇredevˇs´ım na stabilitu spotˇrebitelsk´ ych cen. V praxi se stabilitou nerozum´ı nemˇennost cen, ale jejich m´ırn´ y r˚ ust. Pˇri plnˇen´ı své základn´ı mˇenovˇepolitické u ´lohy, tedy cenové stabiˇ lity, vol´ı centráln´ı banka jeden z nˇekolika moˇzn´ ych mˇenovˇepolitick´ ych reˇzim˚ u. CNB pouˇz´ıvá poˇc´ınaje rokem 1998 reˇzim c´ılován´ı inflace, proto jsou také v naˇsem modelu pouˇzita data aˇz od tohoto roku. Hlavn´ımi rysy c´ılován´ı inflace jsou stˇrednˇedobost této strategie, vyuˇz´ıván´ı prognózy inflace a veˇrejné explicitn´ı vyhláˇsen´ı inflaˇcn´ıho ˇ c´ıle, ˇci posloupnosti c´ıl˚ u. Pro dosaˇzen´ı c´ıl˚ u má CNB k dispozici nˇekolik mˇenovˇepolitick´ ych nástroj˚ u, kter´ ymi jsou operace na volném trhu, automatické facility a repo sazba a povinné minimáln´ı rezervy. Pro téma této diplomové práce je d˚ uleˇzit´ y prvn´ı z nich. O operaci na volném trhu hovoˇr´ıme, je-li jednou ze dvou obchoduj´ıc´ıch stran centráln´ı banka a druhou domác´ı bankovn´ı subjekt. Operace na volném trhu jsou ˇ vˇetˇsinou provádˇeny ve formˇe tzv. repo operac´ı, pˇri kter´ ych CNB pˇrij´ımá od bank 1 pˇrebyteˇcnou likviditu a bankám pˇredává dohodnuté cenné pap´ıry. Obˇe strany se zároveˇ n zavazuj´ı, ˇze po uplynut´ı doby splatnosti probˇehne zpˇetná transakce, v n´ıˇz ˇ CNB jako dluˇzn´ık vrát´ı vˇeˇritelské bance zap˚ ujˇcenou jistinu zv´ yˇsenou o dohodnut´ y ˇ u ´rok a vˇeˇritelská banka vrát´ı CNB poskytnuté cenné pap´ıry. Tento hlavn´ı mˇenov´ y ˇ nástroj se provád´ı ve formˇe tendr˚ u2 , kde CNB pˇrij´ımá nab´ıdky aˇz do v´ yˇse predikovaného pˇrebytku likvidity. Horizont predikce je jeden den, coˇz souvis´ı s t´ım, ˇze se vyhláˇsen´ı repo tendru provád´ı kaˇzd´ y pracovn´ı den. 1

Obecnˇe mohou centráln´ı banky v reˇzimu c´ılován´ı inflace pen´ıze z obˇehu nejen stahovat, ale ˇ také do obˇehu dodávat. Stav, kdy pro udrˇzen´ı u ´rokov´ ych sazeb potˇrebuje CNB likviditu stahovat, ˇ je dán v´ yvojem poˇcátku 90. let, kdy CNB pustila do ekonomiky mnoˇzstv´ı penˇez za u ´ˇcelem udrˇzen´ı stability devizového kurzu 2 Tendr je zp˚ usob operace, kdy banky neznaj´ı podm´ınky (objem, popˇr. ani cenu) nab´ızen´ ych, resp. poptávan´ ych cenn´ ych pap´ır˚ u dalˇs´ıch bank a nemohou tak vlastn´ı nab´ıdky upravovat. Jedná ˇ se o tzv. obálkovou metodu. Nab´ıdky bank vypoˇrádá CNB podle americké aukˇcn´ı procedury, tj. pˇrijme pˇrednostnˇe nab´ıdky poˇzaduj´ıc´ı nejniˇzˇs´ı u ´rokovou sazbu a to aˇz do v´ yˇse predikovaného pˇrebytku likvidity na dan´ y den. V pˇr´ıpadˇe, ˇze objem objednan´ y bankami pˇresáhne predikovan´ y ˇ pˇrebytek likvidity, CNB nab´ıdky za nejvyˇsˇs´ı sazby bud’ zcela odm´ıtne, nebo proporcionálnˇe zkr´ at´ı.

9

2.3

Modelov´ an´ı objemu obˇ eˇ ziva

V souvislosti s ˇr´ızen´ım likvidity se setkáváme s pojmem mˇenová báze. Mˇenová báze je tvoˇrena hotovostn´ımi penˇezi a bezhotovostn´ımi penˇezi. Jej´ı zdroje se dˇel´ı na autonomn´ı faktory a mˇenovˇepolitické faktory. Mezi autonomn´ı faktory patˇr´ı napˇr´ıklad vnˇejˇs´ı sektor, vládn´ı sektor a pohledávky za nebankovn´ımi subjekty. Mezi mˇenovˇepolitické faktory ˇrad´ıme závazky z repo operac´ı, automatické facility a krátkodobé u ´vˇery na zachován´ı likvidity. Jak je ˇr´ızena likvidita? Komerˇcn´ı ˇ banky maj´ı u CNB tzv. clearingové u ´ˇcty (CU), jejichˇz jednu ˇca´st tvoˇr´ı pevnˇe urˇcená v´ yˇse povinn´ ych minimáln´ıch rezerv. Mnoˇzstv´ı penˇez pˇresahuj´ıc´ı tuto hodnotu, které nen´ı nijak u ´roˇceno, naz´ yváme pˇrebytek likvidity. Pˇrebytek likvidity ˇ se do ekonomiky m˚ uˇze dostat napˇr. zv´ yˇsen´ım devizov´ ych rezerv CNB nebo státn´ı emis´ı. V ideáln´ım pˇr´ıpadˇe by byl pˇrebytek likvidity nulov´ y, coˇz ale v praxi neplat´ı. Proto je z d˚ uvodu udrˇzen´ı krátkodob´ ych u ´rokov´ ych sazeb tˇreba tuto pˇrebyteˇcnou likviditu stáhnout, coˇz se ˇreˇs´ı pomoc´ı repo tendr˚ u zm´ınˇen´ ych v pˇredchoz´ı kapitole. Objem obˇeˇziva je nejvˇetˇs´ı autonomn´ı veliˇcinou mˇenové báze a ovlivˇ nuje mnoˇzstv´ı penˇez na clearingov´ ych u ´ˇctech. Konkrétnˇe tak, ˇze o kolik se zv´ yˇs´ı objem obˇeˇziva, o tolik se sn´ıˇz´ı mnoˇzstv´ı bezhotovostn´ıch penˇez na CU a naopak. Predikce objemu obˇeˇziva, jeˇz budou poptávat komerˇcn´ı banky, je tedy nutná k urˇcen´ı v´ yˇse repo tendru na staˇzen´ı pˇrebyteˇcné likvidity. Tato poptávka vykazuje sezónn´ı charakter (Obr. 2.1) souvisej´ıc´ı napˇr´ıklad s vˇetˇs´ı poptávkou po obˇeˇzivu v obdob´ı Vánoc. Dalˇs´ı vlivy jsou popsány v následuj´ıc´ıch kapitolách. ˇ V souˇcasné dobˇe jsou denn´ı predikce objem˚ u obˇeˇziva v CNB realizovány pomoc´ı expertn´ıch odhad˚ u, tj. expertn´ıho sledován´ı pˇredchoz´ıch hodnot a zmˇen v mˇenové bázi a následného vyhodnocen´ı. V této práci je popsáno ˇreˇsen´ı tohoto u ´kolu aplikac´ı ekonometrick´ ych a statistick´ ych teori´ı, pˇriˇcemˇz denn´ı objemy obˇeˇziva jsou chápány jako prvky ˇcasové ˇrady. Jak jiˇz bylo naznaˇceno, základem jsou data z let 1998 aˇz 2003, pˇriˇcemˇz data do roku 2002 jsou pouˇzita pro konstrukci model˚ u a data za rok ’ 2003 jsou ponechána zvláˇst na testován´ı jejich pˇredpov´ıdac´ı schopnosti.

10

280

260

240

220

m

200

180

160

140

120

200

400

600

800

1000

1200

1400

d

Obrázek 2.1: Poptávka po obˇeˇzivu v letech 1998 - 2003 • m - mld. Kˇc • d - denn´ı index

11

Kapitola 3 ARIMA model Prvn´ım modelem, kter´ y je pouˇzit pro modelován´ı poptávky po obˇeˇzivu, je ARIMA model. Jeho základem je kombinace stochastické ˇcásti vycházej´ıc´ı z Box-Jenkinsovy (BJ) metodologie s deterministickou regresn´ı komponentou. Tento model dává základn´ı pˇredstavu o chován´ı ˇcasové ˇrady obˇeˇziva a poskytuje v´ ychoz´ı bod pro dalˇs´ı modely. V prvn´ı ˇcásti jsou nejdˇr´ıve struˇcnˇe nast´ınˇeny principy BJ metodologie a následnˇe popsán model. V tˇret´ı ˇca´sti je model ovˇeˇren a jsou popsány k tomu pouˇzité statistiky. Nakonec je provedena predikce pomoc´ı sestaveného modelu.

3.1

Z´ aklady BJ metodologie

Box-Jenkinsova metodologie pracuje se stacionárn´ımi náhodn´ ymi procesy. Podle Woldovy vˇety lze totiˇz kaˇzd´ y nedeterministick´ y stacionárn´ı proces {yt } rozloˇzit na ˇcistˇe náhodnou sloˇzku a deterministickou ˇca´st yt =

∞ X

ψt εt−j + Dt ,

j=0

P 2 avisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych kde ψ0 = 1, ∞ j=0 ψt < ∞, {εt } je proces nez´ 2 veliˇcin s rozdˇelen´ım N(0, σ ) a pro vˇsechna s a t plat´ı Cov(εs , Dt ) = 0, εt = Pt εt , Dt = Ps Dt a Dt jsou deterministické (Dt − Pt−1 Dt = 0 pro vˇsechny t). Pn znaˇc´ı operátor nejlepˇs´ı lineárn´ı predikce pˇri znalosti pozorován´ı yt , t ≤ n. Náhodnou sloˇzku Woldovy dekompozice oznaˇcovanou jako MA(∞) proces, lze pˇri dostateˇcnˇe

12

velk´ ych p a q dobˇre aproximovat procesem y˜t =

p X

φi y˜t−i +

i=1

q X

θj εt−j + εt ,

j=1

ekvivalentnˇe zapsáno φ(B)˜ yt = θ(B)εt , kde B je operátor zpˇetného posunut´ı. {˜ yt } se naz´ yvá ARMA(p,q) proces a je základn´ım kamenem BJ metodologie. ARMA proces má nulovou stˇredn´ı hodnotu. Podm´ınka stacionarity je, ˇze vˇsechny koˇreny polynomu φ(B) leˇz´ı vnˇe jednotkového kruhu a podm´ınka invertibility je totoˇzná pro polynom θ(B). Podrobnosti lze nalézt v [1].

3.2

ARIMA model popt´ avky po obˇ eˇ zivu

Pˇri konstrukci modelu je tˇreba nejdˇr´ıve ˇradu stabilizovat, pak nalézt vhodnou strukturu stochastické ˇca´sti a relevantn´ı vysvˇetluj´ıc´ı promˇenné a následnˇe správnˇe odhadnout jejich vliv na poptávku po obˇeˇzivu. Jednotlivé kroky jsou popsány podrobnˇeji v podkapitolách. Model vycház´ı z model˚ u popsan´ ych v [2] a [12].

3.2.1

Stabilizace ˇ rady

Pˇri pohledu na graf ˇcasové ˇrady obˇeˇziva (Obr. 2.1) je zˇrejmé, ˇze se jedná o ˇradu silnˇe nestacionárn´ı, která obsahuje nejen trend, ale i ˇradu sezónn´ıch vliv˚ u. Nestacionaritu potvrzuj´ı rovnˇeˇz korelogramy (Obr. 3.1 a 4.1). Protoˇze chceme pˇri modelován´ı vyuˇz´ıt BJ metodologie, která pracuje se stacionárn´ımi ˇradami, je tˇreba ˇradu vhodn´ ym zp˚ usobem stabilizovat. Jako nejvhodnˇejˇs´ı transformaci jsme vyhodnotili prvn´ı diferenci ytd = ∆yt = yt − yt−1 = (1 − B)yt ,

t = 2, . . . , n.

ˇ Radu ytd (Obr. 3.2) dále jiˇz neupravujeme, stále pˇr´ıtomnou nestacionaritu odstran´ıme pˇridán´ım relevantn´ıch vysvˇetluj´ıc´ıch promˇenn´ ych. Po odeˇcten´ı jejich vlivu je ˇcasová ˇrada jiˇz témˇeˇr stacionárn´ı. Dalˇs´ı moˇznou transformac´ı by bylo dodateˇcné pouˇzit´ı diference se zpoˇzdˇen´ım jednoho roku, tj. 252 dn´ı, protoˇze ˇcasová ˇrada obˇeˇziva obsahuje urˇcitou meziroˇcn´ı závislost. Pak by ale bylo nutné pouˇz´ıt sezónn´ı SARIMA model, coˇz software M atlab neumoˇznuje. Tato varianta byla zpracována v [11].

13

Pˇri testován´ı stacionarity jsme také pouˇzili Goldfeld-Quandt˚ uv test pro testován´ı konstantnosti rozptylu. V´ ysledné p-hodnoty (p13 = 0.5834 a p31 = 0.4166) hypotézu homoskedascity nezam´ıtaj´ı. Pˇresto lze pozorovat m´ırn´ y nár˚ ust variance ˇrady v pr˚ ubˇehu ˇcasu, kter´ y by se dal odstranit logaritmován´ım ˇrady, pˇr´ıpadnˇe ˇ adná z uveden´ Box-Coxovou transformac´ı. Z´ ych variant ale nepˇrinesla zlepˇsen´ı pˇredpov´ıdac´ı schopnosti modelu.

3.2.2

Struktura modelu

Zvolen´ y model poptávky po obˇeˇzivu je kombinac´ı ARMA procesu diferencované ˇcasové ˇrady a regresn´ıch promˇenn´ ych a je moˇzno jej zapsat ve tvaru φ(B)δ(B)yt = dt + θ(B)εt , P kde regresn´ı komponenta je sumou jednotliv´ ych kalendáˇrn´ıch vliv˚ u dt = ki=1 dt,i a dt,i jsou funkc´ı pevného vektoru nezávisl´ ych ˇcasovˇe závisl´ ych“ promˇenn´ ych. Poly” nomy φ(B) a θ(B) urˇcuj´ı strukturu a ˇra´d stochastické ARMA komponenty, δ(B) je diferenˇcn´ı operátor a {εt } je proces nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin s normáln´ım rozdˇelen´ım s nulovou stˇredn´ı hodnotou a rozptylem σ 2 . ARMA(p, q) proces d-krát diferencované ˇcasové ˇrady naz´ yváme ARIMA(p, d, q) proces. V naˇsem pˇr´ıpadˇe je p = 11, d = 1, q = 1 a φ(B) = 1 − φ1 B − φ3 B 3 − φ4 B 4 − φ10 B 10 − φ11 B 11 , θ(B) = 1 − θ1 B 1 , δ(B) = 1 − B. ˇ ıselné hodnoty koeficient˚ C´ u polynom˚ u jsou uvedeny v CD pˇr´ıloze. Urˇcen´ı ˇra´d˚ up a q ARMA modelu jsme provedli na základˇe pr˚ ubˇeh˚ u autokorelaˇcn´ı funkce (ACF) a parciáln´ı autokorelaˇcn´ı funkce (PACF) diferencované ˇcasové ˇrady a následn´ ym porovnán´ım v´ ysledk˚ uvu ´vahu pˇricházej´ıc´ıch model˚ u. Do modelu jsme zahrnuli pˇet r˚ uzn´ ych kalendáˇrn´ıch vliv˚ u a konstantn´ı ˇclen ˇ ıselné hodnoty pro zachycen´ı nenulové stˇredn´ı hodnoty ˇrady prvn´ıch diferenc´ı. C´ regresn´ıch koeficient˚ u w jsou opˇet uvedeny v CD pˇr´ıloze. V´ ybˇer a formulace modelovan´ ych kalendáˇrn´ıch vliv˚ u byla provedena na základˇe zkuˇsenost´ı s modelován´ım ˇ poptávky po obˇeˇzivu v CNB a ECB a dále anal´ yzou signifikantnosti parametr˚ ua struktury rezidu´ı.

ustu ˇcasové ˇrady obˇeˇziva je Konstatn´ı ˇ clen Konstatn´ı ˇclen udávaj´ıc´ı rychlost r˚ modelován regresn´ı promˇennou, která je vektorem jedniˇcek, tedy dt,1 = w1 pro vˇsechna t.

14

ˇ Vliv dn˚ u v t´ ydnu Casov´ a ˇrada poptávky po obˇeˇzivu vykazuje urˇcit´ y t´ ydenn´ı cyklus, napˇr´ıklad ve ˇctvrtek je nejvyˇsˇs´ı poptávka po obˇeˇzivu z d˚ uvodu plnˇen´ı bankomat˚ u na v´ıkend. Pro zachycen´ı tohoto vlivu jsme pouˇzili ˇctyˇri v podut st ct statˇe indikátorové promˇenné hpo ı z nich nab´ yvá hodnoty 1, t , ht , ht , ht . Prvn´ je-li pondˇel´ı, -1 pro pátek a 0 jinak. Dalˇs´ı promˇenné maj´ı analogické hodnoty pro u ´ter´ y, stˇredu a ˇctvrtek. Takto definované promˇenné zajist´ı ˇcistou t´ ydenn´ı sezonalitu, tj. aby souˇcet vliv˚ u vˇsech dn˚ u v t´ ydnu byl nulov´ y. ut ut st st ct ct dt,2 = wpo hpo t + w ht + w ht + w ht .

Vliv poˇ rad´ı dne v mˇ es´ıci Závislost hodnot ˇrady na poˇrad´ı dne v mˇes´ıci je ovlivnˇena napˇr´ıklad term´ıny vyplácen´ı mezd. Pˇri modelován´ı mˇes´ıˇcn´ıho vlivu jsme stejnˇe jako v ECB (viz [2]) pouˇzili trigonometrické funkce ve tvaru ¶ p µ X 2πjmt 2πjmt sin cos wj sin + wj cos , dt,3 = Mt Mt j=1 kde mk udává poˇrad´ı dne v mˇes´ıci a Mt je poˇcet dn´ı v daném mˇes´ıci. Hodnota p mus´ı b´ yt dostateˇcnˇe velká, aby dobˇre zachytila modelovanou sezónnost. Bylo ukázáno v [13], ˇze p = 8 je dostaˇcuj´ıc´ı. Vliv st´ atn´ıch sv´ atk˚ u a dn˚ u volna Pr˚ ubˇeh poptávky po obˇeˇzivu ovlivˇ nuj´ı také státn´ı svátky a dny volna. Jsou to Nov´ y rok, Velikonoce, Svátek práce, Den v´ıtˇezstv´ı nad faˇsismem, Den slovansk´ ych vˇerozvˇest˚ u Cyrila a Metodˇeje, Den upálen´ı mistra Jana Husa, Den ˇceské státnosti, Den vzniku samostatného ˇceskoslovenského státu, Den boje za svobodu a demokracii a Vánoce. Komponenty svátk˚ u maj´ı tvar dt = w(B)ht , kde w(B) je polynom udávaj´ıc´ı zpoˇzdˇen´ı, tj. o kolik ˇcasov´ ych index˚ u posuneme promˇennou vpˇred a vzad, a ht je indikátorová promˇenná konkrétn´ıho svátku. Tvary polynom˚ u w(B) jsou uvedeny v tabulce 3.1. Celkov´ y vliv vˇsech svátk˚ u a dn˚ u volna dostaneme seˇcten´ım d´ılˇc´ıch komponent. Vliv mimoˇ r´ adn´ ych ud´ alost´ı V historii v´ yznamnˇe ovlivnili ˇcasovou ˇradu poptávky po obˇeˇzivu dvˇe mimoˇrádné události. Jsou to krach IPB a problém Y2K (problém roku 2000, jednalo se o obavu z nefunkˇcnosti nˇekter´ ych poˇc´ıtaˇcov´ ych program˚ u pˇri pˇrechodu z 31.12.1999 na 1.1.2000). Jejich formulace v modelu je stejná jako u státn´ıch svátk˚ u a dn˚ u volna, tvary polynom˚ u w(B) jsou taktéˇz v tabulce 3.1.

3.2.3

Odhad parametr˚ u modelu

Protoˇze parametry modelu jsou neznámé, je tˇreba je odhadnout. Nejvhodnˇejˇs´ı je pouˇz´ıt metodu maximáln´ı vˇerohodnosti. Parametry ARMA sloˇzky a regresn´ı

15

8

6

0.8

4 0.6

2 ACF

m

0.4

0 0.2

−2

0

−4

−0.2

0

50

100

150 l

200

250

300

−6

200

400

600

800

1000

1200

1400

d

Obrázek 3.1: ACF ˇcasové ˇrady obˇeˇziva Obrázek 3.2: Rada ˇ prvn´ıch diferenc´ı (ˇcervenˇe), hranice v´ yznamnosti • m - mld. Kˇc • d - denn´ı index (modˇre) • l - zpoˇzdˇen´ı 0.5

2

1.5

0.4

1 0.3 0.5 0.2 m

m

0

0.1 −0.5 0 −1

−0.1

−0.2

−1.5

1

2

3 d

4

−2

5

2

4

6

8

10

12

14

16

18

20

22

24

d

(a) Vliv dn˚ u v t´ ydnu • d - index dne v t´ ydnu (b) Vliv poˇrad´ı dne v mˇes´ıci • d - index dne v mˇes´ıci 1

12

10

0.5

8

m

m

0

6

4

2

−0.5 −3

−2

−1

0

1

2

3

4

d

0 −4

−3

−2

−1

0 d

1

2

3

4

(c) Vliv Dne v´ıtˇezstv´ı nad faˇsismem • d - index (d) Vliv krachu IPB • d - index dne od krachu dne od svátku IPB

Obrázek 3.3: Pr˚ ubˇehy sezónn´ıch vliv˚ u • m - mld. Kˇc

16

státn´ı svátky a dny volna krach IPB problém Y2K

Dn˚ u pˇred Dn˚ u po 4 4 4 4 10 4

Polynom posunut´ı (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + . . . + w14 B 14 )B −10

B je operátor zpˇetného posunut´ı a wi jsou r˚ uzn´ a pro jednotlivé svátky

Tabulka 3.1: Poˇcty dn˚ u modelovan´ ych pˇred a po svátc´ıch a událostech

koeficienty se odhaduj´ı souˇcasnˇe. Vzhledem k pˇredpokladu normáln´ıho rozdˇelen´ı lze vˇerohodnostn´ı funkci vyjádˇrit jako ½ ¾ 1 2 −n/2 −1/2 0 −1 L(w, φ, θ, σ ) = (2π) | Γn | exp − (Y − Yˆ ) Γn (Y − Yˆ ) , 2 kde Γn (φ, θ, σ2 ) je kovarianˇcn´ı matice procesu W = Y − Yˆ , Y znaˇc´ı pozorován´ı ˇcasové ˇrady a Yˆ odhad ˇrady modelem. Vzhledem k jej´ı sloˇzitosti je nutné vˇerohodnostn´ı funkci maximalizovat numericky. K tomu je pouˇzita funkce armax softwaru Matlab. Po odhadu parametr˚ u modelu je ovˇeˇrováno, zda nen´ı nˇekter´ y odhadnut´ y koeficient nesignifikantn´ı. K testován´ı je pouˇzita statistika ti =

ˆbi , si

kde ˆbi je odhad i-tého parametru a si jeho smˇerodatná odchylka, i = 1, 2, . . . , k. Nulová hypotéza je, ˇze skuteˇcn´ y parametr je nulov´ y, a pˇri jej´ı platnosti má statistika ti t-rozdˇelen´ı o n−k stupn´ıch volnosti. Vybrán je vˇzdy ten nejv´ıce nesignifikantn´ı parametr, fixován na 0 a znovu je proveden odhad ostatn´ıch parametr˚ u. Postupuje se stejnˇe, aˇz nen´ı ˇzádn´ y koeficient nesignifikantn´ı. Zvolená hladina v´ yznamnosti je α = 10%.

3.3

Verifikace modelu

V této kapitole je ovˇeˇrována oprávnˇenost pouˇzit´ı naˇseho modelu a testována kvalita modelu. Pouˇzita je ˇrada statistik a test˚ u, bl´ıˇze jsou popsány jen ty, které nejsou standardnˇe známé. Hladinu v´ yznamnosti test˚ u vol´ıme α = 5%. Pˇredpokladem ARIMA modelu je, ˇze {εt } je posloupnost nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin s normáln´ım rozdˇelen´ım. Pro testován´ı normality rezidu´ı jsme pouˇzili Kolmogorov-Smirnov˚ uv test, Lilliefors˚ uv test a Jarque-Ber˚ uv test. Souˇcasnˇe jsme analyzovali histogram (Obr. 3.7) a QQ-graf (Obr. 3.6). Pro pˇripomenut´ı uved’me, ˇze prvn´ı dva testy jsou zaloˇzeny na empirické distribuˇcn´ı

17

funkci, tˇret´ı na ˇsikmosti a ˇspiˇcatosti normáln´ıho rozloˇzen´ı. Vˇsechny testy hypotézu normáln´ıho rozloˇzen´ı zam´ıtaj´ı, rezidua maj´ı pˇredevˇs´ım v´ yraznou ˇspiˇcatost. Konstantnost rozptylu jsme ovˇeˇrovali Goldfeld-Quandtovým testem, kter´ y nulovou hypotézu konstantnosti rozptylu zam´ıtá. Za nejd˚ uleˇzitˇejˇs´ı ukazatel povaˇzujeme vzájemnou nezávislost rezidu´ı. Anal´ yza korelogram˚ u (Obr. 3.4 a 3.5) dává uspokojivé v´ ysledky, drobnou korelovanost lze pozorovat pouze pro roˇcn´ı frekvence, coˇz potvrzuje i Ljung-Box˚ uv test. Pro u ´plnost uvád´ıme i v´ ysledky Ljung-Boxova testu pro kvadráty rezidu´ı a ovˇeˇren´ı stacionarity a invertibility procesu. Vˇsechny hodnoty statistik jsou uvedeny v tabulce 4.1. Principem Ljung-Boxova testu je, ˇze pˇri platnosti hypotézy nekorelovanosti jsou hodnoty autokorelaˇcn´ı funkce normálnˇe rozloˇzené. Statistika Ljung-Boxova testu má pak tvar L X Q = n(n + 2) ρˆ2 (j)/(n − j), j=1

kde n je délka ˇcasové ˇrady, L poˇcet frekvenc´ı autokorelaˇcn´ı funkce zahrnut´ ych do statistiky a ρˆ2 (j) v´ ybˇerová autokorelaˇcn´ı funkce na frekvenci j. Pˇri platnosti H0 má Q statistika ch´ı-kvadrát rozloˇzen´ı o L stupn´ıch volnosti. Pro srovnáván´ı kvality alternativn´ıch model˚ u jsme pouˇzili vˇerohodnost modelu, AIC a BIC kritérium. Zároveˇ n jsme provádˇeli anal´ yzu nenormovan´ ych rezidu´ı (Obr. 3.9), zvláˇstˇe odchylek pˇresahuj´ıc´ıch 1 mld. Kˇc. D˚ uleˇzitá vlastnost je pˇredpov´ıdac´ı schopnost, mˇeˇrená stˇredn´ı kvadratickou chybou pˇredpovˇedi. Pˇri jednokrokové pˇredpovˇedi za rok 2003 (Obr. 3.8) jsme zaznamenali pˇet odchylek v´ yraznˇe vˇetˇs´ıch neˇz 1 mld. Kˇc. Jedná se o 14.3.2003, 15.5.2003, 5.6.2003, 18.11.2003 a 12.12.2003. Nejzaj´ımavˇejˇs´ı z nich je 14.3.2003, protoˇze velká chyba predikce v tento den nastává i u ostatn´ıch model˚ u. AIC (Akaike’s information criterion) a BIC (Bayes information criterion) jsou kritéria urˇcená pro správnou identifikaci modelu, hlavnˇe pro urˇcen´ı poˇctu parametr˚ u modelu, napˇr´ıklad ˇrád˚ u p a q u ARMA modelu. Vycházej´ı z vˇerohodnosti modelu, ale penalizuj´ı vzr˚ ustaj´ıc´ı poˇcet parametr˚ u. BIC je historicky pozdˇejˇs´ı u ´pravou kritéria AIC, u kterého docház´ı k tendenci pˇreurˇcen´ı modelu oproti skuteˇcnosti. Spoˇc´ıtaj´ı se jako AIC = −2 log L + 2 ppar a BIC = −2 log L + log(n) ppar, kde log L je logaritmická vˇerohodnost a ppar poˇcet parametr˚ u. C´ılem je minimáln´ı hodnota tˇechto kritéri´ı.

18

0.8

0.6

ACF

0.4

0.2

0

−0.2

0

50

100

150

200

250

300

l

Obrázek 3.4: ACF rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı

0.8

0.6

PACF

0.4

0.2

0

−0.2

0

50

100

150 l

200

250

300

Obrázek 3.5: PACF rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı 0.999 0.997 0.99 0.98 0.95 0.90 0.75 p

0.50 0.25 0.10 0.05 0.02 0.01 0.003 0.001 −5

−4

−3

−2

−1

0

1

data

Obrázek 3.6: QQ-graf rezidu´ı

19

2

3

4

140

120

100

80

60

40

20

0 −6

−4

−2

0

2

4

6

Obrázek 3.7: Histogram rezidu´ı

3.4

Predikce

Vzhledem k lineárn´ı struktuˇre modelu lze ARIMA modely snadno pˇredpov´ıdat. Mˇejme ARIMA proces {yt }, kter´ y lze diferencován´ım pˇrevést na ARMA proces {ytd } (1 − B)d yt = ytd ,

t = 1, . . . , n.

Chceme-li vyjádˇrit z pˇredchoz´ı rovnice yt v závislosti na ytd a minulosti ys , s < t dostáváme d µ ¶ X d d (−1)j yt−j , t = 1, . . . , n, yt = yt − j j=1 z ˇcehoˇz m˚ uˇzeme napsat nejlepˇs´ı lineárn´ı pˇredpovˇed’ yt+h pˇri znalosti pozorován´ı ˇcasové ˇrady do ˇcasu t. Znaˇc´ıme ji Pt yt+h a máme Pt yt+h =

d Pt yt+h

−

d µ ¶ X d j=1

j

(−1)j Pt yt+h−j ,

t = 1, . . . , n,

h > 0.

V naˇsem pˇr´ıpadˇe d = 1 a d Pt yt+h = Pt yt+h + Pt yt+h−1 .

Následnˇe vzhledem k tomu, ˇze {ytd } je ARMA proces s regresn´ı komponentou plat´ı d Pt yt+h = dt+h +

p X

d φi Pt yt+h−i +

i=1

q X

θj Pt εt+h−j .

j=1

Nyn´ı staˇc´ı rekurzivnˇe pokraˇcovat pro ˇcasové okamˇziky t + h, t + h − 1, . . . , t + 1, n Pt εs = 0, s > t a Pt εs = εˆs , s ≤ t, pˇriˇcemˇz vyuˇzijeme Pt ysd = ysd , s ≤ t a zároveˇ

20

kde εˆs =

ysd

− ds −

p X

d φi ys−i

−

i=1

q X

θj εˆs−j ,

s = 1, . . . , t

j=1

a εˆs = 0, ysd = 0 pro s ≤ 0. ˇ Poptávka po obˇeˇzivu se v souˇcasné dobˇe v CNB pˇredpov´ıdá vˇzdy na následuj´ıc´ı den, tedy o jeden krok dopˇredu. V souvislosti s pˇrechodem na jednotnou evropskou ˇ e republiky do Evropské mˇenové unie, mˇenu euro, a tedy uˇzˇs´ımu zapojen´ı Cesk´ ale bude nutno provádˇet pˇredpovˇed’ vˇzdy o jeden a dva t´ ydny dopˇredu, tedy o ’ 5, respektive 10 dn´ı vpˇred. Proto je provedena pˇredpovˇed s krokem h = 1, 5, 10. ˇ Kvalita pˇredpovˇedi je mˇeˇrena kritériem pouˇz´ıvan´ ym v CNB, tj. odmocninou stˇredn´ı kvadratické chyby pˇredpovˇedi (root mean squared error), za testovac´ı rok 2003, rP εˆ2t RMSE = , m kde εˆt jsou chyby predikce a m = 252 je poˇcet pozorován´ı ˇcasové ˇrady (pracovn´ıch dn´ı) v testovaném roce. Stˇredn´ı kvadratická chyba se s velikost´ı kroku zvyˇsuje, ˇc´ıselné hodnoty jsou uvedeny v tabulce 4.1.

21

2

1.5

1

0.5

m

0

−0.5

−1

−1.5

−2

50

100

150

200

250

d

Obrázek 3.8: Rezidua jednokrokové predikce v roce 2003 (ˇcernˇe), hranice ±1 mld. Kˇc (ˇcervenˇe) • m - mld. Kˇc • d - denn´ı index 2

1.5

1

0.5

0 m −0.5

−1

−1.5

−2

−2.5

100

200

300

400

500

600

700

d

Obrázek 3.9: Rezidua transformované ˇrady v letech 2000-2002 • m - mld. Kˇc • d - denn´ı index

22

Kapitola 4 GARCH model ˇ Casto pouˇz´ıvan´ ymi modely pro modelován´ı finanˇcn´ıch ˇcasov´ ych ˇrad, napˇr´ıklad cen akci´ı, jsou modely typu GARCH. Aˇckoliv ˇcasová ˇrada poptávky po obˇeˇzivu nen´ı klasickou finanˇcn´ı ˇradou, obsahuje nˇekteré jej´ı typické znaky, a proto jsme se rozhodli tento typ modelu pouˇz´ıt. Náˇs model se skládá z ARMA sloˇzky, která je z velké ˇcásti totoˇzná s prvn´ım modelem poptávky po obˇeˇzivu, a GARCH sloˇzky, jej´ıˇz principy jsou vysvˇetleny v prvn´ı ˇca´sti této kapitoly. Na prvn´ı pohled se m˚ uˇze zdát, ˇze se jedná pouze o rozˇs´ıˇren´ı ARIMA modelu. Protoˇze je vˇsak moˇzné GARCH sloˇzku pouˇz´ıt samostatnˇe a ˇcasová ˇrada se i jinak transformuje, jak je ukázáno v druhé ˇcásti této kapitoly, povaˇzujeme tento model za samostatn´ y model. V tˇret´ı ˇcásti kapitoly je model verifikován a ve ˇctvrté ˇca´sti je zm´ınˇen postup pˇri predikci.

4.1

Z´ akladn´ı principy

Jak jiˇz bylo naznaˇceno, finanˇcn´ı ˇcasové ˇrady maj´ı nˇekolik typick´ ych vlastnost´ı. Definujme nejdˇr´ıve výnos (return) od ˇcasu t − 1 do ˇcasu t jako rt =

st − 1, st−1

kde st je hodnota ˇcasové ˇrady v ˇcase t, napˇr´ıklad cena akcie, u které se nevyplác´ı dividenda. Potom typické znaky pro finanˇcn´ı ˇcasové ˇrady jsou: ˇ a ˇrada nen´ı stacionárn´ı, ale vykazuje lokáln´ı trendy. • Casov´

23

• Hladina v´ ynos˚ u je pˇribliˇznˇe stabiln´ı. Extrémn´ı hodnoty se vyskytuj´ı ve shluc´ıch nerovnomˇernˇe rozloˇzen´ ych v ˇcase. • V´ ynosy rt se zdaj´ı b´ yt nekorelované, ale |rt | a rt2 vykazuj´ı v´ yznamnou korelovanost. • Rozdˇelen´ı rt je leptokurtické, tj. má ˇsikmost K > 3. Jestliˇze pˇredpokládáme hypotézu efektivn´ıho trhu (viz [7]), dostaneme, ˇze v´ ynosy nelze pˇredpov´ıdat. Napˇr´ıklad pokud bychom chtˇeli aplikovat AR(p) model na modelován´ı ceny akci´ı dostali bychom vˇsechny koeficienty φi = 0, i = 1, . . . , p. V praxi je vˇsak napˇr´ıklad kromˇe odhadu hodnoty tˇreba modelovat volatilitu ceny akci´ı. K tomu se daj´ı dobˇre vyuˇz´ıt modely typu ARCH a GARCH. Náhodn´ y proces {yt } naz´ yváme ARCH (autoregressive conditional heteroscedastic process) proces ˇrádu p, jestliˇze E(yt | Ft−1 ) = µ, Var(yt | Ft−1 ) =

σt2

= ω+

p X

2 αk yt−k ,

k=1

εt = σytt je proces nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin, ω, α1 , . . . , αp > 0 a Ft−1 reprezentuje minulost, Ft−1 = σ{yk , k ≤ t − 1}. Bez omezen´ı obecnosti m˚ uˇzeme uvaˇzovat µ = 0. Jestliˇze druhou rovnici v pˇredchoz´ı definici nahrad´ıme rovnic´ı Var(yt | Ft−1 ) =

σt2

=ω+

p X

2 αk yt−k

+

q X

2 βk σt−k ,

k=1

k=1

dostaneme GARCH (generalized ARCH) proces ˇrádu p, q. Znaˇc´ıme krátce ARCH(p), respektive GARCH(p,q). V následuj´ıc´ıch dvou tvrzen´ıch uved’me nˇekolik základn´ıch vlastnost´ı tˇechto proces˚ u, které potvrzuj´ı správnost jejich aplikace na finanˇcn´ı ˇcasové ˇrady. Tvrzen´ı 1. Necht’ {yt } je stacionárn´ı GARCH(p,q) proces s σ 2 = Var(yt ) < ∞, potom ω σ2 = p q X X 1− αk − βk k=1

a je-li nav´ıc p = q = 1, L(εt ) = K =3+

N(0, 1), E(yt4 )

k=1

< ∞ a K ˇsikmost {yt }, plat´ı

6α2 > 3, 1 − 3α2 − 2αβ − β 2

24

pro

α > 0.

Tvrzen´ı 2. Necht’ {yt } je stacionárn´ı GARCH(p,q) proces s E(yt4 ) < ∞ a mˇejme d´ ale m = max(p, q), αj = 0, j > p a βj = 0, j > q, potom je wt = yt2 − σt2 proces nezávislých stejnˇe N(0, konst.) rozdˇelených náhodných veliˇcin a yt2

=ω+

m X

2 γk yt−k

−

k=1

q X

βk wt−k + wt

k=1

ARMA(m,q) proces s γk = αk + βk . D˚ ukazy obou tvrzen´ı lze nalézt v [7]. Z prvn´ıho tvrzen´ı je moˇzno odvodit nutnou a zároveˇ n postaˇcuj´ıc´ı podm´ınku pro stacionaritu GARCH procesu, která zn´ı p X k=1

αk +

q X

βk < 1

a

ω, α1 , . . . , αp , β1 , . . . , βq > 0.

k=1

Protoˇze vˇsak v praxi pˇredpoklad nezávislosti hodnoty ˇcasové ˇrady na minulosti ˇcasto neplat´ı a také se objevuj´ı r˚ uzné vnˇejˇs´ı vlivy, je tˇreba model vhodnˇe rozˇs´ıˇrit, napˇr´ıklad na ARMA(r,m)-GARCH(p,q) proces {yt } s nenulovou stˇredn´ı hodnotou µ a vysvˇetluj´ıc´ı exogenn´ı promˇennou {xt }, kter´ y má tvar yt = µ +

l X

γk xt−k +

r X

k=0

k=1

k=1

k=1

φk yt−k +

ηt = σt εt , q p X X 2 2 2 βk σt−k , σt = ω + αk yt−k +

m X

θk ηt−k + ηt ,

k=1

(4.1)

kde {εt } je proces nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin.

4.2

GARCH model popt´ avky po obˇ eˇ zivu

Podobnˇe jako u ARIMA modelu je tˇreba ˇcasovou ˇradu poptávky po obˇeˇzivu {yt } nejdˇr´ıve stabilizovat, poté nalézt vhodnou strukturu modelu a následnˇe odhadnout jeho parametry, ˇcemuˇz obsahovˇe odpov´ıdaj´ı jednotlivé podkapitoly.

4.2.1

Stabilizace ˇ rady

Jak jiˇz bylo naznaˇceno, u model˚ u typu GARCH se pouˇz´ıvá pro stabilizaci nestacionárn´ı ˇrady jej´ı prvn´ı relativn´ı diference, kterou jsme na stranˇe 23 u finanˇcn´ıch

25

ˇcasov´ ych ˇrad oznaˇcili jako v´ ynos. Alternativn´ım a námi pouˇzit´ ym zp˚ usobem je ale transformace yt y˜t = log . yt−1 D˚ uvodem pouˇzit´ı této metody u finanˇcn´ıch ˇcasov´ ych ˇrad oproti klasickému diferencován´ı je sloˇzené u ´roˇcen´ı, a tedy hypotetick´ y souˇcinov´ y tvar ˇcasové ˇrady, kter´ y se pˇri pouˇzit´ı logaritmu pˇrevede do souˇctového tvaru. Podrobnosti lze nalézt v [7]. 0.04

1

0.8

0.03

0.6

0.02

0.4

0.01 0.2

m

PACF

0

0

−0.01

−0.2

−0.02

−0.4

−0.6

−0.8

−0.03

0

50

100

150 l

200

250

300

−0.04

200

400

600

800

1000

1200

1400

d

Obrázek 4.1: PACF ˇcasové ˇrady Obrázek 4.2: Transformovaná ˇrada obˇeˇziva (ˇcervenˇe), hranice v´ yznamnosti • m - mld. Kˇc • d - denn´ı index (modˇre) • l - zpoˇzdˇen´ı

4.2.2

Struktura modelu

Námi zkonstruovan´ y model je model ARMA(11,1)-GARCH(1,1) s nˇekolika vysvˇetluj´ıc´ımi exogenn´ımi promˇenn´ ymi, jehoˇz matematick´ y zápis odpov´ıdá rovnic´ım (4.1). Nav´ıc jsme pˇredpokládali, ˇze veliˇciny {εt } maj´ı normované normáln´ı rozdˇelen´ı. Pˇresné hodnoty vˇsech koeficient˚ u jsou uvedeny v CD pˇr´ıloze, zde pouze ’ uved me, ˇze v naˇsem modelu je µ 6= 0,

ω 6= 0

a AR polynom má tvar φ(B) = 1 − φ1 B − φ2 B 2 − φ4 B 4 − φ7 B 7 − φ10 B 10 − φ11 B 11 . Urˇcen´ı ˇra´d˚ u ARMA sloˇzky jsme provedli na základˇe stejn´ ych kritéri´ı jako u ARIMA ˇ ad GARCH komponenty odpov´ıdá bˇeˇznˇe pouˇz´ıvanému ˇrádu GARCH modelu. R´ model˚ u ˇcasov´ ych ˇrad.

26

Deterministickou ˇcást tvoˇr´ı exogenn´ı vysvˇetluj´ıc´ı promˇenné, které jsme stejnˇe jako u ARIMA modelu vyuˇzili k zachycen´ı r˚ uzn´ ych sezónn´ıch vliv˚ u, a také jejich tvar z˚ ustal v podstatˇe stejn´ y. Jsou to Vliv dn˚ u v t´ ydnu, Vliv poˇrad´ı dne v mˇes´ıci, Vliv státn´ıch svátk˚ u a dn˚ u volna a Vliv mimoˇra´dn´ ych událost´ı. Jedinou zmˇenou je slouˇcen´ı Nového roku, Svátku práce, Dne v´ıtˇezstv´ı nad faˇsismem, Dne slovansk´ ych vˇerozvˇest˚ u Cyrila a Metodˇeje, Dne upálen´ı mistra Jana Husa, Dne ˇceské státnosti, Dne vzniku samostatného ˇceskoslovenského státu a Dne boje za svobodu a demokracii a jejich implementace jako v´ıcenásobné opakován´ı jednoho pevného svátku. D˚ uvodem pro slouˇcen´ı byla velká v´ ypoˇcetn´ı nároˇcnost modelu.

4.2.3


Parametry modelu jsme odhadli metodou maximáln´ı vˇerohodnosti s postupn´ ym vyˇrazován´ım nesignifikantn´ıch parametr˚ u stejn´ ym zp˚ usobem jako u ARIMA mode’ lu. Pro u ´plnost uved me vˇerohodnostn´ı funkci modelu GARCH(1,1), která se pak zkombinuje s vˇerohodnostn´ı funkc´ı ARMA modelu, ˇc´ımˇz dostaneme vˇerohodnostn´ı funkci naˇseho modelu poptávky po obˇeˇzivu. Tuto funkci je tˇreba numericky maximalizovat, my jsme pouˇzili funkci garchfit programu Matlab. Pˇri pˇredpokladu normovaného normáln´ıho rozloˇzen´ı veliˇcin {εt } má podm´ınˇená vˇerohodnostn´ı funkce procesu GARCH(1,1) {yt } tvar µ ¶ n Y 1 yt2 p L(ω, α, β| ε1 ) = exp − 2 , 2σt 2πσt2 t=2 pˇriˇcemˇz podm´ınˇená a nepodm´ınˇená vˇerohodnostn´ı funkce se liˇs´ı pˇri dostateˇcnˇe dlouhé ˇcasové ˇradˇe zanedbatelnˇe. V´ıce o vˇerohodnostn´ıch funkc´ıch ARIMA a GARCH proces˚ u lze nalézt napˇr´ıklad v [7].

4.3

Verifikace modelu

V této kapitole je uvedeno ovˇeˇren´ı platnosti pˇrijat´ ych pˇredpoklad˚ u. Základn´ım pˇredpokladem je, ˇze veliˇciny ηt εt = σt jsou nezávislé stejnˇe rozloˇzené s normovan´ ym normáln´ım rozdˇelen´ım. Nazvˇeme {εt } normovan´ ymi reziduii. Aplikujeme testy na normalitu Kolmogorov-Smirnov˚ uv test, Lillefors˚ uv test a Jarque-Ber˚ uv test a analyzujeme histogram (Obr. 4.5) a QQ - graf. I kdyˇz dva testy hypotézu normáln´ıho rozloˇzen´ı stále zam´ıtaj´ı, je zde zˇrejmá vˇetˇs´ı evidence pro oprávnˇenost pˇredpokladu normality oproti ARIMA modelu. Konstantnost rozptylu Goldfeld-Quandt˚ uv test zam´ıtá, coˇz je zp˚ usobeno

27

0.8

0.6

ACF

0.4

0.2

0

−0.2

0

50

100

150 l

200

250

300

Obrázek 4.3: ACF normovan´ ych rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı

0.8

0.6

PACF 0.4

0.2

0

−0.2

0

50

100

150 l

200

250

300

Obrázek 4.4: PACF normovan´ ych rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı pˇredevˇs´ım vˇetˇs´ım rozptylem samotné transformované ˇrady v roce 1998. Moˇzn´ ym ˇreˇsen´ım tohoto problému by bylo zkrácen´ı ˇcasové ˇrady právˇe o tento rok. Pˇredpoklad nezávislosti je pˇribliˇznˇe splnˇen, coˇz ukazuj´ı v´ ysledky Ljung-Boxova testu a korelogramy (Obr. 4.3 a 4.4). Ljung-Box˚ uv test pro kvadráty rezidu´ı poukazuje na typické vlastnosti finanˇcn´ıch ˇcasov´ ych ˇrad. Samozˇrejmˇe jsme také ovˇeˇrili stacionaˇ ıselné v´ ritu modelu. C´ ysledky jsou uvedeny v tabulce 4.1. Kvalitu model˚ u jsme porovnávali vˇerohodnost´ı modelu, AIC a BIC kritériem. Zároveˇ n jsme analyzovali graf nenormovan´ ych rezidu´ı (Obr. 4.7). D˚ uleˇzitá vlastnost je samozˇrejmˇe pˇredpov´ıdac´ı schopnost, mˇeˇrená stˇredn´ı kvadratickou chybou pˇredpovˇedi. Pˇri jednokrokové pˇredpovˇedi za rok 2003 (Obr. 4.6) jsme zaznamenali ˇctyˇri odchylky v´ yraznˇe vˇetˇs´ı neˇz 1 mld. Kˇc. Jedná se o 14.3.2003, 5.6.2003, 18.11.2003 a 12.12.2003.

28

140

120

100

80

60

40

20

0 −6

−5

−4

−3

−2

−1

0

1

2

3

4

Obrázek 4.5: Histogram rezidu´ı

4.4

Predikce

Pˇredpovˇed’ budouc´ı hodnoty ˇcasové ˇrady pomoc´ı GARCH modelu závis´ı pouze na jeho ARMA sloˇzce, nebot’ podle pˇredpokladu je {εt } proces nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin s normovan´ ym normáln´ım rozloˇzen´ım. Proto oznaˇcujeme ˇcisté GARCH modely za nepˇredpov´ıdatelné, nebot’ u nich plat´ı Pt−1 yt = E(yt |Ft−1 ) = Eyt = µ (= 0). Naopak moˇznost vyuˇzit´ı GARCH komponenty by byla pˇri odhadu volatility. Postup pˇredpov´ıdán´ı ARMA procesu byl jiˇz popsán na stranˇe 20. Opˇet jsme provedli pˇredpovˇed’ s krokem h = 1, 5, 10. Velikost stˇredn´ı kvadratické chyby se s rostouc´ım krokem znaˇcnˇe zvyˇsuje, ˇc´ıselné hodnoty jsou uvedeny ˇ sen´ım by mohlo b´ v tabulce 4.1. Reˇ yt zˇrejmˇe zv´ yˇsen´ı ˇra´du AR komponenty, coˇz plat´ı i u ARIMA modelu, zmˇena by ale zároveˇ n vedla ke zhorˇsen´ı chyby jednokrokové pˇredpovˇedi, kterou jsme zvolili za hlavn´ı kritérium.

29

2

1.5

1

0.5

m

0

−0.5

−1

−1.5

−2

50

100

150

200

250

d

Obrázek 4.6: Rezidua jednokrokové predikce v roce 2003 (ˇcernˇe), hranice ±1 mld. Kˇc (ˇcervenˇe) • m - mld. Kˇc • d - denn´ı index 0.01

0.008

0.006

0.004

0.002 m

0

−0.002

−0.004

−0.006

−0.008

−0.01

100

200

300

400

500

600

700

d

Obrázek 4.7: Rezidua transformované ˇrady v letech 2000-2002 • m - mld. Kˇc • d - denn´ı index

30

model Nenormovaná rezidua Stˇredn´ı hodnota Rozptyl Normovaná rezidua Stˇredn´ı hodnota Rozptyl ˇ Sikmost ˇ catost Spiˇ Kolmogorov-Smirnov test p-hodnota Lilliefors test p-hodnota Jarque-Bera test p-hodnota Goldfeld-Quandt test p13-hodnota p31-hodnota Ljung-Box test Q(5) p-hodnota Q(10) p-hodnota Q(22) p-hodnota Q(261) p-hodnota Ljung-Box test Q(5) p-hodnota Q(10) p-hodnota Q(22) p-hodnota Q(261) p-hodnota Log. vˇerohodnost AIC BIC Chyby pˇredpovˇed´ı RMSE krok 1 ˇ RMSE krok 1 v CNB RMSE krok 5 RMSE krok 10 Stacionarita Invertibilita

ARIMA

GARCH

0.00038733 0.20951

9.6617x10−5 9.2609x10−6

-6.185x10−17 1 -0.14759 5.7854 Zam´ıtá H0 0.014344 Zam´ıtá H0 NaN1 Zam´ıtá H0 0 Zam´ıtá H0 0.025101 0.9749 rezidua 0.72279 0.74724 0.63773 0.27431 kvadrát rezidu´ı 0 0 0 2.7641x1010 75.7097 -1.4194 384.0551 [mld. Kˇc] 0.48066 0.43787 1.2867 1.8952 Ano Ano

0.012291 0.99186 -0.28308 5.2802 Nezam´ıtá H0 0.14907 Zam´ıtá H0 NaN1 Zam´ıtá H0 0 Zam´ıtá H0 3.2898x10−10 1 0.09926 0.13377 0.001668 5.0804x108 0 0 0 0 5645.7958 -11163.5917 -10834.6534

1 Procedura

0.50111 0.43787 1.5057 2.2858 Ano Ano

lillietest v Matlab v nˇekter´ ych pˇr´ıpadech nevrac´ı p-hodnotu, ale pouze udává v´ ysledek testu.

Tabulka 4.1: Tabulka v´ ysledk˚ u model˚ u

31

Kapitola 5 STS model Obsahem této kapitoly je tˇret´ı model pouˇzit´ y pro modelován´ı ˇcasové ˇrady obˇeˇziva. Jedná se o strukturovan´ y model ˇcasové ˇrady (STS - structural time series model), tedy o jej´ı rozklad na trend, sezónn´ı a nepravidelnou sloˇzku. Tento model se vyjádˇr´ı jako stavov´ y model (state space model). Pomoc´ı Kalmanov´ ych rekurz´ı a metody maximáln´ı vˇerohodnosti se odhadnou jeho parametry a provede pˇredpovˇed’. V prvn´ı ˇcásti této kapitoly je popsána obecná teorie stavov´ ych model˚ u ˇcasov´ ych ˇrad, která vycház´ı z [6]. V druhé je prezentován konkrétn´ı pouˇzit´ y model pro modelován´ı ˇcasové ˇrady obˇeˇziva. Tˇret´ı a ˇctvrtá ˇca´st je vˇenována jeho verifikaci a pˇredpovˇedi ˇrady pomoc´ı modelu.

5.1 5.1.1

Teorie stavov´ ych model˚ uˇ casov´ ych ˇ rad Z´ akladn´ı model

Obecn´ y lineárn´ı Gaussovsk´ y stavov´ y model lze zapsat jako yt = Zt αt + εt , εt ∼ N(0, Ht ), αt+1 = Tt αt + Rt ηt , ηt ∼ N(0, Qt ), t = 1, . . . , n,

(5.1)

kde yt je p × 1 vektor pozorován´ı a αt je nepozorovan´ y m × 1 vektor naz´ yvan´ y stavový vektor. Principem modelu je, ˇze v´ yvoj modelu v ˇcase je urˇcen vektorem αt podle druhé rovnice (5.1). Protoˇze ale αt nelze pozorovat pˇr´ımo, mus´ıme zaloˇzit anal´ yzu na pozorovaném yt . Prvn´ı rovnice (5.1) se naz´ yvá rovnice pozorov´ an´ı, druhá se naz´ yvá stavová rovnice. Matice Zt , Tt , Rt , Ht a Qt jsou povaˇzovány za známé a chybové ˇcleny εt a ηt za seriálnˇe nezávislé a nezávislé na sobˇe pro vˇsechna t.

32

Matice Zt a Tt−1 mohou záviset na y1 , . . . , yt−1 . Dále pˇredpokládáme, ˇze poˇca´teˇcn´ı stavov´ y vektor α1 je N(a1 , P1 ) rozdˇelen´ y, nezávisl´ y na ε1 , . . . , εn a η1 , . . . , ηn . Pro zaˇca´tek povaˇzujme a1 a P1 za známé. Jeden z nejjednoduˇsˇs´ıch stavov´ ych model˚ u je model lokáln´ıho lineárn´ıho trendu (local linear trend model) yt = µt + εt , εt ∼ N(0, σε2 ), µt+1 = µt + νt + ξt , ξt ∼ N(0, σξ2 ), νt+1 = νt + ζt , ζt ∼ N(0, σζ2 ),

(5.2)

Kdyˇz σξ2 = σζ2 = 0, potom νt+1 = νt = ν0 a µt+1 = µt + ν0 . Trend je tedy pˇresnˇe lineárn´ı a model se redukuje do modelu deterministického lineárn´ıho trendu s náhodn´ ym ˇsumem. Tvar (5.2) se σξ2 > 0 a σζ2 > 0 oproti tomu umoˇznuje zmˇeny hladiny a trendu v ˇcase. Rovnice (5.2) lze pˇrepsat do tvaru µ ¶ ¡ ¢ µt 1 0 yt = + εt , µ ¶ · ¸ µ νt ¶ µ ¶ µt+1 1 1 µt ξt = + , νt+1 0 1 νt ζt coˇz je speciáln´ı pˇr´ıpad (5.1). V´ yhodou stavov´ ych model˚ u je moˇznost jejich spojován´ı v komplexnˇejˇs´ı modely. Jako stavov´ y model lze kromˇe lokáln´ıho lineárn´ıho trendu totiˇz zapsat napˇr´ıklad mˇes´ıˇcn´ı sezónnost modelovanou pomoc´ı kubick´ ych splin˚ u (kapitola 5.2.2). Tyto jednotlivé stavové modely lze za pˇredpokladu vzájemné nezávislosti α1i , {εt } a {ηti } , i = 1, 2 spojit do jednoho stavového modelu ¡ 1 ¢0 (5.3) αt αt2 αt = a

Zt Tt Rt Qt

¡ 1 ¢ Zt ¡ Zt2 , = = diag ¡ Tt1 , Tt2 = diag ¡ Rt1 , Rt2 = diag Q1t , Q2t

¢ ¢, ¢, ,

(5.4)

kde αti je stavov´ y vektor, {ηti } náhodn´ y ˇclen a Zti , Tti , Rti a Qit jsou matice stavového modelu lokáln´ıho lineárn´ıho trendu (i = 1) a mˇes´ıˇcn´ı sezónnosti (i = 2). Znaˇcen´ı odpov´ıdá definici (5.1). Postup lze analogicky rozˇs´ıˇrit na libovoln´ y poˇcet stavov´ ych model˚ u, coˇz pouˇzijeme pˇri konstrukci STS modelu poptávky po obˇeˇzivu.

5.1.2

Kalmanovy rekurze

Základn´ı tˇri u ´lohy anal´ yzy ˇcasové ˇrady stavov´ ymi modely jsou filtrov´ an´ı, vyhlazov´ an´ı a predikce. Jde o problém nalezen´ı nejlepˇs´ıho (ve smyslu nejmenˇs´ı stˇredn´ı kvadratické chyby) lineárn´ıho odhadu stavového vektoru αt pˇri znalosti pozorován´ı

33

• y1 , . . . , yt−1 • y1 , . . . , y t

pˇri predikci, pˇri filtrován´ı,

• y1 , . . . , y n , n > t

pˇri vyhlazován´ı.

Ve stavov´ ych modelech se tyto u ´lohy ˇreˇs´ı r˚ uzn´ ymi skupinami Kalmanov´ ych rekurz´ı. Ty vyuˇz´ıvaj´ı markovské vlastnosti stavov´ ych model˚ u a jsou v´ ypoˇcetnˇe efektivn´ı. Základem je Kalman˚ uv filtr. Následuje jeho odvozen´ı ([6]), znaˇcen´ı je stejné jako v základn´ım tvaru (5.1). C´ılem bude odvodit podm´ınˇené rozdˇelen´ı αt+1 za podm´ınky Yt pro t = 1, . . . , n, kde Yt = {y1 , . . . , yt }, coˇz je jednokroková predikce. Na závˇer drobn´ ym pˇreformulován´ım dostaneme zároveˇ n i v´ ysledek u ´lohy filtrován´ı. Jelikoˇz rozdˇelen´ı vˇsech náhodn´ ych veliˇcin je normáln´ı, podm´ınˇené rozdˇelen´ı podmnoˇziny veliˇcin podm´ınˇené jinou podmnoˇzinou tˇechto veliˇcin je také normáln´ı. Poˇzadované rozdˇelen´ı je tedy urˇceno znalost´ı at+1 = E(αt+1 | Yt ) a Pt+1 = Var(αt+1 | Yt ). Pˇredpokládejme, ˇze αt za podm´ınky Yt−1 má rozdˇelen´ı N(at , Pt ). Jelikoˇz αt+1 = Tt αt + Rt ηt , dostáváme at+1 = E(Tt αt + Rt ηt | Yt ) = Tt E(αt | Yt ),

(5.5)

Pt+1 = Var(Tt αt + Rt ηt | Yt ) = Tt Var(αt | Yt )Tt0 + Rt Qt Rt0 ,

(5.6)

pro t = 1, . . . , n. Oznaˇcme vt = yt − E(yt | Yt−1 ) = yt − E(Zt αt + εt | Yt−1 ) = yt − Zt at . Potom je vt chyba jednokrokové pˇredpovˇedi yt za podm´ınky Yt−1 . Kdyˇz Yt−1 a vt jsou známé, pak také Yt je známé. Tedy E(αt | Yt ) = E(αt | Yt−1 , vt ). Následnˇe E(vt ) = 0 a Cov(yj , vt ) = E[yj E(vt | Yt−1 )0 ] = 0, pro j = 1, . . . , t−1. Z v´ıcerozmˇerné normáln´ı regrese dostáváme E(αt |Yt ) = E(αt | Yt−1 , vt ) = E(αt | Yt−1 ) + Cov(αt , vt )[Var(vt )]−1 vt = at + Mt Ft−1 vt , kde Mt = Cov(αt , vt ) a Ft = Var(vt ). Dále Mt = Cov(αt , vt ) = E[E{αt (Zt αt + εt − Zt at )0 | Yt−1 }] = E[E{αt (αt − at )0 Zt0 | Yt−1 }] = Pt Zt0 a Ft = Var(Zt αt + εt − Zt at ) = Zt Pt Zt0 + Ht .

34

(5.7)

Pˇredpokládejme, ˇze Ft je nesingulárn´ı. Tato podm´ınka je bˇeˇznˇe splnˇena pro vˇsechny dobˇre formulované modely, v´ıce v [6]. Substituc´ı v (5.5) a (5.7) vyjde at+1 = Tt at + Tt Mt Ft−1 vt = Tt at + Kt vt , t = 1, . . . , n,

(5.8)

kde Kt = Tt Mt Ft−1 = Tt Pt Zt0 Ft−1 . Tedy at+1 je lineárn´ı funkce at a vt . Podobnˇe máme Var(αt | Yt ) = = = =

Var(αt | Yt−1 , vt ) Var(αt | Yt−1 ) − Cov(αt , vt )[Var(vt )]−1 Cov(αt , vt )0 Pt − Mt Ft−1 Mt0 Pt − Pt Zt0 Ft−1 Zt Pt .

Substituc´ı do (5.6) dostaneme Pt+1 = Tt Pt L0t + Rt Qt Rt0 ,

t = 1, . . . , n,

(5.9)

kde Lt = Tt − Kt Zt . Rekurze (5.8) a (5.9) jsou základem Kalmanova filtru pro model (5.1). Pro pˇrehlednost vˇsechny rovnice jsou vt = yt − Zt at , Ft = Zt Pt Zt0 + Ht , 0 −1 Kt = Tt Pt Zt Ft , L t = T t − K t Zt , t = 1, . . . , n, 0 0 at+1 = Tt at + Kt vt , Pt+1 = Tt Pt Lt + Rt Qt Rt ,

(5.10)

kde a1 je vektor stˇredn´ıch hodnot a P1 je varianˇcn´ı matice poˇcáteˇcn´ıho stavového vektoru α1 . Rekurze (5.10) naz´ yváme (standardn´ı) Kalman˚ uv filtr. Pˇreformulován´ım Kalmanova filtru lze nyn´ı lehce odvodit rekurentn´ı rovnice pro samotné filtrován´ı. Oznaˇcme odhad stavového vektoru at|t = E(αt | Yt ) a jemu pˇr´ısluˇsné varianˇcn´ı matice Pt|t = Var(αt | Yt ). Dostaneme soubˇeˇzný Kalman˚ uv filtr vt = yt − Zt at , at|t at+1

Ft Mt −1 = at + Mt Ft vt , Pt|t = Tt at|t , Pt+1

= = = =

Zt Pt Zt0 + Ht , Pt Zt0 Pt − Mt Ft−1 Mt0 , t = 1, . . . , n, Tt Pt|t Tt0 + Rt Qt Rt0 .

Dalˇs´ı skupinou rekurentn´ıch rovnic je Kalman˚ uv filtr pro vyhlazován´ı. Jejich v´ ysledkem je odhad αt pˇri znalosti celé ˇcasové ˇrady. Oznaˇcme y = (y10 , . . . , yn0 )0 . Necht’ α bt = E(αt | y) a Vt = Var(αt − α bt ) = Var(αt | y) pro t = 1, . . . , n. Za podm´ınky známého rozdˇelen´ı poˇca´teˇcn´ıho stavového vektoru α1 ∼ N(a1 , P1 ) má Kalman˚ uv filtr pro vyhlazován´ı tvar rt−1 = Zt0 Ft−1 vt + L0t rt , Nt−1 = Zt0 Ft−1 Zt + L0t Nt Lt , α bt = at + Pt rt−1 , Vt = Pt − Pt Nt−1 Pt , 35

(5.11)

pro t = n, . . . , 1, s inicializac´ı rn = 0 a Nn = 0. Pˇri samotném v´ ypoˇctu postupujeme následovnˇe. Nejdˇr´ıve projdeme ˇradu dopˇrednˇe standardn´ım Kalmanov´ ym filtrem, následnˇe zpˇetnˇe Kalmanov´ ym filtrem pro vyhlazován´ı a t´ım z´ıskáme odhady α bt a Vt . Pˇri dopˇredném pr˚ uchodu mus´ıme uloˇzit veliˇciny vt , Ft , Kt , at a Pt pro t = 1, . . . , n. Alternativnˇe lze uloˇzit pouze at a Pt a ostatn´ı veliˇciny dopoˇc´ıtat pˇri zpˇetném pr˚ uchodu, coˇz je pouˇzit´ y zp˚ usob v naˇsem modelu ˇcasové ˇrady obˇeˇziva. Tˇret´ı u ´loha, a to predikce pomoc´ı Kalmanova filtru je popsána v kapitole 5.4. Znaˇcnou v´ yhodou tohoto pˇr´ıstupu je moˇznost anal´ yzy ˇcasové ˇrady v pˇr´ıpadˇe, ˇze nen´ı kompletn´ı a jej´ı pozorován´ı v nˇekter´ ych ˇcasov´ ych okamˇzic´ıch chyb´ı. Pˇredpokládejme, ˇze máme obecn´ y stavov´ y model, kde pozorován´ı yj , j = τ, . . . , τ ∗ − 1, ∗ chyb´ı, pro 1 < τ < τ ≤ n. Pro filtrován´ı v ˇcasech t = τ, . . . , τ ∗ − 1 staˇc´ı pouˇz´ıt v Kalmanovˇe filtru vektor vt = 0 a matici Kt = 0 pro tyto ˇcasové okamˇziky. Dostaneme at+1 = Tt at , Pt+1 = Tt Pt Tt0 + Rt Qt Rt0 , t = τ, . . . , τ ∗ − 1. Podobnˇe pro zpˇetné vyhlazovac´ı rekurze máme rt−1 = Tt0 rt , Nt−1 = Tt0 Nt Tt , t = τ ∗ − 1, . . . , τ. Ostatn´ı rekurzivn´ı rovnice z˚ ustanou nezmˇenˇeny. Ovˇeˇren´ı tˇechto v´ ysledk˚ u lze nalézt v [6]. V naˇsem modelu byl vzhledem k jiˇz hotové implementaci Kalmanova filtru pouˇzit alternativn´ı pˇr´ıstup, kter´ y je s popsan´ ym ekvivalentn´ı. Byl pouˇzit pˇreformulovan´ y model ½ 0 t ∈ {τ, . . . , τ ∗ − 1}, ∗ yt = yt jinak, ½ ½ ∗ 0 t ∈ {τ, . . . , τ − 1}, δt t ∈ {τ, . . . , τ ∗ − 1}, Zt∗ = ε∗t = Zt jinak, εt jinak, kde {δt } je proces nezávisl´ ych stejnˇe rozdˇelen´ ych náhodn´ ych veliˇcin s δt ∼ N(0, I), nezávisl´ ych na α1 , εs , ηs pro t, s = 1, . . . , n. Ostatn´ı ˇcásti modelu z˚ ustávaj´ı nezmˇenˇeny. V pˇr´ıpadˇe v´ıce chybˇej´ıc´ıch ˇca´st´ı pozorován´ı ˇcasové ˇrady se tyto pˇr´ıstupy analogicky rozˇs´ıˇr´ı.

5.1.3

Inicializace Kalmanova filtru

V pˇredcházej´ıc´ıch dvou kapitolách jsme povaˇzovali rozdˇelen´ı poˇca´teˇcn´ıho stavu za známé, α1 ∼ N(a1 , P1 ), a1 a P1 známé. To vˇsak ve vˇetˇsinˇe praktick´ ych aplikac´ı, vˇcetnˇe modelu poptávky po obˇeˇzivu, neplat´ı. Proces ˇreˇs´ıc´ı tento problém se naz´ yvá inicializace. V obecném pˇr´ıpadˇe má poˇcáteˇcn´ı stav modelu tvar α1 = a + A δ + R0 η0 ,

36

η0 ∼ N(0, Q0 ),

(5.12)

kde m × 1 vektor a je znám´ y, δ je q × 1 vektor neznám´ ych veliˇcin, m × q matice A a m × (m − q) matice R0 jsou selekˇcn´ı matice. Selekˇcn´ı matice jsou matice, jejichˇz sloupce jsou sloupce jednotkové matice Im a dané dohromady daj´ı mnoˇzinu g sloupc˚ u matice Im tak, ˇze g ≤ m a A0 R0 = 0. Matice Q0 je pozitivnˇe semidefinitn´ı a známá. Reprezentace konstantn´ı ˇcásti je a, A δ pˇredstavuje nestacionárn´ı ˇca´st a R0 η0 stacionárn´ı ˇca´st. V modelu poptávky po obˇeˇzivu je cel´ y poˇcáteˇcn´ı stav neznám´ y a veliˇciny stavového vektoru nestacionárn´ı, tedy A = Im a α1 = δ. Vektor δ m˚ uˇze b´ yt povaˇzován za pevn´ y vektor neznám´ ych parametr˚ u, tedy δ ∼ N(konst., 0),

(5.13)

nebo jako vektor náhodn´ ych veliˇcin s normáln´ım rozdˇelen´ım s nulov´ ymi stˇredn´ımi hodnotami a nekoneˇcn´ ymi rozptyly. Prvn´ı zp˚ usob odvodil Rosenberg a je v´ıce diskutován v [6]. Druh´ ym zp˚ usobem se mimo jiné zab´ yval De Jong v [4] a naz´ yvá se difusn´ı inicializace. Vycház´ı z pˇredpokladu δ ∼ N(0, κ Iq ),

(5.14)

kde κ → ∞. Pˇri modelován´ı poptávky po obˇeˇzivu byly pouˇzity obˇe varianty. Ve fináln´ım modelu byla zvolena difusn´ı inicializace, aˇckoliv v´ ysledky obou variant byly v podstatˇe rovnocenné. Jednou moˇznost´ı ˇreˇsen´ı je rozˇs´ıˇrený Kalman˚ uv filtr (augmented Kalman filter), jehoˇz základem je rozˇs´ıˇren´ı pozorovaného vektoru ˇcasové ˇrady následuj´ıc´ım zp˚ usobem. Pro dané δ aplikujme Kalman˚ uv filtr s a1 = E(α1 ) = a + A δ, 0 P1 = Var(α1 ) = P∗ = R0 Q0 R0 a oznaˇcme v´ yslednou hodnotu at jako aδ,t . Jelikoˇz aδ,t je lineárn´ı funkce pozorován´ı a a1 = a + A δ m˚ uˇzeme psát aδ,t = aa,t + AA,t δ, kde aa,t je hodnota at v´ ystupu Kalmanova filtru, u kterého je a1 = a, P1 = P∗ a kde j-t´ y sloupec AA,t je hodnota at v´ ystupu fitru, u kterého je vektor pozorován´ı nulov´ y a a1 = Aj , P1 = P∗ , kde Aj je j-t´ y sloupec matice A. Oznaˇcme hodnotu vt v´ ystupu Kalmanova filtru, u kterého je a1 = a + A δ, P1 = P∗ jako vδ,t . Podobnˇe jako pro aδ,t m˚ uˇzeme psát vδ,t = va,t + VA,t δ, kde va,t a VA,t jsou dané v´ ysledkem stejn´ ych filtr˚ u jako aa,t a AA,t . Matice (aa,t , AA,t ) a (va,t , VA,t ) mohou b´ yt spoˇcteny jedin´ ym pr˚ uchodem Kalmanova filtru, do nˇehoˇz vstupuje vektor pozorován´ı yt rozˇs´ıˇren´ y nulov´ ymi hodnotami. Toto je moˇzné, protoˇze v kaˇzdém filtru produkuj´ıc´ım jednotliv´ y sloupec (aa,t , AA,t ) je pouˇzita stejná inicializaˇcn´ı varianˇcn´ı matice P1 = P∗ , a tedy v´ ystup t´ ykaj´ıc´ı se rozptylu, kter´ y oznaˇc´ıme jako Fδ,t , Kδ,t a Pδ,t+1 je stejn´ y. Nahrazen´ım rovnic Kalmanov´ ych filtr˚ u pro vt a at odpov´ıdaj´ıc´ımi rovnicemi pro matice (aa,t , AA,t ) a (va,t , VA,t ) dostaneme (va,t , VA,t ) = (yt , 0) − Zt (aa,t , AA,t ), (aa,t+1 , AA,t+1 ) = Tt (aa,t , AA,t ) + Kδ,t (va,t , VA,t ), 37

(5.15)

kde (aa,1 , AA,1 ) = (a, A). Rekurze pro Ft , Kt a Pt+1 z˚ ustanou stejné jako v klasickém Kalmanovˇe filtru, tedy Fδ,t = Zt Pδ,t Zt0 + Ht , −1 Kδ,t = Tt Pδ,t Zt0 Fδ,t , Lδ,t = Tt − Kδ,t Zt , 0 0 Pδ,t+1 = Tt Pδ,t Lδ,t + Rt Qt Rt ,

(5.16)

pro t = 1, . . . , n a Pδ,1 = P∗ . Index δ je pˇridán proto, ˇze daná veliˇcina je vypoˇc´ıtána za pˇredpokladu daného δ a nikoliv, ˇze na δ pˇr´ımo matematicky závis´ı. Jako rozˇs´ıˇren´ y Kalman˚ uv filtr budeme oznaˇcovat rovnice (5.15) a (5.16). Necht’ nyn´ı plat´ı (5.14). Pro dané κ máme at+1 = E(αt+1 | Yt ) = aa,t+1 + AA,t+1 δ¯t , kde δ¯t = E(δ| Yt ). Dále log p(δ| Yt ) = log p(δ) + log p(Yt | δ) − log p(Yt ) t X = log p(δ) + log p(vδ,j ) − log p(Yt ) j=1

1 1 = − δ 0 δ − b0t δ − δ 0 SA,t δ + ˇcleny nezávislé na δ, 2κ 2 kde bt =

t X

0 VA,j

−1 Fδ,j

va,j ,

SA,t =

j=1

t X

−1 0 VA,j Fδ,j VA,j .

(5.17)

j=1

Diferencován´ım vzhledem k δ a poloˇzen´ım rovno 0 dostaneme maximálnˇe vˇerohodn´ y odhad a pˇri normáln´ımu rozloˇzen´ı veliˇcin máme µ ¶−1 1 ¯ δt = − SA,t + Iq bt . (5.18) κ Tedy Pt+1 = E[(at+1 − αt+1 )(at+1 − αt+1 )0 ] = E[{aa,t+1 − αt+1 + AA,t+1 (δ − δ¯t )}{aa,t+1 − αt+1 + AA,t+1 (δ − δ¯t )}0 ] = Pδ,t+1 + AA,t+1 Var(δ| Yt )A0A,t+1 µ ¶−1 1 = Pδ,t+1 + AA,t+1 SA,t + Iq A0A,t+1 , κ jelikoˇz Var(δ| Yt ) = (SA,t + κ1 Iq )−1 . Pro κ → ∞ dostáváme −1 δ¯t = −SA,t bt , −1 Var(δ| Yt ) = SA,t ,

v pˇr´ıpadˇe, ˇze SA,t je nesingulárn´ı. V´ ypoˇcet bt a SA,t lze snadno pˇridat k rozˇs´ıˇrenému Kalmanovu filtru. Potom −1 bt , at+1 = aa,t+1 − AA,t+1 SA,t −1 0 AA,t+1 , Pt+1 = Pδ,t+1 + AA,t+1 SA,t

38

pˇri κ → ∞. Pro nedegenerativn´ı modely existuje d takové, ˇze pro t < d je SA,t singulárn´ı, takˇze hodnoty at+1 a Pt+1 dané pˇredchoz´ımi rovnicemi neexistuj´ı. Naopak pro t ≥ d je SA,t nesingulárn´ı, hodnoty ad+1 a Pd+1 existuj´ı a následnˇe hodnoty at+1 a Pt+1 pro t > d mohou b´ yt vypoˇc´ıtány standardn´ım Kalmanov´ ym filtrem. Obvykle je d menˇs´ı neˇz délka sezóny modelované ˇcasovou ˇradou, a tedy mnohem menˇs´ı neˇz délka pozorované ˇcasové ˇrady. V modelu poptávky po obˇeˇzivu je vˇsak potˇreba d = 828, coˇz jsou tˇri sezóny (roky) a polovina pozorované ˇcasové ˇrady. D˚ uvodem je zˇrejmˇe soubˇeh tˇr´ı modelovan´ ych vliv˚ u ve stejném obdob´ı, a to na konci kalendáˇrn´ıho roku. Jedná se o roˇcn´ı sezónnost, Nov´ y rok a problém Y2K. Takto dlouhá inicializace znemoˇzn ˇuje pouˇzit´ı jinak v´ ypoˇcetnˇe efektivnˇejˇs´ıho pˇresného inicializaˇcn´ıho Kalmanova filtru (exact initial Kalman filter), protoˇze u nˇej docház´ı k velk´ ym numerick´ ym chybám. Eventuáln´ı moˇznost postupu by byla pouˇzit´ı jeho 0 verze s rozkladem varianˇcn´ıch matic na mocninn´ y souˇcin Pt = P˜t P˜t (square root filter), a zaruˇcit tak jejich pozitivn´ı semidefinitnost, podrobnosti lze nalézt v [6].

5.1.4


Model poptávky po obˇeˇzivu, stejnˇe jako vˇsechny praktické modely, závis´ı na neznám´ ych parametrech. V této kapitole ukáˇzeme, jak lze tyto parametry odhadnout metodou maximáln´ı vˇerohodnosti. Pro zˇrejmost uved’me, ˇze napˇr´ıklad u lokáln´ıho lineárn´ıho modelu (5.2) jsou neznámé parametry velikosti rozptyl˚ u σε2 , σξ2 a σζ2 . Nejdˇr´ıve odvod’me vˇerohodnostn´ı funkci v pˇr´ıpadˇe, ˇze rozdˇelen´ı poˇcáteˇcn´ıho stavového vektoru je známé, α1 ∼ E(a1 , P1 ), a1 a P1 jsou známé. Pak vˇerohodnost modelu je n Y L(y) = p(y1 , . . . , yn ) = p(y1 ) p(yt | Yt−1 ), t=2

kde Yt = {y1 , . . . , yt }. Pro praktick´ y v´ ypoˇcet pak pouˇzijeme logaritmickou vˇerohodnostn´ı funkci n X log L(y) = log p(yt | Yt−1 ), t=1

kde p(y1 | Y0 ) = p(y1 ). Pro model (5.1) máme E(yt | Yt−1 ) = Zt at , vt = yt − Zt at a Ft = Var(yt | Yt−1 ). Substituován´ım N(Zt at , Ft ) za p(yt | Yt−1 ) dostaneme n

np 1X log L(y) = − log 2π − (log |Ft | + vt0 Ft−1 vt ). 2 2 t=1

(5.19)

Hodnoty vt , Ft , a tedy celá log L(y), se vypoˇc´ıtaj´ı pomoc´ı Kalmanova filtru. Pˇredpokládáme jen, ˇze Ft je nesingulárn´ı pro t = 1, . . . , n. Jestliˇze tato podm´ınka

39

nen´ı splnˇena automaticky, je obvykle moˇzné pˇredefinovat model tak, aby splnˇena byla. V naˇsem modelu poptávky po obˇeˇzivu se ale tento problém nevyskytl. V difusn´ım pˇr´ıpadˇe ˇreˇseném rozˇs´ıˇren´ ym Kalmanov´ ym filtrem je sdruˇzená hustota δ a y pro dané κ p(δ, y) = p(δ)p(y| δ) n X = p(δ) p(vδ,t ) t=1

n Y = (2π)(np+q)/2 κ−q/2 | Fδ,t |−1/2 × t=1 µ ¶¸ · 1 δ0δ 0 0 + Sa,n + 2bn δ + δ SA,n δ , × exp − 2 κ

kde vδ,t , bn , SA,n jsou definovány v pˇredcházej´ıc´ı kapitole a Sa,n =

n X

−1 0 va,t Fδ,t va,t .

t=1 −1 Z (5.18) známe δ¯n = − (SA,n + κ−1 Iq ) bn a exponent pˇredcházej´ıc´ıho v´ yrazu m˚ uˇze b´ yt pˇrepsán jako

1 − [Sa,n + (δ + δ¯n )0 (SA,n + κ−1 Iq )(δ − δ¯n ) − δ¯n0 (SA,n + κ−1 Iq )δ¯n ]. 2 Odintegrován´ım δ dostaneme margináln´ı hustotu y. Logaritmická vˇerohodnostn´ı funkce má tvar np q 1 log L(y) = − log 2π − log κ − log | SA,n + κ−1 Iq | − 2 n 2 2 1X 1 − log |Fδ,t | − [Sa,n − δ¯n0 (SA,n + κ−1 Iq )δ¯n ]. 2 t=1 2 Pˇriˇcten´ım funkci

q 2

log κ a pro κ → ∞ dostaneme difusn´ı logaritmickou vˇerohodnostn´ı n

np 1 1X log Ld (y) = − log 2π − log | SA,n | − log | Fδ,t | 2 2 2 t=1 (5.20) 1 0 −1 = − (Sa,n − bn SA,n bn ). 2 V pˇredchoz´ı kapitole jsme ukázali, ˇze rozˇs´ıˇren´ y Kalman˚ uv filtr m˚ uˇze b´ yt pˇreruˇsen v ˇcasovém okamˇziku t = d. M˚ uˇzeme tedy pouˇz´ıt parciáln´ı vˇerohodnost zaloˇzenou na Yd a následnˇe pˇridat pˇr´ıspˇevky inovac´ı vd+1 , . . . , vn ze standardn´ıho Kalmanova filtru. Odhad parametr˚ u pak z´ıskáme numerickou maximalizac´ı této vˇerohodnostn´ı funkce. V naˇsem modelu poptávky po obˇeˇzivu jsme pouˇzili funkci NMinimize programu Mathematica. Numerická maximalizace je znaˇcnˇe v´ ypoˇcetnˇe a ˇcasovˇe nároˇcná, moˇzn´ ym zrychlen´ım by bylo napˇr´ıklad pouˇzit´ı EM algoritmu, nebo pˇrevod problému do spektr´ aln´ı domény.

40

5.2

STS model popt´ avky po obˇ eˇ zivu

Pro modelován´ı poptávky po obˇeˇzivu jsme pouˇzili STS model yt = µt + γt + εt , kde µt je model lokáln´ıho lineárn´ıho trendu (5.2), γt je stochastická sezónn´ı komponenta a nepravidelná sloˇzka εt je proces nezávisl´ ych náhodn´ ych veliˇcin s normáln´ım 2 rozdˇelen´ım s nulovou stˇredn´ı hodnotou a rozptylem σε . Sezónn´ı komponenta je P sumou pˇeti subkomponent γt = 5i=1 γti , které modeluj´ı rozd´ılné sezónn´ı vlivy. Jedná se o sezónnost dn˚ u v týdnu, mˇes´ıˇcn´ı sezónnost, roˇcn´ı sezónnost, státn´ı svátky a dny volna a mimoˇr´ adné události. Kaˇzdá subkomponenta má tvar γti = zti δti , kde zti je pevn´ y 1 × gi vektor nezávisl´ ych ˇcasovˇe závisl´ ych“ promˇenn´ ych. Hodnoty ” u pro zti jsou definovány vˇzdy pˇres celé obdob´ı které modeluj´ı, napˇr´ıklad 5 dn˚ sezónnost dne v t´ ydnu, a jejich hodnoty se periodicky opakuj´ı pro celou ˇcasovou ˇradu. Aby se zajistilo, ˇze γti je ˇcistˇe sezónn´ı komponenta, mus´ı b´ yt souˇcet zti pˇres P i celé sezónn´ı obdob´ı nulov´ y vektor, tedy u sezónnosti dne v t´ ydnu t+4 k=t zk = 0. δti je ˇcasovˇe promˇenliv´ y parametr, kter´ y sleduje náhodnou procházku i δti = δt−1 + χit ,

kde χit je gi × 1 vektor sériálnˇe nezávisl´ ych normálnˇe rozdˇelen´ ych náhodn´ ych proces˚ u s nulovou stˇredn´ı hodnotou a kovarianˇcn´ı matic´ı E(χit χit ‘) = σi2 I. Ne vˇsechny sezónn´ı komponenty v modelu jsou stochastické, nˇekteré jsou uvaˇzovány jako deterministické, a potom je σi2 = 0. Jednotlivé komponenty se spoj´ı podle (5.3) a (5.4) a matice stavového zápisu modelu pak jsou Zt = (1 0 zt1  µ Tt = 

zt2 1 1 0 1 0

zt3

zt4

¶

zt5 )

Ht = σε2

 0  I

Rt = I

Qt = diag(σξ2 , σζ2 , σ12 , . . . , σ12 , σ22 , . . . , σ22 , σ32 , . . . , σ32 , 0, . . . , 0). Aby bylo moˇzné modelovat vˇsechny mˇes´ıce a roky, je nutné uváˇzit mˇes´ıc s 23 (pracovn´ımi) dny a tedy rok s 276 dny. Jestliˇze má nˇekter´ y mˇes´ıc ménˇe neˇz 23 dn˚ u, povaˇzujeme tyto dny za chybˇej´ıc´ı pozorován´ı. Pro ilustraci uved’me 31. duben, nebo 29. u ńor v nepˇrestupn´ ych letech, pokud se samozˇrejmˇe jedná o pracovn´ı dny.

41

5.2.1

Modelovan´ e sez´ onnosti

Sez´ onnost dn˚ u v t´ ydnu K zachycen´ı této závislosti jsme pouˇzili model 1 1 1 γt = zt δt , kde zt1 je vektor 1 × 4 a jeho prvn´ı prvek nab´ yvá 1, je-li den t pondˇel´ı, −1 pátek a 0 pro ostatn´ı dny. Prvky 2, 3 a 4 nab´ yvaj´ı stejn´ ych hodPt+4 1 not pro u ´ter´ y, stˇredu a ˇctvrtek. Z definice je zˇrejmé, ˇze k=t zk je nulov´ y 1 vektor. Sezónnost dn˚ u v t´ ydnu je uvaˇzována jako stochastická, a tedy δt je vektor sleduj´ıc´ı náhodnou procházku s rozptylem σ12 . Mˇ es´ıˇ cn´ı sez´ onnost Pro modelován´ı jsme zvolili kubické spliny promˇenlivé v ˇcase, jejichˇz pˇresná konstrukce je popsána v dalˇs´ı kapitole. Nejd˚ uleˇzitˇejˇs´ım faktorem je zde urˇcen´ı poˇctu uzl˚ u a jejich správné um´ıstˇen´ı v rámci modelovaného obdob´ı. Na základˇe korelogram˚ u, kvality proloˇzen´ı a pˇredpov´ıdac´ı schopnosti modelu jsme vybrali poˇcet ˇsesti uzl˚ u a rozm´ıstˇen´ı (0), 3, 9, 11, 15, 19, 23. Nula je zde uvedena pouze z d˚ uvodu objasnˇen´ı navázán´ı dvou obdob´ı, jedná se vlastnˇe o uzel na pozici 23 pˇredcházej´ıc´ı sezóny. Vzhledem k periodicitˇe kubick´ ych splin˚ u mus´ı m´ıt spliny v deterministickém pˇr´ıpadˇe v prvn´ım (0) a posledn´ım (23) uzlu stejnou hodnotu. Z tohoto d˚ uvodu je d˚ uleˇzité, aby ˇcasová ˇrada na zaˇca´tku modelovaného obdob´ı mˇela pˇribliˇznˇe stejné hodnoty a tedy ˇza´dn´ y prudk´ y trend. Proto prvn´ı den v modelovaném mˇes´ıci neodpov´ıdá prvn´ımu dni v kalendáˇrn´ım mˇes´ıci, ale je posunut aˇz na 5. den kalendáˇrn´ıho mˇes´ıce. Um´ıstˇen´ı chybˇej´ıc´ıch dn˚ u do celkového poˇctu 23 dn˚ u za mˇes´ıc je zaloˇzeno na stejném principu. V modelu poptávky po obˇeˇzivu jsme um´ıstili tyto chybˇej´ıc´ı pozorován´ı mezi 4. a 5. den v mˇes´ıci. Mˇes´ıˇcn´ı komponenta je rovnˇeˇz stochastická s parametrem σ22 . Roˇ cn´ı sez´ onnost Jako pˇr´ıklad závislosti poptávky po obˇeˇzivu na roˇcn´ım obdob´ı uved’me letn´ı prázdniny, kde docház´ı k zv´ yˇsen´ı poptávky z d˚ uvodu turistického ruchu. Zde jsme stejnˇe jako u mˇes´ıˇcn´ı závislosti vyuˇzili kubické spliny. Rozm´ıstˇen´ı i poˇcet uzl˚ u jsme pouˇzili stejné jako v modelu ECB, tj. 16 uzl˚ u s um´ıstˇen´ım (0), 1, 23, 99, 121, 140, 166, 202, 207, 213, 220, 226, 231, 234, 239, 246 a 276. S poˇcátkem modelovaného roku v posledn´ım u ńorovém dnu. D˚ uleˇzitá je zde velká koncentrace uzl˚ u okolo Vánoc a konce roku, aby se zachytil prudk´ y nár˚ ust poptávky po obˇeˇzivu v tomto obdob´ı. Oproti tomu je jen nˇekolik uzl˚ u ve zb´ yvaj´ıc´ı ˇcásti roku, kde kol´ısán´ı ˇrady je menˇs´ı. Zvˇetˇsen´ım poˇctu uzl˚ u se dá sice dosáhnout lepˇs´ı kvality proloˇzen´ı, ale zároveˇ n se zhorˇs´ı pˇredpov´ıdac´ı schopnost modelu. Parametr rozptylu stochastické ˇca´sti oznaˇc´ıme jako σ32 . St´ atn´ı sv´ atky a dny volna Vliv státn´ıch svátk˚ u a dn˚ u volna je modelován deterministick´ ymi vysvˇetluj´ıc´ımi promˇenn´ ymi. Tato komponenta má tvar γt = zt4 δt4 , kde zt4 je v podstatˇe indikátor daného svátku a δt4 je konstantn´ı vektor urˇcuj´ıc´ı samotnou velikost vlivu svátku na poptávku po obˇeˇzivu. ych dn˚ u volna a modelovan´ ych Rozmˇer vektoru zt4 záleˇz´ı na poˇctu modelovan´

42

4 0.4

3

0.2

2 1

m

0

m

0 -1

-0.2

-2 -0.4 -3 1

2

4

3 d

5

(a) Sezónnost dn˚ u v t´ ydnu • d - index dne v t´ ydnu

1

3

5

7

9 11 13 15 17 19 21 23 d

(b) Mˇes´ıˇcn´ı sezónnost • d - index dne v mˇes´ıci

10 0.25 8 0.2 6 0.15 m

m 4 2

0.1

0

0.05

-2

0

-4 -4

1 25 50 75 100 125 150 175 200 225 250 275 d

(c) Roˇcn´ı sezónnost • d - index dne v roce

-3

-2

-1 d

0

1

2

(d) Vliv pevného svátku • d - index dne od svátku 15

15 12.5

10

10 m

m

7.5

5

5

0

2.5

-5

0 -4

-3

-2

-1

0 d

1

2

3

4

1 25 50 75 100 125 150 175 200 225 250 275 d

(e) Vliv krachu IPB • d - index dne od krachu (f) Celkov´ y vliv sezónnost´ı • d - index dne IPB v roce

Obrázek 5.1: Pr˚ ubˇehy sezónn´ıch komponent • m - mld. Kˇc

43

dn˚ u pˇred, respektive po jednotlivém svátku, v naˇsem pˇr´ıpadˇe 1 × 34. Hodnota na vˇsech pozic´ıch vektoru zt4 nab´ yvá bud’ 1, jedná-li se o modelovan´ y den, nebo poˇcet v´ yskyt˚ u modelovaného dne v roce − , poˇcet zb´ yvaj´ıc´ıch dn˚ u v roce P 4 aby bylo zajiˇstˇeno t+275 ze zv´ yˇsená poptávka pˇri k=t zk = 0. Interpretace je, ˇ svátku se rovnomˇernˇe rozloˇz´ı u ´bytkem poptávky v pr˚ ubˇehu zbytku roku. Vzhledem k v´ ypoˇcetn´ı nároˇcnosti modelu jsou Svátek práce, Den osvobozen´ı, Den slovansk´ ych vˇerozvˇest˚ u Cyrila a Metodˇeje a Den upálen´ı mistra Jana Husa, Den ˇceské státnosti, Den vzniku samostatného ˇceskoslovenského státu a Den boje za svobodu a demokracii uvaˇzovány jako v´ıcenásobné opakován´ı jednoho pevného svátku. Svátek Nov´ y rok Velikonoce Vánoce pevn´ y svátek

Dn˚ u pˇred 4 4 4 4

Dn˚ u po 4 4 4 2

Polynom posunut´ı (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + w1 B + . . . + w6 B 6 )B −4


Tabulka 5.1: Poˇcty dn˚ u modelovan´ ych pˇred a po jednotliv´ ych svátc´ıch

Mimoˇ r´ adn´ e ud´ alosti Stejnˇe jako v ARIMA a GARCH modelu jsme zahrnuli dvˇe mimoˇrádné události, krach IPB a problém Y2K, které mˇely vliv na v´ yˇsi poptávky po obˇeˇzivu. Zde ale vzhledem k rekurzivn´ı povaze modelu maj´ı tyto události roku 2000 pouze minimáln´ı vliv na souˇcasnost. Jejich forma v modelu je analogická státn´ım svátk˚ um a dn˚ um volna, pouze s rozd´ılem, ˇze 4 záporn´ y v´ yraz sloˇzek vektoru zt urˇcen´ y k vynulován´ı vlivu bˇehem roku se zde nahrad´ı 0, protoˇze nemá ˇzádn´ y reáln´ y smysl. Oznaˇc´ıme zt5 a δt5 . Událost IPB Y2K

Dn˚ u pˇred 4 10

Dn˚ u po 4 4

Polynom posunut´ı (w0 + w1 B + . . . + w8 B 8 )B −4 (w0 + w1 B + . . . + w14 B 14 )B −10


Tabulka 5.2: Poˇcty dn˚ u modelovan´ ych pˇred a po mimoˇra´dn´ ych událostech

5.2.2

Konstrukce kubick´ ych splin˚ u

C´ılem kubick´ ych splin˚ u je sestrojit funkci g(x), která aproximuje neznámou funkci f (x) tak, ˇze je sloˇzen´ım polynomiáln´ıch funkc´ı diferencovateln´ ych do ˇrádu 3 a

44

má prvn´ı dvˇe derivace spojité. Funkce f (x) je zadána souˇradnicemi uzl˚ u (x†i , yi† ), i = 0, 1, . . . , k. Mnoˇzinu (x†0 , . . . , x†k ) nazveme rozm´ıstˇen´ım uzl˚ u. 4 3 2 1 m

0 -1 -2 -3 3

9

11

15

19

23

d

Obrázek 5.2: Pouˇzité rozm´ıstˇen´ı uzl˚ u mˇes´ıˇcn´ı komponenty • m - mld. Kˇc • d - denn´ı index Nejdˇr´ıve pˇredpokládejme, ˇze souˇradnice uzl˚ u jsou známé, a tedy pevnˇe dané a rozm´ıstˇen´ı je v rostouc´ım poˇrad´ı x†0 < x†1 < · · · < x†k . Oznaˇcme hj = x†j − x†j−1 , dj (x) prvn´ı a d2j (x) druhou derivaci funkce splin˚ u g(x) na † † 2 † intervalu [xj−1 , xj ] a aj = dj (xj ), j = 1, . . . , k. Kubické spliny jsou smˇes´ı polynom˚ u † † tˇret´ıho ˇra´du, a tedy jejich druhá derivace na intervalu [xj−1 , xj ] je lineárn´ı funkce d2j (x)

x − x†j−1 x†j − x aj−1 + aj , = hj hj

j = 1, . . . , k.

V´ yraz pro prvn´ı derivaci a samotnou funkci splin˚ u dostaneme integrac´ı s podm´ın† † † † kami g(xj ) = yj a g(xj−1 ) = yj−1 . Plat´ı " # " # † † 1 (xj − x)2 hj 1 (x − xj−1 )2 hj dj (x) = − − aj−1 + − aj 2 hj 6 2 hj 6 a

x†j − x † [(x†j − x)2 − h2j ] aj−1 + yj−1 gj (x) = − x) 6hj hj x − x†j−1 † [(x − x†j−1 )2 − h2j ] † aj + yj , +(x − xj−1 ) 6hj hj (x†j

(5.21)

kde x†j−1 ≤ x ≤ x†j , j = 1, . . . , k. Posledn´ı rovnici lze zapsat ve vektorovém zápisu gj (x) = rj (x)0 y † + sj (x)0 a,

45

kde a = (a0 , . . . , ak )0 , y † = (y0† , . . . , yk† ) a vektory rj (x) a sj (x) maj´ı na vˇsech pozic´ıch 0, kromˇe j-té a (j + 1)-vé pozice, které odpov´ıdaj´ı pˇr´ısluˇsn´ ym vahám † v rovnici (5.21). Vyuˇzijeme podm´ınku spojitosti prvn´ı derivace dj (xj ) = dj+1 (x†j ) a dostaneme soustavu k − 1 rovnic hj hj+1 aj−1 + 2aj + aj+1 = hj + hj+1 hj + hj+1 † † 6yj−1 6yj† 6yj+1 − + , hj (hj + hj+1 ) hj hj+1 hj+1 (hj + hj+1 )

j = 1, 2, . . . , k − 1.

K soustavˇe k − 1 rovnic o k + 1 neznám´ ych je tˇreba pˇridat dalˇs´ı dvˇe podm´ınky, aby byla ˇreˇsitelná. Jelikoˇz chceme periodickou funkci splin˚ u, mus´ı b´ yt y0† = yk† ,

d21 (x†0 ) = d2k (x†k ).

d1 (x†0 ) = dk (x†k ),

Z druhé podm´ınky plyne † 6yk−1 hk h1 6yk† 6y1† ak−1 + 2ak + a1 = + − . h1 + hk h1 + hk hk (h1 + hk ) h1 hk h1 (h1 + hk )

Posledn´ı podm´ınku pˇrep´ıˇseme jako a0 = ak , a t´ım máme soustavu k rovnic o k neznám´ ych, jej´ıˇz ˇreˇsen´ı v maticovém tvaru je a = P −1 Q y † , kde

       P =              Q=      

h2 h1 +h2

0 0 0 .. .

h3 h3 +h4

0 0 .. .

h4 h4 +h5

0 .. .

··· ··· ··· ··· ··· ...

0

0 0

0 0

··· ···

6 h2 (h1 +h2 ) − h26h3 6 h3 (h3 +h4 )

0

2 h2 h2 +h3

2

h1 h1 +hk

− h16h2

6 h2 (h2 +h3 )

0 0 0 .. . 0

6 h1 (h1 +hk )

0 0 .. .

0 0

0 h3 h2 +h3

2

0 0 0 0 0 .. .

h1 h1 +h2

2

hk hk−1 +hk

0 0 0 0 .. .

hk h1 +hk

             

2 0 0 0 0 0 .. .

6 h1 (h1 +h2 )

0 .. .

··· ··· ··· ··· ··· ...

0 0

··· ···

6 − hk−1 hk

6 hk (hk−1 +hk ) − h16hk

6 h3 (h2 +h3 ) − h36h4 6 h4 (h4 +h5 )

6 hk (h1 +hk )

0 0 0 0 .. .

             

Nyn´ı m˚ uˇzeme vyjádˇrit funkci splin˚ u jako lineárn´ı funkci y † gj (x) = wj (x)0 y˜† ,

46

(5.22)

kde x†j−1 ≤ x ≤ x†j ,

wj (x)0 = r˜j (x)0 + s˜j (x)0 P −1 Q,

j = 1, . . . , k

a k × 1 vektory y˜† , r˜j a s˜j jsou upravené vektory y † , rj a sj vzhledem k y0† = yk† . Parametr x je bod, ve kterém chceme funkci g(x) vypoˇc´ıtat, a j znaˇc´ı, mezi kter´ ymi dvˇema uzly tento bod leˇz´ı. Rovnici(5.22) pˇrep´ıˇseme do tvaru sezónn´ı komponenty stavového modelu γt = wt0 γ † , kde t je ˇcasov´ y okamˇzik, ve kterém se model nacház´ı. Parametr t urˇcuje, ve které ˇcásti obdob´ı (napˇr. kter´ y den v roce) se ˇcasová ˇrada nacház´ı, a j dále nep´ıˇseme, protoˇze je zˇrejmé, mezi kter´ ymi uzly t leˇz´ı. Dále potˇrebujeme zajistit, ˇze celkov´ y † vliv sezónnosti za obdob´ı bude nula. Proto vypust´ıme posledn´ı prvek vektoru γ a v´ ysledn´ y vektor oznaˇc´ıme δ. Následnˇe i-t´ y prvek 1 × (k − 1) vektoru zt vypoˇcteme jako zti = wti − wtk w∗i /w∗k , i = 1, . . . , k − 1, pro vˇsechny dny t v modelovaném obdob´ı, a kde w∗i je i-t´ y prvek vektoru w∗ =

X

wl ,

napˇr. pro rok

w∗ =

276 X

wl .

l=1

obdob´ı

Nakonec pˇrestaneme uvaˇzovat souˇradnice uzl˚ u za pevné, ale umoˇzn´ıme jejich promˇenlivost v ˇcase. Pˇredpokládejme, ˇze sleduj´ı náhodnou procházku † γt† = γt−1 + χt ,

kde χt je vektor seriálnˇe nekorelovan´ ych náhodn´ ych veliˇcin s nulovou stˇredn´ı hodnotou a kovarianˇcn´ı matic´ı Var(χt ) = σχ2 (I − (1/w∗0 w∗ )w∗ w∗0 ). Parametr σχ2 udává rychlost s jakou se spline m˚ uˇze mˇenit a celá kovarianˇcn´ı matice 0 zaruˇcuje, ˇze w∗ χt = 0. Stejnˇe jako u splinu nemˇenného v ˇcase oznaˇc´ıme vektor γt† zkrácen´ y o posledn´ı prvek δt a dostaneme definitivn´ı tvar periodického kubického splinu promˇenného v ˇcase γt = zt δt .

5.2.3

Alternativn´ı modelov´ an´ı mˇ es´ıˇ cn´ı sez´ onnosti

Jako alternativn´ı zp˚ usob zachycen´ı vlivu dne v mˇes´ıci jsme m´ısto proloˇzen´ı kubick´ ych splin˚ u pouˇzili trigonometrické funkce. Sezónn´ı komponenta pˇri délce sezóny s má pak tvar bs/2c X γt = γjt , j=1

47

kde

∗ γj,t+1 = γjt cos λj + γjt sin λj + ωjt , ∗ ∗ ∗ γj,t+1 = −γjt sin λj + γjt cos λj + ωjt ,

j = 1 . . . , bs/2c,

∗ λj = 2πj a ωjt , ωjt jsou nezávislé náhodné veliˇciny s rozdˇelen´ım N(0, σω2 ). V modelu s poptávky po obˇeˇzivu je s = 23 a bylo by tedy tˇreba k = 11 subkomponent γjt , ale lze ukázat (viz [13]), ˇze k = 8 plnˇe postaˇcuje. Pˇredchoz´ı rovnice zap´ıˇseme do stavového modelu se stavov´ ym vektorem ∗ ∗ αt = (γ1t , γ1t , γ2t , . . . , γ8t )

a maticemi Zt = (1, 0, 1, 0, . . . , 1, 0), Tt = diag(C1 , C2 , . . . , C8 ), Rt = I16 , Qt = σω2 I16 , kde

µ Cj =

cos λj sin λj − sin λj cos λj

¶ ,

j = 1, . . . , 8.

Pˇripojen´ı komponenty do modelu poptávky po obˇeˇzivu se provede znám´ ym zp˚ usobem. 4 3 2 1 m

0 -1 -2 -3 1

3

5

7

9

11 13 15 17 19 21 23 d

Obrázek 5.3: Mˇes´ıˇcn´ı sezónnost - trigonometrické funkce • m - mld. Kˇc • d - denn´ı index

5.3

Verifikace modelu

Pˇri testován´ı správnosti a kvality modelu jsme pouˇzili ˇradu nástroj˚ u, jejichˇz ˇc´ıselné hodnoty jsou uvedeny v souhrnné tabulce na konci kapitoly vˇenované STS modelu.

48

Protoˇze se ve velké ˇca´sti jedná o jiˇz popsaná kritéria u ARIMA a GARCH model˚ u, je uveden pouze jejich v´ yˇcet. Hladinu v´ yznamnosti u test˚ u mˇejme α = 5%. Pˇredpokladem stavového modelu je, ˇze chybové ˇcleny εt a ηt jsou normálnˇe rozdˇelené a seriálnˇe nezávislé s konstantn´ımi rozptyly. Vzhledem k tˇemto pˇredpoklad˚ um mus´ı b´ yt standardizované chyby jednokrokové pˇredpovˇedi vt et = √ , Ft

t = 1, . . . , n,

také normálnˇe rozdˇelené a seriálnˇe nezávislé s jednotkov´ ym rozptylem. Normalitu jsme testovali pomoc´ı Kolmogorov-Smirnovova testu, Lillieforsova testu, JarqueBerova testu a vykreslen´ım histogramu a QQ - grafu (Obr. 5.7 a 5.6). Vˇsechny testy, podobnˇe jako u ostatn´ıch model˚ u, hypotézu normáln´ıho rozloˇzen´ı zam´ıtaj´ı. Pro ovˇeˇren´ı nemˇennosti rozptylu jsme pouˇzili Goldfeld-Quandt˚ uv test, kter´ y nulovou hypotézu nezam´ıtá. Seriáln´ı nekorelovanost jsme analyzovali pomoc´ı korelogram˚ u (Obr. 5.4 a 5.5) a Ljung-Boxova testu pro rezidua a také pro jejich druhé mocniny. V´ ysledky lze charakterizovat jako uspokojivé, v´ yznamnˇejˇs´ı korelovanost se objevuje jen u roˇcn´ıch frekvenc´ı. Pˇri srovnáván´ı kvality jednotliv´ ych model˚ u jsme pouˇzili vˇerohodnost modelu, AIC a BIC kritérium. Zároveˇ n jsme provádˇeli anal´ yzu grafu nenormovan´ ych rezidu´ı jednokrokové predikce (Obr. 5.8 a 5.9), zvláˇstˇe odchylek pˇresahuj´ıc´ıch 1 mld. Kˇc. D˚ uleˇzitá vlastnost je samozˇrejmˇe pˇredpov´ıdac´ı schopnost, mˇeˇrená stˇredn´ı kvadratickou chybou pˇredpovˇedi. Pˇri jednokrokové pˇredpovˇedi za rok 2003 jsme zaznamenali tˇri odchylky v´ yraznˇe vˇetˇs´ı neˇz 1 mld. Kˇc. Jedná se o 14.3.2003, 9.5.2003 a 7.11.2003. Posledn´ı dvˇe zˇrejmˇe souvis´ı s urˇcitou ztrátou informace zp˚ usobenou vloˇzen´ım chybˇej´ıc´ıch dn´ı, prvn´ı z nich se nepodaˇrilo uspokojivˇe vysvˇetlit. Vˇsechny testy, histogram a QQ-graf byly provedeny na normovan´ ych rezidu´ıch {et }. Normovaná rezidua uvaˇzujeme bez nulov´ ych rezidu´ı v ˇcasov´ ych okamˇzic´ıch vloˇzen´ ych chybˇej´ıc´ıch dn´ı, protoˇze by doˇslo k neˇza´douc´ımu zkreslen´ı v´ ysledk˚ u. Ostatn´ı statistiky provád´ıme na nenormovan´ ych rezidu´ıch {vt } vˇcetnˇe nul na m´ıstech vloˇzen´ ych chybˇej´ıc´ıch dn´ı, protoˇze zde jde bud’ pouze o relativn´ı vztah mezi veliˇcinami, nebo tato varianta nab´ız´ı vˇetˇs´ı pˇrehlednost v´ ystup˚ u.

5.4

Predikce

Po zkonstruován´ı a ovˇeˇren´ı modelu nyn´ı chceme ˇcasovou ˇradu poptávky po obˇeˇzivu pˇredpov´ıdat. Pˇresnˇeji formulováno naj´ıt nejlepˇs´ı odhad yn+j , j > 0, za pˇredpokladu znalosti y1 , . . . , yn , ve smyslu stˇredn´ı kvadratické chyby. Oznaˇcme jej y n+j a ze standardn´ıch znalost´ı v´ıme, ˇze y n+j = E(yn+j | y). Pro j = 1 se jedná o v´ ystup

49

0.8

0.6

ACF

0.4

0.2

0

−0.2

0

50

100

150 l

200

250

300

Obrázek 5.4: ACF rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı

0.8

0.6

PACF

0.4

0.2

0

−0.2

0

20

40

60

80

100

120

140

160

180

200

l

Obrázek 5.5: PACF rezidu´ı (ˇcervenˇe), hranice v´ yznamnosti (modˇre) • l - zpoˇzdˇen´ı 0.999 0.997 0.99 0.98 0.95 0.90 0.75 p

0.50

0.25 0.10 0.05 0.02 0.01 0.003 0.001 −3

−2

−1

0

1

2

data

Obrázek 5.6: QQ-graf rezidu´ı

50

3

4

70

60

50

40

30

20

10

0 −4

−3

−2

−1

0

1

2

3

4

5

Obrázek 5.7: Histogram rezidu´ı Kalmanova filtru, uvaˇzme tedy j > 1. Jelikoˇz yn+j = Zn+j αn+j + εn+j , dostáváme y n+j = Zn+j E(αn+j | y). Oznaˇcme an+j = E(αn+j | y) a P n+j = E[(an+j − αn+j )(an+j − αn+j )0 | y]. Souˇcasnˇe odvod´ıme F n+j = E[(y n+j − yn+j )(y n+j − yn+j )0 | y] = E[{Zn+j (an+j − αn+j ) − εn+j }{Zn+j (an+j − αn+j ) − εn+j }0 ] 0 = Zn+j P n+j Zn+j + Hn+j . Dále staˇc´ı urˇcit rekurze pro v´ ypoˇcet an+j a P n+j . Z αn+j = Tn+j−1 αn+j−1 + Rn+j−1 ηn+j−1 plyne an+j = Tn+j−1 E(αn+j−1 | y) = Tn+j−1 an+j−1 a

P n+j = E[(an+j − αn+j )(an+j − αn+j )0 | y] 0 = Tn+j−1 E[(an+j−1 − αn+j−1 )(an+j−1 − αn+j−1 )0 | y]Tn+j−1 0 0 +Rn+j−1 E(ηn+j−1 ηn+j−1 )Rn+j−1 0 0 = Tn+j−1 Pn+j−1 Tn+j−1 + Rn+j−1 Qn+j−1 Rn+j−1 .

Takto m˚ uˇzeme postupovat pro n + j, n + j − 1, . . . , n + 2, aˇz an+1 = an+1 a P n+1 = Pn+1 . Vˇsimneme si, ˇze odvozené rovnice jsou rovnice Kalmanova filtru, jestliˇze vk = 0 a Kk = 0 pro k = n + 1, . . . , j − 1. Pro predikci lze tedy pouˇz´ıt standardn´ı Kalman˚ uv filtr, u kterého povaˇzujeme pozorován´ı yt , t > n za chybˇej´ıc´ı. V naˇsem modelu poptávky po obˇeˇzivu jsme provedli pˇredpovˇedi s krokem j = 1, 5, 10, tj. stejnˇe jako u ostatn´ıch model˚ u. Stˇredn´ı kvadratická chyba pˇredpovˇedi se s velikost´ı kroku zvyˇsuje, ˇc´ıselné hodnoty jsou uvedeny v tabulce 5.3. Vzhledem k lepˇs´ım v´ ysledk˚ um uvaˇzujeme v dalˇs´ıch kapitolách jen STS model s mˇes´ıˇcn´ı sezónnost´ı modelovanou kubick´ ymi spliny.

51

2 1 m

0 -1 -2 50

100

150 d

200

250

Obrázek 5.8: Rezidua jednokrokové predikce v roce 2003 (ˇcernˇe), hranice ±1 mld. Kˇc (ˇcervenˇe) • m - mld. Kˇc • d - denn´ı index

2 1 m

0 -1 -2 200

400 d

600

800

Obrázek 5.9: Rezidua jednokrokové predikce v letech 2001-2003 (ˇcernˇe), hranice ±1 mld. Kˇc (ˇcervenˇe) • m - mld. Kˇc • d - denn´ı index

52

Model Nenormovaná rezidua Stˇredn´ı hodnota Rozptyl Normovaná rezidua Stˇredn´ı hodnota Rozptyl ˇ Sikmost ˇ catost Spiˇ Kolmogorov-Smirnov test p-hodnota Lilliefors test p-hodnota Jarque-Bera test p-hodnota Goldfeld-Quandt test p13-hodnota p31-hodnota Ljung-Box test Q(5) p-hodnota Q(10) p-hodnota Q(22) p-hodnota Q(261) p-hodnota Ljung-Box test Q(5) p-hodnota Q(10) p-hodnota Q(22) p-hodnota Q(261) p-hodnota Log. vˇerohodnost AIC BIC Chyby pˇredpovˇed´ı RMSE krok 1 ˇ RMSE krok 1 v CNB RMSE krok 5 RMSE krok 10

s kubick´ ymi spliny

s trigon. funkcemi

-0.01422 0.25276

-0.01317 0.27983

0.02114 0.76767 0.10879 4.73427 Zam´ıtá H0 0.00857 Zam´ıtá H0 0.048 Zam´ıtá H0 2.6312x10−14 Nezam´ıtá H0 0.64623 0.35377 rezidua 0.08078 0.03555 0.01281 0 kvadrát rezidu´ı 7.2217x10−6 0.0004 0.01707 0.66782 -1.79053x106 6487.77 6488.47 [mld. Kˇc] 0.47018 0.43787 1.0122 1.45305

0.00041 0.59438 -0.01544 5.39046 Zam´ıtá H0 0.00053 Zam´ıtá H0 NaN1 Zam´ıtá H0 0 Nezam´ıtá H0 0.79926 0.20074 0.00024 0.00039 0.00011 0 0.00032 0.00578 0.13595 0.96007 730.341 -2.3201 -1.6168 0.50443 0.43787 1.22207 1.56971

1 Procedura

lillietest v Matlab v nˇekter´ ych pˇr´ıpadech nevrac´ı p-hodnotu, ale pouze udává v´ ysledek testu.

Tabulka 5.3: Tabulka v´ ysledk˚ u model˚ u

53

Kapitola 6 Kombinace model˚ u V pˇredcházej´ıc´ıch kapitolách jsme zkonstruovali tˇri r˚ uzné modely pro predikci poptávky po obˇeˇzivu. Nyn´ı bychom mohli podle urˇcitého kritéria, napˇr. velikosti stˇredn´ı kvadratické chyby pˇredpovˇedi, nebo kombinac´ı v´ıce kritéri´ı rozhodnout, kter´ y model povaˇzujeme za nejlepˇs´ı a ten dále pouˇz´ıvat pro predikci ˇrady. T´ım bychom ale ˇza´dn´ ym zp˚ usobem nevyuˇzili dodateˇcnou informaci obsaˇzenou v ostatn´ıch modelech, danou jinou strukturou tˇechto model˚ u, ˇci pouˇzit´ım jin´ ych vysvˇetluj´ıc´ıch exogenn´ıch promˇenn´ ych. Zvláˇstˇe pokud by chyby pˇredpovˇed´ı dvou model˚ u byly negativnˇe korelované, mohli bychom jejich kombinac´ı dosáhnout v´ yraznˇe lepˇs´ı pˇredpovˇedi ˇcasové ˇrady oproti pˇredpovˇedi obou jednotliv´ ych model˚ u. Námi pouˇzitá metoda kombinac´ı model˚ u je zaloˇzená na regresn´ı metodˇe navrˇzené v [8]. Pro pozorovanou ˇcasovou ˇradu yt , t = 1, . . . , T mˇejme dvˇe alternativn´ı pˇredpovˇedi fta a ftb , t = 1, . . . , T . Kombinované pˇredpovˇedi konstruujeme ve tvaru yˆt = αfta + (1 − α)ftb , pro t = 1, . . . , T . Váhy alternativn´ıch model˚ u nalezneme pomoc´ı metody nejmenˇs´ıch ˇctverc˚ u. Vzhledem ke skuteˇcnosti dlouhé inicializace Kalmanova filtru u STS modelu a jeho následné adaptace po pˇrechodu z rozˇs´ıˇreného na standardn´ı Kalman˚ uv filtr hledáme optimáln´ı váhy pouze na ˇcásti ˇcasové ˇrady pˇripadaj´ıc´ı na rok 2002. Model/Pˇredpovˇed’ 1 krok vpˇred 5 krok˚ u vpˇred sts - arima 0.434918 0.904458 sts - garch 0.448697 0.901731 arima - garch 0.476364 1.24903

10 krok˚ u vpˇred 1.23999 1.29987 1.84795

Tabulka 6.1: RMSE pˇredpovˇed´ı kombinovan´ ych model˚ u v mld. Kˇc

54

RMSE krok 1

sts - arima 0.435519

sts - garch 0.449314

arima - garch 0.480947

Tabulka 6.2: RMSE pˇredpovˇed´ı kombinovan´ ych model˚ u po ˇctvrtlet´ıch v mld. Kˇc

Vzhledem k tomu, ˇze ˇcasová ˇrada vykazuje sezonn´ı vlivy, nab´ız´ı se moˇznost nestanovovat váhy pˇredpovˇed´ı model˚ u stejné pro vˇsechna t = 1, . . . , T , ale promˇenlivé po obdob´ıch, a t´ım vyuˇz´ıt pˇr´ıpadné lepˇs´ı vlastnosti jednotliv´ ych model˚ u v tˇechto obdob´ıch. My jsme provedli dˇelen´ı na ˇctvrtlet´ı, a tedy urˇcili váhy α1Q , . . . , α4Q . Kvalita predikce se ovˇsem nezv´ yˇsila. Podrobnˇejˇs´ı dˇelen´ı na mˇes´ıce jako v ECB jsme neaplikovali, protoˇze pro anal´ yzu pouˇz´ıváme pouze hodnoty roku 2002 a v jednotliv´ ych mˇes´ıc´ıch nemáme dostatek pozorován´ı. Kvalita pˇredpovˇed´ı mˇeˇrená odmocninou stˇredn´ı kvadratické chyby se kombinac´ı jednotliv´ ych model˚ u ˇ zlepˇsila, v´ ysledky uvád´ıme v tabulce 6.1 a 6.2. C´ıselné hodnoty koeficient˚ u α a dalˇs´ı statistiky regresn´ı anal´ yzy lze nalézt v CD pˇr´ıloze.

2 1.5 1 0.5 m

0 -0.5 -1 -1.5 50

100

150

200

250

d Obrázek 6.1: Rezidua jednokrokové predikce sts-arima modelu v roce 2003 (ˇcernˇe), hranice ±1 mld. Kˇc (ˇcervenˇe) • m - mld. Kˇc • d - denn´ı index

55

260

250

m 240

230

220 50

100

150

200

250

d

Obrázek 6.2: Jednokroková pˇredpovˇed’ sts-arima modelu (ˇcervenˇe) a skuteˇcnost (ˇcernˇe) v roce 2003 • m - mld. Kˇc • d - denn´ı index

265

260

255

m 250

245

240

10

20

30 d

40

50

60

Obrázek 6.3: Jednokroková pˇredpovˇed’ sts-arima modelu (ˇcervenˇe) a skuteˇcnost (ˇcernˇe) v 4. ˇctvrtlet´ı 2003 • m - mld. Kˇc • d - denn´ı index

56

Kapitola 7 Porovn´ an´ı model˚ u Jak jiˇz bylo naznaˇceno v pˇredcházej´ıc´ıch ˇcástech práce, hlavn´ım kritériem pro porovnán´ı model˚ u bude jejich pˇredpov´ıdac´ı schopnost. Pouˇzit´ ymi statistikami jsou odmocnina stˇredn´ı kvadratické chyby pˇredpovˇedi, Diebold-Marian˚ uv test pro porovnán´ı pˇresnosti pˇredpovˇed´ı a poˇcet odchylek pˇredpovˇedi vˇetˇs´ıch neˇz 1 mld. Kˇc. Nelze nav´ıc ˇr´ıci, ˇze by nˇekter´ y z model˚ u byl statisticky v´ıce správn´ y. Problémy s nenormalitou rezidu´ı nastaly u vˇsech model˚ u, ˇreˇsen´ım by mohl b´ yt pˇredpoklad jiného rozdˇelen´ı rezidu´ı, napˇr´ıklad t-rozdˇelen´ı. Zam´ıtnut´ı konstantnosti rozptylu rezidu´ı u ARIMA a GARCH modelu nepovaˇzujeme za velk´ y problém, nebot’ k nˇemu docház´ı v d˚ usledku vˇetˇs´ıho rozptylu ˇcasové ˇrady obˇeˇziva v roce 1998, kdy ˇ pravdˇepodobnˇe docházelo k adaptaci na nov´ y mˇenovˇepolitick´ y reˇzim CNB. Zkrácen´ı ˇcasové ˇrady jsme neprovedli z d˚ uvodu zachován´ı konzistence podklad˚ u jednotliv´ ych model˚ u. Stˇredn´ı kvadratická chyba, vˇcetnˇe v´ ysledn´ ych hodnot, jiˇz byla uvedena u jednotliv´ ych model˚ u, zb´ yvá tedy objasnit Diebold-Marian˚ uv test. Tento test je urˇcen pro testován´ı nulové hypotézy stejné pˇredpov´ıdac´ı schopnosti dvou model˚ u a vycház´ı z centráln´ı limitn´ı vˇety pro stacionárn´ı náhodné procesy. Diebold a Mariano pracuj´ı s obecnou chybovou funkc´ı (loss function), my jsme pouˇzili kvadratickou chybovou funkci. Testová statistika má pak tvar SDM = s

d¯ 2π fˆd (0) T

d

→ N(0, 1),

P kde d¯ = T1 Tt=1 dî , dî = eˆ2A,i − eˆ2B,i , i = 1, . . . , T , eÂ,i jsou chyby predikce modelu A, eˆB,i jsou chyby predikce modelu B, T je délka predikovaného intervalu a fˆd (0) je konzistentn´ı odhad spektráln´ı hustoty procesu {di } na frekvenci 0. Pˇredpokladem je zde stacionarita procesu {di }. Odhad 2π fˆd (0) dostaneme jako váˇzenou sumu

57

Modely/Pˇredpovˇed’ STS vs ARIMA GARCH sts-arima sts-garch arima-garch ARIMA vs GARCH sts-arima sts-garch arima-garch GARCH vs sts-arima sts-garch arima-garch sts-arima vs sts-garch arima-garch sts-garch vs arima-garch

1 krok vpˇred 0.362 0.000 0.977 0.984 0.421 0.000 0.998 0.891 0.938 1.000 1.000 1.000 0.158 0.006 0.143

5 krok˚ u vpˇred 10 krok˚ u vpˇred 0.063 0.063 0.001 0.000 0.995 0.992 0.990 0.979 0.094 0.070 0.006 0.011 0.980 0.974 0.980 0.966 0.842 0.891 1.000 1.000 1.000 1.000 1.000 0.999 0.532 0.049 0.030 0.026 0.028 0.034

Tabulka 7.1: V´ ysledky Diebold-Marianova testu

dostupn´ ych v´ ybˇerov´ ych autokovarianc´ı γˆd (τ ), µ ¶ (T −1) X τ l 2π fˆd (0) = γˆd (τ ), S(T ) τ =−(T −1)

τ kde funkce l( S(T ) pˇriˇrazuje váhy autokovarianˇcn´ı funkci na jednotliv´ ych zpoˇzdˇe) n´ıch a je naz´ yvána jádrem (kernel) a St udává zpoˇzdˇen´ı, po kterém autokovariance zanedbáváme (truncation lag). Stejnˇe jako Diebold a Mariano jsme pouˇzili jednotkové (obdéln´ıkové) jádro ¶ ½ µ τ 1 pro | S(T |≤1 τ ) = l 0 jinak S(T )

a S(T ) = (k − 1) v pˇr´ıpadˇe k-krokové pˇredpovˇedi, nebot’ dle teoretick´ ych v´ ysledk˚ u jsou chyby optimáln´ı k-krokové pˇredpovˇedi nejv´ yˇse (k − 1)-závislé. V praxi toto sice z mnoha d˚ uvod˚ u nemus´ı platit, pˇresto pouˇzit´ı této hodnoty S(T ) dává dobré v´ ysledky, podrobnosti lze naj´ıt v [5]. Vypoˇctené hodnoty Diebold-Marianova testu, pˇrevedené na pravdˇepodobnosti, jsou uvedeny v tabulce 7.1. Jestliˇze je hodnota malá, napˇr´ıklad pro hladinu testu α = 5% menˇs´ı neˇz 0.025, potom zam´ıtáme nulovou hypotézu ve prospˇech modelu A, tj. pˇredpovˇed’ modelem A je statisticky pˇresnˇejˇs´ı neˇz modelem B. Je-li naopak hodnota velká, zam´ıtáme nulovou hypotézu ve prospˇech modelu B. V tabulce 7.2 následuj´ı ˇcetnosti pˇrekroˇcen´ı hranice 1 mld. Kˇc u chyb pˇredpovˇed´ı jednotliv´ ymi modely. Z v´ ysledk˚ u je patrné, ˇze nejlepˇs´ım modelem je sts-arima model, a to pro vˇsechny provedené horizonty pˇredpovˇed´ı. Pˇri jednokrokové pˇredpovˇedi je pˇredpovˇed’

58

Model/Pˇredpovˇed’ 1 krok vpˇred 5 krok˚ u vpˇred STS 11 63 ARIMA 9 95 GARCH 44 122 sts-arima 10 59 sts-garch 8 62 arima-garch 10 96

10 krok˚ u vpˇred 98 130 185 89 87 138

Tabulka 7.2: Poˇcet chyb predikce vˇetˇs´ıch neˇz 1 mld. Kˇc

sts-arima modelem statisticky v´ yznamnˇe pˇresnˇejˇs´ı neˇz vˇsechny ostatn´ı modely, s vyj´ımkou sts-garch modelu. Nejhorˇs´ım modelem je samotn´ y model GARCH, kter´ y je signifikantnˇe horˇs´ı neˇz ostatn´ı modely a také poˇcet chyb predikce pˇresahuj´ıc´ıch 1 mld. Kˇc je v´ yraznˇe vˇetˇs´ı neˇz u ostatn´ıch model˚ u. Vˇsechny kombinace model˚ u vykazuj´ı lepˇs´ı v´ ysledky neˇz pˇr´ısluˇsné samostatné modely. Nejvˇetˇs´ı pˇr´ınos z kombinace model˚ u lze pozorovat u sts-garch a sts-arima modelu. Naopak u kombinace arima-garch je zlepˇsen´ı pˇredpovˇedi minimáln´ı, nebot’ struktura obou model˚ u je velmi podobná. Budeme tedy povaˇzovat sts-arima za náˇs fináln´ı model. Pˇri anal´ yze predikce (Obr. 6.1) je zaj´ımavé, ˇze k nejvˇetˇs´ım odchylkám predikce od skuteˇcnosti nedocház´ı v obdob´ı Vánoc, ale v kvˇetnu a ˇcervnu, nebot’ právˇe kvalita zachycen´ı vlivu obdob´ı Vánoc se podaˇrila zv´ yˇsit kombinac´ı model˚ u. Nejvˇetˇs´ı odchylka, témˇeˇr 2 mld. Kˇc, nastala 14.3.2003 a je v této ˇca´sti roku zcela ojedinˇelá. V ostatn´ıch mˇes´ıc´ıch je predikce modelem velmi pˇresná. Srovnáme-li fináln´ı model s triviáln´ı pˇredpovˇed´ı, kdy jako pˇredpovˇed’ yt+1 pˇri znalosti y1 , . . . , yt bereme yt , která má RMSE jednokrokové pˇredpovˇedi 1.09739, sn´ıˇzili jsme chybu predikce na necel´ ych 40% oproti této triviáln´ı pˇredpovˇedi.

59

Kapitola 8 Z´ avˇ er C´ılem této práce bylo vytvoˇren´ı stochastického modelu, kter´ y by mohl b´ yt pouˇzit pro pˇredpov´ıdán´ı nejv´ yznamnˇejˇs´ıho autonomn´ıho faktoru ovlivˇ nuj´ıc´ıho likviditu ˇ T´ımto faktorem je objem obˇeˇziva, neboli objem hotovostn´ıch penˇeˇzn´ıho trhu v CR. penˇez drˇzen´ ych mimo centráln´ı banku. Jak jiˇz bylo zm´ınˇeno, tato ˇcasová ˇrada vykazuje velké mnoˇzstv´ı sezónn´ıch závislost´ı, které se navzájem kombinuj´ı a zp˚ usobuj´ı ˇcetné komplikace pˇri tvorbˇe modelu. Pˇri nalezen´ı relevantn´ıch vysvˇetluj´ıc´ıch promˇenn´ ych jsme vycházeli ze ˇ zkuˇsenost´ı expert˚ u CNB, pˇredchoz´ıch model˚ u ([2] a [12]) a vlastn´ı anal´ yzy ˇcasové ˇrady. Nároˇcnost správného zachycen´ı sezónnost´ı potvrzuje napˇr´ıklad 75 nenulov´ ych parametr˚ u u ARIMA modelu. Zkonstruovali jsme tˇri r˚ uzné modely. Jsou jimi ARIMA model, zaloˇzen´ y na Box-Jenkinsovˇe metodologii, GARCH model, kter´ ym jsme se snaˇzili postihnout nˇekteré rysy rezidu´ı ARIMA modelu a STS model zaloˇzen´ y na zcela odliˇsném pˇr´ıstupu k ˇcasov´ ym ˇradám pomoc´ı stavov´ ych prostor˚ u a Kalmanov´ ych rekurz´ı. Vˇsechny tˇri modely dávaj´ı dobré v´ ysledky. Nejlepˇs´ı z nich je STS model, kter´ y po zkombinován´ı s ARIMA modelem vytváˇr´ı fináln´ı model. Oproti tomu aplikace GARCH modelu nesplnila oˇcekáván´ı. Fináln´ı sts-arima model dosáhl nepatrnˇe ˇ menˇs´ı chyby jednokrokové predikce neˇz expertn´ı odhad provádˇen´ y v CNB. V práci byl zkonstruován stochastick´ y model obˇeˇziva, kter´ y dává v´ ysledky na u ´rovni dosud pouˇz´ıvan´ ych odhad˚ u. Pˇridá-li se k v´ ysledk˚ um modelu expertn´ı ’ zkuˇsenost, zkoriguje se pˇredpovˇed modelem v kritick´ ych obdob´ıch roku a pˇri mimoˇra´dn´ ych událostech. Vznikne tak lepˇs´ı odhad poptávky po penˇez´ıch, kter´ y ˇ umoˇzn´ı efektivnˇejˇs´ı aplikaci mˇenovˇepolitick´ ych nástroj˚ u CNB. Vytvoˇren´ y model nab´ız´ı zároveˇ n predikci s delˇs´ım ˇcasov´ ym horizontem, která bude nutná v souvislosti se zaveden´ım jednotné evropské mˇeny euro.

60

Literatura [1] Brockwell, P.J., Davis, R.A. (2001): Introduction to Time Series and Forecasting. Springer-Verlag, New York. [2] Cabrero, A., Camba-Mendez, G., Hirsch, A. (2001): Modeling the Daily Series of Banknotes in Circulation in the Context of the Liquidity Management Operations of the ECB. ECB Working Paper no. 142. [3] Cipra, T. (1986): Anal´ yza ˇcasov´ ych ˇrad s aplikacemi v ekonomii. SNTL, Praha. [4] De Jong, P. (1991): The Diffuse Kalman Filter. The Annals of Statistics 19, 1073-1083. [5] Diebold, F.X., Mariano, R.S. (1995): Comparing Predictive Accuracy. Journal of Business and Economic Statistics 13, 253-263. [6] Durbin, J., Koopman, S.J. (2001): Time Series Analysis by State Space Methods. Oxford University Press. [7] Franke, J., Härdle, W., Hafner, C. (2001): Einf¨ uhrung in die Statistik der Finanzmärkte. Springer-Verlag, Heidelberg. [8] Granger, C.V.J., Ramanathan, R. (1984): Improved Methods of Forecasting. Journal of Forecasting 3, 197-204. [9] Harvey, A.C. (1989): Forecasting, Structural Time Series Models and the Kalman Filter. Cambridge University Press. [10] Harvey, A.C., Koopman, S.J., Riani M. (1997): The Modeling and Seasonal Adjustment of Weekly Observations. Journal of Business and Economic Statistics 15, 354-368. [11] Hlaváˇcek, M. (2005): The Application of Structured Feedforward Neural Networks to the Modeling of the Daily Series of Currency in Circulation. Lecture Notes In Computer Science, Springer-Verlag, New York.

61

[12] Hlaváˇcek, M. (2004): Modelován´ı objemu obˇeˇziva. Bakaláˇrská práce, FJFI ˇ CVUT, Praha. [13] Pierce, D.A., Grupe, M.R., Cleveland, W.P. (1984): Seasonal Adjustment of the Weekly Monetary Aggregates: A Model-Based Approach. Journal of Business & Economic Statistics 2, 260-270.

62

UNIVERZITA KARLOVA V PRAZE Matematicko fyzikální fakulta

Recommend Documents