Fork me on GitHub

LogisticRegression

logistics回归分类器、最优化理论初步、梯度下降最优化算法、数据中的缺失项处理。

logistics回归

假设现在有一些数据点,我们用一条直线进行拟合(该直线为最佳拟合直线),这个拟合的过程就称作回归。利用logistics回归的主要思想是:根据现有的数据对分类边界线建立回归公式,以此进行分类。这里的“回归”一词源于最佳拟合参数,表示要找到最佳拟合参数集。训练分类器时的做法就是寻找最佳拟合参数,使用的是最优算法。

logistics回归的一般流程:

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

基于logistics回归和sigmoid函数的分类

  • 优点:计算代价不高,易于理解和实现。
  • 缺点:容易欠拟合,分类精度不高。
  • 适用数据类型:数值型和标称型数据

为了实现logistics回归分类器,我们可以在每个特征上都乘以一个回归系数,然后把所有的结果值相加,将这个总和代入sigmoid函数中,进而得到一个0~1之间的数值。任何大于0.5的数据即被分入1类,小于0.5的即被分入0类。所以,logistics回归也可以被看成是一种概率估计。

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

sigmoid函数的输入记为z,即z=w0x0+w1x1+w2x2+…+wnxn,如果用向量表示即为z=wTx,它表示将这两个数值向量对应元素相乘然后累加起来。其中,向量x是分类器的输入数据,w即为我们要拟合的最佳参数,从而使分类器预测更加准确。也就是说,logistic回归最重要的是要找到最佳的拟合参数。

梯度上升法

梯度上升法(等同于我们熟知的梯度下降法,前者是寻找最大值,后者寻找最小值),它的基本思想是:要找到某函数的最大值,最好的方法就是沿着该函数的梯度方向搜寻。

假设有100个样本点,每个样本有两个特征:x1和x2.此外为方便考虑,我们额外添加一个x0=1,将线性函数z=wTx+b转为z=wTx(此时向量w和x的维度均价1).那么梯度上升法的伪代码如下:

  初始化每个回归系数为1

  重复R次:

    计算整个数据集梯度

    使用alpha*gradient更新回归系数的向量

  返回回归系数

logistics回归梯度上升优化算法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
# 加载数据
def loadDataSet():
# 创建初级与标签列表
dataMat = []; labelMat = []
# 打开文本数据集
fr = open('testSet.txt')
# 遍历文本的每一行
for line in fr.readlines():
# 对当前行除去首尾空格之后按空格进行分离
lineArr = line.strip().split()
# 将每一行的两个特征x1,x2,加上x0=1,组成列表并添加到数据集列表中
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+exp(-inX))

#梯度上升法更新最优拟合参数
#@dataMatIn:数据集
#@classLabels:数据标签
def gradAscent(dataMatIn, classLabels):
# 将数据集列表转为Numpy矩阵
dataMatrix = mat(dataMatIn)
# 将数据集标签列表转为Numpy矩阵,并转置
labelMat = mat(classLabels).transpose()
# 获取数据集矩阵的行数和列数
m,n = shape(dataMatrix)
# 学习率
alpha = 0.001
# 最大迭代次数
maxCycles = 500
# 初始化权值参数向量每个维度均为1.0
weights = ones((n,1))
# 循环迭代以及向量化计算
for k in range(maxCycles):
# 求当前的sigmoid函数预测概率
h = sigmoid(dataMatrix*weights)
# 计算真实类别和预测类别的差值
error = (labelMat - h)
# 更新权值参数,更新公式通过求导得到
weights = weights + alpha * dataMatrix.transpose()* error
return weights

分析数据:画出决策边界

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
# 画出数据集合logistics回归最佳拟合直线
def plotBestFit(weights):
import matplotlib.pyplot as plt
dataMat,labelMat=loadDataSet()
dataArr = array(dataMat)
n = 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 = 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()

训练算法:随机梯度上升

我们知道梯度上升法每次更新回归系数都需要遍历整个数据集,当样本数量较小时,该方法尚可,但是当样本数据集非常大且特征非常多时,那么随机梯度下降法的计算复杂度就会特别高。一种改进的方法是一次仅用一个样本点来更新回归系数,即随机梯度上升法。由于可以在新样本到来时对分类器进行增量式更新,因此随机梯度上升法是一个在线学习算法。与“在线学习”相对应,一次处理所有数据被称作“批处理”。随机梯度上升法可以写成如下伪代码:

所有回归系数初始化为1

对数据集每个样本

  计算该样本的梯度

  使用alpha*gradient更新回顾系数值

返回回归系数值

1
2
3
4
5
6
7
8
9
10
11
12
13
14
# 随机梯度上升算法
def stocGradAscent0(dataMatrix, classLabels):
# 获取数据集的行数和列数
m,n = shape(dataMatrix)
# 设置学习率为0.01
alpha = 0.01
# 初始化权值向量各个参数为1.0
weights = ones(n)
# 循环m次,每次选取数据集一个样本更新参数
for i in range(m):
h = sigmoid(sum(dataMatrix[i]*weights))
error = classLabels[i] - h
weights = weights + alpha * error * dataMatrix[i]
return weights

拟合效果尚可,但有相当多的样本被错误分类。但是相比于之前的算法,其计算量少了很多。我们可以加大计算量来改善拟合效果。我们知道,评判一个优化算法的优劣的可靠方法是看其是否收敛,也就是说参数的值是否达到稳定值。此外,当参数值接近稳定时,仍然可能会出现一些小的周期性的波动。这种情况发生的原因是样本集中存在一些不能正确分类的样本点(数据集并非线性可分),所以这些点在每次迭代时会引发系数的剧烈改变,造成周期性的波动。显然我们希望算法能够避免来回波动,从而收敛到某个值,并且收敛速度也要足够快。

改进的随机梯度上升算法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
# 改进的随机梯度上升算法
#@dataMatrix:数据集列表
#@classLabels:标签列表
#@numIter:迭代次数,默认150
def stocGradAscent1(dataMatrix, classLabels, numIter=150):
# 获取数据集的行数和列数
m,n = shape(dataMatrix)
# 初始化权值参数向量每个维度均为1
weights = ones(n)
# 循环每次迭代次数
for j in range(numIter):
# 获取数据集行下标列表
dataIndex = list(range(m))
# 遍历行列表
for i in range(m):
# 每次更新参数时设置动态的步长,保证随着更新次数的增多,步长变小,避免在最小值处徘徊
alpha = 4/(1.0+j+i)+0.0001
# 随机获取样本
randIndex = int(random.uniform(0,len(dataIndex)))
# 计算权值更新
h = sigmoid(sum(dataMatrix[randIndex]*weights))
error = classLabels[randIndex] - h
weights = weights + alpha * error * dataMatrix[randIndex]
# 选取该样本后,将该样本下标删除,确保每次迭代时只使用一次
del(dataIndex[randIndex])
return weights

上述代码中有两处改进的地方:

1 alpha在每次迭代更新是都会调整,这会缓解数据波动或者高频运动。此外,alpha还有一个常数项,目的是为了保证在多次迭代后仍然对新数据具有一定的影响,如果要处理的问题是动态变化的,可以适当加大该常数项,从而确保新的值获得更大的回归系数。

2 第二个改进的地方是选择随机的样本对参数进行更新,由于增加了随机性,这就防止参数值发生周期性的波动。

并且采用上述改进算法后,收敛速度更快

总测试

1
2
3
4
5
6
7
8
9
10
11
12
13
# 测试loadDataSet()/sigmoid(inx)/gradAscent(dataMatIn,classLabels)
dataArr,labelMat = loadDataSet()
weight = gradAscent(dataArr,labelMat)
print(weight)
# 测试画线函数
# getA()函数与mat()函数的功能相反,是将一个numpy矩阵转换为数组
plotBestFit(weight.getA())
# 测试随机梯度上升
weight = stocGradAscent0(array(dataArr),labelMat)
plotBestFit(weight)
# 测试改进的随机梯度上升
weight = stocGradAscent1(array(dataArr),labelMat)
plotBestFit(weight)

示例:从疝气病预测病马死亡率

现有数据集有100个样本和20个特征,但是数据存在一定的问题,即数据集有30%的缺失,因此,我们在对病马进行预测死亡率前,首先要解决数据的缺失问题。

  我们可能会遇到数据缺失的情况,但有时候数据相当昂贵,扔掉和重新获取均不可取,这显然是会带来更多的成本负担,所以我们可以选取一些有效的方法来解决该类问题。比如:

  1 使用可用特征的均值填补缺失值

  2 使用特殊值来填补缺失的特征,如-1

  3 忽略有缺失值的样本

  4 使用相似样本的平均值填补缺失值

  5 使用另外的机器学习算法预测缺失值

  这里我们根据logstic回归的函数特征,选择实数0来替换所有缺失值,而这恰好能适用logistic回归。因此,它在参数更新时不会影响参数的值。即如果某特征对应值为0 ,那么由公式w:w+alpha*gradient,可知w不会发生改变。

  此外,由于sigmoid(0)=0.5,表面该特征对结果的预测不具有任何倾向性,因此不会对误差造成影响。

  当然,如果是发生有样本的类标签缺失的情况,此时我们最好的办法是将该样本舍弃,这是因为标签与特征不同,我们很难确定采用某个合适的值替换掉。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
# 分类决策函数
def classifyVector(inX, weights):
prob = sigmoid(sum(inX*weights))
if prob > 0.5: return 1.0
else: return 0.0

#logistic回归预测算法
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列表
lineArr.append(float(currLine[i]))
# 将该样本的特征向量添加到数据集列表
trainingSet.append(lineArr)
# 将该样本标签存入标签列表
trainingLabels.append(float(currLine[21]))
# 调用随机梯度上升法更新logistic回归的权值参数
trainWeights = stocGradAscent1(array(trainingSet), trainingLabels, 5000)
# 统计测试数据集预测错误样本数量和样本总数
errorCount = 0; numTestVec = 0.0
# 遍历测试数据集的每个样本
for line in frTest.readlines():
# 样本总数加1
numTestVec += 1.0
# 对当前行进行处理,分割出各个特征及样本标签
currLine = line.strip().split('\t')
# 新建特征向量
lineArr =[]
# 将各个特征构成特征向量
for i in range(21):
lineArr.append(float(currLine[i]))
# 利用分类预测函数对该样本进行预测,并与样本标签进行比较
if int(classifyVector(array(lineArr), trainWeights))!= int(currLine[21]):
# 如果预测错误,错误数加1
errorCount += 1
# 计算测试集总的预测错误率
errorRate = (float(errorCount)/numTestVec)
# 打印错误率大小
print("the error rate of this test is: %f" % errorRate)
# 返回错误率
return errorRate

#多次测试算法求取预测误差平均值
def multiTest():
# 设置测试次数为10次,并统计错误率总和
numTests = 10
errorSum=0.0
for k in range(numTests):
errorSum += colicTest()
# 打印出测试10次预测错误率平均值
print("after %d iterations the average error rate is: %f" % (numTests, errorSum/float(numTests)))

测试:

1
2
# 进行分类
multiTest()

the error rate of this test is: 0.328358
the error rate of this test is: 0.432836
the error rate of this test is: 0.402985
the error rate of this test is: 0.343284
the error rate of this test is: 0.402985
the error rate of this test is: 0.313433
the error rate of this test is: 0.343284
the error rate of this test is: 0.402985
the error rate of this test is: 0.313433
the error rate of this test is: 0.328358
after 10 iterations the average error rate is: 0.361194

这样,经过10次的迭代,平均误差为36%,显然这跟30%的数据缺失有一定的关系,当然,如果我们调整合适的步长和迭代次数,相信错误率会有所下降。

小结

logistic回归的目的是寻找一个非线性函数sigmoid的最佳拟合参数,从而来相对准确的预测分类结果。为了找出最佳的函数拟合参数,最常用的优化算法为梯度上升法,当然我们为了节省计算损耗,通常选择随机梯度上升法来迭代更新拟合参数。并且,随机梯度上升法是一种在线学习算法,它可以在新数据到来时完成参数的更新,而不需要重新读取整个数据集来进行批处理运算。

参考

1.《机器学习实战》. Peter Harrington

-------------The End-------------