Learning Poisson Binomial Distributions

Introduction

We begin by considering a somewhat fanciful scenario: You are the manager of an independent weekly newspaper in a city of $n$ people. Each week the $i$-th inhabitant of the city independently picks up a copy of your paper with probability $p_{i}$. Of course you do not know the values $p_{1},\dots,p_{n}$; each week you only see the total number of papers that have been picked up. For many reasons (advertising, production, revenue analysis, etc.) you would like to have a detailed “snapshot” of the probability distribution (pdf) describing how many readers you have each week. Is there an efficient algorithm to construct a high-accuracy approximation of the pdf from a number of observations that is independent of the population $n$? We show that the answer is “yes.”

A Poisson Binomial Distribution of order $n$ is the distribution of a sum

where $X_{1},\dots,X_{n}$ are independent Bernoulli (0/1) random variables. The expectations $({\bf E}[X_{i}]=p_{i})_{i}$ need not all be the same, and thus these distributions generalize the Binomial distribution ${\rm Bin}(n,p)$ and, indeed, comprise a much richer class of distributions. (See Section 1.2 below.) It is believed that Poisson [Poi37] was the first to consider this extension of the Binomial distributionWe thank Yuval Peres and Sam Watson for this information [PW11]. and the distribution is sometimes referred to as “Poisson’s Binomial Distribution” in his honor; we shall simply call these distributions PBDs.

PBDs are one of the most basic classes of discrete distributions; indeed, they are arguably the simplest $n$-parameter probability distribution that has some nontrivial structure. As such they have been intensely studied in probability and statistics (see Section 1.2) and arise in many settings; for example, we note here that tail bounds on PBDs form an important special case of Chernoff/Hoeffding bounds [Che52, Hoe63, DP09]. In application domains, PBDs have many uses in research areas such as survey sampling, case-control studies, and survival analysis, see e.g., [CL97] for a survey of the many uses of these distributions in applications. Given the simplicity and ubiquity of these distributions, it is quite surprising that the problem of density estimation for PBDs (i.e., learning an unknown PBD from independent samples) is not well understood in the statistics or learning theory literature. This is the problem we consider, and essentially settle, in this paper.

Let $X=\sum_{i=1}^{n}X_{i}$ be an unknown PBD.

[Learning PBDs from constantly many samples] There is an algorithm with the following properties: given $n,\epsilon,\delta$ and access to independent draws from $X$, the algorithm uses

[Properly learning PBDs from constantly many samples] There is an algorithm with the following properties: given $n,\epsilon,\delta$ and access to independent draws from $X$, the algorithm uses

We note that, since every sample drawn from $X$ is a $\log(n)$-bit string, for constant $\delta$ the number of bit-operations performed by our first algorithm is quasilinear in the length of its input. Moreover, the sample complexity of both algorithms is close to optimal, since $\Omega(1/\epsilon^{2})$ samples are required even to distinguish the (simpler) Binomial distributions ${\rm Bin}(n,1/2)$ and ${\rm Bin}(n,1/2+\epsilon/\sqrt{n})$, which have total variation distance $\Omega(\epsilon).$ Indeed, in view of this observation, our second algorithm is essentially sample-optimal.

Let $X=\sum_{i=1}^{n}a_{i}X_{i}$ be a weighted sum of unknown independent Bernoullis such that there are at most $k$ different values among $a_{1},\dots,a_{n}.$ Then there is an algorithm with the following properties: given $n,\epsilon,\delta,$ $a_{1},\dots,a_{n}$ and access to independent draws from $X$, it uses

To complement Theorem 2, we also show that if there are many distinct weights in the sum, then even for weights with a very simple structure any learning algorithm must use many samples:

2 Related work.

At a high level, there has been a recent surge of interest in the theoretical computer science community on fundamental algorithmic problems involving basic types of probability distributions, see e.g., [KMV10, MV10, BS10, VV11] and other recent papers; our work may be considered as an extension of this theme. More specifically, there is a broad literature in probability theory studying various properties of PBDs; see [Wan93] for an accessible introduction to some of this work. In particular, many results study approximations to the Poisson Binomial distribution via simpler distributions. In a well-known result, Le Cam [Cam60] shows that for any PBD $X=\sum_{i=1}^{n}X_{i}$ with ${\bf E}[X_{i}]=p_{i}$, it holds that

Our approach, which removes $n$ completely from the sample complexity, requires a refined understanding of the structure of the set of all PBDs on $n$ variables, in fact one that is more refined than the understanding provided by the aforementioned results (approximating a PBD by a Poisson, Normal, or Binomial distribution). We give an outline of the approach in the next section.

3 Our approach.

The general learning result that we use (Lemma 10) is the following: We show that for any class ${\cal S}$ of target distributions, if ${\cal S}$ has an $\epsilon$-cover of size $N$ then there is a generic algorithm for learning an unknown distribution from ${\cal S}$ to accuracy ${O(\epsilon)}$ that uses $O((\log N)/\epsilon^{2})$ samples. Our approach is rather similar to the algorithm of [DL01] for choosing a density estimate (but different in some details); it works by carrying out a tournament that matches every pair of distributions in the cover against each other. Our analysis shows that with high probability some $\epsilon$-accurate distribution in the cover will survive the tournament undefeated, and that any undefeated distribution will with high probability be $O(\epsilon)$-accurate.

Applying this general result to the $O(\epsilon)$-cover of size $(1/\epsilon)^{O(\log^{2}(1/\epsilon))}$ described above, we obtain a PBD that is $O(\epsilon)$-close to the target (this accounts for the increased running time in Part (2) versus Part (1)). We stress that for both the non-proper and proper learning algorithms sketched above, many technical subtleties and challenges arise in implementing the high-level plan given above, requiring a careful and detailed analysis.

We prove Theorem 2 using the general approach of Lemma 10 specialized to weighted sums of independent Bernoullis with constantly many distinct weights. We show how the tournament can be implemented efficiently for the class ${\cal S}$ of weighted sums of independent Bernoullis with constantly many distinct weights, and thus obtain Theorem 2. Finally, the lower bound of Theorem 3 is proved by a direct information-theoretic argument.

4 Preliminaries.

For a distribution $X$ supported on $[n]=\{0,1,\dots,n\}$ we write $X(i)$ to denote the value $\Pr[X=i]$ of the probability density function (pdf) at point $i$, and $X(\leq i)$ to denote the value $\Pr[X\leq i]$ of the cumulative density function (cdf) at point $i$. For $S\subseteq[n]$, we write $X(S)$ to denote $\sum_{i\in S}X(i)$ and $X_{S}$ to denote the conditional distribution of $X$ restricted to $S.$ Sometimes we write $X(I)$ and $X_{I}$ for a subset $I\subseteq[0,n]$, meaning $X(I\cap[n])$ and $X_{I\cap[n]}$ respectively.

Total Variation Distance.

Recall that the total variation distance between two distributions $X$ and $Y$ over a finite domain $D$ is

Covers.

Poisson Binomial Distribution.

A Poisson binomial distribution $D\in{\cal S}_{n}$ can be represented uniquely as a vector $(p_{i})_{i=1}^{n}$ satisfying $0\leq p_{1}\leq p_{2}\leq\ldots\leq p_{n}\leq 1$. To go from $D\in{\cal S}_{n}$ to its corresponding vector, we find a collection $X_{1},\ldots,X_{n}$ of mutually independent Bernoullis such that $\sum_{i=1}^{n}X_{i}$ is distributed according to $D$ and ${\bf E}[X_{1}]\leq\ldots\leq{\bf E}[X_{n}]$. (Such a collection exists by the definition of a Poisson binomial distribution.) Then we set $p_{i}={\bf E}[X_{i}]$ for all $i$. Lemma 1 of [DP13] shows that the resulting vector $(p_{1},\ldots,p_{n})$ is unique.

We denote by ${\rm PBD}(p_{1},\ldots,p_{n})$ the distribution of the sum $\sum_{i=1}^{n}X_{i}$ of mutually independent indicators $X_{1},\ldots,X_{n}$ with expectations $p_{i}={\bf E}[X_{i}]$, for all $i$. Given the above discussion ${\rm PBD}(p_{1},\ldots,p_{n})$ is unique up to permutation of the $p_{i}$’s. We also sometimes write $\{X_{i}\}$ to denote the distribution of $\sum_{i=1}^{n}X_{i}.$ Note the difference between $\{X_{i}\}$, which refers to the distribution of $\sum_{i}X_{i}$, and $\{X_{i}\}_{i}$, which refers to the underlying collection of mutually independent Bernoulli random variables.

Translated Poisson Distribution.

We will make use of the translated Poisson distribution for approximating the Poisson Binomial distribution. We define the translated Poisson distribution, and state a known result on how well it approximates the Poisson Binomial distribution.

We say that an integer random variable $Y$ is distributed according to the translated Poisson distribution with parameters $\mu$ and $\sigma^{2}$, denoted $TP(\mu,\sigma^{2})$, iff $Y$ can be written as

where $Z$ is a random variable distributed according to ${\rm Poisson}(\sigma^{2}+\{\mu-\sigma^{2}\})$, where $\{\mu-\sigma^{2}\}$ represents the fractional part of $\mu-\sigma^{2}$.

The following lemma gives a useful bound on the variation distance between a Poisson Binomial Distribution and a suitable translated Poisson distribution. Note that if the variance of the Poisson Binomial Distribution is large, then the lemma gives a strong bound.

Let $J_{1},\ldots,J_{n}$ be independent random indicators with ${\bf E}[J_{i}]=p_{i}$. Then

where $\mu=\sum_{i=1}^{n}p_{i}$ and $\sigma^{2}=\sum_{i=1}^{n}p_{i}(1-p_{i})$.

The following bound on the total variation distance between translated Poisson distributions will be useful.

Running Times, and Bit Complexity.

Throughout this paper, we measure the running times of our algorithms in numbers of bit operations. For a positive integer $n$, we denote by ${\langle{n}\rangle}$ its description complexity in binary, namely ${\langle{n}\rangle}=\lceil\log_{2}n\rceil$. Moreover, we represent a positive rational number $q$ as ${q_{1}\over q_{2}}$, where $q_{1}$ and $q_{2}$ are relatively prime positive integers. The description complexity of $q$ is defined to be ${\langle{q}\rangle}={\langle{q_{1}}\rangle}+{\langle{q_{2}}\rangle}$. We will assume that all $\epsilon$’s and $\delta$’s input to our algorithms are rational numbers.

In this section, we prove Theorem 1 by providing a sample- and time-efficient algorithm for learning an unknown PBD $X=\sum_{i=1}^{n}X_{i}$. We start with an important ingredient in our analysis.

A cover for PBDs. We make use of the following theorem, which provides a cover of the set ${\cal S}={\cal S}_{n}$ of all PBDs of order-$n$. The theorem was given implicitly in [DP11] and explicitly as Theorem 1 in [DP13].

For all $\epsilon>0$, there exists an $\epsilon$-cover ${\cal S}_{\epsilon}\subseteq{\cal S}$ of ${\cal S}$ such that

$|{\cal S}_{\epsilon}|\leq n^{2}+n\cdot\left({1\over\epsilon}\right)^{O(\log^{2}{1/\epsilon})}$; and

${\cal S}_{\epsilon}$ can be constructed in time linear in its representation size, i.e., ${O}(n^{2}\log n)+{O}(n\log n)\cdot\left({1\over\epsilon}\right)^{O(\log^{2}{1/\epsilon})}$.

Moreover, if $\{Y_{i}\}\in{\cal S}_{\epsilon}$, then the collection of $n$ Bernoulli random variables $\{Y_{i}\}_{i=1,\dots,n}$ has one of the following forms, where $k=k(\epsilon)\leq C/\epsilon$ is a positive integer, for some absolute constant $C>0$:

Finally, for every $\{X_{i}\}\in{\cal S}$ for which there is no $\epsilon$-neighbor in ${\cal S}_{\epsilon}$ that is in sparse form, there exists some $\{Y_{i}\}\in{\cal S}_{\epsilon}$ in $k$-heavy Binomial form such that

We remark that the cover theorem as stated in [DP13] does not include the part of the above statement following “finally.” We provide a proof of this extension in Appendix A.

The Basic Learning Algorithm. The high-level structure of our learning algorithms which give Theorem 1 is provided in Algorithm Learn-PBD of Figure 1. We instantiate this high-level structure, with appropriate technical modifications, in Section 2.4, where we give more detailed descriptions of the non-proper and proper algorithms that give parts (1) and (2) of Theorem 1.

At a high level, the subroutine Learn-Sparse is given sample access to $X$ and is designed to find an $\epsilon$-accurate hypothesis $H_{S}$ with probability at least $1-\delta/3$, if the unknown PBD $X$ is $\epsilon$-close to some sparse form PBD inside the cover ${\cal S}_{\epsilon}$. Similarly, Learn-Poisson is designed to find an $\epsilon$-accurate hypothesis $H_{P}$, if $X$ is not $\epsilon$-close to a sparse form PBD (in this case, Theorem 4 implies that $X$ must be $\epsilon$-close to some $k(\epsilon)$-heavy Binomial form PBD). Finally, Choose-Hypothesis is designed to choose one of the two hypotheses $H_{S},H_{P}$ as being $\epsilon$-close to $X.$ The following subsections specify these subroutines, as well as how the algorithm can be used to establish Theorem 1. We note that Learn-Sparse and Learn-Poisson do not return the distributions $H_{S}$ and $H_{P}$ as a list of probabilities for every point in $[n]$. They return instead a succinct description of these distributions in order to keep the running time of the algorithm logarithmic in $n$. Similarly, Choose-Hypothesis operates with succinct descriptions of these distributions.

Our starting point here is the simple observation that any PBD is a unimodal distribution over the domain $\{0,1,\dots,n\}$. (There is a simple inductive proof of this, or see Section 2 of [KG71].) This enables us to use the algorithm of Birgé [Bir97] for learning unimodal distributions. We recall Birgé’s result, and refer the reader to Appendix B for an explanation of how Theorem 5 as stated below follows from [Bir97].

For all $n,\epsilon,\delta>0$, there is an algorithm that draws

samples from an unknown unimodal distribution $X$ over $[n]$, does

The main result of this subsection is the following:

For all $n,\epsilon^{\prime},\delta^{\prime}>0$, there is an algorithm Learn-Sparse${}^{X}(n,\epsilon^{\prime},\delta^{\prime})$ that draws

samples from a target PBD $X$ over $[n]$, does

The high-level idea of Lemma 3 is quite simple. We truncate $O(\epsilon^{\prime})$ of the probability mass from each end of $X$ to obtain a conditional distribution $X_{[\hat{a},\hat{b}]}$; since $X$ is unimodal so is $X_{[\hat{a},\hat{b}]}$. If $\hat{b}-\hat{a}$ is larger than $O(1/\epsilon^{\prime 3})$ then the algorithm outputs “fail” (and $X$ could not have been close to a sparse-form distribution in the cover). Otherwise, we use Birgé’s algorithm to learn the unimodal distribution $X_{[\hat{a},\hat{b}]}$. A detailed description of the algorithm is given in Figure 2 below.

Proof of Lemma 3: As described in Figure 2, algorithm Learn-Sparse${}^{X}(n,\epsilon^{\prime},\delta^{\prime})$ first draws $M=32\log(8/\delta^{\prime})/\epsilon^{\prime 2}$ samples from $X$ and sorts them to obtain a list of values $0\leq s_{1}\leq\cdots\leq s_{M}\leq n.$ We claim the following about the values $\hat{a}$ and $\hat{b}$ defined in Step 2 of the algorithm:

With probability at least $1-\delta^{\prime}/2$, we have $X(\leq\hat{a})\in[3\epsilon^{\prime}/2,5\epsilon^{\prime}/2]$ and $X(\leq\hat{b})\in[1-5\epsilon^{\prime}/2,1-3\epsilon^{\prime}/2]$.

We only show that $X(\leq\hat{a})\geq 3\epsilon^{\prime}/2$ with probability at least $1-\delta^{\prime}/8$, since the arguments for $X(\leq\hat{a})\leq 5\epsilon^{\prime}/2$, $X(\leq\hat{b})\leq 1-3\epsilon^{\prime}/2$ and $X(\leq\hat{b})\geq 1-5\epsilon^{\prime}/2$ are identical. Given that each of these conditions is met with probability at least $1-\delta^{\prime}/8$, the union bound establishes our claim.

To show that $X(\leq\hat{a})\geq 3\epsilon^{\prime}/2$ is satisfied with probability at least $1-\delta^{\prime}/8$ we argue as follows: Let $\alpha^{\prime}=\max\{i~{}|~{}X(\leq i)<3\epsilon^{\prime}/2\}$. Clearly, $X(\leq\alpha^{\prime})<3\epsilon^{\prime}/2$ while $X(\leq\alpha^{\prime}+1)\geq 3\epsilon^{\prime}/2$. Given this, if $M$ samples are drawn from $X$ then the expected number of them that are $\leq\alpha^{\prime}$ is at most $3\epsilon^{\prime}M/2$. It follows then from the Chernoff bound that the probability that more than ${7\over 4}\epsilon^{\prime}M$ samples are $\leq\alpha^{\prime}$ is at most $e^{-(\epsilon^{\prime}/4)^{2}M/2}\leq\delta^{\prime}/8$. Hence except with this failure probability, we have $\hat{a}\geq\alpha^{\prime}+1$, which implies that $X(\leq\hat{a})\geq 3\epsilon^{\prime}/2$. ∎

As specified in Steps 3 and 4, if $\hat{b}-\hat{a}>(C/\epsilon^{\prime})^{3}$, where $C$ is the constant in the statement of Theorem 4, the algorithm outputs “fail”, returning the trivial hypothesis which puts probability mass $1$ on the point . Otherwise, the algorithm runs Birgé’s unimodal distribution learner (Theorem 5) on the conditional distribution $X_{[\hat{a},\hat{b}]}$, and outputs the result of Birgé’s algorithm. Since $X$ is unimodal, it follows that $X_{[\hat{a},\hat{b}]}$ is also unimodal, hence Birgé’s algorithm is appropriate for learning it. The way we apply Birgé’s algorithm to learn $X_{[\hat{a},\hat{b}]}$ given samples from the original distribution $X$ is the obvious one: we draw samples from $X$, ignoring all samples that fall outside of $[\hat{a},\hat{b}]$, until the right $O(\log(1/\delta^{\prime})\log(1/\epsilon^{\prime})/\epsilon^{\prime 3})$ number of samples fall inside $[\hat{a},\hat{b}]$, as required by Birgé’s algorithm for learning a distribution of support of size $(C/\epsilon^{\prime})^{3}$ with probability at least $1-\delta^{\prime}/4$. Once we have the right number of samples in $[\hat{a},\hat{b}]$, we run Birgé’s algorithm to learn the conditional distribution $X_{[\hat{a},\hat{b}]}$. Note that the number of samples we need to draw from $X$ until the right $O(\log(1/\delta^{\prime})\log(1/\epsilon^{\prime})/\epsilon^{\prime 3})$ number of samples fall inside $[\hat{a},\hat{b}]$ is still $O(\log(1/\delta^{\prime})\log(1/\epsilon^{\prime})/\epsilon^{\prime 3})$, with probability at least $1-\delta^{\prime}/4$. Indeed, since $X([\hat{a},\hat{b}])=1-O(\epsilon^{\prime})$, it follows from the Chernoff bound that with probability at least $1-\delta^{\prime}/4$, if $K=\Theta(\log(1/\delta^{\prime})\log(1/\epsilon^{\prime})/\epsilon^{\prime 3})$ samples are drawn from $X$, at least $K(1-O(\epsilon^{\prime}))$ fall inside $[\hat{a},\hat{b}]$.

Analysis: It is easy to see that the sample complexity of our algorithm is as promised. For the running time, notice that, if Birgé’s algorithm is invoked, it will return two lists of numbers $a_{1}$ through $a_{k}$ and $b_{1}$ through $b_{k}$, as well as a list of probability masses $q_{1},\ldots,q_{k}$ assigned to each subinterval $[a_{i},b_{i}]$, $i=1,\ldots,k$, by the hypothesis distribution $H_{S}$, where $k=O(\log(1/\epsilon^{\prime})/\epsilon^{\prime})$. In linear time, we can compute a list of probabilities $\hat{q}_{1},\ldots,\hat{q}_{k}$, representing the probability assigned by $H_{S}$ to every point of subinterval $[a_{i},b_{i}]$, for $i=1,\ldots,k$. So we can represent our output hypothesis $H_{S}$ via a data structure that maintains $O(1/\epsilon^{\prime 3})$ pointers, having one pointer per point inside $[a,b]$. The pointers map points to probabilities assigned by $H_{S}$ to these points. Thus turning the output of Birgé’s algorithm into an explicit distribution over $[a,b]$ incurs linear overhead in our running time, and hence the running time of our algorithm is also as promised. (See Appendix B for an explanation of the running time of Birgé’s algorithm.) Moreover, we also note that the output distribution has the promised structure, since in one case it has a single atom at and in the other case it is the output of Birgé’s algorithm on a distribution of support of size $(C/\epsilon^{\prime})^{3}$.

2 Learning when X𝑋X is close to a k𝑘k-heavy Binomial Form PBD.

For all $n,\epsilon^{\prime},\delta^{\prime}>0$, there is an algorithm Learn-Poisson${}^{X}(n,\epsilon^{\prime},\delta^{\prime})$ that draws

samples from a target PBD $X$ over $[n]$, does

Our proof plan is to exploit the structure of the cover of Theorem 4. In particular, if $X$ is not $\epsilon^{\prime}$-close to any sparse form PBD in the cover, it must be $\epsilon^{\prime}$-close to a PBD in heavy Binomial form with approximately the same mean and variance as $X$, as specified by the final part of the cover theorem. Hence, a natural strategy is to obtain estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ of the mean and variance of the unknown PBD $X$, and output as a hypothesis a translated Poisson distribution with parameters $\hat{\mu}$ and $\hat{\sigma}^{2}$. We show that this strategy is a successful one. Before providing the details, we highlight two facts that we will establish in the subsequent analysis and that will be used later. The first is that, assuming $X$ is not $\epsilon^{\prime}$-close to any sparse form PBD in the cover ${\cal S}_{\epsilon^{\prime}}$, its variance $\sigma^{2}$ satisfies

The second is that under the same assumption, the estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ of the mean $\mu$ and variance $\sigma^{2}$ of $X$ that we obtain satisfy the following bounds with probability at least $1-\delta$:

See Figure 3 and the associated Figure 4 for a detailed description of the Learn-Poisson${}^{X}(n,\epsilon^{\prime},\delta^{\prime})$ algorithm.

Proof of Lemma 5: We start by showing that we can estimate the mean and variance of the target PBD $X$.

We treat the estimation of $\mu$ and $\sigma^{2}$ separately. For both estimation problems we show how to use $O(1/\epsilon^{2})$ samples to obtain estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ achieving the required guarantees with probability at least $2/3$ (we refer to these as “weak estimators”). Then a routine procedure allows us to boost the success probability to $1-\delta$ at the expense of a multiplicative factor $O(\log 1/\delta)$ on the number of samples. While we omit the details of the routine boosting argument, we remind the reader that it involves running the weak estimator $O(\log 1/\delta)$ times to obtain estimates $\hat{\mu}_{1},\ldots,\hat{\mu}_{O(\log{1/\delta})}$ and outputting the median of these estimates, and similarly for estimating $\sigma^{2}$.

We proceed to specify and analyze the weak estimators for $\mu$ and $\sigma^{2}$ separately:

Weak estimator for $\mu$: Let $Z_{1},\ldots,Z_{m}$ be independent samples from $X$, and let $\hat{\mu}={\sum_{i}Z_{i}\over m}$. Then

Choosing $t=\sqrt{3}$ and $m=\lceil 3/\epsilon^{2}\rceil$, the above imply that $|\hat{\mu}-\mu|\leq\epsilon\sigma$ with probability at least $2/3$.

Weak estimator for $\sigma^{2}$: Let $Z_{1},\ldots,Z_{m}$ be independent samples from $X$, and let $\hat{\sigma}^{2}={\sum_{i}(Z_{i}-{1\over m}\sum_{i}Z_{i})^{2}\over m-1}$ be the unbiased sample variance. (Note the use of Bessel’s correction.) Then it can be checked [Joh03] that

where $\kappa$ is the excess kurtosis of the distribution of $X$ (i.e. $\kappa={\frac{{\bf E}[(X-\mu)^{4}]}{\sigma^{4}}}-3$). To bound $\kappa$ in terms of $\sigma^{2}$ suppose that $X=\sum_{i=1}^{n}X_{i}$, where ${\bf E}[X_{i}]=p_{i}$ for all $i$. Then

Choosing $t=\sqrt{3}$ and $m=\lceil 3/\epsilon^{2}\rceil$, the above imply that $|\hat{\sigma}^{2}-\sigma^{2}|\leq\epsilon\sigma^{2}\sqrt{4+{1\over\sigma^{2}}}$ with probability at least $2/3$.

We proceed to prove Lemma 5. Learn-Poisson${}^{X}(n,\epsilon^{\prime},\delta^{\prime})$ runs ${\cal A}(n,\epsilon,\delta)$ from Lemma 6 with appropriately chosen $\epsilon=\epsilon(\epsilon^{\prime})$ and $\delta=\delta(\delta^{\prime})$, given below, and then outputs the translated Poisson distribution $TP(\hat{\mu},\hat{\sigma}^{2})$, where $\hat{\mu}$ and $\hat{\sigma}^{2}$ are the estimated mean and variance of $X$ output by ${\cal A}$. Next, we show how to choose $\epsilon$ and $\delta$, as well as why the desired guarantees are satisfied by the output distribution.

In particular, conditions (b) and (d) above imply that

for some universal constant $\theta$, establishing (1). In terms of this $\theta$, we choose $\epsilon=\epsilon^{\prime}/\sqrt{4+{1\over\theta^{2}}}$ and $\delta=\delta^{\prime}$ for the application of Lemma 6 to obtain—from $O(\log(1/\delta^{\prime})/\epsilon^{\prime 2})$ samples—estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ of $\mu$ and $\sigma^{2}$.

From our choice of parameters and the guarantees of Lemma 6, it follows that, if $X$ is not $\epsilon^{\prime}$-close to any PBD in sparse form inside the cover ${\cal S}_{\epsilon^{\prime}}$, then with probability at least $1-\delta^{\prime}$ the estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ satisfy:

establishing (2). Moreover, if $Y$ is a random variable distributed according to the translated Poisson distribution $TP(\hat{\mu},\hat{\sigma}^{2})$, we show that $X$ and $Y$ are within $O(\epsilon^{\prime})$ in total variation distance, concluding the proof of Lemma 5.

We make use of Lemma 1. Suppose that $X=\sum_{i=1}^{n}X_{i}$, where ${\bf E}[X_{i}]=p_{i}$ for all $i$. Lemma 1 implies that

It remains to bound the total variation distance between the translated Poisson distributions $TP(\mu,\sigma^{2})$ and $TP(\hat{\mu},\hat{\sigma}^{2})$. For this we use Lemma 2. Lemma 2 implies

The claim follows from (3), (4) and the triangle inequality. ∎

The proof of Lemma 5 is concluded. We remark that the algorithm described above does not need to know a priori whether or not $X$ is $\epsilon^{\prime}$-close to a PBD in sparse form inside the cover ${\cal S}_{\epsilon^{\prime}}$ of Theorem 4. The algorithm simply runs the estimator of Lemma 6 with $\epsilon=\epsilon^{\prime}/\sqrt{4+{1\over\theta^{2}}}$ and $\delta^{\prime}=\delta$ and outputs whatever estimates $\hat{\mu}$ and $\hat{\sigma}^{2}$ the algorithm of Lemma 6 produces. $\square$

3 Hypothesis testing.

Our hypothesis testing routine Choose-HypothesisX uses samples from the unknown distribution $X$ to run a “competition” between two candidate hypothesis distributions $H_{1}$ and $H_{2}$ over $[n]$ that are given in the input. We show that if at least one of the two candidate hypotheses is close to the unknown distribution $X$, then with high probability over the samples drawn from $X$ the routine selects as winner a candidate that is close to $X$. This basic approach of running a competition between candidate hypotheses is quite similar to the “Scheffé estimate” proposed by Devroye and Lugosi (see [DL96b, DL96a] and Chapter 6 of [DL01], as well as [Yat85]), but our notion of competition here is different.

We obtain the following lemma, postponing all running-time analysis to the next section.

There is an algorithm Choose-Hypothesis${}^{X}({H}_{1},{H}_{2},\epsilon^{\prime},\delta^{\prime})$ which is given sample access to distribution $X$, two hypothesis distributions $H_{1},H_{2}$ for $X$, an accuracy parameter $\epsilon^{\prime}>0$, and a confidence parameter $\delta^{\prime}>0.$ It makes

Proof of Lemma 8: Figure 5 describes how the competition between $H_{1}$ and $H_{2}$ is carried out.

The correctness of Choose-Hypothesis is an immediate consequence of the following claim. (In fact for Lemma 8 we only need item (i) below, but item (ii) will be handy later in the proof of Lemma 10.)

For part (ii): If $p_{1}-p_{2}\leq 5\epsilon^{\prime}$ then the competition declares a draw, hence $H_{3-i}$ is not the winner. Otherwise we have $p_{1}-p_{2}>5\epsilon^{\prime}$ and the above arguments imply that the competition between $H_{1}$ and $H_{2}$ will declare $H_{3-i}$ as the winner with probability at most $2e^{-{m\epsilon^{\prime 2}/2}}$.

In view of Claim 9, the proof of Lemma 8 is concluded. $\square$

Our Choose-Hypothesis algorithm implies a generic learning algorithm of independent interest.

Let ${\cal S}$ be an arbitrary set of distributions over a finite domain. Moreover, let ${\cal S}_{\epsilon}\subseteq{\cal S}$ be an $\epsilon$-cover of ${\cal S}$ of size $N$, for some $\epsilon>0$. For all $\delta>0$, there is an algorithm that uses

The algorithm performs a tournament, by running Choose-Hypothesis${}^{X}(H_{i},H_{j},\epsilon,$ $\delta/(4N))$ for every pair $(H_{i},H_{j})$, $i<j$, of distributions in ${\cal S}_{\epsilon}$. Then it outputs any distribution $Y^{\star}\in{\cal S}_{\epsilon}$ that was never a loser (i.e., won or tied against all other distributions in the cover). If no such distribution exists in ${\cal S}_{\epsilon}$ then the algorithm says “failure,” and outputs an arbitrary distribution from ${\cal S}_{\epsilon}$.

We note that Devroye and Lugosi (Chapter 7 of [DL01]) prove a similar result, but there are some differences. They also have all pairs of distributions in the cover compete against each other, but they use a different notion of competition between every pair. Moreover, their approach chooses a distribution in the cover that wins the maximum number of competitions, whereas our algorithm chooses a distribution that is never defeated (i.e., won or tied against all other distributions in the cover).

Recent work [DK14, AJOS14, SOAJ14] improves the running time of the tournament approaches of Lemma 10, Devroye-Lugosi and other related tournaments to have a quasilinear dependence of $O(N\log N)$ on the size $N=|{\cal S}_{\epsilon}|$ . In particular, they avoid running Choose-Hypothesis for all pairs of distributions in ${\cal S}_{\epsilon}$.

4 Proof of Theorem 1.

We first show Part (1) of the theorem, where the learning algorithm may output any distribution over $[n]$ and not necessarily a PBD. The algorithm for this part of the theorem, Non-Proper-Learn-PBD, is given in Figure 6. This algorithm follows the high-level structure outlined in Figure 1 with the following modifications: (a) first, if the total variation distance to within which we want to learn $X$ is $\epsilon$, the second argument of both Learn-Sparse and Learn-Poisson is set to ${\epsilon\over 12\max\{c_{1},c_{2}\}}$, where $c_{1}$ and $c_{2}$ are respectively the constants from Lemmas 3 and 5; (b) the third step of Learn-PBD is replaced by Choose-Hypothesis${}^{X}(H_{S},{\widehat{H_{P}}},{\epsilon/8},\delta/3)$, where ${\widehat{H_{P}}}$ is defined in terms of $H_{P}$ as described in Definition 2 below; and (c) if Choose-Hypothesis returns $H_{S}$, then Learn-PBD also returns $H_{S}$, while if Choose-Hypothesis returns $\widehat{H_{P}}$, then Learn-PBD returns $H_{P}$.

(Definition of ${\widehat{H_{P}}}$:) ${\widehat{H_{P}}}$ is defined in terms of $H_{P}$ and the support of $H_{S}$ in three steps:

for all points $i$ such that $H_{S}(i)=0$, we let ${\widehat{H_{P}}}(i)={{H_{P}}}(i)$;

for all points $i$ such that $H_{S}(i)\neq 0$, we describe in Appendix C an efficient deterministic algorithm that numerically approximates $H_{P}(i)$ to within an additive error of $\pm\epsilon/48s$, where $s=O(1/\epsilon^{3})$ is the cardinality of the support of $H_{S}$. If $\widehat{H_{P,i}}$ is the approximation to $H_{P}(i)$ output by the algorithm, we set $\widehat{H_{P}}(i)={\max\{0,\widehat{H_{P,i}}-\epsilon/48s\}}$; notice then that $H_{P}(i)-\epsilon/24s\leq\widehat{H_{P}}(i)\leq H_{P}(i)$; finally,

for an arbitrary point $i$ such that $H_{S}(i)=0$, we set $\widehat{H_{P}}(i)=1-\sum_{j\neq i}\widehat{H_{P}}(j)$, to make sure that $\widehat{H_{P}}$ is a probability distribution.

We remark that the reason why we do not wish to use $H_{P}$ directly in Choose-Hypothesis is purely computational. In particular, since $H_{P}$ is a translated Poisson distribution, we cannot compute its probabilities $H_{P}(i)$ exactly, and we need to approximate them. On the other hand, we need to make sure that using approximate values will not cause Choose-Hypothesis to make a mistake. Our $\widehat{H_{P}}$ is carefully defined so as to make sure that Choose-Hypothesis selects a probability distribution that is close to the unknown $X$, and that all probabilities that Choose-Hypothesis needs to compute can be computed without much overhead. In particular, we remark that, in running Choose-Hypothesis, we do not a priori compute the value of $\widehat{H_{P}}$ at every point; we do instead a lazy evaluation of $\widehat{H_{P}}$, as explained in the running-time analysis below.

Now we turn to Part (2) of Theorem 1, the proper learning result. The algorithm for this part of the theorem, Proper-Learn-PBD, is given in Figure 7. The algorithm is essentially the same as Non-Proper-Learn-PBD but with the following modifications, to produce a PBD that is within $O(\epsilon)$ of the unknown $X$: First, we replace Learn-Sparse with a different learning algorithm, Proper-Learn-Sparse, which is based on Lemma 10, and always outputs a PBD. Second, we add a post-processing step to Learn-Poisson that converts the translated Poisson distribution $H_{P}$ output by this procedure to a PBD (in fact, to a Binomial distribution). After we describe these new ingredients in detail, we explain and analyze our proper learning algorithm.

The procedure Proper-Learn-Sparse${}^{X}(n,\epsilon,\delta)$ is given in Figure 8; we explain the procedure in tandem with a proof of correctness. As in Learn-Sparse, we start by truncating $\Theta(\epsilon)$ of the probability mass from each end of $X$ to obtain a conditional distribution $X_{[\hat{a},\hat{b}]}$. In particular, we compute $\hat{a}$ and $\hat{b}$ as described in the beginning of the proof of Lemma 3 (setting $\epsilon^{\prime}=\epsilon$ and $\delta^{\prime}=\delta$). Claim 4 implies that, with probability at least $1-\delta/2$, $X(\leq\hat{a}),1-X(\leq\hat{b})\in[3\epsilon/2,5\epsilon/2]$. (Let us denote this event by ${\cal G}$.) We distinguish the following cases:

If $\hat{b}-\hat{a}>\omega=(C/\epsilon)^{3}$, where $C$ is the constant in the statement of Theorem 4, the algorithm outputs “fail,” returning the trivial hypothesis that puts probability mass $1$ on the point . Observe that, if $\hat{b}-\hat{a}>\omega$ and $X(\leq\hat{a}),1-X(\leq\hat{b})\in[3\epsilon/2,5\epsilon/2]$, then $X$ cannot be $\epsilon$-close to a sparse-form distribution in the cover.

Locate-Binomial$(\hat{\mu},\hat{\sigma}^{2},n)$: This routine takes as input the output $(\hat{\mu},\hat{\sigma}^{2})$ of Learn-Poisson${}^{X}(n,\epsilon,\delta)$ and computes a Binomial distribution $H_{B}$, without any additional samples from $X$. The guarantee is that, if $X$ is not $\epsilon$-close to any sparse form distribution in the cover $S_{\epsilon}$ of Theorem 4, then, with probability at least $1-\delta$ (over the randomness in the output of Learn-Poisson), $H_{B}$ will be $O(\epsilon)$-close to $X$.

Locate-Binomial is presented in Figure 9; we proceed to explain the algorithm and establish its correctness. This routine has three steps. The first two eliminate corner-cases in the values of $\hat{\mu}$ and $\hat{\sigma}^{2}$, while the last step defines a Binomial distribution $H_{B}\equiv{\rm Bin}(\hat{n},\hat{p})$ with $\hat{n}\leq n$ that is $O(\epsilon)$-close to $H_{P}\equiv TP(\hat{\mu},\hat{\sigma}^{2})$ and hence to $X$ under our working assumptions. (We note that a significant portion of the work below is to ensure that $\hat{n}\leq n$, which does not seem to follow from a more direct approach. Getting $\hat{n}\leq n$ is necessary in order for our learning algorithm for order-$n$ PBDs to be truly proper.) Throughout (a), (b) and (c) below we assume that our working assumptions hold. In particular, our assumptions are used every time we employ the bounds (1) and (2) of Section 2.2.

Tweaking $\hat{\sigma}^{2}$: If $\hat{\sigma}^{2}\leq{n\over 4}$, we set $\sigma_{1}^{2}=\hat{\sigma}^{2}$; otherwise, we set $\sigma_{1}^{2}={n\over 4}$. (As intuition for this tweak, observe that the largest possible variance of a Binomial distribution ${\rm Bin}(n,\cdot)$ is $n/4.$) We note for future reference that in both cases (2) gives

where the lower bound follows from (2) and the fact that any PBD satisfies $\sigma^{2}\leq{n\over 4}$.

where we used the fact that $\sigma^{2}=\Omega(1/\epsilon^{2})$ from (1).

Observe first that $\sigma_{1}^{2}>\sigma_{2}^{2}$ and $\sigma_{2}^{2}\geq 0$, where the last assertion follows from the fact that $\hat{\mu}\leq n$ by construction.

Next, suppose that $X=PBD(p_{1},\ldots,p_{n})$. Then from Cauchy-Schwarz we get that

where the first inequality follows from (2), the second inequality follows from (7) and the fact that any PBD over $n$ variables satisfies $\mu\leq n,$ and the last one from (1).

where we used that $\sigma^{2}=\Omega(1/\epsilon^{2})$ from (1).

Note that, from the way that $\sigma_{2}^{2}$ is set in Step (b) above, we have that $\hat{n}\leq n$ and $\hat{p}\in$, as required for ${\rm Bin}(\hat{n},\hat{p})$ to be a valid Binomial distribution and a valid output for Part 2 of Theorem 1.

Let us bound the total variation distance between ${\rm Bin}(\hat{n},\hat{p})$ and $TP(\hat{\mu},\sigma_{2}^{2})$. First, using Lemma 1 we have:

where the second inequality uses (8) (or (5) depending on which case of Step (b) we fell into) and the last one uses the fact that $\sigma^{2}=\Omega(1/\epsilon^{2})$ from (1). So plugging this into (10) we get:

The next step is to compare $TP(\hat{n}\hat{p},\hat{n}\hat{p}(1-\hat{p}))$ and $TP(\hat{\mu},\sigma_{2}^{2})$. Lemma 2 gives:

Learning weighted sums of independent Bernoullis

Theorem 2. Let $X=\sum_{i=1}^{n}a_{i}X_{i}$ be a weighted sum of unknown independent Bernoulli random variables such that there are at most $k$ different values in the set $\{a_{1},\dots,a_{n}\}.$ Then there is an algorithm with the following properties: given $n,$ $a_{1},\dots,a_{n}$ and access to independent draws from $X$, it uses

samples from the target distribution $X$, runs in time

Given a vector $\overline{a}=(a_{1},\dots,a_{n})$ of weights, we refer to a distribution $X=\sum_{i=1}^{n}a_{i}X_{i}$ (where $X_{1},\dots,X_{n}$ are independent Bernoullis which may have arbitrary means) as an $\overline{a}$-weighted sum of Bernoullis, and we write ${\cal S}_{\overline{a}}$ to denote the space of all such distributions.

To prove Theorem 2 we first show that ${\cal S}_{\overline{a}}$ has an $\epsilon$-cover that is not too large. We then show that by running a “tournament” between all pairs of distributions in the cover, using the hypothesis testing subroutine from Section 2.3, it is possible to identify a distribution in the cover that is close to the target $\overline{a}$-weighted sum of Bernoullis.

Let $\{b_{j}\}_{j=1}^{k}$ denote the set of distinct weights in $a_{1},\dots,a_{n}$, and let n_{j}=\big{|}\{i\in[n]\mid a_{i}=b_{j}\}\big{|}. With this notation, we can write $X=\mathop{{\textstyle\sum}}_{j=1}^{k}b_{j}S_{j}=g(S)$, where $S=(S_{1},\ldots,S_{k})$ with each $S_{j}$ a sum of $n_{j}$ many independent Bernoulli random variables and $g(y_{1},\ldots,y_{k})=\sum_{j=1}^{k}b_{j}y_{j}$. Clearly we have $\mathop{{\textstyle\sum}}_{j=1}^{k}n_{j}=n$. By Theorem 4, for each $j\in\{1,\dots,k\}$ the space of all possible $S_{j}$’s has an explicit $(\epsilon/k)$-cover $\mathcal{S}^{j}_{\epsilon/k}$ of size $|\mathcal{S}^{j}_{\epsilon/k}|\leq n_{j}^{{2}}+n\cdot(k/\epsilon)^{O(\log^{2}(k/\epsilon))}$. By independence across $S_{j}$’s, the product $\mathcal{Q}=\mathop{{\textstyle\prod}}_{j=1}^{k}\mathcal{S}^{j}_{\epsilon/k}$ is an $\epsilon$-cover for the space of all possible $S$’s, and hence the set

is an $\epsilon$-cover for ${\cal S}_{\overline{a}}.$ So ${\cal S}_{\overline{a}}$ has an explicit $\epsilon$-cover of size $|\mathcal{Q}|=\mathop{{\textstyle\prod}}_{j=1}^{k}|\mathcal{S}^{j}_{\epsilon/k}|\leq(n/k)^{{2}k}\cdot(k/\epsilon)^{k\cdot O(\log^{2}(k/\epsilon))}$. ∎

where the sum is over all $k$-tuples $(m_{1},\dots,m_{k})$ such that $0\leq m_{j}\leq n_{j}$ for all $j$ and $b_{1}m_{1}+\cdots+b_{k}m_{k}=w$ (as noted above there are at most $O((n/k)^{k})$ such $k$-tuples). To complete the proof of Theorem 2 we note that $\Pr_{H}[S_{j}=m_{j}]$ can be computed in $O(n_{j}^{2})$ time by standard dynamic programming. $\square$

2 Sample complexity lower bound for learning sums of weighted independent Bernoulli random variables

The intuition underlying this lower bound is straightforward: Suppose there are $n/100$ variables $X_{i}$, chosen uniformly at random, which have $p_{i}=100/n$ (call these the “relevant variables”), and the rest of the $p_{i}$’s are zero. Given at most $c\cdot n$ draws from $X$ for a small constant $c$, with high probability some constant fraction of the relevant $X_{i}$’s will not have been “revealed” as relevant, and from this it is not difficult to show that any hypothesis must have constant error. A detailed argument follows.

Proof of Theorem 3: We define a probability distribution over possible target probability distributions $X$ as follows: A subset $S\subset\{n/2+1,\dots,n\}$ of size $|S|=n/100$ is drawn uniformly at random from all ${n/2\choose n/100}$ possible outcomes.. The vector $\overline{p}=(p_{1},\dots,p_{n})$ is defined as follows: for each $i\in S$ the value $p_{i}$ equals $100/n=1/|S|,$ and for all other $i$ the value $p_{i}$ equals 0. The $i$-th Bernoulli random variable $X_{i}$ has ${\bf E}[X_{i}]=p_{i}$, and the target distribution is $X=X_{\overline{p}}=\sum_{i=1}^{n}iX_{i}.$

Fix any $S,\overline{p}$ as described above. For any $j\in\{n/2+1,\dots,n\}$ we have $X_{\overline{p}}(j)\neq 0$ if and only if $j\in S$. For any $j\in S$ the value $X_{\overline{p}}(j)$ is exactly $(100/n)(1-100/n)^{n/100-1}>35/n$ (for $n$ sufficiently large), and hence $X_{\overline{p}}(\{n/2+1,\dots,n\})>0.35$ (again for $n$ sufficiently large).

The first claim of the lemma holds because any set of $c\geq 2$ numbers from $\{n/2+1,\dots,n\}$ must sum to more than $n$. The second claim holds because the only way a draw $x$ from $X_{\overline{p}}$ can have $x=j$ is if $X_{j}=1$ and all other $X_{i}$ are 0 (here we are using $\lim_{x\rightarrow\infty}(1-1/x)^{x}=1/e$).

The next lemma is an easy consequence of Chernoff bounds:

Fix any $\overline{p}$ as defined above, and consider a sequence of $n/2000$ independent draws from $X_{\overline{p}}=\sum_{i}iX_{i}$. With probability $1-e^{-\Omega(n)}$ the total number of indices $j\in[n]$ such that $X_{j}$ is ever 1 in any of the $n/2000$ draws is at most $n/1000$.

We are now ready to prove Theorem 3. Let $L$ be a learning algorithm that receives $n/2000$ samples. Let $S\subset\{n/2+1,\dots,n\}$ and $\overline{p}$ be chosen randomly as defined above, and set the target to $X=X_{\overline{p}}.$

We consider an augmented learner $L^{\prime}$ that is given “extra information.” For each point in the sample, instead of receiving the value of that draw from $X$ the learner $L^{\prime}$ is given the entire vector $(X_{1},\dots,X_{n})\in\{0,1\}^{n}$. Let $T$ denote the set of elements $j\in\{n/2+1,\dots,n\}$ for which the learner is ever given a vector $(X_{1},\dots,X_{n})$ that has $X_{j}=1.$ By Lemma 16 we have $|T|\leq n/1000$ with probability at least $1-e^{-\Omega(n)}$; we condition on the event $|T|\leq n/1000$ going forth.

Conclusion and open problems

References

Appendix A Extension of the Cover Theorem: Proof of Theorem 4

Theorem 4 is restating the main cover theorem (Theorem 1) of [DP13], except that it claims an additional property, namely what follows the word “finally” in the statement of the theorem. (We will sometimes refer to this property as the last part of Theorem 4 in the following discussion.) Our goal is to show that the cover of [DP13] already satisfies this property without any modifications, thereby establishing Theorem 4. To avoid reproducing the involved constructions of [DP13], we will assume that the reader has some familiarity with them. Still, our proof here will be self-contained.

First, we note that the $\epsilon$-cover ${\cal S}_{\epsilon}$ of Theorem 1 of [DP13] is a subset of a larger ${\epsilon\over 2}$-cover ${\cal S}^{\prime}_{\epsilon/2}$ of size $n^{2}+n\cdot(1/\epsilon)^{O(1/\epsilon^{2})}$, which includes all the $k$-sparse and all the $k$-heavy Binomial PBDs (up to permutations of the underlying $p_{i}$’s), for some $k=O(1/\epsilon)$. Let us call ${\cal S}^{\prime}_{\epsilon/2}$ the “large ${\epsilon\over 2}$-cover” to distinguish it from ${\cal S}_{\epsilon}$, which we will call the “small $\epsilon$-cover.” The reader is referred to Theorem 2 in [DP13] (and the discussion following that theorem) for a description of the large $\epsilon\over 2$-cover, and to Section 3.2 of [DP13] for how this cover is used to construct the small $\epsilon$-cover. In particular, the small $\epsilon$-cover is a subset of the large ${\epsilon/2}$-cover, including only a subset of the sparse form distributions in the large ${\epsilon/2}$-cover. Moreover, for every sparse form distribution in the large $\epsilon/2$-cover, the small $\epsilon$-cover includes at least one sparse form distribution that is $\epsilon/2$-close in total variation distance. Hence, if the large $\epsilon/2$-cover satisfies the last part of Theorem 4 (with $\epsilon/2$ instead of $\epsilon$ and ${\cal S}_{\epsilon/2}^{\prime}$ instead of ${\cal S}_{\epsilon}$), it follows that the small $\epsilon$-cover also satisfies the last part of Theorem 4.