彭陽陽， 張 進(jìn)，2，3??

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

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

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

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

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

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

(1)

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

(2)

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

(3)

(4)

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

(5)

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

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

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

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

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

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

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

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

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

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

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

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

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

(6)

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

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

(7)

(8)

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

(9)

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

(10)

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

4 結(jié)論

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

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

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