#  Beyond Parallel Pancakes: Quasi-Polynomial Time Guarantees for Non-
Spherical Gaussian Mixtures††thanks: This project has received funding from
the European Research Council (ERC) under the European Union’s Horizon 2020
research and innovation programme (grant agreement No 815464).

Rares-Darius Buhai ETH Zürich.  David Steurer ETH Zürich.

###### Abstract

We consider mixtures of $k\geqslant 2$ Gaussian components with unknown means
and unknown covariance (identical for all components) that are well-separated,
i.e., distinct components have statistical overlap at most $k^{-C}$ for a
large enough constant $C\geqslant 1$. Previous statistical-query lower bounds
[DKS17] give formal evidence that even distinguishing such mixtures from
(pure) Gaussians may be exponentially hard (in $k$).

We show that this kind of hardness can only appear if mixing weights are
allowed to be exponentially small, and that for polynomially lower bounded
mixing weights non-trivial algorithmic guarantees are possible in quasi-
polynomial time.

Concretely, we develop an algorithm based on the sum-of-squares method with
running time quasi-polynomial in the minimum mixing weight. The algorithm can
reliably distinguish between a mixture of $k\geqslant 2$ well-separated
Gaussian components and a (pure) Gaussian distribution. As a certificate, the
algorithm computes a bipartition of the input sample that separates some pairs
of mixture components, i.e., both sides of the bipartition contain most of the
sample points of at least one component.

For the special case of colinear means, our algorithm outputs a $k$ clustering
of the input sample that is approximately consistent with the components of
the underlying mixture.

A significant challenge for our results is that they appear to be inherently
sensitive to small fractions of adversarial outliers unlike most previous
algorithmic results for Gaussian mixtures. The reason is that such outliers
can simulate exponentially small mixing weights even for mixtures with
polynomially lower bounded mixing weights.

A key technical ingredient of our algorithms is a characterization of
separating directions for well-separated Gaussian components in terms of
ratios of polynomials that correspond to moments of two carefully chosen
orders logarithmic in the minimum mixing weight.

###### Contents

  1. 1 Introduction
    1. 1.1 Results
  2. 2 Techniques
  3. 3 Preliminaries
  4. 4 Separating polynomial
    1. 4.1 Exact moment results
    2. 4.2 Finite sample bounds
    3. 4.3 Proof of Theorem 4.1
  5. 5 Parallel pancakes
    1. 5.1 Isotropic position properties
    2. 5.2 Exact moment direction recovery
      1. 5.2.1 Sum-of-squares identifiability (proof of Theorem 5.5)
      2. 5.2.2 Sum-of-squares sampling (proof of Theorem 5.6)
    3. 5.3 Finite sample bounds
    4. 5.4 Proof of Theorem 5.1
  6. 6 Colinear means
    1. 6.1 Isotropic position transformation
    2. 6.2 Finite sample isotropic position transformation
    3. 6.3 Finite sample bounds
    4. 6.4 Proof of Theorem 6.1
  7. Acknowledgement
  8. A Appendix
    1. A.1 Bipartition from distance function
    2. A.2 Sum-of-squares lemmas
    3. A.3 Finite sample lemmas
    4. A.4 Isotropic position transformation lemmas
    5. A.5 Miscellaneous lemmas
    6. A.6 Proofs deferred from Section 4
    7. A.7 Proofs deferred from Section 5
    8. A.8 Proofs deferred from Section 6

##  1 Introduction

Gaussian mixture models (GMMs) are among the most extensively studied
statistical models, spanning a diverse range of scientific disciplines
[Pea94]. Over the course of the last two decades, a major research effort has
been aimed at understanding the kinds of algorithmic guarantees that are
feasible for GMMs [Das99, VW02, KMV10, MV10, HK13].

Recent years have seen major algorithmic advances along two dimensions.

The first kinds of advances concern mixtures of a large number of spherical
Gaussians, i.e., Gaussians with identity $I_{d}$ as covariance. Here, we now
have strong guarantees with running times depending gracefully on the number
$k$ of components. Concretely, it is possible to cluster such mixtures under a
separation requirement $O(\sqrt{\log k})$, which up to a constant factor
matches the minimum separation that guarantees clusterability of the mixture
[HL18, DKS18, KSS18]. The running time $(kd)^{(\log k)^{O(1)}}$ is quasi-
polynomial in the number $k$ of components. Even without any separation
requirement, it is possible to compute quasi-polynomially sized covers of the
set of means [DK20].

The second kinds of advances concern mixtures of a small number of Gaussians
with unknown covariances. Here, many previous algorithmic guarantees are now
possible also robustly, i.e., in the presence of a small constant fraction of
adversarially chosen outliers. Concretely, it is possible to estimate the
parameters of an arbitrary mixture of Gaussian components in the presence of
such outliers [BDH+20, LM21]. The running time is polynomial in the ambient
dimension but exponential in the number of components.

So far, mixtures of a large number of Gaussian components with unknown
covariances have defied comparable algorithmic progress.111 A notable
exception is a particular smoothed model for such mixtures when the ambient
dimension is large enough [GHK15]. Indeed, there is formal evidence to suggest
that this setting is computationally inherently harder than the spherical
setting [DKS17] and that even for $k$ components with tiny statistical
overlaps, say at most $2^{-k}$, approximately clustering the components may
require time exponential in $k$ despite the sample complexity being polynomial
in $d$ and $k$. Underlying this evidence is the well-known parallel pancakes
construction.

We show that this kind of hardness can appear only in the arguably
pathological case that mixing weights are allowed to be exponentially small.
Indeed, we develop algorithms with substantial statistical guarantees that run
in quasi-polynomial time whenever the mixing weights are bounded from below by
a polynomial. We hope that our work opens up a new direction of research on
efficient algorithms for mixtures of well-separated Gaussians with
polynomially lower bounded mixing weights. Within this new direction of
research we identify the following appealing open question:

> Consider a mixture of $k\geqslant 2$ Gaussian components with unknown means
> $\mu_{1},\ldots,\mu_{k}\in\mathbb R^{d}$ and unknown covariance
> $\Sigma\in\mathbb R^{d\times d}$ (identical for all components) and with
> minimum mixing weight $p_{\mathrm{min}}>0$. Suppose the components are well-
> separated in the sense that any two distinct components have statistical
> overlap at most $p_{\mathrm{min}}^{C}$ for a large enough constant
> $C\geqslant 1$.
>
> Given a sample of size $n\geqslant d^{O(\log(1/p_{\mathrm{min}}))}$, can we
> compute in time polynomial in $n$ a $k$-clustering of the sample that is
> consistent with the mixture components on all but at most a
> $p_{\mathrm{min}}^{10}$ fraction of the sample?

We conjecture that such an algorithm does exists. Indeed, we confirm the
conjecture for the special case that the means are colinear. In the general
case, our algorithm provides a somewhat weaker guarantee and computes only a
bipartition of the sample that separates at least some pairs of mixture
components.

###  1.1 Results

##### Separating bipartition

Suppose we are given a quasi-polynomial size sample of a uniform mixture of
$k$ Gaussian components with unknown means $\mu_{1},\ldots,\mu_{k}\in\mathbb
R^{d}$ and unknown covariance $\Sigma\in\mathbb R^{d\times d}$ such that the
maximum mean separation satisfies $\max_{a\neq
b}\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})\rVert\gg\sqrt{\log k}$. Then, as
Theorem 1.1 shows, it is possible to compute in quasi-polynomial time a
bipartition of the samples such that, for each side of the bipartition, there
exists a component with $0.99$ of its samples assigned to it.

######  Theorem 1.1.

Given a sample of size $n\geqslant(kd)^{O(\log k)}$ from a uniform mixture of
$k$ Gaussian components $N(\mu_{1},\Sigma),\ldots,N(\mu_{k},\Sigma)$ with
$\max_{a\neq b}\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})\rVert\gg\sqrt{\log k}$,
there exists an algorithm that runs in time $n\cdot d^{O(\log k)}$ and returns
with probability $0.99$ a partition of $[n]$ into two sets $C_{1}$ and $C_{2}$
such that, if true clustering of the samples is $S_{1},...,S_{k}$, then

$\max_{i}\frac{|C_{1}\cap S_{i}|}{|S_{i}|}\geqslant
0.99\quad\text{and}\quad\max_{i}\frac{|C_{2}\cap S_{i}|}{|S_{i}|}\geqslant
0.99.$

The same result holds in the case of non-uniform mixtures of Gaussians, with
$k$ replaced by $1/p_{\mathrm{min}}$ in all guarantees. We note, in
particular, that the separation requirement also increases to
$\sqrt{\log(1/p_{\mathrm{min}})}$.

We obtain the bipartition by running a greedy algorithm with a distance
function that satisfies the following two properties: (1) the distance between
most samples from the same component is small, and (2) there exist two
components $a$ and $b$ such that the distance between most samples from $a$
and $b$ is large. In turn, we extract this distance function from a polynomial
computed in terms of a pseudo-expectation of degree $O(\log k)$.

##### Colinear means

Suppose, in addition, that the mixture of Guassians is well-separated, i.e.,
the minimum mean separation satisfies $\min_{a\neq
b}\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})\rVert\gg\sqrt{\log k}$, and that the
unknown means $\mu_{1},...,\mu_{k}$ are colinear. Given a quasi-polynomial
number of samples from the mixture, Theorem 1.2 shows that it is possible to
compute in quasi-polynomial time a partition of the samples into $k$ clusters
such that the fraction of samples assigned to incorrect clusters is
polynomially small (in $k$).

In the theorem statement below, we assumed for simplicity that all eigenvalues
of $\Sigma$ and all eigenvalues of the covariance matrix of the mixture are
polynomially lower and upper bounded in $k$ and $d$. See Theorem 6.1 for the
full result.

######  Theorem 1.2.

Given a sample of size $n\geqslant(kd)^{O(\log k)}$ from a uniform mixture of
$k$ Gaussian components $N(\mu_{1},\Sigma),\ldots,N(\mu_{k},\Sigma)$ with
$\min_{a\neq b}\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})\rVert\gg\sqrt{\log k}$ and
with $\mu_{1},...,\mu_{k}$ colinear, there exists an algorithm that runs in
time $n^{O(\log k)}$ and returns with high probability a partition of $[n]$
into $k$ sets $C_{1},...,C_{k}$ such that, if true clustering of the samples
is $S_{1},...,S_{k}$, then there exists a permutation $\pi$ of $[k]$ such that

$1-\frac{1}{n}\sum_{i=1}^{k}|C_{i}\cap S_{\pi(i)}|\leqslant k^{-O(1)}.$

The same result holds in the case of non-uniform mixtures of Gaussians, with
$k$ replaced by $1/p_{\mathrm{min}}$ in all guarantees. We note, in
particular, that the separation requirement also increases to
$\sqrt{\log(1/p_{\mathrm{min}})}$.

Given such a clustering, we can also recover the means and the covariance of
the components using robust Gaussian estimation algorithms [DKK+19] or robust
moment estimation algorithms [KSS18]. For example, via [KSS18], we obtain a
multiplicative approximation to the covariance
$\|\hat{\Sigma}-\Sigma\|\leqslant k^{-O(1)}\|\Sigma\|$ and a "covariance-
aware" approximation to the means $\|\hat{\mu}_{i}-\mu_{i}\|\leqslant
k^{-O(1)}\|\Sigma\|^{1/2}$.

##  2 Techniques

We consider uniform222In this section we restrict ourselves for ease of
explanation to uniform mixtures. Our technical sections state all results for
non-uniform mixtures. mixtures of $k\geqslant 2$ well-separated Gaussian
components $N(\mu_{1},\Sigma),\ldots,N(\mu_{k},\Sigma)$ with unknown means
$\mu_{1},\ldots,\mu_{k}\in\mathbb R^{d}$ and unknown covariance
$\Sigma\in\mathbb R^{d\times d}$ (identical for all components). Here, we say
components are _well-separated_ 333 The term ”clusterable mixture” is
sometimes used in the literature to refer to mixtures with well-separated
components. if the maximum affinity444The affinity of two probability measures
is defined to be $1$ minus their statistical distance [Pol02]. (also called
overlap) between two distinct components is bounded by $1/k^{C}$ for a large
enough constant $C\geqslant 1$. For Gaussian components, this notion of well-
separatedness means

$\min_{a\neq
b}\mathopen{}\mathclose{{}\left\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})}\right\rVert\gg\sqrt{\log
k\,}\,.$ (2.1)

##### Distinguishing well-separated mixtures from (pure) Gaussians

Our algorithms are informed by investigating the parallel pancakes
construction underlying Statistical Query lower bounds for such mixtures. This
construction provides a mixture of $k$ well-separated Gaussian components that
appears to be exponentially hard to distinguish from the standard Gaussian
distribution $N(0,I_{d})$. In particular, this mixture matches the first
$\Omega(k)$ moments of $N(0,1)$.

The starting point of our algorithms is the following observation: In order
for a mixture with $k\geqslant 2$ well-separated Gaussian components to match
the first $t$ moments of $N(0,I_{d})$, the minimum mixing weight is
necessarily smaller than $2^{-\Omega(t)}$. In particular, if the mixture has
uniform mixing weights $\frac{1}{k}$, then always one of its first $O(\log k)$
moments distinguishes it from a standard Gaussian.

Underlying this observation is the following simple fact: _A distribution
uniform over $k$ real values can match no more than $O(\log k)$ moments of
$N(0,1)$._ To verify this fact, let $\bm{A}$ be a random variables uniformly
distributed over $k$ (not necessarily distinct) real values. Then, for all
even integers $s\leqslant t$, the ratio of the normalized order-$s$ and
order-$t$ moments of $\bm{A}$ is sandwiched in the following way,

$k^{-1/s}\leqslant\frac{\mathopen{}\mathclose{{}\left(\operatorname*{\mathbb{E}}\bm{A}^{s}}\right)^{1/s}}{\mathopen{}\mathclose{{}\left(\operatorname*{\mathbb{E}}\bm{A}^{t}}\right)^{1/t}}\leqslant
1\,.$ (2.2)

(This ratio is maximized if $\bm{A}$ is constant and minimized if
$\operatorname*{\mathbb{P}}\\{\bm{A}\neq 0\\}=1/k$.) In particular for
$s=\log_{2}k$, this ratio is lower bounded by $1/2$. On the other hand, for
$\bm{B}\sim N(0,1)$, the normalized moments satisfy
$\mathopen{}\mathclose{{}\left(\operatorname*{\mathbb{E}}\bm{B}^{r}}\right)^{1/r}=\Theta(r)^{1/2}$
and thus the ratio of normalized order-$s$ and order-$t$ moments is
$\Theta(s/t)^{1/2}$. In particular, for some choice $t=\Theta(s)$, this ratio
is smaller than $1/2$. It follows that for this choice of $s$ and $t$, the
ratios of normalized moments differ for $\bm{A}$ and for $\bm{B}$, which means
that either their order-$s$ or their order-$t$ moments differ.555 This proof
shows that in order to distinguish a uniform distribution over $k$ values from
$N(0,1)$, it is enough to compare two moments of order logarithmic in $k$
where the choice of orders depends only on $k$ but not on the particular
distribution.

This observation about uniform mixtures of $k\geqslant 2$ well-separated
Gaussian components raises two questions: (1) do the first $O(\log k)$ moments
also allow us to identify parameters of the mixture that are useful for
clustering (in addition to allowing us to distinguish the mixture from
$N(0,1)$), and (2) can we make make these results computationally efficient?

At a high level, we address question (1) by investigating ratios akin to Eq.
2.2 between multivariate polynomials of degree $\Theta(\log k)$ derived from
moments of the underlying mixture. To address question (2), we employ the
proofs-to-algorithms paradigm (cf. [BS14, RSS18]) and translate our arguments
to syntactic proofs captured by the sum-of-squares proof system. These proofs
then allow us to derive efficient algorithms (with running time $(kd)^{O(\log
k)}$ or $(kd)^{(\log k)^{O(1)}}$) in a black-box way.

##### From decision to search: separating directions and ratios of moments

In order to address question (1), we consider the goal of finding a direction
$v\in\mathbb R^{d}$ that may be useful for clustering in the sense that along
direction $v$, two of the components are significantly further apart than
their standard deviation in this direction. More formally, we say that $v$ is
a _separating direction_ for a mixture of $k$ Gaussian components with unknown
means $\mu_{1},\ldots,\mu_{k}\in\mathbb R^{d}$ and unknown covariance
$\Sigma\in\mathbb R^{d\times d}$ if there exist two means $\mu_{a}$ and
$\mu_{b}$ such that

$\lvert\langle\mu_{a}-\mu_{b},v\rangle\rvert\gg\sqrt{\log
k}\cdot\lVert\Sigma^{1/2}v\rVert\,.$ (2.3)

We note that this direction $v$ witnesses that the overlap of the components
$N(\mu_{a},\Sigma)$ and $N(\mu_{b},\Sigma)$ is $k^{-\omega(1)}$. Conversely,
whenever the overlap of two components is that small, there exists a vector
$v$ as above.

We aim to identify separating directions as solutions to inequalities between
the following kind of moment polynomials: For $r\in\mathbb N$, we denote the
degree-$2r$ _moment polynomial_ $p_{2r}\in\mathbb R_{2r}[v]$ by

$p_{2r}(v)\mathrel{\mathop{:}}=\operatorname*{\mathbb{E}}\langle\bm{y}-\bm{y}^{\prime},v\rangle^{2r}\,,$
(2.4)

where $\bm{y},\bm{y}^{\prime}$ are two independent random vectors identically
distributed according to a uniform mixture of $k$ Gaussian components
$N(\mu_{1},\Sigma),\ldots,N(\mu_{k},\Sigma)$.

Using the fact that $\bm{y}-\bm{y}^{\prime}$ can be expressed as sum of two
independent random vectors, one distributed uniformly over
$\\{\mu_{a}-\mu_{b}\\}_{a,b\in[k]}$ and one distributed according to
$N(0,2\Sigma)$, these polynomials turn out to admit the following kind of
approximation,

$p_{2r}(v)=\mathopen{}\mathclose{{}\left(k^{-2/r}\cdot\mathopen{}\mathclose{{}\left\lVert
Mv}\right\rVert_{2r}^{2}+\Theta(r)\cdot\mathopen{}\mathclose{{}\left\lVert
Av}\right\rVert_{2}^{2}}\right)^{r}\,.$ (2.5)

Here, $A=\sqrt{2}\cdot\Sigma^{1/2}$, $M\in\mathbb R^{k^{2}\times d}$ consists
of the differences of means $(\mu_{a}-\mu_{b})_{a,b\in[k]}\subseteq\mathbb
R^{d}$ as rows and $\Theta(r)$ hides a nonnegative function upper bounded
$O(r)$ and lower bounded by $\Omega(r/k^{2/r})$. (Since we will only consider
$r\geqslant\log k$, we have $k^{-2/r}\geqslant\Omega(1)$.) Note that the first
term $k^{-2}\cdot\mathopen{}\mathclose{{}\left\lVert
Mv}\right\rVert_{2r}^{2r}$ in (the binomial expansion of) Eq. 2.5 corresponds
to the order-$r$ moment of the uniform distribution over
$\\{\mu_{a}-\mu_{b}\\}_{a,b\in[k]}$ and the last term
$\Theta(r)^{r}\cdot\langle v,2\Sigma v\rangle^{r}$ to the order-$r$ moment of
$N(0,2\Sigma)$.

We claim that for an appropriate choice $s\leqslant t$ with
$s=\Theta(t)=\Theta(\log k)$, a direction $v$ is separating in the sense of
Eq. 2.3 if and only $p_{2s}(v)^{1/s}\gtrsim p_{2t}(v)^{1/t}$. (Note that by
convexity, $p_{2s}(v)^{1/s}\leqslant p_{2t}(v)^{1/t}$ holds for all directions
$v$.) Underlying this claim is the familiar fact that for all $r\geqslant\log
k$, the norm $\mathopen{}\mathclose{{}\left\lVert Mv}\right\rVert_{2r}$ equals
up to constant factors the maximum entry of $Mv$, i.e., $\max_{a\neq
b}\lvert\langle\mu_{a}-\mu_{b},v\rangle\rvert$.

Indeed, suppose that $v$ is a separating direction. Then, $p_{2s}(v)$
satisfies the lower bound,

$p_{2s}(v)^{1/s}\geqslant k^{-2/s}\cdot\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}\,.$ (2.6)

Since $s=\Theta(\log k)$, we have $p_{2s}(v)^{1/s}\gtrsim\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}$. At the same time, $p_{2t}(v)$
satisfies the upper bound,

$p_{2t}(v)^{1/t}\leqslant\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}+O(t)\cdot\lVert Av\rVert_{2}^{2}\,.$
(2.7)

Since $v$ is a separating direction and $t=\Theta(\log k)$, the upper bound is
dominated by the first term $\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}$. Taking together both bounds, it
follows that $p_{2s}(v)^{1/s}\gtrsim p_{2t}(v)^{1/t}$ for every separating
direction $v$.

Conversely, suppose that $p_{2s}(v)^{1/s}\gtrsim p_{2t}(v)^{1/t}$ and our goal
is to show that $v$ is a separating direction. We lower bound $p_{2t}(v)$
using the last term in the approximation Eq. 2.5 and apply the upper bound
from Eq. 2.7 to $p_{2s}(v)$. In this way, we obtain the inequality

$\Omega(t)\cdot\lVert Av\rVert_{2}^{2}\leqslant\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}+O(s)\cdot\lVert Av\rVert_{2}^{2}\,.$
(2.8)

By choosing $t$ to be a large enough constant multiplied by $s$, we can ensure
that the second term on the right-hand side is negligible. In this case, $v$
satisfies $\max_{a\neq
b}\langle\mu_{a}-\mu_{b},v\rangle^{2}\geqslant\Omega(t)\cdot\lVert
Av\rVert_{2}^{2}$ for $t=\Theta(\log k)$, which means that $v$ is a separating
direction.

##### Challenges toward efficient algorithms for clustering

Disregarding computational efficiency, the above characterization of
separating directions in terms of ratios of moment polynomials suggests the
following simple strategy for clustering uniform mixtures of Gaussian
components $N(\mu_{1},\Sigma),\dots,N(\mu_{k},\Sigma)$: we find an
$\varepsilon$-cover of all separating directions by brute-force searching for
an $\varepsilon$-cover of all solutions to an explicit polynomial system of
the form $\\{p_{2s}(v)=1,\leavevmode\nobreak\ p_{2t}(v)\leqslant O(1)^{t}\\}$.
Each separating direction gives us some information about what pairs of sample
points belong to different components. For large enough mean separation
$\min_{a\neq b}\lVert\Sigma^{-1/2}(\mu_{a}-\mu_{b})\rVert\gg\sqrt{\log k\,}$,
we can hope that by considering all such directions, we collect enough
information to be able to extract a clustering of the sample that is
approximately consistent with the components of the mixture.

This naive approach would require access only to moments of order $O(\log k)$
(which could be accurately estimated from a sample of size $d^{O(\log k)}$)
but the running time is exponentially large (due to brute-force searching for
solutions to a polynomial system).

A natural strategy to make this approach computationally efficient is the sum-
of-squares (sos) hierarchy of semidefinite programming relaxations for systems
of polynomial inequalities. Indeed, we can show that the above
characterization of separating directions is faithfully captured by the sos
proof system underlying the sos hierarchy. Unfortunately, it appears to be
challenging to carry out the rounding step in full generality, i.e.,
extracting from the sos hierarchy enough separating directions to separate all
pairs of components and obtain a complete clustering of the sample.666 In the
context of estimation problems, rounding procedures for sos hierarchies tend
to work well if there is a unique target solution (e.g., a planted sparse
vector in a random subspace) or if there is a small number of target solutions
(e.g., the components of a low-rank tensor). One could try to simplify the
structure of the set of separating directions, e.g., by focusing on ”extreme”
separating directions of the form $v=\Sigma^{-1}(\mu_{a}-\mu_{b})$.
Unfortunately, we do not know the same kind of characterization in terms of
polynomial inequalities for such a simplified set of separating directions.

However, we can show that using the sos hierarchy, it is possible to separate
at least _some_ pairs of components of the mixture by what we call a
separating polynomial. Furthermore, for the special case of well-separated
components with colinear means, we provide a more careful analysis and show
that in this case the sos hierarchy does offer enough information to extract a
complete clustering.

##### Efficiently computing a separating polynomial

As discussed above, we consider the goal of separating some pairs of
components of a mixture (as opposed to the stronger goal of separating all
pairs of components as would be required for a complete clustering). One way
to achieve this goal is by finding a separating direction in the sense of Eq.
2.3. In light of our previous characterization of separating directions, a
natural starting point is a sum-of-squares relaxation for a polynomial system
$\mathcal{A}$ of the form $\\{p_{2s}(v)=1,\leavevmode\nobreak\
p_{2t}(v)\leqslant O(1)^{t}\\}$ for appropriate $s\leqslant t$ satisfying
$s=\Theta(t)=\Theta(\log k)$.

Unfortunately, the structure of the set of separating directions does not
appear to be amenable to the usual kind of rounding techniques for sos
relaxations, and it appears to be challenging to extract a single separating
direction. To overcome this obstacle, we allow our rounding procedure to
output a more general object, called a _separating polynomial_ , that still
allows us to separate some pairs of mixture components.

Recall that a solution to a sos relaxation for a polynomial system
$\mathcal{A}$ can be interpreted as _pseudo-distribution_ $D$ that behaves in
certain ways like a distribution supported on vectors satisfying
$\mathcal{A}$. More concretely, the pseudo-distribution $D$ satsifies (in
expectation) all polynomial inequalities that can be derived syntactically
from $\mathcal{A}$ by a low-degree sum-of-squares proof (see Section 3,
especially Definition 3.1). The previously discussed characterization of
separating directions in terms of the polynomial system $\mathcal{A}$ turns
out to be captured by low-degree sum-of-squares proofs. Concretely, we can
derive from $\mathcal{A}$ via low-degree sum-of-squares proof the polynomial
inequality777 Here, we reuse the notation introduced in the context of Eq.
2.5. $\lVert Mv\rVert_{2s}^{2s}\geqslant(C\log k)^{s}\cdot\lVert
Av\rVert_{2}^{2s}$ (corresponding to Eq. 2.3). Here, $C\geqslant 1$ is an
absolute constant that we can choose as large as we like. Consequently, the
pseudo-distribution $D$ satisfies this inequality in expectation,
$\tilde{\operatorname*{\mathbb{E}}}_{D(v)}\lVert
Mv\rVert_{2s}^{2s}\geqslant(C\log
k)^{s}\cdot\tilde{\operatorname*{\mathbb{E}}}_{D(v)}\lVert
Av\rVert^{2s}_{2}\,.$ By linearity of (pseudo-)expectation, there exist
distinct components $a\neq b$ such that

$\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\langle\mu_{a}-\mu_{b},v\rangle^{2s}\geqslant
k^{-2}\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\lVert
Mv\rVert_{2s}^{2s}\geqslant(k^{-2/s}\cdot C\log
k)^{s}\cdot\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\lVert
Av\rVert^{2s}_{2}\,.$ (2.9)

We extract the following polynomial from this pseudo-distribution,

$q(u)\mathrel{\mathop{:}}=\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\langle
u,v\rangle^{2s}\,.$ (2.10)

By construction, $q(\mu_{a}-\mu_{b})$ equals the left-hand side of Eq. 2.9. At
the same time, letting $\bm{y},\bm{y}^{\prime}\sim N(\mu_{c},\Sigma)$ and
$\bm{w}\sim N(0,I_{d})$, we have

$\displaystyle\operatorname*{\mathbb{E}}q(\bm{y}-\bm{y}^{\prime})$
$\displaystyle=\operatorname*{\mathbb{E}}q(A\bm{w})$
$\displaystyle=\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\operatorname*{\mathbb{E}}\langle
v,A\bm{w}\rangle^{2s}$
$\displaystyle=(2s-1)!!\cdot\tilde{\operatorname*{\mathbb{E}}_{D(v)}}\lVert
A\bm{w}\rVert^{2s}_{2}$

Consequently, since $s=\Theta(\log k)$ and $(2s-1)!!\leqslant O(s)^{s}$,

$\frac{q(\mu_{a}-\mu_{b})}{\operatorname*{\mathbb{E}}q(A\bm{w})}\geqslant\frac{(k^{-2/s}\cdot
C\log k)^{s}}{(2s-1)!!}\geqslant\Omega(C)^{s}\,.$ (2.11)

For an appropriate choice of $C\geqslant 1$, the right-hand side above is at
least $8^{s}$. Since by convexity
$\operatorname*{\mathbb{E}}q(\mu_{a}-\mu_{b}+A\bm{w})\geqslant
q(\mu_{a}-\mu_{b})$, we obtain the following inequality,

$\frac{\operatorname*{\mathbb{E}}q(\mu_{a}-\mu_{b}+A\bm{w})}{\operatorname*{\mathbb{E}}q(A\bm{w})}\geqslant
8^{s}\,.$ (2.12)

This inequality shows that the polynomial $q(u)$ separates the components
$N(\mu_{a},\Sigma)$ and $N(\mu_{b},\Sigma)$ in the following sense: The
numerator of Eq. 2.12 is the typical value of $q(\bm{y}-\bm{y}^{\prime})$ for
$\bm{y}\sim N(\mu_{a},\Sigma)$ and $\bm{y}^{\prime}\sim N(\mu_{b},\Sigma)$.
The denominator of Eq. 2.12 is the typical value of
$q(\bm{y}-\bm{y}^{\prime})$ for $\bm{y},\bm{y}^{\prime}\sim N(\mu_{c},\Sigma)$
and all $c\in[k]$. Equation 2.12 asserts that the gap between these values is
at least $8^{s}$.

The polynomial $q(u)$ can be used to compute a bipartition of the sample that
separates at least one pair of components. Note that
$q(u)^{1/2s}=(\tilde{\operatorname*{\mathbb{E}}}_{D(v)}\langle
u,v\rangle^{2s})^{1/2s}$ satisfies the triangle inequality (see Lemma 4.5 in
[BS14]). Then we can define the distance function $d(x,y)=q(x-y)^{1/2s}$ and
use it in a greedy algorithm in order to obtain the bipartition. This argument
is presented in Lemma A.1.

##### Efficiently computing a clustering for colinear means

For the case that the means are colinear, we consider a strengthening of our
previous approach. Instead of trying to solve a polynomial system of the form
$\\{p_{2s}(v)=1,\leavevmode\nobreak\ p_{2t}(v)\leqslant O(1)^{t}\\}$, we aim
to solve the following related optimization problem:

$\text{minimize}\quad\frac{p_{2t}(v)^{1/t}}{p_{2s}(v)^{1/s}}\quad\text{subject
to}\quad v\in\mathbb R^{d}.$ (2.13)

Algorithmically, we again employ an appropriate sum-of-squares formulation.

To simplify some of our arguments, it is useful to preprocess the mixture and
bring $\bm{y}-\bm{y}^{\prime}$ in isotropic position so that
$\tfrac{1}{k^{2}}\sum_{a,b=1}^{k}(\mu_{a}-\mu_{b})(\mu_{a}-\mu_{b}){}^{\mkern-1.5mu\mathsf{T}}+2\Sigma=I_{d}$.
(Here, $\bm{y},\bm{y}^{\prime}$ are two independent random vectors distributed
according to the mixture.) For every vector $v$, we denote by $v^{\parallel}$
its orthogonal projection into the span of $\\{\mu_{a}-\mu_{b}\\}_{a,b\in[k]}$
and by $v^{\perp}=v-v^{\parallel}$ its projection into the orthogonal
complement.

Every optimizer $v$ of Eq. 2.13 necessarily satisfies,

$\frac{p_{2t}(v)^{1/t}}{p_{2s}(v)^{1/s}}\leqslant\frac{p_{2t}(v^{\parallel})^{1/t}}{p_{2s}(v^{\parallel})^{1/s}}\leqslant
O(1)\,.$ (2.14)

Here, the upper bound $O(1)$ hides an absolute constant whenever we have well-
separated components and $\log k\leqslant s\leqslant t$. The argument for this
upper bound is similar to our discussion for the characterization of
separating directions.

We can also use the decomposition $v=v^{\parallel}+v^{\perp}$ for our previous
approximation Eq. 2.5 of moment polynomials,

$p_{2r}(v)=\mathopen{}\mathclose{{}\left(k^{-2/r}\cdot\mathopen{}\mathclose{{}\left\lVert
Mv^{\parallel}}\right\rVert_{2r}^{2}+\Theta(r)\cdot\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left\lVert
Av^{\parallel}}\right\rVert_{2}^{2}+\mathopen{}\mathclose{{}\left\lVert
v^{\perp}}\right\rVert_{2}^{2}}\right)}\right)^{r}\,.$ (2.15)

Here, we use that after bringing $\bm{y}-\bm{y}^{\prime}$ in isotropic
position, the covariance $\Sigma$ acts as identity orthogonal to the span of
$\\{\mu_{a}-\mu_{b}\\}_{a,b\in[k]}$. In particular,
$Av=Av^{\parallel}+v^{\perp}$ and $\lVert Av\rVert^{2}_{2}=\lVert
Av^{\parallel}\rVert_{2}^{2}+\lVert v^{\perp}\rVert_{2}^{2}$.

An immediate consequence of Eq. 2.15 is the following representation of the
ratio we seek to minimize,

$\frac{p_{2t}(v)^{1/t}}{p_{2s}(v)^{1/s}}=\frac{p_{2t}(v^{\parallel})^{1/t}+\Theta(t)\cdot\lVert
v^{\perp}\rVert_{2}^{2}\pm\Theta(t)\cdot\lVert
Av^{\parallel}\rVert_{2}^{2}}{p_{2s}(v^{\parallel})^{1/s}+\Theta(s)\cdot\lVert
v^{\perp}\rVert_{2}^{2}\pm\Theta(s)\cdot\lVert
Av^{\parallel}\rVert_{2}^{2}}\,.$ (2.16)

We claim that for an appropriate choice of $s$ and $t$ Eq. 2.16 and Eq. 2.14
together imply that $\lVert v^{{\perp}}\rVert\lesssim\lVert
Av^{\parallel}\rVert$. Indeed, for the sake of a contradiction, suppose
$\lVert v^{{\perp}}\rVert\gg\lVert Av^{\parallel}\rVert$ so that the terms
involving $\lVert Av^{\parallel}\rVert$ in Eq. 2.16 are negligible. But then,
if we choose $t$ as $s$ times a sufficiently larger constant factor, the
remaining ratio $\frac{p_{2t}(v^{\parallel})^{1/t}+\Theta(t)\cdot\lVert
v^{\perp}\rVert_{2}^{2}}{p_{2s}(v^{\parallel})^{1/s}+\Theta(s)\cdot\lVert
v^{\perp}\rVert_{2}^{2}}$ is strictly bigger than
$\frac{p_{2t}(v^{\parallel})^{1/t}}{p_{2s}(v^{\parallel})^{1/s}}$, which
contradicts our optimality condition Eq. 2.14. (For this argument, we are also
using the previous upper bound
$\frac{p_{2t}(v^{\parallel})^{1/t}}{p_{2s}(v^{\parallel})^{1/s}}\leqslant
O(1)$ from Eq. 2.14.)

It turns out that in order to compute a clustering for colinear means, it
suffices to find a vector satisfying $\lVert v^{{\perp}}\rVert\lesssim\lVert
Av^{\parallel}\rVert$.

We note that the algorithm we present in Sections 5 and 6 to find such a
direction $v$ follows a somewhat different strategy and minimizes ratios of
the form ${\lVert v\rVert_{2}^{2}}/{p_{2s}(v)^{1/s}}$ or
${p_{2t}(v)^{1/t}}/{\lVert v\rVert_{2}^{2}}$ via appropriate sum-of-squares
formulations.

##  3 Preliminaries

In this section we introduce sum-of-squares proofs and their duals, pseudo-
distributions and pseudo-expectations.

##### Sum-of-squares proofs.

######  Definition 3.1 (Sum-of-squares proofs).

Let $p(x)$ and $q_{1}(x),...,q_{m}(x)$ be polynomials over
$x\in\mathbb{R}^{n}$ and let $\mathcal{A}=\\{q_{1}(x)\geqslant
0,...,q_{m}(x)\geqslant 0\\}$ be a system of polynomial inequalities. A sum-
of-squares proof of degree $t$ that $p(x)\geqslant 0$ under $\mathcal{A}$ is
an identity of the form

$p(x)=\sum_{S\subseteq[m]}\mathopen{}\mathclose{{}\left(\sum_{i=1}^{m_{S}}r_{S,i}(x)^{2}}\right)\prod_{j\in
S}q_{j}(x)$ (3.1)

for polynomials $r_{S,i}(x)$, such that
$\max_{S,i}\operatorname{deg}(r_{S,i}(x)^{2}\prod_{j\in S}q_{j}(x))\leqslant
t$.

If there exists a sum-of-squares proof of degree $t$ that $p(x)\geqslant 0$
under $\mathcal{A}$, we write $\mathcal{A}\sststile{t}{x}p(x)\geqslant 0$. We
also use the notation $\mathcal{A}\sststile{t}{x}p(x)\geqslant q(x)$ if
$\mathcal{A}\sststile{t}{x}p(x)-q(x)\geqslant 0$ and
$\mathcal{A}\sststile{t}{x}p(x)\leqslant q(x)$ if
$\mathcal{A}\sststile{t}{x}q(x)-p(x)\geqslant 0$. If $\mathcal{A}=\emptyset$,
we omit it altogether and write $\sststile{t}{x}p(x)\geqslant 0$. We also
sometimes omit $\mathcal{A}$ if it is clear from context what axioms are
assumed. We note that if $\mathcal{A}\sststile{t}{x}p(x)\geqslant q(x)$ and
$\mathcal{A}\sststile{t}{s}q(x)\geqslant r(x)$, then
$\mathcal{A}\sststile{t}{x}p(x)\geqslant r(x)$, which allows writing chains of
inequalities of the form $\mathcal{A}\sststile{t}{x}p(x)\geqslant
s(x)\geqslant r(x)$.

##### Pseudo-distributions and pseudo-expectations.

We begin by defining pseudo-distributions and pseudo-expectations.

######  Definition 3.2 (Pseudo-distributions).

A pseudo-distribution $D$ of degree $t$ is a function from $\mathbb{R}^{n}$ to
$\mathbb{R}$ with finite support such that
$\sum_{x\in\operatorname{supp}(D)}D(x)=1$ and
$\sum_{x\in\operatorname{supp}(D)}D(x)p(x)^{2}\geqslant 0$ for all polynomials
$p(x)$ with $\operatorname{deg}(p(x)^{2})\leqslant t$.

######  Definition 3.3 (Pseudo-expectations).

Given a pseudo-distribution $D$ of degree $t$, the associated pseudo-
expectation $\tilde{\mathbb{E}}_{D(x)}$ is defined by
$\tilde{\mathbb{E}}_{D(x)}f(x)=\sum_{x\in\operatorname{supp}(D)}D(x)f(x)$ for
a function $f(x)$.

We now define the notion of a pseudo-distribution that satisfies a set of
polynomial inequalities.

######  Definition 3.4 (Constrained pseudo-distributions).

A pseudo-distribution $D$ of degree $t$ satisfies the set of polynomial
inequalities $\mathcal{A}=\\{q_{1}(x)\geqslant 0,...,q_{m}(x)\geqslant 0\\}$
if, for all $S\subseteq[m]$, $\tilde{\mathbb{E}}_{D(x)}r(x)^{2}\prod_{j\in
S}q_{j}(x)\geqslant 0$ for all polynomials $r(x)$ such that
$\operatorname{deg}(r(x)^{2}\prod_{j\in S}q_{j}(x))\leqslant t$.

$D$ approximately satisfies $\mathcal{A}$ up to error $\eta$ if, under the
same conditions as in the previous case,
$\tilde{\mathbb{E}}_{D(x)}r(x)^{2}\prod_{j\in
S}q_{j}(x)\geqslant-\eta\|r(x)^{2}\|_{2}\prod_{j\in S}\|q_{j}(x)\|_{2}$, where
$\|p(x)\|_{2}$ denotes the $2$-norm of the vector of coefficients of the
polynomial $p(x)$.

The connection between pseudo-distributions and sum-of-squares proofs is made
in Fact 3.5, which shows that if a pseudo-distribution satisfies a set of
polynomial inequalities, then it also satisfies any other polynomial
inequalities derived from this set through sum-of-squares proofs.

######  Fact 3.5.

If $D$ is a pseudo-distribution of degree $t$ that satisfies $\mathcal{A}$ and
if $\mathcal{A}\sststile{s}{x}p(x)\geqslant 0$, then
${\tilde{\mathbb{E}}_{D(x)}r(x)^{2}p(x)\geqslant 0}$ for all polynomials
$r(x)$ such that $\operatorname{deg}(r(x)^{2}p(x))\leqslant t$. If $D$
approximately satisfies $\mathcal{A}$ up to error $\eta$, then, under the same
conditions as in the previous case,
${\tilde{\mathbb{E}}_{D(x)}r(x)^{2}p(x)\geqslant-\eta\|r(x)^{2}\|_{2}\|p(x)\|_{2}}$.

Finally, Fact 3.6 shows that there exists an algorithm with time complexity
$(n+m)^{O(t)}$ to compute a pseudo-distribution of degree $t$ that
approximately satisfies $\mathcal{A}$ up to error $2^{-n^{\Theta(t)}}$.

######  Fact 3.6.

For $x\in\mathbb{R}^{n}$, if $\mathcal{A}=\\{q_{1}(x)\geqslant
0,...,q_{m}(x)\geqslant 0\\}$ is feasible and explicitly bounded888Explicit
boundedness means that $\mathcal{A}$ contains a constraint of the form
$x_{1}^{2}+...+x_{n}^{2}\leqslant B$. In our applications it is possible to
add such a constraint with $B$ large enough such that the constraint is always
satisfied by the intended solution., then there exists an algorithm that runs
in time $(n+m)^{O(t)}$ and computes a pseudo-distribution of degree $t$ that
approximately satisfies $\mathcal{A}$ up to error $2^{-n^{\Theta(t)}}$999In
our applications this error is negligible..

##  4 Separating polynomial

##### Setting.

We consider a mixture of $k$ Gaussian distributions $N(\mu_{i},\Sigma)$ with
mixing weights $p_{i}$ for $i=1,...,k$, where $\mu_{i}\in\mathbb{R}^{d}$,
$\Sigma\in\mathbb{R}^{d\times d}$ is positive definite, and $p_{i}\geqslant 0$
and $\sum_{i=1}^{k}p_{i}=1$. Let $p_{\mathrm{min}}=\min_{i}p_{i}$.

The distribution satisfies mean separation for at least one pair of means: for
some $C_{sep}>0$, there exist $a,b\in[k]$ such that

$\mathopen{}\mathclose{{}\left\|\Sigma^{-1/2}(\mu_{a}-\mu_{b})}\right\|^{2}\geqslant
C_{sep}\log p_{\mathrm{min}}^{-1}.$

######  Theorem 4.1 (Separating polynomial algorithm).

Consider the Gaussian mixture model defined above, with $C_{sep}$ larger than
some universal constant. Let $n_{0}=(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}$. Given a sample of size $n\geqslant n_{0}$ from the
mixture, there exists an algorithm that computes in time $n\cdot d^{O(\log
p_{\mathrm{min}}^{-1})}$ a $d$-variate degree-$O(\log p_{\mathrm{min}}^{-1})$
polynomial $q$ such that with high probability the following two properties
hold. Let $s=\lceil\log p_{\mathrm{min}}^{-1}\rceil$. Then:

  * •

There exist distinct $a,b\in[k]$ such that the independent random vectors
$\bm{y}\sim N(\mu_{a},\Sigma)$ and $\bm{y}^{\prime}\sim N(\mu_{b},\Sigma)$
satisfy

$\operatorname*{\mathbb{P}}\mathopen{}\mathclose{{}\left\\{q(\bm{y}-\bm{y}^{\prime})\geqslant\frac{1}{8^{s}}}\right\\}\geqslant
0.99999.\,$

  * •

For all $a\in[k]$, the independent random vectors $\bm{y},\bm{y}^{\prime}\sim
N(\mu_{a},\Sigma)$ satisfy

$\operatorname*{\mathbb{P}}\mathopen{}\mathclose{{}\left\\{q(\bm{y}-\bm{y}^{\prime})\leqslant\frac{1}{80^{s}}}\right\\}\geqslant
0.99999.\,$

We introduce some further notation for this section. Let the random variable
$\bm{z}\in\mathbb{R}^{d}$ be distributed according to the difference of two
independent samples from the mixture. Then $\bm{z}$ is distributed according
to a mixture of Gaussians $N(\mu_{i}-\mu_{j},2\Sigma)$ with mixing weights
$p_{i}p_{j}$ for all $i,j\in[k]$. Let $\Sigma_{z}=2\Sigma$, let $\bm{\mu}_{z}$
be $\mu_{i}-\mu_{j}$ with probability $p_{i}p_{j}$, and let $\bm{w}_{z}\sim
N(0,\Sigma_{z})$. Then we also have that $\bm{z}=\bm{\mu}_{z}+\bm{w}_{z}$,
with $\bm{\mu}_{z}$ and $\bm{w}_{z}$ independent of each other.

###  4.1 Exact moment results

The main ingredient of the algorithm is Lemma 4.2, stated below. This lemma
shows that, given a pseudo-expectation that satisfies the moment lower bound
$\mathbb{E}\langle\bm{y},v\rangle^{2s}\geqslant c^{s}$ and the moment upper
bound $\mathbb{E}\langle\bm{y},v\rangle^{2t}\leqslant C^{t}\\}$ for $s\ll t$,
it is possible to construct a separating polynomial. Note that the constraints
that the pseudo-expectation satisfies are expressed in terms of exact moments
of the distribution, to which we do not have access. Finite sample
considerations are discussed starting with Section 4.2.

######  Lemma 4.2 (Separating polynomial from pseudo-expectation).

Let $c>0$ and $C\geqslant 0$. Let $s\geqslant 1$ and $t\geqslant 50000Cs/c$
integers. Given a pseudo-expectation $\tilde{\mathbb{E}}$ of degree at least
$2t$ over a variable $v\in\mathbb{R}^{d}$ that satisfies
$\\{\mathbb{E}\langle\bm{z},v\rangle^{2s}\geqslant
c^{s},\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant C^{t}\\}$, let
$q(u)=\langle\tilde{\mathbb{E}}v^{\otimes 2s},u^{\otimes 2s}\rangle$. Then:

  * •

There exist distinct $a,b\in[k]$ such that the independent random vectors
$\bm{y}\sim N(\mu_{a},\Sigma)$ and $\bm{y}^{\prime}\sim N(\mu_{b},\Sigma)$
satisfy

$\operatorname*{\mathbb{P}}\mathopen{}\mathclose{{}\left\\{q(\bm{y}-\bm{y}^{\prime})\geqslant\frac{1}{2}\mathopen{}\mathclose{{}\left(\frac{c}{4}}\right)^{s}}\right\\}\geqslant
0.99999\,.$

  * •

For all $a\in[k]$, the independent random vectors $\bm{y},\bm{y}^{\prime}\sim
N(\mu_{a},\Sigma)$ satisfy

$\operatorname*{\mathbb{P}}\mathopen{}\mathclose{{}\left\\{q(\bm{y}-\bm{y}^{\prime})\leqslant
11\mathopen{}\mathclose{{}\left(\frac{4Cs}{t}}\right)^{s}}\right\\}\geqslant
0.99999\,.$

In what follows, we prove a number of supporting lemmas, after which we prove
Lemma 4.2. Then, we state and prove Lemma 4.7, which shows that there exists a
vector $v\in\mathbb{R}^{d}$ which satisfies the constraints required by Lemma
4.2.

We proceed with the supporting lemmas. Lemma 4.3 and Lemma 4.4 give sum-of-
squares bounds on the moments of the mixture.

######  Lemma 4.3 (Moment upper bound).

For $t\geqslant 1$ integer,

$\sststile{2t}{}\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant
2^{2t-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+2^{2t-1}(v^{\top}\Sigma_{z}v)^{t}t^{t}.$

###### Proof.

$\displaystyle\sststile{2t}{}\mathbb{E}\langle\bm{z},v\rangle^{2t}$
$\displaystyle=\mathbb{E}(\langle\bm{\mu}_{z},v\rangle+\langle\bm{w}_{z},v\rangle)^{2t}\stackrel{{\scriptstyle(1)}}{{\leqslant}}2^{2t-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+2^{2t-1}\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t}$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{\leqslant}}2^{2t-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+2^{2t-1}(v^{\top}\Sigma_{z}v)^{t}t^{t}$

where in (1) we used Lemma A.3 and in (2) we used that
$\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t}=(v^{\top}\Sigma_{z}v)^{t}(2t-1)!!\leqslant(v^{\top}\Sigma_{z}v)^{t}t^{t}$.
∎

######  Lemma 4.4 (Moment lower bound).

For $t\geqslant 1$ integer,

$\sststile{2t}{}\mathbb{E}\langle\bm{z},v\rangle^{2t}\geqslant\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+(v^{\top}\Sigma_{z}v)^{t}\frac{t^{t}}{2^{t}}.$

###### Proof.

$\displaystyle\sststile{2t}{}\mathbb{E}\langle\bm{z},v\rangle^{2t}$
$\displaystyle=\mathbb{E}(\langle\bm{\mu}_{z},v\rangle+\langle\bm{w}_{z},v\rangle)^{2t}=\sum_{j=0}^{2t}\binom{2t}{j}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{j}\langle\bm{w}_{z},v\rangle^{2t-j}$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{=}}\sum_{j=0}^{2t}\binom{2t}{j}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{j}\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t-j}\stackrel{{\scriptstyle(2)}}{{=}}\sum_{s=0}^{t}\binom{2t}{2s}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t-2s}$
$\displaystyle\geqslant\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t}\stackrel{{\scriptstyle(3)}}{{\geqslant}}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+(v^{\top}\Sigma_{z}v)^{t}\frac{t^{t}}{2^{t}}$

where in (1) we used that $\bm{\mu}_{z}$ and $\bm{w}_{z}$ are independent, in
(2) we used that $\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t-j}=0$ for $2t-j$
odd, and in (3) we used that
$\mathbb{E}\langle\bm{w}_{z},v\rangle^{2t}=(v^{\top}\Sigma_{z}v)^{t}(2t-1)!!\geqslant(v^{\top}\Sigma_{z}v)^{t}\frac{t^{t}}{2^{t}}$.
∎

Going forward, Lemma 4.5 proves that, if the moments of $\bm{z}$ are small in
direction $v$, then the variance of the components of the mixture is also
small in direction $v$. Given in addition an upper bound on the moments of
$\bm{z}$ in direction $v$ for a sufficiently large moment, Lemma 4.6 proves
that the contribution of the means in direction $v$ is large.

######  Lemma 4.5.

Let $C\geqslant 0$. For $t\geqslant 1$ integer,

$\mathopen{}\mathclose{{}\left\\{\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant
C^{t}}\right\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{v^{\top}\Sigma_{z}v\leqslant\frac{2C}{t}}\right\\}.$

###### Proof.

Substitute the lower bound of Lemma 4.4 into the axiom:

$\sststile{2t}{}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+(v^{\top}\Sigma_{z}v)^{t}t^{t}/2^{t}\leqslant
C^{t}.$

Use that $\sststile{2t}{}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}\geqslant
0$ to drop the first term and then divide by $t^{t}/2^{t}$. This proves that
$\sststile{2t}{}(v^{\top}\Sigma_{z}v)^{t}\leqslant 2^{t}C^{t}/t^{t}.$ Finally,
by Lemma A.7, this implies that $\sststile{2t}{}v^{\top}\Sigma_{z}v\leqslant
2C/t$. ∎

######  Lemma 4.6.

Let $c>0$ and $C\geqslant 0$. For $s\geqslant 1$ and $t\geqslant 16Cs/c$
integers,

$\mathopen{}\mathclose{{}\left\\{\mathbb{E}\langle\bm{z},v\rangle^{2s}\geqslant
c^{s},\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant
C^{t}}\right\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}\geqslant\mathopen{}\mathclose{{}\left(\frac{c}{4}}\right)^{s}}\right\\}.$

###### Proof.

Substitute the upper bound of Lemma 4.3 into the first axiom:

$\sststile{2s}{}2^{2s-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}+2^{2s-1}(v^{\top}\Sigma_{z}v)^{s}s^{s}\geqslant
c^{s}.$

Use, by Lemma 4.5 and Lemma A.6, that
$\sststile{2t}{}(v^{\top}\Sigma_{z}v)^{s}\leqslant\mathopen{}\mathclose{{}\left(2C/t}\right)^{s}$:

$\sststile{2s}{}2^{2s-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}+\mathopen{}\mathclose{{}\left(8Cs/t}\right)^{s}\geqslant
c^{s}.$

Then use that $t\geqslant 16Cs/c$ to obtain that
$\sststile{2s}{}2^{2s-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}\geqslant
c^{s}/2$. Finally, divide by $2^{2s-1}$ to obtain that
$\sststile{2s}{}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}\geqslant(c/4)^{s}$.
∎

Now we prove Lemma 4.2.

###### Proof of Lemma 4.2.

Let $\bm{w}\sim N(0,I_{d})$. Note that, for $\bm{y}\sim N(\mu_{a},\Sigma)$ and
$\bm{y}^{\prime}\sim N(\mu_{b},\Sigma)$ we have that
$\bm{y}-\bm{y}^{\prime}\sim\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w}$, so
$q(\bm{y}-\bm{y}^{\prime})=q(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})$.
Similarly, for $\bm{y},\bm{y}^{\prime}\sim N(\mu_{a},\Sigma)$ we have that
$\bm{y}-\bm{y}^{\prime}\sim\Sigma_{z}^{1/2}\bm{w}$, so
$q(\bm{y}-\bm{y}^{\prime})=q(\Sigma_{z}^{1/2}\bm{w})$. Therefore, we want to
show (1) that there exist distinct $a,b\in[k]$ such that
$q(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})$ is large and (2) that
$q(\Sigma_{z}^{1/2}\bm{w})$ is small.

By Lemma 4.6, we have
$\tilde{\mathbb{E}}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}\geqslant\mathopen{}\mathclose{{}\left(c/4}\right)^{s}$.
By linearity, we also have
$\mathbb{E}\tilde{\mathbb{E}}\langle\bm{\mu}_{z},v\rangle^{2s}\geqslant\mathopen{}\mathclose{{}\left(c/4}\right)^{s}$,
so $\mathbb{E}q(\bm{\mu}_{z})\geqslant(c/4)^{s}$. Therefore there exists some
$\mu_{z}$ in the support of $\bm{\mu}_{z}$ such that
$q(\mu_{z})\geqslant(c/4)^{s}$. Therefore, there exist $a,b\in[k]$ such that
$q(\mu_{a}-\mu_{b})\geqslant(c/4)^{s}$. Furthermore, $a$ and $b$ are distinct,
because otherwise $q(\mu_{a}-\mu_{b})=0$.

We attempt to lower bound $q(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})$:

$\displaystyle q(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})$
$\displaystyle=\langle\tilde{\mathbb{E}}v^{\otimes
2s},(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})^{\otimes 2s}\rangle$
$\displaystyle=\tilde{\mathbb{E}}\langle
v,\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w}\rangle^{2s}$
$\displaystyle\geqslant\frac{1}{2^{2s-1}}\tilde{\mathbb{E}}\langle
v,\mu_{a}-\mu_{b}\rangle^{2s}-\tilde{\mathbb{E}}\langle
v,\Sigma_{z}^{1/2}\bm{w}\rangle^{2s}$
$\displaystyle=\frac{1}{2^{2s-1}}q(\mu_{a}-\mu_{b})-q(\Sigma_{z}^{1/2}\bm{w})$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(c/16}\right)^{s}-q(\Sigma_{z}^{1/2}\bm{w}).$

We want to show that $q(\Sigma_{z}^{1/2}\bm{w})$ is small with high
probability. We start by analyzing the mean and second moment of
$q(\Sigma_{z}^{1/2}\bm{w})$. For $\ell\in\\{1,2\\}$, we have

$\displaystyle\mathbb{E}q(\Sigma_{z}^{1/2}\bm{w})^{\ell}$
$\displaystyle=\mathbb{E}\mathopen{}\mathclose{{}\left(\langle\tilde{\mathbb{E}}v^{\otimes
2s},(\Sigma_{z}^{1/2}\bm{w})^{\otimes
2s}\rangle}\right)^{\ell}=\mathbb{E}\mathopen{}\mathclose{{}\left(\tilde{\mathbb{E}}\langle
v,\Sigma_{z}^{1/2}\bm{w}\rangle^{2s}}\right)^{\ell}$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{\leqslant}}\mathbb{E}\tilde{\mathbb{E}}\langle
v,\Sigma_{z}^{1/2}\bm{w}\rangle^{2s\ell}=\tilde{\mathbb{E}}\mathbb{E}\langle
v,\Sigma_{z}^{1/2}\bm{w}\rangle^{2s\ell}$
$\displaystyle=\tilde{\mathbb{E}}\mathbb{E}\langle\Sigma_{z}^{1/2}v,\bm{w}\rangle^{2s\ell}\leqslant(s\ell)^{s\ell}\tilde{\mathbb{E}}(v^{\top}\Sigma_{z}v)^{s\ell}$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{\leqslant}}\mathopen{}\mathclose{{}\left(2CS\ell/t}\right)^{s\ell},$

where in (1) for $\ell=2$ we used Lemma A.12 and in (2) we used that, by Lemma
4.5 and Lemma A.6, $\tilde{\mathbb{E}}(v^{\top}\Sigma_{z}v)^{s\ell}\leqslant
O\mathopen{}\mathclose{{}\left(2C/t}\right)^{s\ell}$. Then
$\mathbb{E}q(\Sigma_{z}{1/2}\bm{w})\leqslant\mathopen{}\mathclose{{}\left(2Cs/t}\right)^{s}$
and
$\mathbb{E}q(\Sigma_{z}^{1/2}\bm{w})^{2}\leqslant\mathopen{}\mathclose{{}\left(4Cs/t}\right)^{2s}$.
Therefore, by Chebyshev’s inequality, with probability $0.99999$,

$q(\Sigma_{z}^{1/2}\bm{w})\leqslant\mathopen{}\mathclose{{}\left(2Cs/t}\right)^{s}+\sqrt{100000}\mathopen{}\mathclose{{}\left(4Cs/t}\right)^{s}\leqslant
320\mathopen{}\mathclose{{}\left(4Cs/t}\right)^{s}.$

In this case, we also have

$q(\mu_{a}-\mu_{b}+\Sigma_{z}^{1/2}\bm{w})\geqslant\mathopen{}\mathclose{{}\left(c/16}\right)^{s}-320\mathopen{}\mathclose{{}\left(4Cs/t}\right)^{s}\geqslant\mathopen{}\mathclose{{}\left(c/16}\right)^{s}/2,$

where in the last inequality we used that $t\geqslant 50000Cs/c$. This
concludes the proof. ∎

######  Lemma 4.7 (Existence of vector that satisfies moment contraints).

Let $s,t\geqslant\lceil\log p_{\mathrm{min}}^{-1}\rceil$ integers. If
$t\leqslant\max_{i,j}\|\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j})\|^{2}$, there exists
some $v\in\mathbb{R}^{d}$ that satisfies
$\mathbb{E}\langle\bm{z},v\rangle^{2s}=1$ and
$\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant 30^{t}$. Furthermore,
$\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant 8$.

###### Proof.

Let $(a,b)=\arg\max_{(i,j)}\|\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j})\|$ and let
$v=\Sigma_{z}^{-1}(\mu_{a}-\mu_{b})$. The vector for which we will guarantee
the stated properties is
$v^{*}=\frac{v}{(\mathbb{E}\langle\bm{z},v\rangle^{2s})^{1/2s}}$.

We begin by proving that
$\max_{i,j}\langle\mu_{i}-\mu_{j},v\rangle^{2}=\langle\mu_{a}-\mu_{b},v\rangle^{2}$,
which will be used later. We have

$\displaystyle\max_{i,j}\langle\mu_{i}-\mu_{j},v\rangle^{2}$
$\displaystyle=\max_{i,j}\langle\mu_{i}-\mu_{j},\Sigma_{z}^{-1}(\mu_{a}-\mu_{b})\rangle^{2}$
$\displaystyle=\max_{i,j}\langle\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j}),\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\rangle^{2}$
$\displaystyle\leqslant\|\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j})\|^{2}\cdot\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2}$
$\displaystyle\leqslant\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4}$
$\displaystyle=\langle\mu_{a}-\mu_{b},v\rangle^{2}.$

We also have for the variance in direction $v$ that

$(v^{\top}\Sigma_{z}v)^{t}=((\mu_{a}-\mu_{b})^{\top}\Sigma_{z}^{-1}(\mu_{a}-\mu_{b}))^{t}=\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2t}.$

We now derive upper bounds for the $2t$ moments of $\langle\bm{z},v\rangle$ in
direction $v$ and lower bounds for the $2s$ moments in direction $v$. Recall
that we assume
$t\leqslant\max_{i,j}\|\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j})\|^{2}$. For the
upper bound, using Lemma 4.3 we have

$\displaystyle\mathbb{E}\langle\bm{z},v\rangle^{2t}$ $\displaystyle\leqslant
2^{2t-1}\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2t}+2^{2t-1}(v^{\top}\Sigma_{z}v)^{t}t^{t}$
$\displaystyle\leqslant
2^{2t-1}\max_{i,j}\langle\mu_{i}-\mu_{j},v\rangle^{2t}+2^{2t-1}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2t}t^{t}$
$\displaystyle=2^{2t-1}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}+2^{2t-1}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2t}t^{t}$
$\displaystyle\leqslant
2^{2t-1}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}+2^{2t-1}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}$
$\displaystyle\leqslant 2^{2t}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}.$

For the lower bound, using Lemma 4.4 we have

$\displaystyle\mathbb{E}\langle\bm{z},v\rangle^{2s}$
$\displaystyle\geqslant\mathbb{E}\langle\bm{\mu}_{z},v\rangle^{2s}+(v^{\top}\Sigma_{z}v)^{s}\frac{s^{s}}{2^{s}}$
$\displaystyle\geqslant
p_{\mathrm{min}}^{2}\max_{i,j}\langle\mu_{i}-\mu_{j},v\rangle^{2s}+\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2s}\frac{s^{s}}{2^{s}}$
$\displaystyle=p_{\mathrm{min}}^{2}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4s}+\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2s}\frac{s^{s}}{2^{s}}$
$\displaystyle\geqslant
p_{\mathrm{min}}^{2}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4s}.$

Recall that $v^{*}=\frac{v}{(\mathbb{E}\langle\bm{z},v\rangle^{2s})^{1/2s}}$.
Clearly, $\mathbb{E}\langle\bm{z},v^{*}\rangle^{2s}=1$. Furthermore,

$\displaystyle\mathbb{E}\langle\bm{z},v^{*}\rangle^{2t}$
$\displaystyle=\frac{\mathbb{E}\langle\bm{z},v\rangle^{2t}}{(\mathbb{E}\langle\bm{z},v\rangle^{2s})^{t/s}}\leqslant\frac{2^{2t}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}}{p_{\mathrm{min}}^{2t/s}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4t}}=\mathopen{}\mathclose{{}\left(\frac{4}{p_{\mathrm{min}}^{2/s}}}\right)^{t}\leqslant\mathopen{}\mathclose{{}\left(4e^{2}}\right)^{t}\leqslant
30^{t}$

where in the last inequality we used that $p_{\mathrm{min}}^{1/s}\geqslant
e^{-1}$. Therefore $v^{*}$ satisfies the desired moment constraints.

Finally, we prove that $\|\operatorname{cov}(\bm{z})^{1/2}v^{*}\|^{2}\leqslant
8$. Note that
$\operatorname{cov}(\bm{z})=\operatorname{cov}(\bm{\mu}_{z})+\Sigma_{z}$. We
have then

$\displaystyle\|\operatorname{cov}(\bm{z})^{1/2}v^{*}\|^{2}$
$\displaystyle=(v^{*})^{\top}\operatorname{cov}(\bm{z})v^{*}$
$\displaystyle=(v^{*})^{\top}\operatorname{cov}(\bm{\mu}_{z})v^{*}+(v^{*})^{\top}\Sigma_{z}v^{*}$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{\leqslant}}\max_{i,j}\langle\mu_{i}-\mu_{j},v^{*}\rangle^{2}+\|\Sigma_{z}^{1/2}v^{*}\|^{2}$
$\displaystyle=\frac{\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4}}{(\mathbb{E}\langle\bm{z},v\rangle^{2s})^{1/s}}+\frac{\|\Sigma_{z}^{1/2}\Sigma_{z}^{-1}(\mu_{a}-\mu_{b})\|^{2}}{(\mathbb{E}\langle\bm{z},v\rangle^{2s})^{1/s}}$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{\leqslant}}\frac{\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4}}{p_{\mathrm{min}}^{2/s}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4}}+\frac{\|\Sigma_{z}^{1/2}\Sigma_{z}^{-1}(\mu_{a}-\mu_{b})\|^{2}}{p_{\mathrm{min}}^{2/s}\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{4}}$
$\displaystyle\stackrel{{\scriptstyle(3)}}{{\leqslant}}e^{2}+e^{2}\frac{1}{\|\Sigma_{z}^{-1/2}(\mu_{a}-\mu_{b})\|^{2}}$
$\displaystyle\stackrel{{\scriptstyle(4)}}{{=}}e^{2}+o(1)\leqslant 8,$

where in (1) we used that
$(v^{*})^{\top}\operatorname{cov}(\bm{\mu}_{z})v^{*}=(v^{*})^{\top}\mathbb{E}\bm{\mu}_{z}\bm{\mu}_{z}^{\top}v^{*}\leqslant\max_{\mu_{z}}(v^{*})^{\top}\mu_{z}\mu_{z}^{\top}v^{*}$
for $\mu_{z}$ in the support of $\bm{\mu}_{z}$, in (2) we used the lower bound
that we derived above on $\mathbb{E}\langle\bm{z},v\rangle^{2s}$, in (3) we
used that $p_{\mathrm{min}}^{1/s}\geqslant e^{-1}$, and in (4) we used mean
separation.

∎

###  4.2 Finite sample bounds

Recall that, to apply Lemma 4.2, we need to find a pseudo-expectation that
satisfies the moment lower bound
$\mathbb{E}\langle\bm{z},v\rangle^{2s}\geqslant c^{s}$ and the moment upper
bound $\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant C^{t}\\}$ for $s\ll t$.
Lemma 4.8 shows that it suffices to find a pseudo-expectation that satisfies
$\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\lesssim 8$,
$\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2s}\succsim c^{s}$, and
$\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2t}\lesssim C^{t}$. Without the bound
on the norm of $\widehat{\operatorname{cov}}(\bm{z})^{1/2}v$ the errors may be
arbitrarily large.

The result is supported by Lemma 4.9, which shows that quadratics in the
empirical covariance matrix are close to quadratics in the population
covariance matrix of the components, and by Lemma 4.10, which shows that the
empirical moments are close to the population moments. The proofs of these two
lemmas are deferred to the appendix.

######  Lemma 4.8 (Moment constraints from empirical moment constraints).

Let $\eta<0.001$. Let $c,C\geqslant 0$ and let $t\geqslant 1$ integer. For
$n\geqslant(p_{\mathrm{min}}^{-1}d)^{O(t)}\eta^{-2}\varepsilon^{-1}$, with
probability $1-\varepsilon$,

$\displaystyle\mathopen{}\mathclose{{}\left\\{\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\leqslant(1+\eta)\cdot
8,\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2s}\geqslant
c^{s}+\eta,\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2t}\leqslant
C^{t}-\eta}\right\\}$
$\displaystyle\qquad\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\mathbb{E}\langle\bm{z},v\rangle^{2s}\geqslant
c^{s},\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant C^{t}}\right\\}.$

Furthermore, with probability $1-\varepsilon$, the axiom is satisfied with
$s$, $t$, and $v$ as in Lemma 4.7 and with $c=(1-\eta)^{1/s}$ and
$C=(30^{t}+\eta)^{1/t}$.

###### Proof.

By Lemma 4.9,
$\sststile{2}{}\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant(1+\eta)^{2}\cdot
8\leqslant 9$. Then, by Lemma 4.10,
$\sststile{2s}{}\mathbb{E}\langle\bm{z},v\rangle^{2s}\geqslant c^{s}$ and
$\sststile{2t}{}\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant C^{t}$. With
$n$ as given, this holds with probability at least $1-\varepsilon$.

For the second claim of the lemma, we have by Lemma 4.7 that there exists some
$v$ with $\mathbb{E}\langle\bm{z},v\rangle^{2s}=1$,
$\mathbb{E}\langle\bm{z},v\rangle^{2t}\leqslant 30^{t}$, and
$\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant 8$. By Lemma 4.9, we also
have that
$\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\leqslant(1+\eta)^{2}\cdot
8\leqslant 9$, and then by Lemma 4.10, we also have that
$\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant 1-\eta$ and
$\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant 30^{t}+\eta$. Again,
with $n$ as given, this holds with probability at least $1-\varepsilon$. ∎

######  Lemma 4.9.

Let $C\geqslant 0$. For $n\geqslant
kd^{2}\log^{2}(d/\varepsilon)O(\eta^{-2})$, with probability $1-\varepsilon$,

$\\{\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\leqslant
C\\}\sststile{2}{}\\{\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant(1+\eta)C\\},$
$\\{\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant
C\\}\sststile{2}{}\\{\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\leqslant(1+\eta)C\\}.$

###### Proof.

See Section A.6. ∎

######  Lemma 4.10.

Let $C\geqslant 0$ and let $t\geqslant 1$ integer. For
$n\geqslant(Cp_{\mathrm{min}}^{-1}d)^{O(t)}\eta^{-2}\varepsilon^{-1}$, with
probability $1-\varepsilon$,

$\\{\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant
C\\}\sststile{O(t)}{}\\{\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2t}\leqslant\mathbb{E}\langle\bm{z},v\rangle^{2t}+\eta\\},$
$\\{\|\operatorname{cov}(\bm{z})^{1/2}v\|^{2}\leqslant
C\\}\sststile{O(t)}{}\\{\hat{\mathbb{E}}\langle\bm{z},v\rangle^{2t}\geqslant\mathbb{E}\langle\bm{z},v\rangle^{2t}-\eta\\}.$

###### Proof.

See Section A.6. ∎

###  4.3 Proof of Theorem 4.1

###### Proof.

Let $s=\lceil\log p_{\mathrm{min}}^{-1}\rceil$ and $t=10000000s$. The
algorithm is:

  1. 1.

Compute a pseudo-expectation $\tilde{\mathbb{E}}$ of degree $2t$ over
$v\in\mathbb{R}^{d}$ such that
$\|\widehat{\operatorname{cov}}(\bm{z})^{1/2}v\|^{2}\leqslant 1.01\cdot 8$,
$\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant 1-0.005$, and
$\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant 30^{t}+0.005\\}$.

  2. 2.

Construct a separating polynomial $q$ based on $\tilde{\mathbb{E}}$ as in
Lemma 4.2.

  3. 3.

Return $q$.

We now analyze the algorithm. First, we argue that there exists a pseudo-
expectation $\tilde{\mathbb{E}}$ that satisfies the given constraints. Note
that $t\leqslant\max_{i,j}\|\Sigma_{z}^{-1/2}(\mu_{i}-\mu_{j})\|^{2}$ if
$C_{sep}$ is a large enough constant. Therefore, the conditions of Lemma 4.7
for $s$ and $t$ are satisfied. Then, by Lemma 4.8, for $n\geqslant n_{0}$,
there exists a vector $v\in\mathbb{R}^{d}$ that satisfies the given
constraints. Then, there also exists a pseudo-expectation that satisfies the
constraints.

Second, we argue that $q$ has the desired properties. By Lemma 4.8,
$\tilde{\mathbb{E}}$ also sastisfies with high probability that
$\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant 1-0.01\geqslant
0.99^{s}$ and $\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant
30^{t}+0.01\leqslant 31^{t}$. Then the conditions of Lemma 4.2 are satisfied
with $c=0.99$ and $C=31$. Then, we are guaranteed to return a separating
polynomial $q$ with the following properties:

  * •

For independent random vectors $\bm{y}$ and $\bm{y}^{\prime}$ sampled from
different components, we have with probability at least $0.99999$ that

$q(\bm{y}-\bm{y}^{\prime})\geqslant\frac{1}{2}\mathopen{}\mathclose{{}\left(\frac{c}{4}}\right)^{s}\geqslant\frac{1}{2}\mathopen{}\mathclose{{}\left(\frac{0.99}{4}}\right)^{s}\geqslant\frac{1}{8^{s}}.$

  * •

For independent random vectors $\bm{y}$ and $\bm{y}^{\prime}$ sampled from the
same component, we have with probability at least $0.99999$ that

$q(\bm{y}-\bm{y}^{\prime})\leqslant
101\mathopen{}\mathclose{{}\left(\frac{4Cs}{t}}\right)^{s}\leqslant
11\mathopen{}\mathclose{{}\left(\frac{4\cdot 31\cdot
s}{t}}\right)^{s}\leqslant\frac{1}{80^{s}}.$

The time complexity of the algorithm is dominated by the time to compute the
pseudo-expectation. The pseudo-expectation is of degree $O(\log
p_{\mathrm{min}}^{-1})$ over $d$ variables, and each constraint requires
summing over the $n$ samples. Therefore, the time to compute the pseudo-
expectation is $n\cdot d^{O(\log p_{\mathrm{min}}^{-1})}$. ∎

##  5 Parallel pancakes

The model studied in this section is a well-separated mixture of Gaussians
with colinear means that is in isotropic position. As shown in Section 5.1,
the isotropic position property makes this model similar to the parallel
pancakes construction, in the sense that the only direction in which the
components of the mixture have variance different from $1$ is the direction of
the means. In Section 6 we study the same model without the isotropic position
assumption.

##### Setting.

We consider a mixture of $k$ Gaussian distributions $N(\mu_{i},\Sigma)$ with
mixing weights $p_{i}$ for $i=1,...,k$, where $\mu_{i}\in\mathbb{R}^{d}$,
$\Sigma\in\mathbb{R}^{d\times d}$ is positive definite, and $p_{i}\geqslant 0$
and $\sum_{i=1}^{k}p_{i}=1$. Let $p_{\mathrm{min}}=\min_{i}p_{i}$.

The distribution is in isotropic position: for $\bm{y}$ distributed according
to the mixture, we have $\mathbb{E}\bm{y}=0$ and
$\operatorname{cov}(\bm{y})=I_{d}$.

The distribution also satisfies mean separation and mean colinearity:

  * •

Mean separation: for some $C_{sep}>0$ and for all $i\neq j$,

$\mathopen{}\mathclose{{}\left\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})}\right\|^{2}\geqslant
C_{sep}\log p_{\mathrm{min}}^{-1}.$

  * •

Mean colinearity: for some unit vector $u\in\mathbb{R}^{d}$ and for all $i$,

$\mu_{i}=\langle\mu_{i},u\rangle u.$

Also define $\sigma^{2}=u^{\top}\Sigma u$, which is the variance of the
components in the direction of the means.

######  Theorem 5.1 (Parallel pancakes algorithm).

Consider the Gaussian mixture model defined above, with $C_{sep}$ larger than
some universal constant. Let

$n_{0}=\mathopen{}\mathclose{{}\left(\frac{1}{\sigma^{2}}}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}.$

Given a sample of size $n\geqslant n_{0}$ from the mixture, there exists an
algorithm that runs in time $n^{O(\log p_{\mathrm{min}}^{-1})}$ and returns
with high probability a partition of $[n]$ into $k$ sets $C_{1},...,C_{k}$
such that, if the true clustering of the samples is $S_{1},...,S_{k}$, then
there exists a permutation $\pi$ of $[k]$ such that

$1-\frac{1}{n}\sum_{i=1}^{k}|C_{i}\cap
S_{\pi(i)}|\leqslant\mathopen{}\mathclose{{}\left(\frac{p_{\mathrm{min}}}{k}}\right)^{O(1)}.$

We introduce some further notation for this section. Let $\bm{y}$ be
distributed according to the mixture. We specify the model as
$\bm{y}=\bm{\mu}+\bm{w}$, where $\bm{\mu}$ takes value $\mu_{i}$ with
probability $p_{i}$ and $\bm{w}\sim N(0,\Sigma)$, with $\bm{\mu}$ and $\bm{w}$
independent of each other.

###  5.1 Isotropic position properties

In this section we prove some consequences of the fact that $\bm{y}$ is in
isotropic position. Lemma 5.2 shows that
$\Sigma=I_{d}-(1-\sigma^{2})uu^{\top}$ with $0<\sigma^{2}\leqslant 1$. This
means that $\Sigma$ can have at most one eigenvalue less than $1$ and that the
eigenvectors corresponding to this eigenvalue are parallel to the direction of
the means $u$. Then Lemma 5.3 uses this form of $\Sigma$ to quantify the
separation of the means along direction $u$ in terms of $\sigma^{2}$.

######  Lemma 5.2 (Isotropic position component covariance matrix).

We have that (1)
$\sigma^{2}=1-\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2}$, (2)
$\Sigma=I_{d}-(1-\sigma^{2})uu^{\top}$, and (3) $0<\sigma^{2}\leqslant 1$.

###### Proof.

We have that
$\mathbb{E}\bm{y}=\mathbb{E}\bm{\mu}+\mathbb{E}\bm{w}=\mathbb{E}\bm{\mu}$.
Because the distribution is in isotropoic position, we also have that
$\mathbb{E}\bm{y}=0$, so the equation above implies that
$\mathbb{E}\bm{\mu}=0$. Then
$\operatorname{cov}(\bm{\mu})=\sum_{i=1}^{k}p_{i}\mu_{i}\mu_{i}^{\top}$.

Furthermore, since $\bm{\mu}$ and $\bm{w}$ are independent, we have that
$\operatorname{cov}(\bm{y})=\operatorname{cov}(\bm{\mu})+\operatorname{cov}(\bm{w})=\sum_{i=1}^{k}p_{i}\mu_{i}\mu_{i}^{\top}+\Sigma$.
Because the distribution is in isotropic position, we also have that
$\operatorname{cov}(\bm{y})=I_{d}$, so the equation above implies
$\Sigma=I_{d}-\sum_{i=1}^{k}p_{i}\mu_{i}\mu_{i}^{\top}$. Plugging in
$\mu_{i}=\langle\mu_{i},u\rangle u$, we have
$\Sigma=I_{d}-\mathopen{}\mathclose{{}\left(\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2}}\right)uu^{\top}$.

Then, it follows that $\sigma^{2}=u^{\top}\Sigma
u=1-\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2}$. This proves (1). The fact
that $1-\sigma^{2}=\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2}$ also proves
(2). For (3), $\sigma^{2}>0$ follows by the definition using that $\Sigma$ is
positive definite and $\sigma^{2}\leqslant 1$ follows by (1). ∎

######  Lemma 5.3 (Isotropic position mean separation).

For all $i\neq j$,

$\langle\mu_{i}-\mu_{j},u\rangle^{2}\geqslant C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}.$

###### Proof.

By Lemma 5.2, $\Sigma=I_{d}-(1-\sigma^{2})uu^{\top}$. This implies that
$\Sigma^{-1/2}=I_{d}+\mathopen{}\mathclose{{}\left(1/\sqrt{\sigma^{2}}-1}\right)uu^{\top}$.
Then, using that $\mu_{i}=\langle\mu_{i},u\rangle u$, we have that

$\displaystyle\Sigma^{-1/2}(\mu_{i}-\mu_{j})$
$\displaystyle=\mathopen{}\mathclose{{}\left(I_{d}+\mathopen{}\mathclose{{}\left(\frac{1}{\sqrt{\sigma^{2}}}-1}\right)uu^{\top}}\right)\mathopen{}\mathclose{{}\left(\langle\mu_{i}-\mu_{j},u\rangle
u}\right)=\frac{1}{\sqrt{\sigma^{2}}}\langle\mu_{i}-\mu_{j},u\rangle u.$

Therefore, the separation condition
$\mathopen{}\mathclose{{}\left\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})}\right\|^{2}\geqslant
C_{sep}\log p_{\mathrm{min}}^{-1}$ is equivalent to
$\frac{1}{\sigma^{2}}\langle\mu_{i}-\mu_{j},u\rangle^{2}\geqslant C_{sep}\log
p_{\mathrm{min}}^{-1}$. The conclusion follows by multiplying both sides by
$\sigma^{2}$. ∎

###  5.2 Exact moment direction recovery

In this section we discuss how to recover a direction close to the direction
of the means $u$, assuming oracle access to moments $\mathbb{E}\bm{y}^{\otimes
t}$ for any positive integer $t$. Access to these moments allows us to
calculate exactly directional moments of the form
$\mathbb{E}\langle\bm{y},v\rangle^{t}$, which simplifies the analysis. Finite
sample considerations are discussed starting with Section 5.3.

Theorem 5.4 shows that there exists an algorithm that computes a unit vector
$\hat{u}$ with correlation $1-O(\min(\sigma^{2},\frac{1}{k}))$ with the
direction of the means $u$. We remark that it is necessary for $\hat{u}$ to
have a correlation of at least $1-O(\sigma^{2})$ with $u$ in order for the
components of the mixture to be separated along direction $\hat{u}$.

######  Theorem 5.4 (Direction recovery with exact moments).

Assume oracle access to $\mathbb{E}\bm{y}^{\otimes t}$ for any positive
integer $t$. Then there exists an algorithm with time complexity
$\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot d^{O(\log
p_{\mathrm{min}}^{-1})}$ that outputs a unit vector $\hat{u}\in\mathbb{R}^{d}$
such that $\langle u,\hat{u}\rangle^{2}\geqslant
1-320\min(\sigma^{2},\frac{1}{k})$.

The two main ingredients for Theorem 5.4 are Theorem 5.5, which gives an
algorithm to compute pseudo-expectations over unit vectors correlated with
$u$, and Theorem 5.6, which gives an algorithm to sample from such pseudo-
expectations. We note that Theorem 5.5 can be interpreted as a collection of
sum-of-squares identifiability proofs for the direction of the means $u$.

We state these two supporting theorems and then prove Theorem 5.4. After that,
we work toward proving the supporting theorems.

######  Theorem 5.5 (Direction sum-of-squares identifiability).

Assume oracle access to $\mathbb{E}\bm{y}^{\otimes t}$ for any positive
integer $t$. Then there exists an algorithm with time complexity
$\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot d^{O(\log
p_{\mathrm{min}}^{-1})}$ that computes two peseudo-expectations
$\tilde{\mathbb{E}}_{U}$ and $\tilde{\mathbb{E}}_{L}$ of degree $O(\log
p_{\mathrm{min}}^{-1})$ over a variable $v\in\mathbb{R}^{d}$ such that the
following holds. Let $s=\lceil\log p_{\mathrm{min}}^{-1}\rceil$, let
$t=5000s$, and let $\tau=\frac{800e}{C_{sep}k^{2}}$. Then
$\tilde{\mathbb{E}}_{U}\|v\|^{2}=1$, $\tilde{\mathbb{E}}_{L}\|v\|^{2}=1$, and:

  * •

If $\sigma^{2}\geqslant\tau$, then $\tilde{\mathbb{E}}_{U}\langle
u,v\rangle^{2s}\geqslant(1-\tau)^{s}$.

  * •

If $\sigma^{2}<\tau$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\geqslant(4es)^{s}$, then
$\tilde{\mathbb{E}}_{U}\langle u,v\rangle^{2s}\geqslant(1-\sigma^{2})^{s}$.

  * •

If $\sigma^{2}<0.001$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\leqslant(100s)^{s}$, then
$\tilde{\mathbb{E}}_{L}\langle u,v\rangle^{2t}\geqslant(1-20\sigma^{2})^{t}$.

######  Theorem 5.6 (Direction sum-of-squares sampling).

Let $t\in\mathbb{N}$ and $\varepsilon\in\mathbb{R}$ such that $t\geqslant 1$
and $0\leqslant\varepsilon\leqslant 1/(3t^{2})$. Let $u\in\mathbb{R}^{d}$ be a
unit vector. Given a pseudo-expectation $\tilde{\mathbb{E}}$ of degree $2t$
over a variable $v\in\mathbb{R}^{d}$ that satisfies
$\tilde{\mathbb{E}}\|v\|^{2}=1$ and $\tilde{\mathbb{E}}\langle
u,v\rangle^{2t}\geqslant(1-\varepsilon)^{t}$, there exists an algorithm with
time complexity $d^{O(t)}$ that returns a unit vector
$\hat{u}\in\mathbb{R}^{d}$ such that $\langle u,\hat{u}\rangle^{2}\geqslant
1-16\varepsilon$.

###### Proof of Theorem 5.4.

Let $\tau=\frac{800e}{C_{sep}k^{2}}$. The algorithm is:

  1. 1.

Run the algorithm from Theorem 5.5 to obtain pseudo-expectations
$\tilde{\mathbb{E}}_{U}$ and $\tilde{\mathbb{E}}_{L}$.

  2. 2.

Run the algorithm from Theorem 5.6 for pseudo-expectations
$\tilde{\mathbb{E}}_{U}$ and $\tilde{\mathbb{E}}_{L}$ to obtain unit vectors
$\hat{u}_{U}\in\mathbb{R}^{d}$ and $\hat{u}_{L}\in\mathbb{R}^{d}$,
respectively.

  3. 3.

If $\sigma^{2}\geqslant\tau$, return $\hat{u}_{U}$. Else, for $s=\lceil\log
p_{\mathrm{min}}^{-1}\rceil$, if
$\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s}\geqslant(50s)^{s}$, return
$\hat{u}_{U}$. Else, return $\hat{u}_{L}$.

We now analyze the algorithm. We consider the three possible cases in step (3)
of the algorithm:

  * •

Suppose $\sigma^{2}\geqslant\tau$. Then Theorem 5.5 guarantees that
$\tilde{\mathbb{E}}_{U}\langle u,v\rangle^{2s}\geqslant(1-\tau)^{s}$, so by
Theorem 5.6 we have $\langle u,\hat{u}_{U}\rangle^{2}\geqslant
1-16\tau\geqslant 1-16\min(\sigma^{2},\tau)$.

  * •

Suppose $\sigma^{2}<\tau$ and
$(\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s})^{1/s}\geqslant 50s$. We
have by Lemma 5.8 that
$(\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s})^{1/s}\leqslant(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}+es,$
so it must be the case that
$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\geqslant(50-e)s\geqslant
4es$. Then Theorem 5.5 guarantees that $\tilde{\mathbb{E}}_{U}\langle
u,v\rangle^{2s}\geqslant(1-\sigma^{2})^{s}$, so by Theorem 5.6 we have
$\langle u,\hat{u}_{U}\rangle^{2}\geqslant 1-16\sigma^{2}\geqslant
1-16\min(\sigma^{2},\tau)$.

  * •

Suppose $\sigma^{2}<\tau$ and
$\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s}<(50s)^{s}$. We have by Lemma
5.9 that
$(\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s})^{1/s}\geqslant(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s},$
so it must be the case that
$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}<50s\leqslant 100s$. Then
Theorem 5.5 guarantees that $\tilde{\mathbb{E}}_{L}\langle
u,v\rangle^{2t}\geqslant(1-20\sigma^{2})^{t}$, so by Theorem 5.6 we have
$\langle u,\hat{u}_{L}\rangle^{2}\geqslant 1-16\cdot 20\sigma^{2}\geqslant
1-320\min(\sigma^{2},\tau)$.

Let $\hat{u}$ be the unit vector returned by step (3) of the algorithm. Then
we are guaranteed that in all cases $\langle u,\hat{u}\rangle^{2}\geqslant
1-320\min(\sigma^{2},\tau)\geqslant 1-320\min(\sigma^{2},\frac{1}{k})$, where
we used the loose upper bound $\tau\leqslant\frac{1}{k}$.

The time complexity of the algorithm is dominated by the time to run the
algorithm from Theorem 5.5. ∎

####  5.2.1 Sum-of-squares identifiability (proof of Theorem 5.5)

We prove a number of supporting lemmas and then prove Theorem 5.5. The most
important components are Lemma 5.11 and Lemma 5.12, which give sum-of-squares
proofs that, for suitably chosen $s,t=O(\log p_{\mathrm{min}}^{-1})$, either
the maximizer of $\mathbb{E}\langle\bm{y},v\rangle^{2s}$ or the minimizer of
$\mathbb{E}\langle\bm{y},v\rangle^{2t}$ over unit vectors $v$ must be close to
$u$.

We start with Lemma 5.7, Lemma 5.8 and Lemma 5.9, which give sum-of-squares
bounds on the moments of the mixture. Informally, for $t=\Omega(\log
p_{\mathrm{min}}^{-1})$, these bounds correspond to the following
decomposition of the directional $2t$ moments:

$(\mathbb{E}\langle\bm{y},v\rangle^{2t})^{1/t}=\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+\Theta(1)\cdot
t(v^{\top}\Sigma v).$ (5.1)

######  Lemma 5.7 (Moment equality).

For $t\geqslant 1$ integer,

$\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}=\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathbb{E}\langle\bm{\mu},u\rangle^{2s}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!.$

###### Proof.

$\displaystyle\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}$
$\displaystyle=\mathbb{E}(\langle\bm{\mu},v\rangle+\langle\bm{w},v\rangle)^{2t}=\sum_{j=0}^{2t}\binom{2t}{j}\mathbb{E}\langle\bm{\mu},v\rangle^{j}\langle\bm{w},v\rangle^{2t-j}$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{=}}\sum_{j=0}^{2t}\binom{2t}{j}\mathbb{E}\langle\bm{\mu},v\rangle^{j}\mathbb{E}\langle\bm{w},v\rangle^{2t-j}\stackrel{{\scriptstyle(2)}}{{=}}\sum_{s=0}^{t}\binom{2t}{2s}\mathbb{E}\langle\bm{\mu},v\rangle^{2s}\mathbb{E}\langle\bm{w},v\rangle^{2t-2s}$
$\displaystyle\stackrel{{\scriptstyle(3)}}{{=}}\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathbb{E}\langle\bm{\mu},u\rangle^{2s}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$

where in (1) we used that $\bm{\mu}$ and $\bm{w}$ are independent, in (2) we
used that $\mathbb{E}\langle\bm{w},v\rangle^{2t-j}=0$ for $2t-j$ odd, and in
(3) we used that $\bm{\mu}=\langle\bm{\mu},u\rangle u$ and that
$\mathbb{E}\langle\bm{w},v\rangle^{2t-2s}=(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$. ∎

######  Lemma 5.8 (Moment upper bound).

For $t\geqslant 1$ integer,

$\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}\leqslant\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+et(v^{\top}\Sigma
v)}\right)^{t}.$

###### Proof.

Starting with the result in Lemma 5.7,

$\displaystyle\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}$
$\displaystyle=\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathbb{E}\langle\bm{\mu},u\rangle^{2s}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{\leqslant}}\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{\leqslant}}\sum_{s=0}^{t}\binom{t}{s}\langle
u,v\rangle^{2s}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}(v^{\top}\Sigma
v)^{t-s}(et)^{t-s}$ $\displaystyle=\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+et(v^{\top}\Sigma
v)}\right)^{t}.$

In (1) we used that $s\leqslant t$ and Jensen’s inequality as follows:

$\displaystyle\mathbb{E}\langle\bm{\mu},u\rangle^{2s}$
$\displaystyle=\mathbb{E}\langle\bm{\mu},u\rangle^{2t(s/t)}\leqslant\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}.$

In (2) we used that $\binom{2t}{2s}(2t-2s-1)!!\leqslant\binom{t}{s}(et)^{t-s}$
for $0\leqslant s\leqslant t$ integers, which is proved in Lemma A.25. ∎

######  Lemma 5.9 (Moment lower bound).

For $t\geqslant 1$ integer,

$\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}\geqslant\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+p_{\mathrm{min}}^{1/t}t/2(v^{\top}\Sigma
v)}\right)^{t}.$

###### Proof.

Starting with the result in Lemma 5.7,

$\displaystyle\sststile{2t}{}\mathbb{E}\langle\bm{y},v\rangle^{2t}$
$\displaystyle=\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathbb{E}\langle\bm{\mu},u\rangle^{2s}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$
$\displaystyle\stackrel{{\scriptstyle(1)}}{{\geqslant}}\sum_{s=0}^{t}\binom{2t}{2s}\langle
u,v\rangle^{2s}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}\mathopen{}\mathclose{{}\left(p_{\mathrm{min}}^{1/t}}\right)^{t-s}(v^{\top}\Sigma
v)^{t-s}(2t-2s-1)!!$
$\displaystyle\stackrel{{\scriptstyle(2)}}{{\geqslant}}\sum_{s=0}^{t}\binom{t}{s}\langle
u,v\rangle^{2s}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}\mathopen{}\mathclose{{}\left(p_{\mathrm{min}}^{1/t}}\right)^{t-s}(v^{\top}\Sigma
v)^{t-s}(t/2)^{t-s}$ $\displaystyle=\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+p_{\mathrm{min}}^{1/t}t/2(v^{\top}\Sigma
v)}\right)^{t}.$

In (1) we used that $s\leqslant t$ and the fact that the $s$-norm is greater
than or equal to the $t$-norm as follows:

$\displaystyle\mathbb{E}\langle\bm{\mu},u\rangle^{2s}$
$\displaystyle=\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2s}=\sum_{i=1}^{k}(p_{i}^{1/s}\langle\mu_{i},u\rangle^{2})^{s}$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(\sum_{i=1}^{k}(p_{i}^{1/s}\langle\mu_{i},u\rangle^{2})^{t}}\right)^{s/t}=\mathopen{}\mathclose{{}\left(\sum_{i=1}^{k}p_{i}^{t/s}\langle\mu_{i},u\rangle^{2t}}\right)^{s/t}$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(p_{\mathrm{min}}^{t/s-1}\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2t}}\right)^{s/t}=p_{\mathrm{min}}^{1-s/t}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{s/t}.$

In (2) we used that
$\binom{2t}{2s}(2t-2s-1)!!\geqslant\binom{t}{s}(t/2)^{t-s}$ for $0\leqslant
s\leqslant t$ integers, which is proved in Lemma A.25.

∎

Lemma 5.10 shows that the contribution of the means to the $\Omega(\log
p_{\mathrm{min}}^{-1})$ moments in direction $u$ is lower bounded by
$\Omega(k^{2}\sigma^{2}\log p_{\mathrm{min}}^{-1})$. This result is used in
some of the later proofs to argue that if the mean contribution is small, then
$\sigma^{2}$ is small, and conversely, that if $\sigma^{2}$ is large, then the
mean contribution is large.

######  Lemma 5.10.

For $2s\geqslant\lceil\log p_{\mathrm{min}}^{-1}\rceil$ integer,

$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\geqslant\frac{C_{sep}}{100}k^{2}\sigma^{2}\log
p_{\mathrm{min}}^{-1}.$

###### Proof.

By Lemma 5.3, for all $i\neq j$,
$|\langle\mu_{i}-\mu_{j},u\rangle|\geqslant\sqrt{C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}}$. Then there exist $a,b\in[k]$ such that
$|\langle\mu_{a}-\mu_{b},u\rangle|\geqslant(k-1)\sqrt{C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}}$. Hence, there exists $a\in[k]$ such that
$|\langle\mu_{a},u\rangle|\geqslant\frac{k-1}{2}\sqrt{C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}}$. Then

$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\geqslant
p_{\mathrm{min}}^{1/s}\max_{i}\langle\mu_{i},u\rangle^{2}\geqslant
p_{\mathrm{min}}^{1/s}\mathopen{}\mathclose{{}\left(\frac{k-1}{2}}\right)^{2}C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}\geqslant\frac{C_{sep}}{100}k^{2}\sigma^{2}\log
p_{\mathrm{min}}^{-1},$

where we used that $p_{\mathrm{min}}^{1/s}\geqslant e^{-2}$. ∎

We now state and prove the sum-of-squares identifiability proofs of Lemma 5.11
and Lemma 5.12. Let $s,t=O(\log p_{\mathrm{min}}^{-1})$ with $s\ll t$. Lemma
5.11 proves that, in the case
$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\geqslant\Theta(s)$, if
$\mathbb{E}\langle\bm{y},v\rangle^{2s}$ is close to its maximum value over
unit vectors $v$, then $\langle u,v\rangle^{2s}$ is close to $1$. Lemma 5.12
proves that, in the opposite case
$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\leqslant\Theta(s)$, if
$\mathbb{E}\langle\bm{y},v\rangle^{2t}$ is close to its minimum value over
unit vectors $v$, then $\langle u,v\rangle^{2t}$ is close to $1$.

######  Lemma 5.11 (Direction sum-of-squares identifiability from moment
maximization).

Let $M\geqslant 2$. Let $s$ be an integer such that $2s\geqslant\lceil\log
p_{\mathrm{min}}^{-1}\rceil$. Suppose that
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\geqslant(Mes)^{s}$. Then, for
$\varepsilon\leqslant\sigma^{2}/M$,

$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1,\mathbb{E}\langle\bm{y},v\rangle^{2s}\geqslant(1-\varepsilon)\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}-\varepsilon}\right)}\right\\}\sststile{2s}{}\mathopen{}\mathclose{{}\left\\{\langle
u,v\rangle^{2s}\geqslant(1-4\sigma^{2}/M)^{s}}\right\\}.$

Furthermore, $v=u$ satisfies the axiom with $\varepsilon=0$.

###### Proof.

Substitute the upper bound of Lemma 3 into the axiom:

$\sststile{2s}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}}\right)^{1/s}+es(v^{\top}\Sigma
v)}\right)^{s}\geqslant(1-\varepsilon)\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}-\varepsilon}\right).$

Divide by $\mathbb{E}\langle\bm{\mu},u\rangle^{2s}$:

$\sststile{2s}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+(v^{\top}\Sigma
v)\frac{es}{\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}}\right)^{1/s}}}\right)^{s}\geqslant(1-\varepsilon)\mathopen{}\mathclose{{}\left(1-\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2s}}}\right).$

Recall that
$\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}}\right)^{1/s}\geqslant
Mes$, and substitute the lower bound on both sides:

$\displaystyle\sststile{2s}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+\frac{1}{M}(v^{\top}\Sigma v)}\right)^{s}$
$\displaystyle\geqslant(1-\varepsilon)\mathopen{}\mathclose{{}\left(1-\frac{\varepsilon}{(Mes)^{s}}}\right)$

Use that $v^{\top}\Sigma v=1-(1-\sigma^{2})\langle u,v\rangle^{2}$:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+\frac{1}{M}(1-(1-\sigma^{2})\langle
u,v\rangle^{2})}\right)^{s}\geqslant(1-\varepsilon)\mathopen{}\mathclose{{}\left(1-\frac{\varepsilon}{(Mes)^{s}}}\right).$

We simplify now the right-hand side. Use the loose bound
$1-\varepsilon/(Mes)^{s}\geqslant 1-\varepsilon\geqslant(1-\varepsilon)^{s-1}$
to obtain

$\sststile{2s}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+\frac{1}{M}(1-(1-\sigma^{2})\langle
u,v\rangle^{2})}\right)^{s}\geqslant(1-\varepsilon)^{s}.$

Finally, apply Lemma A.15 with $x=\langle u,v\rangle$ and
$\gamma=\frac{1}{1-\varepsilon}$ to obtain that

$\displaystyle\sststile{2s}{}\langle u,v\rangle^{2s}$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(\frac{M-\frac{1}{1-\varepsilon}}{\frac{1}{1-\varepsilon}}\frac{1}{M-1+\sigma^{2}}}\right)^{s}=\mathopen{}\mathclose{{}\left(\frac{M-1-M\varepsilon}{M-1+\sigma^{2}}}\right)^{s}$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(\frac{M-1-\sigma^{2}}{M-1+\sigma^{2}}}\right)^{s}\geqslant\mathopen{}\mathclose{{}\left(1-4\sigma^{2}/M}\right)^{s}.$

To show that $v=u$ satisfies the axiom, simply note that Lemma 5.12 implies
that
$\mathbb{E}\langle\bm{y},u\rangle^{2s}\geqslant\mathbb{E}\langle\bm{\mu},u\rangle^{2s}$.

∎

######  Lemma 5.12 (Direction sum-of-squares identifiability from moment
minimization).

Suppose $\sigma^{2}<0.001$. Let $s$ be an integer such that
$2s\geqslant\lceil\log p_{\mathrm{min}}^{-1}\rceil$. Suppose that
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\leqslant(100s)^{s}$. Let $t$ be an
integer such that $t\geqslant 5000s$. Then, for
$\varepsilon\leqslant\sigma^{2}/100$,

$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1,\mathbb{E}\langle\bm{y},v\rangle^{2t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+et\sigma^{2}}\right)^{t}+\varepsilon}\right)}\right\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\langle
u,v\rangle^{2t}\geqslant(1-20\sigma^{2})^{t}}\right\\}.$

Furthemore, $v=u$ satisfies the axiom with $\varepsilon=0$.

###### Proof.

We start by proving that, for $t\geqslant s$,
$\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\leqslant(100e^{2}s)^{t}$. We have
that

$p_{\mathrm{min}}\cdot\max_{i}\langle\mu_{i},u\rangle^{2s}\leqslant\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\leqslant(100s)^{s}.$

Taking the $s$-th root and using that $p_{\mathrm{min}}^{-1/s}\leqslant
e^{2}$, we obtain that $\max_{i}\langle\mu_{i},u\rangle^{2}\leqslant
100e^{2}s$. Therefore,
$\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\leqslant\max_{i}\langle\mu_{i},u\rangle^{2t}=(100e^{2}s)^{t}$.

We now proceed with the main claim. Substitute the lower bound of Lemma 4 into
the axiom:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+p_{\mathrm{min}}^{1/t}t/2(v^{\top}\Sigma
v)}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+\sigma^{2}et}\right)^{t}+\varepsilon}\right).$

Divide by $\mathbb{E}\langle\bm{\mu},u\rangle^{2t}$:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+(v^{\top}\Sigma
v)\frac{p_{\mathrm{min}}^{1/t}t/2}{\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}}}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(1+\sigma^{2}\frac{et}{\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}}}\right)^{t}+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}}\right).$

Let
$\Delta=\frac{p_{\mathrm{min}}^{1/t}t/2}{\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}}$.
Then

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+\Delta(v^{\top}\Sigma
v)}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(1+\sigma^{2}\frac{2e}{p_{\mathrm{min}}^{1/t}}\Delta}\right)^{t}+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}}\right).$

Note that
$\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}\leqslant
100e^{2}s$ and $p_{\mathrm{min}}^{1/t}=e^{-s/t}\geqslant e^{-1}$. Then
$\Delta\geqslant\frac{e^{-1}t/2}{100e^{2}s}$. For $t\geqslant 5000s$ we have
then $\Delta\geqslant 10$ and $p_{\mathrm{min}}^{-1/t}\leqslant 1.4$. Then
$\frac{2e}{p_{\mathrm{min}}^{1/t}}\Delta\leqslant 8\Delta$. Then:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\langle
u,v\rangle^{2}+\Delta(v^{\top}\Sigma
v)}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(1+8\Delta\sigma^{2}}\right)^{t}+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}}\right).$

Divide by $\mathopen{}\mathclose{{}\left(1+8\Delta\sigma^{2}}\right)^{t}$:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\frac{\langle
u,v\rangle^{2}+\Delta(v^{\top}\Sigma
v)}{1+8\Delta\sigma^{2}}}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(1+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\mathopen{}\mathclose{{}\left(1+8\Delta\sigma^{2}}\right)^{t}}}\right).$

Use that $v^{\top}\Sigma v=1-(1-\sigma^{2})\langle u,v\rangle^{2}$:

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\frac{\langle
u,v\rangle^{2}+\Delta(1-(1-\sigma^{2})\langle
u,v\rangle^{2})}{1+8\Delta\sigma^{2}}}\right)^{t}\leqslant(1+\varepsilon)\mathopen{}\mathclose{{}\left(1+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\mathopen{}\mathclose{{}\left(1+8\Delta\sigma^{2}}\right)^{t}}}\right).$

We simplify now the term involving $\varepsilon$. Note that, by Jensen’s
inequality,
$\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\geqslant(\mathbb{E}\langle\bm{\mu},u\rangle^{2})^{t}=(1-\sigma^{2})^{t}$.
Also note that $(1-\sigma^{2})(1+8\Delta\sigma^{2})\geqslant 1$ for
$\sigma^{2}\leqslant 1/2$ and $\Delta\geqslant 10$. Then use the loose bound

$1+\frac{\varepsilon}{\mathbb{E}\langle\bm{\mu},u\rangle^{2t}\mathopen{}\mathclose{{}\left(1+8\Delta\sigma^{2}}\right)^{t}}\leqslant
1+\frac{\varepsilon}{((1-\sigma^{2})(1+8\Delta\sigma^{2}))^{t}}\leqslant
1+\varepsilon\leqslant(1+\varepsilon)^{t-1}$

to obtain

$\sststile{2t}{}\mathopen{}\mathclose{{}\left(\frac{\langle
u,v\rangle^{2}+\Delta(1-(1-\sigma^{2})\langle
u,v\rangle^{2})}{1+8\Delta\sigma^{2}}}\right)^{t}\leqslant(1+\varepsilon)^{t}.$

Finally, apply Lemma A.16 with $x=\langle u,v\rangle$ and
$\gamma=\frac{1}{1+\varepsilon}$ to obtain that

$\displaystyle\sststile{2t}{}\langle u,v\rangle^{2t}$
$\displaystyle\geqslant\mathopen{}\mathclose{{}\left(\frac{\frac{1}{1+\varepsilon}\Delta-1}{\frac{1}{1+\varepsilon}(\Delta-1)}(1-10\sigma^{2})}\right)^{t}=\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(1-\frac{\varepsilon}{\Delta-1}}\right)(1-10\sigma^{2})}\right)^{t}$
$\displaystyle=\mathopen{}\mathclose{{}\left((1-\sigma^{2}/100)(1-10\sigma^{2})}\right)^{t}\geqslant\mathopen{}\mathclose{{}\left(1-20\sigma^{2}}\right)^{t}.$

To show that $v=u$ satisfies the axiom, simply note that Lemma 5.11 implies
that
${\mathbb{E}\langle\bm{y},u\rangle^{2t}\leqslant\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+et\sigma^{2}}\right)^{t}}$.

∎

We now prove Theorem 5.5.

###### Proof of Theorem 5.5..

Let $s=\lceil\log p_{\mathrm{min}}^{-1}\rceil$, $t=5000s$, and
$\tau=800e/(C_{sep}k^{2})$. The algorithm is:

  1. 1.

If $\sigma^{2}\geqslant\tau$, then let $M=C_{sep}k^{2}\sigma^{2}/(200e)$.
Else, let $M=4$.

  2. 2.

Binary search up to resolution $\sigma^{2}/(100M)$ the largest $T_{U}$ in the
interval $[0,(p_{\mathrm{min}}^{-1})^{s}]$ such that there exists a
degree-$2s$ pseudo-expectation that satisfies
$\\{\|v\|^{2}=1,\mathbb{E}\langle\bm{y},v\rangle^{2s}\geqslant T_{U}\\}$. Let
$\tilde{\mathbb{E}}_{U}$ be the resulting pseudo-expectation for this $T_{U}$.

  3. 3.

Binary search up to resolution $\sigma^{2}/10000$ the smallest $T_{L}$ in the
interval $[0,(p_{\mathrm{min}}^{-1}+et)^{t}]$ such that there exists a
degree-$2t$ pseudo-expectation that satisfies
$\\{\|v\|^{2}=1,\mathbb{E}\langle\bm{y},v\rangle^{2t}\leqslant T_{L}\\}$. Let
$\tilde{\mathbb{E}}_{L}$ be the resulting pseudo-expectation for this $T_{L}$.

  4. 4.

Return $\tilde{\mathbb{E}}_{U}$ and $\tilde{\mathbb{E}}_{L}$.

We now analyze the algorithm. To begin with, suppose that the $T_{U}$ found is
at least the maximum value of $\mathbb{E}\langle\bm{y},v\rangle^{2s}$ and that
the $T_{L}$ found is at most the minimum value of
$\mathbb{E}\langle\bm{y},v\rangle^{2t}$. In this case $\tilde{\mathbb{E}}_{U}$
and $\tilde{\mathbb{E}}_{L}$ satisfy the axioms of Lemma 5.11 and Lemma 5.12,
respectively. Then our algorithm achieves the stated guarantees:

  * •

Suppose $\sigma^{2}\geqslant\tau$. Note that, in this case,
$M=C_{sep}k^{2}\sigma^{2}/(200e)\geqslant 2$. By Lemma 5.10, we have that

$(\mathbb{E}\langle\bm{\mu},u\rangle^{2s})^{1/s}\geqslant\frac{C_{sep}}{100}k^{2}\sigma^{2}\log
p_{\mathrm{min}}^{-1}=2eM\log p_{\mathrm{min}}^{-1}\geqslant Mes.$

Then the conditions of Lemma 5.11 are satisfied, and $\tilde{\mathbb{E}}_{U}$
satisfies $\langle
u,v\rangle^{2s}\geqslant(1-4\sigma^{2}/M)^{s}=(1-\tau)^{s}$.

  * •

Suppose $\sigma^{2}<\tau$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\geqslant(4es)^{s}$. Note that, in
this case, $M=4$. Then the conditions of Lemma 5.11 are satisfied, and
$\tilde{\mathbb{E}}_{U}$ satisfies $\langle
u,v\rangle^{2s}\geqslant(1-4\sigma^{2}/M)^{s}=(1-\sigma^{2})^{s}$.

  * •

Suppose $\sigma^{2}<0.001$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\leqslant(100s)^{s}$. Then the
condition of Lemma 5.12 are satisfied, and $\tilde{\mathbb{E}}_{L}$ satisfies
$\langle u,v\rangle^{2s}\geqslant(1-20\sigma^{2})^{s}$.

We argue now that $T_{U}$ is large enough and that $T_{L}$ is small enough in
order for the pseudo-expectations to satisfy the axioms of the lemmas. For
that, we need

$T_{U}\geqslant(1-\sigma^{2}/M)\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2s}-\sigma^{2}/M}\right),$
$T_{L}\leqslant(1+\sigma^{2}/100)\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+\sigma^{2}et}\right)^{t}+\sigma^{2}/100}\right).$

We prove that the intervals in which we binary search $T_{U}$ and $T_{L}$
contain $\mathbb{E}\langle\bm{\mu},u\rangle^{2s}$ and
$((\mathbb{E}\langle\bm{\mu},u\rangle^{2t})^{1/t}+\sigma^{2}et)^{t}$,
respectively. Then, binary search with the proposed resolutions is guaranteed
to find $T_{U}$ and $T_{L}$ that satisfy the bounds stated above.

Using that $\mathbb{E}\langle\bm{\mu},u\rangle^{2}=1-\sigma^{2}$, we have that

$\displaystyle\mathbb{E}\langle\bm{\mu},u\rangle^{2t}$
$\displaystyle=\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2t}=\sum_{i=1}^{k}(p_{i}^{1/t}\langle\mu_{i},u\rangle^{2})^{t}\leqslant\mathopen{}\mathclose{{}\left(\sum_{i=1}^{k}p_{i}^{1/t}\langle\mu_{i},u\rangle^{2}}\right)^{t}$
$\displaystyle\leqslant\mathopen{}\mathclose{{}\left(p_{\mathrm{min}}^{1/t-1}\sum_{i=1}^{k}p_{i}\langle\mu_{i},u\rangle^{2}}\right)^{t}=p_{\mathrm{min}}^{-(t-1)}\mathbb{E}\langle\bm{\mu},u\rangle^{2}=p_{\mathrm{min}}^{-(t-1)}(1-\sigma^{2})$
$\displaystyle\leqslant p_{\mathrm{min}}^{-t}$

and

$\displaystyle\mathbb{E}\langle\bm{\mu},u\rangle^{2t}$
$\displaystyle\geqslant(\mathbb{E}\langle\bm{\mu},u\rangle^{2})^{t}=(1-\sigma^{2})^{t}.$

Therefore,

$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\in[(1-\sigma^{2})^{s},(p_{\mathrm{min}}^{-1})^{s}],$
$\mathopen{}\mathclose{{}\left(\mathopen{}\mathclose{{}\left(\mathbb{E}\langle\bm{\mu},u\rangle^{2t}}\right)^{1/t}+\sigma^{2}et}\right)^{t}\in[(1-\sigma^{2})^{t},(p_{\mathrm{min}}^{-1}+et)^{t}].$

Then the intervals in which we binary search are wide enough and binary search
is guaranteed to succeed.

The time complexity of the algorithm is given by the number of steps in the
binary search multiplied by the time to compute each of the pseudo-
expectations. The number of steps in the binary search is

$O\mathopen{}\mathclose{{}\left(\max\mathopen{}\mathclose{{}\left\\{\log((p_{\mathrm{min}}^{-1})^{s}k^{2}),\log\frac{(p_{\mathrm{min}}^{-1})^{s}}{\sigma^{2}},\log\frac{(p_{\mathrm{min}}^{-1}+et)^{t}}{\sigma^{2}}}\right\\}}\right)=O\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}+\log^{2}p_{\mathrm{min}}^{-1}}\right).$

For each step, we compute a pseudo-expectation of degree $O(\log
p_{\mathrm{min}}^{-1})$ over $d$ variables, which requires time $d^{O(\log
p_{\mathrm{min}}^{-1})}$. Therefore the time complexity is

$O\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}+\log^{2}p_{\mathrm{min}}^{-1}}\right)\cdot
d^{O(\log
p_{\mathrm{min}}^{-1})}=\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot
d^{O(\log p_{\mathrm{min}}^{-1})}.$

∎

####  5.2.2 Sum-of-squares sampling (proof of Theorem 5.6)

We state and prove Lemma 5.13, which is used in the proof of Theorem 5.6. This
lemma shows that, given a symmetric postivie definite matrix $M$ correlated
with a rank-$1$ matrix $uu^{\top}$ for a unit vector $u$, there exists an
algorithm to recover a unit vector correlated with $u$. After that, we proceed
to prove the theorem.

######  Lemma 5.13 (Matrix rank-$1$ approximation).

Let $0\leqslant\varepsilon<\frac{1}{8}$. Let $u\in\mathbb{R}^{d}$ be a unit
vector. Given a symmetric positive semi-definite matrix
$M\in\mathbb{R}^{d\times d}$ with $\|M\|_{F}\leqslant 1$ such that $\langle
uu^{\top},M\rangle_{F}\geqslant 1-\varepsilon$, there exists a polynomial-time
algorithm that finds a unit vector $\hat{u}\in\mathbb{R}^{d}$ such that
$\langle u,\hat{u}\rangle^{2}\geqslant 1-8\varepsilon$.

###### Proof.

The algorithm is to compute $vv^{\top}$ as the best rank-$1$ approximation of
$M$ and return $\frac{v}{\|v\|}$, which is uniquely defined up to a sign flip.

We now analyze the accuracy of the algorithm. We have that $\langle
uu^{\top},M\rangle_{F}\geqslant 1-\varepsilon$, so
$\|uu^{\top}-M\|_{F}^{2}\leqslant 2-2\langle uu^{\top},M\rangle_{F}\leqslant
2\varepsilon$. For $vv^{\top}$ the best rank-$1$ approximation of $M$, we have
then that

$\|uu^{\top}-vv^{\top}\|_{F}\leqslant\mathopen{}\mathclose{{}\left\|uu^{\top}-M}\right\|_{F}+\mathopen{}\mathclose{{}\left\|M-vv^{\top}}\right\|_{F}\leqslant
2\sqrt{2\varepsilon},$

so $\|uu^{\top}-vv^{\top}\|_{F}^{2}\leqslant 8\varepsilon$. Let
$\hat{u}=\frac{v}{\|v\|}$. To analyze the error of $\hat{u}$, note that

$\displaystyle\|uu^{\top}-vv^{\top}\|_{F}^{2}$
$\displaystyle=1+\|v\|^{4}-2\|v\|^{2}\mathopen{}\mathclose{{}\left\langle
uu^{\top},\hat{u}\hat{u}^{\top}}\right\rangle_{F}$ $\displaystyle\geqslant
1-\mathopen{}\mathclose{{}\left\langle
uu^{\top},\hat{u}\hat{u}^{\top}}\right\rangle_{F}$
$\displaystyle=\frac{1}{2}\mathopen{}\mathclose{{}\left\|uu^{\top}-\hat{u}\hat{u}^{\top}}\right\|_{F}^{2},$

where in the inequality we used that $1+x^{4}-2x^{2}y\geqslant 1-y$ for
$x\in\mathbb{R}$ and $0\leqslant y\leqslant 1$, with $x=\|v\|$ and $y=\langle
uu^{\top},\hat{u}\hat{u}^{\top}\rangle$. Then
${\|uu^{\top}-\hat{u}\hat{u}^{\top}\|_{F}^{2}\leqslant 16\varepsilon}$.
Therefore,

$\langle u,\hat{u}\rangle^{2}=\langle
uu^{\top},\hat{u}\hat{u}^{\top}\rangle_{F}=1-\frac{1}{2}\|uu^{\top}-\hat{u}\hat{u}^{\top}\|_{F}^{2}\geqslant
1-8\varepsilon.$

∎

###### Proof of Theorem 5.6.

The algorithm is to compute $M=\tilde{\mathbb{E}}vv^{\top}$, apply the
algorithm from Lemma 5.13 to $M$ in order to obtain a unit vector $\hat{u}$,
and return $\hat{u}$.

We now analyze the algorithm. We start by analyzing the properties of
$\tilde{\mathbb{E}}$ in more detail. Our first goal is to obtain the lower
bound $\tilde{\mathbb{E}}\langle u,v\rangle^{2}\geqslant 1-2\varepsilon$. We
start by proving the much weaker lower bound $\tilde{\mathbb{E}}\langle
u,v\rangle^{2}\geqslant 1-t\varepsilon$. Then, we use this lower bound to
prove an upper bound $\tilde{\mathbb{E}}\langle u,v\rangle^{2t}\leqslant
1-t(1-\tilde{\mathbb{E}}\langle u,v\rangle^{2})/2$. Comparing this result with
the given lower bound $\tilde{\mathbb{E}}\langle
u,v\rangle^{2t}\geqslant(1-\varepsilon)^{t}\geqslant 1-t\varepsilon$ leads to
the conclusion that $\tilde{\mathbb{E}}\langle u,v\rangle^{2}\geqslant
1-2\varepsilon$.

We proceed with the detailed proof of this fact. Recall that
$\tilde{\mathbb{E}}$ satisfies $\|v\|^{2}=1$ and $\langle
u,v\rangle^{2t}\geqslant(1-\varepsilon)^{t}$. We have that
$\\{\|v\|^{2}=1\\}\sststile{2}{}\\{0\leqslant\langle u,v\rangle^{2}\leqslant
1\\}$, where the lower bound is trivial and the upper bound is by Lemma A.8.
Therefore, $\tilde{\mathbb{E}}$ also satisfies $0\leqslant\langle
u,v\rangle^{2}\leqslant 1$.

By Lemma A.9 and using that $(1-\varepsilon)^{t}\geqslant 1-t\varepsilon$, we
have that

$\\{0\leqslant\langle u,v\rangle^{2}\leqslant 1,\langle
u,v\rangle^{2t}\geqslant(1-\varepsilon)^{t}\\}\sststile{2t}{}\\{\langle
u,v\rangle^{2}\geqslant 1-t\varepsilon\\}.$

By Lemma A.10 applied to $1-\langle u,v\rangle^{2}$ with
$C=\frac{1}{t^{2}\varepsilon}$, we also have that

$\\{1-t\varepsilon\leqslant\langle u,v\rangle^{2}\leqslant
1\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\langle
u,v\rangle^{2t}\leqslant 1-t(1-\langle u,v\rangle^{2})/2}\right\\}.$

Then

$1-t\varepsilon\leqslant\tilde{\mathbb{E}}\langle u,v\rangle^{2t}\leqslant
1-t(1-\tilde{\mathbb{E}}\langle u,v\rangle^{2})/2,$

so by rearranging, $\tilde{\mathbb{E}}\langle u,v\rangle^{2}\geqslant
1-2\varepsilon$.

Then $M=\tilde{\mathbb{E}}vv^{\top}$ satisfies

$\langle uu^{\top},M\rangle_{F}=u^{\top}Mu=\tilde{\mathbb{E}}\langle
u,v\rangle^{2}\geqslant 1-2\varepsilon.$

In addition, $M$ is symmetric positive-definite and

$\|M\|_{F}\leqslant\operatorname{Tr}(M)=\tilde{\mathbb{E}}\operatorname{Tr}(vv^{\top})=\tilde{\mathbb{E}}\|v\|^{2}=1.$

Therefore, $M$ satisfies the conditions of Lemma 5.13, and we are guaranteed
that $\hat{u}$ satisfies $\langle u,\hat{u}\rangle^{2}\geqslant
1-16\varepsilon$.

The given pseudo-expectation is of degree $O(t)$ over $d$ variables, so
representing it requires $d^{O(t)}$ space. Then we simply bound the time
complexity by $d^{O(t)}$, which dominates the other steps of the algorithm. ∎

###  5.3 Finite sample bounds

In Section 5.2 we assumed oracle access to $\mathbb{E}\bm{y}^{\otimes t}$.
However, our algorithm only has access to empirical moments. Lemma 5.14 gives
a sum-of-squares proof that that the empirical moments are in fact close to
the population moments. We defer the proof of the lemma to the appendix.

######  Lemma 5.14 (Closeness of empirical moments and population moments).

For $n\geqslant(p_{\mathrm{min}}^{-1}d)^{O(t)}\eta^{-2}\varepsilon^{-1}$, with
probability $1-\varepsilon$,

$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1}\right\\}\sststile{O(t)}{}\mathopen{}\mathclose{{}\left\\{\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant\mathbb{E}\langle\bm{y},v\rangle^{2t}+\eta}\right\\},$
$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1}\right\\}\sststile{O(t)}{}\mathopen{}\mathclose{{}\left\\{\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\geqslant\mathbb{E}\langle\bm{y},v\rangle^{2t}-\eta}\right\\}.$

###### Proof.

See Section A.7. ∎

###  5.4 Proof of Theorem 5.1

In the setting of Theorem 5.1 we only have access to empirical moments.
Theorem 5.15 and Theorem 5.16 adapt Theorem 5.5 and Theorem 5.4 to this
setting, respectively. Also recall that the goal of Theorem 5.1 is to return a
clustering, not only a unit vector close to $u$. Toward that goal, Theorem
5.17 shows that there exists an algorithm that, given a unit vector close to
$u$, computes such a clustering. We state and prove all of these theorems and
then combine them to prove Theorem 5.1.

######  Theorem 5.15 (Finite sample equivalent of Theorem 5.5).

Let

$n_{0}=\mathopen{}\mathclose{{}\left(\frac{1}{\sigma^{2}}}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}.$

Given a sample of size $n\geqslant n_{0}$ from the mixture, there exists an
algorithm that runs in time
$\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot n\cdot
d^{O(\log p_{\mathrm{min}}^{-1})}$ that computes with high probability two
peseudo-expectations $\tilde{\mathbb{E}}_{U}$ and $\tilde{\mathbb{E}}_{L}$ of
degree $O(\log p_{\mathrm{min}}^{-1})$ over a variable $v\in\mathbb{R}^{d}$
such that the following holds. Let $s=\lceil\log p_{\mathrm{min}}^{-1}\rceil$,
let $t=5000s$, and let $\tau=\frac{800e}{C_{sep}k^{2}}$. Then
$\tilde{\mathbb{E}}_{U}\|v\|^{2}=1$, $\tilde{\mathbb{E}}_{L}\|v\|^{2}=1$, and:

  * •

If $\sigma^{2}\geqslant\tau$, then $\tilde{\mathbb{E}}_{U}\langle
u,v\rangle^{2s}\geqslant(1-\tau)^{s}$.

  * •

If $\sigma^{2}<\tau$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\geqslant(4es)^{s}$, then
$\tilde{\mathbb{E}}_{U}\langle u,v\rangle^{2s}\geqslant(1-\sigma^{2})^{s}$.

  * •

If $\sigma^{2}<0.001$ and
$\mathbb{E}\langle\bm{\mu},u\rangle^{2s}\leqslant(100s)^{s}$, then
$\tilde{\mathbb{E}}_{L}\langle u,v\rangle^{2t}\geqslant(1-20\sigma^{2})^{t}$.

###### Proof.

The algorithm is the same as that in the proof of Theorem 5.5, except that in
step (2) and step (3) of the algorithm the constraints that the pseudo-
expectations are required to satisfy are
$\\{\|v\|^{2}=1,\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant T_{U}\\}$
and $\\{\|v\|^{2}=1,\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant
T_{L}\\}$, respectively.

By Lemma 5.14, for $n\geqslant n_{0}$, we have that with high probability

$\\{\|v\|^{2}=1,\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant
T_{U}\\}\sststile{2s}{}\\{\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2s}\geqslant
T_{U}-\sigma^{2}/(100M)\\},$
$\\{\|v\|^{2}=1,\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant
T_{L}\\}\sststile{2s}{}\\{\hat{\mathbb{E}}\langle\bm{y},v\rangle^{2t}\leqslant
T_{L}+\sigma^{2}/10000\\}.$

These errors, combined with the errors from the binary search resolution, are
still within the amount tolerated by Lemma 5.11 and Lemma 5.12, respectively,
so the same guarantees hold.

The number of steps required by the binary search is the same as in Theorem
5.5. For each step of the binary search, we compute a pseudo-expectation of
degree $O(\log p_{\mathrm{min}}^{-1})$ over $d$ variables, and each constraint
requires summing over the $n$ samples, so the time required is $n\cdot
d^{O(\log p_{\mathrm{min}}^{-1})}$. Therefore the time complexity is

$O\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}+\log^{2}p_{\mathrm{min}}^{-1}}\right)\cdot
n\cdot d^{O(\log
p_{\mathrm{min}}^{-1})}=\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot
n\cdot d^{O(\log p_{\mathrm{min}}^{-1})}.$

∎

######  Theorem 5.16 (Finite sample equivalent of Theorem 5.4).

Let

$n_{0}=\mathopen{}\mathclose{{}\left(\frac{1}{\sigma^{2}}}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}.$

Given a sample of size $n\geqslant n_{0}$ from the mixture, there exists an
algorithm with time complexity
$\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot n\cdot
d^{O(\log p_{\mathrm{min}}^{-1})}$ that outputs with high probability a unit
vector $\hat{u}\in\mathbb{R}^{d}$ such that $\langle
u,\hat{u}\rangle^{2}=1-320\min(\sigma^{2},\frac{1}{k})$.

###### Proof.

The algorithm is the same as that in the proof of Theorem 5.4, with two
exceptions:

  * •

In step (1) of the algorithm, we run the algorithm from Theorem 5.15 instead
of the algorithm from Theorem 5.5. The pseudo-expectations
$\tilde{\mathbb{E}}_{L}$ and $\tilde{\mathbb{E}}_{U}$ satisfy the same
guarantees.

  * •

In step (3) of the algorithm, we check if
$\hat{\mathbb{E}}\langle\bm{y},\hat{u}_{U}\rangle^{2s}\geqslant(50s)^{s}$
instead of
$\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s}\geqslant(50s)^{s}$. By Lemma
5.14, for $n\geqslant n_{0}$, with high probability the difference between the
two moments is less than $1$. Then, if
$\hat{\mathbb{E}}\langle\bm{y},\hat{u}_{U}\rangle^{2s}\geqslant(50s)^{s}$, we
also have
$\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s}\geqslant(50s)^{s}-1$, and if
$\hat{\mathbb{E}}\langle\bm{y},\hat{u}_{U}\rangle^{2s}<(50s)^{s}$, we also
have $\mathbb{E}\langle\bm{y},\hat{u}_{U}\rangle^{2s}<(50s)^{s}+1$. It is easy
to verify that the analysis in Theorem 5.4 is still correct with these
slightly weaker bounds.

Therefore, the same guarantees hold as in Theorem 5.4. The time complexity of
the algorithm is dominated by the time to run the algorithm from Theorem 5.15.

∎

######  Theorem 5.17 (Clustering algorithm).

For some $C>0$, suppose that a unit vector $\hat{u}\in\mathbb{R}^{d}$ is known
such that $\langle u,\hat{u}\rangle^{2}\geqslant
1-C\min(\sigma^{2},\frac{1}{k})$. Suppose that $C_{sep}/C$ is larger than some
universal constant. Then, given a sample of size
$n\geqslant(p_{\mathrm{min}}^{-1})^{O(1)}$ from the mixture, there exists an
algorithm that runs in time $n^{O(\log p_{\mathrm{min}}^{-1})}$ and returns
with high probability a partition of $[n]$ into $k$ sets $C_{1},...,C_{k}$
such that, if the true clustering of the samples is $S_{1},...,S_{k}$, then
there exists a permutation $\pi$ of $[k]$ such that

$1-\frac{1}{n}\sum_{i=1}^{k}|C_{i}\cap
S_{\pi(i)}|\leqslant\mathopen{}\mathclose{{}\left(\frac{p_{\mathrm{min}}}{k}}\right)^{O(1)}.$

###### Proof.

Our algorithm runs the algorithm from Theorem 5.1 of [HL18] with some
$t=O(\log p_{\mathrm{min}}^{-1})$ large enough on input samples $\langle
y_{1},\hat{u}\rangle/(\sqrt{2(C+1)\sigma^{2}})$, …, $\langle
y_{n},\hat{u}\rangle/(\sqrt{2(C+1)\sigma^{2}})$, and returns the clustering
that this algorithm computes as an intermediate step.

We now analyze the algorithm. Note that
$\langle\bm{y},\hat{u}\rangle/\sqrt{2(C+1)\sigma^{2}}$ is distributed
according to a one-dimensional mixture of Gaussians in which all the
components have the same variance
$\hat{u}^{\top}\Sigma\hat{u}/(2(C+1)\sigma^{2})$. We have that

$\hat{u}^{\top}\Sigma\hat{u}=1-(1-\sigma^{2})\langle
u,\hat{u}\rangle^{2}=1-\langle u,\hat{u}\rangle^{2}+\sigma^{2}\langle
u,\hat{u}\rangle^{2}\leqslant C\sigma^{2}+\sigma^{2}=(C+1)\sigma^{2}.$

Therefore, the variance $\hat{u}^{\top}\Sigma\hat{u}/(2(C+1)\sigma^{2})$ is
upper bounded by $1/2$. For the guarantees of the algorithm from [HL18] to
hold, we further need to show that the mixture has large separation between
the means of the components. Note that the mean corresponding to $\mu_{i}$ in
the original mixture becomes
$\langle\mu_{i},\hat{u}\rangle/\sqrt{2(C+1)\sigma^{2}}$ in the new mixture.
For $i\neq j$, we have that

$\displaystyle(\langle\mu_{i}-\mu_{j},\hat{u}\rangle)^{2}$
$\displaystyle=\langle\langle\mu_{i},u\rangle u-\langle\mu_{j},u\rangle
u,\hat{u}\rangle^{2}=\langle
u,\hat{u}\rangle^{2}\langle\mu_{i}-\mu_{j},u\rangle^{2}$
$\displaystyle\geqslant\langle u,\hat{u}\rangle^{2}C_{sep}(u^{\top}\Sigma
u)\log p_{\mathrm{min}}^{-1}=\langle u,\hat{u}\rangle^{2}C_{sep}\sigma^{2}\log
p_{\mathrm{min}}^{-1}$ $\displaystyle\geqslant\frac{C_{sep}}{2}\sigma^{2}\log
p_{\mathrm{min}}^{-1}$

where in the last inequality we used that $\langle
u,\hat{u}\rangle^{2}\geqslant 1-C/k\geqslant 1/2$. Then

$\displaystyle\mathopen{}\mathclose{{}\left(\frac{\langle\mu_{i}-\mu_{j},\hat{u}\rangle}{\sqrt{2(C+1)\sigma^{2}}}}\right)^{2}$
$\displaystyle\geqslant\frac{C_{sep}}{4(C+1)}\log p_{\mathrm{min}}^{-1}.$

For $C_{sep}/C$ larger than some universal constant, the separation
coefficient $C_{sep}/(4(C+1))$ is large enough for the guarantees of Theorem
5.1 of [HL18] to hold meaningfully with $t=O(\log p_{\mathrm{min}}^{-1})$.
Then this algorithm computes a clustering with the stated guarantees. The
algorithm requries $n\geqslant(p_{\mathrm{min}}^{-1})^{O(1)}$ and runs in time
$n^{O(\log p_{\mathrm{min}}^{-1})}$. ∎

###### Proof of Theorem 5.1.

Run the algorithm from Theorem 5.16 to obtain a unit vector $\hat{u}$ that
satisfies $\langle u,\hat{u}\rangle^{2}\geqslant
1-320\min(\sigma^{2},\frac{1}{k})$. Then run the clustering algorithm from
Theorem 5.17 using this unit vector $\hat{u}$. For $C_{sep}/320$ larger than
some universal constant, this algorithm is guaranteed to return a clustering
with the stated guarantees.

The time complexity from Theorem 5.16 is
$\mathopen{}\mathclose{{}\left(\log\frac{1}{\sigma^{2}}}\right)\cdot n\cdot
d^{O(\log p_{\mathrm{min}}^{-1})}$ and the time complexity from Theorem 6.1 is
$n^{O(\log p_{\mathrm{min}}^{-1})}$. We assume
$n\geqslant\mathopen{}\mathclose{{}\left(\frac{1}{\sigma^{2}}}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}$. Therefore, the time complexity is dominated by the
time to run the clustering algorithm from Theorem 6.1. ∎

##  6 Colinear means

In this section we remove the isotropic position assumption from the model in
Section 5. Our strategy is straightfoward: we first put the mixture in
isotropic position and then run the algorithm from Theorem 5.1. The technical
challenge is that we can only put the mixture in approximate isotropic
position. Then, we show that the guarantees of Theorem 5.1 continue to hold
with approxiamte isotropic position, albeit with a sample complexity that
depends on the condition number of the covariance matrix of the mixture.

##### Setting.

We consider a mixture of $k$ Gaussian distributions
$N(\mu_{i}^{0},\Sigma^{0})$ with mixing weights $p_{i}$ for $i=1,...,k$, where
$\mu_{i}^{0}\in\mathbb{R}^{d}$, $\Sigma^{0}\in\mathbb{R}^{d\times d}$ is
positive definite, and $p_{i}\geqslant 0$ and $\sum_{i=1}^{k}p_{i}=1$. Let
$p_{\mathrm{min}}=\min_{i}p_{i}$.

The distribution also satisfies mean separation and mean colinearity:

  * •

Mean separation: for some $C_{sep}>0$ and for all $i\neq j$,

$\mathopen{}\mathclose{{}\left\|\mathopen{}\mathclose{{}\left(\Sigma^{0}}\right)^{-1/2}(\mu_{i}^{0}-\mu_{j}^{0})}\right\|^{2}\geqslant
C_{sep}\log p_{\mathrm{min}}^{-1}.$

  * •

Mean colinearity: for some vector $\mu_{base}^{0}\in\mathbb{R}^{d}$ and some
unit vector $u^{0}\in\mathbb{R}^{d}$ and for all $i$,

$\mu_{i}^{0}=\mu_{base}^{0}+\langle\mu_{i}^{0}-\mu_{base}^{0},u^{0}\rangle
u^{0}.$

Also define, for $\bm{y}^{0}$ distributed according to the mixture,

$\sigma^{2}=\frac{(u^{0})^{\top}\operatorname{cov}(\bm{y}^{0})^{-1}u^{0}}{(u^{0})^{\top}(\Sigma^{0})^{-1}u^{0}}.$

As shown in Lemma 6.3, $\sigma^{2}$ has the same interpretation as in Section
5: it is equal to the variance of the components in the direction of the means
after an isotropic position transformation.

######  Theorem 6.1 (Colinear means algorithm).

Consider the Gaussian mixture model defined above, with $C_{sep}$ larger than
some universal constant. For $\bm{y}^{0}$ distributed according to the
mixture, let

$n_{0}=\mathopen{}\mathclose{{}\left(\frac{1}{\sigma^{2}}\cdot\|\operatorname{cov}(\bm{y}^{0})\|\cdot\|\operatorname{cov}(\bm{y}^{0})^{-1}\|}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(\log
p_{\mathrm{min}}^{-1})}.$

Given a sample of size $n\geqslant n_{0}$ from the mixture, there exists an
algorithm that runs in time $n^{O(\log p_{\mathrm{min}}^{-1})}$ and returns
with high probability a partition of $[n]$ into $k$ sets $C_{1},...,C_{k}$
such that, if the true clustering of the samples is $S_{1},...,S_{k}$, then
there exists a permutation $\pi$ of $[k]$ such that

$1-\frac{1}{n}\sum_{i=1}^{k}|C_{i}\cap
S_{\pi(i)}|\leqslant\mathopen{}\mathclose{{}\left(\frac{p_{\mathrm{min}}}{k}}\right)^{O(1)}.$

We introduce some further notation for this section. Let $\bm{y}^{0}$ be
distributed according to the mixture. We specify the model as
$\bm{y}^{0}=\bm{\mu}^{0}+\bm{w}^{0}$, where $\bm{\mu}^{0}$ takes value
$\mu_{i}$ with probability $p_{i}$ and $\bm{w}^{0}\sim N(0,\Sigma^{0})$, with
$\bm{\mu}^{0}$ and $\bm{w}^{0}$ independent of each other.

###  6.1 Isotropic position transformation

In this section we argue that, if we put the mixture in exact isotropic
position, the conditions of Theorem 5.1 are satisfied and we can simply run
that algorithm.

Assume acces to $\mathbb{E}\bm{y}^{0}$ and to an invertible matrix
$W\in\mathbb{R}^{d\times d}$ such that
$W\operatorname{cov}(\bm{y}^{0})W^{\top}=I_{d}$. Then, define the random
variable $\bm{y}$ by the affine transformation
$\bm{y}=W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0})$. The distribution of $\bm{y}$ is
in isotropic position: it has mean $0$ and covariance matrix $I_{d}$.
Furthermore, Lemma 6.2 shows that $\bm{y}$ is a mixture of Gaussians in which
the components are affine transformed versions of the original components, and
that the mixture continues to satisfy mean separation and mean colinearity.
Then, the conditions of Theorem 5.1 are satisfied. Therefore, if the original
input samples are $y^{0}_{1},...,y^{0}_{n}$, we can simply run that algorithm
on input samples $W(y^{0}_{1}-\mathbb{E}y^{0})$, …,
$W(y^{0}_{n}-\mathbb{E}y^{0})$101010It is straightforward that if sample
$y^{0}_{i}$ comes from the $i$-th component in the original mixture then
$W(y^{0}_{i}-\mathbb{E}y^{0})$ continues to come from the $i$-th component in
the affine transformed mixture..

We define now some variables used to state Lemma 6.2. Recall that
$\bm{y}=W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0})$. Define
$\mu_{i}=W(\mu_{i}^{0}-\mathbb{E}\bm{y}^{0})$, $\Sigma=W\Sigma^{0}W^{\top}$,
and $u=\frac{Wu^{0}}{\|Wu^{0}\|v}$. Also define random variables
$\bm{\mu}=W(\bm{\mu}^{0}-\mathbb{E}\bm{y}^{0})$ and $\bm{w}=W\bm{w}^{0}$.

######  Lemma 6.2 (Model after isotropic position transformation).

The random variable $\bm{y}$ is distributed according to a mixture of $k$
Gaussian distributions $N(\mu_{i},\Sigma)$ with mixing weights $p_{i}$ for
$i=1,...,k$, where $\Sigma$ is positive definite. Alternatively, we specify
the model as $\bm{y}=\bm{\mu}+\bm{w}$, with $\bm{\mu}$ and $\bm{w}$
independent of each other. Furthermore, for all $i\neq j$ we have mean
separation

$\mathopen{}\mathclose{{}\left\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})}\right\|^{2}\geqslant
C_{sep}\log p_{\mathrm{min}}^{-1}$

and for all $i$ we have mean colinearity

$\mu_{i}=\langle\mu_{i},u\rangle u.$

###### Proof.

Recall that $\bm{y}^{0}=\bm{\mu}^{0}+\bm{w}^{0}$. Then
$W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0})=W(\bm{\mu}^{0}-\mathbb{E}\bm{y}^{0})+W\bm{w}^{0}$,
so $\bm{y}=\bm{\mu}+\bm{w}$. Also note that $\bm{\mu}$ takes value
$W(\mu_{i}^{0}-\mathbb{E}\bm{y}^{0})=\mu_{i}$ with probability $p_{i}$ and
$\bm{w}\sim N(0,W\Sigma^{0}W^{\top})=N(0,\Sigma)$. Therefore, $\bm{y}$ is
distributed according to a mixture of $k$ Gaussian distributions
$N(\mu_{i},\Sigma)$ with mixing weights $p_{i}$.

To show that $\Sigma$ is positive definite, we note that for any vector
$v\in\mathbb{R}^{d}$ with $v\neq 0$ we have that

$v^{\top}\Sigma
v=v^{\top}W\Sigma^{0}W^{\top}v=(W^{\top}v)^{\top}\Sigma^{0}W^{\top}v>0,$

where we used that $W^{\top}v\neq 0$ because $W$ is invertible, after which we
used that $\Sigma^{0}$ is positive definite.

We prove now mean colinearity and mean separation. We start with mean
colinearity. Recall that
$\mu_{i}^{0}=\mu_{base}^{0}+\langle\mu_{i}^{0}-\mu_{base}^{0},u^{0}\rangle
u^{0}$. Then, using that $\mathbb{E}\bm{y}^{0}=\mathbb{E}\bm{\mu}^{0}$,

$\mathbb{E}\bm{y}^{0}=\mu_{base}^{0}+\langle\mathbb{E}\bm{y}^{0}-\mu_{base}^{0},u^{0}\rangle
u^{0},$

so

$\mu_{i}=W(\mu_{i}^{0}-\mathbb{E}\bm{y}^{0})=W(\langle\mu_{i}^{0}-\mathbb{E}\bm{y}^{0},u^{0}\rangle
u^{0})=\langle\mu_{i}^{0}-\mathbb{E}\bm{y}^{0},u^{0}\rangle Wu^{0}.$

Then, using that $u=\frac{Wu^{0}}{\|Wu^{0}\|}$,

$\displaystyle\langle\mu_{i},u\rangle u$
$\displaystyle=\mathopen{}\mathclose{{}\left\langle\langle\mu_{i}^{0}-\mathbb{E}\bm{y}^{0},u^{0}\rangle
Wu^{0},\frac{Wu^{0}}{\|Wu^{0}\|}}\right\rangle\frac{Wu^{0}}{\|Wu^{0}\|}$
$\displaystyle=\langle\mu_{i}^{0}-\mathbb{E}\bm{y}^{0},u^{0}\rangle\mathopen{}\mathclose{{}\left\langle\frac{Wu^{0}}{\|Wu^{0}\|},\frac{Wu^{0}}{\|Wu^{0}\|}}\right\rangle
Wu^{0}$ $\displaystyle=\langle\mu_{i}^{0}-\mathbb{E}\bm{y}^{0},u^{0}\rangle
Wu^{0}$ $\displaystyle=\mu_{i},$

which proves mean colinearity. For mean separation, we have that

$\displaystyle\mathopen{}\mathclose{{}\left\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})}\right\|^{2}$
$\displaystyle=\mathopen{}\mathclose{{}\left\|(W\Sigma^{0}W^{\top})^{-1/2}W(\mu_{i}^{0}-\mu_{j}^{0})}\right\|^{2}$
$\displaystyle=(\mu_{i}^{0}-\mu_{j}^{0})^{\top}W^{\top}(W\Sigma^{0}W^{\top})^{-1}W(\mu_{i}^{0}-\mu_{j}^{0})$
$\displaystyle=(\mu_{i}^{0}-\mu_{j}^{0})^{\top}W^{\top}(W^{\top})^{-1}(\Sigma^{0})^{-1}W^{-1}W(\mu_{i}^{0}-\mu_{j}^{0}))$
$\displaystyle=(\mu_{i}^{0}-\mu_{j}^{0})^{\top}(\Sigma^{0})^{-1}(\mu_{i}^{0}-\mu_{j}^{0})$
$\displaystyle=\mathopen{}\mathclose{{}\left\|(\Sigma^{0})^{-1/2}(\mu_{i}^{0}-\mu_{j}^{0})}\right\|^{2}$
$\displaystyle\geqslant C_{sep}\log p_{\mathrm{min}}^{-1}.$

∎

######  Lemma 6.3.

$\sigma^{2}=\frac{(u^{0})^{\top}\operatorname{cov}(\bm{y}^{0})^{-1}u^{0}}{(u^{0})^{\top}(\Sigma^{0})^{-1}u^{0}}=u^{\top}\Sigma
u.$

###### Proof.

For the purposes of this proof, we define $\sigma^{2}=u^{\top}\Sigma u$ as in
Section 5.1 and prove that it also matches the definition in this section.

We have, as in the proof of Lemma 6.2, that

$\|(\Sigma^{0})^{-1/2}(\mu_{i}^{0}-\mu_{j}^{0})\|^{2}=\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})\|^{2}.$

For the left-hand side, we have that

$\displaystyle\|(\Sigma^{0})^{-1/2}(\mu_{i}^{0}-\mu_{j}^{0})\|^{2}$
$\displaystyle=\|(\Sigma^{0})^{-1/2}\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle
u^{0}\|^{2}$
$\displaystyle=\|(\Sigma^{0})^{-1/2}u^{0}\|^{2}\cdot\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle^{2}.$

For the right-hand side, using from Lemma 5.2 that
$\Sigma=I_{d}-(1-\sigma^{2})uu^{\top}$, we have that

$\displaystyle\|\Sigma^{-1/2}(\mu_{i}-\mu_{j})\|^{2}$
$\displaystyle=\mathopen{}\mathclose{{}\left\|\mathopen{}\mathclose{{}\left(I_{d}-(1-\sigma^{2})uu^{\top}}\right)^{-1/2}(\mu_{i}-\mu_{j}))}\right\|^{2}$
$\displaystyle=\mathopen{}\mathclose{{}\left\|\mathopen{}\mathclose{{}\left(I_{d}-(1-\sigma^{2})\frac{(Wu^{0})(Wu^{0})^{\top}}{\|Wu^{0}\|^{2}}}\right)^{-1/2}\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle
Wu^{0}}\right\|^{2}$
$\displaystyle=\mathopen{}\mathclose{{}\left\|\mathopen{}\mathclose{{}\left(I_{d}+\mathopen{}\mathclose{{}\left(\frac{1}{\sigma}-1}\right)\frac{(Wu^{0})(Wu^{0})^{\top}}{\|Wu^{0}\|^{2}}}\right)Wu^{0}}\right\|^{2}\cdot\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle^{2}$
$\displaystyle=\frac{1}{\sigma^{2}}\cdot\|Wu^{0}\|^{2}\cdot\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle^{2}.$

Therefore

$\|(\Sigma^{0})^{-1/2}u^{0}\|^{2}\cdot\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle^{2}=\frac{1}{\sigma^{2}}\cdot\|Wu^{0}\|^{2}\cdot\langle\mu_{i}^{0}-\mu_{j}^{0},u^{0}\rangle^{2},$

so

$\sigma^{2}=\frac{\|Wu^{0}\|^{2}}{\|(\Sigma^{0})^{-1/2}u^{0}\|^{2}}=\frac{(u^{0})^{\top}\operatorname{cov}(\bm{y}^{0})^{-1}u^{0}}{(u^{0})^{\top}(\Sigma^{0})^{-1}u^{0}},$

where we used that, by Lemma A.24, $W=Q\operatorname{cov}(\bm{y}^{0})^{-1/2}$
for an orthogonal matrix $Q$, so
$\|Wu^{0}\|=\|\operatorname{cov}(\bm{y}^{0})^{-1/2}u^{0}\|$. ∎

###  6.2 Finite sample isotropic position transformation

Without access to $\mathbb{E}\bm{y}^{0}$ and to $W$, we apply the isotropic
position transformation with $\hat{\mathbb{E}}\bm{y}^{0}$ and some matrix
$\hat{W}\in\mathbb{R}^{d\times d}$ defined as follows. Let the singular value
decomposition of $\widehat{\operatorname{cov}}(\bm{y}^{0})$ be
$\hat{U}\hat{\Lambda}\hat{U}^{\top}$. Then define $\hat{W}$ and $W$ as

$\hat{W}=(\hat{U}^{\top}\widehat{\operatorname{cov}}(\bm{y}^{0})\hat{U})^{-1/2}\hat{U}^{\top},\quad
W=(\hat{W}\operatorname{cov}(\bm{y}^{0})\hat{W}^{\top})^{-1/2}\hat{W}.$ (6.1)

This choice is analogous to that in Appendix C in [HK13]. By Lemma A.22, we
have that
$\hat{W}\widehat{\operatorname{cov}}(\bm{y}^{0})\hat{W}^{\top}=I_{d}$ and
$W\operatorname{cov}(\bm{y}^{0})W^{\top}=I_{d}$. Hence, $\hat{W}$ corresponds
to an isotropic position transformation for the empirical covariance matrix,
and $W$ to one for the population covariance matrix. In our algorithm, we will
apply the approximate isotropic position transformation to input samples
$y_{1}^{0},...,y_{n}^{0}$ as
$\hat{W}(y_{1}^{0}-\hat{\mathbb{E}}y^{0}),...,\hat{W}(y_{n}^{0}-\hat{\mathbb{E}}y^{0})$.

###  6.3 Finite sample bounds

Lemma 6.4 gives a sum-of-squares proof that the empirical moments of the
mixture with approximate isotropic position transformation are close to the
population moments of the mixture with exact isotropic position
transformation. This lemma is supported by Lemma 6.5 and Lemma 6.6, which
prove that the moments do not change much due to the use of
$\hat{\mathbb{E}}\bm{y}^{0}$ and $\hat{W}$, respectively.

Additionally, for arbitrary unit vectors $v$, Lemma 6.7 proves that
$\langle\hat{W}(\mu_{i}^{0}-\mu_{j}^{0}),v\rangle$ is close to $\langle
W(\mu_{i}^{0}-\mu_{j}^{0}),v\rangle$ and Lemma 6.8 proves that
$v^{\top}\hat{W}(\Sigma^{0})^{1/2}$ is close to $v^{\top}W(\Sigma^{0})^{1/2}$.
These facts are used in the proof of Theorem 6.1 to argue that the clustering
algorithm is correct.

######  Lemma 6.4 (Closeness of empirical approximate isotropic position
moments and population exact isotropic position moments).

Let $\eta<0.001$. For

$n\geqslant\mathopen{}\mathclose{{}\left(\|\operatorname{cov}(\bm{y}^{0})\|\cdot\|\operatorname{cov}(\bm{y}^{0})^{-1}\|}\right)^{O(1)}\cdot(p_{\mathrm{min}}^{-1}d)^{O(t)}\eta^{-2}\varepsilon^{-1},$

with probability $1-\varepsilon$,

$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1}\right\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\hat{\mathbb{E}}\bm{y}^{0}),v\rangle^{2t}\leqslant(1+\eta)\cdot\mathbb{E}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}+\eta}\right\\},$
$\mathopen{}\mathclose{{}\left\\{\|v\|^{2}=1}\right\\}\sststile{2t}{}\mathopen{}\mathclose{{}\left\\{\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\hat{\mathbb{E}}\bm{y}^{0}),v\rangle^{2t}\geqslant(1-\eta)\cdot\mathbb{E}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}-\eta}\right\\}.$

###### Proof.

Select $n$ such that the results of Lemma 6.5, Lemma 6.6, and Lemma 5.14 hold
each with probability $1-\varepsilon/3$. Then, overall, all three results hold
with probability $1-\varepsilon$. Then we have with probability
$1-\varepsilon$ that

$\displaystyle\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\hat{\mathbb{E}}\bm{y}^{0}),v\rangle^{2t}$
$\displaystyle\quad\leqslant(1+\eta/10)\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}+\eta/10$
$\displaystyle\quad\leqslant(1+\eta/10)\mathopen{}\mathclose{{}\left((1+\eta/10)\hat{\mathbb{E}}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}+\eta/10}\right)+\eta/10$
$\displaystyle\quad\leqslant(1+\eta/10)\mathopen{}\mathclose{{}\left((1+\eta/10)\mathopen{}\mathclose{{}\left((1+\eta/10)\mathbb{E}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}+\eta/10}\right)+\eta/10}\right)+\eta/10$
$\displaystyle\quad\leqslant(1+\eta)\mathbb{E}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}+\eta.$

and

$\displaystyle\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\hat{\mathbb{E}}\bm{y}^{0}),v\rangle^{2t}$
$\displaystyle\quad\geqslant(1-\eta/10)\hat{\mathbb{E}}\langle\hat{W}(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}-\eta/10$
$\displaystyle\quad\geqslant(1-\eta/10)\mathopen{}\mathclose{{}\left((1-\eta/10)\hat{\mathbb{E}}\langle
W(\bm{y}^{0}-\mathbb{E}\bm{y}^{0}),v\rangle^{2t}-\eta/10}\right)-\eta/10$
$\displaystyle\quad\geqslant(1-\eta/10)\mathopen{}\mathclose{{}\left((1-\eta/10)\mathopen{}\mathclose{{}\left((1-\eta/10)\mathbb{E}\langle