Урок 2 - основы оптимизации, линейные методы и одномерная оптимизация¶

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

Синтетически сгенерируем игрушечный датасет с продажами мороженого¶

In [ ]:
import numpy as np

samples_count = 20 # кол-во записей в датасете
X_variance = 6.5 # разброс по температуре
y_variance = 15 # разброс по кол-ву продаж
# явно укажем линейную зависимость, выразив её через наклон и смещение.
# наша идеальная модель должна иметь близкие коэффициенты.
k = 5 # угол наклона
b = 12  # смещение
get_noise = lambda s: (np.random.rand(s) - 0.5) * 2


X = np.linspace(5, 45, samples_count) + get_noise(samples_count) * y_variance
y = (X * k + b) + get_noise(samples_count) * y_variance

Посмотрим, как выглядят наши данные¶
In [ ]:
import matplotlib.pyplot as plt

plt.scatter(X, y)
Out[ ]:
<matplotlib.collections.PathCollection at 0x7e068c182950>
No description has been provided for this image

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

Вручную реализуем алгоритм линейной регресси для функции от одной переменной¶

Вначале напишем функцию ошибки, оценивающую, как сильно ошибается наша модель. В рамках решаемой задачи этого будет стандартная средне-квадратичная ошибка MSE.
$MSE=\frac{1}{N}\sum_{i=1}^{N}(y_{i}-ŷ_{i})^{2}$

In [ ]:
def MSE(y_true, y_pred):
  # ваш код здесь
  loss = None
  return loss

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

In [ ]:
def calculate_linear_regression(X, y):
  n = len(X)
  w = (np.sum(np.sum(y) * X) / n - np.sum(X*y)) / (np.sum(np.sum(X) * X) / n - np.sum(X**2))
  b = np.sum(y) / n - w * (np.sum(X) / n)
  return w, b
Посмотрим на результат¶
In [ ]:
w, b = calculate_linear_regression(X, y)
In [ ]:
# предскажем значения с помощью нашей модели

y_pred = w * X + b
In [ ]:
print("Ошибка модели по метрике MSE составила:", MSE(y, y_pred))
Ошибка модели по метрике MSE составила: 51.62708287611527

Обычно само по себе значение средне-квадратичной ошибки даёт мало информации о качестве работы модели. Но мы можем использовать другие метрики для оценки точности модели. Например, MAPE:

In [ ]:
def MAPE(y_true, y_pred):
  return (np.mean(np.abs(y_true - y_pred) / y_true)) * 100

метрика MAPE показывает процент различия предсказаний модели от целевых. Соответственно, чем меньше, тем лучше.

In [ ]:
print("Ошибка модели по метрике MAPE:", MAPE(y, y_pred))
Ошибка модели по метрике MAPE: 6.103235881570222

Если мы сделали вссё правильно, то значение MAPE должно быть не более 7%.

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

In [ ]:
plt.scatter(X, y)
plt.plot([min(X), max(X)], [min(X)*w+b, max(X)*w+b], color="red")
Out[ ]:
[<matplotlib.lines.Line2D at 0x79dcb0eb8760>]
No description has been provided for this image

Результат превосходный!

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

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


Попробуем смоделировать более сложный случай¶

In [ ]:
X = np.array([0, 1, 3, 4, 5, 6, 9, 10, 12, 13, 15, 20, 22, 23, 32, 41, 42, 43]) / 10
y = np.array([1.2, 1.01, 2.1, 1.5, 3.1, 2.8, 4, 4.1, 3, 3.2, 5, 7, 7.5, 5, 5, 4, 4, 5]) / 2

X_variance = 3 / 10
y_variance = 1 / 2
X = np.concatenate([get_noise(len(X))*X_variance+X, get_noise(len(X))*X_variance+X, X])
y = np.concatenate([get_noise(len(y))*y_variance+y, get_noise(len(y))*y_variance+y, y])
In [ ]:
plt.scatter(X, y)
Out[ ]:
<matplotlib.collections.PathCollection at 0x7e06508c2920>
No description has been provided for this image

Инициализируем и обучим полиномиальную модель¶

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

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

In [ ]:
# определим модель
def polynomial_regressor(x, c):
    return c[0] + c[1]*x + c[2]*x**2 + c[3]*x**3
In [ ]:
# Градиент функции ошибки по предсказаниям

def mse_loss_gradient(y_true, y_pred):
    # ващш код здесь

    loss = None
    return loss
In [ ]:
def polynomial_loss_gradient(X, y_true, c):
    y_pred = polynomial_regressor(X, c)
    gradient = np.zeros_like(c)
    for i in range(len(c)):
        gradient[i] = mse_loss_gradient(y_true, y_pred) * np.mean(X**i)
    return gradient

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

In [ ]:
learning_rate = 0.001
max_iterations = 3000

weights = np.random.random(4) # инициализируем наши веса

for iter_num in range(max_iterations):
    # ваш код здесь
iteration: 0 | loss: 27.195408124004448
iteration: 300 | loss: 0.543129138626655
iteration: 600 | loss: 0.543129138626655
iteration: 900 | loss: 0.543129138626655
iteration: 1200 | loss: 0.543129138626655
iteration: 1500 | loss: 0.543129138626655
iteration: 1800 | loss: 0.543129138626655
iteration: 2100 | loss: 0.543129138626655
iteration: 2400 | loss: 0.543129138626655
iteration: 2700 | loss: 0.543129138626655
In [ ]:
xs = np.linspace(0, 4.5, 30)
plt.plot(xs, polynomial_regressor(xs, weights), color="red")
plt.scatter(X, y)
plt.show()
No description has been provided for this image

А теперь попробуем обучить на этих данных модель линейной регрессии

In [ ]:
w, b = calculate_linear_regression(X, y)

plt.scatter(X, y)
plt.plot([min(X), max(X)], [min(X)*w+b, max(X)*w+b], color="red")
Out[ ]:
[<matplotlib.lines.Line2D at 0x7e065054fe20>]
No description has been provided for this image
In [ ]:
y_pred = w*X+b
print("MSE loss:", MSE(y, y_pred))
MSE loss: 0.5171860236690043

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

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

In [ ]:
X_new = np.array([4.5, 4.7, 4.8, 5])
y_new = np.array([1, 1.1, 0.8, 0.75])
In [ ]:
y_pred = w*X_new+b
print("Linear model MSE loss:", MSE(y_new, y_pred))

y_pred = polynomial_regressor(X_new, weights)
print("Polynomial model MSE loss:", MSE(y_new, y_pred))
Linear model MSE loss: 4.963922794003007
Polynomial model MSE loss: 0.5096049283662739

Разница линейной модели и полиномиальной колоссальна. Визуализируем предсказаниях обеих моделей:

In [ ]:
xs = np.linspace(0, 5, 30)
plt.plot(xs, polynomial_regressor(xs, weights), color="red")
plt.scatter(X, y)
plt.scatter(X_new, y_new, color="green")
Out[ ]:
<matplotlib.collections.PathCollection at 0x7e0650bd7f70>
No description has been provided for this image
In [ ]:
X_concated = np.concatenate([X, X_new])
y_concated = np.concatenate([y, y_new])

w, b = calculate_linear_regression(X_concated, y_concated)

plt.scatter(X_new, y_new, color="green")
plt.scatter(X, y)
plt.plot([min(X_concated), max(X_concated)], [min(X_concated)*w+b, max(X_concated)*w+b], color="red")
Out[ ]:
[<matplotlib.lines.Line2D at 0x7e06505a6020>]
No description has been provided for this image

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


Как можно было заметить, в рамках этого занятия рассматривалась только одномерная оптимизация. Но используя теорию вычисления частных производных, метод градиентного спуска, очевидно, можно применить и к многомерным функциям. Самый простой случай - линейная регрессия от нескольких переменных вида $w_1 * x_1 + w_2 * x_2 + ... + w_n * x_n + b$. Но понятно, что подобная теория также накладывается и на любые другие многомерные функции. По-сути, мы имеем только два необходимых условия - выпусклость и дифференцируемость (относительно функции ошибки).
В следующем занятии мы реализуем оптимизацию многомерных функций, уже записывая их в матричном виде и используя матричное дифференцирование. Но в качестве упражнение рекомендуется попробовать реализовать алгоритм оптимизации многомерной линейной регрессии или более сложной функции самостоятельно.

In [ ]:
 
In [ ]:
 
In [ ]: