统计学习方法之 adaBoost

转载：http://www.hankcs.com/ml/adaboost.html

本文是《统计学习方法》第8章提升方法的笔记，整合了《机器学习实战》中的提升树Python代码，并添加了注解和PR值计算代码。《方法》重理论，但不易理解，《实战》重实践，但缺乏理论基础，特别是AdaBoost算法的解释、提升树与加法模型的关系等。两相结合，应该能获得较为全面的知识。

提升方法AdaBoost算法

提升方法的思路是综合多个分类器，得到更准确的分类结果。    

AdaBoost算法的归类

《统计学习方法》称AdaBoost是提升算法的代表，所谓提升算法，指的是一种常用的统计学习方法，应用广泛且有效。在分类问题中，它通过改变训练样本的权重，学习多个分类器，并将这些分类器进行线性组合，提髙分类的性能。

《机器学习实战》称AdaBoost是最流行的元算法，所谓元算法，指的是“学习算法的算法”。

AdaBoost算法的基本思想

1. 多轮训练，多个分类器

2. 每轮训练增加错误分类样本的权值，降低正确分类样本的权值

3. 降低错误率高的分类器的权值，增加正确率高的分类器的权值

AdaBoost算法

给定一个二类分类的训练数据集

其中，每个样本点由实例与标记组成。实例，标记，是实例空间，是标记集合。AdaBoost利用以下算法，从训练数据中学习一系列弱分类器或基本分类器，并将这些弱分类器线性组合成为一个强分类器。

AdaBoost算法

输入：训练数据集:弱学习算法；

输出：最终分类器

(1)初始化训练数据的权值分布

每个w的下标由两部分构成，前一个数表示当前迭代次数，与D的下标保持一致，后一个数表示第几个权值，与位置保持一致。初始情况下，每个权值都是均等的。

(2)对（这里的M原著未做解释，其实是表示训练的迭代次数，是由用户指定的。每轮迭代产生一个分类器，最终就有M个分类器）：

(a)使用具有权值分布的训练数据集学习，得到基本分类器

(b)计算在训练数据集上的分类误差率

分类误差率这个名字可能产生误解，这里其实是个加权和。

(c）计算的系数

这里的对数是自然对数。表示在最终分类器中的重要性。由式可知，当时，,并且随着的减小而增大，所以分类误差率越小的基本分类器在最终分类器中的作用越大。

为什么一定要用这个式子呢？这与前向分步算法的推导有关，在后面的章节会介绍。

(d)更新训练数据集的权值分布

y只有正负一两种取值，所以上式可以写作：

这里，是规范化因子

它使成为一个概率分布。

由此可知，被基本分类器误分类样本的权值得以扩大，而被正确分类样本的权值却得以缩小。两相比较，误分类样本的权值被放大倍。因此，误分类样本在下一轮学习中起更大的作用。不改变所给的训练数据，而不断改变训练数据权值的分布，使得训练数据在基本分类器的学习中起不同的作用，这是AdaBoost的一个特点。

(3)构建基本分类器的线性组合

得到最终分类器

AdaBoost算法的解释

AdaBoost算法还有另一个解释，即可以认为AdaBoost算法是模型为加法模型、损失函数为指数函数、学习算法为前向分步算法时的二类分类学习方法。

为什么还要学习前向分步算法呢？直接给我AdaBoost的代码不就好了吗？因为只有理解了前向分步算法，才能理解AdaBoost为什么能跟决策树组合起来。

前向分步算法

考虑加法模型（additive model)

其中，为基函数，为基函数的参数，为基函数的系数。显然，是一个加法模型。

在给定训练数据及损失函数的条件下，学习加法模型成为经验风险极小化即损失函数极小化问题：

通常这是一个复杂的优化问题。前向分步算法（forward stage wise algorithm)求解这一优化问题的想法是：因为学习的是加法模型，如果能够从前向后，每一步只学习一个基函数及其系数，逐步逼近优化目标函数式（L应该是loss的缩写，表示一个损失函数，输入正确答案yi和模型预测值，输出损失值），那么就可以简化优化的复杂度。具体地，每步只需优化如下损失函数：

也就是说，原来有M个分类器，现在只专注优化一个。

给定训练数据集。损失函数和基函数的集合,学习加法模型的前向分步算法如下：

算法(前向分步算法）

输入：训练数据集，损失函数，基函数集;

输出：加法模型

(1)初始化

(2)对

（a）极小化损失函数

得到参数

 (b)更新

(3)得到加法模型

这样，前向分步算法将同时求解从m=1到M所有参数的优化问题简化为逐次求解各个的优化问题。

前向分步算法与AdaBoost

由前向分步算法可以推导出AdaBoost，用定理叙述这一关系。

定理 AdaBoost算法是前向分歩加法算法的特例。这时，模型是由基本分类器组成的加法模型，损失函数是指数函数。

证明 前向分步算法学习的是加法模型，当基函数为基本分类器时，该加法模型等价于AdaBoost的最终分类器

由基本分类器及其系数组成，m=1,2,…，M。前向分步算法逐一学习基函数，这一过程与AdaBoost算法逐一学习基本分类器的过程一致。下面证明前向分步算法的损失函数是指数损失函数（exponential loss function)

时，其学习的具体操作等价于AdaBoost算法学习的具体操作。

假设经过m-1轮迭代前向分步算法已经得到:

在第m轮迭代得到和。

目标是使前向分步算法得到的使在训练数据集T上的指数损失最小，即

上式可以表示为

其中，（指数中的加法可以拿出来做乘法）。因为既不依赖α也不依赖于G，所以与最小化无关。但依赖于,随着每一轮迭代而发生改变。

现证使式达到最小的就是AdaBoost算法所得到的。

求解式可分两步：

首先，求。对任意a>0,使式最小的由下式得到：

其中，。

此分类器即为AdaBoost算法的基本分类器，因为它是使第m轮加权训练数据分类误差率最小的基本分类器。

之后，求。中

这个转换很简单，当y和G一致时，指数为负，反之为正，第二个等号也是利用这个原理，只不过换成了用指示变量I表述。

将已求得的代入式，对α求导并使导数为0,即得到使式最小的a。

其中，是分类误差率：

这里的与AdaBoost算法第2(c)步的完全一致。

最后来看每一轮样本权值的更新。由

以及，可得

这与AdaBoost算法第2(d)步的样本权值的更新，只相差规范化因子，因而等价。

提升树

提升树是以分类树或回归树为基本分类器的提升方法。提升树被认为是统计学习中性能最好的方法之一。

提升方法实际采用加法模型（即基函数的线性组合）与前向分步算法。以决策树为基函数的提升方法称为提升树（boosting tree)。对分类问题决策树是二叉分类树，对回归问题决策树是二叉回归树。在原著例题中看到的基本分类器，可以看作是由一个根结点直接连接两个叶结点的简单决策树，即所谓的决策树桩（decision stump)。提升树模型可以表示为决策树的加法模型：

其中，表示决策树；为决策树的参数；M为树的个数。

提升树算法

提升树算法采用前向分步算法。首先确定初始提升树/e(x)=0,第m歩的模型是

其中，为当前模型，通过经验风险极小化确定下一棵决策树的参数

由于树的线性组合可以很好地拟合训练数据，即使数据中的输入与输出之间的关系很复杂也是如此，所以提升树是一个髙功能的学习算法。

不同问题有大同小异的提升树学习算法，其主要区别在于使用的损失函数不同。包括用平方误差损失函数的回归问题，用指数损失函数的分类问题，以及用一般损失函数的一般决策问题。

对于二类分类问题，提升树算法只需将AdaBoost算法中的基本分类器限制为二类分类树即可，可以说这时的提升树算法是AdaBoost算法的特殊情况，接下来通过《机器学习实战》中的代码学习其应用。

提升树的Python实现

AdaBoost+决策树=提升树，来看看具体用Python怎么实现。

测试数据

老规矩，写代码前先看数据，有什么样的数据写什么样的代码。

考虑到学习代码的最佳方式是单步，而单步的时候，测试数据越简单越好。所以这里先上一份简单的测试数据：

1. def loadSimpData():
2.     """
3. 加载简单数据集
4.     :return:
5.     """
6.     datMat = matrix([[1., 2.1],
7.                      [2., 1.1],
8.                      [1.3, 1.],
9.                      [1., 1.],
10.                      [2., 1.]])
11.     classLabels = [1.0, 1.0, -1.0, -1.0, 1.0]
12.     return datMat, classLabels

数据很简单，就是一个二维平面上的5个点。

写一段可视化代码：

1. def plotData(datMat, classLabels):
2.     xcord0 = []
3.     ycord0 = []
4.     xcord1 = []
5.     ycord1 = []
6.  
7.     for i in range(len(classLabels)):
8.         if classLabels[i]==1.0:
9.             xcord1.append(datMat[i,0]), ycord1.append(datMat[i,1])
10.         else:
11.             xcord0.append(datMat[i,0]), ycord0.append(datMat[i,1])
12.     fig = plt.figure()
13.     ax = fig.add_subplot(111)
14.     ax.scatter(xcord0,ycord0, marker='s', s=90)
15.     ax.scatter(xcord1,ycord1, marker='o', s=50, c='red')
16.     plt.title('decision stump test data')
17.     plt.show()

画出来是这种效果：

单层决策树分类

单层决策树就是一个树桩，只能利用一个维度的特征进行分类。

1. def stumpClassify(dataMatrix, dimen, threshVal, threshIneq):  # just classify the data
2.     """
3. 用只有一层的树桩决策树对数据进行分类
4.     :param dataMatrix: 数据
5.     :param dimen: 特征的下标
6.     :param threshVal: 阈值
7.     :param threshIneq: 大于或小于
8.     :return: 分类结果
9.     """
10.     retArray = ones((shape(dataMatrix)[0], 1))
11.     if threshIneq == 'lt':
12.         retArray[dataMatrix[:, dimen] <= threshVal] = -1.0
13.     else:
14.         retArray[dataMatrix[:, dimen] > threshVal] = -1.0
15.     return retArray

dimen决定了这个树桩能用的唯一特征。

构建决策树桩

上个函数到底取什么参数才好呢？给定训练数据集的权重向量，利用该向量计算错误率加权和，取最低的参数作为训练结果。

1. def buildStump(dataArr, classLabels, D):
2.     """
3. 构建决策树(一个树桩)
4.     :param dataArr: 数据特征矩阵
5.     :param classLabels: 标签向量
6.     :param D: 训练数据的权重向量
7.     :return: 最佳决策树,最小的错误率加权和,最优预测结果
8.     """
9.     dataMatrix = mat(dataArr)
10.     labelMat = mat(classLabels).T
11.     m, n = shape(dataMatrix)
12.     numSteps = 10.0
13.     bestStump = {}
14.     bestClasEst = mat(zeros((m, 1)))
15.     minError = inf  # 将错误率之和设为正无穷
16.     for i in range(n):  # 遍历所有维度
17.         rangeMin = dataMatrix[:, i].min()   #该维的最小最大值
18.         rangeMax = dataMatrix[:, i].max()
19.         stepSize = (rangeMax - rangeMin) / numSteps
20.         for j in range(-1, int(numSteps) + 1):  # 遍历这个区间
21.             for inequal in ['lt', 'gt']:  # 遍历大于和小于
22.                 threshVal = (rangeMin + float(j) * stepSize)
23.                 predictedVals = stumpClassify(dataMatrix, i, threshVal,
24.                                               inequal)  # 使用参数 i, j, lessThan 调用树桩决策树分类
25.                 errArr = mat(ones((m, 1)))
26.                 errArr[predictedVals == labelMat] = 0 # 预测正确的样本对应的错误率为0,否则为1
27.                 weightedError = D.T * errArr  # 计算错误率加权和
28.                 # print "split: dim %d, thresh %.2f, thresh ineqal: %s, the weighted error is %.3f" % (i, threshVal, inequal, weightedError)
29.                 if weightedError < minError:    # 记录最优树桩决策树分类器
30.                     minError = weightedError
31.                     bestClasEst = predictedVals.copy()
32.                     bestStump['dim'] = i
33.                     bestStump['thresh'] = threshVal
34.                     bestStump['ineq'] = inequal
35.     return bestStump, minError, bestClasEst

上面的这一句

1. weightedError = D.T * errArr  # 计算错误率加权和

其实对应着这个公式：

AdaBoost训练

有了上面的基础，就不难看懂完整的训练代码了：

1. def adaBoostTrainDS(dataArr, classLabels, numIt=40):
2.     """
3. 基于单层决策树的ada训练
4.     :param dataArr: 样本特征矩阵
5.     :param classLabels: 样本分类向量
6.     :param numIt: 迭代次数
7.     :return: 一系列弱分类器及其权重,样本分类结果
8.     """
9.     weakClassArr = []
10.     m = shape(dataArr)[0]
11.     D = mat(ones((m, 1)) / m)  # 将每个样本的权重初始化为均等
12.     aggClassEst = mat(zeros((m, 1)))
13.     for i in range(numIt):
14.         bestStump, error, classEst = buildStump(dataArr, classLabels, D)  # 构建树桩决策树,这是一个若分类器,只能利用一个维度做决策
15.         # print "D:",D.T
16.         alpha = float(
17.             0.5 * log((1.0 - error) / max(error, 1e-16)))  # 计算 alpha, 防止发生除零错误
18.         bestStump['alpha'] = alpha
19.         weakClassArr.append(bestStump)  # 保存树桩决策树
20.         # print "classEst: ",classEst.T
21.         expon = multiply(-1 * alpha * mat(classLabels).T, classEst)  # 每个样本对应的指数,当预测值等于y的时候,恰好为-alpha,否则为alpha
22.         D = multiply(D, exp(expon))  # 计算下一个迭代的D向量
23.         D = D / D.sum() # 归一化
24.         # 计算所有分类器的误差,如果为0则终止训练
25.         aggClassEst += alpha * classEst
26.         # print "aggClassEst: ",aggClassEst.T
27.         aggErrors = multiply(sign(aggClassEst) != mat(classLabels).T, ones((m, 1))) # aggClassEst每个元素的符号代表分类结果,如果与y不等则表示错误
28.         errorRate = aggErrors.sum() / m
29.         print "total error: ", errorRate
30.         if errorRate == 0.0: break
31.     return weakClassArr, aggClassEst

其中alpha的计算：

1.         alpha = float(
2.             0.5 * log((1.0 - error) / max(error, 1e-16)))  # 计算 alpha, 防止发生除零错误

对应着

之后记录下这个alpha和对应的决策树桩，然后算法更新训练数据的权重：

1.         expon = multiply(-1 * alpha * mat(classLabels).T, classEst)  # 每个样本对应的指数,当预测值等于y的时候,恰好为-alpha,否则为alpha
2.         D = multiply(D, exp(expon))  # 计算下一个迭代的D向量
3.         D = D / D.sum() # 归一化

对应着

更新完毕后就可以利用这个权值分布训练下一个决策树桩了，当然，作为一个节省时间的策略，可以检查一下当前错误率是否满足要求，如果满足，则终止训练。

提升树分类

分类很简单，就是一个多数表决的过程：

1. def adaClassify(datToClass, classifierArr):
2.     dataMatrix = mat(datToClass)  
3.     m = shape(dataMatrix)[0]
4.     aggClassEst = mat(zeros((m, 1)))
5.     for i in range(len(classifierArr)):
6.         classEst = stumpClassify(dataMatrix, classifierArr[i]['dim'], \
7.                                  classifierArr[i]['thresh'], \
8.                                  classifierArr[i]['ineq'])  
9.         aggClassEst += classifierArr[i]['alpha'] * classEst
10.         print aggClassEst
11.     return sign(aggClassEst)

调用方法

利用上述简单数据集进行测试：

1. datArr, labelArr = loadSimpData()
2. classifierArr, aggClassEst = adaBoostTrainDS(datArr, labelArr, 30)
3. print adaClassify([0, 0], classifierArr)

输出

1. total error:  0.2
2. total error:  0.2
3. total error:  0.0
4. [[-0.69314718]]
5. [[-1.66610226]]
6. [[-2.56198199]]
7. [[-1.]]

可见前向分步训练过程中的误差率是逐渐降低的：

1. total error:  0.2
2. total error:  0.2
3. total error:  0.0

分类过程中，随着加法模型的不断叠加，对于（0，0）这个点，其累加结果是“负”得越来越厉害的，最终取符号输出-1类别。对应训练数据：

的确应该属于-1。

复杂数据集

《机器学习实战》还给出了一个马疝病数据集：horseColicTraining2.txt

加载代码：

1. def loadDataSet(fileName):
2.     """
3.  从文件中加载以制表符分割的数据集(浮点数格式)
4.     :param fileName:
5.     :return:
6.     """
7.     numFeat = len(open(fileName).readline().split('\t'))
8.     dataMat = []
9.     labelMat = []
10.     fr = open(fileName)
11.     for line in fr.readlines():
12.         lineArr = []
13.         curLine = line.strip().split('\t')
14.         for i in range(numFeat - 1):
15.             lineArr.append(float(curLine[i]))
16.         dataMat.append(lineArr)
17.         labelMat.append(float(curLine[-1]))
18.     return dataMat, labelMat

训练代码：

1. datArr, labelArr = loadDataSet("horseColicTraining2.txt")
2. classifierArr, aggClassEst = adaBoostTrainDS(datArr, labelArr, 10)

我们可以直接将aggClassEst打印出来观察对训练集的预测结果，但这样太抽象，而且我们希望计算准确率和召回率，特别是在测试集上的准确率和召回率。

准确率与召回率

对于二分类问题：

准确率的计算方法为

1. TP / (TP + FP)

召回率的计算方法为

1. TP / (TP + FN)

测试集：horseColicTest2.txt

写一段代码计算在不同的数据集下的准确率与召回率：

1. def evaluate(aggClassEst, classLabels):
2.     """
3. 计算准确率与召回率
4.     :param aggClassEst:
5.     :param classLabels:
6.     :return: P, R
7.     """
8.     TP = 0.
9.     FP = 0.
10.     TN = 0.
11.     FN = 0.
12.     for i in range(len(classLabels)):
13.         if classLabels[i] == 1.0:
14.             if (sign(aggClassEst[i]) == classLabels[i]):
15.                 TP += 1.0
16.             else:
17.                 FP += 1.0
18.         else:
19.             if (sign(aggClassEst[i]) == classLabels[i]):
20.                 TN += 1.0
21.             else:
22.                 FN += 1.0
23.  
24.     return TP / (TP + FP), TP / (TP + FN)
25.  
26.  
27. def train_test(datArr, labelArr, datArrTest, labelArrTest, num):
28.     classifierArr, aggClassEst = adaBoostTrainDS(datArr, labelArr, num)
29.     prTrain = evaluate(aggClassEst, labelArr)
30.     aggClassEst = adaClassify(datArrTest, classifierArr)
31.     prTest = evaluate(aggClassEst, labelArrTest)
32.     return prTrain, prTest
33.  
34.  
35. datArr, labelArr = loadDataSet("horseColicTraining2.txt")
36. datArrTest, labelArrTest = loadDataSet("horseColicTest2.txt")
37. prTrain, prTest = train_test(datArr, labelArr, datArrTest, labelArrTest, 10)
38. print prTrain, prTest

输出

1. (0.8595505617977528, 0.7766497461928934) (0.7872340425531915, 0.8604651162790697)

事实上，通过调整分类器的数量（train_test的最后一个参数），可以得到不同的性能。《机器学习实战》验证了1到10000个分类器的错误率：

上图说明分类器并不是越多越好，50是个最好的值，过了这个值，模型发生过拟合，性能越来越低。

ROC曲线

给定分类器在某个数据集上的输出，我们就能计算出假阳率与真阳率，这两个值决定一个坐标点。以假阳率作为x轴，真阳率作为y轴，制定坐标系。

完美的分类器应该是左上角那个点，分类器越接近左上角，就越完美。这样我们就可以较为直观地比较两个不同的分类器了。

在决策函数中，我们使用的是sign来二值化预测强度，以得到最终分类。也就是说，我们取的阈值是0。分类强度离阈值越远，分类的可信度越高。如果我们改变阈值，我们就能得到同一个分类器在坐标系中不同的点，将它们按照阈值大小顺序连成线，就能得到一条ROC曲线。完美的分类器的ROC曲线应该是正方形的左上角对应的两条边，而随机猜测的ROC曲线则是连接左下角与右上角的一条直线。ROC曲线下的面积（AUC）反映了分类器的平均性能。

绘制代码：

1. def plotROC(predStrengths, classLabels):
2.     import matplotlib.pyplot as plt
3.     cur = (1.0, 1.0)  # 中间变量,初始状态为右上角
4.     ySum = 0.0  # variable to calculate AUC
5.     numPosClas = sum(array(classLabels) == 1.0)  # TP
6.     yStep = 1 / float(numPosClas)
7.     xStep = 1 / float(len(classLabels) - numPosClas)
8.     sortedIndicies = predStrengths.argsort()  # 按元素值排序后的下标,逆序
9.     fig = plt.figure()
10.     fig.clf()
11.     ax = plt.subplot(111)
12.     # 遍历所有的值
13.     for index in sortedIndicies.tolist()[0]:
14.         if classLabels[index] == 1.0:  # 预测正确
15.             delX = 0  # 真阳率不变
16.             delY = yStep  # 假阳率减小
17.         else:
18.             delX = xStep
19.             delY = 0
20.             ySum += cur[1]  # ROC面积的一个小长条
21.         # 从 cur 到 (cur[0]-delX,cur[1]-delY) 画一条线
22.         ax.plot([cur[0], cur[0] - delX], [cur[1], cur[1] - delY], c='b')
23.         cur = (cur[0] - delX, cur[1] - delY)
24.     ax.plot([0, 1], [0, 1], 'b--')  # 随机猜测的ROC线
25.     plt.xlabel('False positive rate')
26.     plt.ylabel('True positive rate')
27.     plt.title('ROC curve for AdaBoost horse colic detection system')
28.     ax.axis([0, 1, 0, 1])
29.     plt.show()
30.     print "the Area Under the Curve is: ", ySum * xStep
31. plotROC(aggClassEst.T, labelArr)

效果：

随机猜测的ROC曲线是一条y=x直线，这是因为对真阳率和假阳率相等，意味着分类器猜对和猜错的比率相等，说明该分类器就是在随机猜。