跟我一起学scikit-learn17：逻辑回归算法

本文主要是介绍跟我一起学scikit-learn17：逻辑回归算法，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

逻辑回归算法的名字虽然带有“回归”二字，但实际上逻辑回归算法是用来解决分类问题的算法。本章首先从二元分类入手，介绍了逻辑回归算法的预测函数、成本函数和梯度下降算法公式；然后介绍了怎样由二元分类延伸至多元分类的问题；接着介绍了正则化，即通过模型的手段来解决模型的过拟合问题；针对正则化，介绍了L1范数和L2范数的含义及区别；最后用一个乳腺癌检测的实例及其模型性能优化来结束本章内容。

1.逻辑回归算法的原理

假设有一场足球赛，我们有两支球队的所有出场球员信息、历史交锋成绩、比赛时间、主客场、裁判和天气等信息，根据这些信息预测球队的输赢。假设比赛结果记为y，赢球标记为1，输球标记为0，这个就是典型的二元分类问题，可以用逻辑回归算法来解决。

从这个例子里可以看出，逻辑回归算法的输出是个离散值，这是与线性回归算法的最大区别。

1.预测函数

需要找出一个预测函数模型，使其值输出在[0,1]之间。然后选择一个基准值，如0.5，如果算出来的预测值大于0.5，就认为其预测值为1，反之，则其预测值为0.

选择
$g(z)=\frac{1}{1+{e}^{-z}}$
来作为预测函数，其中e是自然对数的底数。函数g(z)称为Sigmoid函数，也称为Logistic函数。以z为横坐标，以g(z)为纵坐标，画出的图形如下所示：

从图中可以看出，当z=0时，g(z)=0.5；当z>0时，g(z)>0.5，当z越来越大时，g(z)无限接近于1；当z<0时，g(z)<0.5，当z越来越小时，g(z)无限接近于0。这正是我们想要的针对二元分类算法的预测函数。

问题来了，怎样将输入特征和预测函数结合起来呢？
结合线性回归算法的预测函数 $h_{\theta}(x)={\theta}^Tx$ ，令 $z(x)={\theta}^Tx$ ，则逻辑回归算法的预测函数如下：
$h_{\theta}(x)=g(z)=g({\theta}^Tx)=\frac{1}{1+{e}^{-{\theta}^Tx}}$
下面来解读预测函数：
$h_{\theta}(x)$ 表示在输入值为x，参数为 $\theta$ 的前提条件下y=1的概率。用概率论的公式可以写成：
$h_{\theta}(x)=P(y=1|x,\theta)$
上面的概率公式可以读成：在输入x及参数 $\theta$ 的条件下y=1的概率，这是个条件概率公式。由概率论的知识可以推导出：
$P(y=1|x,\theta)+P(y=0|x,\theta)=1$
对二元分类来说，这是个非黑即白的世界。

2.判定边界

逻辑回归算法的预测函数由下面两个公式给出：
$h_{\theta}(x)=g(\theta^Tx)$
$g(z)=\frac{1}{1+e^{-z}}$
假定y=1的判定条件是 $h_{\theta(x)}>0.5$ ，y=0的判定条件是 $h_{\theta(x)}<0.5$ ，则可以推导出y=1的判定条件是 $\theta^Tx>0$ ，y=0的判定条件就是 $\theta^Tx<0$ 。所以， $\theta^Tx=0$ 即是我们的判定边界。

下面给出两个判定边界的例子。
假设有两个变量x1，x2，其逻辑回归预测函数是 $h_{\theta}(x)=g(\theta_{0}+\theta_{1}x_{1}+\theta_{2}x_{2})$ 。假设给定参数：
$\theta=\begin{bmatrix} -3\\ 1\\ 1\\ \end{bmatrix}$
那么，可以得到判定边界 $3+x_{1}+x_{2}=0$ ，即 $x_{1}+x_{2}=3$ ，如果以 $x_{1}$ 为横坐标， $x_{2}$ 为纵坐标，则这个函数画出来就是一条通过(0,3)和(3,0)两点的斜线。这条线就是判定边界，如图所示：

其中，直线左下方为y=0，直线右上方为y=1。横坐标为 $x_{1}$ ，纵坐标为 $x_{2}$ 。

如果预测函数是多项式 $h_{\theta}(x)=g(\theta_{0}+\theta_{1}x_{1}+\theta_{2}x_{2}+\theta_{3}x_{1}^{2}+\theta_{4}x_{2}^{2})$ ，且给定
$\theta=\begin{bmatrix} -1\\ 0\\ 0\\ 1\\ 1 \end{bmatrix}$
则可以得到判定边界函数 $x_{1}^{2}+x_{2}^{2}=1$ 。还是以 $x_{1}$ 为横坐标，以 $x_{2}$ 为纵坐标，则这是一个半径为1的圆。圆内部是y=0，圆外部是y=1，如上图所示。
这是二阶多项式的情况，更一般的多阶多项式可以表达出更复杂的判定边界。

3.成本函数

我们不能使用线性回归模型的成本函数来推导逻辑回归的成本函数，因为那样的成本函数太复杂，最终很可能会导致无法通过迭代找到成本函数值最小的点。

为了容易地求出成本函数的最小值，我们分成y=1和y=0两种情况来分别考虑其预测值和真实值的误差。我们先考虑最简单的情况，即计算某个样本x，y=1和y=0两种情况下的预测值与真实值的误差，我们选择的成本公式如下：
$Cost(h_{\theta}(x),y)=\left\{\begin{matrix} -log(h_{\theta}(x)), & if & y=1\\ -log(1-h_{\theta}(x)),& if & y=0 \end{matrix}\right.$
其中， $h_{\theta}(x)$ 表示预测为1的概率，log(x)为自然对数。我们以 $h_{\theta}(x)$ 为横坐标，以成本值 $Cost(h_{\theta}(x),y)$ 为纵坐标，把上述两个公式分别画在二维平面上，如图所示：

回顾成本函数的定义，成本是预测值与真实值的差异。当差异越大时，成本越大，模型受到的“惩罚”也越严重。
在左图中，当y=1时，随着 $h_{\theta}(x)$ 的值（预测为1的概率）越来越大，预测值越来越接近真实值，其成本越来越小；在右图中，当y=0时，随着 $h_{\theta}(x)$ 的值（预测为1的概率）越来越大，预测值越来越偏离真实值，其成本越来越大。
思考：符合上述规律的函数模型有很多，为什么我们要选择自然对数函数来作为成本函数呢？
逻辑回归模型的预测函数是Sigmoid函数，而Sigmoid函数里有e的n次方运算，自然对数刚好是其逆运算，比如 $log(e^n)=n$ 。选择自然对数，最终会推导出形式优美的逻辑回归模型参数的迭代函数，而不需要去涉及对数运算和指数函数运算。这就是我们选择自然对数函数来作为成本函数的原因。更进一步，把输入值x从负无穷大到正无穷大映射到[0,1]区间的模型有很多，逻辑回归算法为什么要选择Sigmoid函数作为预测函数的模型呢？严格地讲，不一定非要选择Sigmoid函数作为预测函数。但如果不选择Sigmoid函数作为预测函数，就需要重新选择性质接近的成本函数，这样才能在数学上得到既方便表达、效率又高的成本函数。
下面来看成本函数的统一写法问题。分开表述的成本计算公式始终不方便，能不能合并成一个公式呢？考虑下面的公式：
$Cost(h_{\theta}(x),y)= -ylog(h_{\theta}(x))-(1-y)log(1-h_{\theta}(x))$
由于y只能取0或1，所以上述式子和分开表达的公式是等价的。
介绍到这里，成本函数就要隆重登场了。根据一个样本的成本计算公式，很容易写出所有样本的成本平均值，即成本函数：
$J(\theta)=-\frac{1}{m}\begin{bmatrix}\sum_{i=1}^{m} -y^{(i)}log(h_{\theta}(x^{(i)}))-(1-y^{(i)})log(1-h_{\theta}(x^{(i)}))\end{bmatrix}$

4.梯度下降算法

和线性回归类似，我们使用梯度下降算法来求解逻辑回归模型参数。根据梯度下降算法的定义，可以得出：
$\theta_{j}=\theta_{j}-\alpha\frac{\partial}{\partial\theta_{j}}J(\theta)$
这里的关键是求解成本函数的偏导数。最终推导出来的梯度下降算法公式为：
$\theta_{j}=\theta_{j}-\alpha\frac{1}{m}\sum_{i=1}^{m}(h_{\theta}(x^{(i)})-y^{(i)})x_{j}^{(i)}$
如果对公式推导过程感兴趣的，可以参考本章扩展阅读的内容。
这个公式的形式和线性回归算法的参数迭代公式是一样的。当然，由于这里 $h_{\theta}(x)=\frac{1}{1+e^{-\theta^Tx}}$ ，而线性回归算法里的 $h_{\theta}(x)=\theta^Tx$ 。所以，两者的形式一样，而数值计算方法则完全不同。
至此，我们就把逻辑回归算法的相关原理解释清楚了。

2.多元分类

逻辑回归模型可以解决二元分类问题，即y={0,1}，能不能解决多元分类问题呢？答案是肯定的。针对多元分类问题，y={0,1,2,3,…,n}，总共有n+1个类别。其解决思路是：首先把问题转换为二元分类问题，即y=0是一个类别，y={1,2,3,…,n}作为另外一个类别，然后计算这两个类别的概率；接着，把y=1作为一个类别，把y={0,2,3,…,n}作为另外一个类别，再计算这两个类别的概率。由此推广开，总共需要n+1个预测函数：
$y\in{\{0,1,...,n\}}\\ h_{\theta}^{(0)}(x)=P(y=0|x,\theta)\\ h_{\theta}^{(1)}(x)=P(y=0|x,\theta)\\ ...\\ h_{\theta}^{(n)}(x)=P(y=0|x,\theta)\\ prediction=max_{i}(h_{\theta}^{(i)}(x))$
预测出来的概率最高的那个类别就是样本所属的类别。

3.正则化

回忆之前学过的知识：过拟合是指模型很好地拟合了训练样本，但对新数据预测的准确性很差，这是因为模型太复杂了。解决办法是减少输入特征的个数，或者获取更多的训练样本。这里介绍的正则化也是用来解决过拟合问题的一个方法。

保留所有的特征，减少特征的权重 $\theta_{j}$ 的值。确保所有的特征对预测值都有少量的贡献。
当每个特征 $x_j$ 对预测值y都有少量的贡献时，这样的模型可以良好的工作，这正是正则化的目的可以用它来解决特征过多时的过拟合问题。

1.线性回归模型正则化

我们先来看线性回归模型的成本函数是如何正则化的：
$J(\theta)=\frac{1}{2m}\sum_{i=1}^{m}(h_{\theta}(x^{(i)})-y^{(i)})^2+\lambda\sum_{j=1}^{n}\theta_{j}^{2}$
公式中前半部分就是原来的线性回归模型的成本函数，也称为预测值与实际值的误差。后半部分为加入的正则项。其中 $\lambda$ 的值有两个目的，即要维持对训练样本的拟合，又要避免对训练样本的过拟合。如果 $\lambda$ 的值太大，则能确保不出现过拟合，但可能会导致对现有训练样本出现欠拟合。
思考：怎样从数学上理解正则化后的成本函数解决了过拟合问题呢？
从数学角度看，成本函数增加了一个正则项后，成本函数不再唯一地由预测值和真实值的误差所决定，还和参数 $\theta$ 的大小有关。有了这个限制之后，要实现成本函数最小的目的， $\theta$ 就不能随便取值了。比如某个比较大的 $\theta$ 值可能会让预测值与真实值的误差 $(h_{\theta}(x^{(i)})-y^{(i)})^2$ 值很小，但会导致 $\theta_{j}^{2}$ 很大，最终的结果是成本函数很大。这样通过调节参数 $\lambda$ ，就可以控制正则项的权重，从而避免线性回归算法过拟合。
利用正则化的成本函数，可以推导出正则化后的参数迭代函数：
$\theta_j=\theta_j-\alpha\frac{1}{m}\sum_{i=1}^{m}\begin{bmatrix}(h(x^{(i)})-y^{(i)})x_j^{(i)}+\frac{\lambda}{m}\theta_j\end{bmatrix}$
$\theta_j=\theta_j(1-\alpha\frac{\lambda}{m})-\alpha\frac{1}{m}\sum_{i=1}^{m}(h(x^{(i)})-y^{(i)})x_j^{(i)}$
$(1-\alpha\frac{\lambda}{m})$ 因子在每次迭代时，都将把 $\theta_j$ 收缩一点点。因为 $\alpha$ 和 $\lambda$ 是正数，而m是训练样本的个数，是个比较大的正整数。为什么要对 $\theta_j$ 进行收缩呢？因为加入正则项的成本函数和 $\theta_j^2$ 成正比，所以迭代时需要不断地试图减小 $\theta_j$ 的值。

2.逻辑回归模型正则化

使用相同的思路，我们可以对逻辑回归模型的成本函数进行正则化，其方法也是在原来的成本函数的基础上加上正则项：
$J(\theta)=-\frac{1}{m}\begin{bmatrix} \sum_{i=1}^{m}y^{(i)}log(h_{\theta}(x^{(i)}))+(1-y^{(i)})log(1-h_{\theta}(x^{(i)})) \end{bmatrix}+\frac{\lambda}{2m}\sum_{j=1}^{n}\theta_j^2$
相应地，正则化后的参数迭代公式为：
$\theta_{j}=\theta_{j}-\alpha\frac{\partial}{\partial\theta_{j}}J(\theta)$
$\theta_j=\theta_j-\alpha \begin{bmatrix} \frac{1}{m}\sum_{i=1}^{m}(h(x^{(i)})-y^{(i)})x_j^{(i)}+\frac{\lambda}{m}\theta_j\end{bmatrix}$
$\theta_j=\theta_j(1-\alpha\frac{\lambda}{m})-\alpha\frac{1}{m}\sum_{i=1}^{m}(h(x^{(i)})-y^{(i)})x_j^{(i)}$
需要注意的是，上式中 $j\geqslant1$ ，因为 $\theta_0$ 没有参与正则化。另外，需要注意的是，逻辑回归和线性回归看起来形式是一样的，但其实它们的算法不一样，因为两个式子的预测函数 $h_{\theta}(x)$ 不一样。针对线性回归 $h_{\theta}(x)=\theta^Tx$ ，而针对逻辑回归 $h_{\theta}(x)=\frac{1}{1+e^{-\theta^Tx}}$ 。

4.算法参数

在scikit-learn里，逻辑回归模型由类sklearn.linear_model.LogisticRegression实现。

1.正则项权重

我们上面介绍的正则项权重 $\lambda$ ，在LogisticRegression里有个参数C与此对应，但成反比。即C值越大，正则项权重越小，模型容易出现过拟合；C值越小，正则项权重越大，模型容易出现欠拟合。

2.L1/L2范数

创建逻辑回归模型时，有个参数penalty，其取值有“l1”或“l2”，这是什么意思呢？
这个实际上就是指定我们前面介绍的正则项的形式。回顾之前介绍的内容，在成本函数里添加的正则项 $\sum_{j=1}^{n}\theta_j^2$ ，这个实际上是个L2正则项，即把L2范数作为正则项。当然，也可以添加L1范数来作为正则项。
L1范数作为正则项，会让模型参数 $\theta$ 稀疏化，即让模型参数向量里的0元素尽可能多，只保留模型参数向量中重要特征的贡献。而L2范数作为正则项，则让模型参数尽量小，但不会为0，即尽量让每个特征对应预测值都有一些小的贡献。
问题来了，为什么会造成上述不同的结果呢？
要解释清楚原因，就需要先了解一下L1范数和L2范数的概念，它们都是针对向量的一种运算。为了简单起见，假设模型只有两个参数，它们构成一个二维向量 $\theta=[\theta_1,\theta_2]$ ，则L1范数为：
$||\theta||_1=|\theta|_1+|\theta|_2$
即L1范数是向量里元素的绝对值之和。L2范数为元素平方和的平方根：
$||\theta||_2=\sqrt{\theta_1^2+\theta_2^2}$
定义清楚了之后，我们来介绍它们作为正则项的效果有什么不同。回顾之前介绍的内容，梯度下降算法在参数迭代的过程中，实际上是在成本函数的等高线上跳跃，并最终收敛在误差最小的点上（为了避免误解，此处成未加正则项前的成本为误差）。我们先思考一下，正则项的本质是什么？正则项的本质是惩罚。模型在训练的过程中，如果没有遵循正则项所表达的规则，那么其成本会变大，即受到了惩罚，从而往正则项所表达的规则处收敛。成本函数在这两项规则的综合作用下，正则化后的模型参数应该收敛在误差等高线与正则项等高线相切的点上。
我们把L1范数和L2范数在二维坐标轴上画出图形，即可直观地看到它们所表达的规则的不同。

%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
def L1(x):return 1 - np.abs(x)def L2(x):return np.sqrt(1 - np.power(x, 2))def contour(v, x):return 5 - np.sqrt(v - np.power(x + 2, 2))    # 4x1^2 + 9x2^2 = vdef format_spines(title):    ax = plt.gca()                                  # gca 代表当前坐标轴，即 'get current axis'ax.spines['right'].set_color('none')            # 隐藏坐标轴ax.spines['top'].set_color('none')ax.xaxis.set_ticks_position('bottom')           # 设置刻度显示位置ax.spines['bottom'].set_position(('data',0))    # 设置下方坐标轴位置ax.yaxis.set_ticks_position('left')ax.spines['left'].set_position(('data',0))      # 设置左侧坐标轴位置plt.title(title)plt.xlim(-4, 4)plt.ylim(-4, 4)plt.figure(figsize=(8.4, 4), dpi=144)x = np.linspace(-1, 1, 100)
cx = np.linspace(-3, 1, 100)plt.subplot(1, 2, 1)
format_spines('L1 norm')
plt.plot(x, L1(x), 'r-', x, -L1(x), 'r-')
plt.plot(cx, contour(20, cx), 'r--', cx, contour(15, cx), 'r--', cx, contour(10, cx), 'r--')plt.subplot(1, 2, 2)
format_spines('L2 norm')
plt.plot(x, L2(x), 'b-', x, -L2(x), 'b-')
plt.plot(cx, contour(19, cx), 'b--', cx, contour(15, cx), 'b--', cx, contour(10, cx), 'b--')

在左图中，我们用的是L1范数来作为正则项，L1范数表示的元素的绝对值之和，图中L1范数的值为1，其中 $\theta_{1}$ ， $\theta_{2}$ 坐标轴上的等值线是个正方形，虚线表示的是误差等值线。可以看到，误差等值线和L1范数等值线相切的点位于坐标轴上。
在右图中，我们用的是L2范数来作为正则项，图中L2范数的值为1，在 $\theta_{1}$ ， $\theta_{2}$ 坐标轴上，它的等值线是一个圆。它和模型误差等值线相切的点，一般不在坐标轴上。
至此，我们就清楚了，L1范数作为正则项，会让模型参数稀疏化；而L2范数作为正则项，则会使模型的特征对预测值都有少量的贡献，避免模型过拟合。
作为推论，L1范数作为正则项，有以下几个用途：

特征选择：它会让模型参数向量里的元素为0的点尽量多。因此可以排除掉那些对预测值没有什么影响的特征。从而简化问题。所以L1范数解决过拟合的措施，实际上是减少特征数量。
可解释性：模型参数向量稀疏化后，只会留下那些对预测值有重要影响的特征。这样我们就容易解释模型的因果关系。比如，针对某种癌症的筛查，如果有100个特征，那么我们无从解释到底哪些特征对阳性呈关键作用。稀疏化后，只留下几个关键的特征，就容易看到因果关系。
由此可见，L1范数作为正则项，更多的是一个分析工具，而适合用来对模型求解。因为它会把不重要的特征直接去除。大部分的情况下，我们解决过拟合问题，还是选择L2范数作为正则项，这也是scikit-learn里的默认值。

5.示例：乳腺癌检测

本节来看一个实例，使用逻辑回归算法解决乳腺癌检测问题。我们需要先采集肿瘤病灶造影图片，然后对图片进行分析，从图片中提取特征，再根据特征来训练模型。最终使用模型来检测新采集到的肿瘤病灶造影，一遍判断肿瘤是良性的还是恶性的。这是个典型的二元分类问题。

1.数据采集及特征提取

在工程应用中，数据采集和特征提取工作往往决定着项目的成败。我们可以思考一下，我们要获取足够多数量的良性和恶性分布合理的肿瘤病灶造影图片，需要多大的线下工作？当然，这一工作一般会以科研中心和医院合作的形式来获取数据。其次，拿到病灶造影图片后，要分析图片、提取特征，这也是个费心伤神的工作？这都是值得深入思考和实践的问题。最后，决定要提取的特征集合，还需要编写图片处理程序，以便从病灶造影图片中提取出我们需要的特征。
为了简单起见，我们直接加载scikit-learn自带的一个乳腺癌数据集。这个数据集是已经采集后的数据：

from sklearn.datasets import load_breast_cancer
cancer = load_breast_cancer()
X = cancer.data
y = cancer.target
print('data shape: {0}; no. positive: {1}; no. negative: {2}'.format(X.shape,y[y==1].shape[0],y[y==0].shape[0]))
print(cancer.data[0])

上述代码输出结果如下：

data shape: (569, 30); no. positive: 357; no. negative: 212
[1.799e+01 1.038e+01 1.228e+02 1.001e+03 1.184e-01 2.776e-01 3.001e-011.471e-01 2.419e-01 7.871e-02 1.095e+00 9.053e-01 8.589e+00 1.534e+026.399e-03 4.904e-02 5.373e-02 1.587e-02 3.003e-02 6.193e-03 2.538e+011.733e+01 1.846e+02 2.019e+03 1.622e-01 6.656e-01 7.119e-01 2.654e-014.601e-01 1.189e-01]

我们可以看到，数据集中总共有569个样本，每个样本有30个特征，其中357个阳性（y=1）样本，212个阴性（y=0）样本。同时，我们还打印出一个样本数据，以便直观地进行观察。
为了强调特征提取工作的重要性，这里介绍一下这些特征值的物理意义。这个数据集总共从病灶造影图片中提取了以下10个关键属性：

radius：半径，即病灶中心点离边界的平均距离。
texture：纹理，灰度值的标准偏差。
perimeter：周长，即病灶的大小。
area：面积，也是反映病灶大小的一个指标。
smoothness：平滑度，即半径的变化幅度。
compactness：密实度，周长的平方除以面积，再减去1， $\frac{perimeter^2}{area}-1$
concavity：凹度，凹陷部分轮廓的严重程度。
concave points：凹点，凹陷轮廓的数量。
symmetry：对称性。
fractal demension：分形维度。

从这些指标里可以看出，有些指标是“复合”指标，即由其他的指标经过计算得到的。比如密实度，是由周长和面积计算出来的。不要小看这种运算构建出来的新特征，这是事物内在逻辑关系的体现。
举个例子，我们需要监控数据中心中每台物理主机的运行情况，其中CPU占用率、内存占用率、网络吞吐量十几个重要的指标。问：有台主机CPU占用率80%，这个主机状态是否正常？要不要发布告警？答案：看情况。仅从CPU占用率来看还不能判断主机是否正常，还要看内存占用情况和网络吞吐量情况。如果此时内存占用也成比例上升，且网络吞吐量也在合理的水平，那么造成这一状态的可能是用户访问的流量过大，导致主机负荷增加，不需要告警。但如果内存占用、网络吞吐量和CPU占用不再同一量级，那么这台主机就可能处于不正常的状态。所以，我们需要构建一个复合特征，如CPU占用率和内存占用率的比值，以及CPU占用率和网络吞吐量的比值，这样构造出来的特征更真实地体现出了现实问题中的内在规则。
所以，提取特征时，不妨从事物的内在逻辑关系入手，分析已有特征之间的关系，从而构造出新的特征。这一方法在实际工程应用中是常用的特征提取手段。
回到我们讨论的乳腺癌数据集的特征问题中，实际上它只关注10个特征，然后又构造出了每个特征的标准差及最大值，这样每个特征就衍生出了两个特征，所以总共就有了30个特征。可以通过cancer.feature_names变量来查看这些特征的名称。

2.模型训练

首先，我们数据集分成训练数据集和测试数据集：

from sklearn.model_selection import train_test_split
X_train,X_test,y_train,y_test = train_test_split(X,y,test_size=0.2)

然后使用LogisticRegression模型来训练，并计算训练数据集的评分数据和测试数据集的评分数据：

from sklearn.linear_model import LogisticRegression
model = LogisticRegression()
model.fit(X_train,y_train)
train_score = model.score(X_train,y_train)
test_score = model.score(X_test,y_test)
print('train score: {train_score:.6f}; test_score:{test_score:.6f}'.format(train_score=train_score,test_score=test_score))

输出如下：

train score: 0.947253; test_score:0.973684

我们还可以看一下测试样本中，有几个是预测正确的：

import numpy as np
y_pred = model.predict(X_test)
print('matchs: {0}/{1}'.format(np.equal(y_pred,y_test).shape[0],y_test.shape[0]))

输出如下：

matchs: 114/114

总共114个测试样本，全部预测正确。这里有个疑问，为什么全部都预测正确，而testscore却只有0.973684，而不是1呢？答案是，scikit-learn不是使用这个数据来计算分数，因为这个数据不能完全反映误差情况，而是使用预测概率数据计算模型评分。
针对二元分类问题，LogisticRegression模型会对每个样本输出两个概率，即为0的概率和为1的概率，哪个概率高就预测为哪个类别。
我们可以找出针对测试数据集，模型预测的“自信度”低于90%的样本。怎样找出这些样本呢？我们先计算出测试数据集里的每个样本的预测概率数据，针对每个样本，它会有两个数据，一是预测其为阳性的概率，另外一个是预测其为阴性的概率。接着找出预测为阴性的概率大于0.1的样本，然后在结果集里，找出预测为阳性的概率也大于0.1的样本，这样就找出了模型预测“自信度”低于90%的样本。这是因为所有类别的预测概率加起来，一定是100%，两个都大于0.1，则其最大值一定是小于90%，即“自信度”不足90%。我们可以看下概率数据：

# 预测概率：找出预测概率低于 90% 的样本
y_pred_proba = model.predict_proba(X_test)     # 计算每个测试样本的预测概率
# 打印出第一个样本的数据，以便读者了解数据形式
print('sample of predict probability: {0}'.format(y_pred_proba[0]))# 找出第一列，即预测为阴性的概率大于 0.1 的样本，保存在 result 里
y_pred_proba_0 = y_pred_proba[:, 0] > 0.1
result = y_pred_proba[y_pred_proba_0]# 在 result 结果集里，找出第二列，即预测为阳性的概率大于 0.1 的样本
y_pred_proba_1 = result[:, 1] > 0.1
print(result[y_pred_proba_1])

输出结果如下：

sample of predict probability: [0.29623162 0.70376838]
[[0.29623162 0.70376838][0.54660262 0.45339738][0.17874247 0.82125753][0.20917573 0.79082427][0.10943452 0.89056548][0.35503614 0.64496386][0.23849987 0.76150013][0.13634228 0.86365772][0.80171734 0.19828266][0.21744759 0.78255241][0.81346356 0.18653644][0.2225791  0.7774209 ][0.10788007 0.89211993][0.88068005 0.11931995][0.18189724 0.81810276]]

我们使用model.predict_proba()来计算概率，同时找出那些预测“自信度”低于90%的样本。可以看到，最没有把握的样本是[0.35342587 0.64657413]，即只有64.66%的概率是阳性。大家如果运行这个实例，输出结果可能会略有差异，因为训练样本和测试样本是随机分配的。

3.模型优化

我们使用LogisticRegression模型的默认参数训练出来的模型，准确性看起来还是挺高的。问题是，有没有优化空间呢？如果有，往哪个方向优化呢？

我们先尝试增加多项式特征，实际上，多项式特征和上文介绍的人为添加的复合特征类似，都是从已有特征经过数学运算得来的。只是这里的逻辑关系没那么明显。所幸，虽然我们不能直观地理解多项式特征的逻辑关系，但是有一些方法和工具可以用来过滤出那些对模型准确性有帮助的特征。

首先，我们使用Pipeline来增加多项式特征，就像在前面章节介绍的那样：

from sklearn.linear_model import LogisticRegression
from sklearn.preprocessing import PolynomialFeatures
from sklearn.pipeline import Pipeline# 增加多项式预处理
def polynomial_model(degree=1, **kwarg):polynomial_features = PolynomialFeatures(degree=degree,include_bias=False)logistic_regression = LogisticRegression(**kwarg)pipeline = Pipeline([("polynomial_features", polynomial_features),("logistic_regression", logistic_regression)])return pipeline

接着，增加二阶多项式特征，创建并训练模型：

import time
model = polynomial_model(degree=2, penalty='l1', solver='liblinear')start = time.clock()
model.fit(X_train, y_train)train_score = model.score(X_train, y_train)
cv_score = model.score(X_test, y_test)
print('elaspe: {0:.6f}; train_score: {1:0.6f}; cv_score: {2:.6f}'.format(time.clock()-start, train_score, cv_score))

我们使用L1范数作为正则项（参数penalty=‘l1’），输出如下：

elaspe: 0.109572; train_score: 0.995604; cv_score: 0.991228

可以看到，训练数据集评分和测试数据集评分都增加了。为什么使用L1范数作为正则项呢？前面介绍过，L1范数作为正则项，可以实现参数的稀疏化，即自动帮助我们选择出那些对模型有关联的特征。我们可以观察一下有多少个特征没有被丢弃，即其对应的模型参数 $\theta_{j}$ 非0：

logistic_regression = model.named_steps['logistic_regression']
print('model parameters shape: {0}; count of non-zero element: {1}'.format(logistic_regression.coef_.shape, np.count_nonzero(logistic_regression.coef_)))

输出结果如下：

model parameters shape: (1, 495); count of non-zero element: 110

逻辑回归模型的coef_属性里保存的就是模型参数。从输出结果可以看到，增加二阶多项式特征后，输入特征由原来的30个增加到了495个，最终大多数特征都被丢弃，只保留了110个有效特征。

4.学习曲线

怎么知道使用L1范数作为正则项能提高算法的准确性？答案是：画出学习曲线。学习曲线是模型最有效的诊断工具之一，这也是之前一直强调的内容。

首先画出使用L1范数作为正则项所对应的一阶和二阶多项式的学习曲线：

from common.utils import plot_learning_curve
from sklearn.model_selection import ShuffleSplitcv = ShuffleSplit(n_splits=10, test_size=0.2, random_state=0)
title = 'Learning Curves (degree={0}, penalty={1})'
degrees = [1, 2]
penalty = 'l1'start = time.process_time()
plt.figure(figsize=(12, 4), dpi=144)
for i in range(len(degrees)):plt.subplot(1, len(degrees), i + 1)plot_learning_curve(plt, polynomial_model(degree=degrees[i], penalty=penalty, solver='liblinear', max_iter=300), title.format(degrees[i], penalty), X, y, ylim=(0.8, 1.01), cv=cv)print('elaspe: {0:.6f}'.format(time.process_time()-start))

输出的结果如下：

elaspe: 42.046875

L1范数学习曲线如下图所示：

接着画出使用L2范数作为正则项所对应的一阶和二阶多项式的学习曲线：

import warnings
warnings.filterwarnings("ignore")penalty = 'l2'start = time.clock()
plt.figure(figsize=(12, 4), dpi=144)
for i in range(len(degrees)):plt.subplot(1, len(degrees), i + 1)plot_learning_curve(plt, polynomial_model(degree=degrees[i], penalty=penalty, solver='lbfgs'), title.format(degrees[i], penalty), X, y, ylim=(0.8, 1.01), cv=cv)print('elaspe: {0:.6f}'.format(time.clock()-start))

输出的结果如下：

elaspe: 2.009868

L2范数学习曲线如下图所示：

可以明显地看出，使用二阶多项式并使用L1范数作为正则项的模型最优，因为它的训练样本评分最高，交叉验证样本评分也最高。从图中还可以看出，训练样本评分和交叉验证样本评分之间的间隙还比较大，我们可以采集更多的数据来训练模型，以便进一步优化模型。

L1范数对应的学习曲线，需要花比较长的时间，原因是，scikit-learn的learning_curve()函数在画学习曲线的过程中，要对模型进行多次训练，并计算交叉验证样本评分。同时，为了使曲线更平滑，针对每个点还会进行多次计算求平均值。这个就是ShuffleSplit类的作用。在我们这个实例里，只有569个训练样本，这是个很小的数据集。如果数据集增加100倍，甚至1000倍，拿出来画学习曲线将是场灾难。

问题来了，针对大数据集，怎样高效地画学习曲线？答案很简单，我们可以从大数据集里选择一小部分数据来画学习曲线，待选择好最优的模型之后，再使用全部的数据集来训练模型。有个地方需要注意，我们要尽量保持选择出来的这部分数据的标签分布与大数据集的标签分布相同，如针对二元分类，阳性和阴性比例要一致。

6.拓展阅读

1.梯度下降公式推导

关于逻辑回归模型的梯度下降公式的推导过程，可以参阅博客http://blog.kamidox.com/logistic-regression.html。

2.向量形式

实际上，我们的预测函数就是写成向量形式的：
$h_{\theta}(x)=g(z)=g(\theta^Tx)=\frac{1}{1+e^{-\theta^Tx}}$
这个预测函数一次只计算一个训练样本的预测值，我们要怎样一次性计算出所有样本的预测值呢？
$h=g(X\theta)$
上述公式即可达到目的。其中g(x)为Sigmoid函数。X为m×n的矩阵，即数据集的矩阵表达。
成本函数也有对应的矩阵形式：
$J(\theta)=\frac{1}{m}(-y^Tlog(h)-(1-y)^Tlog(1-h))$
其中，y为目标值向量，h为一次性计算出来的所有样本的预测值。