Time-Series Demand Forecasting for Smart Retail

MEng Capstone · UC Berkeley · 2025

A systematic benchmark of deep learning and classical ML architectures for retail demand forecasting, comparing six models on Walmart's M5 Kaggle dataset (3,049 products × 10 stores, 42,840 hierarchical time series in total). The Encoder-Decoder Transformer achieves the best MAE of 122.01 while simpler models remain surprisingly competitive.

Every component is built from scratch in PyTorch: the Transformer encoder-decoder with cross-attention and teacher forcing, the sliding-window dataset with temporal train/val/test splits, and a comprehensive feature engineering pipeline covering lag features, rolling statistics, calendar events, and price momentum. An XGBoost/LightGBM gradient boosting baseline provides a strong classical ML comparison point.

Model Comparison

MAE comparison across all six models — lower is better. Encoder-Decoder Transformer achieves best result.

Training and validation MAE over 100 epochs. Encoder-Decoder shows smoothest convergence.

28-day forecast on a sample product. Encoder-Decoder (green) tracks actual sales most closely.

Encoder-Decoder Architecture

Input Sequence (28 days)                    Target Sequence (28 days)
        │                                           │
        v                                           v
TimeSeriesEmbedding                    Decoder Embedding
  Linear + PosEnc                        Linear + PosEnc
        │                                           │
        v                                           v
Encoder (N layers)                     Decoder (N layers)
  Self-Attention                         Cross-Attention
  + Feed-Forward   →→ Memory →→      + Causal Self-Attn
  + LayerNorm                            + Feed-Forward
        │                                           │
   Encoder Output                              v
                                       Linear → Forecast (B, 28)

Feature Engineering Pipeline

M5 Raw Data
  │
  └─ sales_train_evaluation.csv → Sales Matrix (42,840 × 1,969)
  │
  └─ calendar.csv → Calendar: day_of_week, month, is_weekend, has_event, SNAP
  │
  ┌─ sell_prices.csv → Price: normalized_price, price_momentum (Δ7d)

Derived Features
  Sales Matrix → Lag features (7d, 14d, 28d)
               → Rolling stats (7d/14d/28d mean, std)

Final Tensor: (42,840 × 1,969 × 15) → Sliding Window → Dataset

Key Features

• Six model architectures — Linear baseline, Encoder-Only Transformer, Encoder-Decoder Transformer, DeepAR, Temporal Fusion Transformer, and XGBoost/LightGBM gradient boosting
• Comprehensive feature engineering — Lag features (7/14/28d), rolling statistics (mean/std over 7/14/28d windows), calendar events & SNAP indicators, and price momentum
• Temporal train/val/test split — Data split along the time axis (80/10/10) to prevent future data leakage, critical for time series validation
• Teacher forcing & autoregressive decoding — Encoder-Decoder uses teacher forcing during training (shifted target input) and step-by-step generation during inference
• M5 competition metrics — MAE, RMSE, RMSSE, and Weighted RMSSE with per-series and aggregate computation matching the official Kaggle evaluation
• YAML-based experiment configs — Frozen dataclass configuration system loads from YAML files, enabling reproducible experiments with different architectures and hyperparameters
• Direct multi-step gradient boosting — One model per forecast step converts time series into tabular regression with lag/rolling/calendar features as columns

Code Highlights

class PositionalEncoding(nn.Module):
    def __init__(self, d_model, max_len=5000, dropout=0.1):
        super().__init__()
        self.dropout = nn.Dropout(p=dropout)
        pe = torch.zeros(max_len, d_model)
        position = torch.arange(0, max_len).unsqueeze(1).float()
        div_term = torch.exp(
            torch.arange(0, d_model, 2).float() * (-math.log(10000.0) / d_model)
        )
        pe[:, 0::2] = torch.sin(position * div_term)
        pe[:, 1::2] = torch.cos(position * div_term)
        self.register_buffer("pe", pe.unsqueeze(0))

    def forward(self, x):
        return self.dropout(x + self.pe[:, :x.size(1)])

def forward(self, x_sales, x_features=None, y=None):
    memory = self._encode(x_sales, x_features)

    if y is not None:  # Teacher forcing
        decoder_input = torch.cat(
            [x_sales[:, -1:], y[:, :-1]], dim=1
        ).unsqueeze(-1)  # shift right, prepend last known
        tgt = self.decoder_embedding(decoder_input)
        tgt_mask = generate_causal_mask(tgt.size(1), tgt.device)
        decoded = self.decoder(tgt, memory, tgt_mask=tgt_mask)
        return self.output_layer(decoded).squeeze(-1)
    else:  # Autoregressive inference
        return self._autoregressive_decode(memory, x_sales[:, -1])

def build_rolling_features(self, sales_matrix):
    """Compute rolling mean & std for each window size."""
    n_series, n_days = sales_matrix.shape
    n_feats = len(self.ROLLING_WINDOWS) * 2  # mean + std
    rolling = np.zeros((n_series, n_days, n_feats), dtype=np.float32)

    for w_idx, window in enumerate(self.ROLLING_WINDOWS):
        for t in range(window, n_days):
            segment = sales_matrix[:, t - window : t]
            rolling[:, t, w_idx * 2]     = segment.mean(axis=1)
            rolling[:, t, w_idx * 2 + 1] = segment.std(axis=1)

    return rolling  # (n_series, n_days, n_windows * 2)

Design Decisions

• Encoder-Decoder over Encoder-Only — The decoder's cross-attention to encoded history produces 11% lower MAE than the encoder-only projection head, at only 1.7× the training cost. The autoregressive decoding captures step-to-step forecast dependencies
• Pre-norm Transformer architecture — LayerNorm before attention/FFN (norm_first=True) provides more stable training gradients than post-norm, critical for convergence on noisy retail time series
• Teacher forcing with shifted targets — During training, the decoder receives the ground truth shifted right by one step, accelerating convergence vs. pure autoregressive training while maintaining generation quality at inference
• Temporal split over random split — Splitting along the time axis (80/10/10) prevents future information leakage that would inflate validation metrics. A random split would leak future sales patterns into the training set
• Direct multi-step for gradient boosting — One model per forecast horizon step avoids error accumulation from recursive single-step prediction, at the cost of training 28 independent models
• Frozen dataclass configs — Immutable configuration objects prevent accidental mutation during training and enable reproducible experiment tracking via YAML serialization

How It Works

Data Pipeline: The M5 dataset contains 3,049 products × 10 stores, 42,840 hierarchical time series in total, each spanning ~1,969 days. Three CSV files (sales, calendar, prices) are merged into a unified feature tensor. The FeatureEngineer class extracts calendar features (weekday, month, events, SNAP), computes lag features at 7/14/28-day offsets, and generates rolling statistics (mean/std) over multiple window sizes. Price features include normalized sell price and 7-day price momentum.

Feature Engineering: Each time series sample consists of a 28-day input window and 28-day forecast horizon. The sliding window M5Dataset generates all valid (input, target) pairs, indexing into a precomputed sales matrix and feature tensor. Lag features are zero-padded for the first lag days; rolling statistics use vectorized numpy operations across all 42,840 series simultaneously for efficiency.

Model Training: The Trainer class manages the training loop with Adam optimizer, learning rate scheduling (Cosine Annealing or ReduceLROnPlateau), gradient clipping at norm 1.0, and best-model checkpointing based on validation MAE. The Encoder-Decoder Transformer uses teacher forcing during training — the decoder receives the shifted ground truth — and switches to autoregressive step-by-step generation at inference time.

Evaluation: Models are evaluated using MAE (primary metric), RMSE, RMSSE (M5 competition metric that normalizes by naive forecast variance), and Weighted RMSSE. The RMSSE normalizes per-series errors by the variance of a naive random-walk forecast on the training history, making errors comparable across products with vastly different sales volumes.

Challenges & Solutions

• Intermittent demand (many zeros) — Many products have long stretches of zero sales. Standard MSE loss would optimize for predicting zero. Using MAE as both loss and metric prevents over-penalizing rare spikes, and the RMSSE metric normalizes per-series to avoid domination by high-volume items
• Complex models underperformed — DeepAR (MAE 144.98) and TFT (MAE 135.07) performed worse than the Linear baseline despite higher model capacity. Root cause: insufficient hyperparameter tuning within the project timeline. These architectures require extensive grid search across learning rates, hidden dimensions, and attention heads
• Temporal data leakage risk — Naive random train/test splits would leak future patterns into training. Solved with strict temporal splitting: train on days 1-1575, validate on 1575-1772, test on 1772-1969. The M5Dataset sliding window enforces this boundary
• 42,840-series scale — Processing all series simultaneously for rolling statistics and lag features requires careful memory management. Vectorized numpy operations across the full matrix (42,840 × 1,969) keep preprocessing under 30 seconds on standard hardware

Tech Stack

References

• Vaswani, A. et al. (2017). “Attention Is All You Need.” NeurIPS. arXiv:1706.03762
• Zeng, A. et al. (2023). “Are Transformers Effective for Time Series Forecasting?” AAAI. arXiv:2205.13504
• Makridakis, S. et al. (2022). “M5 accuracy competition: Results, findings, and conclusions.” International Journal of Forecasting.