Poisson binomial distribution

The probability of having k successful trials out of a total of n can be written as the sum [1]

${\displaystyle \Pr(K=k)=\sum \limits _{A\in F_{k}}\prod \limits _{i\in A}p_{i}\prod \limits _{j\in A^{c}}(1-p_{j})}$

where ${\displaystyle F_{k}}$ is the set of all subsets of k integers that can be selected from {1,2,3,...,n}. For example, if n = 3, then ${\displaystyle F_{2}=\left\{\{1,2\},\{1,3\},\{2,3\}\right\}}$. ${\displaystyle A^{c}}$ is the complement of ${\displaystyle A}$, i.e. ${\displaystyle A^{c}=\{1,2,3,\dots ,n\}\setminus A}$.

${\displaystyle F_{k}}$ will contain ${\displaystyle n!/((n-k)!k!)}$ elements, the sum over which is infeasible to compute in practice unless the number of trials n is small (e.g. if n = 30, ${\displaystyle F_{15}}$ contains over 10²⁰ elements). However, there are other, more efficient ways to calculate ${\displaystyle \Pr(K=k)}$.

As long as none of the success probabilities are equal to one, one can calculate the probability of k successes using the recursive formula [2] [3]

${\displaystyle \Pr(K=k)={\begin{cases}\prod \limits _{i=1}^{n}(1-p_{i})&k=0\\{\frac {1}{k}}\sum \limits _{i=1}^{k}(-1)^{i-1}\Pr(K=k-i)T(i)&k>0\\\end{cases}}}$

where

${\displaystyle T(i)=\sum \limits _{j=1}^{n}\left({\frac {p_{j}}{1-p_{j}}}\right)^{i}.}$

The recursive formula is not numerically stable, and should be avoided if ${\displaystyle n}$ is greater than approximately 20. An alternative is to use a divide-and-conquer algorithm: if we assume ${\displaystyle n=2^{b}}$ is a power of two, denoting by ${\displaystyle f(p_{i:j})}$ the Poisson binomial of ${\displaystyle p_{i},\dots ,p_{j}}$ and ${\displaystyle *}$ the convolution operator, we have ${\displaystyle f(p_{1:2^{b}})=f(p_{1:2^{b-1}})*f(p_{2^{b-1}+1:2^{b}})}$.

Another possibility is using the discrete Fourier transform.[4]

${\displaystyle \Pr(K=k)={\frac {1}{n+1}}\sum \limits _{l=0}^{n}C^{-lk}\prod \limits _{m=1}^{n}\left(1+(C^{l}-1)p_{m}\right)}$

where ${\displaystyle C=\exp \left({\frac {2i\pi }{n+1}}\right)}$ and ${\displaystyle i={\sqrt {-1}}}$.

Still other methods are described in "Statistical Applications of the Poisson-Binomial and conditional Bernoulli distributions" by Chen and Liu.[5]

Since a Poisson binomial distributed variable is a sum of n independent Bernoulli distributed variables, its mean and variance will simply be sums of the mean and variance of the n Bernoulli distributions:

${\displaystyle \mu =\sum \limits _{i=1}^{n}p_{i}}$

${\displaystyle \sigma ^{2}=\sum \limits _{i=1}^{n}(1-p_{i})p_{i}}$

For fixed values of the mean (${\displaystyle \mu }$) and size (n), the variance is maximal when all success probabilities are equal and we have a binomial distribution. When the mean is fixed, the variance is bounded from above by the variance of the Poisson distribution with the same mean which is attained asymptotically as n tends to infinity.

There is no simple formula for the entropy of a Poisson binomial distribution, but the entropy is bounded above by the entropy of a binomial distribution with the same number parameter and the same mean. Therefore, the entropy is also bounded above by the entropy of a Poisson distribution with the same mean.[6]

The Shepp–Olkin concavity conjecture, due to Lawrence Shepp and Ingram Olkin in 1981, states that the entropy of a Poisson binomial distribution is a concave function of the success probabilities ${\displaystyle p_{1},p_{2},\dots ,p_{n}}$.[7] This conjecture was proved by Erwan Hillion and Oliver Johnson in 2015.[8] The Shepp-Olkin monotonicity conjecture, also from the same 1981 paper, is that the entropy is monotone increasing in ${\displaystyle p_{i}}$, if all ${\displaystyle p_{i}\leq 1/2}$. This conjecture was also proved by Hillion and Johnson, in 2019 [9]

The probability that a Poisson binomial distribution gets large, can be bounded using its moment generating function as follows (valid when ${\displaystyle s\geq \mu }$ and for any ${\displaystyle t>0}$):

${\displaystyle {\begin{aligned}\Pr[S\geq s]&\leq \exp(-st)\operatorname {E} \left[\exp \left[t\sum _{i}X_{i}\right]\right]\\&=\exp(-st)\prod _{i}(1-p_{i}+e^{t}p_{i})\\&=\exp \left(-st+\sum _{i}\log \left(p_{i}(e^{t}-1)+1\right)\right)\\&\leq \exp \left(-st+\sum _{i}\log \left(\exp(p_{i}(e^{t}-1))\right)\right)\\&=\exp \left(-st+\sum _{i}p_{i}(e^{t}-1)\right)\\&=\exp \left(s-\mu -s\log {\frac {s}{\mu }}\right),\end{aligned}}}$

where we took ${\textstyle t=\log \left(s/\mu \right)}$. This is similar to the tail bounds of a binomial distribution.