Участник:Каледа Александр/Строгий алгоритм С средних (Hard C-Means)

Строгий алгоритм С средних (Hard C-Means)
Последовательный алгоритм
Последовательная сложность	[math]O(n)[/math]
Объём входных данных	[math]nm+1[/math]
Объём выходных данных	[math]n+Cm[/math]
Параллельный алгоритм
Высота ярусно-параллельной формы	[math]O(1)[/math]
Ширина ярусно-параллельной формы	[math]O(n)[/math]

Авторы описания: Каледа Александр, Шапулин Андрей.

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

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

Строгий алгоритм C средних (Hard C-Means, HCM; Jang, Sun and Mizutani, 1997) один из наиболее широко используемых алгоритмов кластеризации данных. Алгоритм разбивает большие наборы векторов в многомерном пространстве по заранее заданному количеству кластеров (т.е. является управляемым).

На вход алгоритму подается набор векторов и количество кластеров. Дополнительно может задаваться пороговая величина объектной функции, при достижении которой алгоритм останавливается. После окончания работы алгоритма каждый вектор будет отнесет к определенному кластерному центру. Алгоритм, являясь итерационным, зависит от начальных значений кластерных центров, а следовательно может не сходиться к локальному минимуму.

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

Вычисления, выполняющиеся на очередной итерации алгоритма:

1. Расчёт рядовой матрицы (сопоставление вектора соответствующему кластерному центру).
2. Расчёт объектной функции.
3. Пересчёт кластерных центров.

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

Исходные данные: матрица [math]U\in \mathbb{R}^{n\times m}[/math], представляющая собой набор векторов [math]\{u_i\}_{i=1}^{n}, u_i \in \mathbb{R}^m[/math]; количество кластеров [math]C \in \mathbb{N}[/math].

Вычисляемые данные: набор векторов - кластерных центров [math]\{c_i\}_{i=1}^{C}, c_i \in \mathbb{R}^m[/math]; вектор [math]y \in \mathbb{N}^m[/math], содержащий номер найденного кластера для каждого объекта [math]u_i[/math].

Алгоритм:

Шаг 1. Инициализация кластерных центров [math]c_i, i = 1, 2, \dots, C.[/math] Это можно сделать выбрав случайным образом C векторов из входного набора [math]\{u_i\}_{i=1}^{n}[/math]

Шаг 2. Вычисление рядовой матрицы M.

Матрица [math]M[/math] состоит из элементов [math]M = \|m_{ik}\|[/math].

[math] m_{ik} = \begin{cases} 1, & \|u_k-c_i\|^2 \leqslant \|u_k-c_j\|^2 \;\; \text{for all } i \neq j; \;\; (i,j=1,2,\ldots,C; \; k = 1,2,\ldots,n)\\ 0, & \text{otherwise}. \end{cases} [/math]

Шаг 3. Расчет объектной функции:

[math] J = \sum\limits_{i=1}^C J_i = \sum\limits_{i=1}^C \left( \sum\limits_{k:\, u_k \in C_i} \|u_k -c_i \|^2 \right)[/math], где [math]C_i = \{u_k\}, m_{ik} = 1[/math]. То есть [math]C_i[/math] содержит все объекты, принадлежащие кластеру [math]i[/math].

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

Шаг 4. Пересчет кластерных центров.

Пересчет кластерных центров выполняется в соответствии со следующим уравнением:

[math] c_i = \frac{1}{|C_i|} \sum\limits_{k:\, u_k \in C_i} u_k, [/math]

где [math]|C_i|[/math] -- количество элементов в [math]i[/math]-ом кластере.

Шаг 5. Переход на шаг 2.

Описание взято из .

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

Вычислительное ядро алгоритма Hard C-Means можно составить из множественных (всего их [math]nC[/math]) вычислений квадратов нормы разности векторов [math]u_k[/math] и кластерных центров [math]c_i[/math]: [math]\|u_k-c_i\|^2[/math].

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

Как записано и в описании ядра алгоритма, основную часть метода составляют множественные ([math]Cn[/math]) вычисления квадратов норм:

[math]\|u_k-c_i\|^2[/math].

Также присутствуют данные вычислительные операции:

• Суммирование квадратов норм векторов:

[math] J = \sum\limits_{i=1}^C J_i = \sum\limits_{i=1}^C \left( \sum\limits_{k:\, u_k \in C_i} \|u_k -c_i \|^2 \right) [/math]

• Поэлементное суммирование векторов:

[math] c_i = \frac{1}{|C_i|} \sum\limits_{k:\, u_k \in C_i} u_k [/math]

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

Итерационный блок алгоритма можно описать на языке Python следующим образом:

for i in range(n):
    idx_min = 0
    dist_min  = ((u[i] - c[0]) ** 2).sum()

    for j in range(1, C):
        dist  = ((u[i] - c[j]) ** 2).sum()
        if dist < dist_min:
            idx_min = j
            dist_min = dist

    M[i, idx_min] = 1
    c_next[idx_min] += u[i]
    c_next_cnt[idx_min] += 1
    J += dist_min

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

Предполагаем, что предельная точность алгоритма достигаема менее, чем за [math]n[/math] итераций. Иначе использование алгоритма можно считать нерациональным. Также предположим, что размерность [math]m[/math] и количество кластерных центров [math]C[/math] существенно меньше [math]n[/math].

Операция вычисления квадрата нормы разности векторов требует выполнения [math]m-1[/math] арифметических операций сложения, [math]m[/math] операций вычитания и [math]m[/math] операций умножения.

Каждая итерация алгоритма требует выполнения:

• [math]nC[/math] операций вычисления квадрата нормы разности векторов;
• [math]nC[/math] операций сравнения квадрата нормы;
• [math]3n[/math] операций сложения;

Таким образом, на каждой итерации имеем:

• [math]mnC[/math] операций умножения;
• [math](2mC + 3)n[/math] арифметических операций сложения и вычитания;

В итоге имеем сложность последовательной реализации алгоритма [math]O(n)[/math].

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

Основной ресурс параллелизма алгоритма составляет итерационный блок. Эта последовательность итераций может быть выполнена для каждого входного вектора [math]u_i[/math] независимо. Однако, в конце выполнения каждого итерационного блока имеют место последовательные вычисления - сложение значений объектной функции, полученных независимо на каждом процессоре, а также перерасчет кластерных центров.

Высота ЯПФ оценивается как [math]O(1)[/math] (при условии, что число итераций существенно меньше объёма выборки [math]n[/math]). Ширина ЯПФ равна [math]n[/math].

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

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

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

• Последовательная реализация алгоритма на MatLab
• Последовательная реализация алгоритма на С++