I’m trying to do a post-hoc on my GLMM but I’m not sure I’m doing it right.

I test the effect of continuous variables X and Y on a binary response Z with my IDs as a random effect with this code:

model <- glmer(Z ~ ns(X, 3) * ns(Y, 2) + (1 | id),
           data = df,
           family = binomial(link = "logit"),
           control = glmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 1e+05)),
           nAGQ = 1)

After doing the ‘classic’ analyses, I’d like to know if there’s a difference when Y = 0. Here’s what I did:

emm <- emtrends(model, ~ X * Y, var = "X",
            max.degree = 3,
            at = list(Y = c(0)))
test(emm, adjust = "bonferroni")

That give me:

    degree = linear:
         Y        X        X.trend      SE  df z.ratio p.value
        0    -5.79       -0.18800 0.06186 Inf  -3.039  0.0024

degree = quadratic:
         Y        X        X.trend      SE  df z.ratio p.value
        0    -5.79        0.01569 0.00518 Inf   3.030  0.0024

degree = cubic:
         Y        X        X.trend      SE  df z.ratio p.value
        0    -5.79        0.00295 0.00143 Inf   2.067  0.0388

I deduce from this that there does seem to be a difference in the Y = 0 condition. To go further, I’ve selected some X values that I want to compare two by two in this condition:

emm <- emtrends(model, pairwise ~ X * Y, var = "X",
                max.degree = 3,
                at = list(X = c(-16, -14, -12, -10, -8, -6, -4, -2, 0, 2),
                          Y = c(0)))

test(emm, adjust = "bonferroni")

Is it any good?