`emmeans` integration

[mattansb]

It would be great to be able to estimate contrasts / simple effects / simple slopes and their credible intervals from the posterior using emmeans (instead of doing this manually by summing up the correct combination of parameters)

This would require is adding (and exporting) two new functions: recover_data.BFBayesFactor and emm_basis.BFBayesFactor. Some more info can be found here.

Thanks!

[richarddmorey]

I'm not sure what I'd do for the dffun part of the emm_basis; it seems geared toward classical tests. Do you have any suggestions?

[mattansb]

Looking at the emm_basis.stanreg and emm_basis.brmsfit, dffun should be set to function(k, dfargs) Inf (e.g., here). It seems that when the output from emm_basis contains the argument post.beta, that represents sampling from the posterior distribution, emmeans changes its handling of the object from a Frequentist to a Bayesian framework.

[mattansb]

Since BayesFactor::posterior can be coerced into an mcmc object with BayesFactor:::as.mcmc.BFmcmc, it might be possible to pass the object to emmeans:::emm_basis.mcmc?

I've tried playing around with this, but was unsuccessful - but maybe you'll find this useful / inspiring?

[mattansb]

(sorry for spamming you - I seem to have found myself intensely working at this) This issue can be set as resolved (as can #133 ). I will add the code to my BFEffects package, unless you think it should be part of the main BayesFactor.

[mattansb]

Here is my attempt at this integration.

It's not perfect, but it's pretty close. Feel free to make any (or no) use of this going forward.

[richarddmorey]

Do you have some example code that uses this?

On Sun, May 12, 2019 at 9:09 AM Mattan S. Ben-Shachar < notifications@github.com> wrote:

Here is my attempt at this integration. https://gist.github.com/mattansb/f383af8c0dfe88922fb5c75e8572d03e

It's not perfect, but it's pretty close. Feel free to make any (or no) use of this going forward.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/richarddmorey/BayesFactor/issues/128#issuecomment-491574669, or mute the thread https://github.com/notifications/unsubscribe-auth/AAJZVWVJTH5KIGIRNQOPOOLPU7GB7ANCNFSM4GOL3ARA .

-- Richard D. Morey Reader, School of Psychology Cardiff University http://psych.cf.ac.uk/morey

[mattansb]

library(BayesFactor)
#> Loading required package: coda
#> Loading required package: Matrix
#> ************
#> Welcome to BayesFactor 0.9.12-4.2. If you have questions, please contact Richard Morey (richarddmorey@gmail.com).
#> 
#> Type BFManual() to open the manual.
#> ************
library(emmeans)
library(magrittr)
library(purrr)
#> 
#> Attaching package: 'purrr'
#> The following object is masked from 'package:magrittr':
#> 
#>     set_names

source("https://gist.githubusercontent.com/mattansb/f383af8c0dfe88922fb5c75e8572d03e/raw/9616f4a2c959716a1e82bb4a20fd60a8fe8d0f3e/BFBayesFactor_2_emmeans.R")
options(contrasts=c('contr.sum', 'contr.poly'))

data(puzzles)
fit_freq <- aov(RT ~ shape*color + Error(ID/(shape*color)), data = puzzles)

fit_BF <- anovaBF(RT ~ shape*color + ID, data = puzzles, whichRandom = "ID")

emmip(fit_freq,color~shape, CIs = TRUE)
#> Registered S3 methods overwritten by 'ggplot2':
#>   method         from 
#>   [.quosures     rlang
#>   c.quosures     rlang
#>   print.quosures rlang

emmip(fit_BF[4],color~shape, CIs = TRUE)

emmeans(fit_freq,pairwise~shape|color)
#> $emmeans
#> color = color:
#>  shape  emmean    SE   df lower.CL upper.CL
#>  round      45 0.733 16.9     43.5     46.5
#>  square     44 0.733 16.9     42.5     45.5
#> 
#> color = monochromatic:
#>  shape  emmean    SE   df lower.CL upper.CL
#>  round      46 0.733 16.9     44.5     47.5
#>  square     45 0.733 16.9     43.5     46.5
#> 
#> Warning: EMMs are biased unless design is perfectly balanced 
#> Confidence level used: 0.95 
#> 
#> $contrasts
#> color = color:
#>  contrast       estimate    SE   df t.ratio p.value
#>  round - square        1 0.603 20.5 1.658   0.1125 
#> 
#> color = monochromatic:
#>  contrast       estimate    SE   df t.ratio p.value
#>  round - square        1 0.603 20.5 1.658   0.1125
emmeans(fit_BF[4],pairwise~shape|color)
#> $emmeans
#> color = color:
#>  shape  emmean lower.HPD upper.HPD
#>  round    45.0      43.6      46.5
#>  square   44.1      42.7      45.7
#> 
#> color = monochromatic:
#>  shape  emmean lower.HPD upper.HPD
#>  round    45.9      44.4      47.4
#>  square   45.0      43.6      46.6
#> 
#> Results are given on the : (not the response) scale. 
#> HPD interval probability: 0.95 
#> 
#> $contrasts
#> color = color:
#>  contrast       estimate lower.HPD upper.HPD
#>  round - square    0.866    -0.133      1.80
#> 
#> color = monochromatic:
#>  contrast       estimate lower.HPD upper.HPD
#>  round - square    0.847    -0.160      1.84
#> 
#> HPD interval probability: 0.95

^{Created on 2019-05-12 by the reprex package (v0.2.1)}

As you can see, it basically works, but things would sometimes get tricky with categorical-by-continuous interactions, that I just couldn't figure out.

[JanaJarecki]

Hi, nice to see the link to emmeans. Is the issue with categorical-by-continuous interactions still a thing? Cheers Jana

[mattansb]

No change on this end. It seems that the implementation in JASP has some method of posterior estimates, but I'm not sure.