library(readxl)
Spot <- read_excel("~/Nemecko data.xlsx")
mydata<-ts(Spot[,-1],start=2001,frequency = 4)
mdata<-cbind("cons"=log(mydata[,1]),"di"=log(mydata[,2]))Na vysvetlenie postupu tvorby modelu s korekčným členom ECM a následne vytvorenie prognózy na jeho základe využijeme modelovanie spotreby v stálych cenách z roku 2010 v Nemecku. Na modelovanie použijeme dáta z roku 2001. V modeli budeme pracovať so štvrťročnými údajmi, aby sme zvýšili počet pozorovaní. Údaje nebudú očistené od sezónnych vplyvov.
Pre prvotnú špecifikáciu premenných využíva ekonometria výber na základe ekonomickej teórie. Teória rozlišuje niekoľko dôležitých významných faktorov, ktoré nám vplývajú na spotrebu. Tieto údaje si preto adekvátne upravíme:
Takisto na odhadnutie modelu použijeme zlogaritmovanú spotrebu, ktorá sa bude rovnať \(log(CON)\) https://ec.europa.eu/eurostat/databrowser/view/namq_10_gdp__custom_21384760/default/table
Všeobecný tvar modelu, s ktorými budeme pracovať:
\[CON_t = f(DI_t)\]
Dátová množina obsahuje súbor 100 údajov a 2 premenných získaných z Eurostatu od 1. kvartálu roku 2001 do 4. kvartálu 2025 vrátane. Jednotlivé premenné sú označené následovne:
cons - spotreba domácností (v mil. Eur)
di - disponibilný príjem (v mil. Eur)
Prvý krok pri vytvárani modelu s korekčným členom je analýza stacionarity respektíve určenie rádu použitých premenných. Otestujeme premenné CON, DI
#testovanie kointegracie
library(urca)
summary(ur.df(mdata[,"cons"],type="trend",lags=10,selectlags = "BIC"))
summary(ur.df(mdata[,"cons"],type="trend",lags=7))
summary(ur.df(mdata[,"cons"],type="drift",lags=10,selectlags = "BIC"))
summary(ur.df(mdata[,"cons"],type="drift",lags=7))
summary(ur.df(diff(mdata[,"cons"]),type="trend",lags=10,selectlags = "BIC"))
summary(ur.df(diff(mdata[,"cons"]),type="trend",lags=6))Na zistenie rádu integrácie premennej spotreba použijeme rozšírený Dickeyho-Fullerov test. Postup začína nájdením najvhodnejšieho tvaru modelu pomocou bayesovského informačného kritéria (BIC), ktoré vyberie 7 oneskorení. Začíname modelom s trendom (a 7 členmi riešiacimi autokoreláciu podľa BIC).
Výsledkom je model, v ktorom \(\phi_3\) = 4.2409 je menšie ako \(\phi_{3crit}\) = 6.49. Nulovú hypotézu \(H_0:\delta=\alpha_1=0\) nemôžeme na päťpercentnej hladine významnosti zamietnuť a model s trendom nie je vhodný na testovanie.
Pokračujeme modelom s driftom (a 7 členmi riešiacimi autokoreláciu podľa BIC).
Výsledkom je model, v ktorom \(\phi_1\) = 2.5529 je menšie ako \(\phi_{1crit}\) = 4.71. Nulovú hypotézu \(H_0:\delta=\alpha_0=0\) príjimame a na päťpercentnej hladine významnosti je model s driftom nevhodný na testovanie.
Záver, ktorý z toho vyplýva je, že skúmaný časový rad spotreby je random walk bez driftu, teda rad nestacionárny. Na zistenie rádu integrácie pokračujeme opakovaním postupu pre rad prvých diferencií spotreby.
Model prvej diferencie s trendom (a 6 členmi riešiacimi autokoreláciu podľa BIC) má tvar:
\[ \Delta^2 cons = \alpha_0 + \alpha_1t + \delta\Delta cons_{t-1} + \beta_1 \Delta^2 cons_{t-1} + \beta_2 \Delta^2 cons_{t-2} + \beta_3 \Delta^2 cons_{t-3} + \beta_4 \Delta^2 cons_{t-4} + \beta_5 \Delta^2 cons_{t-5} + \beta_6 \Delta^2 cons_{t-6} + u_i \]
Z výsledku vyplýva, že \(\phi_3\) = 17.6038 je väčšie ako \(\phi_{3crit}\) = 6.49. Nulovú hypotézu \(H_0:\delta=\alpha_1=0\) môžeme na päťpercentnej hladine významnosti zamietnuť a model s trendom je vhodný na testovanie.
Hodnota \(\tau_δ\) sa rovná -5.933 a je menšia ako kritická hodnota \(\tau_{crit}\) = -3.45, môžeme zamietnuť nulovú hypotézu \(H_0:\delta=\alpha_0=0\) na päťpercentnej hladine významnosti.
Záver, ktorý z toho vyplýva je, že skúmaný časový rad prvých diferencií spotreby je stacionárny.
Na základe použitej sekvencie DF testu vyplýva, že časový rad spotreby je nestacionárny, ale časový rad prvých diferencií spotreby je už stacionárny. Preto môžeme konštatovať, že časový rad spotreby je integrovaný rádu jedna, teda zapisujeme \[ log(cons)∼I(1) \]
library(urca)
summary(ur.df(mdata[,"di"],type="trend",lags=10,selectlags = "BIC"))
summary(ur.df(mdata[,"di"],type="trend",lags=10))
summary(ur.df(diff(mdata[,"di"]),type="trend",lags=10,selectlags = "BIC"))
summary(ur.df(diff(mdata[,"di"]),type="trend",lags=5))Začíname modelom s trendom (a 10 členmi riešiacimi autokoreláciu podľa BIC).
Výsledkom je model, v ktorom \(\phi_3=8.0449\) je väčšie ako \(\phi_{3,crit}=6.49\). Nulovú hypotézu \(H_0:\delta=\alpha_1=0\) môžeme na päťpercentnej hladine významnosti zamietnuť a model s trendom je vhodný na testovanie stacionarity.
Hodnota \(\tau_3=1.7452\). Keďže táto hodnota nie je menšia ako kritická hodnota \(\tau_{crit}=-3.45\), nemôžeme zamietnuť nulovú hypotézu o existencii jednotkového koreňa na päťpercentnej hladine významnosti.
Záver, ktorý z toho vyplýva je, že skúmaný časový rad disponibilného príjmu je nestacionárny a obsahuje jednotkový koreň. Keďže rad nie je stacionárny, na zistenie rádu integrácie pokračujeme opakovaním postupu pre rad prvých diferencií disponibilného príjmu.
Model prvej diferencie disponibilného príjmu s trendom (a 5 členmi riešiacimi autokoreláciu podľa BIC) má tvar:
\[ \Delta^2 di_t=\alpha_0+\alpha_1t+\delta\Delta di_{t-1}+\beta_1\Delta^2di_{t-1}+\beta_2\Delta^2di_{t-2}+\beta_3\Delta^2di_{t-3}+\beta_4\Delta^2di_{t-4}+\beta_5\Delta^2di_{t-5}+u_t \]
Z výsledku vyplýva, že \(\phi_3=14.7408\) je väčšie ako \(\phi_{3,crit}=6.49\). Nulovú hypotézu \(H_0:\delta=\alpha_1=0\) môžeme na päťpercentnej hladine významnosti zamietnuť a model s trendom je vhodný na testovanie.
Hodnota \(\tau_3\) sa rovná \(-5.4297\) a je menšia ako kritická hodnota \(\tau_{crit}=-3.45\), preto môžeme zamietnuť nulovú hypotézu o existencii jednotkového koreňa na päťpercentnej hladine významnosti.
Záver, ktorý z toho vyplýva je, že skúmaný časový rad prvých diferencií disponibilného príjmu je stacionárny.
Na základe použitej sekvencie DF testu vyplýva, že časový rad disponibilného príjmu je nestacionárny, ale časový rad prvých diferencií disponibilného príjmu už je stacionárny. Preto môžeme konštatovať, že časový rad disponibilného príjmu je integrovaný rádu jedna, teda zapisujeme:
\[ \log(di)\sim I(1) \]
Z uvedených výsledkov vidíme, že všetky rady tvoriace modely sú nestacionárne a integrované rádu jedna. Preto je splnená nevyhnutná podmienka pre využitie koncepcie kointegrácie, lebo rád integrácie všetkých premenných je rovnaký, pričom sa nejedná o prípad integrácie rádu nula – prípad stacionarity.
Ďalším krokom pri tvorbe modelu s korekčným členom je analýza kointegrácie. Využijeme procedúru Englea a Grangera. Na začiatok budeme pracovať so 2 premennými a predpokladať, že spotreba zavisí od disponibilného príjmu.
Odhadneme dlhodobý vzťah medzi spotebou a disponibilným príjmom.
Tvar rovnice kointegrácie:
\[COINT_1 = \beta_1 + \beta_2 di_t + u_t\] Oba parametre nám v rovnici vyšli štatisticky významné, a preto budeme pokračovať na základe rovnice dlhodobej rovnováhy stacionaritu reziduálov. Ak budú reziduály kointegračnej rovnice stacionárne, potom medzi spotrebou a disponibilným príjmom existuje dlhodobý rovnovážny stav a môžeme pristúpiť k tvorbe modelu s korekčným členom ECM.
summary(ur.df(resid(coint),type="none",lags=10, selectlags = "BIC"))
summary(ur.df(resid(coint),type="none",lags=7)) # -2.431 > -3.39 <- nezamietam H0, medzi premennymi nie je kointegracia na 5 % urovniTestujeme reziduály z rovnice \(coint1\) prostredníctvom Dickeyovho-Fullerovho testu. Odhadneme najvhodnejší model podľa BIC so siedmimi členmi riešiacimi autokoreláciu a dostaneme:
\[ \begin{aligned}\Delta coint1_t = \delta coint1_{t-1} + \beta_1 \Delta coint1_{t-1} + \beta_2 \Delta coint1_{t-2} + \beta_3 \Delta coint1_{t-3} + \beta_4 \Delta coint1_{t-4} + \beta_5 \Delta coint1_{t-5} + \beta_6 \Delta coint1_{t-6} + \beta_7 \Delta coint1_{t-7} + u_i\end{aligned} \] Hodnota \(\tau_δ\) sa rovná \(-2.431\) a je menšia ako kritická hodnota \(\tau_{crit}\) = \(-3.39\) pri dvoch premenných, nulovú hypotézu \(H_0:\delta=\alpha_0=0\) na päťpercentnej hladine významnosti zamietame.
Do modelu skúsime pridať umelú premennú \(tt\) - trend.
Tvar novej rovnice kointegrácie:
\[COINT2_1 = \beta_1 + \beta_2 di_t +tt + u_t\] v modeli nám vyšiel štatisticky nevýznamný len parameter di, ale model ako celok je štatisticky významný, preto pokračujeme testovaním reziduálov.
summary(ur.df(resid(coint2),type="none",lags=6, selectlags = "BIC"))
summary(ur.df(resid(coint2),type="none",lags=1))Odhadneme najvhodnejší model podľa BIC s jedným členom riešiacim autokoreláciu a dostaneme:
\[ \Delta coint2_t = \delta coint2_{t-1} + \beta_1 \Delta coint2_{t-1} + u_t \]
Hodnota \(\tau_\delta\) sa rovná \(-7.9997\) a je menšia ako kritická hodnota \(\tau_{crit}=-3.83\). Nulovú hypotézu o existencii jednotkového koreňa v rezíduách zamietame na päťpercentnej hladine významnosti.
Keďže premenné \(log(cons)\) a \(log(di)\) sú integrované rádu jedna a rezíduá z ich dlhodobej rovnováhy sú stacionárne, môžeme konštatovať, že medzi premennými existuje kointegrácia.
Potvrdili sme tak existenciu dlhodobého rovnovážneho vzťahu medzi spotrebou domácností a disponibilným príjmom. Na základe potvrdenej kointegrácie môžeme pristúpiť k odhadu modelu s korekčným členom (ECM).
Po potvrdení existencie kointegrácie medzi premennými spotreba a disponibilný príjem pristúpime k zostaveniu modelu s korekčným členom ECM. Na jeho zostavenie využijeme reziduály z dlhodobého vzťahu rovnováhy \(coint2\), ktoré predstavujú odchýlku od dlhodobej rovnováhy.
Do dátového súboru vložíme vektor reziduálov, z ktorých následne zobrazíme priebeh a vyberieme tie reziduály, ktoré predstavujú šoky. Preto do modelu zaradíme umelú premennú \(u21\), ktorá bude predstavovať šok z 3.Q 2021.
mdata <- ts.union(mdata, c(resid(coint2,)))
colnames(mdata) <- c("cons","di","e")
library(dynlm)
pomtt <- time(mdata[,"cons"])
u21 <- ifelse(pomtt==2021.75,1,0)
# odhad ECM
ecm3 <- dynlm(d(cons) ~ d(di) + L(e,1) +
L(d(cons),1:4) +
L(d(di),1:4) + u21,
data=mdata)
summary(ecm3)Výsledný model korekčného člena ECM má tvar:
\[ \begin{aligned} \widehat{\Delta cons_t} = &\beta_0 +\beta_1\Delta di_t +\gamma_1 e_{t-1} +\beta_2\Delta cons_{t-1} +\beta_3\Delta cons_{t-2} +\beta_4\Delta cons_{t-3} +\beta_5\Delta cons_{t-4}\\ &+\beta_6\Delta di_{t-1} +\beta_7\Delta di_{t-2} +\beta_8\Delta di_{t-3} +\beta_9\Delta di_{t-4} +\gamma_2u21 +u_t \end{aligned} \]
V modeli použijeme oneskorený korekčný člen \(e_{t-1}\), pretože pracujeme so štvrťročnými dátami a odchýlka od dlhodobej rovnováhy z predchádzajúceho obdobia ovplyvňuje aktuálnu zmenu spotreby. Zároveň v modeli využívame oneskorené hodnoty diferencií spotreby a disponibilného príjmu, ktoré zachytávajú krátkodobú dynamiku medzi premennými.
Výsledky odhadu modelu ukazujú, že korekčný člen \(e_{t-1}\) je štatisticky významný a nadobúda zápornú hodnotu \(\gamma_1=-0.6974\). To znamená, že pri odchýlke od dlhodobej rovnováhy dochádza v nasledujúcom období ku korekcii smerom späť k rovnovážnemu stavu. Štatisticky významná je taktiež premenná \(\Delta cons_{t-4}\) a umelá premenná \(u21\).
Model ako celok je štatisticky významný. Hodnota F-štatistiky = 33.08 a p-hodnota < 2.2e-16 potvrdzuje, že vysvetľované premenné ako celok významne vplývajú na krátkodobú dynamiku spotreby.
Následne otestujeme autokoreláciu reziduálov pomocou Ljung-Box testu.
Autokorelácia prvého rádu, \(Q_m\) = \(0.0204\) a p-hodnota = \(0.8864\) je väčšia ako päťpercentná hladina významnosti. Nulovú hypotézu o existencii autokorelácie prvého rádu preto nezamietame. V modeli nie je prítomná autokorelácia prvého rádu.
Ďalším krokom bude otestovanie autokorelácie štvrtého rádu. Hodnota testovacej štatistiky \(Q_m\) = \(7.2627\) a p-hodnota = \(0.1226\) je väčšia ako päťpercentná hladina významnosti. Nulovú hypotézu o existencii autokorelácie štvrtého rádu preto nezamietame.
Na základe výsledkov Ljung-Box testu môžeme konštatovať, že v modeli ECM3 nie je prítomná autokorelácia reziduálov do štvrtého rádu. Model je preto vhodný na ďalšie využitie pri tvorbe prognózy.
Na určenie prognózy disponibilného príjmu bolo využité priemerné kvartálne tempo rastu posledných období. Na základe vypočítaného tempa rastu sme zistili:
\[ \Delta di = 0,007505 \]
Následne boli dopočítané hodnoty disponibilného príjmu pre jednotlivé prognózované obdobia.
Výsledné prognózované hodnoty spotreby a disponibilného príjmu sú uvedené v tabuľke:
| Obdobie | cons | di |
|---|---|---|
| 2026 Q1 | 419 600 | 728 500 |
| 2026 Q2 | 417 600 | 734 000 |
| 2026 Q3 | 416 700 | 739 500 |
| 2026 Q4 | 416 300 | 745 100 |
Z výsledkov prognózy vyplýva, že model predpokladá mierny pokles spotreby v prvých obdobiach prognózy napriek pokračujúcemu rastu disponibilného príjmu. Tento výsledok je spôsobený najmä pôsobením korekčného člena ECM, ktorý zabezpečuje návrat spotreby späť k dlhodobej rovnováhe.