GAM factor smooth interaction--include main effect smooth?

Question

I am working with a dataset it which I am interested in modeling an age*sex interaction in the GAMM framework. From examples I have seen and documentation I have read, this is typically accomplished by setting sex as an ordered factor and testing the following model:

1. y ~ ordered_sex + s(age, by = ordered_sex)

When done this way, the smooth term is significant and the effect is strong. However, there is a very strong effect of age itself, and I am not sure if this is driving the factor smooth effect as well, since the variance explained by age is not accounted for anywhere else (y ~ s(age) is wildly significant).

When I run the model this way, accounting for the age effect,

2. y ~ ordered_sex + s(age) + s(age, by = ordered_sex)

the factor smooth term is no longer significant.

Is this second model a valid way of accounting for variance explained by age? Given the difference in the significance of the factor smooth term in these two models, is it correct to assume that the age effect was driving the factor smooth effect in model 1?

Note all models are being run with gamm (there are repeated measures) in R.

Answer 1

You need to be careful with ordered factors here in mgcv as they aren't doing what I think you want to be fitting.

If you pass an ordered factor to by, then gam() etc set up a smooth for all the levels except the reference level, and further more they are set up as smooth differences between the reference level and the level for a specific smooth. What is happening in your first model is that the reference level of age is modelled as a constant term (it is the intercept), with the effect of age for the other levels being smooth differences from this constant.

In the second model, you add s(age), which then models the smooth effect of age in the reference level. Now, the by smooths model smooth differences from this no-longer-constant reference smooth.

I suspect that in the second model, all the levels of sex respond similarly to age hence there are no large deviations from the smooth for the reference level of sex and hence the terms are not significant. In the first model, the effect of age for the reference level was constant, so the difference smooths picked up the actual non-linear effect of age and hence were significantly different from zero.

If you just want to estimate a model with separate smooth function of age for each level of sex I would use an unordered factor (factor(..., ordered = FALSE), not ordered() or factor(..., ordered = TRUE). The the model would be:

y ~ fsex + s(age, by = fsex)

where fsex is the unordered factor.

If you want the model to be explicitly set up like ANOVA contrasts (estimate an effect for the reference level then have differences between individual levels and the reference), then you need to fit the model as per your second example with and ordered factor

y ~ osex + s(age) + s(age, by = osex)

where osex is the ordered factor. But note that in this model, s(age) is not the main smooth effect of age. It is the smooth effect of age in the reference level of osex.