假定u(n),x(n)是平稳随机信号,u(n)是白噪声,方差是,现在希望能够建立AR模型的参数和x(n)的自相关函数的关系。关于现代谱估计的参数模型,请查看博客:http://blog.csdn.net/u010866505/article/details/78226020。
----------(8)
----------(9)
----------(10)
将上述方程(8),两边同时乘以x(n+m),并且求均值,有:
----------(11)
于是有:
----------(12)
于是有:
-----------(13)
由于u(n)是方差为的白噪声。根据博客http://blog.csdn.net/u010866505/article/details/78226020中的公式(3),有:
= -----------(14)
即有:
------------(15)
有Z变换的定义,,根据(9)式,可以知道当时,h(0)=1.于是:
-----------(16)
在上面推导中,应用了自相关函数的偶对称性,即,于是,上述可写成矩阵形式如下:
-----------(17)
上述(16)(17)就是AR模型的正则方程,又称yule-walker方程。上述矩阵不但是对称的,而且沿着和主对角线平行的任一对角线上的元素都相等,这样的矩阵称为Toeplitz矩阵。
可以看出,一个p阶AR模型有P+1个参数,,只要知道x(n)的前p+1个自相关函数,就可以求解线性方程组,即那p+1个参数,从而可以求出x(n)的功率谱。高斯消元罚求解线性方程组,复杂度是O() ,Levinson-Durbin根据Toeplitz的对称性,给出了一个搞笑的地推算法,需要计算量O().
首先看一下AR模型和线性预测的关系。
设x(n)在n时刻之前的p个数据已知,我们希望利用这p个数据来预测n时刻的值x(n),预测的方法很多,我们用线性预测来实现,是真实值x(n)的预测,那么有:
---------(18)
令真实值和预测值之间的误差是e(n),则有:
--------(19)
因此,总的预测误差功率为:
--------------(20)
为了(20)达到最小,应该是x(n-p)……x(n-1)和预测误差序列e(n)正交<参看《数字信号处理程序》胡光p531>,即:
-------------------(21)
由此可以有:
--------------(22)
根据参看文献[1]中的(10.5.14),有
-----------(23)
(22)和(23)式称为线性预测的wiener-hopf方程。拿(22)和(23)式和AR的正则方程相比较会发现,这俩非常相似,于是乎,假设x(n)是同一个随机信号,如果线性预测器和AR模型的阶数一致,那么就有:
-----------(24)
上面两个说明,一个p阶AR模型的p+1个参数,同样可以用来构建一个p阶的线性预测器。该预测期的最小均方差和AR模型的激励白噪声的能量()相等 .反过来也成立:一AR模型的输出是同阶线性预测器的输出x(n),那么AR模型的系数就是线性预测期的系数,输入的白噪声的能量即方差应该等于。所以,AR模型和线性预测器是等价的,由此可以看出,AR模型是在最小平方意义上对数据的拟合。
如果说x(n)是一个AR(p)过程,是指x(n)是由u(n)激励一个p阶的AR模型所产生的。
现在回归到求解AR正则方程的levison-durbin算法。
定义是p阶AR模型在阶次为m时的第k个参数k=1,2…..m,而m=1,2…..p。为m阶时最小的误差功率(由上面可知)。由(17)可以知道,当m=1时,有:
------(25)
解出上述的解:
----------------(26)
定义初始条件:
--------(27)
那么有:
------(28)
定义第m阶时的第m个系数,即为km,km称为反射系数,根据Toeplitz矩阵的性质,得到levison-durbinde 递推算法:
----------(29)
至此,就可以将所得系数带入到http://blog.csdn.net/u010866505/article/details/78226020中的公式(10)就可以计算功率谱。
参看文献:《数字信号处理程序》胡广书