Properties of the Minimum of Two Poisson Random Variables

Question

I stumbled upon the following problem in my research. We are trying to analyze $Z=\min(X,Y)$ where $X \sim Pois(p\lambda)$ and $Y\sim Pois((1-p)\lambda)$. Note that the RVs expectation is related yet not identical but are independent.

What we are most interested in is a closed form expression for $\mathbb{E}Z$. Or, alternatively, an expression simple enough to prove with that the expectation $\mathbb{E}Z$ is attained at $p=\frac{1}{2}$

I managed to find very little literature on the subject. I saw that in some places this scenario is called a "Poisson Race", but couldn't find anything that is relevant to me.

I tried to go the manual way: \begin{equation} \begin{split} \mathbb{E} Z & = \sum_{n\geq 1} \Pr(min(X,Y) \geq n) \\ & = \sum_{n\geq 1} \Pr(X\geq n\ \text{and}\ Y\geq n) \\ & = \sum_{n\geq 1} \Pr(X\geq n)\cdot \Pr(Y\geq n) \\ & = \sum_{n\geq 1}\Bigg[\Bigg(\sum_{i\geq n} \frac{(p \lambda)^i e^{-p\lambda}}{i!} \Bigg)\Bigg(\sum_{i\geq n} \frac{((1-p) \lambda)^i e^{-(1-p)\lambda}}{i!} \Bigg)\Bigg] \\ & = e^{-\lambda}\sum_{n\geq 1}\Bigg[\Bigg(\sum_{i\geq n} \frac{(p \lambda)^i}{i!} \Bigg)\Bigg(\sum_{i\geq n} \frac{((1-p) \lambda)^i }{i!} \Bigg)\Bigg] \\ & = e^{-\lambda}\sum_{n\geq 1}\Bigg[\Bigg(e^x-e_{n-1}(p\lambda) \Bigg)\Bigg(e^x - e_{n-1}((1-p)\lambda) \Bigg)\Bigg] \\ \end{split} \end{equation}

But this didn't lead to any relatively simple terms. Tried looking into Gamma Taylor partial sums of $e^x$ and Gamma functions $\Gamma (x)$ but again, with no result.

What is obvious, due to the symmetry of the function is that the max is attained at $p=\frac{1}{2}$. Does one see any way to prove so without having to derive once and twice and do all the dirty work?

---

$e_n(x)$ is the Exponential Sum Function

Answer 1

I couldn't find (by math or by Google) any closed form for the expectation, but here is some partial progress.

If $X, Y$ are Poisson with rates $\lambda p$ and $\lambda(1-p)$, then $X+Y$ is Poisson with rate $\lambda$, and when we condition on the value of $X+Y$, we have: $$\Pr[X = k \mid X+Y = n] = \binom nk p^k (1-p)^{n-k}$$ (In other words, $X$ and $Y$ are binomial when $X+Y$ is fixed to $n$.)

It will probably be easier to show that $\mathbb E[Z \mid X+Y=n]$ is maximized when $p = \frac12$, and conclude that $\mathbb E[Z]$ is also maximized when $p = \frac12$, than to deal with $\mathbb E[Z]$ directly. But even this isn't as easy as I thought when I wrote this answer...

I don't think the binomial observation will help you get exact values, but it can help with asymptotics, because for large $n$, you have the Chernoff bound to estimate how often the variable with the larger rate "loses the race" and becomes the minimum.

Answer 2

$\newcommand{\bbx}[1]{\,\bbox[8px,border:1px groove navy]{\displaystyle{#1}}\,} \newcommand{\braces}[1]{\left\lbrace\,{#1}\,\right\rbrace} \newcommand{\bracks}[1]{\left\lbrack\,{#1}\,\right\rbrack} \newcommand{\dd}{\mathrm{d}} \newcommand{\ds}[1]{\displaystyle{#1}} \newcommand{\expo}[1]{\,\mathrm{e}^{#1}\,} \newcommand{\ic}{\mathrm{i}} \newcommand{\mc}[1]{\mathcal{#1}} \newcommand{\mrm}[1]{\mathrm{#1}} \newcommand{\pars}[1]{\left(\,{#1}\,\right)} \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\root}[2][]{\,\sqrt[#1]{\,{#2}\,}\,} \newcommand{\totald}[3][]{\frac{\mathrm{d}^{#1} #2}{\mathrm{d} #3^{#1}}} \newcommand{\verts}[1]{\left\vert\,{#1}\,\right\vert}$ \begin{align} \mathbb{E}\bracks{Z} & = \mathbb{E}\bracks{\min\braces{X,Y}} = \mathbb{E}\bracks{X + Y - \verts{X - Y} \over 2} = {1 \over 2}\,\mathbb{E}\bracks{X} + {1 \over 2}\,\mathbb{E}\bracks{Y} - {1 \over 2}\,\mathbb{E}\bracks{\verts{X - Y}} \\[5mm] & = {1 \over 2}\,p\lambda + {1 \over 2}\,\pars{1 - p}\lambda - {1 \over 2}\,\mathbb{E}\bracks{\verts{X - Y}} = {1 \over 2}\lambda - {1 \over 2}\,\color{#66f}{\mathbb{E}\bracks{\verts{X - Y}}} \end{align}

---

With $\ds{x \equiv p\lambda}$ and $\ds{y \equiv \pars{1 - p}\lambda}$: \begin{align} \color{#66f}{\mathbb{E}\bracks{\verts{X - Y}}} & = \sum_{m = 0}^{\infty}{x^{m}\expo{-p\lambda} \over m!} \sum_{n = 0}^{\infty}{y^{n}\expo{-\pars{1 - p}\lambda} \over n!}\,\verts{m - n} \\[5mm] & = \expo{-\lambda} \sum_{m = 0}^{\infty}\sum_{n = 0}^{m}{x^{m}\,y^{n} \over m!\,n!}\pars{m - n} + \expo{-\lambda} \sum_{m = 0}^{\infty}\sum_{n = m}^{\infty}{x^{m}\,y^{n} \over m!\,n!} \pars{n - m} \\[5mm] & = \expo{-\lambda} \sum_{n = 0}^{\infty}\sum_{m = n}^{\infty}{x^{m}\,y^{n} \over m!\,n!} \pars{m - n} + \expo{-\lambda} \sum_{m = 0}^{\infty}\sum_{n = m}^{\infty}{x^{m}\,y^{n} \over m!\,n!} \pars{n - m} \\[5mm] & = \expo{-\lambda}\sum_{n = 0}^{\infty}\sum_{m = n}^{\infty} {x^{m}\,y^{n} + x^{n}\,y^{m} \over m!\,n!}\pars{m - n} = \expo{-\lambda}\sum_{n = 0}^{\infty}\sum_{m = 0}^{\infty} {x^{m + n}\,y^{n} + x^{n}\,y^{m + n} \over \pars{m + n}!\,n!}m \\[5mm] & = \sum_{m = 0}^{\infty}m\pars{x^{m} + y^{m}} \sum_{n = 0}^{\infty}{\pars{xy}^{n} \over \pars{m + n}!\,n!} \\[5mm] & = \sum_{m = 0}^{\infty}m\bracks{\pars{x \over y}^{m/2} + \pars{y \over x}^{m/2}} \,\mrm{I}_{m}\pars{2\root{xy}} \end{align}

where $\ds{\,\mrm{I}_{\nu}}$ is the Modified Bessel Function of the First Kind.

Our result, '$so\ far$', is given by \begin{align} \mathbb{E}\bracks{Z} & = \mathbb{E}\bracks{\min\braces{X,Y}} \\[5mm] & = {1 \over 2}\lambda - {1 \over 2} \sum_{m = 0}^{\infty}m \bracks{\pars{p \over 1 - p}^{m/2} + \pars{1 - p \over p}^{m/2}} \,\mrm{I}_{m}\pars{2\root{p\bracks{1 - p}}\lambda} \end{align}

Answer 3

Comment. I played with this without getting anything nearly as elegant as @FelixMartin's Answer (+1). I did a quick simulation and found that the relationship between $\mu = E(Z)$ and $p$ depends on the value of $\lambda.$ (In view of @Misha's result, I had initial hopes $\lambda$ might not be crucial, but that seemed counter-intuitive.) For what they may be worth, I post graphs of $\mu/\lambda$ against $p$ for six values of $\lambda.$ (The simulated values should be accurate within the resolution of the plots.)

Addendum. Crude R code is provided below, as requested in Comment. There are two alternative lines beginning z = replicate.... The one with pmin was my initial method. The one with abs was to verify that @FelixMartin's formula gives the same results as mine. Put # at the beginning of the line you want to omit. (Increase 10^3 to 10^4 and 5000 to 10000 for smaller simulation error; slower and not necessary for graphs.) Of course, simulation is for visualization and verification only.

par(mfrow=c(2,3))  # enables six panels per plot
lamb = c(.5, 1, 10, 25, 100, 1000); m=6
for(j in 1:m)      # outer loop for 6 values of lambda
  {
  lam = lamb[j]
  p=seq(.0, 1, by=.05); B = length(p); mu=numeric(B)
  for(i in 1:B)    # inner loop for B values of p
    {   
    pp=p[i]
    z = replicate( 10^3, lam/2 - .5*mean(abs(rpois(5000,pp*lam)-rpois(5000,(1-pp)*lam))) )
    # z = replicate( 10^3, mean(pmin(rpois(5000,pp*lam),rpois(5000,(1-pp)*lam))) )
    # 2nd 'replicate' for z can be substituted for first
    mu[i] = mean(z) }
                   # end inner loop
  plot(p, mu/lam, pch=19, ylim=c(0,.5), main=paste("lambda =",lamb[j]))  }                
                   # end outer loop
par(mfrow=c(1,1))  # returns to default single-panel plot