Module 7 — Calendar features and external regressors

The M4 and M5 competitions ended with the same conclusion: on messy real-world data with many related series, gradient boosting on a carefully engineered supervised table wins. This module turns our pharmacy series into that table, discusses the one design decision that trips up every team — what to do with regressors that are unknown at inference time — and adds a strong row to our results.

From a series to a supervised table

The trick is to reshape the daily observations into rows where each row represents a single day and its features are computed only from the past. For a horizon $h = 28$, we predict $y_{t+h}$ from information known at time $t$.

import pandas as pd

df = pd.read_parquet("pharmacy_daily.parquet")
df = df.set_index("date").sort_index()

def make_features(df: pd.DataFrame, horizon: int = 28) -> pd.DataFrame:
    out = df.copy()
    # Lag features
    for lag in [1, 7, 14, 28, 365]:
        out[f"lag_{lag}"] = out["units"].shift(lag + horizon - 1)
    # Rolling means computed only from the past
    for window in [7, 28]:
        out[f"rmean_{window}"] = (
            out["units"].shift(horizon).rolling(window).mean()
        )
    # Calendar features (safe: known for any date)
    out["dow"] = out.index.dayofweek
    out["dom"] = out.index.day
    out["month"] = out.index.month
    out["is_weekend"] = out["dow"].isin([5, 6]).astype(int)
    return out.dropna()

feat = make_features(df, horizon=28)

Two design choices deserve a comment.

The shift(lag + horizon - 1) pattern. To predict day $t + 28$, we can only use information available at day $t$. So lag_1 at row $t + 28$ must be $y_{t}$, which means shifting by 28 rows. lag_7 uses $y_{t - 6}$ (still available at day $t$), so we shift by $6 + 28 = 34$. Simplifying: to predict $t + h$, feature $\text{lag}_k$ must equal $y_{t + h - k - (h - 1)} = y_{t - k + 1}$, which is shift(h + k - 1). Any different offset leaks the horizon into training.

The rolling mean shift. df["units"].shift(28).rolling(7).mean() computes a 7-day mean ending 28 days ago — safe for a 28-day forecast. Rolling without an explicit shift silently uses the current day and produces the training-time cheating we warned about in module 1.

Cyclic encoding of calendar features

Day-of-week 6 and day-of-week 0 are one day apart, but a tree does not know that. Encoding as sine and cosine on a period of 7 makes the cyclic structure explicit and often helps linear models more than trees, but is still cheap for trees:

import numpy as np
feat["dow_sin"] = np.sin(2 * np.pi * feat["dow"] / 7)
feat["dow_cos"] = np.cos(2 * np.pi * feat["dow"] / 7)
feat["month_sin"] = np.sin(2 * np.pi * feat["month"] / 12)
feat["month_cos"] = np.cos(2 * np.pi * feat["month"] / 12)

For trees, the raw dow integer is usually enough; the sine/cosine pair matters when you switch to linear models or neural networks.

Holidays and promotions

Static holiday flags are easy: pre-compute the calendar for the whole period and merge on date.

from holidays import France

fr = France(years=range(2022, 2027))
feat["is_holiday"] = feat.index.isin(fr).astype(int)
feat["days_to_holiday"] = ...

days_to_holiday is often more useful than the flag alone, because behavior changes in the two days before and the day after (people stock up on the eve of Toussaint, then shop less the next morning).

Promotions are trickier because they are partly known in the future: the team plans them one month ahead. We treat them as an exogenous regressor and pass their planned future values to the model at prediction time.

Fitting a gradient-boosting model

LightGBM is our default here for speed on many series; XGBoost and CatBoost work identically.

import lightgbm as lgb
from sklearn.metrics import mean_absolute_error

feat = feat.dropna()
cutoff = feat.index[-28]
train = feat[feat.index < cutoff]
test = feat[feat.index >= cutoff]

feature_cols = [c for c in feat.columns if c not in ("units",)]

model = lgb.LGBMRegressor(
    n_estimators=1000, learning_rate=0.02, num_leaves=31,
    min_child_samples=20, subsample=0.8, colsample_bytree=0.8,
)
model.fit(
    train[feature_cols], train["units"],
    eval_set=[(test[feature_cols], test["units"])],
    callbacks=[lgb.early_stopping(50)],
)

pred = model.predict(test[feature_cols])
print("LightGBM MAE:", mean_absolute_error(test["units"], pred))

On our data, LightGBM lands at MAE 6.9 units per day — the best result yet.

The critical distinction: regressors known in the future

A regressor helps the model only if we can pass its value at inference time for the horizon we forecast. That splits regressors into two camps.

Known in the future. Day-of-week, holiday flag, planned promotion, planned school closures, planned marketing campaigns. Their future values are on a calendar; we can pass them to the model when producing the 28-day forecast. These are the most valuable features because they encode leading information the past cannot express.

Unknown in the future. Same-day weather, competitor price, real-time footfall, current-day promotions decided by the store manager. Their future values must themselves be forecast, and the resulting error compounds. The safe approach is to not use them as raw features but to use their past (weather one week ago) or to build an intermediate model that predicts them first.

A common mistake: train a gradient boosting model on weather_temperature_today because it improves training-set MAE by 10 %. Then at inference we do not know tomorrow's temperature at inference minute, and we plug in the last observed value — the model was trained to expect near-perfect weather, and its production predictions collapse.

Recursive versus direct multi-step forecasting

We have a 28-day horizon. Two strategies to fill it.

Recursive. Train one model to predict $y_{t+1}$ from features known at $t$. At inference, predict day 1, feed that prediction as lag_1 for day 2, and repeat. Simple, but errors compound.

Direct. Train 28 separate models, each predicting $y_{t+h}$ for one $h \in \{1, \dots, 28\}$, each using only features known at $t$. Slower to train, no error compounding, and each per-horizon model can specialize.

Direct with multi-target. LightGBM does not support multi-output natively; scikit-learn's MultiOutputRegressor wraps it, or you train 28 separate models. Alternatively, produce a single row per (date, horizon) pair and let the horizon be a feature — the "recursive-like" behavior of a single model with better statistical efficiency.

On our pharmacy data, direct with MultiOutputRegressor(LGBMRegressor) produces a marginal 0.3-unit improvement in MAE over the recursive strategy, mostly at the tail of the horizon.

Result on the pharmacy series

Model	MAE (units/day)
naive	24.6
seasonal_naive_7	9.8
SARIMA(1,1,1)(0,1,1)$_7$	8.4
Holt-Winters	8.9
Prophet with FR holidays	7.6
LightGBM with lags + calendar	6.9
LightGBM + planned promotions	6.3

Adding planned promotions as a known-future regressor shaves another 0.6 units. This is the row that will pressure any deep model in module 8 to earn its cost.

Global versus local models

If we have 200 stores, we can either fit one model per store (local) or a single global model trained on all stores with a store-id feature (global). The M5 confirmed that global models usually win on comparable series: they benefit from cross-store learning, share statistical strength on stores with short history, and are dramatically cheaper to operate.

The pattern is: for a chain of pharmacies, one LightGBM model with 200 store ids and 4 years of history vastly outperforms 200 tiny per-store models — even though each store's series individually looks noisier.

A regressor unknown at inference time is not a feature

Any feature that requires knowing information from beyond time $t$ to predict $y_{t+h}$ is either a leak or a liability. If the value is truly known ahead (a calendar entry, a planned campaign), use it happily. If not, either use its lagged version, forecast it first, or leave it out. There is no "in between" that is safe.

Summary

• Reshape the series into a supervised table with lag features and rolling means, all computed only from the past.
• Encode calendar features cyclically (sine, cosine) and add holidays; days_to_holiday often beats the flag alone.
• The critical distinction is regressors known in the future (calendar, planned promotions) versus unknown (weather, competitor moves); using the latter is a common cause of production failure.
• Gradient boosting on this table lands at MAE 6.9 on our pharmacy data — the current best model. Global models across stores beat one-per-store setups in practice.

Next module: LSTM and temporal Transformers. We ask when the deep approach earns its cost and when it does not — and add the last two rows to the results table.