基于直接支持向量機的蒸發波導參數反演研究

呂雍正，芮國勝，田文飚

(海軍航空工程學院，山東煙臺 264001)

大氣波導是電磁波在大氣傳播中的一種異常現象，電磁波會改變原來的傳播路徑，被陷獲在一定高度和范圍內的大氣波導層中以很小的損耗進行超視距傳播，進而影響到雷達和通信系統的使用。隨著我國海洋事業的發展，對海洋大氣信息的全面掌握尤為緊迫。大氣波導影響無線電設備的使用，獲得實時大氣波導信息對于艦船通信導航、雷達偵察探測和武器精確制導具有重要軍事意義。

對海洋大氣波導的反演研究由來已久，Krolik［1］最早提出利用雷達海雜波進行大氣波導的反演研究，后期的反演研究多由此展開;Yardim 等［2］利用馬爾科夫蒙特卡洛方法反演大氣折射率;Caglar Yardim 等［3，4］分別采用遺傳算法和馬爾科夫蒙特卡洛混合算法、卡爾曼濾波和粒子濾波方法估計了波導大氣折射率;這些方法對于大氣折射率剖面的反演研究都提供了一定思路，但對多參數反演仍存在穩定性差、精度不高和收斂速度慢等問題。

傳統支持向量機(SVM)針對有限樣本，能夠克服局部尋優、泛化能力差和多維數問題，Douvenot［5］采用SVM 進行了波導參數的反演，但存在訓練速度慢、算法復雜以及檢測階段運算量大等缺點;楊超［6］提出利用最小二乘支持向量機(LSSVM)反演蒸發波導，LSSVM 把最小二乘法引入到SVM中，簡化了計算復雜性，一定程度上提高了反演效率，線性方程組的系數矩陣對稱但不正定，影響一些快速算法的使用;本文采用直接支持向量機算法(DSVM)，將超平面的閾值平方項加入到非線性最小二乘支持向量機的目標函數中，求解算法更簡單，訓練速度快，泛化能力也沒有明顯減弱。

1 蒸發波導模型選擇

蒸發波導建模過程就是構建合適的折射率剖面模型。蒸發波導模型選取的準確與否關系到能否精確有效地反映真實的電波傳播環境。環境模型中描述折射率剖面的參數越多，即反演參數向量的維數越高，反演的剖面與實際大氣折射率剖面就越接近。實際傳播中折射率剖面隨距離而變化，但通常認為在水平方向上的變化是緩慢的，因此可以假定在水平距離的某特定間隔內其剖面是不變的。

在電磁波的大氣傳播中，大氣折射率可以由大氣溫度、相對濕度和壓強3 個影響參數表示

其中:P 是大氣壓強(hPa);T 是大氣溫度(K);e 是水汽分壓(hPa)。

在實際計算分析中，通常將地球曲面近似為平面處理，用修正折射率M 來表示

其中:z 表示距離地表的高度(m);re表示地球半徑(m)，為常數。

當修正折射率M 和折射率N 滿足如式(3)條件下會出現大氣波導現象

大氣波導現象通常分為以下幾類:蒸發波導、表面波導、抬升波導以及由3 種波導相互作用形成的混合波導。實驗統計，在近海面大氣層中，尤以蒸發波導出現的概率最大，通常發生在海面40 m 以下，如我國的南海和東海海域其發生概率高達80%以上。據此，本文主要研究蒸發波導，從美國海軍物理實驗室提出的五參數經驗模型［7］中選取雙參數模型［8］作為蒸發波導的模型，能夠克服單參數模型下波導高度與強度不一一對應的情況。模型表示如圖1 所示。

圖1 雙參數蒸發波導模型

式(4)中:M0為底層修正折射率指數，一般取330 units/m;c0為中性層結折射率剖面有關的經驗值，取0.13;c1為波導底層斜率，取值區間為［-1，0. 4］;z0為粗糙因子，取1. 5 ×10-4;zb為蒸發波導頂層高度，取值范圍為［0，500］，單位m;d 為蒸發波導高度，當d=0 m 時，對應的是標準大氣情況下的折射率剖面

M1=c0dln (zd/z0)+zd(c1-c0)。一旦確定了波導頂層高度zb，就剩下了d 和c12 個參數。

2 蒸發波導傳播的拋物方程模型

拋物方程是在求解波動方程時得到的一種近似形式，相比于射線追蹤和波導模理論，拋物方程方法對遠距離的電磁波傳播可以提供更加可靠、穩定的計算。

在直角坐標系中，二維標量波動方程在折射率為n(x，z )的各向同性介質中的場量ψ 滿足下面的亥姆霍茲方程［9］

式中:k0=w/c 是真空中的波數;n=n (x，z)是介質的折射指數;x 和z 分別表示地表面的水平距離和距離地表的高度。引入x 軸向衰減函數u (x，z) =ψ(x，z) e-ik0x，得到關于u 的標量方程

通常傳播仰角很小，只考慮前向傳播的電磁波，因此所對應的拋物方程的解

式(7)的求解是一個步進的過程。只要給定一個初始值u (x0，z)，就可以遞推得到其他的區域 (x0+Δx，z )內的場分布。

在一定大氣波導環境下，接收到的雜波功率可以用雷達方程計算得到

其中:L 為傳播損耗，是大氣波導參數M 與雷達到海面的距離r 的函數，L 可以由拋物方程來求解，表達式為L = Lf-20lgF，F 為傳播因子，Lf為自由空間傳播損耗;σ0是海面散射系數，處理時假設它是區域獨立性的，當成常數;C 為與雷達發射功率、天線增益等系統參數有關的常數。

3 直接支持向量機

支持向量機(support vector machine)是由統計學習理論發展而來。它以結構風險最小化原則和VC 維(vapnik chervonenkis dimension)理論為基礎，取代傳統學習算法中的經驗風險最小化原則。SVM 基本思想是:通過非線性映像把輸入向量映像到高維特征空間，通過引入核函數，從而巧妙地將非線性估計問題轉化為高維特征空間的線性估計問題。由于SVM 建立在嚴格的統計學習理論基礎之上，不會陷入局部尋優，有更強的泛化能力，較好地解決了小樣本、非線性、高維數和過學習等實際問題。DSVM 問題的求解規劃［10］

式中:w 為權向量;b 為閾值;ξi∈RN×1為誤差向量;C ＞0 為正則化因子。

引入Lagrange 乘子αi，αi∈RN×1，定義其Lagrange 函數為

根據KKT(Karush-Kuhn-Tucker)優化條件，對w，b，ξ，α求偏導數并令其為0，可得

帶入消去w、ξi，得到方程

使用核函數K=K (A，A) =φ(A)φ(A)T，有方程

得到回歸函數為

4 仿真分析

蒸發波導的參數反演是復雜的非線性過程，將DSVM 引入到反演蒸發波導的過程中。選取目標函數f 來判別仿真信號功率Pcal與計算得到的信號功率Pobs的符合程度，用于確定構建的波導參數M 與實際波導參數的誤差程度。這里選取是Pcal和Pobs的平均值。

4.1 DSVM 對電磁傳播損耗的訓練分析

給定一組維數為m 的蒸發波導高度z，就可以相應得到m 個M 剖面。利用拋物方程法和雷達方程得到m ×N 的正向功率矩陣。

蒸發波導的高度一般不超過40 m，選取波導頂層高度d=20 m，發射天線高度10 m;電磁波頻率區間［1，2］GHz，采樣頻率0.1 GHz，共11 個采樣頻率;接收高度區間［10，28］m，采樣頻率2m，共10 個采樣高度;在每個接收高度處檢測的傳播距離從10 ～60 km 以200 m 為間隔，取251 個點的傳播功率損耗作為訓練數據，共2 510 個訓練樣本。輸出矩陣對應不同頻率和傳播距離處的傳播損耗組成。并以此輸入矩陣和輸出矩陣作為訓練數據來訓練DSVM。

為了驗證預測結果的準確性，給出拋物方程(PE)的計算結果與訓練模型的傳播損耗值的對比圖。選取頻率為1.5 GHz，接收高度分別為10 m，20 m 時DSVM 的預測損耗值與LSSVM、拋物方程計算值的對比圖2、圖3。

圖2 10 m 高度傳播損耗

圖3 20 m 高度傳播損耗

從訓練模型的速度上看，LSSVM 用時183 s，DSVM 用時179s，兩者所用時間相當，沒有明顯優勢。從對比圖中看出，2 種方法的訓練結果都與拋物方程給出的傳播損耗相近，但DSVM 的訓練結果曲線與拋物方程方法計算值吻合結果更精確一些，預測更趨近于拋物方程計算值。據此仿真結果表明對DSVM 在傳播損耗的精度訓練方面有更好的表現。

4.2 反演折射率剖面準確性驗證

折射率剖面與蒸發波導的高度之間的關系式在前面的模型中已經給出，只要反演出高度值和斜率，相應的折射率剖面就可以得到。與前面訓練類似，取發射天線高度10 m;對于蒸發波導高度的訓練數據，取值區間為［10，40］m，采樣頻率取3 m，每個高度處對應的接收功率值共11 組訓練數據;在每個接收高度處檢測的傳播距離從10 ～60 km 以200 m 為間隔，取251 個點的傳播功率損耗作為訓練數據;組成由11 ×251 的輸入矩陣，輸出矩陣為對應每個波導高度值組成的11 ×1 矩陣。以此輸入輸出矩陣作為訓練數據。

同樣的，取蒸發波導高度為18 m，拋物方程計算給出的各對應點處的接收功率值進行反演。由第1 節給出的參數模型作為真實值，對比DSVM 和LSSVM 方法得到的反演折射率剖面圖4 可以看出，2 種優化方法給出的反演剖面與真實值接近，訓練結果如表1 所示，可以看到，在訓練用時方面，DSVM 比LSSVM 要節省1/4 左右，反演得到的蒸發波導高度，DSVM 為19.1 m，LSSVM 為20.3 m，LSSVM 反演結果誤差相對較大，DSVM 給出的反演結果更貼近真實值，適用于后續的反演。

圖4 反演折射率剖面對比圖

表1 反演用時和誤差比較

5 結論

在反演海上蒸發波導的背景下，本文基于直接支持向量機算法進行大氣折射率參數的反演。以雙參數模型作為本次反演的蒸發波導模型，以拋物方程法計算相應傳播損耗，利用改進的算法對大氣折射率參數進行優化反演，仿真實驗結果表明，新的學習算法在模型的訓練速度方面有更為優越的體現，反演結果也更為精確，在蒸發波導反演中有更有效的應用。

［1］Krolik J L，Tabrikian J.Tropospheric refractivity estimation using radar clutter from the sea surface［J］.Proceedings of the 1997 Battlespace Atmospherics Conference，1998:635-642.

［2］Yardim C，Gerstoft P，Hodgkiss W S.Estimation of radio refractivity from radar clutter using Bayesian Monte Carlo analysis［J］. IEEE Transactions on Antennas and Propagation，2006，54(4):1318-1327.

［3］Yardim C，Gerstoft P，Hodgkiss W S.Statistical maritime radar duct estimation using a hybrid genetic algorithms-Markov chain Monte Carlo method［J］.Radio science，2007，42:1-13.

［4］Yardim C，Gerstoft P，Hodgkiss W S. Tracking refractivity from clutter using Kalman and particle filters［J］. IEEE Transactions on Antennas and Propagation，2008，56(4):1058-1070.

［5］Douvenot R，Fabbro V，Gerstoft P，et al. A duct mapping method using least squares support vector machines［J］.Radio science，2008(43):1-12.

［6］楊超，郭立新，吳振森.最小二乘支持向量機在蒸發波導預測中的應用［J］.電波科學學報，2010(8):632-637.

［7］Peter G，Rogers L T，Krolik J L. Inversion for refractivity parameters from radar sea clutter［J］. Radio Sci，38(3):1-18.

［8］王紹班，孫小東，田文飚，等.基于乘子法的海上蒸發波導反演研究［J］. 中國電子科學研究院學報，2013(6):605-610.

［9］M F Levy. Parabolic equation methods for electromagnetic wave propagation［M］.London:IEE Press，2000.

［10］杜喆，劉三陽.最小二乘支持向量機變形算法研究［J］.西安電子科技大學學報，2009(4):331-337.

［11］程煥，謝洪森，孫大軍，等.用雷達海雜波反演蒸發波導的蟻群算法［J］.四川兵工學報，2013(1):91-93.