基于連通性函數(shù)的多點(diǎn)地質(zhì)統(tǒng)計(jì)建模參數(shù)敏感性分析方法

杜 川, 喻思羽, 李少華, 方 紅

(1.長江大學(xué) 地球科學(xué)學(xué)院,武漢 430100；2.中國石油 天然氣遼河油田分公司 勘探開發(fā)研究院,盤錦 124010)

0 引言

隨著多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)(multi-point statistics， MPS)受到廣泛的關(guān)注，涌現(xiàn)出許多相對(duì)成熟的MPS建模算法。相對(duì)于傳統(tǒng)的兩點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)，多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)利用訓(xùn)練圖像代替變差函數(shù)表達(dá)地質(zhì)變量的空間結(jié)構(gòu)性[1]，從中獲取更高階的樣本統(tǒng)計(jì)量[2],MPS建模優(yōu)點(diǎn)在于既能夠模擬地質(zhì)體復(fù)雜的幾何結(jié)構(gòu)，也能忠實(shí)于條件數(shù)據(jù)。Guardiano等[2]提出的首個(gè)MPS算法每模擬一個(gè)網(wǎng)格節(jié)點(diǎn)都需要遍歷整個(gè)訓(xùn)練圖像，計(jì)算成本非常高，使得該算法難以實(shí)際應(yīng)用。之后，更好的算法被開發(fā)出來改善這些缺點(diǎn)。多點(diǎn)統(tǒng)計(jì)技術(shù)包括基于象元的建模算法SNESIM[3]、基于樣式的建模算法FILTERSIM[4]、SIMPAT[5]、DISPAT[6]、基于直接采樣的算法DS[7]、IMPACA[8]、HOSIM[9]、CCSIM[10]、基于圖像縫合技術(shù)的Image Quilting[11]。這些MPS算法使用不同的方法來讀取訓(xùn)練圖像、存儲(chǔ)信息和復(fù)制模式。Mariethoz[12]對(duì)這些多點(diǎn)地質(zhì)建模方法進(jìn)行了定性比較。多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模已用于許多應(yīng)用之中，例如常見的在儲(chǔ)層地質(zhì)建模研究領(lǐng)域，利用Simpat模擬河流相儲(chǔ)層分布[13]；多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)在河流相儲(chǔ)層建模中的應(yīng)用[14]；在三角洲相儲(chǔ)層建模領(lǐng)域的應(yīng)用，如多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模方法在復(fù)雜疊置樣式砂體表征中的應(yīng)用[15]。

在兩點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模中，通過變差函數(shù)考慮統(tǒng)計(jì)數(shù)據(jù)之間的空間相關(guān)性[16]。類似的在多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模過程，需要通過給定的數(shù)據(jù)事件對(duì)訓(xùn)練圖像進(jìn)行掃描來獲取一個(gè)統(tǒng)計(jì)先驗(yàn)?zāi)Ｐ?，相關(guān)的統(tǒng)計(jì)數(shù)據(jù)取決于給定數(shù)據(jù)事件的參數(shù)，通過人為的參數(shù)調(diào)整，找到合適的參數(shù)設(shè)置是獲取良好建模效果的前提，不同的參數(shù)不僅影響建模的質(zhì)量，也涉及CPU及內(nèi)存開銷問題[17]。為此前人也進(jìn)行了大量的探究，常見的參數(shù)敏感性分析，如Liu[18]提供了關(guān)于SNESIM參數(shù)以及每個(gè)參數(shù)對(duì)模型再現(xiàn)質(zhì)量和CPU成本的完整指南；Meerschman[19]也為DS直接采樣算法做了類似的工作；Honarkhah[6]提出了DISPAT，它對(duì)用戶提供的參數(shù)不太敏感。在該方法中，通過計(jì)算不同樣板尺寸的圖案的平均熵，然后使用所謂的肘形圖的拐點(diǎn)，自動(dòng)選擇最佳樣板尺寸。Strebelle[20]提出了一種通過優(yōu)化數(shù)據(jù)樣板的大小來最小化SNESIM的計(jì)算成本，同時(shí)保持模式再現(xiàn)質(zhì)量的方法。在該方法中，最小可接受的數(shù)據(jù)樣板大小是通過閾值來選擇的，超過該閾值，額外的條件數(shù)據(jù)不會(huì)提高估計(jì)的條件概率。Kolbj?rnsen[21]開發(fā)了一種自動(dòng)方法來定量確定馬爾可夫網(wǎng)格建模中使用的多個(gè)網(wǎng)格的最佳數(shù)量，該方法使用從訓(xùn)練圖像獲得的方向相關(guān)函數(shù);Bai[22]還提出了一種量化方法來估計(jì)SNESIM算法的多個(gè)網(wǎng)格的最小數(shù)量。它使用由改進(jìn)的聯(lián)合計(jì)數(shù)統(tǒng)計(jì)量化的空間關(guān)聯(lián)度來尋找最粗略的尺度，并因此基于模擬尺度應(yīng)不小于目標(biāo)尺度的標(biāo)準(zhǔn)來估計(jì)多個(gè)網(wǎng)格的最小數(shù)量。

傳統(tǒng)建模參數(shù)的評(píng)價(jià)主要依靠人工識(shí)別，即首先給定一組(有序)參數(shù)集，使用該參數(shù)集里的每個(gè)參數(shù)模擬一組隨機(jī)模型，然后通過人工判別優(yōu)選與訓(xùn)練圖像相似度較為合適的模型，以該模型對(duì)應(yīng)的參數(shù)作為最優(yōu)參數(shù)。人工識(shí)別的精確性取決于建模工作者的經(jīng)驗(yàn)，具有較強(qiáng)主觀性，同時(shí)人工識(shí)別效率低，不適應(yīng)于現(xiàn)代自動(dòng)化生產(chǎn)的需要。為了解決這一問題，我們引入了一種基于連通性函數(shù)[23-24]的多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模參數(shù)敏感性分析方法，通過計(jì)算條件數(shù)據(jù)下生成的模擬圖像的連通性概率曲線，并把訓(xùn)練圖像與模擬圖像連通性概率曲線之間的相似性用余弦相似度來量化，得到基于連通性函數(shù)的空間相關(guān)性評(píng)價(jià)指標(biāo)與建模參數(shù)的關(guān)系曲線，通過對(duì)量化數(shù)據(jù)進(jìn)行分析，找到相似度曲線趨于平緩的拐點(diǎn)值，以此拐點(diǎn)值作為最優(yōu)參數(shù)，大于最優(yōu)值參數(shù)將變得不敏感，模型效果不會(huì)隨著參數(shù)值增加而得到明顯改善。

1 連通性計(jì)算原理

連通性函數(shù)被定義為網(wǎng)格上的兩個(gè)點(diǎn)屬于同一個(gè)連接區(qū)域的概率[24]，利用連通性函數(shù)計(jì)算區(qū)域化變量指示圖的連通性，可以反映區(qū)域化變量在某一方向某一距離上連通性概率的變化程度。理論上隨機(jī)模型與訓(xùn)練圖像的連通性差異越小，則表明二者沉積相的空間相關(guān)性及位置分布特征越相近，也說明了用于生成該隨機(jī)模型的建模參數(shù)越可靠。Eulogio[25]詳細(xì)地描述了連通性概率計(jì)算原理，其計(jì)算原理如下：

圖1定義工區(qū)網(wǎng)格為矩形G，G的網(wǎng)格節(jié)點(diǎn)(或象元)屬于子集S(砂巖相)或補(bǔ)集Sc(泥巖相)

G=SUSc

(1)

G的象元取值函數(shù)定義為指示函數(shù)I(u)

(2)

如果u和u′=u+h屬于相同的S，用u?u′表示，連通性函數(shù)τ(h)定義為屬于S的兩點(diǎn)連通的概率：

τ(h)=P(u?u+h│u,u+h∈S)

(3)

其近似計(jì)算公式為式(4)。

(4)

以砂巖相連通性為例，如圖1所示u和u′有三種情況，①u和u+h有一個(gè)或者兩個(gè)都屬于泥巖相，例如u3與u3+h或者u4與u4+h，則不參與砂巖相連通計(jì)算；②u和u+h都是砂巖，但不屬于相同砂巖體，例如u2與u2+h，則增加N(u,u+h∈S)數(shù)值；③u和u+h都是砂巖，屬于相同砂巖體，例如u1與u1+h則增加#N(u?u+h│u,u+h∈S)和#N(u,u+h∈S)值。

圖1 連通性示意圖Fig.1 Schematic diagram of connectivity

2 參數(shù)敏感性分析方法

2.1 實(shí)現(xiàn)步驟

基于連通性函數(shù)的多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模參數(shù)優(yōu)選方法是建立在一個(gè)前提條件下：以經(jīng)典多點(diǎn)算法Snesim的一個(gè)重要建模參數(shù)--“搜索節(jié)點(diǎn)數(shù)”為例，基于連通性函數(shù)優(yōu)選搜索節(jié)點(diǎn)數(shù)；隨著搜索節(jié)點(diǎn)的增加，隨機(jī)模型與訓(xùn)練圖像的空間相關(guān)性、幾何形態(tài)越來越相似。

在此前提條件下采用有序排列的建模參數(shù)集建立一組隨機(jī)模型集，然后計(jì)算基于連通性函數(shù)和余弦相似度CosSim函數(shù),評(píng)價(jià)建模參數(shù)集的多點(diǎn)地質(zhì)統(tǒng)計(jì)隨機(jī)模型與訓(xùn)練圖像的空間相關(guān)性及結(jié)構(gòu)特征相似性，進(jìn)而建立基于連通性函數(shù)的空間相關(guān)性評(píng)價(jià)指標(biāo)與建模參數(shù)的關(guān)系曲線，選取評(píng)價(jià)指標(biāo)開始平穩(wěn)、進(jìn)入平臺(tái)區(qū)域時(shí)的拐點(diǎn)所對(duì)應(yīng)的參數(shù)值作為最優(yōu)參數(shù)。

以建模參數(shù)“搜索節(jié)點(diǎn)數(shù)”為例，基于連通性函數(shù)的建模參數(shù)優(yōu)選方法具體流程及流程圖(圖2)如下：

圖2 優(yōu)選方法流程圖Fig.2 Flow chart of optimization method

1)輸入訓(xùn)練圖像TI，給定MPS的關(guān)鍵建模參數(shù)——搜索節(jié)點(diǎn)的K個(gè)連續(xù)單調(diào)遞增數(shù)值(例如：1，2，…，K)組成的集合S={1，2，…，K}。

2)采用MPS，使用TI和S生成隨機(jī)模型集M，其中第k個(gè)參數(shù)Sk相應(yīng)的隨機(jī)模型集定義為Mk，每個(gè)Mk包含n個(gè)隨機(jī)模型，則M中共計(jì)有N=K×n個(gè)模型。

(5)

(6)

計(jì)算Mk所有模型與訓(xùn)練圖像TI關(guān)于連通性函數(shù)參數(shù)組的余弦相似度，然后取其數(shù)學(xué)期望作為建模參數(shù)Sk對(duì)應(yīng)的空間相關(guān)性評(píng)價(jià)指標(biāo)ξk

ξk=E[cosθk(Ai,Bi)]

(7)

5)k增加1，如果k小于K-1，則進(jìn)入步驟4)，否則進(jìn)入步驟6)。

6)建立參數(shù)k與相似度的關(guān)系曲線，分析曲線變化趨勢(shì)，把曲線變化趨于平緩時(shí)的拐點(diǎn)建模參數(shù)k(i)作為最優(yōu)參數(shù)。隨著建模參數(shù)值的遞增，基于建模參數(shù)值的隨機(jī)模型與訓(xùn)練圖像之間相似度逐步增加，并最終達(dá)到難以區(qū)分的程度，即距離也來越小，相似度越來越高，同時(shí)對(duì)應(yīng)的參數(shù)越不敏感。

2.2 Snesim建模參數(shù)敏感性分析實(shí)例

訓(xùn)練圖像是地質(zhì)先驗(yàn)?zāi)Ｐ?，其樣式、大小等都?huì)對(duì)模擬結(jié)果產(chǎn)生重大影響，以一個(gè)網(wǎng)狀河訓(xùn)練圖像channel為例，圖3訓(xùn)練圖像大小為101×101的網(wǎng)格，黑色部分代表河道砂，白色部分代表背景泥，訓(xùn)練圖像數(shù)字化的結(jié)果是一個(gè)101×101的矩陣，“1”代表河道，“0”代表背景泥。其中背景泥占比0.7，河道占比0.3。

圖3 訓(xùn)練圖像channelFig.3 Training image channel

2.2.1 搜索節(jié)點(diǎn)數(shù)

利用Snesim基于像元的多點(diǎn)地質(zhì)建模對(duì)訓(xùn)練圖像進(jìn)行非條件模擬，為了找出最優(yōu)建模參數(shù)值，這個(gè)過程中其他參數(shù)保持不變，只對(duì)Snesim重要建模參數(shù)“搜索節(jié)點(diǎn)”數(shù)進(jìn)行更改，設(shè)置搜索節(jié)點(diǎn)數(shù)5到75，步長10。其他參數(shù)分別設(shè)置為：分類變量比例0.7、0.3，伺服器系統(tǒng)參數(shù)0.5;搜索半徑hx=50，hy=50，搜索角度angel_1=0, angel_2=0，多重網(wǎng)格層數(shù)為3。對(duì)應(yīng)每個(gè)搜素節(jié)點(diǎn)數(shù)下進(jìn)行模擬生成300個(gè)隨機(jī)模型。圖4是不同搜索節(jié)點(diǎn)數(shù)所模擬出來的一個(gè)模型(300個(gè)實(shí)現(xiàn)模型中的一個(gè))，隨著搜索節(jié)點(diǎn)數(shù)量的增加，所模擬出來的河道形態(tài)越來越好，與訓(xùn)練圖像越來越相似，當(dāng)增加到25和35以后，設(shè)置不同搜索節(jié)點(diǎn)數(shù)所模擬出來的模型的空間特征相似用肉眼就無法直接區(qū)分開來，同時(shí)隨著搜索節(jié)點(diǎn)增加，該參數(shù)敏感性降低，所生成的隨機(jī)模型效果改善不明顯甚至幾乎難以區(qū)分，往往在考慮計(jì)算效率以及實(shí)現(xiàn)效果上，通過人為肉眼觀察給定一個(gè)最優(yōu)參數(shù)值往往具有誤差，會(huì)給建模中帶來不確定性。

圖4 不同搜索節(jié)點(diǎn)數(shù)對(duì)應(yīng)的實(shí)現(xiàn)模型Fig.4 Implementation models corresponding to different node numbers(a)Node=5;(b)Node=15;(c)Node=25;(d)Node=35;(e)Node=45;(f)Node=55;(g)Node=65;(h)Node=75

根據(jù)圖2提出的優(yōu)選方法流程，計(jì)算隨機(jī)模型與訓(xùn)練圖像砂巖泥巖相在x、y方向上不同距離下的連通概率，基于連通性概率函數(shù)計(jì)算公式(4)計(jì)算訓(xùn)練圖像channel砂巖相泥巖相在x、y方向上的連通概率，圖5是訓(xùn)練圖像在x和y方向上砂巖泥巖在不同滯后距下的連通概率，隨著滯后距的增加砂巖泥巖相在x、y方向上的連通性總體上不斷趨于減小，其中砂巖的變化有周期化的上下波動(dòng)，并且其變化規(guī)律能大致反映砂巖的空間分布特征。

圖5 訓(xùn)練圖像channel在x、y方向上的連通性概率曲線Fig.5 Connectivity probability curve of training image channel in x and y directions

不同搜索節(jié)點(diǎn)所實(shí)現(xiàn)的隨機(jī)模型與訓(xùn)練圖像在空間型態(tài)結(jié)構(gòu)特征上存在差異，隨著搜索節(jié)點(diǎn)數(shù)量的增加，隨機(jī)模型的形態(tài)結(jié)構(gòu)等特征與訓(xùn)練圖像越來越相似，差異會(huì)越來越小，這種差異也同時(shí)體現(xiàn)在了連通性概率曲線上如圖6所示，其中mod1代表搜索節(jié)點(diǎn)數(shù)等于5的隨機(jī)模型(綠色)，Mod2是搜索節(jié)點(diǎn)數(shù)等于75的隨機(jī)模型(藍(lán)色)，實(shí)線與虛線分別代表泥巖相與砂巖相在x或y方向上的連通概率，由圖6可知，模型與訓(xùn)練圖像的形態(tài)特征和空間結(jié)構(gòu)越相似，二者之間的連通性概率曲線就越相近相似。(圖6曲線紅色和藍(lán)色相近，紅色與綠色差異大)

圖6 訓(xùn)練圖像及其隨機(jī)模型的連通性概率曲線對(duì)比圖Fig.6 Comparison of the connectivity probability curves of the training image and its random model(a)x方向；(a)y方向

利用余弦相似度cosSimilarity函數(shù)計(jì)算訓(xùn)練圖像與隨機(jī)模型連通性概率曲線的相似度，得到隨機(jī)模型(搜索節(jié)點(diǎn)數(shù)從5～75)與訓(xùn)練圖像關(guān)于連通性概率的余弦相似度，然后取其數(shù)學(xué)期望作為建模參數(shù)“搜索節(jié)點(diǎn)”對(duì)應(yīng)的空間相關(guān)性評(píng)價(jià)指標(biāo)ξk建立相似度與搜索節(jié)點(diǎn)數(shù)的關(guān)系曲線(圖7)，隨著搜索節(jié)點(diǎn)數(shù)從5開始增加，300個(gè)隨機(jī)模型與訓(xùn)練圖像的相似度期望不斷增加，表明模型與訓(xùn)練圖像越來越相似(根據(jù)余弦相似度定義，余弦值越接近1，表示越相似)，當(dāng)搜索節(jié)點(diǎn)數(shù)達(dá)到35時(shí)，空間相關(guān)性評(píng)價(jià)指標(biāo)進(jìn)入一個(gè)平臺(tái)區(qū)域，參數(shù)變得不敏感，模型質(zhì)量將不再有明顯提升，是本例中的優(yōu)選參數(shù)。

圖7 建模參數(shù)(搜索節(jié)點(diǎn)數(shù))基于連通性函數(shù)的余弦相似度變化曲線Fig.7 Cosine similarity change curve based on connectivity function for modeling parameters (node numbers)

2.2.2 多重網(wǎng)格層數(shù)

多重網(wǎng)格層數(shù)對(duì)模擬的效果影響很大，但是隨之而來的是消耗更大的計(jì)算性能，利用上述分析方法，可以找到一個(gè)最優(yōu)值，既能達(dá)到良好建模效果，同時(shí)也節(jié)約計(jì)算性能。對(duì)“多重網(wǎng)格重?cái)?shù)”這一參數(shù)利用相同方法測(cè)試，設(shè)置多重網(wǎng)格重?cái)?shù)從1到8。

保持其他參數(shù)不變?cè)O(shè)置搜索節(jié)點(diǎn)數(shù)35、分類變量比例0.7、0.3、伺服系統(tǒng)參數(shù)0.5，搜索半徑hx=50，hy=50、搜索角度angel_1=0, angel_2=0。利用snesim多點(diǎn)地質(zhì)建模對(duì)訓(xùn)練圖像進(jìn)行非條件模擬，在每個(gè)多重網(wǎng)格重?cái)?shù)下進(jìn)行模擬生成300個(gè)隨機(jī)模型(圖8顯示對(duì)應(yīng)網(wǎng)格層數(shù)下所模擬的300個(gè)實(shí)現(xiàn)模型中的一個(gè))。

圖8 不同網(wǎng)格重?cái)?shù)對(duì)應(yīng)的實(shí)現(xiàn)模型Fig.8 Implementation models corresponding to different multigrids(a)Multigrids=1;(b)Multigrids=2;(c)Multigrids=3;(d)Multigrids=4;(e)Multigrids=5;(f)Multigrids=6;(g)Multigrids=7;(h)Multigrids=8

利用余弦相似度cosSimilarity函數(shù)計(jì)算訓(xùn)練圖像與隨機(jī)模型連通性概率曲線的相似度，得到隨機(jī)模型(多重網(wǎng)格重?cái)?shù)1-8)與訓(xùn)練圖像關(guān)于連通性概率的余弦相似度，然后取其相似度均值建立其對(duì)應(yīng)多重網(wǎng)格重?cái)?shù)變化的曲線(圖 9)。

圖9 建模參數(shù)(網(wǎng)格重?cái)?shù))基于連通性函數(shù)的余弦相似度變化曲線Fig.9 Cosine similarity change curve based on connectivity function for modeling parameters (multigrids)

隨著多重網(wǎng)格數(shù)從1開始增加，300個(gè)隨機(jī)模型與訓(xùn)練圖像的相似度期望不斷增加，表明模型與訓(xùn)練圖像越來越相似(根據(jù)余弦相似度定義，余弦值越接近1，表示越相似)，當(dāng)網(wǎng)格重?cái)?shù)到3時(shí)，空間相關(guān)性評(píng)價(jià)指標(biāo)進(jìn)入一個(gè)平臺(tái)區(qū)域，模型質(zhì)量將不再有明顯提升，參數(shù)敏感性降低。這是因?yàn)楫?dāng)網(wǎng)格很粗時(shí)，相應(yīng)的用來掃描訓(xùn)練圖像的樣板規(guī)模會(huì)大于訓(xùn)練圖像的規(guī)模，導(dǎo)致大的數(shù)據(jù)事件很難被找到。此時(shí)數(shù)據(jù)事件就會(huì)相應(yīng)縮小，而這就等同于使用較小的數(shù)據(jù)樣板進(jìn)行掃描，一般來說，網(wǎng)格重?cái)?shù)設(shè)置多少取決于最大尺度結(jié)構(gòu)信息的大小，在本例中101×101的網(wǎng)格訓(xùn)練圖像中，利用連通性方法定量化地給出最佳網(wǎng)格重?cái)?shù)為3，超過3以后，模擬效果不會(huì)有明顯改善，該實(shí)驗(yàn)結(jié)果也符合前人所提出的認(rèn)識(shí)。

3 結(jié)論

采用基于連通性函數(shù)的多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)建模參數(shù)敏感性分析方法，分別對(duì)兩個(gè)關(guān)鍵建模參數(shù)(搜索節(jié)點(diǎn)數(shù)、多重網(wǎng)格層數(shù))進(jìn)行敏感性測(cè)試，得出以下認(rèn)識(shí)：

1)隨著搜索節(jié)點(diǎn)數(shù)和多重網(wǎng)格層數(shù)增加，模型與訓(xùn)練圖像的形態(tài)視覺特征越來越相似，利用連通性函數(shù)來量化這一空間相關(guān)性特征，結(jié)合余弦相似度度量方法對(duì)其進(jìn)行定量的分析，找出價(jià)指標(biāo)開始平穩(wěn)、進(jìn)入平臺(tái)區(qū)域時(shí)的拐點(diǎn)所對(duì)應(yīng)的參數(shù)值作為最佳參數(shù)。

2)針對(duì)本例101×101的河道訓(xùn)練圖像，在snesim模擬過程中搜索節(jié)點(diǎn)數(shù)不能過少，過少捕捉不到河道的結(jié)構(gòu)特征，過多并不能明顯改善建模效果，參數(shù)在到達(dá)一定值后變得不敏感，反而會(huì)增加計(jì)算量，降低效率，多重網(wǎng)格層數(shù)也是如此。利用連通性的方法定量化的分析參數(shù)敏感性，找到搜索節(jié)點(diǎn)數(shù)設(shè)置在35左右，多重網(wǎng)格層數(shù)設(shè)置為3層，能獲取良好模擬效果，也能減少不必要的運(yùn)算，提高效率。