How can I use "brms" to compute a Bayes Factor representing the probability that H0 is true?

[David Sidhu]

I have a simple model in brms:

prior<-c(set_prior("cauchy(0, 2.5)",class = "b", coef = ""),set_prior("cauchy(0, 10)", class = "Intercept", coef = ""))

m1<-brm(DV~IV+(1|Subject), data = data, family="gaussian", prior = prior, warmup = 1000, iter = 2000, chains = 20

The DV is a continuous variable; the IV is a dichotomous variable.

I would like to compute the Bayes Factor that the IV has no effect on the DV (i.e., that its coefficient is equal to 0).

I thought I might do this like this:

h1 <- hypothesis(data, "IV = 0")

However it returns an NA for the evidence ratio:

Hypothesis Tests for class b:

             Estimate Est.Error l-95% CI u-95% CI Evid.Ratio Star

(IV) =       0     -0.2      0.46    -1.12     0.72         NA     

---

'*': The expected value under the hypothesis lies outside the 95% CI.

Is my best option to compute a model without the variable of interest, and then do:

bayes_factor(fit1, fit2)

Is this the same as asking: what is the Bayes Factor that IV = 0?

I'm curious if this is still "Bayesian" as it seems not that different from regular model comparisons, would love to hear about the difference.

Is there another way?

Thanks! 

[Paul Buerkner]

Yes, this is the same (up to numerical error). Moreover, bayes_factor() is more reliable than hypothesis() in case of diffiult priors. For instance, the cauchy() prior may cause problems for hypothesis().

Please see ?hypothesis for explanation why Evid.Ratio is NA by default.
 

--
You received this message because you are subscribed to the Google Groups "brms-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to brms-users+unsubscribe@googlegroups.com.
To post to this group, send email to brms-...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/brms-users/8c44da10-0bd2-4e40-89dd-2b9bc93eaf23%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[David Sidhu]

Perfect, thanks so much! So it sounds like the best way to test the null hypothesis is to compare a model with and without a variable of interest, for instance:

m1<-brm(DV~IV+(1|Subject), data = data, family="gaussian", prior = prior, warmup = 1000, iter = 2000, chains = 20, save_all_pars = TRUE)

m0<-brm(DV~(1|Subject), data = data, family="gaussian", prior = prior, warmup = 1000, iter = 2000, chains = 20, save_all_pars = TRUE)

test<-bayes_factor(h1, h0)

Also, as a side note is this a good guideline for interpretation? https://www.r-bloggers.com/what-does-a-bayes-factor-feel-like/

To unsubscribe from this group and stop receiving emails from it, send an email to brms-users+...@googlegroups.com.

[Paul Buerkner]

I think the blog post you mention is a good introduction to understanding Bayes factors.

To unsubscribe from this group and stop receiving emails from it, send an email to brms-users+unsubscribe@googlegroups.com.

To post to this group, send email to brms-...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/brms-users/c32b98de-f250-437c-8b04-b1aa1711de13%40googlegroups.com.

[David Sidhu]

Hi Paul

One more quick question. I seem to have found that Bayes Factors traditionally compares a posterior distribution to a prior distribution within the same model. However, we are talking about comparing two models. Do you mind explaining how the latter approach approximates the former?

Thanks!

[Paul Buerkner]

See https://www.sciencedirect.com/science/article/pii/S0010028509000826 for an explanation of the Savage Dickey Ratio.

To unsubscribe from this group and stop receiving emails from it, send an email to brms-users+unsubscribe@googlegroups.com.

To post to this group, send email to brms-...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/brms-users/4a01a15c-ff76-4035-91f0-b389cd7fa388%40googlegroups.com.

[David Sidhu]

Thanks very much!

[David Sidhu]

Am I correct that running this:

m1<-brm(DV~IV+(1|Subject), data = data, family="gaussian", prior = prior, warmup = 1000, iter = 2000, chains = 20, save_all_pars = TRUE)

m0<-brm(DV~(1|Subject), data = data, family="gaussian", prior = prior, warmup = 1000, iter = 2000, chains = 20, save_all_pars = TRUE)

test<-bayes_factor(h1, h0)

Will provide P(D|H1) / P(D|H0)? (In this case H1: b_IV ≠ 0; H0: b_IV = 0).

Is there a way to instead compute P(H0|D)? Wagenmakers et al. (2010) suggest that calculating BFH0/(BFH0 + 1) provides the probability of P(H0|D). However, the results don't seem to line up with my results when BF is calculated with the above method. For instance, take the following posterior distribution:

The above method gives P(H0|D) = .70. Just based on an eye test, this doesn't seem to be true; however, maybe I am making a mistake and that is a reasonable probability?

On Tuesday, 6 March 2018 10:57:12 UTC-7, Paul Buerkner wrote:

[Paul Buerkner]

> Will provide P(D|H1) / P(D|H0)? (In this case H1: b_IV ≠ 0; H0: b_IV = 0).

This is correct. If you are interested in posterior probabilities of models / hypothesis, try out method post_prob()

To unsubscribe from this group and stop receiving emails from it, send an email to brms-users+unsubscribe@googlegroups.com.

To post to this group, send email to brms-...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/brms-users/98a13378-bea7-4427-adf8-8a8098268fda%40googlegroups.com.