机器学习基础—集成学习Task3

导语：
今天是集成学习的第三次任务，在上次的任务中，针对Boston数据集，我们给出了三种可用于回归的基础模型：线性回归、决策树、支持向量回归。但都使用的是默认参数，没有任何优化处理。那么，本次内容的重点是，如何优化这些基础模型？怎样评价模型的泛化能力？怎样选择特征？首先，附上学习链接：
集成学习: EnsembleLearning项目-github.

1.优化基础模型

了解机器学习基础的同学都知道，我们的数据集一般可分为训练集和测试集。有时候我们会发现，我们的模型在训练集上能表现得很好，误差很小，但在测试集上的误差却很大！为什么？我们在模型训练过程中，是否只关注训练集，尽量使得其误差最小就ok？
当然不是！我们建立机器学习的目的并不是为了在已有的数据集，也就是训练集上效果表现非常优异，我们希望建立的机器学习模型在未知且情况复杂的测试数据上表现优异，我们称这样的未出现在训练集的未知数据集成为测试数据集，简称测试集。我们希望模型在测试集上表现优异！这也正是我们优化这些基础模型的目的，提高泛化能力！

1.1 训练均方误差与测试均方误差

在上一任务中，我们都知道了“均方误差MSE”这一评价指标。显然，我们训练的过程就是让训练集的MSE越来越小，即模型在训练集的表现越来越好，从拟合的角度上来说，就是使得模型越来越接近真实数据。然而，前面也提到，一般而言，我们并不关心模型在训练集上的训练均方误差，我们关心的是模型面对未知的样本集，即测试集上的测试误差，我们的目标是使得我们建立的模型在测试集上的测试误差最小。
可是，测试集往往是我们在训练时不知道的未知数据，要想使模型的测试MSE最小不是件易事。是不是训练集MSE最小，测试集的MSE也最小呢？事实上却不是如此！

如上图所示，灰色线是训练集，红色线是测试集。
可以看到：当我们的模型的训练均方误差达到很小时，测试均方误差反而很大，但是我们寻找的最优的模型是测试均方误差达到最小时对应的模型，因此基于训练均方误差达到最小选择模型本质上是行不同的。正如上右图所示：模型在训练误差很小，但是测试均方误差很大时，我们称这种情况叫模型的过拟合。

1.2 偏差-方差的权衡

从上图的测试均方误差曲线可以看到：测试均方误差曲线呈现U型曲线，这表明了在测试误差曲线中有两种力量在互相博弈。可以证明：

也就是说，我们的测试均方误差的期望值可以分解为$f(x_0)$ 的方差、$f(x_0)$的偏差平方和误差项𝜖 的方差。为了使得模型的测试均方误差达到最小值，也就是同时最小化偏差的平方和方差。由于我们知道偏差平方和方差本身是非负的，因此测试均方误差的期望不可能会低于误差的方差，因此我们称 Var(𝜀) 为建模任务的难度，这个量在我们的任务确定后是无法改变的，也叫做不可约误差。
所谓模型的方差就是：用不同的数据集去估计 𝑓 时，估计函数的改变量。模型的偏差是指：为了选择一个简单的模型去估计真实函数所带入的误差。
偏差度量了学习算法的期望预测与真实结果的偏离程度，即刻画了学习算法本身的拟合能力。偏差度量的是单个模型的学习能力，而方差度量的是同一个模型在不同数据集上的稳定性。“偏差-方差分解”说明：泛化性能是由学习算法的能力、数据的充分性以及学习任务本身的难度所共同决定的。给定学习任务，为了取得好的泛化性能，则需使偏差较小，即能够充分拟合数据，并且使方差较小，即使得数据扰动产生的影响小。
一般而言，增加模型的复杂度，会增加模型的方差，但是会减少模型的偏差，我们要找到一个方差–偏差的权衡，使得测试均方误差最小。

1.3 特征提取

我们要选择一个测试误差达到最小的模型。但是实际上我们很难对实际的测试误差做精确的计算，因此我们要对测试误差进行估计，估计的方式有两种：训练误差修正与交叉验证。
（1）训练误差修正
模型越复杂，训练误差越小，测试误差先减后增。因此，我们先构造一个特征较多的模型使其过拟合，此时训练误差很小而测试误差很大，那这时我们加入关于特征个数的惩罚。因此，当我们的训练误差随着特征个数的增加而减少时，惩罚项因为特征数量的增加而增大，抑制了训练误差随着特征个数的增加而无休止地减小。具体的数学量如下：

（2）交叉验证
训练误差修正得到测试误差的估计是间接方法，这种方法的桥梁是训练误差，而交叉验证则是对测试误差的直接估计。交叉验证比训练误差修正的优势在于：能够给出测试误差的一个直接估计。在这里只介绍K折交叉验证：我们把训练样本分成K等分，然后用K-1个样本集当做训练集，剩下的一份样本集为验证集去估计由K-1个样本集得到的模型的精度，这个过程重复K次取平均值得到测试误差的一个估计。
以决策树回归为例，进行10折交叉验证：

#k折交叉验证
from sklearn.tree import DecisionTreeRegressor    
from sklearn.model_selection import cross_val_score

reg_tree = DecisionTreeRegressor(criterion = "mse",min_samples_leaf = 5)
scores = cross_val_score(reg_tree, X_train,y_train, cv=10)
print(scores)                                             
print("Accuracy: %0.2f (+/- %0.2f)" % (scores.mean(), scores.std() * 2))

输出结果：

可以看到，10折交叉验证的平均值为0.71，这和测试集的表现是很接近的！
ps：注意这里的数据是训练集的数据，经过随机打乱的，如果使用全集数据且有明显的排序，那么交叉验证的结果会很差！甚至可能会出现负分数！原因不言自明。

在测试误差能够被合理的估计出来以后，我们做特征选择的目标就是：从p个特征中选择m个特征，使得对应的模型的测试误差的估计最小。对应的方法有：最优子集选择和向前逐步选择。

1.4 压缩估计(正则化)

除了刚刚讨论的直接对特征自身进行选择以外，我们还可以对回归的系数进行约束或者加罚的技巧对p个特征的模型进行拟合，显著降低模型方差，这样也会提高模型的拟合效果。具体来说，就是将回归系数往零的方向压缩，这也就是为什么叫压缩估计的原因了。
（1）岭回归(L2正则化的例子)：

调节参数 𝜆 的大小是影响压缩估计的关键， 𝜆 越大，惩罚的力度越大，系数则越趋近于0，反之，选择合适的 𝜆 对模型精度来说十分重要。岭回归通过牺牲线性回归的无偏性降低方差，有可能使得模型整体的测试误差较小，提高模型的泛化能力。
（2）Lasso回归(L1正则化的例子)：
岭回归的一个很显著的特点是：将模型的系数往零的方向压缩，但是岭回归的系数只能呢个趋于0但无法等于0，换句话说，就是无法做特征选择。能否使用压缩估计的思想做到像特征最优子集选择那样提取出重要的特征呢？答案是肯定的！我们只需要对岭回归的优化函数做小小的调整就行了，我们使用系数向量的L1范数替换岭回归中的L2范数：

下图说明了为什么Losso能做到特征选择而岭回归却不能(如图：左边lasso，右边岭回归)

椭圆形曲线为RSS等高线，菱形和圆形区域分别代表了L1和L2约束，Lsaao回归和岭回归都是在约束下的回归，因此最优的参数为椭圆形曲线与菱形和圆形区域相切的点。但是Lasso回归的约束在每个坐标轴上都有拐角，因此当RSS曲线与坐标轴相交时恰好回归系数中的某一个为0，这样就实现了特征提取。反观岭回归的约束是一个圆域，没有尖点，因此与RSS曲线相交的地方一般不会出现在坐标轴上，因此无法让某个特征的系数为0，因此无法做到特征提取。

1.5 降维

到目前为止，我们所讨论的方法对方差的控制有两种方式：一种是使用原始变量的子集，另一种是将变量系数压缩至零。但是这些方法都是基于原始特征 𝑥1,…,𝑥𝑝 得到的，现在我们探讨一类新的方法：将原始的特征空间投影到一个低维的空间实现变量的数量变少，如：将二维的平面投影至一维空间。通过降维,我们希望减少冗余信息 所造成的误差,提高识别（或其他应用）的精度。又或者希望通过降维算法来寻找数据内部的本质结构特征。在很多算法中，降维算法成为了数据预处理的一部分，如PCA。事实上，有一些算法如果没有降维预处理，其实是很难得到很好的效果的。
主成分分析(PCA)：
主成分分析的思想：通过最大投影方差 将原始空间进行重构，即由特征相关重构为无关，即落在某个方向上的点(投影)的方差最大。在进行下一步推导之前，我们先把样本均值和样本协方差矩阵推广至矩阵形式：
样本均值Mean:$\bar{x}=\frac{1}{N}\sum _{i=1}^N{x}_i=\frac{1}{N}{X}^T{1}_N,其中1_N=(1,1,...,1{)}_N^T$
样本协方差矩阵$S^2=\frac{1}{N}\sum _{i=1}^N(x_i-\bar{x})(x_i-\bar{x}{)}^T=\frac{1}{N}{X}^{T}HX,其中，H=I_N-\frac{1}{N}{1}_N{1}_N^T$
最大投影方差的步骤：
(i) 中心化：$x_i-\bar{x}$
(ii) 计算每个点$x_1,...,x_N$至${\vec{u}}_1$方向上的投影：$(x_i-\bar{x}){\vec{u}}_1,||{\vec{u}}_1||=1$
(iii) 计算投影方差：$J=\frac{1}{N}\sum _{i=1}^N[(x_i-\bar{x}{)}^T{\vec{u}}_1{]}^2,||{\vec{u}}_1||=1$
(iv) 最大化投影方差求${\vec{u}}_1$：
$$\begin{gathered} {\bar{u}}_1=argma{x}_{u_1}\frac{1}{N}\sum _{i=1}^N[(x_i-\bar{x}{)}^T{\vec{u}}_1{]}^2 \\ s.t.{\vec{u}}_1^T{\vec{u}}_1=1({\vec{u}}_1往后不带向量符号) \end{gathered}$$
得到：
$$\begin{gathered} J=\frac{1}{N}\sum _{i=1}^N[(x_i-\bar{x}{)}^T{\vec{u}}_1{]}^2=\frac{1}{N}\sum _{i=1}^N[u_1^T(x_i-\bar{x})(x_i-\bar{x}{)}^T{u}_1] \\ =u_1^T[\frac{1}{N}\sum _{i=1}^N(x_i-\bar{x})(x_i-\bar{x}{)}^T]u_1=u_1^T{S}^2{u}_1 \end{gathered}$$
即：
$$\begin{gathered} {\hat{u}}_1=argma{x}_{u_1}{u}_1^T{S}^2{u}_1,s.t.u_1^T{u}_1=1 \\ L(u_1,\lambda )=u_1^T{S}^2{u}_1+\lambda (1-u_1^T{u}_1) \\ \frac{\partial L}{\partial u_1}=2S^2{u}_1-2\lambda u_1=0 \\ 即：S^2{u}_1=\lambda u_1 \end{gathered}$$
可以看到：$\lambda$为$S^2$的特征值，$u_1$为$S^2$的特征向量。因此我们只需要对中心化后的协方差矩阵进行特征值分解，得到的特征向量即为投影方向。如果需要进行降维，那么只需要取p的前M个特征向量即可。

2.向前逐步回归

特征提取的实例：向前逐步回归
案例来源：向前逐步回归
根据AIC准则定义向前逐步回归进行变量筛选

#定义向前逐步回归函数
def forward_select(data,target):
    variate=set(data.columns)  #将字段名转换成字典类型
    variate.remove(target)  #去掉因变量的字段名
    selected=[]
    current_score,best_new_score=float('inf'),float('inf')  #目前的分数和最好分数初始值都为无穷大（因为AIC越小越好）
    #循环筛选变量
    while variate:
        aic_with_variate=[]
        for candidate in variate:  #逐个遍历自变量
            formula="{}~{}".format(target,"+".join(selected+[candidate]))  #将自变量名连接起来
            aic=ols(formula=formula,data=data).fit().aic  #利用ols训练模型得出aic值
            aic_with_variate.append((aic,candidate))  #将第每一次的aic值放进空列表
        aic_with_variate.sort(reverse=True)  #降序排序aic值
        best_new_score,best_candidate=aic_with_variate.pop()  #最好的aic值等于删除列表的最后一个值，以及最好的自变量等于列表最后一个自变量
        if current_score>best_new_score:  #如果目前的aic值大于最好的aic值
            variate.remove(best_candidate)  #移除加进来的变量名，即第二次循环时，不考虑此自变量了
            selected.append(best_candidate)  #将此自变量作为加进模型中的自变量
            current_score=best_new_score  #最新的分数等于最好的分数
            print("aic is {},continuing!".format(current_score))  #输出最小的aic值
        else:
            print("for selection over!")
            break
    formula="{}~{}".format(target,"+".join(selected))  #最终的模型式子
    print("final formula is {}".format(formula))
    model=ols(formula=formula,data=data).fit()
    return(model)


import statsmodels
import statsmodels.api as sm #最小二乘
from statsmodels.formula.api import ols #加载ols模型
forward_select(data=boston_data,target="Price")

lm=ols("Price~LSTAT+RM+PTRATIO+DIS+NOX+CHAS+B+ZN+CRIM+RAD+TAX",data=boston_data).fit()
lm.summary()

结果如下：


可以看出，通过向前逐步回归挑选出了11个特征，而原数据集特征共有13个，其中，特征AGE和INDUS被剔除。在去除这两个特征之外，模型仍能保持0.741的$R^2$分数！

3.岭回归与Lasso回归

3.1 岭回归

岭回归实例分享：
*sklearn.linear_model.ridge_regression(X, y, alpha, , sample_weight=None, solver=‘auto’, max_iter=None, tol=0.001, verbose=0, random_state=None, return_n_iter=False, return_intercept=False, check_input=True)
https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.ridge_regression.html?highlight=rid#sklearn.linear_model.ridge_regression

参数：
alpha：较大的值表示更强的正则化。浮点数
sample_weight：样本权重，默认无。
solver：求解方法，{‘auto’, ‘svd’, ‘cholesky’, ‘lsqr’, ‘sparse_cg’, ‘sag’, ‘saga’}, 默认=’auto’。
“ svd”使用X的奇异值分解来计算Ridge系数。'cholesky’使用标准的scipy.linalg.solve函数通过dot（XT，X）的Cholesky分解获得封闭形式的解。'sparse_cg’使用scipy.sparse.linalg.cg中的共轭梯度求解器。作为一种迭代算法，对于大规模数据（可能设置tol和max_iter），此求解器比“ Cholesky”更合适。 lsqr”使用专用的正则化最小二乘例程scipy.sparse.linalg.lsqr。它是最快的，并且使用迭代过程。“ sag”使用随机平均梯度下降，“ saga”使用其改进的无偏版本SAGA。两种方法都使用迭代过程，并且当n_samples和n_features都很大时，通常比其他求解器更快。请注意，只有在比例大致相同的要素上才能确保“ sag”和“ saga”快速收敛。您可以使用sklearn.preprocessing中的缩放器对数据进行预处理。最后五个求解器均支持密集和稀疏数据。但是，当fit_intercept为True时，仅’sag’和’sparse_cg’支持稀疏输入。

from sklearn import linear_model
reg_rid = linear_model.Ridge(alpha=.5)
reg_rid.fit(X,y)
reg_rid.score(X,y)

输出结果：0.739957023371629
此时，去除特征AGE和INDUS：

from sklearn import linear_model
reg_rid = linear_model.Ridge(alpha=.5)
reg_rid.fit(boston_data[['CRIM', 'ZN', 'CHAS', 'NOX', 'RM', 'DIS', 'RAD',
       'TAX', 'PTRATIO', 'B', 'LSTAT']],y)
reg_rid.score(boston_data[['CRIM', 'ZN', 'CHAS', 'NOX', 'RM', 'DIS', 'RAD',
       'TAX', 'PTRATIO', 'B', 'LSTAT']],y)

输出结果：0.7399738198389612
可以看出，在特征减少的情况下，精度反而有所上升！

3.2 Lasso回归

Lasso实例分享：
class sklearn.linear_model.Lasso(alpha=1.0, *, fit_intercept=True, normalize=False, precompute=False, copy_X=True, max_iter=1000, tol=0.0001, warm_start=False, positive=False, random_state=None, selection=‘cyclic’)
https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.Lasso.html?highlight=lasso#sklearn.linear_model.Lasso

参数：
alpha：正则化强度，1.0代表标准最小二乘。
fit_intercept：是否计算模型截距。默认true。
normalize：是否标准化，默认false。
positive：是否强制系数为正，默认false。

from sklearn import linear_model
reg_lasso = linear_model.Lasso(alpha = 0.5)
reg_lasso.fit(X,y)
reg_lasso.score(X,y)

输出结果：0.7140164719858566
此时，去除特征AGE和INDUS：

from sklearn import linear_model
reg_lasso = linear_model.Lasso(alpha = 0.5)
reg_lasso.fit(boston_data[['CRIM', 'ZN', 'CHAS', 'NOX', 'RM', 'DIS', 'RAD',
       'TAX', 'PTRATIO', 'B', 'LSTAT']],y)
reg_lasso.score(boston_data[['CRIM', 'ZN', 'CHAS', 'NOX', 'RM', 'DIS', 'RAD',
       'TAX', 'PTRATIO', 'B', 'LSTAT']],y)

输出结果：0.7143260947435388
同样地，在特征减少的情况下，精度反而有所上升！
由此可见特征选择的重要性！