Практика применения градиентного бустинга (Gradient Boosting)

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

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

Задание 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.

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

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

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



In [22]:

    
import numpy as np
import pandas as pd
import scipy
import sklearn
from matplotlib import pylab as plt
from sklearn import datasets
from sklearn import ensemble,cross_validation, datasets, metrics, tree, learning_curve, linear_model

%pylab inline









    



Populating the interactive namespace from numpy and matplotlib






    



WARNING: pylab import has clobbered these variables: ['clf', 'plt']
`%matplotlib` prevents importing * from pylab and numpy



In [3]:

    
data=datasets.load_boston()



In [4]:

    
X=np.array(data['data'])
y=np.array(data['target'])



In [5]:

    
print X.shape
print y.shape









    



(506L, 13L)
(506L,)



In [6]:

    
#Делим обучающую и тестовую выборки (примерно 25%)
X_train=X[:380]
X_test=X[380:]
y_train=y[:380]
y_test=y[380:]



In [7]:

    
X_train.shape









    Out[7]:





(380L, 13L)



In [8]:

    
n_trees = range(50)
base_algorithms_list=[tree.DecisionTreeRegressor(max_depth=5, random_state=42) for i in range(50)]
coefficients_list=[0.9 for i in range(50)]

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]

def write_answer_to_file(answer, filename):
    with open(filename, 'w') as f_out:
        f_out.write(str(answer))



In [9]:

    
base_algorithms_list  = []#список базовых алгоритмов
coefficients_list = []#список коэффициентов
mean_sq_err=[]#MSE - среднеквадратичная ошибка
gbm_predict_lst=[]
eta = 0.9
#начальный алгоритм
a_0 = tree.DecisionTreeRegressor(max_depth=5 , random_state=42)
a_0.fit(X_train, y_train)
base_algorithms_list.append(a_0)
coefficients_list.append(eta)
s=0
for k in range(49):
    s = y_train-gbm_predict(X_train)
    clf = tree.DecisionTreeRegressor(max_depth=5, random_state=42)
    clf.fit(X_train, s)
    base_algorithms_list.append(clf)
    coefficients_list.append(eta)
    m = (metrics.mean_squared_error(y_test, gbm_predict(X_test)))**0.5
    mean_sq_err.append(m)
print 'MSE на тестовой выборке:', m

#print 'MSE какие-были ',mean_sq_err









    



MSE на тестовой выборке: 5.45547207453



In [17]:

    
write_answer_to_file(m,'ans_week_4_2.txt')



In [14]:

    
base_algorithms_list  = []
coefficients_list = []
mean_sq_err=[]
gbm_predict_lst=[]
eta = 0.9/(1.0 + 0)
#начальный алгоритм
a_0 = tree.DecisionTreeRegressor(max_depth=5 , random_state=42)
a_0.fit(X_train, y_train)
base_algorithms_list.append(a_0)
coefficients_list.append(eta)
s=0
for k in range(49):
    s = y_train-gbm_predict(X_train)
    clf = tree.DecisionTreeRegressor(max_depth=5, random_state=42)
    clf.fit(X_train, s)
    base_algorithms_list.append(clf)
    coefficients_list.append((eta/(1.0 + (k+1))))
    m2 = (metrics.mean_squared_error(y_test, gbm_predict(X_test)))**0.5
    mean_sq_err.append(m2)
print 'MSE на тестовой выборке:', m2

#print 'MSE какие-были ',mean_sq_err









    



MSE на тестовой выборке: 5.24074258415



In [17]:

    
write_answer_to_file(m2,'ans_week_4_3.txt')



In [18]:

    
base_algorithms_list  = []
coefficients_list = []
mean_sq_err=[]
gbm_predict_lst=[]
eta = 0.9/(1.0 + 0)
#начальный алгоритм
a_0 = tree.DecisionTreeRegressor(max_depth=5 , random_state=42)
a_0.fit(X_train, y_train)
base_algorithms_list.append(a_0)
coefficients_list.append(eta)
s=0
for k in range(199):
    s = y_train-gbm_predict(X_train)
    clf = tree.DecisionTreeRegressor(max_depth=5, random_state=42)
    clf.fit(X_train, s)
    base_algorithms_list.append(clf)
    coefficients_list.append((eta/(1.0 + (k+1))))
    m2 = (metrics.mean_squared_error(y_test, gbm_predict(X_test)))**0.5
    mean_sq_err.append(m2)
print 'MSE на тестовой выборке:', m2









    



MSE на тестовой выборке: 5.29645701182

Согласно эксперименту выше с ростом числа деревьев (в данном случае композиция составила 200 деревьев в сравнении с 50 из предыдущего) видно, что RMSE начала расти (хоть и не существенно(на 0.05), но всё-таки).



In [19]:

    
base_algorithms_list  = []
coefficients_list = []
mean_sq_err=[]
gbm_predict_lst=[]
eta = 0.9/(1.0 + 0)
#начальный алгоритм
a_0 = tree.DecisionTreeRegressor(max_depth=15 , random_state=42)
a_0.fit(X_train, y_train)
base_algorithms_list.append(a_0)
coefficients_list.append(eta)
s=0
for k in range(49):
    s = y_train-gbm_predict(X_train)
    clf = tree.DecisionTreeRegressor(max_depth=15, random_state=42)
    clf.fit(X_train, s)
    base_algorithms_list.append(clf)
    coefficients_list.append((eta/(1.0 + (k+1))))
    m2 = (metrics.mean_squared_error(y_test, gbm_predict(X_test)))**0.5
    mean_sq_err.append(m2)
print 'MSE на тестовой выборке:', m2









    



MSE на тестовой выборке: 6.1315820954

Судя по данным выше увеличение глубины деревьев приводит к более серьезному ухудшению композиции, чем рост их числа, переобучение увеличивается гораздо заметнее (в алгоритме в глубиной 5 RMSE = 5.24074258415, с глубиной 15 - 6.1315820954).



In [20]:

    
write_answer_to_file('2 3','ans_week_4_4.txt')



In [23]:

    
linear_regressor = linear_model.LinearRegression()
linear_regressor.fit(X_train, y_train)
m_linear = (metrics.mean_squared_error(y_test, linear_regressor.predict(X_test)))**0.5
print 'RMSE на тестовой выборке с линейной регрессией:', m_linear









    



RMSE на тестовой выборке с линейной регрессией: 7.84812179648



In [24]:

    
write_answer_to_file(m_linear,'ans_week_4_5.txt')