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回归的一般过程

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

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=w0x0+w1x1+w1x1+...+wnxnz=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+α∇wf(w)w := w+α∇_wf(w) w:=w+α∇wf(w)
该公式将一直被迭代执行，直至达到某个迭代停止条件为止，如：

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

2.梯度下降法

梯度下降法的迭代公式如下：
w:=w−α∇wf(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，赋初值为1labelMat.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 * 3alpha = 0.001maxCycles = 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 * 1error = (labelMat - h)   # 计算真实类别与预测类别的差值  100 * 1weights = weights + alpha * dataMatrix.transpose() * error  # 按照上面计算的差值的方向调整回归系数   3*1 + 3*100 * 100*1 = 3*1return weights  # 返回训练好的回归系数

使用数据集testSet.txt样式：

为了方便计算，上面函数将X0的值设为1.0，此时，w0相当于截距。

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

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

'''分析数据：画出决策的边界'''
def plotBestFit(weights):dataMat, labelMat = loadDataSet()dataArr = np.array(dataMat)  # 将列表转化为数组n = np.shape(dataArr)[0]  # 100xcord1 = []; 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的关系：

分类效果可以，但是数据集虽小，就进行了300次乘积，计算量很大，需要改进。

5.随机梯度上升算法

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

随机梯度上升算法伪代码：

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

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

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

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

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

如下图，画出了x1和x2的关系：效果不是很好，但是判断优化算法优劣的可靠方法是看它是否收敛。此方法收敛效果可以，但是波动较大。

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] - hweights = weights + alpha * error * dataMatrix[randIndex]del(dataIndex[randIndex])return weights

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

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

1.一般过程

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

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

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

缺失值分为：

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

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

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

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

'''从疝气病症预测病马的死亡率'''
def classifyVector(inX, weights):   # 回归系数和特征向量作为参数prob = sigmoid(sum(inX*weights))if prob > 0.5:return 1.0else:return 0.0def 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.0for line in frTest.readlines():numTestVec += 1.0currLine = 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 += 1errorRate = (float(errorCount)/numTestVec)print('the error rate of this test is: %f' % errorRate)return errorRatedef multiTest():   # 调用函数colicTest()十次，并计算平均值numTests = 10; errorSum = 0.0for 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 randomimport 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，赋初值为1labelMat.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 * 3alpha = 0.001maxCycles = 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 * 1error = (labelMat - h)   # 计算真实类别与预测类别的差值  100 * 1weights = weights + alpha * dataMatrix.transpose() * error  # 按照上面计算的差值的方向调整回归系数   3*1 + 3*100 * 100*1 = 3*1return weights'''分析数据：画出决策的边界'''
def plotBestFit(weights):dataMat, labelMat = loadDataSet()dataArr = np.array(dataMat)  # 将列表转化为数组n = np.shape(dataArr)[0]  # 100xcord1 = []; 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.01weights = np.ones(n)# print(weights)   #[1. 1. 1.]for i in range(m):h = sigmoid(sum(dataMatrix[i]*weights))error = classLabels[i] - hweights = 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] - hweights = 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.0else:return 0.0def 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.0for line in frTest.readlines():numTestVec += 1.0currLine = 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 += 1errorRate = (float(errorCount)/numTestVec)print('the error rate of this test is: %f' % errorRate)return errorRatedef multiTest():   # 调用函数colicTest()十次，并计算平均值numTests = 10; errorSum = 0.0for 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()