Home / MLP / Week 12

Week 12 · Unsupervised learning in sklearn

The one-minute story

By the end of this week you can:

• run KMeans, read cluster_centers_, labels_, inertia_, and choose $k$ with the elbow and silhouette;
• explain when DBSCAN beats k-means and what eps, min_samples and the label $-1$ mean;
• fit GaussianMixture, get soft assignments with predict_proba, and choose components with bic;
• apply PCA: n_components, explained_variance_ratio_, inverse_transform, and choosing the number of components.

Words you need

Idea 1: KMeans

from sklearn.cluster import KMeans
km = KMeans(n_clusters=3, n_init=10, random_state=0).fit(X)     # X: scaled features
km.labels_            # the cluster of each training point, shape (n,)
km.cluster_centers_   # the centroids, shape (3, n_features)
km.inertia_           # sum of squared distances to the nearest centre
km.predict(new_X)     # nearest-centre label for new points
km.transform(X)       # distances to every centre, shape (n, 3)

Defaults worth knowing: n_clusters=8, init="k-means++" (smart seeding, MLT week 3), max_iter=300, tol=1e-4. n_init (the number of restarts, keeping the best) was 10 in older versions and is "auto" in recent ones; set it explicitly to avoid surprises. random_state makes the result repeatable. Scale your features first: k-means uses distances. KMeans has no partial_fit; MiniBatchKMeans does (for huge data).

fit_predict(X) fits and returns the labels in one call. inertia_ equals $\sum_i\|\mathbf x_i-\boldsymbol\mu_{z_i}\|^2$ (we checked it by hand).

Choosing $k$

We made 300 points from 3 blobs. The inertia for $k=1,\dots,6$:

It collapses until $k=3$, then barely moves: the elbow is at $k=3$. (Inertia always falls as $k$ rises, so you cannot just pick the smallest.) The silhouette score gives a second opinion: $0.753,\;\mathbf{0.811},\;0.659,\;0.479,\;0.320$ for $k=2,\dots,6$: best at $k=3$ again. Use silhouette_score(X, km.labels_) (needs $k\ge2$).

At $k=3$ the clusters are exactly the three blobs (100 points each, adjusted Rand index $=1.0$ against the truth). The centres were about $(-8.54,\,5.57)$, $(-1.26,\,4.39)$ and $(9.62,\,0.56)$.

Where k-means fails

K-means draws round, equal-sized clusters (it assigns by distance to a centre). On the "two moons" data (300 points, two interleaved crescents) it cuts straight across both (adjusted Rand index only 0.23). It is also sensitive to outliers and needs $k$ in advance.

Idea 2: DBSCAN

from sklearn.cluster import DBSCAN
db = DBSCAN(eps=0.2, min_samples=5).fit(X)       # defaults: eps=0.5, min_samples=5
db.labels_                                         # cluster ids; -1 means NOISE
db.core_sample_indices_

The idea: a point is a core point if at least min_samples points (itself included) lie within distance eps of it. Clusters are chains of core points and their neighbours; points reachable from no core point are noise (label $-1$).

• No n_clusters: DBSCAN finds the number itself.
• Finds arbitrarily shaped clusters: on the moons it recovered both crescents perfectly (adjusted Rand index 1.0, 0 noise points) where k-means failed.
• Handles outliers by labelling them $-1$.
• Weak points: needs a sensible eps (a global density), struggles when clusters have very different densities, and it has no predict (it cannot label new points directly).

Idea 3: GaussianMixture

from sklearn.mixture import GaussianMixture
gm = GaussianMixture(n_components=3, covariance_type="full", n_init=3, random_state=0).fit(X)
gm.weights_       # mixing weights pi_k (sum to 1): about 0.33 each here
gm.means_         # the component means (same as the three blobs)
gm.covariances_   # one covariance per component
gm.predict(X)           # hard label: the most probable component
gm.predict_proba(X)     # SOFT assignments (responsibilities): each row sums to 1
gm.score_samples(X)     # log-likelihood of each point
gm.bic(X), gm.aic(X)

This is exactly the EM algorithm of MLT week 4. Defaults: n_components=1, covariance_type="full" (others: "tied", "diag", "spherical"). A point deep inside one blob has responsibilities like $[0,\,0,\,1]$; a point between two blobs would show a split such as $[0.5,\,0.5,\,0]$.

Parameter count (full covariance, $d=2$, $K=3$): $K\big(d+\frac{d(d+1)}2\big)+(K-1)=3(2+3)+2=\mathbf{17}$, which is what _n_parameters() reports.

Choosing the number of components with BIC (lower is better):

The minimum is at $K=3$: the same answer as the elbow and the silhouette. BIC = $-2\ln L+(\#\text{params})\ln n$: more components always raise the likelihood but also the penalty.

GMM vs k-means: k-means = hard assignment, equal round clusters; GMM = soft assignment, clusters can be elliptical and differently sized, and it is a generative model (gm.sample(n) draws new points).

Idea 4: PCA

from sklearn.decomposition import PCA
pca = PCA(n_components=2).fit(X)        # centres the data for you, but does NOT scale it
Z = pca.transform(X)                      # the new coordinates, shape (n, 2)
pca.components_                           # the principal directions, shape (2, n_features)
pca.explained_variance_ratio_             # share of variance per component
pca.inverse_transform(Z)                  # back to the original space (an approximation)

On the digits data (1797 images, 64 pixels):

• The first two components explain only 28.5% of the variance (0.15 and 0.14), the first ten 73.8%.
• PCA(n_components=0.95) (a number between 0 and 1) keeps as many components as needed to reach 95% of the variance: it chose 29 components, compressing 64 numbers to 29 with little loss.
• Reconstruction error (mean squared error between the original and inverse_transform): 13.4 with 2 components, 0.77 with 30. More components, smaller error: exactly the "discarded eigenvalues" of MLT week 1.
• The explained variances add up to the total variance (1202.15 for the digits), as the theory says.

To plot the cumulative variance: np.cumsum(pca.explained_variance_ratio_) and pick the "elbow" or the 95% mark.

Other facts:

• PCA subtracts the mean for you but does not standardise: do that first when features have different units.
• IncrementalPCA has partial_fit for huge data; KernelPCA(kernel="rbf") for curved structure; n_components="mle" guesses the dimension.
• Use PCA inside a pipeline before a classifier to reduce dimensions and speed up training or reduce noise; for 2-D pictures of high-dimensional data you can also use t-SNE (visualisation only).

Solved examples

Your turn

Exam radar

Quick quiz

Recap

Practise with real exam questions

These come straight from past papers (with the real figures and the official answer key). Open 58 real questions on this week   All End Term papers