1、前向分步算法:
考虑加法模型:
f(x)=∑m=1Mβmb(x;γm)
f
(
x
)
=
∑
m
=
1
M
β
m
b
(
x
;
γ
m
)
其中 b(x;γm) b ( x ; γ m ) 为基函数, γm γ m 是基函数的参数, βm β m 是基函数的系数。
在给定训练数据集及损失函数 L(y,f(x)) L ( y , f ( x ) ) 的条件下,学习加法模型 f(x) f ( x ) 已经成为经验风险极小化即损失函数极小化问题:
minβm,γm∑i=1NL(yi,∑m=1Mβmb(x;γm))
min
β
m
,
γ
m
∑
i
=
1
N
L
(
y
i
,
∑
m
=
1
M
β
m
b
(
x
;
γ
m
)
)
通常这是一个复杂的优化问题,前向分步算法求解这一优化问题的想法是:因为学习的是加法模型,如果能够从前向后每一步只学习一个基函数及其系数,逐步逼近优化目标函数式,那么就可以简化优化的复杂度,具体的,每步只需要优化如下损失函数:
minβ,γ∑i=1NL(yi,βb(x;γ))
min
β
,
γ
∑
i
=
1
N
L
(
y
i
,
β
b
(
x
;
γ
)
)
前向分步算法:
输入:训练数据集 T={(x1,y1),(x2,y2),...,(xn,yn)},xi∈Rn,yi∈{−1,+1} T = { ( x 1 , y 1 ) , ( x 2 , y 2 ) , . . . , ( x n , y n ) } , x i ∈ R n , y i ∈ { − 1 , + 1 } 损失函数 L(y,f(x)) L ( y , f ( x ) ) ,基函数 {b(x;γ)} { b ( x ; γ ) }
输出:加法模型 f(x) f ( x )
(1)初始化 f0(x)=0 f 0 ( x ) = 0
(2)对于 m=1,2,3,...,M m = 1 , 2 , 3 , . . . , M
(a)极小化损失函数:
(βm,γm)=argminβ,γ∑i=1NL(yi,fm−1(xi)+βb(xi;γ))
(
β
m
,
γ
m
)
=
a
r
g
min
β
,
γ
∑
i
=
1
N
L
(
y
i
,
f
m
−
1
(
x
i
)
+
β
b
(
x
i
;
γ
)
)
得到参数 β,γm β , γ m
(b)跟新
fm(x)=fm−1(x)+βmb(x;γm)
f
m
(
x
)
=
f
m
−
1
(
x
)
+
β
m
b
(
x
;
γ
m
)
(3)得到加法模型
f(x)=fM(x)=∑m=1Mβmb(x;γm)
f
(
x
)
=
f
M
(
x
)
=
∑
m
=
1
M
β
m
b
(
x
;
γ
m
)
这样,前向分步算法将同时求解从m=1到M的所有参数 βm,γm β m , γ m 的优化问题简化为逐次求解各个 βm,γm β m , γ m 的优化问题。Adaboost算法是前向分步算法的特例,当前向分步算法的损失函数是指数损失函数时:
L(y,f(x))=exp[−yf(x)]
L
(
y
,
f
(
x
)
)
=
exp
[
−
y
f
(
x
)
]
其学习的具体操作等价于Adaboost算法学习的具体操作,这里就不给出证明了。
2、提升树(Boosting Tree)
提升方法实际采用加法模型与前向分步算法,以决策树为基函数的提升方法称为提升树,提升树模型可以表示为决策树加法模型:
fM(x)=∑m=1MT(x;Θm)
f
M
(
x
)
=
∑
m
=
1
M
T
(
x
;
Θ
m
)
对于二分类问题,提升树的算法只需要将Adaboost算法的中基本分类器限制为二类分类树即可,下面是回归问题的提升树。对于训练数据集 T={(x1,y1),(x2,y2)...,(xn,yn)} T = { ( x 1 , y 1 ) , ( x 2 , y 2 ) . . . , ( x n , y n ) } ,在前面CART树的回归问题中,如果将输入空间划分为J个互不相交的区域 R1,R2,...,RJ R 1 , R 2 , . . . , R J ,并且在每个区域上确定输出常量 cj c j ,那么这棵树可以表示为:
T(x;Θ)=∑j=1JcjI(x∈Rj)
T
(
x
;
Θ
)
=
∑
j
=
1
J
c
j
I
(
x
∈
R
j
)
其中参数 Θ={(R1,c1),(R2,c2),(R3,c3),...,(RJ,cJ)} Θ = { ( R 1 , c 1 ) , ( R 2 , c 2 ) , ( R 3 , c 3 ) , . . . , ( R J , c J ) } 表示树的区域划分和各区域的常数,J是回归树的复杂度即叶节点个数。回归问题提升树使用以下前向分步算法:
f0(x)=0
f
0
(
x
)
=
0
fm(x)=fm−1(x)+T(x;Θm),m=1,2,3,...,M
f
m
(
x
)
=
f
m
−
1
(
x
)
+
T
(
x
;
Θ
m
)
,
m
=
1
,
2
,
3
,
.
.
.
,
M
fM(x)=∑m=1MT(x;Θm)
f
M
(
x
)
=
∑
m
=
1
M
T
(
x
;
Θ
m
)
前向分步算法的第m步,给定当前模型 fm−1(x) f m − 1 ( x ) ,需求解:
Θ^m=argminΘm∑i=1NL(yi,fm−1(xi)+T(xi;Θm))
Θ
^
m
=
a
r
g
min
Θ
m
∑
i
=
1
N
L
(
y
i
,
f
m
−
1
(
x
i
)
+
T
(
x
i
;
Θ
m
)
)
得到的 Θ^m Θ ^ m 为第m棵树的参数。
当采用平方误差损失函数时:
L(y,f(x))=(y−f(x))2=L(y,fm−1(x)+T(x;Θm))=[y−fm−1(x)−T(x;Θm)]2=[r−T(x;Θm)]2
L
(
y
,
f
(
x
)
)
=
(
y
−
f
(
x
)
)
2
=
L
(
y
,
f
m
−
1
(
x
)
+
T
(
x
;
Θ
m
)
)
=
[
y
−
f
m
−
1
(
x
)
−
T
(
x
;
Θ
m
)
]
2
=
[
r
−
T
(
x
;
Θ
m
)
]
2
这里: r=y−fm−1(x) r = y − f m − 1 ( x ) 是当前模型拟合数据的残差,所以对回归问题的提升树算法来说,只需要简单的拟合当前模型的残差。算法描述如下:
输入:训练数据集 T={(x1,y1),(x2,y2)...,(xn,yn)} T = { ( x 1 , y 1 ) , ( x 2 , y 2 ) . . . , ( x n , y n ) }
输出:提升树 fM(x) f M ( x )
(1)初始化 f0(x)=0 f 0 ( x ) = 0
(2)对于 m=1,2,...,M m = 1 , 2 , . . . , M
(a)计算残差:
rmi=yi−fm−1(xi),i=1,2,...,N
r
m
i
=
y
i
−
f
m
−
1
(
x
i
)
,
i
=
1
,
2
,
.
.
.
,
N
(b)拟合残差 rmi r m i 学习的一个回归树 T(x;Θm) T ( x ; Θ m )
(c)更新
fm(x)=fm−1(x)+T(x;Θm)
f
m
(
x
)
=
f
m
−
1
(
x
)
+
T
(
x
;
Θ
m
)
(3)得到提升树
fM(x)=∑m=1MT(x;Θm)
f
M
(
x
)
=
∑
m
=
1
M
T
(
x
;
Θ
m
)
3、梯度提升(gradient boosting)
当损失函数是平方误差和指数损失函数时,提升树的每一步的优化算法是很简单的,但是对于一般损失函数而言,每一步并不是那么容易优化,针对这一问题,Freidman提出了梯度提升算法,这是利用最速下降法的近似方法,其关键是利用损失函数的负梯度在当前模型的值:
−[∂L(y,f(xi))∂f(xi)]f(x)=fm−1(x)
−
[
∂
L
(
y
,
f
(
x
i
)
)
∂
f
(
x
i
)
]
f
(
x
)
=
f
m
−
1
(
x
)
作为回归问题提升树算法中的残差的近似值,拟合一个回归树。具体的算法如下:
输入:训练数据集 T={(x1,y1),(x2,y2)...,(xn,yn)} T = { ( x 1 , y 1 ) , ( x 2 , y 2 ) . . . , ( x n , y n ) }
输出:提升树 fM(x) f M ( x )
(1)初始化
f0(x)=argminc∑i=1NL(yi,c)
f
0
(
x
)
=
a
r
g
min
c
∑
i
=
1
N
L
(
y
i
,
c
)
估计使损失函数极小化的常数值,它是只有一个根节点的树。
(2)对于 m=1,2,...,M m = 1 , 2 , . . . , M
(a)对 i=1,2,...,N i = 1 , 2 , . . . , N 计算:
rmi=−[∂L(yi,f(xi))∂f(xi)]f(x)=fm−1(x)
r
m
i
=
−
[
∂
L
(
y
i
,
f
(
x
i
)
)
∂
f
(
x
i
)
]
f
(
x
)
=
f
m
−
1
(
x
)
该公式计算损失函数的负梯度在当前模型的值,并将它作为残差的估计,对于平方损失函数,它就是通常所说的残差;对于一般的损失函数,它就是残差的近似值。
(b)对 rmi r m i 拟合一个回归树,得到第m棵树的叶节点区域 Rmj,j=1,2,...,J R m j , j = 1 , 2 , . . . , J
(c)对 j=1,2,...,J j = 1 , 2 , . . . , J ,计算
cmj=argminc∑xi∈RmjL(yi,fm−1(xi)+c)
c
m
j
=
a
r
g
min
c
∑
x
i
∈
R
m
j
L
(
y
i
,
f
m
−
1
(
x
i
)
+
c
)
(d)跟新 fm(x)=fm−1(x)+∑Jj=1cmjI(x∈Rmj) f m ( x ) = f m − 1 ( x ) + ∑ j = 1 J c m j I ( x ∈ R m j )
(3)得到回归树:
f^(x)=fM(x)=∑m=1M∑j=1JcmjI(x∈Rmj)
f
^
(
x
)
=
f
M
(
x
)
=
∑
m
=
1
M
∑
j
=
1
J
c
m
j
I
(
x
∈
R
m
j
)
4、提升树(Boosting Tree)的实现
#coding=utf-8
import numpy as np
# 加载数据
def loadDataSet(fileName): #general function to parse tab -delimited floats
dataMat = [] #assume last column is target value
fr = open(fileName)
for line in fr.readlines():
curLine = line.strip().split('\t')
fltLine = map(float,curLine) #map all elements to float()
dataMat.append(fltLine)
return dataMat
# 根据特征和特征值把数据集划分为两个数据集,一个比特征值大,一个比特征值小
def splitData(data_array,col,value):
array_1 = data_array[data_array[:,col] >= value,:]
array_2 = data_array[data_array[:,col] < value,:]
return array_1,array_2
# 计算平方误差
def getErr(data_array):
return np.var(data_array[:,-1]) * data_array.shape[0]
# 返回叶子节点
def regLeaf(data_array):
return np.mean(data_array[:,-1])
#选择最佳的切分点,使用的方法为CART回归,为二叉树
def get_best_split(data_array,ops = (1,4)):
tolS = ops[0] # 误差阈值,用于控制迭代次数
tolN = ops[1] # 每个划分的最小样本个数
if len(set(data_array[:,-1])) == 1:
return None,regLeaf(data_array)
m,n = data_array.shape
best_S = np.inf; # 保存最优的平方误差
best_col = 0; # 最优特征
best_value = 0; # 最优切分点
S = getErr(data_array) # 计算原始平方误差
for col in xrange(n-1):
values = set(data_array[:,col])
for value in values:
array_1,array_2 = splitData(data_array,col,value)
if (array_1.shape[0] < tolN) or (array_2.shape[0] < tolN):
continue
total_error = getErr(array_1) + getErr(array_2)
if total_error < best_S:
best_col = col
best_value = value
best_S = total_error
if (S - best_S) < tolS:
return None,regLeaf(data_array)
array_1,array_2 = splitData(data_array,best_col,best_value)
if (array_1.shape[0] < tolN) or (array_2.shape[0] < tolN):
return None,regLeaf(data_array)
return best_col,best_value
# 用一个对象来保存树的每个节点值(特征列,特征值,结果,左分支,右分支)
class node:
def __init__(self,col = -1,value = None, results = None, gb = None, lb = None):
self.col = col
self.value = value
self.results = results
self.gb = gb
self.lb = lb
# 建立GBDT决策树
def buildTree(data_array,ops=(1,4)):
col, val = get_best_split(data_array, ops)
if col == None:
return node(results = val)
else:
array_1,array_2 = splitData(data_array,col,val)
greater_branch = buildTree(array_1,ops)
less_branch = buildTree(array_2,ops)
return node(col = col,value = val,gb = greater_branch ,lb = less_branch )
# 输入一个样本的所有特征,返回样本回归结果
def treeForeCast(tree, inData):
if tree.results != None:
return tree.results
#print 'tree.col:',tree.col
if inData[tree.col] > tree.value:
return treeForeCast(tree.gb, inData)
else:
return treeForeCast(tree.lb, inData)
# 返回一个测试集的所有回归结果值列表
def createForeCast(tree, testData):
m=len(testData)
yHat = np.mat(np.zeros((m,1)))
for i in range(m):
yHat[i,0] = treeForeCast(tree, testData[i])
return yHat
# 生成提升树
def boostTree(data_array,num_iter,ops = (1,4)):
m,n = data_array.shape
x = data_array[:,0:-1] # 保存所有样本特征列
y = data_array[:,-1].reshape((m,1)) # 保存所有样本结果真实值
list_trees = []
for i in xrange(num_iter):
print 'i: ',i
if i == 0:
tree = buildTree(data_array,ops)
list_trees.append(tree)
yHat = createForeCast(tree,x)
else:
r = y - np.array(yHat) # 计算残差
data_array = np.hstack((x,r)) # hstack()函数水平(按列顺序)的把数组给堆叠起来
tree = buildTree(data_array,ops)
list_trees.append(tree)
rHat = createForeCast(tree, x ) # 用残差拟合的结果
yHat = yHat + rHat
return list_trees
#打印树的节点信息
def printtree(tree,indent=''):
# Is this a leaf node?
if tree.results!=None:
print str(tree.results)
else:
# Print the criteria
print str(tree.col)+':'+str(tree.value)+'? '
# Print the branches
print indent+'T->',
printtree(tree.gb,indent+' ')
print indent+'F->',
printtree(tree.lb,indent+' ')
if __name__ == '__main__':
data = loadDataSet('ex0.txt')
data_array = np.array(data)
# tree = buildTree(data_array)
# printtree(tree)
gbdt_results = boostTree(data_array,10)
5、GBDT(Gradient Boosting Decision Tree)的实现
# coding=utf-8
import numpy as np
from sklearn.ensemble import GradientBoostingRegressor
gbdt=GradientBoostingRegressor(
loss='ls' #损失函数,GBDT回归器可选'ls', 'lad', 'huber', 'quantile'。
, learning_rate=0.1 #学习率/步长。
, n_estimators=100 #迭代次数,和learning_rate存在trade-off关系。
, subsample=1 # 样本采样比例。
, min_samples_split=2 #最大特征数或比例。
, min_samples_leaf=1 #最小特征数或比例。
, max_depth=3 # 树的深度,以下参数多数用来设定决策树分裂停止条件。
, init=None
, random_state=None
, max_features=None
, alpha=0.9
, verbose=0
, max_leaf_nodes=None
, warm_start=False
)
train_feat=np.array([[0.00598802,0.569231,0.647059,0.95122,-0.225434,0.837989,0.357258,-0.0030581,-0.383475],
[0.161677,0.743195,0.682353,0.960976,-0.0867052,0.780527,0.282945,0.149847,-0.0529661],
[0.113772,0.744379,0.541176,0.990244,-0.00578035,0.721468,0.43411,-0.318043,0.288136],
[0.0538922,0.608284,0.764706,0.95122,-0.248555,0.821229,0.848604,-0.0030581,0.239407],
[0.173653,0.866272,0.682353,0.95122,0.017341,0.704709,-0.0210016,-0.195719,0.150424]])
train_id=np.array([320,361,364,336,358])
test_feat=np.array([[0.00598802,0.569231,0.647059,0.95122,-0.225434,0.837989,0.357258,-0.0030581,-0.383475],
[0.161677,0.743195,0.682353,0.960976,-0.0867052,0.780527,0.282945,0.149847,-0.0529661],
[0.113772,0.744379,0.541176,0.990244,-0.00578035,0.721468,0.43411,-0.318043,0.288136],
[0.0538922,0.608284,0.764706,0.95122,-0.248555,0.821229,0.848604,-0.0030581,0.239407],
[0.173653,0.866272,0.682353,0.95122,0.017341,0.704709,-0.0210016,-0.195719,0.150424]])
test_id=np.array([320,361,364,336,358])
gbdt.fit(train_feat,train_id) # 第一个参数为样本特征,第二个参数为样本标签
pred=gbdt.predict(test_feat) # 参数为测试样本特征
total_err=0
for i in range(pred.shape[0]):
print pred[i],test_id[i]
err=(pred[i]-test_id[i])/test_id[i]
total_err+=err*err
print total_err/pred.shape[0]