List-Decodable Robust Mean Estimation and Learning Mixtures of Spherical Gaussians

Introduction

This paper is concerned with the problem of efficiently learning high-dimensional spherical Gaussians in the presence of a large fraction of corrupted data, and in the related problem of parameter estimation for mixtures of high-dimensional spherical Gaussians (henceforth, spherical GMMs). Before we state our main results, we describe and motivate these two fundamental learning problems.

The first problem we study is the following:

A few remarks are in order: We first note that we make no assumptions on the remaining $(1-\alpha)$-fraction of the points in $T$. These points can be arbitrary and may be chosen by an adversary that is computationally unbounded and is allowed to inspect the set of good points. We will henceforth call such a set of points $\alpha$-corrupted. Ideally, we would like to output a single hypothesis vector ${\widehat{\mu}}$ that is close to $\mu$ (with high probability). Unfortunately, this goal is information-theoretically impossible when the fraction $\alpha$ of good samples is less than $1/2$. For example, if the input distribution is a uniform mixture of $1/\alpha$ many Gaussians whose means are pairwise far from each other, there are $\Theta(1/\alpha)$ different valid answers and the list must by definition contain approximations to each of them. It turns out that the information-theoretically best possible size of the candidates list is $s=\Theta(1/\alpha)$. Therefore, the feasible goal is to design an efficient algorithm that minimizes the Euclidean distance between the unknown $\mu$ and its closest ${\widehat{\mu}}_{i}$.

Before we proceed with a detailed background and motivation, we point out the connection between these two problems. Intuitively, Problem 2 can be reduced to Problem 1, as follows: We can think of the samples drawn from a spherical GMM as a set of corrupted samples from a single Gaussian — where the Gaussian in question can be any of the mixture components. The output of the list decoding algorithm will produce a list of hypotheses with the guarantee that every mean vector in the mixture is relatively close to some hypothesis. If in addition the distances between the means and their closest hypotheses are substantially smaller than the distances between the means of different components, this will allow us to reliably cluster our sample points based on which hypothesis they are closest to. We can thus cluster points based on which component they came from, and then we can learn each component independently.

2 List-Decodable Robust Learning

The vast majority of efficient high-dimensional learning algorithms with provable guarantees make strong assumptions about the input data. In the context of unsupervised learning (which is the focus of this paper), the standard assumption is that the input points are independent samples drawn from a known family of generative models (e.g., a mixture of Gaussians). However, this simplifying assumption is rarely true in practice and it is important to design estimators that are robust to deviations from their model assumptions.

This phenomenon had raised the following question: Can we reconcile computational efficiency and robustness in high dimensions? Recent work in the TCS community made the first algorithmic progress on this front: Two contemporaneous works [DKK+16, LRV16] gave the first computationally efficient robust algorithms for learning high-dimensional Gaussians (and many other high-dimensional models) with error close to the information-theoretic optimum. Specifically, for the problem of robustly learning an unknown mean Gaussian $N(\mu,I)$ from an $\eta$-corrupted set of samples, $\eta<1/2$, we now know a polynomial-time algorithm that achieves the information-theoretically optimal error of $O(\eta)$ [DKK+17b].

The aforementioned literature studies the setting where the fraction of corrupted data is relatively small (smaller than $1/2$), therefore the real data is the majority of the input points. A related setting of interest focuses on the regime when the fraction $\alpha$ of real data is small — strictly smaller than $1/2$. From a practical standpoint, this “large error regime” is well-motivated by a number of pressing machine learning applications (see, e.g., [CSV17, SVC16, SKL17]). From a theoretical standpoint, understanding this regime is of fundamental interest and merits investigation in its own right. A specific motivation comes from a previously observed connection to learning mixture models: Suppose we are given samples from the mixture $\alpha\cdot N(\mu,I)+(1-\alpha)E$, i.e., $\alpha$-fraction of the samples are drawn from an unknown Gaussian, while the rest of the data comes from several other populations for which we have limited (or no) information. Can we approximate this “good” Gaussian component, independent of the structure of the remaining components?

The main focus of this work is on the fundamental setting where the good data comes from a Gaussian distribution. Specifically, we ask the following question:

What is the best possible error guarantee (information-theoretically) achievable for list-decodable mean estimation, when the true distribution is an unknown $N(\mu,I)$? More importantly, what is the best error guarantee that we can achieve with a computationally efficient algorithm?

As our first main result, we essentially resolve Question 1.1.

3 Learning Mixtures of Separated Spherical Gaussians

In this paper, we focus on the natural and important case where each of the components is spherical, i.e., each covariance matrix is an unknown multiple of the identity. The majority of prior algorithmic work on this problem studied the setting where there is a minimum separation between the means of the componentsWithout any separation assumptions, it is known that the sample complexity of the problem becomes exponential in the number of components [MV10, HP15].. For the simplicity of this discussion, let us consider the case that the mixing weights are uniform (i.e., equal to $1/k$, where $k$ is the number of components) and each component has identity covariance. (We emphasize that the positive results of this paper hold for the general case of an arbitrary mixture of high-dimensional spherical Gaussians, and apply even in the presence of a small dimension-independent fraction of corrupted data.) The problem of learning separated spherical GMMs was first studied by Dasgupta [Das99], followed by a long line of works that obtained efficient algorithms under weaker separation assumptions.

Interestingly enough, until very recently, even the information-theoretic aspect of this problem was not understood. Specifically, what is the minimum separation that allows the problem to be solvable with $\operatorname{poly}(n,k,1/\delta)$ samples? Recent work by Regev and Vijayraghavan [RV17] characterized this aspect of the problem: Specifically, [RV17] showed that the problem of learning spherical $k$-GMMs (with equal weights and identity covariances) can be solved with $\operatorname{poly}(n,k,1/\delta)$ samples if and only if the means are pairwise separated by at least $\Theta(\sqrt{\log k})$. Unfortunately, the approach of [RV17] is non-constructive in high dimensions. Specifically, they gave a sample-efficient learning algorithm whose running time is exponential in the dimension. This motivates the following question:

Is there a $\operatorname{poly}(n,k)$ time algorithm for learning spherical $k$-GMMs with separation $o(k^{1/4})$, or better $O(k^{\varepsilon})$, for any fixed $\varepsilon>0$? More ambitiously, is there an efficient algorithm that succeeds under the information-theoretically optimal separation?

As our second main result, we make substantial progress towards the resolution of Question 1.2.

4 Our Contributions

In this paper, we develop a set of techniques that yield new efficient algorithms with significantly better guarantees for Problems 1 and 2. Our algorithms depend in an essential way on the analysis of high degree multivariate polynomials. We obtain a detailed structural understanding of the behavior of high degree polynomials under the standard multivariate Gaussian distribution, and leverage this understanding to design our learning algorithms. More concretely, our main technical contribution is a new technique, using degree-$d$ multivariate polynomials, to remove outliers from high-dimensional datasets where the majority of the points are corrupted.

Our main result is an efficient algorithm for list-decodable Gaussian mean estimation with a significantly improved error guarantee:

A major complication that occurs in the regime of $\alpha<1/2$ is that since fewer than half of our samples are good, the values of such a polynomial $p$ might concentrate in several clusters. As a consequence, we will not necessarily be able to identify which cluster contains the good samples. In order to deal with this issue, we need to develop new techniques for outlier removal that handle the setting that the good data is a small fraction of our dataset. Roughly speaking, we achieve this by performing a suitable clustering of points based on the values of $p$, and returning multiple (potentially overlapping) subsets of our original dataset $T$ with the guarantee that at least one of them will be a cleaner version of $T$. This new paradigm for performing outlier removal in the large error regime may prove useful in other contexts as well.

A crucial technical contribution of our approach is the use of degree more than one polynomials for outlier removal in this setting. The intuitive reason for using polynomials of higher degree is this: A small fraction of points that are far from the true mean in some particular direction will have a more pronounced effect on higher degree moments. Therefore, taking advantage of the information contained in higher moments should allow us to discern smaller errors in the distance from the true mean. The difficulty is that it is not clear how to algorithmically exploit the structure of higher degree moments in this setting.

The major obstacle is the following: Since we do not know the mean $\mu$ of $G$ — this is exactly the quantity we are trying to approximate! — we are also not able to evaluate the variance $\operatorname*{Var}[p(G)]$ of $p(G)$. If $p$ was a degree-$1$ polynomial, this would not be a problem, as the variance $\operatorname*{Var}[p(G)]$ does not depend on $\mu$. But for degree at least $2$ polynomials, the dependence of $\operatorname*{Var}[p(G)]$ on $\mu$ becomes a fundamental difficulty. Thus, although we can potentially find polynomials with unexpectedly large empirical variance, we will have no way of knowing whether this is due to corrupted points $x\in T$ (on which $p(x)$ is abnormally far from its true mean), or due to errors in our estimation of the mean of $G$ causing us to underestimate the variance $\operatorname*{Var}[p(G)]$.

In order to circumvent this difficulty, we require a number of new ideas, culminating in an algorithm that allows us to either verify that the variance of $p(G)$ is close to what we are expecting, or to find some other polynomial that allows us to remove outliers.

We leverage the connection between list-decodable learning and learning mixture models to obtain an efficient algorithm for learning spherical GMMs under much weaker separation assumptions. Specifically, by using the algorithm of Theorem 1.3 combined with additional algorithmic ideas, we obtain our second main result:

for all $i\neq j$, for $C>0$ a sufficiently large constant, the algorithm draws $\operatorname{poly}(n,{(dk/\delta)^{d}})$ samples from $F$, runs in time $\operatorname{poly}(n,{(dk/\delta)^{d}})$, and with high probability returns a list $\{(u_{i},\nu_{i},s_{i}),i\in[k]\}$, such that the following conditions hold (up to a permutation): $|u_{i}-w_{i}|=O(\delta)$, $\|\mu_{i}-\nu_{i}\|_{2}/\sigma_{i}=O(\delta/w_{i})$, and $|s_{i}-\sigma_{i}|/\sigma_{i}=O(\delta/w_{i})/\sqrt{n}$.

The reader is also referred to Proposition 4.3 for a more detailed statement that also allows a small, dimension-independent fraction of adversarial noise in the input samples.

We now provide an intuitive explanation of our spherical GMM learning algorithm. First, we note that we can reduce the dimension of the problem from $n$ down to some function of $k$. When the covariance matrices of the components are nearly identical, this can be done with a twist of standard techniques. For the case of arbitrary covariances, we need to employ a few additional ingredients.

When each component has the same covariance matrix, the learning algorithm is quite simple: We start by running our list-decoding algorithm (Theorem 1.3) with appropriate parameters to get a small list of hypothesis means. We then associate each sample with the closest element of our list. At this point, we can cluster the points based on which means they are associated to and use this clustering to accurately learn the correct components.

The general case, when the covariances of the components are arbitrary, is significantly more complex. In this case, we can recover a list $H$ of candidate means only after first guessing the radius of the component that we are looking for. Without too much difficulty, we can find a large list of guesses and thereby produce a list of hypotheses of size $\operatorname{poly}(n/\alpha)$. However, clustering based on this list now becomes somewhat more difficult, as we do not know the radius at which to cluster. We address this issue by performing a secondary test to determine whether or not the cluster that we have found contains many points at approximately the correct distance from each other.

As mentioned in Section 1.2, even the following information-theoretic aspect of list-decodable mean estimation is open: Ignoring sample complexity and running time, how small a distance from the true mean can be achieved with $\operatorname{poly}(1/\alpha)$ many hypotheses or number of hypotheses that is only a function of $\alpha$, i.e., independent of the dimension $n$?

Let $0<\alpha<1/2$. There exists an (inefficient) algorithm that given a set of $\alpha$-corrupted samples from a distribution $D$, where (a) $D$ is subgaussian with bounded variance in each direction, or (b) $D$ has bounded first $k$ moments, for even $k$, outputs a list of $O(1/\alpha)$ vectors one of which is within distance $g(\alpha)$ from the mean $\mu$ of $D$, and $g(\alpha)=O(\sqrt{\log(1/\alpha)})$ in case (a) and $g(\alpha)=O_{k}(\alpha^{-1/k})$ in case (b). Moreover, these error bounds are optimal, up to constant factors. Specifically, the error bound of (a) cannot be asymptotically improved even if $D=N(\mu,I)$, as long as the list size is $\operatorname{poly}(1/\alpha)$. The error bound of (b) cannot be asymptotically improved as long as the list size is only a function of $\alpha$.

For the detailed statements, the reader is referred to Section 5.

We now turn to our computational lower bounds. Given Theorem 1.5, the following natural question arises: For the case of Gaussians, can we achieve the minimax bound in polynomial time? We provide evidence that this may not be possible, by proving a Statistical Query (SQ) lower bound for this problem. Recall that a Statistical Query (SQ) algorithm [Kea98] 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 is a restricted but broad class of algorithms, encompassing many algorithmic techniques in machine learning. A recent line of work [FGR+13, FPV15, FGV17, Fel17] developed a framework of proving unconditional lower bounds on the complexity of SQ algorithms for search problems over distributions.

By leveraging this framework, using the techniques of our previous work [DKS16b], we show that any SQ algorithm for list-decodable Gaussian mean estimation that guarantees error $\alpha^{-1/d}$, for some $d\geq 2$, requires either high accuracy queries or exponentially many queries:

Any SQ list-decodable mean estimation algorithm for $G\sim N(\mu,I)$ that returns a list of sub-exponential size so that some element in the list is within distance $O(\alpha^{-1/d})$ of the mean $\mu$ of $G$ requires either queries of accuracy $2^{O((1/\alpha)^{2/d})}\cdot n^{-\Omega(d)}$ or $2^{n^{\Omega(1)}}$ queries.

The reader is referred to Section 5.2 for the formal statement and proof.

5 Related Work

Robust Estimation. The field of robust statistics [Tuk60, Hub64, HR09, HRRS86, RL05] studies the design of estimators that are stable to model misspecification. After several decades of investigation, the statistics community has discovered a number of estimators that are provably robust in the sense that they can tolerate a constant (less than $1/2$) fraction of corruptions, independent of the dimension. While the information-theoretic aspects of robust estimation have been understood, the central algorithmic question — that of designing robust and computationally efficient estimators in high-dimensions — had remained open.

Recent work in computer science [DKK+16, LRV16] shed light to this question by providing the first efficient robust learning algorithms for a variety of high-dimensional distributions. Specifically, [DKK+16] gave the first robust learning algorithms that can tolerate a constant fraction of corruptions, independent of the dimension. Subsequently, there has been a flurry of research activity on algorithmic robust high-dimensional estimation. This includes robust estimation of graphical models [DKS16a], handling a large fraction of corruptions in the list-decodable model [CSV17, SCV17], developing robust algorithms under sparsity assumptions [BDLS17], obtaining optimal error guarantees [DKK+17b], establishing computational lower bounds for robust estimation [DKS16b], establishing connections with robust supervised learning [DKS17], and designing practical algorithms for data analysis applications [DKK+17a].

Learning GMMs. A long line of work initiated by Dasgupta [Das99], see, e.g., [AK01, VW02, AM05, KSV08, BV08], provides computationally efficient algorithms for recovering the parameters of a GMM under various separation assumptions between the mixture components. More recently, efficient parameter learning algorithms were obtained [MV10, BS10, HP15] under minimal information-theoretic separation assumptions. Without separation conditions, the sample complexity of parameter estimation is known to scale exponentially with the number of components, even in one dimension [MV10, HP15]. To circumvent this information-theoretic bottleneck of parameter learning, a related line of work has studied parameter learning in a smoothed setting [HK13, GVX14, BCMV14, ABG+14, GHK15]. The related problems of density estimation and proper learning for GMMs have also been extensively studied [FOS06, SOAJ14, MV10, HP15, ADLS17, LS17]. 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.

Most relevant to the current work is the classical work of Vempala and Wang [VW02] and the very recent work by Regev and Vijayraghavan [RV17]. Specifically, [VW02] gave an efficient algorithm that learns the parameters of spherical GMMs under the weakest separation conditions known to date. On the other hand, [RV17] characterize the separation conditions under which parameter learning for spherical GMMs can be solved with $\operatorname{poly}(n,k,1/\delta)$ samples. Whether such a separation can be achieved with an efficient algorithm was left open in [RV17]. Our work makes substantial progress in this direction.

6 Detailed Overview of Techniques

We start by reviewing the framework of [DKK+16] for robust mean estimation in the small error regime, followed by an explanation of the main difficulties that arise in the large error regime of the current paper.

In the small error regime, the “filtering” algorithm of [DKK+16] for robust Gaussian mean estimation works by iteratively detecting and removing outliers (corrupted samples) until the empirical variance in every direction is not much larger than expected. If every direction has small empirical variance, then the true mean and the empirical mean are close to each other [DKK+16]. Otherwise, the [DKK+16] algorithm projects the input points in a direction of maximum variance and throws away those points whose projections lie unexpectedly far from the empirical median in this direction. While this iterative spectral technique for outlier removal is by now well-understood for the small error regime (and has been applied to various settings), there are two major obstacles that arise if one wants to generalize it to the large error regime, i.e., where only a small fraction $\alpha$ of samples are good.

The first difficulty is that even the one-dimensional version of the problem in the large error regime is non-trivial. Specifically, consider a direction $v$ of large empirical variance. The [DKK+16] algorithm exploits the fact that the empirical median is a robust estimator of the mean in the one-dimensional setting. In contrast, in the large error regime, it is not clear how to approximate the true mean of a one-dimensional projection. This holds for the following reason: The input distribution can simulate a mixture of $1/\alpha$ many Gaussians whose means are far from each other, and the algorithm will have no way of knowing which is the real one. In order to get around this obstacle, we construct more elaborate outlier-removal algorithms, which we call multifilters. Roughly speaking, a multifilter can return several (potentially overlapping) subsets of the original dataset $T$ with the guarantee that at least one of these subsets is substantially “cleaner” than $T$.

The second difficulty is somewhat harder to deal with. As already mentioned, the filtering algorithm of [DKK+16] iteratively removes outliers by looking for directions in which the empirical distribution has a substantially larger variance than it should. In the low error regime, this approach does a good job of detecting are removing the corrupted points that can move the empirical mean far from the true mean. In the large error regime, the situation is substantially different. In particular, it is entirely possible that the empirical distribution does not have abnormally large variance in any direction, while still the empirical mean is $\Omega(\sqrt{1/\alpha})$-far from the true mean. That is, considering the variance of one-dimensional projections of our dataset in various directions seems inadequate in order to improve the $O(\sqrt{1/\alpha})$ error bound. This obstacle is inherent: the variance of linear polynomials (projections) is not a sufficiently accurate method of detecting a small fraction of good samples being substantially displaced from the mean of the bad samples. To circumvent this obstacle, we will use higher degree polynomials, which are much more sensitive to a small fraction of points being far away from the others. In particular, our algorithms will search for degree-$d$ polynomials that have abnormally large expectation or variance, and use such polynomials to construct our multifilters.

We now sketch how to exploit the existence of a large variance polynomial $p$ to construct a multifilter. Intuitively, the existence of such a polynomial $p$ suggests that there are many points that are far away from other points, and therefore separating these points into (potentially overlapping) clusters should guarantee that almost all good points are in the same cluster. Unfortunately, for this idea to work, we need to know that the variance of $p$ on the good set of points $S$ is not too large. For degree-$1$ polynomials $p$ this condition holds automatically. If $S$ is a sufficiently large set of samples from $G\sim N(\mu,I)$ and $p$ is a normalized linear form, then $\operatorname*{Var}[p(S)]\approx\operatorname*{Var}[p(G)]=1$. But if $p$ has degree at least $2$, the variance $\operatorname*{Var}[p(S)]$ depends on the true mean $\mu$, which unfortunately is unknown. Fortunately, there is a way to circumvent this obstacle by either producing a multifilter or verifying that the variance $\operatorname*{Var}[p(G)]$ is not too large.

We do this as follows: Firstly, we show that the variance $\operatorname*{Var}[p(G)]$, $G\sim N(\mu,I)$, can be expressed as an average of $p_{i}^{2}(\mu)$ for some explicitly computable, normalized, homogeneous polynomials $p_{i}$ (see Lemma 3.24). We then need to algorithmically verify that the polynomials $p_{i}(\mu)$ are not too large. This is difficult to do directly, so instead we replace each $p_{i}$ by the corresponding multilinear polynomial $q_{i}$, and note that $p_{i}(\mu)$ is the average value of $q_{i}$ at many independent copies of $G$. If this is large, then it means that evaluating $q_{i}$ at a random tuple of samples will often have larger than expected size.

This idea will allow us to produce a multifilter for the following reason: Since each $q_{i}$ is multilinear, this essentially allows us to write it as a composition of linear functions. More rigorously, we use the following iterative process: We iteratively plug-in variables one at a time to $q_{i}$. If at any step the size of the resulting polynomial jumps substantially, then the fact that this size is not well-concentrated as we try different samples will allow us to produce a multifilter. The details of this argument are given in Lemma 3.27 of Section 3.7.

6.2 Learning Spherical GMMs

Since a Gaussian mixture model can simultaneously be thought of as a mixture of any one of its components with some error distribution, applying our list-decoding algorithm to samples from a GMM will return a list of hypotheses so that every mean in the mixture is close to some hypothesis in the list. We can then use this list to cluster our samples by component.

In particular, given samples from a Gaussian $G=N(\mu,I)$ and many possible means $h_{1},\ldots,h_{m}$, we consider the process of associating a sample $x$ from $G$ with the nearest $h_{i}$. We note that $x$ is closer to $h_{j}$ than $h_{i}$ if and only if its projection onto the line between them is. Now if $h_{i}$ is substantially closer to $\mu$ than $h_{j}$ is, then this requires that this projection (which is Gaussian distributed) be far from its mean, which happens with tiny probability. Thus, by a union bound, as long as our list contains some $h_{i}$ that is close to $\mu$, the closest hypothesis to $x$ with high probability is not much further. If the separation between the means in our mixture is much larger than the separation between the means and the closest hypotheses, this implies that almost all samples are associated with one of the hypotheses near its component mean, and this will allow us to cluster samples by component. This idea of clustering points based on which of a finite set they are close to is an important idea that shows up in several related contexts in this paper.

The above idea works more or less as stated for mixtures of identity covariance Gaussians, but when dealing with more general mixtures of spherical Gaussians several complications arise. Firstly, in order to run out list-decoding algorithm, we need to know (a good approximation to) the covariance matrix of each component. The other difficulty is that, in order to cluster points, we will take a set of all nearby hypotheses that have reasonable numbers of samples associated with them. The issue is that we no longer know what “nearby” means, as it should depend on the covariance matrix of the associated Gaussian.

To solve the first of these problems we use a trick that will be reused several times. We note that two samples from the same Gaussian $N(\mu,\sigma^{2}I)$ have distance approximately $\Theta(\sigma\sqrt{n})$, and that even one sample from $N(\mu,\sigma^{2}I)$ is unlikely to be much closer than this to samples from different components. Therefore, by simply looking at the distance to the closest other sample gives us a constant factor approximation to the standard deviation of the corresponding component. This allows us to write down a polynomial-size list of viable hypothesis standard deviations. Running our list decoding algorithm for each standard deviation, gives us a polynomial-size list of hypothesis means.

To solve the second problem, we use the above idea to approximate the standard deviations associated to our sample points. When clustering them, we look for collections of sample points with standard deviations approximately the same $\sigma$, whose closest hypotheses are within some reasonable multiple of $\sigma$ of each other. Since we are able to approximate the size of the component that our samples are coming from, we can guarantee that we aren’t accidentally merging several smaller clusters together by using the wrong radius.

One slight wrinkle with the above sketched learning algorithm is that since the number of candidate hypotheses is polynomial in $n$, the separation between the components will be required to be at least $\sqrt{\log(n)}$. This bound is suboptimal, when $n$ is very large. Another issue is that the overall runtime of the learning algorithm would not be a fixed polynomial in $n$, but would scale as $n^{d}$. There is a way around both these issues, by reducing to a lower dimensional problem.

In particular, standard techniques involve looking at the $k$ largest principle values that allow one to project onto a subspace of dimension $k$ without losing too much. Unfortunately, these ideas require that all of the Gaussians involved have roughly the same covariance. Fortunately, if $n$ is large, our ability to approximate the covariance associated to a sample by looking at its distances to other samples becomes more accurate. Using a slightly modification of this idea, we can actually break our samples into subsets so that each subset is a mixture of Gaussians of approximately the save covariance. By projecting each of these in turn, we can reduce the original problem to a $\operatorname{poly}(k)$ number of dimensions and eliminate this extra term.

6.3 Minimax Error Bounds

We now explain our approach to pin down the information-theoretic optimal error for the list-decodable mean estimation problem. Concretely, for the identity covariance Gaussian case we show that there is an (inefficient) algorithm that guarantees that some hypothesis is within $O(\sqrt{\log(1/\alpha)})$ of the true mean. The basic idea is that the true mean must have the property that there is an $\alpha$-fraction of samples that are well-concentrated (in the sense of having good tail bounds in every direction) about the point. The goal of our (inefficient) algorithm will be to find a small number of balls of radius $O(\sqrt{\log(1/\alpha)})$ that covers the set of all such points. We show that such a set exists using the covering/packing duality. In particular, we note that if there are a large number of such sets with means far apart, we get a contradiction since the sets must be individually large but their overlaps must be pairwise small (due to concentration bounds).

Regarding lower bounds, [CSV17] showed an $\Omega(\sqrt{\log(1/\alpha)})$ error lower bound for $N(\mu,\Sigma)$, where $\Sigma$ is unknown and $\Sigma\preceq I$. We strengthen this result by showing that the $\Omega(\sqrt{\log(1/\alpha)})$ lower bound holds even for $N(\mu,I)$. We also prove matching lower bounds of $\Omega_{k}(\alpha^{-1/k})$ for distributions with bounded moments. Our proofs proceed by exhibiting distributions $X$, so that $X$ can be written as $X=\alpha X_{i}+(1-\alpha)E_{i}$ for many different $X_{i}$ satisfying the necessary hypotheses. Then any list-decoding algorithm must return a list of hypotheses close to the mean of every $X_{i}$. If there are many such $X_{i}$’s with means pairwise separated, then the list-decoding algorithm must either return many hypotheses or have large error.

6.4 SQ Lower Bounds

Finally, we prove lower bounds for list-decoding algorithms in the Statistical Query (SQ) model. Roughly speaking, we show that any SQ algorithm must either spend $n^{d}$ time or have accuracy higher than $\alpha^{-1/2d}$, suggesting that our list-decoding algorithm is qualitatively tight in its tradeoff between runtime and sample complexity.

We prove these bounds using the technology developed in [DKS16b]. This basically reduces to finding a $1$-dimensional distribution whose first many moments agree with the corresponding moments of a standard Gaussian. In our case, this amounts to constructing a one-dimensional distribution $A=\alpha N(\alpha^{-1/d},1)+(1-\alpha)E$, so that $A$’s first $d$ moments agree with those of a standard Gaussian. This can be done essentially because the $\alpha N(\alpha^{-1/d},1)$ part of the distribution only contributes at most a constant to any of the first $d$ moments. This allows us to take $E$ approximately Gaussian but slightly tweaked near in order to fix these first few moments.

We note however, that if we move the error component much further from , its contribution to the $d^{th}$ moment becomes super-constant and thus impossible to hide. This corresponds to the fact that degree-$d$ moments are sufficient (and necessary) in order to detect errors of size $\alpha^{-1/d}$.

7 Organization

The structure of the paper is as follows: In Section 2, we provide the necessary definitions and technical facts. In Section 3, we present our list-decoding algorithm. Section 4 gives our algorithm for GMMs. Finally, our minimax error and SQ lower bounds are given in Section 5.

Definitions and Preliminaries

Our basic objects of study are the Gaussian distribution and finite mixtures of spherical Gaussians:

2 Formal Problem Definitions

We record here the formal definitions of the problems that we study. Our first problem is robust mean estimation in the list-decodable learning model. We start by defining the list-decodable model:

We say that a learning problem is $(m,\varepsilon)$-list decodably solvable if there exists an efficient algorithm that can output a set of $m$ hypotheses with the guarantee that at least one is accurate to within error $\varepsilon$ with high probability.

Our notion of robust estimation relies on the following model of corruptions:

Given $0<\alpha\leq 1$ and a distribution family $\mathcal{D}$, an $\alpha$-corrupted set of samples $T$ of size $m$ is generated as follows: First, a set $S$ of $\alpha\cdot m$ many samples are drawn independently from some unknown $D\in\mathcal{D}$. Then an omniscient adversary, that is allowed to inspect the set $S$, adds an arbitrary set of $(1-\alpha)\cdot m$ many points to the set $S$ to obtain the set $T$.

We are now ready to define the problem of list-decodable robust mean estimation:

Our main algorithmic result is for the important special case that $\mathcal{D}$ is the family of unknown mean known covariance Gaussian distributions. We also establish minimax bounds that apply for more general distribution families.

Our second problem is that of learning mixtures of separated spherical Gaussians:

Given a positive integer $k$ and samples from a spherical $k$-GMM $F(x)=\sum_{i=1}^{k}w_{i}N(\mu_{i},\sigma_{i}^{2}\cdot I)$, we want to estimate the parameters $\{(w_{i},\mu_{i},\sigma_{i}),i\in[k]\}$ up to a required accuracy $\delta$. More specifically, we would like to return a list $\{(u_{i},\nu_{i},s_{i}),i\in[k]\}$ so that with high probability the following holds: For some permutation $\pi\in\mathbf{S}_{k}$ we have that for all $i\in[k]$: $|w_{i}-u_{\pi(i)}|\leq\delta$, $\|\mu_{i}-\nu_{\pi(i)}\|_{2}/\sigma_{i}\leq\delta/w_{i}$, and $|\sigma_{i}-s_{\pi(i)}|/\sigma_{i}\leq(\delta/w_{i})/\sqrt{n}$.

3 Basics of Hermite Analysis and Concentration

We will need the following standard concentration bound for degree-$d$ polynomials over independent Gaussians (see, e.g., [Jan97]):

List-Decodable Robust Mean Estimation Algorithm

In this section, we prove our main algorithmic result on list-decodable mean estimation:

The key idea procedure behind our algorithm is a subroutine that given a set of samples either cleans it up producing one or two subsets at least one of which has substantially fewer errors than the original, or certifies that the mean of $G$ must be close to the empirical mean (Proposition 3.6). Using this subroutine, our final algorithm can be obtained by repeatedly applying the subroutine recursively to the returned sets until they produce vectors. The details of this analysis are in Section 3.3.

Before we can get into the detailed overview of this proof, it is necessary to lay out some technical groundwork. First, we will want to have a deterministic condition under which our algorithm will succeed. To that end, we introduce two important definitions. We say that a set $S$ is representative of $G$ if it behaves like a set of independent samples of $G$, in particular in the sense that it is a PRG against low-degree polynomial threshold functions for $G$. We also say that a larger set $T$ is good if (roughly speaking) an $\alpha$-fraction of the elements of $T$ are a representative set for $G$. For technical reasons, will will also want the points of $T$ to be not too far apart from each other. In Section 3.1, we discuss the definitions of representative and good sets and provide some basic results.

In Section 3.2, we show that given a large set of points that contain an $\alpha$-fraction of good points from, one can algorithmically find $O(1/\alpha)$ many subsets so that with high probability at least one of them is good (and thus can be fed into the rest of our algorithm). This would be immediate if t were it not for the requirement that the points in a good set be not too far apart. As it stands, this will require that we perform some very basic clustering algorithms.

The actual design of our multifilter involves working with several types of “pure” degree-$d$ polynomials and their appropriate tensors. In particular, we need to pay attention to harmonic polynomials (which behave well with respect to $L^{2}$-norms), homogeneous polynomials, and multilinear polynomials. In Section 3.5, we introduce these and give several algebraic results relating them that will be required later.

The multifilter at its base level requires a routine that given a polynomial $p$, where $p(T)$ behaves very differently from $p(G)$, allows us to use the values of $p(x)$ to separate the points coming from $G$ from the errors. The basic idea of the technique is to cluster the points $p(x)$, for $x\in T$, and throw away points that are too far from any cluster large enough to correspond to the bulk of the values of $p(G)$ (which must be well-concentrated), or to divide $T$ into two subsets with enough overlap to guarantee that any such cluster could be entirely contained on one side or the other. The details of this basic multifilter algorithm are covered in Section 3.4.

Unfortunately, the application of the multifilter in the application above has a slight catch to it. Our basic multifilter will only apply to $p$ if we can verify that $\operatorname*{Var}[p(G)]$ is not too large. This would be easy to verify if we knew the mean of $G$, but unfortunately, we do not and errors in our approximation may lead to $\operatorname*{Var}[p(G)]$ being much larger than anticipated, and in fact, potentially too large to apply our filter productively. In order to correct this, we will need new techniques to either prove that $\operatorname*{Var}[p(G)]$ is small or to find a filter in the process. Using analytic techniques, in Section 3.5 we show that the $\operatorname*{Var}[p(G)]$ is a weighted average of squares of $q_{i}(\mu-\mu_{T})$ for some normalized, homogeneous polynomials $q_{i}$. Thus, it suffices to verify that each $q_{i}(\mu-\mu_{T})$ is small.

To deal with this issue, it is actually much easier to work with multilinear polynomials, and so instead we deal with multilinear polynomials $r_{i}$ so that $r_{i}(x,\ldots,,x)=q_{i}(x)$. We thus need to verify that $r_{i}(\mu-\mu_{T},\mu-\mu_{T},\ldots,\mu-\mu_{T})$ is small. The discussion of the reduction to this problem is in Section 3.8, while the techniques for verifying that the $r_{i}$ have small mean is in Section 3.7.

In summary, the structure of this section is as follows: In Section 3.1, we define the important notions of a representative set and a good set and prove some basic properties. In Section 3.2, we do some basic clustering and show that a good set can be extracted from a set of corrupted samples. In Section 3.3, we present the main subroutine and show how it can be used to produce our final algorithm.

The remainder of Section 3 will be spent building this subroutine. In Section 3.4, we produce our most basic tool for creating multifilters given a single polynomial where $p(T)$ and $p(G)$ behave substantially differently. In Section 3.5, we present some basic background on harmonic, homogeneous and multilinear polynomials and their associated tensors. In Section 3.6, we use these to produce routines to find multifilters given degree-$1$ polynomials or degree-$2$ polynomials with bounded trace norm for which $p(\mu)$ is too large. In Section 3.7, we leverage these results to extend to produce a similar multifilter for arbitrary multilinear polynomials. In Section 3.8, we use this to get a multifilter for degree-$d$ harmonic polynomials whose $L^{2}$ norms are substantially larger than expected, and in Section 3.9, we combine this with spectral techniques to get the full version of our filtering procedure, thus finishing our algorithm.

1 Representative Sets and Good Sets

We note that even though the definition of “representativeness” of a set $S$ depends on the parameters $\alpha$ and $d$, these quantities will be fixed for the representative set $S$ throughout the execution of our algorithm, and thus the dependence will be implicit. Note that as $\alpha$ increases or $d$ decreases, the representativeness condition becomes weaker. Thus, $S$ being representative with parameters $\alpha$ and $d$ implies that it is also representative with parameter $\beta$ and $d^{\prime}$, for any $\beta\geq\alpha$ and $d^{\prime}\leq d$.

We start by showing that a sufficiently large set of samples drawn from $G\sim N(\mu,I)$ is representative with high probability. This fact follows from standard arguments using the VC-inequality:

for all such polynomials $p$ with probability at least $1-\tau$, if we take

Our list-decodable learning algorithm and its analysis require an appropriate notion of goodness for the corrupted set $T$. At a high-level, throughout its execution, our algorithm will produce several different subsets of the original corrupted set $T$ it starts with, in its attempt to remove outliers. Intuitively, we want a good set $T^{\prime}\subseteq T$ to have the property that a $\beta$-fraction of the points in $T^{\prime}$, with $\beta\geq\alpha$, come from a representative set $S$. However, we also need to account for the possibility that — in the process of removing outliers — our algorithm may also remove a small number of the points in the original representative set $S$. Moreover, for technical reasons, we will require that all points in a good set are contained in a ball of radius $O(\sqrt{n})$. This is formalized in the following definition:

All points in $T$ are within distance $O(\sqrt{n})$ of each other.

2 Naive Clustering

We note that the definition of “goodness” of a set $T$ depends on the parameters $\alpha$ and $d$. The parameter $\alpha$ will change multiple times during the execution of the algorithm (it will increase), while the parameter $d$ does not increase. This justifies our choice of making $\alpha$ explicit in the definition of “$\alpha$-good”, while keeping the dependence on $d$ implicit.

The additional constraint that all points in a good set are contained in a not-too-large ball means that our original set of corrupted samples $T$ may not form a good set. To rectify this issue, we start by performing a very basic clustering step as follows: Since an $\alpha$-fraction of the points in $T$ are concentrated within distance $O(\sqrt{n})$ from the true mean $\mu$, by taking a maximal set of non-overlapping, not-too-large balls each with a large fraction of points, we are guaranteed that at least one of them will contain a good set. This is formalized in the following lemma:

Let $S$ be the subset of $T$ containing the good samples, i.e., the points that are independent samples from $G$. We have that $S\subset T$ and $|S|={2}\alpha\cdot|T|$. By Lemma 3.3, the set $S$ is representative with probability at least $1-\tau$. We henceforth condition on this event.

If all points in $T$ are contained in a ball of radius $O(\sqrt{n})$, there is nothing to prove. If this is not the case, we will show how to efficiently find a collection of at most $1/\alpha$ many balls of radius $O(\sqrt{n})$, so that at least one of the balls contains at least a $(1-\alpha/{6}-e^{-{3}n})$-fraction of the points in $S$.

First note that, by the degree-$2$ Chernoff bound, we have that $\Pr[\|G-\mu\|_{2}^{2}\geq C^{2}n]\leq\exp(-{3}n)$, for a sufficiently large universal constant $C$. Since $\|G-\mu\|_{2}^{2}$ is a degree-$2$ polynomial and $S$ is representative, Definition 3.2 implies that at least an $(1-\exp(-{3}n)-\alpha/{6})$ fraction of points in $S$ are within distance $C\sqrt{n}$ of $\mu$.

Our clustering scheme works as follows: We consider a maximal set of disjoint balls of radius $R=2C\cdot\sqrt{n}$ centered at points of $T$, such that each ball contains at least an $\alpha$-fraction of the points in $T$. Note that this set is non-empty: Since $S$ is representative, the ball of radius $R$ centered at any point $x\in S$ with $\|x-\mu\|_{2}\leq C\sqrt{n}$ contains at least a $(1-e^{-3n}-\alpha/6)$-fraction of the points in $S$, and therefore at least an $\alpha$-fraction of the points in $T$.

Let $B_{1},\ldots,B_{s}$ be the maximal set of disjoint balls described above. Since each $B_{i}$ contains an $\alpha$-fraction of the points in $T$, there are at most $s\leq 1/\alpha$ many such balls. Let $B^{\prime}_{i}$, $i\in[s]$, be the ball with the same center as $B_{i}$ and radius $3R=6C\cdot\sqrt{n}$. Consider the subsets $T_{i}=T\cap B^{\prime}_{i}$, $i\in[s]$. We claim that at least one of the $T_{i}$’s is $\alpha$-good.

The pseudo-code for this clustering algorithm is given below.

We now prove correctness. By definition, all the points of $T_{i}$ are contained in a ball of radius $3R=6C\sqrt{n}$. Let $x$ be any point of $S$ within distance $C\sqrt{n}$ of the mean $\mu$ of $G$. (Recall that, since $S$ is representative, most of the points of $S$ satisfy this property.) Since the initial set of balls $B_{1},\ldots,B_{s}$ was constructed to be maximal, at least one $B_{i}$ must intersect the ball $B_{x}$ of radius $R=2C\sqrt{n}$ centered at $x$. This implies that the ball $B_{x}$ is contained in $B^{\prime}_{i}$. Therefore, $T_{i}$ contains all of the points of $S$ within $C\sqrt{n}$ of the mean $\mu$ of $G$. This completes the proof of Lemma 3.5. ∎

3 Main Multifilter Algorithm and Proof of Theorem 3.1

In this section, we describe the main subroutine that forms the core of our final algorithm. Let $0<\alpha<1/2$ be the fraction of good samples, and let $S$ be the representative set contained in the initial corrupted set $T$. Ideally, starting from the corrupted set $T$, we would like to either certify that the true mean $\mu$ of $G$ is close to the empirical mean $\mu_{T}$ of $T$, or we have an efficient method to “clean up” the set $T$, i.e., obtain a set $T^{\prime}\subset T$ containing a higher fraction of clean samples.

Unfortunately, such a goal seems unattainable in the setting where the fraction $\alpha$ of clean samples is smaller than $1/2$ for the following reason: Suppose that $\alpha=1/k$ for some integer $k\geq 2$. It is quite possible that the original set of corrupted points $T$ is a set of samples drawn from a mixture $\sum_{i=1}^{k}(1/k)\cdot N(\mu_{i},I)$ of $k$ spherical Gaussians with separated means $\mu_{i}$. In this case, the best that we can hope to achieve is break $T$ into several (potentially overlapping) subsets $T_{i}$, with the guarantee that at least one of the subsets $T_{i}$ has a higher fraction of points from $S$. As our final algorithm will then have to recurse on the $T_{i}$’s, we will need to ensure that they are not too large in order to obtain a bound on the total runtime of the algorithm.

Our main sub-routine satisfies the guarantees of the following proposition:

Moreover, with probability at least $1-\tau$, the output of the algorithm satisfies the following guarantees:

In case (1), if $T$ is $\alpha$-good with respect to $G\sim N(\mu,I)$, then we have that

In case (2), the set $T$ is not $\alpha$-good.

In case (3), if $T$ is $\alpha$-good then at least one $T_{i}$ is $\alpha_{i}$-good with respect to $G$.

The algorithm runs in time $O\left((d+\ln(1/\tau))/\alpha\right)^{d}\cdot\operatorname{poly}(|T|,n^{d})$.

The proof of Proposition 3.6 requires several ingredients and is given in the following subsections. We start by showing how Theorem 3.1 follows from this proposition.

The overall algorithm proceeds as follows: Using Lemma 3.5, we produce a list $\{T_{1},\ldots,T_{s}\}$ of subsets of $T$ at least one of which is $\alpha/2$-good with respect to $G$.

The algorithm maintains a list of pairs $\{(T_{i},\alpha_{i})\}_{i=1}^{m}$ with the guarantee that with high probability at least one of the $T_{i}$ is $\alpha_{i}$-good with respect to $G$. It then iteratively applies the routine from Proposition 3.6, to each pair $(T_{i},\alpha_{i})$ in the list, and when its output is not a vector, it replaces the pair $(T_{i},\alpha_{i})$ by either zero, one, or two such pairs in the list. Then it outputs the list of all vectors that Proposition 3.6 ever returned.

We will show that this process produces a list of $O(1/\alpha^{3})$ many candidate hypotheses. We provide an efficient black-box way in which we can reduce the size of this list to $O(1/\alpha)$ without losing more than a constant factor in the final error guarantee. See Proposition B.1.

The pseudo-code for the algorithm is given below:

We now proceed with the analysis. First, it is easy to see that $\sum_{i}\alpha_{i}^{-2}$ can only decrease. Therefore, since all $\alpha_{i}\leq 1$, we can never have more than $O(\alpha^{-{3}})$ many terms in the list. Each iteration of the algorithm either increases the number of $T_{i}$’s or decreases the size of some $T_{i}$ by at least $1$. Therefore, the algorithm terminates after at most $O(|T|/\alpha^{{3}})$ many calls to MainSubroutine.

We claim that with probability at least $1-\tau$, at every step, we either have that one of the $T_{i}$’s is $\alpha_{i}$-good, or we have already produced a good approximation to $\mu$. We can view the set of all $T_{i}$’s ever produced by the algorithm as a tree, where the children of each $T_{i}$ are all subsets produced by MainSubroutine on $T_{i}$. We say that $T_{i}$ belongs to the $j$-th generation if the path from $T$ to $T_{i}$ in this tree has length $j$. Then, it is easy to see inductively that either at least one $T_{i}$ in the $j$-th generation is $\alpha_{i}$-good or a good approximation to the mean was output in an earlier generation, with probability at least $j\tau/|T|$. This is because given a $T_{i}$ that is $\alpha_{i}$-good in the $j^{th}$ generation, with probability at least $1-\tau/|T|$ it either produces a good approximation to the mean or has at least one descendant $T^{\prime}$ that is $\alpha^{\prime}$-good. Since there can only be $|T|$ generations, the list of vectors produced contains a good approximation to the mean with probability at least $1-\tau$. This completes the proof.

To analyze the runtime, note that there are at most $|T|$ generations each involving running MainSubroutine on $O(\alpha^{-3})$ many $T_{i}$’s. This clearly dominates the runtime from the initial use of NaiveClustering, the final use of ListReduction and the other bookkeeping steps. ∎

In the remaining sections, we describe a number of subroutines that either verify some desired properties of low-degree polynomials or produce an output list of subsets satisfying the following condition:

We say that a list of pairs $\{(T_{i},\alpha_{i})\}$, where $T_{i}\subset T$ and $\alpha_{i}\in(0,1)$, satisfies the multifilter condition for $(T,\alpha)$ if the following hold:

$\sum_{i}1/\alpha_{i}^{2}\leq 1/\alpha^{2}$, and

If $T$ is $\alpha$-good (with respect to $G$), then at least one of the $T_{i}$’s is $\alpha_{i}$-good (with respect to $G$).

4 Basic Multifilter Routine

How do we produce the multifilters for corrupted sets $T$ required in Proposition 3.6? The basic idea is the following: If we find a degree-$d$ $n$-variable polynomial $p$ such that the mean and variance of $p(T)$ behave substantially differently from the mean and variance of $p(G)$, we will be able to split $T$ into (overlapping) subsets — based on the values of $p(x)$, $x\in T$ — with the property that at least one of these subsets is “cleaner” than $T$.

The obvious difficulty is that the mean $\mu$ of $G$ is unknown, and therefore so is the mean and variance of $p(G)$. Suppose for now that we have a low-degree polynomial $p$ whose variance $\operatorname*{Var}[p(G)]$ is small. (As we argue in the following sections, we can guarantee that this condition is satisfied with high probability for the polynomials $p$ that we will use.) By concentration, it then follows that at least a $(1-\alpha/2)$ fraction of the points of $\{p(x),x\in S\}$ will lie in an interval of length $\log(1/\alpha)^{O(d)}$. If a decent number of points of $\{p(x),x\in T\}$ fall outside an interval of this length, then we know that the distribution of $p(T)$ is very different. In such a case, be removing points of $T$ beyond some threshold, we can get rid of substantially more bad points than good points.

Recall that we assumed that $\operatorname*{Var}[p(G)]$ is small. However, the mean of $p(G)$ is still unknown. This creates the following complication: There might be several intervals of length $\log(1/\alpha)^{O(d)}$ each of them containing an $\alpha$-fraction of $T$. Any of these intervals would make a reasonable hypothesis and we have no way of distinguishing between them. In such a setting, we split our set $T$ into several overlapping subsets. We make the overlap between them large enough so that no matter where the (unknown) true mean of $p(G)$ is, the vast majority of the points of $\{p(x),x\in S\}$ will lie in the same subset. At the same time, because there is a substantial number of points in $T$ on either side of the divide, we can guarantee that not too many points are kept on either side.

Formally, we establish the following proposition:

“YES”. If the algorithm returns “YES”, then we have that:

$\operatorname*{Var}[p(T)]\leq O(d+\log(1/\alpha))^{d}(\log(2+\log(1/\alpha)))^{2}$, and

“NO”. If the algorithm returns “NO”, then the set $T$ is not $\alpha$-good.

A list of pairs $\{(T_{i},\alpha_{i})\}$, $T_{i}\subset T$, of length one or two satisfying the multifilter condition for $(T,\alpha)$.

The rest of this subsection is devoted to the proof of Proposition 3.8.

We start by describing the basic multifilter routine in pseudocode:

The algorithm begins by defining some appropriate constants and in particular choosing the parameter $R$ to be polylogarithmic in $1/\alpha$, so that if $\operatorname*{Var}[p(G)]\leq 1$ then all but a $\operatorname{poly}(\alpha)$-fraction of the points of $p(G)$ must lie in an interval of length $R$.

We begin in Step 3 by doing a basic sanity check, ensuring that the total range of values of $p(T)$ is not too large. On the one hand, we are guaranteed by Claim 3.9 in the next section that this will be true if $\operatorname*{Var}[p(G)]\leq 1$ and $T$ is good. This assumption will become technically useful as it will imply that if $p(T)$ has too large a mean or variance, then this cannot be due to a negligible number of points at extreme distances from the mean, and instead must be due to a relatively significant number of errors.

We split our analysis into cases roughly corresponding to the one cluster versus two clusters cases. Step 4 corresponds to the case where there is a single interval of length not much more than $R$ that contains all but an $\alpha/2$-fraction of the points of $\{p(x),x\in T\}$. This implies that since an $\alpha$-fraction of points come from $p(G)$, the mean of $p(G)$ cannot be far outside of this interval. Thus, we can construct a slightly larger interval that contains all but a tiny fraction of the points of $p(G)$. From there, it is a relatively easy to show that since $p(T)$ has substantially different behavior from $p(G)$, there must be some (potentially larger interval), where the fraction of points of $p(T)$ outside of this interval is much larger than the fraction of points of $p(G)$ outside of this interval. Thus, throwing away the points outside of the latter interval will serve to clean up our dataset.

Step 5 covers the other case where no short interval contains all but an $\alpha/2$-fraction of the points. In this case, we try to partition $T$ into sets based on whether $p(x)$ lies in $(\infty,t+R)$ or $(t-R,\infty)$, for some appropriate threshold $t$. Notice that since these intervals overlap on an interval of length $2R$, no matter where the mean of $p(G)$ lies, almost all of the points of $p(G)$ will lie in one of the two intervals (if it’s in the middle it could even be the case that almost all points lie in both). We now need to pick a threshold $t$ so that the number of elements of $p(T)$ not on either side is substantial. We are aided in this endeavor by the fact that there must be some gap of length $\Omega(R\log(2+\log(1/\alpha)))$ with at least an $\Omega(\alpha)$-fraction of the points lying on either side of it. By trying many possible values if $t$ within this interval, it is not hard to show that there must be some threshold where neither side keeps too many.

4.2 First Stage of Correctness Proof

In this subsection, we prove the correctness conditions of our algorithm that do not establish goodness of the output.

We start by showing that if the algorithm returns “NO” in Step 3, then the set $T$ is not $\alpha$-good:

Suppose that $T$ is $\alpha$-good and $\operatorname*{Var}[p(G)]\leq 1$. Then, it holds that $\max_{x,y\in T}|p(x)-p(y)|\leq O(1+n/d)^{d/2}$.

Taking the product over $i$ from $1$ to $n$ yields:

If $u\geq 0$, we have that the left hand side is at least

Letting $u=\min(1/2,d/n)$, we find that if $d\ll n$, then the above sum is at most $(n/d)^{d}\exp(O(d^{2}/n)+O(d\|y\|_{2}^{2}/n))=O(n/d)^{d}$. If $d\gg n$, this is $O(1)^{d+n+\|y\|_{2}^{2}}$. In either case, it is $O(1+n/d)^{d}$, and thus $|p(y)|$ is at most $O(1+n/d)^{d/2}$. This completes the proof of Claim 3.9. ∎

The correctness of the other step in which the algorithm can output “NO”, Step 4(b)iii, is deferred to after Lemma 3.18 below.

We next show that the algorithm finds an appropriate threshold $t$ in Steps 4(b)i and 5a.

If the algorithm reaches Step 4(b)i, there exists a threshold $t>2R$ satisfying (2).

The main idea of the proof is quite simple: If there is no $t>2R$ satisfying (2), the random variable $p(T)$ is well-concentrated, which contradicts the lower bound on its variance that holds in Step 4b. Suppose, for the sake of contradiction, that for all $t>2R$, it holds:

By a change of variable, this is equivalent to the following statement. For all $t>2R+(b-a)/2$, it holds:

The proof will require some basic properties of the Gamma function. We will use the following notation: Let $\Gamma(s,x)=\int_{x}^{\infty}\exp(-t)t^{s-1}dt$ be the incomplete gamma function and $\Gamma(s)=\Gamma(s,0)$ be the Gamma function. We will need the following technical claim about the incomplete Gamma function, whose proof can be found in Appendix A:

For $s\geq 1,x\geq 0$, we have that $\Gamma(s,x)\leq\exp(-x)(x+s)^{s-1}$.

We can now bound from below the variance of $p(T)$ as follows:

By the construction of the algorithm, when the variance $\operatorname*{Var}[p(T)]$ is this small, the algorithm does not reach Step 4(b)i. This gives the desired contradiction and completes the proof of Lemma 3.10. ∎

We next show that the algorithm finds an appropriate threshold in Step 5a:

If the algorithm reaches Step 5a, there exists a threshold $t$ satisfying (3).

We will show that if no $t$ satisfies (3), then at least a $1-\alpha/2$ fraction of points $x\in T$ have $p(x)$ in an interval of length $O(\log(2+\log(1/\alpha))\cdot R)$. This would be a contradiction, as no such interval exists if the algorithm reaches Step 5. This is essentially because if there is no such $t$, then if we consider the fraction of $x\in T$ with $p(x)\leq t$, this must grow very quickly with $t$ after it first exceeds $\alpha/4$, and thus there is very little distance between the $\alpha/4$-tails of the distribution and the median. This means that the $\alpha/4$-tails on the left and on the right must be close to each other, implying that the interval from Step 4 must exist.

Let $t_{m}=2m\cdot R+\textrm{median}(p(T))$, for integer $m$. If no $t$ satisfies (3), then $t_{m}$ is not a suitable choice of $t$, and therefore

where $T_{1,m}=\{x\in T:p(x){>}t_{m}-R\}$ and $T_{2,m}=\{x\in T:p(x){<}t_{m}+R\}$ for all $m$ with $|T_{i,m}|\geq\alpha|T|$.

Let $a_{m}=|T_{1,m}|/|T|$. Notice that $|T_{2,m}|/|T|=(1-a_{m+1})$. If Equation (4) holds for a given $m$, then we must have $a_{m}^{2}+(1-a_{m+1})^{2}>(1-\alpha^{2}/100)^{2}.$ This means that $a_{m+1}<(a_{m}^{2}+\alpha^{2}/16)$. If $a_{m}\geq\alpha/4$, then $a_{m+1}\leq 2a_{m}^{2}$. Since $a_{1}<1/2$, it is easy to see that for some $m=\log_{2}\log(1/\alpha)+O(1)$ that $a_{m}\leq\alpha/4$. Thus, $|\{x\in T:p(x)>t+(2m+1)R\}|\leq|T|\alpha/4$. Similarly, we can see that $|\{x\in T:p(x)<t-(2m+1)R\}|\leq|T|\alpha/4$. This gives us an interval of length $O(R\log(2+\log(1/\alpha)))$ that contains all but an $\alpha/2$-fraction of the points of $p(T)$.

This provides the desired contradiction and gives the proof of Lemma 3.12. ∎

We now show that if the algorithm outputs a list of length one, the parameter $\alpha$ improves:

If the algorithm outputs a single pair $(T^{\prime},\alpha^{\prime})$, then $\alpha^{\prime}\geq\alpha$.

Note that $\alpha^{\prime}-\alpha>0$ when $|T^{\prime}|<|T|$. Moreover, $|T^{\prime}|<|T|$ holds true because the algorithm is guaranteed to remove at least one point, since

is non-zero. This completes the proof of Lemma 3.13. ∎

The next lemma shows that whenever the algorithm outputs two sets $T_{1},T_{2}$, the associated parameters $\alpha_{1},\alpha_{2}$ satisfy $1/\alpha_{1}^{2}+1/\alpha_{2}^{2}\leq 1/\alpha^{2}$. This condition will be crucial for bounding the runtime of the overall algorithm. We also note that this condition holds even if the set $T$ was not $\alpha$-good.

If the algorithm outputs a list of two pairs $\{(T_{1},\alpha_{1}),(T_{2},\alpha_{2})\}$, then we have that $1/\alpha_{1}^{2}+1/\alpha_{2}^{2}\leq 1/\alpha^{2}$.

4.3 Second Stage of Correctness Proof

In this subsection, we establish the correctness conditions of the algorithm that hold under the assumption that the set $T$ is $\alpha$-good. Specifically, we will show that if $T$ is $\alpha$-good the following hold: (1) If the algorithm returns “YES”, the means of $p(G)$ and $p(T)$ are close to each other. (2) If the algorithm returns a subset $T^{\prime}$, then $T^{\prime}$ is good. (3) If the algorithm returns two subsets $T_{1},T_{2}$, then at least one of the returned sets is good.

Since $T$ is $\alpha$-good, we have that $\operatorname*{Var}[p(G)]\leq 1$. By the degree-$d$ Chernoff bound and the definition of $R$ we get

We can now apply the bounds given by the definition of $T$ being $\alpha$-good, which give the claim. ∎

We next show that if we find an appropriate interval $[a,b]$ in Step 4, then the mean of $p(G)$ is contained in the interval $[a-R,b+R]$.

By construction, if the algorithm returns “YES”, then

and $[a,b]$ is an interval of length at most $O((C\log(1/\alpha))^{d/2}\log(2+\log(1/\alpha)))$ that contains all except an $\alpha/2$-fraction of points of $p(x)$, $x\in T$. We first consider the contribution of the points in $[a,b]$ to the variance to obtain that

We now prove (2). This proof requires a careful accounting of the number of points from $S$ and points from $T$ that are removed by our algorithm.

If $T$ is $\alpha$-good and the algorithm returns a single pair $(T^{\prime},\alpha^{\prime})$, then $T^{\prime}$ is $\alpha^{\prime}$-good.

Since $T$ is $\alpha$-good, it follows that $S\cap T$ contains a $1-\alpha/6+\exp(-3n)\geq 1/2$ fraction of the points in $S$, hence

But note that the probability of the above event is at least $8/\alpha$ times bigger for $X\in_{u}T$. Since $T^{\prime}$ has the filtered samples removed, we have that

and so $|S\cap T^{\prime}|/|S\cap T|\geq 1-(\alpha/8)(1-|T^{\prime}|/|T|)$. Now we can bound from below the fraction of good samples in $T^{\prime}$:

It remains to show that $|S\cap T^{\prime}|/|S|\geq 1-\alpha^{\prime}/6-\exp(-3n)$. Noting that $(1-\alpha/8(1-|T^{\prime}|/|T|))\alpha=\alpha^{\prime}|T^{\prime}|/|T|$, we get that:

Note that this proof shows that if $T$ is $\alpha$-good, then $T^{\prime}$ is $\alpha^{\prime}$-good. If $\alpha^{\prime}>1$, $T$ cannot be $\alpha^{\prime}$-good, as this would require that $|S\cap T^{\prime}|>|T^{\prime}|$. Thus, if the algorithm outputs “NO” in Step 4(b)iii, then $\alpha^{\prime}>1$ implies that $T$ cannot be $\alpha$-good.

Our last lemma establishes (3), handling the case when the algorithm outputs a list of two sets:

If $T$ is $\alpha$-good and the algorithm outputs a list $\{(T_{1},\alpha_{1}),(T_{2},\alpha_{2})\}$, then for some $i\in\{1,2\}$ $T_{i}$ is $\alpha_{i}$-good.

Note that $|S\cap T|\geq|S|(1-\alpha/6-e^{-3n})$. We have that $|S\cap(T\backslash T^{\prime})|\leq|S|\alpha^{2}/100$. Therefore, as long as $\alpha_{i}/6\geq\alpha/6+\alpha^{2}/100$, this set $T_{i}$ contains enough points of $S$. These points consist of at least $|T|(\alpha-\alpha^{2}/100)/|T_{i}|$ fraction of $T_{i}$. Thus, when

the set $T_{i}$ is $\alpha_{i}$-good. It remains to show that $\alpha_{i}=|T|\alpha(1-\alpha/100)/|T_{i}|\geq\alpha(1+3\alpha/50)$.

Recall that $|T|-\max\{|T_{1}|,|T_{2}|\}\geq\alpha\cdot|T|/4$. Therefore, $|T|/|T_{i}|\geq(1-\alpha/4)^{-1}$, and thus $\alpha_{i}\geq\alpha(1-\alpha/100)/(1-\alpha/4)\geq\alpha(1+3\alpha/50)$. This completes our proof. ∎

In the following subsections, we appropriately use the basic multifilter algorithm to handle increasingly broad families of degree-$d$ polynomials.

5 Useful Results on Polynomials and Tensors

Our algorithm and its analysis will make essential use of degree-$d$ multivariate polynomials. Thus, we will need to establish here some basic facts and terminology. A basic theme throughout will be that “pure” degree-$d$ multivariate polynomials are in $1-1$ correspondence with symmetric order-$d$ tensors. On the other hand, we will require three different notions of what it means for a polynomial to be “pure” degree-$d$: homogeneous polynomials, symmetric multilinear polynomials, and harmonic polynomials (which are useful when considering norms with respect to the Gaussian distribution).

The structure of this section is as follows: In Section 3.5.1, we give the definitions of the various types of polynomials we will require and their associated tensors. In Section 3.5.2, we establish a key lemma (Lemma 3.24) relating the first two moments of such polynomials under an identity covariance Gaussian. This lemma is essential for our algorithm and its analysis.

Any degree-$d$ harmonic polynomial $p(x)$ satisfies $p(x)=h_{A}(x)$, where

We call a homogeneous, harmonic or symmetric multilinear polynomial normalized if the associated tensor $A$ has $\|A\|_{2}=1$.

5.2 Useful Lemma on Polynomials