Джулия собирается убить Питона? Кто знает; возможно нет. Но это, безусловно, забавный, современный и интересный язык, и он особенно привлекателен для ученых с опытом работы с MATLAB, таких как я. У него уникальное сочетание функций: он индексирован 1, динамически типизирован, но (точно вовремя) скомпилирован, имеет прекрасный встроенный менеджер пакетов и прост в использовании в интерактивном режиме (одновременное использование сценариев и REPL ощущается как очень похоже на разработку MATLAB). Мне очень хотелось погрузиться в знакомство с родной Джулией, и, несмотря на то, что я регулярно использую PyTorch и Tensorflow, я никогда раньше не писал код обратного распространения «вручную» — отсюда и этот проект и пост. Два зайца, один камень.

В этой статье мы создадим простой многослойный персептрон (MLP) и реализуем обратное распространение в родном Julia, чтобы классифицировать цифры в наборе данных MNIST (да, это еще один средний пост MNIST — наденьте на меня наручники). Репозиторий mercury, сопровождающий этот пост, доступен здесь, просто откройте и запустите main.jl. Попутно мы рассмотрим некоторые ключевые функции / парадигмы Julia, такие как некоторые важные структуры данных, динамическая диспетчеризация функций, управление пакетами и построение графиков.

Конечно, уже есть подходящие, оптимизированные библиотеки нейронных сетей/машинного обучения (см., например, flux и отважного, малоизвестного аутсайдера tensorflow). Я надеюсь, что эта статья будет интересна для таких людей, как я, которые хотят увидеть гайки и болты обратного распространения ошибки, а также для Джулии в действии.

Часть 1 посвящена созданию с нуля нейронной сети, способной выполнять простой прямой проход. Мы обучим его с помощью алгоритма обратного распространения в Части 2 и оценим его в Части 3.

Часть 1 — пас вперед

Установка и номенклатура

Мы создадим простой трехслойный MLP, взяв в качестве входных данных изображения MNIST [28 x 28] и выбрав один из 10 классов.

Функциональная форма этой модели:

y = g( W₃ f ( W₂ f( W₁ X+ b₁ ) + b₂ ) + b₃)

В которой Wₙ — наши весовые матрицы, а Bₙ — наши векторы смещения. Функции активации f, g, которые мы используем в этой модели, — это ReLU и Softmax соответственно. y(наш прогноз) аппроксимирует Y (истинное значение) как функцию наблюдений (значений пикселей) X. Сторонники чистоты статистики, пожалуйста, простите это обозначение.

Здесь мы не слишком беспокоимся о том, чтобы получить последний 1% точности, поэтому размерность скрытых слоев является довольно произвольным выбором. Мы пойдем [784, 32, 32, 10] для простоты. W₃ — это матрица 32 x 784, b₂ — вектор 32 x 1, и т. д..

Строительные блоки

Julia на самом деле не предназначен для использования в качестве объектно-ориентированного языка, такого как Python. Хотя вы можете определить struct, это «составной тип данных», а не Python object, и в его определении не определены методы. Например, нет смысла вызывать neural_net.forward(), в котором мы можем получить доступ к neural_net.attribute_1.

Вместо этого в Julia мы полагаемся на динамическую диспетчеризацию, при которой функции определяются уникальным образом для различных типов входных данных. Таким образом, forward(x:Array{Float64,2}) может вызывать совершенно другой код для forward(x:Int). Концепция аналогична перегрузке функций, например. Java, но, как указывает этот интересный пост StackOverflow, перегрузка соответствует решению о вызове функции во время компиляции, в то время как в Julia мы принимаем решение во время выполнения.

Итак, наш ключевой строительный блок MLP — плотный слой:

Мы реализуем это как struct с weight, bias и anact_fn. Обратите внимание, что важно, что это mutable, его атрибуты (параметры MLP) разрешено изменять. Ради качества жизни мы не хотим создавать экземпляр этого слоя с использованием синтаксиса: layer_1 = Dense(W_1, b_1, ReLU()), так как это потребует определения матрицы и вектора в другом месте, а это запутанно. Вместо этого мы перегружаем Dense(), чтобы принимать целые числа размерности в качестве входных данных, и генерируем на месте новую матрицу Kaiming init’d для весов и нули для смещения.

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

Для instancelayer из Dense теперь мы можем вызывать layer.weight , layer.bias , layer.act_fn по мере необходимости.

Нелинейность

Для функций активации мы также используем struct с динамически управляемыми методами forward() и gradient(). Каждый бинарный оператор в Julia имеет векторизованную версию, обозначаемую префиксом .. Это похоже на поведение в MATLAB и является одной из моих любимых особенностей этого языка. Это очень полезно для реализации точечных функций активации.

Прямой метод функции ReLU просто возвращает входное значение, если оно больше нуля, иначе ноль. Градиент равен единице, если вход больше нуля, иначе ноль.

Функция softmax, действующая на вектор, softmax(z), возвращает нормализованную поточечную экспоненту элементов вектора, так что сумма выходных элементов равна единице. Его производная по входу z принимает форму grad_softmax(z) = softmax(z) * (1 — softmax(z)) , где 1 — соответствующий вектор единиц, а * — произведение Адамара ( .* в Джулии). В этой реализации мы используем так называемый трюк log-sum-exp для численной стабильности при вычислении softmax. Отличное объяснение этой техники здесь.

Функция softmax гарантирует, что все элементы выходного вектора имеют значение › 0, действительные значения и их сумма равна единице. Эти свойства являются синонимами допустимого распределения вероятностей по элементам вектора, которое мы решили использовать в качестве модели вероятности отнесения входного вектора x к одному из классов, представленных выходным вектором. элементы y. В этой модели входными данными для функции softmax являются логарифмические шансы, связанные с каждым классом. По сути, это все, что нам нужно для создания функциональной формы нашего MLP.

Собираем это вместе

При обучении нейронной сети нам нужно будет отслеживать промежуточные активации, а также значения самих параметров, и здесь мы просто используем dict (чтобы разрешить O(1) доступ к любому элементу) для этого:

Итак, здесь мы определяем словарь net с полями «Слои», «А» и «Z». «Z» и «A» необходимы для алгоритма обратного распространения, как вы увидите позже.

  • «Слои» — это список отдельных Dense слоев.
  • «A» отслеживает значения активаций слоев посредством прямого прохода после применения функции активации.
  • «Z» отслеживает активации до применения функции активации.

Конструктивно, последнее, что нам нужно сделать, это определить метод forward() для нашего MLP, снова используя родную динамическую диспетчеризацию Julia. Мы просто перебираем сохраненные слои в нашем dict, умножая входные данные на матрицу весов, добавляя смещение и применяя нашу точечную функцию активации. Мы отслеживаем активации как до, так и после применения функции активации.

Так что мы в хорошем месте. Мы можем вызвать y = forward(x, net), чтобы сгенерировать набор прогнозов для класса x, хотя на данный момент это явно случайно, так как нам еще предстоит обучить сеть. Время обратного распространения.

Часть 2 — Обучение

Мы построили функциональную форму нашего MLP: y = forward(x, net) выводит вектор 10 x 1 y из вектора 784 x 1 x. Вскоре мы рассмотрим алгоритм обратного распространения для обновления параметров сети, но сначала нам нужно обсудить данные.

Загрузка и потеря

Хотя ранее я говорил, что для этого проекта мы будем использовать только нативную версию Julia, возможно, это было неправильным направлением. В mercury я использовал пакет MLDatasets для эффективной загрузки обучающих и тестовых данных в виде массивов Float64.

На мой взгляд, управление пакетами в Julia намного проще, чем в Python, и это функция, правильно встроенная в язык, а не запоздалая мысль, как в Python. Среды обрабатываются для каждого проекта отдельно, с помощью файлов Project.toml и Manifest.toml. Они записывают зависимости пакетов и позволяют пользователю точно воссоздать среду Julia, в которой разрабатывался проект. Установка пакетов так же проста, как ввод ] для входа в раздел Pkg REPL, а затем использование, например, команды add MLDatasets.

После установки мы можем импортировать с помощью using MLDatasets. Затем мы можем вызватьtrain_x, train_y = MNIST.traindata(Float64). Данные x хороши для дальнейшего использования, но нам нужно использовать «горячее кодирование» для данных y, поскольку train_y по умолчанию выводит целое число от 0 до 9 в качестве метки класса. Делаем это с помощью функции zi_one_hot_encode(train_y) :

Как обсуждалось ранее, наш слой Softmax выводит вектор, представляющий распределение вероятностей по классам, т.е. его элементы действительны, › 0, и сумма равна единице. Мы хотели бы, чтобы выходные данные нашего MLP, вектора y, представляющего распределение вероятностей по классам, были как можно ближе к горячему кодированию обучающих данных. Затем мы можем минимизировать перекрестную энтропию между этими распределениями вероятностей для обучения нашей нейронной сети.

Для более глубокого изучения кросс-энтропийной потери и ее (близкой) связи с дивергенцией Кульбака-Лейблера (относительная энтропия) см. здесь. Это функция максимального правдоподобия для распределения Бернулли по классам (соответствует нашему горячему кодированию), так же как и правильная функция потерь для классификации с одной меткой. Вот реализация Julia вместе с ее производной, которая понадобится нам для обратного распространения:

Итак, мы, наконец, готовы заняться большим зверем — обратным распространением. Я нашел эту запись в блоге Петра Скальского чрезвычайно полезной при решении этой проблемы, и большая часть моей реализации этого алгоритма вдохновлена ​​​​его работой. Я использовал тот же синтаксис, так что ознакомьтесь с его объяснением, если мое вам не подходит!

Обратное распространение

Алгоритм обратного распространения позволяет эффективно вычислять частные производные функции потерь по параметрам нейронной сети послойно. Мы используем цепное правило, проходя по сети в обратном направлении, начиная с Aₙ = выходы, y и заканчивая A₀ = входы, х.

Мы получаем dWₙ и dbₙ для n = 1: глубина сети. Это градиенты функции потерь по отношению к весам и смещениям слоя n. Мы встроили отслеживание A и Z в определение нашей сети (a dict), поэтому определение обновлений градиента для данного слоя — это всего лишь случай кодирования уравнения на рисунке выше:

Это обобщенная функция, которая принимает в качестве входных данных градиент потерь по отношению к предыдущим (ближе к выходам сети) активациям. Так как окончательная активация является выходом сети, то при первом вызове этой функции
dA = xe_loss_derivative(y, Y). Для последующих (ближе к входам) вызовов calculate_gradient сохраняем dA. Здесь мы снова используем динамическую диспетчеризацию, так как работа функции gradient() зависит от типа данных act_fn, который мы определили/перегрузили выше.

Таким образом, у нас остался набор обновлений параметров dW и db.

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

Когда мы хотим выполнить обновление параметров, это так же просто, как:

Создайте тренировочный цикл

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

Как я уже упоминал выше, обычно при обучении нейронной сети мы используем мини-пакетный градиентный спуск, накапливая средние градиенты для каждого параметра (относительно потерь) по ряду записей (без обновления промежуточных параметров) и обновляя в конце партия. В Python (PyTorch, TensorFlow, etc) это обычно делается путем векторизации операций и использования LAPACK или аналогичных сильно оптимизированных функций линейной алгебры для вашей архитектуры CPU/GPU. Здесь мы можем применить более простой подход, используя эффективность Джулии, в которой гораздо более синтаксически очевидно, как работает мини-пакетная обработка. Мы накапливаем градиент (как описано выше) для каждой записи в мини-пакете, а затем обновляем параметры в конце.

Итак, мы можем идти. Мы можем вызвать train(net, mb_size, lr, epochs, train_x, train_y), и наш MLP будет обучаться, сводя к минимуму потерю перекрестной энтропии. Давайте посмотрим, как это работает.

Часть 3 — Результаты

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

Для этого тренировочного прогона я использовал следующие параметры:

  • mb_size = 32
  • epochs = 50
  • lr = 0.01

Мы можем легко построить результаты обучения, используя пакетPlots:

Первая кросс-энтропийная потеря:

А теперь точность:

Я показал вам здесь хороший выбор гиперпараметров (неудивительно, что это был не первый набор, который я пробовал…). Например, установка слишком высокой скорости обучения может привести к колебанию потерь при обучении вокруг своего минимума, а обучение будет иметь недостаточную точность для параметров, чтобы установить хорошие значения.

Выводы

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

Репозиторий mercury сопровождает этот пост и доступен здесь, просто откройте и запустите main.jl после установки MLDatasets, как описано в этом посте.

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