郭書生 朱江梅 王世越 趙前華 王曉飛

(1. 中海石油(中國)有限公司海南分公司 海南?？?570311； 2. 中海油田服務(wù)股份有限公司油田技術(shù)事業(yè)部 河北廊坊 065201)

地震波速度是在地震資料的采集、處理、解釋和評價研究中的關(guān)鍵參數(shù)，地震波速度求取的準(zhǔn)確與否直接影響著地震勘探資料在應(yīng)用中的可信度。隨著油氣勘探程度的不斷提高，勘探區(qū)域地質(zhì)條件變得越來越復(fù)雜、速度橫縱向變化大等問題，對高精度速度預(yù)測提出了挑戰(zhàn)。為此，各種各樣的地震速度預(yù)測技術(shù)應(yīng)運(yùn)而生。這些地震速度預(yù)測方法按照數(shù)據(jù)基礎(chǔ)分為三類：①地震資料為數(shù)據(jù)基礎(chǔ)的疊加速度分析[1-3]、偏移速度分析[4-6]和層析速度反演[7-8]等速度預(yù)測方法；②以測井資料為數(shù)據(jù)基礎(chǔ)的克里金估計[9-10]、隨機(jī)模擬[11]和隨機(jī)反演[12-13]等速度預(yù)測方法；③以井震資料結(jié)合為數(shù)據(jù)基礎(chǔ)的井震匹配[14]、井震約束[15]等速度預(yù)測方法。

本文借鑒“ look ahead VSP ”思路[16]，也就是基于VSP資料進(jìn)行速度預(yù)測。由于VSP技術(shù)是將井、震融于一體的一種可靠地球物理勘探手段，可為預(yù)測高精度的垂向地震速度提供更有力支撐。而基于VSP資料的鉆前預(yù)測研究主要集中在鉆前壓力預(yù)測與深度預(yù)測方面[17-18]，而這些研究歸根還是有賴于高精度鉆前速度的預(yù)測，國內(nèi)外許多學(xué)者對此進(jìn)行了有關(guān)研究。如Khaitan等提出的基于前瞻VSP調(diào)查理論的前瞻VSP反演方法[19]；Esmersoy等提出的通過反演地震波走時估計速度的鉆前速度實(shí)時預(yù)測方法[20]；張科 等提出的根據(jù)經(jīng)驗(yàn)公式、局部構(gòu)造形態(tài)做控制結(jié)合實(shí)際鉆井計算鉆前速度的綜合速度分析方法[21]和一些基于射線追蹤理論和VSP走廊疊加反演波阻抗估算速度的常規(guī)方法[22-23]。但是這些方法都是側(cè)重于利用單一因素來預(yù)測鉆頭前方地層速度，沒有關(guān)注VSP地震資料的多屬性高頻先驗(yàn)信息。

基于以上分析，本文將從以下兩個思路展開研究：一方面，利用高分辨率VSP走廊疊加，作為高精度預(yù)測的基礎(chǔ)；另一方面，針對現(xiàn)有預(yù)測技術(shù)主要關(guān)注單一因素來預(yù)測鉆頭前方地層速度的缺陷，提出了使用多元回歸模型量化刻畫客觀因素對地層速度的影響，從而為高精度預(yù)測提供更好的算法基礎(chǔ)。具體工作步驟可歸納為：首先，對基于VSP走廊疊加資料的多屬性信息構(gòu)建多元回歸模型；其次，反演預(yù)測中頻地層速度，最后再將中頻地層速度等因素特征作為輸入，再通過粒子濾波法[24-25]來求解該多元回歸模型的狀態(tài)參數(shù)，從而有效地預(yù)測鉆頭前方地層速度。以南海鶯歌海盆地L工區(qū)L1井為目標(biāo)研究區(qū)塊，應(yīng)用提出的方法對鉆頭前200 m的地層速度進(jìn)行預(yù)測，以驗(yàn)證本文方法的有效性和穩(wěn)健性。

1 算法的提出

1.1 問題描述與模型

基于以上數(shù)據(jù)條件，我們將鉆頭前方地層速度預(yù)測本質(zhì)問題轉(zhuǎn)化為有監(jiān)督機(jī)器學(xué)習(xí)的回歸問題，已鉆地層速度數(shù)據(jù)作為回歸模型的標(biāo)簽數(shù)據(jù)。由于涉及多個特征量的輸入，所以該回歸問題又屬于多元回歸問題的范疇，如圖1所示，三維坐標(biāo)空間表示由輸入屬性張成的空間，曲面則表示擬合得到的回歸平面。

圖1 多元回歸示意圖Fig.1 Schematic diagram of the multiple regression model

對于多元回歸模型，核心是用系統(tǒng)狀態(tài)方程對問題進(jìn)行建模。根據(jù)系統(tǒng)狀態(tài)方程，這里有兩方面內(nèi)容：一是將待估計的鉆頭前方地層速度建模為一階馬爾科夫鏈，也就是待估計速度序列之間為時間相關(guān)，可以用逼近函數(shù)f(·)表示；二是輸入的多個特征序列和待預(yù)測速度之間為空間相關(guān)，可以用逼近函數(shù)g(·)表示。時間相關(guān)和空間相關(guān)可以用以下系統(tǒng)方程和狀態(tài)測量方程進(jìn)行刻畫，具體表達(dá)式如下：

狀態(tài)轉(zhuǎn)移方程為

xk=f(xk-1)+rk

(1)

狀態(tài)測量方程為

yk=g(xk)+vk

(2)

式(1)、(2)中：xk代表系統(tǒng)是tk時刻的地層速度，值得注意對于時刻tk而言此處為向量；yk代表疊加走廊和各種屬性，包括反演結(jié)果;rk、vk分別代表過程噪聲和測量噪聲。對應(yīng)于本研究，考慮建模過程是存在誤差的，rk是時變建模過程中存在的誤差，主要是線性或者非線性的方式逼近導(dǎo)致產(chǎn)生的誤差。vk是空間關(guān)系建模所導(dǎo)致的誤差，因?yàn)榉蔷€性函數(shù)g(·)需要泰勒展開線性化求解，而線性化則會省略多次項引入誤差，已鉆井段地層速度和VSP走廊疊加數(shù)據(jù)測量過程中存在的誤差。

1.2 基于反演的中頻地層速度預(yù)測

鉆頭前方未知的地層速度因其本身的復(fù)雜性，所以，我們需要分階段來逼近預(yù)測。我們先通過反演獲得鉆頭前方地層速度的中頻信息，本質(zhì)上中頻信息帶寬與走廊疊加資料的帶寬接近。具體來說VSP走廊疊加反演需要在正演基礎(chǔ)上才能夠進(jìn)行，當(dāng)然正演模型在進(jìn)行輸入數(shù)據(jù)井震匹配時也是需要的。對于VSP走廊疊加正演模型，本文采用了常用的褶積模型進(jìn)行構(gòu)建。該模型定義如下：單道VSP走廊疊加S(t)可以建模為地震子波W(t)與反射系數(shù)R(t)卷積的結(jié)果：

S(t)=W(t)*R(t)

(3)

由于ΔZi遠(yuǎn)小于阻抗Zi，因此反射系數(shù)可近似用式(4)表示：

(4)

式(4)中：Zi=ρV，ΔZi為不同深度之間阻抗變化。在實(shí)際砂泥巖的預(yù)測過程中，因?yàn)槊芏圈训淖兓秶^小，所以可以將密度ρ定義為常值，這樣就可以利用Zi求取速度地層速度V。因此，式(3)可以表示為

(5)

式(5)中：St表示時刻t的走廊疊加資料;wt表示時刻t的子波。在實(shí)踐中，反演問題的目的是根據(jù)S和G的知識給出X的估計。然而，由于這個問題的不適定性，方程大部分是使用最小平方正則化，其可表示如下：

(6)

f(x)=‖x‖1

這里‖·‖1為L1范數(shù)，該目標(biāo)函數(shù)的求解有許多求解方法，比如共軛梯度[26]等，鑒于篇幅不在這里贅述。

1.3 完整預(yù)測算法

為了進(jìn)一步提升地層速度預(yù)測效果，需要將反演的中頻地層速度和屬性數(shù)據(jù)一起輸入多元回歸模型，再利用粒子濾波算法來求解多元回歸模型。作為典型的回歸模型，它可以具體分成2個階段，第1個階段是訓(xùn)練階段，也就是已知地層速度情況下，訓(xùn)練1個多元回歸模型；第2個階段是預(yù)測階段，輸入新的反演速度和屬性數(shù)據(jù)，預(yù)測新的地層速度。

在粒子濾波算法中，對于一個動態(tài)非線性離散時間系統(tǒng)，其狀態(tài)空間模型通常寫為：

狀態(tài)轉(zhuǎn)移方程為

xk,1=Axk-1,1+rk,1

(7)

狀態(tài)測量方程為

yk=g1(xk,1)+vk,1

(8)

式(7)、(8)中：A表示將公式(1)中f(·)線性化，具體通過歷史數(shù)據(jù)擬合得到，所謂歷史數(shù)據(jù)是指已鉆探段采集得到的速度；xk,1=[dk,hk]T代表在第k個時刻地層速度預(yù)測所對應(yīng)的模型參數(shù)；rk,1、vk,1分別代表參數(shù)狀態(tài)轉(zhuǎn)移的過程噪聲和地層速度的測量噪聲。

在粒子濾波算法中，系統(tǒng)狀態(tài)的概率分布估計是通過一組具有權(quán)重的隨機(jī)樣本(稱為粒子)給出的。對于一個粒子數(shù)是ND的粒子濾波算法，具體計算步驟如下:

圖2 粒子濾波器示意圖[24]Fig.2 Schematic diagram of particle filter[24]

在利用粒子濾波估計出速度預(yù)測的多元回歸模型的狀態(tài)參數(shù)后，將根據(jù)式(7)、(8)可直接對鉆頭前方地層速度進(jìn)行預(yù)測。

當(dāng)然，在不同的地質(zhì)環(huán)境下，影響鉆頭前方地層速度的客觀因素并不僅僅只有反演結(jié)果和疊加走廊差分衍生屬性，本文提出的方法作為一個開放結(jié)構(gòu)，在未來研究工作中，更多樣和密切的客觀因素可以以系統(tǒng)狀態(tài)方程的特征輸入到我們的方法中，進(jìn)一步提升鉆頭前方地層速度的預(yù)測精度。

1.4 具體算法步驟

本文提出的算法具體步驟如下：

輸入：已鉆井段地層速度X，VSP走廊疊加Y

2) 求解公式(6)的中頻地層速度x；

3) 將以上地震屬性和中頻地層速度輸入回歸模型，利用粒子濾波器進(jìn)行求解；

②根據(jù)粒子的重要度權(quán)值重采樣；

輸出：yk

2 應(yīng)用實(shí)例分析

本文以鶯歌海盆地L工區(qū)L1井為目標(biāo)研究區(qū)塊，該井為高溫超高壓井，目的層溫度為179～188 ℃。該井在鉆到3 834.2 m后，根據(jù)鄰井的信息，認(rèn)為前方將會鉆遇超高壓地層的臺階面，希望將當(dāng)前井段的中完在超高壓地層的臺階面之上，越接近超高壓地層的臺階面越好，因此決定測量中途VSP以標(biāo)定層位深度及預(yù)測下部井段的地層信息。VSP資料采集的最大深度為3 826 m。圖3顯示了基于多元回歸模型的鉆頭前方地層速度預(yù)測的結(jié)果。圖中第2列為錄井巖性剖面(后續(xù)實(shí)鉆)及已知井段通過VSP測量時深關(guān)系轉(zhuǎn)換的地層速度，第3列是垂直地震走廊疊加道；第4列是常規(guī)地震的井旁道；第5列是后續(xù)實(shí)鉆后測量的伽馬測井曲線。第6列是后續(xù)實(shí)鉆后測量的聲波測井曲線(洋紅色)及利用本文方法預(yù)測的獲得速度(黑色)和常規(guī)波阻抗反演獲得的速度(綠色)。通過圖3可以看出，首先垂直地震走廊疊加道(第3列)的分辨率明顯比常規(guī)地震的井旁道(第4列)高，有助于鉆頭前方的砂體的識別。其次，在第6列反演預(yù)測地層速度和已鉆井段實(shí)測地層速度在已鉆井段(3 826 m以上)非常一致。未鉆井段(3 826 m以下)經(jīng)過后續(xù)實(shí)鉆表明，此次預(yù)測的速度趨勢與實(shí)測數(shù)據(jù)相近。

圖3 垂直地震走廊疊加預(yù)測速度與實(shí)鉆結(jié)果對比Fig.3 Vertical seismic data corridor stack section predicted velocity vs. actual drilling

根據(jù)此結(jié)果對目的層進(jìn)行預(yù)測，見表1鉆前預(yù)測部分。本井實(shí)時鉆前深度預(yù)測結(jié)果與后續(xù)實(shí)鉆結(jié)果有很好的一致性，對目的層T31，T31-A1深度的預(yù)測和實(shí)鉆深度誤差0 m，對目的層T31-B1的預(yù)測和實(shí)鉆深度差為-3 m，對目的層T31-C1的預(yù)測和實(shí)鉆深度差為12 m。對于目的層T31-C1，從錄井資料和伽瑪曲線分析，該層為砂泥巖互層，波阻抗變化比較頻繁，而VSP資料的縱向分辨率不及電纜測井，所以對急速變化的波阻抗特征反應(yīng)不及時，加上預(yù)測距離增大等因素的影響，導(dǎo)致預(yù)測誤差較大。

表1 目的層預(yù)測深度與實(shí)鉆結(jié)果對比Table 1 Predicted results of the depth of the target layer

3 結(jié)論

本文提出了基于多元回歸模型的鉆頭前方地層速度預(yù)測方法，該方法通過建立非線性多元回歸模型，充分利用了VSP資料中的多屬性的信息，也即疊加走廊差分衍生屬性(比如振幅導(dǎo)數(shù)、瞬時振幅導(dǎo)數(shù)等)，并使用非線性粒子濾波法對模型進(jìn)行求解。通過在鶯歌海盆地L工區(qū)L1井目標(biāo)研究區(qū)塊的實(shí)際應(yīng)用表明，本文方法反演預(yù)測的速度比常規(guī)波阻抗反演獲得的速度頻率獲得了提高，同時速度變化的趨勢與實(shí)測數(shù)據(jù)更為相近，驗(yàn)證了本文方法的有效性。