机器学习 之 牛顿法和梯度下降法原理与实现

泰勒定理

f(x)在x0点可展开为幂级数f(x)=∑∞i=0aif(x−x0)i，则f(x)在x0的N(x0,σ)邻域内有任意阶导数，且系数an=f(n)(x0)n!。因此

f(x)=∑i=0∞aif(x−x0)i=∑i=0∞f(i)(x0)i!f(x−x0)i

称为f(x)在x0的泰勒级数，系数称为泰勒系数。当x0=0时，称为麦克劳林级数。

牛顿法求根

牛顿法的求根方法其实就是泰勒公式的一阶展开。

首先平方根的函数y=x−−√构造以y为自变量的函数 

f(y)=y2−x

将其按照泰勒公式进行一阶展开后，得到 

f(y)=f(y0)+f′(y0)(y−y0)=0

移项后得到 

y=y0−f(y0)f′(y0)

因此得到一阶展开后的通项公式 

yn+1=yn−f(yn)f′(yn)

将公式带入后，得到 

yn+1=yn−y2n−x2yn=12(yn−xyn)

<code class="hljs cpp has-numbering" style="display: block; padding: 0px; color: inherit; box-sizing: border-box; white-space: pre; border-radius: 0px; word-wrap: normal; background: transparent;"><span style="font-family:Times New Roman;font-size:18px;"><span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">double</span> <span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">sqrt</span>(<span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">double</span> x){    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">double</span> y = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1.0</span>;    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">while</span>(<span class="hljs-built_in" style="color: rgb(102, 0, 102); box-sizing: border-box;">fabs</span>(y * y - x) >= <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1e-9</span>){        y = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">0.5</span> * (y - x / y ); <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">// 对应上面的公式</span>    }     <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">return</span> y;}</span></code><ul class="pre-numbering" style="box-sizing: border-box; position: absolute; width: 50px; top: 0px; left: 0px; margin: 0px; padding: 6px 0px 40px; border-right-width: 1px; border-right-style: solid; border-right-color: rgb(221, 221, 221); list-style: none; text-align: right; background-color: rgb(238, 238, 238);"><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">1</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">2</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">3</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">4</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">5</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">6</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">7</span></li></ul>

牛顿法求平方根倒数

魔数的平方根倒数算法。

<code class="hljs cpp has-numbering" style="display: block; padding: 0px; color: inherit; box-sizing: border-box; white-space: pre; border-radius: 0px; word-wrap: normal; background: transparent;"><span style="font-family:Times New Roman;font-size:18px;"><span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">float</span> Q_rsqrt( <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">float</span> number ){    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">long</span> i;    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">float</span> x2, y;    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">const</span> <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">float</span> threehalfs = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1.5F</span>;    x2 = number * <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">0.5F</span>;    y  = number;    i  = * ( <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">long</span> * ) &y;                       <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">// evil floating point bit level hacking</span>    i  = <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">0x5f3759df</span> - ( i >> <span class="hljs-number" style="color: rgb(0, 102, 102); box-sizing: border-box;">1</span> );               <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">// what the fuck?</span>    y  = * ( <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">float</span> * ) &i;    y  = y * ( threehalfs - ( x2 * y * y ) );   <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">// 1st iteration </span>    <span class="hljs-comment" style="color: rgb(136, 0, 0); box-sizing: border-box;">//y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed</span>    <span class="hljs-keyword" style="color: rgb(0, 0, 136); box-sizing: border-box;">return</span> y;}</span></code><ul class="pre-numbering" style="box-sizing: border-box; position: absolute; width: 50px; top: 0px; left: 0px; margin: 0px; padding: 6px 0px 40px; border-right-width: 1px; border-right-style: solid; border-right-color: rgb(221, 221, 221); list-style: none; text-align: right; background-color: rgb(238, 238, 238);"><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">1</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">2</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">3</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">4</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">5</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">6</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">7</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">8</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">9</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">10</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">11</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">12</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">13</span></li><li style="box-sizing: border-box; padding: 0px 5px;"><span style="font-family:Times New Roman;font-size:18px;">14</span></li></ul>

首先求平方根的倒数y=1x√可以构造以y为自变量的函数，由如下函数表示 

f(y)=1y2−x=0

其一阶泰勒展开的通项公式为 

yn+1=yn−f(yn)f′(yn)=yn(3−xy2n)2=yn(1.5−xy2n2)

该公式对应的代码为y = y * ( threehalfs - ( x2 * y * y ) )，但是由于那个魔数的存在，它只需要迭代一轮。

梯度

假设函数u=f(x,y,z)在p(x,y,z)可微，那么

1. fx(x,y,z)是函数在x轴上的变化率
2. fy(x,y,z)是函数在y轴上的变化率
3. fz(x,y,z)是函数在z轴上的变化率

方向导数就是函数在任意一个方向上的变化率。于是产生一个问题：函数沿着哪个方向变化的时候能够取得最大值？

求解：函数u=f(x,y,z)在沿着向量l⃗ 的方向导数为 

φuφl=φuφxcosα+φuφycosβ+φuφzcosγ=(φuφx,φuφy,φuφz)(cosα,cosβ,cosγ)=g⃗ ∗k⃗ =|g⃗ |∗|k⃗ |cosθ

其中θ表示向量k和向量g的夹角，已知向量k的模为1，所以当θ=0，即向量k和向量g的方向一致时，方向导数取得最大值，最大为|g⃗ |

因此，梯度是一个矢量，它表示函数沿着该矢量的方向导数能够取得最大值，最大值为该矢量的模。 

grad={φuφx,φuφy,φuφz}

梯度下降

设损失函数为 

L=−∑xyi(wxi+b)

沿着梯度的方向，函数可以取得极大值，因此梯度下降的公式为 

w←w−αφLφw=w+η∑xyixib←b−αφLφb=b+η∑xyi

其中的η表示步长，或者叫学习率，用来表示每次更新的程度。

缺点：每次迭代都需要计算全部的数据集，当数据集非常大的时候，计算效率就非常低。

随机梯度下降

即不是一次使用所有误分类点，而是随机选择一个误分类点进行梯度下降。

已知梯度下降公式为 

gradL={φLφw,φLφb}={−∑xyixi,−∑xyi}

不再使用全部的误分类点，而是随机选择一个误分类点的子集对w和b进行更新，极端的情况就是每次只选择一个误分类点。

w←w+ηyixib←b+ηyi

这样，由于算法会由于在初始时选择初始值的不同，或者选择误分类点的顺序上的不同，导致最后得到的模型也会不一样。可以证明，误分类的次数是有限的，当训练数据集线性可分的时候，学习算法是收敛的。但当数据集线性不可分的时候，算法就不收敛了，迭代结果会发生震荡。

牛顿法与凸优化

其实就是泰勒公式的二阶展开。 

f(x)=f(x0)+f′(x0)(x−x0)+12f′′(x0)(x−x0)2

求函数f(x)的极小值点，就是求函数f′(x)=0的点，带入到上面的展开式，得到 

f′(x)=f′(x0)+f′′(x0)(x−x0)=0

移项后得到 

x=x0−f′(x0)f′′(x0)

因此得到迭代公式 

xn+1=xn−f′(xn)f′′(xn)

当x表示的是高维向量时，可将公式表示为 

Xn+1=Xn−∇f(Xn)∇2f(Xn)

使用牛顿法使得目标函数能够二次收敛，相比梯度下降法有更快的收敛速度。然而，牛顿法需要计算目标函数的二阶偏导数，且要求目标Hesse矩阵正定，难度比较大。

拟牛顿法

由于牛顿法需要计算目标函数的二阶偏导，且要求目标Hesse矩阵正定，难度比较大，因此有了拟合牛顿法。不过，不会用！ 
http://en.wikipedia.org/wiki/Quasi-Newton_method