{ "cells": [ { "cell_type": "markdown", "id": "1b33fa45", "metadata": {}, "source": [ "# BitBullet:Lessons: Forecasting Dataset Curation\n", "\n", "A forecast is a claim about the future made from a table, and the table decides\n", "whether the claim can be trusted. A certified forecasting workflow verifies\n", "more about its dataset than any other kind of tabular modelling: one series, an\n", "explicit time axis with one unambiguous period, the observed and unshifted\n", "target, and features that were genuinely knowable when each forecast would have\n", "been made.\n", "\n", "None of that verification builds the table. That part is yours, by design: you\n", "own the semantic correctness of the dataset; BitBullet owns the procedural\n", "correctness of modelling and validation. Whatever prepares your features at\n", "training time has to prepare them identically wherever the model later runs, so\n", "the preparation lives in your pipeline, where it is auditable and reusable —\n", "not hidden inside a training service where drift would be invisible.\n", "\n", "This lesson builds that table end to end from the raw Rossmann panel — the same\n", "data lesson 07 models with the SDK's forecasting toolkit. Here the product is\n", "the dataset itself:\n", "\n", "- a **base layer** of horizon-agnostic features any forecast horizon can use;\n", "- an additive **horizon block** tuned to one horizon;\n", "- a mechanical proof that neither layer can see the future; and\n", "- a controlled comparison measuring exactly what the horizon block is worth.\n", "\n", "### What You Will Build\n", "\n", "| Step | Concept |\n", "| --- | --- |\n", "| 1–2 | Select one store from a 1,115-store panel, and why one series per table |\n", "| 3–4 | A complete daily calendar, verified before any feature exists |\n", "| 5 | Two kinds of features: the base layer and the horizon block |\n", "| 6–7 | The base layer: SDK feature machinery plus the retail semantics you own |\n", "| 8 | Assemble, trim and save a contract-clean file |\n", "| 9 | A mechanical leakage audit: erase the future, prove nothing changes |\n", "| 10–11 | The horizon block, and the additivity property that makes a fair test |\n", "| 12–13 | Model both layers under the certified split protocol; judge against baselines |\n", "| 14 | A second opinion from the evaluation module |\n", "| 15 | Ship both models with a self-describing contract |\n", "| 16 | Advanced features: a field guide to the constructions beyond this recipe |\n", "\n", "**Dataset:** Rossmann Store Sales\n", "\n", "**Source:** [Kaggle – Rossmann Store Sales](https://www.kaggle.com/competitions/rossmann-store-sales/data)\n", "\n", "**Files:** `train.csv`, `store.csv`, `test.csv`\n", "\n", "> Before running this notebook, download the three files above and place them\n", "> in one directory. Point the `ROSSMANN_DIR` environment variable at that\n", "> directory (it defaults to the notebook's working directory). The estimator is\n", "> LightGBM: `pip install lightgbm`." ] }, { "cell_type": "markdown", "id": "9e972eae", "metadata": {}, "source": [ "> **Move from a reviewed dataset to a guided modelling lifecycle without wiring it together in Python.**\n", "> [BitBullet Platform](https://bitbullet.co.uk/platform/datasets) centralises datasets, profiles, managed compute, storage, and modelling workflows in one place. Use guided controls or AI-assisted draft preparation to configure supported experiments, compare completed evidence, and export the trained result with its fitted preprocessing, metadata, and generated inference code. This lesson teaches the forecasting data-curation discipline directly with the BitBullet SDK." ] }, { "cell_type": "markdown", "id": "319454af", "metadata": {}, "source": [ "## 1. Environment and Imports\n", "\n", "Everything for the calendar arithmetic — period inference, the calendar grid,\n", "the target shift, the holdout split, the baseline forecasts — comes from\n", "`bitbullet.forecast`. These are the same primitives a production training flow\n", "should call, which is what makes the protocol in section 12 a faithful\n", "rehearsal rather than an approximation." ] }, { "cell_type": "code", "execution_count": null, "id": "cd2a317a", "metadata": {}, "outputs": [], "source": [ "import os\n", "import sys\n", "from pathlib import Path\n", "from tempfile import TemporaryDirectory\n", "\n", "import lightgbm as lgb\n", "import numpy as np\n", "import pandas as pd\n", "\n", "# Prefer the local SDK when this notebook runs from the lessons directory.\n", "sdk_root = os.path.abspath('..')\n", "if sdk_root not in sys.path:\n", " sys.path.insert(0, sdk_root)\n", "\n", "from bitbullet.evaluate import evaluate_regression\n", "from bitbullet.forecast import (\n", " ForecastFeatureBuilder,\n", " ForecastFrame,\n", " ForecastSchema,\n", " GroupedRollingFeature,\n", " RollingFeature,\n", " grid_holdout_split,\n", " infer_period,\n", " naive_forecasts_on_grid,\n", " naive_scale_on_grid,\n", " place_on_grid,\n", " seasonal_period_for,\n", " shift_on_grid,\n", ")\n", "from bitbullet.model import ModelMetadata, ModelSerializer\n", "\n", "pd.set_option('display.max_columns', 30)\n", "print('Imports are ready.')" ] }, { "cell_type": "markdown", "id": "1f0338ad", "metadata": {}, "source": [ "## 2. One Series Per Table\n", "\n", "`train.csv` stacks 1,115 stores onto one shared daily calendar. A forecasting\n", "table cannot be built from that shape directly, and the reason is worth\n", "internalising: the target shift walks forward along one period grid, so on\n", "stacked entities \"seven periods ahead\" of one store's row would attach a\n", "*different store's* future as the label. One series per table is a correctness\n", "requirement, not a preference.\n", "\n", "Not every store qualifies. We require:\n", "\n", "- **complete history** — all 942 days present, so the calendar has no holes to\n", " reason about (181 stores lose six months of 2014 to refurbishments);\n", "- **a published forward plan** — the store appears in `test.csv`, which is what\n", " lets the final training origins read the operating plan for the period they\n", " predict; and\n", "- **no recording anomalies** — no days marked open with zero sales.\n", "\n", "We take store 1114: a large, Sunday-closed store of the most common type. The\n", "closures are the pedagogically interesting part — they force every decision\n", "this lesson is about." ] }, { "cell_type": "code", "execution_count": null, "id": "17a83499", "metadata": {}, "outputs": [], "source": [ "DATA_DIR = Path(os.environ.get('ROSSMANN_DIR', '.'))\n", "\n", "train_raw = pd.read_csv(DATA_DIR / 'train.csv', parse_dates=['Date'], dtype={'StateHoliday': str})\n", "future_raw = pd.read_csv(DATA_DIR / 'test.csv', parse_dates=['Date'], dtype={'StateHoliday': str})\n", "store_meta = pd.read_csv(DATA_DIR / 'store.csv')\n", "\n", "FULL_DAYS = train_raw['Date'].nunique()\n", "rows_per_store = train_raw.groupby('Store').size()\n", "complete = set(rows_per_store[rows_per_store.eq(FULL_DAYS)].index)\n", "anomalous = set(train_raw.loc[train_raw['Open'].eq(1) & train_raw['Sales'].eq(0), 'Store'])\n", "with_plan = set(future_raw['Store'])\n", "eligible = sorted((complete & with_plan) - anomalous)\n", "\n", "print(f'{rows_per_store.size} stores | {len(complete)} with all {FULL_DAYS} days | {len(eligible)} eligible')\n", "\n", "STORE = 1114\n", "SHIFT = 7 # forecast one week ahead: the horizon declared at training time\n", "\n", "display(store_meta.loc[store_meta['Store'].eq(STORE)])" ] }, { "cell_type": "markdown", "id": "21e51cf5", "metadata": {}, "source": [ "## 3. The Calendar Is the Backbone\n", "\n", "Every temporal quantity in a forecasting job — the horizon, the withheld rows,\n", "the holdout, the validation folds, and every lag you build — is counted in\n", "*periods* of the ordering axis. So the first artefact is not a feature: it is a\n", "complete daily calendar for the store.\n", "\n", "Two details matter here.\n", "\n", "**Closed days stay in.** Sunday sales are exactly zero, and it is tempting to\n", "delete those rows. Resist it: the axis period is inferred from coverage, and a\n", "daily series with every Sunday missing no longer reads as daily. Section 4\n", "shows the refusal. Zeros caused by a published closure are structure, not\n", "noise.\n", "\n", "**The calendar extends past the last observation.** `test.csv` publishes the\n", "store's operating plan — `Open`, `Promo`, `StateHoliday`, `SchoolHoliday` —\n", "for six future weeks. Appending that plan (with no sales) lets the final\n", "training origins read the plan for the period `SHIFT` days later. Without the\n", "extension, the horizon block of section 10 would be blank exactly where the\n", "freshest training rows are.\n", "\n", "Two derived columns join the observed span now, because features in section 6\n", "will group by them: the day of week, and the basket value (sales per customer,\n", "carried forward from the most recent trading day — a forward fill only ever\n", "copies the past)." ] }, { "cell_type": "code", "execution_count": null, "id": "13c2c038", "metadata": {}, "outputs": [], "source": [ "SERIES_START, SERIES_END = pd.Timestamp('2013-01-01'), pd.Timestamp('2015-07-31')\n", "COLS = ['Sales', 'Customers', 'Open', 'Promo', 'StateHoliday', 'SchoolHoliday']\n", "\n", "hist = (\n", " train_raw.loc[train_raw['Store'].eq(STORE), ['Date'] + COLS]\n", " .set_index('Date').sort_index()\n", ")\n", "plan = (\n", " future_raw.loc[future_raw['Store'].eq(STORE), ['Date', 'Open', 'Promo', 'StateHoliday', 'SchoolHoliday']]\n", " .set_index('Date').sort_index()\n", ")\n", "assert hist.index.equals(pd.date_range(SERIES_START, SERIES_END, freq='D'))\n", "\n", "cal = pd.concat([hist, plan.reindex(columns=COLS)]).sort_index()\n", "cal = cal[~cal.index.duplicated(keep='first')]\n", "cal = cal.reindex(pd.date_range(cal.index.min(), cal.index.max(), freq='D'))\n", "for column in ('Open', 'Promo', 'SchoolHoliday'):\n", " cal[column] = cal[column].fillna(0).astype(int)\n", "cal['StateHoliday'] = cal['StateHoliday'].fillna('0').astype(str)\n", "\n", "observed = cal.loc[:SERIES_END].copy()\n", "observed['DayOfWeek'] = observed.index.dayofweek\n", "open_spc = (observed['Sales'] / observed['Customers'].replace(0, np.nan))[observed['Open'].eq(1)]\n", "observed['SPC'] = open_spc.reindex(observed.index).ffill()\n", "\n", "print(f'Observed span: {observed.index[0].date()} -> {observed.index[-1].date()} ({len(observed)} days)')\n", "print(f'Plan extension: {len(cal) - len(observed)} further days of operating plan, no sales')\n", "display(cal.loc[SERIES_END - pd.Timedelta(days=1):SERIES_END + pd.Timedelta(days=2)])" ] }, { "cell_type": "markdown", "id": "7587e12f", "metadata": {}, "source": [ "## 4. Verify the Axis Before Trusting It\n", "\n", "`infer_period` is the check to run the moment an ordering axis is chosen, and\n", "before any feature exists: score candidate calendars by grid coverage, demand\n", "one unambiguous winner, and refuse the axis otherwise. Verifying during\n", "curation means the verdict at training time is never a surprise.\n", "\n", "The second call shows the anti-pattern. Delete the closed days and the best\n", "candidate calendar covers too little of its own grid to be trusted — the axis\n", "is refused, with the observed coverage in the diagnostics. The dataset with\n", "zeros passes; the \"cleaned\" dataset does not." ] }, { "cell_type": "code", "execution_count": null, "id": "aca3a585", "metadata": {}, "outputs": [], "source": [ "full_axis = pd.Series(observed.index)\n", "verdict = infer_period(full_axis)\n", "print(f'Full calendar : {verdict.status} | period {verdict.period.alias} | '\n", " f'{verdict.period.missing_slots} missing periods')\n", "\n", "open_only_axis = pd.Series(observed.index[observed['Open'].eq(1)])\n", "cleaned = infer_period(open_only_axis)\n", "best = (cleaned.period or (cleaned.alternatives or [None])[0])\n", "print(f'Open days only: {cleaned.status}'\n", " + (f' | best candidate covers {best.coverage:.1%} of its grid' if best is not None else ''))\n", "print(f' {cleaned.reason}')" ] }, { "cell_type": "markdown", "id": "d5e7dd60", "metadata": {}, "source": [ "## 5. Two Kinds of Features\n", "\n", "Every column we are about to build falls into one of two categories, and the\n", "whole design of this lesson follows from the split.\n", "\n", "**The base layer is horizon-agnostic.** Each feature describes the origin\n", "period or earlier: lags, rolling levels, momentum, footfall, and the store's\n", "own calendar and operating state on that day. Because nothing references a\n", "later period, the same file is valid at **any** declared horizon — one\n", "artefact, every horizon.\n", "\n", "**The horizon block is horizon-dependent.** Each feature describes the period\n", "being predicted: its weekday, whether the store will be open, whether a\n", "promotion is planned, the school-holiday calendar. These are legitimate\n", "*known-ahead covariates* — the evidence is that Kaggle publishes exactly these\n", "operating columns for the future test span, before any of those days had\n", "happened. The price of the sharper information is that the block is only\n", "correct at the horizon it was built for.\n", "\n", "One boundary keeps the comparison honest, and it is a rule worth adopting in\n", "your own curation: **nothing that reveals the store's forward plan goes into\n", "the base layer.** \"Days until the next promotion\" is knowable at the origin,\n", "but at a one-day horizon it would let a base-only model reconstruct tomorrow's\n", "promo flag, and the comparison between the layers would measure nothing.\n", "Forward-plan features live in the horizon block, where declaring the horizon\n", "is the point." ] }, { "cell_type": "markdown", "id": "b6415217", "metadata": {}, "source": [ "## 6. Base Layer I: The SDK's Standard Machinery\n", "\n", "`ForecastFeatureBuilder` produces the structural families every forecasting\n", "dataset needs, under one convention worth memorising: **lags are one-based and\n", "origin-inclusive**. `lag 1` is the value observed *at* the origin, the freshest\n", "information a forecast can use. That is legitimate because the declared horizon\n", "is at least 1, so every label lies strictly after its origin. `Sales__lag_1`\n", "is also, by definition, the naive forecast — the single most informative raw\n", "input a forecasting model gets.\n", "\n", "Four families, one declaration each:\n", "\n", "- **Lags** of the target and of footfall (`target_lags`, `exogenous_lags`).\n", "- **Rolling statistics** — level, spread and extremes over trailing windows.\n", "- **Calendar components** — deterministic date arithmetic.\n", "- **Grouped rolling statistics** (`GroupedRollingFeature`) — a rolling\n", " statistic over one *group's* history, read as of each row. Grouped by the\n", " day of week, `window=8` gives the store's trailing eight-Thursday mean on\n", " every Thursday: the weekday profile that plain lags and windows cannot\n", " express. The statistic is computed on each group's own past and carried\n", " forward, never backwards, so it obeys the same leakage discipline as\n", " everything else here — section 9 will prove that mechanically rather than\n", " take the API's word for it." ] }, { "cell_type": "code", "execution_count": null, "id": "08c06c8a", "metadata": {}, "outputs": [], "source": [ "SCHEMA = ForecastSchema(time='Date', targets='Sales', cadence='D')\n", "\n", "BUILDER = ForecastFeatureBuilder(\n", " target_lags=(1, 2, 3, 4, 5, 6, 7, 8, 15, 22, 29, 36),\n", " exogenous_lags={'Customers': (1, 2, 8)},\n", " rolling_features=[\n", " RollingFeature(column='Sales', window=w, statistic=s)\n", " for w, s in [\n", " (7, 'mean'), (14, 'mean'), (28, 'mean'), (91, 'mean'),\n", " (7, 'median'), (28, 'median'),\n", " (7, 'std'), (28, 'std'), (91, 'std'),\n", " (28, 'min'), (28, 'max'),\n", " ]\n", " ] + [\n", " RollingFeature(column='Customers', window=7, statistic='mean'),\n", " RollingFeature(column='Customers', window=28, statistic='mean'),\n", " ],\n", " grouped_rolling_features=[\n", " GroupedRollingFeature(column='Sales', by='DayOfWeek', window=4, statistic='mean'),\n", " GroupedRollingFeature(column='Sales', by='DayOfWeek', window=8, statistic='mean'),\n", " GroupedRollingFeature(column='Sales', by='DayOfWeek', window=13, statistic='mean'),\n", " GroupedRollingFeature(column='Sales', by='DayOfWeek', window=8, statistic='std'),\n", " GroupedRollingFeature(column='Customers', by='DayOfWeek', window=4, statistic='mean'),\n", " GroupedRollingFeature(column='SPC', by='DayOfWeek', window=4, statistic='mean', min_periods=1),\n", " ],\n", " calendar_features=(\n", " 'quarter', 'month', 'week', 'day', 'day_of_week', 'day_of_year',\n", " 'is_weekend', 'is_month_start', 'is_month_end',\n", " 'month_sin', 'month_cos', 'day_of_week_sin', 'day_of_week_cos',\n", " ),\n", ")\n", "\n", "preview = BUILDER.transform(ForecastFrame(observed.reset_index(names='Date'), SCHEMA))\n", "sdk_columns = [column for column in preview.columns if '__' in column]\n", "print(f'{len(sdk_columns)} builder features, e.g.:')\n", "print('\\n'.join(' ' + column for column in sdk_columns[:5]))\n", "print('\\n'.join(' ' + column for column in sdk_columns if 'grouped' in column))" ] }, { "cell_type": "markdown", "id": "bca8bdb2", "metadata": {}, "source": [ "## 7. Base Layer II: The Retail Semantics You Own\n", "\n", "The builder covers what is generic. What it deliberately does not decide is\n", "what *this* series means — and that is where the accuracy lives. Four families\n", "encode how a Sunday-closed German drugstore actually behaves:\n", "\n", "| Family | Question it answers |\n", "| --- | --- |\n", "| Trading-day windows | What is the store's true recent level, with structural closures excluded? A calendar-window mean is dragged down by every Sunday. |\n", "| Conditioned levels | What does a promo day look like here, versus a plain day? Computed as trading-day means inside each state. |\n", "| Momentum | Is the store above or below its own recent behaviour, and by how much? Differences, ratios and a z-score against the 28-day window. |\n", "| Event clocks | How long has the current promotion run? How long since a closure? The day after a closure carries pent-up demand. |\n", "\n", "One implementation rule makes these leak-free by construction: every window is\n", "*inclusive of the origin and looks only backwards*. Nothing here — and nothing\n", "in the builder's output — reads a row later than its own date." ] }, { "cell_type": "code", "execution_count": null, "id": "84215635", "metadata": {}, "outputs": [], "source": [ "STATE_HOLIDAY_LABELS = {'0': 'none', 'a': 'public', 'b': 'easter', 'c': 'christmas'}\n", "\n", "\n", "def run_length(flag):\n", " \"\"\"Length of the consecutive run of True ending at (and including) each row.\"\"\"\n", " block = (flag != flag.shift()).cumsum()\n", " return flag.groupby(block).cumcount().add(1).where(flag, 0).astype(int)\n", "\n", "\n", "def days_since(flag, cap):\n", " \"\"\"Days since the most recent True at or before each row (0 on the event).\"\"\"\n", " position = pd.Series(np.arange(len(flag)), index=flag.index, dtype=float)\n", " marker = position.where(flag).ffill()\n", " return (position - marker).fillna(cap).clip(upper=cap).astype(int)\n", "\n", "\n", "def days_until(flag, cap):\n", " \"\"\"Days until the next True at or after each row (0 on the event).\n", "\n", " Reads forward, so it belongs only in the horizon block (section 10).\n", " \"\"\"\n", " position = pd.Series(np.arange(len(flag)), index=flag.index, dtype=float)\n", " marker = position.where(flag).bfill()\n", " return (marker - position).fillna(cap).clip(upper=cap).astype(int)\n", "\n", "\n", "def conditioned_mean(sales, flag, trading, window, state_of_interest):\n", " \"\"\"Trading-day mean of sales inside/outside a binary state, chosen per row.\n", "\n", " Each row reads the mean matching ``state_of_interest`` — the state of the\n", " period the feature describes (the origin for the base layer, the target\n", " period for the horizon block). Windows holding no trading day in the needed\n", " state fall back to the overall trading-day level rather than leaving holes.\n", " Returns the conditioned mean and the on/off uplift ratio.\n", " \"\"\"\n", " live = sales.where(trading, 0.0)\n", " on = flag.where(trading, 0)\n", " off = (1 - flag).where(trading, 0)\n", " mean_on = (live * on).rolling(window).sum() / on.rolling(window).sum().replace(0, np.nan)\n", " mean_off = (live * off).rolling(window).sum() / off.rolling(window).sum().replace(0, np.nan)\n", " level = (\n", " sales[trading].rolling(window, min_periods=max(2, window // 2)).mean()\n", " .reindex(sales.index).ffill()\n", " )\n", " chosen = pd.Series(\n", " np.where(state_of_interest.eq(1), mean_on, mean_off), index=sales.index\n", " )\n", " uplift = (mean_on / mean_off).replace([np.inf, -np.inf], np.nan).fillna(1.0)\n", " return chosen.fillna(level), uplift\n", "\n", "print('Curation helpers are defined.')" ] }, { "cell_type": "code", "execution_count": null, "id": "3281a40b", "metadata": {}, "outputs": [], "source": [ "def build_base(observed_cal):\n", " \"\"\"Every base-layer feature for one store, indexed by origin date.\n", "\n", " The contract of this function is the contract of the base layer: nothing\n", " may read a row later than its own index. Section 9 enforces that claim\n", " mechanically rather than taking the code's word for it.\n", " \"\"\"\n", " idx = observed_cal.index\n", " sales = observed_cal['Sales'].astype(float)\n", " customers = observed_cal['Customers'].astype(float)\n", " trading = observed_cal['Open'].eq(1)\n", " promo = observed_cal['Promo'].astype(int)\n", " school = observed_cal['SchoolHoliday'].astype(int)\n", " holiday = observed_cal['StateHoliday'].map(STATE_HOLIDAY_LABELS).fillna('none')\n", "\n", " # -- the SDK's standard machinery: lags, rolling, grouped, calendar ------\n", " enriched = BUILDER.transform(\n", " ForecastFrame(observed_cal.reset_index(names='Date'), SCHEMA)\n", " ).set_index('Date')\n", " f = enriched[[column for column in enriched.columns if '__' in column]].copy()\n", "\n", " # -- the origin day's own operating state --------------------------------\n", " f['o_is_open'] = trading.astype(int)\n", " f['o_promo'] = promo\n", " f['o_state_holiday'] = holiday\n", " f['o_is_state_holiday'] = holiday.ne('none').astype(int)\n", " f['o_school_holiday'] = school\n", "\n", " # -- calendar the builder does not cover ---------------------------------\n", " f['o_trend_index'] = (idx - idx[0]).days\n", " f['o_days_to_month_end'] = ((idx + pd.offsets.MonthEnd(0)) - idx).days\n", " f['o_is_payday_window'] = ((idx.day <= 5) | (idx.day >= 25)).astype(int)\n", " christmas = pd.to_datetime([f'{year}-12-24' for year in idx.year])\n", " to_christmas = (christmas - idx).days\n", " f['o_days_to_christmas'] = np.clip(np.where(to_christmas < 0, 365 + to_christmas, to_christmas), 0, 60)\n", "\n", " # -- event clocks (all backward-looking) ---------------------------------\n", " f['o_promo_run_len'] = run_length(promo.eq(1))\n", " f['o_days_since_last_promo'] = days_since(promo.eq(1), cap=60)\n", " f['o_days_since_prev_closed'] = days_since(~trading, cap=30)\n", " f['o_is_day_after_closed'] = (~trading).shift(1, fill_value=False).astype(int)\n", " f['o_days_since_prev_state_holiday'] = days_since(holiday.ne('none'), cap=120)\n", " f['o_school_holiday_run_len'] = run_length(school.eq(1))\n", "\n", " # -- trading-day windows: the store's true level, closures excluded ------\n", " open_sales = sales[trading]\n", " for window in (7, 28, 91):\n", " f[f'sales_open_mean_{window}'] = (\n", " open_sales.rolling(window, min_periods=max(2, int(window * 0.7))).mean()\n", " .reindex(idx).ffill()\n", " )\n", " f['sales_open_std_28'] = open_sales.rolling(28, min_periods=20).std().reindex(idx).ffill()\n", " f['customers_open_mean_28'] = (\n", " customers[trading].rolling(28, min_periods=20).mean().reindex(idx).ffill()\n", " )\n", " f['spc_open_mean_28'] = (\n", " observed_cal['SPC'][trading].rolling(28, min_periods=20).mean().reindex(idx).ffill()\n", " )\n", "\n", " # -- exponential smoothing and momentum ----------------------------------\n", " f['sales_ewm_7'] = sales.ewm(span=7, adjust=False).mean()\n", " f['sales_ewm_28'] = sales.ewm(span=28, adjust=False).mean()\n", " for delta in (1, 7, 28):\n", " f[f'sales_diff_{delta}'] = sales.diff(delta)\n", " f['sales_vs_open_mean_28'] = (sales / f['sales_open_mean_28'].replace(0, np.nan)).fillna(1.0)\n", " f['sales_zscore_28'] = (\n", " (sales - f['Sales__rolling_mean_lag_1_window_28'])\n", " / f['Sales__rolling_std_lag_1_window_28'].replace(0, np.nan)\n", " ).fillna(0.0)\n", " f['sales_open7_vs_open91'] = (\n", " f['sales_open_mean_7'] / f['sales_open_mean_91'].replace(0, np.nan)\n", " ).fillna(1.0)\n", "\n", " # -- history conditioned on the origin day's promotion / holiday state ---\n", " f['sales_opromo_mean_28'], _ = conditioned_mean(sales, promo, trading, 28, promo)\n", " f['sales_opromo_mean_91'], f['opromo_uplift_91'] = conditioned_mean(sales, promo, trading, 91, promo)\n", " f['sales_oschool_mean_91'], f['oschool_uplift_91'] = conditioned_mean(sales, school, trading, 91, school)\n", "\n", " # -- how unusual the recent stretch has been -----------------------------\n", " f['n_promo_last_28'] = promo.rolling(28).sum()\n", " f['n_closed_last_28'] = (~trading).astype(int).rolling(28).sum()\n", " f['n_school_last_28'] = school.rolling(28).sum()\n", " return f\n", "\n", "\n", "base_features = build_base(observed)\n", "print(f'Base layer: {base_features.shape[1]} features over {base_features.shape[0]} days')" ] }, { "cell_type": "markdown", "id": "a35c63bd", "metadata": {}, "source": [ "## 8. Assemble, Trim, and Save the Base Dataset\n", "\n", "Three finishing rules turn a feature table into a contract-clean file.\n", "\n", "**Trim the warm-up.** The first 98 days are discarded: every window of 91 days\n", "or less, and every same-weekday statistic of up to 13 occurrences, is fully\n", "populated from day 98 onward. Discarding rows is honest; imputing a 91-day mean\n", "from 12 days of history is not.\n", "\n", "**Drop dead columns.** Constant columns carry nothing. Exact duplicates arise\n", "naturally — this store never trades on a state holiday, so two clocks can\n", "coincide — and keeping one of each pair is enough.\n", "\n", "**Refuse missing values.** The finished file must have none. A hole this late\n", "means the recipe is wrong somewhere, and an assertion is the correct response,\n", "not a fill.\n", "\n", "The file keeps the observed, **unshifted** target. The horizon is declared once\n", "at training time and applied on the verified calendar grid — section 12 does\n", "exactly that with `shift_on_grid` — and `horizon − 1` rows are withheld at\n", "every chronological boundary. A target you shift yourself in the file has no\n", "correct horizon declaration left: any training flow that shifts on the grid\n", "would double-shift it, and one that does not cannot state what the prediction\n", "means. Ship the observed target; declare the horizon once." ] }, { "cell_type": "code", "execution_count": null, "id": "99a955ad", "metadata": {}, "outputs": [], "source": [ "WARMUP = 98\n", "OUT_DIR = Path('curated_datasets')\n", "OUT_DIR.mkdir(exist_ok=True)\n", "\n", "\n", "def finalise(features, filename):\n", " \"\"\"Trim the warm-up, drop dead columns, refuse missing values, save.\"\"\"\n", " frame = features.iloc[WARMUP:].copy()\n", " kept, dropped = {}, []\n", " for name, series in frame.items():\n", " if series.dtype.kind in 'iu':\n", " series = series.astype(np.int64)\n", " elif series.dtype.kind == 'f':\n", " series = series.astype(np.float64)\n", " if series.nunique(dropna=False) <= 1:\n", " dropped.append(f'{name} (constant)')\n", " continue\n", " twin = next((k for k, v in kept.items() if v.dtype == series.dtype and v.equals(series)), None)\n", " if twin is not None:\n", " dropped.append(f'{name} (duplicate of {twin})')\n", " continue\n", " kept[name] = series\n", " frame = pd.DataFrame(kept)\n", "\n", " out = frame.copy()\n", " out.insert(0, 'sales', observed.loc[frame.index, 'Sales'].astype(float))\n", " out.insert(0, 'date', frame.index.strftime('%Y-%m-%d'))\n", " missing = out.isna().sum()\n", " assert missing.eq(0).all(), f'missing values: {missing[missing.gt(0)].to_dict()}'\n", "\n", " path = OUT_DIR / filename\n", " out.reset_index(drop=True).to_csv(path, index=False)\n", " print(f'{filename}: {out.shape[0]} rows x {out.shape[1]} columns '\n", " f'({out.shape[1] - 2} features); dropped {len(dropped)} dead columns')\n", " for entry in dropped:\n", " print(' -', entry)\n", " return path\n", "\n", "\n", "BASE_PATH = finalise(base_features, f'rossmann_store{STORE}_daily_base.csv')" ] }, { "cell_type": "markdown", "id": "9b3651a2", "metadata": {}, "source": [ "## 9. Prove It Cannot See the Future\n", "\n", "Code review can argue a feature is leak-free; an audit can prove it. The test\n", "is brutal and simple: rebuild the entire base layer from a history *truncated\n", "at the origin*, and compare the origin's row against the same row built from\n", "the full history. A feature that used any observation after the origin will\n", "change; a feature that changes is a leak. Identical rows are proof — not\n", "assertion — that every value was computable on the day it claims to describe.\n", "\n", "Note what is being audited: not just the hand-written helpers but the SDK's\n", "lag, rolling and grouped rolling output too. An audit that exempts library\n", "code is a courtesy, not an audit." ] }, { "cell_type": "code", "execution_count": null, "id": "6d1634ad", "metadata": {}, "outputs": [], "source": [ "audit_index = base_features.index[WARMUP:]\n", "probes = audit_index[np.linspace(0, len(audit_index) - 1, 6).astype(int)]\n", "\n", "leaks = 0\n", "for origin in probes:\n", " partial = build_base(observed.loc[:origin])\n", " full_row, partial_row = base_features.loc[origin], partial.loc[origin]\n", " for column in base_features.columns:\n", " a, b = full_row[column], partial_row[column]\n", " if isinstance(a, str):\n", " same = a == b\n", " else:\n", " same = (pd.isna(a) and pd.isna(b)) or np.isclose(\n", " float(a), float(b), rtol=1e-9, atol=1e-9\n", " )\n", " if not same:\n", " leaks += 1\n", " print(f'LEAK {column} @ {origin.date()}: {a!r} != {b!r}')\n", "\n", "print(f'Probed {len(probes)} origins x {base_features.shape[1]} features: '\n", " + ('clean — no feature changes when the future is erased.' if leaks == 0\n", " else f'{leaks} leaking values found.'))" ] }, { "cell_type": "markdown", "id": "e7175b41", "metadata": {}, "source": [ "## 10. The Horizon Block\n", "\n", "Now fix the horizon and buy the sharper information. Every column here\n", "describes the period `SHIFT` days after the row — prefix `t_` — and comes from\n", "sources that are published in advance: the civil calendar and the store's own\n", "operating plan.\n", "\n", "Three constructions, in increasing order of subtlety:\n", "\n", "- **Plan and calendar at the target period.** Build the features on the\n", " plan-extended calendar, then read them at `date + SHIFT`. This is where the\n", " forward-plan clocks (`days_until` and friends) finally belong: declaring the\n", " horizon is the point of this block.\n", "- **History conditioned on the target period's state.** The store's trading\n", " level *under promotion* versus *without* is computed from history at the\n", " origin — but which of the two a row reads is chosen by the **target**\n", " period's planned promo state. Known-ahead conditioning, backward-looking\n", " statistics: sharp and leak-free.\n", "- **The weekly special case.** At `SHIFT = 7` the target weekday equals the\n", " origin weekday, so the same-weekday statistics of section 6 are already\n", " aligned to the target — a seven-day horizon gets the strongest alignment for\n", " free. At any other shift, redirect them with the builder's `key` argument: a\n", " column holding the *target* period's weekday (deterministic calendar\n", " arithmetic, known ahead) makes each row read the weekday profile of the\n", " period it predicts. Section 16 shows the construction." ] }, { "cell_type": "code", "execution_count": null, "id": "dfb9988c", "metadata": {}, "outputs": [], "source": [ "def plan_features(cal_ext):\n", " \"\"\"Calendar and operating plan for every period, with forward-plan clocks.\"\"\"\n", " idx = cal_ext.index\n", " trading = cal_ext['Open'].eq(1)\n", " promo = cal_ext['Promo'].astype(int)\n", " school = cal_ext['SchoolHoliday'].astype(int)\n", " holiday = cal_ext['StateHoliday'].map(STATE_HOLIDAY_LABELS).fillna('none')\n", "\n", " p = pd.DataFrame(index=idx)\n", " dow = idx.dayofweek\n", " p['day_of_week'] = dow.astype(np.int64)\n", " p['dow_sin'] = np.sin(2 * np.pi * dow / 7)\n", " p['dow_cos'] = np.cos(2 * np.pi * dow / 7)\n", " p['day_of_month'] = idx.day\n", " p['days_to_month_end'] = ((idx + pd.offsets.MonthEnd(0)) - idx).days\n", " p['is_payday_window'] = ((idx.day <= 5) | (idx.day >= 25)).astype(int)\n", " p['month'] = idx.month\n", " p['month_sin'] = np.sin(2 * np.pi * (idx.month - 1) / 12)\n", " p['month_cos'] = np.cos(2 * np.pi * (idx.month - 1) / 12)\n", " p['doy_sin'] = np.sin(2 * np.pi * idx.dayofyear / 365)\n", " p['doy_cos'] = np.cos(2 * np.pi * idx.dayofyear / 365)\n", " christmas = pd.to_datetime([f'{year}-12-24' for year in idx.year])\n", " to_christmas = (christmas - idx).days\n", " p['days_to_christmas'] = np.clip(np.where(to_christmas < 0, 365 + to_christmas, to_christmas), 0, 60)\n", "\n", " p['is_open'] = trading.astype(int)\n", " p['promo'] = promo\n", " p['state_holiday'] = holiday\n", " p['is_state_holiday'] = holiday.ne('none').astype(int)\n", " p['school_holiday'] = school\n", "\n", " p['promo_run_len'] = run_length(promo.eq(1))\n", " p['promo_prev_day'] = promo.shift(1).fillna(0).astype(int)\n", " p['promo_next_day'] = promo.shift(-1).fillna(0).astype(int)\n", " p['days_to_next_promo'] = days_until(promo.eq(1), cap=60)\n", " p['days_since_last_promo'] = days_since(promo.eq(1), cap=60)\n", " p['days_to_next_state_holiday'] = days_until(holiday.ne('none'), cap=120)\n", " p['days_to_next_closed'] = days_until(~trading, cap=30)\n", " p['is_day_after_closed'] = (~trading).shift(1, fill_value=False).astype(int)\n", " p['is_day_before_closed'] = (~trading).shift(-1, fill_value=False).astype(int)\n", " return p\n", "\n", "\n", "origins = observed.index\n", "target_dates = origins + pd.Timedelta(days=SHIFT)\n", "\n", "t_block = plan_features(cal).reindex(target_dates)\n", "t_block.index = origins\n", "t_block = t_block.add_prefix('t_')\n", "\n", "sales_series = observed['Sales'].astype(float)\n", "trading_series = observed['Open'].eq(1)\n", "promo_series = observed['Promo'].astype(int)\n", "school_series = observed['SchoolHoliday'].astype(int)\n", "target_promo = pd.Series(cal['Promo'].reindex(target_dates).to_numpy(), index=origins).astype(int)\n", "target_school = pd.Series(cal['SchoolHoliday'].reindex(target_dates).to_numpy(), index=origins).astype(int)\n", "\n", "t_cond = pd.DataFrame(index=origins)\n", "t_cond['sales_tpromo_mean_28'], _ = conditioned_mean(\n", " sales_series, promo_series, trading_series, 28, target_promo)\n", "t_cond['sales_tpromo_mean_91'], t_cond['tpromo_uplift_91'] = conditioned_mean(\n", " sales_series, promo_series, trading_series, 91, target_promo)\n", "t_cond['sales_tschool_mean_91'], _ = conditioned_mean(\n", " sales_series, school_series, trading_series, 91, target_school)\n", "\n", "horizon_features = pd.concat([base_features, t_block, t_cond], axis=1)\n", "HORIZON_PATH = finalise(horizon_features, f'rossmann_store{STORE}_daily_h{SHIFT}.csv')" ] }, { "cell_type": "markdown", "id": "3f9093e5", "metadata": {}, "source": [ "## 11. Additivity: One Controlled Variable\n", "\n", "The horizon file must carry the base layer *unchanged* — the same columns with\n", "the same values, plus the `t_` block. That is not tidiness; it is experimental\n", "design. Because the two files now differ in exactly one thing, the difference\n", "between their holdout scores in section 13 measures the horizon block and\n", "nothing else. If curation had quietly altered a base column along the way, the\n", "comparison would measure two things at once and answer neither." ] }, { "cell_type": "code", "execution_count": null, "id": "71aa4315", "metadata": {}, "outputs": [], "source": [ "base_file = pd.read_csv(BASE_PATH)\n", "horizon_file = pd.read_csv(HORIZON_PATH)\n", "\n", "missing = [column for column in base_file.columns if column not in horizon_file.columns]\n", "altered = [column for column in base_file.columns\n", " if column not in missing and not base_file[column].equals(horizon_file[column])]\n", "assert not missing, f'horizon file lost base columns: {missing}'\n", "assert not altered, f'horizon file altered base columns: {altered}'\n", "extra = horizon_file.shape[1] - base_file.shape[1]\n", "print(f'All {base_file.shape[1]} base columns are byte-identical in the horizon file; '\n", " f'{extra} horizon columns added.')" ] }, { "cell_type": "markdown", "id": "cf5b5fd8", "metadata": {}, "source": [ "## 12. Model Both Under the Certified Split Protocol\n", "\n", "The split protocol every certified forecasting evaluation should follow,\n", "rehearsed with the SDK primitives built for it:\n", "\n", "1. `infer_period` verifies the axis and names the period;\n", "2. `place_on_grid` puts every row on the calendar grid;\n", "3. `shift_on_grid` attaches the label from `SHIFT` periods ahead — rows whose\n", " label slot is empty simply stop being origins;\n", "4. `grid_holdout_split` reserves the final 15% of *periods* as the\n", " chronological holdout, withholding `SHIFT − 1 = 6` periods before the\n", " boundary, because those origins carry labels from after the cutoff.\n", "\n", "A production training flow would tune each model with expanding time-series\n", "folds inside the training span. Here both models get identical fixed\n", "parameters — deliberately, because the question is about the datasets, and the\n", "answer should not depend on who won a lottery of trials." ] }, { "cell_type": "code", "execution_count": null, "id": "387547dd", "metadata": {}, "outputs": [], "source": [ "HOLDOUT = 0.15\n", "PARAMS = dict(\n", " objective='regression_l1', n_estimators=700, learning_rate=0.03,\n", " num_leaves=31, min_child_samples=20, subsample=0.9, subsample_freq=1,\n", " colsample_bytree=0.8, reg_lambda=1.0, random_state=42, verbose=-1,\n", ")\n", "\n", "\n", "def certified_protocol(path, shift, holdout=HOLDOUT):\n", " \"\"\"Prepare one curated file exactly as a certified forecasting job should.\"\"\"\n", " frame = pd.read_csv(path, parse_dates=['date'])\n", " diagnosis = infer_period(frame['date'])\n", " assert diagnosis.status == 'regular', diagnosis.reason\n", " grid = place_on_grid(frame['date'], diagnosis.period)\n", "\n", " labels = shift_on_grid(frame['sales'].to_numpy(float), grid, shift)\n", " is_origin = np.isfinite(labels)\n", "\n", " X = frame.drop(columns=['date', 'sales']).loc[is_origin].reset_index(drop=True)\n", " for column in X.columns:\n", " if X[column].dtype == object:\n", " X[column] = X[column].astype('category')\n", " y = labels[is_origin]\n", "\n", " origin_grid = place_on_grid(frame['date'][is_origin], diagnosis.period)\n", " split = grid_holdout_split(\n", " origin_grid, test_size=holdout, gap=shift - 1, min_train_rows=2, min_test_rows=1\n", " )\n", " return {\n", " 'frame': frame, 'grid': grid, 'diagnosis': diagnosis, 'is_origin': is_origin,\n", " 'X': X, 'y': y,\n", " 'train': np.asarray(split.rows.train_indices),\n", " 'test': np.asarray(split.rows.test_indices),\n", " 'dates': frame['date'][is_origin].reset_index(drop=True),\n", " }\n", "\n", "\n", "runs = {}\n", "for label, path in [('base', BASE_PATH), ('horizon', HORIZON_PATH)]:\n", " prepared = certified_protocol(path, SHIFT)\n", " model = lgb.LGBMRegressor(**PARAMS).fit(\n", " prepared['X'].iloc[prepared['train']], prepared['y'][prepared['train']]\n", " )\n", " predictions = model.predict(prepared['X'].iloc[prepared['test']])\n", " runs[label] = {\n", " 'prepared': prepared, 'model': model,\n", " 'pred': predictions, 'truth': prepared['y'][prepared['test']],\n", " }\n", " span = (prepared['dates'].iloc[prepared['test'][0]].date(),\n", " prepared['dates'].iloc[prepared['test'][-1]].date())\n", " print(f\"{label:>8}: {len(prepared['train'])} training origins | \"\n", " f\"{SHIFT - 1} withheld | {len(prepared['test'])} holdout \"\n", " f\"({span[0]} -> {span[1]}) | {prepared['X'].shape[1]} features\")" ] }, { "cell_type": "markdown", "id": "bd17310b", "metadata": {}, "source": [ "## 13. Judge Against the Baselines Every Forecast Must Beat\n", "\n", "Accuracy numbers mean nothing without a yardstick. The observed unshifted\n", "target is what makes the two canonical baselines computable with no extra\n", "input: the **naive** forecast (each origin's own last observed value) and the\n", "**seasonal naive** forecast (one weekly cycle back). MASE scales the model's\n", "error by the in-sample naive error over the training periods — Hyndman &\n", "Koehler's scale-free yardstick — and the skill score states the headline\n", "plainly: the fraction of baseline error the model removed. A model that does\n", "not beat seasonal naive is not yet earning its complexity.\n", "\n", "At a weekly horizon the two baselines coincide: the value seven periods before\n", "the target *is* the most recent same-weekday value, so naive and seasonal\n", "naive make the same prediction. The table therefore shows one baseline row.\n", "At any horizon that is not a multiple of the seasonal cycle they separate, and\n", "seasonal naive is usually the one to beat.\n", "\n", "One honesty rule for closed-day series: a closed Sunday is a structural zero\n", "that any sane model predicts exactly, so overall R² flatters everyone. The\n", "trading-day columns restate the result on the days that were actually open —\n", "quote those." ] }, { "cell_type": "code", "execution_count": null, "id": "caed4dde", "metadata": {}, "outputs": [], "source": [ "shared = runs['base']['prepared']\n", "seasonal_period = seasonal_period_for(shared['diagnosis'].period)\n", "\n", "baseline_paths = naive_forecasts_on_grid(\n", " shared['frame']['sales'].to_numpy(float), shared['grid'],\n", " horizon=SHIFT, seasonal_period=seasonal_period,\n", ")\n", "test_rows = shared['test']\n", "naive_pred = baseline_paths['naive'][shared['is_origin']][test_rows]\n", "seasonal_pred = baseline_paths['seasonal_naive'][shared['is_origin']][test_rows]\n", "\n", "train_dates = shared['dates'].iloc[shared['train']]\n", "train_values = shared['frame']['sales'].to_numpy(float)[shared['is_origin']][shared['train']]\n", "train_grid = place_on_grid(train_dates.reset_index(drop=True), shared['diagnosis'].period)\n", "seasonal_scale = naive_scale_on_grid(train_values, train_grid, seasonal_period=seasonal_period)\n", "\n", "truth = runs['base']['truth']\n", "\n", "\n", "def score(pred):\n", " error = truth - pred\n", " open_days = truth > 0\n", " open_error = error[open_days]\n", " open_truth = truth[open_days]\n", " return {\n", " 'MAE': np.mean(np.abs(error)),\n", " 'RMSE': np.sqrt(np.mean(error ** 2)),\n", " 'Seasonal MASE': np.mean(np.abs(error)) / seasonal_scale,\n", " 'Trading-day MAE': np.mean(np.abs(open_error)),\n", " 'Trading-day R2': 1 - np.sum(open_error ** 2) / np.sum((open_truth - open_truth.mean()) ** 2),\n", " }\n", "\n", "\n", "rows = {'base model': score(runs['base']['pred']),\n", " 'horizon model': score(runs['horizon']['pred'])}\n", "if np.allclose(naive_pred, seasonal_pred):\n", " rows['naive = seasonal naive'] = score(naive_pred)\n", "else:\n", " rows['naive'] = score(naive_pred)\n", " rows['seasonal naive'] = score(seasonal_pred)\n", "comparison = pd.DataFrame(rows).T.round(3)\n", "display(comparison)\n", "\n", "for label in ('base', 'horizon'):\n", " mae = comparison.loc[f'{label} model', 'MAE']\n", " skill = 1 - mae / comparison.iloc[-1]['MAE']\n", " print(f'{label:>8} model: {skill:.1%} of seasonal-naive error removed')\n", "gain = 1 - comparison.loc['horizon model', 'MAE'] / comparison.loc['base model', 'MAE']\n", "print(f'\\nHorizon block: {gain:.1%} lower holdout MAE than the base layer alone —')\n", "print('the measured worth of knowing the operating plan for the predicted period.')" ] }, { "cell_type": "markdown", "id": "26d821c2", "metadata": {}, "source": [ "## 14. A Second Opinion from the Evaluation Module\n", "\n", "`evaluate_regression` produces the fuller report a regression evaluation\n", "deserves — error metrics plus residual, target and prediction summaries. The\n", "residual mean is worth a glance on every forecasting run: a model can post a\n", "respectable MAE while sitting persistently above or below the truth, and a\n", "biased forecast quietly misallocates stock every single week.\n", "\n", "Expect to see that here. This store trends strongly upward, and tree ensembles\n", "do not extrapolate — they predict within the range they were trained on — so\n", "both models sit below the truth across the holdout, and the residual mean is\n", "where that shows first. The momentum and trading-day-level features shrink the\n", "bias by keeping \"where the store is now\" in view; the residual summary tells\n", "you how much of it remains." ] }, { "cell_type": "code", "execution_count": null, "id": "7ad703fd", "metadata": {}, "outputs": [], "source": [ "reports = {\n", " label: evaluate_regression(run['truth'], run['pred'], include_values=False)\n", " for label, run in runs.items()\n", "}\n", "\n", "metric_table = pd.DataFrame({label: report.metrics for label, report in reports.items()}).round(3)\n", "display(metric_table)\n", "\n", "residual_table = pd.DataFrame(\n", " {label: report.residual_summary for label, report in reports.items()}\n", ").round(1)\n", "display(residual_table)" ] }, { "cell_type": "markdown", "id": "33aff2dd", "metadata": {}, "source": [ "## 15. Ship the Models with Their Contract\n", "\n", "A forecasting artefact has to explain itself: what a prediction *means* — the\n", "target value `SHIFT` periods after the supplied row — is not recoverable from\n", "a pickle. The metadata carries it: the horizon, the ordering column, the\n", "dataset the model was fitted on, its feature list and its holdout scores.\n", "`ModelSerializer` packages model and metadata together, and the reload-replay\n", "check proves the shipped artefact reproduces the exact holdout predictions it\n", "was scored on." ] }, { "cell_type": "code", "execution_count": null, "id": "b4245a5b", "metadata": {}, "outputs": [], "source": [ "with TemporaryDirectory() as artifact_dir:\n", " for label, run in runs.items():\n", " prepared = run['prepared']\n", " source = BASE_PATH if label == 'base' else HORIZON_PATH\n", " metadata = ModelMetadata(\n", " name=f'rossmann_{STORE}_{label}_h{SHIFT}',\n", " version='1.0.0',\n", " model_type='LGBMRegressor',\n", " framework='lightgbm',\n", " task='regression',\n", " hyperparameters=PARAMS,\n", " feature_names=list(prepared['X'].columns),\n", " n_features=prepared['X'].shape[1],\n", " metrics={\n", " 'holdout_mae': float(np.mean(np.abs(run['truth'] - run['pred']))),\n", " 'holdout_seasonal_mase': float(\n", " np.mean(np.abs(run['truth'] - run['pred'])) / seasonal_scale),\n", " },\n", " notes=(\n", " f'Single-series forecasting; horizon {SHIFT} periods (D); order by date; '\n", " f'observed unshifted target; dataset={source.name}; '\n", " f'final {HOLDOUT:.0%} of periods held out chronologically.'\n", " ),\n", " )\n", " artifact = Path(artifact_dir) / f'{metadata.name}.pkl'\n", " ModelSerializer.save(run['model'], artifact, metadata=metadata, include_datasets=False)\n", "\n", " package = ModelSerializer.load(artifact)\n", " replayed = package.model.predict(prepared['X'].iloc[prepared['test']])\n", " assert np.allclose(replayed, run['pred']), 'reloaded model diverged from its record'\n", " print(f\"{artifact.name}: saved, reloaded, replayed — \"\n", " f\"holdout MAE {package.metadata.metrics['holdout_mae']:.1f} confirmed\")" ] }, { "cell_type": "markdown", "id": "156df8fa", "metadata": {}, "source": [ "## 16. Advanced Features: A Field Guide\n", "\n", "The recipe above is deliberately complete without anything in this section.\n", "What follows is vocabulary, not a shopping list: six constructions that come\n", "up in serious forecasting work, each introduced, built, and shown — **and none\n", "of them added to the shipped datasets**. The discipline is the point. Every\n", "candidate feature must re-earn its place under the audit of section 9 and the\n", "controlled comparison of section 13; we measured the first and last families\n", "here on this very store under the section 12 protocol, and neither produced a\n", "consistent improvement over the recipe you already have. Features are not free\n", "— each one buys either signal or noise, and only a measurement tells you\n", "which.\n", "\n", "**Fourier harmonic stacks.** A single sine/cosine pair can only draw one\n", "smooth annual bump. Stacking `K` harmonics — periods of a year, half a year, a\n", "third of a year — lets a model represent multi-modal seasonal shapes (a\n", "Christmas peak *and* a summer dip) from `2K` columns. This is the basis of\n", "dynamic harmonic regression (Hyndman & Athanasopoulos, *FPP3*, ch. 10); tree\n", "ensembles benefit less than linear models, because splits already carve the\n", "calendar, but the columns are cheap and smooth where day-of-year integers are\n", "jagged.\n", "\n", "**Radial basis calendar bumps.** A Gaussian bump centred on an event turns\n", "\"days to Christmas\" into a smooth, bounded response with a chosen width —\n", "narrow for the frantic final week, wide for the season. Unlike a capped linear\n", "clock, the feature is symmetric around the event and dies away on its own.\n", "\n", "**Rolling trend slope.** Differences say *whether* the level moved; a rolling\n", "least-squares slope says *how fast it is moving per day*. For tree models this\n", "matters doubly, because they cannot extrapolate a trend they can only split\n", "on — handing the local slope to the model as a value (section 14's residual\n", "bias is the symptom it treats) is the standard remedy.\n", "\n", "**Rolling quantiles.** The p10–p90 band of the last 28 days describes the\n", "*shape* of recent variation, not just its size: an asymmetric band means the\n", "store's bad days and good days are not mirror images, which a standard\n", "deviation cannot say.\n", "\n", "**Composite-key and redirected grouped statistics.** The\n", "`GroupedRollingFeature` from section 6 generalises twice. A composite key\n", "(weekday × promo state) tracks the trailing level of each *cell* of the\n", "store's week. And the `key` argument redirects reading: a column holding the\n", "*target* period's weekday makes every row carry the weekday profile of the\n", "period it predicts — the general form of the weekly free ride in section 10,\n", "valid at any horizon.\n", "\n", "**Year-over-year alignment.** The ratio of the current 28-day level to the\n", "same window one year earlier isolates annual growth from seasonality. It\n", "demands more than two full years of history (the first year has no\n", "counterpart), and on this store it bought nothing — which is exactly the kind\n", "of thing you want to find out on the holdout, not in production." ] }, { "cell_type": "code", "execution_count": null, "id": "f2b435ea", "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "\n", "year = pd.date_range('2014-01-01', '2014-12-31', freq='D')\n", "doy = year.dayofyear.to_numpy()\n", "\n", "# -- Fourier harmonic stack: 2K columns describing smooth annual shape -------\n", "K = 3\n", "fourier = pd.DataFrame(\n", " {f'annual_{fn}_{k}': getattr(np, fn)(2 * np.pi * k * doy / 365.25)\n", " for k in range(1, K + 1) for fn in ('sin', 'cos')},\n", " index=year,\n", ")\n", "\n", "# -- Radial basis bumps around Christmas, three widths ------------------------\n", "def calendar_bump(index, month, day, width_days):\n", " \"\"\"Gaussian response to the nearest occurrence of one calendar event.\"\"\"\n", " event = pd.to_datetime([f'{y}-{month:02d}-{day:02d}' for y in index.year])\n", " delta = np.abs((index - event).days.to_numpy())\n", " delta = np.minimum(delta, 365 - delta) # wrap around the year boundary\n", " return np.exp(-0.5 * (delta / width_days) ** 2)\n", "\n", "bumps = pd.DataFrame(\n", " {f'christmas_bump_w{w}': calendar_bump(year, 12, 24, w) for w in (5, 15, 30)},\n", " index=year,\n", ")\n", "\n", "fig, axes = plt.subplots(1, 2, figsize=(12, 3.5))\n", "fourier.plot(ax=axes[0], legend=False, linewidth=1)\n", "axes[0].set_title(f'Fourier stack: {K} annual harmonics ({2 * K} columns)')\n", "bumps.plot(ax=axes[1], linewidth=1.5)\n", "axes[1].set_title('Radial basis bumps on Christmas Eve')\n", "for ax in axes:\n", " ax.margins(x=0)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "id": "d0cedc33", "metadata": {}, "outputs": [], "source": [ "sales = observed['Sales'].astype(float)\n", "\n", "# -- Rolling trend slope: local least-squares gradient, in sales per day -----\n", "def rolling_slope(series, window):\n", " \"\"\"Slope of an ordinary least-squares line over each trailing window.\"\"\"\n", " x = np.arange(window) - (window - 1) / 2\n", " denominator = float((x ** 2).sum())\n", " return series.rolling(window).apply(lambda v: float(np.dot(x, v)) / denominator, raw=True)\n", "\n", "slope_28 = rolling_slope(sales, 28)\n", "\n", "# -- Rolling quantile band: the shape of recent variation --------------------\n", "p10 = sales.rolling(28).quantile(0.10)\n", "p90 = sales.rolling(28).quantile(0.90)\n", "\n", "window_view = slice('2014-07-01', '2015-07-31')\n", "fig, axes = plt.subplots(2, 1, figsize=(12, 6), sharex=True)\n", "axes[0].plot(sales.loc[window_view], linewidth=0.6, color='grey', label='daily sales')\n", "axes[0].fill_between(sales.loc[window_view].index, p10.loc[window_view], p90.loc[window_view],\n", " alpha=0.3, label='rolling p10-p90 band')\n", "axes[0].legend(loc='upper left')\n", "axes[0].set_title('Rolling quantile band (28 days)')\n", "axes[1].plot(slope_28.loc[window_view], linewidth=1)\n", "axes[1].axhline(0, color='grey', linewidth=0.5)\n", "axes[1].set_title('Rolling 28-day trend slope (sales per day)')\n", "for ax in axes:\n", " ax.margins(x=0)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": null, "id": "4f86fa41", "metadata": {}, "outputs": [], "source": [ "# -- Composite-key and redirected grouped statistics --------------------------\n", "advanced = observed.copy()\n", "advanced['weekday_promo'] = advanced['DayOfWeek'].astype(str) + '|' + advanced['Promo'].astype(str)\n", "advanced['TargetDayOfWeek'] = (advanced.index + pd.Timedelta(days=3)).dayofweek\n", "\n", "showcase = ForecastFeatureBuilder(\n", " grouped_rolling_features=[\n", " # Trailing mean of each weekday x promo cell of the store's week.\n", " GroupedRollingFeature(column='Sales', by='weekday_promo', window=4,\n", " statistic='mean', min_periods=2),\n", " # The weekday profile of the period three days ahead, read at the\n", " # origin: `by` groups observations by their own weekday, `key`\n", " # redirects each row to the target period's weekday.\n", " GroupedRollingFeature(column='Sales', by='DayOfWeek', window=8,\n", " statistic='mean', key='TargetDayOfWeek'),\n", " ],\n", ")\n", "demo = showcase.transform(ForecastFrame(advanced.reset_index(names='Date'), SCHEMA)).set_index('Date')\n", "\n", "# -- Year-over-year alignment -------------------------------------------------\n", "level_28 = sales.rolling(28, min_periods=14).mean()\n", "yoy_ratio = (level_28 / level_28.shift(364)).replace([np.inf, -np.inf], np.nan)\n", "\n", "summary = pd.DataFrame({\n", " 'weekday_promo_cell_mean': demo['Sales__grouped_mean_by_weekday_promo_lag_1_window_4'],\n", " 'target_weekday_profile': demo['Sales__grouped_mean_by_DayOfWeek_lag_1_window_8_key_TargetDayOfWeek'],\n", " 'yoy_level_ratio': yoy_ratio,\n", "})\n", "display(summary.loc['2015-07-25':'2015-07-31'].round(2))\n", "print('Built and inspected — and deliberately not added to the shipped files:')\n", "print('each of these must first beat the section 13 comparison to earn a column.')" ] }, { "cell_type": "markdown", "id": "f49b5020", "metadata": {}, "source": [ "## What You Built\n", "\n", "Two contract-clean forecasting datasets for one store, and the evidence that\n", "they are trustworthy:\n", "\n", "| Artefact | Property |\n", "| --- | --- |\n", "| `..._daily_base.csv` | Horizon-agnostic. Every feature describes the origin or earlier; valid at any declared horizon. |\n", "| `..._daily_h7.csv` | The identical base columns plus a `t_` block describing the period seven days ahead. Valid only at a seven-day horizon. |\n", "| Leakage audit | Every feature — including the builder's output — rebuilt from truncated history; identical values prove origin-only dependence. |\n", "| Additivity check | The horizon file carries the base layer byte-identically, so the A/B measures one variable. |\n", "| Protocol benchmark | Both files modelled under the certified split arithmetic, judged against naive and seasonal-naive baselines. |\n", "\n", "Wherever these files are trained, the contract travels with them: order by\n", "`date`, keep the observed unshifted target, declare the horizon once at\n", "training time — seven for the `_h7` file, anything for the base file — and\n", "evaluate on a chronological suffix with `horizon − 1` withheld periods. On a\n", "horizon file the declared horizon must match the filename: its `t_` columns\n", "describe one specific future period, and no verification can read your intent\n", "from a column of numbers.\n", "\n", "The layered design is a general recipe, not a Rossmann trick. Build the base\n", "layer first: it is one artefact that serves every horizon, and it sets the\n", "bar. Add a horizon block when the operating plan genuinely drives the target,\n", "and keep it additive so the comparison of section 13 stays a controlled\n", "experiment. Then *measure* whether the block earns its inflexibility — on\n", "stores whose plan is uneventful, the base file can match the tuned one, and\n", "knowing that is worth as much as the uplift elsewhere.\n", "\n", "### Further Reading\n", "\n", "- Hyndman & Athanasopoulos, *Forecasting: Principles and Practice* (3rd ed.) —\n", " benchmark forecasts, evaluation discipline, and dynamic harmonic regression.\n", "- Hyndman & Koehler (2006), *Another look at measures of forecast accuracy* —\n", " the case for MASE.\n", "- Makridakis, Spiliotis & Assimakopoulos (2022), *M5 accuracy competition:\n", " Results, findings, and conclusions*, IJF — gradient-boosted trees on\n", " engineered lag and calendar features as the state of the practice.\n", "\n", "Lesson 07 continues from here: multi-horizon and recursive strategies, panel\n", "models, conformal intervals and hierarchy reconciliation with the SDK's\n", "forecasting toolkit." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "name": "python", "version": "3.11" } }, "nbformat": 4, "nbformat_minor": 5 }