MLE for Poisson-binomial distribution

Question

I am looking for the maximum likelihood estimator (MLE) for the Poisson-binomial distribution. I understand the derivation of the MLE for a Poisson distribution and a binomial distribution, but I am unable to derive the MLE equation for the Poisson-binomial distribution.

Answer 1

Consider a Poisson-binomial distribution with probability parameter $\boldsymbol{\theta} \equiv (\theta_1, ..., \theta_m) \in [0,1]^m$. The probability density function (PDF) for this distribution is given by:

$$p(k | \boldsymbol{\theta}) = \sum_{\mathcal{A} \in \mathfrak{P}_k(m)} \Bigg( \prod_{j \in \mathcal{A}} \theta_j \Bigg) \Bigg( \prod_{j \notin \mathcal{A}} (1-\theta_j) \Bigg),$$

where $\mathfrak{P}_k(m)$ is the class of all subsets of $k$ labels from $\{ 1, 2, ..., m \}$. The sampling density is clearly invariant to permutations of the probability vector $\boldsymbol{\theta}$, and so this parameter is unidentifiable in the distribution. The minimal sufficient parameter for this distribution is the vector of ordered probability values $\boldsymbol{\theta}_*= (\theta_{(1)} \geqslant ... \geqslant \theta_{(m)})$. This means that the MLE for $\boldsymbol{\theta}$ cannot be estimated uniquely, and for any MLE vector, any permutation of that vector will also be an MLE.

---

Setting aside this identifiability issue, we can obtain equations for the MLE. Letting $\boldsymbol{\theta}_{-a} \in [0,1]^{m-1}$ be the probability parameter excluding the element $\theta_a$, we have the recursive formula:

$$p(k | \boldsymbol{\theta}) = \theta_a p(k-1 | \boldsymbol{\theta}_{-a}) + (1-\theta_a) p(k | \boldsymbol{\theta}_{-a}).$$

We therefore have the useful preliminary result:

$$\frac{\partial p}{\partial \theta_a} (k | \boldsymbol{\theta}) = p(k-1 | \boldsymbol{\theta}_{-a}) - p(k | \boldsymbol{\theta}_{-a}).$$

Now, given the observed vector $\boldsymbol{k} \equiv (k_1, ..., k_n)$ taken from IID draws from the Poisson-binomial distribution, we have log-likelihood function $l_{\boldsymbol{k}} (\boldsymbol{\theta}) = \sum_{i=1}^n \ln p(k_i | \boldsymbol{\theta})$ which gives us the score:

$$\frac{\partial l_\boldsymbol{k}}{\partial \theta_a} (\boldsymbol{\theta}) = \sum_{i=1}^n \frac{\partial}{\partial \theta_a} \ln p(k_i | \boldsymbol{\theta}) = \sum_{i=1}^n \frac{p(k_i-1 | \boldsymbol{\theta}_{-a}) - p(k | \boldsymbol{\theta}_{-a})}{p(k_i | \boldsymbol{\theta})}.$$

The MLE occurs at any point $\hat{\boldsymbol{\theta}}$ satisfying:

$$\sum_{i=1}^n \frac{p(k_i-1 | \hat{\boldsymbol{\theta}}_{-a}) - p(k | \hat{\boldsymbol{\theta}}_{-a})}{p(k_i | \hat{\boldsymbol{\theta}})} = 0 \text{ } \text{ } \text{ } \text{ for all } a=1,2,...,m.$$

We have already noted that the MLE is invariant to permutations, and therefore we only require a "representative" of the class of permutations. We can obtain this by imposing the additional order constraint $\theta_{1} \geqslant ... \geqslant \theta_{m}$. For $n \geqslant m$ this set of $m$ score equations --plus the ordering constraint-- should yield a unique solution for $\hat{\boldsymbol{\theta}}$. This solution represents the class of $m!$ permutations of this vector that are all MLEs.

Answer 2

@Ben The formula you write

p(k|θ)=θap(k−1|θ−a)+(1−θa)p(k|θ−a)

says that the likelihood is a linear function of the probability \theta_a. A linear function has CONSTANT derivative: just differentiate the RHS to see what it is. The derivative is only in exceptional cases equal to zero (depends on the configuration of other parameters). That is to say, setting the derivative to zero makes no sense.

The obvious maximiser of p(k|θ) is the profile of success probabilities (1,...,1,0,...,0) with k 1's. For this profile p(k|θ)=1. To set up the maximisation problem in a more reasonable way one must impose constraints on the profile, for instance by requiring the mean to be some constant. This scenario goes back to the work of Chebyshev in the first half of 19th century, leading to a special case of Hoeffding inequalities.