作者：桂。

時間：2017-04-02  08:08:31

宣告：歡迎被轉載，不過記得註明出處哦~

 【讀書筆記08】

前言

西蒙.赫金的《自適應濾波器原理》第四版第五、六章：最小均方自適應濾波器（LMS，Least Mean Square）以及歸一化最小均方自適應濾波器（NLMS，Normalized Least Mean Square）。全文包括：

　　1）LMS與維納濾波器（Wiener Filter）的區別；

　　2）LMS原理及推導；

　　3）NLMS推導；

　　4）應用例項；

內容為自己的讀書記錄，其中錯誤之處，還請各位幫忙指出！

一、LMS與維納濾波器（Wiener Filter）的區別

• 這裡介紹的LMS/NLMS，通常逐點處理，對應思路是：隨機梯度下降；
• 對於Wiener Filter，給定準則函式J，隨機/批量梯度都可以得出最優解；
• LMS雖然基於梯度下降，但準則僅僅是統計意義且通常引入誤差，可以定義為$J_0$，簡而言之$J$通常不等於$J_0$，得出的最優解$w_o$自然也通常不等於維納最優解；
• 分析LMS通常會分析穩定性，穩定性是基於Wiener解，之前已給出分析。但LMS是Wiener解的近似，所以：迭代步長的穩定性，嚴格適用於Wiener解，對於LMS只是一種近似參考，並沒有充分的理論依據。

下文的分析仍然隨機梯度下降的思路進行。

二、LMS原理及推導

LMS是時間換空間的應用，如果迭代步長過大，仍然有不收斂的問題；如果迭代步長過小，對於不平穩訊號，還沒有實現尋優就又引入了新的誤差，屋漏偏逢連夜雨！所以LMS系統是脆弱的，訊號儘量平穩、哪怕短時平穩也湊合呢。

給出框圖：

關於隨機梯度下降，可以參考之前的文章。這裡直接給出定義式：

利用梯度下降：

$- \nabla J = {\bf{x}}{\left( {{{\bf{w}}^T}{\bf{x}} - {d}} \right)^T}$

給出LMS演算法步驟：

1）給定$\bf{w}(0)$，且$1<\mu<1/\lambda_{max}$；

2）計算輸出值：$y\left( k \right) = {\bf{w}}{\left( k \right)^T}{\bf{x}}\left( k \right)$;

3）計算估計誤差：$e\left( k \right) = d\left( k \right) - y\left( k \right)$;

4）權重更新：${\bf{w}}\left( {k + 1} \right) = {\bf{w}}\left( k \right) + \mu e\left( k \right){\bf{x}}\left( k \right)$

三、NLMS推導

 看到Normalized，與之聯絡的通常是約束條件，看到約束不免想起拉格朗日乘子。思路有了，現在開始分析：

假設${\bf{w}}\left( k \right) \Rightarrow {\bf{w}}\left( {k + 1} \right)$得到最優權重，即：

$d\left( k \right) = {\bf{w}}\left( {k + 1} \right){\bf{x}}\left( k \right)$

我們希望在得到期望權重的附近，迭代不要過大以免錯過最優值：

寫出準則函式：

這裡僅僅分析基於歐式距離$p = 2$的情形，其它範數類似。求解得出:

通常為了防止分母為零迭代方程需要修正，而修正後步長存在偏差，故新增調節因子$\mu$：

給出NLMS演算法步驟：

1）給定$\bf{w}(0)$；

2）計算輸出值：$y\left( k \right) = {\bf{w}}{\left( k \right)^T}{\bf{x}}\left( k \right)$;

3）計算估計誤差：$e\left( k \right) = d\left( k \right) - y\left( k \right)$;

4）權重更新：${\bf{w}}\left( {k + 1} \right) = {\bf{w}}\left( k \right) + \frac{\mu }{{\alpha  + {{\left| {{\bf{x}}\left( k \right)} \right|}^2}}}{\bf{x}}\left( k \right){e^*}\left( k \right)$

四、應用例項

　　A-自適應噪聲濾波

這個場景可以簡化為：一個房間兩個麥克風，一個放在遠處採集房間噪聲，一個放在說話人附近採集帶噪語音，認為兩個音訊檔案的噪聲相似。

這裡噪聲直接用白噪聲，對應實際場景可以認為是採集的噪聲資料，給出主要程式碼：

[s, fs, bits] = wavread(filename);          
s=s-mean(s);                           
s=s/max(abs(s));                       
N=length(s);                            
time=(0:N-1)/fs;                       
%%生成帶噪訊號
clean=s';
ref_noise=0.1*randn(1,length(s));
mixed = clean+ref_noise
%NLMS
mu=0.05;M=32;espon=1e-4;
% [en,wn,yn]=lmsFunc(mu,M,ref_noise,mixed);%
[en,wn,yn]=nlmsFunc(mu,M,ref_noise,mixed,espon);

LMS程式碼：

function [e,w,ee]=lmsFunc(mu,M,u,d)
% Normalized LMS
% Call:
% [e,w]=nlms(mu,M,u,d,a);
%
% Input arguments:
% mu = step size, dim 1x1
% M = filter length, dim 1x1
% u = input signal, dim Nx1
% d = desired signal, dim Nx1
% a = constant, dim 1x1
%
% Output arguments:
% e = estimation error, dim Nx1
% w = final filter coefficients, dim Mx1
%intial value 0

w=zeros(M,1); %This is a vertical column

%input signal length
N=length(u);
%make sure that u and d are colon vectors
u=u(:);
d=d(:);
%NLMS
ee=zeros(1,N);
for n=M:N %Start at M (Filter Length) and Loop to N (Length of Sample)
    uvec=u(n:-1:n-M+1); %Array, start at n, decrement to n-m+1
    e(n)=d(n)-w'*uvec;
    w=w+2*mu*uvec*e(n);
    % y(n) = w'*uvec; %In ALE, this will be the narrowband noise.
end

NLMS程式碼：

function [e,w,ee]=nlmsFunc(mu,M,u,d,a)
% Normalized LMS
% Call:
% [e,w]=nlms(mu,M,u,d,a);
%
% Input arguments:
% mu = step size, dim 1x1
% M = filter length, dim 1x1
% u = input signal, dim Nx1
% d = desired signal, dim Nx1
% a = constant, dim 1x1
%
% Output arguments:
% e = estimation error, dim Nx1
% w = final filter coefficients, dim Mx1
%intial value 0

w=zeros(M,1); %This is a vertical column

%input signal length
N=length(u);
%make sure that u and d are colon vectors
u=u(:);
d=d(:);
%NLMS
ee=zeros(1,N);
for n=M:N %Start at M (Filter Length) and Loop to N (Length of Sample)
    uvec=u(n:-1:n-M+1); %Array, start at n, decrement to n-m+1
    e(n)=d(n)-w'*uvec;
    w=w+mu/(a+uvec'*uvec)*uvec*e(n);
    % y(n) = w'*uvec; %In ALE, this will be the narrowband noise.
end

對應結果圖：

可以看出LMS/NLMS在最開始都有一個自適應的過程。

NLMS基於訊號$x$的能量實現變步長，訊號大步長小，訊號小則步長大：目標訊號明顯，則迭代細緻，不明顯，則一帶而過，呵呵，跟平時看書還挺像，聰明的孩子。

再來看一組訊號：

這裡在中間令噪聲突變，可以看到濾波器又需要重新自適應，因此對於短時平穩LMS勉強使用，如果不斷變呢？非平穩LMS自然無效了，這個時候就需要Kalman Filter來搭把手。

　　B-工頻噪聲濾波

現在有一個音訊訊號，分析頻譜：

可以看到訊號帶有明顯的$50Hz$噪聲，我們知道$50Hz$的正弦與餘弦可以組合成任意相位的$50Hz$頻率訊號，基於這個思路，進行自適應濾波：

給出主要的程式碼：

x1=cos(2*pi*50*time);               
x2=sin(2*pi*50*time);
w1=0.1;                               
w2=0.1;
e=zeros(1, N);                         
y=zeros(1, N);
mu=0.05;                            
for i=1: N                           
  y(i)=w1 * x1(i)+ w2 * x2(i);          
  e(i) =x(i)-y(i);                      
  w1=w1+mu * e(i) * x1(i);             
  w2=w2+mu * e(i) * x2(i);
end

結果圖可以看出，工頻50Hz濾除：

基於LMS的應用還有很多，不一一說啦。

參考：

• Simon Haykin 《Adaptive Filter Theory Fourth Edition》.