大氣噪聲模型參數的非線性回歸估計

1 引言

在超低頻通信系統中，大氣噪聲的干擾一直困擾著通信工作者。大氣噪聲是由雷電等噪聲源引起，是一種典型的非高斯噪聲[1-3]。然而目前通信系統的最優接收機都是在高斯噪聲背景假設的前提下設計的，在實際中往往達不到預期的性能，甚至有時會導致接收機性能的惡化。隨著電子技術的蓬勃發展以及對通信質量的需求的不斷提升，對于大氣噪聲的研究越來越引起超低頻工作者的興趣[4,5]。大氣噪聲的研究，早在上個世紀50年代就已開展，但那時主要的工作集中于對大氣噪聲模型的建立。模型流派主要有兩類，一類是經驗模型，另一類是分析模型。其中分析模型雖然比經驗模型表達式更復雜，但因為有實際的物理過程支撐，且同實測數據也有更好的一致性，因此得到了長足的發展。分析模型中最著名也應用最廣的是Middleton的統計物理模型，該模型通過對噪聲源泊松分布的假設，詳細分析了噪聲傳播過程，推導了噪聲的概率模型[6,7]。該模型主要可分為Class A 噪聲模型和Class B噪聲模型。Class A 噪聲帶寬小于接收機帶寬，一般應用于寬帶系統；Class B噪聲帶寬大于接收機帶寬，較多應用于窄帶系統。由于在實際中，超低頻接收機一般設計為窄帶接收機，因此本文研究的大氣噪聲模型為Class B噪聲模型。

自Class B模型的誕生之后，由于其概率密度函數表達式復雜，信號處理難度較大，因此有些學者將α穩定過程作為Class B模型的一種簡化形式，以此來研究非高斯噪聲[8,9]。然而，在實際中，接收機不僅受到非高斯噪聲的影響，而且同時會受到高斯噪聲影響。這部分高斯噪聲不僅會來自于接收機系統內部，也會來自于外部。因此，在實踐中，α穩定過程作為大氣噪聲模型在性能比Class B模型要差[7]。然而，由于Class B模型的復雜，對Class B噪聲的研究甚少，特別在Class B參數估計上，只有Middleton的一些成果。而且這些成果一般都基于經驗分析，在實際中性能不理想，特別是難以應用于實時信號處理中，特別是接收機的設計中。本文從特征函數出發，采用非線性回歸估計，并設計了初始值估計方案和特殊序列，快速估計了Class B噪聲模型參數。實驗表明，該方法精度高，迭代收斂快，有很高的實用價值。

2 Class B模型

Middleton通過對噪聲產生和傳播機制的建模，推導出Class B噪聲模型。Middleton假設噪聲由高斯背景噪聲N(t)和脈沖噪聲X(t)構成。X(t)由一系列接收到的干涉波形組成。

其中Uj(t,θ)代表第j個干擾源產生的波形，θ代表相關的各種參數。假設在時間和空間上波形產生服從泊松分布，最終可推得概率密度函數[7]

其中A是“重疊系數”(overlap index)，也稱為“脈沖系數”(impulsive index)，指干擾源在單位時間平均發射的數目。α是空間傳輸密度因子，一般取0＜α＜2。Ω為歸一化因子。Γ(x)代表Γ函數，而Φ(a;b;x)指的是合流超幾何分布。從式(2)可以看出Class B噪聲模型的概率密度表達式是級數形式，表達式復雜，且是正負項交替，收斂速度慢。也正因為這，使得對Class B噪聲模型參數估計造成了很大難度。雖然其概率密度表達式復雜，但特征函數卻相對簡單，且是閉合形式[7]。

其中b1α是α的函數，是高斯噪聲分量。如果單獨取式(3)中的第1項作為近似，則Class B噪聲模型就退化成對稱α穩定分布(Sα S)。其實高斯分布本質上也是Sα S分布的特例，即α=2。從另一個方面說，Class B模型可以看成對稱α穩定過程(Sα S)和獨立高斯過程(G)的混和模型。下一節將從式(3)出發，推導非線性回歸模型。

3 非線性回歸估計

不妨將式(3)中，用γ替換Ab1α，指代Sα S過程中的離差(dispersion)，于是式(3)等價于

根據特征函數的定義有

若設觀測到的數據為{xk},k= 1 ,… ,N，則式(5)的估計為

結合式(4)和式(6)有

因此通過特定的序列{λk}，可以利用非線性回歸方法求解式(7)。

設θ= {α,γ,δG} = {θ1,θ2,θ3}為F(iλ)的參數空間。將F(iλi,θ)在θ0展開有

對于式(9)，可以通過線性回歸方法求解，即

設平方誤差為

因此，為了使在迭代過程中滿足式(12)不斷收斂，優化后的非線性回歸算法為：

步驟1 估計初始參數θ0;

4 初始值估計和{λk}序列的生成

由于對估計Class B的參數的文獻較少，沒有已知的較好的經驗估計方案，因此本文提出一種新的思路粗略估計Class B參數作為初始值。可以看到，Class B噪聲模型本質上是兩個α穩定過程的混和模型，其中一個是α=2的高斯過程，另一個是α＜ 2 的穩定過程。而在大氣噪聲中，往往是α＜2的穩定過程表達得多，占了主要的部分。因此，不妨先假設接收到的數據{xk}滿足α＜2的穩定過程，首先估計參數α和γ，然后在估計高斯分量的δG。

觀察大氣噪聲數據中，其非高斯特性的表達主要由高幅度的脈沖體現出來。若將高幅度脈沖數據剔除掉，那么其非高斯特征就會減弱。因此，設f0.75(xk)為75%分位點，取xk＜f0.75(xk)的所有值，并計算其方差，即可得到估計值σ2。由于在其中包含了非高斯部分的表達，σ2的值會高于。進一步扣除掉非高斯噪聲的能量，有 =σ2-2γ(α穩定過程γ相當于高斯過程中的/2)，這樣就能較好地初步估計出了高斯分量的δG。

解決了迭代的初始值估計，還有存在{λk}序列選取的問題。由于實際中一般可以假設收到的大氣噪聲數據服從對稱分布，因此有

下面我們引入定理[10]：

圖1 N (i λ )的零點擺動現象

定理1對于任意固定值N＜∞，當k→∞時，有TN,k→AN。(證明從略)

因此，只要構造足夠長的序列{TN,k}，就可以無限逼近AN。在實際中，為了計算方便，可以取α= 1 。只要獲得了AN的值，就可通過將區間(0,AN)劃分為T等分，來構造序列{λk|λk=k(AN/T),k=1,…,T} 。

5 仿真及實際結果

在本文中統計的參數為θ= {α,γ,δG} ={θ1,θ2,θ3} 。由于Class B噪聲模型可以看成是對稱α穩定過程(Sα S)和獨立高斯過程(G)的混和模型。因此噪聲的產生由兩個噪聲源相加，一個Sα S噪聲源，利用文獻[11]中的方法產生；另一個是高斯噪聲源。一般地，以均平方相對錯誤率(MSNRE)作為標準，即

其中此時的N為仿真實驗次數，θi為第i個參數實際值，為第j次實驗第i個參數的估計值。

在仿真實驗中，設置α=0.5,γ=2,= 5 ,{λk}的長度為T為40。每次實驗運行100次。當觀測數據長度K為2500時，圖2為經驗曲線、理論曲線以及估計曲線三者關系的示意圖，X軸為λ的值，Y軸為lnF(iλ)，其中經驗曲線通過FN(iλ)的估計獲得。可以看到，本文提出的方案估計出的參數與實際參數有很高的一致性。如圖3所示，當觀測數據長度為1000, 2500, 5000時，迭代次數同MSNRE的關系圖。數據表明，本方案迭代次數少，一般都能在10次迭代中得到收斂，同時從圖中容易觀察到隨著數據長度的增加，迭代的精度也隨之增加。

表1 α的初步估計值

圖2 經驗曲線，估計曲線及理論曲線關系

圖3 不同數據長度對MSNRE的影響

圖4 大氣噪聲信號

圖5 實測數據及估計結果的APD

6 結論

通過對Class B噪聲模型的研究，本為提出了一種利用非線性回歸模型估計Class B噪聲模型的方法。由于概率密度函數表達式復雜且收斂慢等因素，本文從特征函數出發，推導了非線性回歸模型，并優化了迭代算法，使其加快了收斂速度并減少迭代次數。特別地，利用Class B模型的特殊性，設計了初始值估計方案。同時，分析實際中情況，設計了{λk}選取的方案。最后仿真結果及實測數據的估計結果表明，本文提出的算法收斂快，精度高，估計出的參數能很好地反映實際情況。本文的研究結果對超低頻非高斯噪聲信道的檢測和估計具有較重要的理論和實際意義。

[1]Abraham D A. Detection threshold approximation for non-Gaussian backgrounds[J].IEEE Journal of OceanicEngineering, 2010, 35(8): 335-341.

[2]許玲, 靳致文, 王西乾. SLF/ULF大氣噪聲譜密度全球分布預測[J]. 電波科學學報, 2010, 25(4): 773-778.

Xu Ling, Jin Zhi-wen, and Wang Xi-qian. Estimation of SLF/ULF noise global distribution[J].Chinese Journal of Radio Science, 2010, 25(4): 773-778.

[3]Goken C, Gezici S, and Arikan O. Optimal signaling and detector design for power-constrained binary communications systems over non-Gaussian channels[J].IEEE Communications Letters, 2010, 14(2): 100-102.

[4]蔣宇中. 非高斯噪聲中的 Turbo碼的性能改進研究[J]. 應用科學學報, 2006, 24(4): 336-340.

Jiang Yu-zhong. Performance improvement of Turbo code in non-Gaussian noise[J].Journal of Applied Science, 2006,24(4): 336-340.

[5]Rosado-Munoz A. FPGA implementation of an adaptive filter robust to impulsive noise: two approaches [J].IEEE Transactions on Industrial Electronics, 2011, 58(3): 860-870.

[6]Middleton D. Statistical-physical models of electromagnetic interference[J].IEEE Transactions on Electromagnet Compact, 1977, 19(5): 106-127.

[7]Middleton D. Non-Gaussian noise models in signal processing for telecommunications: new methods and results for Class A and Class B noise models[J].IEEE Transactions on Information Theory, 1999, 45(4): 1129-1149.

[8]Nikias C L and Shao M. Signal Processing with Alpha-Stable Distribution and Applications[M]. New York: Wiley, 1995.

[9]Khalil H K. Performance of an optimal receiver in the presence of alpha-stable and Gaussian noises[J].IEEE Statistical Signal Processing Workshop(SSP), 2011, 56(8):573-576.

[10]Welsh A H. Implementing empirical characteristic function procedures[J].Statistics and Probabilitys Letters, 1986, 4(7):65-67.

[11]Chambers J M, Mallows C L, and Stuck B W. A method for simulation stable random [J].Variables, Journal of American Statistical Association, 1976, 71(5): 340-344.