The Fourier Transform of Poisson Multinomial Distributions and its Algorithmic Applications

Introduction

PMDs comprise a broad class of discrete distributions of fundamental importance in computer science, probability, and statistics. A large body of work in the probability and statistics literature has been devoted to the study of the behavior of PMDs under various structural conditions [Bar88, Loh92, BHJ92, Ben03, Roo99, Roo10]. PMDs generalize the familiar binomial and multinomial distributions, and describe many distributions commonly encountered in computer science (see, e.g., [DP07, DP08, Val08, VV11]). The $k=2$ case corresponds to the Poisson binomial distribution (PBD), introduced by Poisson [Poi37] as a non-trivial generalization of the binomial distribution.

Recent years have witnessed a flurry of research activity on PMDs and related distributions, from several perspectives of theoretical computer science, including learning [DDS12, DDO+13, DKS15a, DKT15, DKS15b], property testing [Val08, VV10, VV11], computational game theory [DP07, DP08, BCI+08, DP09, DP14, GT14], and derandomization [GMRZ11, BDS12, De15, GKM15]. More specifically, the following questions have been of interest to the TCS community:

Is there a statistically and computationally efficient algorithm for learning PMDs from independent samples in total variation distance?

How fast can we compute approximate Nash equilibria in anonymous games with many players and a small number of strategies per player?

How well can a PMD be approximated, in total variation distance, by a discretized Gaussian with the same mean and covariance matrix?

The first question is a fundamental problem in unsupervised learning that has received considerable recent attention in TCS [DDS12, DDO+13, DKS15a, DKT15, DKS15b]. The aforementioned works have studied the learnability of PMDs, and related distribution families, in particular PBDs (i.e., $(n,2)$-PMDs) and sums of independent integer random variables. Prior to this work, no computationally efficient learning algorithm for PMDs was known, even for the case of $k=3.$

The second question concerns an important class of succinct games previously studied in the economics literature [Mil96, Blo99, Blo05], whose (exact) Nash equilibrium computation was recently shown to be intractable [CDO15]. The formal connection between computing Nash equilibria in these games and PMDs was established in a sequence of papers by Daskalakis and Papadimitriou [DP07, DP08, DP09, DP14], who leveraged it to gave the first PTAS for the problem. Prior to this work, no efficient PTAS was known, even for anonymous games with $3$ strategies per player.

The third question refers to the design of Central Limit Theorems (CLTs) for PMDs with respect to the total variation distance. Despite substantial amount of work in probability theory, the first strong CLT of this form appears to have been shown by Valiant and Valiant [VV10, VV11], motivated by applications in distribution property testing. [VV10, VV11] leveraged their CLT to obtain tight lower bounds for several fundamental problems in property testing. We remark that the error bound of the [VV10] CLT has a logarithmic dependence on the size $n$ of the PMD (number of summands), and it was conjectured in [VV10] that this dependence is unnecessary.

2 Our Results

The main technical contribution of this work is the use of Fourier analytic techniques to obtain a refined understanding of the structure of PMDs. As our core structural result, we prove that the Fourier transform of PMDs is approximately sparse, i.e., roughly speaking, its $L_{1}$-norm is small outside a small set. By building on this property, we are able to obtain various new structural results about PMDs, and make progress on the three questions stated in the previous subsection. In this subsection, we describe our algorithmic and structural contributions in detail.

We start by stating our algorithmic results in learning and computational game theory, followed by an informal description of our structural results and the connections between them.

As our main learning result, we obtain the first statistically and computationally efficient learning algorithm for PMDs with respect to the total variation distance. In particular, we show:

We remark that our learning algorithm outputs a succinct description of its hypothesis $\mathbf{H},$ via its Discrete Fourier Transform (DFT), ${\widehat{\mathbf{H}}},$ which is supported on a small size set. We show that the DFT gives both an efficient $\epsilon$-sampler and an efficient $\epsilon$-evaluation oracle for $\mathbf{P}.$

Our learning algorithm and its analysis are described in Section 3.

As our second algorithmic contribution, we give the first efficient polynomial-time approximation scheme (EPTAS) for computing Nash equilibria in anonymous games with many players and a small number of strategies. In anonymous games, all players have the same set of strategies, and the payoff of a player depends on the strategy played by the player and the number of other players who play each of the strategies. In particular, we show:

There is an EPTAS for the mixed Nash equilibrium problem for normalized anonymous games with a constant number of strategies. More precisely, there exists an algorithm with the following performance guarantee: for all $\epsilon>0$, and any normalized anonymous game $\cal G$ of $n$ players and $k$ strategies, the algorithm runs in time $(kn)^{O(k^{3})}(1/\epsilon)^{O(k^{3}\log(k/\epsilon)/\log\log(k/\epsilon))^{k-1}},$ and outputs a (well-supported) $\epsilon$-Nash equilibrium of $\cal G.$

Similarly to [DP08, DP14], our algorithm proceeds by constructing a proper $\epsilon$-cover, in total variation distance, for the space of PMDs. A proper $\epsilon$-cover for $\mathcal{M}_{n,k},$ the set of all $(n,k)$-PMDs, is a subset $C$ of $\mathcal{M}_{n,k}$ such that any distribution in $\mathcal{M}_{n,k}$ is within total variation distance $\epsilon$ from some distribution in $C.$ Our main technical contribution is the efficient construction of a proper $\epsilon$-cover of near-minimum size (see Theorem 1.4). We note that, as follows from Theorem 1.5, the quasi-polynomial dependence on $1/\epsilon$ and the doubly exponential dependence on $k$ in the runtime are unavoidable for any cover-based algorithm. Our cover upper and lower bounds and our Nash approximation algorithm are given in Section 4.

Using our Fourier-based machinery, we prove a strong “size-free” CLT relating the total variation distance between a PMD and an appropriately discretized Gaussian with the same mean and covariance matrix. In particular, we show:

Let $X$ be an $(n,k)$-PMD with covariance matrix $\Sigma.$ Suppose that $\Sigma$ has no eigenvectors other than $\mathbf{1}=(1,1,\ldots,1)$ with eigenvalue less than $\sigma.$ Then, there exists a discrete Gaussian $G$ so that

We remark that our techniques for proving Theorem 1.3 are orthogonal to those of [VV10, VV11]. While Valiant and Valiant use Stein’s method, we prove our strengthened CLT using the Fourier techniques that underly this paper. We view Fourier analysis as the right technical tool to analyze sums of independent random variables. An additional ingredient that we require is the saddlepoint method from complex analysis. We hope that our new CLT will be of broader use as an analytic tool to the TCS community. Our CLT is proved in Section 5.

We now provide a brief intuitive overview of our new structural results for PMDs, the relation between them, and their connection to our algorithmic results mentioned above. The unifying theme of our work is a refined analysis of the structure of PMDs, based on their Fourier transform. The Fourier transform is one of the most natural technical tools to consider for analyzing sums of independent random variables, and indeed one of the classical proofs of the (asymptotic) central limit theorem is based on Fourier methods. The basis of our results, both algorithmic and structural, is the following statement:

Informal Lemma (Sparsity of the Fourier Transform of PMDs.) For any $(n,k)$-PMD $\mathbf{P}$, and any $\epsilon>0$ there exists a “small” set $T=T(\mathbf{P},\epsilon),$ such that the $L_{1}$-norm of its Fourier transform, ${\widehat{\mathbf{P}}},$ outside the set $T$ is at most $\epsilon.$

We will need two different versions of the above statement for our applications, and therefore we do not provide a formal statement at this stage. The precise meaning of the term “small” depends on the setting: For the continuous Fourier transform, we essentially prove that the product of the volume of the effective support of the Fourier transform times the number of points in the effective support of our distribution is small. In particular, the set $T$ is a scaled version of the dual ellipsoid to the ellipsoid defined by the covariance matrix of $\mathbf{P}.$ Hence, roughly speaking, ${\widehat{\mathbf{P}}}$ has an effective support that is the dual of the effective support of $\mathbf{P}.$ (See Lemma 4.2 in Section 4 for the precise statement.)

In the case of the Discrete Fourier Transform (DFT), we show that there exists a discrete set with small cardinality, such that $L_{1}$-norm of the DFT outside this set is small. At a high-level, to prove this statement, we need the appropriate definition of the (multidimensional) DFT, which turns out to be non-trivial, and is crucial for the computational efficiency of our learning algorithm. More specifically, we chose the period of the DFT to reflect the shape of the effective support of our PMD. (See Proposition 3.8 in Section 3 for the statement.)

With Fourier sparsity as our starting point, we obtain new structural results of independent interest for PMDs. The first is a “robust” moment-matching lemma, which we now informally state:

Informal Lemma (Parameter Moment Closeness Implies Closeness in Distribution.) For any pair of $(n,k)$-PMDs $\mathbf{P},\mathbf{Q}$, if the “low-degree” parameter moment profiles of $\mathbf{P}$ and $\mathbf{Q}$ are close, then $\mathbf{P},\mathbf{Q}$ are close in total variation distance.

See Definition 2.2 for the definition of parameter moments of a PMD. The formal statement of the aforementioned lemma appears as Lemma 4.6 in Section 4.1. Our robust moment-matching lemma is the basis for our proper cover algorithm and our EPTAS for Nash equilibria in anonymous games. Our constructive cover upper bound is the following:

A sparse proper cover quantifies the “size” of the space of PMDs and provides useful structural information that can be exploited in a variety of applications. In addition to Nash equilibria in anonymous games, our efficient proper cover construction provides a smaller search space for approximately solving essentially any optimization problem over PMDs. As another corollary of our cover construction, we obtain the first EPTAS for computing threat points in anonymous games.

Perhaps surprisingly, we also prove that our above upper bound is essentially tight:

For any $k>2$, $\epsilon>0$ sufficiently small as a function of $k,$ and $n=\Omega_{k}(\log(1/\epsilon)/\log\log(1/\epsilon))^{k-1}$, any $\epsilon$-cover for ${\cal M}_{n,k}$ has size at least $n^{\Omega(k)}\cdot(1/\epsilon)^{\Omega_{k}(\log(1/\epsilon)/\log\log(1/\epsilon))^{k-1}}.$

We remark that, in previous work [DKS15a], the authors proved a tight cover size bound of $n\cdot(1/\epsilon)^{\Theta(k\log(1/\epsilon))}$ for $(n,k)$-SIIRVs, i.e., sums of $n$ independent scalar random variables each supported on $[k].$ While a cover size lower bound for $(n,k)$-SIIRVs directly implies the same lower bound for $(n,k)$-PMDs, the opposite is not true. Indeed, Theorems 1.4 and 1.5 show that covers for $(n,k)$-PMDs are inherently larger, requiring a doubly exponential dependence on $k.$

3 Our Approach and Techniques

At a high-level, the Fourier techniques of this paper can be viewed as a highly non-trivial generalization of the techniques in our recent paper [DKS15a] on sums of independent scalar random variables. We would like to emphasize that a number of new conceptual and technical ideas are required to overcome the various obstacles arising in the multi-dimensional setting.

We start with an intuitive explanation of two key ideas that form the basis of our approach.

Since the Fourier Transform (FT) of a PMD is the product of the FTs of its component CRVs, its magnitude is the product of terms each bounded from above by $1.$ Note that each term in the product is strictly less than $1$ except in a small region, unless the component CRV is trivial (i.e., essentially deterministic). Roughly speaking, to establish the sparsity of the FT of PMDs, we proceed as follows: We bound from above the magnitude of the FT by the FT of a Gaussian with the same covariance matrix as our PMD. (See, for example, Lemma 3.10.) This gives us tail bounds for the FT of the PMD in terms of the FT of this Gaussian, and when combined with the concentration of the PMD itself, yields the desired property.

A key ingredient in our proofs is the approximation of the logarithm of the Fourier Transform (log FT) of PMDs by low-degree polynomials. Observe that the log FT is a sum of terms, which is convenient for the analysis. We focus on approximating the log FT by a low-degree Taylor polynomial within the effective support of the FT. (Note that outside the effective support the log FT can be infinity.) Morally speaking, the log FT is smooth, i.e., it is approximated by the first several terms of its Taylor series. Formally however, this statement is in general not true and requires various technical conditions, depending on the setting. One important point to note is that the sparsity of the FT controls the domain in which this approximation will need to hold, and thus help us bound the Taylor error. We will need to ensure that the sizes of the Taylor coefficients are not too large given the location of the effective support, which turns out to be a non-trivial technical hurdle. To ensure this, we need to be very careful about how we perform this Taylor expansion. In particular, the correct choice of the point that we Taylor expand around will be critical for our applications. We elaborate on these difficulties in the relevant technical sections. Finally, we remark that the degree of polynomial approximation we will require depends on the setting: In our cover upper bounds, we will require (nearly) logarithmic degree, while for our CLT degree-$2$ approximation suffices.

We are now ready to give an overview of the ideas in the proofs of each of our results.

The high-level structure of our learning algorithm relies on the sparsity of the Fourier transform, and is similar to the algorithm in our previous work [DKS15a] for learning sums of independent integer random variables. More specifically, our learning algorithm estimates the effective support of the DFT, and then computes the empirical DFT in this effective support. This high-level description would perhaps suffice, if we were only interested in bounding the sample complexity. In order to obtain a computationally efficient algorithm, it is crucial to use the appropriate definition of the DFT and its inverse.

The main structural property needed for the analysis of our algorithm is that there exists an explicit set $T$ with integer coordinates and cardinality $(k\log(1/\epsilon))^{O(k)}$ that contains all but $O(\epsilon)$ of the $L_{1}$ mass of ${\widehat{\mathbf{P}}}.$ Given this property, our algorithm draws an additional set of samples of size $(k\log(1/\epsilon))^{O(k)}/\epsilon^{2}$ from the PMD, and computes the empirical DFT (modulo $L$) on its effective support $T.$ Using these ingredients, we are able to show that the inverse of the empirical DFT defines a pseudo-distribution that is $\epsilon$-close to our unknown PMD in total variation distance.

Observe that the support of the inverse DFT can be large, namely $\Omega(n^{k-1}).$ Our algorithm does not explicitly evaluate the inverse DFT at all these points, but outputs a succinct description of its hypothesis $\mathbf{H}$, via its DFT ${\widehat{\mathbf{H}}}.$ We emphasize that this succinct description suffices to efficiently obtain both an approximate evaluation oracle and an approximate sampler for our target PMD $\mathbf{P}.$ Indeed, it is clear that computing the inverse DFT at a single point can be done in time $O(|T|)=(k\log(1/\epsilon))^{O(k)},$ and gives an approximate oracle for the probability mass function of $\mathbf{P}.$ By using additional algorithmic ingredients, we show how to use an oracle for the DFT, ${\widehat{\mathbf{H}}}$, as a black-box to obtain a computationally efficient approximate sampler for $\mathbf{P}.$

Our learning algorithm and its analysis are given in Section 3.

The correctness of our learning algorithm easily implies (see Section 3.3) an algorithm to construct a non-proper $\epsilon$-cover for PMDs of size $n^{O(k^{2})}\cdot(1/\epsilon)^{\log(1/\epsilon))^{O(k)}}.$ While this upper bound is close to being best possible (see Section 4.5), it does not suffice for our algorithmic applications in anonymous games. For these applications, it is crucial to obtain an efficient algorithm that constructs a proper $\epsilon$-cover, and in fact one that works in a certain stylized way.

To construct a proper cover, we rely on the sparsity of the continuous Fourier Transform of PMDs. Namely, we show that for any PMD $\mathbf{P},$ with effective support $S\subseteq[n]^{k},$ there exists an appropriately defined set $T\subseteq^{k}$ such that the contribution of $\overline{T}$ to the $L_{1}$-norm of $|{\widehat{\mathbf{P}}}|$ is at most $\epsilon/|S|.$ By using this property, we show that any two PMDs, with approximately the same variance in each direction, that have continuous Fourier transforms close to each other in the set $T,$ are close in total variation distance. We build on this lemma to prove our robust moment-matching result. Roughly speaking, we show that two PMDs, with approximately the same variance in each direction, that are “close” to each other in their low-degree parameter moments are also close in total variation distance. We emphasize that the meaning of the term “close” here is quite subtle: we need to appropriately partition the component CRVs into groups, and approximate the parameter moments of the PMDs formed by each group within a different degree and different accuracy for each degree. (See Lemma 4.6 in Section 4.1.)

Our algorithm to construct a proper cover, and our EPTAS for Nash equilibria in anonymous games proceed by a careful dynamic programming approach, that is based on our aforementioned robust moment-matching result.

Finally, we note that combining our moment-matching lemma with a recent result in algebraic geometry gives us the following structural result of independent interest: Every PMD is $\epsilon$-close to another PMD that is a sum of at most $O(k+\log(1/\epsilon))^{k}$ distinct $k$-CRVs.

The aforementioned algorithmic and structural results are given in Section 4.

As mentioned above, a crucial ingredient of our cover upper bound is a robust moment-matching lemma, which translates closeness between the low-degree parameter moments of two PMDs to closeness between their Fourier Transforms, and in turn to closeness in total variation distance. To prove our cover lower bound, we follow the opposite direction. We construct an explicit set of PMDs with the property that any pair of distinct PMDs in our set have a non-trivial difference in (at least) one of their low-degree parameter moments. We then show that difference in one of the parameter moments implies that there exists a point where the probability generating functions have a non-trivial difference. Notably, our proof for this step is non-constructive making essential use of Cauchy’s integral formula. Finally, we can easily translate a pointwise difference between the probability generating functions to a non-trivial total variation distance error. We present our cover lower bound construction in Section 4.5.

The basic idea of the proof of our CLT will be to compare the Fourier transform of our PMD $X$ to that of the discrete Gaussian $G$ with the same mean and covariance. By taking the inverse Fourier transform, we will be able to conclude that these distributions are pointwise close. A careful analysis using a Taylor approximation and the fact that both ${\widehat{X}}$ and ${\widehat{G}}$ have small effective support, gives us a total variation distance error independent of the size $n.$ Alas, this approach results in an error dependence that is exponential in $k.$ To obtain an error bound that scales polynomially with $k,$ we require stronger bounds between $X$ and $G$ at points away from the mean. Intuitively, we need to take advantage of cancellation in the inverse Fourier transform integrals. To achieve this, we will use the saddlepoint method from complex analysis. The full proof of our CLT is given in Section 5.

4 Related and Prior Work

There is extensive literature on distribution learning and computation of approximate Nash equilibria in various classes of games. We have already mentioned the most relevant references in the introduction.

Daskalakis et al. [DKT15] studied the structure and learnability of PMDs. They obtained a non-proper $\epsilon$-cover of size $n^{k^{2}}\cdot 2^{O(k^{5k}\log(1/\epsilon)^{k+2})},$ and an information-theoretic upper bound on the learning sample complexity of $O(k^{5k}\log(1/\epsilon)^{k+2}/\epsilon^{2}).$ The dependence on $1/\epsilon$ in their cover size is also quasi-polynomial, but is suboptimal as follows from our upper and lower bounds. Importantly, the [DKT15] construction yields a non-proper cover. As previously mentioned, a proper cover construction is necessary for our algorithmic applications. We note that the learning algorithm of [DKT15] relies on enumeration over a cover, hence runs in time quasi-polynomial in $1/\epsilon,$ even for $k=2.$ The techniques of [DKT15] are orthogonal to ours. Their cover upper bound is obtained by a clever black-box application of the CLT of [VV10], combined with a non-robust moment-matching lemma that they deduce from a result of Roos [Roo02]. We remind the reader that our Fourier techniques strengthen both these technical tools: Theorem 1.3 strengthens the CLT of [VV10], and we prove a robust and quantitatively essentially optimal moment-matching lemma.

In recent work [DKS15a], the authors used Fourier analytic techniques to study the structure and learnability of sums of independent integer random variables (SIIRVs). The techniques of this paper can be viewed as a (highly nontrivial) generalization of those in [DKS15a]. We also note that the upper bounds we obtain in this paper for learning and covering PMDs do not subsume the ones in [DKS15a]. In fact, our cover upper and lower bounds in this work show that optimal covers for PMDs are inherently larger than optimal covers for SIIRVs. Moreover, the sample complexity of our SIIRV learning algorithm [DKS15a] is significantly better than that of our PMD learning algorithm in this paper.

5 Concurrent and Independent Work

Concurrently and independently to our work, [DDKT16] obtained qualitatively similar results using different techniques. We now provide a statement of the [DDKT16] results in tandem with a comparison to our work.

[DDKT16] give a learning algorithm for PMDs with sample complexity $(k\log(1/\epsilon)^{O(k)}/\epsilon^{2})$ and runtime $(k/\epsilon)^{O(k^{2})}.$ The [DDKT16] algorithm uses the continuous Fourier transform, exploiting its sparsity property, plus additional structural and algorithmic ingredients. The aforementioned runtime is not polynomial in the sample size, unless $k$ is fixed. In contrast, our learning algorithm runs in sample–polynomial time, and, for fixed $k$, in nearly-linear time. The [DDKT16] learning algorithm outputs an explicit hypothesis, which can be easily sampled. On the other hand, our algorithm outputs a succinct description of its hypothesis (via its DFT), and we show how to efficiently sample from it.

[DDKT16] also prove a size-free CLT, analogous to our Theorem 1.3, with error polynomial in $k$ and $1/\sigma.$ Their CLT is obtained by bootstrapping the CLT of [VV10, VV11] using techniques from [DKT15]. As previously mentioned, our proof is technically orthogonal to [VV10, VV11, DDKT16], making use of the sparsity of the Fourier transform combined with tools from complex analysis. It is worth noting that our CLT also achieves a near-optimal dependence in the error as a function of $1/\sigma$ (up to log factors).

Finally, [DDKT16] prove analogues of Theorems 1.2, 1.4, and 1.5 with qualitatively similar bounds to ours. We note that [DDKT16] improve the dependence on $n$ in the cover size to an optimal $n^{O(k)},$ while the dependence on $\epsilon$ in their cover upper bound is the same as in [DKT15]. The cover size lower bound of [DDKT16] is qualitatively of the right form, though slightly suboptimal as a function of $\epsilon.$ The algorithms to construct proper covers and the corresponding EPTAS for anonymous games in both works have running time roughly comparable to the PMD cover size.

6 Organization

In Section 3, we describe and analyze our learning algorithm for PMDs. Section 4 contains our proper cover upper bound construction, our cover size lower bound, and the related approximation algorithm for Nash equilibria in anonymous games. Finally, Section 5 contains the proof of our CLT.

Preliminaries

In this section, we record the necessary definitions and terminology that will be used throughout the technical sections of this paper.

We start by defining our basic object of study:

We will require the following notion of a parameter moment for a PMD:

We now define the notion of distribution learning we use in this paper. Note that an explicit description of a discrete distribution via its probability mass function scales linearly with the support size. Since we are interested in the computational complexity of distribution learning, our algorithms will need to use a succinct description of their output hypothesis. A simple succinct representation of a discrete distribution is via an evaluation oracle for the probability mass function:

One of the most general ways to succinctly specify a distribution is to give the code of an efficient algorithm that takes “pure” randomness and transforms it into a sample from the distribution. This is the standard notion of a sampler:

We can now give a formal definition of distribution learning:

Let ${\cal D}$ be a family of distributions. A randomized algorithm $A^{\cal D}$ is a distribution learning algorithm for class $\cal D,$ if for any $\epsilon>0,$ and any $\mathbf{P}\in\cal D,$ on input $\epsilon$ and sample access to $\mathbf{P},$ with probability $9/10,$ algorithm $A^{\cal D}$ outputs an $\epsilon$-sampler (or an $\epsilon$-evaluation oracle) for $\mathbf{P}.$

We emphasize that our learning algorithm in Section 3 outputs both an $\epsilon$-sampler and an $\epsilon$-evaluation oracle for the target distribution.

Let $X=\sum_{i=1}^{n}X_{i}$ be an $(n,k)$-PMD such that for $1\leq i\leq n$ and $1\leq j\leq k$ we denote $p_{i,j}=\Pr[X_{i}=e_{j}]$, where $\sum_{j=1}^{k}p_{i,j}=1.$ To avoid clutter in the notation, we will sometimes use the symbol $X$ to denote the corresponding probability mass function. With this convention, we can write that ${\widehat{X}}(\xi)=\prod_{i=1}^{n}{\widehat{X_{i}}}(\xi)=\prod_{i=1}^{n}\sum_{j=1}^{k}e(\xi_{j})p_{i,j}.$

Efficiently Learning PMDs

In this subsection, we give an algorithm Efficient-Learn-PMD establishing the following theorem:

Our learning algorithm is described in the following pseudo-code:

Let $X$ be the unknown target $(n,k)$-PMD. We will denote by $\mathbf{P}$ the probability mass function of $X$, i.e., $X\sim\mathbf{P}$. Throughout this analysis, we will denote by $\mu$ and $\Sigma$ the mean vector and covariance matrix of $X.$

First, note that the algorithm Efficient-Learn-PMD is easily seen to have the desired sample and time complexity. Indeed, the algorithm draws $m_{0}$ samples in Step 1 and $m$ samples in Step 5, for a total sample complexity of $O(k^{4k}\log^{2k}(k/\epsilon)/\epsilon^{2}).$ The runtime of the algorithm is dominated by computing the DFT in Step 5 which takes time $O(m|T|)=O(k^{6k}\log^{3k}(k/\epsilon)/\epsilon^{2}).$ Computing an approximate eigendecomposition can be done in time $O(k^{4}\log\log n)$(see, e.g., [PC99]). The remaining part of this section is devoted to proving the correctness of our algorithm.

We begin with a brief overview of the analysis. First, we show (Lemma 3.3) that at least $1-O(\epsilon)$ of the probability mass of $X$ lies in the ellipsoid with center $\mu$ and covariance matrix $\widetilde{\Sigma}=O(k\log(k/\epsilon))\Sigma+O(k\log(k/\epsilon))^{2}I.$ Moreover, with high probability over the samples drawn in Step 1 of the algorithm, the estimates ${\widehat{\Sigma}}$ and ${\widehat{\mu}}$ will be good approximations of $\Sigma$ and $\mu$ (Lemma 3.4). By combining these two lemmas, we obtain (Corollary 3.5) that at least $1-O(\epsilon)$ of the probability mass of $X$ lies in the ellipsoid with center ${\widehat{\mu}}$ and covariance matrix $\Sigma^{\prime}=O(k\log(k/\epsilon)){\widehat{\Sigma}}+O(k\log(k/\epsilon))^{2}I.$

Given the above, it it fairly easy to complete the analysis of correctness. For every point in $T$ we can learn the DFT up to absolute error $O(1/\sqrt{m}).$ Since the cardinality of $T$ is appropriately small, this implies that the total error over $T$ is small. The sparsity property of the DFT (Lemma 3.14) completes the proof.

We now proceed with the detailed analysis of our algorithm. We start by showing that PMDs are concentrated with high probability. More specifically, the following lemma shows that an unknown PMD $X$, with mean vector $\mu$ and covariance matrix $\Sigma$, is effectively supported in an ellipsoid centered at $\mu$, whose principal axes are determined by the eigenvectors and eigenvalues of $\Sigma$ and the desired concentration probability:

Let $X$ be an $(n,k)$-PMD with mean vector $\mu$ and covariance matrix $\Sigma.$ For any $0<\epsilon<1,$ consider the positive-definite matrix $\widetilde{\Sigma}=k\ln(k/\epsilon)\Sigma+k^{2}\ln^{2}(k/\epsilon)I$. Then, with probability at least $1-\epsilon/10$ over $X,$ we have that $(X-\mu)^{T}\cdot\widetilde{\Sigma}^{-1}\cdot(X-\mu)=O(1).$

where we used the Cauchy-Schwartz inequality twice, the triangle inequality, and the fact that a $k$-CRV $X_{i}$ with mean $\mu_{i}$ by definition satisfy $\|X_{i}\|_{2}=1,$ and $\|\mu_{i}\|_{1}=1.$

Let $\nu$ be the variance of $u\cdot(X-\mu).$ By Bernstein’s inequality, we obtain that for $t=\sqrt{2\nu\ln(10k/\epsilon)}+2\sqrt{k}\ln(10k/\epsilon)$ it holds

Applying (1) for $u_{j}\cdot(X-\mu),$ with $t_{j}=\sqrt{2\nu_{j}\ln(10k/\epsilon)}+2\sqrt{k}\ln(10k/\epsilon),$ yields that for all $j\in[k]$ we have

Note that this ellipsoid depends on the mean vector $\mu$ and covariance matrix $\Sigma,$ that are unknown to the algorithm. To obtain a bounding ellipsoid that is known to the algorithm, we will use the following lemma (see Appendix A for the simple proof) showing that ${\widehat{\mu}}$ and ${\widehat{\Sigma}}$ are good approximations to $\mu$ and $\Sigma$ respectively.

With probability at least $19/20$ over the samples drawn in Step 1 of the algorithm, we have that $({\widehat{\mu}}-\mu)^{T}\cdot(\Sigma+I)^{-1}\cdot({\widehat{\mu}}-\mu)=O(1)$, and $2(\Sigma+I)\succeq{\widehat{\Sigma}}+I\succeq(\Sigma+I)/2.$

We also need to deal with the error introduced in the eigendecomposition of ${\widehat{\Sigma}}$. Concretely, we factorize ${\widehat{\Sigma}}$ as $V^{T}\Lambda V,$ for an orthogonal matrix $V$ and diagonal matrix $\Lambda.$ This factorization is necessarily inexact. By increasing the precision to which we learn ${\widehat{\Sigma}}$ by a constant factor, we can still have $2(\Sigma+I)\succeq V^{T}\Lambda V+I\succeq(\Sigma+I)/2.$ We could redefine ${\widehat{\Sigma}}$ in terms of our computed orthonormal eigenbasis, i.e., ${\widehat{\Sigma}}:=V^{T}\Lambda V$. Thus, we may henceforth assume that the decomposition ${\widehat{\Sigma}}=V^{T}\Lambda V$ is exact.

For the rest of this section, we will condition on the event that the statements of Lemma 3.4 are satisfied. By combining Lemmas 3.3 and 3.4 , we show that we can get a known ellipsoid containing the effective support of $X,$ by replacing $\mu$ and $\Sigma$ in the definition of $\cal E$ by their sample versions. More specifically, we have the following corollary:

Let $\Sigma^{\prime}=k\ln(k/\epsilon){\widehat{\Sigma}}+k^{2}\ln^{2}(k/\epsilon)I$. Then, with probability at least $1-\epsilon/10$ over $X,$ we have that $(X-{\widehat{\mu}})^{T}\cdot(\Sigma^{\prime})^{-1}\cdot(X-{\widehat{\mu}})=O(1).$

By Lemma 3.4, it holds that $2(\Sigma+I)\succeq{\widehat{\Sigma}}+I\succeq(\Sigma+I)/2$. Hence, we have that

In terms of $\Sigma^{\prime}$ and $\widetilde{\Sigma}$, this is $2\widetilde{\Sigma}\succeq\Sigma^{\prime}\succeq\frac{1}{2}\widetilde{\Sigma}$. By standard results, taking inverses reverses the positive semi-definite ordering (see e.g., Corollary 7.7.4 (a) in [HJ85]). Hence,

Combining the above with Lemma 3.3, with probability at least $1-\epsilon/10$ over $X$ we have that

Since $\Sigma^{\prime}\succeq\frac{1}{2}\widetilde{\Sigma}\succeq\Sigma+I$, and therefore $(\Sigma^{\prime})^{-1}\preceq(\Sigma+I)^{-1}$, Lemma 3.4 gives that

where the last equality follows from (2) and (3). This completes the proof of Corollary 3.5. ∎

With probability at least $1-\epsilon/10$ over $X,$ we have that $X\in{\widehat{\mu}}+M(-1/2,1/2]^{k}.$

For a large enough constant $C$, Corollary 3.5 implies that with probability at least $1-\epsilon/10,$

Note that the above is an equivalent description of the ellipsoid $\mathcal{E}^{\prime}.$ Our lemma will follow from the following claim:

Similarly, we get $\|M^{T}x\|_{2}\leq\|(M^{\prime})^{T}x\|_{2}+\|(M-M^{\prime})^{T}\|_{2}\|x\|_{2}<2\|(M^{\prime})^{T}x\|_{2}.$ In terms of the PSD ordering, we have:

Since $M^{\prime}(M^{\prime})^{T}\succeq I$, both $M^{\prime}(M^{\prime})^{T}$ and $MM^{T}$ are positive-definite, and so $M$ and $M^{\prime}$ are invertible. Taking inverses in Equation (7) reverses the ordering, that is:

Hence, Claim 3.7 implies that with probability at least $1-\epsilon/10$, we have:

where the last inequality follows from (5). In other words, with probability at least $1-\epsilon/10$, $X$ lies in ${\widehat{\mu}}+M(-1/2,1/2]^{k}$, which was to be proved. ∎

The main component of the analysis is the following proposition, establishing that the total contribution to the above sum coming from points $\xi\not\in T$ is small. In particular, we prove the following:

Consider the $k$ fractional parts of the coordinates of $\xi$, i.e., $\xi_{i}-\lfloor\xi_{i}\rfloor$, for $1\leq i\leq k$. Now consider the $k+1$ intervals $I_{a^{\prime}}=\left(\frac{a^{\prime}-1}{k+1},\frac{a^{\prime}}{k+1}\right)$, for $1\leq a^{\prime}\leq k+1$. By the pigeonhole principle, there is an $a^{\prime}$ such that $\xi_{i}-\lfloor\xi_{i}\rfloor\notin I_{a^{\prime}}$, for all $i$, $1\leq i\leq k.$ We define $a=a^{\prime}$ when $a^{\prime}<k+1$, and $a=0$ when $a^{\prime}=k+1$.

For any $i$, with $1\leq i\leq k$, since $\xi_{i}-\lfloor\xi_{i}\rfloor\notin I_{a^{\prime}}$, we have that $\xi_{i}-\lfloor\xi_{i}\rfloor\in\left[0,\frac{a-1}{k+1}\right]\cup\left[\frac{a}{k+1},1\right]$ (taking the first interval to be empty if $a=0$). Hence, by setting one of $b_{i}=\lfloor\xi_{i}\rfloor$, or $b_{i}=\lfloor\xi_{i}\rfloor-1$, we get $\xi_{i}-b_{i}\in\left[\frac{a}{k+1},\frac{a+k}{k+1}\right].$ This completes the proof. ∎

The following lemma gives a “Gaussian decay” upper bound on the magnitude of the DFT, at points $\xi$ whose coordinates lie in an interval of length less than $1.$ Roughly speaking, the proof of Proposition 3.8 proceeds by applying this lemma for all $\xi\notin T.$

Since $\mathbf{P}$ is a PMD, we have $X=\sum_{i=1}^{n}X_{i}$, where $X_{i}\sim\mathbf{P}_{i}$ for independent $k$-CRV’s $\mathbf{P}_{i}$, we have that $|{\widehat{\mathbf{P}}}(\xi)|=\prod_{i=1}^{n}|{\widehat{\mathbf{P}_{i}}}(\xi)|.$ Note also that $\xi^{T}\cdot\Sigma\cdot\xi=\text{Var}[\xi\cdot X]=\sum_{i=1}^{n}\text{Var}[\xi\cdot X_{i}].$ It therefore suffices to show that for each $i\in[n]$ it holds

We will need the following technical claim:

When $0\leq|x|\leq 1/4$, $\sin(2\pi x)$ is concave, since its second derivative is $-4\pi^{2}\sin(2\pi x)\leq 0$. So, we have $\sin(2\pi x)\geq(1-4x)\sin(0)+4x\sin(\pi/2)=4x.$ Integrating the latter inequality, we obtain that $(1-\cos(2\pi x))/2\pi\geq 2x^{2}$, i.e., $(1-\cos(2\pi x))\geq 4\pi x^{2}$. Thus, for $0\leq|x|\leq 1/4$, we have $1-\cos(2\pi x))\geq 4\pi x^{2}\geq\delta^{2}x^{2}$.

When $1/4\leq|x|\leq 3/4$, we have $(1-\cos(2\pi x))\geq 1\geq\delta^{2}x^{2}$. Finally, when $3/4\leq|x|\leq 1-\delta$, we have $0\leq 1-|x|\leq\delta\leq 1/4$, and therefore $1-\cos(2\pi x)=1-\cos(2\pi(1-|x|))\geq 1-\cos(2\pi\delta)\geq 4\pi\delta^{2}\geq\delta^{2}x^{2}$. This establishes the proof of the claim. ∎

Since $\xi\cdot Y_{i}$ is by assumption supported on the interval $[-1+\delta,1-\delta]$, we have that $|{\widehat{\mathbf{P}_{i}}}(\xi)|^{2}$ is

This completes the proof of Lemma 3.10. ∎

We are now ready to prove the following crucial lemma, which shows that the DFT of $\mathbf{P}$ is effectively supported on the set $T.$

For integers $0\leq a\leq k$, we have that

Although the integer points in the above sum are not in the sphere $\|v\|_{2}\leq C^{2}k^{2}\ln(k/\epsilon),$ they lie in some sphere $\|v\|_{2}\leq 2^{t+1}C^{2}k^{2}\ln(k/\epsilon)$, for some integer $t>0$. The number of integral points in one of these spheres is less than that of the appropriate enclosing cube. Namely, we have that

Inequality (8) is obtained by bounding the LHS from above as follows:

This completes the proof of Lemma 3.12. ∎

We are now prepared to prove Proposition 3.8.

Our next simple lemma states that the empirical DFT is a good approximation to the true DFT on the set $T.$

Letting $m=(C^{5}k^{4}\ln^{2}(k/\epsilon))^{k}/\epsilon^{2}$, with $19/20$ probability over the choice of $m$ samples in Step 5, we have that $\sum_{\xi\in T}|{\widehat{\mathbf{H}}}(\xi)-{\widehat{\mathbf{P}}}(\xi)|<\epsilon/10.$

For any given $\xi\in T,$ we note that ${\widehat{\mathbf{H}}}(\xi)$ is the average of $m$ samples from $e(\xi\cdot X)$, a random variable whose distribution has mean ${\widehat{\mathbf{P}}}(\xi)$ and variance at most $O(1)$. Therefore, we have that

Summing over $\xi\in T,$ and noting that $|T|\leq O(C^{2}k^{2}\log(k/\epsilon))^{k}$, we get that the expectation of the quantity in question is less than $\epsilon/400.$ Markov’s inequality completes the argument. ∎

Finally, we bound from above the total variation distance between $\mathbf{P}$ and $\mathbf{H}.$

where the last line follows from Proposition 3.8 and Lemma 3.13. ∎

This completes the analysis and the proof of Theorem 3.1.

2 An Efficient Sampler for our Hypothesis

The learning algorithm of Section 3.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 PMD $\mathbf{P}.$ The running time of this oracle is linear in the size of $T,$ 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 PMD $\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 subsection is not restricted to PMDs, but is more general, applying to all discrete distributions with an approximately sparse DFT (over any dimension) for which an efficient oracle for the DFT is available.

In particular, we prove the following theorem:

We remark that the $\epsilon$-sampler in the above theorem statement can be described as a randomized algorithm that takes as input $M$, $T$, $\widehat{\mathbf{H}}(\xi),$ for $\xi\in T,$ and the Smith normal form decomposition of $M$(see Lemma 3.21).

This section is devoted to the proof of Theorem 3.15. We first handle the case of one-dimensional distributions, and then appropriately reduce the high-dimensional case to the one-dimensional.

We start by providing some high-level intuition. Roughly speaking, we obtain the desired sampler by considering an appropriate definition of a Cumulative Distribution Function (CDF) corresponding to $\mathbf{H}.$ For the $1$-dimensional case (i.e., the case $k=1$ in our theorem statement), the definition of the CDF is clear, and our sampler proceeds as follows: 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, by using the aforementioned remark.

Our first lemma handles the $1$-dimensional case, assuming the existence of an efficient oracle for the CDF:

We begin our analysis by producing an algorithm that works when we are able to exactly sample $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. For the $1$-dimensional case, this follows easily via a closed form expression. For the high-dimensional case, we first obtain a closed form expression for the case that the matrix $M$ is diagonal. We then reduce the general case to the diagonal case, by using a Smith normal form decomposition.

For $\mathbf{H}$ as in Theorem 3.15, we have the following:

In cases (ii) and (iii), we can also compute the embedding of the corresponding ordering onto the integers $[|\det M|]=\{1,2,...,|\det M|\}$, i.e., we can give a monotone bijection $f:S\rightarrow[|\det M|]$ for which we can efficiently compute $f$ and $f^{-1}$ (i.e., with the same running time bound we give for computing $c_{\mathbf{H}}(x)$).

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\cdot x)$ 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)=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(|T|\log(1/\delta)).$

To avoid clutter in the notation, we define $c_{\mathbf{H},i}(x)$ to be one of these sums, i.e.,

where $s_{i}(a^{\prime}_{i},b^{\prime}_{i}):=\sum_{y_{i}=a^{\prime}_{i}}^{b^{\prime}_{i}}e(-\xi_{i}y_{i}).$ As before, this is a geometric series, so either $\xi_{i}=0$, when we have $s_{i}=b^{\prime}_{i}+1-a^{\prime}_{i}$, or $s_{i}=\frac{e(-\xi_{i}a^{\prime}_{i})-e(-\xi(b^{\prime}_{i}+1))}{1-e(-\xi)}$.

We can thus evaluate $c_{\mathbf{H},i}$ in $O(|T|k)$ arithmetic operations and so compute $c_{\mathbf{H}}(x)$ to desired accuracy in $O(|T|k^{2}\log(1/\delta))$ time. We also note that $f:S\rightarrow\{0,1,...,(\prod_{i}|M_{i}|)-1\}$, defined by $f(x):=\sum_{i}(x_{i}-a_{i})\prod_{j=1}^{i}M_{j}$ is a strictly monotone bijection, and that $f$ and $f^{-1}$ can be computed in time $O(k).$ Now, $c_{\mathbf{H}}(f^{-1}(y)$ is the CDF of the distribution on $y\in\{0,1,...,(\prod_{i}|M_{i}|)-1\}$ whose PMF is given by $\mathbf{H}(f^{-1}(y)).$

We will reduce (iii) to (ii). To do this, we use Smith normal form, a canonical factorization of integer matrices:

Note that the Smith normal form satisfies additional conditions on $D$ than those in Lemma 3.21, but we are only interested in finding such a decomposition where $D$ is a diagonal integer matrix.

So, if $y=g(x),$ since $|\det(M)|=|\det(U)|\cdot|\det(D)|\cdot|\det(V)|=|\det(D)|$, we have:

Now we can prove the main theorem of this subsection.

3 Using our Learning Algorithm to Obtain a Cover

As an application of our learning algorithm in Section 3.1, we provide a simple proof that the space $\mathcal{M}_{n,k}$ of all $(n,k)$-PMDs has an $\epsilon$-cover under the total variation distance of size $n^{O(k^{2})}\cdot 2^{O\left(k\log(1/\epsilon)\right)^{O(k)}}.$ Our argument is constructive, yielding an efficient algorithm to construct a non-proper $\epsilon$-cover of this size.

Two remarks are in order: (i) the non-proper cover construction in this subsection does not suffice for our algorithmic applications of Section 4. These applications require the efficient construction of a proper $\epsilon$-cover plus additional algorithmic ingredients. (ii) The upper bound on the cover size obtained here is nearly optimal, as follows from our lower bound in Section 4.5.

The idea behind using our algorithm to obtain a cover is quite simple. In order to determine its hypothesis, $\mathbf{H}$, our algorithm Efficient-Learn-PMD requires the following quantities:

Given this information, the analysis in Section 3.1 carries over immediately. The algorithm Efficient-Learn-PMD works by estimating the mean and covariance using samples, and then taking ${\widehat{\mathbf{H}}}(\xi)$ to be the sample Fourier transform. If we instead guess the values of these quantities using an appropriate discretization, we obtain an $\epsilon$-cover for $\mathcal{M}_{n,k}.$ More specifically, we discretize the above quantities as follows:

We claim that for any $(n,k)$-PMD $\mathbf{P}$ there exists a choice of parameters, so that the returned distribution $\mathbf{H}$ is within total variation distance $\epsilon$ of $\mathbf{P}.$ We show this as follows: Let $\mu$ and $\Sigma$ be the true mean and covariance matrix of $\mathbf{P}.$ We have that $\mu\in\mathcal{Y}$ and that $\Sigma\in\mathcal{S}.$ Therefore, there exist ${\widehat{\mu}}\in\mathcal{Y}_{1}$ and ${\widehat{\Sigma}}\in\mathcal{S}_{1/2}$ so that $|\mu-{\widehat{\mu}}|_{2}\leq 1$ and $I/2\succeq\Sigma-{\widehat{\Sigma}}\succeq-I/2.$ It is easy to see that these conditions imply the conclusions of Lemma 3.4. Additionally, we can pick elements of $\mathcal{C}_{\delta}$ in order to make $|{\widehat{\mathbf{H}}}(\xi)-{\widehat{\mathbf{P}}}(\xi)|<\epsilon/(10|T|)$ for each $\xi\in T.$ This will give that $\sum_{\xi\in T}|{\widehat{\mathbf{H}}}(\xi)-{\widehat{\mathbf{P}}}(\xi)|<\epsilon/10.$ In particular, the hypothesis $\mathbf{H}$ indexed by this collection of parameters will be within variation distance $\epsilon$ of $\mathbf{P}.$ Hence, the set we have constructed is an $\epsilon$-cover, and our proof is complete.

Efficient Proper Covers and Nash Equilibria in Anonymous Games

In this section, we give our efficient proper cover construction for PMDs, and our EPTAS for computing Nash equilibria in anonymous games. These algorithmic results are based on new structural results for PMDs that we establish. The structure of this section is as follows: In Section 4.1, we show the desired sparsity property of the continuous Fourier transform of PMDs, and use it to prove our robust moment-matching lemma. Our dynamic-programming algorithm for efficiently constructing a proper cover relies on this lemma, and is given in Section 4.2. By building on the proper cover construction, in Section 4.3 we give our EPTAS for Nash equilibria in anonymous games. In Section 4.4, we combine our moment-matching lemma with recent results in algebraic geometry, to show that any PMD is close to another PMD with few distinct CRV components. Finally, in Section 4.5 we prove out cover size lower bound.

In this subsection, we establish the sparsity of the continuous Fourier transform of PMDs, and use it to prove our robust moment-matching lemma, translating closeness in the low-degree parameter moments to closeness in total variation distance.

At a high-level, our robust moment-matching lemma (Lemma 4.6) is proved by combining the sparsity of the continuous Fourier transform of PMDs (Lemma 4.2) with very careful Taylor approximations of the logarithm of the Fourier transform (log FT) of our PMDs. For technical reasons related to the convergence of the log FT, we will need one additional property from our PMDs. In particular, we require that each component $k$-CRV has the same most likely outcome. This assumption is essentially without loss of generality. There exist at most $k$ such outcomes, and we can express an arbitrary PMD as a sum of $k$ independent component PMDs whose $k$-CRV components satisfy this property. Formally, we have the following definition:

Any $(n,k)$-PMD $X$ can be written as $X=\sum_{i=1}^{k}X^{i}$, where $X^{i}$ is an $i$-maximal $(n_{i},k)$-PMD, with $\sum_{i}n_{i}=n.$ For the rest of this intuitive explanation, we focus on two $(n,k)$-PMDs $X,Y$ that are promised to be $i$-maximal, for some $i\in[k].$

To guarantee that ${\widehat{X}}$, ${\widehat{Y}}$ have roughly the same effective support, we also assume that they have roughly the same variance in each direction. We will show that if the low-degree parameter moments of $X$ and $Y$ are close to each other, then $X$ and $Y$ are close in total variation distance. We proceed by partitioning the $k$-CRV components of our PMDs into groups, based on their maximum probability element $e_{j}$, with $j\neq i.$ The maximum probability of a $k$-CRV quantifies its maximum contribution to the variance of the PMD in some direction. Roughly speaking, the smaller this contribution is, the fewer terms in the Taylor approximation are needed to achieve a given error. More specifically, we consider three different groups, partitioning the component $k$-CRVs into ones with small, medium, and large contribution to the variance in some direction. For the PMD (defined by the CRVs) of the first group, we only need to approximate the first $2$ parameter moments. For the PMD of the second group, we approximate the low-degree parameter moments up to degree $O_{k}(\log(1/\epsilon)/\log\log(1/\epsilon)).$ Finally, the third group is guaranteed to have very few component $k$-CRVS, hence we can afford to approximate the individual parameters.

To quantify the above, we need some more notation and definitions. To avoid clutter in the notation, we focus without loss of generality on the case $i=k$, i.e., our PMDs are $k$-maximal. For a $k$-maximal $(n,k)$-PMD, $X$, let $X=\sum_{i=1}^{n}X_{i}$, where the $X_{i}$ is a $k$-CRV with $p_{i,j}=\Pr[X_{i}=e_{j}]$ for $1\leq i\leq n$ and $1\leq j\leq k.$ Observe that $\sum_{j=1}^{k}p_{i,j}=1,$ for $1\leq i\leq n,$ hence the definition of $k$-maximality implies that $p_{i,k}\geq 1/k$ for all $i.$ Note that the $j^{th}$ component of the random vector $X$ is a PBD with parameters $p_{i,j}$, $1\leq i\leq n.$ Let $s_{j}(X)=\sum_{i=1}^{n}p_{i,j}$ be the expected value of the $j^{th}$ component of $X$. We can assume that $s_{j}(X)\geq\epsilon/k$, for all $1\leq j\leq k-1$; otherwise, we can remove the corresponding coordinates and introduce an error of at most $\epsilon$ in variation distance.

Note that, for $j\neq k$, the variance of the $j^{th}$ coordinate of $X$ is in $[s_{j}(X)/2,s_{j}(X)]$. Indeed, the aforementioned variance equals $\sum_{i=1}^{n}p_{i,j}(1-p_{i,j})$, which is clearly at most $s_{j}(X)$. The other direction follows by observing that, for all $j\neq k$, we have $1\geq p_{i,k}+p_{i,j}\geq 2p_{i,j}$, or $p_{i,j}\leq 1/2$, where we again used the $k$-maximality of $X$. Therefore, by Bernstein’s inequality and a union bound, there is a set $S\subseteq[n]^{k}$ of size

so that $X$ lies in $S$ with probability at least $1-\epsilon$.

We start by showing that the continuous Fourier transform of a PMD is approximately sparse, namely it is effectively supported on a small set $T.$ More precisely, we prove that there exists a set $T$ in the Fourier domain such that the integral of the absolute value of the Fourier transform outside $T$ multiplied by the size of the effective support $|S|$ of our PMD is small.

Let $X$ be $k$-maximal $(k,n)$-PMD with effective support $S.$ Let

where $[x]$ is the distance between $x$ and the nearest integer, and $C>0$ is a sufficiently large universal constant. Then, we have that

To prove the lemma, we will need the following technical claim:

For all $\xi=(\xi_{1},\ldots,\xi_{k})\in^{k}$, for all $1\leq i\leq n$ and $1\leq j\leq k-1$, it holds:

The claim follows from the following sequence of (in-)equalities:

where the last lines uses the fact $p_{i,k}\geq 1/k$, which follows from $k$-maximality. ∎

As a consequence of Claim 4.3, we have that

We now use the sparsity of the Fourier transform to show that if two $k$-maximal PMDs, with similar variances in each direction, have Fourier transforms that are pointwise sufficiently close to each other in this effective support, then they are close to each other in total variation distance.

We start with an intuitive explanation of the proof. Since $1+s_{j}(X)$, $1+s_{j}(Y)$ are within a factor of $2$ for $1\leq j\leq k-1$, it follows from the above that $X$ and $Y$ are both effectively supported on a set $S\subseteq[n]^{k}$ of size $|S|\leq O\left(\log(k/\epsilon)\right)^{k-1}\cdot\prod_{j=1}^{k-1}\left(1+12s_{j}(X)^{1/2}\right)$. Therefore, to prove the lemma, it is sufficient to establish that $\|X-Y\|_{\infty}\leq O(\epsilon/|S|$).

We prove this statement in two steps by analyzing the continuous Fourier transforms ${\widehat{X}}$ and ${\widehat{Y}}$. The first step of the proof exploits the fact that the Fourier transforms of $X$ and $Y$ are each essentially supported on the set $T$ of the lemma statement. Recalling the assumption that $(1+s_{j}(X))/(1+s_{j}(Y))\in[1/2,2]$, $1\leq j\leq k-1$, an application of Lemma 4.2 yields that $\int_{\overline{T}}|{\widehat{X}}|$ and $\int_{\overline{T}}|{\widehat{Y}}|$ are both at most $\epsilon/|S|.$ Thus, we have that