您當前的位置:首頁 > 攝影

自適應濾波器(二)NLMS自適應濾波器

作者:由 張大俠 發表于 攝影時間:2021-03-27

前一篇文章我們講了LMS自適應濾波器,我們先回顧一下LMS演算法流程:

 \begin{array} yy(n)=\boldsymbol{w}^{T}(n) \boldsymbol{x}(n)  \\ \boldsymbol{e}(n)=d(n)-y(n)  \\ \boldsymbol{w}(n+1)=\boldsymbol{w}(n)+2 \mu e(n) \boldsymbol{x}(n) \end{array}

影響LMS效能的因素,也就是最後一個公式的三個因素:

步長

u

,它是由我們事先指定

輸入向量

\mathbf{x}(n)

估計誤差

e(n)

如果

\mathbf{x}(n)

過大,那麼

2ue(n)\mathbf{x}(n)

的結果中,

\mathbf{x}(n)

佔的權重就太大了,而

\mathbf{x}(n)

是帶噪訊號,這樣

梯度噪聲

就被放大了。為了克服這個問題,可使用歸一化LMS濾波器。在迭代時,對輸入向量

\mathbf{x}(n)

歐式範數(就是模值)

的平方進行歸一化(Normalized LMS)。

歸一化LMS濾波器是最小化干擾原理的一種表現形式,這個原理可以表述如下:

從一次迭代到下一次中,自適應濾波器的權向量應當以最小方式改變,而且受到更新的濾波器輸出所施加的約束。

\hat{\mathbf{w}}(n)

表示第n次迭代濾波器的權向量,

\hat{\mathbf{w}}(n+1)

表示第n+1次迭代濾波器的權向量,那麼NLMS設計準則可表述為約束最佳化問題:給定輸入向量

\mathbf{x}(n)

和目標響應

d(n)

,確定更新抽頭向量

\hat{\mathbf{w}}(n+1)

,以使如下增量

 \delta \hat{\mathbf{w}}(n+1)=\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n)

的歐式範數最小化,並受制於以下約束條件

 \hat{\mathbf{w}}^{\mathrm{H}}(n+1) \mathbf{x}(n)=d(n)

使用拉格朗日乘子法來解決這個約束問題,那麼

代價函式

為:

 J(n)=\|\delta \hat{\mathbf{w}}(n+1)\|^{2}+\operatorname{Re}\left[\lambda^{*}\left(d(n)-\hat{\mathbf{w}}^{\mathrm{H}}(n+1) \mathbf{x}(n)\right)\right]

其中,

\lambda

為複數拉格朗日乘子,

*

表示複共軛,

Re[\cdot]

表示取實部運算,約束對代價函式的貢獻是實值的;

\|\delta \hat{\mathbf{w}}(n+1)\|^{2}

表示歐式範數的平方運算,其結果也是實數。因此,代價函式

J(n)

是實值的

二次函式

,且可表示為:

 J(n)=(\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n))^{\mathrm{H}}(\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n))+\operatorname{Re}\left[\lambda^{*}\left(d(n)-\hat{\mathbf{w}}^{\mathrm{H}}(n+1) \mathbf{x}(n)\right)\right]

採用如下步驟尋找最小的最優更新權向量:

代價函式

J(n)

\hat{\mathbf{w}}(n+1)

求導,可得:

 \frac{\partial J(n)}{\partial \hat{\mathbf{w}}^{\mathrm{H}}(n+1)}=2(\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n))-\lambda^{*} \mathbf{x}(n)

令其為零,即得

最優解

為:

 \hat{\mathbf{w}}(n+1)=\hat{\mathbf{w}}(n)+\frac{1}{2} \lambda^{*} \mathbf{x}(n)

將上式帶入約束條件,求解未知

乘子

\lambda

 \begin{aligned} d(n) &=\hat{\mathbf{w}}^{\mathrm{H}}(n+1) \mathbf{x}(n) \\ &=\left(\hat{\mathbf{w}}(n)+\frac{1}{2} \lambda^{*} \mathbf{x}(n)\right)^{\mathrm{H}} \mathbf{x}(n) \\ &=\hat{\mathbf{w}}^{\mathrm{H}}(n) \mathbf{x}(n)+\frac{1}{2} \lambda \mathbf{u}^{\mathrm{H}}(n) \mathbf{x}(n) \\ &=\hat{\mathbf{w}}^{\mathrm{H}}(n) \mathbf{x}(n)+\frac{1}{2} \lambda\|\mathbf{x}(n)\|^{2} \end{aligned}

可求得:

 \lambda=\frac{2 e(n)}{\|\mathbf{u}(n)\|^{2}}

其中,

 e(n)=d(n)-\hat{\mathbf{w}}^{\mathrm{H}}(n) \mathbf{x}(n)

是誤差訊號。

結合前兩步的結果,可得:

 \begin{aligned} \delta \hat{\mathbf{w}}(n+1) &=\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n) \\ &=\frac{1}{\|\mathbf{x}(n)\|^{2}} \mathbf{x}(n) e^{*}(n) \end{aligned}

為了對一次迭代到下一次迭代抽頭權向量的增量變化進行控制而不改變向量的方向,引入一個正的

實數標度因子

u

,該增量可以寫為:

 \begin{aligned} \delta \hat{\mathbf{w}}(n+1) &=\hat{\mathbf{w}}(n+1)-\hat{\mathbf{w}}(n) \\ &=\frac{\widetilde{\mu}}{\|\mathbf{x}(n)\|^{2}} \mathbf{x}(n) e^{*}(n) \end{aligned}

等價的,我們可以寫出:

 \hat{\mathbf{w}}(n+1)=\hat{\mathbf{w}}(n)+\frac{\mu}{\|\mathbf{x}(n)\|^{2}} \mathbf{x}(n) e^{*}(n)

這個公式就是歸一化LMS演算法抽頭權向量的

遞迴公式

,為什麼叫歸一化呢?因為公式中對輸入向量

 \mathbf{x}(n)

歐式範數的平方就行了歸一化。

當輸入向量

\mathbf{x}(n)

較小時,

\frac{\mu}{\|\mathbf{x}(n)\|^{2}}

的值過小,可能導致

數值

計算困難的情況,為了克服這個情況,將上面的表示式改為:

 \hat{\mathbf{w}}(n+1)=\hat{\mathbf{w}}(n)+\frac{\widetilde{\mu}}{\delta+\|\mathbf{x}(n)\|^{2}} \mathbf{x}(n) e^{*}(n)

其中,

\delta>0

我們總結NLMS演算法的步驟如下:

如果有抽頭權向量

\hat{\mathbf{w}}(n)

的先驗知識,則用它來做初值,否則令

\hat{\mathbf{w}}(n)

為0

計算輸出:

y(n)=\boldsymbol{w}^{T}(n) \boldsymbol{x}(n)

計算誤差:

{e}(n)=d(n)-y(n)

權重更新:

\hat{\mathbf{w}}(n+1)=\hat{\mathbf{w}}(n)+\frac{\widetilde{\mu}}{\delta+\|\mathbf{x}(n)\|^{2}} \mathbf{x}(n) e^{*}(n)

自適應濾波器(二)NLMS自適應濾波器

% 輸入引數:

% xn 輸入的訊號,列向量

% dn 所期望的響應

% M 濾波器的階數

% mu 收斂因子(步長)

% 輸出引數:

% W 濾波器係數矩陣

% en 誤差序列

% yn 濾波器輸出

function

[yn, W, en]

=

nlmsFunc

xn, dn, M, mu, delta

itr

=

length

xn

);

en

=

zeros

itr

1

);

W

=

zeros

M

itr

);

% 每一列代表-次迭代,初始為0

% 迭代計算

for

k

=

M

itr

% 第k次迭代

x

=

xn

k

-

1

k

-

M

+

1

);

% 濾波器M個抽頭的輸入

y

=

W

(:,

k

-

1

)。

*

x

% 濾波器的輸出

en

k

=

dn

k

-

y

% 第k次迭代的誤差

% 濾波器權值計算的迭代式

W

(:,

k

=

W

(:,

k

-

1

+

mu

*

en

k

*

x

/

delta

+

x

*

x

);

end

yn

=

inf

*

ones

size

xn

));

% 初值為無窮大是繪圖使用,無窮大處不會繪圖

for

k

=

M

length

xn

x

=

xn

k

-

1

k

-

M

+

1

);

yn

k

=

W

(:,

end

)。

‘*

x

% 最終輸出結果

end

輸入測試訊號如下:

%% 產生測試訊號

fs

=

1

f0

=

0。02

n

=

1000

t

=

0

n

-

1

’/

fs

xs

=

cos

2

*

pi

*

f0

*

t

);

ws

=

awgn

xs

20

‘measured’

);

M

=

20

% 濾波器的階數

xn

=

ws

dn

=

xs

% rho_max = max(eig(ws*ws。‘)); % 輸入訊號相關矩陣的最大特徵值

% mu = (1/rho_max) ; % 收斂因子 0 < mu < 1/rho

mu1

=

0。01

mu2

=

1

delta

=

1e-3

yn

W

en

=

lmsFunc

xn

dn

M

mu1

);

yn2

W2

en2

=

nlmsFunc

xn

dn

M

mu2

delta

);

歡迎關注微信公眾號:Quant_Times

歡迎大家學習我的課程:

標簽: 濾波器  向量  迭代  xn  LMS