Кластеризация - одна из наиболее распространенных проблем машинного обучения без учителя. Сходство между наблюдениями определяется с помощью некоторых мер расстояния между наблюдениями или мер расстояния на основе корреляции.
Существует 5 классов методов кластеризации:
+ Иерархическая кластеризация
+ Методы разделения (k-means, PAM, CLARA)
+ Кластеризация на основе плотности
+ Кластеризация на основе модели
+ Нечеткая кластеризация
Мое желание написать этот пост возникло в основном из-за того, что я прочитал о пакете clustree, документации dendextend и книге Практическое руководство по кластерному анализу в R. Альбукадель Кассамбара, автор пакета factoextra.
Набор данных
Я буду использовать менее известный набор данных из пакета cluster: all.mammals.milk.1956, который я раньше не рассматривал.
Этот небольшой набор данных содержит список из 25 млекопитающих и компонентов их молока (процентное содержание воды, белка, жира, лактозы, золы) из John Hartigan, Clustering Algorithms, Wiley, 1975.
Сначала загрузим необходимые пакеты.
library(tidyverse) library(magrittr) library(cluster) library(cluster.datasets) library(cowplot) library(NbClust) library(clValid) library(ggfortify) library(clustree) library(dendextend) library(factoextra) library(FactoMineR) library(corrplot) library(GGally) library(ggiraphExtra) library(knitr) library(kableExtra)
Теперь загрузите данные.
data("all.mammals.milk.1956")
raw_mammals <- all.mammals.milk.1956
# subset dataset
mammals <- raw_mammals %>% select(-name) # set rownames
mammals <- as_tibble(mammals)
Давайте изучим и визуализируем данные.
# Glimpse the data set glimpse(mammals) Observations: 25 Variables: 5 $ water <dbl> 90.1, 88.5, 88.4, 90.3, 90.4, 87.7, 86.9, 82.1, 81.9, 81.6, 81.6, 86.5, 90.0,... $ protein <dbl> 2.6, 1.4, 2.2, 1.7, 0.6, 3.5, 4.8, 5.9, 7.4, 10.1, 6.6, 3.9, 2.0, 7.1, 3.0, 5... $ fat <dbl> 1.0, 3.5, 2.7, 1.4, 4.5, 3.4, 1.7, 7.9, 7.2, 6.3, 5.9, 3.2, 1.8, 5.1, 4.8, 6.... $ lactose <dbl> 6.9, 6.0, 6.4, 6.2, 4.4, 4.8, 5.7, 4.7, 2.7, 4.4, 4.9, 5.6, 5.5, 3.7, 5.3, 4.... $ ash <dbl> 0.35, 0.24, 0.18, 0.40, 0.10, 0.71, 0.90, 0.78, 0.85, 0.75, 0.93, 0.80, 0.47,...
Все переменные выражаются числовыми значениями. А как насчет статистического распределения?
# Summary of data set summary(mammals) %>% kable() %>% kable_styling()
# Historgram for each attribute mammals %>% gather(Attributes, value, 1:5) %>% ggplot(aes(x=value)) + geom_histogram(fill = "lightblue2", color = "black") + facet_wrap(~Attributes, scales = "free_x") + labs(x = "Value", y = "Frequency")
Какая связь между различными атрибутами? Используйте corrplot () для создания корреляционной матрицы.
corrplot(cor(mammals), type = "upper", method = "ellipse", tl.cex = 0.9)
Когда у вас есть переменные, которые измеряются в разных масштабах, полезно масштабировать данные.
mammals_scaled <- scale(mammals) rownames(mammals_scaled) <- raw_mammals$name
Снижение размерности может помочь в визуализации данных (например, метод PCA).
res.pca <- PCA(mammals_scaled, graph = FALSE) # Visualize eigenvalues/variances fviz_screeplot(res.pca, addlabels = TRUE, ylim = c(0, 50))
Это 5 компьютеров, на которые приходится 80% дисперсии. График осыпи показывает, что PC1 захватил ~ 75% дисперсии.
# Extract the results for variables
var <- get_pca_var(res.pca)
# Contributions of variables to PC1
fviz_contrib(res.pca, choice = "var", axes = 1, top = 10)
# Contributions of variables to PC2
fviz_contrib(res.pca, choice = "var", axes = 2, top = 10)
# Control variable colors using their contributions to the principle axis
fviz_pca_var(res.pca, col.var="contrib",
gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
repel = TRUE # Avoid text overlapping
) + theme_minimal() + ggtitle("Variables - PCA")
Из этих визуализаций очевидно, что вода и лактоза имеют тенденцию увеличиваться вместе, а белок, зола и жир увеличиваются вместе; две группы обратно связаны.
Наивный (средний) подход
Методы секционирования кластеризации, такие как k-means и Partitioning Around Medoids (PAM), требуют, чтобы вы указали количество кластеров, которые должны быть сгенерированы.
Кластеры k-средних - это, вероятно, один из самых известных методов разделения. Идея кластеризации k-средних состоит в том, чтобы определить кластеры, в которых общая вариация внутри кластера, которая измеряет компактность кластеров, сводится к минимуму.
Мы можем вычислить k-среднее в R с помощью функции kmeans ():
km2 <- kmeans(mammals_scaled, centers = 2, nstart = 30)
В приведенном выше примере данные будут сгруппированы в два кластера, center = 2, и будет предпринята попытка нескольких начальных конфигураций, в результате чего будет получена наилучшая из них. Например, поскольку этот алгоритм чувствителен к начальным положениям центроидов кластера, добавление nstart = 30 сгенерирует 30 начальных конфигураций, а затем усреднит все результаты центроидов.
Поскольку количество кластеров (k) необходимо установить перед тем, как мы начнем, может быть полезно изучить несколько различных значений k.
kmean_calc <- function(df, ...){
kmeans(df, scaled = ..., nstart = 30)
}
km2 <- kmean_calc(mammals_scaled, 2)
km3 <- kmean_calc(mammals_scaled, 3)
km4 <- kmeans(mammals_scaled, 4)
km5 <- kmeans(mammals_scaled, 5)
km6 <- kmeans(mammals_scaled, 6)
km7 <- kmeans(mammals_scaled, 7)
km8 <- kmeans(mammals_scaled, 8)
km9 <- kmeans(mammals_scaled, 9)
km10 <- kmeans(mammals_scaled, 10)
km11 <- kmeans(mammals_scaled, 11)
p1 <- fviz_cluster(km2, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 2")
p2 <- fviz_cluster(km3, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 3")
p3 <- fviz_cluster(km4, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 4")
p4 <- fviz_cluster(km5, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 5")
p5 <- fviz_cluster(km6, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 6")
p6 <- fviz_cluster(km7, data = mammals_scaled, frame.type = "convex") + theme_minimal() + ggtitle("k = 7")
plot_grid(p1, p2, p3, p4, p5, p6, labels = c("k2", "k3", "k4", "k5", "k6", "k7"))
Хотя эта визуальная оценка говорит нам, где происходит разграничение между кластерами, она не говорит нам об оптимальном количестве кластеров.
Определение оптимального количества кластеров
В литературе предлагалось множество мер для оценки результатов кластеризации. Термин проверка кластеризации используется для разработки процедуры оценки результатов алгоритма кластеризации. Существует более тридцати индексов и методов для определения оптимального количества кластеров, поэтому я остановлюсь только на некоторых из них, включая очень удобный пакет clustree.
Метод «локтя»
Вероятно, наиболее известный метод, метод изгиба, в котором сумма квадратов для каждого количества кластеров вычисляется и отображается на графике, и пользователь ищет изменение наклона с крутого на пологий (изгиб), чтобы определить оптимальное количество кластеры. Этот метод неточен, но все же потенциально полезен.
set.seed(31)
# function to compute total within-cluster sum of squares
fviz_nbclust(mammals_scaled, kmeans, method = "wss", k.max = 24) + theme_minimal() + ggtitle("the Elbow Method")
Метод кривой локтя полезен, потому что он показывает, как увеличение числа кластеров способствует разделению кластеров значимым, а не маргинальным образом. Изгиб указывает на то, что дополнительные кластеры помимо третьего имеют небольшую ценность (см. [Здесь] для более математически строгой интерпретации и реализации этого метода). Метод локтя довольно ясный, если не наивное решение, основанное на внутрикластерной дисперсии. Статистика разрыва - более сложный метод работы с данными, имеющими распределение без явной кластеризации (можно найти правильное число k для глобулярных, распределенных по Гауссу, слегка несвязанных распределений данных).
Статистика разрыва
Статистика разрыва сравнивает общую вариацию внутри кластера для различных значений k с их ожидаемыми значениями при нулевом эталонном распределении данных. Оценка оптимальных кластеров будет представлять собой значение, которое максимизирует статистику разрыва (т.е., которая дает статистику наибольшего разрыва). Это означает, что структура кластеризации далека от случайного равномерного распределения точек.
gap_stat <- clusGap(mammals_scaled, FUN = kmeans, nstart = 30, K.max = 24, B = 50)
fviz_gap_stat(gap_stat) + theme_minimal() + ggtitle("fviz_gap_stat: Gap Statistic")
График статистики пробелов показывает статистику по количеству кластеров (k) со стандартными ошибками, нарисованными вертикальными сегментами, и оптимальным значением k, отмеченным вертикальной пунктирной синей линией. Согласно этому наблюдению k = 2 - оптимальное количество кластеров в данных.
Силуэтный метод
Другая визуализация, которая может помочь определить оптимальное количество кластеров, называется методом силуэта. Метод среднего силуэта вычисляет средний силуэт наблюдений для различных значений k. Оптимальное количество кластеров k - это то, которое максимизирует средний силуэт в диапазоне возможных значений для k.
fviz_nbclust(mammals_scaled, kmeans, method = "silhouette", k.max = 24) + theme_minimal() + ggtitle("The Silhouette Plot")
Это также предполагает оптимальный вариант из 2 кластеров.
Метод суммы квадратов
Другой метод проверки кластеризации заключается в выборе оптимального количества кластеров путем минимизации суммы квадратов внутри кластера (мера того, насколько плотен каждый кластер) и максимизации суммы квадратов между кластерами (мера того, насколько разделен каждый кластер). от остальных).
ssc <- data.frame(
kmeans = c(2,3,4,5,6,7,8),
within_ss = c(mean(km2$withinss), mean(km3$withinss), mean(km4$withinss), mean(km5$withinss), mean(km6$withinss), mean(km7$withinss), mean(km8$withinss)),
between_ss = c(km2$betweenss, km3$betweenss, km4$betweenss, km5$betweenss, km6$betweenss, km7$betweenss, km8$betweenss)
)
ssc %<>% gather(., key = "measurement", value = value, -kmeans)
#ssc$value <- log10(ssc$value)
ssc %>% ggplot(., aes(x=kmeans, y=log10(value), fill = measurement)) + geom_bar(stat = "identity", position = "dodge") + ggtitle("Cluster Model Comparison") + xlab("Number of Clusters") + ylab("Log10 Total Sum of Squares") + scale_x_discrete(name = "Number of Clusters", limits = c("0", "2", "3", "4", "5", "6", "7", "8"))
Из этого измерения кажется, что 7 кластеров были бы подходящим выбором.
NbClust
Пакет NbClust предоставляет 30 индексов для определения соответствующего количества кластеров и предлагает пользователям лучшую схему кластеризации из различных результатов, полученных путем варьирования всех комбинаций количества кластеров, мер расстояния и методов кластеризации.
res.nbclust <- NbClust(mammals_scaled, distance = "euclidean",
min.nc = 2, max.nc = 9,
method = "complete", index ="all")
factoextra::fviz_nbclust(res.nbclust) + theme_minimal() + ggtitle("NbClust's optimal number of clusters")
Это говорит о том, что оптимальное количество кластеров - 3.
Clustree
Статистический метод, описанный выше, дает единый балл, который одновременно учитывает только один набор кластеров. Пакет R clustree использует альтернативный подход, рассматривая, как образцы изменяют группировку по мере увеличения числа кластеров. Это полезно для того, чтобы показать, какие кластеры различны, а какие нестабильны. В нем не указывается явно, какой из оптимальных кластеров выбрать, но он полезен для изучения возможных вариантов.
Давайте рассмотрим от 1 до 11 кластеров.
tmp <- NULL
for (k in 1:11){
tmp[k] <- kmeans(mammals_scaled, k, nstart = 30)
}
df <- data.frame(tmp)
# add a prefix to the column names
colnames(df) <- seq(1:11)
colnames(df) <- paste0("k",colnames(df))
# get individual PCA
df.pca <- prcomp(df, center = TRUE, scale. = FALSE)
ind.coord <- df.pca$x
ind.coord <- ind.coord[,1:2]
df <- bind_cols(as.data.frame(df), as.data.frame(ind.coord))
clustree(df, prefix = "k")
На этом рисунке размер каждого узла соответствует количеству выборок в каждом кластере, а стрелки окрашены в соответствии с количеством выборок, которые получает каждый кластер. Отдельный набор прозрачных стрелок, называемый пропорцией входящих узлов, также окрашен и показывает, как образцы из одной группы попадают в другую группу - индикатор нестабильности кластера.
На этом графике мы видим, что при переходе от k = 2 к k = 3 ряд видов из кластера наблюдателей слева переназначаются в третий кластер справа. . При переходе от k = 8 к k = 9 мы видим один узел с несколькими входящими ребрами, что свидетельствует о чрезмерной кластеризации данных.
Также может быть полезно наложить это измерение на другие измерения данных, особенно те, которые получены с помощью методов уменьшения размерности. Мы можем сделать это с помощью функции clustree_overlay ():
df_subset <- df %>% select(1:8,12:13) clustree_overlay(df_subset, prefix = "k", x_value = "PC1", y_value = "PC2")
Я предпочитаю видеть его сбоку, показывая одно из измерений x или y в сравнении с размером разрешения.
overlay_list <- clustree_overlay(df_subset, prefix = "k", x_value = "PC1",
y_value = "PC2", plot_sides = TRUE)
overlay_list$x_side
overlay_list$y_side
Это показывает, что мы можем указать правильное разрешение кластеризации, исследуя края, и мы можем получить избыточную информацию для оценки качества кластеризации.
Выбор подходящего алгоритма
А как насчет выбора подходящего алгоритма кластеризации? Пакет cValid можно использовать для одновременного сравнения нескольких алгоритмов кластеризации, чтобы определить лучший подход к кластеризации и оптимальное количество кластеров. Мы сравним k-средних, иерархическую кластеризацию и кластеризацию PAM.
intern <- clValid(mammals_scaled, nClust = 2:24,
clMethods = c("hierarchical","kmeans","pam"), validation = "internal")
# Summary
summary(intern) %>% kable() %>% kable_styling()
Clustering Methods:
hierarchical kmeans pam
Cluster sizes:
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
Validation Measures:
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
hierarchical Connectivity 4.1829 10.5746 13.2579 20.1579 22.8508 25.8258 32.6270 35.3032 38.2905 39.2405 41.2405 45.7742 47.2742 50.6075 52.6075 55.8575 58.7242 60.7242 63.2242 65.2242 67.2242 69.2242 71.2242
Dunn 0.3595 0.3086 0.3282 0.2978 0.3430 0.3430 0.4390 0.4390 0.5804 0.5938 0.5938 0.8497 0.8497 0.5848 0.5848 0.4926 0.9138 0.9138 0.8892 0.9049 0.9335 1.0558 2.1253
Silhouette 0.5098 0.5091 0.4592 0.4077 0.4077 0.3664 0.3484 0.4060 0.3801 0.3749 0.3322 0.3646 0.3418 0.2650 0.2317 0.2166 0.2469 0.2213 0.1659 0.1207 0.1050 0.0832 0.0691
kmeans Connectivity 7.2385 10.5746 15.8159 20.1579 22.8508 25.8258 33.5198 35.3032 38.2905 39.2405 41.2405 45.7742 47.2742 51.8909 53.8909 57.1409 58.7242 60.7242 63.2242 65.2242 67.2242 69.2242 71.2242
Dunn 0.2070 0.3086 0.2884 0.2978 0.3430 0.3430 0.3861 0.4390 0.5804 0.5938 0.5938 0.8497 0.8497 0.5866 0.5866 0.5725 0.9138 0.9138 0.8892 0.9049 0.9335 1.0558 2.1253
Silhouette 0.5122 0.5091 0.4260 0.4077 0.4077 0.3664 0.3676 0.4060 0.3801 0.3749 0.3322 0.3646 0.3418 0.2811 0.2478 0.2402 0.2469 0.2213 0.1659 0.1207 0.1050 0.0832 0.0691
pam Connectivity 7.2385 14.1385 17.4746 24.0024 26.6857 32.0413 33.8913 36.0579 38.6607 40.6607 42.7869 45.7742 47.2742 51.7242 53.7242 56.9742 58.7242 60.7242 62.7242 64.7242 66.7242 69.2242 71.2242
Dunn 0.2070 0.1462 0.2180 0.2180 0.2978 0.2980 0.4390 0.4390 0.4390 0.4390 0.4390 0.8497 0.8497 0.5314 0.5314 0.4782 0.9138 0.9138 0.8333 0.8189 0.7937 1.0558 2.1253
Silhouette 0.5122 0.3716 0.4250 0.3581 0.3587 0.3318 0.3606 0.3592 0.3664 0.3237 0.3665 0.3646 0.3418 0.2830 0.2497 0.2389 0.2469 0.2213 0.1758 0.1598 0.1380 0.0832 0.0691
Optimal Scores:
Score Method Clusters
Connectivity 4.1829 hierarchical 2
Dunn 2.1253 hierarchical 24
Silhouette 0.5122 kmeans 2
Связность и силуэт являются измерениями связности, в то время как индекс Данна - это отношение наименьшего расстояния между наблюдениями, находящимися не в одном кластере, к наибольшему расстоянию внутри кластера.
Извлечение характеристик кластеров
Как упоминалось ранее, сложно оценить качество результатов кластеризации. У нас нет настоящих ярлыков, поэтому непонятно, как можно измерить, «насколько хорошо это работает» с точки зрения внутренней проверки. Однако кластеризация - отличная отправная точка EDA для более подробного изучения различий между кластерами. Думайте о кластеризации, как о производстве размеров рубашек. Мы могли выбрать только три размера: маленький, средний и большой. мы обязательно сократим расходы, но не всем подойдет. Подумайте о размерах брюк (или о брендах рубашек с множеством размеров (XS, XL, XXL, и т. Д.)), где у вас гораздо больше категорий (или групп). Для некоторых полей выбор оптимального кластера может зависеть от некоторых внешних знаний, таких как стоимость производства k кластеров, чтобы удовлетворить клиентов с наилучшим возможным соответствием. В других областях, таких как биология, где вы пытаетесь определить точное количество клеток, потребуется более глубокий подход. Например, многие из вышеперечисленных эвристик противоречили друг другу в отношении оптимального количества кластеров. Имейте в виду, что все они оценивали алгоритм k-средних при разных значениях k. Это потенциально может означать, что алгоритм k-means не работает, а k не подходит. Алгоритм k-средних не очень надежный алгоритм, чувствительный к выбросам, и этот набор данных очень мал. Лучше всего было бы изучить вышеперечисленные методы на выходе других алгоритмов (например, иерархическая кластеризация, предложенная clValid), собрать больше данных или потратить некоторое время на маркировку образцов для других методов машинного обучения, если возможный.
В конечном итоге мы хотели бы ответить на такие вопросы, как «что отличает этот кластер от других?» и «какие кластеры похожи друг на друга». Выберем пять кластеров и рассмотрим особенности этих кластеров.
# Compute dissimilarity matrix with euclidean distances d <- dist(mammals_scaled, method = "euclidean") # Hierarchical clustering using Ward's method res.hc <- hclust(d, method = "ward.D2" ) # Cut tree into 5 groups grp <- cutree(res.hc, k = 5) # Visualize plot(res.hc, cex = 0.6) # plot tree rect.hclust(res.hc, k = 5, border = 2:5) # add rectangle
# Execution of k-means with k=5
final <- kmeans(mammals_scaled, 5, nstart = 30)
fviz_cluster(final, data = mammals_scaled) + theme_minimal() + ggtitle("k = 5")
Давайте извлечем кластеры и добавим их обратно к нашим исходным данным, чтобы выполнить некоторую описательную статистику на уровне кластера:
as.data.frame(mammals_scaled) %>% mutate(Cluster = final$cluster) %>% group_by(Cluster) %>% summarise_all("mean") %>% kable() %>% kable_styling()
Мы видим, что кластер 2, состоящий исключительно из кролика, имеет высокую зольность. Группа 3, состоящая из тюленей и дельфинов, с высоким содержанием жира, что имеет смысл, учитывая суровые требования такого холодного климата, в то время как группа 4 имеет высокое содержание лактозы.
mammals_df <- as.data.frame(mammals_scaled) %>% rownames_to_column()
cluster_pos <- as.data.frame(final$cluster) %>% rownames_to_column()
colnames(cluster_pos) <- c("rowname", "cluster")
mammals_final <- inner_join(cluster_pos, mammals_df)
ggRadar(mammals_final[-1], aes(group = cluster), rescale = FALSE, legend.position = "none", size = 1, interactive = FALSE, use.label = TRUE) + facet_wrap(~cluster) + scale_y_discrete(breaks = NULL) + # don't show ticks
theme(axis.text.x = element_text(size = 10)) + scale_fill_manual(values = rep("#1c6193", nrow(mammals_final))) +
scale_color_manual(values = rep("#1c6193", nrow(mammals_final))) +
ggtitle("Mammals Milk Attributes")
mammals_df <- as.data.frame(mammals_scaled)
mammals_df$cluster <- final$cluster
mammals_df$cluster <- as.character(mammals_df$cluster)
ggpairs(mammals_df, 1:5, mapping = ggplot2::aes(color = cluster, alpha = 0.5),
diag = list(continuous = wrap("densityDiag")),
lower=list(continuous = wrap("points", alpha=0.9)))
# plot specific graphs from previous matrix with scatterplot
g <- ggplot(mammals_df, aes(x = water, y = lactose, color = cluster)) +
geom_point() +
theme(legend.position = "bottom")
ggExtra::ggMarginal(g, type = "histogram", bins = 20, color = "grey", fill = "blue")
b <- ggplot(mammals_df, aes(x = protein, y = fat, color = cluster)) +
geom_point() +
theme(legend.position = "bottom")
ggExtra::ggMarginal(b, type = "histogram", bins = 20, color = "grey", fill = "blue")
ggplot(mammals_df, aes(x = cluster, y = protein)) +
geom_boxplot(aes(fill = cluster))
ggplot(mammals_df, aes(x = cluster, y = fat)) +
geom_boxplot(aes(fill = cluster))
ggplot(mammals_df, aes(x = cluster, y = lactose)) +
geom_boxplot(aes(fill = cluster))
ggplot(mammals_df, aes(x = cluster, y = ash)) +
geom_boxplot(aes(fill = cluster))
ggplot(mammals_df, aes(x = cluster, y = water)) +
geom_boxplot(aes(fill = cluster))
# Parallel coordiante plots allow us to put each feature on seperate column and lines connecting each column ggparcoord(data = mammals_df, columns = 1:5, groupColumn = 6, alphaLines = 0.4, title = "Parallel Coordinate Plot for the Mammals Milk Data", scale = "globalminmax", showPoints = TRUE) + theme(legend.position = "bottom")
Если вы найдете эту статью полезной, не стесняйтесь поделиться ею с другими или порекомендовать эту статью! 😃
Как всегда, если у вас есть какие-либо вопросы или комментарии, не стесняйтесь оставлять свои отзывы ниже или вы всегда можете связаться со мной в LinkedIn. А пока до встречи в следующем посте! 😄
