Introduction

Quarterly Real Gross Domestic Product

library(Quandl)
## Warning: package 'Quandl' was built under R version 3.3.3
## 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(zoo)
library(xts)
library(dygraphs)
library(knitr)
library(forecast)
## Loading required package: timeDate
## This is forecast 7.3
Quandl.api_key('rxVQhZ8_nxxeo2yy4Uz4')
GDP <- Quandl("FRED/GDPC1", type="ts")

DATA

the time series with log changes in Real Gross Domestic Product \(y_t=\delta\log y_t=log GDP_t - log GDP_{t-1}\) where \(c_t\) is the original quarterly Real Personal Consumption Expenditures.

y.all <- 100* diff(log(GDP),lag=1)

To be able to evaluate the out of sample forecasting performance of the model, we restrict the estimation sample to 1947q1-2008q4.

y <- window(y.all, start=c(), end=c(2008,4))
plot(y, main="Quarter-over-quarter Log Change in GDP")

Methodology

The time series of first differences looks stationary. Now we need to check both ACF and PACF to suggest an ARIMA for analysis.

library(forecast)
acf(y, type="correlation", xlab="Lag", ylab="ACF", main="")

According to ACF graph, we can assume that our time series might be MA(2).

library(forecast)
acf(y, type="partial", xlab="Lag", ylab="PACF", main="")

PACF has decaying pattern. So our model most likely should include both AR and MA components.We will use “auto.arima” function using AICc.

arima.aicc <- auto.arima(y, seasonal = FALSE, stationary = TRUE, stepwise = FALSE, ic="aicc")
arima.aicc
## Series: y 
## ARIMA(2,0,2) with non-zero mean 
## 
## Coefficients:
##          ar1      ar2      ma1     ma2  intercept
##       1.2812  -0.6710  -0.9733  0.4965     0.8122
## s.e.  0.1716   0.1645   0.2080  0.1793     0.0771
## 
## sigma^2 estimated as 0.8294:  log likelihood=-325
## AIC=661.99   AICc=662.34   BIC=683.05

Now we check stationarity and invertibility

plot.Arima(arima.aicc)

The model specifications display stationarity and invertibility of their AR and MA roots.

tsdiag(arima.aicc)

As we can see on the graph, there is no autocorrelation between residuals and they are normally distributed.

Forecasting

For forecasting we restrict the sample to 1947q1-2008q4 and use to create a forecast for the 2009q9-2016q4 period. As shown in the figure below, the forecast is reasonably accurate and the actual quarterly GDP after 2009q1 lies in the 90% confidence interval.

library(Quandl)
rGDP <- Quandl("FRED/GDPC1", type="zoo")
dlrGDP <- diff(log(rGDP))
fstQ <- 1947.25 # 1947Q2
lstQ <- 2008.75 # 2008Q4
dlrGDPp1 <- window(dlrGDP, end=lstQ)
dlrGDPp2 <- window(dlrGDP, start=lstQ+0.25)
library(forecast)
arima202 <- Arima(dlrGDPp1, order=c(2,0,2))
arma22.f.h <- forecast(arima202, length(dlrGDPp2))
plot(arma22.f.h, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03),
main="ARIMA(2.0.2) Model Multistep Forecast - U.S. Real GDP Growth Rate")
lines(arma22.f.h$mean, type="p", pch=16, lty="dashed", col="blue")
lines(dlrGDP, type="o", pch=16, lty="dashed")

Below is the same forecast, but with a retrospective to 1947.

GDP.fcst <- forecast(arima.aicc, h=32)
plot(GDP.fcst, xlim=c(1947,2016))

Rolling scheme forecast

library(forecast)
m1.f.rol <- zoo()
for(i in 1:length(dlrGDPp2))
{
  y.rsc <- window( dlrGDP, start=fstQ+(i-1)/4, end=lstQ+(i-1)/4 )
m1.updt <- arima(y.rsc, order=c(2,0,2))
m1.f.rol <- c(m1.f.rol, forecast(m1.updt, 1)$mean)
}
m1.f.rol <- as.ts(m1.f.rol)
accuracy(m1.f.rol, dlrGDPp2)
##                    ME        RMSE         MAE       MPE     MAPE
## Test set -0.002186505 0.004989311 0.003885289 -96.31497 164.1338

The forecast of AR(1) model

library(forecast)
ar1 <- Arima(dlrGDPp1, order=c(1,0,0))
ar1.f.h <- forecast(ar1, length(dlrGDPp2))
plot(ar1.f.h, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03),
main="AR(1) Model Multistep Forecasts - U.S. Real GDP Growth Rate")
lines(ar1.f.h$mean, type="p", pch=16, lty="dashed", col="blue")
lines(dlrGDP, type="o", pch=16, lty="dashed")

Out-of-Sample Evaluation of Accuracy

Now we compare the out of sample accuracy of the ARIMA(2.0.2) multistep forecast, rolling scheme and compare them with the accuracy of the AR(1) model

accuracy(arma22.f.h)
##                        ME        RMSE         MAE       MPE     MAPE
## Training set 2.372213e-05 0.009014554 0.006668165 -55.80138 174.3994
##                   MASE       ACF1
## Training set 0.6264638 0.01057882
accuracy(ar1.f.h)
##                        ME        RMSE        MAE       MPE     MAPE
## Training set 1.597568e-05 0.009160185 0.00676911 -46.93084 168.4821
##                   MASE        ACF1
## Training set 0.6359474 -0.03047387
accuracy(m1.f.rol, dlrGDPp2)
##                    ME        RMSE         MAE       MPE     MAPE
## Test set -0.002186505 0.004989311 0.003885289 -96.31497 164.1338

Rolling scheme forecast is the best by MAPE, MAE, RMSE, ME. However, Ar(1) has the least Mean absolute percentage error.

Now we plot the two sets of ARIMA forecasts to compare them visually.

plot(arma22.f.h, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03),
main="ARIMA(2.0.2) vs. Rolling Scheme forecasts - U.S. Real GDP Growth Rate")
lines(arma22.f.h$mean, type="p", pch=16, lty="dashed", col="blue")
lines(m1.f.rol, type="p", pch=16, lty="dashed", col="red")
lines(dlrGDP, type="o", pch=16, lty="dashed")

Conclusion

We used auto.arima function to find the best specification for the model of first differences in Real Gross Domestic Product. The suggested specification is (2.0.2), which is also coherent with ACF and PACF. We conclude that rolling scheme ARIMA(2.0.2) has the best forecast accuracy comparing to AR(1) and ARIMA(2.0.2) multistep forecasting.