牛顿法、梯度下降法原理及Python编程应用

牛顿法、梯度下降法原理及Python编程应用

一、项目概述

无论是在学习还是工作中，我们都会遇到很多最优化问题。最优化方法是一种数学方法，它是研究在给定约束之下如何寻求某些因素(的量)，以使某一(或某些)指标达到最优的一些学科的总称。最优化算法在学习和工作中是很重要的，我们学习和工作中遇到的大多问题都可以建模成一种最优化模型进行求解。常见的最优化方法有梯度下降法、牛顿法和拟牛顿法、共轭梯度法等等，本文主要介绍牛顿法和梯度下降法原理以及使用Python编程应用问题。

二、应用领域

使用牛顿法或梯度下降法，在最优化算法所解决的问题中，有着很广阔的应用领域，下面主要列举两点：

1. 企业利润分析，如何使得企业所获得的利润最大化。
2. 构建机器学习算法，大部分的机器学习算法的本质都是建立优化模型，通过最优化方法对目标函数（或损失函数）进行优化，从而训练出最好的模型。

三、算法原理

（一）牛顿法原理及优缺点

原理：

牛顿法是一种在实数域和复数域上近似求解方程的方法。方法使用函数f (x)的泰勒级数的前面几项来寻找方程f (x)=0的根。牛顿法最大的特点就在于它的收敛速度很快。牛顿法是一种迭代算法，每一步都需要求解目标函数的Hessian矩阵的逆矩阵，计算比较复杂。

优点：二阶收敛、收敛速度快。

缺点：计算复杂，局部收敛，当初始点选择不当时，往往会导致不收敛，二阶海塞矩阵必须可逆，否则算法进行困难。

计算过程如下：

g定义为雅克比矩阵，矩阵中各元素对应一阶偏导数，Hessian矩阵中各元素对应二阶偏导数。Hessian矩阵的逆同比于梯度下降法的学习率alpha，Hessian的逆在迭代中不断减小，起到逐渐缩小步长的效果。

（二）梯度下降法原理及优缺点

原理：

梯度下降法实现简单，当目标函数是凸函数时，梯度下降法的解是全局解。一般情况下，其解不保证是全局最优解，梯度下降法的速度也未必是最快的。梯度下降法的优化思想：用当前位置负梯度方向作为搜索方向，因为该方向为当前位置的最快下降方向，所以也被称为是”最速下降法“。最速下降法越接近目标值，步长越小，前进越慢。

优点：实现相对简单。

缺点：靠近极小值时收敛速度减慢，求解需要很多次的迭代；直线搜索时可能会产生一些问题，可能会“之字形”地下降。

计算过程如下：

四、Python编程设计
（一）牛顿法Python编程
用牛顿法求解一元函数：

from sympy import *
# step为迭代步数，x0为初始位置，obj为要求极值的函数
def newtons(step, x0, obj):
    i = 1 # 记录迭代次数的变量
    x0 = float(x0) # 浮点数计算更快
    obj_deri = diff(obj, x) # 定义一阶导数，对应上述公式
    obj_sec_deri = diff(obj, x, 2) # 定义二阶导数，对应上述公式
    while i <= step:
        if i == 1:
            # 第一次迭代的更新公式
            xnew = x0 - (obj_deri.subs(x, x0)/obj_sec_deri.subs(x, x0))
            print('迭代第%d次：%.5f' %(i, xnew))
            i = i + 1
        else:
            #后续迭代的更新公式
            xnew = xnew - (obj_deri.subs(x, xnew)/obj_sec_deri.subs(x, xnew))
            print('迭代第%d次：%.5f' % (i, xnew))
            i = i + 1
    return xnew
x = symbols("x") # x为字符变量
result = newtons(50, 10, x**6+x)
print('最佳迭代的位置：%.5f' %result)

运行结果：

（二）梯度下降法Python编程
用梯度法求解二次函数：

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
 
def Fun(x1,x2):#原函数
    return x1*x1+2*x2*x2-4*x1-2*x1*x2
 
def PxFun(x1,x2):#偏x导
    return 2*x1-2*x2-4
 
def PyFun(x1,x2):#偏y导
    return 4*x2-2*x1
 
#初始化
i=0     #迭代次数
fig=plt.figure()#figure对象
ax=Axes3D(fig)#Axes3D对象
X1,X2=np.mgrid[-2:2:40j,-2:2:40j]#取样并作满射联合
Z=Fun(X1,X2)#取样点Z坐标打表
ax.plot_surface(X1,X2,Z,rstride=1,cstride=1,cmap="rainbow")
ax.set_xlabel('x1')
ax.set_ylabel('x2')
ax.set_zlabel('z')
 
#梯度下降
step=0.1   #下降系数
x1=0
x2=0#初始选取一个点
tag_x1=[x1]
tag_x2=[x2]
tag_z=[Fun(x1,x2)]#三个坐标分别打入表中，该表用于绘制点
new_x1=x1
new_x2=x2
Over=False
while Over==False:
    new_x1-=step*PxFun(x1,x2)
    new_x2-=step*PyFun(x1,x2)#分别作梯度下降
    if Fun(x1,x2)-Fun(new_x1,new_x2)<7e-9:#精度
        Over=True
    x1=new_x1
    x2=new_x2#更新旧点
    tag_x1.append(x1)
    tag_x2.append(x2)
    tag_z.append(Fun(x1,x2))#新点三个坐标打入表中
    i=i+1
 
#绘制点/输出坐标
ax.plot(tag_x1,tag_x2,tag_z,'r.')
plt.title('(x1,x2)~('+str(x1)+","+str(x2)+')')
plt.show()
print("迭代次数：",i)

运行结果：

五、总结
梯度下降法和最小二乘法相比，梯度下降法需要选择步长，而最小二乘法不需要。梯度下降法是迭代求解，最小二乘法是计算解析解。如果样本量不算很大，且存在解析解，最小二乘法比起梯度下降法要有优势，计算速度很快。但是如果样本量很大，用最小二乘法由于需要求一个超级大的逆矩阵，这时就很难或者很慢才能求解解析解了，使用迭代的梯度下降法比较有优势。

参考文献

[1]. https://www.jianshu.com/p/e3552dc9ef09

[2]. https://blog.csdn.net/lsgqjh/article/details/79168095

[3]. https://www.jianshu.com/p/d892d0d13b6d

[4]. http://www.elecfans.com/d/722244.html
p/d892d0d13b6d>

[4]. http://www.elecfans.com/d/722244.html