Regresjon med tidsrekke-data II, fasit
Johannes Mauritzen
Høsten 2018
Noen pakker som vi vil bruke
library(foreign) #for å laste inn data i stata format (read.dta)
library(dynlm) #for dynamiske regresjonsmodeller## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(stargazer) #for fine tabeller##
## Please cite as:## Hlavac, Marek (2018). stargazer: Well-Formatted Regression and Summary Statistics Tables.## R package version 5.2.2. https://CRAN.R-project.org/package=stargazerlibrary(ggplot2) #for elegante plotter
library(zoo)
library(reshape2) #for å omformattere data.frame til long-formatt (melt)Vi begynner med å laste inn data om gjennomsnitts-priser og -avkastning på New York Stock Exchange:
nyse <- read.dta("http://fmwww.bc.edu/ec-p/data/wooldridge/nyse.dta")
head(nyse)## price return return_1 t price_1 price_2 cprice cprice_1
## 1 49.75 NA NA 1 NA NA NA NA
## 2 51.44 3.3969820 NA 2 49.75 NA 1.6899990 NA
## 3 52.05 1.1858490 3.3969820 3 51.44 49.75 0.6100006 1.6899990
## 4 52.28 0.4418819 1.1858490 4 52.05 51.44 0.2299995 0.6100006
## 5 54.24 3.7490489 0.4418819 5 52.28 52.05 1.9600030 0.2299995
## 6 53.76 -0.8849619 3.7490489 6 54.24 52.28 -0.4800034 1.9600030Vi formatterer datasette som zoo
nyse_zoo <- zoo(nyse, order.by = nyse$t)
head(nyse_zoo)## price return return_1 t price_1 price_2 cprice cprice_1
## 1 49.75 NA NA 1 NA NA NA NA
## 2 51.44 3.3969820 NA 2 49.75 NA 1.6899990 NA
## 3 52.05 1.1858490 3.3969820 3 51.44 49.75 0.6100006 1.6899990
## 4 52.28 0.4418819 1.1858490 4 52.05 51.44 0.2299995 0.6100006
## 5 54.24 3.7490489 0.4418819 5 52.28 52.05 1.9600030 0.2299995
## 6 53.76 -0.8849619 3.7490489 6 54.24 52.28 -0.4800034 1.9600030og kjører noen dynamiske regresjoner:
return_mod1 <- dynlm(return ~ L(return), data=nyse_zoo)
return_mod2 <- dynlm(return ~ L(return) + L(return, 2), data= nyse_zoo)
return_mod3 <- dynlm(return ~L(return) + L(return, 2) + L(return, 3), data=nyse_zoo)
stargazer(return_mod1, return_mod2, return_mod3, type="text")##
## ===============================================================================
## Dependent variable:
## -----------------------------------------------------------
## return
## (1) (2) (3)
## -------------------------------------------------------------------------------
## L(return) 0.059 0.060 0.061
## (0.038) (0.038) (0.038)
##
## L(return, 2) -0.038 -0.040
## (0.038) (0.038)
##
## L(return, 3) 0.031
## (0.038)
##
## Constant 0.180** 0.186** 0.179**
## (0.081) (0.081) (0.082)
##
## -------------------------------------------------------------------------------
## Observations 689 688 687
## R2 0.003 0.005 0.006
## Adjusted R2 0.002 0.002 0.001
## Residual Std. Error 2.110 (df = 687) 2.112 (df = 685) 2.114 (df = 683)
## F Statistic 2.399 (df = 1; 687) 1.659 (df = 2; 685) 1.322 (df = 3; 683)
## ===============================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01Dette er noen enkle autoregressive regresjoner. “L” står for lag. Det man tester med denne regresjonen er om avkastning på New York Stock Exchange en dag kan si noe om avkastning neste dag og dagen etter det. Dette er definisjonen av persistent. Hvis en serie ikke er persistent, så er det ingen sammenheng - avkastningen en dag er uavhengig av avkastningen dagen før og dagen før det.
Her er det ingen lagget variabler som er signifikant (hvordan vet man det?), og vi kan dermed ikke si at serien er persistent.
Nå laster vi inn pris-data på GE-aksjen
ge_stock <- read.csv("https://jmaurit.github.io/anvendt_macro/data/GE.csv")
head(ge_stock)## Date Open High Low Close Adj.Close Volume
## 1 2012-09-17 21.93 22.05 21.90 22.05 18.52329 79283200
## 2 2012-09-18 21.99 22.24 21.96 22.24 18.68290 38831200
## 3 2012-09-19 22.30 22.49 22.25 22.43 18.84251 40462600
## 4 2012-09-20 22.19 22.47 22.12 22.43 18.98641 43677000
## 5 2012-09-21 22.55 22.69 22.46 22.53 19.07106 66551700
## 6 2012-09-24 22.40 22.45 22.30 22.36 18.92716 36799800ge_stock["Date"]<- as.Date(ge_stock$Date, format="%Y-%m-%d")
ggplot(ge_stock, aes(Date, Close)) +
geom_line()
Vi kan se at tidsrekken av prisene er “highly persistent” - det vil si at prisen i dag er en god indikasjon på hva prisen imorgen kommer til å være i morgen. Tidsrekken er preget av høy auto-korrelasjon og variansen endrer seg over tid. Dette gjør det sånt at man ikke kan bruke OLS på en konsistent måte. Men hvis vi transformerer datarekken til daglig avkastning - det vil si prosent endring per dag, så får vi en datarekke som vi kan bruke.
ge_stock_zoo <- zoo(ge_stock, order.by=ge_stock$Date)
ge_ret <- diff(log(as.numeric(ge_stock_zoo$Close)))
ge_ret_zoo <- zoo(ge_ret, order.by=ge_stock$Date)
plot(ge_ret_zoo)
Da kan vi kjøre noen dynamiske regresjoner. Her er flere AR modeller.
ret_model1 <- dynlm(ge_ret_zoo ~ L(ge_ret_zoo))
ret_model2 <- dynlm(ge_ret_zoo~L(ge_ret_zoo) + L(ge_ret_zoo,2))
ret_model3 <- dynlm(ge_ret_zoo~L(ge_ret_zoo) + L(ge_ret_zoo,2) + L(ge_ret_zoo,3))
stargazer(ret_model1, ret_model2, ret_model3, type="text")##
## =====================================================================================
## Dependent variable:
## -----------------------------------------------------------------
## ge_ret_zoo
## (1) (2) (3)
## -------------------------------------------------------------------------------------
## L(ge_ret_zoo) 0.041 0.043 0.042
## (0.028) (0.028) (0.028)
##
## L(ge_ret_zoo, 2) -0.057** -0.056**
## (0.028) (0.028)
##
## L(ge_ret_zoo, 3) -0.013
## (0.028)
##
## Constant 0.0001 0.0001 0.0001
## (0.0003) (0.0003) (0.0003)
##
## -------------------------------------------------------------------------------------
## Observations 1,257 1,256 1,255
## R2 0.002 0.005 0.005
## Adjusted R2 0.001 0.003 0.003
## Residual Std. Error 0.011 (df = 1255) 0.011 (df = 1253) 0.011 (df = 1251)
## F Statistic 2.095 (df = 1; 1255) 3.057** (df = 2; 1253) 2.102* (df = 3; 1251)
## =====================================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01Her er L en gyldig forkortelse for Lag
Her har vi kjørt tre regresjoner, en AR(1), AR(2) og AR(3) model av avkastninger på GE aksjer. Det vi finner er at koeffisienten på den første laggen er positiv, men ikke signifikant. Den andre laggen er negative og signifikant på 0,05 prosent nivå. F-testen forteller oss at lag1 og lag2 er signifikant sammen. Vi kan tolke koeffisientene på den måten at hvis avkastningen er +100 en dag, så vil det predikere at det blir +4,3 neste dag, men følgende dag går det tilbake -5,7. Serien er kanskje litt persistent: en høy avkastning en dag, kan bety noe høyere avkastning følgende dag, men det blir også en form av “regression to the mean”, at serien deretter går litt ned igjen.
Random Walk
Den enkleste form av en persistent tidsrekke er en random walk. Dette kan skrives:
\[y_t = y_{t-1} + e_t\]
Intuitivt - prisen i dag er lik prisen i går pluss en tilfeldig (“stokastisk”) hopp, som vi representerer med \(e_t\). Random walk og noen litt mer kompliserte modeller blir brukt til å modellere oljepriser, aksjepriser, og også noen makroøkonomiske variabler.
Vi kan simulere random walk enkelt med r:
#lager en matrise som er 50x30
y<-matrix(,nrow=50, ncol=30)
#fyller den ved bruke av en for-loop
for(s in 1:30){
e <- rnorm(50) #her er vår stokastisk variabel - vi simulerer 50 trekk fra en normal-fordeling
y[,s]<-cumsum(e) #y er da modellert som summen av alle "sjokkene", e.
}
y = data.frame(y)
y["t"]<-1:50
y_l=melt(data.frame(y), id.var="t")
ggplot(y_l, aes(x=t, y=value, group=variable)) +
geom_line(alpha=.5)
I figuren som viser 30 simulerte tidsrekker av 50 perioder, kan man se at variansen øker over tid og det beveger seg ikke mot en bestemt gjennomsnittsverdi - noe som gjør OLS ugyldig. Regresjonen er også ikke-stasjonært (som vi kommer tilbake til), som også gjør det vanskelig å jobbe med tidsrekken uten å transformere den.
I beregningen, så regner vi ut \(y_t\) som en kumulativt sum av \(e_t\) ’ene. det vil si \(y_t = e_0 + e_1 + ... + e_{t-1} + e_t\). Kan du forklare hvorfor?
En annen form av random walk er random walk med drift. Det vil si at random walk mønsteret har en bestemt retning. Dette kan vi skrive som
\[y_t = \alpha_0 + y_{t-1} + e_t\]
#lager en matrise som er 50x30
y<-matrix(,nrow=50, ncol=30)
#fyller det ved bruke av en for-loop
for(s in 1:30){
e <- rnorm(50) #her er vår stokastisk variabel - vi simulerer 50 trekk fra en normal-fordeling
y[,s]<-cumsum(2+e) #y er da modellert som summen av 2 og alle "sjokkene", e.
}
y = data.frame(y)
y["t"]<-1:50
y_l=melt(data.frame(y), id.var="t")
ggplot(y_l, aes(x=t, y=value, group=variable)) +
geom_line(alpha=.5)
Man kunne bruke random walk med drift til å modellere aksjepriser over tid, for eksempel. Vi ser at serien har en tydelig trend oppover, og at den har den vanlige egenskapen til en en random walk - variansen som øker over tid og en ikke-definert gjennomsnitt (forventning).
Differanser av highly persistent data
Highly persistent data (random walk, osv) er, som nevnt, vanskelig å jobbe med i regresjons-sammenheng. Heldigvis er det ofte en enkel løsning - differensiering (first difference).
Tidsrekken:
\[y_t = \alpha_0 + y_{t-1} + e_t\]
blir:
\[\Delta y_t = y_t- y_{t-1} = \alpha + e_t\]
Vi kan simulere:
for(s in 1:30){
e <- rnorm(50) #her er vår stokastisk variabel - vi simulerer 50 trekk fra en normal-fordeling
y[,s]<-cumsum(2+e) #y er da modellert som summen av 2 og alle "sjokkene", e.
}
y = data.frame(y)
y["t"]<-1:50
y <- zoo(y, order.by = y$t)
Dy <- diff(y)
Dy <- data.frame(Dy)
Dy["t"]<-1:49
Dy_l=melt(Dy, id.var="t")
ggplot(Dy_l, aes(x=t, y=value, group=variable)) +
geom_line(alpha=.2) +
geom_hline(yintercept=2)
Nå kan vi se at variansen av datasettet ser omtrent konstant gjennom tid. Vi kunne klart å beregne en OLS regresjon som ville gitt oss omtrent en vannrett linje.
Regresjon med first-difference
Vi bruker fruktbarhetsdataen igjen for å vise hvordan å bruke first-difference:
fertil3 <- read.dta("http://fmwww.bc.edu/ec-p/data/wooldridge/fertil3.dta")
fertil3_zoo <- zoo(fertil3, order.by=fertil3$year)
plot(fertil3_zoo$gfr)
plot(diff(fertil3_zoo$gfr))
Datarekken for fruktbarhet (gfr) ser ut som den er persistent, men vi kan ta en first difference.
Da kan vi kjøre regresjon av fruktbarhet på skattefradrag (pe).
fert_d_mod1 <- dynlm(diff(gfr) ~ diff(pe), data = fertil3_zoo)
fert_d_mod2 <- dynlm(diff(gfr) ~ diff(pe) + L(diff(pe)) + L(d(pe), 2) , data = fertil3_zoo)
stargazer(fert_d_mod1, fert_d_mod2, type="text")##
## ============================================================
## Dependent variable:
## ----------------------------------------
## diff(gfr)
## (1) (2)
## ------------------------------------------------------------
## diff(pe) -0.043 -0.036
## (0.028) (0.027)
##
## L(diff(pe)) -0.014
## (0.028)
##
## L(d(pe), 2) 0.110***
## (0.027)
##
## Constant -0.785 -0.964**
## (0.502) (0.468)
##
## ------------------------------------------------------------
## Observations 71 69
## R2 0.032 0.232
## Adjusted R2 0.018 0.197
## Residual Std. Error 4.221 (df = 69) 3.859 (df = 65)
## F Statistic 2.263 (df = 1; 69) 6.563*** (df = 3; 65)
## ============================================================
## Note: *p<0.1; **p<0.05; ***p<0.01Oppgave: Hvordan kan man tolke resultatene? Hva blir resultatene om man ikke bruker first-diff?
Svar:
Når vi gjør en first-difference, så spør vi – Hva er effekten av en endring i skattefradrag på endring i fruktbarhet? Eller – Hvis vi endrer skattefradrag, hva blir endringen i fruktbarhet? Koeffisienten på skattefradrag og laggen er negativ, men koeffisienten er ikke signifikant annerledes enn 0. Vi kan, med andre ord, ikke konkludere med noe.
Hvis vi ikke brukte first-difference, ville vi gjort som følgende:
fert_mod2 <- dynlm(gfr ~ pe + L(pe) + L(pe, 2) , data = fertil3_zoo)
summary(fert_mod2)##
## Time series regression with "zoo" data:
## Start = 1915, End = 1984
##
## Call:
## dynlm(formula = gfr ~ pe + L(pe) + L(pe, 2), data = fertil3_zoo)
##
## Residuals:
## Min 1Q Median 3Q Max
## -30.146 -16.062 -0.127 17.562 31.842
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 93.15791 4.49965 20.703 <2e-16 ***
## pe -0.01584 0.14026 -0.113 0.910
## L(pe) -0.02134 0.21523 -0.099 0.921
## L(pe, 2) 0.05390 0.13811 0.390 0.698
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.78 on 66 degrees of freedom
## Multiple R-squared: 0.006135, Adjusted R-squared: -0.03904
## F-statistic: 0.1358 on 3 and 66 DF, p-value: 0.9383I denne tilfellen er det ikke store kvalitative forskjeller mellom de to regresjonene – vi får negative koeffisienter, men disse er ikke signifikante annerledes enn 0. Dette kan tyde på at vi ikke behvøde å differensiere de opprinnelige seriene for å få valide inferens.
Seriekorrelasjon
Vi har sett på seriekorrelasjon/autokorrelasjon før - at feilleddet er korrelert over tid i en regresjon. Hvis alt annet er bra (konstant varians, osv) så påvirker ikke dette om estimatet er biased, men standardavvik er ikke lenger “efficient” - det vil si at vi ikke kan stole på standard avvik og p-verdier. Det er best å prøve å kontrollere for seriekorrelasjon.
Testing for seriekorrelasjon - Phillips curve
Phillipskurven er den observerte forholdet mellom inflasjon og arbeidsledighet. Vi tar inn data for å prøve å måle dette forholdet
phillips <- read.dta("http://fmwww.bc.edu/ec-p/data/wooldridge/phillips.dta")
phillips_zoo <- zoo(phillips, order.by=phillips$year)
#scatter plot
plot(phillips[c("unem", "inf")])
#line plots
plot(phillips_zoo$unem)
plot(phillips_zoo$inf)
Vi kjører en enkel regresjon og da får vi ut residualene (husker du hvordan vi definerer residualer?). Vi skal plotte residualene i mot “fitted” verdiene (det vil si, prediksjonene som modellen gir for inflasjon gitt arbeidslediget.) Hvis modellen er grei, burde vi ikke se noe mønster her. Hvis vi ser en tydelig mønster, så er det noe vi burde prøve å modellere. Dette inkluderer seriekorrelasjon – ser det ut som at en høy verdi er følgt av en annen høy verdi, osv?
phillips_mod <- dynlm(inf ~ unem, data=phillips_zoo)
phill_resid <- resid(phillips_mod)
phill_fitted <- fitted(phillips_mod)
#fitted vs. residual plott
plot(phill_fitted, phill_resid)
Fitted vs. residual plotten ser egentlig ganske bra ut (det burde se ganske tilfeldig - ingen bestemte trender.)
Men vi kan teste for seriekorrelasjon (autokorrelasjon) også på en mer formelt måte:
ar_test <- dynlm(phill_resid ~ L(phill_resid))
summary(ar_test)##
## Time series regression with "zoo" data:
## Start = 1949, End = 1996
##
## Call:
## dynlm(formula = phill_resid ~ L(phill_resid))
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.0764 -1.2577 -0.2327 0.9657 6.5375
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.1134 0.3594 -0.316 0.754
## L(phill_resid) 0.5730 0.1161 4.934 1.1e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.49 on 46 degrees of freedom
## Multiple R-squared: 0.346, Adjusted R-squared: 0.3318
## F-statistic: 24.34 on 1 and 46 DF, p-value: 1.098e-05Her kjører vi en regresjon av residualene på sin lagg. Siden laggen er signifikant, har vi grunn til å tro at residualene er seriekorrelert. Dette er en enkel form av det som er kalt Breusch-Godfrey testen
Vi kan gjøre det samme med en model som er “expectations augmented”, det vil si at forholdet er mellom endring i inflasjon fra år til år og arbeidsledighet.
Oppgave: Kjør en regresjon av en “expectations augmented” phillips-kurve, tolk resultatene. Kan du knytte det til vår diskusjon om phillips-kurven tidligere i semesteret? Test om det finnes seriekorrelasjon i feiltermin fra denne modellen.
Svar
phillips_mod2 <- dynlm(d(inf) ~ unem, data=phillips_zoo)
summary(phillips_mod2)##
## Time series regression with "zoo" data:
## Start = 1949, End = 1996
##
## Call:
## dynlm(formula = d(inf) ~ unem, data = phillips_zoo)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.1293 -0.7820 0.0611 0.9887 5.3600
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.0306 1.3768 2.201 0.0328 *
## unem -0.5426 0.2302 -2.357 0.0227 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.451 on 46 degrees of freedom
## Multiple R-squared: 0.1078, Adjusted R-squared: 0.0884
## F-statistic: 5.558 on 1 and 46 DF, p-value: 0.02271Her ser vi at endringer i inflasjon er signifikant og negativt relatert til arbeidsledighet, som teorien om den “expectations-augmented” Phillipskurven ville forutsi.
Vi kan videre teste for seriekorrelasjon:
phill_resid2 <- resid(phillips_mod2)
ar_test2 <- dynlm(phill_resid2 ~ L(phill_resid2))
summary(ar_test2)##
## Time series regression with "zoo" data:
## Start = 1950, End = 1996
##
## Call:
## dynlm(formula = phill_resid2 ~ L(phill_resid2))
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.4062 -0.8323 -0.0874 0.7605 5.2493
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.19417 0.30038 0.646 0.521
## L(phill_resid2) -0.03559 0.12389 -0.287 0.775
##
## Residual standard error: 2.059 on 45 degrees of freedom
## Multiple R-squared: 0.001831, Adjusted R-squared: -0.02035
## F-statistic: 0.08254 on 1 and 45 DF, p-value: 0.7752Ved bruk av pakken lmtest kan vi teste for autokorrelasjon (AR(1)) med Breusch-Godfrey-testen med bare en linje med kode.
library(lmtest)
bgtest(inf ~ unem, data=phillips_zoo)##
## Breusch-Godfrey test for serial correlation of order up to 1
##
## data: inf ~ unem
## LM test = 18.472, df = 1, p-value = 1.725e-05Vi kan også teste for seriekorrelasjon med lengre virkning
\[e_t = \alpha_1 e_{t-1} + \alpha_2 e_{t-2} + ... + \alpha_n e_{n} \]
Vi bruker barium datasettet til å se på dette:
barium <- read.dta("http://fmwww.bc.edu/ec-p/data/wooldridge/barium.dta")
barium_zoo <- zoo(barium, order.by = barium$t, frequency=12)
barium_mod <- dynlm(log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + befile6 + affile6 + afdec6, data=barium_zoo)
bgtest(barium_mod, order=3, type="F")##
## Breusch-Godfrey test for serial correlation of order up to 3
##
## data: barium_mod
## LM test = 5.1247, df1 = 3, df2 = 121, p-value = 0.002264En annen test som brukes ofte for seriekorrelasjon er durbin-watson. Vi viser det med phillips-kurv datasettet:
#funksjonen dwtest kommer med pakken lmtest
dwtest(phillips_mod) ##
## Durbin-Watson test
##
## data: phillips_mod
## DW = 0.8027, p-value = 7.552e-07
## alternative hypothesis: true autocorrelation is greater than 0dwtest(phillips_mod2)##
## Durbin-Watson test
##
## data: phillips_mod2
## DW = 1.7696, p-value = 0.1783
## alternative hypothesis: true autocorrelation is greater than 0Sist men ikke minst, kan vi se på autokorrelasjonsfunksjoner (ACF) og “partial” autokorrelasjonsfunksjoner (pacf)
Vi kan se på inflasjon serien i vår phillipskurv datasett og bruke en autokorrelasjonsfunksjon (som visuelt heter “Correlelogram”) (Se Hill, et al. k. 16).
acf(phillips_zoo$inf)
Du kan tenke på en ACF som om hvis det er en type “sjokk” i t=0 av 1, hvor lenge vedvarer denne sjokken. Kanskje vi får en liten sjokk i inflasjon – høyere bensinpriser på grunn av orkanen i USA – hvor lenge påvirker det inflasjonsserien? Vi ser at inflasjon har noe seriekorrelasjon, men at en sjokk dør ut etter 4-5 år.
Vi kan også se på en pACF:
pacf(phillips_zoo$inf)
Vi kan tolke en pACF mer direkt som hvor mange lagger ville være signifikant i en regresjon. Her ser det ut som det er hovedsakelig en lag som signifikant, og at vi kunne modellere inflasjon som en AR(1) modell.
Regresjon som er robust til seriekorrelasjon
Vi begynner med å laste inn data om Brent oljepriser fra FRED data-basen:
brent<-read.csv("https://fred.stlouisfed.org/graph/fredgraph.csv?chart_type=line&recession_bars=on&log_scales=&bgcolor=%23e1e9f0&graph_bgcolor=%23ffffff&fo=Open+Sans&ts=12&tts=12&txtcolor=%23444444&show_legend=yes&show_axis_titles=yes&drp=0&cosd=1987-05-01&coed=2016-08-01&height=450&stacking=&range=&mode=fred&id=MCOILBRENTEU&transformation=lin&nd=1987-05-01&ost=-99999&oet=99999&lsv=&lev=&mma=0&fml=a&fgst=lin&fgsnd=2009-06-01&fq=Monthly&fam=avg&vintage_date=&revision_date=&line_color=%234572a7&line_style=solid&lw=2&scale=left&mark_type=none&mw=2&width=1168")
head(brent)## DATE MCOILBRENTEU
## 1 1987-05-01 18.58
## 2 1987-06-01 18.86
## 3 1987-07-01 19.86
## 4 1987-08-01 18.98
## 5 1987-09-01 18.31
## 6 1987-10-01 18.76Litt rengjøring…
brent["DATE"] <- as.Date(brent$DATE)
brent_ts <- ts(brent$MCOILBRENTEU, start=c(1987,1), end=c(2016,8), frequency=12)
plot(brent_ts)
Oppgave
Hva kan man si om denne serien bare ved å se på den?
Bruk ACF og pACF til å analysere hvor autokorrelert denne serien er. Hva sier dette intuitivt?
Hvis det er en høy grad av seriekorrelasjon, kan serien være “ikke-stasjonært”, som gjør det vanseklig å analysere med OLS. En løsning er å ta first-difference, eller å se på avkastning. Gjør dette med oljeprisdata og kjør en ACF igjen. Hva har skjedd?
Hvis det er fortsatt autokorrelasjon i serien, kjør en dynamisk regresjon som modellerer denne seriekorrelasjonen eksplisitt.
Svar
En ting vi kan se er variansen virker som at det øker over tid – større svingninger i oljeprisen etter hvert. Vi kan også se at serien ikke går tilbake til en viss gjennomsnittsverdi. Vi må dermed være forsiktig med å kjøre OLS.
acf(brent_ts)
Om det ikke var significant autokorrelasjon, ville de vertikale linjene vært nær 0. Isteden ser vi at det er mye autokorrelasjon i datasettet. Fra det vi så tidligere, så kan oljeprisen minne på en random walk - den er highly persistent.
La oss transformere datasette til månedlig avkastning og sjekke acf igjen.
brent_ret <- diff(log(brent_ts))
plot(brent_ret)
acf(brent_ret)
Nå ser det bedre ut, men det ser ut som avkastning har fortsatt autokorrelasjon. Vi kan modellere dette med dynlm pakken
brent_mod <- dynlm(brent_ret~ L(brent_ret) + L(brent_ret, 2))
summary(brent_mod)##
## Time series regression with "ts" data:
## Start = 1987(4), End = 2016(8)
##
## Call:
## dynlm(formula = brent_ret ~ L(brent_ret) + L(brent_ret, 2))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.90760 -0.05116 0.00785 0.05934 0.43185
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0001582 0.0053258 -0.030 0.976313
## L(brent_ret) 0.2065738 0.0534559 3.864 0.000133 ***
## L(brent_ret, 2) -0.0090583 0.0534737 -0.169 0.865581
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1001 on 350 degrees of freedom
## Multiple R-squared: 0.04199, Adjusted R-squared: 0.03652
## F-statistic: 7.671 on 2 and 350 DF, p-value: 0.000549Den første laggen er positiv og signifikant. Vi kan dermed tolke det til å bety at en endring i oljeprisen en periode vil ofte føre til, i gjennomsnitt, en økning i neste periode også.