珠江口河網(wǎng)水體非光學(xué)活性水質(zhì)參數(shù)高光譜反演

蔡建楠 , 劉海龍 , 姜波 , 陳吟暉 , 李杰鴻 , 吳思曉 , 梁建霞 , 黃華, 邢前國

1. 中山市環(huán)境監(jiān)測站, 廣東 中山 528403;

2. 中山市生態(tài)環(huán)境局, 廣東 中山 528403;

3. 中國科學(xué)院煙臺(tái)海岸帶研究所, 中國科學(xué)院海岸帶環(huán)境過程與生態(tài)修復(fù)重點(diǎn)實(shí)驗(yàn)室, 山東 煙臺(tái) 264003;

4. 廣東環(huán)境保護(hù)工程職業(yè)學(xué)院, 廣東省無人機(jī)環(huán)保應(yīng)用工程技術(shù)研究中心, 廣東 佛山 528216

珠江口地區(qū)河網(wǎng)密布, 區(qū)域內(nèi)河流眾多, 隨著城市化過程的不斷推進(jìn), 水質(zhì)問題越來越突出, 直接影響到居民生活質(zhì)量, 與粵港澳大灣區(qū)生態(tài)環(huán)境建設(shè)的目標(biāo)不符。為此, 珠江口各地投入了巨大的社會(huì)資源, 廣泛開展了河流整治工作, 其中長期、密集的河流水質(zhì)監(jiān)測工作對于城市未達(dá)標(biāo)水體的精準(zhǔn)識別和對河流整治效果的綜合評估具有重要意義(王曉 等, 2017; 張雪 等, 2019)。但傳統(tǒng)的水質(zhì)監(jiān)測分析方法存在成本高昂, 采樣及分析過程耗時(shí)長, 且僅代表特定區(qū)域和時(shí)段等缺點(diǎn), 監(jiān)測結(jié)果往往難以準(zhǔn)確反映監(jiān)測對象整體的水質(zhì)狀況, 難以滿足當(dāng)前形勢下水環(huán)境質(zhì)量評價(jià)的需要(Wang et al, 2004)。

高光譜水質(zhì)遙感監(jiān)測技術(shù)主要基于水體實(shí)測光譜, 通過建立光譜特征和水質(zhì)指標(biāo)的關(guān)系模型, 實(shí)現(xiàn)對水質(zhì)指標(biāo)的反演(王桂芬 等, 2010; 彭建 等, 2018)。高光譜遙感數(shù)據(jù)的光譜譜段豐富, 在水質(zhì)參數(shù)定量化監(jiān)測中具有良好的應(yīng)用前景, 可作為傳統(tǒng)水質(zhì)監(jiān)測方法的有效補(bǔ)充, 近年來得到迅速的發(fā)展(林劍遠(yuǎn) 等, 2019), 其中濁度、懸浮物、葉綠素a 等光學(xué)活性水質(zhì)參數(shù)由于具有機(jī)理明確的光學(xué)響應(yīng)特征而得到廣泛的研究(Wang et al, 2019)。但水環(huán)境管理與整治工作中主要監(jiān)控的指標(biāo)為碳組分、營養(yǎng)鹽等非光學(xué)活性水質(zhì)參數(shù), 該類參數(shù)具有光學(xué)響應(yīng)波段不顯著和復(fù)雜多樣的特點(diǎn), 因而開展該類參數(shù)的光譜特征研究及水質(zhì)反演具有較高的科學(xué)意義和現(xiàn)實(shí)意義(邢前國, 2007)。近年來, 國內(nèi)外學(xué)者采取了間接反演、半經(jīng)驗(yàn)?zāi)Ｐ偷确椒ㄩ_展了水體中總磷(TP)、總氮(TN)、溶解性有機(jī)碳(DOC)、氨氮(NH3-N)等非光學(xué)活性水質(zhì)的高光譜反演研究(Song et al, 2012; Sudduth et al, 2015; Fichot et al, 2016; 張海威 等, 2017), 然而相關(guān)研究主要集中于湖泊水庫, 應(yīng)用于河流特別是水質(zhì)復(fù)雜河網(wǎng)地區(qū)的研究報(bào)道較少。本文通過珠江口中山市各類型河流的高光譜監(jiān)測和同步水質(zhì)分析, 研究受測水體的高光譜特征及其與河流化學(xué)需氧量(CODCr)和總磷(TP)這兩項(xiàng)水質(zhì)指標(biāo)濃度的關(guān)系, 并嘗試建立高光譜數(shù)據(jù)與水質(zhì)指標(biāo)濃度的反演模型, 進(jìn)而為河流非光學(xué)活性水質(zhì)參數(shù)高光譜監(jiān)測的實(shí)施及相關(guān)高光譜傳感器的設(shè)計(jì)提供參考。

1 研究區(qū)概況

中山市位于珠江口西岸, 地處粵港澳大灣區(qū)的幾何中心, 北連廣州, 毗鄰港澳, 陸域總面積約1783km2, 多年平均降雨量為1748.7mm, 多年平均水資源總量17.38 億m3。中山市轄區(qū)內(nèi)河流眾多, 除磨刀門水道、雞鴉水道、小欖水道、橫門水道和石岐河五大主要河流外, 還有各種規(guī)模的河流超過1000 條(圖1)。中山市經(jīng)濟(jì)發(fā)達(dá), 產(chǎn)業(yè)集群密布, 近年來隨著城市的快速發(fā)展和人口的不斷增加, 部分河流水質(zhì)受到了一定的污染。

圖1 研究區(qū)河流水系示意圖 Fig.1 Schematic diagram of river system in the study area

2 數(shù)據(jù)與方法

2.1 光譜測量和水樣分析

2019 年7 月至9 月在中山市實(shí)施了多次河流現(xiàn)場光譜測量。測定點(diǎn)位的選取以圖1 所示的中山市主要河流及其支流為主要對象, 兼顧考慮入河污染特征、流經(jīng)區(qū)域類型、周邊環(huán)境特點(diǎn)的多樣性和代表性, 以盡可能覆蓋中山市各種水質(zhì)類型的河流, 最終獲得151 個(gè)現(xiàn)場實(shí)測光譜數(shù)據(jù)和水質(zhì)樣品(代表性河流照片見圖2)。使用的測量設(shè)備為美國Ocean Optics 公司的USB4000 光譜儀, 測量光譜范圍為360～1000nm, 光譜分辨率<0.2nm。測量時(shí)選擇晴天無云或少云天氣狀況下距正午3h 的時(shí)段內(nèi)進(jìn)行。在每個(gè)點(diǎn)位離河流水面高度0.5m 處, 測量河流水體的上行輻射亮度(Lu)、天空輻射亮度(Ls)和太陽下行輻照度(Ed)。每個(gè)點(diǎn)位測量前對儀器進(jìn)行暗電流校正, 并使用標(biāo)準(zhǔn)灰板進(jìn)行調(diào)試, 使其信號基本保持在80%量程范圍內(nèi)。進(jìn)行3 次測量后取平均值作為光 譜測定結(jié)果, 利用所測數(shù)據(jù)代入公式(1), 即可計(jì)算測點(diǎn)的光譜反射率R(無量綱)。

式中: ρ為天空光在水體表面的菲涅爾反射率, 通常取經(jīng)驗(yàn)值0.028(Mobley, 1999)。

監(jiān)測點(diǎn)的水樣采集與光譜測量同步進(jìn)行, 水樣分析項(xiàng)目為CODCr、TP。水樣取河流表層0～20cm處的混合樣, 按照相應(yīng)的國家技術(shù)規(guī)范要求保存并送至實(shí)驗(yàn)室進(jìn)行分析。實(shí)驗(yàn)室分析過程采用平行樣、標(biāo)準(zhǔn)樣品、加標(biāo)回收等質(zhì)量控制方法, 以確保分析的準(zhǔn)確度和精密度。全部樣品測得的水質(zhì)指標(biāo)特征描述見表1。

圖2 代表性監(jiān)測河流照片 黃色箭頭標(biāo)識處為監(jiān)測點(diǎn)位置 Fig.2 Photos of typical rivers (The position of the arrow indicates the monitoring points)

表1 樣品水質(zhì)特征描述 Tab. 1 Description of characteristics of water quality

2.2 高光譜數(shù)據(jù)預(yù)處理

由于本研究獲得的原始高光譜數(shù)據(jù)的波段數(shù)較多, 為提高計(jì)算效率和反演精度, 且有利于與其他研究成果相比較, 數(shù)據(jù)處理過程中僅保留400～900nm 波段范圍內(nèi)的原始數(shù)據(jù), 并將其光譜分辨率重采樣為1nm, 以降低信息冗余, 重采樣后的數(shù)據(jù)將作為本次研究的數(shù)據(jù)。光譜分辨率重采樣采用分段三次 Hermite 插值多項(xiàng)式(Piecewise Cubic Hermite Interpolating Polynomial, PCHIP)法, 該方法可保持原始數(shù)據(jù)的形狀和相應(yīng)的單調(diào)性(Iqbal et al, 2013), 符合本研究后續(xù)數(shù)據(jù)分析需要。

2.3 相關(guān)分析與模型評價(jià)

2.3.1 水質(zhì)指標(biāo)與光譜反射率相關(guān)系數(shù)的計(jì)算

各水質(zhì)指標(biāo)濃度值與水體不同波段高光譜反射率的相關(guān)性分析有助于識別對水質(zhì)指標(biāo)變化敏感的特征波段, 從而進(jìn)一步篩選可用于水質(zhì)反演的最優(yōu)特征波段及其組合。本研究中各水質(zhì)指標(biāo)濃度值與水體高光譜原始反射率相關(guān)系數(shù)的計(jì)算方法為:

2.3.2 反演模型的建立與評價(jià)

以相關(guān)分析所確定的特征波段反射率及其組合作為自變量, 以水質(zhì)指標(biāo)實(shí)測濃度值作為因變量, 建立相關(guān)水質(zhì)反演模型。本研究從151 個(gè)光譜及水質(zhì)樣品中隨機(jī)選擇130 個(gè)樣本作為建模數(shù)據(jù)集, 并計(jì)算回歸模型的決定系數(shù) R2和均方根誤差 RMSE(Root Mean Squared Error), 以評價(jià)擬合精度并確定最優(yōu)模型。在確定最優(yōu)模型后, 利用未參與建模的21 個(gè)樣本作為驗(yàn)證數(shù)據(jù)集, 進(jìn)行估算結(jié)果與實(shí)測值的對比分析, 以絕對誤差 AE(Absolute Error)和平均相對誤差MRE(Mean Relative Error)評價(jià)模型的泛化能力。R2、RMSE、AE 和MRE 的計(jì)算公式分別為:

3 結(jié)果與討論

3.1 河流水體高光譜反射率特征

本研究中的實(shí)測城市河流水體在400～900nm 波段范圍內(nèi)各波段的光譜反射率均值及第10、第90百分位數(shù)如圖3 所示。所測河流水體高光譜總體特征表現(xiàn)為: 400～500nm 范圍內(nèi)的水體光譜反射率總體較低; 隨后水體光譜反射率隨波長增加而上升, 并呈現(xiàn)“三峰二谷”的光譜形態(tài)特征, 即在585nm、710nm、810nm 附近有明顯的反射峰, 在670nm、760nm 附近有明顯的反射低谷。根據(jù)相關(guān)研究, 地表水在585nm、710nm 附近的反射峰主要由后向散射和在該波段處較弱的吸收共同作用所致; 810nm附近的反射峰則是由于該處水分子吸收系數(shù)較小, 水體懸浮顆粒物散射程度加強(qiáng)而進(jìn)一步提高了水體的反射能力(陳雯揚(yáng) 等, 2010); 670nm、760nm 附近的反射低谷分別由藻類和水分子吸收所致(陳楚群 等, 2001; 劉志國 等, 2007; 楊錦坤 等, 2007)。上述分析表明各河流水體光譜反射率主要受懸浮顆粒物和葉綠素a 影響(劉志國 等, 2007)。

圖3 實(shí)測水體高光譜反射率 Fig.3 Hyperspectral reflectance of sampled water

3.2 水質(zhì)指標(biāo)與水體光譜相關(guān)性分析

3.2.1 水質(zhì)指標(biāo)與高光譜原始反射率的相關(guān)性

根據(jù)圖4 顯示, 各水體樣品中CODCr、TP 濃度的自然對數(shù)值與水體光譜反射率的相關(guān)性特征基本趨于一致, 即在400～900nm 波段上總體呈負(fù)相關(guān), 該負(fù)相關(guān)性先增大后減少; 在550～680nm 波段范圍內(nèi)具有較高的相關(guān)性, CODCr、TP 濃度對數(shù)值與光譜反射率相關(guān)性最高, 分別可達(dá)–0.55 和–0.46; 相關(guān)程度在680nm 之后的波段中迅速下降, 在720nm以上波長范圍內(nèi)相關(guān)系數(shù)均不到–0.15。

圖4 不同波段光譜反射率與CODCr 和TP 濃度對數(shù)值的相關(guān)性 Fig.4 Relationships between the logarithm of CODCr, TP concentration and hyperspectral reflectance

3.2.2 水質(zhì)指標(biāo)與光譜反射率比值的相關(guān)性

對任意兩個(gè)波段的原始光譜反射率數(shù)據(jù)進(jìn)行比值計(jì)算, 形成新的光譜變量, 并與水質(zhì)指標(biāo)進(jìn)行回歸分析, 進(jìn)而突出光譜的吸收與反射特征, 這一方法已被廣泛應(yīng)用于水質(zhì)高光譜反演經(jīng)驗(yàn)?zāi)Ｐ偷臉?gòu)建中(Song et al, 2012)。本研究將每個(gè)樣本光譜數(shù)據(jù)中501 個(gè)波段的反射率進(jìn)行兩兩比值計(jì)算, 形成151個(gè)501×501 的比值數(shù)據(jù)矩陣后, 再按照公式(2)分別計(jì)算矩陣中對應(yīng)光譜反射率比值與CODCr和TP 濃度對數(shù)值的相關(guān)系數(shù), 以識別出與水質(zhì)指標(biāo)相關(guān)性最密切的波段組合, 相關(guān)性矩陣計(jì)算結(jié)果見圖5。

CODCr、TP 濃度對數(shù)值與反射率比值的相關(guān)系數(shù)最高分別為0.795 和0.731, 對應(yīng)的波段組合分別為R675/R794、R690/R815。結(jié)合圖3 的光譜形態(tài)分析, 該最優(yōu)波段組合基本上與光譜反射率曲線吸收峰和近紅外反射峰分別對應(yīng)。這表明與直接利用單波段原始數(shù)據(jù)相比, 反射率比值處理更加強(qiáng)化了光譜的差異特征, 有助于進(jìn)一步篩選與CODCr、TP 濃度值相關(guān)性高的波段組合。

根據(jù)圖5 所示, 與CODCr、TP 濃度對數(shù)值相關(guān)性較高的主要是紅波段與近紅外反射率間的比值, 如R675/R794、R690/R815, 其波長與葉綠素的特征吸收、泥沙的反射峰位置一致, 說明非光學(xué)活性參數(shù)CODCr、TP 與光學(xué)敏感性參數(shù)中的葉綠素及懸浮物有一定的相關(guān)性: 由葉綠素表征的浮游植物本身就是CODCr的重要來源, TP 則可能主要以懸浮顆粒形式出現(xiàn), 且浮游植物本身也對懸浮物濃度有貢獻(xiàn)。

圖5 CODCr(a)和TP(b)濃度對數(shù)值與光譜反射率比值的相關(guān)性矩陣 Fig.5 Correlation matrix of relationships between the logarithm of CODCr (a) and TP (b) concentrations and spectral reflection ratio indexes

3.3 基于反射率比值的水質(zhì)反演模型

3.3.1 模型的建立

基于上述反射率比值的篩選結(jié)果, 選擇R675/R794和R690/R815作為CODCr、TP 反演模型的自變量, 分別采用線性、指數(shù)、冪函數(shù)三種擬合方式建立回歸模型。建立的水質(zhì)指標(biāo)光譜反演回歸模型及結(jié)果評價(jià)見表2、圖6。各擬合方程參數(shù)均通過P<0.05 的顯著性檢驗(yàn)。其中, CODCr、TP 最優(yōu)擬合模型為冪函數(shù)模型, 模型決定系數(shù) R2分別達(dá)到0.661 和0.569; 總體上CODCr擬合效果略優(yōu)于TP。

表2 CODCr 和TP 濃度反演模型擬合結(jié)果比較 Tab. 2 Fitting results of CODCr and TP concentration retrieval models

圖6 CODCr(a)、TP(b)濃度最優(yōu)反演模型擬合結(jié)果 Fig.6 Fitting results of the best CODCr (a) and TP (b) concentration retrieval models

3.3.2 模型的精度驗(yàn)證及評價(jià)

基于上述CODCr、TP 濃度高光譜反演最優(yōu)擬合模型, 利用驗(yàn)證數(shù)據(jù)集進(jìn)行模型估算結(jié)果與實(shí)測值的對比分析。結(jié)果(表 3、圖 7)顯示, CODCr、TP 濃度反演模型的估算精度良好, 兩者的MRE 分別為27.2%和32.1%, 表明采用最優(yōu)反射率比值作為自變量構(gòu)建的冪函數(shù)模型具有良好的適用性, 可較好地反演水體中的CODCr、 TP 濃度值。

表3 CODCr 和TP 濃度估算結(jié)果評價(jià) Tab. 3 Evaluation of estimation results for CODCr and TP concentrations

圖7 CODCr(a)、TP(b)濃度實(shí)測值與模型估算值比較 Fig.7 Comparison of measured and estimated values of CODCr (a) and TP (b) concentrations

4 結(jié)論

本文通過對廣東省中山市151 個(gè)河網(wǎng)水體監(jiān)測點(diǎn)的高光譜監(jiān)測和水質(zhì)分析, 獲得了水體的高光譜反射率特征, 并構(gòu)建了CODCr和TP 濃度反演模型, 得到以下結(jié)論:

1) 各河流水體光譜反射率主要受懸浮顆粒物和葉綠素a 影響。在500～680nm 波段范圍內(nèi), CODCr和TP 濃度對數(shù)值與單波段光譜反射率存在較強(qiáng)的 負(fù)相關(guān)關(guān)系。與單波段相比, 特定波段的反射率比值與CODCr、TP 濃度值的相關(guān)性較高, 與CODCr和 TP 濃度值相關(guān)性最大的波段組合分別為R675/R794、R690/R815。

2) 基于上述反射率波段比值組合建立的CODCr和TP 濃度冪函數(shù)反演模型的決定系數(shù)R2分別為 0.661 和 0.569, 估算平均相對誤差分別為27.2%和32.1%, 表明高光譜技術(shù)在珠江口河網(wǎng)水體CODCr、TP 濃度的水質(zhì)反演中具有一定的應(yīng)用潛景。