滤波辨识(6): 输出误差自回归滑动平均系统的滤波辅助模型广义增广迭代参数估计

丁锋 ,栾小丽 ,徐玲 ,张霄 ,周怡红 ,刘喜梅

(1.江南大学 物联网工程学院,江苏 无锡 214122;2.青岛科技大学 自动化与电子工程学院,山东 青岛 266061)

在辨识研究进程中,辅助模型辨识思想、多新息辨识理论、递阶辨识原理、耦合辨识概念等给我们提供了新的研究方法[1-7],与梯度方法、最小二乘方法、牛顿方法相结合用于研究不同对象,可以提出许多不同辨识方法,如辅助模型辨识方法、多新息辨识方法、递阶辨识方法、耦合辨识方法等。类似地,滤波辨识理念与梯度方法、最小二乘方法、牛顿方法相结合,可以推导出滤波梯度辨识方法、滤波最小二乘辨识方法、滤波牛顿辨识方法等。《青岛科技大学学报(自然科学版)》的连载论文[8-10]研究了递阶递推辨识方法和递阶迭代辨识方法。最近的连载论文应用滤波辨识理念,研究了有限脉冲响应滑动平均系统和方程误差自回归系统的滤波递推辨识方法和滤波迭代辨识方法[11-14],以及输出误差自回归滑动平均系统的滤波辅助模型递推辨识方法[15-16],本研究利用滤波辨识理念,研究输出误差自回归滑动平均系统的滤波辅助模型广义增广迭代辨识方法。

1 系统描述与滤波辨识模型

输出误差自回归滑动平均模型(output-error autoregressive moving average model,OEARMA模型)又称为Box-Jenkins 模型(BJ 模型)。考虑OEARMA 模型描述的动态随机系统[15-16]:

其中u(t)和y(t)分别是系统的输入和输出序列,v(t)是零均值、方差为σ2的随机白噪声序列,A(z),B(z),C(z)和D(z)是单位后移算子z-1的多项式:

设阶次na,nb,nc和nd已知。记n=na+nb+nc+nd。且设t≤0 时,所有变量的初值均为零,如y(t)=0,u(t)=0,v(t)=0。

置有关参数向量和信息向量如下:

其中系统模型输出(即无噪输出)定义为

定义滤波输出信息向量φaf(t)和φdf(t),滤波输入信息向量ϕbf(t)和ϕdf(t)如下:

uf(t)和yf(t)也可按照下列递推式计算,

式(1)两边乘以L(z)得到

这个滤波后的模型是一个输出误差模型(output-error model,OE 模型)。定义滤波系统模型输出(即滤波无噪输出):

利用式(14)～(15),式(13)可以写为向量形式:

将式(12)代入上式得到

在辨识模型(17)中,参数向量θ包含了系统的所有参数ai,bi,ci和di。这里的辨识目标是基于系统的输入输出数据u(t)和y(t),利用滤波辨识理念和辅助模型辨识思想,研究和提出估计参数向量θ的滤波辅助模型广义增广迭代辨识方法。

2 滤波辅助模型广义增广梯度迭代辨识方法

设l为辨识数据起点,L为数据长度。根据式(17),使用t=l到t=l+L-1的L组数据定义堆积输出向量Y(l,L)和堆积信息矩阵Φ(l,L),堆积滤波无噪输出信息矩阵Φaf(l,L),堆积滤波输入信息矩阵Φbf(l,L),堆积滤波输出信息矩阵Φdf(l,L)和堆积输出信息矩阵Φcf(l,L)如下:

2.1 滤波辅助模型广义增广梯度迭代辨识算法

根据辨识模型(17),定义关于参数向量θ的二次准则函数:

设λmax[X]为实对称阵X的最大特征值。令k=1,2,3,…是迭代变量,设∈ℝn是参数向量θ的第k次迭代估计,μk≥0是迭代步长(收敛因子),使用负梯度搜索,极小化J1(θ),得到梯度迭代关系:

故收敛因子μk的保守选择是

由于特征值计算很复杂,故收敛因子也可简单更保守地取为

梯度迭代算法(21)～(22)不可实现,因为右边堆积信息矩阵Φ(l,L)是未知的,故需要先构造Φ(l,L)的估计(l,L)。

2)采集观测数据u(t)和y(t),用式(44)和(47)构造信息向量φc(t)和ϕc(t),

t=1,2,…,l+L。用式(34)构造堆积输出向量Y(l,L),用式(39)造堆积输出信息矩阵Φc(l,L)。

2.2 滤波加权辅助模型广义增广梯度迭代算法

设M(L)∈ℝL×L为非负定加权矩阵。将准则函数J1(θ)修改为加权准则函数:

使用梯度搜索,极小化准则函数J′1(θ)可以得到

那么式(56)～(58)和(34)～(55)构成了加权滤波辅助模型广义增广梯度迭代算法(W-F-AM-GEGI算法),或称为滤波加权辅助模型广义增广梯度迭代算法(F-W-AM-GEGI算法)。

就得到遗忘因子滤波辅助模型广义增广梯度迭代算法(FF-F-AM-GEGI算法),或称为滤波遗忘因子辅助模型广义增广梯度迭代算法(F-FF-AM-GEGI算法)。

3 滤波辅助模型广义增广投影迭代辨识方法

令k=1,2,3,…是迭代变量,μk(L)≥0是数据长度为L时的迭代步长(收敛因子),设(L)∈ℝn是数据长度为L时参数向量θ的第k次迭代估计,使用最速下降法(负梯度搜索),极小化准则函数J1(θ),可以得到梯度递推关系:

投影算法(59)～(61)不可实现,因为右边涉及未知信息矩阵Φ(l,L)。解决的办法是用(l,L)代替式(59)～(61)右边的未知信息矩阵Φ(l,L),得到式(62),联立式(59)～(61),可以得到辨识Box-Jenkins系统参数向量θ的基于滤波的辅助模型广义增广投影迭代辨识算法(filtering-based auxiliary model generalized extended projection iterative identification algorithm),简称为滤波辅助模型广义增广投影迭代算法(filtered auxiliary model generalized extended projection-based iterative algorithm,F-AM-MI-GEPI算法):

因为计算矩阵的迹(trace)比计算特征值简单,所以收敛因子可以更保守取为

因此,可将式(62)修改为

则式(63)和(65)～(84)构成了简化的滤波辅助模型广义增广投影迭代算法(F-AM-GEPI算法),它等同于F-AM-MI-GEGI算法(117)～(137)。

4 滤波辅助模型广义增广最小二乘迭代辨识方法

4.1 滤波辅助模型广义增广最小二乘迭代辨识算法

令准则函数J1(θ)对参数向量θ的偏导数为零:

假设φ(t)是持续激励的,即ΦT(l,L)Φ(l,L)是可逆矩阵,那么从上式可以求得参数向量θ的最小二乘估计:

因为堆积信息矩阵Φ(l,L)是未知的,所以上式无法求解出参数估计。解决的办法是用(l,L)代替式(85)中未知的Φ(l,L),得到式(86)的最小二乘迭代估计,联立式(34)～(55),就得到估计Box-Jenkins系统参数向量θ的基于滤波的辅助模型广义增广最小二乘迭代辨识算法(filtering-based auxiliary model generalized extended least squares iterative identification algorithm),简称为滤波辅助模型广义增广最小二乘迭代算法(filtered auxiliary model generalized extended least squares-based iterative algorithm,F-AM-GELSI算法):

F-AM-GELSI算法(86)～(104)每迭代计算一次参数估计的步骤如下。

2)采集观测数据u(t)和y(t),用式(97)和(100)构造信息向量φc(t)和ϕc(t),

t=1,2,…,l+L。用式(87)构造堆积输出向量Y(l,L),用式(92)造堆积输出信息矩阵Φc(l,L)。

4.2 滤波加权辅助模型广义增广最小二乘迭代算法

极小化准则函数J′1(θ),令其对θ的梯度为零,可以得到

那么式(87)～(105)构成了加权滤波辅助模型广义增广最小二乘迭代辨识算法(W-F-AM-GELSI算法),或称为滤波加权辅助模型广义增广最小二乘迭代算法(F-W-AM-GELSI算法)。进一步,设0<λ≤1为遗忘因子,取加权矩阵为

就得到遗忘因子滤波辅助模型广义增广最小二乘迭代算法(FF-F-AM-GELSI算法),或称为滤波遗忘因子辅助模型广义增广最小二乘迭代算法(F-FF-AMGELSI算法)。

5 滤波辅助模型多新息广义增广梯度迭代辨识方法

5.1 滤波辅助模型多新息广义增广梯度迭代辨识算法

设p为数据窗长度(p≫n)。根据式(17),定义堆积输出向量Y(p,t)和堆积信息矩阵Φ(p,t)如下:

根据辨识模型(17),定义关于参数向量θ的二次准则函数:

令k=1,2,3,…是迭代变量,μk(t)≥0是迭代步长(收敛因子),(t)∈ℝn是参数向量θ在当前时刻t的第k次迭代估计。使用负梯度搜索,极小化准则函数J2(θ),得到梯度迭代关系:

故收敛因子μk(t)的保守选择是

由于特征值计算很复杂,故收敛因子也可简单取为

梯度迭代算法(107)～(108)不可实现,因为堆积信息矩阵Φ(p,t)是未知的,故需要先构造Φ(p,t)的估计(p,t)。

用噪声模型参数估计

构造多项式C(z)和D(z)的估计:

滤波辅助模型多新息广义增广梯度迭代算法又称滤波辅助模型移动数据窗广义增广梯度迭代算法(filtered auxiliary model moving-data-window generalized extended gradient-based iterative algorithm,F-AM-MDW-GEGI算法)。F-AM-MI-GEGI算法(117)～(137)计算参数估计(t)的步骤如下。

1)初始化:令t=1,给定移动数据窗长度p(通常取p≫n),参数估计精度ε,最大迭代次数kmax。置参数估计初值(t)=1n/p0,p0=106。

3)采集观测数据u(t)和y(t),用式(130)和(133)构造信息向量φc(t)和ϕc(t)。用式(120)构造堆积输出向量Y(p,t),用式(125)构造堆积输出信息矩阵Φc(p,t)。

8)如果k

5.2 滤波加权辅助模型多新息广义增广梯度迭代算法

设Λ(t)∈ℝp×p为非负定加权矩阵。将准则函数J2(θ)修改为加权准则函数:

使用梯度搜索,极小化准则函数J′2(θ)可以得到

那么式(138)～(140)和(120)～(137)构成了加权滤波辅助模型多新息广义增广梯度迭代算法(weighted F-AM-MI-GEGI algorithm,W-F-AM-MI-GEGI算法),或称为滤波加权辅助模型多新息广义增广梯度迭代算法(F-W-AM-MI-GEGI算法)。

6 滤波辅助模型多新息广义增广投影迭代辨识方法

令k=1,2,3,…是迭代变量,是μk(t)≥0迭代步长(收敛因子),设(t)∈ℝn是参数向量θ在当前时刻t的第k次迭代估计,使用最速下降法(负梯度搜索),极小化准则函数J2(θ),可以得到梯度递推关系:

投影算法(141)～(143)不可实现,因为右边涉及未知信息矩阵Φ(p,t)。解决的办法是用(p,t)代替式(141)～(143)右边的未知信息矩阵(p,t),得到式(144)～(146),联立式(120)～(137),可以得到辨识Box-Jenkins系统参数向量θ的基于滤波的辅助模型多新息广义增广投影迭代辨识算法 (filtering-based auxiliary model multi-innovation generalized extended projection iterative identification algorithm),简称为滤波辅助模型多新息广义增广投影迭代算法(filtered auxiliary model multi-innovation generalized extended projection-based iterative algorithm,F-AM-MI-GEPI算法):

因为计算矩阵的迹(trace)比计算特征值简单,所以收敛因子可以更保守取为

因此,可将式(144)修改为

则式(145)和(147)～(166)构成了简化的滤波辅助模型多新息广义增广投影迭代算法(F-AM-MI-GEPI算法),它等同于F-AM-MI-GEGI算法(141)～(143)。

7 滤波辅助模型多新息广义增广最小二乘迭代辨识方法

7.1 滤波辅助模型多新息广义增广最小二乘迭代辨识算法

由于最小二乘迭代算法涉及矩阵逆,所以至少要求p≥n,以及信息向量是持续激励的。令准则函数J2(θ)对参数向量θ的偏导数为零:

假设φ(t)是持续激励的,即ΦT(p,t)Φ(p,t)是可逆矩阵,那么从上式可以求得参数向量θ的最小二乘估计:

因为堆积信息矩阵Φ(p,t)是未知的,所以解决的方法是用(p,t)代替式(167)中未知的Φ(p,t),得到式(168)的最小二乘迭代估计(t),联立式(120)～(137),就得到估计Box-Jenkins系统参数向量θ的基于滤波的辅助模型多新息广义增广最小二乘迭代辨识算法(filtering-based auxiliary model multi-innovation generalized extended least squares iterative identification algorithm),简称为滤波辅助模型多新息广义增广最小二乘迭代算法(filtered auxiliary model multi-innovation generalized extended least squares-based iterative algorithm,F-AM-MI-GELSI算法):

滤波辅助模型多新息广义增广最小二乘迭代算法又称滤波辅助模型移动数据窗广义增广最小二乘迭代算法(filtered auxiliary model moving-datawindow generalized extended least squares-based iterative algorithm,F-AM-MDW-GELSI算法)。FAM-MI-GELSI算法(168)～(186)计算参数估计(t)的步骤如下。

1)初始化:令t=1,给定移动数据窗长度p(通常取p≫n),参数估计精度ε,最大迭代次数kmax(t)=1n/p0,p0=106。采集观测数据y(j)和u(j),j=1,2,…,t-1,这里l≥p为一个正整数(请读者思考l的意义)。

2)令k=1。置(t-j)=随机向 量,(t-j)=随机向量,(t-j)=随机向量,j=1,2,…,p-1。随机数初值是为保证式(168)中的矩阵可逆。

3)采集观测数据u(t)和y(t),用式(179)和(182)构造信息向量φc(t)和ϕc(t)。用式(169)构造堆积输出向量Y(p,t),用式(174)构造堆积输出信息矩阵(p,t)。

8)如果k

7.2 滤波加权辅助模型多新息广义增广最小二乘迭代算法

极小化准则函数J′2(θ),令其对θ的梯度为零,可以推导出

那么式(169)～(187)构成了加权滤波辅助模型多新息广义增广最小二乘迭代算法(weighted F-AMMI-GELSI algorithm,W-F-AM-MI-GELSI 算法),或称为滤波加权辅助模型多新息广义增广最小二乘迭代算法(F-W-AM-MI-GELSI算法)。

8 结语

针对输出误差自回归滑动平均(OEARMA)系统,即Box-Jenkins系统,研究和提出了滤波辅助模型(多新息)广义增广梯度迭代辨识方法、滤波辅助模型多新息广义增广投影迭代辨识方法、滤波辅助模型(多新息)递推广义增广梯度迭代辨识方法、滤波辅助模型(多新息)广义增广最小二乘迭代辨识方法。这些滤波辅助模型广义增广迭代辨识方法可以推广到其它有色噪声干扰下的线性和非线性多变量随机系统中[17-28]。