馬天壽，向國(guó)富，石榆帆，桂俊川，張東洋

1 西南石油大學(xué)“油氣藏地質(zhì)及開發(fā)工程”國(guó)家重點(diǎn)實(shí)驗(yàn)室, 成都 610500

2 西南石油大學(xué)工程學(xué)院, 南充 637001

3 中國(guó)石油西南油氣田公司頁(yè)巖氣研究院, 成都 610051

0 引言

準(zhǔn)確預(yù)測(cè)地應(yīng)力對(duì)于油氣勘探開發(fā)具有重要作用和意義[1]。地應(yīng)力是井壁穩(wěn)定分析的邊界條件，其大小決定了維護(hù)井壁穩(wěn)定所需的當(dāng)量泥漿密度[2-4]，同時(shí)地應(yīng)力也是水力壓裂設(shè)計(jì)的關(guān)鍵參數(shù)，其大小直接約束了水力裂縫的起裂和擴(kuò)展[5-6]。地應(yīng)力由垂向地應(yīng)力、最大和最小水平地應(yīng)力構(gòu)成，其中垂向地應(yīng)力可直接由密度測(cè)井資料計(jì)算得到，而最大和最小水平地應(yīng)力由于地質(zhì)環(huán)境和技術(shù)水平的限制難以實(shí)現(xiàn)高效預(yù)測(cè)。目前，確定最大和最小水平地應(yīng)力的方法包括直接測(cè)量法、測(cè)井解釋法、數(shù)值模擬法及地震預(yù)測(cè)法等[7]。其中，直接測(cè)量法和測(cè)井解釋法精度較高且技術(shù)相對(duì)成熟因此被廣泛采用。直接測(cè)量法包括聲發(fā)射測(cè)試、波速各向異性測(cè)試、水力壓裂測(cè)試、診斷裂縫注入測(cè)試等[8-10]，直接測(cè)量法精度最好，但是測(cè)試成本高、耗時(shí)長(zhǎng)，而且只能獲得取芯深度點(diǎn)的地應(yīng)力，所測(cè)得的數(shù)據(jù)十分有限；測(cè)井解釋法可以利用聲波測(cè)井、成像測(cè)井及鉆井誘導(dǎo)裂縫等資料確定地應(yīng)力[11-12]，相較于直接測(cè)量法，測(cè)井解釋法精度偏低，但是具有縱向分辨率高、相對(duì)連續(xù)的特點(diǎn)，可以解釋得到相對(duì)連續(xù)的地應(yīng)力剖面，因而被現(xiàn)場(chǎng)廣泛應(yīng)用。然而，目前的測(cè)井解釋方法依賴于陣列聲波測(cè)井資料和實(shí)測(cè)地應(yīng)力數(shù)據(jù)，嚴(yán)重限制了地應(yīng)力的測(cè)井解釋。此外，由于地下情況復(fù)雜和非均質(zhì)性較強(qiáng)，測(cè)井?dāng)?shù)據(jù)之間表現(xiàn)出極強(qiáng)的非線性關(guān)系[13]，使得傳統(tǒng)的解析模型和神經(jīng)網(wǎng)絡(luò)模型難以表達(dá)測(cè)井?dāng)?shù)據(jù)與地應(yīng)力之間的空間相關(guān)性。因此，亟需尋求一種簡(jiǎn)單、高效的替代模型用于解決地應(yīng)力預(yù)測(cè)問(wèn)題。

近年來(lái)，隨著機(jī)器學(xué)習(xí)在科學(xué)和工程領(lǐng)域的廣泛應(yīng)用，許多學(xué)者提出利用機(jī)器學(xué)習(xí)來(lái)解決地球科學(xué)問(wèn)題[14-15]。張東曉等[15]提出了一種串級(jí)長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)，用于補(bǔ)全和生成人工測(cè)井曲線；Singh[16]提出了一種模糊邏輯模型進(jìn)行地層識(shí)別；Mehrgini等[17]利用Elman神經(jīng)網(wǎng)絡(luò)，通過(guò)訓(xùn)練自然伽馬、電阻率、中子孔隙度、巖石體積密度及含水飽和度等參數(shù)來(lái)預(yù)測(cè)橫波速度；Korjani等[18]基于鄰井地質(zhì)資料，構(gòu)建了一種深度學(xué)習(xí)神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)巖石物理特征的油藏建模方法，用于構(gòu)建缺少測(cè)井和巖心資料的虛擬測(cè)井?dāng)?shù)據(jù)；韓玉嬌等[19]采用AdaBoost機(jī)器學(xué)習(xí)算法，開展了大牛地氣田儲(chǔ)層流體智能識(shí)別研究。在地應(yīng)力預(yù)測(cè)方面：馮鵬等[20]篩選了井徑、補(bǔ)償中子、自然伽馬、密度、深淺側(cè)向電阻率等測(cè)井參數(shù)作為訓(xùn)練因子，利用支持向量回歸方法建立了煤層最小水平地應(yīng)力預(yù)測(cè)模型；尚福華等[21]基于改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)模型對(duì)巖石力學(xué)參數(shù)進(jìn)行反演預(yù)測(cè)，進(jìn)而進(jìn)行地應(yīng)力預(yù)測(cè)；HAN和YIN[22]利用人工神經(jīng)網(wǎng)絡(luò)(ANN)，將井眼形變數(shù)據(jù)作為特征來(lái)預(yù)測(cè)地應(yīng)力；IBRAHIM等[23]基于測(cè)井參數(shù)利用隨機(jī)森林、功能網(wǎng)絡(luò)和模糊推理系統(tǒng)等機(jī)器學(xué)習(xí)技術(shù)來(lái)預(yù)測(cè)地應(yīng)力；GOWIDA等[24]利用ANN，通過(guò)訓(xùn)練鉆壓、扭矩及機(jī)械鉆速等實(shí)鉆數(shù)據(jù)來(lái)預(yù)測(cè)地應(yīng)力；TUNG[25]利用ANN，通過(guò)訓(xùn)練垂深、孔隙壓力、垂向地應(yīng)力等數(shù)據(jù)來(lái)預(yù)測(cè)最小水平地應(yīng)力；PU等[26]通過(guò)決策樹，綜合考慮巖石彈性模量、泊松比及其他參數(shù)的影響來(lái)預(yù)測(cè)地應(yīng)力；MA等[27]建立了一種混合神經(jīng)網(wǎng)絡(luò)模型，用于單井地應(yīng)力的預(yù)測(cè)。

國(guó)內(nèi)外學(xué)者已經(jīng)利用機(jī)器學(xué)方法對(duì)巖性識(shí)別、測(cè)井曲線重構(gòu)、地應(yīng)力預(yù)測(cè)等開展了大量研究，但地應(yīng)力的機(jī)器學(xué)習(xí)方法大都采用傳統(tǒng)全連接神經(jīng)網(wǎng)絡(luò)(FCNN)，這類方法只能構(gòu)造點(diǎn)對(duì)點(diǎn)的映射，無(wú)法有效捕捉測(cè)井資料深度上的變異信息，而且，這類方法通常假設(shè)數(shù)據(jù)是低維和線性的，因而難以有效準(zhǔn)確地預(yù)測(cè)地應(yīng)力。此外，這類方法的泛化能力較差，很難用于新井地應(yīng)力的預(yù)測(cè)。為此，采用一種先進(jìn)的循環(huán)神經(jīng)網(wǎng)絡(luò)，即雙向長(zhǎng)短記憶神經(jīng)網(wǎng)絡(luò)(BiLSTM)，利用常規(guī)測(cè)井參數(shù)來(lái)預(yù)測(cè)地應(yīng)力。首先，選取工區(qū)兩口直井作為研究對(duì)象，分別對(duì)兩口井的測(cè)井資料進(jìn)行地應(yīng)力解釋，形成訓(xùn)練和測(cè)試數(shù)據(jù)集；然后，研究數(shù)據(jù)預(yù)處理與樣本構(gòu)造算法，分析測(cè)井參數(shù)與地應(yīng)力相關(guān)性，建立基于BiLSTM的水平地應(yīng)力預(yù)測(cè)模型，并闡述模型的核心算法和流程；最后，結(jié)合測(cè)井參數(shù)相關(guān)性與地質(zhì)學(xué)含義，分析不同測(cè)井參數(shù)組合模式下水平地應(yīng)力的預(yù)測(cè)效果。本文研究對(duì)于準(zhǔn)確預(yù)測(cè)地應(yīng)力、推廣機(jī)器學(xué)習(xí)在石油工程的應(yīng)用具有重要的作用和意義。

1 地應(yīng)力測(cè)井解釋

本文數(shù)據(jù)來(lái)自于四川盆地CL氣田的2口直井，工區(qū)位于川中隆起區(qū)的川西南低陡褶帶，該工區(qū)整體表現(xiàn)為由北西向南、東方向傾斜的寬緩單斜構(gòu)造，且地層平緩，傾角小，斷裂不發(fā)育。為了獲得2口井的地應(yīng)力測(cè)井解釋數(shù)據(jù)集，需要根據(jù)工區(qū)地層特征建立地質(zhì)力學(xué)模型。

1.1 地質(zhì)力學(xué)模型

地質(zhì)力學(xué)模型是根據(jù)室內(nèi)巖石力學(xué)實(shí)驗(yàn)、礦場(chǎng)實(shí)測(cè)資料和測(cè)井解釋建立地層巖石力學(xué)參數(shù)、孔隙壓力、地應(yīng)力一維剖面的常用方法[27-28]，其建模流程如圖1所示。首先，確定巖石力學(xué)參數(shù)，包括泊松比、彈性模量和Biot系數(shù)；其次，通過(guò)密度測(cè)井資料計(jì)算垂向地應(yīng)力；然后，利用聲波測(cè)井資料確定孔隙壓力；進(jìn)一步，根據(jù)室內(nèi)測(cè)試或小型地破試驗(yàn)結(jié)果，反演構(gòu)造應(yīng)力系數(shù)；最后，利用巖石力學(xué)參數(shù)、孔隙壓力、垂向地應(yīng)力解釋結(jié)果和構(gòu)造應(yīng)力系數(shù)，確定最大和最小水平地應(yīng)力。

圖1 地質(zhì)力學(xué)模型建模流程圖[27]Fig. 1 Flowchart of geo-mechanical modeling[27]

大部分油氣測(cè)井作業(yè)中，由于成本限制，一般只進(jìn)行縱波時(shí)差測(cè)井，缺少橫波時(shí)差測(cè)井資料[29]。因此，需要采用橫波重構(gòu)來(lái)估算橫波時(shí)差，即利用工區(qū)已測(cè)井的縱橫波時(shí)差關(guān)系，建立橫波時(shí)差預(yù)測(cè)的經(jīng)驗(yàn)公式。對(duì)于目標(biāo)工區(qū)，縱橫波時(shí)差之間的經(jīng)驗(yàn)關(guān)系為[27]：

式中：Δtp為縱波時(shí)差，μs/ft；Δts為橫波時(shí)差，μs/ft。

進(jìn)一步，根據(jù)彈性介質(zhì)中縱橫波傳播理論，建立動(dòng)態(tài)彈性參數(shù)與縱波、橫波時(shí)差的關(guān)系[29]：

式中：Ed為動(dòng)態(tài)彈性模量，GPa；μd為動(dòng)態(tài)泊松比，無(wú)因次；ρb為巖石密度，g/cm3。

在實(shí)際工程應(yīng)用中，需要采用靜態(tài)彈性模量和泊松比。因此，需要對(duì)彈性模量和泊松比進(jìn)行動(dòng)靜態(tài)轉(zhuǎn)換。對(duì)于目標(biāo)工區(qū)，彈性模量和泊松比的動(dòng)靜態(tài)轉(zhuǎn)換關(guān)系為[27]：

式中：Es為靜態(tài)彈性模量，GPa；μs為靜態(tài)泊松比，無(wú)因次。

而地層巖石Biot系數(shù)可表示為[29]：

式中：α為Biot系數(shù)，無(wú)因次；ρma為巖石骨架密度，g/cm3；Δtmap和Δtmas分別為巖石骨架的縱波時(shí)差和橫波時(shí)差，μs/ft。

垂向地應(yīng)力也稱上覆巖層壓力，可由密度測(cè)井曲線直接積分求得：

式中：σv為垂向地應(yīng)力，MPa；ρ為缺少密度測(cè)井資料井段的上覆巖層平均密度，g/cm3；TVD0為測(cè)井曲線起始點(diǎn)垂深，m；TVD為給定計(jì)算測(cè)點(diǎn)所對(duì)應(yīng)的垂深，m。

通過(guò)對(duì)比孔隙壓力預(yù)測(cè)方法及預(yù)測(cè)效果，優(yōu)選Eaton法計(jì)算目標(biāo)工區(qū)地層孔隙壓力。Eaton法以泥巖正常壓實(shí)理論為基礎(chǔ)，因此，先要建立正常壓實(shí)趨勢(shì)線，然后才能計(jì)算地層孔隙壓力。對(duì)于目標(biāo)工區(qū)，其正常壓實(shí)趨勢(shì)線可表示為[27]：

而Eaton法計(jì)算地層孔隙壓力的公式可表示為[30]：

式中：pp為地層孔隙壓力，MPa；pw為地層流體靜液柱壓力，MPa；Δtnp為正常趨勢(shì)線上對(duì)應(yīng)的縱波時(shí)差，μs/ft；c為Eaton常數(shù)，無(wú)因次，目標(biāo)工區(qū)Eaton常數(shù)取值為0.914。

最后，計(jì)算最大和最小水平地應(yīng)力，水平地應(yīng)力的計(jì)算模式主要分為兩種[31-32]：一種是單軸應(yīng)變模式，認(rèn)為水平方向無(wú)限大，地層在沉積過(guò)程中只發(fā)生垂向變形，水平方向的應(yīng)力是由上覆巖層重量引起的，主要包括Terzaghi模型、Matthews-Kelly模型、Anderson模型和Newberry模型，這類模式?jīng)]有考慮構(gòu)造應(yīng)力，不符合實(shí)際地質(zhì)條件；另一類是分層計(jì)算模式，主要包括組合彈簧模型、斯倫貝謝模型和黃氏模型，該模式考慮了構(gòu)造應(yīng)力對(duì)水平地應(yīng)力的影響，因此被廣泛應(yīng)用于油田地應(yīng)力預(yù)測(cè)。由于目標(biāo)工區(qū)整體呈低緩構(gòu)造特征，因此，采用黃氏模型計(jì)算水平地應(yīng)力，黃氏模型可表示為[32]：

式中：σH為最大水平地應(yīng)力，MPa；σh為最小水平地應(yīng)力，MPa；ω1和ω2分別為最大和最小水平地應(yīng)力方向的構(gòu)造應(yīng)力系數(shù)，無(wú)因次。

在式(8)所示的黃氏模型中，構(gòu)造應(yīng)力系數(shù)的確定非常關(guān)鍵。對(duì)確定的工區(qū)，其構(gòu)造應(yīng)力系數(shù)可視為常數(shù)。利用室內(nèi)地應(yīng)力差應(yīng)變測(cè)試結(jié)果，并結(jié)合孔隙壓力(pp)、靜態(tài)泊松比(μs)、Biot系數(shù)(α)測(cè)井解釋結(jié)果，可以利用式(8)反演出工區(qū)的構(gòu)造應(yīng)力系數(shù)，目標(biāo)工區(qū)構(gòu)造應(yīng)力系數(shù)取值為ω1=0.458、ω2=0.364，進(jìn)而再利用式(8)計(jì)算連續(xù)的水平地應(yīng)力測(cè)井解釋剖面。

1.2 測(cè)井解釋結(jié)果

通過(guò)測(cè)井解釋得到了工區(qū)2口直井(X-1井和X-2井)的地應(yīng)力剖面，包括垂向地應(yīng)力、最大水平地應(yīng)力和最小水平地應(yīng)力，結(jié)果如圖2所示。由圖可知，每口井有10條測(cè)井曲線，包括自然伽馬(GR)、縱波時(shí)差(DTC)、井徑(CAL)、密度(DEN)、補(bǔ)償中子(CNL)、橫波時(shí)差(DTS)、測(cè)深(MD)、方位角(AZIM)、井斜角(DEVI)和垂深(TVD)。其中，X-1井測(cè)深為2809～3664.3 m，測(cè)段總長(zhǎng)為855.5 m，；X-2井測(cè)深為3218.1～3693.0 m，測(cè)段總長(zhǎng)為474.9 m。不難看出：該工區(qū)垂向應(yīng)力＞最大水平地應(yīng)力＞最小水平地應(yīng)力，屬于潛在正斷層應(yīng)力狀態(tài)，而且，地應(yīng)力隨井深的增加呈增大趨勢(shì)。此外，為了驗(yàn)證黃氏模型的準(zhǔn)確性，將X-1井地應(yīng)力解釋結(jié)果與3663.2 m、3657.3 m和3669.1 m處巖心差應(yīng)變實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比發(fā)現(xiàn)，X-1井垂向地應(yīng)力誤差為0.39%，最大水平地應(yīng)力誤差為0.18%～0.64%，最小水平地應(yīng)力誤差為0.29%。因此，測(cè)井解釋得到的地應(yīng)力與實(shí)測(cè)地應(yīng)力吻合較好，可以將測(cè)井解釋得到的地應(yīng)力作為真實(shí)地應(yīng)力。

圖2 地應(yīng)力測(cè)井解釋成果圖Fig. 2 Logging interpretation results of in-situ stress

2 數(shù)據(jù)預(yù)處理與樣本構(gòu)造

數(shù)據(jù)預(yù)處理和樣本構(gòu)造是機(jī)器學(xué)習(xí)建模的重要環(huán)節(jié)，前者決定了機(jī)器學(xué)習(xí)模型預(yù)測(cè)和泛化能力的好壞，而后者確定了數(shù)據(jù)輸入模型的形式。本文分別將X-1井和X-2井作為訓(xùn)練井和測(cè)試井，其中，訓(xùn)練井用于開發(fā)機(jī)器學(xué)習(xí)模型，測(cè)試井用于驗(yàn)證模型的預(yù)測(cè)和泛化能力。

2.1 數(shù)據(jù)處理

測(cè)井參數(shù)是對(duì)井下巖石物性特征的綜合響應(yīng)，不同測(cè)井參數(shù)能夠反映地層巖性、泥質(zhì)含量、含水飽和度、滲透率等特征[33]。而地應(yīng)力取決于地質(zhì)構(gòu)造、沉積環(huán)境、地層巖性等因素，因此，測(cè)井參數(shù)與地應(yīng)力之間存在著必然聯(lián)系，這也是利用測(cè)井?dāng)?shù)據(jù)預(yù)測(cè)地應(yīng)力的理論依據(jù)。同時(shí)，測(cè)井參數(shù)與地應(yīng)力之間關(guān)系復(fù)雜，不是所有測(cè)井參數(shù)都與地應(yīng)力之間存在關(guān)聯(lián)關(guān)系。因此，需要對(duì)測(cè)井參數(shù)進(jìn)行合理篩選，以消除冗余參數(shù)對(duì)地應(yīng)力預(yù)測(cè)的影響。經(jīng)典機(jī)器學(xué)習(xí)方法往往采用相關(guān)性分析來(lái)挖掘數(shù)據(jù)之間的關(guān)聯(lián)關(guān)系。

考慮到測(cè)井參數(shù)之間一般呈非正態(tài)分布，因此，采用Spearman秩相關(guān)系數(shù)來(lái)表征不同測(cè)井參數(shù)之間的相關(guān)性。Spearman秩相關(guān)系數(shù)是一種非參數(shù)的相關(guān)性度量準(zhǔn)則，通過(guò)對(duì)變量X和Y求秩后再進(jìn)行相關(guān)性分析，能夠滿足非正態(tài)分布數(shù)據(jù)的相關(guān)性計(jì)算，Spearman秩相關(guān)系數(shù)計(jì)算公式為[28]：

式中：RS為Spearman秩相關(guān)系數(shù)；n為樣本的數(shù)量；di為變量X和Y的秩次之差。

同時(shí)，根據(jù)式(9)計(jì)算得到的相關(guān)性系數(shù)，可以按照以下取值范圍判斷相關(guān)性強(qiáng)弱[34]：0.8≤|RS|≤1.0表 示 極 強(qiáng) 相 關(guān)；0.6≤|RS|≤0.79表 示 強(qiáng) 相 關(guān)；0.4≤|RS|≤0.59表示中等相關(guān)；0.2≤|RS|≤0.39表示弱相關(guān)；0.0≤|RS|≤0.19表示不相關(guān)。

圖3所示為X-1井測(cè)井參數(shù)相關(guān)系分析結(jié)果，不難看出：各測(cè)井參數(shù)與水平地應(yīng)力之間具有極其相似的相關(guān)性系數(shù)，但不同測(cè)井參數(shù)與水平地應(yīng)力之間的相關(guān)性存在較大差異：水平地應(yīng)力與垂深、密度及自然伽馬等測(cè)井參數(shù)之間表現(xiàn)為正相關(guān)，與縱波時(shí)差、井徑、補(bǔ)償中子及橫波時(shí)差等測(cè)井參數(shù)之間表現(xiàn)為負(fù)相關(guān)；最大和最小水平地應(yīng)力與垂深之間的相關(guān)性系數(shù)最大且均為0.98，其次是井徑、密度、補(bǔ)償中子和自然伽馬，分別為-0.38、0.37、-0.36和0.34，縱波時(shí)差和橫波時(shí)差相關(guān)性系數(shù)相同且最小為-0.27；根據(jù)相關(guān)性強(qiáng)度判別標(biāo)準(zhǔn)，水平地應(yīng)力與垂深之間相關(guān)性極強(qiáng)，與井徑、密度、補(bǔ)償中子、自然伽馬、縱波時(shí)差、橫波時(shí)差之間相關(guān)性相對(duì)較弱。整體來(lái)看，X-1井各測(cè)井參數(shù)與水平地應(yīng)力之間均存在一定程度的相關(guān)性。

圖3 X-1井測(cè)井參數(shù)相關(guān)性分析圖Fig. 3 Correlation analysis of logging parameters in X-1 well

原始測(cè)井?dāng)?shù)據(jù)受外部環(huán)境以及測(cè)井儀器自身性能的影響，使得原始測(cè)井?dāng)?shù)據(jù)具有較大噪聲[35]。考慮到卡爾曼濾波具有無(wú)偏、穩(wěn)定和概率意義上最優(yōu)的特點(diǎn)[36]，采用卡爾曼濾波算法對(duì)原始測(cè)井?dāng)?shù)據(jù)進(jìn)行去噪。常規(guī)卡爾曼濾波算法包含估計(jì)和修正兩個(gè)過(guò)程，估計(jì)過(guò)程主要依靠時(shí)間更新方程，即首先根據(jù)式(10)建立當(dāng)前狀態(tài)的先驗(yàn)估計(jì)值êk|k-1，然后根據(jù)式(11)推算當(dāng)前狀態(tài)的誤差協(xié)方差估計(jì)值Pk|k-1：

式中：êk|k-1為先驗(yàn)估計(jì)狀態(tài)向量；A為狀態(tài)轉(zhuǎn)移矩陣；Pk為后驗(yàn)估計(jì)誤差協(xié)方差矩陣；Pk|k-1為先驗(yàn)估計(jì)誤差協(xié)方差矩陣；Q為過(guò)程噪聲協(xié)方差矩陣。

然后，在獲得先驗(yàn)估計(jì)的前提下，利用觀測(cè)值對(duì)當(dāng)前狀態(tài)進(jìn)行后驗(yàn)估計(jì)：首先，根據(jù)式(12)計(jì)算當(dāng)前狀態(tài)的卡爾曼增益矩陣Kk；然后，結(jié)合觀測(cè)向量zk，根據(jù)式(13)對(duì)先驗(yàn)估計(jì)狀態(tài)向量êk|k-1進(jìn)行修正，得到狀態(tài)估計(jì)向量êk，即為最優(yōu)估計(jì)；最后，根據(jù)式(14)將誤差協(xié)方差Pk|k-1更新為Pk：

式中：êk為后驗(yàn)狀態(tài)估計(jì)向量；zk為觀測(cè)向量；R為測(cè)量噪聲協(xié)方差矩陣；S為測(cè)量協(xié)方差矩陣；Kk為卡爾曼增益。

采用卡爾曼濾波對(duì)X-1井和X-2井的井徑、密度、自然伽馬、縱波時(shí)差、橫波時(shí)差以及最大和最小水平地應(yīng)力等測(cè)井曲線進(jìn)行了去噪處理，去噪前后的對(duì)比曲線如圖4所示。不難看出：卡爾曼濾波可以有效消除測(cè)井?dāng)?shù)據(jù)中儀器自身噪聲和隨機(jī)性噪聲信息，有助于避免后續(xù)訓(xùn)練過(guò)程中可能會(huì)出現(xiàn)的噪聲擬合，從而避免噪聲引起的過(guò)擬合或欠擬合。因此，卡爾曼濾波去噪也是一種正則化手段。

圖4 測(cè)井參數(shù)去噪前后對(duì)比Fig. 4 Comparison of logging parameters before and after denoising

完成數(shù)據(jù)去噪后，需要對(duì)測(cè)井?dāng)?shù)據(jù)進(jìn)行歸一化處理。歸一化處理的目的是將不同類型數(shù)據(jù)映射到一定范圍內(nèi)，在不影響數(shù)據(jù)特性前提下，消除量綱對(duì)數(shù)據(jù)的影響，同時(shí)，加快機(jī)器學(xué)習(xí)模型的處理過(guò)程，從而提高訓(xùn)練速度。本文采用最大最小歸一化，將測(cè)井?dāng)?shù)據(jù)歸一化到0和1之間，最大最小歸一化表示為[37]：

式中：x*為歸一化之后的數(shù)據(jù)；x為原始數(shù)據(jù)；xmax和xmin分別為數(shù)據(jù)的最大值和最小值。

所有測(cè)井參數(shù)在輸入模型前都要進(jìn)行歸一化處理，值得注意的是最大和最小水平地應(yīng)力作為輸出數(shù)據(jù)，即目標(biāo)數(shù)據(jù)不需要進(jìn)行歸一化，但為了使模型能夠快速收斂，將最大和最小水平地應(yīng)力值縮小100倍。

2.2 樣本構(gòu)造

為有效提取測(cè)井?dāng)?shù)據(jù)序列的深度特征信息，采用了一種滑動(dòng)窗口的樣本構(gòu)造方式[38]。如圖5所示，將歸一化的測(cè)井?dāng)?shù)據(jù)組成二維數(shù)據(jù)矩陣形式，其中，特征表示相關(guān)測(cè)井參數(shù)，標(biāo)簽表示測(cè)井?dāng)?shù)據(jù)點(diǎn)所對(duì)應(yīng)的地應(yīng)力，n和l分別表示滑動(dòng)窗口的步長(zhǎng)和滑窗的大小，每一個(gè)樣本由輸入和輸出樣本組成。在樣本構(gòu)造過(guò)程中，將左側(cè)彩色實(shí)線框中的數(shù)據(jù)作為輸入樣本，右側(cè)彩色實(shí)線框中的數(shù)據(jù)作為輸出樣本?；懊看窝刂疃确较蚧瑒?dòng)一個(gè)步長(zhǎng)，便可構(gòu)造一個(gè)樣本，滑窗滑至最深處，即可完成所有樣本的構(gòu)造。考慮到本研究為一個(gè)深度序列預(yù)測(cè)任務(wù)，為了保證預(yù)測(cè)參數(shù)在深度上連續(xù)，取滑窗步長(zhǎng)為1。

圖5 樣本構(gòu)造方法Fig. 5 Sample construction method

3 雙向長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)模型

基于測(cè)井參數(shù)預(yù)測(cè)最大和最小水平地應(yīng)力，可描述為已知某一井段測(cè)井參數(shù)X={x1, x2,…, xn}以及對(duì)應(yīng)的地應(yīng)力參數(shù)Y={y1, y2,…, yn}，其中，n表示樣本數(shù)量，通過(guò)構(gòu)建模型，由訓(xùn)練得到測(cè)井參數(shù)X與地應(yīng)力參數(shù)Y之間的非線性映射關(guān)系。構(gòu)建模型的目的在于，當(dāng)給出新的且未進(jìn)行地應(yīng)力測(cè)井解釋的測(cè)井段時(shí)，模型能夠預(yù)測(cè)出最優(yōu)地應(yīng)力預(yù)測(cè)結(jié)果Y?={?1, ?2,…, ?n}。

3.1 模型理論

根據(jù)上述水平地應(yīng)力預(yù)測(cè)問(wèn)題的描述，本文使用一種特殊的循環(huán)神經(jīng)網(wǎng)絡(luò)，即雙向長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)(BiLSTM)[39]，BiLSTM網(wǎng)絡(luò)結(jié)構(gòu)如圖6所示。一個(gè)完整的BiLSTM網(wǎng)絡(luò)包括數(shù)據(jù)輸入層、前向LSTM、后向LSTM、激活函數(shù)層及數(shù)據(jù)輸出層。

圖6 BiLSTM網(wǎng)絡(luò)結(jié)構(gòu)圖Fig. 6 Network structure diagram of BiLSTM

BiLSTM主要由前向和后向長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)(LSTM)組成[40]，通過(guò)兩個(gè)平行的隱藏單元層來(lái)捕捉前后兩個(gè)方向上的序列關(guān)系，進(jìn)一步提高了特征提取的完整性和全局性，使得地應(yīng)力預(yù)測(cè)更加準(zhǔn)確。LSTM是BiLSTM的基本模塊，不同于普通的循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)，LSTM由于特殊的結(jié)構(gòu)設(shè)計(jì)使得其能夠有效避免時(shí)間序列模型在反向傳播過(guò)程中因?yàn)閷?duì)序列產(chǎn)生長(zhǎng)期依賴所導(dǎo)致的梯度消失和梯度爆炸問(wèn)題[41]。相較于傳統(tǒng)的RNN，LSTM結(jié)構(gòu)更為復(fù)雜，其包含了四個(gè)交互作用門層(遺忘門層、Tanh門層、輸入門層、輸出門層)[42-43]，從而形成了獨(dú)有的反饋循環(huán)結(jié)構(gòu)，LSTM的網(wǎng)絡(luò)結(jié)構(gòu)如圖7所示。

圖7 LSTM網(wǎng)絡(luò)結(jié)構(gòu)圖Fig. 7 Network structure diagram of LSTM

LSTM的核心思想是細(xì)胞態(tài)，其貫穿整個(gè)循環(huán)過(guò)程，在LSTM的處理過(guò)程中，通過(guò)精心設(shè)計(jì)的交互作用門層，可根據(jù)上一時(shí)刻的隱狀態(tài)和當(dāng)前時(shí)刻的輸入，選擇性刪除上一時(shí)刻細(xì)胞態(tài)中的信息或者添加新信息到當(dāng)前時(shí)刻的細(xì)胞態(tài)中。具體而言，第一個(gè)交互門層稱為遺忘門層，它決定了哪些信息該從細(xì)胞態(tài)中刪除，遺忘門層可表示為[39]：

式中：ft為遺忘門；Wf為遺忘門權(quán)重矩陣；ht-1為上一時(shí)刻的隱狀態(tài)；xt為當(dāng)前時(shí)刻的輸入；ht為當(dāng)前時(shí)刻的隱狀態(tài)；bf為遺忘門的偏置向量；σ為Sigmoid激活函數(shù)。

第二個(gè)交互作用門層稱為輸入門層，它決定了哪些信息將會(huì)被細(xì)胞態(tài)存貯，輸入門層可表示為式(17)。接下來(lái)便利用第三個(gè)交互門層，即Tanh門層生成一個(gè)新的候選值，用于描述當(dāng)前時(shí)刻的細(xì)胞態(tài)，Tanh門層可表示為式(18)[40]：

式中：it為輸入門；t為新的候選值；Wi為輸入門權(quán)重矩陣；Wc為Tanh門權(quán)重矩陣；bi為輸入門偏置向量；bc為Tanh門偏置向量。

當(dāng)經(jīng)過(guò)上述步驟后，舊的細(xì)胞態(tài)Ct-1將被更新成為新的細(xì)胞態(tài)Ct。首先，將遺忘門的值ft與舊的細(xì)胞態(tài)Ct-1進(jìn)行點(diǎn)乘以忘記部分信息，保留更為重要的信息；然后，將輸入門的值it與新的候選值t進(jìn)行點(diǎn)乘以決定哪些輸入信息將會(huì)被保留；最后，將兩者的計(jì)算結(jié)果相加完成對(duì)細(xì)胞態(tài)的更新；該過(guò)程可表示為[41]：

式中：Ct為當(dāng)前時(shí)刻的細(xì)胞態(tài)；Ct-1為上一時(shí)刻的細(xì)胞態(tài)。

第四個(gè)交互門層稱為輸出門層，決定細(xì)胞態(tài)的哪些信息將會(huì)被輸出。在這一階段，首先會(huì)將細(xì)胞態(tài)的值通過(guò)Tanh函數(shù)映射到-1和1之間，然后與輸出門層的計(jì)算結(jié)果進(jìn)行點(diǎn)乘，輸出當(dāng)前時(shí)刻的隱狀態(tài)，該過(guò)程可表示為：

式中：Ot為輸出門；Wo為輸出門權(quán)重矩陣；bo為輸出門偏置向量；ht為當(dāng)前時(shí)刻隱狀態(tài)。

對(duì)于BiLSTM模型，其隱狀態(tài)包含兩部分，一部分來(lái)自于前向LSTM，一部分來(lái)自于后向LSTM，分別用和表示。因此，BiLSTM的隱狀態(tài)可以被定義為[41]：

式中：yt為BiLSTM的輸出；為當(dāng)前時(shí)刻前向LSTM的隱狀態(tài)；為當(dāng)前時(shí)刻后向LSTM的隱狀態(tài)；c1和c2為常數(shù)，且c1+c2=1。

3.2 模型構(gòu)建

由BiLSTM構(gòu)建的預(yù)測(cè)模型具有良好的適應(yīng)性，在解決復(fù)雜地質(zhì)問(wèn)題方面表現(xiàn)出了卓越的預(yù)測(cè)能力[40]。本文所構(gòu)建的地應(yīng)力預(yù)測(cè)模型包括2個(gè)BiLSTM網(wǎng)絡(luò)層、1個(gè)丟棄層(Dropout)、1個(gè)拉直層(Flatten)和1個(gè)全連接層(Linear)，并分別在BiLSTM網(wǎng)絡(luò)層和Linear層后面加上雙曲正切(Tanh)激活函數(shù)層和Sigmoid激活函數(shù)層，模型結(jié)構(gòu)和訓(xùn)練流程如圖8所示。其中，BiLSTM網(wǎng)絡(luò)層是模型的核心層，用于捕捉測(cè)井參數(shù)的序列特征以及不同測(cè)井參數(shù)之間的內(nèi)在聯(lián)系；Dropout層用于避免模型過(guò)擬合；Flatten層用于將高維數(shù)據(jù)特征變換成1維，實(shí)現(xiàn)BiLSTM網(wǎng)絡(luò)層到全連接層的過(guò)渡；Linear層則將所學(xué)習(xí)的特征映射到樣本標(biāo)記空間形成預(yù)測(cè)值；Tanh和Sigmoid激活函數(shù)則完成數(shù)據(jù)的非線性變換。需要注意的是，在模型構(gòu)建初期，暫時(shí)按照經(jīng)驗(yàn)將樣本構(gòu)造中的窗口大小設(shè)置為30，數(shù)據(jù)批尺寸設(shè)置為16，學(xué)習(xí)率設(shè)置為0.01，迭代次數(shù)為50次；模型的訓(xùn)練流程可分為以下幾個(gè)階段：

圖8 模型結(jié)構(gòu)和訓(xùn)練流程圖Fig. 8 Model structure and training flowchart

1) 將測(cè)井?dāng)?shù)據(jù)進(jìn)行去噪和歸一化處理后按照2.2節(jié)所述樣本構(gòu)造方式構(gòu)造出特征和標(biāo)簽樣本，并將樣本進(jìn)行打亂，然后分批次喂入模型；

2) 樣本數(shù)據(jù)進(jìn)入模型后首先經(jīng)過(guò)BiLSTM網(wǎng)絡(luò)，隨著數(shù)據(jù)在網(wǎng)絡(luò)中的自由流動(dòng)，預(yù)測(cè)地應(yīng)力的關(guān)鍵特征被提取出來(lái)，在這一階段BiLSTM網(wǎng)絡(luò)的隱狀態(tài)設(shè)置為12維，并使用Tanh函數(shù)作為非線性激活函數(shù)；

3) 從BiLSTM網(wǎng)絡(luò)層提取的特征被送入Dropout層，Dropout層按照一定的概率對(duì)部分隱藏神經(jīng)元進(jìn)行丟棄操作，丟棄率設(shè)置為20%；

4) 經(jīng)過(guò)Dropout層后的特征進(jìn)入Flatten層，這時(shí)多維特征變被拉直成1維，使得特征能過(guò)渡到Linear層，最后特征經(jīng)過(guò)Linear層，使用Sigmoid函數(shù)進(jìn)行激活后映射成為預(yù)測(cè)值；

5) 計(jì)算預(yù)測(cè)值與標(biāo)簽之間的誤差，并選擇自適應(yīng)矩估計(jì)(Adam)算法來(lái)調(diào)整模型的權(quán)值參數(shù)，使得訓(xùn)練過(guò)程中的誤差可以不斷減小直到達(dá)到迭代終止條件，達(dá)到迭代終止條件后保存模型直接用于后續(xù)測(cè)試集的預(yù)測(cè)。

4 地應(yīng)力預(yù)測(cè)結(jié)果及討論

4.1 結(jié)果分析

將X-1井常規(guī)測(cè)井?dāng)?shù)據(jù)和測(cè)井解釋得到的水平地應(yīng)力數(shù)據(jù)作為訓(xùn)練集，在訓(xùn)練過(guò)程中隨機(jī)選取訓(xùn)練集30%的數(shù)據(jù)作為驗(yàn)證集，將X-2井常規(guī)測(cè)井參數(shù)和測(cè)井解釋得到的水平地應(yīng)力數(shù)據(jù)作為測(cè)試集，以驗(yàn)證該方法的泛化能力。其中，訓(xùn)練集和驗(yàn)證集用于建立和評(píng)估模型，測(cè)試集用于測(cè)試BiLSTM模型預(yù)測(cè)精度。為了分析不同測(cè)井參數(shù)組合對(duì)水平地應(yīng)力預(yù)測(cè)的影響，根據(jù)相關(guān)性分析結(jié)果，按照相關(guān)性系數(shù)的大小順序，設(shè)計(jì)了7種不同的測(cè)井參數(shù)組合方式：組合1僅使用1種測(cè)井參數(shù)：TVD；組合2使用2種測(cè)井參數(shù)：TVD和CAL；組合3使用3種測(cè)井參數(shù)：TVD、CAL和DEN；組合4使用4種測(cè)井參數(shù)：TVD、CAL、DEN和CNL；組合5使用5種測(cè)井參數(shù)：TVD、CAL、DEN、CNL和GR；組合6使用6種測(cè)井參數(shù)：TVD、CAL、DEN、CNL、GR和DTC；組合7使用7種測(cè)井參數(shù)：TVD、CAL、DEN、CNL、GR、DTC和DTS。

為了評(píng)價(jià)模型預(yù)測(cè)效果，采用均方根誤差(RMSE)、平均絕對(duì)誤差(MAE)、平均絕對(duì)百分比誤差(MAPE)作為評(píng)價(jià)指標(biāo)，其計(jì)算公式如下[33,38]：

式中：N為樣本數(shù)量；yi為真實(shí)地應(yīng)力值，MPa；?i為預(yù)測(cè)地應(yīng)力值，MPa；RMSE為均方根誤差，MPa；MAE為平均絕對(duì)誤差，MPa；MAPE為平均絕對(duì)百分比誤差，%。

利用X-1井訓(xùn)練機(jī)器學(xué)習(xí)模型，然后按照7種測(cè)井參數(shù)組合方式，分別對(duì)X-2井3218.1～3693.0 m井段地層的最大和最小水平地應(yīng)力進(jìn)行預(yù)測(cè)，同時(shí)進(jìn)行了超參數(shù)優(yōu)化，優(yōu)化前后預(yù)測(cè)結(jié)果如圖9所示。圖9a為實(shí)際和預(yù)測(cè)最大水平地應(yīng)力對(duì)比，圖9b為實(shí)際和預(yù)測(cè)最小水平地應(yīng)力對(duì)比，其中，紅色實(shí)線代表實(shí)際水平地應(yīng)力，綠色實(shí)線代表預(yù)測(cè)得到的水平地應(yīng)力，藍(lán)色實(shí)線代表超參數(shù)優(yōu)化后得到的水平地應(yīng)力。

不同測(cè)井參數(shù)組合下預(yù)測(cè)結(jié)果的評(píng)價(jià)指標(biāo)如表1所示，結(jié)合圖9分析不難發(fā)現(xiàn)：1) X-2井預(yù)測(cè)井段實(shí)際最大水平地應(yīng)力介于78.5～98.3 MPa，平均值為84.4 MPa；最小水平地應(yīng)力介于74.1～93.7 MPa，平均值為79.7 MPa。2) 7種不同測(cè)井參數(shù)組合模式下，組合6的預(yù)測(cè)結(jié)果與實(shí)際水平地應(yīng)力最為接近，最大和最小水平地應(yīng)力的平均絕對(duì)百分比誤差(MAPE)分別為1.39‰和1.76‰；組合2的預(yù)測(cè)效果最差，最大和最小水平地應(yīng)力的平均絕對(duì)百分比誤差(MAPE)分別為4.26%和4.33%。3) 組合4～7的預(yù)測(cè)效果明顯優(yōu)于組合1～3；從組合3到組合4，在增加了補(bǔ)償中子作為輸入?yún)?shù)之后，預(yù)測(cè)誤差大幅下降，最大和最小水平地應(yīng)力的平均絕對(duì)百分比誤差(MAPE)分別由38.48‰和40.93‰下降至10.37‰和11.35‰，說(shuō)明在增加了補(bǔ)償中子后，模型學(xué)習(xí)到了預(yù)測(cè)地應(yīng)力的關(guān)鍵信息，分配給了補(bǔ)償中子更大的權(quán)重；另一方面，補(bǔ)償中子用于指示巖石的孔隙度，對(duì)地層孔隙壓力極為敏感，說(shuō)明水平地應(yīng)力與孔隙壓力關(guān)系較為密切。

圖9 不同參數(shù)組合下超參數(shù)優(yōu)化前后X-2井最大和最小水平地應(yīng)力預(yù)測(cè)結(jié)果對(duì)比Fig. 9 Comparison of prediction results of maximum and minimum horizontal in-situ stress in X-2 well before and after hyperparameter optimization under different parameter combinations

表1 不同參數(shù)組合下X-2井水平地應(yīng)力預(yù)測(cè)效果評(píng)價(jià)指標(biāo)Table 1 Evaluation indices of horizontal in-situ stresses prediction for X-2 well under different parameter combinations

由于機(jī)器學(xué)習(xí)模型超參數(shù)的選擇與估計(jì)對(duì)模型的預(yù)測(cè)與泛化性能具有重要影響[45]，很多研究學(xué)者在模型超參數(shù)優(yōu)化過(guò)程中將模型學(xué)習(xí)率、批尺寸、迭代次數(shù)等超參數(shù)作為主要影響因素，利用網(wǎng)格搜索、隨機(jī)搜索和貝葉斯優(yōu)化等方法對(duì)模型進(jìn)行優(yōu)化[46]。為此，將模型學(xué)習(xí)率、批尺寸、迭代次數(shù)三種超參數(shù)作為主要影響因素來(lái)研究其對(duì)模型性能的影響，此外，考慮到樣本構(gòu)造方式對(duì)的特殊性，將滑窗的窗口大小也作為主要影響因素進(jìn)行研究。針對(duì)多因素的交互影響，普通實(shí)驗(yàn)設(shè)計(jì)方案需要大量重復(fù)實(shí)驗(yàn)且耗時(shí)長(zhǎng)，為此，采用正交設(shè)計(jì)實(shí)驗(yàn)方法，該方法適合于多因素、多水平的研究，能夠從眾多實(shí)驗(yàn)數(shù)據(jù)中選取適當(dāng)具有代表性的數(shù)據(jù)進(jìn)行測(cè)試，這些實(shí)驗(yàn)數(shù)據(jù)具有分布均勻、齊整可比的特點(diǎn)[44]。此處以學(xué)習(xí)率(A)、批尺寸(B)、迭代次數(shù)(C)以及窗口大小(D)等四個(gè)超參數(shù)為因子，選取各因子的典型值作為水平值，每個(gè)因子選取四個(gè)水平值，其中，學(xué)習(xí)率分為0.005、0.01、0.015、0.02四個(gè)水平值，批尺寸分為16、32、64、128四個(gè)水平值，迭代次數(shù)分為50、100、150、200四個(gè)水平值，窗口大小分為10、20、30、40四個(gè)水平值，從而形成了四因子四水平的16組正交實(shí)驗(yàn)方案，如表2所示。此外，正交實(shí)驗(yàn)采用的輸入特征參數(shù)組合為組合6，計(jì)算了不同實(shí)驗(yàn)方案預(yù)測(cè)得到的最大和最小水平地應(yīng)力的評(píng)價(jià)指標(biāo)RMSE、MAE和MAPE，結(jié)果如表2和圖10所示。

表2 超參數(shù)正交設(shè)計(jì)實(shí)驗(yàn)表及其對(duì)應(yīng)的誤差結(jié)果Table 2 Orthogonal design experiment table of hyperparameters and the corresponding errors

圖10 不同實(shí)驗(yàn)方案預(yù)測(cè)結(jié)果的誤差對(duì)比Fig. 10 Comparison of prediction errors in different experimental schemes

結(jié)合表2和圖10不難看出：1) 不同超參數(shù)條件預(yù)測(cè)的最大和最小水平地應(yīng)力誤差總體相似，且最小水平地應(yīng)力的預(yù)測(cè)誤差略大于最大水平地應(yīng)力；2) 除了第4組實(shí)驗(yàn)外，其余實(shí)驗(yàn)方案均能夠?qū)⑺降貞?yīng)力預(yù)測(cè)誤差控制在較小的范圍內(nèi)，其中，最大水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE介于0.06～0.16 MPa、MAE介于0.04～0.13 MPa、MAPE介于0.48～1.49‰，最小水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE介于0.06～0.20 MPa、MAE介于0.04～0.16 MPa、MAPE介于0.50～1.86‰，說(shuō)明超參數(shù)對(duì)水平地應(yīng)力的影響相對(duì)較?。?) 第13組實(shí)驗(yàn)預(yù)測(cè)得到的水平地應(yīng)力誤差最小，其中，最大水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE、MAE、MAPE分別為0.06 MPa、0.04 MPa和0.48‰，最小水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE、MAE、MAPE分 別 為0.06 MPa、0.04 MPa和0.50‰，因此，推薦采用學(xué)習(xí)率、批尺寸、窗口大小和迭代次數(shù)分別為0.02、16、20和150。

根據(jù)正交設(shè)計(jì)實(shí)驗(yàn)計(jì)算結(jié)果，分析了不同超參數(shù)對(duì)水平地應(yīng)力預(yù)測(cè)誤差的影響規(guī)律，結(jié)果如圖10子圖所示，不難看出：1) 隨著學(xué)習(xí)率的增加，MAPE逐漸降低，當(dāng)學(xué)習(xí)率為0.02時(shí)候，MAPE小幅增加；2) 隨著批尺寸的增加，MAPE逐漸增加，與之相反，隨著迭代次數(shù)的增加，MAPE逐漸下降；3) 窗口大小對(duì)MAPE的影響無(wú)明確規(guī)律，在模擬參數(shù)條件下窗口大小10或30的MAPE更低；4) 不同超參數(shù)影響下最大和最小水平地應(yīng)力預(yù)測(cè)結(jié)果的MAPE變化趨勢(shì)基本一致，但最小水平地應(yīng)力預(yù)測(cè)結(jié)果的MAPE略高于最大水平地應(yīng)力。

此外，按照方案13的超參數(shù)組合，對(duì)表1所述的7種輸入特征參數(shù)組合的水平地應(yīng)力進(jìn)行了預(yù)測(cè)，預(yù)測(cè)結(jié)果如圖9所示。圖9a為參數(shù)優(yōu)化前后最大水平地應(yīng)力預(yù)測(cè)結(jié)果與測(cè)井解釋結(jié)果的對(duì)比，圖9b為參數(shù)優(yōu)化前后最小水平地應(yīng)力預(yù)測(cè)結(jié)果與測(cè)井解釋結(jié)果的對(duì)比，其中，紅色實(shí)線代表實(shí)際水平地應(yīng)力，綠色實(shí)線代表預(yù)測(cè)得到的水平地應(yīng)力，藍(lán)色實(shí)線代表超參數(shù)優(yōu)化后得到的水平地應(yīng)力，通過(guò)分析不難發(fā)現(xiàn)：1) 按照第13組超參數(shù)對(duì)模型進(jìn)行優(yōu)化后，不同輸入特征組合方式預(yù)測(cè)得到的最大和最小水平地應(yīng)力的精度均得到了不同程度的提升，尤其是輸入特征參數(shù)組合1～4提升最為明顯，說(shuō)明第13組超參數(shù)優(yōu)化效果好且非常穩(wěn)定；2) 超參數(shù)優(yōu)化后，輸入?yún)?shù)組合6預(yù)測(cè)的水平地應(yīng)力誤差最小，即輸入?yún)?shù)組合6的預(yù)測(cè)結(jié)果仍然是最優(yōu)的，而且，優(yōu)化前后最大水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE、MAE、MAPE分別降低了0.13 MPa、0.08 MPa和0.91‰，最小水平地應(yīng)力預(yù)測(cè)結(jié)果的RMSE、MAE、MAPE分別降低了0.23 MPa、0.11 MPa和1.26‰，說(shuō)明最小水平地應(yīng)力的預(yù)測(cè)精度提升更加明顯；3) 不同輸入?yún)?shù)組合方式下，在增加或減少某一測(cè)井參數(shù)后的預(yù)測(cè)效果與超參數(shù)優(yōu)化前的結(jié)論基本一致，也進(jìn)一步說(shuō)明輸入?yún)?shù)組合6是最優(yōu)的。

4.2 討論

首先，分析了7種不同測(cè)井參數(shù)組合對(duì)預(yù)測(cè)水平地應(yīng)力的影響。對(duì)于組合1，只加入了垂深作為輸入?yún)?shù)，其預(yù)測(cè)結(jié)果表現(xiàn)為1條斜線，這是可以預(yù)見的，因?yàn)榇股畋旧沓示€性，所以模型也只學(xué)到了這種線性關(guān)系；對(duì)于組合2，在增加了井徑作為輸入?yún)?shù)之后預(yù)測(cè)結(jié)果反而變差，常規(guī)測(cè)井解釋中井徑用于指示地層巖性和判別地應(yīng)力方向，但對(duì)于最大和最小水平地應(yīng)力的大小而言并無(wú)直接關(guān)系，因此，井徑的增加給予了模型錯(cuò)誤的權(quán)重配比，反而使預(yù)測(cè)效果變差；對(duì)于組合3，在增加了密度測(cè)井曲線作為輸入?yún)?shù)之后，模型預(yù)測(cè)效果好轉(zhuǎn)，因?yàn)槊芏葴y(cè)井曲線與最大和最小水平地應(yīng)力的計(jì)算有直接關(guān)系，所以增加密度測(cè)井曲線有利于模型進(jìn)行預(yù)測(cè)；對(duì)于組合5，在增加了自然伽馬作為輸出參數(shù)之后，預(yù)測(cè)誤差進(jìn)一步減小，說(shuō)明自然伽馬對(duì)水平地應(yīng)力的預(yù)測(cè)也起到了一定作用，通常自然伽馬用于識(shí)別地層巖性和確定泥質(zhì)含量，顯然也構(gòu)成了計(jì)算地應(yīng)力的積極要素；對(duì)于組合6，在增加了縱波時(shí)差之后，預(yù)測(cè)效果有了進(jìn)一步提升，事實(shí)上縱波時(shí)差與地層構(gòu)造應(yīng)力之間有直接關(guān)系，通常當(dāng)?shù)貙邮艿筋~外的構(gòu)造應(yīng)力時(shí)，縱波時(shí)差會(huì)偏離正常趨勢(shì)線，因此，有學(xué)者利用這種關(guān)系建立起了預(yù)測(cè)水平地應(yīng)力的經(jīng)驗(yàn)關(guān)系式[47]；對(duì)于組合7，當(dāng)增加橫波時(shí)差之后誤差反而增加，可能是因?yàn)闄M波時(shí)差主要反映地層的巖性、孔隙度特征，與補(bǔ)償中子、自然伽馬、縱波時(shí)差等測(cè)井參數(shù)互為同質(zhì)參數(shù)，說(shuō)明過(guò)多的同質(zhì)參數(shù)會(huì)使數(shù)據(jù)維度增加，一定程度上降低模型的學(xué)習(xí)效果，從而增加預(yù)測(cè)誤差。綜上所述，結(jié)合測(cè)井參數(shù)實(shí)際含義，當(dāng)增加或減少某些測(cè)井參數(shù)后，由于不同參數(shù)對(duì)地應(yīng)力的指示作用不同，使模型在學(xué)習(xí)過(guò)程中發(fā)生權(quán)重偏移，進(jìn)而對(duì)不同參數(shù)組合模式下的預(yù)測(cè)結(jié)果表現(xiàn)出差異性。因此，綜合考慮測(cè)井參數(shù)的實(shí)際含義和參數(shù)獲取的難易程度，推薦使用垂深、井徑、密度、補(bǔ)償中子、自然伽馬和縱波時(shí)差6種測(cè)井參數(shù)作為水平地應(yīng)力的輸入?yún)?shù)組合。

其次，采用正交設(shè)計(jì)實(shí)驗(yàn)方法設(shè)計(jì)了16組實(shí)驗(yàn)，研究了不同超參數(shù)對(duì)BiLSTM模型預(yù)測(cè)性能的影響。結(jié)果發(fā)現(xiàn)在輸入?yún)?shù)組合6條件下16組超參數(shù)組合的精度都比較高，一方面是因?yàn)檩斎雲(yún)?shù)組合6為最優(yōu)輸入組合，可以確保模型訓(xùn)練效果和預(yù)測(cè)精度；另一方面，超參數(shù)取值范圍通過(guò)經(jīng)驗(yàn)設(shè)置在常用范圍內(nèi)，避免了因?yàn)樘荻缺ɑ蛘咛荻认Ф鸬木冉档?。將輸入?yún)?shù)組合6條件下優(yōu)化獲得得超參數(shù)用于其它6種輸入?yún)?shù)組合，結(jié)果發(fā)現(xiàn)不同輸入特征組合方式預(yù)測(cè)得到的最大和最小水平地應(yīng)力的精度均得到了不同程度的提升，而且輸入?yún)?shù)組合6的預(yù)測(cè)結(jié)果仍然是最優(yōu)的。從預(yù)測(cè)結(jié)果來(lái)看，適當(dāng)減小批尺寸和增加迭代次數(shù)可以增加模型預(yù)測(cè)精度，但減小批尺寸將導(dǎo)致模型學(xué)習(xí)效率的下降，增加迭代次數(shù)不僅增加訓(xùn)練成本而且增大了過(guò)擬合的可能性，因此，需要針對(duì)不同數(shù)據(jù)集，按照實(shí)際情況進(jìn)行合理調(diào)參。對(duì)于本研究推薦學(xué)習(xí)率、批尺寸、窗口大小和迭代次數(shù)等超參數(shù)分別取值為0.02、16、20和150。

訓(xùn)練數(shù)據(jù)集數(shù)量限制，數(shù)據(jù)集樣本獲取途徑單一，都是造成預(yù)測(cè)精度低的重要原因。考慮到工程中實(shí)際條件的限制，通常無(wú)法獲取足夠的測(cè)井?dāng)?shù)據(jù)集進(jìn)行訓(xùn)練，因此，采用滑窗采樣的樣本構(gòu)造方式和特征組合篩選實(shí)驗(yàn)來(lái)優(yōu)化輸入特征，并引入BiLSTM這一特殊循環(huán)神經(jīng)網(wǎng)絡(luò)來(lái)優(yōu)化模型結(jié)構(gòu)，另外設(shè)計(jì)了正交實(shí)驗(yàn)對(duì)超參數(shù)取值進(jìn)行優(yōu)化，以提高模型的預(yù)測(cè)精度和泛化能力。但是，本研究尚存在一定的局限性：首先，本研究?jī)H對(duì)CL氣田兩口直井的測(cè)井?dāng)?shù)據(jù)進(jìn)行測(cè)井解釋、訓(xùn)練和驗(yàn)證，所建立的模型數(shù)據(jù)來(lái)源較為單一，因此有必要加入地震、巖心實(shí)驗(yàn)、礦場(chǎng)實(shí)測(cè)等多源、多態(tài)、多尺度的數(shù)據(jù)擴(kuò)充訓(xùn)練樣本集，從而構(gòu)建出普適性和穩(wěn)健性更好的地應(yīng)力預(yù)測(cè)模型。其次，神經(jīng)網(wǎng)絡(luò)模型是基于經(jīng)驗(yàn)風(fēng)險(xiǎn)最小化原則，在小樣本預(yù)測(cè)問(wèn)題中存在過(guò)度擬合的風(fēng)險(xiǎn)，因此在后續(xù)建模過(guò)程中可以嘗試以結(jié)構(gòu)風(fēng)險(xiǎn)為原則或以結(jié)構(gòu)和經(jīng)驗(yàn)風(fēng)險(xiǎn)相結(jié)合為原則的機(jī)器學(xué)習(xí)模型，強(qiáng)化模型在理論方面的引導(dǎo)，避免因純數(shù)據(jù)驅(qū)動(dòng)造成的模型可解釋性和泛化能力降低的問(wèn)題。

5 結(jié)論

1) X-1和X-2井地應(yīng)力測(cè)井解釋結(jié)果表明，地應(yīng)力排序?yàn)榇瓜虻貞?yīng)力＞最大水平地應(yīng)力＞最小水平地應(yīng)力，屬于潛在正斷層應(yīng)力狀態(tài)；對(duì)比測(cè)井解釋和巖心差應(yīng)變測(cè)試結(jié)果發(fā)現(xiàn)，垂向地應(yīng)力解釋誤差為0.39%，最大水平地應(yīng)力解釋誤差為0.18%～0.64%，最小水平地應(yīng)力解釋誤差為0.29%，說(shuō)明測(cè)井解釋得到的地應(yīng)力與實(shí)際地應(yīng)力吻合較好。

2) 水平地應(yīng)力與垂深、密度及自然伽馬等測(cè)井參數(shù)之間表現(xiàn)為正相關(guān)，與縱波時(shí)差、井徑、補(bǔ)償中子及橫波時(shí)差等測(cè)井參數(shù)之間表現(xiàn)為負(fù)相關(guān)；水平地應(yīng)力與垂深之間相關(guān)性極強(qiáng)，與井徑、密度、補(bǔ)償中子、自然伽馬、縱波時(shí)差、橫波時(shí)差之間相關(guān)性相對(duì)較弱。

3) 不同測(cè)井參數(shù)蘊(yùn)含了不同的地質(zhì)信息，不同測(cè)井參數(shù)對(duì)地應(yīng)力預(yù)測(cè)起到了不同的指示作用。通過(guò)對(duì)比7種不同測(cè)井參數(shù)組合模式發(fā)現(xiàn)，組合4～7的預(yù)測(cè)效果明顯優(yōu)于組合1～3，最優(yōu)輸入?yún)?shù)組合為組合6，即包括垂深、井徑、密度、補(bǔ)償中子、自然伽馬、縱波時(shí)差；進(jìn)一步，通過(guò)正交設(shè)計(jì)實(shí)驗(yàn)對(duì)超參數(shù)取值進(jìn)行調(diào)優(yōu)，最佳的超參數(shù)取值為學(xué)習(xí)率0.02、批尺寸16、窗口大小20，迭代次數(shù)150，預(yù)測(cè)得到的最大和最小水平地應(yīng)力平均絕對(duì)百分比誤差分別為0.48‰和0.50‰。

4) BiLSTM模型是一種可以實(shí)現(xiàn)水平地應(yīng)力精準(zhǔn)預(yù)測(cè)的機(jī)器學(xué)習(xí)方法，得益于BiLSTM對(duì)時(shí)間序列出色的處理能力，使得其能夠有效捕捉測(cè)井參數(shù)隨深度的變化趨勢(shì)以及測(cè)井參數(shù)的前后關(guān)聯(lián)信息。因此，推薦使用BiLSTM模型進(jìn)行水平地應(yīng)力的預(yù)測(cè)。