### Environment
# Clear the workspace
rm(list = ls()) # Clear environment
gc() # Clear unused memory
## used (Mb) gc trigger (Mb) limit (Mb) max used (Mb)
## Ncells 521355 27.9 1154940 61.7 NA 669291 35.8
## Vcells 959623 7.4 8388608 64.0 16384 1840004 14.1
cat("\f") # Clear the console
dev.off # Clear the charts
## function (which = dev.cur())
## {
## if (which == 1)
## stop("cannot shut down device 1 (the null device)")
## .External(C_devoff, as.integer(which))
## dev.cur()
## }
## <bytecode: 0x7fb4e1ad25d0>
## <environment: namespace:grDevices>
#make vector containing the packages
packages <- c("ggplot2",
"knitr",
"fpp3",
"forecast",
"dplyr",
"lubridate",
"grid",
"gridExtra")
#Install packages
for (i in 1:length(packages)) {# Iterate through each package in the 'packages' vector.
if (!packages[i] %in% rownames(installed.packages())) {#check if package is installed
install.packages(packages[i], #install package if it isn't already
repos = "http://cran.rstudio.com/",
dependencies = TRUE)}
library(packages[i], character.only = TRUE)} # Load the package into the R session.
## ── Attaching packages ────────────────────────────────────────────── fpp3 0.5 ──
## ✔ tibble 3.2.1 ✔ tsibbledata 0.4.1
## ✔ dplyr 1.1.1 ✔ feasts 0.3.2
## ✔ tidyr 1.3.0 ✔ fable 0.3.3
## ✔ lubridate 1.9.2 ✔ fabletools 0.4.1
## ✔ tsibble 1.1.4
## ── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
## ✖ lubridate::date() masks base::date()
## ✖ dplyr::filter() masks stats::filter()
## ✖ tsibble::intersect() masks base::intersect()
## ✖ tsibble::interval() masks lubridate::interval()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tsibble::setdiff() masks base::setdiff()
## ✖ tsibble::union() masks base::union()
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
# Remove unnecessary values from environment
rm(i,packages)
Source: FRED Economic Data
The core dataset includes quarterly data from the following (as RGDP did not have monthly data:
Unemployment Rate, Percent, Not Seasonally Adjusted (UNRATENSA)
Personal Consumption Expenditures: Chain-type Price Index, Index 2017=100, Seasonally Adjusted (PCEPI)
Real Gross Domestic Product, Billions of Chained 2017 Dollars, Seasonally Adjusted Annual Rate (GDPC1)
### Manipulate Data
# Download data to use for analysis
core_df <-
read.csv("/Users/cam/Desktop/PREDICTIVE ANALYTICS/Discussions/Data/Neural_Net_Data.csv")
# Filter down to last 10 years
core_df$DATE <- as.Date(core_df$DATE) # converts to date format
# Rename variables to simpler terms
core_df <- core_df %>%
rename(UNEM = UNRATENSA,
PCE = PCEPI_PC1,
RGDP = GDPC1_PC1)
# Convert character vars to numeric
core_df$PCE <- as.numeric(core_df$PCE)
## Warning: NAs introduced by coercion
core_df$RGDP <- as.numeric(core_df$RGDP)
## Warning: NAs introduced by coercion
# Drop 2024 quarter 1 data (contains null values)
core_df <- core_df %>%
filter(!(year(DATE) == 2024 & quarter(DATE) == 1))
# Filter to past 10 years
core_df <- core_df %>%
filter(DATE >= max(core_df$DATE) - years(10))
# Split into training and testing data
train <- subset(core_df, DATE <= max(core_df$DATE) - years(2))
test <- subset(core_df, DATE > max(core_df$DATE) - years(2))
## Convert to time series objects
unem_train <- ts(train$UNEM,
start = c(2013,4),
frequency = 4)
unem_test <- ts(test$UNEM,
start = c(2022,1),
frequency = 4)
PCE_train <- ts(train$PCE,
start = c(2013,4),
frequency = 4)
PCE_test <- ts(test$PCE,
start = c(2022,1),
frequency = 4)
RGDP_train <- ts(train$RGDP,
start = c(2013,4),
frequency = 4)
RGDP_test <- ts(test$RGDP,
start = c(2022,1),
frequency = 4)
### Display the data
ggplot() +
geom_line(data = core_df, aes(x = DATE, y = UNEM, color = "Unemployment Rate"), linewidth = 1) +
geom_line(data = core_df, aes(x = DATE, y = PCE, color = "PCE")) +
geom_line(data = core_df, aes(x = DATE, y = RGDP, color = "RGDP Growth Rate")) +
labs(title = "Quarterly Data (Past 10 Years)",
x = "Date",
y = "Percent, Percent Change from 1 Year Ago",
color = "Variable") +
scale_color_manual(values = c("darkgreen", "red", "darkblue"))
# Create Neural Net model
nnet_mod <- nnetar(unem_train,
xreg = cbind(PCE_train,RGDP_train))
# forecast ARIMA model
nnet_forecast <- forecast(nnet_mod,
xreg = cbind(PCE_test,RGDP_test),
h = 8)
## Warning in forecast.nnetar(nnet_mod, xreg = cbind(PCE_test, RGDP_test), : xreg
## contains different column names from the xreg used in training. Please check
## that the regressors are in the same order.
# Display stats for model on training data
summary(nnet_mod)
## Length Class Mode
## x 33 ts numeric
## m 1 -none- numeric
## p 1 -none- numeric
## P 1 -none- numeric
## scalex 2 -none- list
## scalexreg 2 -none- list
## size 1 -none- numeric
## xreg 66 mts numeric
## subset 33 -none- numeric
## model 20 nnetarmodels list
## nnetargs 0 -none- list
## fitted 33 ts numeric
## residuals 33 ts numeric
## lags 2 -none- numeric
## series 1 -none- character
## method 1 -none- character
## call 3 -none- call
# Display Neural Net model
autoplot(unem_train)+
autolayer(nnet_forecast)+
labs(x = "Date",
y = "Unemployment Rate %",
title = "Neural Net Forecast") +
autolayer(unem_test,
series = "True Values")
# Evaluate accuracy
# Make combined df
nnet_data <- nnet_forecast$mean
nnet_data <- data.frame(unem_test,nnet_data)
# Rename variables to proper terms
nnet_data <- nnet_data %>%
rename(True_Vals = unem_test,
Forecast_Vals = nnet_data)
# Residuals
nnet_data$residuals <- nnet_data$True_Vals - nnet_data$Forecast_Vals
# Create ARIMA model
arima_mod <- auto.arima(unem_train)
# forecast ARIMA model
arima_forecast <- forecast(arima_mod, h = 8)
# Display stats for model on training data
summary(arima_mod)
## Series: unem_train
## ARIMA(1,0,0) with non-zero mean
##
## Coefficients:
## ar1 mean
## 0.4914 5.2618
## s.e. 0.1499 0.5206
##
## sigma^2 = 2.607: log likelihood = -61.74
## AIC=129.48 AICc=130.31 BIC=133.97
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.02704315 1.564911 0.7665598 -5.976724 13.60213 0.4965801
## ACF1
## Training set 0.04187761
# Display ARIMA model
autoplot(unem_train)+
autolayer(arima_forecast)+
labs(x = "Date",
y = "Unemployment Rate %",
title = "ARIMA Forecast") +
autolayer(unem_test,
series = "True Values")
# Evaluate accuracy
# Make combined df
arima_data <- arima_forecast$mean
arima_data <- data.frame(unem_test,arima_data)
# Rename variables to proper terms
arima_data <- arima_data %>%
rename(True_Vals = unem_test,
Forecast_Vals = arima_data)
# Residuals
arima_data$residuals <- arima_data$True_Vals - arima_data$Forecast_Vals
# Create ets model
ets_mod <- ets(unem_train)
# forecast ets model
ets_forecast <- forecast(ets_mod, h = 8)
# Display stats for model on training data
summary(ets_mod)
## ETS(A,N,N)
##
## Call:
## ets(y = unem_train)
##
## Smoothing parameters:
## alpha = 0.6072
##
## Initial states:
## l = 6.6523
##
## sigma: 1.779
##
## AIC AICc BIC
## 157.3425 158.1701 161.8320
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.1020116 1.724288 0.7741032 -5.251423 12.22742 0.5014667
## ACF1
## Training set 0.08687058
# Display ETS model
autoplot(unem_train)+
autolayer(ets_forecast)+
labs(x = "Date",
y = "Unemployment Rate %",
title = "ETS Forecast") +
autolayer(unem_test,
series = "True Values")
# Evaluate accuracy
# Make combined df
ets_data <- ets_forecast$mean
ets_data <- data.frame(unem_test,ets_data)
# Rename variables to proper terms
ets_data <- ets_data %>%
rename(True_Vals = unem_test,
Forecast_Vals = ets_data)
# Residuals
ets_data$residuals <- ets_data$True_Vals - ets_data$Forecast_Vals
# Display graph of all forecasts
autoplot(unem_train,
ylab = "Percent, Percent change from a year ago",
size = 0.75) +
ggtitle("Inflation (PCE) Forecasts") +
autolayer(unem_test,
series = "True Values",
size = 0.75) +
autolayer(nnet_forecast$mean,
series = "Neural Network") +
autolayer(arima_forecast$mean,
series = "ARIMA") +
autolayer(ets_forecast$mean,
series = "ETS") +
scale_colour_manual(values = c("red", "blue", "purple", "darkgreen"))
## Create stats table describing forecast performances
# Stats for nnet
nnet_table <- data.frame(
RMSE = sqrt(mean(nnet_data$residuals^2)),
MAE = mean(abs(nnet_data$residuals)),
MAPE = mean(abs((nnet_data$residuals / nnet_data$True_Vals) * 100)),
MASE = mean(abs(nnet_data$residuals) / mean(abs(diff(nnet_data$True_Vals)))),
ME = mean(nnet_data$residuals),
MPE = mean((nnet_data$residuals / nnet_data$True_Vals) * 100)
)
rownames(nnet_table) <- "Neural Network"
#Stats for ARIMA
arima_table <- data.frame(
RMSE = sqrt(mean(arima_data$residuals^2)),
MAE = mean(abs(arima_data$residuals)),
MAPE = mean(abs((arima_data$residuals / arima_data$True_Vals) * 100)),
MASE = mean(abs(arima_data$residuals) / mean(abs(diff(arima_data$True_Vals)))),
ME = mean(arima_data$residuals),
MPE = mean((arima_data$residuals / arima_data$True_Vals) * 100)
)
rownames(arima_table) <- "ARIMA"
#Stats for ETS
ets_table <- data.frame(
RMSE = sqrt(mean(ets_data$residuals^2)),
MAE = mean(abs(ets_data$residuals)),
MAPE = mean(abs((ets_data$residuals / ets_data$True_Vals) * 100)),
MASE = mean(abs(ets_data$residuals) / mean(abs(diff(ets_data$True_Vals)))),
ME = mean(ets_data$residuals),
MPE = mean((ets_data$residuals / ets_data$True_Vals) * 100)
)
rownames(ets_table) <- "ETS"
#Merge
stats_table <- rbind(nnet_table,arima_table,ets_table)
kable(stats_table,
caption = "**Stats for Testing Data**")
| RMSE | MAE | MAPE | MASE | ME | MPE | |
|---|---|---|---|---|---|---|
| Neural Network | 1.0861828 | 0.8792349 | 24.61737 | 2.600554 | -0.8792349 | -24.61737 |
| ARIMA | 1.5144492 | 1.4642001 | 40.91382 | 4.330733 | -1.4642001 | -40.91382 |
| ETS | 0.9921087 | 0.9664736 | 27.00043 | 2.858584 | -0.9664736 | -27.00043 |
Looking at the results from each of the models, it is apparent that the best performing model was the ETS (A,N,N) model. The reason for this may exist outside the predictive power of a neural network in that the unemployment rate appears to have stabilized around its natural rate. The neural network appears to pick up on the trend of unemployment spiking and gradually coming down. Thus, the shape of the predicted neural net forecast mirrors said trend, Given that the economic activity during covid was fairly abnormal, it may have thrown off the forecast.