利用聲壓場的淺海地聲參數貝葉斯反演方法

鄭廣學，樸勝春1,2,，祝捍皓，張海剛1,2,，李楠松1,2,

(1.哈爾濱工程大學 水聲技術重點實驗室，哈爾濱 150001；2.海洋信息獲取與安全工信部重點實驗室(哈爾濱工程大學)，工業和信息化部，哈爾濱 150001；3.哈爾濱工程大學 水聲工程學院，哈爾濱 150001；4.浙江海洋大學 海洋科學與技術學院,浙江 舟山 316022；5.浙江大學 浙江省海洋觀測-成像試驗區重點實驗室，浙江 舟山 316021；6.中國科學院聲學研究所 聲場聲信息國家重點實驗室，北京 100190)

海洋水體環境參數與地聲參數共同構成了水下聲傳播計算最重要的環境參數，相比水體環境參數，地聲參數更難直接測量，因而對其的有效獲取研究一直是水聲學中的熱點課題[1]。由于地聲參數的變化會對聲壓場的分布產生重要影響，因此可以考慮以水聽器接收到的復聲壓作為研究對象來反演地聲參數。相比于直接測量的方法，利用聲學方法反演可以快速、高效的獲取大面積海區的底質參數，避免了人力、物力的浪費，因此受到了廣泛的關注[2-5]。近年來，利用聲學方法反演地聲參數已取得了長足的發展。各種反演方法不斷涌現，如：利用傳播損失的地聲參數反演方法[6]；利用艦船輻射噪聲等各類海洋噪聲信息的地聲參數反演方法[7-10]；利用波導頻散特性的反演方法等。但上述地聲參數反演方法均著重研究反演問題中正演模型的選擇，多仍采用傳統尋優算法，如下降單純形法(downhill simplex，DS)、遺傳算法(genetic algorithm，GA)等對目標函數進行求解。但此類算法一般只能給出待反演參數的點估計，無法對參數的不確定性進行進一步描述，且在尋優過程中對算法初始模型依賴很強，容易陷入局部最優解。貝葉斯反演方法有效地估計最大后驗概率(maximum a posteriori，MAP)模型參數，進而從統計角度定性和定量分析參數反演結果的不確定性。

基于上述原因，本文以單水聽器接收到各距離點上的聲壓場為研究對象，基于貝葉斯理論，建立了針對淺海考慮彈性海底下的地聲參數反演方法，反演對象包括聲速、聲速衰減和密度3類海底聲學參數。

1 反演方法

根據聲場信息估計未知海底參數的聲學反演方法一般由4部分組成：1)構建聲傳播模型；2)建立合適的目標函數；3)全局優化算法；4)反演結果的不確定性分析。在本文的反演方法研究中，同樣針對上述4部分展開工作。

1.1 淺海聲場正演模型

淺海波導中，海底橫波聲速對低頻聲信號的傳播影響不可忽略[11]，本文研究中選擇了如圖1所示的淺海參數化模型。其中，z=0表征海面，z軸表深度，向下為正，r正軸表示聲場向外傳播方向。設水深為H、聲源位于水深zs處、密度ρ1、聲速c1；海底為半無限彈性介質，各向同性。縱波聲速、橫波聲速、密度、縱波聲速衰減和橫波聲速衰減分別用cp、cs、ρb、αp和αs表示，該5項海底參數也是本文反演時的研究對象。

圖1 單層彈性海底的淺海參數化模型Fig.1 Parametric model in shallow sea with elastic bottom

圖1模型下，流體層中的各位置聲壓p(r,z)為：

(1)

(2)

式中：ξ為水平波數；r為水平距離；J0為貝塞爾函數；文獻[12]給出了A、B、C、β等參數的具體推導過程及表示式。對式(1)中聲壓場p(r,z)的數值計算，通常采用簡正波方法(normal mode method，NMM)及快速場方法(fast field method，FFM)[13]。考慮此模型為淺海環境，FFM更適合用于聲場的快速計算，故本文選擇FFM完成對上述參數化模型下聲壓場p(r,z)的數值計算。

在完成對聲壓場p(r,z)數值計算的基礎上，需建立恰當的目標函數來衡量實測聲壓p(r,z)與理論預報聲壓p′(r,z)間的適配程度，并通過對其尋優求解便可實現該模型下各海底參數的獲取。

1.2 目標函數與優化算法

貝葉斯方法用于反演問題時，觀測數據向量d和模型參數向量m被視為隨機的。根據貝葉斯理論有:

P(m|d)=P(d|m)P(m)/P(d)

(3)

其中，P(m|d)為后驗概率密度(posterior probability density，PPD)。條件概率密度P(d|m)定義為似然函數L(m)，為某一(固定)測量數據d下m的函數，代表數據提供的信息，用于衡量在模型參數m下模型預報結果和數據的適配程度。P(m)為先驗概率密度，表示模型參數先驗信息，P(d)為測量數據的概率密度函數[14-15]。因此，后驗概率密度既包含數據信息，也包含模型參數的先驗信息。且P(d)與m無關，可以看作1個常數，所以式(3)可化為：

P(m|d)∝P(d|m)P(m)

(4)

根據貝葉斯理論，目標函數是在高斯數據誤差的假設下基于似然函數建立的，似然函數作為某一指定測量數據下模型參數的函數，代表了數據不確定性分布[14-15]。本文使用的目標函數為：

E(m)=ln[B(m)|p|2]

(5)

式中B(m)表示歸一化的Bartlett失配器：

(6)

式中：p表示單個水聽器接收到的聲壓場(實測場)，p′(m)為待反演參數向量m的預報聲壓(拷貝場)，上標“+”代表共軛轉置。

從式(6)可以看出，B(m)取得最優匹配時的條件為p(r,z)=p′(r,z)，同時，E(m)取得極小值。目標函數的建立，將對m的反演問題轉化為在尋優范圍內對式(5)的極小值求解問題。根據貝葉斯理論，上述待反演參數的PPD將依據Metropolis準則采樣[16]給出，并從統計角度對其進行不確定性分析。

2 聲場正演模型仿真分析

本文在完成了聲場建模、確定了目標函數及尋優算法的基礎上，利用數值模擬，對比圖1模型下各海底參數與其反演結果，對研究方法進行驗證。

2.1 聲傳播特性分析

圖2給出了利用FFM計算所得圖1所示淺海環境下，1～1 500 m 100 Hz聲信號聲壓場的仿真結果，以傳播損失(transmission loss，TL)曲線的形式給出。

圖2 5項參數對TL的影響程度分析Fig.2 The impact of five parameters on transmission loss

仿真中，設聲源深度為zs=20 m，接收點zr位于水下10 m處，表1給出了仿真中各海洋環境參數。圖2(a)～(e)中的實線對應利用表1中各參數真值的計算結果，5項地聲參數取討論值時的計算結果分別用點線(較小值)和虛線(較大值)表示。圖2(f)中應用“距平”的定義[17]，視設定的環境參數真值計算所得傳播損失為均值，視各參數的討論值計算所得傳播損失為離散值來計算距平，用5項地聲參數討論值的距平大小來衡量各參數變化對傳播損失的影響程度。

表1中各參數的討論值均設定為真值偏移±10%，在此基礎上討論各參數變化對聲壓場的影響程度。從圖2(a)～(e)的中可以看出，5項參數的變化對TL曲線影響程度各不相同，從圖2(f)中5項參數討論值的距平大小可以看出，在各參數均偏移±10%的情況下，改變某一參數對TL曲線的影響程度從大到小依次為cp、cs、ρb、αp、αs。因此，5項海底參數對p(r,z)的影響程度可初步定性為：cp、cs>ρb>αp、αs。

2.2 目標函數對各參數的敏感度分析

在完成了正演建模、聲壓場數值計算的基礎上，本文通過數值模擬仿真分析本文所給目標函數對5項地聲參數的敏感度，并討論該目標函數在反演中的可靠性。仿真中，海洋環境參數真值與各參數尋優區間仍取表1中數值，取各參數設定真值下模擬計算得到p(r,z)為參考真值。本文采用固定變量法來分析目標函數對單個地聲參數的敏感度時，即其他參數固定，只在討論區間內改變某一參數計算聲壓場理論預報值p′(r,z)，由式(5)計算p(r,z)與p′(r,z)間的目標函數值E(m) 作為敏感度評價標準。

表1 海洋環境仿真參數Table 1 The simulation parameters of ocean environment

圖3(a)～(e)分別給出了目標函數E(m)隨表1中各單項參數變化時的數值變化規律。從圖3可以看出，在5項地聲參數討論范圍內，目標函數E(m)均只在其真值處取得極小值，且無偽峰或第2極小值，避免了局部最優解的陷阱，證明了該目標函數的可行性。但從圖3(f)可知目標函數隨各參數變化時的波動程度不一致，在5項參數的討論范圍內，目標函數值波動范圍從大到小依次為cp、cs、ρb、αp、αs。因此，5項海底參數對p(r,z)的影響程度可定性為：cp>cs>ρb>αp>αs，進一步證明了圖2的討論結果。

圖3 單參數敏感度Fig.3 The sensitivity of single parameter

2.3 仿真反演結果分析

在完成目標函數對5項海底聲學參數敏感度的研究基礎上，本節通過數值仿真，分析本文所研究的反演方法對目標函數尋優結果的可靠性。海洋環境參數如表1所設，當尋優結束時，各參數尋優結果與仿真真值的對比如圖4給出。圖中水平線段表示了各參數PPD的均值位置及其方差大小，垂直于橫坐標直線段標注仿真真值所在位置，各參數反演結果在表2中給出。

圖4 5項參數及其概率密度分布Fig.4 Five parameters and their probability density distribution

表2 尋優結果Table 2 The results of optimization

從圖4和表2中的結果可以看出，各參數的仿真真值與反演結果吻合程度較高，均處于其尋優結果的均值與方差的和或差范圍內。且各參數的概率密度峰值的尖銳程度也反映了目標函數對各參數的敏感程度，如圖4所示，目標函數對各參數的敏感程度為：cp、cs>ρb>αp、αs，這也與2.1節、2.2節的討論結果相吻合。

為驗證數值仿真中反演結果的可靠性，圖5給出了仿真環境參數下TL曲線與仿真反演結果計算所得TL曲線的對比圖。等間距取12處繪制誤差棒，并標注反演結果計算所得TL曲線的誤差范圍。從圖中可以看出，仿真中TL曲線均處于反推TL曲線誤差范圍之內，證明了該反演方法在數值仿真應用中的可靠性。

圖5 TL曲線對比Fig.5 The comparison of TL

3 消聲水池縮比實驗驗證

3.1 水池實驗

本文結合消聲水池縮比實驗[18]，利用實驗數據對本文所研究反演方法進行討論。實驗中選取質地均勻、高硬度的塑料板模擬“半無限彈性海底”。

表3給出了實驗中各設備的布放參數，其中聲源深度zs為87 mm，接收深度zr為84 mm，水深z為182 mm，水中聲速c1為1 450.212 m/s由測量水中溫度(11.15 ℃)后通過經驗公式求得。測量過程中，聲源固定不動，發射信號為135 kHz的脈沖信號；接收水聽器固定在可移動走架上，采集卡采樣率fs為20 MHz，單個測量點分別測量10次以減少隨機擾動帶來的誤差；完成一點測量后，走架帶動水聽器向遠離聲源方向移動2 mm(誤差不超過20 μm)至下1個測量點，實驗中測量點共720個。測量中TL曲線如圖6(b)所示，圖6(a) 給出了測量中所接收到的第1～100道接收信號時域圖。

圖6 實驗中接收到的信號Fig.6 The signals received in experiment

3.2 反演結果驗證與分析

表3列出了各地聲參數的尋優范圍，以及利用本文所研究反演方法對上述實驗數據反演所得到的“海底”地聲參數值，同時還給出了利用GA對該實驗數據尋優的結果作為對比。

表3 尋優結果Table 3 The results of optimization

圖7給出了算法終止時，各參數在尋優區間上的分布及其概率密度分布，與圖4相仿，圖7中垂直于橫坐標直線段標注GA尋優結果所在位置，水平線段表示了各參數取值PPD的均值位置及其方差大小。圖8定量給出了參數間的相關系數矩陣圖。

從表4、圖7和圖8的結果中可以看到：首先，各待反演參數的GA尋優結果均處于均值與方差的和或差范圍內，證明了利用貝葉斯反演方法獲得的結果與傳統反演方法得到結果相吻合。另外，方差的相對大小從一方面也反映出該參數反演結果不確定性的大小。其次敏感度較小的參數對目標函數貢獻也較小，在圖8中表現為該參數與其他參數的相關性偏小，如αp和αs；換言之，較敏感參數間則表現出較明顯的相關性，如cp、cs、ρb，綜上，可以看出5項待反演參數的敏感度大小分別為cp、cs>ρb>αp、αs，這也與第2節結論一致。

圖7 算法終止時5項參數及其概率密度分布Fig.7 The five parameters and their probability density distribution at the termination of the algorithm

圖8 參數間相關矩陣Fig.8 The correlation matrix of parameters

為進一步驗證表4反演結果的可信度，圖9對比了利用反演結果計算得到的TL曲線和實測TL曲線；圖10對比了同一測量點(第470號接收點)，反演所得時域波形與該點實際接收信號的時域波形圖(均進行歸一化處理)。圖9中，實線為實測TL曲線，點線為應用表4反演結果計算所得的TL曲線并等間距取12處繪制誤差棒，從圖中對比結果可以看到實驗中TL曲線基本處于反推TL曲線誤差范圍之內，證明了該反演方法在實際應用中的可行性；圖10中，實線為第470號接收點實測時域波形，點線為應用表4反演結果“預報”的第470號接收點時域波形，從圖中對比結果可以看到時域波形計算預報結果與實測信號波形匹配較好，若不考慮信號“拖尾”，相關系數可達0.89，證明了本反演方法在實際應用中的可靠性。

圖9 實測TL與反演所得TL對比Fig.9 The comparison of TL between measured and inversed

圖10 時域波形對比Fig.10 The comparison of time domain waveform

4 結論

1)貝葉斯反演能夠從統計角度給出反映地聲參數向量綜合信息的后驗概率來表征反演結果，同時還能夠分析各參數反演結果的不確定性，相比傳統只能給出待反演參數點估計的反演方法，反演信息更加全面、科學。

2)根據貝葉斯反演理論，利用本文給出的目標函數，運用Metropolis準則采樣，可以較為準確得到待反演參數的PPD。數值模擬和水池實驗的結果均表明：PPD分布特性與GA尋優結果吻合程度較好，利用反演參數計算得到的聲場特性與測量場基本一致，時域波形相關系數達0.89。

3)從5項地聲參數反演結果的不確定性大小來看，在仿真模擬及水池實驗的環境下，縱波聲速cp、橫波聲速cs以及密度ρb的不確定性相對較小，對聲壓場的變化更敏感，5項參數依據敏感程度大小依次為cp>cs>ρb>αp>αs，聲場研究結果與反演結果均證實了這一研究結果。