机器学习算法四——Logistic回归（2）（示例：从疝气病症预测病马的死亡率）

示例：从疝气病症预测病马的死亡率
数据包含368个样本和28个特征。疝病是描述马胃肠痛的术语。然而，这种病不一定源自马的胃肠问题，其他问题也可能引发马疝病。
该数据集中包含了医院检测马疝病的一些指标，**有的指标比较主观，有的指标难以测量，例如马的疼痛级别。**除此之外，还存在一个问题，数据集中有30%的值是缺失的。
下面将首先介绍如何处理数据集中的数据缺失问题，然后再利用Logistic回归和随机梯度上升算法来预测病马的生死。

一、准备数据：处理数据中的缺失值
可选做法：
●使用可用特征的均值来填补缺失值
●使用特殊值来填补缺失值，如-1
●忽略有缺失值的样本
●使用相似样本的均值添补缺失值
●使用另外的机器学习算法预测缺失值

数据预处理阶段：1、所有缺失值必须用一个实数值来替换，因为numpy数据类型不允许包含缺失值（这里选择0来替换所有缺失值，恰好能适用于Logistic回归）。这样做的直觉在于，我们需要的是一个在更新时不会影响系数的值。
回归系数的更新公式如下：
weights = weights + alpha * error * dataMatrix[randIndex]
如果dataMatrix的某特征对应值为0，那么该特征的系数将不做更新，即：weights = weights
另，由于sigmoid(0)=0.5，即它对结果的预测不具任何倾向性，因此不会对误差项造成任何影响。
2、如果在测试集中发现了一条数据的类别标签已经缺失，那么简单做法是将该条数据丢弃。

“干净”可用的数据集 + 一个不错的优化算法————>分类器————> 预测病马的生死问题

二、测试算法：用Logistic回归进行分类
Logistic回归方法：所需的只是把测试集上每个特征向量乘以最优化方法得来的回归系数，再将该乘积结果求和，最后输入到Sigmoid函数中。
如果对应的Sigmoid值大于0.5，就预测类别标签为1，否则为0。

Logistic回归分类函数代码如下：

#Logistic回归分类函数
def classifyVector(inX,weights):
    prob = sigmoid(sum(inX*weights)) #计算logistic回归预测概率
    if prob>0.5: #大于0.5预测为1,否则预测为0
        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.append(float(currLine[i])) #将该样本的特征存入lineArr列表
        trainingSet.append(lineArr) #将该样本的特征向量添加到数据集列表
        trainingLabels.append(float(currLine[21])) #将该样本标签存入标签列表      

    trainWeights=stocGradAscent1(np.array(trainingSet),trainingLabels,500) #调用随机梯度上升法更新logistic回归的权值参数 
    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 #如果预测错误，错误数加1   
    errorRate= (float(errorCount)/numTestVec) #计算测试集总的预测错误率
    print('the error rate of this test is: %f' %(errorRate)) #打印错误率大小    
    
    return errorRate #返回错误率

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

def main():
    multiTest()

数据最初有三个类别标签，分别代表马的三种情况：“仍成活”、“已经死亡” 和 “已经安乐死”。这里为了方便，将 “已经死亡” 和 “已经安乐死” 合并成 “未能存活” 这个标签。
整体来看，colicTest()具有完全独立的功能，多次运行得到的结果可能稍有不同，这时因为其中有随机的成分在里面。如果stocGradAscent1()函数中回归系数已经完全收敛，那么结果才将是确定的。

直接这样运行之后会出现一个警告：

>>>Warning (from warnings module):
  File "E:\test\\LogisticRegression\logRegres.py", line 15
    return 1.0/(1+ np.exp(-inx))
RuntimeWarning: overflow encountered in exp

可对sigmoid()函数进行如下优化：

#定义sigmoid函数
def sigmoid(inx):
    #return 1.0/(1+ np.exp(-inx))
    #优化20190310
    if inx>0:
        return 1.0/(1+ np.exp(-inx))
    else:
        return np.exp(inx)/(1+ np.exp(inx))

运行结果：

>>>the error rate of this test is: 0.358209
the error rate of this test is: 0.343284
the error rate of this test is: 0.328358
the error rate of this test is: 0.328358
the error rate of this test is: 0.328358
the error rate of this test is: 0.298507
the error rate of this test is: 0.358209
the error rate of this test is: 0.343284
the error rate of this test is: 0.417910
the error rate of this test is: 0.373134
after 10 iterations the average error rate is: 0.347761

从结果看到，10次迭代后的平均错误率为35%。事实上，这个结果并不差，因为有30%数据缺失。
当然，如果调整colicTest()中的迭代次数和stocGradAscent1()中的步长，平均错误率可以降到20%左右。（调试无果……）

总结

Logistic回归的目的：找一个非线性函数Sigmoid的最佳拟合参数，求解过程可由最优化算法来完成。
其中最常用的就是梯度上升算法，其又可简化为随机梯度上升算法（效果相当，占用更少资源，新数据到来时就完成在线参数更新，不需要重新读取整个数据集）。