# Same or different? The Bayesian way

Say I have the following model:

$$\text{Poisson}(\lambda) \sim \begin{cases} \lambda_1 & \text{if } t \lt \tau \\ \lambda_2 & \text{if } t \geq \tau \end{cases}$$

And I infer the posteriors for $\lambda_1$ and $\lambda_2$ shown below from my data. Is there a Bayesian way of telling (or quantifying) if $\lambda_1$ and $\lambda_2$ are the same or different?

Perhaps measuring the probability that $\lambda_1$ is different from $\lambda_2$? Or perhaps using KL divergences?

For example, how can I measure $p(\lambda_2 \neq \lambda_1)$, or at least, $p(\lambda_2 \gt \lambda_1)$?

In general, once you have the posteriors shown below (assume non-zero PDF values everywhere for both), what is a good way of answering this question?

## Update

It seems that this question can be answered in two ways:

1. If we have samples of the posteriors, we could look at the fraction of the samples where $\lambda_1 \neq \lambda_2$ (or equivalently $\lambda_2 > \lambda_1$). @Cam.Davidson.Pilon included an answer that would address this problem using such samples.

2. Integrating some sort of difference of the posteriors. And that's an important part of my question. What would that integration look like? Presumably the sampling approach would approximate this integral, but I would like to know the formulation of this integral.

Note: The plots above come from this material.

-
You can just calculate the variance of both distributions and add them. That's the variance of the difference in the means. Then calculate the difference in the means and see how many standard deviations it is. You can approximate both distributions with normal to start out and use the usual confidence intervals for a normal distribution. They are clearly different means. – Dave31415 Mar 26 '14 at 19:30
Intrinsic hypothesis testing is an answer – Stéphane Laurent Mar 26 '14 at 19:35
All required calculations are provided in my paper but I have not studied the case of $H_0:\{\phi=1\}$ ($\phi$ is the ratio of the two Poisson rates) – Stéphane Laurent Mar 26 '14 at 19:38
Thanks @StéphaneLaurent. Your paper is a great pointer, but it seems to be specific to Poisson processes. What is the comparison, at a high level, that a Bayesian can do to estimate if $\lambda_2$ is the same or different from $\lambda_1$? Does the analysis have to be distribution specific? – Amelio Vazquez-Reina Mar 26 '14 at 19:56
Sorry @user023472 I don't have energy these days. See Bernardo's papers cited in my paper. "Intrinsic" means that the method is derived from and only from the model. – Stéphane Laurent Mar 26 '14 at 21:08

I think a better question is, are they significantly different?

To answer this, we need to compute $P(\lambda_2 > \lambda_1)$. Call this quantity $p$. If $p \approx 0.50$, then there's equal chance one is larger than the other. On the other hand, If $p$ is really close to 1, then we can be confident that yes $\lambda_2$ is larger (read: different) than $\lambda_1$.

How do we compute $p$? It's trivial in a Bayesian MCMC framework. We have samples from the posterior, so lets just compute the chace that samples from $\lambda_2$ are larger than $\lambda_1$:

 p = np.mean( lambda_2_samples > lambda_1_samples )
print p


I apologize for not including this in the book, I'll definetly add it as I think it's one of the most useful ideas in Bayesian inference

-
The probability is 1.0 they are different, as they are both continuous random variables. Consider: what is your prior guess that $\lambda_1 = \lambda_2$? Do you really think they are actually equal? (Ignore hypothesis testing: we are living in the real world where variables are never actually equal). See this post by my hero, Gelman. Computationally, you can test this by computing np.mean( lambda_2_samples != lambda_1_samples). – Cam.Davidson.Pilon Mar 26 '14 at 20:17
You could define how 'not equal' is meaningful. For example, if, in your example, any difference less than one isn't practically meaningful, then you can look at $P(|\lambda_1-\lambda_2| > 1)$ and that would give you a meaningful statistic for $P(\lambda_1 \ne \lambda_2)$ – Sam Dickson Mar 26 '14 at 20:18
To add to my previous comment, if the variables $\lambda_1$ and $\lambda_2$ were discrete, then there is a chance that $\lambda_2$ = $\lambda_1$. – Cam.Davidson.Pilon Mar 26 '14 at 20:19
oh god, I'd hate to be in that situation! It involves nasty integrals. For most models, you can't actually derive the posteriors. Even if you could, it might still be better to use a computer, just for the sake of getting samples. In summary, samples > formulas for computations like this. – Cam.Davidson.Pilon Mar 26 '14 at 22:08
You're not measuring "sufficiently larger". Consider a distribution with a peak at zero and another with equal masses at peaks -10, 10. Your statistic — the expected value of the indicator that one sample is larger than the other — gives 0.5, but the distributions are clearly totally different. – Neil G Mar 27 '14 at 0:58

As stated, this question is trivial. Assuming $\lambda_1$ and $\lambda_2$ are continuous random variables, $\Pr(\lambda_1=\lambda_2)=0$.

I suspect you are interested in the probability that $\lambda_1$ and $\lambda_2$ are within some $\epsilon$ of each other. In that case, the area of the the difference in the two posterior densities on the interval $[-\epsilon/2, \epsilon/2]$ is your answer. Larger values of overlap indicate that the two posteriors are more similar.

If you would prefer to work with simulated results (and for most problems, we don't have the luxury of choice), simply take the proportion of the results where $\lambda_2>\lambda_1$ as an approximation.

-
Thanks. How does your answer relate to some of the ideas discussed in the comments of the OP? – Amelio Vazquez-Reina Mar 26 '14 at 19:38
Apologies, but I'm not familiar with either of those methods so I cannot meaningfully comment. @Stéphane_Laurent is pretty smart, though, so I'd recommend looking through the link, at a minimum. – General Abrial Mar 26 '14 at 19:41
@user023472 Sorry I don't have the energy today to make an answer about intrinsic discrepancy approach. It is based on the Kullback-Leibler divergence. – Stéphane Laurent Mar 26 '14 at 19:43
@user777 This requires fixing $\epsilon$. What if I just want to see the probability $p(\lambda_2 \gt \lambda_1)$ or $p(\lambda_2 \neq \lambda_1)$? – Amelio Vazquez-Reina Mar 26 '14 at 19:53
Thanks @user777. I am interested in the case when we don't have access to the samples. You had an integral in your post earlier, but you seem to have deleted it. What would that integral look like? – Amelio Vazquez-Reina Mar 26 '14 at 21:56