1、第二章 现代谱估计,现代谱估计概述 AR模型谱估计 线性预测 Burg算法 ARMA模型谱估计 扩展Prony方法 多重信号分类(MUSIC)法,2.1 现代谱估计概述,经典谱估计的主要问题基于信号参数模型的谱估计方法,2.2 AR模型谱估计,AR模型的正则方程 Levinson-Durbin算法 AR谱估计的自相关法 AR模型阶次的选择 AR模型谱估计的性质,2.2.1 AR模型的正则方程,假定 、 都是实平稳的随机信号， 为白 噪声，方差为 ， 为服从AR过程的因果信号。 由AR模型的差分方程，有将上式两边同乘以 ，并求均值，得,2.2.1 AR模型的正则方程,(a)式中， 为AR模型的单
2、位取样响应。由z变换 的性质 ，当 时，有 。 将之代入上式，有 (b),2.2.1 AR模型的正则方程,综合式(a)与式(b)，有 在上述推导中，应用了实信号自相关函数的偶对 称性，即 。由上式可得 个方 程，写成矩阵形式为,2.2.1 AR模型的正则方程,上述两式即为AR模型的正则方程，又称Yule- Walker方程。,2.2.1 AR模型的正则方程,需要指出的是，上式中的自相关矩阵为Toeplitz 矩阵；若 是复过程，那么 ，则 其自相关矩阵是Hermitian对称的Toeplitz矩阵。 这类矩阵具有一系列好的性质，利用这些性质， 可以找到快速求解AR模型参数的高效算法。,2.2.
3、2 Levinson-Durbin算法,Levinson-Durbin递推算法是求解Yule-Walker 方程的快速有效算法，这种算法利用了方程组系 数矩阵（自相关矩阵）所具有的一系列好的性， 使运算量大大减少。其推导的方法有多种，这里只介绍一种较为简便 的推导方法。 设已求得 阶Yule-Walker方程,2.2.2 Levinson-Durbin算法,(a),2.2.2 Levinson-Durbin算法,的参数 ，要求解的m阶Yule-Walker方程为,2.2.2 Levinson-Durbin算法,(b),2.2.2 Levinson-Durbin算法,为此，将式(a)的系数矩阵增
4、加一行和增加一列，成为下式：,2.2.2 Levinson-Durbin算法,(c),2.2.2 Levinson-Durbin算法,式中利用前述的系数矩阵的特点，将式(c)的行倒序， 同时列也倒序，得到,2.2.2 Levinson-Durbin算法,(d),2.2.2 Levinson-Durbin算法,将待求解的m阶Yule-Walker方程表示成式(c)和 式(d)的线性组合形式，即(e),2.2.2 Levinson-Durbin算法,或式中， 是待定系数，称为反射系数。式(e)两 边各右乘以m阶系数矩阵，得到(f),2.2.2 Levinson-Durbin算法,由式(f)可求出由
5、式(c)的第一个方程可求出从上面的推导中可归纳出由m-1阶模型参数求m 阶模型参数的计算公式如下：,2.2.2 Levinson-Durbin算法,对于AR(p)模型，递推计算直到p阶为止。,2.2.3 AR谱估计的自相关法,已知N点观测数据 和AR的阶 数p，则AR谱估计可按下述步骤来进行：由已知的 估计令,2.2.3 AR谱估计的自相关法,用 代替L-D递推算法式中的 ，对于 ，重新求解Yule-Walker方程，这时求出的AR模型参数是真实参数的估计值，即 和将这些参数代入下式，得到 的功率谱的估计，即,2.2.3 AR谱估计的自相关法,若在（0，2）内对 进行N点均匀抽样，则得 到离散
6、谱式中， 。,2.2.4 AR模型阶次的选择,FPE准则（最终预测误差准则） 随着m的增加，使 达到最小值时的 。AIC准则（信息论准则）前者表征 将随着m的增加而单调下降，后者 表示计算误差将随着m的增加而增长。,2.2.5 AR模型谱估计的性质,平滑特性,2.2.5 AR模型谱估计的性质,频率分辨率AR谱估计的频率分辨率，要优于经典谱估计 方法。其原因在于求解AR模型参数的过程，实 际上意味着将根据 估计的 按一定准则进行了外推。 AR谱的匹配性质,2.3 线性预测,前向线性预测后向线性预测格形滤波器,2.3.1 前向线性预测,已知n时刻以前的m个信号数据 ，用这m个数据来线性预测 n时刻
7、信号 的值，如图所示，预测值为式中，上标f表示前向预测。,2.3.1 前向线性预测,其预测误差为 (a)称此预测器为m阶前向线性预测器。 令 ，由此解得 将式(a)代入上式，得,2.3.1 前向线性预测,(b) 由最小均方误差的表达式及正交性原理可知(c) 联立式(b)与式(c)得,2.3.1 前向线性预测,(d)前向线性预测的Wiener-Hopf方程 解此方程则得m阶线性预测器的最佳参数 及 。,2.3.1 前向线性预测,上式与AR模型参数的正则方程式极其相似，若 令 ，则有 ， 成立。这说 明，对于同一个p阶的AR随机信号 ，其AR模 型和同阶的最佳线性预测器模型是等价的。所以 有 (f
8、 )即p阶线性预测器的输出是一个白噪声序列。,2.3.1 前向线性预测,结论：对于给定的随机信号 ，若其最佳前向 线性预测器的阶次等于 的AR模型阶次时，其 前向线性预测误差为白噪声序列。所以阶次等于 AR模型阶次的最佳前向预测误差滤波器实际上 是AR模型的逆系统，即白化滤波器。,2.3.2 后向线性预测,与前向线性预测对应的还有后向线性预测器。即 由n时刻以后的p个数据 来预测 ，即 式中，上标b代表后向预测。在实际工作中，总是用同一组数据来同时实现前向和后向预测，这样上式可改写为,2.3.2 后向线性预测,预测误差 但习惯上常将写成，即仿照前向预测器的推导方法，同样可导出下列公 式： 最佳
9、均方误差及,2.3.2 后向线性预测,对于实序列有 及 若 为复数序列，则,2.3.3 格形滤波器,对于m阶的前、后向预测误差，有如下递推公式 成立式中， 称为反射系数， 。且,2.3.3 格形滤波器,2.4 Burg算法,Burg算法的基本概念 Burg算法存在的问题 改进的协方差算法,2.4.1 Burg算法的基本概念,基本思想自相关法进行AR谱估计时，是遵循以下思路进行的： 由观测的信号数据 先估计自相关函数 。根据估计的自相关函数，利用Levinson-Durbin递推算法求解模型参数 、 。由得出的AR模型参数计算信号的功率谱 。,2.4.1 Burg算法的基本概念,1967年提出的
10、Burg算法在一定程度上改善了这种状况。它所遵循的计算思路是： 由观测的信号数据 先估计反射系数 。根据估计出的反射系数，利用Levinson-Durbin算法递推出AR模型参数 、 。由得出的AR模型参数计算信号的功率谱 。,2.4.1 Burg算法的基本概念,由于在计算中避开了估计自相关函数 ，而 直接从输入数据计算AR模型参数所以减小了计算误差，从而改善了 的频率分辨率。 Burg算法的另一特点是：使用前向、后向预测 误差平方和最小的原则来估计 ，而不是象自相 关法那样只按前向预测误差 的方差最小的原 则导出其正则方程。,2.4.1 Burg算法的基本概念,算法推导令 应满足在上式中代入
11、格形滤波器公式，可得,2.4.1 Burg算法的基本概念,估计出后 ，阶次m时的AR模型参数系数仍然 由Levinson-Durbin算法递推求出，即有,2.4.1 Burg算法的基本概念,计算步骤 由初始条件 ，再由式(2-4-2)求出 ； 由 得 时的参数 由 及式(2.3.15)求出 和 ，再由式(2.4.2) 估计 ； 依照式(2.4.3) 式(2.4.5)，求 出时的参数 、 及 ； 重复上述过程，直到 ，求出所有阶次时的AR参数。,2.4.1 Burg算法的基本概念,若定义式(2.4.2)的分母为那么可以证明， 可以由 和 递推计算，即这样，可以有效地提高计算速度。,2.4.2 B
12、urg算法存在的问题,Burg算法由于避免了直接计算自相关函数，所 以估计误差、频率分辨率要高于自相关法。但由 于它仍然用Levinson-Durbin算法来求AR模型参 数，因此，仍存在谱峰分裂与偏移问题。 当使用单频正弦信号（周期为T） 加白噪声 作试验数据时，有如下规律： 信噪比高时容易发生。由于此时谱峰突出，因此 已不再是AR过程，而是退化的ARMA过程。,2.4.2 Burg算法存在的问题,当信号长度 （ 为采样周期）时， 对任何整数k，任何 的初始相位，谱峰无分裂。当信号长度 ， 的初始相位为0或 时，谱峰无分裂；当信号长度 ， 的初始相位为 的奇数倍时，谱峰分裂严重；N充分大时，
13、分裂现象逐渐减弱。,2.4.3 改进的协方差算法,基本思想：由观测的数据 直接估计 ，而不需要估计 ，从而无需使用 Levinson-Durbin算法。具体说来，其思路是首 先定义：然后，令也就是说，在令 为最小时，不是仅令其相对为最小，而是令其相对 都为最小（m由1至p）。,2.5 ARMA模型谱估计,噪声对AR模型谱估计的影响 MA模型谱估计 ARMA模型谱估计,2.5.1 噪声对AR谱估计的影响,设 是一个p阶AR过程，它被测量噪声 污 染后成为 。即如果 是方差为 的 白噪声，且与 不相关，则有,2.5.2 MA谱估计的计算,(2.5.3)(2.5.4)(2.5.5) 由式(2.5.3
14、)得(2.5.6),2.5.2 MA谱估计的计算,将上式两边同乘以 ，并求均值，得(2.5.7) 式中， 。,2.5.2 MA谱估计的计算,将式(1.6.6)代入 ，并注意到 是方差 为 的白噪声，有(2.5.8) 对MA(q)模型，由式(2.5.4) ，得,2.5.2 MA谱估计的计算,所以，可以求出MA(q)模型的正则方程，即有MA(q)的功率谱为 等效于经典谱估计中的自相关法，即MA谱估 计等效为信号长度为q+1的自相关法谱估计。,2.5.3 ARMA模型谱估计,ARMA(p, q)模型的差分方程式中， 。类似地，可导出其正则方程如下：,2.5.3 ARMA模型谱估计,式中， 是系数 和
15、 的函数，前q+1个方程是 高度非线性的。从第q+1个方程开始是线性的， 可以解出AR部分的系数，将上式中的第二个方 程写成如下展开形式：,2.5.3 ARMA模型谱估计,上式虽然可解出AR部分的系数，但存在以下两 个问题： 由于式中的真实自相关函数 是未知的，因此只能使用估计值 来代替，且要用到大延迟的估计值（最大延迟是p+q），而对于给定的信号长度，这将造成 估计很不准确。因而，也就不能得到AR部分系数的准确估计。,2.5.3 ARMA模型谱估计,式中阶次p和q都是未知的，需要事先指定。而p和q的不正确指定有可能导致式(2.5.13)的系数矩阵奇异。因此，在实际应用中，对式(2.5.13)
16、采用更一般 的形式，即取L个方程，这里 ，即式中，,2.5.3 ARMA模型谱估计,由此得到 的最小二乘解为 求得ARMA(p, q)模型中的AR参数，余下的任 务就是求解MA部分的参数。,2.5.3 ARMA模型谱估计,利用求得的AR系数先得到一个FIR系统为 序列 经此FIR系统滤波，得到一个输出序列ARMA(p, q)模型与FIR系统 级联，近似于模 型 。,2.5.3 ARMA模型谱估计,因此，可以利用输出序列 估计自相关序列 并按MA( q)模型谱估计公式来得到MA谱，即得到MA谱估计 后，利用下式即可求得 ARMA谱估计,2.6 扩展Prony方法,扩展Prony方法对于N点复数序
17、列 采用 的估计模型是一组p个具有任意幅值、相位、 频率与衰减因子的指数函数，即式中， 为幅值； 为相位（单位：弧度）；为衰减因子； 表示振荡频率； 代表 采样间隔。,2.6 扩展Prony方法,为叙述方便起见，先将上式写成如下形式(a)式中， ， 。并定义特 征多项式(b) 式中， 。显然，式(a)中的 即为此多项 式的根。,2.6 扩展Prony方法,由（a）式可得两边同乘 ，并求和，得上式等于零是因为第二个求和项恰好是式 (b)位于根 处的特征多项式 。 上式意味着， 满足递推差分方程式,2.6 扩展Prony方法,令 ，则该式表明：此过程是一个特殊的ARMA(p,p)过程。,2.6 扩
18、展Prony方法,不妨令可得这就将 参数估计问题转化为AR模型 参数估计问题了。当找到 之后，利 用特征方程式求出 。然后，利用式(b) 可构成以下方程：,2.6 扩展Prony方法,式中， ；。这里是一 维Vandermonde矩阵。,2.6 扩展Prony方法,使 为最小的解为从而，正弦信号的振幅、相位、衰减因子与频率 可由以下算法得出： 振幅 相位 衰减因子 频率,2.7 多重信号分类(MUSIC)法,相关矩阵的特征分解基于信号子空间的频率估计基于噪声子空间的频率估计改进的MUSIC法,2.7.1 相关矩阵的特征分解,基本思想：将自相关矩阵中的信息空间分解成两个子空间，即信号子空间和噪声
19、子空间。 设式中,2.7.1 相关矩阵的特征分解,可得自相关函数：式中， 是第i个复正弦的功率。 若定义矢量 是由复正弦波 的N个取样值构成的矢量， 则 可写为式中 称为信号矢量,2.7.1 相关矩阵的特征分解,令 是由 的N个取样数据构成的矢量，即是由白噪声 的N个取样值构成的矢量，即 由前述式可得的自相关矩阵为,2.7.1 相关矩阵的特征分解,若再定义信号自相关矩阵及噪声自相关矩阵如下：则 、 都是N阶方阵，其秩分别为M和N。 而 可写为若 ，显然， 是奇异的，但 由于 的存在而正定,2.7.1 相关矩阵的特征分解,现对方阵 做特征分解，得 ：式中， 为特征向量 所对应的特征值，即并且且特
20、征向量之间是正交的，即,2.7.1 相关矩阵的特征分解,由于 的秩为M，故其特征值中必有 个0。 因此， 可写成如下形式 称 为主特征向量。它所张的空 间称为信号子空间，其维数为M； 而由 所张的空间称为噪声 子空间，其维数为 。显然，这两个子空 间互为正交补空间。,2.7.1 相关矩阵的特征分解,对于单位阵I，也可表示为特征向量 的外积， 即因此， 可表示为,2.7.2 基于信号子空间的频率估计,若舍去特征向量 ，仅保留 信号子空间，那么我们将用秩为M的相关阵 来近似 ，这样可大大提高信号 的信噪比。 基于矩阵 ，再用前面所介绍的方法来估计 的功率谱，将得到好的频率估计和功率谱估计。,2.7
21、.3 基于噪声子空间的频率估计,对于 的非零特征值 所对应的特 征矢量 ，有成立。从而可得 因 为标量，所以有,2.7.3 基于噪声子空间的频率估计,即信号子空间的基底 可以表示为另一组正交基的线性组合。 所以 与 同为信 号子空间的不同正交基组。 又因为信号子空间与噪声子空间互为正交补空 间，所以必有 与噪声子空间中的 任一矢量正交，即,2.7.3 基于噪声子空间的频率估计,Pisarenko法 若 ，则噪声子空间只含有一个基矢 量 ，其特征值为 。且有写成分量形式，得令 ，则有,2.7.3 基于噪声子空间的频率估计,在已知 即 时， 在单位圆上的零点即为待估计的频率 ， 这就是Pisare
22、nko频率估计的基本思想。当 时，将产生虚假的“寄生频率 点”，其原因在于 和 的维数都是N维的， 即,2.7.3 基于噪声子空间的频率估计,这样， 在单位圆上将有个零点，可计算 个频率点。显然，有 点频率为虚假的寄生频率。为避免这一问题，可使用MUSIC法。,2.7.4 改进的MUSIC法,Pisarenko谐波分解法的适用性有一个前题， 即观测样本的长度N应等于M+1（M为复正弦信 号的个数）。MUSIC算法只需满足 ，就可完成 对于样本中所有正弦信号的频率估计。,2.7.4 改进的MUSIC法,对于 ，构成式中， ，显然对于 时，恒有 成立。那么,2.7.4 改进的MUSIC法,在 处应为无限大。但由于 是根据估计的自相关函数计算出来的，其计算误差虽使上式在处为有限值，但仍为充分大的有限值。因此，可以依据 的峰值来估计的频率点。,2.7.4 改进的MUSIC法,若令 ，其中 ，则所得估计为MUSIC估计，即若令 ，其中 ，则所得功率谱称特征向量（eigenvector, EV）估计，即,
