機器學習理論與實戰（七）Adaboost
Adaboost也是一種原理簡單，但很實用的有監督機器學習演算法，它是daptive boosting的簡稱。說到boosting演算法，就不得提一提bagging演算法，他們兩個都是把一些弱分類器組合起來來進行分類的方法，統稱為整合方法（ensemble method）,類似於投資，“不把雞蛋放在一個籃子”，雖然每個弱分類器分類的不那麼準確，但是如果把多個弱分類器組合起來可以得到相當不錯的結果，另外要說的是整合方法還可以組合不同的分類器，而Adaboost和boosting演算法的每個弱分類器的型別都一樣的。他們兩個不同的地方是：boosting的每個弱分類器組合起來的權重不一樣，本節的Adaboost就是一個例子，而bagging的每個弱分類器的組合權重是相等，代表的例子就是random forest。Random forest的每個弱分類器是決策樹，輸出的類別有多個決策樹分類的類別的眾數決定。今天的主題是Adaboost,下面來看看Adaboost的工作原理：
既然Adaboost的每個弱分類器的型別都一樣，那麼怎麼組織安排每個分類器呢？如（圖一）所示：
（圖一）
（圖一）是Adaboost的原理示意圖，左邊矩形表示資料集，中間表示根據特徵閾值來做分類，這樣每一個弱分類器都類似於一個單節點的決策樹，其實就是閾值判斷而已，右邊的三角形對每個弱分類器賦予一個權重，最後根據每個弱分類器的加權組合來判斷總體類別。要注意一下資料集從上到下三個矩形內的直方圖不一樣，這表示每個樣本的權重也發生了變化，樣本權重的一開始初始化成相等的權重，然後根據弱分類器的錯誤率來調整每個弱分類器的全總alpha,如（圖一）中的三角形所示，alpha 的計算如（公式一）所示：
（公式一）
從（公式一）中也能感覺出來，弱分類器權重alpha和弱分類器分類錯誤率epsilon成反比，如果不能看出反比關係，分子分母同時除以epsilon就可以了，而ln是單調函式。這很make sense,當然分類器的錯誤率越高，越不能器重它，它的權重就應該低。同樣的道理，樣本也要區分對待，樣本的權重要用弱分類器權重來計算，其實也是間接靠分類錯誤率，如（公式二）所示：
(公式二)
其中D表示樣本權重向量，有多少個樣本就有多少個權重，下標i表示樣本索引，而上標t表示上一次分類器訓練迭代次數。這樣一直更新迭代，一直到最大迭代次數或者整個分類器錯誤率為0或者不變時停止迭代，就完成了Adaboost的訓練。但是這樣就可以把樣本分開了嗎？下面從一組圖解答這個問題，如（圖二）所示：
（圖二）
這樣Adaboost的原理基本分析完畢，下面進入程式碼實戰階段：
首先來準備個簡單資料集：
from numpy import *
def loadSimpData():
datMat = matrix([[ 1. , 2.1],
[ 2. , 1.1],
[ 1.3, 1. ],
[ 1. , 1. ],
[ 2. , 1. ]])
classLabels = [1.0, 1.0, -1.0, -1.0, 1.0]
return datMat,classLabels上面有5個樣本，接下來就是初始化每個樣本的權重，剛開始相等的：
D = mat(ones((5,1))/5)有了樣本和初始化權重，接下來的任務就是構建一個弱分類器，其實就是一個單節點決策樹，找到決策樹每個特徵維度上對應的最佳閾值以及表示是大於閾值還是小於閾值為正樣本的識別符號。程式碼如下：
def buildStump(dataArr,classLabels,D):
dataMatrix = mat(dataArr); labelMat = mat(classLabels).T
m,n = shape(dataMatrix)
numSteps = 10.0; bestStump = {}; bestClasEst = mat(zeros((m,1)))
minError = inf #init error sum, to +infinity
for i in range(n):#loop over all dimensions
rangeMin = dataMatrix[:,i].min(); rangeMax = dataMatrix[:,i].max();
stepSize = (rangeMax-rangeMin)/numSteps
for j in range(-1,int(numSteps)+1):#loop over all range in current dimension
for inequal in ['lt', 'gt']: #go over less than and greater than
threshVal = (rangeMin + float(j) * stepSize)
predictedVals = stumpClassify(dataMatrix,i,threshVal,inequal)#call stump classify with i, j, lessThan
errArr = mat(ones((m,1)))
errArr[predictedVals == labelMat] = 0
weightedError = D.T*errArr #calc total error multiplied by D
#print "split: dim %d, thresh %.2f, thresh ineqal: %s, the weighted error is %.3f" % (i, threshVal, inequal, weightedError)
if weightedError < minError:
minError = weightedError
bestClasEst = predictedVals.copy()
bestStump['dim'] = i
bestStump['thresh'] = threshVal
bestStump['ineq'] = inequal
return bestStump,minError,bestClasEst注意程式碼中有三個for迴圈，這三個for迴圈其實就是為了完成決策樹的每個特徵維度上對應的最佳閾值以及表示是大於閾值還是小於閾值為正樣本的識別符號，這三個要素。其中it,gt分別表示大於和小於，閾值的選擇是靠增加步長來需找，最終三者的確定是靠決策樹分類錯誤率最小者決定，每個決策樹的分類程式碼如下，很簡單，就是靠閾值判斷：
def stumpClassify(dataMatrix,dimen,threshVal,threshIneq):#just classify the data
retArray = ones((shape(dataMatrix)[0],1))
if threshIneq == 'lt':
retArray[dataMatrix[:,dimen] <= threshVal] = -1.0
else:
retArray[dataMatrix[:,dimen] > threshVal] = -1.0
return retArray有了弱分類器的構造程式碼，下面來看Adaboost的訓練程式碼：
def adaBoostTrainDS(dataArr,classLabels,numIt=40):
weakClassArr = []
m = shape(dataArr)[0]
D = mat(ones((m,1))/m) #init D to all equal
aggClassEst = mat(zeros((m,1)))
for i in range(numIt):
bestStump,error,classEst = buildStump(dataArr,classLabels,D)#build Stump
#print "D:",D.T
alpha = float(0.5*log((1.0-error)/max(error,1e-16)))#calc alpha, throw in max(error,eps) to account for error=0
bestStump['alpha'] = alpha
weakClassArr.append(bestStump) #store Stump Params in Array
#print "classEst: ",classEst.T
expon = multiply(-1*alpha*mat(classLabels).T,classEst) #exponent for D calc, getting messy
D = multiply(D,exp(expon)) #Calc New D for next iteration
D = D/D.sum()
#calc training error of all classifiers, if this is 0 quit for loop early (use break)
aggClassEst += alpha*classEst
#print "aggClassEst: ",aggClassEst.T
aggErrors = multiply(sign(aggClassEst) != mat(classLabels).T,ones((m,1)))
errorRate = aggErrors.sum()/m
print "total error: ",errorRate
if errorRate == 0.0: break
return weakClassArr,aggClassEst上面的程式碼中訓練過程主要任務就是完成（公式二）中的樣本權重D和弱分類器權重alpha的更新，另外還要注意一下，程式碼中迭代了40次，每次都呼叫了buildStump，這就意味著建立了40個弱分類器。當模型收斂後，有了樣本權重和弱弱弱分類器權重，最後就是對測試樣本進行分類，分類程式碼如下：
def adaClassify(datToClass,classifierArr):
dataMatrix = mat(datToClass)#do stuff similar to last aggClassEst in adaBoostTrainDS
m = shape(dataMatrix)[0]
aggClassEst = mat(zeros((m,1)))
for i in range(len(classifierArr)):
classEst = stumpClassify(dataMatrix,classifierArr[i]['dim'],\
classifierArr[i]['thresh'],\
classifierArr[i]['ineq'])#call stump classify
aggClassEst += classifierArr[i]['alpha']*classEst
print aggClassEst
return sign(aggClassEst)考慮到有些做學術的為了比較不同機器學習演算法的好壞，常常需要畫ROC曲線，這裡也給出畫ROC的程式碼：
def plotROC(predStrengths, classLabels):
import matplotlib.pyplot as plt
cur = (1.0,1.0) #cursor
ySum = 0.0 #variable to calculate AUC
numPosClas = sum(array(classLabels)==1.0)
yStep = 1/float(numPosClas); xStep = 1/float(len(classLabels)-numPosClas)
sortedIndicies = predStrengths.argsort()#get sorted index, it's reverse
fig = plt.figure()
fig.clf()
ax = plt.subplot(111)
#loop through all the values, drawing a line segment at each point
for index in sortedIndicies.tolist()[0]:
if classLabels[index] == 1.0:
delX = 0; delY = yStep;
else:
delX = xStep; delY = 0;
ySum += cur[1]
#draw line from cur to (cur[0]-delX,cur[1]-delY)
ax.plot([cur[0],cur[0]-delX],[cur[1],cur[1]-delY], c='b')
cur = (cur[0]-delX,cur[1]-delY)
ax.plot([0,1],[0,1],'b--')
plt.xlabel('False positive rate'); plt.ylabel('True positive rate')
plt.title('ROC curve for AdaBoost horse colic detection system')
ax.axis([0,1,0,1])
plt.show()
print "the Area Under the Curve is: ",ySum*xStep到此位置，Adaboost的程式碼也介紹完了，最終程式的執行結果如（圖三）所示：
（圖三）
而Adaboost的模型ROC執行曲線如（圖四）所示：
（圖四）
最近MIT的幾個人證明了Adaboost可以用一階梯度的角度來解釋，詳見連結
參考文獻：
[1] machinelearning in action.
[2] http://www.robots.ox.ac.uk/~az/lectures/cv/adaboost_matas.pdf