Градиентный бустинг своими руками

Внимание: в тексте задания произошли изменения - поменялось число деревьев (теперь 50), правило изменения величины шага в задании 3 и добавился параметр random_state у решающего дерева. Правильные ответы не поменялись, но теперь их проще получить. Также исправлена опечатка в функции gbm_predict.

В этом задании будет использоваться датасет boston из sklearn.datasets. Оставьте последние 25% объектов для контроля качества, разделив X и y на X_train, y_train и X_test, y_test.

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


In [1]:
import sklearn.datasets as datasets
import numpy as np
from matplotlib import pyplot as plt
import pandas as pd

In [2]:
%matplotlib inline

In [3]:
data = datasets.load_boston()

In [4]:
target = data['target']
data = data['data']
data.shape


Out[4]:
(506, 13)

In [5]:
x_train, x_test, y_train, y_test = data[:380], data[380:], target[:380], target[380:]

Задание 1

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

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

Воспользуйтесь формулой из лекций, задающей ответы на обучающей выборке, на которые нужно обучать новый алгоритм (фактически это лишь чуть более подробно расписанный градиент от ошибки), и получите частный ее случай, если функция потерь L - квадрат отклонения ответа композиции a(x) от правильного ответа y на данном x.

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

Задание 2

Заведите массив для объектов DecisionTreeRegressor (будем их использовать в качестве базовых алгоритмов) и для вещественных чисел (это будут коэффициенты перед базовыми алгоритмами).

В цикле от обучите последовательно 50 решающих деревьев с параметрами max_depth=5 и random_state=42 (остальные параметры - по умолчанию). В бустинге зачастую используются сотни и тысячи деревьев, но мы ограничимся 50, чтобы алгоритм работал быстрее, и его было проще отлаживать (т.к. цель задания разобраться, как работает метод). Каждое дерево должно обучаться на одном и том же множестве объектов, но ответы, которые учится прогнозировать дерево, будут меняться в соответствие с полученным в задании 1 правилом.

Попробуйте для начала всегда брать коэффициент равным 0.9. Обычно оправдано выбирать коэффициент значительно меньшим - порядка 0.05 или 0.1, но т.к. в нашем учебном примере на стандартном датасете будет всего 50 деревьев, возьмем для начала шаг побольше.

В процессе реализации обучения вам потребуется функция, которая будет вычислять прогноз построенной на данный момент композиции деревьев на выборке X:

def gbm_predict(X):
    return [sum([coeff * algo.predict([x])[0] for algo, coeff in zip(base_algorithms_list, coefficients_list)]) for x in X]
(считаем, что base_algorithms_list - список с базовыми алгоритмами, coefficients_list - список с коэффициентами перед алгоритмами)

Эта же функция поможет вам получить прогноз на контрольной выборке и оценить качество работы вашего алгоритма с помощью mean_squared_error в sklearn.metrics.

Возведите результат в степень 0.5, чтобы получить RMSE. Полученное значение RMSEответ в пункте 2.


In [72]:
from sklearn.tree import DecisionTreeRegressor

In [107]:
num_of_estimators = 50
coef = [0.9 for _ in range(num_of_estimators)]

In [108]:
def gbm_predict(X, estimators):
    return np.asarray([sum([coeff * algo.predict([x])[0] for algo, coeff in zip(estimators, coef)]) for x in X], 
                      dtype=np.float32)
def mse(X, Y, estimators):
    error = np.sum(np.power(gbm_predict(X, estimators) - Y, 2)) / X.shape[0]
    return error
def rmse(X, Y, estimators):
    return np.sqrt(mse(X, Y, estimators))

In [109]:
def gradient_boost_fit(x, y):
    estimators = [DecisionTreeRegressor(max_depth=5, random_state=42) for _ in range(num_of_estimators)]
    
    gradient = y
    for i in range(num_of_estimators):
        estimators[i].fit(x, gradient)
        gradient = y - gbm_predict(x, estimators[:(i + 1)])
    return estimators

In [110]:
estimators = gradient_boost_fit(x_train, y_train)

In [111]:
rmse(x_train, y_train, estimators)


Out[111]:
0.0029837646972503428

In [112]:
rmse(x_test, y_test, estimators)


Out[112]:
5.4554720437013398

Задание 3

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

Попробуйте уменьшать вес перед каждым алгоритмом с каждой следующей итерацией по формуле 0.9 / (1.0 + i), где i - номер итерации (от 0 до 49). Используйте качество работы алгоритма как ответ в пункте 3.

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


In [113]:
coef = [0.9 / (1. + i) for i in range(num_of_estimators)]

In [114]:
estimators = gradient_boost_fit(x_train, y_train)

In [115]:
rmse(x_test, y_test, estimators)


Out[115]:
5.2356049450121418

Задание 4

Реализованный вами метод - градиентный бустинг над деревьями - очень популярен в машинном обучении. Он представлен как в самой библиотеке sklearn, так и в сторонней библиотеке XGBoost, которая имеет свой питоновский интерфейс. На практике XGBoost работает заметно лучше GradientBoostingRegressor из sklearn, но для этого задания вы можете использовать любую реализацию.

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

1. С увеличением числа деревьев, начиная с некоторого момента, качество работы градиентного бустинга не меняется существенно.

2. С увеличением числа деревьев, начиная с некоторого момента, градиентный бустинг начинает переобучаться.

3. С ростом глубины деревьев, начиная с некоторого момента, качество работы градиентного бустинга на тестовой выборке начинает ухудшаться.

4. С ростом глубины деревьев, начиная с некоторого момента, качество работы градиентного бустинга перестает существенно изменяться

In [30]:
2 3


  File "<ipython-input-30-847cf70929d6>", line 1
    2 3
      ^
SyntaxError: invalid syntax

Задание 5

Сравните получаемое с помощью градиентного бустинга качество с качеством работы линейной регрессии.

Для этого обучите LinearRegression из sklearn.linear_model (с параметрами по умолчанию) на обучающей выборке и оцените для прогнозов полученного алгоритма на тестовой выборке RMSE. Полученное качество - ответ в пункте 5.

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


In [ ]:
from sklearn.linear_model import LinearRegression
estimator = LinearRegression()
estimator.fit(x_train, y_train)

In [ ]:
def mse(X, Y, estimator):
    error = np.sum(np.power(estimator.predict(X) - Y, 2)) / X.shape[0]
    return error
def rmse(X, Y, estimator):
    return np.sqrt(mse(X, Y, estimator))
rmse(x_test, y_test, estimator)

In [ ]: