CN1790902A - 一种自适应滤波方法及装置 - Google Patents

Abstract

本发明涉及一种自适应滤波方法及装置，包括，将输入信号序列、期望信号分别进行子带变换转化到频域，然后定义一个子带加权递归最小二乘方法的代价函数，得到滤波器系数矢量(频域)的迭代算式，实现系数的自适应调整。本发明的特点是，通过在频域调节子带的加权值，实现对各子带信号在整个方法迭代过程中的份量控制。

Description

一种自适应滤波方法及装置

技术领域

本发明涉及滤波方法及装置，尤其是自适应滤波方法及装置，属于信息技术领域。

背景技术

设计自适应滤波器，主要是在滤波器模块中植入不同的迭代方法而构成。目的是找到一组最优的滤波器系数使残余误差e(n)无限小。所谓递归最小二乘自适应滤波器(RLSAF，RecursiveLeast-Squares Adaptive Filters)方法就是在n-1时刻利用最小二乘方法估计出滤波器抽头系数，在n时刻，利用接收到的新数据信息对滤波器抽头系数进行更新。递归最小二乘自适应滤波器RLS方法由于利用了数据相关矩阵的逆矩阵，白化了输入数据，因此具有比最小均方(LMS，Least-Mean-Square)方法更快的收敛速度，且收敛性能与信号的频谱特性无关，但也因需要对每一个采样时刻的数据进行一次子带变换以及滤波器系数的迭代更新，计算复杂度高，所需存储量大，不利于实时实现，计算量为O(ML²)，O表示量级。

子带加权自适应滤波(WASF，Weighted Adaptive Subband Filtering)方法的主要思想是通过在频域调节子带的加权值，实现对子带信号在整个方法迭代过程中的份量的控制。在低速移动通信的慢时变衰落信道估计中，可以降低或避开功率较强的窄带干扰信号(NBI，Narrow BandInterference)影响，并可较高精度估计出信道的时域响应。另外，加权值的调节过程是自适应的，稳定后的权值分布将与各子带信号信噪比有关，信噪比大的子带信号拥有大的权值。子带加权自适应滤波WASF方法所基于的子带变换的子带数将根据应用的要求主观选定，与滤波器的抽头数无关。同时，任一子带变换，或者空间变换，根据需要，均可以应用到子带加权自适应滤波WASF方法中，也就是说，子带加权自适应滤波WASF方法对不同子带变换的适应能力很强。

发明内容

本发明的目的是提供一种自适应滤波方法及装置，将递归最小二乘自适应滤波器方法的具有较快收敛速度，且收敛性与信号频谱特性无关，即滤波器对输入信号谱分布具有很好的鲁棒性的优点融入子带加权自适应滤波WASF方法中，二者结合在一起，进一步推广。

本发明的实现原理如下：

假设滤波器输入信号为x(n)，输出的期望信号为y(n)，采用横向滤波器结构，逼近估计的滤波器系数矢量 $\hat{H}=[\hat{h}\left(0\right),\hat{h}\left(1\right),L,\hat{h}(L-1)],L$ 为滤波器抽头数目，估计的残余误差为e(n)。应用子带变换，将输入信号x(n)、期望信号y(n)、残余误差e(n)分别向M个不同的子带上映射，即对0≤i≤M-1，k∈Z⁺，有

$X_{i,k}=\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}X\left(n\right){\mathrm{\Ψ}}_{i,k}\left(n\right)---\left(1\right)$

$y_{i,k}=\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}y\left(n\right){\mathrm{\Ψ}}_{i,k}\left(n\right)---\left(2\right)$

$e_{i,k}=\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}e\left(n\right){\mathrm{\Ψ}}_{i,k}\left(n\right)---\left(3\right)$

这里，X(n)＝[x(n)，x(n-1)，L，x(n-L+1)]^T，{Ψ_i，k(n)}代表子带变换的基函数，下标i，k分别表示子带基函数的序号与时间位移，时间位移间隔选择为信号的采样间隔，子带数目M的选择，依据实际应用的要求与方法的复杂度而定。在给出子带加权递归最小二乘自适应滤波RLS-WASF训练方法结构后，各子带上的残余误差表示为

$e_{i,k}=y_{i,k}-{\hat{H}}^H{X}_{i,k}---\left(4\right)$

定义子带加权递归最小二乘自适应滤波RLS-WASF方法的代价函数WLS为

$\mathrm{WLS}=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{m-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,k}{\left|e_{i,k}\right|}^2=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,k}{|y_{i,k}-{\hat{H}}^H{X}_{i,k}|}^2---\left(5\right)$

当代价函数WLS达到最小值，即

$\mathrm{MWLS}=\underset{\hat{H}}{\min }\mathrm{WLS}---\left(6\right)$

时，子带加权递归最小二乘自适应滤波器RLS-WASF方法取得滤波器系数

的最佳估计。为了在各个子带上应用递归最小二乘自适应滤波RLS方法的递归迭代结构，同时又能顾及到在不同子带上配置不同的加权值，定义如下各子带的权值{ω_i，k}

ω_i，k＝ω_i，tη^t-k，0＜η≤1，0≤i≤M-1                  (7)

式中，η为设定的遗忘因子。针对式(5)，计算其梯度

${\▿}_{{\hat{H}}^{*}}=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,l}{\mathrm{\η}}^{t-k}{X}_{i,k}{X}_{i,k}^H\hat{H}-\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-k}{y}_{i,k}^{*}{X}_{i,k}---\left(8\right)$

显然，当式(6)成立时，梯度应当为0，即

$\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-k}{X}_{i,k}{X}_{i,k}^H\hat{H}=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{m-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-k}{y}_{i,k}^{*}{X}_{i,k}---\left(9\right)$

为了简化表达式，令加权相关系数矩阵 $\mathrm{\Φ}\left(t\right)=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-k}{X}_{i,k}{X}_{i,k}^H,$ 加权互相关系数矩阵

$Z\left(t\right)=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-k}{y}_{i,k}^{*}{X}_{i,k},$ 于是(9)式变为

$\mathrm{\Φ}\left(t\right)\hat{H}=Z\left(t\right)---\left(10\right)$

也就是

$\hat{H}={\mathrm{\Φ}}^{-1}\left(t\right)Z\left(t\right)---\left(11\right)$

每一t时刻，直接计算Φ^-1(t)的复杂度非常高，为了简化复杂度，将Φ(t)写为

$\mathrm{\Φ}\left(t\right)=\mathrm{\η}\underset{k=1}{\overset{t-1}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-1-k}{X}_{i,k}{X}_{i,k}^H+\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H$

并假设ω_i，t＝ω_i，t-1，则有

$\mathrm{\Φ}\left(t\right)=\mathrm{\η\Φ}(t-1)+\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H---\left(12\right)$

将上式右侧的第二项，记为A，展开为下式：

$A=\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H=\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\left(\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}X\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\right)\left(\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}X\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\right)}^H$

$=\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}\left[\begin{array}{c}\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\\ \underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x(n-1){\mathrm{\Ψ}}_{i,t}\left(n\right)\\ M\\ \underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x(n-L+1){\mathrm{\Ψ}}_{i,t}\left(n\right)\end{array}\right]{\left[\begin{array}{c}\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\\ \underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x(n-1){\mathrm{\Ψ}}_{i,t}\left(n\right)\\ M\\ \underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}(n-L+1){\mathrm{\Ψ}}_{i,t}\left(n\right)\end{array}\right]}^H$

在实际计算过程中，L＜＜+∞，故有

$\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x(n-a){\mathrm{\Ψ}}_{i,t}\left(n\right)\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x(n-b){\mathrm{\Ψ}}_{i,t}\left(n\right)B{\left(\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\right)}^2$

矩阵A可简记为

$\mathrm{AB}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\left(\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\right)}^2{\left[\begin{array}{cccc}1& 1& L& 1\\ 1& 1& L& 1\\ L& L& L& L\\ 1& 1& L& 1\end{array}\right]}_{L*L}$

记 $U\left(t\right)={\left|\underset{1i=0}{\overset{m-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\left(\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}x\left(n\right){\mathrm{\Ψ}}_{i,t}\left(n\right)\right)}^2\right|}^{1/2}{\left[\begin{array}{cccc}1& 1& L& 1\end{array}\right]}^T,$ 则矩阵A可近似分解为

$A=\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H\mathrm{BU}\left(t\right)U{\left(t\right)}^H---\left(13\right)$

式(12)又可写成

Φ(t)＝ηΦ(t-1)+U(t)U^H(t)                  (14)

根据矩阵理论中的逆矩阵定理，相关系数矩阵的逆矩阵Φ^-1(t)的递归计算公式为

${\mathrm{\Φ}}^{-1}\left(t\right)={\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)-\frac{{\mathrm{\η}}^{-2}{\mathrm{\Φ}}^{-1}(t-1)U\left(t\right)U^H\left(t\right){\mathrm{\Φ}}^{-1}(t-1)}{{1+\mathrm{\η}}^{-1}{U}^H\left(t\right){\mathrm{\Φ}}^{-1}(t-1)U\left(t\right)}---\left(15\right)$

其中，由于U^H(t)Φ^-1(t-1)U(t)是一个数值，故有

U^H(t)Φ^-1(t-1)U(t)＝tr[U^H(t)Φ^-1(t-1)U(t)]

又根据矩阵理论，当U(t)为列向量

tr[U^H(t)Φ^-1(t-1)U(t)]＝tr[Φ^-1(t-1)U(t)U^H(t)]

则有

U^H(t)Φ^-1(t-1)U(t)＝tr[Φ^-1(t-1)U(t)U^H(t)]

其次，根据矩阵理论数值 $\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t},$ 同样有

$\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}=\mathrm{tr}\left[{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{X}_{i,t}\right]$

由前述分析， $\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H\mathrm{BU}\left(t\right)U{\left(t\right)}^H,$ 即为

$\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}=\mathrm{tr}\left[{\mathrm{\Φ}}^{-1}(t-1)U\left(t\right)U^H\left(t\right)\right]$

最终有

$U^H\left(t\right){\mathrm{\Φ}}^{-1}(t-1)U\left(t\right)=\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}---\left(16\right)$

将式(16)带入式(15)得到最终加权相关系数逆矩阵的迭代方程

${\mathrm{\Φ}}^{-1}\left(t\right)={\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\frac{{\mathrm{\η}}^{-2}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}}---\left(17\right)$

再利用互相关系数矩阵的内在联系，将Z(t)写成如下形式

$Z\left(t\right)=\mathrm{\η}\underset{k=1}{\overset{t-1}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{\mathrm{\η}}^{t-1-k}{y}_{i,k}^{*}{X}_{i,k}+\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t}$

并假设ω_i，t＝ω_i，t-1，则有

$Z\left(t\right)=\mathrm{\ηZ}(t-1)+\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t}---\left(18\right)$

将式(18)代入式(11)，至此，得到子带加权递归最小二乘自适应滤波RLS-WASF方法关于横向滤波器系数矢量迭代的关键方程

${\hat{H}}_t=(1-\frac{{\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}}){\hat{H}}_{t-1}+{\mathrm{\Φ}}^{-1}\left(t\right)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t}---\left(19\right)$

其中，I为L*L的单位矩阵，并且假设，Φ^-1(0)＝I。至于各个子带上的权值迭代，简单地按以下方案求出，即

令 $e_{i,k}=y_{i,k}-{\hat{H}}^H{X}_{i,k},$ 则ω_i，k＝1/E|e_i，k|²。考虑到收敛过程中残余误差的可变性，实际操作时，应当用以下公式估计期望值E|e_i，k|²，即

$E{\left|e_{i,k}\right|}^2=f{\left|e_{i,k}\right|}^2+(1-f)E{\left|e_{i,k-1}\right|}^2$

其中，0＜f≤1为调解子带残余误差的遗忘因子，并且当k＜0时，E|e_i，k|²＝0。

综上所述，子带加权递归最小二乘自适应滤波器设计可通过以下步骤得到自适应调整的滤波器系数矢量(频域)

并由滤波器的抽头系数矢量 构成滤波器。本发明自适应滤波方法步骤如下：

第一步，将输入信号序列X(n)＝[x(n)，x(n-1)，L，x(n-L+1)]^T以及期望信号y(n)进行子带变换，分别得到其在完备频域空间中M个不同子空间上的投影X_0，k:X_M-1，k以及y_0，t:y_M-1，t；

第二步，求出t时刻M个子带上的加权值，当前t时刻第k个子带的权值ω_i，t，并用E|e_i，t|²的倒数来估计，即ω_i，t＝1/E|e_i，t|²；

第三步，利用前一个时刻的加权相关系数逆矩阵Φ^-1(t-1)以及第一、二步的结果，计算当前时刻的加权相关系数逆矩阵Φ^-1(t)，迭代方程如下：

${\mathrm{\Φ}}^{-1}\left(t\right)={\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)-\frac{{\mathrm{\η}}^{-2}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}};$

第四步，利用前一个时刻的加权相关系数逆矩阵Φ^-1(t-1)及滤波器系数矢量 等，求得当前t时刻滤波器系数

迭代方程如下：

${\hat{H}}_t=(1-\frac{{\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}}){\hat{H}}_{t-1}+{\mathrm{\Φ}}^{-1}\left(t\right)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t},$

从而得到自适应调整后的滤波器系数矢量(频域)

本发明子带加权递归最小二乘自适应滤波器RLS-WASF方法的良好性能分析如下：假设：(1)遗忘因子η＝1；

(2)不同子带的加权输入信号 $\{\sqrt{{\mathrm{\ω}}_{i,k}}{X}_{i,k},i=0,L,M-1\}$ 独立同分布；

(3)加权输入信号 $\{\sqrt{{\mathrm{\ω}}_{i,k}}{X}_{i,k},k=1,L,t\}$ 是一均值为0，相关矩阵为R的多维高斯分布的随机变量；

(4)式(1)定义的逼近误差e_i，k是一均值为0、方差为σ²的高斯型随机变量。则有：

性能1： 是关于H的无偏估计；

性能2：横向滤波器 的估计方差可以达到克拉美拉下界(Cramer-Rao Lower Bound)；

性能3：子带加权递归最小二乘自适应滤波RLS-WASF方法能够在有限步内收敛。

实际上，传统的自适应滤波器方法，如递归最小二乘RLS方法，是子带加权递归最小二乘自适应滤波RLS-WASF方法的一种特例，因为递归最小二乘RLS方法实际上应用了一种正交变换——单脉冲函数，即

$\mathrm{\δ}\left(n\right)=\left\{\begin{array}{cc}1& n=0\\ 0& n\≠0\end{array}\right.$

而单脉冲函数在时域上保持正交： $\underset{n=0}{\overset{+\∞}{\mathrm{\Σ}}}\mathrm{\δ}\left(n\right)\mathrm{\δ}(n-l)=\mathrm{\δ}\left(l\right)\·$

附图说明

图1是本发明方法的训练结构示意图。

图2是本发明装置在通信系统中一种实施例的原理示意图。其中，201-H为信道冲激响应函数，202-子带加权递归最小二乘自适应滤波器。

图3是本发明一种实施例中信道估计器部分在通信系统中的原理示意图。其中，301-子带加权递归最小二乘自适应滤波器中的信道估计器(RLS-WASF Channel Estimation)。

图4是本发明一种实施例信道估计器模块组成示意图。其中，401-输入信号x(n)子带变换模块，402-期望信号y(n)子带变换模块，403-Φ^-1(t-1)寄存器模块， 寄存器模块，405-权值计算模块，406-求加权相关矩阵逆模块，407-计算滤波器系数模块。

图5是本发明一种实施例信道估计器中模块关系示意图。

具体实施方式

本发明一般用于线性建模的参数估计，在通信系统中可以用来进行信道估计、线性均衡等等。下文以信道估计器为例，基于横向滤波器结构，对本发明方法作进一步详细说明。

请参阅图2-5，本发明的一种实施例，基于子带加权递归最小二乘自适应滤波方法的信道估计器302，主要是将输入信号序列、期望信号分别进行子带变换转化到频域，然后定义一个子带加权递归最小二乘方法的代价函数，推导出滤波器系数矢量(频域)的迭代算式，实现系数的自适应调整。在白适应滤波器202中，RLS-WASF信道估计器302看成是滤波器202的一部分。滤波器202设计时，先用输入信号x(n)和期望信号y(n)进行信道估计，以自适应方法得到子带加权递归最小二乘自适应滤波器202的抽头系数矢量 即信道估计器302的输出——最优的滤波器系数矢量

再由滤波器系数矢量构成横向滤波器202。输入信号x(n)通过横向滤波器202滤波后得到期望信号少y(n)的估计值

也就是，将输入信号序列X(n)＝[x(n)，x(n-1)，L，x(n-L+1)]^T以及期望信号y(n)进行子带变换，经过横向滤波器 $\hat{H}=\left[\begin{array}{cccc}\hat{h}\left(0\right)& \hat{h}\left(1\right)& L& h(n-M+1)\end{array}\right]$ 后，再将其进行子带反变换，最终得到对期望信号y(n)的估计值 。其中，由于滤波器系数矢量

受误差 $e\left(n\right)=y\left(n\right)-\stackrel{o}{y}\left(\mathrm{on}\right)$ 控制，通过对

不断调节，可找到一组最优的滤波器系数矢量

使残余误差e(n)无限小。

下面以估计未知的固定/低速移动通信信道的冲激响应为例，详细叙述基于子带加权递归最小二乘自适应滤波RLS-WASF方法的信道估计器302的实现过程。

本发明植入子带加权递归最小二乘自适应滤波RLS-WASF方法信道估计器，由输入信号x(n)子带变换模块401、期望信号y(n)子带变换模块402、Φ^-1(t-1)寄存器模块403、 寄存器模块404、权值计算模块405、求加权相关矩阵逆模块406和计算滤波器系数矢量模块407组成。

设d(n)为测试信号，η(t)为加性高斯白噪声，H301为信道冲激响应， 为估计出的信道冲激响应的倒数(也就是滤波器系数矢量，用来滤除信道对测试信号的影响)，x(n)为信道估计器输入信号，y(n)为期望信号。测试信号d(n)经过信道，加上高斯白噪声η(t)后得到信道估计器输入信号x(n)，信道估计器期望信号y(n)即为测试信号d(n)，则有表达式如下：

x(n)＝d(n)*h(n)+η(n)

y(n)＝d(n)

设模块已知参数：子带数目M、子带变换基函数的总样点数N、滤波器系数长度L(也即输入信号数据序列长)、遗忘因子f、高斯白噪声η(n)。各模块作用如下：

在输入信号x(n)子带变换模块401中：首先将输入信号x(n)输入移位寄存器，每个信号采样时刻移位一次，得到长为L的数据序列X(n)＝[x(n)，x(n-1)，L，x(n-L+1)]^T，然后将序列按照下式

$X_{i,k}=\underset{n=0}{\overset{N}{\mathrm{\Σ}}}X\left(n\right){\mathrm{\Ψ}}_{i,k}\left(n\right)$

进行子带变换，求得其在完备频域空间中M个不同子空间上的投影X_0，k:X_M-1，k，其中Ψ_i，k(n)代表子带变换的基函数，i表示子带基函数的序号(i＝0，1，…M-1)，k表示时间位移，参数n也为时间位移，N表示子带变换基函数的总样点数(一般，N很大，且N＞＞L)。这里对向量的子带变换就是对向量中的每个元素进行子带变换。

在期望信号y(n)子带变换模块402中，同样也将期望信号y(n)进行子带变换，如下式：

$y_{i,k}=\underset{n=0}{\overset{N}{\mathrm{\Σ}}}y\left(n\right){\mathrm{\Ψ}}_{i,k}\left(n\right)$

得到M个不同子带上的映射y_0，t:y_M-1，t。

在Φ^-1(t-1)寄存器模块403中：Φ^-1(t-1)寄存器用来存放前一时刻的加权相关系数矩阵Φ(t-1)的逆矩阵，寄存器阵列为L*L，按相关矩阵排列。在t＝0时刻，令Φ^-1(0)＝I，I为单位矩阵，即将Φ^-1(t-1)寄存器阵列值按单位阵I初始化。在每个采样时刻之初，模块将寄存器中所存矩阵Φ^-1(t-1)分别送给求加权相关矩阵逆模块406、计算滤波器系数矢量模块407，当计算滤波器系数矢量模块407计算出该时刻滤波器更新系数

后，模块将接收求加权相关矩阵逆模块406的输出矩阵Φ^-1(t)，并用以更新寄存器为当前时刻的加权相关系数矩阵的逆矩阵。

在 寄存器模块404中：

寄存器用来存放前一时刻的滤波器系数矢量 ${\hat{H}}_{t-1}={\left[\begin{array}{ccc}{\hat{H}}_{0,t-1},& L,& {\hat{H}}_{L-1,t-1}\end{array}\right]}^T\·$ 在t＝0时刻，令 ${\hat{H}}_0=0$ 这里0表示零向量。在每个采样时刻之初，该模块将前一时刻滤波器系数 分别送给权值计算模块405、计算滤波器系数矢量模块407，当计算滤波器系数矢量模块407计算出该时刻滤波器更新系数 后，模块接收来自计算滤波器系数矢量模块407的输出 使 寄存器更新为当前时刻的滤波器迭代系数。

在权值计算模块405中：权值计算模块内部除计算单元外还包含一个寄存器，用来存储前一时刻子带残余误差e_i，t-1的模的平方的期望值E|e_i，t-1|²。在每个采样时刻，权值计算模块405接收输入信号x(n)子带变换模块401的输出X_i，t、期望信号y(n)子带变换模块402的输出y_i，t以及 寄存器模块404中寄存器的输出

并由以下步骤计算加权值ω_i，t。这样的处理模块共有M个，对应计算出M个各个子带上的加权值，权值ω_i，t计算步骤如下：

步骤(1)：若t＝0，则初始化内部寄存器E|e_i，t-1|²＝1；

步骤(2)：根据子带残余误差定义求出当前t时刻对应第i各子带上的残余误差e_i，t。子带残余误差定义为期望信号y_i，t与估计期望信号 间的误差，表达式如下：

$e_{i,t}=y_{i,t}-{\hat{H}}_{t-1}{X}_{i,t}$

这里用前一个时刻滤波器系数

对输入信号X_i，t滤波得到当前时刻的期望信号的估计值

步骤(3)：当前时刻子带残余误差模的平方的期望值E|e_i，t|²用当前误差模的平方|e_i，t|²乘上遗忘因子f，加上E|e_i，t-1|²乘上1与遗忘因子差的和来表示，表达式如下：

E|e_i，t|²＝f|e_i，t|²+(1-f)E|e_i，t-1|²

将求得的期望值E|e_i，t|²存入内部寄存器；

步骤(4)：当前t时刻第k个子带的权值ω_i，t用E|e_i，t|²的倒数来估计，即ω_i，t＝1/E|e_i，t|²。

在求加权相关矩阵逆模块406中：每个采样时刻初，求加权相关矩阵逆模块406接收Φ^-1(t-1)寄存器模块403的输出Φ^-1(t-1)、输入信号x(n)子带变换模块401的输出X_0，t:X_M-1，t以及权值计算模块405的输出ω_0，t:ω_M-1，t，并用来计算当前t时刻加权相关系数矩阵的逆矩阵Φ^-1(t)。利用当前时刻输入信号信息在原有逆矩阵的基础上求得的迭代方程如下：

${\mathrm{\Φ}}^{-1}\left(t\right)={\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)-\frac{{\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}}$

在计算滤波器系数矢量模块407中：每个采样时刻初，计算滤波器系数模块接收

寄存器模块404的输出

Φ^-1(t-1)寄存器模块403的输出Φ^-1(t-1)，求加权相关矩阵逆模块406的输出Φ^-1(t)，输入信号x(n)子带变换模块401的输出X_0，t:X_M-1，t，期望信号y(n)子带变换模块402的输出y_0，t:y_M-1，t以及权值计算模块405的输出ω_0，t:ω_M-1，t，并用来计算当前t时刻滤波器系数

迭代方程如下：

${\hat{H}}_t=(I-\frac{{\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{\mathrm{it}}}){\hat{H}}_{t-1}+{\mathrm{\Φ}}^{-1}\left(t\right)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t}\·$

以上模块连接后，信道估计器中信号数据流流程如下：在t时刻信号输入到输入信号x(n)子带变换模块401以及期望信号y(n)子带变换模块402，进行子带变换，变换结果输出到权值计算模块405计算权值。计算出的权值联合Φ^-1(t-1)寄存器模块403以及输入信号x(n)子带变换模块401的输出，输入到求加权相关矩阵逆模块406，进行计算加权相关系数矩阵的逆矩阵。计算出的逆矩阵输出到滤波器系数矢量计算模块407，计算出滤波器系数后输入到

寄存器模块404寄存；同时，求加权相关矩阵逆模块406计算的加权相关系数矩阵逆矩阵，也输入移位寄存器模块403寄存。更新为下一时刻的输入信号数据后，重复以上流程，直至获得稳定的最优的滤波器系数矢量

本发明实施例由滤波器抽头系数矢量

构成的横向滤波器202中，将输入信号序列X(n＝[x(n)，x(n-1)，L，x(n-L+1)]^T以及期望信号y(n)进行子带变换，经过横向滤波器 $\hat{H}=\left[\begin{array}{cccc}\hat{h}\left(0\right)& \hat{h}\left(1\right)& L& h(n-M+1)\end{array}\right]$ 后，再将其进行子带反变换，得到对期望信号y(n)的估计值

。其中，由于滤波器系数矢量 受误差 $e\left(n\right)=y\left(n\right)-\stackrel{o}{y}\left(\mathrm{on}\right)$ 控制，可通过自适应调节，找到一组最优的滤波器系数矢量使残余误差e(n)无限小。

上述的对具体实施例的描述可以使该技术领域的普通技术人员能理解和使用本发明。熟悉本技术的人员显然可以很快地依据这些实施例的揭示做出各种修改，并把在此说明的一般原理应用到其他实施例中而不必经过创造性的劳动。因此，本发明的内容不限于这里所公开的实施例，本领域技术人员根据本发明的揭示，对于本发明做出的改进和修改都应该在本发明的保护范围之内。

Claims (9)

1、一种自适应滤波方法，其特征在于采用以下步骤：

第一步，将输入信号序列以及期望信号进行子带变换，分别得到其在完备频域空间中M个不同子空间上的投影；

第二步，求出t时刻M个子带上的加权值，当前t时刻第k个子带上的权值；

第三步，利用前一个时刻的加权相关系数逆矩阵迭代算式，计算当前时刻的加权相关系数逆矩阵；

第四步，利用前一个时刻的加权相关系数逆矩阵及滤波器系数矢量，求得当前t时刻滤波器系数矢量。

2、根据权利要求1所述的自适应滤波方法，其特征在于：输入信号序列以X(n)＝[x(n)，x(n-1)，L，x(n-L+1)]^T表示，期望信号以y(n)表示，其在完备频域空间中M个不同子空间上的投影分别以X_0，k：X_M-1，k及y_0，t：y_M-1，t表示；当前t时刻第k个子带的权值以ω_i，t表示，并用E|e_i，t|²的倒数来估计，即ω_i，t＝1/E|e_i，t|²；前一个时刻的加权相关系数逆矩阵以Φ^-1(t-1)表示，当前时刻的加权相关系数逆矩阵以Φ^-1(t)表示，则

${\mathrm{\Φ}}^{-1}\left(t\right)={\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)-\frac{{\mathrm{\η}}^{-2}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}};$

前一个时刻的滤波器系数矢量以 表示，当前t时刻滤波器系数以 表示，则；

${\hat{H}}_t=(1-\frac{{\mathrm{\η}}^{-1}{\mathrm{\Φ}}^{-1}(t-1)\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}{X}_{i,t}^H}{1+{\mathrm{\η}}^{-1}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{X}_{i,t}^H{\mathrm{\Φ}}^{-1}(t-1)X_{i,t}}){\hat{H}}_{t-1}+{\mathrm{\Φ}}^{-1}\underset{i-0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}{y}_{i,t}^{*}{X}_{i,t}.$

3、根据权利要求2所述的自适应滤波方法，其特征在于：在得到自适应调整滤波器系数矢量 过程中，各子带上的残余误差为

$e_{i,k}=y_{i,k}-{\hat{H}}^H{X}_{i,k}.$

4、根据权利要求2或3所述的自适应滤波方法，其特征在于：所述方法的代价函数WLS为

$\mathrm{WLS}=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,k}{\left|e_{i,k}\right|}^2=\underset{k=1}{\overset{t}{\mathrm{\Σ}}}\underset{i=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,k}{|y_{i,k}-{\hat{H}}^H{X}_{i,k}|}^2.$

5、根据权利要求4所述的自适应滤波方法，其特征在于：当代价函数WLS达到最小值时，即

$\mathrm{MWLS}=\underset{\hat{H}}{\min }\mathrm{WLS}$ 时

得到频域 的最佳估计值。

6、一种自适应滤波装置，其特征在于：采用权利要求1至5中任一所述的方法进行滤波。

7、根据权利要求6所述的自适应滤波装置，其特征在于：包括：输入信号x(n)子带变换模块、期望信号y(n)子带变换模块、Φ^-1(t-1)寄存器模块、 寄存器模块、权值计算模块、求加权相关矩阵逆模块和计算滤波器系数矢量模块。

8、根据权利要求7所述的自适应滤波装置，其特征在于：该自适应滤波装置中信号数据流的流程如下：在t时刻信号输入到输入信号x(n)子带变换模块401以及期望信号y(n)子带变换模块402，进行子带变换，变换结果输出到权值计算模块405计算权值；计算出的权值联合Φ^-1(t-1)寄存器模块403以及输入信号x(n)子带变换模块401的输出，输入到求加权相关矩阵逆模块406，进行计算加权相关系数矩阵的逆矩阵；计算出的逆矩阵输出到滤波器系数矢量计算模块407，计算出滤波器系数后输入到

寄存器模块404寄存；同时，求加权相关矩阵逆模块406计算的加权相关系数矩阵逆矩阵，也输入移位寄存器模块403寄存；更新为下一时刻的输入信号数据后，重复以上流程，直至获得稳定的最优的滤波器系数矢量

9、根据权利要求7所述的自适应滤波装置，其特征在于：该权值计算模块中，权值ω_i，t的取得步骤如下：

步骤一：若t＝0，则初始化内部寄存器E|e_i，t-1|²＝1；

步骤二：根据子带残余误差定义求出当前t时刻对应第i各子带上的残余误差e_i，t，子带残余误差定义为期望信号y_i，t与估计期望信号

间的误差，表达式如下：

$e_{i,t}=y_{i,t}-{\hat{H}}_{t-1}{X}_{i,t}$

这里用前一个时刻滤波器系数

对输入信号X_i，t滤波得到当前时刻的期望信号的估计值

步骤三：当前时刻子带残余误差模的平方的期望值E|e_i，t|²用当前误差模的平方|e_i，t|²乘上遗忘因子f，加上E|e_i，t-1|²乘上1与遗忘因子差的和来表示，表达式如下：

               E|e_i，t|²＝f|e_i，t|²+(1-f)E|e_i，t-1|²

将求得的期望值E|e_i，t|²存入内部寄存器；

步骤四：当前t时刻第k个子带的权值ω_i，t用E|e_i，t|²的倒数来估计，即ω_i，t＝1/E|e_i，t|²。

Images