# How to learn new clusters on residuals of KMeans

I have 800k long-lat locations that I've clustered into 40k locations. I reverse-geocoded these using a free API.

As I get a new "budget" (40k) next month, I'll be able to process another 40k locations. While there will be a few new points, I would like to retrain the model to decrease the error.

I could for example take the 40k points with the highest error when projected the coordinates back using the KMeans algorithm.

But... this is going to be suboptimal. It could be that there will be 10 exactly the same locations with the same high error. I would like to kind of run a second KMeans, trying to maximize the explainable variance.

So what can I do to make sure I explain as much variance using the second time 40k budget? How should I do this technically?

Summarizing the situation:

Month 1: 800k locations -> 40k locations
Month 2: 810k locations -> 40k existing locations + 40k new budget


That is, I would like to retrain KMeans with fixed 40k clusters and finding 40k new clusters that would be optimal given the 810k locations. Maybe I should take the residuals from the first model using the 810k locations, train a new one, and apply multiple models in a row?

I'm curious to both a technical way to solve it (using something close to a KMeans implementation) as well as a scikit-learn specific solution.

• So, the next month you will get another 40k points to assign to clusters, or you will try to find 40k new clusters in addition to 40k that you have found this month? – Milos Sep 15 '18 at 9:30
• @Milos The next month you get another 40k free API calls, and you get 10k new locations. – PascalVKooten Sep 15 '18 at 9:46
• @Milos So your second thought is accurate. – PascalVKooten Sep 15 '18 at 9:53

Depending on what you mean by "fixing" the old clusters, you could proceed as follows:

(1) If you want to fix the old clusters so as to make the old points stay assigned to their clusters, then use the assignments and centroids from this month to initialize the old clusters, randomly selecting new centroids for the new 40k clusters and randomly populating them with new 10k points. After that, skip updating the old centroids and their points' assignments.

(2). On the other hand, if you want to fix only the centroids of the old clusters, then do the initialization step as in (1), skip updating the means in the iterations of the algorithm, but allow the points to change clusters.

(3) Finally, you could run KMeans as usual, skipping no update steps, but initializing the algorithm as in the first and second approaches. This way, the new results will be influenced by the old clustering, but it will not be guaranteed that the old clusters would stay the same (as in (1), or partially in (2)).

• I think indeed I'm thinking about (2), but feel like the API of scikit-learn does not help me with this. There's no easy way to just add something on top of the interface, particularly initializing the centers and freezing 40k/80k in the second month. – PascalVKooten Sep 15 '18 at 19:42
• @PascalVKooten You could code KMeans from scratch on your own, or look at the source of the scikit's implementation and customize it to fit your needs. – Milos Sep 15 '18 at 20:12
• Yea, I started doing that already, but it's hideous. I'll post it soon. – PascalVKooten Sep 15 '18 at 20:20

I needed to copy some code from kmeans module as it refers to a function defined there, called _mini_batch_step.

Here is the action:

budget = 400
old_centroids=None
for i in range(3):
m, new_clusters = do_mini_batch(X, budget, old_centroids)
print("n_clusters", m.n_clusters)
centroids = m.cluster_centers_
print("first center (shows the centers do not get overwritten per epoch)")
print(centroids[:1])
error = np.mean(np.abs(X - m.cluster_centers_[m.predict(X)]))
print("error", error)
print()
old_centroids = centroids

# prints

# n_clusters 400
# first center (shows the centers do not get overwritten per epoch)
# [[0.28595093 0.56208682]]
# error 0.011584586126928233

# n_clusters 800
# first center (shows the centers do not get overwritten per epoch)
# [[0.28595093 0.56208682]]
# error 0.00873081143320314

# n_clusters 1200
# first center (shows the centers do not get overwritten per epoch)
# [[0.28595093 0.56208682]]
# error 0.006874430587994398


Note that instead of doing this 3 times, you should make sure you just store the new_clusters per month. Concatenating these will make it easy to keep them as history.

Here is the setup that I'm using:

import numpy as np
from sklearn.cluster import MiniBatchKMeans
import sklearn.cluster.k_means_

X = np.random.random((8000, 2))

def patched_mini_batch_step(previous_centers):
if previous_centers is None:
frozen_n = 0
else:
frozen_n = len(previous_centers)
def _mini_batch_step(X, x_squared_norms, centers, counts,
old_center_buffer, compute_squared_diff,
distances, random_reassign=False,
random_state=None, reassignment_ratio=.01,
verbose=False):
centers[:frozen_n] = previous_centers
"""Incremental update of the centers for the Minibatch K-Means algorithm.

Parameters
----------

X : array, shape (n_samples, n_features)
The original data array.

x_squared_norms : array, shape (n_samples,)
Squared euclidean norm of each data point.

centers : array, shape (k, n_features)
The cluster centers. This array is MODIFIED IN PLACE

counts : array, shape (k,)
The vector in which we keep track of the numbers of elements in a
cluster. This array is MODIFIED IN PLACE

distances : array, dtype float, shape (n_samples), optional
If not None, should be a pre-allocated array that will be used to store
the distances of each sample to its closest center.
May not be None when random_reassign is True.

random_state : int, RandomState instance or None, optional, default: None
If int, random_state is the seed used by the random number generator;
If RandomState instance, random_state is the random number generator;
If None, the random number generator is the RandomState instance used
by np.random.

random_reassign : boolean, optional
If True, centers with very low counts are randomly reassigned
to observations.

reassignment_ratio : float, optional
Control the fraction of the maximum number of counts for a
center to be reassigned. A higher value means that low count
centers are more likely to be reassigned, which means that the
model will take longer to converge, but should converge in a
better clustering.

verbose : bool, optional, default False
Controls the verbosity.

compute_squared_diff : bool
If set to False, the squared diff computation is skipped.

old_center_buffer : int
Copy of old centers for monitoring convergence.

Returns
-------
inertia : float
Sum of squared distances of samples to their closest cluster center.

squared_diff : numpy array, shape (n_clusters,)
Squared distances between previous and updated cluster centers.

"""
# Perform label assignment to nearest centers
nearest_center, inertia = sklearn.cluster.k_means_._labels_inertia(X, x_squared_norms, centers,
distances=distances)

if random_reassign and reassignment_ratio > 0:
random_state = sklearn.cluster.k_means_.check_random_state(random_state)
# Reassign clusters that have very low counts
to_reassign = counts < reassignment_ratio * counts.max()
to_reassign[:frozen_n] = False
# pick at most .5 * batch_size samples as new centers
if to_reassign.sum() > .5 * X.shape[0]:
indices_dont_reassign = np.argsort(counts)[int(.5 * X.shape[0]):]
to_reassign[indices_dont_reassign] = False
n_reassigns = to_reassign.sum()
if n_reassigns:
# Pick new clusters amongst observations with uniform probability
new_centers = random_state.choice(X.shape[0], replace=False,
size=n_reassigns)
if verbose:
print("[MiniBatchKMeans] Reassigning %i cluster centers."
% n_reassigns)

if sklearn.cluster.k_means_.sp.issparse(X) and not sklearn.cluster.k_means_.sp.issparse(centers):
sklearn.cluster.k_means_.assign_rows_csr(X, new_centers.astype(np.intp),
np.where(to_reassign)[0].astype(np.intp),
centers)
else:
centers[to_reassign] = X[new_centers]
# reset counts of reassigned centers, but don't reset them too small
# to avoid instant reassignment. This is a pretty dirty hack as it
# also modifies the learning rates.
counts[to_reassign] = np.min(counts[~to_reassign])

# implementation for the sparse CSR representation completely written in
# cython
if sklearn.cluster.k_means_.sp.issparse(X):
return inertia, sklearn.cluster.k_means_._k_means._mini_batch_update_csr(
X, x_squared_norms, centers, counts, nearest_center,
old_center_buffer, compute_squared_diff)

# dense variant in mostly numpy (not as memory efficient though)
k = centers.shape[0]
squared_diff = 0.0
for center_idx in range(k):
# find points from minibatch that are assigned to this center

if count > 0:
if compute_squared_diff:
old_center_buffer[:] = centers[center_idx]

# inplace remove previous count scaling
centers[center_idx] *= counts[center_idx]

# inplace sum with new points members of this cluster

# update the count statistics for this center
counts[center_idx] += count

# inplace rescale to compute mean of all points (old and new)
# Note: numpy >= 1.10 does not support '/=' for the following
# expression for a mixture of int and float (see numpy issue #6464)
centers[center_idx] = centers[center_idx] / counts[center_idx]

# update the squared diff if necessary
if compute_squared_diff:
diff = centers[center_idx].ravel() - old_center_buffer.ravel()
squared_diff += np.dot(diff, diff)

return inertia, squared_diff
return _mini_batch_step

def do_mini_batch(X, num_new_centroids, old_centroids=None, **km_kwargs):
num_old_centroids = len(old_centroids) if old_centroids is not None else 0
print(num_old_centroids)
n_clusters = num_old_centroids + num_new_centroids
try:
sklearn.cluster.k_means_._mini_batch_step = patched_mini_batch_step(old_centroids)
m = MiniBatchKMeans(n_clusters, init_size=3*n_clusters, **km_kwargs)
m.fit(X)
finally:
sklearn.cluster.k_means_._mini_batch_step=old_mini_batch_step
new_clusters = m.cluster_centers_[-num_new_centroids:]
return m, new_clusters