【導(dǎo)讀】最小二乘法是一種典型的有效的數(shù)據(jù)處理方法。由著名學(xué)者高斯在1795年提出，他認(rèn)為。觀測(cè)值和計(jì)算值之間的差的平方乘以度量其精度的數(shù)值以后的和為最小。自適應(yīng)算法的目標(biāo)在于，使濾波器輸出與需。要信號(hào)的誤差的平方的統(tǒng)計(jì)平均值最小。這個(gè)準(zhǔn)則根據(jù)輸入數(shù)據(jù)的長(zhǎng)期統(tǒng)計(jì)特性尋求最佳濾波。小二乘法對(duì)不同的數(shù)據(jù)組導(dǎo)出不同的“最佳”濾波器。遞推最小二乘法是最小二乘法的一類快速算法。本文基于matlab進(jìn)行仿真，研究不同輸入和信噪比濾波器的性能。通過對(duì)比實(shí)驗(yàn)，分析收斂性。對(duì)不同類型濾波器輸入信號(hào)獨(dú)立實(shí)驗(yàn)，周期型號(hào)和非周期信號(hào)的濾波器效果都較好，但周期信。號(hào)的濾波效果更好。但從程序的運(yùn)行時(shí)間上看，其復(fù)雜度高，計(jì)算量比較大。

　　

【正文】 自適應(yīng)濾波器仿真 MATLAB 是由美國(guó) mathworks 公司發(fā)布的主要面對(duì)科學(xué)計(jì)算、可視化以及交互式 程序設(shè)計(jì) 的高科技計(jì)算環(huán)境。它將 數(shù)值分析 、 矩陣計(jì)算 、科學(xué)數(shù)據(jù)可視化以及非 線性 動(dòng)態(tài)系統(tǒng)的 建模和仿真等諸多強(qiáng)大功能集成在一個(gè)易于使用的視窗環(huán)境中，為科學(xué)研究、工程設(shè)計(jì)以及必須進(jìn)行有效 數(shù)值計(jì)算 的眾多科學(xué) 領(lǐng)域 提供了一種全面的解決方案，并在很大程度上擺脫了傳統(tǒng)非交互式 程序設(shè)計(jì)語言 （如 C、 Fortran）的編輯模式，代表了當(dāng)今國(guó)際科學(xué)計(jì)算軟件的先進(jìn)水平。 MATLAB 和 Mathematica、 Maple、 MathCAD 并稱為四大 數(shù)學(xué) 軟件。它在 數(shù)學(xué) 類科技應(yīng)用軟件中在 數(shù)值計(jì)算 方面首屈 一指。 MATLAB 可以進(jìn)行 矩陣 運(yùn)算、繪制 函數(shù) 和數(shù)據(jù)、實(shí)現(xiàn) 算法 、創(chuàng)建用戶界面、連 接其他 編程語言 的 程序等，主要應(yīng)用于工程計(jì)算、控制設(shè)計(jì)、 信號(hào)處理 與通訊、 圖像處理 、信號(hào)檢測(cè) 、 金融建模 設(shè)計(jì) 與分析等 領(lǐng)域 。 本節(jié)中我們利用 MATLAB 仿真軟件，仿真系統(tǒng)框圖如下圖所示。 _+)(0 nw )(1 nw )(2 nw )(1 nw N ?… …… …… .)(nx )1( ?nx )2( ?nx )1( ?? Nnx)( ny)( ne )( nd 圖 仿真框圖 正弦波去噪實(shí)驗(yàn) 原始信號(hào)為單一正弦波 )(xrt ，加入方差為 )(tb ，濾波器輸入信號(hào)為自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) )(tx ，輸出信號(hào)為 )(ty 。 以下為相關(guān)程序： clear。 clc。 randn(39。seed39。, 0) 。 rand(39。seed39。, 0) 。 N=1000。 n=200。 k=12。 b=*randn(1,N)。 x=(zeros(1,N))。 for i=1:N xr(1,i)=sin(k*2*pi*i/N)。 x(1,i)=xr(1,i)+b(i)。 end Delta=。 p=(1/Delta)*eye(n)。 K=zeros(N,n)。 w=zeros(N,n)。 e=zeros(1,N)。 y=zeros(1,N)。 for i=n+1:N K(i,:)=(p*x(in+1:i)39。./(1+x(in+1:i)*p*x(in+1:i)39。))39。 e(1,i)=xr(i)w(i1,:)*x(in+1:i)39。 w(i,:)=w(i1,:)+e(1,i)*K(i,:)。 p=pK(i,:)39。*x(in+1:i)*p。 y(1,i)=w(i,:)*x(in+1:i)39。 end figure(1) subplot(2,1,1)。 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) plot(0:Nn,xr(1,n:N),39。r39。) ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 title(39。原始輸入信號(hào) xr時(shí)域圖 39。) grid on subplot(2,1,2)。 plot(abs(fft(xr)))。grid。 ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 title(39。原始輸入信號(hào) xr頻域圖 39。)。 figure(2) subplot(2,1,1)。 plot(0:Nn,x(1,n:N)) title(39。濾波器輸入信號(hào) x時(shí)域圖 39。) ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 grid subplot(2,1,2)。 plot(abs(fft(x)))。grid。 ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 title(39。濾波器輸入信號(hào) x頻域圖 39。)。 figure(3) plot(0:Nn,e(1,n:N))。hold on plot(0:Nn,y(1,n:N),39。r39。)。hold on grid on。 title(39。誤差信號(hào) e(t)和輸出信號(hào) y(t)39。) 。 gtext(39。e(k)39。),gtext(39。y(k)39。) figure(4) subplot(2,1,1)。 plot(abs(fft(e)))。 ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 grid。 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 title(39。誤差信號(hào) e頻域圖 39。)。 subplot(2,1,2)。 plot(abs(fft(y)))。grid。 ylabel(39。幅值 39。)。xlabel(39。時(shí)間 39。)。 title(39。濾波器輸出信號(hào) y頻域圖 39。)。 運(yùn)行以上程序，得到圖 。 圖 原始輸入信號(hào)時(shí)域及頻域圖 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 圖 濾波器輸入信號(hào)時(shí)域及頻域圖 圖 誤差信號(hào)及輸 出信號(hào)時(shí)域圖 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 圖 濾波器輸誤差信號(hào)和輸出信號(hào)頻域圖 圖 為原始輸入信號(hào) )(txr ，為單一頻率正弦波，幅值為 1，圖 為 的隨機(jī)高斯噪聲，得到濾波器輸入信號(hào) )(tx ，從頻譜圖可以看出，包含原始信號(hào)有很多頻率。圖 )(ty 為經(jīng)過 RLS 濾波器后的輸出信號(hào)，藍(lán)色部分 )(te 為誤差信號(hào)。可以看出誤差信號(hào)較小。圖 ，可以看出輸出信號(hào)頻率較集中。由上可知，濾波效果較好。 濾波器正則化參數(shù)的確定 RLS 算法有兩個(gè)特殊的可變參數(shù)：抽頭輸入數(shù)據(jù)的信噪比，這個(gè)量由運(yùn)行的條件決定；正則化參數(shù) ? 。 以下做了對(duì)比實(shí)驗(yàn)，在輸入為正弦波，控制加性高斯噪聲的均值，相同均 值的隨機(jī)的噪聲都保持不變。來實(shí)現(xiàn)不同信噪比，分析正則化參數(shù)的在最佳選擇。 高信噪比 原始輸入信號(hào)是幅值為 1的正弦波，噪聲方差為 。 實(shí)驗(yàn)組 1： ?? 實(shí)驗(yàn)組 2： ?? 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 程序如下： 實(shí)驗(yàn)組 1 randn(39。seed39。, 0) 。 rand(39。seed39。, 0) 。 N=5000。 n=200。 k=12。 b=*randn(1,N)。 x=(zeros(1,N))。 for i=1:N xr(1,i)=sin(k*2*pi*i/N)。 x(1,i)=xr(1,i)+b(i)。 end Delta=。 p=(1/Delta)*eye(n)。 K=zeros(N,n)。 w=zeros(N,n)。 e=zeros(1,N)。 y=zeros(1,N)。 tr=zeros(1,N)。 for i=n+1:N K(i,:)=(p*x(in+1:i)39。./(1+x(in+1:i)*p*x(in+1:i)39。))39。 e(1,i)=xr(i)w(i1,:)*x(in+1:i)39。 w(i,:)=w(i1,:)+e(1,i)*K(i,:)。 p=pK(i,:)39。*x(in+1:i)*p。 y(1,i)=w(i,:)*x(in+1:i)39。 end figure(1) subplot(2,1,1) plot(0:Nn,e(1,n:N))。hold on 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) plot(0:Nn,y(1,n:N),39。r39。)。hold on grid on。 title(39。誤差信號(hào) e(t)和輸出信號(hào) y(t)39。) 。 ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 subplot(2,1,2)。 plot(0:Nn,w(n:N,1)),hold on plot(0:Nn,w(n:N,2),39。r39。),hold off grid title(39。抽頭權(quán)值 w的變化 39。) 。 ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 運(yùn)行以上程序，得到圖 。 運(yùn)行結(jié)果： 圖 誤差信號(hào)輸出信號(hào)和抽頭權(quán)值 w變化圖 實(shí)驗(yàn)組 2程序與實(shí)驗(yàn)組 1大部分相同，改變 Delta= 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 運(yùn)行得到圖 。 圖 誤差信號(hào)輸出信號(hào)和抽頭權(quán)值 w變化圖 當(dāng)信噪比較高時(shí)，實(shí)驗(yàn)組 1 抽頭系數(shù)在 4000 開始收斂，實(shí)驗(yàn)組 2 抽頭系數(shù)在 500 開始收斂。兩組誤差信號(hào)都很小。不難看出，實(shí)驗(yàn)組 2 收斂速度更快。 低信噪比 原始輸入信號(hào)是幅值為 1 的正弦波，噪聲方差為 。 實(shí)驗(yàn)組 1： ?? 實(shí)驗(yàn)組 2： ?? 實(shí)驗(yàn)程序與高信噪比時(shí)大部分一樣，變動(dòng)的語句： b=*randn(1,N)。其它部分 這里就不再重復(fù)。 運(yùn)行得到圖 ，圖 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 圖 誤差信號(hào)輸出信號(hào)和抽頭權(quán)值 w變化圖 圖 誤差信號(hào)輸出信號(hào)和抽頭權(quán)值 w變化圖 當(dāng)信噪比較低時(shí)，實(shí)驗(yàn)組 1 抽頭系數(shù)在 1800 開始收斂，誤差信號(hào)較小。實(shí)驗(yàn)組 2 抽頭系數(shù)在 500 開始收斂，但是誤差信號(hào)波動(dòng)較大。綜合考慮。實(shí)驗(yàn)組 1 的全局性能較好。 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) 結(jié)論 綜合上述兩組對(duì)比實(shí)驗(yàn)，當(dāng)抽頭輸入噪聲電平低（即輸入信噪比較高）時(shí)， RLS 算法的正則 化系數(shù)取較小的正常數(shù)收斂速度更快，濾波效果也較好。當(dāng)抽頭輸入噪聲電平高（即輸入信噪比較低）時(shí)， RLS 算法的正則化系數(shù)取較大的正常數(shù)收斂速度較快，濾波效果也較好，全局性能更好。 輸入信號(hào)不同對(duì)濾波效果的影響 輸入信號(hào)為周期信號(hào) 原始信號(hào)輸入為周期三角波，噪聲方差為 2。 clear。 clc。 randn(39。seed39。, 0) 。 rand(39。seed39。, 0) 。 N=1000。 n=200。 k=12。 b=2*randn(1,N)。 for i=1:N xr(1,i)=sin(k*2*pi*i/N)。 x(1,i)=xr(1,i)+b(i)。 end t=0::。 xr=10*sawtooth(2*pi*5*t,)。 x=(zeros(1,N))。 for j=1:N x(1,j)=xr(1,j)+b(j)。 end Delta=。 p=(1/Delta)*eye(n)。 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) K=zeros(N,n)。 w=zeros(N,n)。 e=zeros(1,N)。 y=zeros(1,N)。 for i=n+1:N K(i,:)=(p*x(in+1:i)39。./(1+x(in+1:i)*p*x(in+1:i)39。))39。 e(1,i)=xr(i)w(i1,:)*x(in+1:i)39。 w(i,:)=w(i1,:)+e(1,i)*K(i,:)。 p=pK(i,:)39。*x(in+1:i)*p。 y(1,i)=w(i,:)*x(in+1:i)39。 end figure(1) subplot(2,1,1)。 plot(0:Nn,xr(1,n:N),39。r39。) ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 title(39。原始輸入信號(hào) xr 時(shí)域圖 39。) grid on subplot(2,1,2)。 plot(abs(fft(xr)))。grid。 ylabel(39。幅值 39。)。 xlabel(39。頻率 39。)。 title(39。原始輸入信號(hào) xr 頻域圖 39。)。 figure(2) subplot(2,1,1)。 plot(0:Nn,x(1,n:N)) title(39。濾波器輸入信號(hào) x時(shí)域圖 39。) ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 grid 自適應(yīng)濾波器的設(shè)計(jì)與實(shí)現(xiàn) subplot(2,1,2)。 plot(abs(fft(x)))。grid。 ylabel(39。幅值 39。)。 xlabel(39。頻率 39。)。 title(39。濾波器輸入信號(hào) x頻域圖 39。)。 figure(3) plot(0:Nn,e(1,n:N))。hold on plot(0:Nn,y(1,n:N),39。r39。)。hold on grid on。 title(39。誤差信號(hào) e(t)和輸出信號(hào) y(t)39。) 。 ylabel(39。幅值 39。)。 xlabel(39。時(shí)間 39。)。 gtext(39。e(k)39。),gtext(39。y(k)39。) figure(4) subplot(2,1,1)。 plot(