library(Quandl)
## 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)
library(urca)
## Warning: package 'urca' was built under R version 3.2.4
library(vars)
## Warning: package 'vars' was built under R version 3.2.4
## Loading required package: MASS
## Loading required package: strucchange
## Warning: package 'strucchange' was built under R version 3.2.4
## Loading required package: sandwich
## Warning: package 'sandwich' was built under R version 3.2.4
## Loading required package: lmtest
## Warning: package 'lmtest' was built under R version 3.2.4
library(gdata)
## gdata: Unable to locate valid perl interpreter
## gdata:
## gdata: read.xls() will be unable to read Excel XLS and XLSX files
## gdata: unless the 'perl=' argument is used to specify the location
## gdata: of a valid perl intrpreter.
## gdata:
## gdata: (To avoid display of this message in the future, please
## gdata: ensure perl is installed and available on the executable
## gdata: search path.)
## gdata: Unable to load perl libaries needed by read.xls()
## gdata: to support 'XLX' (Excel 97-2004) files.
##
## gdata: Unable to load perl libaries needed by read.xls()
## gdata: to support 'XLSX' (Excel 2007+) files.
##
## gdata: Run the function 'installXLSXsupport()'
## gdata: to automatically download and install the perl
## gdata: libaries needed to support Excel XLS and XLSX formats.
##
## Attaching package: 'gdata'
## The following objects are masked from 'package:xts':
##
## first, last
## The following object is masked from 'package:stats':
##
## nobs
## The following object is masked from 'package:utils':
##
## object.size
fund<-Quandl("FRED/FEDFUNDS",type="zoo")
aaa<-Quandl("FED/RIMLPAAAR_N_M",type="zoo")
y<-cbind(fund,aaa)
y1<-window(y,start="Jul 1954", end = "Jan 2016")
plot(y1,main="Year(1954-2016)")
dfund<-diff(fund)
daaa<-diff(aaa)
yD<-cbind(dfund,daaa)
yD1<-cbind(yD,start="Jul 1954", end = "Jan 2016")
## Warning in merge.zoo(..., all = all, fill = fill, suffixes = suffixes,
## retclass = "zoo", : Index vectors are of different classes: yearmon integer
## integer
Unit root test (On level dataset)
summary( ur.ers(fund, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 14.3817
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary( ur.ers(aaa, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 27.783
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
Unit roor test (On differenced dataset)
summary( ur.ers(dfund, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # 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.8562
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary( ur.ers(daaa, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # 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.1118
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
From the unit root tests, it can bee seen that both Fed Fund and AAA are only stationary in 1st difference (with 1% level of significance).
Johansen’s test for Cointegration
fund_aaa<- ca.jo(y1, ecdet="const", type="trace", K=2, spec="transitory")
summary(fund_aaa)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , without linear trend and constant in cointegration
##
## Eigenvalues (lambda):
## [1] 5.656992e-02 6.236329e-03 1.476527e-18
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 4.61 7.52 9.24 12.97
## r = 0 | 47.53 17.85 19.96 24.60
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## fund.l1 aaa.l1 constant
## fund.l1 1.000000 1.00000 1.00000
## aaa.l1 -1.156911 66.64856 -1.39404
## constant 3.118329 -471.27749 -54.94366
##
## Weights W:
## (This is the loading matrix)
##
## fund.l1 aaa.l1 constant
## fund.d -0.03650126 -1.586843e-04 7.825856e-19
## aaa.d 0.01595198 -7.063239e-05 -2.116698e-19
fund_aaa <- ca.jo(y1, ecdet="const", type="eigen", K=2, spec="transitory")
summary(fund_aaa)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: maximal eigenvalue statistic (lambda max) , without linear trend and constant in cointegration
##
## Eigenvalues (lambda):
## [1] 5.656992e-02 6.236329e-03 1.476527e-18
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 4.61 7.52 9.24 12.97
## r = 0 | 42.92 13.75 15.67 20.20
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## fund.l1 aaa.l1 constant
## fund.l1 1.000000 1.00000 1.00000
## aaa.l1 -1.156911 66.64856 -1.39404
## constant 3.118329 -471.27749 -54.94366
##
## Weights W:
## (This is the loading matrix)
##
## fund.l1 aaa.l1 constant
## fund.d -0.03650126 -1.586843e-04 7.825856e-19
## aaa.d 0.01595198 -7.063239e-05 -2.116698e-19
Constant term test:
lttest(fund_aaa, 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 0.98
Restriction when betta2= -1 :
rest.betta2 <- matrix(c(0,-1,0,0,0,1), c(3,2))
rbetta2 <- blrtest(fund_aaa, H=rest.betta2, r=1)
summary(rbetta2)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Estimation and testing under linear restrictions on beta
##
## The VECM has been estimated subject to:
## beta=H*phi and/or alpha=A*psi
##
## [,1] [,2]
## [1,] 0 0
## [2,] -1 0
## [3,] 0 1
##
## Eigenvalues of restricted VAR (lambda):
## [1] 0.0062 0.0001
##
## The value of the likelihood ratio test statistic:
## 38.3 distributed as chi square with 1 df.
## The p-value of the test statistic is: 0
##
## Eigenvectors, normalised to first column
## of the restricted VAR:
##
## [,1] [,2]
## [1,] NaN NaN
## [2,] -Inf Inf
## [3,] Inf Inf
##
## Weights W of the restricted VAR:
##
## [,1] [,2]
## fund.d NaN NaN
## aaa.d NaN NaN
plotres(fund_aaa)
dpi<-Quandl("FRED/DPI",type="zoo")
pce<-Quandl("FRED/PCEC",type="zoo")
y1<-log(dpi)
y2<-log(pce)
z1<-cbind(y1,y2)
plot(z1, main="Quarterly(1947-2016)")
dy1<-diff(y1)
dy2<-diff(y2)
par(mfrow=c(1,2))
plot(dy1,xlab = "year(Quarterly)",main="Quarterly(1947-2016)")
plot(dy2,xlab = "year(Quarterly)",main="Quarterly(1947-2016)")
Unit root test on level:
summary( ur.ers(y1, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 83.7733
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary( ur.ers(y2, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 63.1974
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
Unit root test on differenced data:
summary( ur.ers(dy1, type="P-test", lag.max=12, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 2.2919
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary( ur.ers(dy2, type="P-test", lag.max=8, model="trend") )
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 1.7594
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
From the unit root tests, it can bee seen that both datasets are only stationary in 1st difference (with 1% level of significance).
Johansen’s test for Cointegration:
coi <- ca.jo(z1, ecdet="trend", type="trace", K=2, spec="transitory")
summary(coi)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , with linear trend in cointegration
##
## Eigenvalues (lambda):
## [1] 1.169941e-01 3.065947e-02 9.504492e-19
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 8.56 10.49 12.25 16.26
## r = 0 | 42.78 22.76 25.32 30.45
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## y1.l1 y2.l1 trend.l1
## y1.l1 1.0000000000 1.000000000 1.00000000
## y2.l1 -0.9467561512 -1.120563179 -1.80151062
## trend.l1 -0.0008356929 0.003153558 0.01371773
##
## Weights W:
## (This is the loading matrix)
##
## y1.l1 y2.l1 trend.l1
## y1.d 0.01960847 -0.024034218 3.575926e-15
## y2.d 0.12974403 -0.007937361 7.403759e-15
coi <- ca.jo(z1, ecdet="trend", type="eigen", K=2, spec="transitory")
summary(coi)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: maximal eigenvalue statistic (lambda max) , with linear trend in cointegration
##
## Eigenvalues (lambda):
## [1] 1.169941e-01 3.065947e-02 9.504492e-19
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 8.56 10.49 12.25 16.26
## r = 0 | 34.22 16.85 18.96 23.65
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## y1.l1 y2.l1 trend.l1
## y1.l1 1.0000000000 1.000000000 1.00000000
## y2.l1 -0.9467561512 -1.120563179 -1.80151062
## trend.l1 -0.0008356929 0.003153558 0.01371773
##
## Weights W:
## (This is the loading matrix)
##
## y1.l1 y2.l1 trend.l1
## y1.d 0.01960847 -0.024034218 3.575926e-15
## y2.d 0.12974403 -0.007937361 7.403759e-15
Constant term test:
lttest(coi, 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 9.9 0
Vector Error Correction (VEC):
coi.VEC <- cajorls(coi, r=1)
coi.VEC
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## y1.d y2.d
## ect1 0.019608 0.129744
## constant 0.001245 -0.039420
## y1.dl1 -0.029967 0.175551
## y2.dl1 0.461941 0.042985
##
##
## $beta
## ect1
## y1.l1 1.0000000000
## y2.l1 -0.9467561512
## trend.l1 -0.0008356929
t-stat and p-value
summary(coi.VEC$rlm)
## Response y1.d :
##
## Call:
## lm(formula = y1.d ~ ect1 + constant + y1.dl1 + y2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.051913 -0.004746 -0.000611 0.004612 0.058281
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.019608 0.027290 0.719 0.473
## constant 0.001245 0.010533 0.118 0.906
## y1.dl1 -0.029967 0.067041 -0.447 0.655
## y2.dl1 0.461941 0.071736 6.439 5.44e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01011 on 271 degrees of freedom
## Multiple R-squared: 0.7314, Adjusted R-squared: 0.7274
## F-statistic: 184.5 on 4 and 271 DF, p-value: < 2.2e-16
##
##
## Response y2.d :
##
## Call:
## lm(formula = y2.d ~ ect1 + constant + y1.dl1 + y2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.036571 -0.004066 -0.000287 0.004549 0.059174
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.129744 0.023428 5.538 7.23e-08 ***
## constant -0.039420 0.009042 -4.360 1.85e-05 ***
## y1.dl1 0.175551 0.057555 3.050 0.00251 **
## y2.dl1 0.042985 0.061585 0.698 0.48579
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.008681 on 271 degrees of freedom
## Multiple R-squared: 0.7853, Adjusted R-squared: 0.7821
## F-statistic: 247.8 on 4 and 271 DF, p-value: < 2.2e-16
Restriction when betta = -1 :
betta <- matrix(c(0,-1,0,0,0,1), c(3,2))
rbetta <- blrtest(coi, H=betta, r=1)
summary(rbetta)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Estimation and testing under linear restrictions on beta
##
## The VECM has been estimated subject to:
## beta=H*phi and/or alpha=A*psi
##
## [,1] [,2]
## [1,] 0 0
## [2,] -1 0
## [3,] 0 1
##
## Eigenvalues of restricted VAR (lambda):
## [1] 0.0385 0.0024
##
## The value of the likelihood ratio test statistic:
## 23.41 distributed as chi square with 1 df.
## The p-value of the test statistic is: 0
##
## Eigenvectors, normalised to first column
## of the restricted VAR:
##
## [,1] [,2]
## [1,] NaN NaN
## [2,] -Inf -Inf
## [3,] Inf Inf
##
## Weights W of the restricted VAR:
##
## [,1] [,2]
## y1.d NaN NaN
## y2.d NaN NaN
Test for restricted alpha:
rest.alpha <- matrix(c(1,0), c(2,1))
ralpha <- alrtest(coi, A=rest.alpha, r=1)
summary(ralpha)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Estimation and testing under linear restrictions on beta
##
## The VECM has been estimated subject to:
## beta=H*phi and/or alpha=A*psi
##
## [,1]
## [1,] 1
## [2,] 0
##
## Eigenvalues of restricted VAR (lambda):
## [1] 0.0429 0.0000 0.0000
##
## The value of the likelihood ratio test statistic:
## 22.15 distributed as chi square with 1 df.
## The p-value of the test statistic is: 0
##
## Eigenvectors, normalised to first column
## of the restricted VAR:
##
## [,1]
## RK.y1.l1 1.0000
## RK.y2.l1 -0.9937
## RK.trend.l1 0.0002
##
## Weights W of the restricted VAR:
##
## [,1]
## [1,] -0.0726
## [2,] 0.0000
Joint test for restricted alpha and betta:
rab <- ablrtest(coi, A=rest.alpha, H=betta, r=1)
summary(rab)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Estimation and testing under linear restrictions on alpha and beta
##
## The VECM has been estimated subject to:
## beta=H*phi and/or alpha=A*psi
##
## [,1] [,2]
## [1,] 0 0
## [2,] -1 0
## [3,] 0 1
##
##
## [,1]
## [1,] 1
## [2,] 0
##
## Eigenvalues of restricted VAR (lambda):
## [1] 0.0153 0.0000
##
## The value of the likelihood ratio test statistic:
## 29.97 distributed as chi square with 1 df.
## The p-value of the test statistic is: 0
##
## Eigenvectors, normalised to first column
## of the restricted VAR:
##
## [,1]
## [1,] NaN
## [2,] Inf
## [3,] -Inf
##
## Weights W of the restricted VAR:
##
## [,1]
## [1,] 0
## [2,] 0
var <- vec2var(coi, r=1)
Forecast:
var.fcst <- predict(var, n.ahead=8)
par( mar=c(4,4,2,1), cex=0.75)
plot(var.fcst)