load data ..
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)
c<-Quandl("FRED/CPILFENS",type="zoo")
str(c)
## 'zooreg' series from Jan 1957 to Nov 2015
## Data: num [1:707] 28.5 28.5 28.7 28.8 28.8 28.9 28.9 29 29.1 29.2 ...
## Index: Class 'yearmon' num [1:707] 1957 1957 1957 1957 1957 ...
## Frequency: 12
construct inflation
A<-diff(c)
B<-(A/c)*100
Inflation<-ts(B,start=c(1957,2),frequency=12)
plot(Inflation,xlab="",ylab="INF",main="INF")
DInflation=diff(Inflation)
plot(DInflation,xlab="",ylab=" INF",main=" InF")
adf.test(Inflation)
##
## Augmented Dickey-Fuller Test
##
## data: Inflation
## Dickey-Fuller = -3.8702, Lag order = 8, p-value = 0.01549
## alternative hypothesis: stationary
adf.test(DInflation)
## Warning in adf.test(DInflation): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: DInflation
## Dickey-Fuller = -12.627, Lag order = 8, p-value = 0.01
## alternative hypothesis: stationary
kpss.test(Inflation)
## Warning in kpss.test(Inflation): p-value smaller than printed p-value
##
## KPSS Test for Level Stationarity
##
## data: Inflation
## KPSS Level = 2.0216, Truncation lag parameter = 6, p-value = 0.01
kpss.test(DInflation)
## Warning in kpss.test(DInflation): p-value greater than printed p-value
##
## KPSS Test for Level Stationarity
##
## data: DInflation
## KPSS Level = 0.020071, Truncation lag parameter = 6, p-value = 0.1
acf(coredata(DInflation),type="correlation",lag=40,ylab="",main="ACF")
acf(coredata(DInflation),type="partial",lag=40,ylab="",main="PACF")
#estimate model
ARMA1<-arima(DInflation,order=c(2,0,2))
ARMA1
##
## Call:
## arima(x = DInflation, order = c(2, 0, 2))
##
## Coefficients:
## ar1 ar2 ma1 ma2 intercept
## 0.8334 -0.4056 -1.3861 0.4866 -0.0001
## s.e. 0.0902 0.0422 0.0935 0.0845 0.0014
##
## sigma^2 estimated as 0.04195: log likelihood = 116.84, aic = -221.69
tsdiag(ARMA1,gof.lag=12)
AIC(ARMA1)
## [1] -221.6853
BIC(ARMA1)
## [1] -194.3361
ARMA2<-arima(DInflation,order=c(4,0,4))
ARMA2
##
## Call:
## arima(x = DInflation, order = c(4, 0, 4))
##
## Coefficients:
## Warning in sqrt(diag(x$var.coef)): NaNs produced
## ar1 ar2 ar3 ar4 ma1 ma2 ma3 ma4
## 0.2235 -0.0639 -0.9355 0.1648 -0.8892 0.0683 0.8796 -0.7862
## s.e. NaN NaN NaN NaN NaN NaN NaN NaN
## intercept
## -0.0001
## s.e. 0.0012
##
## sigma^2 estimated as 0.03437: log likelihood = 185.6, aic = -351.2
tsdiag(ARMA2,gof.lag=12)
AIC(ARMA2)
## [1] -351.2006
BIC(ARMA2)
## [1] -305.6186
ARMA3<-arima(DInflation,order=c(4,0,2))
ARMA3
##
## Call:
## arima(x = DInflation, order = c(4, 0, 2))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ma1 ma2 intercept
## -0.0779 -0.1430 -0.2862 -0.1847 -0.5052 -0.1853 -0.0001
## s.e. 0.2349 0.0788 0.0531 0.0645 0.2358 0.1995 0.0014
##
## sigma^2 estimated as 0.04086: log likelihood = 126.04, aic = -236.08
tsdiag(ARMA3,gof.lag=12)
AIC(ARMA3)
## [1] -236.0777
BIC(ARMA3)
## [1] -199.6122
ARMA4<-arima(DInflation,order=c(3,0,3))
ARMA4
##
## Call:
## arima(x = DInflation, order = c(3, 0, 3))
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 ma3 intercept
## 1.1692 -1.1706 0.1733 -1.8337 1.8012 -0.8276 -0.0001
## s.e. 0.0493 0.0491 0.0489 0.0292 0.0340 0.0276 0.0012
##
## sigma^2 estimated as 0.03439: log likelihood = 185.55, aic = -355.1
tsdiag(ARMA4,gof.lag=12)
AIC(ARMA4)
## [1] -355.1039
BIC(ARMA4)
## [1] -318.6383
# #based on AIC values and Ljunb-Box tests we choose model 3
#forecasting
library("forecast")
## Loading required package: timeDate
## This is forecast 6.2
ARMA3.fcast<-forecast(ARMA3,h=12)
plot(ARMA3.fcast,type="o",pch=16,xlim=c(2005,2017),ylim=c(-0.5,0.5))