Time series

Problem

Sales, server load, energy demand, sensor readings and prices all arrive as one value per time step, and the question is almost always the same: what will the next few values be? Two things make this different from the rows-and-columns problems on the rest of this site. The rows are not independent, since neighbouring values are similar, and the future is not drawn from the same shuffled pool as the past. A model can only use what was known at the moment of the forecast, and it is judged on a stretch of time it never saw.

Ignoring either fact produces models that look excellent in a notebook and fail in production. The standard cross-validation recipe shuffles rows, which lets the model peek at the days on either side of each test day. Features computed over the whole series, such as a mean or a scaler, carry information backwards from the future, which is a form of data leakage.

Intuition

Most business and physical series are a few simple pieces added together: a slow trend, a repeating seasonal pattern with a known period (12 for monthly data with a yearly cycle, 7 for daily data with a weekly one) and a remainder that is noisy but not independent, because a shock this month fades over the next few. Forecasting is mostly about finding those pieces and carrying them forward. A model that captures the trend and the season leaves only the remainder, and the best it can do there is use how quickly shocks fade.

The validation rule follows from how the model will be used: train on everything up to a point in time, test on what comes after, and repeat for several cut-off points. That is all TimeSeriesSplit does.

Mechanics

Classical decomposition

The additive model writes each value as y_t = T_t + S_t + R_t, trend plus season plus remainder. With period m, classical decomposition estimates them in three steps:

1. Trend. A centred moving average over one full period cancels the season. For even m the window has m + 1 points with the two ends weighted ½, the "2×m-MA": T̂_t = (½y_t−6 + y_t−5 + … + y_t+5 + ½y_t+6) / 12 It is undefined for the first and last m/2 points.
2. Season. Average the detrended values y_t − T̂_t separately for each position in the cycle (every January, every February, …), then shift the m averages so they sum to zero.
3. Remainder. R̂_t = y_t − T̂_t − Ŝ_t.

If the seasonal swing grows with the level, the multiplicative model y = T × S × R applies instead, or the additive one on log y. STL (Cleveland et al., 1990) replaces the moving averages with loess smoothers, which handles the ends of the series and lets the season change slowly over time.

Autocorrelation

The sample autocorrelation at lag k measures how much the series resembles itself shifted by k steps:

r_k = Σ_t=1..n−k (y_t − ȳ)(y_t+k − ȳ) / Σ_t=1..n (y_t − ȳ)²

For white noise, about 95% of the r_k fall inside ±1.96/√n. A slowly decaying ACF means a trend; a bump at lag m means a season; a remainder whose ACF is still outside the band at small lags means there is short-term structure a model can still use.

Turning a series into a table

Most machine learning models need rows of features and a target. A forecast at time t for horizon h may use only values up to t, so the usual features are:

• Lags: y_t, y_t−1, …, and the seasonal lag y_t+h−m, whenever h ≤ m.
• Rolling statistics over past windows (mean, max, standard deviation), computed so that the window ends at t, never after it.
• Calendar features: a time index for the trend and the month, weekday or holiday flags for the season. Fourier terms sin(2πkt/m), cos(2πkt/m) for k = 1…K encode a smooth season with 2K columns instead of m − 1 dummies.
• Known future inputs such as planned promotions or the weather forecast, which are available at forecast time.

For several steps ahead there are two main strategies. A recursive model predicts one step, feeds that prediction back in as a lag, and repeats, so its errors compound. A direct model is trained separately for each horizon h, using only features known at t. The "trend + lag 12" model below is direct for every h ≤ 12, because y_t−12 is always in the past.

Classical models

An autoregressive model AR(p) predicts y_t = c + φ₁y_t−1 + … + φ_py_t−p + ε_t. A moving-average model MA(q) uses past noise terms instead: y_t = c + ε_t + θ₁ε_t−1 + … + θ_qε_t−q. Both assume the series is stationary, meaning its mean and autocorrelations do not drift over time. ARIMA(p, d, q) differences the series d times (y_t − y_t−1) to remove a trend first, and seasonal ARIMA adds seasonal differences y_t − y_t−m and seasonal AR and MA terms. Box and Jenkins set out the iterative recipe still used for these models: identify p, d and q from the autocorrelations, estimate the coefficients, then check that the residuals look like noise. Exponential smoothing (simple, Holt's trend, Holt-Winters seasonal) is the other classical family: forecasts are weighted averages of the past with weights that decay geometrically. Neither is in PyTorch or scikit-learn; they live in statistics packages.

Time-ordered validation

TimeSeriesSplit(n_splits, test_size, gap, max_train_size) places n_splits test blocks of test_size consecutive samples at the end of the data. The first starts at n − n_splits × test_size, and each training set is everything before its test block. gap drops that many samples between training and test, for features whose windows would otherwise straddle the boundary. max_train_size turns the growing window into a sliding one. Unlike KFold, there is no shuffle option. The class assumes the samples are equally spaced and already in time order.

The scikit-learn docstring's own example, 12 samples with n_splits=3 and test_size=2, gives training sets ending at 5, 7, 9 and test blocks 6–7, 8–9, 10–11. With gap=2 the training sets shrink to end at 3, 5, 7. This page's time_series_split helper reproduces both and is tested against them.

Worked example

Implementation

import pandas as pd
from sklearn.ensemble import HistGradientBoostingRegressor
from sklearn.model_selection import TimeSeriesSplit, cross_val_score

df = pd.DataFrame({"y": y}, index=pd.period_range("2017-01", periods=len(y), freq="M"))
for lag in (1, 2, 12):
    df[f"lag_{lag}"] = df["y"].shift(lag)            # only past values
df["roll_mean_3"] = df["y"].shift(1).rolling(3).mean()   # window ends at t - 1
df["month"] = df.index.month
df = df.dropna()

X, target = df.drop(columns="y"), df["y"]
cv = TimeSeriesSplit(n_splits=4, test_size=12, gap=0)
scores = cross_val_score(HistGradientBoostingRegressor(), X, target, cv=cv,
                         scoring="neg_mean_absolute_error")
# lag_1 makes this a one-step-ahead model. For a 12-month horizon use only lags >= 12,
# or forecast recursively and feed predictions back into lag_1.


import torch
from torch import nn

class WindowForecaster(nn.Module):
    def __init__(self, hidden=32, horizon=12):
        super().__init__()
        self.lstm = nn.LSTM(input_size=1, hidden_size=hidden, batch_first=True)
        self.head = nn.Linear(hidden, horizon)       # direct: all horizons at once

    def forward(self, window):                       # window: (N, L, 1), scaled
        output, (h_n, c_n) = self.lstm(window)       # output: (N, L, hidden)
        return self.head(output[:, -1])              # (N, horizon)

The pandas shift calls are where leakage is usually introduced or prevented: every feature for row t must be computable at t. The scaler for the LSTM's inputs, like any preprocessing, is fitted on the training fold only. Tree ensembles such as gradient boosting cannot extrapolate a trend beyond the range of targets they were trained on, so a common pattern is to difference or detrend first and model what remains. An LSTM reads a fixed window of L past steps; with batch_first=True the batch of windows is (N, L, features).

Tradeoffs

• Simple baselines first. Naive and seasonal naive forecasts cost nothing and are hard to beat on short or noisy series. MASE against a naive benchmark makes "better than trivial" a number.
• Classical versus learned. ARIMA and exponential smoothing fit one series with a handful of parameters and give prediction intervals. Gradient boosting and neural models need many related series or long histories to pay off, and they can share what they learn across series.
• Growing versus sliding windows. A growing window (the default) uses all history. max_train_size keeps only recent data, which helps when the process changes over time and hurts when it does not.
• The gap. If a feature uses a window of w past steps, or the target is aggregated over a future window, set gap so no training row overlaps the test period.
• Few, correlated folds. Later folds share most of their training data and test on adjacent periods, so fold scores are not independent and their spread understates the true uncertainty.
• Is K-fold ever fine? Bergmeir and Benítez (2012) compared validation schemes on synthetic and real series and recommended blocked cross-validation (contiguous folds, no shuffle) as the standard. Shuffled folds remain wrong for models like the nearest-month one above, which can read neighbours of the test points.
• Metrics. MAE is in the series' units, RMSE punishes large misses, and percentage errors such as MAPE explode near zero. MASE (Hyndman and Koehler, 2006) is scale-free and comparable across series.

Connections

References

• Box and Jenkins. Time Series Analysis: Forecasting and Control. Holden-Day, 1970.
• Cleveland, Cleveland, McRae and Terpenning. STL: A Seasonal-Trend Decomposition Procedure Based on Loess. Journal of Official Statistics 6, 3–73, 1990.
• Hyndman and Koehler. Another look at measures of forecast accuracy. International Journal of Forecasting 22(4), 679–688, 2006. doi:10.1016/j.ijforecast.2006.03.001
• Bergmeir and Benítez. On the use of cross-validation for time series predictor evaluation. Information Sciences, 2012. doi:10.1016/j.ins.2011.12.028
• Hyndman and Athanasopoulos. Forecasting: Principles and Practice, 3rd edition. OTexts, 2021. otexts.com/fpp3
• scikit-learn 1.9.1: sklearn/model_selection/_split.py. PyTorch v2.15.0: torch/nn/modules/rnn.py.