基于支持向量機由海洋地震勘探數(shù)據(jù)對海洋水體溫度的預(yù)測?

彭陽陽， 張 進，2，3??

(1. 中國海洋大學(xué)海洋地球科學(xué)學(xué)院， 山東 青島 266100； 2. 海洋礦產(chǎn)資源評價與探測技術(shù)功能實驗室， 山東 青島 266237；3. 中國海洋大學(xué)海底科學(xué)與探測技術(shù)教育部重點實驗室， 山東 青島 266100)

海水的溫度結(jié)構(gòu)的精細劃分有利于研究水體的溫度結(jié)構(gòu)和海水運動等問題。常規(guī)的海水溫度測量多采用CTD,XBT等方式，雖然這些方式所測得的海洋水體垂向溫度分辨率高，但是橫向連續(xù)性較差。地震海洋學(xué)通過反射地震勘探方法來研究物理海洋現(xiàn)象，相對于傳統(tǒng)的接觸式溫鹽深測量方法，該方法具有較高的橫向分辨率和能快速對整個海水剖面成像等優(yōu)點[1]。

Holbrook等[2-3]首次用地震學(xué)方法獲得海水的地震疊加剖面，并揭示了海洋鋒處的溫鹽結(jié)構(gòu)。隨后Nandi等將常規(guī)反射地震與XCTD數(shù)據(jù)有效結(jié)合，驗證出溫度剖面和地震剖面有很好的一致性，并成功刻畫水團邊界[4]，更加驗證了將地震學(xué)應(yīng)用到海洋溫度反演是可行的。Páramo等[5]利用AVO技術(shù)分析海洋地震數(shù)據(jù)，反演出海水溫度參數(shù)，結(jié)果表明利用AVO技術(shù)雖然能得到海水聲速和溫度的差異，但對地震資料要求較高，并且反演得出海水密度參數(shù)不準(zhǔn)確。為了提高反演準(zhǔn)確度，Wood等[6]利用一維全波形反演方法盡可能的利用豐富的波場信息求取海水溫度參數(shù)，但Kormann等人分析出初始模型對一維全波形反演海水溫度參數(shù)結(jié)果有影響[7]并且運算速度低。為提高溫度剖面反演速度同時考慮到海水的溫度參數(shù)與地震數(shù)據(jù)的非線性關(guān)系，Padhi等人分別利用遺傳算法與非線性最小二乘反演方法結(jié)合反演海水的溫度，并且應(yīng)用效果顯著[8]。然而這些方法在反演海洋水體溫度參數(shù)方面，需要大量的溫度標(biāo)簽數(shù)據(jù)，限制了地震海洋學(xué)的發(fā)展。

支持向量機(Support Vector Machine，SVM)可以在小樣本數(shù)據(jù)情況下得到較好的反演結(jié)果[9]。通過文獻的調(diào)研，前人利用SVM在儲層預(yù)測和海浪預(yù)測等方面行之有效，例如：金權(quán)等人利用SVM對海洋波高進行預(yù)測，證明支持向量機在海浪預(yù)測的可行性[10]。高偉等人提出一種基于支持向量機對海底聲學(xué)參數(shù)快速統(tǒng)計反演的方法[11]。周凡等[12]利用支持向量機對流體進行識別，其預(yù)測準(zhǔn)確率達90.48%。唐小彪提取了地震屬性并且基于支持向量機對地震儲層進行預(yù)測，預(yù)測結(jié)果與設(shè)計的儲層地質(zhì)模型十分吻合[13]。

支持向量機可以在小樣本情況下能夠充分的利用樣本信息，該算法不僅易實現(xiàn)還具有較好的“魯棒”性[14]。本文將一種基于支持向量回歸機的方法應(yīng)用于海水溫度反演研究中。首先，提取海洋地震數(shù)據(jù)的屬性參數(shù)，將主成分分析方法優(yōu)化后的地震屬性作為輸入，再利用網(wǎng)格交叉驗證的手段訓(xùn)練支持向量回歸機，結(jié)合地震數(shù)據(jù)和CTD數(shù)據(jù)來預(yù)測海水溫度。結(jié)果表明支持向量回歸機在少量海水溫度樣本的情況下，對海水溫度預(yù)測的效果顯著，能夠精準(zhǔn)的劃分水層，擬合優(yōu)度達到98%，為后續(xù)的海洋研究奠定基礎(chǔ)。

1 支持向量回歸機(ε-SVR)模型

(1)

式中：b為偏置量；Φ為映射核函數(shù)；w為權(quán)向量。通過求解以下原始的優(yōu)化問題(如(2)式所示)獲得最佳的逼近函數(shù)：

(2)

對原始的優(yōu)化問題求解復(fù)雜，將求解原始問題通過拉格朗日乘子法和強對偶關(guān)系轉(zhuǎn)換為求其對偶問題，如(3)所示：

(3)

(4)

式中：Xj為核函數(shù)中心；σ為函數(shù)的寬度參數(shù)。

(5)

2 基于SVR的海水溫度實例預(yù)測

本文數(shù)據(jù)來自渤海海區(qū)，渤海海區(qū)是深入中國大陸的近封閉性淺海，其平均水深只有18 m，最深處也只有83 m[15]。本文截取1 000道渤海淺剖數(shù)據(jù)，地震采樣間隔為0.5 ms，截取每道86個樣點，記錄時長43.0 ms,深度32 m。8月份在渤海測量CTD為共3道(此時海水溫躍層最強)，截取深度32 m，如圖1所示，研究海區(qū)的海水溫度隨深度總體上呈現(xiàn)下降趨勢。根據(jù)CTD測量的海水深度和速度轉(zhuǎn)換到地震數(shù)據(jù)時間刻度如圖2所示。

(從左向右依次：CTD1,CTD2和CTD3。From left to right: CTD1, CTD2 and CTD3.)

(位于85、471和854地震道的三條黑線表示編號為CTD1、CTD2和CTD3的3個CTD。CTD1, CTD2 and CTD3 are located at CDP85, 741 and 854, respectively.)

利用主成分分析的方法優(yōu)選出對溫度敏感的地震屬性，將CTD標(biāo)簽數(shù)據(jù)和優(yōu)選后的地震屬性作為輸入向量，利用網(wǎng)格交叉驗證方式訓(xùn)練SVR，預(yù)測海水溫度并分析效果。海洋溫度參數(shù)與地震屬性存在非線性關(guān)系，通過ε-SVR構(gòu)造非線性逼近函數(shù)，具體流程包括以下步驟：

第一步：分析并處理原始地震數(shù)據(jù)與CTD數(shù)據(jù)，建立二者時空對應(yīng)關(guān)系。

第二步：從地震數(shù)據(jù)中提取并優(yōu)化地震屬性。

輸入：數(shù)據(jù)集P={p1,p2,……,pn}，需要降到k維。

(1)去平均值(即去中心化)，即每一維特征減去各自的平均值。

(2)計算協(xié)方差矩陣PPT，其中n為樣本數(shù)目。

(3)用特征值分解方法求協(xié)方差矩陣PPT的特征值。

(4)觀察特征值，選擇其中占比較重的特征值對應(yīng)的地震屬性。然后將其對應(yīng)地震屬性向量組成支持向量機輸入集x。

第三步：歸一化處理。對地震屬性和CTD樣本溫度標(biāo)簽進行歸一化處理，其轉(zhuǎn)換函數(shù)如式(6)所示：

(6)

第四步：網(wǎng)格交叉驗證尋找SVR中C和g最優(yōu)參數(shù)組合。選取高斯核函數(shù)進行模型訓(xùn)練，網(wǎng)格離散區(qū)間為[-5,5]，枚舉其中C和g的組合，選擇最優(yōu)的組合，使得SVR誤差最小。

第六步：反歸一化處理并對模型進行評價。對ε-SVR輸出值進行反歸一化處理求取溫度。本文利用擬合優(yōu)度R2和均方誤差函數(shù)(MSE)參數(shù)評價模型的表現(xiàn)性能，其具體公式如式(7)和(8)所示：

(7)

(8)

第七步：預(yù)測擬合CTD3的溫度值并對模型評價。將CTD3對應(yīng)的地震屬性輸入ε-SVR，預(yù)測函數(shù)如式(9)所示，對預(yù)測模型評價不僅采用擬合優(yōu)度(R2)和均方誤差函數(shù)(MSE)的方法，還分析了預(yù)測溫度值偏離真實值百分比，公式如(10)所示。

(9)

式中：L為CTD3對應(yīng)的86個樣本數(shù)目；x為一個對應(yīng)地震屬性向量的樣本。

(10)

第八步:預(yù)測整個海水溫度剖面。將所有道地震屬性向量作為樣本輸入，預(yù)測整個剖面的海水溫度，并分析效果。

本文提取地震屬性共13種：方差(var)、瞬時頻率(Freq)、均方差振幅(Rmms Ampl)、瞬時帶寬(Int)、瞬時相位(Phase)、一階導(dǎo)數(shù)(D1)、二階導(dǎo)數(shù)(D2)、時間(Time)、原始振幅(Ampl)、時間增益(Time Gain)、相對聲波阻抗(Rel Ampl)、道梯度(Trace Grad)、道增益(Iterative AGC)，對地震屬性進行主成分分析。

如圖3主成分分析表明7種屬性值(時間、原始振幅、道梯度、道增益、時間增益、一階導(dǎo)數(shù)、二階導(dǎo)數(shù))占比較重且可表達總屬性95%的信息，因此提取該7種屬性作為輸入數(shù)據(jù)集。

圖3 13種地震屬性主成分分析

使用RBF作為核函數(shù)的SVR需要考慮參數(shù)C和g，確定好C,g使得分類器能正確的預(yù)測未知數(shù)據(jù)，參數(shù)的選擇沒有先驗知識，必須做某種類型的模型選擇(參數(shù)搜索)。本次選擇網(wǎng)格交叉驗證方法，5折交叉驗證充分地利用已知信息提高預(yù)測精度，遍歷網(wǎng)格所有參數(shù)C和g，求取最優(yōu)參數(shù)C和g作為SVR的參數(shù)，其尋優(yōu)步驟為:

(1)確定實際應(yīng)用中參數(shù)的范圍：

-5≤log(C)≤5,-5≤log(g)≤5。

(2)訓(xùn)練支持向量回歸機并預(yù)測。將訓(xùn)練集樣本分為5組，支持向量回歸機模型遍歷網(wǎng)格中的C和g取值并計算每個樣本的預(yù)測值與實際標(biāo)簽的溫度誤差。

(3)尋找溫度誤差最小值，輸出最佳的組合參數(shù)。

利用網(wǎng)格交叉驗證方法編程計算得出最優(yōu)參數(shù)組合參數(shù)C和g分別為4，0.088 4。構(gòu)建支持向量回歸機模型，優(yōu)化后的地震屬性作為輸入，溫度作為輸出。本次將CTD1和CTD2進行訓(xùn)練，其訓(xùn)練效果如圖4所示。

圖4 CTD1和CTD2預(yù)測結(jié)果對比

其中紅色代表CTD溫度的真實值，藍色代表預(yù)測值。通過訓(xùn)練集可以看出其均方誤差為0.7 ℃，擬合優(yōu)度為97%，訓(xùn)練效果較好。抽取其CTD3對應(yīng)處的溫度信息與樣本溫度標(biāo)簽比對，其擬合優(yōu)度達到98%，均方誤差為1.1 ℃，吻合較好，如圖5所示。

圖5 CTD3預(yù)測結(jié)果和標(biāo)簽對比

將CTD3預(yù)測值和真實值做差，誤差如圖6所示，溫度最大差值為1.46 ℃，此時對應(yīng)的標(biāo)簽溫度為23.43 ℃，預(yù)測的溫度為21.96 ℃，其誤差最大百分比為5.92%。誤差主要來源于海水溫度差異小，CTD曲線之間差異較大以及支持向量回歸機在數(shù)值上存在四舍五入的誤差。

圖6 CTD3溫度預(yù)測誤差百分比

將CTD數(shù)據(jù)組合，對數(shù)據(jù)集劃分測試集和預(yù)測集，對組合效果進行測試，其擬合優(yōu)度平均值為0.976，誤差函數(shù)平均值為0.010，可見支持向量機預(yù)測效果良好，其CTD數(shù)據(jù)集組合效果測試如表1所示。

表1 CTD數(shù)據(jù)集組合效果測試

將CTD組合方式1訓(xùn)練好的支持向量機推廣到整個地震數(shù)據(jù)，反演得到的溫度剖面如圖7所示。

圖7 海水溫度預(yù)測剖面

可以看出，反演的海水溫度水平分層明顯，其溫度范圍大致在14.50～25.50 ℃，主要分為5層，溫度變化明顯的層對應(yīng)深度大致在3.04，9.12，15.96和26.22 m，其對應(yīng)層的溫度范圍為：25.5～24.48，24.48～22.21，22.21～19.10，19.10～16.20和16.20～14.50 ℃，計算其溫度隨深度變化的梯度，發(fā)現(xiàn)在第三層，即深度為9.12～15.96 m之間溫度變化最明顯，判斷其為溫躍層，厚度大致為6.84 m，該結(jié)果與文獻[15-16]相符。從整個剖面可以清晰的看出海水溫度分層明顯，橫向連續(xù)性較好，分辨率高。表層海水由于受太陽輻射，水溫25 ℃左右，隨著深度的增加溫度逐漸降低，由于海浪和海水流動的動壓力致使層與層之間有些許起伏，底層海水層區(qū)分明顯。

4 結(jié)論

本文利用主成分分析方法進行地震屬性優(yōu)化，利用CTD溫度標(biāo)簽數(shù)據(jù)和優(yōu)選后的地震屬性訓(xùn)練支持向量回歸機，成功的反演出海水溫度剖面。得到以下結(jié)論：

(1)利用主成分分析方法在13種地震屬性中優(yōu)選出了7種相對海水溫度敏感的地震屬性，應(yīng)用到海水溫度預(yù)測中效果明顯。

(2)支持向量回歸機在海水溫度標(biāo)簽少的情況下可以預(yù)測海水溫度，但是海水溫度本身差異小,再加上標(biāo)簽在差異性較大以及支持向量機在實際應(yīng)用中帶有一定的舍入誤差導(dǎo)致預(yù)測的海水溫度存在一定的誤差。