Machine Leanring Assignment 2: Analysis of Clustering K-Means, DIANA, AGNES, SOM
Sushant Gote 17PHD1054
This R Notebook demonstrates application ?f different Clusterring algorithms as Kmeans, K-medoids, Diana, Ages, SOM on USArrests data availed by R.
`` # Install packages
install.packages("cluster")
install.packages("dendextend")
install.packages('factoextra')
#Load the packages :?``{r} ### Installling addtional Packages if required
install.packages(‘kohonen’) library(‘kohonen’) ``` # Kmeans Clustering Clustering is a broad set of techniques for finding subgroups of observations within a data set. When we cluster observations, we w?nt observations in the same group to be similar and observations in different groups to be dissimilar. Because there isn’t a response variable, this is an unsupervised method, which implies that it seeks to find relationships between the n observations w?thout being trained by a response variable. Clustering allows us to identify which observations are alike, and potentially categorize them therein. K-means clustering is the simplest and the most commonly used clustering method for splitting a dataset into?a set of k groups.
K-means algorithm can be summarized as follows:
Specify the number of clusters (K) to be created (by the analyst) Select randomly k objects from the data set as the initial cluster centers or means Assigns each observation to their ?losest centroid, based on the Euclidean distance between the object and the centroid For each of the k clusters update the cluster centroid by calculating the new mean values of all the data points in the cluster. The centroid of a Kth cluster is a vector ?f length p containing the means of all variables for the observations in the kth cluster; p is the number of variables. Iteratively minimize the total within sum of square (Eq. 7). That is, iterate steps 3 and 4 until the cluster assignments stop changing ?r the maximum number of iterations is reached. By default, the R software uses 10 as the default value for the maximum number of iterations.
Optimal Number of clusters
Algorithm to define the optimal clusters:
Compute clustering algorithm (e.?., k-means clustering) for different values of k. For instance, by varying k from 1 to 10 clusters
For each k, calculate the total within-cluster sum of square (wss)
Plot the curve of wss according to the number of clusters k.
The location o? a bend (knee) in the plot is generally considered as an indicator of the appropriate number of clusters.
The results suggest that 4 is the optimal number of clusters as it appears to be the bend in the knee (or elbow).
set.seed(123)
# functi?n to compute total within-cluster sum of square
wss <- function(k) {
kmeans(df, k, nstart = 10 )$tot.withinss
}
# Compute and plot wss for k = 1 to k = 15
k.values <- 1:15
# extract wss for 2-15 clusters
wss_values <- map_dbl(k.values, wss)
plot(k.va?ues, wss_values,
type="b", pch = 19, frame = FALSE,
xlab="Number of clusters K",
ylab="Total within-clusters sum of squares")set.seed(123)
# function to compute total within-cluster sum of square
wss <- function(k) {
kmeans(df, k, nstart = 10 )$tot.withinss
}
# Compute and plot wss for k = 1 to k = 15
k.values <- 1:15
# extract wss for 2-15 clusters
wss_values <- map_dbl(k.values, wss)
plot(k.values, wss_values,
type="b", pch = 19, frame = FALSE,
xlab="Number of clusters K",
ylab="Total within-clusters sum of squares")
Apply?ng Kmeans for K=4
library('cluster')
library('dendextend')
library('factoextra')
print('Done')[1] "Done" Visualizing Kmeans Clusters
set.seed(123)
final <- kmeans(df, 4, nstart = 25)
print(final)K-means clustering with 4 clusters of sizes 13, 16, 13, 8
Cluster means:
Murder Assault UrbanPop Rape
1 -0.9615407 -1.1066010 -0.9301069 -0.96676331
2 -0.4894375 -0.3826001 0.5758298 -0.26165379
3 0.6950701 1.0394414 0.7226370 1.27693964
4 1.4118898 0.8743346 -0.8145211 0.01927104
Clustering vector:
Alabama Alaska Arizona Arkansas
4 3 3 4
California Colorado Connecticut Delaware
3 3 2 2
Florida Georgia Hawaii Idaho
3 4 2 1
Illinois Indiana Iowa Kansas
3 2 1 2
Kentucky Louisiana Maine Maryland
1 4 1 3
Massachusetts Michigan Minnesota Mississippi
2 3 1 4
Missouri Montana Nebraska Nevada
3 1 1 3
New Hampshire New Jersey New Mexico New York
1 2 3 3
North Carolina North Dakota Ohio Oklahoma
4 1 2 2
Oregon Pennsylvania Rhode Island South Carolina
2 2 2 4
South Dakota Tennessee Texas Utah
1 4 3 2
Vermont Virginia Washington West Virginia
1 2 2 1
Wisconsin Wyoming
1 2
Within cluster sum of squares by cluster:
[1] 11.952463 16.212213 19.922437 8.316061
(between_SS / total_SS = 71.2 %)
Available components:
[1] "cluster" "centers" "totss" "withinss"
[5] "tot.withinss" "betweenss" "size" "iter"
[9] "ifault" K Medoid Clustering
fviz_cluster(final, data = df)
Vi?ualizing K-medioids Clusters
NAHierachical Clustering
Hierarchical clustering can be divided into two main types: agglomerative and divisive.
Agglomerative clustering: It’s also known as AGNES (Agglomerative N?sting). It works in a bottom-up manner. That is, each object is initially considered as a single-element cluster (leaf). At each step of the algorithm, the two clusters that are the most similar are combined into a new bigger cluster (nodes). This procedur? is iterated until all points are member of just one single big cluster (root) (see figure below). The result is a tree which can be plotted as a dendrogram.
Divisive hierarchical clustering: It’s also known as DIANA (Divise Analysis) and it works in a to?-down manner. The algorithm is an inverse order of AGNES. It begins with the root, in which all objects are included in a single cluster. At each step of iteration, the most heterogeneous cluster is divided into two. The process is iterated until all objec?s are in their own cluster .
Note that agglomerative clustering is good at identifying small clusters. Divisive hierarchical clustering is good at identifying large clusters.
The merging or the division of clusters is performed according some (dis)simil?rity measure. In R softwrare, the Euclidean distance is used by default to measure the dissimilarity between each pair of observations.
Tt’s easy to compute dissimilarity measure between two pairs of observations. It’s mentioned above that two clusters th?t are most similar are fused into a new big cluster.
A natural question is : How to measure the dissimilarity between two clusters of observations?
A number of different cluster agglomeration methods (i.e, linkage methods) has been developed to answer to?this question. The most common types methods are:
Maximum or complete linkage clustering: It computes all pairwise dissimilarities between the elements in cluster 1 and the elements in cluster 2, and considers the largest value (i.e., maximum value) of t?ese dissimilarities as the distance between the two clusters. It tends to produce more compact clusters.
Minimum or single linkage clustering: It computes all pairwise dissimilarities between the elements in cluster 1 and the elements in cluster 2, and c?nsiders the smallest of these dissimilarities as a linkage criterion. It tends to produce long, “loose” clusters.
Mean or average linkage clustering: It computes all pairwise dissimilarities between the elements in cluster 1 and the elements in cluster 2,?and considers the average of these dissimilarities as the distance between the two clusters.
Complete linkage and Ward’s method are generally preferred.
Centroid linkage clustering: It computes the dissimilarity between the centroid for cluster 1 (a mean?vector of length p variables) and the centroid for cluster 2.
Ward’s minimum variance method: It minimizes the total within-cluster variance. At each step the pair of clusters with minimum between-cluster distance are merged.
Data preparation and descr?ptive statistics
R dataset USArrest which contains statistics, in arrests per 100,000 residents for assault, murder, and rape in each of the 50 US states in 1973. It includes also the percent of the population living in urban areas.
It contains 50 observa?ions on 4 variables:
[,1] Murder numeric Murder arrests (per 100,000) [,2] Assault numeric Assault arrests (per 100,000) [,3] UrbanPop numeric Percent urban population [,4] Rape numeric Rape arrests (per 100,000)
fviz_cluster(pam.res, data = df)
Before hierarchical clustering, compute some descriptive statistics?
# Load the data set
data("USArrests")
# Remove any missing value (i.e, NA values for not available)
# That might be present in the data
df <- na.omit(USArrests)
# View the firt 6 rows of the data
head(df, n = 6)Note that the variables have a large different means and variances. This is explained by the fact that the variables are measured in different units; Murder, Rape, and Assault are measured as the number of occurrences per?100 000 people, and UrbanPop is the percentage of the state’s population that lives in an urban area.
They must be standardized (i.e., scaled) to make them comparable. Recall that, standardization consists of transforming the variables such that they have?mean zero and standard deviation one.
desc_stats <- data.frame(
Min = apply(df, 2, min), # minimum
Med = apply(df, 2, median), # median
Mean = apply(df, 2, mean), # mean
SD = apply(df, 2, sd), # Standard deviation
Max = apply(df, 2, max) # Maximum
)
desc_stats <- round(desc_stats, 1)
head(desc_stats)R functions for hierarchical clustering
There are different functions available in R for computing hierarchical clustering. The commonly used function? are:
hclust() [in stats package] and agnes() [in cluster package] for agglomerative hierarchical clustering (HC) diana() [in cluster package] for divisive HC #### hclust() function hclust() is the built-in R function [in stats package] for computing hier?rchical clustering.
The simplified format is:
hclust(d, method = “complete”)
d a dissimilarity structure as produced by the dist() function. method: The agglomeration method to be used. Allowed values is one of “ward.D”, “ward.D2”, “single”, “complete”,?“average”, “mcquitty”, “median” or “centroid”.
The dist() function is used to compute the Euclidean distance between observations. Finally, observations are clustered using Ward’s method.
df <- scale(df)
head(df) Murder Assault UrbanPop Rape
Alabama 1.24256408 0.7828393 -0.5209066 -0.003416473
Alaska 0.50786248 1.1068225 -1.2117642 2.484202941
Arizona 0.07163341 1.4788032 0.9989801 1.042878388
Arkansas 0.23234938 0.2308680 -1.0735927 -0.184916602
California 0.27826823 1.2628144 1.7589234 2.067820292
Colorado 0.02571456 0.3988593 0.8608085 1.864967207print('Look @ this head')[1] "Look @ this head"agnes() and diana() functions
The R function agnes() [in cluster package] can be also used to ?ompute agglomerative hierarchical clustering. The R function diana() [ in cluster package ] is an example of divisive hierarchical clustering.
Agglomerative Nesting (Hierarchical Clustering)
agnes(x, metric = “euclidean”, stand = FALSE, method = “aver?ge”)
DIvisive ANAlysis Clustering
diana(x, metric = “euclidean”, stand = FALSE)
x: data matrix or data frame or dissimilarity matrix. In case of matrix and data frame, rows are observations and columns are variables. In case of a dissimilarity matrix? x is typically the output of daisy() or dist(). metric: the metric to be used for calculating dissimilarities between observations. Possible values are “euclidean” and “manhattan”. stand: if TRUE, then the measurements in x are standardized before calcula?ing the dissimilarities. Measurements are standardized for each variable (column), by subtracting the variable’s mean value and dividing by the variable’s mean absolute deviation method: The clustering method. Possible values includes “average”, “single”, ???complete“,”ward".
The function agnes() returns an object of class “agnes” (see ?agnes.object) which has methods for the functions: print(), summary(), plot(), pltree(), as.dendrogram(), as.hclust() and cutree(). The function diana() returns an object of?class “diana” (see ?diana.object) which has also methods for the functions: print(), summary(), plot(), pltree(), as.dendrogram(), as.hclust() and cutree(). Compared to other agglomerative clustering methods such as hclust(), agnes() has the following feat?res:
It yields the agglomerative coefficient (see agnes.object) which measures the amount of clustering structure found Apart from the usual tree it also provides the banner, a novel graphical display (see plot.agnes).
# Dissimilarity matrix
d <- dist(df, method = "euclidean")
# Hierarchical clustering using Ward's method
res.hc <- hclust(d, method = "ward.D2" )
# Plot the obtained dendrogram
plot(res.hc, cex = 0.6, hang = -1)
library("cluster")
# Compute agnes()
res.agnes <- agnes(df, method = "ward")
# Agglomerative coefficient
res.agnes$ac[1] 0.934621# Plot the tree using pltree()
pltree(res.agnes, cex = 0.6, hang = -1, main = "Dendrogram of Agnes") 
print('Done')[1] "Done"It’s also possible to draw AGNES dendrogram using the function plot.hclust() and the function plot.dendrogram() as follow:
pltree(res.agnes, cex = 0.6, hang = -1, main = "Dendrogram of Agnes")
plot.?endrogram()
# plot.hclust()
plot(as.hclust(res.agnes), cex = 0.6, hang = -1)
# plot.dendrogram()
plot(as.dendrogram(res.agnes), cex = 0.6,
horiz = TRUE)
# Cut tree into 4 groups
grp <- cutree(res.hc, k = 4)
# Number of members in each cluster
table(grp)grp
1 2 3 4
7 12 19 12 plot(res.hc, cex = 0.6)
rect.hclust(res.hc, k = 4, border = 2:5)
Diana Clustering for USArrest dataset
library(factoextra)
fviz_cluster(list(data = df, cluster = grp))
# Compute diana()
res.diana <- diana(df)
print('Diana Computed')[1] "Diana Computed"# Plot the tree
pltree(res.diana, cex = 0.6, hang = -1,
main = "Dendrogram of Diana")
# Plot the tree
pltree(res.diana, cex = 0.6, hang = -1,
main = "Second time Dendrogram of Diana")
As for plotting AGNES dendrogram, the functions plot.hclust() and plot.dendrogram() can be used as follow:
# Divise coefficient; amount of clustering structure found
res.diana$dc[1] 0.8514345## [1] 0.8514345pl?t.dendrogram()
# plot.hclust()
plot(as.hclust(res.diana), cex = 0.6, hang = -1)
Interpretation of the dendrogram
In the dendrogram displayed above, each leaf corresponds to one observation. moving up the tree, observations that are similar to e?ch other are combined into branches, which are themselves fused at a higher height.
The height of the fusion, provided on the vertical axis, indicates the (dis)similarity between two observations. The higher the height of the fusion, the less similar the ?bservations are.
Note that, conclusions about the proximity of two observations can be drawn only based on the height where branches containing those two observations first are fused. The proximity of two observations along the horizontal axis can not be ?sed as a criteria of their similarity.
In order to identify sub-groups (i.e. clusters), Cuttingt the dendrogram at a certain height as described in the next section.
Cut the dendrogram into different groups
The height of the cut to the dendrogram co?trols the number of clusters obtained. It plays the same role as the k in k-means clustering.
The function cutree() is used and it returns a vector containing the cluster number of each observation
plot(as.dendrogram(res.diana), cex = 0.6,
horiz = TRUE)
It’s also possible to draw the dendrogram with a border around the 4 clusters. The argument border is used to specify the border ?olors for the rectangles:
# Cut tree into 4 groups
grp <- cutree(res.hc, k = 4)
# Number of members in each cluster
table(grp)grp
1 2 3 4
7 12 19 12
# Get the names for the members of cluster 1
rownames(df)[grp == 1][1] "Alabama" "Georgia" "Louisiana" "Mississippi" "North Carolina"
[6] "South Carolina" "Tennessee" Look at leaves
Observe clusters
plot(res.hc, cex = 0.6)
rect.hclust(res.hc, k = 4, border = 2:5)
Using the function fviz_cluster() [in factoextra], visualizing the result in a scatter plot. Observations are represented by points in the plot, using principal components. A frame is drawn around?each cluster.
plot(res.hc, cex = 0.6)
rect.hclust(res.hc, k = 4, border = 2:5)
plot(res.hc, cex = 0.6)
rect.hclust(res.hc, k = 4, border = 2:5)
library(factoextra)
fviz_cluster(list(data = df, cluster = grp))
The function cutree() c?n be used also to cut the tree generated with agnes() and diana() as follow:
library(factoextra)
fviz_cluster(list(data = df, cluster = grp))
# Cut diana() tree into 4 groups
cutree(as.hclust(res.diana), k = 4) Alabama Alaska Arizona Arkansas California Colorado
1 2 2 3 2 2
Connecticut Delaware Florida Georgia Hawaii Idaho
3 3 2 1 3 4
Illinois Indiana Iowa Kansas Kentucky Louisiana
2 3 4 3 4 1
Maine Maryland Massachusetts Michigan Minnesota Mississippi
4 2 3 2 4 1
Missouri Montana Nebraska Nevada New Hampshire New Jersey
2 4 4 2 4 3
New Mexico New York North Carolina North Dakota Ohio Oklahoma
2 2 1 4 3 3
Oregon Pennsylvania Rhode Island South Carolina South Dakota Tennessee
3 3 3 1 4 1
Texas Utah Vermont Virginia Washington West Virginia
2 3 4 3 3 4
Wisconsin Wyoming
4 3 Self Organizing Map SOM Clustering
# Cut diana() tree into 4 groups
cutree(as.hclust(res.diana), k = 4) Alabama Alaska Arizona Arkansas California Colorado
1 2 2 3 2 2
Connecticut Delaware Florida Georgia Hawaii Idaho
3 3 2 1 3 4
Illinois Indiana Iowa Kansas Kentucky Louisiana
2 3 4 3 4 1
Maine Maryland Massachusetts Michigan Minnesota Mississippi
4 2 3 2 4 1
Missouri Montana Nebraska Nevada New Hampshire New Jersey
2 4 4 2 4 3
New Mexico New York North Carolina North Dakota Ohio Oklahoma
2 2 1 4 3 3
Oregon Pennsylvania Rhode Island South Carolina South Dakota Tennessee
3 3 3 1 4 1
Texas Utah Vermont Virginia Washington West Virginia
2 3 4 3 3 4
Wisconsin Wyoming
4 3 Training progress for SOM
plot(som_model, type="changes")Node Counts
The Kohonen packages allows us to visualise the count of how many samples are mapped to each node on the map. This metric can be used as a measure of map quality - ideally the sample distribution is re?atively uniform. Large values in some map areas suggests that a larger map would be benificial. Empty nodes indicate that your map size is too big for the number of samples. Aim for at least 5-10 samples per node when choosing map size.
# Change the data frame with training data to a matrix
# Also center and scale all variables to give them equal importance during
# the SOM training process.
data_train_matrix <- as.matrix(scale(df))
# Create the SOM Grid - you generally have to specify the size of the
# training grid prior to training the SOM. Hexagonal and Circular
# topologies are possible
som_grid <- somgrid(xdim = 20, ydim=20, topo="hexagonal")
# Finally, train the SOM, options for the number of iterations,
# the learning rates, and the neighbourhood are available
som_model <- som(data_train_matrix, somgrid(5, 5, "hexagonal"))plot(som_model, type="changes")
Neighbour Distance
Often referred to as the “U-Matrix”, this visualisati?n is of the distance between each node and its neighbours. Typically viewed with a grayscale palette, areas of low neighbour distance indicate groups of nodes that are similar. Areas with large distances indicate the nodes are much more dissimilar - and in?icate natural boundaries between node clusters. The U-Matrix can be used to identify clusters within the SOM map. #### U-matrix visualisation plot(som_model, type=“dist.neighbours”, main = “SOM neighbour distances”)
# Node count plot
plot(som_model, type="count", main="Node Counts")
print('many samples are mapped to each node on the map. This metric can be used as a measure of map quality')[1] "many samples are mapped to each node on the map. This metric can be used as a measure of map quality"Nodes / Weight vectors
The node weight vectors, or “codes”, are made up of normalised values of the original variables used to generate th? SOM. Each node’s weight vector is representative / similar of the samples mapped to that node. By visualising the weight vectors across the map, analyzing patterns in the distribution of samples and variables. The default visualisation of the weight vecto?s is a “fan diagram”, where individual fan representations of the magnitude of each variable in the weight vector is shown for each node.
#### U-matrix visualisation
plot(som_model, type="dist.neighbours", main = "SOM neighbour distances")
print('SOM Neighbour Distance')[1] "SOM Neighbour Distance"Heatmaps
H?atmaps are perhaps the most important visualisation possible for Self-Organising Maps. The use of a weight space view as in (4) that tries to view all dimensions on the one diagram is unsuitable for a high-dimensional (>7 variable) SOM. A SOM heatmap allow? the visualisation of the distribution of a single variable across the map. Typically, a SOM investigative process involves the creation of multiple heatmaps, and then the comparison of these heatmaps to identify interesting areas on the map. It is importa?t to remember that the individual sample positions do not move from one visualisation to another, the map is simply coloured by different variables.
# Weight Vector View
plot(som_model, type="codes")
print('Weight vector View')[1] "Weight vector View"# Kohonen Heatmap creation
plot(som_model, type = "property", property = getCodes(som_model)[,4], main=colnames(getCodes(som_model))[4], palette.name=coolBlueHotRed)
plot(som_model, type = "property", property = getCodes(som_model)[,4], main=colnames(getCodes(som_model))[3], palette.name=coolBlueHotRed)
plot(som_model, type = "property", property = getCodes(som_model)[,4], main=colnames(getCodes(som_model))[2], palette.name=coolBlueHotRed)
A more intuitive and useful visualisation is of the variable prior to scaling, which involves some R trickery - using the aggregate function to regenerate the variable from the original training set and th? SOM node/sample mappings. The result is scaled to the real values of the training variable (in this case, unemployment percent).
print('Kohenen HeatMap')[1] "Kohenen HeatMap"Clustering
Clustering can be performed on the SOM nodes to i?olate groups of samples with similar metrics. Manual identification of clusters is completed by exploring the heatmaps for a number of variables and drawing up a “story” about the different areas on the map. An estimate of the number of clusters that would?be suitable can be ascertained using a kmeans algorithm and examing for an “elbow-point” in the plot of “within cluster sum of squares”. The Kohonen package documentation shows how a map can be clustered using hierachical clustering. The results of the cl?stering can be visualised using the SOM plot function again.
var <- 3 #define the variable to plot
var_unscaled <- aggregate(as.numeric(my_data_matrix[,var]), by=list(som_model$unit.classif), FUN=mean, simplify=TRUE)[,2] Conclusion
. The cluster formed by Alabama, Alaska, Arizona, California, Dela?are, Florida, Illinois, Louisiana, Maryland, Michigan, Mississippi, Nevada, New Mexico, New York, North Carolina, South Carolina has the highest Murder, Assault and Rape arests (per 100,00) and, not least, the largest population.
. The cluster formed by A?kansas, Colorado, Georgia, Massachusetts, Missouri, New Jersey, Oklahoma, Oregon, Rhode Island, Tennessee, Texas, Virginia, Washington, Wyoming has the intermediate Murder, Assault and Rape arests (per 100,00) and, not least, the largest population.
. Th? cluster formed by Connecticut, Hawaii, Idaho, Indiana, Iowa, Kansas, Kentucky, Maine, Minnesota, Montana, Nebraska, New Hampshire, North Dakota, Ohio, Pennsylvania, South Dakota, Utah, Vermont, West Virginia, Wisconsin has the lowest Murder, Assault and R?pe arests (per 100,00) and, not least, the largest population.
Remarks for Alogorithms
Kmeans , Kedoid, Hierarchical
K-means clustering is a very simple and fast algorithm. Furthermore, it can efficiently deal with very large data sets? However, there are some weaknesses of the k-means approach.
One potential disadvantage of K-means clustering is that it requires us to pre-specify the number of clusters. Hierarchical clustering is an alternative approach which does not require that we c?mmit to a particular choice of clusters. Hierarchical clustering has an added advantage over K-means clustering in that it results in an attractive tree-based representation of the observations, called a dendrogram.
An additional disadvantage of K-means?is that it’s sensitive to outliers and different results can occur if you change the ordering of your data. The Partitioning Around Medoids (PAM) clustering approach is less sensititive to outliers and provides a robust alternative to k-means to deal with ?hese situations.
References
- https://cran.r-project.org/web/packages/kohonen/kohonen.pdf
- https://www.superdatascience.com/blogs/the-ultimate-guide-to-self-organizing-maps-soms
- https://en.wikibooks.org/wiki/Data_Mini?g_Algorithms_In_R/Clustering/Self-Organizing_Maps_(SOM)#Introduction
- https://www.webpages.uidaho.edu/~stevel/517/Data%20Mining%20Algorithms%20In%20R.pdf
- http://girke.bioinformatics.ucr.edu/GEN242/pages/mydoc/Rclustering.html
- https://clark?atalabs.github.io/soms/SOM_NBA
- https://www.r-bloggers.com/self-organising-maps-for-customer-segmentation-using-r/
- http://www.sthda.com/english/wiki/
- https://uc-r.github.io/kmeans_clustering
- https://www.shanelynn.ie/self-organising-ma?s-for-customer-segmentation-using-r/
- https://www.datanovia.com/en/courses/partitional-clustering-in-r-the-essentials/