Neural models: MLP, RNN, LSTM and GRU

1. Why neural networks for forecasting

The claim the deck opens with is strong and precisely worded: neural networks are very effective in identifying non linear models. Even with input data deriving from complex, nonlinear processes, NN are able to combine linear and nonlinear models to get whichever accuracy level is requested — in prediction. That parenthesis is doing a lot of work, and it is the honest one: the guarantee is about fitting, not about forecasting.

Three capabilities are listed:

• NN can learn patterns in linear time series;
• NN can learn patterns in nonlinear time series;
• NN can generalize linear and nonlinear patterns.

Two structural properties follow. NN are nonparametric: they do not make hypotheses on noise distribution (Gaussian or otherwise) — contrast this with the Kalman filter of chapter 3, whose optimality is conditional on Gaussian noise, and with the AIC approximation of chapter 4, which assumes Gaussian residuals. And NN identify a model from input data, essentially defining a model of the generating process.

A vocabulary you must be able to translate

Given the loosely neurobiological inspiration of the model, terminology is peculiar. The deck supplies the dictionary, and exam questions have been known to consist of nothing but this table.

2. Time-lagged networks and the sliding window

The architecture is described in one sentence: the net is structured according to a feedforward scheme, where inputs are associated to past data and the output is the forecast. Input nodes take x_t-1, x_t-2, …, x_t-K plus a bias; the output node produces x_t.

The slides add a note on convolutional nets that is easy to misread, so here it is verbatim: if the input of hidden and output neurons sees only pre-defined subsets of the output of the previous layer and applies them specific operators, the net is named convolutional. These are networks most effective for classification, less so for regression. For forecasting, then, CNNs are mentioned and set aside.

How learning proceeds

The sliding window is the training loop, in five steps:

1. A window on past data is shown in input.
2. The corresponding output is computed.
3. The output is compared with its corresponding actual value.
4. Backpropagation: weights are changed in order to reduce prediction error.
5. The window is shifted forward, presenting the following input.

3. From series to supervised dataset

The task is stated as: transform a data series into a learning dataset where the next datapoint of each record is the unknown. The deck gives three implementations, differing only in the container they expect — array, series, or dataframe.

# converts a 1D array of values into two np arrays
def create_dataset(arrdata, look_back=1):
    dataX, dataY = [], []
    for i in range(len(arrdata) - look_back):
        a = arrdata[i:(i + look_back)]
        dataX.append(a)
        dataY.append(arrdata[i + look_back])
    return np.array(dataX), np.array(dataY)

# converts a series into a windowed dataframe
def create_dataset2(df, look_back=1, colname='value'):
    cols = ['x{}'.format(x) for x in range(look_back)] + ['y']
    df_win = []
    for i in range(look_back, df.shape[0]):
        df_win.append(df.loc[i-look_back:i, colname].tolist())
    return pd.DataFrame(df_win, columns=cols)

# works directly on the dataframe, building it by columns
def create_dataset3(df, look_back=1, colname='value'):
    dataset = pd.DataFrame()
    for i in range(lookback, 0, -1):
        dataset['t-' + str(i)] = df.Passengers.shift(i)
    dataset['t'] = df.values
    dataset = dataset[lookback:]   # removes the first lookback rows

4. What each architecture is, statistically

The deck states, without ceremony, that NN are able to approximate any function, and that forecasting can be thought of as the identification of a specific (nonlinear) function. It then gives a three-row table that maps each architecture onto its statistical twin — the single most useful table in the chapter, because it lets you carry every intuition from chapter 4 across the divide.

Read the first row and the reason becomes obvious: a single neuron with a linear activation computes a weighted sum of its inputs plus a bias, and its inputs are the last p observations. That is an AR(p). Add a hidden layer with nonlinear activations and you get a nonlinear combination of such sums. Add a feedback loop and the network carries forward an internal state which plays the role of the error memory in an MA term.

The slides illustrate the point with two small examples: a linear regression including exogenous variables — which they label ARX — and a nonlinear (logistic) regression of the same process, showing that the same network topology moves between the two by changing the activation function alone.

5. An MLP in torch, end to end

The worked example is the Box-Jenkins airline series, with a network described as (12-8-1): 12 input nodes, 8 hidden units, 1 output node. The 12 input values are t, t-1, …, t-11 — the last twelve observations — and the prediction value is t+1. The benchmark has 132 input values over 13 time periods of monthly data.

Preliminaries and the split

import numpy as np, pandas as pd, matplotlib.pyplot as plt
import torch, torch.nn as nn, torch.optim as optim
from sklearn.metrics import root_mean_squared_error

df = pd.read_csv('BoxJenkins.csv', usecols=[1])
dataset = df.values.astype('float32')          # COLUMN VECTOR !!!
train_size = int(len(dataset) - 12)            # hold out the last year
train, test = dataset[0:train_size, :], dataset[train_size:len(dataset), :]

look_back = 12
testdata = np.concatenate((train[-look_back:], test))
trainX, trainy = create_dataset(train, look_back)

The line building testdata is worth a second: to predict the first test point the network needs the twelve values that precede it, and those live at the end of the training set. Stitching them on is not leakage — it is giving the model the context it would have in production.

The model and the training loop

X = torch.FloatTensor(trainX)
y = torch.FloatTensor(trainy)

model = nn.Sequential(nn.Linear(12, 10), nn.ReLU(),
                      nn.Linear(10, 8),  nn.ReLU(),
                      nn.Linear(8, 1)   # no activation, allowing negative outputs
                     )

loss_fn   = nn.MSELoss()
optimizer = optim.Adam(model.parameters(), lr=0.001)   # mind the learning rate!
n_epochs, batch_size = 200, 1

for epoch in range(n_epochs):
   for i in range(0, len(X), batch_size):
      Xbatch = X[i:i + batch_size]
      y_pred = model(Xbatch)
      ybatch = y[i:i + batch_size]
      loss   = loss_fn(y_pred, ybatch)
      optimizer.zero_grad()      # clear old gradients
      loss.backward()            # here the learning
      optimizer.step()

Two annotations in the original are the ones that bite in practice. No activation function on the output layer, allowing for negative outputs — essential when the target is a differenced or standardised series, which is centred on zero. And mind the learning rate: with lr too large the loss oscillates and never settles, with lr too small two hundred epochs are not enough. Chapter 9 gives that parameter its proper name, λ, and its proper theory.

Evaluating: RMSE or AIC

The deck offers two alternatives, and the second is a nice bridge back to chapter 4.

# alternative 1: model evaluation as RMSE
model.eval()            # tells the model we are evaluating, not training
with torch.no_grad():   # stop tracking gradient computations
   y_pred = model(X)
   rmse = root_mean_squared_error(trainy, y_pred.numpy())

# alternative 2: model evaluation as AIC
k = sum(p.numel() for p in model.parameters() if p.requires_grad)   # count params
model.eval()
with torch.no_grad():
   y_pred = model(X)
   rss    = torch.sum((y - y_pred) ** 2).item()
   n      = trainX.shape[0]
   aic    = n * np.log(rss / n) + 2 * k

The AIC of a neural network counts every weight and bias as a parameter. For the (12-10-8-1) network above that is 12×10 + 10 + 10×8 + 8 + 8×1 + 1 = 227 parameters, against 132 monthly observations. Chapter 4 warned that with a small sample AIC prefers overfits; here the penalty term alone is 454.

Recursive forecasting

To produce more than one step ahead, the network is fed its own output. This is recursive forecasting, and it is the same pattern XGBoost and random forest will use in chapter 6.

# start from the last window of train (or the first of testX)
input_seq = trainX[-1].tolist()
testForecast = []
model.eval()
with torch.no_grad():
   for _ in range(len(test)):
      input_tensor = torch.tensor(input_seq, dtype=torch.float32).unsqueeze(0)
      pred = model(input_tensor).item()      # scalar prediction
      testForecast.append(pred)
      input_seq = input_seq[1:] + [pred]     # roll the window forward by 1

6. The modelling process, in five phases

The deck lays out the phases for the definition of a neural forecast model. Note that phase 1 is a pointer straight back to chapter 3 — the same as statistical models.

Preprocessing is mandatory here, not optional

The deck lists the procedures needed to improve forecast effectiveness: verification, correction and editing (errors in the data); re-encoding of variables (from categorical to numeric, one-hot); scaling of variables (e.g. linear interval scaling); normalization; selection of independent variables (PCA, …); removal of the outliers; insertion of missing values (means, regression, default, …). All of chapter 3, invoked as a prerequisite.

And then the sentence that turns a recommendation into a requirement: data preprocessing of the input vector is a mandatory requirement for the application of MLPs. The reason is stated precisely: as the sigmoid activation functions in the hidden nodes are only defined in the interval of ]-1, 1[ for the hyperbolic tangent or ]0, 1[ for the logistic function, input data must be scaled to facilitate learning. Some authors recommend linear scaling of data into smaller intervals, e.g. [0.2, 0.8], to avoid saturation effects at the asymptotic bounds of the activation functions.

7. Heuristic choices and sizing rules

How many hidden neurons? The literature disagrees, and the deck simply lists who says what, for a network with n input nodes.

And the activation function:

Three concluding suggestions

• The size of the test set should be about 10% to 30% of that of the training set.
• To avoid overfitting, the training set size should be at least 5 times the number of network weights.
• In theory, a single hidden layer is enough to approximate any continuous function, but in practice a further hidden layer often helps. More than 4 layers (input, output and two hidden) can work wonders but usually only in very controlled contexts.

8. Recurrent networks and BPTT

Recurrent neural networks are introduced as a very effective neural model for forecasting, with five defining properties:

• They have at least one backward connection (feedback loop).
• They can keep activations even with no input → they have an internal state.
• Even with one single neuron, we can input a time series and get a series in the output.
• They can approximate any dynamical system.
• Very complex mathematical analysis, and very complex learning algorithms, hard to control — the most used being Backpropagation Through Time (BTT / BPTT).

The consequence for modelling is stated compactly: RNNs compile preceding data in their context layer, they do not need an input neuron for each past value. That is the structural difference from the time-lagged MLP of section 2, which needs K input nodes for K lags. And then a sentence worth underlining twice: their effect is not function approximation, but process modeling. The prediction model is structured as

outputₜ₊₁ ≈ RNN(NetworkState, inputₜ, outputₜ)

where NetworkState is represented by the activation of the neurons in the context layer.

Activation functions

We can use any activation function we like in the recurrent neural network. Common choices: the sigmoid function, the tanh function, and the ReLU function, max(0, x).

The BPTT(k) algorithm

BPTT(k)      // pₜ is the target value at time step t
1. while not (termination condition)
2.     set all h to zero                     // hidden neurons states
3.     repeat for t = 0 to n-k
4.         forward propagate the network over the unfolded network
           for k time steps to compute all h and y.
5.         compute the error as: e = pₜ - yₜ
6.         backpropagate the error across the unfolded network and
           update the weights.

Line 2 is the key one: the hidden state is reset at the start of each sweep, and then the network is unrolled for k steps so that ordinary backpropagation applies to what is now a deep feedforward network with shared weights.

Advantages and shortcomings

9. LSTM and GRU

Long Short-Term Memory (LSTM) networks are RNN trained by Backpropagation Through Time. The structural novelty: they do not have proper neurons, but memory blocks connecting layers. A block is more complex than a normal neuron — it contains "gates" to maintain internal state and output, and each gate works on the input sequence and uses an activation function that determines its contribution to internal state and output. The weights of all gates must be defined during training.

An LSTM in torch

class LSTM(nn.Module):
   def __init__(self, input_size=1, hidden_size=24, num_layers=1):
      super().__init__()
      self.lstm = nn.LSTM(input_size, hidden_size, num_layers, batch_first=True)
      self.fc   = nn.Linear(hidden_size, 1)

   def forward(self, x):
      out, _ = self.lstm(x)
      return self.fc(out[:, -1, :])      # take the LAST time step only

df   = pd.read_csv("BoxJenkins.csv", usecols=[1])
data = df.values.astype(float)

scaler = MinMaxScaler()          # Normalize, needed
data_scaled = scaler.fit_transform(data)

lag = 12
X, y = create_dataset(data_scaled, lag)
X_train = torch.tensor(X, dtype=torch.float32)
y_train = torch.tensor(y, dtype=torch.float32)

model     = LSTM()
loss_fn   = nn.MSELoss()
optimizer = torch.optim.Adam(model.parameters(), lr=0.01)

for epoch in range(100):
   model.train()
   output = model(X_train)
   loss   = loss_fn(output, y_train)
   optimizer.zero_grad(); loss.backward(); optimizer.step()

# Forecast next 12 values, recursively
model.eval()
fore = []
seq = torch.tensor(X[-1], dtype=torch.float32).unsqueeze(0)
for _ in range(12):
   with torch.no_grad():
      pred = model(seq)
      fore.append(pred.item())
      seq = torch.cat([seq[:, 1:], pred.unsqueeze(1)], dim=1)   # roll forward

forecast = scaler.inverse_transform(np.array(fore).reshape(-1, 1))   # inverse scale

Three details in that block carry the whole chapter. MinMaxScaler with the comment needed — section 6's mandatory preprocessing. out[:, -1, :] — the LSTM emits an output at every step, but only the last one is the forecast. And scaler.inverse_transform — the transform is undone before anyone reads the numbers, exactly as in chapters 3 and 4.

10. The fil rouge scoreboard

The deck evaluates a simple MLP model on the same experiment used for SARIMA in chapter 4: train 2004-2007, test 2008.

The slides annotate this table with two words: easy to improve. And they are right — but the honest reading is the comparison, and it is not flattering for the network.

Check your understanding