A Chebyshev-Gauss Pseudospectral Method for Solving Optimal Control Problems
-
摘要: 提出了一种求解最优控制问题的Chebyshev-Gauss伪谱法, 配点选择为Chebyshev-Gauss点. 通过比较非线性规划问题的Kaursh-Kuhn-Tucker条件和伪谱离散化的最优性条件, 导出了协态和Lagrange乘子的估计公式. 在状态逼近中, 采用了重心Lagrange插值公式, 并提出了一种简单有效的计算状态伪谱微分矩阵的方法. 该法的独特优势是具有良好的数值稳定性和计算效率. 仿真结果表明, 该法能够高精度地求解带有约束的复杂最优控制问题.
-
关键词:
- 最优控制 /
- 伪谱法 /
- 协态估计 /
- Chebyshev-Gauss点
Abstract: A pseudospectral method is presented for direct trajectory optimization of optimal control problems using collocation at Chebyshev-Gauss points, and therefore, it is called Chebyshev-Gauss pseudospectral method. The costate and constraint multiplier estimates for the proposed method are rigorously derived by comparing the discretized optimality conditions of an optimal control problem with the Karush-Kuhn-Tucker conditions of the resulting nonlinear programming problem from collocation. The distinctive advantages of the proposed method over other pseudopsectral methods are the good numerical stability and computational efficiency. In order to achieve this goal, the barycentric Lagrange interpolation is substituted for the classic Lagrange interpolation in the state approximation. Furthermore, a simple yet efficient method is presented to alleviate the numerical errors of state differential matrix using the trigonometric identity especially when the number of Chebyshev-Gauss points is large. The method presented in this paper has been taken to two optimal control problems from the open literature, and the results have indicated its ability to obtain accurate solutions to complex constrained optimal control problems.-
Key words:
- Optimal control /
- pseudospectral methods /
- costate estimation /
- Chebyshev-Gauss points
-
自2006年Hinton等[1]提出第一个深度置信网络开始,经过十年的发展,深度学习已逐渐成为机器学习研究领域的前沿热点.深度置信网络[2]、深度卷积神经网络[3]、深度自动编码器[4]等深度网络也广泛应用于机器学习的各个领域,如图像识别、语音分析、文本分析等[5-7].相对于传统的机器学习网络,深度网络取得了更好的效果,极大地推动了技术发展水平(State-of-the-art)[8].尤其在大数据背景下,针对海量无标签数据的学习,深度网络具有明显的优势[9].
受限玻尔兹曼机(Restricted Boltzmann machine,RBM)[10]是深度学习领域中的一个重要模型,也是构成诸多深度网络的基本单元之一.由于RBM较难训练,所以在很多大数据量任务上使用较少.但相对于其他基本模型,RBM具备较强的理论分析优势和可解释性,是帮助我们理解深度网络和其他基本模型内在机理的重要模型,而且在某些特殊数据集上,RBM可以获得更好的学习效果. 所以,研究RBM仍然很有意义. RBM具有两层结构,在无监督学习下,隐层单元可以对输入层单元进行抽象,提取输入层数据的抽象特征.当多个RBM或RBM与其他基本单元以堆栈的方式构成深度网络时,RBM隐层单元提取到的抽象特征可以作为其他单元的输入,继续进行特征提取.通过这种方式,深度网络可以提取到抽象度非常高的数据特征.当采用逐层贪婪(Greedy layer-wise)[11]训练方法对深度网络进行训练时,各个基本单元是逐一被训练的.因此,RBM训练的优劣将直接影响整个深度网络的性能.
2006年,Hinton等提出了对比散度(Contrastive divergence,CD)算法[12]用以训练RBM网络.在每次训练迭代时,CD算法以数据样本为初始值,通过多步Gibbs迭代获得目标分布的近似采样,然后通过该近似采样来近似目标梯度,取得了较好的效果,是目前RBM训练的标准算法.但研究表明,CD算法对目标梯度的估计是有偏估计[13],而且每次迭代时都需要重新启动Gibbs采样链,这降低了CD算法的训练性能.为此,Tieleman等以CD算法为基础,于2008年提出了持续对比散度(Persistent contrastive divergence,PCD)算法[14].在学习率足够小的前提下,每次参数更新后,RBM模型的变化不大,可以认为RBM网络分布基本不变.基于此假设,PCD算法只运行一条独立的采样链,以上次采样迭代的采样值作为下次采样迭代的初值继续迭代,而不是像CD算法那样每次采样都以样本数据为采样初值,取得了比CD算法更好的训练效果.为了加速PCD算法,Tieleman又于2009年提出了加速持续对比散度(Fast persistent contrastive divergence,FPCD)算法[15],引入了额外的加速参数来提高训练速度. PCD算法和FPCD算法虽然训练性能较CD算法有所提高,但并没有从本质上提高CD算法的混合率[16].不管是CD算法,还是以CD算法为基础的PCD算法、FPCD算法,都是通过一条Gibbs采样链来逼近目标分布,对于目标分布较简单的数据,可以取得较好的效果.但当数据分布复杂,尤其为多模分布时,即目标分布函数存在多个峰值,Gibbs采样链很容易陷入局部极小域,导致样本不能描述数据分布的整体结构[17].为克服这个问题,Desjardins (2010)等[18]、Cho (2010)等[19]、 Brakel (2012)等[20]等分别提出应用并行回火算法(Parallel tempering,PT)来训练RBM. PT算法并行化多条温度链,每条温度链上进行多步Gibbs迭代.高温链采样目标总体分布的结构信息,低温链采样目标局部分布的精确信息.不同温度链之间以一定的交换概率进行交换,不断迭代,最后低温链就可以精确获得目标分布的总体信息.对于多模分布数据,PT算法的训练效果要明显优于CD算法[21].
通过以上描述可知,不管是CD算法还是PT算法,本质上都是以Gibbs采样来获得关于目标分布的采样样本.因此,Gibbs采样性能的优劣将直接影响以上算法的训练效果.本文研究发现,当采用多步Gibbs采样时,在训练初期会发生采样发散现象,严重影响网络收敛速度,而且算法运行速度较慢;当采用单步Gibbs采样时,前期网络收敛性质较好,且算法运行速度较快,但后期采样精度不高.如何在前期保证良好的收敛性质,同时在后期保证网络训练精度并提高算法运行速度,是目前基于Gibbs采样的RBM训练算法亟需解决的问题,从现有文献来看,尚无人对以上问题进行研究.因此,本文将从马尔科夫采样理论的角度对以上问题进行分析,并提出了动态Gibbs采样算法,最后给出了仿真验证.
1. 问题描述
受限玻尔兹曼机是一个马尔科夫随机场模型[22],它具有两层结构,如图 1所示.下层为输入层,包含m个输入单元vi ,用来表示输入数据,每个输入单元包含一个实值偏置量ai ;上层为隐层,包含n个隐层单元hj ,表示受限玻尔兹曼机提取到的输入数据的特征,每个隐层单元包含一个实值偏置bj .受限玻尔兹曼机具有层内无连接,层间全连接的特点.即同层内各节点之间没有连线,每个节点与相邻层所有节点全连接,连线上有实值权重矩阵wij .这一性质保证了各层之间的条件独立性.
本文研究二值受限玻尔兹曼机[23],即随机变量(V,H)取值(v,h)∈{0,1}.由二值受限玻尔兹曼机定义的联合分布满足Gibbs分布 $P(v,h)=\frac{1}{{{Z}_{\theta }}}{{e}^{-{{E}_{\theta }}(v,h)}}$ ,其中θ为网络参数 $\theta =\{{{a}_{i}},{{b}_{j}},{{w}_{ij}}\}$ ,Eθ (v,h)为网络的能量函数:
${{E}_{\theta }}(v,h)=-\sum\limits_{i=1}^{n}{\sum\limits_{j=1}^{m}{{{w}_{ij}}{{v}_{i}}{{h}_{j}}}-\sum\limits_{i=1}^{m}{{{a}_{i}}{{v}_{i}}}-\sum\limits_{j=1}^{n}{{{b}_{j}}{{h}_{j}}}}$
(1) Zθ为配分函数: ${{Z}_{\theta }}=\sum\limits_{v,h}{{{e}^{-{{E}_{\theta }}(v,h)}}}$ .输入层节点v的概率分布P(v)为: $P(v)=\frac{1}{{{Z}_{\theta }}}\sum\limits_{h}{{{e}^{-{{E}_{\theta }}(v,h)}}}$ .由受限玻尔兹曼机各层之间的条件独立性可知,当给定输入层数据时,输出层节点取值满足如下条件概率:
$\begin{array}{*{35}{r}} P({{h}_{k}}=1|v)=\frac{1}{1+\exp (-{{b}_{j}}-\sum\limits_{i=1}^{n}{{{w}_{ij}}{{v}_{i}}})}= \\ sigmoid\left( {{b}_{j}}+\sum\limits_{i=1}^{n}{{{w}_{ij}}{{v}_{i}}} \right)\ \\ \end{array}$
(2) 相应地,当输出层数据确定后,输入层节点取值的条件概率为
$\begin{array}{*{35}{r}} P({{h}_{k}}=1|v)=\frac{1}{1+\exp (-{{a}_{i}}-\sum\limits_{i=1}^{n}{{{w}_{ij}}{{h}_{j}}})}= \\ sigmoid\left( {{a}_{i}}+\sum\limits_{i=1}^{n}{{{w}_{ij}}{{h}_{j}}} \right)\ \\ \end{array}$
(3) 给定一组训练样本 $S=\{{{v}^{1}},{{v}^{2}},\cdots ,{{v}^{n}}\}$ ,训练RBM意味着调整参数θ,以拟合给定的训练样本,使得该参数下由相应RBM表示的概率分布尽可能地与训练数据的经验分布相符合.本文应用最大似然估计的方法对网络参数进行估计.这样,训练RBM的目标就是最大化网络的似然函数: ${{L}_{\theta ,w}}=\prod\limits_{i=1}^{n}{P({{v}^{i}})}$ .为简化计算,将其改写为对数形式: $\ln {{L}_{\theta ,w}}=\sum\limits_{i=1}^{n}{\ln P({{v}^{i}})}$ .进一步推导对数似然函数的参数梯度
$\begin{align} & \frac{\partial \ln P(v)}{\partial {{a}_{i}}}=-\sum\limits_{h}{P(h|v)}\frac{\partial E(v,h)}{\partial {{a}_{i}}}+ \\ & \sum P(v,h)\frac{\partial E(v,h)}{\partial {{a}_{i}}}={{v}_{i}}-\sum\limits_{v}{P(v){{v}_{i}}} \\ & \ \frac{\partial \ln P(v)}{\partial {{b}_{j}}}=-\sum\limits_{h}{P(h|v)}\frac{\partial E(v,h)}{\partial {{b}_{j}}}+ \\ & \sum P(v,h)\frac{\partial E(v,h)}{\partial {{b}_{j}}}=P({{h}_{i}}=1|v)-\ \\ & \sum\limits_{v}{P(v)P({{h}_{i=1}}|v)} \\ & \ \frac{\partial \ln P(v)}{\partial {{w}_{ij}}}=-\sum\limits_{h}{P(h|v)}\frac{\partial E(v,h)}{\partial {{w}_{ij}}}+ \\ & \sum P(v,h)\frac{\partial E(v,h)}{\partial {{w}_{ij}}}=P({{h}_{j}}=1|v){{v}_{i}}-\ \\ & \sum\limits_{v}{P(v)P({{h}_{j=1}}|v){{v}_{i}}}\ \\ \end{align}$
(4) 得到对数似然函数的参数梯度后,可以由梯度上升法求解其最大值.但由于数据分布P(v)未知,且包含配分函数Zθ ,因此,无法给出梯度的解析解.现有训练算法主要是基于采样的方法,首先,构造以P(v) 为平稳分布的马尔科夫链,获得满足P(v)分布的样本;然后,通过蒙特卡洛迭代来近似梯度:
$\begin{matrix} \begin{align} & \nabla {{a}_{i}}=v_{i}^{(0)}-v_{i}^{(k)} \\ & \ \nabla {{b}_{j}}=P({{h}_{j}}=1|{{v}^{(0)}})-P({{h}_{j}}=1|{{v}^{(k)}})\ \\ & \nabla {{w}_{ij}}=P({{h}_{j}}=1|{{v}^{(0)}})v_{i}^{(0)}-P({{h}_{j}}=1|{{v}^{(k)}})v_{i}^{(k)}\ \\ \end{align} \\ \end{matrix}$
(5) 其中,vi(0) 为样本值,vi(k)为通过采样获得的满足P(v)分布的样本.最后,参数更新方程如下:
$\begin{matrix} \begin{align} & {{a}_{i}}={{a}_{i}}+\nabla {{a}_{i}} \\ & \ {{b}_{i}}={{b}_{i}}+\nabla {{b}_{i}}\ \\ & {{w}_{ij}}={{w}_{ij}}+\nabla {{w}_{ij}}\ \\ \end{align} \\ \end{matrix}$
(6) 现有RBM训练算法,包括CD_k算法、并行回火(PT)算法,这两类算法都是以Gibbs采样为基础的,都是通过多步Gibbs采样获得一定精度的目标采样,然后分别通过其他后续操作获得最终的目标梯度近似值. CD_k算法是RBM训练的主流算法,因此,本节以CD_k算法为例,通过仿真的方式,揭示了作为以上算法基本操作单元的Gibbs采样在网络训练过程中出现的问题,研究了它对网络收敛速度和训练精度的影响.
首先给出 CD_k 算法的操作步骤:
步骤 1. 设定网络参数初值.
步骤 2.将训练数据输入到输入层节点,由式(2)对隐层节点值进行采样,
步骤 3.根据式(3)对输入层节点进行采样.再以此采样值作为输入层节点的值重复步骤2,这样就完成了一步Gibbs采样.
步骤 4. 步骤2和步骤3重复k次,完成k步Gibbs采样,即CD_k.
步骤 5.将步骤4获得的采样值带入式(5)中,计算参数梯度.
步骤 6. 将步骤5中获得的参数梯度带入式(6)中,对参数进行更新.
步骤 7.更新训练数据,重复步骤 2～6,直到达到额定迭代次数.
相应的伪代码如算法1所示:
算法 1. CD_k 算法伪代码
Input: RBM(V1,…,Vm,H1,…,H_n ),training batch S
Output: wij,aj and bi for i=1,…,n,j=1,…,m
1: Init ▽wij =▽aj=▽bi=0 for i=1,…,n,j=1,…,m
2: For all the v∈S do
3: v(0)←v
4: for t=0,…,k-1 do
5: for i=1,…,n do sample hi(t)～ p(hi|v(t))
6: for j=1,…,m do sample vj(t+1)～ p(vj|h(t))
7: for i=1,…,n,j=1,…,m do
8: ▽ wij =p(Hi =1|v(0))· vj (0)-p(Hi =1|v(k)· vj(k)
9: ▽aj=vj(0) -vj(k)
10: ▽bi=p(Hi =1|v(0))-p(Hi =0|v(k))
11: wij =wij +η▽wij
12: ai =ai+η▽ai
13: bi=bi+η▽bi
14: End For
其中,a为可见层偏置向量,b为隐层偏置向量,w为网络权值矩阵,η为学习率.
1.1 问题实验描述
1) 实验设计
本文采用的数据集是MNIST数据集,它是二值手写数据集,也是目前训练RBM网络的标准数据集.它总共包含60 000个训练样本和10 000个测试样本,每个样本是一幅28像素×28像素的灰度图.所采用的RBM网络有784×500个节点,输入层有784个可见单元,对应灰度图的784个像素点;输出层有500个隐层节点,这是目前实验显示的训练效果较好的隐层节点数目.具体的网络参数初始值设定如表 1.
表 1 网络参数初值Table 1 Initial value of parameters网络参数 初始值 a zeros(1, 784) b zeros(1, 500) w 0.1 × randn(784, 500) η 0.1 本文设计了6组对比实验,用60000个训练样本对RBM进行训练,分别迭代1000次,如表 2所示.其中CD_k表示进行k步Gibbs迭代.用于显示的样本数据的原始图片如图 2所示.实验结束后,我们比较了各组实验的重构误差,并给出了最终的误差图.
表 2 实验分组Table 2 Experimental grouping数据集 算法 迭代次数 MNIST CD_1 1000 MNIST CD_5 1000 MNIST CD_10 1000 MNIST CD_100 1000 MNIST CD_500 1000 MNIST CD_1000 1000 2) 仿真结果图 3表示整个迭代过程中各组CD算法的重构误差图,图 4给出了各组实验的训练时间,图 5～图 10分别给出了各组实验的采样灰度图.
1.2 问题归纳描述
上节实验给出了CD算法在不同Gibbs采样步数下的仿真图,可以看出,当RBM网络采用多步Gibbs算法进行采样迭代时,会出现如下问题:
问题 1.训练初始阶段,得到的每幅重构采样图几乎完全相同.
如图 11、图 12所示,在训练初始阶段,多步Gibbs采样出现了各组采样数据同分布的现象,这表明各组样本几乎完全相同,这与事实相左.在训练初期,大约0～100次迭代期间,这种现象持续存在.
问题 2. 采样误差分布集中,在批量训练时,存在全0全1现象.
如图 13、图 14所示,当进行多步Gibbs采样时,出现了误差分布集中的现象:有些样本采样几乎全为1,而其他的样本采样几乎全为0.由仿真实验可知,在0～100次迭代期间,这种现象在迭代初期持续存在.
问题3.一步Gibbs采样初始误差小,训练速度快,但后期训练精度低;多步Gibbs采样初始误差大,训练速度慢,但后期训练精高.
如图 15、图 16所示,只进行一步Gibbs采样的CD_1算法在开始时训练误差较小,很快便收敛到较好值,但训练后期精度不如CD_10等进行多步Gibbs迭代的CD算法;进行多步Gibbs采样的CD_k迭代算法,在训练初期误差较大,且不断振荡,而且训练时间较慢,但到训练后期,它们可以达到极高的精度.
以上实验表明,CD算法虽然对RBM具有良好的训练能力,但Gibbs采样的步数对训练性能造成了明显的影响.我们将在下节研究这种影响,并对以上问题给出理论分析.
图 15 采样误差局部放大图Fig. 15 Local enlarged drawing of reconstruction error in initial phase
图 16 CD_1000 采样误差局部放大图Fig. 16 Local enlarged drawing of reconstruction error in later stage2. Gibbs 采样误差的理论分析
Gibbs采样是马尔科夫链蒙特卡洛(Markov chain Monte Carlo,MCMC)采样算法的一种.在RBM训练中,它的转移核是Sigmoid函数.隐层节点和输入层节点交替采样,公式如下:
$\begin{array}{*{35}{l}} \begin{align} & P({{h}_{j}}=1|V)=sigmoid({{b}_{j}}+\sum\limits_{i=1}^{n}{{{w}_{i,j}}{{v}_{i}}})\ \\ & P({{v}_{i}}=1|H)=sigmoid({{a}_{i}}+\sum\limits_{i=1}^{n}{{{w}_{i,j}}{{h}_{j}}})\ \\ \end{align} \\ \end{array}$
(7) 由马尔科夫链收敛定理可知,当n→ +∞ 时,Gibbs采样链会收敛到平衡分布,即:
${{\pi }_{i}}(x)={{\pi }_{i-1}}(x)P={{\pi }_{0}}{{P}^{n}}$
(8) 其中,π(x)为样本x的平衡分布.同时,由细致平衡准则可得:
$\pi (i){{P}_{ij}}=\pi (j){{P}_{ji}},\ \ \forall i,j$
(9) 即Gibbs采样的平稳分布与迭代初始值无关,只与转移概率有关.由上面给出的RBM交替采样概率公式可知,当用Gibbs采样对RBM进行采样训练时,其平稳分布是网络参数的函数:
$\pi (x)=f(a,b,w)$
(10) 从这个角度讲,训练RBM的目的就是调节网络参数,使由网络参数确定的平稳分布等于样本的真实分布.
基于以上描述,下面对第2节中提出的问题给出理论解释.
问题1. 训练初始阶段,得到的每幅重构采样图几乎完全相同.
初始时刻,网络参数初值相同,在早期迭代过程中,网络参数值的变动也不大,满足如下公式:
$\left\{ \begin{array}{*{35}{l}} \begin{align} & {{a}_{i}}-{{a}_{j}}﹤\varepsilon \\ & \ {{b}_{i}}-{{b}_{j}}﹤\varepsilon \ \\ & {{w}_{i}}={{w}_{j}}\ \\ \end{align} \\ \end{array} \right.\ \ \ \ $
(11) ε为一极小正值.由网络参数决定的平稳分布也近乎相同:
$\begin{align} & f({{a}_{i}},{{b}_{i}},w)\approx f({{a}_{j}},{{b}_{j}},w)\Rightarrow \\ & \ \pi ({{x}_{i}})\approx \pi ({{x}_{j}})\ \\ \end{align}$
(12) 即各样本的平稳分布相等.因此,当进行多步Gibbs采样时,各训练样本的采样样本逐渐收敛到相同的平稳分布,这时就出现了问题1描述的现象,各样本的重构采样图几乎完全相同.
问题2. 采样误差分布集中,在批量训练时,存在全0全1现象.
由上一部分分析可知,在训练初期,网络参数改变不大,由RBM参数决定的平衡分布几乎同构,即各采样概率收敛到相同平衡分布值.上述对比实验中,网络参数的初始值为θ =(a,b,w)=(0,0,0.1),此时网络平衡分布收敛在0.5附近,样本数据的收敛概率将在0.5附近浮动,即一部分样本的采样概率略小于0.5,另一部分样本的采样概率略大于0.5,即:
$\begin{align} & \pi (\theta )\to 0.5\ \\ & p({{v}_{i}}|H)=0.5+\varepsilon \\ & \ p({{v}_{n-i}}|H)=0.5-\varepsilon \ \\ \end{align}$
(13) 其中,ε为一极小正值.这时基于随机数对样本进行采样,一部分样本的采样值将全为0,另一部分的采样值将全为1,即全0全1现象.
问题3.一步Gibbs采样初始误差小,训练速度快,但后期训练精度低;多步Gibbs采样初始误差大,训练速度慢,但后期训练精高.
设网络参数期望值为 $\hat{\theta }=\left( \hat{a},\hat{b},\hat{w} \right)$ ,它代表参数的真实值:设网络参数实际值为θ =(a,b,w),这是我们在训练网络过程中,网络参数的实际值,训练的目标就是使网络参数实际值逐渐逼近其真实值.定义网络参数差△θ(△a,△b,△w):
$\begin{array}{*{35}{l}} \begin{align} & \Delta a=\hat{a}-a\ \\ & \Delta b=\hat{b}-b\ \\ & \Delta w=\hat{w}-w\ \\ \end{align} \\ \end{array}$
(14) 在网络训练早期,网络参数差较大,由网络参数定义的平稳分布与真实分布相差也较大,即 $\Delta \pi =\left| {{\pi }_{{\hat{\theta }}}}(x)-{{\pi }_{\theta }}(x) \right|\gg 0$ .此时,如果对样本进行多步迭代采样,采样样本将偏离真实分布,从而不能收敛到真实分布,而是收敛到与真实分布相差较大的其他分布.因此,在迭代初期,CD_1000、CD_500等算法的采样误差非常大,而且运行时间较长.而CD_1算法由于只进行了一次采样迭代,不仅运行速度加快,而且由于采样样本的分布没有偏离真实分布太多,使得这时候的CD_1算法的采样误差非常小.由实验可知,此时采样误差的大小关系为:CD_1﹤CD_5﹤CD_10﹤CD_100﹤CD_500﹤CD_1000.到了网络训练后期,由于网络参数差非常小,网络参数的实际值已经非常接近真实值,这时候进行多步Gibbs迭代能很好地逼近样本真实分布,所以这一阶段,CD_k算法的采样精度要比CD_1高.但由于网络参数差一直存在,所以,Gibbs迭代步数也不宜过高,如实验所示,CD_1000在采样到最后,采样误差仍高于CD_10.
3. 动态Gibbs采样
在现有以Gibbs采样为基础的RBM训练算法中,Gibbs采样的采样步数多为固定值,即在整个训练过程中,每次迭代采样时都进行固定步数的Gibbs采样,这样就难以兼顾训练精度和训练速度这两个训练指标.当进行多步Gibbs采样时,容易在训练前期发生误差发散的现象,且算法运行时间较长;一步Gibbs采样算法运行较快,但后期训练精度不高,基于此,本文提出了动态Gibbs采样(Dynamic Gibbs sampling,DGS)算法.
定义1.动态Gibbs采样是指在迭代训练过程中的不同阶段,根据网络的训练误差,动态地调整Gibbs采样的步数,以达到最优训练效果.
通过上节分析可知,在网络训练初期,网络参数几乎相等,各样本的平稳分布也近乎相等,而且网络参数差较大,样本的平稳分布与真实分布相差也较大,因此,这一阶段应尽量减少采样次数,克服多步Gibbs采样引起的误差发散,提高训练速度,使网络参数尽快逼近真实值;当网络参数逼近真实值时,此时应加大采样迭代次数,提高训练精度.
基于以上定义和描述,DGS算法的操作步骤如下:
步骤 1. 设定网络参数初值和动态策略M.
步骤 2. 在1～m1 迭代范围内,设置Gibbs采样步数k1 =Gibbs_N1 .
步骤 3. 将训练数据输入到输入层节点,由式(2)对隐层节点值进行采样.
步骤 4. 根据式(3)对输入层节点进行采样.再以此采样值作为输入层节点的值重复步骤3,这样就完成了一步Gibbs采样.
步骤 5. 步骤3和步骤4重复k1 次,完成k1步Gibbs采样.
步骤 6. 将步骤5获得的采样值带入式(5)中,计算参数梯度.
步骤 7. 将步骤6中获得的参数梯度带入式(6)中,对参数进行更新.
步骤 8. 更新训练数据,重复步骤3到步骤7,直到迭代次数达到m1 .
步骤 9. 在m1～m2 迭代范围内,设置Gibbs采样步数k2 =Gibbs_N2 .
步骤 10. 重复步骤3到步骤8,直到迭代次数达到m2 .
步骤 11. 在m2～Iter迭代范围内,设置Gibbs采样步数k3 =Gibbs_N3 .
步骤 12. 重复步骤3到步骤8,直到迭代次数达到最大迭代次数Iter.
相应的伪代码如算法2所示.
算法2. DGS 算法伪代码
1: Input: RBM(v1 ,v2 ,v3 ,…,vn ,h1 ,h2 ,h3 ,…,hm ),training batch S
2: Output: wij ,ai andbj for i=1,…,n,j=1,…,m
3: Init:▽wij =▽aj =▽bi =0 for i=1,…,n,j=1,…,m
4: For all the S do
5: for iter=1:m1 do
6: for t=0,…,k-1 do Gibbs_N1
7: for i=1,…,n do sample hi(t)～p(hi|v(t))
8: for j=1,…,m do sample vj(t+1)～p(vj|h(t))
9: for iter=m1 :m2 do
10: for t=0,…,k-1 do Gibbs_N2
11: for i=1,…,n do sample hi(t)～p(hi|v(t))
12: for j=1,…,m do sample vj(t+1)～p(vj|h(t))
13: for iter=m2:Iter do
14: for t=0,…,k-1 do Gibbs_N3
15: for i=1,…,n do sample hi(t)～p(hi|v(t))
16: for j=1,…,m do sample vj(t+1)～p(vj|h(t))
17: for i=1,…,n,j=1,…,m do
18:▽wij =p(Hi =1| v(0))· vj (0)-p(Hi=1|v(k))· vj(k)
19:▽aj =vj(0) -vj(k)
20:▽bi=p(Hi=1|v(0))-p(Hi=0|v(k))
21: wij =wij +η▽wij
22: ai =ai +η▽ai
23:bi=bi+η▽bi
24: End For}
其中,M=(m1 ,m2 )为动态策略,且满足m2 >m1 . Iter为总的迭代次数,iter为当前迭代次数. Gibbs_Ni为Gibbs采样,Ni 表示采样次数,且满足Nn >N_n-1 .其中Gibbs采样次数N与网络训练迭代次数M之间的大致关系如下:
$\begin{align} & Gibbs\_{{N}_{1}}=1\ \ \ \ \ \ \ \ \ \ 若\ iter\in (1\tilde{\ }{{m}_{1}}>) \\ & \ Gibbs\_{{N}_{2}}=2\tilde{\ }10\ \ \ 若\ iter\in ({{m}_{1}}\tilde{\ }{{m}_{2}})\ \\ & Gibbs\_{{N}_{3}}>10\ \ \ \ \ \ \ \ \ 若\ iter\in ({{m}_{2}}\tilde{\ }Iter) \\ & \\ \end{align}$
(15) 4. 仿真实验
本节设计了7组对比实验,第1～6组实验采用固定Gibbs采样步数的CD_k算法进行训练仿真,第6组实验用DGS算法对网络进行训练仿真,如表 3所示.两组实验使用相同的数据集MNIST,网络结构相同,网络参数初始值相同,如表 4所示.本文设计的动态采样策略如表 5所示.下面给出仿真实验结果和分析.
表 3 实验分组Table 3 Experimental grouping数据集 训练算法 Iter MNIST CD_1 1000 MNIST CD_5 1000 MNIST CD_10 1000 MNIST CD_100 1000 MNIST CD_500 1000 MNIST CD_1000 1000 MNIST DGS 1000 表 4 网络参数初值Table 4 Initial values of parameters算法参数 CD_k DGS a zeros(1,784) zeros(1,784) b zeros(1,500) zeros(1,500) w 0.1 × randn(784,500) 0.1 × randn(784,500) η 0.1 0.1 V 784 784 H 500 500 表 5 DGS 迭代策略Table 5 Iterative strategy of DGSM Gibbs_N (1:m1) = (1:300) Gibbs_N1 = 1 (m1:m2) = (300:900) Gibbs_N2 = 5 (m2:Iter) = (900:1000) Gibbs_N3 = 10 4.1 重构误差对比分析
图 17给出了所有算法的重构误差对比图.对比结果显示,本文设计的DGS算法可以很好地训练RBM网络,从而证明了本文算法的有效性.
在迭代初期,DGS算法只进行一次Gibbs采样迭代,避免了采样发散,从而迅速收敛到较好的值,由误差对比图初始阶段的局部放大图(图 18)可以看出,此时误差满足:
图 18 训练初期局部放大图Fig. 18 Local enlarged drawing of reconstruction error in initial phase$\begin{align} & \text{DGS= CD }\!\!\_\!\!\text{ 1}\succ \text{CD5}\succ \text{CD }\!\!\_\!\!\text{ 10}\succ \text{CD }\!\!\_\!\!\text{ 100}\succ \\ & \text{CD }\!\!\_\!\!\text{ 500}\succ \text{CD }\!\!\_\!\!\text{ 1000} \\ & \\ \end{align}$
(16) 在迭代后期,网络参数值已非常接近真实值,此时DGS逐步增大了Gibbs采样的迭代步数,获得了采样精度更高的目标样本,最终获得了更高的训练精度,即:
$\begin{align} & \text{DGS= CD }\!\!\_\!\!\text{ 10}\succ \text{CD }\!\!\_\!\!\text{ 5}\succ \text{CD }\!\!\_\!\!\text{ 1}\succ \text{CD }\!\!\_\!\!\text{ 100}\succ \\ & \text{CD }\!\!\_\!\!\text{ 500}\succ \text{CD }\!\!\_\!\!\text{ 1000} \\ & \\ \end{align}$
(17) 如图 19所示.
4.2 运行时间对比分析
图 20给出了所有算法的运行时间对比图.从图中可以看出,在整个训练过程中,DGS算法、CD_1算法、CD_5算法和CD_10算法的运行速度都明显比其他算法快.因此,下面根据本文设计的动态策略,对各个迭代区间内这4种算法的运行速度进行分析:
在1～300迭代范围内,DGS算法的Gibbs采样步数k设为1,与CD_1算法相同.所以,此时的DGS算法的运行速度与CD_1相同,且快于其他两种算法,如图 21所示.
在300～900迭代范围内,DGS算法的Gibbs采样步数k设为5.由图 22可以看出,此时DGS算法的运行速度逐渐放缓,运行时间明显上升,逐渐大于CD_1算法.
在900～1000迭代范围内,DGS算法的Gibbs采样步数k设为10.所以,这个时期的DGS运行时间持续放缓.但从图 23中可以看出,即便到了训练后期,DGS算法的运行时间仍然小于CD_5算法和其他CD_k (k>5)算法.这说明,DGS算法在后期提高训练精度的同时,只付出了微小的时间代价.
4.3 采样效果图
图 24～图 28分别给出了DGS算法在不同迭代次数下的采样重构图.对比图 11、图 12,可以看出,DGS在训练迭代50次以内就可以很好地重构输入样本,而且没有出现全0全1现象和采样图同构现象,从而克服了第2.2节问题1和问题2中描述的问题.
图 29显示了DGS训练结束后的重构灰度图,图中几乎没有噪点.可见,采用DGS算法训练网络可以获得更高的训练精度,从而解决了第2.2节中问题3描述的问题.
综上所述,本文设计的DGS算法在训练初期克服了多步Gibbs采样发散的缺点,在训练后期获得更高的精度,而且在保证收敛精度的情况下大幅度提高了训练速度,获得了较好的效果.
5. 总结
本文首先通过仿真实验,给出了现有基于Gibbs采样的RBM训练算法在训练初期误差发散和后期训练精度不高等问题的具体描述,然后从马尔科夫采样理论的角度对Gibbs采样误差进行理论分析.证明在RBM网络下,多步Gibbs采样较差的收敛性质是导致前期采样发散和算法运行速度较低的主要原因;单步Gibbs采样是造成后期训练精度不高的主要原因.基于此,本文提出了动态Gibbs采样算法,并给出了验证实验.实验表明,本文提出的动态Gibbs采样算法在训练初期克服了多步Gibbs采样引起的误差发散,后期克服了单步Gibbs采样带来的训练精度低的问题,同时提高了训练速度,以上特点可以弥补现有以Gibbs采样为基础的RBM训练算法的不足.
关于Gibbs采样步数、训练迭代次数与训练精度之间的关系,本文在理论分析部分只给出了定性分析;在动态Gibbs采样算法设计阶段,本文只是根据实验分析,给出Gibbs采样步数和训练迭代次数之间的经验区间. Gibbs采样步数、训练迭代次数以及网络训练精度之间是否存在精确的数学关系,如果存在,其数学模型如何构建.以上问题仍有待进一步研究.
-
[1] Benson D A. A Gauss Pseudospectral Transcription for Optimal Control [Ph.D. dissertation], Massachusetts Institute of Technology, USA, 2004 [2] [2] Benson D A, Huntington G T, Thorvaldsen T P, Rao A V. Direct trajectory optimization and costate estimation via an orthogonal collocation method. Journal of Guidance, Control, and Dynamics, 2006, 29(6): 1435-1440 [3] [3] Huntington G T. Advancement and Analysis of a Gauss Pseudospectral Transcription for Optimal Control Problems [Ph.D. dissertation], Massachusetts Institute of Technology, USA, 2007 [4] [4] Garg D, Patterson M A, Francolin C, Darby C L, Huntington G T, Hager W W, Rao A V. Direct trajectory optimization and costate estimation of finite-horizon and infinite-horizon optimal control problems using a Radau pseudospectral method. Computational Optimization and Applications, 2011, 49(2): 335-358 [5] [5] Garg D, Hager W W, Rao A V. Pseudospectral methods for solving infinite-horizon optimal control problems. Automatica, 2011, 47(4): 829-837 [6] [6] Garg D. Advances in Global Pseudospectral Methods for Optimal Control [Ph.D. dissertation], Massachusetts Institute of Technology, USA, 2011 [7] [7] Elnagar G, Kazemi M A, Rzaazghi M. The pseudospectral Legendre method for discretizing optimal control problems. IEEE Transactions on Automatic Control, 1995, 40(10): 1793-1796 [8] [8] Elnagar G N, Rzaazghi M. A collocation-type method for linear quadratic optimal control problems. Optimal Control Applications and Methods, 1997, 18(3): 227-235 [9] [9] Fahroo F, Ross I M. Costate estimation by a Legendre pseudospectral method. Journal of Guidance, Control, and Dynamics, 2001, 24(2): 270-277 [10] Elnagar G N, Kazemi M A. Pseudospectral Chebyshev optimal control of constrained nonlinear dynamical systems. Computational Optimization and Applications, 1998, 11(2): 195-217 [11] Fahroo F, Ross I M. Direct trajectory optimization by a Chebyshev pseudospectral method. Journal of Guidance, Control, and Dynamics, 2002, 25(1): 160-166 [12] Gong Q, Ross I M, Fahroo F. Costate computation by a Chebyshev pseudospectral method. Journal of Guidance, Control, and Dynamics, 2010, 33(2): 623-628 [13] Fornberg B. A Practical Guide to Pseudospectral Methods. New York: Cambridge University Press, 1998 [14] Weideman J A C, Trefethen L N. The kink phenomenon in Fejr and Clenshaw-Curtis quadrature. Numerische Mathematik, 2007, 107(4): 707-727 [15] Berrut J P, Trefethen L N. Barycentric Lagrange interpolation. SIAM Review, 2004, 46(3): 501-517 [16] Gill P E, Murray W, Saunders M A. SNOPT: an SQP algorithm for large-scale constrained optimization. SIAM Review, 2005, 47(1): 99-131 [17] Biegler L T, Zavala V M. Large-scale nonlinear programming using IPOPT: an integrating framework for enterprise-wide dynamic optimization. Computers and Chemical Engineering, 2009, 33(3): 575-582 [18] Waldvogel J. Fast construction of the Fejr and Clenshaw-Curtis quadrature rules. BIT Numerical Mathematics, 2006, 46(1): 195-202 [19] Costa B, Don W S. On the computation of high order pseudospectral derivatives. Applied Numerical Mathematics, 2000, 33(1-4): 151-159 [20] Gong Q, Kang W, Ross I M. A pseudospectral method for the optimal control of constrained feedback linearizable systems. IEEE Transactions on Automatic Control, 2006, 51(7): 1115-1129 [21] Betts J T. Practical Methods for Optimal Control and Estimation Using Nonlinear Programming (Second edition). Philadelphia: Society for Industrial and Applied Mathematics, 2010 期刊类型引用(14)
1. 王子建. 基于深度聚合网络的单幅图像超分辨率重建. 广播电视信息. 2023(05): 53-58 . 
百度学术2. 董国强，卜丽静，赵瑞山，张正鹏. 多先验约束正则化模型图像复原方法. 测绘工程. 2023(03): 19-26 . 百度学术3. 白勇强，禹晶，李一秾，肖创柏. 基于深度先验的盲图像去模糊算法. 电子学报. 2023(04): 1050-1067 . 百度学术4. 李公平，陆耀，王子建，吴紫薇，汪顺舟. 基于模糊核估计的图像盲超分辨率神经网络. 自动化学报. 2023(10): 2109-2121 . 本站查看5. 彭天奇，禹晶，肖创柏. 基于跨尺度低秩约束的图像盲解卷积算法. 自动化学报. 2022(10): 2508-2525 . 本站查看6. 彭天奇，禹晶，郭乐宁，肖创柏. 基于跨尺度字典学习的图像盲解卷积算法. 光学精密工程. 2021(02): 338-348 . 百度学术7. 兰诚栋，林宇鹏，方大锐，陈建. 多视点稀疏测量的图像绘制方法. 自动化学报. 2021(04): 882-890 . 本站查看8. 鱼轮，李晖晖. 基于学习字典和稀疏约束的湍流退化图像盲复原. 电子测量技术. 2020(17): 105-109 . 百度学术9. 兰叶深，周明安，刘文军，毛建辉. 字典学习在轴承表面缺陷分类中的研究. 内燃机与配件. 2020(22): 147-148 . 百度学术10. 龚军辉，刘小燕，周建松，孙刚. 基于时空双稀疏表示的成人ADHD脑网络检测与分析. 自动化学报. 2019(10): 1903-1914 . 本站查看11. 陈红. 稀疏正则化模糊图像盲复原模型仿真. 计算机仿真. 2019(12): 323-326 . 百度学术12. 李俊山，杨亚威，张姣，李建军. 退化图像复原方法研究进展. 液晶与显示. 2018(08): 676-689 . 百度学术13. 陈春雷，叶东毅，陈昭炯. 多局部模糊核融合的图像盲去模糊算法. 光子学报. 2018(10): 205-215 . 百度学术14. 李喆，李建增，胡永江，张岩. 基于频谱预处理与改进霍夫变换的离焦模糊盲复原算法. 图学学报. 2018(05): 909-916 . 百度学术其他类型引用(13)
-
下载:
计量
- 文章访问数: 2931
- HTML全文浏览量: 196
- PDF下载量: 1352
- 被引次数: 27
下载:
