模型解释
集成(Ensemble)是合并多个机器学习模型来构建更强大模型的方法。在机器学习文献中有许多模型都属于这一类,但已证明有两种集成模型对大量分类和回归的数据集都是有效的,二者都以决策树为基础,分别是随机森林(Random Forest, RF)和提升树 (Boosted Tree, BT)。梯度提升树(Gradient Boosting Tree, GBT)或称梯度提升决策树(Gradient Boosted Decision Tree, GBDT)是一种优化的提升回归树,在损失函数拟合方面,提升回归树损失函数拟合用的是平方损失,而梯度提升树则是使用损失函数的负梯度来拟合本轮损失的近似值,进而拟合一个回归树。梯度提升树在工业上用途广泛,属于最流行、最实用的算法之一。
提升回归树(Boosted Regression Tree, BRT)特指梯度提升树,这是最主流的方法,在一些R函数中GBM就是指BRT。
BRT理解与区别*
BRT 是一种拟合统计模型,与传统的拟合吝啬模型有很大不同,其将回归树(regression tree)和增长(boosting)两种方法结合了起来。回归树(regression tree)是一种用递归二进制拆分(recursive binary splits)将因变量与其预测因子联合起来的模型。增长(boosting),则可以结合多个简单模型并提升预测能力。二者结合而成的 BRT 模型,则可以看做是一种加持的回归模型,将一个一个的简单树(simple tree)向前逐步拟合。
增长回归树模型(Boosted Regression Trees)
在gbm包中,采用的是决策树作为基学习器,重要的参数设置如下:
损失函数的形式(distribution)
迭代次数(n.trees)
学习速率(shrinkage)
再抽样比率(bag.fraction)
决策树的深度(interaction.depth)
损失函数的形式容易设定,分类问题一般选择bernoulli分布,而回归问题可以选择gaussian分布。学习速率方面,我们都知道步子迈得太大容易扯着,所以学习速率是越小越好,但是步子太小的话,步数就得增加,也就是训练的迭代次数需要加大才能使模型达到最优,这样训练所需时间和计算资源也相应加大了。gbm作者的经验法则是设置shrinkage参数在0.01-0.001之间,而n.trees参数在3000-10000之间。
BRT
2. BRT在R语言中的应用
1. gbm:gbm的原始R实现
2. xgboost:一个快速高效的梯度提升框架(C++后端)。
3.h2o:一个强大的基于 Java 的接口,提供并行分布式算法和高效的生产化。
4. dismo:一个用于物种分布建模的包,拓展了gbm包的基础功能。
5. caret:堪称机器学习的百宝箱,可以使用gbm方法进行建模。
2.1 GBM
基本构建
gbm(formula = formula(data),
distribution = "bernoulli", #这里需要改成gaussian
data = list(), #我们的数据
weights, #不用管
var.monotone = NULL, #不用管
n.trees = 100, #一般来说先越大越好,然后选择合适的数目
interaction.depth = 1,
n.minobsinnode = 10,
shrinkage = 0.001, #0.001可能太慢了,可以先试试0.01
bag.fraction = 0.5,#从训练样本选择建立下一个树的比例,不用管
train.fraction = 1.0, #从训练样本选择样本数量建立模型的比例,不用管
cv.folds=0, #交叉验证的则数,可以用来提取最适的回归树数目
keep.data = TRUE, #不用管
verbose = FALSE, #如果为TRUE,将输出进度和性能指标
class.stratify.cv=NULL, #不用管
n.cores = NULL) #使用CPU核心的数量
重要参数理解
**(1)distribution:**模型计算损失函数时,需要对输出变量的数据分布做出假设。一般来说,对于分类问题,选择bernoulli或者adaboost,前者更为推荐;对于连续因变量,选择gaussian或者laplace。此外,gbm包还为一些具体问题提供了不少其它选择。
**(2)shrinkage:**学习速率,即每一步迭代中向梯度下降方向前进的速率。一般来说学习速率越小,模型表现越好。令shrinkage=0.001得出的模型几乎一定比shrinkage=0.01的模型好,然而代价是前者运算所需的时间和所耗内存将是后者的10倍。所以选择的准则是在可以承受的时间和内存范围内,shrinkage越小越好。
(3)n.trees:即number of iteration—迭代次数。迭代次数的选择与学习速率密切相关,下图展示了模型表现、学习速率和迭代次数之间的关系
迭代次数可以设得稍微大一点,因为模型训练完后,gbm中的gbm.perf可以估计出最佳迭代次数以供预测阶段使用。在模型训练阶段,gbm作者的经验法则是:3000-10000之间的迭代次数搭配0.01-0.001之间的学习速率。
(4)interaction.depth和n.minobsinnode:子决策树即基础学习器的深度和决策树叶节点包含的最小观测树,若基础学习器训练得过于复杂,将提升模型对于样本的拟合能力而导致过拟合问题,因此子决策树深度不宜过大,叶节点可包含的最小观测书不宜过小。
2.3 h2o.gbm
2.4 dismo
gbm.step(data, gbm.x, gbm.y, offset = NULL, fold.vector = NULL, tree.complexity = 1,
learning.rate = 0.01, bag.fraction = 0.75, site.weights = rep(1, nrow(data)),
var.monotone = rep(0, length(gbm.x)), n.folds = 10, prev.stratify = TRUE,
family = "bernoulli", n.trees = 50, step.size = n.trees, max.trees = 10000,
tolerance.method = "auto", tolerance = 0.001, plot.main = TRUE, plot.folds = FALSE,
verbose = TRUE, silent = FALSE, keep.fold.models = FALSE, keep.fold.vector = FALSE,
keep.fold.fit = FALSE, ...)
参数说明:
data : 输入数据框
gbm.x : 数据中预测变量的索引或名称
gbm.y : 数据中响应变量的索引或名称
offset : 抵消
fold.vector : 要读入以进行偏移交叉验证的折叠向量
tree.complexity : 设置单个树的复杂性
learning.rate : 设置应用于子树的权重,权重越小,说明需要树的数量就会越多,模型也越准确
bag.fraction : 设置用于选择变量的观察值的比例,默认为0.75
site.weights : 允许站点的不同权重
var.monotone : 将对单个预测值的响应限制为单调
n.folds : 折叠次数
prev.stratify : 患病率分层褶皱-仅用于存在/不存在数据
family : 族-伯努利(=二项式)、泊松、拉普拉斯或高斯
对于分类问题,选择bernoulli或者adaboost,前者更为推荐;对于连续因变量,选择gaussian或者laplace
n.trees : 要拟合的初始树数
step.size : 每个循环中要添加的树数
max.trees : 停止前要适应的最大树数
实例参考
#加载包
library(dismo)
#导入数据
data(Anguilla_train) #建模数据集
data(Anguilla_test) #测试数据集
LST_dismo <- gbm.step(data = Anguilla_train,
gbm.x = 3:13, gbm.y = 2, family = "gaussian",
tree.complexity = 5, learning.rate = 0.005, bag.fraction = 0.5)
LST_dismo$gbm.call$best.trees # 树的最佳数量950
#树数量越多,意味着模型预测能力相对越强。如果想要更多的树数量,则要将 gbm_step() 中的 learning.rate 改到更小:
summary(LST_dismo)#查看各变量的影响程度
#使用pROC包计算AUC值
library(pROC)
auc(test0$Spring_LST,preds)#Area under the
curve:
贡献因子
contributions<-LST_dismo$contributions
#查看各变量的影响程度$rel.inf
a<-summary(LST_dismo)
#然后计算不同类型的因子的贡献比例,画图
library(ggplot2)
ggplot(contributions)+
geom_col(aes(x=rel.inf,y=var,fill=var))+
xlab("Ralative Influence")+
ylab("")+
guides(fill=F)
#或者直接summary(angaus.tc5.lr01)
barplot(a$rel.inf,
names.arg = a$var, # x轴标签
main = "Relative influence", # 主标题
xlab = "factor", # x轴标签
ylab = "importance", # y轴标签
col = "skyblue", # 柱子颜色
ylim = c(0, 70) # y轴范围
)
ggplot(a, aes(x = a$rel.inf, y = a$var)) +
geom_bar(stat = "identity", fill = "skyblue") +
labs(title = "Relative influence",
x = "importance", # 将x轴标签设为Values
y = "factor") + # 将y轴标签设为Categories
coord_flip()+ # 使用coord_flip()函数进行轴互换
theme_classic()# 去除背景
3. 评价
1.求R方
使用模型预测值和实际观测值之间的平方和来计算总变异(Total Sum of Squares,TSS)和残差平方和(Residual Sum of Squares,RSS)。然后,R方可以通过以下公式计算:
R方 = 1 - (RSS / TSS)
其中,TSS表示总变异,RSS表示残差平方和。你可以使用模型的预测值和实际观测值来计算这些值。
#你还可以考虑使用caret包中的R2函数来计算R方值。R2函数可以计算R方值,它通常适用于各种回归模型,包括Gradient Boosting
library(caret)
r_squared <- R2(pred = predict(LST_dismo, newdata = test_data), obs = test_data$your_response_variable)
实例
library(ggplot2)
library(lattice)
library(caret)
install.packages("caret")
r_squared <- R2(pred=preds, obs = test0$Spring_LST)#0.3559662
2.RMSE
rmse <- sqrt(mean((preds-test0$Spring_LST)^2))#1.456581