Normalnie gdy chcemy sprawdzić czy jedna zmienna wpływa na drugą, budujemy model regresji, w którym zmienna objaśniona jest opóźniona o jakiś okres, i testujemy go. Co jednak zrobić, gdy występuje autokorelacja, jak to ma miejsce często w danych ekonomicznych? Dodawanie do takiego modelu autoregresji jest teoretycznie błędne, bo wiadomo, że w normalnej regresji zmienne objaśniające powinny być niezależne od siebie.
Z pomocą przychodzi nam test Grangera. Test sprawdza przyczynowość w sensie Grangera. C. Granger zdefiniował przyczynowość w następujący sposób [1]:
Mówimy, że Xt powoduje Yt, jeśli jesteśmy w stanie lepiej przewidzieć Yt przy użyciu wszystkich dostępnych informacji, niż gdyby wykorzystano wszystkie informacje oprócz Xt.
Jest to definicja bardzo ogólna i abstrakcyjna. Zauważmy, że w tej definicji nie ma opóźnienia Y względem X. Dlaczego? Łatwiejsze zrozumienie przyczynowości Grangera w kontekście procesu ARMA podają autorzy w [2]:
X wywołuje Y, jeśli wcześniejsze informacje o X przewidują zachowanie Y lepiej niż same wcześniejsze informacje o Y.
Innymi słowy Xt i Yt są zmiennymi, które wpływają na siebie w teraźniejszości - natychmiast, ale Granger oddziela zbiór informacji o tych zmiennych, który może mieć przesunięcie w czasie lub ściślej mówiąc zbiór z okresu t zawiera informacje z okresu t-1. Taki podział pozwala rozróżnić zwykłą korelację od przyczynowości. Powiedzmy, że mamy rozkład jazdy: autobus X przejeżdża część trasy od godz. 8 do 9, a autobus Y kolejną część trasy od 9 do 10. Gdy X dojeżdża do swojego końca, Y rusza. Jednak X nie jest przyczyną Y, a jedynie z nim koreluje, ponieważ Y może ruszyć niezależnie od tego czy X dojeżdża czy nie. A teraz wyobraźmy sobie, że jazda obydwu ma się odbywać w tym samym czasie. Gdyby X bardzo przyspieszył i dojeżdżał już o 8:45, to Y musiałby ruszyć o 8:45, a nie 9:00, żeby nie zostać uderzonym przez X. W ten sposób X(t) staje się przyczyną Y(t). Ale ocenić to możemy dopiero, gdy mamy informację, że X wcześniej z jakiegoś powodu przyspieszył, np. nie było nikogo na drodze.
W świecie gospodarki sytuacja wygląda podobnie, choć jest dużo bardziej złożona. Produkcja w jednym kraju będzie wpływać na produkcję w innym kraju, bo np. wzrost produkcji w kraju X to wzrost dochodów, który można przeznaczyć na nowe inwestycje w kraju Y.
Zajmijmy się teraz stroną praktyczną, czyli zastosujmy test Grangera w języku R. Na stronie OECD możemy ściągnąć statystyki produkcji miesięcznej dla krajów OECD. Wybieramy Time -> Select Date Range -> Monthly -> np. od 2001 r. Eksportujemy text file (CSV). Interesują mnie ogólnie 4 obszary:
1. Polska
2. UE
3. USA
4. całe OECD
Niestety wskaźniki na OECD są słabo opisane i trzeba się trochę dokopać, żeby zrozumieć, co jest czym. Przykładowo produkcji nie można oddzielać od sprzedaży, bo zgodnie z zasadami rachunkowości tylko sprzedana produkcja się liczy do rachunków. Podstawową miarą jest indeks produkcji. Żeby odkryć co wchodzi w skład produkcji, możemy poniżej kliknąć w Structural Analysis (STAN) Databases. Pojawi się tablica ze składem industry. Możemy kliknąć w Industry, żeby obejrzeć cały skład. Okazuje się, że mieści się tu cała gospodarka, w tym usługi.
Jeśli chodzi o same dane, to są tu dwie istotne cechy. Po pierwsze jest to indeks w odniesieniu do określonego roku. Nie jest ważne do jakiego, ale należy zmienną przekształcić w stopę zmian (tutaj m/m). Po drugie nie są to dane surowe, ale poddane obróbce - usunięto z nich sezonowość za pomocą X12-ARIMA i TRAMO-SEATS , stąd oznaczenie sa (ang. seasonal adjustment). To drugie może stanowić podstawę różnicy między danymi OECD a GUS.
No dobrze, na początek mamy dwie zmienne:
- stopy zmian produkcji w Polsce, m/m: stopa_pl
- stopy zmian produkcji w USA, m/m: stopa_us
Kod:
plot(stopa_pl, lwd=2, col="blue", main = "Industry, change% m/m (OECD)")
lines(stopa_us, lwd=2, col="red")
legend("topleft", legend=c("Poland", "USA"), col=c("blue", "red"), lty = 1, lwd = 2)
Wykonajmy korelację krzyżową między nimi:
ccf(stopa_pl, stopa_us)
Występuje silna autokorelacja bieżąca (0,6), ale też dość mocna korelacja z pospieszeniem zmiennej stopa_pl o 1 okres. Czyli to produkcja USA jest opóźniona o 1 okres w stos. do polskiej. Wydawałoby się, że powinno być na odwrót - że Polska podąża za USA. Tutaj dostajemy coś mało intuicyjnego. Może problemem są autokorelacje? Sprawdźmy PACF (autokorelacje cząstkowe):
pacf(stopa_pl)
Okazuje się, że autokorelacji 1 rzędu brak. Interesujące, że występuje ujemna korelacja co drugi miesiąc, ale na razie zostawiamy to.
Zastosujmy test Grangera. Użyjemy do tego funkcji grangertest w pakiecie lmtest. Pakiet ten nie jest niczym nowym, bo używamy go zawsze do testów homoskedastyczności. Trzeba nadmienić, że test Grangera dotyczy zmiennych stacjonarnych, więc warto na początku wykonać test stacjonarności, ale to pomijam.
Przetestujmy to co jest mało intuicyjne - czy Polska (x) wpływa na USA (y).
Hipoteza zerowa: Polska nie wpływa na USA.
Hipoteza alternatywna: Polska wpływa na USA.
Jak widzieliśmy w ccf występowało tylko jedno opóźnienie w korelacji, więc przypuszczamy, że w teście wystarczy wskazać order=1:
grangertest(x=stopa_pl, y=stopa_us, order=1)
Granger causality test
Model 1: stopa_us ~ Lags(stopa_us, 1:1) + Lags(stopa_pl, 1:1)
Model 2: stopa_us ~ Lags(stopa_us, 1:1)
Res.Df Df F Pr(>F)
1 265
2 266 -1 17.5 0.000038 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Hipoteza zerowa zostaje odrzucona. Nasze wątpliwości potwierdzają się - zmiany produkcji w Polsce powodują w sensie Grangera zmiany w produkcji USA. Sprawdźmy w drugą stronę - czy produkcja USA (x) wpływa na produkcję Polski (y).
Hipoteza zerowa: USA nie wpływa po miesiącu na Polskę.
Hipoteza alternatywna: USA wpływa po miesiącu na Polskę.
grangertest(x=stopa_us, y=stopa_pl, order=1)
Efekt:
Granger causality test
Model 1: stopa_pl ~ Lags(stopa_pl, 1:1) + Lags(stopa_us, 1:1)
Model 2: stopa_pl ~ Lags(stopa_pl, 1:1)
Res.Df Df F Pr(>F)
1 265
2 266 -1 0.08 0.78
Nie można odrzucić hipotezy zerowej. Okazuje się, że USA nie wpływa w sensie Grangera na Polskę. To wydaje się pozbawione sensu, ale pamiętajmy o silnej bieżącej korelacji między krajami. Możliwe, że opóźnienia występują wewnątrz miesiąca.
Pozostaje zagadka, dlaczego Polska wyprzedza produkcję USA. Zwróćmy uwagę, że Polska jest częścią Unii Europejskiej, która z kolei może wpływać na gospodarkę USA. Tak więc Polska może mieć pośredni wpływ na USA. Sprawdźmy tę hipotezę.
Porównajmy ze sobą:
- stopy zmian produkcji w UE, m/m: stopa_eu
- stopy zmian produkcji w USA, m/m: stopa_us
Nie pokazuję już efektu krzyżowej korelacji, bo jest bardzo podobny do tej pierwszej - występuje głównie korelacja 1 rzędu.
Hipoteza zerowa: UE nie wpływa po miesiącu na USA.
Hipoteza alternatywna: UE wpływa po miesiącu na USA.
Wykonujemy test Grangera. Ponieważ obecnie nie ma jeszcze danych dla UE za czerwiec, odejmuję ostatnią wartość dla USA.
grangertest(stopa_eu, stopa_us[-length(stopa_us)], order=1)
Granger causality test
Model 1: stopa_us[-length(stopa_us)] ~ Lags(stopa_us[-length(stopa_us)], 1:1) + Lags(stopa_eu, 1:1)
Model 2: stopa_us[-length(stopa_us)] ~ Lags(stopa_us[-length(stopa_us)], 1:1)
Res.Df Df F Pr(>F)
1 264
2 265 -1 43.5 0.00000000023 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Hipoteza zerowa zostaje odrzucona, zatem UE istotnie wpływa na USA. Sprawdźmy to samo na odwrót:
Hipoteza zerowa: USA nie wpływa po miesiącu na UE.
Hipoteza alternatywna: USA wpływa po miesiącu na UE.
grangertest(stopa_us[-length(stopa_us)], stopa_eu, order=1)
Granger causality test
Model 1: stopa_eu ~ Lags(stopa_eu, 1:1) + Lags(stopa_us[-length(stopa_us)], 1:1)
Model 2: stopa_eu ~ Lags(stopa_eu, 1:1)
Res.Df Df F Pr(>F)
1 264
2 265 -1 3.48 0.063 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Hipoteza o tym, że USA wpływa na UE musi zostać odrzucona na poziomie 5% istotności, ale już na poziomie 10% nie. Można przypuszczać, że jakiś niewielki feedback, tzn. sprzężenie zwrotne między USA a UE także występuje.
Zaproponowana hipoteza, że Polska wpływa pośrednio na USA poprzez UE jest już częściowo potwierdzona. Żeby w pełni to potwierdzić przetestujmy związek między UE a Polską.
Hipoteza zerowa: UE nie wpływa po miesiącu na Polskę.
Hipoteza alternatywna: UE wpływa po miesiącu na Polskę.
Podobnie jak z USA, odejmuję ostatnią wartość dla Polski i robię test:
grangertest(stopa_eu, stopa_pl[-length(stopa_pl)], order=1)
Granger causality test
Model 1: stopa_pl[-length(stopa_pl)] ~ Lags(stopa_pl[-length(stopa_pl)], 1:1) + Lags(stopa_eu, 1:1)
Model 2: stopa_pl[-length(stopa_pl)] ~ Lags(stopa_pl[-length(stopa_pl)], 1:1)
Res.Df Df F Pr(>F)
1 264
2 265 -1 16.5 0.000065 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Hipoteza zerowa - odrzucona. Produkcja UE wpływa na produkcję Polski. Odwrotnie:
Hipoteza zerowa: Polska nie wpływa po miesiącu na UE.
Hipoteza alternatywna: Polska wpływa po miesiącu na UE.
grangertest(stopa_pl[-length(stopa_pl)], stopa_eu, order=1)
Granger causality test
Model 1: stopa_eu ~ Lags(stopa_eu, 1:1) + Lags(stopa_pl[-length(stopa_pl)], 1:1)
Model 2: stopa_eu ~ Lags(stopa_eu, 1:1)
Res.Df Df F Pr(>F)
1 264
2 265 -1 1.09 0.3
Hipoteza zerowa - zostaje. Produkcja Polski nie wpływa na produkcję UE z opóźnieniem miesięcznym.
W sumie potwierdziliśmy hipotezę, że Polska wpływa na USA pośrednio przez UE.
Zagadka rozwiązana, chociaż pozostaje jeszcze jedna ciekawa zależność, mianowicie między całym OECD a USA, UE i Polską. Podobnie jak w przypadku UE, dane dla całego OECD pojawiają się z opóźnieniem, stąd podobnie trzeba odjąć ostatnie dane USA i Polski. Ponieważ pojawi się tu sporo kombinacji, nie będę już opisywał każdej hipotezy, a podam sam kod, efekt wraz z wnioskami:
OECD -> USA:
grangertest(stopa_oecd, stopa_us[-length(stopa_us)], order=1)
OECD powoduje zmiany w gospodarce USA po miesiącu (duża istotność, p = 0).
USA -> OECD:
grangertest(stopa_us[-length(stopa_us)], stopa_oecd, order=1)
USA powoduje zmiany w gospodarce OECD po miesiącu na poziomie 5% istotności (p = 0.022).
OECD -> UE:
grangertest(x=stopa_oecd, y=stopa_eu, order=1)
OECD nie powoduje zmiany w gospodarce UE po miesiącu (p = 0.36).
UE -> OECD:
grangertest(x=stopa_eu, y=stopa_oecd, order=1)
UE powoduje zmiany w gospodarce OECD po miesiącu (p = 0.0027).
OECD -> Polska:
grangertest(stopa_oecd, stopa_pl[-length(stopa_pl)], order=1)
OECD nie powoduje zmiany w gospodarce Polski po miesiącu na poziomie istotności 5% (p = 0.1).
Polska -> OECD:
grangertest(stopa_pl[-length(stopa_pl)], stopa_oecd, order=1)
Polska nie powoduje zmiany w gospodarce OECD po miesiącu (p = 0.86).
Dostajemy interesujący wniosek. Zmiany w UE wpływają po miesiącu na OECD, w tym na Polskę. Polska nie wpływa na OECD, mimo że wpływa pośrednio na USA, a USA stanowi dużą część OECD. Jednocześnie istnieje sprzężenie zwrotne między USA a OECD. W końcu, UE indukuje zmiany w OECD, ale nie na odwrót.
Wszystkie odkryte zależności w tym artykule ilustruje poniższy schemat:
Strzałka od USA do UE jest przerywana, bo jest to słaba zależność (sprzężenie zwrotne) oraz ma inny kształt, bo omija Polskę, a wskazuje na Unię jako całość.
Literatura:
[1] Granger C.W.J. Investigating causal relations by econometric models and cross spectral methods. Econometrica. 1969;37:424–438.
[2] Amornbunchornvej C., Zheleva E., Berger-Wolf T.Y. Variable-lag Granger Causality for Time Series Analysis. Machine Learning (cs.LG); Econometrics (econ.EM); Quantitative Methods, 2019.
---------------------------
Cały kod:
if (require(lmtest)==FALSE) {
install.packages("lmtest")
library("lmtest")
}
if (require("zoo")==FALSE) {
install.packages("zoo")
library("zoo")
}
if (require("tidyverse")==FALSE) {
install.packages("tidyverse")
library(tidyverse)
}
#zamieniam "\" na "/"
sciezka = r"(C:\Users\Documents\R\testy)"
sciezka = gsub("\\", "/", sciezka, fixed=T)
# albo
# sciezka = gsub("\\\\", "/", sciezka)
# ustawiamy folder roboczy
setwd(sciezka)
fr = 12
nazwaPliku = "MEI_REAL.csv"
if (grepl(";", readLines(nazwaPliku, n=1))) {
plik = read.csv(nazwaPliku, sep=";", dec=",")
} else if (grepl(",", readLines(nazwaPliku, n=1))) {
plik = read.csv(nazwaPliku, sep=",", dec=".")
}
produkcja = "Production of total industry sa, Index"
kraj1 = "Poland"
kraj2 = "United States"
kraj3 = "European Union – 27 countries (from 01/02/2020)"
kraj4 = "OECD - Total"
kraje = c(kraj1, kraj2, kraj3, kraj4)
moje_formaty = c("%d.%m.%Y", "%Y.%m.%d", "%d-%m-%Y", "%Y-%m-%d", "%d/%m/%Y", "%Y/%m/%d", "%b-%Y", "%B-%Y")
for (kraj in kraje) {
filtered <- subset(plik, plik$Subject == produkcja & plik$Country == kraj)
daty <- filtered[["Time"]]
daty = parse_date_time(daty, orders=moje_formaty)
rok = as.numeric(format(daty, "%Y"))
mc = as.numeric(format(daty, "%m"))
wartosci <- filtered[["Value"]]
stopy1m <- diff(wartosci)/wartosci[-length(wartosci)]
if (kraj==kraje[1]) {
stopa1 <- stats::lag(ts(stopy1m, start = c(rok[1], mc[1]), frequency= fr), -1)
} else if (kraj==kraje[2]) {
stopa2 <- stats::lag(ts(stopy1m, start = c(rok[1], mc[1]), frequency= fr), -1)
} else if (kraj==kraje[3]) {
stopa3 <- stats::lag(ts(stopy1m, start = c(rok[1], mc[1]), frequency= fr), -1)
} else if (kraj==kraje[4]) {
stopa4 <- stats::lag(ts(stopy1m, start = c(rok[1], mc[1]), frequency= fr), -1)
}
}
stopa_pl = as.numeric(stopa1)
stopa_us = as.numeric(stopa2)
stopa_eu = as.numeric(stopa3)
stopa_oecd = as.numeric(stopa4)
# test grangera
grangertest(x=stopa_pl, y=stopa_us, order=1)
grangertest(x=stopa_us, y=stopa_pl, order=1)
grangertest(x=stopa_us[-length(stopa_us)], y=stopa_eu, order=1)
grangertest(x=stopa_eu, y=stopa_us[-length(stopa_us)], order=1)
grangertest(x=stopa_eu, y=stopa_pl[-length(stopa_pl)], order=1)
grangertest(x=stopa_pl[-length(stopa_pl)], y=stopa_eu, order=1)
grangertest(x=stopa_oecd, y=stopa_us[-length(stopa_us)], order=1)
grangertest(x=stopa_us[-length(stopa_us)], y=stopa_oecd, order=1)
grangertest(x=stopa_oecd, y=stopa_eu, order=1)
grangertest(x=stopa_eu, y=stopa_oecd, order=1)
grangertest(x=stopa_oecd, y=stopa_pl[-length(stopa_pl)], order=1)
grangertest(x=stopa_pl[-length(stopa_pl)], y=stopa_oecd, order=1)
