Machine Learning in Action 读书笔记
第5章 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
=
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)
该公式将一直被迭代执行,直至达到某个迭代停止条件为止,如:
- 迭代次数达到某个指定值
- 算法达到某个可以允许的误差范围
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样式:
为了方便计算,上面函数将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的关系:
分类效果可以,但是数据集虽小,就进行了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
随机梯度上升算法与梯度上升算法区别:
- 后者的变量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] - h
weights = 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.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()