Inconsistent EMMEANS results in R and Rj Editor

General help and assistance with jamovi. Bug reports can be made at our issues page: https://github.com/jamovi/jamovi/issues . (If you're unsure feel free to discuss it here)

by kris_ariya » Tue Mar 12, 2019 7:27 pm

Hi,

The emmeans for a mixed ANOVA model in the Rj Editor produces different results from those run in R. Here are the codes.
Code: Select all
id <- factor(c(1:30,1:30))
cond <- c("A","A","A","A","A","A","A","A","A","A",
          "B","B","B","B","B","B","B","B","B","B",
          "C","C","C","C","C","C","C","C","C","C",
          "A","A","A","A","A","A","A","A","A","A",
          "B","B","B","B","B","B","B","B","B","B",
          "C","C","C","C","C","C","C","C","C","C")
time <- c(rep("baseline",30), rep("train",30))
y <- c( 8,5,3,5,2,6,5,6,4,4,
        3,5,8,2,5,6,6,4,3,5,
        3,5,8,5,5,6,6,6,3,4,
        9,7,2,7,9,7,8,5,7,9,
        5,5,10,5,3,10,9,5,7,5,
        4,5,6,6,4,7,7,3,2,2)
a <- data.frame(id,cond,time,y)

# Fit Mixed ANOVA
fit <- aov(y ~ cond * time + Error(id/time), data = a)
summary(fit)

# EMMEANS
library("emmeans")
fit.emm <- emmeans(fit, ~ cond | time)
summary(fit.emm)


Here is the result produced in RStudio
Code: Select all
##Error: id
##          Df Sum Sq Mean Sq F value Pr(>F)
##cond       2  11.43   5.717   1.041  0.367
##Residuals 27 148.30   5.493               
##
##Error: id:time
##          Df Sum Sq Mean Sq F value  Pr(>F)   
##time       1  19.27  19.267   9.441 0.00481 **
##cond:time  2  20.63  10.317   5.055 0.01365 *
##Residuals 27  55.10   2.041                   
##---
##Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
##
##time = baseline:
## cond emmean        SE    df lower.CL upper.CL
## A       4.8 0.6137318 44.63 3.563598 6.036402
## B       4.7 0.6137318 44.63 3.463598 5.936402
## C       5.1 0.6137318 44.63 3.863598 6.336402
##
##time = train:
## cond emmean        SE    df lower.CL upper.CL
## A       7.0 0.6137318 44.63 5.763598 8.236402
## B       6.4 0.6137318 44.63 5.163598 7.636402
## C       4.6 0.6137318 44.63 3.363598 5.836402
##
##Confidence level used: 0.95


The result produced by Rj Editor is attached. The anova models agree but the emmeans estimates were all off.

Best,

Kris
Attachments
Rj.png
Rj.png (155.79 KiB) Viewed 270 times
kris_ariya
 
Posts: 5
Joined: Sat Mar 09, 2019 8:27 pm

by jonathon » Wed Mar 13, 2019 1:19 am

hi,

i'd say there are some options different between your local R, and the Rj R.

`options('contrasts')` seems like a likely candidate.

let me know what you find.

cheers

jonathon
User avatar
jonathon
 
Posts: 815
Joined: Fri Jan 27, 2017 10:04 am

by kris_ariya » Wed Mar 13, 2019 8:35 am

Hi Jonathon,

Both local R and Rj R return the following. There is no difference between their results.
Code: Select all
> options("contrasts")
##$contrasts
##        unordered           ordered
##"contr.treatment"      "contr.poly"


However, I did some digging and found that 'fit.emm@linfct' returns different results for the local R and Rj R.

From local R
Code: Select all
$linfct
##     (Intercept) cond1 cond2 time1 cond1:time1 cond2:time1
##[1,]           1     1     0     1           1           0
##[2,]           1     0     1     1           0           1
##[3,]           1    -1    -1     1          -1          -1
##[4,]           1     1     0    -1          -1           0
##[5,]           1     0     1    -1           0          -1
##[6,]           1    -1    -1    -1           1           1


From Rj R
Code: Select all
##$linfct
##     (Intercept) condB condC timetrain condB:timetrain condC:timetrain
##[1,]           1     0     0         0               0               0
##[2,]           1     1     0         0               0               0
##[3,]           1     0     1         0               0               0
##[4,]           1     0     0         1               0               0
##[5,]           1     1     0         1               1               0
##[6,]           1     0     1         1               0               1


The difference in these matrices seems to result in the difference in estimates. However, I have no idea which options lead to this difference.
kris_ariya
 
Posts: 5
Joined: Sat Mar 09, 2019 8:27 pm

by MAgojam » Wed Mar 13, 2019 3:02 pm

Hi, @kris_ariya.
It seems to me a problem related to the version of the emmeans package used.
In jamovi there is the version "1.3.0", used by Rj, while R local uses the updated version to "1.3.3".

Cheers.
Maurizio
User avatar
MAgojam
 
Posts: 46
Joined: Thu Jun 08, 2017 2:33 pm
Location: Parma (Italy)

by kris_ariya » Wed Mar 13, 2019 7:08 pm

Hi @MAgojam,

Thanks for the clue. The way different versions of 'emmeans' handles option('contrasts') seems to be a culprit here.

I think I found out what happen.

a) Using options('contrasts'), I found that both Rj's and local R's defaults are 'contr.treatment.'

b) The 'emmeans' 1.3.0 in a local R will produce a warning (telling us to use 'contr.sum' instead) but it still goes ahead and produces results with 'contr.treatment' anyway. The results are the same as those produced from Rj. I did not catch this at first because Rj does not show the warning message. Also, changing to 'option(contrasts = c("contr.sum","contr.ply"))' will eventually give me the correct results.

c) In 'emmeans' 1.3.3, the program seems to automatically switch the contrasts to "contr.sum" and, therefore, produces the correct estimates right away.

Thank you all for the help. I hope Rj will get a newer version of 'emmeans' soon.

Best,

Kris
kris_ariya
 
Posts: 5
Joined: Sat Mar 09, 2019 8:27 pm


Return to Help