9现代信号处理-功率谱估计(2)

七、最大熵谱估计1、利用最大熵的原则外推自相关函数2、最大熵谱估计与AR模型谱估计的等价性八、最大似然谱估计1、最小方差谱估计2、最大熵谱与最大似然谱估计的关系九、特征分解法谱估计1、正弦波用退化AR模型表示2、白噪声中正弦波组合用一特殊的ARMA模型表示3、特征分解法谱估计功率谱估计十、Prony谱分析法1、利用最大熵的原则外推自相关函数2、最大熵谱估计与AR模型谱估计的等价性十一、多重信号分类MUSIC1、最小方差谱估计2、最大熵谱与最大似然谱估计的关系十二、特征分解法谱估计1、波束形成器2、特征子空间分析3、MUSIC算法及其改进功率谱估计一、最大熵谱估计1.利用最大熵的原则外推自相关函数按照Shannon对熵的定义，当随机变量X取离散值时，熵的定义为(4.6.1)式中pi是出现状态i的概率。当X取连续值时，熵的定义为(4.6.2)iiippHln)]([lnd)(ln)(xpExxpxpH式中,p(x)是X的概率密度函数，对于离散随机序列,概率密度函数用联合概率密度函数代替。显然,熵代表一种不确定性，最大熵代表最大的不确定性，或者说最大的随机性。下面我们研究对于有限的自相关函数值不作任何改变，对于未知自相关函数用最大熵原则外推，即不作任何附加条件的外推方法。假设x(n)是零均值正态分布的平稳随机序列，它的N维高斯概率密度函数为XNRXNRxxxpxxHxxNN12/12/21))((21exp))((det)π2(),,,(式中H21],,,[NxxxX)0()1()()1()0()1()()1()0()(xxxxxxxxxxxxxxxxxxxxrNrNrNrrrNrrrNR按照（4.6.2）式，x(n)信号的熵为（4.6.3）式中det(Rxx(N))表示矩阵Rxx(N)的行列式，由上式表明为使熵最大，要求det(Rxx(N)最大。])))((det()πe2log[(2/12/NRHxxN若已知N+1个自相关函数值rxx(0),rxx(1),…,rxx(N)，下面用最大熵方法外推rxx(N+1)。设rxx(N+1)确实是信号自相关函数的第N+2个值，根据自相关函数的性质，由N+2个自相关函数组成的矩阵为)0()1()()1()()1()0()1()1()()1()0()1(xxxxxxxxxxxxxxxxxxxxxxxxxxrrNrNrNrNrrrNrNrrrNR（4.6.4）它必须是非负定的矩阵，即0)]1(det[NRxx（4.6.5）0)](det[,),1(det[),0(det[NRRRxxxxxx将行列式展开，det(Rxx(N+1))是rxx(N+1)的二次函数，该二次函数系数的符号是：(-1)1+N+2(-1)1+N+1=-1，且det(Rxx(N+1))对rxx(N+1)的二次导数是-2det［Rxx(N-1)］，它是负值，负值表示det(Rxx(N+1))对rxx(N+1)的一次导数是减函数，det(Rxx·(N+1))作为rxx(N+1)的函数，凹口向下，那么只有一个最大值。为选择rxx(N+1)使det(Rxx(N+1)最大，解下列方程：0)]1(det[)1(ddNRNrxxxx(4.6.6)用数学归纳法，得到（4.6.7）上式是rxx(N+1)的一次函数，可以解出rxx(N+1)。继续再将rxx(N+1)代入Rxx(N+2)和det(Rxx(N+2))中，求det(Rxx(N+2))对rxx(N+2)的最大值，得到rxx(N+2);以此类推，可推出任意多个其它自相关函数值，而不必假设它们为零，这就是最大熵谱估计的基本思想。0)1()()1()2()1()2()1()0()1(xxxxxxxxxxxxxxxxxxrNrNrNrrrNrrr2.最大熵谱估计与AR模型谱估计的等价性我们已经知道ARYule-Walker方程，即将m≥1的情况写成矩阵形式：NkwxxkNkxxkxxkmrakmramr121)()()(m0m=0001)0()1()()1()0()1()()1()0(21NxxxxxxxxxxxxxxxxxxaarNrNrNrrrNrrr式中ai是AR模型系数,i=1,2,3,…,N,。在AR模型中,列写齐次方程式,可得（4.6.8）0)0()1()(0)2()1()2(0)1()0()1(111xxNxxxxxxNxxxxxxNxxxxraNraNrNrararNrarar及0)1()()1(1xxNxxxxraNraNr利用N个参数,由齐次方程组即可解得a1,a2,…,aN值,再将得到的参数值代入(4.6.8)式,并将它整理成行列式:0)1()()1()2()1()2()1()0()1(xxxxxxxxxxxxxxxxxxrNrNrNrrrNrrr可以看出AR模型得到的结果与按最大熵外推rxx(N+1)得到的结果一致，这就证明了当x(n)为高斯分布时的最大熵谱估计与AR模型法是等价的。上式（4.6.8)是rxx(N+1)的一次函数，由此可解得rxx(N+1)。再用类似的方法求得rxx(N+2)，rxx(N+3)，┄，然后确定功率谱估计。最大熵谱估计用下式计算信号功率谱：NkkjkwjxxaP12e1)e(（4.6.9）二、最大似然谱估计１、最小方差谱估计最大似然谱估计是用一个FIR滤波器实现，该滤波器对所关心频率的正弦信号，可以无失真地通过，而对于其它频率的信号，让其频响尽可能地小，亦即将它们尽可能地滤除。此时,滤波器输出的均方值，就作为信号的功率谱估计。设实信号用x(n)表示，FIR滤波器系统函数用A(z)表示：pnnznazA0)()(pkpAXknxkananxny0T)()()()()(输出y(n)为(4.6.10)式中TT)](,),1(),([)](,),1(),0([pnxnxnxXpaaaAp输出信号的均方值为ppHppTHppTHpyARAAXXEAAXXAEnyEP][][]|)([|ˆ2(4.6.11)上式中T表示转置,H表示共轭转置,Rp=E［XXT］是Toeplith自相关矩阵,为求,必须先求FIR滤波器的系数。求这些系数的原则是：在所关心频率ωi处，信号x(n)无失真地通过，即在ωi处的传输函数为1：yPˆpiHpnpnAEAnaAiii)()(e1e)()e(jj-0j式中Hjj2j]eee1[pωωωHpiiiE(4.6.12)另外一个原则是在ωi附近的频率分量尽量衰减掉，即ω≠ωi处,滤波器输出y(n)的均方差最小,即(4.6.11)式最小,此时作为信号x(n)的功率谱估计。因此,最大似然谱估计称为最小方差谱估计更为合适，但由于习惯也可以仍称为最大似然谱估计。在以上原则下，使方差最小的滤波器系数和分别为［30］、［31］yPˆyPˆxPˆyPˆxPˆ应该指出,此时并不是真正意义上的信号功率谱,只是描述了信号功率谱的相对强度。yPˆ)()()(1ˆˆ)()()(111ippiHpyxipiHpipppEREPPREERANiMEMMLMiPNnP11),(/11),(2、最大熵谱与最大似然谱估计的关系伯格证明了最大熵谱PMEM与最大似然谱PMLM估计的关系从上式可知最大似然谱估计相当于从最大熵谱估计的最低分辨率到最高分辨率的平均，所以最大熵谱估计的分辨率要比最大似然谱估计的分辨率高。但最大似然谱具有更大的统计稳定性，对模型阶数的依赖性要小于最大熵谱估计。另外在最大熵谱估计中提到，它的最大缺点就是求得最佳频率成分后，其相应的振幅值并不代表原来的振幅值，尚须用其他办法来近似确定。通过其他两位同学的介绍，我们知道，频谱估计中，振幅谱常用傅立叶变换（传统法）求得，功率谱可通过振幅谱的平方求得，另外也可通过自相关函数的傅立叶变换求得。随机信号一般只作功率谱估计，所以功率谱估计在谱估计中占有重要地位。它的主要缺点是失去了相位信息，因此光靠功率谱是无法恢复信号的。三、特征分解法谱估计4.7.1正弦波用退化AR模型表示无论是实正弦波还是复正弦波，都可以用一个退化AR模型表示，设P个实正弦波组成的信号用下式表示：)(sin)(1iiPiinqnx（4.7.1）式中，初相位θi是在区间（-π，π）均匀分布的随机变量，首先分析下面的三角恒等式：])1(sin[cos2])2(sin[)sin(nnn-π＜ω＜π令x(n)=sin(ωn+θ),则上式变为)2()1(cos2)(nxnxnx（4.7.2）将上式进行Z变换，得到0)cos21)((21zzzX（4.7.3）这样(4.7.2)式的特征多项式为0cos2121zz（4.7.4）上式的两个根分别是:z1=ejω,z2=e-jω,它们共轭成对，且模为1。由这两个根可以确定正弦波的频率。对比AR模型的系统函数，可以把正弦波信号用一个特殊的AR（2）模型表示，括弧中的2表示模型是二阶的。该AR模型的激励白噪声方差趋于0，极点趋于单位圆。通常称为退化的AR模型。这一模型系数有两个，即2cosω和1，（4.7.2）式是模型的差分方程。对于P个实正弦波,特征多项式是0)cos21(121piizz（4.7.5）上式是z-1的2P阶多项式，可以表示为10020azaPkkk（4.7.6）注意上式中的系数ak(k=1,2,3,…,2P),必须保证它的根共轭成对。考虑到根共轭成对，也可表示为PiiiPkkkzzzzza1*20))((（4.7.7）这样由(4.7.6)式，P个正弦波组合的模型用下面2P阶差分方程描述Pkkknxanx21)()(（4.7.8）对于复正弦波情况，P个复正弦波组成的信号是Piniiiqnx1)j(e)(（4.7.9）用一个退化的AR（p）模型表示的差分方程为pikn-kxanx1)()(（4.7.10）其特征多项式为1000azaPkkk(4.7.11)其根为ijize1≤i≤P注意这里的根不是共轭成对出现的。总结以上P个正弦波组合是一个退化的AR(2P)过程，独立参量个数为P个；P个复正弦波的组合是退化的AR（P）过程，独立参量个数仍为P个。实正弦过程相应的退化AR过程的阶数比复正弦情况的阶数高1倍。4.7.2白噪声中正弦波组合用一特殊的ARMA模型表示白噪声中正弦波组合的信号为)()sin()()()(1nwnqnwnxnyPiiii（4.7.12）式中，w(n)为白噪声,且0)]()([,)]()([,0)]([,2lwnxElwnwEnwElnw将(4.7.12)式中x(n)的用AR(2P)表示，即将(4.7.8)式带入(4.7.12)式中，得到)()()(21nwinxanyPii（4.7.13）将（4.7.12）式中的n用n-i代替,x(n-i)=y(n-i)-w(n-i)再将上式带入(4.7.13)式,得到PiiPiiinwanwinyany2121)()()()((4.7.14)上式可以看成一个特殊的ARMA(2P,2P)模型，括弧中的两个2P分别表示ARMA模型系统函数分子和分母的阶次。它与一般的ARMA模型比较，有三方面不同：(1)它的AR部分和MA部分具有相同