【机器学习笔记】——Bagging、Boosting、Stacking（RF / Adaboost / Boosting Tree / GBM / GBDT / XGBoost / LightGBM）

时间 2021-01-02 标签集成学习 Bagging Boosting GBDT XGBoost

集成学习

概念

首先，让我们先来了解一下，什么是集成学习（Ensemble Learning，写作银桑，读作昂桑）。集成学习由训练数据构建一组基分类器，然后通过对每个基分类器的预测进行投票来进行分类，对每个个体学习器的预测进行（加权）平均来进行回归，有时也被称为多分类器系统、基于委员会的学习等。如果把单个分类器比作一个决策者的话，集成学习的方法就相当于多个决策者共同进行一项决策。严格来说，集成学习并不算是一种分类器，而是一种分类器结合的方法。

集成学习将多个分类方法聚集在一起，以提高分类的准确率。通常一个集成分类器的分类性能会好于单个分类器

对于二分类问题，假设基分类器的错误率为 $\epsilon$ ：
$P(h_i(\mathcal{x} \ne \mathcal{y})) = \epsilon$

通过投票我们得到集成学习的分类预测：
$H(\mathcal{x}) = sign\left(\sum_{i = 1}^{T} h_i(\mathcal{x})\right)$

假设分类器的错误率相互独立，则由Hoeffding不等式有：
$P(H(\mathcal{x} \ne \mathcal{y})) = \sum_{k = 0}^{\lfloor T/2 \rfloor} \dbinom{T}{k}(1 - \epsilon)^k \epsilon^{T - k} \le exp\left(-\frac{1}{2}T(1-2\epsilon)^2 \right)$

即随着集成中基分类器数目 $T$ 增大，集成的错误率会指数级下降，最终趋于0。当然实际中个体学习器错误率独立的假设是不成立的

集成学习的很多理论研究都是针对弱分类器（精度略高于50%）进行的，但在实践中往往会使用较强的分类器
个体学习器要有一定的“准确性”，并且要有“多样性”。一般的，准确性很高之后，要增加多样性就需牺牲准确性。如何产生并结合“好而不同”的个体学习器，恰是集成学习研究的核心

下面这个例子，图(a)中每个学习器只有67%的精度，但集成学习却达到了100%；图(b)中三个学习器完全一样，集成之后性能没有提高；图©中每个学习器只有33%的精度集成学习的效果更差。

思维导图

Bagging算法

概念

Bagging算法（Bootstrap AGGregatING）是一种并行式集成算法。它基于自助采样法（Bootstrap sampling）采样出T个含有m个训练样本的采样集，然后基于每个采样集训练出一个基学习器，再将这些基学习器进行结合。

算法的伪代码如下：

从偏差-方差分解的角度看，降低一个估计的方差的方式是把多个估计平均起来，所以Bagging主要关注降低方差，因此它在不剪枝决策树、神经网络等易受样本扰动的学习器上效用更为明显。同时，Bagging并不能减小模型的偏差，所以应尽量选择偏差较小的基分类器，如未剪枝的决策树。

编程（分类）

参考：https://github.com/vsmolyakov/experiments_with_python/blob/master/chp01/ensemble_methods.ipynb

我们用iris数据集来看一下不同基分类器Bagging的效果

import itertools
import numpy as np

import seaborn as sns
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec

from sklearn import datasets

from sklearn.tree import DecisionTreeClassifier
from sklearn.neighbors import KNeighborsClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import RandomForestClassifier

from sklearn.ensemble import BaggingClassifier
from sklearn.model_selection import cross_val_score, train_test_split

from mlxtend.plotting import plot_learning_curves
from mlxtend.plotting import plot_decision_regions

np.random.seed(0)

导入数据

iris = datasets.load_iris()
X, y = iris.data[:, 0:2], iris.target

X.shape, y.shape

((150, 2), (150,))

设置基分类器和集成方法

# 配置基分类器的参数
clf1 = DecisionTreeClassifier(criterion='entropy', max_depth=1)    # 深度为1的决策树
clf2 = KNeighborsClassifier(n_neighbors=1)    # k=1的k近邻

# 每个Bagging有10个基分类器
bagging1 = BaggingClassifier(base_estimator=clf1, n_estimators=10, max_samples=0.8, max_features=0.8)
bagging2 = BaggingClassifier(base_estimator=clf2, n_estimators=10, max_samples=0.8, max_features=0.8)

查看单个学习器和Bagging后的效果

label = ['Decision Tree(depth = 1)', 'K-NN(k = 1)', 'Bagging Tree', 'Bagging K-NN']
clf_list = [clf1, clf2, bagging1, bagging2]

fig = plt.figure(figsize=(10, 8))
gs = gridspec.GridSpec(2, 2)    # 创建 2x2 的画布
grid = itertools.product([0,1],repeat=2)    # 返回参数的笛卡尔积的元组，每个元素可以用 2 次

for clf, label, grd in zip(clf_list, label, grid):        
    scores = cross_val_score(clf, X, y, cv=3, scoring='accuracy')    # Array of scores of the estimator for each run of the cross validation.
    print("Accuracy: %.2f (+/- %.2f) [%s]" %(scores.mean(), scores.std(), label))
        
    clf.fit(X, y)
    ax = plt.subplot(gs[grd[0], grd[1]])
    fig = plot_decision_regions(X=X, y=y, clf=clf, legend=2)    # 绘制决策边界
    plt.title(label)

plt.show()

Accuracy: 0.63 (+/- 0.02) [Decision Tree(depth = 1)]
Accuracy: 0.70 (+/- 0.02) [K-NN(k = 1)]
Accuracy: 0.66 (+/- 0.02) [Bagging Tree]
Accuracy: 0.61 (+/- 0.02) [Bagging K-NN]

可以发现决策树的集成有一定效果，而k近邻的效果反而变差了。这是因为k近邻是稳定学习器，不易受样本扰动影响。

再来看看数据集的划分对训练结果的影响

#plot learning curves
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)
    
plt.figure()
plot_learning_curves(X_train, y_train, X_test, y_test, bagging1, print_model=False, style='ggplot')
plt.show()

虽然每次结果有较大的差异，但是我们还是能得到一些结论：对于给定的划分比例，模型在训练集上的错分率基本在30%上下，而在测试集上，当训练集较小时，模型在测试集上表现很差，当训练集比例逐渐增加时，模型在测试集上的表现慢慢变好，一般在训练集占到80%左右时，模型在训练集和测试集上都能取得较好的表现

接下来我们看一下学习器的个数对bagging效果的影响

#Ensemble Size
num_est = np.concatenate(([1],np.arange(5,105,5)))
bg_clf_cv_mean = []
bg_clf_cv_std = []

for n_est in num_est:    
    bg_clf = BaggingClassifier(base_estimator=clf1, n_estimators=n_est, max_samples=0.8, max_features=0.8)
    scores = cross_val_score(bg_clf, X, y, cv=3, scoring='accuracy')
    bg_clf_cv_mean.append(scores.mean())
    bg_clf_cv_std.append(scores.std())
    

plt.figure()
# 绘制误差棒图，横轴为num_est，纵轴为scores.mean()，误差为scores.std()
(_, caps, _) = plt.errorbar(num_est, bg_clf_cv_mean, yerr=bg_clf_cv_std, c='blue', fmt='-o', capsize=5)
for cap in caps:
    cap.set_markeredgewidth(1)                                                           
plt.ylabel('Accuracy'); plt.xlabel('Ensemble Size'); plt.title('Bagging Tree Ensemble');
plt.show()

基于交叉验证的结果，我们可以看到在基分类器数量达到20后，集成的效果不再有明显的提升

随机森林

随机森林（Random Forest, RF）可以看成是改进的Bagging算法，它可以更好地解决决策树的过拟合问题。随机森林以CART树作为基分类器，对比决策树的Bagging，它有一个明显的特点——随机选择特征，对于普通的决策树，我们会在每个节点在当前的所有特征（假设有n个）中寻找最优特征，而随机森林的决策树会在每个节点在随机选取的特征集合（推荐 $n_{rf} = \log_2 n$ ）中寻找最优特征。这样做的好处是如果训练集中的几个特征对输出的结果有很强的预测性，那么这些特征会被每个决策树所应用，这样会导致树之间具有相关性，这样并不会减小模型的方差，但同时随机森林会不可避免的增加偏差。

RF的主要优点有：

训练可以高度并行化，对于大数据时代的大样本训练速度有优势。个人觉得这是的最主要的优点。
由于可以随机选择决策树节点划分特征，这样在样本特征维度很高的时候，仍然能高效的训练模型。
在训练后，可以给出各个特征对于输出的重要性
由于采用了随机采样，训练出的模型的方差小，泛化能力强。
相对于Boosting系列的Adaboost和GBDT， RF实现比较简单。
对部分特征缺失不敏感。

RF的主要缺点有：

在某些噪音比较大的样本集上，RF模型容易陷入过拟合。
取值划分比较多的特征容易对RF的决策产生更大的影响，从而影响拟合的模型的效果。

扩展

RF还有很多变种，其目的都是为了进一步降低方差。

Extremely randomized Trees

Extremely randomized Trees（以下简称Extra Trees）放弃了Boostrap sampling的方式，而选择使用原始数据集，而且在特征选择时更加暴力的随机选择一个特征直接作为最优特征

*Totally Random Trees Embedding

参考：http://www.cnblogs.com/pinard/p/6156009.html

Totally Random Trees Embedding(以下简称 TRTE)是一种非监督学习的数据转化方法。它将低维的数据集映射到高维，从而让映射到高维的数据更好的运用于分类回归模型。我们知道，在支持向量机中运用了核方法来将低维的数据集映射到高维，此处TRTE提供了另外一种方法。

TRTE在数据转化的过程也使用了类似于RF的方法，建立T个决策树来拟合数据。当决策树建立完毕以后，数据集里的每个数据在T个决策树中叶子节点的位置也定下来了。比如我们有3颗决策树，每个决策树有5个叶子节点，某个数据特征xx划分到第一个决策树的第2个叶子节点，第二个决策树的第3个叶子节点，第三个决策树的第5个叶子节点。则x映射后的特征编码为(0,1,0,0,0, 0,0,1,0,0, 0,0,0,0,1), 有15维的高维特征。这里特征维度之间加上空格是为了强调三颗决策树各自的子编码。

映射到高维特征后，可以继续使用监督学习的各种分类回归算法了。

*Isolation Forest

Isolation Forest（以下简称IForest）是一种异常点检测的方法。它也使用了类似于RF的方法来检测异常点。

对于在T个决策树的样本集，IForest也会对训练集进行随机采样,但是采样个数不需要和RF一样，对于RF，需要采样到采样集样本个数等于训练集个数。但是IForest不需要采样这么多，一般来说，采样个数要远远小于训练集个数？为什么呢？因为我们的目的是异常点检测，只需要部分的样本我们一般就可以将异常点区别出来了。

对于每一个决策树的建立， IForest采用随机选择一个划分特征，对划分特征随机选择一个划分阈值。这点也和RF不同。

另外，IForest一般会选择一个比较小的最大决策树深度max_depth,原因同样本采集，用少量的异常点检测一般不需要这么大规模的决策树。

对于异常点的判断，则是将测试样本点 $x$ 拟合到 $T$ 颗决策树。计算在每颗决策树上该样本的叶子节点的深度 $h_t(x)$ ，从而可以计算出平均高度 $h(x)$ 。此时我们用下面的公式计算样本点xx的异常概率:

$s(x,m)=2^{−\frac{h(x)}{c(m)}}$

其中，m为样本个数。 $c(m)$ 的表达式为：

$c(m)=2\ln (m−1)+\xi−2\frac{m−1}{m}$

$\xi$ 为欧拉常数， $s(x,m)$ 的取值范围是 $[0,1]$ ，取值越接近于1，则是异常点的概率也越大。

编程（分类）

import time

# DecisionTree
clf3 = DecisionTreeClassifier(criterion='entropy', max_depth=2)

# Tree Bagging
tree_bagging = BaggingClassifier(base_estimator=clf3, n_estimators=100)

# Random Forest
random_forest = RandomForestClassifier(n_estimators=100, criterion='entropy', max_depth = 2)

label = ['Decision Tree', 'Tree Bagging', 'Random Forest']
clf_list = [clf3, tree_bagging, random_forest]

fig = plt.figure(figsize=(15, 4))
gs = gridspec.GridSpec(1, 3)
grid = itertools.product([0], [0, 1, 2])

for clf, label, grd in zip(clf_list, label, grid):
    b = time.time()
    scores = cross_val_score(clf, X, y, cv=3, scoring='accuracy')
    runtime = time.time() - b
    print("Accuracy: %.2f (+/- %.2f) [%s] Runtime: %.2f s" %(scores.mean(), scores.std(), label, runtime))
        
    clf.fit(X, y)
    ax = plt.subplot(gs[grd[0], grd[1]])
    fig = plot_decision_regions(X=X, y=y, clf=clf, legend=2)
    plt.title(label)

plt.show()

Accuracy: 0.63 (+/- 0.04) [Decision Tree] Runtime: 0.01 s
Accuracy: 0.77 (+/- 0.03) [Tree Bagging] Runtime: 0.40 s
Accuracy: 0.76 (+/- 0.04) [Random Forest] Runtime: 0.36 s

从结果可以看到随机森林的运行速度比Tree Bagging要快一点，这是因为在选择最优特征时，随机森林并没有从全部特征中进行选择

rf_label = ['DecisionTree','1 Trees','5 Trees','10 Trees','20 Trees','100 Trees']

random_forest1 = RandomForestClassifier(n_estimators=1, criterion='entropy', max_depth = 2)
random_forest5 = RandomForestClassifier(n_estimators=5, criterion='entropy', max_depth = 2)
random_forest10 = RandomForestClassifier(n_estimators=10, criterion='entropy', max_depth = 2)
random_forest20 = RandomForestClassifier(n_estimators=20, criterion='entropy', max_depth = 2)
random_forest100 = RandomForestClassifier(n_estimators=100, criterion='entropy', max_depth = 2)

clf_list = [clf3, random_forest1, random_forest5, random_forest10, random_forest20, random_forest100]

fig = plt.figure(figsize=(15, 8))
gs = gridspec.GridSpec(2, 3)
grid = itertools.product([0, 1], [0, 1, 2])

for clf, label, grd in zip(clf_list, rf_label, grid):
    scores = cross_val_score(clf, X, y, cv=3, scoring='accuracy')
    print("Accuracy: %.2f (+/- %.2f) [%s]" %(scores.mean(), scores.std(), label))
        
    clf.fit(X, y)
    ax = plt.subplot(gs[grd[0], grd[1]])
    fig = plot_decision_regions(X=X, y=y, clf=clf, legend=2)
    plt.title(label)

plt.show()

Accuracy: 0.63 (+/- 0.04) [DecisionTree]
Accuracy: 0.59 (+/- 0.08) [1 Trees]
Accuracy: 0.69 (+/- 0.02) [5 Trees]
Accuracy: 0.71 (+/- 0.06) [10 Trees]
Accuracy: 0.74 (+/- 0.04) [20 Trees]
Accuracy: 0.79 (+/- 0.02) [100 Trees]

可以看到随着树的增加，随机森林的决策边界越来越光滑，精度越来越高。而且当树为1个时随机森林色精度是比决策树要差的，这很容易理解，有可能重要的特征被排除在外了，所以同样是一棵树，随机森林分类效果更差。我们再通过误差棒图看一下决策树数量对精度的影响

num_est = np.concatenate(([1],np.arange(5,105,5)))
rf_clf_cv_mean = []
rf_clf_cv_std = []

for n_est in num_est:    
    rf_clf = RandomForestClassifier(n_estimators=n_est, criterion='entropy', max_depth = 2)
    scores = cross_val_score(rf_clf, X, y, cv=3, scoring='accuracy')
    rf_clf_cv_mean.append(scores.mean())
    rf_clf_cv_std.append(scores.std())
    

plt.figure()
(_, caps, _) = plt.errorbar(num_est, rf_clf_cv_mean, yerr=rf_clf_cv_std, c='blue', fmt='-o', capsize=5)
for cap in caps:
    cap.set_markeredgewidth(1)                                                           
plt.ylabel('Accuracy')
plt.xlabel('Trees Number')
plt.show()

可以看到随着树的数量增加，随机森林一开始有明显的提升，但是当树到达10个后就基本趋于稳定了

下面我们比较一下随着学习器的增加，Tree Bagging和Random Forest的差异

num_est = np.concatenate(([1],np.arange(10,1005,10)))
bg_arr = []
rf_arr = []

for n_est in num_est:    
    rf_clf = RandomForestClassifier(n_estimators=n_est, criterion='entropy', max_depth = 2)
    rf_scores = cross_val_score(rf_clf, X, y, cv=3, scoring='accuracy')
    rf_arr.append(rf_scores.mean())
    
    bg_clf = BaggingClassifier(base_estimator=clf3, n_estimators=n_est, max_samples=0.8, max_features=0.8)
    bg_scores = cross_val_score(bg_clf, X, y, cv=3, scoring='accuracy')
    bg_arr.append(bg_scores.mean())
    

plt.figure()
plt.plot(num_est, bg_arr,  color='blue', label='Bagging')
plt.plot(num_est, rf_arr, color='red', label='Random Forest')
plt.legend()
plt.ylabel('Accuracy')
plt.xlabel('Learner Number')
plt.show()

可以看到Bagging和随机森林有相似的收敛性，但是随着个体学习器数量的增加，随机森林的泛化性能要更好。

为什么说Bagging通过减小方差来提升精度

因为Bagging是同质集成并且数据集也相同，因此各个学习器的偏差 $bias(h_t (x))$ 和方差 $Var(h_t (x))$ 近似相等，但是各个学习器并不是独立的，于是
$bias(H(x)) = bias \left(\frac{1}{T}\sum_{t = 1}^{T} h_t (x)\right) = bias(h_t (x))$

$Var(H(x)) = Var \left(\frac{1}{T}\sum_{t = 1}^{T} h_t (x)\right) = \frac{1}{T}Var(h_t (x))$

所以Bagging不能降低偏差，而是会降低方差（因为学习器之间的相关性，方差不会下降到 $\frac{1}{T}Var(h_t (x))$ ），如果各个学习器完全相同，即有 $\frac{1}{T}\sum_{t = 1}^{T} h_t (x) = h_t (x)$ ，那么方差也不会降低，这也就是为什么要求学习器之间要有“多样性”。随机森林通过随机选取特征的方式降低了各决策树的相关性，使得方差进一步降低。

Boosting

不同于bagging方法，boosting方法通过分步迭代（stage-wise）的方式来构建模型，在迭代的每一步构建的弱学习器都是为了弥补已有模型的不足。Boosting算法先从初始训练集训练出一个基学习器，再根据基学习器的表现对训练样本进行调整，使得先前基学习器做错的训练样本在后续受到更多关注，然后基于调整后的样本分布来训练下一个基学习器；如此重复进行，直至基学习器数目达到事先指定的值T，最终将这T个基学习器进行加权。从偏差—方差分解的角度看，Boosting主要关注降低偏差，因此Boosting能基于泛化性能相当弱的学习器构建出很强的集成。下图以Adaboost分类为例简要说明了Boosting的机制

Adaboost

Boosting算法最著名的代表是Adaboost（Adaptive Boosting）算法。Adaboost不改变所给的训练数据，而不断改变训练数据的权值分布，使得训练数据在基本分类器的学习中起不同的作用，然后利用基本分类器的线性组合构建最终的分类器。Adaboost的二类分类算法的描述如下图所示，其中 $y_i\in\{-1,+1\}$ ， $i = 1,2,...,m$ ， $f$ 是真实函数， $\mathcal{D}_t$ 是调整后用于进行第 $t$ 次训练的样本分布， $h_t$ 是基于分布 $\mathcal{D}_t$ 从数据集 $D$ 中训练出的分类器， $\epsilon_t$ 是 $h_t$ 误差的估计， $\alpha_t$ 是分类器 $h_t$ 的权重， $Z_t$ 是规范化因子，以确保 $\mathcal{D}_{t+1}$ 是一个分布， $t=1,2,...,T$ 。

Adaboost是一种比较有特点的算法，可以总结如下：

每次迭代改变的是样本的分布，而不是重复采样（re weight)。
样本分布的改变取决于样本是否被正确分类：总是分类正确的样本权值低，总是分类错误的样本权值高（通常是边界附近的样本）。
最终的结果是弱分类器的加权组合，权值表示该弱分类器的性能。

我们通过一个例子来感受一下Adaboost的机制

实例

我们有如下数据集，其中有红色圆形和绿色三角形共10个样本，希望建立一个分类器，能判断出任意位置样本属于哪一类，我们指定基学习器个数 $T$ 为3

我们初始化样本权值分布 $\mathcal{D}_1 = {0.1, 0.1, \cdots, 0.1}$ ，并基于该分布从训练集中训练出分类器 $h_1$ （理论上 $h_1$ 应使分类错误率达到最小，这里仅用一条线代表分类结果，并未给出这条线是怎么得到的，但是不代表这条线是随意做出的），之后估计 $h_1$ 的误差 $\epsilon_1$ ，确定 $h_1$ 的权重 $\alpha_1$ ，并更新样本的权值分布

重复上面操作依次训练出基学习器 $h_2$ 、 $h_3$ ，此时基学习器个数达到我们预先指定的值。可以按照基学习器的权重给出输出公式。并据此预测其他样本的分类

每个区域是属于哪个类别，由这个区域所在分类器的权值综合决定。比如左下角的区域的样本，属于绿色分类区的权重为 $h_1$ 中的0.42和 $h_2$ 中的0.65，其和为1.07；属于红色分类区域的权重为 $h_3$ 中的0.92；属于红色分类区的权重小于属于绿色分类区的权值，因此左下角属于绿

从结果中看，即使简单的分类器，组合起来也能获得很好的分类效果。

模型推导

Adaboost算法有多种推导方式，一种比较容易理解的是基于“加性模型”，即基学习器的线性组合

$H(x) = \sum_{t = 1}^{T} \alpha_t h_t(x)$

来最小化指数损失函数

$Loss(H | \mathcal{D}) = \mathbb{E}_{\mathcal{D}}[e^{-yH(x)}] = e^{-H(x)}P(y = 1 | x) + e^{H(x)}P(y = -1 | x)$

●什么是前向分步算法（forward stagewise algorithm）？

前向分步算法是求解加法模型损失函数（ $\sum_{i = 1}^{N}Loss(y_i, \alpha h(x_i))$ ）最小化的一种算法，简单来说就是在每一次学习一个基学习器，然后逐步逼近目标函数：

$\begin{array}{lcl} f_0(x) & = & 0 \\ f_1(x) & = & f_0(x) + \alpha_1 h_1(x) & = & \alpha_1 h_1(x) \\ f_2(x) & = & f_1(x) + \alpha_2 h_2(x) & = & \alpha_1 h_1(x) + \alpha_2 h_2(x) \\ \cdots \\ f_M(x) & = & f_{M-1}(x) + \alpha_M h_M(x) & = & \sum_{m = 1}^{M}\alpha_m h_m(x) \end{array}$

●为什么是指数损失函数？

指数损失函数对 $H(x)$ 求偏导得

$\frac{\partial{Loss(H | \mathcal{D})}}{\partial{H(x)}} = -e^{-H(x)}P(y = 1 | x) + e^{H(x)}P(y = -1 | x)$

令偏导数为0有

$H(x) = \frac{1}{2} \ln \frac{P(y = 1 |x)}{P(y = -1 |x)}$

于是有

$\begin{aligned} sign\left(H(x)\right) & = sign\left(\frac{1}{2} \ln \frac{P(y = 1 |x)}{P(y = -1 |x)}\right) \\ & = \begin{cases} 0& \text{if } P(y = 1 |x) \gt P(y = -1 |x)\\ 1& \text{if } P(y = 1 |x) \lt P(y = -1 |x) \end{cases} \\ & = \mathop{\arg\max}_{y \in \{+1,-1\}} P(y | x) \end{aligned}$

过程中忽略了 $P(y = 1 |x) = P(y = -1 |x)$ 的情况，另外西瓜书中对于概率的写法是 $P(f(x) = y|x)$ ，其中 $f$ 是真实函数。从结果可知， $sign\left(H(x)\right)$ 达到了贝叶斯最优错误率（ $if P(\omega_1|x) \gt P(\omega_1|x),\; then\ x\in \omega_1,\;else\ x\in \omega_2$ ），即最小化指数损失函数等价于分类错误率最小化。因此我们可以用它替代0/1损失函数作为最优目标。有了以上结论，我们便可以进行算法的推导。

●推导过程

在算法的第 $t$ 次循环中，我们基于分布 $\mathcal{D}_t$ 产生了基分类器 $h_t(x)$ 和分类器权重 $\alpha_t$ ，那么它们应该使得 $H_t(x) = \sum_{k = 1}^{t}\alpha_k h_k(x)$ 最小化损失函数

$\begin{aligned} Loss(H_t(x) | \mathcal{D}_t) & = Loss(\sum_{k = 1}^{t}\alpha_k h_k(x) \ |\ \mathcal{D}_k) \\ & = Loss(H_{t-1}(x) \ +\ \alpha_t h_t(x)\ |\ \mathcal{D}_k) \\ & = \mathbb{E}_{\mathcal{D}_t} e^{-y (H_{t-1}(x) \ +\ \alpha_t h_t(x))} \\ & = \sum_{i = 1}^{N} e^{-y_i H_{t-1}(x_i)} e^{-y_i\alpha_t h_t(x_i)} \end{aligned}$

我们令 $w_{t,i} = e^{-y_i H_{t-1}(x_i)}$ ，显然， $w_{t,i}$ 与 $\alpha_t$ 和 $h_t$ 无关，所以与最小化损失也无关，于是有

$Loss(H_t(x) | \mathcal{D}_t) = \sum_{i = 1}^{N} w_{t,i} e^{-y_i\alpha_t h_t(x_i)}$

设 $h_t^*$ 和 $\alpha_t^*$ 是满足要求的分类器和分类器权重，即

$(h_t^*,\alpha_t^*) = \mathop{\arg \min}_{\alpha_t, h_t}Loss(H_t(x) | \mathcal{D}_t) = \mathop{\arg \min}_{\alpha_t, h_t}\sum_{i = 1}^{N} w_{t,i} e^{-y_i\alpha_t h_t(x_i)}$

对任意 $\alpha_t \gt 0$ ，考虑

$h_t^*(x) = \mathop{\arg \min}_{h_t} \sum_{i = 1}^{N}w_{t,i}I(h_t(x_i) \ne y_i)$

这意味着 $h_t^*(x)$ 是使第t轮加权训练数据分类误差最小的基本分类器，即为所求。再求 $\alpha_t^*$

$Loss(H_t(x) | \mathcal{D}_t) = \sum_{i = 1}^{N} w_{t,i} e^{-y_i\alpha_t h_t^*(x_i)} = \sum_{h_t^*(x_i) \ne y_i}w_{t,i}e^{\alpha_t}+\sum_{h_t^*(x_i) = y_i}w_{t,i}e^{-\alpha_t}$

损失函数对 $\alpha_t$ 求偏导得

$\frac{\partial{Loss(H_t(x) | \mathcal{D}_t)}}{\partial{\alpha_t}} = -\sum_{h_t^*(x_i) \ne y_i}w_{t,i}e^{\alpha_t}+\sum_{h_t^*(x_i) = y_i}w_{t,i}e^{-\alpha_t}$

令偏导数为0有

$\alpha_t^* = \frac{1}{2}\ln \frac{1-\sum_{h_t^*(x_i) \ne y_i}w_{t,i}}{\sum_{h_t^*(x_i) \ne y_i}w_{t,i}} = \frac{1}{2}\ln \frac{1-\epsilon_t}{\epsilon_t}$

再来看样本权值的更新， $H_t(x_i) = H_{t-1}(x_i) + \alpha_t h_t(x_i)$ 两边同时乘 $-y_i$ 并进行指数运算得

$w_{t+1,i} = e^{-y_iH_t(x_i)} = e^{-y_iH_{t-1}(x_i)}e^{-y_i\alpha_t h_t(x_i)} = w_{t,i}e^{-y_i\alpha_t h_t(x_i)}$

之后为确保 $\mathcal{D}_{t+1}$ 是一个分布，再对权值进行规范化即可

为什么说Boosting通过减小偏差来提升精度

Boosting算法等价于用前向分步算法来最小化损失函数，以Adaboost为例，算法用这样的方式顺序地最小化指数损失函数 $Loss(y,H_t(x))$ ，偏差自然降低。但也正是因为这样造成了子模型之间的强相关性，因此不能显著降低方差。所以所Boosting主要通过降低偏差来提升精度

*多分类任务

参考：http://www.noobyard.com/article/p-kxwistfq-nz.html

目前比较常用的多分类方法有三种

adaboost M1方法

该方法的主要思路是利用多分类的基学习器，最终选择最有可能的分类结果

$H(x) = \mathop{\arg \max}_{y \in \mathcal{Y}} \sum_{t:h_t(x) = y}\alpha_t$

adaboost MH方法

应用“one-vs-all”的原理重构样本空间，最终选择最有可能的结果

对多分类输出进行二进制编码

类似adaboost MH方法，不过不需要重构样本空间，而是把Label映射到 $K$ 位的二进制码上（ $K$ 位分类数），然后训练 $K$ 个二分类分类器，在解码时生成 $K$ 位的二进制数。从而对应到一个label上。

编程（分类）

仍然使用Iris数据集，决策树算法作为个体学习器算法

import itertools
import numpy as np

import seaborn as sns
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec

from sklearn import datasets

from sklearn.tree import DecisionTreeClassifier

from sklearn.ensemble import BaggingClassifier
from sklearn.ensemble import RandomForestClassifier
from sklearn.ensemble import AdaBoostClassifier
from sklearn.model_selection import cross_val_score, train_test_split

from mlxtend.plotting import plot_learning_curves
from mlxtend.plotting import plot_decision_regions

iris = datasets.load_iris()
X, y = iris.data[:, 0:2], iris.target
    
clf = DecisionTreeClassifier(criterion='entropy', max_depth=2)

num_est = [1, 2, 3, 10]
label = ['AdaBoost (n_est=1)', 'AdaBoost (n_est=2)', 'AdaBoost (n_est=3)', 'AdaBoost (n_est=10)']

fig = plt.figure(figsize=(10, 8))
gs = gridspec.GridSpec(2, 2)
grid = itertools.product([0,1],repeat=2)

for n_est, label, grd in zip(num_est, label, grid):     
    boosting = AdaBoostClassifier(base_estimator=clf, n_estimators=n_est)  
    
    scores = cross_val_score(boosting, X, y, cv=3, scoring='accuracy')
    print("Accuracy: %.2f (+/- %.2f) [%s]" %(scores.mean(), scores.std(), label))
    
    boosting.fit(X, y)
    ax = plt.subplot(gs[grd[0], grd[1]])
    fig = plot_decision_regions(X=X, y=y, clf=boosting, legend=2)
    plt.title(label)

plt.show()

Accuracy: 0.63 (+/- 0.04) [AdaBoost (n_est=1)]
Accuracy: 0.56 (+/- 0.11) [AdaBoost (n_est=2)]
Accuracy: 0.73 (+/- 0.04) [AdaBoost (n_est=3)]
Accuracy: 0.71 (+/- 0.08) [AdaBoost (n_est=10)]

再来看看不同训练计划分下Adaboost的表现

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)

boosting = AdaBoostClassifier(base_estimator=clf, n_estimators=10)
        
plt.figure()
plot_learning_curves(X_train, y_train, X_test, y_test, boosting, print_model=False, style='ggplot')
plt.show()

下面是不同基学习器个数下的Adaboost的精度

#Ensemble Size
num_est = np.concatenate(([1],np.arange(5,105,5)))
ada_clf_cv_mean = []
ada_clf_cv_std = []

for n_est in num_est:
    ada_clf = AdaBoostClassifier(base_estimator=clf, n_estimators=n_est)
    scores = cross_val_score(ada_clf, X, y, cv=3, scoring='accuracy')
    ada_clf_cv_mean.append(scores.mean())
    ada_clf_cv_std.append(scores.std())
    


plt.figure()
(_, caps, _) = plt.errorbar(num_est, ada_clf_cv_mean, yerr=ada_clf_cv_std, c='blue', fmt='-o', capsize=5)
for cap in caps:
    cap.set_markeredgewidth(1)                                                                        
plt.ylabel('Accuracy'); plt.xlabel('Ensemble Size'); plt.title('AdaBoost Ensemble');
plt.show()

运行过程中我们可以明显的发现其和Bagging与RF的不同，就是一旦参数确定了结果就是一定的，没有随机性。这是因为Bagging采样的随机性和RF特征选择的随机性。我们再来比较一下这三个模型的精度随着基学习数目增加精度的变化。

num_est = np.concatenate(([1],np.arange(10,100,10),np.arange(100,1001,50)))
bg_arr = []
rf_arr = []
ada_arr = []

for n_est in num_est:    
    rf_clf = RandomForestClassifier(n_estimators=n_est, criterion='entropy', max_depth = 2)
    rf_scores = cross_val_score(rf_clf, X, y, cv=3, scoring='accuracy')
    rf_arr.append(rf_scores.mean())
    
    bg_clf = BaggingClassifier(base_estimator=clf, n_estimators=n_est, max_samples=0.8, max_features=0.8)
    bg_scores = cross_val_score(bg_clf, X, y, cv=3, scoring='accuracy')
    bg_arr.append(bg_scores.mean())
    
    ada_clf = AdaBoostClassifier(base_estimator=clf, n_estimators=n_est)
    scores = cross_val_score(ada_clf, X, y, cv=3, scoring='accuracy')
    ada_arr.append(scores.mean())
    

plt.figure()
plt.plot(num_est, bg_arr,  color='blue', label='Bagging')
plt.plot(num_est, rf_arr, color='red', label='Random Forest')
plt.plot(num_est, ada_arr, color='green', label='Adaboost')
plt.legend()
plt.ylabel('Accuracy')
plt.xlabel('Learner Number')
plt.show()

可见以决策树为基学习器的Adaboost模型表现并不是很好

*训练误差分析

Adaboost最基本的性质是它能在学习过程中不断减少训练误差，下面我们来证明这一性质。我们先令 $g(x) = \sum_{t = 1}{T} \alpha_t h_t(x)$ ，令 $w_{t,i}$ 表示第 $t$ 轮中样本 $i$ 的权值，显然有 $w_{1,i} = \frac{1}{N}$ ，则Adaboost最终分类器 $H(x)$ 的训练误差界为：

$bias(H(x)) \ =\ \frac{1}{N} \sum_{i = 1}^{N} I(H(x_i) \ne y_i) \ \le\ \frac{1}{N} \sum_{i = 1}^{N} e^{-y_i g(x_i)} \ =\ \prod_{t = 1}^{T} Z_t$

证明：

当 $H(x_i) \ne y_i$ 时， $y_i g(x_i) \lt 0$ ，所以 $e^{-y_i g(x_i)} \gt 1 = I(H(x_i) \ne y_i)$ ，于是不等式的前半部分得证，下面只需证明 $\frac{1}{N} \sum_{i = 1}^{N} e^{-y_i g(x_i)} \ =\ \prod_{t = 1}^{T} Z_t$ ：

$\begin{aligned} \frac{1}{N} \sum_{i = 1}^{N} e^{-y_i g(x_i)} & = \frac{1}{N} \sum_{i = 1}^{N} e^{-y_i \sum_{t = 1}^{T} \alpha_t h_t(x_i)} \\ & = w_{1,i} \sum_{i = 1}^{N} e^{-\sum_{t = 1}^{T} \alpha_t h_t(x_i)y_i } \\ & = \sum_{i = 1}^{N} w_{1,i} \prod_{t = 1}^{T} e^{- \alpha_t h_t(x_i)y_i } \\ & = Z_1 \sum_{i = 1}^{N} \frac{w_{1,i} \prod_{t = 1}^{T} e^{- \alpha_t h_t(x_i)y_i }}{Z_1} \\ & = Z_1 \sum_{i = 1}^{N} w_{2,i} \prod_{t = 2}^{T} e^{- \alpha_t h_t(x_i)y_i } \\ & = Z_1 Z_2 \sum_{i = 1}^{N} w_{3,i} \prod_{t = 3}^{T} e^{- \alpha_t h_t(x_i)y_i } \\ & = \cdots \\ & = Z_1 Z_2 \cdots Z_{t-1} \sum_{i = 1}^{N} w_{T,i} e^{- \alpha_T h_T(x_i)y_i } \\ & = \prod_{t = 1}^{T} Z_t \end{aligned}$ N1i=1∑Ne−yig(xi)=N1i=1∑Ne−yi∑t=1Tαtht(xi)=w1,ii=1∑Ne−∑t=1Tαtht(xi)yi=i=1∑Nw1,it=1∏Te−αtht(xi)yi=Z1i=1∑NZ1w1,i∏t=1Te−αtht(xi)yi=Z1i=1∑Nw2,it=2∏Te−αtht(xi)yi=Z1Z2i=1∑Nw3,it=3∏Te−αtht(xi)yi=⋯=Z1Z2⋯Zt−1i=1∑NwT,ie−αThT(xi)yi=t=1∏TZt

■

通过这个定理，我们可以在每一轮都适当修改 $H_t$ 使得 $Z_t$ 最小，从而使训练误差下降最快。特别的，对于二分类问题，我们有以下结论：

$\prod_{t =1}^{T} Z_t \ =\ \prod_{t = 1}^{T}[2\sqrt{\epsilon_t(1-\epsilon_t)} \ = \ \prod_{t = 1}^{T} \sqrt{(1 - 4\gamma_t^2)} \ \le \ exp\left(-2 \sum_{t = 1}^{T} \gamma_t^2 \right)$

其中， $\gamma_t = \frac{1}{2} - \epsilon_t \gt 0$

证明：

先证明等式的部分：

$\begin{aligned} Z_t & = \sum_{i=1}^{N} w_{t,i} exp(-\alpha_t y_i h_t(x_i)) \\ & = \sum_{y_i = h_t(x_i)}w_{t, i} e^{-\alpha_t} \ +\ \sum_{y_i \ne h_t(x_i)} w_{t, i} e^{\alpha_t} \\ & = e^{-\alpha_t}\sum_{y_i = h_t(x_i)}w_{t, i} \ +\ e^{\alpha_t}\sum_{y_i \ne h_t(x_i)} w_{t, i} \\ & = e^{-\alpha_t}(1 - \epsilon_t) \ +\ e^{\alpha_t}\epsilon_t \\ & = {\left(\frac{1-\epsilon_t}{\epsilon_t}\right)}^{-\frac{1}{2}}(1 - \epsilon_t) \ +\ {\left(\frac{1-\epsilon_t}{\epsilon_t}\right)}^{\frac{1}{2}}\epsilon_t \\ & = 2\sqrt{\epsilon_t(1-\epsilon_t)} \\ & = \prod_{t = 1}^{T} \sqrt{(1 - 4\gamma_t^2)} \end{aligned}$ Zt=i=1∑Nwt,iexp(−αtyiht(xi))=yi=ht(xi)∑wt,ie−αt + yi̸=ht(xi)∑wt,ieαt=e−αtyi=ht(xi)∑wt,i + eαtyi̸=ht(xi)∑wt,i=e−αt(1−ϵt) + eαtϵt=(ϵt1−ϵt)−21(1−ϵt) + (ϵt1−ϵt)21ϵt=2ϵt(1−ϵt) =t=1∏T(1−4γt2)

不等式部分可由 $\sqrt{1-x}$ 和 $e^x$ 在 $x = 0$ 处的泰勒展开式推导得出

■

进一步，若存在 $\gamma \gt 0$ 使得对所有 $t$ 有 $\gamma_t \ge \gamma$ ，那么有：

$bias(H(x)) \le exp\left( -2T\gamma^2 \right)$

这表明Adaboost在这种条件下的训练误差是以指数速率下降的

*过拟合分析

这里仅给出结论——Adaboost不会发生过拟合，这是Adaboost的另一个性质。“margin theory”可以比较直观地解释这一性质。

总结

●优点

adaboost是一种有很高精度的分类器。
可以使用各种方法构建子分类器，adaboost算法提供的是框架。
当使用简单分类器时，计算出的结果是可以理解的。而且弱分类器构造极其简单。
简单，不用做特征筛选。
不用担心overfitting！
随着迭代次数的增加，实际上错误率上界在下降。

●缺点

对异常值非常敏感
对损失函数形式有要求

●应用

用于二分类或多分类的应用场景
用于做分类任务的baseline
用于特征选择（feature selection)
用于修正badcase，

提升树

提升树（Boosting Tree）是以二叉分类树或二叉回归树为基本分类器的提升方法，采用加法模型（即基函数的线性组合）与前向分步算法。提升树有很多不同的算法，其主要区别在于使用的损失函数不同，包括用平方误差损失函数的回归问题，用指数函数损失的分类问题，以及用一般损失函数的一般决策问题。

提升树模型为：

$f_M(x) = \sum_{m = 1}^{M}T(x; \Theta_m)$

其中 $T(x; \Theta_m)$ 为决策树； $\Theta_m$ 为决策树的参数； $M$ 为树的个数。首先确定初始模型：

$f_0(x) = 0$

第 $m$ 步的模型是：

$f_m(x) = f_{m-1}(x) + T(x; \Theta_m),\quad m = 1, 2, \cdots, M$

通过最小化损失函数确定决策树 $T(x; \Theta_m)$ 的参数：

$\hat{\Theta}_m = \mathop{\arg \min}_{\Theta_m}Loss(y, f_m(x)) = \mathop{\arg \min}_{\Theta_m}Loss(y, f_{m-1}(x) + T(x; \Theta_m))$

以回归问题的提升树为例，其损失函数是平方误差损失

$\begin{aligned} Loss(y, f_m(x)) & = Loss(y, f_{m-1}(x) + T(x; \Theta_m)) \\ & = {(y - f_{m-1}(x) - T(x; \Theta_m))}^2 \\ & = {(r_m - T(x; \Theta_m))}^2 \end{aligned}$

其中， $r_m$ 是当前模型拟合数据的残差，因此我们只需根据 $(x, r_m)$ 来学习一个回归树 $T(x; \Theta_m)$ ，进而确定 $f_m(x) = f_{m-1}(x) + T(x; \Theta_m)$ ，最终学得 $f_M(x) = \sum_{m = 1}^{M}T(x; \Theta_m)$

实例

考虑用如下数据集学习一个提升树模型（深度为1）

$x_i$	1	2	3	4	5	6	7
$y_i$	4.5	4.7	4.9	5.3	5.8	7.0	7.9

解第1步求 $f_1(x)$ 即回归树 $T(x;\Theta_1)$

根据已知的数据，我们可以考虑如下切分点：

$1.5, \;2.5, \;3.5, \;4.5, \;5.5, \;6.5$

回忆二叉回归树的解法，容易求出各切分点 $s$ 下

$g(s) = \mathop{\min}_{c_1}\sum_{x\in R_1}{(y_i - c_1)}^2 + \mathop{\min}_{c_2}\sum_{x\in R_2}{(y_i - c_2)}^2$

其中 $R_1 = \{x \le s\}$ ， $R_2 = \{x \gt s\}$ ， $c_1$ 、 $c_2$ 易知为 $c_1 = \bar{y}I(x \in R_1)$ ， $c_2 = \bar{y}I(x \in R_2)$

当 $s = 1.5$ 时， $c_1 = 4.5$ ， $c_2 = (4.7+4.9+5.3+5.8+7+7.9)/6 = 5.93$ ， $g(s) = 0 + (4.7-5.93)^2+ (4.9-5.93)^2+ (5.3-5.93)^2+ (5.8-5.93)^2+ (7-5.93)^2+ (7.9-5.93)^2 = 8.0134$

当 $s = 2.5$ 时， $c_1 = (4.5+4.7)/2 = 4.6$ ， $c_2 = (4.9+5.3+5.8+7+7.9)/5 = 6.18$ ， $g(s) = (4.5-4.6)^2 + (4.7-4.6)^2+ (4.9-6.18)^2+ (5.3-6.18)^2+ (5.8-6.18)^2+ (7-6.18)^2+ (7.9-6.18)^2 = 6.208$

当 $s = 3.5$ 时， $c_1 = (4.5+4.7+4.9)/3 = 4.7$ ， $c_2 = (5.3+5.8+7+7.9)/4 = 6.5$ ， $g(s) = (4.5-4.7)^2 + (4.7-4.7)^2+ (4.9-4.7)^2+ (5.3-6.5)^2+ (5.8-6.5)^2+ (7-6.5)^2+ (7.9-6.5)^2 = 4.22$

当 $s = 4.5$ 时， $c_1 = (4.5+4.7+4.9+5.3)/4 = 4.85$ ， $c_2 = (5.8+7+7.9)/3 = 6.90$ ， $g(s) = (4.5-4.85)^2 + (4.7-4.85)^2+ (4.9-4.85)^2+ (5.3-4.85)^2+ (5.8-6.9)^2+ (7-6.9)^2+ (7.9-6.9)^2 = 2.57$

当 $s = 5.5$ 时， $c_1 = (4.5+4.7+4.9+5.3+5.8)/5 = 5.04$ ， $c_2 = (7+7.9)/2 = 7.45$ ， $g(s) = (4.5-5.04)^2 + (4.7-5.04)^2+ (4.9-5.04)^2+ (5.3-5.04)^2+ (5.8-5.04)^2+ (7-7.45)^2+ (7.9-7.45)^2 = 1.477$

当 $s = 6.5$ 时， $c_1 = (4.5+4.7+4.9+5.3+5.8+7)/6 = 5.37$ ， $c_2 = 7.9$ ， $g(s) = (4.5-5.37)^2 + (4.7-5.37)^2+ (4.9-5.37)^2+ (5.3-5.37)^2+ (5.8-5.37)^2+ (7-5.37)^2+ 0 = 4.2734$

$s$	1.5	2.5	3.5	4.5	5.5	6.5
$g(s)$	8.01	6.21	4.22	2.57	1.48	4.27

所以切分点 $s$ 应为5.5，即

$T_1(x) = \begin{cases} 5.04, & x \le 5.5 \\ 7.45, & x \gt 5.5 \end{cases}$

于是 $f_1(x) = T_1(x)$

用 $f_1(x)$ 拟合训练数据的平方损失为

$Loss(y, f_1(x)) = \sum_{i = 1}^{7}{(y_i - f_1(x_i))}^2 = 1.477$

下面用 $f_1(x)$ 拟合残差 $r_{2,i} = y_i - f_1(x_i)$ ：

$x_i$	1	2	3	4	5	6	7
$r_{2,i}$	-0.54	-0.34	-0.14	0.26	0.76	-0.45	0.45

接着求 $T_2(x)$ ，方法同求 $T_1(x)$ ，得

$s$	1.5	2.5	3.5	4.5	5.5	6.5
$g(s)$	1.14	0.92	0.87	1.14	1.48	1.24

所以切分点为3.5，即

$T_2(x) = \begin{cases} -0.34, & x \le 3.5 \\ 0.255, & x \gt 3.5 \end{cases}$

于是 $f_2(x) = f_1(x) + T_2(x)$ ，有

$f_2(x) = \begin{cases} 4.7, & x \le 3.5 \\ 5.295, & 3.5 \lt x \le 5.5 \\ 7.705, & x \gt 5.5 \end{cases}$

用 $f_2(x)$ 拟合训练数据的平方损失为

$Loss(y, f_2(x)) = \sum_{i = 1}^{7}{(y_i - f_2(x_i))}^2 = 0.870$

假设我们提前设定决策树的数量就是2，那么 $f_2(x)$ 即为所求，否则继续重复前面的运算直到决策树数量达到预先设定值 $M$ 或者 $f_m(x)$ 拟合训练数据的平方误差损失满足要求

梯度提升

梯度提升（Gradient boosting）是一种用于回归、分类和排序任务的机器学习技术，属于Boosting算法族的一部分。梯度提升通过集成多个弱学习器，通常是决策树，来构建最终的预测模型，在提升树的基础上，当损失函数不能使用便于计算的平方损失或者指数损失时，利用损失函数的负梯度

$-{\left[\frac{\partial{Loss(y, f(x_i))}}{\partial{f(x_i}}\right]}_{f(x) = f_{m-1}(x)}$

作为回归问题提升树算法中的残差的近似值，即

$r_{m,i} = -{\left[\frac{\partial{Loss(y, f(x_i))}}{\partial{f(x_i}}\right]}_{f(x) = f_{m-1}(x)}$

GBDT

参考：https://www.zybuluo.com/yxd/note/611571#gbdt算法

基于梯度提升算法的学习器叫做GBM(Gradient Boosting Machine)。理论上，GBM可以选择各种不同的学习算法作为基学习器。现实中，用得最多的基学习器是决策树。为什么梯度提升方法倾向于选择决策树（通常是CART树）作为基学习器呢？这与决策树算法自身的优点有很大的关系。决策树可以认为是if-then规则的集合，易于理解，可解释性强，预测速度快。同时，决策树算法相比于其他的算法需要更少的特征工程，比如可以不用做特征标准化，可以很好的处理字段缺失的数据，也可以不用关心特征间是否相互依赖等。决策树能够自动组合多个特征，学习特征之间更高级别的相互关系，它可以毫无压力地处理特征间的交互关系并且是非参数化的，因此你不必担心异常值或者数据是否线性可分（举个例子，决策树能轻松处理好类别A在某个特征维度x的末端，类别B在中间，然后类别A又出现在特征维度x前端的情况）。此外，模型很容易实现扩展。不过，单独使用决策树算法时，有容易过拟合缺点。所幸的是，通过各种方法，抑制决策树的复杂性，降低单颗决策树的拟合能力，再通过梯度提升的方法集成多个决策树，最终能够很好的解决过拟合的问题。由此可见，梯度提升方法和决策树学习算法可以互相取长补短，是一对完美的搭档。至于抑制单颗决策树的复杂度的方法有很多，比如限制树的最大深度、限制叶子节点的最少样本数量、限制节点分裂时的最少样本数量、吸收bagging的思想对训练样本采样（subsample），在学习单颗决策树时只使用一部分训练样本、借鉴随机森林的思路在学习单颗决策树时只采样一部分特征、在目标函数中添加正则项惩罚复杂的树结构（XGboost）等。现在主流的GBDT算法实现中这些方法基本上都有实现，因此GBDT算法的超参数还是比较多的，应用过程中需要精心调参，并用交叉验证的方法选择最佳参数。

具体的算法如下：

输入：训练数据集 $T = \{(x_1,y_1),(x_2,y_2),\cdots,(x_N,y_N)\}$ ， $x_i \in \mathcal{X}\subseteq \mathbb{R}^n$ ， $y_i \in \mathcal{Y}\subseteq \mathbb{R}$ ；损失函数 $Loss(y,f(x))$ ；

输出：回归树 $\hat{f}(x)$

(1) 初始化

$f_0(x) = \mathop{\arg \min}_{c} \sum_{i = 1}^{N} Loss(y_i, c)$

(2) For $m = 1, 2, \cdots, M$ do:

(a) 计算

$r_{m,i} = -{\left[\frac{\partial{Loss(y, f(x_i))}}{\partial{f(x_i}}\right]}_{f(x) = f_{m-1}(x)}$

(b) 对(x_i, r_{m,i})拟合一个回归树，得到第m棵树的叶节点区域 $R_{m,j}$ ， $j = 1, 2, \cdots, J$

$c_{m,j} = \mathop{\arg\min}_{c} \sum_{x_i \in R_{m,j}}Loss(y_i, f_{m-1}+c)$

(d) 更新 $f_m(x) = f_{m-1}(x)+\sum_{j = 1}^{J}c_{m,j}I(x \in R_{m,j})$

(3) 得到回归树

$\hat{f}(x) = f_M(x) = \sum_{m = 1}^{M}\sum_{j = 1}^{J}c_{m,j}I(x \in R_{m,j})$

GBDT实例

我们仍以前面的例子为例，损失函数也仍然选择平方误差损失，我们看看结果是否会有不同

解第一步先初始化 $f_0(x)$ ，易知使平方误差损失达到最小的 $c$ 值为 $c_0 = \frac{1}{7}\sum_{i = 1}^{7}y_i = 5.73$ ，于是 $f_0(x) = c_0 = 5.73$

接着拟合第一轮的残差

$r_{1,i} = -{\left[\frac{\partial{Loss(y_i, f(x_i))}}{\partial{f(x_i)}}\right]}_{f_0(x)} = 2{[y_i-f(x_i)]}_{f_0(x)}$

于是有

$x_i$	1	2	3	4	5	6	7
$r_{1,i}$	-2.45	-2.06	-1.66	-0.86	0.14	2.54	4.34

对 $r_{1,i}$ 拟合回归树，易求得各切分点下的 $g(s)$ 为

$s$	1.5	2.5	3.5	4.5	5.5	6.5
$g(s)$	32.05	24.80	16.87	10.26	5.89	17.06

于是确定切分点 $s = 5.5$ ，叶节点区域 $R_{1,1} = \{x | x \le 5.5\}$ ， $R_{1,2} = \{x | x \gt 5.5\}$ ，有回归树

$T_1(x) = \begin{cases} -1.38, & x \le 5.5 \\ 3.44, & x \gt 5.5 \end{cases}$

进而有

$c_{1,1} = \frac{1}{5}\sum_{i = 1}^{5}[y_i - f_0(x_i)] = \frac{1}{5}\sum_{i = 1}^{5}y_i - f_0(x) = -0.69$

$c_{1,2} = \frac{1}{2}\sum_{i = 6,7}[y_i - f_0(x_i)] = \frac{1}{5}\sum_{i = 6,7}y_i - f_0(x) = 1.72$

最后学得

$f_1(x) = f_0(x) + \sum_{j = 1}^{2}c_{1,j}I(x \in R_{1,j}) = \begin{cases} 5.04, & x \le 5.5 \\ 7.45, & x \gt 5.5 \end{cases}$

用 $f_1(x)$ 拟合训练数据的平方损失误差：

$Loss(y,f_1(x)) = \sum_{i = 1}^{7}{[y_i-f_1(x_i)]}^2 = 1.477$

至此我们可以看到虽然学习过程不同，但是学得的 $f_1(x)$ 和提升树是相同的，那么继续下一轮的学习，计算

$r_{2,i} = -{\left[\frac{\partial{Loss(y_i, f(x_i))}}{\partial{f(x_i)}}\right]}_{f_1(x)} = 2{[y_i-f(x_i)]}_{f_1(x)}$

于是有

$x_i$	1	2	3	4	5	6	7
$r_{2,i}$	-1.08	-0.68	-0.28	0.52	1.52	-0.9	0.9

对 $r_{2,i}$ 拟合回归树，易求得各切分点下的 $g(s)$ 为

$s$	1.5	2.5	3.5	4.5	5.5	6.5
$g(s)$	4.55	3.74	3.48	4.56	5.89	4.96

于是确定切分点 $s = 3.5$ ，叶节点区域 $R_{2,1} = \{x | x \le 3.5\}$ ， $R_{2,2} = \{x | x \gt 3.5\}$ ，有回归树

$T_1(x) = \begin{cases} -0.68, & x \le 3.5 \\ 0.51, & x \gt 3.5 \end{cases}$

进而有

$c_{2,1} = \frac{1}{3}\sum_{i = 1}^{3}[y_i - f_1(x_i)] = \frac{1}{3}\sum_{i = 1}^{3}y_i - f_1(x) = -0.34$

$c_{2,2} = \frac{1}{4}\sum_{i = 4}^{7}[y_i - f_1(x_i)] = \frac{1}{4}\left[ \sum_{i = 4}^{7}y_i - \sum_{i = 4}^{7}f_1(x_i) \right] = 0.255$

最后学得

$f_2(x) = f_1(x) + \sum_{j = 1}^{2}c_{2,j}I(x \in R_{1,j}) = \begin{cases} 4.70, & x \le 3.5 \\ 5.30, & 3.5 \lt x \le 5.5 \\ 7.70, & x \gt 5.5 \end{cases}$

用 $f_2(x)$ 拟合训练数据的平方损失误差：

$Loss(y,f_2(x)) = \sum_{i = 1}^{7}{[y_i-f_2(x_i)]}^2 = 0.870$

最终结果仍和前面是一样的，这说明提升树用平方误差损失实质上也是利用了损失函数的负梯度，提升树用残差进行学习本质上仍然是计算负梯度（因为平方误差损失的负梯度恰好是残差），如果使用其他损失函数就不对了。也就是说残差学习是负梯度学习的一个特例，或者说提升树是梯度提升树的一个特例，但是显然梯度提升适用范围更广，在实际中梯度提升可以使用任意可微的损失函数

GBDT为什么要用负梯度来近似残差？

首先我觉得这个问题本身就是错误的，我们并没有这样做过。GBDT的本来就是在每一轮迭代中对损失函数的负梯度构造回归树，可以证明，这样能够使学习器收敛为一个强学习器（损失函数最小）。而当损失函数为平方误差损失时，恰好损失函数的负梯度就是残差。同时对于损失函数为平方误差的例子，直接用残差来进行解释反而更容易理解，并且由于提出的更早所以变成了负梯度近似残差这一说法。就好像我们做一道数学选择题，我们可能用特殊值法很容易得到了答案，但是并没有接触这道题本质的解法，如果题干变化了（损失函数不再是平方误差损失）或者题型改为解答题，那么特殊值法也不再适用。

GBDT编程（回归）

参考：https://scikit-learn.org/stable/auto_examples/ensemble/plot_gradient_boosting_regression.html#sphx-glr-auto-examples-ensemble-plot-gradient-boosting-regression-py

用波士顿房价数据集做GBDT回归

import numpy as np
import matplotlib.pyplot as plt

from sklearn.ensemble import GradientBoostingRegressor
from sklearn import datasets
from sklearn.utils import shuffle
from sklearn.metrics import mean_squared_error

加载数据

boston = datasets.load_boston()
X, y = shuffle(boston.data, boston.target, random_state=13)
X = X.astype(np.float32)
offset = int(X.shape[0] * 0.9)
X_train, y_train = X[:offset], y[:offset]
X_test, y_test = X[offset:], y[offset:]

训练模型

params = {'n_estimators': 500, 'max_depth': 4, 'min_samples_split': 2,
          'learning_rate': 0.01, 'loss': 'ls'}
clf = GradientBoostingRegressor(**params)

clf.fit(X_train, y_train)
mse = mean_squared_error(y_test, clf.predict(X_test))
print("MSE: %.4f" % mse)

MSE: 6.5401

绘制偏差图

# 计算测试集偏差
test_score = np.zeros((params['n_estimators'],), dtype=np.float64)

for i, y_pred in enumerate(clf.staged_predict(X_test)):
    test_score[i] = clf.loss_(y_test, y_pred)

plt.figure(figsize=(12, 6))
plt.subplot(1, 2, 1)
plt.title('Deviance')
plt.plot(np.arange(params['n_estimators']) + 1, clf.train_score_, 'b-',
         label='Training Set Deviance')
plt.plot(np.arange(params['n_estimators']) + 1, test_score, 'r-',
         label='Test Set Deviance')
plt.legend(loc='upper right')
plt.xlabel('Boosting Iterations')
plt.ylabel('Deviance')

Text(0, 0.5, ‘Deviance’)

绘制特征重要性图

feature_importance = clf.feature_importances_
# make importances relative to max importance
feature_importance = 100.0 * (feature_importance / feature_importance.max())
sorted_idx = np.argsort(feature_importance)
pos = np.arange(sorted_idx.shape[0]) + .5

plt.barh(pos, feature_importance[sorted_idx], align='center')
plt.yticks(pos, boston.feature_names[sorted_idx])
plt.xlabel('Relative Importance')
plt.title('Variable Importance')
plt.show()

梯度提升 Vs Adaboost

参考：https://www.zybuluo.com/yxd/note/611571#gbdt算法

与AdaBoost算法不同，梯度提升方法在迭代的每一步构建一个能够沿着梯度最陡的方向降低损失（steepest-descent）的学习器来弥补已有模型的不足，而传统的boosting算法只能是尽量向梯度方向减小（在Adaboost实例中能够看到随着学习器数量增加精度有起伏的情况），这是GBDT与传统boosting算法最大的区别，这也是为什么GBDT相比传统boosting算法可以用更少的树个数与深度达到更好的效果。经典的AdaBoost算法只能处理采用指数损失函数的二分类学习任务等，而梯度提升方法可以使用更多种类的可微损失函数来处理各类学习任务（多分类、回归、排序等），应用范围大大扩展。另一方面，AdaBoost算法对异常点（outlier）比较敏感（作为错分值而受到较大关注），而梯度提升算法通过引入bagging思想、加入正则项等方法能够有效地抵御训练数据中的噪音，具有更好的健壮性。这也是为什么梯度提升算法（尤其是采用决策树作为弱学习器的GBDT算法）如此流行的原因。

XGboost

参考：https://arxiv.org/pdf/1603.02754v1.pdf

XGboost 的全称是Extreme Gradient Boosting，由华盛顿大学的陈天奇博士提出

●准备知识

在推导算法之前，先回忆一下目标函数（objective function）的概念，目标函数通常定义为如下形式：

$Obj(\Theta) = L(\Theta) + \Omega(\Theta)$

其中， $L(\Theta)$ 是损失函数，用来衡量模型拟合训练数据的好坏程度； $\Omega(\Theta)$ 称之为正则项，用来衡量学习到的模型的复杂度，按照奥卡姆剃刀的原则，因此模型越简单越好。常用的损失函数有平方损失（square loss）： $L(y,\hat{y}) = \sum_{i = 1}^{N}{(y_i - \hat{y}_i)}^2$ ；Logistic损失： $L(y,\hat{y}) = - \sum_{i = 1}^{N} [y_i \log sigmoid(\hat{y}_i) + (1 - y_i) \log sigmoid(1 - \hat{y}_i)]$ 。常用的正则项有L1范数： $\Omega(w) = \lambda ||w||_1$ 和L2范数： $\Omega(w) = \lambda ||w||^2$ 。Ridge regression就是指使用平方损失和L2范数正则项的线性回归模型；Lasso regression就是指使用平方损失和L1范数正则项的线性回归模型；逻辑回归（Logistic regression）指使用logistic损失和L2范数或L1范数正则项的线性模型。

再回顾一下泰勒公式的一些知识，根据泰勒公式把函数 $f(x + \triangle x)$ 在 $x$ 点处二阶展开，可得到如下等式：

$f(x + \triangle x) \approx f(x) + f'(x)\triangle x + \frac{1}{2}f''(x){\triangle x}^2$

那么对于一棵决策树，我们如何描述它的复杂度呢，假设 $f_m(x)$ 是一棵已经建好的决策树，它有 $T$ 个叶子结点， $w_j$ 是叶子结点的权重， $j = 1, 2, \cdots, T$ ，用 $q$ 表示树的结构，它把特征维度为 $d$ 的数据映射到叶子结点的索引： $q: \mathbb{R}^d \to \{1, 2, \cdots, T\}$ ，即 $q(x)$ 表示数据 $x$ 在树中的位置

于是 $f(x)$ 可以表示为：

$f_m(x) = w_{q(x)}, \quad, w \in \mathbb{R}^T, \quad, q: \mathbb{R}^d \to \{1, 2, \cdots, T\}$

有了公式化的决策树，我们便可以用树的叶子数量和叶子权重的平滑程度（用权重的L2范数表示）来定义树的复杂度

$\Omega(f_m(x)) = \gamma T + \frac{1}{2}\lambda \sum_{j = 1}^{T}w_j^2$

●模型

有了前面的准备，我们可以进入XGboost的模型推导，在前向分步算法中，第m步的目标函数为

$\begin{aligned} Obj(\Theta)^{(m)} & = L(\Theta) + \Omega(\Theta) \\ & = \sum_{i = 1}^{N}l(y_i,{\hat{y}_i}^{(m)}) + \sum_{j = 1}^{m} \Omega (f_j) \\ & = \sum_{i = 1}^{N}l(y_i,{\hat{y}_i}^{(m-1)} + f_m(x_i)) + \Omega (f_m) +Constant1 \end{aligned}$

其中 $Constant1 = \sum_{j = 1}^{m-1} \Omega (f_j)$ ，因其不参与第 $m$ 轮的计算，所以可视为常数。

由泰勒公式

$\begin{aligned} Obj(\Theta)^{(t)} & \approx \sum_{i = 1}^{N}[l(y_i,{\hat{y}_i}^{(m-1)}) + l'(y_i,{\hat{y}_i}^{(m-1)})f_m(x_i) + \frac{1}{2}l''(y_i,{\hat{y}_i}^{(m-1)})f_m^2(x_i)] + \Omega (f_m) +Constant1 \\ & \triangleq \sum_{i = 1}^{N}[l(y_i,{\hat{y}_i}^{(m-1)}) + g_i f_m(x_i) + \frac{1}{2}h_i f_m^2(x_i)] + \Omega (f_m) +Constant1 \\ & = \sum_{i = 1}^{N}[g_i f_m(x_i) + \frac{1}{2}h_i f_m^2(x_i)] + \Omega (f_m) + Constant1 + Constant2 \end{aligned}$

其中 $Constant2 = l(y_i,{\hat{y}_i}^{(m-1)})$ 同样对于目标函数求最优值无影响，所以目标函数可以进一步简化为

$\begin{aligned} Obj(\Theta)^{(m)} & \approx \sum_{i = 1}^{N}[g_i f_m(x_i) + \frac{1}{2}h_i f_m^2(x_i)] + \Omega (f_m) \\ & = \sum_{i = 1}^{N}[g_i w_{q(x_i)} + \frac{1}{2}h_i w_{q(x_i)}^2] + \gamma T + \frac{1}{2}\lambda \sum_{j = 1}^{T}w_j^2 \end{aligned}$

定义 $I_j = \{i | q(x_i) = j\}$ 表示第 $j$ 个叶子结点里的样本序号的集合，于是目标函数又可表示为

$\begin{aligned} Obj(\Theta)^{(m)} & = \sum_{j = 1}^{T}\left[\sum_{i \in I_j}\left(g_i w_j + \frac{1}{2}h_i w_j^2\right)\right] + \gamma T + \frac{1}{2}\lambda \sum_{j = 1}^{T}w_j^2 \\ & = \sum_{j = 1}^{T} \left[\left(\sum_{i \in I_j}g_i \right)w_j + \frac{1}{2}\left(\sum_{i \in I_j}h_i + \lambda\right)w_j^2\right] + \gamma T \\ & \triangleq \sum_{j = 1}^{T} [G_j w_j + \frac{1}{2}(H_j + \lambda)w_j^2] + \gamma T \end{aligned}$

假设我们的树的结构 $q(x)$ 已经确定，为求目标函数最小值，令目标函数对 $w_j$ 求偏导，令偏导数为0得

$w_j^* = - \frac{G_j}{H_j + \lambda}$

于是

$\min_{w}Obj(\Theta)^{(m)} = - \frac{1}{2}\sum_{j = 1}^{T} \frac{G_j^2}{H_j + \lambda} + \gamma T$

那么现在现在只要知道树的结构，就能得到一个该结构下的最好分数（ $Obj$ ）。可是树的结构应该怎么确定呢？没法用枚举，毕竟可能的状态基本属于无穷种。这种情况，贪婪算法是个好方法。从树的深度为0开始，每一节点都遍历所有的特征。对于某个特征，先按照该特征里的值进行排序（复杂度为 $O(n\log n)$ ），然后线性扫描（复杂度为 $O(n)$ ）该特征来决定最好的分割点，最后在所有特征（复杂度为 $O(d(n \log n +n)))$ ）里选择分割后，分数增加（Gain）最高的那个特征。

$Obj_{split} = − \frac{1}{2}[\frac{G^2_L}{H_L+\lambda}+\frac{G^2_R}{H_R+\lambda}]+\gamma T_{split}$

$Obj_{no-split} = − \frac{1}{2}\frac{{G_L + G_R}^2}{H_L+H_R +\lambda}+\gamma T_{no-split}$

$Gain=Obj_{no-split}−Obj_{split}=\frac{1}{2}[\frac{G^2_L}{H_L+\lambda}+\frac{G^2_R}{H_R+\lambda} - \frac{{G_L + G_R}^2}{H_L+H_R +\lambda}]−\gamma(T_{split}−T_{no-split})$

这时，就有两种后续。一种是当最好的分割的情况下， $Gain$ 为负时就停止树的生长（类似前剪枝），这样的话效率会比较高也简单，但是这样就放弃了未来可能会有更好的情况（因为引入了正则项，所以增加分支并不能保证 $Gain$ 为正）。另外一种就是一直分割到最大深度，然后进行修剪，递归得把划分叶子得到的 $Gain$ 为负的收回（类似后剪枝）。一般来说，后一种要好一些，于是我们采用后一种，完整的（无修剪）算法如下（若树的深度为 $k$ ，那么总的复杂度为 $O(kdn \log n)$ ）

此外还有近似算法，该算法更适合对数据量较大、难以计算的任务，算法描述如下，这里不做过多介绍

●注意

对某个节点的分割时，是需要按某特征的值排序，那么对于无序的类别变量，就需要进行one-hot化。不然的话，假设某特征有1，2，3三种变量。我们进行比较时，就会只比较左子树为1,2或者右子树为2,3，或者不分割，哪个更好。这样的话就没有考虑，左子树为1,3的分割。因为 $Gain$ 于特征的值范围是无关的，它采用的是已经生成的树的结构与权重来计算的。所以不需要对特征进行归一化处理。

●总述

我们回顾一下求解模型的整个过程来看看在具体应用XGboost时需要哪些操作

输出：回归树 $\hat{f}(x)$

(1) 初始化 $f_0(x) = 0$

(2) For $m = 1, 2, \cdots, M$ do:

a) 计算

$g_i = \frac{\partial{l(y_i,{\hat{y}_i}^{(m-1)})}}{\partial{{\hat{y}_i}^{(m-1)}}}$

$h_i = \frac{\partial^2{l(y_i,{\hat{y}_i}^{(m-1)})}}{{\left(\partial{\hat{y}_i}^{(m-1)}\right)}^2}$

b) 计算

$Obj = - \frac{1}{2}\sum_{j = 1}^{T} \frac{G_j^2}{H_j + \lambda} + \gamma T$

c) 采用贪婪算法生成树 $f_m(x)$

d) 更新 $\hat{y}^{m} = \hat{y}^{m-1} + \epsilon_m f_m(x)$ ，其中 $\epsilon_m$ 是一个收缩系数，或者叫做步进。加上它的好处就是每一步我们都不是做一个完全的最优化，留下余地给未来的循环，这样能防止过拟合。

(3) 得到回归树

$\hat{f}(x) = \sum_{m = 1}^{M}\epsilon_m f_m(x)$

XGboost编程

官网：https://xgboost.readthedocs.io/en/latest/python/python_intro.html

数据接口

xgboost可以把不同类型的数据转化为模块内置的DMatrix数据类型：xgboost.DMatrix(data, label=None, missing=None, weight=None, silent=False, feature_names=None, feature_types=None, nthread=None)，各参数的说明如下：

参数名	参数说明	备注
data	数据来源：可以是文本文件/csv文件/二进制文件/np.array/scipy.sparse/pd.DataFrame/dt.Frame	1. 注意对于分类特征，即非数值型特征要转化成one-hot编码。2. 对于有标题的csv文件，不能直接导入，要先转化为pdDataFrame格式
label	标签列
missing	缺失值	如果不设置，默认为np.nan
weight	样本权值	在排序任务中，权值是针对每个分组的，单个样本的权值没有意义
silent	是否打印信息
feature_names	特征名
feature_types	特征类型
nthread	使用线程数量，默认使用最大线程

Tips：把数据转化为二进制可以加快加载速度

dtrain = xgb.DMatrix('train.svm.txt')
dtrain.save_binary('train.buffer')

设置参数

可以通过字典或者列表的方式初始化参数

param = {'max_depth': 2, 'eta': 1, 'silent': 1, 'objective': 'binary:logistic'}
param['nthread'] = 4
param['eval_metric'] = 'auc'

通过字典的方法修改或者增加参数

param['eval_metric'] = ['auc', '[email protected]']

# alternatively:
# plst = param.items()
# plst += [('eval_metric', '[email protected]')]

设置验证集

evallist = [(dtest, 'eval'), (dtrain, 'train')]

很明显这里的param非常重要，那么它有哪些可以设置的参数呢

增加随机性

eta 这个就是学习步进，也就是上面中的 $\epsilon$ 。
subsample 这个就是随机森林的方式，每次不是取出全部样本，而是有放回地取出部分样本。有人把这个称为行抽取，subsample就表示抽取比例
colsample_bytree和colsample_bylevel 这个是模仿随机森林的方式，这是特征抽取。colsample_bytree是每次准备构造一棵新树时，选取部分特征来构造，colsample_bytree就是抽取比例。colsample_bylevel表示的是每次分割节点时，抽取特征的比例。
max_delta_step 这个是构造树时，允许得到 $f_m(x)$ 的最大值。如果为0，表示无限制。也是为了后续构造树留出空间，和eta相似

控制模型复杂度

max_depth 树的最大深度
min_child_weight 如果一个节点的权重和小于这玩意，那就不分了
gamma每次分开一个节点后，造成的最小下降的分数。类似于上面的 $Gain$
alpha和lambda就是目标函数里的表示模型复杂度中的L1范数和L2范数前面的系数

其他参数

booster 表示用哪种模型，一共有gbtree, gbline, dart三种选择。一般用gbtree。
nthread 并行线成数。如果不设置就是能采用的最大线程。
sketch_eps 这个就是近似算法里的 $\epsilon$ 。
scale_pos_weight 这个是针对二分类问题时，正负样例的数量差距过大。

训练模型

训练一个模型需要一个参数列表和一个数据集

bst = xgb.train(param, dtrain, num_round = 10, evallist)

保存模型

bst.save_model('0001.model')

模型及其特征映射也可以转储到文本文件中

# dump model
bst.dump_model('dump.raw.txt')
# dump model with feature map
bst.dump_model('dump.raw.txt', 'featmap.txt')

加载模型

bst = xgb.Booster({'nthread': 4})  # init model
bst.load_model('model.bin')  # load data

提前结束

如果你有一个验证集，你可以使用early_stopping_rounds来找到更合适的boost轮次。提前停止需要至少一组evals（估计指标）。如果有多个，它将使用最后一个。

train(..., evals=evals, early_stopping_rounds=10)

设置该参数后，训练会持续进行直到验证集的得分不再增加（评价标准是每early_stopping_rounds轮，验证集的错分率不再降低），当提前结束发生时，模型会额外输出三个字段：bst.best_score, bst.best_iteration 和 bst.best_ntree_limit，注意返回的模型是最后一次迭代的模型而不是最优的模型。

预测

# 7 entities, each contains 10 features
data = np.random.rand(7, 10)
dtest = xgb.DMatrix(data)
ypred = bst.predict(dtest)

如果生成模型的过程中出发了提前结束，那么在预测的时候可以设置bst.best_ntree_limit参数来使用最优的模型

ypred = bst.predict(dtest, ntree_limit=bst.best_ntree_limit)

可视化

特征重要性可视化，需要导入matplotlib模块

xgb.plot_importance(bst)

指定序号的树可视化，需要导入graphviz和matplotlib模块

xgb.plot_tree(bst, num_trees=2)

如果通过ipython运行，可以把指定的树转化为graphviz实例

xgb.to_graphviz(bst, num_trees=2)

实例

参考：https://github.com/dmlc/xgboost/tree/master/demo/guide-python

基础练习

# import numpy as np
import xgboost as xgb

# 加载数据
dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')

# 指定参数（树的深度为2，学习步进为1，打印信息，目标函数为二分类逻辑回归）
param = {'max_depth':2, 'eta':1, 'silent':1, 'objective':'binary:logistic'}

# 指定验证集并开始训练
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
num_round = 2
bst = xgb.train(param, dtrain, num_round, watchlist)

# 预测并保存模型
preds = bst.predict(dtest)
labels = dtest.get_label()
print('error=%f' % (sum(1 for i in range(len(preds)) if int(preds[i] > 0.5) != labels[i]) / float(len(preds))))
bst.save_model('model/2/0001.model')

# 将模型转储为文本文件
bst.dump_model('model/2/dump.raw.txt')
# 将模型和特征映射转储为文本文件
# bst.dump_model('model/2/dump.raw.txt', 'model/2/featmap.txt') # 这里无法生成featmap文件，我们一会再来解决

[11:02:22] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[11:02:22] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263
error=0.021726

import numpy as np

# 把DMatrix格式数据转化为二进制格式数据
dtest.save_binary('file/2/dtest.buffer')

# 加载训练好的模型和转化好的数据
bst2 = xgb.Booster(model_file='model/2/0001.model')
dtest2 = xgb.DMatrix('file/2/dtest.buffer')

preds2 = bst2.predict(dtest2)

# 查看两次预测结果是否相同
assert np.sum(np.abs(preds2 - preds)) == 0

[11:02:57] 1611x127 matrix with 35442 entries loaded from file/2/dtest.buffer

import pickle   # 可以序列化对象并保存到磁盘中，用于存储训练好的模型

# 序列化模型
pks = pickle.dumps(bst2)

# 加载序列化模型
bst3 = pickle.loads(pks)
preds3 = bst3.predict(dtest2)

assert np.sum(np.abs(preds3 - preds)) == 0

import scipy.sparse

# 通过CSR稀疏矩阵加载数据
print('start running example of build DMatrix from scipy.sparse CSR Matrix')
labels = []
row = []; col = []; dat = []
i = 0

for line in open('file/2/agaricus.txt.train'):
    arr = line.split()
    labels.append(int(arr[0]))
    for item in arr[1:]:
        k,v = item.split(':')
        row.append(i); col.append(int(k)); dat.append(float(v))
    i += 1
    
csr = scipy.sparse.csr_matrix((dat, (row, col)))    # 创建稀疏矩阵
dtrain = xgb.DMatrix(csr, label=labels)    # 加载稀疏矩阵
watchlist = [(dtest, 'eval'), (dtrain, 'train')]    # 指定训练集，测试集
bst = xgb.train(param, dtrain, num_round, watchlist)    # 训练模型


# 通过CSC稀疏矩阵加载数据
print('='*30, '\n', 'start running example of build DMatrix from scipy.sparse CSC Matrix')

csc = scipy.sparse.csc_matrix((dat, (row, col)))
dtrain = xgb.DMatrix(csc, label=labels)
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
bst = xgb.train(param, dtrain, num_round, watchlist)

start running example of build DMatrix from scipy.sparse CSR Matrix
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263
start running example of build DMatrix from scipy.sparse CSC Matrix
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263

# 通过numpy矩阵加载数据
print('start running example of build DMatrix from numpy array')

npymat = csr.todense()
dtrain = xgb.DMatrix(npymat, label=labels)
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
bst = xgb.train(param, dtrain, num_round, watchlist)

start running example of build DMatrix from numpy array
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263

定义损失函数和评估指标

注意：

自定义的评估函数需要两个返回值，分别是评估指标的名称和评估结果，并且指标名称中不能含有冒号和空格
注意模型给出的预测是margin，而在自定义损失函数有时候需要对预测进行处理，如用sigmoid函数

import numpy as np
import xgboost as xgb

# 自定义损失函数
print('start running example to used customized objective function')

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')

param = {'max_depth': 2, 'eta': 1, 'silent': 1}
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
num_round = 2

start running example to used customized objective function
[13:18:20] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[13:18:20] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test

我们以Logistic损失为例

$l(y_i, \hat{y}_i) = y_i\log (1 + e^{-\hat{y}_i}) +(1-y_i)\log (1 + e^{\hat{y}_i})$

于是

$\begin{aligned} g_i & = \frac{ \partial l(y_i,{\hat{y}_i}^{(m-1)}) } {\partial {\hat{y}_i}^{(m-1)} } \\ & = -y_i \frac{e^{-{\hat{y}_i}^{(m-1)}}} {1 + e^{-{\hat{y}_i}^{(m-1)}}} + (1 - y_i) \frac { e^{ {\hat{y}_i}^{(m-1)} } } {1 + e^{{\hat{y}_i}^{(m-1)}}} \\ & = -y_i \left( 1-\frac{1}{1 + e^{-{\hat{y}_i}^{(m-1)}}}\right) + (1 - y_i)\frac{1}{1 + e^{-{\hat{y}_i}^{(m-1)}}} \\ & = \frac{1}{1 + e^{-{\hat{y}_i}^{(m-1)}}} - y_i \\ & = Pred - Label \\ \\ h_i & = \frac{ \partial g_i } {\partial {\hat{y}_i}^{(m-1)} } \\ & = \frac{ e^{-{\hat{y}_i}^{(m-1)}} }{ {(1+e^{-{\hat{y}_i}^{(m-1)}})}^2 } \\ & = \frac{1}{1 + e^{-{\hat{y}_i}^{(m-1)}}} \left( \frac{1}{1 + e^{-{\hat{y}_i}^{(m-1)}}} -1 \right) \\ & = Pred * (Pred - 1) \end{aligned}$ gihi=∂y^i(m−1)∂l(yi,y^i(m−1))=−yi1+e−y^i(m−1)e−y^i(m−1)+(1−yi)1+ey^i(m−1)ey^i(m−1)=−yi(1−1+e−y^i(m−1)1)+(1−yi)1+e−y^i(m−1)1=1+e−y^i(m−1)1−yi=Pred−Label=∂y^i(m−1)∂gi=(1+e−y^i(m−1))2e−y^i(m−1)=1+e−y^i(m−1)1(1+e−y^i(m−1)1−1)=Pred∗(Pred−1)

# 自定义目标函数——对数似然损失，计算g_i,h_i
def logregobj(preds, dtrain):
    labels = dtrain.get_label()
    preds = 1.0 / (1.0 + np.exp(-preds))
    grad = preds - labels
    hess = preds * (1.0 - preds)
    return grad, hess


# 自定义评估函数
# 注意在自定义损失函数时，默认的预测值（Pred）是margin，如logistic回归的默认预测值是未经过logistic变换的值
# 所以在损失函数定义中要对Pred做一下siagmoid处理。而内置的评估函数是假定pred经过logistic变换的
def evalerror(preds, dtrain):
    labels = dtrain.get_label()
    # since preds are margin(before logistic transformation, cutoff at 0)
    return 'my-error', float(sum(labels != (preds > 0.0))) / len(labels)

bst = xgb.train(param, dtrain, num_round, watchlist, obj=logregobj, feval=evalerror)    # obj参数传入目标函数，feval传入评估函数

[0] eval-rmse:1.59229 train-rmse:1.59597 eval-my-error:0.042831 train-my-error:0.046522
[1] eval-rmse:2.40519 train-rmse:2.40977 eval-my-error:0.021726 train-my-error:0.022263

分步训练

注意

因为模型内置的评估函数是假定对预测进行处理的，因为它需要给我们直接的结果，而在分步训练中，我们要使用的上一轮未经处理的预测值，所以需要设置参数output_margin=True使我们获得上一轮未经处理的结果。以logistic损失为例，我们要的是上一轮的预测得分而不是预测标签

import numpy as np
import xgboost as xgb

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')
watchlist = [(dtest, 'eval'), (dtrain, 'train')]

print ('='*30, '\n', 'start running example to start from a initial prediction')
param = {'max_depth':2, 'eta':1, 'silent':1, 'objective':'binary:logistic'}

# 第一轮训练
bst = xgb.train(param, dtrain, 1, watchlist)

# 获取第一轮训练结果（margin），而不是预测标签
ptrain = bst.predict(dtrain, output_margin=True)
ptest = bst.predict(dtest, output_margin=True)

# 更新训练集和测试集的margin
dtrain.set_base_margin(ptrain)
dtest.set_base_margin(ptest)


print('='*30, '\n', 'this is result of running from initial prediction')
bst = xgb.train(param, dtrain, 1, watchlist)

[15:39:34] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[15:39:34] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test
start running example to start from a initial prediction
[0] eval-error:0.042831 train-error:0.046522
this is result of running from initial prediction
[0] eval-error:0.021726 train-error:0.022263

使用最初生成的n棵树进行预测

import numpy as np
import xgboost as xgb

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')
param = {'max_depth':2, 'eta':1, 'silent':1, 'objective':'binary:logistic'}
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
num_round = 3
bst = xgb.train(param, dtrain, num_round, watchlist)

print('='*30, '\n', 'start testing prediction from first n trees')

# 指定用1棵树进行预测
label = dtest.get_label()
ypred1 = bst.predict(dtest, ntree_limit=1)

# 指定用2棵树进行预测
ypred2 = bst.predict(dtest, ntree_limit=2)

# 默认情况下使用所有树进行预测
ypred3 = bst.predict(dtest)
print('error of ypred from first 1 tree =%f' % (np.sum((ypred1 > 0.5) != label) / float(len(label))))
print('error of ypred from first 2 trees =%f' % (np.sum((ypred2 > 0.5) != label) / float(len(label))))
print('error of ypred from all trees =%f' % (np.sum((ypred3 > 0.5) != label) / float(len(label))))

[15:44:55] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[15:44:55] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263
[2] eval-error:0.006207 train-error:0.007063
start testing prediction from first n trees
error of ypred from first 1 tree =0.042831
error of ypred from first 2 trees =0.021726
error of ypred from all trees =0.006207

使用广义线性模型作为基学习器

import xgboost as xgb

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')

# change booster to gblinear, so that we are fitting a linear model
# alpha is the L1 regularizer
# lambda is the L2 regularizer
# you can also set lambda_bias which is L2 regularizer on the bias term
param = {'silent':1, 'objective':'binary:logistic', 'booster':'gblinear',
         'alpha': 0.0001, 'lambda': 1}

# normally, you do not need to set eta (step_size)
# XGBoost uses a parallel coordinate descent algorithm (shotgun)（平行坐标下降算法）
# there could be affection on convergence（收敛） with parallelization on certain cases
# setting eta to be smaller value, e.g 0.5 can make the optimization more stable
# param['eta'] = 1

watchlist = [(dtest, 'eval'), (dtrain, 'train')]
num_round = 4
bst = xgb.train(param, dtrain, num_round, watchlist)
preds = bst.predict(dtest)
labels = dtest.get_label()
print('error=%f' % (sum(1 for i in range(len(preds)) if int(preds[i] > 0.5) != labels[i]) / float(len(preds))))

[15:56:57] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[15:56:57] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test
[0] eval-error:0.113594 train-error:0.104099
[1] eval-error:0.117939 train-error:0.105788
[2] eval-error:0.121043 train-error:0.107631
[3] eval-error:0.121043 train-error:0.107477
error=0.121043

交叉验证

import numpy as np
import xgboost as xgb

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
param = {'max_depth':2, 'eta':1, 'silent':1, 'objective':'binary:logistic'}
num_round = 2

print('='*30, '\n', 'running cross validation')
# 结果：[轮次] 指标名:指标均值+指标标准差
# 设置5折交叉验证
xgb.cv(param, dtrain, num_round, nfold=5,
       metrics={'error'}, seed=0,
       callbacks=[xgb.callback.print_evaluation(show_stdv=True)])

print('='*30, '\n', 'running cross validation, disable standard deviation display')
# num_boost_round是最大迭代次数，early_stopping_rounds是提前停止，若最近3轮内指标不再提高，停止并输出最好轮次，结果不显示指标的标准差
res = xgb.cv(param, dtrain, num_boost_round=10, nfold=5,
             metrics={'error'}, seed=0,
             callbacks=[xgb.callback.print_evaluation(show_stdv=False),
                        xgb.callback.early_stop(3)])

print(res)


print('='*30, '\n', 'running cross validation, with preprocessing function')
# 定义预处理函数对训练集、测试集和参数进行处理，可以用来重塑权重
def fpreproc(dtrain, dtest, param):
    # 重塑权重
    label = dtrain.get_label()
    ratio = float(np.sum(label == 0)) / np.sum(label == 1)
    param['scale_pos_weight'] = ratio
    return (dtrain, dtest, param)

# 进行每折验证前先对数据和参数进行预处理
res = xgb.cv(param, dtrain, num_round, nfold=5,
             metrics={'auc'}, seed=0, fpreproc=fpreproc)

print(res)


# 用自定义的损失函数和评估函数进行交叉验证
print('='*30, '\n', 'running cross validation, with cutomsized loss function')
def logregobj(preds, dtrain):
    labels = dtrain.get_label()
    preds = 1.0 / (1.0 + np.exp(-preds))
    grad = preds - labels
    hess = preds * (1.0 - preds)
    return grad, hess
def evalerror(preds, dtrain):
    labels = dtrain.get_label()
    return 'error', float(sum(labels != (preds > 0.0))) / len(labels)

param = {'max_depth':2, 'eta':1, 'silent':1}

res = xgb.cv(param, dtrain, num_round, nfold=5, seed=0,
             obj=logregobj, feval=evalerror)

print(res)

[16:14:25] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train

==============================
running cross validation
[0] train-error:0.0506682+0.009201 test-error:0.0557316+0.0158887
[1] train-error:0.0213034+0.00205561 test-error:0.0211884+0.00365323

==============================
running cross validation, disable standard deviation display
[0] train-error:0.0506682 test-error:0.0557316
Multiple eval metrics have been passed: ‘test-error’ will be used for early stopping.

Will train until test-error hasn’t improved in 3 rounds.
[1] train-error:0.0213034 test-error:0.0211884
[2] train-error:0.0099418 test-error:0.0099786
[3] train-error:0.0141256 test-error:0.0144336
[4] train-error:0.0059878 test-error:0.0062948
[5] train-error:0.0020344 test-error:0.0016886
[6] train-error:0.0012284 test-error:0.001228
[7] train-error:0.0012284 test-error:0.001228
[8] train-error:0.0009212 test-error:0.001228
[9] train-error:0.0006142 test-error:0.001228
Stopping. Best iteration:
[6] train-error:0.0012284+0.000260265 test-error:0.001228+0.00104094

train-error-mean train-error-std test-error-mean test-error-std
0 0.050668 0.009201 0.055732 0.015889
1 0.021303 0.002056 0.021188 0.003653
2 0.009942 0.006076 0.009979 0.004828
3 0.014126 0.001706 0.014434 0.003517
4 0.005988 0.001878 0.006295 0.003123
5 0.002034 0.001470 0.001689 0.000574
6 0.001228 0.000260 0.001228 0.001041

==============================
running cross validation, with preprocessing function
train-auc-mean train-auc-std test-auc-mean test-auc-std
0 0.958228 0.001442 0.958232 0.005778
1 0.981414 0.000647 0.981431 0.002595

==============================
running cross validation, with cutomsized loss function
train-error-mean train-error-std train-rmse-mean train-rmse-std
0 0.050668 0.009201 1.595072 0.003868
1 0.021303 0.002056 2.442600 0.076834

test-error-mean test-error-std test-rmse-mean test-rmse-std
0 0.055732 0.015889 1.598043 0.012826
1 0.021188 0.003653 2.449282 0.080900

预测叶子索引

import xgboost as xgb

dtrain = xgb.DMatrix('file/2/agaricus.txt.train')
dtest = xgb.DMatrix('file/2/agaricus.txt.test')
param = {'max_depth':2, 'eta':1, 'silent':1, 'objective':'binary:logistic'}
watchlist = [(dtest, 'eval'), (dtrain, 'train')]
num_round = 3
bst = xgb.train(param, dtrain, num_round, watchlist)

print ('start testing predict the leaf indices')

# 使用最初的2棵树
print('='*30, '\n', 'using first 2 trees')
leafindex = bst.predict(dtest, ntree_limit=2, pred_leaf=True)
print(leafindex.shape)
print(leafindex)

# 使用所有树
print('='*30, '\n', 'using all tree')
leafindex = bst.predict(dtest, pred_leaf=True)
print(leafindex.shape)
print(leafindex)

[16:49:55] 6513x127 matrix with 143286 entries loaded from file/2/agaricus.txt.train
[16:49:55] 1611x127 matrix with 35442 entries loaded from file/2/agaricus.txt.test
[0] eval-error:0.042831 train-error:0.046522
[1] eval-error:0.021726 train-error:0.022263
[2] eval-error:0.006207 train-error:0.007063
start testing predict the leaf indices
using first 2 trees
(1611, 2)
[[4 3]
[3 3]
[4 3]
…
[3 3]
[5 4]
[3 3]]
using all tree
(1611, 3)
[[4 3 5]
[3 3 5]
[4 3 5]
…
[3 3 3]
[5 4 5]
[3 3 3]]

结果显示1611个测试样本在每棵树的叶子索引

利用sklearn中的数据学习xgboost

iris：鸢尾花数据集：多分类（3）数据集
digits：手写数字数据集：多分类（10）数据集
boston：波士顿房价数据集：回归任务数据集

import pickle
import xgboost as xgb

import numpy as np
from sklearn.model_selection import KFold, train_test_split, GridSearchCV
from sklearn.metrics import confusion_matrix, mean_squared_error
from sklearn.datasets import load_iris, load_digits, load_boston

rng = np.random.RandomState(31337)    # 用于产生为随机数


# 使用手写数字数据集中标签为0和1的数据做二分类任务
print("Zeros and Ones from the Digits dataset: binary classification")
digits = load_digits(2)
y = digits['target']
X = digits['data']
kf = KFold(n_splits=2, shuffle=True, random_state=rng)    # 交叉验证
for train_index, test_index in kf.split(X):
    xgb_model = xgb.XGBClassifier().fit(X[train_index], y[train_index])
    predictions = xgb_model.predict(X[test_index])
    actuals = y[test_index]
    print(confusion_matrix(actuals, predictions))

    
# 使用鸢尾花数据集做多分类任务
print('='*30, '\n', "Iris: multiclass classification")
iris = load_iris()
y = iris['target']
X = iris['data']
kf = KFold(n_splits=2, shuffle=True, random_state=rng)
for train_index, test_index in kf.split(X):
    xgb_model = xgb.XGBClassifier().fit(X[train_index], y[train_index])
    predictions = xgb_model.predict(X[test_index])
    actuals = y[test_index]
    print(confusion_matrix(actuals, predictions))


# 使用波士顿房价数据集做回归任务
print('='*30, '\n', "Boston Housing: regression")
boston = load_boston()
y = boston['target']
X = boston['data']
kf = KFold(n_splits=2, shuffle=True, random_state=rng)
for train_index, test_index in kf.split(X):
    xgb_model = xgb.XGBRegressor().fit(X[train_index], y[train_index])
    predictions = xgb_model.predict(X[test_index])
    actuals = y[test_index]
    print(mean_squared_error(actuals, predictions))

    
# 参数优化
print('='*30, '\n', "Parameter optimization")
y = boston['target']
X = boston['data']
xgb_model = xgb.XGBRegressor()
clf = GridSearchCV(xgb_model,
                   {'max_depth': [2,4,6],
                    'n_estimators': [50,100,200]}, verbose=1, cv=5)    # 5折交叉验证，需要对max_depth和n_estimators进行优化，偶尔输出日志
clf.fit(X,y)
print(clf.best_score_)
print(clf.best_params_)


# 对sklearn模型进行序列化
print('='*30, '\n', "Pickling sklearn API models")
# 必须用二进制格式打开
pickle.dump(clf, open("best_boston.pkl", "wb"))
clf2 = pickle.load(open("best_boston.pkl", "rb"))
print(np.allclose(clf.predict(X), clf2.predict(X)))


# 设置提前停止
print('='*30, '\n', "With early stopping")
X = digits['data']
y = digits['target']
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
clf = xgb.XGBClassifier()
clf.fit(X_train, y_train, early_stopping_rounds=10, eval_metric="auc",
        eval_set=[(X_test, y_test)])

Zeros and Ones from the Digits dataset: binary classification
[[87 0]
[ 1 92]]
[[91 0]
[ 3 86]]

==============================
Iris: multiclass classification
[[19 0 0]
[ 0 31 3]
[ 0 1 21]]
[[31 0 0]
[ 0 16 0]
[ 0 3 25]]

==============================
Boston Housing: regression
9.860776812557337
15.942418468446029

==============================
Parameter optimization
Fitting 5 folds for each of 9 candidates, totalling 45 fits
[Parallel(n_jobs=1)]: Using backend SequentialBackend with 1 concurrent workers.
0.6699572097100618
{‘max_depth’: 2, ‘n_estimators’: 100}

==============================
Pickling sklearn API models
True

==============================
With early stopping
[0] validation_0-auc:0.999497
Will train until validation_0-auc hasn’t improved in 10 rounds.
[1] validation_0-auc:0.999497
[2] validation_0-auc:0.999497
[3] validation_0-auc:0.999749
[4] validation_0-auc:0.999749
[5] validation_0-auc:0.999749
[6] validation_0-auc:0.999749
[7] validation_0-auc:0.999749
[8] validation_0-auc:0.999749
[9] validation_0-auc:0.999749
[10] validation_0-auc:1
[11] validation_0-auc:1
[12] validation_0-auc:1
[13] validation_0-auc:1
[14] validation_0-auc:1
[Parallel(n_jobs=1)]: Done 45 out of 45 | elapsed: 8.4s finished

[15] validation_0-auc:1
[16] validation_0-auc:1
[17] validation_0-auc:1
[18] validation_0-auc:1
[19] validation_0-auc:1
[20] validation_0-auc:1
Stopping. Best iteration:
[10] validation_0-auc:1

XGBClassifier(base_score=0.5, booster=‘gbtree’, colsample_bylevel=1,
colsample_bytree=1, gamma=0, learning_rate=0.1, max_delta_step=0,
max_depth=3, min_child_weight=1, missing=None, n_estimators=100,
n_jobs=1, nthread=None, objective=‘binary:logistic’, random_state=0,
reg_alpha=0, reg_lambda=1, scale_pos_weight=1, seed=None,
silent=True, subsample=1)

*使用klearn进行并行

import os

if __name__ == "__main__":
    # NOTE: on posix systems, this *has* to be here and in the
    # `__name__ == "__main__"` clause to run XGBoost in parallel processes
    # using fork, if XGBoost was built with OpenMP support. Otherwise, if you
    # build XGBoost without OpenMP support, you can use fork, which is the
    # default backend for joblib, and omit this.
    try:
        from multiprocessing import set_start_method
    except ImportError:
        raise ImportError("Unable to import multiprocessing.set_start_method."
                          " This example only runs on Python 3.4")
    set_start_method("forkserver")

    import numpy as np
    from sklearn.model_selection import GridSearchCV
    from sklearn.datasets import load_boston
    import xgboost as xgb

    rng = np.random.RandomState(31337)

    print("Parallel Parameter optimization")
    boston = load_boston()

    os.environ["OMP_NUM_THREADS"] = "2"  # or to whatever you want
    y = boston['target']
    X = boston['data']
    xgb_model = xgb.XGBRegressor()
    clf = GridSearchCV(xgb_model, {'max_depth': [2, 4, 6],
                                   'n_estimators': [50, 100, 200]}, verbose=1,
                       n_jobs=2)
    clf.fit(X, y)
    print(clf.best_score_)
    print(clf.best_params_)

使用sklearn评估结果

import xgboost as xgb
import numpy as np
from sklearn.datasets import make_hastie_10_2

X, y = make_hastie_10_2(n_samples=2000, random_state=42)

# Map labels from {-1, 1} to {0, 1}
labels, y = np.unique(y, return_inverse=True)

X_train, X_test = X[:1600], X[1600:]
y_train, y_test = y[:1600], y[1600:]

param_dist = {'objective':'binary:logistic', 'n_estimators':2}

clf = xgb.XGBModel(**param_dist)
# Or you can use: clf = xgb.XGBClassifier(**param_dist)

clf.fit(X_train, y_train,
        eval_set=[(X_train, y_train), (X_test, y_test)], 
        eval_metric='logloss',
        verbose=True)

# Load evals result by calling the evals_result() function
evals_result = clf.evals_result()

print('Access logloss metric directly from validation_0:')
print(evals_result['validation_0']['logloss'])

print('')
print('Access metrics through a loop:')
for e_name, e_mtrs in evals_result.items():
    print('- {}'.format(e_name))
    for e_mtr_name, e_mtr_vals in e_mtrs.items():
        print(' - {}'.format(e_mtr_name))
        print(' - {}'.format(e_mtr_vals))
 

print

【机器学习笔记】——Bagging、Boosting、Stacking（RF / Adaboost / Boosting Tree / GBM / GBDT / XGBoost / LightGBM）

目 录

集成学习

概念

思维导图

Bagging算法

概念

编程（分类）

随机森林

扩展

Extremely randomized Trees

*Totally Random Trees Embedding

*Isolation Forest

编程（分类）

为什么说Bagging通过减小方差来提升精度

Boosting

Adaboost

实例

模型推导

为什么说Boosting通过减小偏差来提升精度

*多分类任务

adaboost M1方法

adaboost MH方法

对多分类输出进行二进制编码

编程（分类）

*训练误差分析

*过拟合分析

总结

提升树

实例

梯度提升

GBDT

GBDT实例

GBDT为什么要用负梯度来近似残差？

GBDT编程（回归）

梯度提升 Vs Adaboost

XGboost

XGboost编程

数据接口

设置参数

训练模型

提前结束

预测

可视化

实例

基础练习

定义损失函数和评估指标

分步训练

使用最初生成的n棵树进行预测

使用广义线性模型作为基学习器

交叉验证

预测叶子索引

利用sklearn中的数据学习xgboost

*使用klearn进行并行

使用sklearn评估结果

目录