岁虚山

行有不得，反求诸己

0%

Logistic 回归

发表于 更新于 分类于

前言

假设现在平面上有一些数据点，我们需要使用一条直线对这些点进行拟合，也就是寻找最佳拟合直线，这个拟合过程就称作回归。

而利用 Logistic 回归进行分类的主要思想就是：根据现有数据对分类边界线建立回归公式，以此进行分类。

“回归”一词实际上就是源于最佳拟合。而我们训练的过程就是找到最佳拟合参数，从而得到我们的回归公式。

接下来，我们将对 Logistic 回归分类器的推导和实现过程进行展示。

原理

所谓拟合，其实就是根据现有的样本，寻找某种统一的规律，使得绝大部分样本的分布都尽可能地满足这种这种规律。

如下图：

01-拟合示意图

上图中，样本大致分布在 这条直线附近，那么我们就称方程 为对样本的坐标分布的拟合。

其中，各变量的系数 $(2,-1,-1)$ 就叫做拟合参数。

当然，上图中的例子，我们所做的拟合严格来讲应该属于线性回归（Linear Regression）的范围。这是因为上例中的因变量 的值并不是一个二类值（即非 0 即 1）。

Logistic 回归是基于线性回归的基础上，通过引入 Sigmoid 函数，将线性问题转换为非线性问题，进而达到分类的效果。

线性回归

在了解 Logistic 回归之前，我们首先来了解一下线性回归。

所谓线性，就是一个函数，或者称为映射，且要求同时满足两个条件：可加性和齐次性。

这两个概念解释起来比较麻烦，在这里我们可以将其简单的理解为对于多元变量函数：

即函数 的值是其各个自变量通过加减运算得到的，各自变量之间的幂数相同。

而对于类似：

的函数就不属于线性函数，因为函数中的各项自变量的幂次不同。

建立模型

对于数据集 $D$，共有 个样本 ，其中， 是样本 的标签，每个样本 有 个维度：

线性模型的目的是为了找到一个线性组合使得对于 中的所有样本：

用向量形式表示为：

其中：

我们将每个样本的标签 也写成矩阵形式：

通常来讲，带有 的形式不太好求解，因此可以根据线性表达式：

将矩阵 和 进行改造，将 放进 矩阵中，将 矩阵增加一列，即令：

则矩阵表达可简化为：

我们的目标就是找到合适的 ，使得 尽量趋向于 ，为了表示这两者之间的差异，通常可以用以下几种方式进行衡量：

1、误差平方和：

2、欧氏距离

3、曼哈顿距离

当然，还有其他的距离表示方式，在这里不一一例举。

对于样本 ，若其特征是一维的，则称其回归问题为一元线性回归。一元线性回归的问题的本质就是在 平面上求的一条直线，来尽量拟合在该平面内的所有样本点。其中 就是直线的斜率， 就是直线在 轴上的截距。

若 的特征是二维的，即二元线性回归，则其本质是在 空间中求的一个平面，来尽量拟合在该三维空间中的所有样本点。

对于多元线性回归问题，则求的就是一个尽量拟合所有样本的超平面。

求解

我们以误差平方和作为度量（也可以称之为损失函数）：

则目标就是将损失降到最小：

矩阵形式表示为：

对于求解，通常有两种方式：

最小二乘法

由于目标函数是关于 的凸函数，因此当其关于 和 的导数为 0 时，得到最优解，即：

上式中，系数 2 并不影响求导结果，为了求导方便，一般将误差平方和写作以下形式：

这样，求导后就不会出现 2 的系数了。

因此：

用矩阵的形式来表示就是：

通过极值点求最优解：

通常情况下 往往不是满秩矩阵，所以无法求解矩阵的逆，故无法求的唯一的阶。因此通常会引入正则化的方式，将矩阵补成满秩。

L2正则：

求得：

通过正则化，还可以限制模型复杂度。

当然，除了 L2正则，还有岭回归、lasso回归、弹性网络等正则化方式。

梯度下降法

虽然使用最小二乘法可以一步就能求得最优解，但是其中会涉及到大量的矩阵的运算，因此在一般情况下会使用梯度下降法来进行求解。

我们仍然使用误差平方和来作为损失函数：

损失函数 在任意一点的梯度为：

梯度下降法通过随机化初始的 并通过迭代的方式来逐渐逼近最优解。

在迭代的过程中， 向负梯度的方向增长。过程如下：

将 初始化为随机数，作为迭代开始的初始值：

在第 次迭代时，通过第 次迭代时 的值计算梯度：

计算完梯度后，需要判断梯度是否接近于 ，若接近于 则认为已经逼近了最优解，可以跳出迭代，否则继续迭代。

通过学习率，来计算 移动的距离，方向为负梯度方向：

其中， 为学习率，也可以理解为每次移动的步长。在迭代刚开始的时候，学习率可以设置大一点，前期移动的距离可以大一点；迭代后期可以设置小一点，这样移动距离小，可以更精确。

重复迭代过程，直到达到最大迭代次数或梯度趋近于 。

伪代码如下：

为了简化计算过程，仍然可以采用最小二乘法中的思想，在样本 中添加一个值为 1 的特征，将 维特征变为 维，即：

然后同样在 添加一个元素 ，将其视作 ，即：

这样，在迭代过程中就只用计算 ，并对 进行迭代即可。

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
import numpy as np

def create_samples(size: int, dimension: int):
    x = np.mat(np.random.rand(size, dimension))
    alphas = np.mat(np.random.rand(dimension, 1)) * 10
    noise = np.mat(np.random.rand(size, 1))
    y = x * alphas + noise
    return x, y, alphas

def linear_regression(data_mat: np.matrix, value: np.matrix):
    m, n = data_mat.shape
    data_mat = np.insert(data_mat, 0, np.ones(m), axis=1)
    weights = np.ones((n + 1, 1))

    max_cycles = 10000

    step = 0.001

    for i in range(max_cycles):
        partial_w = data_mat.transpose() * (data_mat * weights - value)
        if np.fabs(np.max(partial_w)) < 1e-7:
            break
        weights = weights - step * partial_w

    return weights[0, 0], weights[1:, 0]

def predict(x: np.matrix, w: np.matrix, b: float):
    return x * w + b

x, y, alpha = create_samples(400, 4)
b, weight = linear_regression(x, y)
print(alpha)
print(weight)

create_samples 函数：

代码的第 3-8 行，定义了一个 create_samples 函数，用来创建样本数据。该函数接收两个 int 型的参数 size 和 dimension，分别代表样本的数量和样本的特征的维度。

代码的第 4 行，创建了一个 size 行 dimension 列 的随机数矩阵，作为样本集。

代码的第 5 行，创建了一个 dimension 行 1 列 的随机数矩阵作为样本特征的权重。

代码的第 6 行，创建了一个 size 行 1 列 的随机数矩阵作为噪音。

代码的第 7 行，创建通过样本 x 和权重 alpha 来计算出样本的标签值，并加上噪音 noise 数据。

linear_regression 函数：

代码的第 10-25 行，定义了一个 linear_regression 函数，用来使用线性回归来计算回归系数。

该函数接收两个 np.matrix 类型的参数 data_mat 和 value，分别代表样本和其标签值。

代码的第 11 行，计算出样本矩阵的行和列，并分别保存在变量 m 和 n 中。

代码的第 12 行，通过对输入的样本数据矩阵进行扩展，添加了一列全 1 的列向量，形成了一个新的 m 行 n+1 列 的矩阵。

代码的第 13 行，创建了一个全 1 的 n+1 列向量，作为初始的回归系数。

代码的第 15 行，定义了最大的迭代次数。

代码的第 17 行，定义了每次移动的步长，也可以称为学习率。

代码的第 19-23 行，开始进行迭代计算。

代码的第 20 行，计算当前的梯度值（偏导）。

代码的第 21-22 行，通过梯度值来判断是否满足停止迭代的条件（偏导接近于 0）。

代码的第 23 行，通过梯度和学习率来更新回归系数。

代码的第 25 行，返回回归系数向量的第一个元素，作为 ，返回剩下的元素作为回归系数 。

predict 函数：

代码的第 27-28 行，定义了一个 predict 函数，用来对输入的样本的标签值进行预测。

该函数接收三个参数，np.matrix 型的参数 x 为待预测的样本，是一个 1 行 n 列 的行向量；np.matrix 型的参数 w 为和 float 型的参数 b 为回归系数。

代码的第 28 行，计算 作为预测值进行返回。

测试：

代码的第 30-33 行，为测试语句。

代码的第 30 行，通过 create_samples 函数创建了 400 个 4 维特征的样本集。

代码的第 31 行，通过 linear_regression 函数计算回归系数 和 。

代码的第 32、33 行，即打印创建样本时所使用的权重和线性通过线性回归计算出的回归系数。

下图为迭代过程中的误差平方和变化：

误差

可以看出，在迭代刚开始时，误差的下降速度非常快，最终在第 7600 左右次迭代时，偏导已经趋向于 0 。

为了更好地展示线性回归的效果，我们使用 1 维特征的样本集来对回归效果进行展示：

线性回归测试

Logistic 回归

Logistic 回归于线性回归原理一致，但与线性回归的用途不同。

一般来讲，线性回归用于数值型的拟合问题，即 ，而 Logistic 运用于分类问题，即 。

如果说，能够想办法将数值问题转化为分类问题，那么就能够使用线性回归的思想来解决分类问题。

我们假想这样一种函数 ，能够接收所有的输入值，而其输出只有两个值，如下：

其函数图像如下：

越阶函数图像

那么只需要将线性回归 的结果通过 映射到 上，即：，即可将线性回归转化为 Logistic 回归。

但是上面的越阶函数存在一定的问题，其在 点处是瞬间跳跃的，这个跳跃过程有时不好处理。因此我们采用另一个具有类似性质的连续函数 Sigmoid 函数来将其替代：

其函数图像如下所示：

Sigmoid

可以看出，该函数是处处光滑的，在 点的值为 ，当我们将 轴的范围扩大时：

Sigmoid2

或者在 前加上系数时：

Sigmoid3

建立模型

在线性回归中，我们使用： 来表示样本 的预测值，现在，我们要将线性回归应用到 Logistic 回归上，因此可以借助 Sigmoid 函数，将 映射为 ，即：

其中，。

为了方便，我们将 记作 ，此时， 不再代表对样本的回归预测的值，而是代表对样本的预测类别。

我们将式子进行转化：

其中，我们可以认为 为类别为 1 的可能性，即 ； 为类别为 0 的可能性，即 。

因此二者的比值表示一种几率，即样本 为 1 类别的相对可能性。取对数表示“对数几率”。

Logistic 回归就是使用线性回归的方式，去预测这个对数几率，从而根据对数几率来实现分类的效果。因此， Logistic 回归的实质含义应该是：对数几率回归。

即：

也可以理解为将 作为样本，将对数几率作为样本值，使用线性回归去训练回归系数。

那么可以得到：

而 ，则可得：

将上面两式合并可得：

因为 只能取 0 或 1，因此，上式的含义为：

对于样本集而言，其中共有 个样本，样本分别为 ，各个样本对应的样本值（标签）为 。

设 为 的概率，即：，则 的概率为 ，则：

又因为各样本之间相互独立，那么他们的联合分布为各边缘分布的乘积，那么可以得到似然函数：

我们的目标就是求出使这一似然函数的值最大的参数估计（）。

对函数取对数可得：

则目标参数为：

将其转化为最小化负的对数似然函数：

其中：

如此，就得到了 Logistic 回归的损失函数，即机器学习中的二元交叉墒（Binary crossentropy）：

使用梯度下降法进行优化，首先对 求导：

为了简化运算，我们将样本和回归系数进行了改造，去掉了 ，简化为：

因此得到：

伪代码

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
from sklearn.datasets import make_blobs
import numpy as np

def create_samples(size: int, dimension: int):
    samples, labels = make_blobs(n_samples = size, n_features = dimension, centers = 2, cluster_std = 3.5)
    return np.mat(samples), np.mat(labels).transpose()

def sigmoid(x: np.matrix):
    return 1.0 / (1 + np.exp(-x))

def logistic_regression(data_mat: np.matrix, value: np.matrix):
    m, n = data_mat.shape
    data_mat = np.insert(data_mat, 0, np.ones(m), axis=1)
    weights = np.ones((n + 1, 1))

    max_cycles = 500

    step = 0.01

    for i in range(max_cycles):
        partial_w = data_mat.transpose() * (sigmoid(data_mat * weights) - value)
        weights = weights - step * partial_w

    return weights[0, 0], weights[1:, 0]

def predict(x: np.matrix, w: np.matrix, b: float):
    return 0 if x * w.transpose() + b < 0.5 else 1


x, y = create_samples(100, 2)
b, weight = logistic_regression(x, y)

create_samples 函数

代码的第 4-6 行，定义了一个 create_samples 函数，用来创建样本数据。该函数接收两个 int 型的参数 size 和 dimension，分别代表样本的数量和样本的特征的维度。

代码的第 4 行，调用 make_blobs 函数创建样本和标签。

sigmoid 函数：

代码的第 8-9 行，定义了越阶函数 。

logistic_regression 函数：

代码的第 11-24 行，定义了一个 logistic_regression 函数，用来使用线性回归来计算回归系数。

该函数接收两个 np.matrix 类型的参数 data_mat 和 value，分别代表样本和其标签值。

代码的第 12 行，计算出样本矩阵的行和列，并分别保存在变量 m 和 n 中。

代码的第 13 行，通过对输入的样本数据矩阵进行扩展，添加了一列全 1 的列向量，形成了一个新的 m 行 n+1 列 的矩阵。

代码的第 14 行，创建了一个全 1 的 n+1 列向量，作为初始的回归系数。

代码的第 16 行，定义了最大的迭代次数。

代码的第 18 行，定义了每次移动的步长，也可以称为学习率。

代码的第 20-24 行，开始进行迭代计算。

代码的第 21 行，计算当前的梯度值（偏导）。

代码的第 22 行，通过梯度和学习率来更新回归系数。

代码的第 24 行，返回回归系数向量的第一个元素，作为 ，返回剩下的元素作为回归系数 。

predict 函数：

代码的第 26-27 行，定义了一个 predict 函数，用来对输入的样本的标签值进行预测。

该函数接收三个参数，np.matrix 型的参数 x 为待预测的样本，是一个 1 行 n 列 的行向量；np.matrix 型的参数 w 为和 float 型的参数 b 为回归系数。

代码的第 27 行，计算 与 0.5 进行比较，返回类别。

下图为迭代过程中的误差变化：

Logistic 回归误差

分类效果：

Logistic 分类效果

优化

目前，我们已经成功使用 Logistic 回归实现了分类模型。但是在计算过程中，使用的梯度下降法在每次更新回归系数 时都需要遍历整个数据集。当数据集中的样本数量较少时，计算速度尚在可接受的范围内。但是当数据集的规模达到一定的量级或样本的特征较多时，那么使用这种方式进行回归系数的更新的计算复杂度就太高了。

为了减少因样本数量太多或样本特征数量太大所造成的计算开销，可以改进回归系数的计算方式，将迭代更新回归系数的方式改为遍历样本集，每次选择一个样本来更新回归系数，我们暂且称之为随机梯度下降法。

采用这种形式，可以在新的样本加入时，对分类器进行增量式更新，达到在线学习的效果。

伪代码如下：

只需要将 logistic_regression 函数修改为：

1
2
3
4
5
6
7
8
9
10
11
12
def logistic_regression(data_mat: np.matrix, value: np.matrix):
    m, n = data_mat.shape
    data_mat = np.insert(data_mat, 0, np.ones(m), axis=1)
    weights = np.ones((n + 1, 1))

    step = 0.01

    for i in range(m):
        partial_w = data_mat[i].transpose() * (sigmoid(sum(data_mat[i] * weights)) - value[i])[0, 0]
        weights = weights - step * partial_w

    return weights[0, 0], weights[1:, 0]

分类效果如下：

随机Logistic 分类效果

我们将原版的 Logistic 回归的效果与采用随机梯度下降进行比较：

Logistic 随机梯度下降效果对比

可以看出，随机梯度下降的方式虽然降低了运算量，但是分裂效果却并不好。

因此，我们可以对算法进行改造，增加迭代次数：

1
2
3
4
5
6
7
8
9
10
11
12
13
def logistic_regression(data_mat: np.matrix, value: np.matrix):
    m, n = data_mat.shape
    data_mat = np.insert(data_mat, 0, np.ones(m), axis=1)
    weights = np.ones((n + 1, 1))

    step = 0.01
    max_cycle = 200
    for it in range(max_cycle):
        for i in range(m):
            partial_w = data_mat[i].transpose() * (sigmoid(sum(data_mat[i] * weights)) - value[i])[0, 0]
            weights = weights - step * partial_w
            
    return weights[0, 0], weights[1:, 0]

效果对比如下：

我们看一下使用随机梯度下降时，回归系数的变化：

Logistic 随机梯度下降系数变化

可以看到，回归系数 $w_0$ 很快便收敛，而 $w_1、w_2$ 则在迭代程中则不断波动。衡量一个算法的优劣的方法，是看其是否收敛，很明显，$w_1、w_2$ 未能很快的达到收敛的效果。因此，我们还需要对算法进行改进，使其尽量避免产生波动，从而收敛于某个值，同时，也想使其收敛速度加快。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
import random

def logistic_regression(data_mat: np.matrix, value: np.matrix):
    m, n = data_mat.shape
    data_mat = np.insert(data_mat, 0, np.ones(m), axis=1)
    weights = np.ones((n + 1, 1))

    step = 0.01
    max_cycle = 100
    for it in range(max_cycle):
        indexs = [index for index in range(m)]
        random.shuffle(indexs)
        #indexs = random.sample([i for i in range(0, m)], m)
        for i in range(m):
            step = 4 / (1.0 + it + i) + 0.01
            rand_index = indexs[i]
            partial_w = data_mat[rand_index].transpose() * (sigmoid(sum(data_mat[rand_index] * weights)) - value[rand_index])[0, 0]
            weights = weights - step * partial_w
            
    return weights[0, 0], weights[1:, 0]

新的函数中，在第 11 行，根据样本的数量 $m$ 生成了下标列表，并在代码的第 12 行将列表打乱，得到了一个打乱后的下标列表。这两行代码的效果与第 13 行相同，也可以使用第 13 行代码来生成这样的一个乱序的下标列表。

在 14-18 行的循环中，从乱序的下标列表中选取一个值，作为更新回归系数的样本下标（相当于遍历样本集时，将样本集随机打乱后再进行遍历），以此来消除每次迭代时因为样本顺序而产生的波动。

代码的第15 行，对学习步长（学习率）进行更新，使其在迭代开始时的值较大，而随着迭代次数的增加，步长越来越小，使得目标能够快速收敛。

程序分类效果如下：

Logistic 随机梯度下降效果对比2

系数变化如下：

Logistic 随机梯度下降系数变化2

可以看出，随机选择样本样本避免了波动，而随着迭代变化的学习步长，明显能够加快算法的收敛速度。

-------------本文结束 感谢阅读-------------
打赏一包辣条