Logistic回归实战篇之预测病马死亡率(一)

作 者:崔家华 编 辑:李文臣

Python版本: Python3.x

运行平台 : Windows

IDE : Sublime text3

一、前言

本系列文章对梯度上升算法和改进的随机梯度上升算法进行了对比,总结了各自的优缺点,并对sklearn.linear_model.LogisticRegression进行了详细介绍。

二、改进的随机梯度上升算法

梯度上升算法在每次更新回归系数(最优参数)时,都需要遍历整个数据集。可以看一下我们之前写的梯度上升算法:

def gradAscent(dataMatIn, classLabels):
    dataMatrix = np.mat(dataMatIn) 
           #转换成numpy的mat
    labelMat = np.mat(classLabels).transpose() 
           #转换成numpy的mat,并进行转置
    m, n = np.shape(dataMatrix)
      #返回dataMatrix的大小。m为行数,n为列数。
    alpha = 0.01 
      #移动步长,也就是学习速率,控制更新的幅度。
    maxCycles = 500 #最大迭代次数
    weights = np.ones((n,1))
    for k in range(maxCycles):
        h = sigmoid(dataMatrix * weights) 
           #梯度上升矢量化公式
        error = labelMat - h
        weights = weights + alpha * 
            dataMatrix.transpose() * error
    return weights.getA(),weights_array 
           #将矩阵转换为数组,返回权重数组

假设,我们使用的数据集一共有100个样本。那么,dataMatrix就是一个100*3的矩阵。每次计算h的时候,都要计算dataMatrix*weights这个矩阵乘法运算,要进行100*3次乘法运算和100*2次加法运算。同理,更新回归系数(最优参数)weights时,也需要用到整个数据集,要进行矩阵乘法运算。总而言之,该方法处理100个左右的数据集时尚可,但如果有数十亿样本和成千上万的特征,那么该方法的计算复杂度就太高了。因此,需要对算法进行改进,我们每次更新回归系数(最优参数)的时候,能不能不用所有样本呢?一次只用一个样本点去更新回归系数(最优参数)?这样就可以有效减少计算量了,这种方法就叫做随机梯度上升算法。

1、随机梯度上升算法

让我们直接看代码:

def stocGradAscent1(dataMatrix, classLabels, numIter=150):
    m,n =np.shape(dataMatrix)
            #返回dataMatrix的大小,m为行数,n为列数。
    weights = np.ones(n)  #参数初始化
    for j in range(numIter): 
        dataIndex = list(range(m)
        for i in range(m): 
               #降低alpha的大小,每次减小1/(j+i)。
            randIndex = int(random.uniform(0,len(dataIndex))) 
                #随机选取样本
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
                #选择随机选取的一个样本,计算h
            error = classLabels[randIndex] - h 
                #计算误差
            weights = weights + alpha 
                   * error * dataMatrix[randIndex] 
               #更新回归系数
            del(dataIndex[randIndex]) 
                     #删除已经使用的样本
    return weights #返回

该算法第一个改进之处在于,alpha在每次迭代的时候都会调整,并且,虽然alpha会随着迭代次数不断减小,但永远不会减小到0,因为这里还存在一个常数项。必须这样做的原因是为了保证在多次迭代之后新数据仍然具有一定的影响。如果需要处理的问题是动态变化的,那么可以适当加大上述常数项,来确保新的值获得更大的回归系数。另一点值得注意的是,在降低alpha的函数中,alpha每次减少1/(j+i),其中j是迭代次数,i是样本点的下标。第二个改进的地方在于跟新回归系数(最优参数)时,只使用一个样本点,并且选择的样本点是随机的,每次迭代不使用已经用过的样本点。这样的方法,就有效地减少了计算量,并保证了回归效果。

编写代码如下,看下改进的随机梯度上升算法分类效果如何:

# -*- coding:UTF-8 -*-from matplotlib.font_manager import FontPropertiesimport matplotlib.pyplot as pltimport numpy as npimport random"""函数说明:加载数据Parameters:    无Returns:    dataMat - 数据列表    labelMat - 标签列表Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762 Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""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])]) #添加数据
        labelMat.append(int(lineArr[2])) #添加标签
    fr.close()        #关闭文件    return dataMat, labelMat   #返回"""函数说明:sigmoid函数Parameters:    inX - 数据Returns:    sigmoid函数Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""def sigmoid(inX):    return 1.0 / (1 + np.exp(-inX))"""函数说明:绘制数据集Parameters:    weights - 权重参数数组Returns:    无Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/ Modify:    2017-08-30"""def plotBestFit(weights):    dataMat, labelMat = loadDataSet() #加载数据集
    dataArr = np.array(dataMat)
                   #转换成numpy的array数组
    n = np.shape(dataMat)[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])
                                        #1为正样本
        else:
            xcord2.append(dataArr[i,1])
            ycord2.append(dataArr[i,2])
                                       #0为负样本
    fig = plt.figure()
    ax = fig.add_subplot(111)     #添加subplot
    ax.scatter(xcord1, ycord1, s = 20, c = 'red', 
           marker = 's',alpha=.5) #绘制正样本
    ax.scatter(xcord2, ycord2, s = 20,
           c = 'green',alpha=.5) #绘制负样本
    x = np.arange(-3.0, 3.0, 0.1)
    y = (-weights[0] - weights[1] * x) / weights[2]
    ax.plot(x, y)
    plt.title('BestFit')       #绘制title    plt.xlabel('X1');     plt.ylabel('X2')         #绘制label    plt.show()"""函数说明:改进的随机梯度上升算法Parameters:    dataMatrix - 数据数组    classLabels - 数据标签    numIter - 迭代次数Returns:    weights - 求得的回归系数数组(最优参数)Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-31"""def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n= np.shape(dataMatrix)
           #返回dataMatrix的大小m为行数,n为列数。
    weights = np.ones(n)       #参数初始化
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):            alpha = 4/(1.0+j+i)+0.01                #降低alpha的大小,每次减小1/(j+i)。            randIndex = int(random.uniform(0,len(dataIndex))) 
            #随机选取样本
            h = sigmoid(sum(dataMatrix[randIndex]*weights)) 
            #选择随机选取的一个样本,计算h
            error = classLabels[randIndex] - h #计算误差
            weights=weights+alpha*error*dataMatrix[randIndex] 
             #更新回归系数
            del(dataIndex[randIndex])#删除已经使用的样本    return weights    #返回if __name__ == '__main__':    dataMat, labelMat = loadDataSet()
    weights = stocGradAscent1(np.array(dataMat), labelMat)
    plotBestFit(weights)

代码运行结果:

2、回归系数与迭代次数的关系

可以看到分类效果也是不错的。不过,从这个分类结果中,我们不好看出迭代次数和回归系数的关系,也就不能直观的看到每个回归方法的收敛情况。因此,我们编写程序,绘制出回归系数和迭代次数的关系曲线:

# -*- coding:UTF-8 -*-from matplotlib.font_manager import FontPropertiesimport matplotlib.pyplot as pltimport numpy as npimport random"""函数说明:加载数据Parameters:    无Returns:    dataMat - 数据列表    labelMat - 标签列表Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""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])]) #添加数据
        labelMat.append(int(lineArr[2]))  #添加标签
    fr.close()           #关闭文件    return dataMat, labelMat       #返回"""函数说明:sigmoid函数Parameters:    inX - 数据Returns:    sigmoid函数Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""def sigmoid(inX):    return 1.0 / (1 + np.exp(-inX))"""函数说明:梯度上升算法Parameters:    dataMatIn - 数据集    classLabels - 数据标签Returns:    weights.getA() - 求得的权重数组(最优参数)    weights_array - 每次更新的回归系数Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-28"""def gradAscent(dataMatIn, classLabels):    dataMatrix = np.mat(dataMatIn)#转换成numpy的mat
    labelMat = np.mat(classLabels).transpose() 
               #转换成numpy的mat,并进行转置
    m, n = np.shape(dataMatrix)   
          #返回dataMatrix的大小。m为行数,n为列数。
    alpha = 0.01 
        #移动步长,也就是学习速率,控制更新的幅度。
    maxCycles = 500   #最大迭代次数
    weights = np.ones((n,1))
    weights_array = np.array([])
    for k in range(maxCycles):
        h = sigmoid(dataMatrix * weights)
                    #梯度上升矢量化公式
        error = labelMat - h
        weights = weights + alpha 
                   * dataMatrix.transpose() * error
        weights_array = np.append(weights_array,weights)
    weights_array = weights_array.reshape(maxCycles,n)    return weights.getA(),weights_array 
                 #将矩阵转换为数组,并返回"""函数说明:改进的随机梯度上升算法Parameters:    dataMatrix - 数据数组    classLabels - 数据标签    numIter - 迭代次数Returns:    weights - 求得的回归系数数组(最优参数)    weights_array - 每次更新的回归系数Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-31"""def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n= np.shape(dataMatrix)
       #返回dataMatrix的大小,m为行数,n为列数。
    weights = np.ones(n)    #参数初始化
    weights_array = np.array([]) 
               #存储每次更新的回归系数
    for j in range(numIter):                                      
        dataIndex = list(range(m))
        for i in range(m):           
            alpha = 4/(1.0+j+i)+0.01
             #降低alpha的大小,每次减小1/(j+i)。
            randIndex =int(random.uniform(0,len(dataIndex)))
                  #随机选取样本
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
               #选择随机选取的一个样本,计算h
            error = classLabels[randIndex] - h
                 #计算误差
            weights = weights + alpha * error 
           * dataMatrix[randIndex] #更新回归系数
            weights_array = np.append(weights_array,
             weights,axis=0) #添加回归系数到数组中
            del(dataIndex[randIndex]) 
                 #删除已经使用的样本
    weights_array = weights_array.reshape(numIter*m,n) 
                 #改变维度    return weights,weights_array      #返回"""函数说明:绘制回归系数与迭代次数的关系Parameters:    weights_array1 - 回归系数数组1    weights_array2 - 回归系数数组2Returns:    无Author:    Jack CuiBlog:    http://blog.csdn.net/c406495762Zhihu:    https://www.zhihu.com/people/Jack--Cui/Modify:    2017-08-30"""def plotWeights(weights_array1,weights_array2):    #设置汉字格式
    font =
    FontProperties(fname=r"c:\windows\fonts\simsun.ttc",size=14)
    #将fig画布分隔成1行1列,不共享x轴和y轴
    #fig画布的大小为(13,8)当nrow=3,nclos=2时
    #代表fig画布被分为个区域
    #axs[0][0]表示第一行第一列

    fig, axs = 
      plt.subplots(nrows=3, ncols=2,sharex=False,
                  sharey=False, figsize=(20,10))
    x1 = np.arange(0, len(weights_array1), 1)
    #绘制w0与迭代次数的关系
    axs[0][0].plot(x1,weights_array1[:,0])
    axs0_title_text = 
       axs[0][0].set_title(u'梯度上升算法:
         回归系数与迭代次数关系',FontProperties=font)
    axs0_ylabel_text = axs[0][0].
             set_ylabel(u'W0', FontProperties=font)
    plt.setp(axs0_title_text, 
                size=20,weight='bold',color='black') 
    plt.setp(axs0_ylabel_text,
                  size=20,weight='bold',color='black')
    #绘制w1与迭代次数的关系
    axs[1][0].plot(x1,weights_array1[:,1])
    axs1_ylabel_text = axs[1][0].
            set_ylabel(u'W1',FontProperties=font)
    plt.setp(axs1_ylabel_text,
                size=20,weight='bold',color='black')
                        #绘制w2与迭代次数的关系
    axs[2][0].plot(x1,weights_array1[:,2])
    axs2_xlabel_text = axs[2][0].
        set_xlabel(u'迭代次数',FontProperties=font)
    axs2_ylabel_text = axs[2][0].
         set_ylabel(u'W1', FontProperties=font)
    plt.setp(axs2_xlabel_text,
              size=20,weight='bold',color='black') 
    plt.setp(axs2_ylabel_text,
              size=20,weight='bold',color='black')


    x2 = np.arange(0, len(weights_array2), 1)
    #绘制w0与迭代次数的关系
    axs[0][1].plot(x2,weights_array2[:,0])
    axs0_title_text = 
       axs[0][1].set_title(u'改进的随机梯度上升算法:
         回归系数与迭代次数关系',FontProperties=font)
    axs0_ylabel_text = axs[0][1].
             set_ylabel(u'W0',FontProperties=font)
    plt.setp(axs0_title_text,
            size=20, weight='bold',color='black') 
    plt.setp(axs0_ylabel_text,
            size=20,weight='bold',color='black')
    #绘制w1与迭代次数的关系
    axs[1][1].plot(x2,weights_array2[:,1])
    axs1_ylabel_text = axs[1][1].
             set_ylabel(u'W1',FontProperties=font)
    plt.setp(axs1_ylabel_text,
            size=20,weight='bold',color='black')
    #绘制w2与迭代次数的关系
    axs[2][1].plot(x2,weights_array2[:,2])
    axs2_xlabel_text = axs[2][1].
       set_xlabel(u'迭代次数', FontProperties=font)
    axs2_ylabel_text = axs[2][1].
          set_ylabel(u'W1', FontProperties=font)
    plt.setp(axs2_xlabel_text,
            size=20,weight='bold',color='black') 
    plt.setp(axs2_ylabel_text,
            size=20,weight='bold',color='black')

    plt.show()

if __name__ == '__main__':
    dataMat, labelMat = loadDataSet()           
    weights1,weights_array1 =
       stocGradAscent1(np.array(dataMat), labelMat)
    weights2,weights_array2 = 
            gradAscent(dataMat, labelMat)
    plotWeights(weights_array1, weights_array2)

运行结果如下:

由于改进的随机梯度上升算法,随机选取样本点,所以每次的运行结果是不同的。但是大体趋势是一样的。我们改进的随机梯度上升算法收敛效果更好。为什么这么说呢?让我们分析一下。我们一共有100个样本点,改进的随机梯度上升算法迭代次数为150。而上图显示15000次迭代次数的原因是,使用一次样本就更新一下回归系数。因此,迭代150次,相当于更新回归系数150*100=15000次。简而言之,迭代150次,更新1.5万次回归参数。从上图左侧的改进随机梯度上升算法回归效果中可以看出,其实在更新2000次回归系数的时候,已经收敛了。相当于遍历整个数据集20次的时候,回归系数已收敛。训练已完成。

再让我们看看上图右侧的梯度上升算法回归效果,梯度上升算法每次更新回归系数都要遍历整个数据集。从图中可以看出,当迭代次数为300多次的时候,回归系数才收敛。凑个整,就当它在遍历整个数据集300次的时候已经收敛好了。

没有对比就没有伤害,改进的随机梯度上升算法,在遍历数据集的第20次开始收敛。而梯度上升算法,在遍历数据集的第300次才开始收敛。想像一下,大量数据的情况下,谁更牛逼?

本系列篇章:

Logistic回归实战篇之预测病马死亡率(一)

Logistic回归实战篇之预测病马死亡率(二)

Logistic回归实战篇之预测病马死亡率(三)

原文发布于微信公众号 - 机器学习算法全栈工程师(Jeemy110)

原文发表时间:2017-10-22

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏fangyangcoder

python机器学习实战(四)

http://www.cnblogs.com/fydeblog/p/7364317.html

982
来自专栏鸿的学习笔记

写给开发者的机器学习指南(十三)

在我们实际使用支持向量机(SVM)之前,我先简要介绍一下SVM是什么。 基本SVM是一个二元分类器,它通过选取代表数据点之间最大间隔的超平面将数据集分成2部分。...

421
来自专栏iOSDevLog

决策树模型概述

2005
来自专栏深度学习与计算机视觉

理解Logistic回归算法原理与Python实现

一般的机器学习的实现大致都是这样的步骤: 1.准备数据,包括数据的收集,整理等等 2.定义一个学习模型(learning function model)...

3368
来自专栏超然的博客

MIT-线性代数笔记(1-6)

  对方程组中某个方程进行时的那个的数乘和加减,将某一未知系数变为零,来削弱未知数个数

1002
来自专栏ml

调参过程中的参数 学习率,权重衰减,冲量(learning_rate , weight_decay , momentum)

无论是深度学习还是机器学习,大多情况下训练中都会遇到这几个参数,今天依据我自己的理解具体的总结一下,可能会存在错误,还请指正. learning_rate , ...

5138
来自专栏机器学习算法全栈工程师

Logistic回归实战篇之预测病马死亡率(一)

作 者:崔家华 编 辑:李文臣 Python版本: Python3.x 运行平台 : Windows IDE : Sublime...

2696
来自专栏书山有路勤为径

特征向量(Feature Vectors)

已经可以从图像中提取基于形状的特征,如何使用这一组特征来检测整个对象,以山峰图像角点检测举例:

1002
来自专栏开心的学习之路

贝叶斯决策理论(数学部分)

概率质量函数(Probability Mass Function)是针对离散值而言的,通常用大写字母P表示。假设某个事

1333
来自专栏机器学习算法全栈工程师

Logistic回归实战篇之预测病马死亡率(二)

作 者:崔家华 编 辑:李文臣 三、从疝气病症状预测病马的死亡率 1、实战背景 本次实战内容,将使用Logistic回归来预测患疝气病的马的存活问题。原始数据集...

3174

扫码关注云+社区