极大似然辨识及其MTLAB实现.doc
《极大似然辨识及其MTLAB实现.doc》由会员分享，可在线阅读，更多相关《极大似然辨识及其MTLAB实现.doc（10页珍藏版）》请在冰豆网上搜索。
极大似然辨识及其MATLAB实现
摘要：
极大似然参数估计方法是以观测值的出现概率为最大作为准则的，这是一种很普遍的参数估计方法，在系统辨识中有着广泛的应用。
本文主要探讨了极大似然参数估计方法以及动态模型参数的极大似然辨识并且对其进行了MATLAB实现。
关键词：
极大似然辨识MATLAB仿真迭代计算
1极大似然原理
设有离散随机过程与未知参数有关，假定已知概率分布密度。
如果我们得到n个独立的观测值…，则可得分布密度，，…,。
要求根据这些观测值来估计未知参数，估计的准则是观测值{}的出现概率为最大。
为此，定义一个似然函数
（1.1）
上式的右边是n个概率密度函数的连乘，似然函数L是的函数。
如果L达到极大值，的出现概率为最大。
因此，极大似然法的实质就是求出使L达到极大值的的估值。
为了便于求，对式（1.1）等号两边取对数，则把连乘变成连加，即
（1.2）
由于对数函数是单调递增函数，当L取极大值时，lnL也同时取极大值。
求式（1.2）对的偏导数，令偏导数为0，可得
（1.3）
解上式可得的极大似然估计。
2系统参数的极大似然估计
设系统的差分方程为
（2.1）
式中
因为是相关随机向量，故（2.1）可写成
（2.2）
式中
（2.3）
（2.4）
是均值为0的高斯分布白噪声序列。
多项式，和中的系数和序列的均方差都是未知参数。
设待估参数
（2.5）
并设的预测值为
（2.6）
式中为预测误差；，，为，，的估值。
预测误差可表示为
（2.7）
或者
=
（2.8）
因此预测误差满足关系式
（2.9）
式中
假定预测误差服从均值为0的高斯分布，并设序列具有相同的方差。
因为与，和有关，所以是被估参数的函数。
为了书写方便，把式（2.9）写成
（2.10）
（2.11）
或写成
（2.12）
令k=n+1,n+2,…,n+N,可得的N个方程式，把这N个方程式写成向量-矩阵形式
（2.13）
式中
,
因为已假定是均值为0的高斯噪声序列，高斯噪声序列的概率密度函数为
（2.14）
式中y为观测值，和m为y的方差和均值，那么
（2.15）
对于符合高斯噪声序列的极大似然函数为
（2.16）
或
（2.17）
对上式（2.17）等号两边取对数得
（2.18）
或写为
（2.19）
求对的偏导数，令其等于0，可得
（2.20）
则
（2.21）
式中
（2.22）
越小越好，因为当方差最小时，最小，即残差最小。
因此希望的估值取最小
（2.23）
因为式（2.10）可理解为预测模型，而e（k）可看做预测误差。
因此使式（2.22）最小就是使误差的平方之和最小，即使对概率密度不作任何假设，这样的准则也是有意义的。
因此可按J最小来求的估计值。
由于e（k）式参数的线性函数，因此J是这些参数的二次型函数。
求使最大的，等价于在式（2.10）的约束条件下求使J为最小。
由于J对是非线性的，因而求J的极小值问题并不好解，只能用迭代方法求解。
求J极小值的常用迭代算法有拉格朗日乘子法和牛顿-拉卜森法。
下面介绍牛顿-拉卜森法。
整个迭代计算步骤如下：
（1）确定初始的值。
对于中的可按模型
（2.24）
用最小二乘法来求，而对于中的可先假定一些值。
（2）计算预测误差
（2.25）
给出
并计算
（2.26）
（3）计算J的梯度和海赛矩阵,有
（2.27）
式中
（2.28）
即
（2.29）
同理可得
（2.30）
（2.31）
将式（2.29）移项化简，有
（2.32）
因为
（2.33）
由求偏导，故
（2.34）
将（2.34）代入（2.32），所以
（2.35）
所以得
（2.36）
同理可得（2.30）和（2.31）为
（2.37）
（2.38）
根据（2.36）构造公式
（2.39）
将其代入（2.36），可得
（2.40）
消除可得
（2.41）
同理可得（2.37）和（2.38）式
（2.42）
（2.43）
式（2.29）、式（2.30）和式（2.31）均为差分方程，这些差分方程的初始条件为0，可通过求解这些差分方程，分别求出e（k）关于的全部偏导数，而这些偏导数分别为，和的线性函数。
下面求关于的二阶偏导数，即
（2.44）
当接近于真值时，e（k）接近于0。
在这种情况下，式（2.44）等号右边第2项接近于0，可近似表示为
（2.45）
则利用式（2.45）计算比较简单。
（4）按牛顿-拉卜森计算的新估值，有
（2.46）
重复
（2）至（4）的计算步骤，经过r次迭代计算之后可得，近一步迭代计算可得
（2.47）
如果
（2.48）
则可停止计算，否则继续迭代计算。
式（2.48）表明，当残差方差的计算误差小于时就停止计算。
这一方法即使在噪声比较大的情况也能得到较好的估计值。
3动态模型参数极大似然辨识及其MATLAB实现
设动态系统的模型表示为
式中，是均值为0，方差为2，服从正态分布的不相关随机噪声；u（k）和z（k）表示系统的输入输出变量。
现给出一系统模型为
z（k）-1.2z（k-1）+0.6z（k-2）=u（k-1）+0.5（k-2）+e（k）
e（k）=v（k）-v（k-1）+0.2v（k-2）
其中v（k）为随机信号，输入信号是幅值为1的M系列或随机信号，试用递推的极大似然法求系统辨识的参数Φ。
程序如下：
clear
a
（1）=1;b
（1）=0;d
（1）=0;u
（1）=d
（1）;z
（1）=0;z
（2）=0;%初始化
fori=2:
1200%产生m序列u（i）
a（i）=xor（c（i-1）,d（i-1））;
b（i）=a（i-i）;
c（i）=b（i-1）;
d（i）=c（i-1）;
u（i）=d（i）;
end
u;
v=randn（1200,1）;%产生正态分布随机数
V=0;%计算噪声方差
fori=1:
1200
V=V+v（i）*v（i）;
end
V1=V/1200;
fork=3:
1200%根据v和u计算z
z（k）=1.2*z（k-1）-0.6*z（k-2）+u（k-1）+0.5*u（k-2）+v（k）-v（k-1）+0.2*v（k-2）;
end
o1=0.001*ones（6,1）;p0=eye（6,6）;%幅初值
zf
（1）=0.1;zf
（2）=0.1;vf
（2）=0.1;vf
（1）=0.1;uf
（2）=0.1;uf
（1）=0.1;
%迭代计算参数值和误差值
fork=3:
1200
h=[-z（k
升级会员 