极大似然参数辨识方法

《极大似然参数辨识方法》由会员分享，可在线阅读，更多相关《极大似然参数辨识方法（8页珍藏版）》请在金锄头文库上搜索。

1、2.1.3)(2.2.1)(2.2.2)(2.2.3)(2.2.4)2 极大似然参数辨识方法极大似然参数估计方法是以观测值的出现概率为最大作为准则的，这是一种很普遍的参 数估计方法，在系统辨识中有着广泛的应用。2.1 极大似然原理设有离散随机过程V与未知参数有关，假定已知概率分布密度f(V p)。如果我们 得到n个独立的观测值V ,V,V，贝y可得分布密度f (Vp )，f (V p)，f (V p)。12,n12n要求根据这些观测值来估计未知参数p，估计的准则是观测值V 的出现概率为最大。k 为此，定义一个似然函数L(V ,V，,V )二 f (V ) f (V p ”f(V )(2 1 1

2、)12n12n(211丿上式的右边是n个概率密度函数的连乘，似然函数L是的函数。如果L达到极大值，V k 的出现概率为最大。因此，极大似然法的实质就是求出使L达到极大值的的估值。为了便于求，对式(2.1.1)等号两边取对数，贝把连乘变成连加，即InL = Inf (V )(2.1.2)ii=1 由于对数函数是单调递增函数，当L取极大值时，lnL也同时取极大值。求式(2.1.2) 对 的偏导数，令偏导数为0，可得d In L =0 d解上式可得的极大似然估计ML。2.2 系统参数的极大似然估计 设系统的差分方程为a(z-1)y(k) = b(z-Qu(k) + g (k) 式中a(z-1) =

3、1 + a z-1 +. + a z-n1nb(z-1)=b +bz-1+.+b z-n01n因为乍(k)|是相关随机向量，故(2.2.1)可写成a(z-1)y(k) = b(z-1)u(k) + c(z-1)e (k) 式中c (z-1)8 (k) = E (k)c(z-1) = 1 + c z-1 + c z-n1n8 (k)是均值为0的高斯分布白噪声序列。多项式 a(z-1)，b(z-1)和c(z-1)中的系数 a , a, b，b , c,c和序列8 (k)的均方差b都是未知参数。1,0n 1n设待估参数0= a a b b c c 】T(2.2.5)1 n 0 n 1 n并设 y(k

4、) 的预测值为y (k) = - a y (k -1)a y (k - n) + b u (k) HF b u (k - n) +1n0nc e(k -1) HF c e(k - n)12.2.6)式中e(k -i)为预测误差；a.，b.，c.为a，b， iiiiic 的估值。预测误差可表示为 i(k) = y (k) y(k) = y (k) 工 a y (k i) + 工 b u(k i) +i=1i=0或者工 c e(ki) = (1+a z-1 + ai1n z -n) y (k) (b 0+哲 z -1 + + b ” z -n )u (k)i=1(c z-1 + c z-2 + c

5、 z-n )e(k)12n2.2.7)(1 + c z -i + c” z -n )e(k) = (1 + z -i + a” z-n) y (k)(b + b z-i hf b z-n )u (k)01n因此预测误差b(k)满足关系式2.2.8)式中c (z -i)e(k) = a (z -i) y (k) 一 b (z -1)u (k)2.2.9)ccca(z-1) =1+ a z-1 + a z-n1nbc(z-1) = bc + bc z-1 + bc z-n01ncccc(z-1) =1+ c z-1 + c z-n1n假定预测误差e(k)服从均值为0的高斯分布，并设序列t?(k)具

6、有相同的方差e。因为4(k)与c(z-1)，a(z-1 )和b(z-1)有关，所以e是被估参数0的函数。为了书写方便,把式(2.2.9)写成c(z-1)e(k) = a(z-1)y(k) -b(z-1)u(k)2.2.10)或写成e(k) = y(k) + a y(k -1) + + a y(k - n) - b u(k -1) - b u(k -1) - -1n01b u(k - n) - c e(k -1)c (k - n), k = n +1, n + 2,n12.2.11)e(k) = y(k) + 工a y(k - i) -工bu(k - i) -工ce(k - i)iiii=1i=

7、0i=1令k=n+l,n+2,n+N,可得e(k)的N个方程式，把这N个方程式写成向量-矩阵形式e 二 Y 0(2.2.13)N N N式中y (n +1)e(n +1)Y =Ny (n + 2),eN =e(n + 2)_ y (n + N)_e(n + N),0 y (n)yu (n +1) uy (n +1)一 y(2)u (n + 2)u (2)y (n + N 1)一 y (N)u (n + N)u (N)e(n) e(1)e(n +1)e(2)e(n + N 1)e(N)因为已假定伙是均值为0的高斯噪声序列，高斯噪声序列的概率密度函数为exp1 (2kq 2)212c 2( y m

8、)22.2.14)式中y为观测值，c 2和m为y的方差和均值，那么exp1 (2kq 2)212c 2e2(k)2.2.15)对于e(k)符合高斯噪声序列的极大似然函数为L(S 0，c) = Le(n +1), e(n + 2),e(n + N) 0 = f e( n +1) 0 f e(n + 2) 0 f e(n + N) 0 1expe2(n +1) + e 2(n + 2) Hb e2(n + N) =1 1 N2c 2N-(2 兀c 2)2(2 兀c 2)2exp(eTe )2c 2 N N2.2.16)1L(Y 0 ,c)二expNN(2兀c 2)22c22.2.17)对上式(2.

9、2.17)等号两边取对数得1lnL(Y 0,c) = ln+ lnexp(一ereNN2c 2 NN(2 兀c 2)2)一Nln2K Nlnc2 22 一eTe2c 2 N N或写为2.2.19)2.2.20)则式中Ac2=丄艺e2(k)k = n +1討艺e2仗)=吕Jk = n +12.2.21)In L(Y p, c)二一In 2兀一In c 2 一 艺 e2(k)N222c 2k =n+1求lnL(Ynp,c)对c2的偏导数，令其等于0,可得空LY四=一丄+丄垃Ne2(k) = 0 dc 22c 2 2c 4k = n +1k = n +1= min J2.2.22)c 2越小越好，因

10、为当方差c2最小时，e 2( k )最小，即残差最小。因此希望c 2的估值取最2.2.23)因为式(2.2.10)可理解为预测模型，而e(k)可看做预测误差。因此使式(2.2.22) 最小就是使误差的平方之和最小，即使对概率密度不作任何假设，这样的准则也是有意义的 因此可按J最小来求a , a, b，b , c,c的估计值。1,0n 1n由于e(k)式参数a ,a,b，b ,c ,c的线性函数，因此J是这些参数的二次型函数。1,0n 1n求使lnL(Y P,c)最大的(A，等价于在式(2.2.10)的约束条件下求$使J为最小。由于J对c 是非线性的，因而求J的极小值问题并不好解，只能用迭代方法

11、求解。求J极小值的常 i用迭代算法有拉格朗日乘子法和牛顿-拉卜森法。下面介绍牛顿-拉卜森法。整个迭代计算步骤如下：(1) 确定初始的p0值。对于p0中的a ,a,b，b可按模型001,0ne(k) = a (z -1) y (k) - b( z-1)u (k)(2.2.24)用最小二乘法来求，而对于 0中的1，cn可先假定一些值。(2) 计算预测误差e(k) = y(k) y(k)(2 2 25)给出J = 1 龙e2( k)2并计算2 二龙e2(k)Nk-n+1JQ 2 J(3)计算J的梯度60和海赛矩阵,有Q0 2QJ垃n (k)Qe(k)6060k - n+1式中2.2.26)2.2.2

12、7)6e(k)Qe(k)60 Qa16e(k) Qe(k)Qe(k) Qe(k)QaQb Qb Qcn0n1Qe(k)QcnQe(k)Q y(k) + a y(k -1) + a y(k -n) -b u(k) -bu(k -1)b u(k -n)-01QaiQaic e(k -1)c e(k - n)1n-y(k - i) - c Qe(k -1)1Qe(k -2)-c Qa2iQaiQe(k - n)Qai2.2.28)Qe(k) _Qaiy(k-i)-Ycjj-1Qe(k j)Qai2.2.29)同理可得詈-u (k - i)丄 cjij -1Qe(k j)Qbi2.2.30)込-e(k

13、-i)-Yc Qe(k-j)Qcj Qcij -1i2.2.31)将式(2.2.29)移项化简，有y(k - i)-沁 + Zc 吨-j)-Yc 吨-j) QajQaij -1ijQ aj -0i2.2.32)因为由e(k - j)求偏导，故2.2.34)de(k j)de(k) z - j将（2.2.34）代入（2.2.32），所以y(k-i) _c j) _c de(k)z-j _空2工jj_0jdaj _0ic z-jdaji j _ 02.2.35)c ( z T)= 1 + c z-1 +1所以得2.2.36)(xde(k)c( z-1)_ y (k-i)dai同理可得（2.2.30）和（2.2.31）为2.2.37)c( z-1)讐 _-U (k - i)i(门 de(k)2.2.38)c( z -1)_ -e(k-i)dci根据（2.2.36）构造公式c(zj)如l(ij)二 yk - (i - j) - j二 y(k - i)daj2.2.39)将其