Skip to content

Time Series Forecasting with RNNs

What a time series is

A time series is a sequence of one or more values recorded over time — active users per hour, daily temperature, quarterly revenue. If there’s a single value per time step it’s univariate; if there are several (revenue, debt, and so on) it’s multivariate. The main task is forecasting: predicting future values. A related task is imputation: filling in missing values from the past.

diagram Diagram mermaid

Generating a toy time series

The book generates a synthetic univariate series — the sum of two sine waves with random frequency and phase, plus a little noise — so you can focus on the model instead of data cleaning:

generate_series.py
import numpy as np
 
def generate_time_series(batch_size, n_steps):
    freq1, freq2, offsets1, offsets2 = np.random.rand(4, batch_size, 1)
    time = np.linspace(0, 1, n_steps)
    series = 0.5 * np.sin((time - offsets1) * (freq1 * 10 + 10))       # wave 1
    series += 0.2 * np.sin((time - offsets2) * (freq2 * 20 + 20))      # + wave 2
    series += 0.1 * (np.random.rand(batch_size, n_steps) - 0.5)        # + noise
    return series[..., np.newaxis].astype(np.float32)
 
n_steps = 50
series = generate_time_series(10000, n_steps + 1)
X_train, y_train = series[:7000, :n_steps], series[:7000, -1]
X_valid, y_valid = series[7000:9000, :n_steps], series[7000:9000, -1]
X_test, y_test = series[9000:, :n_steps], series[9000:, -1]
 
print(X_train.shape, y_train.shape)   # (7000, 50, 1) (7000, 1)
generate_series.py
import numpy as np
 
def generate_time_series(batch_size, n_steps):
    freq1, freq2, offsets1, offsets2 = np.random.rand(4, batch_size, 1)
    time = np.linspace(0, 1, n_steps)
    series = 0.5 * np.sin((time - offsets1) * (freq1 * 10 + 10))       # wave 1
    series += 0.2 * np.sin((time - offsets2) * (freq2 * 20 + 20))      # + wave 2
    series += 0.1 * (np.random.rand(batch_size, n_steps) - 0.5)        # + noise
    return series[..., np.newaxis].astype(np.float32)
 
n_steps = 50
series = generate_time_series(10000, n_steps + 1)
X_train, y_train = series[:7000, :n_steps], series[:7000, -1]
X_valid, y_valid = series[7000:9000, :n_steps], series[7000:9000, -1]
X_test, y_test = series[9000:, :n_steps], series[9000:, -1]
 
print(X_train.shape, y_train.shape)   # (7000, 50, 1) (7000, 1)

Notice the shape: RNN inputs are always [batch_size, time_steps, dimensionality][batch_size, time_steps, dimensionality]dimensionalitydimensionality is 1 for a univariate series, more for a multivariate one.

Baseline metrics first

Before reaching for an RNN, always check a couple of simple baselines — otherwise you might think a fancy model works great when it’s actually worse than doing almost nothing.

Naive forecasting — just predict that the series repeats its last value:

naive_baseline.py
import numpy as np
from tensorflow import keras
 
y_pred_naive = X_valid[:, -1]
mse_naive = np.mean(keras.losses.mean_squared_error(y_valid, y_pred_naive))
print("Naive MSE:", mse_naive)   # ~0.020
naive_baseline.py
import numpy as np
from tensorflow import keras
 
y_pred_naive = X_valid[:, -1]
mse_naive = np.mean(keras.losses.mean_squared_error(y_valid, y_pred_naive))
print("Naive MSE:", mse_naive)   # ~0.020

A plain linear model (a DenseDense layer on the flattened input) does noticeably better:

linear_baseline.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.Flatten(input_shape=[50, 1]),
    keras.layers.Dense(1),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("Linear MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.004
linear_baseline.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.Flatten(input_shape=[50, 1]),
    keras.layers.Dense(1),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("Linear MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.004

That linear model is your target to beat with an RNN.

A simple RNN

simple_rnn_forecast.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(1, input_shape=[None, 1]),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("SimpleRNN MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.014
simple_rnn_forecast.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(1, input_shape=[None, 1]),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("SimpleRNN MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.014

A single-neuron SimpleRNNSimpleRNN beats the naive baseline but doesn’t beat the linear model — it only has 3 parameters total, which just isn’t much capacity.

Deep RNN

Stacking recurrent layers gives the model more capacity — this is called a deep RNN. The key rule: every recurrent layer feeding another recurrent layer needs return_sequences=Truereturn_sequences=True, so it hands over a full 3D sequence instead of just the last step’s output:

deep_rnn.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20),
    keras.layers.Dense(1),   # a Dense output layer (not the last recurrent layer)
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("Deep RNN MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.003
deep_rnn.py
from tensorflow import keras
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20),
    keras.layers.Dense(1),   # a Dense output layer (not the last recurrent layer)
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, y_train, epochs=20,
          validation_data=(X_valid, y_valid), verbose=0)
print("Deep RNN MSE:", model.evaluate(X_valid, y_valid, verbose=0))   # ~0.003

Using a Dense(1)Dense(1) output instead of ending on SimpleRNN(1)SimpleRNN(1) is usually better: it trains a little faster, performs about the same, and lets you pick any output activation — a final recurrent layer with only 1 unit is stuck with a tiny hidden state and (by default) a tanhtanh output squeezed into [-1, 1][-1, 1].

Forecasting several steps ahead

What if you need the next 10 values, not just the next one? There are two common strategies.

Strategy 1 — repeat one-step forecasts. Predict the next value, append it to the input, predict again, and so on:

iterative_forecast.py
import numpy as np
 
series = generate_time_series(1, n_steps + 10)
X_new, Y_new = series[:, :n_steps], series[:, n_steps:]
X = X_new
for step_ahead in range(10):
    y_pred_one = model.predict(X[:, step_ahead:], verbose=0)[:, np.newaxis, :]
    X = np.concatenate([X, y_pred_one], axis=1)
 
Y_pred = X[:, n_steps:]
iterative_forecast.py
import numpy as np
 
series = generate_time_series(1, n_steps + 10)
X_new, Y_new = series[:, :n_steps], series[:, n_steps:]
X = X_new
for step_ahead in range(10):
    y_pred_one = model.predict(X[:, step_ahead:], verbose=0)[:, np.newaxis, :]
    X = np.concatenate([X, y_pred_one], axis=1)
 
Y_pred = X[:, n_steps:]

Errors accumulate with each extra step, so this works best when you only need a handful of steps ahead.

Strategy 2 — predict all 10 steps at once. Change the targets to be 10-value vectors and give the output layer 10 units:

multi_step_forecast.py
from tensorflow import keras
 
series = generate_time_series(10000, n_steps + 10)
X_train, Y_train = series[:7000, :n_steps], series[:7000, -10:, 0]
X_valid, Y_valid = series[7000:9000, :n_steps], series[7000:9000, -10:, 0]
X_test, Y_test = series[9000:, :n_steps], series[9000:, -10:, 0]
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20),
    keras.layers.Dense(10),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, Y_train, epochs=20,
          validation_data=(X_valid, Y_valid), verbose=0)
multi_step_forecast.py
from tensorflow import keras
 
series = generate_time_series(10000, n_steps + 10)
X_train, Y_train = series[:7000, :n_steps], series[:7000, -10:, 0]
X_valid, Y_valid = series[7000:9000, :n_steps], series[7000:9000, -10:, 0]
X_test, Y_test = series[9000:, :n_steps], series[9000:, -10:, 0]
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20),
    keras.layers.Dense(10),
])
model.compile(loss="mse", optimizer="adam")
model.fit(X_train, Y_train, epochs=20,
          validation_data=(X_valid, Y_valid), verbose=0)

This is more accurate than the iterative approach, but it only trains on the error at the last time step. Better still: turn it into a true sequence-to-sequence model, forecasting the next 10 values at every time step (not only the last one) with return_sequences=Truereturn_sequences=True on every recurrent layer and a TimeDistributed(Dense(10))TimeDistributed(Dense(10)) output. This feeds far more error gradients through the network — and stabilizes and speeds up training:

seq_to_seq_forecast.py
import numpy as np
from tensorflow import keras
 
Y = np.empty((10000, n_steps, 10))
for step_ahead in range(1, 11):
    Y[:, :, step_ahead - 1] = series[:, step_ahead:step_ahead + n_steps, 0]
Y_train, Y_valid, Y_test = Y[:7000], Y[7000:9000], Y[9000:]
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20, return_sequences=True),
    keras.layers.TimeDistributed(keras.layers.Dense(10)),
])
 
def last_time_step_mse(Y_true, Y_pred):
    return keras.metrics.mean_squared_error(Y_true[:, -1], Y_pred[:, -1])
 
model.compile(loss="mse", optimizer=keras.optimizers.Adam(learning_rate=0.01),
              metrics=[last_time_step_mse])
model.fit(X_train, Y_train, epochs=20, validation_data=(X_valid, Y_valid), verbose=0)
seq_to_seq_forecast.py
import numpy as np
from tensorflow import keras
 
Y = np.empty((10000, n_steps, 10))
for step_ahead in range(1, 11):
    Y[:, :, step_ahead - 1] = series[:, step_ahead:step_ahead + n_steps, 0]
Y_train, Y_valid, Y_test = Y[:7000], Y[7000:9000], Y[9000:]
 
model = keras.models.Sequential([
    keras.layers.SimpleRNN(20, return_sequences=True, input_shape=[None, 1]),
    keras.layers.SimpleRNN(20, return_sequences=True),
    keras.layers.TimeDistributed(keras.layers.Dense(10)),
])
 
def last_time_step_mse(Y_true, Y_pred):
    return keras.metrics.mean_squared_error(Y_true[:, -1], Y_pred[:, -1])
 
model.compile(loss="mse", optimizer=keras.optimizers.Adam(learning_rate=0.01),
              metrics=[last_time_step_mse])
model.fit(X_train, Y_train, epochs=20, validation_data=(X_valid, Y_valid), verbose=0)

All outputs are used for the training loss, but only the last time step matters for evaluation — hence the custom last_time_step_mselast_time_step_mse metric. TimeDistributedTimeDistributed applies its wrapped layer (here, Dense(10)Dense(10)) independently to every time step.

Comparing the approaches

ModelApprox. MSE (1-step)Notes
Naive (repeat last value)0.020free, always compute this first
Linear (FlattenFlatten + DenseDense)0.004strong, cheap baseline
SimpleRNN(1)SimpleRNN(1)0.014worse than linear — too little capacity
Deep RNN (3 SimpleRNNSimpleRNN layers)0.003finally beats the linear model

For 10-step-ahead forecasts, sequence-to-vector (~0.008 MSE) beats the iterative one-step-at-a-time approach (~0.029 MSE), and a full sequence-to-sequence model (~0.006 MSE) beats both.

A real-world example: the Jena temperature dataset (Chollet)

The synthetic sine-wave series above is great for learning the mechanics, but real sensor data is messier. Deep Learning with Python (Chollet) walks through a concrete problem: given hourly weather measurements — temperature, pressure, humidity, and 11 other quantities — from the past 5 days, predict the temperature 24 hours from now. The raw file has 420,551 rows (one every 10 minutes, from 2009-2016) and 14 columns.

diagram Diagram mermaid

Windowing data with timeseries_dataset_from_arraytimeseries_dataset_from_array

Instead of writing your own generator for overlapping windows, Keras ships a utility that slices a raw array into fixed-length sequences (and lines up a matching target for each one) on the fly, without duplicating the data in memory. A tiny example makes the idea concrete — feed in the numbers 0-9, ask for windows of length 3, and target each window with the value 3 steps ahead:

windowing_demo.py
import numpy as np
from tensorflow import keras
 
int_sequence = np.arange(10)
dummy_dataset = keras.utils.timeseries_dataset_from_array(
    data=int_sequence[:-3],
    targets=int_sequence[3:],
    sequence_length=3,
    batch_size=2,
)
for inputs, targets in dummy_dataset:
    for i in range(inputs.shape[0]):
        print([int(x) for x in inputs[i]], int(targets[i]))
 
# ── Output ───────────────────────────
# [0, 1, 2] 3
# [1, 2, 3] 4
# [2, 3, 4] 5
# [3, 4, 5] 6
# [4, 5, 6] 7
# ─────────────────────────────────────
windowing_demo.py
import numpy as np
from tensorflow import keras
 
int_sequence = np.arange(10)
dummy_dataset = keras.utils.timeseries_dataset_from_array(
    data=int_sequence[:-3],
    targets=int_sequence[3:],
    sequence_length=3,
    batch_size=2,
)
for inputs, targets in dummy_dataset:
    for i in range(inputs.shape[0]):
        print([int(x) for x in inputs[i]], int(targets[i]))
 
# ── Output ───────────────────────────
# [0, 1, 2] 3
# [1, 2, 3] 4
# [2, 3, 4] 5
# [3, 4, 5] 6
# [4, 5, 6] 7
# ─────────────────────────────────────

For the real temperature task, the book uses sampling_rate=6sampling_rate=6 (keep one point per hour instead of every 10 minutes), sequence_length=120sequence_length=120 (5 days of hourly data), and delay = sampling_rate * (sequence_length + 24 - 1)delay = sampling_rate * (sequence_length + 24 - 1) (so each target lands exactly 24 hours past the end of its window):

jena_windows.py
from tensorflow import keras
 
sampling_rate = 6
sequence_length = 120
delay = sampling_rate * (sequence_length + 24 - 1)
batch_size = 256
 
train_dataset = keras.utils.timeseries_dataset_from_array(
    raw_data[:-delay],
    targets=temperature[delay:],
    sampling_rate=sampling_rate,
    sequence_length=sequence_length,
    shuffle=True,
    batch_size=batch_size,
    start_index=0,
    end_index=num_train_samples,
)
# val_dataset / test_dataset follow the same call with shifted start/end_index
 
for samples, targets in train_dataset:
    print("samples shape:", samples.shape)   # (256, 120, 14)
    print("targets shape:", targets.shape)   # (256,)
    break
jena_windows.py
from tensorflow import keras
 
sampling_rate = 6
sequence_length = 120
delay = sampling_rate * (sequence_length + 24 - 1)
batch_size = 256
 
train_dataset = keras.utils.timeseries_dataset_from_array(
    raw_data[:-delay],
    targets=temperature[delay:],
    sampling_rate=sampling_rate,
    sequence_length=sequence_length,
    shuffle=True,
    batch_size=batch_size,
    start_index=0,
    end_index=num_train_samples,
)
# val_dataset / test_dataset follow the same call with shifted start/end_index
 
for samples, targets in train_dataset:
    print("samples shape:", samples.shape)   # (256, 120, 14)
    print("targets shape:", targets.shape)   # (256,)
    break

A common-sense baseline (before touching a neural net)

Chollet’s rule: always write down the dumbest possible heuristic and measure it before reaching for deep learning. Here, temperature is both continuous and roughly periodic day-to-day, so a reasonable guess is “24 hours from now, the temperature will be about the same as it is right now”:

common_sense_baseline.py
import numpy as np
 
def evaluate_naive_method(dataset, temperature_mean, temperature_std):
    total_abs_err = 0.0
    samples_seen = 0
    for samples, targets in dataset:
        # column 1 is temperature; un-normalize it back to degrees Celsius
        preds = samples[:, -1, 1] * temperature_std + temperature_mean
        total_abs_err += np.sum(np.abs(preds - targets))
        samples_seen += samples.shape[0]
    return total_abs_err / samples_seen
 
# Validation MAE: 2.44 degrees Celsius
# Test MAE:       2.62 degrees Celsius
common_sense_baseline.py
import numpy as np
 
def evaluate_naive_method(dataset, temperature_mean, temperature_std):
    total_abs_err = 0.0
    samples_seen = 0
    for samples, targets in dataset:
        # column 1 is temperature; un-normalize it back to degrees Celsius
        preds = samples[:, -1, 1] * temperature_std + temperature_mean
        total_abs_err += np.sum(np.abs(preds - targets))
        samples_seen += samples.shape[0]
    return total_abs_err / samples_seen
 
# Validation MAE: 2.44 degrees Celsius
# Test MAE:       2.62 degrees Celsius

This uses mean absolute error (MAE) rather than MSE, since MAE stays in the original units (degrees) and is easy to reason about: “off by 2.44 degrees on average.”

Dense and Conv1D models don’t beat the baseline

It’s tempting to assume “more complex = better,” so it’s worth checking a plain DenseDense model and a Conv1DConv1D model against that 2.44-degree baseline first:

dense_vs_conv_jena.py
from tensorflow import keras
from tensorflow.keras import layers
 
# A flattened, densely-connected model
inputs = keras.Input(shape=(sequence_length, 14))
x = layers.Flatten()(inputs)
x = layers.Dense(16, activation="relu")(x)
outputs = layers.Dense(1)(x)
dense_model = keras.Model(inputs, outputs)
 
# A 1D convnet
inputs = keras.Input(shape=(sequence_length, 14))
x = layers.Conv1D(8, 24, activation="relu")(inputs)
x = layers.MaxPooling1D(2)(x)
x = layers.Conv1D(8, 12, activation="relu")(x)
x = layers.MaxPooling1D(2)(x)
x = layers.Conv1D(8, 6, activation="relu")(x)
x = layers.GlobalAveragePooling1D()(x)
outputs = layers.Dense(1)(x)
conv_model = keras.Model(inputs, outputs)
dense_vs_conv_jena.py
from tensorflow import keras
from tensorflow.keras import layers
 
# A flattened, densely-connected model
inputs = keras.Input(shape=(sequence_length, 14))
x = layers.Flatten()(inputs)
x = layers.Dense(16, activation="relu")(x)
outputs = layers.Dense(1)(x)
dense_model = keras.Model(inputs, outputs)
 
# A 1D convnet
inputs = keras.Input(shape=(sequence_length, 14))
x = layers.Conv1D(8, 24, activation="relu")(inputs)
x = layers.MaxPooling1D(2)(x)
x = layers.Conv1D(8, 12, activation="relu")(x)
x = layers.MaxPooling1D(2)(x)
x = layers.Conv1D(8, 6, activation="relu")(x)
x = layers.GlobalAveragePooling1D()(x)
outputs = layers.Dense(1)(x)
conv_model = keras.Model(inputs, outputs)
ModelVal MAEvs. baseline (2.44)
Common-sense baseline2.44
FlattenFlatten + DenseDense~2.5-3.0rarely beats it — flattening throws away the notion of time
Conv1DConv1D + pooling~2.9worse — pooling destroys order, and weather isn’t shift-invariant across a day
LSTM(16)LSTM(16)2.36finally beats the baseline

Flattening the window into one long vector destroys any notion of “before” and “after.” A Conv1DConv1D at least respects local order within its window, but pooling still throws away where in the day a pattern happened — morning weather behaves differently than evening weather, so the “same pattern anywhere in the window” assumption behind convolutions doesn’t hold as well as it does for images. An RNN, which walks through the sequence in order and never destroys that order, is the natural fit — exactly the deep RNNs built earlier on this page.

Watch a window slide across the raw series

Every sample an RNN sees is a fixed-length window cut from the raw series. As the window slides forward one step at a time, the point it’s predicting (24 hours past the window’s end) slides forward with it:

sketch Sliding window over a temperature series p5.js
The amber box is one training window (past values), sliding across the series; the connected blue dot is its 24-hour-ahead target.

Mini-checkpoint

Why must you set return_sequences=Truereturn_sequences=True on every recurrent layer except possibly the last one, when stacking recurrent layers?

  • because each layer needs to receive a full 3D sequence [batch, time, features][batch, time, features] as input, not just the 2D output of the previous layer’s last time step.

Visualize it

Watch a noisy signal (blue, the actual values) and a forecast curve (amber, what the model predicts) — the forecast tracks the true signal for the steps it has already seen, then keeps extrapolating for the future steps marked past the dashed line:

sketch Time series forecast p5.js
The blue curve is the true signal, traced out over time; the amber curve is the model's forecast, which stays close for known steps and drifts a little past the forecast boundary.

🧪 Try It Yourself

Exercise 1 – Compute a Naive Baseline

Exercise 2 – Stack a Deep RNN

Exercise 3 – Prepare Multi-Step Targets

Exercise 4 – Window a Series with timeseries_dataset_from_array

Next

Continue to Advanced Recurrent Layers (Dropout, Stacking, Bidirectional) — refine these RNNs with recurrent dropout, deeper stacks, and bidirectional wrappers before moving on to text data in Phase 5.

If this helped you, consider buying me a coffee ☕

Buy me a coffee

Was this page helpful?

Let us know how we did