Уровень алгоритма

Участник:Бротиковская Данута/Алгоритм k-means: различия между версиями

Материал из Алговики
Перейти к навигации Перейти к поиску
 
(не показано 379 промежуточных версий 5 участников)
Строка 1: Строка 1:
 +
{{Assignment|Zhum|Konshin}}
 +
 
{{algorithm
 
{{algorithm
 
| name              = Алгоритм k средних (k means)
 
| name              = Алгоритм k средних (k means)
| serial_complexity = <math>O(n^3)</math>
+
| serial_complexity = <math>O(ikdn)</math>
| pf_height        = <math>O(n)</math>
+
| input_data        = <math> dn </math>
| pf_width          = <math>O(n^2)</math>
+
| output_data      = <math> n </math>
| input_data        = <math>\frac{n (n + 1)}{2}</math>
 
| output_data      = <math>\frac{n (n + 1)}{2}</math>
 
 
}}
 
}}
  
Авторы страницы <b>Данута Бротиковская</b> и <b>Денис Зобнин</b>
+
Авторы страницы <b>[[Участник:Бротиковская Данута|Данута Бротиковская]]</b> и <b>[[Участник:DennZo1993|Денис Зобнин]]</b>.
 +
 
 +
Оба автора внесли одинаковый вклад в написание статьи (поиск материалов, изучение статей, заполнение контента статьи), четких разграничений по пунктам статьи сделать невозможно.
  
 
== Свойства и структура алгоритма ==
 
== Свойства и структура алгоритма ==
  
 
=== Общее описание алгоритма ===
 
=== Общее описание алгоритма ===
Алгоритм <b>k средних</b> (<b>k means</b>) -- наиболее популярный метод кластеризации. Был изобретен в 1950-х годах математиком <i>Гуго Штейнгаузом</i> и почти одновременно Стюартом Ллойдом. Особую популярность приобрел после публикации работы МакКуина в 1967.
+
Алгоритм <b><i>k средних</i></b> (<b><i>k means</i></b>) &ndash; неиерархический<ref>"https://ru.wikipedia.org/wiki/Иерархическая_кластеризация"</ref>, итерационный метод кластеризации<ref>"https://ru.wikipedia.org/wiki/Кластерный_анализ"</ref>, получивший большую популярность благодаря своей простоте, наглядности реализации и достаточно высокому качеству работы. Был изобретен в 1950-х годах математиком <i>Гуго Штейнгаузом</i><ref>Steinhaus, Hugo. "Sur la division des corp materiels en parties." Bull. Acad. Polon. Sci 1.804 (1956): 801.</ref> и почти одновременно <i>Стюартом Ллойдом</i><ref>Lloyd, S. P. "Least square quantization in PCM. Bell Telephone Laboratories Paper. Published in journal much later: Lloyd, SP: Least squares quantization in PCM." IEEE Trans. Inform. Theor.(1957/1982).</ref>. Особую популярность приобрел после публикации работы <i>МакКуина</i><ref>MacQueen, James. "Some methods for classification and analysis of multivariate observations." Proceedings of the fifth Berkeley symposium on mathematical statistics and probability. Vol. 1. No. 14. 1967.</ref> в 1967.
Цель алгоритма заключается в разделении N наблюдений на K кластеров таким образом, чтобы каждое наблюдение придележало ровно одному кластеру, расположенному на наименьшем расстоянии от наблюдения.
+
 
 +
Алгоритм представляет собой версию EM-алгоритма<ref>"https://ru.wikipedia.org/wiki/EM-алгоритм"</ref>, применяемого также для разделения смеси гауссиан. Основная идея <i>k means</i> заключается в том, что данные произвольно разбиваются на кластеры, после чего итеративно перевычисляется центр масс для каждого кластера, полученного на предыдущем шаге, затем векторы разбиваются на кластеры вновь в соответствии с тем, какой из новых центров оказался ближе по выбранной метрике.
 +
 
 +
Цель алгоритма заключается в разделении <math>n</math> наблюдений на <math>k</math> кластеров таким образом, чтобы каждое наблюдение принадлежало ровно одному кластеру, расположенному на наименьшем расстоянии от наблюдения.
  
 
=== Математическое описание алгоритма ===
 
=== Математическое описание алгоритма ===
Дан набор из <math>n</math> d-мерных векторов <math>X</math>&nbsp;=&nbsp;<math>\{x_1, x_2, ..., x_n\}</math>. Алгоритм k средних разбивает набор <math>X</math> на <math>k, k<=n</math> наборов <math>S=\{S_1, S_2, ..., S_k\}</math> таким образом, чтобы минимизировать сумму квадратов расстояний от каждой точки кластера до его центра. Другими словами:
+
<b><u>Дано:</u></b>
 +
* набор из <math>n</math> наблюдений <math>X=\{\mathbf{x}_1, \mathbf{x}_2, ..., \mathbf{x}_n\}, \mathbf{x}_i \in \mathbb{R}^d, i=\overline{1, n}</math>;
 +
* <math>k</math> - требуемое число кластеров, <math>k \in \mathbb{N}, k \leq n</math>.
  
<math>\arg\min_{s} \sum\limits_{i=1}^k \sum\limits x \in S_i \lVert x-\mu \rVert^2</math>
+
<b><u>Требуется:</u></b>
 +
 
 +
Разделить множество наблюдений <math>X</math> на <math>k</math> кластеров <math>S_1, S_2, ..., S_k</math>:
 +
*  <math>S_i \cap S_j= \varnothing, \quad i \ne j</math>
 +
 
 +
*  <math>\bigcup_{i=1}^{k} S_i = X</math>
 +
 
 +
 
 +
<b><u>Действие алгоритма:</u></b>
 +
<p>Алгоритм <i>k means</i> разбивает набор <math>X</math> на <math>k</math> наборов <math>S_1, S_2, ..., S_k,</math> таким образом, чтобы минимизировать сумму квадратов расстояний от каждой точки кластера до его центра (центр масс кластера). Введем обозначение, <math>S=\{S_1, S_2, ..., S_k\}</math>. Тогда действие алгоритма <i>k means</i> равносильно поиску:</p>
 +
<table style="width:100%">
 +
  <tr>
 +
    <td align="center" width="90%"><math>\arg\min_{S} \sum\limits_{i=1}^k \sum\limits_{\mathbf{x} \in S_i} \rho(\mathbf{x}, \mathbf{\mu}_i )^2,</math></td>
 +
    <td align="right"><math>(1)</math></td>
 +
  </tr>
 +
</table>
 +
где <math>\mathbf{\mu}_i</math> &ndash; центры кластеров, <math>i=\overline{1,k}, \quad \rho(\mathbf{x}, \mathbf{\mu}_i)</math> &ndash; функция расстояния между <math>\mathbf{x}</math> и <math>\mu_i</math>
  
  
 
<b><u>Шаги алгоритма: </u></b>
 
<b><u>Шаги алгоритма: </u></b>
 +
<ol>
 +
  <li>
 +
    <b>Начальный шаг: инициализация кластеров</b>
 +
    <p>Выбирается произвольное множество точек <math>\mu_i, i=\overline{1,k}</math>, рассматриваемых как начальные центры кластеров: <math>\mu_i^{(0)} = \mu_i, \quad i=\overline{1, k}</math></p>
 +
  </li>
 +
  <li>
 +
    <b>Распределение векторов по кластерам</b>
 +
    <p><b>Шаг</b> <math>t: \forall \mathbf{x}_i \in X, i=\overline{1,n}: \mathbf{x}_i \in S_j \iff  j=\arg\min_{k}\rho(\mathbf{x}_i,\mathbf{\mu}_k^{(t-1)})^2</math></p>
 +
  </li>
 +
  <li>
 +
    <b>Пересчет центров кластеров</b>
 +
    <p><b>Шаг </b> <math>t: \forall i=\overline{1,k}:  \mu_i^{(t)} = \cfrac{1}{|S_i|}\sum_{\mathbf{x}\in S_i}\mathbf{x}</math></p>
 +
  </li>
 +
  <li>
 +
    <b>Проверка условия останова:</b>
 +
    <p></p>
 +
    '''if''' <math>\exist i\in \overline{1,k}: \mu_i^{(t)} \ne \mu_i^{(t-1)}</math> '''then'''
 +
        t = t + 1;
 +
        goto 2;
 +
    '''else'''
 +
        '''stop'''
 +
  </li>
 +
</ol>
 +
 +
=== Вычислительное ядро алгоритма ===
  
<b>Начальный шаг: </b>Инициализация кластеров. Выбирается произвольное множество k точек, рассматриваемых как начальные центры кластеров.
+
Вычислительным ядром являются шаги 2 и 3 приведенного выше алгоритма: <b><i>распределение векторов по кластерам</i></b> и <b><i>пересчет центров кластеров</i></b>.  
  
<b>
+
<p>
 +
<b><i>Распределение векторов</i></b> по кластерам предполагает вычисление расстояний между каждым вектором <math>\mathbf{x}_i \in X, i= \overline{1,n}</math> и центрами кластера <math>\mathbf{\mu}_j, j= \overline{1,k}</math>. Таким образом, данный шаг предполагает <math>kn</math> вычислений расстояний между d-мерными векторами.
 +
</p>
 +
<p>
 +
<b><i>Пересчет центров кластеров</i></b> предполагает <math>k</math> вычислений центров масс <math>\mathbf{\mu}_i</math> множеств <math>S_i, i=\overline{1,k}</math> представленных выражением в шаге 3 представленного выше алгоритма.
 +
</p>
  
=== Вычислительное ядро алгоритма ===
 
 
=== Макроструктура алгоритма ===
 
=== Макроструктура алгоритма ===
 +
 +
<b>Инициализация центров масс <math>\mu_1, ..., \mu_k</math></b>.
 +
 +
Наиболее распространенными являются следующие стратегии:
 +
<ul>
 +
  <li>
 +
    <b>Метод Forgy</b><br>В качестве начальных значений <math>\mu_1, ..., \mu_k</math> берутся случайно выбранные векторы.
 +
  </li>
 +
  <li>
 +
    <b>Метод случайно разделения (Random Partitioning)</b><br>Для каждого вектора <math>\mathbf{x}_i \in X, i=\overline{1, n}</math> выбирается случайным образом кластер <math>S_1, ..., S_k</math>, после чего для каждого полученного кластера вычисляются значения <math>\mu_1, ..., \mu_k</math>.
 +
  </li>
 +
</ul>
 +
 +
 +
 +
<b>Распределение векторов по кластерам</b>
 +
 +
Для этого шага алгоритма между векторами <math>\mathbf{x}_i \in X, i=\overline{1, n}</math> и центрами кластеров <math>\mu_1, ..., \mu_k</math> вычисляются <b>расстояния</b> по формуле (как правило, используется Евлидово расстояние):
 +
<table style="width:100%">
 +
  <tr>
 +
    <td align="center" width="90%"><math> \mathbf{v}_1, \mathbf{v}_2 \in \mathbb{R}^d,  \quad \rho(\mathbf{v}_1, \mathbf{v}_2) = \lVert \mathbf{v}_1- \mathbf{v}_2 \rVert= \sqrt{\sum_{i=1}^{d}(\mathbf{v}_{1,i} - \mathbf{v}_{2,i})^2}</math></td>
 +
    <td align="right"><math>(2)</math></td>
 +
  </tr>
 +
</table>
 +
 +
 +
 +
<b>Пересчет центров кластеров</b>
 +
 +
Для этого шага алгоритма производится пересчет центров кластера по <b>формуле вычисления центра масс</b>:
 +
<table style="width:100%">
 +
  <tr>
 +
    <td align="center" width="90%"><math> \mu = \cfrac{1}{|S|}\sum_{\mathbf{x}\in S}\mathbf{x}</math></td>
 +
    <td align="right"><math>(3)</math></td>
 +
  </tr>
 +
</table>
 +
 
=== Схема реализации последовательного алгоритма ===
 
=== Схема реализации последовательного алгоритма ===
=== Последовательная сложность алгоритма ===
+
 
 +
'''algorithm''' k-means '''is'''
 +
    '''1.''' Инициализировать центры кластеров <math>\mathbf{\mu}_i^{(1)}, i=\overline{1,k}</math>
 +
    '''2.''' <math>t \leftarrow 1</math>
 +
    '''3.''' '''Распределение по кластерам'''<br>
 +
        <math>S_i^{(t)}=\{\mathbf{x}_p: \lVert\mathbf{x}_p-\mathbf{\mu}_i^{(t)}\rVert^2 \leq \lVert\mathbf{x}_p-\mathbf{\mu}_j^{(t)}\rVert^2 \quad \forall j=\overline{1,k}\},</math><br>
 +
        где каждый вектор <math>\mathbf{x}_p</math> соотносится единственному кластеру <math>S^{(t)}</math>
 +
    '''4.''' '''Обновление центров кластеров'''<br>
 +
        <math>\mathbf{\mu}_i^{(t+1)} = \frac{1}{|S^{(t)}_i|} \sum_{\mathbf{x}_j \in S^{(t)}_i} \mathbf{x}_j </math><br>
 +
    '''5.''' '''if''' <math>\exist i \in \overline{1,k}: \mathbf{\mu}_i^{(t+1)} \ne \mathbf{\mu}_i^{(t)}</math> '''then'''
 +
            t = t + 1;
 +
            goto '''3''';
 +
        '''else'''
 +
            '''stop'''
 +
 
 +
=== Последовательная сложность алгоритма===
 +
<div style="padding-bottom: 20px">
 +
<p>Обозначим <math>\Theta_{\rm centroid}^{d, m}</math> временную сложность вычисления центорида кластера, число элементов которого равна <math>m</math>,  в d-мерном пространстве.</p>
 +
<p>Аналогично  <math>\Theta_{\rm distance}^d</math> &ndash; временная сложность вычисления расстояния между двумя d-мерными векторами.</p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<b>Сложность шага инициализации <math>k</math> кластеров мощности <math>m</math> в d-мерном пространстве &ndash; <math>\Theta_{\rm init}^{k, d, m}</math></b>
 +
<ul>
 +
  <li><i>Стратерия Forgy</i>: вычисления не требуются, <math>\Theta_{\rm init}^{k, d, m} = 0</math></li>
 +
  <li><i>Стратегия случайного разбиения</i>: вычисление центров <math>k</math> кластеров,  <math>\Theta_{\rm init}^{k, d, m} = k \cdot \Theta_{\rm centroid}^{d, m}, m \le n</math></li>
 +
</ul>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p><b>Cложность шага распределения d мерных векторов по <math>k</math> кластерам &ndash; <math>\Theta_{\rm distribute}^{k, d}</math></b></p>
 +
<p>
 +
На этом шаге для каждого вектора <math>\mathbf{x}_i \in X, i=\overline{1, n}</math> вычисляется k расстояний до центров кластеров <math>\mathbf{\mu}_1, ...\mathbf{\mu}_k</math>
 +
</p>
 +
<p align="center">
 +
  <math>\Theta_{\rm distribute}^{k, d} = n \cdot k \cdot \Theta_{\rm distance}^d</math>
 +
</p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p>
 +
<b>Сложность шага пересчета центров <math>k</math> кластеров размера <math>m</math> в d-мерном пространстве &ndash; <math>\Theta_{\rm recenter}^{k, d, m}</math></b>
 +
</p>
 +
<p>На этом шаге вычисляется <math>k</math> центров кластеров <math>\mathbf{\mu}_1, ...\mathbf{\mu}_k</math></p>
 +
<p align="center"><math>\Theta_{\rm recenter}^{k, d, m} = k \cdot \Theta_{\rm centroid}^{d, m}</math> </p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p><b>Рассчитаем <math>\Theta_{\rm centroid}^{d, m}</math> для кластера, число элементов которого равно <math>m</math></b></p>
 +
<p align="center"><math>\Theta_{\rm centroid}^{d, m}</math> = <math>m \cdot d</math> сложений + <math>d</math> делений </p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p><b>Рассчитаем <math>\Theta_{\rm distance}^d</math> в соответствие с формулой <math>(2)</math></b></p>
 +
<p align="center"><math>\Theta_{\rm distance}^d</math> = <math>d</math> вычитаний + <math>d</math> умножений + <math>(d-1)</math> сложение</p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
  <p>Предположим, что алгоритм сошелся за <math>i</math> итераций, тогда временная сложность алгоритма <math>\Theta_{\rm k-means}^{d, n}</math></p>
 +
  <p align="center"><b><math>\Theta_{\rm k-means}^{d, n} \le \Theta_{\rm init}^{k, d, n} + i(\Theta_{\rm distribute}^{k, d} + \Theta_{\rm recenter}^{k, d, n})</math></b></p>
 +
</div>
 +
<div style="padding-bottom: 5px">
 +
<p>Операции сложения/вычитания:</p>
 +
<p align="center"><b><math>\Theta_{\rm k-means}^{d, n} \le knd+ i(kn(2d-1) + knd) = knd+ i(kn(3d-1)) \thicksim O(ikdn)</math></b></p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p>Операции умножения/деления:</p>
 +
<p align="center"><b><math>\Theta_{\rm k-means}^{d, n} \le kd + i(knd + kd) = kd + ikd(n+1) \thicksim O(ikdn) </math></b></p>
 +
</div>
 +
<div style="padding-bottom: 20px">
 +
<p>Получаем, что <b>временная сложность</b> алгоритма <b>k means</b> кластеризации <math>n</math> <b>d-мерных</b> векторов на <math>k</math> кластеров за <math>i</math> итераций:</p>
 +
<p align="center"><b><math> \Theta_{\rm k-means}^{d, n} \thicksim O(ikdn) </math></b></p>
 +
</div>
 +
 
 
=== Информационный граф ===
 
=== Информационный граф ===
 +
 +
Рассмотрим информационный граф алгоритма. Алгоритм <i>k means</i> начинается с этапа инициализации, после которого следуют итерации, на каждой из которых выполняется два последовательных шага (см. [[#Схема реализации последовательного алгоритма|"Схема реализации последовательного алгоритма"]]):
 +
* распределение векторов по кластерам
 +
* перерасчет центров кластеров
 +
 +
Поскольку основная часть вычислений приходится на шаги итераций, распишем информационные графы данных шагов.
 +
 +
 
 +
<p><b>Распределение векторов по кластерам</b></p>
 +
Информационный граф шага распределения векторов по кластерам представлен на ''рисунке 1''. Исходами данного графа является исходные векторы <math>\mathbf{x}_1, ..., \mathbf{x}_n</math>, а также центры кластеров <math>\mathbf{\mu}_1, ... \mathbf{\mu}_k</math>, вычисленные ранее (на шаге инициализации, если рассматривается первая итерация алгоритма, или на шаге пересчета центров кластеров предыдущей итерации в противном случае). Каждая пара векторов данных <math>\mathbf{x}_i, i=\overline{1,n}</math> и центров кластера <math>\mathbf{\mu}_j, j=\overline{1, k}</math> : (<math>\mathbf{x}_i</math>, <math>\mathbf{\mu}_j</math>) подаются на независимые узлы <i>"d"</i> вычисления расстояния между векторами (более подробная схема вычисления расстояния представлена далее, ''рисунок 2''). Далее узлы вычисления расстояния <i>"d"</i>, соответствующие одному и тому же исходному вектору <math>\mathbf{x}_i</math> передаются на один узел <i>"m"</i>, где далее происходит вычисление новой метки кластера для каждого вектора <math>\mathbf{x}_i</math> (берется кластер с минимальным результатом вычисления расстояния). На выходе графа выдаются метки кластеров , <math>L_1, ..., L_n</math>, такие что <math>\forall \mathbf{x}_i, i=\overline{1, n}, \mathbf{x}_i \in S_j \Leftrightarrow L_i = j</math>.
 +
 +
[[file:Clusterization.png|thumb|center|800px|Рис. 1. Схема распределения векторов по кластерам. ''d'' &ndash; вычисление расстояния между векторами; ''m'' &ndash; вычисление минимума.]]
 +
 +
 +
<p><b>Вычисление расстояния между векторами</b></p>
 +
Подробная схема вычисления расстояния между векторами <math>\mathbf{x}_i, \mathbf{\mu}_j</math> представлена на ''рисунке 2''. Как показано на графе, узел вычисления расстояния между векторами <i>"d"</i> состоит из шага взятия разности между векторами (узел "<math>-</math>") и взятия нормы получившегося вектора разности (узел "<math>||\cdot||^2</math>"). Более подробно, вычисление расстояния между векторами <math>\mathbf{x}_i = {x_{i1,}, ...,{x_{in}}}, \mathbf{\mu}_j = {\mu_{j1}, ...,\mu_{jn}}</math> может быть представлено как вычисление разности между каждой парой компонент <math>(x_{iz}, \mu_{jz}), z=\overline{1, d}</math> (узел "<math>-</math>"), далее возведение в квадрат для каждого узла "<math>-</math>" (узел "<math>()^2</math>") и суммирования выходов всех узлов "<math>()^2</math>" (узел "<math>+</math>"). 
 +
 +
[[file:dist_calc.png|thumb|center|800px|Рис. 2. Схема вычисления расстояния между вектором и центром кластера.]]
 +
 +
<p><b>Пересчет центров кластеров</b></p>
 +
Информационный граф шага пересчета центров кластеров представлен на ''рисунке 3''. Исходами данного графа является исходные векторы <math>\mathbf{x}_1, ..., \mathbf{x}_n</math>, а также им соответствующие метки кластера, <math>L_1, ..., L_n</math>, такие что <math>\forall x_i, i=\overline{1, n}, \mathbf{x}_i \in S_j \Leftrightarrow L_i = j</math>, вычисленные на этапе распределения векторов по кластерам. Все векторы <math>\mathbf{x}_1, ..., \mathbf{x}_n</math> подаются в узлы <math>+_1, ...+_k</math>, каждый узел <math>+_m, m = \overline{1,k}</math> соответствует операции сложения векторов кластера с номером <math>m</math>. Метки кластера <math>L_1, ..., L_n</math> также совместно передаются на узлы <math>S_m, m = \overline{1,k}</math>, на каждом из которых вычисляется количество векторов в соответствующем кластере (количество меток с соответствующим значением). Далее каждая пара выходов узлов <math>+_m</math> и <math>S_m</math> подается на узел "<math>/</math>", где производится деление суммы векторов кластера на количество элементов в нем. Значения, вычисленные на узлах "<math>/</math>", присваиваются новым центрам кластеров (выходные значения графа).
 +
 +
[[file:Recluster.png|thumb|center|800px|Рис. 3. Схема пересчета центров кластеров]]
 +
 
=== Ресурс параллелизма алгоритма ===
 
=== Ресурс параллелизма алгоритма ===
 +
<p>
 +
Работа алгоритма состоит из <math>i</math> итераций, в каждой из которых происходит <b>распределение d-мерных векторов по <math>k</math> кластерам</b>, а также <b>пересчет центров кластеров в d-мерном пространстве</b>. В шаге <b>распределения d-мерных векторов по <math>k</math> кластерам</b> расстояния между вектором и центрами кластеров вычисляются независимо (отсутствуют информационные зависимости). <b>Центры масс кластеров</b> также пересчитываются независимо друг от друга. Таким образом, имеет место [[глоссарий#Массовый параллелизм|''массовый параллелизм'']]. Вычислим параллельную сложность <math>\Psi_*</math> каждого из шагов, а также параллельную сложность всего алгоритма, <math>\Psi_{\rm k-means}</math>. Будем исходить из предположения, что может быть использовано любое необходимое число потоков.
 +
</p>
 +
 +
 +
<p><b>Распределение d-мерных векторов по <math>k</math> кластерам</b></p> 
 +
<p>
 +
Поскольку на данном шаге для каждой пары векторов <math>\mathbf{x}_i, i=\overline{1,n}</math> и <math>\mathbf{\mu}_j, j=\overline{1,k}</math> операции вычисления расстояния не зависят друг от друга, они могут выполняться параллельно. Тогда, разделив все вычисление расстояний на <math>n</math> потоков, получим, что в каждом потоке будет выполняться только одна операция вычисления расстояния между векторами размерности d. При этом каждому вычислительному потоку передаются координаты центров всех кластеров <math>\mathbf{\mu}_1, ..., \mathbf{\mu}_k</math>. Таким образом, параллельная сложность данного шага определяется <i> сложностью параллельной операции вычисления расстояния между d-мерными векторами</i>, <math>\Psi_{\rm distance}^d</math>  и <i>сложностью определения наиболее близкого кластера</i> (паралельное взятие минимума по расстояниям), <math>\Psi_{\rm min}^k</math>. Для оценки <math>\Psi_{\rm distance}^d</math> воспользуемся  [[Нахождение_частных_сумм_элементов_массива_сдваиванием | параллельной реализацией нахождения частичной суммы элементов массива путем сдваивания]]. Аналогично, <math>\Psi_{\rm min}^k = log(k)</math>. В результате, <math>\Psi_{\rm distance}^d = O(log(d))</math>.  Таким образом:
 +
</p>
 +
<p align="center"><math>\Psi_{\rm distribute}^{k, d}  = \Psi_{\rm distance}^d + \Psi_{\rm min}^k = O(log(d))+O(log(k)) = O(log(kd))</math></p>
 +
 +
 +
<p><b>Пересчет центров кластеров в d-мерном пространстве</b></p>
 +
<p>
 +
Поскольку на данном шаге для каждого из <math>k</math> кластеров центр масс может быть вычислен независимо, данные операции могут быть выполнены в отдельных потоках. Таким образом, параллельная сложность данного шага, <math>\Psi_{\rm recenter}^{k, d}</math>,  будет определяться <i>параллельной сложностью вычисления одного центра масс кластера размера <math>m</math></i>, <math>\Psi_{\rm recenter}^{k, d}</math>,  а так как <math>m \le n \Rightarrow \Psi_{\rm recenter}^{d, m} \le \Psi_{\rm recenter}^{d, n}</math>. Сложность вычисления центра масс кластера d-мерных векторов размера n аналогично предыдущим вычислениям равна <math>O(log(n))</math>. Тогда: 
 +
</p>
 +
<p align="center"><math>\Psi_{\rm recenter}^{k, d}  \le \Psi_{\rm recenter}^{d, n} = O(log(n))</math></p>
 +
 +
 +
<p><b>Общая параллельная сложность алгоритма</b></p>
 +
<p>
 +
На каждой итерации необходимо обновление центров кластеров, которые будут использованы на следующей итерации. Таким образом, итерационный процесс выполняется последовательно<ref>Zhao, Weizhong, Huifang Ma, and Qing He. "Parallel k-means clustering based on mapreduce." IEEE International Conference on Cloud Computing. Springer Berlin Heidelberg, 2009.</ref>. Тогда, поскольку сложность каждой итерации  определяется <math>\Psi_{\rm distribute}^{k, d}</math> и <math>\Psi_{\rm recenter}</math>, сложность всего алгоритма, <math>\Psi_{\rm k-means}</math> в предположении, что было сделано <math>i</math> операций определяется выражением
 +
</p>
 +
<p align="center"><math>\Psi_{\rm k-means} \approx i \cdot (\Psi_{\rm distribute}^{k, d} + \Psi_{\rm recenter}^{k, d}) \le i \cdot O(log(kdn))</math></p>
 +
 
=== Входные и выходные данные алгоритма ===
 
=== Входные и выходные данные алгоритма ===
 +
<p>
 +
<b>Входные данные</b><br>
 +
<ul>
 +
  <li>Матрица из <math>n \cdot d</math> элементов <math>x_{i, j} \in \mathbb{R}, i=\overline{1, n}, j=\overline{1, d}</math>&ndash; координат векторов (наблюдений).</li>
 +
  <li>Целое положительное число <math>k, k \le n</math> &ndash; количество кластеров.</li>
 +
</ul>
 +
</p>
 +
<p>
 +
<b>Объем входных данных</b><br>
 +
<math>1</math> целое число + <math>n \cdot d</math> вещественных чисел (при условии, что координаты &ndash; вещественные числа).
 +
</p>
 +
<p>
 +
<b>Выходные данные</b><br>
 +
<math>n</math> целых положительных чисел <math>L_1, ..., L_n</math>&ndash; номера кластеров, соотвествующие каждому вектору (при условии, что нумерация кластеров начинается с <math>1</math>).
 +
</p>
 +
<p>
 +
<b>Объем выходных данных</b><br>
 +
<math>n</math> целых положительных чисел.
 +
</p>
 +
 
=== Свойства алгоритма ===
 
=== Свойства алгоритма ===
 +
<p><b>[[глоссарий#Вычислительная мощность|''Вычислительная мощность'']]</b></p>
 +
Вычислительная мощность алгоритма <i>k means</i> равна <math>\frac{ikdn}{nd} = ki </math>, где <math>k</math> &ndash; число кластеров, <math>i</math> &ndash; число итераций алгоритма.
 +
 +
<p><b>[[глоссарий#Детерминированность |''Детерминированность'']] и [[глоссарий#Устойчивость |''Устойчивость'']] </b></p>
 +
Алгоритм <i>k means</i> является итерационным. Количество итераций алгоритма в общем случае не фиксируется и зависит от начального расположения объектов в пространстве, параметра <math>k</math>, а также от начального приближения центров кластеров, <math>\mu_1, ..., \mu_k</math>. В результате этого может варьироваться результат работы алгоритма. При неудачном выборе начальных параметров итерационный процесс может сойтись к локальному оптимуму<ref>Von Luxburg, Ulrike. Clustering Stability. Now Publishers Inc, 2010.</ref>. По этим причинам алгоритм не является ни <b>детермирированным</b>, ни <b>устойчивым</b>.
 +
 +
<p><b>Соотношение последовательной и параллельной сложности алгоритма</b></p>
 +
 +
<math>\frac{\Theta_{\rm k-means}}{\Psi_{\rm k-means}} = \frac{O(ikdn)}{O(i \cdot log(kdn))}</math>
 +
 +
 +
 +
<b>Сильные стороны алгоритма</b>:
 +
<ul>
 +
  <li><i>Сравнительно высокая эффективность при простоте реализации</i></li>
 +
  <li><i>Высокое качество кластеризации</i></li>
 +
  <li><i>Возможность распараллеливания</i></li>
 +
  <li><i>Существование множества модификаций</i></li>
 +
</ul>
 +
<b>Недостатки алгоритма</b><ref>Ortega, Joaquín Pérez, Ma Del Rocío Boone Rojas, and María J. Somodevilla. "Research issues on K-means Algorithm: An Experimental Trial Using Matlab."</ref>:
 +
<ul>
 +
  <li><i>Количество кластеров является параметром алгоритма</i></li>
 +
  <li>
 +
    <p><i>Чувствительность к начальным условиям</i></p>
 +
    <p>Инициализация центров кластеров в значительной степени влияет на результат кластеризации.</p>
 +
  </li>
 +
  <li>
 +
    <p><i>Чувствительность к выбросам и шумам</i></p>
 +
    <p>Выбросы, далекие от центров настоящих кластеров, все равно учитываются при вычислении их центров.</p>
 +
  </li>
 +
  <li>
 +
    <p><i>Возможность сходимости к локальному оптимуму</i></p>
 +
    <p>Итеративный подход не дает гарантии сходимости к оптимальному решению.</p>
 +
  </li>
 +
  <li>
 +
    <p><i>Использование понятия "среднего"</i></p>
 +
    <p>Алгоритм неприменим к данным, для которых не определено понятие "среднего", например, категориальным данным.</p>
 +
  </li>
 +
</ul>
 +
 
== Программная реализация алгоритма ==
 
== Программная реализация алгоритма ==
 
=== Особенности реализации последовательного алгоритма ===
 
=== Особенности реализации последовательного алгоритма ===
Строка 44: Строка 309:
 
=== Возможные способы и особенности параллельной реализации алгоритма ===
 
=== Возможные способы и особенности параллельной реализации алгоритма ===
 
=== Масштабируемость алгоритма и его реализации ===
 
=== Масштабируемость алгоритма и его реализации ===
==== Масштабируемость алгоритма ====
+
 
==== Масштабируемость реализации алгоритма ====
+
Исследование масштабируемости параллельной реализации алгоритма k means проводилось на суперкомпьютере "Ломоносов"<ref name="Lom">Воеводин Вл., Жуматий С., Соболев С., Антонов А., Брызгалов П., Никитенко Д., Стефанов К., Воеводин Вад. Практика суперкомпьютера «Ломоносов» // Открытые системы, 2012, N 7, С. 36-39.</ref> [http://parallel.ru/cluster Суперкомпьютерного комплекса Московского университета]. Алгоритм реализован на языке C с использованием средств MPI.
 +
Для исследования масштабируемости проводилось множество запусков программы с разным значением параметра (количество векторов для кластеризации), а также с различным числом процессоров. Фиксировались результаты запусков &ndash; время работы <math>t</math> и количество произведенных итераций алгоритма <math>i</math>.
 +
 
 +
Параметры запусков для экспериментальной оценки:
 +
<ul>
 +
  <li>Значения <b>количества векторов</b> <math>n</math>: 20'000, 30'000, 50'000, 100'000, 200'000, 300'000, 500'000, 700'000, 1'000'000, 1'500'000, 2'000'000.</li>
 +
  <li>Значения <b>количества процессоров</b> <math>p</math>: 1, 8, 16, 32, 64, 128, 160, 192, 224, 256, 288, 320, 352, 384, 416, 448, 480, 512.</li>
 +
  <li>Значение <b>количества кластеров</b> <math>k</math>: 100.</li>
 +
  <li>Значение <b>размерности векторов</b> <math>d</math>: 10.</li>
 +
</ul>
 +
 
 +
 
 +
 
 +
Для проведения экспериментов были сгенерированы нормально распределенные псевдослучайные данные (с использованием Python библиотеки [http://scikit-learn.org/ scikit-learn]):
 +
<syntaxhighlight lang="python">
 +
from sklearn.datasets import make_classification
 +
X1, Y1 = make_classification(n_features=10, n_redundant=2, n_informative=8,
 +
                            n_clusters_per_class=1, n_classes=100,
 +
                            n_samples=2000000)
 +
</syntaxhighlight>
 +
 
 +
 
 +
 
 +
Для заданной конфигурации эксперимента (<math>n, d, p, k</math>) и полученных результатов (<math>t, i</math>)  [[Глоссарий#Производительность|производительность]]  и эффективность реализации расчитывались по формулам:
 +
<p>
 +
<ul><li><math>{\rm Performance} = \frac{N_{\rm k-means}^{d, n}}{t}\ \ (FLOPS),</math></li></ul>
 +
где <math>N_{\rm k-means}^{d, n}</math> &ndash; точное число операций с плавающей точкой (операции с памятью, а также целочисленные операции не учитывались), вычисленное в соответствие с разделом [[#Последовательная сложность алгоритма| "Последовательная сложность алгоритма"]];
 +
</p>
 +
<p>
 +
<ul><li><math>{\rm Efficiency} = \frac{100 \cdot {\rm Performance}}{{\rm Performance}_{\rm Peak}^{p}}\ \ (%),</math></li></ul>
 +
где <math>{\rm Performance}_{\rm Peak}^{p}</math> &ndash; пиковая производительность суперкомпьютера при <math>p</math> процессорах, вычисленная согласно спецификациям Intel<sup>&reg;</sup> XEON<sup>&reg;</sup> X5670<ref>"http://ark.intel.com/ru/products/47920/Intel-Xeon-Processor-X5670-12M-Cache-2_93-GHz-6_40-GTs-Intel-QPI"</ref>.
 +
</p>
 +
 
 +
 
 +
 
 +
Графики зависимости производительности и эффективности параллельной реализации k means от числа векторов для кластеризации (<math>n</math>) и числа процессоров (<math>p</math>) представлены на рисунках 4 и 5 соответственно.
 +
 
 +
[[file:Flops.png|thumb|center|800px|Рис. 4. Параллельная реализация k means. График зависимости производительности реализации алгоритма от числа векторов для кластеризации (<math>n</math>) и числа процессоров (<math>p</math>).]]
 +
[[file:Efficiency.png|thumb|center|800px|Рис. 5. Параллельная реализация k means. График зависимости эффективности реализации алгоритма от числа векторов для кластеризации (<math>n</math>) и числа процессоров (<math>p</math>).]]
 +
 
 +
<p>
 +
<b>В результате</b> экспериментальной оценки были получены следующие оценки эффективности реализации:
 +
<ul>
 +
  <li><b>Минимальное</b> значение: 0.000409 % достигается при <math>n=20'000, p=480</math></li>
 +
  <li><b>Максимальное</b> значение: 0.741119 % достигается при <math>n=300'000, p=1</math></li>
 +
</ul>
 +
</p>
 +
 
 +
<p>
 +
Оценки масштабируемости реализации алгоритма k means:
 +
<ul>
 +
  <li><b>По числу процессоров:</b> -0.002683 &ndash; эффективность убывает с ростом числа процессоров. Данный результат вызван ростом накладных расходов для обеспечения параллельного выполнения алгоритма.</li>
 +
  <li><b>По размеру задачи:</b> 0.002779 &ndash; эффективность растет с ростом числа векторов. Данный результат вызван тем, что при увеличении размера задачи, количество вычислений растет по сравнению с временем, затрачиваемым на пересылку данных.</li>
 +
  <li><b>Общая оценка:</b> -0.000058 &ndash; можно сделать вывод, что в целом эффективность реализации незначительно уменьшается с ростом размера задачи и числа процессоров.</li>
 +
</ul>
 +
</p>
 +
 
 +
 
 +
[https://github.com/serban/kmeans Использованная параллельная реализация алгоритма k means]
 +
 
 
=== Динамические характеристики и эффективность реализации алгоритма ===
 
=== Динамические характеристики и эффективность реализации алгоритма ===
 
=== Выводы для классов архитектур ===
 
=== Выводы для классов архитектур ===
 +
 +
В однопоточном режиме на наборах данных, представляющих практический интерес (порядка нескольких десятков тысяч векторов и выше), время работы алгоритма неприемлемо велико. Благодаря свойству массового параллелизма должно наблюдаться значительное ускорение алгоритма на многоядерных архитектурах (Intel Xeon), а также на графических процессорах, даже на мобильных вычислительных системах (ноутбуках), оснащенных видеокартой. Также алгоритм k means будет демонстрировать значительное ускорение на сверхмощных вычислительных комплексах (суперкомпьютерах, системах облачных вычислений<ref>"Issa. "Performance characterization and analysis for Hadoop K-means iteration". Journal of Cloud Computing, 2016"</ref>).
 +
 +
На сегодняшний день существует множество реализаций алгоритма k means, в частности, направленных на оптимизацию параллельной работы на различных архитектурах<ref>"Raghavan R. A fast and scalable hardware architecture for K-means clustering for big data analysis : дис. – University of Colorado Colorado Springs. Kraemer Family Library, 2016."</ref><ref>"Yang, Luobin, et al. "High performance data clustering: a comparative analysis of performance for GPU, RASC, MPI, and OpenMP implementations." The Journal of supercomputing 70.1 (2014): 284-300."</ref><ref>"Li, You, et al. "Speeding up k-means algorithm by GPUs." Computer and Information Technology (CIT), 2010 IEEE 10th International Conference on. IEEE, 2010."</ref>. Предлагается множество адаптаций алгоритма под конкретные архитектуры. Например, авторы работы<ref>"Kanan, Gebali, Ibrahim. "Fast and Area-Efficient Hardware Implementation of the K-means
 +
Clustering Algorithm". WSEAS Transactions on circuits and systems. Vol. 15. 2016"</ref> производят перерасчет центров кластеров на этапе распределения векторов по кластерам.
 +
 
=== Существующие реализации алгоритма ===
 
=== Существующие реализации алгоритма ===
 +
==== Открытое программное обеспечение ====
 +
<ol>
 +
  <li>[https://www.icpsr.umich.edu/CrimeStat/ CrimeStat]</li> Программное обеспечение, созданное для операционных систем Windows, предоставляющее инструменты статистического  и пространственного анализа для решения задачи картирования преступности.
 +
  <li>[http://juliastats.github.io Julia]</li> Высокоуровневый высокопроизводительный свободный язык программирования с динамической типизацией, созданный для математических вычислений, содержит реализацию k-means.
 +
  <li>[https://mahout.apache.org Mahout]</li> Apache Mahout - Java библиотека для работы с алгоритмами машинного обучения с использованием MapReduce. Содержит реализацию k-means.
 +
  <li>[https://www.gnu.org/software/octave/ Octave]</li> Написанная на C++  свободная система для математических вычислений, использующая совместимый с MATLAB язык высокого уровня, содержит реализацию k-means.
 +
  <li>[http://spark.apache.org/docs/latest/mllib-clustering.html Spark]</li> Распределенная реализация k-means содержится в библиотеке Mlib для работы с алгоритмами машинного обучения, взаимодействующая с Python библиотекой NumPy и библиотека R.
 +
  <li>[http://torch.ch Torch]</li>  MATLAB-подобная библиотека для языка программирования Lua с открытым исходным кодом, предоставляет большое количество алгоритмов для глубинного обучения и научных расчётов. Ядро написано на Си, прикладная часть выполняется на LuaJIT, поддерживается распараллеливание вычислений средствами CUDA и OpenMP. Существуют реализации k-means.
 +
  <li>[http://www.cs.waikato.ac.nz/ml/weka/ Weka]</li>  Cвободное программное обеспечение для анализа данных, написанное на Java. Содержит k-means и x-means.
 +
  <li>[http://accord-framework.net/docs/html/T_Accord_MachineLearning_KMeans.htm Accord.NET]</li> C# реализация алгоритмов k-means, k-means++, k-modes.
 +
  <li>[http://docs.opencv.org/3.0-beta/doc/py_tutorials/py_ml/py_kmeans/py_kmeans_opencv/py_kmeans_opencv.html OpenCV]</li> Написанная на С++ библиотека, направленная в основном на решение задач компьютерного зрения. Содержит реализацию k-means.
 +
  <li>[http://mlpack.org/ MLPACK]</li> Масштабируемая С++ библиотека для работы с алгоритмами машинного обучения, содержит реализацию k-means.
 +
  <li>[https://www.scipy.org/ SciPy]</li> Библиотека Python, содержит множество реализаций k-means.
 +
  <li>[http://scikit-learn.org/ scikit-learn]</li> Библиотека Python, содержит множество реализаций k-means.
 +
  <li>[https://www.r-bloggers.com/k-means-clustering-in-r/ R]</li> Язык программирования для статистической обработки данных и работы с графикой, а также свободная программная среда вычислений с открытым исходным кодом в рамках проекта GNU, содержит три реализации k-means.
 +
  <li>[http://elki.dbs.ifi.lmu.de ELKI]</li> Java фреймворк, содержащий реализацию k-means, а также множество других алгоритмов кластеризации.
 +
</ol>
 +
==== Проприетарное программное обеспечение  ====
 +
 +
<ol>
 +
  <li>[https://www.ayasdi.com/blog/bigdata/topological-data-analysis-of-oil-and-gas-petrophysical-data/ Ayasdi]</li>
 +
  <li>[http://www.stata.com/features/cluster-analysis/ Stata]</li>
 +
  <li>[http://mathworld.wolfram.com/K-MeansClusteringAlgorithm.html Mathematica]</li>
 +
  <li>[http://www.mathworks.com/help/stats/kmeans.html?requestedDomain=www.mathworks.com MATLAB]</li>
 +
  <li>[https://support.sas.com/rnd/app/stat/procedures/fastclus.html SAS]</li>
 +
  <li>[http://docs.rapidminer.com/studio/operators/modeling/segmentation/k_means.html RapidMiner]</li>
 +
  <li>[https://blogs.sap.com/2013/03/28/sap-hana-pal-k-means-algorithm-or-how-to-do-customer-segmentation-for-the-telecommunications-industry/ SAP HANA]</li>
 +
</ol>
 
== Литература ==
 
== Литература ==
 +
<references />

Текущая версия на 17:18, 19 декабря 2016

Symbol confirmed.svgЭта работа успешно выполнена
Преподавателю: в основное пространство, в подстраницу

Данное задание было проверено и зачтено.
Проверено Konshin и Zhum.



Алгоритм k средних (k means)
Последовательный алгоритм
Последовательная сложность [math]O(ikdn)[/math]
Объём входных данных [math] dn [/math]
Объём выходных данных [math] n [/math]


Авторы страницы Данута Бротиковская и Денис Зобнин.

Оба автора внесли одинаковый вклад в написание статьи (поиск материалов, изучение статей, заполнение контента статьи), четких разграничений по пунктам статьи сделать невозможно.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Алгоритм k средних (k means) – неиерархический[1], итерационный метод кластеризации[2], получивший большую популярность благодаря своей простоте, наглядности реализации и достаточно высокому качеству работы. Был изобретен в 1950-х годах математиком Гуго Штейнгаузом[3] и почти одновременно Стюартом Ллойдом[4]. Особую популярность приобрел после публикации работы МакКуина[5] в 1967.

Алгоритм представляет собой версию EM-алгоритма[6], применяемого также для разделения смеси гауссиан. Основная идея k means заключается в том, что данные произвольно разбиваются на кластеры, после чего итеративно перевычисляется центр масс для каждого кластера, полученного на предыдущем шаге, затем векторы разбиваются на кластеры вновь в соответствии с тем, какой из новых центров оказался ближе по выбранной метрике.

Цель алгоритма заключается в разделении [math]n[/math] наблюдений на [math]k[/math] кластеров таким образом, чтобы каждое наблюдение принадлежало ровно одному кластеру, расположенному на наименьшем расстоянии от наблюдения.

1.2 Математическое описание алгоритма

Дано:

  • набор из [math]n[/math] наблюдений [math]X=\{\mathbf{x}_1, \mathbf{x}_2, ..., \mathbf{x}_n\}, \mathbf{x}_i \in \mathbb{R}^d, i=\overline{1, n}[/math];
  • [math]k[/math] - требуемое число кластеров, [math]k \in \mathbb{N}, k \leq n[/math].

Требуется:

Разделить множество наблюдений [math]X[/math] на [math]k[/math] кластеров [math]S_1, S_2, ..., S_k[/math]:

  • [math]S_i \cap S_j= \varnothing, \quad i \ne j[/math]
  • [math]\bigcup_{i=1}^{k} S_i = X[/math]


Действие алгоритма:

Алгоритм k means разбивает набор [math]X[/math] на [math]k[/math] наборов [math]S_1, S_2, ..., S_k,[/math] таким образом, чтобы минимизировать сумму квадратов расстояний от каждой точки кластера до его центра (центр масс кластера). Введем обозначение, [math]S=\{S_1, S_2, ..., S_k\}[/math]. Тогда действие алгоритма k means равносильно поиску:

[math]\arg\min_{S} \sum\limits_{i=1}^k \sum\limits_{\mathbf{x} \in S_i} \rho(\mathbf{x}, \mathbf{\mu}_i )^2,[/math] [math](1)[/math]

где [math]\mathbf{\mu}_i[/math] – центры кластеров, [math]i=\overline{1,k}, \quad \rho(\mathbf{x}, \mathbf{\mu}_i)[/math] – функция расстояния между [math]\mathbf{x}[/math] и [math]\mu_i[/math]


Шаги алгоритма:

  1. Начальный шаг: инициализация кластеров

    Выбирается произвольное множество точек [math]\mu_i, i=\overline{1,k}[/math], рассматриваемых как начальные центры кластеров: [math]\mu_i^{(0)} = \mu_i, \quad i=\overline{1, k}[/math]

  2. Распределение векторов по кластерам

    Шаг [math]t: \forall \mathbf{x}_i \in X, i=\overline{1,n}: \mathbf{x}_i \in S_j \iff j=\arg\min_{k}\rho(\mathbf{x}_i,\mathbf{\mu}_k^{(t-1)})^2[/math]

  3. Пересчет центров кластеров

    Шаг [math]t: \forall i=\overline{1,k}: \mu_i^{(t)} = \cfrac{1}{|S_i|}\sum_{\mathbf{x}\in S_i}\mathbf{x}[/math]

  4. Проверка условия останова:

       if [math]\exist i\in \overline{1,k}: \mu_i^{(t)} \ne \mu_i^{(t-1)}[/math] then
           t = t + 1;
           goto 2;
       else
           stop
    

1.3 Вычислительное ядро алгоритма

Вычислительным ядром являются шаги 2 и 3 приведенного выше алгоритма: распределение векторов по кластерам и пересчет центров кластеров.

Распределение векторов по кластерам предполагает вычисление расстояний между каждым вектором [math]\mathbf{x}_i \in X, i= \overline{1,n}[/math] и центрами кластера [math]\mathbf{\mu}_j, j= \overline{1,k}[/math]. Таким образом, данный шаг предполагает [math]kn[/math] вычислений расстояний между d-мерными векторами.

Пересчет центров кластеров предполагает [math]k[/math] вычислений центров масс [math]\mathbf{\mu}_i[/math] множеств [math]S_i, i=\overline{1,k}[/math] представленных выражением в шаге 3 представленного выше алгоритма.

1.4 Макроструктура алгоритма

Инициализация центров масс [math]\mu_1, ..., \mu_k[/math].

Наиболее распространенными являются следующие стратегии:

  • Метод Forgy
    В качестве начальных значений [math]\mu_1, ..., \mu_k[/math] берутся случайно выбранные векторы.
  • Метод случайно разделения (Random Partitioning)
    Для каждого вектора [math]\mathbf{x}_i \in X, i=\overline{1, n}[/math] выбирается случайным образом кластер [math]S_1, ..., S_k[/math], после чего для каждого полученного кластера вычисляются значения [math]\mu_1, ..., \mu_k[/math].


Распределение векторов по кластерам

Для этого шага алгоритма между векторами [math]\mathbf{x}_i \in X, i=\overline{1, n}[/math] и центрами кластеров [math]\mu_1, ..., \mu_k[/math] вычисляются расстояния по формуле (как правило, используется Евлидово расстояние):

[math] \mathbf{v}_1, \mathbf{v}_2 \in \mathbb{R}^d, \quad \rho(\mathbf{v}_1, \mathbf{v}_2) = \lVert \mathbf{v}_1- \mathbf{v}_2 \rVert= \sqrt{\sum_{i=1}^{d}(\mathbf{v}_{1,i} - \mathbf{v}_{2,i})^2}[/math] [math](2)[/math]


Пересчет центров кластеров

Для этого шага алгоритма производится пересчет центров кластера по формуле вычисления центра масс:

[math] \mu = \cfrac{1}{|S|}\sum_{\mathbf{x}\in S}\mathbf{x}[/math] [math](3)[/math]

1.5 Схема реализации последовательного алгоритма

algorithm k-means is
    1. Инициализировать центры кластеров [math]\mathbf{\mu}_i^{(1)}, i=\overline{1,k}[/math]
    2. [math]t \leftarrow 1[/math]
    3. Распределение по кластерам
[math]S_i^{(t)}=\{\mathbf{x}_p: \lVert\mathbf{x}_p-\mathbf{\mu}_i^{(t)}\rVert^2 \leq \lVert\mathbf{x}_p-\mathbf{\mu}_j^{(t)}\rVert^2 \quad \forall j=\overline{1,k}\},[/math]
где каждый вектор [math]\mathbf{x}_p[/math] соотносится единственному кластеру [math]S^{(t)}[/math] 4. Обновление центров кластеров
[math]\mathbf{\mu}_i^{(t+1)} = \frac{1}{|S^{(t)}_i|} \sum_{\mathbf{x}_j \in S^{(t)}_i} \mathbf{x}_j [/math]
5. if [math]\exist i \in \overline{1,k}: \mathbf{\mu}_i^{(t+1)} \ne \mathbf{\mu}_i^{(t)}[/math] then t = t + 1; goto 3; else stop

1.6 Последовательная сложность алгоритма

Обозначим [math]\Theta_{\rm centroid}^{d, m}[/math] временную сложность вычисления центорида кластера, число элементов которого равна [math]m[/math], в d-мерном пространстве.

Аналогично [math]\Theta_{\rm distance}^d[/math] – временная сложность вычисления расстояния между двумя d-мерными векторами.

Сложность шага инициализации [math]k[/math] кластеров мощности [math]m[/math] в d-мерном пространстве – [math]\Theta_{\rm init}^{k, d, m}[/math]

  • Стратерия Forgy: вычисления не требуются, [math]\Theta_{\rm init}^{k, d, m} = 0[/math]
  • Стратегия случайного разбиения: вычисление центров [math]k[/math] кластеров, [math]\Theta_{\rm init}^{k, d, m} = k \cdot \Theta_{\rm centroid}^{d, m}, m \le n[/math]

Cложность шага распределения d мерных векторов по [math]k[/math] кластерам – [math]\Theta_{\rm distribute}^{k, d}[/math]

На этом шаге для каждого вектора [math]\mathbf{x}_i \in X, i=\overline{1, n}[/math] вычисляется k расстояний до центров кластеров [math]\mathbf{\mu}_1, ...\mathbf{\mu}_k[/math]

[math]\Theta_{\rm distribute}^{k, d} = n \cdot k \cdot \Theta_{\rm distance}^d[/math]

Сложность шага пересчета центров [math]k[/math] кластеров размера [math]m[/math] в d-мерном пространстве – [math]\Theta_{\rm recenter}^{k, d, m}[/math]

На этом шаге вычисляется [math]k[/math] центров кластеров [math]\mathbf{\mu}_1, ...\mathbf{\mu}_k[/math]

[math]\Theta_{\rm recenter}^{k, d, m} = k \cdot \Theta_{\rm centroid}^{d, m}[/math]

Рассчитаем [math]\Theta_{\rm centroid}^{d, m}[/math] для кластера, число элементов которого равно [math]m[/math]

[math]\Theta_{\rm centroid}^{d, m}[/math] = [math]m \cdot d[/math] сложений + [math]d[/math] делений

Рассчитаем [math]\Theta_{\rm distance}^d[/math] в соответствие с формулой [math](2)[/math]

[math]\Theta_{\rm distance}^d[/math] = [math]d[/math] вычитаний + [math]d[/math] умножений + [math](d-1)[/math] сложение

Предположим, что алгоритм сошелся за [math]i[/math] итераций, тогда временная сложность алгоритма [math]\Theta_{\rm k-means}^{d, n}[/math]

[math]\Theta_{\rm k-means}^{d, n} \le \Theta_{\rm init}^{k, d, n} + i(\Theta_{\rm distribute}^{k, d} + \Theta_{\rm recenter}^{k, d, n})[/math]

Операции сложения/вычитания:

[math]\Theta_{\rm k-means}^{d, n} \le knd+ i(kn(2d-1) + knd) = knd+ i(kn(3d-1)) \thicksim O(ikdn)[/math]

Операции умножения/деления:

[math]\Theta_{\rm k-means}^{d, n} \le kd + i(knd + kd) = kd + ikd(n+1) \thicksim O(ikdn) [/math]

Получаем, что временная сложность алгоритма k means кластеризации [math]n[/math] d-мерных векторов на [math]k[/math] кластеров за [math]i[/math] итераций:

[math] \Theta_{\rm k-means}^{d, n} \thicksim O(ikdn) [/math]

1.7 Информационный граф

Рассмотрим информационный граф алгоритма. Алгоритм k means начинается с этапа инициализации, после которого следуют итерации, на каждой из которых выполняется два последовательных шага (см. "Схема реализации последовательного алгоритма"):

  • распределение векторов по кластерам
  • перерасчет центров кластеров

Поскольку основная часть вычислений приходится на шаги итераций, распишем информационные графы данных шагов.


Распределение векторов по кластерам

Информационный граф шага распределения векторов по кластерам представлен на рисунке 1. Исходами данного графа является исходные векторы [math]\mathbf{x}_1, ..., \mathbf{x}_n[/math], а также центры кластеров [math]\mathbf{\mu}_1, ... \mathbf{\mu}_k[/math], вычисленные ранее (на шаге инициализации, если рассматривается первая итерация алгоритма, или на шаге пересчета центров кластеров предыдущей итерации в противном случае). Каждая пара векторов данных [math]\mathbf{x}_i, i=\overline{1,n}[/math] и центров кластера [math]\mathbf{\mu}_j, j=\overline{1, k}[/math] : ([math]\mathbf{x}_i[/math], [math]\mathbf{\mu}_j[/math]) подаются на независимые узлы "d" вычисления расстояния между векторами (более подробная схема вычисления расстояния представлена далее, рисунок 2). Далее узлы вычисления расстояния "d", соответствующие одному и тому же исходному вектору [math]\mathbf{x}_i[/math] передаются на один узел "m", где далее происходит вычисление новой метки кластера для каждого вектора [math]\mathbf{x}_i[/math] (берется кластер с минимальным результатом вычисления расстояния). На выходе графа выдаются метки кластеров , [math]L_1, ..., L_n[/math], такие что [math]\forall \mathbf{x}_i, i=\overline{1, n}, \mathbf{x}_i \in S_j \Leftrightarrow L_i = j[/math].

Рис. 1. Схема распределения векторов по кластерам. d – вычисление расстояния между векторами; m – вычисление минимума.


Вычисление расстояния между векторами

Подробная схема вычисления расстояния между векторами [math]\mathbf{x}_i, \mathbf{\mu}_j[/math] представлена на рисунке 2. Как показано на графе, узел вычисления расстояния между векторами "d" состоит из шага взятия разности между векторами (узел "[math]-[/math]") и взятия нормы получившегося вектора разности (узел "[math]||\cdot||^2[/math]"). Более подробно, вычисление расстояния между векторами [math]\mathbf{x}_i = {x_{i1,}, ...,{x_{in}}}, \mathbf{\mu}_j = {\mu_{j1}, ...,\mu_{jn}}[/math] может быть представлено как вычисление разности между каждой парой компонент [math](x_{iz}, \mu_{jz}), z=\overline{1, d}[/math] (узел "[math]-[/math]"), далее возведение в квадрат для каждого узла "[math]-[/math]" (узел "[math]()^2[/math]") и суммирования выходов всех узлов "[math]()^2[/math]" (узел "[math]+[/math]").

Рис. 2. Схема вычисления расстояния между вектором и центром кластера.

Пересчет центров кластеров

Информационный граф шага пересчета центров кластеров представлен на рисунке 3. Исходами данного графа является исходные векторы [math]\mathbf{x}_1, ..., \mathbf{x}_n[/math], а также им соответствующие метки кластера, [math]L_1, ..., L_n[/math], такие что [math]\forall x_i, i=\overline{1, n}, \mathbf{x}_i \in S_j \Leftrightarrow L_i = j[/math], вычисленные на этапе распределения векторов по кластерам. Все векторы [math]\mathbf{x}_1, ..., \mathbf{x}_n[/math] подаются в узлы [math]+_1, ...+_k[/math], каждый узел [math]+_m, m = \overline{1,k}[/math] соответствует операции сложения векторов кластера с номером [math]m[/math]. Метки кластера [math]L_1, ..., L_n[/math] также совместно передаются на узлы [math]S_m, m = \overline{1,k}[/math], на каждом из которых вычисляется количество векторов в соответствующем кластере (количество меток с соответствующим значением). Далее каждая пара выходов узлов [math]+_m[/math] и [math]S_m[/math] подается на узел "[math]/[/math]", где производится деление суммы векторов кластера на количество элементов в нем. Значения, вычисленные на узлах "[math]/[/math]", присваиваются новым центрам кластеров (выходные значения графа).

Рис. 3. Схема пересчета центров кластеров

1.8 Ресурс параллелизма алгоритма

Работа алгоритма состоит из [math]i[/math] итераций, в каждой из которых происходит распределение d-мерных векторов по [math]k[/math] кластерам, а также пересчет центров кластеров в d-мерном пространстве. В шаге распределения d-мерных векторов по [math]k[/math] кластерам расстояния между вектором и центрами кластеров вычисляются независимо (отсутствуют информационные зависимости). Центры масс кластеров также пересчитываются независимо друг от друга. Таким образом, имеет место массовый параллелизм. Вычислим параллельную сложность [math]\Psi_*[/math] каждого из шагов, а также параллельную сложность всего алгоритма, [math]\Psi_{\rm k-means}[/math]. Будем исходить из предположения, что может быть использовано любое необходимое число потоков.


Распределение d-мерных векторов по [math]k[/math] кластерам

Поскольку на данном шаге для каждой пары векторов [math]\mathbf{x}_i, i=\overline{1,n}[/math] и [math]\mathbf{\mu}_j, j=\overline{1,k}[/math] операции вычисления расстояния не зависят друг от друга, они могут выполняться параллельно. Тогда, разделив все вычисление расстояний на [math]n[/math] потоков, получим, что в каждом потоке будет выполняться только одна операция вычисления расстояния между векторами размерности d. При этом каждому вычислительному потоку передаются координаты центров всех кластеров [math]\mathbf{\mu}_1, ..., \mathbf{\mu}_k[/math]. Таким образом, параллельная сложность данного шага определяется сложностью параллельной операции вычисления расстояния между d-мерными векторами, [math]\Psi_{\rm distance}^d[/math] и сложностью определения наиболее близкого кластера (паралельное взятие минимума по расстояниям), [math]\Psi_{\rm min}^k[/math]. Для оценки [math]\Psi_{\rm distance}^d[/math] воспользуемся параллельной реализацией нахождения частичной суммы элементов массива путем сдваивания. Аналогично, [math]\Psi_{\rm min}^k = log(k)[/math]. В результате, [math]\Psi_{\rm distance}^d = O(log(d))[/math]. Таким образом:

[math]\Psi_{\rm distribute}^{k, d} = \Psi_{\rm distance}^d + \Psi_{\rm min}^k = O(log(d))+O(log(k)) = O(log(kd))[/math]


Пересчет центров кластеров в d-мерном пространстве

Поскольку на данном шаге для каждого из [math]k[/math] кластеров центр масс может быть вычислен независимо, данные операции могут быть выполнены в отдельных потоках. Таким образом, параллельная сложность данного шага, [math]\Psi_{\rm recenter}^{k, d}[/math], будет определяться параллельной сложностью вычисления одного центра масс кластера размера [math]m[/math], [math]\Psi_{\rm recenter}^{k, d}[/math], а так как [math]m \le n \Rightarrow \Psi_{\rm recenter}^{d, m} \le \Psi_{\rm recenter}^{d, n}[/math]. Сложность вычисления центра масс кластера d-мерных векторов размера n аналогично предыдущим вычислениям равна [math]O(log(n))[/math]. Тогда:

[math]\Psi_{\rm recenter}^{k, d} \le \Psi_{\rm recenter}^{d, n} = O(log(n))[/math]


Общая параллельная сложность алгоритма

На каждой итерации необходимо обновление центров кластеров, которые будут использованы на следующей итерации. Таким образом, итерационный процесс выполняется последовательно[7]. Тогда, поскольку сложность каждой итерации определяется [math]\Psi_{\rm distribute}^{k, d}[/math] и [math]\Psi_{\rm recenter}[/math], сложность всего алгоритма, [math]\Psi_{\rm k-means}[/math] в предположении, что было сделано [math]i[/math] операций определяется выражением

[math]\Psi_{\rm k-means} \approx i \cdot (\Psi_{\rm distribute}^{k, d} + \Psi_{\rm recenter}^{k, d}) \le i \cdot O(log(kdn))[/math]

1.9 Входные и выходные данные алгоритма

Входные данные

  • Матрица из [math]n \cdot d[/math] элементов [math]x_{i, j} \in \mathbb{R}, i=\overline{1, n}, j=\overline{1, d}[/math]– координат векторов (наблюдений).
  • Целое положительное число [math]k, k \le n[/math] – количество кластеров.

Объем входных данных
[math]1[/math] целое число + [math]n \cdot d[/math] вещественных чисел (при условии, что координаты – вещественные числа).

Выходные данные
[math]n[/math] целых положительных чисел [math]L_1, ..., L_n[/math]– номера кластеров, соотвествующие каждому вектору (при условии, что нумерация кластеров начинается с [math]1[/math]).

Объем выходных данных
[math]n[/math] целых положительных чисел.

1.10 Свойства алгоритма

Вычислительная мощность

Вычислительная мощность алгоритма k means равна [math]\frac{ikdn}{nd} = ki [/math], где [math]k[/math] – число кластеров, [math]i[/math] – число итераций алгоритма.

Детерминированность и Устойчивость

Алгоритм k means является итерационным. Количество итераций алгоритма в общем случае не фиксируется и зависит от начального расположения объектов в пространстве, параметра [math]k[/math], а также от начального приближения центров кластеров, [math]\mu_1, ..., \mu_k[/math]. В результате этого может варьироваться результат работы алгоритма. При неудачном выборе начальных параметров итерационный процесс может сойтись к локальному оптимуму[8]. По этим причинам алгоритм не является ни детермирированным, ни устойчивым.

Соотношение последовательной и параллельной сложности алгоритма

[math]\frac{\Theta_{\rm k-means}}{\Psi_{\rm k-means}} = \frac{O(ikdn)}{O(i \cdot log(kdn))}[/math]


Сильные стороны алгоритма:

  • Сравнительно высокая эффективность при простоте реализации
  • Высокое качество кластеризации
  • Возможность распараллеливания
  • Существование множества модификаций

Недостатки алгоритма[9]:

  • Количество кластеров является параметром алгоритма
  • Чувствительность к начальным условиям

    Инициализация центров кластеров в значительной степени влияет на результат кластеризации.

  • Чувствительность к выбросам и шумам

    Выбросы, далекие от центров настоящих кластеров, все равно учитываются при вычислении их центров.

  • Возможность сходимости к локальному оптимуму

    Итеративный подход не дает гарантии сходимости к оптимальному решению.

  • Использование понятия "среднего"

    Алгоритм неприменим к данным, для которых не определено понятие "среднего", например, категориальным данным.

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

2.2 Локальность данных и вычислений

2.2.1 Локальность реализации алгоритма

2.2.1.1 Структура обращений в память и качественная оценка локальности
2.2.1.2 Количественная оценка локальности

2.3 Возможные способы и особенности параллельной реализации алгоритма

2.4 Масштабируемость алгоритма и его реализации

Исследование масштабируемости параллельной реализации алгоритма k means проводилось на суперкомпьютере "Ломоносов"[10] Суперкомпьютерного комплекса Московского университета. Алгоритм реализован на языке C с использованием средств MPI. Для исследования масштабируемости проводилось множество запусков программы с разным значением параметра (количество векторов для кластеризации), а также с различным числом процессоров. Фиксировались результаты запусков – время работы [math]t[/math] и количество произведенных итераций алгоритма [math]i[/math].

Параметры запусков для экспериментальной оценки:

  • Значения количества векторов [math]n[/math]: 20'000, 30'000, 50'000, 100'000, 200'000, 300'000, 500'000, 700'000, 1'000'000, 1'500'000, 2'000'000.
  • Значения количества процессоров [math]p[/math]: 1, 8, 16, 32, 64, 128, 160, 192, 224, 256, 288, 320, 352, 384, 416, 448, 480, 512.
  • Значение количества кластеров [math]k[/math]: 100.
  • Значение размерности векторов [math]d[/math]: 10.


Для проведения экспериментов были сгенерированы нормально распределенные псевдослучайные данные (с использованием Python библиотеки scikit-learn):

from sklearn.datasets import make_classification
X1, Y1 = make_classification(n_features=10, n_redundant=2, n_informative=8,
                             n_clusters_per_class=1, n_classes=100,
                             n_samples=2000000)


Для заданной конфигурации эксперимента ([math]n, d, p, k[/math]) и полученных результатов ([math]t, i[/math]) производительность и эффективность реализации расчитывались по формулам:

  • [math]{\rm Performance} = \frac{N_{\rm k-means}^{d, n}}{t}\ \ (FLOPS),[/math]

где [math]N_{\rm k-means}^{d, n}[/math] – точное число операций с плавающей точкой (операции с памятью, а также целочисленные операции не учитывались), вычисленное в соответствие с разделом "Последовательная сложность алгоритма";

  • [math]{\rm Efficiency} = \frac{100 \cdot {\rm Performance}}{{\rm Performance}_{\rm Peak}^{p}}\ \ (%),[/math]

где [math]{\rm Performance}_{\rm Peak}^{p}[/math] – пиковая производительность суперкомпьютера при [math]p[/math] процессорах, вычисленная согласно спецификациям Intel® XEON® X5670[11].


Графики зависимости производительности и эффективности параллельной реализации k means от числа векторов для кластеризации ([math]n[/math]) и числа процессоров ([math]p[/math]) представлены на рисунках 4 и 5 соответственно.

Рис. 4. Параллельная реализация k means. График зависимости производительности реализации алгоритма от числа векторов для кластеризации ([math]n[/math]) и числа процессоров ([math]p[/math]).
Рис. 5. Параллельная реализация k means. График зависимости эффективности реализации алгоритма от числа векторов для кластеризации ([math]n[/math]) и числа процессоров ([math]p[/math]).

В результате экспериментальной оценки были получены следующие оценки эффективности реализации:

  • Минимальное значение: 0.000409 % достигается при [math]n=20'000, p=480[/math]
  • Максимальное значение: 0.741119 % достигается при [math]n=300'000, p=1[/math]

Оценки масштабируемости реализации алгоритма k means:

  • По числу процессоров: -0.002683 – эффективность убывает с ростом числа процессоров. Данный результат вызван ростом накладных расходов для обеспечения параллельного выполнения алгоритма.
  • По размеру задачи: 0.002779 – эффективность растет с ростом числа векторов. Данный результат вызван тем, что при увеличении размера задачи, количество вычислений растет по сравнению с временем, затрачиваемым на пересылку данных.
  • Общая оценка: -0.000058 – можно сделать вывод, что в целом эффективность реализации незначительно уменьшается с ростом размера задачи и числа процессоров.


Использованная параллельная реализация алгоритма k means

2.5 Динамические характеристики и эффективность реализации алгоритма

2.6 Выводы для классов архитектур

В однопоточном режиме на наборах данных, представляющих практический интерес (порядка нескольких десятков тысяч векторов и выше), время работы алгоритма неприемлемо велико. Благодаря свойству массового параллелизма должно наблюдаться значительное ускорение алгоритма на многоядерных архитектурах (Intel Xeon), а также на графических процессорах, даже на мобильных вычислительных системах (ноутбуках), оснащенных видеокартой. Также алгоритм k means будет демонстрировать значительное ускорение на сверхмощных вычислительных комплексах (суперкомпьютерах, системах облачных вычислений[12]).

На сегодняшний день существует множество реализаций алгоритма k means, в частности, направленных на оптимизацию параллельной работы на различных архитектурах[13][14][15]. Предлагается множество адаптаций алгоритма под конкретные архитектуры. Например, авторы работы[16] производят перерасчет центров кластеров на этапе распределения векторов по кластерам.

2.7 Существующие реализации алгоритма

2.7.1 Открытое программное обеспечение

  1. CrimeStat
  2. Программное обеспечение, созданное для операционных систем Windows, предоставляющее инструменты статистического и пространственного анализа для решения задачи картирования преступности.
  3. Julia
  4. Высокоуровневый высокопроизводительный свободный язык программирования с динамической типизацией, созданный для математических вычислений, содержит реализацию k-means.
  5. Mahout
  6. Apache Mahout - Java библиотека для работы с алгоритмами машинного обучения с использованием MapReduce. Содержит реализацию k-means.
  7. Octave
  8. Написанная на C++ свободная система для математических вычислений, использующая совместимый с MATLAB язык высокого уровня, содержит реализацию k-means.
  9. Spark
  10. Распределенная реализация k-means содержится в библиотеке Mlib для работы с алгоритмами машинного обучения, взаимодействующая с Python библиотекой NumPy и библиотека R.
  11. Torch
  12. MATLAB-подобная библиотека для языка программирования Lua с открытым исходным кодом, предоставляет большое количество алгоритмов для глубинного обучения и научных расчётов. Ядро написано на Си, прикладная часть выполняется на LuaJIT, поддерживается распараллеливание вычислений средствами CUDA и OpenMP. Существуют реализации k-means.
  13. Weka
  14. Cвободное программное обеспечение для анализа данных, написанное на Java. Содержит k-means и x-means.
  15. Accord.NET
  16. C# реализация алгоритмов k-means, k-means++, k-modes.
  17. OpenCV
  18. Написанная на С++ библиотека, направленная в основном на решение задач компьютерного зрения. Содержит реализацию k-means.
  19. MLPACK
  20. Масштабируемая С++ библиотека для работы с алгоритмами машинного обучения, содержит реализацию k-means.
  21. SciPy
  22. Библиотека Python, содержит множество реализаций k-means.
  23. scikit-learn
  24. Библиотека Python, содержит множество реализаций k-means.
  25. R
  26. Язык программирования для статистической обработки данных и работы с графикой, а также свободная программная среда вычислений с открытым исходным кодом в рамках проекта GNU, содержит три реализации k-means.
  27. ELKI
  28. Java фреймворк, содержащий реализацию k-means, а также множество других алгоритмов кластеризации.

2.7.2 Проприетарное программное обеспечение

  1. Ayasdi
  2. Stata
  3. Mathematica
  4. MATLAB
  5. SAS
  6. RapidMiner
  7. SAP HANA

3 Литература

  1. "https://ru.wikipedia.org/wiki/Иерархическая_кластеризация"
  2. "https://ru.wikipedia.org/wiki/Кластерный_анализ"
  3. Steinhaus, Hugo. "Sur la division des corp materiels en parties." Bull. Acad. Polon. Sci 1.804 (1956): 801.
  4. Lloyd, S. P. "Least square quantization in PCM. Bell Telephone Laboratories Paper. Published in journal much later: Lloyd, SP: Least squares quantization in PCM." IEEE Trans. Inform. Theor.(1957/1982).
  5. MacQueen, James. "Some methods for classification and analysis of multivariate observations." Proceedings of the fifth Berkeley symposium on mathematical statistics and probability. Vol. 1. No. 14. 1967.
  6. "https://ru.wikipedia.org/wiki/EM-алгоритм"
  7. Zhao, Weizhong, Huifang Ma, and Qing He. "Parallel k-means clustering based on mapreduce." IEEE International Conference on Cloud Computing. Springer Berlin Heidelberg, 2009.
  8. Von Luxburg, Ulrike. Clustering Stability. Now Publishers Inc, 2010.
  9. Ortega, Joaquín Pérez, Ma Del Rocío Boone Rojas, and María J. Somodevilla. "Research issues on K-means Algorithm: An Experimental Trial Using Matlab."
  10. Воеводин Вл., Жуматий С., Соболев С., Антонов А., Брызгалов П., Никитенко Д., Стефанов К., Воеводин Вад. Практика суперкомпьютера «Ломоносов» // Открытые системы, 2012, N 7, С. 36-39.
  11. "http://ark.intel.com/ru/products/47920/Intel-Xeon-Processor-X5670-12M-Cache-2_93-GHz-6_40-GTs-Intel-QPI"
  12. "Issa. "Performance characterization and analysis for Hadoop K-means iteration". Journal of Cloud Computing, 2016"
  13. "Raghavan R. A fast and scalable hardware architecture for K-means clustering for big data analysis : дис. – University of Colorado Colorado Springs. Kraemer Family Library, 2016."
  14. "Yang, Luobin, et al. "High performance data clustering: a comparative analysis of performance for GPU, RASC, MPI, and OpenMP implementations." The Journal of supercomputing 70.1 (2014): 284-300."
  15. "Li, You, et al. "Speeding up k-means algorithm by GPUs." Computer and Information Technology (CIT), 2010 IEEE 10th International Conference on. IEEE, 2010."
  16. "Kanan, Gebali, Ibrahim. "Fast and Area-Efficient Hardware Implementation of the K-means Clustering Algorithm". WSEAS Transactions on circuits and systems. Vol. 15. 2016"