Участник:Astandrik/Алгоритм кластеризации, основанный на минимальном остовном дереве

Алгоритм кластеризации, основанный на минимальном остовном дереве
Последовательный алгоритм
Последовательная сложность	[math] O(N^2 log(N)) [/math]
Объём входных данных	[math] N [/math]
Объём выходных данных	[math] N [/math]
Параллельный алгоритм
Высота ярусно-параллельной формы	[math] O(log(N)) [/math]
Ширина ярусно-параллельной формы	[math] O(N) [/math]

Основные авторы описания: Стандрик Антон Сергеевич

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

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

Данный алгоритм предназначен для решения задач кластеризации. Говоря в общем, кластеризация — многомерная статистическая процедура, выполняющая сбор данных, содержащих информацию о выборке объектов, и затем упорядочивающая объекты в сравнительно однородные группы .

Впервые идея использования минимального остовного дерева для кластеризации данных была преждожена C. T. Zah в 1971 году. Этот метод широко используется исследователями при анализе изображений , а также при анализе биологических данных .

В общем описании работы алгоритма можно выделить следующие шаги:

1. Разбиение входных векторов на пары и вычисление расстояний в них. Для этого небходим выбор метрики пространства, в котором осуществляется кластерный анализ. В качестве примеров метрик можно привести:
   • Метрика Евклида: [math]d_{ij} = \sqrt{\sum^n_{k=1}(x_{ik} - x_{jk})^2}[/math]
   • Метрика city-block (манхэттенская метрика): [math]d_{ij} = \sum^n_{k=1}|x_{ik} - x_{jk}| [/math]
   • Метрика Брея-Картиса: [math]d_{ij} = \sum^n_{k=1}\frac{|x_{ik} - x_{jk}|}{|x_{ik} + x_{jk}| }[/math]
2. Построение минимального остовного дерева по заданному множеству. Существует три основных алгоритма построения минимального остовного дерева:
   • Алгоритм Прима
   • Алгоритм Крускала
   • Алгоритм Борувки
3. Разделение полученного дерева на несвязное множество поддеревьев в соответствии с некоторой стратегией. Можно выделить три популярные стратегии кластеризации:
   • SEMST (standard Euclidean minimum spanning tree clustering algorithm) - предполагается, что нам заранее известно число кластеров. В таком случае, если число кластеров равно [math]k[/math], то из полученного в первом пункте остовного дерева удаляются [math]k-1[/math] связей с наибольшим весом.
   • CEMST(The maximum cost minimum spanning tree) - Работает так же как и SEMST, однако вместо веса вершины используется её стоимость, которая вычисляется следующим образом: Пусть ребро [math]e = (x_{i}, x_{j})[/math], тогда [math]cost(e) = w(e) * degree(x_{i}) * degree(x_{j})[/math]. Далее, зная число кластеров [math]k[/math] и стоимости рёбер, удаляются [math]k-1[/math] рёбер с максимальной стоимостью.
   • ZEMST (The Zahn’s minimum spanning tree) - от двух предыдущих эта стратегия отличается тем, что здесь число кластеров не известно заранее, а вычисляется в процессе кластеризации. Его суть в том, что удаляются те рёбра, веса которых значительно превосходят средний вес ребра в соседних поддеревьях.

1.1.1 Достоинства

Основным достоинством алгоритма является его способность выделять кластеры произвольной формы. На следующем примере можно видеть, как алгоритм даже в случае SEMST алгоритма корректно разбивает вложенные окружности на три кластера:

Слева направо: Рисунок 1: Начальные данные. Рисунок 2: построенное минимальное остовное дерево. Рисунок 3: разбиение на кластеры путём удаления самых длинных связей.

1.1.2 Недостатки

Недостатком является трудность определения четких критериев, при которых осуществляется разбиение данных на разные кластеры.

Кроме того, существуют наборы данных, на которых популярные алгоритмы (SEMST, CEMST, ZEMST) не дают правильного результата, и для которых требуется применение более сложных алгоритмов (например HEMST с временем работы [math]O(n^2 * log(n))[/math])

В качестве иллюстрации можно привести два похожих визуально, но различных с точки зрения разбиения на кластеры с помощью описываемого алгоритма, примера. (В приведенных примерах первый рисунок обозначает начальный набор данных. Второй - построенное на основе данных минимальное остовное дерево. Третий рисунок - разбиение на кластеры.):

Правильное разбиение на два кластера.

Неверное разбиение на один кластер.

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

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

Исходные данные: набор векторов размерности [math]n[/math] : [math](a_1, a_2, ..., a_n)[/math], где [math] a_j \in \mathbb {R} \quad \forall j \in [1,n] [/math].

Вычисляемые данные: Сгруппированные по кластерам данные.

1.2.1 Основные обозначения

• [math]E[/math] - множество рёбер.
• [math]V[/math] - множество вершин.
• [math]G = (V,E) [/math] - граф, заданный множеством вершин [math]V[/math] и множеством рёбер [math]E[/math]
• [math]MST(G)[/math] (Minimum Spanning Tree) - минимальное остовное дерево графа [math][/math]. Остовное дерево — ациклический связный подграф данного связного неориентированного графа, в который входят все его вершины. Минимальное остовное дерево в связанном взвешенном неориентированном графе — это остовное дерево этого графа, имеющее минимальный возможный вес, где под весом дерева понимается сумма весов входящих в него рёбер.
• [math]MSF(G)[/math] (Minimum Spanning Forest) - набор минимальных остовных деревьев для всех компонент связности графа [math] G [/math]

1.2.2 Алгоритмы кластеризации

SEMST(standard Euclidean minimum spanning tree) - простейший алгоритм, основанный на заранее известном числе кластеров [math] k [/math]. Путём удаления [math] k - 1 [/math] связей с наибольшим весом предполагается разбиение остовного дерева на нужное количество кластеров. Иными словами минимальное остовное дерево разбивается на [math] k [/math] поддеревьев: [math] G = \{G_1,G_2,...,G_k\} [/math] путём минимизации суммы всех весов: [math] W = \sum ^k_{l=1} \sum _{x_i,x_j \in V_l} w_{x_i,x_j} -\gt min [/math], где [math] V_l [/math] - множество вершин в поддереве [math] G_l [/math]. В конце такой процедуры каждое множество [math] V_l [/math] становится кластером [math] C_l [/math]. Поскольку (в случае использования алгоритма Крускала для построения остовного дерева) рёбра перед началом кластеризации уже отсортированы, для удаления наибольших [math] k [/math] рёбер требуется [math] O(k) [/math] операций, и для распределения получившихся вершин по кластерам требуется время [math] O(n) [/math].

CEMST TODO ZEMST TODO

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

Процесс решения задачи кластеризации при помощи минимального остовного дерева может быть разделен на три независимых этапа:

1. Приведение задачи наборе векторов данных к задаче о поиске минимального остовного дерева в графе (то есть нахождение расстояний между всеми парами данных векторов с использованием выбранной метрики, и запись этого расстояния как веса ребра, соединяющего эти вершины). Последовательная сложность: [math] O(n^2) [/math]
2. Поиск минимального остовного дерева. Последовательная сложность: [math] O(E*log(E)) [/math]
3. Разбиение остовного дерева на кластеры. Последовательная сложность зависит от реализации. (от [math] O(1) [/math] в случае SEMST до [math] O(n^2) [/math] в случае ZEMST).

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

int N;        // Размер входных данных
int number_of_clusters;        // Число кластеров
int size = N* (N-1) / 2;        // Количество связей в полном графе
Node nodes[N];        // Вершины(точки), заданные начальными условиями
Edge edges[size];        // Все расстояния между вершинами 
vector<Edge> MST;        // Рёбра минимального остовного дерева
read_data(&nodes, file_name);        // Считывание начальных данных

//----------------Часть 1--------------//
for(int i = 0; i < N-1; i++) {
  for(int j = i+1; j < N; j++) {
	double distance = getDistance(node[i], node[j]);        // Вычисление расстояния между вершинами
	edges.[i*N + j] = (Edge(i,j,distance));        // Запись значения расстояния в массив
  }
}


//----------------Часть 2--------------//
sort(edges, size);        // Сортировка рёбер графа в порядке возрастания
for(int i = 0; i < size; i++) {        // Обходим список рёбер по возрастанию 
	if(find_set(edges[i].first) != find_set(edges[i].second)) {        // Если вершины принадлежат разным подграфам - связываем их
      Union(edges[i].first, edges[i].second)        // Если вершины принадлежат разным подграфам - связываем их
	  MST.push_back(edges[i])        // Запись использованного ребра в минимальное остовное дерево
	}
}


//----------------Часть 3--------------//
for(int i = MST.size() - 1; i > MST.size() - 1 - number_of_clusters; i--) {
	remove_edge(MST, i);        // Удаление самых длинных рёбер из MST - формирование кластеров
}
write_data(MST);        // Запись выходных данных

Как описано в пункте 1.3, весь описываемый алгоритм (начиная от момента получения данных и заканчивая моментом, когда данные отданы на вывод), может быть условно разделён на три части. В то время как первая часть может быть достаточно эффективно решена "в лоб", вторая и третьи части, как указано в пункте 1.1, могут иметь альтернативы. В данном фрагменте кода подробно описано использование Алгоритма Крускала для построения MST, и алгоритма SEMST для разбиения на кластеры.

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

Как и в остальных пунктах, здесь сложность алгоритма состоит их суммы сложностей трёх его составных частей. Если оценивать сложность асимптотически в зависимости от входных данных, то:

[math] T(n) = O(n^2) + O(\frac{n (n-1)}{2} log(\frac{n (n-1)}{2})) + O(1) [/math]

На первый взгляд второе слагаемое может показаться слишком большим. Однако, верны и следующие равенства:

[math] O(\frac{n (n-1)}{2} log(\frac{n (n-1)}{2}) = O(n^2 log(n^2)) = O(2 n^2 log(n)) = O(n^2 log(n)) [/math]

Отсюда видно, что для входных данных этого алгоритма асимптотические сложности алгоритмов Крускала, Прима и Борувки совпадают.

Таким образом, последовательная сложность алгоритма [math] O(n^2 log(n)) [/math]

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

Макро-описание параллельного алгоритма кластеризации

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

Данная задача хорошо распараллеливается благодаря тому, что все "тяжелые" вычислительные задачи, начиная с момента считывания входных данных, можно распределить по процессам практически равномерно. Это достигается за счет следующей реализации:

1. Мастер-процесс считывает входные данные и отправляет их остальным процессам с помощью MPI_Broadcast.
2. Все пары точек разбиваются так, что каждому процессу достается одинаковое количество пар.
3. Процесс-рабочий считает расстояния в своих парах
4. Процесс-рабочий сортирует пары по расстояниям
5. Процесс-рабочий строит MST из своих пар точек.

......

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

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

В качестве входных данных программе дается множество векторов [math] x_j = (a_{1j}, a_{2j}, ..., a_{nj})[/math] размерности [math]n[/math]. Для удобства изучения свойств алгоритма и отображения результатов можно брать массив векторов размерности 2. В этом случае наглядность будет обеспечена естественной возможностью отобразить вектора в плоскости.

Кроме того может оказаться полезным особым образом подготовить данные перед анализом (удаление "шумов", выбросов и масштабирование данных по тем или иным координатам).

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

Выходными данными программы будет множество [math] M = ((x_1,x_2,..,x_k),(x_{k+1},...),...) [/math], состоящее из сгруппированных по кластерам векторов из входных данных.