locus
locus

Reputation: 431

custom contrasts in emtrends

With this this example data set

dat <- data.frame(sample=c(1,2,3,4,5,6, 7,8,9,10,11,12,13,14), treatment=c(1,2,3,4,5,6, 1,2,3,4,5,6,7,8), condition=factor(c("A","A","A","A","A","A", "B","B","B","B","B","B","B","B")), scores=c(2,4,3,5,6,13, 41,30,30,23,24,24,10,8))
mod <- lm(scores~treatment*condition, data=dat)

I set up custom contrasts to look at a linear increase for condition A across treatment, and linear decrease for condition B across treatment variable:

con <- matrix(c(contr.poly(6)[,1],rep(0,8),rep(0,6),contr.poly(8)[,1]*-1), ncol=2)
con
            [,1]        [,2]
 [1,] -0.5976143  0.00000000
 [2,] -0.3585686  0.00000000
 [3,] -0.1195229  0.00000000
 [4,]  0.1195229  0.00000000
 [5,]  0.3585686  0.00000000
 [6,]  0.5976143  0.00000000
 [7,]  0.0000000  0.54006172
 [8,]  0.0000000  0.38575837
 [9,]  0.0000000  0.23145502
[10,]  0.0000000  0.07715167
[11,]  0.0000000 -0.07715167
[12,]  0.0000000 -0.23145502
[13,]  0.0000000 -0.38575837
[14,]  0.0000000 -0.54006172

Is it possible to compute these contrasts with emmeans?

I know I can do emm <- emtrends(mod, ~ condition, var="treatment"), which fits a linear function, but what if I'd like to use my custom contrasts, can they be incorporated into emtrends?

Upvotes: 0

Views: 664

Answers (2)

Russ Lenth
Russ Lenth

Reputation: 6830

I am trying an entirely different answer now that I know more from a few comments.

In the OP, we have a test dataset dat to which we have fitted the model below, and obtained estimated trends

mod <- lm(scores ~ treatment*condition, data = dat)
emt <- emtrends(mod, ~ condition, var = "treatment")   # I renamed it

Here are the estimates, and a direct comparison of those estimates:

> emt
 condition treatment.trend    SE df lower.CL upper.CL
 A                    1.80 0.805 10  0.00604     3.59
 B                   -4.14 0.520 10 -5.30085    -2.98

Confidence level used: 0.95 

> pairs(emt)
 contrast estimate    SE df t.ratio p.value
 A - B        5.94 0.958 10   6.201  0.0001

From a comment:

let's say I wish to get the estimate of a linear increase for condition A and linear decrease for condition B (across treatment)...

I believe that the above comparison of trends does exactly that. The difference of slopes is estimated to be 5.94. By the way, you can get this estimate without using anything in the emmeans package:

> summary(mod) $ coefficients
                      Estimate Std. Error    t value     Pr(>|t|)
(Intercept)          -0.800000  3.1355565 -0.2551381 8.037871e-01
treatment             1.800000  0.8051366  2.2356456 4.936755e-02
conditionB           43.192857  4.0889261 10.5633745 9.597218e-07
treatment:conditionB -5.942857  0.9583042 -6.2014308 1.012643e-04

Note that the estimate and test for treatment:condition is exactly the same.

But continuing the comment:

... and have the coefficients in the same contrast: con <- list(c1=c(contr.poly(6)[,1], contr.poly(8)[,1]*-1))

First of all, this is a contrast of the means, not the slopes. Second, I definitely do not think this is what you want, because you are putting the slope for treatments 1--6 on one scale, and for treatments 7--14 on another.

If you want it all in one contrast and on the same scale, use emmeans (not emtrends) to estimate the means at two different treatment values (say 5 and 6) and two conditions, and construct the appropriate contrast:

> EMM <- emmeans(mod, ~ treatment*condition, at = list(treatment = 5:6))
> contrast(EMM, list(con = c(c(-1, 1), -c(-1, 1))))
 contrast estimate    SE df t.ratio p.value
 con          5.94 0.958 10   6.201  0.0001

This works because the model fits straight lines which have the same slopes regardless of the value of treatment, and the slope is equal to the difference in predictions for two treatments spaced one apart.

Upvotes: 1

Russ Lenth
Russ Lenth

Reputation: 6830

You can do it using

contrast(emm, con)

That's because emtrends(), like emmeans(), returns an object of class emmGrid, which is what contrast() can work with.

By the way, polynomial contrasts are available via a built-in function, so they don't require a custom setup. You can do something like

contrast(emm, "poly", max.degree = 3)

addendum

Oops, I overlooked the fact that con is a matrix. You need to first re-create it as a named list with each element being the desired contrast coefficients. Or a data.frame with the contrast coefficients as columns.

Upvotes: 1

Related Questions