机器学习实战_线性回归&逻辑回归（二）

摘要：岭回归正则岭回归也称为正则化是线性回归的正则化版注意到这个正则项只有在训练过程中才会被加到代价函数。如果那此时的岭回归便变为了线性回归。

正如我们在第一和第二章看到的那样，降低模型的过拟合的好方法是正则化这个模型（即限制它）：模型有越少的自由度，就越难以拟合数据。例如，正则化一个多项式模型，一个简单的方法就是减少多项式的阶数。

对于一个线性模型，正则化的典型实现就是约束模型中参数的权重。 接下来我们将介绍三种不同约束权重的方法：Ridge回归，Lasso回归和Elastic Net。

岭回归（也称为Tikhonov正则化）是线性回归的正则化版：注意到这个正则项只有在训练过程中才会被加到代价函数。当得到完成训练的模型后，我们应该使用没有正则化的测量方法去评价模型的表现

一般情况下，训练过程使用的代价函数和测试过程使用的评价函数不一样样的。除了正则化，还有一个不同：训练时的代价函数应该在优化过程中易于求导，而在测试过程中，评价函数更应该接近最后的客观表现。一个好的例子：在分类训练中我们使用对数损失（马上我们会讨论它）作为代价函数，但是我们却使用精确率/召回率来作为它的评价函数。

岭回归(L2正则)代价函数:
$$J( heta)=MSE( heta)+alphafrac{1}{2}sumlimits_{i=1}^n heta_i^2$$

超参数$alpha$ 决定了你想正则化这个模型的强度,正则化强度越大，模型会越简单。如果$alpha$=0 那此时的岭回归便变为了线性回归。如果alpha非常的大，所有的权重最后都接近与零，最后结果将是一条穿过数据平均值的水平直线

(公式4-6：代价函数的梯度向量，加上$alpha w$是$1/2*alpha*w^2$求偏导的结果 )

在使用岭回归前，对数据进行放缩（可以使用StandardScaler）是非常重要的,算法对于输入特征的数值尺度（scale）非常敏感。大多数的正则化模型都是这样的。

下图展示了在相同线性数据上使用不同$alpha$ 值的岭回归模型最后的表现。左图中，使用简单的岭回归模型，最后得到了线性的预测。右图中的数据首先使用10 阶的PolynomialFearuresj进行扩展，然后使用StandardScaler进行缩放，最后将岭模型应用在处理过后的特征上。这就是带有岭正则项的多项式回归。注意当$alpha$ 增大的时候，导致预测曲线变得扁平（即少了极端值，多了一般值），这样减少了模型的方差，去增加了模型的偏差

对线性回归来说，对于岭回归，我们可以使用封闭方程去计算，也可以使用梯度下降去处理。它们的缺点和优点是一样的。公式4-9表示封闭方程的解（矩阵$mathbf{A}$ 是一个除了左上角有一个0的$n imes n$的单位矩，这个0代表偏差项。译者注：偏差$ heta_0$ 不被正则化的）。

岭回归的封闭方程的解(公式4-6：代价函数的梯度向量，加上$alpha w$是$1/2*alpha*w^2$求偏导的结果 ):
$$hat{ heta} =({mathbf{X}}^Tcdotmathbf{X}+alphamathbf{A})^{-1}cdot{mathbf{X}}^Tcdotmathbf{y}$$

下面是如何使用Scikit-Learn来进行封闭方程的求解（使用Cholesky法进行矩阵分解对上面公式行变形）

>>> from sklearn.linear_model import Ridge
>>> ridge_reg = Ridge(alpha=1, solver="cholesky")
>>> ridge_reg.fit(X, y)
>>> ridge_reg.predict([[1.5]])
array([[ 1.55071465]]

使用随机梯度法进行求解：

>>> sgd_reg = SGDRegressor(penalty="l2")
>>> sgd_reg.fit(X, y.ravel())
>>> sgd_reg.predict([[1.5]])
array([[ 1.13500145]])

penalty参数指的是正则项的惩罚类型。指定“l2”表明你要在代价函数上添加一项：权重向量$ell_2$ 范数平方的一半，这就是简单的岭回归。

Lasso回归（也称Least Absolute Shrinkage，或者Selection Operator Regression）是另一种正则化版的线性回归：就像岭回归那样，它也在代价函数上添加了一个正则化项，但是它使用权重向量的$ell_1$ 范数而不是权重向量$ell_2$ 范数平方的一半。

Lasso回归的代价函数:
$$J( heta)=MSE( heta)+alphasumlimits_{i=1}^nleft| heta_i ight|$$

下图展示了和之前相同的事情，仅仅是用Lasso模型代替了Ridge模型，同时调小了$alpha$ 的值

Lasso回归的一个重要特征是它倾向于完全消除最不重要的特征的权重（即将它们设置为零）。例如，右图中的虚线所示（$alpha=10^{-7}$ ），曲线看起来像一条二次曲线，而且几乎是线性的，这是因为所有的高阶多项特征都被设置为零

下面是一个使用Lasso类的小Scikit-Learn示例。你也可以使用SGDRegressor(penalty="l1")来代替它

>>> from sklearn.linear_model import Lasso
>>> lasso_reg = Lasso(alpha=0.1)
>>> lasso_reg.fit(X, y)
>>> lasso_reg.predict([[1.5]])
array([ 1.53788174]

弹性网络介于Ridge回归和Lasso回归之间。它的正则项是Ridge回归和Lasso回归正则项的简单混合，同时你可以控制它们的混合率r，当r=0时，弹性网络就是Ridge回归，当r=1时，其就是Lasso回归

弹性网络代价函数：

$$J( heta)=MSE( heta)+ralphasumlimits_{i=1}^nleft| heta_i ight|+frac{1-r}{2}alphasumlimits_{i=1}^n heta_i^2$$

那么我们该如何选择线性回归，岭回归，Lasso回归，弹性网络呢？一般来说有一点正则项的表现更好，因此通常你应该避免使用简单的线性回归。岭回归是一个很好的首选项，但是如果你的特征仅有少数是真正有用的，你应该选择Lasso和弹性网络。就像我们讨论的那样，它两能够将无用特征的权重降为零。一般来说，弹性网络的表现要比Lasso好，因为当特征数量比样例的数量大的时候，或者特征之间有很强的相关性时，Lasso可能会表现的不规律。下面是一个使用Scikit-Learn 弹性网络ElasticNet（l1_ratio指的就是混合率r）的简单样例:

>>> from sklearn.linear_model import ElasticNet
>>> elastic_net = ElasticNet(alpha=0.1, l1_ratio=0.5)
>>> elastic_net.fit(X, y)
>>> elastic_net.predict([[1.5]])
array([ 1.54333232])

随着训练的进行，算法一直学习，它在训练集上的预测误差（RMSE）自然而然的下降。然而一段时间后，验证误差停止下降，并开始上升。这意味着模型在训练集上开始出现过拟合。一旦验证错误达到最小值，便提早停止训练.

随机梯度和小批量梯度下降不是平滑曲线，你可能很难知道它是否达到最小值。 一种解决方案是，只有在验证误差高于最小值一段时间后（你确信该模型不会变得更好了），才停止，之后将模型参数回滚到验证误差最小值。

下面是一个早期停止法的基础应用

from sklearn.base import clone
sgd_reg = SGDRegressor(n_iter=1, warm_start=True, penalty=None,learning_rate="constant", eta0=0.0005)

minimum_val_error = float("inf")
best_epoch = None
best_model = None
for epoch in range(1000):
    sgd_reg.fit(X_train_poly_scaled, y_train)    # 训练多项式的新特征，拟合非线性
    y_val_predict = sgd_reg.predict(X_val_poly_scaled)
    val_error = mean_squared_error(y_val_predict, y_val)
    if val_error < minimum_val_error:
        minimum_val_error = val_error
        best_epoch = epoch
        best_model = clone(sgd_reg)

注意：当warm_start=True时，调用fit()方法后，训练会从停下来的地方继续，而不是从头重新开始

Logistic回归模型计算输入特征的加权和（加上偏差项），但它不像线性回归模型那样直接输出结果，而是把结果输入logistic()函数进行二次加工后进行输出

逻辑回归模型的概率估计（向量形式）：

$$hat{p}=h_ heta(mathbf{x})=sigma( heta^T cdot mathbf{x})$$

Logistic函数（也称为logit），用$sigma()$ 表示，其是一个sigmoid函数（图像呈S型），它的输出是一个介于0和1之间的数字
逻辑函数(S函数)

$$sigma(t)=frac{1}{1+exp(-t)}$$

逻辑回归预测模型($sigma()$ 概率输出以0.5作为二分类门槛):

单个样例的代价函数:

这个代价函数是合理的，因为当t接近0时，-log(t)变得非常大，所以如果模型估计一个正例概率接近于0，那么代价函数将会很大，同时如果模型估计一个负例的概率接近1，那么代价函数同样会很大。 另一方面，当t接近于1时， -log(t)接近0，所以如果模型估计一个正例概率接近于0，那么代价函数接近于0，同时如果模型估计一个负例的概率接近0，那么代价函数同样会接近于0， 这正是我们想的.（简单来说,y=1时，概率p越接近1损失越小；相反y=0时，概率p越接近0时损失越小）

整个训练集的代价函数只是所有训练实例的平均值。可以用一个表达式（你可以很容易证明）来统一表示，称为对数损失

逻辑回归的代价函数（对数损失）：

$$J( heta)=-frac{1}{m}sumlimits_{i=1}^mleft[y^{(i)}logleft(hat{p}^{(i)} ight)+left(1-y^{(i)} ight)logleft(1-hat{p}^{(i)} ight) ight]$$

但是这个代价函数对于求解最小化代价函数的$ heta$ 是没有公式解的（没有等价的正态方程）。 但好消息是，这个代价函数是凸的，所以梯度下降（或任何其他优化算法）一定能够找到全局最小值（如果学习速率不是太大，并且你等待足够长的时间）。下面公式给出了代价函数关于第j个模型参数$ heta_j$ 的偏导数

逻辑回归代价函数的偏导数:

$$frac{partial}{partial heta_j}J( heta_j)=frac{1}{m} sumlimits_{i=1}^m{left(sigmaleft( heta^T cdot mathbf{x}^{(i)} ight)-y^{(i)} ight)}{x_j}^{(i)}$$

这个公式首先计算每个样例的预测误差，然后误差项乘以第j项特征值，最后求出所有训练样例的平均值。 一旦你有了包含所有的偏导数的梯度向量，你便可以在梯度向量上使用批量梯度下降算法。 也就是说：你已经知道如何训练Logistic回归模型。 对于随机梯度下降，你当然只需要每一次使用一个实例，对于小批量梯度下降，你将每一次使用一个小型实例集。

我们使用鸢尾花数据集来分析Logistic回归。 这是一个著名的数据集，其中包含150朵三种不同的鸢尾花的萼片和花瓣的长度和宽度。这三种鸢尾花为：Setosa，Versicolor，Virginica

让我们尝试建立一个分类器，仅仅使用花瓣的宽度特征来识别Virginica，首先让我们加载数据：

>>> from sklearn import datasets
>>> iris = datasets.load_iris()
>>> list(iris.keys())
["data", "target_names", "feature_names", "target", "DESCR"]
>>> X = iris["data"][:, 3:] # petal width
>>> y = (iris["target"] == 2).astype(np.int)

接下来，我们训练一个逻辑回归模型：

from sklearn.linear_model import LogisticRegression

log_reg = LogisticRegression()
log_reg.fit(X, y) # 训练模型

我们来看看模型估计的花瓣宽度从0到3厘米的概率估计

X_new = np.linspace(0, 3, 1000).reshape(-1, 1)    # 构造花瓣宽度从0到3厘米的所有特征
y_proba = log_reg.predict_proba(X_new)    # 预测概率
plt.plot(X_new, y_proba[:, 1], "g-", label="Iris-Virginica")
plt.plot(X_new, y_proba[:, 0], "b--", label="Not Iris-Virginica"

Virginica花的花瓣宽度（用三角形表示）在1.4厘米到2.5厘米之间，而其他种类的花（由正方形表示）通常具有较小的花瓣宽度，范围从0.1厘米到1.8厘米。注意，它们之间会有一些重叠。在大约2厘米以上时，分类器非常肯定这朵花是Virginica花（分类器此时输出一个非常高的概率值），而在1厘米以下时，它非常肯定这朵花不是Virginica花（不是Virginica花有非常高的概率）。在这两个极端之间，分类器是不确定的。但是，如果你使用它进行预测（使用predict()方法而不是predict_proba()方法），它将返回一个最可能的结果。因此，在1.6厘米左右存在一个决策边界，这时两类情况出现的概率都等于50％：如果花瓣宽度大于1.6厘米，则分类器将预测该花是Virginica，否则预测它不是（即使它有可能错了）：

>>> log_reg.predict([[1.7], [1.5]])
array([1, 0])

下图的线性决策边界表示相同的数据集，但是这次使用了两个特征进行判断：花瓣的宽度和长度。 一旦训练完毕，Logistic回归分类器就可以根据这两个特征来估计一朵花是Virginica的可能性。 虚线表示这时两类情况出现的概率都等于50％：这是模型的决策边界。 请注意，它是一个线性边界。每条平行线都代表一个分类标准下的两两个不同类的概率，从15％（左下角）到90％（右上角）。越过右上角分界线的点都有超过90％的概率是Virginica花

就像其他线性模型，逻辑回归模型也可以$ell_1$或者$ell_2$ 惩罚使用进行正则化。Scikit-Learn默认添加了$ell_2$ 惩罚

在Scikit-Learn的LogisticRegression模型中控制正则化强度的超参数不是$alpha$ （与其他线性模型一样），而是是它的逆：C. C的值越大，模型正则化强度越低

Logistic回归模型可以直接推广到支持多类别分类，不必组合和训练多个二分类器， 其称为Softmax回归或多类别Logistic回归.

这个想法很简单：当给定一个实例$mathbf{x}$ 时，Softmax回归模型首先计算k类的分数$s_k(mathbf{x})$ ，然后将分数应用在Softmax函数（也称为归一化指数）上，估计出每类的概率。 计算$s_k(mathbf{x})$ 的公式看起来很熟悉，因为它就像线性回归预测的公式一样

k类的Softmax得分: $s_k(mathbf{x})= 	heta^T  cdot mathbf{x}$

注意，每个类都有自己独一无二的参数向量$ heta_k$ 。 所有这些向量通常作为行放在参数矩阵$Theta$ 中

一旦你计算了样例$mathbf{x}$ 的每一类的得分，你便可以通过Softmax函数估计出样例属于第k类的概率$hat{p}_k$ ：通过计算e的$s_k(mathbf{x})$ 次方，然后对它们进行归一化（除以所有分子的总和）。

和Logistic回归分类器一样，Softmax回归分类器将估计概率最高（它只是得分最高的类）的那类作为预测结果，如公式4-21所示

Softmax回归分类器一次只能预测一个类（即它是多类的，但不是多输出的），因此它只能用于判断互斥的类别，如不同类型的植物。 你不能用它来识别一张照片中的多个人。

现在我们知道这个模型如何估计概率并进行预测，接下来将介绍如何训练。我们的目标是建立一个模型在目标类别上有着较高的概率（因此其他类别的概率较低），最小化公式4-22可以达到这个目标，其表示了当前模型的代价函数，称为交叉熵，当模型对目标类得出了一个较低的概率，其会惩罚这个模型。 交叉熵通常用于衡量待测类别与目标类别的匹配程度（我们将在后面的章节中多次使用它）

交叉熵源于信息论。假设你想要高效地传输每天的天气信息。如果有八个选项（晴天，雨天等），则可以使用3位对每个选项进行编码，因为2^3=8。但是，如果你认为几乎每天都是晴天，更高效的编码“晴天”的方式是：只用一位（0）。剩下的七项使用四位（从1开始）。交叉熵度量每个选项实际发送的平均比特数。 如果你对天气的假设是完美的，交叉熵就等于天气本身的熵（即其内部的不确定性）。 但是，如果你的假设是错误的（例如，如果经常下雨）交叉熵将会更大，称为Kullback-Leibler散度（KL散度）。

现在你可以计算每一类的梯度向量，然后使用梯度下降（或者其他的优化算法）找到使得代价函数达到最小值的参数矩阵$Theta$ 。

让我们使用Softmax回归对三种鸢尾花进行分类。当你使用LogisticRregression对模型进行训练时，Scikit_Learn默认使用的是一对多模型，但是你可以设置multi_class参数为“multinomial”来把它改变为Softmax回归。你还必须指定一个支持Softmax回归的求解器，例如“lbfgs”求解器（有关更多详细信息，请参阅Scikit-Learn的文档）。其默认使用$ell_12$ 正则化，你可以使用超参数C控制它。

X = iris["data"][:, (2, 3)] # petal length, petal width
y = iris["target"]

softmax_reg = LogisticRegression(multi_class="multinomial",solver="lbfgs", C=10)
softmax_reg.fit(X, y)

所以下次你发现一个花瓣长为5厘米，宽为2厘米的鸢尾花时，你可以问你的模型你它是哪一类鸢尾花，它会回答94.2％是Virginica花（第二类），或者5.8％是其他鸢尾花

>>> softmax_reg.predict([[5, 2]])
array([2])
>>> softmax_reg.predict_proba([[5, 2]])
array([[ 6.33134078e-07, 5.75276067e-02, 9.42471760e-01]])是

图4-25用不同背景色表示了结果的决策边界。注意，任何两个类之间的决策边界是线性的。 该图的曲线表示Versicolor类的概率（例如，用0.450标记的曲线表示45％的概率边界）。注意模型也可以预测一个概率低于50％的类。 例如，在所有决策边界相遇的地方，所有类的估计概率相等，分别为33％。