首页 » 机器学习实战 » 机器学习实战全文在线阅读

《机器学习实战》5.2 基于最优化方法的最佳回归系数确定

关灯直达底部

Sigmoid函数的输入记为z,由下面公式得出:

如果采用向量的写法,上述公式可以写成z = wTx。它表示将这两个数值向量的对应元素相乘然后全部加起来即得到z值。其中的向量x是分类器的输入数据,向量w也就是我们要找到的最佳参数(系数),从而使得分类尽可能地精确。为了寻找该最佳参数,需要用到最优化理论的一些知识。

下面首先介绍梯度上升的最优化方法,我们将学习到如何使用该方法求得数据集的最佳参数。接下来,展示如何绘制梯度上升法产生的决策边界图,该图能将梯度上升法的分类效果可视化地呈现出来。最后我们将学习随机梯度上升算法,以及如何对其进行修改以获得更好的结果。

5.2.1 梯度上升法

我们介绍的第一个最优化算法叫做梯度上升法。梯度上升法基于的思想是:要找到某函数的最大值,最好的方法是沿着该函数的梯度方向探寻。如果梯度记为∇,则函数f(x,y)的梯度由下式表示:

这是机器学习中最易造成混淆的一个地方,但在数学上并不难,需要做的只是牢记这些符号的意义。这个梯度意味着要沿x的方向移动,沿y的方向移动。其中,函数f(x,y)必须要在待计算的点上有定义并且可微。一个具体的函数例子见图5-2。

图5-2 梯度上升算法到达每个点后都会重新估计移动的方向。从P0开始,计算完该点的梯度,函数就根据梯度移动到下一点P1。在P1点,梯度再次被重新计算,并沿新的梯度方向移动到P2。如此循环迭代,直到满足停止条件。迭代的过程中,梯度算子总是保证我们能选取到最佳的移动方向

图5-2中的梯度上升算法沿梯度方向移动了一步。可以看到,梯度算子总是指向函数值增长最快的方向。这里所说的是移动方向,而未提到移动量的大小。该量值称为步长,记做。用向量来表示的话,梯度上升算法的迭代公式如下:

该公式将一直被迭代执行,直至达到某个停止条件为止,比如迭代次数达到某个指定值或算法达到某个可以允许的误差范围。

梯度下降算法

你最经常听到的应该是梯度下降算法,它与这里的梯度上升算法是一样的,只是公式中的加法需要变成减法。因此,对应的公式可以写成

梯度上升算法用来求函数的最大值,而梯度下降算法用来求函数的最小值。

基于上面的内容,我们来看一个Logistic回归分类器的应用例子,从图5-3可以看到我们采用的数据集。

图5-3 一个简单数据集,下面将采用梯度上升法找到Logistic回归分类器在此数据集上的最佳回归系数

5.2.2 训练算法:使用梯度上升找到最佳参数

图5-3中有100个样本点,每个点包含两个数值型特征:X1和X2。在此数据集上,我们将通过使用梯度上升法找到最佳回归系数,也就是拟合出Logistic回归模型的最佳参数。

梯度上升法的伪代码如下:

每个回归系数初始化为1重复R次:    计算整个数据集的梯度    使用`alpha × gradient`更新回归系数的向量    返回回归系数  

下面的代码是梯度上升算法的具体实现。为了解实际效果,打开文本编辑器并创建一个名为logRegres.py的文件,输入下列代码:

程序清单5-1 logistic回归梯度上升优化算法

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]))    return dataMat,labelMatdef sigmoid(inX):    return 1.0/(1+exp(-inX))def gradAscent(dataMatIn, classLabels):   # ❶(以下两行)转换为NumPy矩阵数据类型     dataMatrix = mat(dataMatIn)    labelMat = mat(classLabels).transpose    m,n = shape(dataMatrix)    alpha = 0.001    maxCycles = 500    weights = ones((n,1))    for k in range(maxCycles):        #❷(以下三行)矩阵相乘        h = sigmoid(dataMatrix*weights)        error = (labelMat - h)        weights = weights + alpha * dataMatrix.transpose* error    return weights 

程序清单5-1的代码在开头提供了一个便利函数loadDataSet,它的主要功能是打开文本文件testSet.txt并逐行读取。每行前两个值分别是X1和X2,第三个值是数据对应的类别标签。此外,为了方便计算,该函数还将X0的值设为1.0。接下来的函数是5.2节提到的函数sigmoid

梯度上升算法的实际工作是在函数gradAscent里完成的,该函数有两个参数。第一个参数是dataMatIn,它是一个2维NumPy数组,每列分别代表每个不同的特征,每行则代表每个训练样本。我们现在采用的是100个样本的简单数据集,它包含了两个特征X1和X2,再加上第0维特征X0,所以dataMath里存放的将是100×3的矩阵。在❶处,我们获得输入数据并将它们转换成NumPy矩阵。这是本书首次使用NumPy矩阵,如果你对矩阵数学不太熟悉,那么一些运算可能就会不易理解。比如,NumPy对2维数组和矩阵都提供一些操作支持,如果混淆了数据类型和对应的操作,执行结果将与预期截然不同。对此,本书附录A给出了对NumPy矩阵的介绍。第二个参数是类别标签,它是一个1×100的行向量。为了便于矩阵运算,需要将该行向量转换为列向量,做法是将原向量转置,再将它赋值给labelMat。接下来的代码是得到矩阵大小,再设置一些梯度上升算法所需的参数。

变量alpha是向目标移动的步长,maxCycles是迭代次数。在for循环迭代完成后,将返回训练好的回归系数。需要强调的是,在❷处的运算是矩阵运算。变量h不是一个数而是一个列向量,列向量的元素个数等于样本个数,这里是100。对应地,运算dataMatrix * weights代表的不是一次乘积计算,事实上该运算包含了300次的乘积。

最后还需说明一点,你可能对❷中公式的前两行觉得陌生。此处本书略去了一个简单的数学推导,我把它留给有兴趣的读者。定性地说,这里是在计算真实类别与预测类别的差值,接下来就是按照该差值的方向调整回归系数。

接下来看看实际效果,打开文本编辑器,添加程序清单5-1的代码。

在Python提示符下,敲入下面的代码:

>>> import logRegres>>> dataArr,labelMat=logRegres.loadDataSet>>> logRegres.gradAscent(dataArr,labelMat)matrix([[ 4.12414349],        [ 0.48007329],        [-0.6168482 ]])  

5.2.3 分析数据:画出决策边界

上面已经解出了一组回归系数,它确定了不同类别数据之间的分隔线。那么怎样能画出该分隔线,从而使得优化的过程便于理解呢?下面将解决这个问题,打开logRegres.py并添加如下代码。

程序清单5-2 画出数据集和logistic回归最佳拟合直线的函数

def plotBestFit(wei):    import matplotlib.pyplot as plt    weights = wei.getA    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  

程序清单5-2中的代码是直接用Matplotlib画出来的。唯一要指出的是,❶处设置了sigmoid函数为0。回忆5.2节,0是两个类别(类别1和类别0)的分界处。因此,我们设定 0 = w0x0 + w1x1 + w2x2,然后解出X2和X1的关系式(即分隔线的方程,注意X0=1)。

运行程序清单5-2的代码,在Python提示符下输入:

>>> reload(logRegres)<module /'logRegres/' from /'logRegres.py/'>>>> logRegres.plotBestFit(weights.getA)   

输出的结果如图5-4所示。

图5-4 梯度上升算法在500次迭代后得到的Logistic回归最佳拟合直线

这个分类结果相当不错,从图上看只错分了两到四个点。但是,尽管例子简单且数据集很小,这个方法却需要大量的计算(300次乘法)。因此下一节将对该算法稍作改进,从而使它可以用在真实数据集上。

5.2.4 训练算法:随机梯度上升

梯度上升算法在每次更新回归系数时都需要遍历整个数据集,该方法在处理100个左右的数据集时尚可,但如果有数十亿样本和成千上万的特征,那么该方法的计算复杂度就太高了。一种改进方法是一次仅用一个样本点来更新回归系数,该方法称为随机梯度上升算法。由于可以在新样本到来时对分类器进行增量式更新,因而随机梯度上升算法是一个在线学习算法。与“在线学习”相对应,一次处理所有数据被称作是“批处理”。

随机梯度上升算法可以写成如下的伪代码:

所有回归系数初始化为1  对数据集中每个样本        计算该样本的梯度          使用alpha × gradient更新回归系数值  返回回归系数值      

伪代码最后一行应顶格,原文错误,已修改

以下是随机梯度上升算法的实现代码。

程序清单5-3 随机梯度上升算法

def stocGradAscent0(dataMatrix, classLabels):    m,n = shape(dataMatrix)    alpha = 0.01    weights = ones(n)    for i in range(m):        h = sigmoid(sum(dataMatrix[i]*weights))        error = classLabels[i] - h        weights = weights + alpha * error * dataMatrix[i]    return weights    

可以看到,随机梯度上升算法与梯度上升算法在代码上很相似,但也有一些区别:第一,后者的变量h和误差error都是向量,而前者则全是数值;第二,前者没有矩阵的转换过程,所有变量的数据类型都是NumPy数组。

为了验证该方法的结果,我们将程序清单5-3的代码添加到logRegres.py中,并在Python提示符下输入如下命令:

>>> from numpy import *>>> reload(logRegres)<module /'logRegres/' from /'logRegres.py/'>>>> dataArr,labelMat=logRegres.loadDataSet>>> weights=logRegres.stocGradAscent0(array(dataArr),labelMat)>>> logRegres.plotBestFit(weights)    

执行完毕后将得到图5-5所示的最佳拟合直线图,该图与图5-4有一些相似之处。可以看到,拟合出来的直线效果还不错,但并不像图5-4那样完美。这里的分类器错分了三分之一的样本。

直接比较程序清单5-3和程序清单5-1的代码结果是不公平的,后者的结果是在整个数据集上迭代了500次才得到的。一个判断优化算法优劣的可靠方法是看它是否收敛,也就是说参数是否达到了稳定值,是否还会不断地变化?对此,我们在程序清单5-3中随机梯度上升算法上做了些修改,使其在整个数据集上运行200次。最终绘制的三个回归系数的变化情况如图5-6所示。

图5-5 随机梯度上升算法在上述数据集上的执行结果,最佳拟合直线并非最佳分类线

图5-6 运行随机梯度上升算法,在数据集的一次遍历中回归系数与迭代次数的关系图。回归系数经过大量迭代才能达到稳定值,并且仍然有局部的波动现象

图5-6展示了随机梯度上升算法在200次迭代过程中回归系数的变化情况。其中的系数2,也就是图5-5中的X2只经过了50次迭代就达到了稳定值,但系数1和0则需要更多次的迭代。另外值得注意的是,在大的波动停止后,还有一些小的周期性波动。不难理解,产生这种现象的原因是存在一些不能正确分类的样本点(数据集并非线性可分),在每次迭代时会引发系数的剧烈改变。我们期望算法能避免来回波动,从而收敛到某个值。另外,收敛速度也需要加快。

对于图5-6存在的问题,可以通过修改程序清单5-3的随机梯度上升算法来解决,具体代码如下:

程序清单5-4  改进的随机梯度上升算法

def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n = shape(dataMatrix)    weights = ones(n)    for j in range(numIter):          dataIndex = range(m)        for i in range(m):            #❶  alpha每次迭代时需要调整            alpha = 4/(1.0+j+i)+0.01            #❷  随机选取更新             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  

程序清单5-4与5-3类似,但增加了两处代码来进行改进。

第一处改进在❶处。一方面,alpha在每次迭代的时候都会调整,这会缓解图5-6上的数据波动或者高频波动。另外,虽然alpha会随着迭代次数不断减小,但永远不会减小到0,这是因为❶中还存在一个常数项。必须这样做的原因是为了保证在多次迭代之后新数据仍然具有一定的影响。如果要处理的问题是动态变化的,那么可以适当加大上述常数项,来确保新的值获得更大的回归系数。另一点值得注意的是,在降低alpha的函数中,alpha每次减少 1/(j+i) ,其中j是迭代次数,i是样本点的下标1 。这样当j<<max(i)时,alpha就不是严格下降的。避免参数的严格下降也常见于模拟退火算法等其他优化算法中。

1. 要注意区分这里的下标与样本编号,编号表示了样本在矩阵中的位置(代码中为randIndex),而这里的下标i表示本次迭代中第i个选出来的样本。——译者注

程序清单5-4第二个改进的地方在❷处,这里通过随机选取样本来更新回归系数。这种方法将减少周期性的波动(如图5-6中的波动)。具体实现方法与第3章类似,这种方法每次随机从列表中选出一个值,然后从列表中删掉该值(再进行下一次迭代)。

此外,改进算法还增加了一个迭代次数作为第3个参数。如果该参数没有给定的话,算法将默认迭代150次。如果给定,那么算法将按照新的参数值进行迭代。

stocGradAscent1类似,图5-7显示了每次迭代时各个回归系数的变化情况。

图5-7 使用样本随机选择和alpha动态减少机制的随机梯度上升算法stocGradAscent1所生成的系数收敛示意图。该方法比采用固定alpha的方法收敛速度更快

比较图5-7和图5-6可以看到两点不同。第一点是,图5-7中的系数没有像图5-6里那样出现周期性的波动,这归功于stocGradAscent1里的样本随机选择机制;第二点是,图5-7的水平轴比图5-6短了很多,这是由于stocGradAscent1可以收敛得更快。这次我们仅仅对数据集做了20次遍历,而之前的方法是500次。

下面看看在同一个数据集上的分类效果。将程序清单5-4的代码添加到logRegres.py文件中,并在Python提示符下输入:

>>> reload(logRegres)<module /'logRegres/' from /'logRegres.py/'>>>> dataArr,labelMat=logRegres.loadDataSet>>> weights=logRegres.stocGradAscent1(array(dataArr),labelMat)>>> logRegres.plotBestFit(weights)  

程序运行结束之后应该可以看到与图5-8类似的结果图。该分隔线达到了与GradientAscent差不多的效果,但是所使用的计算量更少。

图5-8 使用改进的随机梯度上升算法得到的系数

默认迭代次数是150,可以通过stocGradAscent的第3个参数来对此进行修改,例如:

>>> weights=logRegres.stocGradAscent1(array(dataArr),labelMat, 500)     

目前,我们已经学习了几个优化算法,但还有很多优化算法值得探讨,所幸这方面已有大量的文献可供参考。另外再说明一下,针对给定的数据集,读者完全可以对算法的各种参数进行调整,从而达到更好的效果。

迄今为止我们分析了回归系数的变化情况,但还没有达到本章的最终目标,即完成具体的分类任务。下一节将使用随机梯度上升算法来解决病马的生死预测问题。