一種新型跟蹤離散1/f噪聲信號遞歸RLS 算法

曹昌勇

皖西學院 機械與電子工程學院，安徽 六安237012

1 引言

1/f 噪聲是一種普遍存在于物理、生物、社會、醫學、天文、地理等系統中的自然現象，1925年約翰遜（Johnson）在電子管中首次觀察到1/f噪聲。在以后的研究中，1/f噪聲一直是電子元器件和電子線路低頻噪聲的主要形式[1]。20 世紀80 年代以來的研究發現，電子元器件中的1/f噪聲的大小與電路的可靠性有密切的關系。一個1/f 噪聲很大的元器件或電路，很可能存在較嚴重的可靠性問題，其壽命則會較短，而且抗惡劣環境的能力也就較差。因此，人們通過檢測和分析電子系統1/f噪聲的性質，就可以間接地預測電路或系統的質量和可靠性[1-2]。自從發現1/f 噪聲的近百年時間里，人們對它的產生機理進行了大量的研究，提出了各種不同的模型，但仍然還沒有一個統一的理論。人們發現1/f 噪聲與系統的缺陷有關，是一個復雜的動力系統[3]，滿足長期統計規律，并且是長程相關的。對如何生成1/f噪聲的研究，也有不少的模型，比如分數階布朗運動過程（fBm）；通過混沌模型[4]來生成1/f 噪聲，該方法用一個logistic 混沌序列輸入到hosking 分數階差分濾波器中，并能產生具有混沌特性的1/f噪聲序列；通過小波變換[5]，該方法根據白噪聲的小波變換后的系數仍服從白噪聲規律，然后正交歸一化小波系數找到滿足1/f信號生成定理的正交小波系數集，通過Karhunen-Loeve 展開式生成1/f 信號。對1/f 噪聲的參數估計，主要采用的辦法有小波變換、最小二乘、時頻分析等；對于跟蹤1/f噪聲，普遍采用自適應算法，比如最小均方誤差（LMS）算法，文獻[6]提出了一種變步長LMS 算法跟蹤fBm 過程。由于LMS 算法收斂速度慢，因此本文采用改進的最小遞歸二乘法（RLS）來跟蹤1/f噪聲，但是這種算法存在計算復雜，迭代步數增加到一定程度濾波器系數必定發散等缺點。

自適應算法被廣泛應用在系統辨識、預測建模、噪聲對消以及信道估計等領域。它與濾波器一起構成自適應濾波器，可以在一定程度上對噪聲進行過濾，基本原理如圖1 所示。

圖1 自適應濾波器原理圖

自適應算法有很多，用的最多的要數LMS 和RLS 算法。為了更好地對非平穩信號進行處理，需要對算法進行改進。人們在改進這兩種算法的過程中做了大量的工作。對LMS 算法的改進一般采用的辦法是變步長[6-7]，或者是考慮引入梯度算法[8]或其他算法實現對快速變化的信號進行跟蹤。對于RLS 算法來說，它的收斂速度比LMS 要快得多，但它是以增加計算量為代價的。影響計算復雜度的關鍵在于該濾波器算法的結構上。文獻[9]提出了一種分級實現的RLS 算法，它把大的濾波器階數分解成幾個級聯的階數小的濾波器；文獻[10]通過輸入濾波器的信號的自相關函數是服從遞歸規律建立了擴展算法（ERLS）；文獻[11-12]通過誤差反饋的辦法，在文獻[9]的基礎上增加一個誤差判斷開關，通過開關切換去修改變化比較大的級聯濾波器的系數，達到部分更新濾波器系數的效果，進而減少運算量。對RLS 算法的收斂速度上的改進主要是實時修正RLS 算法的記憶因子[13]。對于非平穩信號來說，主要是期望能快速根據輸入信號的變化來調整濾波器系數。

正規或者快速RLS 算法是建立在迭代計算基礎上的，因此在有限的計算機精度的條件下，必定要發散。為了解決濾波器權系數發散的問題以及快速準確的跟蹤信號變化的問題，本文根據輸入濾波器的信號的自相關矩陣應該滿足慢變化的規律，引入誤差的非線性函數來控制這個變化以實現實時快速跟蹤1/f信號。該算法繼承了正規RLS算法的快收斂特性，并且不增加正規的RLS 算法的計算量，誤差也比正規RLS 算法有所減小。通過該方法對混沌1/f噪聲過程與fBm 過程進行跟蹤仿真，實驗表明它具有良好的穩定性，算法收斂速度較正規RLS 算法沒有降低，在小的濾波器階數的條件下也具有良好的跟蹤能力。該算法在遇到突變輸入時會快速調整濾波器系數，因此在短時間內濾波器系數的改變比正規RLS 算法要快得多，這就實現了快速跟蹤目標，從而減小了誤差。

RLS 算法中記憶因子對算法的跟蹤速度有比較大的影響。當記憶因子比較小時，算法響應快變信號的能力強，但抗干擾能力弱，適合于跟蹤非平穩信號；當記憶因子較大時，算法響應快變信號能力比較差，但抗干擾能力強，適合跟蹤平穩信號。記憶因子的大小與抗干擾強弱是矛盾的。因此記憶因子的調整量太大和太小都不好，本文給出了記憶因子與輸入信號自相關矩陣的特征值（fBm 噪聲信號）的一個表達式，為動態調整記憶因子提供了理論依據。

2 1/f噪聲的兩個模型

所謂1/f過程是指隨機過程x（t）的功率譜是冪函數的形式，即S(f)=σ2/f a，其中σ2是常數，由于其產生機理不明，一般都把它當做是一種隨機過程，它在低頻部分的功率特別大，特別是頻率靠近0 的時刻。

2.1 1/f混沌模型

文獻[4]通過采用logistic 方程生成混沌序列，然后把混沌序列輸入到hosking 分數階濾波器，表示如下。

其中，H(Z)為hosking 分數階濾波器的傳遞函數，其沖擊響應為：

其中，u(n)是單位階躍函數，式（1b）的離散遞推表達為：

通過求式（1）與式（1c）的卷積即可得到混沌1/f噪聲信號。式（1）產生的混沌序列的均值為0，自相關函數與白噪聲的自相關函數相似。

2.2 分數維布朗運動過程

文獻[6]提出了一種變步長的LMS 算法跟蹤fBm 過程。如果一個時間連續的隨機過程的統計特性是尺度不變的就叫自相似過程，自相似可用公式表示如下：

z(t)即為一個具有自相似指數H（Hurst 指數）的自相似隨機過程，它對任意的標量c＞0 都具有自相似性。式（2）中的不等號僅僅在統計意義上對所有有限分布滿足自相似性。具有平穩增加的非平穩的，且具有自相似性的高斯隨機過程就是分數維布朗運動過程，計為BH(t)，一個fBm 過程的一階導數稱為分數維高斯噪聲（fGn）。離散的fBm 可表示如下：

Ts為采樣周期，因為分數維布朗運動具有自相似性，不失一般性，取Ts=1。離散fBm 的均值、方差與自相關函數表示如下：

引理1 Hurst 指數為0 ＜H ＜1 的離散第一階分數維布朗運動（1-fBm）的自相關函數R(n)（M×M 階矩陣）能被對角化為R(n)=QΛ(n)QT，如果它的近似估計在n 比較大時成立，這里Q(n)≈Q 是一個常數正交矩陣，Λ(n)=diag{λ1,λ2,…,λM-1,λM(n)}。 R(n)的前面M-1個特征值是與時間無關的，特征值λM(n)是時間n 的函數。所有特征值都與Hurst指數H 相關，指數H 刻畫了離散分數維布朗運動，具體見公式（2a）。其最大的特征值是時間n 和Hurst指數H 的函數。該特征值可用下式表示：

3 RLS 算法及修正

3.1 基本RLS 算法

假如濾波器階數為M，n 時刻的權系數為W(n)=[W1(n),W2(n),…,WM(n)]T，輸入數據為向量U(n)=[u(n),u(n-1),…,u(n-M+1)]T，則輸出設為y(n)其表達如下：

期望輸出為d(n)，表達如下：

其中W0為最優濾波器系數，假設η(n)是均值為0，方差為的高斯白噪聲，文獻[14]考慮了η(n)服從獨立同分布的噪聲，采用類A 噪聲模型進行分析，本文僅考慮高斯白噪聲的情況。n 時刻的估計誤差為：

基于最小二乘法的RLS 算法的最小化代價函數如下：

λ 為遞歸最小二乘算法的記憶因子，滿足0 ＜λ ＜1，它對當前的誤差記憶強，對過去時刻的誤差記憶弱。把式（5）代入式（6）并兩邊對W 求導數，令其導數為0，即

對式（9）利用矩陣求逆引理可得到其逆矩陣的遞推表達式：

令P(n)=R-1(n)為輸入信號的自相關矩陣的逆矩陣：

將式（11）與式（12）代入到式（8）得到權系數的遞推公式（13b）：

上面的四個表達式給出了常規的RLS 算法的一個迭代表達，要滿足迭代，需要初始化的項為P(n)和W(n)。一般W(0)初始化為長度為M 的0 向量，P(0)=δ-1IM，IM為M 階的單位矩陣，δ 當輸入信號的信噪比大時取小的正數，反之取大的正常數[15]。

在式（12）中，P(n)是輸入向量的時間平均相關陣的逆矩陣。

3.2 修正RLS 算法

本文對式（13c）P(n)進行修正，用一個非線性的誤差函數來修正它，并用式（15a）取代式（13c）。

公式（14）所對應的關系曲線如圖2 所示，其中a 取0.8，b 取3，可以根據程序設置不同的值來討論曲線的變化情況；公式（15）不需要贅畫，因為IM是一個單位矩陣。 f(e(n))是關于誤差e(n)的非負非線性函數，把絕對值與a 去掉即為參考文獻[3]中函數的形式。式（14）中的a 與b 是正數，一般與被跟蹤的信號有關，可以通過圖2 曲線圖討論來確定。誤差一般不為0，當誤差增大時P(n)增大，即表示信號在快速變化，信號發生了突然的變化，這需要在下一次比較大的調整權系數W(n)，通過觀察式（13），由于P(n)增加，導致k(n)增加比較大，從而W(n)得到了大的調整，說明信號是快變的；相反，當誤差很小時，表示濾波器權系數不需要進行大的調整，從式（14）可見小的誤差不會對原來的P(n) 產生大的改變，因此k(n) 的變化也比較小，從而W(n)的調整很小，這表明輸入信號在這一小段時間內是平穩的。另一方面，由于P(n)的改變，可以說是對它進行了一個變化的初始化工作。

圖2 F（n）同e（n）的關系曲線

對于RLS 算法，P(n)是集平均輸入向量自相關陣R(n)在n 時刻的逆矩陣，在正規的RLS 算法中，該項是一直迭代下去的，由于計算機精度問題，必定會導致P(n)不可逆，導致R(n)的特征值發生突變。特征值發生突變的根本原因是有部分非平穩數據的輸入，造成誤差增大，導致權系數調整過大，最終算法發散。因此一種辦法是通過關于誤差的非線性函數來修正權系數的調整[14]，即非線性RLS 算法。本文用非線性函數來微量調整P(n)矩陣，這樣可以避免特征值的突變，但算法仍然是線性的RLS 算法。

3.3 與λ 相關的一個表達式

下面推導記憶因子λ 與輸入信號（1-fBm）的自相關矩陣的特征值的一個重要關系。

把式（4）、式（16）代入式（13a）得到：

由式（12）、式（13b）、式（16）和（17）推出：

對式（18）兩邊分別取期望得到：

根據P(n)滿足慢變規律，可知它相對于輸入U(n)與誤差e(n)來說，它的變化顯得很小，因此可以把它看做是與輸入和誤差相互獨立的，故可表達如下：

由式（16）和式（13b）知，V(n)與U(n)獨立，假設觀測噪聲η(n)與U(n)相互統計獨立，且η(n)的均值為0，將式（17）代入到式（20）得到：

其中，RU(n)=E[U(n)UT(n)]為n 時刻輸入濾波器數據的自相關函數，在式（21）中，根據式（12）得到，且RLS 算法需要初始化P(0)=δ-1IM，即R(0)=δIM，故：

把式（22）、式（23）代入式（21），并考慮到輸入數據的自相關矩陣是緩變的，把當前時刻自相關陣與之前時刻的相關陣看做是可以近似統一的（記憶因子的存在也為這種近似統一提供了依據），因此得到：

由引理1 可知，式（24）RU(n)=QΛ(n)QT，其中Q 為正交的矩陣，Λ(n)是由RU(n)的特征值構成的對角陣，設其特征值為λ1,λ2,…,λM-1,λM(n)，其中λM(n)在式（2c）中給出。因此式（24）轉化為：

設B=I-Λ2(n)-1，由于，因此B 的特征值可以用Λ(n)的特征值來表示。式（26）給出了一個迭代過程。RLS 算法收斂意味著當n 趨于無窮時：

即當n 趨于無窮時，W(n+1)→W0達到最優值。

其中λΛ2(j,i)為Λ2在第i 次迭代過程中第j 個特征值。由可以把Λ2表示出來：

最后由式（28）得到收斂的條件為：

公式（30）給出了記憶因子λ 與輸入信號的自相關矩陣的特征值λΛ(j,n)在n 時刻的一個表達式，這為在跟蹤過程中動態調整RLS 算法的記憶因子λ 給出了理論上的依據。

公式（30）說明調整過程中算法的記憶因子是受到自相關矩陣的特征值λΛ(j,n)控制的。該公式是記憶因子滿足的一個關系表達式，可以根據這個公式來設置記憶因子，當然，如果記憶因子在每次迭代過程中的設置不滿足這個條件時，算法肯定會發散，從這個角度看，屬于收斂條件范疇的討論。

4 仿真結果

4.1 對混沌1/f信號跟蹤的仿真

根據第2章中的第一個模型產生1/f時間序列（初始值取0.3，hosking 分數階濾波器d=0.4），取3 000 個離散數據。設修正的RLS 算法濾波器階數為M=10，觀測噪聲方差為0.4，得到如圖3 的跟蹤結果（誤差均值與方差分別為0.046 5、0.467 4）。

圖3 修正RLS 算法對混沌1/f噪聲序列的跟蹤（a=5，b=0.005）

標準RLS 算法的跟蹤結果如圖4 所示（誤差均值為0.047 7，方差為0.480 4）。

圖4 標準RLS 算法跟蹤混沌1/f噪聲信號

從對混沌1/f噪聲信號的跟蹤情況看，修正RLS 算法比標準RLS 算法弱占優勢，其均值和方差均比標準算法小。

4.2 對1-fBm 噪聲信號的跟蹤

取H=0.7，濾波器階數為10，修正算法與標準算法結果如圖5（誤差均值與方差分別為0.001 7、0.104 0）和圖6 所示（均值為-0.005 6，方差為0.274 2）。

圖5 修正RLS 算法對混沌1/f噪聲序列的跟蹤（a=5，b=0.5）

圖6 標準算法跟蹤1-fBm 信號

從對混沌1/f 信號或1-fBm 的結果看，修正的算法比標準算法要好。從圖7 與圖8 在對跟蹤1-fBm 過程中濾波器系數的調整過程可以看出，修正算法濾波器的系數可以發生快速改變，而標準算法只是很慢地變化，這也導致它的誤差比較大，跟蹤速度達不到的原因。

圖7 修正算法對濾波器系數的調整過程

圖8 標準RLS 算法對濾波器系數的調整過程

圖7、圖8 中不同的顏色曲線表示的是RLS 算法中各個濾波器系數的動態曲線，由于濾波器系數一般都不止一個，因此為了區別每一個系數的動態變化，作圖時采用了不同的顏色。

從圖7 可見，迭代時會出現濾波器系數的調整過程；對比圖8 的系數調整可見，算法不會影響收斂速度。

5 結束語

RLS 算法是一種快速收斂的算法，但是在跟蹤非平穩信號時無法快速調整濾波器系數以便實時跟蹤，并且在迭代次數達到很大時算法必將發散，本文通過引入一個非線性函數來對輸入信號相關逆矩陣進行修正，這樣可以快速調整濾波器系數，以適應信號的快速變化。對兩個1/f噪聲模型產生的序列進行仿真實驗，結果表明修正方法比正規的RLS 算法要好。此外，本文還推導了RLS 算法的記憶因子與輸入信號自相關矩陣特征值的一個關系，這為動態調整記憶因子提供了理論上的依據。

[1] 莊奕琪，馬中發，杜磊.1/f 噪聲之謎[J].世界科技研究與發展，1999，21（4）：69-72.

[2] 楊冬云，王司.1/f 分形噪聲理論及其在信號處理中的應用研究綜述[J].黑龍江工程學院學報：自然科學版，2004（3）：31-35.

[3] 許生龍.1/f噪聲的動力學統計理論[J].紅外技術，2003（5）：63-67.

[4] 張強，馬潤年，許進.1/f過程的混沌模型[J].信號處理，2001（5）：392-394.

[5] 李光，于盛林.基于小波變換的擬1/f 信號生成[J].數據采集與處理，2003（1）：119-122.

[6] Gupta A，Joshi S D.A novel least mean squares algorithm for tracking a discrete-time fBm process[C]//Proceedings of 2006 Annual India Conference，2006.

[7] 覃景繁，歐陽景正.一種新的變步長LMS 自應用濾波算法[J].數據采集與處理，1997（3）：171-174.

[8] Ang W P，Farhang-Boroujeny B.A new class of gradient adaptive step-size LMS algorithms[J].IEEE Trans on Signal Processing，2001，49：805-809.

[9] Wu An-Yeu，Liu K J R.Split recursive least-squares-algorithms，architectures，and applications[J].IEEE Trans on Circuits and Systems Part II：Analog and Digital Signal Processing，1996，43：645-658.

[10] Haykin S，Sayed A H，Zeidler J R，et al.Adaptive tracking of linear time-variant systems by extended RLS algorithms[J].IEEE Transactions on Signal Processing，1997，45（5）.

[11] Jin Minglia.Partial updating RLS algorithm[C]//Proceedings of the 7th International Conference on Signal Processing（ICSP’04），Beijing，China，2004：392-395.

[12] Niu Wei，Wang Minxi，Chen Kaya.A novel reduced order RLS predistortion[C]//Proceedings of Microwave Conference，Asia-Pacific，2005.

[13] Park D，Jun B E，Kim J H.Fast tracking RLS algorithm using novel variable forgetting factor with unity zone[J].Electronics Letters，1991，27：2150-2151.

[14] Weng J F，Leung S H.Nonlinear RLS algorithm for amplitude estimation in class A noise[J].IEEE Proceedings Communications，2000，147：81-86.

[15] Wei P C，Han Jun，Zeidler J R，et al.Comparative tracking performance of the LMS and RLS algorithms for chirped narrowband signal recovery[J].IEEE Transactions on Signal Processing，2002，50：1602-1609.

[16] Chun B，Kim B，Lee Y H.Generalization of exponentially weighted RLS algorithm based on a state-space model[C]//Proceedings of the IEEE International Symposium on Circuits and Systems.[S.l.]：IEEE Press，1998：198-201.

[17] Yeo H K，Sharif B S，Hinton O R，et al.Improved RLS algorithm for time-variant underwater acoustic communications[J].Electronics Letters，2000，36：191-192.