逻辑回归算法的名字虽然带有“回归”二字，但实际上是用来解决分类问题的算法。

1.逻辑回归算法的原理

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

与线性回归算法的最大区别是，逻辑回归算法的输出是个离散值。

1.预测函数

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

选择Sigmoid函数（也称为Logistic函数，逻辑回归的名字由此而来）

来作为预测函数，其中e是自然对数的底数。以z为横坐标，以g(z)为纵坐标，画出的图形如下所示：

figure6_1.png

问题来了，怎样将输入特征和预测函数结合起来呢？
结合线性回归算法的预测函数，令，则逻辑回归算法的预测函数如下：

表示在输入值为x，参数为的前提条件下y=1的概率。用概率论的公式可以写成：

在输入x及参数的条件下y=1的概率，这是个条件概率公式，

对二元分类来说，这是个非黑即白的世界。

2.判定边界

逻辑回归算法的预测函数由下面两个公式给出：


假定y=1的判定条件是，y=0的判定条件是，则可以推导出y=1的判定条件是，y=0的判定条件就是。所以，即是我们的判定边界。

下面给出两个判定边界的例子。
假设有两个变量x1，x2，其逻辑回归预测函数是。假设给定参数：

那么，可以得到判定边界，即，如果以为横坐标，为纵坐标，则这个函数画出来就是一条通过(0,3)和(3,0)两点的直线。这条线就是判定边界，如图所示：

figure6_2.png

如果预测函数是多项式，且给定

则可以得到判定边界函数。还是以为横坐标，以为纵坐标，则这是一个半径为1的圆。圆内部是y=0，圆外部是y=1，如上图所示。
这是二阶多项式的情况，更一般的多阶多项式可以表达出更复杂的判定边界。

3.成本函数

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

为了容易地求出成本函数的最小值，我们分成y=1和y=0两种情况来分别考虑其预测值和真实值的误差。我们先考虑最简单的情况，即计算某个样本x，y=1和y=0两种情况下的预测值与真实值的误差，我们选择的成本公式如下：

其中，表示预测为1的概率，log(x)为自然对数。如图所示：

figure6_3.png

思考：符合上述规律的函数有很多，逻辑回归算法为什么要选择自然对数函数来作为成本函数呢？
因为逻辑回归模型的预测函数是Sigmoid函数，里面有e的n次方运算，自然对数刚好是其逆运算：。选择自然对数，可以推导出形式简单的逻辑回归模型参数的迭代公式，而不需要去涉及对数运算和指数运算。这就是我们选择自然对数函数来作为逻辑回归的成本函数的原因。

更进一步，把输入值x从映射到区间的函数有很多，逻辑回归算法为什么要选择Sigmoid函数作为预测函数呢？
严格地讲，不一定非要选择Sigmoid函数作为预测函数。但如果不选择Sigmoid函数，就需要重新选择性质接近的成本函数，使用现成的Sigmoid函数既方便表达、效率又高。

分开表述的成本计算公式始终不方便，能不能合并成一个公式呢？考虑下面的公式：

由于y只能取0或1，所以上述式子和分开表达的公式是等价的。
根据一个样本的成本计算公式，很容易写出所有样本的成本平均值，即成本函数：

4.梯度下降算法

和线性回归类似，这里使用梯度下降算法来求解逻辑回归模型参数。根据梯度下降算法的定义，可以得出：

推导出来的梯度下降算法参数迭代公式为：

如果对公式推导过程感兴趣，可以参考后面扩展阅读的内容。
这个公式的形式和线性回归算法的参数迭代公式是一样的。当然，由于这里，而线性回归算法里的。所以，两者的形式一样，而数值计算方法则完全不同。

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个预测函数：

预测出来的概率最高的那个类别就是样本所属的类别。

3.正则化

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

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

1.线性回归模型正则化

我们先来看线性回归模型的成本函数是如何正则化的：

公式中前半部分就是原来的线性回归模型的成本函数，也称为预测值与实际值的误差。后半部分为加入的正则项。其中的值有两个目的，即要维持对训练样本的拟合，又要避免对训练样本的过拟合。如果的值太大，则能确保不出现过拟合，但可能会导致对现有训练样本出现欠拟合。

思考：怎样从数学上理解正则化后的成本函数解决了过拟合问题呢？
从数学角度看，成本函数增加了一个正则项后，成本函数不再唯一地由预测值和真实值的误差所决定，还和参数的大小有关。有了这个限制之后，要实现成本函数最小的目的，就不能随便取值了。比如某个比较大的值可能会让预测值与真实值的误差值很小，但会导致很大，最终的结果是成本函数很大。这样通过调节参数，就可以控制正则项的权重，从而避免线性回归算法过拟合。
利用正则化的成本函数，可以推导出正则化后的参数迭代公式：


因子在每次迭代时，都将把收缩一点点。因为和是正数，而m是训练样本的个数，是个比较大的正整数。为什么要对进行收缩呢？因为加入正则项的成本函数和成正比，所以迭代时需要不断地试图减小的值。

2.逻辑回归模型正则化

同样，可以对逻辑回归模型的成本函数进行正则化，其方法也是在原来的成本函数的基础上加上正则项：

相应地，正则化后的参数迭代公式为：



上式中，因为没有参与正则化。另外，逻辑回归和线性回归看起来形式是一样的，但其实它们的算法不一样，因为两个式子的预测函数不一样：线性回归，逻辑回归。

4.算法参数

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

1.正则项权重

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

2.L1/L2范数

创建逻辑回归模型时，有个参数penalty（惩罚），其取值有“l1”或“l2”，即指定正则项的形式。在成本函数里添加的正则项，这个实际上是个L2正则项，即把L2范数作为正则项。当然，也可以添加L1范数来作为正则项。

L1范数作为正则项，会让模型参数稀疏化，即让模型参数向量里的0元素尽可能多，只保留模型参数向量中重要特征的贡献。而L2范数作为正则项，则让模型参数尽量小，但不会为0，即尽量让每个特征对应预测值都有一些小的贡献。

问题来了，为什么会造成上述不同的结果呢？这就需要从L1范数和L2范数的定义开始讨论了。
为了简单起见，假设模型只有两个参数，它们构成一个二维向量，则L1范数为：

即L1范数是向量里元素的绝对值之和。L2范数为向量里所有元素的平方和的算术平方根：

我们知道，梯度下降算法在参数迭代的过程中，实际上是在成本函数的等高线上跳跃，并最终收敛在误差最小的点上。那么正则项的本质是什么？正则项的本质是惩罚。在参数迭代的过程中，如果没有遵循正则项所表达的规则，那么其成本会变大，即受到了惩罚，从而往正则项所表达的规则处收敛。正则化后的模型参数应该收敛在误差等高线与正则项等高线相切的点上。

我们把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 = v

def 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--')

figure6_4.png

在右图中，使用L2范数来作为正则项，L2范数表示的是元素平方和的算术平方根，这里L2范数的值为1，在，坐标轴上，它的等值线是一个圆。它和模型误差等值线相切的点，一般不在坐标轴上（切点坐标没有0元素）。

因此L1范数作为正则项，会让模型参数稀疏化；而L2范数作为正则项，则会使模型的特征对预测值都有少量的贡献，避免模型过拟合。

作为推论，L1范数作为正则项，有以下几个用途：

• 选择重要特征：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-01
 1.471e-01 2.419e-01 7.871e-02 1.095e+00 9.053e-01 8.589e+00 1.534e+02
 6.399e-03 4.904e-02 5.373e-02 1.587e-02 3.003e-02 6.193e-03 2.538e+01
 1.733e+01 1.846e+02 2.019e+03 1.622e-01 6.656e-01 7.119e-01 2.654e-01
 4.601e-01 1.189e-01]

数据集中总共有569个样本，每个样本有30个特征，其中357个阳性（y=1）样本，212个阴性（y=0）样本。同时，还打印出一个样本数据，以便直观地进行观察。

这30个特征是怎么来的呢？这个数据集总共从病灶造影图片中提取了以下10个关键属性：

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

实际上它只关注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.9的样本（同时也是预测为阳性的概率大于0.1小于0.9），这些样本就是“自信度”不足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()函数，而计算预测分类用model.predict()函数。

3.模型优化

首先，使用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范数作为正则项可以实现参数的稀疏化，即自动选择出那些对模型有关联的重要特征。

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范数作为正则项所对应的一阶和二阶多项式的学习曲线：

from common.utils import plot_learning_curve
from sklearn.model_selection import ShuffleSplit

cv = 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范数学习曲线如下图所示：

figure6_5.png

接着画出使用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范数学习曲线如下图所示：

figure6_6.png

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

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

那么，针对大数据集，怎样高效地画学习曲线？答案很简单，可以从大数据集里选择一小部分数据来画学习曲线，待选择好最优的模型之后，再使用全部的数据集来训练模型。但是要尽量保持选择出来的这部分数据的标签分布与大数据集的标签分布相同，如针对二元分类，阳性和阴性比例要一致。更直观的说就是，抽取出来的样本集为原来数据集的一个缩影，尽可能相似。

6.拓展阅读

1.梯度下降公式推导

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

2.向量形式

实际上，我们的预测函数就是写成向量形式的：

这个预测函数一次只计算一个训练样本的预测值，怎样一次性计算出所有样本的预测值呢？答案是把预测函数的参数写成向量的形式：

其中g(x)为Sigmoid函数。X为m×n的矩阵，即数据集的矩阵表达。
成本函数也有对应的矩阵形式：

其中，y为目标值向量，h为一次性计算出来的所有样本的预测值。

3.算法性能优化

梯度下降算法的效率比较低，优化的梯度下降算法有Conjugate Gradient、BFGS、L-BFGS等。这些算法比较复杂，实现这些算法是数值计算专家的工作，一般工程人员只需要了解这些算法且会使用即可。