From c3ee3a05f708868d246f38ab59b496b53d6e4648 Mon Sep 17 00:00:00 2001 From: Nathan Walker Date: Thu, 24 Sep 2026 20:35:11 -0400 Subject: [PATCH 1/4] Refresh slice 1, pass 1: store training rows compactly (#131) TrainingRows (chimeraboost/training_rows.py) keeps the rows a fitted booster's leaves came from: plain numeric columns as bins of the fitted binner, cross-feature parents as raw floats, categoricals as codes plus their category values, y and raw weights. rebuild_X turns the store back into a raw-equivalent X that the existing replay refit consumes unchanged, so there is no second preprocessing path: every bin maps back to a value that bins identically, and only cross parents are read raw. GradientBoosting.replay_kwargs() returns the constructor kwargs of an equivalent booster, read off the base signature. Tests prove the round trip exact against the raw rows: the binned matrix, the fitted preprocessor state, and the replayed booster's leaf values and predictions, with and without appended rows. The public store_training_data / refresh() API follows in pass 2. Co-Authored-By: Claude Opus 5.5 --- benchmarks/CAMPAIGN_PLAN.md | 54 ++++++ benchmarks/REFRESH_PLAN.md | 116 +++++++++++++ chimeraboost/booster.py | 25 +++ chimeraboost/training_rows.py | 191 +++++++++++++++++++++ tests/test_refresh.py | 306 ++++++++++++++++++++++++++++++++++ 5 files changed, 692 insertions(+) create mode 100644 benchmarks/REFRESH_PLAN.md create mode 100644 chimeraboost/training_rows.py create mode 100644 tests/test_refresh.py diff --git a/benchmarks/CAMPAIGN_PLAN.md b/benchmarks/CAMPAIGN_PLAN.md index 17d6962..a76d53b 100644 --- a/benchmarks/CAMPAIGN_PLAN.md +++ b/benchmarks/CAMPAIGN_PLAN.md @@ -452,6 +452,60 @@ Recommended pick: **R1, R2, R3, R4 + H(1)(4)(5)**. R1 and R2 have free probes an ## Iteration log (append-only) +#### I066 2026-09-24 issue #131 slice 1 (`refresh(X, y)` on an opt-in stored training set; LIBRARY feature, opt-in, pre-registered) +why now: the focus rule's third issue. PR #169 merged (5b27b06), #81 +closed. Program file `benchmarks/REFRESH_PLAN.md` (design, decisions, +later slices). Branch `campaign/issue131-refresh-slice1` from main 5b27b06; +two muse passes: `20260924-issue131-slice1a-internals.md`, then +`...-slice1b-api.md` after review. +change: slice 1 as REFRESH_PLAN.md specifies (regressor + binary +classifier, single model; bagging, multiclass, `loss="Quantile"` and +random effects raise at fit when `store_training_data=True`). The default +path changes only by a bit-identical refactor (`_fit_gdiff` arithmetic +shared with the stored-rows path). +barriers: B13 (replay is a screening, not a selection instrument) and B2 +(the refit amplifies a bad audition) matched on "replay". Refresh selects +nothing and changes no default or audition; B13's bit-identical round trip +is the invariant this rung relies on. +forecast: (1) pass 1a: the stored-rows preprocessing path reproduces the +replay path's binned matrix bit for bit on the same rows, and with new rows +equals the replay path on the concatenated raw rows, over every column +block; (2) identity snapshot 186/186 and the goldens unchanged after 1a +and after 1b; (3) pass 1b: refresh with zero new rows is bit-identical to +the fitted model in every configuration of REFRESH_PLAN.md's test (a); +refresh(A) then refresh(B) equals refresh(A+B); (4) full suite green. Both +Pareto axes untouched (opt-in, default path bit-identical). +pass 1a (muse exit 0): the twin design as first written +(`fit_transform_stored`, `Binner.transform_block`, shared `_gdiff_means`, +`stored=` booster plumbing): forecast (1) HIT, every equality exact (binned +matrix, fitted preprocessor state, leaf values, predict_raw), 5 tests, full +suite green in the sandbox. DISCARDED at review: +551 library lines, ~400 of +them a second implementation of the preprocessing to be kept in step by +hand. Redesign (REFRESH_PLAN.md, "the ONE preprocessing path"): keep bins +for plain numerics and RAW floats only for cross parents, map each bin back +to a value that bins identically (bin = number of borders <= v, binning.py: +82-93), and feed the rebuilt X to the existing replay refit unchanged. +Linear leaves read the binned matrix (booster.py:888), so they are +unaffected. Same storage, one code path, no booster plumbing. +forecast for pass 1a' (`20260924-issue131-slice1a2-rows.md`): (1') every +stored bin round-trips exactly (every feature, every bin index incl. the +missing slot); `fit_transform` on the rebuilt X equals `fit_transform` on +the raw rows (matrix and fitted state), with and without new rows; the +booster replay refit on the rebuilt X equals it on the raw rows (leaf +values, predict_raw); library diff under ~200 lines. +pass 1a' (muse exit 0): `chimeraboost/training_rows.py` (173 lines) + +`GradientBoosting.replay_kwargs()` (25 lines, read off +`inspect.signature(_BaseBooster.__init__)`): 198 library lines against the +twin's 551. Forecast (1') HIT: tests (a)-(f) all bit for bit on a forced +block of all three cross kinds (12 diff/prod, 8 gdiff; plain [2, 3], +parents [0, 1, 4, 5]), a count column, one combo pair, linear leaves (99 of +100 trees with `lin_coef`), 5% NaN numerics, NaN categories, zero weights, +300 appended rows with unseen categories; the missing category maps back +through `"__nan__"` -> `np.nan`. Full suite 1230 passed, 1 skipped (conda +python); identity snapshot 186/186 bit-identical. Committed on the branch +as the pass-1 checkpoint. +verdict: PENDING(muse 1b, `20260924-issue131-slice1b-api.md`) + #### I065 2026-09-24 issue #81 (research cascade: dead self-test anchor, stale `ideas.py` flags; BENCH tooling + test, pre-registered) why now: the focus rule's second issue. PR #168 merged (3e02ab1), #84 closed. Branch `campaign/issue81-research-ideas` from main 3e02ab1; muse diff --git a/benchmarks/REFRESH_PLAN.md b/benchmarks/REFRESH_PLAN.md new file mode 100644 index 0000000..818c0ee --- /dev/null +++ b/benchmarks/REFRESH_PLAN.md @@ -0,0 +1,116 @@ +# REFRESH_PLAN — `refresh(X, y)` on an opt-in stored training set (issue #131) + +Started 2026-09-24. Program file for GitHub issue #131 (the maintainer's spec). +Design read from the code by a planning pass on 2026-09-24; file:line refs are +at main 3e02ab1. + +## What refresh is + +Tree structure, `lr_`, binner borders, count/cross selections and +`cat_combinations` stay pinned; leaf values, linear-leaf coefficients, ordered +target statistics and gdiff group means are recomputed by replaying every tree +on (stored rows + new rows). Mechanism = the existing structure-transfer refit +(`replay_donor`, booster.py:934-1011, tree.py:2126-2169), fed from stored rows +instead of raw X. + +## Findings that shape the design + +- Invariant holds: replaying a fitted model's own trees on its own training + rows reproduces it bit for bit, for `refit_full` "replay", True and False, on + the regressor and the binary classifier (diff/prod/gdiff crosses, count + column, linear leaves, weights, NaNs). Grown and replayed leaves share the + float-gradient kernels (quantization touches split search only, + tree.py:2433-2437, 2501-2507); replay consumes the random stream in the same + order (booster.py:901-903, 1002-1003). +- The rows to store are the rows the FINAL booster's leaves came from: all + rows after a refit (default), the 80% split in splitter order when the refit + is skipped (`refit_full=False`, quality 1-2, `loss="Quantile"`, + sklearn_api.py:890), X as given with an explicit `eval_set` or no early + stopping. +- Storage (REVISED 2026-09-24, see Log): numeric columns as BINS of the + pinned binner (uint16 today, uint8 when every column has <= 256 bins), except + the parents of any cross column (diff, prod or gdiff), which are kept as RAW + float64 because crosses are computed from raw values and gdiff group means + are refit on replay rows (preprocessing.py:420-465, 506-510, 685); + categoricals as int32 codes plus each column's categories in code order (new + rows append unseen categories in first-appearance order); y as the booster + saw it; raw weights. Row order is the booster's order (first-appearance + codes, positional TS permutations). Count, cross and TS columns are not + stored: they are rebuilt. +- The ONE preprocessing path: the store is turned back into a raw-equivalent X + and fed to the existing replay refit unchanged. A bin maps back to a value + that bins identically, exactly: the binner puts v in the bin equal to the + number of borders <= v, and non-finite v in the missing slot + (binning.py:82-93), so bin 0 -> just below the first border, bin k (1..m) -> + border k-1, the missing slot -> NaN. Linear leaves read the binned matrix + (booster.py:888), not raw values, so they are unaffected; only the cross + parents need their raw values. +- The spec's 10-20x storage saving is optimistic for numeric data: about 3x + when crosses engage (up to ~6 parent columns stay float64), about 6.5x when + they do not; far more for object-dtype categoricals. +- Exact reproduction needs an integer `random_state` (TS permutations are + redrawn from `random_state + t`). + +## Slice 1 (regressor + binary classifier, single model) + +- `store_training_data=False` (last ctor param, `_SKLEARN_ONLY`, validated in + `_check_flag_params`); `n_samples_trained_`; private `_training_data_` + (a `TrainingRows`), not a public `X_train_` (in sklearn that means raw X). +- `refresh(X, y, sample_weight=None) -> self`; zero rows allowed; replays the + FITTED configuration (ignores later `set_params`). +- Raises at fit: bagging (`n_ensembles > 1`, quality 4/5), multiclass, + `loss="Quantile"`, `random_effects=True`. Raises at refresh: no store, + unseen class label, weighted fit refreshed without weights. +- `temperature_` (binary) stays frozen. +- Files: a new `chimeraboost/training_rows.py` (`TrainingRows`: capture from + raw rows with a fitted preprocessor, append new rows, rebuild a + raw-equivalent X), booster.py (`replay_kwargs()` only), sklearn_api.py + (capture at the end of both `_fit_single`, `refresh` via a shared + `_refresh_single` that rebuilds X and runs the existing replay refit). +- Tests (`tests/test_refresh.py`): (a) zero new rows = identity, bit for bit, + over refit modes, eval_set, no ES, 256+ category column, all cross kinds, + linear leaves, weights, NaNs, ordered boosting, MAE/Poisson/custom + objective, subsample/colsample, binary; (b) refresh = a manual replay refit + on the concatenated raw rows; (c) nothing pinned moves, `n_samples_trained_` + grows; (c2) refresh(A) then refresh(B) = refresh(A+B); (d) pickle round + trip; (e) every error. +- Risk: a future preprocessing input that reads raw numeric values beyond + binning (as crosses do) would silently lose exactness for non-parent + columns. The refresh tests (a)/(b) over every column kind are the guard; + `training_rows.py` says so where it decides which columns stay raw. + +## Later slices (one GitHub issue each when slice 1 ships) + +2 bagging (per-member stores) · 3 multiclass (vector replay; keep multiclass +`refit_full` as is) · 4 quantiles (`loss="Quantile"`, the multi-quantile head) +· 5 random effects · 6 reservoir `store_training_data=int` · 7 `strict_ids`. + +## Decisions and open questions + +Slice 1 (stated to the maintainer 2026-09-24 as "unless you say otherwise", +no objection when he merged #169; revisit if he says so): +1. DECIDED: private `_training_data_` + public `n_samples_trained_`, not the + spec's public `X_train_` / `y_train_` / `sample_weight_` (in sklearn + `X_train_` means raw X; the store is bins, codes and a few raw columns). +2. DECIDED: a weighted fit refreshed without `sample_weight` raises. +3. DECIDED: `temperature_` stays frozen on refresh. + +Still open, for their slices: +4. (slice 2) bag members: every member gets all new rows, or a member-seeded + draw? +5. (slice 6) reservoir rows up-weighted by seen/capacity, or plain rows? + +## Log + +- 2026-09-24: design pass; slice 1 started as campaign rung I066 (branch + `campaign/issue131-refresh-slice1`), in two muse passes: internals, then + the public API and its tests. Docs by Claude. +- 2026-09-24: pass 1a built the design as first written (a stored-rows twin + of `fit_transform`: `fit_transform_stored`, `Binner.transform_block`, a + shared `_gdiff_means`, `stored=` plumbing through the booster). Every + equality was exact (binned matrix, fitted state, leaf values), but it added + ~550 library lines, ~400 of them a second implementation of the + preprocessing that would have to be kept in step by hand. DISCARDED at + review for the reconstruct design above: the same storage, one code path, + no booster plumbing. The twin's diff and tests are kept outside the repo + only as a reference. diff --git a/chimeraboost/booster.py b/chimeraboost/booster.py index f73f8a2..afa7ac0 100644 --- a/chimeraboost/booster.py +++ b/chimeraboost/booster.py @@ -824,6 +824,31 @@ def __init__(self, loss="RMSE", loss_kwargs=None, **kw): self.loss_name = loss self.loss_kwargs = loss_kwargs or {} + def replay_kwargs(self): + """Constructor kwargs that rebuild an equivalent booster. + + Every ``_BaseBooster.__init__`` parameter (each stored under its own + name, so read off ``inspect.signature`` rather than a hand-kept + list) plus ``loss`` (from ``loss_name``) and ``loss_kwargs``; lists + and dicts are copied. + """ + import inspect + + kw = {} + params = inspect.signature(_BaseBooster.__init__).parameters + for name in params: + if name == "self": + continue + v = getattr(self, name) + if isinstance(v, list): + v = list(v) + elif isinstance(v, dict): + v = dict(v) + kw[name] = v + kw["loss"] = self.loss_name + kw["loss_kwargs"] = dict(self.loss_kwargs) + return kw + def _fit_impl(self, X, y, cat_features=None, eval_set=None, sample_weight=None, callbacks=None, prep_cache=None): """Fit the additive model. diff --git a/chimeraboost/training_rows.py b/chimeraboost/training_rows.py new file mode 100644 index 0000000..f8822a9 --- /dev/null +++ b/chimeraboost/training_rows.py @@ -0,0 +1,191 @@ +"""Compact storage of the rows a fitted booster's leaves came from. + +Plain numeric columns are kept as bins of the fitted binner (feature-major), +cross parents as raw float64, categoricals as int32 codes plus per-column +categories in code order, y as float64 and the raw weights. ``rebuild_X`` +turns the store back into a raw-equivalent array that the existing replay +refit consumes unchanged -- there is deliberately no second preprocessing +path here. +""" + +import numpy as np + +from .binning import BIN_DTYPE +from .preprocessing import as_model_array +from .target_encoding import factorize + +# What factorize() turns a missing categorical (None, NaN, or refusing +# self-comparison) into; rebuild_X maps this category back to np.nan. +_MISSING = "__nan__" + + +def _split_plain_parent(prep): + """Split ``prep.num_features_`` into plain (binned) and parent (raw). + + Every numeric column appearing in any of ``prep.cross_pairs`` stays raw, + whatever the op: diff/prod crosses are computed from raw floats and gdiff + group means are refit on the replay rows. Any future preprocessing input + that reads raw numeric values beyond binning (as crosses do) must be + added to the raw set here, or refresh loses exactness. + """ + num = list(prep.num_features_) + inset = set(num) + parents = set() + for i, j, _op in prep.cross_pairs: + if i in inset: + parents.add(i) + if j in inset: + parents.add(j) + return [f for f in num if f not in parents], [f for f in num if f in parents] + + +def _values_for_bins(bins, borders): + """A value per bin that bins back identically (binning.py:82-93).""" + m = len(borders) + out = np.empty(len(bins), dtype=np.float64) + if m == 0: + out[:] = 0.0 + out[bins == 1] = np.nan + return out + miss = bins == m + 1 + zero = bins == 0 + out[zero] = np.nextafter(borders[0], -np.inf) + mid = ~(miss | zero) + out[mid] = borders[bins[mid].astype(np.int64) - 1] + out[miss] = np.nan + return out + + +def _extend_codes(stored_vals, new_col, stored_codes): + """Codes for [stored; new] with unseen categories appended in order.""" + ext = list(stored_vals) + idx = {v: i for i, v in enumerate(ext)} + ncodes, ncats = factorize(new_col) + loc = np.empty(len(ncats), dtype=np.int64) + for li, nc in enumerate(ncats): + if nc in idx: + loc[li] = idx[nc] + else: + idx[nc] = len(ext) + ext.append(nc) + loc[li] = len(ext) - 1 + new_codes = loc[ncodes].astype(np.int32) + full = np.concatenate([stored_codes, new_codes]) + return full, np.asarray(ext, dtype=object) + + +class TrainingRows: + """The rows a fitted booster's leaves came from, in its row order.""" + + def __init__(self, n_features, cat_features, num_features, plain_features, + parent_features, plain_bins, parent_raw, cat_codes, + cat_values, y, sample_weight): + self.n_features = int(n_features) + self.cat_features = list(cat_features) + self.num_features = list(num_features) + self.plain_features = list(plain_features) + self.parent_features = list(parent_features) + self.plain_bins = plain_bins + self.parent_raw = parent_raw + self.cat_codes = cat_codes + self.cat_values = cat_values + self.y = y + self.sample_weight = sample_weight + + @property + def n_rows(self): + return int(self.y.shape[0]) + + @classmethod + def capture(cls, prep, X, y, sample_weight=None): + """Capture raw rows with a FITTED ``FeaturePreprocessor``.""" + X = as_model_array(X, bool(prep.cat_features_)) + y = np.asarray(y, dtype=np.float64) + w = (None if sample_weight is None + else np.asarray(sample_weight, dtype=np.float64)) + plain, parents = _split_plain_parent(prep) + pos = {f: k for k, f in enumerate(prep.num_features_)} + cols = [pos[f] for f in plain] + nb = prep.binner_.n_bins_ + dt = (np.uint8 if all(nb[k] <= 256 for k in cols) else BIN_DTYPE) + if cols: + bins = np.ascontiguousarray(prep.transform(X)[:, cols].T, dtype=dt) + else: + bins = np.empty((0, X.shape[0]), dtype=dt) + pcols = [pos[f] for f in parents] + num = prep._numeric_block(X) + if pcols: + raw = np.ascontiguousarray(num[:, pcols], dtype=np.float64) + else: + raw = np.empty((X.shape[0], 0), dtype=np.float64) + codes = prep._codes_for_transform(X).astype(np.int32) + vals = [] + for m in prep.cat_maps_: + arr = np.empty(len(m), dtype=object) + for v, c in m.items(): + arr[c] = v + vals.append(arr) + return cls(X.shape[1], prep.cat_features_, prep.num_features_, + plain, parents, bins, raw, codes, vals, y, w) + + def append(self, prep, X_new, y_new, sample_weight_new=None): + """A NEW ``TrainingRows`` for [stored; new]; neither mutates.""" + X_new = as_model_array(X_new, bool(self.cat_features)) + y_new = np.asarray(y_new, dtype=np.float64) + w_new = (None if sample_weight_new is None + else np.asarray(sample_weight_new, dtype=np.float64)) + if X_new.shape[1] != self.n_features: + raise ValueError("X_new has %d columns, stored %d" + % (X_new.shape[1], self.n_features)) + if (self.sample_weight is None) != (w_new is None): + raise ValueError("weighted/unweighted append mismatch") + pos = {f: k for k, f in enumerate(prep.num_features_)} + cols = [pos[f] for f in self.plain_features] + if cols: + nb = np.ascontiguousarray( + prep.transform(X_new)[:, cols].T, + dtype=self.plain_bins.dtype) + bins = np.concatenate([self.plain_bins, nb], axis=1) + else: + bins = np.empty((0, self.n_rows + len(y_new)), + dtype=self.plain_bins.dtype) + pcols = [pos[f] for f in self.parent_features] + num = prep._numeric_block(X_new) + praw = np.concatenate( + [self.parent_raw, + np.ascontiguousarray(num[:, pcols], dtype=np.float64) + if pcols else np.empty((len(y_new), 0))]) + codes, vals = [], [] + for j, f in enumerate(self.cat_features): + c, v = _extend_codes(self.cat_values[j], X_new[:, f], + self.cat_codes[:, j]) + codes.append(c) + vals.append(v) + stacked = (np.column_stack(codes).astype(np.int32) if codes + else np.empty((self.n_rows + len(y_new), 0), + dtype=np.int32)) + y = np.concatenate([self.y, y_new]) + w = (None if w_new is None + else np.concatenate([self.sample_weight, w_new])) + return TrainingRows(self.n_features, self.cat_features, + self.num_features, self.plain_features, + self.parent_features, bins, praw, stacked, + vals, y, w) + + def rebuild_X(self, prep): + """The raw-equivalent array in the original column layout.""" + n = self.n_rows + X = (np.empty((n, self.n_features), dtype=object) + if self.cat_features + else np.empty((n, self.n_features), dtype=np.float64)) + pos = {f: k for k, f in enumerate(prep.num_features_)} + for p, f in enumerate(self.plain_features): + X[:, f] = _values_for_bins( + self.plain_bins[p], prep.binner_.borders_[pos[f]]) + for j, f in enumerate(self.parent_features): + X[:, f] = self.parent_raw[:, j] + for j, f in enumerate(self.cat_features): + col = self.cat_values[j][self.cat_codes[:, j].astype(np.int64)] + col[col == _MISSING] = np.nan + X[:, f] = col + return X diff --git a/tests/test_refresh.py b/tests/test_refresh.py new file mode 100644 index 0000000..c9d442b --- /dev/null +++ b/tests/test_refresh.py @@ -0,0 +1,306 @@ +"""Stored training rows round-trip (issue #131 slice 1a2). + +A fitted booster's rows are kept as binner bins (plain numerics), raw +floats (cross parents) and int codes (categoricals); ``rebuild_X`` turns +them back into a raw-equivalent array the existing replay refit consumes +unchanged. Every comparison here is bit for bit. +""" + +import inspect +import pickle + +import numpy as np +import pytest + +from chimeraboost import ChimeraBoostRegressor +from chimeraboost.booster import GradientBoosting, _BaseBooster +from chimeraboost.preprocessing import FeaturePreprocessor +from chimeraboost.target_encoding import factorize +from chimeraboost.training_rows import TrainingRows + +CAT = [6, 7] +N_NUM = 6 + + +def _cats(rng, n, levels, p_nan): + c = np.array(rng.choice(levels, n), dtype=object) + c[rng.random(n) < p_nan] = np.nan + return c + + +def _target(rng, num, hi, lo): + x = np.where(np.isfinite(num), num, 0.0) + base = (3.0 * (x[:, 0] > x[:, 1]) + x[:, 2] * x[:, 3] + + 0.5 * x[:, 4] - 0.3 * x[:, 5]) + lo_eff = np.array([{"a": 0.0, "b": 0.5, "c": -0.5, "d": 1.0, "e": -1.0}.get( + v, 0.0) for v in lo]) + hi_eff = np.zeros(len(lo)) + for k, v in enumerate(hi): + if isinstance(v, str) and v.startswith("L"): + try: + hi_eff[k] = (int(v[1:]) % 10) * 0.2 + except ValueError: + pass + return base + lo_eff + hi_eff + 0.1 * rng.standard_normal(len(lo)) + + +def _weights(rng, n): + w = rng.uniform(0.5, 1.5, n) + w[rng.random(n) < 0.01] = 0.0 + return w + + +def _frame(num, hi, lo): + X = np.empty((len(hi), N_NUM + 2), dtype=object) + X[:, :N_NUM] = num + X[:, 6] = hi + X[:, 7] = lo + return X + + +def _make_data(seed=0): + rng = np.random.default_rng(seed) + n_full, n_held, n_new = 2600, 400, 300 + num = rng.standard_normal((n_full + n_held + n_new, N_NUM)) + num[rng.random(num.shape) < 0.05] = np.nan + hi_levels = [f"L{i}" for i in range(350)] + lo_levels = ["a", "b", "c", "d", "e"] + hi_full = _cats(rng, n_full, hi_levels, 0.02) + lo_full = _cats(rng, n_full, lo_levels, 0.02) + hi_held = _cats(rng, n_held, hi_levels, 0.02) + lo_held = _cats(rng, n_held, lo_levels, 0.02) + hi_new = _cats(rng, n_new, hi_levels, 0.05) + lo_new = _cats(rng, n_new, lo_levels, 0.05) + unseen = rng.random(n_new) < 0.25 + hi_new[unseen] = rng.choice([f"NEW{i}" for i in range(10)], unseen.sum()) + lo_new[rng.random(n_new) < 0.15] = "z" + n0, n1 = n_full, n_full + n_held + d = {} + d["X_full"] = _frame(num[:n_full], hi_full, lo_full) + d["X_held"] = _frame(num[n0:n1], hi_held, lo_held) + d["X_new"] = _frame(num[n1:], hi_new, lo_new) + d["y_full"] = _target(rng, num[:n_full], hi_full, lo_full) + d["y_held"] = _target(rng, num[n0:n1], hi_held, lo_held) + d["y_new"] = _target(rng, num[n1:], hi_new, lo_new) + d["w_full"] = _weights(rng, n_full) + d["w_new"] = _weights(rng, n_new) + return d + + +@pytest.fixture(scope="module") +def donor(): + """A realistic fitted donor: forced cross block, count column, combos.""" + d = _make_data() + est = ChimeraBoostRegressor(random_state=0, n_estimators=100, + cross_features="always", cat_combinations=True, + linear_leaves=True) + est.fit(d["X_full"], d["y_full"], cat_features=CAT, + sample_weight=d["w_full"]) + b, prep = est.model_, est.model_.prep_ + assert len(b.trees_) > 0 + assert any(t.lin_coef is not None for t in b.trees_) + assert prep.count_features_ == [6] + assert len(prep.cat_maps_[0]) >= 256 + ops = {op for _, _, op in b.cross_pairs} + assert ops == {"diff", "prod", "gdiff"}, b.cross_pairs + # Default refit_full="replay" retrains on every row: the RMSE init is the + # FULL-data weighted mean, not the 80% split's. + assert b.init_ == pytest.approx( + float(np.average(d["y_full"], weights=d["w_full"]))) + d["est"], d["booster"], d["prep"] = est, b, prep + d["rows"] = TrainingRows.capture(prep, d["X_full"], d["y_full"], + d["w_full"]) + assert d["rows"].plain_features and d["rows"].parent_features + return d + + +def _pinned_prep(donor_booster, donor_prep): + """A fresh prep with the replay-refit pinning (_prep_or_replay_matrices).""" + fresh = FeaturePreprocessor( + donor_booster.max_bins, donor_booster.cat_smoothing, + donor_booster.random_state, donor_booster.cat_n_permutations, + bool(donor_prep.combo_pairs_), donor_booster.cross_pairs, + donor_booster.cat_count_features) + fresh._pinned_count_features = list(donor_prep.count_features_) + fresh._pinned_cat_counts = (donor_prep.cat_maps_, donor_prep.cat_counts_) + return fresh + + +def _fit_pinned(donor, X, y, w): + prep = _pinned_prep(donor["booster"], donor["prep"]) + mat = prep.fit_transform(X, [y], CAT, w, + binner=donor["prep"].binner_) + return prep, mat + + +def _assert_prep_equal(p1, p2): + assert p1.cat_maps_ == p2.cat_maps_ + assert p1.combo_pairs_ == p2.combo_pairs_ + assert p1.combo_maps_ == p2.combo_maps_ + assert len(p1.gdiff_maps_) == len(p2.gdiff_maps_) + for (d1, g1), (d2, g2) in zip(p1.gdiff_maps_, p2.gdiff_maps_): + assert d1 == d2 + assert g1 == g2 + assert p1.count_features_ == p2.count_features_ + assert len(p1.cat_counts_) == len(p2.cat_counts_) + for c1, c2 in zip(p1.cat_counts_, p2.cat_counts_): + np.testing.assert_array_equal(c1, c2) + assert len(p1.encoders_) == len(p2.encoders_) + for e1, e2 in zip(p1.encoders_, p2.encoders_): + assert e1.prior_ == e2.prior_ + assert e1.n_cat_ == e2.n_cat_ + for a1, a2 in zip(e1.sums_, e2.sums_): + np.testing.assert_array_equal(a1, a2) + for a1, a2 in zip(e1.counts_, e2.counts_): + np.testing.assert_array_equal(a1, a2) + np.testing.assert_array_equal(p1.feature_map_, p2.feature_map_) + np.testing.assert_array_equal(p1.is_numeric_binned_, p2.is_numeric_binned_) + np.testing.assert_array_equal(p1.n_bins_, p2.n_bins_) + + +def test_round_trip_same_rows(donor): + """(a) Rebuilt X refits exactly like the raw rows, same rows.""" + rows = donor["rows"] + assert rows.n_rows == len(donor["y_full"]) + X_rb = rows.rebuild_X(donor["prep"]) + assert X_rb.shape == donor["X_full"].shape + assert X_rb.dtype == object + w = GradientBoosting._normalize_weights(donor["w_full"], + len(donor["y_full"])) + p1, m1 = _fit_pinned(donor, X_rb, donor["y_full"], w) + p2, m2 = _fit_pinned(donor, donor["X_full"], donor["y_full"], w) + assert m1.dtype == m2.dtype + np.testing.assert_array_equal(m1, m2) + _assert_prep_equal(p1, p2) + np.testing.assert_array_equal(p1.transform(donor["X_held"]), + p2.transform(donor["X_held"])) + + +def test_round_trip_appended_rows(donor): + """(b) The same with new rows appended through ``append``.""" + prep = donor["prep"] + before = ([dict(m) for m in prep.cat_maps_], + [(dict(d), g) for d, g in prep.gdiff_maps_], + [b.copy() for b in prep.binner_.borders_]) + rows = donor["rows"] + n_old = rows.n_rows + rows2 = rows.append(prep, donor["X_new"], donor["y_new"], donor["w_new"]) + assert rows.n_rows == n_old + assert rows2.n_rows == n_old + len(donor["y_new"]) + assert [dict(m) for m in prep.cat_maps_] == before[0] + assert [(dict(d), g) for d, g in prep.gdiff_maps_] == before[1] + for b, b0 in zip(prep.binner_.borders_, before[2]): + np.testing.assert_array_equal(b, b0) + X_cat = np.concatenate([donor["X_full"], donor["X_new"]]) + y_cat = np.concatenate([donor["y_full"], donor["y_new"]]) + w_cat = np.concatenate([donor["w_full"], donor["w_new"]]) + for j, f in enumerate(CAT): + codes, cats = factorize(X_cat[:, f]) + np.testing.assert_array_equal(rows2.cat_codes[:, j], codes) + assert list(rows2.cat_values[j]) == list(cats) + X_rb = rows2.rebuild_X(prep) + assert X_rb.shape == X_cat.shape + w = GradientBoosting._normalize_weights(w_cat, len(y_cat)) + p1, m1 = _fit_pinned(donor, X_rb, y_cat, w) + p2, m2 = _fit_pinned(donor, X_cat, y_cat, w) + np.testing.assert_array_equal(m1, m2) + _assert_prep_equal(p1, p2) + np.testing.assert_array_equal(p1.transform(donor["X_held"]), + p2.transform(donor["X_held"])) + + +def _fit_replay(donor, X, linear): + kw = donor["booster"].replay_kwargs() + kw["n_estimators"] = len(donor["booster"].trees_) + kw["learning_rate"] = float(donor["booster"].lr_) + kw["early_stopping_rounds"] = None + kw["replay_donor"] = (donor["booster"].trees_, donor["prep"]) + if linear: + kw["linear_leaves"] = True + b = GradientBoosting(**kw) + b.fit(X, donor["y_full"], cat_features=CAT, + sample_weight=donor["w_full"]) + return b + + +def _assert_booster_equal(b1, b2): + assert b1.init_ == b2.init_ + assert len(b1.trees_) == len(b2.trees_) + for t1, t2 in zip(b1.trees_, b2.trees_): + np.testing.assert_array_equal(t1.values, t2.values) + if t1.lin_coef is None or t2.lin_coef is None: + assert t1.lin_coef is None and t2.lin_coef is None + else: + np.testing.assert_array_equal(t1.lin_coef, t2.lin_coef) + + +@pytest.mark.parametrize("linear", [False, True]) +def test_booster_replay_on_rebuilt(donor, linear): + """(c) Replay refit on rebuilt X equals the one on raw rows. + + False is the donor as fitted (linear leaves), True forces + ``linear_leaves=True``; replay refits ``lin_coef`` off a linear donor. + """ + X_rb = donor["rows"].rebuild_X(donor["prep"]) + b1 = _fit_replay(donor, X_rb, linear) + b2 = _fit_replay(donor, donor["X_full"], linear) + assert any(t.lin_coef is not None for t in b1.trees_) + _assert_booster_equal(b1, b2) + np.testing.assert_array_equal(b1.predict_raw(donor["X_held"]), + b2.predict_raw(donor["X_held"])) + + +def test_exhaustive_bin_round_trip(donor): + """(d) Every plain-numeric bin maps to a value that bins back to it.""" + binner = donor["prep"].binner_ + pos = {f: k for k, f in enumerate(donor["prep"].num_features_)} + for f in donor["rows"].plain_features: + borders = binner.borders_[pos[f]] + m = len(borders) + vals = np.empty(m + 2, dtype=np.float64) + if m == 0: + vals[0], vals[1] = 0.0, np.nan + else: + vals[0] = np.nextafter(borders[0], -np.inf) + vals[1:m + 1] = borders + vals[m + 1] = np.nan + feat = np.zeros((m + 2, len(binner.borders_))) + feat[:, pos[f]] = vals + got = binner.transform(feat)[:, pos[f]] + np.testing.assert_array_equal(got, np.arange(m + 2)) + + +def test_replay_kwargs_covers_init(donor): + """(e) ``replay_kwargs`` covers every booster init parameter.""" + kw = donor["booster"].replay_kwargs() + base = set(inspect.signature(_BaseBooster.__init__).parameters) - {"self"} + assert set(kw) == base | {"loss", "loss_kwargs"} + assert kw["loss"] == donor["booster"].loss_name + assert kw["loss_kwargs"] == donor["booster"].loss_kwargs + assert kw["cross_pairs"] == donor["booster"].cross_pairs + assert kw["cross_pairs"] is not donor["booster"].cross_pairs + assert kw["loss_kwargs"] is not donor["booster"].loss_kwargs + b2 = GradientBoosting(**kw) + assert b2.loss_name == donor["booster"].loss_name + assert b2.cross_pairs == donor["booster"].cross_pairs + + +def test_training_rows_pickle_and_dtype(donor): + """(f) ``TrainingRows`` pickles exactly; plain bins are uint8.""" + rows = donor["rows"] + assert rows.plain_bins.dtype == np.uint8 + rt = pickle.loads(pickle.dumps(rows)) + assert rt.n_features == rows.n_features + assert rt.cat_features == rows.cat_features + assert rt.num_features == rows.num_features + assert rt.plain_features == rows.plain_features + assert rt.parent_features == rows.parent_features + assert rt.plain_bins.dtype == rows.plain_bins.dtype + np.testing.assert_array_equal(rt.plain_bins, rows.plain_bins) + np.testing.assert_array_equal(rt.parent_raw, rows.parent_raw) + np.testing.assert_array_equal(rt.cat_codes, rows.cat_codes) + for v1, v2 in zip(rt.cat_values, rows.cat_values): + np.testing.assert_array_equal(v1, v2) + np.testing.assert_array_equal(rt.y, rows.y) + np.testing.assert_array_equal(rt.sample_weight, rows.sample_weight) + assert rt.n_rows == rows.n_rows From 8b4cbcc13f2cb1d5cac09ce99c8a612eb3ef2759 Mon Sep 17 00:00:00 2001 From: Nathan Walker Date: Thu, 24 Sep 2026 21:00:27 -0400 Subject: [PATCH 2/4] Refresh slice 1, pass 2: store_training_data and refresh(X, y) (#131) ChimeraBoostRegressor and the binary ChimeraBoostClassifier gain store_training_data=False. When True, the fit keeps the rows the final booster's leaves came from (TrainingRows, from pass 1), and refresh(X, y, sample_weight=None) appends new rows, rebuilds a raw-equivalent X from the store, and replays the fitted booster's configuration on it: the same structure-transfer replay the default full-data refit uses, with every tree's splits, the tree count and the learning rate pinned, and no size-adaptive setting re-resolved. n_samples_trained_ reports the rows the leaves come from. Bagging, multiclass, loss="Quantile" and random_effects=True raise a clear NotImplementedError at fit when the flag is on. refresh refuses a model without a store, unseen class labels, and a weights mismatch. Tests: a refresh with no new rows is bit-identical to the fitted model in 15 configurations; a refresh equals a manual replay on the stacked raw rows; two refreshes equal one combined refresh; pinned state never moves; pickling and every error are covered. The default fit is unchanged (identity snapshot 186/186). Docs: a parameters row, a "Refreshing with new rows" recipe, and a CHANGELOG entry. Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 13 + benchmarks/CAMPAIGN_PLAN.md | 27 +- benchmarks/REFRESH_PLAN.md | 5 + chimeraboost/sklearn_api.py | 299 ++++++++++++++++++++++- docs/parameters.md | 1 + docs/recipes.md | 39 +++ tests/test_refresh.py | 474 +++++++++++++++++++++++++++++++++++- 7 files changed, 852 insertions(+), 6 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 2b4897a..c94bd0d 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,6 +4,19 @@ All notable changes to ChimeraBoost are documented here. The format follows [Keep a Changelog](https://keepachangelog.com/). ## [Unreleased] +### Added +- **`store_training_data` and `refresh(X, y)`** on the regressor and the + binary classifier (#131, first slice). With `store_training_data=True` a + model keeps its training rows in compact form (bins for plain numeric + columns, raw values only for the columns that feed cross features, + category codes), and `refresh` folds in new rows by replaying every + tree's splits on the old and new rows together and refitting only the + leaf values. It is the same replay the default full-data refit uses, + with the tree count, learning rate and splits pinned; refreshing with no + new rows reproduces the model bit for bit. Bagged models, multiclass, + `loss="Quantile"` and `random_effects=True` are not supported yet. The + default fit is unchanged. + ### Changed - **The quantile benchmark's NGBoost opponent now uses the RoNGBa settings** (Ren, Sun and Wu 2019; #163): trees of up to 31 leaves, a diff --git a/benchmarks/CAMPAIGN_PLAN.md b/benchmarks/CAMPAIGN_PLAN.md index a76d53b..54a16f7 100644 --- a/benchmarks/CAMPAIGN_PLAN.md +++ b/benchmarks/CAMPAIGN_PLAN.md @@ -504,7 +504,32 @@ parents [0, 1, 4, 5]), a count column, one combo pair, linear leaves (99 of through `"__nan__"` -> `np.nan`. Full suite 1230 passed, 1 skipped (conda python); identity snapshot 186/186 bit-identical. Committed on the branch as the pass-1 checkpoint. -verdict: PENDING(muse 1b, `20260924-issue131-slice1b-api.md`) +pass 1b (muse exit 0): `sklearn_api.py` +291 (the parameter, the slice-1 +gate inside `_quality_applied` plus the multiclass check once classes are +known, capture at the end of both `_fit_single`, `refresh` through +`_refresh_single` and helpers, each C901 <= 6); `training_rows.py` +unchanged; 29 new tests. Slice 1 library total 489 lines (training_rows +173, booster 25, sklearn_api 291). +result: (2) HIT: identity snapshot 186/186 after both passes. (3) HIT: +zero-row refresh bit-identical in all 15 configurations (refit_full +"replay"/True/False, eval_set, no early stopping, a 350-level categorical, +all three cross kinds, linear leaves, weights, NaNs, ordered boosting, MAE, +Poisson, a custom objective, subsample + colsample, a weighted binary +classifier with string labels); refresh equals a manual replay on the +stacked raw rows; refresh(A) then refresh(B) equals refresh(A + B); pinned +state unchanged; pickle; every error. (4) HIT: full suite 1259 passed, 1 +skipped (conda python); ruff clean on `chimeraboost/`. Review: a DataFrame +with reordered columns is refused by `_check_feature_names_match` before +any append; SHAP after a refresh adds up to `predict_raw` within 4e-13. +Usefulness smoke (not a gate): test RMSE 0.8543 on the 60% model, 0.7975 +after refreshing with the next 30%, 0.7698 for a full refit on 90%: +refresh recovers ~67% of the full refit's gain without growing a tree. +docs (Claude): `docs/parameters.md` (the row), `docs/recipes.md` +("Refreshing with new rows"), CHANGELOG (Unreleased, Added); the API pages +render the `refresh` docstring. +verdict: **PASS → PR for the maintainer** (library, tests, docs). Issue +#131 stays open for slices 2-7 (REFRESH_PLAN.md). +next: issue #113 (random effects, slice 2) per the focus rule. #### I065 2026-09-24 issue #81 (research cascade: dead self-test anchor, stale `ideas.py` flags; BENCH tooling + test, pre-registered) why now: the focus rule's second issue. PR #168 merged (3e02ab1), #84 diff --git a/benchmarks/REFRESH_PLAN.md b/benchmarks/REFRESH_PLAN.md index 818c0ee..421cd71 100644 --- a/benchmarks/REFRESH_PLAN.md +++ b/benchmarks/REFRESH_PLAN.md @@ -114,3 +114,8 @@ Still open, for their slices: review for the reconstruct design above: the same storage, one code path, no booster plumbing. The twin's diff and tests are kept outside the repo only as a reference. +- 2026-09-24: slice 1 done (I066): pass 1a' `training_rows.py` + `replay_kwargs` + (198 lines), pass 1b the public API (291 lines); zero-row refresh exact in + 15 configurations, chained refreshes compose, identity snapshot unchanged. + Smoke: a 60% model refreshed with 30% more rows recovered ~67% of a full + 90% refit's RMSE gain. Awaiting the maintainer's merge. diff --git a/chimeraboost/sklearn_api.py b/chimeraboost/sklearn_api.py index c348789..02e74f8 100644 --- a/chimeraboost/sklearn_api.py +++ b/chimeraboost/sklearn_api.py @@ -12,6 +12,7 @@ from .random_effects import (codes_for_labels, estimate_ratio_reml, solve_intercepts) from .target_encoding import factorize +from .training_rows import TrainingRows from sklearn.base import BaseEstimator, RegressorMixin, ClassifierMixin @@ -60,7 +61,7 @@ def loss(T): "cat_features", "cross_features", "cross_top_columns", "selection_rounds", "refit_full", "refit_members", - "quality", "random_effects"}) + "quality", "random_effects", "store_training_data"}) # --- quality: named operating points on the strength/slowdown Pareto -------- # Evidence: benchmarks/SELECT_PLAN.md. Every recipe only pins parameters that @@ -235,7 +236,7 @@ def _check_depth(p, name): def _check_flag_params(estimator, p): - """The tri-state string/bool flags: refit_full, cross_features, quality.""" + """The string/bool flags: refit_full, cross_features, quality, store.""" v = p.get("refit_full") if v is not None and v != "replay" and not isinstance(v, (bool, np.bool_)): raise ValueError( @@ -257,6 +258,12 @@ def _check_flag_params(estimator, p): + ", ".join(f"{k} ({n})" for k, n in QUALITY_NAMES.items()) + f"; got {v!r}.") + if "store_training_data" in p: + v = p["store_training_data"] + if not isinstance(v, (bool, np.bool_)): + raise ValueError( + f"store_training_data must be True or False; got {v!r}.") + def _check_loss_family(p): # Regressor-only loss / alpha (the classifier picks its loss automatically). @@ -1870,6 +1877,174 @@ def _add_callback(callbacks, extra): return base + [extra] +def _raise_store_unsupported(what): + """The slice-1 coverage error for ``store_training_data=True`` fits.""" + raise NotImplementedError( + f"store_training_data=True does not support {what} yet; refresh() " + "currently covers single-model regression and binary classification.") + + +def _check_store_fit(est): + """Raise unless this ``store_training_data=True`` fit is slice-1 shaped. + + Called inside ``_quality_applied`` so quality 4/5 already resolved to + ``n_ensembles``; the multiclass arm runs in the classifier's + ``_fit_single`` once the classes are known. + """ + if not est.store_training_data: + return + if est.n_ensembles and est.n_ensembles > 1: + _raise_store_unsupported("n_ensembles > 1 (including quality=4 and 5)") + if getattr(est, "loss", None) == "Quantile": + _raise_store_unsupported("loss='Quantile'") + if getattr(est, "random_effects", False): + _raise_store_unsupported("random_effects=True") + + +def _check_store_multiclass(est): + """The multiclass arm of the slice-1 gate (classes are known by now).""" + if est.store_training_data and est.n_classes_ > 2: + _raise_store_unsupported("multiclass classification") + + +def _capture_training_store(est, cap_full, cap_train, refit, *, + classification): + """Set ``n_samples_trained_`` and ``_training_data_`` at the end of fit. + + The captured rows are ``cap_full`` (every row) when the full-data refit + replaced the model, else ``cap_train`` (the post-split training rows -- + X as given with an explicit eval_set or ``early_stopping=False``), each + an ``(X, y, sample_weight)`` triple in the booster's row order. The + captured ``y`` is exactly as the booster saw it (0/1 floats for binary), + the weights raw and un-normalized. + """ + cap_X, cap_y, cap_w = cap_full if refit else cap_train + est.n_samples_trained_ = int(len(cap_y)) + if not est.store_training_data: + est._training_data_ = None + return + if classification: + # Binary only (multiclass raised before any fitting work). + cap_y = (cap_y == est.classes_[1]).astype(np.float64) + est._training_data_ = TrainingRows.capture( + est.model_.prep_, cap_X, cap_y, cap_w) + + +def _refresh_check_y(y, n_new, classification): + """Validate refresh targets; the empty case skips the target-type probe.""" + if n_new: + return _check_y_target(y, n_new, classification) + y_arr = np.asarray(y) + if y_arr.shape[0] != 0: + raise ValueError( + "X and y have inconsistent lengths: X has 0 samples, " + f"y has {y_arr.shape[0]}.") + if y_arr.ndim == 2: + if y_arr.shape[1] == 1: + return y_arr.ravel() + raise ValueError( + "Multi-output y is not supported; pass a 1D y of shape " + "(n_samples,).") + return y_arr + + +def _refresh_check_w(sample_weight, n_new): + """Validate refresh weights; the empty case only checks the shape.""" + if sample_weight is None: + return None + if n_new: + _check_sample_weight_arr(sample_weight, n_new) + return np.asarray(sample_weight, dtype=np.float64) + sw = np.asarray(sample_weight, dtype=np.float64) + if sw.ndim != 1 or sw.shape[0] != 0: + raise ValueError( + f"sample_weight must be 1D of length 0; got shape {sw.shape}.") + return sw + + +def _refresh_check_store_weights(rows0, sample_weight): + """The weighted/unweighted refresh presence checks.""" + had_weights = rows0.sample_weight is not None + if had_weights and sample_weight is None: + raise ValueError( + "This model was fit with sample_weight, so refresh() needs " + "sample_weight for the new rows too, in the same units.") + if not had_weights and sample_weight is not None: + raise ValueError( + "This model was fit without sample_weight, so refresh() cannot " + "take sample_weight for the new rows; refit with sample_weight " + "to use weights.") + + +def _refresh_map_labels(est, y_arr): + """Map refresh labels to the fitted 0/1 encoding, rejecting unseen ones.""" + yv = np.asarray(y_arr) + if yv.size == 0: + return np.empty(0, dtype=np.float64) + try: + unseen = np.setdiff1d(np.unique(yv), + np.unique(np.asarray(est.classes_))) + except TypeError: # mixed un-orderable label types + seen = set(np.asarray(est.classes_).tolist()) + unseen = np.array([v for v in dict.fromkeys(yv.tolist()) + if v not in seen], dtype=object) + if unseen.size: + raise ValueError( + f"y contains class label(s) {unseen.tolist()} not seen at fit; " + "refresh() cannot add classes -- refit instead.") + return (yv == est.classes_[1]).astype(np.float64) + + +def _refresh_replay(est, rows, X_all): + """Replay the fitted booster's configuration on the combined rows. + + Everything is pinned at its fitted value -- the size-adaptive autos are + NOT re-resolved the way ``_refit_on_full`` does -- and the donor's + histories are carried over so ``validation_history_`` still reports the + curve that chose the budget. + """ + donor = est.model_ + kw = donor.replay_kwargs() + kw["n_estimators"] = len(donor.trees_) + kw["learning_rate"] = float(donor.lr_) + kw["early_stopping_rounds"] = None + kw["replay_donor"] = (donor.trees_, donor.prep_) + b = GradientBoosting(**kw) + b.fit(X_all, rows.y, cat_features=list(donor.prep_.cat_features_), + sample_weight=rows.sample_weight) + b.train_history_ = donor.train_history_ + b.valid_history_ = donor.valid_history_ + est.model_ = b + est._training_data_ = rows + est.n_samples_trained_ = rows.n_rows + est._shap_importances_cache_ = None + if hasattr(est, "expected_value_"): + del est.expected_value_ + return est + + +def _refresh_single(est, X, y, sample_weight, *, classification): + """Shared body of the estimators' ``refresh`` (issue #131 slice 1).""" + Xv = _check_predict_input(est, X) + rows0 = getattr(est, "_training_data_", None) + if rows0 is None: + raise ValueError( + "refresh() needs the rows this model was trained on, but it was " + "fit with store_training_data=False. Refit with " + "store_training_data=True to enable refresh().") + n_new = len(Xv) if Xv is not None else len(X) + y_arr = _refresh_check_y(y, n_new, classification) + w_new = _refresh_check_w(sample_weight, n_new) + _refresh_check_store_weights(rows0, sample_weight) + if classification: + y_arr = _refresh_map_labels(est, y_arr) + else: + y_arr = np.asarray(y_arr, dtype=np.float64) + prep = est.model_.prep_ + rows = rows0.append(prep, X, y_arr, w_new) + return _refresh_replay(est, rows, rows.rebuild_X(prep)) + + class _RegBoosterFactory: """Value-capturing replacement for the regressor `_fit_single`'s old `_fit_booster` / `_screen` closures. @@ -2174,6 +2349,17 @@ class ChimeraBoostRegressor(RegressorMixin, BaseEstimator): ``predict`` takes the row groups and adds the fitted intercepts; unseen groups get exactly 0. Slice 1: ``loss="RMSE"`` single models only (not ``n_ensembles > 1``). + store_training_data : bool, default False + Keep the rows the final booster's leaves came from, so ``refresh`` + can fold in new rows later without a from-scratch refit. What is + stored is compact, not raw X: plain numeric columns as binner bins, + the raw values of only the cross-feature parent columns, + categoricals as integer codes plus their per-column categories, the + target as the booster saw it, and the raw sample weights. The store + rides along in pickles and grows with every ``refresh`` call. + Roughly 3x smaller than float64 X for numeric data when cross + features engage, about 6.5x when they do not; far smaller for + object-dtype categoricals. Attributes ---------- @@ -2212,6 +2398,10 @@ class ChimeraBoostRegressor(RegressorMixin, BaseEstimator): group_ratio_ : float or None Fitted noise-to-group variance ratio (``inf`` means no group signal was found). ``None`` unless fit with ``random_effects=True``. + n_samples_trained_ : int + Number of rows the final booster's leaves came from: all rows when + the full-data refit ran, else the post-split training rows. Grown + by ``refresh``. """ # Both estimators accept cross_features="always" (the unrefereed forced @@ -2239,7 +2429,8 @@ def __init__(self, n_estimators=2000, learning_rate=None, depth=None, cat_features=None, quantize_gradients=True, eval_metric=None, delta=1.0, tweedie_variance_power=1.5, refit_full="replay", refit_members=False, quality=None, - adaptive_learning_rate=True, random_effects=False): + adaptive_learning_rate=True, random_effects=False, + store_training_data=False): self.n_estimators = n_estimators self.learning_rate = learning_rate self.depth = depth @@ -2282,6 +2473,7 @@ def __init__(self, n_estimators=2000, learning_rate=None, depth=None, # consulted when learning_rate is None. False == the historical flat 0.1. self.adaptive_learning_rate = adaptive_learning_rate self.random_effects = random_effects + self.store_training_data = store_training_data def fit(self, X, y, cat_features=None, eval_set=None, groups=None, sample_weight=None, callbacks=None): @@ -2344,6 +2536,7 @@ def fit(self, X, y, cat_features=None, eval_set=None, groups=None, _check_feature_names_match(self, eval_set[0]) with _quality_applied(self): + _check_store_fit(self) if self.random_effects: if self.loss != "RMSE": raise ValueError( @@ -2587,6 +2780,7 @@ def _fit_single(self, X, y, cat_features, eval_set, groups, sample_weight, n_re_groups) y_full_for_refit = (np.asarray(y_full, dtype=np.float64) - b_pre[group_codes_full]) + pre_refit = self.model_ self._dispatch_reg_refit(kw, loss_kwargs, X_full, y_full_for_refit, sw_full, cat_features, auto_split, cat_ctx=full_ctx) @@ -2600,6 +2794,9 @@ def _fit_single(self, X, y, cat_features, eval_set, groups, sample_weight, self.model_, X_full, y_full, sw_full, group_codes_full, n_re_groups) + _capture_training_store( + self, (X_full, y_full, sw_full), (X, y, sample_weight), + self.model_ is not pre_refit, classification=False) return self def _arm_forced_cross(self, fb, select_ll, ll): @@ -2810,6 +3007,40 @@ def _transform_raw(self, raw): tf = getattr(self.model_.loss_, "transform", None) return raw if tf is None else tf(raw) + def refresh(self, X, y, sample_weight=None): + """Fold new rows into this fitted model without a from-scratch refit. + + Appends ``(X, y[, sample_weight])`` to the rows stored at fit time + (``store_training_data=True``), rebuilds a raw-equivalent training + matrix from the store, and replays the FITTED booster's configuration + on it -- the same structure-transfer replay the default full-data + refit uses. Every tree's structure, the round count and the learning + rate stay pinned; leaf values and linear-leaf coefficients are refit + against the combined gradients. + + ``refresh`` replays the fitted configuration and ignores any + ``set_params`` made after fit. The binner borders, the count and + cross selections, the validation history and the fitted selections + are unchanged; ``n_samples_trained_`` grows by the number of new + rows. + + Parameters + ---------- + X, y : array-like + New rows with the same features the model was fit on. Zero rows + are allowed, in which case predictions are unchanged. + sample_weight : array-like of shape (n_samples,) or None + Weights for the new rows, in the same units as the fit weights. + Required when the model was fit with weights, rejected when it + was not. + + Returns + ------- + self + """ + return _refresh_single(self, X, y, sample_weight, + classification=False) + def predict_raw(self, X): """Raw additive score before the loss link and conformal quantile offset. @@ -3235,6 +3466,17 @@ class ChimeraBoostClassifier(ClassifierMixin, BaseEstimator): called without its own ``cat_features`` (the fit argument overrides). Provided as a constructor argument so ``GridSearchCV``/``Pipeline`` can carry it. + store_training_data : bool, default False + Keep the rows the final booster's leaves came from, so ``refresh`` + can fold in new rows later without a from-scratch refit. What is + stored is compact, not raw X: plain numeric columns as binner bins, + the raw values of only the cross-feature parent columns, + categoricals as integer codes plus their per-column categories, the + target as the booster saw it, and the raw sample weights. The store + rides along in pickles and grows with every ``refresh`` call. + Roughly 3x smaller than float64 X for numeric data when cross + features engage, about 6.5x when they do not; far smaller for + object-dtype categoricals. Attributes ---------- @@ -3254,6 +3496,10 @@ class ChimeraBoostClassifier(ClassifierMixin, BaseEstimator): Bagged-mode member defaults that were auto-applied (params the user left on auto resolve to tuned member values inside a bag; explicit values always win). Set only when ``n_ensembles > 1``. + n_samples_trained_ : int + Number of rows the final booster's leaves came from: all rows when + the full-data refit ran, else the post-split training rows. Grown + by ``refresh``. """ # Not pinned on the classifier: linear_leaves=None is already an auto rule @@ -3281,7 +3527,8 @@ def __init__(self, n_estimators=2000, learning_rate=None, depth=6, n_ensembles=None, ensemble_n_jobs=-1, max_samples=0.8, cat_features=None, quantize_gradients=True, eval_metric=None, refit_full="replay", refit_members=False, - quality=None, adaptive_learning_rate=True): + quality=None, adaptive_learning_rate=True, + store_training_data=False): self.n_estimators = n_estimators self.learning_rate = learning_rate self.depth = depth @@ -3319,6 +3566,7 @@ def __init__(self, n_estimators=2000, learning_rate=None, depth=6, # Size fade for the auto learning rate, default-on since 0.30.0; only # consulted when learning_rate is None. False == the historical flat 0.1. self.adaptive_learning_rate = adaptive_learning_rate + self.store_training_data = store_training_data def fit(self, X, y, cat_features=None, eval_set=None, groups=None, sample_weight=None, callbacks=None): @@ -3375,6 +3623,7 @@ def fit(self, X, y, cat_features=None, eval_set=None, groups=None, _check_eval_labels(eval_set, y) with _quality_applied(self): + _check_store_fit(self) if self.n_ensembles and self.n_ensembles > 1: if callbacks is not None: raise ValueError( @@ -3668,6 +3917,7 @@ def _fit_single(self, X, y, cat_features, eval_set, groups, sample_weight, (X, y, sample_weight, eval_set, es_active, auto_split, X_full, y_full, sw_full, split_idx) = self._resolve_classes_and_split( X, y, cat_features, eval_set, groups, sample_weight) + _check_store_multiclass(self) full_ctx, train_ctx, val_ctx = _shared_cat_ctxs( X_full, split_idx, cat_features) @@ -3754,11 +4004,52 @@ def _fit_single(self, X, y, cat_features, eval_set, groups, sample_weight, self.temperature_ = _fit_temperature(raw, cal_y, self._multiclass, sample_weight=cal_w) + pre_refit = self.model_ self._dispatch_cls_refit(kw, X_full, y_full, sw_full, cat_features, auto_split, cat_ctx=full_ctx) + _capture_training_store( + self, (X_full, y_full, sw_full), (X, y, sample_weight), + self.model_ is not pre_refit, classification=True) return self + def refresh(self, X, y, sample_weight=None): + """Fold new rows into this fitted binary model without a from-scratch + refit. + + Appends ``(X, y[, sample_weight])`` to the rows stored at fit time + (``store_training_data=True``), rebuilds a raw-equivalent training + matrix from the store, and replays the FITTED booster's configuration + on it -- the same structure-transfer replay the default full-data + refit uses. Every tree's structure, the round count and the learning + rate stay pinned; leaf values and linear-leaf coefficients are refit + against the combined gradients. + + ``refresh`` replays the fitted configuration and ignores any + ``set_params`` made after fit. The binner borders, the count and + cross selections, the validation history, the fitted selections and + the frozen ``temperature_`` are unchanged; ``n_samples_trained_`` + grows by the number of new rows. + + Parameters + ---------- + X, y : array-like + New rows with the same features the model was fit on. ``y`` + carries original class labels, mapped through the fitted + ``classes_`` -- labels never seen at fit raise. Zero rows are + allowed, in which case predictions are unchanged. + sample_weight : array-like of shape (n_samples,) or None + Weights for the new rows, in the same units as the fit weights. + Required when the model was fit with weights, rejected when it + was not. + + Returns + ------- + self + """ + return _refresh_single(self, X, y, sample_weight, + classification=True) + def predict_proba(self, X): Xv = _check_predict_input(self, X) X = X if Xv is None else Xv diff --git a/docs/parameters.md b/docs/parameters.md index fcc89e8..bfc75b4 100644 --- a/docs/parameters.md +++ b/docs/parameters.md @@ -137,6 +137,7 @@ See the User Guide: [early stopping](recipes.md#early-stopping) for `eval_set` a | Parameter | Default | Effect | |---|---|---| | `random_effects` | `False` | Fit one shrunk intercept per `groups=` label on top of the trees (regressor, RMSE single-model only). Small groups pool toward zero; unseen groups predict trees-only. See the User Guide: [grouped data](recipes.md#grouped-data-random-intercepts). | +| `store_training_data` | `False` | Keep the rows the final model was fit on, in compact form, so `refresh(X, y)` can fold in new rows later by refitting only the leaf values. Regressor and binary classifier, single models only. The stored rows travel with the pickled model and grow with every refresh. See the User Guide: [refreshing with new rows](recipes.md#refreshing-with-new-rows). | ## System diff --git a/docs/recipes.md b/docs/recipes.md index 523e875..1157800 100644 --- a/docs/recipes.md +++ b/docs/recipes.md @@ -361,6 +361,45 @@ reg = ChimeraBoostRegressor(refit_full=False, random_state=0).fit(X_train, y_tra reg = ChimeraBoostRegressor(quality=2, random_state=0).fit(X_train, y_train) ``` +## Refreshing with new rows + +When new labelled rows arrive after a fit, `refresh` folds them in without starting +over. It keeps every tree's splits and refits only the leaf values on the old rows plus +the new ones, which is the same replay the default full-data refit uses, so it costs +about what that refit costs: roughly a third of growing the model again. The model has +to keep its training rows for this, so ask for it at fit time: + +```python +reg = ChimeraBoostRegressor(store_training_data=True, random_state=0) +reg.fit(X_old, y_old) +reg.refresh(X_new, y_new) # leaves refit on the old and new rows together +reg.n_samples_trained_ # how many rows the leaves now come from +``` + +What a refresh leaves alone: + +- The trees' splits, the number of trees and the learning rate. A refresh cannot add a + split the original data never called for, so a shift that needs new splits still + needs a refit. +- The binary classifier's probability calibration (`temperature_`). +- Every setting chosen at fit time. `set_params` after the fit does not change what a + refresh does. + +Worth knowing: + +- The stored rows take much less room than X itself. Plain numeric columns are kept as + bin indices, only the few columns that feed cross features keep their raw values, and + categories are kept as integer codes. For numeric data expect about a third of the size + of X as 64-bit floats. The store travels with the pickled model and grows with every + refresh. +- If the model was fit with `sample_weight`, pass weights for the new rows too, in the + same units. +- The classifier cannot learn new classes this way. Refit instead. +- Fit with an integer `random_state` if you need refreshes to be exactly reproducible. +- Bagged models (`n_ensembles > 1`, including `quality=4` and `5`), multiclass + classification, `loss="Quantile"` and `random_effects=True` do not support + `store_training_data` yet. + ## Early stopping Early stopping is on by default. With no `eval_set`, the estimator holds out a diff --git a/tests/test_refresh.py b/tests/test_refresh.py index c9d442b..674b0dc 100644 --- a/tests/test_refresh.py +++ b/tests/test_refresh.py @@ -4,6 +4,9 @@ floats (cross parents) and int codes (categoricals); ``rebuild_X`` turns them back into a raw-equivalent array the existing replay refit consumes unchanged. Every comparison here is bit for bit. + +Slice 1b (below the line) covers the public ``store_training_data`` / +``refresh`` API; every comparison there is bit for bit too, never allclose. """ import inspect @@ -11,8 +14,11 @@ import numpy as np import pytest +from sklearn.base import clone +from sklearn.exceptions import NotFittedError -from chimeraboost import ChimeraBoostRegressor +from chimeraboost import (ChimeraBoostClassifier, ChimeraBoostRegressor, + CustomObjective) from chimeraboost.booster import GradientBoosting, _BaseBooster from chimeraboost.preprocessing import FeaturePreprocessor from chimeraboost.target_encoding import factorize @@ -304,3 +310,469 @@ def test_training_rows_pickle_and_dtype(donor): np.testing.assert_array_equal(rt.y, rows.y) np.testing.assert_array_equal(rt.sample_weight, rows.sample_weight) assert rt.n_rows == rows.n_rows + + +# --------------------------------------------------------------------------- +# Slice 1b: the public store_training_data / refresh API (issue #131). +# --------------------------------------------------------------------------- + + +class _SqErr(CustomObjective): + """Squared error through the public custom-objective hook. + + Module level, so instances are picklable. + """ + + def init(self, y, sample_weight=None): + return float(np.average(y, weights=sample_weight)) + + def grad_hess(self, y, raw): + return raw - y, np.ones_like(raw) + + def eval(self, y, raw, sample_weight=None): + return float(np.sqrt(np.average((raw - y) ** 2, + weights=sample_weight))) + + +def _binary_labels(y): + return np.where(np.asarray(y) > np.median(y), "pos", "neg") + + +def _snapshot(est, Xh): + """Predictions, init and per-tree leaves to compare a refresh against.""" + snap = {"pred": est.predict(Xh), "init": est.model_.init_, + "n": est.n_samples_trained_, + "trees": [(t.values.copy(), + None if t.lin_coef is None else t.lin_coef.copy()) + for t in est.model_.trees_]} + if hasattr(est, "predict_proba"): + snap["proba"] = est.predict_proba(Xh) + return snap + + +def _assert_snapshot_equal(est, snap, Xh): + np.testing.assert_array_equal(est.predict(Xh), snap["pred"]) + if "proba" in snap: + np.testing.assert_array_equal(est.predict_proba(Xh), snap["proba"]) + assert est.model_.init_ == snap["init"] + assert len(est.model_.trees_) == len(snap["trees"]) + for t, (v, lc) in zip(est.model_.trees_, snap["trees"]): + np.testing.assert_array_equal(t.values, v) + if lc is None or t.lin_coef is None: + assert lc is None and t.lin_coef is None + else: + np.testing.assert_array_equal(t.lin_coef, lc) + assert est.n_samples_trained_ == snap["n"] + + +# (a) estimator-kwarg overrides per configuration. Every fit also pins an +# integer random_state, store_training_data=True, and (regressor) the +# single-fit linear_leaves=False / cross_features=False unless the +# configuration says otherwise. NaNs ride in every configuration: _make_data +# salts 5% of the numerics and 2-5% of the categoricals with NaN. +_A_EST = { + "refit-replay": dict(refit_full="replay"), + "refit-true": dict(refit_full=True), + "refit-false": dict(refit_full=False), + "eval-set": {}, + "no-es": dict(early_stopping=False, n_estimators=20), + "count": {}, + "cross-always": dict(cross_features="always"), + "linear": dict(linear_leaves=True), + "weighted": {}, + "ordered": dict(ordered_boosting=True), + "mae": dict(loss="MAE"), + "poisson": dict(loss="Poisson"), + "custom": dict(loss=_SqErr()), + "subsample": dict(subsample=0.8, colsample=0.7), + "binary": {}, +} + + +def _cfg_parts(cfg, d): + """The (a) estimator, fit kwargs and fit target/weights for one config.""" + base = dict(random_state=0, n_estimators=50, store_training_data=True, + linear_leaves=False, cross_features=False) + base.update(_A_EST[cfg]) + fit_kw = dict(cat_features=CAT) + y_fit, w_fit = d["y_full"], None + cls = ChimeraBoostRegressor + if cfg == "eval-set": + fit_kw["eval_set"] = (d["X_held"], d["y_held"]) + if cfg in ("weighted", "binary"): + fit_kw["sample_weight"] = d["w_full"] + w_fit = d["w_full"] + if cfg == "poisson": + y_fit = np.exp((y_fit - y_fit.mean()) / y_fit.std()) + if cfg == "binary": + cls = ChimeraBoostClassifier + del base["linear_leaves"] + y_fit = _binary_labels(d["y_full"]) + return cls(**base), fit_kw, y_fit, w_fit + + +@pytest.mark.parametrize("cfg", ["refit-replay", "refit-true", "refit-false", + "eval-set", "no-es", "count", "cross-always", + "linear", "weighted", "ordered", "mae", + "poisson", "custom", "subsample", "binary"]) +def test_zero_rows_is_identity(cfg): + """(a) refresh() with zero new rows changes nothing, bit for bit.""" + d = _make_data() + est, fit_kw, y_fit, w_fit = _cfg_parts(cfg, d) + est.fit(d["X_full"], y_fit, **fit_kw) + if cfg == "count": + assert est.model_.prep_.count_features_ == [6] + if cfg == "cross-always": + ops = {op for _, _, op in est.model_.cross_pairs} + assert ops == {"diff", "prod", "gdiff"}, est.model_.cross_pairs + if cfg == "linear": + assert any(t.lin_coef is not None for t in est.model_.trees_) + snap = _snapshot(est, d["X_held"]) + w0 = None if w_fit is None else w_fit[:0] + out = est.refresh(d["X_full"][:0], y_fit[:0], sample_weight=w0) + assert out is est + _assert_snapshot_equal(est, snap, d["X_held"]) + + +def test_refresh_equals_manual_replay(): + """(b) refresh() equals a manual replay refit on the concatenated rows. + + The new rows carry unseen categories, NaN categories and NaN numerics + (see _make_data). + """ + d = _make_data() + est = ChimeraBoostRegressor(random_state=0, n_estimators=60, + cross_features="always", linear_leaves=True, + store_training_data=True) + est.fit(d["X_full"], d["y_full"], cat_features=CAT, + sample_weight=d["w_full"]) + donor = est.model_ + assert donor.prep_.count_features_ == [6] + est.refresh(d["X_new"], d["y_new"], sample_weight=d["w_new"]) + assert est.n_samples_trained_ == len(d["y_full"]) + len(d["y_new"]) + + kw = donor.replay_kwargs() + kw["n_estimators"] = len(donor.trees_) + kw["learning_rate"] = float(donor.lr_) + kw["early_stopping_rounds"] = None + kw["replay_donor"] = (donor.trees_, donor.prep_) + man = GradientBoosting(**kw) + man.fit(np.concatenate([d["X_full"], d["X_new"]]), + np.concatenate([d["y_full"], d["y_new"]]), + cat_features=CAT, + sample_weight=np.concatenate([d["w_full"], d["w_new"]])) + _assert_booster_equal(est.model_, man) + np.testing.assert_array_equal(est.predict_raw(d["X_held"]), + man.predict_raw(d["X_held"])) + + +def _pinned_snapshot(est): + """Every refresh-pinned quantity the (c) test compares.""" + m = est.model_ + return { + "n_estimators": est.n_estimators, + "best": est.best_iteration_, + "lr": m.lr_, + "depths": [t.depth for t in m.trees_], + "feats": [np.asarray(t.splits_feat).copy() for t in m.trees_], + "thrs": [np.asarray(t.splits_thr).copy() for t in m.trees_], + "lin_feats": [(None if t.lin_feats is None + else np.asarray(t.lin_feats).copy()) + for t in m.trees_], + "borders": [b.copy() for b in m.prep_.binner_.borders_], + "count": list(m.prep_.count_features_), + "cross_pairs": list(m.cross_pairs), + "valid": list(m.valid_history_), + "n_features": est.n_features_in_, + "names": getattr(est, "feature_names_in_", None), + "ll_sel": getattr(est, "linear_leaves_selected_", None), + "cx_sel": est.cross_features_selected_, + "cx_pairs": est.cross_pairs_, + "temp": getattr(est, "temperature_", None), + "classes": getattr(est, "classes_", None), + } + + +def _assert_pinned_equal(est, snap): + m = est.model_ + assert est.n_estimators == snap["n_estimators"] + assert est.best_iteration_ == snap["best"] + assert m.lr_ == snap["lr"] + assert [t.depth for t in m.trees_] == snap["depths"] + assert len(m.trees_) == len(snap["feats"]) + for t, f, th, lf in zip(m.trees_, snap["feats"], snap["thrs"], + snap["lin_feats"]): + np.testing.assert_array_equal(np.asarray(t.splits_feat), f) + np.testing.assert_array_equal(np.asarray(t.splits_thr), th) + if lf is None or t.lin_feats is None: + assert lf is None and t.lin_feats is None + else: + np.testing.assert_array_equal(np.asarray(t.lin_feats), lf) + for b, b0 in zip(m.prep_.binner_.borders_, snap["borders"]): + np.testing.assert_array_equal(b, b0) + assert list(m.prep_.count_features_) == snap["count"] + assert list(m.cross_pairs) == snap["cross_pairs"] + assert list(m.valid_history_) == snap["valid"] + assert est.n_features_in_ == snap["n_features"] + if snap["names"] is None: + assert getattr(est, "feature_names_in_", None) is None + else: + np.testing.assert_array_equal(est.feature_names_in_, snap["names"]) + assert getattr(est, "linear_leaves_selected_", None) == snap["ll_sel"] + assert est.cross_features_selected_ == snap["cx_sel"] + assert est.cross_pairs_ == snap["cx_pairs"] + if snap["temp"] is not None: + assert est.temperature_ == snap["temp"] + if snap["classes"] is not None: + np.testing.assert_array_equal(est.classes_, snap["classes"]) + + +@pytest.mark.parametrize("kind", ["regression", "binary"]) +def test_refresh_keeps_pinned_state(kind): + """(c) Nothing pinned moves; n_samples_trained_ grows by the new rows.""" + d = _make_data() + if kind == "regression": + est = ChimeraBoostRegressor(random_state=0, n_estimators=60, + cross_features="always", + linear_leaves=True, + store_training_data=True) + y_fit, y_new = d["y_full"], d["y_new"] + else: + est = ChimeraBoostClassifier(random_state=0, n_estimators=60, + cross_features="always", + store_training_data=True) + y_fit = _binary_labels(d["y_full"]) + y_new = _binary_labels(d["y_new"]) + est.fit(d["X_full"], y_fit, cat_features=CAT, sample_weight=d["w_full"]) + assert est.model_.prep_.count_features_ == [6] + est.shap_importances(d["X_full"][:20]) + assert est._shap_importances_cache_ is not None + assert hasattr(est, "expected_value_") + snap = _pinned_snapshot(est) + n0 = est.n_samples_trained_ + est.refresh(d["X_new"], y_new, sample_weight=d["w_new"]) + _assert_pinned_equal(est, snap) + assert est.n_samples_trained_ == n0 + len(y_new) + assert est._shap_importances_cache_ is None + assert not hasattr(est, "expected_value_") + + +def test_refresh_preserves_feature_names(): + """(c) feature_names_in_ survives a refresh (DataFrame fit).""" + pd = pytest.importorskip("pandas") + d = _make_data() + num_cols = [f"n{i}" for i in range(N_NUM)] + + def frame(X): + df = pd.DataFrame(np.asarray(X[:, :N_NUM], dtype=np.float64), + columns=num_cols) + df["hi"] = list(X[:, 6]) + df["lo"] = list(X[:, 7]) + return df + + cols = num_cols + ["hi", "lo"] + est = ChimeraBoostRegressor(random_state=0, n_estimators=20, + store_training_data=True) + est.fit(frame(d["X_full"]), d["y_full"], cat_features=["hi", "lo"]) + assert list(est.feature_names_in_) == cols + est.refresh(frame(d["X_new"]), d["y_new"]) + assert list(est.feature_names_in_) == cols + assert est.n_samples_trained_ == len(d["y_full"]) + len(d["y_new"]) + + +def test_refresh_ignores_post_fit_set_params(): + """refresh() replays the fitted configuration, not later set_params.""" + d = _make_data() + est = ChimeraBoostRegressor(random_state=0, n_estimators=40, + store_training_data=True) + est.fit(d["X_full"][:1500], d["y_full"][:1500], cat_features=CAT) + snap = _snapshot(est, d["X_held"]) + est.set_params(n_estimators=5, learning_rate=0.5, max_bins=16, depth=2, + subsample=0.5, l2_leaf_reg=99.0) + est.refresh(d["X_full"][:0], d["y_full"][:0]) + _assert_snapshot_equal(est, snap, d["X_held"]) + + +def test_chained_refresh_equals_single(): + """(c2) refresh(A) then refresh(B) equals one refresh(A+B), bit for bit.""" + d = _make_data() + + def fresh(): + e = ChimeraBoostRegressor(random_state=0, n_estimators=50, + cross_features="always", linear_leaves=True, + store_training_data=True) + e.fit(d["X_full"], d["y_full"], cat_features=CAT, + sample_weight=d["w_full"]) + return e + + n_a = 120 + e1 = fresh() + e1.refresh(d["X_new"][:n_a], d["y_new"][:n_a], + sample_weight=d["w_new"][:n_a]) + e1.refresh(d["X_new"][n_a:], d["y_new"][n_a:], + sample_weight=d["w_new"][n_a:]) + e2 = fresh() + e2.refresh(d["X_new"], d["y_new"], sample_weight=d["w_new"]) + assert e1.n_samples_trained_ == e2.n_samples_trained_ + _assert_booster_equal(e1.model_, e2.model_) + np.testing.assert_array_equal(e1.predict(d["X_held"]), + e2.predict(d["X_held"])) + + +def test_pickle_round_trip(): + """(d) A refreshed-after-round-trip model equals the refreshed original.""" + d = _make_data() + est = ChimeraBoostRegressor(random_state=0, n_estimators=50, + cross_features="always", linear_leaves=True, + store_training_data=True) + est.fit(d["X_full"], d["y_full"], cat_features=CAT, + sample_weight=d["w_full"]) + rt = pickle.loads(pickle.dumps(est)) + assert rt._training_data_.n_rows == est._training_data_.n_rows + est.refresh(d["X_new"], d["y_new"], sample_weight=d["w_new"]) + rt.refresh(d["X_new"], d["y_new"], sample_weight=d["w_new"]) + _assert_booster_equal(rt.model_, est.model_) + np.testing.assert_array_equal(rt.predict(d["X_held"]), + est.predict(d["X_held"])) + assert rt.n_samples_trained_ == est.n_samples_trained_ + assert est.model_.replay_donor is None + assert rt.model_.replay_donor is None + + plain = ChimeraBoostRegressor(random_state=0, n_estimators=10).fit( + d["X_full"][:500], d["y_full"][:500], cat_features=CAT) + assert plain._training_data_ is None + + +def test_store_training_data_errors(): + """(e) Every slice-1 error message.""" + d = _make_data() + X, y = d["X_full"][:800], d["y_full"][:800] + + for bad in ("yes", 1, None): + with pytest.raises(ValueError, match="store_training_data must be"): + ChimeraBoostRegressor(random_state=0, + store_training_data=bad).fit( + X, y, cat_features=CAT) + + with pytest.raises(NotImplementedError, match="n_ensembles > 1"): + ChimeraBoostRegressor(random_state=0, n_ensembles=2, + store_training_data=True).fit( + X, y, cat_features=CAT) + + with pytest.raises(NotImplementedError, match="n_ensembles > 1"): + ChimeraBoostRegressor(random_state=0, quality=4, + store_training_data=True).fit( + X, y, cat_features=CAT) + + y3 = np.where(y > np.quantile(y, 2 / 3), "c", + np.where(y > np.quantile(y, 1 / 3), "b", "a")) + with pytest.raises(NotImplementedError, match="multiclass classification"): + ChimeraBoostClassifier(random_state=0, + store_training_data=True).fit( + X, y3, cat_features=CAT) + + with pytest.raises(NotImplementedError, match="loss='Quantile'"): + ChimeraBoostRegressor(random_state=0, loss="Quantile", + store_training_data=True).fit( + X, y, cat_features=CAT) + + groups = np.arange(len(y)) % 20 + with pytest.raises(NotImplementedError, match="random_effects=True"): + ChimeraBoostRegressor(random_state=0, random_effects=True, + store_training_data=True).fit( + X, y, cat_features=CAT, groups=groups) + + with pytest.raises(NotFittedError): + ChimeraBoostRegressor().refresh(X[:10], y[:10]) + + plain = ChimeraBoostRegressor(random_state=0, n_estimators=10).fit( + X, y, cat_features=CAT) + with pytest.raises(ValueError, match="store_training_data=False"): + plain.refresh(X[:10], y[:10]) + + +def test_refresh_weight_and_label_errors(): + """(e) Weighted/unweighted mismatches and unseen refresh labels.""" + d = _make_data() + X, y = d["X_full"][:800], d["y_full"][:800] + + w = ChimeraBoostRegressor(random_state=0, n_estimators=10, + store_training_data=True).fit( + X, y, cat_features=CAT, sample_weight=d["w_full"][:800]) + with pytest.raises(ValueError, match="needs sample_weight"): + w.refresh(X[:5], y[:5]) + + u = ChimeraBoostRegressor(random_state=0, n_estimators=10, + store_training_data=True).fit( + X, y, cat_features=CAT) + with pytest.raises(ValueError, match="cannot take sample_weight"): + u.refresh(X[:5], y[:5], sample_weight=np.ones(5)) + + yb = _binary_labels(y) + clf = ChimeraBoostClassifier(random_state=0, n_estimators=10, + store_training_data=True).fit( + X, yb, cat_features=CAT) + bad_y = np.array(["pos", "nope", "neg"], dtype=object) + with pytest.raises(ValueError, match="cannot add classes"): + clf.refresh(X[:3], bad_y) + + +@pytest.mark.parametrize("cls", [ChimeraBoostRegressor, ChimeraBoostClassifier]) +def test_store_param_round_trip(cls): + """(e) clone and get_params round-trip the new parameter.""" + est = cls(store_training_data=True) + assert est.get_params()["store_training_data"] is True + assert clone(est).store_training_data is True + est.set_params(store_training_data=False) + assert est.get_params()["store_training_data"] is False + assert cls().get_params()["store_training_data"] is False + + +@pytest.mark.parametrize("kind", ["regression", "binary"]) +def test_default_fit_unaffected(kind): + """(e) The default fit is unaffected: flag False vs True agree exactly.""" + d = _make_data() + if kind == "regression": + cls, y = ChimeraBoostRegressor, d["y_full"] + else: + cls, y = ChimeraBoostClassifier, _binary_labels(d["y_full"]) + kw = dict(random_state=0, n_estimators=50) + a = cls(store_training_data=False, **kw).fit( + d["X_full"], y, cat_features=CAT, sample_weight=d["w_full"]) + b = cls(store_training_data=True, **kw).fit( + d["X_full"], y, cat_features=CAT, sample_weight=d["w_full"]) + np.testing.assert_array_equal(a.predict(d["X_held"]), b.predict(d["X_held"])) + if kind == "binary": + np.testing.assert_array_equal(a.predict_proba(d["X_held"]), + b.predict_proba(d["X_held"])) + assert a._training_data_ is None + assert b._training_data_ is not None + assert a.n_samples_trained_ == b.n_samples_trained_ + + +def test_refresh_usefulness_smoke(): + """(f) Usefulness smoke (not a gate): refresh lands between 60% and 90%. + + Fit on the first 60% of rows, refresh with the next 30%, and compare + test RMSE against the 60% model and a full refit on 90%. Refresh is + expected between the two, but that is not asserted. + """ + d = _make_data() + X60, y60 = d["X_full"][:1560], d["y_full"][:1560] + X30, y30 = d["X_full"][1560:2340], d["y_full"][1560:2340] + + def rmse(est): + return float(np.sqrt(np.mean((est.predict(d["X_held"]) - d["y_held"]) ** 2))) + + base = ChimeraBoostRegressor(random_state=0, n_estimators=200, + store_training_data=True).fit( + X60, y60, cat_features=CAT) + rmse_base = rmse(base) + base.refresh(X30, y30) + rmse_ref = rmse(base) + full = ChimeraBoostRegressor(random_state=0, n_estimators=200).fit( + d["X_full"][:2340], d["y_full"][:2340], cat_features=CAT) + rmse_full = rmse(full) + print(f"(f) RMSE: 60%={rmse_base:.4f} refreshed={rmse_ref:.4f} " + f"90%={rmse_full:.4f}") + assert np.isfinite([rmse_base, rmse_ref, rmse_full]).all() From 435e866e323ba1b3ed9670a8d4014d995db4bab2 Mon Sep 17 00:00:00 2001 From: Nathan Walker Date: Thu, 24 Sep 2026 22:37:26 -0400 Subject: [PATCH 3/4] Faster linear-leaf fit in replay, bit-identical (#131) The per-leaf ridge sums in _linear_leaf_fit were bound by the largest leaf, which often holds 40-60% of the rows and ran on one thread, after a serial counting sort. Above _SMALL_N rows the kernel now sorts in parallel, gathers each leaf's rows into contiguous buffers once, and spreads the work over (leaf, accumulator) tasks, so a big leaf's sums run on 1 + k threads. Every sum still adds its leaf's rows in increasing row order with the same expressions, so results are bit-identical: the identity snapshot is unchanged (186/186), and new tests compare the kernel against a verbatim copy of the old one. replay_oblivious_tree also uses the parallel leaf assignment. On a 517k-row store: the linear-leaf fit drops from 9.9 to 7.8 ms per tree, a daily refresh from 5.0 to 3.6 s, and a 500k-row fit with linear leaves from 15.9 to 13.4 s, since the grow path shares the kernel. Co-Authored-By: Claude Opus 5.5 --- benchmarks/CAMPAIGN_PLAN.md | 33 ++++ chimeraboost/tree.py | 312 +++++++++++++++++++++++++++--------- tests/test_replay_kernel.py | 210 ++++++++++++++++++++++++ 3 files changed, 483 insertions(+), 72 deletions(-) create mode 100644 tests/test_replay_kernel.py diff --git a/benchmarks/CAMPAIGN_PLAN.md b/benchmarks/CAMPAIGN_PLAN.md index 54a16f7..1c5a64b 100644 --- a/benchmarks/CAMPAIGN_PLAN.md +++ b/benchmarks/CAMPAIGN_PLAN.md @@ -530,6 +530,39 @@ render the `refresh` docstring. verdict: **PASS → PR for the maintainer** (library, tests, docs). Issue #131 stays open for slices 2-7 (REFRESH_PLAN.md). next: issue #113 (random effects, slice 2) per the focus rule. +follow-up (the maintainer, 2026-09-24, on PR #170): "Use a larger dataset +to prove refresh's purpose", then "the prime use case is more of a 'daily +refit' after a day of data comes in", then "Let's put the fast kernel in. +Pretend you're like a maintainer adjusting an initial PR". Measured: +Zurich delays at scale (fit 3.3M rows 259 s, refresh +1.6M 141 s, full +refit on 4.9M 459 s; RMSE 3.0406 / 3.0402 / 3.0402); a daily refresh +(517k-row store + 17k new rows) 5.0 s against a 15.6 s fit and a 0.56 s +predict pass. Per tree 11.0 ms, of which `_linear_leaf_fit` 9.2 ms: a +2.3 ms serial counting sort, then per-leaf sums bound by the largest leaf +(median 38%, max 60% of rows). Muse task `20260924-issue131-fast-replay.md` +on the PR branch: parallel stable sort, contiguous leaf-sorted gather, +parallel over (leaf, accumulator) so one big leaf no longer serializes, +every sum in the same row order. +forecast (kernel): identity snapshot 186/186 and the goldens unchanged; +`_linear_leaf_fit` 9.2 -> <= 3 ms at 517k rows; the daily refresh 5.0 -> +<= 2.5 s; linear-leaf default fits faster by the kernel's share. +muse (exit 0), `tree.py` only: `_linear_leaf_fit` keeps today's code as +the arm for `n <= _SMALL_N`; above it, a parallel stable counting sort, one +parallel gather into leaf-sorted contiguous buffers (`gs`, `hs`, `Xd` as +(k, n)), and parallel (leaf, accumulator-group) tasks, accumulator-major, +each sum in increasing row order with today's expressions; +`replay_oblivious_tree` uses the parallel `donor.apply`. 21 new tests +against a verbatim copy of the old kernel (sizes, empty and tiny leaves, a +60% leaf, NaN bins, k 1 and 6, 1 and 12 threads). +result: bit-identical HIT (identity snapshot 186/186; full suite 1280 +passed, 1 skipped; ruff clean). Speed MISS against the targets: at 517k +rows `_linear_leaf_fit` 9.93 -> 7.81 ms (1.27x, target <= 3), replay per +tree 11.74 -> 9.17 ms, the daily refresh 5.0 -> 3.6 s (target <= 2.5), a +500k-row linear-leaf fit 15.9 -> 13.4 s (1.19x, the grow path shares the +kernel). The gather alone moves ~76 MB per call (3.55 ms, RAM-bound). +Untried lead: gather uint16 bins (6 MB) instead of float64 design values +(24 MB), looking centres up in L1. +verdict (kernel): PASS as a bit-identical speedup, shipped in PR #170. #### I065 2026-09-24 issue #81 (research cascade: dead self-test anchor, stale `ideas.py` flags; BENCH tooling + test, pre-registered) why now: the focus rule's second issue. PR #168 merged (3e02ab1), #84 diff --git a/chimeraboost/tree.py b/chimeraboost/tree.py index 565e31c..99f6845 100644 --- a/chimeraboost/tree.py +++ b/chimeraboost/tree.py @@ -868,95 +868,263 @@ def _linear_leaf_fit(leaf, grad, hess, n_leaves, lin_feats, centers_std, Xb, # Notes ----- - Parallel over leaves and bit-identical to the old serial global scan. A - stable counting sort groups sample indices by leaf in original order, so each - leaf accumulates in exactly the float-add sequence the serial version used -- - a leaf only ever saw its own samples, in increasing i. Thread-count invariant - for the same reason. - - Design values are gathered per sample inside the leaf loop: no (k, n) scratch - matrix, and one parallel region holds the JIT compile cost down. + Two bit-identical arms on `n <= _SMALL_N`: the small arm uses a serial + stable counting sort and parallelizes over leaves; the large arm sorts in + parallel, gathers leaf-sorted contiguous buffers once, and parallelizes + over (leaf, accumulator-group) tasks so a 60%-of-rows leaf spreads over + 1+k threads. Every accumulator sums in increasing row index with the same + expressions, so both arms match the old serial global scan bit for bit. """ n = leaf.shape[0] k = lin_feats.shape[0] d = 1 + k - coef = np.zeros((n_leaves, d)) + if n <= _SMALL_N: + coef = np.zeros((n_leaves, d)) - # Per-leaf grad/hess totals (for the constant fallback) and counts. - counts = np.zeros(n_leaves, dtype=np.int64) - Gtot = np.zeros(n_leaves) - Htot = np.zeros(n_leaves) - for i in range(n): - l = leaf[i] - counts[l] += 1 - Gtot[l] += grad[i] - Htot[l] += hess[i] + # Per-leaf grad/hess totals (for the constant fallback) and counts. + counts = np.zeros(n_leaves, dtype=np.int64) + Gtot = np.zeros(n_leaves) + Htot = np.zeros(n_leaves) + for i in range(n): + l = leaf[i] + counts[l] += 1 + Gtot[l] += grad[i] + Htot[l] += hess[i] - # Stable counting sort: order[start[l]:start[l+1]] = leaf-l samples in - # increasing original index. - start = np.zeros(n_leaves + 1, dtype=np.int64) - for l in range(n_leaves): - start[l + 1] = start[l] + counts[l] + # Stable counting sort: order[start[l]:start[l+1]] = leaf-l samples in + # increasing original index. + start = np.zeros(n_leaves + 1, dtype=np.int64) + for l in range(n_leaves): + start[l + 1] = start[l] + counts[l] - pos = start[:n_leaves].copy() - order = np.empty(n, dtype=np.int64) - for i in range(n): - l = leaf[i] - order[pos[l]] = i - pos[l] += 1 + pos = start[:n_leaves].copy() + order = np.empty(n, dtype=np.int64) + for i in range(n): + l = leaf[i] + order[pos[l]] = i + pos[l] += 1 - # Per-leaf normal equations + solve; leaves are independent. - for l in prange(n_leaves): - if counts[l] == 0: - continue + # Per-leaf normal equations + solve; leaves are independent. + for l in prange(n_leaves): + if counts[l] == 0: + continue - if counts[l] < 2 * d or k == 0: - if Htot[l] > 0.0: - coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) - continue + if counts[l] < 2 * d or k == 0: + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue - Ml = np.zeros((d, d)) - rl = np.zeros(d) - xrow = np.empty(k) + Ml = np.zeros((d, d)) + rl = np.zeros(d) + xrow = np.empty(k) - for q in range(start[l], start[l + 1]): - i = order[q] - h = hess[i] - g = grad[i] + for q in range(start[l], start[l + 1]): + i = order[q] + h = hess[i] + g = grad[i] + + # Standardized design values for this sample; missing bins -> 0. + for j in range(k): + f = lin_feats[j] + v = centers_std[f, Xb[f, i]] + xrow[j] = v if np.isfinite(v) else 0.0 + + Ml[0, 0] += h + rl[0] += -g + for j in range(k): + xj = xrow[j] + Ml[0, 1 + j] += h * xj + Ml[1 + j, 0] += h * xj + rl[1 + j] += -g * xj + for jj in range(k): + Ml[1 + j, 1 + jj] += h * xj * xrow[jj] + + Ml[0, 0] += l2_intercept + for j in range(1, d): + Ml[j, j] += lin_lambda + for j in range(d): + Ml[j, j] += 1e-9 # jitter: keep the solve well-posed + + beta = _solve_small(Ml, rl) + if np.isnan(beta[0]): + # Singular pivot (unreachable given the diagonal ridge + jitter): + # keep the plain constant Newton value rather than a broken slope. + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue + + for j in range(d): + coef[l, j] = lr * beta[j] + else: - # Standardized design values for this sample; missing bins -> 0. + # Parallel stable counting sort: per-chunk leaf counts, prefix sums in + # (leaf, chunk) order, each chunk scatters its rows in increasing index. + n_chunks = (n + 4095) // 4096 + if n_chunks < 1: + n_chunks = 1 + if n_chunks > 64: + n_chunks = 64 + cnt = np.zeros((n_chunks, n_leaves), dtype=np.int64) + for c in prange(n_chunks): + lo = c * n // n_chunks + hi = (c + 1) * n // n_chunks + for i in range(lo, hi): + cnt[c, leaf[i]] += 1 + counts = np.zeros(n_leaves, dtype=np.int64) + for l in range(n_leaves): + s = 0 + for c in range(n_chunks): + s += cnt[c, l] + counts[l] = s + start = np.zeros(n_leaves + 1, dtype=np.int64) + for l in range(n_leaves): + start[l + 1] = start[l] + counts[l] + cursor = np.empty((n_chunks, n_leaves), dtype=np.int64) + for l in range(n_leaves): + s = start[l] + for c in range(n_chunks): + cursor[c, l] = s + s += cnt[c, l] + order = np.empty(n, dtype=np.int64) + for c in prange(n_chunks): + lo = c * n // n_chunks + hi = (c + 1) * n // n_chunks + for i in range(lo, hi): + l = leaf[i] + order[cursor[c, l]] = i + cursor[c, l] += 1 + + # One parallel gather into leaf-sorted contiguous buffers. + gs = np.empty(n, dtype=np.float64) + hs = np.empty(n, dtype=np.float64) + if k > 0: + Xd = np.empty((k, n), dtype=np.float64) + else: + Xd = np.empty((1, n), dtype=np.float64) + for q in prange(n): + i = order[q] + gs[q] = grad[i] + hs[q] = hess[i] for j in range(k): f = lin_feats[j] v = centers_std[f, Xb[f, i]] - xrow[j] = v if np.isfinite(v) else 0.0 + Xd[j, q] = v if np.isfinite(v) else 0.0 - Ml[0, 0] += h - rl[0] += -g - for j in range(k): - xj = xrow[j] - Ml[0, 1 + j] += h * xj - Ml[1 + j, 0] += h * xj - rl[1 + j] += -g * xj + # Task map, accumulator-major so a big leaf's groups spread over threads: + # acc 0 = fallback (Gtot, Htot), acc 1 = full intercept group, + # acc 2+j = full slope group for feature j. + n_full = 0 + n_fallback = 0 + for l in range(n_leaves): + if counts[l] == 0: + continue + if counts[l] >= 2 * d and k > 0: + n_full += 1 + else: + n_fallback += 1 + T = n_fallback + n_full * (1 + k) + task_leaf = np.empty(T, dtype=np.int64) + task_acc = np.empty(T, dtype=np.int64) + p = 0 + for l in range(n_leaves): + if counts[l] > 0 and not (counts[l] >= 2 * d and k > 0): + task_leaf[p] = l + task_acc[p] = 0 + p += 1 + for l in range(n_leaves): + if counts[l] >= 2 * d and k > 0: + task_leaf[p] = l + task_acc[p] = 1 + p += 1 + for j in range(k): + for l in range(n_leaves): + if counts[l] >= 2 * d and k > 0: + task_leaf[p] = l + task_acc[p] = 2 + j + p += 1 + + Gtot = np.zeros(n_leaves) + Htot = np.zeros(n_leaves) + Ml_all = np.zeros((n_leaves, d, d)) + rl_all = np.zeros((n_leaves, d)) + for t in prange(T): + l = task_leaf[t] + acc = task_acc[t] + lo = start[l] + hi = start[l + 1] + if acc == 0: + sG = 0.0 + sH = 0.0 + for q in range(lo, hi): + sG += gs[q] + sH += hs[q] + Gtot[l] = sG + Htot[l] = sH + elif acc == 1: + sG = 0.0 + sH = 0.0 + sM00 = 0.0 + sr0 = 0.0 + for q in range(lo, hi): + g = gs[q] + h = hs[q] + sG += g + sH += h + sM00 += h + sr0 += -g + Gtot[l] = sG + Htot[l] = sH + Ml_all[l, 0, 0] = sM00 + rl_all[l, 0] = sr0 + else: + j = acc - 2 + s01 = 0.0 + s10 = 0.0 + sr = 0.0 + sblk = np.zeros(k) + for q in range(lo, hi): + h = hs[q] + g = gs[q] + xj = Xd[j, q] + hx = h * xj + s01 += hx + s10 += hx + sr += -g * xj + for jj in range(k): + sblk[jj] += hx * Xd[jj, q] + Ml_all[l, 0, 1 + j] = s01 + Ml_all[l, 1 + j, 0] = s10 + rl_all[l, 1 + j] = sr for jj in range(k): - Ml[1 + j, 1 + jj] += h * xj * xrow[jj] - - Ml[0, 0] += l2_intercept - for j in range(1, d): - Ml[j, j] += lin_lambda - for j in range(d): - Ml[j, j] += 1e-9 # jitter: keep the solve well-posed - - beta = _solve_small(Ml, rl) - if np.isnan(beta[0]): - # Singular pivot (unreachable given the diagonal ridge + jitter): - # keep the plain constant Newton value rather than a broken slope. - if Htot[l] > 0.0: - coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) - continue - - for j in range(d): - coef[l, j] = lr * beta[j] + Ml_all[l, 1 + j, 1 + jj] = sblk[jj] + # Same regularization, solve and fallbacks as the small-n arm. + coef = np.zeros((n_leaves, d)) + for l in range(n_leaves): + if counts[l] == 0: + continue + if counts[l] < 2 * d or k == 0: + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue + Ml = np.empty((d, d)) + rl = np.empty(d) + for a in range(d): + rl[a] = rl_all[l, a] + for b in range(d): + Ml[a, b] = Ml_all[l, a, b] + Ml[0, 0] += l2_intercept + for j in range(1, d): + Ml[j, j] += lin_lambda + for j in range(d): + Ml[j, j] += 1e-9 + beta = _solve_small(Ml, rl) + if np.isnan(beta[0]): + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue + for j in range(d): + coef[l, j] = lr * beta[j] return coef @@ -2146,7 +2314,7 @@ def replay_oblivious_tree(donor, Xb, grad, hess, l2, lr, linear_leaves=False, return (ObliviousTree(sf, st, np.zeros(1), np.zeros(0)), np.zeros(Xb.shape[1], dtype=np.int64)) - leaf = _assign_leaves(Xb, sf, st) + leaf = donor.apply(Xb) n_leaves = 1 << len(sf) values = _leaf_values(leaf, grad, hess, n_leaves, l2, lr) diff --git a/tests/test_replay_kernel.py b/tests/test_replay_kernel.py new file mode 100644 index 0000000..ade237a --- /dev/null +++ b/tests/test_replay_kernel.py @@ -0,0 +1,210 @@ +"""Bit-identity guards for the fast replay kernel (issue131, I066). + +``_linear_leaf_fit_ref`` is a frozen copy of ``tree._linear_leaf_fit`` as of +2026-09-24, before the parallel rewrite. Every test asserts exact +``np.array_equal`` against it: any summation-order or expression change fails +in the last ulp. +""" +import numba +import numpy as np +import pytest +from numba import njit, prange + +from chimeraboost.tree import ( + ObliviousTree, + _assign_leaves, + _leaf_values, + _linear_leaf_fit, + _solve_small, + replay_oblivious_tree, +) + +DEFAULT_THREADS = numba.get_num_threads() + + +@njit(cache=False, parallel=True) +def _linear_leaf_fit_ref(leaf, grad, hess, n_leaves, lin_feats, centers_std, Xb, # noqa: C901 -- frozen reference kernel + l2_intercept, lin_lambda, lr): + """Frozen copy of tree._linear_leaf_fit (2026-09-24).""" + n = leaf.shape[0] + k = lin_feats.shape[0] + d = 1 + k + coef = np.zeros((n_leaves, d)) + + counts = np.zeros(n_leaves, dtype=np.int64) + Gtot = np.zeros(n_leaves) + Htot = np.zeros(n_leaves) + for i in range(n): + l = leaf[i] + counts[l] += 1 + Gtot[l] += grad[i] + Htot[l] += hess[i] + + start = np.zeros(n_leaves + 1, dtype=np.int64) + for l in range(n_leaves): + start[l + 1] = start[l] + counts[l] + + pos = start[:n_leaves].copy() + order = np.empty(n, dtype=np.int64) + for i in range(n): + l = leaf[i] + order[pos[l]] = i + pos[l] += 1 + + for l in prange(n_leaves): + if counts[l] == 0: + continue + + if counts[l] < 2 * d or k == 0: + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue + + Ml = np.zeros((d, d)) + rl = np.zeros(d) + xrow = np.empty(k) + + for q in range(start[l], start[l + 1]): + i = order[q] + h = hess[i] + g = grad[i] + + for j in range(k): + f = lin_feats[j] + v = centers_std[f, Xb[f, i]] + xrow[j] = v if np.isfinite(v) else 0.0 + + Ml[0, 0] += h + rl[0] += -g + for j in range(k): + xj = xrow[j] + Ml[0, 1 + j] += h * xj + Ml[1 + j, 0] += h * xj + rl[1 + j] += -g * xj + for jj in range(k): + Ml[1 + j, 1 + jj] += h * xj * xrow[jj] + + Ml[0, 0] += l2_intercept + for j in range(1, d): + Ml[j, j] += lin_lambda + for j in range(d): + Ml[j, j] += 1e-9 + + beta = _solve_small(Ml, rl) + if np.isnan(beta[0]): + if Htot[l] > 0.0: + coef[l, 0] = -lr * Gtot[l] / (Htot[l] + l2_intercept) + continue + + for j in range(d): + coef[l, j] = lr * beta[j] + + return coef + + +def _make_inputs(n, n_leaves, k, seed): + """Random kernel inputs covering empty, tiny (<2d) and 60%+ leaves.""" + rng = np.random.default_rng(seed) + d = 1 + k + n_features, max_bins = 8, 32 + + leaf = np.empty(n, dtype=np.int64) + n_big = max(int(0.6 * n), 2 * d) + if n_big + 5 <= n: + n_s1, n_s2 = 2, 3 + elif n_big + 2 <= n: + n_s1, n_s2 = 1, 1 + else: + n_big = n - 2 + n_s1, n_s2 = 1, 1 + leaf[:n_big] = 0 + leaf[n_big:n_big + n_s1] = 1 + leaf[n_big + n_s1:n_big + n_s1 + n_s2] = 2 + rest = n - (n_big + n_s1 + n_s2) + if rest > 0: + leaf[n_big + n_s1 + n_s2:] = rng.integers(3, n_leaves - 2, rest) + + grad = rng.standard_normal(n) + hess = rng.random(n) + 0.1 + hess[rng.random(n) < 0.05] = 0.0 + lin_feats = np.arange(k, dtype=np.int64) + centers_std = rng.standard_normal((n_features, max_bins)) + centers_std[rng.random((n_features, max_bins)) < 0.05] = np.nan + Xb = rng.integers(0, max_bins, size=(n_features, n)).astype(np.uint16) + return leaf, grad, hess, lin_feats, centers_std, Xb + + +@pytest.mark.parametrize("n,n_leaves", [(16, 8), (2000, 16), (200_000, 64)]) +@pytest.mark.parametrize("k", [1, 6]) +@pytest.mark.parametrize("n_threads", [1, DEFAULT_THREADS]) +def test_linear_leaf_fit_matches_ref(n, n_leaves, k, n_threads): + leaf, grad, hess, lin_feats, cs, Xb = _make_inputs(n, n_leaves, k, 1000 + n + k) + prev = numba.get_num_threads() + numba.set_num_threads(n_threads) + try: + ref = _linear_leaf_fit_ref(leaf, grad, hess, n_leaves, lin_feats, cs, + Xb, 1.0, 1.0, 0.1) + got = _linear_leaf_fit(leaf, grad, hess, n_leaves, lin_feats, cs, + Xb, 1.0, 1.0, 0.1) + finally: + numba.set_num_threads(prev) + assert got.shape == ref.shape + assert np.array_equal(got, ref) + + +def _make_donor(n_features, max_bins, depth, k, seed): + rng = np.random.default_rng(seed) + sf = rng.integers(0, n_features, depth).astype(np.int64) + st = rng.integers(0, max_bins - 1, depth).astype(np.int64) + n_leaves = 1 << depth + vals = rng.standard_normal(n_leaves) + lin_feats = np.arange(k, dtype=np.int64) + lin_coef = np.zeros((n_leaves, 1 + k)) + return ObliviousTree(sf, st, vals, np.zeros(depth), + lin_feats=lin_feats, lin_coef=lin_coef, + centers_std=None), n_leaves + + +@pytest.mark.parametrize("n,depth", [(2000, 4), (200_000, 6)]) +@pytest.mark.parametrize("k", [1, 6]) +@pytest.mark.parametrize("n_threads", [1, DEFAULT_THREADS]) +def test_replay_matches_ref(n, depth, k, n_threads): + rng = np.random.default_rng(77 + n + k) + n_features, max_bins = 8, 32 + donor, n_leaves = _make_donor(n_features, max_bins, depth, k, 5) + Xb = rng.integers(0, max_bins, size=(n_features, n)).astype(np.uint16) + grad = rng.standard_normal(n) + hess = rng.random(n) + 0.1 + cs = rng.standard_normal((n_features, max_bins)) + cs[rng.random((n_features, max_bins)) < 0.05] = np.nan + + prev = numba.get_num_threads() + numba.set_num_threads(n_threads) + try: + tree, leaf = replay_oblivious_tree(donor, Xb, grad, hess, 1.0, 0.1, + linear_leaves=True, centers_std=cs, + linear_lambda=1.0) + leaf_ref = _assign_leaves(Xb, donor.splits_feat, donor.splits_thr) + values_ref = _leaf_values(leaf_ref, grad, hess, n_leaves, 1.0, 0.1) + coef_ref = _linear_leaf_fit_ref(leaf_ref, grad, hess, n_leaves, + donor.lin_feats, cs, Xb, 1.0, 1.0, 0.1) + finally: + numba.set_num_threads(prev) + assert np.array_equal(leaf, leaf_ref) + assert np.array_equal(tree.values, values_ref) + assert np.array_equal(tree.lin_coef, coef_ref) + + +def test_replay_constant_matches_ref(): + """Linear off: replay still matches (covers the assign swap alone).""" + rng = np.random.default_rng(9) + n, n_features, max_bins = 200_000, 8, 32 + donor, n_leaves = _make_donor(n_features, max_bins, 6, 2, 6) + Xb = rng.integers(0, max_bins, size=(n_features, n)).astype(np.uint16) + grad = rng.standard_normal(n) + hess = rng.random(n) + 0.1 + tree, leaf = replay_oblivious_tree(donor, Xb, grad, hess, 1.0, 0.1) + leaf_ref = _assign_leaves(Xb, donor.splits_feat, donor.splits_thr) + values_ref = _leaf_values(leaf_ref, grad, hess, n_leaves, 1.0, 0.1) + assert np.array_equal(leaf, leaf_ref) + assert np.array_equal(tree.values, values_ref) From c139c77567444b0556abf54d7be28f1fd84f172a Mon Sep 17 00:00:00 2001 From: Nathan Walker Date: Thu, 24 Sep 2026 22:46:45 -0400 Subject: [PATCH 4/4] Changelog and plan: the faster linear-leaf fit, and refresh at scale (#131) Co-Authored-By: Claude Opus 5.5 --- CHANGELOG.md | 9 ++++++++- benchmarks/REFRESH_PLAN.md | 8 ++++++++ 2 files changed, 16 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index c94bd0d..521d6c9 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -15,9 +15,16 @@ The format follows [Keep a Changelog](https://keepachangelog.com/). with the tree count, learning rate and splits pinned; refreshing with no new rows reproduces the model bit for bit. Bagged models, multiclass, `loss="Quantile"` and `random_effects=True` are not supported yet. The - default fit is unchanged. + default fit's results are unchanged. On a 3.3M-row regression, folding in + a day of new rows took 70 s against 232 s for a full refit, at the same + accuracy. ### Changed +- **Fits with linear leaves are faster, with identical results.** The + per-leaf linear fit used to sum each leaf on one thread, and one leaf + often holds 40-60% of the rows. Large leaves now spread their sums over + several threads, each sum still adding the rows in the same order. A + 500k-row fit went from 15.9 to 13.4 s, and a refresh from 5.0 to 3.6 s. - **The quantile benchmark's NGBoost opponent now uses the RoNGBa settings** (Ren, Sun and Wu 2019; #163): trees of up to 31 leaves, a learning rate of 0.04 and at most 500 rounds, in place of NGBoost's stock diff --git a/benchmarks/REFRESH_PLAN.md b/benchmarks/REFRESH_PLAN.md index 421cd71..7058a2d 100644 --- a/benchmarks/REFRESH_PLAN.md +++ b/benchmarks/REFRESH_PLAN.md @@ -119,3 +119,11 @@ Still open, for their slices: 15 configurations, chained refreshes compose, identity snapshot unchanged. Smoke: a 60% model refreshed with 30% more rows recovered ~67% of a full 90% refit's RMSE gain. Awaiting the maintainer's merge. +- 2026-09-24 (maintainer follow-up on PR #170): at scale, a daily update + (Zurich delays, 3.3M-row model, 50k new rows) refreshes in 70 s against + a 232 s full refit and a 216 s original fit, same RMSE (3.04146 vs + 3.04145); a predict pass over the store is 9.0 s, so a refresh is ~8 + predict passes, not the spec's ~1. Most of it was the linear-leaf fit, + now faster and bit-identical (PR #170, commit 435e866): 9.9 -> 7.8 ms per + tree at 517k rows. Lead for more: gather uint16 bins instead of float64 + design values in `_linear_leaf_fit` (the gather is ~76 MB per call).