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
Quandl.api_key('4-KG5x_Vo7rXzmZNAHch')
library("tseries")
library("urca")
library("vars")
## Loading required package: MASS
## Loading required package: strucchange
## Loading required package: sandwich
## Loading required package: lmtest
library("gdata")
## gdata: read.xls support for 'XLS' (Excel 97-2004) files ENABLED.
##
## gdata: read.xls support for 'XLSX' (Excel 2007+) files ENABLED.
##
## 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
dpi<-Quandl("FRED/DPI",type="zoo")
pces <- Quandl("FRED/PCES", type="zoo")
str(dpi)
## 'zooreg' series from 1947 Q1 to 2016 Q1
## Data: num [1:277] 171 170 178 180 185 ...
## Index: Class 'yearqtr' num [1:277] 1947 1947 1948 1948 1948 ...
## Frequency: 4
str(pces)
## 'zooreg' series from Jan 1959 to Mar 2016
## Data: num [1:687] 139 140 141 142 143 ...
## Index: Class 'yearmon' num [1:687] 1959 1959 1959 1959 1959 ...
## Frequency: 12
x1=log(dpi)
x2=log(pces)
v1 <- cbind(x1,x2)
## Warning in merge.zoo(..., all = all, fill = fill, suffixes = suffixes,
## retclass = "zoo", : Index vectors are of different classes: yearqtr yearmon
plot(v1,xlab="year", main="1947-Q1-2016-Q1")
dx1<-diff(x1)
dx2<-diff(x2)
par(mfrow=c(1,2))
plot(dx1,xlab = "year", main="1947-2016")
plot(dx2,xlab = "year", main="1947-2016")
summary( ur.ers(x1, 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: 83.7733
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary(ur.ers(x2, 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: 362.5563
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary (ur.ers(dx1, 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: 2.2919
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
summary( ur.ers(dx2, 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: 3.4123
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
A.ca <- ca.jo(v1, ecdet="trend", type="trace", K=2, spec="transitory")
summary(A.ca)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , with linear trend in cointegration
##
## Eigenvalues (lambda):
## [1] 2.276300e-01 1.592440e-02 9.782889e-17
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 3.64 10.49 12.25 16.26
## r = 0 | 62.28 22.76 25.32 30.45
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## x1.l1 x2.l1 trend.l1
## x1.l1 1.000000000 1.000000000 1.000000000
## x2.l1 -0.827565825 -0.860233917 -1.222375450
## trend.l1 -0.001526432 -0.000175721 0.006445846
##
## Weights W:
## (This is the loading matrix)
##
## x1.l1 x2.l1 trend.l1
## x1.d 0.03960338 -0.033694408 4.935164e-14
## x2.d 0.11497726 -0.003844646 6.509465e-13
A.ca <- ca.jo(v1, ecdet="trend", type="eigen", K=2, spec="transitory")
lttest(A.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 2.07 0.15
A.ca.vec <- cajorls(A.ca, r=1)
A.ca.vec
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## x1.d x2.d
## ect1 0.03960 0.11498
## constant -0.05948 -0.18114
## x1.dl1 -0.11663 0.05330
## x2.dl1 0.59846 0.28171
##
##
## $beta
## ect1
## x1.l1 1.000000000
## x2.l1 -0.827565825
## trend.l1 -0.001526432
summary(A.ca.vec$rlm)
## Response x1.d :
##
## Call:
## lm(formula = x1.d ~ ect1 + constant + x1.dl1 + x2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.046299 -0.004782 -0.000343 0.004465 0.033381
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.03960 0.02627 1.508 0.133
## constant -0.05948 0.04259 -1.397 0.164
## x1.dl1 -0.11663 0.07364 -1.584 0.115
## x2.dl1 0.59846 0.11354 5.271 3.2e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.008679 on 223 degrees of freedom
## Multiple R-squared: 0.7918, Adjusted R-squared: 0.788
## F-statistic: 212 on 4 and 223 DF, p-value: < 2.2e-16
##
##
## Response x2.d :
##
## Call:
## lm(formula = x2.d ~ ect1 + constant + x1.dl1 + x2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0130967 -0.0026803 0.0004772 0.0025095 0.0156722
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.11498 0.01449 7.938 9.99e-14 ***
## constant -0.18114 0.02348 -7.713 4.05e-13 ***
## x1.dl1 0.05330 0.04061 1.313 0.191
## x2.dl1 0.28171 0.06261 4.499 1.10e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.004786 on 223 degrees of freedom
## Multiple R-squared: 0.9413, Adjusted R-squared: 0.9402
## F-statistic: 893.2 on 4 and 223 DF, p-value: < 2.2e-16
rest.bettal <- matrix(c(0,-1,0,0,0,1), c(3,2))
a.rbettal <- blrtest(A.ca, H=rest.bettal, r=1)
summary(a.rbettal)
##
## ######################
## # 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.0929 0.0027
##
## The value of the likelihood ratio test statistic:
## 36.51 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]
## x1.d NaN NaN
## x2.d NaN NaN
r.alpha <- matrix(c(1,0), c(2,1))
a.r.alpha<- alrtest(A.ca, A=r.alpha, r=1)
summary(a.r.alpha)
##
## ######################
## # 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.0245 0.0000 0.0000
##
## The value of the likelihood ratio test statistic:
## 52.99 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.x1.l1 1.0000
## RK.x2.l1 -0.8420
## RK.trend.l1 -0.0009
##
## Weights W of the restricted VAR:
##
## [,1]
## [1,] -0.0705
## [2,] 0.0000
a.rboth1 <- ablrtest(A.ca, A=r.alpha, H=rest.bettal, r=1)
summary(a.rboth1)
##
## ######################
## # 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.0044 0.0000
##
## The value of the likelihood ratio test statistic:
## 57.64 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
a.var <- vec2var(A.ca, r=1)
a.var.fcst <- predict(a.var, n.ahead = 8)
par( mar=c(4,4,2,1), cex=0.75)
plot(a.var.fcst)
