When you fit a GAM, the smoothing parameters \lambda control the trade-off between fidelity to the data and smoothness of the fitted curves. REML, GCV, and ML each choose a particular set of \lambda values, but the choice is never exact (it is itself an estimate, subject to uncertainty). A natural question is: how much would the predictions change if the smoothing parameters were somewhat different?
The smoothing_sensitivity() method answers this by re-fitting the model across a grid of multiplier values applied to the estimated smoothing parameters. If predictions are stable across a wide range of multipliers, you can be confident that the conclusions do not hinge on the exact \lambda values chosen. If they change substantially, that signals the data do not strongly constrain the smoothness of the fit, and results should be interpreted more cautiously.
Fit a model, then call smoothing_sensitivity(). The method scales all smoothing parameters by each multiplier in the grid, re-fits with those fixed values, and collects the predictions and fit statistics at each step.
import numpy as np
import whittaker as wk
rng = np.random.default_rng(23)
n = 200
x = np.linspace(0, 2 * np.pi, n)
y = np.sin(x) + rng.normal(0, 0.3, n)
model = wk.GAM("y ~ s(x)").fit({"x": x, "y": y})
sens = model.smoothing_sensitivity()
print(sens)SensitivityResult(11 steps, 200 observations)
Multiplier range: [0.01, 100] (baseline=1)
EDF range: [3.0, 9.9]
Max |prediction change|: 0.699The printed summary reports the multiplier range, the span of EDF values across the grid, and the largest absolute prediction change at any observation. A small maximum change relative to the response scale is a sign of robustness.
The SensitivityResult contains arrays indexed by the multiplier step. Each row corresponds to one multiplier value, applied uniformly to all smoothing parameters in the model.
multipliers: the grid of scaling factors. By default, 11 values log-spaced from 0.01 to 100 (i.e., the estimated \lambda divided by 100 up to multiplied by 100). The baseline (the original fit) sits at multiplier 1.0 in the middle of the grid.
predictions: the fitted values at each step, shape (n_steps, n_obs). Use baseline_predictions to get the row corresponding to the original fit.
edf_total: total effective degrees of freedom at each step. Small multipliers (less smoothing) yield higher EDF; large multipliers (more smoothing) shrink the EDF toward 1.
deviance_explained, gcv_scores, aic_values: fit-quality metrics at each step. These help you see whether the chosen \lambda sits near the optimum.
# The baseline index identifies which row matches the original fit
print(f"Baseline multiplier: {sens.multipliers[sens.baseline_idx]:.2f}")
print(f"Baseline EDF: {sens.edf_total[sens.baseline_idx]:.1f}")
print(f"EDF range: [{sens.edf_total.min():.1f}, {sens.edf_total.max():.1f}]")Baseline multiplier: 1.00
Baseline EDF: 6.9
EDF range: [3.0, 9.9]A prediction envelope shows how much the fitted curve varies across the multiplier grid. A narrow envelope means the predictions are robust; a wide one reveals regions where the data do not strongly constrain the fit.
import altair as alt
x_plot = np.linspace(0, 2 * np.pi, 200)
sens_plot = model.smoothing_sensitivity(new_data={"x": x_plot}, n_steps=21)
# Build envelope: min and max prediction at each x across all multipliers
pred_min = sens_plot.predictions.min(axis=0)
pred_max = sens_plot.predictions.max(axis=0)
baseline = sens_plot.baseline_predictions
envelope_data = [
{
"x": float(x_plot[i]),
"lower": float(pred_min[i]),
"upper": float(pred_max[i]),
"baseline": float(baseline[i]),
}
for i in range(len(x_plot))
]
band = alt.Chart({"values": envelope_data}).mark_area(
opacity=0.2, color="steelblue"
).encode(
x=alt.X("x:Q"),
y=alt.Y("lower:Q", title="Prediction"),
y2="upper:Q",
)
line = alt.Chart({"values": envelope_data}).mark_line(
color="steelblue"
).encode(x="x:Q", y="baseline:Q")
obs_data = [{"x": float(x[i]), "y": float(y[i])} for i in range(n)]
points = alt.Chart({"values": obs_data}).mark_circle(
size=15, opacity=0.3, color="gray"
).encode(x="x:Q", y="y:Q")
(band + line + points).properties(
width=500, height=300,
title="Prediction envelope across smoothing parameter multipliers (0.01x–100x)"
)The shaded region covers all predictions from the most wiggly (0.01x) to the most smooth (100x) setting. The solid line is the baseline fit at the estimated \lambda.
Plotting EDF, GCV, or AIC against the multiplier reveals where the optimum sits and how flat or peaked the criterion surface is. A flat minimum means the data tolerate a range of smoothing levels. A sharp minimum means the criterion strongly prefers one setting.
metric_data = [
{"multiplier": float(m), "GCV": float(g), "AIC": float(a), "EDF": float(e)}
for m, g, a, e in zip(
sens.multipliers, sens.gcv_scores, sens.aic_values, sens.edf_total
)
]
gcv_chart = alt.Chart({"values": metric_data}).mark_line(
point=True, color="steelblue"
).encode(
x=alt.X("multiplier:Q", scale=alt.Scale(type="log"), title="Smoothing multiplier"),
y=alt.Y("GCV:Q", title="GCV score"),
).properties(width=350, height=200, title="GCV across multipliers")
aic_chart = alt.Chart({"values": metric_data}).mark_line(
point=True, color="firebrick"
).encode(
x=alt.X("multiplier:Q", scale=alt.Scale(type="log"), title="Smoothing multiplier"),
y=alt.Y("AIC:Q", title="AIC"),
).properties(width=350, height=200, title="AIC across multipliers")
gcv_chart | aic_chartBoth GCV and AIC reach their minimum near multiplier 1.0, confirming that the automatic selection found a good setting. The curves are relatively flat near the minimum, which means moderate changes to \lambda would not substantially affect the fit.
By default, smoothing_sensitivity() uses 11 log-spaced multipliers from 0.01 to 100. You can customize this in several ways:
# Fewer steps, narrower range
sens_narrow = model.smoothing_sensitivity(n_steps=5, log_range=(-1.0, 1.0))
print(f"Multipliers: {sens_narrow.multipliers.round(2)}")Multipliers: [ 0.1 0.32 1. 3.16 10. ]# Explicit multiplier values
sens_custom = model.smoothing_sensitivity(multipliers=[0.1, 0.5, 1.0, 2.0, 10.0])
print(f"Multipliers: {sens_custom.multipliers}")Multipliers: [ 0.1 0.5 1. 2. 10. ]A narrower range is useful when you only care about local sensitivity (e.g., “what if \lambda were half or double its current value?”). A wider range shows the full spectrum from severely underfitting (very large \lambda) to overfitting (very small \lambda).
If the prediction envelope is wide, the fit is sensitive to the choice of \lambda. This can happen when:
k gives the optimizer more freedom, and the criterion surface can become flatter. Reducing k to a value supported by wk.check(model) can sharpen the optimum.In these cases, unconditional=True in predict() already inflates confidence intervals to account for smoothing-parameter uncertainty. The sensitivity analysis complements this by showing the full range of plausible fitted curves, not just the interval at the estimated \lambda.
unconditional=True for smoothing-parameter uncertainty.