Logistic回归算法实战

实战——预测病马的死亡率

前面，我简单介绍了「Logistic分类器」。

Logistic回归算法

今天，我们将进行实战操练：使用Logistic分类器来预测患有疝病的马的存活率。

疝病是描述马胃肠痛的术语。引发该病的因素（这里就是样本特征值）有很多。

一、数据集

原始数据集来源2010年1月11日的UCI机器学习数据库。该数据集包含368个样本和28个特征。

由于该数据集中有30%的值是缺失的，所以我们需要处理数据集中的缺失值。

几种可选的做法如下：

1.使用可用特征的均值来填补缺失值。
2.使用特殊值来填补缺失值。
3.忽略有缺失值的样本。
4.使用相似样本的均值填补缺失值。
5.使用另外的机器学习算法预测缺失值。

我们对原始数据做如下预处理：

1.所有特征缺失值必须用实数0来替换.

因为我们采用Logistic分类器，sigmod函数在x=0处的取值不影响最后的概率估计和更新回归系数。所以，我们恰好能采用实数0来进行填补缺失值。

2.若某条样本的类别标签缺失，则直接将该条数据丢弃。

因为，我们实在很难通过某种做法来判断该标签是什么。

预处理完成后的数据集有两个：horseColicTest.txt和horseColicTraining.txt。

处理后的数据集包含21个特征值和1个类别标签向量。

细心的朋友可能会问，前面不是28个特征值吗？为什么处理后变成21个了呢？这里我也没想通，可能数据缺失，书上确实是这样，但是这并不妨碍我们使用Logistic分类器。

二、随机梯度下降算法

数据集已经准备好了，那么我们该说说「随机梯度下降算法」了。

咦？我记得前面好像是梯度下降算法吧，这里怎么成随机的了？

至于为什么，咱们接着往下讲。

【一点小问题】

没错，前面我们确实讲的是梯度下降算法。但是，你们有没有发现，这个算法的计算量是非常大的。

梯度下降算法在每次更新回归系数时都需要遍历整个数据集。

回归系数更新公式：weights = weights - alpha * dataMatrix.transpose()* error

这里的dataMatrix是整个数据集，其计算量大致为数据集的行数x列数。

如果样本和特征值足够大，那么计算量是非常大的。

那有没有什么改进方法？有！终于轮到「随机梯度下降算法」出场了。

【随机梯度下降算法】

随机梯度下降算法：一次只用一个样本点来更新回归系数。

其代码如下：

'''
Created on Oct 27, 2010
Logistic Regression Working Module
@author: Peter
@modified by Albert on July 18,2020
'''def stocGradDecline0(dataArray, classLabels):"""随机梯度下降算法"""m,n = np.shape(dataArray)# 设置步长alpha = 0.01# 初始化回归系数weights = np.ones(n)   #initialize to all onesfor i in range(m):# 每次一个样本，注意这里是数组运算，需要手动sumh = sigmoid(sum(dataArray[i]*weights))# 预测值-真实值error = h - classLabels[i]# 更新回归系数，注意这里都是数值计算，不同于梯度下降算法的矩阵运算weights = weights - alpha * error * dataArray[i]return weights

那么，该随机梯度下降算法的分类效果究竟怎么样呢？

我们运行程序：

# 返回的数据集是矩阵类型
dataMat,labelMat = loadDataSet()
weights = stocGradDecline0(np.array(dataMat),labelMat)
# 注意，使用随机梯度下降算法时，weights是数组类型，不需要调用getA()
plotBestFit(weights)

效果如下：

诶呀，一看这结果，傻眼了，这效果明显没有前面的梯度下降算法好啊!

各位别急，请听我说。

你们是不是把迭代次数给忘了？

# 迭代次数 maxCycles = 500

上次的结果是在整个数据集上迭代500次后的结果，而这次仅仅是只是迭代了一次，或者压根没迭代。

判断一个优化算法是否优劣的可靠性方法，是看它最终的参数是否收敛。

对此，我们假定迭代次数为500，然后研究随机梯度下降算法在数据集上遍历500次时，回归系数与迭代次数的关系

核心代码如下：

'''
Created on Oct 27, 2010
Logistic Regression Working Module
@author: Peter
@modified by Albert on July 18,2020
'''
import numpy as np
import matplotlib.pyplot as plt
import logRegresdef stocGradDecline0(dataArr, classLabels):m,n = np.shape(dataArr)alpha = 0.5weights = np.ones(n)   #initialize to all ones# 设置500次迭代（遍历整个数据集500次），并记录（500x样本个数）次回归系数weightsHistory=np.zeros((500*m,n))for j in range(500):for i in range(m):h = logRegres.sigmoid(sum(dataArr[i]*weights))error = h - classLabels[i] weights = weights - alpha * error * dataArr[i]# 记录回归系数，共（500x样本个数）weightsHistory[j*m + i,:] = weightsreturn weightsHistory# 运行主程序
dataMat,labelMat=logRegres.loadDataSet()
myHist = stocGradDecline0(np.array(dataMat),labelMat)# 绘图
n = np.shape(np.array(dataMat))[0] #number of points to create
xcord1 = []; ycord1 = []
xcord2 = []; ycord2 = []markers =[]
colors =[]fig = plt.figure()
ax = fig.add_subplot(311)
type1 = ax.plot(myHist[:,0])
plt.ylabel('X0')
ax = fig.add_subplot(312)
type1 = ax.plot(myHist[:,1])
plt.ylabel('X1')
ax = fig.add_subplot(313)
type1 = ax.plot(myHist[:,2])
plt.xlabel('iteration')
plt.ylabel('X2')
plt.show()

运行结果如下：

可见，在500次迭代中，回归系数都有着不同程度的周期性波动（数据集不是线性可分的，每次迭代都会引发系数剧烈改变），且三个回归系数的收敛速度也不一样。

对此，我们希望算法能够避免周期性波动，并且尽快收敛到某一个值。

升级版核心代码如下，其余部分除调用处外，都相同：


def stocGradDecline1(dataArr, classLabels,numIter=100):"""改进的随机梯度下降算法"""m,n = np.shape(dataArr)weights = np.ones(n)   #initialize to all ones# weightsHistory 这里是为了画图需要，实际算法时，相关代码可以去掉weightsHistory=np.zeros((numIter*m,n))for j in range(numIter):dataIndex = list(range(m))for i in range(m):alpha = 4/(1.0+j+i)+0.01# 随机选择一个样本randIndex = int(np.random.uniform(0,len(dataIndex)))h = logRegres.sigmoid(sum(dataArr[randIndex]*weights))error = h -classLabels[randIndex]weights = weights - alpha * error * dataArr[randIndex]weightsHistory[j*m + i,:] = weights# 删掉已经使用过的样本索引值del(dataIndex[randIndex])print (weights)return weightsHistory

上述代码改进的地方有三处：

由于alpha的值会影响回归系数，所以我们让alpha在每次迭代时都进行调整。
随机选择某个样本来更新回归系数，这种做法可以减少回归系数的周期性波动。
添加了一个迭代参数，默认100。

运行效果如下：

可见，回归系数没有出现明显的周期性波动。另外，各参数的收敛速度也明显加快了。

我们设置迭代次数为3次。

# 返回的数据集是矩阵类型
dataMat,labelMat = loadDataSet()
weights = stocGradDecline1(np.array(dataMat),labelMat,3)
# 注意，使用随机梯度下降算法时，weights是数组类型，不需要调用getA()
plotBestFit(weights)

最后，我们来看看分类效果。

怎么样，效果不错吧，这可只是迭代了三次的结果啊，和迭代500次的梯度下降算法，随机梯度下降算法的计算量大大减少了。

另外，如果，你运行同样的代码，发现分类效果有一些不一样，不用担心，因为这是随机梯度下降算法，它每次选择的样本都是具有随机性的，不过其分类结果大致是差不多的。

多次运行上面的代码，分类结果如下：

三、病马预测

数据集准备好了，算法也优化好了，那么现在就可以进行我们的最终目标——病马预测了。

直接上代码：

1.sigmod 分类函数

'''
Created on Oct 27, 2010
Logistic Regression Working Module
@author: Peter
@modified by Albert on July 18,2020
'''def classifyVector(inX, weights):"""Sigmod分类函数"""prob = sigmoid(sum(inX*weights))# 预测结果大于0.5， 归为1类；反之，归为0类if prob > 0.5:return 1.0else: return 0.0

前面已经对sigmod函数进行了充分的讲解，这里就不再重复了。

2.疝气病预测

'''
Created on Oct 27, 2010
Logistic Regression Working Module
@author: Peter
@modified by Albert on July 18,2020
'''
def colicTest(numIter=100):"""病马预测函数"""frTrain = open('horseColicTraining.txt')frTest = open('horseColicTest.txt')trainingSet = []; trainingLabels = []for line in frTrain.readlines():currLine = line.strip().split('\t')lineArr =[]# 21个特征值for i in range(21):lineArr.append(float(currLine[i]))trainingSet.append(lineArr)# 类别标签trainingLabels.append(float(currLine[21]))trainWeights = stocGradDecline1(np.array(trainingSet), trainingLabels, numIter)# 定义错误率errorCount = 0numTestVec = 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 errorRate

上面这个函数，打开测试集和训练集，来训练和测试我们构建的分类器，并计算错误率。另外，改变alpha和迭代次数，错误率也是会跟着改变的。我们我们设置迭代次数为400次，alpha同上，计算10次预测的平均错误率。

运行代码：multiTest(400)

结果如下：

可见，平均错误率差不多35%左右，这个结果还行。另外，这还不是最优的结果，因为数据有30%的丢失，alpha和迭代次数也有待调试。最终，分类器的错误率还会下降很多。

四、总结

本文主要介绍了以下三部分内容：

1.有缺失值的数据集的相关处理方法。
2.一种优化的随机梯度下降算法。
3.Logsitic实战项目编程。

Logistic回归的目的是寻找一个非线性函数Sigmod的最佳拟合系数，求解过程可以由最优化算法来完成，其中最常用的算法就是随机梯度下降算法。

最后，Logistic回归算法系列就先介绍到这里。

参考资料：

《机器学习实战》， [美] ，Peter Harriington，人民邮电出版社.

Python 3 代码下载：
链接：https://pan.baidu.com/s/10yIYk_-j1pcfZOErVyB-TA
提取码：dwps