Sample-Efficient Private Learning of Mixtures of Gaussians

Authors: hassan  ashtiani
mahbod  majid
shyam  narayanan
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 08:52:47

Abstract: We study the problem of learning mixtures of Gaussians with approximate differential privacy. We prove that roughly $kd^2 + k^{1.5} d^{1.75} + k^2 d$ samples suffice to learn a mixture of $k$ arbitrary $d$-dimensional Gaussians up to low total variation distance, with differential privacy. Our work improves over the previous best result (which required roughly $k^2 d^4$ samples) and is provably optimal when $d$ is much larger than $k^2$. Moreover, we give the first optimal bound for privately learning mixtures of $k$ univariate (i.e., $1$-dimensional) Gaussians. Importantly, we show that the sample complexity for learning mixtures of univariate Gaussians is linear in the number of components $k$, whereas the previous best sample complexity was quadratic in $k$. Our algorithms utilize various techniques, including the inverse sensitivity mechanism, sample compression for distributions, and methods for bounding volumes of sumsets.

Main paper: 
 introduction: 
 Learning Gaussian Mixture Models (GMMs) is one of the most fundamental problems in algorithmic statistics. Gaussianity is a common data assumption, and the setting of Gaussian mixture models is motivated by heterogeneous data that can be split into numerous clusters, where each cluster follows a Gaussian distribution. Learning mixture models is among the most important problems in machine learning [Bis06] , and is at the heart of several unsupervised and semi-supervised machine learning models. The study of Gaussian mixture models has had numerous scientific applications dating back to the 1890s [Pea94] , and is a crucial tool in modern data analysis techniques in a variety of fields, including bioinformatics [LKWB22] , anomaly detection [ZSM + 18], and handwriting analysis [Bis06] .In this work, we study the problem of learning a GMM from samples. We focus on the density estimation setting, where the goal is to learn the overall mixture distribution up to low total variation distance. Unlike the parameter estimation setting for GMMs, density estimation can be done even without any boundedness or separation assumptions on the parameters of the components. In fact, it is known that mixtures of k Gaussians in d-dimensions can be learned up to total variation distance α using O(kd 2 /α 2 ) samples [ABH + 18].Ensuring data privacy has emerged as an increasingly important challenge in modern data analysis and statistics. Differential privacy (DP) [DMNS06] is a rigorous way of defining privacy, and is considered to be the gold standard both in theory and practice, with deployments by Apple [Tea17] , Google [EPK14] , Microsoft [DKY17] , and the US Census Bureau [DLS + 17]. As is the case for many data analysis tasks, standard algorithms for learning GMMs leak potentially sensitive information about the individuals who contributed data. This raises the question of whether we can do density estimation for GMMs under the constraint of differential privacy.Private density estimation for GMMs with unrestricted Gaussian components is a challenging task. In fact, privately learning a single unrestricted Gaussian has been the subject of multiple recent studies [AAK21, KMS + 22b, AL22, KMV22, AKT + 23, HKMN23]. Private learning of GMMs is significantly more challenging, because even without privacy constraints, parameter estimation for GMMs requires exponentially many samples in terms of the number of components [MV10] . Therefore, it is not clear how to use the typical recipe of "adding noise" to the estimated parameters or "privately choosing" from the finite-dimensional space of parameters. Consequently, the only known sample complexity bounds for privately learning unrestricted GMMs are loose [AAL24b, AAL21] .Let us first formally define the problem of learning GMMs. We represent a GMM D = k i=1 w i N (µ i , Σ i ) by its parameters, namely {(w i , µ i , Σ i )} k i=1 , where w i ≥ 0, i w i = 1, µ i ∈ R d , and Σ i is a positive definite matrix. In the following, a GMM learning algorithm A receives a set of data points in R d and outputs a (representation of) a GMM. The total variation distance between two distributions is d TV ( D, D) = 12 R d |D(x) -D(x)|dx 2 . Definition 1.1 (Learning GMMs). For α, β ∈ (0, 1), we say A learns GMMs with n samples up to accuracy α and failure probability β if for every GMM D, given samples X 1 , . . . , X n i.i.d.∼ D, it outputs (a representation of) a GMM D such that d TV ( D, D) ≤ α with probability at least 1 -β. α and β are called the accuracy and failure probability, respectively. For clarity of presentation, we will typically fix the value of β (e.g., β = 1/3). The above definition does not enforce the constraint of differential privacy. The following definitions formalizes (approximate) differential privacy. Definition 1.2 (Differential Privacy (DP) [DMNS06, DKM + 06]). Let ε, δ ≥ 0. A randomized algorithm A : X n → O is said to be (ε, δ)-differentially private ((ε, δ)-DP) if for any two neighboring datasets X, X ′ ∈ X n and any measurable subset O ⊂ O,P[A(X ′ ) ∈ O] ≤ e ε • P[A(X) ∈ O] + δ.If the GMM learner A of Definition 1.1 is (ε, δ)-DP, we say that A privately learns GMMs. Formally, we have the following definition. Definition 1.3 (Privately learning GMMs). Fix the number of samples n, dimension d, and number of mixture components k. For α, β ∈ (0, 1) and ε, δ ≥ 0, a randomized algorithm A, that takes as input X 1 , . . . , X n ∈ R d , (ε, δ)-privately learns GMMs up to accuracy α and failure probability β, if:1. For any GMM D that is a mixture of up to k Gaussians in d dimensions, if X = {X 1 , . . . , X n } i.i.d.∼ D, A(X 1 , . . . , X n ) outputs a GMM D such that d TV ( D, D) ≤ α with probability at least 1 -β (over the randomness of the data X and the algorithm A).2. For any neighboring datasets X, X ′ ∈ (R d ) n (not necessarily drawn from any GMM) and any measurable subsetO ⊂ O, P[A(X ′ ) ∈ O] ≤ e ε • P[A(X) ∈ O] + δ.Finally, we assume a default value for β of 1/3, meaning that if not stated, the failure probability β is assumed to equal 1/3.Our main goal in this paper is to understand the number of samples (as a function of the dimension d, the number of mixture components k, the accuracy α, and the privacy parameters ε, δ) that are needed to privately and accurately learn the GMM up to low total variation distance.In this work, we provide improved sample complexity bounds for privately learning mixtures of arbitrary Gaussians, improving over previous work of [AAL21, AAL24b] . Moreover, our sample complexity bounds are optimal in certain regimes, when the dimension is either 1 or a sufficiently large polynomial in k and log 1 δ . For general dimension d, we prove the following theorem. Theorem 1.4. For any α, ε, δ ∈ (0, 1), k, d ∈ N, there exists an inefficient (ε, δ)-DP algorithm that can learn a mixture of k arbitrary full-dimensional Gaussians in d dimensions up to accuracy α, using the following number of samples:n = O kd 2 α 2 +kd 2 + d 1.75 k 1.5 log 0.5 (1/δ) + k 1.5 log 1. 5 (1/δ) αε + k 2 d α .Notably, the mixing weights and the means can be arbitrary and the covariances of the Gaussians can be arbitrarily poorly conditioned, as long as the covariances are non-singular 3 .We remark that we omit the dependence on β (and assume by default a failure probability of 1/3). However, it is well-known that one can obtain failure probability β with only a multiplicative O(log 1/β) blowup in sample complexity, in a black-box fashion 4 . In fact, our analysis can yield even better dependencies on β in some regimes, though to avoid too much complication, we do not analyze this.For reasonably large dimension, i.e., d ≥ k 2 log 2 (1/δ), this can be simplified to Õ kd 2 α 2 + kd 2 αε , which is in fact optimal (see Theorem 1.6). Hence, we obtain the optimal sample complexity for sufficiently large dimension. Theorem 1.4 also improves over the previous best sample complexity upper bound of [AAL24b] , which usesO k 2 d 4 + kd 2 log(1/δ) α 2 ε + kd log(1/δ) α 3 ε + k 2 d 2 α 4 ε samples.Our results provide a polynomial improvement in all parameters, but to simplify the comparison, if we ignore dependencies in the error parameter α and privacy parameters ε, δ, we improve the sample complexity from k 2 d 4 to kd 2 + k 2 d + k 1.5 d 1.75 : note that our result is quadratic in the dimension whereas [AAL24b] is quartic.When the dimension is d = 1, we can provide an improved result, which is optimal for learning mixtures of univariate Gaussians (see Theorem 1.6 for a matching lower bound). Theorem 1.5. For any α, ε, δ ∈ (0, 1), k ∈ N, there exists an inefficient (ε, δ)-DP algorithm that can learn a mixture of k arbitrary univariate Gaussians (of nonzero variance) up to accuracy α, using the following number of samples:n = O k α 2 + k log(1/δ) αε .For privately learning mixtures of univariate Gaussians, the previous best-known result for arbitrary Gaussians required O k 2 log 3/2 (1/δ)α 2 εsamples [AAL21] . Importantly, we are the first paper to show that the sample complexity can be linear in the number of components.Our work purely focuses on sample complexity, and as noted in Theorems 1.4 and 1.5, they do not have polynomial time algorithms. We note that the previous works of [AAL21, AAL24b] also do not run in polynomial time. Indeed, there is reason to believe that even non-privately, it is impossible to learn GMMs in polynomial time (in terms of the optimal sample complexity) [DKS17, BRST21, GVV22] .Finally, we prove the following lower bound for learning GMMs in any fixed dimension d. Theorem 1.6. Fix any dimension d ≥ 1 number of components k ≥ 2, any α, ε at most a sufficiently small constant c * , and δ ≤ (αε/d) O(1) . Then, any (ε, δ)-DP algorithm that can learn a mixture of k arbitrary full-dimensional Gaussians in d dimensions up to total variation distance α, with probability at least 2/3, requires at least the following number of samples:Ω kd 2 α 2 + kd 2 αε + k log(1/δ) αε .Note that for d = 1, this matches the upper bound of Theorem 1.5, thus showing that our univariate result is near-optimal in all parameters α, ε, δ. Moreover, our lower bound refutes the conjecture of [AAL21] , which conjectures that only Θ k α 2 + k αε + log(1/δ) ε samples are needed in the univariate case and Θ kd 2 α 2 + kd 2 αε + log(1/δ) ε samples are needed in the d-dimensional case. However, we note that our lower bound asymptotically differs from the conjectured bound in [AAL21] only when δ is extremely small.In the non-private setting, the sample complexity of learning unrestricted GMMs with respect to total variation distance (a.k.a. density estimation) is known to be Θ(kd 2 /α 2 ) [ABM18, ABH + 18], where the upper bound is obtained by the so-called distributional compression schemes.In the private setting, the only known sample complexity upper bound for unrestricted GMMs [AAL24b] is roughly k 2 d 4 log(1/δ)/(α 4 ε), which exhibits sub-optimal dependence on various parameters 5 . This bound is achieved by running multiple non-private list-decoders and then privately aggregating the results. For the special case of axis-aligned GMMs, an upper bound of k 2 d log(1/δ) 3/2 /(α 2 ε) is known [AAL21] . These are the only known results even for privately learning (unbounded) univariate GMMs. In other words, the best known upper bound for sample complexity of privately learning univariate GMMs has quadratic dependence on k.In the related public-private setting [BKS22, BBC + 23], it is assumed that the learner has access to some public data. In this setting, [BBC + 23] show that unrestricted GMMs can be learned with a moderate amount of public and private data.Assuming the parameters of the Gaussian components (and the condition numbers of the covariance matrices) are bounded, one can create a cover for GMMs and use private hypothesis selection [BSKW19] or the private minimum distance estimator [AAK21] to learn the GMM. On the flip side, [ASZ21] prove a lower bound on the sample complexity of learning GMMs, though their lower bound is weaker than ours and is only against pure-DP algorithms.The focus of our work is on density estimation. A related problem is learning the parameters a GMM, which has received extensive attention in the (non-private) literature (e.g., [Das99, MV10, BS10, LM21, BDJ + 22, LL22] among many other papers). To avoid identifiability issues, one has to assume that the Gaussian components are sufficiently separated and have large-enough weights. In the private setting, the early work of [NRS07] demonstrated a privatized version of [VW04] for learning GMMs with fixed (known) covariance matrices. The strong separation assumption (of Ω(k 1/4 )) between the Gaussian components in [NRS07] was later relaxed to a weaker separation assumption [CCAd + 23]. A substantially more general result for privately learning GMMs with unknown covariance matrices was established in [KSSU19] , based on a privatized version of [AM05] . Yet, this approach also requires a polynomial separation (in terms of k) between the components, as well as a bound on the spectrum of the covariance matrices. [CKM + 21] weakened the separation assumption of [KSSU19] and improved over their sample complexity. This result is based on a generic method that learns a GMM using a private learner for Gaussians and a non-private clustering method for GMMs. Finally, [AAL23] designed an efficient reduction from private learning of GMMs to its non-private counterpart, removing the boundedness assumptions on the parameters and achieving minimal separation (e.g., by reducing to [MV10] ). Nevertheless, unlike density estimation, parameter estimation for unrestricted GMMs requires exponentially many samples in terms of k [MV10] .A final important question is that of efficient algorithms for learning GMMs. Much of the work on learning GMM parameters focuses on computational efficiency (e.g,. [MV10, BS10, LM21, BDJ + 22, LL22]), as does some work on density estimation (e.g., [CDSS14, ADLS17] ). However, under some standard hardness assumptions, it is known that even non-privately learning mixtures of k d-dimensional Gaussians with respect to total variation distance cannot be done in polynomial time as a function of k and d [DKS17, BRST21, GVV22].Addendum. In a concurrent submission, [AAL24a] extended the result of [AAL24b] for learning unrestricted GMMs to the agnostic (i.e., robust) setting. In contrast, our algorithm works only in the realizable (non-robust) setting. Moreover, [AAL24a] slightly improved the sample complexity result of [AAL24b] fromO(log(1/δ)k 2 d 4 /(εα 4 )) to O(log(1/δ)k 2 d 4 /(εα 2 )).The sample complexity of our approach is still significantly better than [AAL24a] in terms of all parameters-similar to the way it improved over [AAL24b] . 5 More precisely, the upper bound isO k 2 d 4 α 2 ε + kd 2 log(1/δ) α 2 ε + kd log(1/δ) α 3 ε + k 2 d 2 α 4 εX i ∈ X, ∥X i ∥ 2 ≤ L.2. For all subsets Y ⊂ X of size n ′ , there does not exist any real number r and any unit vectorv such that |⟨X i , v⟩ -r| ≤ 1 L for all X i ∈ Y. Proof. If ∥X i ∥ 2 ≥ L, then X ⊤ M X ≥ L 2 for M the d × d identity matrix I. However, E[X ⊤ M X] = d. So, by Lemma A.7, the probability of this event is at most 2e -c•min(n 4 /d,n 2 ) ≤ 2e -c•n 2 ≤ n -Ω(n ′ ) .Next, we bound the probability that the first item holds but the second item doesn't hold. First, we can create a 1 L 2 -net of unit vectors v ′ , of size (3L 2 ) d ≤ L 3d , and a 1 L -net of real numbers r ′ from [-L, L], of size 2L 2 . Now, suppose that the first item holds, but there is some unit vector v and real number r with |⟨X i , v⟩ -r| ≤ 1 L for all X i ∈ Y, for some Y ⊂ X of size n ′ . In this case, ⟨X i , v⟩ ≤ ∥X i ∥ 2 ≤ L for all X i ∈ X, so we can assume that |r| ≤ L. Thus, if v ′ is the closest unit vector to v in the net and r ′ is the closest real number to r in the net, then|⟨X i , v ′ ⟩ -r ′ | ≤ |⟨X i , v⟩ -r| + |⟨X i , v ′ -v⟩| + |r -r ′ | ≤ 1 L + ∥X i ∥ 2 • 1 L 2 + 1 L ≤ 3 L .In other words, if the first event holds and the second does not, there exists v ′ , r ′ in the net andY ⊂ X of size n ′ , such that |⟨X i , v ′ ⟩ -r ′ | ≤ 3L for all X i ∈ Y. We can now bound this via a union bound. To perform the union bound, we first bound the probability of this event holding for some fixed v ′ , r ′ , Y.For any fixed v ′ , r ′ , Y, note that ⟨X i , v ′ ⟩ Xi∈Y is just n ′ i.i.d. copies of N (0, 1). Let's call these values z 1 , . . . , z n ′ . If there is some r ′ such that |z i -r ′ | ≤ 3 L for all i, then |z i -z 1 | ≤ 6L for all i. Since the PDF of a N (0, 1) is uniformly bounded by at most 1 2 , for any fixed z 1 , the probablity that any other z i is within 6 L of z 1 is at most 6 L , so the overall probability is at most (6/L) n ′ -1 . Now, the union bound is done over n n ′ choices of Y, at most L 3d choices of v ′ , and at most 2L 2 choices of Y. Thus, the overall probability is at most(6/L) n ′ -1 • n n ′ • L 3d • 2L 2 ≤ (6n) n ′ /L n ′ -3-3d ≥ (6n/L 0.7 ) n ′. Thus, as long as L ≥ 100n 2 , this is much smaller than n -Ω(n ′ ) .By shifting and scaling appropriately, we have the following corollary.Corollary B.5. Let n ≥ n ′ ≥ 20d, and L ≥ C 1 n 2 , for C 1 the same constant as in Proposition B.4. Suppose that some data points X = {X 1 , . . . , X n } are drawn from N (µ, Σ). Then, with probability at least 1 -n -Ω(n ′ ) , the following both hold.X i ∈ X, ∥Σ -1/2 (X i -µ)∥ 2 ≤ L.2. For all subsets Y ⊂ X of size n ′ , there does not exist any real number r and any nonzero vector v such that |⟨X i , v⟩ -r| ≤ 1 L ∥Σ 1/2 v∥ 2 for all X i ∈ Y.We can modify the robust algorithm of Theorem B.3 to have the following guarantees.Lemma B.6. Let c 0 , C 0 be as in Theorem B.3, and C 1 be as in Proposition B.4. Let L := C 1 • n 2 . Also, fix any η ≤ γ ≤ c 0 , any e -d ≤ β < 1, and and anyρ such that O(η) ≤ ρ ≤ c 0 √ d. Finally, suppose n ≥ m ≥ O d γ 2 + d 2 ρ 2 .Then, there exists an (inefficient, deterministic) algorithm A that, on a dataset Y of size m, outputs either some (μ, Σ) or ⊥, with the following properties.Y, Y ′ with d H (Y, Y ′ ) ≤ η • m, if A(Y) = (μ, Σ) ̸ =⊥ and A(Y ′ ) = (μ ′ , Σ′ ) ̸ =⊥, then (μ ′ , Σ′ ) ≈ 8C0γ,8C0ρ,8C0γ (μ, Σ).2. For any fixed µ, Σ, with probability at least1 -O(β) over Y 1 , . . . , Y m ∼ N (µ, Σ), A(Y) ̸ =⊥ and A(Y) ≈ C0γ,C0ρ,C0γ (µ, Σ).3. Fix an integer k ≥ 1, and additionally suppose that m ≥ 40d • k. Suppose that X = {X 1 , . . . , X n } are drawn i.i.d. from some mixture of k Gaussians (µ i , Σ i ). Then, with probability at least 1 -O(β) over X, for every Y ′ ⊂ X of size m and every Y withd H (Y, Y ′ ) ≤ m/2, either A(Y) =⊥, or A(Y) = (μ, Σ), where there is some i ∈ [k] such that 1 9L 4 • Σ i ≼ Σ ≼ 9L 4 • Σ i and ∥Σ -1/2 i (μ -µ i )∥ 2 ≤ 10L 3 .As in Theorem B.3, A has knowledge of η, γ, ρ, β, but not µ or Σ.Proof. The algorithm A works as follows. Start off by computing A 0 (Y) = (μ, Σ), where A 0 is the algorithm of Theorem B.3. We can run the algorithm on any arbitrary dataset, even if it didn't come from Gaussian samples, and we can assume that the algorithm A 0 always outputs some mean-covariance pair (for instance, by having it output (0, I) by default instead of ⊥). Next, A, on a dataset Y, checks if any of the following three conditions hold.(a) There existsY ′ such that d H (Y, Y ′ ) ≤ η • m, and A 0 (Y ′ ) ̸ ≈ 8C0γ,8C0ρ,8C0γ A 0 (Y). (b) There exists Y i ∈ Y such that ∥ Σ-1/2 (Y i -μ)∥ 2 > 3L.(c) There exists a subset Z ⊂ Y of size 20d, and a unit vector v and real number r, such that for allY i ∈ Z, |⟨Y i , v⟩ -r| < 1 2L • ∥ Σ1/2 v∥ 2 .If any of these conditions hold, the output isA(Y) =⊥. Otherwise, the output is A(Y) = A 0 (Y).We now verify that the algorithm satisfies the required properties. Property 1 clearly holds, be- cause if there exist datasets Y, Y ′ with d H (Y, Y ′ ) ≤ η • m, A(Y) ̸ =⊥ and A(Y ′ ) ̸ =⊥,-Ω(20d) ≤ β, ∥Σ -1/2 (Y i -µ)∥ 2 ≤ L for all Y i ∈ Y,and for all Z ⊂ Y of size 20d, there does not exist a real r and a nonzero vectorv such that |⟨Y i , v⟩ -r| ≤ 1 L • ∥Σ 1/2 v∥ 2 for all Y i ∈ Y. Now, we know that by Theorem B.3, (μ, Σ) ≈ 0.5,0.5 (µ, Σ). So, ∥ Σ-1/2 (Y i -μ)∥ 2 ≤ 2∥Σ -1/2 (Y i -μ)∥ 2 ≤ 2(∥Σ -1/2 (Y i -µ)∥ + 0.5) = 2L + 1 ≤ 3L,and if |⟨Y i , v⟩ -r| ≤ 1 2L • ∥ Σ1/2 v∥ 2 , then |⟨Y i , v⟩ -r| ≤ 1 L • ∥Σ 1/2 v∥ 2 , for all Y i ∈ Y.Thus, conditions b) and c) are also not met, so Property 2 holds.For Property 3, note that if m ≥ 40d • k, then |Y ∩ X| ≥ 20d • k, i.e., Y contains at least 20dk uncorrupted points. So, by the Pigeonhole Principle, for every possible Y, there is some index i ∈ [k] such that at least 20d points in Y ∩ X come from the i th Gaussian component. Now, let us condition on the event that Corollary B.5 holds for X, where n ′ = 20d, which happens with at least 1 -n -Ω(n ′ ) ≥ 1 -β probability. We claim that for every possible Y, and for an i∈ [k] such that 20d points Z ⊂ Y ∩ X came from the i th Gaussian component, then either A(Y) =⊥ or if A(Y) = (μ, Σ) then 1 9L 4 • Σ i ≼ Σ ≼ 9L 4 • Σ i and ∥Σ -1/2 i (μ -µ i )∥ 2 ≤ 10L 3 .First, we verify that1 9L 4 • Σ i ≼ Σ ≼ 9L 4 • Σ i .Otherwise, there exists a unit vector v such that eitherv ⊤ Σv ≥ 9L 4 • v ⊤ Σ i v or v ⊤ Σv ≤ 1 9L 4 • v ⊤ Σ i v.In the former case, by Part 1 of Corollary B.5, for allY i ∈ Z, ∥Σ -1/2 i (Y i -µ i )∥ 2 ≤ L, so |⟨v, Y i -µ i ⟩| = |⟨Σ 1/2 i v, Σ -1/2 i (Y i -µ i )⟩| ≤ L•∥Σ 1/2 i v∥ 2 = L• v ⊤ Σ i v ≤ L• v ⊤ Σv 9L 4 = 1 3L •∥ Σ1/2 v∥ 2 .Thus, if we set r = ⟨v, µ i ⟩, then |⟨Y i , v⟩ -r| < 1 2L • ∥ Σ1/2 v∥ 2 for all Y i ∈ Z, so the algorithm A would have output ⊥ due to condition c). Alternatively, if there exists a unit vector v such thatv ⊤ Σv ≤ 1 9L 4 • v ⊤ Σ i v, then by condition b), ∥ Σ-1/2 (Y i -μ)∥ 2 ≤ 3L for all Y i ∈ Z, so |⟨v, Y i -μ⟩| = |⟨ Σ1/2 v, Σ-1/2 (Y i -µ)⟩| ≤ 3L • ∥ Σ1/2 v∥ 2 ≤ 1 L • ∥Σ 1/2 i v∥ 2 .So, there exists r = ⟨v, μ⟩ such that |⟨v,Y i ⟩ -r| ≤ 1 L • v ⊤ Σ i v which is impossible by Part 2 of Corollary B.5 Next, we verify that ∥Σ -1/2 i (μ -µ i )∥ 2 ≤ 10L 3 . By Triangle inequality, for every Y i ∈ Z, ∥Σ -1/2 i (Y i -μ)∥ 2 + ∥Σ -1/2 i (Y i -µ i )∥ 2 ≥ ∥Σ -1/2 i (µ i -μ)∥ 2 . But, ∥Σ -1/2 i (Y i -µ i )∥ 2 ≤ L by Part 1 of Corollary B.5, and ∥Σ -1/2 i (Y i -μ)∥ 2 ≤ 3L 2 • ∥ Σ-1/2 (Y i -μ)∥ 2 ≤ 3L 2 • 3L ≤ 9L 3 ,because we just proved that 9L 4 Σ i ≽ Σ and assuming we do not reject due to condition b). Thus,∥Σ -1/2 i (µ i -μ)∥ 2 ≤ 9L 3 + L ≤ 10L 3 .In this section, we explain the robustness-to-privacy conversion [HKMN23] that we will utilize in proving Theorem 1.4. We will also need some relevant results about computing volume, which will be important in the robustness-to-privacy conversion.Given a pair (µ, Σ), where µ ∈ R d and Σ ∈ R d×d is positive definite, we will define Proj(µ, Σ) ∈ R d•(d+3)/2 to represent the coordinates of µ along with the upper-diagonal coordinates of Σ. For any set Ω of mean-covariance pairs (µ, Σ), define Proj(Ω) := {Proj(µ, Σ) : (µ, Σ) ∈ Ω}. Because Σ is symmetric, Proj(µ, Σ) fully encodes the information about µ and Σ. We also define vol(Ω) to be the Lebesgue measure of Proj(Ω), i.e., the d(d+3)-dimensional measure of all points Proj(µ, Σ) for (µ, Σ) ∈ Ω. Next, we will define the normalized volumeEQUATIONwhere θ = Proj(µ, Σ), and we take a Lebesgue integral over θ ∈ Proj(Ω).To motivate this choice of normalized volume, we will see that the volume is invariant under some basic transformations.Given a function f : R D → R D , for some integer D ≥ 1 we recall that the Jacobian J at some x is the matrix with J ij = ∂fi ∂xj (x). For a function that takes a symmetric matrix Σ and outputs the symmetric matrix AΣA ⊤ , we view D = d(d+1) 2 and the function Σ → AΣA ⊤ as a function f from R D → R D by taking the upper triangular part of both Σ and AΣA ⊤ . Note that f is a linear function (thus f (x) = J • x where J ∈ R D×D , and moreover, the following fact is well-known.Lemma C.1. [MH08, Theorem 11.1.5] For J as defined above, det(J) = det(A) d+1 . Now, for any fixed µ ∈ R d and positive definite Σ ∈ R d×d , consider the transformation (μ, Σ) → (Σ 1/2 ΣΣ 1/2 , Σ 1/2 μ + µ), viewed as a linear map g from Proj(μ, Σ) → Proj(Σ 1/2 μ + µ, Σ 1/2 ΣΣ 1/2 ). By setting A := Σ 1/2 , note that the map g behaves like f on the last d(d+1) 2 coordinates, and on the first d coordinates, it is simply an affine map μ → Σ 1/2 μ + µ. Therefore, the overall linear map g has determinant det(Σ1/2 ) • det(Σ 1/2 ) d+1 = (det Σ) (d+2)/2 .From this, we can infer the following.Lemma C.2. Fix any µ ∈ R d and positive definite Σ ∈ R d×d . Let h be the map (μ, Σ) → (Σ 1/2 μ + µ, Σ 1/2 ΣΣ 1/2 ). Then, for any set S of mean-covariance pairs, the normalized volume of S equals the normalized volume of h(S).Proof. Let g be the corresponding map from Proj(μ, Σ) to Proj(Σ 1/2 μ + µ, Σ 1/2 ΣΣ 1/2 ). Using a simple integration by substitution, we havevol n (h(S)) = θ=Proj(μ, Σ)∈Proj(h(S)) 1 (det Σ) (d+2)/2 d θ = θ∈g(Proj(S)) 1 (det Σ) (d+2)/2 d θ = θ=Proj(μ, Σ)∈Proj(S) 1 (det(Σ 1/2 ΣΣ 1/2 )) (d+2)/2 • (det g)d θ = θ∈Proj(S) 1 (det Σ) (d+2)/2 d θ = vol n (S).We note the following restatement of (the inefficient, approx-DP version of) the main robustness-toprivacy conversion of [HKMN23] .In the following theorem, we will think of the parameter space as lying in R D , with some normalized volume vol n (Ω) = θ∈Ω p(θ)dθ, where p ≥ 0 is some nonnegative Lebesgue-measurable function and the integration is Lebesgue integration. In our application, we will think of θ = Proj(µ, Σ), and p(θ) = (det Σ) -(d+2)/2 , to match with (1).Theorem C.3. Restatement of [HKMN23, Lemma 4.2] Let 0 < η < 0.1 and 10η ≤ η * ≤ 1 be fixed parameters. Also, fix privacy parameters ε 0 , δ 0 < 1 and confidence parameter β 0 < 1. Let S = S(θ, X) ∈ R ≥0 be a score function that takes as input a dataset X = {X 1 , . . . , X n } and a parameter θ ∈ Θ ⊂ R D . For any dataset X and any0 ≤ η ′ ≤ 1, let V η ′ (X) be the D-dimensional normalized volume of points θ ∈ Θ with score at most η ′ •n, i.e., V η ′ (X) = vol n ({θ ∈ Θ : S(θ, X) ≤ η ′ • n}).Suppose the following properties hold:• (Bounded Sensitivity) For any two adjacent datasets X, X ′ and any θ ∈ Θ, |S(θ, X) -S(θ, X ′ )| ≤ 1.• (Volume) For some universal constant C, and for any X of size n such that there exists θ with S(θ, X)≤ 0.7η * n, n ≥ C • log(V η * (X)/V 0.8η * (X))+log(1/δ0) ε0•η * .Then, there exists an (ε 0 , δ 0 )-DP algorithm M, that takes as input X and outputs either some θ ∈ Θ or ⊥, such that for any datasetX, if n ≥ C • max η ′ :η≤η ′ ≤η * log(V η ′ (X)/Vη(X))+log(1/(β0•η)) ε0•η ′ , then M(X)outputs some θ ∈ Θ of score at most 2ηn with probability 1 -β 0 .We remark that the algorithm M is allowed prior knowledge of n, D, η, η * , ε 0 , δ 0 , β 0 , as well as the domain Θ, function p(θ) that dictates the normalized volume, and score function S.We remark that the original result in [HKMN23] assumes the volumes are unnormalized, but the result immediately generalizes to normalized volumes. To see why, consider a modified domainΘ ′ ∈ R D+1, where θ ′ = (θ, z) ∈ Θ ′ if and only if θ ∈ Θ and 0 ≤ z ≤ p(θ). Also, consider a modified score function S ′ acting on Θ ′ , where S ′ ((θ, z), X) := S(θ, X) for for any (θ, z) ∈ Θ ′ .Then, note that the unnormalized volume of Θ ′ , by Fubini's theorem, is precisely (θ,z)∈Θ ′ 1dzdθ = θ∈Θ p(θ)dθ, which is precisely the normalized volume of corresponding to the new Θ ′ precisely match the normalized volumes corresponding to the old Θ. A similar calculation will give us that the unnormalized volume of points (θ, z) in Θ ′ with S ′ ((θ, z), X) ≤ t equals the normalized volume of points θ ∈ Θ with S(θ, X) ≤ t, for any t.Thus, to privately find a parameter θ ∈ Θ of low score with respect to S (assuming the normalized volume constraints hold), can create Θ ′ and S ′ , and apply the unnormalized version of Theorem C.3 to obtain some (θ, z). Also, if S ′ ((θ, z), X) ≤ 2ηn, then S(θ) = S ′ ((θ, z), X) ≤ 2ηn, and if outputting (θ, z) is (ε, δ)-DP, then so is simply outputting θ.Lemma C.4. Let M ∈ R d×d be symmetric with ∥M ∥ op ≤ 0.1, and let 1 ≤ ν ≤ 2 be a scaling parameter. Then,ν -d ≤ det(I + νM ) det(I + M ) ≤ ν d .Proof. If the eigenvalues of M are λ 1 , . . . , λ d ∈ [-0.1, 0.1], thendet(I + νM ) det(I + M ) = d i=1 1 + ν • λ i 1 + λ i .Note that for x ∈ [-0.1, 0.1] and ν ≥ 1, 1+ν•x 1+x is an increasing function. Therefore,d i=1 1 + ν • λ i 1 + λ i ≤ d i=1 1 + 0.1ν 1 + 0.1 ≤ d i=1 ν = ν d and d i=1 1 + ν • λ i 1 + λ i ≥ d i=1 1 -0.1ν 1 -0.1 ≥ d i=1 1 ν = ν -d ,where we used the fact that 1 ≤ ν ≤ 2.We note an important lemmas about the normalized volumes of certain sets.Lemma C.5. Fix some γ, τ ≤ 0.1 and some scaling parameters 1 ≤ ν 1 ≤ 2 and 1 ≤ ν 2 . Let R 1 be the set of (µ, Σ) ≈ γ,τ (0, I) and R 2 be the set of (µ, Σ) ≈ ν1•γ,ν2•τ (0, I). Then, the ratio of normalized volumesvol n (R 2 ) vol n (R 1 ) ≤ ν 2d 2 1 • ν d 2Proof. By definition of normalized volume, we havevol n (R 2 ) = ∥M ∥op≤ν1•γ,∥µ∥2≤ν2•τ 1 det(I + M ) (d+2)/2 dµdM,where we are slightly abusing notation because we are truly integrating along the upper-triangular part of M . (Overall, the integral is d(d+3) 2 -dimensional.) We have a similar expression for vol n (R 1 ).Let us consider the map sending (µ, I + M ) → (ν 2 • µ, I + ν 1 • M ). This is a linear map, and for ∥µ∥ 2 ≤ γ and symmetric ∥M ∥ op ≤ τ , this is a bijective map from R 1 to R 2 . Therefore, by an integration by substitution, we havevol n (R 2 ) = ∥M ∥op≤γ,|µ∥2≤τ 1 det(I + ν 1 • M ) (d+2)/2 • ν d(d+1)/2 1 ν d 2 dµdM ≤ ν d(d+1)/2 1 ν d 2 • (ν -d ) -(d+2)/2 • ∥M ∥op≤γ,|µ∥2≤τ 1 det(I + M ) (d+2)/2 dµdM = ν 2d 2 1 ν d 2 • vol n (R 1). Above, the first line is integration by substitution, and the second line uses Lemma C.4.Next, we note the following bound on the size of a net of matrices with small Frobenius norm. Proof. First, note that the unit sphere in R d has a 1 d 10 -net (in Euclidean distance) of size e O(d log d) . Let B 0 be this net. The net B will be the set of matrices ℓ i=1 κ i w i w ⊤ i , where every w i ∈ B 0 , every κ i is an integral multiple of 1 d 10 , and every|κ i | ≤ γ. The cardinality of B is at most |B 0 | ℓ • (d 10 ) ℓ ≤ e O(ℓ•d log d) .Now, we show that B is actually a net. For any matrix M ∈ R 3 , we can write M = d i=1 λ i v i v ⊤ i , where λ i , v i are the eigenvalues and eigenvectors, respectively, of M . Assume the eigenvalues are sorted so that |λ i | are in decreasing order. Now, sinced i=1 λ 2 i ≤ ρ 2 , which means that |λ i | ≤ ρ √ ℓfor all i > ℓ.Now, let B 1 = i>ℓ λ i v i v ⊤ i : since the v i 's are orthogonal and |λ i | ≤ ρ √ ℓ , this means that ∥B 1 ∥ op ≤ α √ ℓ . Also, note that M = B 1 + ℓ i=1 λ i v i v ⊤ i .Suppose we replace each λ i with κ i by rounding to the nearest multiple of 1/d 10 , and replace each v i with w i ∈ B 0 such that ∥v i -w i ∥ 2 ≤ 1/d 10 . Then,ℓ i=1 κ i w i w ⊤i is in B, and by Triangle inequality, we haveℓ i=1 λ i v i v ⊤ i - ℓ i=1 κ i w i w ⊤ i op ≤ ℓ • (|λ i -κ i | + 2 • λ i • ∥v i -w i ∥ 2 ) ≤ O 1 d 9 ≤ ρ √ ℓ .Thus, ∥M -ℓ i=1 κ i w i w ⊤ i ∥ op ≤ 2 • ρ √ ℓ .We also note the following lemma. We defer the proof to Appendix H. Lemma C.7. Fix any µ and positive definite Σ. Then, for anyµ 1 , µ 2 , Σ 1 , Σ 2 , we have that (µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ) if and only if (Σ 1/2 µ 1 + µ, Σ 1/2 Σ 1 Σ 1/2 ) ≈ γ,ρ,τ (Σ 1/2 µ 2 + µ, Σ 1/2 Σ 2 Σ 1/2 ).Our main lemma in this subsection is the following, which roughly bounds the normalized volume of a Frobenius norm ball "fattened" by an operator norm ball. Lemma C.8. Let c 1 ∈ (0, 0.01) be a sufficiently small universal constant. Fix any µ and positive definite Σ. Fix some parameters γ 1 , γ 2 ∈ ( 1 d , c 1 ) and ρ 2 ∈ ( 1 d ,γ 3 1 1000 • √ d).Let T 1 (µ, Σ) represent the set of {(µ 1 , Σ 1 )} such that (µ 1 , Σ 1 ) ≈ γ1,γ1 (µ, Σ). Let T 2 (µ, Σ) represent the set of {(µ 2 , Σ 2 )} such that (µ 2 , Σ 2 ) ≈ γ1,γ1 (µ 1 , Σ 1 ) for some (µ 1 , Σ 1 ) ≈ γ2,ρ2,γ2 (µ, Σ). Then,vol n (T 2 (µ, Σ)) vol n (T 1 (µ, Σ)) ≤ e O(d 5/3 log d•ρ 2/3 2 /γ 2 1 ) .Furthermore, both vol n (T 1 (µ, Σ)) and vol n (T 2 (µ, Σ)) do not change even if we change µ, Σ.Proof. For now, we also assume that µ = 0 and Σ = I.Let 1 ≤ ℓ ≤ d be decided later, and let B be the net from Lemma C.6, where we set γ = γ 2 and ρ = ρ 2 . Let B 1 be a 1 d 10 -net (in Euclidean distance) over the d-dimensional unit ball, i.e., over µ with ∥µ∥ 2 ≤ 1. Note that |B 1 | ≤ e O(d log d) , and recall that |B| ≤ e O(ℓ•d log d) . Now, for any (µ 2 , Σ 2 ) ∈ T 2 (0, I), we can write(µ 2 , Σ 2 ) ≈ γ1,γ1 (µ 1 , Σ 1 ), where ∥Σ 1 -I∥ op ≤ γ 2 , ∥Σ 1 -I∥ F ≤ ρ 2 , and ∥µ 1 ∥ 2 ≤ γ 2 . Now, we can choose some µ ′ ∈ B 1 with ∥µ 1 -µ ′ ∥ 2 ≤ d -10 and B ∈ B such that ∥(Σ 1 -I) -B∥ op ≤ 2ρ2 √ℓ . What can we say about the relationship between (µ 2 , Σ 2 ) and (µ ′ , I + B)? First, note that because(µ 2 , Σ 2 ) ≈ γ1,γ1 (µ 1 , Σ 1 ), we have (1-γ 1 )•Σ 1 ≼ Σ 2 ≼ (1+γ 1 )•Σ 1 . Next, since ∥Σ 1 -(I +B)∥ op = ∥(Σ 1 -I) -B∥ op ≤ 2ρ2 √ℓ , and since Σ 1 has all eigenvalues between 1 -γ 2 ≥ 1 2 and 1 + γ 2 ≤ 2, this means 1 -4ρ2 √ ℓ • Σ 1 ≼ (I + B) ≼ 1 + 4ρ2 √ ℓ • Σ 1 . Thus, if 4ρ < √ ℓ, we have that 1 -γ 1 1 + (4ρ 2 / √ ℓ) • (I + B) ≼ Σ 1 ≼ 1 + γ 1 1 -(4ρ 2 / √ ℓ) • (I + B).Next, because (µ 2 , Σ 2 ) ≈ γ1,γ1 (µ 1 , Σ 1 ), we have ∥Σ-1/2 1 (µ 2 -µ 1 )∥ 2 ≤ γ 1 , which means ∥µ 2 - µ 1 ∥ 2 ≤ 2γ 1 . Thus, because ∥µ 1 -µ ′ ∥ 2 ≤ d -10 ≤ γ 1 , this means that ∥µ 2 -µ ′ ∥ 2 ≤ 3γ 1 . Moreover, assuming 10ρ ≤ √ ℓ, the eigenvalues of I + B are at least (1 -γ 1 ) • 1 -4ρ2 √ ℓ ≥ 1 4 , which means that ∥(I + B) -1/2 (µ 2 -µ ′ )∥ 2 ≤ 6γ 1 . In summary, if 10ρ ≤ √ ℓ, then 1 -γ 1 1 + (4ρ 2 / √ ℓ) • (I + B) ≼ Σ 2 ≼ 1 + γ 1 1 -(4ρ 2 / √ ℓ) • (I + B), ∥(I + B) -1/2 (µ 2 -µ ′ )∥ 2 ≤ 6γ 1 .Note that if γ 1 ≤ 0.1 and 10ρ ≤ √ ℓ, then by a simple calculation,1+γ1 1-(4ρ2/ √ ℓ) ≤ 1 + γ 1 + 10ρ2 √ ℓ and 1-γ1 1+(4ρ2/ √ ℓ) ≥ 1 -γ 1 -4ρ2 √ ℓ . This implies that (µ 2 , Σ 2 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (µ ′ , I + B), for some µ ′ ∈ B 1 and B ∈ B. Note that this holds for any (µ 2 , Σ 2 ) ∈ T 2 (0, I). (µ 2 , Σ 2 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (µ ′ , Σ ′ ). Now, by Lemma C.7, (µ 2 , Σ 2 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (µ ′ , Σ ′ ) if and only if µ 2 = Σ ′ 1/2 µ 0 + µ ′ and Σ 2 = Σ ′ 1/2 Σ 0 Σ ′ 1/2 , where (µ 0 , Σ 0 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (0, I). By Lemma C.2, this means that the volume of {(µ 2 , Σ 2 ) : (µ 2 , Σ 2 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (µ ′ , Σ ′ )} equals the volume of {(µ 0 , Σ 0 ) : (µ 0 , Σ 0 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (0, I)}. Thus, if 10ρ ≤ √ℓ, then T 2 (0, I) can be covered by e O(ℓ•d log d) regions, each of which has the same normalized volume as the volume of {(µ 0 , Σ 0 ) :(µ 0 , Σ 0 ) ≈ γ1+10ρ2/ √ ℓ,6γ1 (0, I)}. Moreover, if we further have 10ρ2 √ ℓ ≤ γ 1 , then by Lemma C.8, this is at most 1 + 10ρ2 √ ℓ•γ1 2d 2 • 6 d • vol n (T 1 (0, I)).We will set ℓ = 100(dρ2) 2/3γ 2 1 . Since ρ 2 ≥ 1 d , ℓ ≥ 100 as long as γ 1 ≤ c 1 ≤ 1. Also, ℓ ≤ 100 • γ 3 1 1000 • d 3/2 2/3 • 1 γ 2 1 = d. Finally, 10ρ2 √ ℓ = γ1•ρ 2/3 2 d 1/3 ≤ γ 1 , since ρ 2 ≤ √ d.Overall, the ratiovol n (T 2 (0, I)) vol n (T 1 (0, I)) ≤ e O(ℓ•d log d) • 1 + 10ρ 2 √ ℓ • γ 1 2d 2 • 6 d ≤ e O(d 5/3 log d•ρ 2/3 2 /γ 2 1 ).Finally, we show how to remove the assumption that µ = 0 and Σ = I. First, note that for any general µ, Σ, (µ 1 , Σ 1 ) ∈ T 1 (0, I) if and only if(Σ 1/2 µ 1 + µ, Σ 1/2 Σ 1 Σ 1/2 ) ∈ T 1 (µ, Σ), by Lemma C.7. So, by Lemma C.2, vol n (T 1 (µ, Σ)) = vol n (T 1 (0, I)). Next, (µ 2 , Σ 2 ) ≈ γ1,γ1 (µ 1 , Σ 1 ) and (µ 1 , Σ 1 ) ≈ γ2,ρ2,γ2 (0, I) if and only if (Σ 1/2 µ 2 + µ, Σ 1/2 Σ 2 Σ 1/2 ) ≈ γ1,γ1 (Σ 1/2 µ 1 + µ, Σ 1/2 Σ 1 Σ 1/2 ) and (Σ 1/2 µ 1 + µ, Σ 1/2 Σ 1 Σ 1/2 ) ≈ γ2,ρ2,γ2 (µ, Σ), by two applications of Lemma C.7. So, by Lemma C.2, vol n (T 2 (µ, Σ)) = vol n (T 2 (0, I)). Thus, the volume ratios stay the same as well.In this section, we prove that if we know a very crude approximation to all of the Gaussian mixture components with sufficiently large weight, we can privately learn the full density of the Gaussian mixture model up to low total variation distance. This will be useful for proving both Theorem 1.4 and Theorem 1.5.We start with the following auxiliary lemma. Lemma D.1. Let G ≥ 1 and ζ ≤ 1 be some parameters. For some d ≥ 1, let μ ∈ R d , and let Σ ∈ R d×d be positive definite. Let U μ, Σ be the set of (µ, Σ)such that 1 G • Σ ≼ Σ ≼ G • Σ and ∥ Σ-1/2 (µ -μ)∥ 2 ≤ G. Then, there exists a net B μ, Σ of size O(G √ d/ζ) d(d+3) such that for every (µ, Σ) ∈ U μ, Σ, there exists (µ ′ , Σ ′ ) ∈ B μ, Σ such that d TV (N (µ, Σ), N (µ ′ , Σ ′ )) ≤ ζ.Proof. First, assume that Σ = I and μ = 0. Now, let us consider the set U := U 0,I = {(µ, Σ) :∥µ∥ 2 ≤ G, 1 G • I ≼ Σ ≼ G • I}.Let G ′ ≥ G be a parameter that we will set later. Now, we can look at the un-normalized volume of U viewed in R d(d+3)/2 (i.e., we are projecting to the upper-triangular part of Σ). First, we claim (using a standard volume argument) that there is a cover B ⊂ U of size (2G ′ ) d(d+3) , such that for every (µ, Σ) ∈ U, there is(µ ′ , Σ ′ ) ∈ B with ∥µ -µ ′ ∥ 2 ≤ 1 G ′ and ∥Σ -Σ ′ ∥ op ≤ 1 G ′ . To see why, consider a maximal packing of (µ ′ i , Σ ′ i ) ∈ U such that for every i ̸ = j, either ∥µ ′ i -µ ′ j ∥ 2 > 1 G ′ or ∥Σ ′ i -Σ ′ j ∥ op > 1 G ′ .Then, this set is clearly a cover of U (or else we could have increased the size of the packing, breaking maximality). Then, the setsS i := {(µ, Σ) : ∥µ -µ ′ i ∥ ≤ 1 2G ′ , ∥Σ -Σ ′ i ∥ op ≤ 1 2G ′ } are disjoint, and are all contained in the set S := {(µ, Σ) : ∥µ∥ 2 ≤ 2G ′ , ∥Σ∥ op ≤ 2G ′ }. S is just a shifting and scaling of S i by a factor of 4(G ′ ) 2 , and thus vol(S) = (4(G ′ ) 2 ) d(d+3)/2 •vol(S i ) = (2G ′ ) d(d+3) • vol(S i ).Because every S i is disjoint and contained in S, the number of such indices i is at most(2G ′ ) d(d+3) .Next, consider any (µ, Σ) ∈ U and(µ ′ , Σ ′ ) ∈ B with ∥µ -µ ′ ∥ 2 ≤ 1 G ′ , ∥Σ -Σ ′ ∥ op ≤ 1 G ′ . Then, ∥Σ -1/2 Σ ′ Σ -1/2 -I∥ op = ∥Σ -1/2 (Σ ′ -Σ)Σ -1/2 ∥ op ≤ ∥Σ -1 ∥ op • ∥Σ ′ -Σ∥ op ≤ G G ′ . Also, ∥Σ -1/2 (µ -µ ′ )∥ 2 ≤ ∥Σ -1/2 ∥ op • ∥µ -µ ′ ∥ 2 ≤ √ G G ′ .In other words, we have found a net B of size (2G ′ ) d(d+3) such that for every (µ, Σ) ∈ U, there exists(µ ′ , Σ ′ ) ∈ B with (µ ′ , Σ ′ ) ≈ G/G ′ ,G/G ′ (µ, Σ).Next, consider general μ, Σ. Note that (µ, Σ) ∈ U μ, Σ is equivalent to µ = Σ1/2 µ 0 + μ, where ∥µ 0 ∥ 2 ≤ G, and Σ = Σ1/2 Σ 0 Σ1/2 , where 1G • I ≼ Σ 0 ≼ G • I. So, if we consider the map f : (µ 0 , Σ 0 ) → ( Σ1/2 µ 0 + μ, Σ1/2 Σ 0 Σ1/2 ), this bijectively maps U to U μ, Σ.We can also consider the cover B μ, Σ = f (B), where B is the cover constructed for U. Then, by Lemma C.7,(µ ′ 0 , Σ ′ 0 ) ≈ G/G ′ ,G/G ′ (µ 0 , Σ 0 ) if and only if f (µ ′ 0 , Σ ′ 0 ) ≈ G/G ′ ,G/G ′ f (µ 0 , Σ 0 ).Hence, regardless of the choice of μ, Σ, for every (µ, Σ) ∈ U μ, Σ, there exists d+3) and for every (µ, Σ) ∈ U μ, Σ, there exists(µ ′ , Σ ′ ) ∈ B μ, Σ with (µ ′ , Σ ′ ) ≈ G/G ′ ,G/G ′ (µ, Σ). Moreover, |B μ, Σ| = |B| ≤ (2G ′ ) d(d+3) . Moreover, if (µ ′ , Σ ′ ) ≈ G/G ′ ,G/G ′ (µ, Σ), then ∥Σ -1/2 (µ ′ -µ)∥ 2 ≤ G G ′ and ∥Σ -1/2 ΣΣ -1/2 -I∥ F ≤ √ d • ∥Σ -1/2 ΣΣ -1/2 -I∥ op ≤ √ d•G G ′ . So, by Lemma A.8, d TV (N (µ, Σ), N (µ ′ , Σ ′ )) ≤ √ d•G G ′ . Hence, if we set G ′ = G √ d ζ , then the size of the net B μ, Σ is O(G √ d/ζ) d((µ ′ , Σ ′ ) ∈ B μ, Σ such that d TV (N (µ, Σ), N (µ ′ , Σ ′ )) ≤ ζ.Next, we note the following result about differentially private hypothesis selection.Let H = {H 1 , . . . , H M } be a set of probability distributions over some domain D. There exists an ε-differentially private algorithm (with respect to a dataset X = {X 1 , . . . , X n }) which has following guarantees.Let D be an unknown probability distribution over D, and suppose there exists a distributionH * ∈ H such that d TV (D, H * ) ≤ α. If n ≥ O log M α 2 + log M αε, and if X 1 , . . . , X n are samples drawn independently from D, then the algorithm will output a distribution Ĥ ∈ H such that d TV (D, Ĥ) ≤ 4α with probability at least 9/10.Given Lemma D.1 and Theorem D.2, we are in a position to convert any crude approximation into a fine approximation. Namely, we prove the following result.Lemma D.3. Let D represent an unknown GMM with representation {(w i , µ i , Σ i )} k i=1 . Let G ≥ 1 be some fixed parameter. Then there exists an ε-differentially algorithm, onn ≥ O d 2 k•log(G•k•k ′ • √ d/α) α 2 + d 2 k•log(G•k•k ′ • √ d/α) αεsamples, with the following property.Suppose the algorithm is given as input a set {(μ j , Σj )} k ′ j=1 for some k ′ ≥ 1, such that for every (µ i , Σ i ) with weight w i ≥ α k , there exists j ≤ k ′ such that 1 G • Σj ≼ Σ i ≼ G • Σj and ∥ Σ-1/2 j (µ i - μj )∥ 2 ≤ G.M = |H| is at most | B| O(k) • (1/ζ) O(k) = (k ′ • G • √ d/ζ) O(d 2 •k) . Since we set ζ = α/k, M ≤ (G • kk ′ √ d α ) O(d 2 •k) . Next, let us consider the true GMM D, with representation {(w i , µ i , Σ i )} k i=1 . Because every (µ i , Σ i ) with w i ≥ α k satisfies 1 G • Σj ≼ Σ i ≼ G • Σj and ∥ Σ-1/2 j (µ i -μj )∥ 2 ≤ G, by Lemma D.1, there exists (µ ′ i , Σ ′ i ) ∈ B such that d TV (N (µ i , Σ i ), (µ ′ i , Σ ′ i )) ≤ ζ. Moreover, we can round each w i ≥ α k to some w ′ i which is an integral multiple of ζ, such that |w i -w ′ i | ≤ ζ.Finally, for each w i < α k , we can choose an arbitrary Gaussian in B and round w i to some w ′ i . Then, the total variation distance between D and the GMM with representation{(w ′ i , µ ′ i , Σ ′ i )} k i=1 is at most i≤k:wi>α/k α k + |w i -w ′ i | + i≤k:wi≥α/k (ζ + |w i -w ′ i |) ≤ α + O(k • ζ).Hence, if we set ζ = α k , there exists a distributionH * ∈ H with d TV (D, H * ) ≤ O(α).Therefore, the algorithm of Theorem D.2, using n≥ O d 2 k•log(G•k•k ′ • √ d/α) α 2 + d 2 k•log(G•k•k ′ • √ d/α) αε samples, will find H ∈ H such that d TV (D, H) ≤ O(α), and is ε-differentially private.Pseudocode: We give a simple pseudocode for the algorithm of Lemma D.3, in Algorithm 1.Algorithm 1: FINE-ESTIMATE(X 1 , X 2 , . . . , X n ∈ R d , d, k, k ′ , ε, α, G, {(μ j , Σj )} j≤k ′ )Input: Samples X 1 , . . . , X n , and crude predictions (μ j , Σj ), where 1 ≤ j ≤ k ′ for some k ′ . Output: An (ε, 0)-DP prediction H, which is a mixture over at most k Gaussians. 6 Let H be the set of mixtures {(w ′ i , µ ′ i , Σ ′ i )} of k or fewer components, where every (µ ′ i , Σ ′ i ) ∈ B and every w ′ i is an integral multiple of ζ. 7 Run the ε-DP algorithm of Theorem D.2 on X 1 , . . . , X n with respect to H.In this section, we provide the algorithm and analysis for Theorem 1.4. We first describe and provide pseudocode for the algorithm, and then we prove that it can privately learn mixtures of d-dimensional Gaussians with low sample complexity.High-Level Approach: Suppose that the unknown distribution is a GMM with representation {(w i , µ i , Σ i )} k i=1 . Define Θ to be the set of all feasible mean-covariance pairs (µ, Σ), i.e., where µ ∈ R d and Σ ∈ R d×d is positive definite. We also start off with a set Ω = Θ, which will roughly characterize the region of "remaining" mean-covariance pairs. At a high level, our algorithm proceeds as follows. We will first learn a very crude approximation of each (µ i , Σ i ), one at a time. Namely, using roughlyε √ k log(1/δ) , δk -DP, we will learn some (μ 1 , Σ1 ) ∈ Ω which is "vaguely close" to some (µ i , Σ i ). We will then remove every (µ, Σ) that is vaguely close to this (μ, Σ) from our Ω, and then attempt to repeat this process. We will repeat it up to k times, in hopes that every (µ i , Σ i ) is vaguely close to some (μ j , Σj ). By advanced composition, the full set of {(μ j , Σj )} is still (ε, δ)-DP.At this point, the remainder of the algorithm is quite simple. Namely, we have some crude estimate for every (µ i , Σ i ) (it is "vaguely close" to some (μ j , Σj )). Moreover, one can create a fine net of roughly e Õ(d 2 ) (µ, Σ)-pairs that cover all mean-covariance pairs vaguely close to each (μ j , Σj ), since the dimension of (µ, Σ) is O(d 2 ). Moreover, the weight vector (w 1 , . . . , w k ) has a fine net of size roughly e Õ(k) . As a result, we can reduce the problem to a small set of hypotheses: there are roughly e Õ(d 2 ) choices for each (µ i , Σ i ), and e Õ(k) choices for w, for a total of e Õ(k•d 2 ) choices for the mixture of Gaussians. We can then apply known results on private hypothesis selection [BSKW19] , which will suffice.The main difficulty in the algorithm is in privately learning some (μ, Σ) which is "vaguely close" to some (µ i , Σ i ). We accomplish this task by applying the robustness-to-privacy conversion of [HKMN23] , along with a carefully constructed score function, which we will describe later in this section.Algorithm Description: Let c 0 ≤ 0.01 and C 0 > 1 be the constants of Lemma B.6, and c 1 ≤ 0.01 be the constant of Lemma C.8. Let c 2 ≤ 0.1 • min(c 0 , c 1 ) be a sufficiently small constant. We set parametersη * = c2•α 8k , η = η * 10 , ε 0 = ε √ 4k log(1/δ) , δ 0 = δ 2k , β 0 = e -d . We also set m = O(d 1.75 ), N = O m • √ k log(1/δ) ε + √ k•log 3/2 (1/δ) ε+ kd , and n = N • 2k α . Finally, we define parameters γ = τ = c 2 and ρ = c 2 • d 1/8 . (See Lines 1-4 of Algorithm 2 for a more precise description of some of the parameters.)We now define the main score function. To do so, we first set some auxiliary parameter η ′ = γ ′ = c2 8C0 and ρ ′ = ρ 8C0 . We will consider the (deterministic) algorithm A of Lemma B.6, where A is given parameters η ′ , γ ′ , ρ ′ , β 0 , that acts on a dataset Z and outputs either ⊥ or some (μ, Σ). Next, for some domain Ω ⊂ Θ of "feasible" mean-covariance pairs (µ, Σ), we will define a function f Ω , which takes as input a dataset Y of size N and a mean-covariance pair (µ, Σ) ∈ Θ, and outputs (2) Finally, given a dataset X of size n and ( µ, Σ) ∈ Θ, we define the score functionf Ω (Y, µ, Σ) =    1 if (µ, Σ) ∈ Ω , A(Y) ≈ γ,S Ω (( µ, Σ), X) = min t : ∃X ′ , Y ′ , µ, Σ such that d H (X, X ′ ) = t, Y ′ ⊂ X ′ , |Y ′ | = N, ( µ, Σ) ≈ γ,τ (µ, Σ), f Ω (Y ′ , µ, Σ) = 1 . (3)Note that ( µ, Σ) does not need to be in Ω, but it must satisfy ( µ, Σ) ≈ γ,τ (µ, Σ) for some (µ, Σ) ∈ Ω. Also, note that it is possible that no matter how one changes the data points, the conditions are never met (for instance if ( µ, Σ) ̸ ≈ γ,τ (µ, Σ) for any (µ, Σ) ∈ Ω). This is not a problem: we will simply set the score to be +∞ if this occurs.With this definition of score function S Ω , we can let M be an (ε 0 , δ 0 )-DP algorithm based on Theorem C.3, with the settings of n, η, η * , ε 0 , δ 0 , β 0 as above, and with D = n(n+3). Here, we recall that (μ, Σ) is viewed as n(n+3) 2 -dimensional by only considering the upper-triangular part of Σ. We also assume the domain is Θ and the normalized volume is as in (1).Given this algorithm M (which implicitly depends on Ω), the algorithm works as follows. We will actually draw a total of n + n ′ samples (for n as above, and n ′ to be defined later), though we start by just looking at the first n samples X = {X 1 , . . . , X n }. We initially set Ω = Θ, and run the following procedure up to k times, or until the algorithm outputs ⊥. For the j th iteration, we use M(X) to compute some prediction (μ j , Σj ) (or possibly we output ⊥). Assuming we do not output ⊥, we remove from Ω every (µ, Σ) such that n-12 • Σj ≼ Σ ≼ n 12 • Σj and ∥ Σ-1/2 j (µ -μj )∥ 2 ≤ n 12 .At the end, we have computed some {(μ j , Σj )} k ′ j=1 where 0 ≤ k ′ ≤ k. If k ′ = 0 we simply output ⊥. Otherwise, we will draw n ′ = Õ d 2 k α 2 + d 2 k αε fresh samples. We run the ε-DP hypothesis selection based algorithm, Algorithm 1 on the fresh samples, using parameters G = n 12 and {(μ j , Σj )} k ′ j=1 .Pseudocode: We give pseudocode for the algorithm in Algorithm 2.Algorithm 2: ESTIMATE(X 1 , X 2 , . . . , X n+n ′ ∈ R d , k, ε, δ, α) Input: Samples X 1 , . . . , X n , X n+1 , . . . , X n+n ′ . Output: An (ε, δ)-DP prediction {( µ i , Σ i ), w i }./* Set parameters */ 1 Let c 2 be a sufficiently small constant and K = poly log(d, k, 1 α , 1 ε , log 1 δ ) be sufficiently large.2 η * ← c2•α 8k , η ← η * 10 , ε 0 ← ε √ 4k log(1/δ) , δ 0 ← δ 2k , β 0 ← e -d . 3 m = K • d 1.75 , N ← K • m+log(1/δ0) ε0 + kd , n ← N • 2k α . 4 γ, τ ← c 2 ρ ← 10 • d 1/8 . 5 Θ, Ω ← {(µ, Σ) : µ ∈ R d , Σ ∈ R d×d Symmetric, Positive Definite}. /* Get X */ 6 Obtain samples X ← {X 1 , X 2 , . . . , X n }./* Learn a crude approximation (μ j , Σj ) of the mean-covariance pairs, one at a time */ 7 for j = 1 to k do 8 Define f Ω (Y, µ, Σ) and S Ω ((µ, Σ), X) as in (2) and (3), respectively. 9 Let M be the (ε 0 , δ 0 )-DP algorithm obtained by Theorem C.3 with the score function S Ω , with parameters n, η, η * , ε 0 , δ 0 , β 0 as defined above, D = n(n+3), domain Θ, and p(µ, Σ) = (det Σ) -(d+2)/2 . 10 A ← M(X) 11 if A ̸ =⊥ then 12 μj , Σj ← A 13 Ω ← Ω\{(µ, Σ) : n -12 • Σj ≼ Σ ≼ n 12 • Σj and ∥ Σ-1/2 j (µ -μj )∥ 2 ≤Set n ′ ← K • d 2 k α 2 + d 2 k αε . 19 Sample X n+1 , . . . , X n+n ′ , and redefine X ← {X n+1 , . . . , X n+n ′ }. 20 Run Algorithm 1 on X, {(μ j , Σj )}, with k ′ = #{(μ j , Σj )} and G = n 12 .In this section, we analyze Lines 7-17 of Algorithm 2. The main goal is to show that the algorithm is private, and that for samples drawn from a mixture of Gaussians, every component (µ i , Σ i ) with large enough weight w i is "vaguely close" to some (μ j , Σj ) computed by the algorithm.First, we show that when the samples are actually drawn as a Gaussian mixture, then under some reasonable conditions, any (µ, Σ) close to a true mean-covariance pair has low score. Proposition E.1. Suppose X = {X 1 , . . . , X n } is drawn from a Gaussian mixture model, with representation {(w i , µ i , Σ i )} k i=1 . Then, with 1 -O(k • β 0 ) probability over X, for any set Ω of mean-covariance pairs, for all i ∈ [k] such that (µ i , Σ i ) ∈ Ω and w i ≥ α/k, and for all ( µ, Σ) ≈ γ,τ (µ i , Σ i ), S Ω (( µ, Σ), X) = 0.Proof. Fix any i ∈ [k] with w i ≥ α/k. Set µ = µ i and Σ = Σ i . Also, let X ′ = X, and Y be a random subset of size N among the samples in X actually drawn from N (µ i , Σ i ). Note that the number of such samples has distribution Bin(n, w i ), which for w i ≥ α k and n = 2k α • N , is at least N with e -Ω(N ) ≤ β 0 failure probability.Then, Y is just N i.i.d. samples from N (µ i , Σ i ). So, if we draw a random subset Z of size m of Y, it has the same distribution as m i.i.d. samples from N i , Σ i ). We apply Part 2 of Lemma B.6 (wherewe use parameters η ′ , γ ′ , ρ ′ in the application). Note that m ≥ O d γ ′ 2 + d 2 ρ ′ 2 , so with probability at least 1 -O(β 0 ) over Z, A(Z) ≈ 8C0γ ′ ,8C0ρ ′ ,8C0γ ′ (µ i , Σ i ). By our setting of γ ′ , ρ ′ , this means that A(Z) ≈ γ,ρ,γ (µ i , Σ i ).In other words, for each index i ∈ N m and corresponding subset Z of Y, if we let W i be the indicator that A(Z) ̸ ≈ γ,ρ,γ (µ i , Σ i ), then P(W i = 1) ≤ O(β 0 ). While the values of W i are not necessarily independent, by linearity of expectation we have thatE[ W i ] ≤ N m • O(β 0 ), so the probability that E[ W i ] ≥ 1 3 • N m is at most O(β 0 ) by Markov's inequality. Moreover, because N ≥ m ≥ O d γ ′2 + d 2 ρ ′2, we can apply A on Y and we again obtain that with probability1-O(β 0 ), A(Y) ≈ γ,ρ,γ (µ i , Σ i ).i ∈ [k] with w i ≥ α k , with probability at least 1 -O(β 0 ) over X, there exists Y ⊂ X of size N , such that A(Y) ≈ γ,ρ,γ (µ i , Σ i ) and at least 2/3 of the subsets Z ⊂ Y of size m have A(Z) ≈ γ,ρ,γ (µ i , Σ i ). Thus, for any set Ω, if (µ i , Σ i ) ∈ Ω then S Ω (( µ, Σ), X) = 0 for all ( µ, Σ) ≈ γ,τ (µ i , Σ i ). The proof follows by a union bound over all i ∈ [k].Next, we show that for any dataset X, if there is even a single (µ * , Σ * ) with low score, there must be a region of ( µ, Σ) which all has low score. Proposition E.2. Suppose that t = min µ * ,Σ * S Ω ((µ * , Σ * ), X). Then, there exists some µ, Σ such that for all ( µ, Σ) ≈ γ,τ (µ, Σ), we have that S Ω (( µ, Σ), X) = t.Proof. Fix µ * , Σ * so that S Ω ((µ * , Σ * ), X) = t. Then, there is some (µ, Σ) and some X ′ , Y such that d H (X, X ′ ) = t, Y ⊂ X ′ , |Y| = N , and f Ω (Y, µ, Σ) = 1. Thus, for any ( µ, Σ) ≈ γ,τ (µ, Σ), by definition, we have that S Ω (( µ, Σ), X) ≤ t. But since t is the minimum possible score, we in fact have equality: S Ω (( µ, Σ), X) = t for all such ( µ, Σ).Next, we want to show that regardless of the dataset (even if not drawn from any Gaussian mixture distribution), the set of points of low score can't have that much volume.Proposition E.3. Fix a dataset X of size n. Then, the set of µ, Σ such that S Ω (( µ, Σ), X) ≤ η * • n can be partitioned into at most n m regions S i , which is indexed by some (µ ′ i , Σ ′ i ). Moreover, for all ( µ, Σ) ∈ S i , there exists (µ, Σ) such that ( µ, Σ) ≈ γ,τ (µ, Σ) and (µ, Σ) ≈ 8γ,8ρ,8τ (µ ′ i , Σ ′ i ). Proof. Pick any ( µ, Σ) with score at most η * • n. Let X ′ , Y ′ , µ, Σ be such that d H (X, X ′ ) ≤ η * • n, Y ′ ⊂ X ′ , |Y ′ | = N , ( µ, Σ) ≈ γ,τ (µ, Σ), and f Ω (Y ′ , µ, Σ) = 1. If we define Y ⊂ X of size N to be the corresponding subset as Y ′ is to X ′ , then d H (Y, Y ′ ) ≤ η * • n = c2 4 • N . Now, ifwe take a random subset Z ′ of size m in Y ′ , and look at the corresponding subset Z of size m in Y, by a Chernoff bound, with at least 0.99 probability d H (Z, Z ′ ) ≤ c2 2 • m. Moreover, with at least 2/3 probability over the random subset Z ′ , A(Z ′ ) ≈ γ,ρ,τ (µ, Σ). Therefore, there always exists a subset Z ⊂ X of size m and a setZ ′ of size m such that d H (Z, Z ′ ) ≤ c2 2 • m and A(Z ′ ) ≈ γ,ρ,τ (µ, Σ).Now, for any fixed Z ⊂ X of size m, if we look at any sets Z ′ , Z ′′ of size m and Hamming distance at most c2 2 • m from Z, then d H (Z ′ , Z ′′ ) ≤ c 2 • m. So, by Property 1 of Lemma B.6, for every such Z ′ and Z ′′ with A(Z ′ ), A(Z ′′ ) ̸ =⊥, we must have A(Z ′′ ) ≈ γ,ρ,τ A(Z ′ ).To complete the proof, we order the subsets Z 1 , Z 2 , . . . , Z ( n m ) of size m in X. For each i ≤ n m , if there exists someZ ′ such that d H (Z i , Z ′ ) ≤ c2 2 • m and A(Z ′ ) ̸ =⊥, choose an arbitrary such Z ′ and let (µ ′ i , Σ ′ i ) := A(Z ′ ). Otherwise, we do not define (µ ′ i , Σ ′ i ).Then, for any ( µ, Σ) of score at most η * • n, there exists (µ, Σ), a subset Z i ⊂ X of size m, and a setZ ′ i of size m such that ( µ, Σ) ≈ γ,τ (µ, Σ), d H (Z i , Z ′ i ) ≤ c2 2 • m, and A(Z ′ i ) ≈ γ,ρ,τ (µ, Σ). Because A(Z ′ i ) ̸ =⊥ and Z ′ i has Hamming distance at most c2 2 m from Z i , this also implies that A(Z ′ i ) ≈ γ,ρ,τ (µ ′ i , Σ ′ i ). By Proposition B.2, this means that (µ, Σ) ≈ 8γ,8ρ,8τ (µ ′ i , Σ ′ i ), which completes the proof.Proposition E.4. Suppose X = {X 1 , . . . , X n } is drawn from a Gaussian mixture model, with representation {(w i , µ i , Σ i )} k i=1 . Then, with probability at least 1 -O(β 0 ) over the randomness of X, for any set Ω of mean-covariance pairs, and for every ( µ, Σ) with S Ω (( µ, Σ), X) ≤ N 2 , there exists an index i ∈[k] such that 1 O(n 8 ) • Σ i ≼ Σ ≼ O(n 8 ) • Σ i and ∥Σ -1/2 i ( µ -µ i )∥ 2 ≤ O(n 6 ).Proof. We apply Part 3 of Lemma B.6, though we use N to replace the value m in Lemma B.6. We are assuming N ≥ 40d • k. So, by definition of score, ( µ, Σ) ≈ γ,τ (µ, Σ) where f Ω (Y ′ , µ, Σ) = 1, for some Y ′ which can be generated by taking a subset of X of size at most N and altering at most N/2 elements. Sincef Ω (Y ′ , µ, Σ) = 1, this means that if (μ, Σ) = A(Y ′ ) then (μ, Σ) ≈ γ,ρ,τ (µ, Σ). So, by Proposition B.2, ( µ, Σ) ≈ 8γ,8τ (μ, Σ). Assuming γ = τ = c 2 is sufficiently small, this means that ( µ, Σ) ≈ 1/2,1/2 (μ, Σ).By Part 3 of Lemma B.6, with probability at least 1-O(β 0 ), there exists i∈ [k] such that 1 O(n 8 ) •Σ i ≼ Σ ≼ O(n 8 ) • Σ i and ∥Σ -1/2 i (μ -µ i )∥ 2 ≤ O(n 6 ). However, we know that 1 2 Σ ≼ Σ ≼ 2 Σ, which means that 1 O(n 8 ) •Σ i ≼ Σ ≼ O(n 8 )•Σ i . Moreover, ∥ Σ-1/2 (μ-µ)∥ 2 ≤ 0.5, so ∥Σ -1/2 i (μ-µ)∥ 2 ≤ O(n 4). Thus, by triangle inequality, we have that ∥Σ-1/2 i ( µ -µ i )∥ 2 ≤ O(n 6 ).Now, given a set X = {X 1 , . . . , X n } drawn from a GMM with representation {(w i , µ i , Σ i )} k i=1 , we say X is regular if it satisfies Propositions E.1 and E.4. In other words: 1. for any set Ω of mean-covariance pairs, for all i ∈ [k] such that (µ i , Σ i ) ∈ Ω and w i ≥ α/k, and for all ( µ, Σ) ≈ γ,τ (µ i , Σ i ), S Ω (( µ, Σ), X) = 0.2. for any set Ω of mean-covariance pairs, and for every ( µ, Σ) with S Ω (( µ, Σ), X) ≤ N 2 , there exists an index i ∈[k] such that 1 O(n 8 ) • Σ i ≼ Σ ≼ O(n 8 ) • Σ i and ∥Σ -1/2 i ( µ -µ i )∥ 2 ≤ O(n 6 ).When we say X is regular, the components (µ i , Σ i ) and weights w i are implicit.We now show that every step of the crude approximation (Lines 8-16 in Algorithm 2) will, with 1-O(β 0 ) probability, find some crude approximation (μ j , Σj ) to some (µ i , Σ i ), as long as(µ i , Σ i ) ∈ Ω.Lemma E.5. For any Ω ⊂ Θ, the corresponding algorithm M (which depends on Ω) is (ε 0 , δ 0 )-DP on datasets X.Suppose additionally that X is regular, and that there exists i ∈ [k] such that w i ≥ α k and (µ i , Σ i ) ∈ Ω. Then, with 1 -O(β 0 ) probability (just over the mechanism M), M(X) = ( µ, Σ) and• There exists i ∈ [k] such that 1 O(n 8 ) • Σ i ≼ Σ ≼ O(n 8 ) • Σ i and ∥Σ -1/2 i ( µ -µ i )∥ 2 ≤ O(n 6 ).• There exists (µ, Σ) ∈ Ω such that ( µ, Σ) ≈ γ,τ (µ, Σ).Proof. We will apply Theorem C.3. We first need to check the necessary conditions. It follows almost immediately from the definition of S Ω that S Ω (( µ, Σ), •) has the bounded sensitivity property with respect to neighboring datasets, for any fixed Ω, µ, Σ.To check the volume condition, note that for any dataset X of size n, if min µ,Σ S Ω ((µ, Σ), X) ≤ 0.7η * n, then by Proposition E.2, there exists µ, Σ such that S Ω (( µ, Σ), X) ≤ 0.8η * n for all ( µ, Σ) ≈ γ,τ (µ, Σ). Conversely, for any X, by Proposition E.3, the set of ( µ, Σ) of score at most η * n can be partitioned into n m regions S i , indexed by (µ ′ i , Σ ′ i ). Moreover, each S i is contained in the set of ( µ, Σ) such that there exists (µ, Σ) such that ( µ, Σ) ≈ γ,τ (µ, Σ) ≈ 8γ,8ρ,8τ (µ ′ i , Σ ′ i ). Let us use the notation of Lemma C.8, where γ 1 = γ, γ 2 = 8γ, and ρ 2 = 8ρ. Then, we have that the set of points with score at most 0.8η * n contains T 1 (µ, Σ) some (µ, Σ), but the set of points with score at most η * n is contained in the union of T 2 (µ ′ i , Σ ′ i ) for at most n m choices of (µ ′ i , Σ ′ i ). So, as long as min µ,Σ S Ω ((µ, Σ), X) ≤ 0.7η * n, we can apply Lemma C.8 to obtainV η * (X) V 0.8η * (X) ≤ n m • exp O d 5/3 • log d • ρ 2/3 2 γ 2 1 = exp O d 7/4 log d + m log n , since ρ 2 = 8ρ = 8c 2 • d 1/8and γ 1 = γ = c 2 , and since c 2 is a constant. Thus, ifn ≥ C • O(d 7/4 log d + m log n) + log(1/δ 0 ) ε 0 • η * ,we have that M is (ε 0 , δ 0 )-DP, by Theorem C.3. This holds by our parameter settings, assuming K is sufficiently large.Next, we prove accuracy. Assume that X is regular. We again adopt the notation of Lemma C.8, (with γ 1 = γ, γ 2 = 8γ, and ρ 2 = 8ρ). By the first condition of regularity, for all i ∈ [k] such that (µ i , Σ i ) ∈ Ω and w i ≥ α/k, every ( µ, Σ) ∈ T 1 (µ i , Σ i ) satisfies S Ω (( µ, Σ), X) = 0. We still have that the set of points of score at most η * n is contained in the union of at most n m sets T2 (µ ′ i , Σ ′ i ). Thus, V η * (X) V η (X) ≤ exp O d 7/4 log d + m log n ,where we recall that we set η = η * 10 . So, by Theorem C.3, as long asn ≥ C • O(d 7/4 log d + m log n) + log(1/(β 0 • η)) ε 0 • η ,which indeed holds, M(X) outputs some ( µ, Σ) of score at most 2ηn, with 1 -β 0 probability. By the second condition of regularity" there exists an index i ∈[k] such that 1 O(n 8 ) • Σ i ≼ Σ ≼ O(n 8 ) • Σ i and ∥Σ -1/2 i ( µ -µ i )∥ 2 ≤ O(n 6). Moreover, because ( µ, Σ) has finite score, ( µ, Σ) ≈ γ,τ (µ, Σ) for some (µ, Σ) ∈ Ω.Lemma E.6. Lines 7-17 of Algorithm 2 are (ε, δ)-DP. Moreover, for every regular set X, with probability at least 1 -O(k • β 0 ) over the randomness of Algorithm 2, for every i ∈[k] such that w i ≥ α k , there exists j such that n -12 • Σj ≼ Σ i ≼ n 12 • Σj and ∥ Σ-1/2 j (µ -μj )∥ 2 ≤ n 12 .Proof. The privacy guarantee is immediate by (adaptive) advanced composition. Indeed, the algorithm M at each step is (ε 0 , δ 0 )-DP, and only depends on X and Ω, which is determined only by the output of all previous runs of M. Now, let's say that X = {X 1 , . . . , X n } is regular, and we run Algorithm 2 on X. Now, after j ≥ 0 steps of the loop, define P j to be the set of indices i ∈ [k] such that w i ≥ α k and (µ i , Σ i ) ∈ Ω j (where Ω j refers to the set Ω after j steps of the loop are completed). Likewise, define Q j to be the set of indices i such that there exists (μ, Σ) ∈ Ω j with 1n 10 • Σ i ≼ Σ ≼ n 10 • Σ i and ∥Σ -1/2 i (μ -µ i )∥ 2 ≤ n 10 . Note that P 0 ⊂ Q 0 = [k],and that P j ⊂ Q j always. Moreover, because Ω only shrinks, P j+1 ⊂ P j and Q j+1 ⊂ Q j always. Now, after some step j < k, we claim that if P j ̸ = ∅, then |Q j+1 | ≤ |Q j | -1 with failure probability at most O(β 0 ), i.e., Q decreases in size for the next step if P isn't currently empty. The failure probability is only over the randomness of M at each step, and holds for any regular dataset X. A union bound says the total failure probability is O(k • β 0 ).To see why this holds, note that if some index i ∈ P j , then (µ i , Σ i ) ∈ Ω j . So, we can apply Lemma E.5. At step j + 1, we find some (μ j+1 , Σj+1 ), with the following properties. First, there exists an index i′ ∈ [k] with 1 O(n 8 ) • Σ i ′ ≼ Σj+1 ≼ O(n 8 ) • Σ i ′ and ∥Σ -1/2 i ′ (μ j+1 -µ i ′ )∥ 2 ≤ O(n 6 ).Second, (μ j+1 , Σj+1 ) ≈ γ,τ (µ, Σ), for some (µ, Σ) ∈ Ω j .We claim this means i ′ ∈ Q j . To see why, it suffices to show that 1n 10 • Σ i ′ ≼ Σ ≼ n 10 • Σ i ′ and ∥Σ -1/2 i ′ (µ -µ i ′ )∥ 2 ≤ n 10, by definition of Q j and because (µ, Σ) ∈ Ω j . However, we know that1 O(n 8 ) •Σ i ′ ≼ Σj+1 ≼ O(n 8 )•Σ i ′ and 1 2 Σj+1 ≼ Σ ≼ 2 Σj+1 , which implies 1 n 10 •Σ i ′ ≼ Σ ≼ n 10 •Σ i ′ . Also, we know that ∥Σ -1/2 i ′ (μ j+1 -µ i ′ )∥ 2 ≤ O(n 6 ) and ∥Σ -1/2 (µ -μj+1 )∥ 2 ≤ τ ≤ 1 The latter inequality along with the fact that Σ ≼ n 10 • Σ i ′ implies that ∥Σ -1/2 i ′ (µ -μj+1 )∥ 2 ≤ n 5 . So, by triangle inequality, ∥Σ -1/2 i ′ (µ -µ i ′ )∥ 2 ≤ O(n 6 ) ≤ n 10 . Next, we claim that i ′ ̸ ∈ P j+1 , so i ′ ̸ ∈ Q j+1 . Indeed, we have that 1 O(n 8 ) • Σ i ′ ≼ Σj+1 ≼ O(n 8 ) • Σ i ′ , which immediately means 1 n 10 • Σj+1 ≼ 1 O(n 8 ) • Σj+1 ≼ Σ i ′ ≼ O(n 8 ) • Σj+1 ≼ n 10 • Σj+1 . Also, ∥Σ -1/2 i ′ (μ j+1 -µ i ′ )∥ 2 ≤ O(n 6 ), and since Σ i ′ ≼ n 10 • Σj+1 , this means ∥ Σ-1/2 j+1 (µ i ′ -μj+1 )∥ 2 ≤ O(n 11 ) ≤ n 12. Thus, the algorithm will remove(µ i ′ , Σ i ′ ) from Ω at step j + 1, so i ′ ̸ ∈ Q j+1 .Thus, either P j is empty, or Q j decreases in size to Q j+1 for each 0 ≤ j ≤ k + 1. This implies that P k is empty, i.e., for every i ∈ [k] such that w i ≥ α k , (µ i , Σ i ) ̸ ∈ Ω j . This can only happen if each such (µ i , Σ i ) was removed at some point. So, for every such i, there was some index j such thatn -12 • Σj ≼ Σ i ≼ n 12 • Σj and ∥ Σ-1/2 j (µ -μj )∥ 2 ≤ n 12 .We first quickly summarize how to put everything together to prove Theorem 1.4. To prove our main result, need to ensure that our algorithm is private, uses few enough samples, and accurately learns the mixture. Privacy will be simple, as we have shown in Lemma E.6 that the crude approximation is private, and the hypothesis selection is private as shown in Appendix D. The sample complexity will come from the settings of n and n ′ in lines 3 and 18 of the algorithm, and from our setting of parameters. Finally, we showed in Lemma E.6 that we have found a set of (μ j , Σj ), of at most k meancovariance pairs, such that every true (µ i , Σ i ) of sufficiently large weight is crudely approximated by some (μ j , Σj ). Indeed, this is exactly the situation for which we can apply the result on private hypothesis selection (Lemma D.3, based on [BSKW19] ).We are now ready to complete the analysis and prove Theorem 1.4.Proof of Theorem 1.4. By Lemma E.6, note that the algorithm up to Line 17 is (ε, δ)-DP with respect to X 1 , . . . , X n , and does not depend on X n+1 , . . . , X n+n ′ . Assuming {(μ j , Σj )} from these lines is fixed, Lines 18-20 are (ε, δ)-DP with respect to X n+1 , . . . , X n+n ′ , by Lemma D.3, and do not depend on X 1 , . . . , X n . So, the overall algorithm is (ε, δ)-DP.Next, we verify accuracy. Note that by Lemma E.6, and for G = n 12 , the sets (μ j , Σj ) that we find satisfy the conditions for Lemma D.3, with failure probability at most O(k • β 0 ). This probability is at most 0.1, assuming e d is significantly larger than k. So, it suffices for n ′ , the number of samples used in Line 18 of Algorithm 2, to satisfyn ′ ≥ O d 2 k•log(G•k √ d/α) α 2 + d 2 k•log(G•k √ d/α) αε = O d 2 k α 2 + d 2 kαε , where the last part of the bound holds by our assumptions on G and n. Thus, the total sample complexity isn + n ′ = O kd 2 α 2 + kd 2 + d 1.75 k 1.5 log 0.5 (1/δ) + k 1.5 log 1.5 (1/δ) αε + k 2 d α .By Lemma D.3, with failure probability at most 0.1, Lines 18-20 of the algorithm will successfully output a hypothesis which is a mixture of k Gaussians, with total variation distance at most O(α)from the right answer. So, the overall success probability is at least 0.8.Finally, we note that the assumption that e d is much larger than k can be made WLOG, by padding d to be a sufficiently large multiple of log k. Namely, if given samples X i = (X i,1 , . . . , X i,d ), we can add coordinates X i,d+1 , . . . , X O(log k) ∼ N (0, 1). Then, if the original samples were a mixture of kGaussians N (µ i , Σ i ), the new distribution is a mixture of k Gaussians N µ i 0 , Σ i 0 0 I , with the same mixing weights. So, we can learn the mixture distribution up to total variation distance α in the larger space, and then remove all except the first d coordinates, to output our final answer. This will not affect the sample complexity by more than a poly(log k) factor.In this section, we provide the algorithm and analysis for Theorem 1.5. We first describe and provide pseudocode for the algorithm, and then we prove that it can privately learn mixtures of Gaussians with low sample complexity.We will use σ i := √ Σ i to denote the standard deviation of a univariate Gaussian N (µ i , Σ i ).As in the algorithm in Appendix E, we will actually draw a total of n + n ′ samples. We start off by only considering the first n data points X = {X 1 , . . . , X n }, where eachX j ∈ R. Now, let Y = {Y 1 , . . . , Y n } be the sorted version of X, i.e., Y 1 ≤ Y 2 ≤ • • • ≤ Y n ,and (X 1 , . . . , X n ) and (Y 1 , . . . , Y n ) are the same up to permutation. Finally, for each j ≤ n -1, define Z j to be the ordered pair (Y j , Y j+1 -Y j ), and define Z = Z(X) to be the unordered multiset of Z j 's. (Note that Z depends deterministically on X.)The algorithm now works as follows. Suppose we have data X = {X 1 , . . . , X n }. After sorting to obtain Y 1 , . . . , Y n , let r j = Y j , s j = Y j+1 -Y j , and Z j = (r j , s j ). So, Z(X) = {(r j , s j )} 1≤j≤n-1 , viewed as an unordered set. Note that every s j ≥ 0, and the r j 's are in nondecreasing order. We will create a set of buckets that can be bijected onto Z 2 . For each Z j = (r j , s j ), if s j > 0 we assign Z j to the bucket labeled (a, b) ∈ Z 2 if 2 a ≤ s j < 2 a+1 and b • n 5 • 2 a ≤ r j < (b + 1) • n 5 • 2 a . If s j = 0, we do not assign Z j to any bucket.For each element e ∈ Z 2 , we keep track of the number of indices j ∈ [n -1] with Z j sent to e. In other words, for e = (a, b), we define the count c e = #{j :2 a ≤ s j < 2 a+1 , b • n 5 • 2 a ≤ r j < (b + 1) • n 5 • 2 a }.For each c e , we sample an independent draw g e ∼ TLap(1, ε/10, δ/10), and define ce = c e + g e . Finally, we let S = {(μ i , Σi )} be the the set of pairs (b • n 5 • 2 a , 2 2a ) where e = (a, b) satisfies ce > 100 ε log 1 δ . Hence, we have computed some {(μ i , Σi )} k ′ i=1 , for some k ′ ≥ 0. (We will show in the analysis that k ′ ≤ n -1 always, so k ′ is finite.) If k ′ = 0 we simply output ⊥. Otherwise, we will draw Õ k α 2 + k αε fresh samples. We run the ε-DP hypothesis selection based algorithm, Algorithm 1 on the fresh samples, using parameters G = n 10 and {(μ i , Σi )} k ′ i=1 .Pseudocode. We give pseudocode for the algorithm in Algorithm 3.We start by analyzing Lines 1-19, which generates the set S = {(μ i , Σi )} of candidate meancovariance pairs.First, we note a very basic proposition. with probability 1. For ε, δ ≤ 1, this is less than 100ε • log 1 δ . probability at most 1 1000k , at least α 4k • n points X r : r ∈ T i are in the interval [µ i -σ i , µ i + σ i ]. Next, note that for any X r , X r ′ ∼ N (µ i , σ 2 i ), X r -X r ′ ∼ N (0, 2σ 2 i ), and thus has magnitude at least σi 10 4 n 3 with at most 1 1000n 3 failure probability. So, by a union bound over all r ̸ = r ′ ∈ T i , with at most 1 1000n ≤ 1 1000k failure probability, all data points X r , X r ′ for r ̸ = r ′ ∈ T i are separated by at least σi 10 4 n 3 . So, with at least 3 1000k failure probability, there is a subset U i ⊂ T i of size at least α 4k • n, such that for every r, r ′ ∈ U i , |X r -X r ′ |, and for every r ∈ U i , X r ∈ [µ i -σ i , µ i + σ i ]. By a union bound, this holds simultaneously for all i ∈ [k] with w j ≥ α k , with probability at least 0.997. Conditioned on this event, if we consider the elements X in sorted order (i.e., Y), between every consecutive pair X r , X r ′ : r, r ′ ∈ U i (after sorting), we know X r ′ -X r ≥ σi 10 4 •n 3 . So, there exist someX r ≤ Y j , Y j+1 ≤ X r ′ such that Y j+1 -Y j ≥ σi 10 4 •n 4 . Hence, for every i ∈ [k] with w i ≥ α k , there exist at least |U i | -1 ≥ α 4k -1 ≥ α 8k indices j such that Y j+1 -Y j ≥ σi 10 4 •n 4 . Moreover, note that Y j , Y j+1 ∈ [µ i -σ i , µ i + σ i ], so Y j+1 -Y j ≤ 2σ i .Hence, conditioned on the event from the previous paragraph, the lemma holds.Given the above lemma, we can prove that for every (µ i , Σ i ) in the Gaussian mixture with at least α k weight, there is some corresponding crude approximation (μ, Σ) ∈ S.Lemma F.4. Assume that the conditions and result of Lemma F.3 holds. Then, there exists (μ, Σ) ∈ S such that n-10 • Σ ≼ Σ i ≼ n 10 • Σ and ∥ Σ-1/2 (µ i -μ)∥ 2 ≤ n 6 .Proof. Fix any i with w i ≥ α k . Then, for any j with Y j , Y j+1 ∈ [µ i -σ i , µ i +σ i ] and σi 10 4 •n 4 ≤ Y j+1 -Y j ≤ 2σ i , the index i contributes to the count of some bucket e = (a, b), where ⌊log 2 σi 10 5 •n 4 ⌋ ≤ a ≤ ⌊log 2 (2σ i )⌋. Moreover, Y j ∈ [b • n 5 • 2 a , (b + 1) • n 5 • 2 a ). Therefore, if we consider the set V i of e = (a, b) such that ⌊log 2 σi 10 5 •n 4 ⌋ ≤ a ≤ ⌊log 2 (2σ i )⌋ and [b • n 5 • 2 a , (b + 1) • n 5 • 2 a ) ∩ [µ i -σ i , µ i + σ i ] isV i such that (μ, Σ) := (b • n 5 • 2 a , 2 2a ) in S. By definition of V i , σi n 5 ≤ 2 a ≤ 2σ i , so 1 n 10 • Σ i ≤ 2 2a ≤ 4Σ i ≤ n 10 • Σ i . Moreover, | Σ-1/2 (µ i -μ)| = 2 -a • |µ i -μ|. But the definition of V i means μ = b • n 5 • 2 a satisfies μ ≤ µ i + σ i and μ ≥ µ i -σ i -n 5 • 2 a ≥ µ i -(2n 5 + 1) • σ i . So, |μ -µ i | ≤ (2n 5 + 1) • σ i ≤ n 6 • 2 a . Thus, | Σ-1/2 (μ -µ i )| ≤ n 6 .We now prove that the mechanism (until Line 19) is differentially private. First, we note the following auxiliary claim. Lemma F.5. Let X, X ′ be adjacent datasets of size n (i.e., only differing on a single element). Then, the corresponding sets Z = Z(X) and Z ′ = Z(X ′ ) differ in distance at most 3, where by distance we mean that there exists a permutation of the elements in Z and Z ′ , respectively, such that at most three indicesi ≤ n -1 satisfy Z i ̸ = Z ′ i .Proof. Note that for the sorted versions Y, Y ′ of X, X ′ , respectively, we can convert from Y to Y ′ by removing one data point and adding one more data point, without affecting the order of any other data points.Suppose we remove some Y j from Y . If j = 1, then this just removes Z 1 , and if j = n, then this just removes Z n-1 . If j ≥ 2, this modifies Z j-1 and removes Z j . Likewise, if we add some new Y j ′ , this will either add one new ordered pair to Z (if Y j ′ is either the smallest or largest element), or replace one ordered pair in Z with two new pairs. Therefore, if we remove a Y j and then add a Y j ′ , this will change at most 3 of the ordered pairs in Z.We now prove privacy.making learning the full mixture equivalent to learning each component (on average). Namely, even if we are given the information of which Gaussian each sample comes from, we will need to learn at least 2/3 of the Gaussians up to total variation distance O(α), to learn the full mixture up to total variation distance α. Hence, we will need at least k times as many samples as for learning a single Gaussian, which means we need Ω kd 2 αε total samples.The lower bound of k log(1/δ) αε will follow from Theorem G.1. Note that it suffices to prove the lower bound in the univariate case. We plant k distinct Gaussians N (µ i , 1), where the µ i are very far away from each other, i.e., pairwise |µ i -µ j | ≫ (1/δ) 10C . We also assume that µ 1 = 0 is known, and the remaining µ i are unknown but we are promised the value of each µ i up to error (1/δ) C . The overall mixture will have the first Gaussian N (0, 1) of weight w 1 = 1 -α/c * , and the remaining Gaussians N (µ i , 1) each have weight w i = α/(c * • (k -1)). Even if we are given the information of which Gaussian component each sample comes from, to learn the overall mixture up to error α, we need to learn at least 2/3 of the small-weight components, each up to total variation distance O(c * ). Hence, we will need log(1/δ) ε samples from most of the small-weight components. Since the small weight components have weight Θ(α/k), we need Ω k log(1/δ) αε total samples. We now give a formal proof of Theorem 1.6.First, we note two lemmas that will be helpful in proving the Theorem. Lemma G.3. Let α ∈ [0, 1). Let f be a probability density function over R d . Assume g : R d → R ≥0 exists such thatR d |f (x) -g(x)| dx ≤ α.Then there exists h such that h is a probability density function over R d , andR d |f (x) -h(x)| dx ≤ 2α. Proof. We know that R d f (x) dx - R d g(x) dx ≤ R d |f (x) -g(x)| dx ≤ α.Therefore, G := R d g(x) dx = 1 ± α. Now two cases are possible either G < 1, or G > 1, otherwise we are done. If G < 1, let h(x) be the density function corresponding to the following distribution: with probability G take a sample from g(x)/G, and with probability 1 -G select 0. ThenR d |f (x) -h(x)| dx ≤ 1 -G + R d |f (x) -g(x)| dx ≤ 2α, as desired. If G > 1, let h(x)be a density function as follows: ∀x : 0 ≤ h(x) ≤ g(x), and R d h(x) dx = 1. It is easy to see such an h exists by greedily picking h. Then we haveR d |g(x) -h(x)| dx = R d g(x) -h(x) dx ≤ G -1 ≤ α. Therefore, we can write R d |f (x) -h(x)| dx ≤ R d |f (x) -g(x)| dx + R d |g(x) -h(x)| dx ≤ 2α,as desired.Lemma G.4. Suppose w ∈ (0, 1) is fixed. Suppose an (ε, δ) differentially private algorithm exists that takes n samples from the mixture D = wD 1 + (1 -w)D 2 , and learns D 1 in total variation distance up to error α, with success probability 1 -β. Moreover, assume while sampling it is known which component the sample is sampled from. Then there exists an (ε, δ) differentially private algorithm A that takes nw/γ samples from D 1 and outputs an estimate of D 1 up to total variation distance α, with success probability 1 -β -γ.where the last inequality comes from taking M ≥ 30 log(k/100α). We show that using the distribution D, we can construct an answer to the problem of learning the Gaussian D i . To do so first take g i to be equal to k f (x) over B i and 0 everywhere else. We haveR d |f i (x) -g i (x)| dx = Bi f i (x) -k f (x) dx + R d \Bi f i (x) dx ≤ 300α + 100α/k ≤ 400α.Now we may apply Lemma G.3, and deduce that given D we can construct probability density functions and distributions Di 's such that d TV ( Di , D i ) ≤ 800α, for all i ∈ G. To recap, so far we have shown that given an (ε, δ) differentially private algorithm that takes as inputs samples from our constructed mixture of Gaussians and outputs a density D that has total variation distance at most α, from the ground truth distribution with success probability 2/3, we can use D to construct densities Di such that d TV ( Di , D i ) ≤ 800α, for 0.99k of the indices i. This implies that there exists a fixed index i for which the component D i is learned up to error 800α with success probability 2/3 -0.01 ≥ 0.65. Applying Lemma G.4, implies that there must exist an (ε, δ) differentially private algorithm that takes 100n/k samples from D i and estimates its density up to total variation distance 800α, with success probability at least 0.6. Therefore, applying Theorem G.2, we conclude that n = Ω( kd 2 αε ). Now let's prove the last term k log(1/δ) αε . We aim to apply Theorem G.1. Let µ i 's be k distinct values in R, each having distance M from each other, for M ≫ (1/δ) 10C to be set later, where µ 1 = 0. It is easy to see such a set exists. Now consider the following set of Gaussians: D i = N (µ i , 1), where µ i 's are known up to log(1/δ) C , and µ 1 = 0 is also known. Consider the mixture D over these Gaussians, with weights w 1 = 1 -α/c * , and w i = α/(c * (k -1)). We also assume that when sampling from the mixture we know which component each sample comes from. Consider an (ε, δ) differentially private algorithm A that takes n samples from D and outputs a distribution D such that d TV (D, D) ≤ α with probability 2/3. Now consider a sample from D i , from standard Gaussian tail bounds we know that at least 1exp(-M 2 /200) fraction of the mass of D i is contained within a ball of radius M/10, around µ i . Let B i denote this ball, and note that B i 's are disjoint.Let f, f i , f be the probability density functions corresponding to D, D i , D respectively. Assuming, we are in the success regime, similar to the proof of the previous term, we can show that there exists a set G of indices such that |G| ≥ 0.99k, and ∀i ∈ G : Bi f (x) -f (x) dx ≤ 200α/k. Moreover, with a similar argument as the previous term for i ∈ G we can say as long as M ≥ 20 log(k/100α), Bi w i f i (x) -f (x) dx ≤ 300α/k. We show that using the distribution D, we can construct an answer to the problem of learning the Gaussian D i , for i ̸ = 1. To do so take g i to be equal to f (x)/w i , over B i and 0 everywhere else. We haveR d |f i (x) -g i (x)| dx = Bi f i (x) -f (x)/w i dx+ R d \Bi f i (x) dx ≤ 300α kw i + 100α k ≤ 400c * .Now we may apply Lemma G.3, and deduce that given D, we can construct probability density functions and distributions Di 's such that d TV ( Di , D i ) ≤ 800c * , for all i ∈ G. To recap, so far we have shown that given an (ε, δ) differentially private algorithm that takes as inputs samples from our constructed mixture of Gaussians and outputs density D that has total variation distance at most α from the ground truth distribution with success probability 2/3, we can use D to construct densities Di such that d TV ( Di , D i ) ≤ 800c * , for 0.99k of the indices i. This implies that there exists a fixed index i ̸ = 1, for which the component D i is learned up to error 800c * , with success probability 2/3 -0.01 ≥ 0.65. Applying Lemma G.4, implies that there must exist an (ε, δ) differentially private algorithm that takes 100nw i samples from D i and estimates its density up to total variation distance 800c * , with success probability 0.6. Therefore, applying Theorem G.1, and noting that w i = α/(c * (k -1)) we conclude that n = Ω( k log(1/δ)).In this section, we prove Proposition B.2, Theorem B.3, and Lemma C.7.H.1 Proof of Proposition B.2First, we note a basic fact.Fact H.1. For any x, y ∈ [0.9, 1.1], we have that (xy -1) 2 ≤ 4((x -1) 2 + (y -1) 2 ). First, assume (µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ). This means ∥J 1 J ⊤ 1 -I∥ op ≤ γ and ∥J 1 J ⊤ 1 -I∥ F ≤ ρ. We then have that Σ-1/2 1 Σ 2 Σ -1/2 1 = J -11 (J -1 1 ) ⊤ = (J ⊤ 1 J 1 ) -1 . Now, note that J 1 J ⊤ 1 and J ⊤ 1 J 1 are both symmetric and have the same eigenvalues. If we call these eigenvalues λ 1 , . . . , λ d , then the eigenvalues of Σ -1/2 1 Σ 2 Σ -1/2 1 are λ -1 1 , . . . , λ -1 d . Now, our assumption (µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ) implies that 1 -γ ≤ λ i ≤ 1 + γ and (1 -λ i ) 2 ≤ ρ 2 . This means that, assuming γ ≤ 0.1, 1 -2γ ≤ 1 1+γ ≤ λ -1 i ≤ 1 1-γ ≤ 1 + 2γ, and1 -λ -1 i 2 ≤ (1 -λ i ) 2 • λ -2 i ≤ 2 • (1 -λ i ) 2 = 2ρ 2 .This means that ∥Σ-1/2 1 Σ 2 Σ -1/2 1 ∥ op ≤ 2γ and ∥Σ -1/2 1 Σ 2 Σ -1/2 1 ∥ F ≤ 2ρ.Finally, Σ -1/2 1 (µ 1 -µ 2 ) = J -1 1 Σ -1/2 2 (µ 1 -µ 2 ). Because Σ -1/2 2 (µ 1 -µ 2 ) has magnitude at most τ by our assumption, J -1 1 Σ -1/2 2 (µ 1 -µ 2 ) has magnitude at most the maximum singular value of J -1 1 times τ . But every singular value of J -1 1 is some λ -1/2 i which is at most 2, so ∥Σ-1/2 1 (µ 2 -µ 1 )∥ 2 = ∥J -1 1 Σ -1/2 2 (µ 1 -µ 2 )∥ 2 ≤ 2τ .Next, assume (µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ) and (µ 2 , Σ 2 ) ≈ γ,ρ,τ (µ 3 , Σ 3 ). First, note that Σ -1/2 3 Σ 1 Σ -1/2 3 = J 2 J 1 J ⊤ 1 J ⊤ 2 = (J 2 J 1 )(J 2 J 1 ) ⊤ . If the eigenvalues of J 1 J ⊤ 1 are {λ i } and the eigenvalues of J 2 J ⊤ 2 are {λ ′ i }, then the singular values of J 1 and J 2 are { √ λ i } and { λ ′ i }, respectively. By our assumption, 1-γ ≤ λ i , λ ′ i ≤ 1+γ, which means that√ 1 -γ ≤ √ λ i , λ ′ i ≤ √ 1 + γ.Thus, the singular values of J 2 J 1 are between 1 -γ and 1 + γ, which means that the eigenvalues of (J 2 J 1 )(J 2 J 1 ) ⊤ are between (1 -γ) 2 and (1 + γ) 2 . Hence, for γ ≤ 0.1, ∥J 2 J 1 J ⊤ 1 J ⊤ 2 -I∥ op ≤ 4γ. Assume that λ i , λ ′ i are in decreasing order. We now consider the k th largest singular value of J 2 J 1 . If σ k := √ λ k is the k th largest singular value of J 1 and σ ′ k := λ ′ k is the k th largest singular value of J 2 , by Corollary A.5 there exist subspaces V k , V ′ k of dimension d -k + 1 such that ∥J 1 v∥ 2 ≤ σ k ∥v∥ 2 for all v ∈ V k and ∥J 2 v∥ 2 ≤ σ ′ k ∥v∥ 2 for all v ∈ V ′ k . Therefore, for every v ∈ V k ∩ J -1 1 V ′ k (note that J 1 is invertible since J 1 J ⊤ 1 has all eigenvalues between 1 -γ and 1 + γ) we have that ∥J 2 J 1 v∥ 2 ≤ σ ′ k ∥J 1 v∥ 2 ≤ σ k σ ′ k ∥v∥ 2 . Because V k and J -1 1 V ′ k both have dimension d -k + 1, their intersection has dimension at least d -2k + 2. So, there is a subspace of dimension at least d -2k + 2 such that every v in the subspace has ∥J 2 J 1 v∥ 2 ≤ σ k σ ′ k • ∥v∥ 2 . Thus, the (2k -1) th largest singular value of J 2 J 1 is at most σ k σ ′ k , so the (2k -1) th largest eigenvalue of J 2 J 1 J ⊤ 1 J ⊤ 2 is at most λ k λ ′ k . The same argument, looking at the smallest singular values, tells us that the (2k -1) th smallest eigenvalue of J 2 J 1 J ⊤ 1 J ⊤ 2 is at least λ d-k+1 λ ′ d-k+1 . Thus, for any t, the t th largest eigenvalue of J 2 J 1 J ⊤ 1 J ⊤ 2 is at most λ ⌊(t+1)/2⌋ λ ′ ⌊(t+1)/2⌋ and at least λ ⌈(d+t)/2⌉ λ ′ ⌈(d+t)/2⌉ . For P , we can choose a 1/n 10 -sized net over the Frobenius norm metric (i.e., the distance between two matrices P 1 , P 2 is ∥P 1 -P 2 ∥ F ) for each of the k nonzero eigenvalues and eigenvectors in the unit d-dimensional sphere, which has size at most n 100d•k . Therefore, by a union bound, the probability that 1 n i∈S ⟨X i X ⊤ i -I, P ′ ⟩ ≤ 0.01 for every |S| ≥ (1 -c)n and every P ′ in the net is at least 1 -(2e/c) cn • n 100d•k • e -c1n/2 . Finally, we consider P outside of the net. For any symmetric P of rank k and Frobenius norm 1, it has Frobenius distance at most 1/n 8 from some P ′ in the net. Let us consider the event that every ∥X i ∥ 2 2 ≤ 10n, which for n ≥ d occurs with failure probability at most 2n • e -c1n by Hanson-Wright. Under this event, ⟨X i X ⊤ i -I, P -P ′ ⟩ ≤ ∥X i X ⊤ i -I∥ F • ∥P -P ′ ∥ F ≤ (10n + √ d)/n 8 ≤ 0.01.As a result, with failure probability at most (2e/c) cn •n 100d•k •e -c1n/2 +2n•e -c1n ≤ e -cn (assuming c is sufficiently small), we have both properties. Namely, 1 n i∈S ⟨X i X ⊤ i -I, P ′ ⟩ ≤ 0.01 for every |S| ≥ (1 -c)n and every P ′ in the net, and for any P , ⟨X i X ⊤ i -I, P -P ′ ⟩0.01 for the closest P ′ in the net to P . Overall, this means that for all such P and S, 1 n i∈S ⟨X i X ⊤ i -I, P ⟩ ≤ 0.1.We prove another lemma which contrasts with Lemma H.3. Lemma H.4. Fix any 2 ≤ ρ ≤ √ d. Let k = 4d/ρ 2 , and let n ≥ O(d • k) and c ∈ (0, 0.1) be as in Lemma H.3. Fix any covariance matrix Σ and let X = {X 1 . . . , X n } ∼ N (0, Σ). Then, with probability at least 1 -e -cn , for every Σ such that 0.95 • Σ ≼ Σ ≼ 1.05 • Σ and ∥Σ -1/2 ΣΣ -1/2 -I∥ F ≥ ρ, for every symmetric matrix P ∈ R d×d of rank at most k and Frobenius norm 1, and for every S ⊂ [n] of size at least (1 -c)n, 1 n i∈S ⟨ Σ -1/2 X i X ⊤ i Σ -1/2 -I, P ⟩ ≥ 0.7.Proof. Let J = Σ -1/2 Σ 1/2 , and let Y i = Σ -1/2 X i . Then, Σ -1/2 X i = JY i , which means that Σ -1/2 X i X ⊤ i Σ -1/2 -I = JY i Y ⊤ i J ⊤ -I = J(Y i Y ⊤ i -I)J ⊤ + (JJ ⊤ -I).We are slightly abusing the notation and using D(x) as the pdf of D at points x.For clarity of presentation, we assume the covariance matrices are not singular. However, extending our results to degenerate matrices is straightforward.To obtain success probability β with O(n • log 1/β) samples, we repeat the procedure T = O(log 1/β) times on independent groups of n samples each, to get T estimates D1, . . . , DT , and by a Chernoff bound, at least51% of the estimates are within total variation distance α of the true mixture D. So, by choosing an estimate that is within 2α of at least 51% of the estimates, it is still within 3α total variation distance of D. 
technical overview and roadmap: 
 We highlight some of our conceptual and technical contributions. We mainly focus on the highdimensional upper bound, and discuss the univariate upper bound at the end.Suppose we are promised a bound on the means and covariances of the Gaussian components, i.e.,1 R • I ≼ Σ i ≼ R • I and ∥µ i ∥ 2 ≤ R for all i ∈ [k].i ∈ [k] we know some Σi such that 1 R • Σ i ≼ Σi ≼ R • Σ i , then it suffices to have n = O kd 2 log R α 2+ kd 2 log R αε samples to learn the full GMM in total variation distance.Our main goal will be to learn every covariance Σ i with such an approximation, for R = poly k, d, 1 α , 1 ε , so that log R can be hidden in the Õ factor. To explain why this goal is sufficient, suppose we can crudely learn every covariance Σ i with approximation ratio R, using n ′ samples. We then needO kd 2 log R α 2 + kd 2 log R αε = Õ kd 2 α 2 + kd 2αε additional samples to learn the full distribution using hypothesis selection, so the total sample complexity isÕ n ′ + kd 2 log R α 2 + kd 2 log R αε. Hence, we will aim for this easier goal of crudely learning each covariance, for both Theorem 1.4 and Theorem 1.5, using as few samples n ′ as possible. We will also need to approximate each mean µ i , though for simplicity we will just focus on covariances in this section.The main goal will be to simply provide a rough approximation to the set of standard deviations σ i = √ Σ i , as we can finish the procedure with hypothesis selection, as discussed above.Bird's eye view: Say we are given a dataset X = {X 1 , . . . , X n }: note that every X j ∈ R since we are dealing with univariate Gaussians. The main insight is to sort the data in increasing order (i.e., reorder so thatX 1 ≤ X 2 ≤ • • • ≤ X n ) and consider the unordered multiset of successive differences {X 2 -X 1 , X 3 -X 2 , . . . , X n -X n-1 }.One can show that if a single datapoint X j is changed, then the set of consecutive differences (up to permutation) does not change in more than 3 locations (see Lemma F.5 for a formal proof).We then apply a standard private histogram approach. Namely, for each integer a ∈ Z, we create a corresponding bucket B a , and map eachX j+1 -X j into B a if 2 a ≤ X j+1 -X j < 2 a+1 . If some mixture component i had variance Σ i = σ 2 i ,we should expect a significant number of X j+1 -X j to at least be crudely close to σ i , such as for the X j drawn from the i th mixture component. So, some corresponding bucket should be reasonably large. Finally, by adding noise to the count of each bucket and taking the largest noisy counts, we will successfully find an approximation to all variances.In more detail: For each (possibly negative) integer a let c a be the number of indices i such that 2 a ≤ X j+1 -X j < 2 a+1 . We will prove that, if the weight of the i th component in the mixture is w i and there are n points, then we should expect at least Ω(w i • n) indices j to be in a bucket a with 2 a within a poly(n) multiplicative factor of the standard deviation σ i (see Lemma F.3). The point of this observation is that there are at most O(log n) buckets B a with 2 a between σi poly(n) and O(σ i ), but there are at least Ω(w i • n) indices mapping to one of these buckets. So by the Pigeonhole principle, one of these buckets has at least Ω wi•n log n indices, i.e., c a ≥ Ω wi•n log n for some a with σipoly(n) ≤ 2 a ≤ O(σ i ).Moreover, we know that if we change a single data point X j , the set of consecutive differences {X j+1 -X j } after sorting changes by at most 3. So, if we change a single X j , at most 3 of the counts c a can change, each by at most 3.Now, for every integer a, draw a noise value from the Truncated Laplace distribution (see Definition A.2 and Lemma A.3), and add it to c a to get a noisy count ca . The details of the noise distribution are not important, but the idea is that this distribution is always bounded by O 1 ε log 1 δ . Moreover, the Truncated Laplace distribution preserves (ε, δ)-DP. This means that the counts {c a } a∈Z will have (O(ε), O(δ))-DP, because the true counts c a only change minimally across adjacent datasets.Our crude approximation to the set of standard deviations will be the set of 2 a such that ca exceeds some threshold which is a large multiple of1 ε log 1 δ . If n ≥ Õ k log(1/δ) αεand the weight w i ≥ α/k, it is not hard to verify that wi•n log n exceeds a large multiple of 1 ε log 1 δ . So, for each i ≤ k with weight at least α/k, some corresponding a with σi poly(n) ≤ 2 a ≤ O(σ i ) will have count c a significantly exceeding the threshold, and thus noisy count ca exceeding the threshold. This will be enough to crudely approximate the values σ i coming from large enough weight. We can ignore any component with weight less than α/k, as even if all but one of components have such small weight, together they only contribute α weight. So, we can ignore them and it will only cost us α in total variation distance, which we can afford.Putting everything together: In summary, we needed Õ k log(1/δ) αε samples to approximate each covariance (of sufficient weight) up to a poly(n) multiplicative factor. By setting R = poly(n) and using the reduction described in Section 2.1, we then need an additional O k log n α 2+ k log n αε samples. If we set n = Õ k α 2 + k log(1/δ) αε , we will obtain that n ≥ Õ k log(1/δ) αε + O k log n α 2 + k log n αε, which is sufficient to solve our desired problem in the univariate setting.Note that this proof relies heavily on the use of private histograms and the order of the data points in the real line. Therefore, it cannot be extended to the high-dimensional setting. We will use a completely different approach to prove Theorem 1.4.As in the univariate case, we only need rough approximations of the covariances. We will learn the covariances one at a time: in each iteration, we privately identify a single covariance Σi which crudely approximates some covariance Σ i in the mixture, with (ε/ k log(1/δ), δ/k)-DP. Using the well-known advanced composition theorem (see Theorem A.1), we will get an overall privacy guarantee of (ε, δ)-DP. However, to keep things simple, we will usually aim for (ε, δ)-DP when learning each covariance, and we can replace ε with ε/ k log(1/δ) and δ with δ/k at the end.A natural approach for parameter estimation, rather than learn Σ i one at a time, is to learn all of the covariances together. However, we believe that this approach would cause the sample complexity to multiply by a factor of k, compared to learning a single covariance. The advantage of learning the covariances one at a time is that we can apply advanced composition: this will cause the sample complexity to multiply by roughly k log(1/δ) instead.In the rest of this subsection, our main focus is to identify a single crude approximation Σi .Applying robustness-to-privacy: The first main insight is to use the robustness-to-privacy conversion of Hopkins et al. [HKMN23] (see also [AUZ23] ). Hopkins et al. prove a black-box (but not computationally efficient) approach that can convert robust algorithms into differentially private algorithms, using the exponential mechanism and a well-designed score function. This reduction only works for finite dimensional parameter estimation and therefore cannot be applied directly to density estimation. On the other hand, parameter estimation for arbitrary GMMs requires exponentially many samples in the number of components [MV10] . However, we will demonstrate that this lower bound does not hold when we only need a very crude estimation of the parameters.The idea is the following. For a dataset X = {X 1 , . . . , X n }, let S = S( Σ, X) be a score function, which takes in a dataset X of size n and some candidate covariance Σ, and outputs some non-negative integer. At a high level, the score function S( Σ, X) will be the smallest number of data points t that we should change in X to get to some new data set X ′ with a specific desired property: X ′ should "look like" a sample generated from a mixture distribution with Σ being the covariance of one of its components -namely, a component with a significant mixing weight. One can define "looks like" in different ways, and we will adjust the precise definition later. We remark that this notion of score roughly characterizes robustness, because if the samples in X were truly drawn from a Gaussian mixture model with covariances {Σ i } k i=1 , we should expect S(Σ i , X) to be 0 (since X should already satisfy the desired property), but if we altered k data points, the score should be at most k. The high-level choice of score function is somewhat inspired by a version of the exponential mechanism called the inverse sensitivity mechanism [AD20b, AD20a] , though the precise way of defining the score function requires significant care and is an important contribution of this paper.The robustness-to-privacy framework of [HKMN23] , tailored to learning a single covariance, implies the following general result, which holds for any score function S following the blueprint above. For now, we state an informal (and slightly incorrect) version. Theorem 2.1 (informal -see Theorem C.3 for the formal statement). For any η ∈ [0, 1), and any dataset X of size n, define the value V η (X) to be the volume (i.e., Lebesgue measure) of the set of covariance matrices Σ (where the covariance can be viewed as a vector by flattening), such that S( Σ, X) ≤ η • n.Fix a parameter η < 0.1, and suppose that for any dataset X of size n such that V η/2 (X) is strictly positive, the ratio of volumes V η (X)/V η/2 (X) is at most some K (which doesn't depend on X). Then, if n ≥ log K ε•η , there is a differentially private algorithm that can find a covariance Σ of low score (i.e., where S( Σ, X) ≤ η • n) using n samples.Note that Theorem 2.1 does not seem to say anything about whether X comes from a mixture of Gaussians. However, we aim to instantiate this theorem with a score function that is carefully designed for GMMs. Recall that we want S( Σ, X) to capture the number of points in X that need to be altered to make it look like a data set that was generated from a mixture, with Σ being the covariance of one of the Gaussian components. In other words, S( Σ, X) should be small if (and hopefully only if) a "mildly corrupted" version of X includes a subset of points that are generated from a Gaussian with covariance Σ. At the heart of designing such a score function, one needs to come up with a form of "robust Gaussianity tester" that tells whether a given set of data points are generated from a Gaussian distribution. Aside from this challenge, the volume ratio associated with the chosen score function needs to be small for every dataset X otherwise the above theorem would require a large n (i.e., number of samples). These two challenges are, however, related. If the robust Gaussianity tester has high specificity-i.e., rejects most of the sets that are not generated from a (corrupted) Gaussian-then the volume ratio is likely to be small (i.e., a smaller number of candidates Σ would receive a good/low score).First Attempt: We first try an approach which resembles that of [HKMN23] for privately learning a single Gaussian. We "define" S( Σ, X) as the smallest integer t satisfying the following property: there exists a subset Y ⊂ X of size n/k, such that we can change t data points from Y to get to Y ′ , where Y ′ "looks like" i.i.d. samples from a Gaussian with some covariance Σ that is "similar to" Σ. The choice of Y having size n/k is motivated by the fact that each mixture component, on average, has n/k data points in X.The notions of "looks like" and "similar to" of course need to be formally defined. We say Σ is "similar to" Σ (or Σ ≈ Σ) if they are spectrally close, i.e., 0.5Σ ≼ Σ ≼ 2Σ. We say that Y ′ "looks like" samples from a Gaussian with covariance Σ if some covariance estimation algorithm predicts that Y ′ came from a Gaussian with covariance Σ. The choice of covariance estimation algorithm will be quite nontrivial and ends up being a key ingredient in proving Theorem 1.4.To apply Theorem 2.1, we first set η = c/k for some small constant c. We cannot set a larger value η, because if we change t ≈ n/k data points, we could in fact create a totally arbitrary new Gaussian component with large weight. Since there is no bound on the eigenvalues of the covariance matrix, this could cause the volume V η (X) to be infinite. The main question we must answer is how to bound the volume ratio V η (X)/V η/2 (X). To answer this question, we first need to understand what it means for S( Σ, X) ≤ η • n. If S( Σ, X) ≤ η • n, then there exists a corresponding set Y ⊂ X of size n/k, and we can change η • n = c • |Y| points from Y to get to some Y ′ which looks like samples from a Gaussian with covariance Σ ≈ Σ. Thus, Y looks like c-corrupted samples from such a Gaussian (i.e., a c fraction of the data is corrupted). This motivates using a robust covariance estimation algorithm: indeed, robust algorithms can still approximately learn Σ even if a small constant fraction of data is corrupted, so for any Y ⊂ X, we expect that no matter how we change a c fraction of Y to obtain Y ′ , the robust algorithm's covariance estimate should not change much. So, for any fixed Y, the set of possible Σ, and thus the set of possible Σ, should not be that large.In summary, to bound V η (X) versus V η/2 (X), there are at most n n/k choices for Y ⊂ X in the former case, and at least 1 choice in the latter case (since we assume V η/2 (X) > 0 in Theorem 2.1). Moreover, for any such Y, the volume of corresponding Σ should be exponential in d 2 (either for V η (X) or V η/2 (X)), since the dimensionality of the covariance is roughly d 2 . So, this suggests that the overall volume ratio is at mostn n/k • e O(d 2 ) . Since log n n/k ≈ (n/k) • log k, if we plug into Theorem 2.1 it suffices to have n ≥ (n/k) log k+d 2 ε•(c/k). Unfortunately this is impossible unless ε ≥ log k. These ideas will serve as a good starting point, though we need to improve the volume ratio analysis. To do so, we also modify the robust algorithm, by strengthening the assumptions on what it means for samples to "look like" they came from a Gaussian.To improve the volume ratio, we draw inspiration from a technique called sample compression, which has been used in previous work on non-private density estimation for mixtures of Gaussians [ABH + 18, ABDH + 20]. The idea behind sample compression is that one does not need the full set Y to do robust covariance estimation; instead, we look at a smaller set of samples. For instance, if Y ⊂ X looks like c-corrupted samples from a Gaussian of covariance Σ ≈ Σ, we expect that a random subset Z of Y also looks like c-corrupted samples from such a Gaussian. Moreover, as long as one uses m ≥ O(d) corrupted samples from a Gaussian, we can still (inefficiently) approximate the covariance. This motivates us to modify the robust algorithm as follows: rather than just checking whether Y looks like c-corrupted samples from a Gaussian of covariance roughly Σ, we also test whether an average subset Z ⊂ Y of size m does as well. Therefore, if Σ has low score, there exists a corresponding set Z ⊂ X of size m, and there are only n m ≤ e m log n choices for Z. So, now it suffices to have n ≥ m log n+d 2 ε•(c/k) , which will give us a bound of Õ(d 2 k/ε), as long as m ≤ O(d 2 ). Importantly, we still check the robust algorithm on Y of size roughly n/k, which allows us to keep the robustness threshold η at roughly c/k.There is one important caveat that for each subset Z, there is a distinct corresponding covariance Σ, and the volume of Σ ≈ Σ can change drastically as Σ changes. (For instance, the volume of Σ ≈ T • Σ is T Θ(d 2 ) times as large as the volume of Σ ≈ Σ. Since we have no bounds on the possible covariances, T could be unbounded.) For our volume ratio to actually be bounded by about n m • e O(d 2 ) , we want the volume of Σ ≈ Σ to stay invariant with respect to Σ. This can be done by defining a "normalized volume" where the normalization is inversely proportional to the determinant (see Appendix C.1 for more details). The robustness-to-privacy conversion (Theorem 2.1) will still hold.While the bound of Õ(d 2 k/ε) seems good, we recall that this bound is merely the sample complexity for (ε, δ)-DP crude approximation of a single Gaussian component. As discussed at the beginning of this subsection, to learn all k Gaussians, we actually need (ε/ k log(1/δ), δ/k)-DP, rather than (ε, δ)-DP, when crudely approximating a single component. This will still end up leading to a significant improvement over previous work [AAL24b] , but we can improve the volume ratio even further and thus obtain even better sample complexity.Improving Dimension Dependence: Previously, we used the fact that the volume of candidate Σ (corresponding to a fixed Z) was roughly exponential in d 2 for either V η/2 (X) or V η (X), so the ratio should also be roughly exponential in d 2 . Here, we improve this ratio, which will improve the overall volume ratio.First, we return to understanding the guarantees of the robust algorithm. It is known that, given m ≥ O(d) samples from a Gaussian of covariance Σ, we can provide an estimate Σ such that(1 -O( d/m))Σ ≼ Σ ≼ (1 -O( d/m))Σ.As above, we need to solve this even if a c fraction of samples are corrupted. While this can cause the relative spectral error to increase from 1 ± O( d/m) to 1 ± O(c + d/m), for now let us ignore the additional c factor. If V η/2 (X) > 0, then there is some covariance Σ and some set Y of size n/k, where the robust algorithm thinks Y looks like (possibly corrupted) Gaussian samples of covariance Σ. So, every Σ such that 0.5Σ ≼ Σ ≼ 2Σ has score of at most η/2 • n. This gives us a lower bound on V η/2 (X). We now want to upper bound V η (X). If S( Σ, X) < ηn, we still have that the robust algorithm thinks some Y looks like Gaussian samples of covariance Σ, where 0.5Σ ≼ Σ ≼ 2Σ. But now, we use the additional fact that for some Z ⊂ Y of size m, the robust algorithm on Z finds a covariance Σ. By the accuracy of the robust algorithm, Σ and Σ should be similar, i.e., (1-O( d/m))Σ ≼ Σ ≼ (1 -O( d/m))Σ (where we ignored the c factor). Thus, there exists some Z ⊂ X of size m and aΣ corresponding to Z, such that 0.5(1 -O( d/m)) Σ ≼ Σ ≼ 2(1 + O( d/m)) Σ.Therefore, from η/2 to η, we have dilated the candidate set of Σ by a factor of 1 + O( d/m) in the worst case, and we have at least 1 choice in the η/2 case but at most n m choices in the η case. Thus, the overall volume ratio is at mostn m • (1 + O( d/m)) d 2 = e O(m log n+d 2 • √ d/m) , since the dimension of Σ is roughly d 2 . Consequently, it now suffices to have n ≥ m log n+d 2 • √ d/m ε•(c/k) : setting m = d 5/3 gives us an improved bound of Õ(d 5/3 k/ε) for learning a single Σ i .There are some issues with this approach: most notably, we ignored the fact that the spectral error is really c + d/m rather than d/m. However, the robust algorithm can do better than just estimating up to spectral error c + d/m: it can also get an improved Frobenius error. While we will not formally state the guarantees on the robust algorithm here (see Theorem B.3 for the formal statement), the main high-level observation is that if the robust estimator Σ can be 1 ± c times as large as the true covariance Σ in only O(1) directions then for an average direction the ratio of Σ to Σ will be 1 ± O( d/m). We can utilize this observation to bound the volume ratio, using some careful ε-net arguments (this is executed in Appendix C.3). Our dimension dependence of d 5/3 will increase to d 7/4 , though this still improves over the previous d 2 bound.We will formally define the score function S( Σ, X) in Appendix E.1 and fully analyze the application of the robustness-to-privacy conversion, as outlined here, in Appendix E.2.Accuracy: One important final step that we have so far neglected is ensuring that any Σ of low score must be crudely close to some Σ i , if X is actually drawn from a GMM. We will just focus on the case where S( Σ, X) = 0, so some Y ⊂ X of size n/k looks like a set of samples from a Gaussian with covariance Σ. If the samples Y all came from the i-th component of the GMM, then it will not be difficult to show Σ is similar to Σ i . The more difficult case is when data point in Y are generated from several components. However, if n ≥ 20k 2 d, then |Y| ≥ 20kd, which means that by the Pigeonhole Principle, at least 20d points in Y come from the same mixture component (µ i , Σ i ). We are able to prove that, with high probability over samples drawn from a single Gaussian component N (µ i , Σ i ), that the empirical covariance of any subset of size at least 20d is crudely close to Σ i (see Corollary B.5). As a result, when verifying that a subset Y "looks like" i.i.d. Gaussian samples with covariance Σ, we can ensure that the empirical covariance of every subset Z ⊂ Y of size 20d is crudely close to Σ. Thus, if the score S( Σ, X) = 0, Σ is close to Σ, which is crudely close to some Σ i .We also formally analyze the accuracy in Appendix E.2.Putting everything together: To crudely learn a single Gaussian component with (ε, δ)-DP, we will need n ≥ Õ(d 7/4 k/ε) samples to find some covariance Σ with low score, and we also need n ≥ O(k 2 d) so that a covariance Σ of low score is actually a crude approximation to one of the real mixture components. To crudely approximate all components, we learn each Gaussian component with (ε/ k log(1/δ), δ/k)-DP. The advanced composition theorem will imply that repeating this procedure k times on the same data (to learn all k components) will be (ε, δ)-DP. Hence, by replacingε with ε/ k log(1/δ), we get that it suffices for n ≥ Õ d 7/4 k 3/2 √ log(1/δ) ε + k 2 d , if we need tocrudely learn all of the covariances. Finally, we can apply the private hypothesis selection technique (recall Section 2.1), which requires an additional Õ d 2 k α 2 + d 2 k αε . Combining these terms will give the final sample complexity.We remark that the sample complexity obtained above is actually better than the complexity in Theorem 1.4. There are two reasons for this. The first is that we have been assuming each component has weight 1/k, meaning it contributes about n/k data points. In reality, the weights may be arbitrary and thus some components may have much fewer data points. However, it turns out that one can actually ignore any component with less than α/k weight, if we want to solve density estimation up to total variation distance α. This will multiply the sample complexity termsÕ d 7/4 k 3/2 √ log(1/δ) ε + k 2 d ,needed for crude approximation, by a factor of 1/α. Finally, the informal Theorem 2.1 is slightly inaccurate, and the accurate version of the theorem will end up adding the additional term ofk 3/2 log 3/2 (1/δ) αε . Along with the final term Õ d 2 k α 2 + d 2 kαε from the private hypothesis selection, these terms will exactly match Theorem 1.4.In Appendix A, we note some general preliminary results. In Appendix B, we note some additional preliminaries on robust learning of a single Gaussian. In Appendix C, we discuss the robustnessto-privacy conversion and prove some volume arguments needed for Theorem 1.4. In Appendix D, we explain how to reduce to the crude approximation setting, using private hypothesis selection. In Appendix E, we design and analyze the algorithm for multivariate Gaussians, and prove Theorem 1.4. In Appendix F, we design and analyze the algorithm for univariate Gaussians, and prove Theorem 1.5. In Appendix G, we prove Theorem 1.6. Finally, Appendix H proves some auxiliary results that we state in Appendices B and C.Our results are on theoretical guarantees on the sample complexity of privately learning Mixtures of Gaussians. We do not provide any efficient or practical algorithms, and focus on statistical guarantees. We also do not discuss how to set the parameters and accuracy guarantees for practical applications, this is a question best left to practitioners. Finally, we assume each sample is i.i.d. drawn from a Gaussian Mixture Model distribution, though we remark that we use a "robustness-to-privacy" framework that will automatically make our algorithm robust to a roughly α/k fraction of corruptions.We state the advanced composition theorem of differential privacy. Theorem A.1 (Advanced Composition Theorem [DR14, Theorem 3.20]). Let ε, δ, δ ′ ≥ 0 be arbitrary parameters. Let A 1 , . . . , A k be algorithms on a dataset X, where eachA i is (ε, δ)-DP. Then, the concatenation A(X) = (A 1 (X), . . . , A k (X)) is 2k log 1 δ ′ • ε + kε • (e ε -1), k • δ + δ ′ -DP.Moreover, this holds even if the algorithms A i are adaptive. By this, we mean that for all i ≥ 2 the algorithm A i is allowed to depend on A 1 (X), . . . , A i-1 (X). However, privacy must still hold for A i , conditioned on the previous outputs A 1 (X), . . . , A i-1 (X).Next, we note the properties of the Truncated Laplace distribution and mechanism. Definition A.2 (Truncated Laplace Distribution). For ∆, ε, δ > 0, the Truncated Laplace Distribution TLap(∆, ε, δ) is the distribution with PDF proportional to e -|x|•ε/∆ on the region [-A, A], whereEQUATION, and PDF 0 outside the region[-A, A]. Lemma A.3 (Truncated Laplace Mechanism [GDGK20, Theorem 1]). Let f : X n → R be a real- valued function, and let ∆ > 0, such that for any neighboring datasets X, X ′ , |f (X) -f (X ′ )| ≤ ∆. Then, the mechanism A that outputs f (X) + TLap(∆, ε, δ) is (ε, δ)-DP.In this section, we note some standard but useful lemmas. We first note the Courant-Fischer theorem. Theorem A.4 (Courant-Fischer). Let A ∈ R d×d be a real symmetric matrix, with eigenvaluesλ 1 ≥ λ 2 ≥ • • • ≥ λ d . Then, for each k, λ k = max V ⊂R d dim(V )=k min x∈V ∥x∥2=1 x ⊤ Ax = min V ⊂R d dim(V )=d-k+1 max x∈V ∥x∥2=1 x ⊤ Ax,where V refers to a linear subspace of R d .By setting A = J ⊤ J, we have the following corollary. Corollary A.5. Let J ∈ R d×d be a real (possibly asymmetric) matrix with singular valuesσ 1 ≥ σ 2 ≥ • • • ≥ σ d . Then, for each k, σ k = max V ⊂R d dim(V )=k min x∈V ∥x∥2=1 ∥Jx∥ 2 = min V ⊂R d dim(V )=d-k+1 max x∈V ∥x∥2=1 ∥Jx∥ 2 .Next, we note a basic proposition. Proposition A.6. [HKMN23, Lemma 6.8] Let M ∈ R d×d be a real symmetric matrix, and J ∈ R d×d be any real-valued (but not necessarily symmetric) matrix such that ∥JJ ⊤ -I∥ op ≤ ϕ ≤ 1, for some parameter ϕ. Then,∥J ⊤ M J∥ 2 F ≤ (1 + 3ϕ) • ∥M ∥ 2 F .We also note the Hanson-Wright inequality. Lemma A.7 (Hanson-Wright). Given a d × d matrix M ∈ R d×d and a d-dimensional Gaussian vector X ∼ N (0, I), for any t ≥ 0,P X ⊤ M X -E[X ⊤ M X] ≥ t ≤ 2 exp -c min t 2 ∥M ∥ 2 F , t ∥M ∥ op ,for some universal constant c > 0.Finally, we note a folklore simple characterization of total variation distance (see also [AAL23, Theorem 1.8]). Lemma A.8. For anyµ 1 , µ 2 ∈ R d and positive definite Σ 1 , Σ 2 ∈ R d×d , d TV (N (µ 1 , Σ 1 ), N (µ 2 , Σ 2 )) ≤ 1 √ 2 • max ∥Σ -1/2 1 Σ 2 Σ -1/2 1 -I∥ F , ∥Σ -1/2 1 (µ 1 -µ 2 )∥ 2 .In this section, we establish what a robust algorithm can do for learning a single high-dimensional Gaussian. We will state the accuracy guarantees in terms of the following definition, which combines Mahalanobis distance, spectral distance, and Frobenius distance.Definition B.1. Given mean-covariance pairs (µ, Σ) and (μ, Σ), and given parameters γ, ρ, τ > 0, we say that (μ, Σ)≈ γ,ρ,τ (µ, Σ) if • ∥Σ -1/2 ΣΣ -1/2 -I∥ op ≤ γ. (Spectral distance) • ∥Σ -1/2 ΣΣ -1/2 -I∥ F ≤ ρ. (Frobenius distance) • ∥Σ -1/2 (μ -µ)∥ 2 ≤ τ . (Mahalanobis distance)We will also use ≈ γ,τ as a shorthand for ≈ γ,√ d•γ,τ , i.e., there is no additional Frobenius norm condition beyond what is already imposed by the operator norm condition.Although ≈ γ,ρ,τ is neither symmetric nor transitive, symmetry and transitivity happen in an approximate sense. Namely, the following result holds: we defer the proof to Appendix H. We remark that similar results are known (e.g., [AL22, Lemma 3.2], [AAL23, Lemma 4.1]). Proposition B.2 (Approximate Symmetry and Transitivity). Fix any γ, ρ, τ > 0 such that γ ≤ 0.1. Then, for any (µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ), we have that (µ 2 , Σ 2 ) ≈ 2γ,2ρ,2τ (µ 1 , Σ 1 ), and for any(µ 1 , Σ 1 ) ≈ γ,ρ,τ (µ 2 , Σ 2 ) and (µ 2 , Σ 2 ) ≈ γ,ρ,τ (µ 3 , Σ 3 ), we have that (µ 1 , Σ 1 ) ≈ 4γ,4ρ,4τ (µ 3 , Σ 3 ).We note the following theorem about the accuracy of robustly learning a single Gaussian. While this result will roughly follow from known results and techniques in the literature, we were unable to find a formal statement of the theorem, so we prove it in Appendix H. Theorem B.3. For some universal constants c 0 ∈ (0, 0.01) and C 0 > 1, the following holds. Fix anyη ≤ γ ≤ c 0 , β < 1, and ρ such that O(η) ≤ ρ ≤ c 0 √ d. Let n ≥ O d+log(1/β) γ 2 + (d+log(1/β)) 2 ρ 2. Then, there exists an (inefficient, deterministic) algorithm A 0 with the following property. For any µ, Σ, with probability at least 1 -β over X = {X 1 , . . . , X n } ∼ N (µ, Σ), and for all datasetsX ′ with d H (X, X ′ ) ≤ η • n, we have A 0 (X ′ ) ≈ C0γ,C0ρ,C0γ (µ, Σ).We remark that A 0 may have knowledge of η, γ, ρ, β, but does not have knowledge of µ, Σ, or the uncorrupted data X.In this section, we note that, when given i.i.d. samples from a Gaussian, with high probability one cannot choose a reasonably large subset with an extremely different empirical mean or covariance.First, we note the following proposition, for which the proof follows by an ε-net type argument. Proposition B.4. Let n ≥ n ′ ≥ 20d, and L ≥ C 1 n 2 for a sufficiently large constant C 1 . Suppose that some data points X = {X 1 , . . . , X n } are i.i.d. sampled from N (0, I). Then, with probability at least 1 -n -Ω(n ′ ) , the following both hold. 
 
 Appendix:  


Decision: The paper considers the problem of estimating k-component and d-dimensional Gaussian mixtures with central differential privacy. Previous algorithms required k^2 d^4 samples and the proposed work shows that kd^2 + k^1.5 d^1.75 + k^2 d samples suffice. In the high-dimensional setting, when d >>k^2, they further show that this bound is optimal. They also show that sample complexity of learning one dimensional Gaussian mixtures is linear in k. The results are interesting and I recommend acceptance.

As reviewers mention, the writing is dense and based on my own read of the paper, I completely agree.  I strongly urge authors to utilize the full page limit and improve the writing in the final version of the paper by addressing reviewer concerns.  It would also be great if authors can add  additional intuition on how each of the terms in the sample complexity arise and also highlight the technical novelty that enabled this result.

Reviews: review: paper_summary: paper_summary: The paper studies the problem of learning a mixture of $k$ $d$-dimensional Gaussians using a differentially private mechanism with respect to the samples. It provides an improved sample complexity which has asymptotically optimal dependence on the dimension $d$ for small $k$, the lower bound is also given by the paper. Additionally, for $d=1$ the paper provides optimal sample complexity.  
The paper follows a high level approach of obtaining crude approximations of the mean and covariances of the component Gaussians, which can be used to reduce – using a net based argument – the problem to private hypothesis selection for which existing algorithm (needing only log the size of the net number of samples) can be used. To obtain the crude approximations respecting the DP guarantee, the paper uses a variant of the exponential mechanism with a carefully constructed scoring function measuring the distance between the a candidate Gaussian and any “heavy” component of Gaussian mixture. This uses a number of samples depending on the size of the sample set $n$ and the dimensionality the hypotheses $d^2$. A naive approach fails due to the blowup incurred, and the paper utilizes sample compression to reduce the dependence to $O(d\\log n)$ and improves the dimensionality dependence to $O(d)$ via an estimation argument.  
The above is outlined in more detail in the main paper and proved formally in the appendices. While all the details have not been verified by the reviewer, the approach taken by the paper seems correct.
strength_weakness: strength_weakness: 1. The sample complexity does not match the lower bound in the dependence on $k$.
2. The degradation on the dependence on $k$ from the univariate to the multi-variate case is not explained in the main paper.
3. The details of the various technical parts are tedious and could be alleviated by better presentation, e.g. by listing all the parameters and their dependencies in a table.
questions: questions: It will be useful to combine definitions 1.1 and 1.2 into a precise definition of privately learning GMMs.
limitations: limitations: Yes, addressed.
paper_summary: paper_summary: This works focuses on the task of density estimation of a mixture of Gaussians under the restriction of differential privacy. Unlike parameter estimation, density estimation does not require accurately estimating the mixture's parameters, but instead bounds the total variation distance between estimated and true distributions. This task can be achieved even without any boundedness or separation assumptions on the parameters of the components.

This task is known to require sample size of $O(k d^{2}/\alpha^{2})$ even in the non-private setting, where $k$ is the number of components, $d$ is the dimension, $\alpha$ is the bound on the TV distance, and the $O$ notation ignores logarithmic factors. In the private setting, previous work [1] achieved $\alpha$ accuracy guarantee under $(\epsilon, \delta)$-DP with $O(k^{2} d^{4}/\alpha^{2} \epsilon)$ (the exact bound includes several other terms that depend on $\log(1/\delta)$ as well).

This work provides an improved bound, reducing the dependence on the parameters to $O \left(\frac{k d^{2}}{\alpha^{2}} + \frac{k d^{2}}{\alpha \epsilon} \right)$ in the high dimensional regime, $O \left(\frac{k}{\alpha^{2}} + \frac{k \cdot \log(1/\delta)}{\alpha \epsilon} \right)$ in the univariate setting, and lower bounds which asymptotically match the upper bound in the univariate case and nearly match it in the multivariate one.

The proposed algorithm relays on first achieving crude estimation of the parameters and then using hypothesis selection to improve the estimation. The crude estimation is achieved using robust GMM estimation techniques and robusteness-to-privacy conversion, based on an inverse sensitivity-like instantiation of the exponential mechanism. This method is computationally inefficient.

[1] Mohammad Afzali, Hassan Ashtiani, and Christopher Liaw. Mixtures of gaussians are privately learnable with a polynomial number of samples. In Algorithmic Learning Theory, pages 1–27. PMLR, 2024.
strength_weakness: strength_weakness: Despite the great work done by the authors and my best efforts, I was not able to follow up all the proof structure. In particular, I could not find the justification for some of terms that appeared in the bound presented on Theorem 1.3. It seems to me like a section that "puts everything together" will be of great benefit. I will try to describe my current understanding and existing gaps.

To the best of my understanding:
* The $\frac{k d^{2}}{\alpha^{2}}$ and $\frac{k d^{2}}{\alpha \epsilon}$ were explained at the "Reducing to “crude approximation”" section, and they represent the sample size required to accurately and privately learn the GMM given some crude estimation of its components, using hypothesis selection method.
* The $\frac{d^{1.75} k^{1.5} \sqrt{\log(1/\delta)}}{\alpha \epsilon}$ was explained at the "Improving Dimension Dependence" section, where $\frac{d^{1.75} k}{\alpha \epsilon}$ results from the fact $O\left(\frac{d^{1.75} k}{c \epsilon} \right)$ points are required to get the crude estimation of the parameters of a single component under $c$-robustness, which then can be translated to a DP estimation with $\alpha$ accuracy using rubousteness-to-privacy conversion techniques (Theorem 2.1), and the additional $\sqrt{k \cdot \log(1/\delta)}$ term results from advanced composition over $k$ components.
* I failed to understand where the additional two terms ($\frac{(k \cdot \log(1/\delta))^{1.5}}{\alpha \epsilon}$ and $\frac{k^{2} d}{\alpha}$) come from, and I suspect they were accumulated at some point during the rubousteness-to-privacy transformation.
questions: questions: As I mentioned before, I will find an additional section that brings all the components together to outline the final proof method, focusing on the final achieved bound, to be very useful.
limitations: limitations: N/A
paper_summary: paper_summary: This paper investigates the sample complexity of privately learning mixtures of Gaussians. The authors achieve a sample complexity of approximately $O(kd^2+k^{1.5}d^{1.75}+k^2d)\log R$ where $R$ is an upper bound on the condition number of the covariance matrix and the norm of the mean. This result improves upon the previous best bound of $O(k^2d^4)$. Additionally, they constructed a lower bound of $\Omega(kd^2+k\log(1/\delta)/\epsilon)$, which refutes a conjecture from prior work, and they achieve optimal sample complexity when $d=1$.
strength_weakness: strength_weakness: The paper is technically dense, making it challenging to understand and verify all the details.
The authors do not fully utilize the page limit, which could have been used to explain the techniques more clearly.
The technique for establishing the lower bound is not discussed in the main body of the paper.
Including formal definitions and crucial lemmas in the main text could help readers better understand the key ideas. For instance, the definition of volume was unclear until I consulted the appendix.
questions: questions: 1. Where does the $k^2$ term come from in the upper bound?
2. Are there alternative methods to address this problem that do not rely on the robustness-to-privacy conversion technique?
limitations: limitations: Yes
paper_summary: paper_summary: The paper studies the important problem of private learning of the Gaussian mixture model to estimate the underlying distribution within a desired total variation distance. By combining different techniques, the authors succeed in deriving bounds that are of quadratic dimension, thus significantly improving the existing results.
strength_weakness: strength_weakness: - The paper could benefit from some numerical verifications of the results/algorithms used.

- The appendix section is poorly organized. Without an outline and proof, it is difficult to follow such a long appendix.
questions: questions: 1. It would be better to explain more the difference between density estimation and parameter estimation.

2. How do the bounds of Theorems 1.3, 1.4., and 1.5 behave with respect to the failure probability $\beta$ order-wise?

3. Can the authors discuss the assumption $R=n^{100}$ (line 160) a bit more? Why does it depend on $n$ (and not $k$ or $d$?)? Also, is the polynomial dependence restrictive?

4. The first sentence of the informal statement of Theorem 2.1. does not seem rigorous enough. Is such a $\tilde{\Sigma}$ unique? I don't think so. Then, if not, which choice of $\tilde{\Sigma}$ is taken into account when computing $V_{\eta}(\mathbf{X})$? Probably maximum volume?

5. The appendix sections are very difficult to follow. There should be a detailed outline at the beginning to guide the reader as to what each appendix deals with. Also, there should be sufficient references in the main text. For example, after the informal statement of Theorem 2.1., it should be clearly stated where the full statement and proof can be found. Similarly, after each paragraph (e.g., Section 2.1) or after each mentioned previously established result (e.g., advanced composition theorem), there should be an exact pointer to where the complete version can be found in the appendices.

6. Since the submission was allowed for 9 pages, the authors had 1.5 more pages. I think some material from the appendices, for example the algorithms used to privately learn GMM, could be moved to the main text.

7. Finally, what is missing the most is the experimental section. If I'm not mistaken, this should be easy to do. Although not all theoretical work requires experimental verification, I believe that if the experiments verify the theoretical finding, this would greatly increase the importance of the results.

Minor comments:

- Line 52: it is better to add "The parameters $\alpha$ and $\beta$...".

- The constant $c^*$ in Theorem 1.5. does not appear in the bound. 

- In the first paragraph of Section 2.1, it is mentioned that "... as we can finish the procedure with hypothesis selection, as discussed above". However, hypothesis selection has not been sufficiently discussed.

- Is it true that the intuition given in Section 2.1. holds if $ n \gg k^2$?

- Typo line 178 (a a)

- Is $k$ in lines 203-204 the number of components in the Gaussian mixture? Why does this statement ("if we altered $k$ data points...") only hold for $k$ changes?
limitations: limitations: NA

Comments: comment: Thanks for the response. I maintain my score.
I acknowledge the rebuttal by the authors. Just to clarify, in my weakness no. 3 comment, I had suggested listing the parameters used in the different proofs and their dependencies in table(s). This would make the proofs a bit more accessible in my opinion and I hope the authors will look into it.
Overall my rating of Accept remains unchanged.
I thank the authors for their explanation, and hope it will be reflected in the final version, so that all readers will have the opportunity to fully comprehend this important result.
I thank the authors for the clarifications. I believe the paper merits the publication, and I maintain my initial rating.