馬豐魁，姜群鷗,2*，徐藜丹，梁勇，王榮臣，蘇帥

1.北京林業(yè)大學(xué)水土保持學(xué)院，100083 北京；2.北京林業(yè)大學(xué)/水土保持與荒漠化防治教育部重點(diǎn)實(shí)驗(yàn)室，100083 北京；3.北京市密云水庫管理處，101512 北京

水是人類生存發(fā)展必需的自然資源，是社會(huì)經(jīng)濟(jì)發(fā)展的重要物質(zhì)保障，是生態(tài)環(huán)境的控制性因素。然而隨著中國經(jīng)濟(jì)社會(huì)的發(fā)展，水資源需求量日益增加，水污染整體呈現(xiàn)惡化態(tài)勢。水資源短缺和水污染已經(jīng)成為制約中國可持續(xù)發(fā)展的瓶頸因素（曾昭等，2013）。北京是一座重度資源型缺水、水污染嚴(yán)重的城市，經(jīng)濟(jì)快速發(fā)展，城市擴(kuò)張以及人口持續(xù)增長導(dǎo)致北京市的水質(zhì)不斷惡化（高陽等，2018）。在這種背景下，密云水庫作為北京市最大的飲用水源供應(yīng)地，其水質(zhì)安全不僅關(guān)乎北京市兩千多萬市民的生命健康和生活用水問題，而且對(duì)于北京市發(fā)展和首都各項(xiàng)職能的正常運(yùn)轉(zhuǎn)具有重大影響，引起了國家和政府的高度重視。保障城市水源安全、提高城鎮(zhèn)飲用水源的安全保障水平，是城市生態(tài)文明建設(shè)的重要內(nèi)容，也是生態(tài)文明建設(shè)的重要組成和基礎(chǔ)保障（馬巍等，2016）。因而探究監(jiān)測密云水庫主要水質(zhì)參數(shù)變化的有效方法，實(shí)現(xiàn)對(duì)密云水庫主要水質(zhì)指標(biāo)的快速反演尤為重要。

傳統(tǒng)的水質(zhì)監(jiān)測一般以實(shí)地采樣分析為主，可以較為準(zhǔn)確地測得多種水質(zhì)參數(shù)，但是需要耗費(fèi)大量的時(shí)間和成本，且只能獲得點(diǎn)狀數(shù)據(jù)，所得數(shù)據(jù)在時(shí)間和空間上都不連續(xù)，無法滿足大范圍、實(shí)時(shí)的水質(zhì)監(jiān)測要求。（趙玉芹等，2009）。水質(zhì)遙感技術(shù)作為一種具有區(qū)域化監(jiān)測能力的方法，可以快速、實(shí)時(shí)的獲取水質(zhì)的時(shí)空分布，在水質(zhì)監(jiān)測中日益受到重視（曹引等，2019）。常規(guī)的遙感水質(zhì)參數(shù)反演方法主要有經(jīng)驗(yàn)方法、半經(jīng)驗(yàn)法和半分析法（張海威等，2017）。翟召坤等（2018）基于高分一號(hào)衛(wèi)星數(shù)據(jù)和實(shí)測數(shù)據(jù)，采用典型的經(jīng)驗(yàn)方法對(duì)潘家口水庫水質(zhì)參數(shù)和敏感波段進(jìn)行了線性擬合，得到水質(zhì)參數(shù)遙感反演模型。林劍遠(yuǎn)等（2019）基于航空和水表高光譜遙感數(shù)據(jù)，利用半經(jīng)驗(yàn)法對(duì)浙江省嘉興市主城區(qū)河網(wǎng)水質(zhì)進(jìn)行了監(jiān)測。石亮亮（2019）基于實(shí)測數(shù)據(jù)，開發(fā)了水體IOPs及CDOM反演的半分析算法，對(duì)東海、千島湖水域進(jìn)行了反演。這些方法本質(zhì)上都是通過線性回歸來實(shí)現(xiàn)水質(zhì)反演，而水環(huán)境系統(tǒng)具有強(qiáng)烈的非線性和不確定性特征，傳統(tǒng)的線性回歸不能完全反映其變化規(guī)律，通過線性回歸分析會(huì)導(dǎo)致結(jié)果不準(zhǔn)確（Chen et al.，2018）。

近年來，人工神經(jīng)網(wǎng)絡(luò)理論的迅速發(fā)展為水質(zhì)評(píng)價(jià)和預(yù)測提供了一種新的有效方法。人工神經(jīng)網(wǎng)絡(luò)具有獨(dú)特的分布并行處理、非線性映射和自適應(yīng)學(xué)習(xí)等優(yōu)點(diǎn)，BP神經(jīng)網(wǎng)絡(luò)是其中一種較為成熟的非線性函數(shù)逼近的方法，所得結(jié)果更為準(zhǔn)確（李輝東等，2015；黃俊等，2015；易湘生等，2012）。國內(nèi)外許多專家學(xué)者運(yùn)用神經(jīng)網(wǎng)絡(luò)進(jìn)行了水質(zhì)參數(shù)的研究。Gazzaz et al.（2012）通過建立神經(jīng)網(wǎng)絡(luò)模型完成了對(duì)馬來西亞金塔河的水質(zhì)參數(shù)的預(yù)測，研究結(jié)果表明，該模型預(yù)測結(jié)果可以解釋95.4%左右的實(shí)測水質(zhì)參數(shù)的變化。Wu et al.（2009）基于MODIS數(shù)據(jù)集，采用BP神經(jīng)網(wǎng)絡(luò)、多元線性回歸等方法對(duì)巢湖的水質(zhì)參數(shù)進(jìn)行了反演，結(jié)果表明，BP神經(jīng)網(wǎng)絡(luò)模型優(yōu)于多元線性回歸，相對(duì)誤差低于35%。Song et al.（2011）使用Landsat TM數(shù)據(jù)和實(shí)測數(shù)據(jù)，建立了經(jīng)驗(yàn)回歸和神經(jīng)網(wǎng)絡(luò)，對(duì)查干湖的水質(zhì)參數(shù)進(jìn)行反演，結(jié)果表明神經(jīng)網(wǎng)絡(luò)的均方根誤差小于10%，而回歸的均方根誤差小于25%。Xiao et al.（2017）以漢江中下游為研究區(qū)，建立了VIP-BP模型估計(jì)葉綠素含量，結(jié)果表明，VIP-BP模型可以準(zhǔn)確真實(shí)的反映葉綠素含量變化，且誤差明顯低于3波段模型。

上述研究表明，神經(jīng)網(wǎng)絡(luò)與遙感技術(shù)的結(jié)合可以很好地實(shí)現(xiàn)水質(zhì)參數(shù)的反演。在此背景下，本研究基于遙感和GIS技術(shù)，采用BP神經(jīng)網(wǎng)絡(luò)算法，結(jié)合地面監(jiān)測數(shù)據(jù)和Landsat 8遙感影像，分別建立了總磷、總氮、氨氮、COD（化學(xué)需氧量）4個(gè)主要水質(zhì)參數(shù)的神經(jīng)網(wǎng)絡(luò)模型，實(shí)現(xiàn)了密云水庫主要水質(zhì)參數(shù)的反演。旨在為密云水庫水質(zhì)監(jiān)測探索一種便捷、可靠、高效的方法，從而對(duì)密云水庫主要水質(zhì)參數(shù)變化做出科學(xué)評(píng)估，為密云水庫水質(zhì)監(jiān)測和保護(hù)工作提出重要的科學(xué)建議。

1 材料與方法

1.1 研究區(qū)概況

圖1 密云水庫區(qū)位圖Fig.1 Location of Miyun Reservoir Area

密云水庫位于北京市東北部密云縣城北13 km處（116°07′—117°30′E，40°14′—41°05′N），地處燕山群山之中（圖1）。庫區(qū)建成于1960年9月，流域面積15788 km2，最大水面面積為188 km2，最大庫容43.75 ×108m3，最大水深43.5 m，是華北地區(qū)最大的水庫。密云水庫有2大入庫河流，分別是白河和潮河，分為白河、潮河和內(nèi)湖3個(gè)庫區(qū)（張敏等，2019）。該區(qū)屬于暖溫帶半濕潤季風(fēng)氣候，年均溫 10.5 ℃，年均最低氣溫-18 ℃，最高氣溫38 ℃，降水主要集中在6—8月（李東青等，2015）。庫區(qū)主要功能是防洪、灌溉、發(fā)電等，為北京市、天津市和河北省服務(wù)。隨著北京城市化和經(jīng)濟(jì)的快速發(fā)展，水資源短缺的局面日益加重，密云水庫于1981年開始專為北京市供水，且生活用水比重逐年上升，當(dāng)前已成為北京市最大的飲用水源供應(yīng)地，被譽(yù)為“北京生命之水”。自20世紀(jì)90年代以來，水體的富營養(yǎng)化已成為密云水庫最主要的水質(zhì)污染問題（許爾琪等，2018）。因此開展對(duì)密云水庫的水質(zhì)參數(shù)的反演研究十分必要。

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

1.2.1 地面監(jiān)測數(shù)據(jù)

研究團(tuán)隊(duì)于2018年10月17日乘船進(jìn)入密云水庫，在規(guī)劃采樣區(qū)域利用取水器距離水面30 cm深處共采集了80個(gè)樣品，同時(shí)記錄了采樣點(diǎn)的GPS坐標(biāo)。采樣點(diǎn)的布設(shè)原則是盡可能均勻分布，敏感區(qū)域多布設(shè)一些取樣點(diǎn)。采集回樣品之后，當(dāng)天將樣品送檢相關(guān)部門，在實(shí)驗(yàn)室對(duì)每個(gè)樣品進(jìn)行測試：水質(zhì)總氮采用化學(xué)發(fā)光器測量基態(tài)NO2，將結(jié)果換算為總氮；水質(zhì)氨態(tài)氮采用SEALAA3連續(xù)流動(dòng)分析儀測試；水質(zhì)總磷采用鉬酸銨分光光度法測定；水質(zhì)化學(xué)需氧量采用快速消解分光光度法測定。（陳穎等，2018；廖蕾等，2016；GB/T 11893—1989；HJ/T 399—2007）最終得到80采樣點(diǎn)的總磷、氨氮、總磷、COD（化學(xué)需氧量）4個(gè)指標(biāo)的質(zhì)量濃度。

1.2.2 Landsat 8遙感影像

本研究采用的遙感數(shù)據(jù)為Landsat 8，遙感影像下載自中國科學(xué)院遙感與數(shù)字地球研究所（http://eds.ceode.ac.cn/）和地理空間數(shù)據(jù)云（http://www.gscloud.cn/sources），空間分辨率為 30 m×30 m。選取研究區(qū)的影像時(shí)間跨度為 2013—2018年，其中用于構(gòu)建神經(jīng)網(wǎng)絡(luò)模型的遙感影像時(shí)間為2018年10月17日，與水質(zhì)樣品采集時(shí)間一致。由于密云水庫從12月至次年2月有冰期，該時(shí)段遙感影像的反射率無法真實(shí)反映水質(zhì)參數(shù)的變化，再加上一些時(shí)段的遙感影像受到云層的影響，最終選取2013—2018年非冰期共36期云量低于30%的遙感影像，代表非冰期主要水質(zhì)指標(biāo)的變化。原始影像含5個(gè)數(shù)據(jù)集：多光譜數(shù)據(jù)集（1—7波段）、全色波段數(shù)據(jù)（8波段）、卷云波段數(shù)據(jù)（9波段）、熱紅外數(shù)據(jù)（10、11波段）和質(zhì)量波段數(shù)據(jù)（12波段）以及含有整個(gè)圖像內(nèi)容的MTL文件，共13個(gè)文件。本研究選取多光譜數(shù)據(jù)（1—7波段）作為反演密云水庫水質(zhì)的原始數(shù)據(jù)。

1.3 數(shù)據(jù)處理

1.3.1 Landsat 8遙感影像預(yù)處理

本研究首先對(duì)Landsat 8影像進(jìn)行輻射定標(biāo)和大氣校正，使遙感影像表達(dá)真實(shí)的反射率（秦雁等，2011）。處理完成后將Landsat 8多光譜數(shù)據(jù)的7個(gè)波段分別導(dǎo)出為tif圖像，便于后期建模和預(yù)測。再通過歸一化水體指數(shù)NDWI，對(duì)研究區(qū)水域邊界進(jìn)行提取。然后使用提取出的密云水庫的水域邊界對(duì)遙感影像進(jìn)行裁剪，得到研究區(qū)各波段的遙感影像。最后將遙感影像進(jìn)行轉(zhuǎn)換，得到研究區(qū)全區(qū)的反射率值。

1.3.2 采樣點(diǎn)反射率提取

本研究將記錄的采樣點(diǎn)的 GPS坐標(biāo)離散到空間中，并將采樣點(diǎn)與Landsat 8影像疊加，用ArcGIS提取Landsat 8各波段在采樣點(diǎn)的反射率值，將提取的采樣點(diǎn)反射率數(shù)據(jù)進(jìn)行整理，作為建立神經(jīng)網(wǎng)絡(luò)模型的輸入數(shù)據(jù)。

1.3.3 敏感因子提取

水體的光譜特征由其中的各種物質(zhì)對(duì)光輻射的吸收和散射性質(zhì)所決定，當(dāng)水體受到污染時(shí)，光譜特征會(huì)發(fā)生變化，所以遙感影像各波段與水質(zhì)參數(shù)之間存在一定的相關(guān)性。使用 SPSS軟件，對(duì)Landsat 8的1—7波段的反射率和4個(gè)水質(zhì)參數(shù)（總氮、氨氮、總磷、COD）的實(shí)測值進(jìn)行皮爾遜相關(guān)性分析，以相關(guān)性分析為輔助，在保留最大信息量的原則下剔除干擾波段，得到水質(zhì)參數(shù)的敏感波段組合（表1），使得BP神經(jīng)網(wǎng)絡(luò)的擬合效果最優(yōu)。

1.4 研究方法

1.4.1 BP神經(jīng)網(wǎng)絡(luò)

BP（Back Propagation）神經(jīng)網(wǎng)絡(luò)是按誤差反向傳播算法訓(xùn)練的多層前饋式網(wǎng)絡(luò)。由Rumelhart和Mccelland提出，其運(yùn)算原理為使用最小二乘法使網(wǎng)絡(luò)的誤差平方和達(dá)到最小。網(wǎng)絡(luò)在訓(xùn)練過程中利用錯(cuò)誤反饋傳播不斷調(diào)整閾值和權(quán)值，使輸出結(jié)果達(dá)到最佳值，實(shí)現(xiàn)最佳的擬合效果（盛夏等，2013）。典型的BP神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)包括輸入層、隱含層和輸出層，每一層包含若干神經(jīng)元，上下層神經(jīng)元之間全連接，而同層神經(jīng)元之間無連接，輸入層與隱含層神經(jīng)元的連接值是網(wǎng)絡(luò)的權(quán)值，學(xué)習(xí)過程的目標(biāo)是找到最優(yōu)的權(quán)值集，該權(quán)值集可以產(chǎn)生正確的輸出（李麗等，2008；張斌等，2012）。模型的基本建立過程是：給定神經(jīng)網(wǎng)絡(luò)輸入和期望值，由輸入層傳遞到隱含層，經(jīng)隱含層處理后，再傳遞到輸出層，由輸出層處理產(chǎn)生結(jié)果，此為信息前向傳播；計(jì)算實(shí)際輸出與期望輸出之間的誤差，將誤差值沿網(wǎng)絡(luò)反向傳播，并修正連接權(quán)值，此為誤差反向傳播；給定另一個(gè)輸入，重復(fù)上述過程，直到網(wǎng)絡(luò)輸出的誤差達(dá)到允許的范圍或達(dá)到設(shè)定的訓(xùn)練次數(shù)為止，得到研究所需的神經(jīng)網(wǎng)絡(luò)模型（張青等，2016；劉淑梅等，2015；劉曉芬等，2011）。

表1 水質(zhì)參數(shù)敏感波段Table 1 Water quality parameter sensitive band

1.4.2 納什系數(shù)

納什效率系數(shù)（Nash-Sutcliffeefficiency coefficient），用于驗(yàn)證水文模型模擬結(jié)果的好壞，納什效率系數(shù)（E）的計(jì)算公式為：

式中：Qo指觀測值；Qm指模擬值；Qt表示第t時(shí)刻的某個(gè)值；Qo表示觀測值的平均值。E取值為負(fù)無窮至1，表示模式質(zhì)量好，模型可信度高；E接近 0，表示模擬結(jié)果接近觀測值的平均值水平，即總體結(jié)果可信，但過程模擬誤差大；E遠(yuǎn)遠(yuǎn)小于0，則模型不可信。

2 結(jié)果與分析

2.1 BP神經(jīng)網(wǎng)絡(luò)模型的構(gòu)建

本研究采用的是3層結(jié)構(gòu)的BP神經(jīng)網(wǎng)絡(luò)，在MATLAB中使用feedforwardnet函數(shù)，默認(rèn)的訓(xùn)練函數(shù)為 trainlm（Levenberg-Marquardt算法），選擇傳遞函數(shù)tansig，輸出函數(shù)pureline，隱含層節(jié)點(diǎn)數(shù)通過試選法確定。以2018年10月17日的Landsat 8遙感影像提取的 80個(gè)采樣點(diǎn)敏感波段反射率為輸入數(shù)據(jù)，以同天對(duì)應(yīng)的地面監(jiān)測水質(zhì)參數(shù)作為目標(biāo)數(shù)據(jù)，構(gòu)建一個(gè)3層的BP神經(jīng)網(wǎng)絡(luò)，包括一個(gè)輸入層、隱含層和輸出層。feedforwardnet函數(shù)將輸入數(shù)據(jù)分為3份，70%用于訓(xùn)練，15%用于檢驗(yàn)，15%用于測試，防止過度擬合。feedforwardnet默認(rèn)隱含層節(jié)點(diǎn)數(shù)為 10，學(xué)習(xí)步長為 0.05，迭代次數(shù)為50000，學(xué)習(xí)誤差為 1e-10。通過調(diào)整隱含層節(jié)點(diǎn)數(shù)、學(xué)習(xí)步長和迭代次數(shù)，對(duì)模型進(jìn)行反復(fù)訓(xùn)練，直到BP神經(jīng)網(wǎng)絡(luò)模型的r值達(dá)到最優(yōu)（表2），即r值不再上升為止。得到4個(gè)參數(shù)的最優(yōu)神經(jīng)網(wǎng)絡(luò)模型，網(wǎng)絡(luò)結(jié)構(gòu)見（表3）。

表2 BP神經(jīng)網(wǎng)絡(luò)模型相關(guān)系數(shù)rTable 2 The correlation coefficient of BP neural network model r

表3 各水質(zhì)參數(shù)BP-ANN結(jié)構(gòu)Table 3 BP-ANN structure of each water quality parameter

2.2 BP神經(jīng)網(wǎng)絡(luò)模型的驗(yàn)證

在訓(xùn)練區(qū)中隨機(jī)挑選 10個(gè)樣點(diǎn)，以訓(xùn)練好的BP神經(jīng)網(wǎng)絡(luò)分別對(duì)各采樣點(diǎn)的總氮、氨氮、總磷、COD含量進(jìn)行預(yù)測，得到預(yù)測結(jié)果對(duì)BP神經(jīng)網(wǎng)絡(luò)模型做進(jìn)一步驗(yàn)證（表4）。

從表4可以看出，驗(yàn)證點(diǎn)BP神經(jīng)網(wǎng)絡(luò)的預(yù)測值和實(shí)測值比較接近，特別是在平均值水平上預(yù)測值和實(shí)測值非常接近，總氮在誤差0.1的情況下可以認(rèn)為預(yù)測值和實(shí)測值的平均值相等，氨氮和總磷在誤差 0.01的情況下可以認(rèn)為預(yù)測值與實(shí)測值的平均值相等，COD的預(yù)測值與實(shí)測值的差距相對(duì)較大，但是 COD含量實(shí)測值的變化幅度較大，而且COD含量的國家水質(zhì)標(biāo)準(zhǔn)分類區(qū)間大，反演結(jié)果雖然與預(yù)測值有差距，仍可以在一定程度上反映密云水庫的COD含量變化情況。

為了更直觀的表現(xiàn)預(yù)測和實(shí)測值的變化情況，將實(shí)測值與預(yù)測值繪制成折線圖（圖2）。從圖2可以看出，總氮、氨氮、總磷含量的預(yù)測值與實(shí)測值比較接近，變化趨勢基本一致，預(yù)測值相對(duì)于實(shí)測值的浮動(dòng)較小，預(yù)測效果較好。COD含量的預(yù)測值與實(shí)測值的變化趨勢大致相同，但是預(yù)測值相對(duì)于實(shí)測值的偏差相對(duì)較大，COD的預(yù)測效果相對(duì)較差。因此本研究所建立的4個(gè)BP神經(jīng)網(wǎng)絡(luò)模型中，總氮、氨氮、總磷具有較好的反演能力，COD的反演能力相對(duì)較差。

為進(jìn)一步驗(yàn)證BP神經(jīng)網(wǎng)絡(luò)模型質(zhì)量，采用納什效率系數(shù)和均方根誤差，以 10個(gè)驗(yàn)證點(diǎn)的觀測值和預(yù)測值進(jìn)行計(jì)算，結(jié)果見（表5）。

從計(jì)算結(jié)果來看，氨氮、總磷和 COD的納什系數(shù)接近1，說明這幾個(gè)水質(zhì)參數(shù)的模型質(zhì)量較好，模型可信度高，其中是氨氮和總磷納什系數(shù)較高，而總氮的納什系數(shù)更接近于 0，說明總氮的神經(jīng)網(wǎng)絡(luò)模型其模擬結(jié)果接近觀測值的平均值水平，即總氮模擬結(jié)果總體結(jié)果可信，但模擬過程誤差大。均方根誤差顯示，總磷、總氮和氨氮的均方根誤差均小于0.03，其中最小的為總磷，其次是氨氮，再次是總氮，說明總磷、總氮和氨氮的預(yù)測值和實(shí)測值較為接近。COD的均方根誤差大于1，說明 COD的真實(shí)值和預(yù)測值有偏差。綜合來看，本研究建立的氨氮和總磷的BP神經(jīng)網(wǎng)絡(luò)模型預(yù)測效果較好，其次是氨氮，COD的神經(jīng)網(wǎng)絡(luò)模型預(yù)測效果較差。

表4 BP神經(jīng)網(wǎng)絡(luò)預(yù)測值與實(shí)測值的比較Table 4 Comparison between measured values and predicted values of BP neural network mg·L-1

圖2 BP神經(jīng)網(wǎng)絡(luò)模型結(jié)果驗(yàn)證Fig.2 Verification of BP neural network Model

表5 納什效率系數(shù)、均方根誤差Table 5 NSE (Nash-Sutcliffe efficiency coefficient),RMSE(Root Mean Square Error)

2.3 密云水庫主要水質(zhì)參數(shù)的變化分析

2.3.1 主要水質(zhì)指標(biāo)的年際變化

本研究基于構(gòu)建的 BP神經(jīng)網(wǎng)絡(luò)模型反演了2013—2018年密云水庫非結(jié)冰期每月的主要水質(zhì)參數(shù)，估算了每年主要水質(zhì)參數(shù)總氮、氨氮、總磷、COD的年平均值（圖3），分析2013—2018年密云水庫水質(zhì)的年際變化情況。

結(jié)果顯示，2013—2018年密云水庫的主要水質(zhì)指標(biāo)總氮、氨氮、總磷、COD整體均呈下降趨勢?？偟|(zhì)量濃度在2013年最高，達(dá)到0.932 mg·L-1，2017 年最低，為0.868 mg·L-1，2015—2018 年年際波動(dòng)較小，含量基本穩(wěn)定；氨氮質(zhì)量濃度在2013—2018年變化比較平緩，2015達(dá)到最高值0.018 mg·L-1，2015—2017年連續(xù)下降，2017年為最低值 0.0099 mg·L-1；總磷質(zhì)量濃度變化為先呈上升趨勢，2015年達(dá)到最高值0.027 mg·L-1，2015—2018年總磷含量表現(xiàn)為逐年小幅下降且趨勢變緩；COD變化趨勢與氨氮相似，2013—2015年呈上升的趨勢，2015年達(dá)到最高值 12.974 mg·L-1，2015—2017年快速下降，2017年為最低值6.17 mg·L-1，年際波動(dòng)比較大。

從 2013—2018年間密云水庫年均水質(zhì)指標(biāo)變化來看，氨氮和 COD的含量一直處于Ⅰ類水質(zhì)標(biāo)準(zhǔn)，且在Ⅰ類類別中含量較小，指標(biāo)優(yōu)良；而總磷近年為Ⅱ類水質(zhì)標(biāo)準(zhǔn)，總氮為Ⅲ類水質(zhì)?？傮w來看，密云水庫水質(zhì)較好，且年際變化顯示近年來密云水庫各項(xiàng)水質(zhì)指標(biāo)均呈下降趨勢，水質(zhì)整體呈上升趨勢。同時(shí)由年際變化可以看出，2018年總氮、氨氮、COD比上年有所增長。因此，密云水庫的水質(zhì)風(fēng)險(xiǎn)仍然存在，保護(hù)與監(jiān)測工作不容忽視。

2.3.2 主要水質(zhì)參數(shù)的季節(jié)變化分析

本研究選取的是 2013—2018非結(jié)冰期遙感影像，因冬季密云水庫區(qū)域人為活動(dòng)較少，生物生長緩慢，庫區(qū)水質(zhì)比較穩(wěn)定，故只研究春、夏、秋 3個(gè)季節(jié)主要水質(zhì)指標(biāo)的變化狀況。為了剔除偶然年份的影響，得到研究時(shí)段密云水庫水質(zhì)參數(shù)的季節(jié)性變化特征，本研究將反演結(jié)果按季節(jié)求得多年平均值，獲取了4個(gè)水質(zhì)指標(biāo)不同季節(jié)的變化分布圖（圖4）。

圖3 主要水質(zhì)參數(shù)的年際變化Fig.3 Interannual variation of major water quality parameters

由圖4可以看出，密云水庫主要水質(zhì)參數(shù)總氮、氨氮、總磷、COD的季節(jié)變化趨勢基本一致，均表現(xiàn)為春季濃度最高，秋季次之，夏季濃度最低，即密云水庫夏季水質(zhì)最優(yōu)，后因水資源的補(bǔ)充不足，水質(zhì)參數(shù)含量升高，導(dǎo)致秋春季水質(zhì)較差，春季達(dá)到最差。

2.3.3 密云水庫主要水質(zhì)參數(shù)的空間變化特征分析

本研究在空間柵格尺度上估算了 2013—2018年密云水庫水質(zhì)參數(shù)總氮、氨氮、總磷、COD的含量，揭示了主要水質(zhì)參數(shù)的空間分異與變化特征。

2.3.3.1 水質(zhì)參數(shù)總氮的空間變化特征

密云水庫總氮含量按照國家水質(zhì)分類標(biāo)準(zhǔn)總體以0.7 5mg·L-1和1 mg·L-1為分界點(diǎn)劃分為三類區(qū)域。從反演結(jié)果來看，總氮質(zhì)量濃度比較集中分布于0.75—1 mg·L-1這個(gè)區(qū)段。為了更加顯著突出總氮含量的空間變化，將這一區(qū)段再以 0.83 mg·L-1和0.91mg·L-1為分界點(diǎn)細(xì)分，共分為五級(jí)（圖5）。

從分級(jí)結(jié)果來看，北部區(qū)域的總氮含量高于其他區(qū)域。其中，2013年水庫的總氮質(zhì)量濃度主要集中在0.91—1 mg·L-1區(qū)域，中部和北部地區(qū)含量高于其他區(qū)域。2014年總氮含量明顯下降，集中區(qū)域轉(zhuǎn)變?yōu)?.83—0.91 mg·L-1區(qū)域，西部和北部岸邊含量顯著高于其他區(qū)域。2015年水域面積縮小，總氮含量明顯降低，東北部含量高于其他區(qū)域。到2016年水域面積擴(kuò)大，雖然總氮含量仍然集中在0.83—0.91 mg·L-1區(qū)段，但是相比2014年整體降低了，北部水域總氮含量高于其他區(qū)域。2017年水域面積進(jìn)一步擴(kuò)大，總氮含量在0.75—0.83 mg·L-1區(qū)段分布面積擴(kuò)大。2018年水域面積達(dá)到近年最大，總氮含量相比 2017年有所上升，主要分布區(qū)域位于東南部地區(qū)，這些區(qū)域遙感影像的云量偏高，這會(huì)造成估算值比實(shí)際值偏高?？傮w來看，總氮含量較高的水域面積呈減小趨勢，且集中分布于密云水庫北部和東南部。

2.3.3.2 水質(zhì)參數(shù)氨氮的空間變化特征

2013—2018年 密云水庫的氨氮質(zhì)量濃度小于0.03 mg·L-1。為了顯示密云水庫氨氮含量的空間差異性，以 0.006、0.012、0.018、0.024 mg·L-1為分界點(diǎn)將氨氮含量分為五級(jí)（圖6）。2013—2015年密云水庫的氨氮含量以0.012—0.018 mg·L-1為主，北部和東北部含量較高。2014年西部水域氨氮含量明顯降低，2015年東北部0.018—0.024 mg·L-1的氨氮水域面積不斷擴(kuò)大。2016年之后密云水庫氨氮含量明顯下降，以0.012—0.018 mg·L-1為主；2017年氨氮含量進(jìn)一步下降，以0.006—0.012 mg·L-1為主；2018年東部地區(qū)氨氮含量較上年有所上升。從密云水庫氨氮的空間分布上看，2015年較氨氮含量較高的區(qū)域分布面積最大，2015年后氨氮含量呈下降趨勢，2017年氨氮含量較低的區(qū)域分布面積最大。從水質(zhì)參數(shù)氨氮標(biāo)準(zhǔn)來看，2013—2018年密云水庫的水質(zhì)良好，氨氮含量處于Ⅰ類水質(zhì)標(biāo)準(zhǔn)，且遠(yuǎn)低于Ⅰ類閾值。

圖4 主要水質(zhì)參數(shù)季節(jié)變化Fig.4 Seasonal changes in major water quality parameters

圖5 2013—2018密云水庫水質(zhì)參數(shù)總氮空間分布圖Fig.5 Spatial distribution of total nitrogen in Miyun Reservoir during the period of 2013-2018

2.3.3.3 水質(zhì)參數(shù)總磷的空間變化特征

圖6 2013—2018密云水庫水質(zhì)參數(shù)氨氮空間分布圖Fig.6 Spatial distribution of ammonia nitrogen in Miyun Reservoir during the period of 2013-2018

圖7 2013—2018密云水庫水質(zhì)參數(shù)總磷空間分布圖Fig.7 Spatial distribution of total phosphorus in Miyun Reservoir during the period of 2013-2018

總磷按國家水質(zhì)分類標(biāo)準(zhǔn)標(biāo)分為四級(jí)（圖7）。2013年密云水庫的總磷質(zhì)量濃度集中在 0.01—0.025 mg·L-1區(qū)域，整體含量較低。2014年總磷含量略有升高，集中在 0.025—0.05 mg·L-1區(qū)域，北部和南部區(qū)域含量較高。2015年總磷含量較高區(qū)域主要分布于東北部。2016年密云水庫北部0.025—0.05 mg·L-1總磷分布面積擴(kuò)大。而到2017年之后北部區(qū)域總磷含量明顯降低；2018年南部地區(qū)總磷含量升高?？傮w來看，密云水庫的總磷含量處于Ⅱ類和Ⅲ類水質(zhì)標(biāo)準(zhǔn)。因此，密云水庫仍需加強(qiáng)上中游地區(qū)污染物總磷的排放和監(jiān)管。

2.3.3.4 水質(zhì)參數(shù)COD的空間變化特征

為了體現(xiàn)COD質(zhì)量濃度的空間變化，以7.5、15、20、30 mg·L-1為分界點(diǎn)將COD含量劃分為五級(jí)（圖8）。2013年密云水庫的COD質(zhì)量濃度以7.5—15 mg·L-1為主，東北部區(qū)域COD含量較高。東北部COD含量較高的區(qū)域隨著時(shí)間推移而不斷擴(kuò)大，2015年COD含量較高的區(qū)域擴(kuò)到最大。2016年東北部COD含量最高等級(jí)的區(qū)域迅速縮小，而20—30 mg·L-1的COD面積不斷擴(kuò)大。到2017年COD含量全面下降，COD含量主要集中在0—7.5 mg·L-1區(qū)段。2018年東部COD含量上升為 7.5—15 mg·L-1。總體而言，密云水庫的COD含量處于Ⅰ類和Ⅱ類水質(zhì)標(biāo)準(zhǔn)，符合飲用水源要求。

圖8 2013—2018密云水庫水質(zhì)參數(shù)COD空間分布特征Fig.8 Spatial distribution of COD in Miyun Reservoir during the period of 2013-2018

3 討論

分析結(jié)果表明，密云水庫的主要水質(zhì)參數(shù)總氮、氨氮、總磷、COD含量雖然年際之間有波動(dòng)，但整體呈下降趨勢，表明密云水庫的水質(zhì)在不斷改善，水質(zhì)保護(hù)工作取得了顯著成果。從各項(xiàng)水質(zhì)參數(shù)的季節(jié)平均值來看，夏季含量最低，秋季次之，春季最高，可能是由于夏季降雨較多，密云水庫的水源得到補(bǔ)給，水體得到更新，導(dǎo)致水質(zhì)相比其他季節(jié)更佳。密云水庫主要水質(zhì)指標(biāo)的空間分布顯示，庫區(qū)北部、和東部水質(zhì)指標(biāo)偏高，這種分布可能與北部和東部村莊密集，以及密云水庫兩大入庫河流：潮河和白河，從北部入庫有關(guān)。秦麗歡等（2017）研究表明，深水區(qū)和朝河及白河的過度區(qū)域總磷含量較高，這與本研究所得空間分布吻合。此外，班靜雅等（2019）研究表明南水北調(diào)工程使得庫區(qū)水位大幅上升，庫區(qū)周邊大面積土壤被永久淹沒，產(chǎn)生了土壤污物溶出，模擬結(jié)果表明對(duì)密云水庫水質(zhì)影響范圍小、程度低。筆者認(rèn)為，這也可能是導(dǎo)致北部和東部水域指標(biāo)含量較高的另一個(gè)原因?？傮w來看，密云水庫水質(zhì)較好，氨氮和COD兩個(gè)參數(shù)為Ⅰ類水質(zhì)，總磷為Ⅱ類，總氮為Ⅲ類，各項(xiàng)水質(zhì)參數(shù)指標(biāo)滿足飲用水源要求。然而密云水庫水質(zhì)監(jiān)測和保護(hù)工作不能放松，總磷和總氮在個(gè)別年份也存在小幅超標(biāo)的現(xiàn)象。

本研究開展的基于BP神經(jīng)網(wǎng)絡(luò)算法反演密云水庫水質(zhì)參數(shù)的研究達(dá)到了一定精度，但也存在一些不足。由于Landsat 8數(shù)據(jù)有限，再加上密云地區(qū)特殊的氣象條件，一些時(shí)段的遙感影像受到霧霾和云層的影響無法使用，導(dǎo)致遙感影像無法被充分利用，以后的研究需要進(jìn)一步探索針對(duì)水質(zhì)反演的大氣校正模型和遙感影像的去云技術(shù)，從而減小云量對(duì)于輻射值的影響，使遙感影像得到充分利用。另外，本研究采用的多光譜遙感影像分辨率為 30 m×30 m，而記錄采樣點(diǎn)坐標(biāo)的手持GPS的精度為15 m，因此兩者的配準(zhǔn)會(huì)產(chǎn)生一定的誤差，期待以后的研究通過技術(shù)手段的進(jìn)步，提高遙感影像與采樣點(diǎn)的配準(zhǔn)程度，使數(shù)據(jù)能夠更加準(zhǔn)確的反映水庫的水質(zhì)指標(biāo)的變化。

4 結(jié)論

本研究基于Landsa 8遙感數(shù)據(jù)與地面監(jiān)測水質(zhì)數(shù)據(jù)，采用遙感和GIS技術(shù)，通過構(gòu)建BP神經(jīng)網(wǎng)絡(luò)模型，實(shí)現(xiàn)了密云水庫2013—2018總磷、總氮、氨氮和COD（化學(xué)需氧量）4個(gè)水質(zhì)參數(shù)的反演，并對(duì)水質(zhì)參數(shù)的時(shí)空變化做了分析，得到如下結(jié)論：

（1）主要水質(zhì)參數(shù)的Landsat 8敏感波段為：總氮為1、4波段，氨氮為1—7波段，總磷為1、3—7波段，COD為2—5波段。對(duì)BP神經(jīng)網(wǎng)絡(luò)模型的預(yù)測結(jié)果驗(yàn)證顯示，密云水庫水質(zhì)指標(biāo)的預(yù)測值與實(shí)測值在平均值水平上非常接近，納什系數(shù)和均方根誤差評(píng)價(jià)結(jié)果表明，研究構(gòu)建的總磷和氨氮的BP神經(jīng)網(wǎng)絡(luò)模型質(zhì)量較好，可以很好實(shí)現(xiàn)水質(zhì)反演，而氨氮和COD模型質(zhì)量不佳，但在平均值水平上仍然具有準(zhǔn)確性。

（2）密云水庫主要水質(zhì)參數(shù)總氮、氨氮、總磷、COD的含量在2013—2018年總體呈下降趨勢，年均值顯示，氨氮和 COD兩個(gè)參數(shù)處于Ⅰ類水質(zhì)，總磷為Ⅱ類，總氮為Ⅲ類，主要各項(xiàng)水質(zhì)指標(biāo)滿足飲用水源要求，需注意加強(qiáng)總磷和總氮的控制。

（3）從季節(jié)變化特征分析顯示，4個(gè)水質(zhì)均表現(xiàn)為春季最高、秋季次之、夏季最低。

（4）密云水庫主要水質(zhì)參數(shù)的空間差異性較大，總體來看，水庫北部和東部水域的4個(gè)水質(zhì)參數(shù)含量相對(duì)較高，這與北部和東部村莊分布，以及密云水庫的兩大入庫河流：潮河和白河，從北部入庫有關(guān)。