StreamingGAM(
formula,
*,
family=None,
decay=1.0,
smoothing_params=None,
)
Fits a GAM incrementally by accumulating weighted sufficient statistics (X'WX and X'Wz) across data batches, rather than storing and re-fitting on the full dataset each time. This makes it suitable for data arriving continuously or in a stream too large to hold in memory at once: the memory footprint depends only on the number of basis coefficients p (via a p x p matrix), not on the number of observations seen. The model structure (formula, basis dimensions, penalties) is fixed at initialisation from a small pilot batch, and subsequent partial_fit() calls add data without storing the raw observations.
Two update modes are supported: accumulate mode (decay=1.0, the default), where every batch contributes equally regardless of when it arrived, and sliding-window mode (decay < 1.0), where older batches’ contributions to the accumulated statistics are exponentially downweighted, allowing the model to adapt to a slowly drifting data-generating process.
Use StreamingGAM for large or continuously-arriving datasets where holding the full data in memory (as ordinary GAM does) is impractical, or where the underlying relationship may drift over time and old data should be forgotten.
Parameters
formula: str | Formula
-
Model formula (e.g. "y ~ s(x1) + s(x2)").
family: Family | None = None
-
Response distribution family. Defaults to Gaussian().
decay: float = 1.0
-
Exponential decay factor for sliding window. 1.0 (default) means no decay (accumulate all data equally). Values less than 1.0 downweight older batches’ sufficient statistics by a factor of decay every time a new batch arrives.
smoothing_params: list[float] | None = None
-
Fixed smoothing parameters. If
None, estimated from the pilot batch (via a one-off ordinary GAM fit with REML) and optionally re-estimated later via solve(reestimate_smoothing=True).
Notes
Each partial_fit() call treats the incoming batch as one step of iteratively reweighted least squares: given the current coefficients, it computes working responses z and IRLS weights W for the batch, forms the batch’s contribution to the weighted normal equations,
X_{\text{batch}}' W X_{\text{batch}}, \qquad X_{\text{batch}}' W z_{\text{batch}},
and adds these to the running totals (after applying the decay factor, if any, to the existing totals). Calling solve() then solves the penalized normal equations
\left(\sum_{\text{batches}} X'WX + S_\lambda\right) \beta = \sum_{\text{batches}} X'Wz
via a Cholesky factorization, where S_\lambda = \sum_j \lambda_j S_j is the weighted sum of penalty matrices. Because only the accumulated p x p matrix X'WX and length-p vector X'Wz are retained, this scales to arbitrarily many observations at fixed memory cost in p.
Examples
import numpy as np
from whittaker.streaming import StreamingGAM
rng = np.random.default_rng(0)
model = StreamingGAM("y ~ s(x)")
for _ in range(5):
x = rng.uniform(0, 1, 200)
y = np.sin(2 * np.pi * x) + rng.normal(scale=0.2, size=200)
model.partial_fit({"x": x, "y": y})
model.solve()
print(model.summary())
StreamingGAM summary
============================================================
Formula: y ~ s(x)
Family: Gaussian
Decay: 1.0
N obs: 1000
N batches: 5
EDF total: 9.6
Scale: 0.3044
Deviance: 301.5
Smooth terms:
s(x): edf=8.6
Snapshots: 1
Attributes
|
Name
|
Description
|
|
coefficients
|
Current coefficient vector from the most recent solve() (or the pilot fit).
|
|
edf_total
|
Total effective degrees of freedom from the most recent solve().
|
|
family
|
Response distribution family used for the link function and variance model.
|
|
formula
|
Parsed model formula this streaming GAM was constructed with.
|
|
is_initialised
|
Whether the model structure has been built from a pilot batch.
|
|
is_solved
|
Whether coefficients are available from an initial fit or a solve() call.
|
|
n_batches
|
Total number of batches passed to partial_fit() so far.
|
|
n_obs
|
Accumulated observation count across all ingested batches.
|
|
scale
|
Estimated dispersion (scale) parameter from the most recent solve().
|
|
smoothing_params
|
Smoothing parameters currently in use, one per penalized smooth term.
|
coefficients
Current coefficient vector from the most recent solve() (or the pilot fit).
Returns a copy, so mutating the returned array does not affect the model.
edf_total
Total effective degrees of freedom from the most recent solve().
Computed as the trace of the hat matrix implied by the accumulated X'WX and the current smoothing parameters.
family
Response distribution family used for the link function and variance model.
is_initialised
Whether the model structure has been built from a pilot batch.
Becomes True after the first partial_fit() call, once the model matrix, basis dimensions, and penalties have been established.
is_solved
Whether coefficients are available from an initial fit or a solve() call.
predict(), coefficients, and other solved-state accessors raise RuntimeError when this is False.
n_batches
Unlike n_obs, this count is never decayed.
n_obs
Accumulated observation count across all ingested batches.
Under sliding-window mode (decay < 1.0), this count is itself decayed at each partial_fit() call, so it reflects an effective rather than a raw cumulative count.
scale
Estimated dispersion (scale) parameter from the most recent solve().
For families with a known scale (e.g. Poisson, binomial), this is fixed at 1.0; otherwise it is estimated from the accumulated deviance and edf_total.
smoothing_params
Smoothing parameters currently in use, one per penalized smooth term.
smoothing_params: list[float]
These come from the constructor’s smoothing_params argument if fixed, or from the pilot fit / most recent solve(reestimate_smoothing=True) call otherwise. Returns a copy, so mutating the returned list does not affect the model.
Methods
|
Name
|
Description
|
|
partial_fit()
|
Ingest a batch of data.
|
|
predict()
|
Predict on new data.
|
|
reset()
|
Reset accumulated statistics (keep model structure).
|
|
should_refit()
|
Check if a refit (solve) is recommended.
|
|
smoothing_history()
|
Return history of snapshots taken at each solve() call.
|
|
solve()
|
Solve for coefficients from accumulated sufficient statistics.
|
|
summary()
|
Build a human-readable text summary of the streaming GAM’s current state.
|
partial_fit()
On the first call, builds the model matrix structure (basis dimensions, penalties, knot locations) from this batch via a pilot GAM fit, and initializes the coefficients and smoothing parameters from that fit (unless fixed smoothing parameters were supplied). Subsequent calls reuse this fixed structure: the batch’s working response and IRLS weights are computed using the current coefficients (from the most recent solve(), or the pilot fit), and its contribution to the running X'WX / X'Wz sufficient statistics is added (after applying the decay factor to existing totals, if decay < 1.0).
Parameters
data: InputData
-
Column-oriented batch data containing the response and all covariates in
formula.
Returns
StreamingGAM
-
Returns
self for method chaining.
predict()
predict(
new_data,
*,
se=False,
)
Builds the prediction design matrix using the fixed basis structure established at initialisation, and forms the linear predictor and (via the family’s inverse link) fitted values from the current coefficients. Requires that solve() has been called at least once.
Parameters
new_data: InputData
-
Column-oriented covariate data.
se: bool = False
-
If
True, compute standard errors on the linear predictor scale from the Bayesian posterior covariance implied by the current accumulated X'WX and smoothing parameters.
reset()
Reset accumulated statistics (keep model structure).
should_refit()
Check if a refit (solve) is recommended.
should_refit(
*,
min_batches=10,
)
Returns True when enough new data has accumulated since the last solve to justify re-estimating coefficients.
Parameters
min_batches: int = 10
-
Minimum batches since last solve before recommending refit.
smoothing_history()
Return history of snapshots taken at each solve() call.
solve()
Solve for coefficients from accumulated sufficient statistics.
solve(
*,
reestimate_smoothing=False,
)
Forms the penalized normal equations (X'WX + S_lambda) beta = X'Wz from the currently accumulated statistics and the current smoothing parameters, and solves them via a Cholesky factorization (with a small ridge added for numerical stability). Also updates effective degrees of freedom, the scale estimate, and appends a StreamingSnapshot to the fit history.
Parameters
reestimate_smoothing: bool = False
-
If
True, re-estimate smoothing parameters from the current sufficient statistics using a GCV line search over each smoothing parameter in turn (coordinate-wise). Default False uses the current (pilot or fixed) smoothing parameters. Ignored (has no effect) if the model was constructed with fixed smoothing_params.
Returns
StreamingGAM
-
Returns
self for method chaining.
summary()
Build a human-readable text summary of the streaming GAM’s current state.
Reports the formula, family, decay factor, accumulated observation and batch counts, total EDF, scale, and accumulated deviance, followed by a per-smooth-term EDF breakdown and the number of solve() snapshots recorded so far.
Returns
str
-
Multi-line summary text.