Choosing the basis dimension k in mgcv

GAMs

smoothing

model selection

The k in a GAM smooth is a ceiling on wiggliness, not the fit itself. Set it generously, let the penalty work, and check with k.check that the basis is large enough.

Author

Tidy Ecology

Published

2026-08-07

Every s(x) in an mgcv model carries a k, the basis dimension, and it invites a familiar worry: what if I pick it wrong? The worry is misplaced, because k does not set how wiggly the fitted smooth is. The penalty does that. What k sets is a ceiling: the most wiggliness the smooth is allowed to reach. Provided the ceiling is high enough to contain the structure in the data, the penalty pulls the effective flexibility down to whatever the data support, and the exact value of k above that point barely matters.

This is the same idea as mesh convergence in a numerical solver. You refine the grid until the answer stops moving, then stop worrying about the grid. Here you raise k until the fit stops changing, confirm the basis is not the limiting factor, and move on. This post shows the ceiling in action and demonstrates the k.check diagnostic that tells you when the ceiling is too low.

A fit that needs some wiggliness

The truth combines a slow cycle with a faster ripple, so it needs a fair amount of flexibility to represent, roughly ten to twelve effective degrees of freedom. We fit the same smooth at a ladder of k values and record two things: the effective degrees of freedom the penalty settles on, and how close the fit gets to the known truth.

set.seed(4207)
n <- 300
x <- sort(runif(n))
f <- function(x) sin(2 * pi * x) + 0.6 * sin(6 * pi * x)
sigma <- 0.3
y <- f(x) + rnorm(n, 0, sigma)
d <- data.frame(x = x, y = y)
truef <- f(x)
rmse  <- function(fit) sqrt(mean((fit - truef)^2))

ks <- c(4, 6, 8, 12, 16, 24, 32)
res <- lapply(ks, function(k) {
  m  <- gam(y ~ s(x, bs = "tp", k = k), data = d, method = "REML")
  kc <- k.check(m)
  data.frame(k = k, edf = summary(m)$edf, rmse = rmse(fitted(m)),
             ceiling = k - 1, kindex = kc[1, "k-index"], kp = kc[1, "p-value"])
})
res <- do.call(rbind, res)

At k equal to 4 the basis is too small to hold the ripple. The fit uses 3.0 effective degrees of freedom (near its ceiling of three) and lands 0.376 from the truth. By k equal to 8 the error has dropped to 0.075 and stays there: at k of 24 it is 0.073, and at k of 32 it is 0.073. Adding basis functions past the point of sufficiency does not change the fit, because the penalty absorbs them.

The effective degrees of freedom tell the complementary story. They keep climbing, from 7.0 at k of 8 to 16.0 at 24 and 16.7 at 32, but the growth decelerates and stays far below the ceiling of 31. The gap between the effective degrees of freedom and the ceiling widens as k grows: that gap is the penalty working harder to hold the extra flexibility in check.

Scatter of points with three curves: a dotted true curve with a clear ripple, an orange curve at low k that misses the ripple entirely, and a green curve at high k that follows the ripple closely. — Figure 1: At k equal to 4 (orange) the basis cannot represent the ripple and the fit is smoothed past the real structure. At k equal to 24 (green) the fit recovers the truth (dotted). More basis functions do not overfit; the penalty holds them back.

The k.check diagnostic

If the fit is insensitive to k once the basis is large enough, how do you know you have cleared that threshold? You check whether the residuals still carry pattern the smooth could not reach. The k.check function in mgcv does this by ordering the residuals along the covariate and comparing the variance of differences between neighbours to the overall residual variance. If the basis was too small, neighbouring residuals stay on the same side of zero across the stretch the smooth missed, so their differences are unusually small and the ratio, the k-index, falls below one.

hand_kindex <- function(m) {
  e  <- residuals(m, type = "deviance")
  eo <- e[order(d$x)]
  sqrt(mean(diff(eo)^2) / (2 * var(e)))
}
p7_hk4  <- hand_kindex(m4)
p7_hk24 <- hand_kindex(m24)

The hand version captures the principle: the neighbour-difference index is 0.60 at k of 4 and 1.02 at k of 24, below one where structure remains and near one where it does not. The k.check implementation refines this and adds a permutation p-value, which is the number to read in practice. It flags k of 4 firmly, with a k-index of 0.36 and a p-value below 0.001, and clears k of 8, with a k-index of 1.00 and a p-value of 0.50. At k of 24 the k-index is 1.05, comfortably clear. A low k-index with a small p-value is the signal to raise k and refit.

Two stacked panels against k on the x axis: the top panel shows RMSE dropping steeply then flattening near a low value; the bottom panel shows effective degrees of freedom rising and levelling off, with a dashed line for the k minus 1 ceiling rising well above them. — Figure 2: Error against basis dimension falls sharply and then plateaus; effective degrees of freedom climb but stay well under the k minus 1 ceiling (dashed). Beyond the threshold, k costs computation, not accuracy.

Setting k in practice

The practical rule follows from the ceiling picture. Start k generously, larger than you expect the effective degrees of freedom to be, and let restricted maximum likelihood set the smoothness. Run k.check. If a smooth shows a low k-index and a small p-value, its basis is the binding constraint, so raise its k and refit. If the check is clear, the basis is not limiting the fit and there is nothing to tune.

Two limits are worth stating plainly. The check is a guide, not a guarantee: it detects leftover pattern along the covariate, and a basis that misses structure in some other way can still pass. And the true wiggliness is never known; k only needs to exceed the effective degrees of freedom the data can support, which is itself an estimate from one sample. Raising k past the plateau is cheap insurance against the first problem and costs only computation. Under-setting it is the error that biases the fit, as the k equal to 4 curve showed.

References

Wood SN 2003. Journal of the Royal Statistical Society Series B 65(1):95-114 (10.1111/1467-9868.00374).

Wood SN 2011. Journal of the Royal Statistical Society Series B 73(1):3-36 (10.1111/j.1467-9868.2010.00749.x).

Simpson GL 2018. Frontiers in Ecology and Evolution 6:149 (10.3389/fevo.2018.00149).

Zuur AF, Ieno EN, Walker NJ, Saveliev AA, Smith GM 2009. Springer. ISBN 978-0-387-87457-9.

Wood SN 2017. 2nd edn. CRC Press. ISBN 978-1-4987-2833-1.

--- title: "Choosing the basis dimension k in mgcv" description: "The k in a GAM smooth is a ceiling on wiggliness, not the fit itself. Set it generously, let the penalty work, and check with k.check that the basis is large enough." date: "2026-08-07 10:00" categories: [GAMs, smoothing, model selection] image: thumbnail.png image-alt: "Two panels showing effective degrees of freedom and error against basis dimension k, both flattening once k is large enough." --- Every `s(x)` in an `mgcv` model carries a `k`, the basis dimension, and it invites a familiar worry: what if I pick it wrong? The worry is misplaced, because `k` does not set how wiggly the fitted smooth is. The penalty does that. What `k` sets is a ceiling: the most wiggliness the smooth is allowed to reach. Provided the ceiling is high enough to contain the structure in the data, the penalty pulls the effective flexibility down to whatever the data support, and the exact value of `k` above that point barely matters. This is the same idea as mesh convergence in a numerical solver. You refine the grid until the answer stops moving, then stop worrying about the grid. Here you raise `k` until the fit stops changing, confirm the basis is not the limiting factor, and move on. This post shows the ceiling in action and demonstrates the `k.check` diagnostic that tells you when the ceiling is too low. ## A fit that needs some wiggliness The truth combines a slow cycle with a faster ripple, so it needs a fair amount of flexibility to represent, roughly ten to twelve effective degrees of freedom. We fit the same smooth at a ladder of `k` values and record two things: the effective degrees of freedom the penalty settles on, and how close the fit gets to the known truth. ```{r} #| label: setup #| echo: false #| message: false #| warning: false library(mgcv) library(ggplot2) theme_te <- function() { theme_minimal(base_size = 12) + theme( panel.grid.minor = element_blank(), panel.grid.major = element_line(colour = "#e5e4dc", linewidth = 0.3), plot.background = element_rect(fill = "#f5f4ee", colour = NA), panel.background = element_rect(fill = "#f5f4ee", colour = NA), strip.text = element_text(colour = "#16241d", face = "bold"), text = element_text(colour = "#16241d"), axis.text = element_text(colour = "#16241d"), plot.title = element_text(face = "bold") ) } col_truth <- "#16241d"; col_low <- "#c0622d" col_ok <- "#275139"; col_sage <- "#93a87f" ``` ```{r} #| label: k-ladder #| echo: true set.seed(4207) n <- 300 x <- sort(runif(n)) f <- function(x) sin(2 * pi * x) + 0.6 * sin(6 * pi * x) sigma <- 0.3 y <- f(x) + rnorm(n, 0, sigma) d <- data.frame(x = x, y = y) truef <- f(x) rmse <- function(fit) sqrt(mean((fit - truef)^2)) ks <- c(4, 6, 8, 12, 16, 24, 32) res <- lapply(ks, function(k) { m <- gam(y ~ s(x, bs = "tp", k = k), data = d, method = "REML") kc <- k.check(m) data.frame(k = k, edf = summary(m)$edf, rmse = rmse(fitted(m)), ceiling = k - 1, kindex = kc[1, "k-index"], kp = kc[1, "p-value"]) }) res <- do.call(rbind, res) ``` ```{r} #| label: pull-numbers #| echo: false g <- function(kk, col) res[[col]][res$k == kk] p7_edf4 <- g(4,"edf"); p7_rmse4 <- g(4,"rmse"); p7_ki4 <- g(4,"kindex") p7_edf8 <- g(8,"edf"); p7_rmse8 <- g(8,"rmse"); p7_ki8 <- g(8,"kindex"); p7_p8 <- g(8,"kp") p7_edf24 <- g(24,"edf"); p7_rmse24<- g(24,"rmse"); p7_ki24 <- g(24,"kindex"); p7_p24<- g(24,"kp") p7_edf32 <- g(32,"edf"); p7_rmse32<- g(32,"rmse") ``` At `k` equal to 4 the basis is too small to hold the ripple. The fit uses `r sprintf("%.1f", p7_edf4)` effective degrees of freedom (near its ceiling of three) and lands `r sprintf("%.3f", p7_rmse4)` from the truth. By `k` equal to 8 the error has dropped to `r sprintf("%.3f", p7_rmse8)` and stays there: at `k` of 24 it is `r sprintf("%.3f", p7_rmse24)`, and at `k` of 32 it is `r sprintf("%.3f", p7_rmse32)`. Adding basis functions past the point of sufficiency does not change the fit, because the penalty absorbs them. The effective degrees of freedom tell the complementary story. They keep climbing, from `r sprintf("%.1f", p7_edf8)` at `k` of 8 to `r sprintf("%.1f", p7_edf24)` at 24 and `r sprintf("%.1f", p7_edf32)` at 32, but the growth decelerates and stays far below the ceiling of `r round(res$ceiling[res$k==32])`. The gap between the effective degrees of freedom and the ceiling widens as `k` grows: that gap is the penalty working harder to hold the extra flexibility in check. ```{r} #| label: fig-fits #| echo: false #| fig-width: 7 #| fig-height: 4.5 #| fig-cap: "At k equal to 4 (orange) the basis cannot represent the ripple and the fit is smoothed past the real structure. At k equal to 24 (green) the fit recovers the truth (dotted). More basis functions do not overfit; the penalty holds them back." #| fig-alt: "Scatter of points with three curves: a dotted true curve with a clear ripple, an orange curve at low k that misses the ripple entirely, and a green curve at high k that follows the ripple closely." m4 <- gam(y ~ s(x, bs = "tp", k = 4), data = d, method = "REML") m24 <- gam(y ~ s(x, bs = "tp", k = 24), data = d, method = "REML") gr <- data.frame(x = seq(min(x), max(x), length = 400)) gr$k4 <- predict(m4, gr); gr$k24 <- predict(m24, gr); gr$tru <- f(gr$x) ggplot() + geom_point(data = d, aes(x, y), colour = col_sage, alpha = 0.5, size = 1.3) + geom_line(data = gr, aes(x, tru), colour = col_truth, linetype = "dotted", linewidth = 0.8) + geom_line(data = gr, aes(x, k4), colour = col_low, linewidth = 1) + geom_line(data = gr, aes(x, k24), colour = col_ok, linewidth = 1) + labs(x = "x", y = "y", title = "Too small a basis smooths away real structure") + theme_te() ``` ## The k.check diagnostic If the fit is insensitive to `k` once the basis is large enough, how do you know you have cleared that threshold? You check whether the residuals still carry pattern the smooth could not reach. The `k.check` function in `mgcv` does this by ordering the residuals along the covariate and comparing the variance of differences between neighbours to the overall residual variance. If the basis was too small, neighbouring residuals stay on the same side of zero across the stretch the smooth missed, so their differences are unusually small and the ratio, the k-index, falls below one. ```{r} #| label: hand-kindex #| echo: true hand_kindex <- function(m) { e <- residuals(m, type = "deviance") eo <- e[order(d$x)] sqrt(mean(diff(eo)^2) / (2 * var(e))) } p7_hk4 <- hand_kindex(m4) p7_hk24 <- hand_kindex(m24) ``` The hand version captures the principle: the neighbour-difference index is `r sprintf("%.2f", p7_hk4)` at `k` of 4 and `r sprintf("%.2f", p7_hk24)` at `k` of 24, below one where structure remains and near one where it does not. The `k.check` implementation refines this and adds a permutation p-value, which is the number to read in practice. It flags `k` of 4 firmly, with a k-index of `r sprintf("%.2f", p7_ki4)` and a p-value below 0.001, and clears `k` of 8, with a k-index of `r sprintf("%.2f", p7_ki8)` and a p-value of `r sprintf("%.2f", p7_p8)`. At `k` of 24 the k-index is `r sprintf("%.2f", p7_ki24)`, comfortably clear. A low k-index with a small p-value is the signal to raise `k` and refit. ```{r} #| label: fig-plateau #| echo: false #| fig-width: 7 #| fig-height: 4.5 #| fig-cap: "Error against basis dimension falls sharply and then plateaus; effective degrees of freedom climb but stay well under the k minus 1 ceiling (dashed). Beyond the threshold, k costs computation, not accuracy." #| fig-alt: "Two stacked panels against k on the x axis: the top panel shows RMSE dropping steeply then flattening near a low value; the bottom panel shows effective degrees of freedom rising and levelling off, with a dashed line for the k minus 1 ceiling rising well above them." long <- rbind( data.frame(k = res$k, value = res$rmse, panel = "error to truth"), data.frame(k = res$k, value = res$edf, panel = "effective df")) ceil <- data.frame(k = res$k, value = res$ceiling, panel = "effective df") ggplot(long, aes(k, value)) + geom_line(colour = col_ok, linewidth = 1) + geom_point(colour = col_ok, size = 1.8) + geom_line(data = ceil, aes(k, value), colour = col_low, linetype = "dashed", linewidth = 0.7) + facet_wrap(~panel, ncol = 1, scales = "free_y") + labs(x = "basis dimension k", y = NULL, title = "Fit plateaus; effective df stays under the ceiling") + theme_te() ``` ## Setting k in practice The practical rule follows from the ceiling picture. Start `k` generously, larger than you expect the effective degrees of freedom to be, and let restricted maximum likelihood set the smoothness. Run `k.check`. If a smooth shows a low k-index and a small p-value, its basis is the binding constraint, so raise its `k` and refit. If the check is clear, the basis is not limiting the fit and there is nothing to tune. Two limits are worth stating plainly. The check is a guide, not a guarantee: it detects leftover pattern along the covariate, and a basis that misses structure in some other way can still pass. And the true wiggliness is never known; `k` only needs to exceed the effective degrees of freedom the data can support, which is itself an estimate from one sample. Raising `k` past the plateau is cheap insurance against the first problem and costs only computation. Under-setting it is the error that biases the fit, as the `k` equal to 4 curve showed. ## References Wood SN 2003. Journal of the Royal Statistical Society Series B 65(1):95-114 (10.1111/1467-9868.00374). Wood SN 2011. Journal of the Royal Statistical Society Series B 73(1):3-36 (10.1111/j.1467-9868.2010.00749.x). Simpson GL 2018. Frontiers in Ecology and Evolution 6:149 (10.3389/fevo.2018.00149). Zuur AF, Ieno EN, Walker NJ, Saveliev AA, Smith GM 2009. Springer. ISBN 978-0-387-87457-9. Wood SN 2017. 2nd edn. CRC Press. ISBN 978-1-4987-2833-1. ## Related tutorials - [Penalised regression splines from scratch](../penalised-regression-splines/) - [Concurvity in additive models](../concurvity-in-additive-models/) - [Checking a generalised additive model](../checking-a-gam/) - [Model selection with AIC](../model-selection-aic/)