Activity 4.1 - Hierarchical clustering

SUBMISSION INSTRUCTIONS

Render to html
Publish your html to RPubs
Submit a link to your published solutions

#Loading packages needed for this assignment
library(tidyverse)
library(cluster)
library(factoextra)
library(patchwork)

Question 1

Consider three data sets below. Each data set contains three “clusters” in two dimensions. In the first, the clusters are three convex spheres. (A convex cluster is one where all points in the cluster can be connected with a straight line that does not leave the cluster.) In the second, one cluster is a sphere; one is a ring; and one is a half-moon. In the third, two clusters are spirals and one is a sphere. Our goal is to compare the performance of various hierarchical methods in clustering different cluster shapes.

three_spheres <- read.csv('Data/cluster_data1.csv')
ring_moon_sphere <- read.csv('Data/cluster_data2.csv')
two_spirals_sphere <- read.csv('Data/cluster_data3.csv')

Perform agglomerative clustering with single, complete, average, and Ward linkages. Cut each tree to produce three clusters. Produce 12 scatterplots, one per data set/linkage combination, showing the 3-cluster solution. Title each graph with the linkage used, as well as the average silhouette width (\(\bar s\)) for that clustering solution. Use patchwork to create a nice 4x3 grid of your plots.

#Scaling Data
ts_scaled <- scale(three_spheres)
rms_scaled <- scale(ring_moon_sphere)
tss_scaled <- scale(two_spirals_sphere)

#Three Spheres
ts_single <- agnes(ts_scaled, metric='euclidean', method='single')
ts_complete <- agnes(ts_scaled, metric='euclidean', method='complete')
ts_average <- agnes(ts_scaled, metric='euclidean', method='average')
ts_ward <- agnes(ts_scaled, metric='euclidean', method='ward')

ts_single3 <- cutree(ts_single, k = 3)
ts_complete3 <- cutree(ts_complete, k = 3)
ts_average3 <- cutree(ts_average, k = 3)
ts_ward3 <- cutree(ts_ward, k = 3)

#Three Spheres Silhouette
sil_ts_single <- silhouette(ts_single3, dist(ts_scaled))
mean_sil_ts_single <- mean(sil_ts_single[, "sil_width"])
mean_sil_ts_single

[1] 0.7203813

sil_ts_complete <- silhouette(ts_complete3, dist(ts_scaled))
mean_sil_ts_complete <- mean(sil_ts_complete[, "sil_width"])
mean_sil_ts_complete

[1] 0.7198019

sil_ts_average <- silhouette(ts_average3, dist(ts_scaled))
mean_sil_ts_average <- mean(sil_ts_average[, "sil_width"])
mean_sil_ts_average

[1] 0.7198019

sil_ts_ward <- silhouette(ts_ward3, dist(ts_scaled))
mean_sil_ts_ward <- mean(sil_ts_ward[, "sil_width"])
mean_sil_ts_ward

[1] 0.7203813

#Ring Moon Sphere
rms_single <- agnes(rms_scaled, metric='euclidean', method='single')
rms_complete <- agnes(rms_scaled, metric='euclidean', method='complete')
rms_average <- agnes(rms_scaled, metric='euclidean', method='average')
rms_ward <- agnes(rms_scaled, metric='euclidean', method='ward')

rms_single3 <- cutree(rms_single, k = 3)
rms_complete3 <- cutree(rms_complete, k = 3)
rms_average3 <- cutree(rms_average, k = 3)
rms_ward3 <- cutree(rms_ward, k = 3)

#RMS Silhouette
sil_rms_single <- silhouette(rms_single3, dist(rms_scaled))
mean_sil_rms_single <- mean(sil_rms_single[, "sil_width"])
mean_sil_rms_single

[1] 0.1447351

sil_rms_complete <- silhouette(rms_complete3, dist(rms_scaled))
mean_sil_rms_complete <- mean(sil_rms_complete[, "sil_width"])
mean_sil_rms_complete

[1] 0.4876003

sil_rms_average <- silhouette(rms_average3, dist(rms_scaled))
mean_sil_rms_average <- mean(sil_rms_average[, "sil_width"])
mean_sil_rms_average

[1] 0.5866697

sil_rms_ward <- silhouette(rms_ward3, dist(rms_scaled))
mean_sil_rms_ward <- mean(sil_rms_ward[, "sil_width"])
mean_sil_rms_ward

[1] 0.5859922

#Two Spiral Sphere
tss_single <- agnes(tss_scaled, metric='euclidean', method='single')
tss_complete <- agnes(tss_scaled, metric='euclidean', method='complete')
tss_average <- agnes(tss_scaled, metric='euclidean', method='average')
tss_ward <- agnes(tss_scaled, metric='euclidean', method='ward')

tss_single3 <- cutree(tss_single, k = 3)
tss_complete3 <- cutree(tss_complete, k = 3)
tss_average3 <- cutree(tss_average, k = 3)
tss_ward3 <- cutree(tss_ward, k = 3)

#Two Spirals and Sphere Silhouette
sil_tss_single <- silhouette(tss_single3, dist(tss_scaled))
mean_sil_tss_single <- mean(sil_tss_single[, "sil_width"])
mean_sil_tss_single

[1] 0.1581627

sil_tss_complete <- silhouette(tss_complete3, dist(tss_scaled))
mean_sil_tss_complete <- mean(sil_tss_complete[, "sil_width"])
mean_sil_tss_complete

[1] 0.3856684

sil_tss_average <- silhouette(tss_average3, dist(tss_scaled))
mean_sil_tss_average <- mean(sil_tss_average[, "sil_width"])
mean_sil_tss_average

[1] 0.4224457

sil_tss_ward <- silhouette(tss_ward3, dist(tss_scaled))
mean_sil_tss_ward <- mean(sil_tss_ward[, "sil_width"])
mean_sil_tss_ward

[1] 0.4234588

make_cluster_plot <- function(data, clusters, title_text) {
  df <- data.frame(
    x = data[,1],
    y = data[,2],
    cluster = factor(clusters)
  )
  
  ggplot(df, aes(x, y, color = cluster)) +
    geom_point(size = 1.8, alpha = 0.8) +
    theme_minimal(base_size = 12) +
    ggtitle(title_text) +
    theme(legend.position = "none")
}

title_ts_single   <- paste0("TS:S (s̄ = ", round(mean_sil_ts_single, 3), ")")
title_ts_complete <- paste0("TS:C (s̄ = ", round(mean_sil_ts_complete, 3), ")")
title_ts_average  <- paste0("TS:A (s̄ = ", round(mean_sil_ts_average, 3), ")")
title_ts_ward     <- paste0("TS:W (s̄ = ", round(mean_sil_ts_ward, 3), ")")

title_rms_single   <- paste0("RMS:S (s̄ = ", round(mean_sil_rms_single, 3), ")")
title_rms_complete <- paste0("RMS:C (s̄ = ", round(mean_sil_rms_complete, 3),")")
title_rms_average  <- paste0("RMS:A (s̄ = ", round(mean_sil_rms_average, 3), ")")
title_rms_ward     <- paste0("RMS:W (s̄ = ", round(mean_sil_rms_ward, 3), ")")


title_tss_single   <- paste0("TSS:S (s̄ = ", round(mean_sil_tss_single, 3), ")")
title_tss_complete <- paste0("TSS:C (s̄ = ", round(mean_sil_tss_complete, 3), ")")
title_tss_average  <- paste0("TSS:A (s̄ = ", round(mean_sil_tss_average, 3),")")
title_tss_ward     <- paste0("TSS:W (s̄ = ", round(mean_sil_tss_ward, 3), ")")

p_ts_single   <- make_cluster_plot(ts_scaled, ts_single3, title_ts_single)
p_ts_complete <- make_cluster_plot(ts_scaled, ts_complete3, title_ts_complete)
p_ts_average  <- make_cluster_plot(ts_scaled, ts_average3, title_ts_average)
p_ts_ward     <- make_cluster_plot(ts_scaled, ts_ward3, title_ts_ward)

p_rms_single   <- make_cluster_plot(rms_scaled, rms_single3, title_rms_single)
p_rms_complete <- make_cluster_plot(rms_scaled, rms_complete3, title_rms_complete)
p_rms_average  <- make_cluster_plot(rms_scaled, rms_average3, title_rms_average)
p_rms_ward     <- make_cluster_plot(rms_scaled, rms_ward3, title_rms_ward)

p_tss_single   <- make_cluster_plot(tss_scaled, tss_single3, title_tss_single)
p_tss_complete <- make_cluster_plot(tss_scaled, tss_complete3, title_tss_complete)
p_tss_average  <- make_cluster_plot(tss_scaled, tss_average3, title_tss_average)
p_tss_ward     <- make_cluster_plot(tss_scaled, tss_ward3, title_tss_ward)


(p_ts_single | p_ts_complete | p_ts_average | p_ts_ward) /
(p_rms_single | p_rms_complete | p_rms_average | p_rms_ward) /
(p_tss_single | p_tss_complete | p_tss_average | p_tss_ward)

Discuss the following:

Which linkage works best for which scenario?

Three Spheres: All three are equal
Ring, Moon, Sphere: Average is the best with a average silhouette of 0.587
Two Spiral + Sphere: Ward is the best with a average silhouette of 0.423

Does the average silhouette width always do a good job of measuring the quality of the clustering solution?

No the average silhouette width does not always do a good job at measuring the quality of the clustering because looking at the ring and and the spirals becasue each one of the linkages give very different clustering.

(Hint: you have a lot of repetitive code to write. You may find it helpful to write a function that takes a data set and a linkage method as arguments, does the clustering and computes average silhouette width, and produces the desired plot.)

Question 2

Consider the data set below on milk content of 25 mammals. The variables have been pre-scaled to z-scores, hence no additional standardizing is necessary. (Data source: Everitt et al. Cluster analysis 4ed)

mammals <- read.csv('Data/mammal_milk.csv') %>% 
  column_to_rownames('Mammal')

A)

Perform agglomerative clustering with single, complete, average, and Ward linkages. Which has the best agglomerative coefficient?

single <- agnes(mammals, metric='euclidean', method='single')
complete <- agnes(mammals, metric='euclidean', method='complete')
average <- agnes(mammals, metric='euclidean', method='average')
ward <- agnes(mammals, metric='euclidean', method='ward')

single$ac

[1] 0.7875718

complete$ac

[1] 0.8985539

average$ac

[1] 0.8706571

ward$ac

[1] 0.9413994

The ward linkage has the best agglomerative coefficient.

B)

Plot a dendrogram of the method with the highest AC. Which mammals cluster together first?

fviz_dend(ward) + 
  labs(title='Ward linkage')

Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
ℹ Please use tidy evaluation idioms with `aes()`.
ℹ See also `vignette("ggplot2-in-packages")` for more information.
ℹ The deprecated feature was likely used in the factoextra package.
  Please report the issue at <https://github.com/kassambara/factoextra/issues>.

Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.
ℹ The deprecated feature was likely used in the factoextra package.
  Please report the issue at <https://github.com/kassambara/factoextra/issues>.

Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
of ggplot2 3.3.4.
ℹ The deprecated feature was likely used in the factoextra package.
  Please report the issue at <https://github.com/kassambara/factoextra/issues>.

Fox and Sheep, donkey and mule, dog and rat, deer and reindeer, camel and zebra and a couple of the mamals that are first clustered together

C)

If the tree is cut at a height of 4, how many clusters will form? Which cluster will have the fewest mammals, and which mammals will they be?

fviz_dend(ward, h = 4) + 
  labs(title='Ward linkage')

There is 4 clusters when the height is 4, the purple cluster has the fewest mammals with only dolphin and seal in the cluster.

D)

Use WSS and average silhouette method to suggest the optimal number of clusters. Re-create the dendrogram with the cluster memberships indicated.

fviz_nbclust(mammals,
             FUNcluster = hcut,
             k.max = 10,
             method='wss',
             hc_method='ward',
             hc_func = 'agnes')

fviz_nbclust(mammals,
             FUNcluster = hcut,
             k.max = 10,
             method='silhouette',
             hc_method='ward',
             hc_func = 'agnes')

Both the WSS and the average silhouette seem to suggest the optimal number of clusters is 3 for the mammals data set

fviz_dend(ward, k = 3) + 
  labs(title='Ward linkage')

E)

Use suitable visualizations, including dimension reduction techniques, to explore the different milk characteristics of the assigned clusters. Discuss.

(ward_clusters <- cutree(ward, k = 3) 
    %>% data.frame()
    %>% setNames('k3')
) %>% head

mammal_pca <- prcomp(mammals, center=TRUE, scale. = TRUE)

k2_biplot <- fviz_pca(mammal_pca, 
         habillage = factor(ward_clusters$k3),
         repel = TRUE) + 
      ggtitle('3-cluster solution') + 
      guides(color='none',shape='none')
k2_biplot

Looking at the PCA bi plot we can see the distinct regions that these mammals are in and what variables they have large amounts of. The green group with rabbit, whale, rat, and others have large amounts of Ash and Protein in their milk. The blue group with only seal and dolphin have large amounts of fat in their milk. Then finally the red grouping with the most amount of mammals within it have large amounts of water and lactose within their milk.