Homework 3 Problem 3
Abeer A lthawad
February 14, 2017
Import data:
time series data for Real Gross Domestic Product
library(Quandl)
library(forecast)
Quandl.api_key("sfE3eBMHtyGwpQTDwkWQ")
rGDP <- Quandl("FRED/GDPC1", type="zoo")Plots
plot(rGDP, type= "l",xlab= "Years", ylab="rGDP", main="Real Gross Domestic Product", major.format="%Y Q%q")
Two Part of The Sample
fstQ <- 1947.25 # 1947Q2
lstQ <- 2008.75 # 2008Q4Stationary Time Series Data
dlrGDP <- diff(log(rGDP))
dlrGDPp1 <- window(dlrGDP, end="2008 Q4")
dlrGDPp2 <- window(dlrGDP, start="2009 Q1")plot(dlrGDP, type= "l", xlab= "Years", ylab="dlrGDP", main="log-change in Real Gross Domestic Product", major.format="%Y Q%q")
ACF
Acf(dlrGDPp1 ,type='correlation',lag=24, main="ACF")
PACF
Acf(dlrGDPp1,type='partial',lag=24, main="PACF")
Find the best model
m1 <- auto.arima(dlrGDPp1, seasonal=FALSE, stationary=TRUE, stepwise=FALSE, ic="aicc")
m1## Series: dlrGDPp1
## ARIMA(2,0,2) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 ma2 intercept
## 1.2808 -0.6706 -0.9728 0.4962 0.0081
## s.e. 0.1719 0.1646 0.2083 0.1792 0.0008
##
## sigma^2 estimated as 8.294e-05: log likelihood=812.48
## AIC=-1612.96 AICc=-1612.61 BIC=-1591.91m1.f <- forecast(m1, length(dlrGDPp2))
plot(m1.f, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03),
main="ARMA(2,2) Model")
lines(m1.f$mean, type="p", pch=16, lty="dashed", col="blue")
lines(dlrGDP, type="o", pch=16, lty="dashed")
AR(1)
m2 <- Arima(dlrGDPp1, order=c(1,0,0))
m2## Series: dlrGDPp1
## ARIMA(1,0,0) with non-zero mean
##
## Coefficients:
## ar1 intercept
## 0.3607 0.0081
## s.e. 0.0605 0.0009
##
## sigma^2 estimated as 8.459e-05: log likelihood=808.6
## AIC=-1611.19 AICc=-1611.09 BIC=-1600.66m2.f<- forecast(m2, length(dlrGDPp2))
plot(m2.f, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03),
main="AR(1)")
lines(m1.f$mean, type="p", pch=16, lty="dashed", col="blue")
lines(dlrGDP, type="o", pch=16, lty="dashed")
par(mar = rep(2, 4))
tsdiag(m1, gof.lag = 20)
Inverted AR & MA
plot.Arima(m1)
Thus ARMA(2,2) is stationary and invertible.
Generate 1 to 32 step ahead forecast(Recursive Scheme)
m1.f.rec <- zoo()
for(i in 1:length(dlrGDPp2))
{
m <- window( dlrGDP, end=lstQ+(i-1)/4 )
m1.updt <- arima(m, order=c(2,0,2))
m1.f.rec <- c(m1.f.rec, forecast(m1.updt, 1)$mean)
}
m1.f.rec <- as.ts(m1.f.rec)
accuracy(m1.f.rec, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.002300581 0.005038677 0.003921268 -99.50147 166.961m2.f.rec <- zoo()
for(i in 1:length(dlrGDPp2))
{
n <- window( dlrGDP, end=lstQ+(i-1)/4 )
m2.updt <- arima(n, order=c(1,0,0))
m2.f.rec <- c(m2.f.rec, forecast(m2.updt, 1)$mean)
}
m2.f.rec <- as.ts(m2.f.rec)
accuracy(m2.f.rec, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.001862434 0.004853304 0.003965669 -94.9437 170.2363Rolling scheme
m1.f.rol <- zoo()
for(i in 1:length(dlrGDPp2))
{
y <- window( dlrGDP, start=fstQ+(i-1)/4, end=lstQ+(i-1)/4 )
m1.updt <- arima(y, 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.002202307 0.004997696 0.003901092 -96.81412 164.633m2.f.rol <- zoo()
for(i in 1:length(dlrGDPp2))
{
x <- window( dlrGDP, start=fstQ+(i-1)/4, end=lstQ+(i-1)/4 )
m2.updt <- arima(x, order=c(1,0,0))
m2.f.rol <- c(m2.f.rol, forecast(m2.updt, 1)$mean)
}
m2.f.rol <- as.ts(m2.f.rol)
accuracy(m2.f.rol, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.001829009 0.004835928 0.003944857 -94.16649 169.6169plot(m1.f, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03), main="ARMA(2,2) Forecasts")
lines(m1.f.rec, type="p", pch=16, lty="dotdash", col="blue2")
lines(m1.f.rol, type="p", pch=16, lty="twodash", col="brown2")
lines(dlrGDP, type="o", pch=16, lty="dashed", col="darkolivegreen4")
legend("topleft", legend=c("Recursive", "Rolling", "In-Sample"), col=c("blue2", "brown2", "darkolivegreen4"), lty=c("dotdash", "twodash", "dashed"))
plot(m2.f, type="o", pch=16, xlim=c(2005,2016), ylim=c(-0.03,0.03), main="AR(1) Forecasts")
lines(m2.f.rec, type="p", pch=16, lty="dotdash", col="blue2")
lines(m2.f.rol, type="p", pch=16, lty="twodash", col="brown2")
lines(dlrGDP, type="o", pch=16, lty="dashed", col="darkolivegreen4")
legend("topleft", legend=c("Recursive", "Rolling", "In-Sample"), col=c("blue2", "brown2", "darkolivegreen4"), lty=c("dotdash", "twodash", "dashed")) 
AR(1) is more stable than ARMA(2,2)
Accuracy
In-sample
ARMA(2,2)
accuracy(m1)## 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.01057882AR(1) model
accuracy(m2)## 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.03047387Out-of-Sample
Recursive Scheme
m1.f.rec <- zoo()
for(i in 1:length(dlrGDPp2))
{
m <- window( dlrGDP, end=lstQ+(i-1)/4 )
m1.updt <- arima(m, order=c(2,0,2))
m1.f.rec <- c(m1.f.rec, forecast(m1.updt, 1)$mean)
}
m1.f.rec <- as.ts(m1.f.rec)
accuracy(m1.f.rec, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.002300581 0.005038677 0.003921268 -99.50147 166.961m2.f.rec <- zoo()
for(i in 1:length(dlrGDPp2))
{
n <- window( dlrGDP, end=lstQ+(i-1)/4 )
m2.updt <- arima(n, order=c(1,0,0))
m2.f.rec <- c(m2.f.rec, forecast(m2.updt, 1)$mean)
}
m2.f.rec <- as.ts(m2.f.rec)
accuracy(m2.f.rec, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.001862434 0.004853304 0.003965669 -94.9437 170.2363Rolling Scheme
accuracy(m1.f.rol, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.002202307 0.004997696 0.003901092 -96.81412 164.633accuracy(m2.f.rol, dlrGDPp2)## ME RMSE MAE MPE MAPE
## Test set -0.001829009 0.004835928 0.003944857 -94.16649 169.6169AR(1) has less error than ARMA(2,2).
Conclusion
AR(1) has an accurate estimation, but ARMA(2,2) overestimating the model, thus AR(1)is the best model.