In this problem we test “Civilian Unemployment Rate” for nonstationarity using ADF and KPPS tests. Let’s first import and plot the 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
Quandl.api_key("t1BuHRnYBbU8N1_qwu5b")
CUR <- Quandl("FRED/UNRATE", type="zoo")
str(CUR)
## 'zooreg' series from Jan 1948 to Jan 2016
## Data: num [1:817] 3.4 3.8 4 3.9 3.5 3.6 3.6 3.9 3.8 3.7 ...
## Index: Class 'yearmon' num [1:817] 1948 1948 1948 1948 1948 ...
## Frequency: 12
head(CUR)
## Jan 1948 Feb 1948 Mar 1948 Apr 1948 May 1948 Jun 1948
## 3.4 3.8 4.0 3.9 3.5 3.6
tail(CUR)
## Aug 2015 Sep 2015 Oct 2015 Nov 2015 Dec 2015 Jan 2016
## 5.1 5.1 5.0 5.0 5.0 4.9
plot(CUR, xlab="", ylab="", main="Civilian Unemployment Rate (CUR)")
And let’s estimate and plot the difference and log difference of series:
dCUR<-diff(CUR)
lCUR<-log(CUR)
dlCUR<-diff(lCUR)
par(mfrow=c(1,2))
plot(dlCUR,xlab="",ylab="",main="First Differance of CUR (dCUR)")
plot(dlCUR,xlab="",ylab="",main="Log Differance of CUR (dlCUR)")
And ADF and KPSS tests results are:
ADF and KPSS tests are:
library(tseries)
adf.test(dCUR)
## Warning in adf.test(dCUR): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: dCUR
## Dickey-Fuller = -8.2434, Lag order = 9, p-value = 0.01
## alternative hypothesis: stationary
kpss.test(dCUR)
## Warning in kpss.test(dCUR): p-value greater than printed p-value
##
## KPSS Test for Level Stationarity
##
## data: dCUR
## KPSS Level = 0.064885, Truncation lag parameter = 6, p-value = 0.1
adf.test(dlCUR)
## Warning in adf.test(dlCUR): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: dlCUR
## Dickey-Fuller = -8.3101, Lag order = 9, p-value = 0.01
## alternative hypothesis: stationary
kpss.test(dlCUR)
## Warning in kpss.test(dlCUR): p-value greater than printed p-value
##
## KPSS Test for Level Stationarity
##
## data: dlCUR
## KPSS Level = 0.067206, Truncation lag parameter = 6, p-value = 0.1
At 10% level we can reject the null hypothesis of non-stationarity. Both first difference and log difference are giving similar results. Now let’s derive ACF and PACF, necessary for forecasting:
par(mfrow=c(2,2))
acf(as.data.frame(dCUR),type='correlation',lag=24, main="dCUR")
acf(as.data.frame(dCUR),type='partial',lag=24, main="dCUR")
acf(as.data.frame(dlCUR),type='correlation',lag=24, main="dlCUR")
acf(as.data.frame(dlCUR),type='partial',lag=24, main="dlCUR")
Based on ACF and PACF, I’d like to estimate ARiMA model of:
m1 <- arima(dCUR,order=c(3,1,3),seasonal=list(order=c(1,0,1),period=12))
m1
##
## Call:
## arima(x = dCUR, order = c(3, 1, 3), seasonal = list(order = c(1, 0, 1), period = 12))
##
## Coefficients:
## ar1 ar2 ar3 ma1 ma2 ma3 sar1 sma1
## -1.8954 -1.3042 -0.1897 0.9374 -0.3106 -0.7162 0.4866 -0.8045
## s.e. 0.0848 0.1316 0.0534 0.0816 0.0766 0.0760 0.0694 0.0491
##
## sigma^2 estimated as 0.0373: log likelihood = 181.35, aic = -344.69
tsdiag(m1,gof.lag=12)
Forecasting
library(forecast)
## Loading required package: timeDate
## This is forecast 6.2
m1.fcast <- forecast(m1, h=24)
plot(m1.fcast, xlim=c(2014, 2017))