Commit efb7bad1 authored by Paul Rougieux's avatar Paul Rougieux
Browse files

Test TimesFM 3 on agriculture and forestry time series #2

parent 1e4f5914
Loading
Loading
Loading
Loading
+6 −0
Original line number Diff line number Diff line
data/raw/*
data/processed/*
!data/**/.gitkeep
output/*
!output/.gitkeep
__pycache__/
+360 −0
Original line number Diff line number Diff line
---
title: "TimesFM 3.0 on German agricultural and forestry price indices"
date: 2026-09-17
issue: https://aidaho-edu.uni-hohenheim.de/gitlab/modellierung/modellierungs-fitnessclub/-/work_items/2
---

# 1. Objective

Test whether the zero-shot multivariate foundation model **TimesFM 3.0**
(`google/timesfm-3.0-pytorch`, released August 2026) forecasts monthly German
producer price indices for agricultural and forestry products better than
classical econometric models that are re-estimated on the same information set.

"Better" is decided on an out-of-sample rolling-origin backtest, not on a single
train/test cut, because a single 12-month test window on four series gives an
effective sample size close to one.

Three questions from the work item, made operational:

| Question | Operationalisation |
|---|---|
| Which models to compare? | Random walk, seasonal naive, drift, damped-trend ETS, SARIMA, VAR (multivariate) vs. TimesFM 3.0 univariate and TimesFM 3.0 multivariate |
| Which metrics? | MASE, MAE, RMSE, sMAPE for the point forecast; pinball loss and 80 % interval coverage for the predictive distribution; Diebold-Mariano test against the random walk |
| Non-stationarity, seasonality, regime shifts? | Section 5: log transform + differencing for the econometric models, raw levels for TimesFM (it normalises internally); seasonality is explicit in SARIMA/ETS, implicit in TimesFM; regime shifts get their own named-origin stress set |

A secondary objective is methodological: build a backtest harness where any new
model is one function, so the club can re-use it on other datasets in later
sessions.

# 2. Input data

## 2.1 Primary source

Destatis GENESIS-Online, price indices in agriculture and forestry:
<https://www.destatis.de/EN/Themes/Economy/Prices/Price-Indices-In-Agriculture-And-Forestry/_node.html>

| Table | Content | Frequency | Base |
|---|---|---|---|
| `61231` | Producer price index for forestry products — oak, beech, spruce, pine stemwood | monthly | 2020 = 100 |
| `61211` | Producer price index for agricultural products — wheat, barley, sugar beet, potatoes | monthly | 2020 = 100 |

Both are index series, so they are unit-free and directly comparable across
species and crops — convenient for a multivariate model that shares one
normalisation across channels.

## 2.2 Access

GENESIS closed its anonymous `GAST` account. Verified on 2026-09-17: the REST
API at `https://www-genesis.destatis.de/genesisWS/rest/2020/` answers
`helloworld/whoami` without credentials, but `catalogue/tables` and
`data/tablefile` return `401 Code 15`. Registration is free.

`src/get_data.py` therefore supports three input paths:

1. `--source genesis` — REST API, credentials from `GENESIS_USERNAME` /
   `GENESIS_PASSWORD` (or `GENESIS_TOKEN`).
2. `--source csv --csv-file <path>` — a flat CSV (`ffcsv`) downloaded by hand
   from the GENESIS web interface. **No account needed.** This is the fallback
   for the session if nobody has registered.
3. `--source synthetic` — a generated stand-in panel with the same shape,
   seasonality and two deliberate regime shifts. Lets the whole pipeline run
   offline, and is the right thing to develop the harness against.

## 2.3 Structure after ingestion

`data/processed/prices_monthly.csv`, tidy long format:

| column | type | meaning |
|---|---|---|
| `date` | ISO date, first of month | observation month |
| `series_id` | string | e.g. `spruce_stemwood`, `wheat` |
| `group` | string | `forestry` or `agriculture` |
| `value` | float | price index, base 2020 = 100 |

Expected coverage: roughly 2005-01 to the most recent published month
(2026-07 or 2026-08 given the one-to-two month publication lag), i.e. about
250 monthly observations per series. The loader prints the coverage it actually
found; every date in the design below is derived from the last observation `T`,
not hard-coded.

# 3. Output data

Everything lands in `output/`.

| File | Schema |
|---|---|
| `forecasts.csv` | `origin, series_id, model, horizon, target_date, y_true, y_pred, q10 … q90` — one row per (origin, series, model, horizon) |
| `metrics_by_horizon.csv` | `model, series_id, horizon, n, mae, rmse, smape, mase, pinball, coverage_80` |
| `metrics_overall.csv` | `model, series_id, split` aggregated over horizons, plus `dm_stat`, `dm_pvalue` against the random walk |
| `fanchart_<series_id>.png` | last origin, median + 80 % band, TimesFM vs SARIMA vs truth |
| `mase_by_horizon.png` | MASE against horizon, one line per model |
| `summary.md` | the table the group reads out loud in the wrap-up phase |

`forecasts.csv` is the atomic artefact: metrics and plots are pure functions of
it, so a metric can be added after the fact without re-running the model.

## 3.1 Visualising the forecasts against the test set

`src/plots.py` reads `forecasts.csv` and replots without touching a model. Note
that with rolling origins there is **no single out-of-sample path to draw**:
each series carries 36 overlapping 12-month trajectories. That is why there are
several figures rather than one, each answering a different question.

| File | Question it answers |
|---|---|
| `fixed_horizon_h<h>.png` | The literal one. Fix h, and each point comes from a different origin, so it collapses to a genuine out-of-sample series: what each model said h months ahead against what happened. One panel per series. |
| `trajectories_<series_id>.png` | All 36 forecast paths laid over the realised series, one panel per model. Shows *how* a model fails - lagging a turning point looks like a comb of paths peeling off in the wrong direction, which no summary metric conveys. |
| `scatter_pred_actual.png` | Predicted against realised, coloured by horizon. Systematic bias is a cloud off the diagonal; a model that only fails far ahead shows as colour separation. |
| `calibration.png` | Reliability over the 9 quantiles. A curve sitting on one side of the diagonal is bias; a curve *crossing* it is the wrong interval width. |
| `coverage_by_horizon.png` | Where the 80 % interval breaks down as the horizon grows. |
| `error_distribution.png` | Boxplot of absolute scaled error. The mean MASE in the summary hides its own tail, and for a price forecast the worst month is often what matters. |

Read `fixed_horizon` and `trajectories` together: the first says how big the
errors are, the second says what kind of error it is.

# 4. Training set and test set

TimesFM is **zero-shot** — it is never fitted on these data. So "training set"
means two different things and the split must be fair to both:

- for TimesFM: the **context window** handed to the model at an origin;
- for the econometric models: the **estimation sample**, re-fitted from scratch
  at every origin.

Both see exactly the same information: all observations up to and including the
origin month, nothing after. This is enforced in one place, `iter_origins()` in
`src/backtest.py`.

## 4.1 The split, anchored on the last observation `T`

With `T = 2026-07` (adjust automatically if the loader finds a later month):

| Split | Origins | Count | Purpose |
|---|---|---|---|
| **Context / estimation pool** | 2005-01 … origin | expanding | never evaluated, only consumed |
| **Development** | `T-95``T-48` → 2018-08 … 2022-07 | 48 | choose context length, log vs. level, covariate set, SARIMA order |
| **Test (hold-out)** | `T-47``T-12` → 2022-08 … 2025-07 | 36 | reported once, after development is frozen |

Horizon `h = 1 … 12` at every origin. The test window stops at `T-12` so every
origin has a full 12 months of realised values — no ragged edge, no horizon
where the model count differs.

36 origins × 12 horizons × 8 series ≈ 3 456 forecast points per model. Enough
for a Diebold-Mariano test at short horizons; at `h = 12` the overlapping
windows make the errors strongly autocorrelated, so the DM test uses a
Newey-West correction with lag `h-1`.

Rolling origins are **expanding**, not sliding: a 2005 observation still helps a
2025 SARIMA. For TimesFM the context is truncated to the last `--context-length`
months (default 512, well inside the model's 15 360-step limit).

## 4.2 Two extra evaluation sets

**Regime-shift stress set.** Six named origins, each reported on its own rather
than averaged away:

| Origin | Event |
|---|---|
| 2018-06 | onset of the spruce bark-beetle calamity — supply shock, price collapse |
| 2020-03 | COVID-19 |
| 2021-09 | construction-timber price surge |
| 2022-02 | Russian invasion of Ukraine — grain and energy |
| 2022-08 | peak of the energy-price spike |
| 2023-06 | disinflation / mean reversion |

A foundation model that has seen millions of series may extrapolate a shock
better — or may smooth it into a trend. Either result is interesting, and both
are invisible in an average over 36 origins.

**Leakage-aware honest set.** TimesFM 3.0 was released in August 2026 and its
pretraining corpus is not documented at series level. Destatis indices are
public and widely mirrored, so origins before ~2026 may be **contaminated**: the
model may have memorised the realised path. Origins from `T-11` onward
(2026-01 …) with `h = 1 … 6` are almost certainly after the pretraining cutoff.
Few origins, weak statistics, but it is the only honest zero-shot comparison and
it must be reported next to the main result, not instead of it.

This is the single biggest threat to the validity of the whole exercise and
should be stated in the wrap-up.

## 4.3 What counts as a win

TimesFM beats a baseline if its MASE is lower at the same horizon **and** the DM
test rejects equal predictive accuracy at 5 %. On monthly price indices the
random walk is famously hard to beat, so "TimesFM ties the random walk at
`h = 1` and beats seasonal naive at `h = 12`" is a plausible, publishable and
entirely acceptable outcome.

# 5. Non-stationarity, seasonality, regime shifts

| Issue | Econometric models | TimesFM 3.0 |
|---|---|---|
| Non-stationarity | log transform, then `d=1`; unit-root behaviour is assumed rather than tested per origin (an ADF test at every origin would be a pre-test bias) | raw levels; the model applies reversible instance normalisation and optional linear detrending internally (`use_linear_detrending=True`) |
| Seasonality | explicit: `D=1, s=12` in SARIMA, additive seasonal in ETS, and the seasonal naive baseline | implicit; no frequency indicator since 2.5, the model infers periodicity from the context |
| Regime shifts | none of the baselines model breaks; that is the point — they are the honest status quo | unknown behaviour, which is what the stress set measures |

Forecasts from log-space models are converted back with `exp()` on the median,
which is the median of the level — appropriate here, since MASE and pinball loss
are both median-oriented. No smearing correction is applied; if the group wants
mean forecasts this has to change and should be flagged.

# 6. Models

| Key | Type | Specification |
|---|---|---|
| `rw` | baseline | last observed value, carried forward |
| `snaive` | baseline | value 12 months earlier |
| `drift` | baseline | random walk with drift estimated over the whole estimation sample |
| `ets` | econometric | `statsmodels` Holt-Winters, damped additive trend, additive seasonal, on logs |
| `sarima` | econometric | `statsmodels` SARIMAX(1,1,1)(0,1,1)₁₂ on logs |
| `var` | econometric, multivariate | VAR on log first differences, lag order by AIC (max 12) — the fair multivariate comparator |
| `timesfm_uni` | foundation | one call per series, levels, 9 quantiles |
| `timesfm_mv` | foundation | one call per group, all four channels jointly, levels, 9 quantiles |

`timesfm_uni` vs `timesfm_mv` isolates what the multivariate machinery of 3.0
actually buys on these data — the headline feature of the release.

Prediction intervals: TimesFM returns quantiles natively. SARIMA and ETS get
theirs from `get_forecast().conf_int()`. The naive baselines get an empirical
interval from the in-sample residual quantiles at each horizon, so that pinball
loss and coverage are defined for every model.

# 7. Environment

`environment.yml` builds a conda environment `timesfm-fitness` on Python 3.11
(PyTorch wheels are most reliable there). TimesFM itself is pip-only:

```bash
conda env create -f environment.yml
conda activate timesfm-fitness
python -c "import timesfm3, torch; print(torch.__version__)"
```

CPU is enough. The 3.0 checkpoint is ~330 M parameters; a 512-step context with
horizon 12 takes well under a second per series on a laptop, and the whole
backtest is about 300 model calls.

**Licence warning.** The TimesFM *source code* is Apache-2.0, but the **3.0
pretrained weights are released under `timesfm-non-commercial-license-v1.0`**
and are restricted to non-commercial, non-production use. A university workshop
is fine. Anything that ends up in a contract-funded deliverable is not. Weights
up to 2.5 remain Apache-2.0, which is the fallback if this becomes a problem.

# 8. How to run

```bash
conda activate timesfm-fitness
cd events/20260917

# 1. Data. Pick one:
python src/get_data.py --source synthetic                      # works offline, right now
python src/get_data.py --source csv --csv-file data/raw/61231.csv
GENESIS_USERNAME=... GENESIS_PASSWORD=... python src/get_data.py --source genesis

# 2. Backtest. Start without TimesFM to check the harness (seconds):
python src/backtest.py --models rw,snaive,drift,ets,sarima,var --split test

# 3. Add the foundation model (downloads ~1.3 GB on first run):
python src/backtest.py --models all --split test
python src/backtest.py --models all --split stress
python src/backtest.py --models all --split honest --horizon 6
```

```bash
# 4. Replot from the saved forecasts - no model is re-run, takes seconds:
python src/plots.py --split test
python src/plots.py --split test --models timesfm_mv,sarima,rw --horizons 1,3,12
python src/plots.py --split stress --series spruce_stemwood,wheat
```

`--split dev` during the specification phase, `--split test` once and once only.

## 8.1 Measured runtime

The whole pipeline was run end to end on the synthetic panel (259 months,
8 series, 36 origins, h = 1…12) on CPU, Python 3.11, single-threaded BLAS:

| Step | Time |
|---|---|
| `rw`, `snaive`, `drift` | 0.6 s each, both groups |
| `ets` | ~20 s per group |
| `sarima` | ~90 s per group (warm-started from the previous origin) |
| `var` | ~1.5 s per group |
| TimesFM checkpoint load | ~30 s, once |
| `timesfm_uni` + `timesfm_mv` | ~2 min per group |
| **full `--models all --split test`** | **~8 min** |
| `--models baselines --split test` | ~4 min |

So a full re-run fits inside the implementation hour, and the baseline-only run
is fast enough to iterate on. SARIMA dominates the cost: it is 288 maximum-
likelihood fits per group. It is warm-started from the previous origin's
estimates, which roughly triples its speed, and the warm start is discarded
whenever a fit fails to converge so a bad optimum cannot propagate down the
backtest.

Two behaviours already visible on the synthetic panel, both worth watching for
on the real data:

- **ETS is over-confident in the wrong direction**: 80 % coverage came out at
  0.92 against SARIMA's 0.88. That is the `sqrt(h)` interval approximation in
  `forecast_ets` being too generous, exactly as flagged in the code comment.
- **Every model's intervals collapse on the stress set**: coverage fell to
  0.37-0.63 across the six regime-shift origins, against a 0.80 target. If that
  reproduces on the real series it is the most decision-relevant result of the
  session — a forecast interval that fails precisely when prices move is worse
  than no interval at all.

# 9. Open questions for the specification phase

1. **Data access** — does anyone have a GENESIS account, or do we download the
   two CSVs by hand at the start of the session?
2. **Scope** — forestry (61231) only, agriculture (61211) only, or both? Both is
   eight series and a richer multivariate test; one table is faster.
3. **Covariates** — TimesFM 3.0 supports past-only and past-and-future
   covariates. Obvious candidates: energy price index, diesel, fertiliser,
   HICP, harvest volumes. Worth the extra ingestion effort today, or a follow-up
   session?
4. **Horizon** — is 12 months the decision-relevant horizon for the people in
   the room, or is 1-3 months (contract negotiation) what actually matters?
5. **Leakage** — do we accept the contaminated main result with a caveat, or
   make the post-cutoff honest set the headline?
6. **Higher frequency** — the work item welcomes higher-frequency data. Weekly
   timber auction or grain spot prices would make the foundation model much more
   interesting, but needs a source.

# 10. Status

The pipeline runs end to end on the synthetic panel: all three splits
(`test`, `stress`, `honest`) produce forecasts, metrics, plots and a summary.
The test split yields exactly 8 models x 8 series x 36 origins x 12 horizons =
27 648 forecast rows.

The TimesFM 3.0 calls are verified by execution, not only by reading the
package: `TimesFM3Forecaster.from_pretrained("google/timesfm-3.0-pytorch",
device="cpu")` loads in ~30 s, univariate `predict_batch` returns
`forecast (12,)` and `quantiles (12, 9)`, multivariate returns
`forecast (k, 12)` and `quantiles (k, 12, 9)` — the shapes the harness assumes.

**Not yet verified:** `parse_genesis_csv` has never seen a real GENESIS export.
The column sniffing is deliberately tolerant and prints what it found, but
budget five minutes at the start of the session to adjust it.

**Meaningless so far:** every number in `output/` comes from the synthetic
panel. It demonstrates that the harness works, and nothing about German price
formation.

# 11. Files

```
events/20260917/
├── README.md              this document
├── environment.yml        conda environment
├── src/
│   ├── get_data.py        GENESIS API / CSV / synthetic → tidy panel
│   ├── metrics.py         MASE, pinball, coverage, Diebold-Mariano
│   ├── backtest.py        rolling-origin harness, all models
│   └── plots.py           replots forecasts.csv, no model re-run
├── data/                  gitignored, except this structure
└── output/                gitignored
```
+0 −0

Empty file added.

+0 −0

Empty file added.

+53 −0
Original line number Diff line number Diff line
# Modellierungs-Fitnessclub, session 2026-09-17
# TimesFM 3.0 on German agricultural and forestry producer price indices
#
#   conda env create -f environment.yml
#   conda activate timesfm-fitness
#
# Python 3.11: the version with the least friction for PyTorch wheels.
name: timesfm-fitness

channels:
  - conda-forge

dependencies:
  - python=3.11

  # data handling
  - numpy>=1.26.4
  - pandas>=2.2
  - pyarrow          # parquet, and a faster csv reader
  - requests         # GENESIS REST API

  # econometric baselines
  - statsmodels>=0.14   # SARIMAX, ExponentialSmoothing, VAR
  - scipy>=1.11
  - scikit-learn>=1.4

  # plots
  - matplotlib>=3.8

  # notebooks, optional but handy in a workshop
  - jupyterlab
  - ipykernel

  - pip>=24
  - pip:
      # TimesFM is not on conda-forge. The [torch] extra pulls a CPU-capable
      # torch wheel; on a machine with CUDA, install torch from the PyTorch
      # index first and pip will keep it.
      - timesfm[torch]==3.0.2

      # checkpoint download, pinned by timesfm but listed for clarity
      - huggingface_hub>=0.28.0
      - safetensors>=0.5.3

# Checkpoint (~1.3 GB, cached in ~/.cache/huggingface after the first run):
#     google/timesfm-3.0-pytorch
#
# Licence: the TimesFM source is Apache-2.0, but the 3.0 *weights* are under
# timesfm-non-commercial-license-v1.0 — non-commercial, non-production use only.
# Weights up to 2.5 are Apache-2.0 if that ever becomes a constraint.
#
# Smoke test after creating the environment:
#     python -c "import torch, statsmodels, timesfm3; print(torch.__version__)"
Loading