生物医学信号处理-42-估计理论课件.ppt

1、最简单且常见的情况是取得一组随机信号的样本后估计信号的主要统计特性，如均值、均方、方差、相关函数、功率谱等。对于前三者可以采用“直接估计”法，即按照定义，用有限样本来估计：估计的任务就是从含噪信号中估计出信号的某个参数或某些参数估计的任务就是从含噪信号中估计出信号的某个参数或某些参数1221211()1()1()()limlimlimNnNnNnNnNnNnE xxNE xxNVar xxE xN12121111()()NxnnNxnnNnnmxNDxNVar xxE xN均值：均方：方差：均值估计：均方估计：方差估计：4.2 4.2 估计理论估计理论 N=5000;N=5000;y=rand

2、n(1,N);y_junzhi=mean(y)y_junzhi=0.0116 sum(y)/N ans=0.0116 y_junfang=y*y/N y_junfang=0.9988 sum(y.2)/N ans=0.9988 y_fangcha=sum(y-y_junzhi).2)/N y_fangcha=0.9986 N=50000;N=50000;y=randn(1,N);y_junzhi=mean(y)y_junzhi=-0.0017 y_junfang=y*y/Ny_junfang=1.0056 y_fangcha=sum(y-y_junzhi).2)/Ny_fangcha=1.00

3、56更复杂的情况是从含噪的观察中估计信号的参数，如下图：非线性估计：非线性估计：以待估参数a的先验概率知识f(a)和条件先验概率 f(x|a)为基础，依据某些最优判据，通过非线性的 数理统计算法做出估计。线性估计：线性估计：限定估计必须是观察值x的线性函数，并在最小均方 误差意义下做估计。估计质量：估计质量：估计的偏差偏差：若待估参数a和其估计值的均值E()相等，就称无偏估计。偏差b=E()-a估计的方差方差：表示各次估计值相对于估计值均值E（）的 分散程度，即2=E-E()2实际估计中偏差和方差往往有矛盾，两者结合起来才能较全的表示估计质量。通常用D=2+b2表示估计的好坏，可以证明，D就是

4、估计值的均方误差均方误差：D=E(a-)2有效估计：若用某方法得出的估计值的方差小于任何其他估计 方法得出的方差值，则称此估计为有效估计。一致估计：如果随着样本数目的增加，估计的均方误差趋于0 （要求当N+时，偏差和方差都趋于0），则称 此估计为一致估计。2211()()kkE aE aE aE a2()0limNE aa问题提出：设观察x=s+n，s是随机变量，f(s)，f(x|s)都已知，现做单次观察x=x1，试根据一定的估计准则由x1做出s的估计值。极大似然估计贝叶斯估计：对不同估计结果给予不同代价，并要求平均代价最小。通常代价c定义为估计误差(s-)的单变量函数。mincE c,ccf

5、 x s dsdx ,|f x sf s x f x|ccf s x ds f x dx|Rcf s x dsc为x和s的联合函数内积分和f(x)都非负条件风险,条件平均代价贝叶斯估计贝叶斯估计|Rcf s x ds2css2|Rssf s x ds2|0dRss f s x dsds|sf s x dssf s x ds|sf s x dssf s x dss|MSssf s x ds(|)MSsE s x|f s f x sf s xf x|1|MSsf s f x s dsssf s f x s dsf xf s f x s ds（1）均方误差估计条件均值|Rcf s x ds|css|

6、ssRss f s x dsss f s x ds|ABSABSssf s x dsf s x ds|f s f x sf s xf x|ABSABSssf s f x s dsf s f x s ds对求导，并令结果为0（2）绝对值估计ABS2222|1|ssssRf s x dsf s x dsf s x ds ln|0MAPs sf s xs|f s f x sf s xf x ln|lnln|lnf s xf sf x sf x lnln|0MAPs sf sf x sss（3）均匀估计 MAP为了使上式最小,显然估计值应取在后验概率密度函数最大之处，极大后验概率估计(也叫MAP估计)

7、极大似然估计极大似然估计1）在对s没有先验知识的情况下，只能用其代替极大后验概率估计2）若s是未知的确定性量，则谈不上先验概率f(s)，此时则用ML估计|0MLs sf x ssln|0MLs sf x ss似然函数就是条件先验概率密度函数f(x|s)。极大似然估计的准则是,以使似然函数f(x|s)最大的s值作为估计值 ML(也叫ML估计)|ln|0|0MLMLs ss sfsfsss 或 XX极大似然估计 lnln|0MAPs sf sfsssX极大后验概率估计 MAP|MSsf s fs dssf s fs dsXX最小均方误差估计12,.,Mx xxX观察是矢量情况观察是矢量情况222e

8、xp()(),020,ssssf ss其他221()(|)exp()22nnxsf x s 22lnln|10snf sf x ssxssss220nxssaa221nsa（）224(11)2nMAPaxsax例1 x=s+n,n为噪声,做正态分布N(0,n2),s为随机变量,呈瑞利分布：则f(x|s)是均值为s的高斯分布：2ln|0nf x sxssMLsx|MSsf s f x s dssf s f x s dsABSMAPMSssssm如要做ML估计,则如果要做MS估计,则如果有理由判定后验概率密度函数具有高斯分布的形式,如左图所示,则必有因为高斯函数的中值、峰值和条件均值都处在同一横坐

9、标上221()exp()22sssf s221221()11(|)exp()()exp()2222MkMMkkknnnnxsxsfsX 22122|11|()exp()222MkMknsnxsf s fssf sffXXXX222211|()exp()2Mpkkpnf slsxXX22222snpsnM 例2 由M次观察xk=s+nk(k=1,2,M)做MS、ABS、MAP、ML各估计。式中,s为待估计的随机变量,做正态分布N(0,s2);nk为互相独立的观察噪声,都做正态分布N(0,n2)展开此式的指数,必定是s的二次函数。因此可对s配平方,并将只与xk有关而与s无关的剩余部分并入系数l(X

10、 X)中,得式中221MpABSMAPMSkknsssx2121(|)()exp()22MkMknnxsfsX11MMLkksxM在给定了观察值X后,l(X)便是一个给定常数。因此后验概率密度函数f(s|X)显然是高斯型的。这样便可直接得出结论:再求ML,因为取对数并对s求导,令其等于0,得 可见,样本均值便是极大似然估计。值得一提的是,做极大似然估计时不考虑f(s)。因此上面所得结论的适用范围很广,不管信号本身的概率分布如何,只要加性噪声做髙斯分布则样本均值总是信号的极大似然估计。(做MAP估计时要注意概率密度函数的定义域)(|)MSaaf a x da()MSMSE aaf x dx(|)

11、()af a x f x dxda(,)af a x dxda()af a da E a均方估计的无偏性下面证明均方估计是无偏估计。按定义,MS是a在后验概率密度函数意义下的均值：但MS由观察值x决定,所以其均值概述概述不需要关于概率密度的知识只要知道信号和噪声的一二阶统计知识估计算法是观察值的线性函数估计的判据一般采用二次代价函数，即均方误差最小。2,cssf x s dsdx2,csaxbf x s dsdx=s axb,0,0 x saxb f x s dsdxsaxb f x s dsdx 22,(),()xf x s dsdxxf x dxE xx f x s dsdxE x,(),

12、()sf x s dsdxE sxsf x s dsdxE xs2()()()()()aE xbE xE xsaE xbE s22()()()()()()()E xsE s E xaE xExbE saE x当采用误差平方为代价函数时对a，b求导，并令导数为0线性均方估计设观察的逐次样本是xi=s+ni，i=1,2,N。现在用观察值(x1,x2,xN)的线性组合来得到s的估计值：使得与信号s间的均方误差最小:求定各系数hj。线性最小均方(LMS)估计,以下记作 LMSE(e2)对各hj求导,并令导数等于01=Njjjsh x2minE ss221()()NjjjE eE sh x221()()

13、2 2()0()0,1,2,.,iiiiNjjijE eeeEE eE exhhhE sh x xiN1212111212122212.NNNNNNNNHhhhGgggrrrrrrRrrr1HR G1()(),1,2,.,ijijiiNj ijijrE x xgE sxh rg iN1()0,1,2,.,NjjijE sh xxiN推导得出的方程组选择处理器系数的原则是使估计误差e和所有xi(i=1N)正交。这一原则称为正交原理由于估计是xi的线性组合,因此e和各xi正交便导致e和也正交,即采用LMS估计后得到的最小均方误差是可见,正交条件下估计的均方误差是估计误差和信号乘积的均值()0,1,

14、2,.,iiE exE ss xiN()0E esE ss s2minmin()E ssE ss sE ss sE ss sE es递归线性最小均方估计当观察数据数目逐次增加时,总是由上一次估计结果结合本次新输入的数据,把上次估计值做修改得出本次估计。不论平稳或非平稳过程都可应用。(1)问题提法观察值为xj=s+nj，j=1,2,.先验统计知识为要求用递归线性最小均方估计算法由上一次估计k结合本次新观察xk+1,做出本次估计 k+1。估计的判据是均方误差最小,即这一步的关键是求系数ak+1、bk+1。只要系数求得,k+1便可估出。22()()0()()jjkjnkjsE nE snE n nE

15、 sD 1111kkkkksasbx211()minkkE ss(2)基本公式的推导先推导计算ak+1、bk+1的有关公式。先证明线性均方误差估计是无偏估计。再利用这一特点推证系数和的内在联系是ak+1+bk+1=1证明:先证明LMS估计是无偏估计。为简单起见,只用单样本来估计的情况来做证明。由该节所述,此时要用 k=ax+b来做估计。前面已求得因为 E(k)=a E(x)+b化简得出 E(k)=E(s)对更复杂的情况(观察是矢量)也可以类似地证明。22()()()()()()()E xsE s E xaE xExbE saE x再导出ak+1、bk+1间的关系等号两边取均值化简便得再由正交原

16、理证明此时的最小均方误差1111kkkkksasbx1111()()()kkkkkE saE sbE x1()()()kkE sE sE s11()()()()kkE xE sE nE s111kkab221111()()kkkknE ssE sssb因为11()0kkE ssx11111()()kkkkkE sssE ssnE sn 1111111()kkkkkkkkE snE asbxn11121111()()kkkkkkknbE xnbE snnb11kkkbbb2112112112()()()()kknskknksknksknE sssbDE ssbE ssDbE s sDb证明：由正

17、交原理将xk+1=s+nk+1带入得即再经定义证明即或同理1111()()kkkkkE ssE asbxs11111()()()kkkkkkksaE s sbE xsaE s sbD22111()sknksknksDbaDbbD211(1)()ksknksbDbbD11kkkbbb又结合前面结果有化简即初始化可由 0=E(s)，结合正交原理得出：例如，设E(s)=0，则 0=0，因而 1=a1 0+b1x1=b1x1。系数b1可由正交原理决定111 11112221()0()0()1()1snsnsE ss xE sb x xDE sxbE xDD1)估计值 k+1是否会趋于一个极限?因为b1

18、=1/(1+n2/Ds)是一个小于的正数,因此b2=b1/(1+b1)也必然是一个小于1的正数,而且其值比b1更小。类似地推论下去,bk必定愈来愈小,终至趋于0。此时ak=1-bk将趋于1。因此,k很大后,将有可见,在信号平稳的情况下,递归过程最终会收敛。2)由b1、bk+1及k+1各表示式可见,各次估计所用的系数和估计的均方误差都和观察xk无关,因此它们完全可以在开始输入数据之前先算好,存储在计算机内,不必临时计算。这样可以缩短每次递归计算所需时间,但以增加内存为代价。在复杂的递归估计问题中,这一特点对实时工作颇为有利。1111kkkkkksasbxs最小二乘估计为了估计未知纯量s,对它进行

19、N次线性测量,得观察值 式中,cks为用线性测量方法测得的理论值;ck为由测量方法决定的已知系数(由于各次测量方法和仪器未必相同,所以各次测量的ck值未必相同);vk为各次测量的误差。最小二乘估计准则是选s的估计值,使各次测量误差的平方和最小,即估计值可以直接用对求导得出:以上便是最小二乘法的基本思路,1kkkxc sv kN2211()minkNNkkkkvxc s1()0Nkkkkxc s c121kkNkkLSNkc xsc则 最小二乘法不仅用来做估计,在信号处理技术中还有许多其他应用,其中曲线拟合最为典型。事实上高斯最早提出最小二乘法时正是针对天体运动轨道的拟合问题。下面以直线拟合为例

20、加以说明。设有一组观察值(x1,y1),(x2,y2),(xN,yN)。现在要用一条直线来拟合这组数据。设直线方程为 y=kx+m 要求选择k、m之值,使各点和此直线间的估计误差ek的平方和最小2211,1()minkkkkkkNNkkkkeyyykxm kNeykxm解答可从微分求极值的概念入手,将对k与m求导,并令导数等于0,得到一组关于k、m的线性方程组,从而解出k、m,即测得某克山病区10名健康儿童头发与全血中的硒含量，见下表，试用最小二乘法建立由发硒x推算血硒y的经验公式。发硒发硒x x74668869917366965873血硒血硒y y13101311169714510 x=74

21、 66 88 69 91 73 66 96 58 73;y=13 10 13 11 16 9 7 14 5 10;A=sortrows(x;y)A=58 5 66 7 66 10 69 11 73 9 73 10 74 13 88 13 91 16 96 14 x=A(:,1)y=A(:,2)p=polyfit(x,y,1)p=0.2358 -6.9803 pp=polyfit(x,y,2)pp=-0.0062 1.1969 -43.6255 yg=polyval(p,x)yg2=polyval(pp,x)yg=6.6969 8.5834 8.5834 9.2908 10.2340 10.23

22、40 10.4699 13.7712 14.4787 15.6577yg2=5.0949 8.5656 8.5656 9.6641 10.9564 10.9564 11.2487 14.0485 14.3346 14.5653y=y=5 5 7 7 10 10 11 11 9 9 10 10 13 13 13 13 16 16 14 14 figure plot(x,y,k*,x,yg,ko-,x,yg2,k-.)legend(实测值,拟合直线,拟合二次曲线)err1=sum(yg-y).2)err2=sum(yg2-y).2)text(70,8,strcat(y=,num2str(p(1),x,num2str(p(2),err1,num2str(err1)text(60,6,strcat(y=,num2str(pp(1),x2+,num2str(pp(2),x,num2str(pp(3),err2=,num2str(err2)err1=23.9512err2=18.3039