The Difference Of Two Independent Poisson Variables

Why a Difference Has Its Own Distribution

Subtracting two independent counts produces an integer-valued random variable that can be positive, negative, or zero. Its distribution is not generally Poisson: even when both counts have simple Poisson laws, their difference has probability on every integer when both rates are positive. Characteristic functions are especially useful here because independence turns the transform of the difference into a product, and the resulting Fourier series reveals each probability mass.

Work on a probability space \((\Omega,\mathcal{F},\mathbb{P})\). Let \(X\) and \(Y\) be independent random variables with Poisson distributions of rates \(\lambda\) and \(\mu\), where \(\lambda,\mu\ge0\). Thus, for \(r=0,1,2,\ldots\), \(\mathbb{P}(X=r)=e^{-\lambda}\lambda^r/r!\), with the rate-zero distribution interpreted as mass one at zero; the same convention applies to \(Y\). Set \(D=X-Y\), whose state space is the integers \(\mathbb{Z}\). The characteristic function of a real-valued random variable \(Z\) is \(\phi_Z(t)=\mathbb{E}[e^{itZ}]\), for real \(t\).

The probability mass formula comes from taking Fourier coefficients of the characteristic function. This approach both derives the distribution and makes the role of the difference’s negative sign explicit.

Proof. For a Poisson random variable \(X\) of rate \(\lambda\), summing its defining probabilities gives \(\mathbb{E}[e^{itX}]=e^{-\lambda}\sum_{a=0}^{\infty}(\lambda e^{it})^a/a! =\exp(\lambda(e^{it}-1))\). Since \(D=X-Y\), we have \(e^{itD}=e^{itX}e^{-itY}\). Independence therefore gives \(\phi_D(t)=\mathbb{E}[e^{itX}]\mathbb{E}[e^{-itY}] =\exp(\lambda(e^{it}-1)+\mu(e^{-it}-1))\).

Expand both exponentials in their power series: \[ \phi_D(t)=e^{-(\lambda+\mu)} \sum_{a=0}^{\infty}\sum_{b=0}^{\infty} \frac{\lambda^a\mu^b}{a!b!}e^{i(a-b)t}. \] The sum of the absolute values of the coefficients before the constant factor is \(e^\lambda e^\mu\), so this double series is absolutely convergent. We may group its terms by the integer \(a-b\). For each integer \(k\), the coefficient of \(e^{ikt}\) is the sum over \(a-b=k\). If \(k\ge0\), put \(a=r+k\) and \(b=r\), giving \[ e^{-(\lambda+\mu)} \sum_{r=0}^{\infty} \frac{\lambda^{r+k}\mu^r}{(r+k)!r!}. \] If \(k<0\), put \(q=|k|\), \(a=r\), and \(b=r+q\), giving \[ e^{-(\lambda+\mu)} \sum_{r=0}^{\infty} \frac{\lambda^r\mu^{r+q}}{r!(r+q)!}. \] These are the probability masses as well as the Fourier coefficients: for an integer-valued \(D\), \(\phi_D(t)=\sum_{j\in\mathbb{Z}}\mathbb{P}(D=j)e^{ijt}\), and integrating \(\phi_D(t)e^{-ikt}\) over \([-\pi,\pi]\), divided by \(2\pi\), selects the \(j=k\) term because the integral of \(e^{i(j-k)t}\) is zero for \(j\ne k\) and \(2\pi\) for \(j=k\). The exchange of sum and integral is justified by \(\sum_j\mathbb{P}(D=j)=1\).

For \(k\ge0\), factor \((\lambda/\mu)^{k/2}\) from the coefficient sum; its remaining terms are \((\lambda\mu)^{r+k/2}/(r!(r+k)!)\), precisely the series for \(I_k(2\sqrt{\lambda\mu})\). For \(k<0\), with \(q=|k|\), the corresponding factor is \((\mu/\lambda)^{q/2}=(\lambda/\mu)^{k/2}\), and the remaining series is \(I_q(2\sqrt{\lambda\mu})\). This proves the stated mass formula in both cases. Finally, independence gives \(\mathbb{E}[X-Y]=\lambda-\mu\) and \(\operatorname{Var}(X-Y)=\operatorname{Var}(X)+\operatorname{Var}(Y)=\lambda+\mu\), using the Poisson mean and variance. \(\square\)

Reading the Formula and Checking Moments

The factor \(I_{|k|}\) depends on the magnitude of the integer, while \((\lambda/\mu)^{k/2}\) records the direction of the difference. When the rates are equal, this directional factor is one and the mass function is symmetric about zero. When the rates differ, the mean is shifted toward the count with the larger rate. A useful independent check comes directly from differentiating the characteristic function; the sign in the second derivative matters.

Write \(a=\lambda-\mu\) and \(s=\lambda+\mu\). Differentiating the displayed characteristic function gives \(\phi_D'(0)=ia\) and \(\phi_D''(0)=-a^2-s\). For any random variable with finite second moment, differentiating under the expectation yields \(\phi_D'(0)=i\mathbb{E}[D]\) and \(\phi_D''(0)=-\mathbb{E}[D^2]\). Here the moments are finite because \(X\) and \(Y\) are Poisson. Thus \(\mathbb{E}[D]=a\), \(\mathbb{E}[D^2]=a^2+s\), and \(\operatorname{Var}(D)=\mathbb{E}[D^2]-(\mathbb{E}[D])^2=s\). In particular, the second derivative is the negative, not the positive, of the second moment.

Common Errors and a Reliable Check

A frequent mistake is to identify the second derivative of a characteristic function with the second moment without its minus sign. For the rates \(2\) and \(1\), that incorrect calculation would say \(\mathbb{E}[D^2]=\phi_D''(0)=-4\), which is impossible for a square. The correction is \(\mathbb{E}[D^2]=-\phi_D''(0)=4\); subtracting the squared mean gives \(4-1^2=3\).

A second error is to use \(e^{it}\) for both Poisson factors. The difference is \(X-Y\), so the \(Y\) contribution is \(\mathbb{E}[e^{-itY}]=\exp(\mu(e^{-it}-1))\). Replacing it by \(\exp(\mu(e^{it}-1))\) instead derives the characteristic function of the sum \(X+Y\), not the difference. A final check is to inspect the mean: the correct difference has mean \(\lambda-\mu\), while the mistaken sum has mean \(\lambda+\mu\).

The distribution is often called the Skellam distribution. Its mass formula is useful whenever a model compares two independent Poisson counts, such as counts of events in two categories or arrivals from two sources. The characteristic function compactly records the distribution, while coefficient extraction gives individual probabilities and the moment calculation provides a quick consistency check.