石綏祥，王蕾,2*，余璇，徐凌宇

( 1. 國家海洋信息中心 數(shù)字海洋實驗室，天津 300171；2. 國家海洋局東海信息中心，上海 200136；3. 上海大學(xué) 計算機(jī)工程與科學(xué)學(xué)院，上海 200444)

1 引言

葉綠素a 濃度作為水質(zhì)狀況及富營養(yǎng)化程度的衡量指標(biāo)，一直是水質(zhì)監(jiān)測的重要參數(shù)[1–3]，是判斷赤潮發(fā)生與否的重要因素，因此越來越多的科研人員投入到對葉綠素a 濃度的預(yù)測研究中。而限于早期葉綠素a 濃度相關(guān)數(shù)據(jù)的匱乏，對葉綠素a 濃度的預(yù)測大致可以分為兩類方法，其中一類方法為統(tǒng)計學(xué)方法，這一類的方法最先由加拿大專家Vollenweider[4]提出，他利用了統(tǒng)計模型來預(yù)測富營養(yǎng)化問題，但是這一類傳統(tǒng)的統(tǒng)計方法只能求解某一要素的平均濃度分布，無法模擬相關(guān)因素與葉綠素a 的影響關(guān)系；另一類方法主要是以生態(tài)動力學(xué)為理論依據(jù)，基于對流?擴(kuò)散方程建立模型[5]來預(yù)測葉綠素a 濃度，此類方法的優(yōu)點是考慮了自然界中多種因素間的相互作用，但生態(tài)動力學(xué)模型有一個共同的缺陷是其包含太多參數(shù)，而這些參數(shù)的設(shè)定十分依賴對于特定問題的經(jīng)驗，并且難以確定合適的參數(shù)。

隨著數(shù)據(jù)量和種類的增多，收集的數(shù)據(jù)也越來越多，傳統(tǒng)的統(tǒng)計學(xué)習(xí)方法越來越難以適應(yīng)，鑒于神經(jīng)網(wǎng)絡(luò)及深度學(xué)習(xí)在大數(shù)據(jù)量及復(fù)雜場景中，對于特征的表達(dá)能力，越來越多的學(xué)者開始利用神經(jīng)網(wǎng)絡(luò)及深度學(xué)習(xí)對葉綠素a 進(jìn)行預(yù)測。人工神經(jīng)網(wǎng)絡(luò)是一種模擬生物神經(jīng)網(wǎng)絡(luò)信息處理功能的信息處理模型，也是一個高度復(fù)雜的非線性動力學(xué)系統(tǒng)[6]。在此類研究中，最具代表性的人工神經(jīng)網(wǎng)絡(luò)模型是反向傳播神經(jīng)網(wǎng)絡(luò)(Back Propagation Neural Network, BPNN)，趙玉芹等[7]利用BP 神經(jīng)網(wǎng)絡(luò)成功地對渭河水質(zhì)參數(shù)進(jìn)行了遙感反演；盧志娟等[8]利用BP 神經(jīng)網(wǎng)絡(luò)實現(xiàn)了對西湖湖心區(qū)葉綠素a 濃度的周預(yù)測；周露洪等[9]通過對2006?2008 年的常規(guī)水質(zhì)參數(shù)進(jìn)行主成分分析，建立BP 神經(jīng)網(wǎng)絡(luò)模型對葉綠素a 濃度進(jìn)行月預(yù)測；楊柳等[10]利用BP 神經(jīng)網(wǎng)絡(luò)對溫榆河進(jìn)行了水質(zhì)參數(shù)反演，反演結(jié)果優(yōu)于傳統(tǒng)線性回歸模型。應(yīng)用廣泛的BP 網(wǎng)絡(luò)大多只考慮相關(guān)因素對葉綠素a 的影響關(guān)系，而忽略了所有影響葉綠素a 的因素在不同時延上影響關(guān)系的差異性，無法動態(tài)的依賴歷史時序信息進(jìn)行葉綠素a 濃度預(yù)測，并存在局部極小、收斂時間過長從而泛化能力差的問題[11]。

Dekker 等[12]建立了葉綠素a 濃度與TM 影像的線性和指數(shù)回歸模型，指出指數(shù)模式要優(yōu)于線性模式，但TM 的分辨率較低，不利于水質(zhì)參數(shù)監(jiān)測，且T 的波段組合缺乏物理解釋；柴永強(qiáng)等[13]利用決策樹模型智能檢測赤潮現(xiàn)象，采用機(jī)器學(xué)習(xí)技術(shù)訓(xùn)練檢測赤潮的決策樹分類模型，此模型對渤海等海域取得了較滿意的結(jié)果，但是決策樹模型對缺失數(shù)據(jù)處理比較困難，也容易出現(xiàn)過擬合現(xiàn)像；李修竹等[14]采用支持向量機(jī)的方法，以溫度、鹽度等8 種要素作為輸入，葉綠素a 濃度作為輸出對葉綠素a 濃度進(jìn)行預(yù)測取得了較好的結(jié)果，此種方法雖然避免了神經(jīng)網(wǎng)絡(luò)方法的局部極小值問題，但是對核函數(shù)的高維映射解釋力欠缺。

基于目前的研究成果，本文在海洋數(shù)據(jù)多要素關(guān)聯(lián)關(guān)系的基礎(chǔ)上，提出了一種結(jié)合時序方法的遞歸神經(jīng)網(wǎng)絡(luò)智能預(yù)測模型，對判斷藻類赤潮的重要指標(biāo)葉綠素a 濃度進(jìn)行預(yù)測，以解決傳統(tǒng)的遞歸神經(jīng)網(wǎng)絡(luò)的梯度消失或爆炸問題。

本文的創(chuàng)新點主要集中在兩個方面：（1） 量化不同要素對于葉綠素a 濃度在不同時延上的關(guān)系，依據(jù)長短期影響關(guān)系將其分類；（2） 對于不同的類別要素，構(gòu)建差異性的子網(wǎng)絡(luò)對其分別建模，每個子網(wǎng)絡(luò)采用不同的訓(xùn)練方式，利用融合層將不同子網(wǎng)絡(luò)的特征進(jìn)行融合以得到更為穩(wěn)定準(zhǔn)確的結(jié)果。

2 數(shù)據(jù)來源與研究方法

2.1 數(shù)據(jù)來源

近年來，三都澳水產(chǎn)養(yǎng)殖的開發(fā)利用幾乎是直線上升，宜養(yǎng)面積的利用超過100%，由于三都澳是一個口小腹大的內(nèi)灣，一般海水進(jìn)行一次完全循環(huán)要1 周，垃圾在海面上漂流很長時間，沉積慢，嚴(yán)重影響了海區(qū)環(huán)境衛(wèi)生，影響海區(qū)水交換，使水域污染日益嚴(yán)重，超過水域的自凈能力，導(dǎo)致海水富營養(yǎng)化加劇，病害發(fā)生頻繁，因此本文選取三都澳的連續(xù)監(jiān)測資料作為本文的實驗數(shù)據(jù)，其中三都澳站位數(shù)據(jù)時間段是從2015 年5 月12 日至2015 年7 月2 日，監(jiān)測要素包括表層溫度（WD）、PH、電導(dǎo)率（DDL）、葉綠素（YLS）、濁度（ZD）、溶解氧（RJY）、鹽度（YD）、電壓（DY）共8 個，監(jiān)測頻率的時間間隔為5 min 一次，站位分布情況如圖1 所示，圖中紅色位置是三都澳站點所在位置。

2.2 研究方法

本文利用長短期記憶（Long Short-Term Memory,LSTM）模型強(qiáng)大的長短期記憶能力分析三都澳站位水質(zhì)監(jiān)測數(shù)據(jù)，選擇對三都澳水質(zhì)影響較大的因子網(wǎng)絡(luò)模型，預(yù)測三都澳站位的葉綠素a 濃度的變化趨勢，從而為有效及時地控制水質(zhì)提供科學(xué)依據(jù)，總體研究方法由圖2 所示。

從圖2 中可以看出本文提出的模型一共分為3 層：數(shù)據(jù)分析層、模型層和決策層。其中數(shù)據(jù)分析層為原始數(shù)據(jù)的預(yù)處理，第二層模型層包含自相關(guān)性分析、多要素與葉綠素a 的關(guān)聯(lián)關(guān)系分析、長短期依賴時間分析以及構(gòu)建融合的LSTM 預(yù)測模型，第三層通過對預(yù)測分析結(jié)果的反饋來提供的有效分析決策。

3 數(shù)據(jù)預(yù)處理

3.1 數(shù)據(jù)歸一化方法

由于各要素原始數(shù)據(jù)的量綱及數(shù)量級不同，在進(jìn)行數(shù)據(jù)相關(guān)性分析和網(wǎng)絡(luò)訓(xùn)練前要先對數(shù)據(jù)進(jìn)行歸一化處理。

圖 1 三都澳站位分布Fig. 1 Station location distribution map of Sandu Ao

經(jīng)過標(biāo)準(zhǔn)化處理后的數(shù)據(jù)一般在0 到1 之間，這樣有利于網(wǎng)絡(luò)訓(xùn)練，在量綱統(tǒng)一的基礎(chǔ)上為下一步相關(guān)性分析和網(wǎng)絡(luò)訓(xùn)練奠定基礎(chǔ)。

3.2 自相關(guān)分析方法

在進(jìn)行葉綠素a 濃度預(yù)測前，首先需要知道葉綠素的自相關(guān)情況，以便在后面的網(wǎng)絡(luò)訓(xùn)練和預(yù)測過程中如何考慮其自相關(guān)情況。

相隔 k 期的兩個變量 xt與 xt?k的協(xié)方差即滯后 k期的自協(xié)方差可以定義為：程的自協(xié)方差函數(shù)。當(dāng) k =0 時， γ0=Var(xt)=自相關(guān)系數(shù)的定義如下：

4 構(gòu)建融合的（Merged）LSTM 模型

4.1 多要素間關(guān)聯(lián)分析方法

m×n

散布矩陣為 的半正定矩陣如下：

其中 T 表示矩陣的轉(zhuǎn)置，散布矩陣可以簡要的表示為S=XCnXTCnXT，在此 Cn定 義為定心矩陣，其中 Cn公式如下：

式中，O 表示所有元素都是1 的矩陣；在最大似然估計中，給定 n個樣本，一個多元正太分布的協(xié)方差可以表示為歸一化的散度矩陣：

4.2 長短期依賴時間分析法

X 、 Y （也可以看作兩個集合），它們的元素個數(shù)均為N （ N 表示時間序列的長度），兩個要素取的第i(1 ≤i ≤N)個值分別用 Xi、 Yi表示。對 X 、 Y進(jìn)行排序（同時為升序或降序），得到兩個元素排行集合 x 、 y ， 其中元素 xi、

圖 2 研究方法圖Fig. 2 Research method map

本文對不同海洋要素間關(guān)系采用相關(guān)系數(shù)方法計算出不同時間延遲下的相關(guān)系數(shù)大小，確定最優(yōu)時間間隔，相同時間間隔的作為一個整體劃分。對于葉綠素a 要素與其他任意一個要素的時間序列定義為yi分別為 Xi在 X 中的排行以及 Yi在 Y中的排行。將集合x 、 y中 的元素對應(yīng)相減得到一個排行差分集合d，其中di=xi?yi， 1 ≤i ≤N 。隨機(jī)變量 X 、 Y之間的斯皮爾曼等級相關(guān)系數(shù)可以由x 、 y或 者 d計算得到，其計算方式如下所示

相關(guān)關(guān)系是一種非確定性的關(guān)系，相關(guān)系數(shù)是研究變量之間線性相關(guān)程度的量。我們將相關(guān)系數(shù)定義為：

4.3 LSTM 神經(jīng)元

LSTM 由Hochreiter 和 Schmidhuber(1997)提出[15–16]，隨后Alex Graves 對算法進(jìn)行了改良和推廣[17]，基于改進(jìn)的LSTM 相關(guān)網(wǎng)絡(luò)結(jié)構(gòu)得到了廣泛的應(yīng)用。例如，文獻(xiàn)[18]將降雨量的預(yù)報作為一個時空序列預(yù)報問題，應(yīng)用LSTM 結(jié)合CNN 來預(yù)測一個地區(qū)在相對較短時間內(nèi)的未來降雨強(qiáng)度；文獻(xiàn)[19]將海水海表面溫度（Sea Surface Temperature, SST）預(yù)測作為一個序列預(yù)測問題，建立了一個端到端可訓(xùn)練的LSTM 神經(jīng)網(wǎng)絡(luò)模型，從時空角度利用歷史數(shù)據(jù)預(yù)測未來的SST 值，各像素的局部相關(guān)性和整體性可以通過固定尺寸的板塊來表達(dá)和保持。LSTM 通過遺忘門和輸出激活功能的設(shè)計來處理信息的長短期記憶，其神經(jīng)元結(jié)構(gòu)如圖3 所示，它包含一個動態(tài)的門機(jī)制，由輸入門、輸出門、遺忘門和記憶單元組成[20–22]。

圖4 詳細(xì)的描述了LSTM 內(nèi)部的數(shù)據(jù)流，其中遺忘門讀取上一個狀態(tài) h _(t ?1) 和 當(dāng)前輸入狀態(tài) x_t的信息，通過 S igmoid層輸出一個在0 到1 之間的數(shù)值給每個細(xì)胞狀態(tài) C_(t ?1)，C _(t ?1)中的數(shù)字決定從細(xì)胞狀態(tài)中丟棄什么信息，1 代表“完全保留”，0 代表“完全舍棄”。我們通過輸入層來決定什么樣的新信息將被更新并且放在細(xì)胞狀態(tài)中，首先將 h_(t ?1) 與 x_t輸入Sigmoid 函數(shù)確定將要更新的值，然后通過 t anh 層創(chuàng)建候選值向量隨后將舊狀態(tài)與 ft相乘，確定我們需要遺忘的信息，加上與的乘積產(chǎn)生新的候選值，最終我們根據(jù)新的細(xì)胞狀態(tài)來決定輸出什么值，通過Sigmoid 層決定輸出的細(xì)胞狀態(tài)，將細(xì)胞狀態(tài)通過tanh層進(jìn)行處理并將其與 Sigmoid的輸出相乘得到這一時間的輸出，可以形式化的描述如下：

4.4 Merged-LSTM 時間序列學(xué)習(xí)結(jié)構(gòu)

根據(jù)前面數(shù)據(jù)分析結(jié)果構(gòu)建融合的LSTM 時間序列學(xué)習(xí)模型，模型結(jié)構(gòu)如圖5 所示。模型的輸入數(shù)據(jù)分為兩類，其中第一類是根據(jù)4.2 節(jié)長短期依賴時間分析方法得到的短期依賴要素，標(biāo)記為OS，即在5 d 內(nèi)相關(guān)系數(shù)相似的因素類，另一類是根據(jù)長短期依賴時間分析方法得到的長期依賴要素，記為BS，即在5 d 到15 d 時間內(nèi)相關(guān)系數(shù)相似的因素類。如圖5所示，從模型的第二層開始，我們分別對BS 和OS 兩類不同的依賴關(guān)系要素采用不同的結(jié)構(gòu)進(jìn)行訓(xùn)練，其中一個用于訓(xùn)練OS 序列，使用一層的LSTM 來進(jìn)行短期依賴記憶；另一個用于專門訓(xùn)練BS 的時序依賴關(guān)系，使用兩層的LSTM 來進(jìn)行長期依賴記憶，之后使用一層Merge Layer 將兩類經(jīng)過LSTM 訓(xùn)練之后的數(shù)據(jù)合并起來，最后將輸出數(shù)據(jù)送入Rectified Linear Units 進(jìn)行激活，以便達(dá)到更快的收斂速度，最后通過輸出層Linear 線性輸出[23]得到最終預(yù)測結(jié)果。

圖 3 LSTM 神經(jīng)元結(jié)構(gòu)Fig. 3 Neuron Structure of LSTM

圖 4 LSTM 內(nèi)部計算流程圖Fig. 4 LSTM internal calculation flow chart

圖 5 Merged-LSTM 結(jié)構(gòu)圖Fig. 5 The structure of Merged-LSTM

網(wǎng)絡(luò)采用如下交叉熵公式：

其中 y 是我們預(yù)測的分布， y′是真實的分布，交叉熵作為指標(biāo)來衡量模型的預(yù)測用于描述真相的好壞，在這里我們通過梯度下降方法最小化交叉熵來使模型的輸出更符合真實分布。

5 實驗結(jié)果與分析

本文采用2.1 節(jié)的數(shù)據(jù)源去除無效數(shù)據(jù)，然后對數(shù)據(jù)進(jìn)行歸一化和自相關(guān)性分析，然后利用多要素間關(guān)聯(lián)分析方法得到要素之間的關(guān)聯(lián)關(guān)系如圖6 所示。

圖 6 要素間相似矩陣圖Fig. 6 Similarity matrix diagram between elements

圖6 的每一個小圖都表示任意兩個要素之間的相關(guān)程度，其中對角線所示的子圖為每個要素的自相關(guān)曲線，我們可以觀察葉綠素a 與其他各要素之間的關(guān)系。從圖中可以看到，葉綠素a 與其他要素間的相關(guān)散布矩陣多呈現(xiàn)有規(guī)則的點云形狀，形狀越規(guī)則，相關(guān)性越強(qiáng)。從分析結(jié)果可以看出各個屬性與葉綠素a 之間都存在一定的相關(guān)性，因此本文選取全部8 個要素作為預(yù)測葉綠素a 濃度的基礎(chǔ)數(shù)據(jù)。

根據(jù)選取的所有要素，利用長短期依賴時間分析法計算出其他要素與葉綠素a 濃度在不同時延下的相關(guān)系數(shù)，見表1 所示。

表1 為葉綠素a 濃度與各個要素在不同時延下的相關(guān)系數(shù)，每一列代表5 個數(shù)據(jù)的時延關(guān)系，根據(jù)計算出的相關(guān)系數(shù)，我們可以把采集到的不同維度進(jìn)行分類，在本文中把相關(guān)系數(shù)絕對值0.2 以下的定義為短期依賴關(guān)系，0.2 以上的定義為長期依賴關(guān)系，這樣可以把這8 個要素劃分為兩類，一類是與葉綠素a 濃度具有長期依賴關(guān)系，另一類與葉綠素a 濃度具有短期依賴關(guān)系。不同的依賴關(guān)系分別用不同的神經(jīng)網(wǎng)絡(luò)單元及層數(shù)訓(xùn)練，然后在融合層進(jìn)行信息融合，為模型建立提供依據(jù)。

本文模型中OS 子網(wǎng)絡(luò)中，每個隱藏層LSTM 神經(jīng)元個數(shù)設(shè)置為64，BS 子網(wǎng)絡(luò)中，每個隱藏層LSTM神經(jīng)元個數(shù)設(shè)置為128，參數(shù)利用隨機(jī)正態(tài)分布進(jìn)行初始化，初始學(xué)習(xí)率設(shè)置為0.01，學(xué)習(xí)率按照lrt=lr0/(1+kt)進(jìn)行遞減，其中 k為控制減緩幅度， t 為訓(xùn)練次數(shù)，本文中 k 設(shè)置為0.005。訓(xùn)練集數(shù)據(jù)為樣本數(shù)據(jù)的70%，驗證數(shù)據(jù)集為樣本數(shù)據(jù)的30%，由圖7 可知誤差下降速度較快，在訓(xùn)練次數(shù)為300 次時基本達(dá)到收斂，并且在下降過程中沒有出現(xiàn)劇烈的抖動，表明本文所提出的模型結(jié)構(gòu)具有收斂速度快和訓(xùn)練穩(wěn)定性高的優(yōu)點。

我們對比了本文所提出的模型與傳統(tǒng)的多元線性回歸(Regression)，多層感知器(Multi-layer Perceptron, MLP)以及遞歸神經(jīng)網(wǎng)絡(luò)(Recursive Neural Network, RNN)對葉綠素a 濃度的預(yù)測結(jié)果，圖8 為實驗對比分析圖，其中橫坐標(biāo)為時間點的個數(shù)，共30 個時間點，每個時間點相差5 min，縱坐標(biāo)為歸一化后的葉綠素a 濃度值。圖9 為實際葉綠素a 濃度值和預(yù)測葉綠素a 濃度值的真實誤差對比，從圖中可以看出這組實驗在預(yù)測誤差方面本文方法比以上所提的4 種方法低很多。

我們通過計算不同方法的平均誤差發(fā)現(xiàn)，本文所提出的模型在實驗中誤差結(jié)果別為0.11，取得了最好的結(jié)果。傳統(tǒng)的多元線性回歸的誤差結(jié)果為0.48，優(yōu)于多層感知器模型MLP 的誤差結(jié)果0.49，遞歸神經(jīng)網(wǎng)絡(luò)的誤差結(jié)果0.16，優(yōu)于傳統(tǒng)的多元線性回歸的平均誤差結(jié)果。綜合以上分析，我們發(fā)現(xiàn)考慮不同要素在不同時間點之間的內(nèi)部關(guān)系，對葉綠素a 濃度的預(yù)測結(jié)果準(zhǔn)確性可以得到大幅度的提升，其次由于遞歸神經(jīng)網(wǎng)絡(luò)將要素間的長短期依賴關(guān)系看做相同的整體，因此其預(yù)測誤差比本文所提出的模型大。同樣的，這也表明我們不能簡單的將多要素的信息看做一個整體，而是需要先對各要素間的長短期依賴關(guān)系進(jìn)行分析，再依據(jù)分析結(jié)果設(shè)計合理的網(wǎng)絡(luò)結(jié)構(gòu)。

表 1 不同時延下相關(guān)系數(shù)Table 1 Correlation coefficient under different time delays

圖 7 模型訓(xùn)練誤差下降情況Fig. 7 Decline of model training error

圖 8 實際葉綠素a 濃度值和預(yù)測葉綠素a 濃度值Fig. 8 Actual and predicted chlorophyll a concentration

6 結(jié)論

本文考慮到多要素與葉綠素a 濃度之間的關(guān)聯(lián)性，再加上多要素與葉綠素a 濃度的長短期依賴關(guān)系，并結(jié)合長短期記憶的人工神經(jīng)網(wǎng)絡(luò)預(yù)測葉綠素a 濃度取得了較高的預(yù)測精度。

（1）本文分析了與葉綠素a 濃度相關(guān)要素的長短期依賴關(guān)系，發(fā)現(xiàn)在不同的時間點，各要素間依據(jù)相關(guān)系數(shù)可以分為兩種依賴關(guān)系，從查詢文獻(xiàn)中發(fā)現(xiàn)，本文作為第一次提出對于原始要素的分析研究。

（2）本文依據(jù)要素間不同的依賴關(guān)系所得的分類結(jié)果，我們將原始輸入要素分為兩個不同的子網(wǎng)絡(luò)進(jìn)行訓(xùn)練，不同的子網(wǎng)絡(luò)所采用的神經(jīng)元和激活函數(shù)不同，最后在融合層進(jìn)行特征的融合，該網(wǎng)絡(luò)結(jié)構(gòu)具有一定創(chuàng)新性。本文所提出的網(wǎng)絡(luò)模型結(jié)構(gòu)與Regression、MLP、RNN 在葉綠素a 濃度的預(yù)測結(jié)果中，本文的模型誤差均有大幅度的下降（圖8，圖9）。

（3）本文所提出的模型具有收斂速度快，訓(xùn)練過程平穩(wěn)的特點。從圖7 中可以發(fā)現(xiàn)，本文所提模型結(jié)構(gòu)訓(xùn)練誤差下降較快，在300 步后，誤差趨于收斂狀態(tài)，其次我們發(fā)現(xiàn)誤差下降過程平穩(wěn)，沒有出現(xiàn)劇烈的抖動情況，這一現(xiàn)象說明了模型在收斂速度及訓(xùn)練過程中的優(yōu)點，另一方面也證明了本文所提模型的可行性。

圖 9 實驗結(jié)果圖Fig. 9 Experimental results