Quantum Entropy Scoring for Fast Robust Mean Estimation and Improved Outlier Detection

Introduction

We study outlier-robust statistics in high dimensions, focusing on the question: can theoretically sound outlier robust algorithms have practical running times for large, high-dimensional data sets? We address two related problems: robust mean estimation, which is primarily theoretical, and an applied counterpart, outlier detection.

Robust mean estimation Our main theoretical contribution is the first nearly-linear time algorithm for robust mean estimation with nearly-optimal error. Here the goal is to estimate the mean $\mu\in\R^{d}$ of a $d$-dimensional distribution $D$ given $\epsilon$-corrupted samples $X_{1},\ldots,X_{n}$ – that is, i.i.d. samples, an unknown $\epsilon$-fraction of which have been maliciously corrupted. Under (for instance) the assumption that the covariance of $D$ is bounded by $\Id$, it has been long known to be possible in exponential time to estimate $\mu$ by $\hat{\mu}$ having $\|\mu-\hat{\mu}\|_{2}\leq O(\sqrt{\epsilon})$. In particular, this rate of error is independent of $d$.

Polynomial-time algorithms provably achieving such $d$-independent error became known only recently, starting with the works . Until our work, the running time of algorithms with provably $d$-independent error remained suboptimal by polynomial factors in $d$ or $\epsilon$: the fastest running time achieved before this work was $\widetilde{O}(\min(nd^{2},nd/\epsilon^{6}))$ . (Here $\widetilde{O}(\cdot)$ notation hides logarithmic factors in $n$ and $d$). While these running times represent a dramatic improvement over previous exponential-time algorithms, there are still many interesting regimes where the additional runtime overheads these algorithms incur render them impractically slow. We give the first algorithm for robust mean estimation with running time $\widetilde{O}(nd)$ which achieves error $\|\mu-\hat{\mu}\|_{2}\leq O(\sqrt{\epsilon})$. Note that this running time is nearly-linear in the input size $nd$. Similar to prior works, our algorithm has information-theoretically optimal sample complexity and nearly-optimal error rates in both the bounded-covariance and sub-Gaussian regimes.

We compare our method to baselines based on PCA and Euclidean distances, as well as more sophisticated algorithms from existing literature based on nearest-neighbor distances. Our algorithm has nearly-linear running time in theory, and simple implementations in practice incur minimal overhead beyond standard spectral methods, allowing us to run on $1024$-dimensional data with no special optimizations and $8192$-dimensional data with a fast approximate implementation. It can therefore be used in practice to complement existing approaches to outlier detection in exploratory data analysis.

For us, an outlier is an element of a data set which was generated according to a different process than the majority of the data. For instance, we may imagine that our samples $X_{1},\ldots,X_{n}$ were sampled i.i.d. from a distribution $(1-\epsilon)D+\epsilon N$ over $\R^{d}$, where $D$ is the distribution of inliers, $N$ is the distribution of outliers, and $\epsilon>0$ is a small number – that is, we imagine that a constant fraction of our data may be outliers.

For this discussion, we also informally imagine that $N$ is sufficiently distinct from $D$ that the set of outliers could be approximately identified by brute-force search over subsets of $(1-\epsilon)n$ samples, if given unlimited computational resources. Otherwise, outlier detection is not a meaningful problem, and robust mean estimation is easy (because the empirical mean will be a good estimator). Under these circumstances, what makes identifying outliers and estimating the mean in their presence difficult? Chiefly:

Outliers may not be identifiable in isolation. On its own, a typical outlier $X_{i}\sim N$ may look much like a typical inlier $X_{j}\sim D$. For instance, it could be $\|X_{i}\|_{2}\approx\|X_{j}\|_{2}$, and $X_{i},X_{j}$ may have similar distance to the nearest few neighboring samples, especially in high dimensions where samples are far apart.

Outliers still introduce bias, collectively Even if individual outliers look innocuous, the collective effect a modified $\epsilon$-fraction of samples $X_{i}$ can still substantially change the empirical distribution of $X_{1},\ldots,X_{n}$. As a result, even simple statistical tasks like estimating the mean or covariance of $D$ require sophisticated estimators: naively pruning individual outliers and then employing standard empirical estimators typically leads to far-suboptimal error rates. For example, an $\epsilon$-fraction of $X_{1},\ldots,X_{n}$ which are all slightly biased in a single direction may shift the empirical mean of $X_{1},\ldots,X_{n}$, but this bias will be difficult to detect by looking at small numbers of samples at once. This also demonstrates that successful outlier detection can require global geometric information about a high-dimensional dataset, such as whether or not a direction exists in which many (say, $\epsilon n$) samples are unusually biased.

Outliers may be inhomogeneous. Outliers need not exhibit unusual bias in only one direction, or all have the same norm, or lie in a single cluster. Rather, if a dataset exhibits several forms of corruption, there may be as many different-looking kinds of outliers. In the theoretical robust mean estimation setting, the adversary producing $\epsilon$-corrupted samples may corrupt $\epsilon n/10$ samples by biasing them in some direction, another $\epsilon n/10$ samples by unusually enlarging their norms, and so forth.

Since robust mean estimation involves a malicious adversary, all of the above phenomena must be addressed by our robust mean estimation algorithm. In the empirical section of this paper, we focus on designing an outlier detection method suited to situations where at least one of them occurs – in other cases, existing methods (such as those based on Euclidean norms or local neighborhoods of individual samples ) may be more appropriate.

2 QUE: Quantum Entropy Scoring

Recent innovations in robust mean estimation rely on the following crucial observation about $\epsilon$-corrupted samples $X_{1},\ldots,X_{n}$ from a distribution $D$ with covariance $\Sigma\preceq\Id$. Namely: any subset $S\subseteq\{X_{1},\ldots,X_{n}\}$ of samples which shift the empirical mean by distance more than $\sqrt{\epsilon}$ in some direction $v$ also introduce an eigenvalue of magnitude greater than $1$ to the empirical covariance.

In robust mean estimation, this leads to (amongst others) the filter algorithm of , one of the first to achieve dimension-independent error rates. Roughly speaking, the algorithm iterates the following until the empirical covariance $\overline{\Sigma}$ has small spectral norm: (1) compute the top eigenvector $v$ of the empirical covariance of $\overline{\Sigma}$, then (2) throw out samples $X_{i}$ whose projections $|\langle X_{i}-\overline{\mu},v\rangle|\gg 1$ is unusually large, where $\overline{\mu}$ is the empirical mean of the corrupted dataset. For outlier detection this suggests a natural scoring rule – let the outlier score $\tau_{i}$ of sample $X_{i}$ be proportional to $|\langle X_{i}-\overline{\mu},v\rangle|$.

The main drawback of these algorithms is that they do not adequately account for inhomogeneity of outliers. For the filter, this leads to a worst-case running time of $\widetilde{O}(nd^{2})$, because the filter operation (which can be implemented in $\widetilde{O}(nd)$ time) may have to be repeated as many as $d$ times if the adversary introduces outliers lying in $d$ orthogonal directions. The rule $\tau_{i}=|\langle X_{i}-\overline{\mu},v\rangle|$ may miss outliers causing a large eigenvalue of $\overline{\Sigma}$, but in a direction orthogonal to the top eigenvector $v$.

Our first observation is that any eigenvalue/eigenvector $\lambda,v$ – not just the top ones – of the empirical covariance with $\lambda\gg 1$ must be due to outliers. We therefore consider the intermediate goal of finding a distribution over directions $v\in\R^{d}$ containing information about as many outlier directions as possible. We formalize this as the following entropy-regularized convex program over $d\times d$ positive semidefinite matrices:

QUE scores are also appealing from a computational perspective: we show that a list of approximate QUE scores $\tau_{i}^{\prime}=(1\pm 0.01)\tau_{i}$ can be computed from $X_{1},\ldots,X_{n}$ in nearly-linear time, by appropriate use of Johnson-Lindenstrauss sketching and efficient computation of the matrix exponential by series expansion. This is crucial to both the nearly-linear running time of our algorithm for robust mean estimation and to the scalability of our outlier detection method.

In Section 1.4 we describe refinements of QUE scoring which fit it into the matrix multiplicative weights framework , leading to our nearly-linear time algorithms for robust mean estimation. We give two very similar algorithms, one for when the distribution of inliers is only assumed to have bounded covariance, and one when the inliers are assumed to be subgaussian. The resulting algorithms are conceptually similar to the following modification of the filter mentioned above: until $\|\overline{\Sigma}\|_{2}\leq O(1)$, compute QUE scores $\tau_{i}$, throw out data points $X_{i}$ with $\tau_{i}\gg 1$, and repeat. (To obtain provable guarantees, our final algorithms are somewhat more complex: in some iterations we use QUE scores based on certain reweightings of the data learned in previous iterations.)

In Section 1.5 we describe experiments validating the QUE scoring rule on both synthetic and real data sets. We show that it performs especially well by comparison to local-neighborhood methods and to scoring based on only the top eigenvector in data sets where the inliers are close to isotropic (or can be made so by applying data whitening procedures) and in which there are heterogeneous outliers.

3 Related work

Robust mean estimation: The study of robust statistics and in particular robust mean estimation began with major works by Anscombe, Huber, Tukey and others in the 1960s . The literature on polynomial-time algorithms for robust statistics has exploded in recent years, following works by Diakonikolas et al and Lai, Rao and Vempala giving the first polynomial-time algorithms for robust mean estimation with dimension-independent (or nearly dimension-indepedent) error . A full survey is beyond our scope here – see e.g. the recent theses for a thorough account. Particularly relevant to our work is the recent work of Cheng, Diakonikolas, and Ge who design an algorithm for robust mean estimatin with running time $\widetilde{O}(nd/\epsilon^{6})$ – the first to achieve nearly linear time for constant $\epsilon$ – by appeal to nearly linear time solvers for packing and covering semidefinite programs . Our algorithms carry two advantages over this prior work: first, our algorithm runs in nearly linear time for any choice of $\epsilon=\epsilon(n,d)$, and second, because we avoid the $1/\epsilon^{6}$ scaling and appeal to semidefinite programming, our theoretical ideas lead to a practical method for outlier detection. The techniques of Diakonikolas et al. were later extended to robust covariance estimation ; it remains an interesting direction to extend our techniques to covariance estimation.

Concurrent work: After this manuscript was initially submitted, we became aware of the concurrent work , which also obtains a nearly-linear time algorithm for robust mean estimation of distributions with bounded covariance. The algorithm of also obtains subgaussian confidence intervals (see e.g. ), which the algorithm in this work does not. By contrast, the algorithms in our work also obtain improved rates of error with respect to $\epsilon$ when the underlying distribution is sub-Gaussian, and our method is sufficiently practical that we are able to implement parts of it to run our experiments on outlier detection. (The method of relies on nearly-linear time solvers for packing/covering semidefinite programs, which are not yet practical.) Finally, implicit in the work is a reduction from arbitrary $\epsilon$ to the case $\epsilon=1/100$; we describe this reduction and some consequences in Appendix B.1.

Outlier detection Detection of outliers goes back nearly to the beginning of statistics itself . Even restricting to the high dimensional case it has a literature too broad to survey here. Much recent work has focused on so-called local outlier factor-based methods, which assign outlier scores based on the local density of other samples near each $X_{i}$ – see e.g. and further references in . We find that QUE scoring compares favorably to such local methods in high-dimensional datasets like we describe in Section 1.1 – see Sections 1.5 and 7 for details.

4 Robust mean estimation: results and algorithm overview

We turn to our algorithm for robust mean estimation, deferring details to Sections 2—4.

Let $D$ be a distribution on $\R^{d}$. We say that $X_{1},\ldots,X_{n}$ are an $\epsilon$-corrupted set of samples from $D$ if they are first drawn i.i.d. from $D$, then modified by an adversary who may adaptively inspect all the samples, remove $\epsilon n$ of them, and replace them with arbitrary vectors in $\R^{d}$.

Note that $\epsilon$-corruption is a stronger outlier model than the $(1-\epsilon)D+\epsilon N$ mixture model we described in Section 1; our algorithms also work in this milder mixture model. Our main theoretical result is:

For every $n,d\in\N$ and $\epsilon>0$ there are algorithms QUEScoreFilter ,s.g.-QUEScoreFilter with running time $\widetilde{O}(nd)$, such that for every distribution $D$ on $\R^{d}$ with mean $\mu$ and covariance $\Sigma$, given $n$ $\epsilon$-corrupted samples from $D$, QUEScoreFilter produces $\hat{\mu}$ such that $\|\hat{\mu}-\mu\|_{2}\leq O(\sqrt{\epsilon})+\widetilde{O}(\sqrt{d/n})$ if $\Sigma\preceq\Id$, and s.g.-QUEScoreFilter produces $\hat{\mu}$ such that $\|\hat{\mu}-\mu\|_{2}\leq O(\epsilon\sqrt{\log(1/\epsilon)}+\sqrt{d/n})$ if $D$ is sub-Gaussian with $\Sigma=\Id$, all with probability at least $0.99$.

For the bounded covariance case, the $O(\sqrt{\epsilon})$ term information-theoretically optimal up to constant factors. The other term, $\widetilde{O}(\sqrt{d/n})$, is information-theoretically optimal up to the logarithmic factors in the $\widetilde{O}(\cdot)$ even without corruptions. For the sub-Gaussian case, the $O(\varepsilon\sqrt{\log 1/\epsilon})$ term is believed to be necessary for computationally efficient algorithms (see e.g the statistical-query lower bound ), although that term can be made $O(\epsilon)$ by using computationally-intractable estimators such as Tukey median, and the latter is information-theoretically optimal . The $\sqrt{d/n}$ term is information-theoretically optimal even without corruptions.

In this section we discuss our algorithm for the bounded-covariance case $\Sigma\preceq\Id$ in the setting that the adversary may not remove samples, leaving technical details and the modifications necessary to handle removed samples and sub-Gaussian $D$ to Section 4.

Let $S=\{X_{1},\ldots,X_{n}\}\subseteq\R^{d}$ be a dataset with the property that $S$ partitions into $S=S_{g}\cup S_{b}$ with $|S_{b}|\leq\epsilon n$ and $\E_{i\sim S_{g}}(X_{i}-\mu_{g})(X_{i}-\mu_{g})^{\top}\preceq\Id$, where $\mu_{g}=\E_{i\sim S_{g}}X_{i}$. Given $S$, the goal is to find a vector $\hat{\mu}$ with $\|\mu_{g}-\hat{\mu}\|_{2}\leq O(\sqrt{\epsilon})$.

Like prior algorithms for robust mean estimation, ours maintains a weight vector $w_{1},\ldots,w_{n}\geq 0$ with $\sum w_{i}\leq 1$, initialized to $w_{i}=1/n$. The algorithm iteratively decreases the weight of points suspected to be outliers that are causing $\|\mu(w)-\mu_{g}\|_{2}$ to be large. Some prior algorithms, e.g. the filter of instead iteratively throw out points suspected to be outliers. However, since those algorithms are (necessarily) randomized, they can also be viewed as weighting points, where the weight of $X_{i}$ is the probability it has not been thrown out. The algorithm we present here can also be implemented by throwing out points in a randomized fashion – we discuss further in the appendix. A key insight of recent work on robust mean estimation is that it suffices to find weights $w$ which place almost as much mass on $S_{g}$ as does the uniform weighting and whose empirical covariance is small. This is formalized in the following lemma. For a weight vector $w$, let $|w|=\sum w_{i}$, $\mu(w)=\tfrac{1}{|w|}\sum w_{i}X_{i}$, and $M(w)=\tfrac{1}{|w|}\sum w_{i}(X_{i}-\mu(w))(X_{i}-\mu(w))^{\top}$. Let $\|M\|_{2}$ be the spectral norm of a matrix $M$.

Let $S=\{X_{1},\ldots,X_{n}\}$ be as in Definition 1.3. Suppose that $w$ is a weight vector such that $\|M(w)\|_{2}\leq O(1)$ and $w$ is mostly good, by which we mean $|\tfrac{1}{n}\textbf{1}_{S_{g}}-w_{g}|\leq|\tfrac{1}{n}\textbf{1}_{S_{b}}-w_{b}|$, where $\textbf{1}_{S_{g}},\textbf{1}_{S_{b}}$ are the indicators of $S_{g},S_{b}$ and $w_{g},w_{b}$ are $w$ restricted to $S_{g},S_{b}$ respectively. (Intuitively, $w$ is mostly good if it results by removing from the uniform weighting $\textbf{1}_{S}/n$ more weight from $S_{b}$ than from $S_{g}$.) Then $\|\mu(w)-\mu_{g}\|_{2}\leq O(\sqrt{\epsilon})$.

Lemma 1.2 captures the following geometric intuition: if the bad points $S_{b}$ receive enough weight in $w$ to cause $\|\mu(w)-\mu_{g}\|_{2}\gg\sqrt{\epsilon}$, then an $O(\epsilon)$-fraction of the mass of $w$ is on $X_{i}$ which are unusually correlated with the vector $\mu(w)-\mu_{g}$, which leads to a large maximum eigenvalue in $M(w)$. Prior works employ a variety of methods to find a mostly good weight vector $w$ with $\|M(w)\|_{2}\leq O(1)$. Perhaps the simplest is the filter of , which iterates: While $\|M(w)\|_{2}\gg 1$, compute its top eigenvector $v$ and naive spectral scores $\tau_{i}=\langle X_{i}-\mu(w),v\rangle^{2}$. Throw out $X_{i}$ with large $\tau_{i}$ and repeat.

The filter ensures that the weight vector it maintains is mostly good because (in an averaged sense) $\tau_{i}$ can be large only for $X_{i}$ which are corrupted. This is because the (weighted) sum of all scores $\sum w_{i}\tau_{i}=\langle M(w),vv^{\top}\rangle\gg 1$, while the contribution to this sum from $S_{g}$ has $\sum_{i\in S_{g}}w_{i}\tau_{i}\approx\langle\tfrac{1}{n}\sum_{i\in S_{g}}(X_{i}-\mu_{g})(X_{i}-\mu_{g})^{\top},vv^{\top}\rangle\leq 1$. (Here we ignore some details about centering $X_{i}$ at $\mu_{g}$ rather than $\mu(w)$.) Thus, the $\tau_{i}$ from $S_{b}$ must make up almost all of $\sum w_{i}\tau_{i}$. Simple approaches to removing or downweighting $X_{i}$ with large $\tau_{i}$ then remove strictly more weight from $S_{b}$ than from $S_{g}$.

However, filtering based on naive spectral scores alone faces a barrier to achieving nearly-linear running-time. If the corruptions $S_{b}$ are split among many orthogonal directions, the naive spectral filter will have to find those directions one at a time. Thus, it may require $\Omega(d)$ iterations (leading to $\Omega(nd^{2})$ running time) to arrive at $w$ with $\|M(w)\|_{2}\leq O(1)$.

Our main idea is that by replacing naive spectral scores with slightly modified QUE scores, each iteration of the filter can take into account projections of each sample onto many large eigenvectors of $M(w)$. We show that our modified QUE scores $\tau_{i}$ maintain the property that $\sum_{i\in S_{b}}w_{i}\tau_{i}\gg\sum_{i\in S_{g}}w_{i}\tau_{i}$, and so downweighting according to $\tau_{i}$ removes more mass from $S_{b}$ than $S_{g}$. However, filtering with QUE scores makes faster progress than with naive spectral scores: roughly speaking, we show that only $O(\log d)^{2}$ rounds of filtering according to QUE scores are required to find a mostly-good weight vector $w$ with $\|M(w)\|_{2}\leq O(1)$.

The core of our algorithm is a subroutine, DecreaseSpectralNorm, to take a mostly good weight vector $w$ with $\|M(w)\|_{2}\gg 1$ and in $O(\log d)$ rounds of QUE filtering produce another mostly good $w^{\prime}$ with $\|M(w^{\prime})\|_{2}\leq\tfrac{3}{4}\|M(w)\|_{2}$. Repeating this subroutine $O(\log d)$ times and then outputting the resulting $\mu(w)$ yields our main algorithm. An outline of this subroutine is presented as Algorithm 1. We first establish a rigorous sense in which downweighting according to outlier scores $\tau_{i}$ makes progress: it decreases the weighted average of the scores while removing more weight from bad points than good.

There is a downweighting algorithm which takes a density matrix $U$ and a mostly good weight vector $w$ and produces a mostly good weight vector $w^{\prime}$ by downweighting points with large score $\tau_{i}=\langle X_{i}-\mu(w),U(X_{i}-\mu(w)\rangle$ such that $\sum w_{i}^{\prime}\tau_{i}\leq\tfrac{1}{3}\sum w_{i}\tau_{i}$ so long as $\sum w_{i}\tau_{i}\gg 1$. Furthermore, $M(w^{\prime})\preceq M(w)$.

Let us give a geometric interpretation to Lemma 1.3: it establishes that if $\sum w_{i}\tau_{i}=\langle U,M(w)\rangle\gg 1$ then the quadratic form of $M(w^{\prime})$ decreases in the directions defined by $U$, since

This guarantee becomes more meaningful as the entropy $S(U)$ increases, because it suggests the quadratic form of $M(w)$ has decreased in more directions. To make this formal, we appeal to the matrix multiplicative weights framework. DecreaseSpectralNorm applies downweighting iteratively using a sequence of entropy-maximizing density matrices $U_{1},\ldots,U_{T}$ chosen according to the matrix multiplicactive weights update rule, leading to a series of mostly good weight vectors $w_{1},\ldots,w_{T}$ such that $\|M(w_{T})\|_{2}\leq\tfrac{3}{4}\|M(w_{0})\|_{2}$. We choose

where $w_{0}=w$ is the input weight vector, $U_{0}=\Id$, and $w_{t}$ results from applying the downweighting of Lemma 1.3 to $w_{t-1}$ using $U_{t}$ (if $\langle M(w_{t-1}),U_{t}\rangle\gg 1$). The following lemma is a special case of the standard (local norm) regret bound for matrix multiplicative weights.

For any $w_{0},\ldots,w_{T}$, if $\alpha\leq 1/\|M(w_{t})\|_{2}$ for all $t\leq T$, then

Now we sketch the analysis of DecreaseSpectralNorm.

If $w=w_{0}$ is mostly good, with $\|M(w_{0})\|_{2}\geq 100$, then DecreaseSpectralNorm produces mostly good $w_{T}$ with $\|M(w_{T})\|_{2}\leq\tfrac{3}{4}\|M(w)\|_{2}$.

Since $M(w_{t})\preceq M(w_{t+1})$ by Lemma 1.3, we have $\|M(w_{t})\|_{2}\leq\|M(w_{0})\|_{2}$ for all $t$, and hence $\alpha=1/\|M(w_{0})\|_{2}\leq 1/\|M(w_{t})\|_{2}$ for all $t$, so $w_{0},\ldots,w_{T}$ and $U_{0},\ldots,U_{T-1}$ satisfy the hypotheses of Lemma 1.4. By our choice of $\alpha$ and $M(w_{T})\preceq M(w_{t})$ for all $t$, (4) implies

If $\langle U_{t},M(w_{t-1})\rangle\geq\|M(w_{0})\|_{2}/3\gg 1$, then DecreaseSpectralNorm performs downweighting, and by Lemma 1.3 and (2) (which we establish rigorously in supplemental material), $\langle M(w_{t}),U_{t}\rangle\leq\tfrac{1}{3}\langle M(w_{t-1}),U_{t}\rangle\leq\tfrac{1}{3}\|M(w_{0})\|$. Otherwise, by hypothesis $\langle M(w_{t}),U_{t}\rangle=\langle M(w_{t-1}),U_{t}\rangle\leq\|M(w_{0})\|_{2}/3$. Using this bound and dividing by $T$, we obtain $\|M(w_{T})\|_{2}\leq(\tfrac{2}{3}+\tfrac{\log d}{T})\|M(w_{0})\|_{2}$. Choosing $T\geq 20\log d$ completes the proof sketch. ∎

Running time: Our overall algorithm only requires $\log(nd)^{O(1)}$ iterations of DecreaseSpectralNorm, and the latter only requires $O(\log(d))$ iterations of downweighting, so we just have to implement downweighting in nearly-linear time. We show in supplemental material that this can be done by avoiding representing any of the matrices $U_{t}$ explicitly in memory: instead, we maintain only low-rank sketches of them. This leads to some approximation error in computing the QUE scores, but we show that approximations to the QUE scores suffice for all arguments above.

For remaining technical details and full proofs, see Sections 5-9 of supplemental materials.

5 Outlier detection: algorithm and experimental results

In this section, we empirically evaluate outlier detection using QUE scoring. We must work with data containing well-defined and known inliers and outliers so that we can compare our results to ground-truth. We generate such data sets in three distinct ways, leading to three main experiments. (In supplemental material we also study some outlier-detection data sets appearing in prior work .)

Synthetic: We create synthetic data sets in $128$ dimensions and $10^{3}-10^{4}$ samples with an $\epsilon$-fraction of inhomogeneous outliers in $k$ directions by sampling from a mixture of $k+1$ Gaussians $(1-\epsilon)\mathcal{N}(0,\Id)+\sum_{i=1}^{k}\epsilon_{i}[\tfrac{1}{2}\mathcal{N}(C\sqrt{k/\epsilon}\cdot e_{i},\sigma^{2}\Id)+\tfrac{1}{2}\mathcal{N}(-C\sqrt{k/\epsilon}\cdot e_{i},\sigma^{2}\Id)]$, where $e_{1},\ldots,e_{k}$ are standard basis vectors, with $C\approx 1$ and $\sigma\ll 1$. The outliers are the samples from $\mathcal{N}(\pm C\sqrt{k/\epsilon}e_{i},\sigma^{2}\Id)$. By varying $\epsilon,k$ and the distribution $\epsilon_{1},\ldots,\epsilon_{k}$ of outlier weights, we demostrate in this simplified model how max-entropy outlier scoring improves on baseline algorithms in the presence of inhomogeneous outliers. We choose the scaling $\sqrt{k/\epsilon}\cdot e_{i}$ because then standard calculations predict that if $\epsilon_{i}\approx\epsilon/k$ the outliers from $\mathcal{N}(\pm C\sqrt{k/\epsilon}e_{i},\sigma^{2}\Id)$ will contribute an eigenvalue greater than $1$ to the overall empirical covariance.

Mixed – word embeddings: We create a data set consisting of word embeddings drawn from several sources. Inliers are the $100$-dimensional GloVe embeddings () of the words in a random $\approx 10^{3}$ word long section of a novel (we use Sherlock Holmes) and outliers are embeddings of the first paragraphs of $k$ featured Wikipedia articles from May 2019 .

Perturbed – images: We create a data set consisting of CIFAR10 images some of which have artificially-introduced dead pixels. Inliers are $\approx 4500$ random CIFAR images $X\in\{1,\ldots,256\}^{1024}$ (restricted to the red color channel). Outliers are $\approx 500$ random CIFAR images, partitioned into groups $S_{1},\ldots,S_{k}$, such that for each group $i$ a random coordinate $p_{i}\in\{1,\ldots,1024\}$ and a random value $c_{i}\in\{1,\ldots,256\}$ is chosen and for each $X\in S_{i}$ we set $X_{p_{i}}=c_{i}$.

Metric: All the methods we evaluate produce a vector of scores $\tau_{1},\ldots,\tau_{n}\in\R$. We use the standard ROCAUC metric to compare these scores to a ground-truth partition $S=S_{g}\cup S_{b}$ into inlier and outlier sets. $\text{ROCAUC}(\tau_{1},\ldots,\tau_{n},S_{b},S_{g})=\Pr_{i\sim S_{b},j\sim S_{g}}(\tau_{i}\geq\tau_{j})$ is simply the probability that a randomly chosen outlier is scored higher than a random inlier.

Whitening: Scoring methods based on the projection of data points $X_{i}$ onto large eigenvectors of the empirical covariance work best when those eigenvectors correspond to directions in which many outliers lie. In particular, if $\Sigma_{g}$, the covariance of $S_{g}$, itself has large eigenvalues then such spectral methods perform poorly. We assume access to a whitening transformation $W\in\R^{d\times d}$, which captures a small amount of prior knowledge about the distribution of inliers $S_{g}$. For best performance $W$ should approximate $W^{*}=(\Sigma_{g})^{-1/2}$ since $W^{*}X_{i}$ form an isotropic set of vectors. Of course, to compute $W^{*}$ exactly would require knowing which points are inliers, but we find that relatively naive approximations suffice. In particular, if a clean dataset $Y_{1},\ldots,Y_{m}$ whose distribution is similar to the distribution of inliers is available, its empirical covariance can be used to find a good whitening transformation $W$. In our synthetic data we use $W=\Id$. In our word embeddings experiment, we obtain $W$ using the empirical covariance of the embedding of another random section of Sherlock Holmes. In our CIFAR-10 experiment, we obtain $W$ from the empirical covariance of a fresh sample of $\approx 5000$ randomly chosen images from CIFAR-10.

High-dimensional scaling: Implementing Algorithm 2 by explicitly forming the matrix $\overline{\Sigma}$ and performing a singular value decomposition (SVD) to compute $\exp(\alpha\overline{\Sigma})$ is feasible on relatively low-dimensional data ($d\approx 100$). See Section 7 for discussion and results of a nearly-linear time implementation.

Robust mean estimation: results and preliminaries

for all $k$ even. We say a distribution $D$ over $\R^{d}$ and mean $\mu$ is sub-gaussian with variance proxy $\Sigma\preceq I$, if for all unit vectors $v$, the distribution of $\left\langle v,X-\mu\right\rangle$ is sub-Gaussian with variance proxy $v^{\top}\Sigma v$. Intuitively, a sub-Gaussian distribution is simply any distribution which concentrates as well as a Gaussian.

With this terminology in place, we are now ready to state our main results on robust mean estimation. Our first result is for robust mean estimation under the assumption of bounded covariance:

Let $D$ be a distribution on $\R^{d}$ with mean $\mu$ and covariance $\Sigma\preceq I$. Let $\varepsilon>0$ be sufficiently small, and let $\delta>0$. Let $S$ be an $\varepsilon$-corrupted set of samples from $D$ of size $n$. There is an algorithm $\textsc{QUEScoreFilter}(S,\delta,\epsilon)$ which takes $S$ and $\delta$, and outputs $\widehat{\mu}$ so that with probability $1-\delta-\exp(-\varepsilon n)$, we have

Moreover, the algorithm runs in time $\widetilde{O}(nd\log 1/\delta)$.

We make two observations about this problem. First, it is well-understood (see e.g. ) that $\Omega(\sqrt{\varepsilon})$ error is unavoidable for this problem, no matter how many samples are given. Second, observe that the rate $O(\sqrt{d/n})$ is necessary for this problem even without corruptions. Thus, up to log factors, and the dependence on $\delta$, this guarantee is information-theoretically optimal. We note also that the algorithm does not need to know $\epsilon$ exactly; any upper bound will suffice (with commensurate weakening in the error rate).

We also prove a strong statement for the case of robust mean estimation for sub-gaussian distributions:

Let $D$ be an isotropic sub-gaussian distribution with variance proxy $I$ and mean $\mu$. Let $\varepsilon>0$ be sufficiently small, and let $\delta>0$. Let $S$ be an $\varepsilon$-corrupted set of samples from $D$ of size $n$. There is an algorithm $\textsc{s.g.-QUEScoreFilter}(S,\delta,\varepsilon)$ which takes $S,\delta$, and $\varepsilon$, and outputs $\widehat{\mu}$ so that with probability $1-\delta$, we have

Moreover, the algorithm runs in time $\widetilde{O}(nd\log 1/\delta)$.

It is suspected, based on statistical-query lower bounds, that $\Omega(\varepsilon\sqrt{\log 1/\varepsilon})$ error is incurred by any computationally-efficient algorithm in this setting, although $\Theta(\epsilon)$ is the minimax optimal dependence of the error rate on $\epsilon$ . Moreover, $\Omega\left(\sqrt{(d+\log 1/\delta)/n}\right)$ is the minimax rate for mean estimation for Gaussians without noise. Thus, the error guarantees of this algorithm are minimax optimal up to constants and the factor $\sqrt{\log(1/\epsilon)}$.

This algorithm assumes that the distribution is isotropic. There is evidence that such an assumption is necessary to get error beyond $\sqrt{\varepsilon}$ using computationally efficient (i.e. poly-time) algorithms . In our theorem statement above we also assume that the variance proxy is at most $I$. This is done for simplicity: it is easily verifiable that our algorithm works (with an appropriate scaling in front of the error guarantee) if the variance proxy is PSD upper bounded by $\sigma^{2}I$ for any $\sigma^{2}$.

Before we describe our techniques, we require a few additional algorithmic tools, which we describe below.

2 Soft selection of subsets of points

In our presentation of our filtering algorithm for robust mean estimation, it will be convenient for us to work with a “soft” version of the filter. Instead of wholly removing points that we deem suspicious, we will maintain a set of weights for each point, and downweight those that we find suspicious. In this section, we establish notation for dealing with such operations. However, we briefly remark that, as we will explain later in Appendix A.3, the same results (up to log factors) can be established using “hard” filtering more akin to the algorithms presented in prior work, e.g. .

Throughout this paper, we will let $\Delta_{n}$ denote the simplex in $n$ dimensions, and we let

Typically when the set $S$ is understood, we will omit the dependence on $S$ in the notation. For any set $T\subseteq S$, we let $\mu(T)=\mu(w)$ where $w\in{\Gamma}_{n}$ is the vector $w_{i}=1/n$ for $i\in T$ and $w_{i}=0$ otherwise, and similarly we let $M(T)=M(w)$.

3 Naive pruning

One primitive we will require will be the ability to removes points which are “obviously” outliers. It is well-known that there exist randomized nearly-linear time algorithms for achieving this. For completeness we prove this lemma in Appendix A.

There is an algorithm NaivePrune with the following guarantees. Let $\varepsilon>1/2$, and let $\delta>0$. Let $S\subset\R^{d}$ be a set of $n$ points so that there exists a ball $B$ of radius $r$ and a subset $S^{\prime}\subseteq S$ so that $|S^{\prime}|\geq(1-\varepsilon)n$, and $S^{\prime}\subset B$. Then $\textsc{NaivePrune}(S,r,\delta)$ runs in time $O(nd\log 1/\delta)$ and with probability $1-\delta$ outputs a set of points $T\subseteq S$ so that $S^{\prime}\subseteq T$, and $T$ is contained in a ball of radius $4r$.

In the case where the output of NaivePrune satisfies the conditions of the lemma, we say that NaivePrune succeeds.

4 The one-dimensional filter

An important algorithmic primitive for us will be an univariate soft outlier removal step. The sub-problem considered here is as follows: we are given a set of nonegative scores $\tau_{1},\ldots,\tau_{m}$, with the guarantee that there is a small subset $S\subseteq[m]$ so that $\sum_{i\in S}\tau_{i}>\frac{1}{2}\sum_{i=1}^{m}\tau_{i}$, that is, they contribute a majority of the mass of the points. The goal is to then either downweight (or remove) the overall set of scores in such a way so that more mass from $S$ is removed than from outside of $S$, or alternatively, more points are removed from $S$ than from outside $S$. An algorithm for achieving this via downweighting has already been described in , and a randomized algorithm that achieves the same sorts of guarantees with high probability by removing points is implicit in the filtering algorithm of (e.g. in Algorithm 3 in Appendix A of ). In this paper, we will require a slight strengthening of these algorithms. We require that not only do we remove more weight from the bad points than the good points, but we also decrease the overall sum by a constant factor. We observe that while we will present a method for acheving this via downweighting, one can achieve the same guarantee (with high probability) by removing points. In the main text, we choose to present the soft downweighting method for robust mean estimation for simplicity. See Appendix A.3 for details.

Formally, we describe an algorithm 1DFilter and prove the following guarantee for the algorithm. The algorithm and its analysis are fairly straightforward so we defer the formal descriptions and proofs to Appendix A.2.

Let $\eta\in(0,1/2)$, let $b\geq 2\eta$, and let $w_{1},\ldots,w_{m}$ and $\tau_{1},\ldots,\tau_{m}$ be non-negative numbers so that $\sum_{i=1}^{m}w_{i}\leq 1$. Let $\tau_{\max}=\max_{i\in[m]}\tau_{i}$. Suppose there exist two disjoint sets $S_{g},S_{b}$ so that $S_{g}\cup S_{b}=[m]$, and moreover,

Then $\textsc{1DFilter}(w,\tau,b)$ runs in time $O\left(\left(1+\log\frac{\tau_{\max}}{b\sigma}\right)m\right)$ and outputs $0\leq w^{\prime}\leq w$ so that:

more weight is removed from $S_{b}$ than $S_{g}$, i.e. $\sum_{i\in S_{g}}w_{i}-w_{i}^{\prime}\leq\sum_{i\in S_{b}}w_{i}-w_{i}^{\prime}$, and

the weighted sum of the $\tau$ has decreased, i.e. $w^{\prime}$ satisfies

In particular, note that if $b=\Omega(1)$ and $\tau_{\max}/\sigma\leq m^{O(1)}$, then this algorithm runs in nearly linear time.

As mentioned previously, there is also a randomized strategy that avoids downweighting and achieves the same guarantee with high probability (up to logarithmic factors in runtime). Our overall robust mean estimation algorithm (for both settings presented in the paper) can be instantiated using this algorithm rather than 1DFilter. While as far as we know this yields no theoretical improvements for robust mean estimation (indeed, our analysis of it proves bound which are worse by logarithmic factors than our analysis of soft downweighting), it is much closer to the practical outlier detection method used in the experiments in Section 1.5 and also to prior algorithms presented in , and may be of instructive value. For this reason we describe this algorithm in Appendix A.3.

5 Matrix Multiplicative weights

We will use the following form of the MMW update, which is essentially the same as presented in . In each iteration $t=0,\ldots,T$, the player chooses an action $U_{k}\in\Delta_{d\times d}$, receives a gain matrix $F_{t}\in\R^{d\times d}$, and receives reward $\langle F_{k},U_{k}\rangle$. Then the player sees $F_{k}$. In , they demonstrate that if the player plays according to the entropy regularizer (or equivalently, matrix multiplicative weights), namely,

Here, for any symmetric matrix $A=\sum_{i=1}^{d}\lambda_{i}v_{i}v_{i}^{\top}$, we let $|A|$ denote $|A|=\sum_{i=1}^{d}|\lambda_{i}|v_{i}v_{i}^{\top}$. Equivalently, by rearranging terms, and taking a supremum over $U$ of (8), we obtain that the update satisfies

An MMW algorithm for robust mean estimation with bounded covariance

In this section we describe the algorithm which achieves Theorem 2.1. We first identify a deterministic condition on the set of inliers under which our algorithm is guaranteed to be correct. It is a very mild condition: at a high level, it simply states that the empirical mean of the samples is converging to the true mean, and the empirical covariance is bounded.

We say a set of points $S$ is $(\gamma_{1},\gamma_{2})$-good with respect to a distribution $D$ with mean $\mu$ and covariance $\Sigma\preceq\Id$ if the following two properties hold:

The following is a generalization of Lemma A.18 in , which states that, with high probability, any set of i.i.d. points from a distribution with bounded covariance will contain a large set which is good with respect to that distribution. For completeness we prove this lemma in Appendix B.

Let $\varepsilon\in[0,1/2)$, and let $n$ be a positive integer. Let $D$ be a distribution with mean $\mu$ and covariance $\Sigma\preceq\Id$, and let $X_{1},\ldots,X_{n}$ be independent draws from $D$. Then, with probability $1-\delta-\exp(-\varepsilon n)$, there exists a set $S\subseteq\{X_{1},\ldots,X_{n}\}$ so that the following two conditions are simultaneously satisfied:

$S$ is $(\gamma_{1},\gamma_{2})$-good with respect to $D$, where

for some universal constants $c,c^{\prime}>0$.

In particular, observe that for constant $\delta$ and for $n=\Omega(d\log d/\varepsilon)$, we have that $\gamma_{1}=O(\varepsilon)$ and $\gamma_{2}=O(1)$.

Throughout we will let $S=S_{g}\cup S_{b}\setminus S_{r}$, where $S_{g}$ is $(\gamma_{1},\gamma_{2})$-good with respect to $D$, and $|S_{b}|,|S_{r}|\leq\varepsilon|S|$. For any $w\in{\Gamma}_{n}$, we let $w_{g}\in{\Gamma}_{|S_{g}|}$ denote the restriction of $w$ to the indices in $S_{g}$, and similarly define $w_{b}\in{\Gamma}_{|S_{b}|}$.

Our set of weights of interest will be slightly different than those considered in prior papers, but morally captures the same concept, up to issues of reweighting. We will always guarantee that the weights we consider lie within the following set:

Intuitively, weights in the set $\mathfrak{S}_{n,\varepsilon}$ are what happens when we start with the uniform weighting $\tfrac{1}{n}\cdot\mathbf{1}$ and remove weight from points in the data set, always removing at least as much mass from the bad set as we do from the good set.

2 Geometric lemmata

We first prove the following sequence of structural lemmata. The first, which is implicit in the earlier work, and which is in some sense the fundamental geometric fact which guides our algorithmic design, gives an upper bound on the deviation between the weighted empirical mean of the data set and the true mean of the distribution in terms of the spectral norm of the weighted covariance of the dataset.

Let $S,\gamma_{1},\gamma_{2}$ be as above, and let $w\in\mathfrak{S}_{n,\varepsilon}$. Then

Let $\rho=\mu(w)-\mu$. We have the following sequence of identities:

We now upper bound each term separately. By Cauchy-Schwarz, we have

We now turn our attention to $W_{2}$ and $W_{3}$. Both bounds will follow from the following claim:

Let $w^{\prime},\alpha\in{\Gamma}_{n}$ be so that $|w^{\prime}|\leq\varepsilon$, $w^{\prime}_{i}\leq\frac{1}{n}$ for all $i\in[n]$, and $w^{\prime}\leq\alpha$. Then, for any $v\in\R^{d}$, we have

Here (a) follows since $w_{i}^{\prime}\leq\alpha_{i}$, and (b) follows from the definition of spectral norm, and the assumption on $|w^{\prime}|$. Thus, by taking square roots and combining terms, we have

With this claim, we can now bound $W_{2}$ and $W_{3}$. To bound $W_{2}$, let $w^{\prime}_{i}=\frac{1}{n}-w_{i}$ for $i\in S_{g}$ and $w_{i}=0$ otherwise, and let $\alpha_{i}=\frac{1}{n}$ if $i\in S_{g}$, and $0$ otherwise. Then, applying the claim with $v=\rho$ yields that

Similarly, to bound $W_{3}$, let $w^{\prime}_{i}=w_{i}$ if $i\in S_{b}$ and $w^{\prime}_{i}=0$ otherwise, and let $\alpha_{i}=w_{i}$ for all $i\in S$. Again, letting $v=\rho$, we get that

as well. Combining these three bounds, and using the fact that $|w|\leq 1$, yields that

Simplifying this expression then yields the desired bound on $\|\rho\|_{2}$. ∎

We also require the following linear algebraic fact:

Let $w^{\prime},w\in{\Gamma}_{n}$ so that $w^{\prime}\leq w$. Then $|w^{\prime}|M(w^{\prime})\preceq|w|M(w)$.

We first observe that if $w^{\prime}\leq w$, then

and so $M(w^{\prime})\preceq M(w)$, as claimed. ∎

3 General algorithm description, bounded second moment

In this section we describe the algorithm we would like to run via matrix multiplicative weights, and demonstrate that it will terminate in a small number of iterations, when given good approximations to the entropic scores.

The algorithm, which we call QUEScoreFilter , proceeds in epochs, and takes as input a corrupted dataset $S$, and a score oracle $\mathscr{O}$. Initially, in epoch $s=0$, we let $w^{\left(0\right)}=\frac{1}{n}\mathbf{1}_{n}$. Then, in epoch $s$, the algorithm proceeds iteratively as follows. First, approximately compute $\lambda^{\left(s\right)}\approx_{0.1}\|M(w^{\left(s\right)})\|_{2}$, and if $\lambda^{\left(s\right)}\leq 100\gamma_{2}$, then we terminate and output $\mu(w^{\left(s\right)})$.

where $U_{t}^{\left(s\right)}$ is as in Algorithm 3, namely,

for all $t,i,s$, where $\tau_{t,i}^{\left(s\right)}$ is defined in (10) (the choice of $0.1$ in the approximation here as well as for $\lambda^{(s)}_{t}$ is arbitrary; any constant sufficiently small will suffice). In Section 5 we will demonstrate how to implement such an approximate score oracle in (randomized) nearly-linear time by sketching.

4 Correctness of QUEScoreFilter

We first prove correctness. The main lemma is the following per-epoch guarantee:

The following invariants always hold. For all epochs $s$, we have:

$w^{s}\in\mathfrak{S}_{n,\varepsilon}$, and

If $\|M\left(w^{\left(s\right)}\right)\|_{2}>\frac{100}{1.1}\gamma_{2}$ and $\|M\left(w^{\left(s\right)}\right)\|_{2}>\frac{100}{1.1}\gamma_{1}^{2}$, then epoch $s$ finishes after $O(\log d)$ iterations, and outputs $w^{\left(s+1\right)}$ so that $\|M\left(w^{\left(s+1\right)}\right)\|_{2}\leq\frac{2}{3}\|M\left(w^{\left(s\right)}\right)\|_{2}$.

We first show how the lemma implies the theorem.

We now prove the runtime bound. In every iteration, besides the call the the oracle, the only costly operations are the approximate top eigenvalue computations and running 1DFilter. However, the approximate top eigenvalue computations can be done in time $\widetilde{O}(nd)$ via power method since we only ask for a constant multiplicative approximation, and the bound on $\left\lVert X_{i}\right\rVert_{2}$ implies that 1DFilter runs in $O(n\log\kappa)$ time. ∎

We first show that $\sum_{i\in S_{g}}w_{i}\tau_{i}\leq c\sum_{i=1}^{n}w_{i}\tau_{i}$ for some universal constant $c\leq 0.11$. Let ${\widetilde{w}}_{i}=\frac{1}{n}$ if $i\in S_{g}$ and ${\widetilde{w}}_{i}=0$ otherwise. Then, we have

Here (a) follows from $M(w^{\prime})\preceq\sum_{i=1}^{n}w_{i}^{\prime}(X_{i}-\mu(w))(X_{i}-\mu(w))^{\top}$. This completes the proof of the claim. ∎

Notice that Claim 3.7 immediately implies that the first invariant $w\in\mathfrak{S}_{n,\varepsilon}$ always holds. We now turn to proving the second invariant. Let $T$ be the number of iterations that the epoch runs for. Observe that for all $t=0,\ldots,T$, we have that $M(w_{t})\preceq M(w_{0})\preceq\frac{1}{\alpha}I$, and so we are indeed in the setting of the guarantee in Section 2.5. Thus, by (8), since $F_{t}=M(w_{t+1})$, we obtain the following regret bound:

where the second inequality follows by our choice of $\alpha$. We claim that for all $t>0$, we must have $\left\langle M(w_{t}),U_{t}\right\rangle\leq 0.31\|M(w_{0})\|_{2}$. There are two cases. If we enter the if statement in Line 18, then

and $M(w_{t+1})=M(w_{t})$, so this is clearly satisfied. Otherwise, the desired bound follows by Claim 3.7. Thus overall, by (14) we have that

By Lemma 3.4, we further have that $M(w_{t+1})\preceq M(w_{t})$ for all $t=0,\ldots,T-1$, and so this implies that

Simplifying both sides yields that if $T=C\log n$ for some sufficiently large constant $C$, then $\|M(w_{T})\|_{2}\leq\frac{2}{3}\|M(w_{0})\|_{2}$. Thus after $O(\log d)$ iterations, we must terminate. ∎

An MMW algorithm for robust mean estimation for sub-gaussian distributions

In this section we give an analog of the result in Section 3 but but in the setting where the distribution $D$ is subgaussian, and has identity covariance. We again first identify a deterministic condition for the inlers under which our algorithms will succeed. In this case, we need a stricter analog of $\varepsilon$-goodness. Specifically, we will require:

Let $D$ be a distribution with covariance $\Id$ and mean $\mu$. We say a set of points $S$ is $(\varepsilon,\gamma_{1},\gamma_{2},\beta_{1},\beta_{2})$-subgaussian good (or $(\varepsilon,\gamma_{1},\gamma_{2},\beta_{1},\beta_{2})$-s.g. good for short) with respect to $D$ if there exist universal constants $C_{1},C_{2}$ so that the following inequalities are satisfied:

$\left\lVert\mu(S)-\mu\right\rVert\leq\gamma_{1}$ and $\left\lVert\frac{1}{|S|}\sum_{i\in S}\left(X_{i}-\mu(S)\right)\left(X_{i}-\mu(S)\right)^{\top}-\Id\right\rVert_{2}\leq\gamma_{2}$, and

For any subset $T\subset S$ so that $|T|=2\varepsilon|S|$, we have

For conciseness, when the parameters $\varepsilon,\gamma_{1},\gamma_{2},\beta_{1},\beta_{2}$ are understood, we will omit them and refer to the data set as s.g.-good.

We have the following concentration inequality. The proof is very similar to that of Lemmata 2.1.8 and 2.1.9 in . For completeness we include a proof of this lemma in Appendix C.

Let $X_{1},\ldots,X_{n}\sim D$, where $D$ is subgaussian with variance proxy $1$. Then, for any $\varepsilon>0$ sufficiently small, we have that $S=\{X_{1},\ldots,X_{n}\}$ is $\left(\varepsilon,\gamma_{1},\gamma_{2}\right)$-s.g. good with probability $1-\delta$, where

In particular, note that when $n=\Omega\left(\frac{d+\log 1/\delta}{\varepsilon^{2}\log 1/\varepsilon}\right)$, then Lemma 4.1 implies that $n$ i.i.d. samples from an isotropic sub-gaussian distribution is $(\varepsilon,O(\varepsilon\sqrt{\log 1/\varepsilon}),O(\varepsilon\log 1/\varepsilon),O(\sqrt{\log 1/\varepsilon}),O(\log 1/\varepsilon))$-s.g. good with probability $1-\delta$. We will also require the following simple consequences of subgaussian goodness.

Let $D$ be an isotropic distribution. Let $S$ be $(\varepsilon,\gamma_{1},\gamma_{2},\beta_{1},\beta_{2})$-s.g. good w.r.t. $D$. Then:

for all $w^{\prime}\in{\Gamma}_{n}$ with $w^{\prime}\leq\frac{1}{n}\mathbf{1}_{n}$ and $|w^{\prime}|\leq 2\varepsilon$, and for all unit vectors $v\in\R^{d}$, we have

if $w\in{\Gamma}_{n}$ satisfies $w\leq\frac{1}{n}\mathbf{1}_{n}$ and $\left|\frac{1}{n}\mathbf{1}_{n}-w\right|\leq 2\varepsilon$, then

We first prove the first claim. Let $w^{\prime\prime}\in\Gamma_{n}$ be anything so that $w^{\prime}\leq w^{\prime\prime}\leq\frac{1}{n}\mathbf{1}_{n}$ and $|w^{\prime\prime}|=\varepsilon$. Then since all quantities on the LHS of the expression are nonnegative, we have that

Now let $A=\{w^{\prime\prime}\in\Gamma_{n}:|w^{\prime\prime}|=\varepsilon\}$. This set is clearly convex, and moreover, by inspection, the vertices of $A$ are exactly given by $\frac{1}{n}\mathbf{1}_{T}$ where $|T|=2\varepsilon n$. Thus, by convexity, the maximum of the RHS of (20) over $w^{\prime\prime}\in A$ is obtained by $w^{\prime\prime}=\frac{1}{n}\mathbf{1}_{T}$ for some $T$ with $|T|=2\varepsilon n$. But then we have

by the s.g.-goodness of $S$. This completes the proof of the first bullet point.

We now turn our attention to the second claim. We have that

by the first claim. Further expanding, we have

Putting it all together yields (18). To prove (19), simply observe that

In the first inequality we have used the definition of subgaussian goodness and convexity. ∎

As a result, we also have the following tail bound on mean shifts caused by small subsets of points:

Let $D$ be an isotropic distribution. Let $S$ be $(\varepsilon,\gamma_{1},\gamma_{2},\beta_{1},\beta_{2})$-s.g. good w.r.t. $D$. Then:

for all $w^{\prime}\in{\Gamma}_{n}$ with $w^{\prime}\leq\frac{1}{n}\mathbf{1}_{n}$ and $|w^{\prime}|\leq 2\varepsilon$, we have

if $w\in{\Gamma}_{n}$ satisfies $w\leq\frac{1}{n}\mathbf{1}_{n}$ and $\left|\frac{1}{n}\mathbf{1}_{n}-w\right|\leq 2\varepsilon$, then