K-means and Single-Linkage Hierarchical Clustering

The four observation example from K-means by hand

1 Aim and how to run this file

We will solve the example in the supplied slides using two methods:

1. K-means, with two clusters and starting centres at A and B.
2. Hierarchical clustering with single linkage, followed by a cut that gives two clusters.

The data are included below. Open this file in RStudio and click Knit to produce an HTML document containing the R commands, their output, and the plots. You can also run the code chunks in order.

The clustering and plotting functions come with R. Knitting requires the rmarkdown package and its dependencies. If necessary, run the following once in the R console:

install.packages("rmarkdown")

The final answer is {A, B} and {C, D} for both methods. The sections below explain how each method reaches that answer.

2 Enter and plot the data

The four observations are A = (1, 1), B = (2, 1), C = (4, 3), and D = (5, 4), exactly as in the slides.

X <- data.frame(
  x = c(1, 2, 4, 5),
  y = c(1, 1, 3, 4),
  row.names = c("A", "B", "C", "D")
)

X

#>   x y
#> A 1 1
#> B 2 1
#> C 4 3
#> D 5 4

Each row is one observation. The columns x and y are the two numerical features used for clustering. A, B, C, and D are observation names. They are not known class labels.

This is an unsupervised learning problem: we use the numerical features to find groups without being given the correct group for each observation.

plot(
  X$x, X$y,
  pch = 19, cex = 1.6, col = "#526477",
  xlim = c(0, 6), ylim = c(0, 5), asp = 1,
  xlab = "x", ylab = "y",
  main = "Four observations"
)
text(X$x, X$y, labels = rownames(X), pos = 3, cex = 1.1)

The original observations, before assigning any clusters.

Look at the plot: A and B are close to each other, and C and D are close to each other.

We use the original coordinates throughout so the distances match the slides. For real data, consider the units and scales of the variables: changing a variable’s scale changes its contribution to distance.

3 Solve the example with k-means

3.1 Choose the same starting centres as the slides

We choose \(K = 2\) and start with

\[ c_1 = (1,1) \quad\text{and}\quad c_2 = (2,1). \]

These are the coordinates of A and B.

initial_centres <- as.matrix(X[c("A", "B"), ])
initial_centres

#>   x y
#> A 1 1
#> B 2 1

X[c("A", "B"), ] selects rows A and B and both columns. The result supplies the actual starting centres to kmeans().

3.2 Run kmeans()

km <- kmeans(
  X,
  centers = initial_centres,
  algorithm = "Lloyd",
  iter.max = 100
)

The arguments mean:

• X: the numerical data to cluster.
• centers = initial_centres: start at A and B; the two rows specify two clusters.
• algorithm = "Lloyd": alternate between assigning observations to their nearest centre and updating the centres to the group means, as in the slides.
• iter.max = 100: allow up to 100 iterations; this is a limit, not a request to perform 100 iterations.

R’s default k-means algorithm is Hartigan-Wong. We explicitly select Lloyd’s algorithm to match the hand calculation. Supplying the starting centres makes this example reproducible without random initialisation. See the official kmeans documentation.

3.3 Read the cluster assignments and centres

km$cluster

#> A B C D 
#> 1 1 2 2

data.frame(X, K_means_cluster = km$cluster)

#>   x y K_means_cluster
#> A 1 1               1
#> B 2 1               1
#> C 4 3               2
#> D 5 4               2

The assignment is A and B in one cluster, and C and D in the other.

km$centers

#>     x   y
#> 1 1.5 1.0
#> 2 4.5 3.5

km$size

#> [1] 2 2

The final centres are

\[ c_1 = \left(\frac{1+2}{2},\frac{1+1}{2}\right) = (1.5,1), \qquad c_2 = \left(\frac{4+5}{2},\frac{3+4}{2}\right) = (4.5,3.5). \]

Each cluster contains two observations. The $ symbol extracts a named part of the saved result: for example, km$centers extracts the centres.

3.4 (Optional) Plot the k-means solution

cluster_colours <- c("#0072B2", "#D55E00")

plot(
  X$x, X$y,
  col = cluster_colours[km$cluster],
  pch = 19, cex = 1.6,
  xlim = c(0, 6), ylim = c(0, 5), asp = 1,
  xlab = "x", ylab = "y",
  main = "K-means: two clusters"
)
text(X$x, X$y, labels = rownames(X), pos = 3, cex = 1.1)
points(
  km$centers[, "x"], km$centers[, "y"],
  pch = 4, cex = 2, lwd = 2,
  col = cluster_colours
)
text(
  km$centers[, "x"], km$centers[, "y"],
  labels = c("c1", "c2"), pos = 1,
  col = cluster_colours
)
legend(
  "topleft",
  legend = c("Cluster 1", "Cluster 2", "Centre"),
  col = c(cluster_colours, "black"),
  pch = c(19, 19, 4), bty = "n"
)

K-means assigns A and B to one cluster and C and D to the other. Crosses show the final centres.

4 (Optional) Connect the R solution to the hand calculations

This section reproduces the assignment and update steps in the slides. It can be used to explain what the kmeans() call does.

For an observation \(P=(x,y)\) and a centre \(c=(a,b)\), the squared Euclidean distance is

\[ d^2(P,c) = (x-a)^2 + (y-b)^2. \]

Comparing squared distances gives the same nearest centre as comparing their square roots.

4.1 Round 1: assign points to the initial centres

d2_c1 <- (X$x - initial_centres[1, "x"])^2 +
         (X$y - initial_centres[1, "y"])^2

d2_c2 <- (X$x - initial_centres[2, "x"])^2 +
         (X$y - initial_centres[2, "y"])^2

round1 <- data.frame(
  Squared_to_c1 = d2_c1,
  Squared_to_c2 = d2_c2,
  Cluster = ifelse(d2_c1 <= d2_c2, 1, 2),
  row.names = rownames(X)
)

round1

#>   Squared_to_c1 Squared_to_c2 Cluster
#> A             0             1       1
#> B             1             0       2
#> C            13             8       2
#> D            25            18       2

ifelse(d2_c1 <= d2_c2, 1, 2) assigns a point to cluster 1 when its first squared distance is smaller, and to cluster 2 otherwise. This code assigns an exact tie to cluster 1; no ties occur in this example.

The first assignment is {A} and {B, C, D}. In particular, C joins cluster 2 because \(8 < 13\).

4.2 Round 1: update the centres

colMeans() calculates the mean of each coordinate within a cluster.

centres_round1 <- rbind(
  c1 = colMeans(X[round1$Cluster == 1, , drop = FALSE]),
  c2 = colMeans(X[round1$Cluster == 2, , drop = FALSE])
)

centres_round1

#>          x       y
#> c1 1.00000 1.00000
#> c2 3.66667 2.66667

The updated centres are \(c_1=(1,1)\) and \(c_2=(11/3,8/3)\). The condition round1$Cluster == 1 selects the observations assigned to cluster 1. drop = FALSE keeps the selected data in a two-column table, including when a cluster has just one observation.

Keep the stored values at full precision for the next calculation.

4.3 Round 2: assign points again

d2_c1 <- (X$x - centres_round1[1, "x"])^2 +
         (X$y - centres_round1[1, "y"])^2

d2_c2 <- (X$x - centres_round1[2, "x"])^2 +
         (X$y - centres_round1[2, "y"])^2

round2 <- data.frame(
  Squared_to_c1 = d2_c1,
  Squared_to_c2 = d2_c2,
  Cluster = ifelse(d2_c1 <= d2_c2, 1, 2),
  row.names = rownames(X)
)

round(round2, 4)

#>   Squared_to_c1 Squared_to_c2 Cluster
#> A             0        9.8889       1
#> B             1        5.5556       1
#> C            13        0.2222       2
#> D            25        3.5556       2

The second squared-distance column is \(89/9\), \(50/9\), \(2/9\), and \(32/9\). B changes to cluster 1 because \(1 < 50/9\). The groups are now {A, B} and {C, D}.

round(round2, 4) rounds the displayed table; it does not replace the stored, unrounded values.

4.4 Round 2: update the centres

centres_round2 <- rbind(
  c1 = colMeans(X[round2$Cluster == 1, , drop = FALSE]),
  c2 = colMeans(X[round2$Cluster == 2, , drop = FALSE])
)

centres_round2

#>      x   y
#> c1 1.5 1.0
#> c2 4.5 3.5

The new centres are (1.5, 1) and (4.5, 3.5), matching the kmeans() solution.

4.5 Final check: do any assignments change?

d2_c1 <- (X$x - centres_round2[1, "x"])^2 +
         (X$y - centres_round2[1, "y"])^2

d2_c2 <- (X$x - centres_round2[2, "x"])^2 +
         (X$y - centres_round2[2, "y"])^2

final_check <- data.frame(
  Squared_to_c1 = d2_c1,
  Squared_to_c2 = d2_c2,
  Previous_cluster = round2$Cluster,
  New_cluster = ifelse(d2_c1 <= d2_c2, 1, 2),
  row.names = rownames(X)
)

final_check

#>   Squared_to_c1 Squared_to_c2 Previous_cluster New_cluster
#> A          0.25          18.5                1           1
#> B          0.25          12.5                1           1
#> C         10.25           0.5                2           2
#> D         21.25           0.5                2           2

No assignments change. Therefore, recalculating the means would leave the centres unchanged, and the algorithm stops. There are two rounds that change the centres, followed by a final assignment check.

5 Solve the example with single-linkage hierarchical clustering

5.1 Calculate Euclidean distances between observations

Hierarchical clustering starts with each observation in its own cluster. At each step it merges two clusters.

We first calculate the ordinary Euclidean distances between all pairs of observations:

\[ d(P,Q) = \sqrt{(x_P-x_Q)^2 + (y_P-y_Q)^2}. \]

d <- dist(X, method = "euclidean")
distance_matrix <- as.matrix(d)

round(distance_matrix, 4)

#>        A      B      C      D
#> A 0.0000 1.0000 3.6056 5.0000
#> B 1.0000 0.0000 2.8284 4.2426
#> C 3.6056 2.8284 0.0000 1.4142
#> D 5.0000 4.2426 1.4142 0.0000

For example,

\[ d(A,B)=1,\qquad d(C,D)=\sqrt{2},\qquad d(B,C)=\sqrt{8}. \]

dist() creates a distance object. as.matrix() displays it as a full table. The diagonal entries are zero, and the table is symmetric. These are Euclidean distances, while the earlier k-means assignment tables show squared distances. See the official dist documentation.

5.2 Apply single linkage

The single-linkage distance between two clusters is the distance between their closest pair of observations, with one observation in each cluster:

\[ d_{\mathrm{single}}(G_1,G_2) = \min_{P\in G_1,\ Q\in G_2} d(P,Q). \]

hc_single <- hclust(d, method = "single")

The argument must be method = "single" to request single linkage. hclust() uses complete linkage by default. See the official hclust documentation.

5.3 Understand the three merges

There are four observations, so there are three merges before everything belongs to one cluster.

Step	Clusters merged	Merge height	Closest pair
1	{A} and {B}	\(1\)	A and B
2	{C} and {D}	\(\sqrt{2}\approx 1.4142\)	C and D
3	{A, B} and {C, D}	\(\sqrt{8}\approx 2.8284\)	B and C

R stores the merge heights in the fitted object:

data.frame(
  Step = 1:3,
  Height = hc_single$height
)

#>   Step  Height
#> 1    1 1.00000
#> 2    2 1.41421
#> 3    3 2.82843

After the first merge, the distances between the remaining clusters are

\[ \begin{aligned} d_{\mathrm{single}}(\{A,B\},\{C\}) &= \min(\sqrt{13},\sqrt{8})=\sqrt{8},\\ d_{\mathrm{single}}(\{A,B\},\{D\}) &= \min(5,\sqrt{18})=\sqrt{18},\\ d(\{C\},\{D\}) &= \sqrt{2}. \end{aligned} \]

Thus C and D merge next. For the final merge, examine all four distances between the two remaining clusters:

between_cluster_distances <- distance_matrix[
  c("A", "B"), c("C", "D")
]

round(between_cluster_distances, 4)

#>        C      D
#> A 3.6056 5.0000
#> B 2.8284 4.2426

min(between_cluster_distances)

#> [1] 2.82843

The minimum is the distance between B and C, \(\sqrt{8}\). Single linkage uses the nearest pair across the two clusters. It does not calculate a distance between their means.

5.4 Draw the dendrogram

plot(
  hc_single,
  hang = -1,
  main = "Single-linkage hierarchical clustering",
  xlab = "Observations",
  ylab = "Merge height (Euclidean distance)",
  sub = ""
)
rect.hclust(hc_single, k = 2, border = cluster_colours)
abline(h = 2, lty = 2, col = "#666666")
legend(
  "topleft", legend = "Cut at height 2",
  lty = 2, col = "#666666", bty = "n", cex = 0.9
)

Single linkage merges A with B at height 1 and C with D at height sqrt(2). Cutting at height 2 gives two clusters.

Read the dendrogram from the bottom upward. The height of each horizontal joining line is the distance at which that merge occurs. With hang = -1, the observation labels share a common baseline.

The dashed line is at height 2. Since \(\sqrt{2}<2<\sqrt{8}\), the two lower merges have occurred and the final merge has not. The result is two clusters.

5.5 Cut the tree into two clusters

hc_groups <- cutree(hc_single, k = 2)
hc_groups

#> A B C D 
#> 1 1 2 2

data.frame(X, Single_linkage_cluster = hc_groups)

#>   x y Single_linkage_cluster
#> A 1 1                      1
#> B 2 1                      1
#> C 4 3                      2
#> D 5 4                      2

We can also specify the cut height:

cutree(hc_single, h = 2)

#> A B C D 
#> 1 1 2 2

Both commands give {A, B} and {C, D}. hclust() builds the full hierarchy; cutree() extracts a partition with a chosen number of groups or a chosen cut height. See the official cutree documentation.

6 Compare the two methods

6.1 Compare memberships

comparison <- data.frame(
  X,
  K_means = km$cluster,
  Single_linkage = hc_groups
)

comparison

#>   x y K_means Single_linkage
#> A 1 1       1              1
#> B 2 1       1              1
#> C 4 3       2              2
#> D 5 4       2              2

table(K_means = km$cluster, Single_linkage = hc_groups)

#>        Single_linkage
#> K_means 1 2
#>       1 2 0
#>       2 0 2

The table shows the same two pairs for this example: A with B, and C with D. The cross-tabulation counts how many observations belong to each combination of clusters.

Cluster numbers are names. The assignments (1, 1, 2, 2) and (2, 2, 1, 1) describe the same grouping. Agreement concerns which observations stay together, even when the numerical labels are different.

6.2 Plot both partitions

old_par <- par(mfrow = c(1, 2), mar = c(4, 4, 3, 1))

plot(
  X$x, X$y,
  col = cluster_colours[km$cluster], pch = 19, cex = 1.6,
  xlim = c(0, 6), ylim = c(0, 5), asp = 1,
  xlab = "x", ylab = "y", main = "K-means"
)
text(X$x, X$y, labels = rownames(X), pos = 3)
points(
  km$centers[, "x"], km$centers[, "y"],
  pch = 4, cex = 1.8, lwd = 2, col = cluster_colours
)

plot(
  X$x, X$y,
  col = cluster_colours[hc_groups], pch = 19, cex = 1.6,
  xlim = c(0, 6), ylim = c(0, 5), asp = 1,
  xlab = "x", ylab = "y", main = "Single linkage: cut at K = 2"
)
text(X$x, X$y, labels = rownames(X), pos = 3)

Both methods identify the same two pairs when two clusters are requested.

par(old_par)

The crosses in the left panel are k-means centres. The right panel shows the partition obtained by cutting the hierarchical tree.

6.3 What is different about the methods?

Feature	K-means used here	Single linkage used here
Starting point	Two supplied centres, A and B	Each observation in its own cluster
Main step	Assign to nearest centre, then update means	Merge the clusters with the smallest single-linkage distance
Distance used in the step	Observation to centre	Closest pair of observations across clusters
Number of clusters	Choose \(K=2\) before fitting	Build the tree, then cut it into two clusters
Main result	Cluster assignments and centres	A hierarchy, represented by a dendrogram

These methods give the same partition for this small example. They can produce different partitions for other data. This example does not establish that two clusters are optimal for every problem.

7 Optional: use the alternative start from the slides

The slides also ask what happens when the starting centres are A and D.

alternative_centres <- as.matrix(X[c("A", "D"), ])

km_AD <- kmeans(
  X,
  centers = alternative_centres,
  algorithm = "Lloyd",
  iter.max = 100
)

km_AD$cluster

#> A B C D 
#> 1 1 2 2

km_AD$centers

#>     x   y
#> 1 1.5 1.0
#> 2 4.5 3.5

km_AD$tot.withinss

#> [1] 1.5

This start produces the same two groups, the same centres, and total WSS 1.5. The initial assignment already gives {A, B} and {C, D}; one centre update reaches the final centres, followed by a check that the assignments remain unchanged. Different starting centres can affect k-means results on other datasets.

8 Sources and reproducibility

• Example data and k-means hand calculations: Adel M., K-means by hand, supplied file KMeans_By_Hand-1.pdf, slides 2-13.
• R function references: kmeans, dist, hclust, and cutree.

All calculations and figures in executable chunks are generated from the four observations entered at the beginning. No external data files are required. Values in the explanatory text refer to this fixed four-point example.

To render from the R console after saving this file in your current working directory:

rmarkdown::render("KMeans_and_Single_Linkage.Rmd")

The following records the R version used when you knit the document:

R.version.string

#> [1] "R version 4.4.1 (2024-06-14 ucrt)"