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

Автор: Гой Антон, 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] b = \arg \min_{\forall j \in \{1, \dots, L\}} \rho(\mathbf{x}_i, \mathbf{w}_{j})[/math], наиболее близкий к текущему объекту [math]\mathbf{x}_i[/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}}[/math], [math]\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} \ MAX\_ITERATIONS \ \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.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].