馮春健 劉漢露 劉錦昆 賈永剛 侯 方 薛 涼 權(quán)永崢

(①中石化石油工程設(shè)計有限公司， 東營 257020， 中國) (②中國海洋大學(xué)， 青島 266100， 中國)

0 引 言

海底沉積物發(fā)生侵蝕再懸浮的主要驅(qū)動力是波浪和海流。一種侵蝕再懸浮作用機(jī)制是波浪產(chǎn)生的軌道切應(yīng)力和流致剪切應(yīng)力超過了海床表層沉積物的強(qiáng)度，使其發(fā)生再懸浮(Green et al.，2014)。另一種侵蝕再懸浮作用機(jī)制是：在波浪循環(huán)荷載的作用下，滲透性差的海床會發(fā)生孔隙水壓力的累積甚至液化(張少同等， 2016)，當(dāng)液化作用產(chǎn)生時，海床一定深度范圍內(nèi)細(xì)顆粒物被滲流輸運(yùn)到海水中(Zhang et al.，2017)，這種“瞬態(tài)泵送”過程對侵蝕再懸浮的貢獻(xiàn)率達(dá)到20%～60%(Zhang et al.，2018)。進(jìn)一步準(zhǔn)確預(yù)測不同海況下瞬態(tài)泵送對再懸浮的貢獻(xiàn)率對近岸侵蝕量的定量計算和侵蝕災(zāi)害的有效評估具有重要意義。

海底侵蝕再懸浮研究的一個重要指標(biāo)是近底海水的懸沙濃度(Suspended Sediment Concentration， SSC)。海水懸沙濃度預(yù)測的經(jīng)典方法是通過現(xiàn)場采集沉積物和海水樣品進(jìn)行室內(nèi)分析，建立物理模型和統(tǒng)計模型(Dey et al.，2018)。這些方法針對實(shí)時短期預(yù)測難以適用，尤其是實(shí)時評估預(yù)測波致瞬態(tài)泵送對再懸浮的貢獻(xiàn)率。隨著計算機(jī)和人工智能的發(fā)展，機(jī)器學(xué)習(xí)和深度學(xué)習(xí)方法在地質(zhì)災(zāi)害(鄢好等， 2019； 杜星等， 2020； 劉艷輝等， 2021)和海洋泥沙動力學(xué)研究中逐漸發(fā)揮重要作用(Mohammad et al.，2020； Huang et al.，2021)。深度學(xué)習(xí)方法可以從大量數(shù)據(jù)集中提取信息，這些信息可以作為理解深層次物理規(guī)律的佐證(Goldstein et al.，2019)。人工神經(jīng)網(wǎng)絡(luò)(ANN)已廣泛應(yīng)用于海岸泥沙動力學(xué)研究，在懸沙濃度(時間序列數(shù)據(jù))的預(yù)測(Mohammad et al.，2016)和沿岸泥沙輸移速率計算(Maanen et al.，2010)等方面表現(xiàn)突出。神經(jīng)網(wǎng)絡(luò)強(qiáng)大的學(xué)習(xí)能力也伴隨著訓(xùn)練時容易過擬合和誤差疊加的缺點(diǎn)，訓(xùn)練出來的模型只適用于當(dāng)前的數(shù)據(jù)。此外，它們需要大量的超參數(shù)調(diào)整，這大大增加了訓(xùn)練所需的時間。

長短時記憶(Long Short-Term Memory， LSTM)循環(huán)神經(jīng)網(wǎng)絡(luò)由一個輸入門、遺忘門和輸出門組成，具有不同于一般神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)，可以解決上述問題(Hochreiter et al.， 1997； Kaveh et al.，2020)。LSTM神經(jīng)網(wǎng)絡(luò)可以用特殊的結(jié)構(gòu)將無意義的信息舍棄，從而防止了梯度下降引起的梯度消失或梯度爆炸的現(xiàn)象，解決神經(jīng)網(wǎng)絡(luò)應(yīng)用于時間序列預(yù)測中所帶來的誤差疊加問題。本文應(yīng)用2016～2017年冬季在中國黃河水下三角洲觀測到的懸沙濃度、瞬態(tài)泵送和水動力學(xué)數(shù)據(jù)集訓(xùn)練和測試LSTM深度學(xué)習(xí)預(yù)測模型，并與兩種常用的時間序列預(yù)測模型(支持向量回歸和人工神經(jīng)網(wǎng)絡(luò))進(jìn)行性能比較。

1 研究區(qū)和數(shù)據(jù)集

為了訓(xùn)練和測試波致海床瞬態(tài)液化對再懸浮貢獻(xiàn)率的深度學(xué)習(xí)預(yù)測模型，于2016～2017年冬季在黃河水下三角洲地區(qū)的埕島海域進(jìn)行了海水懸沙濃度和水動力條件的原位長期觀測，觀測點(diǎn)(118.900 255 2°E， 38.166 741 14°N)位置如圖 1 中所示。

圖 1 研究區(qū)地理位置：(a)現(xiàn)代黃河水下三角洲地理位置(b)觀測點(diǎn)位置(c)觀測所用到的儀器設(shè)備Fig. 1 Location of the study area：(a) Geographical location of the underwater delta of the modern Yellow River(b) Position of observation point(c) Instruments and equipment used for observation

1.1 研究區(qū)概況

現(xiàn)代黃河三角洲形成于1855年黃河改道后，研究區(qū)所在的埕島海域位于中國半封閉式渤海西南部，在渤海灣和萊州灣之間( 圖 1b) 。研究區(qū)位于濟(jì)陽斷陷區(qū)，地形較平坦，表層沉積物由粉砂(62.3%)和黏土(37.4%)組成，平均粒徑為0.043 mm。冬季風(fēng)暴事件的經(jīng)常發(fā)生導(dǎo)致海底侵蝕、液化等地質(zhì)災(zāi)害頻繁(楊作升， 1993)。1976年黃河改道至清水溝后，刁口葉瓣廢棄，河流泥沙供應(yīng)減少(Xue et al.，1993)，北岸發(fā)生了強(qiáng)烈的風(fēng)浪侵蝕再懸浮和泥沙輸運(yùn)，岸線不斷后退。

影響再懸浮量(懸沙濃度)變化的機(jī)制主要包括河流輸入、波浪或水流引起的局部再懸浮和沉降過程(Maa et al.， 2007)。在本研究區(qū)，由于刁口葉瓣廢棄，河流輸入導(dǎo)致懸沙濃度變化可以忽略，波浪是造成黃河水下三角洲發(fā)生大量侵蝕再懸浮的控制因素(張少同等， 2016)。波浪對海床沉積物的作用不僅包括軌道切應(yīng)力，海平面的波峰波谷循環(huán)交替變化對海床施加垂向循環(huán)荷載，會導(dǎo)致瞬態(tài)和殘余兩種孔壓響應(yīng)(Jeng， 2003)，對應(yīng)著兩種不同的滲流模式——瞬態(tài)滲流和累積滲流。波峰波谷的垂向振蕩過程中，波峰時在海床中引起下壓力，波谷時在海床中產(chǎn)生上吸力，導(dǎo)致海床和海水的界面處發(fā)生垂向的瞬態(tài)滲流。波致海床孔壓所引發(fā)的垂向滲流，在海床液化之前，強(qiáng)度往往較小，但持續(xù)的常態(tài)波浪作用或驟增的風(fēng)暴浪引起海床瞬態(tài)液化后，會發(fā)生“瞬態(tài)泵送”現(xiàn)象，即海床中的細(xì)顆粒會被輸運(yùn)到海床表面，進(jìn)一步被水流搬運(yùn)、懸浮，這一過程對再懸浮的貢獻(xiàn)量不可忽視，可達(dá)20%～60%(Zhang et al.， 2018)。

1.2 數(shù)據(jù)集

在觀測點(diǎn)位置(118.9002552°E， 38.16674114°N)的海床表面布放了四腳架平臺和海底實(shí)驗艙(圖 1)，四腳架平臺上和海底實(shí)驗艙內(nèi)各裝有一臺波潮儀(virtuoso， RBR，加拿大，內(nèi)置壓力傳感器)和一臺濁度計(Concerto RBR，加拿大，內(nèi)置光學(xué)后向散射傳感器)。海水濁度和波浪參數(shù)的記錄分別以1次/3 s的連續(xù)采集模式和5 min·h-1(6 Hz)的burst采樣模式進(jìn)行。連續(xù)、同步測量2016年12月18日至2017年1月11日的水深、波浪、海水濁度以及海底實(shí)驗艙內(nèi)的波浪、海水濁度參數(shù)。

為了分析瞬態(tài)泵送對再懸浮的貢獻(xiàn)率，需要得到實(shí)驗艙內(nèi)外同一高度的海水懸沙濃度。四腳架搭載的海底邊界層懸浮物濃度剖面儀(Argus Surface Meter IV，Argus公司，德國)，通過內(nèi)置的OBS傳感器陣列(間距1 cm)以1次/30 min的burst模式來記錄近底邊界層內(nèi)(距底0～1.8 m)的懸沙濃度剖面。

本文數(shù)據(jù)集包括的參數(shù)有：水深(D)、有效波高(Hs)、有效波周期(Ts)、實(shí)驗艙內(nèi)懸沙濃度(SSCin)、實(shí)驗艙外懸沙濃度(SSCout)。為了使所有時間序列的長度都相同，同時消除測量噪聲，對所有數(shù)據(jù)每半小時取平均。上述5個參數(shù)作為自變量，懸沙濃度貢獻(xiàn)率(C)為目標(biāo)變量，定義為因變量。建立基于LSTM的深度學(xué)習(xí)預(yù)測模型需要用到的時間序列數(shù)據(jù)如圖 2 所示。

圖 2 LSTM的訓(xùn)練和預(yù)測數(shù)據(jù)Fig. 2 Training and prediction data of LSTM

2 研究方法

2.1 LSTM神經(jīng)網(wǎng)絡(luò)

循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network， RNN)是一類用于處理序列數(shù)據(jù)的神經(jīng)網(wǎng)絡(luò)(Rumelhart et al.，1986)。相比于ANN的優(yōu)點(diǎn)在于，在獲取每個神經(jīng)元的輸出的同時，也將其作為輸入再次反饋給神經(jīng)元。這種方式不僅在每個時間步長中接收新的信息，而且還向這些新信息中添加先前輸出的加權(quán)值。這樣RNN的神經(jīng)元就具備了先前輸入的“記憶”。長短時記憶(Long Short-Term Memory， LSTM)循環(huán)神經(jīng)網(wǎng)絡(luò)是一種特殊類型的RNN，可以規(guī)避常規(guī)RNN中梯度爆炸和梯度消失的問題(Hochreiter et al.， 1997)。LSTM網(wǎng)絡(luò)由輸入層、隱含層和輸出層組成，其結(jié)構(gòu)類似于RNN。但RNN與LSTM的不同之處在于，LSTM用一個記憶單元取代了常規(guī)RNN的基本單元(Graves et al.，2013)。如圖 3 所示，該記憶單元包含3個門函數(shù)(輸入門、遺忘門和輸出門)，而且定義了“細(xì)胞狀態(tài)”(ct)。門函數(shù)可以在訓(xùn)練中記錄過去的關(guān)鍵特征，并根據(jù)權(quán)重選擇不重要的特征來遺忘，可以輕松地讀取和更新長時間序列的信息。

圖 3 用于數(shù)據(jù)預(yù)測的長短時記憶循環(huán)神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig. 3 The structure of short-and long-term memory cyclic neural network for data prediction

在t時刻，xt為LSTM細(xì)胞的輸入數(shù)據(jù)(多參數(shù)時則是多維度矩陣)；ht-1為LSTM單元在t-1時刻的輸出；ct為t時刻細(xì)胞狀態(tài)值，用來保存重要信息；ht是t時刻LSTM單元的輸出。LSTM單元的計算過程包括以下步驟：

首先，細(xì)胞狀態(tài)候選值

(1)

式中：Wc為權(quán)重矩陣；bc為偏差。

輸入門的值it用式(2)計算，輸入門控制將當(dāng)前輸入數(shù)據(jù)更新為細(xì)胞的狀態(tài)值。

it=σ(Wi· [ht-1，xt]+bi)

(2)

式中： σ為sigmoid函數(shù)；Wi為權(quán)重矩陣；bi為偏差。

遺忘門用來控制遺忘上一層細(xì)胞狀態(tài)的內(nèi)容，根據(jù)上一序列的ht-1和本序列的xt為輸入，通過sigmoid激活函數(shù)，判斷上一層需要保留或去除的細(xì)胞狀態(tài)內(nèi)容。遺忘門的值ft用式(3)計算。

ft=σ(Wf· [ht-1，xt]+bf)

(3)

式中：Wf為權(quán)重矩陣；bf為偏差。

用遺忘門和輸入門將t-1時刻細(xì)胞狀態(tài)ct-1更新為t時刻細(xì)胞狀態(tài)ct：

(4)

最后要基于細(xì)胞狀態(tài)保存的內(nèi)容來確定輸出什么內(nèi)容。即選擇地的輸出細(xì)胞狀態(tài)保存的內(nèi)容。輸出門的值ot用式(5)計算，輸出門控制存儲單元狀態(tài)值的輸出。

ot=σ(Wo·[ht-1，xt]+bo

(5)

式中：Wo為權(quán)重矩陣；bo為偏差。

LSTM單元的最終輸出：

ht=ot×tanh(ct)

(6)

由于上述這種LSTM內(nèi)部參數(shù)的共享機(jī)制，可以通過設(shè)置權(quán)重矩陣的維數(shù)來控制輸出的維數(shù)。在LSTM神經(jīng)單元中，輸入之間有一個很長的延遲和反饋。梯度既不會爆炸也不會消失，因為LSTM體系結(jié)構(gòu)中存儲單元的內(nèi)部狀態(tài)保持著一個恒定的錯誤流。

在神經(jīng)網(wǎng)絡(luò)訓(xùn)練中，防止過擬合是至關(guān)重要的。Srivastava et al. (2014)提出了Dropout方法來防止神經(jīng)網(wǎng)絡(luò)訓(xùn)練中的過擬合。在神經(jīng)網(wǎng)絡(luò)的層與層之間的訓(xùn)練過程中，一些神經(jīng)元以一定的概率隨機(jī)從網(wǎng)絡(luò)中脫落。

2.2 LSTM 預(yù)測模型

2.2.1 數(shù)據(jù)預(yù)處理——?dú)w一化

因為LSTM遞歸神經(jīng)網(wǎng)絡(luò)的轉(zhuǎn)換函數(shù)是雙曲正切函數(shù)(tanh函數(shù)：ex-e-x/(ex+e-x)其函數(shù)值范圍為- 1～1)，訓(xùn)練集和測試集數(shù)據(jù)歸一化。歸一化方程如下：

(7)

式中：xmax和xmin分別為數(shù)據(jù)集(訓(xùn)練集和測試集)的最大值和最小值；x為輸入值；z為x的轉(zhuǎn)換值。

2.2.2 分割訓(xùn)練集和測試集

將數(shù)據(jù)集拆分為訓(xùn)練集和測試集，這里將前9天的數(shù)據(jù)作為訓(xùn)練集，后3天的數(shù)據(jù)作為測試集。然后將訓(xùn)練集和測試集分別拆分為輸入和輸出變量。最后將輸入變量(X)轉(zhuǎn)變成LSTM需要的三維格式，即[samples， timesteps， features]。

2.2.3 定義和擬合LSTM模型

為了建立預(yù)報模型，本研究采用了兩層LSTM單元，預(yù)測值在第3層輸出。LSTM的超參數(shù)包括隱含層數(shù)、神經(jīng)元數(shù)、學(xué)習(xí)率、激活函數(shù)、批處理規(guī)模、時間點(diǎn)數(shù)(epoch)和損失函數(shù)。神經(jīng)網(wǎng)絡(luò)的超參數(shù)決定了網(wǎng)絡(luò)如何運(yùn)行，進(jìn)而決定了網(wǎng)絡(luò)的準(zhǔn)確性和有效性。

在本研究中，根據(jù)調(diào)整后的超參數(shù)值對模型性能的影響，手動調(diào)整LSTM超參數(shù)進(jìn)行模型優(yōu)化。超參數(shù)設(shè)置如下，學(xué)習(xí)率為0.0001，兩層LSTM神經(jīng)元數(shù)量均為100個，輸出層設(shè)置一個神經(jīng)元，也作為預(yù)測值yt的輸出。Dropout設(shè)置為0.5。所有的激活函數(shù)都使用了整流線性單元(ReLus)，因為ReLus除了計算過程簡單之外，還可以有效地進(jìn)行梯度下降和反向傳遞，從而防止梯度爆炸和消失。

2.3 機(jī)器學(xué)習(xí)

2.3.1 支持向量回歸(Support Vector Regression， SVR)

支持向量回歸是由Vapnik et al. (1997)提出的。支持向量回歸算法包括不敏感損失函數(shù)和懲罰因子函數(shù)； 因此，它比支持向量機(jī)算法更具有魯棒性。SVR算法是將數(shù)據(jù)投影到一個高維超平面后，計算每個點(diǎn)到超平面的總距離。如果一個超平面是用最小總距離確定的，那么這個超平面就是解。SVR是預(yù)測整個時間序列最常用的方法。

2.3.2 人工神經(jīng)網(wǎng)絡(luò)(Artificial Neural Network， ANN)

反向傳播網(wǎng)絡(luò)(Back Propagation Network， BPN)是目前最常用的監(jiān)督學(xué)習(xí)ANN模型(Rumelhart et al.，1986)。BPN是一種結(jié)合了逆向傳遞、梯度下降和鏈?zhǔn)椒▌t的優(yōu)化算法。采用梯度下降法從參數(shù)的初始位置推進(jìn)到最陡的下坡方向，并更新參數(shù)位置。斜率信息通過導(dǎo)數(shù)函數(shù)(計算函數(shù)斜率)得到。梯度下降利用這一特性來優(yōu)化損失函數(shù)。

2.4 性能評價標(biāo)準(zhǔn)

采用均方根誤差(Root Mean Squared Error， RMSE)、平均絕對百分比誤差(Mean Absolute Percentage Error， MAPE)、平均平方誤差-標(biāo)準(zhǔn)偏差(Root Mean Squared error-Standard Deviation Ratio， RSR)這3個常用的統(tǒng)計變量來比較LSTM的真實(shí)值與預(yù)測值的偏差，以評價LSTM的預(yù)測性能。RMSE、MAPE、RSR分別表示為式(8)～式(10)。

(8)

(9)

(10)

平均平方誤差-標(biāo)準(zhǔn)偏差在0～1之間，越接近0的情況下模型性能越好，一般認(rèn)為RSR在0.7以下時模型能很好地預(yù)測實(shí)測值。

3 結(jié) 果

3.1 預(yù)測結(jié)果

圖 4和圖 5 用不同的形式展示了SVR模型、ANN模型、LSTM模型做的瞬態(tài)泵送貢獻(xiàn)率實(shí)際值和預(yù)測值的對比。從圖 4 中，散點(diǎn)圖的擬合直線斜率和R2越接近1，說明預(yù)測的誤差越小。LSTM模型和ANN模型具有相似的精度，SVR模型的預(yù)測比較分散，誤差更大。

圖 4 83.5 h內(nèi)實(shí)測與預(yù)報貢獻(xiàn)率的散點(diǎn)圖Fig. 4 Scatter plot of the measured and forecast contribution rate within 83.5 hoursa. SVR; b. ANN; c. LSTM

圖 5 各模型預(yù)測結(jié)果Fig. 5 Prediction results of each modela. SVR; b. ANN; c. LSTM

應(yīng)用MAPE值、RMSE值和RSR值評價3種模型的性能。由表 1 可以看出，LSTM模型、ANN模型的誤差值要比SVR模型小得多。LSTM模型的MAPE值分別比SVR模型和ANN模型的MAPE值低85.31%、70.11%； LSTM模型的RMSE值分別比SVR模型和ANN模型的RMSE值低83.07%、63.97%； LSTM模型的RSR值分別比SVR模型和ANN模型的RSR值低79.09%、55.49%。這些結(jié)果表明，LSTM模型顯著優(yōu)于SVR模型和ANN模型。

表 1 SVR，ANN，LSTM 模型的MAPE值、RMSE值、RSR值Table 1 MAPE values，RMSE values and RSR values of SVR，ANN and LSTM models

3.2 SVR，ANN，LSTM模型不同預(yù)測時間效果對比

從誤差變化的趨勢上來看，圖 6中顯示的是MAPE和RMSE隨不同預(yù)測時間的變化。紅色實(shí)線代表RMSE，黑色實(shí)線代表MAPE。SVR和LSTM的預(yù)測時間為36 h的條件下，誤差明顯增大； 而后SVR和LSTM的MAPE逐漸減小，但ANN和LSTM的RMSE逐漸增大。ANN的RMSE和MAPE逐漸減小。從數(shù)值大小來看，如圖 7 和圖 8，SVR和ANN的兩種誤差在12 h和24 h相差不大，但遠(yuǎn)大于LSTM； 后面預(yù)測時間增加，誤差遠(yuǎn)均大于LSTM和ANN。ANN和LSTM相比，LSTM的兩種誤差均小于ANN。此外，RMSE雖然在48 h后有下降趨勢，但是變化量和絕對值小于ANN。由此也可以證明，LSTM的性能優(yōu)于ANN。

圖 6 各模型不同預(yù)測時長結(jié)果對比Fig. 6 Comparison of the results of different prediction durations of different modelsa. SVR; b. ANN; c. LSTM

圖 7 各模型不同預(yù)測時間MAPE的值Fig. 7 MAPE values of various models at different prediction times

圖 8 各模型不同預(yù)測時間RMSE的值Fig. 8 RMSE values of each model at different prediction times

4 結(jié) 論

本文建立了基于LSTM的瞬態(tài)液化對再懸浮貢獻(xiàn)率的深度學(xué)習(xí)預(yù)測模型，利用黃河水下三角洲埕島海域的現(xiàn)場觀測數(shù)據(jù)集對該模型進(jìn)行了訓(xùn)練和測試，并與SVR、ANN時間序列預(yù)測模型的性能進(jìn)行了比較。結(jié)果表明：

(1)LSTM模型顯著優(yōu)于SVR模型和ANN模型。與其他模型相比，LSTM模型在訓(xùn)練過程中除了誤差值更小外，還表現(xiàn)出更高的穩(wěn)定性和更快的收斂速度。

(2)LSTM模型對3.5 d以內(nèi)的瞬態(tài)泵送再懸浮貢獻(xiàn)率預(yù)報誤差最小，其MAPE、RMSE和RSR均值分別為5.87%、1.6730、0.1574，均明顯低于其他預(yù)測模型。在所有的預(yù)測區(qū)間內(nèi)，LSTM模型的R2值在所有比較模型中均最高。