graphtoolbox.aggregation.aggregation

Online and batch expert aggregation (Opera-style robust prediction).

This module provides Aggregation, a self-contained implementation of sequential expert aggregation in the spirit of the R package opera (Gaillard & Goude). Given a matrix of expert forecasts and the corresponding observations, the class combines the experts into a single forecast whose weights adapt over time so as to track, and provably compete with, the best expert (or the best fixed convex combination) in hindsight.

Implemented aggregation rules

• "EWA"Exponentially Weighted Average (Hedge). Fixed learning rate,

  calibrated automatically when observations are available.

• "MLpol"ML-Poly with the polynomial potential of Gaillard, Stoltz &

  van Erven (2014). Fully parameter-free (per-expert adaptive learning rate); the robust default.

• "BOA"Bernstein Online Aggregation (Wintenberger, 2017). Second-order

  update with an exponential potential.

• "uniform": Static uniform mean of the experts (baseline).

• "best"Oracle constant weight on the single best expert over the whole

  horizon (baseline; requires observations).

References

• Cesa-Bianchi, N. & Lugosi, G. (2006). Prediction, Learning, and Games.

• Gaillard, P., Stoltz, G. & van Erven, T. (2014). A second-order bound with excess losses. COLT.

• Wintenberger, O. (2017). Optimal learning with Bernstein online aggregation. Machine Learning.

• Gaillard, P. & Goude, Y. opera: Online Prediction by Expert Aggregation (R package).

class graphtoolbox.aggregation.aggregation.Aggregation(model: str = 'MLpol', loss: str = 'square', learning_rate: float | None = None, gradient: bool = True, prior: ndarray | Tensor | None = None)

Bases: object

Sequential aggregation of expert forecasts.

The estimator consumes a matrix of expert predictions of shape [T, K] (T time steps, K experts) and, when observations y of shape [T] are available, produces an aggregated forecast whose weights are updated online. The full trajectory of weights is exposed in weights_, and the final convex weight vector in coefficients_, so that the fitted mixture can be reused as a static combiner on out-of-sample experts through predict().

Parameters:

• model ({"MLpol", "EWA", "BOA", "uniform", "best"}, default "MLpol") – Aggregation rule. "MLpol" is parameter-free and is a robust default.

• loss ({"square", "absolute", "percentage"}, default "square") – Loss used both to evaluate experts and to drive the weight updates.

• learning_rate (float or None, default None) – Learning rate eta for "EWA" and "BOA". When None and observations are provided, it is calibrated automatically by replaying the sequence over a logarithmic grid and keeping the value with the lowest cumulative loss (see learning_rate_). Ignored by "MLpol", "uniform" and "best".

• gradient (bool, default True) – If True, updates use the linearised (gradient) pseudo-loss g_t * x_{t,k}, as in opera. This guarantees bounded, convex regrets. If False, updates use the raw per-expert losses.

• prior (numpy.ndarray or None, default None) – Optional prior weights of shape [K] (non-negative, summing to one). Defaults to the uniform prior.

weights_

Weight trajectory of shape [T, K]; row t holds the weights used to form the forecast at time t (before observing y_t).

Type:

numpy.ndarray

coefficients_

Final weight vector of shape [K] (the weights that would be applied to the next, unseen time step).

Type:

numpy.ndarray

prediction_

Aggregated forecast of shape [T] from the last online run.

Type:

numpy.ndarray

experts_loss_

Per-expert cumulative loss of shape [K].

Type:

numpy.ndarray

loss_

Cumulative loss of the aggregated forecast.

Type:

float

learning_rate_

Learning rate actually used (relevant for "EWA" / "BOA").

Type:

float

Examples

>>> import numpy as np
>>> rng = np.random.default_rng(0)
>>> T, K = 500, 4
>>> y = np.sin(np.linspace(0, 20, T))
>>> experts = y[:, None] + rng.normal(0, [0.1, 0.3, 0.6, 1.0], size=(T, K))
>>> agg = Aggregation(model="MLpol").run(experts, y)
>>> agg.loss_ < (experts - y[:, None]).__pow__(2).mean(0).min() * T
True

run(experts: ndarray | Tensor, y: ndarray | Tensor, block_size: int = 1) → Aggregation

Run the online aggregation over a full sequence.

Parameters:

• experts (array-like of shape [T, K]) – Expert forecasts (numpy array or torch tensor).

• y (array-like of shape [T]) – Observations.

• block_size (int, default 1) – Forecast-commitment granularity, following opera’s predict / update split. With block_size=1 the weights are updated at every step. With block_size=h the forecast for a whole block of h steps is committed with the weights known at the block’s start; once the block is observed the internal weights still advance one step at a time within it, as opera loops over the block inside update. For a day-ahead forecast served in daily batches of 48 half-hours, set block_size=48.

Returns:

self, with weights_, prediction_, coefficients_, experts_loss_ and loss_ set.

Return type:

Aggregation

fit(experts: ndarray | Tensor, y: ndarray | Tensor, block_size: int = 1) → Aggregation

Alias of run() that returns self (scikit-style).

predict(experts: ndarray | Tensor) → ndarray

Apply the fitted static weights (coefficients_) to new experts.

Use this for out-of-sample combination once the mixture has been fitted on a history via run(). For time-varying online weights on a fresh sequence, call run() again with the new observations.

Parameters:

experts (array-like of shape [T, K]) – New expert forecasts.

Returns:

Aggregated forecast of shape [T].

Return type:

numpy.ndarray

reset(n_experts: int) → Aggregation

Initialise streaming state for n_experts experts.

After reset, alternate partial_fit() (to get the forecast for the current step) and update() (to feed back the observation).

Parameters:

n_experts (int) – Number of experts K.

Returns:

self.

Return type:

Aggregation

partial_fit(expert_row: ndarray | Tensor) → float

Return the aggregated forecast for the current step.

Parameters:

expert_row (array-like of shape [K]) – The experts’ forecasts for the current time step.

Returns:

Aggregated forecast.

Return type:

float

update(y_t: float) → None

Feed back the observation for the step served by partial_fit().

Parameters:

y_t (float) – Observation for the current time step.

property streaming_weights_: ndarray

Weight trajectory accumulated through the streaming interface.

summary() → Dict[str, float]

Return a compact performance summary of the fitted mixture.

Returns:

Keys: model, loss (mean aggregated loss), best_expert_loss (mean loss of the best single expert), mean_expert_loss (average over experts of their mean loss) and learning_rate.

Return type:

dict