基于PLSR的珠江口城市河流水質高光譜反演

黃 華，李茂億，陳吟暉，陳 耿，劉海龍，邢前國，蔡建楠

(1.廣東環境保護工程職業學院，廣東 佛山 528216；2.中山市環境監測站，廣東 中山 528403； 3.廣東省環境監測中心，廣東 廣州 510308；4.中國科學院煙臺海岸帶研究所，山東 煙臺 264003)

珠江口地區河網密布，河流水質狀況直接影響城市的整體風貌和居民生活質量[1]。近年來隨著國家《水污染防治行動計劃》的深入推進，城市河流的綜合整治和水質改善已成為珠江口各市的重要工作任務。為實現對河流水質狀況監控及河流整治效果的評估，往往需要配套開展密集、長期的水質監測工作，而采用常規的水樣現場采集和實驗室分析的方法存在耗時長、花費高、采樣點位有限等局限，無法滿足新形勢下水環境管理工作的需要[2]。

水質遙感技術基于水體水質與光譜特征之間的響應關系實現對水質指標的反演[3]，獲得的水質信息具有時間和空間上的相對連續性，可以作為傳統水質監測方法的有效補充，近年來得到迅速發展。其中，高光譜水質遙感技術由于具有納米級的光譜分辨率和連續的光譜譜段，極大地增加了水體光譜信息量，有助于精準建立各種水質參數與光譜響應關系以提高反演模型的精度，成為水質遙感研究領域的熱點。目前，相關研究主要面向的監測對象為水庫和海洋水體，而珠江口河網水體受陸源物質大量輸入和海洋潮汐因素的共同作用，有著與內陸二類水體和開闊大洋水體不同的水體光學特征，開展珠江口水體水質光學遙感研究具有重要的意義[4]。目前城市河流水質遙感有關研究中的目標水質參數主要為水體濁度、懸浮物、葉綠素a、有色可溶性有機物(CDOM)等光學活性參數[5-7]，對于河流水質綜合指標的遙感反演研究較少。

本文以珠江口城市中山市為研究區，通過開展河流高光譜監測和同步水質分析，研究受測水體高光譜特征與城市河流水質指數(city river water quality index，CWQI)的關系，并基于偏最小二乘回歸(partial least squares regression， PLSR)建立高光譜數據與河流水質指數的反演模型，同時探索水體高光譜監測的最佳光譜分辨率，以期為城市河流水質光譜監測及有關高光譜傳感器的設計和應用提供參考。

1 研究方法

中山市地處珠江口西岸，北接廣州，與深圳隔海相望，緊鄰港澳，處在粵港澳大灣區的地理中心位置。中山市境內河流密布，不同規模的內河涌、排洪渠高達1 000余條。改革開放以來，中山市經濟迅猛發展，但隨著人口和產業聚集程度不斷提高，河流水體受到了一定程度的污染。

1.1 河流光譜監測和數據處理

1.1.1地面高光譜監測

于2019年7—9月選擇中山市境內具有代表性的河流水體進行水質高光譜監測，共獲得151個點位的現場實測數據。高光譜監測使用美國Ocean Optics公司USB4000光譜儀，光譜分辨率約為 0.20 nm，光譜范圍為345.34～1 046.12 nm，視場角為15°。每次測量前均使用標準白板進行儀器定標。測量時將光譜儀探頭垂直置于離河流水面 0.5 m 處測量河流水體上行輻射亮度Lu、天空輻射亮度Ls和太陽下行輻照度Ed，每個點位重復進行3次測量后取平均值作為測定結果。測點的光譜反射率Rrs計算公式[8]為

(1)

式中γ為天空光在水體表面的菲涅爾反射率，通常取經驗值0.028?？紤]到400～900 nm波段范圍以外的光譜反射率更易受外界條件及儀器設計誤差等因素的影響，可能對監測結果帶來信號噪聲[9]，因此僅選擇400～900 nm范圍內共2 534個波段的原始光譜反射率數據用于構建反演模型。

1.1.2無人機載高光譜監測

無人機遙感監測技術具有靈活性大、影像分辨率高、時效性強、成本低等優點，是水環境遙感研究的一個重要發展方向[10]。在利用地面實測高光譜構建CQWI反演模型后，將其推廣應用于無人機高光譜影像中，以評價模型應用于無人機載高光譜監測的適用性。無人機載高光譜監測采用北京智科遠達數據技術有限公司ZK-VNIR-FPG480型高光譜成像掃描儀，其有效譜段范圍為400～1 000 nm，光譜通道數大于或等于270，光譜分辨率小于或等于 2.8 nm，空間分辨率為9 cm(100 m飛行高度)，視場角為26°。高光譜成像掃描儀搭載在大疆M600 PRO無人機平臺上，在晴朗天氣條件下，根據設定的航飛路線采集研究水體的光譜，并同步采集靶標的光譜數據，對水體數據進行定標運算得到水體反射率數據。無人機航拍區域主要為內河與外江交匯處、主要水工設施周邊、碼頭作業區等潛在的水質突變點區域，以盡可能同時獲取不同水質水體的光譜信息。

1.2 CWQI構建和水質監測

目前，國內外對于地表水環境質量評價方法主要有單因子水質類別評價法、水質指數法和不確定性方法三大類[11]，其中水質指數法可綜合多項水質指標的監測結果，計算過程簡單，評價結果便于比較。本文采用原環境保護部《城市地表水環境質量排名技術規定(試行)》中規定的河流水質指數作為反映河流水質狀況的綜合指標，單項指標(除溶解氧外)和溶解氧的CWQI計算公式分別為

(2)

(3)

式中：Ii、IDO分別為水質指標i和溶解氧的CWQI；ρi、ρ(DO)分別為水質指標i和溶解氧的質量濃度測定值；ρsi、ρs(DO)分別為水質指標i和溶解氧Ⅲ類地表水對應的質量濃度標準限值。所有水質指標的CWQI值相加即為該水樣的CWQI。

《城市地表水環境質量排名技術規定(試行)》推薦CWQI計算的基本水質指標為GB 3838—2002《地表水環境質量標準》表1中除水溫、糞大腸菌群和總氮以外的21項指標。為優化監測項目，對中山市各河流監測斷面2018—2019年的數據(共計水質監測樣本2 669個)進行了分析，結果顯示中山市各河流的主要超標因子為化學需氧量(COD)、總磷(TP)、氨氮(NH3-N)和溶解氧(DO)，這4項水質監測指標計算得到的CWQI值(I4)與21項水質監測指標計算得到的CWQI值(I21)之比的均值和標準差分別為53.7%和0.08%，相關系數r=0.94(P<0.05)，兩者具有顯著的相關性(圖1)，表明I4可較好地反映河流水質狀況，故本文中選取上述4項指標作為與光譜測量同步進行的水質監測指標，并計算I4值作為河流水質綜合評價結果。

表1 樣品水質特征

圖1 4項和21項水質監測指標CWQI值的關系

各河流監測點的水樣采集與地面高光譜測量同步進行，水質監測指標為COD、NH3-N、TP和DO。其中DO的質量濃度使用美國YSI公司Pro 2030多參數水質分析儀現場測定，其他指標取河流表層 0～20 cm處混合樣，裝入潔凈聚乙烯塑料瓶中并加入相應保存劑后置于低溫保存箱及時送至實驗室分析。COD質量濃度采用重鉻酸鉀法測定，NH3-N質量濃度采用納氏試劑光度法測定，TP質量濃度采用鉬銻抗分光光度法測定[12]，并使用空白樣、平行樣、標準樣品和加標回收等措施對樣品監測全過程進行了質量控制。水質樣品各特征指標見表1，可見多數樣品水質屬于Ⅳ或Ⅴ類，主要定類水質指標為TP和NH3-N。

1.3 PLSR反演模型建立與評價

由于水體反射率數據波段數眾多，相鄰波段間隔緊密，而不同波段光譜反射率數據間最大方差膨脹因子大于80 000，表明不同波段光譜反射率存在極嚴重的多重相關性[13]。為提高反演模型精度，采用PLSR建立水體不同波段光譜反射率與CWQI的反演模型。PLSR基于高維投影思想，通過提取主成分的方法對系統中的數據信息進行分解和篩選，實現對系統解釋性最強的綜合變量的有效抽取[13-15]；PLSR常用于樣本量小、變量間存在多重相關性的數據處理，適用于中高光譜遙感數據反演模型的建立[16-17]。

PLSR建模的基本步驟[18]為：①標準化處理水體光譜反射率數據矩陣X和水體CWQI數據矩陣Y，得到的標準化矩陣分別記為E0和F0；②分別提取E0和F0第一對成分T1和U1，T1和U1分別為E0和F0的線性組合并使之相關性最大；③分別求E0和F0在T1上的回歸，并得到回歸殘差矩陣E1和F1；④用殘差矩陣E1和F1代替E0和F0重復上述步驟；⑤檢查收斂性以確定提取的主成分數。

將151個樣本原始數據隨機劃分為建模集(120個樣本)和驗證集(31個樣本)兩部分，其中建模集用建立PLSR反演模型；驗證集用于反演模型對CWQI反演效果的驗證。模型驗證效果分別采用相關系數(r)、均方根誤差(RMSE)、平均相對誤差(MRE)進行評價。

2 結果與分析

2.1 光譜反射率與CWQI的關系

根據CWQI計算結果，151個水質樣品的CWQI值范圍為1.36～10.75，平均值為3.60，標準差為1.91。將不同CWQI值范圍的水樣進行分類并計算其高光譜反射率均值，結果(圖2)顯示，在400～700 nm波段范圍內光譜反射率對不同CWQI值的水樣有較好的區分度，水體光譜反射率均值隨CWQI的增大依次遞減，其中I4≤2的水體反射率均值在整個波段范圍內均高于其他類別的水體。但在 720～900 nm波段范圍內，出現I4>5的較差水質水體的光譜分辨率高于2

圖2 不同CWQI范圍水樣對應的高光譜反射率均值

2.2 不同光譜分辨率對PLSR反演模型的影響

為分析光譜分辨率對水體PLSR反演模型估算效果的影響，需找出水體CWQI高光譜遙感監測最佳光譜分辨率。采用分段三次Hermite插值多項式法對原始數據進行重采樣，采樣間隔依次為約 0.20 nm、1 nm、10 nm、20 nm、30 nm、40 nm、50 nm、60 nm、70 nm和80 nm，連同兩端波段的總波段數分別為2 534、900、51、26、18、14、11、10、9和8。在此基礎上建立不同光譜分辨率樣本集的PLSR反演模型，并根據相應的建模與驗證精度評價結果研究模型估算精度隨光譜分辨率的變化。結果(表2)表明，在各光譜分辨率的數據集中，原始光譜數據集的PLSR反演模型的估算效果最差；對原始數據進行不同光譜分辨率的重采樣后，均有助于提升PLSR反演模型的估算效果。PLSR反演模型的驗證集RMSE隨著光譜分辨率的降低呈現先減小后增大的趨勢。光譜分辨率為50 nm、波段數為11個時，驗證集RMSE最小達0.790，MRE為18.63%，模型估算效果最優。但總體上重采樣后各驗證集的PLSR反演模型估算效果均處于接近的水平。該現象表明，在CWQI的遙感監測中，過高的光譜分辨率可能會造成數據冗余和引入環境噪聲，造成模型復雜化并影響模型穩定性。反之，在較低的光譜分辨率下，數量較少的波段可在基本涵蓋水體主要污染物特征譜段的同時，增強光譜信噪比，從而提高估算精度。

表2 不同波段數的PLSR建模效果評價

2.3 提取的主成分數對PLSR反演模型的影響

確定提取的主成分數是PLSR反演模型建立和優化的另一重要環節。目前PLSR建模確定主成分數的方法主要有舍一交叉驗證方法、分批交叉驗證方法、分裂樣本交叉驗證方法、隨機樣本交叉驗證方法等[17]，其中舍一交叉驗證方法的應用較為廣泛，但該方法得到的建模主成分數可能不是最優值，主成分數過多或過少可能會導致欠擬合或過擬合，影響模型的估算準確度[23]。本文采用50 nm分辨率的水體光譜反射率數據，分別提取2～10個主成分進行PLSR建模，以驗證集RMSE最小時對應的主成分數為最佳的主成分數。圖3表明，模型的RMSE隨主成分數的增加呈先下降后增加的趨勢，當主成分數為8時，模型的RMSE達到最低值。

圖3 不同主成分數與驗證集RMSE的關系

因此，確定建模最佳主成分數為8，并建立PLSR反演模型如下：

I=4.84-23.33R400+22.62R450+53.53R500+

4.69R550-24.19R600-102.32R650+

75.96R700+25.38R750-30.93R800+

46.15R850-19.99R900

(4)

式中Rj為不同波長處的水體光譜反射率。

2.4 PLSR反演模型估算效果評價

以31個驗證集樣本的CWQI實測值為縱坐標、CWQI高光譜反演模型(式(4))計算值為橫坐標，繪制其實測值與計算值的對比散點圖(圖4)。結果表明，該反演模型的估算精度良好，RMSE、MRE分別為0.768和18.1%，說明該反演模型在有效縮減了高光譜輸入參數的前提下，仍然保持了良好的泛化能力和反演精度。

圖4 驗證集樣本CWQI實測值與反演模型計算值比較

2.5 反演結果

由于無人機航拍覆蓋范圍較大，本文未進行同步的水質分析，主要對基于無人機高光譜影像的CWQI反演效果做初步定性分析，并重點關注河流水質突變區域CWQI反演結果的合理性。從應用結果(圖5)上看，反演得到的CWQI值可較合理地反映出河流水質的空間差異。如圖5(a)反演結果突出了地理位置接近，但因水閘控制導致連通性差的水體水質的顯著差異；圖5(b)反演結果反映了河流交匯處水體混合過程中的局部水質差異；圖5(c)反演結果突出了碼頭船舶作業引起的水體擾動和水質變化。該結果表明，基于地面高光譜測點數據建立的CWQI反演模型應用于空間影像數據反演CWQI具有較好的可行性。

(a) 某水閘兩側水體

受限于無人機航程的限制，每次作業時間約在20 min左右，對于河流水質大范圍監測的能力較為不足，僅能重點開展局部區域的水質高光譜反演。同時，現場航飛采用多航帶作業方式，但對于水體流動性強，漲落潮劇烈的水域，水體流動導致水質的時空變化不能被無人機同步采集；這種準同步的水體光譜數據與實際河流的水質分布會有一定差異。但整體上反演模型應用于無人機高光譜影像的CWQI反演可取得良好效果，表明該方法具有良好的業務化運行潛力。

3 結 語

2019年7—9月在珠三角地區中山市選擇代表性河流進行了水體高光譜反射率監測和同步水質監測，并以基于COD、TP、NH3-N和DO 4項水質指標計算得到的CWQI來表征河流水質狀況。CWQI與水體高光譜反射率的相關分析表明，不同波段高光譜反射率對水質狀況有較好的光學敏感性；基于PLSR建立高光譜數據與CWQI的最優反演模型，其驗證數據集RMSE和MRE分別為0.768和18.1%，具有良好的反演精度。將該模型應用于無人機高光譜影像的CWQI反演，可以較合理地反映出河流水質的空間差異。盡管無人機航程、水流變化復雜等因素可能對無人機高光譜水質反演產生影響，但水質高光譜反演模型和無人機結合，可為河流水質綜合評價提供一種大范圍覆蓋、可視化表達的技術手段，具有良好的業務化運行潛力。