李 瀾,田 華,季鐵梅,鞏彩蘭,胡 勇,王歆暉,何志杰
1. 中國科學院上海技術物理研究所,上海 200083 2. 中國科學院紅外探測與成像技術重點實驗室,上海 200083 3. 中國科學院大學,北京 100049 4. 上海市水文總站,上海 200232
城市河道是城市生態(tài)環(huán)境的重要組成部分,改革開放以來我國城市經濟快速發(fā)展,但也造成河道出現不同程度的污染[1-2]。城市河道水污染治理的前提和關鍵是水質狀況的實時監(jiān)測; 目前水質監(jiān)測方法以傳統(tǒng)采樣監(jiān)測為主遙感監(jiān)測為輔,傳統(tǒng)的水質監(jiān)測通過實驗室化驗確定水質指標濃度再根據相關標準確定水質級別; 遙感監(jiān)測方法是通過衛(wèi)星、飛機、無人機等平臺搭載的傳感器獲取河道區(qū)域的遙感影像,再經過一系列處理分析建立水質參數與光譜反射率之間的定量關系,進而確定河道的水質級別。傳統(tǒng)采樣監(jiān)測技術已較為成熟,但存在費時、費力、只能獲得斷面點位水質等不足。隨著遙感和高光譜技術的發(fā)展,可以準確、快速地獲取各種地物的各波段的光譜信息,國內外學者提出了眾多基于遙感反演水質參數進而判斷水質類型的研究方法[3-6]。王麗艷等使用MODIS數據反演了呼倫湖水體COD濃度并確定水質級別[7]。Wang等提出了一種綜合水質評價方法,并利用該方法反演了多個水質參數濃度[8]。但河道水質類別的劃分指標如溶解氧(DO)、總磷(TP)、氨氮(NH3-N)等不具有明顯的光學特性,通過遙感反演各水質指標進而判斷水質類型這類方法的局限性在于最終級別的識別精度受各個水質指標反演精度共同決定,單一指標反演精度不高必然整體識別精度低下,模型應用局限性大、泛化能力不強,適用性不足。
水體反射率是諸多水質參數綜合作用下的結果,水體的光譜特征與水質類型之間存在某種映射關系。提出一種基于光譜二階微分波動指數的光譜分析方法,能有效反映各波段光譜二階微分的波動情況,在此基礎上可實現城市河道水質類型的快速識別,可推廣應用于利用其他遙感平臺獲得的高光譜數據監(jiān)測內陸水體的水質類型。
基于光譜二階微分波動特性的水質光譜特征提取及分類方法,其基本思想在于利用上下包絡線來提取光譜二階微分曲線的振動范圍,具體的數據處理流程如圖1。
圖1 基于光譜二階微分波動指數計算流程圖Fig.1 Flow chart of second-order differentialfluctuation index calculation
(1)光譜二階微分計算。光譜二階微分計算包括差分計算和平滑兩部分,首先利用光譜差分計算原始光譜的二階微分,隨后使用Savitzky-Golay Smoothing法[9]對原始二階光譜進行平滑處理,消除噪點和其他干擾。
(2)局部極大值極小值提取。曲線的波峰和波谷通常對應局部極大值和極小值,而波峰波谷的值可以反映曲線局部的變化范圍,使用3×3的滑動窗口逐點提取局部極值,滿足式(1)或(2)的即為局部極大值或極小值點。
(1)
(2)
(3)無關極值點去除。曲線上的波峰波谷對應著局部極值點,通常光譜二階微分曲線某一波段附近只會對應一個峰值,由于儀器探測性能限制或其他原因的影響,峰值周圍會出現一系列無關極值點,這會對后續(xù)處理產生干擾,為此可設定一個距離閾值m,采用循環(huán)迭代方式去除無關極值點,下面以極大值為例介紹具體實現步驟。①尋找當前所有極大值中最大值以及對應的波長位置; ②去除最大值周圍距離小于閾值m的所有極大值點; ③彈出并保存最大值和波長位置。重復①②③,直到所有的無關極大值點均消除。無關極小值的去除與上述步驟類似,只需每次彈出最小值和位置即可。
(4)雙包絡線生成。在上一步處理后可得到光譜二階微分曲線上的各波峰點和波谷點,使用三次樣條插值法得到光譜二階微分曲線的上下兩條包絡線。圖2是某光譜二階微分及雙包絡線示意圖。
圖2 光譜二階微分曲線與雙包絡線Fig.2 Second-order differential curve and double envelope
(5)光譜二階微分波動指數計算。上下包絡線能表示光譜二階微分的振動范圍,而某一波段處的對應上下包絡線函數值之間的差值能表示該波段附近光譜二階微分變化的劇烈程度。因此可定義一種光譜二階微分波動特性指數來描述這一特征,具體計算公式如式(3)
F(λ)=eA(f(λ)-g(λ))
(3)
式(3)中,F(λ)表示波長為λ處附近的變化特征,f(λ)和g(λ)分別對應該波段上下包絡線的函數值,A是固定的放大系數用于增加區(qū)分度。
先后多次在上海市嘉定區(qū)采集了各水質類型水體的光譜數據,光譜數據獲取采用ASD公司生產的FieldSpec3便攜式地物光譜儀,波長范圍選擇為400~900 nm。測量和處理參照國家標準《水體可見光-短波紅外光譜反射率測量GB/T 36540—2018》進行。圖3是獲取的水體反射光譜曲線圖。
圖3 反射光譜曲線圖Fig.3 Reflectance spectral curves
采集光譜數據的同時,還需要獲取相應的水質數據。將現場采樣的水體樣本密封保存帶回實驗室化驗獲得各水質指標,再根據《地表水環(huán)境質量評價辦法(試行)》等標準劃分類別,共得到Ⅰ和Ⅱ類樣本11個(由于Ⅰ類樣本稀少,將其與Ⅱ類樣本合并)、Ⅲ類14個、Ⅳ類20個、Ⅴ類23,劣Ⅴ類20個。
采用上述的數據處理方法計算各類水體的光譜二階微分波動系數,合適的距離閾值和放大系數可以有效去除虛假極值點,增加不同類別之間的區(qū)分度,經多次試驗最小距離閾值m設置為10,放大系數A設為500效果較好,圖4是計算的各類水體的平均光譜二階微分波動指數曲線圖。
圖4 各類型水體的平均波動指數Fig.4 Average fluctuation index of each type of water body
從圖4可以看出,各類水體波動指數在720~740,750~770以及820~840 nm處均有明顯的峰值,且不同水質類型水體峰值高度不同,隨后統(tǒng)計了各類樣本中這三個波段內的平均波動指數的均值、標準差等特征,結果見表1—表3。
上述統(tǒng)計結果表明,水體在720~740,750~770和820~840 nm內的光譜二階微分平均波動系數與水質級別具有正向相關關系,隨著水質級別的上升,這三個波段內的平均波動指數的平均值和中位數也隨之上升,且水質級別越高上升幅度越大,變化越劇烈,同樣平均波動指數的標準差也在增加。因此可利用這一特性實現河道水質類型的自動識別分類。
表1 各類型水體720~740 nm平均波動指數統(tǒng)計Table 1 Average fluctuation index from 720~740 nm
表2 各類型水體750~770 nm平均波動指數統(tǒng)計Table 2 Average fluctuation index from 750~770 nm
表3 各類型水體820~840 nm平均波動指數統(tǒng)計Table 3 Average fluctuation index from 820~840 nm
首先將各類樣本按照2∶1的比例隨機劃分為測試集和訓練集,隨后采用本方法計算訓練集內各樣本的二階微分波動指數曲線,得到三個特征波段內的平均波動指數,最后以這三個特征波段的平均波動指數為輸入,使用LSSVM[10]構建水質類型識別模型,訓練完成后進行測試,取重復10次實驗的平均值作為最終結果,具體結果見圖5(a)。
圖5(a)顯示平均識別準確度達80.65%,而通過對誤分類的樣本進行分析后發(fā)現,誤分類樣本多集中在兩類樣本特征邊界處,故單獨測試了分類誤差不超過1類的識別準確度,重復測試10次,如圖5(b)所示。結果表明誤差不超過1類的平均識別準確度達96.77%,在31個測試樣本中,僅1個樣本分類誤差超過1類,因此光譜二階微分波動系數結合機器學習模型可實現城市河道水質的快速識別,且識別誤差不超過1類。
針對目前城市水體大面積水質類型調查的需求,提出了基于光譜二階微分波動指數的高光譜數據處理和分析方法,在對各類型水體的波動指數分析后得到如下結論: (1)720~740,750~770和820~840 nm波段內的平均波動指數與水質級別具有正相關關系,隨著水質級別的上升,平均波動指數值升高,對應波段附近的光譜二階微分振幅越大,這表明水質狀況越差的水體,其光譜二階微分波動越大。(2)以三個特征波段的平均波動指數作為輸入特征,采用LSSVM模型構建水質類型識別模型,經測試,模型的平均識別準確度達80.65%,誤差不超過一類的識別準確度達96.77%,能較準確地識別城市河道水質類型。