极大似然辨识及其MTLAB实现
篇一:matlab极大似然估计
极大似然辨识及其MATLAB实现
摘 要:极大似然参数估计方法是以观测值的出现概率为最大作为准则的,这是一种很普遍的参数估计方法,在系统辨识中有着广泛的应用。本文主要探讨了极大似然参数估计方法以及动态模型参数的极大似然辨识并且对其进行了MATLAB实现。 关键词:极大似然辨识 MATLAB仿真 迭代计算
1 极大似然原理
设有离散随机过程{Vk}与未知参数有关,假定已知概率分布密度f(Vk)。如果我们得到n个独立的观测值V1,V2,…,Vn,则可得分布密度f(V1),f(V2),…,f(Vn)。要求根据这些观测值来估计未知参数,估计的准则是观测值{{Vk}}的出现概率为最大。为此,定义一个似然函数
L(V1,V2,,Vn)f(V1)f(V2)f(Vn)
(1.1)
上式的右边是n个概率密度函数的连乘,似然函数L是的函数。如果L达到极大值,{Vk}
的出现概率为最大。因此,极大似然法的实质就是求出使L达到极大值的的估值。为了
便于求,对式(1.1)等号两边取对数,则把连乘变成连加,即
n
lnL
ln
i1
f(Vi) (1.2)
由于对数函数是单调递增函数,当L取极大值时,lnL也同时取极大值。求式(1.2)对的偏导数,令偏导数为0,可得
lnL
0
(1.3)
解上式可得的极大似然估计ML。
2 系统参数的极大似然估计
设系统的差分方程为
11
a(z)y(k)b(z)u(k)(k) (2.1) 式中
a(z)1a1z
11
1
...anz
n
n
b(z)b0b1z
1
...bnz
因为(k)是相关随机向量,故(2.1)可写成
a(z)y(k)b(z)u(k)c(z)(k) (2.2) 式中
c(z)(k)(k) (2.3)
1
111
c(z1)1c1z1cnzn (2.4)
111
(k)是均值为0的高斯分布白噪声序列。多项式a(z),b(z)和c(z)中的系数
a1,,a,b0,bn,c1,cn和序列{(k)}的均方差都是未知参数。
设待估参数
[a1an b0bn c1cn (2.5) 并设y(k)的预测值为
T
y(k)a1y(k1)any(kn)b0u(k)bnu(kn)
c1e(k1)cne(kn) (2.6)
式中e(ki)为预测误差;ai,bi,ci为ai,bi,ci的估值。预测误差可表示为
e(k)y(k)y(k)y(k)
n
n
a
i1
i
n
y(ki)
bu(ki)
i
i0
i1
1n1n
cie(ki)(1a1zanz)y(k)(b0b1zbnz)u(k)
(c1z或者
1
c2z
2
cnz
n
)e(k) (2.7)
(1c1z
1
cnz
n
)e(k)=(1a1z
1
1
anz
n
)y(k)
(b0b1z因此预测误差e(k)满足关系式
bnz
n
)u(k) (2.8)
c(z)e(k)a(z)y(k)b(z)u(k) (2.9) 式中
1
1
1
a(z)1a1z
11
anz
n
2
b(z)b0b1z
1
matlab极大似然估计。
1
bnz
n
c(z)1c1z
1
1
cnz
n
假定预测误差e(k)服从均值为0的高斯分布,并设序列e(k)具有相同的方差。因为e(k)与c(z),a(z)和b(z)有关,所以是被估参数的函数。为了书写方便,
1
1
1
2
把式(2.9)写成
c(z)e(k)a(z)y(k)b(z)u(k) (2.10)
e(k)y(k)a1y(k1)any(kn)b0u(k1)b1u(k1)
1
1
1
bnu(kn)c1e(k1)cn(kn),kn1,n2, (2.11)
或写成
n
n
n
e(k)y(k)aiy(ki)biu(ki)cie(ki) (2.12)
i1
i0
i1
令k=n+1,n+2,…,n+N,可得e(k)的N个方程式,把这N个方程式写成向量-矩阵形式
eNYNN (2.13) 式中
a1e(n1)
ane(n2) , b0e(nN)
bn
y(n1)
y(n2)
,e YNN
y(nN)
y(n)
y(n1)
y(nN1)
e(1)
e(2)e(n1)y(2)u(n2)u(2)
y(N)u(nN)u(N)e(nN1)e(N)y(1)
u(n1)u(1)
e(n)
N
因为已假定e(k)是均值为0的高斯噪声序列,高斯噪声序列的概率密度函数为
f
1
1
exp[
12
2
(ym)]
2
(2
2
2
)
2
(2.14)
式中y为观测值,和m为y的方差和均值,那么
f
1
1
exp[
12
2
e(k)] (2.15)
2
(2
2
)
2
对于e(k)符合高斯噪声序列的极大似然函数为
L(YN,)L[e(n1),e(n2),,e(nN)]f[e(n1)]f[e(n2)]f[e(nN)]
1
N
exp{
12
2
[e(n1)e(n2)e(nN)]}
222
1
N
exp(
12
2
eNeN)
T
(2
2
)2(2
2
)2
(2.16)
或
L(YN,)
1
N
T
exp[
(YN)(YN)
2
2
] (2.17)
(2)2
2
对上式(2.17)等号两边取对数得
lnL(YN,)ln
1
N
lnexp(
12
2
eNeN)
T
N2
ln2
N2
ln
2
12
2
eNeN
T
(2
2
)2
(2.18)
或写为
lnL(YN,)
N2ln2
N2ln
2
nN2
12
e
kn1
2
(k) (2.19)
求lnL(YN,)对2的偏导数,令其等于0,可得 则
lnL(YN,)
2
N2
2
12
4
nN
e
kn1
2
(k)0 (2.20)
式中
2
1N
nN
e(k)
2
21N2
nN
e(k)
2
2N
J (2.21)
kn1kn1
J
2
2
12
nN
e
kn1
2
(k) (2.22)
2
2
越小越好,因为当方差最小时,e(k)最小,即残差最小。因此希望的估值取最
小
2
2N
minJ (2.23)
因为式(2.10)可理解为预测模型,而e(k)可看做预测误差。因此使式(2.22)最小就是使误差的平方之和最小,即使对概率密度不作任何假设,这样的准则也是有意义的。因此可按J最小来求a1,,a,b0,bn,c1,cn的估计值。
由于e(k)式参数a1,,a,b0,bn,c1,cn的线性函数,因此J是这些参数的二次型函数。
求使lnL(YN,)最大的,等价于在式(2.10)的约束条件下求使J为最小。由于J对
ci是非线性的,因而求J的极小值问题并不好解,只能用迭代方法求解。求J极小值的常用
迭代算法有拉格朗日乘子法和牛顿-拉卜森法。下面介绍牛顿-拉卜森法。整个迭代计算步骤如下:
(1)确定初始的0值。对于0中的a1,,a,b0,bn可按模型
e(k)a(z)y(k)b(z)u(k) (2.24)
1
1
用最小二乘法来求,而对于0中的
(2)计算预测误差
c1,cn
可先假定一些值。
e(k)y(k)y(k)
(2.25)
给出
J
12
nN
并计算
e
kn1
2
(k)
2
1N
nN
J
e
kn1
2
(k)
(2.26) J
22
(3)计算J的梯度 和海赛矩阵
J
nN
,有
式中
e(k)
e(k)
(2.27)
kn1
T
e(k)
e(k)e(k)e(k)e(k)e(k)e(k)
aabbcc1n0n1n
e(k)ai
ai
[y(k)a1y(k1)any(kn)b0u(k)b1u(k1)bnu(kn)
c1e(k1)cne(kn)]
y(ki)c1即
同理可得
e(k)bie(k)cie(k)ai
e(k1)ai
c2
e(k2)ai
cn
e(kn)ai
(2.28)
n
y(ki)cj
j1
e(kj)ai
(2.29)
n
u(ki)cj
j1n
e(kj)bie(kj)ci
(2.30)
e(ki)cj
j1
(2.31)
将式(2.29)移项化简,有
y(ki)
e(k)ai
n
j1
cj
e(kj)ai
n
j0
cj
e(kj)ai
(2.32)
极大似然估计
篇二:matlab极大似然估计
极大似然估计(maximum likelihood estimination)
极大似然估计法是求点估计的一种方法,最早由高斯提出,后来费歇尔(Fisher)在1912年重新提出。它属于数理统计的范畴。
大学期间我们都学过概率论和数理统计这门课程。
概率论和数理统计是互逆的过程。概率论可以看成是由因推果,数理统计则是由果溯因。
用两个简单的例子来说明它们之间的区别。
由因推果(概率论)
例1:设有一枚骰子,2面标记的是“正”,4面标记的是“反”。共投掷10次,问:5次“正”面朝上的概率?
解:记“正面”朝上为事件A,正面朝上的次数为x。
有题意可知:P A =3。
115P x=5 =
http://m.myl5520.com/mingrenmingyan/91894.html