逻辑回归实战（Python批量梯度下降，随机梯度下降，从疝气病症预测病马的死亡率）

最新推荐文章于 2024-04-23 19:26:48 发布

阿雅Yage

最新推荐文章于 2024-04-23 19:26:48 发布

阅读量2.6k

点赞数 1

分类专栏： Machine Learning python

本文链接：https://blog.csdn.net/u012369535/article/details/97034272

版权

Machine Learning 同时被 2 个专栏收录

10 篇文章 3 订阅

订阅专栏

python

6 篇文章 0 订阅

订阅专栏

一、逻辑回归代码

逻辑回归其实是求一个分类器，是二分类问题，其利用一个sigmoid函数去定义样本属于正类的概率，sigmoid函数的输入值z利用了线性回归的wx，所得出的数值是在[0,1]内的概率，即样本属于正类的概率值。在实践中，我们需要从训练样本中学习出一个w和b，这就是需要的逻辑回归模型，w和b包含在一个权重系数向量w中，w的更新是依据对损失函数即负对数似然函数求偏导数得梯度，然后每轮迭代中w根据负梯度来更新，梯度是一个向量，表示函数值上升最快的方向，由于这里的损失函数（负对数似然函数）是凸函数，因此在训练中求使其值最小的w值，就利用到梯度下降法（负梯度方向），当迭代到一定次数后所得到的w即为模型，预测时将利用sigmoid函数（此时输入值即wx，w是所得模型，x是预测样本）求得x属于正样本的概率（回归），若值大于阈值，即判定为正类（分类）。

关于回归预测模型和损失函数的定义，以及梯度更新公式的推导，
参考：https://mp.weixin.qq.com/s/BJxdDz7DQg5QIWzzgNQrgA
其梯度公式最后少了一个负号

其中每轮迭代中权重系数矩阵weights的更新运用到训练样本数据矩阵dataMatrix，训练样本标签矩阵labelMat和训练样本预测值矩阵h，这种向量化或矩阵化运算提高计算效率，避免for循环带来的计算开销，其矩阵计算的原理如下：这里样本特征是2，加上偏置项1，所以有3个特征，故其权重向量的维度也是3
在这里插入图片描述

完整的逻辑回归代码如下：

# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
def loadDataSet(file_name):
    """
    :desc: 加载并解析数据
    :param file_name: 文件名称，要解析的文件所在路径
    :return:
        dataMat -- 原始数据的特征
        labelMat -- 原始数据的标签，也就是每条样本对应的类别
    """
    dataMat = []
    labelMat = []
    fr = open(file_name)
    for line in fr.readlines():
        lineArr = line.strip().split()
        if len(lineArr) == 1:
            continue
        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])  # 此处为每个样本添加一个特征1.0
        labelMat.append(int(lineArr[2]))
    return dataMat, labelMat


# sigmoid跃阶函数
def sigmoid(x):

    return 1.0 / (1 + np.exp(-x))

# 批量梯度下降
def gradDescent(dataMatIn, classLabels):
    """
    :param dataMatIn: 是一个2维Numpy数组，每行代表每个训练样本，每列分别代表样本不同特征
    :param classLabels: 是一个List，是训练样本的类别标签，它是一个1*100的行向量。为了便于矩阵计算，需要将该行向量转换为列向量，
                        即将原向量转置后赋值给labelMat
    :return: Numpy数组，权重系数
    """
    # 将numpy数组dataMatIn转化为numpy矩阵[[1 1 2]
    #                                   [1 1 2]
    #                                    ...  ]
    dataMatrix = np.mat(dataMatIn) # 转换为Nupmy矩阵
    # 将标签List转化为行矩阵[[0 1 0 1 0 1 ...]]，并转置成列矩阵[[0]
    #                                                [1]
    #                                                [0]
    #                                                ..]
    labelMat = np.mat(classLabels).transpose()
    # m是训练样本个数，n是特征数
    m, n = np.shape(dataMatrix)
    # alpha代表移动步长
    alpha = 0.001
    # 迭代次数
    maxCycles = 500
    # 初始权重系数矩阵,长度与特征数相同[[1.]
    #                              [1.]
    #                              [1.]]
    weights = np.ones((n, 1))  # numpy数组
    for k in range(maxCycles):
        h = sigmoid(dataMatrix*weights)  # 对应于训练样本的预测值 100*1的列矩阵
        errors = labelMat - h # 100*1的列矩阵
        weights = weights + alpha * dataMatrix.transpose() * errors
    return np.array(weights)


# 随机梯度下降
# 梯度下降算法在每次更新权重系数时都需要遍历整个数据集，计算复杂度较高
# 随机梯度下降每次更新权重只需要一个样本
def stocGradDescent(dataMatrix, classLabels):
    """

    :param dataMatrix: Numpy2维数组
    :param classLabels: List
    :return: numpy数组，权重系数
    """
    m, n = np.shape(dataMatrix) # m是第一维度的元素个数，即样本数，n是第二维度元素个数，即特征数
    dataMat = np.mat(dataMatrix)
    alpha = 0.001
    weights = np.ones((n, 1))
    print(np.shape(weights))
    print(np.shape(dataMat[1]))
    for i in range(m):
        h = sigmoid(np.sum(dataMat[i]*weights))  # 取一条样本与权重系数相乘，这里的h是一个数值而不再是矩阵了
        error = classLabels[i] - h
        weights = weights + alpha * error * dataMat[i].transpose()
    return np.array(weights)

# 改进随机梯度下降
def stocGradDescent1(dataMatrix, classLabels, numIter=150):
    """
    :param dataMatrix: Numpy2维数组
    :param classLabels: List
    :param numIter: 迭代次数
    :return:
    """
    m, n = np.shape(dataMatrix)
    dataMat = np.mat(dataMatrix)
    weights = np.ones((n, 1))
    for j in range(numIter):
        dataIndex = list(range(m)) # 存放样本索引[0 1 2 .. m-1]
        for i in range(m):
            alpha = 4/(1.0+j+i)+0.0001  # alpha 会随着迭代不断减小
            # 随机取一个样本的索引,其值在[0，len(dataIndex))左闭右开区间范围内
            randIndex = int(np.random.uniform(0, len(dataIndex)))
            h = sigmoid(np.sum(dataMat[dataIndex[randIndex]]*weights))
            error = classLabels[dataIndex[randIndex]] - h
            weights = weights + alpha * error * dataMat[dataIndex[randIndex]].transpose()
            # 删除该随机索引
            del dataIndex[randIndex]
    return np.array(weights)

def plotBestFit(dataArr, labelMat, weights):
    """
    Desc: 将训练样本和训练出的权重系数可视化
    :param dataArr: numpy2维数组，样本数据的特征，[[1 1 2] [1 2 1] [1 1 2]]
    :param labelMat: List，样本数据的类别标签
    :param weights: numpy数组，权重系数
    :return: None
    """
    n = np.shape(dataArr)[0]  # 样本个数
    xcord1 = []; ycord1 = []  # 分别存放类别为1的样本特征
    xcord2 = []; ycord2 = []  # 分别存放类别为2的样本特征
    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('X')
    plt.ylabel('Y')
    plt.savefig("D://fig_2.png")
    plt.show()

def testLR():
    #  1.准备训练数据
    dataMat, labelMat = loadDataSet("D:\\testSet.TXT")
    #  2.训练模型
    dataArr = np.array(dataMat)
    # weights = gradDescent(dataArr, labelMat)
    weights = stocGradDescent1(dataArr, labelMat)

    print("The final weight is: \n", weights, "\n", type(weights), np.shape(weights))
    # 数据可视化
    plotBestFit(dataArr, labelMat, weights)


#  预测样本
def classifyVector(x, weights):
    """

    :param x: numpy matrix, 待预测的样本特征向量
    :param weights: numpy array, 根据训练样本学习由梯度下降得到的权重
    :return: 如果prob大于0.5返回1，否则返回0
    """
    prob = sigmoid(np.sum(x*weights))
    if prob > 0.5:
        return 1.0
    else:
        return 0.0

if __name__ == "__main__":
    testLR()

其中用于学习权重系数有三种方法：分别是批量梯度下降gradDescent、随机梯度下降stocGradDescent和改进随机梯度下降stocGradDescent1。批量梯度下降也就是常说的梯度下降法，梯度下降算法在每次更新回归系数时都需要遍历整个数据集，该方法在处理 100 个左右的数据集时尚可，但如果有数十亿样本和成千上万的特征，那么该方法的计算复杂度就太高了。一种改进方法是一次仅用一个样本点来更新回归系数，该方法称为随机梯度下降算法。相比之下，批量梯度下降比随机梯度下降具有更好的学习效果，学习得到的回归系数用作预测的效果更准确，下面是两种算法寻优后的分类超平面：

由于w0*x0+w1*x1+w2*x2=f(x)， x0最开始就设置为1， x2就是我们画图的y值，令f(x)=0，
所以： w0+w1*x+w2*y=0 => y = (-w0-w1*x)/w2 即分类超平面直线

批量梯度下降：
在这里插入图片描述

随机梯度下降：
在这里插入图片描述

随机梯度下降每次迭代只用到一条样本数据，其下降方向并不是全局最优的，也就是不一定向着最低点的方向，因此虽然迭代速度快，但是很容易陷入局部最优，最终的学习效果就没有批量梯度下降好。

因此对随机梯度下降作一次改进，代码里是stocGradDescent1，第一处改进为 alpha 的值。alpha 在每次迭代的时候都会调整，会随着迭代次数不断减少，但永远不会减小到 0，因为在计算公式中添加了一个常数项。

第二处修改为 randIndex 更新，这里通过随机选取样本来更新回归系数。这种方法每次随机从列表中选出一个值，然后从列表中删掉该值（再进行下一次迭代）。。

改进随机梯度下降：
在这里插入图片描述

在机器学习领域，梯度下降分为3类：
1、批量梯度下降BGD：每轮迭代计算所有样本数据，所有样本共同决定最优下降方向
优点：逼近全局最优
缺点：计算量大，批量完成，速度慢，样本数量大时不适用
2、随机梯度下降SGD：每轮迭代从训练样本中抽取一个样本进行计算更新，每次都不用遍历所有数据集
优点：速度快，逐条样本进行计算
缺点：陷入局部最优，下降幅度相对慢，需要迭代更多次数，因为每次选取的方向不一定是最优
3、小批量的梯度下降MBGD：每轮迭代计算所有数据中的子集，是1和2的中和方案，节省时间，并且寻优也准确

二、从疝气病症预测病马的死亡率

使用逻辑回归来预测患有疝病的马的存活问题。疝病是描述马胃肠痛的术语。然而，这种病不一定源自马的胃肠问题，其他问题也可能引发马疝病。这个数据集中包含了医院检测马疝病的一些指标，分为训练集和测试集，每条数据共22个字段，最后一个字段是其类别标签0或1，项目代码如下

病马数据集：
链接：https://pan.baidu.com/s/19qLbQ2TL757r0zyJjbLTUQ
提取码：gw3j

#  预测样本
def classifyVector(x, weights):
    """

    :param x: numpy matrix, 待预测的样本特征向量
    :param weights: numpy array, 根据训练样本学习由梯度下降得到的权重
    :return: 如果prob大于0.5返回1，否则返回0
    """
    prob = sigmoid(np.sum(x*weights))
    if prob > 0.5:
        return 1.0
    else:
        return 0.0


#  从疝气病症预测病马的死亡率
def colicTest():
    """
    Desc：读入测试集和训练集，并对数据进行格式化处理
    :return: errorRate 分类错误率
    """
    frTrain = open('D:\\horseColicTraining.txt')  # 训练样本
    frTest = open('D:\\horseColicTest.txt')  # 测试样本
    trainingSet = [] # [[],[],[],...]
    trainingLabels = []
    for line in frTrain.readlines():
        currLine = line.strip().split('\t') # 得到每一条数据字符串列表（共22个字段，最后一个字段为标签）
        lineArr = [] # 存储每条样本特征
        for i in range(len(currLine)-1):
           lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabels.append(float(currLine[len(currLine)-1]))
    # 使用改进后的随机梯度下降算法求最佳权重系数
    # trainWeights = stocGradDescent1(np.array(trainingSet), trainingLabels, 500) # 改进后随机梯度下降
    # trainWeights = stocGradDescent(np.array(trainingSet), trainingLabels) # 随机梯度下降
    trainWeights = gradDescent(np.array(trainingSet), trainingLabels) # 梯度下降
    errorCount = 0
    numTestVec = 0.0

    # 读取测试样本进行测试，计算分类错误的样本条数和最终的错误率
    for line in frTest.readlines():
        numTestVec += 1.0 # 测试样本数量加1
        currLine = line.strip().split('\t')
        lineArr = [] # 存储每条样本特征
        for i in range(len(currLine)-1):
           lineArr.append(float(currLine[i]))
        if int(classifyVector(np.mat(lineArr), trainWeights)) != int(currLine[len(currLine)-1]):
            errorCount += 1 # 预测错误，错误数加1

    errorRate = float(errorCount) / numTestVec
    print("测试的错误率是: ", errorRate)
    return errorRate
if __name__ == "__main__":
    colicTest()