Neural Network Models for Time Series Prediction

Feed-forward neural networks with a single hidden layer and lagged inputs for forecasting univariate time series.The nnetar function in the forecast package for R fits a neural network model to a time series with lagged values of the time series as inputs (and possibly some other exogenous inputs). So it is a nonlinear autogressive model, and it is not possible to analytically derive prediction intervals. Therefore we use simulation.The neural networks is fit by the function:

nnetar(y, p, P = 1, size, repeats = 20, xreg = NULL, lambda = NULL, model = NULL, subset = NULL, scale.inputs = TRUE, x = y, …)

Suppose we fit a NNETAR model to the famous Canadian lynx data:

library(ggthemes)
library(forecast)
library(tidyverse)
library(tseries)
library(lubridate)
library(timetk)
library(readxl)
library(tidyquant)
library(scales)
library(forecast)   #  forecasting pkg
library(sweep)   # Broom tidiers for forecast pkg
library(broom)
library(tibble)
library(stringr)
library(highcharter)
library(knitr)
theme_set(theme_bw())

modeldat=read_excel("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/Time series/Manpower/actvsfcdatajuly311.xlsx",sheet = "Sheet1")%>%na.omit()%>%dplyr::rename(Date=Week)
modeldat%>%head()

modeldat%>%ggplot(aes(Date,Actual))+geom_line()+
geom_point(alpha = 0.5, color = palette_light()[[1]], shape=20,size=2) +
   labs(title = "Durability Vehicles Forecasting", x = "Time", y = "Number of \n Durability Vehicles",
       subtitle = "data from 2016 & 2017") +
    theme_tq()

#Data has to be numeric for the tsclean function to work
dat=read_excel("/Users/nanaakwasiabayieboateng/Documents/memphisclassesbooks/Time series/Manpower/actvsfcdatajuly311.xlsx",sheet = "Sheet1")%>%dplyr::rename(Date=Week)%>%na.omit()
dat_ts = ts(dat["Actual"])

(fit <- nnetar(dat_ts, lambda=0.5))

Series: dat_ts 
Model:  NNAR(1,1) 
Call:   nnetar(y = dat_ts, lambda = 0.5)

Average of 20 networks, each of which is
a 1-1-1 network with 4 weights
options were - linear output units 

sigma^2 estimated as 0.2705

I’ve used a Box-Cox transformation with \(\lambda=0.5\)

to ensure the residuals will be roughly homoscedastic.

The model can be written as

\(y_{i}=f(y_{t-i})+\epsilon_{t}\)

where \(\textbf{y}_{t-1}=(y_{t-1},y_{t-2},\cdots,y_{t-8})\prime\) is a vector containing lagged values of the series, and \(f\) is a neural network with 4 hidden nodes in a single layer.\

The error series \(\{\epsilon_{t} \}\) is assumed to be homoscedastic (and possibly also normally distributed).

We can simulate future sample paths of this model iteratively, by randomly generating a value for \(\epsilon_{t}\), either from a normal distribution, or by resampling from the historical values.

So if \(\{\epsilon_{T+1}^{*} \}\) is a random draw from the distribution of errors at time \(T+1\) , then

\(y_{T+1}^{*}=f(\textbf{y}_{T+1})+\epsilon_{T+1}^{*}\)

is one possible draw from the forecast distribution for \(\textbf{y}_{T+1}\). Setting

\(\textbf{y}_{T+1}^{*}=(y^{*}_{T+1},y_{T},\cdots,y_{T-6})\prime\) we can then repeat the process to get

\(y_{T+2}=f(\textbf{y}_{T+1}^{*})+\epsilon_{T+1}^{*}\)

In this way, we can iteratively simulate a future sample path. By repeatedly simulating sample paths, we build up knowledge of the distribution for all future values based on the fitted neural network. Here is a simulation of 9 possible future sample paths for the lynx data. Each sample path covers the next 20 years after the observed data.

sim <- ts(matrix(0, nrow=20, ncol=5), start=end(dat_ts)[1]+1)
for(i in seq(5))
  sim[,i] <- simulate(fit, nsim=20)
library(ggplot2)
autoplot(dat_ts) + forecast::autolayer(sim)

For a multivariate timeseries, specify a seriesname for each timeseries. Defaulting to column names.

If we do this a few hundred or thousand times, we can get a very good picture of the forecast distributions. This is how the forecast.nnetar function produces prediction intervals:

fcast <- forecast(fit, PI=TRUE, h=20)
autoplot(fcast)

Because it is a little slow, PI=FALSE is the default, so prediction intervals are not computed unless requested. The npaths argument in forecast.nnetar controls how many simulations are done (default 1000). By default, the errors are drawn from a normal distribution. The bootstrap argument allows the errors to be “bootstrapped” (i.e., randomly drawn from the historical errors).

sweep::sw_tidy(fit)

sweep::sw_glance(fit) 

sweep::sw_augment(fit )%>%head()

#fcast%>%head()%>%kable()
forecast::forecast(fit,h=12,level = c(90, 95))%>%tk_tbl()