О SOM

Представленные финским профессором Теуво Кохоненом в 1980-х годах, самоорганизующиеся карты или SOM предоставляют средства для более низкоразмерного и дискретного представления, называемого картой, наборов данных, сохраняя при этом топологию данных.

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

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

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

Давайте разберем это:

"Самоорганизация"

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

"Карта"

Полученные в результате полученные веса представлены U-матрицей (матрица унифицированных расстояний, подробно описанная ниже), «картой», которая представляет расстояние между нейронами. Идея обучения SOM заключается в том, что сформированная карта аналогичным образом реагирует на аналогичные входные данные, при этом аналогичные входные данные группируются вместе в определенных областях результирующей карты.

Реализация на Python

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



Цель состоит в том, чтобы изучить массив весов W , который в нашем случае является трехмерным массивом, но интерпретируется как двумерный массив нейронов размером length x width, где каждый нейрон представляет собой одномерный вектор длины n_features . Мы изучаем 2D-представление набора данных. Затем W используется для создания U-матрицы — единой 2D-карты всего набора данных.

Каждый входной вектор используется для обновления W . Ближайший нейрон W к точке данных — это наилучшая совпадающая единица (BMU). Расстояния остальных нейронов от BMU используются для обновления функции соседства, которая является основой обновления W. Большие значения в W представляют кластеры похожих входных векторов. Скорость обучения и радиус функции окрестности уменьшаются со временем по мере того, как окрестности становятся меньше, т. Е. Похожие входные данные группируются ближе друг к другу.

Во-первых, импортируйте все необходимые библиотеки.

from sklearn.preprocessing import MinMaxScaler, StandardScaler
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import random
import time
from itertools import product

Для данных мы используем набор данных, используемый для прогнозирования здоровья плода, который классифицирует здоровье плода как Normal, Suspect или Pathological с использованием данных КТГ. У него 22 функции, но мы сосредоточимся только на первых 10:

filein = "fetal_health.csv"
data = np.loadtxt(open(filein, "rb"), delimiter = ",", skiprows = 1)
features = ['baseline value', 'accelerations', 'fetal_movement', 
            'uterine_contractions', 'light_decelerations', 
            'severe_decelerations', 'prolongued_decelerations', 
            'abnormal_short_term_variability', 
            'mean_value_of_short_term_variability', 
            'percentage_of_time_with_abnormal_long_term_variability']
X_data = np.array(data)[:,0:len(features)]

print("Number of datapoints: ", X_data.shape[0])
print("Number of features: ", X_data.shape[1])

И, поскольку он использует метрику на основе расстояния, мы продолжим и масштабируем входные данные:

scaler = MinMaxScaler()
X_data = scaler.fit_transform(X_data) # normalizing, scaling

Теперь углубимся в детали реализации алгоритма.

Во-первых, мы посмотрим на гиперпараметры:

n_features = X_data.shape[1]
max_iter = 1000                                   # max. no. of iterations
length = 50                                       # length of W
width = 50                                        # width of W
W = np.random.rand([length, width, n_features])   # initial weights
aplha0 = 0.9                                      # initial learning rate
beta = [[]]                                       # neighbourhood function
sigma0 = max(length, width)                       # initial radius of beta
TC = max_iter/np.log(sigma0)                      # time constant
  • length, width, n_feautures — размеры матрицы весов, Wгде lengthиwidthзадают размер карты признаков, причем большее значение дает более детализированный результат. W инициализируется случайными значениями.
  • max_iter — максимальное количество итераций, на которых мы хотим запустить алгоритм.
  • TC — постоянная времени, используемая для скорости обучения и затухания радиуса окрестности.
  • alpha0 — начальное значение скорости динамического обучения alpha. Скорость обучения затухает как экспоненциальная функция временного шага или номера итерации t :
def aplhaUpdate(t):
    return aplha0*np.exp(-t/TC)
  • sigma0 — начальное значение динамического (сокращающегося) радиуса функции соседства. Мы начинаем с радиуса, который покрывает всю 2D-карту нейронов в W. Обновляется так же, как alpha .
def sigmaUpdate(t):
    return sigma0*np.exp(-t/TC)
  • beta Двумерный вектор, в котором хранятся значения функции соседства для каждого узла из BMU .
  • BMUB наибольшая M привязка U nit — это узел, ближайший к точке входных данных x_n, т. е. имеющий кратчайшее расстояние от нее, в нашем случае евклидово расстояние, как мы видим ниже.
def winningNode(x_n, W): # Find the winning node and distance to it
    dist = np.linalg.norm(x_n - W, ord =2 , axis = 2)
    d = np.min(dist)
    BMU = np.argwhere(dist == np.min(dist))
    return d, BMU
  • Чтобы обновить beta, нам сначала нужно найти расстояние каждого узла в весовом векторе от BMU. hardDist — это вектор, представляющий евклидово расстояние каждого длинного n_feature вектора в W до любого другого вектора с точки зрения позиций векторов в W, где каждый вектор представлен своими x и y координаты, coords .
hardDist = np.zeros([length, width, length, width])
for i,j,k,l in product(range(length),range(width),range(length),range(width)):
      hardDist[i,j,k,l]= ((i-k)**2 + (j-l)**2)**0.5
  • gridDistances используется для доступа к 2D-массиву, представляющему расстояние или каждый вектор от вектора в позиции coords.
def gridDistances(coords):
    return hardDist[coords[0][0],coords[0][1],:,:]
  • Затем мы можем обновить beta . Для каждого входного вектора вычисляется beta, чтобы дать наименьшие значения точкам, ближайшим к BMU .
def betaUpdate(coords, sigma): # neighbourhood function
    return np.exp(-gridDistances(coords)**2/(2*sigma**2))
  • На рисунках ниже показан массив весов W, который можно представить как двумерный массив векторных нейронов на втором рисунке. Вектор, выделенный красным цветом, — это BMU. Синяя поверхность показывает значение beta для каждого нейрона, которое является самым высоким в BMU. Значения beta представляют собой двумерное нормальное распределение.

  • Наконец, мы можем обновить W. Мы видим, что есть три фактора, влияющие на обновление W : (1) alpha , которое является постоянным в течение одной итерации набора данных; (2) f=x_n-W, отличие каждого вектора в W от входного вектора; и (3) beta . Коэффициент alpha*beta*f добавляется к W , поэтому для векторов в W ближе к входу значение beta больше, что приводит к большему значению вектора в W .
def weightUpdate(W, aplha, beta, x_t):
    f= (x_t-W)
    beta = np.reshape(beta,[beta.shape[0],beta.shape[1],1])
    return W+aplha*beta*f

Теперь, когда мы обсудили предварительные условия, мы можем собрать все воедино:

errors=np.zeros([max_iter])

for t in range(max_iter): 

    t1 = time.time()
    alpha = aplhaUpdate(t)
    sigma = sigmaUpdate(t)
    
    np.random.shuffle(X_data)
    n_sample = len(X_data) # iterating over whole dataset
    QE = 0

    for n in range(n_sample):
        d, BMU = winningNode(X_data[n,:], W)
        QE = QE + d
        beta = betaUpdate(BMU, sigma)
        W = weightUpdate(W, alpha, beta, X_data[n,:])
        W /= np.linalg.norm(W, axis=2).reshape((length, width, 1))

    QE /= n_sample
    errors[t]=QE
    t1 = time.time() - t1

    print("Iteration: ", "{:05d}".format(t), " | Quantization Error: ", "{:10f}".format(QE),  " | Time: ", "{:10f}".format(t1))

Давайте рассмотрим это шаг за шагом.

  • Для каждой итерации alpha и sigma слегка затухают, чтобы приспособиться к постепенному изменению W.
  • Для каждого входного вектора в наборе данных мы находим BMU и расстояние BMU от входного вектора, d.
  • beta обновляется в соответствии с BMU, который затем используется для обновления W.
  • W нормируется так, что W при разных значениях t остаются сопоставимыми.
  • Чтобы отслеживать процесс обучения, мы используем ошибку квантования QE, которая представляет собой среднее значение расстояний каждой точки данных до ее BMU. Если W обновляется правильно, оно должно уменьшаться по мере итераций.

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

Мы можем построить карты отдельных объектов, разрезав W по последнему измерению, которое представляет объекты.

Это завершает первый шаг обучения W получению отдельных карт объектов.

Вторая часть — это просто вычисление расстояний между нейронами для получения U-матрицы (Unified-distance matrix), которая может быть как шестиугольной, так и прямоугольной. На изображении ниже мы видим оба способа просмотра двумерного массива нейронов в W.

U-матрица представляет расстояния между нейронами. На рисунках ниже показаны U-матрицы как для шестиугольного, так и для прямоугольного представления. На обоих изображениях клетки желтого цвета представляют расстояния между соседними нейронами, например, {1,2} — это расстояние между нейронами 1 и 2 на изображениях выше. Оранжевые клетки представляют собой среднее расстояние от нейрона с номером до окружающих его нейронов, т. е. {1}=mean({1,2}, {1,4}) или {5}=mean({2,5}, {4,5}, {5,6}, {5,7}, {5,8}, {5,9}) для шестиугольной U -matrix и {5}=mean({2,5}, {4,5}, {5,6}, {5,8}) для прямоугольного.

В прямоугольной U-матрице есть дополнительный тип ячеек, выделенный синим цветом. Мы не можем считать это расстоянием между оранжевыми ячейками по углам синих ячеек, например, средним ({1}, {5}) или средним ({2}, {4}), поскольку оба эти значения будут справедливо для синей клетки, но не равно, что было бы противоречиво. Таким образом, мы принимаем значение как среднее значение желтых ячеек вокруг него, что обеспечивает более плавный переход между ячейками в окончательной U-матрице. Дополнительные разъяснения см. в этой статье. Расчет обоих типов U-матриц приведен ниже.

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

# CALCULATION OF RECTANGULAR U-MATRIX FROM W
def make_u_rect(W):

    U = np.zeros([W.shape[0]*2-1, W.shape[1]*2-1], dtype=np.float64)

    # YELLOW CELLS
    for i in range(W.shape[0]): # across columns
        k=1
        for j in range(W.shape[1]-1):
            U[2*i, k]= np.linalg.norm(W[i,j]-W[i,j+1], ord=2)
            k += 2

    for j in range(W.shape[1]): # down rows
        k=1
        for i in range(W.shape[0]-1):
            U[k,2*j] = np.linalg.norm(W[i,j]-W[i+1,j], ord=2)
            k+=2

    # ORANGE AND BLUE CELLS - average of cells top, bottom, left, right.
    for (i,j) in product(range(U.shape[0]), range(U.shape[1])):
        if U[i,j] !=0: continue
        all_vals = np.concatenate((
               U[(i-1 if i>0 else i): (i+2 if i<=U.shape[0]-1 else i), j],
               U[i, (j-1 if j>0 else j): (j+2 if j<=U.shape[1]-1 else j)]))
        U[i,j] = all_vals[all_vals!=0].mean()
    
    # Normalizing in [0-1] range for better visualization.
    scaler = MinMaxScaler()
    return scaler.fit_transform(U)

Для шестиугольной матрицы мы создаем временную U-матрицу, используя значения самих нейронов в качестве заполнителей в позициях {1}, {2}, …, чтобы упростить вычисление {1,2}, {2,3} и т. д. в U. Затем {1}, {2}, … заменяются средним значением окружающих их значений.

# CALCULATION OF HEXAGONAL U-MATRIX FROM W
def make_u_hex(W): 
    
    # Creating arrays with extra rows to accommodate hexagonal shape.
    U_temp =  np.zeros([4*W.shape[0]-1, 2*W.shape[1]-1, W.shape[2]])
    U = np.zeros([4*W.shape[0]-1, 2*W.shape[1]-1])
    """
    The U matrix is mapped to a numpy array as shown below.
    U_temp holds neuron at postion 1 in place of {1} for easy computation
    of {1,2}, {2,3} etc. in U. {1}, {2}, {3}, ... are computed later.
    [
      [    (1),        0,       0,       0,    (3)],
      [      0,    (1,2),       0,   (2,3),      0],
      [   (1,4),       0,     (2),       0,  (3,6)],
      [      0,    (2,4),       0 ,  (2,6),      0],
      [    (4),        0,    (2,5),      0,    (6)],
      [      0,    (4,5),        0,  (5,6),      0],
      [   (4,7),       0,      (5),      0,  (6,9)],
      [      0,    (2,4),        0,  (5,9),      0],
      [    (7),        0,    (5,8),      0,    (9)],
      [      0,    (7,8),        0,  (8,9),      0],
      [      0,        0,      (8),      0,      0]
    ]
    """

    # Creating a temporary array placing neuron vectors in 
    # place of orange cells.
    k=0
    indices = []
    for i in range(W.shape[0]):
        l=0
        for j in range(W.shape[1]):
            U_temp[k+2 if l%4!=0 else k,l,:] = W[i,j,:]
            indices.append((k+2 if l%4!=0 else k,l))
            l+=2
        k += 4
    
    # Finding distances for YELLOW cells.
    for (i,j),(k,l) in product(indices, indices):
        if abs(i-k)==2 and abs(j-l)==2:      # Along diagonals
            U[int((i+k)/2), int((j+l)/2)] = 
                np.linalg.norm(U_temp[i,j,:]-U_temp[k,l,:], ord=2)
        
        if abs(i-k)==4 and abs(j-l)==0:       # In vertical direction
            U[int((i+k)/2), int((j+l)/2)] = 
                np.linalg.norm(U_temp[i,j,:]-U_temp[k,l,:], ord=2)
    
    # Calculating ORANGE cells as mean of immediate surrounding cells.
    for (i,j) in indices:
        all_vals = 
          U[(i-2 if i-1>0 else i): (i+3 if i+2<U.shape[0]-1 else i), 
            (j-1 if j>0 else j): (j+2 if j<=U.shape[1]-1 else j)]
        U[i,j] = np.average(all_vals[all_vals!=0])
    
    # To remove extra rows introduced in above function.
    new_U = collapse_hex_u(U)

    # Normalizing in [0-1] range for better visualization.
    scaler = MinMaxScaler()
    return scaler.fit_transform(new_U)


def collapse_hex_u(U):

    new_U = np.zeros([int((U.shape[0]+1)/2), U.shape[1]])
    # Moving up values in every alternate column
    for j in range(1, U.shape[1], 2):
            for i in range(U.shape[0]-1):
                U[i,j]=U[i+1,j]

    # Removing extra rows
    for i in range(new_U.shape[0]): new_U[i,:] = U[2*i,:]

    return new_U

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

def draw_hex_grid(ax, U):

    def make_hex(ax, x, y, colour):
        if x%2==1: y=y+0.5
        xs = np.array([x-0.333,x+0.333,x+0.667,x+0.333,x-0.333,x-0.667])
        ys = np.array([  y+0.5,  y+0.5,      y,  y-0.5,  y-0.5,      y])
        ax.fill(xs, ys, facecolor = colour)
    
    ax.invert_yaxis()
    cmap = matplotlib.cm.get_cmap('jet')

    for (i,j) in product(range(U.shape[0]), range(U.shape[1])):
        if U[i,j]==0: rgba='white'
        else: rgba = cmap(U[i,j])
        make_hex(ax, i, j, rgba)
        
    ax.set_title("U-Matrix (Hexagonal)")

Результирующие U-матрицы для обоих типов с использованием W показаны на рис. 3 это:

Делая снимки W как на рисунке 3 во время обучения, мы можем наблюдать за изучением каждой карты признаков и U-матрицы. Цветовая карта, показанная выше, показывает, что красный представляет самое высокое значение, указывающее на большее количество точек входных данных, сгруппированных в этой области, в то время как области размытия представляют пространства, где точки входных данных встречаются редко.

По мере прохождения итераций обучение W замедляется. В начальных итерациях карты признаков сильно различаются, тогда как в более поздних итерациях они становятся несколько стабильными. Та же тенденция наблюдается и с матрицей U с уменьшением размеров окрестностей. Это также отражено в QE, для которого скорость снижения уменьшается со временем.