参考资料:
李航《统计学习方法》
《机器学习实战》
推荐阅读:
【机器学习】Boosting与AdaBoost分类与回归原理详解与公式推导
集成学习原理小结
集成学习之Adaboost算法原理小结
1. 基本算法思路
提升(boosting)算法实质上是基于“三个臭皮匠赛过诸葛亮”的思想进行的,也就是通过训练多个弱分类器的(有权)组合形成一个强分类器。
对于提升方法来说,有两个问题需要回答:
- 一是在每一轮如何改变训练数据的权值或概率分布
- 二是如何将弱分类器组合成一个强分类器
AdaBoost是adaptive boosting(自适应boosting)的缩写,在训练过程中,首先先训练一系列弱分类器,然后将这些分类器组合起来,就形成了强分类器。关于训练数据的权值或概率分布,AdaBoost的做法是提高前一轮弱分类器错误分类样本的权重,而降低被正确分类样本的权值。这样一来,没有得到正确分类的数据,由于其权值的加大而受到后一轮的弱分类器的更大关注。通俗来说,就是“爱哭的孩子有奶喝”。对于弱分类器的组合,AdaBoost采取加权多数表决的方法,加大分类误差概率小的弱分类器的权值,使其在表决中起较大的作用。
如下图所示:
左边是数据集,其中直方图的不同宽度表示每个样例上的不同权重。在经过一个分类器之后,加权的预测结果会通过三角形中的alpha值进行加权。每个三角形中输出的加权结果在圆形中求和,从而得到最终的输出结果。
假设数据集中有N个训练样本,每个样本都由一个权重\omega_{ki},表示第k次迭代时第i个样本的权重值,这些权重值组成一个向量D_k,表示第k次迭代时的权重向量。
D_k=(\omega_{k1},\omega_{k2},\omega_{k3},\cdots,\omega_{ki},\cdots,\omega_{kN})
关键步骤的计算如下:
(1)初始化样本权重D_1,一开始,这些权重都初始化成相等值。即D_1=(\omega_{11},\omega_{12},\omega_{13},\cdots,\omega_{1i},\cdots,\omega_{1N}),其中\omega_{1i}=\frac{1}{N}。
(2)迭代训练。
对m=1,2,\cdots,M,每个弱分类器,使用具有权值分布D_m的训练数据集学习,得到基本分类器
G_m(x)=\chi \to {-1,+1}
每次训练都会得到一个错误率\varepsilon,
e_m=\sum_{i=1}^N P(G_m(x_i) \neq y_i)=\sum_{i=1}^N \omega_{mi}I(G_m(x_i) \neq y_i)
计算G_m(x)的系数
\alpha_m = \frac{1}{2} \log \frac{1-e_m}{e_m}
然后在同一数据集上再次训练弱分类器。在分类器的第二次训练当中,将会重新调整每个样本的权重,其中第一次分对的样本的权重将会降低,而第一次分错的样本的权重将会提高。训练数据集的权值分布更新:
D_{m+1}=(\omega_{m+1,1},\cdots,\omega_{m+1,i},\cdots,\omega_{m+1,N})
\omega_{m+1,i}=\frac{\omega_{mi}}{Z_m} \exp (-\alpha_m y_i G_m(x_i)), \text{其中} i=1,2,\cdots,N
这里的Z_m为规范化因子
Z_m=\sum_{i=1}^N \omega_{mi} \exp (-\alpha_m y_i G_m(x_i))
这个规范化因子使D_{m+1}成为一个概率分布。
(3)弱分类器的基本组合
设弱分类器的线性组合为
f(x)=\sum_{m=1}^M \alpha_m G_m(x)
则可以得到最终的分类器
G(x)=sign(f(x))=sign(\sum_{m=1}^M \alpha_m G_m(x))
2. AdaBoost算法流程
AdaBoost算法流程如下:
(1)假设训练数据集具有均匀的权值分布,即初始化每个训练样本在基本分类器的学习中的作用相同。
(2)AdaBoost反复学习基本分类器,在每一轮m=1,2,\cdots,M,依次进行如下操作。
(a)使用当前分布D_m加权的训练数据集,学习基本分类器G_m(x);
(b)计算基本分类器G_m(x)在加权训练数据集上的分类误差率:
e_m=\sum_{i=1}^N P(G_m(x_i) \neq y_i)=\sum_{G_m(x_i) \neq y_i} \omega_{}mi
(c)计算基本分类器G_m(x)的系数\alpha_m。
(d)更新训练数据的权值,为下一轮作准备。
(3)得到基本分类器的线性组合f(x)。
3. AdaBoost算例
假设训练数据为
序号 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 |
---|---|---|---|---|---|---|---|---|---|---|
x | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 |
y | 1 | 1 | 1 | -1 | -1 | -1 | 1 | 1 | 1 | -1 |
假设弱分类器由x>v或x{<}v产生,阈值v使该分类器在训练数据集上分类误差率最低。
def basicClassifierTraining(x, y, D): ## 训练基本分类器
v = 0.5
error = 1000000
best_v = None
best_direct = None
best_pred_y = None
while v<11:
pred_y_pos = np.array([1 if xi > v else -1 for xi in x])
emPositive = sum([D[i] for i in range(len(x)) if pred_y_pos[i] != y[i]])
pred_y_neg = np.array([-1 if xi > v else 1 for xi in x])
emNegetive = sum([D[i] for i in range(len(x)) if pred_y_neg[i] != y[i]])
if emPositive<emNegetive:
direct = 'positive'
em = emPositive
pred_y = pred_y_pos
else:
direct = 'negetive'
em = emNegetive
pred_y = pred_y_neg
if em<error:
error = em
best_v = v
best_direct = direct
best_pred_y = pred_y
v = v + 1
return best_v, error, best_direct, best_pred_y
在本算例中,
D_1=(\omega_{11},\omega_{12},\cdots,\omega_{110})
\omega_{1i}=0.1,\text{其中} i=1,2,\cdots,10
进行迭代操作如下。
对m=1:
(a)在权值分布为D_1的训练数据集上,阈值v=2.5时分类误差最低,所以基本分类器为
G_1(x)=
\begin{cases}
1, x{<}2.5 \\
-1, x>2.5 \\
\end{cases}
(b)G_1(x)在训练数据集上的误差率为e_1=P(G_1(x_i) \neq y_i)=0.3。
(c)计算G_1(x)的系数:\alpha_1=\frac{1}{2} \log \frac{1-e_1}{e_1}=0.4236
def calalpha(error): ## 计算权重
alpha = 0.5*np.log((1-error)/error)
return alpha
(d)更新训练数据的权值分布:
D_2=(\omega_{21},\cdots,\omega_{2i},\cdots,\omega_{210})
\omega_{2i}=\frac{\omega_{1i}}{Z_1} \exp (-\alpha_1 y_i G_1(x_i)),\text{其中}i=1,2,\cdots,10
def updateD(D, alpha, v, x, y, best_pred_y): ## 更新D值
newD = []
Z = 0
for i in range(len(D)):
w = D[i]*np.exp(-alpha*y[i]*best_pred_y[i])
newD.append(w)
Z += w
out = []
for d in newD:
out.append(d/Z)
return out
D_2=(0.07143,0.07143,0.07143,0.07143,0.07143,0.07143,0.16667,0.16667,0.16667,0.07143)
f_1(x)=0.4236G_1(x)
训练器sign[f(x)]在训练数据集上有3个误分类点(点7,8,9)
def sign(x):
if x>0:
return 1
if x<0:
return -1
def calerrorpoint(gx, x, y):
## gx为组合分类器,是一个列表
## 列表的每个元素是一个二元组
## 二元组的第一个元素表示该分类器的权重alpha,第二个元素表示阈值v
num = 0
for i in range(len(x)):
pred_y = 0
for classifier in gx:
alpha = classifier[0]
v = classifier[1]
direct = classifier[2]
if direct == 'positive':
if x[i]>v:
pred_y += 1
else:
pred_y += -1
else:
if x[i]>v:
pred_y += -1
else:
pred_y += 1
## 进行sign函数计算]
pred_y = sign(pred_y)
if pred_y != y[i]:
num +=1
return num
对m=2
同样循环上述的操作(注意此时的训练数据权值为D_2),可以得到基本分类器为
G_1(x)=
\begin{cases}
1, x{<}8.5 \\
-1, x>8.5 \\
\end{cases}
可以得到训练器f(x)=0.4236G_1(x)+0.6496G_2(x),在训练数据集上有3个误分类点。
对m=3,
可以得到基本分类器为
G_1(x)=
\begin{cases}
1, x{>}5.5 \\
-1, x{<}5.5 \\
\end{cases}
可以得到训练器f(x)=0.4236G_1(x)+0.6496G_2(x)+0.7514G_3(x),在训练数据集上有0个误分类点。
所以最终分类器为G(x)=sign[f(x)]=sign[0.4236G_1(x)+0.6496G_2(x)+0.7514G_3(x)]
主程序如下:
if __name__ == '__main__':
x = [0,1,2,3,4,5,6,7,8,9]
y = [1,1,1,-1,-1,-1,1,1,1,-1]
D = [0.1]*10
## 初始化组合分类器, 每个元素是一个三元组
## 三元组的第一个元素表示该分类器的权重alpha,第二个元素表示阈值v,第三个元素表示比较的方向
gx = []
num = 10
m = 1
while num>0 and m<4:
print(str(m) + '-'*10)
v, error, direct, best_pred_y = basicClassifierTraining(x, y, D)
alpha = calalpha(error)
D = updateD(D, alpha, v, x, y, best_pred_y)
gx.append((alpha, v, direct))
print(v, error, direct, alpha)
print(D)
num = calerrorpoint(gx, x, y)
print(num)
m = m + 1
最后的输出:
1----------
2.5 0.30000000000000004 negetive 0.4236489301936017
[0.07142857142857142, 0.07142857142857142, 0.07142857142857142, 0.07142857142857142, 0.07142857142857142, 0.07142857142857142, 0.16666666666666663, 0.16666666666666663, 0.16666666666666663, 0.07142857142857142]
3
2----------
8.5 0.21428571428571427 negetive 0.6496414920651304
[0.04545454545454546, 0.04545454545454546, 0.04545454545454546, 0.16666666666666669, 0.16666666666666669, 0.16666666666666669, 0.10606060606060606, 0.10606060606060606, 0.10606060606060606, 0.04545454545454546]
6
3----------
5.5 0.18181818181818185 positive 0.752038698388137
[0.12499999999999996, 0.12499999999999996, 0.12499999999999996, 0.10185185185185185, 0.10185185185185185, 0.10185185185185185, 0.0648148148148148, 0.0648148148148148, 0.0648148148148148, 0.12499999999999996]
0
4. AdaBoost回归
以下算法流程摘自集成学习之Adaboost算法原理小结
Adaboost R2回归算法过程:
输入:样本集T=\{(x_1,y_1),(x_2,y_2),\cdots, (x_m,y_m)\},弱学习器算法G(x), 弱学习器迭代次数K。
输出:最终的强学习器f(x)
(1)初始化样本权重为
D_1=(\omega_{11},\omega_{12},\cdots,\omega_{1m});\quad \text{其中 }\omega_{1i}=\frac{1}{m};i=1,2,\cdots,m
(2)对于k=1,2,\cdots,K
使用具有权重D_k的权重来训练样本,得到弱学习器G_k(x);
计算训练集上的最大误差E_k=\max |y_i-G_k(x_i)|,i=1,2,\cdots,m。
计算每个样本的相对误差:
线性误差——e_{ki}=\frac{|y_i-G_k(x_i)|}{E_k}
平方误差——e_{ki}=\frac{(y_i-G_k(x_i))^2}{E_k^2}
指数误差——e_{ki}=1-\exp (\frac{-|y_i-G_k(x_i)|}{E_k})
之后,计算回归误差率e_k=\sum_{i=1}^m w_{ki} e_{ki}
然后,计算弱学习器的系数\alpha_k = \frac{e_k}{1-e_k}
接着,更新样本集的权重分布为w_{k+1,i}=\frac{w_{ki}}{Z_k} \alpha_k^{1-e_{ki}},其中Z_k = \sum_{i=1}^m w_{ki} \alpha_{k}^{1-e_{ki}}。
(3)构建最后的强学习器为f(x)=G_{k^{\ast}}(x),其中G_{k^{\ast}}(x)是所有\ln \frac{1}{\alpha_k},k=1,2,\cdots,K的中位数值对应序号k^{\ast}对应的弱分类器。
5. 提升树
提升树(boosting tree)是以分类树或者回归树为基本分类器的提升方法,被认为是统计学习中性能最好的方法之一。对分类问决策树是二叉分类树,对回归问题决策树是二叉回归树。
5.1 提升决策树算法
单层决策树(decision stump,也称决策树桩)是一种简单的决策树(关于决策树可以看我之前的文章——【机器学习系列】决策树)。我们前面提到的那个弱分类器——x>v或x{<}v就可以看成是一个决策树桩。
提升树模型可以看成是决策树的加法模型(即基函数的线性组合):
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)
其中,f_{m-1}(x)为当前模型,经过经验风险极小化确定下一棵决策树的参数\Theta_m:
\hat{\Theta}_m = \arg \min_{\Theta_m} \sum_{i=1}^N L(y_i,f_{m-1}(x_i)+T(x_i;\Theta_m))
5.2 基于单层决策树构建弱分类器
首先,我们先来看一下单层决策树这一弱分类器的生成和择优。
import numpy as np
def stumpClassify(dataMatrix, dimen, threshVal, threshIneq):
retArray = np.ones((np.shape(dataMatrix)[0], 1))
if threshIneq == 'lt':
retArray[dataMatrix[:, dimen] <= threshVal] = -1.0
else:
retArray[dataMatrix[:, dimen] > threshVal] = -1.0
return retArray
首先,stumpClassify()
是通过阈值比较对数据进行分类的。所有在阈值一边的数据会分到类别-1,而在另外一边的数据分到类别+1。这个函数一共有四个输入参数,dataMatrix为输入的特征值,dimen表示分类是依据第几维特征进行的,threshVal为分类阈值,threshIneq为分类模式——‘lt’表示小于阈值的归为-1类,‘gt’表示大于阈值的归为-1类。
def buildStump(dataArr, classLabels, D):
dataMatrix = np.mat(dataArr)
labelMat = np.mat(classLabels).T
m, n = np.shape(dataMatrix)
numsteps = 10.0
bestStump = {}
bestClassEst = np.mat(np.zeros((m, 1)))
minError = np.inf
for i in range(n):
rangeMin = dataMatrix[:,i].min()
rangeMax = dataMatrix[:,i].max()
stepSize = (rangeMax-rangeMin)/numsteps
for j in range(-1, int(numsteps)+1):
for inequal in ['lt', 'gt']:
threshVal = (rangeMin + float(j) * stepSize)
predictedVals = stumpClassify(dataMatrix, i, threshVal, inequal)
errArr = np.mat(np.ones((m,1)))
errArr[predictedVals == labelMat] = 0
weightedError = D.T * errArr
print("spit: dim %d, thresh %.2f, thresh ineqal: %s, the weighted error is %.3f"%(i, threshVal, inequal, weightedError))
if weightedError < minError:
minError = weightedError
bestClassEst = predictedVals.copy()
bestStump['dim'] = i
bestStump['thresh'] = threshVal
bestStump['ineq'] = inequal
return bestStump, minError, bestClassEst
buildStump()
函数会遍历stumpClassify()
函数所有的可能输入值,并找到数据集上最佳的单层决策树。这个函数一共有三层循环,第一层主要是遍历所有的特征维度,第二层以一定的步长进行该特征维度上的取值范围,第三层遍历阈值比较符号——‘gt’或者’lt’,从而找到最佳的单层决策树。
最佳决策树的判断标准为权重向量——D,表示该决策树中各个数据的权重。和前面提到的弱分类器一样,初始情况下,每个训练样本数据的权重是一样的:
m, n = np.shape(np.mat(dataArr))
D = np.mat(np.ones((m,1))/m)
在每一次迭代训练中,我们还是一样地更新\alpha值和权重向量D。
def adaBoostTrainDS(dataArr, classLabels, numIt = 40):
weakClassArr = []
m = np.shape(dataArr)[0]
D = np.mat(np.ones((m, 1))/m)
aggClassEst = np.mat(np.zeros((m,1)))
for i in range(numIt):
bestStump, error, classEst = buildStump(dataArr, classLabels, D)
print("D: ", D.T)
alpha = float(0.5*np.log((1.0-error)/max(error, 1e-16)))
bestStump['alpha'] = alpha
weakClassArr.append(bestStump)
print("classEst: ", classEst.T)
## 更新D值
expon = np.multiply(-1*alpha*np.mat(classLabels).T, classEst)
D = np.multiply(D, np.exp(expon))
D = D/D.sum()
## 错误率累加
aggClassEst += alpha * classEst
print("aggClassEst: ", aggClassEst.T)
aggErrors = np.multiply(np.sign(aggClassEst) != np.mat(classLabels).T, np.ones((m,1)))
errorRate = aggErrors.sum()/m
print('total error: ', errorRate, "\n")
if errorRate == 0.0:
break
return weakClassArr
我们可以自己准备一些数据集来测试一下:
dataArr = [[1],[2],[3],[4],[5],[6],[7],[8],[9],[10]]
classLabels = [1, 1, -1, 1, 1, -1, 1, 1, 1, -1]
weakClassArr = adaBoostTrainDS(dataArr, classLabels, numIt = 40)
需要注意的是,虽然都是二分类,但是上面的代码要求数据的标签应该是+1和-1,而不是0和1。
6. 回归问题——梯度提升回归(GBR, gradient boosting regression)
6.1 梯度提升算法
输入:训练数据集T=\{(x_1, y_1),(x_2,y_2),\cdots,(x_N,y_N)\},x_i \in \chi \subseteq R^n,y_i \in \mathcal{Y} \subseteq R^n;损失函数L(y,f(x));
输出:回归树\tilde{f}(x)
(1)初始化
f_0(x)=\arg \min_c \sum_{i=1}^N L(y_i, c)
(2)对m=1,2,\cdots, M
对i=1,2,\cdots,N,计算损失函数的负梯度在当前模型的值,将它作为残差的估计:
r_{mi} = -[\frac{\partial L(y_i,f(x_i))}{\partial f(x_i)}]_{f(x)=f_{m-1}(x)}
对于平方损失函数,它就是通常所说的残差;对于一般损失函数,它就是残差的近似值。
对r_{mi}拟合一个回归树,得到第m棵树的叶结点区域R_{mj},j=1,2,\cdots,J,以拟合残差的近似值。
接下来,利用线性搜索估计叶结点区域的值,使损失函数极小化。即,对j=1,2,\cdots,J,计算
c_{mj}=\arg \min_{c} \sum_{x_i \in R_{mj}} L(y_i,f_{m-1}(x_i)+c)
最后,更新回归树f_m(x)=f_{m-1} (x)+\sum_{j=1}^J c_{mj} I(x \in R_{mj})
(3)得到回归树
\hat{f}(x)=f_M(x)=\sum_{m=1}^M \sum_{j=1}^J c_{mj} I(x \in R_{mj})
6.2 案例分析
案例来源于sklearn的官网教程——Gradient Boosting regression。
数据集共有442个样本,每个输入有10个特征值。
import matplotlib.pyplot as plt
import numpy as np
from sklearn import datasets, ensemble
from sklearn.inspection import permutation_importance
from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split
diabetes = datasets.load_diabetes()
X, y = diabetes.data, diabetes.target
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.1, random_state=13)
首先,对数据集进行预处理,test_size=0.1
表示训练集占数据集的0.1。test_size
若为浮点时,表示测试集占总样本的百分比;若为整数时,表示测试样本样本数;若为None时,test size自动设置成0.25。random_state=13
表示随机数的种子数为13;若为None时,每次生成的数据都是随机,可能不一样;若为整数时,每次生成的数据都相同。
在sklearn中使用GBR模型特别简单,只需要调用一下函数GradientBoostingRegressor()
,设置一下相关参数就好了。关于数据的调参可以看Gradient Boosting Regressor机器学习超参数调整
params = {'n_estimators': 500,
'max_depth': 4,
'min_samples_split': 5,
'learning_rate': 0.01,
'loss': 'ls'}
reg = ensemble.GradientBoostingRegressor(**params)
reg.fit(X_train, y_train)
评论(0)
您还未登录,请登录后发表或查看评论