sklearn与XGBoost

1 在学习XGBoost之前

1.1 xgboost库与XGB的sklearn API

陈天奇创造了XGBoost算法后，很快和一群机器学习爱好者建立了专门调用XGBoost库，名为xgboost。xgboost是一个独立的、开源的，并且专门提供梯度提升树以及XGBoost算法应用的算法库。它和sklearn类似，有一个详细的官方网站可以提供学习资料，并且可以与C、Python、R、Julia等语言连用，但需要单独安装和下载。
xgboost documents：https://xgboost.readthedocs.io/en/latest/index.html
本文基于Pythonl来运行。xgboost库要求必须提供适合的Scipy环境。接下来提供在windows中使用pip安装xgboost的代码：

pip install xgboost  #安装xgboost库
pip install --upgrade xgboost   #更新xgboost库

安装完毕之后，就能够使用这个库中所带的XGB相关的类了。

import xgboost as xgb

现在有两种方式来使用xgboost库。第一种方式，是直接使用xgboost库自己的建模流程。

其中最核心的是DMatrix这个读取数据的类，以及train()这个用于训练的类。与sklearn把所有的参数都写在类中的方式不同，xgboost库中必须先使用字典设定参数集，再使用train来将参数集输入，然后进行训练。会这样设计的原因，是因为XGB所涉及到的参数实在太多，全部写在xgb.train()太长也容易出错。下面给出params可能的取值以及xgboost.train的列表。
params {eta, gamma, max_depth, min_child_weight, max_delta_step, subsample, colsample_bytree, colsample_bylevel, colsample_bynode, lambda, alpha, tree_method_string, sketch_eps, scale_pos_weight, updater, refresh_leaf, process_type, grow_policy, max_leaves, max_bin, predictor, num_parallel_tree}
xgboost.train(params, dtrain, num_boost_round = 10, evals = (), obj = None, feval = None, maximize = False, early_stopping_rounds = None, evals_result = None, verbose_eval = True, xgb_model = None, callbacks = None, learning_rates = None)
或者，也可以选择第二种方法，使用xgboost库中的sklearn中的API。这就是说，可以调用如下的类，并用sklearn当中的惯例的实例化，fit和predict的流程来运行XGB，并且也可以调用属性（比如coef_）等等。当然，这是回归的类，也有用于分类、排序的类。他们与回归的类非常相似，因此了解一个类即可。
xgboost.XGBRegressor(max_depth = 3, learning_rate = 0.1, n_estimators = 100, silent = True, objective = ‘reg:linear’, booster = ‘gbtree’, n_jobs = 1, nthread = None, gamma = 0, min_child_weight = 1, max_delta_step = 0, subsample = 1, colsample_bytree = 1, colsample_bylevel = 1, reg_alpha = 0, reg_lambda = 1, scale_pos_weight = 1, base_score = 0.5, random_state = 0, seed = None, missing = None, importance_type = ‘gain’, $\ast \ast$kwargs)
调用xgboost.train和调用sklearn API中的类XGBRegressor，需要输入的参数是不同的，而且看起来相当的不同。但其实，这些参数只是写法不同，功能是相同的。如params字典中的第一个参数eta，其实就是XGBRegressor里面的参数learning_rate，他们的含义和实现的功能是一模一样的。只不过在sklearn API中，开发团队友好的将参数的名称调节成了与sklearn中其他算法类更相似的样子。
所以，使用xgboost中设定的建模流程来建模，和使用sklearn API中的类来建模，模型效果是比较相似的，但是xgboost库本身的运算速度（尤其是交叉验证）以及调参手段比sklearn要简单。

1.2 XGBoost的三大板块

XGBoost本身的核心是基于梯度提升树实现的集成算法，整体来说可以有三个核心部分：集成算法本身，用于集成的弱评估器，以及应用中的其他过程。三个部分中，前两个部分包含了XGBoost的核心原理以及数学过程，最后的部分主要是在XGBoost应用中占有一席之地。

参数	集成算法	弱评估器	其他过程
n_estimators	$\ast$
learning_rate	$\ast$
silent(verbosity)	$\ast$
subsample	$\ast$
max_depth		$\ast$
objective		$\ast$
booster		$\ast$
gamma		$\ast$
min_child_weight		$\ast$
max_delta_step		$\ast$
colsample_bytree		$\ast$
colsample_bylevel		$\ast$
reg_alpha		$\ast$
reg_lambda		$\ast$
nthread			$\ast$
n_jobs			$\ast$
scale_pos_weight			$\ast$
base_score			$\ast$
seed			$\ast$
random_state			$\ast$
missing			$\ast$
importance_type			$\ast$

2 梯度提升树

2.1 提升集成算法：重要参数n_estimators

XGBoost的基础是梯度提升算法，因此必须先从了解梯度提升算法开始。梯度提升（Gradient boosting）算法是构建预测模型的最强大技术之一，它是集成算法中提升法（Boosting）的代表算法。集成算法通过在数据集上建立多个弱评估器，汇总所有弱评估器的预测模型，以获取比单个模型更好的回归或分类表现。弱评估器被要求至少比随机猜测更好的模型，即预测准确度不低于50%的任意模型。
集成不同弱评估器的方法有很多种，有一次性建立多个平行独立的弱评估器的装袋法（Bagging），也有逐一构建弱评估器，经过多次迭代逐渐累积多个弱评估器的提升法（Boosting）。提升法中最著名的算法包括Adaboost和梯度提升树，XGBoost就是由梯度提升树发展而来的。梯度提升树中可以有回归树也可以有分类树，两者都是以CART树算法作为主流，XGBoost背后也是CART树，这意味着XGBoost中所有的树都是二叉的。以梯度提升回归树为例，了解一下Boosting算法的工作原理。
首先梯度提升回归树是专注于回归的树模型的提升集成模型，其建模过程大致如下：最开始先建立一棵树，然后逐渐迭代，每次迭代过程中都增加一棵树，逐渐形成众多树模型集成的强评估器。

对于决策树而言，每个被放入模型的任意样本$i$最终一定都会落到一个叶子节点上。而对于回归树，每个叶子节点上的值是这个叶子节点上所有样本的均值。

对于梯度提升树回归树来说，每个样本的预测结果可以表示为所有树上的结果的加权求和：$${\hat{y}}_i^{(k)}=\sum _k^K{\gamma }_k{h}_k(x_i)$$
其中，$K$是树的总数量，$k$代表第$k$棵树，${\gamma }_k$是这棵树的权重，$h_k$表示这棵树上的预测结果。
值得注意的是，XGB作为GBDT的改进，在$\hat{y}$上却有所不同。对于XGB来说，每个叶子节点上会有一个预测分数（prediction score），也被称为叶子权重。这个叶子权重就是所有在这个叶子节点上的样本在这一棵树上的回归取值，用$f_k(x_i)$或者$w$来表示，其中$f_k$表示第$k$棵决策树，$x_i$表示样本$i$对应的特征向量。当只有一棵树的时候，$f_1(x_i)$就是提升集成算法返回的结果，但这个结果往往非常糟糕。当有多棵树的时候，集成模型的回归结果就是所有树的预测分数之和，假设这个集成模型中总共有$K$棵决策树，则整个模型在这个样本$i$上给出的预测结果为：$${\hat{y}}_i^{(k)}=\sum _k^K{f}_k(x_i)$$

XGB和GBDT核心区别1：求解预测值$\hat{y}$的方式不同
GBDT中预测值是由所有弱分类器上的预测结果的加权求和，其中每个样本上的预测结果就是样本所在的叶子节点的均值。而GBDT中的预测值是所有弱分类器上的叶子权重直接求和得到，计算叶子权重是一个复杂的过程。

从上式来看，在集成中需要考虑的第一件事是超参数$K$，究竟要建多少棵树？

参数含义	xgb.train()	xgb.XGBRegressor()
集成中弱评估器的数量	num_round，默认为10	n_estimators，默认为100
训练者是否打印每次训练的结果	silent，默认为False	silent，默认True

在最新版本的xgboost中没有参数silent，改为verbosity，可以取值0-3，其中取0时表示silent。
在随机森林中理解n_estimators：n_estimators越大，模型的学习能力就会越强，模型也越容易过拟合。在随机森林中，调整的第一个参数就是n_estimators，这个参数非常强大，常常能够一次性将模型调整到极限。在XGB中，也期待相似的表现，虽然XGB的集成方式与随机森林不同，但使用更多的弱分类器来增强模型整体的学习能力这件事是一致的。

#导入需要的库，模块以及数据
from xgboost import XGBRegressor as XGBR
from sklearn.ensemble import RandomForestRegressor as RFR
from sklearn.linear_model import LinearRegression as LinearR
from sklearn.datasets import load_boston
from sklearn.model_selection import KFold,cross_val_score as CVS,train_test_split as TTS
from sklearn.metrics import mean_squared_error as MSE
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from time import time
import datetime

data = load_boston()#波士顿数据集非常简单，但它涉及到的问题却很多
#data
x = data.data
y = data.target
x.shape
#结果：(506, 13)
y.shape
#结果：(506,)

#建模，查看其他接口和属性
xtrain,xtest,ytrain,ytest = TTS(x,y,test_size = 0.3, random_state = 420)
reg = XGBR(n_estimators = 100).fit(xtrain,ytrain)#训练
reg.predict(xtest)#传统接口predict
'''
结果：
array([ 6.6689262, 22.34918  , 31.052807 , 13.911595 ,  9.467966 ,
       22.658588 , 14.514282 , 15.092699 , 15.293644 , 12.680115 ,
       24.140797 , 35.890083 , 21.573483 , 27.07066  , 19.052658 ,
        9.89033  , 23.386076 , 23.588493 , 23.311466 , 22.401644 ,
       18.98444  , 15.766946 , 25.8352   , 20.193802 , 19.982517 ,
       15.611423 , 22.883228 , 29.838228 , 22.815304 , 16.779037 ,
       37.13194  , 20.133305 , 19.67352  , 23.525528 , 22.845137 ,
       23.87397  , 15.17887  , 23.45934  , 16.685331 , 31.761686 ,
       18.525843 , 22.441063 , 38.48728  , 17.93719  , 15.10122  ,
       28.980541 , 46.363487 , 12.842797 ,  9.618281 , 35.40579  ,
       25.657566 , 20.605602 , 20.800055 , 49.228447 , 31.355848 ,
       29.382515 , 18.911947 , 21.049877 , 16.165182 , 18.098577 ,
       14.659002 , 21.720213 , 19.413454 , 28.932102 , 30.573524 ,
       19.228426 , 20.531511 , 15.666289 , 23.52929  , 19.30554  ,
       28.384985 , 42.83562  , 29.429724 , 23.306015 , 19.741224 ,
       24.202463 , 38.735516 , 26.782232 , 22.168324 , 20.67139  ,
       19.534992 , 47.360317 , 24.008802 ,  8.184814 , 25.613977 ,
       20.691843 , 17.127483 , 21.10027  , 22.279167 ,  7.755469 ,
       20.0476   , 15.406112 , 42.38165  , 33.828186 , 22.788246 ,
        9.302419 , 10.416187 , 20.033014 ,  8.241165 , 12.816478 ,
       30.793974 , 10.078776 , 25.383692 , 21.933594 , 30.53568  ,
       42.866497 , 19.598145 ,  8.321976 , 23.194368 , 12.547767 ,
       46.838867 , 22.961082 , 20.244467 , 23.170694 , 18.948856 ,
       29.682056 , 24.363865 , 19.715958 , 44.975193 , 17.64582  ,
       24.3169   , 24.946495 , 18.23235  , 16.922577 , 14.77714  ,
       21.757038 , 33.83876  , 10.938419 , 20.035763 , 21.085218 ,
       19.331562 , 20.505526 ,  8.285518 , 28.01946  , 29.631227 ,
       19.908175 , 18.30325  , 14.204149 , 10.795732 , 18.232307 ,
       42.266888 , 17.304502 , 22.974121 , 20.946724 , 30.724297 ,
       20.072989 , 12.772859 , 41.472908 , 27.652838 , 22.20476  ,
       14.235871 , 25.88964  ], dtype=float32)
'''
reg.score(xtest,ytest)#R^2，shift+tab可以看函数具体原理
#结果：0.9050988968414799
MSE(ytest,reg.predict(xtest))
#结果：8.830916343629323
reg.feature_importances_#树模型的优势之一：能够查看模型的重要性分数，可以使用嵌入法(SelectFromModel)进行特征选择
'''
结果：
array([0.01902167, 0.0042109 , 0.01478316, 0.00553537, 0.02222196,
       0.37914088, 0.01679686, 0.0469872 , 0.04073574, 0.05491759,
       0.06684221, 0.00869464, 0.3201119 ], dtype=float32)
'''
#交叉验证，与线性回归&随机森林回归进行对比
reg = XGBR(n_estimators = 100)#交叉验证中导入没有经过训练的模型
CVS(reg,xtrain,ytrain,cv = 5).mean()#这里返回的模型评估指标是R^2
#严谨的交叉验证与不严谨的交叉验证之间的讨论：训练集 or 全数据？
#结果：0.7995062821902295
CVS(reg,xtrain,ytrain,cv = 5,scoring = "neg_mean_squared_error").mean()
#结果：-16.215644229762717
#来查看一下sklearn中所有的模型评估指标
import sklearn
sorted(sklearn.metrics.SCORERS.keys())
'''
结果：
['accuracy',
 'adjusted_mutual_info_score',
 'adjusted_rand_score',
 'average_precision',
 'balanced_accuracy',
 'completeness_score',
 'explained_variance',
 'f1',
 'f1_macro',
 'f1_micro',
 'f1_samples',
 'f1_weighted',
 'fowlkes_mallows_score',
 'homogeneity_score',
 'jaccard',
 'jaccard_macro',
 'jaccard_micro',
 'jaccard_samples',
 'jaccard_weighted',
 'max_error',
 'mutual_info_score',
 'neg_brier_score',
 'neg_log_loss',
 'neg_mean_absolute_error',
 'neg_mean_absolute_percentage_error',
 'neg_mean_gamma_deviance',
 'neg_mean_poisson_deviance',
 'neg_mean_squared_error',
 'neg_mean_squared_log_error',
 'neg_median_absolute_error',
 'neg_root_mean_squared_error',
 'normalized_mutual_info_score',
 'precision',
 'precision_macro',
 'precision_micro',
 'precision_samples',
 'precision_weighted',
 'r2',
 'rand_score',
 'recall',
 'recall_macro',
 'recall_micro',
 'recall_samples',
 'recall_weighted',
 'roc_auc',
 'roc_auc_ovo',
 'roc_auc_ovo_weighted',
 'roc_auc_ovr',
 'roc_auc_ovr_weighted',
 'top_k_accuracy',
 'v_measure_score']
'''
#使用随机森林和线性回归进行一个对比
rfr = RFR(n_estimators = 100)
CVS(rfr,xtrain,ytrain,cv = 5).mean()
#结果：0.7922215794897552
CVS(rfr,xtrain,ytrain,cv = 5,scoring = "neg_mean_squared_error").mean()
#结果：-16.71229234925553
lr = LinearR()
CVS(lr,xtrain,ytrain,cv = 5).mean()
#结果：0.6835070597278079
CVS(lr,xtrain,ytrain,cv = 5,scoring = "neg_mean_squared_error").mean()
#结果：-25.349507493648463
#如果开启参数silent：在数据巨大，预料到算法运行会非常缓慢的时候可以使用这个参数来监控模型的训练进度
#最新版本的sklearn中的XGBR的参数中没有silent这一参数，改为verbosity
#verbosity : Optional[int]        The degree of verbosity. Valid values are 0 (silent) - 3 (debug).
reg = XGBR(n_estimators = 10, verbosity = 0)
CVS(reg,xtrain,ytrain,cv = 5,scoring = 'neg_mean_squared_error').mean()
#结果为：-18.633733952333067

#定义绘制以训练样本数为横坐标的学习曲线的函数
def plot_learning_curve(estimator, title, x, y
                       ,ax = None  #选择子图
                       ,ylim = None  #设置纵坐标的取值范围
                       ,cv = None  #交叉验证
                       ,n_jobs = None  #设定索要使用的线程
                       ):
    from sklearn.model_selection import learning_curve
    import matplotlib.pyplot as plt
    import numpy as np
    train_sizes, train_scores, test_scores = learning_curve(estimator, x, y
                                                           ,shuffle = True
                                                           ,cv = cv
                                                           #,random_state = 420
                                                            #random_state固定时绘制的图像不会每次运行都变化
                                                           ,n_jobs = n_jobs)
    if ax == None:
        ax = plt.gca()
    else:
        ax = plt.figure()
    ax.set_title(title)
    if ylim is not None:
        ax.set_ylim(*ylim)
    ax.set_xlabel("Training examples")
    ax.set_ylabel("Score")
    ax.grid()#绘制网格，不是必须
    ax.plot(train_sizes, np.mean(train_scores, axis = 1),"o-", color = "r", label = "Training score")
    ax.plot(train_sizes, np.mean(test_scores,axis = 1), "o-", color = "g", label = "Test score")
    ax.legend(loc = "best")
    return ax
#使用学习曲线观察XGB在波士顿数据集上的潜力
cv = KFold(n_splits = 5, shuffle = True, random_state = 42)#交叉验证模式
plot_learning_curve(XGBR(n_estimators = 100, random_state = 420)
                   ,"XGB", xtrain,ytrain,ax = None,cv = cv)
plt.show()
#多次运行，观察结果，这是怎么造成的？——random_state
#在现在的状况下，如何看数据的潜力，还能调上去吗？

训练集上的表现展示了模型的学习能力，测试集上的表现展示了模型的泛化能力，通常模型在测试集上的表现不太可能超过训练集，因此希望测试集的学习曲线能够努力逼近训练集的学习曲线。观察三种学习曲线组合：希望将模型调整成什么样？能够将模型调整成什么样？

#使用参数学习曲线观察n_estimators对模型的影响
axisx = range(10,1010,50)
rs = []
for i in axisx:
    reg = XGBR(n_estimators = i, random_state = 420)
    rs.append(CVS(reg,xtrain,ytrain,cv = cv).mean())
print(axisx[rs.index(max(rs))],max(rs))
plt.figure(figsize = (20,5))
plt.plot(axisx,rs,c = "red",label = "XGB")
plt.legend()
plt.show()
#选出来的n_estimators非常不寻常，是否要选择准确率最高的n_estimators值呢？
#结果为：160 0.8320776498992342

进化的学习曲线：方差与泛化误差
在随机森林中有方差-偏差困境。在机器学习中，用来衡量模型在未知数据上的准确率的指标，叫做泛化误差（Genelization error）。一个集成模型(f)在未知数据集(D)上的泛化误差$E(f;D)$，由方差(var)，偏差(bias)和噪声(e)共同决定。其中偏差就是训练集上的拟合程度决定，方差是模型的稳定性决定，噪音是不可控的。而泛化误差越小，模型就越理想。
$E(f;D)=bia{s}^2+var+e^2$
过去往往直接取学习曲线获得的分数的最高点，即考虑偏差最小的点，是因为模型极度不稳定，方差很大的情况其实比较少见。但是现在的数据量非常少，模型会相对不稳定，因此应当方差也纳入考虑的范围。在绘制学习曲线时，不仅要考虑偏差的大小，还要考虑方差的大小，更要考虑泛化误差中可控的部分。当然，并不是说可控的部分比较小，整体的泛化误差就一定小，因为误差有时候可能占主导。基于这种思路改进学习曲线：

axisx = range(50,1050,50)
rs = []
var = []
ge = []
for i in axisx:
    reg = XGBR(n_estimators = i, random_state = 420)
    cvresult = CVS(reg, xtrain, ytrain, cv = cv)
    #记录  1-偏差
    rs.append(cvresult.mean())
    #记录  方差
    var.append(cvresult.var())
    #计算泛化误差的可控部分
    ge.append((1-cvresult.mean())**2 + cvresult.var())
#打印R2最高所对应的参数取值，并打印这个参数下的方差
print(axisx[rs.index(max(rs))], max(rs), var[rs.index(max(rs))])
#打印方差最低时对应的参数取值，并打印这个参数下的R2
print(axisx[var.index(min(var))], rs[var.index(min(var))], min(var))
#打印泛化误差可控部分的参数取值，并打印这个参数下的R2，方差以及泛化误差的可控部分
print(axisx[ge.index(min(ge))], rs[ge.index(min(ge))], var[ge.index(min(ge))], min(ge))
plt.figure(figsize = (20,5))
plt.plot(axisx, rs, c = "red", label = "XGB")
plt.legend()
plt.show()
'''
结果：
100 0.8320924293483107 0.005344212126112929
100 0.8320924293483107 0.005344212126112929
100 0.8320924293483107 0.005344212126112929 0.03353716440826495
'''

#细化学习曲线，找出最佳n_estimators
axisx = range(50,250,10)
rs = []
var = []
ge = []
for i in axisx:
    reg = XGBR(n_estimators = i, random_state = 420)
    cvresult = CVS(reg, xtrain, ytrain, cv = cv)
    rs.append(cvresult.mean())
    var.append(cvresult.var())
    ge.append((1-cvresult.mean())**2 + cvresult.var())
print(axisx[rs.index(max(rs))], max(rs), var[rs.index(max(rs))])
print(axisx[var.index(min(var))], rs[var.index(min(var))], min(var))
print(axisx[ge.index(min(ge))], rs[ge.index(min(ge))], var[ge.index(min(ge))], min(ge))
rs = np.array(rs)
var = np.array(var)*0.01
plt.figure(figsize = (20,5))
plt.plot(axisx,rs,c = "black",label = "XGB")
#添加方差线
plt.plot(axisx, rs+var, c = "red", linestyle = "-.")
plt.plot(axisx, rs-var, c = "red", linestyle = "-.")
plt.legend()
plt.show()
'''
结果：
100 0.8320924293483107 0.005344212126112929
100 0.8320924293483107 0.005344212126112929
100 0.8320924293483107 0.005344212126112929 0.03353716440826495
'''

#看看泛化误差的可控部分如何？
plt.figure(figsize = (20,5))
plt.plot(axisx, ge, c = "gray", linestyle = "-.")
plt.show()

#检测模型效果是否提高了？
time0 = time()
print(XGBR(n_estimators = 10, random_state = 420).fit(xtrain, ytrain).score(xtest,ytest))
print(time()-time0)
'''
结果：
0.8894915063570491
0.02552485466003418
'''
time0 = time()
print(XGBR(n_estimators = 100, random_state = 420).fit(xtrain, ytrain).score(xtest,ytest))
print(time()-time0)
'''
结果：
0.9050988968414799
0.17356395721435547
'''

从这个过程中观察n_estimators参数对模型的影响，可以得出以下结论：
首先，XGB中的树的数量决定了模型的学习能力，树的数量越多，模型的学习能力越强。只要XGB中树的数量足够了，即便只有很少的数据，模型也能够学到训练数据100%的信息，所以XGB也是天生过拟合的模型。但在这种情况下，模型会变得非常不稳定。
第二，XGB中树的数量很少的时候，对模型的影响较大，当树的数量已经很多的时候，对模型的影响比较小，只能有微弱的变化。当数据本身就处于过拟合的时候，再使用过多的树能达到的效果甚微，反而浪费计算资源。当唯一指标$R^2$或者准确率给出的n_estimators看起来不太可靠的时候，可以改造学习曲线来帮助选择参数。
第三，树的数量提升对模型的影响有极限，最开始，模型的表现会随着XGB的树的数量一起提升，但到达某个点之后，树的数量越多，模型的效果会逐步下降，这也说明了暴力增加n_estimators不一定有效果。
这些都和随机森林中的参数n_estimators表现出一致的状态。在随机森林中总是先调整n_estimators，当n_estimators的极限已达到，才考虑其他参数，但XGB中的状况明显更加复杂，当数据集不太寻常的时候会更加复杂，这是要给出的第一个超参数，因此还是建议优先调整n_estimators，一般都不会建议一个太大的数目，300以下为佳。

2.2 有放回随机抽样：重要参数subsample

确认了有多少棵树之后，来思考一个问题：建立了众多的树，怎么能够保证模型整体的效果变强呢？集成的目的是为了模型在样本上能表现出更好的效果，所以对于所有的提升集成算法，每构建一个评估器，集成模型的效果都会比之前更好。也就是随着迭代的进行，模型整体的效果必须要逐渐提升，最后要实现集成模型的效果最优，要实现这个目标，可以首先从训练数据上着手。
训练模型之前，必然会有一个巨大的数据集。都知道树模型是天生过拟合的模型，并且如果数据量太大，树模型的计算会非常缓慢。因此，要对原始数据集进行有放回抽样（bootstrap）。有放回的抽样每次只能抽取一个样本，若需要总共N个样本，就需要抽取N次。每次抽取一个样本的过程是独立的，这一次被抽到的样本会被放回数据集中，下一次还可能被抽到，因此抽出的数据集中，可能有一些重复的数据。

在无论是装袋还是提升的集成算法中，有放回抽样都是防止过拟合，让单一弱分类器变得更轻量的必要操作。实际应用中，每次抽取50%左右的数据就能够有不错的效果了。sklearn的随机森林类中也有名为bootstrap的参数来帮助控制这种随机有放回抽样。同时，这样做还可以保证集成算法中的每个弱分类器（每棵树）都是不同的模型，基于不同的数据集建立的自然是不同的模型，而集成一系列一模一样的弱分类器是没有意义的。
在梯度提升树中，每一次迭代都要建立一棵新的树，因此每次迭代中，都要有放回抽取一个新的训练样本。不过，这并不能保证每次建新树后，集成的效果都比之前要好。因此我们规定，这梯度提升树中，每构建一个评估器，都让模型更加集中于数据集中容易被判错的那些样本。

首先有一个巨大的数据集，在建第一棵树时，对数据进行初次有放回抽样，然后建模。建模完毕后，对模型进行一个评估，然后将模型预测错误的样本反馈给数据集，一次迭代就算完成。紧接着，建立第二棵决策树，于是开始进行第二次有放回抽样。但这次有放回抽样和初次的随机有放回抽样不同，在这次的抽样中，加大了被第一棵树判断错误的样本的权重。也就是说，被第一棵树判断错误的样本，更有可能被抽中。
基于这个有权者的训练集来建模，新建的决策树就会更加倾向于这些权重更大的，很容易被判错的样本。建模完毕之后，又将判错的样本反馈给原始数据集，下一次迭代的时候，被判错的样本的权重会更大，新的模型会更加倾向于很难被判断的这些样本，如此反复迭代，越后面建的树，越是之前的树们判错样本上的专家，越专注于攻克那些之前的树们不擅长的数据。对于一个样本而言，它被预测错误的次数越多，被加大权重的次数也就越多。相信只要弱分类器足够强大，随着模型整体不断在被判错的样本上发力，这些样本会渐渐被判断正确。如此就一定程度上实现了每新建一棵树模型的效果都会提升的目标。
在sklearn中，使用参数subsample来控制随机抽样。在xgb和sklearn中，这个参数都默认为1且不能取到0，这说明无法控制模型是否进行随机有放回抽样，只能控制抽样抽出来的样本量大概是多少。

参数含义	xgb.train()	xgb.XGBRegressor()
随机抽样的时候抽取的样本比例，范围(0,1]	subsample，默认为1	subsample，默认为1

那除了让模型更加集中于那些困难样本，采样还对模型造成了什么样的影响？采样会减少样本数量，而从学习曲线来看样本数量越少模型的过拟合会越严重，因为对模型来说，数据量越少模型学习越容易，学到的规则也会越具体越不适用于测试样本。所以subsample参数通常是在样本量本身很大的时候来调整和使用。
模型现在正处于样本量过少并且过拟合的状态，根据学习曲线展现出来的规律，训练样本量在200左右的时候，模型的效果有可能反而比更多训练数据的时候好，但这并不代表模型的泛化能力在更小的训练样本量下会更强。正常来说样本量越大，模型才不容易过拟合，现在展现出来的效果，是由于样本量太小造成的一个巧合。从这个角度来看，subsample参数对模型的影响应该会非常不稳定，大概率应该是无法提升模型的泛化能力的，但也不乏提升模型的可能性。依然使用波士顿房价数据集，来看学习曲线：

axisx = np.linspace(0,1,20)
rs = []
for i in axisx:
    reg = XGBR(n_estimators = 100,subsample = i, random_state = 420)
    rs.append(CVS(reg,xtrain,ytrain,cv = cv).mean())
print(axisx[rs.index(max(rs))],max(rs))
plt.figure(figsize = (20,5))
plt.plot(axisx,rs,c = "green", label = "XGB")
plt.legend()
plt.show()
#结果：1.0 0.8320924293483107

#继续细化学习曲线
axisx = np.linspace(0.05,1,20)
rs = []
var = []
ge = []
for i in axisx:
    reg = XGBR(n_estimators = 100, subsample = i, random_state = 420)
    cvresult = CVS(reg,xtrain,ytrain,cv = cv)
    rs.append(cvresult.mean())
    var.append(cvresult.var())
    ge.append((1-cvresult.mean())**2 + cvresult.var())
print(axisx[rs.index(max(rs))], max(rs), var[rs.index(max(rs))])
print(axisx[var.index(min(var))], rs[var.index(min(var))], min(var))
print(axisx[ge.index(min(ge))], rs[ge.index(min(ge))], var[ge.index(min(ge))], min(ge))
rs = np.array(rs)
var = np.array(var)
plt.figure(figsize = (20,5))
plt.plot(axisx,rs,c = "black", label = "XGB")
plt.plot(axisx,rs+var,c = "red",linestyle = "-.")
plt.plot(axisx,rs-var,c = "red",linestyle = "-.")
plt.legend()
plt.show()
'''
结果：
1.0 0.8320924293483107 0.005344212126112929
0.75 0.8173142984564995 0.0026778019326057657
1.0 0.8320924293483107 0.005344212126112929 0.03353716440826495
'''
#不要盲目找寻泛化误差可控部分的最低值，注意观察结果

#看看泛化误差的情况如何
reg = XGBR(n_estimators = 100
          ,subsample = 1.0
          ,random_state = 420).fit(xtrain,ytrain)
reg.score(xtest,ytest)
#结果：0.9050988968414799
MSE(ytest,reg.predict(xtest))
#结果：MSE(ytest,reg.predict(xtest))

参数选择的结果再预料之中，总体来说选择的这个参数还是使用全部训练数据集，这是由于数据量过少，降低抽样的比例反而会让数据的效果降低。

2.3 迭代决策树：重要参数eta

从数据的角度而言，让模型更加倾向于努力攻克那些难以判断的样本。但是，并不是说只要新建了一棵倾向于困难样本的决策树，它就能够帮助把困难样本判断正确。困难样本被加重权重是因为前面的树没有把它判断正确，所以对于下一棵树来说，它要判断的测试集的难度，是比之前的树所遇到的数据的难度都要高的，那要把这些样本都判断正确，会越来越难。如果新建的树在判断困难样本这件事上还没有前面的树做得好呢？如果新建的树刚好是一棵特别糟糕的树呢？所以，除了保证模型逐渐倾向于困难样本的方向，还必须控制新弱分类器的生成，必须保证每次新添加的树一定得是对这个新数据集预测效果最优的那一棵树。
思考：怎么保证每次新添加的树一定让集成学习的效果提升？
也许可以枚举？也许可以学习sklearn中的决策树构建树时一样随机生成固定数目的树，然后生成最好的那一棵？
平衡算法表现和运算速度是机器学习的艺术，希望能找出一种方法，直接帮我们求解出最优的集成算法结果。求解最优结果，能否把它转化成一个传统的最优化问题呢？
回顾一下最优化问题的老朋友——逻辑回归模型。在逻辑回归中，有方程：
$y(x)=\frac{1}{1+e^{-{\theta }^T\text{x}}}$
目标是求解让逻辑回归的拟合效果最优的参数组合$\theta$。首先找出了逻辑回归的损失函数$J(\theta )$，这个损失函数可以通过带入$\theta$来衡量逻辑回归在训练集上的拟合效果。然后，利用梯度下降来迭代$\theta$：
${\theta }_{k+1}={\theta }_k-\alpha \ast d_{ki}$
让第$k$次迭代中的${\theta }_k$减去通过步长和特征取值$x$计算出来的一个量，以此来得到$k+1$次迭代后的参数向量${\theta }_{k+1}$。可以让这个过程持续下去，知道找到能够让损失函数最小化的参数$\theta$为止。这是一个最典型的最优化过程。这个过程其实和现在希望做的事是相似的。

现在希望求解集成算法的最优结果，那应该可以使用同样的思路：首先找到一个损失函数$Obj$，这个损失函数应该可以通过代入预测结果${\hat{y}}_i$来衡量梯度提升树在样本的预测结果。然后，利用梯度下降来迭代集成算法：
${\hat{y}}_i^{(k+1)}={\hat{y}}_i^{(k)}+f_{k+1}(x_i)$
在$k$次迭代后，集成算法中总共有$k$棵树，而前面讲明了，$k$棵树的集成结果是前面所有树上的叶子权重的累加${\sum }_k^K{f}_k(x_i)$。所有让$k$棵树的集成结果${\hat{y}}_i^{(k)}$加上新建的树上的叶子权重$f_{k+1}(x_i)$，就可以得到第$k+1$次迭代后，总共$k+1$棵树的预测结果${\hat{y}}_i^{(k+1)}$了。让这个过程持续下去，直到找到能够让损失函数最小化的$\hat{y}$，这个$\hat{y}$就是模型的预测结果。参数可迭代，集成的树林也可以迭代。
但要注意，在逻辑回归中参数$\theta$迭代的时候减去的部分是人为规定的步长和梯度相乘的结果。而在GBDT和XGB中，却希望能够求解出让预测结果$\hat{y}$不断迭代的部分$f_{k+1}(x_i)$。但无论如何，现在已经有了最优化的思路了，只要顺着这个思路求解下去，必然能够在每一个数据集上找到最优的$\hat{y}$。
在逻辑回归中，自定义步长$\alpha$来干涉迭代速率，在XGB中看起来却没有这样的设置，但其实不然。在XGB中，完整的迭代决策树的公式应该写作：
${\hat{y}}_i^{(k+1)}={\hat{y}}_i^{(k)}+\eta f_{k+1}(x_i)$
其中$\eta$读作“eta”，是迭代决策树时的步长（shrinkage），又叫做学习率（learning rate）。和逻辑回归中的$\alpha$类似，$\eta$越大，迭代的速度越快，算法的极限很快就被达到，有可能无法收敛到真正的最佳。$\eta$越小，越有可能找到最精确的最佳值，更多的空间被留给了后面建立的树，但迭代的速度会比较缓慢。

在sklearn中，使用参数learning_rate来干涉学习速率：

参数含义	xgb.train()	xgb.XGBRegressor()
集成中的学习率，又称为步长，以控制迭代速率，常用于防止过拟合	eta，默认0.3，取值范围(0,1]	learning_rate，默认0.1，取值范围(0,1]

探索一下参数eta的性质：

#首先先来定义一个评分函数，这个评分函数能够帮助我们直接打印xtrain上的交叉验证结果
def regassess(reg,xtrain,ytrain,cv,scoring = ["r2"],show = True):
    score = []
    for i in range(len(scoring)):
        if show:
            print("{}:{:.2f}".format(scoring[i]  #模型评估指标的名字
                                    ,CVS(reg
                                        ,xtrain,ytrain
                                        ,cv = cv, scoring = scoring[i]).mean()))
        score.append(CVS(reg,xtrain,ytrain,cv = cv,scoring = scoring[i]).mean())
    return score
#运行一下函数来看看效果
regassess(reg,xtrain,ytrain,cv,scoring = ["r2","neg_mean_squared_error"])
'''
结果：
r2:0.83
neg_mean_squared_error:-12.11
[0.8320924293483107, -12.108495815458545]
'''
#关闭打印功能试试看
regassess(reg,xtrain,ytrain,cv,scoring = ["r2","neg_mean_squared_error"],show = False)
#结果：[0.8320924293483107, -12.108495815458545]
#观察一下eta如何影响模型
from time import time
import datetime
for i in [0,0.2,0.5,1]:
    time0 = time()
    reg = XGBR(n_estimators = 100, random_state = 420, learning_rate = i)
    print("learning_rate = {}".format(i))
    regassess(reg,xtrain,ytrain,cv,scoring = ["r2","neg_mean_squared_error"])
    print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))
    print("\t")
'''
结果：
learning_rate = 0
r2:-6.76
neg_mean_squared_error:-567.55
00:02:518612
	
learning_rate = 0.2
r2:0.83
neg_mean_squared_error:-12.30
00:03:455323
	
learning_rate = 0.5
r2:0.82
neg_mean_squared_error:-12.48
00:03:074720
	
learning_rate = 1
r2:0.71
neg_mean_squared_error:-20.06
00:02:524670
'''

除了运行时间，步长还是一个对模型效果影响巨大的参数，如果设置太大模型就无法收敛（可能导致$R^2$很小或者MSE很大的情况），如果设置太小模型速度就会非常缓慢，但最后究竟会收敛到何处很难由经验来判定，在训练集上表现出来的模样和在测试集上相差甚远，很难直接探索出一个泛化误差很低的步长。

axisx = np.arange(0.05,1,0.05)
rs = []
te = []
for i in axisx:
    reg = XGBR(n_estimators = 100,random_state = 420,learning_rate = i)
    score = regassess(reg,xtrain,ytrain,cv,scoring = ["r2","neg_mean_squared_error"],show = False)
    test = reg.fit(xtrain,ytrain).score(xtest,ytest)
    rs.append(score[0])
    te.append(test)
print(axisx[rs.index(max(rs))],max(rs))
plt.figure(figsize = (20,5))
plt.plot(axisx,te,c = "gray",label = "test")
plt.plot(axisx,rs,c = "green",label = "train")
plt.legend()
plt.show()
#结果：0.1 0.8347040420775244

虽然从图上来说，默认的0.1看起来是一个比较理想的情况，并且看起来更小的步长更利于现在的数据，但是无法确定对于其他数据会有怎么样的效果。所以通常不调整$\eta$，即便调整，一般也会在[0.01,0.2]之间变动。如果希望模型的效果更好，更多的可能是从树本身的角度来说，对树进行剪枝，而不会寄希望于调整$\eta$。
梯度提升树是XGB的基础，已经介绍了XGB中与梯度提升树的过程相关的四个参数：n_estimators，learning_rate，silent（verbosity）,subsample。这四个参数的主要目的，其实并不是提升模型表现，更多是了解梯度提升树的原理。现在来看，梯度提升树可以说是由三个重要的部分组成：

1. 一个能够衡量集成算法效果的，能够被最优化的损失函数$Obj$
2. 一个能够实现预测的弱评估器$f_k(x)$
3. 一种能够让弱评估器集成的手段，包括讲解的迭代方法，抽样手段，样本加权等等过程
   XGBoost是在梯度提升树的这三个核心要素上运行，它重新定义了损失函数和弱评估器，并且对提升算法的集成手段进行了改进，实现了运算速度和模型效果的高度平衡。并且，XGBoost将原本的梯度提升树拓展开来，让XGBoost不再是单纯的树的集成模型，也不只是单单的回归模型。只要调节参数，可以选择任何希望集成的算法，以及任何希望实现的功能。

3 XGBoost的智慧

xgboost.XGBRegressor(max_depth = 3, learning_rate = 0.1, n_estimators = 100, silent = True, objective = ‘reg:linear’, booster = ‘gbtree’, n_jobs = 1, nthread = None, gamma = 0, min_child_weight = 1, max_delta_step = 0, subsample = 1, colsample_bytree = 1, colsample_bylevel = 1, reg_alpha = 0, reg_lambda = 1, scale_pos_weight = 1, base_score = 0.5, random_state = 0, seed = None, missing = None, importance_type = ‘gain’, **kwargs)

*args -> 指的是普通的参数
** kwargs -> key-value 参数，指的是那种字典参数

3.1 选择弱评估器：重要参数booster

梯度提升算法中不只有梯度提升树，XGB作为梯度提升算法的进化，自然也不只有树模型一种弱评估器。在XGB中，除了树模型，还可以选用线性模型，如逻辑回归、线性回归，来进行集成。虽然主流的XGB依然是树模型，但必须了解可以使用其他的模型。基于XGB的这种性质，有参数"booster"来控制究竟使用怎样的弱评估器。

xgb.train() & params	xgb.XGBRegressor()
xgb_model	booster
使用哪种弱评估器。可以输入gbtree、gblinear或dart。输入的评估器不同，使用的params参数也不同，每种评估器都有自己的params列表。评估器必须与params参数相匹配，否则报错。	使用哪种弱评估器。可以输入gbtree、gblinear或dart。gbtree代表梯度提升树，dart是Dropouts meet Multiple Additive Regression Trees，可译为抛弃提升树，在建树的过程中会抛弃一部分树，比梯度提升树有更好的防过拟合功能。输入gblinear使用线性模型。

两个参数都默认为“gbtree”，如果不想使用树模型，则可以自行调整。当XGB使用线性模型的时候，它的许多数学过程就与使用普通的Boosting集成非常相似，因此以使用更多，也更加核心的基于树模型的XGBoost为主。
简单看看现有的数据集上，是什么样的弱评估器表现更好：

for booster in ["gbtree","gblinear","dart"]:
    reg = XGBR(n_estimators = 100
              ,learning_rate = 0.1
              ,random_state = 420
              ,booster = booster).fit(xtrain,ytrain)
    print(booster)
    print(reg.score(xtest,ytest))
'''
结果：
gbtree
0.9257336963525182
gblinear
0.5988666782512837
dart
0.9257337087151665
'''

3.2 XGB的目标函数：重要参数objective

梯度提升算法中都存在着损失函数，不同于逻辑回归和SVM等算法中固定的损失函数写法，集成算法中的损失函数是可选的，要选用什么损失函数取决于希望解决什么问题，以及希望使用怎样的模型。比如，如果目标是进行回归预测，那可以选择调节后的均方误差RMSE作为损失函数；如果进行分类预测，那可以选择错误率error或者对数损失log_loss。只要选出的函数是一个可微的，能够代表某种损失函数，它就可以是XGB中的损失函数。
在众多机器学习算法中，损失函数的核心是衡量模型的泛化能力，即模型在未知数据上预测的准确与否，训练模型的核心目标也是希望模型能够预测准确。在XGB中，预测准确自然是非常重要的因素，但之前提到过，XGB是实现了模型表现和运算速度的平衡的算法。普通的损失函数，如错误率、均方误差等等，都只能够衡量模型的表现，无法衡量模型的运算速度。在许多模型中，使用空间复杂度和时间复杂度来衡量模型的运算速率。XGB因此引入了模型复杂度来衡量算法的运算效率。因此XGB的目标函数被写作：传统损失函数+模型复杂度。$$Obj=\sum _{i=1}^{m}l(y_i,{\hat{y}}_i)+\sum _{k=1}^K\Omega (f_k)$$其中$i$代表数据集中的第$i$个样本，$m$表示导入第$k$棵树的数据总量，$K$代表建立的所有树（n_estimators），当只建立了$t$棵树的时候，式子应当为${\sum }_{k=1}^t\Omega (f_k)$。第一项代表传统的损失函数，衡量真实标签$y_i$与预测值${\hat{y}}_i$之间的差异，通常是RMSE，调节后的均方误差。第二项代表模型的复杂度，使用树模型的某种变换$\Omega$表示，这个变化代表了一个从树的结构来衡量树模型的复杂度的式子，可以有多种定义。注意，第二项中没有特征矩阵$x_i$的介入。在迭代每一棵树的过程中，都最小化$Obj$来力求获取最优的$\hat{y}$，因此同时最小化了模型的错误率和模型的复杂度，这种设计目标函数的方法不得不说实在是非常巧妙和聪明。

目标函数：可能的困惑
与其他算法一样，最小化目标函数以求模型效果最佳，并且可以通过限制n_estimators来限制迭代次数，因此可以看到生成的每棵树所对应的目标函数取值。目标函数中的第二项看起来是与$K$棵树都相关，但第一个式子看起来却只和样本量相关，仿佛只与当前迭代到的这棵树相关，看起来很奇怪。

其实并非如此，第一项传统损失函数也是与已经建好的所有树相关的，相关在：$${\hat{y}}_i^{(k)}=\sum _k^K{f}_k(x_i)=\sum _k^K{f}_{k-1}(x_i)+f_k(x_i)，K=k-1$$其中${\hat{y}}_i$已经包含了所有树的迭代结果，因此整个目标函数都与$K$棵树相关。
还可以从另一个角度去理解目标函数。回忆在随机森林中的方差-偏差困境。在机器学习中，用来衡量模型在未知数据上的准确率的指标，叫做泛化误差（Genelization error）。一个集成模型(f)在未知数据集(D)上的泛化误差$E(f;D)$，由方差(var)、偏差(bias)和噪声(e)共同决定，而泛化误差越小，模型就越理想。从下图可以看出，方差和偏差是此消彼长的，并且模型的复杂度越高，方差越大，偏差越小。

方差可以被简单地解释为模型在不同数据集上表现出来的稳定性，而偏差是模型预测的准确度。那方差-偏差困境就可以对应到$Obj$中了：$$Obj=\sum _{i=1}^{m}l(y_i,{\hat{y}}_i)+\sum _{k=1}^K\Omega (f_k)$$第一项是衡量偏差，模型越不准确，第一项就会越大。第二项是衡量方差，模型越复杂，模型的学习就会越具体，到不同数据集上的表现就会差异巨大，方差就会越大。所以求解$Obj$的最小值，其实是在求解方差与偏差的平衡点，以求模型的泛化误差最小，运行速度最快。树模型和树的集成模型都是学习天才，是天生过拟合的模型，因此大多数树模型最初都会出现在图像的右上方，必须通过剪枝来控制模型不要过拟合。现在XGBoost的损失函数中自带限制方差变大的部分，也就是说XGBoost会比其他的树模型更加聪明，不会轻易落到图像的右上方。可见，这个模型在设计的时候的确考虑了方方面面，无关XGBoost会如此强大。
在应用中，使用参数“objective”来确定目标函数的第一部分中的$l(y_i,{\hat{y}}_i)$，也就是衡量损失的部分。

xgb.train()	xgb.XGBRegressor()	xgb.XGBClassifier()
obj：默认binary:logistic	objective：默认reg:linear	objective：默认binary:logistic

常用的选择有：

输入	选用的损失函数
reg:linar	使用线性回归的损失函数，均方误差，回归时使用
binary:logistic	使用逻辑回归的损失函数，对数损失log_loss，二分类时使用
binary:hinge	使用支持向量机的损失函数，Hinge Loss，二分类时使用
multi:softmax	使用softmax损失函数，多分类时使用

还可以选择自定义损失函数，比如，可以选择输入平方损失$l(y_j,{\hat{y}}_j)=(y_j-{\hat{y}}_j{)}^2$，此时XGBoost其实就是算法梯度提升机器（gradient boosted machine）。在xgboost中，被允许自定义损失函数，但通常还是使用类已经设置好的损失函数。回归类中本来使用的就是reg:linear，因此在这里无需做任何调整。注意：分类型的目标函数导入回归类中会直接报错。
讲到这个参数，可以试试看xgb自身的调佣方式：

由于xgb中所有的参数都需要自己输入，并且objective参数的默认值是二分类，因此必须手动调节。试试看在其他参数相同的情况下，xgboost库本身和sklearn比起来，效果如何：

#默认reg:linear
reg = XGBR(n_estimators = 100,random_state = 420).fit(xtrain,ytrain)
reg.score(xtest,ytest)
MSE(ytest,reg.predict(xtest))
#结果：8.830916343629323

#xgb实现法
import xgboost as xgb
#使用类Dmatrix读取数据
dtrain = xgb.DMatrix(xtrain,ytrain)
dtest = xgb.DMatrix(xtest,ytest)
#非常遗憾无法打开来看看，所以通常都是先读到pandas里面查看之后再放到DMatrix中
dtrain
#结果：<xgboost.core.DMatrix at 0x2dcddb1b4f0>
#写明函数
param = {"verbosity":0, "objective":"reg:linear", "eta":0.1}
num_round = 100  #n_estimators
#类train,可以直接导入的参数是训练数据，树的数量，其他参数都需要通过params来导入
bst = xgb.train(param, dtrain, num_round)
#接口predict
preds = bst.predict(dtest)
from sklearn.metrics import r2_score
r2_score(ytest,preds)
#结果：0.9257336963525182
MSE(ytest,preds)
#结果：6.910768082073709

看得出来，无论是从$R^2$还是从MSE的角度来看，都是xgb库本身表现更优秀，这也许是由于底层的代码是由不同的团队创造的缘故。随着样本量的逐渐上升，sklearnAPI中调用的结果与xgboost中直接训练的结果会比较相似，如果希望的话可以分别训练，然后选取泛化误差较小的库。如果可以的话，建议脱离sklearnAPI直接调用xgboost库，因为xgboost库本身的调参要方便许多。

3.3 求解XGB的目标函数

有了如下目标函数，来看看如何求解它。在逻辑回归和支持向量中，通常先将目标函数转化成一个容易求解的方式（如对偶），然后使用梯度下降或者SMO之类的数学方法来执行最优化过程。之前通过逻辑回归的迭代过程理解了在梯度提升树中树时如何迭代的，是否可以使用逻辑回归的参数求解方式来求解XGB的目标函数呢？$$Obj=\sum _{i=1}^{m}l(y_i,{\hat{y}}_i)+\sum _{k=1}^K\Omega (f_k)$$很遗憾，在XGB中无法使用梯度下降，原因是XGB的损失函数没有需要求解的参数。在传统梯度下降中迭代的是参数，而在XGB中迭代的是树，树$f_k$不是数字组成的向量，并且不受特征矩阵$x$的影响，尽管这个迭代过程可以被类比到梯度下降上，但真实的求解过程却是完全不同的。
在求解XGB的目标函数的过程中，考虑的是如何能够将目标函数转化为更简单的、与树的结构直接相关的写法，以此来建立树的结构与模型的效果（包括泛化能力与运行速度）之间的直接联系。也因为这种联系的存在，XGB的目标函数又被称为“结构分数”。$${\hat{y}}_i^{(t)}=\sum _k^t{f}_k(x_i)=\sum _k^{(}t-1)f_k(x_i)+f_t(x_i)={\hat{y}}_i^{(t-1)}+f_t(x_i)$$首先，先来进行第一步转换：

其中$g_i$和$h_i$分别是在损失函数$l(y_i^t,{\hat{y}}_i^{(t-1)})$上对${\hat{y}}_i^{(t-1)}$所求的一阶导数和二阶导数，它们被统称为每个样本的梯度统计量（gradient statistics）。在许多算法的解法推导中，求解导数都是为了让一阶导数等于0来求解极值，而现在求解导数只是为了配合泰勒展开中的形式，仅仅是简化公式的目的罢了。所以GBDT和XGB的区别之中，GBDT求一阶导数，XGB求二阶导数，这两个过程根本是不可类比的。XGB在求解极值为目标的求导中也是求解一阶导数。
**可能的困惑：泰勒展开好像不是长这样？**单独查看泰勒展开的数学公式，会看到这样的公式：
$f(x)=\frac{f(c)}{0!}+\frac{f^{\prime }(c)}{1!}(x-c)+\frac{f^{\prime \prime }(c)}{2!}(x-c{)}^2+\frac{f^{\prime \prime \prime }(c)}{3!}(x-c{)}^3+...$
其中$f^{\prime }(c)$表示$f(x)$上对$x$求导后，令$x$的值等于$c$所取得的值。其中有假设：$c$与$x$非常接近，$(x-c)$非常接近0，于是可以将式子改写成：
$f(x+x-c)\approx \frac{f(c)}{0!}+\frac{f^{\prime }(c)}{1!}(x-c)+\frac{f^{\prime \prime }(c)}{2!}(x-c{)}^2+\frac{f^{\prime \prime \prime }(c)}{3!}(x-c{)}^3+...$
只取前两项，约等于：
$\approx \frac{f(c)}{0!}+\frac{f^{\prime }(c)}{1!}(x-c)+\frac{f^{\prime \prime }(c)}{2!}(x-c{)}^2$
$\approx f(c)+f^{\prime }(c)(x-c)+\frac{f^{\prime \prime }(c)}{2}(x-c{)}^2$
令：$x_1=x$,$x_2=x-c$：
$f(x_1+x_2)\approx f(x_1)+f^{\prime }(x_1)x_2+\frac{f^{\prime \prime }(x_1)}{2}{x}_2^2$
$F({\hat{y}}_i^{(t-1)}+f_t(x_i))\approx F({\hat{y}}_i^{(t-1)})+g_i\ast f_t(x_i)+h_i\ast \frac{1}{2}(f_t(x_i){)}^2$
如刚才所说，$x-c$需要很小，与$x$相比起来越小越好，在式子中，需要很小的这部分就是$f_t(x_i)$。这其实很好理解，对于一个集成算法来说，每次增加的一棵树对模型的影响其实非常小，尤其是当有许多树的时候，如n_estimators = 500的情况，$f_t(x_i)$与$x$相比总是非常小的，因此这个条件可以被满足，泰勒展开可以被使用。如此，目标函数可以被顺利转化为：$$Obj=\sum _{i=1}^m[f_t(x_i)g_i+\frac{1}{2}(f_t(x_i){)}^2{h}_i]+\Omega (f_t)$$这个式子中，$g_i$和$h_I$只与传统损失函数相关，核心的部分时需要决定的树$f_t$。接下来，研究一下$f_t$。

3.4 参数化决策树$f_k(x)$：参数alpha、lambda

xgboost.XGBRegressor(max_depth = 3, learning_rate = 0.1, n_estimators = 100, silent = True, objective = ‘reg:linear’, booster = ‘gbtree’, n_jobs = 1, nthread = None, gamma = 0, min_child_weight = 1, max_delta_step = 0, subsample = 1, colsample_bytree = 1, colsample_bylevel = 1, reg_alpha = 0, reg_lambda = 1, scale_pos_weight = 1, base_score = 0.5, random_state = 0, seed = None, missing = None, importance_type = ‘gain’, **kwargs)
在参数化决策树之前，先来简单复习一下回归树的原理。对于决策树而言，每个被放入模型的任意样本$i$最终都会落到一个叶子节点上。对于回归树，通常来说每个叶子节点上的预测值是这个叶子节点上所有样本的标签的均值。但值得注意的是，XGB作为普通回归树的改进算法，在$\hat{y}$上却有所不同。
对于XGB来说，每个叶子节点上会有一个预测分数（prediction score），也被称为叶子权重。这个叶子权重就是所有在这个叶子节点上的样本在这一棵树上的回归取值，用$f_k(x_i)$或者$w$来表示。

当有多棵树的时候，集成模型的回归结果就是所有树的预测分数之和，假设这个集成模型中总共有$K$棵决策树，则整个模型在这个样本$i$上给出的预测结果为：$${\hat{y}}_i^{(k)}=\sum _k^K{f}_k(x_i)$$
记者这个理解，来考虑每一棵树。对每一棵树，它都有自己独特的结构，这个结构即是叶子节点的数量、树的深度、叶子的位置等等所形成的一个可以定义唯一模型的树结构。在这个结构中，使用$q(x_i)$表示样本$x$所在的叶子节点，并且使用$w_{q(x_i)}$来表示这个样本落在第$k$棵树上的第$q(x_i)$个叶子节点中所获得的分数，于是有：
$f_k(x_i)=w_{q(x_i)}$
**这是对于每一个样本而言的叶子权重，然后在一个叶子节点上所有样本所对应的叶子权重是相同的。**设一棵树上总共包含了$T$个叶子节点，其中每个叶子节点的索引为$j$，则这个叶子节点上的样本权重是$w_j$。依据这个，定义模型的复杂度$\Omega (f)$为（注意这不是唯一可能的定义，当然还可以使用其他的定义，只要满足叶子越多/深度越大，复杂度越大的理论，可以自己决定$\Omega (f)$是一个怎样的式子）：
$\Omega (f)=\gamma T+正则化(Regularization)$
如果使用$L2$正则项：
$=\gamma T+\frac{1}{2}\lambda \mid \mid w\mid {\mid }^2=\gamma T+\frac{1}{2}\lambda {\sum }_{j=1}^T{w}_j^2$
如果使用$L1$正则项：
$=\gamma T+\frac{1}{2}\alpha \mid w\mid =\gamma T+\frac{1}{2}\alpha {\sum }_{j=1}^T\mid w_j\mid$
还可以两个一起使用：
$=\gamma T+\frac{1}{2}\lambda {\sum }_{j=1}^T{w}_j^2+\frac{1}{2}\alpha {\sum }_{j=1}^T\mid w_j\mid$
这个结构中有两部分内容，一部分是控制结构的$\gamma T$，另一部分则是正则项。叶子数量$T$可以代表整个树结构，这是因为在XGBoost中所有的树都是CART树（二叉树），所以可以根据叶子的数量$T$判断出树的深度，而$\gamma$是自定的控制叶子数量的参数。
至于第二部分正则化，类比一下岭回归和Lasso的结构，参数$\alpha$和$\lambda$的作用其实非常容易理解，它们都是控制正则化强度的参数，可以二选一使用，也可以一起使用加大正则化的力度。当$\lambda$和$\alpha$都为0的时候，目标函数就是普通的梯度提升树的目标函数。

XGB vs GBDT核心区别2：正则项的存在
在普通的梯度提升树GBDT中，是不在目标函数中使用正则项的。但GBDT借用正则项来修正树模型天生容易过拟合这个缺陷，在剪枝之前让模型能够尽量不过拟合。

来看正则化系数分别对应的参数：

参数含义	xgb.train()	xgb.XGBRegressor()
L1正则化的参数$\alpha$	alpha，默认为0，取值范围为$[0,+\infty ]$	reg_alpha，默认为0，取值范围为$[0,+\infty ]$
L2正则化的参数$\lambda$	lambda，默认为1，取值范围为$[0,+\infty ]$	reg_lambda，默认为1，取值范围为$[0,+\infty ]$

根据以往的经验，往往认为两种正则化达到的效果是相似的，只不过细节不同。如在逻辑回归中，两种正则化都会压缩$\theta$参数的大小，只不过L1正则化会让$\theta$为0，而L2正则化不会。在XGB中也是如此，当$\lambda$和$\alpha$越大，惩罚越重，正则项所占的比例就越大，在尽全力最小化目标函数的最优化方向下，叶子节点数量就会被压制，模型的复杂度就越来越低，所以对于天生过拟合的XGB来说，正则化可以一定程度上提升模型效果。
对于两种正则化如何选择的问题，从XGB的默认参数来看，优先选择的是L2正则化。当然，如果想尝试L1也不是不可以。两种正则项还可以交互，因此这两个参数的使用其实比较复杂。在实际应用中，正则化参数往往不是我们调参的最优选择，如果真的希望控制模型复杂度，会调整$\gamma$而不是调整这两个正则化参数，因此不必过于在意这两个参数最终如何影响了模型效果。对于树模型来说，还是剪枝参数地位更高更优先。只需要理解这两个参数从数学层面上如何影响模型就足够了。如果希望调整$\lambda$和$\alpha$，往往会使用网格搜索。

#使用网格搜索来查找最佳的参数组合
from sklearn.model_selection import GridSearchCV
param = {"reg_alpha":np.arange(0,5,0.05),"reg_lambda":np.arange(0,2,0.05)}
gscv = GridSearchCV(reg,param_grid = param,scoring = "neg_mean_squared_error",cv = cv)

time0 = time()
gscv.fit(xtrain,ytrain)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))
#结果：30:04:877986
gscv.best_params_
#结果：{'reg_alpha': 1.5, 'reg_lambda': 1.3}
gscv.best_score_
#结果：-11.109242037961907
preds = gscv.predict(xtest)
from sklearn.metrics import r2_score,mean_squared_error as MSE
r2_score(ytest,preds)
#结果：0.9130653431705259
MSE(ytest,preds)
#结果：8.089607562736592

3.5 寻找最佳树结构：求解w与T

上一节中定义了树和树的复杂度的表达式，树使用叶子节点上的预测分数来表达，而树的复杂度则是叶子数目加上正则项：
$f_t(x_i)=w_{q(x_i)}$，$\Omega (f_t)=\gamma T+\frac{1}{2}\lambda {\sum }_{j=1}^T{w}_j^2$
假设现在第$t$棵树的结构已经被确定为$q$，可以将树的结构代入损失函数，来继续转化目标函数。转化目标函数的目的是：建立树的结构（叶子节点的数量）与目标函数的大小之间的直接联系，以求出在第$t$次迭代中需要求解的最优的树$f_t$。注意，假设使用的是L2正则化（这也是参数lambda和alpha的默认设置，lambda为1，alpha为0），接下来的推导也会根据L2正则化来进行。$$\sum _{i=1}^m[f_t(x_i)g_i+\frac{1}{2}{f}_t(x_i{)}^2{h}_i]+\Omega (f_t)$$$$=\sum _{i=1}^m[w_{q(x_i)}{g}_i+\frac{1}{2}{w}_{q(x_i)}^2{h}_i]+\gamma T+\frac{1}{2}\lambda \sum _{j=1}^T{w}_j^2$$$$=\sum _{i=1}^m{w}_{q(x_i)}{g}_i+\sum _{i=1}^m\frac{1}{2}{w}_{q(x_i)}^2{h}_i+\gamma T+\frac{1}{2}\lambda \sum _{j=1}^T{w}_j^2$$其实，每片叶子上的$w_j$是一致的，唯一不同的是每个样本所对应的$g_i$。所有的样本必然会被分到$T$片叶子节点中的任一节点上，定义索引为$j$的叶子上含有的样本的集合是$I_j$。

于是可以有：$$\sum _{i=1}^m{w}_{q(x_i)}\ast g_i=w_{q(x_1)}\ast g_1+w_{q(x_2)}\ast g_2+w_{q(x_3)}\ast g_3=w_1(g_1+g_2)+w_2\ast g_3=\sum _{j=1}^T(w_j\sum _{i\in I_j}{g}_i)$$如此实现了这个转化。
$$=\sum _{j=1}^T(w_j\ast \sum _{i\in I_j}{g}_i)+\frac{1}{2}\sum _{j=1}^T(w_j^2\ast \sum _{i\in I_j}{h}_i)+\gamma T+\frac{1}{2}\lambda \sum _{j=1}^T{w}_j^2$$$$=\sum _{j=1}^T[w_j\sum _{i\in I_j}{g}_i+\frac{1}{2}{w}_j^2(\sum _{i\in I_j}{h}_i+\lambda )]+\gamma T$$
对于最终的式子，定义：$G_j={\sum }_{i\in I_j}{g}_i$，$H_j={\sum }_{i\in I_j}{h}_i$
于是可以有：$$Ob{j}^{(t)}=\sum _{j=1}^T[w_j{G}_j+\frac{1}{2}{w}_j^2(H_j+\lambda )]+\gamma T$$$$F^{\ast }(w_j)=w_j{G}_j+\frac{1}{2}{w}_j^2(H_j+\lambda )$$其中每个$j$取值下都是一个以$w_j$为自变量的二次函数$F^{\ast }$，目标是追求让$Obj$最小，只要单独的每一个叶子$j$取值下的二次函数都最小，那他们的加权和必然也会最小。于是，在$F^{\ast }$上对$w_j$求导，让一阶导数等于0以求极值，可得：
$\frac{\partial F^{\ast }(w_j)}{\partial w_j}=G_j+w_j(H_j+\lambda )$
$0=G_j+w_j(H_j+\lambda )$
$w_j=-\frac{G_j}{H_j+\lambda }$
把这个公式代入目标函数中，则有：$$Ob{j}^{(t)}=\sum _{j=1}^T[-\frac{G_j}{H_j+\lambda }\ast G_j+\frac{1}{2}(-\frac{G_j}{H_j+\lambda }{)}^2(H_j+\lambda )]+\gamma T$$$$=\sum _{j=1}^T[-\frac{G_j^2}{H_j+\lambda }+\frac{1}{2}\ast \frac{G_j^2}{H_j+\lambda }]+\gamma T$$$$=-\frac{1}{2}\sum _{j=1}^T\frac{G_j^2}{H_j+\lambda }+\gamma T$$
注意到，比起最初的损失函数+复杂度的样子，目标函数已经发生了巨大变化。样本量$i$已经被归结到每个叶子当中去，目标函数是基于每个叶子节点，也就是树的结构来计算。所以，目标函数又叫做“结构分数（structure score）”，分数月底，树整体的结构越好。如此，就建立了树的结构（叶子）和模型效果的直接联系。
更具体一点，来看一个例子：

所以在XGB的运行过程中，会根据Obj的表达式直接探索最好的结构，也就是说找寻最佳的树。从式子中可以看出，$\lambda$和$\gamma$是设定好的超参数，$G_j$和$H_j$是由损失函数和这个特征结构下树的预测结果${\hat{y}}_i^{(t-1)}$共同决定，而$T$只由树结构决定。则通过最小化Obj所求解出的其实是$T$，叶子数量。所以本质也就是求解树的结构了。
在这个算式下，可以有一种思路，那就是枚举所有可能的树结构$q$，然后一个个计算$Obj$，待选定了最佳的树结构（最佳的$T$）之后，使用这种树结构下计算出来的$G_j$和$H_j$就可以求解出每个叶子上的权重$w_j$，如此就找到最佳树结构，完成这次迭代。
可能的困惑：求解$w_j$的一些细节

1. 用$w$求解$w$？
   一个可能会感到困惑的点是，$G_j$和$H_j$的本质其实是损失函数上的一阶导数$g_i$和二阶导数$h_i$之和，而一阶和二阶导数本质是：
   $g_i=\frac{\partial l(y_i,{\hat{y}}_i^{(t-1)})}{\partial {\hat{y}}_i^{(t-1)}}$，$h_i=\frac{{\partial }^{2}l(y_i,{\hat{y}}_i^{(t-1)})}{\partial ({\hat{y}}_i^{(t-1)}{)}^2}$
   $y_i$是已知的标签，但有预测值的求解公式：
   ${\hat{y}}_i^{(k)}={\sum }_k^K{f}_k(x_i)={\sum }_k^K{w}_{q(x_i)}$
   这其实是说，$G_j$和$H_j$的计算中带有$w$，那先确定最佳的$T$，再求出$G_j$和$H_j$，结合$\lambda$求出叶子权重$w_j$的思路不觉得有些问题吗？仿佛在代入$w$求解$w$?对于这样的疑惑，要注意${\hat{y}}_i^{(t-1)}$与现在要求解的$w_j$其实不是在同一棵树上的。别忘记一直在迭代，现在求解的$w_j$是第$t$棵树上的结果，而${\hat{y}}_i^{(t-1)}$是前面$(t-1)$棵树的累计$w$，是在前面所有的迭代中已经求解出来的已知的部分。
2. 求解第一棵树时，没有“前面已经迭代完毕的部分”，怎么办？
   在建立第一棵树时，并不存在任何前面的迭代已经计算出来的信息，但根据公式需要使用如下式子来求解$f_1(x_i)$，并且在求解过程中还需要对前面所有的树的结果进行求导。
   ${\hat{y}}_i^{(1)}={\hat{y}}_i^{(0)}+f_1(x_i)$
   这时候，假设${\hat{y}}_i^{(0)}=0$来解决问题。事实是，由于前面没有已经测出来的树的结果，整个集成算法的结果现在也的确为0。
3. 第0棵树的预测值假设为0，那求解第一棵树的$g_i$和$h_i$的过程是在对0求导？
   这个问题其实很简单。在进行求导时，所有的求导过程都和之前推导的过程相一致，之所以能够这么做，是因为其实不是在对0求导，而是对一个变量${\hat{y}}_i^{(t-1)}$求导。只是除了求导之外，还需要在求导后的结果中带入这个变量此时此刻的取值，而这个取值在求解第一棵树时刚好等于0而已。更具体地如下，对0求导，和对变量求导后变量取值为0的区别：
   对常数0进行求导：$\frac{\partial (x^2+x)}{\partial (0)}=0$
   对变量$x$进行求导，但变量$x$等于0：$\frac{\partial (x^2+x)}{\partial (x)}=2x+1=2\ast 0+1=1$
   这些细节都理解了之后，对于先求解$Obj$的最小值来求解树结构$T$，然后以此为基础求解出$w_j$的过程已经没有什么问题了。回忆一下，为了找出最小的$Obj$，需要枚举所有可能的树结构，这似乎又回到了最初的困惑——之所以要使用迭代和最优化的方式，就是因为不希望进行枚举，这样既浪费资源又浪费时间。因此，可以采用贪婪算法寻找最佳结构。

3.6 寻找最佳分枝：结构分数之差

贪婪算法指的是控制局部最优来达到全局最优的算法，决策树算法本身就是一种使用贪婪算法的方法。XGB作为树的集成模型，自然也想到采用这样的方法来进行计算，所以认为，如果每片叶子都是最优的，则整体生成的树结构就是最优的，如此就可以避免去枚举所有可能的树结构。

回忆一下决策树中是如何进行计算的：使用基尼系数或信息熵来衡量分枝之后叶子节点的不纯度，分枝前的信息熵与分枝后的信息熵之差叫做信息增益，信息增益最大的特征上的分枝就被选中，当信息增益低于某个阈值时，就让树停止生长。在XGB中，使用的方式是类似的：首先使用目标函数来衡量树的结构的优劣，然后让树从深度0开始生长，每进行一次分枝，就计算目标函数减少了多少，当目标函数的降低低于设定的某个阈值时，就让树停止生长。
来看具体的例子，下图中，有中间节点“是否是男性”，这个中间节点下面有两个叶子节点，分别是样本弟弟和妹妹，看这个分枝点上，树的结构分数之差如何表示。

对于中间节点这个叶子节点而言，$T=1$，则这个节点上的结构分数为：
$I=1,4$
$G=g_1+g_4$
$H=h_1+h_4$
$Scor{e}_{middle}=-\frac{1}{2}\frac{G^2}{H+\lambda }+\gamma$
对于弟弟和妹妹节点而言，则有：
$Scor{e}_{sis}=-\frac{1}{2}\frac{g_4^2}{h_4+\lambda }+\gamma$
$Scor{e}_{bro}=-\frac{1}{2}\frac{g_1^2}{h_1+\lambda }+\gamma$
则分枝后的结构分数之差为：
$Gain=Scor{e}_{sis}+Scor{e}_{bro}-Scor{e}_{middle}$
$=-\frac{1}{2}\frac{g_4^2}{h_4+\lambda }+\gamma -\frac{1}{2}\frac{g_1^2}{h_1+\lambda }+\gamma -(-\frac{1}{2}\frac{G^2}{H+\lambda }+\gamma )$
$=-\frac{1}{2}\frac{g_4^2}{h_4+\lambda }+\gamma -\frac{1}{2}\frac{g_1^2}{h_1+\lambda }+\gamma +\frac{1}{2}\frac{G^2}{H+\lambda }-\gamma$
$=-\frac{1}{2}[\frac{g_4^2}{h_4+\lambda }+\frac{g_1^2}{h_1+\lambda }-\frac{G^2}{H+\lambda }]+\gamma$
将负号去除：
$=\frac{1}{2}[\frac{g_4^2}{h_4+\lambda }+\frac{g_1^2}{h_1+\lambda }-\frac{(g_1+g_4{)}^2}{(h_1+h_4)+\lambda }]-\gamma$
CART树全部是二叉树，因此这个式子可以推广。从这个式子可以总结出，其实分枝后的结构分数之差为：
$Gain=\frac{1}{2}[\frac{G_L^2}{H_L+\lambda }+\frac{G_R^2}{H_R+\lambda }-\frac{(G_L+G_R{)}^2}{(H_L+H_R)+\lambda }]-\gamma$
其中$G_L$和$H_L$是从左节点（弟弟节点）上计算得出，$G_R$和$H_R$是从右节点（妹妹节点）上计算得出，而$(G_L+G_R)$和$(H_L+H_R)$是从中间节点上计算得出。对于任意分枝，都可以这样来进行计算。在现实中，会对所有特征的所有分枝点进行如上计算，然后选出让目标函数下降最快的节点来进行分枝。对每一棵树的每一层，都进行这样的计算，比起原始的梯度下降，实践证明这样的求解最佳树结构的方法运算最快，并且在大型数据下也能够表现不错。至此，作为XGBoost的使用者，已经将需要理解的XGB的原理理解完毕了。

3.7 让树停止生长：重要参数gamma

在之前所有的推导过程中，都没有提到$\gamma$这个变量。从目标函数和结构分数之差$Gain$的式子中来看，$\gamma$是每增加一片叶子就会被剪去的惩罚项。增加的叶子越多，结构分数之差$Gain$会被惩罚越重，所以$\gamma$又被称之为是“复杂性控制（complexity control）”，所以$\gamma$是用来防止过拟合的重要参数。实践证明，$\gamma$是对梯度提升树影响最大的参数之一，其效果丝毫不逊色于n_estimators和防止过拟合的神器max_depth。同时，$\gamma$还是让树停止生长的重要参数。
在逻辑回归中，使用参数$tol$来设定阈值，并规定如果梯度下降时损失函数减小量小于$tol$下降就会停止。**在XGB中，规定只要结构分数之差$Gain$是大于0的，即只要目标函数还能够继续减小，就允许树继续进行分枝。**也就是说，对于目标函数减小量的要求是：
$\frac{1}{2}[\frac{G_L^2}{H_L+\lambda }+\frac{G_R^2}{H_R+\lambda }-\frac{(G_L+G_R{)}^2}{(H_L+H_R)+\lambda }]-\gamma >0$
$\frac{1}{2}[\frac{G_L^2}{H_L+\lambda }+\frac{G_R^2}{H_R+\lambda }-\frac{(G_L+G_R{)}^2}{(H_L+H_R)+\lambda }]>\gamma$
如此，可以直接通过设定$\gamma$的大小来让XGB中的树停止生长。$\gamma$因此被定义为，在树的叶节点上进行进一步分枝所需的最小目标函数减少量，在决策树和随机森林中也有类似的参数（min_split_loss，min_samples_split）。$\gamma$设定越大，算法就越保守，树的叶子数量就越少，模型的复杂度就越低。

参数含义	xgb.train()	xgb.XGBRegressor()
复杂度的惩罚项$\gamma$	gamma，默认为0，取值范围$[0,+\infty ]$	gamma，默认为0，取值范围$[0,+\infty ]$

如果希望从代码中来观察$\gamma$的作用，使用sklearn中传统的学习曲线等工具就比较苦难。下面这段代码，是一段让参数$\gamma$在0-5之间均匀取值的学习曲线。其运行速度较缓慢并且曲线的效果匪夷所思。

axisx = np.arange(0,5,0.05)
rs = []
var = []
ge = []
for i in axisx:
    reg = XGBR(n_estimators = 100, random_state = 420, gamma = i)
    result = CVS(reg, xtrain, ytrain, cv = cv)
    rs.append(result.mean())
    var.append(result.var())
    ge.append((1-result.mean())**2+result.var())
print(axisx[rs.index(max(rs))], max(rs), var[rs.index(max(rs))])
print(axisx[var.index(min(var))],rs[var.index(min(var))], min(var))
print(axisx[ge.index(min(ge))],rs[ge.index(min(ge))], var[ge.index(min(ge))],min(ge))
rs = np.array(rs)
var = np.array(var)*0.1
plt.figure(figsize = (20,5))
plt.plot(axisx,rs,c = "black",label = "XGB")
plt.plot(axisx,rs+var,c = "red",linestyle = "-.")
plt.plot(axisx,rs-var,c = "red",linestyle = "-.")
plt.legend()
plt.show()
'''
结果：
0.8 0.8408699341307502 0.005891402360996227
4.55 0.8270850743887515 0.004605425470082616
0.6000000000000001 0.8398453472024396 0.004660990778716479 0.030310503591423586
'''

可以看到，完全无法从中看出什么趋势，偏差时高时低，方差时大时小，参数$\gamma$引起的波动远远超过其他参数（其他参数至少还有一个先上升再平稳的过程，而$\gamma$则是仿佛完全无规律）。在sklearn下XGBoost太不稳定，如果这样来调整参数的话，效果就很难保证。因此，为了调整$\gamma$，需要来引入新的工具，xgboost库中的类xgboost.cv。
xgboost.cv(params, dtrain, num_boost_round = 10, nfold = 3, stratified = False, folds = None, metrics = (), obj = None, feval = None, maximize = False, early_stopping_rounds = None, fpreproc = None, as_pandas = True, verbose_eval = None, show_stdv = True, seed = 0, callbacks = None, shuffle = True)

import xgboost as xgb
#为了便捷，使用全数据
dfull = xgb.DMatrix(x,y)
#设定参数
param1 = {"verbosity":0,"obj":"reg:linear","gamma":0}
num_round = 100
n_fold = 5
#使用类xgb.cv
time0 = time()
cvresult1 = xgb.cv(param1,dfull,num_round,n_fold)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))
#结果：00:01:174784
#看看类xgb.cv生成了什么结果
cvresult1

plt.figure(figsize = (20,5))
plt.grid()
plt.plot(range(1,101),cvresult1.iloc[:,0],c = "red", label = "train,gamma=0")
plt.plot(range(1,101),cvresult1.iloc[:,2],c = "orange", label = "test,gamma=0")
plt.legend()
plt.show()
#xgboost中回归模型的默认模型评估指标是什么？——RMSE

为了使用xgboost.cv，必须要熟悉xgboost自带的模型评估指标。xgboost在建库的时候奔着大而全的目标，和sklearn类似，包括了大约20个模型评估指标，然而用于回归和分类的其实只有几个，大部分是用于一些更加高级的功能，如ranking。看看用于回归和分类的指标都有哪些：

指标	含义
rmse	回归用，调整后的均方误差
mae	回归用，绝对平方误差
logloss	二分类用，对数损失
mlogloss	多分类用，对数损失
error	分类用，分类误差，等于1-准确率
auc	分类用，AUC面积

param1 = {"verbosity":0,"obj":"reg:linear","gamma":0, "eval_metrics":"mae"}
cvresult1 = xgb.cv(param1,dfull,num_round,n_fold)
plt.figure(figsize = (20,5))
plt.grid()
plt.plot(range(1,101),cvresult1.iloc[:,0],c = "red", label = "train,gamma=0")
plt.plot(range(1,101),cvresult1.iloc[:,2],c = "orange", label = "test,gamma=0")
plt.legend()
plt.show()
#从这个图中，可以看出什么？
#怎样从图中观察模型的泛化能力？
#从这个图的角度来说，模型的调参目标是什么？

来看看如果调整$\gamma$，会发生怎样的变化：

param1 = {"verbosity":0,"obj":"reg:linear","gamma":0}
param2 = {"verbosity":0,"obj":"reg:linear","gamma":20}
num_round = 100
n_fold = 5

time0 = time()
cvresult1 = xgb.cv(param1,dfull,num_round,n_fold)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

time0 = time()
cvresult2 = xgb.cv(param2,dfull,num_round,n_fold)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

plt.figure(figsize = (20,5))
plt.grid()
plt.plot(range(1,101),cvresult1.iloc[:,0],c = "red", label = "train,gamma=0")
plt.plot(range(1,101),cvresult1.iloc[:,2],c = "orange", label = "test,gamma=0")
plt.plot(range(1,101),cvresult2.iloc[:,0],c = "green", label = "train,gamma=20")
plt.plot(range(1,101),cvresult2.iloc[:,2],c = "blue", label = "test,gamma=20")
plt.legend()
plt.show()
#从这里，看出gamma是如何控制过拟合了吗？——控制训练集上的训练，降低训练集上的表现
'''
结果：
00:01:139495
00:01:287364
'''
#gamma变大，运行时间变长，说明是后剪枝

试一个分类的例子：

from sklearn.datasets import load_breast_cancer
data2 = load_breast_cancer()

x2 = data2.data
y2 = data2.target
dfull2 = xgb.DMatrix(x2,y2)
param1 = {"verbosity":0, "obj":"binary:logistic", "gamma":0, "nfold":5}
param2 = {"verbosity":0, "obj":"binary:logistic", "gamma":2, "nfold":5}
num_round = 100

time0 = time()
cvresult1 = xgb.cv(param1,dfull2,num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

time0 = time()
cvresult2 = xgb.cv(param2,dfull2,num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

plt.figure(figsize = (20,5))
plt.grid()
plt.plot(range(1,101),cvresult1.iloc[:,0],c = "red", label = "train,gamma=0")
plt.plot(range(1,101),cvresult1.iloc[:,2],c = "orange", label = "test,gamma=0")
plt.plot(range(1,101),cvresult2.iloc[:,0],c = "green", label = "train,gamma=2")
plt.plot(range(1,101),cvresult2.iloc[:,2],c = "blue", label = "test,gamma=2")
plt.legend()
plt.show()
'''
结果：
00:00:702175
00:00:835806
'''

有了xgboost.cv这个工具，参数调整就容易多了。这个工具可以让我们直接看到参数如何影响了模型的泛化能力。

4 XGBoost应用中的其他问题

4.1 过拟合：剪枝参数与回归模型调参

xgboost.XGBRegressor(max_depth = 3, learning_rate = 0.1, n_estimators = 100, silent = True, objective = ‘reg:linear’, booster = ‘gbtree’, n_jobs = 1, nthread = None, gamma = 0, min_child_weight = 1, max_delta_step = 0, subsample = 1, colsample_bytree = 1, colsample_bylevel = 1, reg_alpha = 0, reg_lambda = 1, scale_pos_weight = 1, base_score = 0.5, random_state = 0, seed = None, missing = None, importance_type = ‘gain’, **kwargs)
作为天生过拟合的模型，XGBoost应用的核心之一就是减轻过拟合带来的影响。作为树模型，减轻过拟合的方式主要是靠对决策树剪枝来降低模型的复杂度，以求降低方差。已经学习了好几个可以用来防止过拟合的参数，包括复杂度控制$\gamma$，正则化的两个参数$\lambda$和$\alpha$，控制迭代速度的参数$\eta$以及管理每次迭代前进行的随机有放回抽样的参数subsample。所有的这些参数都可以用来减轻过拟合。但除此之外，还有几个影响重大的，专用于剪枝的参数：

参数含义	xgb.train()	xgb.XGBRegressor()
树的最大深度	max_depth，默认为6	max_depth，默认为6
每次生成树时随机抽样特征的比例	colsample_bytree，默认为1	colsample_bytree，默认为1
每次生成树的一层时，随机抽样特征的比例	colsample_bylevel，默认为1	colsample_bylevel，默认为1
每次生成一个叶子节点时，随机抽样特征的比例	colsample_bynode，默认为1	N.A.
一个叶子节点上所需要的最小$h_i$，即叶子节点上的二阶导数之和，类似于样本权重	min_child_weight，默认为1	min_child_weight，默认为1

这些参数中，树的最大深度是决策树中的剪枝法宝，算是最常用的剪枝参数，不过在XGBoost中，最大深度的功能与参数$\gamma$相似，因此如果先调节了$\gamma$，则最大深度可能无法展示出巨大的效果。当然，如果先调整了最大深度，则$\gamma$也有可能无法显示明显的效果。通常来说，这两个参数中只使用一个，不过两个都试试也没有坏处。
三个随机抽样特征的参数中，前两个比较常用。在建立树时对特征进行抽样其实是决策树和随机森林中比较常见的一种方法，但是在XGBoost之前，这种方法并没有被使用到boosting算法当中过，Boosting算法一直以抽取样本（横向抽样）来调整模型过拟合的程度， 而实践证明其实纵向抽样（抽取特征）更能够防止过拟合。
参数min_child_weight不太常用，它是一片叶子上的二阶导数$h_i$之和，当样本所对应的二阶导数很小时，比如说为0.01，min_child_weight若设定为1，则说明一片叶子上至少需要100个样本。本质上来说，这个参数其实是在控制叶子上所需的最小样本量，因此对于样本量很大的数据会比较有效。如果样本量很小（如现在使用的波士顿房价数据集，则这个参数效用不大）。就剪枝的效果来说，这个参数的功能也被$\gamma$替代了一部分，通常来说会试试看这个参数，但这个参数不是优先选择。
通常当获得了一个数据集后，先使用网络搜索找出比较合适的n_estimators和eat组合，然后使用gamma或者max_depth观察模型处于什么样的状态（过拟合还是欠拟合，处于方差-偏差图像的左边还是右边），最后再决定是否要进行剪枝。通常来说，对于XGB模型，大多数时候都是需要剪枝的。接下来就来看看使用xgb.cv这个类来进行剪枝调参，以调整出一组泛化能力很强的参数。
先从最原始的，设定默认参数开始，先观察一下默认参数下，交叉验证曲线的趋势：

dfull = xgb.DMatrix(x,y)
param1 = {"verbosity":0  #并非默认
         ,"obj":"reg:linear"   #并非默认
         ,"subsample":1
         ,"max_depth":6
         ,"eta" :0.3
         ,"gamma":0
         ,"lambda":1
         ,"alpha":0
         ,"colsample_bytree":1
         ,"colsample_bylevel":1
         ,"colsample_bynode":1
         ,"nfold":5}
num_round = 200

time0 = time()
cvresult1 = xgb.cv(param1, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

fig,ax = plt.subplots(1, figsize = (15,10))
#ax.set_ylim(top = 5)
ax.grid()
ax.plot(range(1,201),cvresult1.iloc[:,0],c = "red", label = "train,original")
ax.plot(range(1,201),cvresult1.iloc[:,2],c = "orange", label = "train,original")
ax.legend(fontsize = "xx-large")
plt.show()
#结果：00:01:387368

从曲线上可以看出，模型现在处于过拟合的状态，决定要进行剪枝。目标是：训练集和测试集的曲线尽量接近，如果测试集上的结果不能上升，那训练集上的结果降下来也是不错的选择（让模型不那么具体到训练数据，增加泛化能力）。在这里，要使用三组曲线。一组用于展示原始数据上的结果，一组用于展示上一个参数调节完毕后的结果，最后一组用于展示现在调节的参数的结果。

dfull = xgb.DMatrix(x,y)
param1 = {"verbosity":0  #并非默认
         ,"obj":"reg:linear"   #并非默认
         ,"subsample":1
         ,"max_depth":6
         ,"eta" :0.3
         ,"gamma":0
         ,"lambda":1
         ,"alpha":0
         ,"colsample_bytree":1
         ,"colsample_bylevel":1
         ,"colsample_bynode":1
         ,"nfold":5}
num_round = 200

time0 = time()
cvresult1 = xgb.cv(param1, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

fig,ax = plt.subplots(1, figsize = (15,10))
ax.set_ylim(top = 5)
ax.grid()
ax.plot(range(1,201),cvresult1.iloc[:,0],c = "red", label = "train,original")
ax.plot(range(1,201),cvresult1.iloc[:,2],c = "orange", label = "test,original")

param2 = {"verbosity":0
         ,"max_depth":2
         ,"obj":"reg:linear"
         ,"nfold":5}
param3 = {"verbosity":0
         ,"max_depth":1
         ,"obj":"reg:linear"
         ,"nfold":5}

time0 = time()
cvresult2 = xgb.cv(param2, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

time0 = time()
cvresult3 = xgb.cv(param3, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

ax.plot(range(1,201),cvresult2.iloc[:,0],c = "green", label = "train,last")
ax.plot(range(1,201),cvresult2.iloc[:,2],c = "blue", label = "test,last")
ax.plot(range(1,201),cvresult3.iloc[:,0],c = "gray", label = "train,this")
ax.plot(range(1,201),cvresult3.iloc[:,2],c = "pink", label = "test,this")
ax.legend(fontsize = "xx-large")
plt.show()
'''
结果：
00:01:353010
00:01:305881
00:01:250819
'''

利用xgb.cv进行调参，逐个改变，每个参数的数值，来观察训练集和测试集的曲线变化，以找到训练集和测试集曲线之间距离最小，测试集效果不降低且运行时间更少的模型。同时调参的顺序不同，最终的结果会有所不同。

dfull = xgb.DMatrix(x,y)
param1 = {"verbosity":0  #并非默认
         ,"obj":"reg:linear"   #并非默认
         ,"subsample":1
         ,"max_depth":6
         ,"eta" :0.3
         ,"gamma":0
         ,"lambda":1
         ,"alpha":0
         ,"colsample_bytree":1
         ,"colsample_bylevel":1
         ,"colsample_bynode":1
         ,"nfold":5}
num_round = 200

time0 = time()
cvresult1 = xgb.cv(param1, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

fig,ax = plt.subplots(1, figsize = (15,10))
ax.set_ylim(top = 5)
ax.grid()
ax.plot(range(1,201),cvresult1.iloc[:,0],c = "red", label = "train,original")
ax.plot(range(1,201),cvresult1.iloc[:,2],c = "orange", label = "test,original")

param2 = {"verbosity":0
         ,"obj":"reg:linear"
         ,"max_depth":2
         ,"eta":0.05
         ,"gamma":0
         ,"lambda":1
         ,"alpha":0
         ,"colsample_bytree":1
         ,"colsample_bylevel":0.4
         ,"colsample_bynode":1
         ,"nfold":5}
param3 = {"verbosity":0
         ,"obj":"reg:linear"
         ,"subsample":1
         ,"eta":0.05
         ,"gamma":20
         ,"lambda":3.5
         ,"alpha":0.2
         ,"max_depth":4
         ,"colsample_bytree":0.4
         ,"colsample_bylevel":0.6
         ,"colsample_bynode":1
         ,"nfold":5}

time0 = time()
cvresult2 = xgb.cv(param2, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

time0 = time()
cvresult3 = xgb.cv(param3, dfull, num_round)
print(datetime.datetime.fromtimestamp(time()-time0).strftime("%M:%S:%f"))

ax.plot(range(1,201),cvresult2.iloc[:,0],c = "green", label = "train,last")
ax.plot(range(1,201),cvresult2.iloc[:,2],c = "blue", label = "test,last")
ax.plot(range(1,201),cvresult3.iloc[:,0],c = "gray", label = "train,this")
ax.plot(range(1,201),cvresult3.iloc[:,2],c = "pink", label = "test,this")
ax.legend(fontsize = "xx-large")
plt.show()
'''
结果：
00:01:405866
00:01:491525
00:01:410872
'''

在这个调整过程中，可能会有几个问题：

1. 一个个参数调整太麻烦，可不可以使用网格搜索？
   当然可以！主要电脑有足够的计算资源，并且信任网格搜索，那任何时候都可以使用网格搜索。只是使用的时候要注意，首先XGB的参数非常多，参数可取的范围也很广，究竟是使用np.linspace或者np.arange作为参数的备选值也会影响结果，而且网格搜索的运行速度往往不容乐观，因此建议至少先使用xgboost.cv来确认参数的范围，否则很可能花很长时间做了无用功。
   并且，在使用网格搜索的时候，最好不要一次性将所有的参数都放入进行搜索，最多一次两三个。有一些互相影响的参数需要放在一起使用，比如学习率eta和树的数量n_estimators。
   另外，如果网格搜索的结果与理解相违背，与手动调参的结果相违背，选择模型效果较好的一个。如果两者效果差不多，那选择相信手动调参的结果。网格毕竟是枚举出结果，很多时候得出的结果可能会是具体到数据的巧合，无法去一一解释网格搜索得出的结论为何是这样。如果感觉都无法解释，那就不要在意，直接选择结果较好的一个。
2. 调参的时候参数的顺序会影响调参结果吗？
   会影响，因此在现实中，会优先调整那些对模型影响巨大的参数。在这里，建议的剪枝上的调参顺序是：n_estimators与eta共同调节，gamma或者max_depth，采样和抽样参数（纵向抽样影响更大），最后才是正则化的两个参数。当然，可以根据自己的需求来进行调整。
3. 调参之后测试集上的效果害没有原始设定上的效果好怎么办？
   如果调参之后，交叉验证曲线确实显示测试集和训练集上的模型评估效果是更加接近的，推荐使用调参之后的效果。希望增强模型的泛化能力，然而泛化能力的增强并不代表着在新数据集上模型的结果一定优秀，因为未知数据集并未一定符合全数据的分布，在一组未知数据上表现十分优秀，也不一定就能够在其他的未知数据集上表现优秀。因此不必过于纠结在现有的测试集上是否表现优秀。当然了，在现有数据上如果能够实现训练集和测试集都非常优秀，那模型的泛化能力自然也会是很强的。

4.2 XGBoost模型的保存和调用

在使用Python进行编程时，可能会需要编写较为复杂的程序或者建立复杂的模型。比如XGBoost模型，这个模型的参数复杂繁多，并且调参过程不是太容易，一旦训练完毕，往往希望将训练完毕后的模型保存下来，以便日后用于新的数据集。在Python中，保存模型的方法有许多种。以XGBoost为例，来讲解两种主要的模型保存和调用方法。

4.2.1 使用Pickle保存和调用模型

pickle是python编程中比较标准的一个保存和调用模型的库，可以使用pickle和open函数的连用，来将模型保存到本地。以刚才已经调整好的参数和训练好的模型为例，可以这样来使用pickle：

import pickle
dtrain = xgb.DMatrix(xtrain,ytrain)
#设定参数，对模型进行训练
param = {"verbosity":0
         ,"obj":"reg:linear"
         ,"subsample":1
         ,"eta":0.05
         ,"gamma":20
         ,"lambda":3.5
         ,"alpha":0.2
         ,"max_depth":4
         ,"colsample_bytree":0.4
         ,"colsample_bylevel":0.6
         ,"colsample_bynode":1}
num_round = 200
bst = xgb.train(param, dtrain, num_round)
#保存模型
pickle.dump(bst, open("xgboostonboston.dat","wb"))
#注意，open中往往使用w或者r作为读取的模式，但其实w与r只能用于文本文件(.txt)，当希望导入的不是文本文件，而是模型本身的时候，使用"wb"和"rb"作为读取的模式。其中wb表示以二进制写入，rb表示以二进制读入

#看看模型被保存到了哪里？
import sys
sys.path
'''
结果：
['C:\\Users\\86188\\Desktop\\coding',
 'D:\\Anaconda\\python38.zip',
 'D:\\Anaconda\\DLLs',
 'D:\\Anaconda\\lib',
 'D:\\Anaconda',
 '',
 'D:\\Anaconda\\lib\\site-packages',
 'D:\\Anaconda\\lib\\site-packages\\locket-0.2.1-py3.8.egg',
 'D:\\Anaconda\\lib\\site-packages\\win32',
 'D:\\Anaconda\\lib\\site-packages\\win32\\lib',
 'D:\\Anaconda\\lib\\site-packages\\Pythonwin',
 'D:\\Anaconda\\lib\\site-packages\\IPython\\extensions',
 'C:\\Users\\86188\\.ipython']
'''

#重新打开jupyter notebook
from sklearn.datasets import load_boston
from sklearn.model_selection import train_test_split as TTS
from sklearn.metrics import mean_squared_error as MSE
import pickle
import xgboost as xgb

data = load_boston()
x = data.data
y = data.target

xtrain,xtest,ytrain,ytest = TTS(x,y,test_size = 0.3, random_state = 420)
#注意，如果保存的模型是xgboost库中建立的模型，则导入的数据类型也必须是xgboost库中的数据类型
dtest = xgb.DMatrix(xtest,ytest)

#导入模型
loaded_model = pickle.load(open("xgboostonboston.dat","rb"))
print("Loaded model from: xgboostonboston.dat")
#结果：Loaded model from: xgboostonboston.dat

#做预测，直接调用接口predict
ypreds = loaded_model.predict(dtest)
from sklearn.metrics import mean_squared_error as MSE, r2_score
MSE(ytest,ypreds)
#结果：9.18818551539731
r2_score(ytest,ypreds)
#结果：0.9012595174151566

4.2.2 使用Joblib保存和调用模型

Joblib是Scipy生态系统中的一部分，它为python提供保存和调用管道和对象的功能，处理NumPy结构的数据尤其高效，对于很大的数据集和巨大的模型非常有用。Joblib与pickle API非常相似，来看看代码：

import joblib
#同样可以看看模型被保存到哪里了
joblib.dump(bst,"xgboost-boston.dat")
loaded_model = joblib.load("xgboost-boston.dat")
ypreds = loaded_model.predict(dtest)
MSE(ytest,ypreds)
#结果：9.18818551539731
r2_score(ytest,ypreds)
#结果：0.9012595174151566

#使用sklearn中的模型
from xgboost import XGBRegressor as XGBR
bst = XGBR(n_estimators = 200
          ,eta = 0.05
          ,gamma = 20
          ,reg_lambda = 3.5
          ,reg_alpha = 0.2
          ,max_depth = 4
          ,colsample_bytree = 0.4
          ,colsample_bylevel = 0.6).fit(xtrain,ytrain)#训练完毕
joblib.dump(bst,"xgboost-boston-sklearn.dat")
loaded_model = joblib.load("xgboost-boston-sklearn.dat")
#这里可以直接导入xtest
ypreds = loaded_model.predict(xtest)
MSE(ytest,ypreds)
#结果：9.18818551539731
r2_score(ytest,ypreds)
#结果：0.9012595174151566

在这两种保存方法下，都可以找到保存下来的dat文件，将这些文件移动到任意计算机上的python下的环境变量路径中（使用sys.path查看），则可以使用import来对模型进行调用。注意，模型的保存调用与自写函数的保存调用是两回事，要注意区分。

4.3 分类案例：XGB中的样本不均衡问题

虽然回归是XGB的常用领域，但是作为机器学习中的大头，分类算法也是不可忽视的，XGB作为分类的例子也非常多。存在分类，就会存在样本不平衡问题带来的影响，XGB中存在着调节样本不平衡的参数scale_pos_weight，这个参数非常类似于之前随机森林、支持向量机中都使用过的class_weight参数，通常在参数中输入的是负样本量与正样本量之比$\frac{sum(negativeinstances)}{sum(positiveinstances)}$。

参数含义	xgb.train()	xgb.XGBClassifier()
控制正负样本比例，表示为负/正样本比例，在样本不平衡问题中使用	scale_pos_weight，默认为1	scale_pos_weight，默认为1

看看如何使用这个参数：

#导库，创建样本不平衡的数据集
import numpy as np
import xgboost as xgb
import matplotlib.pyplot as plt
from xgboost import XGBClassifier as XGBC
from sklearn.datasets import make_blobs
from sklearn.model_selection import train_test_split as TTS
from sklearn.metrics import confusion_matrix as cm, recall_score as recall, roc_auc_score as auc

class_1 = 500 #类别1有500个样本
class_2 = 50  #类别2只有50个
centers = [[0.0,0.0],[2.0,2.0]]  #设定两个类别的中心
clusters_std = [1.5,0.5]  #设定两个类别的方差，通常来说，样本量比较大的类别会更加松散
x,y = make_blobs(n_samples = [class_1,class_2]
                ,centers = centers
                ,cluster_std = clusters_std
                ,random_state = 0
                ,shuffle = False)
xtrain,xtest,ytrain,ytest = TTS(x,y,test_size = 0.3, random_state = 420)

(y==1).sum()/y.shape[0]
#结果：0.09090909090909091

#在数据集上建模：sklearn模式
#在sklearn下建模
clf = XGBC().fit(xtrain,ytrain)
ypred = clf.predict(xtest)
clf.score(xtest,ytest)#默认模型评估指标——准确率
#结果：0.9272727272727272
cm(ytest,ypred,labels = [1,0])#少数类写在前面
'''
结果：
array([[  9,   4],
       [  8, 144]], dtype=int64)
'''
recall(ytest,ypred)
#结果：0.6923076923076923
auc(ytest,clf.predict_proba(xtest)[:,1])
#结果：0.9701417004048585
#负/正样本比例
clf_ = XGBC(scale_pos_weight = 10).fit(xtrain,ytrain)
ypred_ = clf_.predict(xtest)
clf_.score(xtest,ytest)
#结果：0.9333333333333333
cm(ytest,ypred_,labels = [1,0])
'''
结果：
array([[ 10,   3],
       [  8, 144]], dtype=int64)
'''
recall(ytest,ypred_)
#结果：0.7692307692307693
auc(ytest,clf_.predict_proba(xtest)[:,1])
#结果：0.9696356275303644

#随着样本权重逐渐增加，模型的recall,auc和准确率如何变化？
for i in [1,5,10,20,30]:
    clf_ = XGBC(scale_pos_weight = i).fit(xtrain,ytrain)
    ypred_ = clf_.predict(xtest)
    print(i)
    print("\tAccuracy:{}".format(clf_.score(xtest,ytest)))
    print("\tRecall:{}".format(recall(ytest,ypred_)))
    print("\tAUC:{}".format(auc(ytest,clf_.predict_proba(xtest)[:,1])))
'''
结果：
1
	Accuracy:0.9272727272727272
	Recall:0.6923076923076923
	AUC:0.9701417004048585
5
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9660931174089069
10
	Accuracy:0.9333333333333333
	Recall:0.7692307692307693
	AUC:0.9696356275303644
20
	Accuracy:0.9333333333333333
	Recall:0.7692307692307693
	AUC:0.9686234817813765
30
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9701417004048583	

'''

#在数据集上建模：xgboost模式
dtrain = xgb.DMatrix(xtrain,ytrain)
dtest = xgb.DMatrix(xtest,ytest)
#看看xgboost库自带的predict接口
param = {"verbosity":0
        ,"objective":"binary:logistic"
        ,"eta":0.1
        ,"scale_pos_weight":1}
num_round = 100
bst = xgb.train(param,dtrain,num_round)
preds = bst.predict(dtest)
#看看preds返回了什么——概率值
preds
'''
结果：
array([0.00110357, 0.00761518, 0.00110357, 0.00110357, 0.93531454,
       0.00466839, 0.00110357, 0.00110357, 0.00110357, 0.00110357,
       0.00110357, 0.00410493, 0.00454478, 0.00571528, 0.00751026,
       0.00110357, 0.00110357, 0.00110357, 0.00110357, 0.00110357,
       0.00110357, 0.00110357, 0.00110357, 0.00110357, 0.00110357,
       0.00712637, 0.00110357, 0.00110357, 0.00110357, 0.00110357,
       0.00110357, 0.00110357, 0.00110357, 0.00793251, 0.00466839,
       0.00110357, 0.00339395, 0.00657186, 0.00110357, 0.00457053,
       0.00571528, 0.0026763 , 0.00110357, 0.00110357, 0.00110357,
       0.00884932, 0.00712637, 0.00110357, 0.00712637, 0.00466839,
       0.00110357, 0.00110357, 0.00712637, 0.00110357, 0.00110357,
       0.00110357, 0.00110357, 0.63748044, 0.00110357, 0.00793251,
       0.00110357, 0.00451971, 0.00644181, 0.00110357, 0.00110357,
       0.00110357, 0.00110357, 0.00751026, 0.00712637, 0.00110357,
       0.00866458, 0.00110357, 0.00110357, 0.00110357, 0.91610426,
       0.00110357, 0.00110357, 0.89246494, 0.0026763 , 0.00501714,
       0.00761518, 0.00884932, 0.00339395, 0.00110357, 0.93531454,
       0.00110357, 0.00110357, 0.00110357, 0.82530665, 0.00751026,
       0.00110357, 0.35174078, 0.00110357, 0.00110357, 0.70393246,
       0.00110357, 0.76804197, 0.00110357, 0.00110357, 0.00110357,
       0.00110357, 0.96656513, 0.00110357, 0.00571528, 0.25400913,
       0.00110357, 0.00110357, 0.00110357, 0.00110357, 0.00457053,
       0.00110357, 0.00110357, 0.00110357, 0.89246494, 0.00110357,
       0.9518535 , 0.0026763 , 0.00712637, 0.00110357, 0.00501714,
       0.00110357, 0.00110357, 0.00571528, 0.00110357, 0.00110357,
       0.00712637, 0.00110357, 0.00110357, 0.00712637, 0.00110357,
       0.25136763, 0.00110357, 0.00110357, 0.00110357, 0.00110357,
       0.00110357, 0.8904051 , 0.3876418 , 0.00110357, 0.00457053,
       0.00657186, 0.9366597 , 0.00866458, 0.00110357, 0.00501714,
       0.00501714, 0.00110357, 0.00110357, 0.00368543, 0.00501714,
       0.9830577 , 0.00110357, 0.00644181, 0.00110357, 0.00571528,
       0.00110357, 0.00110357, 0.00110357, 0.00110357, 0.00466839,
       0.00110357, 0.00110357, 0.92388713, 0.90231985, 0.80084217],
      dtype=float32)

'''
#自己的设定阈值
ypred = preds.copy()
ypred[preds>0.5] = 1
ypred[preds != 1] = 0
#写明参数
scale_pos_weight = [1,5,10]
names = ["negative vs positive: 1"
        ,"negative vs positive: 5
        ,"negative vs positive: 10"]
#导入模型评估指标
from sklearn.metrics import accuracy_score as accuracy, recall_score as recall, roc_auc_score as auc
for name,i in zip(names,scale_pos_weight):
    param = {"verbosity":0
            ,"objective":"binary:logistic"
            ,"eta":0.1
            ,"scale_pos_weight":i}
    clf = xgb.train(param,dtrain,num_round)
    preds = clf.predict(dtest)
    ypred = preds.copy()
    ypred[preds > 0.5] = 1
    ypred[ypred != 1] = 0
    print(name)
    print("\tAccuracy:{}".format(accuracy(ytest,ypred)))
    print("\tRecall:{}".format(recall(ytest,ypred)))
    print("\tAUC:{}".format(auc(ytest,preds)))
'''
结果：
negative vs positive: 1
	Accuracy:0.9272727272727272
	Recall:0.6923076923076923
	AUC:0.9741902834008097
negative vs positive: 5
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9635627530364372
negative vs positive: 10
	Accuracy:0.9515151515151515
	Recall:1.0
	AUC:0.9665991902834008
'''
#当然也可以尝试不同的阈值
for name,i in zip(names,scale_pos_weight):
    for thres in [0.3,0.5,0.7,0.9]:
            param = {"verbosity":0
                    ,"objective":"binary:logistic"
                    ,"eta":0.1
                    ,"scale_pos_weight":i}
            clf = xgb.train(param,dtrain,num_round)
            preds = clf.predict(dtest)
            ypred = preds.copy()
            ypred[preds > thres] = 1
            ypred[ypred != 1] = 0
            print("{},thresholds:{}".format(name,thres))
            print("\tAccuracy:{}".format(accuracy(ytest,ypred)))
            print("\tRecall:{}".format(recall(ytest,ypred)))
            print("\tAUC:{}".format(auc(ytest,preds)))
'''
结果：
negative vs positive: 1,thresholds:0.3
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9741902834008097
negative vs positive: 1,thresholds:0.5
	Accuracy:0.9272727272727272
	Recall:0.6923076923076923
	AUC:0.9741902834008097
negative vs positive: 1,thresholds:0.7
	Accuracy:0.9212121212121213
	Recall:0.6153846153846154
	AUC:0.9741902834008097
negative vs positive: 1,thresholds:0.9
	Accuracy:0.9515151515151515
	Recall:0.5384615384615384
	AUC:0.9741902834008097
negative vs positive: 5,thresholds:0.3
	Accuracy:0.9515151515151515
	Recall:1.0
	AUC:0.9635627530364372
negative vs positive: 5,thresholds:0.5
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9635627530364372
negative vs positive: 5,thresholds:0.7
	Accuracy:0.9272727272727272
	Recall:0.6923076923076923
	AUC:0.9635627530364372
negative vs positive: 5,thresholds:0.9
	Accuracy:0.9212121212121213
	Recall:0.6153846153846154
	AUC:0.9635627530364372
negative vs positive: 10,thresholds:0.3
	Accuracy:0.9515151515151515
	Recall:1.0
	AUC:0.9665991902834008
negative vs positive: 10,thresholds:0.5
	Accuracy:0.9515151515151515
	Recall:1.0
	AUC:0.9665991902834008
negative vs positive: 10,thresholds:0.7
	Accuracy:0.9393939393939394
	Recall:0.8461538461538461
	AUC:0.9665991902834008
negative vs positive: 10,thresholds:0.9
	Accuracy:0.9212121212121213
	Recall:0.6153846153846154
	AUC:0.9665991902834008
'''

4.4 XGBoost类中的其他参数和功能

更多计算资源：n_jobs
nthread和n_jobs都是算法运行所使用的线程，与sklearn中规则一样，输入整数表示使用的线程，输入-1表示使用计算机全部的计算资源。如果数据量很大，则需要这个参数来调用更多线程。
降低学习难度：base_score
base_score是一个比较容易被混淆的参数，它叫做全局偏差，在分类问题中，它是希望关注的分类的先验概率。比如，如果有1000个样本，其中300个正样本，700个负样本，则base_score就是0.3。对于回归来说，这分数默认0.5，但其实这个分数在这种情况下并不有效。许多使用XGBoost的人已经提出，当使用回归的时候base_score的默认应该是标签的均值，不过现在xgboost库尚未对比做出改进。使用这个参数，便是在告诉模型一些我们已经了解但模型不一定能够从数据中学习到的信息。通常不会使用这个参数，但对于严重的样本不平衡问题，设置一个正确的base_score取值是很有必要的。
生成树的随机模式：random_state
在xgb库和sklearn中，都存在控制生成树的随机模式的参数random_state。在剪枝中，提到可以通过随机抽样样本、随机抽样特征来减轻过拟合的影响，可以通过其他参数来影响随机抽样的比例，却无法对随机抽样干涉更多。因此，真正的随机性还是由模型自己生成的。如果希望控制这种随机性，可以在random_state参数中输入固定整数。需要注意的是，xgb库和sklearn库中，在random_state参数中输入同一个整数未必表示同一个随机模式，不一定会得到相同的结果，因此导致模型的feature_importances也会不一致。
自动处理缺失值：missing
XGBoost被设计成是能够自动处理缺失值的模型，这个设计的初衷其实是为了让XGBoost能够处理稀疏矩阵，可以在参数missing中输入一个对象，比如np.nan，或数据的任意取值，表示将所有含有这个对象的数据作为空值处理。XGBoost会将所有的空值当作稀疏矩阵中的0来进行处理，因此在使用XGB

XGBoost结语

作为工程能力强，效果优秀的算法，XGBoost应用广泛并且原理复杂。可以意识到，XGBoost的难点在于它是一个集大成的模型。它所涉及的知识点和模型流程，多半在其他常用的机器学习模型中出现过：如树的迭代过程，其实可以和逻辑回归中的梯度下降过程进行类比；如剪枝的过程，许多都可以与随机森林和决策树中的知识进行对比。当然，XGBoost还有很多可以进行探索的，能够使用XGB算法的库和模块也不止sklearn和xgboost，许多库对xgboost的底层原理进行了更多的优化，让它变得更快更优秀（如lightBGM、使用分布式进行计算等等）。

参考资料：https://www.bilibili.com/video/BV1P7411P78r?p=258&spm_id_from=333.788.top_right_bar_window_history.content.click