机器学习实战教程(六):Logistic回归实战篇之预测病马死亡率(a

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

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

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、随机梯度上升算法

让我们直接看代码:

Python

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[dataIndex[randIndex]]*weights))  #选择随机选取的一个样本,计算h

            error = classLabels[dataIndex[randIndex]] - h                          #计算误差

            weights = weights + alpha * error * dataMatrix[dataIndex[randIndex]]  #更新回归系数

            del(dataIndex[randIndex])                                        #删除已经使用的样本

    return weights                                                      #返回

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[dataIndex[randIndex]]*weights))  #选择随机选取的一个样本,计算h

            error = classLabels[dataIndex[randIndex]] - h                          #计算误差

            weights = weights + alpha * error * dataMatrix[dataIndex[randIndex]]  #更新回归系数

            del(dataIndex[randIndex])                                        #删除已经使用的样本

    return weights                                                      #返回

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

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

# -*- coding:UTF-8 -*-

from matplotlib.font_manager import FontProperties

import matplotlib.pyplot as plt

import numpy as np

import random

"""

函数说明:加载数据

Parameters:

    无

Returns:

    dataMat - 数据列表

    labelMat - 标签列表

"""

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函数

"""

def sigmoid(inX):

    return 1.0 / (1 + np.exp(-inX))

"""

函数说明:绘制数据集

Parameters:

    weights - 权重参数数组

"""

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 - 求得的回归系数数组(最优参数)

"""

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[dataIndex[randIndex]]*weights))        #选择随机选取的一个样本,计算h

            error = classLabels[dataIndex[randIndex]] - h                        #计算误差

            weights = weights + alpha * error * dataMatrix[dataIndex[randIndex]]  #更新回归系数

            del(dataIndex[randIndex])                                        #删除已经使用的样本

    return weights                                                            #返回

if __name__ == '__main__':

    dataMat, labelMat = loadDataSet()

    weights = stocGradAscent1(np.array(dataMat), labelMat)

    plotBestFit(weights)

代码运行结果:

免费视频教程:www.mlxs.top      


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

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

# -*- coding:UTF-8 -*-

from matplotlib.font_manager import FontProperties

import matplotlib.pyplot as plt

import numpy as np

import random

"""

函数说明:加载数据

Parameters:

    无

Returns:

    dataMat - 数据列表

    labelMat - 标签列表

"""

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函数

"""

def sigmoid(inX):

    return 1.0 / (1 + np.exp(-inX))

"""

函数说明:梯度上升算法

Parameters:

    dataMatIn - 数据集

    classLabels - 数据标签

Returns:

    weights.getA() - 求得的权重数组(最优参数)

    weights_array - 每次更新的回归系数

"""

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 - 每次更新的回归系数

"""

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[dataIndex[randIndex]]*weights))          #选择随机选取的一个样本,计算h

            error = classLabels[dataIndex[randIndex]] - h                          #计算误差

            weights = weights + alpha * error * dataMatrix[dataIndex[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 - 回归系数数组2

"""

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'W2',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)

运行结果如下:

免费视频教程:www.mlxs.top      

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

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

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

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 204,293评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,604评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 150,958评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,729评论 1 277
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,719评论 5 366
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,630评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,000评论 3 397
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,665评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,909评论 1 299
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,646评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,726评论 1 330
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,400评论 4 321
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,986评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,959评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,197评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 44,996评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,481评论 2 342