《机器学习实战》分类篇04.Logistic回归

最新推荐文章于 2023-12-04 17:07:03 发布

ZPILOTE

最新推荐文章于 2023-12-04 17:07:03 发布

阅读量499

点赞数

分类专栏：机器学习文章标签： Logistic回归 Sigmoid函数随机梯度上升法疝病马死亡率处理数据中的缺失值

本文链接：https://blog.youkuaiyun.com/MRZHUGH/article/details/102991100

版权

机器学习专栏收录该内容

4 篇文章

订阅专栏

本文深入讲解Logistic回归的原理及应用，包括Sigmoid函数、梯度上升法、随机梯度上升法，以及如何处理缺失数据。通过病马生死预测案例，展示了Logistic回归在实际问题中的实施步骤。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

文章目录

Logistic回归

假设现在有一些数据点，用一条直线对这些点进行拟合（该线称为最佳拟合直线），这个拟合过程就称为回归。利用Logistic回归进行分类的主要思想是：根据现有数据对分类边界线建立回归公式，以此进行分类。

Logistic回归的一般过程 1. 收集数据：采用任意方法收集数据 2. 准备数据：由于需要进行距离计算，因此要求数据类型为数值型。另外，结构化数据格式则最佳 3. 分析数据：采用任意方法对数据进行分析 4. 训练算法：大部分时间用于训练，训练的目的是为了找到最佳的分类回归系数 5. 测试算法：一旦训练步骤完成，分类将会很快 6. 使用算法：首先，我们需要输入一些数据，并将其转换成对应的结构化数值；接着，基于训练好的回归系数就可以对这些数值进行简单的回归计算，判定它们属于哪个类别；在这之后，我们就可以在输出的类别上做一些其它分析工作。

接下来学习最优化算法（梯度上升法和一个改进的随机梯度上升法），并利用它们训练出一个非线性函数用于分类。

1.基于Logistic回归和Sigmoid函数的分类

Logistic回归优点：计算代价不高，易于理解和实现

缺点：容易欠拟合，分类精度可能不高

适用数据类型：数值型和标称型数据

理想的函数应该是：能接受所有的输入然后预测出类别。例如：在两个类的情况下，上述函数输出0或1。数学上有个函数Sigmoid函数。具体的计算公式如下：
$\sigma(z)=\frac{1}{1+e^z}$
下图给出了Sigmoid函数在不同坐标尺度下的两条曲线图。当x为0时，Sigmoid函数值为0.5。随着x的增大，对应的Sigmoid值将逼近于1；而随着x的减小，Sigmoid值将逼近于0。如果横坐标刻度足够大，Sigmoid函数看起来很像一个阶跃函数。

因此，为了实现Logistic回归分类器，可以在每个特征上都乘以一个回归系数，然后把所有的结果值相加，将这个总和代入Sigmoid函数中，进而得到一个范围0~1之间的数值。任何大于0.5的数据被分入1类，小于0.5的被分入0类。

那么问题来了：最佳回归系数是多少？

2.基于最优化方法的最佳回归系数确定

Sigmoid函数的输入记为z，由下面公式得出：
$z = w_0x_0+w_1x_1+w_2x_2+...+w_nx_n$
如果采用向量的写法，上述公式可写出：
$z=w^Tx$
它表示将两个数值向量对应元素相乘然后全部加起来即得到z值。其中的向量x是分类器的输入数据，向量w就是要找到的最佳参数（系数），从而使得分类器尽可能地精确。为了寻找该最佳参数，需要用到最优化理论。

2.1梯度上升法

思想是：要找到某函数的最大值，最好的方法是沿着该函数的梯度方向探寻。

数学原理详解： https://www.cnblogs.com/HongjianChen/p/8718988.html

可以类比吴恩达大佬的机器学习课程中梯度下降法，与这里的梯度上升法一样的，只是加法需要变成减法。

而且Logistic回归的目标函数：代表的是概率，求概率最大值，用梯度上升法；

但线性回归的代价函数：代表的则是误差，求误差最小值，用梯度下降法。

2.2训练算法：适用梯度上升找到最佳参数

给定数据样本点，两个数值型特性：X1和X2。通过使用梯度上升法找到最佳回归系数。

梯度上升法的伪代码：

每个回归系数初始化为1
重复R次：
	计算整个数据集的梯度
	使用alpha x gradient更新回归系数的向量
	返回回归系数

在Pycharm中新建Logistic.py

import numpy as np

# 打开文本并逐行读取
def loadDataSet():
    dataMat = []
    labelMat = []
    fr = open('testSet.txt')
    for line in fr.readlines():
        lineArr = line.strip().split()
        # 每行前两个值分别是X1和X2，
        dataMat.append([1.0,float(lineArr[0]),float(lineArr[1])])
        # 第三个值是数据对应的类别标签
        labelMat.append(int(lineArr[2]))
    return dataMat, labelMat

# sigmoid函数
def sigmoid(inX):
    return 1.0/(1+ np.exp(-inX))

# 转换为Numpy矩阵数据类型
def gradAscent(dataMatIn, classLabels):
    dataMatrix = np.mat(dataMatIn)
    labelMat = np.mat(classLabels).transpose()
    m,n = np.shape(dataMatrix)
    alpha = 0.01 # 向目标移动的步长
    maxCycles = 500 # 迭代次数
    weights = np.ones((n,1))  # 回归系数
    for k in range(maxCycles):
        # 矩阵相乘
        h = sigmoid(dataMatrix * weights)
        error = (labelMat - h)
        weights = weights + alpha * dataMatrix.transpose() * error
    return weights

然后把testSet.txt文本文件放到当前项目目录下；

继续在Pycharm中新建main.py：

import Logistic

if __name__ == '__main__':
    dataArr,labelMat = Logistic.loadDataSet()
    print(Logistic.gradAscent(dataArr,labelMat))

输出结果：

在这里插入图片描述

2.3分析数据：画出决策边界

系数已经求出，为了便于可视化理解，在Logistic.py中添加代码：

# 画出数据集和Logistic回归最佳拟合直线的函数
def plotBestFit(weights):
    weights = weights.getA() # getA()函数将矩阵类型转化为数组，与mat函数正好相反
    dataMat,labelMat = loadDataSet()
    dataArr = np.array(dataMat)
    n = np.shape(dataArr)[0]
    xcord1 = []; ycord1 = []
    xcord2 = []; ycord2 = []
    for i in range(n):
        if int(labelMat[i]) == 1:
            xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])
        else:
            xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')
    ax.scatter(xcord2, ycord2, s=30, c='green')
    x = np.arange(-3.0,3.0,0.1)
    # 最佳直线拟合
    y = (-weights[0]-weights[1]*x)/weights[2]
    ax.plot(x,y)
    plt.xlabel('X1'); plt.ylabel('X2')
    plt.show()

然后再main.py中添加：

    weights = Logistic.gradAscent(dataArr,labelMat)
    print(weights)
    Logistic.plotBestFit(weights)

运行，输出结果：

在这里插入图片描述

从结果图来看，只分错了5个点，但是，这个方法却需要了300次乘法运算，因此，还需要对优化算法进行改进。

2.4训练算法：随机梯度上升

梯度上升算法在每次更新回归系数时都需要遍历整个数据集，如果有数十亿样本和成千上万的特征，该方法的计算复杂度就太高了。

一种改进方法是一次仅用一个样本点来更新回归系数，该方法称为随机梯度上升算法。

由于可以在新样本到来时对分类器进行增量式更新，因而随机梯度上升算法是一个在线学习算法。与“在线学习”相对应，一次处理所有数据被称为“批处理”。

随机梯度的伪代码和代码：

所有回归系数初始化为1
对数据集中每个样本
		计算该样本的梯度
		使用alhpa x gradient更新回归系数值
返回回归系数

def stocGradAscent0(dataMatrix, classLables):
    m,n = np.shape(dataMatrix)
    alpha = 0.01
    weights = np.ones(n)
    for i in range(m):
        h = sigmoid(sum(dataMatrix[i] * weights))
        error = classLables[i] - h
        weights = weights + alpha * error * dataMatrix[i]
    return weights

stocGradAscent0()函数gradAscent()函数很相似，但也有区别：第一，后者的变量h和误差error都是向量，前者则全是数值；第二，前者没有矩阵的转换过程，所有的变量的数据类型都是NumPy数组。

然后在main()添加：

import Logistic
import numpy as np # 用到array,导入numpy

if __name__ == '__main__':
    dataArr,labelMat = Logistic.loadDataSet()
    # weights = Logistic.gradAscent(dataArr,labelMat)
    # print(weights)
    # Logistic.plotBestFit(weights)
    weights = Logistic.stocGradAscent0(np.array(dataArr),labelMat)
    print(weights)
    Logistic.plotBestFit(weights)

由于之前把weights = weights.getA() 写进plotBestFit()函数，而这里dataArr已经转换为数组了，故不需要weights = weights.getA()这一行了，所以先在plotBestFit()函数中注释掉这一行，再运行函数。（也可以把这行挪到main.py中，然后注释掉）

输出结果：

L5.jpy

可以看出拟合的直线效果还不错，但这里分类器错分了三分之一左右的样本。

一个判断优化算法优劣的可靠方法是看它是否收敛，也就是说参数是否达到了稳定值，是否还会不断地变化？这里主要看系数X0,X1和X2，三个系数的收敛速度各不相同，主要是因为存在一些不能正确分类的样本点（数据集并非线性可分）。而步长alpha也会影响该收敛速度。所以也要把alpha设为动态可变的。

对此，对随机梯度上升算法做了修改，使其在整个数据集上运行150次。

改进的随机梯度上升算法的代码：

# 比stocGradAscent0()函数多了迭代次数numIter
def stocGradAscent1(dataMatrix, classLables, numIter = 150):
    m,n = np.shape(dataMatrix)
    weights = np.ones(n)
    for j in range(numIter):
        # python3.x中range返回的是range对象，不返回数组对象，所以用list
        dataIndex = list(range(m))
        for i in range(m):
            # alpha 每次迭代时更新，这里j是迭代次数，i是样本点的下标，alpha每次减少，但不会为0
            alpha = 4/(1.0 + j + i) + 0.01
            # 随机选取样本来更新回归系数，减少系数周期性的波动
            randIndex = int(np.random.uniform(0,len(dataIndex)))
            h = sigmoid(sum(dataMatrix[randIndex] * weights))
            error = classLables[randIndex] - h
            weights = weights + alpha * error * dataMatrix[randIndex]
            del(dataIndex[randIndex])
    return weights

下面来看看同一个数据集上的分类效果，在main.py中添加：

    dataArr,labelMat = Logistic.loadDataSet()
    weights = Logistic.stocGradAscent1(np.array(dataArr),labelMat)
    print(weights)
    Logistic.plotBestFit(weights)

输出结果：

从输出结果来看该分割线达到了与gradAscent()差不多的效果，但所使用的计算量更少了。

当然也可以通过stocGradAscent()的第3个参数来对此进行修改，例如：

    weights = Logistic.stocGradAscent1(np.array(dataArr),labelMat, 500)

输出结果：

目前，已经学习了几个优化算法，当然，还有很多优化算法值得学习，但接下来就先用随机梯度上升算法来解决病马的生死预测问题。

3.示例：从疝气病症预测病马的死亡率

这里的数据包含368个样本和28个特征。影响马死亡的不只是疝气病症，可能还有其它因素。此外，数据集中30%的值是缺失的，所以要首先处理数据缺失问题。

示例：使用Logistic回归估计马疝病的死亡率

收集数据：给定数据文件
准备数据：用Python解析文本文件并填充缺失值
分析数据：可视化并观察数据
训练算法：使用优化算法，找到最佳的系数
测试算法：为了量化回归的效果，需要观察错误率。根据错误率决定是否回退到训练阶段，通过改变迭代的次数和步长等参数来得到更好的回归系数
使用算法：实现一个简单的命令行程序来收集马的症状并输出预测结果并非难事

3.1准备数据：处理数据中的缺失值

数据中的缺失值是个非常棘手的问题，有很多文献都致力于解决这个问题。比如：某个样本20个特征缺少了某个特征，但还有19个特征，那它是否还有用？结果是当然可用，因为有时候数据相当昂贵，扔掉和重新获取都是不可取的，所以要采用一些方法来解决这个问题。

一些可选的做法：

使用可用特征的均值来填补缺失值
使用特殊值来填补缺失值，如-1
忽略有缺失值的样本
使用相似样本的均值添加缺失值
使用另外的机器学习算法预测缺失值

在预处理阶段做两件事：

第一，所有的缺失值必须用一个实数值来替换，因为我们使用的NumPy数据类型不允许包含缺失值。这里选择实数0来替换所有缺失值，恰好能适用于Logistic回归。回归系数的更新公式是：

weights = weights + alpha * error * dataMatrix[randIndex]

如果dataMatrix的某特征对应值为0，那么该特征的系数将不做更新，即：

weights = weights

另外，由于sigmoid(0)=0.5,即它对结果的预测不具有任何倾向性，因此也不会对误差造成任何影响。

第二，如果在测试数据集中发现了一条数据的类别标签已经缺失，那么简单的做法就是将该条数据丢弃，这是因为类别标签于特征不同，很难确定采用某个合适的值来代替。

此外采用类似kNN的方法，这种方法可能不太行。

3.2测试算法：用Logistic回归进行分类

在Logistic.py中添加代码：

# 以回归系数和特征向量作为输入来计算对应的Sigmoid值
def classifyVector(inX, weights):
    prob = sigmoid(sum(inX * weights))
    # 如果Sigmoid值大于0.5，返回1，否则返回0
    if prob > 0.5: return 1.0
    else: return 0.0

# 用来打开测试集和训练集，并对数据进行格式化处理
def colicTest():
    # 导入训练集
    frTrain = open('horseColicTraining.txt')
    frTest = open('horseColicTest.txt')
    trainingSet = []
    trainingLabels = []
    for line in frTrain.readlines():
        currLine = line.strip().split('\t')
        lineArr = []
        for i in range(21):
            lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabels.append(float(currLine[21]))
    # 计算回归系数向量，这里可以自由设定迭代次数
    trainWeights = stocGradAscent1(np.array(trainingSet), trainingLabels, 500)
    errorCount = 0
    numTestVec = 0.0
    for line in frTest.readlines():
        numTestVec += 1.0
        currLine = line.strip().split('\t')
        lineArr = []
        for i in range(21):
            lineArr.append(float(currLine[i]))
        if int(classifyVector(np.array(lineArr),trainWeights)) != int(currLine[21]):
            errorCount += 1
    # 计算分类错误率
    errorRate = (float(errorCount/numTestVec))
    print("这个训练集的分类错误率是: %f" % errorRate)
    return errorRate

# 调用colicTest()10次并求其平均值
def multiTest():
    numTests = 10
    errorSum = 0.0
    for k in range(numTests):
        errorSum += colicTest()
    print("训练 %d 之后，平均分类错误率是： %f" % (numTests,errorSum/float(numTests)))

然后在更新main.py：

import Logistic
import numpy as np

if __name__ == '__main__':
    dataArr,labelMat = Logistic.loadDataSet()
    # weights = Logistic.gradAscent(dataArr,labelMat)
    # print(weights)
    # Logistic.plotBestFit(weights)
    # weights = Logistic.stocGradAscent0(np.array(dataArr),labelMat)
    # print(weights)
    # Logistic.plotBestFit(weights)
    # weights = Logistic.stocGradAscent1(np.array(dataArr),labelMat, 500)
    # print(weights)
    # Logistic.plotBestFit(weights)
    Logistic.multiTest()

输出结果：

有个警告：

RuntimeWarning: overflow encountered in exp return 1.0/(1+ np.exp(-inX))

查询得知优化Sigmoid函数就好了，数据溢出的问题：

# 当x是一个非常小的负数时，exp(-x)会过大，导致溢出，下面进行优化：
def sigmoid(x):
    if x >= 0:
        return 1.0/(1 + np.exp(-x))
    else:
        # 原式分子分母同乘exp(x)这个很小的数，可以防止数据溢出
        return np.exp(x)/(1 + np.exp(x))