Кластеризація з пакетом ClusterR, частина 1

Ця стаття присвячена кластеризації, а точніше, мою нещодавно доданого в CRAN пакету ClusterR. Деталі та приклади нижче в більшості своїй засновані на пакеті Vignette.

Кластерний аналіз або кластеризація — завдання групування набору об'єктів таким чином, щоб об'єкти всередині однієї групи (званої кластером) були більше схожі (в тому чи іншому сенсі) один на одного, ніж на об'єкти в інших групах (кластерах). Це одна з головних завдань дослідницького аналізу даних і стандартна техніка статистичного аналізу, що застосовується в різних сферах, у т. ч. машинному навчанні, розпізнаванні образів, аналізі зображень, пошуку інформації, біоінформатики, стисненні даних, комп'ютерної графіки.

Найбільш відомі приклади алгоритмів кластеризації кластеризація на основі зв'язності (ієрархічна кластеризація), кластеризація на основі центрів (метод k-середніх, метод k-медоидов), кластеризація на основі розподілів (GMM — Gaussian mixture models — Гаусова суміш розподілів) і кластеризація на основі щільності (DBSCAN — Density-based spatial clustering of applications with noise — просторова кластеризація додатків з шумом на основі щільності, OPTICS — Ordering points to identify the clustering structure — впорядкування точок для визначення структури кластеризації, і ін).

Пакет ClusterR складається з алгоритмів кластеризації на основі центрів (метод k-середніх, k-середніх у міні-групах, k-медоидов) і розподілів (GMM). Також пакет пропонує функції для:

  • валідації вихідних даних з допомогою міток істинності;
  • виводу результатів на контурному або двовимірному графіку;
  • передбачення нових результатів спостережень;
  • оцінки оптимального числа кластерів для кожного окремого алгоритму.

Гаусового суміш розподілів (GMM — Gaussian mixture models)

Гаусового суміш розподілів — статистична модель для представлення нормально розподілених субпопуляцій межах загальної популяції. Гаусова суміш розподілів параметризируется двома типами значень — суміш ваг компонентів середніх компонентів або коваріацій (багатовимірного випадку). Якщо кількість компонентів відомо, техніка, найчастіше використовується для оцінки параметрів суміші розподілів — ЕМ-алгоритм.

Функція GMM в пакеті ClusterR — реалізація на R класу для моделювання даних як гауссових суміші розподілів (GMM) з бібліотеки Armadillo з припущенням про діагональних коваріаційних матрицях. Можна настроювати параметри функції, в тому числі gaussian_comps, dist_mode (eucl_dist, maha_dist), seed_mode (static_subset, random_subset, static_spread, random_spread), km_iter em_iter (більше інформації про параметри у документації до пакета). Проілюструю функцію GMM на синтетичних даних dietary_survey_IBS.

library(ClusterR)
data(dietary_survey_IBS)
dim(dietary_survey_IBS)

## [1] 400 43

X = dietary_survey_IBS[, -ncol(dietary_survey_IBS)] # дані (без залежної змінної)
y = dietary_survey_IBS[, ncol(dietary_survey_IBS)] # залежна змінна
dat = center_scale(X, mean_center = T, sd_scale = T) # центрування і масштабування даних
gmm = GMM(dat, 2, dist_mode = "maha_dist", seed_mode = "random_subset", km_iter = 10,
em_iter = 10, verbose = F) 

# передбачення центроїдів, ковариационной матриці і ваг
pr = predict_GMM(dat, gmm$centroids, gmm$covariance_matrices, gmm$weights) 

Початково функція GMM повертає центроїди, матрицю коваріації (де кожен рядок являє діагональну ковариационную матрицю), ваги логарифмічні функції правдоподібності для кожної гауссових компоненти. Потім функція predict_GMM бере вихідні дані моделі GMM і повертає ймовірні кластери.

На додаток до вже згаданих функцій можна використовувати Optimal_Clusters_GMM для оцінки кількості кластерів даних за допомогою інформаційного критерію Акаіке (AIC, Akaike information), або байєсівського інформаційного критерію (BIC, Bayesian information).

opt_gmm = Optimal_Clusters_GMM(dat, max_clusters = 10, criterion = "BIC", 
dist_mode = "maha_dist", seed_mode = "random_subset",
km_iter = 10, em_iter = 10, var_floor = 1e-10, 
plot_data = T)


При виборі моделі з зумовленого безлічі краще вибирати ту, у якої найнижчий BIC, тут це справедливо для кількості кластерів, рівного 2.

Припускаючи, що доступні мітки істинності, для валідації вихідних кластерів можна використовувати методи external_validation (rand_index, adjusted_rand_index, jaccard_index, fowlkes_Mallows_index, mirkin_metric, purity, entropy, nmi (нормалізована взаємна інформація) і var_info (варіація інформації).

res = external_validation(dietary_survey_IBS$class, pr$cluster_labels, 
method = "adjusted_rand_index", summary_stats = T)
res

## 
## ---------------------------------------- 
## purity : 1 
## entropy : 0 
## normalized mutual information : 1 
## variation of information : 0 
## ---------------------------------------- 
## specificity : 1 
## sensitivity : 1 
## precision : 1 
## recall : 1 
## F-measure : 1 
## ---------------------------------------- 
## OR accuracy rand-index : 1 
## adjusted-rand-index : 1 
## jaccard-index : 1 
## fowlkes-mallows-index : 1 
## mirkin-metric : 0 
## ----------------------------------------

І якщо параметр summary_stats встановлено в TRUE, то повертаються також метрики специфічності, чутливості, точності, повноти, F-заходи (specificity, sensitivity, precision, recall, F-measure відповідно).

Метод k-середніх

Кластеризація методом k-середніх — метод векторного квантування, що застосовується в обробці сигналів, часто використовується для кластерного аналізу даних. Мета кластеризації методом k-середніх — розділити n значення k кластерів, в яких кожне значення належить кластеру з найближчим середнім, виступає прототипом кластера. Це призводить до поділу області даних на клітинки Вороного. Найбільш часто застосовуваний алгоритм використовує итеративную уточнюючу техніку. Через повсюдне вживання, його називають алгоритмом k-середніх; зокрема, серед фахівців у сфері комп'ютерних наук він також відомий як алгоритм Ллойда.

Пакет ClusterR надає дві різні функції k-середніх, KMeans_arma, реалізацію на R методу k-середніх з бібліотеки armadillo, і KMeans_rcpp, яка використовує пакет RcppArmadillo. Обидві функції приводять до одних і тих же результатів, однак, повертають різні ознаки (код нижче ілюструє це).

KMeans_arma
KMeans_arma швидше, ніж функція KMeans_rcpp, однак, початково вона виводить центроїди тільки деяких кластерів. Більш того, кількість колонок в даних повинна перевищувати кількість кластерів, інакше функція поверне помилку. Кластеризація буде працювати швидше на багатоядерних машинах з включеним OpenMP (наприклад, -fopenmp в GCC). Алгоритм ініціалізується один раз, і зазвичай 10 ітерацій достатньо для збіжності. Вихідні центроїди розподіляються з допомогою одного з алгоритмів keep_existing, static_subset, random_subset, static_spread або random_spread. Якщо seed_mode одно keep_existing, користувач повинен передати на вхід матрицю центроїдів.

Я зменшу кількість вимірів у даних dietary_survey_IBS з допомогою аналізу головних компонент (PCA — principal component analysis), а саме — функції princomp з пакету stats, щоб можна було вивести двовимірний графік кластерів, побудованих в результаті.

pca_dat = stats::princomp(dat)$scores[, 1:2]
km = KMeans_arma(pca_dat, clusters = 2, n_iter = 10, seed_mode = "random_subset", 
verbose = T, CENTROIDS = NULL)
pr = predict_KMeans(pca_dat, km)
table(dietary_survey_IBS$class, pr)
class(km) = 'matrix'
plot_2d(data = pca_dat, clusters = as.vector(pr), centroids_medoids = as.matrix(km)) 


KMeans_rcpp
Як стверджувалося вище, функція KMeans_rcpp пропонує деякі додаткові можливості в порівнянні з функцією KMeans_arma:

  • вона дозволяє більш однієї ініціалізації (можна параллелизировать з OpenMP)
  • крім инициализаций optimal_init, quantile_init, random або kmeans++, можна задати центроїди в параметрі CENTROIDS
  • час роботи алгоритму і його збіжність можна налаштовувати параметрами num_init, max_iters tol
  • якщо num_init > 1, KMeans_rcpp повертає атрибути найкращої ініціалізації, використовуючи критерій within-cluster-sum-of-squared-помилка (сума квадратів помилок-всередині кластера)
  • алгоритм повертає наступні атрибути: clusters, fuzzy_clusters (якщо fuzzy = TRUE), centroids, total_SSE, best_initialization, WCSS_per_cluster, obs_per_cluster, between.SS_DIV_total.SS
Більше подробиць про KMeans_rcpp є в документації пакета. Проілюструю різні параметри KMeans_rcpp на прикладі векторного квантування і пакета OpenImageR.

library(OpenImageR)
path = 'elephant.jpg'
im = readImage(path)

# спочатку змінимо розмір картинки, щоб зменшити кількість вимірювань
im = resizeImage(im, 75, 75, method = 'bilinear') 
imageShow(im) # вивід вихідного зображення
im2 = apply(im, 3, as.vector) # векторизированный RGB 


# кластеризація з KMeans_rcpp
km_rc = KMeans_rcpp(im2, clusters = 5, num_init = 5, max_iters = 100, 
initializer = 'optimal_init', threads = 1, verbose = F)
km_rc$between.SS_DIV_total.SS

## [1] 0.9873009

Атрибут between.SS_DIV_total.SS дорівнює (total_SSE — sum(WCSS_per_cluster)) / total_SSE. Якщо у кластеризації немає закономірності, проміжна сума квадратів буде дуже невеликою частиною загальної суми квадратів, а якщо атрибут between.SS_DIV_total.SS близький до 1.0, значення кластеризуются досить добре.

getcent = km_rc$centroids
getclust = km_rc$clusters
new_im = getcent[getclust, ] # кожне значення асоціюється з найближчим центроїдом
dim(new_im) = c(nrow(im), ncol(im), 3) # зворотне перетворення до тривимірної картинки
imageShow(new_im)


Крім того, можна скористатися функцією Optimal_Clusters_KMeans (яка неявно використовують KMeans_rcpp), щоб визначити оптимальну кількість кластерів. Доступні такі критерії: variance_explained, WCSSE (within-cluster-sum-of-squared-error — сума квадратів помилок-всередині кластера), dissimilarity, silhouette, distortion_fK, AIC, BIC Adjusted_Rsquared. У документації по пакету є більше інформації по кожному критерію.

У наступному прикладі коду використовується критерій distortion_fK, який повністю описаний в статті “Selection of K in K-means clustering, Pham., Dimov., Nguyen., (2004)" («Вибір До кластеризації методом К-середніх»).

opt = Optimal_Clusters_KMeans(im2, max_clusters = 10, plot_clusters = T, 
verbose = F, criterion = 'distortion_fK', fK_threshold = 0.85)


Значення менше фіксованого кордону (тут fK_threshold = 0.85) можна рекомендувати для кластеризації. Однак, існує одного оптимального розбиття на кластери, тому f(K) слід використовувати тільки для припущення про їх кількість, але остаточне рішення, приймати це значення чи ні — за користувачем.

Метод k-середніх у міні-групах

Метод k-середніх у міні-групах — варіація класичного алгоритму k-середніх. Він особливо корисний для великих наборів даних, оскільки замість використання всього набору (як в k-середніх) беруться міні-групи з випадкових значень даних, щоб оптимізувати цільову функцію.

Параметри алгоритму MiniBatchKmeans майже такі ж, як і функції KMeans_rcpp з пакету ClusterR. Найбільш важлива відмінність — параметр batch_size (розмір міні-груп) і init_fraction (відсоток даних для визначення вихідних центроїдів, який застосовується, якщо initializer дорівнює 'kmeans++' або 'quantile_init').

Скористаюся прикладом векторного квантування, щоб показати різницю в часі обчислення і якості висновку функцій KMeans_rcpp і MiniBatchKmeans.

path_d = 'dog.jpg'
im_d = readImage(path_d)

# спочатку змінимо розмір картинки, щоб зменшити кількість вимірювань
im_d = resizeImage(im_d, 350, 350, method = 'bilinear') 
imageShow(im_d) # вивід вихідного зображення


im3 = apply(im_d, 3, as.vector) # векторизированный RGB
dim(im3) # вихідні вимірювання даних

# 122500 3

Спочатку виконаємо кластеризацію k-середніх.

start = Sys.time()
km_init = KMeans_rcpp(im3, clusters = 5, num_init = 5, max_iters = 100, 
initializer = 'kmeans++', threads = 1, verbose = F)
end = Sys.time()
t = end - start
cat('time to complete :', t, attributes(t)$units, '\n')

# час виконання : 2.44029 secs

getcent_init = km_init$centroids
getclust_init = km_init$clusters
new_im_init = getcent_init[getclust_init, ] # кожне значення асоціюється з найближчим центроїдом
dim(new_im_init) = c(nrow(im_d), ncol(im_d), 3) # зворотне перетворення до тривимірної картинки
imageShow(new_im_init)


Тепер кластеризацію k-середніх у міні-групах.

start = Sys.time()
km_mb = MiniBatchKmeans(im3, clusters = 5, batch_size = 20, num_init = 5, max_iters = 100, 
init_fraction = 0.2, initializer = 'kmeans++', early_stop_iter = 10,
verbose = F)
pr_mb = predict_MBatchKMeans(im3, km_mb$centroids)
end = Sys.time()
t = end - start
cat('time to complete :', t, attributes(t)$units, '\n')

# час виконання : 0.8346727 secs 

getcent_mb = km_mb$centroids
new_im_mb = getcent_mb[pr_mb, ] # кожне значення асоціюється з найближчим центроїдом
dim(new_im_mb) = c(nrow(im_d), ncol(im_d), 3) # зворотне перетворення до тривимірної картинки
imageShow(new_im_mb)


Незважаючи на невелику різницю в якості висновку, k-середніх у міні-групах повертає результат в середньому в два рази швидше, ніж класичний метод k-середніх для цього набору даних.

Щоб реалізувати метод k-середніх у міні-групах, я використовував наступні ресурси:

Джерело: Хабрахабр

0 коментарів

Тільки зареєстровані та авторизовані користувачі можуть залишати коментарі.