机器/深度学习中的参数优化原理初探

1. 最小二乘估计法

0x1：最小二乘估计的发展历史

最小二乘法(Least Squares Method,简记为LSE)是一个比较古老的方法，源于天文学和测地学上的应用需要。

在早期数理统计方法的发展中，这两门科学起了很大的作用。丹麦统计学家霍尔把它们称为“数理统计学的母亲”。此后近三百年来，它广泛应用于科学实验与工程技术中。美国统计史学家斯蒂格勒( S. M. Stigler)指出, 最小二乘方法是19世纪数理统计学的压倒一切的主题。1815年时,这方法已成为法国、意大利和普鲁士在天文和测地学中的标准工具,到1825年时已在英国普遍使用。

追溯到1801年，意大利天文学家朱赛普·皮亚齐发现了第一颗小行星谷神星。经过40天的跟踪观测后，由于谷神星运行至太阳背后，使得皮亚齐失去了谷神星的位置。随后全世界的科学家利用皮亚齐的观测数据开始寻找谷神星，但是根据大多数人计算的结果来寻找谷神星都没有结果。时年24岁的高斯也计算了谷神星的轨道。奥地利天文学家海因里希·奥尔伯斯根据高斯计算出来的轨道重新发现了谷神星。高斯于其1809年的著作《关于绕日行星运动的理论》中。在此书中声称他自1799年以来就使用最小二乘方法,由此爆发了一场与勒让德的优先权之争。

近代学者经过对原始文献的研究,认为两人可能是独立发明了这个方法,但首先见于书面形式的，以勒让德为早。然而,现今教科书和著作中,多把这个发明权归功于高斯。其原因,除了高斯有更大的名气外,主要可能是因为其正态误差理论对这个方法的重要意义。

勒让德在其著作中,对最小二乘方法的优点有所阐述。然而,缺少误差分析。我们不知道,使用这个方法引起的误差如何,就需建立一种误差分析理论。

高斯于1823年在误差e1 ,… , en独立同分布的假定下,证明了最小二乘方法的一个最优性质: 在所有无偏的线性估计类中,最小二乘方法是其中方差最小的！

在德国10马克的钞票上有高斯像，并配了一条正态曲线。在高斯众多伟大的数学成就中挑选了这一条,亦可见这一成就对世界文明的影响。

现行的最小二乘法是勒让德( A. M. Legendre)于1805年在其著作《计算慧星轨道的新方法》中提出的。它的主要思想就是选择未知参数，使得理论值与观测值之差的平方和达到最小：

H=\sum_{0}^{m}{(y-y_{i})^2}

我们现在看来会觉得这个方法似乎平淡无奇,甚至是理所当然的。这正说明了创造性思维之可贵和不易。从一些数学大家未能在这个问题上有所突破,可以看出当时这个问题之困难。欧拉、拉普拉斯在许多很困难的数学问题上有伟大的建树,但在这个问题上未能成功。

在高斯发表其1809年著作之前，约在1780年左右,拉普拉斯已发现了概率论中的“中心极限定理”。

根据这个定理，大量独立的随机变量之和,若每个变量在和中起的作用都比较小，则和的分布必接近于正态。

测量误差正具有这种性质。一般地说,随机(而非系统)的测量误差，是出自大量不显著的来源的叠加。因此，中心极限定理给误差的正态性提供了一种合理的理论解释。这一点对高斯理论的圆满化很有意义,因为高斯原来的假定(平均数天然合理)总难免给人一种不自然的感觉。

耐人寻味的是,无论是中心极限定理的发明者拉普拉斯,还是早就了解这一结果的高斯,都没有从这个结果的启示中去考察误差分布问题。对前者而言,可能是出于思维定势的束缚,这对拉普拉斯来说可算不幸,他因此失掉了把这个重要分布冠以自己名字的机会(正态分布这个形式最早是狄莫弗( De Moiv re) 1733年在研究二项概率的近似计算时得出的。以后也有其他学者使用过,但都没有被冠以他们的名字。高斯之所以获得这一殊荣,无疑是因为他把正态分布与误差理论联系了起来) 。

可以说,没有高斯的正态误差理论配合, 最小二乘方法的意义和重要性可能还不到其现今所具有的十分之一。最小二乘方法方法与高斯误差理论的结合,是数理统计史上最重大的成就之一,其影响直到今日也尚未过时！

0x2：从一个日常的例子来看最小二乘估计的思想

1. 求平均数 - 数学直觉

来看一个生活中的例子。比如说，有五把尺子：

用它们来分别测量一线段的长度，得到的数值分别为（颜色指不同的尺子）：

这5次测量本质上就是代表了5次实验观测得到的观测样本。这种现象在实际的工程项目中是十分常见的，甚至可以说是百分百出现的，即样本噪音或样本偏差问题。

之所以出现不同的值可能因为：

  • 不同厂家的尺子的生产精度不同
  • 尺子材质不同，热胀冷缩不一样
  • 测量的时候心情起伏不定
  • ....

这种情况下，一个直观的解决思路是：取平均值来作为线段的长度。

这种解决思路是最符合”直觉“的，但是有人可能会提出质疑：

  • 这样做有道理吗？
  • 用调和平均数行不行？
  • 用中位数行不行？
  • 用几何平均数行不行？

2. 用最小二乘估计的思路来解这个问题

首先，把测试得到的值画在笛卡尔坐标系中，分别记作y_i ：

其次，把要猜测的线段长度的真实值用平行于横轴的直线来表示（因为是猜测的，所以用虚线来画），记作y ：

每个点都向y 做垂线，垂线的长度就是|y-y_i| ，也可以理解为测量值和真实值之间的误差：

因为误差是长度，取绝对值计算计算比较复杂，所以设计者用平方来代表误差：

|y-y_i|\to (y-y_i)^2\\

误差的平方和就是（\epsilon 代表误差）：

S_{\epsilon^2}=\sum (y-y_i)^2\\

下面的动图表示了y在不断变化下，误差的变化情况。

自然，误差的平方和S_{\epsilon^2} 在不断变化的。

法国数学家，阿德里安-马里·勒让德提出了一个划时代的定理：让总的误差的平方最小的y 就是真值。这是基于，如果误差是随机的，应该围绕真值上下波动。

勒让德的想法变成代数式就是：

这是一个二次函数，对其求导，导数为0的时候取得最小值：

化简得：

可以看到正好是算术平均数。这个巧合在1元1次方程中是刚好成立的，也就是说，算数平均是最小二乘估计的一个特例。

3. 最小二乘估计法

最小二乘估计法包含两个基本的元素：

1. 利用平方作为评价指标。要注意的是，深度学习中的平方损失函数只是一个指标，即用平方的方式进行损失函数的评估，不要和最小二乘法混淆。
2. 利用求导极值的方式，直接得到全局最优值。这点是最小二乘估计和梯度下降这种迭代算法最大的区别。

笔者插入：这里需要注意的是，导数为0并不是得到极值的充要条件，只是一个必要条件，对于多元高阶方程组，通过求导极值并不能保证得到最优值。

0x3：最小二乘估计的适用范围

最小二乘估计是一个通用的未知参数最优化方案，它可以对任何的模型进行参数优化。

从历史的角度来说，最小二乘估计的出现是因为当时回归建模的出现以及对回归模型进行参数估计的需求。

所以我们先花一些时间，来讨论下回归模型。

0. 回归模型

回归模型来源于当时科学和工程领域一个很基本的问题，数学建模。建模的思想在处理统计数据分析的每一个学科中都很有用。

例如，假设给定一个随机变量集，要完成的任务是：找到可能存在于这些变量之间的关系，如果这种关系存在的话。

作为函数逼近的一种特殊形式，回归的典型方案如下：

1. 选择随机变量中的一个变量作为关注变量，这一随机变量被称为依赖变量，或者响应变量（response）
2. 剩下的随机变量称为独立变量，或者回归量（regressor）。它们的作用是用来解释或者预测响应的统计行为。
3. 响应对回归量的依赖还包括一个附加的误差项，用来说明在对依赖程度公式化时候的不确定性。误差项称为期望误差（expectational error）或解释误差（explanational error）

上述这样的模型称为回归模型（regression model）。

回归模型有很多种类，例如我们熟知的线性回归模型就是其中一种：

1. Linear Regression线性回归
2. Logistic Regression逻辑回归
3. Polynomial Regression多项式回归
4. Stepwise Regression逐步回归
5. Ridge Regression岭回归
6. Lasso Regression套索回归
7. ElasticNet回归

在对具体问题建模的时候，要选取哪种模型呢？这个问题似乎没有最佳答案，需要具体问题具体分析，结合专家知识，领域经验来选定最适合的模型。例如我们在量化交易中第一步要做的事情是设计一个符合股市历史波动的回归模型，而公式的形式各家量化交易公司都不一样，这里面考验的就是对股市交易的本质理解了。

我们通过一个简单例子来讨论这个问题，比如温度与冰淇淋的销量（一个自变量、一个因变量）：

1. 线性回归估计

可以假设这种线性关系为：

上图的i,x,y 分别为：

总误差的平方为：

不同的a,b 会导致不同的S_{\epsilon^2} ，根据多元微积分的知识求偏导：

这个时候S_{\epsilon^2} 取最小值。对于a,b 而言，上述方程组为线性方程组，用之前的数据解出来：

也就是这根直线：

2. 非线性回归估计

同一个问题，还可以假设（这其实就是数学建模过程）：

在这个假设下，可以根据最小二乘法（求导极值），算出a,b,c ，得到下面这根红色的二次曲线：

3. 任意模型的最小二乘估计

同一组数据，选择不同的f(x) ，通过最小二乘法可以得到不一样的拟合曲线

 

0x4：最小二乘法和矩阵投影的关系 - 最小二乘估计的的几何解释

从矩阵投影的角度来看最小二乘，最小二乘法的本质是最小化系数矩阵所张成的向量空间到观测向量的欧式误差距离。

我们从几何的角度来帮助我们理解其中的概念含义。

1. 线性方程组的几何意义

为了从几何的角度解释最小二乘法，我们先回顾一下，线性方程组的几何意义。线性方程组可以从行和列两个角度看。 
举例，看如下简单的线性方程组：

1）从行的角度看

线性方程组的解就是直线a和直线b的交点，这一点应该谁都知道，如下图所示：

从图上很容易可以看出该线性方程组的解为x1=1，x2=2

从行的角度看，线性方程组的每一行都代表空间中的几何图形（二维时是直线，三维时是平面，高维时是高维空间中的一个”平面“），线性方程组的解就是这些图形的交点or交线or交面。 

2）从列的角度看

从列的角度看时，先把上面的线性方程组写成矩阵的形式

再把上式拆开，写成下面这个样子

表示向量 a1 的 x1 倍加上向量 a2 的 x2 倍等于向量 b 。这种角度就是从列向量的角度看线性方程组。该方程的解已经从行的角度看出来了，为x1=1，x2=2，于是我们把向量a1,a2,b画到下面这张图上：

向量 a1 的1倍加上向量 a2 的2倍刚好等于向量 b，而倍数1和2,就是我们的解 x1 和 x2

从列的角度看线性方程组Ax=b的解，就是为系数矩阵A里的每一列都寻找一个合适的倍数，使每一列乘上这个倍数后再相加刚好等于向量b，这个倍数就是解。本质上就是找到A里的列向量的一个线性组合使之等于b。

2. 最小二乘法就是解一个无解的线性方程组

这个topic是一个笔者认为非常重要意义的概念。

在机器/深度学习或者统计学任务中，我们基本上不可能找到一个100% matched的model来对目标问题域进行建模（我们没有上帝视角），也即真实物理模型是非常难获取的，像爱因斯坦质能方程那种跨世纪的发现在大多数的业务场景中都难以得到。所以，统计学家和机器学习工作者只能退而求其次，寻找一个“大致能拟合（尽可能拟合）”的模型来描述目标问题域。寻找这个大致拟合的模型的方法之一就是最小二乘法。

我们从一个最简单的例子开始，已知平面上有3个点(1,2)，(0,2)，(2,3)  

我们想用一条直线去拟合它。如上所述，我们没有上帝视角，所以我们并不知道真实物理模型的函数形式。

我们只能寻找一个大致能拟合的模型，这个过程笔者相信大家都经历过，例如在深度学习中编写网络结构的代码（DNN/LSTM/CNN层等），通过建立高斯模型的假设进行异常发现。

回到这个简单的线性例子，我们现在设这条直线的方程为 y=kx+b 。我们希望这条直线可以同时通过这三个点（100%拟合），也就是这条直线的参数要满足:

从图中直观的看，没有一条直线可以同时过这三个点，所以这个方程是无解的。怎么解一个无解的方程组呢？

为了表述方便，我们换一下符号，用x1表示k，用x2表示b。即：

写成矩阵形式：

从列的角度看：

一但化成列的形式，我们就很自然想到把向量a1，a2，b 画到图上:

要找到解，就要找到a1，a2的一个线性组合，使得组合后的向量刚好等于b。可惜的是，根据向量相加的原理，任何的a1和a2线性组合，只可能出现在a1,a2所在的平面S上（这个平面S就是a1和a2组成的向量空间），但是向量b不在平面S上，如下图。不可能找到解：

找不到完美的解，就只能找到一个最接近的解。所以现在问题转换为了哪个解是最接近的解。

我们的思路是在平面S上找一个最接近向量b的向量来代替向量b，记这个替代品向量为P。所以现在问题又转换为了如何评价最接近向量b。

根据向量相关性的内积定理，和向量b相关性最高的向量就是过向量b的终点做平面S 的垂线（也就是做投影），垂足就是代替向量P（最接近向量b的向量）的终点。

P与b之间的误差：e=b−P

原来的方程为Ax=b是无解的，我们用P代替b后，P在a1a2所在的平面上，所以现在方程Ax̂ =P就有解了。

但是这个解和目标真实值之间存在误差，误差为：

此时，误差向量e一定是垂直于平面S的，也就是要垂直于a1和a2。垂直在数学上是有明确公式表示的。也就是e∙a1=0，e∙a2=0，用矩阵表示出来就是，即：

综合上面公式得：，化简一下就是，所以最佳的近似解就是

这样，按公式解出的解出了最近似的解为(1/2,11/6)。 

从列的角度，我们就可以用a1和a2的线性组合来表示P，下图所示

那么最优的直线的斜率和截距就是我们解出的k=1/2，b=11/6=1.8333。如下图

0x5：最小二乘估计和最大似然估计的关系 - 从概率统计的角度来看

最小二乘估计的基本假设是：

这个假设真的合理吗？数学王子高斯（1777－1855）也像我们一样心存怀疑。

高斯换了一个思考框架，通过概率统计的思维框架来思考。

让我们回到最初测量线段长度的问题。高斯想，通过测量得到了这些值：

每次的测量值x_i 都和线段长度的真值x 之间存在一个误差：

根据大数定理，误差本身也属于非系统随机变量，因此这些误差最终会形成一个概率分布，只是现在不知道误差的概率分布是什么。假设概率密度函数为：

假设一个误差的联合概率，把所有的测量数据整合起来：

把x 作为变量的时候，上面就是似然函数了。关于极大似然估计的讨论，可以参阅我的另一篇blog。

L(x) 的图像可能是这样的：

根据极大似然估计的思想，L(x)的极大似然估计最应该是下面这点：

从最优化的角度来看，当下面这个式子成立时，取得最大值：

笔者插入：极大似然估计得到的极大似然点也我们用求导极值求最优点，在数学公式上是一致的，读者在看到这里的时候有没有体会到一种极致的美感。

如果最小二乘法是对的，那么x=\overline{x} 时应该取得最大值，即：

现在可以来解这个微分方程了。最终得到：

可以看到，这就是一个正态分布，并且这还是一个充要条件：

也就是说，如果误差的分布是正态分布，那么最小二乘法得到的就是极大似然的值。

接下来问题来了，那么误差的分布是正态分布吗？

如果误差是由于随机的、无数的、独立的、多个因素造成的，比如之前提到的：

  • 不同厂家的尺子的生产精度不同
  • 尺子材质不同，热胀冷缩不一样
  • 测量的时候心情起伏不定
  • ......

那么根据中心极限定理，误差的分布就应该是正态分布。

Relevant Link: 

https://www.zhihu.com/question/20447622
https://blog.csdn.net/Macer3/article/details/48394239 
https://www.zhihu.com/question/37031188
https://www.zhihu.com/question/20447622
https://blog.csdn.net/Macer3/article/details/48394239
https://www.zhihu.com/question/37031188

 

2. 维纳滤波器 - 一种负反馈自动控制系统

0x1：维纳的生平简略

诺伯特·维纳（Norbert Wiener，1894年11月26日－1964年3月18日），生于美国密苏里州哥伦比亚，美国应用数学家，在电子工程方面贡献良多。他是随机过程和噪声信号处理的先驱，又提出“控制论”一词。
他主要在哥廷根和剑桥，研究布朗运动、傅里叶变换、调和分析、狄利克雷问题和Tauber型定理（Abelian and tauberian theorems）等。1929年维纳指导当时在贝尔电话公司学习的博士生李郁荣研制了李-维纳网络，获得美国专利。1935年维纳应国立清华大学校长梅贻琦和数学系主任熊庆来之聘，到清华大学讲学，主讲傅立叶变换，听讲者包括华罗庚、段学复等。维纳曾推荐华罗庚和徐贤修合写的论文发表在麻省理工学院的《数学物理学报》（上）
二战时，他在枪炮控制方面工作，引发了他对通讯理论和反馈的兴趣，著《控制论：或关于在动物和机器中控制和通信的科学》一书，促成了控制论的诞生。 
http://www.ee.cityu.edu.hk/~gchen/pdf/Wiener_S.pdf

维纳的个人履历实在令人望而生畏，他真的将各种领域的知识融合在了一起。同时，他创造了自动控制理论，深刻地影响了机器/深度学习的发展历程。神经元模型以及GD梯度下降的发明也是基于自动控制理论基础上创立的。

笔者思考：
自动负反馈控制理论，从另一个角度讨论了一个信息处理系统如何达到稳态（即找到最优参数）。
笔者任务其核心思想是：比起一次接受所有历史状态信息，通过最小二乘估计得到最优参数。信息系统其实可以将时间序列切片，通过只获取一个有限时间段内的信息，对通过滤波器后的信息误差进行评估，进而对系统本身进行负反馈调节，这种方式最终可以让系统趋向于稳态。
这也是GD（梯度下降）思想的核心之一，当然，GD本身还解决了对高维多元随机变量求极值的计算复杂度问题，GD本身是易于通过计算机实现的。

0x2：数字滤波器的基本概念以及它和机器学习算法的关系

数字滤波器是数字信号处理中使用得罪广泛的一种线性系统环节，是数字信号处理的重要基础。

数字滤波器的本质是将一组输入的数字序列通过一定的运算后转变为另一组输出的数字序列。

1. 滤波器的定义

实现滤波处理的运算电路、或设备称为滤波器。

2. 滤波的定义

对输入信号通过一定的处理得到输出信号，这个处理通常是滤除输入信号的某些频率成分；保留信号中某些频率范围内的有用信号成分。所以把这种处理的过程称为滤波。

笔者思考：CNN卷积网络的训练过程体现了非常明显的滤波过程，CNN网络在训练样本数据的过程中，会自动地保留样本数据中和target有关的“关键性像素区域”，例如小猫识别任务中，CNN会保留图像样本中各种姿势的猫，而对背景这些“冗余信息”会逐渐通过权重调整过滤掉，从某种程度上来说，这就是一种滤波过程。

3. 自适应滤波器结构 - 一种包含负反馈的神经元模型

几乎所以的机器学习算法的参数优化（训练）过程都包含负反馈，算法通过在训练中不断根据本轮迭代的预测结果和目标结果之间的差距来动态调整自己的负反馈，从而逐渐将权值参数调整到”尽量完美“的状态（即拟合）。

下图是一个未知动态系统的流程方框图，由包含元素的输入向量所刺激，其中 i 是刺激（激励）应用于系统时的瞬时时间，作为时间序列的索引。

作为对刺激的相应，系统产生一个输出 y(i) 作为相应。因此，次系统的外部行为由下述数据集描述：

从数字信号的时空特性角度来看，刺激向量 x(i) 能够以两种根本不同的方式出现，一种是空间的，另一种是空间的：

1. x(i) 的 M 个元素代表空间中的不同点，输入向量 x(i) 的维数 M 称为输入空间的维数，或简称输入维数。在这种情况下我们称 x(i) 为数据的瞬像（snapshot）
2. x(i) 的 M 个元素代表在时间上均匀分布的某个刺激，的现在和（M-1）个过去的值组成的集合。即它是时域累计刺激。

我们现在面对的问题是如何通过建立一个简单线性神经元来设计未知动态系统的一个多输入-单输出模型（即滤波器模型）。

这个神经元模型是在一个算法的影响下运行的，此算法控制对神经元的突触权值的必要调整，同时记住以下要点：

1. 此算法从任意设定的一个神经元突触权值开始。这意味着建模者不需要对目标问题域有先验知识。
2. 为相应系统行为的统计变化，突触权值的调整是建立在连续的基础之上的，即把时间加进算法中，权值的调整是连续进行的。
3. 调整突触权值的计算在长度为一个采样周期的时间段内完成。

这样描述的神经元模型称为“自适应滤波器（adaptive filter）”，而其中负责进行调整的算法理论就是LMS（最小均方算法），LMS我们放到下一个章节来展开讨论，我们这里先关注滤波器系统本身。

虽然是在作为系统辨识（system identification）的任务背景下给出的描述，但自适应滤波器的特征还是有很广泛的应用。

下图是一个自适应滤波器的示意图，它的运行由两个连续过程组成：

1. 过滤过程（滤波过程），涉及两个信号计算
　　1）一个输入，记为 y(i)，它被产生以响应刺激向量 x(i) 的 M 个元素，即 x1(i),x2(i),.....,xM(i)
　　2）一个误差信号，记为 e(i)，它是通过比较输出 y(i) 和未知目标系统（要模型的实际物理系统）的相应输出 d(i) 而获得的。从效果上讲，d(i) 可作为一个期望响应信号（desired response）或者目标（target）信号。
2. 自适应过程，包括根据误差 e(i) 对神经元突触权值的自动调整。

这两个共同运作过程的组合构成了一个围绕神经元运作的反馈环（feedback loop）。

上述的这两个连续过程的产生原理如下：

因为神经元是线性的，输出 y(i) 为局部诱导域：。

神经元的输出 y(i) 要和未知系统在时刻 i 的相应输出 d(i) 作比较，通常 y(i) 与 d(i) 不等，因此它们的比较结果得到了误差信号：

误差信号 e(i) 用来对神经元突触权值的调整进行控制的方式，是由用于导出自适应滤波算法的代价函数决定的。

这个问题与无约束最优化问题密切相关，无约束最优化不仅可以用在线性自适应滤波器上，还可以应用在一般的神经网络上。

为了下一章节讨论LMS作准备，我们这里先讨论下自适应滤波算法中的无约束最优化问题。

0x3：无约束最优化

考虑代价函数，它是一个对未知权值（参数）向量 w， 连续可微（continuously differentiable）的函数。

函数将 w 映射为单个实数，它是一种度量，用来选择自适应滤波算法的权值（参数）向量 w，以使得模型以最优方式运行。

这样，代价函数就成功地将一个学习问题转换为了最优化问题。

我们想找到一个最优解满足条件：

也就是说，需要解决一个无约束的最优化问题，即：

选择适当的权值向量 w 最小化代价函数。

最优性的必要条件（注意不是充要条件）是：

这里是梯度算子（gradient operator）。权值向量的梯度为0，代表了已经优化到了极值，至少是局部极值。

同时是代价函数的梯度向量（gradient vector）：

一类特别适合自适应滤波器设计的无约束最优化算法是以局部迭代下降（iterative descent）思想为基础的：

从一个初始估计值 w(0) 开始，产生一系列权值向量 w(1),w(2),...，使得代价函数在算法的每次迭代中都要下降，即

我们希望算法最终收敛到最优解，但是这里也仅仅是希望，因为局部迭代下降并不能保证一定收敛到全局最优解。

下面我们来讨论几种以迭代下降思想的基本形式或变种形式的无约束最优化方法。

1. 最速下降法

在最速下降法中，对权值向量 w 的连续调整是在醉宿下降的方法进行的，即它是与梯度向量方向相反的，记为：

因此，最速下降法一般表示为：，这里是一个正常数，称为步长（stepsize）或学习速率（learning-rate parameter）。g(n) 是在 w(n) 处的梯度向量值。

为了证明最速下降法的公式满足，我们用 w(n) 附近的一阶泰勒级数展开来逼近，即：

上式对小的是适用的，将代入上式得：

上式表明，对正的学习率参数，代价函数每次迭代都必定是下降的，除非权值向量 w(n) 的梯度已经等于0（达到局部最优极值）。但这里提供的推导是近似的，只有当学习率足够小时才是正确的（泰勒分解成立）。

最速下降法收敛到最优解的速度是很慢的。此外，学习率参数对收敛行为有重要影响：

1. 当 η 很小时：算法的瞬时相应是平缓的（overdamped），这是由于 w(n) 的轨迹是 w 平面的一个光滑曲线，如下图 a）
2. 当 η 很大时：算法的瞬时相应是剧烈的（underdamped），这是由于 w(n) 的轨迹是锯齿（振荡）形的，如下图 b）
3. 当 η 超过某一临界值时，算法是不稳定的（即不收敛）

2. 牛顿法（Newton's method）

牛顿法的基本思想是：最小化代价函数在当前点 w(n) 周围的二次近似值；最小化在算法的每次迭代中都要进行。

具体来说，就是利用代价函数在点 w(n) 周围的二次泰勒级数展开式，我们得到：

式中，g(n) 是代价函数在点 w(n) 处的处的 M x 1 梯度向量。矩阵 H(n) 是在 w(n) 的 m 行 m 列 Hessian矩阵。的Hessian矩阵定义为：

上式需要代价函数关于 w 的元素二阶可微。

对式中进行微分并求极值得：，我们最小化了改变量。

解上式有关的方程得到：。

整合上面公式得：。

一般来说，牛顿法收敛得很快，而且不会出现最速下降法有时会出现的锯齿形情况。但是，应用牛顿法时， Hessian矩阵必须对每个 n 都是正定矩阵。

遗憾的是，一般不能保证在算法的每次迭代中 H(n) 都是正定矩阵。

假如 Hessian矩阵 H(n) 不正定，对牛顿法进行修正就有必要。在很多时候，牛顿法的最主要局限在于其计算复杂度。

3. Gauss-Newton法

为了处理牛顿法的计算复杂度过大，而不对其收敛行为做太严重的让步，可以使用Gauss-Newton法。为了应用这一方法，我们采用表示为误差平方和的代价函数，令：

误差信号 e(i) 是可调权值向量 w 的函数。给定操作点 w(n)，通过引入下面新的项来线性化 e(i) 对 w 的依赖性：

用矩阵记号可以写成等价的形式：

其中，e(n) 是误差向量

J(n) 是 e(n) 的 n x m  Jacobi 矩阵：

Jacobi 矩阵 J(n) 是 m x n 梯度矩阵的转置，这里：

现在更新后的权值向量 w(n+1) 定义为：

将代入上式得：

将以上表达式对 w 求微分并设结果为零，得到：

综合上式，可得：

上式描述了 Gauss-Newton方法的纯粹形式。

不想牛顿法必须知道代价函数的 Hessian矩阵，Gauss-Newton法只需要知道误差向量 e(n) 的 Jacobi矩阵。

但是，为了使 Gauss-Newton 迭代可计算，矩阵乘积必须是非奇异的。

4. GD（梯度下降法）

用负梯度作搜索方向，即令，是一种自然的选择。这种方法被称为梯度方法或梯度下降方法。

停止准则通常取为，其中是小正数。

注意：梯度下降是最速下降在欧式范数下的特殊情况。

Relevant Link:

https://www.zhihu.com/question/267021131

 

3. LMS（最小均方）算法简介

我们从最小二乘估计器引入最小二乘滤波器，这样可以很自然地进入对LMS的讨论中。最小二乘滤波器和最小二乘估计器虽然只有几字之差，但是其整个优化运算过程是不一样的。最小二乘滤波器引入了自适应反馈的思想。

0x1：最小二乘估计器 - 维纳滤波器的一种特定形式

我们在前面的章节中讨论了最小二乘估计器，它利用极小化（求导极值）的传统放来从环境的观测模型中找到最小二乘解。

从这个小节开始，我们将最小二乘估计器放到一个维纳滤波器的框架中进行讨论，我们称之为最小二乘滤波器（least-squares filter）。我们接下来利用 Gauss-Newton法来重新推导这个滤波器公式。

1. 最小二乘滤波器的误差向量

我们定义如下误差向量：

其中 d(n) 是 n x 1 的期望响应向量：

X(n) 是 n x M 的数据矩阵：

误差向量 e(n) 对 w(n) 取微分得到梯度矩阵：

相应的，e(n) 的 Jacobi矩阵是：

因为误差向量对权值向量 w(n) 已经是线性的，Gauss-Newton法在一次迭代后收敛，将上式综合后代入可得：

项被看作是数据矩阵 X(n) 的伪逆，即：

因此，上式可写为：

读者注意！！

这个公式和我们在文章之前推导的最小二乘的几何意义得到的公式是一致的。通俗地说：

Gauss-Newton（以及其他迭代算法）的每一次迭代，本质上就是在这个 n 的时域内，进行最小二乘估计，并根据得到的本次最优解对权值向量进行更新。

这个公式表示了下面所陈述的一个简便途径：

权值向量 w(n+1) 求解定义在持续时间 n 的一个观察区间上的线性最小二乘问题，是如下两项的乘积：伪逆和期望的响应向量 d(n)

2. 遍历环境下线性最小二乘滤波器的极限形式

我们已经知道了，LMS算法在一次迭代中（时间 n 时域区间），本质上是在进行最小二乘估计。接下来继续思考，如果这个过程无限进行下去会得到什么呢？即 n 趋近于无穷。

令表示线性最小二乘滤波器关于观测数 n 的极限形式，允许 n 趋于无穷。

基于公式

得到：

现在假设输入向量 x(i) 和相应的期望响应 d(i) 来自于联合遍历。我们可以用时间均值来代替总体均值。

输入向量 x(i) 的相关矩阵（correlation matrix）的总体平均形式是：

并且，相应地，输入向量 x(i) 和期望响应 d(i) 之间的互相关系（cross-correlation vector）的总体平均形式是：

，其中 E 表示期望算子。

从而，在遍历假设下，有：和

综上，可将式：

改写为：

权值向量称为线性最优滤波问题的维纳解。

因此，我们可以做以下的陈述：对一个遍历过程，当观察样本数趋于无穷时，线性最小二乘滤波器渐进趋于维纳滤波器。

3. 维纳滤波器的局限性

虽然，当样本量趋近于无穷时，线性最小二乘滤波器趋近于维纳滤波器，但是设计维纳滤波器需要二阶统计量的知识：

1）输入向量 x(n) 的相关矩阵

2）x(n) 与期望响应 d(n) 的互相关向量

但是，在实际的情况下，这些信息都是未知的，所以维纳滤波器只是一个理论上的最优滤波器。

在实际工程实践中，我们可以利用线性自适应滤波器（linear adaptive filter）来处理未知的环境，自适应在这里的含义就是滤波器能够调整自己的自由参数来响应环境的统计变化。在连续的时间基础上做这类调整的一个流行的算法就是最小均方算法（LMS）。

接下来，我们进入对LMS的讨论。

0x2： LMS（Least mean square）一种线性自适应滤波算法

LMS最小均方算法是第一个解决如预测和信道均等化等问题的线性自适应滤波算法。

值得注意的是，LSM算法自身不仅可以作为自适应滤波应用机器，它还可以作为其他自适应滤波算法的评价准则，这里面的原因包括：

1. 从计算复杂度来说，对于可调参数而言 LMS算法的复杂度是线性的，这使得算法是计算高效的（computationlly efficient）的，而算法从性能上来说依然是有效的。
2. 算法可以简单地用代码来实现，因而是容易建立的。
3. 对于外部扰动来说，算法是鲁棒的，LMS算法是模型独立的。对任何模型，LMS都可以不辱使命，通过每次迭代的局部最优尽可能地去靠近全局最优。

对工程来说，上述性能都是非常重要的。之所以强调说工程，是因为其实LMS并不是理论上最优的算法，但是却是最实际工程有效的。

因为在实际情况中，我们很难获得全局最优解，甚至说都无法完整按照最速下降的思想进行最优方向的梯度下降，原因大致如下：

1. 神经网络的目标函数多数是非凸函数，不一定有全局最优解
2. 多维函数偏导为0并非是局部极值
3. 解析式法运算复杂，多维函数通过解析式方法求解极值不是一件容易的事
4. 在一个 n 时域区间内，梯度并不是保持不变的，最速下降要求同时考虑所有方向上的梯度，这样得到的梯度矩阵不一定是正定，可逆的。

但是LMS拥有计算简单、鲁棒性等优点，使得LMS在之后的深度学习/BP理论的发展中被不断继承和发扬光大。

0x3：LMS算法公式

最小均方（least mean square，LMS）算法的建立是基于极小化代价函数的瞬时值。注意！是瞬时值。

代价函数为：

这里 e(n) 是 n 时刻的瞬时误差信号。

把对权值向量 w 求微分得到：

和在最小二乘滤波器上一样，LMS算法运行在一个线性神经元上，可以把误差信号表示为：

因此，

综上公式得：

，代价函数对 w(n) 的导数。

上式作为梯度向量的一种瞬时估计（instantaneous estimate），可以记为：

最后，将上式梯度的瞬时估计公式，带入最速下降法作为最速下降法的梯度向量，可以得到LMS的算法公式：

这里值得注意的是：

学习率参数的倒数可以作为度量LMS算法的记忆（memory）：

给赋的值越小，LMS算法将记忆的过去的数据就越多。因此，值小的话，LMS算法执行得更精确，但算法的收敛速度慢

1. LMS用权重向量的瞬时估计代替实时变化

在LMS的公式中，我们用代替来强调这样一个事实：

利用最速下降法可以得到一个权值向量，而LMS算法产生该权值向量的一个瞬时估计。所以，利用LMS算法时我们牺牲掉最速下降法的一个明显特征。

在最速下降法中，对一个给定的，权值向量在权值空间中遵循一个明确定义的轨迹。对比之下，在LMS算法中，权值向量则跟踪一个随机的轨迹。

由于这个原因，LMS算法有时也被称为”随机梯度算法“。当LMS算法的迭代次数趋于无穷时，在维纳解周围随机行走（布朗运动）。

一个重要的事实是，与最速下降法不同，LMS算法不需要知道环境的统计特征。从实际的角度来看，LMS的这一特征是非常重要的。

0x4：LMS算法的信号流图表示

我们可以把LMS算法中的权值向量演变过程表示如下：

这里，I 是单位矩阵。通过运用LMS算法，我们认识到：

这里是单位时间延迟算子（unit-time delay operator），意味着存储。

我们利用信号流图来表示LMS算法，这图揭示了LMS算法是随机反馈系统的一个实例。反馈的出现对LMS算法的收敛行为有重要影响。

0x5：用马尔科夫模型来描述LMS算法和维纳滤波器的偏差

为了给LMS算法提供一种统计分析，我们利用下式定义的权值误差向量（weight-error vector）更加方便。

其中，是最优维纳解，是相应的由LMS算法计算的权值向量的估计。

因此，利用术语，假设其为一个状态（state），可以将式重写为紧凑形式：

这里，我们有

右边的噪音项为，其中，是由维纳滤波器产生的估计误差。

1. LMS算法的马尔科夫模型的特征

1. 模型的更新状态（updated state），由向量定义，依赖于老的状态，其自依赖性由转移矩阵（transition matrix） A(n) 定义。

2. 在时间 n 上状态的演化被内部所产生的噪音 f(n) 所扰动，这一噪声扮演者”驱动力“的角色。

下图给出了表示这一模型的向量值信号流图。标志为的分支表示模型的记忆，作为单位时间延迟算子，由下式定义：

和上面的原始形式相比，这个图中用紧凑形式重点强调了LMS算法中的反馈过程。

需要注意的是！

上面的信号流图以及相应的方程提供了在小学习率参数的假设下LMS算法的收敛性分析的框架。

下面我们开始讨论小学习率以及实现这一目标的两个基础知识：朗之万方程和Kushner直接平均法。

Relevant Link: 

https://www.cnblogs.com/xingshansi/p/6603263.html

 

4. LMS算法收敛的基础知识

0x1：朗之万方程 - 布朗运动的特点

考虑稳定性和收敛性，我们可以说LMS算法（对足够小的）从未达到完美的稳定或者收敛条件。

而且，在大量的迭代时间步 n 之后，算法达到”伪平衡“条件，这从定性上讲，可由算法围绕着维纳解执行布朗运动来描述。

这一类统计行为可通过非平衡热力学的朗之万方程（Langevin equation）来很好地解释。

由 v(t) 定义质量为 m 的宏观粒子陷入粘滞铃的速度。假设粒子足够小以使其由热起伏而来的速度被视为重要的（即小粒子的运动可以由热力学定律描述，因为质量较大的物体主要由牛顿力学和流体力学定律所描述）。

然后，由均分人力学原理（equipartition law of thermodynamics），粒子的平均能量由下式定义：

其中是Boltzmann常数，T为绝对温度（absolute temperature）。粘滞铃中分子作用于粒子的总驱动由两部分组成：

1. 根据Stoke定律（Stoke's lay）的等于- av(t) 的阻尼力（damping force），其中 a 是摩擦系数

2. 涨落力（fluctuating force），其性质是平均指定的。

粒子的运动方程在缺少外部驱动时由下式给定：

两边同除以 m，我们得：

其中，和

项是每单位质量涨落力（fluctuating force per unit mass）

因为它依赖于组成粒子的极为大量的原子数量，所以它是一个统计驱动力，它处于不规则运动状态。

公式被称为朗之万方程（Langevin equation），被称为朗之万力（Langevin force）。

朗之万方程描述了在粘滞铃中所有时间下粒子的运动（如果其初始条件是指定），它是描述非平衡热力学的第一个数学公式。
在后面的章节中，我们会证明LMS算法的一个变换版本具有和朗之万方程的离散时间版本相同的数学形式。但是在证明之前，需要先讨论一个基础知识。

0x2：Kushner直接平均法

公式的马尔科夫模型是非线性随机差分方程（nonlinerar stochastic difference equation）。

这一方程是非线性的原因是因为：

1. 转移矩阵 A(n) 依赖于输入向量 x(n) 的外积。因此，权值误差向量对于 x(n) 的依赖性和叠加原理相冲突，而这一原则是线性的需要。

2. 而且，方程是随机的是因为训练样本是从随机环境中取得的。

基于这2个原因，我们发现对LMS算法作严格的统计分析是困难的任务。

然后，在满足一定的条件下，LMS算法的统计分析能够通过将Kushner直接平均法（kushner direct-averaging method）应用于上式的模型从而显著地简化。

1. 对Kushner直接平均法的正规陈述

考虑由马尔科夫模型所描述的随机学习系统：

其中，对一些输入向量 x(n) ，我们有

而且附加噪声 f(n) 是由学习率参数线性拉伸的。有

1. 学习率参数是充分小的。

2. 附加噪声 f(n) 本质上独立于状态，修正马尔科夫模型的状态演化由下式来描述：

上式确定矩阵是修正马尔科夫模型的转移矩阵，我们用来表示修正马尔科夫模型的状态，用来强调这一模型随时间的演化仅在微笑的学习率参数的有限情况下等同于原始马尔科夫小莫名。

上述的讨论可以充分说明：

1. 学习率参数小的时候，LMS算法具有长记忆（long memory）。因此，更新状态的演化可以通过时间一步一步追踪所有的路径直到初始条件。

2. 当小的时候，可以在的展开式序列中忽略二阶和高阶项。

3. 最后，在上式修正马尔科夫模型中，通过调用遍历性原理，此时总体平均可以被时间平均所替。

0x3：小学习率参数下统计LMS学习理论

我们已经有了Kushner直接平均法，现在可以建立LMS算法的统计分析原则了。我们做三个合理的假设：

假设一：学习率参数是小的

从实际角度，的小的选择也是有意义的。特别是，当小时，LMS算法对于外部扰动是鲁棒的

假设二：维纳滤波器产生的估计误差是白噪

假设三：输入向量 x(n) 和期望响应 d(n) 是联合高斯分布

由物理现象产生的随机过程频繁地出现使得高斯模型是适当的。

更重要的是，朗之万力对于LMS算法的非平衡行为负责，它证明了自身在大量的足够多的迭代次数 n 的时候，算法在最优维纳解的周围进行布朗运动（Brownian motion）。当然，前提是基于学习率参数的情况下。

Relevant Link:

http://chem.xmu.edu.cn/teach/chemistry-net-teaching/wuhua/chapter3/part12/12-1.html

 

posted @ 2018-11-06 14:49  郑瀚Andrew  阅读(1557)  评论(0)  编辑  收藏  举报