7
$\begingroup$

Let $X$ and $Y$ be two Poisson variables with different mean.

Is there a better (as in more concise or numerically faster) way to compute $P(X\leq Y)$ than using

$ P(X\leq Y) = \sum_{y=0}^\infty P(Y=y) P(X\leq y) $

  • 0
    *Independen$t$* Poisson random variables.2011-05-07

2 Answers 2

4

OK, this is how I've done it for now. It's very likely that this is not the best way to do it, but I just decided to post it, for the fun of it. There are probably some mistakes, so go ahead and point them out.

First, for the sake of definiteness, I will posit $X \sim Poiss(\lambda)$ and $Y \sim Poiss(\mu)$. In particular, this means that

$ \mathbb{P}(X=k) = \frac{e^{-\lambda} \lambda^k}{k!} \; .$

Then, notice that

$ \frac{d}{d\lambda} \mathbb{P}(X=k) = \mathbb{P}(X=k-1) - \mathbb{P}(X=k) $

with special case

$ \frac{d}{d\lambda} \mathbb{P}(X=0) = - \mathbb{P}(X=0) \; . $

This implies that

$ \frac{d}{d\lambda} \mathbb{P}(X \leq y) = - \mathbb{P}(X=y) $

by summing our first identity over $k$ from $0$ to $y$.

Now, differentiating the sought after probability with respect to $\lambda$ gives us

$ \frac{d}{d\lambda}\mathbb{P}(X\leq Y) = - \sum_{y=0}^{\infty} \mathbb{P}(X=y)\mathbb{P}(Y=y) $

The sum on the right hand side can be expressed as

$ - \sum_{y=0}^{\infty} \mathbb{P}(X=y)\mathbb{P}(Y=y) = - e^{-(\lambda + \mu)} \sum_{y=0}^{\infty} \frac{\lambda^y \mu^y}{(y!)^2} \; . $

The last part can be rewritten noting that the modified Bessel function $I_0(x)$ has the following Taylor expansion around $0$

$I_0(x) = \sum_{n=0}^{\infty} \frac{(x/2)^{2n}}{(n!)^2} \; .$

We then obtain

$ \frac{d}{d\lambda}\mathbb{P}(X\leq Y) = - e^{-(\lambda + \mu)} I_0(2 \sqrt{\lambda \mu})$

Now, integrating this we get the following expression

$ \mathbb{P}(X\leq Y) - \mathbb{P}(0\leq Y)= - \int_0^{\lambda} e^{-(t + \mu)} I_0(2 \sqrt{t \mu}) dt $

or

$ \mathbb{P}(X\leq Y) = 1 - e^{-\mu}\int_0^{\lambda} e^{-t} I_0(2 \sqrt{\mu t}) dt \; . $

This integral should be in principle more numerically tractable than the previous sum, provided you have modified Bessel functions implemented in whatever program you are using to solve this.

  • 0
    Nice stuff! Thanks much, that will indeed be easier to calculate with the tools I have.2011-02-16
2

Rasholnikov has shown you a method using modified Bessel functions. You could also get that sort of result from a Skellam distribution, which is the difference between two Poisson distributions.

If you want a quick and simple approximation for $X \sim \text{Poiss}(\lambda)$ and $Y \sim \text{Poiss}(\mu)$ then

$\text{Pr}(X \leq Y) \approx \Phi\left( \frac{1/2 - (\lambda - \mu )}{\sqrt{\lambda + \mu}} \right)$

where $\Phi()$ is the cumulative distribution function of a standard normal distribution. The $1/2$ is there as a continuity adjustment to deal with the possibility that $X=Y$.

For example with $\lambda =5$ and $\mu =10$, the true result is about 0.9256 while the approximation gives 0.9222.

  • 0
    Thanks! This approximation is indeed.. simpler =) I'll try if its good enough for my case.2011-02-16