Show Sidebar Hide Sidebar

Empirical evaluation of the impact of k-means initialization in Scikit-learn

Evaluate the ability of k-means initializations strategies to make the algorithm convergence robust as measured by the relative standard deviation of the inertia of the clustering (i.e. the sum of distances to the nearest cluster center).

The first plot shows the best inertia reached for each combination of the model (KMeans or MiniBatchKMeans) and the init method (init="random" or init="kmeans++") for increasing values of the n_init parameter that controls the number of initializations.

The second plot demonstrate one single run of the MiniBatchKMeans estimator using a init="random" and n_init=1. This run leads to a bad convergence (local optimum) with estimated centers stuck between ground truth clusters.

The dataset used for evaluation is a 2D grid of isotropic Gaussian clusters widely spaced.

New to Plotly?

Plotly's Python library is free and open source! Get started by downloading the client and reading the primer.
You can set up Plotly to work in online or offline mode, or in jupyter notebooks.
We also have a quick-reference cheatsheet (new!) to help you get started!

Version

In [1]:
import sklearn
sklearn.__version__
Out[1]:
'0.18'

Imports

This tutorial imports shuffle, check_random_state, MiniBatchKMeans and KMeans.

In [2]:
import plotly.plotly as py
import plotly.graph_objs as go

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import matplotlib

from sklearn.utils import shuffle
from sklearn.utils import check_random_state
from sklearn.cluster import MiniBatchKMeans
from sklearn.cluster import KMeans

Calculations

In [3]:
random_state = np.random.RandomState(0)

# Number of run (with randomly generated dataset) for each strategy so as
# to be able to compute an estimate of the standard deviation
n_runs = 5

# k-means models can do several random inits so as to be able to trade
# CPU time for convergence robustness
n_init_range = np.array([1, 5, 10, 15, 20])

# Datasets generation parameters
n_samples_per_center = 100
grid_size = 3
scale = 0.1
n_clusters = grid_size ** 2


def make_data(random_state, n_samples_per_center, grid_size, scale):
    random_state = check_random_state(random_state)
    centers = np.array([[i, j]
                        for i in range(grid_size)
                        for j in range(grid_size)])
    n_clusters_true, n_features = centers.shape

    noise = random_state.normal(
        scale=scale, size=(n_samples_per_center, centers.shape[1]))

    X = np.concatenate([c + noise for c in centers])
    y = np.concatenate([[i] * n_samples_per_center
                        for i in range(n_clusters_true)])
    return shuffle(X, y, random_state=random_state)

Quantitative evaluation of various init methods

In [4]:
fig = plt.figure()
plots = []
legends = []

cases = [
    (KMeans, 'k-means++', {}),
    (KMeans, 'random', {}),
    (MiniBatchKMeans, 'k-means++', {'max_no_improvement': 3}),
    (MiniBatchKMeans, 'random', {'max_no_improvement': 3, 'init_size': 500}),
]

for factory, init, params in cases:
    print("Evaluation of %s with %s init" % (factory.__name__, init))
    inertia = np.empty((len(n_init_range), n_runs))

    for run_id in range(n_runs):
        X, y = make_data(run_id, n_samples_per_center, grid_size, scale)
        for i, n_init in enumerate(n_init_range):
            km = factory(n_clusters=n_clusters, init=init, random_state=run_id,
                         n_init=n_init, **params).fit(X)
            inertia[i, run_id] = km.inertia_
    p = go.Scatter(y=inertia.mean(axis=1),
                   x=n_init_range,
                   name="%s with %s init" % (factory.__name__, init),
                   error_y=dict(visible=True, arrayminus=inertia.std(axis=1)),
                   line=dict(width=1))
    plots.append(p)

layout=go.Layout(title = "Mean inertia for various k-means init across %d runs" % n_runs,
                 hovermode='closest',
                 xaxis = dict(title='n_init'),
                 yaxis = dict(title='inertia'))

fig = go.Figure(data=plots, layout=layout)
Evaluation of KMeans with k-means++ init
Evaluation of KMeans with random init
Evaluation of MiniBatchKMeans with k-means++ init
Evaluation of MiniBatchKMeans with random init
In [5]:
py.iplot(fig)
Out[5]:

Qualitative visual inspection of the convergence

In [6]:
X, y = make_data(random_state, n_samples_per_center, grid_size, scale)
km = MiniBatchKMeans(n_clusters=n_clusters, init='random', n_init=1,
                     random_state=random_state).fit(X)
plot=[]

for k in range(n_clusters):
    my_members = km.labels_ == k
    c = matplotlib.colors.colorConverter.to_rgb(cm.spectral(float(k) / n_clusters, 1))
    c ='rgb'+str(c)
    cluster = go.Scatter(x=X[my_members, 0], y=X[my_members, 1],
                         showlegend=False,
                         mode='markers', marker=dict(color=c, size=4,
                                           line=dict(color='black', width=1)))
    center = km.cluster_centers_[k]
    
    cluster_center = go.Scatter(x=[center[0]],y=[center[1]],
                                showlegend=False,
                                mode='markers', marker=dict(color=c, size=10,
                                         line=dict(color='black', width=1)))
    plot.append(cluster)
    plot.append(cluster_center)

layout = go.Layout(title="Example cluster allocation with a single random init<br>"
                          "with MiniBatchKMeans",
                   hovermode='closest',
                   xaxis=dict(zeroline = False),
                   yaxis=dict(zeroline = False))

fig = go.Figure(data=plot, layout=layout)
In [7]:
py.iplot(fig)
Out[7]:

License

Author:

    Olivier Grisel <olivier.grisel@ensta.org>

License:

    BSD 3 clause
Still need help?
Contact Us

For guaranteed 24 hour response turnarounds, upgrade to a Developer Support Plan.