袁藝，姚華建*，秦巖

1 中國科學技術(shù)大學地球和空間科學學院，地震與地球內(nèi)部物理實驗室,合肥　230026　2　蒙城地球物理國家野外科學觀測研究站,安徽 亳州　233527

?

基于鄰域算法的瑞利面波垂直-水平振幅比及頻散曲線聯(lián)合反演及應用

袁藝1,2，姚華建1,2*，秦巖1,2

1 中國科學技術(shù)大學地球和空間科學學院，地震與地球內(nèi)部物理實驗室,合肥2300262蒙城地球物理國家野外科學觀測研究站,安徽 亳州233527

摘要瑞利面波垂直-水平振幅比(或ZH振幅比)是一個隨頻率變化的函數(shù)，對于臺站下方淺層地殼結(jié)構(gòu)非常敏感，且具有和頻散資料不同的深度敏感核，是傳統(tǒng)頻散反演方法的一個很好的補充，從而可以將基階瑞利面波的ZH振幅比和面波頻散數(shù)據(jù)聯(lián)合起來更好地反演獲得觀測臺站下方的速度結(jié)構(gòu).本文提出了基于鄰域算法的面波頻散曲線與ZH振幅比聯(lián)合反演方法，我們進行了基于理論模型的模擬測試，證明了聯(lián)合反演是一種更為可靠的反演方法，且能更好地約束淺層地殼結(jié)構(gòu).相比于頻散曲線單獨反演，聯(lián)合反演不僅可以精確反演獲得地殼的Vs結(jié)構(gòu)，對分層地殼的Vp/Vs也能很好地約束.然后我們將聯(lián)合反演算法應用于實際測量數(shù)據(jù)，獲得了中國西南昆明臺(KMI)下方更為準確的地殼橫波速度結(jié)構(gòu)及Vp/Vs模型.

關(guān)鍵詞聯(lián)合反演； 鄰域算法； 瑞利面波； 相速度頻散； ZH振幅比； 地殼結(jié)構(gòu)

1引言

利用地震面波和背景噪聲數(shù)據(jù)測量得到的面波頻散數(shù)據(jù)是反演地下結(jié)構(gòu)一種十分重要且日趨成熟的方法(Shapiro et al.,2005；Yao et al.,2008)，但是這種方法存在一定的缺陷：因為某一周期的面波相速度或群速度對地下一定深度范圍的結(jié)構(gòu)比較敏感，周期越長，敏感深度越大，但敏感深度范圍越寬(徐果明等，2007)，從而導致深部的反演結(jié)果的分辨率降低；且由于高質(zhì)量的短周期頻散數(shù)據(jù)較難獲得，所以通常很難準確獲得地殼淺部幾公里的速度結(jié)構(gòu).如果能提高對地表幾百米到幾千米的速度結(jié)構(gòu)的分辨率，我們將能更好地關(guān)聯(lián)地表的地質(zhì)構(gòu)造特征和地殼深部的結(jié)構(gòu)異常及演化過程，且對認識和預測地震地質(zhì)災害及礦產(chǎn)油氣資源勘查都具有重要的意義.

基階瑞利面波在傳播的過程中徑向(H)和垂向(Z)的振動分量具有90°的相位差，其垂向和徑向分量的振幅比值與地下淺部結(jié)構(gòu)有密切的關(guān)系，我們把這一參數(shù)稱為ZH振幅比.ZH振幅比是一個隨頻率變化的函數(shù)，對于橫向各向同性，縱向隨深度變化的分層地球模型來說，這一函數(shù)只與基階瑞利面波所到達位置的地下介質(zhì)力學參數(shù)有關(guān)，與其他參數(shù)(如震源)無關(guān)(Tanimoto and Alvizuri，2006).隨著寬頻帶地震臺的增加及數(shù)據(jù)測量、處理技術(shù)的提升，關(guān)于ZH比的研究成果也越來越多.Tanimoto和Alvizuri(2006)利用基階瑞利面波徑向和垂向分量之間的90°相位差從連續(xù)的噪聲信號中搜尋基階瑞利面波，計算其ZH振幅比并反演地下淺層橫波速度結(jié)構(gòu).Tanimoto和Rivera(2008)利用數(shù)值方法得出了長周期(20～250 s)面波ZH振幅比的深度敏感核曲線，相比于基階瑞利面波相速度頻散曲線，ZH比對淺層結(jié)構(gòu)更加敏感.Yano等(2009)利用中長周期瑞利面波的ZH振幅比反演了GEOSCOPE臺網(wǎng)不同臺站下方的地殼和上地幔結(jié)構(gòu).Tanimoto等(2013)發(fā)展了從連續(xù)噪聲信號中計算ZH比的新方法.

面波頻散曲線和ZH振幅比都可以用來反演地下速度結(jié)構(gòu)，但是同樣的頻率范圍兩者的敏感深度范圍不同，ZH振幅比正好可以彌補頻散曲線對地殼淺層結(jié)構(gòu)約束不足的缺陷.Boore和Toks?z(1969)很早就做了頻散曲線和ZH振幅比聯(lián)合反演的試驗.Lin等(2012)將中長周期(24～100 s)的ZH振幅比和8～100 s的頻散曲線結(jié)合起來用于反演美國西部地區(qū)地殼上地幔的結(jié)構(gòu).Lin等(2014)利用三分量背景噪聲的互相關(guān)數(shù)據(jù)計算了臺站下方中短周期的ZH振幅比(8～24 s)，和之前從地震數(shù)據(jù)中得到的24～100 s ZH振幅比及8～100 s相速度頻散結(jié)果結(jié)合，可以更好地用于聯(lián)合反演美國西部地區(qū)地殼及上地幔的三維結(jié)構(gòu).上述這些ZH比和頻散曲線的聯(lián)合反演均是基于線性最小二乘的反演方法，需要計算ZH比和頻散數(shù)據(jù)的敏感核.Chong等(2015)提出了基于模擬退火算法的面波頻散曲線及ZH振幅比聯(lián)合反演，直接在模型空間進行隨機搜索從而獲得全局最優(yōu)模型，所以不需要計算敏感核.

本文提出一種基于鄰域算法(Neighborhood Algorithm,或NA)(Sambridge，1999a，1999b)的面波頻散曲線和ZH振幅比的聯(lián)合反演，并通過貝葉斯分析方法獲得模型反演的誤差、分辨率及參數(shù)之間的相關(guān)性.我們首先通過模擬測試數(shù)據(jù)證明基于鄰域算法的聯(lián)合反演的可靠性，接下來將聯(lián)合反演算法應用于實際地震數(shù)據(jù).

2瑞利面波ZH振幅比及頻散的深度敏感核

我們定義ZH振幅比為基階瑞利面波在地表處(z=0時)垂直向與水平向振幅的比值，即

(1)

其中，ω代表頻率，|r1(ω|z)|和|r2(ω|z)|為對應于頻率ω在不同深度z處的瑞利面波本征函數(shù)值(Aki and Richards，1980).

(2)

η為ZH振幅比(或相速度c)；ρ為密度，α為縱波波速，β為橫波波速，δp(p=ρ,α,β)分別為ZH振幅比(或相速度)、密度、縱波波速、橫波波速的擾動，這些參數(shù)均為深度z的函數(shù)；Kp被稱為敏感核，是頻率和深度的函數(shù).

圖1給出了一維平層速度結(jié)構(gòu)下(地殼厚度為45 km)3個不同周期瑞利面波ZH振幅比及相速度對結(jié)構(gòu)參數(shù)的深度敏感核，可以清楚地觀察到ZH振幅比和相速度具有非常不同的深度敏感核，尤其是對橫波速度.與相速度相比，ZH振幅比的敏感深度更淺，約為頻散曲線的一半(Tanimoto and Rivera,2008)，而且ZH振幅比在地表處的敏感性最強，因此在相同頻帶下能反演比頻散曲線更淺深度的橫波速度結(jié)構(gòu).瑞利波相速度頻散對地殼的縱波速度及密度有一定的敏感性，ZH比對地殼(淺部)密度也有一定的敏感性，但對縱波速度的敏感性很低.

圖1　不同周期基階瑞利面波相速度(a)和ZH振幅比(b)對縱橫波速和密度的深度敏感核

3ZH振幅比及頻散的聯(lián)合反演方法

因為ZH振幅比具有和面波頻散不一樣的深度敏感核，且它對淺部結(jié)構(gòu)的約束性更強，所以我們可以將其與面波頻散曲線結(jié)合進行聯(lián)合反演，以更好地約束地殼的速度結(jié)構(gòu)(Lin et al.,2012).本文中，我們發(fā)展了一套基于鄰域算法(NA)(Sambridge,1999a,1999b)的ZH振幅比及頻散曲線聯(lián)合反演算法.

3.1鄰域算法(NA)

鄰域算法是一種基于全局搜索的反演方法，它包括兩個階段：第一階段應用“Voronoi cells”(維諾單元)的特殊幾何分割方法在模型空間中搜索可接受的(“好”的)模型(Sambridge,1999a)，第二階段對第一階段生成的模型進行貝葉斯統(tǒng)計分析.得到的模型參數(shù)的一維邊際后驗概率分布和高斯函數(shù)分布的偏離程度可以作為反演結(jié)果可靠與否的判斷標準，且以一維概率分布的寬度(即標準差)來表征模型的不確定性.二維邊際后驗概率密度函數(shù)可以用于衡量兩個模型參數(shù)之間的相互折衷性(trade-off)，其分布越集中越接近圓形，說明二者的相互折衷性越小.一維和二維邊際概率分布函數(shù)都是估計分辨率和模型參數(shù)不確定性的有力工具(Sambridge,1999b).

3.2理論模型反演

為了驗證頻散曲線與ZH振幅比聯(lián)合反演相對于頻散曲線單獨反演地殼結(jié)構(gòu)的優(yōu)越性，我們先用理論模型正演產(chǎn)生的數(shù)據(jù)來進行NA反演測試.

3.2.1理論模型構(gòu)建及數(shù)據(jù)準備

我們構(gòu)建的理論模型包括地殼和上地幔兩部分，上地幔310 km深度之下設(shè)為半無限空間.地殼厚度設(shè)為45 km，分為四層，第一層的縱橫波速比(Vp/Vs)為1.9，其余三層均為1.7.我們先設(shè)定了地殼中四層的橫波速度Vs，根據(jù)Vp/Vs值計算相應層的縱波波速Vp，再根據(jù)經(jīng)驗公式(Brocher,2005)可得密度值；上地幔部分從45 km到310 km，來源于ak135模型(Kennett et al.,1995)，如表1所示.根據(jù)此模型我們采用Herrmann和Ammon(2004)的方法正演得到瑞利面波相速度頻散數(shù)據(jù)31個(周期范圍為5～150 s)以及ZH振幅比數(shù)據(jù)26個(周期范圍為5～120 s)，這些數(shù)據(jù)(未加噪聲)將用于下面的NA測試反演之中.

3.2.2NA聯(lián)合反演

表1中給出我們采用NA算法擬反演的10個參數(shù)，包括4層地殼和4層上地幔的橫波速度，以及地殼第一層(0～6 km)和其余三層(6～45 km)的縱橫波速比(Vp/Vs).對于橫波速度結(jié)構(gòu)，我們以理論模型為參考值，對每個參數(shù)給定相應的變化范圍(如表1所示).這里上地幔的各層包含ak135模型中的若干厚度較薄的層，在反演中模型參數(shù)值將在參考模型值的基礎(chǔ)上整體變化.由于面波數(shù)據(jù)對間斷面約束較差，且和相鄰層的速度結(jié)構(gòu)有很大的折衷性(例如Yao et al.,2008)，所以我們通常不反演層厚(或間斷面的位置).

NA算法中將產(chǎn)生各層的橫波速度，我們將根據(jù)縱橫波速比(如果也作為模型參數(shù)的話)或經(jīng)驗公式(Brocher,2005)計算出地殼部分相應層的縱波速度，之后再根據(jù)經(jīng)驗公式(Brocher,2005)從縱波速度計算出密度值.但是該經(jīng)驗公式不適用于上地幔部分，所以我們采用如下另一經(jīng)驗關(guān)系式來計算縱波速度和密度相對于橫波速度的變化(Masters et al.,2000；Yao et al.,2008)：

(3)

反演中誤差函數(shù)采用了L1范數(shù)的形式：

(4)

考慮到反演結(jié)果的可靠性以及收斂速度，經(jīng)過反復試驗，我們在每一次迭代中將產(chǎn)生20個新模型(Ns=20)，結(jié)合先前生成的模型，每次迭代時模型空間將被分割生成新的維諾單元圖，接著從所有的維諾單元中重新選出10個最好的單元(Nr=10)再生成20個新模型，如此反復迭代，直至搜索過程收斂.之后我們對生成的模型進行貝葉斯統(tǒng)計分析獲得每個模型參數(shù)的后驗平均值(posterior mean)及其標準誤差(對應于68%的置信區(qū)間).

圖2—4顯示了我們直接反演10個參數(shù)(見表1)的結(jié)果.圖2黑色細實線給出了后驗平均模型，而虛線給出了搜索得到的最優(yōu)模型，最優(yōu)模型和理論 模型符合得非常好，但后驗平均模型和理論模型在部分參數(shù)上存在一定的偏差.由于Moho面附近的第4層(下地殼)和第5層(上地幔頂部)的橫波速度之間存在相互折衷(如圖4的二維后驗概率密度函數(shù)所示)，導致下地殼后驗平均橫波速度較真實模型偏高約0.12 km·s-1，而上地幔頂部后驗橫波速度則偏低約0.18 km·s-1.

表1　鄰域算法反演的10個參數(shù)的參考值及其變化范圍

圖2一步聯(lián)合反演得到地殼及上地幔的橫波速度結(jié)構(gòu)(左)，由統(tǒng)計后驗平均模型和最優(yōu)模型正演得到的相速度頻散曲線(中)和ZH振幅比(右)與理論計算值的比較左圖中粗實線、細實線、點化線分別代筆理論模型、后驗平均模型和最優(yōu)模型，灰色區(qū)域給出后驗平均模型的標準偏差；中間及右圖中三角形代表由后驗平均模型正演出的頻散數(shù)據(jù)及ZH振幅比，星號代表最優(yōu)值，菱形代表理論模型計算值.

圖3　頻散及ZH比數(shù)據(jù)同時反演地殼及上地幔10個參數(shù)的一維后驗概率密度函數(shù)分布(陰影區(qū)域所示)橫軸為橫波速度或縱橫波速比的變化范圍，縱軸代表歸一化的后驗概率密度，每幅圖中的黑色豎直實線分別對應于后驗平均模型.對于橫波速度參數(shù)(上兩行)，橫坐標0值處對應于真實模型的值.

圖4　頻散及ZH比數(shù)據(jù)同時反演地殼及上地幔10個參數(shù)的二維后驗概率密度函數(shù)分布黑線、藍線和紅線分別代表99%、90%和60%置信級的等值線，白色三角形代表后驗平均模型.

由圖3可見，雖然后驗概率密度函數(shù)呈高斯函數(shù)分布，且大部分模型參數(shù)的后驗平均值與真實模型值也比較相符，但模型后驗函數(shù)寬度較大，說明反演模型的標準誤差偏大.由于NA算法的全局搜索參數(shù)較多時，各參數(shù)之間的相互折衷性較大，可能會導致反演結(jié)果的不穩(wěn)定，從而使得模型搜索落入到區(qū)域極小值附近，使得最終得到的后驗平均模型與真實模型存在一定的偏差.

因此，我們嘗試采用兩步反演的模式，第一步僅反演地殼上地幔若干層的橫波速度結(jié)構(gòu)，這里為8個參數(shù)，分別為4層地殼和4層上地幔的橫波速度(見表1的前8個參數(shù))，以求將各層的橫波速度初步確定.獲得較為可靠的上地幔模型后，我們進行第二步反演，共反演7個參數(shù)，分別為4層地殼(同第一步分層)和1層上地幔(45～310 km)的橫波速度Vs，以及地殼第一層(0～6 km)和其余三層(6～45 km)的縱橫波速比(Vp/Vs).第二步的參考模型選取第一步反演獲得的地殼上地幔橫波速度結(jié)構(gòu)的后驗平均模型值，并將上地幔的反演參數(shù)由4個變成1個，即Moho-310 km深度的上地幔橫波速度在第二次反演時僅相對于第一次反演獲得的該深度范圍的橫波模型做整體的擾動變化.通過反演參數(shù)的變化，我們希望更好地約束地殼橫波速度結(jié)構(gòu)及第一步反演中未涉及的Vp/Vs值.

兩步NA聯(lián)合反演結(jié)果如圖5—7所示.圖5給出NA反演所得后驗平均模型值與理論模型的對比圖，兩模型非常接近，且相比于一步反演的結(jié)果，最終反演模型的誤差很小.另外，面波頻散和ZH比的后驗平均模型計算值與理論值幾乎重合(圖5)，這個結(jié)果證明了反演所得后驗平均模型的正確性，驗證了NA聯(lián)合反演算法的可靠性和精確度，結(jié)果也明顯優(yōu)于單步反演法(圖2—4).

圖6的陰影部分給出聯(lián)合反演所得模型的一維邊際后驗概率密度函數(shù)分布，各模型參數(shù)的后驗平均值如圖中豎直黑線所示，最大似然模型可以從這些一維分布圖的峰值得到(Sambridge,1999b).相比于一步反演的結(jié)果(見圖3)，圖6顯示各反演參數(shù)的一維邊際后驗概率密度函數(shù)分布基本都呈高斯分布，且范圍較窄，顯示出模型的標準誤差小，分辨能力強.除了橫波速度Vs，縱橫波速比Vp/Vs也與理論設(shè)定值相符.6 km-Moho深度的Vp/Vs值的標準誤差較0～6 km深度Vp/Vs值的標準誤差明顯要大，顯示出聯(lián)合反演對地殼淺部的Vp/Vs值有 更強的約束.

圖5　兩步聯(lián)合反演得到地殼及上地幔的橫波速度后驗平均模型(實線)及理論模型(虛線)(左)，由后驗平均模型計算得到的相速度頻散曲線(中)和ZH振幅比(右)與理論模型計算值的比較

圖6　頻散與ZH振幅比兩步法聯(lián)合反演(黑色實線及陰影部分)和頻散數(shù)據(jù)單獨反演(虛線)的一維后驗概率密度函數(shù)分布對比.每幅圖中的豎直實線或虛線分別對應于聯(lián)合反演或頻散曲線單獨反演的后驗平均模型

圖7給出了不同模型參數(shù)之間的二維邊際后驗概率密度函數(shù)分布.與圖4中一步反演結(jié)果對比，橫波速度參數(shù)之間的二維邊際后驗概率密度函數(shù)等值線彌散性明顯變小，幾乎呈圓形分布，且半徑很小，顯示出地殼中4層的橫波速度均能很好地確定，且各層速度之間的相關(guān)性也很小.0～6 km這一層的橫波速度與Vp/Vs值幾乎沒有相關(guān)性，但0～6 km的橫波速度與6 km-Moho深度的Vp/Vs值存在較為明顯的負相關(guān)性.0～6 km的Vp/Vs值與6 km-Moho深度的Vp/Vs值沒有明顯的相關(guān)性，但0～6 km的Vp/Vs值的反演誤差要明顯小一些，這與圖6中一維邊際后驗概率密度函數(shù)的分析結(jié)果一致.

經(jīng)過對比，兩步反演的結(jié)果更好，分辨率更高，所以之后我們采用兩步反演法對實際數(shù)據(jù)進行反演.

3.2.3瑞利面波頻散數(shù)據(jù)單獨反演

使用頻散數(shù)據(jù)來反演地下速度結(jié)構(gòu)是比較常見的方法.這里我們執(zhí)行和3.2.2中聯(lián)合反演相同的反演流程，但僅使用31個面波頻散數(shù)據(jù)參與反演.本次反演中得到的后驗平均模型和理論值差距較大(圖6垂直虛線所示)，反演得到的一維后驗概率分布函數(shù)分布較寬，尤其是0～6 km和6～15 km這 兩層的橫波速度；且有的參數(shù)出現(xiàn)明顯的多峰值情況，說明僅用頻散曲線反演得到的模型參數(shù)的非唯一性較強.經(jīng)比較，我們舍棄了反演得到的后驗平均模型而采用了所有模型中的“最優(yōu)”模型(數(shù)據(jù)擬合誤差最小)，結(jié)果如圖8所示.由圖可見，單獨反演的最優(yōu)模型雖然能夠較完美地擬合頻散數(shù)據(jù)，但反演得出的ZH振幅比和理論值的差距卻較大.這說明傳統(tǒng)的頻散數(shù)據(jù)反演方法存在一定的局限性，對一些模型參數(shù)的約束不強，從而會導致部分參數(shù)的反演結(jié)果存在較大的誤差.

圖7　基于理論模型合成數(shù)據(jù)的NA兩步聯(lián)合反演所獲得的不同模型參數(shù)之間的二維后驗概率密度函數(shù)分布

圖8　頻散曲線單獨反演結(jié)果，左圖中黑色實線為NA反演中的最優(yōu)模型，其它同圖6

4應用

4.1數(shù)據(jù)

通常無論是地震面波還是背景噪聲，ZH振幅比的實際測量結(jié)果都存在較大的誤差.Tanimoto和Rivera(2008)認為這種誤差分布是穩(wěn)定的，可以用統(tǒng)計方式進行處理.為了減小誤差，我們在提取地震面波振幅比的數(shù)據(jù)處理過程中設(shè)定了多個篩選條件.如信噪比(SNR)須大于10，瑞利面波徑向和垂向分量之間要近似滿足90°相位差，并且刪除ZH振幅比明顯異常的數(shù)據(jù)點等.圖9為我們得到的昆明臺下方45 s和85 s周期的ZH振幅比分布，可以看到它們大致呈現(xiàn)高斯分布.最終處理所得ZH振幅比曲線與已有的地殼上地幔模型(Yao et al.,2010)產(chǎn)生的ZH振幅比曲線基本吻合，證明本次反演所用ZH振幅比測量數(shù)據(jù)較為可靠.

本次聯(lián)合反演頻散觀測數(shù)據(jù)共39個，來源于Yao等(2010)對青藏高原東南部地區(qū)瑞利面波成像在昆明臺所處經(jīng)緯度處的插值結(jié)果，周期范圍10～150 s.ZH振幅比的觀測值共29個，周期范圍25～110 s.

4.2參考模型

NA聯(lián)合反演參考模型包括地殼及上地幔兩部分，上地幔310 km以下設(shè)為半無限空間.地殼分為三層(0～15 km，15～30 km，30～45.554 km)，從上到下各層橫波速度分別為3.3 km·s-1、3.6 km·s-1、3.8 km·s-1，Moho面深度參考Yao等(2010)，設(shè)為45.554 km，上地幔部分從45.554 km到310 km仍然取自ak135模型.

4.3NA聯(lián)合反演步驟及結(jié)果

我們采用和理論測試中相同的反演流程，Ns和Nr取值也與理論測試相同.

第一步：反演7個參數(shù)，分別為3層地殼和4層上地幔(Moho-80 km，80～120 km，120～210 km，210～310 km)的橫波速度.

第二步：反演6個參數(shù)，分別為3層地殼和1層上地幔(Moho-310 km)的橫波速度Vs，以及地殼第一層(0～15 km)和其余兩層(15 km-Moho)的縱橫波速比Vp/Vs.第二步的參考模型選取第一步反演獲得的后驗平均模型值，在第一步的基礎(chǔ)上對地殼各層的橫波速度進行調(diào)整，并將重點放在縱橫波速比的反演上.

聯(lián)合反演結(jié)果如圖10—12所示.圖10中，黑色粗實線為最終反演模型，由該模型正演得到的相速度頻散曲線及ZH振幅比曲線基本都處于觀測值的誤差范圍內(nèi).但當周期在140 s以上時，計算值與觀測值差距較大，可能是因為長周期的觀測頻散數(shù)據(jù)存在較大的誤差，所以上地幔深部的反演結(jié)果不太可信.圖11顯示各反演參數(shù)的一維邊際后驗概率密度函數(shù)分布基本呈高斯型分布，且范圍較窄，顯示出很好的模型分辨能力.需要指出的是，15 km-Moho范圍內(nèi)反演得到的Vp/Vs值的后驗概率分布有兩個峰值，其中一個明顯高于另一個，此處我們并不像選取其他參數(shù)那樣選取后驗概率平均模型值(位于兩峰之間)，而是取高的峰(最大似然值)對應的Vp/Vs值作為最終反演結(jié)果.另外，我們將地殼分為兩層并分別反演它們的Vp/Vs值，結(jié)果都在1.72左右，這和已有的利用接收函數(shù)方法得到的該地區(qū)的泊松比大致符合(黃建平等，2008；徐強等，2009)，顯示出ZH振幅比和頻散曲線的聯(lián)合反演可以用于獲得分層地殼的Vp/Vs值(或泊松比).圖12為各反演參數(shù)值之間的二維邊際后驗概率密度函數(shù)等值線，和理論測試結(jié)果(圖7)對比，實際數(shù)據(jù)結(jié)果的彌散性更強，但最內(nèi)層的60%等值線(紅色)幾乎呈圓形且半徑很小.

值得注意的是，雖然第一步反演的參考模型的選取比較隨意，但如果參數(shù)范圍合適，最終結(jié)果并未受影響：反演第一步所得模型(即為第二步反演的參考模型)的正演結(jié)果(頻散曲線和ZH比)已經(jīng)和觀測值比較接近了，到反演第二步，其反演模型的正演結(jié)果已經(jīng)基本落在了觀測值的誤差范圍以內(nèi).由圖10可見，最初的參考模型對應的正演曲線(虛線)和最終反演模型對應的曲線(粗實線)相差較大，但經(jīng)過兩步聯(lián)合反演，即使是不太準確的初始模型也可以得到比較準確的結(jié)果，有力驗證了我們所使用的聯(lián)合反演算法的穩(wěn)定性.

圖10　實際聯(lián)合反演得到的昆明臺下方地殼及上地幔的橫波速度結(jié)構(gòu)(左)，實際數(shù)據(jù)聯(lián)合反演得到的頻散數(shù)據(jù)(中)和ZH振幅比(右)與觀測值的對比.左圖中黑色實線為最終反演得到的后驗平均模型，灰色區(qū)域顯示其標準誤差，虛線和點虛線分別為第一次和第二次NA反演的參考模型；中間和右圖中帶誤差棒的星號代表觀測數(shù)據(jù)，虛線和三角分別對應于第一步和第二步NA反演的參考模型計算得到的理論數(shù)據(jù)，實線對應最終反演獲得的后驗平均模型計算得到的理論數(shù)據(jù)

圖11　實際地震數(shù)據(jù)聯(lián)合反演得到的模型參數(shù)的一維后驗概率密度函數(shù)分布(與圖6類似).除右下角的圖外，黑色豎直實線對應于后驗平均模型；右下角的圖(Vp/Vs [15-Moho])的黑色豎直實線對應于模型的最大似然值

5討論及結(jié)論

用相速度頻散曲線來反演地下速度結(jié)構(gòu)的方法應用廣泛，但在某些深度范圍的分辨率不足，尤其對地殼淺部的約束不夠.ZH振幅比有獨立于相速度的深度敏感核，對于同一周期，其敏感深度范圍與相速度的敏感深度范圍顯著不同，更加接近地殼淺層.二者結(jié)合起來，可顯著提高對地殼淺部結(jié)構(gòu)的分辨能力.

本文提出了基于鄰域算法(NA)的頻散曲線和ZH比聯(lián)合反演方法，相比于一步反演，采用兩步反演更穩(wěn)定，能夠更為準確地獲得地殼上地幔橫波速度結(jié)構(gòu)及分層地殼的Vp/Vs值.通過理論測試我們發(fā)現(xiàn)，當一步反演的參數(shù)較多時，即使得到的最優(yōu)模型和理論模型非常接近，統(tǒng)計所得的后驗平均模型卻可能和理論模型相差較遠，模型搜索容易陷入局部極小.另外，聯(lián)合反演明顯優(yōu)于頻散曲線單獨反演，能更好地約束模型參數(shù)并降低其反演誤差.之后我們將兩步聯(lián)合反演應用于實際測量數(shù)據(jù)，獲得了昆明臺下方更為可靠的地殼速度結(jié)構(gòu)模型.此外，基于鄰域算法的聯(lián)合反演可以通過后驗概率密度分布函數(shù)更好地分析模型的誤差和各參數(shù)間的相關(guān)性.

之前的相關(guān)研究或者沒有給出Vp/Vs的反演值，或者發(fā)現(xiàn)其與密度值之間很強的相互折衷性而沒能得到確定的結(jié)果.在本文中，我們采用了更為合理的縱波速度與密度之間的經(jīng)驗公式(Brocher,2005)，并不反演密度，從而更好地約束Vp/Vs值.從昆明臺真實數(shù)據(jù)的反演結(jié)果看，聯(lián)合反演獲得的分層地殼的Vp/Vs值與以前通過接收函數(shù)獲得的地殼的平均泊松比結(jié)果比較吻合，說明頻散曲線和ZH比的聯(lián)合反演可以獨立地提供分層地殼泊松比的約束，這為認識地殼不同深度介質(zhì)的物理和成分屬性提供新的約束.

另外有一點值得注意：在用實際地震數(shù)據(jù)進行反演時，我們將地殼分為了三層，其厚度基本相等.而理論測試時，我們卻分了四層，相當于將第一層(15 km)又分為兩層(0～6 km和6～15 km).這是因為ZH振幅比主要是對淺層的約束比較強，且起始周期越小反演效果越好.理論測試數(shù)據(jù)的周期從5 s開始，所以對所劃分的四層地殼的約束都很好.而實際反演中，相速度的周期從10 s開始，ZH振幅比的周期則從25 s開始，所以很難約束淺部幾千米的結(jié)構(gòu).對于地震面波，我們通常測量得到的ZH振幅比的周期都在20 s以上(Lin et al.，2012)，為了更好地約束地殼淺層的速度結(jié)構(gòu)，我們可以通過連續(xù)背景噪聲信號(Tanimoto et al.，2013)或背景噪聲的互相關(guān)數(shù)據(jù)(Lin et al.，2014)來得到短周期的ZH振幅比.

綜上，頻散曲線和ZH比的聯(lián)合反演確實要優(yōu)于頻散曲線的單獨反演，但聯(lián)合反演結(jié)果的分辨率也受到數(shù)據(jù)質(zhì)量、周期范圍等因素的影響.此外，考慮到面波頻散及ZH振幅比對界面信息的約束不足，本文并沒有反演Moho面或殼內(nèi)其他界面的深度，若在現(xiàn)有基礎(chǔ)上再加入對分界面敏感的接收函數(shù)數(shù)據(jù)，可能會得到更好的效果.面波頻散和接收函數(shù)均對介質(zhì)的S波速度結(jié)構(gòu)敏感，但二者的側(cè)重點不同.面波頻散曲線可以揭示介質(zhì)的橫波速度結(jié)構(gòu)變化，而接收函數(shù)對橫波的絕對速度約束比較差，但在確定轉(zhuǎn)換界面深度方面具有獨特的優(yōu)勢，故聯(lián)合面波頻散和接收函數(shù)資料可以獲得比單一方法更為可靠的橫波速度結(jié)構(gòu)(胡家富等，2005；劉啟元等，2010).近年來，面波頻散和接收函數(shù)的聯(lián)合反演廣泛應用于地殼上地幔速度結(jié)構(gòu)研究中(胡家富等，2005；劉啟元等，2010).開展聯(lián)合反演的目的是為了實現(xiàn)優(yōu)勢互補，將面波頻散曲線、ZH振幅比曲線及接收函數(shù)聯(lián)合來反演地殼上地幔速度結(jié)構(gòu)，可以更為準確地約束模型信息，減小反演結(jié)果的不確定性，這也是目前研究的熱點方向之一.

人才的培養(yǎng)始終圍繞為社會服務而展開。在經(jīng)濟飛速發(fā)展的大背景下，體育旅游消費人群的產(chǎn)品需求也在發(fā)生改變。加強學生基礎(chǔ)能力培養(yǎng)不放松，適時根據(jù)社會需求作出相應調(diào)整是必要的。避免長期的同質(zhì)化培養(yǎng)，高效區(qū)分體育旅游市場，不斷完善多元體育旅游人才的有效供給。

致謝感謝兩位審稿人對本論文提出的修改建議，感謝美國University of Colorado,Boulder的沈偉森博士提供自動的面波時頻分析程序用于地震面波振幅比的測量工作.

References

Aki K,Richards P G.1980.Quantitative Seismology,Vol.1: Theory and Methods.San Francisco:WH Freeman and Company.Boore D M,Toks?z M N.1969.Rayleigh wave particle motion and crustal structure.Bulletin of the Seismological Society of America,59(1): 331-346.

Brocher T M.2005.Empirical relations between elastic wavespeeds and density in the earth′s crust.Bulletin of the Seismological Society of America,95(6): 2081-2092.

Chong J J,Ni S D,Zhao L.2015.Joint inversion of crustal structure with the rayleigh wave phase velocity dispersion and theZHRatio.Pure and Applied Geophysics,172(10): 2585-2600.

Herrmann R B,Ammon C J.2004.Computers programs in seismology,surface waves,receiver functions and crustal structure,version 3.30.Department of Earth and Atmospheric Sciences.Saint Louis University.Hu J F,Zhu X G,Xia J Y,et al.2005.Using surface wave and receiver function to jointly inverse the crust-mantle velocity structure in the West Yunnan area.Chinese J.Geophys.(in Chinese),48(5): 1069-1076.

Huang J P,Chong J J,Ni S D,et al.2008.Inversing the crustal thickness under the stations of China via H-Kappa method.Journal of University fo Science and Technology of China (in Chinese),38(1): 33-40.

Levshin A L,Pisarenko V F,PogrebinskyG A.1972.On a frequency-time analysis of oscillations.Annales Geophysicae,28(2): 211-218.

Lin F C,Schmandt B,Tsai V C.2012.Joint inversion of rayleigh wave phase velocity and ellipticity using usarray: Constraining velocity and density structure in the upper crust.Geophysical Research Letters,39: L12303.

Lin F C,Tsai V C,Schmandt B.2014.3-D crustal structure of the western United States: Application of Rayleigh-wave ellipticity extracted from noise cross-correlations.Geophysical Journal International,198(2): 656-670.

Liu Q Y,Li Y,Chen J H,et al.2010.Joint inversion of receiver function and ambient noise based on Bayesian theory.Chinese J.Geophys.(in Chinese),53(11): 2603-2612.

Masters G,Laske G,Bolton H,et al.2000.The relative behavior of shear velocity,bulk sound speed,and compressional velocity in the mantle: implications for chemical and thermal structure.∥ Earth′s Deep Interior: Mineral Physics and Tomography From the Atomic to the Global Scale.Washington DC,American Geophysical Union,63-87.

Sambridge M.1999a.Geophysical inversion with a neighbourhood algorithm—Ⅰ.Searching a parameter space.Geophysical Journal International,138(2): 479-494.

Sambridge M.1999b.Geophysical inversion with a neighbourhood algorithm—Ⅱ.Appraising the ensemble.Geophysical Journal International,138(3): 727-746.

Shapiro N M,Campillo M,Stehly L,et al.2005.High-resolution surface-wave tomography from ambient seismic noise.Science,307(5715): 1615-1618.

Tanimoto T,Alvizuri C.2006.Inversion of the HZ ratio of microseisms forS-wave velocity in the crust.Geophysical Journal International,165(1): 323-335.

Tanimoto T,Rivera L.2008.The ZH ratio method for long-period seismic data: sensitivity kernels and observational techniques.Geophysical Journal International,172(1): 187-198.

Tanimoto T,Yano T,Hakamata T.2013.An approach to improve rayleigh-wave ellipticity estimates from seismic noise: Application to the Los Angeles Basin.Geophysical Journal International,193(1): 407-420.

Xu G M,Yao H J,Zhu L B,et al.2007.Shear wave velocity structure of the crust and upper mantle in western China and its adjacent area.Chinese J.Geophys.(in Chinese),50(1): 193-208.

Xu Q,Zhao J M,Cui Z X,et al.2009.Structure of the crust and upper mantle beneath the southeastern Tibetan Plateau by P and S receiver functions.Chinese J.Geophys.(in Chinese),52(12): 3001-3008.

Yano T,Tanimoto T,Rivera L.2009.The ZH ratio method for long-period seismic data: Inversion for S-wave velocity structure.Geophysical Journal International,179(1): 413-424.

Yao H J,Beghein C,Van Der Hilst R D.2008.Surface wave array tomography in SE Tibet from ambient seismic noise and two-station analysis—Ⅱ.Crustal and upper-mantle structure.Geophysical Journal International,173(1): 205-219.

Yao H J,Van Der Hilst R D,Montagner JP.2010.Heterogeneity and anisotropy of the lithosphere of Se Tibet from surface wave array tomography.Journal of Geophysical Research,115(B12),doi: 10.1029/2009JB007142.

附中文參考文獻

胡家富,朱雄關(guān),夏靜瑜等.2005.利用面波和接收函數(shù)聯(lián)合反演滇西地區(qū)殼幔速度結(jié)構(gòu).地球物理學報,48(5): 1069-1076.

黃建平,崇加軍,倪四道.2008.利用h-Kappa方法反演中國地區(qū)臺站下地殼厚度.中國科學技術(shù)大學學報,38(1): 33-40.

劉啟元,李昱,陳九輝等.2010.基于貝葉斯理論的接收函數(shù)與環(huán)境噪聲聯(lián)合反演.地球物理學報,53(11): 2603-2612.

徐果明,姚華建,朱良保等.2007.中國西部及其鄰域地殼上地幔橫波速度結(jié)構(gòu).地球物理學報,50(1): 193-208.

徐強,趙俊猛,崔仲雄等.2009.利用接收函數(shù)研究青藏高原東南緣的地殼上地幔結(jié)構(gòu).地球物理學報,52(12): 3001-3008.

(本文編輯胡素芳)

基金項目國家自然科學基金(41222028)，地震行業(yè)專項(201508008)，中央高?；究蒲袠I(yè)務費專項(WK2080000053)資助.

作者簡介袁藝，女，1987年生，2015年碩士畢業(yè)于中國科學技術(shù)大學地球物理專業(yè)，目前在中國地震局第二監(jiān)測中心工作，主要從事面波成像方面的研究工作.E-mail：yuanyigre@126.com *通訊作者姚華建，男，1979年生，現(xiàn)為中國科學技術(shù)大學地球和空間科學學院教授、博士生導師.2001年及2004年分別獲得中國科學技術(shù)大學本科和碩士學位，2009年獲得麻省理工學院地球物理學博士學位.主要從事地震波及背景噪聲成像、大地震破裂過程、臺陣技術(shù)、地球物理聯(lián)合反演等領(lǐng)域的研究.E-mail: hjyao@ustc.edu.cn

doi:10.6038/cjg20160318 中圖分類號P315

收稿日期2015-03-10，2015-10-21收修定稿

Joint inversion of Rayleigh wave vertical-horizontal amplitude ratios and dispersion based on the Neighborhood Algorithm and its application

YUAN Yi1,2,YAO Hua-Jian1,2*,QIN Yan1,2

1LaboratoryofSeismologyandPhysicsofEarth′sInterior,SchoolofEarthandSpaceSciences,UniversityofScienceandTechnologyofChina,Hefei230026,China2NationalGeophysicalObservatoryatMengcheng,AnhuiBozhou233527,China

AbstractRayleigh wave ellipticity (or ZH ratio) is a function of frequency and is particularly sensitive to shallow crustal structure beneath the seismograph station.Since depth sensitivity kernels of ZH ratios are different from those of dispersion data,the ZH ratio provides good complementary information for the dispersion-based inversion method.Therefore,we can combine the ZH ratio and dispersion data of Rayleigh wave fundamental mode to better invert for the velocity structure under a specific seismograph station.In this paper,we propose a joint inversion method using the dispersion and ZH ratio data based on the Neighborhood Algorithm.We conduct synthetic tests based on a theoretical model and prove the robustness of the joint inversion method,which can better constrain the shallow crustal structure.Compared to traditional inversion methods that only use dispersion data,the joint inversion can provide a more accurate crustal Vs model as well as Vp/Vs ratios for the layered crust.Finally,we apply the joint inversion technique to real measurements and obtain a more accurate crust shear velocity and Vp/Vs model beneath the station at Kunming (KMI) in southwest China.

KeywordsJoint inversion； Neighborhood Algorithm； Rayleigh waves； Phase velocity dispersion； ZH ratio； Crustal structure

袁藝，姚華建，秦巖.2016.基于鄰域算法的瑞利面波垂直-水平振幅比及頻散曲線聯(lián)合反演及應用.地球物理學報，59(3):959-971,doi:10.6038/cjg20160318.

Yuan Y,Yao H J,Qin Y.2016.Joint inversion of Rayleigh wave vertical-horizontal amplitude ratios and dispersion based on the Neighborhood Algorithm and its application.Chinese J.Geophys.(in Chinese),59(3):959-971,doi:10.6038/cjg20160318.