Clarification on the use of include or terms in fitted_samples()

[Pembs]

I'm unsure as to what the difference is between using include or terms in fitted_samples, as they result in different .fitted values etc (Option 1 vs Option 2 below).

An example model :

     modGI_all <- gam(OTU237 ~ s(days, k=9) + age + sex + mouthwash + s(samplecode, k=10, bs="re"), 
               data = df, method = "REML" )

Generating new data and using fitted samples:

 new_data_global3 <- tibble(
           days= evenly(df$days, n=280), 
           age=mean(df$age),
                  samplecode=rep("SB3", times=280),
                  sex="Female",
                  mouthwash="No",
          ) |> mutate(.row=row_number())

Option 1:

fv_global  <- gratia::fitted_samples(
                      modGI_all,
                      data = new_data_global3,
                      scale=c("response"),
                      method=c("gaussian"),
                      unconditional=T,
                      n = 1000,
                      include="s(days)"
                       )

Option 2 - I would have thought option 1 and option 2 would have resulted in similar .fitted values, but they are quite different.

  fv_global3  <- gratia::fitted_samples(
                         modGI_all,
                        data = new_data_global2,
                        scale=c("response"),
                        method=c("gaussian"),
                        unconditional=T,
                        n = 1000,
                       terms="s(days)")

I can also replicate the results of option 2,if I put exclude=c("(Intercept)", "age", "sexMale", "mouthwashYes", "s(samplecode)"), instead of using terms.

Could somebody please clarify the difference in using include="s(days)" vs terms="s(days)" when specifying the scale=response? Ideally, I want to show the predicted response of s(days) without including the intercept, covariates or subject info (samplecode). However, there isn't explicit documentation as to what include does in the arguments of fitted_samples https://gavinsimpson.github.io/gratia/reference/fitted_samples.html. Thank you

[gavinsimpson]

I don't know why you think include is a valid argument to predict.gam(), because it isn't. There are two arguments that control which arguments are included in the predictions:

1. terms, and
2. exclude

which is why "there isn't explicit documentation as to what include does in the arguments of fitted_samples", although ?fitted_samples does direct you to ?predict.gam (for the ... argument), which would have alerted you to the fact that you had made a mistake.

terms is how you select which terms need to be included in the predictions, while exclude is the complement of terms and indicates which model terms should not be included in the predictions.

If all you want is posterior draws from a specific smooth, then you would be better off with smooth_samples() as it only requires you specify a smooth's label. You say you want to have things on the response scale - in general this may not be sensible for non-Gaussian data, but as you are fitting a Gaussian GAM with identity link function here there is no distinction between the link scale and the response scale - hence smooth_samples() gives you directly what you want.

As you have already deduced that option 2 with either terms or the complement exclude gives you what you want for fitted_samples() you could also use either form of option 2. Option 1 is just wrong.

[Pembs]

Thank you so much for your reply Gavin- I think I came across an example somewhere online that had used include and when I used it as an option, it ran without flagging up an error which confused me- but you're right I should have just used the help information from ?predict.gam. Greatly appreciate all of your help , all the best