李修竹, 蘇榮國??, 張傳松, 石曉勇,2
(1.中國海洋大學化學化工學院,山東 青島 266100; 2.國家海洋局海洋減災中心, 北京 100194)
隨著沿海地區經濟的飛速發展,大量工業廢水和污染物排放入海,導致近岸海域富營養化嚴重,赤潮頻發。赤潮對生態環境和人類健康造成巨大破壞,打破了海域生態系統平衡,并給水產養殖造成了巨大損失,引起了國內外的廣泛關注[1-2]。因此,對水體富營養化進行早期預測很有必要。眾所周知,海洋生態系統十分復雜,各因素之間的相互作用及其動態變化過程未被完全知曉,制約了傳統的生態水動力模型的發展[3]。
近年來,數據智能技術在預測模型中興起,主要包括遺傳算法[4]、人工神經網絡[5]、多項式回歸[6]、灰色理論[7]等方法。許多學者運用這些方法對能反映近岸海域富營養化狀況的葉綠素a進行預測,并取得了一定成果。其中,BP神經網絡(BP-ANN)在葉綠素a濃度預測中運用最廣,實例最多。但是神經網絡的經驗風險最小化是基于樣本夠多的情況,在處理小樣本和多變量數據時,并不能很好的保證模型的泛化能力,存在過擬合現象,即樣本數較少時,預測結果不能保證。
支持向量機(SVM)是根據Vapnik[8]建立的統計學理論為基礎,遵循結構風險最小化原理建立的一種處理數據方法。該方法通過引用核函數,實現了在小樣本、非線性的條件下提高模型預測的準確性,避免了BP-ANN局部最優和收斂時長等問題,彌補了神經網絡的不足。
長江口海域生態環境復雜多樣[9],隨著人類活動的增加,長江口鄰近海域富營養化嚴重,赤潮頻發,漁業環境質量下降,嚴重影響了江蘇沿岸經濟的發展,因此對該海域進行環境預測技術研究很有必要。葉綠素a是反映海水中藻類生物量的重要指標,也是表征水體富營養化程度的重要指示之一[10]。預測葉綠素a濃度可以為長江口及鄰近海域富營養化的監測和管理提供支持。
本文以長江口鄰近海域作為研究對象,2015年春季(3月)和夏季(7月)在長江口鄰近海域(29°77′N~32°25′N,122°00′E~124°00′E),分別設置了9個斷面,60個站位進行現場調查,共采集了172個表層和次表層的海水樣品,具體站位分布圖見圖1。具體采樣方法:現場用Niskin采水器根據站位水深進行采水,采水后立即用0.7 μm的GF/F膜過濾,并用馬弗爐燒過的錫紙包好冷凍保存,帶回實驗室測定葉綠素a濃度;過濾后的水樣部分裝在灼燒后的棕色玻璃瓶(已預先在400 ℃的馬弗爐中灼燒4 h)中冷凍保存,用于測定有色溶解有機物(CDOM)的特征吸收系數;另一部分水樣冷凍保存于100 mL的聚乙烯瓶中,用于總氮、總磷的測定。

圖1 2015年3和7月采樣站位示意圖Fig.1 Map of sampling station in March and July 2015
溫度(T)、鹽度(S)和溶解氧(DO):利用多參數水質儀CTD獲得。
TN和TP:采用磷鉬藍法和Cr-Cu還原法獲得,即先用含硼酸的堿性過硫酸鉀氧化消解海水樣品,再用Cr-Cu還原法把氧化后的硝酸鹽轉化為亞硝酸鹽,用重氮-偶氮反應顯色測定,無機磷酸鹽則采用磷鉬藍法測定。
CDOM的特征吸收系數:采用紫外可見分光光度法測定,測樣儀器為島津UV-2550紫外可見分光光度計,并用5 cm石英比色皿,以Mill-Q水為參比測定過濾水樣。
葉綠素a:采用分光光度法獲得,先將冷凍的GF/F濾膜放入離心管,并用10 mL 90%的丙酮溶液振蕩萃取得到上清液,將上清液置于1 cm比色皿中,用島津UV-2550紫外可見分光光度計,以丙酮作參比,測定630、647、664、750 nm的其吸光度,再利用Jeffrey-Humphrey的方程計算Chl-a的濃度[11]:
Chl-a(mg·m-3)=[11.85×(A664-A750)-
1.54×(A647-A750)-0.08×(A630-A750)]×
Ve/(L×Vf)。
其中:L為比色皿長度;Ve為萃取液的體積;Vf為過濾水樣的體積;A為吸光度。
支持向量機(SVM)的基本思想就是利用核函數將低維輸入空間中線性不可分的點映射成高維特征空間中線性可分的點,并通過劃分超平面使同類樣本之間相似性盡可能的大,即所有的點到分類超平面的距離最大化,達到最大泛化能力[12-13]。
近海海域富營養化是一個多因素耦合、多維度協同作用的結果,關系復雜且具有多維的非線性特征[14-15],此外,葉綠素a濃度與各影響因素之間也存在典型的非線性關系,而支持向量機回歸(SVR)就是將實際問題通過非線性映射到高維特征空間,并在高維特征空間構建線性回歸,從而得到低維空間的非線性回歸效果[16]。支持向量機回歸(SVR)模型的樣本只有一類,所尋求的最優平面是使所有樣本點離超平面的“總偏差”最小,樣本點都在兩條邊界之間,相當于求最大間隔的問題[17]。
支持向量機回歸模型與人工神經網絡類似,主要包括三個部分:輸入層、中間層和輸出層,具體模型見圖2。

圖2 支持向量機模型的輸入與輸出Fig.2 The inputs and outputs of support vector machine
其中:x1~xn為低維輸入向量,即xi=(x1,x2, …,xn);K(xi,x) =〈φ(xi) ·φ(x)〉為核函數,xi為支持向量,φ向高維空間映射的函數[13],通過對應支持向量的拉格朗日系數(β1, β2, …, βn)連接輸入向量和核函數就能得到線性組合函數f(x)。
(1)
兩個低維空間中的輸入向量經過某種變化后計算出其在高維空間中的向量內積值即為核函數,巧妙的避免了向量由低維向高維空間映射時計算復雜的問題。目前支持向量機回歸模型常用的核函數有線性核函數、多項式核函數、Sigmoid核函數和徑向基核函數(RBF核函數)[18]。與多項式和Sigmoid核函數相比,徑向基核函數參數少,更適合非線性映射,且具有較寬的收斂域,計算更為簡單,因此,徑向基核函數應用更廣[18-19]。本文使用徑向基核函數,公式如下:
(2)
其中g為核參數。
核函數選定后,需要確定相應的最優懲罰參數c和核參數g,其中c主要表示懲罰系數,即對誤差的容忍度,g表示主要影響樣本數據子空間分布的復雜程度。最優懲罰參數c和核參數g通常通過K-折交叉驗證法(K-fold Cross Validation)得到[20]。具體操作是將數據集平均分為K組,輪流將其中的K-1組做訓練,剩余的1組做驗證,在給定參數情況下,K次的結果的均值作為對模型的評價指標,對支持向量機回歸模型表示為預測值和實際值的均方誤差(MSE)。
有色溶解有機物(CDOM)是水體中溶解有機物(DOM)的基本組分,能夠影響控制C、N、P等元素的生物地球化學循環和浮游植物進行光合作用[21]。CDOM的特征吸收系數與水體DOM含量及性質、水體濁度等密切相關,是近海海域生態環境監測主要指標的組成部分[22]。水溫(T)、鹽度(S)、總氮(TN)、總磷(TP)、溶解氧(DO)這些參數是評價海水水質的基本指標,對海洋中藻類產生直接或間接影響,進而影響水體中葉綠素a濃度[23]。因此,本文確定的監測參數包括水溫(T)、鹽度(S)、總氮(TN)、總磷(TP)、溶解氧(DO)等基本參數和有色溶解有機物(CDOM)特征吸收系數aCDOM(355)和aCDOM(455)。
對獲得的172個樣品的溶解氧、鹽度、溫度、TN、TP、aCDOM(355)、aCDOM(455)參數及測得的Chl-a濃度進行統計分析,由表1可知,Chl-a的平均值為1.754 9 μg·L-1,變化范圍在0.013 1~18.954 4 μg·L-1之間。其中,長江口附近海域表層和近岸站位所采集的水樣中Chl-a濃度較高。長江口近岸海域受陸源輸入影響較大,營養物質濃度較高,浮游植物生長旺盛[24]。

表1 各參數數值特征Table 1 The numerical characteristic of each parameter
此外,考慮到數據中各變量存在量綱和數量級的差異,對數據進行了歸一化處理以減小數值差異帶來的影響[25]。具體處理方法如下:
其中:x′為變量x歸一化后的值;xmax和xmin分別是數據的最大值和最小值。
對支持向量機回歸模型的輸入變量進行相關性分析,去除不相關或重復變量,對于保證模型的合理性以及提高模型的準確度具有重要意義。本文以葉綠素a濃度作為輸出變量,以水溫、鹽度、總氮(TN)、總磷(TP)、溶解氧和有色溶解有機物(CDOM)特征吸收系數aCDOM(355)、aCDOM(455)作為候選輸入變量,分別計算了各輸入變量與葉綠素a濃度的Pearson相關系數,結果如表2。
由表2可知,在0.01的顯著性水平下,溶解氧(DO)、溫度和吸收系數aCDOM(355)、aCDOM(455)與Chl-a都具有顯著相關關系,相關系數在-0.204~0.479之間,而鹽度和TN與Chl-a在0.05的顯著性水平下也呈顯著相關關系,其相關系性系數分別為-0.192和0.165,這表明本研究所選取的大多參數與Chl-a之間具有相關性,考慮到TP與其余參數在0.01的顯著性水平下有顯著相關關系,且TP是海水富營養化監測的基本要素,故最終的輸入變量確定為溫度、鹽度、總氮(TN)、總磷(TP)、溶解氧(DO)以及CDOM的紫外特征吸收系數aCDOM(355)和aCDOM(455)。

表2 輸入變量和葉綠素a濃度的Pearson相關系數Table 2 The Pearson correlation coefficient between the input variables and Chl-a
注:**表示P<0.01;*表示P< 0.05。**indicates very significant association;*indicates very significant association.
從172個樣品中隨機抽取112個樣品作為訓練集,剩余的60個樣品作驗證集。以CDOM特征吸收系數aCDOM(355)和aCDOM(455)以及溫度、鹽度、溶解氧等7個參數作為輸入變量,以葉綠素a濃度為因變量。
支持向量機回歸采用臺灣大學林智仁開發設計的LIBSVM-3.1工具包實現,在MATLAB平臺下進行建模[26],設置核函數為徑向基核函數,相應參數為最優懲罰參數c和核參數g。為確定c和g最佳參數值,將c和g分別取以2為底的指數離散值,代入 K-CV交叉驗證的算法中,選取這K個模型中平均驗證準確度最大,即平均驗證均方根誤差(MSE)最小的 那組c、g值作為該模型的參數,該方法被稱為“網格尋優法”(GS)。為減少計算量,把c、g的間隔設置大一點,再通過最佳參數位置范圍逐漸減小其范圍和間隔,進行精細的網格尋優,以此確定最終的參數值。
設置模型參數c∈ {2-10,2-9.5, …,210},g∈ {2-10,2-9.5, …,210},V=10,進行網格尋優搜索,結果如圖3所示。

圖3 網格尋優搜索結果Fig.3 The optimization results of Grid Search
利用網格尋優搜索方法得到的最佳參數值為:c=11.313 4、g=0.5,得到最佳參數后,核函數為,按該參數進行設置后,輸入訓練集數據,就能得到最終的葉綠素a濃度預測模型,該支持向量機回歸模型的函數可表達為:
模型的性能由可決系數(R2)和均方誤差(MSE)決定,R2表示測量值與預測值之間的相關性,R2越接近于1,表示樣本的預測值對實測樣本的擬合度越好,模型的擬合效果越好[27]。而MSE主要用來表征樣本數據之間的變化程度,MSE的數值越小,表明預測模型對實驗數據的分析具有越好的精確度[28],該模型的MSE=0.048 7。
將訓練集和驗證集數據輸入上述預測模型中,對輸出值進行反歸一化,得到模擬的葉綠素a濃度,并將實測葉綠素a濃度值和預測值進行對比,對比結果具體見圖4所示。
由圖4可知,以7個變量構建的GS-SVR模型所輸出的葉綠素a濃度預測值和實測值在變化趨勢上大致相同。在訓練集中,模擬值和實測值在0.01的顯著性水平下,Pearson相關系數為0.886(p<0.01),均方誤差MSE為0.024 0;而驗證集中,模擬值和實測值在0.01的顯著性水平下,Pearson相關系數為0.840(p<0.01),均方誤差MSE為0.041 8。Zhang等[29]利用基于主成分分析(PCA)方法的模糊BP神經網絡模型預測中國東海近岸海域葉綠素a濃度,預測結果與實測值具有良好一致性,MSE為0.109;Rocha等[30]通過多元線性回歸方法預測了巴西帕爾杜河葉綠素a濃度,其預測結果與實測值的Pearson相關系數為0.520;Zheng等[31]利用元胞自動機與支持向量機結合(CA-SVM)建立了渤海灣葉綠素a濃度預測模型,其預測結果與實測值的R2為0.861,均方差MSE為0.190。相較而言,本研究所建立的模型得到的預測值和實測值具有更好的一致性。

圖4 實際值與預測值對比Fig.4 Comparison of measured values and predicted values
根據7個輸入變量建立的支持向量機回歸預測模型,由w=∑s.v.βiφ(xi)求出各變量權重系數,并得到各輸入變量對輸出變量的重要性,將權重系數最大的溫度賦值100,可以得到各輸入參數的相對重要性,具體情況見表3。

表3 輸入變量的權重系數和相對重要性Table 3 Input variables in importance according to their weights and standardized weights
由上述表3中的權重系數可知,支持向量機回歸預測模型中對葉綠素a濃度預測影響顯著的輸入變量是溫度和CDOM 特征吸收系數aCDOM(355)。海水中的浮游植物通過光合作用進行初級生產,而葉綠素a是浮游植物進行光合作用的重要色素。有研究表明,溫度是一切酶促反應的控制因子,水溫與浮游植物的初級生產密切相關[32],浮游植物代謝率和光合作用暗反應都取決于水溫,當光照充足時,光合作用的速度與溫度呈正相關[33]。

除此以外,鹽度和溶氧對長江口鄰近海域葉綠素a濃度預測模型影響也較大。鹽度是反應近岸海域特別是海域陸源輸入的常用指標[41],在這些區域鹽度與營養鹽之間都有較為明顯的負相關關系,水系混合影響著營養鹽的消長[42]。溶解氧是海洋浮游植物光合作用的產物,也是海洋中影響異養生物活動的主要因素[43],是衡量海水水質的基本參數之一。
基于CDOM 特征吸收系數aCDOM(355)和aCDOM(455)以及溫度、鹽度、溶解氧、TP、TN等5個基本水質參數作為輸入變量,利用支持向量機回歸(SVR)建立了長江口鄰近海域葉綠素a濃度預測模型,預測值與實測值具有較好的一致性,且溫度和CDOM 特征吸收系數影響顯著,表明該模型能較好的預測長江口鄰近海域葉綠素a濃度,可為長江口及鄰近海域富營養化監測提供技術支持。