Identification of Linear Time-delay Systems With Unknown Delay Distributions in Its Value Range
-
摘要: 在大多数系统辨识方法中, 通常假设时变时滞在其可能的取值范围内服从均匀分布. 但是这种假设是非常受限的且在实际过程中常常无法得到满足. 因此在时滞取值概率条件未知的情况下, 针对一类线性时变时滞系统提出有效的辨识方法. 利用期望最大化(Expectation maximization, EM)算法将拟研究的辨识问题公式化, 期望最大化算法通过不断地迭代执行期望步骤和最大化步骤得到优化的参数估计. 在期望步骤中, 将未知的时变时滞当作隐含变量来处理并且假设可能的取值范围已知. 在每一个采样时刻, 时滞的变换由一个概率向量控制, 并且该向量中的每一个元素是未知的, 将其当作待估计的未知参数处理. 在算法的每次迭代过程中, 计算时滞的后验概率密度函数(Probability density function, PDF), 并在此基础上构造代价函数(Q-函数). 在最大化步骤中, 通过不断优化(Q-函数)来估计想要的参数, 包括模型参数、噪声参数、控制概率向量中的每一个元素和未知的时滞. 最后通过一个数值例子验证提出算法的有效性.Abstract: In majority system identification methods, the varying system delay is generally assumed to be uniformly distributed in its possible value range. However this assumption is very limited and it is often not satisfied in practical settings. Hence this paper addresses the identification of linear time-delay systems with unknown delay distributions in its value range. The formulation of the above identification problem is realized by using the expectation maximization (EM) algorithm which iteratively performs the expectation step (E-step) and the maximization step (M-step) until the desired optimal parameters obtained. In the E-step, the unknown varying time-delay is processed as the latent variable and its possible value range is assumed to be known as a priori. At each sampling instant, the variant of the time-delay is governed by a probability vector which is processed as the unknown parameter. In each iteration, the posterior probability density function (PDF) of the unknown delay is calculated which facilitates the construction of the cost function (Q-function). In the M-step, the calculated Q-function is optimized in order to estimate the unknown parameters including the model parameters, the noise parameter, each element in the governmental probability vector and the unknown delay. Finally, the verification tests performed on a numerical example are provided to illustrate the effectiveness of the proposed algorithm.
-
系统辨识作为一种有效的系统建模手段, 在过去的一段时间里受到国内外学者的广泛关注[1-5]. 与传统的第一原理建模方法相比, 系统辨识方法直接通过可测的系统数据提炼出与实际系统等价的数学模型, 无需深入窥探系统的内部复杂机理, 因而使得建模过程变得简单、高效[6-11]. 在系统辨识领域, 由于长时间的数据传输、数据获取方式的不同等因素, 常常导致不同的过程变量之间存在未知的时间延迟(时滞). 时滞作为一个常见的现实问题, 很容易降低辨识数据的数据质量并且对辨识算法的设计提出更高的要求[11-13]. 在辨识算法设计的过程中, 忽略时滞因素的影响会导致有偏的参数估计, 甚至错误的辨识结果[14-16].
在实际的工业过程中, 过程变量之间的时滞大体上可以分为两类: 固定时滞和时变时滞[13-18]. 对于包含固定时滞的系统辨识问题, 传统的方法一般是分为两个步骤进行: 1) 对输入−输出数据利用相关分析法, 先求解未知的时滞参数; 2) 基于时滞参数求解未知的模型参数[17-21]. 但是这样的解决方案实质上是将时滞估计和参数估计独立分开执行, 很容易造成误差累计, 进而影响辨识结果的精度[17]. 为了解决误差累计的问题, 文献[22]在统计学的框架下并基于极大似然判据提出了一类有效的迭代辨识算法. 利用期望最大化(Expectation maximization, EM)算法首先将辨识问题公式化, 将未知的时滞变量当作隐含变量来处理. 在辨识过程中, 不断计算时滞的后验概率密度函数(Probability density function, PDF), 通过比较后验概率密度函数的值来判断未知时滞的取值. 该方法的主要优点包括: 1) 给出了显性的时滞估计公式; 2) 能够同时给出时滞参数和模型参数的极大似然估计, 避免了误差累计降低辨识结果的精度. 随后该辨识思想被重新利用在线性变参数时滞系统过程中, 文献[23]采用多模型方法对线性变参数时滞系统设计有效的局部辨识算法. 首先在每一个工作点处, 选取线性自回归时滞模型为局部模型, 每一个局部模型的时滞参数各不相同; 然后利用期望最大化算法同时估计各子模型的时滞参数、模型参数以及有效宽度参数; 最后利用输出插值策略估计系统的全局输出. 在每一个采样时刻, 将系统全局输出写成每个子模型输出数据的加权组合的形式. 此外, 文献[24]结合了冗余法则以及递归最小二乘算法针对一类线性时滞系统进行了辨识研究. 核心思想是扩展原系统的维数, 然后再辨识扩展系统的模型参数. 如果某项系数趋近于零, 则认为该项在原系统中不存在, 并依据于此来判断时滞参数. 这种辨识思路无形中增加了算法的复杂度.
针对时变时滞系统辨识而言, 首先要确定未知时滞的变化机制. 在现有的文献中, 最常用的变化机制是首先假定时滞可能的取值范围, 然后假定时滞在取值范围内服从均匀分布, 即在每个采样时刻, 时滞取每个值的概率相同[3, 14, 16, 18-19]. 文献[25]针对多率采样时滞系统提出了一类有效的辨识方法. 输入、输出变量分别采用快率、慢率采样, 且输入和输出变量之间存在服从均匀分布的时变时滞. 利用期望最大化算法实现上述系统辨识问题, 将未知的时滞变量当作隐含变量来处理, 上述辨识问题可归纳为: 在慢率输出基础上同时估计时滞参数和模型参数. 另外, 基于辨识得到的有限脉冲响应模型还可以进一步促进输出丢失情况下的输出误差模型和传递函数模型的辨识问题. 文献[14]采用多模型局部辨识的思想针对线性变参数输出误差模型辨识进行深入研究. 能够同时给出各子模型的时滞参数和模型参数的估计公式, 由于目标代价函数与模型参数呈非线性关系, 采用阻尼牛顿法来迭代搜索模型参数的最优解. 此外, 文献[18]在输入输出数据双率采样的条件下考虑了线性状态空间时滞模型的辨识问题. 首先利用离散化技术得到待辨识系统的慢率采样模型; 利用数学工具将慢率采样模型转化成能观标准型; 利用卡尔曼滤波算法估计未知的状态变量; 最后利用随机梯度算法和递归最小二乘算法辨识模型参数. 在辨识过程中考虑了状态变量与输出变量之间存在未知的时变时滞, 且假定时滞在可能的取值范围内服从均匀分布.
然而, 对于时变时滞系统而言, 认为时滞在可能的取值范围内服从均匀分布这个假设通常过于严格, 在很多情况下很难得到满足. 为了解决这个问题, 本文在时滞取值概率未知的条件下研究线性时滞系统辨识问题. 仍旧假定时滞可能的取值范围先验已知, 但是时滞取每个值的概率各不相同, 即在时滞的取值范围内, 时滞取每一个值的比重各不相同. 不难看出, 时滞服从均匀分布可以看作本文研究内容的一个特例, 因此本文的研究更加贴近实际应用. 采用极大似然算法搭建辨识框架, 将时滞当作隐含变量来处理. 在每个采样时刻不断计算时滞的后验概率密度函数, 且该函数作为权重自适应地分配给每一个数据点, 因此提高时滞的估计精度也能促进模型参数估计精度的提升. 同时给出模型参数、噪声参数、控制概率向量中的每一个元素和未知的时滞的估计公式, 最后利用一个数值例子验证算法的有效性.
1. 辨识问题阐述
考虑以下线性双率采样时滞系统
$$\left\{ {\begin{aligned} & {{{\boldsymbol{x}}_k} = {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_k}} \\ &{{{\boldsymbol{y}}_{{T_i}}} = {{\boldsymbol{x}}_{{T_i} - {{\boldsymbol{\tau}} _i}}} + {{\boldsymbol{e}}_{{T_i}}}} \end{aligned}} \right.$$ (1) 其中, ${{\boldsymbol{u}}_{1:N}} = {\{ {{\boldsymbol{u}}_k}\} _{k = 1, \cdots ,N}}$表示快速率采样的输入数据且假设采样周期为${\Delta _t}$; ${{\boldsymbol{x}}_k}$表示理想的无噪声快率采样输出; ${\boldsymbol{y} _{{T_1}:{T_M}}} = {{\rm{\{ }}{\boldsymbol{y} _{{T_i}}}{\rm{\} }}_{i = 1, \cdots ,M}}$表示慢速率采样的真实输出数据, 它的采样周期假设为快率采样的整数倍, 即$f{\Delta _t}$, 其中$f$表示一个任意的正整数; ${\boldsymbol{\tau }_{1:M}} = {\{ {\boldsymbol{\tau }_i}\} _{i = 1, \cdots ,M}}$表示未知的时变时滞; ${{\boldsymbol{e}}_{{T_i}}}$表示输出测量噪声且服从零均值的高斯分布, 即${\boldsymbol{e} _{{T_i}}} \sim {\rm{N(}}{\boldsymbol{e} _{{T_i}}}|0,r)$. $\boldsymbol{G} ({z^{ - 1}})$表示系统的多项式传递函数并且有
$${\boldsymbol{G}}({z^{ - 1}}) = {g_1}{z^{ - 1}} + {g_2}{z^{ - 2}}{\rm{ + }} \cdots + {g_m}{z^{ - m}}$$ (2) 其中, $m$表示系统的阶次. 结合式(1)和式(2), 系统的模型可改写为
$${\boldsymbol{x}_k} = \boldsymbol{\phi }_k^{\rm{T}}\boldsymbol{\theta }$$ (3) 其中
$${\boldsymbol{\phi }_k} = {{\rm{[}}{\boldsymbol{u} _{k - 1}}, \;{\boldsymbol{u} _{k - 2}}, \cdots ,{\rm{ }}{\boldsymbol{u} _{k - m}}{\rm{]}}^{\rm{T}}}$$ (4) 并且待辨识的模型参数${\boldsymbol{\theta}} $定义为
$$\boldsymbol{\theta } = {[{g_1}, {g_2}, \cdots ,{\rm{ }}{g_m}]^{\rm{T}}}$$ (5) 此外, 系统输出量测过程则可改写为
$${{\boldsymbol{y}}_{{T_i}}} = {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {\tau _i}}} + {{\boldsymbol{e}}_{{T_i}}}$$ (6) 由式(6)不难看出${{\boldsymbol{y}}_{{T_i}}}$服从以下高斯分布
$${\boldsymbol{y} _{{T_i}}} \sim {\rm{N(}}{\boldsymbol{y} _{{T_i}}}|\boldsymbol{G} ({z^{ - 1}}){\boldsymbol{u} _{{T_i} - {\tau _i}}},r)$$ (7) 在本文中, 将未知的时变时滞${\tau _i}$当作隐含变量来处理. 首先假设时滞的取值范围已知, 即$[1,{\rm{ }}J]$. 在每个采样时刻, 假设时滞的取值为$[1,{\rm{ }}J]$中的任意一个整数, 但取值概率不同且未知, 用数学语言描述如下:
$$p({\boldsymbol{\tau }_i} = j) = {\boldsymbol{\beta }_j},\;\;\;0 \leq {\boldsymbol{\beta }_j} \leq 1$$ (8) 并且
$$\sum\limits_{j = 1}^J {{\boldsymbol{\beta }_j}} = 1$$ (9) 不难看出, 服从均匀分布的时变时滞系统可看作是本文的一个特例, 即
$${\boldsymbol{\beta }_1} = {\boldsymbol{\beta }_2} = \cdots = {\boldsymbol{\beta }_J} = \frac{1}{J}$$ (10) 在本文中, 慢率采样的输出数据和快率采样的输入数据构成了可用辨识数据集${\boldsymbol{C} _{{\rm{obs}}}} = {\rm{\{ }}{\boldsymbol{y} _{{T_1}:{T_M}}}, {\boldsymbol{u} _{1:N}}{\rm{\} }}$; 利用未知的时变时滞构造隐含数据集${\boldsymbol{C} _{{\rm{mis}}}} = {\rm{\{ }}{\boldsymbol{\tau }_{1:M}}{\rm{\} }}$; 待辨识的参数集构造为${\boldsymbol{\Theta }}=\{\boldsymbol{\theta }, {\boldsymbol{\beta }_{1:J}}, r{\rm{\} }}$. 因此本文的主要任务可以进一步总结为: 基于可用数据集${{\boldsymbol{C}}_{{\rm{obs}}}}$设计有效的辨识算法, 同时估计参数集${\boldsymbol{\Theta }}$和未知的时变时滞${\boldsymbol{\tau }_{1:M}}$.
2. 基于期望最大化(EM)算法的辨识算法推导
2.1 EM算法简介
在系统辨识领域, EM算法对于解决数据丢失问题或者隐含变量问题非常有效. EM算法是一个迭代优化算法, 通过不断重复执行期望步骤(E-step)和最大化步骤(M-step)以实现待辨识参数的估计. 在E-step中首先构造系统对数似然函数为
$${L_1} = \log p({{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{C}}_{{\rm{mis}}}}|{\boldsymbol{\Theta }})$$ (11) 由于隐含数据集不可用, 使得计算上述似然函数变得尤为困难. 这里通过计算其相对于隐含变量的期望来构造代价函数如下:
$$\begin{split} & {\rm{ }}Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s}) =\\ &\quad{{\rm{E}}_{{{\boldsymbol{C}}_{{\rm{mis}}}}|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s}}}\{ \log p({{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{C}}_{{\rm{mis}}}}|{\boldsymbol{\Theta }})\} =\\ &\quad\int {p({{\boldsymbol{C}}_{{\rm{mis}}}}|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} \log p({{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{C}}_{{\rm{mis}}}}|{\boldsymbol{\Theta }}){\rm{d}}{{\boldsymbol{C}}_{{\rm{mis}}}} \end{split} $$ (12) 其中, ${{\boldsymbol{\Theta }}^s}$表示在第$s$次迭代中得到的参数估计值. 在M-step中, 通过优化目标代价函数以得到新的模型参数估计值
$${{\boldsymbol{\Theta }}^{s + 1}} = \arg \mathop {\max }\limits_{\boldsymbol{\Theta }} Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s})$$ (13) 重复执行E-step和M-step, 直至参数收敛[26]. 期望最大化算法实质上是通过不断迭代优化${L_1}$的条件期望函数来代替直接优化${L_1}$函数, 最终得到所需的模型参数估计值. 算法1总结了EM算法的具体实施步骤.
算法 1. EM算法
1) 将待辨识的参数初始化${{\boldsymbol{\Theta }}^s} = {{\boldsymbol{\Theta }}^0}$且设置$s = 0$;
2) E-step: 根据式 (12) 计算目标代价函数$Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s})$;
3) M-step: 根据式 (13) 优化代价函数并得到新的参数估计值;
4) $s = s + 1$, 转回到步骤 2), 直至参数收敛[26].
2.2 基于EM算法的辨识过程推导
1) E-${\rm{step}} $
首先构造系统的对数似然函数如下:
$$\begin{split} \log {L_1} =\;& \log p({{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{C}}_{{\rm{mis}}}}|{\boldsymbol{\Theta }})= \\ &\log p({{\boldsymbol{y}}_{{T_1}:{T_M}}}, {{\boldsymbol{\tau }}_{1:M}}, {{\boldsymbol{u}}_{1:N}}{\rm{|}}{\boldsymbol{\Theta }}) =\\ &\log p({{\boldsymbol{y}}_{{T_1}:{T_M}}}{\rm{|}}{{\boldsymbol{\tau }}_{1:M}}, {{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }}) \;+ \\ &\log p({{\boldsymbol{\tau }}_{1:M}}{\rm{|}}{{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }})+ {{\boldsymbol{C}}_1} \end{split} $$ (14) 其中, ${{\boldsymbol{C}}_1} = \log p({{\boldsymbol{u}}_{1:N}}|{\boldsymbol{\Theta }})$表示一个与参数无关的常数项. 基于式 (4)和式(6) 不难看出, 在每个采样时刻, 慢采样输出${{\boldsymbol{y}}_{{T_i}}}$与之前一段时间的输入变量${{\boldsymbol{u}}_{1:{T_i} - 1}}$、未知时滞${{\boldsymbol{\tau }}_i}$以及模型参数${\boldsymbol{\Theta }}$有关, 因此等式右边的概率密度函数$\log p({{\boldsymbol{y}}_{{T_1}:{T_M}}}{\rm{|}}{{\boldsymbol{\tau }}_{1:M}}, {{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }})$可进一步化简为
$$\begin{split} &\log p({{\boldsymbol{y}}_{{T_1}:{T_M}}}{\rm{|}}{{\boldsymbol{\tau }}_{1:M}}, {{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }}) = \\ &\qquad\sum\limits_{i = 1}^M {\log p({{\boldsymbol{y}}_{{T_i}}}|{{\boldsymbol{u}}_{1:{T_i} - 1}},{{\boldsymbol{\tau }}_i},{\boldsymbol{\Theta }})} \end{split} $$ (15) 再结合式 (7), 可以进一步化简为
$$\begin{split} &\log p({{\boldsymbol{y}}_{{T_1}:{T_M}}}{\rm{|}}{{\boldsymbol{\tau }}_{1:M}}, {{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }}) =\\ &\qquad- \frac{M}{2}\log 2{{\pi }} - \frac{M}{2}\log r \;- \\ &\qquad \sum\limits_{i = 1}^M {\frac{{{{[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {{\boldsymbol{\tau }}_i}}}]}^2}}}{{2r}}} \end{split} $$ (16) 此外, 由于时滞${{\boldsymbol{\tau }}_i}$与模型参数以及输入变量无关, 因此$\log p({{\boldsymbol{\tau }}_{1:M}}{\rm{|}}{{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }})$可以进一步化简为
$$\begin{gathered} {\rm{ }}\log p({{\boldsymbol{\tau }}_{1:M}}{\rm{|}}{{\boldsymbol{u}}_{1:N}},{\boldsymbol{\Theta }}) = \sum\limits_{i = 1}^M {\log p({{\boldsymbol{\tau }}_i})} \end{gathered} $$ (17) 结合式 (16)和式(17), 系统的对数似然函数最终化简为
$$\begin{split} \log {L_1} = \;& - \sum\limits_{i = 1}^M {\frac{{{{[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {{\boldsymbol{\tau }}_i}}}]}^2}}}{{2r}}} \; + \\ & \sum\limits_{i = 1}^M {\log p({{\boldsymbol{\tau }}_i})} - \frac{M}{2}\log r + {{\boldsymbol{C}}_2} \end{split} $$ (18) 其中, ${{\boldsymbol{C}}_2} = {{\boldsymbol{C}}_1} - (M/2)\log 2\pi$是参数无关项. 此时再基于式 (12), 可求得
$$\begin{split} &Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s}) = - \frac{M}{2}\log r \;- \\ &\;\; \frac{1}{{2r}}\sum\limits_{i = 1}^M {\int {p({{\boldsymbol{\tau }}_i}|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } {[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {{\boldsymbol{\tau }}_i}}}]^2}{\rm{d}}{{\boldsymbol{\tau }}_i} \; +\\ &\;\;\sum\limits_{i = 1}^M {\int {p({{\boldsymbol{\tau }}_i}|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \log p({{\boldsymbol{\tau }}_i}){\rm{d}}{{\boldsymbol{\tau }}_i} + {{\boldsymbol{C}}_2} \\[-20pt] \end{split} $$ (19) 由于时滞变量${{\boldsymbol{\tau }}_i}$属于离散变量, 因此
$$\begin{split} &Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s}) = - \frac{M}{2}\log r \;- \\ & \frac{1}{{2r}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s}){{[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {{\boldsymbol{\tau }}_i}}}]}^2}} } \;+\\ &\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \log p({{\boldsymbol{\tau }}_i} = j)+{{\boldsymbol{C}}_2} \\[-20pt] \end{split} $$ (20) 观察式 (20) 可得, 想要优化代价函数$Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s})$, 首先需要计算时滞取值的后验概率密度函数$p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})$. 根据文献[17]中给出的计算公式可得
$$\begin{split} &p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s}) = \\ &\qquad\frac{{p({{\boldsymbol{y}}_{{T_i}}}|{{\boldsymbol{u}}_{1:{T_i} - 1}},{{\boldsymbol{\tau }}_i} = j,{{\boldsymbol{\Theta }}^s})p({{\boldsymbol{\tau }}_i} = j)}}{{\sum\limits_{j = 1}^J {p({{\boldsymbol{y}}_{{T_i}}}|{{\boldsymbol{u}}_{1:{T_i} - 1}},{{\boldsymbol{\tau }}_i} = j,{{\boldsymbol{\Theta }}^s})p({{\boldsymbol{\tau }}_i} = j)} }}= \\ &\qquad\frac{{{\rm{N}}({{\boldsymbol{y}}_{{T_i}}}|{{\boldsymbol{G}}^s}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - j}},{r^s})}}{{\sum\limits_{j = 1}^J {{\rm{N}}({{\boldsymbol{y}}_{{T_i}}}|{{\boldsymbol{G}}^s}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - j}},{r^s})} }} \end{split} $$ (21) 为了方便优化, 可将代价函数$Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s})$拆分成以下形式
$$Q({\boldsymbol{\Theta }}|{{\boldsymbol{\Theta }}^s}) = - {Q_1}({\boldsymbol{\theta }},r) + {Q_2}({{\boldsymbol{\beta }}_j}) + {{\boldsymbol{C}}_2}$$ (22) 其中,
$$\begin{split} {Q_1}({\boldsymbol{\theta }},r) =\;& \frac{1}{{2r}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {\Big\{ p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \;\times \\ &{[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{G}}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - {{\boldsymbol{\tau }}_i}}}]^2}\Big\} {\rm{ + }}\frac{M}{2}\log r \end{split} $$ (23) $${Q_2}({{\boldsymbol{\beta }}_j}) = \sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \log {{\boldsymbol{\beta }}_j}$$ (24) 2) M-step
在此步骤中通过优化上述代价函数来得到参数估计. 首先${Q_1}({\boldsymbol{\theta }},r)$属于一个复合函数, 它同时是模型参数${\boldsymbol{\theta }}$以及噪声参数$r$的函数, 只能采取循环优化的方式得到优化的参数估计. 当优化${\boldsymbol{\theta }}$时, 令$r = {r^s}$为定值; 反之, 当优化$r$时, 令${\boldsymbol{\theta }}={{\boldsymbol{\theta }}^s}$为定值. 将函数${Q_1}({\boldsymbol{\theta }},r)$对$r$求偏导并令导数等于零, 得
$$\begin{split} &\frac{{\partial {Q_1}({\boldsymbol{\theta }},r)}}{{\partial r}} = 0{\rm{ }} \Rightarrow {r^{s + 1}} = \\ &\qquad\frac{1}{M}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {\Big\{ p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} }\;\times \\ &\qquad {[{{\boldsymbol{y}}_{{T_i}}} - {{\boldsymbol{G}}^s}({z^{ - 1}}){{\boldsymbol{u}}_{{T_i} - j}}]^2}\Big\} \end{split} $$ (25) 由于模型参数${\boldsymbol{\theta }}$与代价函数${Q_1}({\boldsymbol{\theta }},r)$并非呈线性关系, 因此无法直接利用导数的方式求得${\boldsymbol{\theta }}$的估计公式. 结合式 (3)、(4)和(23), 可将代价函数${Q_1}({\boldsymbol{\theta }},r)$改写为
$$\begin{split} {Q_1}({\boldsymbol{\theta }},r) =\;& \frac{1}{{2r}}\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {\Big\{ p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \;\times\\ & {[{{\boldsymbol{y}}_{{T_i}}} - {\boldsymbol{\phi }}_{_{{T_i} - j}}^{\rm{T}}{\boldsymbol{\theta }}]^2}\Big\} +\frac{M}{2}\log r \end{split} $$ (26) 其中, 向量${{\boldsymbol{\phi }}_{{T_i} - j}}$定义为
$${{\boldsymbol{\phi }}_{{T_i} - j}} = {[{{\boldsymbol{u}}_{{T_i} - j - 1}}, {{\boldsymbol{u}}_{{T_i} - j - 2}}, \cdots ,{\rm{ }}{{\boldsymbol{u}}_{{T_i} - j - m}}]^{\rm{T}}}$$ (27) 此时, 函数${Q_1}({\boldsymbol{\theta }},r)$与模型参数${\boldsymbol{\theta }}$呈线性关系. 再进一步对${Q_1}({\boldsymbol{\theta }},r)$求偏导可得
$$\begin{split} &\frac{{\partial {Q_1}({\boldsymbol{\theta }},r)}}{{\partial {\boldsymbol{\theta }}}} = 0{\rm{ }} \Rightarrow {{\boldsymbol{\theta }}^{s + 1}} = \\ & \quad \frac{{\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s}){{\boldsymbol{\phi }}_{{T_i} - j}}{{\boldsymbol{y}}_{{T_i}}}} } }}{{\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s}){{\boldsymbol{\phi }}_{{T_i} - j}}{\boldsymbol{\phi }}_{{T_i} - j}^{\rm{T}}} } }} \end{split} $$ (28) 此外, 根据函数${Q_2}({{\boldsymbol{\beta }}_j})$求解${{\boldsymbol{\beta }}_j}$的迭代估计公式, 需要构造拉格朗日函数. 首先假设${L_{\boldsymbol{\beta }}}$为拉格朗日因子, 结合式 (9) 可构造所需的拉格朗日函数如下:
$$\begin{split} f({{\boldsymbol{\beta }}_j},{L_{\boldsymbol{\beta }}}) =\;& \sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } \log {{\boldsymbol{\beta }}_j}\;+ \\ & {L_{\boldsymbol{\beta }}}\left(\sum\limits_{j = 1}^J {{{\boldsymbol{\beta }}_j}} - 1\right)\\[-15pt] \end{split} $$ (29) 通过求解上述拉格朗日函数可得
$$\begin{split} &\left\{ {\begin{aligned} & {\frac{{\partial f({{\boldsymbol{\beta }}_j},{L_{\boldsymbol{\beta }}})}}{{\partial {{\boldsymbol{\beta }}_j}}} = 0} \\ & {\frac{{\partial f({{\boldsymbol{\beta }}_j},{L_{\boldsymbol{\beta }}})}}{{\partial {L_{\boldsymbol{\beta }}}}} = 0} \end{aligned}} \right.{\rm{ }} \Rightarrow \\ &\qquad{\boldsymbol{\beta }}_j^{s + 1} = \frac{{\sum\limits_{i = 1}^M {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} }}{{\sum\limits_{i = 1}^M {\sum\limits_{j = 1}^J {p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})} } }} \end{split} $$ (30) 并且在每个采样时刻, 未知时变时滞的迭代估计公式如下:
$${\boldsymbol{\tau }}_i^{s + 1} =\arg \mathop { \max }\limits_{j \in [1,{\rm{ }}J]} p({{\boldsymbol{\tau }}_i} = j|{{\boldsymbol{C}}_{{\rm{obs}}}},{{\boldsymbol{\Theta }}^s})$$ (31) 至此, 基于期望最大化算法的未知参数和未知时变时滞估计公式全部推导完成. 值得注意的是, 在参数循环迭代的过程中, 只用到了前一次的模型参数${{\boldsymbol{\theta }}^s}$和噪声参数${r^s}$. 因此在执行算法初始化时, 只需设置${\boldsymbol{\theta }}$和$r$的初值. 本文提出的辨识算法可总结如算法2所示.
算法 2. 时滞取值概率未知下的线性时滞系统辨识算法
1) 模型参数初始化${{\boldsymbol{\Theta }}^0}{\rm{ = \{ }}{{\boldsymbol{\theta }}^0}, {r^0}{\rm{\} }}$, 并且设定$s = 0$;
2) 根据式 (21) 计算每个采样时刻时滞可能取值的后验概率密度函数;
3) 根据式 (22) ~ (24) 构造目标代价函数;
4) 根据式 (25) 计算噪声方差的迭代估计值;
5) 根据式 (28) 计算模型参数的估计值;
6) 基于拉格朗日方程并且根据式 (30) 计算时滞取值概率的值;
7) 根据式 (31) 不断估计未知的时变时滞;
8) 如果${\rm{||}}\frac{{{{\boldsymbol{\Theta }}^{s + 1}} - \;{{\boldsymbol{\Theta }}^s}}}{{{{\boldsymbol{\Theta }}^s}}}{\rm{|}}{{\rm{|}}_2} \leq {\rm{1}}{{\rm{0}}^{ - 3}}$, 结束此次循环, 否则$s = s + 1$并且回到步骤 2), 最终得到的模型参数估计表示为$\substack{{{\boldsymbol{\Theta }}^*}{\rm{ = \{ }}{{\boldsymbol{\theta }}^ * }, {{\boldsymbol{\beta }}^ * }_{1:J}, {r^ * }{\rm{\} }}}$.
3. 辨识算法验证
3.1 算法性能验证
考虑以下双率采样的线性时滞系统
$$ \begin{split}{{\boldsymbol{x}}_k} =\;& 1.5{{\boldsymbol{u}}_{k - 1}} - 0.7{{\boldsymbol{u}}_{k - 2}} \;+\\ &{{\boldsymbol{u}}_{k - 3}} + 0.5{{\boldsymbol{u}}_{k - 4}} \end{split}$$ (32) $${{\boldsymbol{y}}_{{T_i}}} = {{\boldsymbol{x}}_{{T_i} - {{\boldsymbol{\tau}} _i}}} + {{\boldsymbol{e}}_{{T_i}{\rm{ }}}}\;\;\;\;\qquad$$ (33) 其中, 输入信号选取为幅值在 [−1, 1]之间的高斯白噪声, 即${{\boldsymbol{u}}_k} = - 1 + 2\times{\rm{N}}(0,1)$. 输出数据采用慢率采样且采样周期是输入的5倍, 即$f = 5$. 在每个慢率采样时刻, 时滞取值范围分布在[1, 4]中, 并且按照概率[0.35, 0.2, 0.25, 0.2]取值, 即
$$[{{\boldsymbol{\beta }}_1},{\rm{ }}{{\boldsymbol{\beta }}_2},{\rm{ }}{{\boldsymbol{\beta }}_3},{\rm{ }}{{\boldsymbol{\beta }}_4}] = [0.35,{\rm{ 0}}{\rm{.2, 0}}{\rm{.25, 0}}{\rm{.2}}]$$ (34) 输出测量噪声选取为方差$r = 0.01$的高斯噪声, 且输出噪声与输入信号相互独立. 对于输入数据生成$N = 2\;000$个数据点, 则收集到的辨识数据如图1所示.
在执行辨识算法之前, 首先需要将待辨识的参数初始化. 模型参数${\boldsymbol{\theta}} $的初始值选取为
$${{\boldsymbol{\theta }}^0}= [0{\rm{.7, - 0}}{\rm{.2, 0}}{\rm{.5, 0}}{\rm{.1]}}^{\rm{T}}$$ (35) 噪声方差$r$的初始值选取为
$${r^0} = 0.05$$ (36) 然后执行算法2, 经过50次迭代之后, 得到的时滞分布向量估计值如下:
$$\begin{gathered} {\rm{ }}[{\boldsymbol{\beta }}_1^{\rm{*}}, {\boldsymbol{\beta }}_2^{\rm{*}}, {\boldsymbol{\beta }}_3^{\rm{*}}, {\boldsymbol{\beta }}_4^{\rm{*}}] = [0.3300,{\rm{ 0}}{\rm{.1943, 0}}{\rm{.2625, 0}}{\rm{.2098}}] \end{gathered} $$ (37) 同时, 得到模型参数的收敛曲线如图2所示.
通过图2可以看出, 本文提出的辨识算法具有非常良好的收敛性, 经过少量的迭代之后, 每个模型参数都能精确地收敛到真实值附近. 此时, 未知的噪声参数迭代曲线如图3所示.
此外, 结合式 (31), 得到未知时滞的估计结果如图4所示.
在图4中, 蓝色的实线表示真实的时滞; 红色的虚线表示估计的时滞. 经过统计, 时滞估计的正确率在91% 以上. 由此可见, 本文提出的辨识算法能够同时估计模型参数以及未知的时变时滞.
为了更加系统地验证本算法的有效性, 本文还设计了蒙特卡洛模拟仿真实验对算法进行更加系统的测试. 在接下来的测试中, 运行100次独立的蒙特卡洛仿真, 迭代次数设置为50次. 与文献[1]中类似, 每次独立的蒙特卡洛仿真的初始值都从以真实参数为中心的对称区间中选取未知的初值${{\boldsymbol{\theta }}^0}$. 得到的蒙特卡洛辨识结果如图 5 ~ 7所示.
图 5 参数${g_1}$和${g_2}$的蒙特卡洛辨识结果Fig. 5 The identified${g_1}$and${g_2}$in Monte Carlo simulations
图 6 参数${g_3}$和${g_4}$的蒙特卡洛辨识结果Fig. 6 The identified${g_3}$and${g_4}$in Monte Carlo simulations此外, 本文还在不同的噪声条件下测试了算法的性能. 首先定义信噪比 (Signal-to-noise rate, SNR) 的计算式为
$$SNR = 10\log \frac{{{\rm{var}} ({{\boldsymbol{y}}_k})}}{{{\rm{var}} ({{\boldsymbol{e}}_k})}}$$ (38) 其中, ${\rm{var}} ( \cdot )$表示方差操作符. 分别在信噪比等于$10\;{\rm{dB}}$, $15\;{\rm{dB}}$, $20\;{\rm{dB}}$, $25\;{\rm{dB}}$, $30\;{\rm{dB}}$时执行算法2, 并且得到辨识结果, 如图8、图9和表1所示.
图 8 不同噪声下的${g_1}$和${g_2}$估计结果Fig. 8 The identified${g_1}$and${g_2}$under different noise levels
图 9 不同噪声下的${g_3}$和${g_4}$估计结果Fig. 9 The identified${g_3}$and${g_4}$under different noise levels表 1 不同信噪比下的未知时滞估计精度Table 1 The estimation accuracy of the unknown delays under different SNRs信噪比 (dB) 时滞匹配精度 (%) 10 73.00 15 81.25 20 88.75 25 93.25 30 97.50 从上述辨识结果可以得出以下结论:
1) 从图2 ~ 4可以看出, 本文提出的辨识算法能够同时给出模型参数、噪声方差和时变时滞的估计值. 经过少数迭代之后, 模型参数和噪声参数可以精确地收敛到真实值; 同时时滞的估计值和真实值之间的匹配精度达到91%以上, 证明了所提算法的有效性.
2) 从图5 ~ 7可以看出, 在蒙特卡洛仿真中, 本文提出的辨识算法依然稳定, 在不同的初始值条件下得到的模型参数估计以及时滞参数估计都能精确地收敛到真实值, 这进一步证明了所提算法的稳定性. 对于迭代过程中出现的数值差异, 是由于参数初始值不同而造成的.
3) 从图8、图9和表1可以看出, 当信噪比不断增大时, 意味着噪声不断减少且数据质量不断变好, 模型参数和时滞参数的估计精度不断提高, 并且随着数据质量不断提升, 模型参数收敛到真实值的速度也不断加快. 上述辨识结果也再次验证了所提算法的有效性.
3.2 对比实验验证
为了进一步验证算法2的有效性, 本文还设计了比较实验. 将算法2分别与两个基于递推最小二乘 (Recursive least squares, RLS) 辨识方法进行比较. 在第1种方法中, 假设时滞参数精确已知, 将这种算法命名为 RLS-accurate; 在第2种方法中, 忽略未知的时滞, 将这种算法命名为 RLS-ignore. 输出量测噪声采用信噪比为$20\;{\rm{dB}}$的高斯白噪声, 通过执行上述三种辨识算法得到的结果如表2所示.
表 2 对比实验结果Table 2 The identification results of the comparison tests参数 真实值 RLS-accurate 算法2 RLS-ignore ${g_1}$ 1.5 1.4982 1.4861 1.1553 ${g_2}$ −0.7 −0.7025 −0.6969 −0.4136 ${g_3}$ 1.0 1.0092 1.0146 0.8348 ${g_4}$ 0.5 0.5052 0.4941 0.5256 由上述辨识结果可以看出:
1) 本文提出的辨识算法在参数估计精度上与RLS-accurate方法相当. 再结合表1的辨识结果可以进一步看出, 本文提出的辨识算法在估计时滞参数的同时, 也能够精确地给出模型参数的估计结果.
2) 由RLS-ignore方法的辨识结果可以看出, 如果忽略时滞因素的影响会造成有偏的参数估计, 这从侧面证明了本文提出的算法2的有效性.
4. 结束语
针对时滞取值概率未知的情况, 本文基于期望最大化算法提出了一类有效的线性时滞系统辨识算法. 将时滞的取值概率当成未知的参数来处理, 在辨识过程中所提出的算法能够有效地同时估计模型参数、噪声参数和时滞参数. 仿真结果表明, 本文提出的辨识算法具有快速的收敛速度以及较强的稳定性. 基于本文的研究内容, 未来的研究工作可注重以下两方面:
1) 跳过高斯假设, 在更加复杂的数据条件下, 例如在非高斯噪声、自相关噪声以及量化的输出条件下研究时滞系统辨识方法;
2) 进一步考虑相邻时刻时滞间的相关性, 对于包含相关时滞 (例如: 马尔科夫时滞) 的系统辨识问题进行深入研究. 此方向的研究甚至可以扩展到马尔科夫切换系统辨识方法的研究.
-
图 5 参数${g_1}$和${g_2}$的蒙特卡洛辨识结果
Fig. 5 The identified${g_1}$and${g_2}$in Monte Carlo simulations
图 6 参数${g_3}$和${g_4}$的蒙特卡洛辨识结果
Fig. 6 The identified${g_3}$and${g_4}$in Monte Carlo simulations
图 8 不同噪声下的${g_1}$和${g_2}$估计结果
Fig. 8 The identified${g_1}$and${g_2}$under different noise levels
图 9 不同噪声下的${g_3}$和${g_4}$估计结果
Fig. 9 The identified${g_3}$and${g_4}$under different noise levels
表 1 不同信噪比下的未知时滞估计精度
Table 1 The estimation accuracy of the unknown delays under different SNRs
信噪比 (dB) 时滞匹配精度 (%) 10 73.00 15 81.25 20 88.75 25 93.25 30 97.50 
下载: 导出CSV
表 2 对比实验结果
Table 2 The identification results of the comparison tests
参数 真实值 RLS-accurate 算法2 RLS-ignore ${g_1}$ 1.5 1.4982 1.4861 1.1553 ${g_2}$ −0.7 −0.7025 −0.6969 −0.4136 ${g_3}$ 1.0 1.0092 1.0146 0.8348 ${g_4}$ 0.5 0.5052 0.4941 0.5256
下载: 导出CSV
-
[1] 黄玉龙, 张勇刚, 李宁, 赵琳. 一种带有色量测噪声的非线性系统辨识方法. 自动化学报, 2015, 41(11): 1877--1892.Huang Yu-Long, Zhang Yong-Gang, Li Ning, Zhao Lin. An identification method for nonlinear systems with colored measurement noise. Acta Automatica Sinica, 2015, 41(11): 1877--1892. [2] 段广全, 孙书利. 带未知模型参数和衰减观测率系统自校正分布式融合估计. 自动化学报, 2021, 47(2): 423--431.Duan Guang-Quan, Sun Shu-Li. Self-tuning distributed fusion estimation for systems with unknown model parameters and fading measurement rates. Acta Automatica Sinica, 2021, 47(2): 423--431. [3] Chen J, Huang B, Ding F, Gu Y. Variational Bayesian approach for ARX systems with missing observations and varying time-delays. Automatica, 2018, 94: 194--204. doi: 10.1016/j.automatica.2018.04.003 [4] Ding F, Liu X M, Hayat T. Hierarchical least squares identification for feedback nonlinear equation-error systems. Journal of the Franklin Institute, 2020, 357: 2958--2977. doi: 10.1016/j.jfranklin.2019.12.007 [5] 李沛, 孙书利. 阳春华, 贺建军, 桂卫华. 基于影子趋势对比的矿热炉炉况在线辨识及趋势预测. 自动化学报, 2020, 41(x): 1--12doi: 10.16383/j.aas.c190827.Li Pei, Yang Chun-Hua, He Jian-Jun, Gui Wei-Hua. Smelting Condition Identification and Prediction for Submerged Arc Furnace Based on Shadow trend comparison. Acta Automatica Sinica, 2020, 41(x): 1−12doi: 10.16383/j.aas.c190827. [6] Lindfors M, Chen T S. Regularized LTI system identification in the presence of outliers: A variational EM approach. Automatica, 2020, 121: 1--13. [7] Liu X, Liu X F. An improved identification method for a class of time-delay systems. Signal Processing, 2020, 167: 1--9. [8] Battilotti S, d’Angelo M. Stochastic output delay identification of discrete-time Gaussian systems. Automatica, 2019, 109: 1--6. [9] Sammaknejad N, Zhao Y J, Huang B. A review of the Expectation Maximization algorithm in data-driven process identification. Journal of Process Control, 2019, 73: 123--136. doi: 10.1016/j.jprocont.2018.12.010 [10] Zhao W X, Yin G, Bai E W. Sparse system identification for stochastic systems with general observation sequences. Automatica, 2020, 121: 1--13. [11] Chen J, Shen Q Y, Liu Y J, Wan L J. Expectation maximization identification algorithm for time-delay two-dimensional systems. Journal of the Franklin Institute, 2020, 37: 9992--10009. [12] Xia W F, Zheng W X, Xu S Y. Extended dissipativity analysis of digital filters with time delay and Markovian jumping parameters. Signal Processing, 2018, 152: 247--254. doi: 10.1016/j.sigpro.2018.06.004 [13] 刘青松. 基于嵌套-伪预估器反馈的时滞控制系统输入时滞补偿. 自动化学报, 2020, 45(x): 1--8doi: 10.16383/j.aas.c190830.Liu Qing-Song. Nested-pseudo predictor feedback based input delay compensation for time-delay control systems. Acta Automatica Sinica, 2020, 45(x): 1−18doi: 10.16383/j.aas.c190830. [14] Yang X Q, Yang X B. Local identification of LPV dual-rate system with random measurements delays. IEEE Transactions on Industrial Electronics, 2017, 62(2): 1499--1507. [15] Li X D, Yang X Y, Song S J. Lyapunov conditions for finite-time stability of time-varying time-delay systems. Automatica, 2019, 103: 135--140. doi: 10.1016/j.automatica.2019.01.031 [16] Ding F, Wang X H, Mao L, Xu L. Joint state and multi-innovation parameter estimation for time-delay linear systems and its convergence based on the Kalman filtering. Digital Signal Processing, 2017, 62: 211--223. doi: 10.1016/j.dsp.2016.11.010 [17] Yang X Q, Yin S. Robust global identification and output estimation for LPV dual-rate systems subjected to random output time-delays. IEEE Transactions on Industrial Informatics, 2017, 13(6): 2876--2885. doi: 10.1109/TII.2017.2702754 [18] Chen L, Han L L, Huang B, Liu F. Parameter estimation for a dual-rate system with time delay. ISA Transactions, 2014, 53(5): 1368--1376. doi: 10.1016/j.isatra.2014.01.001 [19] Ma J X, Chen J, Xiong W L, Ding Feng. Expectation maximization estimation algorithm for Hammerstein models with non-Gaussian noise and random time delay from dual-rate sampled-data. Digital Signal Processing, 2018, 73: 135--144. doi: 10.1016/j.dsp.2017.11.009 [20] Xu L, Ding F, Gu Y, Alsaedi A, Hayat T. A multi-innovation state and parameter estimation algorithm for a state space system with d-step state-delay. Signal Processing, 2017, 140: 97--103. doi: 10.1016/j.sigpro.2017.05.006 [21] Kodamana H, Huang B, Ranjan R, Zhao Y J, Tan R M, Sammaknejad N. Approaches to robust process identification: A review and tutorial of probabilistic methods. Journal of Process Control, 2018, 66: 68--83. doi: 10.1016/j.jprocont.2018.02.011 [22] Yang X Q, Karimi H R. Identification of LTI time-delay systems with missing output data using GEM algorithm. Mathematical Problems in Engineering, 2014, 15: 1−8 [23] Yang X Q, Gao H J. Multiple model approach to linear parameter varying time-delay system identification with EM algorithm. Journal of the Franklin Institute, 2014, 351(12): 5565--5581. doi: 10.1016/j.jfranklin.2014.09.015 [24] Chen J, Ma J X, Liu Y J, Ding F. Identification methods for time-delay systems based on the redundant rules. Signal Processing, 2017, 137: 192--198. doi: 10.1016/j.sigpro.2017.02.006 [25] Xie L, Yang H Z, Huang B. FIR model identification of multirate processes with random delays using EM algorithm. AIChE Journal, 2013, 59(11): 4124--4132. doi: 10.1002/aic.14147 [26] Wu C. On the convergence properties of the EM algorithm. The Annals of Statistics, 1983, 11(1): 95--103. doi: 10.1214/aos/1176346060 期刊类型引用(3)
1. 刘切，李俊豪，王浩，曾建学，柴毅. 不确定性环境下维纳模型的随机变分贝叶斯学习. 自动化学报. 2024(06): 1185-1198 . 
本站查看2. 刘鑫，陈强，王兰豪，代伟. 非对称偏斜噪声条件下一种鲁棒概率系统辨识算法研究. 自动化学报. 2024(10): 2022-2035 . 本站查看3. 刘鑫，海洋，代伟. 基于厚尾双学生氏t分布的非线性状态空间系统鲁棒辨识方法. 电子学报. 2024(09): 3052-3064 . 百度学术其他类型引用(7)
-
下载:
计量
- 文章访问数: 1130
- HTML全文浏览量: 311
- PDF下载量: 195
- 被引次数: 10
