python梯度下降
‘壹’ 一文搞懂梯度下降&反向传播
如果把神经网络模型比作一个黑箱,把模型参数比作黑箱上面一个个小旋钮,那么根据通用近似理论(universal approximation theorem),只要黑箱上的旋钮数量足够多,而且每个旋钮都被调节到合适的位置,那这个模型就可以实现近乎任意功能(可以逼近任意的数学模型)。
显然,这些旋钮(参数)不是由人工调节的,所谓的机器学习,就是通过程序来自动调节这些参数。神经网络不仅参数众多(少则十几万,多则上亿),而且网络是由线性层和非线性层交替叠加而成,上层参数的变化会对下层的输出产生非线性的影响,因此,早期的神经网络流派一度无法往多层方向发展,因为他们找不到能用于任意多层网络的、简洁的自动调节参数的方法。
直到上世纪80年代,祖师爷辛顿发明了反向传播算法,用输出误差的均方差(就是loss值)一层一层递进地反馈到各层神经网络,用梯度下降法来调节每层网络的参数。至此,神经网络才得以开始它的深度之旅。
本文用python自己动手实现梯度下降和反向传播算法。 请点击这里 到Github上查看源码。
梯度下降法是一种将输出误差反馈到神经网络并自动调节参数的方法,它通过计算输出误差的loss值( J )对参数 W 的导数,并沿着导数的反方向来调节 W ,经过多次这样的操作,就能将输出误差减小到最小值,即曲线的最低点。
虽然Tensorflow、Pytorch这些框架都实现了自动求导的功能,但为了彻底理解参数调节的过程,还是有必要自己动手实现梯度下降和反向传播算法。我相信你和我一样,已经忘了之前学的微积分知识,因此,到可汗学院复习下 Calculus
和 Multivariable Calculus 是个不错的方法,或是拜读 这篇关于神经网络矩阵微积分的文章 。
Figure2是求导的基本公式,其中最重要的是 Chain Rule ,它通过引入中间变量,将“ y 对 x 求导”的过程转换为“ y 对中间变量 u 求导,再乘以 u 对 x 求导”,这样就将一个复杂的函数链求导简化为多个简单函数求导。
如果你不想涉及这些求导的细节,可以跳过具体的计算,领会其思想就好。
对于神经网络模型: Linear -> ReLu -> Linear -> MSE(Loss function) 来说,反向传播就是根据链式法则对 求导,用输出误差的均方差(MSE)对模型的输出求导,并将导数传回上一层神经网络,用于它们来对 w 、 b 和 x (上上层的输出)求导,再将 x 的导数传回到它的上一层神经网络,由此将输出误差的均方差通过递进的方式反馈到各神经网络层。
对于 求导的第一步是为这个函数链引入中间变量:
接着第二步是对各中间变量求导,最后才是将这些导数乘起来。
首先,反向传播的起点是对loss function求导,即 。 :
mse_grad()之所以用unsqueeze(-1)给导数增加一个维度,是为了让导数的shape和tensor shape保持一致。
linear层的反向传播是对 求导,它也是一个函数链,也要先对中间变量求导再将所有导数相乘:
这些中间变量的导数分别是:
对向量 求导,指的是对向量所有的标量求偏导( ),即: ,这个横向量也称为y的梯度。
这里 ,是一个向量,因此, 求导,指的是y的所有标量(y_1, y_2, ..., y_n)对向量x求偏导,即:
。
这个矩阵称为雅克比矩阵,它是个对角矩阵,因为 ,因此 。
同理, 。
因此,所有中间导数相乘的结果:
lin_grad() 中的inp.g、w.g和b.g分别是求 的导数,以inp.g为例,它等于 ,且需要乘以前面各层的导数,即 outp.g @ w.t() ,之所以要用点积运算符(@)而不是标量相乘,是为了让它的导数shape和tensor shape保持一致。同理,w.g和b.g也是根据相同逻辑来计算的。
ReLu层的求导相对来说就简单多了,当输入 <= 0时,导数为0,当输入 > 0时,导数为1。
求导运算终于结束了,接下来就是验证我们的反向传播是否正确。验证方法是将forward_backward()计算的导数和Pytorch自动微分得到的导数相比较,如果它们相近,就认为我们的反向传播算法是正确的。
首先,将计算好的参数导数保存到w1g、b1g、w2g和b2g中,再用Pytorch的自动微分来求w11、b11、w22和b22的导数。
最后,用np.allclose()来比较导数间的差异,如果有任何一个导数不相近,assert就会报错。结果证明,我们自己动手实现的算法是正确的。
反向传播是遵循链式法则的,它将前向传播的输出作为输入,输入作为输出,通过递进的方式将求导这个动作从后向前传递回各层。神经网络参数的求导需要进行矩阵微积分计算,根据这些导数的反方向来调节参数,就可以让模型的输出误差的优化到最小值。
欢迎关注和点赞,你的鼓励将是我创作的动力
‘贰’ 梯度下降算法的原理是什么
梯度下降是非常常用的优化算法。作为机器学习的基础知识,这是一个必须要掌握的算法。借助本文,让我们来一起详细了解一下这个算法。
前言
本文的代码可以到我的Github上获取:
https://github.com/paulQuei/gradient_descent
本文的算法示例通过Python语言实现,在实现中使用到了numpy和matplotlib。如果你不熟悉这两个工具,请自行在网上搜索教程。
关于优化
大多数学习算法都涉及某种形式的优化。优化指的是改变x以最小化或者最大化某个函数的任务。
我们通常以最小化指代大多数最优化问题。最大化可经由最小化来实现。
我们把要最小化或最大化的函数成为目标函数(objective function)或准则(criterion)。
我们通常使用一个上标*表示最小化或最大化函数的x值,记做这样:
[x^* = arg; min; f(x)]
优化本身是一个非常大的话题。如果有兴趣,可以通过《数值优化》和《运筹学》的书籍进行学习。
模型与假设函数
所有的模型都是错误的,但其中有些是有用的。– George Edward Pelham Box
模型是我们对要分析的数据的一种假设,它是为解决某个具体问题从数据中学习到的,因此它是机器学习最核心的概念。
针对一个问题,通常有大量的模型可以选择。
本文不会深入讨论这方面的内容,关于各种模型请参阅机器学习的相关书籍。本文仅以最简单的线性模型为基础来讨论梯度下降算法。
这里我们先介绍一下在监督学习(supervised learning)中常见的三个符号:
m,描述训练样本的数量
x,描述输入变量或特征
y,描述输出变量或者叫目标值
- 请注意,一个样本可能有很多的特征,因此x和y通常是一个向量。不过在刚开始学习的时候,为了便于理解,你可以暂时理解为这就是一个具体的数值。
- 代价函数也叫损失函数。
- 不同的模型可能会用不同的损失函数。例如,logistic回归的假设函数是这样的:。其代价函数是这样的:
对于一个函数,怎么确定下行的方向?
每一步该往前走多远?
有没有可能停留在半山腰的平台上?
- 这里的下标i表示第i个参数。 上标k指的是第k步的计算结果,而非k次方。在能够理解的基础上,下文的公式中将省略上标k。
收敛是指函数的变化率很小。具体选择多少合适需要根据具体的项目来确定。在演示项目中我们可以选择0.01或者0.001这样的值。不同的值将影响算法的迭代次数,因为在梯度下降的最后,我们会越来越接近平坦的地方,这个时候函数的变化率也越来越小。如果选择一个很小的值,将可能导致算法迭代次数暴增。
公式中的 称作步长,也称作学习率(learning rate)。它决定了每一步往前走多远,关于这个值我们会在下文中详细讲解。你可以暂时人为它是一个类似0.01或0.001的固定值。
在具体的项目,我们不会让算法无休止的运行下去,所以通常会设置一个迭代次数的最大上限。
我们随机选择了 都为10作为起点
设置最多迭代1000次
收敛的范围设为0.001
学习步长设为0.01
如果样本数量较小(例如小于等于2000),选择BGD即可。
如果样本数量很大,选择 来进行MBGD,例如:64,128,256,512。
- 《深度学习》一书中是这样描述的:“与其说是科学,这更像是一门艺术,我们应该谨慎地参考关于这个问题的大部分指导。”。
对于凸函数或者凹函数来说,不存在局部极值的问题。其局部极值一定是全局极值。
最近的一些研究表明,某些局部极值并没有想象中的那么糟糕,它们已经非常的接近全局极值所带来的结果了。
Wikipeida: Gradient descent
Sebastian Ruder: An overview of gradient descent optimization algorithms
吴恩达:机器学习
吴恩达:深度学习
Peter Flach:机器学习
李宏毅 - ML Lecture 3-1: Gradient Descent
PDF: 李宏毅 - Gradient Descent
Intro to optimization in deep learning: Gradient Descent
Intro to optimization in deep learning: Momentum, RMSProp and Adam
Stochastic Gradient Descent – Mini-batch and more
刘建平Pinard - 梯度下降(Gradient Descent)小结
多元函数的偏导数、方向导数、梯度以及微分之间的关系思考
[Machine Learning] 梯度下降法的三种形式BGD、SGD以及MBGD
- 作者:阿Paul https://paul.pub/gradient-descent/
训练集会包含很多的样本,我们用 表示其中第i个样本。
x是数据样本的特征,y是其目标值。例如,在预测房价的模型中,x是房子的各种信息,例如:面积,楼层,位置等等,y是房子的价格。在图像识别的任务中,x是图形的所有像素点数据,y是图像中包含的目标对象。
我们是希望寻找一个函数,将x映射到y,这个函数要足够的好,以至于能够预测对应的y。由于历史原因,这个函数叫做假设函数(hypothesis function)。
学习的过程如下图所示。即:首先根据已有的数据(称之为训练集)训练我们的算法模型,然后根据模型的假设函数来进行新数据的预测。
线性模型(linear model)正如其名称那样:是希望通过一个直线的形式来描述模式。线性模型的假设函数如下所示:
[h_{ heta}(x) = heta_{0} + heta_{1} * x]
这个公式对于大家来说应该都是非常简单的。如果把它绘制出来,其实就是一条直线。
下图是一个具体的例子,即: 的图形:
在实际的机器学习工程中,你会拥有大量的数据。这些数据会来自于某个数据源。它们存储在csv文件中,或者以其他的形式打包。
但是本文作为演示使用,我们通过一些简单的代码自动生成了需要的数据。为了便于计算,演示的数据量也很小。
import numpy as np
max_x = 10
data_size = 10
theta_0 = 5
theta_1 = 2
def get_data:
x = np.linspace(1, max_x, data_size)
noise = np.random.normal(0, 0.2, len(x))
y = theta_0 + theta_1 * x + noise
return x, y
这段代码很简单,我们生成了x范围是 [1, 10] 整数的10条数据。对应的y是以线性模型的形式计算得到,其函数是:。现实中的数据常常受到各种因素的干扰,所以对于y我们故意加上了一些高斯噪声。因此最终的y值为比原先会有轻微的偏离。
最后我们的数据如下所示:
x = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
y = [6.66, 9.11, 11.08, 12.67, 15.12, 16.76, 18.75, 21.35, 22.77, 24.56]
我们可以把这10条数据绘制出来这样就有一个直观的了解了,如下图所示:
虽然演示用的数据是我们通过公式计算得到的。但在实际的工程中,模型的参数是需要我们通过数据学习到的。所以下文我们假设我们不知道这里线性模式的两个参数是什么,而是通过算法的形式求得。
最后再跟已知的参数进行对比以验证我们的算法是否正确。
有了上面的数据,我们可以尝试画一条直线来描述我们的模型。
例如,像下面这样画一条水平的直线:
很显然,这条水平线离数据太远了,非常的不匹配。
那我们可以再画一条斜线。
我们初次画的斜线可能也不贴切,它可能像下面这样:
最后我们通过不断尝试,找到了最终最合适的那条,如下所示:
梯度下降算法的计算过程,就和这种本能式的试探是类似的,它就是不停的迭代,一步步的接近最终的结果。
代价函数
上面我们尝试了几次通过一条直线来拟合(fitting)已有的数据。
二维平面上的一条直线可以通过两个参数唯一的确定,两个参数的确定也即模型的确定。那如何描述模型与数据的拟合程度呢?答案就是代价函数。
代价函数(cost function)描述了学习到的模型与实际结果的偏差程度。以上面的三幅图为例,最后一幅图中的红线相比第一条水平的绿线,其偏离程度(代价)应该是更小的。
很显然,我们希望我们的假设函数与数据尽可能的贴近,也就是说:希望代价函数的结果尽可能的小。这就涉及到结果的优化,而梯度下降就是寻找最小值的方法之一。
对于每一个样本,假设函数会依据计算出一个估算值,我们常常用来表示。即 。
很自然的,我们会想到,通过下面这个公式来描述我们的模型与实际值的偏差程度:
[(h_ heta(x^i) - y^i)^2 = (widehat{y}^{i} - y^i)^2 = ( heta_{0} + heta_{1} * x^{i} - y^{i})^2]
请注意, 是实际数据的值, 是我们的模型的估算值。前者对应了上图中的离散点的y坐标,后者对应了离散点在直线上投影点的y坐标。
每一条数据都会存在一个偏差值,而代价函数就是对所有样本的偏差求平均值,其计算公式如下所示:
[L( heta) = frac {1}{m} sum_{i=1}^{m}(h_ heta(x^i) - y^i)^2 = frac {1}{m} sum_{i=1}^{m}( heta_{0} + heta_{1} * x^{i} - y^{i})^2]
当损失函数的结果越小,则意味着通过我们的假设函数估算出的结果与真实值越接近。这也就是为什么我们要最小化损失函数的原因。
借助上面这个公式,我们可以写一个函数来实现代价函数:
def cost_function(x, y, t0, t1):
cost_sum = 0
for i in range(len(x)):
cost_item = np.power(t0 + t1 * x[i] - y[i], 2)
cost_sum += cost_item
return cost_sum / len(x)
这个函数的代码应该不用多做解释,它就是根据上面的完成计算。
我们可以尝试选取不同的 和 组合来计算代价函数的值,然后将结果绘制出来:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D
theta_0 = 5
theta_1 = 2
def draw_cost(x, y):
fig = plt.figure(figsize=(10, 8))
ax = fig.gca(projection='3d')
scatter_count = 100
radius = 1
t0_range = np.linspace(theta_0 - radius, theta_0 + radius, scatter_count)
t1_range = np.linspace(theta_1 - radius, theta_1 + radius, scatter_count)
cost = np.zeros((len(t0_range), len(t1_range)))
for a in range(len(t0_range)):
for b in range(len(t1_range)):
cost[a][b] = cost_function(x, y, t0_range[a], t1_range[b])
t0, t1 = np.meshgrid(t0_range, t1_range)
ax.set_xlabel('theta_0')
ax.set_ylabel('theta_1')
ax.plot_surface(t0, t1, cost, cmap=cm.hsv)
在这段代码中,我们对 和 各自指定了一个范围进行100次的采样,然后以不同的 组合对来计算代价函数的值。
如果我们将所有点的代价函数值绘制出来,其结果如下图所示:
从这个图形中我们可以看出,当 越接近 [5, 2]时其结果(偏差)越小。相反,离得越远,结果越大。
直观解释
从上面这幅图中我们可以看出,代价函数在不同的位置结果大小不同。
从三维的角度来看,这就和地面的高低起伏一样。最高的地方就好像是山顶。
而我们的目标就是:从任意一点作为起点,能够快速寻找到一条路径并以此到达图形最低点(代价值最小)的位置。
而梯度下降的算法过程就和我们从山顶想要快速下山的做法是一样的。
在生活中,我们很自然会想到沿着最陡峭的路往下行是下山速度最快的。如下面这幅图所示:
针对这幅图,细心的读者可能很快就会有很多的疑问,例如:
这些问题也就是本文接下来要讨论的内容。
算法描述
梯度下降算法最开始的一点就是需要确定下降的方向,即:梯度。
我们常常用 来表示梯度。
对于一个二维空间的曲线来说,梯度就是其切线的方向。如下图所示:
而对于更高维空间的函数来说,梯度由所有变量的偏导数决定。
其表达式如下所示:
[ abla f({ heta}) = ( frac{partial f({ heta})}{partial heta_1} , frac{partial f({ heta})}{partial heta_2} , ... , frac{partial f({ heta})}{partial heta_n} )]
在机器学习中,我们主要是用梯度下降算法来最小化代价函数,记做:
[ heta ^* = arg min L( heta)]
其中,L是代价函数,是参数。
梯度下降算法的主体逻辑很简单,就是沿着梯度的方向一直下降,直到参数收敛为止。
记做:
[ heta ^{k + 1}_i = heta^{k}_i - lambda abla f( heta^{k})]
这里有几点需要说明:
线性回归的梯度下降
有了上面的知识,我们可以回到线性模型代价函数的梯度下降算法实现了。
首先,根据代价函数我们可以得到梯度向量如下:
[ abla f({ heta}) = (frac{partial L( heta)}{ partial heta_{0}}, frac{ partial L( heta)}{ partial heta_{1}}) = (frac {2}{m} sum_{i=1}^{m}( heta_{0} + heta_{1} * x^{i} - y^{i}) , frac {2}{m} sum_{i=1}^{m}( heta_{0} + heta_{1} * x^{i} - y^{i}) x^{i})]
接着,将每个偏导数带入迭代的公式中,得到:
[ heta_{0} := heta_{0} - lambda frac{partial L( heta_{0})}{ partial heta_{0}} = heta_{0} - frac {2 lambda }{m} sum_{i=1}^{m}( heta_{0} + heta_{1} * x^{i} - y^{i}) heta_{1} := heta_{1} - lambda frac{partial L( heta_{1})}{ partial heta_{1}} = heta_{1} - frac {2 lambda }{m} sum_{i=1}^{m}( heta_{0} + heta_{1} * x^{i} - y^{i}) x^{i}]
由此就可以通过代码实现我们的梯度下降算法了,算法逻辑并不复杂:
learning_rate = 0.01
def gradient_descent(x, y):
t0 = 10
t1 = 10
delta = 0.001
for times in range(1000):
sum1 = 0
sum2 = 0
for i in range(len(x)):
sum1 += (t0 + t1 * x[i] - y[i])
sum2 += (t0 + t1 * x[i] - y[i]) * x[i]
t0_ = t0 - 2 * learning_rate * sum1 / len(x)
t1_ = t1 - 2 * learning_rate * sum2 / len(x)
print('Times: {}, gradient: [{}, {}]'.format(times, t0_, t1_))
if (abs(t0 - t0_) < delta and abs(t1 - t1_) < delta):
print('Gradient descent finish')
return t0_, t1_
t0 = t0_
t1 = t1_
print('Gradient descent too many times')
return t0, t1
这段代码说明如下:
如果我们将算法迭代过程中求得的线性模式绘制出来,可以得到下面这幅动态图:
最后算法得到的结果如下:
Times: 657, gradient: [5.196562662718697, 1.952931052920264]
Times: 658, gradient: [5.195558390180733, 1.9530753071808193]
Times: 659, gradient: [5.194558335124868, 1.9532189556399233]
Times: 660, gradient: [5.193562479839619, 1.9533620008416623]
Gradient descent finish
从输出中可以看出,算法迭代了660次就收敛了。这时的结果[5.193562479839619, 1.9533620008416623],这已经比较接近目标值 [5, 2]了。如果需要更高的精度,可以将delta的值调的更小,当然,此时会需要更多的迭代次数。
高维扩展
虽然我们举的例子是二维的,但是对于更高维的情况也是类似的。同样是根据迭代的公式进行运算即可:
[ heta_{i} = heta_{i} - lambda frac {partial L( heta)}{partial heta_i} = heta_{i} - frac{2lambda}{m} sum_{i=1}^{m}(h_ heta(x^{k})-y^k)x_i^k]
这里的下标i表示第i个参数,上标k表示第k个数据。
梯度下降家族BGD
在上面的内容中我们看到,算法的每一次迭代都需要把所有样本进行遍历处理。这种做法称为之Batch Gradient Descent,简称BGD。作为演示示例只有10条数据,这是没有问题的。
但在实际的项目中,数据集的数量可能是几百万几千万条,这时候每一步迭代的计算量就会非常的大了。
于是就有了下面两个变种。
SGD
Stochastic Gradient Descent,简称SGD,这种算法是每次从样本集中仅仅选择一个样本来进行计算。很显然,这样做算法在每一步的计算量一下就少了很多。
其算法公式如下:
[ heta_{i} = heta_{i} - lambda frac {partial L( heta)}{partial heta_i} = heta_{i} - lambda(h_ heta(x^k)-y^k)x_i^k]
当然,减少算法计算量也是有代价的,那就是:算法结果会强依赖于随机取到的数据情况,这可能会导致算法的最终结果不太令人满意。
MBGD
以上两种做法其实是两个极端,一个是每次用到了所有数据,另一个是每次只用一个数据。
我们自然就会想到两者取其中的方法:每次选择一小部分数据进行迭代。这样既避免了数据集过大导致每次迭代计算量过大的问题,也避免了单个数据对算法的影响。
这种算法称之为Mini-batch Gradient Descent,简称MBGD。
其算法公式如下:
[ heta_{i} = heta_{i} - lambda frac {partial L( heta)}{partial heta_i} = heta_{i} - frac{2lambda}{m} sum_{i=a}^{a + b}(h_ heta(x^k)-y^k)x_i^k]
当然,我们可以认为SGD是Mini-batch为1的特例。
针对上面提到的算法变种,该如何选择呢?
下面是Andrew Ng给出的建议:
下表是 Optimization for Deep Learning 中对三种算法的对比
方法准确性更新速度内存占用在线学习BGD好慢高否SGD好(with annealing)快低是MBGD好中等中等是
算法优化
式7是算法的基本形式,在这个基础上有很多人进行了更多的研究。接下来我们介绍几种梯度下降算法的优化方法。
Momentum
Momentum是动量的意思。这个算法的思想就是借助了动力学的模型:每次算法的迭代会使用到上一次的速度作为依据。
算法的公式如下:
[v^t = gamma v^{t - 1} + lambda abla f( heta) heta = heta - v_t]
对比式7可以看出,这个算法的主要区别就是引入了,并且,每个时刻的受前一个时刻的影响。
从形式上看,动量算法引入了变量 v 充当速度角色——它代表参数在参数空间移动的方向和速率。速度被设为负梯度的指数衰减平均。名称动量来自物理类比,根据牛顿运动定律,负梯度是移动参数空间中粒子的力。动量在物理学上定义为质量乘以速度。在动量学习算法中,我们假设是单位质量,因此速度向量 v 也可以看作是粒子的动量。
对于可以取值0,而是一个常量,设为0.9是一个比较好的选择。
下图是momentum算法的效果对比:
对原来的算法稍加修改就可以增加动量效果:
def gradient_descent_with_momentum(x, y):
t0 = 10
t1 = 10
delta = 0.001
v0 = 0
v1 = 0
gamma = 0.9
for times in range(1000):
sum1 = 0
sum2 = 0
for i in range(len(x)):
sum1 += (t0 + t1 * x[i] - y[i])
sum2 += (t0 + t1 * x[i] - y[i]) * x[i]
v0 = gamma * v0 + 2 * learning_rate * sum1 / len(x)
v1 = gamma * v1 + 2 * learning_rate * sum2 / len(x)
t0_ = t0 - v0
t1_ = t1 - v1
print('Times: {}, gradient: [{}, {}]'.format(times, t0_, t1_))
if (abs(t0 - t0_) < delta and abs(t1 - t1_) < delta):
print('Gradient descent finish')
return t0_, t1_
t0 = t0_
t1 = t1_
print('Gradient descent too many times')
return t0, t1
以下是该算法的输出:
Times: 125, gradient: [4.955453758569991, 2.000005017897775]
Times: 126, gradient: [4.955309381126545, 1.9956928964532015]
Times: 127, gradient: [4.9542964317327005, 1.9855674828684156]
Times: 128, gradient: [4.9536358220657, 1.9781180992510465]
Times: 129, gradient: [4.95412496254411, 1.9788858350530971]
Gradient descent finish
从结果可以看出,改进的算法只用了129次迭代就收敛了。速度比原来660次快了很多。
同样的,我们可以把算法计算的过程做成动态图:
对比原始的算法过程可以看出,改进算法最大的区别是:在寻找目标值时会在最终结果上下跳动,但是越往后跳动的幅度越小,这也就是动量所产生的效果。
Learning Rate 优化
至此,你可能还是好奇该如何设定学习率的值。
事实上,这个值的选取需要一定的经验或者反复尝试才能确定。
关键在于,这个值的选取不能过大也不能过小。
如果这个值过小,会导致每一次迭代的步长很小,其结果就是算法需要迭代非常多的次数。
那么,如果这个值过大会怎么样呢?其结果就是:算法可能在结果的周围来回震荡,却落不到目标的点上。下面这幅图描述了这个现象:
事实上,学习率的取值未必一定要是一个常数,关于这个值的设定有很多的研究。
下面是比较常见的一些改进算法。
AdaGrad
AdaGrad是Adaptive Gradient的简写,该算法会为每个参数设定不同的学习率。它使用历史梯度的平方和作为基础来进行计算。
其算法公式如下:
[ heta_i = heta_i - frac{lambda}{sqrt{G_t + epsilon}} abla f( heta)]
对比式7,这里的改动就在于分号下面的根号。
根号中有两个符号,第二个符号比较好理解,它就是为了避免除0而人为引入的一个很小的常数,例如可以设为:0.001。
第一个符号的表达式展开如下:
[G_t = sum_{i = 1}^{t} abla f( heta){i} abla f( heta){i}^{T}]
这个值其实是历史中每次梯度的平方的累加和。
AdaGrad算法能够在训练中自动的对learning rate进行调整,对于出现频率较低参数采用较大的学习率;相反,对于出现频率较高的参数采用较小的学习率。因此,Adagrad非常适合处理稀疏数据。
但该算法的缺点是它可能导致学习率非常小以至于算法收敛非常的慢。
关于这个算法的直观解释可以看李宏毅教授的视频课程:ML Lecture 3-1: Gradient Descent。
RMSProp
RMS是Root Mean Square的简写。RMSProp是AI教父Geoff Hinton提出的一种自适应学习率方法。AdaGrad会累加之前所有的梯度平方,而RMSProp仅仅是计算对应的平均值,因此可缓解Adagrad算法学习率下降较快的问题。
该算法的公式如下:
[E[ abla f( heta_{i})^2]^{t} = gamma E[ abla f( heta_{i})^2]^{t - 1} + (1-gamma)( abla f( heta_{i})^{t})^{2} heta_i = heta_i - frac{lambda}{sqrt{E[g^2]^{t+1} + epsilon}} abla f( heta_{i})]
类似的,是为了避免除0而引入。 是衰退参数,通常设为0.9。
这里的 是t时刻梯度平方的平均值。
Adam
Adam是Adaptive Moment Estimation的简写。它利用梯度的一阶矩估计和二阶矩估计动态调整每个参数的学习率。
Adam的优点主要在于经过偏置校正后,每一次迭代学习率都有个确定范围,使得参数比较平稳。
该算法公式如下:
[m^{t} = eta_{1} m^{t-1} + (1-eta_{1}) abla f( heta) v^{t} = eta_{2} v^{t-1} + (1-eta_{2}) abla f( heta)^2 widehat{m}^{t} = frac{m^{t}}{1 - eta^{t}_1} widehat{v}^{t} = frac{v^{t}}{1 - eta^{t}_2} heta = heta - frac{lambda}{sqrt{widehat{v}^{t}} + epsilon}widehat{m}^{t}]
,分别是对梯度的一阶矩估计和二阶矩估计。, 是对,的校正,这样可以近似为对期望的无偏估计。
Adam算法的提出者建议 默认值为0.9,默认值为0.999,默认值为 。
在实际应用中 ,Adam较为常用,它可以比较快地得到一个预估结果。
优化小结
这里我们列举了几种优化算法。它们很难说哪种最好,不同的算法适合于不同的场景。在实际的工程中,可能需要逐个尝试一下才能确定选择哪一个,这个过程也是目前现阶段AI项目要经历的工序之一。
实际上,该方面的研究远不止于此,如果有兴趣,可以继续阅读 《Sebastian Ruder: An overview of gradient descent optimization algorithms》 这篇论文或者 Optimization for Deep Learning 这个Slides进行更多的研究。
由于篇幅所限,这里不再继续展开了。
算法限制
梯度下降算法存在一定的限制。首先,它要求函数必须是可微分的,对于不可微的函数,无法使用这种方法。
除此之外,在某些情况下,使用梯度下降算法在接近极值点的时候可能收敛速度很慢,或者产生Z字形的震荡。这一点需要通过调整学习率来回避。
另外,梯度下降还会遇到下面两类问题。
局部最小值
局部最小值(Local Minima)指的是,我们找到的最小值仅仅是一个区域内的最小值,而并非全局的。由于算法的起点是随意取的,以下面这个图形为例,我们很容易落到局部最小值的点里面。
这就是好像你从上顶往下走,你第一次走到的平台未必是山脚,它有可能只是半山腰的一个平台的而已。
算法的起点决定了算法收敛的速度以及是否会落到局部最小值上。
坏消息是,目前似乎没有特别好的方法来确定选取那个点作为起点是比较好的,这就有一点看运气的成分了。多次尝试不同的随机点或许是一个比较好的方法,这也就是为什么做算法的优化这项工作是特别消耗时间的了。
但好消息是:
鞍点
除了Local Minima,在梯度下降的过程中,还有可能遇到另外一种情况,即:鞍点(Saddle Point)。鞍点指的是我们找到点某个点确实是梯度为0,但它却不是函数的极值,它的周围既有比它小的值,也有比它大的值。这就好像马鞍一样。
如下图所示:
多类随机函数表现出以下性质:在低维空间中,局部极值很普遍。但在高维空间中,局部极值比较少见,而鞍点则很常见。
不过对于鞍点,可以通过数学方法Hessian矩阵来确定。关于这点,这里就不再展开了,有兴趣的读者可以以这里提供的几个链接继续探索。
参考资料与推荐读物
‘叁’ python逻辑回归怎么求正系数
Python 逻辑回归求正系数的方法可以分为两种:
1. 使用线性模型的求解方法:可以使用sklearn中的LogisticRegression类来求解正系数,调用其中的fit()方法就可以求解出正系数。
2. 使用梯度下降法:可以自己实现梯度下降法,通过不断迭代更新正系数,最终获得最优的正系数。
‘肆’ python gradientboostingregressor可以做预测吗
可以
最近项目中涉及基于Gradient Boosting Regression 算法拟合时间序列曲线的内容,利用python机器学习包scikit-learn 中的GradientBoostingRegressor完成
因此就学习了下Gradient Boosting算法,在这里分享下我的理解
Boosting 算法简介
Boosting算法,我理解的就是两个思想:
1)“三个臭皮匠顶个诸葛亮”,一堆弱分类器的组合就可以成为一个强分类器;
2)“知错能改,善莫大焉”,不断地在错误中学习,迭代来降低犯错概率
当然,要理解好Boosting的思想,首先还是从弱学习算法和强学习算法来引入:
1)强学习算法:存在一个多项式时间的学习算法以识别一组概念,且识别的正确率很高;
2)弱学习算法:识别一组概念的正确率仅比随机猜测略好;
Kearns & Valiant证明了弱学习算法与强学习算法的等价问题,如果两者等价,只需找到一个比随机猜测略好的学习算法,就可以将其提升为强学习算法。
那么是怎么实现“知错就改”的呢?
Boosting算法,通过一系列的迭代来优化分类结果,每迭代一次引入一个弱分类器,来克服现在已经存在的弱分类器组合的shortcomings
在Adaboost算法中,这个shortcomings的表征就是权值高的样本点
而在Gradient Boosting算法中,这个shortcomings的表征就是梯度
无论是Adaboost还是Gradient Boosting,都是通过这个shortcomings来告诉学习器怎么去提升模型,也就是“Boosting”这个名字的由来吧
Adaboost算法
Adaboost是由Freund 和 Schapire在1997年提出的,在整个训练集上维护一个分布权值向量W,用赋予权重的训练集通过弱分类算法产生分类假设(基学习器)y(x),然后计算错误率,用得到的错误率去更新分布权值向量w,对错误分类的样本分配更大的权值,正确分类的样本赋予更小的权值。每次更新后用相同的弱分类算法产生新的分类假设,这些分类假设的序列构成多分类器。对这些多分类器用加权的方法进行联合,最后得到决策结果。
其结构如下图所示:
可以发现,如果要用Gradient Boosting 算法的话,在sklearn包里调用还是非常方便的,几行代码即可完成,大部分的工作应该是在特征提取上。
感觉目前做数据挖掘的工作,特征设计是最重要的,据说现在kaggle竞赛基本是GBDT的天下,优劣其实还是特征上,感觉做项目也是,不断的在研究数据中培养对数据的敏感度。
‘伍’ 如何通过Python进行深度学习
作者 | Vihar Kurama
编译 | 荷叶
来源 | 云栖社区
摘要:深度学习背后的主要原因是人工智能应该从人脑中汲取灵感。本文就用一个小例子无死角的介绍一下深度学习!
人脑模拟
深度学习背后的主要原因是人工智能应该从人脑中汲取灵感。此观点引出了“神经网络”这一术语。人脑中包含数十亿个神经元,它们之间有数万个连接。很多情况下,深度学习算法和人脑相似,因为人脑和深度学习模型都拥有大量的编译单元(神经元),这些编译单元(神经元)在独立的情况下都不太智能,但是当他们相互作用时就会变得智能。
我认为人们需要了解到深度学习正在使得很多幕后的事物变得更好。深度学习已经应用于谷歌搜索和图像搜索,你可以通过它搜索像“拥抱”这样的词语以获得相应的图像。-杰弗里·辛顿
神经元
神经网络的基本构建模块是人工神经元,它模仿了人类大脑的神经元。这些神经元是简单、强大的计算单元,拥有加权输入信号并且使用激活函数产生输出信号。这些神经元分布在神经网络的几个层中。
inputs 输入 outputs 输出 weights 权值 activation 激活
人工神经网络的工作原理是什么?
深度学习由人工神经网络构成,该网络模拟了人脑中类似的网络。当数据穿过这个人工网络时,每一层都会处理这个数据的一方面,过滤掉异常值,辨认出熟悉的实体,并产生最终输出。
输入层:该层由神经元组成,这些神经元只接收输入信息并将它传递到其他层。输入层的图层数应等于数据集里的属性或要素的数量。输出层:输出层具有预测性,其主要取决于你所构建的模型类型。隐含层:隐含层处于输入层和输出层之间,以模型类型为基础。隐含层包含大量的神经元。处于隐含层的神经元会先转化输入信息,再将它们传递出去。随着网络受训练,权重得到更新,从而使其更具前瞻性。
神经元的权重
权重是指两个神经元之间的连接的强度或幅度。你如果熟悉线性回归的话,可以将输入的权重类比为我们在回归方程中用的系数。权重通常被初始化为小的随机数值,比如数值0-1。
前馈深度网络
前馈监督神经网络曾是第一个也是最成功的学习算法。该网络也可被称为深度网络、多层感知机(MLP)或简单神经网络,并且阐明了具有单一隐含层的原始架构。每个神经元通过某个权重和另一个神经元相关联。
该网络处理向前处理输入信息,激活神经元,最终产生输出值。在此网络中,这称为前向传递。
inputlayer 输入层 hidden layer 输出层 output layer 输出层
激活函数
激活函数就是求和加权的输入到神经元的输出的映射。之所以称之为激活函数或传递函数是因为它控制着激活神经元的初始值和输出信号的强度。
用数学表示为:
我们有许多激活函数,其中使用最多的是整流线性单元函数、双曲正切函数和solfPlus函数。
激活函数的速查表如下:
反向传播
在网络中,我们将预测值与预期输出值相比较,并使用函数计算其误差。然后,这个误差会传回这个网络,每次传回一个层,权重也会根绝其导致的误差值进行更新。这个聪明的数学法是反向传播算法。这个步骤会在训练数据的所有样本中反复进行,整个训练数据集的网络更新一轮称为一个时期。一个网络可受训练数十、数百或数千个时期。
prediction error 预测误差
代价函数和梯度下降
代价函数度量了神经网络对给定的训练输入和预期输出“有多好”。该函数可能取决于权重、偏差等属性。
代价函数是单值的,并不是一个向量,因为它从整体上评估神经网络的性能。在运用梯度下降最优算法时,权重在每个时期后都会得到增量式地更新。
兼容代价函数
用数学表述为差值平方和:
target 目标值 output 输出值
权重更新的大小和方向是由在代价梯度的反向上采取步骤计算出的。
其中η 是学习率
其中Δw是包含每个权重系数w的权重更新的向量,其计算方式如下:
target 目标值 output 输出值
图表中会考虑到单系数的代价函数
initial weight 初始权重 gradient 梯度 global cost minimum 代价极小值
在导数达到最小误差值之前,我们会一直计算梯度下降,并且每个步骤都会取决于斜率(梯度)的陡度。
多层感知器(前向传播)
这类网络由多层神经元组成,通常这些神经元以前馈方式(向前传播)相互连接。一层中的每个神经元可以直接连接后续层的神经元。在许多应用中,这些网络的单元会采用S型函数或整流线性单元(整流线性激活)函数作为激活函数。
现在想想看要找出处理次数这个问题,给定的账户和家庭成员作为输入
要解决这个问题,首先,我们需要先创建一个前向传播神经网络。我们的输入层将是家庭成员和账户的数量,隐含层数为1, 输出层将是处理次数。
将图中输入层到输出层的给定权重作为输入:家庭成员数为2、账户数为3。
现在将通过以下步骤使用前向传播来计算隐含层(i,j)和输出层(k)的值。
步骤:
1, 乘法-添加方法。
2, 点积(输入*权重)。
3,一次一个数据点的前向传播。
4, 输出是该数据点的预测。
i的值将从相连接的神经元所对应的输入值和权重中计算出来。
i = (2 * 1) + (3* 1) → i = 5
同样地,j = (2 * -1) + (3 * 1) → j =1
K = (5 * 2) + (1* -1) → k = 9
Python中的多层感知器问题的解决
激活函数的使用
为了使神经网络达到其最大预测能力,我们需要在隐含层应用一个激活函数,以捕捉非线性。我们通过将值代入方程式的方式来在输入层和输出层应用激活函数。
这里我们使用整流线性激活(ReLU):
用Keras开发第一个神经网络
关于Keras:
Keras是一个高级神经网络的应用程序编程接口,由Python编写,能够搭建在TensorFlow,CNTK,或Theano上。
使用PIP在设备上安装Keras,并且运行下列指令。
在keras执行深度学习程序的步骤
1,加载数据;
2,创建模型;
3,编译模型;
4,拟合模型;
5,评估模型。
开发Keras模型
全连接层用Dense表示。我们可以指定层中神经元的数量作为第一参数,指定初始化方法为第二参数,即初始化参数,并且用激活参数确定激活函数。既然模型已经创建,我们就可以编译它。我们在底层库(也称为后端)用高效数字库编译模型,底层库可以用Theano或TensorFlow。目前为止,我们已经完成了创建模型和编译模型,为进行有效计算做好了准备。现在可以在PIMA数据上运行模型了。我们可以在模型上调用拟合函数f(),以在数据上训练或拟合模型。
我们先从KERAS中的程序开始,
神经网络一直训练到150个时期,并返回精确值。
‘陆’ Python实现简单多线程任务队列
Python实现简单多线程任务队列
最近我在用梯度下降算法绘制神经网络的数据时,遇到了一些算法性能的问题。梯度下降算法的代码如下(伪代码):
defgradient_descent(): # the gradient descent code plotly.write(X, Y)
一般来说,当网络请求 plot.ly 绘图时会阻塞等待返回,于是也会影响到其他的梯度下降函数的执行速度。
一种解决办法是每调用一次 plotly.write 函数就开启一个新的线程,但是这种方法感觉不是很好。 我不想用一个像 cerely(一种分布式任务队列)一样大而全的任务队列框架,因为框架对于我的这点需求来说太重了,并且我的绘图也并不需要 redis 来持久化数据。
那用什么办法解决呢?我在 python 中写了一个很小的任务队列,它可以在一个单独的线程中调用 plotly.write函数。下面是程序代码。
classTaskQueue(Queue.Queue):
首先我们继承 Queue.Queue 类。从 Queue.Queue 类可以继承 get 和 put 方法,以及队列的行为。
def__init__(self, num_workers=1): Queue.Queue.__init__(self) self.num_workers=num_workers self.start_workers()
初始化的时候,我们可以不用考虑工作线程的数量。
defadd_task(self, task,*args,**kwargs): args=argsor() kwargs=kwargsor{} self.put((task, args, kwargs))
我们把 task, args, kwargs 以元组的形式存储在队列中。*args 可以传递数量不等的参数,**kwargs 可以传递命名参数。
defstart_workers(self): foriinrange(self.num_workers): t=Thread(target=self.worker) t.daemon=True t.start()
我们为每个 worker 创建一个线程,然后在后台删除。
下面是 worker 函数的代码:
defworker(self): whileTrue: tupl=self.get() item, args, kwargs=self.get() item(*args,**kwargs) self.task_done()
worker 函数获取队列顶端的任务,并根据输入参数运行,除此之外,没有其他的功能。下面是队列的代码:
我们可以通过下面的代码测试:
defblokkah(*args,**kwargs): time.sleep(5) print“Blokkah mofo!” q=TaskQueue(num_workers=5) foriteminrange(1): q.add_task(blokkah) q.join()# wait for all the tasks to finish. print“Alldone!”
Blokkah 是我们要做的任务名称。队列已经缓存在内存中,并且没有执行很多任务。下面的步骤是把主队列当做单独的进程来运行,这样主程序退出以及执行数据库持久化时,队列任务不会停止运行。但是这个例子很好地展示了如何从一个很简单的小任务写成像工作队列这样复杂的程序。
defgradient_descent(): # the gradient descent code queue.add_task(plotly.write, x=X, y=Y)
修改之后,我的梯度下降算法工作效率似乎更高了。如果你很感兴趣的话,可以参考下面的代码。 classTaskQueue(Queue.Queue): def__init__(self, num_workers=1):Queue.Queue.__init__(self)self.num_workers=num_workersself.start_workers() defadd_task(self, task,*args,**kwargs):args=argsor()kwargs=kwargsor{}self.put((task, args, kwargs)) defstart_workers(self):foriinrange(self.num_workers):t=Thread(target=self.worker)t.daemon=Truet.start() defworker(self):whileTrue:tupl=self.get()item, args, kwargs=self.get()item(*args,**kwargs)self.task_done() deftests():defblokkah(*args,**kwargs):time.sleep(5)print"Blokkah mofo!" q=TaskQueue(num_workers=5) foriteminrange(10):q.add_task(blokkah) q.join()# block until all tasks are doneprint"All done!" if__name__=="__main__":tests()
‘柒’ 从零开始用Python构建神经网络
从零开始用Python构建神经网络
动机:为了更加深入的理解深度学习,我们将使用 python 语言从头搭建一个神经网络,而不是使用像 Tensorflow 那样的封装好的框架。我认为理解神经网络的内部工作原理,对数据科学家来说至关重要。
这篇文章的内容是我的所学,希望也能对你有所帮助。
神经网络是什么?
介绍神经网络的文章大多数都会将它和大脑进行类比。如果你没有深入研究过大脑与神经网络的类比,那么将神经网络解释为一种将给定输入映射为期望输出的数学关系会更容易理解。
神经网络包括以下组成部分
? 一个输入层,x
? 任意数量的隐藏层
? 一个输出层,?
? 每层之间有一组权值和偏置,W and b
? 为隐藏层选择一种激活函数,σ。在教程中我们使用 Sigmoid 激活函数
下图展示了 2 层神经网络的结构(注意:我们在计算网络层数时通常排除输入层)
2 层神经网络的结构
用 Python 可以很容易的构建神经网络类
训练神经网络
这个网络的输出 ? 为:
你可能会注意到,在上面的等式中,输出 ? 是 W 和 b 函数。
因此 W 和 b 的值影响预测的准确率. 所以根据输入数据对 W 和 b 调优的过程就被成为训练神经网络。
每步训练迭代包含以下两个部分:
? 计算预测结果 ?,这一步称为前向传播
? 更新 W 和 b,,这一步成为反向传播
下面的顺序图展示了这个过程:
前向传播
正如我们在上图中看到的,前向传播只是简单的计算。对于一个基本的 2 层网络来说,它的输出是这样的:
我们在 NeuralNetwork 类中增加一个计算前向传播的函数。为了简单起见我们假设偏置 b 为0:
但是我们还需要一个方法来评估预测结果的好坏(即预测值和真实值的误差)。这就要用到损失函数。
损失函数
常用的损失函数有很多种,根据模型的需求来选择。在本教程中,我们使用误差平方和作为损失函数。
误差平方和是求每个预测值和真实值之间的误差再求和,这个误差是他们的差值求平方以便我们观察误差的绝对值。
训练的目标是找到一组 W 和 b,使得损失函数最好小,也即预测值和真实值之间的距离最小。
反向传播
我们已经度量出了预测的误差(损失),现在需要找到一种方法来传播误差,并以此更新权值和偏置。
为了知道如何适当的调整权值和偏置,我们需要知道损失函数对权值 W 和偏置 b 的导数。
回想微积分中的概念,函数的导数就是函数的斜率。
梯度下降法
如果我们已经求出了导数,我们就可以通过增加或减少导数值来更新权值 W 和偏置 b(参考上图)。这种方式被称为梯度下降法。
但是我们不能直接计算损失函数对权值和偏置的导数,因为在损失函数的等式中并没有显式的包含他们。因此,我们需要运用链式求导发在来帮助计算导数。
链式法则用于计算损失函数对 W 和 b 的导数。注意,为了简单起见。我们只展示了假设网络只有 1 层的偏导数。
这虽然很简陋,但是我们依然能得到想要的结果—损失函数对权值 W 的导数(斜率),因此我们可以相应的调整权值。
现在我们将反向传播算法的函数添加到 Python 代码中
为了更深入的理解微积分原理和反向传播中的链式求导法则,我强烈推荐 3Blue1Brown 的如下教程:
Youtube:https://youtu.be/tIeHLnjs5U8
整合并完成一个实例
既然我们已经有了包括前向传播和反向传播的完整 Python 代码,那么就将其应用到一个例子上看看它是如何工作的吧。
神经网络可以通过学习得到函数的权重。而我们仅靠观察是不太可能得到函数的权重的。
让我们训练神经网络进行 1500 次迭代,看看会发生什么。 注意观察下面每次迭代的损失函数,我们可以清楚地看到损失函数单调递减到最小值。这与我们之前介绍的梯度下降法一致。
让我们看看经过 1500 次迭代后的神经网络的最终预测结果:
经过 1500 次迭代训练后的预测结果
我们成功了!我们应用前向和方向传播算法成功的训练了神经网络并且预测结果收敛于真实值。
注意预测值和真实值之间存在细微的误差是允许的。这样可以防止模型过拟合并且使得神经网络对于未知数据有着更强的泛化能力。
下一步是什么?
幸运的是我们的学习之旅还没有结束,仍然有很多关于神经网络和深度学习的内容需要学习。例如:
? 除了 Sigmoid 以外,还可以用哪些激活函数
? 在训练网络的时候应用学习率
? 在面对图像分类任务的时候使用卷积神经网络
我很快会写更多关于这个主题的内容,敬请期待!
最后的想法
我自己也从零开始写了很多神经网络的代码
虽然可以使用诸如 Tensorflow 和 Keras 这样的深度学习框架方便的搭建深层网络而不需要完全理解其内部工作原理。但是我觉得对于有追求的数据科学家来说,理解内部原理是非常有益的。
这种练习对我自己来说已成成为重要的时间投入,希望也能对你有所帮助
‘捌’ python怎么实现人工智能
程序学习的过程就是使用梯度下降改变算法模型参数的过程。
比如说f(x) = aX+b; 这里面的参数是a和b,使用数据训练算法模型来改变参数,达到算法模型可以实现人脸识别、语音识别的目的。
实现人工智能的根本是算法,python是实现算法的一种语言,因为python语言的易用性和数据处理的友好性,所以现在很多用python语言做机器学习。其它语言比如java、c++等也也可以实现人工智能相关算法。下图是一个神经网络的示意图。
‘玖’ 怎么看python中逻辑回归输出的解释
以下为python代码,由于训练数据比较少,这边使用了批处理梯度下降法,没有使用增量梯度下降法。
##author:lijiayan##data:2016/10/27
##name:logReg.pyfrom numpy import *import matplotlib.pyplot as pltdef loadData(filename):
data = loadtxt(filename)
m,n = data.shape print 'the number of examples:',m print 'the number of features:',n-1 x = data[:,0:n-1]
y = data[:,n-1:n] return x,y#the sigmoid functiondef sigmoid(z): return 1.0 / (1 + exp(-z))#the cost functiondef costfunction(y,h):
y = array(y)
h = array(h)
J = sum(y*log(h))+sum((1-y)*log(1-h)) return J# the batch gradient descent algrithmdef gradescent(x,y):
m,n = shape(x) #m: number of training example; n: number of features x = c_[ones(m),x] #add x0 x = mat(x) # to matrix y = mat(y)
a = 0.0000025 # learning rate maxcycle = 4000 theta = zeros((n+1,1)) #initial theta J = [] for i in range(maxcycle):
h = sigmoid(x*theta)
theta = theta + a * (x.T)*(y-h)
cost = costfunction(y,h)
J.append(cost)
plt.plot(J)
plt.show() return theta,cost#the stochastic gradient descent (m should be large,if you want the result is good)def stocGraddescent(x,y):
m,n = shape(x) #m: number of training example; n: number of features x = c_[ones(m),x] #add x0 x = mat(x) # to matrix y = mat(y)
a = 0.01 # learning rate theta = ones((n+1,1)) #initial theta J = [] for i in range(m):
h = sigmoid(x[i]*theta)
theta = theta + a * x[i].transpose()*(y[i]-h)
cost = costfunction(y,h)
J.append(cost)
plt.plot(J)
plt.show() return theta,cost#plot the decision boundarydef plotbestfit(x,y,theta):
plt.plot(x[:,0:1][where(y==1)],x[:,1:2][where(y==1)],'ro')
plt.plot(x[:,0:1][where(y!=1)],x[:,1:2][where(y!=1)],'bx')
x1= arange(-4,4,0.1)
x2 =(-float(theta[0])-float(theta[1])*x1) /float(theta[2])
plt.plot(x1,x2)
plt.xlabel('x1')
plt.ylabel(('x2'))
plt.show()def classifyVector(inX,theta):
prob = sigmoid((inX*theta).sum(1)) return where(prob >= 0.5, 1, 0)def accuracy(x, y, theta):
m = shape(y)[0]
x = c_[ones(m),x]
y_p = classifyVector(x,theta)
accuracy = sum(y_p==y)/float(m) return accuracy
调用上面代码:
from logReg import *
x,y = loadData("horseColicTraining.txt")
theta,cost = gradescent(x,y)print 'J:',cost
ac_train = accuracy(x, y, theta)print 'accuracy of the training examples:', ac_train
x_test,y_test = loadData('horseColicTest.txt')
ac_test = accuracy(x_test, y_test, theta)print 'accuracy of the test examples:', ac_test
学习速率=0.0000025,迭代次数=4000时的结果:
似然函数走势(J = sum(y*log(h))+sum((1-y)*log(1-h))),似然函数是求最大值,一般是要稳定了才算最好。
从上面这个例子,我们可以看到对特征进行归一化操作的重要性。