Градиентный спуск и библиотека Theano

Введение в градиентный спуск

Здравствуйте и вновь добро пожаловать на занятия по теме «Машинное обучение без учителя: скрытые марковские модели на языке Python».

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

Суть в следующем. Пусть у вас есть функция, минимум которой вы хотите найти, и пусть вам нужно найти такие входные данные, при которых функция была бы в минимуме. Как правило, мы хотим минимизировать функцию затрат или ошибок. Может потребоваться также и найти максимумы – например, когда мы ищем максимум для функции правдоподобности некоторого распределения вероятностей. Всё, что нужно сделать – это просто поменять местами стороны. Чтобы объяснить это, рассмотрим очень простой одномерный пример. Как правило, в машинном обучении мы используем размерности, гораздо большие единицы, но этот пример позволит наглядно увидеть суть.

Итак, пусть у нас есть простая функция J=w². Мы знаем, что минимум функции при w=0, но, предположим, мы этого не знаем. Наш весовой коэффициент установим случайным образом. Предположим, w=20. Мы знаем, что производная dJ/dw равна 2w. Установим коэффициент обучения равным 0,1. В первом приближении мы имеем:

w – 0,1*2w = 20 – 0,1*40 = 16.

Поэтому установим w=16. Во втором приближении

w – 0,1*2w = 16 – 0,1*2*16 = 12,8.

Это даёт нам новое значение w=12,8. В третьм приближении

w – 0,1*2w = 12,8 – 0,1*2*12,8 = 10,24.

Как вы можете видеть, на каждом шаге мы всё ближе и ближе к нулю, зато каждый шаг становится меньшим, так как при приближении к нулю наклон становится меньше.

Теперь давайте попробуем реализовать это в коде и посмотрим, сможем ли мы полностью дойти до нуля. Импортируем библиотеку NumPy, установим значение w=20 и будем печатать результат на каждом шаге итераций. Количество приближений установим равным 30.

import numpy as np

w = 20

for i in xrange(30):

w = w – 0.1*2*w

print w

Как видим, w достигает значения 0,02, так что, похоже, 30 приближений недостаточно. Попробуем 100 приближений.

w = 20

for i in xrange(100):

w = w – 0.1*2*w

print w

Теперь результат равен 4,07*10^-9 – это очень близко к нулю.

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

Почему этот метод так важен? По мере продвижения далее в глубокое обучение и машинное обучение функции будут становиться всё более сложными. Для нейронных сетей с softmax нахождение производной может занять у вас несколько часов или даже дней. При переходе к свёрточным и рекуррентным нейронным сетям градиенты, конечно, можно найти на бумаге, но нет никакого желания тратить на это время. Куда лучше потратить время на проверку различных архитектур и параметров, не заботясь о градиентах. Тем более, что для вычисления градиентов можно использовать специальные библиотеки, такие как Theano и TensorFlow.

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

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

J(w₁, w₂) = w₁² + w₂⁴.

Введение в Theano Scan

Для начала нам нужно познакомимся с очень важной функцией scan из библиотеки Theano.

Чем она так важна? Задумаемся о том, как работает библиотека Theano. Мы создаём переменные и функционально их связываем, хотя они не имеют никаких значений, пока мы не запустим функцию. Поэтому когда мы создаём матрицу X, мы не указываем её размерность, а просто сообщаем, что это матрица, и Theano понимает, что это объект с двумя размерностями. Что будет, если мы захотим просмотреть все элементы X? Theano «не знает», где конец, поскольку мы ещё не сообщили ей этого. X всегда должен иметь определённую длину N; в таком случае мы можем использовать цикл for от 0 до N. В общем случае мы не можем точно знать длину наших обучающих последовательностей, поэтому если попробовать написать что-то вроде «for i in xrange(X.shape[0])», то ничего не получится, поскольку X.shape[0] пока не имеет значения. Вот тут-то и вступает в игру функция scan. Функция scan позволяет передавать длину конечной переменной как количество выполняемых циклов.

В некоторых случаях цикл for можно использовать, но всё равно лучше пользоваться функцией scan. Вот некоторые преимущества scan перед циклом for:

• во-первых, позволяет количеству проходов цикла быть частью символьного графа Theano;
• во-вторых, минимизирует количество обращений к графическому процессору, если тот задействован в вычислениях;
• в-третьих, она может вычислять градиенты при помощи последовательных этапов вычислений;
• в-четвёртых, она работает значительно быстрее циклa for для скомпилированной функции Theano;
• в-пятых, она может снизить общее количество используемой памяти, определяя необходимый её объём.

Итак, рассмотрим строение функции scan. Вот она в своей простейшей форме:

outputs, updates = theano.scan(

    fn=some_function,

    sequences=thing_to_loop_over,

    n_steps=number_of_times_to_iterate,

)

def some_function(element):

    return do_something_to(element)

Первый аргумент – это некоторая функция, которая будет применяться к каждому элементу последовательности. Второй аргумент – фактическая последовательность. Таким образом, к каждому отдельному элементу последовательности будет применена одна и та же функция fn. n_steps – это количество повторов; как правило, это просто длина последовательности.

Сама функция scan возвращает две вещи – outputs и updates. Updates – это объект Theano типа упорядоченных обновлений. Обычно он просто игнорируется, так что больше не будем о нём говорить. Действительно интересным является outputs. Предполагается, что некоторая функция в результате работы возвращает некоторое значение функции. Так вот, выходом в данном случае будут все значения, объединённые вместе. Например, если в качестве аргумента будут числа 1, 2, 3, а функцией будет квадрат числа, то выходом будут числа 1, 4, 9.

Разумеется, функция scan может делать вещи и посложнее. Включим ещё один аргумент outputs_info, позволяющий вычислять рекуррентные отношения:

outputs, updates = theano.scan(

    fn=reccurence,

    sequences=thing_to_loop_over,

    n_steps=number_of_times_to_iterate,

    outputs_info=[initial_value]

)

def reccurence(element, recurring_variable):

    return do_something_to(element, recurring_variable)

Аргумент outputs_info представляет начальные значения рекуррентных переменных. Обратите внимание, что они включены в список – это потому, что их может быть более одной, поэтому Theano ожидает здесь именно списка. В результате мы можем вычислить, например, последовательность Фибоначчи, в которой текущий член последовательности зависит от двух предыдущих. Заметьте, что при этом, когда рекурсия возвращает более одного элемента, исходящие переменные функции scan возвращают все значения в списке. На самом деле сколько аргументов вы укажете для рекурсии, столько и получите на выходе.

Ещё один аргумент, который можно добавить, – non_sequences. Он нужен для помещения данных в рекурсию, когда мы не хотим их перебирать, а получить в целом:

outputs, updates = theano.scan(

    fn=reccurence,

    sequences=thing_to_loop_over,

    n_steps=number_of_times_to_iterate,

    outputs_info=[initial_value]

    non_sequences=decay

)

def reccurence(x, last, decay):

    return (1-decay)*x + decay*last

Теперь попробуем реализовать это в коде и посмотрим, что получится. Сначала – файл scan1.py, который находится репозитарии github, на случай, если вам не хочется самим писать код.

Импортируем Numpy, Theano и Theano.tensor

import numpy as np

import theano

import theano.tensor as T

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

x = T.vector(‘x’)

Определим функцию square, возвращающую x².

def square(x):

  return x*x

И вызываем функцию scan.

outputs, updates = theano.scan(

  fn=square,

  sequences=x,

  n_steps=x.shape[0],

)

Теперь собственно функция Theano.

square_op = theano.function(

  inputs=[x],

  outputs=[outputs],

)

o_val = square_op(np.array([1, 2, 3, 4, 5]))

print(“output:”, o_val)

Запустим код и получим ожидаемый ответ.

Теперь перейдем к файлу scan2.py. Импорт библиотек остаётся прежним, но на этот раз мы попробуем вычислить последовательность Фибоначчи.

import numpy as np

import theano

import theano.tensor as T

N = T.iscalar(‘N’)

def recurrence(n, fn_1, fn_2):

  return fn_1 + fn_2, fn_1

outputs, updates = theano.scan(

  fn=recurrence,

  sequences=T.arange(N),

  n_steps=N,

  outputs_info=[1., 1.]

)

Создаём функцию Theano

fibonacci = theano.function(

  inputs=[N],

  outputs=outputs,

)

o_val = fibonacci(8)

print(“output:”, o_val)

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

Переходим к следующему примеру. Теперь будем работать с файлом scan3.py. Импортируем библиотеки:

import numpy as np

import matplotlib.pyplot as plt

import theano

import theano.tensor as T

В этом примере мы создадим фильтр нижних частот. Вначале создадим сигнал с большим количеством шума и длиной 300. Основой будет синусоида. Сразу и выведем его на экран, чтобы видеть сигнал с шумом.

X = 2*np.random.randn(300) + np.sin(np.linspace(0, 3*np.pi, 300))

plt.plot(X)

plt.title(“original”)

plt.show()

Определим переменную decay и вектор sequence. Это и будет наш зашумленный сигнал.

decay = T.scalar(‘decay’)

sequence = T.vector(‘sequence’)

Теперь функция reccurence.

def recurrence(x, last, decay):

  return (1-decay)*x + decay*last

В функции scan на этот раз просто пропустим updates:

outputs, _ = theano.scan(

  fn=recurrence,

  sequences=sequence,

  n_steps=sequence.shape[0],

  outputs_info=[np.float64(0)],

  non_sequences=[decay]

)

Теперь создадим функцию Theano с названием lpf.

lpf = theano.function(

  inputs=[sequence, decay],

  outputs=outputs,

)

Y = lpf(X, 0.99)

plt.plot(Y)

plt.title(“filtered”)

plt.show()

Запустим код и посмотрим, что выйдет.

Сначала мы видим первоначальный сигнал, весьма зашумленный. Далее мы видим очищенный сигнал.

Дискретная скрытая марковская модель в Theano

Как вы можете заметить, много кода взято из предыдущих файлов. Это связано с тем, что API остаётся без изменений. Так, у нас по-прежнему остаётся объект HMM, функции fit и get_cost_multi – тут ничего не меняется. Единственное, что претерпевает изменения – реализация fit.

Итак, начнём. Все необходимые импорты библиотек уже есть, поэтому приступим к написанию функции fit.

Кстати, у нас ещё будет функция set, позволяющая устанаваливать π, A и B в качестве переменных Theano.

Итак, функция fit. Мы добавим к ней несколько дополнительных параметров. В частности, V – размер словаря.

    def fit(self, X, learning_rate=0.001, max_iter=10, V=None, p_cost=1.0, print_period=10):

        if V is None:

            V = max(max(x) for x in X) + 1

        N = len(X)

        pi0 = np.ones(self.M) / self.M

        A0 = random_normalized(self.M, self.M)

        B0 = random_normalized(self.M, V)

        thx, cost = self.set(pi0, A0, B0)

Теперь функция set.

    def set(self, pi, A, B):

        self.pi = theano.shared(pi)

        self.A = theano.shared(A)

        self.B = theano.shared(B)

        thx = T.ivector(‘thx’)

Добавляем функцию reccurence. В качестве аргументов ей служат время, старое значение α и X. Мы сразу будем писать код для масштабированной скрытой марковской модели.

        def recurrence(t, old_a, x):

            a = old_a.dot(self.A) * self.B[:, x[t]]

            s = a.sum()

            return (a / s), s

Теперь функция scan библиотеки Theano, с которой вы только что познакомились. Updates, как и ранее, просто игнорируем.

        [alpha, scale], _ = theano.scan(

            fn=recurrence,

            sequences=T.arange(1, thx.shape[0]),

            outputs_info=[self.pi*self.B[:,thx[0]], None],

            n_steps=thx.shape[0]-1,

            non_sequences=thx

        )

        cost = -T.log(scale).sum()

        self.cost_op = theano.function(

            inputs=[thx],

            outputs=cost,

            allow_input_downcast=True,

        )

        return thx, cost

Теперь вернёмся к функции fit. Имея thx и cost, мы можем написать обновления:

        updates = [

            (self.pi, (self.pi – learning_rate*T.grad(cost, self.pi)) / self.pi.sum()),

            (self.A, (self.A – learning_rate*T.grad(cost, self.A)) / self.A.sum(axis=1).dimshuffle(0, ‘x’)),

            (self.B, (self.B – learning_rate*T.grad(cost, self.B)) / self.B.sum(axis=1).dimshuffle(0, ‘x’)),

        ]

Далее функция train_op. Это тоже функция Theano.

        train_op = theano.function(

            inputs=[thx],

            updates=updates,

            allow_input_downcast=True,

        )

Затем переходим к основному циклу.

        costs = []

        for it in xrange(max_iter):

            if it % print_period == 0:

                print(“it:”, it)

            for n in xrange(N):

                c = self.get_cost_multi(X, p_cost).sum()

                costs.append(c)

                train_op(X[n])

Выводим на экран значения π, A и B:

        print(“A:”, self.A.get_value())

        print(“B:”, self.B.get_value())

        print(“pi:”, self.pi.get_value())

        print(“len(costs):”, len(costs))

        plt.plot(costs)

        plt.show()

И теперь остальные функции, которые также надо определеить.

    def get_cost(self, x):

        return self.cost_op(x)

    def get_cost_multi(self, X, p_cost=1.0):

        P = np.random.random(len(X))

        return np.array([self.get_cost(x) for x, p in zip(X, P) if p < p_cost])

    def log_likelihood(self, x):

        return -self.cost_op(x)

Мы готовы запустить программу.

Сходимость вполне хорошая. Логарифм правдоподобия опять чуть лучше, чем при настоящих параметрах.