1、第五讲 处理有色噪声扰动的LS类方法(1/3),第五讲 有色噪声扰动的最小二乘类方法由上一讲中所讨论的LS估计的统计特性可知,在动态系统的参数估计中,LS法仅当随机扰动为零均值的白噪声时,才能保证得到无偏一致的估计值.在实际工程和社会系统的辨识中,随机扰动w(k)只是各种系统内外扰动和结构建模误差等因素的综合反映.因此,并不能保证随机扰动为统计独立的白噪声.故,实际系统辨识时,并不能保证LS估计值为无偏一致估计值.,第五讲 处理有色噪声扰动的LS类方法(2/3),从本讲开始,我们将讨论具有随机扰动为相关性扰动(有色噪声)的动态系统的无偏一致的参数估计问题,具体的LS类方法有增广最小二乘法(Ex
2、tended Least-square Method, ELS)广义最小二乘法(Generalized Least-square Method, GLS)辅助变量法(Instrumental Variables Method, IV)等.其它无偏参数估计算法还有多步最小二乘法偏差补偿最小二乘法,第五讲 处理有色噪声扰动的LS类方法(3/3),这些不同的处理有色噪声的辨识方法主要是针对不同的有色噪声的特性、有色噪声的不同模型表达、以及不同的辨识要求提出的.下面就分别讨论:增广最小二乘法(ELS)广义最小二乘法(GLS)辅助变量法(IV)其它方法可参阅其它教材和文献.,1. ELS法(1/19),
3、1. 增广最小二乘法考虑如下SISO随机离散系统A(z-1)y(k)=B(z-1)u(k)+v(k) (1)其中y(k),u(k)和v(k)分别为系统的输出、输入和相关的随机扰动;,对所考虑的模型的相关随机扰动v(k),一般假定其为平稳相关序列.,1. ELS法(2/19),由第二讲中的关于有色噪声的结论1和假设2可知,平稳的相关扰动v(k)可被建模如下v(k)=C(z-1)w(k) (2)其中w(k)为白噪声序列;C(z-1)为未知的、稳定的、有限阶的线性滤波器,并可表示为如下首一多项式,因此,由系统模型(1)和噪声模型(2),可得如下表示A(z-1)y(k)=B(z-1)u(k)+C(z-
4、1)w(k),1. ELS法(3/19),当假定噪声w(k)可测已知时,上式又可表示为如下自回归方程y(k)=t(k-1)q+w(k) (3)其中,A(z-1)y(k)=B(z-1)u(k)+C(z-1)w(k),1. ELS法(4/19),当上述观测数据向量(k-1)精确已知时,利用前面讨论的成批或RLS法可求得回归参数向量的LS估计值.但是,实际上上述数据向量(k-1)中包含有不可测的噪声量w(k-1),.,w(k-nc),因此对自回归式(3)并不能直接套用LS估计方法.为此,引入通过在递推参数估计过程中在线估计噪声w(k)以实现模型参数在线递推估计的ELS法.,1. ELS法(5/19)
5、,ELS法的思想就是:,利用该参数估计值来在线估计白噪声w(k-i)的值 (k-i)以替代数据向量(k-1)中的白噪声w(k-i),然后进行下一步的参数估计.,在递推估计过程中,假设当前或前一步的在线参数估计值已相当程度可用的前提下,1. ELS法(6/19),噪声w(k)的具体的估计算法是如下的事后估计或事前估计算法:,其中(k-1)为数据向量(k-1)的如下在线估计值,1. ELS法(7/19),因此,基于渐消记忆的RLS估计算法可得如下递推的ELS法,其中加权因子l为1时就为普通的ELS法.上述ELS辨识实际上是的参数估计和噪声wk的估计交替进行,计算顺序为:,1. ELS法(8/19)
6、,上述分析过程表明,ELS法是LS法的一种简单推广.它只是扩充了参数向量q和数据向量(k)的维数,也同时辨识噪声模型.就这种意义上说,可称之为ELS法.值得指出的是,ELS法虽然可同时得到噪声模型的参数估计,但其收敛过程却比过程模型A(z-1)和B(z-1)的估计值的收敛慢许多.从实用角度来说,噪声模型C(z-1)的阶次不宜取得太高.,1. ELS法(9/19)-ELS法计算步骤,综上所述,ELS法的基本计算步骤可总结如下确定被辨识系统模型的结构,以及多项式A(z-1)、B(z-1)和C(z-1)的阶次;设定递推参数初值(0),P(-1),w(0);采样获取新的观测数据y(k)和u(k),并组
7、成观测数据向量(k-1);用式(5)(7)所示的ELS法计算当前参数递推估计值;用(4)式计算白噪声w(k)的事后或事前在线估计值w(k);循环次数k加1,然后转回到第3步继续循环.,1. ELS法(10/19)-ELS法仿真程序,下面给出针对随机线性离散系统(XARMA),的ELS在线辨识仿真程序伪代码.,1. ELS法(11/19)-ELS法仿真程序),1. ELS法(12/19)-ELS法仿真程序,1. ELS法(13/19)-ELS法仿真程序,也可采用事前估计,1. ELS法(14/19)例1,例1 考虑如图下所示的仿真对象,图中,通过控制w值来改变数据的噪信比.辨识中,选择如下模型结
8、构y(k)+a1y(k-1)+a2y(k-2)=b1u(k-1)+b2u(k-2)+w(k) +c1w(k-1)+c2w(k-2) (8)其它条件同第四讲中的例1.结果如下表所示.,w(k)是服从均值为零,方差为1的正态分布的不相关随机噪声;输入信号u(k)采用伪随机二进制序列;,1. ELS法(15/19),表1 计算机仿真结果(噪信比=23%,数据组数1000),1. ELS法(16/19),递推辨识过程的辨识值如下图所示,遗忘因子=1时递推辨识结果(1),噪声估计误差,1. ELS法(17/19),遗忘因子=1时递推辨识结果(2),参数估计误差的平方和,噪声估计误差,1. ELS法(18
9、/19),遗忘因子=0.98时递推辨识结果(1),噪声估计误差,1. ELS法(19/19),遗忘因子=0.98时递推辨识结果(2),噪声估计误差,2. GLS法(1/3),2. 广义最小二乘法上一讲讨论了相关随机扰动v(k)可用C(z-1)w(k)建模情况下的参数估计问题.但有些相关扰动用C(z-1)w(k)来建模的话,线性滤波器C(z-1)的阶次相当高,这加大了参数估计的工作量,也极大地影响了建模的精度和使用上的困难性.在此情况下,相关扰动v(k)可用如下方式建模v(k)=w(k)/D(z-1) (1)其中w(k)为零均值的白噪声;,2. GLS法(2/3),1/D(z-1)为未知的、稳定
10、的、有限阶的线性滤波器,且D(z-1)可表示为如下稳定的、阶次较低的首一多项式,因此,具有该类相关随机扰动的随机离散系统的数学模型为:A(z-1)y(k)=B(z-1)u(k)+v(k)=B(z-1)u(k)+w(k)/D(z-1) (2),2. GLS法(3/3),广义最小二乘(Generalized Least-squares, GLS)法讨论的是可用模型(2)表示的随机离散系统的参数估计问题.下面,将分别讨论成批型和递推型GLS法的思想和算法,以及GLS仿真GLS评述,一、成批型GLS法(1/9),一、成批型GLS法成批型GLS法的基本思想是:先预选定一个线性滤波器Df(z-1)将模型中
11、的输出y(k)、输入u(k)和有色噪声w(k)/D(z-1)白色化,即将模型A(z-1)y(k)=B(z-1)u(k)+w(k)/D(z-1)两边左乘线性滤波器Df(z-1)后,有A(z-1)Df(z-1)y(k)B(z-1)Df(z-1)u(k)+w(k)并记为A(z-1)yf(k)B(z-1)uf(k)+w(k)再对上述白色化后的模型中多项式A(z-1)和B(z-1)的系数进行LS估计;,一、成批型GLS法(2/9),然后基于对A(z-1)和B(z-1)的估计,利用模型A(z-1)y(k)=B(z-1)u(k)+v(k)来计算模型残差(相关噪声)v(k)的估计值再基于噪声模型D(z-1)v
12、(k)=w(k)和相关噪声v(k)的估计值,利用LS法辨识D(z-1),并基于此修正选定的滤波器Df(z-1);如此循环往复地进行Df(z-1)的修正及A(z-1)和B(z-1)的迭代估计,直到模型的估计残差满足给定的精度或者Df(z-1)、A(z-1)和B(z-1)的变化量相对较小为止.,一、成批型GLS法(3/9)计算步骤,上述成批型GLS法基本思想可用如下流程图图示.,综上所述,成批型的GLS法的算法和计算步骤如下:Step 1. 确定被辨识系统模型的结构,以及多项式A(z-1)、B(z-1)和D(z-1)的阶次;Step 2. 选定稳定的初始滤波器Df(z-1);,一、成批型GLS法(
13、4/9)计算步骤,Step 3. 采样获取新的观测数据y(k)和u(k);Step 4. 基于滤波器Df(z-1),进行如下滤波计算yf(k)=Df(z-1)y(k) (3)uf(k)=Df(z-1)u(k) (4)Step 5. 由(2)(4)式,列成如下自回归方程:,其中,一、成批型GLS法(5/9)计算步骤,并由此得如下向量式的自回归方程组Yf=ff+W (6)其中Yf=yf(1), yf(2), ., yf(L)W=w(1), w(2), ., w(L)f=f(0), f(1), ., f(L-1)Step 6. 用如下LS法计算f参数估计值,一、成批型GLS法(6/9)计算步骤,St
14、ep 7. 由(2)式有y(k)=(k-1)f+v(k) (8)其中(k-1)=y(k-1), , y(k-na) u(k-1), , u(k-nb),因此,基于(7)式的参数估计值 f,计算相关扰动v(k)的估计值,即回归式(8)的如下估计残差,一、成批型GLS法(7/9)计算步骤,Step 8. 由(1)式,可得如下关于相关扰动v(k)和白噪声w(k)的自回归方程,其中,并由此得如下向量式的自回归方程组V=vv+W (11)其中V=v(1), v(2), ., v(L)v=v(0), v(1), ., v(L-1),一、成批型GLS法(8/9)计算步骤,Step 9. 选取上一步估计得的D
15、(z-1)的估计值作进行系统白色化处理的线性滤波器Df(z-1).若模型的估计残差满足给定的精度或者Df(z-1)、A(z-1)和B(z-1)的变化量相对较小则循环结束,否则返回第4步继续循环.,因此,将用式(9)计算好的相关扰动v(k)的估计值 (k)来替代自回归方程组(11)的相关扰动v(k),则可得如下回归参数向量v,即噪声模型D(z-1)的系数的LS估计值,一、成批型GLS法(9/9),以上辨识过程表明,GLS法的思想是对输入输出数据先进行一次滤波预处理,然后利用普通LS法对滤波后的数据进行辨识,并反复迭代此过程.可想而知,这种方法受滤波模型的好坏的影响较大,而在迭代过程中滤波模型的好
16、坏也直接与系统模型辨识结果有直接的关系.从优化理论的角度来说,GLS法其实属于非线性优化方法.因此,难以避免出现非线性优化中的局部极值点情况,即该方法并不能保证得到的估计值是一致无偏的.这是GLS法的一个不太令人满意之处.,二、递推型GLS法(1/6),二、递推GLS法递推GLS法的基本思想是与成批型LS法大致相同,所不同的是:由于所考虑的是递推估计,不能像成批型那样作反复迭代.解决的方法是分别对(5)和(10)所示的两个模型的辨识设计两个递推估计算法,并在每一个递推步中,让它们依顺序递推一次.随着递推过程的深入,将不断改进噪声模型D(z-1)的辨识结果,同时亦得到较佳的A(z-1)和B(z-
17、1)的辨识结果.,二、递推型GLS法(2/6),根据上述基本思想,有如下递推GLS估计算法:Step 1. 确定被辨识系统模型的结构,以及多项式A(z-1)、B(z-1)和D(z-1)的阶次;,Step 2. 选定递推参数初值 f(0)和Pf(-1), v(0)和Pv(-1),以及稳定的初始滤波器Df(z-1);,Step 3. 采样获取新的观测数据y(k)和u(k);Step 4. 基于滤波器Df(z-1),分别由下式计算y(k),u(k)的滤波值yf(k)和uf(k);yf(k)=Df(z-1)y(k) uf(k)=Df(z-1)u(k),二、递推型GLS法(3/6),Step 5. 构造
18、观测数据向量f(k-1),并依如下RLS法估计多项式A(z-1)和B(z-1)的系数,即回归参数向量f,二、递推型GLS法(4/6),Step 6. 基于式,计算相关扰动v(k)的估计值,并构造观测数据向量v(k-1).,二、递推型GLS法(5/6),Step 7. 依如下RLS法估计多项式D(z-1)的系数,即回归参数向量v,Step 8. 选取上一步估计得的D(z-1)的估计值作进行系统白色化处理的线性滤波器Df(z-1).将循环次数加1,返回第3步继续循环递推辨识.,三、GLS仿真(1/14),三、 广义最小二乘法仿真考虑如下有色噪声扰动的随机线性离散系统,的GLS在线辨识仿真程序伪代码
19、.,三、GLS仿真(2/14),三、GLS仿真(3/14),三、GLS仿真(4/14),三、GLS仿真(5/14),三、GLS仿真(6/14),三、GLS仿真(7/14),例1 本例中的仿真对象、实验条件与第四讲中例1相同,模型结构选择如下y(k)+a1y(k-1)+a2y(k-2)=b1u(k-1)+b2u(k-2)+v(k) (19)v(k)+d1v(k-1)+d2v(k-2)=w(k) (20)仿真结果如表1所示.,三、GLS仿真(8/14),表1 计算机仿真结果,三、GLS仿真(9/14),递推辨识过程的辨识值如下图所示,噪信比=73%,初始白化滤波器Df=1时递推辨识结果(1),噪声
20、估计误差,三、GLS仿真(10/14),噪信比=73%,初始白化滤波器Df=1时递推辨识结果(2),噪声估计误差,三、GLS仿真(11/14),噪信比=73%,初始白化滤波器Df=1+0.5z-1-0.5z-2时递推辨识结果(1),噪声估计误差,三、GLS仿真(12/14),噪信比=73%,初始白化滤波器Df=1+0.5z-1-0.5z-2时递推辨识结果(2),噪声估计误差,三、GLS仿真(13/14),噪信比=23%,初始白化滤波器Df=1时递推辨识结果(1),噪声估计误差,三、GLS仿真(14/14),噪信比=23%,初始白化滤波器Df=1时递推辨识结果(2),噪声估计误差,三、GLS评述
21、(1/1),三、GLS评述对GLS算法,有如下评述:有色干扰下估计精度较高迭代收敛较快,但是收敛性未得到证明能同时得到过程参数和噪声参数的估计计算量大,费机时,3. IV法(1/2),3. 辅助变量法上两讲分别讨论了同时辨识系统模型和噪声模型的ELS法和GLS法.在有些建模及其模型的运用问题中,并不需要知道噪声模型,即不需要对噪声建模(辨识).此时,若采用ELS法和GLS法来辨识,需花费较多的计算时间,而且辨识的参数越多则辨识的精度和效果越差.这一点,GLS法尤其突出.,3. IV法(2/2),本讲讨论有相关随机扰动时的随机离散系统参数估计的辅助变量(Instrument Variable,
22、IV)法.该方法不需要辨识系统的噪声模型,只要辅助系统选择得恰当,便可获得较高精度的无偏一致估计. 本讲介绍的主要内容为:成批型算法辅助变量的选择递推型算法IV法的评价,一、成批型算法(1/6),一、成批型算法考虑如下自回归方程的一致辨识问题y(k)=(k-1)+w(k), k=1,2,.,L (1)或YL=L+WL (2)其中对扰动量w(k)未加白噪声的假定,其它符号与第三讲一致.下面先考察一般LS法的参数辨识值的一致收敛性.,一、成批型算法(2/6),由第三讲的提出的一般LS估计式,有,一、成批型算法(3/6),由式(3)可知,欲获得无偏一致的参数估计值,若,依据Frechek定理: 若随
23、机序列xk,则f(xk)f() 则由条件(4)与(5),有LS 0即LS估计为无偏估计.,一、成批型算法(4/6),但,当w(k)不为白噪声时,条件,并不能保证成立,相应的LS估计值,亦不能保证为无偏一致的.因此下面的问题是:如何在有色噪声干扰条件下得到无偏估计？,一、成批型算法(4/6),为了对w(k)不为白噪声时,可获得q的无偏一致估计值,定义一个如下辅助观测数据矩阵*L=*(0), *(1), , *(L-1) (6)其中辅助变量矩阵须满足,一、成批型算法(5/6),因此,使得对不定(不相容)方程组YL=L (9)和回归方程YL=L+WL可获得如下IV解,一、成批型算法(6/6),当辅助
24、变量矩阵满足条件(7)和(8)时,则有,可见,只要选择适当的辅助变量,使之满足上述的两个条件,IV估计值就可以是无偏一致的.IV估计法的关键在于如何具体构造出辅助变量矩阵*L.,依据Frechek定理与条件(7)和(8),二、辅助变量的选择(1/9),二、辅助变量的选择IV法的关键是选择适当的辅助变量使得条件(7)和(8)得到满足.如果采用下图所示的输出x(k)作为辅助变量,置*(k-1)=x(k-1), ., x(k-na), u(k-1), ., u(k-nb) (12),那么,当u(k)为高阶持续激励信号时,条件(7)必定满足.又因x(k)只与u(k)有关,也就是说与噪声无关,故条件(8
25、)也满足.,二、辅助变量的选择(2/9),人们已经得出几种不同的方案,要严格证明它们满足上述条件是很难的,但是在实际应用中已经取得很好的结果.常用的辅助变量选择方法有:自适应滤波纯滞后Tally原理,二、辅助变量的选择(3/9)-自适应滤波,A. 自适应滤波当把上图所示的辅助系统模型看成是自适应滤波器时,辅助变量可按如下关系生成或求得x(k)=*(k)(k) (13)或x(k)=*(k)(k) (14a)(k)=(1-)(k-1)+(k-d), 取0.010.1; d取010 (14b)其中(k)表示k时刻的IV参数估计值,它可由递推算法给出.,二、辅助变量的选择(4/9)-自适应滤波,下面简
26、单分析一下上述自适应滤波法所选择的辅助变量与噪声变量的相关性。对*(k)与噪声w(k)的相关性分析如下,因此,经过层层传递(每次传递相关性小于1),*(k)与噪声w(k)的相关性将变得非常小。当噪声w(k)与输入u(k)无关时,且u(k)为高阶持续激励信号时,可以证明选用式(13)和(14)作为辅助变量以较大概率满足上述两个条件Soderstrom,1974.,二、辅助变量的选择(5/9)-纯滞后,B. 纯滞后当图1所示的辅助系统模型选为纯滞后环节时,辅助变量取作x(k)=u(k-d) (15)其中dnb,nb为多项式B(z-1)的阶次.这时辅助变量向量*(k-1)可写成*(k-1)=u(k-
27、d-1), ., u(k-d-na), u(k-1), ., u(k-nb),二、辅助变量的选择(6/9)-纯滞后,值得指出的是,只有u(k)为高阶持续激励信号时,上述选取的*(k-1)的列之间现行无关,收敛条件(7)才能满足。此外,当输入u(k)与噪声w(k)无关,*(k)与噪声w(k)统计无关。因此u(k)为高阶持续激励信号时,上述选取的*(k)满足上述辅助变量可满足上述两个条件.,二、辅助变量的选择(7/9)-Tally原理,C. Tally原理如果相关噪声的模型可以看成为,时,其中dnc,则辅助变量可取作x(k)=y(k-d) (16)即*(k-1)=y(k-d-1), ., y(k-
28、d-na), u(k-1), ., u(k-nb),二、辅助变量的选择(8/9)-Tally原理,下面简单分析一下Tally原理所选择的辅助变量与噪声变量的相关性。对*(k)与噪声w(k)的相关性分析如下,因此,经过层层传递(每次传递相关性小于1),*(k)与噪声w(k)的相关性将变得非常小。当噪声w(k)与输入u(k)无关时,且u(k)为高阶持续激励信号时,可以证明Tally原理所选择的辅助变量以较大概率满足上述两个条件方崇智等,1988.,二、辅助变量的选择(9/9),在上述3种选择辅助变量方法中，自适应滤波法辨识收敛快，精度高，递推辨识结果平滑，波动小。Tally原理法次之。但纯滞后法辨
29、识收敛慢，精度低，数据增多对精度提高效果不显著；再之，递推辨识结果波动大，难于达到控制的要求。,三、递推型算法(1/15),三、递推型算法下面讨论IV参数估计算法(10)的递推算法.类似于RLS法的推导,由成批型IV法可得如下递推IV算法,其中P(k-1)和K(k-1)的定义式为,K(k-1)=P(k-1)*(k-1) (21),三、递推型算法(2/15),综上所述,递推IV法的基本计算步骤可总结如下确定被辨识系统模型的结构,以及多项式A(z-1)和B(z-1)的阶次;确定或设计所采用的辅助变量系统;设定递推参数初值(0),P(-1);采样获取新的观测数据y(k)和u(k),并组成观测数据向量
30、(k-1);计算辅助变量x(k),并组成辅助变量观测数据向量*(k-1);用式(17)(19)所示的递推IV法计算当前参数递推估计值;循环次数k加1,然后转回到第4步继续循环.,三、递推型算法(3/15),下面给出针对随机线性离散系统(XARMA),的ELS在线辨识仿真程序伪代码.,三、递推型算法(4/15),三、递推型算法(5/15),三、递推型算法(6/15),三、递推型算法(7/15),例1 本例的仿真对象及实验条件与第五讲相同.若利用一般LS法是得不到无偏一致估计的.这里采用IV法,选择如下模型结构,仿真结果如表1所示.,三、递推型算法(8/15),表1 计算机仿真结果(噪信比=73%
31、,C(z-1)=1-1.0z-1+0.2z-2),三、递推型算法(9/15),递推辨识过程的辨识值如下图所示,自适应滤波法辅助变量递推辨识结果,三、递推型算法(10/15),d=3时的纯滞后法辅助变量递推辨识结果,三、递推型算法(11/15),d=4时的纯滞后法辅助变量递推辨识结果,三、递推型算法(12/15),d=5时的纯滞后法辅助变量递推辨识结果,三、递推型算法(13/15),d=3时的Tally原理辅助变量递推辨识结果,三、递推型算法(14/15),d=4时的Tally原理辅助变量递推辨识结果,三、递推型算法(15/15),d=5时的Tally原理辅助变量递推辨识结果,四、IV 法评价(1/1),四、IV 法评价对IV法与LS法,有如下评述: IV法计算量比LS估计增加不多,而在相关噪声情况下估计精度有明显改善.IV法与递推IV法思路相似,但不等价.递推IV法计算量与RLS法接近,在相关噪声情况下估计精度优于RLS法.为了提高递推IV法的可靠性,在递推的前几十步最好用RLS法过渡.,
