Guide

How the pieces fit: the Dist object, transforms, conjugation, and ensembles — everything returns a distribution, everything composes.

The formal write-up and the full benchmark study are in the papers.

Architecture: transforms all the way down

Every model in skaters is a chain of bijective transforms with a distributional leaf at the bottom:

If you prefer a football metaphor to a formula, the tiki-taka page explains the same idea as moving the ball into open space and back, in your choice of language and team.

$y \;\xrightarrow{T_1}\; y' \;\xrightarrow{T_2}\; y'' \;\xrightarrow{\cdots}\; \text{leaf} \;\rightarrow\; \hat{D}$

The leaf estimates $\hat{D} = \mathcal{N}(0, \hat\sigma^2)$ from residuals via Welford's algorithm. Predictions in the original space come from inverting the transform chain: $\hat{D}_{\text{original}} = T_1^{-1}\bigl(T_2^{-1}\bigl(\cdots(\hat{D})\bigr)\bigr).$

Every node returns a list[Dist]. There is no separate “point forecast” vs “uncertainty” — both are aspects of the same $\hat{D}$. An EMA doesn't “predict”; it strips off a running level $\ell_t$, leaving simpler residuals $\varepsilon_t = y_t - \ell_t$, and the prediction is whatever the leaf's estimate becomes when run back through the inverse.

Notation: hatted quantities are online estimates using information up to $y_{t-1}$ (the pre-update state as $y_t$ arrives). In drift, $\hat\mu_t$ is the smoothed mean increment (the drift estimate, shrunk toward zero by $\lambda$); in standardize, $\hat\mu_t$ and $\hat\sigma_t$ are a running mean and scale of the series itself; in markov_drift, $\hat\mu^{M}_t$ is the one-step conditional mean of the embedded regime forecaster and $g_t \in [-\tfrac12, \tfrac12]$ its learned gate.

Transforms

TransformForwardUse case
ema_transform(α)$y'_t = y_t - \ell_t$Remove level
difference()$y'_t = y_t - y_{t-1}$Random walk
drift(α, λ)$y'_t = \Delta y_t - \hat\mu_t$Random walk + drift
holt_linear(α, β)$y'_t = y_t - (\ell_t + b_t)$Level + trend (Holt 1957)
ar(p)$y'_t = y_t - \sum \hat\phi_j y_{t-j}$Autoregression (online RLS)
fractional_difference(d)$y'_t = (1-B)^d y_t$Long memory
standardize(α)$y'_t = (y_t - \hat\mu_t)/\hat\sigma_t$Remove scale
garch(ω, α, β)$y'_t = y_t / \hat\sigma_t$Volatility clustering
seasonal_difference(s)$y'_t = y_t - y_{t-s}$Periodicity
yeo_johnson(λ)signed Box–Cox coordinateLearn the coordinate (log / root / linear)
markov_drift(η, λ)$y'_t = y_t - g_t\,\hat\mu^{M}_t$Regime response (machine-generated; opt-in)

markov_drift is opt-in and composes like any other transform. Standalone, or as an extra candidate beside the standard population (the benchmarked form, +61% one-step log-likelihood wins over plain laplace on the FRED daily universe, worst series −0.012 nats):

# standalone: markov signal wrapping the standard leaf
from skaters.conjugate import conjugate
from skaters.leaf import leaf
from skaters.markov import markov_drift

f = conjugate(leaf(k=1), markov_drift(), k=1)

# beside laplace's population (see benchmarks/markov_forecaster.py):
from skaters.api import _build_candidates, _objective_leaf
from skaters.sticky import sticky
from skaters.tails import gpdtails
from skaters.terminal import terminal_leaf_ensemble

candidates, depths, _ = _build_candidates(1)
candidates.append(conjugate(leaf(k=1), markov_drift(), k=1))
depths.append(2)
f = terminal_leaf_ensemble(candidates, k=1, leaf_fn=_objective_leaf("crps", 0.03),
                           learning_rate=0.8, complexity_penalty=0.005,
                           depths=depths, max_components=20, forget=0.99)
f = gpdtails(sticky(f, k=1), k=1)

Conjugation

Transforms compose via conjugation. Given a transform $T$ and a skater $f$, $f_{\text{conj}}(y) = T^{-1}\bigl(f\bigl(T(y)\bigr)\bigr).$ The whole predictive distribution is mapped back through the inverse, so multi-step uncertainty accumulates correctly and the coordinate choice is handled coherently rather than by a biased point back-transform.

from skaters import conjugate, ema, difference, standardize

# diff removes trend, then EMA predicts the differenced series
f = conjugate(ema(alpha=0.1, k=3), difference(), k=3)

# Chain: standardize, then difference, then EMA
f = conjugate(
    conjugate(ema(alpha=0.1, k=3), difference(), k=3),
    standardize(),
    k=3,
)
# canonical name: std|diff|ema_t|leaf

Ensembles

Precision-weighted (MSE)

Weights by $w_i \propto 1/\text{MSE}_i$ where $\text{MSE} = \text{bias}^2 + \text{variance}$.

from skaters import precision_weighted_ensemble, ema

f = precision_weighted_ensemble([
    ema(alpha=0.05, k=1),
    ema(alpha=0.2, k=1),
], k=1)

Bayesian (log-likelihood with XGBoost-inspired regularisation)

Each model accumulates a log-weight updated at every observation: $\log w_i \mathrel{+}= \eta \cdot \log p_i(y_t) - \lambda \cdot d_i,$ where $\eta$ is the learning rate (shrinkage), $\lambda$ is the complexity penalty, and $d_i$ is the model's depth. Predictions are combined via Dist.combine with softmax weights.

from skaters import bayesian_ensemble, ema

f = bayesian_ensemble(
    [ema(alpha=0.05, k=1), ema(alpha=0.2, k=1)],
    k=1,
    learning_rate=0.5,       # eta: prevents over-concentrating
    complexity_penalty=0.02, # lambda: penalises deeper chains
    depths=[1, 1],
)

Adaptive search

Beam search over the transform grammar. Grows the candidate population online: expand top performers with new transforms, replay recent history to warm-start, prune losers.

from skaters import search

f = search(
    k=1,
    expand_interval=100,   # expand top performers every 100 obs
    max_depth=3,           # max transform chain depth
    replay_buffer=500,     # warm-start new candidates on recent history
    max_pool=30,
)

Design