(三峽大學(xué) 水利與環(huán)境學(xué)院,湖北 宜昌 443002)

水華是在溫度、光照、營養(yǎng)鹽等環(huán)境條件適宜時，藻類大量生長繁殖并富集成一定濃度，導(dǎo)致區(qū)域水體變色的現(xiàn)象，是水體富營養(yǎng)化的典型代表。水華暴發(fā)易引起明顯的水質(zhì)變化，嚴(yán)重阻隔生態(tài)能量的有效循環(huán)，破壞水體的生物多樣性，已成為水環(huán)境治理中的難題[1]。葉綠素a含量常被作為特征指標(biāo)用以預(yù)測藻類生長暴發(fā)，在水體理化性質(zhì)和生物存量分析指標(biāo)中占據(jù)重要地位[2]。研究藻類水華的暴發(fā)機制，分析水體中葉綠素a的時空變化規(guī)律，預(yù)測水華的時空分布，對水華及其生態(tài)影響的預(yù)警和防范都具有重大意義[3]。

利用實時監(jiān)測數(shù)據(jù)，采用神經(jīng)網(wǎng)絡(luò)對水體中葉綠素a含量進(jìn)行預(yù)測，是水體中葉綠素a含量預(yù)測的主要手段之一[4]，得到了廣泛應(yīng)用[5]，取得了比較有效的預(yù)測結(jié)果。例如Velo-Suarez等[6]在Andalucia的大西洋沿岸流域?qū)崿F(xiàn)了對藻類的神經(jīng)網(wǎng)絡(luò)預(yù)測模型構(gòu)建；Guallar等[7]在地中海Alfacs灣建立了對雙鞭毛藻和硅藻的神經(jīng)網(wǎng)絡(luò)預(yù)測機制；趙文喜等[8]在中國海河干流完成了基于BP神經(jīng)網(wǎng)絡(luò)的葉綠素a含量預(yù)測短時預(yù)測研究，皆對相應(yīng)流域的藻類生長趨勢進(jìn)行了有效預(yù)測。

輸入數(shù)據(jù)的準(zhǔn)確性對神經(jīng)網(wǎng)絡(luò)的預(yù)測精度有決定性的影響，但實時監(jiān)測數(shù)據(jù)受各種隨機因素影響，會存在異值點、數(shù)據(jù)缺失、數(shù)據(jù)不光滑等不符合物理規(guī)律的情況[9]，需要對原始監(jiān)測數(shù)據(jù)處理，以提高數(shù)據(jù)質(zhì)量及預(yù)測精度。相關(guān)研究表明，不同處理方法對預(yù)測精度影響較大。例如王亞宸等[10]在澳大利亞能源市場研究中，采用小波變換處理噪音數(shù)據(jù)，成功實現(xiàn)了維多利亞州電力負(fù)荷以及電價的高精度預(yù)測；Iliou等[11]在骨質(zhì)疏松癥預(yù)測案例中，應(yīng)用了一種基于MLP分類器的新型數(shù)據(jù)處理方法，達(dá)到神經(jīng)網(wǎng)絡(luò)高效分類的目的。但數(shù)據(jù)處理對基于神經(jīng)網(wǎng)絡(luò)的水體葉綠素a含量預(yù)測精度影響方面尚沒有專門研究。

因此，本文在實時監(jiān)測數(shù)據(jù)基礎(chǔ)上，分別采用3種異值點處理方法與2種數(shù)據(jù)光滑方法組合，討論不同數(shù)據(jù)處理方法對主成分影響，研究神經(jīng)網(wǎng)絡(luò)的輸入?yún)?shù)確定方法，選擇5種神經(jīng)網(wǎng)絡(luò)輸入?yún)?shù)組合，分析不同數(shù)據(jù)處理方案下神經(jīng)網(wǎng)絡(luò)預(yù)測精度，研究樣本數(shù)據(jù)處理策略對基于神經(jīng)網(wǎng)絡(luò)的葉綠素a預(yù)測精度的影響，為提高基于神經(jīng)網(wǎng)絡(luò)的水華暴發(fā)預(yù)警技術(shù)的預(yù)警精度提供支撐。

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

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

研究所采用的數(shù)據(jù)為研究水域2016年4月下旬到12月監(jiān)測時段內(nèi)的pH、氨氮(NH3-N)、電導(dǎo)率(COND)、水溫(WT)、溶解氧(DO)、葉綠素a(Chl-a)、淡水藍(lán)綠藻(AFA)、氧化還原電位(ORP)、氣溫(AT)、氣壓(AP)、相對濕度(RH)、降雨(Rainfall)及光強(Lux)13項實時監(jiān)測數(shù)據(jù)，數(shù)據(jù)頻率為每10 min1次。受多種因素影響，數(shù)據(jù)存在異常值，部分?jǐn)?shù)據(jù)缺失，數(shù)據(jù)連續(xù)性及光滑性不足等問題(圖1(a)中，氨氮監(jiān)測數(shù)據(jù)存在異常值；圖1(b)中，水溫監(jiān)測值連續(xù)性、光滑性不足)。

圖1 監(jiān)測數(shù)據(jù)異常分布示意Fig.1 Unusual distribution of monitoring data

1.2 異常值及缺失值處理

結(jié)合本次監(jiān)測數(shù)據(jù)的采集頻率及實際分布狀態(tài)，本文依次嘗試使用以下3種判據(jù)準(zhǔn)則進(jìn)行異常值處理。

(1) 拉依達(dá)準(zhǔn)則[12](Pauta Criterion)，又稱3σ準(zhǔn)則，假設(shè)監(jiān)測數(shù)據(jù)只存在隨機誤差，對數(shù)據(jù)計算出標(biāo)準(zhǔn)差σ，隨機誤差在指定概率區(qū)間(-3σ,3σ)的分布概率約為99.7%，監(jiān)測誤差超過指定概率區(qū)間就判定為異常值。

(2) 肖維勒準(zhǔn)則[13](Chauvenet Criterion)，在監(jiān)測數(shù)據(jù)中，臨近時段的n次監(jiān)測數(shù)據(jù)，如果某監(jiān)測值xi與平均值x之差的絕對值大于標(biāo)準(zhǔn)偏差與肖維勒系數(shù)之積，則該監(jiān)測數(shù)據(jù)為異常值。

(1)

ωn=1+0.4ln(n)

(2)

針對藻類水華數(shù)據(jù)變幅較大的特性，肖維勒法剔除異常值時樣本數(shù)需靈活取值；本文以樣本數(shù)500為一組，肖維勒系數(shù)ω500=3.20進(jìn)行異常值判斷。

(3) 格拉布斯準(zhǔn)則[14](Grubbs Criterion),監(jiān)測值對應(yīng)殘差的絕對值滿足下式時，判定監(jiān)測值為異常值。

(3)

式中，g(n,a)為格拉布斯臨界系數(shù)，與監(jiān)測樣本數(shù)以及顯著水平有關(guān)，本文臨界系數(shù)取g(100,0.05)=3.17，以樣本數(shù)100為一組進(jìn)行分批檢驗。

本文基于以上3種準(zhǔn)則對監(jiān)測數(shù)據(jù)進(jìn)行異常值篩選處理。對于同一監(jiān)測數(shù)據(jù)，3種準(zhǔn)則選取評判的樣本群體不同(σ=0.398 2，ω500=3.20，g(100,0.05)=3.17)，對應(yīng)的判定結(jié)果亦相異，拉依達(dá)準(zhǔn)則判定標(biāo)準(zhǔn)最為寬松，格拉布斯準(zhǔn)則最為嚴(yán)格，肖維勒準(zhǔn)則居于兩者之間。此外，若該點監(jiān)測數(shù)據(jù)為異常值，則當(dāng)作缺失值進(jìn)行處理。

常用缺失值處理方法有個案剔除法、均值替換法、熱卡填充法、回歸替換法、期望最大化法等。本文采用如下方法處理：同一時間點若單項監(jiān)測數(shù)據(jù)缺失則采用均值替換法，若同時多項監(jiān)測數(shù)據(jù)缺失則采用個案剔除法去除該時間點所有監(jiān)測數(shù)據(jù)。

1.3 數(shù)據(jù)平滑濾波

水體的理化監(jiān)測指標(biāo)理論上應(yīng)是光滑連續(xù)變化的，但受監(jiān)測頻次限制及外界隨機因素干擾，數(shù)據(jù)光滑性常常不滿足要求，需要利用平滑濾波方法處理。本文在異常值處理的基礎(chǔ)上，分別采用臨近加權(quán)平均法、局部多項式回歸法[15]對監(jiān)測數(shù)據(jù)平滑。

臨近加權(quán)平均法中，以監(jiān)測點i為計算中心，計算臨近個點k的加權(quán)平均值作為監(jiān)測數(shù)據(jù)xi對應(yīng)的平滑值：

(4)

圖2 監(jiān)測數(shù)據(jù)平滑效果對比Fig.2 Smoothing effect comparison of monitoring data

局部多項式回歸，是基于最小二乘法原理在該點擬合回歸的多項式代入值與監(jiān)測值之差的平方和最小時，確定局部多項式回歸效果最佳，平滑值計算公式如下：

j=-m,…,0,…,+m

(5)

本文采用納什效率系數(shù)分析數(shù)據(jù)平滑效果：

(6)

對以拉依達(dá)準(zhǔn)則、肖維勒準(zhǔn)則和格拉布斯準(zhǔn)則進(jìn)行數(shù)據(jù)異值處理后得到的3組數(shù)據(jù)，分別進(jìn)行臨近加權(quán)平均、局部多項式回歸處理得到6組數(shù)據(jù)變換預(yù)處理方案，方案中各監(jiān)測指標(biāo)以納什效率系數(shù)不低于0.985作為平滑標(biāo)準(zhǔn)。各方案既有效消除了數(shù)據(jù)中的高頻“噪音”影響，同時也盡可能地保持了數(shù)據(jù)原有的客觀真實性，可以為后續(xù)研究提供數(shù)據(jù)支撐。部分?jǐn)?shù)據(jù)平滑前后的對比如圖2所示。

1.4 主成分分析

表1 不同數(shù)據(jù)處理方案主成分分析結(jié)果Tab.1 Principal component analysis results of different data processing schemes

注：1.“3σ”，拉依達(dá)準(zhǔn)則；“Chauvenet”，肖維勒準(zhǔn)則；“Grubbs”,格拉布斯準(zhǔn)則；2.“AAv”，臨近加權(quán)平均法；“SG”，局部多項式回歸；3.基于特征值0.6；具有Kaiser標(biāo)準(zhǔn)化的正交旋轉(zhuǎn)法。

以pH、氨氮、電導(dǎo)率、水溫、溶解氧、葉綠素a、氧化還原電位、氣溫、氣壓、相對濕度、降雨及光強標(biāo)準(zhǔn)化后數(shù)據(jù)進(jìn)行主成分分析，分析結(jié)果擬合度均達(dá)到85%以上。顯然，不同數(shù)據(jù)處理方案得到的主成分也有區(qū)別，說明數(shù)據(jù)處理方案對主成分分析結(jié)果有一定影響。但總體而言占據(jù)優(yōu)勢的成分依次為：葉綠素a、氣溫、光強、氣壓、降雨、電導(dǎo)率、相對濕度；可將其作為葉綠素a含量預(yù)測輸入?yún)?shù)。

2 基于BP神經(jīng)網(wǎng)絡(luò)葉綠素a預(yù)測模型構(gòu)建

2.1 輸入輸出參數(shù)的選擇

由于葉綠素a的含量與藻類的數(shù)量密切相關(guān)，在一定程度上能夠反映水質(zhì)狀況，是判斷水體富營養(yǎng)化的重要指標(biāo)之一[18]，前文主成分分析表明，葉綠素a需為模型輸入?yún)?shù)，才可保證模型的預(yù)報精度。需要定義新的輸出參數(shù)，以預(yù)測下一時刻葉綠素a含量。鑒于此，定義單位時間內(nèi)葉綠素a含量的變化量為平均生長率GR，用公式表示為

(7)

其中，GR為t1～t2時刻之間的平均生長率。

下一時刻葉綠素a含量，由下式計算：

(8)

2.2 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)設(shè)置及訓(xùn)練

適合的BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)設(shè)置既能加快收斂速度，也能保證預(yù)測精度[19]。本文依據(jù)預(yù)測誤差最小確定葉綠素a含量預(yù)測的BP神經(jīng)網(wǎng)絡(luò)隱含層為2層，隱含層神經(jīng)元個數(shù)為12；隱含層采用learngdm閾值學(xué)習(xí)函數(shù)，輸出層采用purelin傳遞函數(shù)；最大迭代次數(shù)1 000次，設(shè)置0.05的步長。有效數(shù)據(jù)11 000組，其中8 000組作為訓(xùn)練集，2 000組作為測試集，1 000組作為預(yù)測集。葉綠素a含量預(yù)測的BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)圖如圖3所示。

圖3 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig.3 BP neural network structure

研究不同數(shù)量輸入?yún)?shù)組合對神經(jīng)網(wǎng)絡(luò)預(yù)測精度的影響，當(dāng)參數(shù)輸入數(shù)量相同時，以預(yù)測精度最好的組合作為該輸入?yún)?shù)下的最優(yōu)組合。以格拉布斯準(zhǔn)則異常值處理組合局部多項式回歸平滑處理數(shù)據(jù)為例，不同輸入?yún)?shù)下葉綠素a含量預(yù)測值分布對比見圖4。不同數(shù)量的最優(yōu)輸入?yún)?shù)組合，對葉綠素a含量均有較好的預(yù)測精度，說明影響葉綠素a含量變化的因子間不是相互獨立的，各因子間存在交織作用。因此，關(guān)注主要影響因子，精簡預(yù)測模型輸入是提高葉綠素a含量預(yù)測效率的可行途徑。

不同數(shù)量的最優(yōu)參數(shù)組合預(yù)測精度如表2所示，對比主成分分析結(jié)果可發(fā)現(xiàn)，7參數(shù)下輸入?yún)?shù)的最優(yōu)組合與主成分分析結(jié)果相同，當(dāng)輸入?yún)?shù)減少時，去除的輸入?yún)?shù)為主成分分析結(jié)果中權(quán)重最小的因子，主成分分析結(jié)果可用于指導(dǎo)神經(jīng)網(wǎng)絡(luò)輸入?yún)?shù)選擇；神經(jīng)網(wǎng)絡(luò)預(yù)測精度隨著輸入?yún)?shù)的增加而提高，將所有監(jiān)測參數(shù)(13個)當(dāng)作輸入?yún)?shù)時，預(yù)測精度最高，可達(dá)0.994，相比4參數(shù)輸入的0.986，預(yù)測精度有所提高，但提高程度有限。其中，4、13參數(shù)葉綠素預(yù)測效果對比如圖5所示。

圖4 不同參數(shù)輸入下葉綠素a含量分布Fig.4 Distribution of Chlorophyll-a content under different parameters input

表2 不同輸入?yún)?shù)組合的平均生長率預(yù)測精度Tab.2 Average growth rate prediction accuracy of different input parameter combinations

3 數(shù)據(jù)處理對預(yù)測精度影響

將3種異值點處理方法與2種數(shù)據(jù)光滑方法組合，結(jié)合原始數(shù)據(jù)，得到了7組數(shù)據(jù)處理方案，與5種輸入?yún)?shù)數(shù)量方案組合，共可得到35種組合方案。針對每種數(shù)據(jù)處理方案，進(jìn)行主成分分析，根據(jù)各影響因子權(quán)值，由小到大剔除輸入?yún)?shù)，對BP神經(jīng)網(wǎng)絡(luò)訓(xùn)練，得到不同數(shù)據(jù)處理方案及不同輸入?yún)?shù)組合下的預(yù)測精度結(jié)果，如表3所示。

采用原始數(shù)據(jù)進(jìn)行預(yù)測時，預(yù)測精度普遍不高，如圖6所示，預(yù)測精度不隨輸入?yún)?shù)增加而增加，說明原始數(shù)據(jù)存在干擾項，影響神經(jīng)網(wǎng)絡(luò)預(yù)測精度。

對比原始數(shù)據(jù)方案與數(shù)據(jù)處理后方案，在輸入?yún)?shù)數(shù)目相同情況下，對數(shù)據(jù)進(jìn)行處理可顯著提高神經(jīng)網(wǎng)絡(luò)的預(yù)測精度，但不同的數(shù)據(jù)處理方法對其預(yù)測精度的提高程度不同。在采用臨近加權(quán)平均法對數(shù)據(jù)進(jìn)行平滑處理的條件下，不同數(shù)據(jù)異值處理方法對預(yù)測精度提高效果明顯，但預(yù)測精度與輸入?yún)?shù)數(shù)量間相關(guān)性較差，不同數(shù)據(jù)異值處理方法優(yōu)劣無法確定，但以拉依達(dá)準(zhǔn)則處理后得到的預(yù)測精度最高，為0.938。在采用局部多項式回歸方法對數(shù)據(jù)進(jìn)行平滑處理的條件下，不同數(shù)據(jù)異值處理方法對預(yù)測精度提高顯著，預(yù)測精度與輸入?yún)?shù)數(shù)量間相關(guān)性強，以格拉布斯準(zhǔn)則處理效果最好，肖維勒準(zhǔn)則處理效果較差。

圖5 4，13參數(shù)輸入的葉綠素a含量預(yù)測效果Fig.5 Comparison of Chlorophyll-a content prediction results with 4 and 13 input parameters

表3 不同數(shù)據(jù)處理方案的預(yù)測精度Tab.3 Prediction accuracy of different data processing schemes

注：不同數(shù)據(jù)處理方案得到的主影響因素不同，相應(yīng)不同方案的最優(yōu)輸入?yún)?shù)組合不同。

整體而言，采用局部多項式回歸方法對數(shù)據(jù)進(jìn)行平滑處理，對預(yù)測精度改善優(yōu)于臨近加權(quán)平均方法，如圖7、8所示，采用局部多項式回歸的葉綠素a含量預(yù)測誤差，整體波動幅度更小。格拉布斯準(zhǔn)則異值處理組合局部多項式回歸法平滑數(shù)據(jù)，在不同輸入?yún)?shù)數(shù)量下均可達(dá)到最佳預(yù)測結(jié)果。該方案可以將以葉綠素a、氣溫、光強及氣壓4項因素下的葉綠素a預(yù)測精度從原始數(shù)據(jù)的0.800提高至0.986，同比提高23.25%。

圖6 原始數(shù)據(jù)的不同輸入?yún)?shù)預(yù)測效果Fig.6 Different input parameters′ prediction effects with initial data

圖7 臨近加權(quán)平均法的葉綠素a含量預(yù)測效果Fig.7 Prediction effect of chlorophyll-a content by weighted average method

圖8 局部多項式回歸法的葉綠素a含量預(yù)測效果Fig.8 Prediction effect of chlorophyll-a content by local polynomial regression method

4 結(jié) 論

本文對神經(jīng)網(wǎng)絡(luò)的訓(xùn)練數(shù)據(jù)進(jìn)行預(yù)處理，建立了35種組合方案，對基于BP神經(jīng)網(wǎng)絡(luò)葉綠素a含量平均生長率進(jìn)行預(yù)測，對比分析預(yù)測結(jié)果，評估數(shù)據(jù)處理對基于神經(jīng)網(wǎng)絡(luò)的葉綠素a預(yù)測精度的影響。得到如下結(jié)論。

(1) 主成分分析方法可為BP神經(jīng)網(wǎng)絡(luò)輸入?yún)?shù)選擇和簡化提供極為重要的參考。

(2) 利用數(shù)據(jù)處理技術(shù)對基礎(chǔ)監(jiān)測數(shù)據(jù)進(jìn)行處理，可顯著提高基于神經(jīng)網(wǎng)絡(luò)的葉綠素a含量預(yù)測精度。

(3) 不同的數(shù)據(jù)處理方案對基于神經(jīng)網(wǎng)絡(luò)的葉綠素a預(yù)測精度影響幅度不同；以格拉布斯準(zhǔn)則進(jìn)行異值處理，組合局部多項式回歸法進(jìn)行數(shù)據(jù)平滑，是本研究的最佳數(shù)據(jù)處理方案。

(4) 本研究方法的內(nèi)在作用機理，需要在對數(shù)據(jù)處理前后幅值、頻譜等數(shù)據(jù)特征變化進(jìn)行深入研究后明確。