backend/scripts/backtest_rate_sensitivity.py

@@ -99,6 +99,30 @@ CAVEATS (read before trusting the numbers)
       confounded-window flag, the Z-bucket phrase). It answers one question:
       does the fitted slope predict direction out-of-sample?
 
+CANDIDATE-METHOD VARIANTS (#978 Almon / #979 deseasonalize)
+-----------------------------------------------------------
+Two forecast modules were shipped advisory-only and NOT wired into prod, pending
+a backtest proving they recover OOS signal the production ``best_lag``-on-raw
+path missed. This harness evaluates BOTH as opt-in variants, each isolated
+against the always-on raw reference (we never explode into all combinations):
+
+  • ``--deseasonalize`` (#979 ``normalize.py``) — divide month-of-year seasonal
+    factors out of the units series BEFORE ``log_diff`` (cleaner regressand),
+    then score with the SAME ``best_lag`` engine. Factors are fit on TRAIN
+    months ONLY and applied point-in-time onto the test months (a test-month
+    observation must NEVER influence the factors — same leakage discipline as
+    ``--detrend``).
+  • ``--almon`` (#978 ``regression.py``) — replace ``best_lag``'s single winning
+    lag with an Almon polynomial distributed-lag estimator (a smooth β curve over
+    lags 0..6, reported via its long-run multiplier). A NEW OOS evaluator
+    (``evaluate_oos_almon``) fits on TRAIN only and predicts each test month
+    ``Σ_j β_j·Δrate[t−j]`` strictly point-in-time.
+
+A NEGATIVE result is a valid, honest outcome — you cannot extract signal that
+isn't there. The failure mode this harness guards against is a leaky evaluator
+that FALSELY shows signal, so both methods' TRAIN/TEST boundaries are pinned and
+unit-tested on synthetic series with known answers.
+
 USAGE
 -----
     DATABASE_URL=postgresql+psycopg://... \
@@ -109,6 +133,9 @@ USAGE
 
     # the #978b cross-checks: survivorship-free Source A + a detrend control:
     python -m scripts.backtest_rate_sensitivity --source both --detrend
+
+    # evaluate the candidate methods OOS (#978 Almon-ADL + #979 deseasonalize):
+    python -m scripts.backtest_rate_sensitivity --source B --almon --deseasonalize
 """
 
 from __future__ import annotations
@@ -181,6 +208,13 @@ _SOURCES: tuple[str, ...] = (_SOURCE_B, _SOURCE_A)
 # (the difference step then behaves exactly like the raw log_diff path).
 _DETREND_MIN_POINTS: int = 3
 
+# Estimator selector — which OOS evaluator backtest_tier dispatches to.
+#   best_lag = single-winning-lag OLS (the production §9.6 core, evaluate_oos).
+#   almon    = #978 Almon polynomial distributed-lag (evaluate_oos_almon).
+# Default is best_lag so existing variants/tests are byte-identical (back-compat).
+_ESTIMATOR_BEST_LAG: str = "best_lag"
+_ESTIMATOR_ALMON: str = "almon"
+
 
 def _import_engine() -> tuple[Any, Any, Any]:
     """Lazy import of the §9.6 engine's pure funcs + Δln helper.
@@ -214,6 +248,50 @@ def _import_lags() -> tuple[int, ...]:
     return _LAGS
 
 
+def _import_regression() -> tuple[Any, int]:
+    """Lazy import of the #978 Almon distributed-lag estimator + its lag window.
+
+    Returns ``(fit_almon_dl, _MAX_LAG)``. Deferred for the SAME reason as
+    ``_import_engine`` — app.services.forecasting.regression pulls
+    app.core.config.Settings (fail-fasts without DATABASE_URL), and the
+    pure-logic unit tests must import this module cheaply. ``_MAX_LAG`` (the
+    regression module's distributed-lag window upper bound, default 6) is read
+    so the Almon evaluator's point-in-time prediction iterates the SAME lag span
+    the estimator fitted.
+    """
+    try:
+        from app.services.forecasting.regression import _MAX_LAG, fit_almon_dl
+    except ImportError:  # pragma: no cover — fallback for adhoc invocation
+        import sys
+
+        sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+        from app.services.forecasting.regression import _MAX_LAG, fit_almon_dl
+    return fit_almon_dl, _MAX_LAG
+
+
+def _import_normalize() -> tuple[Any, Any]:
+    """Lazy import of the #979 seasonal deseasonalization helpers.
+
+    Returns ``(seasonal_factors, deseasonalize_values)`` from
+    app.services.forecasting.normalize. Deferred for the SAME reason as
+    ``_import_engine`` (the normalize module pulls SalesSeries → Settings).
+    """
+    try:
+        from app.services.forecasting.normalize import (
+            deseasonalize_values,
+            seasonal_factors,
+        )
+    except ImportError:  # pragma: no cover — fallback for adhoc invocation
+        import sys
+
+        sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+        from app.services.forecasting.normalize import (
+            deseasonalize_values,
+            seasonal_factors,
+        )
+    return seasonal_factors, deseasonalize_values
+
+
 def _session() -> Session:
     """Lazy SessionLocal factory — see _import_engine for why it's deferred."""
     try:
@@ -255,12 +333,16 @@ class TierResult:
     full_sample_lag: int | None
     lag_stable: bool
     skipped: str | None
+    deseasonalized: bool = False
+    estimator: str = _ESTIMATOR_BEST_LAG
 
     def as_dict(self) -> dict[str, Any]:
         return {
             "tier": self.tier,
             "source": self.source,
             "detrended": self.detrended,
+            "deseasonalized": self.deseasonalized,
+            "estimator": self.estimator,
             "n_aligned": self.n_aligned,
             "n_train": self.n_train,
             "n_test": self.n_test,
@@ -500,6 +582,136 @@ def evaluate_oos(
     }
 
 
+def evaluate_oos_almon(
+    delta_sales: list[float | None],
+    rate_deltas: list[float | None],
+    *,
+    holdout_frac: float = _HOLDOUT_FRAC,
+) -> dict[str, Any]:
+    """Time-ordered OOS backtest of the #978 Almon distributed-lag fit. PURE (no DB).
+
+    Mirrors ``evaluate_oos``'s contract and return-dict shape exactly so
+    ``backtest_tier`` can wrap it identically — but the estimator is the Almon
+    polynomial distributed-lag model (``regression.fit_almon_dl``) instead of the
+    single-winning-lag OLS (``rate_sensitivity.best_lag``). The Almon model fits
+    ONE smooth coefficient curve β_0..β_K over the whole 0..max_lag window, so a
+    TEST month's prediction sums the contributions of ALL lags, not just one.
+
+    Steps:
+      1. Split the aligned months time-ordered (same boundary ``evaluate_oos``
+         uses): fit on the oldest ``holdout_frac``, test on the newest remainder.
+      2. ``fit_almon_dl(rate_deltas[:n_train], delta_sales[:n_train])`` on the
+         TRAIN slice only — NB the regression module's signature is
+         ``fit_almon_dl(x, y)`` with ``x`` the regressor (Δrate) and ``y`` the
+         regressand (Δln(sales)), so Δrate is the FIRST arg. None (math-infeasible
+         / too thin) → empty result (nothing to validate). The fit's ``best_lag``
+         (peak-|β_j| lag) is the "train_lag"; its ``long_run_coef`` (Σ_j β_j) is
+         reported as "train_beta" (the long-run multiplier); its ``r2`` is the
+         in-sample R².
+      3. For each TEST month at absolute index t, predict
+         ``Σ_{j=0..max_lag} per_lag[j]·Δrate[t−j]`` STRICTLY point-in-time: each
+         lagged regressor is built over the FULL aligned series via
+         ``_shift_for_lag(rate_deltas, j)`` then read at index t, so a test month
+         only ever reads Δrate at indices ≤ t (never the future). Score the SIGN
+         vs the actual Δln(sales); skip a month when the actual OR any required
+         lag value is None (don't fabricate a partial lag profile).
+      4. Refit ``fit_almon_dl`` on the FULL aligned series → full-sample peak lag
+         (for the verdict's lag-stability check).
+
+    Returns the SAME dict keys as ``evaluate_oos`` (n_aligned, n_train, n_test,
+    train_lag, train_beta [=long-run Σβ here], in_sample_r2, oos_hit_rate,
+    oos_signed_mae, full_sample_lag, lag_stable). The in-sample R² is high by
+    construction (honesty block applies identically); only oos_hit_rate is
+    trustworthy.
+    """
+    fit_almon_dl, max_lag = _import_regression()
+
+    n = len(delta_sales)
+    n_train = _time_ordered_split(n, holdout_frac)
+    n_test = n - n_train
+
+    empty: dict[str, Any] = {
+        "n_aligned": n,
+        "n_train": n_train,
+        "n_test": n_test,
+        "train_lag": None,
+        "train_beta": None,
+        "in_sample_r2": None,
+        "oos_hit_rate": None,
+        "oos_signed_mae": None,
+        "full_sample_lag": None,
+        "lag_stable": False,
+    }
+    if n_train < 2 or n_test < 1:
+        return empty
+
+    train_sales = delta_sales[:n_train]
+    train_rate = rate_deltas[:n_train]
+    # fit_almon_dl(x, y): x = Δrate (regressor) FIRST, y = Δln(sales) (regressand).
+    fit = fit_almon_dl(train_rate, train_sales, max_lag=max_lag)
+    if fit is None:
+        # Almon fit math-infeasible on TRAIN → nothing to validate here.
+        return empty
+
+    per_lag = fit["per_lag_coef"]  # tuple β_0..β_max_lag from the Almon form
+    train_lag = int(fit["best_lag"])
+    train_beta = float(fit["long_run_coef"])  # long-run Σβ reported as the "beta" field
+    in_sample_r2 = float(fit["r2"]) if fit["r2"] is not None else None
+
+    # Build one point-in-time-safe shifted view per lag over the FULL aligned
+    # series. shifted_by_lag[j][t] == rate_deltas[t-j] (None where t-j < 0), so a
+    # test month at index t reads only Δrate at indices ≤ t — never the future.
+    shifted_by_lag = [_shift_for_lag(rate_deltas, j) for j in range(max_lag + 1)]
+
+    hits = 0
+    scored = 0
+    abs_err_sum = 0.0
+    for t in range(n_train, n):
+        actual = delta_sales[t]
+        if actual is None:
+            continue
+        # The Almon prediction needs the FULL lag profile at t; a single missing
+        # lag means we can't form Σ per_lag[j]·Δrate[t-j] → skip (no fabrication).
+        lag_vals: list[float] = []
+        complete = True
+        for j in range(max_lag + 1):
+            xv = shifted_by_lag[j][t]
+            if xv is None:
+                complete = False
+                break
+            lag_vals.append(float(xv))
+        if not complete:
+            continue
+        predicted = sum(per_lag[j] * lag_vals[j] for j in range(max_lag + 1))
+        scored += 1
+        abs_err_sum += abs(predicted - float(actual))
+        # Directional hit: same nonzero sign (a flat 0.0 actual can't be
+        # directionally predicted, so it never counts as a hit).
+        if (predicted > 0 and actual > 0) or (predicted < 0 and actual < 0):
+            hits += 1
+
+    oos_hit_rate = (hits / scored) if scored > 0 else None
+    oos_signed_mae = (abs_err_sum / scored) if scored > 0 else None
+
+    # Full-sample refit (x=Δrate first, y=Δln(sales)) for the lag-stability check.
+    full_fit = fit_almon_dl(rate_deltas, delta_sales, max_lag=max_lag)
+    full_lag = int(full_fit["best_lag"]) if full_fit is not None else None
+    lag_stable = full_lag is not None and full_lag == train_lag
+
+    return {
+        "n_aligned": n,
+        "n_train": n_train,
+        "n_test": scored,  # report the number actually SCORED, not the raw span
+        "train_lag": train_lag,
+        "train_beta": train_beta,
+        "in_sample_r2": in_sample_r2,
+        "oos_hit_rate": oos_hit_rate,
+        "oos_signed_mae": oos_signed_mae,
+        "full_sample_lag": full_lag,
+        "lag_stable": lag_stable,
+    }
+
+
 def align_series(
     sales_by_month: dict[date, int],
     rate_by_month: dict[date, float],
@@ -542,6 +754,34 @@ def _delta_sales_series(
     return _rate_first_diff(_detrend_log(units, fit_n=fit_n))
 
 
+def _deseasonalize_units(months: list[date], units: list[int], *, fit_n: int) -> list[float | None]:
+    """Deseasonalize a units series → Δln(deseasonalized units). PURE (deferred import).
+
+    The #979 control. Mirrors ``_delta_sales_series``'s train-only discipline:
+      1. Fit month-of-year seasonal factors on the TRAIN months ONLY
+         (``seasonal_factors(months[:fit_n], units[:fit_n])``). A test-month
+         observation must NEVER influence the factors (look-ahead leakage rule).
+      2. Deseasonalize the FULL series with the TRAIN-fit factors
+         (``deseasonalize_values(months, units, factors)`` → raw_t / factor[m]).
+      3. ``log_diff`` the deseasonalized units → Δln regressand (the same Y-axis
+         ``evaluate_oos`` scores).
+
+    ``log_diff`` itself is point-in-time (a first difference reads only t and
+    t−1), so the ONLY leakage risk is the factor fit — pinned to the TRAIN slice
+    here. The deseasonalized values are floats (raw / factor); ``log_diff`` is
+    float-math throughout, so no int-rounding is applied (unlike the prod
+    ``normalize_demand`` wrapper, which rounds only to honour SalesSeries.units's
+    int contract — irrelevant for the regressand).
+    """
+    seasonal_factors, deseasonalize_values = _import_normalize()
+    _bl, _ols, log_diff = _import_engine()
+    train_months = months[:fit_n]
+    train_units = units[:fit_n]
+    adjustment = seasonal_factors(train_months, train_units)
+    deseasonalized = deseasonalize_values(months, units, adjustment.factors)
+    return log_diff(deseasonalized)
+
+
 def backtest_tier(
     sales_by_month: dict[date, int],
     rate_by_month: dict[date, float],
@@ -549,19 +789,28 @@ def backtest_tier(
     tier: str,
     source: str = _SOURCE_B,
     detrend: bool = False,
+    deseasonalize: bool = False,
+    estimator: str = _ESTIMATOR_BEST_LAG,
     holdout_frac: float = _HOLDOUT_FRAC,
     min_months: int = _MIN_BACKTEST_MONTHS,
 ) -> TierResult:
     """Build Δ-series for one tier, run the OOS backtest, wrap as TierResult.
 
     Aligns the tier's monthly sold-units to the monthly key_rate, computes the
-    regressand (``log_diff`` raw, or Δ of linear-detrended ``ln`` when
+    regressand and Δrate (first diff), then delegates to the OOS evaluator named
-    ``detrend`` — see ``_delta_sales_series``) and Δrate (first diff), then
+    by ``estimator`` (``evaluate_oos`` for ``best_lag``, ``evaluate_oos_almon``
-    delegates to ``evaluate_oos``. ``source`` (B/A) and ``detrend`` are recorded
+    for ``almon``). The regressand is one of three PREPROCESSING variants, all
-    on the result for labelling, not used in the math here. Tiers with fewer than
+    fit point-in-time on TRAIN months only (no look-ahead leakage):
-    ``min_months`` aligned months are SKIPPED (TierResult with ``skipped`` set,
+      • plain ``log_diff(units)`` (default),
-    all metrics None) — no silent drop. PURE aside from the deferred engine
+      • Δ of linear-detrended ``ln(units)`` when ``detrend`` (#978b control), or
-    import.
+      • ``log_diff`` of month-of-year deseasonalized units when ``deseasonalize``
+        (#979 control — see ``_deseasonalize_units``).
+    ``detrend`` and ``deseasonalize`` are mutually exclusive preprocessing paths;
+    if both are set ``deseasonalize`` takes precedence (the planner never emits
+    that combination). ``source``/``detrend``/``deseasonalize``/``estimator`` are
+    recorded on the result for labelling. Tiers with fewer than ``min_months``
+    aligned months are SKIPPED (TierResult with ``skipped`` set, all metrics
+    None) — no silent drop. PURE aside from the deferred engine import.
     """
     months, units, rates = align_series(sales_by_month, rate_by_month)
     n_aligned = len(months)
@@ -570,6 +819,8 @@ def backtest_tier(
             tier=tier,
             source=source,
             detrended=detrend,
+            deseasonalized=deseasonalize,
+            estimator=estimator,
             n_aligned=n_aligned,
             n_train=0,
             n_test=0,
@@ -583,21 +834,31 @@ def backtest_tier(
             skipped=f"only {n_aligned} aligned months (< {min_months})",
         )
 
-    # Detrend (when enabled) must be fit on TRAIN months ONLY, then projected
+    # All preprocessing that ESTIMATES anything (the detrend trend, the seasonal
-    # point-in-time onto the test months — otherwise the held-out TEST data
+    # factors) must be fit on TRAIN months ONLY, then projected point-in-time
-    # shapes the trend and the OOS hit-rate is inflated by look-ahead leakage
+    # onto the test months — otherwise the held-out TEST data shapes the
-    # (#978 Part A). We compute the SAME train boundary evaluate_oos will use
+    # preprocessing and the OOS hit-rate is inflated by look-ahead leakage (#978
-    # (len(delta_sales) == len(units) == n_aligned, so the split index matches)
+    # Part A for detrend; the #979 leakage rule for deseasonalize). We compute the
-    # and pass it as the detrend fit window.
+    # SAME train boundary the evaluator will use (len(delta_sales) == len(units)
+    # == n_aligned, so the split index matches) and pass it as the fit window.
     n_train = _time_ordered_split(n_aligned, holdout_frac)
-    delta_sales = _delta_sales_series(units, detrend=detrend, fit_n=n_train)
+    if deseasonalize:
+        delta_sales = _deseasonalize_units(months, units, fit_n=n_train)
+    else:
+        delta_sales = _delta_sales_series(units, detrend=detrend, fit_n=n_train)
     rate_deltas = _rate_first_diff([float(r) for r in rates])
-    res = evaluate_oos(delta_sales, rate_deltas, holdout_frac=holdout_frac)
+
+    if estimator == _ESTIMATOR_ALMON:
+        res = evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=holdout_frac)
+    else:
+        res = evaluate_oos(delta_sales, rate_deltas, holdout_frac=holdout_frac)
 
     return TierResult(
         tier=tier,
         source=source,
         detrended=detrend,
+        deseasonalized=deseasonalize,
+        estimator=estimator,
         n_aligned=res["n_aligned"],
         n_train=res["n_train"],
         n_test=res["n_test"],
@@ -895,6 +1156,8 @@ def run_backtest(
     district: str | None,
     source: str = _SOURCE_B,
     detrend: bool = False,
+    deseasonalize: bool = False,
+    estimator: str = _ESTIMATOR_BEST_LAG,
     rate_by_month: dict[date, float] | None = None,
 ) -> dict[str, Any]:
     """Drive ONE source/variant of the read-only backtest → results dict. No writes.
@@ -902,12 +1165,14 @@ def run_backtest(
     Loads the monthly key_rate (or reuses ``rate_by_month`` when the caller has
     already loaded it once for several variants), then the EKB-wide and per-class
     sold-units series for ``source``, backtests each tier (``backtest_tier``,
-    with ``detrend`` applied), and assembles the per-source verdict + per-tier
+    with the ``detrend`` / ``deseasonalize`` preprocessing and the ``estimator``
-    OOS lifts.
+    applied), and assembles the per-source verdict + per-tier OOS lifts.
 
     ``classes`` None → auto-discover every class present in the chosen source; an
     empty list → EKB-wide only. ``district`` narrows ALL tiers for Source B only
-    (ignored for Source A).
+    (ignored for Source A). ``deseasonalize`` (#979) and ``estimator`` ("almon"
+    selects the #978 Almon distributed-lag evaluator) are recorded on the run
+    dict + params alongside ``detrend`` so labels/JSON carry the full variant.
     """
     if rate_by_month is None:
         rate_by_month = load_rate_by_month(db, since=since)
@@ -928,12 +1193,17 @@ def run_backtest(
         tier=_EKB_WIDE,
         source=source,
         detrend=detrend,
+        deseasonalize=deseasonalize,
+        estimator=estimator,
         holdout_frac=holdout_frac,
     )
     logger.info(
-        "source=%s detrend=%s EKB-wide: aligned=%d train=%d test=%d lag=%s hit_rate=%s",
+        "source=%s detrend=%s deseasonalize=%s estimator=%s EKB-wide: "
+        "aligned=%d train=%d test=%d lag=%s hit_rate=%s",
         source,
         detrend,
+        deseasonalize,
+        estimator,
         ekb.n_aligned,
         ekb.n_train,
         ekb.n_test,
@@ -953,6 +1223,8 @@ def run_backtest(
             tier=cls,
             source=source,
             detrend=detrend,
+            deseasonalize=deseasonalize,
+            estimator=estimator,
             holdout_frac=holdout_frac,
         )
         tiers.append(res)
@@ -972,6 +1244,8 @@ def run_backtest(
     return {
         "source": source,
         "detrended": detrend,
+        "deseasonalized": deseasonalize,
+        "estimator": estimator,
         "a_district_ignored": a_district_ignored,
         "params": {
             "since": since,
@@ -980,6 +1254,8 @@ def run_backtest(
             "classes": classes,
             "source": source,
             "detrended": detrend,
+            "deseasonalized": deseasonalize,
+            "estimator": estimator,
             "min_backtest_months": _MIN_BACKTEST_MONTHS,
             "lags": list(_import_lags()),
         },
@@ -992,24 +1268,70 @@ def run_backtest(
     }
 
 
-def _variant_label(source: str, detrend: bool) -> str:
+def _variant_label(
-    """Human label for a (source, detrend) run, e.g. 'B raw' / 'B detrended' / 'A raw'."""
+    source: str,
-    return f"{source} {'detrended' if detrend else 'raw'}"
+    detrend: bool,
+    *,
+    deseasonalize: bool = False,
+    estimator: str = _ESTIMATOR_BEST_LAG,
+) -> str:
+    """Human label for a variant run. PURE.
 
-
+    The method descriptor is mutually-exclusive (the planner never combines
-def _plan_variants(sources: list[str], detrend: bool) -> list[tuple[str, bool]]:
+    them), so the label names the ONE active method:
-    """Which (source, detrend) variants to run, in report order. PURE.
+      • ``estimator == 'almon'``   → ``f"{source} Almon-ADL"`` (the #978 estimator),
-
+      • ``deseasonalize``          → ``f"{source} deseasonalized"`` (#979 preprocessing),
-    For each requested source we always run the RAW variant (the reference). When
+      • ``detrend``                → ``f"{source} detrended"`` (#978b control),
-    ``--detrend`` is set we ALSO run the detrended variant of that source, so a
+      • otherwise (best_lag, raw)  → ``f"{source} raw"`` (the always-on reference).
-    single invocation can show ``B raw`` next to ``B detrended`` (the survivorship
-    control) for the verdict's side-by-side comparison.
     """
-    variants: list[tuple[str, bool]] = []
+    if estimator == _ESTIMATOR_ALMON:
+        return f"{source} Almon-ADL"
+    if deseasonalize:
+        return f"{source} deseasonalized"
+    if detrend:
+        return f"{source} detrended"
+    return f"{source} raw"
+
+
+def _variant_label_for_run(run: dict[str, Any]) -> str:
+    """Label a run dict via ``_variant_label`` (reads source/detrend/deseason/estim)."""
+    return _variant_label(
+        run["source"],
+        run["detrended"],
+        deseasonalize=run.get("deseasonalized", False),
+        estimator=run.get("estimator", _ESTIMATOR_BEST_LAG),
+    )
+
+
+def _plan_variants(
+    sources: list[str],
+    detrend: bool,
+    *,
+    deseasonalize: bool = False,
+    almon: bool = False,
+) -> list[tuple[str, bool, bool, str]]:
+    """Which variants to run, in report order. PURE.
+
+    Each entry is ``(source, detrend, deseasonalize, estimator)``. For each
+    requested source we ALWAYS run the RAW reference (best_lag on raw units),
+    then ADD only the explicitly-requested method variants — we do NOT explode
+    into all combinations:
+      • ``--detrend``      → ``(src, True,  False, best_lag)`` (#978b survivorship control),
+      • ``--deseasonalize``→ ``(src, False, True,  best_lag)`` (#979 seasonal preprocessing),
+      • ``--almon``        → ``(src, False, False, almon)``    (#978 Almon distributed-lag).
+    The method flags are independent: passing several adds several variants per
+    source (raw + each requested method), each isolating ONE method against the
+    raw reference for the verdict's side-by-side comparison.
+    """
+    variants: list[tuple[str, bool, bool, str]] = []
     for src in sources:
-        variants.append((src, False))
+        variants.append((src, False, False, _ESTIMATOR_BEST_LAG))
         if detrend:
-            variants.append((src, True))
+            variants.append((src, True, False, _ESTIMATOR_BEST_LAG))
+        if deseasonalize:
+            variants.append((src, False, True, _ESTIMATOR_BEST_LAG))
+        if almon:
+            variants.append((src, False, False, _ESTIMATOR_ALMON))
     return variants
 
 
@@ -1019,19 +1341,25 @@ def cross_source_verdict(
     margin: float = _VERDICT_HITRATE_MARGIN,
     min_months: int = _MIN_BACKTEST_MONTHS,
 ) -> dict[str, Any]:
-    """Compare the EKB-wide OOS verdict across variants (B raw / B detrended / A).
+    """Compare the EKB-wide OOS verdict across every method variant. PURE.
 
-    The #978b question: is Source B's negative OOS verdict a SURVIVORSHIP
+    The §9.6 question, generalised across the variants present: does ANY
-    ARTIFACT or a real ``no signal``? We line up each variant's EKB-wide
+    method/source recover OOS directional signal the raw best_lag baseline
-    OOS hit-rate vs the 0.5 coin-flip baseline and synthesise a conclusion:
+    missed? Variants can be the survivorship controls (B detrended, A) AND the
+    two new candidate methods — ``deseasonalized`` (#979 month-of-year seasonal
+    preprocessing) and ``Almon-ADL`` (#978 polynomial distributed-lag estimator).
+    We line up each variant's EKB-wide OOS hit-rate vs the 0.5 coin-flip baseline
+    (decision rule UNCHANGED: a variant "beats" only if scorable AND
+    ``oos_hit_rate >= 0.5 + margin`` AND ``lag_stable``) and synthesise:
 
-      • If NO variant (B raw, B detrended, or survivorship-free A) clears
+      • If NO variant clears coin-flip+margin → the negative verdict is
-        coin-flip+margin → the negative verdict is corroborated as a real
+        corroborated as a real ``no signal``, not an artifact and not something
-        ``no signal``, not an artifact (the detrend + survivorship-free controls
+        the new methods rescue — neither deseasonalizing the input nor the
-        agree). Source A's thin-data caveat is attached when A drove a verdict.
+        smoother Almon estimator extracts signal that isn't there. Source A's
-      • If the detrended or survivorship-free variant DOES clear the bar while
+        thin-data caveat is attached when an A row drove the comparison.
-        raw B did not → the raw verdict may have been a survivorship artifact;
+      • If a control OR a candidate method DOES clear the bar while raw best_lag
-        flag the variant that shows signal.
+        did not → flag that variant: the raw verdict may be a survivorship
+        artifact, or the method recovered signal worth promoting from advisory.
 
     PURE — operates on already-computed run dicts. Returns a dict with a
     ``lines`` list (rendered as-is) plus structured fields for JSON.
@@ -1041,7 +1369,7 @@ def cross_source_verdict(
     thin_variants: list[str] = []
     for run in runs:
         ekb: TierResult = run["ekb_result"]
-        label = _variant_label(run["source"], run["detrended"])
+        label = _variant_label_for_run(run)
         hr = ekb.oos_hit_rate
         scorable = ekb.skipped is None and hr is not None and ekb.n_test >= 1
         beats = bool(scorable and hr is not None and hr >= 0.5 + margin and ekb.lag_stable)
@@ -1055,6 +1383,8 @@ def cross_source_verdict(
                 "variant": label,
                 "source": run["source"],
                 "detrended": run["detrended"],
+                "deseasonalized": run.get("deseasonalized", False),
+                "estimator": run.get("estimator", _ESTIMATOR_BEST_LAG),
                 "scorable": scorable,
                 "oos_hit_rate": _round_or_none(hr, 4),
                 "n_test": ekb.n_test,
@@ -1064,16 +1394,24 @@ def cross_source_verdict(
             }
         )
 
+    # Width the label column to the widest variant label present (Almon-ADL /
+    # deseasonalized are longer than the original "B detrended").
+    label_w = max((len(r["variant"]) for r in rows), default=13)
+    label_w = max(label_w, 13)
+
     lines: list[str] = []
-    lines.append("CROSS-SOURCE VERDICT (B raw vs B detrended vs A — #978b):")
+    lines.append(
+        "CROSS-SOURCE / CROSS-METHOD VERDICT (raw best_lag vs controls "
+        "[detrended, A] vs candidate methods [deseasonalized #979, Almon-ADL #978]):"
+    )
     for r in rows:
         if not r["scorable"]:
             why = r["skipped"] or "no gated lag / empty test window"
-            lines.append(f"  {r['variant']:<13} → not scorable ({why})")
+            lines.append(f"  {r['variant']:<{label_w}} → not scorable ({why})")
         else:
             tag = "SIGNAL > coin-flip" if r["beats_coin"] else "no signal (≤ coin-flip)"
             lines.append(
-                f"  {r['variant']:<13} → OOS_hit={_fmt_rate(r['oos_hit_rate'])} "
+                f"  {r['variant']:<{label_w}} → OOS_hit={_fmt_rate(r['oos_hit_rate'])} "
                 f"(n_test={r['n_test']}, lag_stable={'yes' if r['lag_stable'] else 'no'}) "
                 f"→ {tag}"
             )
@@ -1082,16 +1420,19 @@ def cross_source_verdict(
         conclusion = (
             "CONCLUSION: OOS signal above coin-flip appears in: "
             + ", ".join(signal_variants)
-            + ". The §9.6 negative verdict on raw Source B may be a SURVIVORSHIP "
+            + ". The §9.6 negative verdict on raw best_lag may be a SURVIVORSHIP "
-            "ARTIFACT — the control(s) above recover directional signal."
+            "ARTIFACT, or a candidate method (deseasonalize / Almon-ADL) recovers "
+            "directional signal worth promoting from advisory — inspect the flagged "
+            "variant(s) above."
         )
         promote_any = True
     else:
         conclusion = (
-            "CONCLUSION: NO variant (raw B, detrended B, or survivorship-free A) "
+            "CONCLUSION: NO variant — neither the survivorship controls (detrended B, "
-            "beats coin-flip+margin out-of-sample. The §9.6 negative verdict is a "
+            "survivorship-free A) NOR the candidate methods (#979 deseasonalize, #978 "
-            "REAL 'no signal', NOT a survivorship artifact — detrending B and the "
+            "Almon-ADL) — beats coin-flip+margin out-of-sample. The §9.6 negative "
-            "survivorship-free Source A both agree. Keep advisory."
+            "verdict is a REAL 'no signal', NOT a survivorship artifact, and the new "
+            "methods do not recover signal that isn't there. Keep advisory."
         )
         promote_any = False
     lines.append("  " + conclusion)
@@ -1125,19 +1466,23 @@ def run_all(
     district: str | None,
     sources: list[str],
     detrend: bool,
+    deseasonalize: bool = False,
+    almon: bool = False,
 ) -> dict[str, Any]:
-    """Run every requested (source, detrend) variant + the cross-source verdict.
+    """Run every requested variant + the cross-source/method verdict. No writes.
 
     Loads the monthly key_rate ONCE and reuses it across variants. ``sources`` is
-    a subset of (B, A); ``detrend`` adds the detrended variant of each. No
+    a subset of (B, A); the always-on RAW reference runs per source, plus one
-    writes. Returns ``{"variants": [run, ...], "cross_verdict": {...}}``.
+    variant per requested method: ``detrend`` (#978b control), ``deseasonalize``
+    (#979 preprocessing), ``almon`` (#978 distributed-lag estimator). Returns
+    ``{"variants": [run, ...], "cross_verdict": {...}}``.
     """
     rate_by_month = load_rate_by_month(db, since=since)
     logger.info("loaded key_rate months: %d (since=%s)", len(rate_by_month), since)
 
-    variants = _plan_variants(sources, detrend)
+    variants = _plan_variants(sources, detrend, deseasonalize=deseasonalize, almon=almon)
     runs: list[dict[str, Any]] = []
-    for src, dt_flag in variants:
+    for src, dt_flag, deseason_flag, estimator in variants:
         runs.append(
             run_backtest(
                 db,
@@ -1147,6 +1492,8 @@ def run_all(
                 district=district,
                 source=src,
                 detrend=dt_flag,
+                deseasonalize=deseason_flag,
+                estimator=estimator,
                 rate_by_month=rate_by_month,
             )
         )
@@ -1177,13 +1524,22 @@ def render_table(results: dict[str, Any]) -> str:
     vd = results["verdict"]
     source = results["source"]
     detrended = results["detrended"]
+    deseasonalized = results.get("deseasonalized", False)
+    estimator = results.get("estimator", _ESTIMATOR_BEST_LAG)
+
+    # Short method tag for the title line (one active method per variant).
+    if estimator == _ESTIMATOR_ALMON:
+        method_tag = " · Almon-ADL"
+    elif deseasonalized:
+        method_tag = " · deseasonalized"
+    elif detrended:
+        method_tag = " · detrended"
+    else:
+        method_tag = ""
 
     lines: list[str] = []
     lines.append("=" * 78)
-    lines.append(
+    lines.append(f"BACKTEST [source {source}{method_tag}]: §9.6 rate-sensitivity OOS validation")
-        f"BACKTEST [source {source}{' · detrended' if detrended else ''}]: "
-        "§9.6 rate-sensitivity OOS validation"
-    )
     lines.append("=" * 78)
     lines.append(
         f"since={params['since']}  holdout_frac={params['holdout_frac']}  "
@@ -1196,6 +1552,19 @@ def render_table(results: dict[str, Any]) -> str:
             "removes a spurious monotone (survivorship) trend so it can't drive β. "
             "Trend fit on TRAIN months only, projected point-in-time onto test (no leakage)."
         )
+    if deseasonalized:
+        lines.append(
+            "DESEASONALIZED (#979): month-of-year seasonal factors divided out of units "
+            "BEFORE log_diff — strips the calendar pattern so a seasonal peak isn't read as "
+            "a rate effect. Factors fit on TRAIN months only, applied point-in-time (no leakage)."
+        )
+    if estimator == _ESTIMATOR_ALMON:
+        lines.append(
+            "ALMON-ADL (#978): Almon polynomial distributed-lag estimator (smooth β over lags "
+            "0..6) replaces single-lag best_lag. 'beta' column = LONG-RUN Σβ multiplier; 'lag' "
+            "= peak-|β_j| lag. Fit on TRAIN months only; test prediction sums all lags "
+            "point-in-time (no leakage)."
+        )
     if results.get("a_district_ignored"):
         lines.append(
             "NOTE: --district was IGNORED for Source A (corp_sum aggregates are not "
@@ -1362,6 +1731,22 @@ def _parse_args(argv: list[str] | None = None) -> argparse.Namespace:
         "ln(units) before differencing (removes a spurious monotone "
         "survivorship trend so it can't drive the regression).",
     )
+    p.add_argument(
+        "--deseasonalize",
+        action="store_true",
+        help="Also run a DESEASONALIZED variant of each source (#979): divide out "
+        "month-of-year seasonal factors before log_diff (strips the calendar "
+        "pattern so a seasonal peak isn't read as a rate effect). Factors fit "
+        "on TRAIN months only, applied point-in-time.",
+    )
+    p.add_argument(
+        "--almon",
+        action="store_true",
+        help="Also run an ALMON-ADL variant of each source (#978): replace the "
+        "single-lag best_lag estimator with an Almon polynomial distributed-lag "
+        "fit (smooth β over lags 0..6, long-run multiplier). Fit on TRAIN "
+        "months only; test prediction sums all lags point-in-time.",
+    )
     p.add_argument(
         "--holdout-frac",
         type=float,
@@ -1409,13 +1794,16 @@ def main(argv: list[str] | None = None) -> int:
     classes = _parse_classes(args.classes)
     sources = _parse_source(args.source)
     logger.info(
-        "backtest start: since=%s holdout_frac=%.2f classes=%s district=%s sources=%s detrend=%s",
+        "backtest start: since=%s holdout_frac=%.2f classes=%s district=%s sources=%s "
+        "detrend=%s deseasonalize=%s almon=%s",
         args.since,
         args.holdout_frac,
         "auto" if classes is None else classes,
         args.district,
         sources,
         args.detrend,
+        args.deseasonalize,
+        args.almon,
     )
 
     db = _session()
@@ -1428,6 +1816,8 @@ def main(argv: list[str] | None = None) -> int:
             district=args.district,
             sources=sources,
             detrend=args.detrend,
+            deseasonalize=args.deseasonalize,
+            almon=args.almon,
         )
     finally:
         db.close()

backend/tests/scripts/test_backtest_rate_sensitivity.py

@@ -9,11 +9,22 @@ Covers the PURE backtest logic on SYNTHETIC series (no live DB):
   - align_series         — inner-join by year-month
   - evaluate_oos         — inject sales=f(rate@lag) → high OOS hit-rate;
                            inject noise → hit-rate ≈ 0.5; point-in-time honesty
-  - backtest_tier        — thin-tier skip; happy path; (#978b) detrended variant
+  - evaluate_oos_almon   — (#978) Almon distributed-lag OOS evaluator: recovers a
-                           recovers an injected signal masked by a trend
+                           known peak lag + negative long-run on a clean signal;
+                           train fit IMMUNE to test-half corruption (no leakage);
+                           predictor never reads a future rate index; same return
+                           keys as evaluate_oos
+  - _deseasonalize_units — (#979) seasonal factors fit on TRAIN months only,
+                           applied point-in-time; recovers a known month pattern;
+                           a TEST-window spike does NOT move the fitted factors
+  - backtest_tier        — thin-tier skip; happy path; (#978b) detrended variant;
+                           (#978) almon estimator path; (#979) deseasonalize path;
+                           BACKWARD-COMPAT: default args == original raw best_lag
   - verdict / tier_lift  — promotion criterion, coin-flip baseline, lag stability
-  - _parse_source / _plan_variants — (#978b) B/A/both selection + variant plan
+  - _variant_label / _plan_variants — raw/detrended/deseasonalized/Almon-ADL
-  - cross_source_verdict — (#978b) B raw vs B detrended vs A conclusion
+                           labels + the per-flag variant plan (no all-combos)
+  - cross_source_verdict — controls (detrended/A) + candidate methods
+                           (deseasonalize #979, Almon-ADL #978) verdict + labels
 
 DB is MOCKED (a fake session) only to assert the Source A/B SQL SHAPE — that it
 uses CAST(:x AS type) and never the psycopg3-incompatible :x::type form, hits the
@@ -141,6 +152,98 @@ def _units_from_rate_with_trend(
     return units
 
 
+# --------------------------------------------------------------------------- #
+# Almon distributed-lag synthetic helpers (#978) — MIRROR the proven
+# construction in tests/services/forecasting/test_regression.py so the Almon
+# evaluator is exercised on a signal the estimator demonstrably recovers. The
+# regressor is a DIRECT LCG-jittered Δrate series (low cross-lag autocorrelation
+# → the per-lag reconstruction is faithful); the regressand is a quadratic-shaped
+# distributed lag the Almon deg-2 polynomial represents exactly.
+# --------------------------------------------------------------------------- #
+
+
+def _aperiodic_rate_deltas(n: int, *, seed: int = 13) -> list[float]:
+    """Δrate series with APERIODIC (LCG) jitter → low autocorrelation across lags.
+
+    Mirrors regression's ``_aperiodic_rate_deltas``: a periodic regressor would let
+    false lags compete with the injected one; LCG jitter keeps successive Δ weakly
+    correlated so the true lag shape wins. out[0] = 0.0 (finite from index 0); the
+    Almon lag-matrix builder drops incomplete leading rows itself.
+    """
+    lvl = 10.0
+    state = seed
+    levels: list[float] = []
+    for _ in range(n):
+        state = (state * 1103515245 + 12345) % 2147483648
+        lvl += 0.3 + (state / 2147483648.0 - 0.5) * 0.8
+        levels.append(lvl)
+    return [0.0] + [levels[i] - levels[i - 1] for i in range(1, n)]
+
+
+def _hump_beta(max_lag: int, *, peak: int, scale: float = 0.06) -> list[float]:
+    """A negative 'hump' lag shape peaking (in magnitude) at ``peak``. Mirror reg.
+
+    |β_j| = scale − 0.012·(j−peak)² (floored at 0.005), all signs negative — the
+    economically expected shape (rate ↑ → demand ↓, response builds then fades),
+    representable by an Almon deg-2 polynomial so the fit recovers the peak.
+    """
+    betas: list[float] = []
+    for j in range(max_lag + 1):
+        mag = scale - 0.012 * (j - peak) ** 2
+        betas.append(-max(0.005, mag))
+    return betas
+
+
+def _delta_sales_from_lag_shape(
+    rate_deltas: list[float], beta: list[float], *, max_lag: int
+) -> list[float | None]:
+    """delta_sales[t] = Σ_j β_j·rate_deltas[t−j]; leading (t<max_lag) → None.
+
+    The clean, noiseless distributed-lag regressand carrying the injected shape
+    exactly. ``evaluate_oos_almon`` fits β on TRAIN and predicts the same Σ form,
+    so on this construction the OOS directional hit-rate is ~1.0.
+    """
+    out: list[float | None] = [None] * max_lag
+    for t in range(max_lag, len(rate_deltas)):
+        out.append(sum(beta[j] * rate_deltas[t - j] for j in range(max_lag + 1)))
+    return out
+
+
+# --------------------------------------------------------------------------- #
+# Seasonal synthetic helpers (#979) — a units series carrying a KNOWN
+# month-of-year multiplicative pattern over ≥2 full years.
+# --------------------------------------------------------------------------- #
+
+# A known month-of-year seasonal pattern (multiplicative): summer peak, winter dip.
+_KNOWN_SEASONAL: dict[int, float] = {
+    1: 0.70,
+    2: 0.80,
+    3: 1.00,
+    4: 1.10,
+    5: 1.20,
+    6: 1.30,
+    7: 1.40,
+    8: 1.20,
+    9: 1.00,
+    10: 0.90,
+    11: 0.80,
+    12: 0.60,
+}
+
+
+def _seasonal_units(
+    months: list[dt.date], *, base: float = 1000.0, factor: dict[int, float] | None = None
+) -> list[float]:
+    """units[t] = base · factor[month_of(t)] — a clean known seasonal pattern.
+
+    Float values (the deseasonalize path is float-math throughout: divide by
+    factor then log_diff). With ≥2 full years the seasonal guard passes and
+    ``seasonal_factors`` recovers ``factor`` up to the overall-mean normalisation.
+    """
+    fac = factor or _KNOWN_SEASONAL
+    return [base * fac[m.month] for m in months]
+
+
 # --------------------------------------------------------------------------- #
 # _time_ordered_split
 # --------------------------------------------------------------------------- #
@@ -462,6 +565,247 @@ class TestEvaluateOos:
         assert shifted[res["n_train"]] is None or isinstance(shifted[res["n_train"]], float)
 
 
+# --------------------------------------------------------------------------- #
+# evaluate_oos_almon (#978) — the new Almon distributed-lag OOS evaluator
+# --------------------------------------------------------------------------- #
+
+
+class TestEvaluateOosAlmon:
+    def test_recovers_known_distributed_lag(self) -> None:
+        # Clean noiseless distributed lag with a quadratic hump peaking at lag 2.
+        # The Almon deg-2 fit on TRAIN must recover that peak and a negative
+        # long-run multiplier, and predict direction OOS ~perfectly (clean signal).
+        max_lag = 6
+        n = 72
+        rate_deltas = _aperiodic_rate_deltas(n)
+        beta = _hump_beta(max_lag, peak=2)
+        delta_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+
+        res = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        # train_lag = the fitted peak-|β_j| lag; matches the injected peak (±0).
+        assert res["train_lag"] == 2
+        # "train_beta" reports the long-run Σβ multiplier — negative here.
+        assert res["train_beta"] is not None and res["train_beta"] < 0
+        # Clean noiseless construction → directional hit-rate clearly beats coin.
+        assert res["oos_hit_rate"] is not None and res["oos_hit_rate"] > 0.5
+        assert res["oos_hit_rate"] >= 0.9
+        # In-sample R² is high by construction (reported, not trusted).
+        assert res["in_sample_r2"] is not None and res["in_sample_r2"] > 0.9
+        # Lag stable: the full-sample refit finds the same peak lag.
+        assert res["full_sample_lag"] == 2
+        assert res["lag_stable"] is True
+
+    def test_recovers_different_peak_lag(self) -> None:
+        # Shift the injected peak to lag 4 → the Almon fit must track it.
+        max_lag = 6
+        n = 80
+        rate_deltas = _aperiodic_rate_deltas(n, seed=29)
+        beta = _hump_beta(max_lag, peak=4)
+        delta_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+        res = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        assert res["train_lag"] == 4
+        assert res["oos_hit_rate"] is not None and res["oos_hit_rate"] >= 0.9
+
+    def test_no_look_ahead_leakage_train_fit_immune_to_test_corruption(self) -> None:
+        # Build a clean signal, then corrupt ONLY the test-half delta_sales (flip
+        # sign + scale + offset). The TRAIN fit cannot see the test window, so
+        # train_lag / train_beta / in_sample_r2 must be byte-identical to the
+        # uncorrupted run; only the OOS score may move.
+        max_lag = 6
+        n = 72
+        rate_deltas = _aperiodic_rate_deltas(n)
+        beta = _hump_beta(max_lag, peak=2)
+        clean_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+
+        n_train = bt._time_ordered_split(n, 0.7)
+        corrupt_sales: list[float | None] = list(clean_sales)
+        for t in range(n_train, n):
+            v = corrupt_sales[t]
+            if v is not None:
+                corrupt_sales[t] = -v * 5.0 + 0.123  # arbitrary test-only corruption
+
+        clean = bt.evaluate_oos_almon(clean_sales, rate_deltas, holdout_frac=0.7)
+        corrupt = bt.evaluate_oos_almon(corrupt_sales, rate_deltas, holdout_frac=0.7)
+
+        # TRAIN fit identical — the corruption is entirely in the held-out window.
+        assert clean["train_lag"] == corrupt["train_lag"]
+        assert clean["train_beta"] is not None and corrupt["train_beta"] is not None
+        assert math.isclose(clean["train_beta"], corrupt["train_beta"], rel_tol=0, abs_tol=1e-12)
+        assert clean["in_sample_r2"] is not None and corrupt["in_sample_r2"] is not None
+        assert math.isclose(
+            clean["in_sample_r2"], corrupt["in_sample_r2"], rel_tol=0, abs_tol=1e-12
+        )
+        # The OOS hit-rate DID respond to the corruption (flipped signs miss) —
+        # proving the test window is actually scored, not ignored.
+        assert clean["oos_hit_rate"] is not None and corrupt["oos_hit_rate"] is not None
+        assert corrupt["oos_hit_rate"] < clean["oos_hit_rate"]
+
+    def test_point_in_time_predictor_never_reads_future_rate(self) -> None:
+        # Structural no-future-leak assertion: the per-lag shifted views the
+        # evaluator reads at a test index t are _shift_for_lag(rate_deltas, j),
+        # whose element at t equals the ORIGINAL rate_deltas[t-j] (≤ t) — never an
+        # index > t. We assert this for every lag j across every test month.
+        max_lag = 6
+        n = 60
+        rate_deltas = _aperiodic_rate_deltas(n)
+        beta = _hump_beta(max_lag, peak=2)
+        delta_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+        res = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        n_train = res["n_train"]
+        for j in range(max_lag + 1):
+            shifted = bt._shift_for_lag(rate_deltas, j)
+            for t in range(n_train, n):
+                # The value the predictor uses at (t, lag j) is rate_deltas[t-j],
+                # which is at or before t (None when t-j < 0). It is NEVER t+k.
+                if shifted[t] is not None:
+                    assert t - j >= 0
+                    assert shifted[t] == rate_deltas[t - j]
+
+    def test_skips_test_month_with_incomplete_lag_profile(self) -> None:
+        # A None in the rate series punches a hole: a test month whose full lag
+        # profile can't be formed is skipped (not fabricated). With one rate hole
+        # near the test boundary, the evaluator still scores the remaining months
+        # and never crashes / never counts the holed month.
+        max_lag = 6
+        n = 72
+        rate_deltas: list[float | None] = list(_aperiodic_rate_deltas(n))
+        beta = _hump_beta(max_lag, peak=2)
+        delta_sales = _delta_sales_from_lag_shape(
+            [r if r is not None else 0.0 for r in rate_deltas], beta, max_lag=max_lag
+        )
+        n_train = bt._time_ordered_split(n, 0.7)
+        # Punch a hole in a test-window rate delta → the months that read it via
+        # any lag j become unscorable.
+        hole = n_train + 2
+        rate_deltas[hole] = None
+        res = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        # Still produced a result, fewer scored months than the raw test span.
+        assert res["oos_hit_rate"] is not None
+        assert res["n_test"] <= n - n_train
+        assert math.isfinite(res["oos_signed_mae"])
+
+    def test_too_few_months_returns_empty(self) -> None:
+        # 1 month → can't split → empty result (all metrics None, not a crash).
+        res = bt.evaluate_oos_almon([None], [None], holdout_frac=0.7)
+        assert res["train_lag"] is None
+        assert res["oos_hit_rate"] is None
+        assert res["n_train"] == 1 and res["n_test"] == 0
+
+    def test_infeasible_fit_returns_empty(self) -> None:
+        # Too few aligned points for the Almon fit (< _MIN_FIT_OBS usable rows) →
+        # fit_almon_dl returns None → empty result, no crash.
+        n = 20  # > min split but Almon needs more usable rows after max_lag drop
+        rate_deltas = _aperiodic_rate_deltas(n)
+        # Flat regressand → zero-variance / infeasible fit on the train slice.
+        delta_sales: list[float | None] = [None] + [0.0] * (n - 1)
+        res = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        assert res["train_lag"] is None
+        assert res["oos_hit_rate"] is None
+
+    def test_return_dict_has_same_keys_as_evaluate_oos(self) -> None:
+        # backtest_tier wraps both evaluators identically → identical key sets.
+        max_lag = 6
+        n = 60
+        rate_deltas = _aperiodic_rate_deltas(n)
+        beta = _hump_beta(max_lag, peak=2)
+        delta_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+        almon = bt.evaluate_oos_almon(delta_sales, rate_deltas, holdout_frac=0.7)
+        # evaluate_oos on the same arrays (best_lag) for a key-set comparison.
+        bl = bt.evaluate_oos(delta_sales, rate_deltas, holdout_frac=0.7)
+        assert set(almon.keys()) == set(bl.keys())
+
+
+# --------------------------------------------------------------------------- #
+# Deseasonalization (#979) — month-of-year factors recovered + TRAIN-only fit
+# --------------------------------------------------------------------------- #
+
+
+class TestDeseasonalize:
+    def test_recovers_known_seasonal_pattern(self) -> None:
+        # units = base · known_factor[month] over 3 full years → seasonal_factors
+        # must recover the known pattern (up to overall-mean normalisation) and
+        # deseasonalize_values must flatten the month-means to ~equal.
+        seasonal_factors, deseasonalize_values = bt._import_normalize()
+        n = 36  # 3 full years
+        ms = _months(n)
+        units = _seasonal_units(ms)
+        adj = seasonal_factors(ms, units)
+        assert adj.applied is True
+        assert adj.n_full_years == 3
+        # Expected normalised factor = known[m] / mean(known).
+        overall = sum(_KNOWN_SEASONAL.values()) / 12.0
+        for m in range(1, 13):
+            expected = _KNOWN_SEASONAL[m] / overall
+            assert math.isclose(adj.factors[m], expected, abs_tol=1e-9)
+        # Deseasonalized month-means collapse to a single value (pattern removed).
+        des = deseasonalize_values(ms, units, adj.factors)
+        by_month: dict[int, list[float]] = {}
+        for d, v in zip(ms, des, strict=False):
+            assert v is not None
+            by_month.setdefault(d.month, []).append(v)
+        means = [sum(vs) / len(vs) for vs in by_month.values()]
+        assert max(means) - min(means) < 1e-6
+
+    def test_factors_fit_on_train_only_immune_to_test_spike(self) -> None:
+        # Insert an EXTREME spike in a TEST-window month and assert the seasonal
+        # factors fit on the TRAIN slice are UNCHANGED vs the no-spike series. The
+        # train/test boundary is _time_ordered_split — exactly what
+        # _deseasonalize_units slices to.
+        seasonal_factors, _deseason = bt._import_normalize()
+        n = 48
+        ms = _months(n)
+        units = _seasonal_units(ms)
+        n_train = bt._time_ordered_split(n, 0.7)
+
+        clean = seasonal_factors(ms[:n_train], units[:n_train])
+        spiked_units = list(units)
+        spiked_units[n - 1] = spiked_units[n - 1] * 100.0  # extreme TEST-window spike
+        spiked = seasonal_factors(ms[:n_train], spiked_units[:n_train])
+
+        for m in range(1, 13):
+            assert math.isclose(clean.factors[m], spiked.factors[m], abs_tol=1e-12)
+        # Sanity: a LEAKY full-series fit WOULD have moved (the spike is real) —
+        # so the train-only slice is what protects us, not a no-op.
+        full_clean = seasonal_factors(ms, units)
+        full_spiked = seasonal_factors(ms, spiked_units)
+        assert any(abs(full_clean.factors[m] - full_spiked.factors[m]) > 1e-9 for m in range(1, 13))
+
+    def test_deseasonalize_units_helper_uses_time_ordered_boundary(self) -> None:
+        # The backtest helper _deseasonalize_units must fit factors on months[:fit_n]
+        # ONLY. We feed fit_n = _time_ordered_split and confirm the regressand it
+        # builds equals a manual TRAIN-fit-then-full-apply-then-log_diff, and is
+        # NOT equal to a leaky full-sample-fit version (when they differ).
+        seasonal_factors, deseasonalize_values = bt._import_normalize()
+        _bl, _ols, log_diff = bt._import_engine()
+        n = 48
+        ms = _months(n)
+        units_f = _seasonal_units(ms)
+        units = [max(1, round(v)) for v in units_f]
+        n_train = bt._time_ordered_split(n, 0.7)
+
+        # Spike a TEST-window month so train-fit and full-fit factors differ.
+        units[n - 1] = units[n - 1] * 50
+
+        got = bt._deseasonalize_units(ms, units, fit_n=n_train)
+
+        train_factors = seasonal_factors(ms[:n_train], units[:n_train]).factors
+        expected = log_diff(deseasonalize_values(ms, units, train_factors))
+        full_factors = seasonal_factors(ms, units).factors
+        leaky = log_diff(deseasonalize_values(ms, units, full_factors))
+
+        # The helper matches the TRAIN-only path exactly.
+        assert len(got) == len(expected)
+        for g, e in zip(got, expected, strict=False):
+            assert (g is None and e is None) or (
+                g is not None and e is not None and math.isclose(g, e, abs_tol=1e-12)
+            )
+        # And the train-only vs leaky paths genuinely differ (fix is observable).
+        assert any(
+            g is not None and lk is not None and abs(g - lk) > 1e-9
+            for g, lk in zip(got, leaky, strict=False)
+        )
+
+
 # --------------------------------------------------------------------------- #
 # backtest_tier — thin-tier skip + happy path
 # --------------------------------------------------------------------------- #
@@ -558,6 +902,124 @@ class TestBacktestTier:
         if raw.oos_hit_rate is not None and detr.oos_hit_rate is not None:
             assert detr.oos_hit_rate >= raw.oos_hit_rate
 
+    def test_records_deseasonalized_and_estimator_flags(self) -> None:
+        # The TierResult carries the new deseasonalize flag and estimator label.
+        n = 48
+        ms = _months(n)
+        rate = _aperiodic_rate_levels(n)
+        units = _units_from_rate(rate, lag=2, beta=-0.05)
+        sales = {ms[i]: units[i] for i in range(n)}
+        rate_by = {ms[i]: rate[i] for i in range(n)}
+        res = bt.backtest_tier(
+            sales, rate_by, tier=bt._EKB_WIDE, deseasonalize=True, estimator=bt._ESTIMATOR_ALMON
+        )
+        assert res.deseasonalized is True
+        assert res.estimator == bt._ESTIMATOR_ALMON
+        d = res.as_dict()
+        assert d["deseasonalized"] is True
+        assert d["estimator"] == bt._ESTIMATOR_ALMON
+
+    def test_almon_estimator_path_runs(self) -> None:
+        # estimator="almon" routes backtest_tier to evaluate_oos_almon. On a clean
+        # distributed-lag series it recovers the peak lag and scores OOS well.
+        max_lag = 6
+        n = 72
+        ms = _months(n)
+        rate_deltas = _aperiodic_rate_deltas(n)
+        # Reconstruct rate LEVELS from the deltas so align_series has a rate series;
+        # the tier re-differences them → the same rate_deltas reach the evaluator.
+        rate_levels = [10.0]
+        for j in range(1, n):
+            rate_levels.append(rate_levels[-1] + rate_deltas[j])
+        beta = _hump_beta(max_lag, peak=2)
+        delta_sales = _delta_sales_from_lag_shape(rate_deltas, beta, max_lag=max_lag)
+        # Turn the Δln signal into a units series (cumulative exp) so the tier's
+        # log_diff(units) reproduces delta_sales on the finite region.
+        ln_u = math.log(1000.0)
+        units: list[int] = [round(math.exp(ln_u))]
+        for t in range(1, n):
+            step = delta_sales[t] if delta_sales[t] is not None else 0.0
+            ln_u += step
+            units.append(max(1, round(math.exp(ln_u))))
+        sales = {ms[i]: units[i] for i in range(n)}
+        rate_by = {ms[i]: rate_levels[i] for i in range(n)}
+        res = bt.backtest_tier(
+            sales, rate_by, tier=bt._EKB_WIDE, estimator=bt._ESTIMATOR_ALMON, holdout_frac=0.7
+        )
+        assert res.skipped is None
+        assert res.estimator == bt._ESTIMATOR_ALMON
+        assert res.train_lag == 2
+        assert res.oos_hit_rate is not None and res.oos_hit_rate >= 0.8
+
+    def test_deseasonalize_path_runs_and_uses_train_only_fit(self) -> None:
+        # deseasonalize=True routes through _deseasonalize_units; the regressand it
+        # builds must equal a TRAIN-only-fit reconstruction (no leakage) and the
+        # tier still produces a scored result.
+        seasonal_factors, deseasonalize_values = bt._import_normalize()
+        _bl, _ols, log_diff = bt._import_engine()
+        n = 48
+        ms = _months(n)
+        # Seasonal units with a mild rate-driven drift so a lag can gate.
+        rate = _aperiodic_rate_levels(n)
+        rate_deltas = bt._rate_first_diff(rate)
+        ln_u = math.log(1000.0)
+        units: list[int] = [round(math.exp(ln_u) * _KNOWN_SEASONAL[ms[0].month])]
+        for t in range(1, n):
+            src = rate_deltas[t - 2] if t - 2 >= 1 and rate_deltas[t - 2] is not None else 0.0
+            ln_u += -0.04 * src
+            units.append(max(1, round(math.exp(ln_u) * _KNOWN_SEASONAL[ms[t].month])))
+        sales = {ms[i]: units[i] for i in range(n)}
+        rate_by = {ms[i]: rate[i] for i in range(n)}
+
+        res = bt.backtest_tier(
+            sales, rate_by, tier=bt._EKB_WIDE, deseasonalize=True, holdout_frac=0.7
+        )
+        assert res.skipped is None
+        assert res.deseasonalized is True
+
+        # The regressand the tier built equals a TRAIN-only-fit reconstruction.
+        months_aligned, units_aligned, _rates = bt.align_series(sales, rate_by)
+        n_train = bt._time_ordered_split(len(months_aligned), 0.7)
+        train_factors = seasonal_factors(months_aligned[:n_train], units_aligned[:n_train]).factors
+        expected = log_diff(deseasonalize_values(months_aligned, units_aligned, train_factors))
+        got = bt._deseasonalize_units(months_aligned, units_aligned, fit_n=n_train)
+        for g, e in zip(got, expected, strict=False):
+            assert (g is None and e is None) or (
+                g is not None and e is not None and math.isclose(g, e, abs_tol=1e-12)
+            )
+
+    def test_backward_compat_defaults_unchanged(self) -> None:
+        # The CRITICAL back-compat check: a default backtest_tier call (no
+        # deseasonalize, estimator=best_lag) must produce the SAME metric fields
+        # as the pre-change raw path. We pin every metric to an explicit raw
+        # best_lag run and confirm the new descriptor fields default correctly.
+        n = 48
+        ms = _months(n)
+        rate = _aperiodic_rate_levels(n)
+        units = _units_from_rate(rate, lag=2, beta=-0.05)
+        sales = {ms[i]: units[i] for i in range(n)}
+        rate_by = {ms[i]: rate[i] for i in range(n)}
+
+        res = bt.backtest_tier(sales, rate_by, tier=bt._EKB_WIDE, holdout_frac=0.7)
+        # New descriptor fields default to the production raw path.
+        assert res.deseasonalized is False
+        assert res.estimator == bt._ESTIMATOR_BEST_LAG
+        assert res.detrended is False
+
+        # Metric fields equal a direct evaluate_oos (best_lag) on the same arrays —
+        # i.e. the default path is byte-identical to the original implementation.
+        n_train = bt._time_ordered_split(n, 0.7)
+        delta_sales = bt._delta_sales_series(units, detrend=False, fit_n=n_train)
+        rate_deltas = bt._rate_first_diff([float(r) for r in rate])
+        direct = bt.evaluate_oos(delta_sales, rate_deltas, holdout_frac=0.7)
+        assert res.train_lag == direct["train_lag"]
+        assert res.train_beta == direct["train_beta"]
+        assert res.in_sample_r2 == direct["in_sample_r2"]
+        assert res.oos_hit_rate == direct["oos_hit_rate"]
+        assert res.oos_signed_mae == direct["oos_signed_mae"]
+        assert res.full_sample_lag == direct["full_sample_lag"]
+        assert res.lag_stable == direct["lag_stable"]
+
 
 # --------------------------------------------------------------------------- #
 # verdict / tier_lift
@@ -689,27 +1151,107 @@ class TestParseSource:
 
 
 class TestPlanVariants:
-    def test_raw_only_without_detrend(self) -> None:
+    # Each entry is (source, detrend, deseasonalize, estimator). The RAW
-        assert bt._plan_variants([bt._SOURCE_B], detrend=False) == [(bt._SOURCE_B, False)]
+    # reference (best_lag on raw units) is always first per source; method flags
+    # ADD one variant each (no all-combinations explosion).
+    _BL = bt._ESTIMATOR_BEST_LAG
+    _AL = bt._ESTIMATOR_ALMON
+
+    def test_raw_only_without_any_flag(self) -> None:
+        assert bt._plan_variants([bt._SOURCE_B], detrend=False) == [
+            (bt._SOURCE_B, False, False, self._BL)
+        ]
 
     def test_detrend_adds_detrended_variant_per_source(self) -> None:
         plan = bt._plan_variants([bt._SOURCE_B, bt._SOURCE_A], detrend=True)
         assert plan == [
-            (bt._SOURCE_B, False),
+            (bt._SOURCE_B, False, False, self._BL),
-            (bt._SOURCE_B, True),
+            (bt._SOURCE_B, True, False, self._BL),
-            (bt._SOURCE_A, False),
+            (bt._SOURCE_A, False, False, self._BL),
-            (bt._SOURCE_A, True),
+            (bt._SOURCE_A, True, False, self._BL),
         ]
 
+    def test_deseasonalize_adds_deseasonalized_variant(self) -> None:
+        plan = bt._plan_variants([bt._SOURCE_B], detrend=False, deseasonalize=True)
+        assert plan == [
+            (bt._SOURCE_B, False, False, self._BL),
+            (bt._SOURCE_B, False, True, self._BL),
+        ]
+
+    def test_almon_adds_almon_variant(self) -> None:
+        plan = bt._plan_variants([bt._SOURCE_B], detrend=False, almon=True)
+        assert plan == [
+            (bt._SOURCE_B, False, False, self._BL),
+            (bt._SOURCE_B, False, False, self._AL),
+        ]
+
+    def test_all_flags_add_one_variant_each_per_source(self) -> None:
+        # raw + detrended + deseasonalized + Almon-ADL, in that order, per source.
+        plan = bt._plan_variants(
+            [bt._SOURCE_B, bt._SOURCE_A], detrend=True, deseasonalize=True, almon=True
+        )
+        assert plan == [
+            (bt._SOURCE_B, False, False, self._BL),
+            (bt._SOURCE_B, True, False, self._BL),
+            (bt._SOURCE_B, False, True, self._BL),
+            (bt._SOURCE_B, False, False, self._AL),
+            (bt._SOURCE_A, False, False, self._BL),
+            (bt._SOURCE_A, True, False, self._BL),
+            (bt._SOURCE_A, False, True, self._BL),
+            (bt._SOURCE_A, False, False, self._AL),
+        ]
+
+    def test_no_all_combinations_explosion(self) -> None:
+        # Two method flags on one source → 1 raw + 2 method variants = 3, NOT the
+        # 2x2x... cross-product of preprocessing x estimator.
+        plan = bt._plan_variants([bt._SOURCE_B], detrend=True, almon=True)
+        assert len(plan) == 3
+        assert plan == [
+            (bt._SOURCE_B, False, False, self._BL),
+            (bt._SOURCE_B, True, False, self._BL),
+            (bt._SOURCE_B, False, False, self._AL),
+        ]
+
+
+class TestVariantLabel:
+    def test_raw_detrended_deseasonalized_almon_labels(self) -> None:
+        assert bt._variant_label(bt._SOURCE_B, False) == "B raw"
+        assert bt._variant_label(bt._SOURCE_B, True) == "B detrended"
+        assert bt._variant_label(bt._SOURCE_B, False, deseasonalize=True) == "B deseasonalized"
+        assert (
+            bt._variant_label(bt._SOURCE_A, False, estimator=bt._ESTIMATOR_ALMON) == "A Almon-ADL"
+        )
+
+    def test_estimator_takes_precedence_in_label(self) -> None:
+        # The planner never combines methods, but if both were set the estimator
+        # (the strongest method signal) names the variant.
+        assert (
+            bt._variant_label(bt._SOURCE_B, True, deseasonalize=True, estimator=bt._ESTIMATOR_ALMON)
+            == "B Almon-ADL"
+        )
+
 
 # --------------------------------------------------------------------------- #
 # cross_source_verdict (#978b) — B raw vs B detrended vs A
 # --------------------------------------------------------------------------- #
 
 
-def _run(source: str, detrended: bool, ekb: bt.TierResult) -> dict:
+def _run(
+    source: str,
+    detrended: bool,
+    ekb: bt.TierResult,
+    *,
+    deseasonalized: bool = False,
+    estimator: str = bt._ESTIMATOR_BEST_LAG,
+) -> dict:
     """Minimal run dict (only the fields cross_source_verdict reads)."""
-    return {"source": source, "detrended": detrended, "ekb_result": ekb}
+    return {
+        "source": source,
+        "detrended": detrended,
+        "deseasonalized": deseasonalized,
+        "estimator": estimator,
+        "ekb_result": ekb,
+    }
 
 
 class TestCrossSourceVerdict:
@@ -759,6 +1301,48 @@ class TestCrossSourceVerdict:
         assert cv["promote_any"] is False
         assert cv["signal_variants"] == []
 
+    def test_candidate_methods_labelled_and_no_signal(self) -> None:
+        # raw + deseasonalized + Almon-ADL all at/below coin-flip → REAL no signal,
+        # the conclusion mentions the candidate methods, and each variant is
+        # labelled by its method (not lumped under "raw"/"detrended").
+        runs = [
+            _run(bt._SOURCE_B, False, _tier(oos_hit_rate=0.48)),
+            _run(bt._SOURCE_B, False, _tier(oos_hit_rate=0.50), deseasonalized=True),
+            _run(
+                bt._SOURCE_B,
+                False,
+                _tier(oos_hit_rate=0.47),
+                estimator=bt._ESTIMATOR_ALMON,
+            ),
+        ]
+        cv = bt.cross_source_verdict(runs)
+        assert cv["promote_any"] is False
+        labels = [r["variant"] for r in cv["rows"]]
+        assert labels == ["B raw", "B deseasonalized", "B Almon-ADL"]
+        # The conclusion is generalised to the candidate methods.
+        assert "deseasonalize" in cv["conclusion"] and "Almon-ADL" in cv["conclusion"]
+        # Row descriptors carry the method so JSON consumers can filter.
+        assert cv["rows"][1]["deseasonalized"] is True
+        assert cv["rows"][2]["estimator"] == bt._ESTIMATOR_ALMON
+
+    def test_candidate_method_recovers_signal_is_flagged(self) -> None:
+        # raw best_lag no signal, but the Almon-ADL variant clears coin-flip+margin
+        # (lag stable) → flagged as a variant recovering signal worth inspecting.
+        runs = [
+            _run(bt._SOURCE_B, False, _tier(oos_hit_rate=0.49)),
+            _run(
+                bt._SOURCE_B,
+                False,
+                _tier(oos_hit_rate=0.82),
+                estimator=bt._ESTIMATOR_ALMON,
+            ),
+        ]
+        cv = bt.cross_source_verdict(runs)
+        assert cv["promote_any"] is True
+        assert "B Almon-ADL" in cv["signal_variants"]
+        # Conclusion offers the candidate-method reading.
+        assert "candidate method" in cv["conclusion"]
+
 
 # --------------------------------------------------------------------------- #
 # DB layer SQL SHAPE — mocked session, asserts CAST not :: and read-only