Сверточная сеть на python. Часть 1. Определение основных параметров модели

qstw6bd1jmstsaadj__koiwuywu.png


Несмотря на то, что можно найти не одну статью, объясняющую принцип метода обратного распространения ошибки в сверточных сетях (раз, два, три, четыре, пять и даже дающих «интуитивное» понимание — шесть), мне, тем не менее, никак не удавалось полностью понять эту тему. Кажется, что авторы недостаточно внимания уделяют обычным примерам либо же опускают какие-то хорошо понятные им, но не очевидные другим особенности, и весь материал по этой причине становится неподъемным. Мне хотелось разложить все по полочкам для самого себя и в итоге конспекты вылились в статью. Я постарался исключить все недостатки существующих объяснений и надеюсь, что эта статья ни у кого не вызовет вопросов или недопониманий. И, может, следующий новичок, который, также как и я, захочет во всем разобраться, потратит уже меньше времени.

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

Чего не будет в этих статьях? Объяснений основ математики и частных производных, подробностей «интуитивного» понимания сути backpropagation (для начала можете прочесть эту отличную статью) или того, как работают сети вообще, в том числе сверточные. Для лучшего понимания материала желательно знать эти вещи и особенно — основы работы нейросетей.
Итак, первая статья.

Конволюция


mqpm5qxvazi8mbb_n3t-ef500vo.jpeg


Иллюстрация выше обозначает основные переменные, которые будут использоваться в дальнейшем изложении.
Давайте рассмотрим формулу конволюции. Но сначала, что мы хотим увидеть в формуле, что она должна отражать? Давайте обратимся к википедии:

«В сверточной нейронной сети в операции свертки используется лишь ограниченная матрица весов небольшого размера, которую «двигают» по всему обрабатываемому слою (в самом начале — непосредственно по входному изображению), формируя после каждого сдвига сигнал активации для нейрона следующего слоя с аналогичной позицией. То есть для различных нейронов выходного слоя используются одна и та же матрица весов, которую также называют ядром свертки… Тогда следующий слой, получившийся в результате операции свертки такой матрицей весов, показывает наличие данного признака в обрабатываемом слое и её координаты, формируя так называемую карту признаков (англ. feature map).»


Значит, формула свертки должна показывать «движение» ядра $w^l$ по входному изображению или карте признаков $y^{l-1}$. Именно это и показывает следующая формула:

$x^l_{ij}=\sum_{a=-\infty}^{+\infty}\sum_{b=-\infty}^{+\infty}w^l_{ab}\cdot y^{l-1}_{(i\cdot s-a)(j\cdot s-b)}+b^l_{ij} \qquad \forall i\in (0,...,N) \enspace \forall j\in (0,...,M)$

Здесь подстрочные индексы $i$, $j$, $a$, $b$ — это индексы элементов в матрицах, а $s$ — величина шага свертки (stride).
Надстрочные индексы $l$ и $l-1$ — это индексы слоев сети.
$x^{l-1}$ — выход какой-то предыдущей функции, либо входное изображение сети
$y^{l-1}$ — это $x^{l-1}$ после прохождения функции активации (например, relu или сигмоида; пункт про функции активации будет немного позже)
$w^l$ — ядро свертки
$b^l$ — bias или смещение (на картинке выше отсутствует)
$x^l$ — результат операции конволюции. То есть операции проходят отдельно для каждого элемента $ij$ матрицы $x^l$, размерность которой $(N,M)$.

Ниже отличная иллюстрация работы формулы конволюции. Синим цветом отображена $y^{l-1}$, зеленым — $x^l$, и серая движущаяся матрица три на три — это ядро свертки $w^l$:

bc055edc77decd9f656eab01d6af5b34.gif

Разобравшись с обозначениями параметров и поняв основной принцип операции свертки (детали формулы конволюции станут понятны ближе к концу статьи, так что, возможно, вам потребуется вернуться сюда позже), можно перейти к важному нюансу этой формулы и самой операции вообще: у ядра существует центральный элемент. Причем центральный элемент может находиться не обязательно в центре (где центр матрицы два на два пикселя?), а вообще в любой ячейке ядра — какой вы захотите. Например, центральный элемент ядра в анимации выше находится в позиции (1,1), и относительного этого элемента и происходит операция свертки. Так, что же за центральный элемент?

Центральный элемент ядра


Итак, индексация элементов ядра происходит в зависимости от расположения центрального элемента. Фактически, центральный элемент определяет начало «координатной оси» ядра свертки. Посмотрите картинку ниже, слева ядро с центральным элементом в нулевой строке и нулевом столбце, справа — в первой строке и первом столбце:

x--m_5swdcjkhaux373ex74oezy.jpeg

Видите, как изменилась индексация элементов? И да, я говорю, что во втором ядре центральный элемент «переместился» в первую строку первого столбца, имея в виду, что это произошло относительно нумерации первого ядра. На самом деле, индексы центрального элемента всегда равны (0,0), и для удобства дальнейшего изложения я буду говорить «старые» координаты, говоря о позиции центрального элемента в левом верхнем углу, и «новые» координаты — о положении центрального элемента в любом другом месте ядра, используя при этом нумерацию «старой» координатной сетки (например, (1,1) для правого ядра свертки на рисунке выше).

Но, возвращаясь к формуле конволюции, что означает — сумма от минус бесконечности до плюс бесконечности? Ведь само ядро имеет вполне определенные размеры и у него нет бесконечного числа элементов. Я находил разные варианты написания формулы, например, $\sum_{u=-a}^{a}\sum_{v=-b}^{b}$ (вот и вот). Также находил и варианты с бесконечностью, как в варианте в начале статьи (вот и вот). Но последний показался мне более «общим» случаем.

Минус в формуле для ядра свертки — это следствие расположения центрального элемента. Нам следует «перебрать» все возможные существующие элементы, и начать можно от минус бесконечности. Или от минус $a$ и $b$. Если элемент по этим индексам не определен для данного ядра, то умножение происходит на ноль, и фактически операции начинаются не от минус бесконечности, а от умноженной на минус один позиции центрального элемента (в нумерации «старых» координат). А закончится операция не на плюс бесконечности, а на разности между количеством элементов ядра по рассматриваемой оси минус индекс центрального элемента (опять-таки в нумерации «старой» координатной оси). Причем последнее значение полученного диапазона не включается (так как индексация от нуля).

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

code_demo_indexes.py
ссылка на git
import numpy as np

size_axis = (3,3)
center_w_l = (1,1)

def create_axis_indexes(size_axis, center_w_l):
	coordinates = []
	for i in range(-center_w_l, size_axis-center_w_l):
		coordinates.append(i)
	return coordinates

def create_indexes(size_axis, center_w_l):
	# расчет координат на осях ядра свертки в зависимости от номера центрального элемента ядра
	coordinates_a = create_axis_indexes(size_axis=size_axis[0], center_w_l=center_w_l[0])
	coordinates_b = create_axis_indexes(size_axis=size_axis[1], center_w_l=center_w_l[1])
	return coordinates_a, coordinates_b

print(create_indexes(size_axis, center_w_l))


На картинке ниже мы объявили центром ядра тот элемент, который находится в позиции (1,1).

vxp3seaugvyvakgxht1v4xuhmvu.jpeg

Но «старые» координаты говорят нам, что позиция центрального элемента должна находиться по индексам (0,0), а значит, необходимо переопределить координатные оси для нового положения центрального элемента.

Если мы подставим в код выше наши значения, то получим заполненный питоновский лист значениями из range (-1, 2), то есть лист будет содержать [-1,0,1]. Еще раз, почему range (-1, 2)? «Минус один» потому, что операция начинается от минус индекса нашего центрального элемента, а «два» получается как длина оси (равная трем) минус индекс центрального элемента в старых координатах (то есть один). Последний элемент range не включается.

Кросс-корреляция


Приведу еще раз формулу конволюции:

$x^l_{ij}=\sum_{a=-\infty}^{+\infty}\sum_{b=-\infty}^{+\infty}w^l_{ab}\cdot y^{l-1}_{(i\cdot s-a)(j\cdot s-b)}+b^l_{ij} \qquad \forall i\in (0,...,N) \enspace \forall j\in (0,...,M)$

И ниже формула кросс-корреляции:

$x^l_{ij}=\sum_{a=-\infty}^{+\infty}\sum_{b=-\infty}^{+\infty}w^l_{ab}\cdot y^{l-1}_{(i\cdot s+a)(j\cdot s+b)}+b^l_{ij}$

Да, отличие только в том, что минусы в расчете подстрочных индексов $y^{l-1}$ заменены на плюсы. На практике, применяя формулу конволюции, мы можем видеть, что ядро при свертке «переворачивается» (причем переворачивается относительно центрального элемента!), в то время как при кросс-корреляции элементы ядра при свертке сохраняют свои позиции. Посмотрите иллюстрацию, чтобы лучше понять, что здесь имеется в виду:

exqvppkt1kufrfnggmwtoouqegu.jpeg


Здесь можно видеть позицию ядра, его расположение во время свертки относительно матрицы $y^{l-1}$. Ниже код, который также выводит на экран демонстрационные примеры, аналогичные тем, что на картинке выше, только уже для всех $i$ и $j$

code_demo_convolution_feed_x_l.py
ссылка на git
import numpy as np

# w_l = np.array([
# 	[1,2,3,4],
# 	[5,6,7,8],
# 	[9,10,11,12],
# 	[13,14,15,16]])

w_l = np.array([
	[1,2],
	[3,4]])

y_l_minus_1 = np.zeros((3,3))

other_parameters={
	'convolution':True,
	'stride':1,
	'center_w_l':(0,0)
}

def convolution_feed_x_l(y_l_minus_1, w_l, conv_params):
	indexes_a, indexes_b = create_indexes(size_axis=w_l.shape, center_w_l=conv_params['center_w_l'])
	stride = conv_params['stride']
	# матрица выхода будет расширяться по мере добавления новых элементов
	x_l = np.zeros((1,1))
	# в зависимости от типа операции меняется основная формула функции
	if conv_params['convolution']:
		g = 1 # операция конволюции
	else:
		g = -1 # операция корреляции
	# итерация по i и j входной матрицы y_l_minus_1 из предположения, что размерность выходной матрицы x_l будет такой же
	for i in range(y_l_minus_1.shape[0]):
		for j in range(y_l_minus_1.shape[1]):
			demo = np.zeros([y_l_minus_1.shape[0], y_l_minus_1.shape[1]]) # матрица для демонстрации конволюции
			result = 0
			element_exists = False
			for a in indexes_a:
				for b in indexes_b:
					# проверка, чтобы значения индексов не выходили за границы
					if i*stride - g*a >= 0 and j*stride - g*b >= 0 \
					and i*stride - g*a < y_l_minus_1.shape[0] and j*stride - g*b < y_l_minus_1.shape[1]:
						result += y_l_minus_1[i*stride - g*a][j*stride - g*b] * w_l[indexes_a.index(a)][indexes_b.index(b)] # перевод индексов в "нормальные" для извлечения элементов из матрицы w_l
						demo[i*stride - g*a][j*stride - g*b] = w_l[indexes_a.index(a)][indexes_b.index(b)]
						element_exists = True
			# запись полученных результатов только в том случае, если для данных i и j были произведены вычисления
			if element_exists:
				if i >= x_l.shape[0]:
					# добавление строки, если не существует
					x_l = np.vstack((x_l, np.zeros(x_l.shape[1])))
				if j >= x_l.shape[1]:
					# добавление столбца, если не существует
					x_l = np.hstack((x_l, np.zeros((x_l.shape[0],1))))
				x_l[i][j] = result
				# вывод матрицы demo для отслеживания хода свертки
				print('i=' + str(i) + '; j=' + str(j) + '\n', demo)
	return x_l

def create_axis_indexes(size_axis, center_w_l):
	coordinates = []
	for i in range(-center_w_l, size_axis-center_w_l):
		coordinates.append(i)
	return coordinates

def create_indexes(size_axis, center_w_l):
	# расчет координат на осях ядра свертки в зависимости от номера центрального элемента ядра
	coordinates_a = create_axis_indexes(size_axis=size_axis[0], center_w_l=center_w_l[0])
	coordinates_b = create_axis_indexes(size_axis=size_axis[1], center_w_l=center_w_l[1])
	return coordinates_a, coordinates_b

print(convolution_feed_x_l(y_l_minus_1, w_l, other_parameters))

Пример выхода работы скрипта
ygyy83-emhed8kexvtkbbxc-uha.jpeg


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

В зависимости от способа свертки — конволюции или кросс-корреляции, различной величины шага и выбора центрального элемента ядра — размерность выходной матрицы $x^l$ может варьироваться. В самом простом случае, при шаге равном единице, размерность матрицы $x^l$ будет равна той, что была у матрицы $y^{l-1}$. Уверен, что общая формула для нахождения этой размерности, то есть именно значений $(N, M)$, учитывая все перечисленные особенности, существует (см., например, документацию tensorflow, но там учтена только возможность различной величины шага), но в представленной выше функции на python изначально размерность матрицы $x^l$ не задана, и к этой матрице добавляются строки и столбцы по мере выполнения расчетов, что, конечно, не является оптимальным решением.

Функции активации


Ниже формулы функций активации, которые можно будет использовать в будущей модели. Фактически, это просто «превращение» $x^l$ в $y^l$ таким образом: $y^l=f(x^l)$. Функция активации позволяет сделать сеть нелинейной, и если бы мы не использовали функции активации (тогда получалось бы, что $y^l=x^l$) или использовали бы линейную функцию, то неважно какое количество слоев тогда было бы в сети: их все можно было бы заменить одним единственным слоем с линейной функцией активации.
Итак, ReLU:

$f_{ReLU}=max(0,x)$

И сигмоида:

$f_{sigmoid} =\frac{1}{1+e^{-x}}$

image

Сигмоида используется только если классов (для задачи классификации) не больше двух: выход модели будет числом от нуля (первый класс) до единицы (второй класс). Для большего же числа классов, чтобы выход модели отражал вероятность этих классов (и сумма вероятностей по выходам сети равнялась единице), используется softmax. Функция выглядит просто, но будут определенные сложности при вычислении формулы для backprop.

$f_{softmax} =y^l_i=f(x^l_i)=\frac{e^{x^l_i}}{\sum_{k=0}^{n}e^{x^l_k}}$

где $n$ — это количество классов.

Слой макспулинга


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

image

Реализация макспулинга в python:

code_demo_maxpool.py
ссылка на git
import numpy as np

y_l = np.array([
	[1,0,2,3],
	[4,6,6,8],
	[3,1,1,0],
	[1,2,2,4]])

other_parameters={
	'convolution':False,
	'stride':2,
	'center_window':(0,0),
	'window_shape':(2,2)
}

def maxpool(y_l, conv_params):
	indexes_a, indexes_b = create_indexes(size_axis=conv_params['window_shape'], center_w_l=conv_params['center_window'])
	stride = conv_params['stride']
	# выходные матрицы будут расширяться по мере добавления новых элементов
	y_l_mp = np.zeros((1,1)) # матрица y_l после операции макспулинга
	y_l_mp_to_y_l = np.zeros((1,1), dtype='= 0 and j*stride - g*b >= 0 \
					and i*stride - g*a < y_l.shape[0] and j*stride - g*b < y_l.shape[1]:
						if y_l[i*stride - g*a][j*stride - g*b] > result:
							result = y_l[i*stride - g*a][j*stride - g*b]
							i_back = i*stride - g*a
							j_back = j*stride - g*b
						element_exists = True
			# запись полученных результатов только в том случае, если для данных i и j были произведены вычисления
			if element_exists:
				if i >= y_l_mp.shape[0]:
					# добавление строки, если не существует
					y_l_mp = np.vstack((y_l_mp, np.zeros(y_l_mp.shape[1])))
					# матрица y_l_mp_to_y_l расширяется соответственно матрице y_l_mp
					y_l_mp_to_y_l = np.vstack((y_l_mp_to_y_l, np.zeros(y_l_mp_to_y_l.shape[1])))
				if j >= y_l_mp.shape[1]:
					# добавление столбца, если не существует
					y_l_mp = np.hstack((y_l_mp, np.zeros((y_l_mp.shape[0],1))))
					y_l_mp_to_y_l = np.hstack((y_l_mp_to_y_l, np.zeros((y_l_mp_to_y_l.shape[0],1))))
				y_l_mp[i][j] = result
				# в матрице y_l_mp_to_y_l хранятся координаты значений,
				# которые соответствуют выбранным в операции максипулинга ячейкам из матрицы y_l
				y_l_mp_to_y_l[i][j] = str(i_back) + ',' + str(j_back)
	return y_l_mp, y_l_mp_to_y_l

def create_axis_indexes(size_axis, center_w_l):
	coordinates = []
	for i in range(-center_w_l, size_axis-center_w_l):
		coordinates.append(i)
	return coordinates

def create_indexes(size_axis, center_w_l):
	# расчет координат на осях ядра свертки в зависимости от номера центрального элемента ядра
	coordinates_a = create_axis_indexes(size_axis=size_axis[0], center_w_l=center_w_l[0])
	coordinates_b = create_axis_indexes(size_axis=size_axis[1], center_w_l=center_w_l[1])
	return coordinates_a, coordinates_b

out_maxpooling = maxpool(y_l, other_parameters)
print('выходная матрица:', '\n', out_maxpooling[0])
print('\n', 'матрица с координатами для backprop:', '\n', out_maxpooling[1])

Пример выхода работы скрипта
unxtwpdhrfd2j3d_6pnzdufnrwc.jpeg


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

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

Слой полносвязной сети


После слоев конволюции мы получим множество карт признаков. Их соединим в один вектор и этот вектор подадим на вход fully connected сети.

qso0uo3uzq2pqjjyeyeplzn6hwa.jpeg

Формула для fc-слоя (fully connected) выглядит так:

$x^l_i=\sum_{k=0}^{m}w^l_{ki}y^{l-1}_k+b^l_i \qquad \forall i\in (0,...,n)$

И в матричном виде (под строкой я написал размерности матриц):

$X^l=Y^{l-1}W^l+B^l_i \\ \tiny (1 \times n) \quad (1 \times m)(m \times n) \quad (1 \times n)$

И вот как эти матрицы выглядят во время вычислений:

$ \begin{aligned} \begin{bmatrix} & x^{l}_0 & x^{l}_1 & x^{l}_2 & ... & x^{l}_n \end{bmatrix} = \begin{bmatrix} & y^{l-1}_0 & y^{l-1}_1 & ... & y^{l-1}_m \end{bmatrix} \begin{bmatrix} & w^{l}_{00} & w^{l}_{01} & ... & w^{l}_{0n} \\ & w^{l}_{10} & w^{l}_{11} & ... & w^{l}_{1n} \\ &... & ... & ... & ... \\ & w^{l}_{m0} & w^{l}_{m1} & ... & w^{l}_{mn} \end{bmatrix} \\ \\ + \begin{bmatrix} & b^{l}_0 & b^{l}_1 & ... & b^{l}_n \end{bmatrix} \end{aligned} $

Функция потерь


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

$ E=\sum_{i=0}^{n}\frac{1}{2}(y^{truth}_i-y^l_i)^2$

где $n$ — это количество классов, $y^l$ — выход модели, а $y^{truth}$ — правильные ответы.
$\frac{1}{2}$ здесь нужна только для сокращения формулы во время обратного распространения ошибки по сети. Можно убрать $\frac{1}{2}$ и ничего принципиально не изменится.
После прочтения этой статьи решил использовать cross-entropy:

$E=-\sum_{i=0}^{n}y^{truth}_iln(y^l_i)$

Структура будущей модели


Теперь, разобрав основные слои сети, мы можем представить примерный вид будущей модели:

  1. функция, которая извлекает из датасета следующее изображение/батч для проведения обучения;
  2. первый слой сверточной сети, который на вход принимает изображение, на выходе отдает карты признаков;
  3. слой макспулинга, который снижает размерность карт признаков;
  4. второй слой сверточной сети принимает полученные на предыдущем шаге карты, и на выходе дает другие карты признаков;
  5. сложение полученных на предыдущем шаге карт в один вектор;
  6. первый слой полносвязной сети принимает вектор, производит вычисления, которые дают значения для скрытого полносвязного слоя;
  7. второй слой полносвязной сети, количество выходных нейронов которого равно количеству классов в используемом датасете;
  8. выход всей модели подается в функцию потерь, которая сравнивает прогнозируемое значение с истинным, и вычисляет разницу между этими значениями.


Итоговая функция потерь является своего рода количественным «штрафом», который можно рассматривать как меру качества прогноза модели. Это значение мы и будем использовать для обучения модели с помощью backpropagation — обратного распространения ошибки. Формулы, которые используют эту ошибку и «протягивают» ее сквозь все слои для обновления параметров и обучения модели, мы рассмотрим в следующей части статьи.

© Habrahabr.ru