In K-means let's assume there are M prototypes denoted by

\(Z = {z_1, z_2, \cdots , z_M}\)

This set is usually smaller than the original data set. If the data points reside in a p-dimensional Euclidean space, the prototypes reside in the same space. They will also be p-dimensional vectors. They may not be samples from the training data set, however, they should well represent the training dataset.

Each training sample is assigned to one of the prototypes. This is denoted using the assignment function by \(A(\cdot)\). Then \(A(x_i) = j\) means the \(i^{th}\) training sample is assigned to the \(j^{th}\) prototype, that is, the\(j^{th}\) prototype is used to approximate or represent point \(x_i\).

In k-means, we need to solve two unknowns. The first is to select a set of prototypes; the second is the assignment function.

The Objective Function in K-Means

In K-means, the optimization criterion is to minimize the total squared error between the training samples and their representative prototypes.  This is equivalent to minimizing the trace of the pooled within covariance matrix. The objective function is:

\(\text{arg } \underset{\mathcal{Z}, A}{\text{ min }}\sum_{i=1}^{N}||x_i-z_{A(x_i)}||^2 \)

For every point \(x_i\), first we find which z_j will be used to represent \(x_i\). As we mentioned above, if \(A(x_i) = j\) then the prototype \(z_j\) will be used. Hence, we compute the distance between \(x_i\) with \(z_{A(x_i)}\) . The norm in the above formula denotes the Euclidean distance.

Ideally, we want every data point to be perfectly represented by a prototype.  Since this is generally impossible, we try to minimize the total discrepancy.

Denote the objective function by:

\(L(\mathcal{Z},A)=\sum_{i=1}^{N}||x_i-z_{A(x_i)}||^2\)

From the perspective of data compression, this seems obvious.   But from the clustering perspective, one may ask why we use such an optimization criterion.

Here is the intuition for this...

Evaluating clustering results is inevitably subjective.  For classification, we have the true labels. Therefore, we can compute the error rate to assess the result.  However, in (unsupervised) clustering, we don’t have the class labels given for the training data. 

One heuristic generally accepted is that points in the same cluster should be tight and points in different groups should be as far apart as possible. The k-means algorithm reflects the heuristic by attempting to minimize the total within-cluster distances between each data point and its corresponding prototype.

Necessary Conditions

How do we optimize the objective function of k-means?

What follows are necessary conditions for an optimal solution. Necessary means that if a solution is indeed optimal, then it has to satisfy these conditions, but not vice versa. This means we might have a solution satisfying the conditions that are not optimal.  We call such solutions locally optimal, in contrast to globally optimal.

Condition 1: If the prototypes, Z, are fixed, the optimal assignment function \(A(\cdot)\) (using Euclidean distance) should follow the nearest neighbor rule, that is,

\(A(x_i) = \text{ arg }\underset{j \in {1,2,...,M}}{ min } || x_i − z_j || \)

This means that we should always find among those prototypes the one closest to the new point. Intuitively this should be clear because we want to minimize the summed distance between every point and its prototype. We have the freedom to choose the prototype for each data point.  Obviously, for a single point, the best prototype to choose is the closest one to it.

Condition 2 : Given a fixed assignment function, \(A(\cdot)\), the prototype, \(z_j\), should be the average (centroid) of all the samples assigned to the \(j^{th}\) prototype:

\(z_j=\frac{\Sigma_{i:A(x_i)=j}x_i}{N_j}\)

where \(N_j\) is the number of samples assigned to prototype j.

This is also easy to prove.  For a given group of data points associated with a single prototype, the total sum of squared distances is a quadratic function in terms of the prototype.  This function is convex and its minimum is unique and is given by setting the first order derivatives to zero, which leads to the arithmetic mean of the data points.

Sometimes we call the assignment function "the partition" and the prototypes "the centroids."

The two necessary conditions on the previous page inspired the k-means algorithm.

It is interesting to note that in data compression for electrical engineers, k-means is better known as the Lloyd algorithm.

The k-means algorithm alternates the two steps:

• For a fixed set of centroids (prototypes), optimize A(•) by assigning each sample to its closest centroid using Euclidean distance.
• Update the centroids by computing the average of all the samples assigned to it.

The algorithm converges since, after each iteration, the objective function is non-increasing.

Based on the necessary conditions, the algorithm improves the prototypes and the assignment function iteratively. It is possible that the steps taken might not lead to any change in either the prototypes or the partition. When either becomes static, the algorithm has converged and there is no need to continue.

How do we decide when to stop?

One criterion for stopping is if we observe the assignment functions in the two iterations are exactly the same. If the assignment function doesn't change anymore, then the prototypes won't change either (and vice versa).

In practice, we often stop when the decrease in the objective function becomes small. We can compute the ratio between the decrease and the value of the objective function in the previous iteration. If the ratio is below a threshold, e.g., 0.0001, the iteration is stopped.

So, how do we initiate this procedure? How does it get started?

We bypassed this issue earlier. We had an example where we showed that k-means can result in different solutions for different initializations. The problem is that we do not have any theoretical guidance on how to pick the initial prototypes. Many times researchers don't know whether the initialization is a good one until they have completed the k-means algorithm and have the results in hand. Different ways have been proposed to initialize the k-means algorithm.

The simplest scheme is just to randomly pick some data points from the training data set and treat these as the initial prototypes. There are some other more sophisticated ways for initialization. Sometimes, people try different initializations and pick one that results in the minimum value for the objective function.

There are also methods that randomly pick the prototypes using another method to start the k-means iteration. For example, one may just pick a set of initial prototypes as random samples of a certain distribution.

Initialization in the above simulation:

• Generated M random vectors with independent dimensions. For each dimension, the feature is uniformly distributed in [-1, 1].

For instance, suppose the dimension is p and the number of prototypes is M. Randomly pick a point in [-1,1] and do this p times to obtain a p-dimensional vector. If the random number generator generates a value in the range [0,1] you can simply multiply the value by 2 and subtract 1 to achieve a uniform distribution on [-1,1].

– Linearly transform the jth variable, \(Z_j , j = 1, 2, \cdots , p\) in each prototype (a vector) by: \(Z_js_j + m_j\), where \(s_j\) is the sample standard deviation of dimension j and \(m_j\) is the sample mean of dimension j, both computed using the training data.

During a k-means iteration, it may happen that a prototype has no data point assigned to it. We call this an “empty” prototype. When this happens, we can remove this prototype to reduce the number of prototypes by 1. If we don’t want to reduce the number of prototypes, we can further divide points in a non-empty prototype arbitrarily into two groups and compute the means for the two groups. The two means will serve as two prototypes to replace the previous prototype.

Agglomerative clustering can be used as long as we have pairwise distances between any two objects. The mathematical representation of the objects is irrelevant when the pairwise distances are given. Hence agglomerative clustering readily applies for non-vector data.

Let's denote the data set as \(A = {x_1, \cdots , x_n}\).

The agglomerative clustering method is also called a bottom-up method as opposed to k-means or k-center methods that are top-down. In a top-down method, a data set is divided into more and more clusters. In a bottom-up approach, all the data points are treated as individual clusters to start with and gradually merged into bigger and bigger clusters.

In agglomerative clustering, clusters are generated hierarchically. We start by taking every data point as a cluster. Then we merge two clusters at a time. In order to decide which two clusters to merge, we compare the pairwise distances between any two clusters and pick a pair with the minimum distance. Once we merge two clusters into a bigger one, a new cluster is created. The distances between this new cluster and the existing ones are not given. Some scheme has to be used to obtain these distances based on the two merged clusters. We call this the update of distances. Various schemes of updating distances will be described shortly.

We can keep merging clusters until all the points are merged into one cluster. A tree can be used to visualize the merging process. This tree is called a dendrogram.

The key technical detail here is how we define the between-cluster distance. At the bottom level of the dendrogram where every cluster only contains a single point, the between-cluster distance is simply the Euclidian distance between the points.  (In some applications, the between-point distances are pre-given.)  However, once we merge two points into a cluster how do we get the distance between the new cluster and the existing clusters?

The idea is to somehow aggregate the distances between the objects contained in the clusters.

• For clusters containing only one data point, the between-cluster distance is the between-object distance.
• For clusters containing multiple data points, the between-cluster distance is an agglomerative version of the between-object distances. There are a number of ways to accomplish this. Examples of these aggregated distances include the minimum or maximum between-object distances for pairs of objects across the two clusters.

For example, suppose we have two clusters and each one has three points. One thing we could do is to find distances between these points. In this case, there would be nine distances between pairs of points across the two clusters. The minimum (or maximum, or average) of the nine distances can be used as the distance between the two clusters.

How do we decide which aggregation scheme to use? Depending on how we update the distances, dramatically different results may come up. Therefore, it is always good practice to look at the results using scatterplots or other visualization methods instead of blindly taking the output of any algorithm. Clustering is inevitably subjective since there is no gold standard.

Normally the agglomerative between-cluster distance can be computed recursively. The aggregation as explained above sounds computationally intensive and seemingly impractical. If we have two very large clusters, we have to check all the between-object distances for all pairs across the two clusters. The good news is that for many of these aggregated distances, we can update them recursively without checking all the pairs of objects. The update approach will be described soon.

Psudo Code Steps:

1. Begin with n clusters, each containing one object and we will number the clusters 1 through n.
2. Compute the between-cluster distance D(r, s) as the between-object distance of the two objects in r and s respectively, r, s =1, 2, ..., n. Let the square matrix D = (D(r, s)). If the objects are represented by quantitative vectors we can use Euclidean distance.
3. Next, find the most similar pair of clusters r and s, such that the distance, D(r, s), is minimum among all the pairwise distances.
4. Merge r and s to a new cluster t and compute the between-cluster distance D(t, k) for any existing cluster k ≠ r, s . Once the distances are obtained,  delete the rows and columns corresponding to the old cluster r and s in the D matrix, because r and s do not exist anymore. Then add a new row and column in D corresponding to cluster t.
5. Repeat Step 3 a total of n − 1 times until there is only one cluster left.