基于Bayes-MCMC 的水聲雙程信道建模及自適應采樣反演

趙 罡 ,孫乃葳 ,申 珅 ,楊益新

(1.西北工業(yè)大學 航海學院,陜西 西安,710072;2.中國船舶集團有限公司 第705 研究所,陜西 西安,710077)

0 引言

迄今為止聲脈沖仍是水下通信及目標探測的主要手段。在水下傳播過程中,受信道多徑及目標閃爍等多種因素的作用,聲脈沖信號會發(fā)生衰落和起伏等變化,這一變化對工作在淺海水域的魚/水雷等主戰(zhàn)武器,以及無人水下航行器(unmanned undersea vehicle,UUV)、微型潛艇、蛙人探測聲吶等水下特種裝備的作用尤其顯著,對水下通信質(zhì)量及目標檢測有效性帶來較大影響[1-3]。為消除這一影響,人們在針對水聲衰落信號的自適應處理方面開展了較多工作[4-5],但對于提供驗證環(huán)境相關(guān)的仿真方法和模型研究卻所做不多。實際水聲信道衰落情況復雜多變[6],外場試驗本身又存在成本高昂難以開展的問題,在包括魚/水雷、UUV及特種作戰(zhàn)裝備在內(nèi)的聲探測武器研發(fā)過程中,使用外場試驗的方式對信號處理方法進行改進和驗證顯然是不現(xiàn)實的?？紤]到信噪比(signal noise ratio,SNR)在目標檢測中的重要性[7],建立與實際情況有較高符合度的水聲雙程信道目標回波時域起伏仿真模型,并提供一種可有效利用實航數(shù)據(jù)的模型參數(shù)反演方法,是水下目標精細回波仿真所必須考慮的問題,在水下聲探測武器的研制過程中具有重要的意義。

1 水聲雙程信道建模

首先以經(jīng)典水聲單程信道模型為基礎進行雙程信道模型的建立: 考慮到一定時空不變條件下信道的對稱性將去程信道與返程信道特性等同,然后在直達和多徑信號強度2 個方向,通過積分計算二者聯(lián)合的全概率表達式,最后建立收發(fā)同置應用背景下的水聲往返雙程信道模型。

1.1 單程信道基本模型

水聲信道的信號傳輸主要受界面、介質(zhì)和障礙物的反射、衍射等作用,其信道物理過程復雜,難以解析描述,常利用經(jīng)驗數(shù)據(jù)擬合概率分布函數(shù)來表示信號幅度變化等重要特征,針對水聲環(huán)境特點,可通過直達信號和多徑信號的加性合成來表示,即

式中:z(t)為直達信號;d(t)為多徑信號。

魚/水雷、UUV 及特種作戰(zhàn)裝備常用高頻短脈沖作為目標聲源探測信號,其頻率高且傳播距離較近,因此開展信道物理過程復雜性建模時可做如下合理假設: 1)信道的直射分量受水體中障礙物遮擋的影響遠小于界面反射的多徑效應;2)傳播衰減造成的短脈沖內(nèi)慢起伏可忽略。此時可認為z(t)保持包絡不變、而d(t)受界面及溫度梯度的影響變化明顯,則r(t)符合Ricean 分布且滿足條件概率分布

式中: 忽略時間參量t;表示d(t)平均散射功率;I0(·)表示第1 類零階修正貝塞爾函數(shù)。

1.2 雙程信道模型

雙程信道模型的建立需聯(lián)合考慮信號去返2 次行程,即分別設定r1(t)、r2(t)為去程信號和返程信號,則二者對應的直達和多徑成分滿足

r2(t)是由水下目標體對去程信號r1(t)反射產(chǎn)生的,故其直達成分z2(t)的包絡與目標反射特性相似呈對數(shù)正態(tài)分布,可用r1(t)的條件概率表達,即

其中,表示去程信號的平均散射功率。

文中所用實航數(shù)據(jù)以聲周期方式保存,相對于聲在水中傳播時1 個聲周期尋的短脈沖時間片段內(nèi)可認為目標回波過程是確定的,其強度并不發(fā)生閃爍,此時雙程水聲信道模型可簡化為

2 水聲雙程信道參數(shù)反演算法

水聲雙程信道參數(shù)的反演問題中,非啟發(fā)式的非線性反演方法如遺傳算法(genetic algorithm,GA)、粒子群算法(particle swarm optimization,PSO)等受算法本身設計限制只能以概率密度函數(shù)(probability density function,PDF)差異作為代價函數(shù)開展反演[10],其對實航數(shù)據(jù)時域真值信息利用不足。為有效利用實航數(shù)據(jù)時域采樣先驗結(jié)果[11],流程中以馬爾科夫鏈蒙特卡洛(Markov chain Monte Carlo,MCMC)方法為主要方法實現(xiàn)水聲雙程信道模型抽樣,并進一步結(jié)合Bayes 方法的最大后驗估計(maximum a posteriori,MAP)完成參數(shù)反演過程,而反演參數(shù)的有效性判斷主要通過全局誤差評價模型對比分析海試數(shù)據(jù)、反演結(jié)果分布特性差異來進行。水聲雙程信道反演算法基本流程如圖1所示。

圖1 水聲雙程信道反演算法基本流程Fig.1 Basic process of underwater acoustic go-back channel inversion algorithm

一般而言,Bayes-MCMC 反演的基本思路是通過采樣構(gòu)造一條在給定真值樣本條件下穩(wěn)態(tài)分布服從所需模型參數(shù)向量后驗概率密度的馬爾科夫鏈,當數(shù)據(jù)鏈收斂時使用采樣數(shù)據(jù)計算所求模型參數(shù)集的各矩特征量,數(shù)據(jù)鏈同時作為以反演結(jié)果為參數(shù)的正演模型數(shù)據(jù)樣本[12]。MCMC 方法的核心在于構(gòu)造轉(zhuǎn)移核接收函數(shù),其設計原則與采樣器直接相關(guān),最為通用的采樣器包括Metropolis-Hastings(M-H)、吉布斯采樣(Gibbs sampling,GS)方法等。GS 方法作為M-H 方法的特例是對參數(shù)集逐一完成單維度正交采樣,其效率雖然較高但要求準確的參數(shù)條件建議分布,而水聲反演問題往往準確的模型參數(shù)向量先驗條件概率密度難于確定;M-H 方法使用任意建議分布在整個參數(shù)集空間整體采樣以生成數(shù)據(jù)樣本,且建議分布取對稱分布時可進一步簡化轉(zhuǎn)移核接受函數(shù)形式,具有較好的實現(xiàn)性,因此文中采用M-H 采樣方法[13]。

由于雙程信道實航數(shù)據(jù)分布特性已知的情況與通常MCMC 采樣先驗概率未知的假設條件不同[14]: 當先驗概率未知時,MCMC 采樣一般是獨立并發(fā)抽取2 組隨機序列,并以2 組序列收斂情況一致性作為停止采樣的判據(jù)[15],或者通過對數(shù)據(jù)采樣過程的設計達到降低計算量提高反演效率的目的[16-17];而雙程信道實航數(shù)據(jù)分布特性已知時,可據(jù)此對每次抽取結(jié)果的有效性進行判斷,并使用抽樣分布與實航數(shù)據(jù)分布的絕對差異或相似程度作為采樣狀態(tài)改變的判定條件。結(jié)合圖1 反演算法特征和M-H 采樣存在參數(shù)多維耦合條件下命中效率低的問題,為保證數(shù)據(jù)的高效抽樣,文中基于MCMC 進行了優(yōu)化,提出一種M-H 自適應單維度串行采樣算法,即通過M-H 方法進行精度自適應調(diào)節(jié)的單維度串行采樣過程。

首先給出如下定義: 定義模型參數(shù)向量m={m1,m2,···,mn}的取值空間為 ?;定義m中mi的補集為?mmi={m1,m2,···,mi-1,mi+1,···,mn}。

具體反演算法步驟如下:

1)設定目標反演偏差 εT、反演偏差增量 Δε、最大反演次數(shù)Rm,初始化反演偏差ε0=∞、反演次數(shù)Rc=0、反演參數(shù)序號i,j=1;

2)在 ?上隨機抽取或采用全局尋優(yōu)方法生成模型參數(shù)向量m={m1,m2,···,mn};

4)對每一個被接受的模型參數(shù)ms,使用計算模型正演PDF;

5)對模型正演PDF 及真值PDF 進行同步離散采樣并計算其面積誤差ε,當不滿足誤差收斂條件即ε ＞ε0時回到步驟3);否則接受,同時令i=i+1 及ε0=ε;

6)判斷 εT與 ε0關(guān)系: 若ε0＜εT,即目標反演偏差已經(jīng)達到,則令m=并停止反演;若εT＜ε0,則令m=,同時判斷若i＜n,即尚有模型參數(shù)分量未執(zhí)行M-H 采樣,則回到步驟3);否則反演次數(shù)Rc加1,同時令i=1;

7)當Rc＞Rm時即達到單輪最大反演次數(shù)時,為避免陷入局部優(yōu)化令εT=εT+Δε進行反演精度調(diào)整,并將反演次數(shù)計數(shù)置0 后回到步驟3);否則直接回到步驟3)。

算法步驟1)中的各初始設定參數(shù)的選擇可根據(jù)反演時預期的結(jié)果精度設定,但精度要求越高、反演次數(shù)越多,雖然會得到更為精確的反演結(jié)果,也會帶來更大的計算壓力,對算法執(zhí)行效率產(chǎn)生影響;算法步驟3)中的對稱分布僅僅是為滿足簡化M-H 采樣接受條件的需求,常用如均勻分布、高斯分布等都可使用?？偠灾?基于MCMC 方法的M-H 自適應單維度串行采樣算法可綜合考慮精度和效率因素,通過初始參數(shù)設置靈活控制參數(shù)采樣過程,且在算法終止條件上通過對全局優(yōu)化代價函數(shù)結(jié)果的判斷以達到有效收斂的目的,可以有效避免因初始參數(shù)設置不當而引起的反演結(jié)果偏差過大的問題。

3 反演算法有效性驗證

3.1 反演模型仿真設計

為驗證反演算法的有效性,首先對經(jīng)典分布模型參數(shù)進行反演驗證。Ricean 分布是描述海洋信道時常用的統(tǒng)計模型,當使用直達信號與多徑信號功率比即K因子進行描述時其表達式為

此時反演模型流程步驟中參數(shù)MAP 估計的主要任務是完成參數(shù)和K的估計。記真值樣本為Si,構(gòu)建Ricean 模型參數(shù)反演轉(zhuǎn)移核接受函數(shù)為

式中,Ri(m1)及Ri(m2)為符合式(7)在設定模型參數(shù)條件下的離散化采樣結(jié)果,其中m1和m2分別代表迭代采樣過程中前后2 次由σ和K構(gòu)成的參數(shù)向量。

由于各采樣獨立,v為各維數(shù)據(jù)誤差協(xié)方差矩陣Cd特征值積乘,通?？墒褂肎A、模擬退火(simulated annealing,SA)、PSO 等全局尋優(yōu)方法獲取極大似然估計?？紤]到文中方法Cd矩陣是由真值和采樣數(shù)據(jù)PDF 函數(shù)值殘差構(gòu)成,同時方法本身存在參數(shù)收斂性判斷,且v值與數(shù)據(jù)矩陣相似性反向相關(guān)在反演過程中逐漸趨小,因此可在收斂性判斷同時進行v的迭代更新,而其初值只需在0～1 范圍任意取值即可,不失一般性可將v設定為1。

3.2 仿真收斂條件設計

常規(guī)MCMC 反演一般采用針對2 組獨立并行采樣序列的非參數(shù)檢驗方法,如Kolmogorov-Smirnov(K-S)檢驗等的結(jié)果作為收斂判斷和算法終止條件,但在文中方法中MCMC 作為參數(shù)采樣過程主要用于提高采樣效率,而對反演過程的終止條件判斷可針對整體反演結(jié)果進行評價?；诖?設計針對模型參數(shù)反演結(jié)果的全局誤差模型如下

式中:p、q分別代表反演及真值樣本數(shù)據(jù)PDF,(p||q)代表2 個PDF 全局差異程度,當其小于εT時,認為反演目標達成,終止采樣,進行上述判斷的參數(shù)區(qū)間取Z2即二西格瑪水平(2-Sigma Level)。

全局誤差模型可對仿真采樣結(jié)果進行定量評價,但對于多組采樣及反演結(jié)果還需要有相似性的定性評價。在機器學習中,KL(Kullback-leibler)散度常用來對2 個概率分布的相似性進行評價,其定義為

KL 散度越小,代表2 個分布越相似,2 個相同分布的KL 散度值為0。文中將直接使用KL 散度指標對多組反演全局誤差模型計算結(jié)果進行定性評價。

3.3 算法驗證

3.3.1 信道模型真值樣本生成及評價

在水深一定的條件下,對于幾公里內(nèi)較近的傳輸距離,由于直達信道的存在以及多途情況不那么嚴重,以Ricean 描述的水聲信道表現(xiàn)為較小的σ值和較大的K值。取一組典型值采用MCMC 采樣方法生成符合由該組參數(shù)定義的Ricean 分布仿真樣本如圖2(a)所示,縱坐標r為采樣信號幅值;其采樣結(jié)果及其與真值PDF 和累積分布函數(shù)(cumulative distribution function,CDF)對比如圖2(b)和圖2(c)所示,其中縱坐標p為某幅值信號出現(xiàn)概率。

圖2 MCMC 采樣結(jié)果及其與真值分布對比Fig.2 MCMC Sampling results and comparison with their true distributions

仿真結(jié)果分析表明,與真值模型設定參數(shù)相比,各組采樣均可較快收斂。σ2參數(shù)仿真結(jié)果誤差0.4%,K參數(shù)仿真結(jié)果誤差0.5%,PDF 全局誤差1.2%。仿真采樣生成的Ricean 分布樣本與設定模型具有較高一致性,可以作為參數(shù)反演真值。

3.3.2 信道模型參數(shù)反演

使用所設計參數(shù)抽樣及結(jié)果評價模型,以上節(jié)所生成數(shù)據(jù)樣本作為真值對Ricean 分布模型參數(shù)開展反演,為比對反演結(jié)果將上節(jié)給定參數(shù)PDF作為擬合結(jié)果。反演設定條件為: 目標反演偏差εT=1%、反演偏差增量Δε=1%、最大反演次數(shù)Rm=5;每趟次參數(shù)采樣數(shù)5 000 個。

參數(shù)迭代反演過程以及結(jié)果對應的Ricean 分布PDF/CDF 與真值分布對比如圖3 所示。

圖3 模型參數(shù)采樣過程及反演結(jié)果與真值分布對比Fig.3 Sampling process of model parameter and comparison of inversion results with true distributions

反演結(jié)果分析表明,與正演模型設定參數(shù)相比,參數(shù) σ2反演結(jié)果誤差為0.8%,參數(shù)K反演結(jié)果誤差為1.5%,PDF 全局誤差為2.6%。各參數(shù)均有較好的收斂一致性,單維度串行采樣模型具有較高的反演精度和效率,可用于水聲雙程信道模型參數(shù)反演。

4 反演模型實測數(shù)據(jù)驗證

4.1 真值樣本預處理

不同于單程信道測量信號,水下自導武器聲尋的脈沖通常為矩形或梯形包絡的脈沖信號,其長度一般在10 ms 量級,而回波脈沖受傳播、吸收衰減及信道衰落等作用,在接收端信號幅值往往較低并常采用時變或自動增益技術(shù)以提升其SNR,這些因素的影響在信道模型參數(shù)反演時都需考慮。一般來說,為了得到原始信號數(shù)據(jù),需要基于補償曲線消除增益作用影響,同時為了得到較為準確的脈沖起始和終止時刻,副本相關(guān)等基本信號處理過程也是必須的。

在一段實際回波信號中選擇出反演所需的脈沖回波后,還需要對信號進行幅值歸一化處理,同時考慮到采樣及計算均在離散點上進行,所選擇PDF 值域及積分計算值域范圍均會對結(jié)果帶來一定的影響。即使經(jīng)過了預處理過程,受水文及界面環(huán)境的影響,弱回波脈沖信號被本底噪聲部分淹沒,仍會對反演結(jié)果帶來一定的影響,下面的反演就是針對SNR 條件不同的幾組回波展開的。以某湖上試驗數(shù)據(jù)為例,采樣時信道傳輸距離為500～2 000 m,2 組不同SNR 條件真值采樣如圖4 所示。

圖4 不同SNR 下脈沖回波時域信號Fig.4 Pulse echo time domain signals under different SNRs

4.2 反演模型設計

反演過程的M-H 抽樣模型轉(zhuǎn)移核接受函數(shù)形式及參數(shù)定義同式(8),不同之處在于符合設定模型參數(shù)的離散化采樣結(jié)果生成時,采取式(6)所定義不考慮目標閃爍條件的水聲往返雙程信道時域起伏模型。

反演結(jié)果評價同樣采用式(9)和式(10)模型開展,但是考慮到真實信號的模糊性,當信號樣本屬于較低SNR 條件(小于3 dB)時反演設定條件及相關(guān)評價準則放寬至目標反演偏差εT=0.03,其他不變;同時為了避免模型陷入單個參數(shù)采樣的局部優(yōu)化,在采樣算法的步驟5)設定單參數(shù)最大無效反演5 輪次后,以新的隨機參數(shù)向量重新開始反演。

4.3 模型參數(shù)反演

基于實際數(shù)據(jù)樣本對無目標閃爍條件的水聲往返雙程信道時域起伏模型參數(shù)開展反演。設定σ2隨機采樣區(qū)間[0.1 1],K隨機采樣區(qū)間[0.5 3],反演最終結(jié)果如下。

4.3.1 不同SNR 信號反演

表1 不同SNR 信號參數(shù)反演結(jié)果對比Table 1 Comparison of parameter inversion results under different SNR signals

2 組反演結(jié)果對應的Ricean 分布PDF/CDF 與真值分布對比如圖5 所示。

圖5 不同SNR 信號反演結(jié)果與真值分布對比Fig.5 Comparison of inversion results with true distributions under different SNR signals

反演結(jié)果表明,SNR 水平對于反演結(jié)果影響較大: 在相同發(fā)射接收水平下,SNR 越高反映信道質(zhì)量越好,代表直達信號強度的K因子越大,而表征多途信號功率的 σ2參數(shù)越小。此外,在較高SNR 條件下信號起伏分布規(guī)律與正演模型相似度更好,表現(xiàn)為用于評價PDF 整體差異的結(jié)果更小,表征分布相似性的KL 散度指標也更小。同時方差分析表明不同信號條件下,各評價參數(shù)均有較好的收斂一致性。最后,反演結(jié)果與真值的總體差異也表明所建立的水聲雙程信道仿真模型是有效的,尤其是在較高SNR 條件下。

4.3.2 多樣本反演

文中最后選擇了來自于某海域?qū)嵑皆囼灥? 組在時空上有一定聯(lián)系的實際信號進行了模型參數(shù)反演。其中,前8 組數(shù)據(jù)信道傳輸距離自1 000～2 000 m 非等間隔分布,SNR 為1～2 dB,第9 組數(shù)據(jù)作為對照SNR 條件約3 dB,信號時域形態(tài)如圖6 所示。

圖6 多樣本脈沖回波時域信號Fig.6 Multiple groups of pulse time domain echo

反演初始設定及過程參數(shù)與4.3.1 節(jié)相同,均進行100 趟反演。1～8 組反演結(jié)果對應的Ricean分布PDF/CDF 與真值分布對比如圖7 所示。

1～8 組反演結(jié)果如表2 所示。

各參數(shù)的變化趨勢如圖8 所示。

作為對照,對第5、9 組數(shù)據(jù)按同樣條件進行了100 趟反演,參數(shù)反演過程結(jié)果對比見圖9。

第5、9 組信號參數(shù)反演結(jié)果的對比如表3所示。

表3 第5、9 組信號參數(shù)反演結(jié)果對比Table 3 Comparison of inversion results at 5th&9thgroups signals

第5、9 組反演結(jié)果對應的Ricean 分布PDF/CDF與真值分布對比如圖10 所示,其中橫縱坐標與前文定義相同。

針對上述反演統(tǒng)計結(jié)果的分析表明:

1)從表1 及圖7 結(jié)果可知,各組反演過程均成功收斂,反演結(jié)果與真值誤差按模型所定義的總體差異均在5%以內(nèi);

圖7 1～8 組信號反演結(jié)果與真值分布對比Fig.7 Comparison of inversion results with true distributions at 1thto 8thgroups signals

圖8 1～8 組反演參數(shù)變化趨勢Fig.8 Change trend of inversion parameters at 1thto 8thgroups signals

表2 1～8 組雙程信道模型參數(shù)反演結(jié)果Table 2 Parameter inversion results of go-back channel model at 1thto 8thgroups signals

圖9 第5、9 組雙程信道模型參數(shù)反演過程結(jié)果對比Fig.9 Comparison of parameter inversion results of go-back channel model at 5th&9thgroups signals

圖10 第5、9 組信號反演結(jié)果與真值分布對比Fig.10 Comparison of inversion results with true distributions at the 5th&9thgroups signals

5 結(jié)束語

針對水下作戰(zhàn)裝備研發(fā)對水聲環(huán)境仿真需求,基于貝葉斯全概率公式構(gòu)建了水下目標回波雙程衰落信道模型,并結(jié)合針對實測時域數(shù)據(jù)的模型參數(shù)反演特點,設計了單維度串行采樣模型。模型通過M-H 核函數(shù)設計,將MCMC 方法引入全局反演算法中,以有效利用有限實測數(shù)據(jù)信息,為水下目標回波仿真模型置信度提升提供了有力而靈活的保證。文中使用仿真雙程信道數(shù)據(jù)對反演算法的有效性進行了驗證,并在最后使用多組實測數(shù)據(jù)開展反演,通過對反演結(jié)果的分析確認了雙程衰落信道模型的正確性。

值得注意的是,在使用實測數(shù)據(jù)開展模型參數(shù)反演時,除了信號SNR 條件會對反演效率和結(jié)果產(chǎn)生較大的影響外,準確的識別接收信號中脈沖回波所在空時區(qū)間是反演有效開展的重要保證,在算法的后續(xù)優(yōu)化中,可針對實測數(shù)據(jù)的自動選取及預處理開展更多工作,以進一步提升反演速度和結(jié)果的準確性。