(a)
library(Quandl)
## Warning: package 'Quandl' was built under R version 3.4.4
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(tseries)
## Warning: package 'tseries' was built under R version 3.4.4
Poil <- Quandl("FRED/MCOILWTICO", type ="ts", collapse = "monthly", start_date = "1995-01-01", end_date = "2017-12-01")
Pgas <- Quandl("FRED/GASREGCOVM", type ="ts", collapse = "monthly", start_date = "1995-01-01", end_date = "2017-12-01")
LPoil <- log(Poil)
LPgas <- log(Pgas)
# Create a single time seris plot
seqplot.ts(LPoil, LPgas, xlab="year", ylab="Price", colx = "black", coly = "red", main = "The log of Oil and Gas Price" )
legend("topleft", legend=c("LPoil","LPgas"), col=1:2, lty=1:1)
# The log of Oil and gas price move together.
(b)
# Unit Root Test - Series in Level
library(urca)
# KPSS test
KPSS.testLPoil <- ur.kpss(LPoil, type = "tau", lags = "short")
KPSS.testLPgas <- ur.kpss(LPgas, type = "tau", lags = "short")
summary(KPSS.testLPoil)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.6529
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
summary(KPSS.testLPgas)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.5874
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
# The KPSS test has a test statistic of 0.6529 for the log of oil price, meaning that LPoil has a unit root.
# The KPSS test has a test statistic of 0.0.5874 for the log of gas price, meaning that LPgas has a unit root.
# ERS test
ERS.testLPoil <- ur.ers(LPoil, type = "P-test", model = "trend")
ERS.testLPgas <- ur.ers(LPgas, type = "P-test", model = "trend")
summary(ERS.testLPoil)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 9.4782
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary(ERS.testLPgas )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 8.2276
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
# The ERS test has a test statistic of 9.4782 for the log of oil price, meaning that LPoil has a unit root.
# The ERS test has a test statistic of 8.2276 for the log of gas price, meaning that LPgas has a unit root.
# Unit Root Test - Series in First Difference
dLPoil <- diff(LPoil)
dLPgas <- diff(LPgas)
KPSS.testdLPoil <- ur.kpss(dLPoil, type = "tau", lags = "short")
KPSS.testdLPgas <- ur.kpss(dLPgas, type = "tau", lags = "short")
summary(KPSS.testdLPoil)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.0444
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
summary(KPSS.testdLPgas)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.0377
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
# The KPSS test has a test statistic of 0.0444 for the log of oil price, meaning that LPoil has a I(1) process.
# The KPSS test has a test statistic of 0.0377 for the log of gas price, meaning that LPgas has a I(1) process.
ERS.testdLPoil <- ur.ers(dLPoil, type = "P-test", lag.max=12, model = "trend")
ERS.testdLPgas <- ur.ers(dLPgas, type = "P-test", lag.max=12, model = "trend")
summary(ERS.testdLPoil)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 0.7873
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary(ERS.testdLPgas)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 0.4781
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
# The ERS test has a test statistic of 0.7873 for the log of oil price, meaning that LPoil has a I(1) process.
# The ERS test has a test statistic of 0.4781 for the log of gas price, meaning that LPgas has a I(1) process.
# As a result, the first difference of the LPoil and LPgas are approximately weakly stationary.
(c)
library(vars)
## Loading required package: MASS
## Warning: package 'MASS' was built under R version 3.4.4
## Loading required package: strucchange
## Loading required package: sandwich
## Loading required package: lmtest
# Cointegration Analysis
v1 <- cbind(LPoil, LPgas)
# Determine the Number of Lags to include in Cointegration Analysis
var <- VARselect(v1, type = "trend")
nlags <- var$selection["SC(n)"]
# nlags is the number of lags to include in cointegration analysis using SIC
# case 2, restricted constant: ecdet="const"
# case 3, unrestricted constant: ecdet="none"
# case 4, restricted trend: ecdet="trend"
# Using time series plots from part (a) as a guide to determine the specification of the
# deterministic components in the cointegration test, we consider Case 4.
# trace test
v1.CA <- ca.jo(v1, ecdet = "trend", type = "trace", K=nlags, spec = "transitory", season = 12)
summary(v1.CA)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , with linear trend in cointegration
##
## Eigenvalues (lambda):
## [1] 1.037680e-01 1.873353e-02 -1.699366e-19
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 5.16 10.49 12.25 16.26
## r = 0 | 35.07 22.76 25.32 30.45
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## LPoil.l1 LPgas.l1 trend.l1
## LPoil.l1 1.0000000000 1.00000000 1.0000000
## LPgas.l1 -1.6267117480 3.75026554 -10.7065612
## trend.l1 0.0009487176 -0.01919425 0.1692144
##
## Weights W:
## (This is the loading matrix)
##
## LPoil.l1 LPgas.l1 trend.l1
## LPoil.d -0.006283109 -0.008367640 -2.545462e-19
## LPgas.d 0.138252677 -0.003301629 8.623684e-20
# The results show that null hypothesis of r=0 is rejected.
# The test statistic is 35.07 and the critical value at the 1% level is 30.45.
# We fail to reject the test for r<=1 since the test statistic is 5.16
# and the critical value at the 10% level is 10.49.
# The cointegration trace test suggests that r should be equal to 1.
# max eigenvalue test
v1.CA <- ca.jo(v1, ecdet = "trend", type = "eigen", K=nlags, spec = "transitory", season = 12)
summary(v1.CA)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: maximal eigenvalue statistic (lambda max) , with linear trend in cointegration
##
## Eigenvalues (lambda):
## [1] 1.037680e-01 1.873353e-02 -1.699366e-19
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 5.16 10.49 12.25 16.26
## r = 0 | 29.91 16.85 18.96 23.65
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## LPoil.l1 LPgas.l1 trend.l1
## LPoil.l1 1.0000000000 1.00000000 1.0000000
## LPgas.l1 -1.6267117480 3.75026554 -10.7065612
## trend.l1 0.0009487176 -0.01919425 0.1692144
##
## Weights W:
## (This is the loading matrix)
##
## LPoil.l1 LPgas.l1 trend.l1
## LPoil.d -0.006283109 -0.008367640 -2.545462e-19
## LPgas.d 0.138252677 -0.003301629 8.623684e-20
# The results show that the null hypothesis of r=0 is rejected.
# The test statistic is 29.91 and the critical value at the 1% level is 23.65.
# However, we fail to reject the test for r<=1 since the test statistic is 5.16
# and the critical value at the 10% level is 10.49.
# The cointegration eigen test suggests that r should be equal to 1.
(d)
# LR-test for no linear trend
lttest(v1.CA, r = 1)
## LR-test for no linear trend
##
## H0: H*2(r<=1)
## H1: H2(r<=1)
##
## Test statistic is distributed as chi-square
## with 1 degress of freedom
## test statistic p-value
## LR test 0.63 0.43
# the results of the test justify the specification used on part (c).
(e)
# estimate a bivariate VEC model for LPoil and LPgas
v1.CA.VEC <- cajorls(v1.CA , r = 1)
v1.CA.VEC
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## LPoil.d LPgas.d
## ect1 -0.006283 0.138253
## constant 0.021126 -0.386041
## sd1 0.040323 0.037469
## sd2 0.042653 0.037282
## sd3 0.082663 0.080786
## sd4 0.057719 0.065159
## sd5 0.046381 0.068482
## sd6 0.036897 0.051379
## sd7 0.040373 0.034379
## sd8 0.029176 0.044451
## sd9 0.036513 0.042695
## sd10 0.018676 0.002836
## sd11 -0.005861 0.002439
## LPoil.dl1 0.344830 0.229270
## LPgas.dl1 -0.308645 0.110101
## LPoil.dl2 0.164034 0.068542
## LPgas.dl2 -0.057460 -0.209944
##
##
## $beta
## ect1
## LPoil.l1 1.0000000000
## LPgas.l1 -1.6267117480
## trend.l1 0.0009487176
(f)
summary(v1.CA.VEC$rlm)
## Response LPoil.d :
##
## Call:
## lm(formula = LPoil.d ~ ect1 + constant + sd1 + sd2 + sd3 + sd4 +
## sd5 + sd6 + sd7 + sd8 + sd9 + sd10 + sd11 + LPoil.dl1 + LPgas.dl1 +
## LPoil.dl2 + LPgas.dl2 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.271652 -0.047691 0.004689 0.057003 0.196982
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 -0.006283 0.062206 -0.101 0.919625
## constant 0.021126 0.174776 0.121 0.903886
## sd1 0.040323 0.024569 1.641 0.101980
## sd2 0.042653 0.025210 1.692 0.091883 .
## sd3 0.082663 0.025428 3.251 0.001305 **
## sd4 0.057719 0.026145 2.208 0.028154 *
## sd5 0.046381 0.026298 1.764 0.078976 .
## sd6 0.036897 0.026260 1.405 0.161219
## sd7 0.040373 0.025707 1.571 0.117533
## sd8 0.029176 0.024979 1.168 0.243881
## sd9 0.036513 0.024803 1.472 0.142222
## sd10 0.018676 0.024992 0.747 0.455592
## sd11 -0.005861 0.024562 -0.239 0.811593
## LPoil.dl1 0.344830 0.091963 3.750 0.000219 ***
## LPgas.dl1 -0.308645 0.143061 -2.157 0.031903 *
## LPoil.dl2 0.164034 0.093158 1.761 0.079466 .
## LPgas.dl2 -0.057460 0.125308 -0.459 0.646947
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.08161 on 256 degrees of freedom
## Multiple R-squared: 0.1367, Adjusted R-squared: 0.07932
## F-statistic: 2.384 on 17 and 256 DF, p-value: 0.001969
##
##
## Response LPgas.d :
##
## Call:
## lm(formula = LPgas.d ~ ect1 + constant + sd1 + sd2 + sd3 + sd4 +
## sd5 + sd6 + sd7 + sd8 + sd9 + sd10 + sd11 + LPoil.dl1 + LPgas.dl1 +
## LPoil.dl2 + LPgas.dl2 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.212348 -0.022546 0.001538 0.029239 0.122802
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.138253 0.035314 3.915 0.000116 ***
## constant -0.386041 0.099220 -3.891 0.000127 ***
## sd1 0.037469 0.013948 2.686 0.007696 **
## sd2 0.037282 0.014312 2.605 0.009727 **
## sd3 0.080786 0.014435 5.596 5.62e-08 ***
## sd4 0.065159 0.014842 4.390 1.66e-05 ***
## sd5 0.068482 0.014929 4.587 7.04e-06 ***
## sd6 0.051379 0.014908 3.446 0.000664 ***
## sd7 0.034379 0.014594 2.356 0.019244 *
## sd8 0.044451 0.014180 3.135 0.001921 **
## sd9 0.042695 0.014081 3.032 0.002678 **
## sd10 0.002836 0.014188 0.200 0.841728
## sd11 0.002439 0.013944 0.175 0.861300
## LPoil.dl1 0.229270 0.052207 4.392 1.65e-05 ***
## LPgas.dl1 0.110101 0.081216 1.356 0.176404
## LPoil.dl2 0.068542 0.052886 1.296 0.196131
## LPgas.dl2 -0.209944 0.071137 -2.951 0.003458 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04633 on 256 degrees of freedom
## Multiple R-squared: 0.4907, Adjusted R-squared: 0.4568
## F-statistic: 14.51 on 17 and 256 DF, p-value: < 2.2e-16
# The results show that a1 is not statistically significant and is negative,
# but a2 is statistically significant at the 0.01% level and the sign is positive.
# It means Whenever there is a disruption in the long run equilibirum,
# the signs of the adjustment parameters are consistent with the error correction mechanism.
(g)
# We now re-estimate the VEC model with the restriction that a2 =0
# test for restricted adjustment parameters alpha
restricted.alpha <- matrix(c(1,0), c(2,1))
v1.CA.restricted.alpha <- alrtest(v1.CA, A = restricted.alpha, r = 1)
v1.VEC.restricted <- cajorls(v1.CA.restricted.alpha, r = 1)
v1.VEC.restricted
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## LPoil.d LPgas.d
## ect1 -0.0557042 0.1004001
## constant 0.1611240 -0.2818895
## sd1 0.0420914 0.0393371
## sd2 0.0444985 0.0395397
## sd3 0.0840561 0.0831298
## sd4 0.0576318 0.0658194
## sd5 0.0447996 0.0677313
## sd6 0.0340551 0.0489595
## sd7 0.0367786 0.0309679
## sd8 0.0259067 0.0410684
## sd9 0.0333962 0.0393111
## sd10 0.0161125 0.0001366
## sd11 -0.0071299 0.0012353
## LPoil.dl1 0.3779774 0.2583171
## LPgas.dl1 -0.3430617 0.0741168
## LPoil.dl2 0.1890390 0.0886464
## LPgas.dl2 -0.0839394 -0.2457155
##
##
## $beta
## ect1
## LPoil.l1 1.0000000000
## LPgas.l1 -1.4901472670
## trend.l1 0.0004371265
summary(v1.VEC.restricted)
## Length Class Mode
## rlm 12 mlm list
## beta 3 -none- numeric
(h)
# The intuition behind the restriction is the idea that innovations in oil price affect
# gas price, but innovations in gas price do not affect oil price.
# Therefore, we should restrict the model to imply that gas price cannot affect oil price.
(i)
v1.VAR <- vec2var(v1.CA, r = 1)
v1.VAR.fcst <- predict(v1.VAR, n.ahead = 180)
library(ggfortify)
## Warning: package 'ggfortify' was built under R version 3.4.4
## Loading required package: ggplot2
## Warning: package 'ggplot2' was built under R version 3.4.4
autoplot(v1.VAR.fcst)
