Training Gaussian Mixture Models at Scale via Coresets

1. Introduction

• Mixture-model inference over latent variables depends on the full data set, which may be distributed across machines or arrive as a memory-limited stream.The paper focuses on mixtures of Gaussians as a statistical mixture-model case.Source: §1 · p. 1
• A coreset is a weighted subset whose fitted models also provide a good fit to the original data, with size independent of data-set size.The construction targets mixtures whose component covariance eigenvalues lie in [λ, 1/λ].Source: §1 · p. 2
• Coresets approximate GMM log-likelihoods to arbitrary accuracy as ε → 0 and can be constructed in parallel or in streaming settings.The stated space and per-point update-time bounds are polynomial in d, k, λ^-1, ε^-1, log n, and log(1/δ).Source: §1 · p. 2

2.1 Fitting Gaussian Mixture Models by Maximum Likelihood Estimation

• GMM fitting seeks maximum-likelihood parameters for a mixture whose component weights sum to one, with covariance eigenvalues constrained to [λ, 1/λ] to avoid degenerate solutions.The section frames this as inference over mixture components and their means and covariance matrices.Source: §2.1 · p. 2

2.2 Approximating the Log-likelihood

• The scale dependence of log-likelihood prevents a general multiplicative-likelihood guarantee, motivating approximation of log-likelihood values over all constrained models.The section decomposes negative log-likelihood into a data-dependent cost and a data-independent normalizing term, then targets the cost with a weighted set.Source: §2.2 · p. 3
• Uniform subsampling can miss a small Gaussian component, making the fitted model's difference from the full-data MLE arbitrarily large as component means separate.In the example, the small component has weight 1/√n; the figure notes that m ∈ Θ(√n) is needed for the uniform sample to include it with high probability.Source: Figure 1 · p. 4

2.3 Coresets for Semi-spherical Gaussian Mixtures

• A (k, ε)-coreset is intended to approximate the data-dependent cost uniformly across k-Gaussian mixtures, which in turn uniformly approximates log-likelihood.The cost for a weighted set uses each point's weight in its contribution.Source: §2.3 · p. 4
• Because a coreset is much smaller than the original data, weighted EM on it can be orders of magnitude faster while its optimal model approximates the full-data optimum in log-likelihood.The passage notes that finding the optimal coreset model is difficult, so weighted EM can be used to seek a good solution.Source: §2.3 · p. 4

3.1 Naive Approach via Uniform Sampling

• Uniform sampling can miss a low-weight Gaussian component unless Ω(√n) points are sampled, producing arbitrarily poor log-likelihood approximations even for two well-separated Gaussians.The problem becomes more pronounced as k and d grow.Source: §3.1 · p. 5

3.2 Better Approximation via Importance Sampling

• Approximate k-means solutions suffice to construct a sampling strategy that biases samples toward less dense regions and yields valid coresets.The k-means approximation partitions the data into Voronoi cells used to estimate density and variance.Source: §3.2 · pp. 5–6
• Theorem 2 guarantees that the constructed set is a (k, ε)-coreset with probability at least 1 − δ, and its size is independent of n and polynomial in k, d, and ε.The theorem’s construction selects the best of p independent bicriteria runs by φ(X, B_i).Source: §3.2 · p. 6
• The construction uses a bicriteria k-means approximation, computes point sensitivities, and samples weighted points according to q(x).The suggested k-means++ solution has size k and an expected O(log^2 k) approximation; the construction steps are parallelizable.Source: §3.2 · pp. 6–7

• Weighted EM accounts for coreset point weights when fitting a Gaussian mixture, whereas exhaustive search over coreset partitions is infeasible for soft-assigned mixture components.The paper uses weighted EM in its experiments.Source: §4 · p. 7
• Because EM runs on a significantly smaller data set, it can use multiple random restarts; weighted EM on the coreset typically achieves comparable log-likelihood performance to full-data EM.The experiments constrain variances to be greater than a small preselected λ to avoid zero-variance estimates with infinite density and log-likelihood.Source: §4 · p. 7
• The weighted EM procedure initializes with weighted k-means, then alternates expectation and maximization steps using the point weights and variance threshold λ.The maximization step updates component weights, means, and covariances.Source: §4 · p. 8

5. Streaming and Parallel Computation

• Coreset composition allows coresets for separate data subsets to be combined and then compressed, supporting both streaming and parallel construction.The composition properties build on earlier work on geometric clustering coresets.Source: §5 · p. 8

5.1 Streaming Computation

• A binary merge-and-compress tree keeps one coreset per level in memory, avoiding the exponential error growth possible with immediate repeated compression.For the illustrated stream, intermediate errors grow from ε to 4ε and then 13ε; the general procedure sets intermediate error to control the final bound.Source: §5.1 · p. 9
• A stream of n points admits a λ-semi-spherical GMM coreset with per-point update time and memory polynomial in d, k, λ^-1, ε^-1, log n, and log(1/δ), with probability at least 1 − δ.The construction uses exponentially increasing batches when n is unknown in advance.Source: §5.1 · p. 9

5.2 Distributed Computation

• Distributed construction partitions the data, computes partition coresets, and performs O(log n) parallel merge-and-compress iterations.The process maps naturally to map-reduce, with map tasks constructing disjoint-partition coresets and reduce tasks merging them.Source: §5.2 · p. 10
• With m machines, a λ-semi-spherical GMM coreset can be computed in time (n/m) · poly(d, k, λ^-1, ε^-1, log(1/δ), log n), with probability at least 1 − δ.The simpler union-based approach can first produce a coreset whose size depends on m, then compress that union to remove the machine-count dependence.Source: §5.2 · p. 10

6. Experimental Evaluation

• The evaluation compared coreset sizes from 1’000 to 10’000, using 80% of each dataset for training and the remaining 20% to compute relative error.Subsampling results are medians over 200 runs; the full-data baseline is the best of 100 runs.Source: §6 · pp. 10–11
• At fixed subsample sizes, coresets had smaller approximation errors than uniform subsamples and achieved significant speedups over full-data training.Figure 4 reports median relative error against the best models trained on the full data set.Source: Figure 4 · p. 11, §6 · p. 11
• For Higgs at m = 10’000, coresets achieved 57.5× speedup and 4% relative error, versus 12.2% error for uniform subsampling.The dataset contains 11’000’000 instances; evaluation used 150-component GMMs on the first two principal components.Source: §6 · p. 11
• For csn at m = 10’000, coresets achieved 7× speedup and 6.6% relative error, versus 58.7% error for uniform subsampling.The GMM had 100 components.Source: §6 · p. 11
• For Kdd at m = 10’000, coresets achieved 8.3× speedup and 2.7% relative error, versus 16.8% error for uniform subsampling.The GMM had 10 components; uniform subsampling may be more sensitive to initialization.Source: §6 · p. 11
• For msyp at m = 5’000, coresets achieved 78× speedup and 2.3% relative error, versus 3.9% error for uniform subsampling.The GMM had 50 components and used the top 25 principal components.Source: §6 · p. 11

7. Related Work

• This extended version improves on prior work with direct negative log-likelihood analysis, a new importance-sampling strategy, a practical linear-time algorithm, and substantially larger empirical datasets.It also reports constructing larger coresets in less time, by approximately two orders of magnitude.Source: §7 · p. 12

7.1 Learning Mixtures of Gaussians

• Earlier Gaussian-mixture learning results variously assume separation, axis-aligned components, non-degeneracy, or identical spherical Gaussians; the cited algorithms do not apply to streaming or parallel settings.The paper contrasts these constraints with its milder component assumptions and streaming and parallel construction.Source: §7.1 · p. 12

7.2 Approximation Algorithms via Coresets

• Coresets have been applied to numerous geometric, graph, streaming, and learning problems, and this paper uses the sensitivity framework of Feldman and Langberg (2011).The related work also notes a reduction of quadratic dependence on total sensitivity to near-linear via Braverman et al. (2016).Source: §7.2 · pp. 12–13

8. Conclusion

• The construction uses a computational-geometry connection and importance sampling, runs linearly in n, and yields coresets whose size is independent of dataset size.The guarantees apply to possibly unstructured data without assumptions on its generating distribution or model.Source: §8 · p. 13
• On several real-world datasets, training time fell by up to two orders of magnitude while hold-out likelihood remained competitive with full-data models.Coresets can be constructed in parallel or in one streaming pass using poly(d, k, λ−1, ε−1, log n, log(1/δ)) space and update time.Source: §8 · p. 13
• Open problems are whether coreset size can be independent of λ for any set of k Gaussian mixtures and independent of d.Source: §8 · p. 13

Appendix A. Bounding the Coreset Size

• Sensitivity is the worst-case contribution of a point to log-likelihood, and a looser uniform upper bound increases the resulting coreset size.For Gaussian mixtures, the query space comprises all k-component GMMs in d dimensions.Source: Appendix A · p. 14

A.1 Total Sensitivity by Reduction to the Euclidean Space

• The proof bounds log-likelihood contribution differences for points using a Voronoi partition induced by the bicriteria set, then applies the bounds uniformly over mixture models.The analysis uses covariance eigendecomposition and the geometric distance to the corresponding Voronoi representative.Source: §A.1 · pp. 14, 16, 17
• A rough approximation to optimal k-means clustering lower-bounds the denominator needed to upper-bound each point’s sensitivity.Source: §A.1 · p. 15
• The resulting total sensitivity is independent of |X|, and its bound can be computed from a bicriteria k-means approximation available in linear time.Using k-Means++ gives α = O(log2 k) and β = k; the bound is tight up to a constant because some datasets have S ∈ Θ(k).Source: §A.1 · p. 17

A.2 Pseudo-dimension of Gaussian Mixtures

• Pseudo-dimension is analyzed by converting real-valued functions into a subgraph class whose VC dimension is equal to the original family’s pseudo-dimension.The paper defines the subgraph indicator as sgn(f(x) − y).Source: §A.2 · pp. 17–18
• The threshold function hθ,r(x) is computed from quadratic-form exponents whose sum is compared with zero.Its computation uses O(m) arithmetic operations and k + 1 exponential evaluations.Source: §A.2 · p. 18, Figure 5 · p. 19
• The pseudo-dimension bound for Gaussian-mixture functions is O(m^2k^2 + mk(t + log mk)).Here m = k(d + 1)(d + 2)/2 − 1, and t counts the operations in the threshold-function computation.Source: §A.2 · p. 18
• Whether the gap between the upper bound and the Ω(kd^2) lower bound can be closed remains open in the general setting.The lower bound is attributed to Akama and Irie (2011).Source: §A.2 · p. 18

A.3 Sufficient Coreset Size

• A bound on coreset size can also be obtained using Feldman and Langberg (2011), by bounding primal shattering dimension instead of pseudo-dimension.The passage refers readers to Bachem et al. (2017a) for discussion of the difference.Source: §A.3 · p. 19
• Theorem 2 establishes uniform approximation of the coreset and full-data log-likelihoods as ε → 0 under its stated assumptions.Its construction uses independent runs of Algorithm 2 and selects B⋆ by minimizing φ(X, Bi).Source: §A.3 · p. 19
• With probability at least 1 − δ, the sampled set C is a (k, ε)-coreset of X.The stated sample-size bound uses total sensitivity and the function-family dimension.Source: §A.3 · pp. 19–20

A.4 Directly Approximating the Log-Likelihood

• Under the condition that each component’s covariance eigenvalue product is at least 1/(2π)^d, the coreset and full-data costs satisfy multiplicative (1 ± ε) bounds.The bounds yield corresponding inequalities for log-likelihood using the stated decomposition.Source: §A.4 · pp. 20–21

Appendix B. Convergence of Weighted EM for Gaussian Mixtures

• Weighted EM differs from the non-weighted version by scaling updates proportionally to each point’s weight, and it converges to a stationary point.The convergence statement cites Dempster et al. (1977).Source: Appendix B · p. 21