4
$\begingroup$

I've fitted 2 different binomial GAMM models, both with a smooth for a predictor Timesec that varies by Condition (a categorical factor). In the first, Condition was a nominal, unordered factor:

bam(cbind(FixTarget, NFixations - FixTarget)
~ 0 + Condition
+ s(Timesec, by=Condition, bs="tp")
+ s(Participant, by=Condition, bs="re")
+ s(Participant, Timesec, by=Condition, bs="re"),
data = by.subject,
family = "binomial",
discrete = TRUE)

In the second, Condition was an ordered factor, and as such, I have also included a 'simple' smooth for Timesec:

bam(cbind(FixTarget, NFixations - FixTarget)
~ 1 + Condition
+ s(Timesec, bs="tp") + s(Timesec, by=Condition, bs="tp")
+ s(Participant, bs="re") + s(Participant, by=Condition, bs="re")
+ s(Participant, Timesec, bs="re") + s(Participant, Timesec, by=Condition, bs="re"),
data = by.subject,
family = "binomial",
discrete = TRUE)

The first model estimates separate smooths for Timesec for the levels of Condition, whereas in the second we estimate a smooth for the reference level of Condition and additional difference smooths, relative to the reference level (see GAM factor smooth interaction--include main effect smooth?).

What I didn't expect was that the predictions from these two models differ (see image): in the model with an ordered factor (bottom left), the predictions for the non-reference level has a much wider CI than in the model with an unordered factor (top left).

Predictions from GAMMs with unordered factor (top row) and ordered factor (bottom row), with random slopes (left column) and without random effects (right column).

I was expecting these two models to be equivalent but with different parametrizations. But it seems the uncertainty about the difference smooth in the ordered factor model is carrying to the non-reference level. Why so and what am I doing wrong?

Also note: In both of these models, there were by-Participant random effects. I've fitted separate models without random slopes, in which I removed s(Participant, Timesec, by=Condition, bs="re"). These models show consistent CIs for ordered vs. unordered models (see right column in the image). This seems to be a relevant piece of the puzzle.

Finally, for completeness, predictions were obtained excluding the random effects, as such (after creating a data frame with new data called predictions.df):

mterms <- sapply(m$smooth, "[[",  "label")
excluded.terms <- mterms[grep("Participant", mterms)]

predictions <- predict(m, newdata = predictions.df,
                       type = "link",
                       se.fit = TRUE,
                       exclude = excluded.terms)

EDIT

Following Gavin's suggestion below, I have run models with factor smooths, instead of with linear adjustments.

The unordered factor model:

bam(cbind(FixTarget, NFixations - FixTarget)
         ~ 0 + Condition
         + s(Timesec, by=Condition, bs="tp")
         + s(Timesec, Participant, by=Condition, bs="fs", m=1),
         data = by.subject,
         family = "binomial",
         discrete = TRUE,
         nthreads = parallel::detectCores() - 1)

The ordered model:

bam(cbind(FixTarget, NFixations - FixTarget)
         ~ 1 + Condition
         + s(Timesec, bs="tp") + s(Timesec, by=Condition, bs="tp")
         + s(Timesec, Participant, bs="fs", m=1)
         + s(Timesec, Participant, by=Condition, bs="fs", m=1),
         data = by.subject,
         family = "binomial",
         discrete = TRUE,
         nthreads = parallel::detectCores() - 1)

The problem remained. In the ordered model, the reference level still yields smooths with smaller CIs, the non-reference level yields larger CIs: Predictions from GAMMs with unordered factor (leftmost) and ordered factor with one and the other factor levels as the reference (middle and rightmost).

$\endgroup$

1 Answer 1

1
$\begingroup$

Isn't this extra uncertainty due to uncertainty in the estimated random slopes and the potential for lack of identification of the linear basis in each smooth and the linear random slopes you introduced to the model.

I suspect the problem is worse in the ordered factor model because now the linear basis in the ordered-factor by smooth is representing deviations from the global (or average) smooth and the random slopes are modelling deviations in another linear function about the linear basis function in the global smooth.

I'm not convinced you should have random slopes in a model where you also include smooth functions of the covariate that is part of the random slopes. Why not allow for random smooths (say via a bs = "fs" smooth)?

$\endgroup$
1
  • $\begingroup$ Many thanks, Gavin. I have now fitted models with factor smooths (see edit), but the problem remains. (I was using these "linear random slopes" because I wanted a more constrained model in which each subject varied little from the average smooth). $\endgroup$
    – JVerissimo
    Mar 3, 2023 at 19:44

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.