An Identification Method for Nonlinear Systems with Colored Measurement Noise
-
摘要: 利用最大似然判据, 本文提出了一种带有色量测噪声的非线性系统辨识方法. 首先, 利用量测差分方法将有色量测噪声白色化, 获得新的量测方程, 从而将带有色量测噪声的非线性系统辨识问题转化成带白色量测噪声和一步延迟状态的非线性系统辨识问题. 其次, 利用期望最大化(Expectation maximization, EM)算法提出了一种新的基于最大似然估计的非线性系统辨识方法, 该算法由期望步骤(Expectation step, E-step)和最大化步骤(Maximization step, M-step)两部分组成. 在期望步骤中, 基于当前估计的参数并利用带有色量测噪声的高斯近似滤波器和平滑器, 近似计算完整的对数似然函数的期望. 在最大化步骤中, 近似计算的似然函数期望值被最大化, 并且通过解析更新获得噪声参数估计, 通过Newton更新方法获得模型参数的估计. 最后, 数值仿真验证了本文提出算法的有效性.Abstract: In this paper, an identification method for nonlinear systems with colored measurement noise is proposed by using the maximum likelihood criterion. Firstly, the colored measurement noise is decorrelated based on the measurement differencing approach, and a new measurement equation is derived. Thus, the nonlinear system identification problem with colored measurement noise is transformed into the nonlinear system identification problem with white measurement noise and one-step delayed state. Secondly, a new nonlinear system identification method with maximum likelihood estimation is proposed based on the expectation maximization (EM) algorithm, which consists of expectation step (E-step) and the maximization step (M-step). In the E-step, the expectation of the complete data log-likelihood function is approximately calculated based on currently estimated parameters and the Gaussian approximated filter and smoother for nonlinear system with colored measurement noise. In the M-step, the approximately calculated expectation value is maximized, and noise parameter estimations are updated analytically and model parameter estimations are updated approximately by using Newton method. Finally, the efficiency of the proposed algorithm is illustrated in numerical simulations.
-
磨浆过程主要为后续造纸过程提供满足相应物理特性的纸浆纤维, 同时也是保证纸品质量的前提.但由于磨机物理结构复杂性以及磨浆运行过程中受外界不确定因素的影响, 导致纤维长度具有较强的随机分布特性, 而传统的均值或者方差并不足以表征整个纤维长度随机分布(Fiber length stochastic distribution, FLSD)特征, 即具有非高斯分布特性[1].而事实上, FLSD形状作为衡量磨浆过程中最主要的生产指标之一, 不仅影响到纸浆脱水效率和后续造纸过程的电耗, 而且直接影响到整个制浆和造纸工业的能耗和最终的纸品质量[2-6].在这种情况下, 就必须寻找一种能够对磨浆过程输出FLSD进行有效控制的方法.
虽然一些学者较早地意识到FLSD在纸浆生产过程中的重要作用, 但目前仍利用离线获得的纤维长度的均值和方差作为衡量纸浆纤维质量的生产指标[3-9].如文献[3$-$6]采用长纤维百分含量作为评价纸浆质量的工艺指标, 而事实上, 这种以统计纤维束长度的均值和方差并足以有效表征整个纤维长度的分布特征, 甚至会导致FLSD信息的缺失, 这主要因为木片等富含纤维的纸浆原料经过盘磨机的横向挤压和纵向帚化后, 虽然使得纤维束逐渐分解成单根纤维, 然而, 也导致纤维长度的分布形状具有很强的随机性和不确定性.此外, 加之纤维长度分布在线检测仪器缺失, 使得通过离线纤维长度的均值和方差来调整过程操作变量, 其检测的滞后性难以使纸浆纤维质量稳定在工艺要求范围内, 致使长纤维百分含量的控制大多过度依赖操作人员的经验, 然而人工调整主观性较强, 常常造成工况波动大, 严重影响磨浆过程生产指标的稳定性.
另一方面, 现有针对随机动态系统的最小方差控制、自校正控制、随机线性二次型控制, 均假设系统服从高斯分布, 主要集中针对过程输出随机变量的均值和方差进行建模和控制, 然而, 针对输出变量服从非高斯分布的有界动态随机系统, 1996年王宏教授提出了直接设计控制器以使输出PDF形状跟踪期望PDF形状的控制策略, 并系统地提出了多种建模和控制方法[10-18].这类控制策略包含了传统以输出均值和方差为目的随机控制方法, 具有更为广泛的应用.近些年, 随机分布控制理论已成功应用于各类具有随机分布动态特性的工业过程, 如造纸过程的絮凝粒径分布[10, 12]、燃烧过程的火焰分布[18-19], 聚合过程的分子量分布[20], 铜粗选过程的泡沫尺寸分布[21]等, 这些过程输出随机变量均不能满足高斯分布特征, 并且具有较强的随机分布动态特性.
另外, 在制浆和造纸工业领域, 目前已有多种用于测量纸浆纤维各种形态参数的在线自动化检测和分析仪器, 这些检测仪器常采用数码CCD摄像机获取的纤维图像经数字化后传输到计算机系统进行处理, 通过实时二维图像分析软件将每根纤维从图像信号中识别出来, 逐一测量纤维的形态参数, 能迅速准确地获得纤维形态参数是测量结果.如加拿大Optest公司FQA-360、芬兰Kajaani公司FS-300、丹麦Fiber-Visions等[8-9], 这些先进测量仪器为研究基于磨浆过程输出FLSD的建模及控制提供了技术支撑.
磨浆过程输出FLSD具有典型的非高斯分布动态特性, 采用传统纤维长度的均值和方差难以有效描述其分布特征, 本文根据随机分布相关控制理论[10], 利用RBF神经网络逼近输出FLSD的PDF, 为了改善传统线性权值模型[14-19]精度不高、泛化能力不强等问题, 采用随机权神经网络(Random vector functional-link networks, RVFLNs)[22-24]建立表征输入变量和权值向量之间的非线性模型, 基于磨浆过程输出FLSD模型提出了一种预测PDF控制方法, 实现了对输出FLSD形状的跟踪控制, 基于工业数据实验表明了所提方法的有效性.
1. 磨浆过程描述
典型的磨浆过程工艺流程如图 1所示, 磨浆过程即是将植物原料经盘磨机反复研磨后, 经汽浆分离后获得造纸所需的纸浆纤维.其主要包括喂料系统, 供水系统, 液压伺服系统和磨盘调速系统.当磨机运行时, 首先将经过将被筛选木片在蒸煮仓里进行高温预处理, 经清除杂质后的木片在螺旋喂料器作用下送入磨室.磨机作为磨浆过程中的核心设备, 主要有定盘、动盘、电液伺服装置和主电机等组成.当预处理后木片和稀释水注入磨区后, 利用电液伺服装置可以实时调节磨盘间隙, 动盘在主电机带动下通过机械摩擦、剪切、撕裂、切割等作用, 最终将预处理后木片分解为单根纤维.然后, 纸浆通过送入旋风分离器实现汽浆分离, 最终获得满足造纸过程所需的纸浆纤维.可以看出, 由于工艺流程长及现场环境恶劣等原因, 若操作变量调节不当, 即便通过磨机的反复研磨, 也难以获得满意的FLSD形状, 这样不但导致纸浆质量无法满足工艺要求, 而且也极易造成原料浪费和过程能耗过大.
磨浆过程的主要操作变量包括磨盘间隙、稀释水流量、动盘转速和螺旋喂料量等.研究表明:长纤维百分含量作为衡量磨浆过程纤维长度均值的工艺指标, 主要与注入磨区的稀释水流量和磨盘间隙密切相关[3, 5-6], 由于磨浆过程主要目的是对纤维进行切断、压溃、吸水膨胀, 最终使纤维束分离为单根纤维, 所以, 磨盘间隙直接影响纸浆纤维被切断和压溃强度.同时, 为使植物纤维能够较好地分离为单根纤维, 需要让纤维获得足够的水份进行膨胀, 纸浆纤维的吸水膨胀程度主要取决于注入磨区的稀释水流量.而植物纤维被切断、压溃强度以及吸水膨胀程度决定了最终获得的纤维长度随机分布形状.另外, 结合实际工程经验, 通常情况下分别通过改变喂料螺旋转速和供水泵转速来调节螺旋喂料量和稀释水流量, 根据工程实践经验, 当产量一定的情况下, 动盘转速是固定不变的, 螺旋喂料量也是恒定的.因此, 稀释水流量和磨盘间隙不但可以看作影响长纤维百分含量的主要变量, 同样也可以作为影响最终FLSD形状的关键变量, 对整个制浆生产流程都起着极为关键的作用.因此, 本文将稀释水流量和磨盘间隙作为影响磨浆过程输出FLSD形状的关键输入变量.
2. 控制策略
根据有界动态随机分布控制相关理论[10]可知, 随机分布系统模型主要由随机变量的PDF输出部分和权值与输入变量之间的动态部分组成.为了表示输入变量和输出PDF之间的动态关系, 常引入一组基函数(如B样条基函数[16-17, 20-21]、RBF基函数[15, 18-19])来逼近输出随机变量的PDF, 通过调节基函数的权值来控制输出PDF形状, 这样将随机分布系统的输出PDF和输入变量之间的动态关系转化为权值向量和输入变量之间的动态关系, 最终通过对相对应权值的控制实现对输出PDF形状的动态调节.为此, 本文针对磨浆过程输出FLSD提出预测PDF控制策略如图 2所示, 具体如下:
图 2 输出纤维长度随机分布预测PDF控制策略图Fig. 2 Strategy diagram of the predictive PDF control for the output FLSD1) 首先, 需寻找一组合适的RBF基函数来逼近输出FLSD的PDF.采用RBF神经网络逼近输出PDF的均方根, 基于迭代学习方法实现RBF基函数参数整定, 并对实际输出PDF相应权值向量进行估计.
2) 其次, 针对常规线性权值向量模型精度不高, 泛化能力不强等缺点, 本文基于随机权神经网络[22-24]方法建立输入变量和前$n-1$个权值向量之间的非线性模型, 进而获得磨浆过程输出FLSD模型.
3) 最后, 基于输出FLSD模型设计预测PDF控制器, 使得输出PDF获得良好的目标跟踪能力.
3. 纤维长度随机分布预测PDF控制
3.1 随机分布系统模型
随着数据采集技术和检测仪器的快速发展, 对随机变量的输出PDF等已经有了较为成熟的检测方法.为了方便描述各种随机过程, 假设为描述动态随机系统输出的一致有界随机过程变量, 为$k$时刻控制随机系统分布形状的输入向量, 这表明在任一采样时刻$k$, $z(k)$就可以通过其概率密度函数来描述, 其定义式如下:
$ P({a<z(k)<\zeta}, {u(k)})=\int_{a}^\zeta{\gamma(y, u(k))}\textrm{d}y $
(1) 式中, $P({a<z(k)<\zeta}, {u(k)})$表示随机系统在$u(k)$作用下输出落在区间$[a, \zeta]$内的概率, 即$z(k)$的输出PDF${\gamma(y, u(k))}$形状由输入变量$u(k)$控制.假设区间$[a, b]$为已知, 并且输出PDF${\gamma(y, u(k))}$连续且有界, 采用如下具有高斯型RBF神经网络来逼近输出的PDF均方根
$ R_l(y)=\textrm{exp}\left(-\frac{(y-\mu_l)^2}{\sigma_l^2} \right), l=1, 2, \cdots, n $
(2) 式中, $l$表示第$l$个网络节点, $n$为网络节点总数, $\mu_l$和$\sigma_l$分别表示第$l$个网络节点函数的中心值和宽度.根据RBF神经网络逼近原理, 此时输出PDF$\gamma(y, u(k))$的均方根可以表示为
$ \sqrt{\gamma(y, u(k))}=\pmb C(y)\pmb V(k)+R_n(y){\omega_n(k)}+e_0(y, k) $
(3) 式中, $\pmb C(y)=[R_1(y), R_2(y), \cdots, R_{n-1}(y)]$, , $\omega_{n}(k)$为第$n$个基函数相对应的权值, $e_0(y, k)$为逼近误差.此外, 由于输出PDF需要满足隐含条件
$ \int_{a}^b{\gamma({y, u(k)})}\textrm{d}y=1 $
(4) 为便于分析忽略逼近误差$e_0(y, k)$, 则第$n$个权值$\omega_n(k)$可用权值向量的非线性函数$h(\pmb V(k))$表示为
$ h(\pmb V(k))=\frac{{\sqrt{\Sigma_2-\pmb V^\textrm{T}(k)\Sigma_3\pmb V(k)}}-\Sigma_1\pmb V(k)}{\Sigma_2} $
(5) 式中, , \textrm{d}\textit{y}, $\Sigma_2=\int_{a}^{b}{{R_n^2(y)}}\textrm{d}y$, $\Sigma_3=\Sigma_2\Sigma_0-\Sigma_1^{\rm{T}}\Sigma_1.$
从式(5)可以看出, 若保证非线性函数$h(\pmb V(k))$存在, 需满足如下约束条件:
$ {\pmb V}^\textrm{T}(k)\Sigma_4 \pmb {V}(k)<1 $
(6) 式中, $\Sigma_4=\Sigma_2^{-1} \Sigma_3>0$, 不等式(6)可以看作在动态权值向量$\pmb V(k)$需满足的约束条件.由式(5)可以看出, 一旦所有的基函数确定, 由于有界区间$[a, b]$已知, $\pmb C(y)$和$R_n(y)$均已知, 这就意味着在$n$个权值中有$n-1$个是相互独立的.
当实际工业过程中输出随机变量的PDF可测量时, 可通过如式(7)$ \sim $(9)方法对相应的权值进行估计.结合式(3)和式(5)所示的输出PDF的均方根可以表示为
$ \sqrt{\gamma(y, u(k))}=\left[ \begin{array}{cc} \pmb C(y)&R_n(y)\\ \end{array} \right]\left[ \begin{array}{c} \pmb V(k)\\ h(\pmb V(k))\\ \end{array} \right] $
(7) 对式(7)两边左乘$[\pmb C(y)R_n(y)]^{\rm T}$, 并对两端在区间$[a, b]$上进行积分, 由此可得到:
$ \left[ \begin{array}{cc} A_1\\ A_2\\ \end{array} \right] =\left[ \begin{array}{cc} \Sigma_0 & \Sigma_1 \\ \Sigma_1^{\rm{T}} & \Sigma_2 \\ \end{array} \right] \left[ \begin{array}{cc} \pmb V(k)\\ h(\pmb V(k))\\ \end{array} \right] $
(8) 式中, , $A_2=\int_{a}^{b}{R_n}(y)$ $\sqrt{\gamma(y, u(k))}\textrm{d}y$.
当矩阵非奇异时, 式(8)所示权值可以通过矩阵求逆获得:
$ \left[ \begin{array}{cc} \pmb V(k)\\ h(\pmb V(k))\\ \end{array} \right] =\left[ \begin{array}{cc} \Sigma_0 & \Sigma_1 \\ \Sigma_1^{\rm{T}} & \Sigma_2 \\ \end{array} \right]^{-1}\left[ \begin{array}{c} A_1\\ A_2\\ \end{array} \right] $
(9) 式(9)揭示了输出随机变量的PDF与权值向量之间的关系, 可以看出当RBF基函数确定后, 只要输出PDF可测量, 便可很容易获得相应的权值向量.通常情况下在获得相应的权值向量之后, 采用最小二乘法或子空间辨识法获得权值和输入变量之间的线性动态模型[15, 18-20]. 从式(9)可以看出, 在RBF基函数已知的情况下, 若要获得理想的输出FLSD模型.首先, 需要对不同时刻实际输出PDF进行权值估计.其次, 建立输入变量和权值向量之间动态模型.然而, 由于实际工业过程高度非线性、机理复杂等原因存在, 采用常规线性模型难以有效描述权值向量的动态特性.为此, 必须采用有效手段, 获取基于磨浆过程输出FLSD模型, 以实现FLSD形状的在线连续估计及控制.
另外, 对于不能或者难以获得机理模型的复杂工业过程, 基于数据驱动建模方法通常被看作一种非常有效的替代手段.目前, 常见的数据驱动建模方法主要有支持向量机方法[25]、模糊推理方法[26]、案例推理方法[27]以及随机权神经网络方法[22]等.其中, 随机权神经网络作为一种简单易用、有效的单隐层前馈神经网络学习算法, 在保证逼近任意连续函数的前提下, 采用随机给定神经元隐含层权值和偏置, 通过计算隐含层输出矩阵的广义逆建立学习网络, 克服了传统单隐层神经网络的缺点, 由于训练速度快, 模型结构简单、易于实现以及泛化能力强等鲜明特点, 在很多领域获得广泛应用[23-24].本文利用随机权神经网络方法建立输出变量和前$n-1$个权值向量之间非线性预测模型.因此, 最终磨浆过程输出FLSD模型可以表示为
$ \left \{\begin{aligned}\pmb V(k+1)=f \left( \pmb V(k), \pmb u(k) \right)\\ \sqrt{\gamma(y, \pmb u(k))}=\pmb C(y)\pmb V(k)+R_n(y)h(\pmb V(k)) \end{aligned}\right. $
(10) 式中, $f(\cdot)$为表示过程输入和权值之间的非线性表达式, 可以看出在控制输入和输出PDF可测量情况下, 同时在RBF基函数已知时, 在利用式(9)获得权值向量之后, 通过随机权神经网络方法很容易获得输入变量和权值向量之间的非线性模型.
从式(10)明显发现要获得理想的输出FLSD模型, 首先需要选择一组合适的RBF基函数, 若RBF基函数的中心值和宽度选择不当, 不但影响PDF输出部分的近似精度, 而且也难以获得满意的输出FLSD模型.因此, RBF基函数的选择对于磨浆过程输出FLSD的模建模精度显得至关重要.
3.1.1 RBF基函数参数整定
为了提高随机分布模型输出PDF对实际输出PDF逼近精度, 文献[15]中提出了一种基于迭代学习控制原理的RBF基函数参数更新算法, 假设给出初始的RBF基函数参数, 结合实际输出PDF数据, 通过式(9)可以获得初始的权值向量, 并将此初始向量乘以相应的初始RBF基函数, 便可获得近似的输出PDF.基于迭代学习控制原理, RBF基函数的中心值和宽度要根据上一迭代学习周期的逼近误差来调节, 通过对RBF基函数参数的反复调节, 最终使得将近似输出PDF与实际输出PDF之间误差最小.因此, 所采用的性能指标主要考察近似输出PDF跟踪期望的输出PDF能力, 所以所用误差来自迭代周期内的每一个采样点, 采用如下性能指标
$ J_{p}(i)=\int_a^b \left({\sqrt{\gamma_{p+1}(y, u_i)}-{\sqrt{g(y)}}}\right)^2\textrm{d}y $
(11) 式中, $J_{p}(i)$可以表示为第$p$个迭代周期后第$i$个采样时刻RBF神经网络近似输出PDF的均方根与期望输出PDF$g(y)$的均方根之间的逼近性能指标.同时, 定义在第$p$个迭代周期后$M$个采样点近似输出PDF分别与期望输出PDF之间构成性能指标向量可表示如下
$ \pmb E_{p}=[J_{p}(1), J_{p}(2), \cdots, J_{p}(M)]^{\rm T} $
(12) 为了实现式(10)所示的RBF基函数参数整定, 在相邻的第$p$次和第$(p+1)$次迭代学习周期内采用如下P型迭代学习控制率
$ \begin{cases} \mu_{l, p+1}=\mu_{l, p}+{\pmb\alpha}_{\mu}{\pmb E_{p}} \\ \sigma_{l, p+1}=\sigma_{l, p}+{\pmb \beta}_{\sigma}{\pmb E_{p}}\\ \end{cases} $
(13) 式中, 学习参数$\alpha_\mu$和$\beta_\sigma$分别定义如下
$ \left\{ \begin{align} & {{\pmb\alpha}_{\mu}}={{\zeta }_{\mu }}[{{\mathit{\lambda }}_{1}},{{\mathit{\lambda }}_{2}},\cdots ,{{\mathit{\lambda }}_{\mathit{M}}}] \\ & {{\pmb \beta}_{\sigma}}={{\zeta }_{\sigma }}[\mathit{\lambda }{{\mathit{'}}_{1}},\mathit{\lambda }{{\mathit{'}}_{2}},\cdots ,\mathit{\lambda }{{\mathit{'}}_{M}}] \\ \end{align} \right. $
其中, $\lambda$和$\lambda'$分别为学习元素, 参数$\zeta_\mu$和$\zeta_{\sigma}$分别为迭代学习率.从式(13)看出所有元素均为非负, 这表明迭代学习率可以为正数也可以为负数, 这就意味着RBF基函数的中心值和宽度随着迭代学习次数呈现出增加或者降低趋势.
3.1.2 基于RVFLNs权值模型
从式(10)可以看出当前时刻输出PDF不但与磨浆过程输入变量有关, 同时与前一时刻的输出PDF形状也密切相关, 而输出PDF通过调节RBF基函数中心值和宽度以及相应的权值, 所以在获得理想RBF基函数参数之后, 通过式(9)可以获得所有时刻输出PDF相对应权值向量, 然而由于权值向量之间相互耦合, 因此, 权值向量模型可以看作是一个多输入多输出的回归建模问题.为此, 采用基于RVFLNs[22-24]建立输入变量与前$n-1$个权值向量之间的非线性模型.
假设分别有$m$个输入变量与权值向量组成的的样本集合, 其中为磨浆过程的输入变量, 表示$n$个权值中的前$n-1$个权值向量.对于一个有$L$个隐层节点, 若以$g(x)$作为激活函数的RVFLNs输出可以表示为
$ f_{\textrm{R}}({\boldsymbol{u}}_i)=\sum\limits_{j=1}^L {\pmb \beta_j g(\pmb \omega_j \cdot \pmb u_i+b_j)} $
(14) 式中, $\pmb u_i$为磨浆过程的输入变量, 为$m$个输入节点连接第$j$个隐含层的输入权重, 为第$j$个隐含层连接输出节点的输出权重, $b_j$是第$j$个隐层单元的偏置, $\pmb \omega_j \cdot \pmb u_j$表示和$\pmb u_i$的内积.
随机权神经网络和其他单隐层神经网络学习目标一样均是使得模型输出与实际输出$\pmb V_i$之间误差最小, 即有.当存在$\pmb \beta_j$, $\pmb \omega_j$和$\pmb b_j$, 使得
$ \pmb V_i=\sum\limits_{j=1}^L {\pmb \beta_j g(\pmb \omega_j \cdot \pmb u_j+b_j)} $
(15) 此时, 将式(15)可以矩阵表示为
$ {H} {\beta}={Y} $
(16) 式中, $H$为隐含层输出矩阵, $\beta$为输出权重, $Y$为预测模型的权值输出.且有
$ \begin{array}{l} H\left( {{\mathit{\boldsymbol{\omega }}_1}, \cdots ,{\mathit{\boldsymbol{\omega }}_L},{b_1}, \cdots ,{b_L},{\mathit{\boldsymbol{u}}_1}, \cdots ,{\mathit{\boldsymbol{u}}_L}} \right) = \\ \;\;\;\;{\left[ \begin{array}{l} g\left( {{\mathit{\boldsymbol{\omega }}_1} \cdot {\mathit{\boldsymbol{u}}_1} + {b_1}} \right)\;\;\; \cdots \;\;\;g\left( {{\mathit{\boldsymbol{\omega }}_L} \cdot {\mathit{\boldsymbol{u}}_1} + {b_L}} \right)\\ \;\;\;\;\;\;\;\;\; \vdots \;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \cdots \\ g\left( {{\mathit{\boldsymbol{\omega }}_1} \cdot {\mathit{\boldsymbol{u}}_N} + {b_1}} \right)\;\; \cdots \;\;\;g\left( {{\mathit{\boldsymbol{\omega }}_L} \cdot {\mathit{\boldsymbol{u}}_N} + {b_L}} \right) \end{array} \right]_{N \times L}} \end{array}\\ \beta = {\left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{\beta }}_1^{\rm{T}}}\\ \vdots \\ {\mathit{\boldsymbol{\beta }}_L^{\rm{T}}} \end{array}} \right]_{L \times (n - 1)}}\quad ,\quad Y = {\left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{V}}_1^{\rm{T}}}\\ \vdots \\ {\mathit{\boldsymbol{V}}_L^{\rm{T}}} \end{array}} \right]_{N \times (n - 1)}} $
从式(16)可以看出, 当输入权重和隐层偏置$b_j$被随机确定后, 只需调整输出层权值就可以使网络具有较好的逼近性能.为了能够较好地训练上述网络, 希望获得最优的输出权重$\hat{\beta}$, 使得
$ \hat{\beta}=\arg\min\limits_{\beta}\|H{\beta}-Y\|^2 $
(17) 此时随机权神经网络的学习问题就转化为式(16)所示的线性系统${H\beta}=Y$的最小二乘求解问题, 为此隐层输出矩阵$H$就能被唯一确定, 此时可以获得输出权重$\beta$可以表示为
$ \hat{\beta }={{H}^{\dagger }}Y $
(18) 式中, $ {{H}^{\dagger }}$是矩阵$H$的\textrm{Moore-Penrose}广义逆.可以看出此算法只需要设置网络的隐层节点个数, 便可以随机初始化输入权重和偏置并得到相应的输出权重.所以该算法在执行过程中不需要调整网络的输入权值以及隐元的偏置, 便可以获得唯一的最优解.
为了更好地反映磨浆过程输出FLSD动态特性, 将当前采样时刻输入变量, $m$为输出变量个数以及当前时刻相对应的权值向量作为非线性模型综合输入, 即建立的权值动态预测模型用于实现如下的非线性动态映射关系:
$ \pmb V_{\textrm{m}}(k+1)=f_{\textrm{R}}\{\pmb V(k), \pmb u(k)\} $
(19) 式中, $\pmb V_{\textrm{m}}(k+1)$为模型输出的前$n-1$个权值向量.此时, 在第$k$时刻预测输出PDF为
$ {\sqrt{\gamma_{\textrm{m}}(y, k+1)}}= \pmb C(y)\pmb V_{\textrm{m}}(k+1)+R_n(y) {}\nonumber \\ {}h(\pmb V_{\textrm{m}}(k+1))= \pmb C(y) f_{\textrm{R}}\{\pmb V(k), \pmb u(k)\}+{}\nonumber \\ {}R_n(y)h(f_{\textrm{R}}\{\pmb V(k), \pmb u(k)\}) $
(20) 3.2 预测PDF控制
由于实际工业过程中非线性、时变、模型失配和随机扰动等不确定性因素的存在, 模型输出PDF难以与实际输出PDF完全一致, 然而在滚动优化过程中, 需要实际输出PDF与模型输出PDF保持一致, 因此, 通常采用反馈校正来降低过程的不确定性对系统性能的影响, 提高系统的控制精度和鲁棒性.假设模型在$k$时刻第$j$步预测输出PDF为
$ {\sqrt{\gamma_{\textrm{m}}(y, k+j)}}=\pmb C(y)\pmb V_{\textrm{m}}(k+j)+ \nonumber \\ \qquad{}R_n(y)h(\pmb V_{\textrm{m}} {}(k+j))={}\nonumber \\ \qquad{}R_n(y\pmb C(y) f_{\textrm{R}}\{\pmb V(k+j-1), \pmb u(k+j-1)\}+{}\nonumber \\ \qquad{}R_n(y)h(f_{\textrm{R}}\{\pmb V(k+j-1), \pmb u(k+j-1)\}) $
(21) 另外, 在第$k$时刻实际输出PDF和预测输出PDF之间的误差为
$ e(y, k)=\sqrt{\gamma(y, k)}-\sqrt{\gamma_{\textrm{m}}(y, k)}={}\pmb C(y)(\pmb V(k)-\nonumber \\ \qquad{}f_{\textrm{R}}\{\pmb V(k-1)\pmb u(k-1)\})+R_n(y)\cdot\nonumber \\ \qquad{}{(h(\pmb V(k)-h(f_{\textrm{R}}\{\pmb V(k-1), \pmb u(k-1)\}))} $
(22) 利用该误差对第$j$步预测输出PDF ${\sqrt{\gamma_{\textrm{m}}(y, k+j)}}$进行反馈修正, 补偿后预测输出PDF为
$ \sqrt{\gamma_{\textrm{p}}(y, k+j)}=\sqrt{\gamma_{\textrm{m}}(y, k+j)}+\beta_j e(y, k) $
(23) 式中, $\beta_j(0<\beta_j<1)$为校正系数.此时结合式(21)$ \sim $(23)在$k$时刻第$j$步期望输出PDF和补偿后的预测输出PDF之间误差为
$ e_{\textrm{p}}(y, k+j)=\sqrt{\gamma_{\textrm{g}}(y, k+j)}- \sqrt{\gamma_{\textrm{p}}(y, k+j)}= {}\nonumber \\ {} \pmb C(y)[\pmb V_{\textrm{g}}(k+j)-\nonumber \\ {}f_{\textrm{R}}\{\pmb V(k+j-1), \pmb u(k+j-1)\}] +{}\nonumber \\ {} R_n(y)[h(\pmb V_{\textrm{g}}(k+j))-\nonumber \\ h(f_{\textrm{R}} \{\pmb V(k+j-1), {}{} \pmb u(k+j-1)\})]+\nonumber \\ {}\beta_j[\pmb C(y)(\pmb V(k)-f_{\textrm{R}} \{(\pmb V(k-1), {}\nonumber \\ {} \pmb u(k-1)\})]+\beta_j R_n(y)[h(\pmb V(k))-\nonumber \\ {}h(f_{\textrm{R}} \{\pmb V(k-1), {} \pmb u(k-1)\})] $
(24) 式中, ${\gamma_{\textrm{g}}(y, k+j)}$和${\gamma_{\textrm{p}}(y, k+j)}$分别为$k$时刻第$j$步的期望输出PDF和预测输出PDF.
此外, 预测控制作为一种优化控制算法, 常通过最小化系统的性能指标函数来确定未来的最优控制序列, 使得未来预测输出尽可能接近期望的目标输出.在实际工业过程中, 为了保证操作的可行性等要求, 普遍存在着输入变量带约束的情形, 同时对控制作用的大小加以约束, 避免控制作用变化过于剧烈.本文设计预测PDF控制器目的是尽可能地使输出PDF尽可能跟踪期望输出PDF, 所以选取如下所示性能指标函数
$ \min\limits_{\triangle u} J=\nonumber\\ \quad\sum\limits_{j=1}^{N_\textrm{p}}\int_a^b\left({\sqrt{\gamma_{\textrm{g}}(y, k+j)}}-{\sqrt{\gamma_{\textrm{p}}(y, k+j)}}\right)^2\textrm{d}y+\nonumber \\ \quad\sum\limits_{j=0}^{N_\textrm{c}-1}[\lambda_j \Delta u(k+j)]^2 {}\nonumber \\ {\rm s.t.}u_{\textrm{min}}<\Delta u(k+j)+u(k+j-1)<u_{\textrm{max}}, {}\nonumber \\ \quad\Delta u_{\textrm{min}}<\Delta u(k+j)<\Delta u_{\textrm{max}} $
(25) 式中, $N_\textrm{p}$和$N_\textrm{u}$分别为预测时域和控制时域, $\lambda_j$为控制增量加权系数, $u_{\textrm{max}}$和$u_{\textrm{min}}$分别为输入变量的上限和下限值, $\Delta u_{\textrm{max}}$和分别为输入变量变化率的上限和下限值.可以看出对输出PDF跟踪控制最终转化为对前$n-1$个权值跟踪控制.
从式(25)明显可以看出上述预测PDF控制器的设计可以看作是一个带有约束条件的非线性优化求解问题.针对上述求解问题通常采用遗传算法、粒子群算法、序列二次规划算法(Sequence quadratic program, SQP)等优化算法获得非线性最优预测控制率, 其中, SQP算法作为一种求解约束非线性优化问题的有效方法之一, 具有收敛性快、计算效率高、边界搜索能力强, 在实际中受到广泛重视和应用.本文采用SQP方法求解式(25)所示的带约束的非线性规划问题设计预测PDF控制器, 使得磨浆过程输出PDF具有良好的目标跟踪能力.
4. 工业数据验证
本文利用某化机浆磨浆过程的稀释水流量、磨盘间隙以及FLSD PDF生产数据对所提方法进行数据验证, 具体如下:
4.1 RBF基函数参数整定
为了获得磨浆过程输出FLSD模型, 首先, 需要确定一组合适的RBF基函数近似输出PDF, 采用基于迭代学习方法研究RBF神经网络对期望输出PDF的逼近效果, 并将得到RBF基函数作为本批次近似输出PDF的基函数, 此外, 通过大量FLSD的PDF数据分析获得期望输出PDF, 本文选择4个RBF基函数来验证对期望输出PDF的近似效果, 这里中心值和宽度参数迭代学习率分别为$\alpha_\mu=0.01$, $\beta_\sigma=0.005$, 另外, 假设中心值和宽度的初始值如下所示
$ \mu_1=0.3, \mu_2=0.8, \mu_3=1.3, \mu_4=1.8, {}\nonumber \\ {}\sigma_1^2=\sigma_2^2=\sigma_3^2=\sigma_4^2=0.06 $
(26) 基于式(26)所示的初始RBF基函数, 首先, 可以通过式(9)获得期望输出PDF权值估计, 然后利用得到的估计权值与初始RBF基函数相乘便得到对应的逼近值, 以此利用迭代学习方法通过调整中心值和宽度, 直到获得理想的逼近效果.经过100次迭代学习后, 获得中心值和宽度分别为
$ \mu_1=0.40, \mu_2=0.91, \mu_3=1.25, \mu_4=1.46, {}\nonumber \\ {}\sigma_1^2=0.068, \sigma_2^2=0.074, \sigma_3^2=0.098, \sigma_4^2=0.027 $
(27) 另外, 图 3为在迭代学习50次和100次后的RBF基函数位置变化趋势图, 可以看出随着迭代次数的增加, 中心值和宽度逐渐向理想位置移动. 图 4为性能指标函数值随迭代学习次数的变化趋势, 可以看出随着迭代次数的增加, 目标性能函数逐渐减小, 并在迭代学习80次时基本不再变化.中心值和宽度随迭代次数变化趋势分别如图 5和图 6所示, 从图 5和图 6看出在迭代学习100次后, 中心值和宽度均趋于平稳. 图 7为在迭代学习100次后, 近似输出PDF与期望输出PDF的逼近结果, 可以看出本文方法对输出PDF具有满意的逼近效果.同时, 利用式(9)对期望输出PDF数据进行权值估计, 此时获得相对应的期望权值为.
4.2 预测PDF控制效果
在完成RBF基函数参数整定同时利用式(4)对不同时刻输出PDF进行权值估计, 然后利用RVFLNs方法建立前三组权值的非线性预测模型, 利用稀释水流量、磨盘间隙和输出PDF数据, 采用所提方法建立磨浆过程输出FLSD模型, 并基于SQP算法优化式(25)设计预测PDF控制器.
本文取预测时域$N_\textrm{p}=3$, 控制时域$N_\textrm{u}=2$, 控制增量加权系数$\lambda_j=0.05$, 反馈校正系数$\beta_j=0.55$.另外, 根据实际操作经验, 输入变量稀释水流量($u_1$)和磨盘间隙($u_2$)分别满足$70\textrm{L/min}<u_1<75\textrm{L/min}$, $0.8\textrm{mm}<u_2<1.2\textrm{mm}$, $|\triangle u_1|<1$, $|\triangle u_2|<1$.根据上述分析可知, 在4个权值中有3个是相互独立的, 因此, 期望PDF相对的应权值, 则.另外, 假设初始输出PDF所对应的权值, 则$h(\pmb V_0)=1.205$, 输入变量稀释水流量($u_1$)和磨盘间隙($u_2$)的初始值分别为$u_1=71\textrm{L/min}$, $u_2=1.0\textrm{mm}$.
图 8和图 9分别在预测PDF控制器下预测权值输出响应曲线、控制输入的动态响应, 从图 8可以看出预测权值输出能够实现对期望权值的跟踪, 但权值动态模型由于非线性存在, 在一定程度上影响到预测权值输出对期望权值跟踪控制性能.此外, 图 9所示输入变量稀释水流量和磨盘间隙均能较好地稳定在可操作区间内. 图 10和图 11分别为输出FLSD的PDF3D图以及初始时刻、目标、最终时刻输出PDF, 从图 10和图 11明显能够看出在预测PDF控制器作用下, 实际输出PDF从初始输出PDF形状具有很明显逼近期望输出PDF趋势, 并最终实现对输出PDF跟踪控制.
图 10 预测PDF控制器下输出PDF3D响应Fig. 10 3D responses of the output PDF with the predictive PDF5. 结论
本文从当前磨浆过程实际控制问题出发, 针对具有典型非高斯分布特征的输出FLSD形状提出了一种预测PDF控制方法.采用迭代学习方法获得理想的RBF基函数基础上对不同时刻输出PDF相对应的权值进行估计, 针对权值之间强耦合、非线性强等特点, 采用RVFLNs建立表征输出变量和权值向量之间关系的预测模型, 最终将输出PDF的控制转化为对权值向量的控制, 基于工业数据实验结果表明了所提方法的有效性.
-
[1] Zhang Yong-Gang, Huang Yu-Long, Zhao Lin. A general framework solution of Gaussian filter with multiple step randomly delayed measurements. Acta Automatica Sinica, 2015, 41(1):122-135(张勇刚, 黄玉龙, 赵琳. 一种带多步随机延迟量测高斯滤波器的一般框架解. 自动化学报, 2015, 41(1):122-135) [2] [2] Lei M, van Wyk B J, Qi Y. Online estimation of the approximate posterior Cramer-Rao lower bound for discrete-time nonlinear filtering. IEEE Transactions on Aerospace and Electronic Systems, 2011, 47(1):37-57 [3] Zhang Yong-Gang, Huang Yu-Long, Li Ning, Zhao Lin. Conditional posterior Cramr-Rao lower bound for nonlinear sequential Bayesian estimation with one-step randomly delayed measurements. Acta Automatica Sinica, 2015, 41(3):559-574(张勇刚, 黄玉龙, 李宁, 赵琳. 带一步随机延迟量测非线性序列贝叶斯估计的条件后验克拉美罗下界. 自动化学报, 2015, 41(3):559-574) [4] [4] Arasaratnam I, Haykin S. Cubature Kalman filters. IEEE Transactions on Automatic Control, 2009, 54(6):1254-1269 [5] [5] Arasaratnam I, Haykin S. Cubature Kalman smoothers. Automatica, 2011, 47(10):2245-2250 [6] [6] Sarmavuori J, Srkk S. Fourier-Hermite Kalman filter. IEEE Transactions on Automatic Control, 2012, 57(6):1511-1515 [7] [7] Sarmavuori J, Srkk S. Fourier-Hermite Rauch-Tung-Striebel smoother. In:Proceedings of the 20th European Signal Processing Conference. Bucharest, Romania:IEEE, 2012. 2109-2113 [8] [8] Jia B, Xin M, Cheng Y. High-degree cubature Kalman filter. Automatica, 2013, 49(2):510-518 [9] [9] Jia B, Xin M. A new class of nonlinear Rauch-Tung-Striebel cubature Kalman smoothers. ISA Transactions, 2015, 55(3):72-80 [10] Zhang Yong-Gang, Huang Yu-Long, Wu Zhe-Min, Li Ning. A high order unscented Kalman filtering method. Acta Automatica Sinica, 2014, 40(5):838-848(张勇刚, 黄玉龙, 武哲民, 李宁. 一种高阶无迹卡尔曼滤波方法. 自动化学报, 2014, 40(5):838-848) [11] Jia B, Xin M, Cheng Y. Sparse-grid quadrature nonlinear filtering. Automatica, 2012, 48(2):327-341 [12] Dunk J, Straka O, imandl M. Stochastic integration filter. IEEE Transactions on Automatic Control, 2013, 58(6):1561-1566 [13] Zhang X C. A novel cubature Kalman filter for nonlinear state estimation. In:Proceedings of the 2013 IEEE 52nd Annual Conference on Decision and Control. Firenze:IEEE, 2013. 7797-7802 [14] Zhang Y G, Huang Y L, Li N, Zhao L. Embedded cubature Kalman filter with adaptive setting of free parameter. Signal Processing, 2015, 114:112-116 [15] Wang S Y, Feng J C, Tse C K. Spherical simplex-radial cubature Kalman filter. IEEE Signal Processing Letters, 2014, 21(1):43-46 [16] Chang L B, Hu B Q, Li A, Qin F J. Transformed unscented Kalman filter. IEEE Transactions on Automatic Control, 2013, 58(1):252-257 [17] Ljung L. Perspectives on system identification. In:Proceedings of the 17th IFAC World Congress. Seoul, Korea:IFAC, 2008. 7172-7184 [18] Gaperin M, Juričić D. Application of unscented transformation in nonlinear system identification. In:Proceedings of the 18th IFAC World Congress. Milano, Italy:IFAC, 2011. 4428-4433 [19] Zhao Lin, Wang Xiao-Xu, Sun Ming, Ding Ji-Cheng, Yan Chao. Adaptive UKF filtering algorithm based on maximum a posterior estimation and exponential weighting. Acta Automatica Sinica, 2010, 36(7):1007-1019(赵琳, 王小旭, 孙明, 丁继成, 闫超. 基于极大后验估计和指数加权的自适应UKF滤波算法. 自动化学报, 2010, 36(7):1007-1019) [20] Zhao L Q, Wang J L, Yu T, Jian H, Liu T J. Design of adaptive robust square-root cubature Kalman filter with noise statistic estimator. Applied Mathematics and Computation, 2015, 256:352-367 [21] Wang Lu, Li Guang-Chun, Qiao Xiang-Wei, Wang Zhao-Long, Ma Tao. An adaptive UKF algorithm based on maximum likelihood principle and expectation maximization algorithm. Acta Automatica Sinica, 2012, 38(7):1200-1210(王璐, 李光春, 乔相伟, 王兆龙, 马涛. 基于极大似然准则和最大期望算法的自适应UKF算法. 自动化学报, 2012, 38(7):1200-1210) [22] Bavdekar V A, Deshpande A P, Patwardhan S C. Identification of process and measurement noise covariance for state and parameter estimation using extended Kalman filter. Journal of Process Control, 2011, 21(4):585-601 [23] Schn T B, Wills A, Ninness B. System identification of nonlinear state-space models. Automatica, 2011, 47(1):39-49 [24] Vnnen V. Gaussian Filtering and Smoothing based Parameter Estimation in Nonlinear Models for Sequential Data[Master dissertation], Aalto University, Finland, 2012 [25] Kokkala J, Solin A, Srkk S. Expectation maximization based parameter estimation by sigma-point and particle smoothing. In:Proceedings of the 17th International Conference on Information Fusion. Salamanca, Spain:IEEE, 2014. 1-8 [26] Wu W R, Chang D C. Maneuvering target tracking with colored noise. IEEE Transactions on Aerospace and Electronic Systems, 1996, 32(4):1311-1320 [27] Srinivasan S, Aichner R, Kleijn W B, Kellermann W. Multi-channel parametric speech enhancement. IEEE Signal Processing Letters, 2006, 13(5):304-307 [28] Jamoos A, Grivel E, Bobillet W, Guidorzi R. Errors-in-variables-based approach for the identification of AR time-varying fading channels. IEEE Signal Processing Letters, 2007, 14(11):793-796 [29] Mahmoudi A, Karimi M, Amindavar H. Parameter estimation of autoregressive signals in presence of colored AR(1) noise as a quadratic eigenvalue problem. Signal Processing, 2012, 92(4):1151-1156 [30] Wang Xiao-Xu, Liang Yan, Pan Quan, Zhao Chun-Hui, Li Han-Zhou. Unscented Kalman filter for nonlinear systems with colored measurement noise. Acta Automatica Sinica, 2012, 38(6):986-998(王小旭, 梁彦, 潘泉, 赵春晖, 李汉舟. 带有色量测噪声的非线性系统Unscented卡尔曼滤波器. 自动化学报, 2012, 38(6):986-998) [31] Wang X X, Pan Q. Nonlinear Gaussian filter with the colored measurement noise. In:Proceedings of the 17th International Conference on Information Fusion. Salamanca, Spain:IEEE, 2014. 1-7 [32] Wang X X, Liang Y, Pan Q, Zhao C H, Yang F. Nonlinear Gaussian smoothers with colored measurement noise. IEEE Transactions on Automatic Control, 2015, 60(3):870-876 [33] Yuan G N, Xie Y J, Song Y, Liang H B. Multipath parameters estimation of weak GPS signal based on new colored noise unscented Kalman filter. In:Proceedings of the 2010 IEEE International Conference on Information and Automation. Harbin, China:IEEE, 2010. 1852-1856 [34] Titterton D H, Weston J L. Strapdown Inertial Navigation Technology. London:Peter Peregrinus, Ltd., 1997. [35] Goto S, Nakamura M, Uosaki K. On-line spectral estimation of nonstationary time series based on AR model parameter estimation and order selection with a forgetting factor. IEEE Transactions on Signal Processing, 1995, 43(6):1519-1522 [36] Giurcăneanu C D, Razavi S A. AR order selection in the case when the model parameters are estimated by forgetting factor least-squares algorithms. Signal Processing, 2010, 90(2):451-466 [37] Kok M, Dahlin J, Schn T B, Wills A. Newton-based maximum likelihood estimation in nonlinear state space models. In:Proceedings of the 17th IFAC Symposium on System Identification. IFAC, 2015 [38] Poyiadjis G, Doucet A, Singh S S. Particle approximations of the score and observed information matrix in state space models with application to parameter estimation. Biometrika, 2011, 98(1):65-80 [39] Li X R. Measure of nonlinearity for stochastic systems. In:Proceedings of the 15th International Conference on Information Fusion. Singapore:IEEE, 2012. 1073-1080 [40] Bucy R S, Senne K D. Digital synthesis of non-linear filters. Automatica, 1971, 7(3):287-298 [41] Dunik J, Simandl M, Straka O. Unscented Kalman filter:aspects and adaptive setting of scaling parameter. IEEE Transactions on Automatic Control, 2012, 57(9):2411-2416 [42] Einicke G A, White L B. Robust extended Kalman filtering. IEEE Transactions on Signal Processing, 1999, 47(9):2596-2599 [43] Li W L, Jia Y M. H-infinity filtering for a class of nonlinear discrete-time systems based on unscented transform. Signal Processing, 2010, 90(12):3301-3307 [44] Jia B, Xin M. Sparse-grid quadrature H_{} filter for discrete-time systems with uncertain noise statistics. IEEE Transactions on Aerospace and Electronic Systems, 2013, 49(3):1626-1636 期刊类型引用(3)
1. 李秀英，尹帅，孙书利. 传感器饱和的非线性网络化系统模糊H_∞滤波. 自动化学报. 2021(05): 1149-1158 . 
本站查看2. 刘秀翀，王占山. 系统H_∞范数计算:Lyapunov函数的直接优化方法. 自动化学报. 2019(08): 1606-1610 . 本站查看3. 张艳玲，刘爱志，孙长银. 间接互惠与合作演化的若干问题研究进展. 自动化学报. 2018(01): 1-12 . 本站查看其他类型引用(3)
-
下载:
计量
- 文章访问数: 2050
- HTML全文浏览量: 76
- PDF下载量: 954
- 被引次数: 6
下载:
