5第五章现代谱估计

第五章现代谱估计主要内容：5.1谱估计的发展5.2AR模型参数法5.3最大熵谱估计5.4最大似然谱功率估计5.5谐波分解模型5.1谱估计的发展19世纪末舒斯特（Schuster）周期图时间序列由角频率的正弦信号与噪声叠加而成。则周期图（寻找数据的隐周期性即频率）在处会出现峰值。通过计算周期图。由各峰值可显示出正弦频率信号。1930年，维纳-辛钦定理,证明自相关函数和功率谱互为傅立叶变换，建立了使用傅氏方法处理随机过程的理论体系。谱分析的第二步。22121...jjjNNSxexexeN001958年，布莱克曼（Blackman）和图基（Tukey）经典论文“由通信工程观点对功率谱的测量”给出用维纳相关法从抽样序列得到功率谱的实现方法——BT法。其性能与窗函数选择有关。周期图和BT法称为经典谱估计方法。（且是线性估计方法）。上述方法的最大问题是由于数据截断（或开窗）带来的频率泄漏。弱信号的主瓣很容易被强信号的旁瓣所淹没。对于短序列这一情况尤为突出。1967年以后为现代谱估计阶段，是基于非线性的谱分析。包括各种模型参数法，熵谱法及最小方差法，与奇异值/特征值分解处理的非参数法等。当数据序列较长时，经典谱估计有较好的性能，但对于短序列，经典谱估计存在分辨率不高等致命弱点。现代谱估计具有优良性能，为经典法所不及。现代谱估计的主要方法有：1、模型参量谱估计——可得到高分辨率的谱估计。而这取决于假定模型对观察数据的适配能力。2、非参量谱估计不用有限参数描述的信号模型，直接由自相关延迟序列得到。高信噪比下不如模型法，但在低信噪比下，模型参量谱估计的分辨率大为下降。1973年，皮萨伦科（Pisarenko）提出特征矢量法，开辟了基于自相关矩阵或数据矩阵进行特征分解的非参量谱估计。3、熵谱估计1967年，伯格提出最大熵谱分析法。其方法是对已知延迟点上的自相关函数不加修改，而是对未知延迟点上的自相关函数按信息论中的最大熵外推而得。1971年，范登博斯证明，最大熵谱分析I法与AR模型参量法等效。因此也将该法列为模型参量法谱估计。1979年，美国海军实验室绍尔（Shore）提出最小交叉熵谱分析I法。以后又有最大熵谱分析II法、最小交叉熵谱分析II法、最大熵分析拓广、最小交叉熵分析拓广、多信号交叉熵谱分析等，已可独立成为现代谱估计的重要组成部分。4、多谱（高阶谱）与多维谱估计通常的谱估计或功率谱估计，只包含振幅信息。但许多实际问题需要相位信息。这种情况可考虑使用多谱（高阶谱）——高阶累积量的多维傅氏变换。5、自适应谱估计和鲁棒估计,lnqxHpqqxdxpx5.2AR模型参量法模型参量法的基本思想是根据所研究信号的先验知识，对信号窗外的数据作出某种比较合理的假设（外推或预测）。具体步骤为：1、选择一个好的模型。输入为或或白噪声情况下，使输出等于所研究的信号，至少也是对该信号的一个近似。AR,MA,ARMA及谐波信号模型。2、利用已知的自相关函数或数据求模型参数。3、利用求出的模型参数估计该信号的功率谱。nt5.2.1尤拉-沃克（Yule-Walker）方程5.2.2莱文森-德宾（Levison-Durbin）算法5.2.3BURG算法5.2.4AR模型阶数的判定5.2.1尤拉-沃克（Yule-Walker）方程AR（p）模型(1)1pnlnlnlxax(2)此式(2)可简写成：*11121,0,0xxnknplnklnknlplxxnknlplxxlplxxlrkExxEaxxarklExarklkarklk1()()pxxlxxlrkarklk取k=0,1,2,…,p，可得矩阵表达为：(4)211(0)(1)()(1)(0)(1)0(3)()(1)(0)0xxxxxxxxxxxxpxxxxxxrrrparrrparprpr221()1exp()ARpkkpfajk5.2.2莱文森-德宾（Levinson-Durbin）算法此算法是以递推求解的方法。221,112,12,22{,},{,,}aaa21,11,21,11{,,,,}pppppaaa设p-1阶AR参量矩阵方程为：211,11,11(0)(1)(1)(1)(0)(2)0(1)(2)(0)0xxxxxxppxxxxxxppxxxxxxrrrparrrparprpr(5)11(1)111,11,121(6)[1,,,][,00]ppppTppppTpRAEAaaE这里：，，即：(1)(0)(1)(1)(1)(0)(2)(1)(2)(0)1,2,,1.xxxxxxxxxxxxpxxxxxxrrrprrrpRrprprkp1,121,10222211,11,10(1)/(0)(1)(1)(0)(0)(1)(0)(1)xxxxxxxxxxxxxxarrrrarrara解得：211,121,111,111(0)(1)(1)(0)0(0)(1)(1)(0)0xxxxxxxxxxxxxxxxkrrarrrarrar当时，上式(5)成为：即：递推求解p阶AR参数：取(7)式的共轭，并利用的关系，得到：2111,1(5)0()()(7)00ppxxpixxinrnarnin将式写成2,1,2,,,,pppppaaa和*()()xxxxrnrn21*11,10()()01,2,,1ppxxpixxinrnarinnpnn排序倒置于是有：*1,1*1,2*21,110(0)(1)(1)0(1)(0)(2)(8)(1)(2)(0)1ppxxxxxxppxxxxxxppxxxxxxarrrparrrparprpr***11,11,21,1211111[,,,,1][0,0,,](9)TTppppppTTpppTTppAaaaERAE（）设：，则有：而：如果令：11()(0)(1)()(1)(10)(1)()(1)(0)()xxxxxxxxpxxppxxxxxxxxxxrprrrpRrRrRrprrrp（）（）（）*,11,11,111*,11,11,1,1100(11)001pppppTpppppppppppaaaAAKKAaaaa那么有：111111111,1()(1)0()(1)(0)(0)(1)()0(1)()()()()xxppppxxxxxxxxxxxxxxxxTpppxxppxxpppxxpkxxkrpRARArrprrrrrprKARrpRArpaKrparpk（）（）（）（）1*,111()pkxxkpTprpkRA（）1*11,111,1121*11,1121,11()()()()()()0000()()ppxxpkxxpkpTxxpkxxpkppxxpkxxkppxxpkxxpkErparpkKrparpkErparpkKrparpk2()0(12)0ppE1211,11*2*11,111,121[()()]0(13)()()(14)[()()](15)pppxxpkxxkpPpxxpkxxkpxxpkxxkppKrparpkKrparpkrparpkK由同样有：由式（13）得：,,1,1,22*2221111~117(1)18ppppipiPppipppppppaKaaKaipKKK*（）由(11)式,（）由(12),(14)式，可得：（）下面，进一步对反射系数和预测误差功率进行分析：1）反射系数为p阶AR模型的第p个系数将H(z)用极点表示为：,pppKa1()1piiiGHzaz11112()(1)(1)...(1)pGHzzzzzzz对于因果稳定系统，极点全在单位圆内，因此有：应等于p个绝对值小于1的zi相乘，其绝对值必小于1，即有：因此，若由于有限字长或相关函数估计误差影响递推过程中出现||1,即|Kp|1,则模型不稳定，递推应立即停止。||1,1,2,...,;izip,||||11,2,...pppKap,ppa,ppa2）预测功率因为所以可见AR模型阶次越高，预测误差功率越小。如果，则有，表明在p阶时，预测误差功率已达到最小，说明AR过程的正确模型应是p阶，也可以此来确定AR模型的适合阶次。综上所述，一个p阶AR过程，可以等效用三组参数来表示：（1）p+1个自相关函数（2）p+1个AR模型系数和G（或，）（3）和p个反射系数这三组系数可相互导出2p222211||0||1ppppKK且110......pp10pK221pp()xn(0)(1)...()xxxRRRp,1,2,,,...ppppaaa2wp22wpG12,...pKKK(0)xR莱文森-德宾(L-D)算法步骤：1、初始化2、递推公式201,111,1(0){()},(1)/(0),xxxxxxrExnarrKa1,,1,1,......1ppipipipaaKaip,pppaK2221(1)pppK1,121()()pxxpkxxkpprparpkK5.2.3BURG算法此算法在1967年被提出，不同于前节Levinson—Durbin算法基于信号的自相关函数进行递推估计，而是直接从已知信号数据中直接进行递推估计，其估计准则为保证前后预测平均功率最小。p阶前后向预测如下图示：p阶前后向预测误差为：,1*,1()()()(5.3.1)()()()(5.3.2)pfppiipbppiienxnaxnienxnpaxnpi