Optimal Learning via the Fourier Transform for Sums of Independent Integer Random Variables

Introduction

We study sums of independent integer random variables:

For convenience, throughout this paper, we will often blur the distinction between a random variable and its distribution. In particular, we will use the term $k$-SIIRV for the random variable or its corresponding distribution, and the distinction will be clear from the context.

Sums of independent integer random variables (SIIRVs) comprise a rich class of distributions that arise in many settings. The special case of $k=2$, ${\cal S}_{n,2}$, was first considered by Poisson [Poi37] as a non-trivial extension of the Binomial distribution, and is known as Poisson binomial distribution (PBD). In application domains, SIIRVs 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 practical uses of these distributions. We remark that these distributions are of fundamental interest and have been extensively studied in probability and statistics. For example, tail bounds on SIIRVs form an important special case of Chernoff/Hoeffding bounds [Che52, Hoe63, DP09b]. Moreover, there is a long line of research on approximate limit theorems for SIIRVs, dating back several decades (see e.g., [Pre83, Kru86, BHJ92]), and [CL10, CGS11] for some recent results.

The main motivation of this work was the problem of learning an unknown $k$-SIIRV given access to independent samples. Understanding this problem is intimately related to obtaining a refined structural understanding of the space of $k$-SIIRVs. The connection between structure and distribution learning is the main thrust of this paper.

While density estimation has been studied for several decades, the number of samples required to learn is not yet well understood, even for surprisingly simple and natural classes of univariate discrete distributions. More specifically, there is no known complexity measure of a distribution family $\mathcal{D}$ that characterizes the sample complexity of learning an unknown distribution from $\mathcal{D}.$ In contrast, the VC dimension of a concept class plays such a role in the PAC model of learning Boolean functions (see, e.g, [BEHW89, KV94]).

Obtaining a computationally efficient learning algorithm with optimal (or near-optimal) sample complexity is an important goal. In many learning settings, achieving this goal turns out to be quite challenging. More specifically, in many scenarios, both supervised and unsupervised, the only computationally efficient learning algorithms known use a (provably) suboptimal sample size. Intuitively, increasing the sample size (e.g., by a polynomial factor) can make the algorithmic task substantially easier. Characterizing the tradeoff between sample complexity and computational complexity is of fundamental importance in learning theory. In this work, we essentially characterize this tradeoff for the unsupervised problem of learning SIIRVs.

2 Our Results

The main technical contribution of this paper is the use of Fourier analytic and geometric tools to obtain a refined structural understanding of the space of $k$-SIIRVs. As a byproduct of our techniques, we characterize the sample complexity of learning $k$-SIIRVs (up to constant factors), and moreover we obtain a computationally efficient learning algorithm with near-optimal sample complexity. Our results answer the main open questions of [DDS12b, DDO+13].

Along the way we prove several new structural results of independent interest about $k$-SIIRVs, including: the approximate sparsity of their Fourier transform; tight upper and lower bounds on $\epsilon$-covers (in total variation distance and Kolmogorov distance); and a novel geometric characterization of the space of $k$-SIIRVs, that is crucial for our sample complexity lower bound. Below, we state our results in detail and elaborate on their context and the connections between them.

As our first result, we give a sample near-optimal and computationally efficient learning algorithm for $k$-SIIRVs:

Our algorithm outputs a succinct description of the hypothesis $\mathbf{H},$ via its Discrete Fourier Transform (DFT) ${\widehat{\mathbf{H}}},$ which is supported on a set of small cardinality. The DFT immediately gives a fast evaluation oracle for $\mathbf{H}.$ We also show how to use the DFT, in a black-box manner, to obtain an efficient approximate sampler for the target distribution $\mathbf{P}.$ Our efficient learning algorithm is described in Section 2.1. In Section 2.3 we give the efficient construction of our sampler.

Given our $\widetilde{O}(k/\epsilon^{2})$ sample upper bound, it would be tempting to conjecture that $\Theta(k/\epsilon^{2})$ is in fact the optimal sample complexity of learning $k$-SIIRVs. If true, this would imply that learning a $k$-SIIRV is as easy as learning a $k$-IRV. Surprisingly, we show that this is not the case:

Theorem 1.2 precisely characterizes the sample complexity of learning $k$-SIIRVs (up to constant factors) by giving an upper bound and a matching information-theoretic sample lower bound. The sharp sample complexity bound of $\Theta((k/\epsilon^{2})\sqrt{\log(1/\epsilon)})$ is surprising, and cannot be obtained using standard information-theoretic tools (e.g., metric entropy). We elaborate on this issue in Section 1.4.

We remark that the upper bound of Theorem 1.2 does not specify the running time of the corresponding algorithm. This is because the simplest such algorithm actually runs in time exponential in $k.$ For the important special case of $k=2,$ we obtain a sample–optimal learning algorithm that runs in sample–linear time:

For any $\epsilon>0,$ there is an algorithm that learns PBDs within variation distance $\epsilon$ using $O((1/\epsilon^{2})\sqrt{\log(1/\epsilon)})$ samples and running in time $O((1/\epsilon^{2})\sqrt{\log(1/\epsilon)}).$

The upper bound of Theorem 1.2 and Theorem 1.3 are established in Section 2.4. Our tight sample complexity lower bound is proved in Section 5.

Our learning upper bounds are obtained via an approach which is novel in this context. Specifically, we show that the Fourier transform of $k$-SIIRVs is approximately sparse, and exploit this property to learn the distribution via learning its Fourier transform in its effective support. The sparsity of the Fourier transform explains why this family of distributions is learnable with sample complexity independent of $n,$ and moreover it yields the sharp sample-complexity bound. The algorithmic idea of exploiting Fourier sparsity for distribution learning is general (see Section 2.2), and was subsequently used by the authors in other related settings [DKS15a, DKS15b].

Our core structural result is the following simple property of the Fourier transform of $k$-SIIRVs:

Any $k$-SIIRV with “large” variance has a Fourier transform with “small” effective support.

One can obtain different versions of the above informal statement depending on the setting and the desired application. See Lemma 2.3 for a formal statement in the context of the DFT. The Fourier sparsity of $k$-SIIRVs forms the basis for our upper bounds in this paper. As previously mentioned, this structural property motivates and enables our learning algorithm. Moreover, it is useful in order to obtain sparse $\epsilon$-covers for $\mathcal{S}_{n,k},$ the space of $k$-SIIRVs, under the total variation distance.

More specifically, using the approximate sparsity of the Fourier transform of SIIRVs combined with analytic arguments, we obtain a computationally efficient algorithm to construct a proper $\epsilon$-cover for $\mathcal{S}_{n,k},$ of near-minimum size. In particular, we show:

For $\epsilon\leq 1/k$, there exists a proper $\epsilon$-cover ${\cal S}_{n,k,\epsilon}\subseteq{\cal S}_{n,k}$ of ${\cal S}_{n,k}$ under the total variation distance of size $|{\cal S}_{n,k,\epsilon}|\leq n\cdot(1/\epsilon)^{O(k\log(1/\epsilon))}$ that can be constructed in polynomial time.

We also prove a matching lower bound on the cover size, showing that our above construction is essentially optimal:

Before our work, no non-trivial lower bound on the cover size was known. We view the inherent quasi-polynomial dependence on $1/\epsilon$ of the cover size established here as a rather surprising fact. Our cover size lower bound proof relies on a new geometric characterization of the space of $k$-SIIRVs that we believe is of independent interest, and may find other applications. Our tight lower bound on the sample complexity of learning $k$-SIIRVs relies critically on this characterization. Our cover size lower bound is proved in Section 4.

3 Preliminaries

We record a few definitions that will be used throughout this paper.

We emphasize that our learning algorithms output both an $\epsilon$-sampler and an $\epsilon$-evaluation oracle for the target distribution.

Let ${\cal F}$ be a family of probability distributions. Given $\delta>0$, a subset ${\cal G}\subseteq{\cal F}$ is said to be a proper $\delta$-cover of ${\cal F}$ with respect to the metric $d(\cdot,\cdot)$ if for every distribution $\mathbf{P}\in{\cal F}$ there exists some $\mathbf{Q}\in{\cal G}$ such that $d(\mathbf{P},\mathbf{Q})\leq\delta.$ If ${\cal G}$ is not a subset of ${\cal F},$ then the cover is called non-proper. The $\delta$-covering number for $({\cal F},d)$ is the minimum cardinality of a $\delta$-cover. The $\delta$-packing number for $({\cal F},d)$ is the maximum number of points (distributions) in $\cal{F}$ at pairwise distance at least $\delta$ from each other.

4 Our Approach and Techniques

The unifying idea of this work is an analysis of the structure of the Fourier Transform (FT) of $k$-SIIRVs. The FT is a natural tool to consider in this context. Recall that the FT of a sum of independent random variables is the product of the FT’s of the individual variables. Moreover, if two random variables have similar FT’s, they also have similar distributions. These two basic facts are the starting point of our analysis. We now provide an overview of the ideas underlying our results, and give a comparison to previous techniques.

Let $\mathcal{D}$ be a family of distributions over a domain of size $N.$ How many samples are required to learn an arbitrary $\mathbf{P}\in\mathcal{D}$ within variation distance $\epsilon$? Without any restrictions on $\mathcal{D},$ it is a folklore fact that the sample complexity learning is $\Theta(N/\epsilon^{2}).$ The optimal learning algorithm in this case is the obvious one: output the empirical distribution. By exploiting the structure of the family $\mathcal{D},$ one may obtain better results.

A very natural type of structure to consider is some sort of “shape constraint” on the probability density function, such as log-concavity or unimodality. There is a long line of work in statistics on this topic (see, e.g., the books [BBBB72, GJ14]), and more recently in TCS [DDS12a, CDSS14a, CDSS14b, ADLS15]. Alas, it turns out that $k$-SIIRVs do not satisfy any of the shape constraints considered in the literature (see [DDO+13] for a discussion).

A different type of structure, based on the notion of metric entropy [Yat85, Bir86, DL01], yields the following implication: If a distribution class $\mathcal{D}$ has an $\epsilon/2$-cover of size $M,$ then it is learnable with $O(\log M/\epsilon^{2})$ samples.We remark that the running time of this method is $\Omega(M/\epsilon^{2}),$ which is not necessarily polynomial in the sample size. In a celebrated paper in information theory [YB99], Yang and Barron show that, for broad families of (continuous) distributions, the metric entropy characterizes the sample complexity of learning. For $k$-SIIRVs, however, this is not the case: Via Theorem 1.4, the metric entropy method implies a sample upper bound of $O((1/\epsilon^{2})\cdot\log n+(k/\epsilon^{2})\cdot\log^{2}(1/\epsilon)).$ Note that, since our cover size upper bound is tight, this sample bound is the limit of the metric entropy method for $k$-SIIRVs. Thus, this method gives a suboptimal sample upper bound for our learning problem, both qualitatively (dependence on $n$), and quantitatively (dependence on $\epsilon$).

Previous work on learning $k$-SIIRVs [DDS12b, DDO+13] relies on a certain “regularity” lemma about the structure of these distributions: Any $k$-SIIRV is either $\epsilon$-close in total variation distance to being $L=\Theta(k^{9}/\epsilon^{4})$- “sparse”, i.e., it is supported on a set of at most $L$ consecutive integers, or $\epsilon$-close to being “Gaussian like”. In the former case, the distribution can be learned using $O(L/\epsilon^{2})$ samples, and in the latter case one can exploit the Gaussian structure to learn with a small number of samples as well. Unfortunately, the sparse case is a bottleneck for this approach, as any algorithm to learn a istribution over support $L$ requires $\Omega(L/\epsilon^{2})$ samples. Hence, one needs to exploit the structure of $k$-SIIRVs beyond the aforementioned.

In this paper, we depart from the aforementioned approaches. We identify a simple condition – the approximate sparsity of the Fourier transform – as the “right” property that determines the sample complexity of our learning problem. The Fourier sparsity explains why the sample complexity of learning $k$-SIIRVs is independent of $n,$ and allows us to obtain the sharp sample bound as a function of both $k$ and $\epsilon.$ We show that this is a more general phenomenon (see Theorem 2.5 in Section 2.2): any univariate distribution that has an $s$-sparse Fourier transform, in a certain well-defined technical sense, is learnable with $\widetilde{O}(s/\epsilon^{2})$ samples.

Our computationally efficient learning algorithm proceeds as follows: It starts by drawing an initial set of samples to determine the effective support of the target distribution and its Fourier transform. This is achieved by estimating the mean and variance of our SIIRV. We remark that, for computational purposes, our algorithm uses the Discrete Fourier Transform (DFT). For the appropriate definition of the DFT, we show (Lemma 2.3) there exists an explicit set $S$ of cardinality $|S|=O(k^{2}\log(k/\epsilon))$ that contains all the “heavy” Fourier coefficientsWe moreover show that there exists a set of cardinality $O(k\log(k/\epsilon))$ that contains all the “heavy” Fourier coefficients, alas this smaller set is not explicitly known a priori.. Our algorithm then draws an additional set of samples to estimate the DFT of the target distribution at the points of the effective support $S,$ and sets the DFT to everywhere else. By exploiting the sparsity in the Fourier domain, we show that the inverse of the empirical DFT achieves total variation distance $\epsilon/2$ after $\widetilde{O}(k/\epsilon^{2})$ samples. Note that an explicit description of an accurate hypothesis for our learning problem can have an effective support of size $\Omega(k\sqrt{n}).$ While we can easily obtain such a description (by explicitly computing the inverse DFT), this would not lead to a computationally efficient algorithm. We instead output a succinct description of our hypothesis (in time that is independent of $n$). In particular, our algorithm outputs the empirical DFT at the points of its effective support. Our learning algorithm is given in Section 2.1.

Finally, we note that our above-described Fourier-learning algorithm achieves a near-optimal sample complexity (up to logarithmic factors). The basic idea to obtain the optimal sample complexity is to smoothly mollify the DFT instead of truncating it. This removes some artifacts caused by a sharp truncation and yields a hypothesis whose error from the true distribution decays rapidly as we move away from the mean. Our sample-optimal upper bound is established in Section 2.4.

We start by commenting on previous approaches for proving cover upper bounds in this context. The main technique for the $2$-SIIRV cover upper bound of [DP09a] is the following lemma (that is deduced in [DP09a] using a result from [Roo00]): If two $2$-SIIRVs agree on their first $\Omega(\log(1/\epsilon))$ moments, then their total variation distance is at most $\epsilon$. First, we show that this moment-matching lemma is quantitatively tight: we give an example of two $2$-SIIRVs over $k+1$ variables that agree on the first $k$ moments and have variation distance $2^{-\Omega(k)}$ (Proposition B.1).

We emphasize however that such a moment-matching technique cannot be generalized to $k$-SIIRVs, even for $k=3.$ Intuitively, this is because knowledge about moments fails to account for potential periodic structure of the probability mass function that comes into play for $k>2$. For example, $\Omega(n)$ moments do not suffice to distinguish between the cases that a $3$-SIIRV of order $n$ is supported on the even versus the odd integers. More specifically, in Proposition B.2 (Appendix B), we give an explicit example of two $3$-SIIRVs of order $n/2$ that agree exactly on the first $n-1$ moments and have disjoint supports.

In conclusion, moment-based approaches fail to detect periodic structure. On the other hand, this type of structure is easily detectable by considering the Fourier transform. Our cover upper bound hinges on showing that the Fourier transform of a $k$-SIIRV is necessarily of low complexity, i.e., it can be succinctly described up to small error. In particular, since the Fourier transform is smooth, we show (Lemma 3.6), roughly, that its logarithm can be well approximated by a low degree Taylor polynomial on intervals of length $O(1/k).$ (Our actual statement is somewhat more complicated as it needs to account for roots of the Fourier transform close to the unit circle.) Therefore, providing approximations to the low-degree Taylor coefficients of the logarithm of the Fourier transform provides a concise approximate description of the distribution.

Our lower bounds take a geometric view of the problem. At a high-level, we consider the function that maps the set of $n(k-1)$ parameters defining a $k$-SIIRV to the corresponding probability mass function. We show that there exists a region of the space of distributions where this function is locally invertible. For $k=2$, we in fact show that the distribution of any $2$-SIIRV with distinct parameters lies in the interior of this region. This structural understanding allows us to use certain appropriately defined expectations to extract the effect of individual parameters on the distribution. In addition, for $n=\Theta(\log(1/\epsilon)),$ we show that near a particular $k$-SIIRV not only is the map from parameters to distribution locally a bijection, but that this map is actually surjective onto a ball of reasonable size. In other words, near this particular distribution, the $\Omega(k\log(1/\epsilon))$ parameters of the output distribution are effectively independent, which intuitively implies the $(1/\epsilon)^{\Omega(k\log(1/\epsilon))}$ lower bound on the cover size.

To prove our sample lower bound, at a high-level, we combine the aforementioned geometric understanding with Assouad’s lemma [Ass83]. We note that one might naively expect that such a situation would lead to a lower bound of $\Omega(k\log(1/\epsilon)/\epsilon^{2})$, but since the distributions under consideration have additional structure, it turns out that the best lower bound that can be obtained is $\Omega(k\sqrt{\log(1/\epsilon)}/\epsilon^{2}).$

5 Related Work

Density estimation is a classical topic in statistics and machine learning with a rich history and extensive literature (see e.g., [BBBB72, DG85, Sil86, Sco92, DL01]). The reader is referred to [Ize91] for a survey of statistical techniques in this context. In recent years, a large body of work in TCS has been studying these questions from a computational perspective; see e.g., [KMR+94, FM99, AK01, CGG02, VW02, FOS05, BS10, KMV10, MV10, DDS12a, DDS12b, DDO+13, CDSS14a, CDSS14b, ADLS15].

Covering numbers (and their logarithms, known as metric entropy numbers) were first defined by A. N. Kolmogorov in the 1950’s and have since played a central role in a number of areas, including approximation theory, geometric functional analysis (see, e.g., [Dud74, Mak86, BGL07] and the books [KT59, Lor66, CS90, ET96]), geometric approximation algorithms [Hp11], information theory, statistics, and machine learning (see, e.g., [Yat85, Bir86, HI90, HO97, YB99, GS13] and the books [vdVW96, DL01, Tsy08]).

Concurrent work by Daskalakis et al. [DKT15], using different techniques, gives upper bounds on the learning sample complexity of Poisson Multinomial Distributions (PMDs). While upper bounds on the sample complexity of PMDs yield similar upper bounds for $k$-SIIRVs, the implied upper bounds for $k$-SIIRVs are quantitatively significantly weaker than ours. Moreover, the [DKT15] learning algorithm has running time exponential in $k$ and super-polynomial in $1/\epsilon.$

6 Organization

In Section 2 we describe and analyze our learning algorithms for $k$-SIIRVs. Section 3 contains our cover upper bound construction. Our cover lower bound is given in Section 4, and our sample lower bound in Section 5.

Learning SIIRVs

In this section, we describe our algorithms for learning $k$-SIIRVs. The structure of this section is as follows: In Section 2.1, we give our sample near-optimal and computationally efficient learning algorithm. As mentioned in the introduction, our algorithm outputs a succinct description of its hypothesis $\mathbf{H}$, via its DFT. In Section 2.2, we provide a simple general algorithm that learns any one-dimensional discrete distribution with a sparse Fourier support. In Section 2.3, we show how to efficiently obtain an $\epsilon$-sampler for our unknown $k$-SIIRV, using the DFT representation of $\mathbf{H}$ as a black-box. Finally, in Section 2.4 we present our more sophisticated Fourier-based learning algorithm with optimal sample complexity.

The main result of this subsection is Theorem 1.1, which we state below in more detail for the sake of completeness.

For computational purposes, our learning algorithm in this section uses the Discrete Fourier Transform, which we now define.

This obstacle can be circumvented by relying on a new structural result that we believe may be of independent interest. We show that for any $k$-SIIRV with large variance, its Fourier Transform will have small effective support. In particular, for any $k$-SIIRV with standard deviation $s$ and $\epsilon>0$ we consider its Discrete Fourier transform modulo $M$, and show the set of points in $[M-1]$ whose Fourier transform is bigger than $\epsilon$ in magnitude has size at most $O(Mks^{-1}\sqrt{\log(1/\epsilon)})$. By choosing $M$ to be approximately $s\sqrt{\log(1/\epsilon)}$, i.e., of the same order as the effective support of $\mathbf{P}$, we conclude that the effective support of $\widehat{\mathbf{P}}$ (modulo $M$) is $O(k\log(1/\epsilon))$.

If the effective support for $\widehat{\mathbf{P}}$ was explicitly known, we could truncate our empirical Discrete Fourier transform $\widehat{\mathbf{Q}}$ (modulo $M$) outside this set and reduce the $L_{2}$ error $\|\widehat{\mathbf{Q}}-\widehat{\mathbf{P}}\|_{2}$ to $N^{-1/2}k^{1/2}s^{-1/2}\log^{1/4}(1/\epsilon)$. This in turn would correspond to an $L_{1}$ error of $O(N^{-1/2}k^{1/2}\sqrt{\log(1/\epsilon)})$. Unfortunately, we do not know exactly where the support of the Fourier transform is, so we will need to approximate it by calculating the empirical DFT where the support might be, and then simply truncating this empirical DFT whenever it is sufficiently small. Fortunately, we do have some idea of where the support is and it is not hard to show that we can truncate at all of the appropriate points with high probability.

The bulk of our analysis will depend on showing that the Fourier transform of $\mathbf{P}$ has appropriately small effective support. To do this we need the following lemma:

At most $4Mks^{-1}\sqrt{\log(1/\delta)}$ many integers $0\leq\xi\leq M-1$ have $|\widehat{\mathbf{P}}(\xi)|>\delta\;.$

Before we proceed with the proof of the lemma some comments are in order. Statement (i) of the lemma exhibits an explicit set $\mathcal{L}$ of cardinality $O(Mk^{2}s^{-1}\sqrt{\log(1/\delta)})$ that contains all the points $\xi\in[M-1]$ such that $|\widehat{\mathbf{P}}(\xi)|>\delta.$ Note that the set $\mathcal{L}$ can be efficiently computed from $M$, $\delta$, $s$, and does not otherwise depend on the particular $k$-SIIRV $\mathbf{P}$. Statement (ii) of the lemma shows that the effective support $\mathcal{L}^{\prime}=\mathcal{L}^{\prime}(\delta)=\{\xi\in[M-1]\mid|\widehat{\mathbf{P}}(\xi)|>\delta\}$ is in fact significantly smaller than $\mathcal{L}$, namely $|\mathcal{L}^{\prime}|=O(Mks^{-1}\sqrt{\log(1/\delta)})$. This part of the lemma is non-constructive in the sense that it does not provide an explicit description for $\mathcal{L}^{\prime}$ (beyond the fact that $\mathcal{L}^{\prime}\subseteq\mathcal{L}$). The upper bound on the size of the effective support is the basis for the analysis of our algorithm.

Since $\mathbf{P}\in\mathcal{S}_{n,k}$, for $X\sim\mathbf{P}$, we have $X=\sum_{i=1}^{n}X_{i}$ where each $X_{i}\sim\mathbf{P}_{i}$ for a $k$-IRV $\mathbf{P}_{i}$. Let $Y_{i}=X_{i}-X^{\prime}_{i}$ be the difference of two independent copies of $X_{i}$. Let $p_{ij}=\Pr\left[|Y_{i}|=j\right].$ Note that $Y_{i}$ is a symmetric random variable. Consider its DFT modulo $M$ which we will write as $\widehat{Y_{i}}$. We have the following sequence of (in)equalities:

Therefore, we have that $|\widehat{\mathbf{P}}(\xi)|^{2}=\prod_{i=1}^{n}|\widehat{\mathbf{P}_{i}}(\xi)|^{2}\leq\exp(-8\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}[\xi j/M]^{2})$. Taking square roots, we obtain

Note that we can relate the variance of $\mathbf{P}$ to the $p_{ij}$’s as follows:

To complete the proof of (i), we will need a simple counting argument given in the following claim:

For each $c$ satisfying $1\leq c\leq j-1$ there are either $\lfloor 2Ma\rfloor$ or $\lfloor 2Ma\rfloor+1$ integers $0\leq\xi\leq M-1$ with $|\frac{\xi}{M}-\frac{c}{j}|<a$. For $c=0$ and $c=j$ there are either $\lfloor Ma\rfloor$ or $\lfloor Ma\rfloor+1$ integers with $|\frac{\xi}{M}-\frac{c}{j}|<a$. Finally, note that $|\frac{\xi}{M}-\frac{c}{j}|<a$ for some $1\leq c\leq j-1$ if and only if $[j\xi/M]<aj$. ∎

An application of the above claim for $a=(1/2s)\sqrt{\ln(1/\delta)}$ implies that there are at most

integers $0\leq\xi\leq M-1$ with \min_{j}\big{(}\frac{[\xi j/M]}{j}\big{)}^{2}<\ln(1/\delta)/(4s^{2}). For all other integers we have $|\widehat{\mathbf{P}}(\xi)|\leq\delta\;,$ which completes the proof of (i).

Note that $[\xi j/M]\geq\sqrt{1-N_{j}}\cdot ks^{-1}\sqrt{\ln(1/\delta)}/2$. Plugging this into (1) gives

Since $s^{2}=\frac{1}{2}\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}j^{2}\leq\frac{k^{2}}{2}\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}$, it follows that $\theta:=\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}\geq 2s^{2}/k^{2}.$ Therefore,

By Markov’s inequality, except with probability $4ks^{-1}\sqrt{\ln(1/\delta)}$, we have that $\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}N_{j}\leq\frac{\theta}{2}.$ In this event, we have $\sum_{i=1}^{n}\sum_{j=1}^{k-1}p_{ij}(1-N_{j})\geq\frac{\theta}{2}$ and hence

Since $\xi$ is uniformly distributed on $[M-1]$, it follows that $|\widehat{\mathbf{P}}(\xi)|>\delta$ for at most $4Mks^{-1}\sqrt{\ln(1/\delta)}$ integers $\xi$ in $[M-1]$. This completes the proof of (ii). ∎

Note that it is straightforward to verify the sample complexity bound. The running time of the algorithm is dominated by computing the DFT $\widehat{\mathbf{Q}}$. Since the support of $\mathbf{Q}$ is at most $N$, for each $\xi\in S$, we sum at most $N$ terms to calculate $\widehat{\mathbf{Q}}(\xi)$. Therefore, the overall running time is $O(N\cdot|S|)=O(k\log^{2}(k/\epsilon)/\epsilon^{2}\cdot k^{2}\log(k/\epsilon))=O(k^{3}\log^{3}(k/\epsilon)/\epsilon^{2})$ as claimed.

To show correctness, we will prove that the expected squared $L_{2}$ norm between $\widehat{\mathbf{H}}$ and $\widehat{\mathbf{P}}$ is small, i.e., that $\|\widehat{\mathbf{H}}-\widehat{\mathbf{P}}\|_{2}^{2}=(1/M)\cdot\sum_{\xi=0}^{M-1}|\widehat{\mathbf{H}}(\xi)-\widehat{\mathbf{P}}(\xi)|^{2}$ has small expected value.

It is easy to see that, after drawing a constant number of samples, the quantities $\widetilde{\mu}$ and $\widetilde{\sigma}$ can be estimated to satisfy the required conditions with probability at least $19/20$. (This follows for example by Lemma 6 of [DDS12b] with $\epsilon=1/2.$) We will henceforth condition on this event.

Since $M=1+2\lceil 6\widetilde{\sigma}\sqrt{\ln(4/\epsilon)})\rceil$, a random variable $X\sim\mathbf{P}$ lies in $[\lfloor\widetilde{\mu}\rfloor-(M-1)/2,\lfloor\widetilde{\mu}\rfloor+(M-1)/2]$ with probability at least $1-\frac{\epsilon}{2}$. Indeed, an application of Bernstein’s inequality for $X$ yields that

where $\mu$ is the mean of $\mathbf{P}$, for any $t>0$. For $t=2s\sqrt{\ln(4/\epsilon)}$, we have $t^{2}=(\ln(4/\epsilon))4s^{2}$ and $2s^{2}+\frac{2}{3}kt=2s^{2}+\frac{4}{3}ks\sqrt{\ln(4/\epsilon)}\leq\frac{8}{3}s^{2}\leq 4s^{2}$. Thus, $\Pr(X>\mu+t)\leq\epsilon/4$. Similarly, it holds $\Pr(X<\mu-t)\leq\epsilon/4.$ Now note that $\lfloor\widetilde{\mu}\rfloor+(M-1)/2\geq(\mu-s)+\lceil 3s\sqrt{\ln(4/\epsilon)})\rceil\geq\mu+t$ and $\lfloor\widetilde{\mu}\rfloor-(M-1)/2\leq\mu-t$. Hence, $X$ is in $[\lfloor\widetilde{\mu}\rfloor-(M-1)/2,\lfloor\widetilde{\mu}\rfloor+(M-1)/2]$ with probability at least $1-\epsilon/2$ as desired.

Fix $T=R/2=C^{-1}\epsilon/(\sqrt{k\ln(k/\epsilon))}$. We analyze separately the contribution to the squared $L_{2}$ norm coming from $\xi$ with $|\widehat{\mathbf{P}}(\xi)|>T$ and with $|\widehat{\mathbf{P}}(\xi)|\leq T$. Let us denote $\mathcal{L^{\prime}}(T)=\{\xi\in[M-1]\mid|\widehat{\mathbf{P}}(\xi)|>T\}$. First consider

We first claim that with high probability $\widehat{\mathbf{H}}(\xi)=0$ for all $\xi\in\overline{\mathcal{L^{\prime}}}(T)$. This happens automatically when $\xi\not\in S$, where the $S$ is defined in the algorithm description. Note that $|S|=O(k^{2}\log(k/\epsilon))$. For $\xi\in S\setminus\mathcal{L^{\prime}}(T)$, we note that $\widehat{\mathbf{Q}}(\xi)$ is an average of $N$ i.i.d. numbers each of absolute value $1$ and mean $\widehat{\mathbf{P}}(\xi)$ (which has absolute value less than $1$). Note that if $|\widehat{\mathbf{Q}}(\xi)-\widehat{\mathbf{P}}(\xi)|\geq R-T$, then either the real or the imaginary part is at least $(R-T)/\sqrt{2}$. By a Chernoff bound, the probability that for a given $\xi\in S\setminus\mathcal{L^{\prime}}(T)$, $\Re(\widehat{\mathbf{Q}}(\xi)-\widehat{\mathbf{P}}(\xi))\geq(R-T)/\sqrt{2}$ is at most $2\exp(-N(R-T)^{2}/4)$. The same is true of the imaginary part so by a union bound the probability that $|\widehat{\mathbf{Q}}(\xi)-\widehat{\mathbf{P}}(\xi)|\geq R-T$ is at most $4\exp(-N(R-T)^{2}/4)$. Again by a union bound we get that the probability that any $\xi\in S\setminus\mathcal{L^{\prime}}(T)$ has $|\widehat{\mathbf{Q}}(\xi)-\widehat{\mathbf{P}}(\xi)|\geq R-T$ is at most $O(k^{2}\log(k/\epsilon)\exp(-N(R-T)^{2}/4))=O(k^{2}\log(k/\epsilon)\exp(-C\ln(k/\epsilon)))=O(\epsilon^{C-1})$. Hence, except with probability $O(\epsilon^{C-1})$, for all $\xi$ in $S\setminus\mathcal{L^{\prime}}(T)$ we have $|\widehat{\mathbf{Q}}(\xi)-\widehat{\mathbf{P}}(\xi)|<R-T$ and so $|\widehat{\mathbf{Q}}(\xi)|\leq R$. In fact, the total expected contribution to the squared $L_{2}$ norm coming from cases when $\widehat{\mathbf{H}}(\xi)$ is not identically on all such $\xi$ is also $O(\epsilon^{C-1})$. Therefore, up to negligible error, the squared $L_{2}$ error coming from this range is at most

Applying Lemma 2.3 (ii) with $\delta:=T2^{-r-1}$ for each $r\geq 0$, this is at most

We now consider the remaining contribution

By Lemma 2.3 (ii) applied with $\delta:=T,$ we have $|\mathcal{L^{\prime}}(T)|\leq 4ks^{-1}\sqrt{\ln(1/T)}$. So, the expected size of the $L_{2}^{2}$ error on $\mathcal{L^{\prime}}(T)$ has

Combining the above results, we find that the expected $L_{2}^{2}$ error between $\widehat{\mathbf{H}}$ and $\widehat{\mathbf{P}}$ is at most

Therefore, if $C$ is sufficiently large, Markov’s inequality yields that, with probability $\frac{2}{3}$, we have $\|\widehat{\mathbf{H}}-\widehat{\mathbf{P}}\|_{2}^{2}<\epsilon^{2}/M.$

Since $X\sim\mathbf{P}$ is in $[\lfloor\widetilde{\mu}\rfloor-(M-1)/2,\lfloor\widetilde{\mu}\rfloor+(M-1)/2]$ with probability at least $1-\epsilon/2$, we have $\|\mathbf{P}-\mathbf{P}^{\prime}\|_{1}\leq\epsilon$ and so $\|\mathbf{P}-\mathbf{H}\|_{1}\leq\|\mathbf{P}-\mathbf{P}^{\prime}\|_{1}+\|\mathbf{H}-\mathbf{P}^{\prime}\|_{1}\leq 2\epsilon$. Since $\widehat{\mathbf{H}}(0)=\widehat{\mathbf{Q}}(0)=1$, it follows that $\sum_{i=0}^{n}\mathbf{H}(i)=1$. Also, by symmetry, all the $\mathbf{H}(i)$’s are real. This completes the proof of Theorem 2.1. ∎

2 A General Fourier Learning Algorithm

The algorithmic approach of the previous subsection is not specialized to $k$-SIIRVs, but is applicable more generally. In essence, the approach really only depended upon two facts:

$\mathbf{P}$ is effectively supported on a small set $T.$

${\widehat{\mathbf{P}}}$ is effectively supported on a small set $S.$

It turns out that by using similar ideas, we can learn any probability distribution with these properties. The following simple theorem provides a generalization for integer-valued random variables. However, the approach can also be generalized to higher dimensions and to continuous distributions.

Then, there exists an algorithm which learns $\mathbf{P}$ to total variational distance $\epsilon$ using $N=O(|T|\mu(S)/\epsilon^{2})$ samples, where $\mu(S)$ is the Lebesgue measure of $S.$

Algorithm Learn-Sparse-FT Input: sample access to a distribution $\mathbf{P}$ over $[n]$ and $\epsilon>0.$ Let $C$ be a sufficiently large universal constant. 1. Take $N=C|T|\mu(S)/\epsilon^{2}$ samples from $\mathbf{P}$ to get an empirical distribution $\mathbf{Q}$. 2. Compute $\widehat{\mathbf{H}}$ which is defined as $\widehat{\mathbf{H}}(\xi)=\widehat{\mathbf{Q}}(\xi),$ if $\xi\in S,$ and $\widehat{\mathbf{H}}(\xi)=0$ otherwise. 3. Output $\mathbf{H},$ where $\mathbf{H}$ is the inverse Fourier transform on ${\widehat{\mathbf{H}}}$ restricted to $T.$ In particular $\mathbf{H}(i)=\int_{\xi\in S}e(-n\xi){\widehat{\mathbf{H}}}(\xi)d\xi$ for $i\in T$ and for $i\not\in T.$

Note that this is exactly the form of the algorithm for learning $k$-SIIRVs, except that the latter algorithm must also learn $T$ (which is done by computing an approximate mean and variance) and $S$ (which is obtained through a thresholding procedure). Also, note that we use the continuous Fourier transform here rather than a discrete Fourier transform. This is mostly for conceptual convenience. In practice, the continuous Fourier transform can be replaced by a sufficiently fine discrete Fourier transform, yielding an algorithm in which the integrals can be replaced by finite sums.

The analysis of the algorithm is not difficult. We begin by bounding that the expected $L_{2}$ difference between ${\widehat{\mathbf{P}}}$ and ${\widehat{\mathbf{H}}}.$ In particular, we note that

Now, for any given value of $\xi$, we note that ${\widehat{\mathbf{P}}}(\xi)-{\widehat{\mathbf{H}}}(\xi)$ has mean and variance at most $1/N$. Therefore, we have that

For $C$ large enough, by the Markov inequality, this is at most $\epsilon^{2}/(9|T|)$ with probability at least $2/3.$ If this holds, then

By Plancherel’s Theorem, this would imply that the squared $L^{2}$ distance between $\mathbf{P}$ and the inverse Fourier transform of $\mathbf{H}$ is at most $\epsilon^{2}/(4|T|).$ Along with Cauchy-Schwartz, this implies that

3 An Efficient Sampler for our Hypothesis

The learning algorithm of Section 2.1 outputs a succinct description of the hypothesis pseudo-distribution $\mathbf{H}$, via its DFT. This immediately provides us with an efficient evaluation oracle for $\mathbf{H},$ i.e., an $\epsilon$-evaluation oracle for our target SIIRV $\mathbf{P}.$ The running time of this oracle is linear in the size of $S,$ the effective support of the DFT.

Note that we can explicitly output the hypothesis $\mathbf{H}$ by computing the inverse DFT at all the points of the support of $\mathbf{H}.$ However, in contrast to the effective support of ${\widehat{\mathbf{H}}},$ the support of $\mathbf{H}$ can be large, and this explicit description would not lead to a computationally efficient algorithm. In this subsection, we show how to efficiently obtain an $\epsilon$-sampler for our unknown $k$-SIIRV $\mathbf{P},$ using the DFT representation of $\mathbf{H}$ as a black-box. In particular, starting with the DFT of an accurate hypothesis $\mathbf{H},$ represented via its DFT, we show how to efficiently obtain an $\epsilon$-sampler for the unknown target distribution. We remark that the efficient procedure of this section is not restricted to $k$-SIIRVs, but is more general, applying to all univariate discrete distributions for which an efficient oracle for the DFT is available.

In particular, we prove the following theorem:

Combining the above with Theorem 2.1, we get:

For the output of algorithm Learn-SIIRV, $M=O((1+\sigma)\sqrt{\log(1/\epsilon)})=O(kn)$ and $|S|\leq|\mathcal{L^{\prime}}(T)|\leq 2Mks^{-1}\sqrt{\ln(1/T)}=O(k\log(k/\epsilon))$. ∎

This section is devoted to the proof of Theorem 2.6. We start by providing some high-level intuition. Roughly speaking, we obtain the desired sampler by the Cumulative Distribution Function (CDF) corresponding to $\mathbf{H}.$ We use the DFT to obtain a closed form expression for the CDF of $\mathbf{H},$ and then we query the CDF using an appropriate binary search procedure to sample from the distribution. One subtle point is that $\mathbf{H}(x)$ is a pseudo-distribution, i.e. it is not necessarily non-negative at all points. Our analysis shows that this does not pose any problems with correctness.

Our first lemma shows that it is sufficient to have an efficient oracle for the CDF:

We begin our analysis by producing an algorithm that works when we are able to exactly compute $c_{\mathbf{H}}(x).$

We have an interval $[a^{\prime},b^{\prime}]$, initially $[a-1,b]$, with $c_{\mathbf{H}}(a^{\prime})\leq y\leq c_{\mathbf{H}}(b^{\prime})$ and $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(b^{\prime}).$

If $b^{\prime}-a^{\prime}=1$, output $d_{\mathbf{H}}(y)=b^{\prime}.$

Otherwise, find the midpoint $c^{\prime}=\lfloor(a^{\prime}+b^{\prime})/2\rfloor.$

If $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(c^{\prime})$ and $y\leq c_{\mathbf{H}}(c^{\prime})$, repeat with $[a^{\prime},c^{\prime}]$; else repeat with $[c^{\prime},b]$.

The function $d_{\mathbf{H}}$ satisfies: For any $y\in$, it holds $c_{\mathbf{H}}(d_{\mathbf{H}}(y)-1)\leq y\leq c_{\mathbf{H}}(d_{\mathbf{H}}(y))$ and $c_{\mathbf{H}}(d_{\mathbf{H}}(y)-1)<c_{\mathbf{H}}(d_{\mathbf{H}}(y)).$

Note that if we don’t have $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(c^{\prime})$ and $y\leq c_{\mathbf{H}}(c^{\prime})$, then $c_{\mathbf{H}}(c^{\prime})<y\leq c_{\mathbf{H}}(b^{\prime})$. So, Step 4 gives an interval $[a^{\prime},b^{\prime}]$ which satisfies $c_{\mathbf{H}}(a^{\prime})\leq y\leq c_{\mathbf{H}}(b^{\prime})$ and $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(b^{\prime})$. The initial interval $[a-1,b]$ satisfies these conditions since $c_{\mathbf{H}}(a-1)=0$ and $c_{\mathbf{H}}(b)=1$. By induction, all $[a^{\prime},b^{\prime}]$ in the execution of the above algorithm have $c_{\mathbf{H}}(a^{\prime})\leq y\leq c_{\mathbf{H}}(b^{\prime})$ and $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(b^{\prime})$. Since this is impossible if $a^{\prime}=b^{\prime}$, and Step 4 always recurses on a shorter interval, we eventually have $b^{\prime}-a^{\prime}=1$. Then, the conditions $c_{\mathbf{H}}(a^{\prime})\leq y\leq c_{\mathbf{H}}(b^{\prime})$ and $c_{\mathbf{H}}(a^{\prime})<c_{\mathbf{H}}(b^{\prime})$ give the claim. ∎

Computing $d_{\mathbf{H}}(y)$ requires $O(\log(b-a+1))$ evaluations of $c_{\mathbf{H}}$, and $O(\log(b-a+1))$ comparisons of $y.$ For the rest of this proof, we will use $n=b-a+1$ to denote the support size.

We now show how to effectively sample from $\mathbf{Q}^{\prime}$. The issue is how to simulate a sample from the uniform distribution on $$with uniform random bits. We do this by flipping coins for the bits of$Y$lazily. We note that we will only need to know more than$m$bits of$Y$if$Y$is within$2^{-m}$of one of the values of$c_{\mathbf{H}}(x)$for some$x$. By a union bound, this happens with probability at most$n2^{-m}$over the choice of$Y$. Therefore, for$m>\log_{2}(10n/\epsilon)$, the probability that this will happen is at most$\epsilon/10$ and can be ignored.

We next show that we can efficiently compute an appropriate CDF, using the DFT.

Recall that the PMF of $\mathbf{H}$ at $x\in S$ is given by the inverse DFT:

When $\xi\not=0$, the term $\sum_{i:a\leq i\leq x}e(-\xi x/M)$ is a geometric series. By standard results on its sum, we have:

When $\xi=0$, $e(-\xi)=1$, and we get $\sum_{a\leq i\leq x}e(-\xi x/M)=i+1-a.$ In this case, we also have $\widehat{\mathbf{H}}(\xi)=1.$ Putting this together we have:

Hence, we obtain a closed form expression for the CDF that can be approximated to desired precision in time $O(|S|\log(1/\delta)).$ ∎

Now we can prove the main theorem of this subsection.

By Proposition 2.10, we can efficiently calculate the CDF of $\mathbf{H}$. So, we can apply Lemma 2.8 to this CDF. This gives us an $\epsilon$-sampler for $\mathbf{H}.$ To find the time it takes to compute each sample, we need to substitute $D=O\left(|S|\log(M/\epsilon)\right)$ from the running time of the CDF into the bound in Lemma 2.8, yielding $O\left(\log M\cdot\log(M/\epsilon)\right)\cdot|S|$ time. This completes the proof. ∎

4 Sample–Optimal Learning Algorithm

In this subsection, we show how to improve the sample complexity of our learning algorithm for $k$-SIIRVs given in Section 2.1, and obtain an algorithm with optimal sample complexity (up to constant factors). The basic idea behind the improvement is as follows: In our previous analysis, we made critical use of the fact that essentially all of the mass of the distribution in question lies in an explicit interval of length $O(s\sqrt{\log(1/\epsilon)}),$ where $s$ is the standard deviation. By using our Fourier learning approach, we were able to learn a distribution that approximated our target on this support. In order to improve this algorithm, we observe that although it is necessary to move $\Omega(\sqrt{\log(1/\epsilon)})$ standard deviations from the mean before the cumulative density function (CDF) drops below $\epsilon,$ the CDF has already begun to drop off exponentially after only a single standard deviation from the mean.

We will warm up in Section 2.4.1, where we describe our algorithm in the case of $2$-SIIRVs. This will exhibit the important new ideas of this technique. Then, in Section 2.4.2, we extend these results to $k$-SIIRVs, which brings with it several technical complications, mostly arising from the fact that we do not know a priori a good effective support for the Fourier transform.

In this subsection, we will prove the following theorem:

There exists an algorithm that given $N=O(\sqrt{\log(1/\epsilon)}/\epsilon^{2})$ independent samples to a $2$-SIIRV $X,$ runs in time $O(N)$ and with probability at least $2/3$ outputs a hypothesis distribution $Y$ that is within $\epsilon$ of $X$ in total variational distance.

Our new algorithm Learn-$2$-SIIRV-Optimal-Sample is described in pseudocode below. We first provide an equivalent alternative interpretation of our algorithm in terms of truncating the Fourier transform. As in our algorithm Learn-SIIRV of Section 2.1, we start by obtain approximations $\widetilde{\sigma}^{2}$ and $\widetilde{\mu}$ for the variance and mean. Similarly, we output the empirical distribution if $\widetilde{\sigma}\leq\Theta(\sqrt{\ln(1/\epsilon)}).$ This allows us to assume that $\widetilde{\sigma}=\Omega(\sqrt{\ln(1/\epsilon)}).$ (Note that this bound is not as strong as that in Learn-SIIRV because in the current setting we aim to use fewer samples.)

Let $I(\xi)$ be the indicator function of the interval $[-C\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},C\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}],$ for $C$ a sufficiently large constant. Let ${\widehat{F}}$ be the convolution of $I$ and ${\widehat{G}},$ i.e., ${\widehat{F}}=I\ast{\widehat{G}}.$ We note that multiplication by ${\widehat{F}}$ is an appropriate method of thresholding. In particular, we start by showing that ${\widehat{F}}$ approximates $I$ in the following way:

${\widehat{F}}(\xi)\geq 1-\epsilon^{2}$ for $|\xi|\leq(C-3)\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}$.

${\widehat{F}}(\xi)\leq\epsilon^{2}$ for $\frac{1}{4}\geq|\xi|\geq(C+3)\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}$.

Note that ${\widehat{F}}$ is the convolution of $I$ and ${\widehat{G}}.$ We can write:

Clearly, this convolution is positive at all points. Thus, we get (i).

When $|\xi|\leq(C-3)\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},$ note that the integral in (5) is over

By standard tail bounds, the Gaussian $\frac{s}{\sqrt{2\pi}}e^{-\widetilde{\sigma}^{2}\nu^{2}/2}$ has all but $1-\epsilon^{2}$ of its mass in the interval $[-3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}],$ and so ${\widehat{F}}(\xi)\geq 1-\epsilon^{2}.$ This gives us (ii).

When $\frac{1}{4}\geq|\xi|\geq(C+3)\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},$ the integral in (5) is over $\nu\in[\xi-C\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},\xi+C\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}],$ which is disjoint from the interval $[-3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}].$ By the same bound, the Gaussian has at most $\epsilon^{2}$ of its mass outside $[-3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)},3\widetilde{\sigma}^{-1}\sqrt{\log(1/\epsilon)}].$ So, we deduce (iii). ∎

At a high-level, our new algorithm involves the following steps:

Let $Z$ be the empirical distribution and ${\widehat{Z}}$ be the (continuous) Fourier transform of $Z.$

Let ${\widehat{Y}}(\xi)={\widehat{Z}}(\xi){\widehat{F}}(\xi).$

Let $Y$ be the truncation of the inverse Fourier transform of ${\widehat{Y}}$ to $[\widetilde{\mu}-C\widetilde{\sigma}\sqrt{\log(1/\epsilon)},\widetilde{\mu}+C\widetilde{\sigma}\sqrt{\log(1/\epsilon)}],$ for $C$ a sufficiently large constant.

We now proceed to show correctness. In the proof of Theorem 2.1, we argued that $O(1)$ samples suffice to get that with high probability $\widetilde{\sigma}$ and $\widetilde{\mu}$ satisfy the desired bounds. We condition on this event. We claim that when $\widetilde{\sigma}$ is small, namely $O(\sqrt{\ln(1/\epsilon)}),$ the empirical distribution suffices. This follows from the fact that the empirical estimate of a discrete distribution $\mathbf{P}$ has expected variation distance $\leq\epsilon$ from $\mathbf{P}$ after $O(\|\mathbf{P}\|_{1/2}/\epsilon^{2})$ samples. By an application of Bernstein’s inequality (see Lemma C.3) it follows that a $2$-SIIRV with standard deviation $\sigma$ has $1/2$-norm bounded from above by $O(\sigma+1).$ This proves our claim.

So, we henceforth assume that $\widetilde{\sigma}$ is $\Omega(\sqrt{\ln(1/\epsilon)})$ and $O(1/\epsilon).$ We now proceed with the main part of the analysis. We start with the following simple claim:

Recall that $Y$ is a the convolution of $Z$ with $F.$ If we consider our samples to be random variables $X_{(1)},\ldots,X_{(N)}$ each of which is an i.i.d. copy of $X,$ we can express $Y(p)$ for a given $p$ as a random variable: $Y(p)=\frac{1}{N}\sum_{i=1}^{N}F(p-X_{(i)})\;,$ for $a\leq p\leq b,$ where $a=\widetilde{\mu}-C\widetilde{\sigma}\sqrt{\log(1/\epsilon)}$ and $b=\widetilde{\mu}+C\widetilde{\sigma}\sqrt{\log(1/\epsilon)}.$ Note that the expectation of $Y(p)$ is

We claim that this quantity will become small as $p$ moves away from $\mu.$ Intuitively, this should be the case because for $p$ far from $\mu,$ then for all $q$ either $|p-q|$ will be large or $|q-\mu|$ will be large. In the former case, $F(p-q)$ is small, and in the latter $X(q)$ is. In order to properly analyze this quantity, we will have to group up these errors for $p$ in blocks of size $\widetilde{\sigma}.$ In particular, we have that

for $N=\sqrt{\log(1/\epsilon)}/\epsilon^{2}.$ This completes the proof. ∎

4.2 Sample Optimal Learning Algorithm for k𝑘k-SIIRVs

The proof of this theorem is somewhat analogous to that of Theorem 2.11. However, it should be noted that the runtime of this algorithm is not given. This is because the runtime of the simplest such algorithm is actually exponential in $k.$ The difficulty is that while in the $2$-SIIRV case we could determine the effective support of the Fourier transform just from the standard deviation, in the case of $k$-SIIRVs this is not the case. In essence, our algorithm will first guess this effective support (of which there are exponentially many possibilities), and then given this guess will run an appropriate algorithm. At the end, we will need to run a standard tournament procedure (e.g., [DL01]) to determine which of these guesses lead to the closest approximation to $X.$ Since the number of possibilities is $2^{O(k)}$ (see Claim 2.19), the sample complexity of this tournament is $O(k/\epsilon^{2}).$

Once again, under these assumptions, we can use Bernstein’s inequality to prove concentration bounds for $X$:

Suppose that $\widetilde{\sigma}\geq Ck\sqrt{\log(1/\epsilon)}$, for $C$ sufficiently large, then for all $t\geq 0$ we have that

We assumed that $|\mu-\widetilde{\mu}|\leq\widetilde{\sigma}$. So, if $|X-\widetilde{\mu}|\geq(2+t)\widetilde{\sigma})$, then $|X-\mu|\geq(1+t)\widetilde{\sigma}$. Bernstein’s inequality gives that

Since $\widetilde{\sigma}=\Omega(k\sqrt{\log(1/\epsilon)})=\Omega(k)$ and $\mathop{\textnormal{Var}}\nolimits[X]=O(\widetilde{\sigma}^{2})$, we have that $\mathop{\textnormal{Var}}\nolimits[X]+\frac{1}{3}k=O(\widetilde{\sigma}^{2})$. ∎