MathOverflow is a question and answer site for professional mathematicians. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

Let $X$ and $Y$ be Poisson random variables with means $\lambda$ and $1$, respectively. The difference of $X$ and $Y$ is a Skellam random variable, with probability density function $$\mathbb P(X - Y = k) = \mathrm e^{-\lambda - 1} \lambda^{k/2} I_k(2\sqrt{\lambda}) =: S(\lambda, k),$$ where $I_k$ denotes the modified Bessel function of the first kind. Let $F(\lambda)$ denote the probability that $X$ is larger than $Y$: $$F(\lambda) := \mathbb P(X > Y) = \sum_{k=1}^{\infty} S(\lambda, k) = \mathrm e^{-\lambda - 1} \sum_{k=1}^\infty \lambda^{k/2} I_k(2\sqrt{\lambda}).$$ According to Mathematica, the graph of the function $F$ looks like

My question:

  • Is there a closed-form expression for the function $F$?
  • If not, what are $\lim_{\lambda \to 0} F'(\lambda)$ and $F'(1)$? What is the asymptotic behavior as $\lambda \to \infty$?

  • share|cite|improve this question
    Sheldon Ross suggests that for Poisson random variables $X$ and $Y$ with parameters $\lambda$ and $\mu$ respectively, $P(X < Y) = \lambda/(\lambda + \mu)$, which is the same for the probability of two Exponential random variables with the same respective intensity parameters. Not sure how he justifies it, although for a Poisson PROCESS, it makes sense. Any thoughts on why this may be? – user39646 Sep 6 '13 at 21:12

    Regarding the asymptotic behavior when $\lambda \to \infty$:

    To get an estimate one simply finds that the dominating element among $$A_k= \mathbb{P}(Y=k+1)\mathbb{P}(X=k) = e^{-1-\lambda} \frac{\lambda^k}{k!(k+1)!}$$ is $A_{\sqrt{\lambda}}$ which gives roughly $e^{2\sqrt{\lambda}-\lambda}$. Probably there's also a Poly($\lambda$) factor here.

    Oh, and $F'(\lambda)$ is simply $\mathbb{P}(X=Y)$ since the probability of adding 1 when increasing the intensity of a Poisson RV by $\epsilon$ is $\epsilon$. In the special case $\lambda=1$ we get $$F'(1)=e^{-2} \sum_{k=0}^\infty \frac{1}{(k!)^2}$$

    share|cite|improve this answer

    Numerically, $\lim_{\lambda \to 0} F^\prime(\lambda) = 1/e$. Heuristically this should be true because to have $X > Y$ when $\lambda$ is very small, the most likely case will be $X = 1, Y = 0$ by far; that occurs with probability $\lambda e^{-\lambda} e^{-1}$.

    Maple gives an explicit formula for $F^\prime(1)$ involving sums; evaluating gives $F^\prime(1) = 0.3085083225$, and the Inverse Symbolic Calculator says this is somehow to $I_0(2) e^{-2}$. I'm not sure how to prove this at all but maybe knowing the answer helps?

    share|cite|improve this answer
    That's a nice application of the law of small numbers, Michael. Thanks! – Tom LaGatta Apr 9 '11 at 5:04
    Wait, what's the "law of small numbers"? In other words, my cheap trick has a name? – Michael Lugo Apr 9 '11 at 5:28

    An alternative way to see the $1/e$ is as follows.

    Let $x=2\sqrt{\lambda}$. Recall that for small argument $0 < x \ll \sqrt{k+1}$ we have

    $$I_k(x) \approx \frac{1}{\Gamma(k+1)}(x/2)^k$$

    Using this, we see that $$\sum_{k \ge 0} \lambda^{k/2}I_k(2\sqrt{\lambda}) \approx \sum_{k\ge 0}\frac{(x/2)^{2k}}{\Gamma(k+1)}.$$

    This sum is nothing but $e^{x^2/4} = e^\lambda$. Multiplying with $e^{-\lambda-1}$ we obtain the said $1/e$ approximation.

    Perhaps better approximations to $F$ can be obtained in a similar vein.

    share|cite|improve this answer

    By simple computations : The definition of the modified Bessel function of the first kind yields $$ I_k(\lambda)=\sum_{n\geq 0}\frac{1}{n!(n+k)!}\left(\frac{\lambda}{2}\right)^{2n+k} $$ so that we get (the sums transpositions are clearly allowed) $$F(\lambda)=e^{-\lambda-1}\sum_{k\geq 1}\sum_{n\geq 0}\frac{\lambda^{k+n}}{n!(k+n)!}=e^{-\lambda-1}\sum_{n\geq 1}a_n\lambda^n \qquad \mbox{where}\qquad a_n=\frac{1}{n!}\sum_{k=0}^{n-1}\frac{1}{k!}.$$ Thus, deriving under the sign sum

    $$ F'(\lambda) = e^{-\lambda-1}\Big(1+\sum_{n\geq 1 }[(n+1)a_{n+1}-a_n)]\lambda^n\Big) = e^{-\lambda-1}\sum_{n\geq 0}\frac{1}{(n!)^2}\lambda^n $$ we obtain the closed form $$ F'(\lambda)=e^{-\lambda-1}I_0(2\sqrt{\lambda}). $$ One finally get $$F'(0)=e^{-1}, \quad F'(1)=e^{-2}I_0(2)=e^{-2}\sum_{n\geq0}\frac{1}{(n!)^2}$$ and, using the asymptotic formula when $\lambda\rightarrow+\infty$ for all $k$ $$ I_k(\lambda)=\frac{e^{\lambda}}{\sqrt{2\pi\lambda}}\Big(1+O(\lambda^{-1})\Big), $$ that $$ F'(\lambda)=\frac{e^{2\sqrt{\lambda}-\lambda-1}}{2\sqrt{\pi\sqrt{\lambda}}}\Big(1+O(\lambda^{-1/2})\Big) $$ when $\lambda\rightarrow+\infty$.

    share|cite|improve this answer

    Your Answer


    By posting your answer, you agree to the privacy policy and terms of service.

    Not the answer you're looking for? Browse other questions tagged or ask your own question.