Участник:Anton goy/Самоорганизующиеся карты Кохонена: различия между версиями

Материал из Алговики
Перейти к навигации Перейти к поиску
 
(не показаны 94 промежуточные версии 4 участников)
Строка 1: Строка 1:
 +
{{Assignment|ASA}}
  
 
+
Автор: [[U:Anton goy|Гой Антон]], 617 группа.
Автор: Гой Антон, 617 группа.
 
  
 
= Свойства и структура алгоритмов =
 
= Свойства и структура алгоритмов =
Строка 33: Строка 33:
 
# Задать начальные приближение весов <math>\mathbf{w}_{j}</math>, <math>\forall j=1, \dots, L</math>.
 
# Задать начальные приближение весов <math>\mathbf{w}_{j}</math>, <math>\forall j=1, \dots, L</math>.
 
# Выбрать случайно номер объекта <math>i \sim Unif[1, N]</math>.
 
# Выбрать случайно номер объекта <math>i \sim Unif[1, N]</math>.
# Найти расстояние между объектом <math>\mathbf{x}_i</math> и всеми векторами весов <math>\mathbf{w}_{j}</math>: <math>\rho(\mathbf{x}_i, \mathbf{w}_{j}) = \left\| \mathbf{x}_i - \mathbf{w}_{j} \right\|^2, \quad \forall j=1, \dots, L</math>
+
# Найти расстояние между объектом <math>\mathbf{x}_i</math> и всеми векторами весов <math>\mathbf{w}_{j}</math>: <math display=block>\rho(\mathbf{x}_i, \mathbf{w}_{j}) = \left\| \mathbf{x}_i - \mathbf{w}_{j} \right\|^2, \quad \forall j=1, \dots, L</math>
# Найти нейрон-победитель <math>n_b</math>, <math> b = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j})</math>, наиболее близкий к текущему объекту <math>\mathbf{x}_i</math>.
+
# Найти нейрон-победитель <math>n_b</math> наиболее близкий к текущему объекту <math>\mathbf{x}_i</math>, где <math display=block> b = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j}).</math>
# Для всех неронов <math>n_j</math> пересчитать их веса по следующей формуле: <math>\mathbf{w}_{j} \leftarrow \mathbf{w}_{j} + \alpha(t)h_{b, j}(t)(\mathbf{x}_i - \mathbf{w}_{j})</math>, где <math>\alpha(t)</math> - темп обучения (learning rate), монотонно убывающая функция от номера итерации <math>t</math>, <math>h_{b, j}(t)</math> - функция определяющая "меру соседства" нейронов <math>n_b</math> и <math>n_j</math>.
+
# Для всех неронов <math>n_j</math> пересчитать их веса по следующей формуле: <math display=block>\mathbf{w}_{j} \leftarrow \mathbf{w}_{j} + \alpha(t)h_{b, j}(t)(\mathbf{x}_i - \mathbf{w}_{j}),</math> где <math>\alpha(t)</math> - темп обучения (learning rate), монотонно убывающая функция от номера итерации <math>t</math>, <math>h_{b, j}(t)</math> - функция определяющая "меру соседства" нейронов <math>n_b</math> и <math>n_j</math>.
 
# Проверить критерий останова. При необходмости перейти на шаг 2.
 
# Проверить критерий останова. При необходмости перейти на шаг 2.
  
Строка 48: Строка 48:
 
** <math>\alpha(t) = \lambda^{\frac{t}{T}}, \lambda \in (0, 1)</math>
 
** <math>\alpha(t) = \lambda^{\frac{t}{T}}, \lambda \in (0, 1)</math>
 
* ''' Задание <math>h_{b,j}(t)</math> '''
 
* ''' Задание <math>h_{b,j}(t)</math> '''
** Наиболее популярный способ: <math>h_{b,j}(t) = \exp \left\{ - \frac{\| \mathbf{r}_b -\mathbf{r}_j \|^2}{2\sigma(t)^2} \right\}</math>, где <math>\sigma(t)</math> - убывающая функция.
+
** Наиболее популярный способ: <math display=block>h_{b,j}(t) = \exp \left\{ - \frac{\| \mathbf{r}_b -\mathbf{r}_j \|^2}{2\sigma(t)^2} \right\},</math> где <math>\sigma(t)</math> - убывающая функция, которую можно интерпретировать как величину пропорциональную радиусу соседства. Концептуально  принцип действия такой меры соседства можно представить следующим образом: пусть в точке <math>\mathbf{r}_b</math> задана окружность с радиусом пропорциональным <math>\sigma(t)</math>, тогда за пределом  этой окружности степень соседства становится пренебрежимо мала. Пример: <math display=block>\sigma(t) = \frac{\sigma_0}{1 + \frac{t}{T}},</math> где <math>\sigma_0 = \max \left\{ \frac{X}{2}, \frac{Y}{2} \right\} </math>, <math>T</math> - максимальное количество итераций. Для достижения лучших результатов на практике процедуру обучения разбивают на два этапа: первый этап - итерации выполняются с очень большим радиусом, что позволяет сойтись в хоррошую область быстро, а второй этап - радиус выбирается маленьким и постепенно уменьшается, что дает возможность весам более точно уловить зависимости из обучающей выборки.
 +
** Существуют более разреженные меры сходства, например, присваивающие ненулевые веса, только внутри некоторого радиуса. Такие меры сходства могут значительно ускорить вычисления, так как не требуется вычислять меру сходства для всех нейронов.
  
 
Кроме онлайн алгоритма обучения часто на практике применяют более эффективную с точки зрения вычислений batch-версию. Пусть выбрано <math>B \in \mathbb{N}</math>. Тогда batch-версия алгоритма описывается следующим образом:
 
Кроме онлайн алгоритма обучения часто на практике применяют более эффективную с точки зрения вычислений batch-версию. Пусть выбрано <math>B \in \mathbb{N}</math>. Тогда batch-версия алгоритма описывается следующим образом:
Строка 55: Строка 56:
  
 
# Задать начальные приближение весов <math>\mathbf{w}_{j}</math>, <math>\forall j=1, \dots, L</math>.
 
# Задать начальные приближение весов <math>\mathbf{w}_{j}</math>, <math>\forall j=1, \dots, L</math>.
# Случайно выбрать множество <math>\mathcal{B} \subset \{1, \dots, N \}</math> размера <math>B</math>.
+
# Случайно выбрать множество <math>\mathcal{B} \subseteq \{1, \dots, N \}</math> размера <math>B</math>.
# Для каждого объекта <math>\mathbf{x}_i</math>, <math>i \in \mathcal{B}</math>, найти нейрон-победитель <math>n_{b(i)}</math>, где <math> b(i) = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j})</math>
+
# Для каждого объекта <math>\mathbf{x}_i</math>, <math>i \in \mathcal{B}</math>, найти нейрон-победитель <math>n_{b(i)}</math>, где <math display=block> b(i) = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j})</math>
# Пересчитать веса для всех нейронов <math>\mathbf{w}_j \leftarrow \frac{\sum_{i \in \mathcal{B}}h_{b(i),i}\mathbf{x}_i}{\sum_{i \in \mathcal{B}}h_{b(i),i}}</math>, <math>\forall j=1, \dots, L</math>.
+
# Пересчитать веса для всех нейронов <math display=block>\mathbf{w}_j \leftarrow \frac{\sum_{i \in \mathcal{B}}h_{b(i),j}\mathbf{x}_i}{\sum_{i \in \mathcal{B}}h_{b(i),j}}, \forall j=1, \dots, L</math>.
 
# Проверить критерий останова. При необходимости перейти на шаг 2.
 
# Проверить критерий останова. При необходимости перейти на шаг 2.
  
Интерпретация шага 2 состоит в следующем: веса нейронов преставляют собой выпуклую комбинацию объектов с номерами из <math>\mathcal{B}</math>.
+
Интерпретация шага 2 состоит в следующем: веса нейронов преставляют собой выпуклую комбинацию объектов с номерами из <math>\mathcal{B}</math>. Часто на практике полагают <math>\mathcal{B} = \{1, \dots, N\}</math>, то есть на каждой итерации мы совершаем полный проход по всей выборке объектов.
 +
 
 +
Выходом алгоритма можно считать веса нейронов, которые затем можно использовать для различных визуализаций и определния кластерной структуры обучающей выборки.
  
 
== Вычислительное ядро алгоритма ==
 
== Вычислительное ядро алгоритма ==
Строка 72: Строка 75:
 
Сложность вычисления <math>\mathbf{D}</math> составляет <math>O\left(L D B \right)</math>, учитывая что <math>\left\|\mathbf{w}_j\right\|^2</math> можно вычислить за один проход по всем нейронам.  
 
Сложность вычисления <math>\mathbf{D}</math> составляет <math>O\left(L D B \right)</math>, учитывая что <math>\left\|\mathbf{w}_j\right\|^2</math> можно вычислить за один проход по всем нейронам.  
  
Далее поиск минимум в каждой строке матрицы <math>\mathbf{D}</math> занимает <math>O(LB)</math>
+
Еще одной вычислительно интенсивной частью является вычисление <math>\mathbf{w}_j</math>. Здесь мы можем вычислить все необходимые <math>h_{b,j}(t)</math>, а затем, используя эти значения, найти новые приблидения для <math>\mathbf{w}_j</math>.
  
 
== Макроструктура алгоритма ==
 
== Макроструктура алгоритма ==
Строка 91: Строка 94:
  
 
Функция <math>ComputeWeightedSumsAcrossObjects</math> вычисляет новые значения весов для каждого нейрона, суммируя строки матрицы <math>\mathbf{X}_{\mathcal{B}}</math> с весами пропорциональными их мере соседства с соответствующим нейроном-победителем. Неявно данная функция также использует вектор положений <math>\mathbf{r}_j</math> нейронов на двумерной плоскости, но так как эти положения остаются неизменными в продолжении всего обучения, мы не будем указывать их явно в аргументах.
 
Функция <math>ComputeWeightedSumsAcrossObjects</math> вычисляет новые значения весов для каждого нейрона, суммируя строки матрицы <math>\mathbf{X}_{\mathcal{B}}</math> с весами пропорциональными их мере соседства с соответствующим нейроном-победителем. Неявно данная функция также использует вектор положений <math>\mathbf{r}_j</math> нейронов на двумерной плоскости, но так как эти положения остаются неизменными в продолжении всего обучения, мы не будем указывать их явно в аргументах.
 +
 +
Если рассмотреть вычисления выше описанных функций, то можно также выделить следующий набор макроопераций более низкого порядка, которые используются выше определенными функциями:
 +
 +
# <math>\mathbf{x}^{T} \mathbf{y}</math> - скалярное произведение двух вектров;
 +
# <math>\alpha \mathbf{x}</math> - умножение вектора на скаляр;
 +
# <math>\mathbf{x} + \mathbf{y}</math> - поэлементное сложение двух векторов;
 +
# <math>\min_i x_i</math> - поиск минимального элемента в векторе <math>\mathbf{x}</math>;
 +
# Вычисление функции <math>h_{b,j}(t)</math>.
  
 
== Схема реализации последовательного алгоритма ==
 
== Схема реализации последовательного алгоритма ==
Строка 96: Строка 107:
 
В этой секции, используя макрооперции из предыдущего раздела, дадим очень краткий псевдокод batch-версии алгоритма, который легко позволяет понимать основные строительные блоки данного алгоритма:
 
В этой секции, используя макрооперции из предыдущего раздела, дадим очень краткий псевдокод batch-версии алгоритма, который легко позволяет понимать основные строительные блоки данного алгоритма:
  
<math display="block">
+
[[File:Kohonen_network_hex_init.png|right|thumb|250px|Рисунок 3. Способ инициализации <math>\mathbf{r}_j</math> как узлов шестиугольных ячеек на плоскости]]
 +
 
 +
<math>
 
\begin{align}
 
\begin{align}
 
\mathbf{W} := &InitializeWeights() \\
 
\mathbf{W} := &InitializeWeights() \\
 
\mathbf{R} := &InitializeNeuronPositionsOnPlane() \\
 
\mathbf{R} := &InitializeNeuronPositionsOnPlane() \\
\mathbf{FOR} & \ t=1 \ \mathbf{TO} \ MAX\_ITERATIONS \ \mathbf{DO} \\
+
\mathbf{FOR} & \ t=1 \ \mathbf{TO} \ T \ \mathbf{DO} \\
 
         & \mathcal{B} \leftarrow RandomSample\left( \{1, \dots, N\} \right) \\
 
         & \mathcal{B} \leftarrow RandomSample\left( \{1, \dots, N\} \right) \\
 
         & \mathbf{D} \leftarrow ComputePairwiseDistances(\mathbf{X}_{\mathcal{B}}, \mathbf{W}) \\
 
         & \mathbf{D} \leftarrow ComputePairwiseDistances(\mathbf{X}_{\mathcal{B}}, \mathbf{W}) \\
Строка 108: Строка 121:
 
\end{align}
 
\end{align}
 
</math>
 
</math>
 +
 +
 +
 +
Здесь мы использовали <math>\mathbf{R} \in \mathbb{R}^{L \times 2}</math> как обозначение для матрицы, строками которой являются <math>\mathbf{r}_j</math>. Способ инициализации этих векторов для решетки с шестиугольными ячейками показан на Рис. 3. Кроме того, стоит отметить, что при нам никто не запрещает посчитать расстояния между положениями нейронов на двумерной сетки заранее, до основного цикла алгоритма. Это позволит на при подсчете меры соседства <math>h_{b,j}(t)</math> использовать уже вычисленные <math>\rho(\mathbf{r}_b, \mathbf{r}_j)</math>.
  
 
Как видно из алгоритма, основной структурой данных, которая нам потребуется, будет обычная матрица, а также вектора. Стоит отметить, что алгоритм в целом неплохо оптимизируем в терминах матричных и веторных операций. Например, вычисление матрицы <math>\mathbf{D}</math> на языке Python с использованием библиотеки оптимизированных матричных вычислений NumPy будет выглядить следующим образом:
 
Как видно из алгоритма, основной структурой данных, которая нам потребуется, будет обычная матрица, а также вектора. Стоит отметить, что алгоритм в целом неплохо оптимизируем в терминах матричных и веторных операций. Например, вычисление матрицы <math>\mathbf{D}</math> на языке Python с использованием библиотеки оптимизированных матричных вычислений NumPy будет выглядить следующим образом:
Строка 125: Строка 142:
 
В разделе [[#Вычислительное ядро алгоритма|1.3]] мы уже попытались дать некоторые оценки сложности алгоритма обучения самоорганизующейся карты Кохонена. Ниже мы рассмотрим последовательную сложность более подробно.
 
В разделе [[#Вычислительное ядро алгоритма|1.3]] мы уже попытались дать некоторые оценки сложности алгоритма обучения самоорганизующейся карты Кохонена. Ниже мы рассмотрим последовательную сложность более подробно.
  
Прежде всего необходимо понять в терминах каких операций нам следует оценить сложность алгоритма. В алгоритме преобладают арифметические операции над вещественными числами
+
Прежде всего необходимо понять в терминах каких операций нам следует оценить сложность алгоритма. В алгоритме преобладают арифметические операции и операции сравнения над вещественными числами, поэтому в их терминах и будем определять сложность.
 +
 
 +
Пусть, как и ранее, <math>N</math> - количество объектов в обучающей выборке, <math>L</math> - количество нейронов в слое Кохонена, <math>D</math> - размерность признакового пространства объектов (размерность векторов <math>\mathbf{x}_i</math>) и <math>B</math> - количество объектов в случайной подвыборке отбираемой на каждой итерации алгоритма.
 +
 
 +
Оценим сложность одной итерации алгоритма:
 +
 
 +
# Функция <math>ComputePairwiseDistances</math>.
 +
## <math>L \times D</math> умножений и <math>L \times (D - 1)</math> сложений для вычисления всех норм <math>\left\|\mathbf{w}_j\right\|^2</math>;
 +
## <math>B \times L \times D</math> умножений и <math>B \times  L \times (D - 1)</math> сложений для вычисления всех скалярных произведений <math>\mathbf{x}_i^{T}\mathbf{w}_j</math>;
 +
## <math>B \times L</math> сложений для вычисления <math>\hat{d}_{ij}</math>;
 +
## '''Итого:''' <math>L \times D + L \times (D - 1) + B \times L \times D + B \times  L \times (D - 1) + B \times L = 2 \times L \times D \times (B + 1) - L</math> арифметических операций над вещественными числами. Тогда асимптотическая сложность равна <math>O(B \times L \times D)</math>;
 +
# Функция <math>FindClosestNeurons</math>.
 +
## В этой функции используются только операции сравнения между вещественными числами и их количество составляет <math>B \times (L - 1)</math> или асимптотически <math>O(B \times L)</math>;
 +
# Функция <math>ComputeWeightedSumsAcrossObjects</math>.
 +
## Будем считать что для вычисления <math>h_{b, j}(t)</math> требуется некоторое константное число операций <math>Const</math>. Тогда для вычисления всех таких  функций требуется <math>B \times L \times Const</math> операций.
 +
## Вычисление <math>\mathbf{w}_j</math> делится на вычисление числителя и выисление знаменателя. В первом случае нам требуется выполнить <math>B \times D</math> умножений и <math>(B - 1) \times D</math> сложений. Во втором случае мы выполняем <math>B</math> сложений.
 +
## '''Итого:''' <math>(2 \times B \times D + B - D + 1) \times L + B \times L \times Const</math> операций. Асимптотическая сложность: <math>O(B \times L \times D)</math>
 +
 
 +
В итоге асимптотическая сложность одной итерации будет составлять <math>O(B \times L \times D)</math>.
 +
 
 +
Отметим, что не смотря на все количество операций, используя векторные и матричные операции процессора можно добится значительной эффективности.
 +
 
 +
== Информационный граф ==
 +
 
 +
Для того чтобы сохранить простоту графа, мы выбрали максимально крупные операции, которые в тоже время хорошо показывают структуру алгоритма. Кроме того, в графе мы оперируем не на уровне отдельных компонент вектора, а на уровне целых векторов-объектов и векторов-весов и векторных операций над ними.
 +
 
 +
Нами были рассмотрены следующие операции:
 +
 
 +
# <tt>sq_norm</tt> - вычисление квадрата нормы <math>\| \mathbf{w}_j \|^2 = \mathbf{w}_j^{T} \mathbf{w}_j</math>
 +
# <tt>distance</tt> - вычисления квазирасстояния между объектом и вектором весов <math>\| \mathbf{w}_j \|^2 - 2 \mathbf{x}_i^{T} \mathbf{w}_j</math>
 +
# <tt>min</tt> - поиск ближайшего нейрона на основе расстояний, которые получены на входе операции
 +
# <tt>h <math>\ast</math> '''x''' </tt> - вычисление веса соседства и умножения этого веса на вектор <math>\mathbf{x}</math> поэлементно
 +
# <tt>update</tt> - суммирование взвешенных векторов <math>\mathbf{x}</math> и обновление весов нейронов
 +
 
 +
На Рис. 4 и Рис. 5 представлен информационный граф алгоритма для трех объектов <math>\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3</math>  и двух нейронов с весами <math>\mathbf{w}_1, \mathbf{w}_2</math> в разных проекциях на плоскость.
 +
 
 +
[[File:Information_graph1.png|thumb|center|700px|Рисунок 4. Информационный граф batch-алгоритма обучения самоорганизующейся карты Кохонена с тремя входными объектами <math>\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3</math> и с двумя нейронами с весами <math>\mathbf{w}_1, \mathbf{w}_2</math>. Проекция 1.]]
 +
 
 +
[[File:Information_graph2.png|thumb|center|700px|Рисунок 5. Информационный граф batch-алгоритма обучения самоорганизующейся карты Кохонена с тремя входными объектами <math>\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3</math> и с двумя нейронами с весами <math>\mathbf{w}_1, \mathbf{w}_2</math>. Проекция 2.]]
 +
 
 +
== Ресурс параллелизма алгоритма ==
 +
 
 +
Наиболее очевидными способами использования параллельных ресурсов для ускорения алгоритма являются:
 +
* распределить объекты обучающей выборки по процессам (распределение данных);
 +
* распределить веса нейронов по процессам (распределение сети);
 +
 
 +
В первом случае в каждом узле вычислительной сети мы храним часть обучающей выборки и всю нейронную сеть целиком (все веса <math>\mathbf{w}_j</math>), а во втором в каждом узле хранится вся обучающая выборка и некоторая подсеть исходной нейронной сети.
 +
 
 +
Допустим неограниченное количество вычислительных узлов.
 +
 
 +
Рассмотрим первый случай. Имея в своем распоряжении неограниченное количество вычислительных узлов мы можем полностью распараллелить поиск минимальных расстояний, который теперь можно сделать за <math>O(L \times D)</math> арифметических операций с вещественными числами. Поиск ближайшего нейрона параллельно можно сделать за <math>O(L)</math> сравнений. Сложность  вычисления <math>\mathbf{w}_j</math>  - <math>O(B \times L \times D)</math>, так как нам необходима пересылка данных одному процессу, который вычислит новые значения для весов. Таким образом асимптотически мы не получили улучшения, но если рассмотреть точное количество операций при, то улучшение есть.
 +
 
 +
Рассмотрим второй случай. Расстояния мы по прежнему можем считать параллельно - <math>O(B \times D)</math>.  Но поиск минимального расстояния нужно выполнять последовательно, так как у одного  рассматриваемого процесса нет расстояний до всех нейронов, поэтому сложность -  <math>O(B \times L)</math>.  Сложность обновления нейронов - <math>O(B \times D)</math>. В итоге получаем <math>O(B \times ( L + D) )</math>.
 +
 
 +
== Входные и выходные данные алгоритма ==
 +
 
 +
Входными данными алгоритма являются:
 +
# множество <math>\mathcal{L} = \left\{ \mathbf{x}_i \right\}_{i=1}^N</math> - обучающая выборка. Наиболее просто можно представить обучающую выборку как матрицу <math>\mathbf{X} \in \mathbb{R}^{N \times D}</math>. Важно понимать, что на практике часто встречаются большие выборки состоящие из сотней тысяч объектов, каждый из которых может иметь десятки или сотни признаков (компонет), поэтомк размеры матрицы могут быть значительными;
 +
# количество нейронов <math>L</math>, кроме того явно конкретизированы размеры выходной карты <math>X, Y</math>;
 +
# положения нейронов на двумерной карте <math>\mathbf{r}_j</math>, которые мы можем упорядочить в матрицу <math>\mathbf{R}</math>. На самом деле задав размеры сетки и форму ячеек, мы неявно уже определили расстояния между положениями нейронов, поэтому, вообще говоря, описывать матрицу <math>\mathbf{R}</math> как входные данные излишне, но для ясности сделаем это;
 +
# количество итераций <math>T</math>;
 +
# начальный радиус <math>\sigma_0</math>.
 +
 
 +
Общая память, которая требуется нам для хранения входных данных, равна <math>O(N \times D)</math>.
 +
 
 +
Выходными данными алгоритма являются:
 +
 
 +
# обученные веса нейронов <math>\mathbf{w}_j</math>, которые удобно хранить в виде матрицы <math>\mathbf{W}</math>.
 +
 
 +
Общая память, которая требуется нам для хранения выходных данных, равна <math>O(L \times D)</math>.
 +
 
 +
== Свойства алгоритма ==
 +
 
 +
Вычислительная мощность алгоритма равна <math>O\left(\frac{BL}{N + L}\right)</math>, что при <math>B=N</math> равно <math>O\left(\frac{NL}{N + L}\right)</math>, что асимптотически при большом размере выборке равно <math>O(1)</math>.
 +
 
 +
Задача, которую решают самоогранизующиеся карты, является чрезвычайно сложной и некорректно поставленной, поэтому важной составляющей алгоритма является случайность. Она проявляется в случайном порядке выбора объектов (или случайных подмножеств объектов) для обучения. Но данное стохастическое поведение алгоритма не имеет существенного влияния на параллельную структуру алгоритма, кроме того если в batch-версии выбрать <math>B=N</math>, то случайность будет устранена.
 +
 
 +
Рассматрим свойства параллельной реализации. Основной проблемой реализации параллельной версии алгоритма является то, что нам необходимо выполнять массово две операции:
 +
* искать ближайший нейрон к текущему объекту;
 +
* обновлять веса каждого нейрона в конце итерации.
 +
 
 +
Так, если мы выполнили распределение данных, то мы имеем часть объектов, но в тоже время у нас есть локально вся сеть, следовательно найти ближайший нейрон для каждого объекта из известной нами подвыборки не требует взаимодействия с другими процессами. Но для того чтобы обновить значения весов нейронов нам необходима информация от других вычислительных узлов, в частности нам требуются <math>h_{b(i), j}</math> и <math>\mathbf{x}_i</math>, чтобы обновить веса нашей сети. Следовательно, нам необходим отдельный процесс-аггрегатор, который, получив все необходимые данные, обновил значения весов и переслал бы обновленную нейронную сеть всем остальным процессам. В целом, такой подход полностью укладывается в парадигму MapReduce.
 +
 
 +
Если же мы выбрали стратегию распределения сети, то теперь мы имеем только часть сети, но зато вся обучающая выборка. Теперь найти ближайший нейрон не так просто, так как локально мы можем проверить только те нейроны, которые хранятся в памяти процесса, поэтому остальные процессы должны переслать нам их минимальное расстояние. Далее зная номер ближайшего нейрона (при условии, что мы храним локально все положения нейронов <math>\mathbf{R}</math>), мы можем без труда найти выпуклую комбинацию объектов выборки и обновить веса нейронов.
 +
 
 +
Ясно, что каждый из подходов имеет свои преимущества: первый более масштабируем для выборок большого размера, когда мы физически не можем поместить все объекты в один узел, но в тоже время размеры самой нейронной сети незначительны; во втором случае мы можем строить сети состоящие из большого количества нейронов при минимальном взаимодействии процессов (процессу требуется только передать два числа - минимальное расстояние и номер нейрона, которому оно соответствует).
 +
 
 +
= Программная реализация алгоритма =
 +
 
 +
== Особенности реализации последовательного алгоритма ==
 +
 
 +
== Локальность данных и вычислений ==
 +
 
 +
== Возможные способы и особенности параллельной реализации алгоритма ==
 +
 
 +
== Масштабируемость алгоритма и его реализации ==
 +
В данном разделе мы рассмотрим паралелльную реализацию алгоритма, которая основана на разделении входной выборки по всем доступным процессам. Каждый процесс работает со своим кусоком данным и в конце каждой итерации все процессы синхронизируют состояние весов нейронов. В качестве допущения мы будем рассматривать только одну итерацию алгоритма (если итераций больше чем одна, то результаты просто будут отличаться на множитель равный количеству итераций). Кроме того, допустим <math>B = N</math>.
 +
 
 +
Все эксперименты были выполнены на суперкомпьютере IBM Blue Gene/P. Выборка была сгенерировна из стандартного нормального распределения для <math>N = D = 1000</math>. Набор изменяемых параметров запуска:
 +
 
 +
# <math>L</math> (число нейронов с решетке) = [ 1024,  4096,  9216, 16384, 25600, 36864, 50176, 65536]
 +
# <math>P</math> (число процессоров) = [  8,  16,  32,  48,  64, 128, 384, 256, 512]
 +
 
 +
Для оценки производительности (GFLOPS) была использована оценка числа операций приведенная ранее.
 +
 
 +
На Рис 6 показана зависимость производительности алгоритма от числа процессоров и количества нейронов в решетке. График подтверждает наше интуитивную догадку, что для больших размеров решетки алгоритм имеет большую производительность в терминах операций с плавающей точкой. Максимальная производитлеьность 24 GFLOPS достигается на 512 процессорах и для 50176 (224 * 224) нейронов. Минимальная же производительность (0.6 GFLOPS) наблюдается на 8 процессорах при любом числе нейронов.
 +
 
 +
[[File:Performance som1.png|1000px|thumb|center|Рисунок 6. Производительность алгоритма SOM]]
 +
 
 +
На Рис 7 показана аналогичная зависимость эффективности реализации. Здесь стоит отметить, что эффективность реализации очень низкая (в среднем 2% от пиковой производительности). Это долольно просто объяснить: в алгоритме требуется передовать больше количество данных (веса нейронов) для синхронизации процессов, поэтому значительную часть времени процессоры простаивают, ожидая передачи или отправки данных. Данная ситуация ухудшается, когда нейронов становится больше, что отлично видно на графике. Максимальная эффектвность (2.7 %) была достигнута для 8 процессов и остается неизменной при любом числе нейрнов, в данном случае более высокая эффективность достигается за счет того, что требуется осуществить синхронизацию только между 8 процессами, что намного дешевле, чем в случае 512 процессоров (минимально достигнутая эффективость - 1%).
 +
 
 +
[[File:Efficiency som1.png|1000px|thumb|center|Рисунок 7. Эффективность алгоритма SOM]]
 +
 
 +
В целом можно отметить, что алгоритм обладает свойством сильной масштабируемости: на Рис 6 можно наблюдать сублинейный рост производительности с ростом числа процессоров при фиксированном числе нейронов. Так же заметим, что алгоритм проявляет свойство слабой масштабируемости: так, например, при <math>W / p = 128</math> отмечаем следующий рост производительности: [0.61, 2.17, 7.92, 21.4] при соответствующих значенияз числа процессоров: [8, 32, 128, 512]. Другой пример: <math>W / p = 64</math>, производительности - [1.21, 4.17, 12.96] при [16, 64, 256].
 +
 
 +
== Динамические характеристики и эффективность реализации алгоритма ==
 +
 
 +
== Выводы для классов архитектур ==
 +
 
 +
== Существующие реализации алгоритма ==
 +
 
 +
Существующие параллельные реализации:
 +
 
 +
* [https://github.com/peterwittek/somoclu Somoclu]: Massively parallel self-organizing maps: accelerate training on multicore CPUs, GPUs, and clusters (GNU GPLv3).
 +
** поддержка OpenMP, MPI и CUDA;
 +
** интерфейсы для Python, R, Julia, и Matlab
 +
* [https://github.com/wimvanderbauwhede/Parallel-SOM Parallel-SOM]: Parallel Implementation of Self Organzing Map in OpenCL
 +
* [https://github.com/andreyto/mr-mpi-som MapReduce-MPI SOM]: Parallel implementation of Self-Organizing Map (SOM) algorithm with MapReduce-MPI (GNU GPLv3).
 +
 
 +
Существующие последовательные реализации:
 +
 
 +
* [https://github.com/ilarinieminen/SOM-Toolbox SOM-Toolbox]: A Matlab toolbox for Self-Organizing Maps (SOM) and others (GNU GPLv2).
 +
** предоставляет инструменты для удобной визуализации полученной карты.
 +
* [https://github.com/JustGlowing/minisom MiniSom]: A minimalistic implementation of the Self Organizing Maps (Creative Commons Attribution 3.0 Unported License).
 +
** реализация на языке Python с использованием библиотеки NumPy.
 +
* [https://github.com/sevamoo/SOMPY SOMPY]: A Python Library for Self Organizing Map (Apache License v2).
 +
* [http://knnl.sourceforge.net/ KNNL]: C++ Kohonen Neural Network Library (BSD License).
 +
* [http://jknnl.sourceforge.net/ JKNNL]: Java Kohonen Neural Network Library (BSD License).
 +
* [https://github.com/dashaub/kohonen4j kohonen4j]: Self-Organizing Maps in Java (GNU GPLv3).
 +
* [https://cran.r-project.org/web/packages/kohonen/ kohonen] Supervised and Unsupervised Self-Organising Maps in R (GNU GPLv3).
 +
 
 +
= Литература =
 +
 +
# Teuvo Kohonen (1997). Self-Organizing Maps.  Springer-Verlag New York, Inc., Secaucus, NJ, USA.
 +
# Peter Wittek, Shi Chao Gao, Ik Soo Lim, Li Zhao (2015). Somoclu: An Efficient Parallel Library for Self-Organizing Maps. arXiv:1305.1422.
 +
# Juha Vesanto (1999).  SOM-based data visualization methods. Intell. Data Anal. 3, 2 (March 1999), 111-126.
 +
# Timo D. Hämäläinen (2001). Parallel implementation of self-organizing maps. In Self-Organizing neural networks, Udo Seiffert and Lakhmi C. Jain (Eds.). Springer Studies In Fuzziness And Soft Computing Series, Vol. 78. Springer-Verlag New York, Inc., New York, NY, USA 245-278.
 +
# Seung-Jin Sul, Andrey Tovchigrechko (2011). Parallelizing BLAST and SOM Algorithms with MapReduce-MPI Library. IEEE International Symposium on Parallel and Distributed Processing.

Текущая версия на 10:54, 6 февраля 2017

Symbol wait.svgЭта работа прошла предварительную проверку
Дата последней правки страницы:
06.02.2017
Данная работа соответствует формальным критериям.
Проверено ASA.


Автор: Гой Антон, 617 группа.

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

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

Рисунок 1. Структура самоогранизующейся карты Кохонена с шестиугольной решеткой

Самоорганизующаяся карта Кохонена (англ. Self-Organizing Map или сокращено SOM) - это разновидность нейронных сетей, относящаяся к алгоритмам обучения без учителя. Основная цель - найти скрытые закономерности в данных по средством снижения размерности исходного пространства. Важным свойством карт Кохонена является то, что они строят отображение в пространство низкой размерности (обычно двумерное) таким образом, что топология исходного пространства сохраняется. Результат данного отображения - правильная решетка из обученных нейронов - называется "картой" исходного пространства. Алгоритм был разработан известным финским учёным, заслуженным академиком Финской Академии Наук Теуво Кохоненом в 1984(2) году. Карты Кохенана находят успешное применение в задачах кластеризации и визуализации, а также для снижения размерности и детектирования аномалий в данных.

Карты Кохонена и по своей архитектуре, и по методу обучения отличаются от обычных нейронных сетей прямого распространения. C точки зрения метода обучения, карты Кохонена не используют градиентные методы для минимизации ошибки (как это делается в сетях прямого распространения), поскольку являются алгоритмом обучения без учителя и никак не могут учитывать информацию и метках классов. Поэтому нейронная сеть обучается через соревнование между нейронами: на каждом шаге алгоритма для случайного объекта из обучающей выборки выбирается нейрон-победитель (best matching unit, BMU), который в определенном смысле наиболее похож на данный объект.

Рисунок 2б. Топология слоя Кохонена с прямоугольными ячейками
Рисунок 2а. Топология слоя Кохонена с шестиугольными ячейками

А архитектура карты Кохонена (Рис. 1) представляет два слоя: первый слой состоит из входных нейронов (их количество равно размерности исходного пространства), второй слой (его называют слоем Кохонена) представляет собой строго упорядоченную решетку, в узлах которой расположены нейроны-сумматоры, соедененные со всеми входами сети. Для визуализации топологии слоя Кохонена (Рис. 2) используют либо шестиугольные (Рис. 2а), либо прямоугольные (Рис. 2б) ячейки, в которых распологаются нейроны. Шестиугольные ячейки часто являются более предпочтительными, так как расстояние от центра выбранной ячейки до ее соседей одинаково.

Каждый нейрон [math]n_j[/math] слоя Кохонена описывается вектором весов [math]\mathbf{w}_j[/math], размерность которого совпадает с размерностью исходного пространства, и координатами нейрона на двумерной плоскости - [math]\mathbf{r}_j[/math].

Процесс обучения состоит в настройке векторов весов [math]\mathbf{w}_j[/math].

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

Пусть [math]\mathcal{L} = \left\{ \mathbf{x}_i \right\}_{i=1}^{N}[/math] - некоторое подмножество точек пространства [math]\mathbb{R}^D[/math], [math]\mathbf{x}_i = \left( x_{i}^{(1)}, \dots, x_{i}^{(D)} \right) \in \mathbb{R}^D[/math]. В машинном обучении множество [math]\mathcal{L}[/math] называют обучающей выборкой.

Пусть задано [math]L = X \times Y[/math] общее количество нейронов в слое Кохонена. И каждому нейрону [math]n_j[/math] поствлен в соответствие вектор весов [math]\mathbf{w}_{j} = \left( w_{j}^{(1)}, \dots, w_{j}^{(D)} \right) \in \mathbb{R}^D[/math] и вектор [math]\mathbf{r}_j \in \mathbb{R}^2[/math], определющий его положение на двумерной плоскости.

Алгоритм обучения (online version):

  1. Задать начальные приближение весов [math]\mathbf{w}_{j}[/math], [math]\forall j=1, \dots, L[/math].
  2. Выбрать случайно номер объекта [math]i \sim Unif[1, N][/math].
  3. Найти расстояние между объектом [math]\mathbf{x}_i[/math] и всеми векторами весов [math]\mathbf{w}_{j}[/math]: [math]\rho(\mathbf{x}_i, \mathbf{w}_{j}) = \left\| \mathbf{x}_i - \mathbf{w}_{j} \right\|^2, \quad \forall j=1, \dots, L[/math]
  4. Найти нейрон-победитель [math]n_b[/math] наиболее близкий к текущему объекту [math]\mathbf{x}_i[/math], где [math] b = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j}).[/math]
  5. Для всех неронов [math]n_j[/math] пересчитать их веса по следующей формуле: [math]\mathbf{w}_{j} \leftarrow \mathbf{w}_{j} + \alpha(t)h_{b, j}(t)(\mathbf{x}_i - \mathbf{w}_{j}),[/math] где [math]\alpha(t)[/math] - темп обучения (learning rate), монотонно убывающая функция от номера итерации [math]t[/math], [math]h_{b, j}(t)[/math] - функция определяющая "меру соседства" нейронов [math]n_b[/math] и [math]n_j[/math].
  6. Проверить критерий останова. При необходмости перейти на шаг 2.

Данное выше описание является довольно общим и не раскрывает некоторых деталей:

  • Начальное приближение для [math]\mathbf{w}_j[/math]
    • Инициализация случайными значениями. Самый простой способ - [math]w_j^{(d)} \sim Unif[0, 1][/math]. Данный подход может потребовать большое количество итераций, так как вектора весов неупорядочены и могут находится на значительном расстоянии от точек обучающей выборки.
    • Инициализация объектами обучающей выборки. Данный способ может позволить уменьшить количество итераций, так как вектора весов будут изначально находятся в области пространства, в которой расположены объекты выборки. Но проблема с неупорядоченностью остается.
    • Инициализация при помощи PCA. При помощи PCA находятся первые две главные компоненты [math]e_1, e_2[/math]. Затем в плоскости натянутой на ветора [math]e_1, e_2[/math] выбирается регулярная сетка точек, которыми и инициализируются вектора [math]\mathbf{w}_j[/math]. Этот способ является наиболее предпочтительным и на практике показывает хорошие результаты.
  • Задание [math]\alpha(t)[/math]
    • [math]\alpha(t) = \alpha_0 \exp \left\{ -\frac{t}{\lambda} \right\}[/math]
    • [math]\alpha(t) = \frac{1}{t}[/math]
    • [math]\alpha(t) = \lambda^{\frac{t}{T}}, \lambda \in (0, 1)[/math]
  • Задание [math]h_{b,j}(t)[/math]
    • Наиболее популярный способ: [math]h_{b,j}(t) = \exp \left\{ - \frac{\| \mathbf{r}_b -\mathbf{r}_j \|^2}{2\sigma(t)^2} \right\},[/math] где [math]\sigma(t)[/math] - убывающая функция, которую можно интерпретировать как величину пропорциональную радиусу соседства. Концептуально принцип действия такой меры соседства можно представить следующим образом: пусть в точке [math]\mathbf{r}_b[/math] задана окружность с радиусом пропорциональным [math]\sigma(t)[/math], тогда за пределом этой окружности степень соседства становится пренебрежимо мала. Пример: [math]\sigma(t) = \frac{\sigma_0}{1 + \frac{t}{T}},[/math] где [math]\sigma_0 = \max \left\{ \frac{X}{2}, \frac{Y}{2} \right\} [/math], [math]T[/math] - максимальное количество итераций. Для достижения лучших результатов на практике процедуру обучения разбивают на два этапа: первый этап - итерации выполняются с очень большим радиусом, что позволяет сойтись в хоррошую область быстро, а второй этап - радиус выбирается маленьким и постепенно уменьшается, что дает возможность весам более точно уловить зависимости из обучающей выборки.
    • Существуют более разреженные меры сходства, например, присваивающие ненулевые веса, только внутри некоторого радиуса. Такие меры сходства могут значительно ускорить вычисления, так как не требуется вычислять меру сходства для всех нейронов.

Кроме онлайн алгоритма обучения часто на практике применяют более эффективную с точки зрения вычислений batch-версию. Пусть выбрано [math]B \in \mathbb{N}[/math]. Тогда batch-версия алгоритма описывается следующим образом:

Алгоритм обучения (batch version):

  1. Задать начальные приближение весов [math]\mathbf{w}_{j}[/math], [math]\forall j=1, \dots, L[/math].
  2. Случайно выбрать множество [math]\mathcal{B} \subseteq \{1, \dots, N \}[/math] размера [math]B[/math].
  3. Для каждого объекта [math]\mathbf{x}_i[/math], [math]i \in \mathcal{B}[/math], найти нейрон-победитель [math]n_{b(i)}[/math], где [math] b(i) = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j})[/math]
  4. Пересчитать веса для всех нейронов [math]\mathbf{w}_j \leftarrow \frac{\sum_{i \in \mathcal{B}}h_{b(i),j}\mathbf{x}_i}{\sum_{i \in \mathcal{B}}h_{b(i),j}}, \forall j=1, \dots, L[/math].
  5. Проверить критерий останова. При необходимости перейти на шаг 2.

Интерпретация шага 2 состоит в следующем: веса нейронов преставляют собой выпуклую комбинацию объектов с номерами из [math]\mathcal{B}[/math]. Часто на практике полагают [math]\mathcal{B} = \{1, \dots, N\}[/math], то есть на каждой итерации мы совершаем полный проход по всей выборке объектов.

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

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

Рассмотрим batch-версию алгоритма. Основные вычисления приходятся на шаг 2, где необходимо вычислять расстояния между всеми [math]\mathbf{x}_i[/math], [math]i \in \mathcal{B}[/math] и весами нейронов [math]\mathbf{w}_j[/math], [math]j \in \{ 1, \dots, L \}[/math] и затем искать минимум:

[math]d_{ij} = \left\| \mathbf{x}_i - \mathbf{w}_j \right\|^2 = \left\|\mathbf{x}_i\right\|^2 - 2 \mathbf{x}_i^{T}\mathbf{w}_j + \left\|\mathbf{w}_j\right\|^2[/math]

Так как нас будет интересовать минимум по [math]j[/math], то первое слагаемое [math]\left\|\mathbf{x}_i\right\|^2[/math] можно не учитывать. Тогда все что нам нужно это вычислить матрицу квази-расстояний [math]\mathbf{D} = \| \hat{d}_{ij} \|_{B \times L}[/math], где [math]\hat{d}_{ij} = \left\|\mathbf{w}_j\right\|^2 - 2 \mathbf{x}_i^{T}\mathbf{w}_j [/math]

Сложность вычисления [math]\mathbf{D}[/math] составляет [math]O\left(L D B \right)[/math], учитывая что [math]\left\|\mathbf{w}_j\right\|^2[/math] можно вычислить за один проход по всем нейронам.

Еще одной вычислительно интенсивной частью является вычисление [math]\mathbf{w}_j[/math]. Здесь мы можем вычислить все необходимые [math]h_{b,j}(t)[/math], а затем, используя эти значения, найти новые приблидения для [math]\mathbf{w}_j[/math].

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

Пусть [math]\mathbf{X}_{\mathcal{B}}[/math] - матрица, строками которой являются объекты [math]\mathbf{x}_i[/math] таких, что [math]i \in \mathcal{B}[/math], а [math]\mathbf{W}[/math] - матрица, строками которой являются веса нейронов.

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

  1. [math]\mathbf{D} \leftarrow ComputePairwiseDistances(\mathbf{X}_{\mathcal{B}}, \mathbf{W})[/math]
  2. [math]\mathbf{b} \leftarrow FindClosestNeurons(\mathbf{D})[/math]
  3. [math]\mathbf{W} \leftarrow ComputeWeightedSumsAcrossObjects(\mathbf{X}_{\mathcal{B}}, \mathbf{b})[/math]

Здесь мы дали описание макроопераций в терминах функций, возвращающих соответствующии значения.

Функция [math]ComputePairwiseDistances[/math] вычисляет квазирасстояния между объектами и нейронами и возвращает матрицу [math]\mathbf{D} \in \mathbb{R}^{B \times L}[/math].

Функция [math]FindClosestNeurons[/math] в каждой строке матрицы [math]\mathbf{D}[/math] находит минимальное значение и возвращает вектор [math]\mathbf{b} \in \mathbb{R}^B[/math], такой что элемент в позиции [math]i[/math] соответствует номеру нейрона-победителя для объекта [math]\mathbf{x}_i[/math].

Функция [math]ComputeWeightedSumsAcrossObjects[/math] вычисляет новые значения весов для каждого нейрона, суммируя строки матрицы [math]\mathbf{X}_{\mathcal{B}}[/math] с весами пропорциональными их мере соседства с соответствующим нейроном-победителем. Неявно данная функция также использует вектор положений [math]\mathbf{r}_j[/math] нейронов на двумерной плоскости, но так как эти положения остаются неизменными в продолжении всего обучения, мы не будем указывать их явно в аргументах.

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

  1. [math]\mathbf{x}^{T} \mathbf{y}[/math] - скалярное произведение двух вектров;
  2. [math]\alpha \mathbf{x}[/math] - умножение вектора на скаляр;
  3. [math]\mathbf{x} + \mathbf{y}[/math] - поэлементное сложение двух векторов;
  4. [math]\min_i x_i[/math] - поиск минимального элемента в векторе [math]\mathbf{x}[/math];
  5. Вычисление функции [math]h_{b,j}(t)[/math].

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

В этой секции, используя макрооперции из предыдущего раздела, дадим очень краткий псевдокод batch-версии алгоритма, который легко позволяет понимать основные строительные блоки данного алгоритма:

Рисунок 3. Способ инициализации [math]\mathbf{r}_j[/math] как узлов шестиугольных ячеек на плоскости

[math] \begin{align} \mathbf{W} := &InitializeWeights() \\ \mathbf{R} := &InitializeNeuronPositionsOnPlane() \\ \mathbf{FOR} & \ t=1 \ \mathbf{TO} \ T \ \mathbf{DO} \\ & \mathcal{B} \leftarrow RandomSample\left( \{1, \dots, N\} \right) \\ & \mathbf{D} \leftarrow ComputePairwiseDistances(\mathbf{X}_{\mathcal{B}}, \mathbf{W}) \\ \quad & \mathbf{b} \leftarrow FindClosestNeurons(\mathbf{D}) \\ \quad & \mathbf{W} \leftarrow ComputeWeightedSumsAcrossObjects(\mathbf{X}_{\mathcal{B}}, \mathbf{b}) \\ \end{align} [/math]


Здесь мы использовали [math]\mathbf{R} \in \mathbb{R}^{L \times 2}[/math] как обозначение для матрицы, строками которой являются [math]\mathbf{r}_j[/math]. Способ инициализации этих векторов для решетки с шестиугольными ячейками показан на Рис. 3. Кроме того, стоит отметить, что при нам никто не запрещает посчитать расстояния между положениями нейронов на двумерной сетки заранее, до основного цикла алгоритма. Это позволит на при подсчете меры соседства [math]h_{b,j}(t)[/math] использовать уже вычисленные [math]\rho(\mathbf{r}_b, \mathbf{r}_j)[/math].

Как видно из алгоритма, основной структурой данных, которая нам потребуется, будет обычная матрица, а также вектора. Стоит отметить, что алгоритм в целом неплохо оптимизируем в терминах матричных и веторных операций. Например, вычисление матрицы [math]\mathbf{D}[/math] на языке Python с использованием библиотеки оптимизированных матричных вычислений NumPy будет выглядить следующим образом:

  W_norms = numpy.linalg.norm(W, axis=1)^2 # Вычисляем квадрат нормы для каждого вектора весов
  D = W_norms.T - 2 * numpy.dot(X, W.T) # Вычисляем матрицу квазирасстояний

А вычисление вектора [math]\mathbf{b}[/math] можно реализовать всего одной строчкой:

  b = D.argmin(axis=1)

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

В разделе 1.3 мы уже попытались дать некоторые оценки сложности алгоритма обучения самоорганизующейся карты Кохонена. Ниже мы рассмотрим последовательную сложность более подробно.

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

Пусть, как и ранее, [math]N[/math] - количество объектов в обучающей выборке, [math]L[/math] - количество нейронов в слое Кохонена, [math]D[/math] - размерность признакового пространства объектов (размерность векторов [math]\mathbf{x}_i[/math]) и [math]B[/math] - количество объектов в случайной подвыборке отбираемой на каждой итерации алгоритма.

Оценим сложность одной итерации алгоритма:

  1. Функция [math]ComputePairwiseDistances[/math].
    1. [math]L \times D[/math] умножений и [math]L \times (D - 1)[/math] сложений для вычисления всех норм [math]\left\|\mathbf{w}_j\right\|^2[/math];
    2. [math]B \times L \times D[/math] умножений и [math]B \times L \times (D - 1)[/math] сложений для вычисления всех скалярных произведений [math]\mathbf{x}_i^{T}\mathbf{w}_j[/math];
    3. [math]B \times L[/math] сложений для вычисления [math]\hat{d}_{ij}[/math];
    4. Итого: [math]L \times D + L \times (D - 1) + B \times L \times D + B \times L \times (D - 1) + B \times L = 2 \times L \times D \times (B + 1) - L[/math] арифметических операций над вещественными числами. Тогда асимптотическая сложность равна [math]O(B \times L \times D)[/math];
  2. Функция [math]FindClosestNeurons[/math].
    1. В этой функции используются только операции сравнения между вещественными числами и их количество составляет [math]B \times (L - 1)[/math] или асимптотически [math]O(B \times L)[/math];
  3. Функция [math]ComputeWeightedSumsAcrossObjects[/math].
    1. Будем считать что для вычисления [math]h_{b, j}(t)[/math] требуется некоторое константное число операций [math]Const[/math]. Тогда для вычисления всех таких функций требуется [math]B \times L \times Const[/math] операций.
    2. Вычисление [math]\mathbf{w}_j[/math] делится на вычисление числителя и выисление знаменателя. В первом случае нам требуется выполнить [math]B \times D[/math] умножений и [math](B - 1) \times D[/math] сложений. Во втором случае мы выполняем [math]B[/math] сложений.
    3. Итого: [math](2 \times B \times D + B - D + 1) \times L + B \times L \times Const[/math] операций. Асимптотическая сложность: [math]O(B \times L \times D)[/math]

В итоге асимптотическая сложность одной итерации будет составлять [math]O(B \times L \times D)[/math].

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

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

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

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

  1. sq_norm - вычисление квадрата нормы [math]\| \mathbf{w}_j \|^2 = \mathbf{w}_j^{T} \mathbf{w}_j[/math]
  2. distance - вычисления квазирасстояния между объектом и вектором весов [math]\| \mathbf{w}_j \|^2 - 2 \mathbf{x}_i^{T} \mathbf{w}_j[/math]
  3. min - поиск ближайшего нейрона на основе расстояний, которые получены на входе операции
  4. h [math]\ast[/math] x - вычисление веса соседства и умножения этого веса на вектор [math]\mathbf{x}[/math] поэлементно
  5. update - суммирование взвешенных векторов [math]\mathbf{x}[/math] и обновление весов нейронов

На Рис. 4 и Рис. 5 представлен информационный граф алгоритма для трех объектов [math]\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3[/math] и двух нейронов с весами [math]\mathbf{w}_1, \mathbf{w}_2[/math] в разных проекциях на плоскость.

Рисунок 4. Информационный граф batch-алгоритма обучения самоорганизующейся карты Кохонена с тремя входными объектами [math]\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3[/math] и с двумя нейронами с весами [math]\mathbf{w}_1, \mathbf{w}_2[/math]. Проекция 1.
Рисунок 5. Информационный граф batch-алгоритма обучения самоорганизующейся карты Кохонена с тремя входными объектами [math]\mathbf{x}_1, \mathbf{x}_2, \mathbf{x}_3[/math] и с двумя нейронами с весами [math]\mathbf{w}_1, \mathbf{w}_2[/math]. Проекция 2.

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

Наиболее очевидными способами использования параллельных ресурсов для ускорения алгоритма являются:

  • распределить объекты обучающей выборки по процессам (распределение данных);
  • распределить веса нейронов по процессам (распределение сети);

В первом случае в каждом узле вычислительной сети мы храним часть обучающей выборки и всю нейронную сеть целиком (все веса [math]\mathbf{w}_j[/math]), а во втором в каждом узле хранится вся обучающая выборка и некоторая подсеть исходной нейронной сети.

Допустим неограниченное количество вычислительных узлов.

Рассмотрим первый случай. Имея в своем распоряжении неограниченное количество вычислительных узлов мы можем полностью распараллелить поиск минимальных расстояний, который теперь можно сделать за [math]O(L \times D)[/math] арифметических операций с вещественными числами. Поиск ближайшего нейрона параллельно можно сделать за [math]O(L)[/math] сравнений. Сложность вычисления [math]\mathbf{w}_j[/math] - [math]O(B \times L \times D)[/math], так как нам необходима пересылка данных одному процессу, который вычислит новые значения для весов. Таким образом асимптотически мы не получили улучшения, но если рассмотреть точное количество операций при, то улучшение есть.

Рассмотрим второй случай. Расстояния мы по прежнему можем считать параллельно - [math]O(B \times D)[/math]. Но поиск минимального расстояния нужно выполнять последовательно, так как у одного рассматриваемого процесса нет расстояний до всех нейронов, поэтому сложность - [math]O(B \times L)[/math]. Сложность обновления нейронов - [math]O(B \times D)[/math]. В итоге получаем [math]O(B \times ( L + D) )[/math].

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

Входными данными алгоритма являются:

  1. множество [math]\mathcal{L} = \left\{ \mathbf{x}_i \right\}_{i=1}^N[/math] - обучающая выборка. Наиболее просто можно представить обучающую выборку как матрицу [math]\mathbf{X} \in \mathbb{R}^{N \times D}[/math]. Важно понимать, что на практике часто встречаются большие выборки состоящие из сотней тысяч объектов, каждый из которых может иметь десятки или сотни признаков (компонет), поэтомк размеры матрицы могут быть значительными;
  2. количество нейронов [math]L[/math], кроме того явно конкретизированы размеры выходной карты [math]X, Y[/math];
  3. положения нейронов на двумерной карте [math]\mathbf{r}_j[/math], которые мы можем упорядочить в матрицу [math]\mathbf{R}[/math]. На самом деле задав размеры сетки и форму ячеек, мы неявно уже определили расстояния между положениями нейронов, поэтому, вообще говоря, описывать матрицу [math]\mathbf{R}[/math] как входные данные излишне, но для ясности сделаем это;
  4. количество итераций [math]T[/math];
  5. начальный радиус [math]\sigma_0[/math].

Общая память, которая требуется нам для хранения входных данных, равна [math]O(N \times D)[/math].

Выходными данными алгоритма являются:

  1. обученные веса нейронов [math]\mathbf{w}_j[/math], которые удобно хранить в виде матрицы [math]\mathbf{W}[/math].

Общая память, которая требуется нам для хранения выходных данных, равна [math]O(L \times D)[/math].

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

Вычислительная мощность алгоритма равна [math]O\left(\frac{BL}{N + L}\right)[/math], что при [math]B=N[/math] равно [math]O\left(\frac{NL}{N + L}\right)[/math], что асимптотически при большом размере выборке равно [math]O(1)[/math].

Задача, которую решают самоогранизующиеся карты, является чрезвычайно сложной и некорректно поставленной, поэтому важной составляющей алгоритма является случайность. Она проявляется в случайном порядке выбора объектов (или случайных подмножеств объектов) для обучения. Но данное стохастическое поведение алгоритма не имеет существенного влияния на параллельную структуру алгоритма, кроме того если в batch-версии выбрать [math]B=N[/math], то случайность будет устранена.

Рассматрим свойства параллельной реализации. Основной проблемой реализации параллельной версии алгоритма является то, что нам необходимо выполнять массово две операции:

  • искать ближайший нейрон к текущему объекту;
  • обновлять веса каждого нейрона в конце итерации.

Так, если мы выполнили распределение данных, то мы имеем часть объектов, но в тоже время у нас есть локально вся сеть, следовательно найти ближайший нейрон для каждого объекта из известной нами подвыборки не требует взаимодействия с другими процессами. Но для того чтобы обновить значения весов нейронов нам необходима информация от других вычислительных узлов, в частности нам требуются [math]h_{b(i), j}[/math] и [math]\mathbf{x}_i[/math], чтобы обновить веса нашей сети. Следовательно, нам необходим отдельный процесс-аггрегатор, который, получив все необходимые данные, обновил значения весов и переслал бы обновленную нейронную сеть всем остальным процессам. В целом, такой подход полностью укладывается в парадигму MapReduce.

Если же мы выбрали стратегию распределения сети, то теперь мы имеем только часть сети, но зато вся обучающая выборка. Теперь найти ближайший нейрон не так просто, так как локально мы можем проверить только те нейроны, которые хранятся в памяти процесса, поэтому остальные процессы должны переслать нам их минимальное расстояние. Далее зная номер ближайшего нейрона (при условии, что мы храним локально все положения нейронов [math]\mathbf{R}[/math]), мы можем без труда найти выпуклую комбинацию объектов выборки и обновить веса нейронов.

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

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

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

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

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

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

В данном разделе мы рассмотрим паралелльную реализацию алгоритма, которая основана на разделении входной выборки по всем доступным процессам. Каждый процесс работает со своим кусоком данным и в конце каждой итерации все процессы синхронизируют состояние весов нейронов. В качестве допущения мы будем рассматривать только одну итерацию алгоритма (если итераций больше чем одна, то результаты просто будут отличаться на множитель равный количеству итераций). Кроме того, допустим [math]B = N[/math].

Все эксперименты были выполнены на суперкомпьютере IBM Blue Gene/P. Выборка была сгенерировна из стандартного нормального распределения для [math]N = D = 1000[/math]. Набор изменяемых параметров запуска:

  1. [math]L[/math] (число нейронов с решетке) = [ 1024, 4096, 9216, 16384, 25600, 36864, 50176, 65536]
  2. [math]P[/math] (число процессоров) = [ 8, 16, 32, 48, 64, 128, 384, 256, 512]

Для оценки производительности (GFLOPS) была использована оценка числа операций приведенная ранее.

На Рис 6 показана зависимость производительности алгоритма от числа процессоров и количества нейронов в решетке. График подтверждает наше интуитивную догадку, что для больших размеров решетки алгоритм имеет большую производительность в терминах операций с плавающей точкой. Максимальная производитлеьность 24 GFLOPS достигается на 512 процессорах и для 50176 (224 * 224) нейронов. Минимальная же производительность (0.6 GFLOPS) наблюдается на 8 процессорах при любом числе нейронов.

Рисунок 6. Производительность алгоритма SOM

На Рис 7 показана аналогичная зависимость эффективности реализации. Здесь стоит отметить, что эффективность реализации очень низкая (в среднем 2% от пиковой производительности). Это долольно просто объяснить: в алгоритме требуется передовать больше количество данных (веса нейронов) для синхронизации процессов, поэтому значительную часть времени процессоры простаивают, ожидая передачи или отправки данных. Данная ситуация ухудшается, когда нейронов становится больше, что отлично видно на графике. Максимальная эффектвность (2.7 %) была достигнута для 8 процессов и остается неизменной при любом числе нейрнов, в данном случае более высокая эффективность достигается за счет того, что требуется осуществить синхронизацию только между 8 процессами, что намного дешевле, чем в случае 512 процессоров (минимально достигнутая эффективость - 1%).

Рисунок 7. Эффективность алгоритма SOM

В целом можно отметить, что алгоритм обладает свойством сильной масштабируемости: на Рис 6 можно наблюдать сублинейный рост производительности с ростом числа процессоров при фиксированном числе нейронов. Так же заметим, что алгоритм проявляет свойство слабой масштабируемости: так, например, при [math]W / p = 128[/math] отмечаем следующий рост производительности: [0.61, 2.17, 7.92, 21.4] при соответствующих значенияз числа процессоров: [8, 32, 128, 512]. Другой пример: [math]W / p = 64[/math], производительности - [1.21, 4.17, 12.96] при [16, 64, 256].

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

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

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

Существующие параллельные реализации:

  • Somoclu: Massively parallel self-organizing maps: accelerate training on multicore CPUs, GPUs, and clusters (GNU GPLv3).
    • поддержка OpenMP, MPI и CUDA;
    • интерфейсы для Python, R, Julia, и Matlab
  • Parallel-SOM: Parallel Implementation of Self Organzing Map in OpenCL
  • MapReduce-MPI SOM: Parallel implementation of Self-Organizing Map (SOM) algorithm with MapReduce-MPI (GNU GPLv3).

Существующие последовательные реализации:

  • SOM-Toolbox: A Matlab toolbox for Self-Organizing Maps (SOM) and others (GNU GPLv2).
    • предоставляет инструменты для удобной визуализации полученной карты.
  • MiniSom: A minimalistic implementation of the Self Organizing Maps (Creative Commons Attribution 3.0 Unported License).
    • реализация на языке Python с использованием библиотеки NumPy.
  • SOMPY: A Python Library for Self Organizing Map (Apache License v2).
  • KNNL: C++ Kohonen Neural Network Library (BSD License).
  • JKNNL: Java Kohonen Neural Network Library (BSD License).
  • kohonen4j: Self-Organizing Maps in Java (GNU GPLv3).
  • kohonen Supervised and Unsupervised Self-Organising Maps in R (GNU GPLv3).

3 Литература

  1. Teuvo Kohonen (1997). Self-Organizing Maps. Springer-Verlag New York, Inc., Secaucus, NJ, USA.
  2. Peter Wittek, Shi Chao Gao, Ik Soo Lim, Li Zhao (2015). Somoclu: An Efficient Parallel Library for Self-Organizing Maps. arXiv:1305.1422.
  3. Juha Vesanto (1999). SOM-based data visualization methods. Intell. Data Anal. 3, 2 (March 1999), 111-126.
  4. Timo D. Hämäläinen (2001). Parallel implementation of self-organizing maps. In Self-Organizing neural networks, Udo Seiffert and Lakhmi C. Jain (Eds.). Springer Studies In Fuzziness And Soft Computing Series, Vol. 78. Springer-Verlag New York, Inc., New York, NY, USA 245-278.
  5. Seung-Jin Sul, Andrey Tovchigrechko (2011). Parallelizing BLAST and SOM Algorithms with MapReduce-MPI Library. IEEE International Symposium on Parallel and Distributed Processing.