系统辨识第四章课件.ppt

1、系统辨识第四章2第一章第一章 模型方法与辨识模型方法与辨识第二章第二章 脉冲响应辨识脉冲响应辨识 第三章第三章 最小二乘辨识最小二乘辨识第四章第四章 极大似然辨识极大似然辨识第五章第五章 时间序列建模与随机时间序列建模与随机逼近辨识逼近辨识第六章第六章 模型阶次的辨识模型阶次的辨识第七章第七章 闭环系统辨识闭环系统辨识3第四章第四章 极大似然辨识极大似然辨识 前言前言 4-1 极大似然原理极大似然原理 4-2 动态系统模型参数的极大似然估计动态系统模型参数的极大似然估计 4-3 极大似然估计的一致性极大似然估计的一致性 4-4 预报误差参数辨识法预报误差参数辨识法4第四章第四章 极大似然辨识极

2、大似然辨识 极大似然法，是一种适用范围非常广泛的传统辨识方法，1906年，由R.A.Fisher提出。极大似然估计方法在随机系统参数估计、故障检测及容错控制等方面，有广泛应用。把这种经典的估计方法用于动态过程或动态系统辨识，可以获得良好的估计性质。极大似然法要求已知输出量的条件概率密度函数，建立随机观测数据与未知参数之间的概率特性和统计关系，通过使条件概率密度函数为极大的准则，求出未知参数的估计值。因而，极大似然辨识法是一种概率性的参数估计方法。54-1 极大似然原理极大似然原理一、似然函数一、似然函数。密度函数为的条件概率条件下，数是一个随机变量，在参设)(zfzz总体：次观测，可得观测数据

3、进行对被估量N TNNzzzz)(),2(),1(件概率密度函数为：为随机序列，其联合条则)(kzzN)(Nzf6以示区别。记作的似然函数，有时叫做称这时的概率密度函数，的函数，不再表示条件参数只是为一组确定的数据，则如果)()()(NNNNzLzfFisherzfz )(Nzf件概率密度函数为：为随机序列，其联合条则)(kzzNTNNzzzz)(),2(),1(7，如果及分别取为时，函数。当的可以认为只是而言，对于具体观测值21)(zzzfzN)()(21zfzf 可见，条件概率密度函数与似然函数有不同的物理含义，但其数学表达形式一致，即)()(NNzfzL8。确定估值然函数达到极大来称为似

4、然函数，并以似或数把验前条件概率密度函大。故为准确值的可能性要为准确值的可能性比则)()(21zLzfFisher，如果及分别取为时，函数。当的可以认为只是而言，对于具体观测值21)(zzzfzN)()(21zfzf9下图所示。为标量，其似然函数如测值两种可能，观和只有例如，设未知参数z21 )(zf)(1zf)(2zf1z2z0z10为极大似然估计故：若1211),()(zfzfzzN为极大似然估计故：若2212),()(zfzfzzN11二、极大似然估计求法二、极大似然估计求法 极大似然估计定义的极大值，即是值，似然函数如果对于所有可能的)()(NNzfzf)(max)(NNzfzf。大似

5、然估计，记作的极是的可能性最大，称是准确值则ML 0012似然方程与对数似然方程似然方程与对数似然方程)(max)(NNzfzfML故可通过似然方程0022MLMLNNzfzf)(,)(。求出极大似然估计ML由于。由于数似然方程确定为了简化运算，常用对ML)()(1)(lnNNNzfzfzf13)()(1)(lnNNNzfzfzf等效。与的值不可能为零，所以的极大值附近，而在0)(0)(ln)()(NNNNzfzfzfzf从而可由值达到极大在相同的与因此由于对数为单调增函数指出与,)()(ln,:.NNzfzfDoobLJWaldA对数似然方程0022MLMLNNzfzf)(ln,)(ln。确

6、定极大似然估计ML14正态独立同分布随机过程均值与方差的极大似然估计正态独立同分布随机过程均值与方差的极大似然估计则其概率分布密度为正态分布随机变量，设iz22)(21exp21),(mzmzfii为均方差。为方差，为均值，其中，2m条件下的似然函数为和在参数布特性，则随机变量机序列，具有独立同分下获得的随是一组在独立观测条件设2)(,),2(),1(mzkzzzkikmizfmkzfmzfmzfmzf122222),)(),)(),)2(),)1(),(15正态分布，故似然函数因为zkiikmzmzf1222)(21exp21),(取对数似然函数ln21ln)(21)(2121ln)(21e

7、xp21ln),(ln1221221222kkmzmzmzmzfkiikiikiikkikmizfmkzfmzfmzfmzf122222),)(),)(),)2(),)1(),(16ln21ln)(21),(ln1222kkmzmzfkiik用对数似然方程：令kikiikmkizmzmmzf121220)(1),(ln有kiizkm1)(117令0ln)(42)(41),(ln222142kmzmzmzfiikik 0因VVV)(ln故得：0)(1213kmzkiiln21ln)(21),(ln1222kkmzmzfkiik180)(1213kmzkiikmzkii212)(1求出212)(1

8、mzkkiikiizkm1)(119验证：4122222222230kiikkmzkmzfkkmmzf)(),(ln)(),(ln为观测次数，恒正ln21ln)(21),(ln1222kkmzmzfkiik023),(ln)(2422222221kkkmzfkmzkkii确是极大似然估计。和表明求出的2m20例4-1的极大似然估计。求参数，试，的概率密度为条件下，随机变量过程，在参数为一个独立同分布随机设0)()()(2 xxexftxkx解个观测向量为的设随机变量Ntx)(TNNxxxx)()2()1(条件下的似然函数为在参数则随机变量)(txNkNkNNkNkNkxkxkxkxkxfxf1

9、12121)(exp)()(exp)()()(21NkNkNNkNkNkxkxkxkxkxfxf112121)(exp)()(exp)()()(相应的对数似然函数0)(2 )(ln2)(ln)()(ln2)(ln1111NkNkNNkNkNkxNkxNxfkxkxNxf22NkMLkxN1)(2且02)(ln222MLNNxf的极大似然估计。是参数极大值，因此所求使对数似然函数达到了表明MLML 02221111NkNkNNkNkNkxNkxNxfkxkxNxf)()(ln)(ln)()(ln)(ln 23例4-2，其概率密度为独立同分布随机序列设)(kx0 x,0 0 x,exp4)(223

10、2xxxf的极大似然估计。，为待估参数，求其中，0解的似然函数条件下，随机变量在参数xNkNkNNNkNkNkxkxkxkxkxfxf112223122321)(1exp)()4()(exp)(4)()(24TNNxxxx)(,),(),(21 其中相应的对数似然函数NkNkNkxkxNNxf11222)(1)(ln3)4ln()(ln，求得令0)(23(123NkNkxNxfNkMLkxN12)(32NkNkNNNkxkxxf112223)(1exp)()4()(254-2 动态系统模型参数的极大似然估计动态系统模型参数的极大似然估计一、第一、第1种模型噪声情况种模型噪声情况设动态系统差分方

11、程为niniiikvikubikyaky1100)()()()(，得系统量测方程：，令Nnnk1)()()(0NVNNY26)()()(0NVNNY式中)()1()()(,)()1()()(,)()1()()(NnnNNvnvnvNVNynynyNYTTTTnnbbbaaa00201002010)()2()1()()2()1()(nkukukunkykykykT27NnnjkjkjkjvkvENnnkkvE,)()(,)(10102 时当，时，当布、零均值，则有同分的统计特性为正态独立设噪声序列)(kv28噪声的联合概率密度函数噪声的联合概率密度函数率密度函数为正态分布，故其联合概因：维向量，

12、其统计分布为噪声)(),()()(NVINNVnNNV201 11220112121nNinNiiYvfVf)(exp)()(20021221exp)2()(YYNVfTnN29向量方程误差的似然函数向量方程误差的似然函数，其中代替如果用0nnbbbaaa2121则向量方程误差为（残差）YN)(联合概率密度函数为一样，属正态分布，其与，故，由于)()()()(0NVNYNYNV 2212212YYNfTnNexp)()(30 2212212YYNfTnNexp)()(的似然函数：即是条件概率密度函数，合概率密度函数，它应条件下得到的联和是在实际上，上述)()(NNf2 22122212YYfT

13、nNexp)(),(31输出观测向量的似然函数输出观测向量的似然函数根据随机向量变换法则，可以导出221222212)()(exp)(),(),(YYfYfTnN32模型参数的极大似然估计模型参数的极大似然估计观测向量Y的对数似然函数为222)()(21)ln2(ln21),(lnYYnNYfT，有）已知，令当量测矩阵022MLMLYf,(ln0212122)()(TTTTTTMLTYYYYYY 33得：0)22(212MLTTMLY整理：)(0)(12aYTMLTML，有再令0),(ln22Yf)(02)()(2142bYYnNMLMLTMLML0 21)()(2122TTTTTTMLTYY

14、YYYY34定义为：因向量方程误差),(N),(),1(),(),(NeneneN),(),(),()()(2MLNnkMLMLTMLTMLkeNNYY而大似然估计：式，分别求得参数的极式及故由)()(ba),(111)()()(221MLNnkMLTMLMLLSTTMLkenNnNYYy35两点注意事项两点注意事项的一致估计。是相同，与最小二乘估计的极大似然估计存在条件下，模型参数阵已知、分布噪声序列、为零均值、正态独立同只有当001MLLSMLTkv)()(。于推导递推极大似然法的估计。这一结果将用取极小时的数相当于下列性能指标函然估计的极大似似然函数取极大所得已知时可见，当噪声方差由对数

15、似然函数02 ,NnKTTTkeYYJYYnNYf),()()()()()()ln(ln),(ln202222122136二、第二、第2种模型噪声情况种模型噪声情况设动态系统的差分方程为：niiniiniiikvckvikubikyaky11010)()()()()(。零均值正态独立同分布噪声序列)(kv令组合噪声其余为简便研究，令,0,01iccniiniiikubikyakykvckvkz101011)()()()()()(3701)()()(,NNYNzNnnk得组合噪声动态方程：令种噪声情况，而表达式同第及、式中，1)()(0NNY)()()()()()()(111111NvcNvnv

16、cnvnvcnvNz因z(k)是v(k)与v(k-1)的线性组合，故z(k)也是零均值正态分布噪声序列，但不再是无关序列。38的相关函数：考察)(kzNnnjkjkcjkcjkcjvkvcjvkvcjvkvcjvkvEjvcjvkvckvEjzkzE,)()()()()()()()()()()()()()()(101111111112121221211111其他种噪声模型不同。为一步相关序列，与第表明1)(kz噪声的统计特性均值NnnkkzE,)(1039方差阵RcccccccNzNzET22111211121100001001)()(联合概率密度函数)()(exp)()(01012121YR

17、YRNzfTnN的行列式。表示式中，RR40观测向量的似然函数观测向量的似然函数数：似然函数及对数似然函的关于可以导出代替种情况，用类似于第Y,1 0)()(21exp)2(1)(11YRYRYfTnN)()(21ln21)2ln(21)(ln1YRYRnNYfT模型参数的极大似然估计模型参数的极大似然估计已知，令及噪声方差阵如果量测矩阵R410)(lnMLYf有0222121111111)(MLTTTTTTTTRYRRYRRYYRYML 存在，得：若11)(RTLSRTTMLYRR1)(111表明：此时模型参数的极大似然估计正好等于模型参数的马尔可夫估计。424-3 极大似然估计的一致性极大

18、似然估计的一致性 动态系统极大似然估计的一致性问题，经常转换为预报误差方程的预报误差估计的一致性问题。因为一大类含噪声的线性系统的差分方程或状态方程，可以转化为预报误差方程，而预报误差方程与似然函数之间可以建立起直接的联系，所得结果具有较宽的适用范围。根据输入和输出数据，可得相应集合),(),()(),(),()(1211kukukUkykykY43构造预报误差方程：)(),(),()(kvkUkYfky01动态系统：为噪声。对于线性为待辨识参数，其中，)(kv0)()()()(kvikubikyakyniinii1010因此，预报误差方程中的)()()()()()()()(),(),(nku

19、bkubkubnkyakyakyaikubikyakUkYfnnniinii0020100201101002121144Tnnbbbaaa00201002010可见，预报误差方程代表一大类含噪声的线性动态系统。预报误差模型：为如下时，预报误差方程演变代替真值当用估值0 ),(),(),()(kekUkYfky1为预报误差。其中：)e(k,假设：统计独立，即和输入序列序列，且与输出序列的新息、方差为为独立同分布、零均值噪声序列)()(v(k)2kUkYo114501)(),()(kUkYkvE噪声的条件期望序列。平稳随机为具有各态历经特性的输出序列)(kyo2取性能指标为NnkNkenNDJ),

20、()()(211(残差平方和的均值)，有种噪声差分方程模型中上节已证：在第的样本协方差。为预报误差序列其中1),()(keD46),()(MLNnkMLLSTTMLkenNY22111。，也是最小二乘估计件下的极大似然估计，就是正态条的预报误差估计因此，LSMLNNNJ)(min)(。从而证得，表为）收敛到其期望值（真值时，依概率在致性，即具有一条件下，预报误差估计不难证明：在上述假设NsaMLNsaNNN)(:1)(47下面给出简要证明。由预报误差模型，预报误差：),(),()(),(kUkYfkyke1由预报误差方程知：)(),(),()(kvkUkYfky1因而预报误差)(),()(),

21、(),(),(),(),(kvfkvkUkYfkUkYfke1148式中),(),(),(),(),(kUkYfkUkYff11表达式代入性能指标：将预报误差),(keNnkNnkNnkNnkNnkNnkNfnNfkvnNkvnNffkvkvnNkvfnNkenNDJ),(),()()(),(),()()()(),(),()()(2222221121111211111149收敛到其期望值，即各项均以概率中统计独立，故与成立，且与因假设条件1)(),()(NooJfkv2101111222),()(),()(),()()()(fEkvEfkvEfkvnNkvEkvnNNnKNsaNnKNsa(噪

22、声方差)0),(),(2211fEfnNNsaNnK方差指 f50收敛到如下期望值：依概率时，性能指标因而，当1)(NJN),()(22fENJNsaN的估计，必有为性能指标由于min)()(NJNN0),(2fE，则而)(N011),(),(),(),(),(kUkYfkUkYff51于是有 NsaN)(从而证明了预报误差估计和正态条件下的极大似然估计具有一致性，都是真实参数 的一致估计。524-4 预报误差参数辨识法预报误差参数辨识法 极大似然法要求已知数据的概率分布，通常都假设数据服从正态(高斯)分布。然而，实际问题中的数据不一定都是正态分布的。当数据的概率分布不知道时，无法应用极大似然

23、估计。预报误差参数辨识法不要求数据概率分布先验知识，是一种更加一般的参数辨识方法，也是极大似然估计的一种推广。业已证明，当数据的概率服从正态分布时，预报误差估计法等价于极大似然法（Goodwin澳大利亚教授,1977）.53一、预报误差准则预报误差模型，构造向量，模型参数，输入向量设输出向量ormRkuRky11)()()()()()(,)()()()(12111211ukukukUykykykY54则预报误差模型：)(),(),()(kekUkYfky11。阵为，协方差为预报误差向量，均值其中)()()()(jekeEkeERkeTem01 预报误差模型表明：k时刻的输出，可以用k时刻以前的

24、数据来“预报”。的条件数学期望值，即的“最好”预报，可取的条件下，对输出和在获得数据集合)()()()(kykykUkY1155),(),()()(11kUkYkyEky这种预报，可使预报误差范数平方的条件期望最小，即min),(),()()(112kUkYkykyE 显然，这种“最好”的输出预报，应是“最好”模型的输出。由最优控制理论知，“最好”输出预报，应是使某一个预报误差准则(即性能指标)为极小而获得。56预报误差准则常用的误差准则有如下两种：)(detlog)()()()2()1(DJWDtrJNN单出系统，有对于单入。，为正定矩阵。若取式中，-)()()(trDJIWWN1NnKNk

25、enNJ),()()(211157而在多入-多出情况下),(),()()()()()(1111kUkYfkykekekenNDNnkT。差阵的协方将收敛于预报误差，上节已证：当ekeDN)()(知识。着数据概率分布的先验用不。显然，确定计称为预报误差估计极小化，得到的参数估或通过使)()()()()2()1(NNJJNN58二、预报误差法与极大似然法之间的关系预报误差模型的似然函数统计独立，则对数据与分布，且统计独立，但不一定同若新息序列)()()(kykuke TTTTnyNyNyNY)()1()()(应用Bayes公式，得条件概率密度函数NnkkUkYkyfnUnYnyfNUNYNyfNU

26、NYNyfNUNYf),(),()(),(),()(),(),()(),(),()(),()(111122111159由预报误差模型)(),(),()(kekUkYfky11，对的条件概率分布。因而决于预报误差的条件概率密度分布取条件下，和知，在)()()(),(kekykUkY11),(),()(),()(111kUkYkyfNUNYfNnk得进行随机变量置换，可中的),(),()(11kUkYkyf)()(det),(),()(),(),()(kykekUkYkefkUkYkyf111160绝对值。因为是雅可比矩阵行列式的其中，)()(detkyke),(),()()(11kUkYfkyk

27、e必有1)()(detkyke故),(),()(),()(111kUkYkefNUNYfNnk61可进一步表示为：的似然函数，则上述预报误差模型为零，协方差阵为独立同分布，且其均值维新息序列如果进一步设),1()()(NUNYfkeme )()(expdet)()()(expdet),()(kekekekeNUNYfeNnkTnNemeTNnkem1211212122121幅值相乘相角相加62预报误差协方差 已知时的似然函数e由似然函数)()(expdet)(),()(kekeNUNYfeNnkTnNem1212121。等价于使，已知时，使似然函数易见，当min)()(max),()(keke

28、NUNYfNnkeTe1163根据矩阵迹的运算性质：.TTAXXtrAXX故nnkTeeNnkTNkeketrkekeJ)()(.)()()(11已令NnKTkekenND)()()(11有)()1()(1DnNtrJeN(样本协方差)64，则有若取正定矩阵1)1(enNW)()()()1(WDtrJJNN表明：法是等价的。极大似然法与预报误差意义下种的预报误差估计。在这准则时，使预报误差就是当的极大似然估计，已知时，当协方差阵,min)()()(111NeeJnNW65预报误差协方差 未知时的似然函数e取负对数似然函数，有NnkeTeekekenNNUNYfJ)()(detlog)log()

29、,()(log),(121212121)n-m(N)()(keketrTe1根据矩阵迹的微分运算法则：66TAAA)(detlog1为对称矩阵，则有注意到eNnkeTeeeekekenNJ1112121)()(),(TBAXXBAXtrX)(.111的极大似然估计：，得预报误差协方差阵令0),(eeJ67e以 代替 ，负对数似然函数为eNnkeTeekekenNnNmJ)()(detlog)log()(),(12121221，于是，因)(detlog)()()2(DJDNeNnkTNkeDkeJnNnNmDJ)()()()()log()()(,()(12212122168：当 为正态分布不相关

30、随机向量时，等价于 取极小，又等价于第2种预报误差准则 取极小。在这种意义下，极大似然法与预报误差法也是等价的。)(,(DJ)()2(NJ结论)(ke预报误差法与极大似然法等价。或者说，极大似然法是预报误差法的特例。在上式中,右端各项为正,必有：等价于min)()2(NJmin)(,(DJ),()(1NUNYf .由于似然函数 取极大,69三、预报误差参数估计方法(Newton-Raphson法)预报误差参数估计法实质 由于预报误差准则 或 一般都是参数)()1(NJ)()2(NJ的非线性函数，故令 极小化求 的方法，)2,1()()(lJlN归纳为极小化 的最优化算法。)()(lNJe若预报

31、误差 的协方差阵 已知，则取)(ke)()1(NJ作为预报误差准则，且取权阵 ；1)1(enNW若 未知，则应选 为预报误差准则。e)()2(NJ70预报误差准则极小化的最优化算法根据NewtonRaphson原理，)2,1(min)()(lJlN的最优化算法归纳为如下迭代方程：nTllnnnJJ)(12)(21)()(111式中：1n 第 次迭代的参数估计值；1n)()(lJ 预报误差准则 关于 的梯度；)()(lJ712)(2)(lJ Hessian矩阵；1n 迭代步长，使min)(1)(1nlJ显然，上述最优化算法的关键是：关于 的梯度)()(lJ及Hessian矩阵的具体计算式；利用一

32、维搜索法求 ，使1nmin)(1)(1nlJ。72梯度与Hessian矩阵的计算设 的协方差矩阵已知，即 已知，取权阵)(kee1)1(enNW，而NnkTDkekenN)()()(11(n为系统阶次)因而：)(1 DW由矩阵运算公式：AXAXXXAXXAXXtrTTTT2;73其中TAA 可得梯度向量 的第i个元素为：)()1(JiNnkTNnkTiNnkTiiikeWkenNkWekenNkekWetrnNWDtrJ)()(12)()(11)()(11)()()1(代入 表达式)(DAXXAXXtrTT)(AXAXXXTT274Hessian矩阵 的第i行，第j列的元素为：2)1(2)(J

33、)()()()()()()()()()()()(ieTNnkjijTiTNnkjNnkiTjijijkekekeWkekeWkenNkeWkenNJJ1211221212代入 表达式iJ)()1(略去75需要指出：计算Hessian矩阵元素时，忽略了e(k)关于 的二阶导数，目的是简化Hessian矩阵计算，并可保证Hessian矩阵正定性。式中 为 的维数，即待辨识的参数nnNji2;1,2,1,),(),()(11kUkYfkeii个数，且 (预报误差模型)()()()()()()()()()()()(ieTNnkjijTiTNnkjNnkiTjijijkekekeWkekeWkenNke

34、WkenNJJ121122121276Newton-Raphson法迭代计算步骤o1设置初始状态 和步长 ，赋 ；1n1n0:1no2计算 和 ，)(,)(11nDken1nke)(Nnnk,1),(),()()(111nkUkYfkykeNnkTkekenND)()()(11)(),(),()(121111nN,ikUkYfkenii 77o3计算梯度 和Hessian矩阵1)1()(nJ12)1(2)(nJiNnkeTikekenNJ)()()()(1112);,()()()()(1212112nNjikekeJieTNnkjij o4利用一维搜索法，求步长 ，使1nmin;)()(111)1(12)1(2)1(nTnnJJJ78o5计算11n;)()(1111)1(12)1(21nTnnnJJo6给定迭代精度正小数 ，若1111nnn则停止迭代；否则，赋 ，返回 。1:11 nno2