Dynamiske Prognosemodeller

Denne labben er delvis basert på k. 9 av Forecasting boka

Det er ofte at man ønsker å modellere en tidsrekke avhengig av andre variabler og størrelser. Det vil si, vi ønsker å kjøre regresjoner som:

\(y_t = \beta_0 + \beta_1 x_t + \eta\)

For eksempel, vi vil kanskje estimere effekten av høyere oljepris på fastlands-BNP. Deretter kan vi snu modellen og lage en prognose eller scenario for hva BNP blir under ulike oljepriser-baner.

Problemet er at hvis det finnes dynamikk/autokorrelasjon i \(y_t\) som ikke blir modellert, så havner det i vår feil-ledd. Og det gjør alt fra t-tester til vår usikkerhetsbånder ugyldige.

En løsning er å modellere dynamikken i \(y_t\) i tillegg til å estimere effekten av \(x_t\). Denne type modellen heter ARMAX, der X står for exogen. Men det er problematisk å tolke resultatene for en sånn modell. Du kan lese mer om hvorfor her.

En alternativ løsning er å modellere dynamikken i feilleddet. Så at man estimerer et system som kan, for eksempel se ut som:

\(y_t = \beta_0 + \beta_1 x_t + \eta\)

\(\eta_t = \theta \eta_{t-1} + \epsilon\)

Her ser vi nå at vi ikke lenger modellerer dynamikken i y-variablen, men isteden i feil-leddet, \(\eta\). Dermed kan vi estimere den ubetinget, fulle effekten av x på y.

Vi laster igjen noen pakker vi vil trenge:

library(fpp2)
library(tidyverse)
library(zoo)

Vi bruker uschange datasetteet fra fpp2 pakken som gir oss endringer i noen makroøkonomiske størrelser fra USA.

head(uschange)

##         Consumption     Income Production   Savings Unemployment
## 1970 Q1   0.6159862  0.9722610 -2.4527003 4.8103115          0.9
## 1970 Q2   0.4603757  1.1690847 -0.5515251 7.2879923          0.5
## 1970 Q3   0.8767914  1.5532705 -0.3587079 7.2890131          0.5
## 1970 Q4  -0.2742451 -0.2552724 -2.1854549 0.9852296          0.7
## 1971 Q1   1.8973708  1.9871536  1.9097341 3.6577706         -0.1
## 1971 Q2   0.9119929  1.4473342  0.9015358 6.0513418         -0.1

La oss modellere effekten av inntekt har på konsum.

konsum = uschange[, 1] #første kollone 
inntekt = uschange[, 2] #andre kollone 
ledighet = uschange[, 5] #femte kollone

konsum %>% ggtsdisplay()

inntekt %>% ggtsdisplay()

library(tseries)
adf.test(konsum)

## Warning in adf.test(konsum): p-value smaller than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  konsum
## Dickey-Fuller = -4.5318, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary

Det vi ser her er at det er noe autokorrelasjon i konsum, men serien er stasjonær. Vi kan gjette kanskje tre AR ledd ut i fra pACF figuren.
Inntekt er både stasjonær og uten særlig mye autokorrelasjon.

Så kan vi modellere

armax_fit1 = konsum %>% Arima(order = c(3,0,0), xreg=inntekt)

summary(armax_fit1)

## Series: . 
## Regression with ARIMA(3,0,0) errors 
## 
## Coefficients:
##          ar1     ar2     ar3  intercept    xreg
##       0.1256  0.2377  0.1595     0.6030  0.1972
## s.e.  0.0782  0.0718  0.0724     0.0915  0.0483
## 
## sigma^2 estimated as 0.3213:  log likelihood=-156.77
## AIC=325.55   AICc=326.02   BIC=344.94
## 
## Training set error measures:
##                       ME      RMSE       MAE      MPE     MAPE      MASE
## Training set 0.001673788 0.5591672 0.4173871 30.20659 163.0622 0.6539603
##                     ACF1
## Training set 0.005555023

vi sjekker residualene:

checkresiduals(armax_fit1)

## 
##  Ljung-Box test
## 
## data:  Residuals from Regression with ARIMA(3,0,0) errors
## Q* = 5.5454, df = 3, p-value = 0.1359
## 
## Model df: 5.   Total lags used: 8

Dette ser ganske ok ut. Nå kan vi skape vår forecast, 8 kvartaler fram i tid. Vi lager to scenariorer. En der inntektsvekst er likt sit gjennomsnitt. Og en der inntekt er 0,5 prosentpoeng lavere

h=8
fcast_mean = forecast(armax_fit1, xreg=rep(mean(inntekt),h))
fcast_low =  forecast(armax_fit1, xreg=rep(mean(inntekt)-.5,h))

autoplot(konsum) +
  ylab("% change in US consumption") +
  autolayer(fcast_mean, PI = FALSE, series = "mean") +
  autolayer(fcast_low, PI = FALSE, series = "low", alpha=.5) +
  guides(colour = guide_legend(title = "Scenario"))

Og med usikkerhets-intervaller.

autoplot(konsum) +
  ylab("% change in US consumption") +
  autolayer(fcast_low, PI = TRUE, series = "low", alpha=.5) +
  guides(colour = guide_legend(title = "Scenario"))

autoplot(konsum) +
  ylab("% change in US consumption") +
  autolayer(fcast_low, PI = TRUE, series = "mean", alpha=.5) +
  guides(colour = guide_legend(title = "Scenario"))

Oppgaver

Det kan også henne at man vil bruke lagget verdier av x-variablen i regresjonen. Kjør en regresjon der man også inkluderer lagget inntekt. Se her for hjelp.

inntektM = cbind(
    inntekt= inntekt,
    inntektL1 = stats::lag(inntekt,-1),
    inntektL2 = stats::lag(inntekt,-2),
    inntektL3 = stats::lag(inntekt,-3))

NROW(inntekt)

## [1] 187

armax_fit2 = konsum[4:187] %>% Arima(order = c(3,0,0), xreg=inntektM[4:187,1:4])

summary(armax_fit2)

## Series: . 
## Regression with ARIMA(3,0,0) errors 
## 
## Coefficients:
##          ar1     ar2     ar3  intercept  inntekt  inntektL1  inntektL2
##       0.0830  0.1576  0.1713     0.3903   0.2375     0.1643     0.0442
## s.e.  0.0762  0.0767  0.0798     0.1051   0.0472     0.0463     0.0465
##       inntektL3
##          0.0550
## s.e.     0.0441
## 
## sigma^2 estimated as 0.309:  log likelihood=-149.05
## AIC=316.1   AICc=317.13   BIC=345.03
## 
## Training set error measures:
##                        ME      RMSE       MAE      MPE     MAPE      MASE
## Training set 0.0009321757 0.5437028 0.4150936 23.95729 150.0177 0.7427494
##                    ACF1
## Training set 0.01690297

La oss også inkludere vekst i ledighetsraten, variablen som vi kalte ledighet tidligere. Er det en viktig faktor i konsum? Skap to ulike scenarior: en med høy (positiv) ledighetsvekst og en med lav (negativ) ledighetsvekst.
Bruk dataseriene fra lab 12 om oljepriser og valutakurs. Skap en dynamisk prognosemodell der man estimerer valutakurset 12 måneder fram i tid under lav og høy oljepris.

Dynamiske Prognosemodeller

Johannes Mauritzen

Høsten 2020

Oppgaver