浅海
参数
变维贝叶斯
反演
吴伟文
浅海地声参数的变维贝叶斯反演*吴伟文1,2任群言1,2鹿力成1,2郭圣明1,2马力1,2(1中国科学院水声环境特性重点实验室(声学研究所)北京100190)(2中国科学院大学北京100049)2022 年 6 月 6 日收到2022 年 11 月 2 日定稿摘要为解决海底沉积层分层结构未知时的反演问题,提出了一种变维粒子滤波方法,利用声场的互谱密度,估计沉积层分层结构以及地声参数。仿真结果表明:变维粒子滤波在沉积层分层结构未知时,能有效反演沉积层层数以及地声参数,粒子的并行计算能使其相较于可逆跳蒙特卡洛马尔可夫链(rjMCMC)更加高效。利用变维粒子滤波,对南海垂直线阵列接收到的线性调频信号进行处理,反演结果与 rjMCMC 反演得到的沉积层层数和地声参数结果相近,说明了此方法能有效估计沉积层层数的同时反演浅海地声参数,得到可靠的参数后验概率密度。关键词地声反演,地声模型,贝叶斯理论,不确定性分析PACS 数数43.30,43.60Trans-dimensional Bayesian geoacoustic inversion in shallow waterWUWeiwen1,2RENQunyan1,2LULicheng1,2GUOShengming1,2MALi1,2(1Key Laboratory of Underwater Acoustic Environment,Institute of Acoustics,Chinese Academy of SciencesBeijing100190)(2University of Chinese Academy of SciencesBeijing100049)ReceivedJun.6,2022RevisedNov.2,2022Abstract Todealwiththeinversionproblemwhenthespatialstructureofseafloorsedimentsisunknown,atrans-dimensionalparticlefiltermethodisproposedinthispaper,wherethecross-spectraldensityofthepressurefieldisusedtoestimatethesedimentlayeringstructureandgeoacousticparameters.Thesimulationresultsshowthatthenumberofsedimentlayersandgeoacousticparameterscanbeeffectivelyestimatedusingtheproposedmethod,andtheparallelparticlecalculationmakethismethodmoreefficientthanreversiblejumpMonteCarloMarkovchain(rjMCMC)sampling.Alinearfrequency-modulatedsignalreceivedbyaverticallinearrayintheSouthChinaSeaisprocessedusingtheproposedtrans-dimensionalparticlefilter.TheinversionresultsofthesedimentlayersandgeoacousticparametersaresimilartothoseacquiredbyrjMCMCinversion.Thenumberofsedimentlayersandtheposteriorprobabilitydensityofparameterscanbeeffectivelyestimatedusingthismethod.KeywordsGeoacousticinversion,Geoacousticmodel,Bayesiantheory,Uncertaintyanalysis引言声信号在水中具有良好的传播性能,常用来进行水下目标探测或信息传递。浅海环境中,声传播以及匹配场定位性能等受地声参数及其空间分布的影响较大,浅海地声参数的估计是海洋声学中的一个基础研究问题1-2。地声参数可以通过直接测量3获得,但此方法时间成本和经济代价都很高昂。匹配场反演(MFI)4是比直接测量法更有效率的方法。MFI 将地声参数的获取作为优化问题,通过目标函数衡量观测数据*国家自然科学基金项目(12274431)和中国科学院声学研究所自主部署项目(MBDX202102)资助通讯作者:任群言,第48卷第3期声学学报Vol.48,No.32023年5月ACTAACUSTICAMay,2023DOI:10.15949/ki.0371-0025.2023.03.017与拷贝场数据之间的误差,目的就是要在选定的模型中,通过搜索算法寻找一组地声参数向量,最小化目标函数。目标函数中待匹配的数据,通常是对地声参数变化较敏感的物理量。声压互谱密度5、声压时频谱、脉冲信号时延等随地声参数的变化较为明显,常作为目标函数中待匹配的物理量。由于参数反演问题具有维度高,前向模型非线性等特点,参数搜索算法一般选用遗传算法(GA)6、模拟退火(SA)6等全局优化算法。近年来,为了量化分析反演结果的不确定性,贝叶斯方法被广泛应用于参数反演中7-10,用后验概率密度(PPD)来表示反演结果及其不确定性。Gerstoft 等7首次利用贝叶斯理论表征多维参数的边缘概率密度,并根据参数的后验概率密度解释反演结果。Bo 等11用自适应粒子滤波(APF)反演浅海地声参数,分析反演结果的不确定性,其利用船自噪声作为声源,拖曳水平阵作为接收,取得了较好的反演效果。实际海域中,海底沉积层分层结构的先验信息可能不足,如沉积层层数未知。这导致了使用固定模型的反演方式不能较好反应海底的真实分层结构。近场声传播往往受海底分层结构影响较大,错误的海底分层模型会导致模型失配,存在过拟合或欠拟合现象,无法得到符合实际的地声反演结果。Battle 等12利用贝叶斯方法研究了贝叶斯反演中模型选择的问题,分析了不同模型的选择对反演结果的影响。Dettmer 等13-14使用 rjMCMC 采样解决了沉积层层数未知时的地声反演问题。但 rjMCMC 计算量巨大且每条马尔可夫链无法并行计算,难以满足实时性要求较高的系统。本文提出一种变维粒子滤波(TDPF)方法,利用垂直线阵接收到的线性调频信号的多个频点,进行浅海地声参数贝叶斯反演。该方法可以在滤波过程中快速实现模型的选择,同时获得参数的后验概率密度。仿真与实验结果表明,TDPF 能较好地反演出沉积层分层信息,并获得可靠的参数后验概率密度。沉积层层数已知时,TDPF 与使用遗传算法的反演结果作比较,两者获得了相近的参数反演结果。分层未知时,对于地声参数的估计,TDPF 取得了与rjMCMC 相近的反演效果,此外,相较于 rjMCMC 采样,TDPF 每一次循环中都可对所有粒子同时计算似然函数的大小,而 rjMCMC 链中,样本需按照采样顺序逐个计算似然函数以判断是否接收新样本。此方法相较于 rjMCMC 具有较好的实时性,能更快地反演地声参数,反演效率较高。1 反演模型 1.1反演参数模型z0rH1H2c1c21212c333浅海环境中,垂直线阵接收声源发射的线性调频信号,包含多条海底海面反射路径。接收阵与声源距离较近,可以假设此段距离内地声参数近似不变。采用图 1 所示的声传播模型进行反演,其中沉积层层数未知导致地声模型未知。模型参数包括声源参数、阵列参数以及地声参数。声源为换能器发射的线性调频信号,深度表示为,与垂直阵的距离表示为。声源深度实验中有较为准确的监测记录。水体参数包括水深 D 和声速剖面(SSP)参数。实验中,水深由船载测深仪大致获得,声速剖面由CTD 测量获得,并认为实验过程中声速剖面较为稳定,声速剖面的微小变化对反演结果的影响不大5,故本文未反演声速剖面。海底分为沉积层和基底层,其中沉积层层数未知,假设两层沉积层时,层厚分别为,每层沉积层声速、密度和衰减恒定,分别为,基底声速、密度和衰减为,。阵列参数包括阵列深度和倾斜角度,其中阵列深度有可靠的测量结果,同时阵列深度不是敏感的物理量,阵列参数简化为只包含倾斜角度。考虑近场声传播存在泄露简正波,本文使用考虑泄露简正波的 Krakenc 计算近场声传播。LFM 声源z0DH1H2c1,1,1c2,2,2c3,3,3r垂直阵海水沉积层基底SSP图1MFI 环境模型 1.2数据模型及似然函数Nx假设垂直阵包含个全向阵元,在参数向量为时接收到的窄带复声压数据可以表示为y(fj)=(fj)d(fj,x)+v(fj),(1)fjxd(fj,x)(fj)v(fj)其中,为第 j 个频率点,为参数向量,为信道响应,为复声压,为噪声,并假设噪声为高斯噪声。标准化的 Bartlett 失配度为11B(fj,x)=1dH(fj,x)?R(fj)d(fj,x)d(fj,x)2Tr(?R(fj),(2)3期吴伟文等:浅海地声参数的变维贝叶斯反演497?R(fj)=y(fj)yH(fj)其中,为互谱密度。J假设不同频率的数据误差不相关,给出个频率的联合似然函数为L(xt)=Jj=1|N?V(fj)Tr?R(fj)exp|B(fj,xt)?V(fj)|N,(3)NJ其中,为传声器阵列阵元个数,为频点个数,?V(fj)=B(fj,?x),?x?V(fj)其中,是参数向量的最优估计,与有效信噪比的关系为?ESNR=10lg|1?V(fj)?V(fj)|.(4)NNNN当较大时,此似然函数较窄,当较小时,似然函数较宽。的大小的改变会影响似然函数的宽度,也会影响粒子滤波的参数后验概率密度。本文采用文献 12 中启发式的。2 反演原理及实现 2.1贝叶斯反演原理声源距离接收阵列较近,在此范围,可将参数近似认为不变的。为了获得更准确的参数,需要一定数量的粒子以计算参数后验概率密度,可建立状态方程和观测方程分别为xt=f(xt1)+nt,(5)yt=g(xt)+vt,(6)xttntvtytxtvt其中,为第 次迭代的环境参数,包含声源深度、声源距离、阵列倾斜角、水深、沉积层声速、密度、衰减、层厚、基底声速、密度和衰减。为状态噪声,为观测噪声。为参数为时计算得到的声场加上观测噪声后的声场。xt=xt1+nt地声反演中,对于参数的变化规律没有足够的先验知识的情况下,状态方程一般可以设置为。贝叶斯方法中,参数后验概率密度等价于先验与似然的乘积:(xt)=L(xt)(xt),(7)(xt)L(xt)(xt)其中,为参数后验概率密度,为似然函数,为参数先验概率密度。2.2粒子滤波原理(xt)p(x1:t1|y1:t1)p(xt|xt1)L(xt)yt粒子滤波(PF)作为一种序贯贝叶斯方法,对所有粒子而言,即为参数的先验概率密度,即为加入新数据后的p(yt|xt)(xt)p(x1:t|y1:t)似然函数,粒子后验概率密度即为。粒子滤波目标分布(即参数后验概率密度)的递推式为p(x1:t|y1:t)=p(x1:t1|y1:t1)p(yt|xt)p(xt|xt1)p(yt|y1:t1),(8)p(yt|y1:t1)其中,为归一化系数。传统粒子滤波目标分布用粒子形式表示5:P(x1:t|y1:t)Ni=1w(i)t(x1:tx(i)1:t),(9)w(i)tNi=1w(i)t=1其中,为粒子权重,。q(x1:t)粒子以重要性概率密度采样,可得:w(i)tp(x(i)1:t1|y1:t1)p(yt|x(i)t)p(x(i)t|x(i)t1)q(x(i)1:t)w(i)t1p(yt|x(i)t)p(x(i)t|x(i)t1)q(x(i)1:t).(10)q(x1:t)令重要性概率密度函数:q(x1:t)=p(xt|xt1).(11)则粒子权重简化为w(i)t w(i)t1 p(yt|x(i)t).(12)2