Statistical Query Lower Bounds for Robust Estimation
of High-Dimensional Gaussians and Gaussian Mixtures

arXiv:1611.03473v2 [cs.LG] 17 May 2017

Ilias Diakonikolas∗
University of Southern California
diakonik@usc.edu

Daniel M. Kane†
University of California, San Diego
dakane@cs.ucsd.edu

Alistair Stewart
University of Southern California
alistais@usc.edu
May 18, 2017
Abstract
We describe a general technique that yields the first Statistical Query lower bounds for a range of
fundamental high-dimensional learning problems involving Gaussian distributions. Our main results are
for the problems of (1) learning Gaussian mixture models (GMMs), and (2) robust (agnostic) learning
of a single unknown Gaussian distribution. For each of these problems, we show a super-polynomial
gap between the (information-theoretic) sample complexity and the computational complexity of any
Statistical Query algorithm for the problem. Statistical Query (SQ) algorithms are a class of algorithms
that are only allowed to query expectations of functions of the distribution rather than directly access
samples. This class of algorithms is quite broad: a wide range of known algorithmic techniques in
machine learning are known to be implementable using SQs. Moreover, for the unsupervised learning
problems studied in this paper, all known algorithms with non-trivial performance guarantees are SQ or
are easily implementable using SQs.
Our SQ lower bound for Problem (1) is qualitatively matched by known learning algorithms for
GMMs. At a conceptual level, this result implies that – as far as SQ algorithms are concerned – the
computational complexity of learning GMMs is inherently exponential in the dimension of the latent
space – even though there is no such information-theoretic barrier. Our lower bound for Problem (2)
implies that the accuracy of the robust learning algorithm in [DKK+ 16] is essentially best possible
among all polynomial-time SQ algorithms. On the positive side, we also give a new (SQ) learning
algorithm for Problem (2) achieving the information-theoretically optimal accuracy, up to a constant
factor, whose running time essentially matches our lower bound. Our algorithm relies on a filtering
technique generalizing [DKK+ 16] that removes outliers based on higher-order tensors.
Our SQ lower bounds are attained via a unified moment-matching technique that is useful in other
contexts and may be of broader interest. Our technique yields nearly-tight lower bounds for a number
of related unsupervised estimation problems. Specifically, for the problems of (3) robust covariance
estimation in spectral norm, and (4) robust sparse mean estimation, we establish a quadratic statistical–
computational tradeoff for SQ algorithms, matching known upper bounds. Finally, our technique can
be used to obtain tight sample complexity lower bounds for high-dimensional testing problems. Specifically, for the classical problem of robustly testing an unknown mean (known covariance) Gaussian, our
technique implies an information-theoretic sample lower bound that scales linearly in the dimension.
Our sample lower bound matches the sample complexity of the corresponding robust learning problem
and separates the sample complexity of robust testing from standard (non-robust) testing. This separation
is surprising because such a gap does not exist for the corresponding learning problem.
∗
†

Supported by NSF Award CCF-1652862 (CAREER) and a Sloan Research Fellowship.
Supported by NSF Award CCF-1553288 (CAREER) and a Sloan Research Fellowship.

1 Introduction
1.1 Background and Overview For the unsupervised estimation problems considered here, the input is
a probability distribution which is accessed via a sampling oracle, i.e., an oracle that provides i.i.d. samples
from the underlying distribution. Statistical Query (SQ) algorithms are a restricted class of algorithms that
are only allowed to query expectations of functions of the distribution rather than directly access samples.
This class of algorithms is quite broad: a wide range of known algorithmic techniques in machine learning
are known to be implementable using SQs. These include spectral techniques, moment and tensor methods,
local search (e.g., Expectation Maximization), and many others (see, e.g., [CKL+ 06, FGR+ 13] for a detailed
discussion). Moreover, for the unsupervised learning problems studied in this paper, all known algorithms
with non-trivial performance guarantees are SQ or are easily implementable using SQs.
A number of techniques have been developed in information theory and statistics to characterize the
sample complexity of inference tasks. These involve both techniques for proving sample complexity upper
bounds (e.g., VC dimension, metric/bracketing entropy) and information-theoretic lower bounds (e.g., Fano
and Le Cam methods). On the other hand, computational lower bounds have been much more scarce in
the unsupervised setting. Perhaps surprisingly, it is possible to prove unconditional lower bounds on the
computational complexity of any SQ algorithm that solves a given learning problem. Given the ubiquity and
generality of SQ algorithms, an SQ lower bound provides strong evidence of the problem’s computational
intractability.
In this paper, we describe a general technique that yields the first Statistical Query lower bounds for a
range of fundamental high-dimensional learning problems involving Gaussian distributions. Such problems
are ubiquitous in applications across the data sciences and have been intensely investigated by different
communities of researchers for several decades. Our main results are for the problems of (1) learning Gaussian mixture models (GMMs), and (2) robust (agnostic) learning of a single unknown Gaussian distribution.
In particular, we show a super-polynomial gap between the (information-theoretic) sample complexity and
the computational complexity of any Statistical Query algorithm for these problems. In more detail, our SQ
lower bound for Problem (1) is qualitatively matched by known learning algorithms for GMMs (all of which
can be implemented as SQ algorithms). For Problem (2), we give a new (SQ) algorithm in this paper whose
running time nearly matches our SQ lower bound.
Our SQ lower bounds are attained via a unified moment-matching technique that is useful in other
contexts and may be of broader interest. Our technique yields nearly-tight lower bounds for a number of
related unsupervised estimation problems. Specifically, for the problems of (3) robust covariance estimation
in spectral norm, and (4) robust sparse mean estimation, we establish a quadratic statistical–computational
tradeoff for SQ algorithms, matching known upper bounds.
Finally, we use our technique to obtain tight sample complexity lower bounds for high-dimensional
testing problems. Specifically, for the classical problem of robustly testing an unknown mean (known covariance) Gaussian, our technique implies an information-theoretic lower bound that scales linearly in the
dimension. This lower bound matches the sample complexity of the corresponding robust learning problem
and separates the sample complexity of robust testing from standard (non-robust) testing. This separation is
surprising because such a gap does not exist for the corresponding learning problem.
Before we discuss our contributions in detail, we provide the necessary background for the Statistical
Query model and the unsupervised estimation problems that we study.
Statistical Query Algorithms. A Statistical Query (SQ) algorithm relies on an oracle that given any
bounded function on a single domain element provides an estimate of the expectation of the function on
a random sample from the input distribution. This computational model was introduced by Kearns [Kea98]
in the context of supervised learning as a natural restriction of the PAC model [Val84]. Subsequently, the
SQ model has been extensively studied in a plethora of contexts (see, e.g., [Fel16b] and references therein).
1

A recent line of work [FGR+ 13, FPV15, FGV15, Fel16a] developed a framework of SQ algorithms
for search problems over distributions – encompassing the distribution estimation problems we study in
this work. It turns out that one can prove unconditional lower bounds on the computational complexity
of SQ algorithms via the notion of Statistical Query dimension. This complexity measure was introduced
in [BFJ+ 94] for PAC learning of Boolean functions and was recently generalized to the unsupervised setting [FGR+ 13, Fel16a]. A lower bound on the SQ dimension of a learning problem provides an unconditional lower bound on the computational complexity of any SQ algorithm for the problem.
Remark. We would like to emphasize here that the SQ lower bounds shown in this paper apply to the
running time of an SQ algorithm and not on its sample complexity (when we simulate the SQ algorithm by
drawing samples to answer its SQ queries). Specifically, for all learning problems considered in this paper,
there exist straightforward SQ algorithms (that can be simulated with sample access to the distribution) with
near-optimal sample complexity, albeit with exponential running time. Specifically, lower bounds on the SQ
dimension of the corresponding problems establish lower bounds on the running time of any SQ algorithm
for the problem – not on its sample complexity.
Learning Gaussian Mixture Models. A mixture model is a convex combination of distributions of known
type. The most commonly studied case is a Gaussian mixture model (GMM). An n-dimensional
k-GMM
P
is a distribution in Rn that is composed of k unknown Gaussian components, i.e., F = ki=1 wi N (µi , Σi ),
where the weights wi , mean vectors µi , and covariance matrices Σi are unknown. The problem of learning
a GMM from samples has received tremendous attention in statistics and, more recently, in TCS. A long
line of work initiated by Dasgupta [Das99, AK01, VW02, AM05, KSV08, BV08] provides computationally efficient algorithms for recovering the parameters of a GMM under separability assumptions. Subsequently, efficient parameter learning algorithms have been obtained [MV10, BS10, HP15] under minimal
information-theoretic separation assumptions. The related problems of density estimation and proper learning have also been extensively studied [FOS06, SOAJ14, DK14, MV10, HP15, LS15]. In density estimation
(resp. proper learning), the goal is to output some hypothesis (resp. GMM) that is close to the unknown
mixture in total variation distance.
The sample complexity of density estimation (and proper learning) for n-dimensional k-GMMs, up
to variation distance ǫ, is easily seen to be poly(n, k, 1/ǫ) – without any assumptions.(In Appendix A,
we describe a simple SQ algorithm for this learning problem with sample complexity poly(n, k, 1/ǫ), albeit
exponential running time). Given that there is no information-theoretic barrier for learnability in this setting,
the following question arises: Is there a poly(n, k, 1/ǫ) time algorithm for density estimation (or proper
learning) of n-dimensional k-GMMs? This question has been raised as an open problem in a number of
settings (see, e.g., [Moi14, Dia16] and references therein).
For parameter learning, the situation is somewhat subtle: In full generality, the sample complexity
is of the form poly(n) · (1/γ)Ω(k) , where the parameter γ > 0 quantifies the “separation” between the
components. Even in one-dimension, a sample complexity lower bound of (1/γ)Ω(k) is known [MV10,
HP15]1 . The corresponding “hard” instances [MV10, HP15] consist of GMMs whose components have
large overlap, so many samples are required to distinguish between them. Is this the only obstacle towards
a poly(n, k) time parameter learning algorithm? Specifically, suppose that we are given an instance of the
problem with the additional promise that the components are “nearly non-overlapping” – so that poly(n, k)
samples suffice for the parameter learning problem as well. (In Appendix B, we show that when the total
variation distance between any pair of components in the given mixture is close to 1, parameter learning
reduces to proper learning; hence, there is a poly(n, k)-sample parameter learning (SQ) algorithm that runs
in exponential time.) Is there a poly(n, k) time parameter learning algorithm for such instances?
1
To circumvent the information-theoretic bottleneck of parameter learning, a related line of work has studied parameter learning
in a smoothed setting [HK13, BCMV14, ABG+ 14, GHK15].

2

In summary, the sample complexity of both versions of the learning problem is poly(n)f (k). On the
other hand, the running time of all known algorithms for either version scales as ng(k) , where g(k) ≥ k.
This runtime is super-polynomial in the sample complexity of the problem for super-constant values of k
and is tight for these algorithms, even for GMMs with almost non-overlapping components. The preceding
discussion is summarized in the following:
Question 1.1. Is there a poly(n, k)-time density estimation algorithm for n-dimensional k-GMMs? Is there
a poly(n, k)-time parameter learning algorithm for nearly non-overlapping n-dimensional k-GMMs?
Robust Learning of a Gaussian. In the preceding paragraphs, we were working under the assumption
that the unknown distribution generating the samples is exactly a mixture of Gaussians. The more general
and realistic setting of robust (or agnostic) learning – when our assumption about the model is approximately
true – turns out to be significantly more challenging. Specifically, until recently, even the most basic setting
of robustly learning an unknown mean Gaussian with identity covariance matrix was poorly understood.
Without corruptions, this problem is straightforward: The empirical mean gives a sample-optimal efficient
estimator. Unfortunately, the empirical estimate is very brittle and fails in the presence of corruptions.
The standard definition of agnostically learning a Gaussian (see, e.g., Definition 2.1 in [DKK+ 16] and
references therein) is the following: Instead of drawing samples from a perfect Gaussian, we have access to a
distribution D that is promised to be close to an unknown Gaussian G – specifically ǫ-close in total variation
distance. This is the only assumption about the distribution D, which may otherwise be arbitrary: the ǫfraction of “errors” can be adversarially selected. The goal of an agnostic learning algorithm is to output
a hypothesis distribution H that is as close as possible to G (or, equivalently, D) in variation distance.
Note that the minimum variation distance, dTV (H, G), information-theoretically achievable under these
assumptions is Θ(ǫ), and we would like to obtain a polynomial-time algorithm with this guarantee.
Agnostically learning a single high-dimensional Gaussian is arguably the prototypical problem in robust
statistics [Hub64, HRRS86, HR09]. Early work in this field [Tuk75, DG92] studied the sample complexity
of robust estimation. Specifically, for the case of an unknown mean and known covariance Gaussian, the
Tukey median [Tuk75] achieves O(ǫ)-error with O(n/ǫ2 ) samples (see, e.g., [CGR15] for a simple proof).
Since Ω(n/ǫ2 ) samples are information-theoretically necessary – even without noise – the robustness requirement does not change the sample complexity of the problem.
The computational complexity of agnostically learning a Gaussian is less understood. Until recently,
√
all known polynomial time estimators could only guarantee error of Θ(ǫ n). Two recent works [DKK+ 16,
LRV16] made a first step in designing robust polynomial-time estimators for this problem. The results
of [DKK+ 16] apply in the standard agnostic model; [LRV16] works in a weaker model – known as Huber’s contamination model [Hub64] – where the noisy distribution D is of the form (1 − ǫ)G + ǫN ,
where N is an unknown “noise” distribution. For the problem of
√ robustly estimating +an unknown mean
log n), while [DKK 16] obtains error
Gaussian
N
(µ,
I),
[LRV16]
obtains
an
error
guarantee
of
O(ǫ
p
O(ǫ log(1/ǫ)), independent of the dimension2 .
A natural and important open problem, put forth by these works [DKK+ 16, LRV16], is the following:
Question 1.2. Is there a poly(n/ǫ)- time agnostic learning algorithm, with error O(ǫ), for an n-dimensional
Gaussian?
Statistical–Computational Tradeoffs. A statistical–computational tradeoff refers to the phenomenon that
there is an inherent gap between the information-theoretic sample complexity of a learning problem and its
computational sample complexity, i.e, the minimum sample complexity attainable by any polynomial time
2

The algorithm
p of [LRV16] can be extended to work in the standard agnostic model at the expense of an increased error
guarantee of O(ǫ log n log(1/ǫ)).

3

algorithm for the problem. The prototypical example is the estimation of a covariance matrix under sparsity
constraints (sparse PCA) [JL09, CMW13, CMW15], where a nearly-quadratic gap between informationtheoretic and computational sample complexity has been established (see [BR13b, WBS16b]) – assuming
the computational hardness of the planted clique problem.
For a number of high-dimensional learning problems (including the problem of robustly learning a
Gaussian under the total variation distance), it is known that the robustness requirement does not change
the information-theoretic sample complexity of the problem. On the other hand, it is an intriguing possibility that injecting noise into a high-dimensional learning problem may change its computational sample
complexity.
Question 1.3. Does robustness create inherent statistical–computational tradeoffs for natural high-dimensional
estimation problems?
In this work, we consider two natural instantiations of the above general question: (i) robust estimation
of the covariance matrix in spectral norm, and (ii) robust sparse mean estimation. We give basic background
for these problems in the following paragraphs.
For (i), suppose we have sample access to a (zero-mean) n-dimensional unknown-covariance Gaussian, and we want to estimate the covariance matrix with respect to the spectral norm. It is known (see, e.g.,
[Ver12]) that O(n/ǫ2 ) samples suffice so that the empirical covariance is within spectral error at most ǫ from
the true covariance; and this bound is information-theoretically optimal, to constant factors, for any estimator. For simplicity, let us assume that the desired accuracy is a small positive constant, e.g., ǫ = 1/10. Now
suppose that we observe samples from a corrupted Gaussian in Huber’s contamination model (the weaker
adversarial model) where the noise rate δ ≪ 1/10. First, it is not hard to see that the injection of noise
does not change the information-theoretic sample complexity of the problem: there exist (computationally
inefficient) robust estimators (see, e.g., [CGR15]) that use O(n) samples. (There is a straightforward SQ
algorithm for this problem as well that uses O(n) samples, but again runs in exponential time.) On the other
hand, if we are willing to use Õ(n2 ) samples, a polynomial-time robust estimator with constant spectral
error guarantee is known [DKK+ 16, DKK+ 17a]3 . The immediate question that follows is this:
Is there a computationally efficient robust covariance estimator in spectral error that uses a
strongly sub-quadratic sample size, i.e., O(n2−c ) for a constant 0 < c < 1?
For (ii), suppose we want to estimate the mean µ ∈ Rn of an identity covariance Gaussian up to
ℓ2 -distance ǫ, under the additional promise that µ is k-sparse, and suppose that k ≪ n1/2 . It is wellknown that the information-theoretic sample complexity of this problem is O(k log n/ǫ2 ), and the truncated empirical mean achieves the optimal bound. For simplicity, let us assume that ǫ = 1/10. Now suppose that we observe samples from a corrupted sparse mean Gaussian (in Huber’s contamination model),
where the noise rate δ ≪ 1/10. As in the setting of the previous paragraph, the injection of noise does
not change the information-theoretic sample complexity of the problem: there exist a (computationally
inefficient) robust SQ algorithm for this problem (see [Li17]) that use O(k log n) samples. Two recent
works [Li17, DBS17] gave polynomial time robust algorithms for robust sparse mean estimation with sample complexity Õ(k2 log n). In summary, in the absence of robustness, the information-theoretically optimal
sample bound is known to be achievable by a computationally efficient algorithm. In contrast, in the presence of robustness, there is a quadratic gap between the information-theoretic optimum and the sample
complexity of known polynomial-time algorithms. The immediate question is whether this gap is inherent:
We note that the robust covariance estimators of [DKK+ 16, DKK+ 17a] provide error guarantees under the Mahalanobis
distance, which is stronger than the spectral norm. Under the stronger metric, Ω(n2 ) samples are information-theoretically required
even without noise.
3

4

Is there a computationally efficient robust k-sparse mean estimator that uses a strongly subquadratic sample size , i.e., O(k2−c ) for a constant 0 < c < 1?
It is conjectured in [Li17] that a quadratic gap is in fact inherent for efficient algorithms.
High-Dimensional Hypothesis Testing. So far, we have discussed the problem of learning an unknown
distribution that is promised to belong (exactly or approximately) in a given family (Gaussians, mixtures
of Gaussians). A related inference problem is that of hypothesis testing [NP33, LR05]: Given samples
from a distribution in a given family, we want to distinguish between a null hypothesis and an alternative
hypothesis. Starting with [GR00, BFR+ 00], this broad question has been extensively investigated in TCS
with a focus on discrete probability distributions. A natural way to solve a distribution testing problem is to
learn the distribution in question to good accuracy and then check if the corresponding hypothesis is close to
one satisfying the null hypothesis. This testing-via-learning approach is typically suboptimal and the main
goal in this area has been to obtain testers with sub-learning sample complexity.
In this paper, we study natural hypothesis testing analogues of the high-dimensional learning problems
discussed in the previous paragraphs. Specifically, we study the sample complexity of (i) robustly testing an
unknown mean Gaussian, and (ii) testing a GMM.
To motivate (i), we consider arguably the most basic high-dimensional testing task: Given samples from
a Gaussian N (µ, I), where µ ∈ Rn is unknown, distinguish between the case that µ = 0 versus kµk2 ≥ ǫ.
(The latter condition is equivalent, up to constant factors, to dTV (N (µ, I), N (0, I)) ≥ ǫ.) The classical test
for this task is Hotelling’s T-squared statistic [Hot31], which is unfortunately not defined when the sample
size is smaller than the dimension [ZB96]. More recently, testers that succeed in the sub-linear regime have
been developed [SD08] (also see [ZB96, CQ10]). In Appendix C, we give a simple and natural tester for
√
this problem that uses O( n/ǫ2 ) samples, and show that this sample bound is information-theoretically
optimal, up to constant factors.
Now suppose that our Gaussianity assumption about the unknown distribution is only approximately
satisfied. Formally, we are given samples from a distribution D on Rn which is promised to be either (a)
a standard Gaussian N (0, I), or (b) a δ-noisy version of N (µ, I), where µ ∈ Rn satisfies kµk2 ≥ ǫ, and
the noise rate δ satisfies δ ≪ ǫ. The robust hypothesis testing problem is to distinguish, with high constant
probability, between these two cases. Note that condition (b) implies that dTV (D, N (0, I)) = Ω(ǫ), and
therefore the two cases are distinguishable. 4
Robust hypothesis testing is of fundamental importance and has been extensively studied in robust
statistics [HR09, HRRS86, Wil97]. Perhaps surprisingly, it is poorly understood in the most basic settings,
even information-theoretically. Specifically, the sample complexity of our aforementioned robust mean
testing problem has remained open. It is easy to see that the tester of Appendix C fails in the robust setting.
On the other hand, the testing-via-learning approach implies a sample upper bound of O(n/ǫ2 ) for our
robust testing problem – by using, e.g., the Tukey median. The following question arises:
Question 1.4. Is there an information-theoretic gap between robust testing and non-robust testing? What is
the sample complexity of robustly testing the mean of a high-dimensional Gaussian?
We conclude with our hypothesis testing problem regarding GMMs: Given samples from a distribution D
on Rn , we want to distinguish between the case that D = N (0, I), or D is a 2-mixture of identity covariance
Gaussians. This is a natural high-dimensional testing problem that we believe merits investigation in its own
right. The obvious open question here is whether there exists a tester for this problem with sub-learning
sample complexity.
4

Robust testing should not be confused with tolerant testing, where the completeness is relaxed. In our context, tolerant testing
corresponds to distinguishing
√ between dTV (D, N (0, I)) ≤ ǫ/2 versus dTV (D, N (0, I)) ≥ ǫ, where D = N (µ, I), and is easily
seen to be solvable with O( n/ǫ2 ) samples as well.

5

1.2 Our Results The main contribution of this paper is a general technique to prove lower bounds for
a range of high-dimensional estimation problems involving Gaussian distributions. We use analytic and
probabilistic ideas to construct explicit families of hard instances for the estimation problems described in
Section 1.1. Using our technique, we prove super-polynomial Statistical Query (SQ) lower bounds that
answer Questions 1.1 and 1.2 in the negative for the class of SQ algorithms. We also show that the observed
quadratic statistical–computational gap for robust sparse mean estimation and robust spectral covariance
estimation is inherent for SQ algorithms. As an additional important application of our technique, we obtain
information-theoretic lower bounds on the sample complexity of the corresponding testing problems. (We
note that our testing lower bounds apply to all algorithms.) Specifically, we answer Question 1.4 in the
affirmative, by showing that the robustness requirement makes the Gaussian testing problem informationtheoretically harder. In the body of this section, we state our results and elaborate on their implications and
the connections between them.
SQ Lower Bound for Learning GMMs. Our first main result is a lower bound of nΩ(k) on the complexity
of any SQ algorithm that learns an arbitrary n-dimensional k-GMM to constant accuracy (see Theorem 4.1
for the formal statement):
Theorem 1.1 (SQ Lower Bound for Learning GMMs). Any SQ algorithm that learns an arbitrary nΩ(1)
≥ nΩ(k) queries to an
dimensional k-GMM to constant accuracy, for all n ≥ poly(k), requires 2n
−O(k)
SQ oracle of precision n
.
Theorem 1.1 establishes a super-polynomial gap between the information-theoretic sample complexity of
learning GMMs and the complexity of any SQ learning algorithm for this problem. It is worth noting that
our hard instance is a family of high-dimensional
GMMs whose components are almost non-overlapping.
P
Specifically, for each GMM F = ki=1 wi N (µi , Σi ) in the family, the total variation distance between any
pair of Gaussian components can be made as large as 1 − 1/poly(n, k). More specifically, for our family of
hard instances, the sample complexity of both density and parameter learning is Θ(k · log n) (the standard
cover-based algorithm that achieves this sample upper bound is SQ). In contrast, any SQ learning algorithm
for this family of instances requires runtime at least nΩ(k) .
At a conceptual level, Theorem 1.1 implies that – as far as SQ algorithms are concerned – the computational complexity of learning high-dimensional GMMs is inherently exponential in the dimension of the
latent space – even though there is no such information-theoretic barrier in general. Our SQ lower bound
identifies a common barrier of the strongest known algorithmic approaches for this learning problem, and
provides a rigorous explanation why a long line of algorithmic research on this front either relied on strong
separation assumptions or resulted in runtimes of the form nΩ(k) .
SQ Lower Bound for Robustly Learning a Gaussian. Our second main result concerns the agnostic
learning of a single n-dimensional Gaussian. We prove two SQ lower bounds with qualitatively similar
guarantees for different versions of this problem. Our first lower bound is for the problem of agnostically
learning a Gaussian with unknown mean and identity covariance. Roughly speaking, we show that any SQ
1/4
algorithm that solves this learning problem to accuracy O(ǫ) requires complexity nΩ(log (1/ǫ)) . We show
(see Theorem 5.1 for a more detailed statement):
Theorem 1.2 (SQ Lower Bound for Robust Learning of Unknown Mean Gaussian). Let ǫ > 0, 0 < c ≤ 1/2,
and n ≥ poly(log(1/ǫ)). Any SQ algorithm that robustly learns an n-dimensional Gaussian N (µ, I),
c/2
Ω(1)
≥ nΩ(log(1/ǫ) ) queries to an SQ oracle
within total variation distance O(ǫ log(1/ǫ)1/2−c ), requires 2n
c/2
of precision n−Ω(log(1/ǫ) ) .
Some comments are in order. First, Theorem 1.2 shows a super-polynomial gap between the sample complexity of agnostically learning an unknown mean Gaussian and the complexity of SQ learning algorithms
6

for this problem. As mentioned in the introduction, O(n/ǫ2 ) samples information-theoretically suffice to
agnostically learn an unknown mean Gaussian to within error O(ǫ). Second, the robust learning algo+
rithm
SQ model, and achieves error
pof [DKK 16] runs in poly(n, 1/ǫ) time, can be implemented in thep
O(ǫ log(1/ǫ)). As a corollary of Theorem 1.2, we obtain that the O(ǫ log(1/ǫ)) error guarantee of
the [DKK+ 16] algorithm is best possible among all polynomial-time SQ algorithms.
Roughly speaking, Theorem 1.2 shows that any SQ algorithm that solves the (unknown mean Gaussian)
1/4
robust learning problem to accuracy O(ǫ) needs to have running time at least nΩ(log (1/ǫ)) , i.e., quasipolynomial in 1/ǫ. It is natural to ask whether this quasi-polynomial lower bound can be improved to, say,
exponential, e.g., nΩ(1/ǫ) . We show
√ that the lower bound of Theorem 1.2 is qualitatively tight. We design
an (SQ) algorithm that uses Oǫ (n log(1/ǫ) ) SQ queries of inverse quasi-polynomial precision. Moreover,
we can turn this SQ algorithm into an algorithm in the sampling oracle model with similar complexity.
Specifically, we show (see Theorem 8.7 and Corollary 8.8):
Theorem 1.3 (SQ Algorithm for Robust Learning of Unknown Mean Gaussian). Let D be a distribution
on
√
O( log(1/ǫ))
n
n
R such that dTV (D, N (µ, I))
)
√ ≤ ǫ for some µ ∈ R . There is an SQ algorithm that uses Oǫ (n
O( log(1/ǫ))
n
, and outputs µ
e ∈ R such that dTV (N (e
µ, I), N (µ, I)) ≤ O(ǫ). The
SQ’s to D of precision ǫ/n
SQ algorithm can be turned into an algorithm√(in the sample model) with the same error guarantee that has
sample complexity and running time Oǫ (nO( log(1/ǫ)) ).
Theorems 1.2 and 1.3 give a qualitatively tight characterization of the complexity of robustly learning
an unknown mean Gaussian in the standard agnostic model, where the noisy distribution D is such that
dTV (D, N (µ, I)) ≤ ǫ. Equivalently, D satisfies (1 − ǫ1 )D + ǫ1 N1 = (1 − ǫ2 )N (µ, I) + ǫ2 N2 , where
N1 , N2 are unknown distributions and ǫ1 + ǫ2 ≤ ǫ. A weaker error model, known as Huber’s contamination
model in the statistics literature [Hub64, HRRS86, HR09], prescribes that the noisy distribution D is of the
form D = (1−ǫ)N (µ, I)+ǫN , where N is an unknown distribution. Intuitively, the difference is that in the
former model the adversary is allowed to subtract good samples and add corrupted ones, while in the latter
the adversary is only allowed to add corrupted ones. We note that the lower bound of Theorem 1.2 does
not apply in Huber’s contamination model. This holds for a reason: Concurrent work [DKK+ 17b] gives
a poly(n/ǫ) time algorithm with O(ǫ) error for robustly learning N (µ, I) in Huber’s model. Hence, as a
corollary, we establish a computational separation between these two models of corruptions. We provide an
intuitive justification in Section 1.3.
Our second super-polynomial SQ lower bound is for the problem of robustly learning a zero-mean
unknown covariance Gaussian with respect to the spectral norm. Specifically, we show (see Theorem 5.12
for a detailed statement):
Theorem 1.4 (SQ Lower Bound for Robust Learning of Unknown Covariance Gaussian). Let ǫ > 0, 0 <
c ≤ 1, and n ≥ poly(log(1/ǫ)). Any SQ algorithm that, given access to an ǫ-corrupted n-dimensional
e with kΣ
e − Σk2 ≤ O(ǫ log(1/ǫ)1−c ), requires at least
Gaussian N (0, Σ), with I/2  Σ  2I, returns Σ
c/4
c/4
Ω(1)
≥ nΩ(log(1/ǫ) ) queries to an SQ oracle of precision n−Ω(log(1/ǫ) ) .
2n

Similarly, Theorem 1.4 shows a super-polynomial gap between the information-theoretic sample complexity and the complexity of any SQ algorithm for this problem. As mentioned in the introduction, O(n/ǫ2 )
samples information-theoretically suffice to agnostically learn the covariance to within spectral error O(ǫ).
Second, the robust learning algorithm of [DKK+ 16] runs in poly(n, 1/ǫ) time, can be implemented in the
SQ model, and achieves error O(ǫ log(1/ǫ)) in Mahalanobis distance (hence, also in spectral norm). Again,
the immediate corollary is that the O(ǫ log(1/ǫ)) error guarantee of the [DKK+ 16] algorithm is best possible among all polynomial-time SQ algorithms. The lower bound of Theorem 1.4 does not apply in Huber’s
contamination model. This holds for a reason: [DKK+ 17b] gives a poly(n) · 2poly log(1/ǫ) time algorithm
with O(ǫ) error in Huber’s model.
7

Statistical–Computational Tradeoffs for SQ algorithms. Our next SQ lower bounds establish nearly
quadratic statistical–computational tradeoffs for robust spectral covariance estimation and robust sparse
mean estimation. We note that both these lower bounds also hold in Huber’s contamination model. For
the former problem, we show (see Theorem 6.1 for the formal statement):
Theorem 1.5. Let 0 < c < 1/6, and n sufficiently large. Any SQ algorithm that, given access to an
ǫ-corrupted N (0, Σ), where ǫ ≤ c/ ln(n) for kΣk2 ≤ poly(n/ǫ), and returns Σ̃ with Σ̃/2  Σ  2Σ̃,
c/3
requires at least 2Ω(n ) queries to an SQ oracle of precision γ = O(n)−(1−5c/2) .
We note that, in order to simulate a single query of the above precision, we need to draw Ω(1/γ 2 ) =
Ω(n2−5c ) samples from our distribution. Roughly speaking, Theorem 1.5 shows that if an SQ algorithm
c/3
uses less than this many samples, then it needs to run in 2Ω(n ) time. This suggests a nearly-quadratic
statistical-computational tradeoff for this problem.
For robust sparse mean estimation we show (see Theorem 6.6 for the detailed statement):
Theorem 1.6. Fix any 0 < c < 1 and let n ≥ 8k2 . Any SQ algorithm that, given access to an ǫ-corrupted
N (µ, I), where ǫ = k−c/4 , and µ ∈ Rn is promised to be k-sparse with kµk2 = 1, and outputs a hypothesis
c
vector µ
b satisfying kb
µ − µk2 ≤ 1/2, requires at least nΩ(ck ) queries to an SQ oracle of precision γ =
O(k)3c/2−1 .

Similarly, to simulate a single query of the above precision, we need to draw Ω(1/γ 2 ) = Ω(k2−3c )
samples from our distribution. Hence, any SQ algorithm that uses this many samples requires runtime at
c
least nΩ(ck ) . This suggests a nearly-quadratic statistical-computational tradeoff for this problem.
Sample Complexity Lower Bounds for High-Dimensional Testing. We now turn to our informationtheoretic lower bounds on the sample complexity of the corresponding high-dimensional testing problems.
For the robust Gaussian mean testing problem in Huber’s contamination model, we show (see Theorem 7.5
for a more detailed statement):

Theorem 1.7 (Sample Complexity Lower Bound for Robust Testing of Unknown Mean Gaussian). Fix
ǫ > 0. Any algorithm with sample access to a distribution D on Rn which satisfies either (a) D = N (0, I)
or (b) D is a δ-noisy N (µ, I), and kµk2 ≥ ǫ, and distinguishes between the two cases with probability 2/3
requires (i) Ω(n) samples if δ = ǫ/100, (ii) Ω(n1−c ) samples if δ = ǫ/nc/4 , for any constant 0 < c < 1.
As stated in the Introduction, without the robustness requirement, for any constant ǫ > 0, the Gaussian
√
mean testing problem can be solved with Oǫ ( n) samples. Hence, the conceptual message of Theorem 1.7
is that robustness makes the Gaussian mean testing problem information-theoretically harder. In particular,
the sample complexity of robust testing is essentially the same as that of the corresponding learning problem.
Theorem 1.7 can be viewed as a surprising fact because it implies that the effect of robustness can be very
different for testing versus learning of the same distribution family. Indeed, recall that the sample complexity
of robustly learning an ǫ-corrupted unknown mean Gaussian, up to error O(ǫ), is O(n/ǫ2 ) – i.e., the same
as in the noiseless case.
As a final application of our techniques, we show a sample complexity lower bound for the problem of
testing whether a spherical GMM is close to a Gaussian (see Theorem 7.6 for the detailed statement):
Theorem 1.8 (Sample Complexity Lower Bound for Testing a GMM). Any algorithm with sample access to
a distribution D on Rn which satisfies either (a) D = N (0, I), or (b) D = (1/2)N (µ1 , I)+ (1/2)N (µ2 , I),
such that dTV (D, N (0, I)) ≥ ǫ, and distinguishes between the two cases with probability at least 2/3
requires Ω(n/ǫ2 ) samples.
Similarly, the sample lower bound of Theorem 1.8 is optimal, up to constant factors, and coincides with
the sample complexity of learning the underlying distribution.
8

1.3 Our Approach and Techniques In this section, we provide a detailed outline of our approach and
techniques. The structure of this section is as follows: We start by describing our Generic Lower Bound
Construction, followed by our main applications to the problems of Learning GMMs and Robustly Learning
an Unknown Gaussian. We continue with our applications to statistical–computational tradeoffs. We then
explain how our generic technique can be used to obtain our Sample Complexity Testing Lower Bounds,
which rely on essentially the same hard instances as our SQ lower bounds. We conclude with a sketch of
our new (SQ) Algorithm for Robustly Learning an Unknown Mean Gaussian to optimal accuracy.
Generic Lower Bound Construction. The main idea of our lower bound construction is quite simple: We
construct a family of distributions D that are standard Gaussians in all but one direction, but are somewhat
different in the remaining direction (Definition 3.1). Effectively, we are hiding the interesting information
about our distributions in this unknown choice of direction. By exploiting the simple fact that it is possible
to find exponentially many nearly-orthogonal directions (Lemma 3.7), we are able to show that any SQ
algorithm with insufficient precision needs many queries in order to learn an unknown distribution from D.
To prove our generic SQ lower bound, we need to bound from below the SQ-dimension of our hard
family of distributions D. Roughly speaking, the SQ-dimension of a distribution family (Definition 2.11)
corresponds to the number of nearly uncorrelated distributions (with respect to some fixed distribution) in
the family (see Definitions 2.9 and 2.10). It is known that a lower bound on the SQ-dimension implies a
corresponding lower bound on the number and precision of queries of any SQ algorithm (see Lemma 2.12).
More concretely, our hard families of distributions are constructed as follows: Given a distribution A on
the real-line, we define a family of high-dimensional distributions Pv (x), for v ∈ Sn a unit n-dimensional
vector. The distribution Pv gives a copy of A in the v-direction, while being an independent standard
Gaussian in the orthogonal directions (Definition 3.1). Our hard family will be the set D = {Pv | v ∈ Sn }.
For the sake of the intuition, we make two observations: (1) If A and N (0, 1) have substantially different
moments of degree at most m, for some m, then Pv and N (0, I) can be easily distinguished by comparing
their mth -order moment tensors. Since these tensors can be approximated in roughly nm queries (and
time), the aforementioned lower bound construction would necessarily fail unless the low-order moments
of A match the corresponding low-order moments of G. We show that, aside from a few mild technical
conditions (see Condition 3.2), this moment-matching condition is essentially sufficient for our purposes.
If the degree at most m tensors agree, we need to approximate tensors of degree m + 1. Intuitively, in
order to extract useful information from these higher degree tensors, one needs to approximate essentially
all of the nm+1 many such tensor entries. (2) A natural approach to distinguish between Pv and N (0, I)
would be via random projections. As a critical component of our proof, we show (see Lemma 3.5) that a
random projection of Pv will be exponentially close to N (0, 1) with high probability. Therefore, a random
projection-based algorithm would require exponentially many random directions until it found a good one.
We now proceed with a somewhat more technical description of our proof. To bound from below the
SQ-dimension of our hard family of distributions, we proceed
R as follows: The definition of the pairwise
correlation (Definition 2.9) implies we need to show that Pv Pv′ /G ≈ 1, where G ∼ N (0, I) is the
Gaussian measure, for any pair of unit vectors v, v ′ that are nearly orthogonal. To prove this fact, we make
essential use of the Gaussian (Ornstein–Uhlenbeck) noise operator and its properties (see, e.g., [O’D14]).
We explain this connection in the following paragraph.
By construction of the distributions Pv , Pv′ , it follows that in the directions perpendicular to both v and
v ′ , the relevant factors integrate to 1. Letting y = v · x and z = v ′ · x and letting y ′ , z ′ be the orthogonal
directions to y and z, we need to consider the integral
Z
A(y)A(z)G(y ′ )G(z ′ )/G(x) .

9

Fixing y and integrating over the orthogonal direction, we get
Z
Z
A(y)/G(y) A(z)G(z ′ )dy ′ .
Now, if v and v ′ are (exactly) orthogonal, z = y ′ and the inner integral equals G(y). When this is not the
case, the A(z) term is not quite vertical and the G(z ′ ) term not quite horizontal, so instead what we get is
only nearly Gaussian. In general, the inner integral is equal to
Uv·v′ A(y) ,
√
where Ut is the member of the Ornstein–Uhlenbeck semigroup, Ut f (z) = E[f (tz + 1 − t2 G)]. We show
that this quantity is close to a Gaussian, when v · v ′ is close to 0 (see Lemma 3.4).
The core idea of the analysis relies on the fact that Ut A is a smeared out version of A. As such, it
only retains the most prominent features of A, namely its low-order moments. In fact, we are able to show
that if A and G agree in their first m moments, then Ut A is Om (tm )-close to a Gaussian (see Lemma 3.5),
and thus the integral in question is Om ((|v · v ′ |)m )-close to 1. This intuition is borne out in a particularly
clean way by writing A/G in the basis of Hermite polynomials. The moment-matching condition implies
that the decomposition involves none of the Hermite polynomials of degrees 1 through m. However, the
Ornstein–Uhlenbeck operator, Ut , is diagonalized by the basis Hi G with eigenvalue ti . Thus, if A − G can
be written in this basis with no terms of degree less than m, applying Ut decreases the size of the function
by a multiple of approximately tm .
So far, we have provided a proof sketch of the following statement (Lemma 3.4): When two unit vectors
′
v, v are nearly orthogonal, then the distributions Pv , Pv′ are nearly uncorrelated. Since, for 0 < c <
c
1/2, we can pack 2Ω(n ) unit vectors v onto the sphere so that their pairwise inner products are at most
nc−1/2 (Lemma 3.7), we obtain an SQ-dimension lower bound of our hard family. In particular, to learn
c
the distribution Pv , for unknown v, any SQ algorithm requires either 2Ω(n ) queries or queries of accuracy
better than O(n)(m+1)(c−1/2) (Proposition 3.3). This completes the proof sketch of our generic construction.
In our main applications, we construct one-dimensional distributions A satisfying the necessary momentmatching conditions for m taken to be super-constant, thus obtaining super-polynomial SQ lower bounds.
For our quadratic statistical–computational tradeoffs, we match a constant number of moments. In the following paragraphs, we explain how we apply our framework to bound the SQ dimension for: (i) learning
k-GMMs to constant accuracy, (ii) robustly learning an ǫ-corrupted Gaussian to accuracy O(ǫ), and (iii)
robustly estimating a Gaussian covariance within constant spectral error and robustly estimating a sparse
Gaussian mean to constant ℓ2 -error. In all cases, we construct a distribution A on the real-line that satisfies the necessary moment-matching conditions such that the family D = {Pv | v ∈ Sn } belongs in the
appropriate class, e.g., is a k-GMM for (i), an ǫ-corrupted Gaussian for (ii), etc.
SQ Lower Bound for Learning k-GMMs. The properties of our one-dimensional distribution A are summarized in Proposition 4.2. Specifically, we construct a distribution A on the real line that is a k-mixture
of one-dimensional “skinny” Gaussians, Ai , that agrees with N (0, 1) on the first m = 2k − 1 moments
(condition (i)). For technical reasons, we require that the chi-squared divergence of A to N (0, 1) is bounded
from above by an appropriate quantity (condition (iv)). The Gaussian components, Ai , have the same variance and appropriately bounded means (condition (ii)). We can also guarantee that the components Ai are
almost non-overlapping (condition (iii)). This implies that the corresponding high-dimensional distributions
Pv , P′v will be at total variation distance close to 1 from each other when the directions v, v ′ are nearly
orthogonal, and moreover their means will be sufficiently separated.
To establish the existence of a distribution A with the above properties, we proceed in two steps: First,
we construct (Lemma 4.3) a discrete one-dimensional distribution B supported on k points, lying in an
10

√
O( k) length interval, that agrees with N (0, 1) on the first k moments. The existence of such a distribution
B essentially follows from standard tools on Gauss-Hermite quadrature. The distribution A is then obtained
(Corollary 4.4) by adding a zero-mean skinny Gaussian to an appropriately rescaled version of B. Additional
technical work (Lemmas 4.5 and 4.6) gives the other conditions.
Our family of hard high-dimensional instances will consist of GMMs that look like almost non-overlapping
“parallel pancakes” and is reminiscent of the family of instances considered in Brubaker and Vempala [BV08].
For the case of k = 2, consider a 2-GMM where both components have the same covariance that is far from
spherical, the vector between the means is parallel to the unit eigenvector with smallest eigenvalue, and
the distance between the means is a large multiple of the standard deviation in this direction (but a small
multiple of that in the orthogonal direction). This family of instances was considered in [BV08], who gave
an efficient spectral algorithm to learn them.
Our lower bound construction can be thought of as k “parallel pancakes” in which the means lie in a
one-dimensional subspace, corresponding to the smallest eigenvalue of the identical covariance matrices of
the components. All n − 1 orthogonal directions will have an eigenvalue of 1, which is much larger than
the smallest eigenvalue. In other words, for each unit vector v, the k-GMM Pv will consist of k “skinny”
Gaussians whose mean vectors all lie in the direction of v. Moreover, each pair of components
will have total
√
variation distance very close to 1 and their mean vectors are separated by Ω(1/ k). We emphasize once
more that our hard family of instances is learnable with O(k log n) samples – both for density estimation
and parameter estimation. On the other hand, any SQ learning algorithm for the family requires nΩ(k) time.
SQ Lower Bounds for Robustly Learning Unknown Gaussian. In the agnostic model, there are two
types of adversarial noise to handle: subtractive noise – corresponding to the good samples removed by
the adversary – and additive noise – corresponding to the bad points added by the adversary. The approach
of [DKK+ 16] does not do anything
p to address subtractive noise, but shows that this type of noise can incur
“small” error, e.g., at most O(ǫ log(1/ǫ)) for the case of unknown mean. For additive noise, [DKK+ 16]
uses an iterative spectral algorithm to filter out outliers.
For concreteness, let us consider the case of robustly learning N (µ, I). Intuitively, achieving error O(ǫ)
in the agnostic model is hard for the following reason: the two types of noise can collude so that the first few
moments of the p
corrupted distribution are indistinguishable from those of a Gaussian whose mean vector
has distance Ω(ǫ log(1/ǫ)) from the true mean.
To formalize this intuition, for our robust SQ learning lower bound, we construct a distribution A on
the real line that agrees with N (0, 1) on the first m = Ω(log1/4 (1/ǫ)) moments and is ǫ/100-close in
total variation distance to G′ = N (ǫ, 1) (see Proposition 5.2). We achieve this by taking A to be the
Gaussian N (ǫ, 1) outside its effective support, while in the effective support we add an appropriate degreem univariate polynomial p satisfying the appropriate moment conditions. By expressing this polynomial as
a linear combination of appropriately scaled Legendre polynomials, we can prove that its L1 and L∞ norms
within the effective support of G′ are much smaller than ǫ (see Lemma 5.6). This result is then used to
bound from above the distance of A from G′ , which gives our SQ lower bound.
We use a similar technique to prove our SQ lower bound for robust covariance estimation in spectral
norm. Specifically, we construct a distribution A that agrees with N (0, 1) on the first m = Ω(log(1/ǫ))
moments and is ǫ/100-close in total variation distance to G′ = N (0, (1 − δ)2 ), for some δ = O(ǫ) (see
Proposition 5.13). We similarly take A to be the Gaussian G′ outside its effective support, while in the
effective support we add an appropriate degree-m univariate polynomial p satisfying the appropriate moment
conditions. The analysis proceeds similarly as above.
Statistical–Computational Tradeoffs for SQ algorithms. For robust covariance estimation in spectral
norm, our one-dimensional distribution is selected to be A = (1 − ǫ)N (0, σ) + ǫN1 , where N1 is a mixture
of 2 unit-variance Gaussians with opposite means. By selecting σ appropriately, we can have A match
the first 3 moments of N (0, 1), see Theorem 6.1. For robust sparse mean estimation, it suffices to take
11

A = (1 − δ)N (ǫ, 1) + δN1 , where N1 is a unit-variance Gaussian selected so that E[A] = 0. An important
aspect of both these constructions is that the chi-squared distance χ2 (A, N (0, 1)) needs to be as small as
possible. Indeed, since we only match a small number of moments, our bound on χ2 (A, N (0, 1)) crucially
affects the accuracy of our SQ queries (Proposition 3.3).
Sample Complexity Testing Lower Bounds. Our sample complexity lower bounds follow from standard
information-theoretic
arguments, and rely on the same lower bound instances and correlation bounds (i.e.,
R
bounds on Pv Pv′ /G) established in our SQ lower bounds. In particular, we consider the problem of distinguishing between the distribution G ∼ N (0, I) and the distribution Pv for a randomly chosen unit vector
v ∈ Sn using N independent samples. Let G⊗N denote the distribution on N independent samples from G,
and P⊗N
the distribution obtained by picking a random v and then taking N independent samples from Pv .
v
If it is possible to reliably distinguish between these cases, itR must be
QNthe case that the chi-squared diver⊗N
⊗N
gence χ(Pv , G ) is substantially larger than 1. This is v,v′ ,xi i=1 Pv (xi )Pv′ (xi )/G(xi )dvdv ′ dxi .
Note that after fixing v and v ′ the above integral separates as a product, giving
Z

v,v′

Z

Pv (x)Pv′ (x)/G(x)dx

N

dvdv ′ .

(1)

Note that the inner integral was bounded from above by roughly (1 + (v · v ′ )m ). A careful analysis of the
distribution of the angle between two random unit vectors allows us to show that, unless N = Ω(n), the
chi-squared divergence is close to 1, and thus that this testing problem is impossible.
Algorithm for Robustly Learning Unknown Mean Gaussian. We give an SQ algorithm with O(ǫ)error for robustly learning an unknown mean Gaussian, showing that our corresponding SQ lower bound is
qualitatively tight. Our algorithm builds on the filter technique of [DKK+ 16], generalizing it to the more
involved setting of higher-order tensors.
As is suggested by our SQ lower bounds, the obstacle to learning the mean robustly, is that there are
ǫ-noisy Gaussians that are Ω(ǫ)-far in variation distance from a target Gaussian G, and yet match G in
all of their first O(log1/4 (1/ǫ)) moments. For our algorithm to circumvent this difficulty, it will need to
approximate all of the tth -order tensors for t ≤ k = Ω(log1/4 (1/ǫ)). Note that this already requires nk SQ
queries.
The first thing we will need to show is that k moments suffice, for an appropriate parameter k. Because of our lower bound construction, we know that k needs to be at least Ω(log1/4 (1/ǫ)). We show
that k = O(log1/2 (1/ǫ)) suffices. Specifically, we prove a one-dimensional moment-matching lemma
(Lemma 8.1) establishing the following: If an ǫ-noisy one-dimensional Gaussian approximately matches a
reference Gaussian G in all of its first k moments, where k = Θ(log1/2 (1/ǫ)) (i.e., quadratically larger than
our lower bound), then it must be O(ǫ)-close to G in variation distance. We note that it suffices to prove this
statement in the one-dimensional case, as we can just project onto the line between the means.
We now proceed to describe ourpalgorithm: Using the basic filter algorithm from [DKK+ 16], we start by
learning the true mean to error O(ǫ log(1/ǫ)). By translating, we can assume that the mean is this close to
0. We need to robustly approximate the low-order moments of our target Gaussian G′ . This is complicated
by the fact that even a small fraction of errors can have a huge impact on the moments of the distribution.
However, any large errors are easily detectable. In particular, if any tth moment tensor differs substantially
from that of the standard Gaussian, it will necessarily imply the presence of errors. In particular, it will
allow us to construct a polynomial p so that E[p(X)] − E[p(G′ )] (where X is a noisy version of G′ ) is much
larger than ǫkp(G′ )k2 . If this is the case, then many of our errors, x, must have p(x) very far from the mean.
By standard concentration inequalities, this will allow us to identify these points as almost certainly being
errors. This in turn lets us build a filter to clean-up our distribution X, making it closer to G′ .

12

Repeatedly applying filters as necessary, we can reduce to the case where the higher-order moments of
X are close to the higher-order moments of G. This will tell us that, in almost all directions, the first k
moments of X match the corresponding moments of G. By our moment-matching lemma, this will imply
that the mean of G′ is close to 0 in these directions. We will then only need to approximate the mean of
the projection of G′ onto the low-dimensional subspace V in which these moments fail to match. This
approximation can be done in a brute-force manner (in time exponential in dim(V ), which is still relatively
small), completing the description of the algorithm.
1.4 Related Work This work studies learning and testing high-dimensional structured distributions. Distribution learning and testing are two of the most fundamental inference tasks in statistics with a rich history
(see, e.g., [NP33, BBBB72, DG85, Sil86, Sco92, DL01, LR05]) that date back to Karl Pearson. The main
criteria to evaluate the performance of an estimator are its sample complexity and its computational complexity. Despite intensive investigation for several decades by different communities, the (sample and/or
computational) complexity of many learning and testing problems is still not well-understood, even for
some surprisingly simple high-dimensional settings. In the past few decades, a long line of work within
TCS [KMR+ 94, Das99, FM99, AK01, VW02, CGG02, MR05, BV08, KMV10, MV10, BS10, DDS12a,
DDS12b, CDSS13, DDO+ 13, CDSS14a, CDSS14b, ADLS17, DDS15, DDKT16, DKS16b, DKS16a] has
focused on designing efficient estimators in a variety of settings. We have already mentioned the most
relevant references for the specific questions we consider in Section 1.1.
With respect to computational lower bounds for unsupervised estimation problems, the most relevant
references are the works [FGR+ 13, FPV15, KV16] that show SQ lower bounds for the planted clique and
related planted-like problems. It should be noted that, beyond the fact that we also use the concept of
SQ dimension, our techniques are entirely different than theirs. Prior work by Feldman, O’Donnell, and
Servedio [FOS08] implicitly showed an SQ lower bound of nΩ(log k) for the problem of learning k-mixtures
of product distributions over {0, 1}n . This was obtained by a straightforward reduction from the problem of
learning k-leaf decision trees over n Boolean variables. Our lower bound construction for learning GMMs
is entirely different from [FOS08] that relied on the obvious combinatorial structure of the discrete setting.
A related line of work gives statistical-computational tradeoffs for sparse PCA [BR13a, BR13b, MW15,
WBS16a], based on various computational hardness assumptions. These results are of similar flavor as our
statistical–computational tradeoffs for SQ algorithms (Theorems 1.5 and 1.6). An important difference between these tradeoffs and the super-polynomial SQ lower bounds we prove in this paper (Theorems 1.1, 1.2,
and 1.4) is that the aforementioned sparse problems are known to be tractable if we increase the sample size
by a quadratic factor beyond the information-theoretic limit. In contrast, our main SQ lower bound results
establish a super-polynomial gap between the information-theoretic limit and the computational complexity
of any SQ algorithm.
Finally, we remark that in the supervised setting of PAC learning Boolean functions, a number of hardness results are known based on various complexity assumptions, see, e.g., [KKMS08, KS06, FGKP06,
KK14, DLS14, Dan16] for the problems of learning halfspaces and learning intersections thereof.
1.5 Discussion and Future Directions The main contribution of this paper is a technique that gives essentially tight SQ lower bounds for a number of fundamental high-dimensional learning problems, including
learning GMMs and robustly learning a single Gaussian. To the best of our knowledge, these are the first
such lower bounds for high-dimensional distribution learning problems in the continuous setting. As a
corollary, we provide a rigorous explanation of the observed (super-polynomial) gap between the sample
complexity of these problems and the runtime of the best known algorithms.
Our work naturally raises a number of interesting future directions. A natural open problem is to extend
our lower bound technique to broader families of high-dimensional distributions. More concretely, is there
a kω(1) poly(n) SQ lower bound for learning k-mixtures of n-dimensional spherical Gaussians? Note that
13

our nΩ(k) lower bound does not apply for the spherical case, as it crucially exploits the structure of the
covariance matrices. In fact, faster learning algorithms for the spherical case are known [SOAJ14], albeit
with exponential dependence on the number k of components. More broadly, can we extend our techniques
to other families of structured high-dimensional distributions (e.g., mixtures of other distribution families)?
1.6 Organization The structure of this paper is as follows: In Section 2, we introduce basic notation,
definitions, and a number of useful facts that will be required throughout the paper. Our SQ lower bounds
are established in Sections 3–6. Specifically, in Section 3, we give our generic high-dimensional SQ lower
bound construction, assuming the existence of a one-dimensional density satisfying the necessary moment
conditions. In Sections 4, 5, and 6, we construct the appropriate one-dimensional densities, thereby establishing our SQ lower bounds. Specifically, Sections 4 and 5 give our super-polynomial SQ lower bounds for
the problems of learning GMMs and robustly learning an unknown Gaussian. Section 6 gives our quadratic
statistical–computational tradeoffs (for SQ algorithms) for the problems of robust covariance estimation in
spectral norm and robust sparse mean estimation. Section 7 gives our (information-theoretic) sample complexity lower bounds for high-dimensional testing. Finally, in Section 8 we present our (SQ) algorithm for
robustly learning an unknown mean Gaussian with optimal accuracy, whose runtime qualitatively matches
our SQ lower bound from Section 5.
Acknowledgements. This project evolved over a number of years. We would like to thank Vitaly Feldman
for answering numerous questions about the Statistical Query model; Andy Drucker for useful discussions
on Question 1.1; Anup B. Rao for asking a question that motivated Theorem 6.1; Weihao Kong and Gregory
Valiant for useful discussions on Question 1.4; and Ankur Moitra and Eric Price for feedback on a previous
version of this paper.

2 Definitions and Preliminaries
2.1 Notation and Basic Definitions For n ∈ Z+ , we denote by [n] the set {1, 2, . . . , n}. We will denote
by Sn the Euclidean unit sphere in Rn . If v is a vector, we will let kvk2 denote its Euclidean norm. If M is
a matrix, we will let kM k2 denote its spectral norm, and kM kF denote its Frobenius norm.
Our basic object of study is the Gaussian (or Normal) distribution and finite mixtures of Gaussians:
Definition 2.1. The n-dimensional Gaussian distribution N (µ, Σ) with mean vector µ ∈ Rn and covariance
matrix Σ ∈ Rn×n is the distribution with probability density function


1
T −1
−n/2
1/2
f (x) = (2π)
det(Σ) exp − (x − µ) Σ (x − µ) .
2
Definition 2.2. An n-dimensional k-mixture
of Gaussians (k-GMM) is a distribution onP
Rn with probability
Pk
density function defined by F (x) = j=1 wj N (µj , Σj ), where wj ≥ 0, for all j, and kj=1 wj = 1.

Throughout the paper, we will make extensive use of the pdf of the standard one-dimensional Gaussian
N (0, 1), which we will denote by G(x).
Definition 2.3. The total variation distance between two distributions (with probability density functions)
R
def
P, Q : Rn → R+ is defined to be dTV (P, Q) = (1/2) · kP − Qk1 = (1/2) · x∈Rn |P(x) − Q(x)|dx. The
R
def R
χ2 -divergence of P, Q is χ2 (P, Q) = x∈Rn (P(x) − Q(x))2 /Q(x)dx = x∈Rn P2 (x)/Q(x)dx − 1.

2.2 Formal Problem Definitions We record here the formal definitions of the problems that we study.
Our first problem of interest is learning a mixture of k arbitrary high-dimensional Gaussians:

14

Definition 2.4 (Density Estimation/Proper Learning of GMMs). Let Gn,k be the family of n-dimensional
k-GMMs. The problem of density estimation for Gn,k is the following: Given ǫ > 0 and sample access to
an unknown P ∈ Gn,k , with probability 9/10, output a hypothesis distribution H such that dTV (H, P) ≤ ǫ.
The problem of proper learning for Gn,k is the same with the additional requirement that H ∈ Gn,k .
We also
Pkconsider the more challenging task of parameter estimation: Given samples from a distribution
P = i=1 Gi , where Gi is a weighted Gaussian, the goal is to learn a distribution Q that can be written as
P
Q = ki=1 Hi , where Hi is a weighted Gaussian and dTV (Hi , Gi ) is small for all i.
Our next question is the problem of robustly learning a Gaussian in the standard agnostic model:

Definition 2.5 (Robust Learning of Unknown Gaussian). The problem of robust (agnostic) learning of an
unknown Gaussian is the following: Given ǫ > 0, and sample access to an unknown distribution D with
dTV (D, N (µ, Σ)) ≤ ǫ, for some µ ∈ Rn and Σ ∈ Rn×n , output a hypothesis distribution H such that, with
probability at least 9/10, it holds that dTV (H, N (µ, Σ)) = f (ǫ), for some f : R+ → R+ .
Our SQ lower bounds apply to two special cases of this problem: when µ is unknown and Σ = I, and
when µ = 0 and Σ is unknown. For the latter case, our lower bound applies even for learning with respect
to the spectral norm (which is weaker than approximation in variation distance).
We now define the problem of robustly testing a Gaussian in Huber’s model. We remind the reader that
our tight sample complexity lower bound applies to this weaker model as well.
Definition 2.6 (Robust Testing of Unknown Mean Gaussian). The problem of robust testing of an unknown
mean Gaussian is the following: Given ǫ > 0, 0 < δ < ǫ/2 and sample access to an unknown distribution
D over Rn with the promise that one of the following two cases is satisfied: (i) D = N (0, I), or (ii)
D = (1 − δ)N (µ, I) + δN1 , where N1 is an unknown noise distribution and kµk2 ≥ ǫ, the goal is to
correctly distinguish between the two cases with confidence probability 2/3.
Finally, our problem of testing GMMs is the following:
Definition 2.7 (Testing Spherical GMMs). The problem of testing of a spherical GMM is the following:
Given ǫ > 0 and sample access to an unknown distribution D with the promise that one of the following two
cases is satisfied: (a) D = N (0, I), or (b) D is a 2-GMM w1 N (µ1 , I) + w2 N (µ2 , I) in Rn which satisfies
dTV (D, N (0, I)) ≥ ǫ, distinguish between the two cases with probability at least 2/3.
2.3 Basics on Statistical Query Algorithms over Distributions We begin by recording the necessary
definitions of Statistical algorithms for problems over distributions. All the definitions and facts in this
section are from [FGR+ 13]. We start by defining a general search problem over distributions.
Definition 2.8 (Search problems over distributions). Let D be a set of distributions over Rn , let F be a set
of solutions and Z : D → 2F be a map from a distribution D ∈ D to a subset of solutions Z(D) ⊆ F that
are defined to be valid solutions for D. The distributional search problem Z over D and F is to find a valid
solution f ∈ Z(D) given access to (an oracle or samples from) an unknown D ∈ D.
For general search problems over a distribution, we define SQ algorithms as algorithms that do not see
samples from the distribution but instead have access to an SQ oracle. We consider two types of SQ oracles
from the literature.
1. STAT(τ ): For a tolerance parameter τ > 0 and any bounded function f : Rn → [−1, 1], STAT(τ )
returns a value v ∈ [Ex∼D [f (x)] − τ, Ex∼D [f (x)] + τ ].

15

2. VSTAT(t): For a sample size parameter t > 0 and any bounded function f : Rn→ [0, 1], VSTAT(t)

q
Varx∼D [f (x)]
1
returns a value v ∈ [Ex∼D [f (x)] − τ, Ex∼D [f (x)] + τ ], where τ = max t ,
,
t
where Varx∼D [f (x)] = Ex∼D [f (x)] (1 − Ex∼D [f (x)]).

The first oracle was defined by Kearns [Kea98] and the second was introduced in [FGR+ 13]. These oracles
are known to be polynomially equivalent [FGR+ 13]. Also note that these oracles can return any value within
the given tolerance, and therefore can make adversarial choices.
The main technical tool that allows us to prove unconditional lower bounds on the complexity of SQ
algorithms is an appropriate notion of Statistical Query (SQ) dimension. Such a notion was defined in the
context of PAC learning of Boolean functions in [BFJ+ 94], and subsequently generalized to search problems
over distributions in [FGR+ 13]. We will require the simpler definition from Section 3 of that work that relies
on pairwise correlations:
Definition 2.9 (Pairwise Correlation). The pairwise correlation of two distributions with probability density
functions D1 , D2 : Rn → R+ with respect to a distribution with density D : Rn → R+ , where the support
def R
of D contains the supports of D1 and D2 , is defined as χD (D1 , D2 ) = Rn D1 (x)D2 (x)/D(x)dx − 1.

We remark that when D1 = D2 in the above definition, the pairwise correlation is identified with the
def R
between D1 and D, i.e., χ2 (D1 , D) = Rn D1 (x)2 /D(x)dx − 1.
We will also need the following definition:

χ2 -divergence

Definition 2.10. We say that a set of m distributions D = {D1 , . . . , Dm } over Rn is (γ, β)-correlated
relative to a distribution D over Rn if
(
γ if i 6= j
|χD (Di , Dj )| ≤
β if i = j.
We are now ready to define our notion of dimension:
Definition 2.11 (Statistical Query Dimension). For β, γ > 0, a search problem Z over a set of solutions
F, and a class of distributions D over Rn , let m be the maximum integer such that there exists a reference
distribution D over Rn and a finite set of distributions DD ⊆ D such that for any solution f ∈ F, Df =
DD \ Z −1 (f ) is (γ, β)-correlated relative to D and |Df | ≥ m. We define the statistical (query) dimension
with pairwise correlations (γ, β) of Z to be m and denote it by SD(Z, γ, β).
Our lower bounds proceed by bounding from below the statistical query dimension of the considered
distribution learning problems. The corresponding lower bounds on the complexity of SQ algorithms for
these problems are a corollary of the following result from [FGR+ 13]:
Lemma 2.12 (Corollary 3.12 in [FGR+ 13]). Let Z be a search problem over a set of solutions F and a
class of distributions D over Rn . For γ, β > 0, let s = SD(Z,
γ, β). For any γ ′ > 0, any SQ algorithm for
√
Z requires at least s · γ ′ /(β − γ) queries to the STAT( γ + γ ′ ) or VSTAT(1/(3(γ + γ ′ ))) oracles.

3 Statistical Query Lower Bounds: From One–Dimension to High–Dimensions
All our statistical query lower bounds are shown in two steps: We first construct a one-dimensional density
A satisfying certain technical conditions, and then use A to construct a high-dimensional distribution which
is Gaussian in all but one directions. The second step is the same for all the problems that we consider. To
formally define it, we require the following construction:

16

Definition 3.1 (High-Dimensional Hidden Direction Distribution). For a distribution A on the real line
with probability density function A(x) and a unit vector v ∈ Rn , consider the distribution over Rn with
probability density function

Pv (x) = A(v · x) exp −kx − (v · x)vk22 /2 /(2π)(n−1)/2 .

That is, Pv is the product distribution whose orthogonal projection onto the direction of v is A, and onto the
subspace perpendicular to v is the standard (n − 1)-dimensional normal distribution.
Suppose that we have constructed a one-dimensional distribution A satisfying the following condition:

Condition 3.2. Let m ∈ Z+ . The distribution A on R is such that (i) the first m moments of A agree with
the first m moments of N (0, 1), and (ii) χ2 (A, N (0, 1)) is finite.
Note that Condition 3.2-(ii) above implies that the distribution A has a probability density function
(pdf), which we will denote by A(x). We will henceforth blur the distinction between a distribution and its
pdf. The main result of this section is the following:
Proposition 3.3. Given a distribution A on R that satisfies Condition 3.2 for some m ∈ Z+ and a constant
0 < c < 1/2, consider the set of distributions {Pv }v∈Sn , for n ≥ mΩ(1/c) . For a given ǫ > 0, suppose that
dTV (Pv , Pv′ ) > 2ǫ whenever |v · v ′ | is at most 1/8. Then, any SQ algorithm which, given access to Pv (x)
c/2
for an unknown v ∈ Sn , outputs a hypothesis
Q with dTV (Q, Pv ) ≤ ǫ needs at least 2Ω(n ) ≥ nm+1 
p
queries to STAT(O(n)−(m+1)(1/4−c/2) χ2 (A, N (0, 1))) or to VSTAT O(n)(m+1)(1/2−c) /χ2 (A, N (0, 1)) .

An intuitive interpretation of the proposition is as follows: If we do not want our SQ algorithm to use
a number of queries exponential in nΩ(1) , then we would need Ω(n)Ω(m+1) samples to simulate a single
statistical query.
The rest of this section is devoted to the proof of Proposition 3.3. In Section 3.1, we prove a correlation
bound which is the main technical ingredient for the proof. In Section 3.2, we show a simple packing for
unit vectors over the sphere and put the pieces together to complete the proof.
3.1

Main Correlation Lemma The main technical result of this section is the following:

Lemma 3.4 (Correlation Lemma). Let m ∈ Z+ . If the distribution A over R agrees with the first m moments
of N (0, 1), then for all v, v ′ ∈ Rn , we have that
|χN (0,I) (Pv , Pv′ )| ≤ |v · v ′ |m+1 χ2 (A, N (0, 1)) .

(2)

Note that we may assume that χ2 (A, N (0, 1)) is finite, otherwise the lemma statement is trivial. Hence,
we can henceforth assume that Condition 3.2 is satisfied. In particular the distributions A, Pv and Pv′ all
have probability density functions.
To prove Lemma 3.4 we proceed as follows: We start by bounding the χ2 -divergence between the onedimensional projection of Pv onto v ′ and N (0, 1). As well as being a critical component towards the proof
of Lemma 3.4, this fact can be used to show that random projections of Pv are close to N (0, 1) with high
probability. Specifically, we show:
Lemma 3.5. Let Q be the distribution of v ′ · X, for X ∼ Pv . Then, we have that
χ2 (Q, N (0, 1)) ≤ (v · v ′ )2(m+1) χ2 (A, N (0, 1)) .

17

Proof. Let θ be the angle between v and v ′ . Let x, y be orthogonal coordinates for the plane spanned by v
and v ′ , with the x-axis in the v ′ direction. Note that Pv is a product of a distribution on this plane and a
standard Gaussian perpendicular to it. On this plane, Pv is a product of A and N (0, 1). Thus, we have that
Z
Pv (x)dx
Q(x) =
′
Zx:v ·x=x
A(x cos θ + y sin θ)G(x sin θ − y cos θ)dy .
=
y∈R

Let Uθ be the linear operator that maps f : R → R to
Z
f (x cos θ + y sin θ)G(x sin θ − y cos θ)dy ,
y∈R

so that Q = Uθ (A). We will show that we can expand A as a linear combination of eigenfunctions of Uθ . √
Let Hei (x) denote the i-th (probabilists’) Hermite polynomial. We note
R that the functions Hei (x)G(x)/ i!,
for i ∈ N, are orthonormal with respect to the inner product hf, gi = x∈R f (x)g(x)/G(x)dx. Indeed, by
p
using the fact that the Hermite functions, which can be written as He(x) RG(x), are a complete orthonormal family for L2 (R), we get that any function f : R → R such that R f (x)g(x)/G(x)dx
< ∞ is
R∞
almost everywhere equal to a linear combination of these Hei (x)G(x)/i!. Since −∞ A(x)2 /G(x)dx =
1 + χ2 (A, N (0, 1)) is finite, we can write
A(x) =

∞
X

√
ai Hei (x)G(x)/ i! .

(3)

i=0

Using the orthogonality of Hei (x) we can extract these coefficients, since
h
√ i Z
ai Hei (x)2 G(x)/i!dx = ai .
EX∼A Hei (X)/ i! =

(4)

R

Since A agrees with the first m moments of the standard Gaussian, for 0 ≤ i ≤ m, we have that
√
√
EX∼A [Hei (X)/ i!] = EX∼N (0,1) [Hei (X)/ i!] = δi,0 .
This implies that a0 = 1 and a1 , . . . , am = 0. Thus, we have
A(x) = G(x) +

∞
X

√
ai Hei (x)G(x)/ i! .

(5)

i=m+1

The orthonormality of these functions with respect to the inner product hf, gi =
allows us to express the χ2 -divergence in terms of these coefficients:
Z ∞
(A(x) − G(x))2 /G(x)dx
χ2 (A, N (0, 1)) =
=
=

Z

−∞
∞

−∞

∞
X

∞
X

i=m+1

R

R f (x)g(x)/G(x)dx

also

!2
√
ai Hei (x)G(x)/ i! /G(x)dx

a2i .

(6)

i=m+1

18

Now we consider the effect of Uθ on this orthogonal family. From the definition of Uθ , we have
Z ∞
Uθ (Hei G)(x) =
Hei (x cos θ + y sin θ)G(x sin θ − y cos θ)dy .
−∞

We will use the well-known fact that Hei (x)G(x) is an eigenfunction of Uθ :
Fact 3.6. We have that: Uθ (Hei G)(x) = cosi (θ)Hei (x)G(x) .

For completeness, we include a proof in Appendix D.
We can now use these eigenfunctions and eigenvalues with Equation (5) to get an expression for Uθ A:
Uθ A(x) = G(x) +

∞
X

√
ai cosi θHei (x)G(x)/ i! ,

(7)

i=m+1

which can be used to express its χ2 -divergence:

∞
X

χ2 (Uθ A, N (0, 1)) =

a2i cos2i θ

i=m+1

≤ cos2(m+1) θ
= cos

2(m+1)

∞
X

a2i

i=m+1
2

θ · χ (A, N (0, 1)) ,

(8)

where the last line uses (6). Recalling that Q = Uθ A and cos θ = v · v ′ , this completes the proof.
Proof of Lemma 3.4. We first show that the correlation between the high-dimensional densities Pv and Pv′
needed for Lemma 3.4 can be reduced to a one-dimensional correlation. Just as in the proof of Lemma
3.5, let θ = arccos(v · v ′ ) and let x, y be coordinates for the plane spanned by v and v ′ with the x-axis in
the v ′ direction. Each of Pv and Pv′ is a product of a distribution on this plane and a standard Gaussian
perpendicular to it. On this plane, they are both products of A and N (0, 1) with different rotations applied.
Thus, we have that
Z
Pv (x)Pv′ (x)/G(x)dx
χN (0,I) (Pv , Pv′ ) + 1 =
n
ZR∞ Z ∞
=
A(x)G(y)A(x cos θ + y sin θ)G(x sin θ − y cos θ)/G(x)G(y)dxdy
−∞ −∞
Z ∞
Z ∞
A(x cos θ + y sin θ)G(x sin θ − y cos θ)dydx
A(x)/G(x) ·
=
−∞

−∞

= 1 + χN (0,1) (A, Uθ A) .

Now we can bound from above this correlation in terms of the χ2 -divergences of both distributions from
N (0, 1), one of which we can bound using Lemma 3.5:
Z
|A(x) − G(x)||Uθ A(x) − G(x)|/G(x)dx
|χN (0,1) (A, Uθ A)| ≤
sZ
sRZ
≤

=

p

R

(A(x) − G(x))2 /G(x)dx ·

χ2 (A, N (0, 1))

≤ | cos

m+1

2

·

χ2 (U

R

(Uθ A(x) − G(x))2 /G(x)dx

θ A, N (0, 1))

(θ)| · χ (A, N (0, 1)) ,

where the first line follows by triangle inequality, the second inequality is Cauchy-Schwarz, and the last line
uses Lemma 3.5. The proof of Lemma 3.4 is now complete.
19

3.2 Proof of Proposition 3.3 We note that our distribution learning problem can be expressed as a search
problem in the sense of [FGR+ 13]. Consider the following search problem Z: given access to (an oracle or
samples from) Pv , for an unknown unit vector v, find a distribution f such that dTV (Pv , f ) ≤ ǫ. Thus, for
us, the set of solutions F is the set of all distributions on Rn , and D ⊆ F is the set of Pv for all unit vectors
v. For a unit vector v, Z(Pv ) is the set of all distributions f such that dTV (Pv , f ) ≤ ǫ. For a distribution f
on Rd , Z −1 (f ) is the set of Pv such that dTV (f, Pv ) ≤ ǫ.
To prove a lower bound on the statistical dimension, we will take D = N (0, I), D = {Pv : v ∈ Sn }
and construct a suitable finite set DD .
c

Lemma 3.7. For any 0 < c < 1/2, there is a set S of at least 2Ω(n ) unit vectors in Rn such that for each
pair of distinct v, v ′ ∈ S, it holds |v · v ′ | ≤ O(nc−1/2 ).
c

The lemma can be shown by a probabilistic argument. Specifically, if we take |S| = 2Ω(n ) unit vectors
uniformly at random from the unit sphere Sn , then the desired event happens with positive probability. The
proof is given in Appendix D.
c

Proof of Proposition 3.3. For a given constant 0 < c < 1/2, fix a set S of 2Ω(n ) unit vectors in Rn
satisfying the statement of Lemma 3.7. Let DD = {Pv | v ∈ S}. Then, by Lemma 3.4, we have that, for
v, v ′ ∈ S with v 6= v ′ , it holds
χN (0,I) (Pv , Pv′ ) ≤ |v · v ′ |m+1 χ2 (A, N (0, 1)) = Ω(n)−(m+1)(1/2−c) χ2 (A, N (0, 1)) .
If v = v ′ , then χN (0,I) (Pv , Pv ) = χ2 (Pv , N (0, I)) = χ2 (A, N (0, 1)). We thus have that, for
def

def

γ = Ω(n)−(m+1)(1/2−c) χ2 (A, N (0, 1)) and β = χ2 (A, N (0, 1)) ,
DD is (γ, β) correlated with respect to D = N (0, I).
Since dTV (Pv , Pv′ ) > 2ǫ for distinct v and v ′ in S, for any distribution f over Rn , we have that
Z −1 (f ) = {Pv : v ∈ S and dTV (Pv , f ) ≤ ǫ} has |Z −1 (f )| ≤ 1 using the triangle inequality for dTV .
We thus have that DD \ Z −1 (f ) is (γ, β) correlated with respect to D = N (0, I), and |DD \ Z −1 (f )| ≥
c
|S| − 1 = 2Ω(n ) . That is,
c
SD (Z, γ, β) ≥ 2Ω(n ) .
def

An application of Lemma 2.12 for γ ′ = γ = Ω(n)−(m+1)(1/2−c) χ2 (A, N (0, 1)), we obtain that any SQ
c
algorithm requires at least 2Ω(n ) n−(m+1)(1/2−c) calls to the




p
STAT O(n)−(m+1)(1/4−c/2) χ2 (A, N (0, 1)) or VSTAT O(n)(m+1)(1/2−c) /χ2 (A, N (0, 1))
oracle to solve Z. From our assumption that n ≥ mΩ(1/c) , it follows that n ≥ Ω((m + 1) log n)2/c , and
c/2
c/2
therefore 2Ω(n ) ≥ nm+1 . Hence, the total number of required queries is at least 2Ω(n ) ≥ nm+1 . This
completes the proof.

4 SQ Lower Bound for Learning Gaussian Mixtures
The main result of this section is the following:
Theorem 4.1. Fix 0 < ǫ < 1. Any SQ algorithm that given SQ access to a k-mixture P of n-dimensional
Gaussians, N (µi , Σi ), i ∈ [k], for n ≥ Ω(k 8 log(1/ǫ)3 ), which are promised to satisfy dTV (Pi , Pj ) ≥ 1−ǫ,
for all i 6= j, and moreover are such that




1/2
1/2
max max kµi − µj k2 , max kΣi k2
k
≤ poly(k, log(1/ǫ)) min 1/kΣ−1
,
2
i
i,j

i

i

20

and outputs a distribution Q with dTV (P, Q) ≤ 1/2, needs at least 2Ω(n

or to VSTAT O(n)k/3 .

1/8 )

≥ n2k calls to STAT O(n)−k/6



Remark. We remark that the well-conditioned assumption in Theorem 4.1 (i.e., that the distances between
the means and the largest and smallest eigenvalues of any covariance matrix are bounded) guarantees that
an SQ algorithm with a bounded number of SQ queries is possible.
The proof of Theorem 4.1 follows by an application of the framework developed in Section 3 and the
following proposition:
Proposition 4.2. For any ǫ > 0, there exists a distribution A on R that is a mixture of k Gaussians Ai ,
i ∈ [k], and satisfies the following conditions:
(i) A agrees with N (0, 1) on the first 2k − 1 moments.


√
1
and mean of magnitude O( k).
(ii) Each Gaussian component Ai has variance Θ k2 log2 (k+1/ǫ)

(iii) It holds dTV (Ai , Aj ) ≥ 1 − ǫ, for all i 6= j.

(iv) We have χ2 (A, N (0, 1)) ≤ exp(O(k)) log(1/ǫ).
(v) In the high-dimensional construction of Definition 3.1, we have that dTV (Pv , Pv′ ) ≥ 1/2 whenever
|v · v ′ | ≤ 1/2.
Given Proposition 4.2, the proof of Theorem 4.1 follows easily.
Proof of Theorem 4.1. First note that the distribution A given by Proposition 4.2, satisfies Condition 3.2 for
m = 2k − 1. By Proposition 3.3, applied for c = 1/4, any algorithm that is given SQ access to Pv , for an
1/8
unknown unit vector v ∈ S, and outputs a distribution Q with dTV (P, Q) ≤ ǫ, needs at least 2Ω(n ) ≥ n2k
calls to




p
STAT O(n)−k/4 · exp(O(k)) log(1/ǫ) or VSTAT O(n)k/2 /(exp(O(k)) log(1/ǫ)) .

For n = Ω(k8 log(1/ǫ)3 ), we have nk/2 ≥ (exp(O(k)) log(1/ǫ))3 , and so we need precision O(n)−k/6 for
STAT or O(n)k/3 for VSTAT.
It remains to show that Pv is a mixture of k Gaussians that satisfies the necessary conditions. Note
that Pv , when expressed in an appropriate basis, is a product of the mixture of k univariate Gaussians
and the standard (n − 1)-dimensional normal distribution.
Recall that
 the product of two Gaussians is

Pk
1
′
a Gaussian. If A =
i=1 wi N (µi , δ), where δ = Θ k 2 log(k+1/ǫ) (by Proposition 4.2 (ii)), then we

Pk
′
T . We can bound the variation distance between two
have that Pv =
i=1 wi N vµi , I − (1 − δ)vv
components by:


dTV N (vµ′i , I − (1 − δ)vv T ), N (vµ′j , I − (1 − δ)vv T ) = dTV N (µ′i , δ), N (µ′j , δ) ≥ 1 − ǫ ,

by Proposition 4.2 (iii). Also, we have that
n
o
1/2
max maxi,j kµi − µj k2 , maxi kΣi k2
√
max {maxi,j kµi − µj k2 , 1}


=
≤ O( k/δ) ≤ poly(k log(1/ǫ)) .
1/2
δ
mini 1/kΣ−1
i k2

This completes the proof.

21

4.1

Proof of Proposition 4.2 We start with the following lemma:

Lemma 4.3. There is a discrete distribution B on the real line, supported on k points,
that agrees with
√
N (0, 1) on the first 2k − 1 moments. All points x in the support of B have |x| = O( k).
Proof. This lemma essentially follows from standard techniques for Gaussian quadrature [AS72]. Given a
(possibly infinite) interval [a, b], a weighting function ω(x), and an integer k > 0, we can find xi and wi for
1 ≤ i ≤ k such that
Z b
k
X
wi p(xi ) ,
ω(x)p(x)dx =
a

i=1

for all polynomials p(x) of degree at most 2k − 1. The Gauss-Hermite quadrature is a standard implemen2
x ’s to be the
tation of this general scheme on the interval (−∞, ∞) with ω(x) = e−x . Here, we take the
√ i
2k−1 k! π
roots of the k-th (physicist’s) Hermite polynomial Hk (x). Then, we have that wi = k2 H (xi )2 .
k−1

2
√1 e−x /2 ,
2π

the pdf of N (0, 1). To do this, we need to rescale
√
def
the above wi and xi , and use the probabilist’s Hermite polynomials
Hek (x) = 2−k/2 Hk (x/ 2). We claim
√
k!
that we can take the x′i ’s to be the roots of Hek (x), i.e., x′i = 2xi and wi′ = k2 Hek−1
(x′ )2 . Indeed, we have
We would like to take ω(x) = G(x) :=

i

k
X

wi′ p(x′i ) =

i=1

k
X

k!

√

k2 Hek−1 (
i=1
k
k−1 √
X

2xi

)2

√
p( 2xi )

√
2 k! π
1
√
2xi )
·
·
p(
k2 Hk−1 (xi )2
π
i=1
Z ∞ √
1
2
p( 2y)e−y dy
=√
π −∞
Z ∞
√
1
2
=√
p(x)e−x /2 (1/ 2)dx
π
Z ∞ −∞
p(x)G(x)dx ,
=
=

−∞

for all polynomials p(x) of degree at most 2k − 1.
R∞
P
Note that all the weights are nonnegative by definition. Also note that ki=1 wi′ = −∞ 1 · G(x) = 1.
We take B to be the probability distribution with probability wi′ of being x′i , for each 1 ≤ i ≤ k. Then we
have
Z ∞
k
X
′ ′j
j
wi x i =
xj G(x)dx = EX∼N (0,1) [X j ] ,
EX∼B [X ] =
i=1

−∞

for all integers 1 ≤ j ≤ 2k − 1.
√
k), and so all roots of
It is known (see, e.g., [Sze89]) that all roots of Hk (x) have
absolute
value
O(
√
Hek . Hence, all points x in the support of B have |x| = O( k). This completes the proof.
On the other hand, if we want χ2 (A, N (0, 1)) to be finite, we need to have a mixture of Gaussians each
with positive variance δ > 0.

Corollary 4.4. For any 0 < δ < 1, there is a distribution A on R that is a mixture of k Gaussians each with
variance δ that agrees
√ with N (0, 1) on the first 2k − 1 moments. The means of all the Gaussian components
have magnitude O( k).

22

Proof. By rescaling the distribution B given by Lemma
4.3, we can find a discrete distribution B ′ supported
√
on k points with absolute value no bigger than O( k) that agrees with the first
√ 2k − 1 moments of N (0, 1 −
δ). The rescaled distribution B ′ assigns probability mass wi′ to the points 1 − δx′i , for 1 ≤ i ≤ k. Let
X ∼ B ′ , X ′ ∼ N (0, 1 − δ), and Y ∼ N (0, δ) that is independent of X, X ′ . We take A to be the distribution
of X + Y . Then, we have
j  
j  
X
X
i
i
E[X ′i ]E[Y j−i ] = E[(X ′ + Y )j ] ,
E[X i ]E[Y j−i ] =
E[(X + Y )j ] =
j
j
i=0

i=0

for all integers 1 ≤ j ≤ 2k − 1. By standard facts about Gaussians, √
X ′ +Y is distributed as N (0, 1). Finally,
note that the distribution of X + Y is a mixture of k Gaussians N ( 1 − δx′i , δ) with weights wi′ .

To appropriately set the parameter δ, we need to consider the high-dimensional construction (Definition 3.1):
 √

Lemma 4.5. For v, v ′ ∈ Sn with |v · v ′ | ≤ 1/2, we have that dTV (Pv , Pv′ ) ≥ 1 − O k δ log(1/δ) csc θ .

Proof. We write Ai , for 1 ≤ i ≤ k, for the Gaussians N (µi , δ) that p
A is a mixture of. Fix ǫ > 0. By a
Chernoff bound, Ai is within the interval [µi − a, µi + a], where a = 2 δ log(1/ǫ) with probability at least
1 − ǫ.
We again consider the plane spanned by v and v ′ . Let x, y be the orthogonal coordinates with v in the
direction of the x-axis. Similarly, let x′ , y ′ be the orthogonal coordinates with v ′ in the direction of the
x′ -axis. Let θ be the angle between v and v ′ . We have that
Z ∞ Z ∞
Z
min{Pv (x), Pv′ (x)}dx =
min{A(x)G(y), A(x′ )G(y ′ )}dxdy
x=−∞ y=−∞
x
Z ∞ Z ∞
=
min{A(x)G(y), A(x′ )G(y ′ )} csc θdxdx′
′
x=−∞ x =−∞
Z ∞ Z ∞
min{Ai (x)G(y), Aj (x′ )G(y ′ )} csc θdxdx′
≤ k max
i,j
x=−∞ x′ =−∞
Z µi +a Z µj +a
min{Ai (x)G(y), Aj (x′ )G(y ′ )} csc θdxdx′
≤ kǫ + k max
i,j

≤ kǫ +

x=µi −a x′ =µj −a
2
k max a csc θ max
min{Ai (x)G(y), Aj (x′ )G(y ′ )}
i,j
x,x′ ∈R

√
≤ kǫ + ka2 csc θ/(2π δ)
√
= kǫ + k δ csc θ log(1/ǫ)/π .
√
√
R
Taking ǫ = δ, we obtain that x min{Pv (x), Pv′ (x)}dx ≤ O(k δ log(1/δ) csc θ). On the other hand,
Z
1
|Pv (x) − Pv′ (x)|dx
dTV (Pv , Pv′ ) =
2 x
Z
1
(max{Pv (x), Pv′ (x)} − min{Pv (x), Pv′ (x)}) dx
=
2 x
Z
1
(Pv (x) + Pv′ (x) − 2 min{Pv (x), Pv′ (x)}) dx
=
2 x
Z
= 1 − min{Pv (x), Pv′ (x)}dx
x
√
≥ 1 − O(k δ log(1/δ) csc θ) .
This completes the proof.
23

This gives an upper bound on δ. We don’t want δ to be too small, because of the following lemma:
√
Lemma 4.6. We have that χ2 (A, N (0, 1)) ≤ exp(O(k))/ δ.
Proof. Each component Ai , for 1 ≤ i ≤ k, satisfies the following:
Z
2
1 + χ (Ai , N (0, 1)) = Ai (x)2 /G(x)dx
x
Z

1
exp −(x − µi )2 /δ + x2 /2 dx
=√
2πδ Zx

1
=√
exp −x2 (1/δ − 1/2) + 2µi x/δ − µ2i /δ dx
2πδ Zx

1
=√
exp −(x − 2µi /(2 − δ))2 ((2 − δ)/2δ) + 2µ2i /(δ(2 − δ)) − µ2i /δ dx
2πδ Zx

1
=√
exp −(x − 2µi /(2 − δ))2 ((2 − δ)/2δ) + 2µ2i /(δ(2 − δ)) − (2 − δ)µ2i /δ(2 − δ) dx
2πδ x
√
Z  p


2 exp(µ2i /(2 − δ))
p
1/ 2π(2δ/(2 − δ)) exp −(x − 2µi /(2 − δ))2 (1/δ − 1/2) dx
=
(2 − δ)δ
x
√

2
2 exp µi /(2 − δ)
p
=
(2 − δ)δ
√
≤ exp(O(k))/ 2δ .
Thus, for the mixture A we have that:
2

1 + χ (A, N (0, 1)) =
≤
≤

XX
i

j

i

j

i

j

XX

XX

wi wj /(1 − δ)

Z

Ai (x)Aj (x)/G(x)dx
x

q
wi wj /(1 − δ) (1 + χ2 (Ai , N (0, 1)) (1 + χ2 (Aj , N (0, 1))
√
wi wj /(1 − δ) · exp(O(k))/ 2δ

XX
√
wi wj
= exp(O(k))/ 2δ ·
i

√
= exp(O(k))/ 2δ · 1 .

j

This completes the proof.
The following simple lemma helps us enforce the condition that the Gaussian components are wellseparated:


1
Lemma 4.7. Given ǫ > 0, if δ ≤ O k log(1/ǫ)
, then dTV (Ai , Aj ) ≥ 1 − ǫ.

√
Proof. It is known (see, e.g., [Sze89]) that the difference between two roots of Hk (x) is√Ω(1/ k). Thus,
the same is true of Hek (x) and by our construction, we have that |µi − µj | ≥ Ω((1 − δ)/pk), for i 6= j. By
standard Chernoff bounds, with probability at least 1 − ǫ/2, Ai lies in the range (µi − p 2δ ln(2/ǫ), µi +
p
p2δ ln(2/ǫ)). Similarly, with probability at least 1 − ǫ/2, Aj lies in the range (µj − 2δ ln(2/ǫ),
√µj +
2δ
ln(2/ǫ)).
If
these
intervals
are
disjoint,
we
have
d
(A
,
A
)
≥
1−ǫ.
This
holds
when
(1−δ)/
k) =
TV
i
j
p
Ω( δ log(1/ǫ)), which is true when δ ≤ O(1/k log(1/ǫ)).
24

We now have all the necessary tools to prove Proposition 4.2. We take δ = C/(k2 log2 (k + 1/ǫ)) for a
sufficiently small constant C > 0. Combined with Corollary 4.4, this gives condition (ii). For condition (i),
note that, by Corollary 4.4, A agrees with N (0, 1) on the first 2k − 1 moments. Since δ was selected to be
smaller than O(1/k log(1/ǫ)), Lemma 4.7 gives condition (iii). Lemma 4.6 gives condition (iv). Finally, by
our choice of δ and Lemma 4.5, we get condition (v). This completes the proof.

5 SQ Lower Bounds for Robust Learning of a Gaussian
In this section, we prove our super-polynomial SQ lower bounds for robustly learning a high-dimensional
Gaussian. In Section 5.1, we show our lower bound for robustly learning an unknown mean spherical
Gaussian. In Section 5.2, we give our lower bound for robustly learning a zero mean unknown covariance
Gaussian with respect to the spectral norm.
5.1 Robust Learning Lower Bound for Unknown Mean Gaussian In this subsection, we use the framework of Section 3 to prove the following theorem:
Theorem 5.1. Let 0 < ǫ < 1 and n ≥ log(1/ǫ)Ω(1) . Fix any M ∈ Z+ such that M = O(log1/2 (1/ǫ)),
where the universal constant in the O(·) is assumed to be sufficiently small. Any algorithm that, given
SQ access to a distribution P on Rn which satisfies dTV (P, N (µ, I)) ≤ ǫ for an unknown µ ∈ Rn with
kµk2 ≤ poly(n/ǫ), and returns a hypothesis distribution Q with dTV (Q, P) ≤ O(ǫ log(1/ǫ)1/2 /M 2 ),


1/12
requires at least 2Ω(n ) ≥ nM calls to STAT O(n)−M/6 or to VSTAT O(n)M/3 .
The theorem will follow from the following proposition:

Proposition 5.2. For any δ > 0 and m ∈ Z+ , there is a distribution A on R satisfying the following
conditions:
(i) A and N (0, 1) agree on the first m moments.
p
(ii) dTV (A, N (δ, 1)) ≤ O(δm2 / log(1/δ)).

(iii) χ2 (A, N (0, 1)) = O(δ).

Before we prove Proposition 5.2, we show how Theorem 5.1 easily follows from it using the machinery
developed in Section 3.
Proof of Theorem 5.1. We can assume without loss of generality that ǫ > 0 is smaller than a sufficiently
small universal constant. We apply Propositions 3.3 and 5.2, with the parameter c set to c = 1/6, m = M ,
and δ = Cǫ ln(1/ǫ)1/2 /M 2 , where C > 0 is a sufficiently large constant. By Proposition
5.2, we have that
p
(i) A and N (0, 1) agree on the first m moments, (ii) dTV (A, N (δ, 1)) ≤ O(δm2 / log(1/δ)) = O(ǫ) and
(iii) χ2 (A, N (0, 1)) = O(δ). Note that for any unit vector v ∈ Sn , it holds that
dTV (Pv , N (vδ, I)) = dTV (A, N (δ, 1)) ≤ O(Cǫ) .
Therefore, for any unit vectors v, v ′ with |v · v ′ | ≤ 1/8 we have that:
dTV (Pv , Pv′ ) ≥ dTV (N (δv, I), N (δv ′ , I)) − dTV (Pv , N (vδ, I)) − dTV (Pv′ , N (v ′ δ, I))
≥ Ω(δkv − v ′ k2 ) − O(Cǫ)
√
= Ω(δ 2 − 2v · v ′ ) − O(Cǫ)
= Ω(δ) − O(Cǫ) = Ω(δ) ,

25

(9)

where the first line is the triangle inequality, the second line uses (9) and that dTV (N (µ1 , I), N (µ2 , I)) =
Ω(kµ1 − µ2 k) when kµ1 − µ2 k is smaller than an absolute constant, the third line uses the assumption that
|v · v ′ | ≤ 1/8, and the last line follows from the definition of δ.
We want to apply Proposition 3.3 with its “local parameter” ǫ taken to be Ω(δ). The assumption on n
in the statement of Theorem 5.1, i.e., n ≥ log(1/ǫ)Ω(1) , implies the condition n ≥ mΩ(1) in the statement
of Proposition 3.3. By our choice of c = 1/6, we conclude
that any
SQ algorithm for our learning problem

√ 

1/12 )
Ω(n
M
−M/6
δ or to VSTAT O(n)M/3 /δ to produce
requires at least 2
≥ n queries to STAT O(n)
a hypothesis distribution Q with
dTV (Q, Pv ) ≤ O(δ) ≤ O(ǫ ln(1/ǫ)1/2 /M 2 ) ,
where we used the assumption that C is sufficiently large. This completes the proof.
Proof of Proposition 5.2 The rest of this section is devoted to the proof of Proposition 5.2. We start
by describing the outline of the proof. We then provide a number of intermediate useful lemmas that we
subsequently combine to complete the proof.
p
The proof plan proceeds as follows. For some C = Θ( log(1/δ)), we define the one-dimensional
distribution A to be:
• For x ∈
/ [−C, C], we define A(x) = G(x − δ).
• For x ∈ [−C, C], we define A(x) = G(x − δ) + p(x), where p(x) is the degree-m polynomial with
RC
−C p(x)dx = 0 and
Z

C

i

p(x)x dx =
−C

Z

∞

−∞

(G(x) − G(x − δ))xi dx ,

(10)

for 1 ≤ i ≤ m. (We note that p is unique after fixing m, C and δ.)
We need to show that we can
p find appropriate values for the parameters m, C, and δ such that the L1 -norm
of p(x) is at most O(δm2 / log(1/δ)) and that A(x) is non-negative. To achieve that, we will express p(x)
as a linear combination of (appropriately scaled) Legendre polynomials, a family of orthogonal polynomials
on [−C, C]. Rather than directly showing that the first m moments agree, we will instead want that the
expectations of the first m scaled Legendre polynomials agree. Bounds on the coefficients of the Legendre
polynomials in p(x) allow us to obtain bounds on the L1 and L∞ norms of p(x) on [−C, C]. Choosing m,
δ, and C appropriately will complete the proof of the proposition.
Properties of Legendre Polynomials
we will need:

We start by recording the properties of Legendre polynomials that

Fact 5.3. [Sze89] The Legendre polynomials, Pk (x), for k ∈ Z+ , satisfy the following properties:
(i) Pk (x) is a degree-k polynomial, P0 (x) = 1, and P1 (x) = x.
R1
(ii) −1 Pi (x)Pj (x)dx = (2/(2i + 1))δi,j for all i, j ≥ 0.

(iii) |Pk (x)| ≤ 1 for all |x| ≤ 1.

(iv) Pk (x) = (−1)k Pk (−x).
P⌊k/2⌋ 
(v) Pk (x) = (1/2k ) i=0 ki

2k−2i k−2i
x
.
k

26

As a simple corollary we obtain the following lemma:
Corollary 5.4. We have:
(i) |Pk (x)| ≤ (4|x|)k for all |x| ≥ 1.
√
R1
(ii) −1 |Pk (x)|dx ≤ O(1/ k).

We are now ready to proceed with the formal proof. The main technical result of this section is the following
lemma:
P
3/2 /C 2 ), for 0 ≤ k ≤ m.
Lemma 5.5. We can write p(x) = m
k=0 ak Pk (x/C), where |ak | = O(δk

Before we give the proof of Lemma 5.5, we deduce two corollaries that will be useful in the proof of
Proposition 5.2. First, we can obtain bounds on the L1 and L∞ norms of p(x) on [−C, C]. As an immediate
corollary of Lemma 5.5 and the aforementioned properties of Legendre polynomials, we deduce:
RC
Corollary 5.6. We have that: −C |p(x)|dx ≤ O(δm2 /C) and |p(x)| ≤ δm5/2 /C 2 , for all x ∈ [−C, C].
We now bound from above the desired χ2 -divergence:


Lemma 5.7. χ2 (A, N (0, 1)) = O δ2 + δm5/2 /C 2 · (C 2 δ2 + max|x|≤C |p(x)|/G(x)) .
Proof. We have the following:
Z ∞
A(x)2 /G(x)dx − 1
χ2 (A, N (0, 1)) =
=

−∞
∞

Z

−∞

Z

G(x − δ)2 /G(x)dx − 1 +

For the first term, we note that:
Z ∞
Z
2
G(x − δ) /G(x)dx =
−∞

∞

−∞

C
−C

2p(x)G(x − δ)/G(x)dx +

Z

C

p(x)2 /G(x)dx .

−C

G(x − 2δ) exp(δ2 )dx = exp(δ2 ) ≤ 1 + 2δ2 .

We bound the second term from above as follows:
Z C
Z
2p(x)G(x − δ)/G(x)dx =
−C

=

C

−C
Z C
−C

2p(x) exp(xδ − δ2 /2)dx
p(x) · (1 + xδ + O(C 2 δ2 ))dx
2 2

≤ C|a0 | + O(Cδ|a1 |) + C δ

Z

≤ 0 + O(δ2 /C) + O(δ3 m5/2 ) ,

C
−C

|p(x)|dx

where the last lines uses Lemma 5.5 and Corollary 5.6. Finally, for the third term we have:
Z

C
−C

p(x)2 /G(x)dx ≤ δm5/2 /C 2 max |p(x)|/G(x) ,
x∈[−C,C]

where the inequality follows from Corollary 5.6. This completes the proof of Lemma 5.7.

27

Proof of Lemma 5.5 We first note that we can express p(x) as a linear combination of scaled Legendre
polynomials whose coefficients are explicitly given by integrals:
RC
P
Claim 5.8. We can write p(x) = m
k=0 ak Pk (x/C), where ak = ((2k + 1)/2C) −C Pk (x/C)p(x)dx.

Proof. Since p(x) has degree at most m and the set of polynomials Pk (x/C),
Pm0 ≤ k ≤ m, contains a
polynomial of each degree from 0 to m, there exists ak ∈ R such that p(x) = k=0 ak Pk (x/C).
RC
It follows from Fact 5.3 (ii) and a change of variables that −C Pi (x/C)Pj (x/C)dx = (2C/(2i+1))δi,j ,
for all i, j ≥ 0. We can use this to extract the ak ’s. For 1 ≤ k ≤ m, we have
RC

Pk (x/C)p(x)dx =

RC

Pk (x/C)

ai Pi (x/C)dx =

i=0

−C

−C

m
P

m
P

ai

i=0

RC

Pk (x/C)Pi (x/C)dx = (2C/(2k+1))ak .

−C

Since the first m moments of p are fixed, via 10, we obtain:
Z

C

p(x)Pk (x/C)dx =
−C

Z

∞
−∞

(G(x) − G(x − δ))Pk (x/C)dx ,

(11)

for any 0 ≤ k ≤ m.
Since we will apply this with the parameter 1/δ exponential in m and C, we will be able to ignore O(δ2 )
terms. We use Taylor’s theorem to expand (G(x) − G(x − δ)) up to second order terms:
Fact 5.9. G(x) − G(x − δ) = xG(x)δ + (ξ(x)2 − 1)/2 · G(ξ(x))δ2 , for some x ≤ ξ(x) ≤ x + δ.
R∞
By (11) and Fact 5.9, to bound the magnitude of the ak ’s, it suffices to bound the terms −∞ Pk (x/C)xG(x)dx
R∞
and −∞ Pk (x)(ξ(x)2 − 1)/2 · G(ξ(x))dx. This is done in the following two lemmas.
√
R∞
Lemma 5.10. For k ≤ 4C, we have that −∞ Pk (x/C)xG(x)dx ≤ O( k/C).

Proof. When k is even, using Fact 5.3 (iv), we have that Pk (x/C)xG(x) = −(Pk (−x/C) · (−x)G(−x)),
and so the integral is zero. When k is odd, we can rewrite Fact 5.3 (v) in ascending order of terms, by using
the change of variables j = (k + 1)/2 − i, as
k

Pk (x) = (1/2 )

(k+1)/2 

X
j=1

k
(k − 1)/2 + j




k + 2j − 1 2j−1
x
.
2j − 1

By standard results about the moments of Gaussians, for all j ≥ 1, we have that
Q
1)!! := ji=1 (2i + 1). Thus, we can write
Z

∞

−∞

Pk (x/C)xG(x)dx =

Z

∞

k

1/2

−∞
k

= (1/2 )

(k+1)/2 

X
j=1

(k+1)/2 

X
j=1

k
(k − 1)/2 + j

k
(k − 1)/2 + j

28





R∞

−∞ x

2j G(x)dx

= (2j +


k + 2j − 1 2j
x G(x)/C 2j−1 dx
2j − 1


k + 2j − 1
(2j + 1)!!/C 2j−1 .
2j − 1

Note that this quantity is non-negative. We can bound it from above as follows:


Z ∞
(k+1)/2 
X
k
k + 2j − 1
k
Pk (x/C)xG(x)dx = (1/2 )
(2j + 1)!!/C 2j−1
(k
−
1)/2
+
j
2j
−
1
−∞
j=1

X



√
k + 2j − 1
(1/ k)
(2j + 1)!!/C 2j−1
2j − 1

X

√
(1/ k)(k + 2j − 1)2j−1 (2j + 1)!!/C 2j−1 (2j − 1)!

X

√
(1/ k)2(k + 2j − 1)2j−1 /C 2j−1

(k+1)/2

≤

j=1

(k+1)/2

≤

j=1

(k+1)/2

≤

j=1

X
√
√ (k+1)/2
(2k/C)2j−1 ≤ 8 k/C .
≤ (2/ k)
j=1

The proof of Lemma 5.10 is now complete.
R∞
Lemma 5.11. For any integer 1 ≤ k ≤ 4C, we have that −∞ Pk (x)(ξ(x)2 − 1)/2 · G(ξ(x))dx ≤ O(1) ,
for any function ξ(x) with x ≤ ξ(x) ≤ x + δ, for all x ∈ R.

Proof. We separate this integral into the interval [−C, C] and the tails. We can use Fact 5.3 (iii) to bound
the integral on [−C, C], as follows:
Z C
Z C
Pk (x)(ξ(x)2 − 1)/2 · G(ξ(x))dx ≤
(ξ(x)2 − 1)/2 · G(ξ(x))dx
−C

≤

Z

−C
C

−C

(|x| + δ + 1)2 /2 · G(min{0, |x| − δ})dx

≤ O(δ) +

Z

C+δ

−C−δ

(|x| + 2δ + 1)2 /2 · G(x)dx


≤ O δ + EX∼G [1] + EX∼G [|X|] + EX∼G [X 2 ]
= O(1) .

For the tails, we need Corollary 5.4(i). For the right tail, [C, ∞), we have
Z ∞
Z ∞
2
Pk (x)(ξ(x) − 1)/2 · G(ξ(x))dx ≤
(4|x|/C)k (ξ(x)2 − 1)/2 · G(ξ(x))dx
C
ZC∞
≤
(4|x|/C)k (x + δ)2 G(x − δ)dx
ZC∞
≤
(4/C)k |x + 2δ|k+2 G(x − δ)dx
C
Z ∞
≤
(4/C)k |x|k+2 · (1 + 2δ/C)k G(x − δ)dx
C−δ
Z ∞
≤2
(4/C)k |x|k+2 G(x − δ)dx
−∞

≤ O((4/C)k (k + 3)!!)
√
≤ O((4 k/C)k ) ≤ O(1) .
29

A similar bound holds for the left tail, which completes the proof.
Putting everything together, gives Lemma 5.5.
To prove Proposition 5.2, we need to set C appropriately and check the bounds on m needed for A(x)
to satisfy the necessary properties.
p
Proof of Proposition 5.2. Note that unless δ is sufficiently small and m2 ≤ O( log(1/δ)),
p taking A =
N (0, 1), instead of using our√construction, satisfies the proposition. We will take C = Θ( log(1/δ)), and
so we can assume that m ≤ C.
Recall that A(x) is defined to be G(x − δ) + p(x) on [−C, C] and G(x − δ) outside of [−C, C]. Firstly,
A(x) needs to be the pdf of a distribution. Since
Z ∞
Z C
(G(x) − G(x − δ))Pk (x/C)dx
Pk (x/C)p(x)dx =
−C

R∞

−∞

for k = 0, when Pk (x) = 1, we have that −∞ A(x)dx = 1. We also need that A(x) is non-negative, i.e.,
that A(x) = G(x − δ) + p(x) ≥ 0 for all x ∈ [−C, C]. Note that
G(x) + p(x) ≥ G(C + δ) − δm5/2 /C 2 ,

p
using Corollary 5.6. Since m2 ≤ C, we need G(C + δ) ≥ δC 3/4 . This holds when C = ln(1/δ) − δ,
p
p
3/4
since then we have G(C +δ) = δ/2π ≥ δ ln(1/δ) for sufficiently small δ. Note that this also implies
that A(x) ≤ 2G(x − δ) for all x, and |p(x)| ≤ G(x) for all −C ≤ x ≤ C.
The second of these and Lemma 5.7 imply (iii). For (i), by construction, we have that the first m
moments agree.
A satisfies (ii), since by Lemma 5.6 ,
Z
p
1 C
|p(x)|dx ≤ O(δm2 /C) = O(δm2 / log(1/δ)).
dTV (A, N (δ, 1)) =
2 −C
The proof of Proposition 5.2 is now compete.

5.2 Robust Learning Lower Bound for Unknown Covariance Gaussian In this subsection, we prove
an SQ lower bound for robustly learning the covariance matrix of a high-dimensional Gaussian with known
mean. We note that our lower bound applies even for spectral norm approximation. In particular, we show:
Theorem 5.12. Let ǫ > 0 and n ≥ Ω(log2 (1/ǫ)). Fix any M ∈ Z+ such that M = O(log1/4 (1/ǫ)),
where the universal constant in the O(·) is assumed to be sufficiently small. Any algorithm that, given SQ
access to a distribution P on Rn which has dTV (P, N (0, Σ)) ≤ ǫ for some positive-definite Σ ∈ Rn×n with
e with kΣ
e − Σk2 ≤ O(ǫ log(1/ǫ)/M 4 ), requires at least 2n2/15 ≥ nM
I/2  Σ  2I, and returnsa matrix Σ

calls to STAT O(n)−M/6 or to VSTAT O(n)M/3 .
The theorem will follow from the following proposition:

Proposition 5.13. For any 1/3 > δ > 0 and integer log(1/δ)1/4 ≫ m > 0, there is a distribution A on R
satisfying the following conditions:
(i) A and N (0, 1) agree on the first m moments.
(ii) dTV (A, N (0, (1 − δ)2 ) ≤ O(δm4 / log(1/δ)).
(iii) χ2 (A, N (0, 1)) = O(1 + m8 δ3/2 /log(1/δ)5/2 ).
30

As in the previous subsection, Theorem 5.12 follows easily from Proposition 3.3 and Proposition 5.13.
Proof of Theorem 5.12. We apply Propositions 3.3 and 5.13 with m = M and δ = 2ǫ ln(1/ǫ)/M 4 . By
Proposition 5.13, we have that (i) A and N (0, 1) agree on the first M moments, (ii) dTV (A, N (0, (1 −
δ)2 )) ≤ O(δm4 / log(1/δ)) = O(ǫ) and (iii) χ2 (A, N (0, 1)) = O(1).
Note that we cannot directly apply Proposition 3.3, since we are not aiming to learn within small total
variation distance. Instead, we are interested in a different search problem, that of finding an approximation
e to the covariance Σ with kΣ
e − Σk2 ≤ ǫ log(1/ǫ)/M 4 , where Σ is the covariance of a mean 0 Gaussian
Σ
within ǫ total variation distance. Note that for Pv , we have that Σ = I − (1 − (1 − δ)2 )vv T . We need to
e the set Z −1 (Σ)
e = {Pv :
argue that this search problem has at most one solution in S, i.e., that for any Σ,
2
T
4
−1
e
e
v ∈ S and kΣ−I −(1−(1−δ) )vv k2 ≤ ǫ log(1/ǫ)/M } has |Z (Σ)| ≤ 1, where S is as in Lemma 3.7.
Lemma 5.14. For S as in Lemma 3.7 with c = 1/6 and with n larger than a sufficiently large constant,
e − I − (1 − (1 − δ)2 )vv T k2 ≤ ǫ log(1/ǫ)/M 4 }| ≤ 1, for all Σ.
e
|{Pv : v ∈ S and kΣ

e and let v, v ′ be distinct
Proof. Suppose for a contradiction that this set has size at least 2 for some Σ
2
T
e − Σv k2 ≤
elements. Let Σv = I − (1 − (1 − δ) )vv and define Σv′ similarly. Then we have kΣ
4
4
e
ǫ log(1/ǫ)/M and kΣ − Σv′ k2 ≤ ǫ log(1/ǫ)/M . By the triangle inequality, we have kΣv − Σv′ k2 ≤
2ǫ log(1/ǫ)/M 4 . However, we also have that |v ·v ′ | ≤ O(n−1/3 ) ≤ 1/2. Now we get that v T Σv = (1− δ)2 ,
but v T Σv′ v = (1 − |v · v ′ |2 ) · 1 + |v · v ′ |2 · (1 − δ)2 ≥ 3/4 + (1/4)(1 − δ)2 . We thus obtain
kΣv − Σv′ k2 ≥ v T (Σv′ − Σv )v ≥ (3/4)(1 − (1 − δ)2 ) ≥ (3/2)δ − O(δ2 ) > δ ,
where the last inequality assumes that ǫ is at most an appropriately small universal constant. Since δ =
2ǫ ln(1/ǫ)/M 4 , this leads to a contradiction.
Now we use the the proof of Proposition 3.3 with c = 1/6 and ǫ taken to be Ω(δ). Our condition on n,
8 ). We conclude that it requires at least 2n2/15 ≥ nM calls to
n ≥ ln(1/ǫ)2 , implies
the condition n ≥ Ω(m


e such that kΣ
e − Σk2 ≤ ǫ log(1/ǫ)/M 4 .
STAT O(n)−M/6 or to VSTAT O(n)M/3 to produce a a Σ
Proof of Proposition 5.13. Similarly to the previous subsection, we choose to define the univariate distribution A to have probability density function given by
A(x) = G(x/(1 − δ))/(1 − δ) − p(x)1[−C,C] ,

p
where C is a sufficiently small multiple of log(1/δ) and p(x) is the unique degree-m polynomial that
causes APand G(x) (the pdf of N (0, 1)) to have the same
R ∞ first m moments. Once again, we may write
m
p(x) = k=0 ak Pk (x/C), where ak = ((2k + 1)/2C) −∞ (G(x) − G(x/(1 − δ))/(1 − δ))Pk (x/C)dx.
The bulk of our proof will now be in bounding the ak ’s.
The first thing to note is that since G(x) − G(x/(1 − δ)) is even, ak is 0 for k odd. For k even, we will

31

need to compute this expression using Fact 5.3 (v). In particular, we have that
Z
2k + 1 ∞
ak =
(G(x) − G(x/(1 − δ))/(1 − δ))Pk (x/C)dx
2C
−∞

Z
⌊k/2⌋  
X k
2k − 2i
2k + 1 ∞
−k
(G(x) − G(x/(1 − δ))/(1 − δ))2
(x/C)k−2i dx
=
2C
i
k
−∞
i=0

2k + 1
=
2C
=

Z

∞
−∞

−k

(G(x) − G(x/(1 − δ))/(1 − δ))2

⌊k/2⌋ 

2k + 1 X
2k+1 C
j=0

=

k
k/2 + j





k + 2j
C −2j
k

Z

∞
−∞

⌊k/2⌋ 

X

k
k/2 + j

j=0




k + 2j
(x/C)2j dx
k

(G(x) − G(x/(1 − δ))/(1 − δ))x2j dx



⌊k/2⌋ 
k
k + 2j
2k + 1 X
C −2j (2j − 1)!!(1 − (1 − δ)2 j)
2k+1 C
k/2 + j
k
j=0


⌊k/2⌋ 
δ(2k + 1) X k + 2j
≤
C −2j (2j − 1)!!j
C
2j
j=0

∞ 
δ(2k + 1) X 2k 2j
≤
C
C
j=1

3

≤ δ10k C

−3

,

where in the last step we assume that k is less than a sufficiently small multiple of C.
It is now clear that A is a pseudo-distribution that matches its first m moments with N (0, 1). Firstly,
in order to check that A is a distribution, it
clear that A(x) > 0 for |x| > C. For |x| ≤ C we have
Pis
m
4 −3 . Since this is smaller than δ 1/2 <
that |A(x) − G(x/(1 − δ))/(1 − δ)| ≤
k=0 |ak | ≤ δ10m C
G(x/(1 − δ))/(1 − δ), we have that A(x) ≥ 0 everywhere.
Next, we need to bound from above dTV (A, N (0, 1 − δ)), i.e., the L1 -distance between A(x) and
G(x/(1 − δ))/(1 − δ). This in turn is at most
m Z
X
k=0

C
−C

|ak Pk (x/C)|dx ≤

m
X

δ20k3 C −2 = O(δm4 / log(1/δ)).

k=0

Finally, we need to bound from above χ2 (A, N (0, 1)). Note that
χ2 (A, N (0, 1)) ≤ O(χ2 (N (0, 1 − δ), N (0, 1)) + χ2 (A − N (0, 1 − δ), N (0, 1))) .
It is easy to see that χ2 (N (0, 1 − δ), N (0, 1)) = O(1 + δ). On the other hand, we have that
χ2 (A − N (0, 1 − δ), N (0, 1))) ≤ (m + 1)

m
X

a2k

k=0

= O(m7 δ2 C −6 )
= O(m8 δ3/2 C
This completes the proof.
32

Z

C

Pk (x/C)2 /G(x)dx

−C

m Z C
X

k=0
−5

).

−C

G(x)−2 dx

6 Statistical and Computational Tradeoffs
In this section, we prove our SQ lower bounds establishing statistical-computational tradeoffs for two natural
robust estimation problems. In Section 6.1, we give a sharp-tradeoff for the problem of robustly estimating
the covariance matrix in spectral norm. In Section 6.2, we show such a tradeoff for robust sparse mean
estimation.
6.1 Robust Estimation of Covariance Matrix in Spectral Norm In this subsection, we establish an SQ
lower bound for robust covariance estimation in spectral norm. Our SQ lower bound provides evidence for
the existence of a statistical-computational tradeoff for this problem. Roughly speaking, we show that, for
any constant c > 0, given samples from a corrupted n-dimensional Gaussian N (0, Σ), any computationally
efficient SQ algorithm that approximates Σ within a factor of 2 requires Ω(n2−c ) samples. Our lower bound
applies even to the weaker Huber contamination model.
We note that the information-theoretic optimum for this problem is known to be Θ(n) samples (and is
achievable by an exponential time SQ algorithm). Hence, our lower bound establishes a nearly-quadratic
gap in the sample complexity between efficient and inefficient SQ algorithms for this problem. Formally,
we show:
Theorem 6.1. Let 0 < c < 1/6, and n > 125. Any algorithm that, given SQ access to a distribution
P on Rn of the form P = (1 − ǫ)N (0, Σ) + ǫN1 , where ǫ ≤ c/ ln(n) and N1 is a noise distribution,
for some covariance Σ with kΣk2 ≤ poly(n/ǫ), and returns a matrix Σ̃ with Σ̃/2  Σ  2Σ̃, requires


c/3
at least 2Ω(n ) calls to STAT O(n)−(1−5c/2) or to VSTAT O(n)2−5c . Furthermore, the result holds
even when the noise distribution N1 is a mixture of 2 Gaussians.
Proof. Let ǫ = c/ ln(n). We consider the following mixture of 3 Gaussians:
p
p
A = (1 − ǫ)N (0, (1/5 − ǫ)/(1 − ǫ)) + (ǫ/2) · N ( 4/(5ǫ), 1) + ǫ/2 · N (− 4/(5ǫ), 1) .

Note that A is symmetric about 0 and so, for X ∼ A, we have EX∼A [X] = EX∼A [X 3 ] = 0. The variance
of A is VarX∼A [X] = EX∼A [X 2 ] = (1/5 − ǫ) + ǫ · 4/(5ǫ) + ǫ = 1. That is, A agrees with N (0, 1) on the
first 3 moments. We need a bound on χ2 (A, N (0, 1)). For this, we use the following three easy facts (see
Appendix D for the simple proofs):
Fact 6.2. For distributions B, C, D and w ∈ [0, 1], we have that χ2 (wB + (1 − w)C, D) = w2 χ2 (B, D)+
(1 − w)2 χ2 (C, D) + 2w(1 − w)χD (B, C).
Fact 6.3. For µ, µ′ ∈ R, we have that χN (0,1) (N (µ′ , 1), N (µ, 1)) = exp(µµ′ ) − 1.
p
Fact 6.4. We have that χ2 (N (0, σ 2 ), N (0, 1)) = 2/σ 4 − 1/σ 2 − 1.

Note that 1/6 ≤ (1/5 − ǫ)/(1 − ǫ) ≤ 1/5. Fact 6.4 now yields
p
√
χ2 (N (0, (1/5 − ǫ)/(1 − ǫ)), N (0, 1)) ≤ 2/(1/6)2 + 1/(1/6) − 1 = 78 − 1 ≤ 8 .

Using Fact 6.2, we can write:

p
χ2 (A, N (0, 1)) ≤ (1 − ǫ)2 χ2 (N (0, (1/5 − ǫ)/(1 − ǫ)), N (0, 1)) + (ǫ2 /2)χ2 (N ( 4/5ǫ, 1), N (0, 1))
p
p
p
+ (ǫ2 /2)χN (0,1) (N ( 4/5ǫ, 1), N (− 4/5ǫ, 1)) + (1 − ǫ)ǫχN (0,1) (N (0, (1/5 − ǫ)), N ( 4/5ǫ, 1))
≤ 8 + (ǫ2 /2)(exp(4/(5ǫ)) − 1) + 0
q
p
+ (1 − ǫ)ǫ χ2 (N ( 4/5ǫ, 1), N (0, 1))χ2 (N (0, (1/5 − ǫ)/(1 − ǫ)), N (0, 1))
√
≤ 8 + ǫ2 exp(4/(5ǫ)) + ǫ 8 exp(2/(5ǫ)) ≤ O(1 + exp(1/ǫ))
≤ O(nc ) .

33

Note that we cannot directly apply Proposition 3.3, since we are not aiming to learn within small variation
distance. Instead, we are interested in a different search problem, that of approximating the covariance Σv
of the (1 − ǫ) weight component of Pv to within a factor of 2. We need to argue that this search problem
has at most one solution in S, i.e., that for any Σ, the set Z −1 (Σ) = {Pv : v ∈ S and Σ  Σv  2Σ} has
|Z −1 (Σ)| ≤ 1, where S is as in Lemma 3.7.
Lemma 6.5. For S as in Lemma 3.7, |{Pv : v ∈ S and Σ  Σv  2Σ}| ≤ 1 for all Σ.
Proof. Suppose for a contradiction that |Z −1 (Σ)| ≥ 2 for some Σ. Then there are distinct v, v ′ ∈ S with
Σ  Σv  2Σ and Σ  Σv′  2Σ. However, we have that |v · v ′ | ≤ O(nc−1/2 ) ≤ n−1/3 . Now
v T Σv v = (1/5 − ǫ)/(1 − ǫ) < 1/5, but
v T Σv′ v = (1 − |v · v ′ |2 ) · 1 + |v · v ′ |2 · (1/5 − ǫ)/(1 − ǫ) ≥ 1 − 5|v · v ′ |2 /6 ≥ 1 − n−1/3 > 4/5 .
Thus, we need v T Σv ≤ 2v T Σv v < 2/5, but v T Σv ≥ v T Σv′ v/2 > 2/5. This is a contradiction and so
|Z −1 (Σ)| < 1.
Now the proof of Proposition 3.3 applies and we obtain that any algorithm that outputspa Σ̃ satisfying the
c/2
desired conditions uses at least 2Ω(n ) ≥ nm+1 queries to STAT(O(n)−(m+1)(1/4−c/2) χ2 (A, N (0, 1)))
or to VSTAT(O(n)(m+1)(1/2−c) /χ2 (A, N (0, 1))), where m = 3. Now substituting χ2 (A, N (0, 1) ≤
c/2
O(nc ), we get that we need at least 2Ω(n ) queries to STAT(O(n)−(1−5c/2) ) or to VSTAT(O(n)(2−5c) ).
The proof is now complete.
6.2 Robust Sparse Mean Estimation In this subsection, we establish an SQ lower bound for robust
sparse mean estimation. Our SQ lower bound gives evidence for the existence of a statistical-computational
tradeoff for this problem. Roughly speaking, we show that, for any constant c > 0, given samples from
a corrupted n-dimensional Gaussian N (µ, I), where the mean vector µ is k-sparse, any computationally
efficient SQ algorithm that approximates the true mean requires Ω(k2−c ) samples. Our lower bound applies
even to the weaker Huber contamination model.
We note that the information-theoretic optimum for this problem is known to be Θ(k log n) (and is
achievable by an exponential time SQ algorithm). Hence, our lower bound establishes a nearly-quadratic
gap in the sample complexity between efficient and inefficient SQ algorithms. Formally, we show:
Theorem 6.6. Fix any constant 0 < c < 1. Let k, n ∈ Z+ be sufficiently large positive integers satisfying
n ≥ 8k2 . Any algorithm which, given SQ access to a distribution P on Rn such that P = (1 − δ)N (µ, I) +
δN1 , where N1 is an arbitrary distribution, δ = ǫ/kc/4 , and µ ∈ Rn is promised to be k-sparse with
c
kµk2 = ǫ, and outputs a hypothesis vector µ
b satisfying kb
µ − µk2 ≤ ǫ/2, requires at least Ω(nck /8 ) queries
to STAT(O(k)3c/2−1 ) or to VSTAT(O(k)2−3c ).

To prove our result, we will use the framework of Section 3: We will construct a suitable one-dimensional
distribution A and consider an appropriate collection of distributions Pv , but this time only for k-sparse unit
vectors v on Rn . We start by showing an analogue of Lemma 3.7 for k-sparse vectors, and then use it to
prove an analogue of Proposition 3.3. Our analogue of Lemma 3.7 is the following:
Lemma 6.7. Fix a constant 0 < c < 1. There exists a set S of k-sparse unit vectors on Rn of cardinality
c
|S| = ⌊nck /8 ⌋ such that for each pair of distinct vectors v, v ′ ∈ S we have that |v · v ′ | ≤ 2kc−1 .

Proof. Let √
D be the uniform distribution over the set of vectors v on Rn that have exactly k coordinates
equal to 1/ k and the rest n − k coordinates equal to zero. Consider the distribution D ′ of the inner product
v · v ′ , with v and v ′ independently drawn from D. Then, we have that v · v ′ = i/k, where i is the number of
34

non-zero coordinates that v ′ and v ′ have in common. Note that the distribution D ′ does not change if we fix
v ′ and therefore we have that
 
  
k
n−k
n
′
Pr[v · v = i/k] =
/
.
i
k−i
k
In other words the random variable k(v ·v ′ ) is distributed as the hypergeometric distribution with parameters
(n, k, k). By standard tail bounds on the hypergeometric distribution, for t > 0, we have


Pr k(v · v ′ ) ≥ k(k/n + t) ≤ exp (−k · KL (t + k/n||k/n)) ,

where KL(a||b) = a ln(a/b) + (1 − a) ln( 1−a
1−b ). We apply this concentration inequality in a regime where
a ≥ 9b, in which case KL(a||b) ≥ a ln(a/b) − 1 ≥ (a/2) ln(a/b).
def

Fix any constant 0 < c < 1. We apply the above concentration bound for t = kc−1 . Recalling the
assumption n ≥ k2 , we have that t + k/n ≥ 9k/n, and therefore


Pr k(v · v ′ ) ≥ k(k/n + t) ≤ exp (−k · KL (t + k/n||k/n))
≤ exp (−kt ln (tn/k))
 −kt 
tn
n −kc
=
=
k
k2−c
≤ n−ck

c /2

.

ckc /4 
c
Now if we let S be a set of ⌊nck /8 ⌋ unit vectors drawn independently from D, there are ⌊n 2 ⌋ <
c
⌊nck /2 ⌋ distinct pairs of v, v ′ ∈ S, and by a union bound the probability there exist distinct v, v ′ ∈ S with
c
c
(v · v ′ ) ≥ 2kc−1 is less than ⌊nck /2 ⌋n−ck /2 < 1. Thus, there exists a set S such that all distinct pairs
v, v ′ ∈ S satisfy |v · v ′ | ≤ 2kc−1 . This completes the proof.

Before we proceed with the proof of Theorem 6.6, we make a useful observation: By following the
proof of Proposition 3.3 using Lemma 6.7 instead of Lemma 3.7, mutatis mutandis, we obtain:
Proposition 6.8. Given a distribution A over R that satisfies Condition 3.2 for some m ∈ Z+ , and any
constant 0 < c < 1, consider the set of distributions Pv for v ∈ Sn that are k-sparse, with n ≥ max{2(m +
1) ln n, 8k2 }. For a given ǫ > 0, suppose that dTV (Pv , Pv′ ) > 2ǫ whenever |v·v ′ | is smaller than 1/8. Then,
any SQ algorithm which, given access to Pv (x) for an unknown k-sparse v ∈ Sn , outputs a p
hypothesis Q
c /8
ck
2(m+1)
−(m+1)(1/2−c)
with dTV (Q, Pv ) ≤ ǫ needs at least ⌊n
⌋≥k
queries to STAT(O(k)
χ2 (A, N (0, 1)))
(m+1)(1−2c)
2
or to VSTAT(O(k)
/χ (A, N (0, 1))).
The above proposition can be used for m = 1 to establish a similar but quantitatively somewhat weaker
SQ lower bound. We can make a crucial improvement to this proposition for the specific A we use in the
proof below.
Proof of Theorem 6.6. We select the one-dimensional distribution A as follows:
A = (1 − δ)N (ǫ, 1) + δN (−(1 − δ)ǫ/δ, 1) ,
where δ = ǫk−c/4 . Note that A has mean 0, i.e., matches m = 1 moments of N (0, 1).
2
2
2 2
We could use Facts 6.2 and 6.3 to obtain
√ χ (A, N (0, 1) ≤ O(ǫ exp(ǫ /δ )). However, this would
require the parameter δ to be equal to ǫ/ c ln k to get the required bounds from Proposition 6.8. The
issue here is that χ2 (A, N (0, 1)) is much bigger than the variance of A, which means that the correlation
inequality |χN (0,1) (Pv , Pv′ )| ≤ (v · v ′ )2 χ2 (A, N (0, 1)) is far from tight for most v and v ′ . For our choice
of A, we prove the following lemma:
35

Lemma 6.9. If A = (1 − δ)N (ǫ, 1) + δN (−(1 − δ)ǫ/δ, 1), then for v, v ′ ∈ Sn , we have

1 + |χN (0,1) (Pv , Pv′ )| ≤ exp ǫ4 (v · v ′ )2 /δ4 .
√
P
Proof. Let θ be the angle between v and v ′ . As in (3), we willP
use the expansion A(x) = √∞
i=0 ai Hei (x)G(x)/ i!.
i
As in (7), we also have the expansion Uθ A(x) = G(x) + ∞
i=2 ai cos θHei (x)G(x)/ i!. Thus, we can
write
χN (0,I) (Pv , Pv′ ) = χN (0,1) (A, Uθ A)
Z ∞
=
(A(x) − G(x))(Uθ A − G(x))/G(x)dx
−∞
! ∞
!
Z ∞ X
∞
X
√
√
ai Hei (x)G(x)/ i!
=
ai cosi θHei (x)G(x)/ i! /G(x)dx
−∞

=

∞
X

i=1

i=2

a2i cosi θ .

i=2

Since A is a distribution with mean zero, we have a0 = 1, a1 = 0. We need to take advantage of the fact
than χ2 (A, N (0, 1)).
that for our selected probability density function A, the coefficient a22 is much smaller
√
We can find the ai explicitly using (4), which gives that ai = EX∼A [Hei (x)/ i!]. We have the following
well-known fact:
Fact 6.10. For µ ∈ R, we have EX∼N (µ,1) [Hei (x)] = µi .
Proof. Note that the i-th derivative of G(x) isP
(−1)i Hei (x)G(x). Using Taylor’s theorem, we can expand
∞
i
G(x − µ) around x to obtain G(x − µ) =
i=0 µ Hei (x)G(x)/i! . Taking the expectation of Hei (x)
extracts the i-th term, establishing the fact.
Thus, we have
√
i!ai = EX∼A [Hei (x)]
= (1 − δ)EX∼N (ǫ,1) [Hei (x)] + δEX∼N (−(1−δ)ǫ/δ,1) [Hei (x)]

= (1 − δ)ǫi + δ(−(1 − δ)ǫ/δ)i .

√
In addition to a0 = 1,a1 = 0, we can derive the bound |ai | ≤ (ǫ/δ)i / i!. Recalling the special case a1 = 0
and summing over i, we have
|χN (0,I) (Pv , Pv′ )| ≤
≤

∞
X

i=1
∞
X
i=2

a2i | cos θ|i

(ǫ/δ)2i | cos θ|i /i!


= exp ǫ2 |v · v ′ |/δ2 − 1 − ǫ2 |v · v ′ |/δ2 .

To complete the proof of the lemma, it is sufficient to show that exp(x)−x ≤ exp(x2 ) for all x ≥ 0. We note
that both expressions are 1 and have derivative 0 at x = 0. It suffices to show that d2 (exp(x) − x)/dx2 ≤
d2 (exp(x2 )/dx2 for x ≥ 0. Note that
2

d2 ex d2 (ex − x)
2
2
2
2
/
= (2ex + 4x2 ex )e−x ≥ 2ex −x = 2e(x−1/2) −1/4 ≥ 2e−1/4 > 1 .
2
2
dx
dx
This completes the proof.
36

We now have all the necessary ingredients to complete the proof of Theorem 6.6. For distinct k-sparse
unit vectors v, v ′ ∈ S, where S is given by Lemma 6.7, we have that

|χN (0,1) (Pv , Pv′ )| ≤ exp ǫ4 (v · v ′ )2 /δ4 − 1

≤ exp ǫ4 4k2c−2 /δ4 − 1

≤ exp 4k2c−2 · kc − 1

= e4 exp k3c−2 − 1
≤ e4 k3c−2 .

Following the proof of Proposition 3.3, we now have that it takes Ω(nck
or to VSTAT(O(k)2−3c ) to learn Pv .

c /8

) queries to STAT(O(k)3c/2−1 )

7 Sample Complexity Lower Bounds for High–Dimensional Testing
In this section, we use our framework to prove information-theoretic lower bounds on the sample complexity
of our two high-dimensional testing problems: (i) robustly testing the mean of a single unknown mean
identity covariance Gaussian in Huber’s contamination model, and (ii) (non-robustly) testing between a
single spherical Gaussian and a mixture of 2 spherical Gaussians.
Both these statements follow from the structural results established in the previous sections using the
following proposition:
Proposition 7.1. Let A be a distribution on R such that A has mean 0 and χ2 (A, N (0, 1)) is finite. Then,
there is no algorithm that, for any n, given N < n/(8χ2 (A, N (0, 1))) samples from a distribution D over
Rn which is either N (0, I) or Pv , for some unit vector v ∈ Rn , correctly distinguishes between the two
cases with probability at least 2/3.
Proof. At a high-level, the proof of the proposition uses the structure of the set of Pv ’s and standard
information-theoretic arguments.
Suppose that, after fixing the dimension n, the algorithm takes at most N samples. We can consider the
testing algorithm as a (possibly randomized) function from N -tuples of samples to its output. For a distribution D, let D ⊗N denote the distribution over independent N -tuples drawn from D. We write f (D ⊗N )
for the Bernoulli distribution that gives the output of the algorithm given a single sample of D ⊗N . Let QN
be the distribution obtained by choosing v uniformly at random over the unit sphere Sn , and then drawing
N samples from Pv . Then, f (QN ) should be “NO” with probability at least 2/3, since the probability that
⊗N ) is “YES” with probability at least
each f (P⊗N
v ) is “NO” is at least 2/3. On the other hand, f (N (0, I)
2/3. By the data processing inequality, it follows that
dTV (QN , N (0, I)⊗N ) ≥ dTV (f (QN ), f (N (0, I)⊗N )) ≥ 1/3 .
Suppose for the sake of contradiction that N < n/(8χ2 (A, N (0, 1))). Then, we claim that
dTV (QN , N (0, I)⊗N ) < 1/3 .

37

Indeed, we have that:
4dTV (QN , N (0, I)⊗N )2 + 1 ≤ χ2 (QN , N (0, I)⊗N ) + 1 =
Z
Z
N
Y
(1)
(N ) 2
G(x(i) )dx(N ) . . . dx(1)
QN (x , . . . x ) /
···
=
=

=
=
≤

Z

Z

Z

x(1) ∈Rn

x(1) ∈Rn

v∈Sn

Z

Z

Zv∈Sn Z
v∈Sn

···

Z

v′ ∈Sn
v′ ∈S
v′ ∈S

x(N) ∈Rn

x(N) ∈Rn

Z

i=1

Z

x(1) ∈Rn

v∈Sn

···

Z

Z

v′ ∈S

(1)
(N )
(1)
(N )
PN
)PN
)/
v (x , . . . , x
v′ (x , . . . , x
n

N
Y

G(x(i) )dv ′ dvdx(N ) . . . dx(1)

i=1

N
Y

x(N) ∈Rn i=1

Pv (x(i) )Pv′ (x(i) )/G(x(i) )dx(N ) . . . dx(1) dv ′ dv

(1 + χN (0,I) (Pv , Pv′ ))N dv ′ dv
n

n

1 + |v · v ′ |2 χ2 (A, N (0, 1))

N

dv ′ dv ,

where the last line follows from Lemma 3.4, since A satisfies Condition 3.2 for m = 1. We will need the
following facts about the Beta function B(x, y):
R π/2
Fact 7.2.
(i) For x > −1, y > −1 we have that: 0 sinx (θ) cosy (θ) = B((x + 1)/2, (y + 1)/2)/2.
(ii) For all x, y ∈ R, we have that B(x, y + 1) = B(x, y) · y/(x + y).

Consider choosing v and v ′ independently uniformly at random over Sn . To obtain our tight sample
complexity results, we will need a more precise analysis for the distribution of the angle θ between v and
v ′ . Specifically, we show the following:
Lemma 7.3. If we choose v and v ′ uniformly at random from Sn , the angle θ between them is distributed
with the probability density function sinn−2 (θ)/B((n − 1)/2, 1/2), where B(x, y) is the Beta function.
Proof. Since a rotation of the sphere moves both v and v ′ , θ is invariant under such rotations. Thus, we
get the same distribution by fixing v ′ = e1 , the unit vector in the x1 -direction, and choosing v uniformly at
random over the sphere. Now we have that cos θ = x1 . Let Sn (r) denote the surface area of the sphere of
radius r in (n + 1) dimensions and note that Sn (r) = r n Sn (1).
For any measurable function f ,we have that
Z π
Z
Z
f (cos θ) sinn−2 (θ)Sn−2 (1)dθ .
f (x1 )Sn−2 (x2 )ds =
f (x1 )dΣ =
0

S2

Sn

We thus have that the pdf of θ is proportional to sinn−2 (θ). Taking f (x1 ) ≡ 1, note that
Z π
Z
sinn−2 (θ)Sn−2 (1)dθ = B((n − 1)/2, 1/2) ,
1dΣ =
Sn−2 (1) =
Sn

0

using Fact 7.2 (i). We thus have that the pdf of θ is
sinn−2 (θ)/Sn−2 (1) = sinn−2 (θ)/B((n − 1)/2, 1/2) .
This completes the proof.

38

We now have that
4d2TV

⊗N

(QN , N (0, I)

)+1≤

Z

Z

v′ ∈S

1 + |v · v ′ |2 χ2 (A, N (0, 1))

N

dv ′ dv

n
n
Zv∈S
π
N
=
1 + | cos θ|2 χ2 (A, N (0, 1)) sinn−2 (θ)/B((n − 1)/2, 1/2)dθ

0

=2

Z

π/2

1 + cos2 (θ)χ2 (A, N (0, 1))

0

N

sinn−2 (θ)/B((n − 1)/2, 1/2)dθ

N   Z π/2
X
N
χ2 (A, N (0, 1))i cos2i (θ) sinn−2 (θ)/B((n − 1)/2, 1/2)dθ
2
=
i
0
i=0
N  
X
N 2
=
χ (A, N (0, 1))i B((n − 1)/2, i + 1/2)/B((n − 1)/2, 1/2) .
i
i=0

PN

We rewrite the above sum as i=0 bi , where
 
N 2
bi =
χ (A, N (0, 1))i B((n − 1)/2, i + 1/2)/B((n − 1)/2, 1/2) .
i
N  N
/ i = (N − i)/(i + 1) and using
Note that b0 = 1. Now consider the ratio bi+1 /bi . Note that that i+1
Fact 7.2, we have that B((n − 1)/2, i + 3/2)/B((n − 1)/2, i + 1/2) = (i + 1/2)/(i + n/2). Therefore, it
follows that
bi+1 /bi = χ2 (A, N (0, 1)) · (N − i)/(i + n/2) · (i + 1/2)/(i + 1) .
For all i, we have

bi+1 /bi ≤ χ2 (A, N (0, 1)) · N/(n/2) = 2N χ2 (A, N (0, 1))/n .

When N ≤ n/(8χ2 (A, N (0, 1))), we have bi+1 /bi ≤ 1/4 for i ≥ 0, and therefore
N
X
i=0

Hence, we have

bi ≤ 1/(1 − 1/4) = 4/3 .

4d2TV (QN , N (0, I)⊗N ) + 1 ≤ 1/3 < 4/9 .

This implies that dTV (QN , N (0, I)⊗N ) < 1/3, which is the desired contradiction. This completes the
proof.
Remark 7.4. It is worth noting that matching m > 1 many moments does not seem to help in the setting
of the previous proposition, as long as χ2 (A, N (0, 1)) ≤ 1. This may seem to some extent unsurprising,
given that O(n) samples suffice for some of the learning problems we consider here. On the other hand,
we consider it somewhat surprising looking at the proof of Proposition 7.1. Specifically, for general m, we
would have that
N  
X
N 2
2
⊗N
4dTV (QN , N (0, I) ) + 1 ≤
χ (A, N (0, 1))i B((n−1)/2, i(m+1)/2+1/2)/B((n−1)/2, 1/2) .
i
i=0

Note that the ratio of one term to the next approximately grows as N χ2 (A, N (0, 1))i(m−1)/2 /n(m+1)/2 . For
this to be less than 1/2, for all 0 ≤ i ≤ N , we need N (m+1)/2 χ2 (A, N (0, 1)) ≤ O(n(m+1)/2 ). Thus, we
need at least N = Ω(n/(χ2 (A, N (0, 1)))2/m ) samples . This suggests that we should be able to obtain a
tighter lower bound if χ2 (A, N (0, 1)) > 1 using this technique. We omit the details here, as we are mainly
interested in the regime χ2 (A, N (0, 1)) ≤ O(1) for our applications in this paper.
39

Using Proposition 7.1, we establish the two main results of this section:
Theorem 7.5 (Sample Complexity Lower Bound for Robustly Testing Unknown Mean Gaussian). There is
no algorithm that, for every ǫ > 0 and positive integer n, given fewer than Ω(n) samples from a distribution
P on Rn which is promised to satisfy either (a) P = N (0, I) or (b) P = (1 − δ)N (µ, I) + δN1 , where
δ = ǫ/100, kµk2 ≥ ǫ and the noise distribution N1 is a spherical Gaussian, can distinguish between the
two cases with probability at least 2/3.
If instead, for any constant 0 < c < 1, we are promised that P = (1 − δ)N (µ, I) + δN1 , where
δ = ǫ/nc/4 in case (b), then no algorithm that takes less than Ω(n1−c ) samples can distinguish between (a)
and (b) with probability at least 2/3.
Proof. Let δ be the noise rate. We will take δ = ǫ/100 or δ = ǫ/nc/4 . In both cases, we select our
one-dimensional distribution to be the following:
A = (1 − δ)N (ǫ, 1) + δN (−(1 − δ)ǫ/δ, 1) .
Now Lemma 6.9 yields that for v, v ′ ∈ Sn , we have


1 + |χN (0,1) (Pv , Pv′ )| ≤ exp ǫ4 (v · v ′ )2 /δ4 .

We will not apply Proposition 7.1 directly but follow its proof using the aforementioned stronger correlation bound. We have:
Z
Z
(1 + χN (0,I) (Pv , Pv′ ))N dv ′ dv
4dTV (QN , N (0, I)⊗N )2 + 1 ≤
′
Zv∈Sn Zv ∈Sn

exp N ǫ4 (v · v ′ )2 /δ4 dv ′ dv
≤
′
n v ∈Sn
Zv∈S
π

exp N ǫ4 (cos θ)2 /δ4 sinn−2 (θ)/B ((n − 1)/2, 1/2) dθ
=
0

=

∞
X
i=0

(N ǫ4 /δ4 )i B((n − 1)/2, i + 1/2)/i!B((n − 1)/2, 1/2) .

Now note that the ratio of the (i + 1)-th term to the i-th term of the corresponding series is
2N ǫ4
N ǫ4 /δ4 i + 1/2
·
≤
.
i + 1 i + n/2
nδ4
When N ≤ nδ4 /2ǫ4 , the 0-th term is 1 and the ratio of the (i + 1)-th to i-th term is less than 1/4. Therefore,
the above sum is less than 1/(1 − 1/4) = 4/3, which implies that
√
dTV (QN , N (0, I)⊗N ) ≤ 1/ 12 < 1/3 .
Following the proof of Proposition 7.1, we conclude that no algorithm satisfying the necessary conditions
exists. To complete the proof, note that for δ = ǫ/100, we need at least Ω(n) samples. And for δ = n−c/4 ǫ,
we need at least Ω(n1−c ) samples.
Theorem 7.6 (Sample Complexity Lower Bound for Testing GMMs). Fix 0 < ǫ < 1 and n ∈ Z+ . There is
no algorithm that, given less than Ω(n/ǫ2 ) samples from a distribution P on Rn that is promised to be either
(a) P = N (0, I), or (b) P is a mixture of two Gaussians each with weight 1/2 and identity covariance, such
that dTV (P, N (0, I)) ≥ ǫ, distinguishes between the two cases with probability at least 2/3.
40

Proof. We choose our one-dimensional distribution as
A = (1/2)N (−δ, 1) + (1/2)N (δ, 1) ,
where we set (with hindsight) δ = Θ(ǫ1/2 ).
Note that A has mean 0. By Claims 6.2 and 6.3, we have that


χ2 (A, N (0, 1)) = (1/2) exp(δ2 )/2 + exp(−δ2 )/2 − 1 = (1/2) cosh(δ2 ) − 1 = Θ(δ4 ) .

Similarly, it can be shown that dTV (A, N (0, 1)) = Θ(δ2 ), which is Ω(ǫ) by our choice of δ. Applying
Proposition 7.1 completes the proof.

8 SQ Algorithms for Robustly Learning and Testing a Gaussian
The structure of this section is as follows: In Section 8.1, we prove a moment–matching structural result that
forms the basis of our algorithms. In Section 8.2, we present our robust testing algorithm, and in Section 8.3
we give our robust learning algorithm.
8.1 One-Dimensional Moment Matching Lemma The main result of this section is the following structural result:
p
Lemma 8.1. Let G ∼ N (0, 1). For δ > ǫ > 0, define k = 2⌈ǫ ln(1/ǫ)/δ⌉. Let G′ be an ǫ-noisy onedimensional Gaussian with unit variance so that for all t ≤ k we have that the tth moments of G and G′
agree to within an additive (t − 1)!(δ/ǫ)t ǫ/t. Then, we have that dTV (G, G′ ) = O(δ).
e = N (µ, 1) be such that dTV (G′ , G)
e ≤ ǫ. We can assume without loss of generality that µ > 0.
Proof. Let G
p
Looking at just the mean andp
variance suffices to get dTV (G, G′ ) = O(δ log(1/δ)) using techniques
similar to the proof of the O(ǫ log(1/ǫ)) filter algorithm from [DKK+ 16]. This allows us to focus on the
case that µ ≥ 1. Formally, we have the following claim:
e is at least 1.
Claim 8.2. Lemma 8.1 holds when the mean µ of G

Proof. We start by noting that we can assume ǫ > 0 is smaller that a sufficiently small universal constant.
Assuming otherwise and recalling that δ > ǫ gives that dTV (G, G′ ) ≤ 1 = O(δ), in which case the lemma
statement is trivial.
Let µ′ be the mean of G′ . We can write
e + ǫ′ E − ǫ′ L ,
G′ = G

(12)

e ≥ ǫ′ L .
G

(13)

e ≤ ǫ. Moreover, it holds
for distributions E, L with disjoint supports where ǫ′ = dTV (G′ , G)

Since k ≥ 2 by definition, the lemma assumptions imply that |µ′ | ≤ δ and EX∼G′ [X 2 ] ≤ δ2 /ǫ. By (12)
we have that µ′ = µ + ǫ′ EX∼E [X] − ǫ′ EX∼L [X], and similarly EX∼G′ [(X − µ′ )2 ] = 1 + (µ − µ′ )2 +
ǫ′ EX∼E [(X − µ′ )2 ] − ǫ′ EX∼L [(X − µ′ )2 ]. Therefore, we get
ǫ′ EX∼E [(X − µ′ )2 ] ≤ δ2 /ǫ − 1 − (µ′ − µ)2 + ǫ′ EX∼L [(X − µ)2 ] .

41

Note that PrX∼L [|X − µ′ | ≥ T ] ≤ PrX∼Ge [|X − µ′ | ≥ T ]/ǫ′ . As in Corollary 8.8 of [DKK+ 16], we
have that
Z ∞
′ 2
PrX∼L [|X − µ′ | ≥ T ]2T dT
EX∼L [(X − µ ) ] =
0
Z ∞

2T min 1, PrX∼N (µ,1) [|X − µ′ | ≥ T ]/ǫ′ dT
≤
Z0 ∞
2T min{1, exp(−((T − |µ − µ′ |)2 /2)2T /ǫ′ }dT
≤
0
Z √2 ln(1/ǫ′ )+|µ−µ′ |
Z ∞
e
=
2T dT + √
exp(−T 2 /2)2(T + |µ − µ′ |)/dTV (G′ , G)dT
2 ln(1/ǫ′ )

0

′

′ 2

= O(ln(1/ǫ ) + 1 + |µ − µ | ) .

By combining the above, we get that
ǫ′ EX∼E [(X − µ′ )2 ] ≤ δ2 /ǫ − 1 − (1 − O(ǫ′ ))|µ − µ′ |2 + O(ǫ′ ln(1/ǫ′ )) ,
and thus EX∼E [(X − µ′ )2 ] ≤ δ2 /(ǫǫ′ ).
′ 2
′ 2
′ 2
However, the means of L and E, µL and µE have
p |µL − µ | ′ ≤′ EX∼L [(X −pµ ) ] and |µE′ − µ | ≤
′
2
′
′
′
EX∼E [(X − µ ) ],p
and therefore ǫ |µL − µ | ≤ O(ǫ ln 1/ǫ′ + ǫ |µ − µ|) ≤ O(ǫ ln 1/ǫ + ǫ|µ − µ|) and
′
′
ǫ |µE − µ| ≤ O(ǫ δ2 /ǫǫ′ ) = O(δ).
p
We thus have that µ′ − µ = ǫ′ µE − ǫ′ µL satisfies |µ′ − µ| ≤ O(δ) + O(ǫ ln 1/ǫ) + O(ǫ|µ − µ′ |). Since
ǫ is sufficiently small, we have |µ′ − µ| ≤ O(δ) + 1/2. Since |µ′ | ≤ δ, We thus have that µ ≤ O(δ) + 1/2.
Since we assumed that µ ≥ 1, we get that δ = Ω(µ − 1/2) = Ω(1/2) and so µ = O(δ). We therefore
conclude that dTV (G, G′ ) = O(µ) = O(δ).
e = Θ(µ). Suppose for the sake
We will henceforth assume that µ ≤ 1 and thus we have that dTV (G, G)
′
e ≥
of contradiction that dTV (G, G ) = Ω(δ), for a sufficiently large constant in the big-Ω. Then, dTV (G, G)
2
e
dTV (G, G) − ǫ ≫ δ. We may assume that µ ≥ C δ for a sufficiently large constant C > 0. Thus, we have
that ǫ ≤ µ/C 2 ≤ 1/C 2 .
The proof will proceed as follows: Let f (x) = sin(xCǫ/µ). We note that f (x) has a simple expectation
e and we can easily get a lower bound on their difference. We will also use the Taylor series
under G or G,
for f (x) and our moment bounds to derive an upper bound on this difference which contradicts this lower
bound.
For x ∈ R+ , we want an expression for EX∼N (x,1) [f (X)]. By standard facts on the Fourier transform,
we have that
EX∼N (x,1) [exp(−iωX)] = exp(−iωx)EX∼N (0,1) [exp(−iωX)] = exp(−ω 2 /2 − iωx) ,
for any ω ∈ R. Since sin(xCǫ/µ) = (exp(−ixCǫ/µ) − exp(ixCǫ/µ))/2i, we obtain that
EX∼N (x,1) [f (X)] = exp(−(Cǫ/µ)2 /2) sin(xCǫ/µ) .
Therefore, EX∼N (0,1) [f (X)] = 0 and EX∼N (µ,1) [f (X)] ≥ exp(1/C) sin(Cǫ) > (C/2)ǫ, and thus
EX∼G′ [f (X)] ≥ EX∼Ge [f (X)] − ǫ > (C/3)ǫ.

(14)

Let h be the degree-(k − 1) Taylor polynomial of f plus the term (Cxǫ/µ)k /k!. By the Lagrange form of
the remainder in Taylor’s theorem, we have that h(x) − (Cxǫ/µ)k /k! − f (x) ≤ f (k) (ξ)xk /k!, for some
42

ξ ∈ [0, x]. Since k is even, we have that the k-th derivative of f , |f (k) (ξ)| = (Cǫ/µ)k | sin(ξCǫ/µ)| ≤
(Cǫ/µ)k . Thus, we get
f (x) ≤ h(x) ≤ f (x) + 2(xCǫ/µ)k /k! .
Our goal will be to show that EX∼G′ [h(X)] is substantially larger than EX∼N (0,1) [h(X)], which will contradict the assumption about approximately matching moments. We start by considering EX∼N (µ,1) [h(X)]
versus EX∼N (0,1) [h(X)]. We can write
EX∼N (µ,1) [h(X)] − EX∼N (0,1) [h(X)]

h
i
h
i

= EX∼N (µ,1) [f (X)] − EX∼N (0,1) [f (X)] + O EX∼N (µ,1) (XCǫ/µ)k /k! − EX∼N (0,1) (XCǫ/µ)k /k!
h
i
h
i

,
≥ (C/3)ǫ + O EX∼N (µ,1) (XCǫ/µ)k /k! − EX∼N (0,1) (XCǫ/µ)k /k!

where the last inequality follows from (14). To bound this latter term, we make the following claim:
Claim 8.3. We have that

√
|G(x − µ) − G(x)| ≤ O(µ)G(x/ 2).

Proof. First, we note that G(x − µ)/G(x) = exp(−2xµ + µ2 /2).
Recalling our assumption that 0 ≤ µ < 1, for |x| ≤ O(1/µ),
√ we have that |G(x − µ) − G(x)| ≤
2
O(|x|µ + µ )G(x) ≤ O((|x| + 1)µ)G(x).
√ Then, since G(x)/G(x/ 2) = O(G(x)) ≤ O(1/(|x| + 1)), we
get |G(x − µ) − G(x)| ≤p
O(µ)G(x/ 2).
For |x| ≥ 5/µ,√since 2 ln(1/µ) + 2 ≤ 2 ln(1/µ) + 3 ≤ 2/µ + 3 ≤ x, we have that G(|x| − 2) ≤ µ.
Thus, G(x)/G(x/ 2) ≤ O(G(x)) ≤ O(µ) and
√
G(x − µ)/G(x/ 2) = O(exp(−x2 /4 + µx − µ2 /2))
√
= O(G((x − 2µ)/ 2) exp(µ2 /2))
≤ O(G(|x| − 2)) ≤ O(µ) ,

√
and hence |G(x − µ) − G(x)| ≤ O(µ)G(x/ 2).
Using the previous claim, we note that
Z

∞

−∞

√
(xCǫ/µ)k /k! |G(x − δ) − G(x)| dx = O(ǫ)EX∼N (0,1/2) [(Cǫ/µ)k /k!] = O(ǫ/(C/ 2)k ) = O(ǫ).

Therefore, we have that
EX∼Ge [h(X)] ≥ EX∼G [h(X)] + (C/3)ǫ.

Next, we wish to compare EX∼G′ [h(X)] to EX∼Ge [h(X)]. Using (12), we have the following for ǫ′ =
e
dTV (G′ , G):
EX∼G′ [h(X)] = EX∼Ge [h(X)] + ǫ′ (EX∼E [h(X)] − EX∼L [h(X)])
h
i

≥ EX∼Ge [h(X)] + ǫ′ EX∼E [f (X)] − EX∼L [f (X)] − EX∼L (Xǫ/δ)k /k!
h
i
= EX∼Ge [h(X)] + O(ǫ) − ǫ′ EX∼L (XCǫ/µ)k /k! .

e ≥ ǫ′ L, it follows that L satisfies the concentration inequality
From (13), i.e., G
Pr [|X − µ| > T ] ≤ 2 exp(−T 2 /2)/ǫ′ .

X∼L

43

We now proceed to bound the subtractive term from above:
h
i
ǫ′ EX∼L (XCǫ/µ)k /k!
≤ ǫ′

∞
X

J=0

Pr

X∼L

h

|X − µ| ≥ 2J

p

log(1/ǫ)

i

k
p
2(J + 1) log(1/ǫ) + µ (Cǫ/µ)k /k!

∞
k
 p
X
p
2ǫJ (2(J + 2) log(1/ǫ)Cǫ/µ)k /k!
≤ ǫ 3 log(1/ǫ)Cǫ/µ /k! +
J=1

≤ǫ·O

p

log(1/ǫ)Cǫ/µ

k

/k! ·

1+

∞
X

ǫJ (J + 2)k

J=1

!

k √
p
log(1/ǫ)ǫ/Cδ /k! ≤ ǫ · O
log(1/ǫ)ǫ/Cδ /( 2πk(k/e)k )
≤ǫ·O
p
√
√
≤ ǫO( log(1/ǫ)ǫ/Cδk)k / k ≤ O(ǫ/ kC k ) ≤ O(ǫ) .
p

k

Therefore, we conclude that

EX∼G′ [h(X)] ≥ EX∼N (µ,1) [h(X)] + O(ǫ) ≥ EX∼N (0,1) [h(X)] + (C/4)ǫ .
On the other hand, P
recalling that h is the degree-(k − 1) Taylor expansion of f (x) = sin(xCǫ/µ), we
k−1
i
i
i
can write h(x) =
i=0 ai x , with ai = O((Cǫ/µ) /i!) = O((ǫ/Cδ) /i!). Therefore, the difference
EX∼G′ [h(X)] − EX∼N (0,1) [h(X)] is the sum over i of ai times the difference in the ith moments, which
by assumption is at most
ǫ

k−1
k−1
X
X
(1/i)2 = O(ǫ) .
(i − 2)!(δ/Cǫ)i O((ǫ/δ)i /i!) = O(ǫ)
i=1

i=0

This contradicts the fact that their difference is at least (C/4)ǫ, and concludes the proof.
8.2 Robust Testing Algorithm In this subsection, we give a robust testing algorithm, i.e., an algorithm
that distinguishes between an ǫ-noisy Gaussian and N (0, I). This algorithm will form the basis for our
robust learning algorithm of the following subsection.
Theorem 8.4 (Robust Testing Algorithm). Let G′ be an ǫ-noisy version of an n-dimensional Gaussian
with identity covariance. Let δ be at least a sufficiently large constant multiple of ǫ. There
p exists an SQ
algorithm that makes O(nk ) queries to STAT(ǫ · O(n log(n/ǫ)2 )−k ) where k = 2⌈O(ǫ log(1/ǫ)/δ)⌉,
and distinguishes between the cases that G′ is the standard normal distribution N (0, I), and the case that
G′ is at least δ-far from N (0, I). The algorithm has running time nO(k) .
By simulating the statistical queries with samples, we obtain:
Corollary 8.5. Given sample access to G′ , an ǫ-noisy version of an n-dimensional Gaussian with identity
covariance and ǫ, δ > 0 with δ be at least a sufficiently large constant multiple of ǫ, there is an algorithm that
with probability 9/10 distinguishes between the cases that G′ is the standard normal distribution N (0, I),
and the case that G′ is at leastp
δ-far from N (0, I) and requires at most (n log(1/ǫ))O(k) /ǫ2 samples and
running time where k = 2⌈O(ǫ log(1/ǫ)/δ)⌉.
e = N (µ, I) be such that dTV (G′ , G)
e ≤ ǫ. We need
Proof. In the case when G′ is δ-far from N (0, I), let G
′
′
to show that we can distinguish between the cases G = N (0, I) and G is an ǫ-noisy version of a Gaussian
44

e = N (µ, I) with dTV (G,
e G) ≥ δ − ǫ. We assume from now on that the completeness case is to show that
G
e
dTV (G, G) ≥ δ, since replacing δ with δ + ǫ does not affect the statement of the theorem.
The algorithm is quite simple. Let C be a sufficiently large universal constant such that the O(δ) total
variation distance bound in Lemma 8.1 is less than Cδ. We assume that δ > 4Cǫ.
Robust Testing Algorithm:
• For each coordinate axis 1 ≤ i ≤ n, use STAT to approximate PrX∼G′ [X ≤ ǫ] and PrX∼G′ [X ≥ ǫ]
to within ǫ/2. If any of these approximations are bigger than 1/2 + ǫ, then output “NO”.
p
• Let k = 2⌈2Cǫ ln(1/ǫ)/δ⌉. Let C ′ be a sufficiently large constant.
p
• Using access to STAT, find all the mixed moments of X ∼ G′ , conditioned on kXk2 ≤ C ′ k n log(n/ǫ),
of order at most k to within error n−k/2 ǫ/2.
• If the difference between any moment of order t ≤ k that we measured and that of N (0, I) is more
than ((t − 1)!(δ/2Cǫ)t /t − 1) · n−k/2 ǫ, then output “NO”.
• Otherwise, output “YES”.
The idea is to use Lemma 8.1 with the approximations the moments. However, we have the issue that
the STAT oracle
p can only be used to approximate the expectation of a bounded function. Using the condition
kXk2 ≤ C ′ k n log(n/ǫ) allows us to avoid this. But we first need to show that conditioning on it does not
affect the moments too much and does not move the distribution far in total variational distance.
Note that the first step of the algorithm will reject if the median of N (µ, I) projected onto any coordinate
axis is outside of the interval [−ǫ, ǫ]. If this occurs, then µ 6= 0. If this step does not reject, then µ projected
√
onto any coordinate axis is O(ǫ), and
so kµk2 ≤ O(ǫ n).
p
The condition kxk2 ≤ C ′ k n log(n/ǫ) ensures that the any degree less than k monomial in X
is at most (C ′ nk2 log(n/ǫ)))k/2 . Thus, we can approximate this expectation to precision n−k/2 ǫ using
√
STAT(ǫ · O(n log(n/ǫ)2 )−k ). However, since kµk2 ≤ O(ǫ n), by standard concentration bounds, the
probability that X ∼ N (µ, I) does
not satisfy this condition is at most ǫ · (C ′ nk2 log(n/ǫ)))−k . Let G′′ be
p
G′ conditioned on kXk2 ≤ C ′ k n log(n/ǫ). If G′ = N (0, I), then we need to show that the moments of
G′′ and G′ are within ln−k/2 (ǫ/2). To show this, we use the following lemma:
Lemma 8.6. For 0 ≤ δ ≤ exp(−k), the difference
between any mixed moment of degree at most k of
p
N (0, I) and N (0, I) conditioned on kxk2 ≤ O( nk log(1/δ)), is at most δ.

Proof. Let X be distributed as N (0, I). Let ǫ be δ/(k + ln 1/δ)(k−1) C for ap
sufficiently large constant
C. Thus,
ln(1/ǫ) = ln(1/δ) + ln C + (k − 1) ln(k + ln(1/δ)). Let T be 2n log 1/ǫ. Then, T =
p
OC ( nk log(1/δ)) and by standard concentration inequalities, we have that Pr[kXk2 ≥QT ] ≤ ǫ.
For a ∈ Nn with kak1 ≤ k, consider the monomial of degree
ma (x) = ni=1 xai i . First we
Qnat most k,
E(Xiai ) = 0, sinceQthe odd moments
consider its mean and variance. If any ai is odd, E(ma )(X) = i=1 Q
n
ai
ai /2 (a /2)! ≤
of Xi ∼ N (0, 1) are zero. If all ai are even, then E(ma )(X) = Q
i
i=1 E(Xi ) =
i2
Q
2a
n
a a ! ≤ 2k k!.
2] =
i
E(X
)
=
2
2k/2 (k/2)!. For the variance,
we
have
Var[m
(X)]
≤
E[m
(X)
i
a
a
i
i=1
i
i
√
Let pa (X) = ma (X)/2k/2 k!. Then we have that 0 ≤ E[pa (X)] ≤ 1 and Var[pa (X)] ≤ 1.
By the standard concentration inequality given in Lemma 8.16 below, we have that for all t > 0,
Pr[|pa (X)| ≥ t + 1] ≤ exp(2 − (t/R)2/k ) for some R > 0. Thus, we have Pr[|pa (X)| ≥ c + 1] ≤ ǫ, for

45

c = R(ln(1/ǫ) − 2)k/2 . Let I(x) be the indicator function of kxk2 ≥ T . Then we have that
|E[I(X)pa (X)]| ≤ E[I(X)|pa (X)|]
Z ∞
Pr[I(X)|pa (X)| ≥ t]dt
=
0
Z ∞
Z c+1
exp(2 − (t − 1/R)2/k )dt
ǫdt +
≤
c+1
0
Z ∞
exp(2 − (t/R)2/k )dt
= (c + 1)ǫ +
Zc ∞
= (c + 1)ǫ +
exp(2 − x)(dt/dx)dx
ln(1/ǫ)
Z ∞
exp(2 − x)xk/2−1 dx
= (c + 1)ǫ + (Rk/2)

(where x = (t/R)2/k )

ln(1/ǫ)

k/2−1

= (c + 1)ǫ + (Rk/2)ǫ ·
k/2

≤ R ln(1/ǫ)
where the integral
ma (x), we have

R∞

2

X
j=0

((k/2)!/(k/2 − j)!) ln(1/ǫ)k/2−j

ǫ + (k /8)ǫ(k + log(1/ǫ))k/2−1

≤ O(k2 ǫ(k + log(1/ǫ))k/2−1 ) ,

ln(1/ǫ) exp(2

− x)xk/2−1 dx is calculated explicitly below in Claim 8.18. In terms of

√
|E[I(X)ma (X)]| ≤ O(k2 ǫ(k + log(1/ǫ))k/2−1 2k/2 k!) ≤ O(ǫ(k + log(1/ǫ)(k−1) ) ≤ O(δ/C) ≤ δ/2 .
Then, for X ′ distributed as N (0, I) conditioned on kX ′ k2 ≤ T , we have
|E[ma (X)] − E[ma (X ′ )]| = |E[ma (X)] − E[ma (X)(1 − I(X))]/(1 − Pr[kXk2 > T ])|

= |(E[I(X)ma (X)] − E[ma (X)] Pr[kXk2 > T ]) /1 − Pr[kXk2 > T ])|

√ 
≤ δ/2 + 2k/2 kǫ /(1 − ǫ)
≤ 2δ/3(1 − ǫ) ≤ δ .

p
p
Applying Lemma 8.6 for δ = n−k/2 ǫ, noting that C ′ k n log(n/ǫ) = Ω(C ′ nk log(1/δ)) yields that
the moments of G′′ and G′ = N (0, I) are within ln−k/2 (ǫ/2). Thus, in this case, the approximations of the
moments of G′′ are within n−k/2 ǫ of the moments of G′ .
For the soundness case, we just note that since (t − 1)!(δ/2Cǫ)t ǫ/t − 1 ≥ (δ/2Cǫ)/2 − 1 ≥ 1, the
bounds on the moments we need to fail are bigger than the precision of the statistical queries we use to
approximate them, and therefore we never output “NO” when G′ = N (0, I).
e Then G′′ is a 2ǫ-noisy
Now suppose that G′ is an ǫ-noisy version of an identity covariance Gaussian G.
e
e
version of G. We will denote µ the mean vector of G and will assume that kµk2 ≥ δ. We need to show that
the algorithm outputs “NO”.
Consider the unit vector v = µ/kµk2 which has v · µ ≥ δ. Consider the one-dimensional distributions
ev of the form v · X, where either X ∼ G′′ or X ∼ G.
e Note that G
ev = N (kµk2 , 1) has mean larger
G′′v and G
′′
ev ) ≤ 2ǫ. We can now apply the contrapositive of Lemma 8.1 with δ/C in place
than δ and that dTV (Gv , G
46

of δ and 2ǫ in place of ǫ, which implies that there is a t ≤ k such that the t-th moment of G′v is more than
(t − 1)!(δ/2Cǫ)t ǫ/t far from that of N (0, 1). That is,
|EX∼G′ [(v · X)t ] − EX∼N (0,I) [(v · X)t ]| ≥ (t − 1)!(δ/2Cǫ)t ǫ/t .
Q ai
n
Now consider the polynomial (v·x)t . Note that the
coefficient
of
a
monomial
i xi , for a ∈ Z>0 , kak1 = t,
Q
n
ai
t
is given by the multinomial theorem as a1 ,...,a
i=1 vi . The L1 -norm of its coefficients is the same as the
n
Q
L1 -norm of entries of the rank-t tensor v ⊗t , which has i-th entry tj=1 vij , for i ∈ {0, . . . , n}t , since there

Q
t
entries in this symmetric tensor which are given by ni=1 viai for any a. The Frobenius norm of
are a1 ,...,a
n
Q
P
P
Q
v ⊗t , the L2 -norm of its entries, is i∈{0,...,n}t tj=1 vi2j = tj=1 ni=1 vi2 = 1. Since there are nt entries,

the L1 -norm of its entries must be at most nt/2 ≤ nk/2 . We can write EX∼G′′ [(v·X)t ]−EX∼N (0,I) [(v·X)t ]
as a linear combination of differences in moments with these coefficients, and we can thus bound this from
above by the product of the L1 -norm of
Qthe coefficients and the L∞ -norm of the differences in moments.
Hence, there must be some moment E[ i Xiai ] with a ∈ Zn>0 , kak ≤ k, such that
"
"
#
#
Y a
Y a
EX∼G′
Xi i − EX∼N (0,I)
Xi i ≥ (t − 1)!(δ/2Cǫ)t ǫ/t · n−k/2 .
i

i

This in turn means that the difference in the approximation of this moment of G′′ and that of N (0, I) is at
most ((t − 1)!(δ/2Cǫ)t /t − 1) · n−k/2 ǫ. Thus, the testing algorithm outputs “NO”.
8.3 Robust Learning Algorithm In this section, we build on the testing algorithm of the previous section
to design our robust learning algorithm. Formally, we prove:

Theorem 8.7. Let G′ be an ǫ-noisy version of an n-dimensional Gaussian with identity covariance matrix,
N (µ, I) with kµk2 ≤ poly(n, ǫ). There is an algorithm that, given statistical query access√to G′ , outO( log(1/ǫ))
+√
puts an approximation
µ
e to the mean µ such that√
kµ − µ
ek2 ≤ O(ǫ). The algorithm uses
√
√ n
O(

2log(1/ǫ)

log(1/ǫ))

calls to STAT(ǫ/(n ln(1/ǫ))O(

log(1/ǫ))

) and has running time nO(

log(1/ǫ))

O(

+2log(1/ǫ)

log(1/ǫ))

By simulating the statistical queries with samples, we obtain:

Corollary 8.8. Given sample access to G′ , an ǫ-noisy version of an n-dimensional Gaussian N (µ, I),√there
is an algorithm that
µ
e with kµ−e
µk2 ≤ O(ǫ) and requires (n log(1/ǫ))O( log(1/ǫ)) /ǫ2
√
√ with probability 9/10 outputs
samples and nO(

log(1/ǫ))

O(

/ǫ2 + 2log(1/ǫ)

log(1/ǫ))

time.

The work [DKK+ 16] gives algorithms which can compute an approximation µ′ with kµ − µ′ k2 ≤
p
O(ǫ log(1/ǫ)). These algorithms can be expressed as Statistical Query algorithms. However, due to the
model of adversary used for robustness in [DKK+ 16], the algorithms were expressed there in terms of operations on sets of samples that were drawn before the execution of the algorithm. The filtering algorithms
work by successively removing samples from this set and then computing expectations of the current set of
remaining samples. The samples that are removed are those that satisfy an explicit condition, we say that
they are rejected by a filter. We can implement these algorithms as SQ algorithms by replacing expectations
of the current set of remaining samples with the conditional expectation of the input distribution, conditioned
on all previous filters accepting. This is similar to the filtering algorithm for learning binary Bayesian networks given in [DKS16c]. Even there, we still used samples to compute the threshold for the filter. We note
that using arguments similar to those we use for the algorithm below, all theses algorithms can be expressed
as SQ algorithms. In particular, this is the case for Algorithm Filter-Gaussian-Unknown-Mean,
p
which we will use as a black box pre-processing step to approximate the mean within O(ǫ log(1/ǫ)).
Instead of dealing with moments, i.e., the expectations of monomials, directly, we will consider expectations of Hermite polynomials, which have a simpler form for normal distributions.
47

.

Definition Q
8.9. We define multi-dimensional normalized
Hermite polynomials as follows: for a ∈ Zn ,
Q
Hea (x) = ni=1 Heai (xi ). We define n(a) = ni=1 ai !.
Thus we have the following:

Fact 8.10. EX∼N (0,I) [Hea (X)Heb (X)] = δab n(a).
We define a linear combination of these Hermite polynomials with degree t associated with a tensor
of rank t. For i ∈ {0, . . . , n}t , we define the count vector c(i) ∈ Zn≥0 such that c(i)j is the number of
coordinates of i that are j.
Definition 8.11. For a rank-t tensor A over Rn , we define
X
√
hA (x) =
Hec(a) (x)/ t! .
i∈{0,...,n}t

We are now ready to describe our learning algorithm.
Robust Learning Algorithm:
p
1. Let k = 2⌈ ln(1/ǫ))⌉.

p
2. Compute an approximation µ′ with kµ′ − µk2 ≤ O(ǫ log(1/ǫ)) by iterating Algorithm
Filter-Gaussian-Unknown-Mean from [DKK+ 16]. We change the origin so that µ′ = 0.
p
3. Let N be the filter that accepts when kxk2 ≤ 2n log(1/ǫ).
√
4. For 1 ≤ t ≤ k, let Pet be the rank-t tensor with i1 , . . . , it entry given by t! times the result of asking
an SQ oracle for EX∼G′ [hc(i) (X)] conditioned on N accepting to within precision ǫ/nt/2 .
5. While kPet kF ≥ ǫΩ(log(1/ǫ))t/2 for some t,

′

• Let t′ be the least t such that kPet′ kF ≥ ǫΩ(log(1/ǫ))t /2 .
√
P
• Let A = Pet′ /kPet′ kF . Let hA (x) =
i1 ,...it′ Ai Hec(i) (x)/ t! For each positive integer T ,
approximate
Pr [|hA (X)| ≥ T + 1]
e
X∼G

until one is found that is at least
′

′

3 exp(2 − Ω(T )2/t )) + ǫ/Cn2t

for a sufficiently large constant C. Let F be the filter that accepts when |hA (x)| ≤ T + 1.
• Recalculate Pet , for all 1 ≤ t ≤ k, where all expectations are conditioned N and the filters F
from all previous iterations.

6. End while.

7. For 1 ≤ t ≤ k, compute the SVD of M (Pet ), Pet considered as an n × nt−1 matrix, and let Vt ⊆ Rn
be the subspace spanned by all right singular vectors of Pt with singular value more than ǫ.

8. Let V be the span of V1 , . . . , Vk .

9. Let S ⊂ V be a set of unit vectors of size | dim(V )|O(dim(V )) such that for any unit vector v ∈ V ,
there is a v ′ ∈ S with kv − v ′ k2 ≤ 1/2.
48

√
10. For each v ∈ S, compute the median mv of v T X, for X ∼ G′ , to within ǫ/ pdim V using bisection
and statistical queries to approximate the Pr[v T X ≤ m] for m = µ′ + O(ǫ log(1/ǫ)). (We don’t
need to condition on any filters here).
11. Find a feasible point µ
eV of the LP µ
eV ∈ V with |v T µ
eV − v T mv | ≤ O(ǫ) for all v ∈ S.

12. Return µ
eV .

Using similar techniques to those used to express this algorithm in terms of Statistical Queries, we can
e
run Algorithm Filter-Gaussian-Unknown-Mean using poly(n/ǫ) time and calls to STAT(O(ǫ/poly(n))).
Note that we can approximate conditional expectations easily as a ratio of expectations approximated
by two SQ queries. Since, as we will show, our filters only throw away at most an O(ǫ) fraction of points,
we will not need to increase the precision beyond a constant factor to
√ do this.
O( log(1/ǫ))
. To show that we can use the
The algorithm needs approximate
expectations to within ǫ/n
√
O( log(1/ǫ))
) to obtain this, we need to note that√the distributions we approxioracle STAT(ǫ/(n ln(1/ǫ))
mate the expectations of are supported p
in an interval of length (n ln(1/ǫ))O( log(1/ǫ)) . Thanks to the naive
pruning of Step 3, only x with kxk2 ≤ 2n log(1/ǫ) contribute to these expectations. This suffices
√ to show
that the maximum value of all polynomials we consider on any such x is at most (n ln(1/ǫ))O(
We need to show the following for the filter step of our algorithm:

log(1/ǫ))

.

Proposition 8.12. The loop in Step 5 takes O(n2k ) iterations and all filters together accept with probability
at least 1 − O(ǫ).
8.3.1 Proof of Proposition 8.12 We now proceed with the proof. By standard concentration bounds, N
accepts with probability 1 − O(ǫ). Let F ′ be the event that N and all filters F from previous iterations
accept. We assume inductively that PrG′ [F ′ ] ≤ O(ǫ), and need to show that the same holds if we include
the filter F produced in the current iteration.
p
e = N (µ, I) be a Gaussian with dTV (G′ , G)
e ≤ ǫ such that kµk2 ≤ O(ǫ log(1/ǫ)). We write
Let G
G′ |F ′ for the distribution obtained by conditioning on F ′ . Since PrG′ [F ′ ] ≤ Cǫ, we have that
e ≤ dTV (G′ , G)
e + dTV (G′ , G′ |F ′ ) ≤ (C + 1)ǫ .
dTV (G′ |F, G)

e + dTV (G′ |F, G)E
e − dTV (G′ |F, G)L,
e
Thus, for the first iteration, we have G′ |F ′ = G
for distributions E
and L with disjoint supports.
e + wE E − wL L, for distributions E and L with disjoint
For any iteration, we will write G′ |F ′ = wGe G
supports, where wE + wL = O(ǫ) and wGe = 1 + O(ǫ). In the first iteration, we will take wGe = 1 and
e
wE = wL = dTV (G′ |F ′ , G).
We will need properties of the polynomials hA (x) for the analysis. In particular, we show the following:

Lemma 8.13. Given a rank-t symmetric tensor A over Rn , let hA (x) be as in Definition 8.11. Then, we
have:
(i) EX∼N (0,I) [hA (X)2 ] = kAk2F .
(ii) If
PB is a rank-t tensor with Bi =
i Ai Bi .

√

t!EX∼P [Hec(i) (X)], for a distribution P, then EX∼P [hA (X)] =

(iii) We can recover A from hA (x) using

√
t!Ai1 ,...,it =

∂
∂xi1

· · · ∂x∂i hA (x).
t

(iv) If O is an orthogonal matrix, then hA (Ox) = hB (x) for a symmetric rank-t tensor B with kBkF =
kAkF .
49

√
(v) If B is a rank-t tensor with Bi = t!EX∼P [Hec(i) (X)], for a distribution P, and j > 0, v ∈ Rn ,
√
then EX∼P [Hej (v · X)]/ t! = B(v, . . . , v).
Proof. For (i), we first need
Pto get an expression for√the coefficients of each Hea (x), since they appear
multiple times in hA (x) = i∈{0,...,n}t Ai Hec(a) (x)/ t!. Let c−1 (a) be a function mapping a ∈ Zt≥0 with
−1
kak1 = t to i ∈ {1, . . . , n}t , with c(c−1 (a)) = a for all a. Since A is
 symmetric, the choice of c does
t
not affect Ac−1 (a) for any a. Note that, for a given a, there are a1 ,...,an = t!/n(a) possible i with c(i) = a.
Thus, we have:
X
√
hA (x) =
Ai Hec(a) (x)/ t!
i∈{0,...,n}t

X

=

kak1 =t

=

√
t!/n(a) · Ac−1 (a) Hea (x)/ t!

X √
t!/n(a) · Ac−1 (a) Hea (x) .

kak1 =t

Now, by orthogonality of Hea (x) with distinct a, we have that, for X ∼ N (0, I), it holds:
X
E[hA (X)2 ] =
t!/n(a)2 · A2c−1 (a) E[Hea (X)2 ]
kak1 =t

=

X

kak1 =t

=

X
i

t!/n(a) · A2c−1 (a)

A2i = kAk2F .

For (ii), we now have that
√ X
EX∼P [hA (X)] = 1/ t! ·
Ai EX∼P [Hec(i) (x)]
=

X

i

Ai Bi .

i

Q
For (iii), note that a with kak1 = t, Hea (x) has only one monomial of degree t, which is xai i . Thus,
given i ∈ {1, . . . , n}t , there is only one a with kak1 = t and ∂x∂i · · · ∂x∂i Hea (x) 6= 0, which is a = c(i)
t
1
and has
∂
∂
∂
∂ Y c(i)j
Hec(i) (x) =
xj
···
···
∂xi1
∂xit
∂xi1
∂xit
j
Y
=
c(i)j ! = n(c(i)) .
j

Thus, we have
√
∂
∂
hA (x) = n(c(i)) · t!/n(a) · Ac−1 (a)
···
∂xi1
∂xit
√
= t!Ai .
For (iv), consider the linear transformation of rank-t tensors O⊗t , which for our purposes can be defined as
the unique function such that (O ⊗t A)(v1 , . . . , vt ) = A(OT v1 , . . . , OT vt ), for all rank-t tensors over Rn , A,
50

and vectors v1 . . . vt ∈ Rn . It is a fact that kO⊗t AkF = kAkF . Now,
P we∂gconsider the t-th order directional
derivatives of a function f , ∇v1 , · · · ∇vt f (x), where ∇v g(x) = i vi ∂xi (x). These can be expressed in
terms of the rank-t tensor F , with Fi = ∂x∂i · · · ∂x∂i f (x), as ∇v1 , . . . ∇vt f (x) = F (v1 , . . . vt ). We have
t
1
that, for all i,
∂
∂
hA (Ox)
···
∂xi1
∂xit
= ∇ei1 · · · ∇eit hA (Ox)

= ∇OT ei · · · ∇OT eit hA (x)
1
√
= t!A(OT ei1 , . . . , O T eit )
√
= t!(O⊗t A)i
∂
∂
h ⊗t (a) .
...
=
∂xi1
∂xit O A

Since hA (Ox) and hO⊗t A (x) are both multivariate polynomials of degree t, that these derivatives agree
means that the coefficients of all monomials of degree-t agree. We thus have hA (Ox) = hO⊗t A (x) +
p(x), where p is a polynomial of degree at most t − 1. Since hO⊗t A (x) is a linear combination of Hermite polynomials of degree t, which are orthogonal to all polynomials of degree smaller than t, we have
EX∼N (0,I) [hO⊗t A (X)p(X)] = 0, and so
1 = kAkF = EX∼N (0,I) [hA (Ox)2 ]

= EX∼N (0,I) [hO⊗t A (X)2 ] + EX∼N (0,I) [p(X)2 ] + 2EX∼N (0,I) [hO⊗t A (X)p(X)]

= kO ⊗t AkF + EX∼N (0,I) [p(X)2 ] + 0
= 1 + EX∼N (0,I) [p(X)2 ] .

Since EX∼N (0,I) [p(X)2 ] = 0, we must have p(x) ≡ 0, and thus
hA (Ox) = hO⊗t A (x) .
Taking B = O⊗t A gives (iv).
For (v), let O be an orthogonal matrix that gives a rotation mapping e1 to v, a = {t, 0, . . . , 0} and T1 be
the rank-t tensor with (1, . . . , 1) entry 1 and every other entry 0. Then, we can rewrite
√
√
Het (v · x) = Het ((Ox)1 ) = Hea ((Ox)1 ) = t!hT1 (Ox) = t!hO⊗t T1 (x) .
For any i, we have
(O⊗t T1 )i = (O ⊗t T1 )(ei1 , . . . , eit )
= T1 (OT ei1 , . . . , O T eit )
=

t
Y

(OT eij )1

j=1

=

t
Y

(Oe1 )ij

j=1

=

t
Y

vi j .

j=1

51

Thus we have O⊗t T1 = v ⊗t , the tensor with entries

Qt

j=1 vij ,

and so

√
EX∼P [Het (v · X)]/ t! = EX∼P [hv⊗t (X)]
=

X

Bi

i

t
Y

vij = B(v, . . . , v) .

(by (ii))

j=1

This completes the proof.
√
We write Pt′ or Gt′ for the rank-t′ tensor with entries t!E[hc(i) (X)], where X is distributed according
e respectively. We know that kPet′ kF ≥ ǫΩ(log(1/ǫ)t′ /2 ).
to G′ |F ′ or G
′

Lemma 8.14. When kPet′ kF ≥ ǫΩ(log(1/ǫ)t /2 ), we have |EX∼G′ [hA (X)]| ≥ ǫ · Ω(log(1/ǫ))t/2 .

Proof. The assumption on the SQ errors imply that the corresponding entries of Pt′ and Pet′ are within
′
ǫ/nt /2 . It follows that
√
kPt′ − Pet′ kF ≤ t′ !ǫ .

Using Lemma 8.13 (iii), we have

EX∼G′ [hA (X)] =

X

Ai (Pt′ )i

i∈[n]t′

= 1/kPet kF ·

X

i∈[n]

t′



(Pet′ )i (Pt′ )i

= 1/kPet kF · kPet′ k2F +

X

i∈[n]

t′



(Pet′ )i (Pt′ − Pet′ )i 

≥ kPet′ kF − kPt′ − Pet′ kF ≥ kPet′ kF −
p
t′
′
≥ ǫΩ(log(1/ǫ)t /2 )) − ǫ log(1/ǫ)

√

t′ !ǫ

′

≥ ǫΩ(log(1/ǫ)t /2 ) .

Lemma 8.15. EX∼Ge [hA (X)] = O(ǫ

p

log(1/ǫ))t and EX∼Ge [hA (X)2 ] = O(1).

Proof. Since for all a ∈ Zk>0 , EX∼N (0,I) [Hea (X)] = 0, we have EX∼N (0,I) [hA (X] = 0. By Lemma 8.13
(i), EX∼N (0,I) [hA (X)2 ] = kAk2F = 1.
e = N (µ, I) instead of N (0, I). Consider a rotation given
We need to take these expectations under G
by an orthogonal matrix O that maps kµk2 e1 to µ. By Lemma 8.13 (i), there is a symmetric rank-t tensor
B with kBkF = 1 such that hA (OX) = hB (X). Now we have that
EX∼G′ [hA (X)] = EX∼N (kµk2 e1 ,I) [hA (OX)] = EX∼N (kµk2 e1 ,I) [hB (X)] ,
and similarly for hA (X)2 .
n−1
for the vector dropping the first coordinate of a, we have
For a ∈ Zn≥0 with a 6= 0, writing a−1 ∈ Z≥0
EX∼N (kµk2 e1 ,I) [Hea (X)] = EX∼N (0,I) [Hea−1 (X)]EX∼N (kµk2 ,1) [Hea1 (X)]
= δa−1 ,0 EX∼N (0,1) [Hea1 (X + kµk2 )] .
52

Note that there is only one index i such that c(i) is zero in all except the first coordinate, and so we have
√
EX∼Ge [hA (X)] = EX∼N (kµk2 e1 ,I) [hB (X)] = B1,...,1 EX∼N (0,1) [Het (X + kµk2 )/ t!] .
Since kBkF = 1, we have

√
EX∼Ge [hB (X)] ≤ |EX∼N (0,1) [Het (X + kµk2 )/ t!]|.

i
By standard results, we have that dHe
dx (x) = iHei−1 (x), and so by Taylor’s theorem we have Hei (x +

Pi
j
kµk2 ) = j=0 ji kµk2 Hei−j (x). Thus,

|EX∼N (0,1) [Het (X + kµk2 )]| =
=
This gives

t  
X
t

i=0
kµkt2

i

kµki2 |EX∼N (0,1) [Het−i (X)]

.

p
√
√
EX∼Ge [hA (X)] ≤ kµkt2 / t! = O(ǫ log(1/ǫ))t / t! ,

as required.
Similarly, we have, for all a, b ∈ Zn>0 ,

EX∼N (kµk2 e1 ,I) [Hea (X)Heb (X)] = EX∼N (0,I) [Hea−1 (X)Heb−1 (X)]EX∼N (kµk2 ,1) [Hea1 (X)Heb1 (X)]
= δa−1 b−1 · EX∼N (0,1) [Hea1 (X + kµk2 )Heb1 (X + kµk2 )] .
When kak1 = kbk1 = t, if a−1 = b−1 , then a = b (since a1 = b1 = t − kak1 ). For 1 ≤ j ≤ t, we have:

!2 
j  
X
j
EX∼N (0,1) [Hej (X + kµk2 )2 ] = EX∼N (0,1) 
kµki2 Hej−i (x) 
i
i=0

=

=

=

≤

j  2
X
j
i=0
j 
X

i=0
j
X

i=0
j
X

j
i

2

2
kµk2i
2 |EX∼N (0,1) [Hej−i (X) ]

kµk2i
2 i!

2
kµk2i
2 (j!/(j − i)!) /i!
2
kµk2i
2 j i/i!

i=0
j
X

≤2

i

2−2i

i=0

(since kµk2 ≤ 1/2k ≤ 1/2j)

≤3.

Putting these together, for a, b with kak1 = kbk1 = t, we have
EX∼N (kµk2 e1 ,I) [Hea (X)Heb (X)] ≤ 3δab .
The sum of squares of coefficients of all Hea (x) in hB (X) is EX∼N (0,1) [hB (X)2 ] = kBkF = 1, and so we
have that EX∼N (kµk2 e1 ,I) [hB (X)2 ]| ≤ 3. Finally, recall that EX∼G [h(A)2 ] = EX∼N (kµk2 e1 ,I) [hB (X)2 ]|,
and so this is O(1), as required.
53

Now consider the equation
EX∼G′ [hA (X)] = wGe EX∼Ge [hA (X)] + wE EX∼E [hA (X)] − wL EX∼L [hA (X)].
We know that the LHS is Ω(ǫ log(1/ǫ)t/2 ) and that the first term on the RHS is smaller. Therefore, one
e we will use standard concentration inequalities to show
of the last two terms is small. Since wL L ≤ G,
t/2
that wL EX∼L [hA (X)] is O(ǫ log(1/ǫ) ). If we cannot find a filter, then wE E ≤ G′ must satisfy similar
concentration inequalities, which would imply that wE EX∼E [hA (X)] is O(ǫ log(1/ǫ)t/2 ). Since some term
on the RHS must be bigger than this, we can find a filter.
Lemma 8.16. For X ∼ N (0, I), if p(x) is a degree-d polynomial with E[p(X)2 ] ≤ 1, we have that
Pr [|p(X)| ≥ T + E[p(X)]] ≤ exp(2 − Ω(T )2/d )) .
′

Lemma 8.17. We have that wL |EX∼L [hA (X)]| ≤ ǫ · O(log(1/ǫ))t /2 ).
Proof. We start with the following claim:
Claim 8.18. For any R > 0, d ∈ Z+ , ǫ > 0, and exp(−(a/R)2/d ) = ǫ, we have
Z ∞
exp(−(T /R)2/d )T dT ≤ (d2 /2)ǫ(d + ln(1/ǫ))d−1 .
a

Proof. Note that (a/R)2/d = ln(1/ǫ). First, we change variables to x = (T /R)2/d to obtain
Z ∞
Z ∞
dT
2/d
exp(−x)xd/2
exp(−(T /R) )T dT =
dx
dx
ln(1/ǫ)
ln(1/ǫ)
Z ∞
exp(−x)xd/2 · (Rd/2)xd/2−1 dx
=
ln(1/ǫ)
Z ∞
exp(−x)xd−1 dx .
= (Rd/2)
ln(1/ǫ)

We can now integrate by parts
Z
Z ∞
d−1
d−1
exp(−x)x dx = ǫ ln(1/ǫ)
+ (d − 1)
ln(1/ǫ)

∞

exp(−x)xd−2 dx .

ln(1/ǫ)

By a simple induction, we have
Z

∞

exp(−x)xd−1 dx = ǫ

ln(1/ǫ)

d−1
X
j=0

d!/(d − j)! ln(1/ǫ)(d−j−1) .

Now we have
Z

∞

2/d

exp(−(T /R)

)T dt = (d/2)ǫ

a

d−1
X
j=0

d!/(d − j)! ln(1/ǫ)(d−j−1)
2/d

≤ (d/2) exp(−(a/R)
2

)

d−1
X
j=0

d−1

≤ (d /2)ǫ(d + ln(1/ǫ))

54

.

dj ln(1/ǫ)(d−j−1)

Let µh = EX∼Ge [hA (X)]. We have the following sequence of inequalities:
Z ∞
T wL Pr [|hA (X)| > T ]dT
wL EX∼L| [hA (X)2 ] =
X∼L
0
Z ∞
T min{wL , Pr [|hA (X)| > T ]}dT
≤
e
X∼G
0
Z ∞
′
T min{wL , exp(2 − Ω(T − µh )2/t ))}dT
≤
Z0 ∞
′
T min{wL , exp(2 − ((T − µh )/R)2/t ))}dT
(for some R > 0)
≤
0
Z ∞
′
(T + µh ) min{wL , exp(2 − (T /R)2/t ))}dT
=
−µ
Z ∞
Z ch
′
exp(2 − (T /R)2/t )(T + µh )dT
(T + µh )wL dT +
=
c

−µh

′

≤ wL (c + µh )2 /2 + e2 (c + µh )/c ·
t′

′2

′

Z

(where c = (R ln 1/wL )t /2 )

∞

′

exp(−(T /R)2/t )T dT

c
′

≤ wL O(ln(1/wL )) + O(t wL (t + ln(1/wL ))t −1 )
′

≤ ǫ · O(ln(1/ǫ))t ,

p
where the last line follows from t′ ≤ k = O( log(1/ǫ)) ≤ O(log(1/ǫ)). Then, by an application of the
Cauchy-Schwarz inequality, we have that
q
t′ /2
2E
2
).
wL |EX∼L [hA (X)]| ≤ wL EX∼L [|hA (X)|] ≤ wL
X∼L [hA (X)] ≤ ǫ · O(log(1/ǫ))
This completes the proof of the lemma.

′

Lemma 8.19. If PrX∼G′ |F ′ [|hA (X)| ≥ T + 1] ≤ O(exp(2 − Ω(T )2/d )) + ǫ/(2n)2t ), for all integers T ,
′
then wE |EX∼E [hA (X)]| ≤ O(ǫ ln(1/ǫ)t /2 ).
′ includes the filter M , we have that the support of G′ |F ′ and the support of E includes only
Proof. Since F p
x with kxk2 ≤ 2n ln(1/ǫ):
p
p
Claim 8.20. When kxk2 ≤ 2n ln(1/ǫ), then |hA (x)| ≤ (2n ln(1/ǫ))t .
p
Proof. Note that t′ ≤ k ≤ 2n log(1/ǫ). Using the explicit formula for the coefficient Hei (x), we can

55

show that for |x| ≤

p

2n ln(1/ǫ), with k ≥ i, the Hei (x) is dominated by its leading coefficient:
|Hei (x)| =
≤
≤
≤

⌊i/2⌋

X
j=0

⌊i/2⌋

i!(−1)j xi−2j /j!(i − 2j)!2j

X

j!(3/2)j |x|i−2j

X

p
j!(3/2)j ( 2n ln(1/ǫ))i−2j

j=0

⌊i/2⌋
j=0

⌊i/2⌋

X

(2n log(1/ǫ))(i−j)/2

j=0

≤ 2 · (n log(1/ǫ))i/2 .
Therefore, for any a ∈ Zn≥0 with

P

i

ai = t, and kxk2 ≤
|ha (x)| =
≤

′

n
Y

i=1
n
Y
i=1

p

2n log(1/ǫ), we have that

p
|Heai (x)/ ai !|

p
ai
2 n log(1/ǫ)

p
≤ (2 n log(1/ǫ))t .

t
the L1 -norm of the entries is at most nt/2 . Thus, we have that
Since kAkF = 1 and
p A has n tentries, p
t/2
|hA (x)| ≤ n · (2 n log(1/ǫ)) = (2n log(1/ǫ))t .

We note that since

′

Pr

[|hA (X)| ≥ T + 1] ≤ O(exp(2 − Ω(T )2/d )) + ǫ/(2n)2t ) ,

Pr

[|hA (X)| ≥ T + 2] ≤ O(exp(2 − Ω(T )2/d )) + ǫ/(2n)2t ) ,

X∼G′ |F ′

for integers T , that
X∼G′ |F ′

′

for all T .

56

Similarly to the proof of Lemma 8.17, we obtain:
Z ∞
T wE Pr [|hA (X)| > T ]dT
wE EX∼E [hA (X)2 ] =
X∼E
0
√
Z
t
(2n

log(1/ǫ))

=

≤
≤
≤
≤

T wE Pr [|hA (X)| > T ]dT

Z

Z

Z

Z

0

(2n

√

X∼E

log(1/ǫ))t

T min{wE , Pr ′ [|hA (X)| > T ]}dT

0
(2n

√

(2n

√

0

−2

(2n

√

X∼G

log(1/ǫ))t

′

T min{wE , O(exp(2 − Ω(T − 2)2/t ) + ǫ/Cn2t )}dT

log(1/ǫ))t −1

′

(T + 2) min{wE , O(exp(2 − (T /R)2/t ) + ǫ/Cn2t )}dT

log(1/ǫ))t

2t′

O(T ǫ/(2n) )dT +
0

Z

c

(T + 2)wE dT
−2

Z ∞
O(exp(−(T /R)2/t ))T dT
(where c = (R ln 1/wE )t/2 )
+ (c + 1)/c ·
c
p
′
′
= (2n log(1/ǫ))2t · ǫ/(2n)2t + O(wE c2 /2) + O(t′2 wE (t′ + ln(1/wE ))t −1 )
≤ ǫ · (4 log(1/ǫ)t + wE · O(ln 1/wE )t

≤ ǫ · O(log(1/ǫ))t .

Then, by the Cauchy-Schwarz inequality, we conclude that
q
t/2
2E
2
.
wE |EX∼E [hA (X)]| ≤ wE EX∼E [|hA (X)|] ≤ wE
X∼E [hA (X)] ≤ ǫ · O(log(1/ǫ))

This completes the proof.

As an immediate consequence, we obtain:
Corollary 8.21. There is an integer 0 ≤ T ≤ O(n
Pr

X∼G′ |F ′

p

log(1/ǫ))t such that
′

[|hA (X)| ≥ T + 1] ≥ 3 exp(2 − Ω(T )2/d )) + 2ǫ/Cn2t .

We can now prove the following crucial lemma:
Lemma 8.22. The algorithm finds a T with
Pr

X∼G′ |F ′

′

[|hA (X)| ≥ T + 1] ≤ 3 exp(2 − Ω(T )2/d )) + ǫ/Cn2t .

Proof.
8.21, such a T exists, and therefore our algorithm will find one after enumerating
p By Corollary
t
O(n log(1/ǫ)) possibilities.

Let F be the event that the new filter accepts. In the next iterations, we will use G′ |F ′ ∩ F instead of
G′ |F ′ . We need to show that the parameters wGe , wE and wL improve in such a way that we only need a
bounded number of iterations:
e + w′ E ′ − w′ L′ , where L′ and E ′ have disjoint supports
Claim 8.23. We can write G′ |F ′ ∩ F = w′e G
E
L
G

′

′ , w ′ > 0 and w ′ + w ′ ≤ w + w − ǫ/Cn2t . The probability that the filter rejects is at most
wE
E
L
L
E
L
′ − w ′ ).
O(wE + wL − wE
L

57

Proof. This proof is very similar to that of Claim 26 from [DKS16c]. Let ¬F be the event that the filter
′
rejects, i.e., that |hA (X)| ≥ T + 1. We have that PrG′ |F ′ [¬F ] ≥ 3 exp(2 − Ω(T )2/d )) + ǫ/Cn2t . On the
other hand, by the concentration inequality, PrGe [¬F ] ≤ exp(2 − Ω(T )2/d )). Thus, we have that
′

Pr [¬F ] ≥ 3 Pr[¬F ] + ǫ/Cn2t .

G′ |F ′

e
G

e E and L yields for the event ¬F that
However, the defining relation between G′ |F ′ and G,
Pr [¬F ] ≥ wGe Pr[¬F ] + wE Pr[¬F ] − wE Pr[¬F ] .

G′ |F ′

E

e
G

E

Since
Pr [¬F ] ≤ wGe Pr′ [¬F ] + wE Pr[¬F ] ,

G′ |F ′

and wGe ≤ 1 + O(ǫ), we must have

E

G

Pr [¬F ] ≤ (2 + O(ǫ))wE Pr[¬F ] ,

G′ |F ′

E

and
Pr′ [¬F ] ≤ (1/3 + O(ǫ))wE Pr[¬F ].
E

G

Then, we get that




′
′
1 − Pr [¬F ] (G |F )(x) = 1 − Pr [¬F ] (G′ |F ′ ∩ ¬F )(x)
G′ |F ′
G′ |F ′




e
e
= wGe G(x) + wE 1 − Pr[¬F ] E(x) + wL 1 − Pr[¬F ] L(x) − wGe Pr′ [¬F ]G(x)
.
E

L

G

Thus, we have

wL (1 − PrL [¬F ]) − wGe PrG′ [¬F ]
1 − PrPe [¬F ]


≤ wL + (1 + O(ǫ)) Pr′ [¬F ] + O ǫ Pr[¬F ]
G
Pe


≤ wL + (1/3 + O(ǫ)) wE Pr[¬F ] + O ǫ Pr [¬F ] .

′
wL
=

E

Also we have

G′ |F ′

wE (1 − PrE [¬F ])
1 − PrG′ |F ′ [¬F ]




≤ wE 1 − Pr[¬F ] + O ǫ Pr [¬F ] .

′
wE
=

E

Thus,
wL + wE −

′
wL

−

′
wE

G′ |F ′


≥ (2/3 − O(ǫ)) wE Pr[¬F ] − O ǫ Pr[¬F ]
E



≥ (1/3 − O(ǫ)) Pr [¬F ] ≥ ǫ/Cn

Pe
2t′

G′ |F ′

Note that the penultimate inequality also gives that
′
′
Pr [¬F ] ≤ (3 + O(ǫ)) (wL + wE − wL
− wE
).

G′ |F ′

This completes the proof.
Proposition 8.12 now follows using induction on the iterations.
58

.

8.3.2

Completing the Proof of Correctness

Lemma 8.24. dim(V ) ≤ O(log(1/ǫ))k .
Proof. After leaving the filter loop, for all 1 ≤ t ≤ k, we have that kPet kF ≤ ǫO(log(1/ǫ))t/2 . M (Pet ) has
the same Frobenius norm, and thus the L2 -norm of its singular values, when considered as a matrix. Thus,
t/2 singular values bigger than ǫ. So, we have that dim(V ) = O(log(1/ǫ))t/2 ,
there are at most O(log(1/ǫ))
t
Pk
and so dim(V ) ≤ t=1 dim Vt ≤ O(log(1/ǫ))k .
Let µV be the projection of µ onto the subspace V . Now we can show using our moment matching
lemma that it suffices to approximate µV .

Lemma 8.25. We have that kµV − µk2 ≤ O(ǫ).
Proof. Let v = (µV − µ)/kµV − µk2 . Note that v is a unit vector perpendicular to V and we need to show
that v T µ ≤ O(ǫ). We will apply Lemma 8.1 to G′′ , the projection G′ conditioned on the event that all the
e G′ |F ′ ) ≤ O(ǫ), and so we have that
filters we produced accept F ′ , onto v. Note that G′ |F ′ has dTV (G,
′′
T
dTV (G , N (0, v µ)) ≤ O(ǫ). We can bound the expectation of the Hermite polynomials as follows, for
1 ≤ t ≤ k:
√
√
|EX∼G′′ [Het (X)/ t!]| = |EX∼G′ [Het (v · X)/ t!]| = |P t (v, . . . , v)|
= |(v ⊗t−1 )T M (P t )v| ≤ kv ⊗t−1 k2 kM (P t )vk2

≤ 1 · O(ǫ) .

√
√
On the other hand, we have EX∼N (0,1) [Het (X)/ t!] = 0, for 1 ≤ t ≤ k. For t = 0, Het (X)/ t! = 1,
which has expectation 1 under both G′′ and N (0, 1). We want to consider the difference in the expectations
of X t , for√1 ≤ t ≤ k. We can write xt as a linear combination of Hermite polynomials, xt =
Pt
we have that EX∼N (0,1) [(X t )2 ] =
Pi=02 ai Hei (x)/ i. Using the orthonormality of these polynomials,
t 2
t
i ai . On the other hand, by standard results, EX∼N (0,1) [(X ) ] = 2 t!. Thus, we have:
|EX∼G′′ [X t ] − EX∼N (0,1) [X t ]| =

t
X
i=0


√
√ 
ai EX∼G′′ [Hei (X)/ i!] − EX∼N (0,1) [Hei (X)/ i!]

≤ O(ǫ) ·
≤ O(ǫ) ·

t
X

ai

i=0

√

t·

√

2t t!

√
√ √
Note that for t ≥ 10, t2t t! ≤ (t − 1)!/t. Thus, there is a constant c > 0 such that this O(ǫ) · t · 2t t!
is smaller than (t − 1)!ct ǫ/t, for all 1 ≤ t ≤ k. Now we can apply Lemma 8.1 with
p δ = cǫ and obtain
T
that |v µ| ≤ O(δ) = O(ǫ). We need to set k to be a sufficiently high multiple of ǫ ln(1/ǫ) to make this
work.
It remains to analyze the rest of the algorithm and show that µ
eV it produces is close to µV .

Lemma 8.26. We can construct a set S ⊂ V of unit vectors of size dim(V )O(dim(V )) such that for any unit
vector v ∈ V , there is a v ′ ∈ S with kv − v ′ k2 ≤ 1/2, in time dim(V )O(dim(V )) .

59

Proof. Let ℓ = dim(V ). We will construct such a cover for Rℓ and translate that to V by using the orthonormal basis for V given by the right singular vectors of M (Pt√
) with singular values bigger than ǫ. We
can divide the cube [−1, 1]ℓ into ℓO(ℓ) cubes of side length 1/(2 ℓ). For each cube, we check if it has
a corner with L2 -norm ≥ 1 and a corner with L2 -norm ≤ 1. If it does not, it does not contain any unit
vectors so we can ignore it. If the cube does contain any unit vectors, then if its center is v, we add the
normalized vector
v/kvk2 to S. Since there is a unit vector v ′ in the cube, and all vectors in the cube have
√
kv ′ − vk2 ≤ ℓkv ′ − vk∞ ≤ 1/4, we have |kvk2 − 1| ≤ 1/4, and so kv − v/kvk2 k2 ≤ 1/4. Thus, for any
unit vector v ′ in this cube, we have kv ′ − v/kvk2 k2 ≤ 1/4 + 1/4 ≤ 1/2. Since every unit v ′ is in some cube
whose normalized center we added to S, we are done.
Firstly, we note that in order to approximate v T µ, it is sufficient to find an x with PrX∼G′ [v T X ≥ x] =
1/2 + O(ǫ):
Lemma 8.27. For all x ∈ R with | PrX∼G′ [v T X ≥ x] − 1/2| ≤ 3ǫ, we have that |v T µ − x| ≤ O(ǫ) .
e projected onto v has G(x − v T µ) ≥ 1/2 for all x with |x − v T µ| ≤ O(ǫ).
Proof. First note that the pdf of G
T
Supposing that |x − v µ| ≥ 8ǫ, we have that | PrX∼Ge [v T X ≥ x] − 1/2| ≥ 4ǫ, and so | PrX∼G′ [v T X ≥
e ≥ 3ǫ.
x] − 1/2| ≥ 4ǫ − dTV (G′ , G)

To show that we can find such a point by bisection, we need to show that there is an interval of such
points of reasonable length where we are looking for them:
Lemma 8.28. Given a unit vector v ∈ Rn , there is an interval [a, b] such that
• for all x ∈ [a, b], we have that | PrX∼G′ [v T X ≥ x] − 1/2| ≤ 2ǫ,
• b − a = Θ(ǫ),
• and |a|, |b| ≤ O(ǫ

p

log 1/ǫ).

Proof. We can take [a, b] to be the set of x with | PrX∼Ge [v T X ≥ x] − 1/2| ≤ ǫ. This is an interval since
e ≤ 2ǫ.
PrX∼Ge [v T X ≥ x] is monotone. All x in it have | PrX∼G′ [v T X ≥ x] − 1/2| ≤ ǫ + dTV (G′ , G)
T
T
T
Thus,
pby the previous lemma, |b − v µ|, |a
p− v µ| ≤ O(ǫ) and thus |b − a| ≤ O(ǫ), and since |v µ| ≤
O(ǫ log 1/ǫ), we have that |a|, |b| ≤ O(ǫ log 1/ǫ).
Thus, we obtain:

Lemma 8.29. Given a unit vector v ∈ Rn , we can find an mv with |v T µ−mv | ≤ O(ǫ) using O(log log 1/ǫ)
statistical queries of precision ǫ/2.
Proof. We use bisection to find a point where our SQ approximation pe to PrX∼G′ [v T X ≥ x] is within 5ǫ/2
of 1/2. If we find such a point, it has |v T µ − mv | ≤ O(ǫ), by Lemma 8.27. Lemma
p8.28 yields that there
is an interval [a, b] of length O(ǫ) containing such points in the interval |x| ≤ O(ǫ log(1/ǫ)). Indeed, if
our test point x has pe > 1/2 + 5ǫ/2, then x > b and if pe < 1/2 − 5ǫ/2, then x < a. Thus, [a, b] remains a
subinterval of the interval we are considering.
We now have that µv is a feasible point of the LP considered in Step 11. The following lemma completes
the proof:
Lemma 8.30. Any feasible point of the LP considered in Step 11, µ
eV has kµV − µ
eV k2 ≤ O(ǫ).
60

Proof. Consider the vector v = (µV − µ
eV )/kµV − µ
eV k2 . Note that v is in V , since µC , µ
eV are. Since
′
′
v is a unit vector in V , there is a v ∈ S with kv − v k2 ≤ 1/2. Since µ
eV is a solution to the LP,
v ′T (µV − µ
eV ) ≤ O(ǫ). Thus, we have that
kµV − µ
eV k2 = v T (µV − µ
eV )

= v ′T (µV − µ
eV ) + (v − v ′ )T (µV − µ
eV )

≤ O(ǫ) + kµV − µ
eV k2 /2 .

Therefore, kµV − µ
eV k2 ≤ O(ǫ), as required.

Proof of Theorem 8.7. Since the LP has a feasible point, we can find such a point µ
eV that has kµV − µ
eV k2 ≤
O(ǫ). By the previous lemma, we have that kµV −µk2 ≤ O(ǫ). Thus, the algorithm is correct. All statistical
queries are of the claimed precision. We need to get bounds on the running time and number of statistical
queries.
Step 2 that uses the algorithm from [DKK+ 16], takes poly(n/ǫ) time and statistical queries. Finding Pet ,
for each 1 ≤ t ≤ k, takes nt statistical queries, giving nO (k) time total. There are at most O(n2k ) iterations
of the loop. Each iteration takes poly(n/ǫ) time and statistical queries to find T , and nt statistical queries
to recompute Pet . It suffices to compute the SVD to within Frobenius norm 1/poly(n/ǫ), which takes time
O(k)
k
. Computing it takes
poly(n/ǫ). The set S has size dim(V )O(dim(V )) = log(1/ǫ)kO(log(1/ǫ)) = 2log(1/ǫ)
O(k)
O(k)
log(1/ǫ)
log(1/ǫ)
statistical queries and time. The LP has
. Approximating the medians takes 2
time 2
O(k)
O(k)
O(k)
log(1/ǫ)
bits, and so with a
constraints and log(1/ǫ)
variables. The size of the LP is 2log(1/ǫ)
2
O(k)
time.
polynomial time LP solver, we can get 2log(1/ǫ)
We thus have that the total time and statistical queries are both at most
√
√
log(1/ǫ)O( log(1/ǫ))
O(k)
log(1/ǫ)O(k)
O( log(1/ǫ))
.
+2
n
poly(1/ǫ) + 2
=n

References
[ABG+ 14]

J. Anderson, M. Belkin, N. Goyal, L. Rademacher, and J. R. Voss. The more, the merrier: the
blessing of dimensionality for learning large gaussian mixtures. In Proceedings of The 27th
Conference on Learning Theory, COLT 2014, pages 1135–1164, 2014.

[ADLS17]

J. Acharya, I. Diakonikolas, J. Li, and L. Schmidt. Sample-optimal density estimation in
nearly-linear time. In Proceedings of the Twenty-Eighth Annual ACM-SIAM Symposium on
Discrete Algorithms, SODA 2017, pages 1278–1289, 2017.

[AK01]

S. Arora and R. Kannan. Learning mixtures of arbitrary Gaussians. In Proceedings of the 33rd
Symposium on Theory of Computing, pages 247–257, 2001.

[AM05]

D. Achlioptas and F. McSherry. On spectral learning of mixtures of distributions. In Proceedings of the Eighteenth Annual Conference on Learning Theory (COLT), pages 458–469,
2005.

[AS72]

M. Abramowitz and I. Stegun. Handbook of Mathematical Functions. Dover, 1972.

[BBBB72]

R.E. Barlow, D.J. Bartholomew, J.M. Bremner, and H.D. Brunk. Statistical Inference under
Order Restrictions. Wiley, New York, 1972.

61

[BCMV14] A. Bhaskara, M. Charikar, A. Moitra, and A. Vijayaraghavan. Smoothed analysis of tensor
decompositions. In Symposium on Theory of Computing, STOC 2014, pages 594–603, 2014.
[BFJ+ 94]

A. Blum, M. Furst, J. Jackson, M. Kearns, Y. Mansour, and S. Rudich. Weakly learning DNF
and characterizing statistical query learning using Fourier analysis. In Proceedings of the
Twenty-Sixth Annual Symposium on Theory of Computing, pages 253–262, 1994.

[BFR+ 00]

T. Batu, L. Fortnow, R. Rubinfeld, W. D. Smith, and P. White. Testing that distributions are
close. In IEEE Symposium on Foundations of Computer Science, pages 259–269, 2000.

[BR13a]

Q. Berthet and P. Rigollet. Complexity theoretic lower bounds for sparse principal component
detection. In COLT 2013 - The 26th Annual Conference on Learning Theory, pages 1046–
1066, 2013.

[BR13b]

Q. Berthet and P. Rigollet. Optimal detection of sparse principal components in high dimension. Ann. Statist., 41(4):1780–1815, 2013.

[BS10]

M. Belkin and K. Sinha. Polynomial learning of distribution families. In FOCS, pages 103–
112, 2010.

[BV08]

S. C. Brubaker and S. Vempala. Isotropic PCA and Affine-Invariant Clustering. In Proc. 49th
IEEE Symposium on Foundations of Computer Science, pages 551–560, 2008.

[CDSS13]

S. Chan, I. Diakonikolas, R. Servedio, and X. Sun. Learning mixtures of structured distributions over discrete domains. In SODA, pages 1380–1394, 2013.

[CDSS14a] S. Chan, I. Diakonikolas, R. Servedio, and X. Sun. Efficient density estimation via piecewise
polynomial approximation. In STOC, pages 604–613, 2014.
[CDSS14b] S. Chan, I. Diakonikolas, R. Servedio, and X. Sun. Near-optimal density estimation in nearlinear time using variable-width histograms. In NIPS, pages 1844–1852, 2014.
[CFJ13]

T. Cai, J. Fan, and T. Jiang. Distributions of angles in random packing on spheres. Journal of
Machine Learning Research, 14(1):1837–1864, 2013.

[CGG02]

M. Cryan, L. Goldberg, and P. Goldberg. Evolutionary trees can be learned in polynomial time
in the two state general Markov model. SIAM Journal on Computing, 31(2):375–397, 2002.

[CGR15]

M. Chen, C. Gao, and Z. Ren. Robust covariance matrix estimation via matrix depth. CoRR,
abs/1506.00691, 2015.

[CKL+ 06]

C.-T. Chu, S. K. Kim, Y. A. Lin, Y. Yu, G. Bradski, A. Y. Ng, and K. Olukotun. Map-reduce
for machine learning on multicore. In Proceedings of the 19th International Conference on
Neural Information Processing Systems, NIPS’06, pages 281–288, Cambridge, MA, USA,
2006. MIT Press.

[CMW13]

T. T. Cai, Z. Ma, and Y. Wu. Sparse pca: Optimal rates and adaptive estimation. Ann. Statist.,
41(6):3074–3110, 12 2013.

[CMW15]

T. Cai, Z. Ma, and Y. Wu. Optimal estimation and rank detection for sparse spiked covariance
matrices. Probability Theory and Related Fields, 161(3):781–815, 2015.

[CQ10]

S. X. Chen and Y. L. Qin. A two-sample test for high-dimensional data with applications to
gene-set testing. Ann. Statist., 38(2):808–835, 04 2010.
62

[Dan16]

A. Daniely. Complexity theoretic limitations on learning halfspaces. In Proceedings of the
48th Annual Symposium on Theory of Computing, STOC 2016, pages 105–117, 2016.

[Das99]

S. Dasgupta. Learning mixtures of Gaussians. In Proceedings of the 40th Annual Symposium
on Foundations of Computer Science, pages 634–644, 1999.

[DBS17]

S. Du, S. Balakrishnan, and A. Singh. Computationally efficient robust estimation of sparse
functionals. CoRR, abs/1702.07709, 2017.

[DDKT16]

C. Daskalakis, A. De, G. Kamath, and C. Tzamos. A size-free CLT for poisson multinomials
and its applications. In Proceedings of STOC’16, 2016.

[DDO+ 13]

C. Daskalakis, I. Diakonikolas, R. O’Donnell, R.A. Servedio, and L. Tan. Learning Sums of
Independent Integer Random Variables. In FOCS, pages 217–226, 2013.

[DDS12a]

C. Daskalakis, I. Diakonikolas, and R.A. Servedio. Learning k-modal distributions via testing.
In SODA, pages 1371–1385, 2012.

[DDS12b]

C. Daskalakis, I. Diakonikolas, and R.A. Servedio. Learning Poisson Binomial Distributions.
In STOC, pages 709–728, 2012.

[DDS15]

A. De, I. Diakonikolas, and R. Servedio. Learning from satisfying assignments. In Proceedings of the 26th Annual ACM-SIAM Symposium on Discrete Algorithms, SODA 2015, pages
478–497, 2015.

[DG85]

L. Devroye and L. Györfi. Nonparametric Density Estimation: The L1 View. John Wiley &
Sons, 1985.

[DG92]

D. L. Donoho and M. Gasko. Breakdown properties of location estimates based on halfspace
depth and projected outlyingness. Ann. Statist., 20(4):1803–1827, 12 1992.

[Dia16]

I. Diakonikolas. Learning structured distributions. In P. Bühlmann, P. Drineas, M. Kane, and
M.J. van Der Laan, editors, Handbook of Big Data, Chapman & Hall/CRC Handbooks of
Modern Statistical Methods, chapter 15, pages 267–284. Taylor & Francis, 2016.

[DK14]

C. Daskalakis and G. Kamath. Faster and sample near-optimal algorithms for proper learning
mixtures of gaussians. In Proceedings of The 27th Conference on Learning Theory, COLT
2014, pages 1183–1213, 2014.

[DKK+ 16]

I. Diakonikolas, G. Kamath, D. M. Kane, J. Li, A. Moitra, and A. Stewart. Robust estimators in
high dimensions without the computational intractability. In Proceedings of FOCS’16, 2016.

[DKK+ 17a] I. Diakonikolas, G. Kamath, D. M. Kane, J. Li, A. Moitra, and A. Stewart. Being robust (in
high dimensions) can be practical. CoRR, abs/1703.00893, 2017.
[DKK+ 17b] I. Diakonikolas, G. Kamath, D. M. Kane, J. Li, A. Moitra, and A. Stewart. Robustly learning
a gaussian: Getting optimal error, efficiently. CoRR, abs/1704.03866, 2017.
[DKS16a]

I. Diakonikolas, D. M. Kane, and A. Stewart. The fourier transform of poisson multinomial
distributions and its algorithmic applications. In Proceedings of STOC’16, 2016.

[DKS16b]

I. Diakonikolas, D. M. Kane, and A. Stewart. Optimal learning via the fourier transform for
sums of independent integer random variables. In Proceedings of COLT 2016, pages 831–849,
2016. Full version available at https://arxiv.org/abs/1505.00662.
63

[DKS16c]

I. Diakonikolas, D. M. Kane, and A. Stewart. Robust learning of fixed-structure bayesian
networks. CoRR, abs/1606.07384, 2016.

[DL01]

L. Devroye and G. Lugosi. Combinatorial methods in density estimation. Springer Series in
Statistics, Springer, 2001.

[DLS14]

A. Daniely, N. Linial, and S. S.-Shwartz. From average case complexity to improper learning
complexity. In Symposium on Theory of Computing, STOC 2014, pages 441–448, 2014.

[Fel16a]

V. Feldman. A general characterization of the statistical query complexity.
abs/1608.02198, 2016.

[Fel16b]

V. Feldman. Statistical query learning. In Encyclopedia of Algorithms, pages 2090–2095.
2016.

[FGKP06]

V. Feldman, P. Gopalan, S. Khot, and A. Ponnuswami. New results for learning noisy parities
and halfspaces. In Proc. 47th IEEE Symposium on Foundations of Computer Science (FOCS),
pages 563–576, 2006.

[FGR+ 13]

V. Feldman, E. Grigorescu, L. Reyzin, S. Vempala, and Y. Xiao. Statistical algorithms and a
lower bound for detecting planted cliques. In Proceedings of STOC’13, pages 655–664, 2013.

[FGV15]

V. Feldman, C. Guzman, and S. Vempala. Statistical query algorithms for stochastic convex
optimization. CoRR, abs/1512.09170, 2015.

[FM99]

Y. Freund and Y. Mansour. Estimating a mixture of two product distributions. In Proceedings
of the 12th Annual COLT, pages 183–192, 1999.

[FOS06]

J. Feldman, R. O’Donnell, and R. Servedio. PAC learning mixtures of Gaussians with no
separation assumption. In Proc. 19th Annual Conference on Learning Theory (COLT), pages
20–34, 2006.

[FOS08]

J. Feldman, R. O’Donnell, and R. A. Servedio. Learning mixtures of product distributions
over discrete domains. SIAM J. Comput., 37(5):1536–1564, 2008.

[FPV15]

V. Feldman, W. Perkins, and S. Vempala. On the complexity of random satisfiability problems
with planted solutions. In Proceedings of the Forty-Seventh Annual ACM on Symposium on
Theory of Computing, STOC, 2015, pages 77–86, 2015.

[GHK15]

R. Ge, Q. Huang, and S. M. Kakade. Learning mixtures of gaussians in high dimensions. In
Proceedings of the Forty-Seventh Annual ACM on Symposium on Theory of Computing, STOC
2015, pages 761–770, 2015.

[GR00]

O. Goldreich and D. Ron. On testing expansion in bounded-degree graphs. Technical Report
TR00-020, Electronic Colloquium in Computational Complexity, 2000.

[HK13]

D. Hsu and S. M. Kakade. Learning mixtures of spherical gaussians: moment methods and
spectral decompositions. In Innovations in Theoretical Computer Science, ITCS ’13, pages
11–20, 2013.

[Hot31]

H. Hotelling. The generalization of student’s ratio. Ann. Math. Statist., 2(3):360–378, 08
1931.

64

CoRR,

[HP15]

M. Hardt and E. Price. Tight bounds for learning a mixture of two gaussians. In Proceedings
of the Forty-Seventh Annual ACM on Symposium on Theory of Computing, STOC 2015, pages
753–760, 2015.

[HR09]

P.J. Huber and E. M. Ronchetti. Robust statistics. Wiley New York, 2009.

[HRRS86]

F. R. Hampel, E. M. Ronchetti, P. J. Rousseeuw, and W. A. Stahel. Robust statistics. The
approach based on influence functions. Wiley New York, 1986.

[Hub64]

P. J. Huber. Robust estimation of a location parameter. Ann. Math. Statist., 35(1):73–101, 03
1964.

[JL09]

I. M. Johnstone and A. Y. Lu. On consistency and sparsity for principal components analysis
in high dimensions. Journal of the American Statistical Association, 104(486):682–693, 2009.

[Kea98]

M. Kearns. Efficient noise-tolerant learning from statistical queries. Journal of the ACM,
45(6):983–1006, 1998.

[KK14]

A. R. Klivans and P. Kothari. Embedding hard learning problems into gaussian space. In
Approximation, Randomization, and Combinatorial Optimization. Algorithms and Techniques,
APPROX/RANDOM 2014, pages 793–809, 2014.

[KKMS08] A. Kalai, A. Klivans, Y. Mansour, and R. Servedio. Agnostically learning halfspaces. SIAM
Journal on Computing, 37(6):1777–1805, 2008.
[KMR+ 94] M. Kearns, Y. Mansour, D. Ron, R. Rubinfeld, R. Schapire, and L. Sellie. On the learnability
of discrete distributions. In Proc. 26th STOC, pages 273–282, 1994.
[KMV10]

A. T. Kalai, A. Moitra, and G. Valiant. Efficiently learning mixtures of two Gaussians. In
STOC, pages 553–562, 2010.

[KS06]

A. Klivans and A. Sherstov. Cryptographic hardness for learning intersections of halfspaces.
In Proc. 47th IEEE Symposium on Foundations of Computer Science (FOCS), pages 553–562,
2006.

[KSV08]

R. Kannan, H. Salmasian, and S. Vempala. The spectral method for general mixture models.
SIAM J. Comput., 38(3):1141–1156, 2008.

[KV16]

R. Kannan and S. Vempala. Beyond spectral: Tight bounds for planted gaussians. CoRR,
abs/1608.03643, 2016.

[Li17]

J. Li. Robust sparse estimation tasks in high dimensions. CoRR, abs/1702.05860, 2017.

[LR05]

E. L. Lehmann and J. P. Romano. Testing statistical hypotheses. Springer Texts in Statistics.
Springer, 2005.

[LRV16]

K. A. Lai, A. B. Rao, and S. Vempala. Agnostic estimation of mean and covariance. In
Proceedings of FOCS’16, 2016.

[LS15]

J. Li and L. Schmidt. A nearly optimal and agnostic algorithm for properly learning a mixture
of k gaussians, for any constant k. CoRR, abs/1506.01367, 2015. To appear in COLT’17.

[Moi14]

A. Moitra. Algorithmic aspects of machine learning.
http://people.csail.mit.edu/moitra/docs/bookex.pdf.
65

2014.

Available online at

[MR05]

E. Mossel and S. Roch. Learning nonsingular phylogenies and Hidden Markov Models. In To
appear in Proceedings of the 37th Annual Symposium on Theory of Computing (STOC), 2005.

[MV10]

A. Moitra and G. Valiant. Settling the polynomial learnability of mixtures of Gaussians. In
FOCS, pages 93–102, 2010.

[MW15]

T. Ma and A. Wigderson. Sum-of-squares lower bounds for sparse pca. In Proceedings of the
28th International Conference on Neural Information Processing Systems, NIPS’15, pages
1612–1620, 2015.

[NP33]

J. Neyman and E. S. Pearson. On the problem of the most efficient tests of statistical hypotheses. Philosophical Transactions of the Royal Society of London. Series A, Containing Papers
of a Mathematical or Physical Character, 231(694-706):289–337, 1933.

[O’D14]

R. O’Donnell. Analysis of Boolean Functions. Cambridge University Press, New York, NY,
USA, 2014.

[Sco92]

D.W. Scott. Multivariate Density Estimation: Theory, Practice and Visualization. Wiley, New
York, 1992.

[SD08]

M. S. Srivastava and M. Du. A test for the mean vector with fewer observations than the
dimension. Journal of Multivariate Analysis, 99(3):386 – 402, 2008.

[Sil86]

B. W. Silverman. Density Estimation. Chapman and Hall, London, 1986.

[SOAJ14]

A. T. Suresh, A. Orlitsky, J. Acharya, and A. Jafarpour. Near-optimal-sample estimators for
spherical gaussian mixtures. In Advances in Neural Information Processing Systems (NIPS),
pages 1395–1403, 2014.

[Sze89]

G. Szegö. Orthogonal Polynomials, volume XXIII of American Mathematical Society Colloquium Publications. A.M.S, Providence, 1989.

[Tuk75]

J.W. Tukey. Mathematics and picturing of data. In Proceedings of ICM, volume 6, pages
523–531, 1975.

[Val84]

L. Valiant. A theory of the learnable. Communications of the ACM, 27(11):1134–1142, 1984.

[Ver12]

R. Vershynin. How close is the sample covariance matrix totheactual covariance matrix? Journal of Theoretical Probability, 25(3):655–686, 2012.

[VW02]

S. Vempala and G. Wang. A spectral algorithm for learning mixtures of distributions. In
Proceedings of the 43rd Annual Symposium on Foundations of Computer Science, pages 113–
122, 2002.

[WBS16a]

T. Wang, Q. Berthet, and R. J. Samworth. Statistical and computational trade-offs in estimation
of sparse principal components. Ann. Statist., 44(5):1896–1930, 10 2016.

[WBS16b]

T. Wang, Q. Berthet, and R.J. Samworth. Statistical and computational trade-offs in estimation
of sparse principal components. Ann. Statist., 44(5):1896–1930, 2016.

[Wil97]

R. R. Wilcox. Introduction to robust estimation and hypothesis testing. Statistical modeling
and decision science. Acad. Press, San Diego, Calif. [u.a.], 1997.

[ZB96]

H. Saranadasa Z. Bai. Effect of high dimension: by an example of a two sample problem.
Statist. Sinica,, 6:311–329, 1996.
66

Appendix
A

Sample Complexity Upper Bound for Learning GMMs

In this section, we show that learning a k-mixture of n-dimensional Gaussians to variation distance error ǫ
is easy information theoretically. In particular, we have:
Theorem A.1. Given ǫ > 0 and positive integers k and n, there exists an algorithm that, given a probability
distribution P which is a k-mixture of n-dimensional Gaussians, takes O(n2 k3 log2 (k)/ǫ5 ) samples from P
and with probability at least 2/3 returns a distribution Q with dTV (P, Q) < ǫ.
Note that the algorithm given in Theorem A.1 will not be computationally efficient.
The basic idea of Theorem A.1 will be to make many guesses as to the mixture, at least one of which
is close, and then run a tournament to find the true answer. We approximate the mixture by first guessing
approximations to the weights and then approximating each individual Gaussian. If we had polynomially
many samples from a single part of the mixture, it would be easy to learn:
Lemma A.2 (Folklore). Let G be an n-dimensional Gaussian, and let δ > 0. There exists a polynomial
time algorithm that, given O(n2 /δ2 ) independent samples from G, returns a probability distribution P so
that, with probability at least 2/3, dTV (G, P) < δ.
Note that we can easily improve the success probability in Lemma A.3 to 1− δ at the cost of multiplying
the sample complexity by log(1/δ). In particular, we have:
Corollary A.3. Let G be an n-dimensional Gaussian, and let ǫ > 0. There exists a polynomial time algorithm that given M independent samples from G returns a probability distribution P so that with probability
at least 1 − exp(−Ω(M ǫ2 /d2 )), we have dTV (G, P) < ǫ.
Unfortunately, we cannot simply run this algorithm for each component in our mixture, since we do not
know which samples come from which component. However, if we manage to correctly guess where each
sample comes from this will not be an issue.
Proposition A.4. Given ǫ > 0 and positive integers k and n, there exists an algorithm that given a probability distribution P, which is a k-mixture of n-dimensional Gaussians, and Θ(n2 k3 log(k)/ǫ3 ) independent
samples from P, returns a set of exp(O(n2 k3 log2 (k)/ǫ3 )) distributions Qi so that with probability at least
2/3 there exists an i so that dTV (Pi , Qi ) < ǫ.
Proof. If our algorithm is given N samples, it will return a qi for each function f : [N ] → [k]. Intuitively, f
encodes our guess as to which sample came from which component of the mixture. Note that there are only
exp(N log(k)) many such f ’s.
The algorithm is quite simple. Let s1 , s2 , . . . , sN be our samples, and let Si = {sj : f (j) = i}. Letting
A be the algorithm from Corollary A.3 with δ taken to be ǫ/(10k), we let
Qf =


k 
X
|Si |
i=1

N

A(Si ).

We claim that at least one of these works
Pkwith probability 2/3.
In particular, let P be the mixture i=1 wi Gi . Consider the case where f correctly guesses which part
of the mixture each sample was taken from. In particular, f (i) = j if and only if si was taken from Gj . We
claim that, with probability at least 2/3, this choice of f leads to dTV (Qf , P) < ǫ. We will henceforth use
Si to denote the set of samples actually taken from the ith component of the mixture.

67

Firstly, note that by standard concentration bounds, we have that |SNi | − wi < ǫ/(10k), for all i, with
probability at least 9/10.
Secondly, note that after conditioning on which samples of P were taken form which part, the samples
themselves are independent samples from the appropriate Gi ’s, that with probability at least 9/10 we have
that dTV (Gi , A(Si )) < ǫ/(10k) for all i with |Si | ≫ n2 k2 log(k)/ǫ2 .
We claim that if both of the above conditions hold (which happens with probability at least 2/3) that
dTV (P, Qf ) < ǫ. Letting S be the set of indices i so that wi < ǫ/(5k), we note that for i 6∈ S that
|Si | ≫ n2 k2 log(k)/ǫ2 . We then have that

k 
k
X
1 X
|Si |
dTV (P, Qf ) =
wi Gi −
A(Si )
2
N
i=1
i=1
1


k
X
|Si |
1
≤
wi Gi −
A(Si )
2
N
2
i=1
 



k
k
|Si |
|Si |
1X
1X
max wi ,
wi −
+
|Gi − A(Si )|1
≤
2
N
2
N
i=1
i=1
X
1X
(wi + ǫ/(10k))|Gi − A(Si )|1 +
dTV (Gi , A(Si ))
≤ ǫ/20 +
2
i∈S

X
1X
≤ ǫ/20 +
3ǫ/(10k) +
ǫ/(10k)
2
i∈S

i6∈S

i6∈S

≤ ǫ/20 + 3ǫ/20 + ǫ/20

<ǫ.
This completes the proof.

Theorem A.1 now follows immediately form a standard tournament argument (see, e.g., [DL01, DDS12b,
DDS15]).

B Sample Complexity Upper Bound for Parameter Estimation of Separated GMMs
Next we consider P
the more complicated task of parameter estimation. In particular, given samples from a
distribution P = ki=1 Gi , where each Gi is a weighted Gaussian, we would like to learn a distribution
P
Q that is not only close to P but that can be written as Q = ki=1 Hi with kHi − Gi k1 small for all i.
Now, in general, this task will require number of samples exponential in k, simply because there are pairs
of mixtures that are ǫΩ(k) -close in variation distance and yet ǫ-far in terms of their individual components.
However, we will show that if the components are separated, this cannot be the case and thus learning the
distribution in variation distance will be sufficient.
Before we begin, we need to clarify
R our notion of separation. Given two pseudo-distributions, p and q,
we define their overlap as V (p, q) := min(dp, dq). We should note that if p and q are honest distributions,
then dTV (p, q) = 1 − V (p, q). We have the following theorem:
P
P
Theorem B.1. There exists a constant C so that if we have two mixtures p = ki=1 Gi and q = ki=1 Hi ,
where p and q are normalized distributions with Hi , Gi weighted Gaussians, such that dTV (p, q) < (δ/k)C
for some sufficiently small δ > 0, and so that for any i 6= j, V (Gi , Gj ), V (Hi , Hj ) < (δ/k)C , then there
exists a permutation π : [k] → [k] so that kGi − Hπ(i) k1 < δ for all i.

68

We begin by producing a proxy for the overlap between distributions. In particular, for pseudo-distributions
p and q, we define
Z

p
h(p, q) = − log
dpdq .
Notice that if p and q are true distributions, this is related to the Hellinger distance by H(p, q) = 2(1 −
e−h(p,q) ). We also note the relationship to the overlap:
Lemma B.2. If p and q are pseudo-distributions with L1 norm at most 1, then
V (p, q) = exp(−Θ(h(p, q)) + O(1)) .
Proof. On the one hand, there is an easy upper bound
Z
Z p
dpdq = exp(−h(p, q)).
V (p, q) = min(dp, dq) ≤

The lower bound is by Cauchy-Schwarz
exp(−h(p, q)) =

Z p

This completes our proof.

dpdq ≤

Z

1/2 Z
1/2
p
≤ 2V (p, q).
min(dp, dq)
max(dp, dq)

Ideally we would like to show that h is nearly a metric for Gaussians. Namely that h(A, C) =
O(h(A, B) + h(B, C)). This would imply that Hi could not have large overlap with more than one Gj ,
since if V (Hi , Ga ) and V (Hi , Gb ) were both large, then h(Hi , Ga ), h(Hi , Gb ) would be small and therefore, h(Ga , Gb ) would be small. This would contradict our assumption that Ga and Gb have small overlap.
Unfortunately, this is not true. In one dimension, a very wide Gaussian may have non-trivial overlap with
two narrow Gaussians with widely separated means, neither of which overlaps the other substantially. We
will need to develop techniques to deal with this circumstance.
To do this, we introduce an intermediate notation. If Gi = wi N (µi , Σi ) are weighted Gaussians, we
define
hΣ (G1 , G2 ) := h(N (0, Σ1 ), N (0, Σ2 )).
This is useful because it does satisfy an approximate triangle inequality.
Proposition B.3. For F, G, H weighted Gaussians, we have that
hΣ (F, H) = O(hΣ (F, G) + hΣ (G, H)).
Before we prove this, we will first need to find an approximation to hΣ .
Lemma B.4. If G and H are weighted Gaussians with covariance matrices A and B respectively, then


hΣ (G, H) = Θ Σλ eigenvalue of B −1/2 AB −1/2 min(| log(λ)|, | log(λ)|2 ) .

Proof. By making an appropriate change of variables, we can assume that H has identity covariance and G
has covariance B −1/2 AB −1/2 . Thus, it suffices to consider the case where B = I. In this case, we may

69

diagonalize A to get A = diag(λi ). We then have that
hΣ (G, H) = h(N (0, A), N (0, I))
= − log (2π)−n/2
= − log
=

n
X

n
Y

−1/4

λi

n
Y

−1/4

λi

i=1

Z

exp −

−1/2
((1 + λ−1
i )/2)

i=1

!

n
X

!

x2i /2(1/(2λi ) + 1/2) dx

i=1

!

log(λi )/4 + log((1 + λ−1
i )/2)/2.

i=1

We claim that
log(λ)/4 + log((1 + λ−1 )/2)/2 = Θ(min(| log(λ)|, | log(λ)|2 )).
To see this note that when λ = 1 + ǫ for small values of ǫ, the left hand side above is
(ǫ/4 − ǫ2 /8 + O(ǫ3 )) + (−ǫ/4 + 3ǫ/16 + O(ǫ3 )) = ǫ2 /16 + O(ǫ3 ) = Θ(ǫ2 ).
On the other hand, when λ ≫ 1, this is asymptotic to | log(λ)|/4, and when λ ≪ 1, it is similarly asymptotic
to − log(λ)/4. Finally, since it is easily verified that log(λ)/4+log((1+λ−1 )/2)/2 is never 0 unless λ = 1,
this proves the claim, from which our lemma follows easily.
We will also need the following fact about eigenvalues of a product of matrices:
Lemma B.5. Let A and B be symmetric matrices with eigenvalues ν1 ≥ ν2 ≥ . . . ≥ νn > 0 and µ1 ≥
µ2 ≥ . . . ≥ µn > 0, respectively. Let λ1 ≥ λ2 ≥ . . . ≥ λ2n > 0 be the sorting of the νi and µi together.
Let M be a matrix with M T M = A. Then, the kth largest eigenvalue of M T BM is at most λ2k .
Proof. We need to show that there is an (n − k + 1)-dimensional subspace V so that for v ∈ V we have
that vA1/2 BA1/2 v ≤ λ2k |v|2 . Suppose that λ1 , . . . , λk−1 contains m of the νi and k − m − 1 of the µi . Let
V be the subspace of vectors v so that v is perpendicular to the top m eigenvectors of A and so that A1/2 v
is perpendicular to the top m − k − 1 eigenvalues of B. Then
vM T BM v ≤ λk |M v|2 = λk vAv ≤ λ2k |v|2 .
This completes the proof.
We are now ready to prove Proposition B.3.
Proof. Let F, G, H have covariance matrices A, B, C respectively. Let Σ1 = A−1/2 BA−1/2 , Σ2 =
B −1/2 CB −1/2 and Σ3 = A−1/2 CA−1/2 = (A−1/2 B 1/2 )Σ2 (B 1/2 A−1/2 ). Let the eigenvalues of Σi be
(i)
(i)
(i)
λ1 ≥ λ2 ≥ . . . ≥ λn > 0. Let f (x) = max(0, min(log(x), log 2 (x))). We have by Lemma B.4 that
!
n
X
(1)
(1)
f (λi ) + f (1/λi ) ,
hΣ (F, G) = Θ
i=1

hΣ (G, H) = Θ

n
X

(2)
f (1/λi )

(3)
f (λi )

(3)
f (1/λi )

i=1

hΣ (F, H) = Θ

n
X

!

(2)
f (λi ) +

i=1

70

+

!

,

.

(3)

(1)

On the other hand, Lemma B.5 says that λi is at most the square of the ith largest of the λj
Therefore,
!
n
n
n
X
X
X
(2)
(1)
(3)
f (λi ) .
f (λi ) +
f (λi ) = O
i=1

i=1

i=1

(2)

and λj .

Similarly, by considering the inverses of these matrices, we find that
n
X

(3)

f (1/λi ) = O

n
X

(1)

f (1/λi ) +

(2)

!

f (1/λi ) .

i=1

i=1

i=1

n
X

Together these complete the proof.
In addition to this, we need to know what else contributes to h(G, H). We define
hµ (G, H) = h(G, H) − hΣ (G, H).
We make the following claim:
Proposition B.6.
−1
hµ (w1 N (µ1 , Σ1 ), w2 N (µ2 , Σ2 )) = −1/2 log(w1 w2 )+inf ((x−µ1 )Σ−1
1 (x−µ1 )+(x−µ2 )Σ2 (x−µ2 ))/4.
x

Proof. We have that
h(w1 N (µ1 , Σ1 ), w2 N (µ2 , Σ2 ))

Z
−1
−1
−n/2
−1/4
(2π)
(det(Σ1 Σ2 ))
exp(−((x − µ1 )Σ1 (x − µ1 ) + (x − µ2 )Σ2 (x − µ2 ))/4)dx .
= − 1/2 log(w1 w2 ) − log
−1
Letting x0 achieve the minimum value of ((x − µ1 )Σ−1
1 (x − µ1 ) + (x − µ2 )Σ2 (x − µ2 ))/4, this is
−1
−1/2 log(w1 w2 ) + ((x0 − µ1 )Σ−1
1 (x0 − µ1 ) + (x0 − µ2 )Σ2 (x0 − µ2 ))/4
Z

−1
−1
−n/2
−1/4
+ log
(2π)
(det(Σ1 Σ2 ))
exp(−(x − x0 )(Σ1 + Σ2 )(x − x0 )/4)dx .

Noting that the term at the end is simply hΣ (w1 N (µ1 , Σ1 ), w2 N (µ2 , Σ2 )) completes the proof.
We need one further proposition from which Theorem B.1 will follow easily.
Proposition B.7. Under
the assumptions of Theorem B.1, for each i there exists at most one j so that
√
C
h(Gi , Hj ) > (δ/k) .
To prove this, we will need one further lemma:
√

Lemma B.8. If h(Gi , Hj ) > (δ/k) C , with ΣG and ΣH the covariance matrices of the corresponding
Gaussians, then for A a sufficiently large constant (independent of C) ΣG ≤ AΣH .
Proof. Suppose for sake of contradiction that this is not the case. By making a change of variables, we can
′
assume that ΣG = I. This
p means that ΣH has some eigenvector v with eigenvalue less than 1/A. Let H be
1/4
log(k/δ) in the direction closer to the mean of Gi . We have that
Hj translated by C
√
h(Hj , H ′ ) = hµ (Hj , H ′ ) = Θ(A C log(k/δ)).
71

Therefore, V (Hj , H ′ ) = (δ/k)Ω(A

√

C) .

On the other hand,

h(Gi , H ′ ) = hΣ (Gi , H ′ ) + hµ (Gi , H ′ )
√
≤ hΣ (Gi , H) + O( C log(k/δ))
√
≤ h(Gi , H) + O( C log(k/δ))
√
= O( C log(k/δ)).
√

This means that V (G, H ′ ) = (δ/k)O( C) .
√
P
P
This means that p = ℓ Gℓ has V (p,√H ′ ) = (δ/k)O( C) , and since p is close to q = ℓ Hℓ , there must
be some ℓ so that V (H ′ , Hℓ ) > (δ/k)O( C) . Note that ℓ here cannot be j. On the other hand, this implies
that
h(Hj , Hℓ ) = hΣ (Hj , Hℓ ) + hµ (Hj , Hℓ )

√
≤ hΣ (H ′ , Hℓ ) + O(hµ (H ′ , Hℓ ) + A C log(k/δ))
√
≤ O(h(H ′ , Hℓ ) + A C log(k/δ))
√
= O(A C log(k/δ)).

Therefore, V (Hℓ , Hj ) = (δ/k)O(A
the proof.

√

C) ,

which for C ≫ A2 contradicts our assumptions. This completes

We are now prepared to prove Proposition B.7.
√

Proof. Suppose for sake of contradiction that V (Gi , Hj ), V (Gi , Hℓ ) > (δ/k) C for some j 6= ℓ. Then, by
Lemma B.8, we have that all of the covariance matrices of Gi , Hj , Hℓ are comparable to each other (namely
each is no more than a constant√multiple of any other). We claim that this implies that hµ (Hj , Hℓ ) =
O(hµ (Hj , Gi ) + hµ (Hℓ , Gi ) + C log(k/δ)). This is because, letting ΣG be the covariance matrix of Gi ,
′
′
and letting
√ wG = |Gi |1 , wH = |Hj |1 , wH = |Hℓ |1 , we have the following: First, each of wG , wH , wH =
exp(O( C log(δ/k))), because each distribution has large overlap with some other. Next, we have that
√
hµ (Gi , Hj ) = O( C log(δ/k)) + inf Θ((x − µGi )ΣG (x − µGi ) + (x − µHj )ΣG (x − µHh )
x
√
= O( C log(δ/k)) + Θ((µGi − µHj )ΣG (µGi − µHj )).
Similarly,

√
hµ (Gi , Hℓ ) = O( C log(δ/k)) + Θ((µGi − µHℓ )ΣG (µGi − µHℓ )),

and

√
hµ (Hj , Hℓ ) = O( C log(δ/k)) + Θ((µHj − µHℓ )ΣG (µHj − µHℓ )).

This implies that
hµ (Hj , Hℓ ) = O(hµ (Hj , G) + hµ (G, Hℓ )).
Therefore, we have that
h(Hj , Hℓ ) = hΣ (Hj , Hℓ ) + hµ (Hj , Hℓ )
= O(hΣ (Hj , G) + hµ (Hj , G) + hΣ (G, Hℓ ) + hµ (G, Hℓ ))
= O(h(Hj , G) + h(G, Hℓ ))
= O(log(1/V (Hj , G)) + log(1/V (G, Hℓ )))
√
= O( C log(k/δ)).
72

√

However, this implies that V (Hj , Hℓ ) = (δ/k)O(
This completes our proof.

C) ,

a contradiction.

We are now ready to prove Theorem B.1
√

Proof. For each Gi that has overlap more than (δ/k) C with some Hj , let π(i) be that j. For other i, define
π(i) arbitrarily subject to π being √
a permutation.
Note that V (Gi , Hj ) < (δ/k) C for any j 6= π(i). Also note that
V (Gi , q) ≥ V (Gi , p) − |p − q|1 = |Gi |1 − 2(δ/k)C .
On the other hand,
V (Gi , q) ≤

X

V (Gi , Hj )

j

≤ V (Gi , Hπ(i) ) +

X

V (Gi , Hj )

j6=π(i)

≤ V (Gi , Hπ(i) ) + δ/3.
Therefore V (Gi , Hπ(i) ) ≥ |Gi |1 − δ/2. It is also at most |Gi |1 − δ/2. On the other hand |Gi − Hπ(i) |1 =
|Gi |1 + |Hπ(i) |1 − 2V (Gi , Hπ(i) ) ≤ δ. This completes the proof.

C

Testing the Mean of a High-Dimensional Gaussian

√
Theorem C.1. There exists an algorithm that given ǫ > 0 and k = O( n/ǫ2 ) samples from an ndimensional Gaussian G = N (µ, I) distinguishes between the cases
• µ=0
• kµk2 > ǫ
with probability at least 2/3.
Proof. The tester is fairly simple. Let Xi be the ith sample, and let
k
1 X
Xi .
Z := √
k i=1

The algorithm returns “YES” if kZk22 < ǫ2 k/2 + n and “NO”
√ otherwise.
To show correctness, note that Z is distributed as N (µ k, I). If µ = 0, then kZk22 has mean n and
variance O(n), and so it is less than n + ǫ2k/2 with probability at least 2/3, assuming that k is a sufficiently
√
large multiple of n/ǫ2 . On the other hand, if kµk2 > ǫ, we note that kZk22 has mean n + kkµk22 and
√
variance O(n) + O(kkµk22 ). Thus, if kkµk22 ≫ n, the algorithm rejects with probability 2/3. Again, this
√
happens if kµk2 > ǫ and k is a sufficiently large multiple of n/ǫ2 . This completes the proof.
We also note that this tester can be implemented in the SQ model simply by verifying that each
√
coordinate-wise median has absolute value less than ǫ/ n, which can be verified by showing that Pr(xi >
√
0) = 1/2 + O(ǫ/ n).
We also show that the tester above is sample-optimal, up to a constant factor:
√
Theorem C.2. There is no algorithm that given k = o( n/ǫ2 ) samples from an n-dimensional Gaussian
G = N (µ, I) distinguishes between the cases
73

• µ=0
• kµk2 > ǫ
with probability at least 2/3.
Proof. Suppose for sake of contradiction that such an algorithm does exist. Consider the following scenario:
√
Let µ be taken from the distribution N (0, (2ǫ/ n)I). Note that kµk2 > ǫ with probability at least 9/10.
Let Y1 , Y2 , . . . , Yk be independent samples taken from N (µ, I). And let Z1 , . . . , Zk be independent samples
from N (0, I). Assuming that our algorithm exists, it can distinguish between a sample from Y1 , . . . , Yk
and a sample from Z1 , . . . , Zk with probability better than 1/2. This means that these distributions must
have constant variational distance. However, note that the vector (Z1 , . . . , Zk ) is simply a standard nkdimensional Gaussian. The vector (Y1 , . . . , Yk ) on the other hand is an nk-dimensional Gaussian with mean
0 and with

2

1 + 2ǫ /n , if ab = cd
Cov(Yab , Ycd ) = 2ǫ2 /n ,
if b = d and a 6= c


0,
otherwise .
By standard results, G′ = N (0, Σ) has constant variation distance from N (0, I) if and only if kΣ−IkF ≫ 1.
Taking Σ to be the covariance matrix for the Y ’s, we have that
√
kΣ − Ik2F = nk(2ǫ/ n)2 = 4kǫ2 /n = o(1) .

This implies that the distribution on Y ’s is close, in total variation distance, to the distribution on Z’s, and
gives a contradiction.

D

Omitted Proofs

D.1 Proof of Fact 3.6 We will need the following claim:
Claim D.1. We have that:
Hei (x cos θ + y sin θ) =

i  
X
i

j

j=1

cosj θ sini−j θHej (x)Hei−j (y) .

Proof. The Hei (x) are monic polynomials: the lead term is xi with coefficient 1. Thus, all the degree-i
terms of Hei (x cos θ + y sin θ) are given by
i  
X
i
(x cos θ + y sin θ)i
cosj θ sini−j θxj y i−j .
j
j=1

It follows that the degree-i terms of the LHS and RHS of the lemma agree. Therefore, we have
i  
X
i
Hei (x cos θ + y sin θ) = p(x, y) +
cosj θ sini−j θHej (x)Hei−j (y) ,
j
j=1

for some polynomial p(x, y) of degree at most i − 1. We need to show that p(x, y) is identically zero. To
show this we consider E[Hei (X cos θ + Y sin θ)2 ], for (X, Y ) ∼ N (0, I). Since the Gaussian is unaltered
by rotations, by a change of coordinates we have that:
Z ∞Z ∞
2
E[Hei (X cos θ + Y sin θ) ] =
Hei (x cos θ + y sin θ)2 G(x)G(y)dxdy
−∞
Z−∞
∞ Z ∞
=
Hei (x′ )2 G(x′ )G(y ′ )dx′ dy ′ = i! .
−∞

−∞

74

However, pairs of distinct Hej (x)Hei−j (y) are orthogonal to each other and they are all orthogonal to the
lower degree polynomial p(x, y). Thus, we have
E[Hei (X cos θ + Y sin θ)2 ] = E[p(X, Y )2 ] +

i  2
X
i
j=1

= E[p(X, Y )2 ] +

j

i  2
X
i
j=1

= E[p(X, Y )2 ] + i!

j

cos2j θ sin2(i−j) θi!(i − j)!

i  
X
i
j=1

2

cos2j θ sin2(i−j) θE[Hej (X)2 Hei−j (Y )2 ]

j

cos2j θ sin2(i−j) θ

= E[p(X, Y ) ] + i!(cos2 θ + sin2 θ)i = E[p(X, Y )2 ] + i! .
We must therefore have that E[p(X, Y )2 ] = 0. Since the Gaussian has positive pdf everywhere, this implies
that p(X, Y ) is identically zero.
We now have:
Uθ (Hei G)(x) =
=
=

Z

∞

Z−∞
∞

Hei (x cos θ + y sin θ)G(x cos θ + y sin θ)G(x sin θ − y cos θ)dy
Hei (x cos θ + y sin θ)G(x)G(y)dy

−∞
i 
X
j=1
i


Z ∞
i
j
i−j
Hej (x)Hei−j (y)G(x sin θ − y cos θ)dy
cos θ sin θ
j
−∞

= cos θHei (x)G(x) ,
since

R∞

−∞ Hei−j (y)G(y)dy

= δij . This completes the proof.

D.2 Proof of Lemma 3.7 We use the following lemma:
Lemma D.2 (Proposition 1 from [CFJ13]). Given any 0 < ǫ < π/2, let θ be the angle between two random
unit vectors uniformly distributed over Sn . Then we have that:
√
Pr[|θ − π/2| ≥ ǫ] ≤ O( n(cos ǫ)n−2 ) .
As a corollary, we have:
Corollary D.3. Let θ be the angle between two random unit vectors uniformly distributed over Sn . Then we
have that:



Pr | cos θ| ≥ Ω(n−α ) ≤ exp −Ω(n1−2α ) ,
for any 0 ≤ α ≤ 1/2.

Proof. We apply Lemma D.2 with ǫ = n−α . If n−α = O(1), the result is trivial, so we may assume that
ǫ ≤ 1/100. Then we have that cos ǫ ≤ 1 − ǫ2 /2 + ǫ2 /24 ≤ 1 − ǫ2 /3 ≤ exp(ǫ2 /4). Lemma D.2 now gives
that



√
Pr |θ − π/2| ≥ n−α ≤ O n exp(−n−2α /4)n−2 ≤ exp(−n1−2α /5) .
Note that if |θ − π/2|≤n−α , it follows that | cos θ| ≤ |θ − π/2| ≤ n−α .
75

Using Corollary D.3 for α = 1/2 − c, and a union bound over all pairs of distinct vectors in S, the
probability that there exist v 6= v ′ ∈ S such that |v · v ′ |≥Ω(nc−1/2 ) is less than
|S|2 2−Ω(n

2c )

<1.

Therefore, the set S will satisfy the statement of Lemma 3.7 with positive probability, as desired.
D.3 Proof of Fact 6.2 We have:
Z
2
1 + χ (wB + (1 − w)C, D) = (wB(x) + (1 − w)C(x))2 /D(x)dx
Z
Z
Z
2
2
2
2
=w
B(x) /D(x)dx + (1 − w)
C(x) /D(x)dx + 2w(1 − w) B(x)C(x)/D(x)dx
= w2 (1 + χ2 (B, D)) + (1 − w)2 (1 + χ2 (C, D)) + 2w(1 − w)(1 + χD (B, C))

= 1 + w2 χ2 (B, D) + (1 − w)2 χ2 (C, D) + 2w(1 − w)χD (B, C) .
This completes the proof.

D.4 Proof of Fact 6.3 By definition, we can write
Z ∞
′
G(x − µ′ )G(x − µ)/G(x)dx
1 + χN (0,1) (N (µ , 1), N (µ, 1)) =
−∞
Z ∞
√

= (1/ 2π) ·
exp −(x − µ′ )2 /2 − (x − µ)2 /2 + x2 /2 dx
Z−∞
∞
√

exp −x2 /2 + (µ′ + µ)x − µ2 /2 − µ′2 /2 dx
= (1/ 2π) ·
−∞
Z ∞

G(x − µ − µ′ ) exp −µ2 /2 + −µ′2 /2 + (µ + µ′ )2 /2 dx
=
−∞

= exp(µ′ µ) .

This completes the proof.
D.5 Proof of Fact 6.4 By definition, we have that
Z ∞
2
2
G(x/σ)2 /G(x)dx
1 + χ (N (0, σ ), N (0, 1)) = (1/σ)
−∞
Z ∞
1
= √
exp(x2 /2 − x2 /σ 2 )dx
2
σ 2π −∞
p
Z ∞
p
2/σ 2 − 1
1
p
G(x/
2/σ 2 − 1)dx
=
·
σ2
2/σ 2 − 1
−∞
p
2/σ 2 − 1 p
=
= 2/σ 4 − 1/σ 2 .
σ2

This completes the proof.

76

