松嫩平原水體葉綠素a、懸浮物及鹽分濃度遙感反演研究

馬 馳

（遼寧省交通高等?？茖W校，沈陽 110122）

水體中葉綠素、懸浮物以及鹽分是內陸水體水質的重要指標[1,2]。近年來，隨著全球氣候的變化、人類活動的加劇、農牧業(yè)生產者不斷向耕地牧場施加化肥與農藥，造成河流湖泊等水體水質下降，主要體現(xiàn)在水體葉綠素、懸浮物及鹽分的濃度升高，嚴重影響了區(qū)域生態(tài)環(huán)境和經濟的可持續(xù)發(fā)展[3-5]。諸多研究表明，水體的光學特征控制著傳感器所成遙感影像的反射率，而葉綠素、懸浮物及鹽分的濃度直接影響著水體的光學特征，從而決定了遙感影像反射率的光譜特征[6,7]。由于遙感技術具有影像覆蓋范圍廣、數(shù)據(jù)獲取經濟方便、影像數(shù)據(jù)更新快等優(yōu)勢，因此，利用遙感技術實時監(jiān)測水體水質潛力巨大。

近年來，國內外學者相繼展開了對水體水質的遙感研究。Dekker 等[8]在利用Landsat TM 遙感數(shù)據(jù)研究荷蘭的Loosdrecht湖水水質時發(fā)現(xiàn)，Landsat TM 遙感影像第2 波段反射率與Loosdrecht 湖水體葉綠素a和總懸浮物濃度均具有較好相關性；周德明等[9]以太湖為研究對象，利用Landsat TM 遙感影像估測水體葉綠素a 和水體懸浮物濃度，發(fā)現(xiàn)波段反射率組合(TM1+TM3)/TM1/TM3與葉綠素a 相關性最好，TM2+TM3與懸浮物相關性最好；汪雨豪等[10]利用高分2號（GF-2）遙感數(shù)據(jù)研究蘇州市區(qū)水體水質時發(fā)現(xiàn)，以近紅外波段、紅光波段估測水體葉綠素a濃度精度較好，以綠波段、近紅外波段、紅光波段估測水體懸浮物濃度精度較好；孫德勇等[11]使用ASD Field Spec Pro野外光譜輻射儀測量太湖水體高光譜數(shù)據(jù)，輔以水體葉綠素a濃度的實驗室化驗值，采用光學分類的方法建立3個類別水體的反演模型，研究結果顯示，四波段模型適合第一類水體，一階微分模型適合第二、三類水體；徐祎凡等[12]利用環(huán)境災害小衛(wèi)星高光譜影像（HJ-1A HSI）研究太湖水體營養(yǎng)化狀態(tài)，經過相關性分析發(fā)現(xiàn)，利用第73、76、80 三個波段建立的葉綠素a濃度的反演模型決定系數(shù)R2達到0.862 6；董舜丹[13]等比較了Landsat 8 與Sentinel-2 兩種多光譜遙感影像對香港近海葉綠素a濃度的反演能力，結果顯示：兩種遙感影像對于反演香港近海海域葉綠素a濃度均具有較強的實用性，且Landsat 8遙感影像的反演精度略高于Sentinel-2。

綜上所述，現(xiàn)階段針對水體水質遙感方面的研究多以單個河流、湖泊等小區(qū)域為研究對象，常用的遙感影像數(shù)據(jù)中，Landsat TM、ETM 傳感器已經停止接收數(shù)據(jù)，無法保證其現(xiàn)勢性；EO-1、GF-2 等商業(yè)數(shù)據(jù)費用較高，限制了其使用范圍；ASD Field Spec Pro 實測光譜數(shù)據(jù)、HJ-1A HSI 等高光譜數(shù)據(jù)，由于其光譜波段太多，數(shù)據(jù)冗余度較高，數(shù)據(jù)處理復雜。哨兵2A （Sentinel-2A）遙感衛(wèi)星于2015年升空，其遙感影像免費提供給世界用戶使用，與Landsat TM 等多光譜遙感影像相比較，Sentinel-2A 遙感影像具有更高的空間分辨率、更精細的波譜分辨率以及更短的重訪周期，必將在地表地物識別等方面發(fā)揮巨大潛力，但利用Sentinel-2A 遙感數(shù)據(jù)估測水體水質等方面的研究仍處于試驗階段。松嫩平原是我國重要的商品糧基地，上世紀中后期以來，受自然條件、地質地貌、人為因素等多方面影響，松嫩平原生態(tài)環(huán)境趨向惡化，主要表現(xiàn)為土壤沙化、鹽堿化嚴重，河流湖泊等水體的葉綠素、懸浮物濃度升高。為此，本文試驗利用Sentinel-2A 遙感影像反演松嫩平原河流、湖泊等水體葉綠素a、懸浮物、鹽分等3 種水質參數(shù)的濃度，為松嫩平原區(qū)域水體水質監(jiān)測提供數(shù)據(jù)支持，為Sentinel-2A 遙感數(shù)據(jù)在水體水質等方面的研究提供參考。

1 材料與方法

1.1 研究區(qū)概況

松嫩平原位于東北平原中部，主要包括黑龍江省中南部、吉林省大部、遼寧省北部，由松花江、嫩江沖積而成，面積約10 萬km2。松嫩平原東部、北部、西部地區(qū)為山前臺地地貌，平均海拔150～300 m，中部、南部地區(qū)為沖積平原地貌，地勢開闊平坦，多河流、湖泊、沼澤，平均海拔110～180 m。

1.2 研究方法與技術路線

利用Sentinel-2A 遙感數(shù)據(jù)進行松嫩平原水質遙感監(jiān)測的技術路線如圖1所示。

圖1 技術路線Fig.1 Technical route

1.3 數(shù)據(jù)獲取與處理

2018年5月5日至5月11日在研究區(qū)內沿交通干線選擇中型以上水庫、湖泊以及嫩江干、支流進行水體采樣，采得嫩江干、支流水樣9 個，天然湖泊水樣15 個，人工水庫水樣8個，共計水樣32 個，同時測量采樣點坐標，采樣點分布見圖2。采樣后在實驗室化驗水樣的葉綠素a、懸浮物、鹽分等3種水質參數(shù)濃度。其中，葉綠素a 濃度的測定采用分光光度法；懸浮物濃度采用稱重法，即：利用玻璃纖維濾膜過濾水樣，以60 ℃烘烤24 h 去除水分，以0.001 g 天平稱重獲得懸浮物重量，再除以樣品體積獲得懸浮物濃度；鹽分濃度測定采用滴定法，測算每升水樣中鹽分陰、陽離子的質量，將鹽分陰、陽離子求和從而計算出水樣鹽分濃度。

圖2 水體采樣點分布Fig.2 Distribution of water sampling points

為了保證水體采樣時間與遙感影像獲取時相具有一致性，本試驗選擇成像時間為2018年5月7日的覆蓋松嫩平原全境的sentinel-2A 遙感影像共12 景。遙感影像的預處理主要包括影像的輻射校正、大氣校正、幾何校正以及影像的裁剪與拼接，數(shù)據(jù)預處理工作采用ENVI 5.3軟件進行。

為了在遙感影像中清晰表現(xiàn)出松嫩平原水體信息，本試驗提取了sentinel-2A 遙感影像中的河流、湖泊等水體，即：通過分析水體、城鎮(zhèn)、裸土等地表地物在sentinel-2A 遙感影像中的光譜特征，利用ENVI 5.3 軟件采用決策樹分類的方法從選取的12 景Sentinel-2A 遙感影像中提取出河流、湖泊等水體。

1.4 水體葉綠素a、懸浮物、鹽分濃度與水體反射率的相關性分析

將采得的32 個水樣分成建模樣本和檢驗樣本兩部分：隨機選取22 個水樣作為建模樣本，剩余的10 個水樣作為模型檢驗樣本。將水樣葉綠素a、懸浮物、鹽分3 種水質參數(shù)的濃度與sentinel-2A 遙感影像反射率逐波段進行相關性分析[見式（1）]，獲取松嫩平原水體葉綠素a、懸浮物、鹽分濃度的敏感波段，用以建立3種水質參數(shù)的反演模型。

式中：r為3 種水質參數(shù)的濃度與sentinel-2A 遙感影像反射率的相關系數(shù)；xi與x-為3 種水質參數(shù)濃度的實測值與平均值；yi與y-為研究區(qū)水體反射率的實測值與平均值。

諸多學者的研究表明[3,14,15]，將遙感影像反射率進行適當?shù)牟ǘ谓M合或數(shù)學變換，可以有效抑制遙感影像中噪聲影響，提高遙感影像反射率與水體葉綠素a、懸浮物、鹽分的相關性。本試驗引入歸一化植被指數(shù)、懸浮物指數(shù)、鹽分指數(shù)等多種波段組合，并與研究區(qū)水體水質參數(shù)濃度進行相關性分析，選取相關性最好的波段組合作為反演因子建立研究區(qū)3種水質參數(shù)的反演模型。

1.5 模型建立及模型精度評定

以sentinel-2A 遙感影像反射率為自變量，研究區(qū)采樣水體葉綠素a、懸浮物、鹽分等水質參數(shù)濃度為因變量，利用SPSS軟件采用回歸分析的方法建立研究區(qū)水體葉綠素a、懸浮物、鹽分濃度的反演模型。模型的精度采用決定系數(shù)（Coefficient of determination，R2）和均方根誤差（Root mean square error，RMSE）予以評定，選取R2較大、RMSE較小的模型作為研究區(qū)水質參數(shù)的反演模型。

式中：yi為水質參數(shù)濃度的實測值；為水質參數(shù)濃度的平均值；為水質參數(shù)濃度的反演值；n為建模樣本個數(shù)。

2 光譜特征分析與建模

水體葉綠素a、懸浮物、鹽分等水質參數(shù)濃度直接影響水體的光譜特征，并分別反映于不同的光譜波段[16,17]。因此，通過分析葉綠素a、懸浮物、鹽分濃度與水體光譜之間的關系，能夠定量反演出3 種水質參數(shù)的濃度。本試驗在分析3 個水質參數(shù)濃度與sentinel-2A 反射率相關性的基礎上，選取葉綠素a、懸浮物、鹽分的敏感波段，建立3種水質參數(shù)的反演模型。

2.1 葉綠素a濃度的反演

2.1.1 相關性分析

總結前人經驗，現(xiàn)階段常用的水質參數(shù)遙感反演的波段類型主要包括：單波段、雙波段求商、多波段組合、比值植被指數(shù)以及歸一化植被指數(shù)等。本試驗將建模水樣葉綠素a濃度的化驗值與水樣對應的sentinel-2A 各波段反射率進行相關性分析，分析結果如表1所示。從相關性分析結果看出，sentinel-2A 第8 波段與松嫩平原水體葉綠素a 濃度相關性最好，相關系數(shù)達到0.688，第10、11、12 等3 個波段與葉綠素a 濃度相關性較差，相關系數(shù)均未達到0.5。雙波段求商中，藍波段、深藍波段與綠波段的反射率求商，與研究區(qū)水體葉綠素a 濃度的相關性較好，分別達到-0.721 和-0.785；多波段組合中，近紅外（第8波段）、窄近紅外（第8A 波段）與紅波段的組合（歸一化植被指數(shù)）及與藍波段、深藍波段的組合，與研究區(qū)葉綠素a 濃度相關性較好，其中，窄近紅外（第8A波段）與深藍波段（第1 波段）的組合(B8A-B1)/(B8A+B1)與葉綠素a 的濃度相關系數(shù)為0.835，窄近紅外（第8A 波段）與紅波段（第4 波段）的組合(B8A-B4)/(B8A+B4)與葉綠素a 的濃度相關系數(shù)達到了0.862。

表1 波段反射率及其組合與葉綠素a濃度的相關系數(shù)Tab.1 Correlation coefficient between band reflectance and combination and chlorophyll a concentration

2.1.2 葉綠素a濃度反演模型的建立與精度評價

參考相關性分析結果，以sentinel-2A 影像多波段組合(B8A-B4)/(B8A+B4)為反演因子，Pearson 為相關系數(shù)，反演因子與水體葉綠素a 濃度的相關性分析結果顯示，其顯著性p＜0.001，表明在超過99%的置信區(qū)間內兩者顯著相關。故將(B8A-B4)/(B8A+B4)作為反演松嫩平原水體葉綠素a 濃度的特征因子，建立研究區(qū)水體葉綠素a濃度的一元一次、一元二次、對數(shù)、倒數(shù)、指數(shù)、指數(shù)、冪等反演模型。

利用模型決定系數(shù)R2和模型均方根誤差RMSE對葉綠素a濃度反演模型進行精度評價。表2顯示，利用一元二次建模方法建立的葉綠素a 濃度反演模型Y=15.13X+19.606X2+5.849，模型決定系數(shù)最大，為R2=0.887，表明其模型精度高于其他4種反演模型，均方根誤差最小，為RMSE=0.523 μg/L，表明一元二次模型的葉綠素a濃度反演值與實測值的偏差均小于其他4種模型，故認為利用一元二次建模方法建立的松嫩平原水體葉綠素a濃度反演模型Y=15.13X+19.606X2+5.849最優(yōu)。

表2 葉綠素a濃度反演模型及精度對比Tab.2 Inversion model and precision comparison of chlorophyll a concentration

2.2 水體懸浮物濃度的反演

2.2.1 相關性分析

松嫩平原水體懸浮物濃度與Sentinel-2A 影像反射率相關性分析結果顯示，研究區(qū)水體懸浮物濃度與水體反射率呈正相關性，且第5 波段相關性最好，為r=0.664，其次為第4 波段，相關系數(shù)為0.651，第10 波段相關性最差，僅為0.104，未達到顯著水平（表3）。

表3 波段反射率及其組合與懸浮物濃度的相關系數(shù)Tab.3 Correlation coefficient between band reflectance and combination and suspended solids concentration

雙波段求商、多波段組合與研究區(qū)水體懸浮物濃度相關性分析結果顯示（圖3），第5波段與第1波段的商B5/B1與研究區(qū)水體懸浮物濃度的相關性最好，為r=0.751，其次為第6 與第2 波段求商B6/B2，與懸浮物濃度相關系數(shù)為0.714，再次為第8 與第9 波段求商B9/B8，相關系數(shù)為-0.709。多波段組合中，第8 波段與第1 波段的組合(B8-B1)/(B8+B1)與研究區(qū)水體懸浮物濃度相關性最好，為r=0.746，其次為第3、第5波段的組合(B5-B3)/(B5+B3)，與懸浮物濃度的相關系數(shù)為r=0.729。

圖3 波段反射率組合與懸浮物相關性分析Fig.3 Correlation analysis between band reflectance combination and suspended solids

2.2.2 水體懸浮物濃度反演模型的建立與精度評價

依據(jù)反射率的單波段、雙波段求商、多波段組合與松嫩平原水體懸浮物濃度相關性分析結果，選取雙波段求商B5/B1作為研究區(qū)水體懸浮物濃度反演的特征因子，建立懸浮物濃度反演的一元一次、一元二次、倒數(shù)、對數(shù)、指數(shù)、冪等反演模型（表4）。

表4 懸浮物濃度反演模型及精度對比Tab.4 Suspended solids concentration estimation model and precision comparison

建模結果顯示，利用一元二次建模方法建立的研究區(qū)水體懸浮物濃度反演模型Y=-0.205X+0.403X2+0.191，模型決定系數(shù)最大，為R2=0.765，表明該模型精度高于其他5 種模型，模型的均方根誤差為RMSE=0.074 g/L，均小于其他5個模型的均方根誤差。因此，將一元二次模型Y=-0.205X+0.403X2+0.191 作為研究區(qū)水體懸浮物濃度的最優(yōu)反演模型，并將其應用于松嫩平原水體懸浮物濃度制圖。

2.3 水體鹽分濃度的反演

2.3.1 相關性分析

松嫩平原水體鹽分濃度與sentinel-2A 各波段反射率的相關性分析結果表明（表5），第4 至第9 波段的相關性大于0.5，其中，第9 波段反射率與水體鹽分的相關性最好，達到r=0.780，其次為第8A 波段，相關系數(shù)為r=0.672。第一波段反射率與水體鹽分濃度相關性最差，僅為r=0.313。

表5 波段反射率及其組合與鹽分濃度的相關系數(shù)Tab.5 Correlation coefficient between band reflectance and combination and salt concentration

為了探究波段組合對遙感影像噪聲的削弱效果，參考前人研究成果，本試驗引入多種鹽分指數(shù)，包括(R-NIR)/(R+NIR)、B/R、(B-R)/(B+R)、(B×R)/B、(B×R)/G、(B×NIR)/G等。多波段組合與研究區(qū)水體鹽分濃度的相關性分析結果顯示，(B2×B8)/B3鹽分指數(shù)與水體鹽分濃度相關性最好，達到r=0.876，其次為(B8-B4)/(B8+B4)鹽分指數(shù)，與水體鹽分濃度的相關系數(shù)為r=0.819。B2/B3與水體鹽分濃度相關性較差，僅為r=-0.615。

2.3.2 水體鹽分濃度反演模型的建立與精度評價

選擇鹽分指數(shù)(B2×B8)/B3作為松嫩平原水體鹽分濃度的特征因子，以倒數(shù)、對數(shù)、指數(shù)、一元線性、一元二次、冪等方法建立研究區(qū)水體鹽分濃度的反演模型（表6）。建模結果表明，以(B2×B8)/B3為特征因子建立的一元二次反演模型Y=-5.592X+29.841X2+0.728，模型決定系數(shù)均大于其他5 個模型，達到R2=0.848，相應的，該模型的均方根誤差最小，為RMSE=0.385 g/L，故將一元二次反演模型Y=-5.592X+29.841X2+0.728作為研究區(qū)水體鹽分濃度的最優(yōu)反演模型，并應用于制圖。

表6 鹽分濃度反演模型及精度對比Tab.6 Salt concentration estimation model and precision comparison

2.4 模型精度檢驗

根據(jù)研究區(qū)水體葉綠素a、懸浮物、鹽分濃度的建模結果，選擇3種水質參數(shù)的最優(yōu)反演模型，利用檢驗樣本所對應sentinel-2A 影像的反射率計算葉綠素a、懸浮物、鹽分濃度的反演值，并利用相對誤差[式（4）]對3種水質參數(shù)反演模型進行精度檢驗（見圖4）。

圖4 水質參數(shù)相對誤差散點圖Fig.4 Scatter diagram of relative error of water quality parameters

水體葉綠素a、懸浮物、鹽分的相對誤差分析結果顯示，分別有8個葉綠素a檢驗樣本、7個懸浮物檢驗樣本、7個鹽分檢驗樣本的相對誤差處于-20%至+20%之間，分別占總檢驗樣本的80%、70%、70%。其中，水體葉綠素a 濃度的反演值與實測值相對誤差的最大值為48.4%，最小值為6.2%；懸浮物濃度的反演值與實測值相對誤差的最大值為41.1%，最小值為7.7%；水體鹽分濃度的反演值與實測值相對誤差的最大值為46.3%，最小值為1.1%，表明本試驗選定的3個水質參數(shù)反演模型具有足夠的精度和穩(wěn)定性，能夠勝任松嫩平原水體3個水質參數(shù)的反演工作。

2.5 研究區(qū)水體水質制圖

參照建模結果，利用松嫩平原水質參數(shù)的3個最優(yōu)反演模型反演研究區(qū)水體葉綠素a、懸浮物、鹽分濃度并制圖（圖5～圖7）。水體葉綠素a 濃度制圖結果顯示，研究區(qū)北部的連環(huán)湖、大慶水庫、青肯泡，中部的向陽湖、嫩江一帶，水體葉綠素a 濃度普遍低于3 μm/L，西部的洋沙泡、中部的南引水庫、南部的花敖泡、工農湖等水域葉綠素a 濃度普遍高于5 μm/L；懸浮物濃度制圖結果顯示，研究區(qū)中部的月亮湖、南引水庫一帶，南部的查干湖、花敖泡、工農湖一帶，西部的洋沙泡、頭道泡一帶水體懸浮物濃度較高；水體鹽分濃度制圖結果顯示，研究區(qū)北部的龍虎泡、西部的洋沙泡一帶、中部的月亮湖與南引水庫一帶、南部的查干湖與工農湖一帶，水體鹽分濃度普遍高于0.5 g/L。

圖5 水體葉綠素a濃度反演結果圖Fig.5 Estimation map of chlorophyll a concentration in water

圖6 水體懸浮物濃度估反演果圖Fig.6 Estimation map of suspended solids concentration in water

圖7 水體鹽分濃度反演結果圖Fig.7 Salt concentration estimation map of water body

3 討 論

水體成分及其濃度顯著影響著水體的光譜特征，并直接體現(xiàn)在遙感影像的水體反射率上。相關研究表明，水體葉綠素a 濃度升高后會使水體的光譜曲線在近紅外波段出現(xiàn)反射峰。本試驗中，sentinel-2A的近紅外光譜區(qū)域的第8、8A兩個波段反射率與葉綠素a濃度相關性最好，與前人的研究結果相同[10,18,19]；純凈水在藍、綠光譜范圍存在較強反射，在紅波范圍存在吸收谷。但是，當水體中的懸浮物濃度增加時，水體顏色變亮，反射峰向紅波區(qū)域移動，從而將提升懸浮物濃度與遙感影像紅波區(qū)域反射率的相關性[20]，本試驗研究結果顯示，松嫩平原水體懸浮物濃度與sentinel-2A 影像反射率的相關系數(shù)峰值出現(xiàn)于紅邊1（第5）波段；隨著水體中鹽分濃度的增加而產生黃色物質，引起水體在可見光與近紅外光譜范圍出現(xiàn)較強的反射，進而提升與水體鹽分濃度的相關性[21]，本試驗中研究區(qū)水體鹽分濃度與sentinel-2A 影像反射率的相關性在第8A、第9波段較好。

將波段反射率進行波段組合、波段求商等變換，尤其在引入某些有針對性的指數(shù)（如歸一化植被指數(shù)、懸浮物指數(shù)、鹽分指數(shù)等）以后，可以有效削弱遙感影像中噪聲對成像質量的影響，顯著提升與3個水質參數(shù)之間的相關性。經過波段組合以后，研究區(qū)水體葉綠素a濃度與sentinel-2A遙感影像反射率的相關系數(shù)從0.688 提升至0.862，水體懸浮物濃度與sentinel-2A 遙感影像反射率的相關系數(shù)從0.664 提升至0.751，水體鹽分濃度與遙感影像反射率的相關系數(shù)從0.780 提升至0.876。本試驗利用回歸分析的方法建立了松嫩平原水體葉綠素a、懸浮物、鹽分3 個水質參數(shù)的一元二次反演模型，模型的判定系數(shù)分別達到0.887、0.765、0.848，用該模型反演松嫩平原3個水質參數(shù)濃度，獲得了較好的反演效果，其原因總結如下：①本試驗所選取的遙感影像成像時間與水體采樣時間同步，影像可以真實反映采樣時刻研究區(qū)水體葉綠素a、懸浮物及鹽分等水質參數(shù)的光譜特征；②本試驗將遙感影像進行了輻射校正與大氣校正，可以有效削弱影像成像過程中噪聲的影響，提高了影像反射率與水體葉綠素a、懸浮物及鹽分等水質參數(shù)的相關性；③構建了歸一化植被指數(shù)、懸浮物指數(shù)、鹽分指數(shù)等適當?shù)牟ǘ谓M合后，顯著提高了sentinel-2A 遙感影像對研究區(qū)水體葉綠素a、懸浮物及鹽分等水質參數(shù)的敏感度，從而改善了水質參數(shù)反演模型的反演精度。

4 結 論

本試驗利用sentinel-2A 遙感影像結合研究區(qū)水體采樣的水質化驗數(shù)據(jù)，反演了松嫩平原水體葉綠素a、懸浮物、鹽分濃度，并得到以下結論：

（1）松嫩平原水體葉綠素a、懸浮物、鹽分濃度與sentinel-2A 遙感影像的反射率均呈正相關性，相關系數(shù)的峰值分別出現(xiàn)于第8、第5和第9波段。

（2）通過引入歸一化植被指數(shù)、懸浮物指數(shù)、鹽分指數(shù)等波段組合后，可以有效改善與水體葉綠素a、懸浮物、鹽分濃度的相關性。

（3）以一元二次建模方法建立的研究區(qū)3種水質參數(shù)反演模型精度較高，可以很好的反演松嫩平原水體葉綠素a、懸浮物、鹽分濃度。制圖結果顯示，松嫩平原的連環(huán)湖、向陽湖、大慶水庫及嫩江干支流水域的葉綠素a、懸浮物及鹽分濃度較低，水質較好。