ROPE method for emmGrid and emm_list objects

[rvlenth]

I note that there is a rope method for objects created by the emmeans package. I was curious about what it does, so I looked at the code:

> bayestestR:::rope.emmGrid
function (x, range = "default", ci = 0.95, ci_method = "ETI", 
    verbose = TRUE, ...) 
{
    xdf <- insight::get_parameters(x)
    dat <- rope(xdf, range = range, ci = ci, ci_method = ci_method, 
        verbose = verbose, ...)
    attr(dat, "object_name") <- insight::safe_deparse(substitute(x))
    dat
}

My understanding is that the default for range is c(-0.1, 0.1) * sigma where sigma is the error SD. And you can't recover that from the emmGrid object via insight::get_parameters(). So what is used for sigma when we have range = "default"?

[rvlenth]

OK, I did more experimenting and now believe that the default is just $\pm 0.1$ on an absolute scale. So to do reasonable equivalence testing, we'd have to supply the range manually -- which is actually what I would recommend anyway.

I'd suggest not having a default at all for range, because this is something users should always think about, not just blindly use the default.

[rvlenth]

Here's an illustration that may be useful for showing how ROPE compares with what emmeans provides.

# simulated data for simple one-way layout
set.seed(22.0831)
true.mu = c(21, 24, 21)  # true means -- two are equal
foo = data.frame(trt = rep(factor(rep(1:3, 90))))
foo$y = true.mu[as.numeric(foo$trt)] + rnorm(270, sd = 2.5)

# model
foo.brm = brm(y ~ trt, data = foo)
(foo.prs = pairs(emmeans(foo.brm, "trt")))

##  contrast    estimate lower.HPD upper.HPD
##  trt1 - trt2  -2.7945     -3.51    -2.018
##  trt1 - trt3  -0.0528     -0.80     0.669
##  trt2 - trt3   2.7344      2.02     3.431
## 
## Point estimate displayed: median 
## HPD interval probability: 0.95 

rope(foo.prs, range = c(-.5,.5))

## # Proportion of samples inside the ROPE [-0.50, 0.50]:
## 
## Parameter   | inside ROPE
## -------------------------
## trt1 - trt2 |      0.00 %
## trt1 - trt3 |     85.71 %
## trt2 - trt3 |      0.00 %

For comparison, emmeans supports a delta argument in summary() and test(), which for Bayesian models, adds a column with the fraction having absolute value less than delta

test(foo.prs, delta = 0.5)

##  contrast    estimate lower.HPD upper.HPD p.equiv
##  trt1 - trt2  -2.7945     -3.51    -2.018   0.000
##  trt1 - trt3  -0.0528     -0.80     0.669   0.814
##  trt2 - trt3   2.7344      2.02     3.431   0.000
## 
## Point estimate displayed: median 
## 'p.equiv' based on posterior P(|lin. pred.| < 0.5) 
## HPD interval probability: 0.95 

I think ROPE is expressed relative to the HPD. In this case, note that 0.814 / 0 .950 = 0.857 which is the ROPE value.

In frequentist equivalence testing by the Schuirmann method (aka TOST), we would conclude equivalence at the 0.05 level if the 90% (not 95%) confidence interval lies entirely in $[-\delta, \delta]$. This is because two one-sided tests are done and we require both tail areas to be at least .05. So the rough equivalent is having my p.equiv at least 0.90, and ROPE at least 0.90 / 0.95 = 0.947.

[strengejacke]

@mattansb we should add this to the docs, or maybe to a vignette?

[mattansb]

I'd suggest not having a default at all for range, because this is something users should always think about, not just blindly use the default.

I definitely agree with this in general (that there shouldn't be a default rope range at all, i.e., deprecate rope_range() completely...).

@mattansb we should add this to the docs, or maybe to a vignette?

Perhaps @rvlenth can add the delta example to the docs over at emmeans and reference back to here? Both are fine.

[rvlenth]

Obviously I agree with having to specify a range. There already is an example for the test() part of this with nonzero delta, but I added a link to this issue, which will be available once I push it up.

[mattansb]

Thanks!