calibrate_sigma()

Find the ELF bandwidth sigma that minimises out-of-sample ELF loss.

Usage

Source

calibrate_sigma(
    formula,
    data,
    tau=0.5,
    *,
    n_folds=5,
    sigma_values=None,
    method="GCV",
    seed=None,
)

The quantile family used by QuantileGAM and GAM with quantile loss replaces the non-differentiable check function with a smooth “extended log-F” (ELF) surrogate controlled by a bandwidth sigma. Too large a sigma over-smooths the check loss and biases the fitted quantile; too small a sigma makes the surrogate nearly non-differentiable again and can destabilize IRLS. calibrate_sigma selects sigma empirically by K-fold cross-validation: for each candidate value, the model is fit on K - 1 folds, predictions are made on the held-out fold, and the true (non-smoothed) pinball loss

\rho_\tau(y - \hat q_\tau) = (y - \hat q_\tau)\,(\tau - \mathbb{1}[y < \hat q_\tau])

is accumulated across folds. The sigma minimizing total out-of-sample pinball loss is refined with a second, finer grid search around the best value from the coarse grid.

Parameters

formula: str

GAM formula string, e.g. "y ~ s(x)".

data: InputData

Column-oriented data dict.

tau: float = 0.5

Target quantile level in (0, 1).

n_folds: int = 5

Number of CV folds.

sigma_values: NDArray | list[float] | None = None

Candidate sigma values to evaluate. If None, a log-spaced grid of 10 values from 0.01 * sd(y) to 2 * sd(y) is used, followed by a refinement grid around the best value. If an explicit grid is passed, no refinement step is performed.

method: str = "GCV"

Smoothing parameter selection method used when fitting each candidate model ("GCV", "REML", "ML").

seed: int | None = None
Random seed for fold assignment.

Returns

float
Calibrated sigma value (the one minimizing CV pinball loss).

Examples

import numpy as np
from whittaker.calibration import calibrate_sigma

rng = np.random.default_rng(0)
n = 400
x = rng.uniform(0, 1, n)
y = np.sin(2 * np.pi * x) + rng.normal(scale=0.2 + 0.3 * x, size=n)

best_sigma = calibrate_sigma("y ~ s(x)", {"x": x, "y": y}, tau=0.9, n_folds=5, seed=0)
print(round(best_sigma, 4))
0.0143