← All tutorials

Gaussian Mixture Models (GMM)

Soft clustering with overlapping bell curves.

On this page

k-Means Clustering gives every point one cluster and every cluster a circle. A Gaussian mixture gives each point a share in every cluster and lets clusters be tilted ellipses.

1 Hard circles or soft ellipses

A Gaussian mixture model (GMM) describes the data as a weighted sum of Gaussian bumps, gives each point a share in every bump, fits the bumps with the EM loop from Maximum Likelihood and EM, and picks the number of bumps with a score that charges for extra parameters.

2 A mixture adds weighted bumps

Each Gaussian, a component, has a center, a covariance and a weight, and EM picks the values that maximize the likelihood of the data. EM runs once over the training set; after that, scoring a new point is one pass over the K components.

3 Each point gets a share, by Bayes rule

Which component made a point? Bayes rule answers with shares: each component’s height at the point, times its weight, rescaled to add to 1.

4 Fitting the mixture with EM

The shares need the components and the components need the shares, the loop EM solves. The E step computes every r_ik; the M step refits each component from the points counted at their shares.

This snippet runs one EM round in NumPy on two toy blobs and prints the new weights and means; after the starting guess, its lines follow the card from top to bottom.

import numpy as np
from scipy.stats import multivariate_normal as mvn
rng = np.random.default_rng(0)
blob_a = rng.normal([-1, -1], 0.5, (100, 2))
blob_b = rng.normal([1, 1], 0.5, (150, 2))
X = np.vstack([blob_a, blob_b])
w = np.array([0.5, 0.5])
mu = np.array([[-0.5, 0.0], [0.5, 0.0]])
Sigma = np.array([np.eye(2), np.eye(2)])
K = len(w)
components = [mvn(mu[k], Sigma[k]) for k in range(K)]
a = np.column_stack([w[k] * components[k].pdf(X) for k in range(K)])
r = a / a.sum(axis=1, keepdims=True)
N_k = r.sum(axis=0)
mu = r.T @ X / N_k[:, None]
covs = [np.cov(X.T, aweights=r[:, k], bias=True) for k in range(K)]
Sigma = np.stack(covs)
w = N_k / len(X)
print(w.round(2), mu.round(2))

1. 1Draws two toy clumps of 100 and 150 points.
2. 2Sets starting weights w, means μ and covariances Σ for K = 2.
3. 3Computes the scores a_ik, eq card in Section 3, line 1, then the E-step shares r_ik, eq card in Section 4, line 1.
4. 4Moves each mean μ_k, eq card in Section 4, line 2.
5. 5Refits each Σ_k, eq card in Section 4, lines 3 and 4.
6. 6Resets the weights w_k, eq card in Section 4, line 5.

5 k-means is the hard-label limit

Give every component the same round covariance, with its standard deviation shrinking toward zero. Each share becomes 0 or 1 for the nearest center, and the center update becomes k-means’. k-means never estimates a covariance, hence its circles.

6 Covariance shapes trade flexibility for parameters

A full covariance lets every component tilt and stretch freely, at three numbers per component in 2D. Cheaper options constrain it; scikit-learn names them with covariance_type.

The long axes of each ellipse are the directions PCA (principal component analysis) would find inside that cluster.

7 How many components? Let BIC decide

More components always raise the likelihood, so the likelihood alone cannot pick K. The Bayesian information criterion (BIC) can.

This snippet fits 1 to 6 components to the same two blobs and prints each BIC; look for the smallest.

import numpy as np
from sklearn.mixture import GaussianMixture
rng = np.random.default_rng(0)
blob_a = rng.normal([-1, -1], 0.5, (100, 2))
blob_b = rng.normal([1, 1], 0.5, (150, 2))
X = np.vstack([blob_a, blob_b])
for K in range(1, 7):
    gm = GaussianMixture(n_components=K, random_state=0).fit(X)
    print(K, gm.bic(X))  # lowest wins

With two separate blobs, expect the smallest BIC at 2.

8 Gaussian mixtures in practice

In scikit-learn a mixture sits after a scaler; fit runs EM once on the training points, and predict_proba or score reads new points.

import numpy as np
from sklearn.mixture import GaussianMixture
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler
rng = np.random.default_rng(0)
short_eruptions = rng.normal([2.0, 55], [0.3, 6], (100, 2))
long_eruptions = rng.normal([4.3, 80], [0.4, 6], (172, 2))
X = np.vstack([short_eruptions, long_eruptions])
gm = GaussianMixture(n_components=2, covariance_type="full", random_state=0)
model = make_pipeline(StandardScaler(), gm)
model.fit(X)
new = [[2.1, 54.0], [4.0, 79.0]]
print(model.predict_proba(new).round(2))  # shares r_ik
print(model.score(new), gm.n_iter_)  # mean log-likelihood, rounds

1. 1Draws a geyser-like toy cloud, eruption length and waiting time both in minutes, 100 short and 172 long eruptions.
2. 2Chains StandardScaler with a two-component GaussianMixture.
3. 3Runs EM rounds until the average log-likelihood per point gains less than tol.
4. 4Gives new points their shares r_ik, eq card in Section 3, line 2.
5. 5Prints the mean log-likelihood and the rounds run.

Defaults are from the GaussianMixture docs. Gotchas are in Section 10.

9 Where GMMs are used

In speaker verification, the classic use, Reynolds et al. (2000) fit a large background GMM to speech from many people and adapt each speaker’s model from it. Each short frame becomes a vector of MFCCs, numbers describing its spectrum shape (see MFCC and Other Audio Features), and both models score the frames.

A fitted mixture scores every new point (score_samples); a very low score flags an anomaly. Classic speech recognizers used a mixture per hidden state to score each frame (Rabiner 1989, eq. 49), as in Hidden Markov Models.

10 Common traps

• A component can collapse onto one point and send the likelihood to infinity; reg_covar prevents it.
• EM finds a local peak; raise n_init to try several starts.
• Unscaled features bend every ellipse; standardize first, shifting each feature to mean 0 and scaling it to standard deviation 1.
• Component names are arbitrary, so K components give K! (K factorial) equal answers (label switching).
• Too many components overfit; let BIC choose.

11 What if the hidden cluster changes over time?

When the hidden cluster follows a chain over time, as words in speech or a casino’s dice do, the model is a hidden Markov model.

12 References

• Bishop, C. M. (2006). Pattern Recognition and Machine Learning. Springer. Ch. 9 (mixtures and EM) and App. A (Old Faithful). PDF.
• Dempster, A. P., et al. (1977). Maximum likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society B, 39(1), 1 to 38.
• Reynolds, D. A., et al. (2000). Speaker verification using adapted Gaussian mixture models. Digital Signal Processing, 10, 19 to 41.
• Rabiner, L. R. (1989). A tutorial on hidden Markov models and selected applications in speech recognition. Proceedings of the IEEE, 77(2), 257 to 286.
• scikit-learn, sklearn.mixture.GaussianMixture documentation.