《数字信号处理题解及电子》第12章-2课件.ppt
对同样一组数据，我们可以实现双向预测：()x np(1)x n(1)x np( )x npkffknxkanx1)()()(Forward Prediction)()()(nxnxneff前向预测误差序列2|( )| ffEen误差功率1 ( )( ) ()pbbkx na k x nk Backward Prediction2|( )| bbEe n后向预测误差功率1 ()( ) ()pbbkxnpak x npk 对同一组数据的后向预测()()()()( )bbbbe npx npxnpe npe n后向预测误差序列( )()()bbe nx npe np令：0,1,2,( )bbkpa k可以得到使 最小的 及 。当然也可使用正交原理得：b(1),( )bbaapminbmin11(0)( ) ( )( )( ) (),1,2,.,pbbxxkpbxxkrak r krmak r mkmp 后向预测的Wiener-Hopef Eqfbminmin)()()()()(*kakakakakaffbfb可以证明：前、后向预测对等关系上述结果表明，使用已知的 p 个数据，我们可以实现前向预测，也可以实现后向预测，两种情况下可各自得到对等的Wiener-Hopf方程。将它们单独使用，所得分辨率都不理想。可以设想，如将二者结合起来，即同时使前向、后向预测误差功率为最小，应能得到更好的分辨率。人们在线性预测方面进行了大量的研究。11*11( )( )(1)( )(1)( )ffbmmmmbbfmmmmenenk enenenk en前、后向预测误差序列有如下的关系：00( )( )( )fbene nx n1,2,mp初始条件11111111(1),( )|(1)|*|(1)|cov(1),( )var( )var(1)bfmmmfbmmbfmmfbmmenenkenenenenenen 反射系数上述关系引出了线性预测中的Lattice结构。这一结构在现代谱估计、语音信号处理中有着重要的应用。22|( )|( )|ffnbbnene n上述的关系还是集总平均。对实际的信号：单个样本有限长，求均值要简化，对( ),0,1,1x n nN1( )( )( )( )( ) ()ffpfkenx nxnx nak x nk22|( )| |( )| ffbbEenEe n取代 的范围n1( )( )( )( )( ) ()ffpfkenx nxnx nak x nk (0)0 x(0)(0)fex(1)(1) (0)fxax (1)(1)(1)fexx(1)x N (1)x Np ( )x p(0)x(1)x p(1)( ) (1)fx Npap x N (1)(1)( ) (1)ffeNpx Npap x N (0)(0)(1)(0)(1)(1)(2)(0)(1)( )(1)(1)(0)( )(1)(2)()(1)(1)(1)(1) ()( )(1)(2)(1)fffffffxexxex px pxepx px pxxepx Nx Nx Npx NpeNx Nx Npx NpeNx Nx NeNpx 1(1)(2)( )(1)fffaaapNN点数据，前向预测误差序列范围3X2X1X(0)(1)(0)(1)(2)(0)( )(1)(1)(0)(1)(2)()(1)(1)(1) ()(1)(2)(1)xxxx px pxx px pxxx Nx Nx Npx Npx Nx Npx Npx Nx Nx N0:X上三角+中间块+下三角：上、下加窗；(0)( )(1)(1)ffffeepeNeNP 0:X1:X( )(1)ffepeN 中间块：上、下不加窗；(0)( )(1)fffeepeN 2:X中间块+上三角：下不加窗、上加窗；( )(1)(1)fffepeNeNp 3:X中间块+下三角：上不加窗、下加窗；12.6 AR模型系数求解算法AR模型系数求解算法很多，人们目前仍在探讨新的求解算法。目前，常用的算法是： 1. 自相关法 2. Burg算法 3. 协方差(covariance)方法； 4. 改进的协方差算法（modified ) , 又称：Marple 算法 5. 最大似然（Maximum Likelihood)估计 3. 递推算法：由 求 ，由 递推，还是直接由 递推)(nx)( mrx)( mrx)(nx各算法之间的主要区别：1. 的取值范围，即nnef),(10,XX23,XX选择那一个？ 2. 仅用前向预测，还是前后向都预测？即 令 最小，还是 最小？ffbTfpfpfpfpfppNepeeee)1(,),(,),1 (),0(一、自相关法(1),(2),( )ffffTpaaaapfpfpaXe10fpHfppNnfpfeene| )(|102pffpHoaXXmin001令：使用0X使用前向预测使最小，得注意：矩阵 的结果，即是对有限长数据求出的自相关函数，因此，上式等效于：pffppoaRmin100HXXN自相关法的特点：1. 只用前向预测，且 等效前、后加窗, 分辨率不好;)(nef2. 用 ，得到的 是Toeplits阵，才可能用Levinson算法求解；00XXH1pR3. 实际上是我们前面讨论过的Yule-Walker 方 程。方法最简单。12121|( )|1|( )|Nffppn pNbbppn penNpenNp11*1100( )( )(1)( )(1)( )1,2,( )( )( )ffbmmm mbbfmmmmfbenenk enenenk enmpene nx n二、Burg算法使用前、后向预测12fbfb前、后都不加窗Lattice 结构，递推算法1211211*11| ) 1(| )(|) 1()(2NmnbmNmnfmNmnbmfmmnenenenek12*11)|1 (1, 2 , 1,)()()()(mmmmmmmmmkmkkmakmakkaka先求： ( )mmkam令：0fbmk得到 的求解公式：mk 再用Levinson 递推求其它递推步骤1. 令： 求出2. 求 时的参数3. 求出 ，再求4. 用Levinson算法，求 时的5. 重复上述过程，直到);()()(00nxnenebf1k1m1111(1),(1 |) (0)xakkr11( )( )fbene n、2k2m22 (1),apm Burg算法：一个公认的较好的算法。Burg 算法的特点：112211|( )| ,|( )|NNffbbppppn pn penenNpNp1. 同时使用前向后后向预测，即使12fbfb最小2. 的选择保证前、后不加 窗，即( ),( )fbppen en3. 在每一级， 仅对 最小，然后套用自相关法的Levinson递推算法，影响分辨率；fbmk4. 直接用数据递推，方法简单。三、改进的协方差法Marple方法11221211|( )| ,|( )|fbfbNNffbbppppn pn penenNpNp同Burg算法0 ( )fbmai1,2,1,2,immp注意：这是Marple 算法和Burg算法的最大区别。Burg算法仅：/0,1,fbmkmp (1,1)(1,2)(1, )(1,0)(1)(2,1)(2,2)(2, )(2,0)(2)( )( ,1)( ,2)( , )( ,0)xxxxxxxxxxxxcccpcacccpcaa pcpcpcp pcp 上述最小化的结果是得到一个协方差方程：注意：该矩阵不是Toeplitz矩阵，因此不能用Levinson算法求解。Marple于1983年给出的求解上式的快速递归算法。所以，该算法称作“改进的协方差法，或Marple算法。该算法的估计性能最好，但计算复杂。（e）Burg算法 Burg算法10, ( )pf13p（g）Marple算法 Marple算法10, ( )ph13p12.7 MA模型qkknukbnunx1)()()()(qkkzkbzH1)(1)(221()|1( )|qjjkxkP eb k e( )( )H zB z)(nu)(nx(0)1b10( )( ) ()( ) ()() ( )( )()xqkqxukx n x nmb k u nmku nm x nb k rr mEEmk0,1,mq220( ) ()( ) ()( )0,qq mk mkxb k b kmb k b kmr mmq再推导一步，有：非线性方程组MA模型的正则方程()( )qjj mBTxmqPer m e222()()( )jjMAqj mxmqPeB er m e从谱估计的角度，MA模型等效于经典法中的间接法，所以分辨率低。因此，MA模型用于谱估计无优势。但，MA模型： 1. 常用于系统辨识； 2. ARMA模型中包含了MA部分。令其等效为 模型求解算法：由于MA模型的正则方程是非线性方程，所以人们提出了很多的求解算法，如谱分解、基于迭代的方法、基于高阶AR模型近似的方法。后者最好用，基础是Wold分解定理。1( )( )1( )qkqkHzHzb k z )MA(q1)(11)(1)(kkzkazAzH对 建立一个无穷阶的AR模型( )x n( ) ( )1Az B z于是有：( )( )a kb k步骤：1. 由 ,建立 得 ;2. 对 建立 阶线性预测器，系数为 ,即建立两次AR模型。1,.1 , 0),(Nnnxqpp),AR（pkkap1),()(kapqqkkb1),.,(1( )( ) ()( )qka mb k a mkm1,m 1,mp 12( )( )()( )| ( )|qppkMAmamb k amke me m近似12.8 ARMA(p,q)模型2101( ) ()( ) ()( )( ) (),pq mxkkxpxka k r mkh k b mkr ma k r mkmq0,1,kqARMA模型的正则方程pm 对第二个式子，( )(1)(1)(1)(1)(1)( )(2)(2)(2)(1)(2)( )()( )xxxxxxxxxxxxr qr qr qpr qar qr qr qpr qar qpr qpr qr qpa p ,RarrRae 可以先 求 ，然后再解第一个方程，求出 ；但这样做的效果不好，一是 的性能不好，二是第一个方程也不好求解。首先，建立一个超定方程（方程个数未知数）： paa 11qbb( )(1)(1)(1)(1)(2)(1)(2)( )()( )(1)(2)()()xxxxxxxxxxxxr qr qr qpr qaar qpr qpr qr qpa pr Mr Mr Mpr M paa 11( )( ) ()( ),1pxxkr ma k r mke mmq 21| ( )|MHn qe n e e1()HH aR RR r用求伪逆的方法可求出 ；注意，伪逆可用奇异值分解（SVD）的方法求解；求出 后，剩下的工作是求b a apaaa,.,212. 用 对 滤波；3. 滤波输出 相当于一 MA（q) 过程，按 上节MA模型的求解方法，可求 出ARMA(p,q)模型 的 参数。( )A z)(nx)(nyARMA 模型系数求解的方法：1. 先求出： ，它们可构成 ；( )A z( )( )B zA z( )u n( )y n)(nx( )A z( )u n( )y n( )B z（a）MA(10) （b）MA(16) （c）ARMA(10,10) （d）ARMA(10,13)12.10 基于矩阵特征分解的功率谱估计假定信号由 M 个复正弦加白噪声组成：1( )exp( )MkkkkX nAjnju n 21( )exp()( )Mxiiwir kAjkk 21( )2()MxiiwiPA 已知：(0),(1),( )xxxrrrp不会奇异*(0)(1)( )(1)(0)(1)( )(1)(0)xxxxxxpxxxrrr prrr pRr pr pr(1) (1)pp可构成目标：1. 由该矩阵估计 个正弦信号的频率和幅度； 2. 估计信号 的功率谱；( )X nM定义：1,exp(),exp(),1,2,TiiijjpiMe为信号向量，它包含了 个复正弦，其频率和原信号的频率相同。M求解的关键是自相关矩阵的分解：信号相关阵的表示1( )exp()( )MxiiwirkAjkk 因为：1MHpiiiwiARe eI所以：1MHpiiiwiARe eI1MHpiiiiSA e epwWIpppRSW相关矩阵的分解：信号部分和噪声部分秩是秩为M秩为1p再定义11pHpiiiiSV V1,0,HijijV Vij特征分解121:,det()0,0,:pMMpNote ifpM thenSandso1MHpiiiiSV V11pHiiiIV V借用特征向量的特点 主特征向量1MVV构成的p+1维空间11pVV构成的M维信号空间1MVV构成的噪声空间11MpVV111()pMHHiwiiwiiii MVVVV信号空间特征值111pMHHpiiiwiiiiRVVVV基于噪声子空间的频率估计和功率谱估计：111MHHpiiiwMMiifMpthenRVVVV噪声空间只有一个特征向量可以证明：1,0,1,2,MiiM Ve 和信号向量正交1M V110( )exp()01,2,MHiMMike Vvkj kiM即：110( )exp()01,2,MHiMMike Vvkj kiM求解上式，可得到 的 个根，它们都在单位圆上，因此可求出 12,M ( )V z实现了频率估计M10()( )0Mjj kMkV evk e10( )( )0MkMkV zvk zM 阶多项式方法： 由 估计 ，由 构成 ，并假定 ；( ),0,1x n nN ( )xr m(0),( )xxrr ppRMp2. 对 作特征分解，找最小的 ，及pR1p1pV3. 代入上式，解出： 实现了频率估计。 12,M 4. 由下式，求12,MA AA12112212(1)exp()exp()exp()(2)exp( 2)exp( 2)exp( 2)()exp()exp()exp()xMxMxMMrjjjArjjjAr MjMjMjMA求出5. 由1(0)MxiwirAw按上述步骤，可求出正弦信号的参数 Pisarenko 谐波分解若噪声空间向量不止一个，估计信号的频率，可应用谱估计的方法。1211( )( )xpHkkk MPeV 111( )( )MUSICpHHkkk MPeV Ve1. 若1,11kkMpMUSIC(Multiple Signal Classification)方法 111( )1( )EVpHHkkk MKPeV Ve2. 若1/,11kkkMpEV(Eigenvector)方法用特征分解求出的功率谱曲线与本章内容有关的MATLAB文件： 1. pyulear.m 用AR模型的自相关法估 计信号的功率谱，其基本调用格式是： Px, F = pyulear(x, order, Nfft, Fs) 2. pburg.m 用AR模型的Burg算法估计信 号的功率谱，其基本调用格式是： Px, F = pburg(x, order, Nfft, Fs)（一） 有关功率谱估计的MATLAB文件3. pcov.m 用AR模型方差方法估计信号的 功率谱，其基本调用格式是： Px, F = pcov(x, order, Nfft, Fs)4. pmcov.m 用AR模型的改进的方差方法估 计信号的功率谱，其基本调用格式是：Px, F = pmcov(x, order, Nfft, Fs) 5. pmem.m 最大熵功率谱估计，其估计 性能类似pyulear, 其基本调用格式是： Px, F = pmem(x, order, Nfft, Fs)6. pmusic.m 用自相关矩阵分解的MUSIC 算法估计信号的功率谱，其基本调用格 式是： Px, F = pmusic(x, order, Nfft, Fs)7. peig.m 用自相关矩阵分解的特征向量 法估计信号的功率谱，其基本调用格式是：Px, F = peig(x, order, Nfft, Fs), Px, F，V, E = peig(x, order, Nfft, Fs),x :信号向量，order:模型的阶次，Fs:抽样频率，Nfft:对x作FFT时的长度。Px:估计出的功率谱，F是频率轴坐标。对peig, 输出的E 是由自相关矩阵的特征值所组成的向量，V是由特征向量组成的矩阵。V的列向量张成了噪声子空间，V的行数减去列数即是信号子空间的维数。 （二）有关（二）有关AR模型参数估计的文件：模型参数估计的文件：包括：aryule, arburg, arcov 及 armcov。8. aryule.m 用自相关法(即Yule-Walker法)估 计AR模型的参数，其基本调用格式是： a, E = aryule(x, order), a, E，k = aryule(x, order)9. arburg.m 用Burg算法估计AR模型的参数， 其基本调用格式是： a, E = arburg(x, order) a, E，k = arburg(x, order)10. arcov.m 用方差方法估计AR模型的参数， 其基本调用格式是： a, E = arcov(x, order) 11. armcov.m 用改进的方差方法估计AR模型 的参数，其基本调用格式是： a, E = armcov(x, order) x ：信号向量；order：模型的阶次； a：AR模型系数向量； E：AR模型输入白噪声的功率，或order阶线 性预测器的最小预测误差。 k： 反射系数向量。 （三）有关线性预测的MATLAB文件 . lpc ： 用来计算线性预测系数。 a=lpc(x, order); 其作用等同于 aryule； . ac2poly ：由自相关函数求线性预测系数 a, E=ac2poly(R);14. Poly2ac： 由线性预测系数求自相关函数。 Rpoly2ac(a, E); . ac2rc 由自相关函数得到反射系数及。 k, R0=ac2rc(R);16. rc2ac 由反射系数及得到自相关函数。 Rrc2ac(k, R0);17. poly2rc 由线性预测系数得到反射系数 k=poly2rc(a), 或 k, R0=poly2rc(a, E); 17. rc2poly 由反射系数及得到线性预测系数。 a=rc2poly(k)，或 a, E= rc2poly(k, R0); 18. Levinson：用Levinson-Durbin 算法求解 Toeplitz 矩阵，该文件是一个C-MEX 内部 文件。以上多个文件都要调用它。 a, E, k= levinson(R, order)。
