陳文浩, 王志章, 董少群, 侯加根

(1.中國石油大學(xué)(北京)地球科學(xué)學(xué)院, 北京 102249;2.中國石油大學(xué)(北京)油氣資源與探測國家重點(diǎn)實(shí)驗(yàn)室, 北京 102249)

0 引 言

常規(guī)砂巖儲(chǔ)層通常依據(jù)聲波線性解釋儲(chǔ)層的孔隙度,其合理解釋是儲(chǔ)量估算及儲(chǔ)層評(píng)價(jià)工作的重點(diǎn)和難點(diǎn)。對(duì)于致密砂巖儲(chǔ)層,聲波與孔隙度不具線性關(guān)系,已有的孔隙度線性解釋方法不適于解釋致密砂巖儲(chǔ)層的孔隙度。王洪輝等[1]認(rèn)為簡單的線性解釋致密砂巖儲(chǔ)層孔隙度難以滿足致密砂巖儲(chǔ)層評(píng)價(jià)的要求。張永浩等[2]模擬致密砂巖儲(chǔ)層條件時(shí),對(duì)聲波速度與孔隙度進(jìn)行實(shí)驗(yàn)分析,發(fā)現(xiàn)孔隙度隨壓力增大呈冪函數(shù)減小,而聲波速度隨壓力增大呈冪函數(shù)增大,所以孔隙度與聲波速度的線性關(guān)系隨壓力增大而變差。常規(guī)砂巖儲(chǔ)層孔隙度線性解釋模型不適于致密砂巖儲(chǔ)層的孔隙度解釋。與核函數(shù)融合形成的核嶺回歸算法[3-6]不僅可以有效處理非線性問題,而且具有參數(shù)少、建模時(shí)間短、運(yùn)算速度快等優(yōu)點(diǎn)[7],具有廣泛的應(yīng)用空間,但在油氣儲(chǔ)層預(yù)測中的應(yīng)用研究較少。本文闡述核嶺回歸建立致密砂巖孔隙度解釋模型,并用來解釋致密砂巖儲(chǔ)層的孔隙度,分析多種方法的解釋結(jié)果,認(rèn)為核嶺回歸是建立致密砂巖儲(chǔ)層孔隙度解釋模型的有效方法。

1 核嶺回歸方法概述

線性回歸中的欠定問題,最小二乘回歸的效果較差。為此,眾多學(xué)者提出了最小二乘回歸的改進(jìn)方法,其中比較經(jīng)典的是Hoerl和Kennard[8-9]提出的嶺回歸算法。之后Saunders等嘗試用核技巧改進(jìn)嶺回歸算法,核嶺回歸(Kernel Ridge Regression,KRR)由Cristianini和Shawe-Taylor正式使用[10],是簡化的支持向量回歸機(jī),相對(duì)于支持向量回歸機(jī),具有引用參數(shù)少、節(jié)省機(jī)時(shí)等優(yōu)點(diǎn)。

1.1 嶺回歸算法

設(shè)有n個(gè)樣品,每個(gè)樣品有m個(gè)自變量和1個(gè)因變量,它們的觀測值構(gòu)成數(shù)據(jù)矩陣X和Y,即

X=x11x12…x1n

x21x22…x2n

???

xm1xm2…xmn

Y=y1y2…yn

(1)

式中,xij(i=1,2,…,m;j=1,2,…,n)為j樣品第i個(gè)自變量的觀測值;yj(j=1,2,…,n)為因變量的第j個(gè)觀測值。

最小二乘法線性回歸依據(jù)式(1)數(shù)據(jù),求解線性回歸模型y=b0+∑mi=1bixi中待定常數(shù)與待定系數(shù)b0、b1、…、bm的最佳估計(jì)值a0、a1、…、am,得線性回歸方程

j=A·Xj+a0(j=1,2,…,n)

(2)

式中,A=a1a2…am為線性回歸方程的系數(shù)向量;Xj為矩陣X的第j元素構(gòu)成的列向量;A·Xj為二者的內(nèi)積;a0為常數(shù)。求解線性回歸方程常數(shù)和系數(shù)的原則是使損失

Q=∑nj=1(yj-j)2

達(dá)到最小。

對(duì)于線性回歸中的欠定問題,在損失最小原則下求解A和a0得到的方程會(huì)是病態(tài)方程。為解決該問題,Hoerl和Kennard引入了損失函數(shù)

L(A,a0)=γA·A+Q

(3)

式中,A·A為正則化項(xiàng);Q為擬合損失項(xiàng);γ為正則化參數(shù),γ>0,它對(duì)A·A和Q起平衡作用,稱這種有權(quán)重因子的回歸算法為嶺回歸算法。若式(3)中γ=0,則為最小二乘意義下的回歸分析。

對(duì)于標(biāo)準(zhǔn)差標(biāo)準(zhǔn)化變量,式(3)中a0等于0。因此。可將損失函數(shù)寫為矩陣形式

L(A)=γAAT+(Y-AX)(Y-AX)T

(4)

對(duì)式(4)中參數(shù)A求偏導(dǎo),并令其為0,可以得到嶺回歸模型的待定系數(shù)

A=(XXT+γI)-1XYT

(5)

式中,I為m×m的單位矩陣。

線性嶺回歸算法中正則化項(xiàng)的引入使回歸方程具有較好的穩(wěn)定性和泛化能力。但是,這種算法不能有效處理自變量間存在非線性相關(guān)的情況,這是線性回歸算法普遍存在的一個(gè)局限性。實(shí)際應(yīng)用中,用核技巧改進(jìn)后的嶺回歸算法,即核嶺回歸算法可以有效解決上述問題。

1.2 核嶺回歸算法

1.2.1 核嶺回歸算法基本原理

將原始空間Rm中的n個(gè)樣品用某一非線性映射函數(shù)φ映射到一個(gè)高維特征空間F中,樣本數(shù)據(jù)變?yōu)棣?Xj),(j=1,2,…,n)。映射后可以在F中進(jìn)行嶺回歸,回歸方程的形式與式(2)類似,即

j=A·φ(Xj)+a0

(6)

相應(yīng)的損失函數(shù)也有類似的形式

L(A,a0)=γA·A+Q

(7)

對(duì)于標(biāo)準(zhǔn)差標(biāo)準(zhǔn)化變量,損失函數(shù)的矩陣形式為

L(A)=γAAT+(Y-AX)(Y-AX)T

(8)

式中,X=φ(X1)φ(X2)…φ(Xn),權(quán)系數(shù)向量為A=a1a2…。對(duì)式(8)中A求偏導(dǎo),并令其為0,得

λA+XTAX=XTY

(9)

因非線性映射φ未知,故X、A未知。在實(shí)際應(yīng)用中一般求不出A的顯式表達(dá)式。將式(9)重新改寫后可得

A=λ-1XT(Y-AX)=XTα

(10)

其中,α=λ-1(Y-AX)。式(10)表明,A可以變換成樣品點(diǎn)的線性組合,α又可改寫為

α=(XXT+λI)-1Y=(K+λI)-1Y

(11)

其中,K稱為核矩陣,kij=φ(Xi)·φ(Xj)為核函數(shù),記為K(Xi,Xj)

對(duì)于新的樣品Xi,y的估計(jì)值

i=α·φ(Xj)=YT(K+λI)-1k

其中,k是n維的列向量,向量的元素為φ(Xi)與φ(Xj)的內(nèi)積。

通過上述變換,可以避免直接顯式地計(jì)算映射φ(xi),解決了非線性回歸的問題。

1.2.2 核嶺回歸方法中的參數(shù)

(1) 核函數(shù)中參數(shù)。應(yīng)用核嶺回歸算法,首先要確定核函數(shù),總體思路是將輸入向量表示為內(nèi)積的形式。目前,應(yīng)用較廣泛的核函數(shù)有3種類型。

①多項(xiàng)式核函數(shù)K(x,y)d,d為常數(shù);

②神經(jīng)網(wǎng)絡(luò)核函數(shù)K(x,y)=tanh(β0(x,y)+β1)

③高斯徑向基核函數(shù)K(x,y)=e-‖x-y‖2σ2,σ為控制著2個(gè)樣品的相似程度的參數(shù)。

高斯徑向基核函數(shù)具有良好的泛化性能[11],本文將使用高斯徑向基核函數(shù),σ過小時(shí)易出現(xiàn)過渡擬合,過大時(shí)又會(huì)使擬合精度下降,因此,σ的取值涉及算法的擬合精度。

(2) KRR中正則化參數(shù)。對(duì)于核嶺回歸,除了核函數(shù),還需要確定正則化參數(shù)γ。從均方根誤差圖(見圖1)可知,正則化參數(shù)γ過大或過小都將減低擬合精度。

圖1 λ對(duì)擬合誤差的影響

2 實(shí)例研究

以紅崗油田H90區(qū)塊致密砂巖儲(chǔ)層為例進(jìn)行孔隙度解釋研究,該儲(chǔ)層主要分布于扶余油層(泉四段)。區(qū)內(nèi)發(fā)育辮狀河三角洲平原和前緣亞相,主要砂體為平原亞相的分流河道砂、溢岸薄層砂及前緣亞相的水下分流河道砂。主要巖性為長石巖屑砂巖,多中砂質(zhì)細(xì)粒。泉四段油層物性較差,孔隙度為3%～13%,最高14.6%,平均8.03%;滲透率多為0.024～1 mD*非法定計(jì)量單位,1 mD=9.87×10-4 μm2; 1 ft=12 in=0.304 8 m; 下同,最高3.6 mD,平均0.294 mD,為低孔隙度超低滲透率儲(chǔ)層。

2.1 不同孔隙度解釋方法的解釋精度

2.1.1 核嶺回歸

(1) 樣本曲線。選擇與孔隙度相關(guān)密切的樣本曲線,以減少輸入變量,提高訓(xùn)練速度和精度。該地區(qū)三孔隙度測井響應(yīng)與孔隙度關(guān)系較為密切,其中聲波時(shí)差A(yù)C與孔隙度的相關(guān)系數(shù)為0.716 2,中子CNL與孔隙度的相關(guān)系數(shù)為0.664 1,密度DEN與孔隙度的相關(guān)系數(shù)為-0.730 1(見圖2)。

(2) 核嶺回歸解釋孔隙度。以三孔隙度測井曲線作為孔隙度解釋模型的輸入。利用核嶺回歸解釋孔隙度,還需要確定模型中的正則化參數(shù)γ和核函數(shù)參數(shù)σ。對(duì)于這2個(gè)參數(shù)的選取,還沒有統(tǒng)一公認(rèn)的最優(yōu)方法,目前常用網(wǎng)格搜索法和交叉驗(yàn)證法進(jìn)行優(yōu)選(見圖3)。

參數(shù)優(yōu)選的實(shí)現(xiàn)過程是給出γ和σ的取值區(qū)間,對(duì)于取定的γ和σ,把訓(xùn)練集作為原始數(shù)據(jù)集利用K-CV方法得到該組γ和σ下訓(xùn)練集預(yù)測結(jié)果,根據(jù)預(yù)測精度選擇γ和σ。網(wǎng)格搜索即在選定的范圍內(nèi),讓?duì)煤挺页手笖?shù)增長,首先選取粗網(wǎng)格進(jìn)行搜索,得到一組最優(yōu)參數(shù);然后逐步減小步長,最終獲得所需參數(shù)。圖3中紅點(diǎn)對(duì)應(yīng)解釋模型的最優(yōu)參數(shù),其中γ=0.064 8,σ=1.393 7。

圖2 三孔隙度參數(shù)與孔隙度交會(huì)圖

圖3 網(wǎng)格搜索法與交叉驗(yàn)證法確定γ和σ

2.1.2 解釋孔隙度精度比較

為了比較不同孔隙度解釋方法的解釋精度,分別采用線性回歸、逐步回歸、BP神經(jīng)網(wǎng)絡(luò)、支持向量回歸機(jī)方法解釋孔隙度,與核嶺回解釋孔隙度精度比較(見表1)。

表1 不同孔隙度解釋方法的解釋精度統(tǒng)計(jì)表

(1) 一元線性回歸解釋孔隙度模型分別為

φ=0.140AC-22.37

φ=0.664CNL+2.077

φ=-17.72DEN+53.94

(2) 多元線性回歸解釋孔隙度解釋模型為

φ=0.0753AC+0.0040CNL-9.8790DEN+12.8709

(3) 逐步回歸解釋孔隙度模型為

φ=0.076AC-10.1711DEN+17.8858

φ=-0.0013AC2-35.1002DEN2+

0.6739AC+0.1996CNL+169.809DEN

(4) BP神經(jīng)網(wǎng)絡(luò)解釋孔隙度

(5) 支持向量回歸機(jī)方法解釋孔隙度

根據(jù)表1各種方法的均方根誤差值可以看出,線性回歸解釋孔隙度方法雖然快捷,但是預(yù)測精度相對(duì)較差,4種非線性方法預(yù)測結(jié)果精度提升明顯,但相對(duì)而言,KRR方法的預(yù)測精度要略高于這些方法。表明針對(duì)致密砂巖儲(chǔ)層孔隙度的測井解釋,采用非線性算法提取的特征比線性算法提取的特征包含更準(zhǔn)確的信息,證明致密砂巖孔隙度解釋呈現(xiàn)非線性趨勢大于呈現(xiàn)線性趨勢。

為了對(duì)比3種非線性模型的執(zhí)行效率,對(duì)訓(xùn)練時(shí)間進(jìn)行了統(tǒng)計(jì)。表1中運(yùn)算時(shí)間為模型訓(xùn)練時(shí)間,統(tǒng)計(jì)結(jié)果在i5 2.0 GHz CPU、1 G RAM、MATLAB R2009a環(huán)境下得到。從表1中運(yùn)算時(shí)間可以看出支持向量回歸算法需要求解二次規(guī)劃問題,運(yùn)算時(shí)間大大高于只需進(jìn)行線性方程組求解的核嶺回歸算法;同樣,神經(jīng)網(wǎng)絡(luò)模型訓(xùn)練耗時(shí)1.9 s,核嶺回歸算法僅需0.5 s。核嶺回歸算法具有參數(shù)少、訓(xùn)練時(shí)間短、計(jì)算效率高的特點(diǎn)。

圖4 不同孔隙度解釋方法對(duì)H90井孔隙度解釋效果對(duì)比

由不同孔隙度解釋方法解釋H90井的孔隙度圖(見圖4)可以看出,單參數(shù)線性回歸和多參數(shù)線性回歸解釋的上部砂體的孔隙度均大于真實(shí)的孔隙度,3種非線性方法解釋的孔隙度與樣品孔隙度較接近。

3 結(jié) 論

(1) 核嶺回歸擬合方法能較好地解決非線性問題。對(duì)于致密砂巖儲(chǔ)層,采用KRR解釋孔隙度具有更好的預(yù)測效果。

(2) 基于核函數(shù)的非線性方法方法建立嶺回歸解釋模型能夠盡可能地提取測井信息與儲(chǔ)層孔隙度之間的非線性映射關(guān)系,對(duì)其他地區(qū)致密砂巖儲(chǔ)層的孔隙度解釋具有一定的參考作用。

參考文獻(xiàn):

[1] 王洪輝, 黎鵬, 段新國. 四川盆地須家河組低孔致密砂巖孔隙度測井解釋研究 [J]. 成都理工大學(xué)學(xué)報(bào): 自然科學(xué)版, 2009(3): 249-252.

[2] 張永浩, 杜環(huán)虹, 李新, 等. 蘇里格氣田致密砂巖儲(chǔ)層條件下聲波速度與孔隙度實(shí)驗(yàn)研究 [J]. 測井技術(shù), 2013, 37(3): 229-234.

[3] Vladimir N V. The Nature of Statistical Learning Theory [J]. New York: Springer-Verlag, 1995: 20-30.

[4] 印興耀, 孔國英, 張廣智. 基于核主成分分析的地震屬性優(yōu)化方法及應(yīng)用 [J]. 石油地球物理勘探, 2008(2): 179-183.

[5] 唐耀華, 張向君, 高靜懷. 基于地震屬性優(yōu)選與支持向量機(jī)的油氣預(yù)測方法 [J]. 石油地球物理勘探, 2009(1): 75-80.

[6] 鐘儀華, 李榕. 基于主成分分析的最小二乘支持向量機(jī)巖性識(shí)別方法 [J]. 測井技術(shù), 2009, 33(5): 425-429.

[7] 李琦, 邵誠. 基于核嶺回歸的非線性系統(tǒng)辨識(shí)及其應(yīng)用 [J]. 系統(tǒng)仿真學(xué)報(bào), 2009(8): 2152-2155.

[8] Arthur E Hoerl, Robert W Kennard. Ridge Regression: Biased Estimation for Nonorthogonal Problems [J]. Technometrics, 1970, 12(1): 55-67.

[9] Arthur E Hoerl, Robert W Kennard. Ridge Regression: Applications to Nonorthogonal Problems [J]. Technometrics, 1970, 12(1): 69-82.

[10] Nello Cristianini, John Shawe-taylor. An Introduction to Support Vector Machines and Other Kernel-based Learning Methods [M]. Cambridge: Cambridge University Press, 2000.

[11] Kristiaan Pelckmans, Jos De Brabanter, J A K Suykens. The Differogram: Non-parametric Noise Variance Estimation and its Use for Model Selection [J]. Neurocomputing (S0925-2312), 2005, 69(3): 100-122.