Can you calculate Bayes Factors for a Bayesian Random-effects meta-analysis?

Question

I have the following problem and I wanted to see if somebody with more experience can help me.

I'm doing a Bayesian random-effects meta-analysis in rJAGS with models like these:

model {
for (i in 1:12) {
P[i] <- 1/S.sqr[i]      # Calculate precision
T[i] ~ dnorm(theta[i], P[i]) # study effects
theta[i] ~ dnorm(mu, prec) # random effects

}
mu ~ dnorm(0, 0.5)  # mean difference prior
tau ~ dunif(0,10)   # Uniform on SE
tau.sqr <- tau*tau   # between-study variance
prec<-1/(tau.sqr) # precision of tau
}

Basically, I have about a dozen studies and for each of them I have calculated an effect size in ms, e.g.:

T<- c(7, 12, 6, 23, 15, 9, 17, 20, 15, 11, 25, 9)

I have a theoretically-motivated prior for a null hypothesis (e.g. μ ~ N(0, 0.1)) and an alternative hypothesis (e.g. μ ~ N(15, 0.1)). My question is: is it possible to calculate Bayes factors to quantify the evidence in support of the alternative hypothesis compared to the null hypothesis?

So far, I haven't been able to find a solution to this. I know that the R package "BayesFactor" has a function for a meta-analysis, but it uses t-values from individual studies, which doesn't work for me (it also assumes that the studies estimate the same effect size, which is also not what I want to do).

Answer 1

If you want to know which of the hypotheses is more likely, why don't you calculate the probability that the a posteriori mean $\mu$ is below a certain threshold close to zero, say, $Pr(\mu \le 0.5\, |\, D)$ (where $D$ is your data)? That would allow you to make direct probability statements about the mean of the effect sizes (e.g. "there's an $6\%$ chance that the effect size is, on average, close to zero").

If you are set on using the Bayes factor, a conceptually simple way to calculate it is to set up a model with a random discrete indicator variable whose value decides which prior ($\mu \sim N(0, 0.1)$ or $\mu \sim N(15, 0.1)$) is used in each step. The posterior and prior distributions of this indicator variable can then be used to calculate the Bayes factor, sidestepping the difficult issue of multi-dimensional numerical integration. Unfortunately, I'm not fluent in R, so I can't help you with the implementation.

Cheers!