# cross_validate()


K-fold cross-validation for a GAM specification.


Usage

``` python
cross_validate(
    formula,
    data,
    *,
    family=None,
    n_folds=10,
    method="GCV",
    metric="deviance",
    select=False,
    seed=None,
)
```


Estimates the out-of-sample predictive performance of a GAM by repeatedly refitting it on `n_folds - 1` folds of the data and scoring the fit on the remaining held-out fold, then aggregating the resulting losses across folds. Because each fold is scored on data that was not used to fit that particular model, [cross_validate()](cross_validate.md#whittaker.cross_validate) gives a more honest estimate of generalization error than simply scoring a single fit on the data it was trained on, which is optimistic (the fit has already adapted to that data's noise). Use it to compare candidate formulas, families, fitting methods, or basis choices on equal footing, or as a sanity check that a model selected by GCV/REML/ML also performs well out of sample.

Two loss metrics are available via `metric`:

- `"deviance"` (default): for each fold, the family's `deviance()` between the held-out responses and the predictions from a model fit on the training folds, divided by the number of test observations in that fold. This matches the deviance-based loss the family itself uses during fitting, so it is comparable across [method](CausalGAM.md#whittaker.CausalGAM.method) choices for the same `family`.
- `"mse"`: the mean squared error, `mean((y_test - pred) ** 2)`, on the response scale. This is family-agnostic and directly interpretable in the response's original units, but does not account for family-specific variance structure the way deviance does.

Folds are constructed by randomly permuting the row indices with `rng.permutation(n)` (where `rng` is seeded from `seed`) and then assigning fold id `i * n_folds // n` to the row that ends up in position `i` of the permutation. This produces a non-stratified partition into `n_folds` contiguous-in-permutation-order, roughly (but not exactly, when `n` is not a multiple of `n_folds`) equal-size groups; no attempt is made to balance the distribution of the response or any covariate across folds.


## Parameters


`formula: str`  
GAM formula string, e.g. `"y ~ s(x1) + s(x2) + x3"`. The response named on the left-hand side is looked up in `data` to build the fold assignment and to compute the loss; the right-hand side is passed unchanged to `~whittaker.gam.GAM` for every fold.

`data: dict[str, numpy.ndarray] or InputData`  
Column-oriented data as `{name: 1-D array}` (or any `InputData`-compatible object, such as a `pandas.DataFrame` or `polars.DataFrame`). Must contain every column referenced by `formula`, all of equal length.

`family: Family = None`  
Response distribution family, e.g. [Gaussian()](Gaussian.md#whittaker.Gaussian), [Binomial()](Binomial.md#whittaker.Binomial), [Poisson()](Poisson.md#whittaker.Poisson), [Gamma()](Gamma.md#whittaker.Gamma), or [Tweedie()](Tweedie.md#whittaker.Tweedie). Used both to fit each fold's `~whittaker.gam.GAM` and, when `metric="deviance"`, to compute each fold's loss via `family.deviance()`. Defaults to [Gaussian()](Gaussian.md#whittaker.Gaussian).

`n_folds: int = ``10`  
Number of folds to split the data into. Must be at least 2 and, for every fold to receive at least one test observation, should not exceed the number of rows in `data`. Defaults to `10`. See the Notes section below for guidance on choosing this value.

`method: str = ``"GCV"`  
Smoothing-parameter selection method passed through to [GAM.fit()](GAM.md#whittaker.GAM.fit) for every fold. One of `"GCV"` (default), `"REML"`, or `"ML"`; see [GAM.fit()](GAM.md#whittaker.GAM.fit) for what each criterion optimizes.

`metric: str = ``"deviance"`  
Loss metric to compute on each held-out fold: `"deviance"` (default) or `"mse"`. See the discussion above for exactly how each is computed.

`select: bool = ``False`  
Whether to add shrinkage penalties for automatic smooth-term selection, forwarded to `GAM.fit(select=...)` for every fold. Defaults to `False`.

`seed: int = None`  
Seed for the `numpy.random.default_rng()` random number generator used to build the fold assignment. Pass a fixed integer to make the fold split (and hence the resulting [CVResult](CVResult.md#whittaker.CVResult)) reproducible across calls; `None` (the default) uses a fresh, non-reproducible seed.


## Returns


`CVResult`  
Cross-validation result holding the mean out-of-sample loss (`cv_score`), the per-fold losses (`cv_scores`), their standard error (`cv_se`), and the number of folds used (`n_folds`).


## Notes

The number of folds controls a bias-variance tradeoff in the CV estimate itself. With a small `n_folds` (e.g. `3`-`5`), each training fold omits a large fraction of the data, so the fitted model is somewhat different from (typically smoother/less flexible than) a model fit on the full dataset; the resulting `cv_score` tends to be pessimistically biased, but because there are only a few, relatively large folds, `cv_scores` tends to have lower variance across repeated runs. With a large `n_folds` (up to the leave-one-out limit, `n_folds = n`), each training fold is nearly the full dataset, so bias shrinks toward the true generalization error of the full-data fit -- but the individual test folds are tiny (a single point at `n_folds = n`), so `cv_scores` becomes noisier (higher variance), and fitting cost grows linearly with `n_folds` since a full [GAM.fit()](GAM.md#whittaker.GAM.fit) is performed once per fold. In practice, `n_folds = 5` or `n_folds = 10` are common compromises between these effects. Leave-one-out cross-validation is rarely used directly for GAMs because of its cost; `method="GCV"` in [GAM.fit()](GAM.md#whittaker.GAM.fit) already computes an efficient analytical approximation to the leave-one-out error from a single fit, without refitting the model `n` times.


## Examples


``` python
import numpy as np
import whittaker as wt

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

result = wt.cross_validate("y ~ s(x)", {"x": x, "y": y}, n_folds=5, seed=0)
print(result.cv_score, result.cv_se)
```


    0.044761796080802364 0.0034832382870951125


``` python
# Compare two candidate formulas on the same folds using the seed.
linear_result = wt.cross_validate("y ~ x", {"x": x, "y": y}, n_folds=5, seed=0)
smooth_result = wt.cross_validate("y ~ s(x)", {"x": x, "y": y}, n_folds=5, seed=0)
linear_result.cv_score, smooth_result.cv_score
```


    (0.2517163130471863, 0.044761796080802364)