How to report the results of anovaBF.

Aram

My apologies if this has already been asked; I couldn't find any answer. I'm new to Bayesian inference, and I've been trying to use BayesFactor. I am reading Bayesian analysis of factorial designs (Rouder, Jeffrey N., et al., 2017) and though it is very useful, I didn't quite understand some parts.

First, on page 30 it says: "Consider the 15th model, which is [...]. The Bayes factor is the comparison of this model against the null model with no age effect, and this age-effect model is less preferable to the null model by 0.69-to-1"

The line on the output looks like this:

Bayes factor analysis

[15] a + s : 1.181082 ±1.16%

And here is an image of the console:

How do I get that .69-to-1?

Second, to report Bayes factors, I read the classification of Harold Jeffreys (1961). He said that BF_{10} > 100 is Extreme evidence for H1. Are the numbers of the second column this BF_{10}? Would this mean that all models from 1 to 14 are strong evidence in favour of the H1?

Could model 1 be interpreted as strong evidence of an effect of age and distance, and is it correct to say that model 18 shows _extreme evidence_of the lack of distance effect?

Third, how can I report this last part? Saying that the best model predicts the data 1.0404717e+18 times better than the model including distance?

Sorry for the very long post. Any reading suggestions on how to report and interpret the results in R using BayesFactor is welcome.

EJ

I will leave this one for Richard to address.

richarddmorey

Hi Aram,

I'm not sure how to reproduce what you did because you did not provide code, but when I load the data and run:

bf = anovaBF(rt ~ a + d + a:d + p + a:p + d:p + a:d:p + s, data = dat, whichRandom = "s")
sort(bf, decreasing = TRUE)

I get

Bayes factor analysis
--------------
[1] a + d + a:d + s                          : 1.335145e+16 ±8.99%
[2] a + d + a:d + p + s                      : 3.495641e+15 ±16.16%
[3] a + d + a:d + p + a:p + s                : 5.291642e+14 ±6.82%
[4] d + s                                    : 3.008386e+14 ±0.52%
[5] a + d + s                                : 1.685939e+14 ±7.62%
[6] a + d + a:d + p + d:p + s                : 1.429879e+14 ±14.98%
[7] d + p + s                                : 6.239307e+13 ±1.17%
[8] a + d + p + s                            : 4.714871e+13 ±26.1%
[9] a + d + a:d + p + a:p + d:p + s          : 1.965467e+13 ±5.88%
[10] a + d + p + a:p + s                     : 8.65209e+12  ±18.49%
[11] a + d + a:d + p + a:p + d:p + a:d:p + s : 5.055364e+12 ±51.9%
[12] d + p + d:p + s                         : 2.338971e+12 ±2.07%
[13] a + d + p + d:p + s                     : 1.173893e+12 ±7.02%
[14] a + d + p + a:p + d:p + s               : 206743920613 ±4.8%
[15] a + s                                   : 0.6343836    ±8.01%
[16] a + p + s                               : 0.311274     ±55.31%
[17] p + s                                   : 0.1833042    ±1.28%
[18] a + p + a:p + s                         : 0.01847103   ±7.28%

Against denominator:
  rt ~ s 
---
Bayes factor type: BFlinearModel, JZS

Notice that model 15 is not exactly .69 to 1 due to the Monte Carlo sampling used to estimate the Bayes factor, but it is close.

As to interpretation, you can think of the numbers in the middle column -- the Bayes factors -- as (roughly) being "improvements" in the denominator model listed at the bottom. If the number is greater than 1 it is an improvement on the denominator model; less than one, and it is not. The rows with very large numbers
are big improvements on the denominator model.

But you probably want to know more; like, are any simpler models improvements on the full (most complex) model? Since model [18] in bf is the full model, you can try:

sort(bf/bf[18], decreasing = TRUE)

which will yield:

Bayes factor analysis
--------------
[1] a + d + a:d + s                          : 2641.046     ±52.67%
[2] a + d + a:d + p + s                      : 691.4715     ±54.35%
[3] a + d + a:d + p + a:p + s                : 104.6738     ±52.34%
[4] d + s                                    : 59.50878     ±51.9%
[5] a + d + s                                : 33.34951     ±52.45%
[6] a + d + a:d + p + d:p + s                : 28.28438     ±54.02%
[7] d + p + s                                : 12.34195     ±51.91%
[8] a + d + p + s                            : 9.32647      ±58.09%
[9] a + d + a:d + p + a:p + d:p + s          : 3.887883     ±52.23%
[10] a + d + p + a:p + s                     : 1.711467     ±55.09%
[11] a + d + a:d + p + a:p + d:p + a:d:p + s : 1            ±0%
[12] d + p + d:p + s                         : 0.462671     ±51.94%
[13] a + d + p + d:p + s                     : 0.2322074    ±52.37%
[14] a + d + p + a:p + d:p + s               : 0.04089595   ±52.12%
[15] a + s                                   : 1.254872e-13 ±52.51%
[16] a + p + s                               : 6.1573e-14   ±75.84%
[17] p + s                                   : 3.625934e-14 ±51.91%
[18] a + p + a:p + s                         : 3.653748e-15 ±52.4%

Against denominator:
  rt ~ a + d + a:d + p + a:p + d:p + a:d:p + s 
---
Bayes factor type: BFlinearModel, JZS

Notice that many of the simpler models appear to be "improvements" on the most complex model, some by a considerable amount.

You could also see how much worse other models are than the "best" model, with:

sort(bf/max(bf), decreasing = TRUE)

if you like.

(I'd also like to point out the high error estimates in the example; before one made any decisions about reporting these data, it would be prudent to try to reduce the error through more iterations; say, 100000 or more)

Aram

Thank you both.

Richard, I did it again and got almost the same number, so maybe I did something wrong the first time:

Now it makes sense. So the normal thing to report is the ratio, the Bayes Factor and the errors of the models one is interested in?

By the way, thanks for BayesFactor! So useful

richarddmorey

I wouldn't bother reporting the error, provided you get it sufficiently low that it seems trustworthy. Try

plot(bf)

and look at the red error bars to see if you have reason for concern.