一維大地電磁貝葉斯反演研究

郭 曼，吳姿穎，張晶晶，陳姝霓，曾志文，丁文偉，郝冬冬

（1.東華理工大學(xué)，江西南昌330013；2.楊凌龍翔數(shù)字科技有限公司，陜西咸陽712000）

在地球物反演過程中解的多解性是必然存在的，貝葉斯反演理論為解決解的多解性提供了思路。王文濤[1]針對地震數(shù)據(jù)中地震參數(shù)的非線性問題，采用貝葉斯反演基于模擬退火算法進(jìn)行采樣，實現(xiàn)了地震儲層的不確定性評估。趙巒嘯等[2]結(jié)合貝葉斯框架中不同的地球物理信息進(jìn)行水庫倒塌實現(xiàn)定量評估。楊迪坤等[3]將數(shù)據(jù)空間的不確定性轉(zhuǎn)移到大地電磁反演中的模型空間，利用概率分布來反映模型參數(shù)的可能性。郭榮文等[4]利用改進(jìn)的自適應(yīng)純模擬退火算法（AS?SA）實現(xiàn)了一維磁層貝葉斯反演，利用Metroplois-Hasing采樣獲得了有關(guān)地基介質(zhì)不確定度的信息。史學(xué)東等[5]通過在一維大地電磁反演中使用貝葉斯方法獲得了模型參數(shù)與參數(shù)之間的相關(guān)性，得出來其實數(shù)據(jù)本身除了含有地電模型信息，還有其他有用的反演信息。陳曉等[6]采用模擬退火算法對貝葉斯框架下的大地電磁和地震數(shù)據(jù)進(jìn)行聯(lián)合反演，利用物性參數(shù)的后驗概率分布來表達(dá)反演結(jié)果，并取得了較好的效果。但在實現(xiàn)過程中沒有對不同物性參數(shù)對貝葉斯反演結(jié)果的影響進(jìn)行系統(tǒng)的討論，鑒于此，本文基于非?？焖俚哪M退火算法討論各參數(shù)對一維MT貝葉斯反演結(jié)果的影響進(jìn)行討論。

1 貝葉斯反演理論

在地球物理反演中，貝葉斯反演理論框架于1987年Tarantola就提出了結(jié)合數(shù)據(jù)信息和模型信息的模型后驗概率的公式：

式中：p(m)——模型m的先驗概率;

p(m|d)——在模型m下的條件概率，也稱似然函數(shù);

p(d)——觀測數(shù)據(jù)模型全空間概率；

σ(m|d)——在觀測數(shù)據(jù)d下模型m的后驗概率。

在計算后驗概率分布時，p(d)可視為一個常數(shù)僅起正則化因子的作用。那么后驗概率分布可進(jìn)一步表述為：

上式P(d|m)P(m)稱為后驗概率分布的核，實際勘探中的觀測數(shù)據(jù)總會存在一定的噪音,同時反演的模型結(jié)果也會有誤差,如果假設(shè)這些噪音和誤差的分布均為高斯分布時,后驗概率分布就可表述為：

計算出樣本模型的后驗概率分布后可以按如下公式來求取模型參數(shù)的均值[1]：

2 相關(guān)參數(shù)對后驗概率密度分布的影響

為了討論各參數(shù)對貝葉斯反演理論的影響，設(shè)計如圖1所示5層地質(zhì)模型：計算頻點在對數(shù)區(qū)間取40個。反演過程中，非常快速地模擬退火算法（VFSA）的初始溫度為10℃,降溫系數(shù)為0.96，每個溫度進(jìn)行10次擾動，進(jìn)行500次降溫。取真實值的±30%作為擾動區(qū)間，最終得到反演結(jié)果。

本文選取一個初始擾動模型作為對比模型，進(jìn)行研究相關(guān)參數(shù)的變化對后驗概率分布提取的影響，物性參數(shù)如表1所示，各參數(shù)反演后驗概率分布如圖2所示。

2.1 先驗概率P（m）的影響

假若先驗概率分布P（m）滿足高斯分布，則所涉及的2個關(guān)鍵參數(shù)為先驗均值和方差。為了討論先驗概 率P（m）的影響，本文以上述模型為例設(shè)計如下方案。

圖1 理論模型

表1 參數(shù)未改變原始參考模型

圖2 各參數(shù)反演后驗概率分布

方案一：保持先驗方差不變，通過改變單個先驗均值和改變所有先驗均值試驗先驗均值影響，該方案的后驗概率密度分布如圖3所示；

方案二：保持先驗均值不變，通過改變單個先驗方差和改變所有先驗方差試驗先驗方差影響；該方案的后驗概率密度分布如圖4所示。

2.1.1 先驗均值的影響

（1）單個均值的改變。以模型的第四層厚度為例，先驗均值為4500Ω·m，先驗方差為850000。將先驗均值從 4500Ω·m 增大至 5500Ω·m 和 7500Ω·m ，先驗方差保持不變。可以看出，先驗均值為5500Ω·m時除在真實值附近有峰值，在5500Ω·m附近也形成了極高的的峰值。先驗均值再變大至7500Ω·m則后驗概率密度分布大都集中于5000Ω·m。這是由于擾動取值范圍是2800～5200Ω·m的緣故，所以后驗概率分布受先驗影響較大。在先驗方差確定的情況下，后驗概率密度分布在先驗均值附近集中。當(dāng)先驗均值偏離真實值較大時，后驗概率密度分布呈現(xiàn)多峰。

圖3 各均值改變后概率分布

（2）改變模型所有參數(shù)的先驗均值。見表2，后驗概率密度分布參見圖3。參考對比模型圖2，若是給定的先驗均值相對真實值越小，其后驗概率分布會在真實值左邊靠近給定的先驗均值范圍分布。相反相對真實值越大，其真實值右邊靠近給定的先驗均值范圍分布?？梢缘贸鱿闰灳禃绊懞篁灨怕史植?，其模型可能取值范圍分布會偏向于先驗區(qū)間。當(dāng)先驗信息取高斯分布時，最好是能均值能接近于真實值的效果才最好。

2.1.2 先驗方差的影響

（1）單個方差的改變。仍以第四層厚度模型為例：方差由850000減至85000和8500，經(jīng)驗證得出，方差越小其后驗概率峰值分布會在真實值右邊偏離。方差增大至85000000和8500000000，經(jīng)驗證得出，峰值區(qū)間依然靠近真實值附近。

（2）改變模型所有參數(shù)的方差。見表3，后驗概率密度分布參見圖4?？梢钥闯觯合闰灧讲畋碚髁讼闰灨怕拭芏绕x先驗均值的偏離度。方差較大時對后驗概率分布的提取基本影響不大；較小時則會壓制均值兩邊分布比重，增大均值附近的分布比重形成密集分布，所以當(dāng)先驗均值的可信性不高時，應(yīng)設(shè)置較大的先驗方差。

表2 先驗均值改變前后對比值

表3 先驗均值改變前后對比值

圖4 方差改變后后驗概率分布

2.2 溫度的影響（初始溫度T=1℃）

從圖5可以看出，溫度增高至T=10℃，大量干擾模型被接受，其分布形態(tài)越來越明顯，峰值區(qū)間基本是在真實值附近。繼續(xù)增高溫度至T=1000℃，接受高能量模型更多，但峰值區(qū)間卻越來越不能凸顯。

圖5 T=10℃

圖6 T=0.1℃

從圖6可以看出，當(dāng)降低溫度至T=0.1℃，許多原本被接受收的模型消失了，雖然可以大致反映出后驗概率分布，但局部值有凸顯出來。當(dāng)溫度降溫繼續(xù)降低至T=0.05℃，模型采樣樣本越來越少，效果加劇。也就是說高溫狀態(tài)雖然采樣比較充分但是會接受較多的高能量解，這會影響后驗概率的提取和定量評價，導(dǎo)致后驗概率密度分布不集中；而在低溫狀態(tài)則會出現(xiàn)采樣不足的現(xiàn)象。結(jié)合模擬退火算法全局尋優(yōu)的特點以及貝葉斯反演采樣的充足性，應(yīng)舍棄高溫時接受的電阻率模型，從中溫開始接受樣本發(fā)揮各自的優(yōu)點。

3 結(jié)論

本文基于模擬退火算法實現(xiàn)一維大地電磁貝葉斯反演。模型試驗表明：貝葉斯反演理論可實現(xiàn)大地電磁反演解的定量評價。除此，相關(guān)參數(shù)的選取對后驗概率密度分布具有一定的影響。具體而言，對于均值和方差，當(dāng)均值和模型真實值存在一定偏差時，方差應(yīng)選取較大值才可對參數(shù)后驗概率密度分布影響較小。對于溫度，高溫采樣充分，但會接受較多的高能量解，這會影響后驗概率的提取和定量評價，導(dǎo)致后驗概率密度分布不集中。低溫則會采樣不足?？紤]到模擬退火算法的全局尋優(yōu)特點以及貝葉會采樣不足，所以應(yīng)舍棄過高溫和過低溫，從中溫開始接受樣本。

[1] 王文濤.地震儲層評價與預(yù)測的貝葉斯反演方法研究[D].中國地質(zhì)大學(xué)（武漢）,2008.

[2] 趙巒嘯.貝葉斯框架下基于疊前地震數(shù)據(jù)的巖性（流體）預(yù)測[D].上海：同濟(jì)大學(xué),2010.

[3] 楊迪坤，胡祥云.含噪聲數(shù)據(jù)反演的概率描述[J].地球物理學(xué)報，2008，51(3):901-907.

[4] 郭榮文.貝葉斯MT反演的非線性和不確定度分析[D].中南大學(xué),2011.

[5] 史學(xué)東,柳建新,郭榮文，等.大地電磁法的1D無偏差貝葉斯反演[J].物探化探計算技術(shù)，2012（4）：371-379.

[6] 陳曉,李文喬,郭曼,等.大地電磁測深和重力數(shù)據(jù)貝葉斯聯(lián)合反演[J].科學(xué)技術(shù)與工程,2016，16(15):30-35.