第五章功率谱估计1-2节

1、1第五章功率谱估计2第一节引言3/1131、信号与系统的分析研究方法n（1）在时域进行（一维信号）n（2）在频域进行（一维信号）n（3）在时频域进行（二维信号）n（4）在空域进行（多维信号）4/113（1）在时域进行n如我们前面所学的维纳滤波、卡尔曼滤波和自适应滤波都属于这种方法。5/113（2）在频域进行n对确定性信号傅里叶变换是在频率域分析研究的理论基础。n对随机性信号其傅里叶变换并不存在，因此研究它的功率谱。6/113（3）在时频域进行（二维信号）n时频分析是一种新的信号分析方法，对于非平稳信号、宽带信号，采用时频分析方法显得特别重要，近年来得到越来越广泛的重视。n由于时频特征分离的固有

2、性质，使局限时间或频谱的一维信号处理扩展到二维的时频平面，信号的时频特征得到了分离，对信号的分析可以具体到信号的特定时间和特定频率，全面反映观测信号的时频联合特征，相应的时频滤波也变得更为容易、更为直观 。n时频分布可分为线性时频表示（如短时傅立叶变换、小波分析等）和非线性时频表示。7/113（4）在空域进行n空域信号处理器是雷达必不可少的最重要分系统。 n阵列信号处理。属于空域信号处理。n空域相关去噪法。是一种新的自适应小波滤噪法。n极化域空域联合谱估计。n分布式移动通信系统。充分利用空域来提高系统容量及频谱利用率。n空时二维自适应信号处理是新一代机载雷达的一项关键技术。 8/1132、随机

3、信号的频域分析n随机信号是无始无终，能量无限的信号，其傅里叶变换不收敛，因而不能确定此类信号的频谱。n随机信号能量无限，其功率未必无限，因而常用功率谱来描述其频率特性。n随机信号自相关函数的傅里叶变换是信号的功率谱密度。9/1133、谱估计定义n谱估计或功率谱估值：根据有限个观测数据，估计平稳随机信号的功率谱。10/1133、谱估计的应用n谱估计的应用遍及雷达、声纳、通信、生物医学、地震勘探等诸多领域。n在雷达信号处理中，由回波信号功率谱密度、谱峰的密度、高度和位置，可以确定运动目标的位置、幅射强度和运动速度等。n在被动式声纳信号处理中，谱峰的位置可提供鱼雷的方向（方位角）。11/113n在生

4、物医学工程中，有关生理电信号的功率谱密度的谱峰可以指示癫痫病的发作周期。n在电子战中，谱分析可用来对目标进行分类，识别等。根据信号、干扰与噪声的功率谱，可设计适当的滤波器，尽量不失真地重现信号，最大限度地抑制干扰与噪声。12/113功率谱估计分为两大类：n（1）经典谱估计（非参数化方法）（线性谱估计）。n（2）现代谱估计（参数化方法）（非线性谱估计）。4、谱估计的分类13/113n实质上仍依赖于传统的傅里叶变换法。n经典的谱估计法又分为两种：n（a）相关图法（又称间接法（BT法）。是由布莱克曼(R.B.Blackman)和图基(J.W.Tukey)于1958年提出的，称为BT法。n（b）直接法

5、（又称周期图periodogram法）。是由舒斯特(Schuster)于1898年提出的。（1）经典谱估计（非参数化方法）14/113（a）间接法（BT法）n直到1965年快速傅里叶变换算法（FFT）问世以前，是最流行的谱估计方法。BT法又称为相关图法对信号序列估计求其自相关函数值对自相关函数的估计进行加权对加权的自相关函数做傅里叶变换获得功率谱估计。15/113（b）直接法（又称周期图（periodogram)法）对观测到的数据样本直接进行傅里叶变换取模的平方，再除以N得到功率谱估计。不用估计自相关函数，且可以用FFT进行计算，在FFT出现以后，周期图法才得到了广泛的应用。16/113（2）

6、现代谱估计n其基本思想是根据已有的观测数据，建立信号所服从的模型，从而在观测不到的区间上，信号的取值服从模型的分布情况，不再认为是零。n主要讨论参数模型（AR、MA、ARMA）法。17/1135、经典谱估计与现代谱估计方法的比较n(1)经典谱估计的优缺点n(2)现代谱估计的优缺点18/113(1)经典谱估计的优缺点n优点：原理简单，便于实现，并可采用FFT等技术而使计算量大为减小。n缺点：它的估计方差大，谱分辨率差（分辨率大约为数据长度的倒数）。使得这种方法难以应用于短数据记录等情况。n原因：一个假设，即除了得到N个数据外，序列的其它值均被认为是零。但实际上序列或其自相关函数在我们未能观测到或

7、未估计的值，实际上并不全是零。19/113(2)现代谱估计的优缺点n优点：现代谱估计质量比经典谱估计的质量有很大提高。n缺点：没有理论能指导选择一个合适的模型，只能根据功率谱的一些先验知识，即重要的谱特性，选择模型。选型合适，可以得到偏差小分辨率高、方差小的谱估计质量。20/1136、功率谱估计的发展n功率谱估计的分辨能力用参量法可以改进。如自回归模型法、最大熵法和最大似然估计等。n在这些方法中，不再认为观测到的N个数据以外的数据全为零，因此克服了经典法的这个缺点，提高了谱估计的分辨率。21/113n然而这些方法在信噪比（SNR）较低时性能并不好，为此，1982年以来，人们陆续提出了多种基于矩

8、阵奇异值分解或特征值分解的改进的谱估计方法，也叫做超分辨方法。22第二节经典谱估计方法23/113一、相关图法-1958BlackmanTukey根据维纳 辛欣定理，年和给出了相关图法的具体实现。(1)( )( ) ( ) ()xxx nmE x n x nm计算随机信号的自相关函数：-(2)( )( )xxmxxmm 若可由的傅里叶变换得到它的功率谱。24/113-(3)( )( )jmwxxxxmPwm e 求功率谱。( )( )xxPwx n是通过的自相关函数间接得到的，所以，也称为间接法。25/113(1)( )( ),0,1,-1( )x nx n nNx nN先估计随机信号的自相关

9、函数。设为实随机序列的一批样本，共有 个值。二、相关图法功率谱估计步骤估 计 的 方 法 有 两 种 ：(a)有 偏 自 相 关 函 数 估 计(b)无 偏 自 相 关 函 数 估 计26/113-101()( ) (),1-()21N mxxnxxmx n x nmmNNmmN它的自相关函数由下式估计：式中的长度为。27/113(2)()( )xxmx n求的离散傅里叶变换，得到的功率谱的估值。-( )( ),-1MjwmxxxxmMPwm eMN()( )xxPwx n是通过的自相关函数间接得到的，所以又称由于为间接法。28/113(3)FFT采用计算上式功率谱。FFT(0 -1)(-)(

10、0 -1)LMML设变换域为，须将求和域，移到。-1-0( )( )LjwmxxxxmPwm e( )0( )01-(- )-(+1)-1xxxxxxmmMmMmL Mm LLMmL 29/113FFT,:经变换 得2-1-0( )( )( )LjkmLxxxxxxmPkFFTmm e0,1,2,-1kL30/113()()()xxxxx nNmm用的个 有 限 值 得 到自 相 关 函 数 的 估 计，它 与 信 号 本 身 自 相 关 函 数有 一 定 程 度 的 差 别 。三、相关图法功率谱估计质量它对相关图法的估计性能有很大影响。31/113(1)()()xxxxmm讨论接近的程度。(

11、 )xxm主要由的估计偏差和估计方差来衡量。对于不同的样本，得出的估计是不同的，且估计值是一个随机变量。通过确定它们的均值、方差等参数，衡量估计方法的好坏。32/113(a)( )xxm先求的均值-10-101( ) ()-1=( )=( )-(N mnN mxxxxxxnEE x n x nmNmmNmmm( )( )xxmx n即自相关函数的估计值的数学期望等于序列的自相关函数的真实值。( )( )xxxxmm因此，是自相关函数的无偏估计。33/113(b)()xxm其次求的方差22var( )( ) -( )xxxxxxmEmEm根据方差定义：2-1-1200( )1( ) () ( )

12、 ()1-xxN mN mnkEmE x n x nm x k x kmNm 式中第一项： （ ）34/113( )x n零均值实正当随态机序列是序列时，它的各高阶矩都可以用其一阶和二阶矩来表示： 1 2 341 2341 324142 3=+E x x x xE x x E x xE x x E x xE x xE x x) () ( ) ()E x n x n m x k x km则：（ ) ()( ) () ( )() () ()() ( )E x n x n m E x k x kmE x n x k E x n m x kmE x n x km E x n m x k（22( )(

13、- )(- )( - ) (2)xxxxxxxxmk nkm nk m n35/1132-1-122200( )1( )( - )(- )( - )-xxN mN mxxxxxxxxnkEmmk nkm nk m nNm 所以：22(1)( ),xxxxEmm将上式代入中，且（得：22-1-12200var)( ) -( )1( - )(- )( - )3-xxxxxxN mN mxxxxxxnkmEmEmk nkm nk m nNm （ ）36/113-1-1200minmax( - )(- )( - )- ,0,- ,-1-(-1), 0, -1-1N mN mxxxxxxnkk nkm

14、nk m nlk n lkln nNmlNmnlk kNmlNm 针对令最大值、最小值为：当取最大，得：当取最大，0()-lknNm即的 情 况 将 出 现次 ，1-1lNm 的 情 况 将 出 现次 。-lNml依 此 类 推 ， 不 同 的 的 情 况 ，出 现 的 次 数 将 为。37/113-122-(-1)var)1(- )( )()( -)-xxN mxxxxxxlN mmNmlllml mNm上式可写成：（-122-(-1)-122-(-1)(1-)( )()( -)-( )()( -)-N mxxxxxxlN mN mxxxxxxlN mmlNllml mNNmNllml mN

15、m38/113(i)1lim var( )0)xxNxxNmNmm当时，上式以趋于零，即故 （满足一致估计的条件。)xxm（为无偏一致估计。39/113(ii)NmN对于固定 ，当越靠近 时，( )xxm的估计方差越大，因而不能得到有用的估计。-10()-1()( ) ()(),-1Nmxxxxxxnmx n x nmmNmmmNNN此时可按下式估计：40/113- ( )( )xxxxNmEmmN其均值：0 ()()xxxxmEmm表 明 ： 当时 ，由 于()xxm相当于真值用三角窗函数加权。()xxxxmm是（的 有 偏 估 计 。41/113(iii)mN 若取有限值，则当时，( )x

16、xm的偏差为零。lim( ) -)0 xxxxNEmm即：（( )xxxxmm因此，是 （的渐近无偏估计。其估计方差为：2-var( )var( )var( )xxxxxxNmmmmN42/113()()xxxxmm是的渐近无偏一致估计，( )xxm且其估计方差小于的估计方差()( )( )xxxxmm所以一般用 有偏估计 渐近无偏估计做为自相关函数的估计，而较少使用无偏一致估计。43/113-10( )1( )( ) ()xxN mxxnmmx n x nmN为了方便，仍用表示，即：N总之：自相关函数的两种估计方法：即估值是渐近无偏或无偏的。当采样数目增大时，其方差均减少，N因此，这两种自相

17、关估计可通过增大 而得到改进。44/113 自相关函数是渐近无偏估计，经过傅里叶变换得其到的功率谱,估计质量?功率谱估计不一定是渐近无偏估计，是有偏估计，是一种不好的估计方法。45/113二、周期图法n1、周期图法定义n2、周期图与BT法的等价关系n3、周期图法的谱估计质量n4、周期图法的改进46/1131、周期图法定义-101()()()Nmxxnmx n x nmN有 偏 自 相 关 函 数 的 定 义 ：10-1( )0nNd n定义矩形窗函数其它-1( )( ) ( ) () ()xxnmd n x n d nm x nmN 则式可写成：47/113( ) ( )( )d n x ny

18、 n令-1()( ) ()()(-)()()(-)xxnxxmy n y nmy mymNmy mym 上式可简化为：即可看成是和的卷积。-1-0()()() ()()(jwmjwmmmNjwmjwNmy my m ed m x m ex m eXe 而的 傅 里 叶 变 换 为 ：）48/113()jwNXe显 然是 周 期 性 的 。()( )jwNXex n即是有限长序列的傅里叶变换。2*11( )()()()jwjwjwxxNNNPwXeXeXeNN功率谱为自相关函数的傅里叶变换，即( )( )x ny n当为实序列时，也为实序列。*(- )()jwNy mXe的傅里叶变换则为。49/

19、113()jwNNXeNIw直接将模的平方除以 求得的功率谱估周期图法，其结果计的方法称为用 （表示。21( )()jwxxNNPwIwXeN即：（FFT由于序列的傅里叶变换可利用计算而提高运算效率，这是周期图法的主要优点。周期图法是用数据的傅里叶变换直接求得，不计算自相关函数，所以又称直接法。50/113( )x n观测数据FFT取模的平方1 N()jwxxP e周期图法计算功率谱框图51/1132、周期图与BT法的等价关系2101()()Njw nNnIwx n eN11*001( )( )NNjwkjwnknx k exn eN11*()001( )( )NNjw knknx k xn

20、eN52/113,mknkmn即，则：BTBTBT证明了有偏自相关函数的法和周期图法的等价关系。即：法和用周期图法估计功率谱是等价的，法估计质量和周期法的估计质量是一样的。11-*-(-1)01( )( ) ()NNmjwmNmNnIwx n x mneN1-*01( )( ) ()Nmxxnmx n x mnN1-(-1)( )( )=( )NjwmNxxxxmNIwm ePw53/1133、周期图法的谱估计质量为了解周期图法的谱估计效果，讨论它的和方差。1（ ）周期图法的估计均值2*11)()()()jwjwjwNNNNE IwEXeE XeXeNN（-1-1-001=( )( )NNjw

21、njwknkEx n ex k eN-1=( ) ( )( ) ( )jwnjwknkEd n x n ed k x k eN -( - )-1=( ) ( ) ( ) ( )jw k nnkEd n d k x n x k eN 54/113- ,+ ,mk nkm n令即：代入上式得：-1)( ) ()( ) ()1( )( )(1)jwmNmnjwmNxxmE Iwd n d nm E x n x nm eNqmm eN （-11-,-1( )( ) ()0NnmmNqmd n d n mNN 对于其它55/113( )Bartlett)Nqm由于为二个矩形函数的卷积，因此是个三角窗函数

22、，称为（窗函数。22-sin112( )( )sin2jwnNnwNQwd n ewNN 其傅里叶变换为56/113(1)式为两个函数相乘的傅里叶变换，等于它们各自傅里叶变换的卷积。-()()()1( )(-)(2)2NNxxNxxE IwQwPwQPwd即：()()()NNxxQwEIwPw由 上 式 可 知 ， 除 非为函 数 ，否 则将 不 等 于，故 周 期 图 作 为 功 率 谱 的 估 计 是 有 偏 的 。57/113lim( )lim 1-1NNNNmqmN 当时，( )( )lim( )( )NNxxNQwwE IwPw故此时N 因此，周期图作为功率谱的估计，当时是无偏的，即

23、渐近无偏。58/1132（ ）周期图法估计的方差2( )(0-1)xx nnN假定序列是零均值、方差为的实高斯白噪声序列。22var()() -()NNNIwEIwEIw按方差的定义有：59/1132( )()( )xxxx nmw 假设是白色的，当代入：-2-( )()()=jwmNNxxxmE Iwqmm e 得：2(a)()NEIw求-()()()jwmNNxxmE Iwqmm e 根据式：24( )NxEIw60/1132(b)( )NE Iw求12( )NIwww先求在两个频率和处的协方差。1212-1-1-1-1-( - )(- )20000()()1( ) ( ) () ( )N

24、NNNNNj wk lwm nklmnE Iw IwE x k x l x m x neN12-1-1-1-120000-( - )(- )1(- )(-)(-)( -)(-)( -)NNNNxxxxklmnj wk lwm nxxxxxxxxEklmnNkmlnknlme 61/1132( )()()xxxx nmw 假设是白色的，故2( ) ( ) () ( ) =0 xklmnkmlnE x k x l x m x nknlm当且或且或且其它62/11322121241212()(-)sinsin221()(-)sinsin22xN wwN wwwwwwNN将上式代入：1212124-1

25、-1-1-1-(- )()-(- )(-)220000()()NNNNNNjm nwwjm nw wxmnmnE Iw IwNeeN63/11312www当时 ， 上 式 为 ：224sin()11sinNxNwEIwNw 24sin2+sinxNwNw64/113(c)()NIw求 周 期 图的 方 差22var()()-()NNNIwEIwEIw244sin2-sinxxNwNw24sin1sinxNwNw65/1134var()0NxNIw 可见，当时，4var()NxNIw无论取何值，与是同一数量级。结论：周期图法不是功率谱的一致估计，它不是功率谱最好的估计。66/11322sinva

26、r()() 1sinNxxNwIwPwNw对于非白噪声，下列关系近似成立( )x n上面我们假设了是白噪声的情况。2var()()NxxNIwPw 当时，指出：本结论是以假设高斯概率密度为根据的，但其在一个相当宽的范围内定性结果也成立。67/113例子2( )0=143264128256256FFT(a) (b) (c) (d)xx nNNNN信号是均值为 ，方差的白噪声，分别观测了 次。，得到了 ， ， 四种不同长度的白噪声。用点作它们的周期图，表示如下图、 、 、 。68/113图(a)69/113程序nN=64;Nfft=64;nsigma2= 1;nwn=sigma2.*randn(1

27、,N);nW = fft(wn,Nfft);nPxx=(abs(W).2)/Nfft;nf=(0:length(Pxx)-1)/length(Pxx)*2*pi;nplot(f(1:Nfft/2),Pxx(1:Nfft/2)70/113图(b)71/113图(c)72/113图(d)73/113结论1信号的真正功率谱是常数 ，但采用周期图法看出它的上下摆动幅度很大，NN当 增加时，摆动的幅度并不减少，只是随 的加大，摆动的频率加快。周期图法估计功率谱误差很大，功率谱不可靠，它的分辨率低（因取只取一段观测数据）。74/1133、周期图法的改进n周期图法作为功率谱的估计不满足一致估计的条件，必须进

28、行改进。n改进的措施主要是将周期图进行平滑（平均是一种主要的平滑方法），使估计方差减小，从而得到一致谱估计。n平滑的方法主要有两种：（1）窗口处理法。目的平滑周期图（2）周期图平均法（3）Welch法：修正的周期图平均法周期法从原理上讲，是基于相关函数值取傅里叶变换的一种方法。75/113（1）窗口处理法()jwW e选用一适应的功率谱窗函数与周期图进行，达到使周期图平卷积滑的目的。FFT它是出现之前常用的窗口处理法。实质是：选择适当的窗口函数作为加权函数进行加权平均。76/113（a）窗口处理法原理-101( )( ) (),-1-N mxxnmx n x nmmNNm已知：m作相关函数估计

29、时，滞后量增大，参与求和的项就会减少，因而平均效果减弱。-1mN如：当时，求和项只有一项。()xxmm即：较大时的那些相关估值不太可靠，估计方差将会增大。77/1130-1mN的全部相关周期图法用了从函数估计值。如果对自相关函数估计值进行加窗，使方差大的自相关函数值加权小，以减小其对谱估计的影响，就有可能提高谱估计的性能。正由于方差大的相关函数估计值造成了谱估计有较大的方差。78/113( )2-1w mM设是长度为的窗函数，-1-(-1)( )( ) ( )MjwmxxxxmMPwm w m e窗口处理法平滑周期图谱估计为：79/113( )(),jww mW e设的傅里叶变换为( )xxP

30、w由于为两个序列相乘的傅里叶变换，( - )1( )( )2j wxxNPwIW ed可表示为其各自傅里叶变换的卷积，即80/113( )w m因此要求窗函数应是一个偶序列，( )xxPww应是 的实、偶、非负函数，( - )0,-j wW ew 且满足条件：(Hanning)(Hamming)三角窗满足此条件，海宁与汉明窗函数不满足此条件。81/113（b)几个谱窗1,-1( )0mM MNw mn（1）矩形滞后窗82/1131-,-1( )0mmM MNw mMn（2）巴特利特窗（BartlettsWindow)83/11311+cos,-1( )20mmM MNw mMn（3）布莱克曼图

31、基窗（Blackman-Tukey Window)84/1133231-662( )2 1-20mmMmNNw mMmmMN其它n（4）帕曾窗（Parzen Window)85/113结论jw( )W(e)w m由于窗函数比原先周期图宽度窄，比原先主瓣宽度宽，利用函数法可进一步平滑周期图，减少估计误差。但是偏移加大了，使分辨率降低。86/113（2）周期图平均法n原理：对一随机变量进行观测，得到L组独立记录数据。用每一组数据求其平均值，然后将L个均值加起来求平均。或是对一随机变量进行观测，得到一组独立记录数据，先将数据分L段，对各段数据分别求其周期图，再将各段周期图进行平均作为功率谱的估计。n

32、又称为Bartlett方法，由于其采用FFT运算可提高运算速度，是当前用得最多的一种平滑方法。87/113L设：对于一个随机变量进行观测，得到 组独立记录数据，L用每一组数据求平均值，然后将 个均值起来求平均，1/ L结论：得到均值相同，方差是用一组数据得到的方差的。88/113122,Lx xx如果是不相关的随机变量，且每一个有均值 及方差，a（ ）定理122()/LxxxxLL则它们的数学平均的均值等于 ，方差为。89/113b（ ）结论1/,L是每个随机变量方差的L具有 个独立同分布随机变量平均的方差 var0Lx 当时，可达到一致估计的目的。因而将若干个独立估计值进行平均可有效降低估计

33、量的方差。90/113根据定理，为降低周期图方法的估计方差，c（ ）应用Bartlett提出了平均周期图方法，( )(0-1)x nnN将序列分段求周期图然后再平均。91/113( )x nLM设将分成 段，每段有个样本，,( ),1,2,3,( )(-) 0-1,1iiNLMix n iLx nx niM MnMiL 则第 段序列的样本可写成：2-1-01( )( )MjwniiniI wx n eM第 段的周期图为：92/11311( )( )LxxiiP wI wL即：,( )xxmMm如果很小，( )iI w则可假定各段的周期图是互相独立的，L总的谱估计可定义为 段周期图的平均，93/

34、113为了分析偏移，对上式求统计平均，( - )-1()(2jj wxxW eP ed11:( )( ) =( )LxxiiiE PwE I wE I wL得到均值21sin(2)()( ) =sin(2)jMwW eFT w mMw其中：94/113,MMN偏移和每一段的数据个数有关。由于平均周期图的偏移比周期图的偏移大。由于三角谱窗主瓣的宽度变宽，分辨率更加降低，因此，偏移的大小反映分辨率的低与高。表现在：三角谱窗主瓣的宽度比周期图主瓣的宽度宽。平均周期图仍是有偏估计。95/113L由于 各段数据互相独立，1var()var( )jwxxiPeI wL因此平均周期图的方差：1 L平均周期图的估计方差是周期图的方差的。L 当时，估计方差趋于零，达到一致估计的目的。96/113LM随着分段数 的增加，方差愈小，功率谱愈平滑。点数减小，分辨率下降，将使得估计变成有偏。L