# How to do bayesian model comparison for control and treatment

I have a model for a biological experiment which has a typical bayesian structure, that is $\lambda \rightarrow X$. Now let's assume for control the parameters are $\lambda_1 \rightarrow X_1$ and $\lambda_2 \rightarrow X_2$. Here $X_1$ is the observed data in one condition, and $X_2$ is the observed data in another condition.

Now I have two hypotheses, $H_0 \sim \lambda_1=\lambda_2$ and $H_1 \sim \lambda_1 \neq \lambda_2$. Now according to the bayesian hypotheses testing the bayes factor will be calculated as follows, \begin{align} BF & = \frac{\int \int P(X_1,X_2|H_1)P(\lambda_1,\lambda_2|H_1)d\lambda_1d\lambda_2}{\int P(X_1,X_2|H_0)P(\lambda_1,\lambda_2|H_0)d\lambda_1d\lambda_2} \\ \end{align} Now I am unable to move further from it, because I don't know how to reflect $\lambda_1$ and $\lambda_2$ are differenty in the integration. I can assume $\lambda' = \lambda_1 -\lambda_2$ then I can do anormal hypothesis testing. But my distribution does not alow me to have a closed form with such a parameter $\lambda'$.

You need to calculate the marginal likelihood for the model with $\lambda_1=\lambda_2=\lambda$ and another one without this assumption, this is:

\begin{align} BF & = \frac{\int \int f(data\vert \lambda_1,\lambda_2)\pi(\lambda_1,\lambda_2)d\lambda_1d\lambda_2}{\int f(data\vert\lambda)\pi(\lambda)d\lambda} ,\\ \end{align}

where $f(data\vert \lambda_1,\lambda_2)$ represents the likelihood for the model without assumptions on $(\lambda_1,\lambda_2)$ and $\pi(\lambda_1,\lambda_2)$ is the prior on these parameters. In the denominator, $f(data\vert\lambda)$ represents the likelihood for the model where you assume $\lambda_1=\lambda_2=\lambda$, and $\pi(\lambda)$ is the prior on $\lambda$. See

https://en.wikipedia.org/wiki/Bayes_factor

If at all possible, I would recommend re-parameterizing the model in terms of $\eta=\lambda_1-\lambda_2$ (alternatively $\eta=\lambda_2-\lambda_1$) and testing $H_0:\eta=0$. This can done using Bayes factors approximated via the Savage Dickey density ratio, which would save you from the integration burden, provided that you can sample from the posterior distribution of $\eta$. See:

https://projecteuclid.org/euclid.ejs/1278682959

• Thanks for making it more clear. I have one doubt, if $\pi(\lambda)$ is further hyper-parametrized, say $\pi(\lambda)=\Gamma(\alpha,\beta)$, and we assume that the $\alpha,\beta$ are estimated on $data$, then don't you think that they yield the same value ? Commented Jul 27, 2016 at 14:36
• @HirakSarkar All priors are hyperparametrized, but I would not recommend estimating the hyper-parameters from the data, unless you meant that these data are from a previous experiment. Data dependent priors are inventions of the devil. Anyway, they should not lead to the same result since $\pi(\lambda_1,\lambda_2)$ is different from $\pi(\lambda)$, and the model is different. tThe model with two $\lambda$s is more flexible. Commented Jul 27, 2016 at 14:46
• I am sorry for my ignorance but what I thought is, I have hyper-priors $\alpha,\beta$ then $\pi(\lambda_1,\lambda_2) = \Gamma(\lambda_1;\alpha,\beta)\Gamma(\lambda_2;\alpha,\beta)$, and $\pi(\lambda) = \Gamma(\alpha,\beta)$. When I integrate out $\lambda$'s they should be different, right ? Here, I am talking about numerical integration, since I don't have a closed form. Thanks again for bearing with me. Commented Jul 27, 2016 at 15:54