3
$\begingroup$

Let's say I have a model with two continuous predictors (nitrogen and temperature) and one categorical variable (variety).

dat <- data.frame(blocks=rep(c(1:15),each=2), variety=rep(c("A","B","C"),each=10),soil=runif(30,0,10), nitro = runif(30, 0, 10), temp= rnorm(30, 10, 3));

mod <- lmer(soil ~ variety*nitro*temp + (1|blocks), data=dat)

Let's pretend the variety:nitro:temp interaction is significant. How can I calculate the trend involving these two continuous predictors? For example maybe I believe that the slope between soil and nitrogen will be steeper at increasingly higher temperatures only in variety A. Or, in other words, how can I test if the slope between soil and nitrogen changes as a function of temperature in variety A?

I thought this would give me that information

> emtrends(mod, pairwise~variety, var=c("nitro","temp"))

but it throws an error

Error in .subset2(x, i, exact = exact) : subscript out of bounds
Improve this question
$\endgroup$

1 Answer 1

Reset to default
4
$\begingroup$

The var argument specifies the variable whose slope you are interested in. Thus it can be only one character string, not a vector; in your case nitro.

The at argument allows you to specify values of other variables for which you wish do to tests and comparisons. at is a named list. In your case, to address the question you ask that would be:

emtrends(mod, pairwise ~ temp, var="nitro", at=list(variety="A", temp=c(20,40))

The above line tests whether the trend of nitro is different between 20° and 40° in variety A. I picked 20° and 40° arbitrarily; a common choice is $mean-SD$ and $mean+SD$, or the medians of the first and third terciles.

Improve this answer
$\endgroup$
5
  • $\begingroup$ Thanks @Ous, your answer was helpful. However can't you have a regression slope for each categorical predictor that reflects the relationship between the continuous predictors and the DV? For example, I can imagine a regression line in which the DV is in y axis, nitro in x axis and temp in z axis. Each variety would then get its own slope in this 3D plane. It would then just be a matter of comparing the slope for variety A with the slope for variety B for example $\endgroup$
    – locus
    Commented Jun 22, 2019 at 23:34
  • $\begingroup$ I think this might be possible to test with the summary() function if the contrasts are done right, but I thought it would be easier with emtrends somehow. $\endgroup$
    – locus
    Commented Jun 22, 2019 at 23:35
  • $\begingroup$ Actually you would get a 2D plane for the interaction, not a line. The slope of it would be of dimension 2, and I have no idea how you might contrasts such slopes. Imagine you had a 3-way interaction with categorical variables. How would you proceed for post hoc tests using emmeans? $\endgroup$
    – Ous
    Commented Jun 23, 2019 at 0:39
  • $\begingroup$ You are right, it would be a plane not a line. But if you look at the output of summary(mod) you get varietyB:nitro:temp ... 1.75 0.09 I read this like: how different is the nitro:temp interaction for varietyB from the nitro:temp interaction for varietyA (the reference). So here you could potentially calculate whether this interaction is significantly different between the two varieties. Or am I misunderstanding this output? $\endgroup$
    – locus
    Commented Jun 23, 2019 at 22:59
  • $\begingroup$ Also, the code you provided returns the same t.ratio and p.value regardless of what contrast I specify. Changing the code from at=list(variety="A", temp=c(20,40)) to at=list(variety="A", temp=c(20,30)) results in the same t.ratio and p.value $\endgroup$
    – locus
    Commented Jun 23, 2019 at 23:01

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.