# Model comparison with WAIC

The **Widely Applicable Information Criterion (WAIC)** (Watanabe, 2010; Gelman, Hwang & Vehtari, 2014) is a Bayesian model-comparison criterion that estimates out-of-sample predictive accuracy from the posterior log-likelihood matrix. Like [PSIS-LOO](loo.md), it returns an estimate of the **expected log predictive density (ELPD)**, but it uses the variance of the log-likelihood across posterior draws instead of importance-sampling corrections. This makes WAIC cheaper to compute and free of the Pareto k diagnostic concerns that can arise with LOO.


# When to use WAIC

Use `model.waic()` when:

- you want a quick Bayesian model-comparison metric without the importance-sampling overhead of LOO
- you are comparing several Bayesian fits (`method="VI"` or `method="MCMC"`) on the same data
- none of your observations are highly influential (if some are, LOO's per-observation diagnostics are more informative)

WAIC and LOO are asymptotically equivalent and typically give very similar results. LOO is generally preferred when Pareto k diagnostics are clean, because it is more robust to influential observations. WAIC is a good first-pass alternative, especially when fitting many candidate models.

For frequentist fits, use `model.aic` and `model.bic` instead.


# Basic usage

Fit two competing models and call `.waic()` on each.


``` python
import numpy as np
import whittaker as wk
from whittaker.families.poisson import Poisson

rng = np.random.default_rng(23)
n = 300
x = np.linspace(0, 2 * np.pi, n)
lam = np.exp(1.5 * np.sin(x))
y = rng.poisson(lam).astype(float)

data = {"x": x, "y": y}

# Model 1: smooth effect of x
m1 = wk.GAM("y ~ s(x)", family=Poisson()).fit(data, method="VI")

# Model 2: linear effect of x only
m2 = wk.GAM("y ~ x", family=Poisson()).fit(data, method="VI")

w1 = m1.waic(seed=0)
w2 = m2.waic(seed=0)
print(w1)
```


    WAICResult
      ELPD_WAIC: -405.40  (SE 17.58)
      p_WAIC:    7.77
      WAIC:      810.80


# Understanding the WAIC result

The [WAICResult](../reference/WAICResult.md#whittaker.WAICResult) object contains:

**`elpd_waic`** is the summed expected log pointwise predictive density. Higher values mean better predictive accuracy.

**`se_elpd_waic`** is the standard error, computed as \sqrt{n \cdot \text{Var}(\text{pointwise})}.

**`p_waic`** is the effective number of parameters, computed as the sum of the per-observation variance of the log-likelihood across posterior draws. It measures how much the posterior predictions vary from observation to observation. When `p_waic` is much larger than the actual parameter count, it suggests the model may be overfit or misspecified.

**[waic](../reference/GAM.md#whittaker.GAM.waic)** is the WAIC on the deviance scale: -2 \cdot \text{ELPD}\_\text{WAIC}. Lower is better. This scale matches the familiar AIC/BIC convention.


``` python
print(f"ELPD_WAIC:     {w1.elpd_waic:.2f} (SE {w1.se_elpd_waic:.2f})")
print(f"p_WAIC:        {w1.p_waic:.2f}")
print(f"WAIC:          {w1.waic:.2f}")
print(f"Observations:  {len(w1.pointwise)}")
```


    ELPD_WAIC:     -405.40 (SE 17.58)
    p_WAIC:        7.77
    WAIC:          810.80
    Observations:  300


# Comparing models

[waic_compare()](../reference/waic_compare.md#whittaker.waic_compare) computes the paired difference in ELPD between two models, using the pointwise values to account for correlation across observations.


``` python
from whittaker import waic_compare

cmp = waic_compare(w1, w2)
print(cmp)
```


    WAICComparison
      ELPD diff: +86.11  (SE 15.57)
      model 1 preferred


A positive `elpd_diff` means the first model is preferred. As with LOO, a difference larger than two standard errors is strong evidence.


``` python
ratio = abs(cmp.elpd_diff) / cmp.se_diff
print(f"|ELPD diff| / SE = {ratio:.1f}")
if ratio > 2:
    preferred = "model 1 (smooth)" if cmp.elpd_diff > 0 else "model 2 (linear)"
    print(f"Strong evidence for {preferred}")
else:
    print("Models are not clearly distinguishable")
```


    |ELPD diff| / SE = 5.5
    Strong evidence for model 1 (smooth)


# WAIC vs LOO

WAIC and LOO estimate the same quantity (ELPD) and are asymptotically equivalent. For well-behaved models they typically agree closely:


``` python
import warnings

with warnings.catch_warnings():
    warnings.simplefilter("ignore")
    loo1 = m1.loo(seed=0)

print(f"ELPD_WAIC: {w1.elpd_waic:.2f}  (SE {w1.se_elpd_waic:.2f})")
print(f"ELPD_LOO:  {loo1.elpd_loo:.2f}  (SE {loo1.se_elpd_loo:.2f})")
print(f"p_WAIC:    {w1.p_waic:.2f}")
print(f"p_LOO:     {loo1.p_loo:.2f}")
```


    ELPD_WAIC: -405.40  (SE 17.58)
    ELPD_LOO:  -405.48  (SE 17.59)
    p_WAIC:    7.77
    p_LOO:     7.44


The key differences:

|  | WAIC | PSIS-LOO |
|----|----|----|
| **Speed** | Fast (just means and variances) | Requires PSIS smoothing per observation |
| **Diagnostics** | No per-observation diagnostics | Pareto k flags unreliable observations |
| **Robustness** | Sensitive to influential observations | PSIS stabilizes extreme weights |
| **Recommendation** | Good first pass | Preferred when Pareto k diagnostics are clean |


# WAIC with MCMC fits

The workflow is the same for MCMC fits. All stored posterior draws are used automatically.

``` python
model = wk.GAM("y ~ s(x)", family=Poisson()).fit(data, method="MCMC")
w = model.waic()
```


# Where to go next

- **[Model averaging with stacking](stacking.md)**: combine multiple models using optimal stacking weights computed from LOO or WAIC results.
- **[Model comparison with LOO](loo.md)**: PSIS-LOO with per-observation Pareto k diagnostics.
- **[Posterior predictive checks](ppc.md)**: assess whether the model generates realistic data.
- **[MCMC sampling](mcmc.md)**: fitting models with the No-U-Turn Sampler.
- **[Variational inference](variational-inference.md)**: a faster Bayesian alternative.
- **[Cross-validation](cross-validation.md)**: K-fold cross-validation for frequentist fits.