rules determining nesting in generalized additive models

Question

I understand the concept of nesting when all of the free parameters are clearly identified in linear or generalized linear models (one model is nested within another if its free parameters are a subset of the free parameters of the other). However, it is less clear how to determine nesting in generalized additive models. I would like to find some clear explanation of how to determine if one generalized additive model (gam) is statistically nested within another.

For example, given two models, and using R coding:

m1<-gam(y~x,family="gaussian",...)

m2<-gam(y~s(x), family="gaussian",...)

What are the rules that tell you that m1 is nested within m2, and so allow you to use a likelihood ratio test of the null hypothesis that y in linearly related to x? My intuitive guess is to think that the smoother function of m2 is a bit like a weighted Taylor expansion of the basis functions; something like m2=a+w_1*bx+w_2*cx^2+..., such that the free parameters in m1 (intercept, slope, residual variance) are a subset of the free parameters in this Talyor expansion. However, I am just guessing...

This is only a simple example, but I would like to know the rules, or a reference that gives such rules, in the general case.

Answer 1

Properly nesting these things is going to be difficult in general but potentially not necessary if the reason for nesting is just to compare the models.

Essentially what you are asking is if the effect of x on y is wiggly (non-linear) or will a linear effect suffice.

You can actually test this given the assumptions of the generalised likelihood ratio test for GAMs (see ?anova.gam) and given the basis you used here, a thin plate spline, there should be correct nesting as the basis includes a linear basis function such that you could set all other parameters for basis functions to 0 and the models would be nested (except for the disappearance of the smoothness parameter in the linear version of the model).

I'm not 100% sure how this would work for the other bases in mgcv; there doesn't appear to be any basis function that is a linear function in the cubic regression spline basis (bs = 'cr') in the sense that you can zero all coefficients other than the one for the linear function.

You can also test this idea directly with a single model by modifying the thin plate spline basis when fitting only the second model.

gam(y ~ x + s(x, bs = "tp", m = c(2,0)), family = "gaussian", ...)

The m bit sets up the thin plate regression spline basis (bs = "tp", the default so your models used this anyway) with a second order penalty (again the default penalty is on the integrated squared second derivative of the smooth evaluated over the range of $x$), but no null space (the functions that are perfectly smooth). This is shown below

The above image shows the bases for your m2 (left) and for my model (right). Notice that the pink linear basis function going from bottom left to top right in the left hand figure is missing from the right hand figure, and the basis functions themselves have been modified accordingly.

As such, you can use the no-null-space model parameterisation to test for a non-linear effect using summary().

If you are interested in this more generally, then nesting is not assured in models with say

y ~ s(x) + s(z)

vs

y ~ te(x, z)

again because there aren't terms in the te() basis that map to terms in the individual s() bases. In such cases, the ti() basis can be used to get proper nesting:

y ~ s(x) + s(z) + ti(x, z)

where in the ti() basis the "main" smooth effects of $x$ and $z$ have been removed from the usual tensor product basis.

I should also note that you can use AIC() to compare the models without needing them to be strictly nested.