One way of fitting Generalised Additive Models (GAM) involves using cubic splines as basis functions. Cubic splines are splines constructed of piecewise third-order polynomials (source). If I understand correctly, given n knots, the data will be subset in n+1 intervals, and the data in each interval will be modelled by a third-order polynomial.
Given a set of data subdivided by n knots, I would expect a cubic spline to be described by at least n + 1 + (n+1)*4 parameters, namely:
- the coordinates of the n knots,
- four parameters for each cubic spline used to model each of the n-1 subsets of data (a, b, c, and d in y = a + bx + cx^2 + dx^3*),
- at least one parameter for the function of the sum of the basis functions that produces the model estimates.
To test this I generated data according to a sinusoid pattern using R (after this example):
set.seed(1) x <- seq(0, pi * 2, 0.1) sin_x <- sin(x) y <- sin_x + rnorm(n = length(x), mean = 0, sd = sd(sin_x / 2)) Sample_data <- data.frame(y,x)
And I modelled them as:
gam_y <- mgcv::gam(y ~ s(x, bs="cr"), knots=list(1,2,3,4,5), method = "REML") # plot data plot(y~x) # plot GAM predictions gam_pred <- predict(gam_y, newdata = data.frame(x = Sample_data[,2])) lines(x, gam_pred, col="red", lwd=2) # display basis functions model_matrix <- predict(gam_y, type = "lpmatrix") matplot(x, model_matrix[,-1], type = "l", lty = 2, add = T)
(GAM fit as a solid line, basis functions as dashed lines)
If I were right, each basis function should be described by 30 parameters (5 + 1 + (5+1)*4). Instead, each basis function is defined by a single parameter:
gam_y$coefficients (Intercept) s(x).1 s(x).2 s(x).3 s(x).4 s(x).5 s(x).6 s(x).7 s(x).8 0.05361059 0.67279976 0.99708672 0.88434236 0.28890966 -0.35159312 -0.82555194 -0.96022992 -0.61303024 s(x).9 0.05663602
I am clearly confused. Can anyone help me understand how each basis function can be defined by a single parameter?