王歆暉，田 華，季鐵梅，鞏彩蘭，胡 勇，李 瀾，何志杰

（1.中國科學院大學 物理科學學院，北京 100049；2.中國科學院 上海技術(shù)物理研究所，上海 200083；3.上海市水文總站，上海 200232；4.中國科學院 紅外成像光譜重點實驗室，上海 200083）

0 引言

內(nèi)陸河流作為重要的資源和載體，是生態(tài)系統(tǒng)的綠色生命線，關(guān)系到人類的生存和發(fā)展。然而，河流水環(huán)境極易受到污染，特別是城市河道，隨著現(xiàn)代城市的飛速發(fā)展與擴張，人口密度增加，給城市河流水環(huán)境帶來了挑戰(zhàn)。從2017 年起，上海市就推行“河長制”，目的在于加強對上海市市管河道的監(jiān)管，對水污染進行防治，保護水資源。因此，保障居民的正常生產(chǎn)生活及推進可持續(xù)發(fā)展，對城市河流水質(zhì)快速、準確的監(jiān)測工作提出了高要求。目前，河流水質(zhì)監(jiān)測主要采用人工實地采樣與實驗室分析的方法，該方法結(jié)果準確，但只能獲取河流截面采樣點數(shù)據(jù)［1-2］。與其相比，內(nèi)陸河流遙感監(jiān)測具有以下優(yōu)勢：1）遙感技術(shù)具有快速、大范圍和周期性的特點，可以彌補常規(guī)定時、定點監(jiān)測的不足，節(jié)省大量人力、物力和財力［3］；2）隨著對水體光譜特性的深入研究以及水質(zhì)參數(shù)反演模型的改進，遙感圖像能夠較為準確地獲取河流水質(zhì)信息，與地面實測水文參數(shù)、地理位置信息等結(jié)合，可有效發(fā)現(xiàn)污染水體以及水質(zhì)變化趨勢，作為水質(zhì)管理的參考。

隨著新的衛(wèi)星傳感器不斷升空，時間分辨率和空間分辨率不斷提高，水質(zhì)遙感監(jiān)測從定性發(fā)展到定量，可實現(xiàn)遙感反演的水質(zhì)指標也由葉綠素a、懸浮物，發(fā)展到總氮、總磷、高錳酸鹽指數(shù)等［4-7］。當前，利用衛(wèi)星遙感技術(shù)對水質(zhì)監(jiān)測的研究主要集中在區(qū)域跨度較大的研究對象上，例如大型湖泊水庫和河口海岸的監(jiān)測，對小尺度和微型態(tài)環(huán)境問題研究相對較少。MATTEWS［8］利用MERIS 數(shù)據(jù)發(fā)展了改進的最大峰值算法，對內(nèi)陸和近岸水體葉綠素a 濃度進行探測，研究了非洲南部50 多個湖庫水體從2002 到2012 年的富營養(yǎng)化及其變化狀況。HOU等［9］利用MODIS 地表反射率數(shù)據(jù)645 nm 波段反演了長江中下游102 個大型湖泊的水體總懸浮濃度，并分析了其變化趨勢。國內(nèi)學者臧友華［10］將遺傳算法和支持向量機組合，提出一種新的反演模型GA-SVM 模型，采用ETM+數(shù)據(jù)對渭河水質(zhì)參數(shù)葉綠素a、透明度、總磷進行了反演。鞏彩蘭等［11］對黃浦江進行了水體光譜測量和同步水質(zhì)采樣，分析了各水質(zhì)參數(shù)與反射率及相關(guān)反射率特征的相關(guān)性，建立了常規(guī)水質(zhì)參數(shù)與反射率之間的關(guān)系模型。以上研究主要是針對部分具有光學特性的遙感水質(zhì)參數(shù)或者是國家《地表水環(huán)境質(zhì)量標準》中規(guī)定的單一水質(zhì)參數(shù)，雖然簡單直觀地反映了單因子污染狀況，但是缺乏對河流水質(zhì)的綜合研究。因此，本文利用5 項主要水質(zhì)參數(shù)建立河流綜合水質(zhì)遙感模型，以綜合水質(zhì)指標來表征河流水質(zhì)狀況。

本文采用Sentinel-2A 衛(wèi)星遙感影像數(shù)據(jù)和同步地面河流斷面實測水質(zhì)參數(shù)數(shù)據(jù)，對5 項水質(zhì)參數(shù)溶解氧、高錳酸鹽指數(shù)、五日生化需氧量、氨氮以及總磷進行因子分析，在此基礎(chǔ)上建立水質(zhì)遙感反演的模型，得到綜合水質(zhì)指標用于確定水質(zhì)類別。

1 實驗數(shù)據(jù)及預處理

1.1 實驗數(shù)據(jù)

本文采用歐洲航空局發(fā)射的Sentinel-2A 衛(wèi)星搭載的多光譜成像儀獲取的遙感影像，影像成像日期為2018 年4 月9 日，區(qū)域為上海市，具體研究選擇上海市青浦區(qū)和松江區(qū)的部分河流流域，研究區(qū)域如圖1所示。影像覆蓋13個光譜波段，幅寬達290 km，重訪周期10 d。從可見光和近紅外到短波紅外，具有不同的空間分辨率。

圖1 研究區(qū)域衛(wèi)星影像band 8、4、3 假彩色合成圖Fig.1 False color composite images of bands 8，4，and 3 in the study area

考慮空間分辨率，本研究用到的波段為band 2、band 3、band 4、band 8（見表1）。

表1 Sentinel-2A 衛(wèi)星部分波段信息Tab.1 Partial band information of Sentinel-2A satellite

水質(zhì)參數(shù)采樣數(shù)據(jù)獲取的是衛(wèi)星成像同一天的上海市河流斷面實地采樣分析數(shù)據(jù)，按照《地表水環(huán)境質(zhì)量標準（GB 3838—2002）》對上海市主要河流94 個斷面采樣點溶解氧（DO）、高錳酸鹽指數(shù)（CODMn）、五日生化需氧量（BOD5）、氨氮（NH3-N）、總磷（TP）共5 項指標進行評價，如圖2 所示。總體上，溶解氧、高錳酸鹽指數(shù)、總磷評價等級較好，基本達到3 類水的標準。超標的主要水質(zhì)參數(shù)是氨氮，通常來源于飼料、水生動物的排泄物、肥料及動植物尸體分解，未處理或處理過的城市生活和工業(yè)廢水、各種浸濾液和地表徑流是氨氮含量升高的主要原因［12］。

圖2 各水質(zhì)參數(shù)所屬水體類型分布情況圖Fig.2 Water grade summary of each water quality parameter

1.2 實驗數(shù)據(jù)預處理

對于經(jīng)過輻射定標等預處理后的反射率影像，運用歸一化差異水體指數(shù)（NDWI）對2018年4月9 日的實驗影像研究區(qū)進行河流水體提取。對NDWI 圖像設(shè)置掩模模板提取水體，閾值往往設(shè)為0，即NDWI 圖像大于0 的像素被認為是水體［13］：

式中：bandgreen為第3 波段（綠光）反射率；bandNIR為第8 波段（近紅外）反射率。

將影像與水體提取掩模模板疊加對比，掩模結(jié)果如圖3 所示，其中白色部分為識別的水體，而黑色部分為非水體。目視解譯結(jié)果表明，歸一化差異水體指數(shù)NDWI 能夠較好地識別出河流。

圖3 2018 年4 月研究區(qū)水體NDWI 提取掩模結(jié)果Fig.3 Mask results of water extraction with NDWI in the study area in April 2018

2 基于因子分析的水質(zhì)遙感反演方法

2.1 數(shù)據(jù)波段篩選

在正式進行建模前對波段進行了篩選，結(jié)合單波段與各主要水質(zhì)參數(shù)的相關(guān)性以及OIF 指數(shù)，來選擇對水質(zhì)參數(shù)變化較為敏感的特征波段或者波段組合。本文分析了各單波段與各水質(zhì)參數(shù)的相關(guān)性，整體相關(guān)性較低，Pearson 系數(shù)整體不超過0.4，所以Sentinel-2 單波段并不適合進行水質(zhì)遙感反演，需要考慮其他波段組合。研究區(qū)主要超標的水質(zhì)參數(shù)為氨氮，各波段與氨氮相關(guān)系數(shù)按照從大到小順序排列，波段順序為band 3、band 4、band 2、band 8。OIF 指數(shù)［14］是進行波段篩選的依據(jù)之一，通常波段數(shù)據(jù)的標準差越大，包含的信息量越多；而波段之間的相關(guān)系數(shù)越小，則波段之間獨立性越高，數(shù)學表達式為

式中：Si為第i個波段的標準差；Rij為第i、j兩波段的相關(guān)系數(shù)。

不同三波段組合的OIF 指數(shù)見表2，按照從大到小的順序排列來選擇最優(yōu)組合方案，理想組合順序為：1）band 2、band 4、band 8；2）band 3、band 4、band 8；3）band 2、band 3、band 8。因組合1 和組合2的OIF 指數(shù)值差異并不明顯，結(jié)合研究區(qū)主要超標的水質(zhì)參數(shù)氨氮與波段的相關(guān)性，最終確定選用組合2 波段3、4、8 組合。

2.2 綜合因子分析理論

因子分析是一種處理多元變量數(shù)據(jù)的統(tǒng)計方法，研究眾多變量之間的內(nèi)在關(guān)系，對可觀測變量進行概括和抽象，用較少的因子變量來最大程度上解釋原始的觀測變量［15-16］，其主要表達形式為

式中：X1，X2，…，Xp為p個存在相關(guān)關(guān)系的可觀測變量；F1，F(xiàn)2，…，F(xiàn)m為m個互相獨立的因子，又稱為公共因子，是不可觀測的變量；ε1，ε2，…，εp為特殊因子，是無法包含在公共因子中的部分；aij(1≤i≤p，1≤j≤m)為第i個變量與第j個因子的相關(guān)系數(shù)。

基于因子分析建立遙感反演模型，首先計算原始變量的相關(guān)系數(shù)矩陣，然后利用相關(guān)系數(shù)矩陣結(jié)合反演模型和先驗知識提取公共因子，通常只取方差大于1 或者特征值大于1 的因子，或者按照因子的累計方差貢獻率來確定，一般要求至少達到70%。因子Fj的方差貢獻是式（3）中對應列相關(guān)系數(shù)的平方和，根據(jù)如下公式計算方差貢獻率，方差貢獻率進行累加即為累計方差貢獻率：

接著進行因子旋轉(zhuǎn)，其目的是通過坐標變換使因子實際意義更易解釋，最后根據(jù)如下公式計算綜合因子得分：

式中：ωj(1≤j≤n)為對應因子Fj的方差貢獻率。

2.3 基于綜合水質(zhì)因子的遙感反演模型

由于水質(zhì)參數(shù)本身數(shù)量級的不一致，為了避免對模型造成影響，對數(shù)據(jù)進行最大最小化標準化，將其映射到-1 至1 范圍內(nèi)。采用因子分析法對研究區(qū)水質(zhì)參數(shù)進行分析，首先KMO 值為0.763，大于閾值0.5，說明了變量之間是存在相關(guān)性的，符合要求；然后是Bartlett 球形檢驗的結(jié)果，Sig 檢驗概率值為0.000，小于閾值0.05，說明結(jié)果顯著，表明數(shù)據(jù)可以進行因子分析。

河流水質(zhì)參數(shù)因子分析的主要因子成分如圖4所示。主要因子1 的特征值為2.8，方差貢獻率為56.74%；主要因子2 的特征值為0.95，方差貢獻率為18.93%；前3 個主要因子的方差貢獻率之和為87.5%。因此，水質(zhì)參數(shù)監(jiān)測數(shù)據(jù)原始5 個變量的信息能夠用變換后的3 個主要因子來表達。

圖4 監(jiān)測數(shù)據(jù)的主要因子成分Fig.4 Major factor components of the monitoring data

由表3 可知，河流水質(zhì)參數(shù)因子分析的旋轉(zhuǎn)載荷矩陣和貢獻率，主要抽取了方差貢獻率較大的前3 個主要因子成分，因子成分F1主要受高錳酸鹽指數(shù)和五日生化需氧量控制，系數(shù)分別為0.928 和0.832，這兩項指標越大，反映水體受有機物的污染越嚴重。因子成分F2主要由總磷和氨氮控制，系數(shù)分別為0.921 和0.649。通常是城市生活洗滌污水、垃圾等包含的無機鹽類，與水體的富營養(yǎng)化有關(guān)。因子成分F3則主要由溶解氧控制，系數(shù)值最大為0.970。

表3 因子旋轉(zhuǎn)載荷矩陣與方差貢獻率Tab.3 Factor loading matrix and variance contribution rate

5 項主要水質(zhì)參數(shù)作為式（1）中各可觀測變量，經(jīng)過因子旋轉(zhuǎn)后借助表3 的因子旋轉(zhuǎn)載荷矩陣，可以得到如下3 個主要因子成分的表達式：

將各個公共因子的方差貢獻率等參數(shù)代入式（5），得到水質(zhì)綜合因子得分模型為

為了得到最終水質(zhì)遙感反演模型，需要建立各個水質(zhì)參數(shù)和遙感影像波段3、4、8 之間的關(guān)系。以影像波段反射率數(shù)據(jù)作為自變量，水質(zhì)參數(shù)數(shù)據(jù)作為因變量，采用最小二乘法確定多元線性回歸模型的待定參數(shù)，當實測值與回歸模型得到的擬合值的差值平方和達到最小時的模型參數(shù)即為所求待定參數(shù)。根據(jù)研究數(shù)據(jù)得到的回歸模型結(jié)果見表4。由回歸分析結(jié)果可見，溶解氧、高錳酸鹽與波段之間的關(guān)系相對較密切。

表4 不同水質(zhì)參數(shù)的線性回歸模型Tab.4 Regression models of different water quality parameters

根據(jù)得出的綜合因子得分模型以及水質(zhì)參數(shù)和波段的回歸分析結(jié)果，計算得到最終水質(zhì)遙感反演模型為

3 水質(zhì)遙感反演結(jié)果與精度分析

對預處理［17］后的2019 年4 月9 日的研究區(qū)影像，應用綜合水質(zhì)因子反演模型進行反演，得到研究區(qū)的綜合水質(zhì)圖，如圖5 所示。一般情況下，綜合水質(zhì)F值越小，代表水體越潔凈，否則水體則存在水質(zhì)參數(shù)超標，或者受到污染的可能性。從結(jié)果圖來看，主干河流的值基本低于支流，是符合實際水質(zhì)情況的。

圖5 2019 年4 月研究區(qū)綜合水質(zhì)圖Fig.5 Comprehensive water quality map of the study area in April 2019

在綜合水質(zhì)圖的基礎(chǔ)上，據(jù)不同水質(zhì)等級的情況分別設(shè)定閾值，最終確定不同情況下綜合水質(zhì)F對應的水質(zhì)等級，繪制綜合水質(zhì)等級圖，如圖6 所示。因衛(wèi)星遙感影像分辨率的限制，圖中小河流存在斷裂、破碎等狀況。從結(jié)果圖看Ⅲ類和Ⅳ類水主要分布在主干河流，相對而言，Ⅴ類水體在支流較多，整體水質(zhì)狀況良好，劣5 類水體只出現(xiàn)在局部區(qū)域。

隨機選取30 個樣本數(shù)據(jù)，將反演的綜合水質(zhì)等級結(jié)果與實際測量采樣的結(jié)果進行對比驗證，如圖7 所示，水質(zhì)類別2、3、4、5、6 分別對應標準水質(zhì)類別Ⅱ類、Ⅲ類、Ⅳ類、Ⅴ類以及劣Ⅴ類。

圖7 實測水質(zhì)與反演綜合水質(zhì)類別對比圖Fig.7 Comparison of the measured water quality and the inversion comprehensive water quality

在大河道區(qū)域準確度能夠達到71.2%，但是細小的支流準確度只能夠達到42.8%，誤差的主要來源如下：1）地面實測采樣水質(zhì)數(shù)據(jù)與衛(wèi)星影像時間和空間的不一致性；2）僅以5 項主要水質(zhì)參數(shù)表征河流水質(zhì)狀況，仍存在其他因素的干擾；3）衛(wèi)星數(shù)據(jù)空間分辨率的限制會帶來的小河流混合像元問題。

4 結(jié)束語

利用2018 年4 月Sentinel-2A 遙感影像以及同步獲取的河流斷面實測水質(zhì)參數(shù)數(shù)據(jù)，結(jié)合5 項主要水質(zhì)參數(shù)溶解氧、高錳酸鹽指數(shù)、五日生化需氧量、氨氮以及總磷，在因子分析的基礎(chǔ)上建立水質(zhì)遙感反演的模型，得到綜合水質(zhì)指標用于確定水質(zhì)類別。實驗結(jié)果表明：1）上海市主要超標的水質(zhì)參數(shù)指標為氨氮，水污染治理應著重處理難降解的有機物以及含氮的可溶性無機物等；2）綜合水質(zhì)指標值越小，內(nèi)陸河流水質(zhì)狀況越良好，水體更清潔；3）本文方法對于大河流的水質(zhì)等級判定的正確率達到了71.2%?；谝蜃臃治龅乃|(zhì)遙感反演模型能夠應用于內(nèi)陸河流水質(zhì)遙感監(jiān)測，為水環(huán)境部門提供參考信息。但是本文僅考慮了5 項水質(zhì)參數(shù)，對于總體水質(zhì)狀況的表達還略有不足，需要考慮其他水質(zhì)遙感監(jiān)測指標，輔助綜合水質(zhì)類別判定；另外，該模型僅適用于一個季節(jié)，獲取長時間序列的影像以及同步水質(zhì)參數(shù)數(shù)據(jù)，能夠?qū)δＰ瓦M行改進，并且利用多時相數(shù)據(jù)對內(nèi)陸河流水質(zhì)狀況進行動態(tài)變化研究。