Machine Learning in Action 读书笔记---第5章 Logistic回归

柔光的暖阳◎ 2022-09-03 00:18 168阅读 0赞

# Machine Learning in Action 读书笔记 #

第5章 Logistic回归

--------------------

### 文章目录 ###

*  Machine Learning in Action 读书笔记
 *  一、Logistic回归
 *   *  1.Logistic回归的一般过程
     *  2.Logistic回归的特点
     *  3.Logistic回归分类器的构建
 *  二、优化算法
 *   *  1.梯度上升法
     *  2.梯度下降法
     *  3.训练算法：利用梯度上升找到最佳参数
     *  4.分析数据：画出决策边界
     *  5.随机梯度上升算法
     *  6.改进的随机梯度算法
 *  三、示例：从疝气病预测病马的死亡率
 *   *  1.一般过程
     *  2. 准备数据：处理数据中的缺失值
     *  3.测试算法：用Logistic回归进行分类
 *  四、所有代码

--------------------

> **本章内容**：介绍几个最优化算法（梯度上升算法、随机梯度上升算法、改进的随机梯度上升算法），并利用它们训练出一个非线性函数用于分类。

# 一、Logistic回归 #

假设现在有一些数据点，用一条直线对这些点进行拟合，这条直线被称为**最佳拟合直线**，这个拟合过程就称为**回归**。  
利用Logistic回归进行分类的**主要思想**：根据现有数据对分类边界线建立回归公式，以此进行分类。  
训练分类器时的做法就是**寻找最佳拟合参数**，使用的是**最优化算法**。

## 1.Logistic回归的一般过程 ##

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

## 2.Logistic回归的特点 ##

*  **优点**：计算代价不高，易于理解和实现
 *  **缺点**：容易欠拟合，分类精度可能不高
 *  **适用数据类型**：数值型和标称型数据

## 3.Logistic回归分类器的构建 ##

首先在每个特征上都乘以一个回归系数，然后把所有的结果值相加，将这个总和代入**Sigmoid函数**中，进而得到一个范围在0~1之间的数值。任何大于0.5的结果被分入1类，小于0.5即被归入0类。所以，Logistic回归也可以被看成是一种**概率估计**。  
  
  
其中，在两个类的情况下，函数输出0或1，这种性质的函数称为海维塞德阶跃函数，或单位阶跃函数。  
Sigmoid函数是一种阶跃函数，sigmoid函数的计算公式如下：  
 σ ( z ) = 1 ( 1 + e − z ) σ(z) =\\frac\{1\}\{(1+e^\{-z\})\} σ(z)=(1\+e−z)1  
Sigmoid函数的输入记为z，由下面公式：  
 z = w 0 x 0 + w 1 x 1 + w 1 x 1 + . . . + w n x n z=w\_0x\_0+w\_1x\_1+w\_1x\_1+...+w\_nx\_n z=w0x0\+w1x1\+w1x1\+...\+wnxn  
向量**x**是分类器的输入数据，向量**w**是需要找到的最佳参数（系数）  
在本例子中，当z=0 时，可以用于二分类回归，代码如下：

y = (-weights[0]-weights[1]*x)/weights[2]  # 其中，y为特征x2,x为特征x1

# 二、优化算法 #

> 一个判断优化算法优劣的可靠方法是看它是否收敛

## 1.梯度上升法 ##

梯度上升法的**基本思想**是：要找到某个函数的最大值，最好的方法是沿着该函数的梯度方向探寻。  
**梯度算子**：也就是所说的移动方向，梯度算子总是指向函数值增长最快的方向。  
**步长**：移动量的大小  
梯度上升法的迭代公式如下：  
 w : = w + α ∇ w f ( w ) w := w+α∇\_wf(w) w:=w\+α∇wf(w)  
该公式将一直被迭代执行，直至达到某个**迭代停止条件**为止，如：

1.  迭代次数达到某个指定值
2.  算法达到某个可以允许的误差范围

## 2.梯度下降法 ##

梯度下降法的迭代公式如下：  
 w : = w − α ∇ w f ( w ) w := w-α∇\_wf(w) w:=w−α∇wf(w)  
梯度上升算法用来求函数的最大值，而梯度下降算法用来求函数的最小值。

## 3.训练算法：利用梯度上升找到最佳参数 ##

梯度上升的伪代码如下：

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

Logistic回归梯度上升优化算法代码：

'''下面用到了三种容易混淆的数据类型：
        列表：python自带的数据类型，特征为：list
        数组：numpy中的数据类型，特征为：array
        矩阵：numpy中的数据类型，特征为：mat 或 matrix
    '''
    
    '''便利函数，用于打开文本文件并逐行读取'''
    def loadDataSet():
        dataMat = []; labelMat = []
        fr = open('testSet.txt')
        for line in fr.readlines():
            lineArr = line.strip().split()
            dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])]) # 为了方便计算，初始化了一个特征X0，赋初值为1
            labelMat.append(int(lineArr[2]))  # 存储样本类型标签，本例子中分为0和1两类
        return dataMat, labelMat
    
    '''sigmoid函数：坐标轴跨度大的情况下，在0附近位置可以看成是一个阶跃函数'''
    def sigmoid(inX):  # 输入的是一个100 * 1维的矩阵
        return 1.0/(1 + np.exp(-inX))
    
    '''梯度上升算法：迭代五百次得到最佳权重'''
    def gradAscent(dataMatIn, classLabels): # dataMatIn:包含三个特征X1，X2,X3的二维列表； classLabels:每个样本的标签 一维列表
        dataMatrix = np.mat(dataMatIn)   # 将二维列表转化为二维矩阵
        # print(dataMatrix)
        labelMat = np.mat(classLabels).transpose() # transpose()在不指定参数默认为转置
        # print(labelMat)
        m, n = np.shape(dataMatrix)  # 100 * 3
        alpha = 0.001
        maxCycles = 500   # 迭代次数
        weights = np.ones((n, 1)) # 初始化每个特征的权重参数为1
        # print(weights)
        for k in range(maxCycles):
            # print(type(dataMatrix * weights))
            h = sigmoid(dataMatrix * weights)  # dataMatrix:100 * 3   weights:3 * 1
            # print(h)  # h:100 * 1
            error = (labelMat - h)   # 计算真实类别与预测类别的差值  100 * 1
            weights = weights + alpha * dataMatrix.transpose() * error  # 按照上面计算的差值的方向调整回归系数   3*1 + 3*100 * 100*1 = 3*1
        return weights  # 返回训练好的回归系数

使用数据集testSet.txt样式：  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70]  
为了方便计算，上面函数将X0的值设为1.0，此时，w0相当于截距。

## 4.分析数据：画出决策边界 ##

画出数据集和Logistic回归最佳拟合直线的函数

'''分析数据：画出决策的边界'''
    def plotBestFit(weights):
        dataMat, labelMat = loadDataSet()
        dataArr = np.array(dataMat)  # 将列表转化为数组
        n = np.shape(dataArr)[0]  # 100
        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)  # arange(start, stop, step)
        y = (-weights[0]-weights[1]*x)/weights[2]
        # print(x)
        # print(y)
        ax.plot(x, y)
        plt.xlabel('X1')
        plt.ylabel('X2')
        plt.show()

如下图，画出了x1和x2的关系：  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 1]  
分类效果可以，但是数据集虽小，就进行了300次乘积，计算量很大，需要改进。

## 5.随机梯度上升算法 ##

梯度上升算法在每次更新回归系数时都需要遍历整个数据集，而**随机梯度上升算法可以一次仅用一个样本点来更新回归系数**，这样在新样本到来时可以对分类器进行增量式更新，因此随机梯度上升算法是一个在线学习算法。与在线学习相对应，一次处理所有数据被称作是批处理。在线算法可以在新数据到来时就完成参数更新，而不需要重新读取整个数据集来进行批处理算法。  
  
  
随机梯度上升算法伪代码：

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

随机梯度上升算法实现代码：

'''随机梯度上升算法'''
    def stocGradAscent0(dataMatrix, classLabels):
        m, n = np.shape(dataMatrix)
        alpha = 0.01
        weights = np.ones(n)
        # print(weights)   #[1. 1. 1.]
        for i in range(m):
            h = sigmoid(sum(dataMatrix[i]*weights))
            error = classLabels[i] - h
            weights = weights + alpha * error * dataMatrix[i]
        return weights

随机梯度上升算法与梯度上升算法区别：

1.  后者的变量h和误差error都是向量，而前者全是数值
2.  前者没有矩阵的转换过程，所有变量的数据类型都是Numpy数组

如下图，画出了x1和x2的关系：效果不是很好，但是判断优化算法优劣的可靠方法是看它是否收敛。此方法收敛效果可以，但是波动较大。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 2]

## 6.改进的随机梯度算法 ##

'''改进的随机梯度上升算法'''
    def stocGradAscent1(dataMatrix, classLabels, numIter=150):
        m, n = np.shape(dataMatrix)
        weights = np.ones(n)
        for j in range(numIter):  # j 为迭代次数
            dataIndex = list(range(m))
            for i in range(m):   # i 为样本点的下标，即本次迭代第 i个选出样本
                # alpha会随着迭代次数不断减少，但是因为有常数项，所以永远不会减小到0，这样可以保证在多次迭代之后新数据仍然具有一定的影响
                alpha = 4/(1.0+j+i)+0.01  # alpha每次迭代时需要调整，可以缓解数据波动或者高频波动
                # 通过随机选取样本来更新回归系数，可以减少周期性的波动
                randIndex = int(random.uniform(0, len(dataIndex)))  # 随机选取更新, randIndex表示样本在矩阵中的位置
                h = sigmoid(sum(dataMatrix[randIndex]*weights))
                error = classLabels[randIndex] - h
                weights = weights + alpha * error * dataMatrix[randIndex]
                del(dataIndex[randIndex])
        return weights

stocGradAscent1中加入了样本随机选择机制，算法的收敛速度更快，波动小，分类图如下：

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 3]

# 三、示例：从疝气病预测病马的死亡率 #

## 1.一般过程 ##

1.  收集数据：给定数据文件
2.  准备数据：用python解析文本文件并**填充缺失值**
3.  分析数据：可视化并观察数据
4.  **训练算法**：使用优化算法，找到最佳系数
5.  **测试算法**：为了量化回归的效果，需观察错误率，根据错误率决定是否回退到训练阶段，通过改变迭代的次数和步长等参数来得到更好的回归系数。
6.  使用算法：收集马的症状并输出预测结果

## 2. 准备数据：处理数据中的缺失值 ##

> 机器学习的一个重要问题就是如何处理缺失数据。

缺失值分为：

*  特征缺失（补）
 *  类标签缺失（删）

下面给出一些可选的处理数据中缺失值的做法：

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

## 3.测试算法：用Logistic回归进行分类 ##

'''从疝气病症预测病马的死亡率'''
    def classifyVector(inX, weights):   # 回归系数和特征向量作为参数
        prob = sigmoid(sum(inX*weights))
        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):  # 一共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('the error rate of this test is: %f' % errorRate)
        return errorRate
    
    def multiTest():   # 调用函数colicTest()十次，并计算平均值
        numTests = 10; errorSum = 0.0
        for k in range(numTests):
            errorSum += colicTest()
        print('after %d iterations the average error rate is: %f' % (numTests, errorSum/float(numTests)))

可以通过调整以上代码中的colicTest()中的迭代次数和stocGradAscent1()中的步长来优化，还可以降低平均错误率。

# 四、所有代码 #

'''logRegres.py'''
    
    import random
    
    import numpy as np
    import matplotlib.pyplot as plt
    
    '''下面用到了三种容易混淆的数据类型：
        列表：python自带的数据类型，特征为：list
        数组：numpy中的数据类型，特征为：array
        矩阵：numpy中的数据类型，特征为：mat 或 matrix
    '''
    
    '''便利函数，用于打开文本文件并逐行读取'''
    def loadDataSet():
        dataMat = []; labelMat = []
        fr = open('testSet.txt')
        for line in fr.readlines():
            lineArr = line.strip().split()
            dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])]) # 为了方便计算，初始化了一个特征X0，赋初值为1
            labelMat.append(int(lineArr[2]))  # 存储样本类型标签，本例子中分为0和1两类
        return dataMat, labelMat
    
    '''sigmoid函数：坐标轴跨度大的情况下，在0附近位置可以看成是一个阶跃函数'''
    def sigmoid(inX):  # 输入的是一个100 * 1维的矩阵
        return 1.0/(1 + np.exp(-inX))
    
    '''梯度上升算法：迭代五百次得到最佳权重'''
    def gradAscent(dataMatIn, classLabels): # dataMatIn:包含三个特征X1，X2,X3的二维列表； classLabels:每个样本的标签 一维列表
        dataMatrix = np.mat(dataMatIn)   # 将二维列表转化为二维矩阵
        # print(dataMatrix)
        labelMat = np.mat(classLabels).transpose() # transpose()在不指定参数默认为转置
        # print(labelMat)
        m, n = np.shape(dataMatrix)  # 100 * 3
        alpha = 0.001
        maxCycles = 500   # 迭代次数
        weights = np.ones((n, 1)) # 初始化每个特征的权重参数为1
        # print(weights)
        for k in range(maxCycles):
            # print(type(dataMatrix * weights))
            h = sigmoid(dataMatrix * weights)  # dataMatrix:100 * 3   weights:3 * 1
            # print(h)  # h:100 * 1
            error = (labelMat - h)   # 计算真实类别与预测类别的差值  100 * 1
            weights = weights + alpha * dataMatrix.transpose() * error  # 按照上面计算的差值的方向调整回归系数   3*1 + 3*100 * 100*1 = 3*1
        return weights
    
    '''分析数据：画出决策的边界'''
    def plotBestFit(weights):
        dataMat, labelMat = loadDataSet()
        dataArr = np.array(dataMat)  # 将列表转化为数组
        n = np.shape(dataArr)[0]  # 100
        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)  # arange(start, stop, step)
        y = (-weights[0]-weights[1]*x)/weights[2]
        # print(x)
        # print(y)
        ax.plot(x, y)
        plt.xlabel('X1')
        plt.ylabel('X2')
        plt.show()
    
    '''随机梯度上升算法'''
    def stocGradAscent0(dataMatrix, classLabels):
        m, n = np.shape(dataMatrix)
        alpha = 0.01
        weights = np.ones(n)
        # print(weights)   #[1. 1. 1.]
        for i in range(m):
            h = sigmoid(sum(dataMatrix[i]*weights))
            error = classLabels[i] - h
            weights = weights + alpha * error * dataMatrix[i]
        return weights
    
    '''改进的随机梯度上升算法'''
    def stocGradAscent1(dataMatrix, classLabels, numIter=150):
        m, n = np.shape(dataMatrix)
        weights = np.ones(n)
        for j in range(numIter):  # j 为迭代次数
            dataIndex = list(range(m))
            for i in range(m):   # i 为样本点的下标，即本次迭代第 i个选出样本
                # alpha会随着迭代次数不断减少，但是因为有常数项，所以永远不会减小到0，这样可以保证在多次迭代之后新数据仍然具有一定的影响
                alpha = 4/(1.0+j+i)+0.01  # alpha每次迭代时需要调整，可以缓解数据波动或者高频波动
                # 通过随机选取样本来更新回归系数，可以减少周期性的波动
                randIndex = int(random.uniform(0, len(dataIndex)))  # 随机选取更新, randIndex表示样本在矩阵中的位置
                h = sigmoid(sum(dataMatrix[randIndex]*weights))
                error = classLabels[randIndex] - h
                weights = weights + alpha * error * dataMatrix[randIndex]
                del(dataIndex[randIndex])
        return weights
    
    '''从疝气病症预测病马的死亡率'''
    def classifyVector(inX, weights):   # 回归系数和特征向量作为参数
        prob = sigmoid(sum(inX*weights))
        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):  # 一共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('the error rate of this test is: %f' % errorRate)
        return errorRate
    
    def multiTest():   # 调用函数colicTest()十次，并计算平均值
        numTests = 10; errorSum = 0.0
        for k in range(numTests):
            errorSum += colicTest()
        print('after %d iterations the average error rate is: %f' % (numTests, errorSum/float(numTests)))
    
    if __name__ == '__main__':
        dataArr, labelMat = loadDataSet()
        # print(dataArr)
        # print(labelMat)
        myWeights = gradAscent(dataArr, labelMat)
        # print(myWeights)
        '''
        [[ 4.12414349]
         [ 0.48007329]
         [-0.6168482 ]]
        '''
    
        # plotBestFit(myWeights.getA())  # matrix.getA()可以将矩阵自身返回成一个n维数组对象，否则会报错：x and y must have same first dimension, but have shapes (60,) and (1, 60)
        '''测试随机梯度上升算法'''
        weights = stocGradAscent0(np.array(dataArr), labelMat)
        # print(weights)  # [ 1.01702007  0.85914348 -0.36579921]
        # plotBestFit(weights)
        '''测试改进的随机梯度上升算法'''
        weights = stocGradAscent1(np.array(dataArr), labelMat)
        # plotBestFit(weights)
    
        '''从疝气病症预测病马的死亡率'''
        multiTest()

[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70]: /images/20220829/827090ac7bda4944a30ab48f5c001bb4.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 1]: /images/20220829/763899aa0dfc4dc3b99c17d810398b5b.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 2]: /images/20220829/9de5d942275f4a338163d3bcca58ba63.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2lkX18zOQ_size_16_color_FFFFFF_t_70 3]: /images/20220829/c85f3cf7453e4b73bbd892f93fb0a52d.png