肖臣稷，王 卿，王 敏，阮俊杰，陳 敏，丁 玲，黃沈發(fā),

(1.東華大學(xué) 環(huán)境科學(xué)與工程學(xué)院，上海 201620；2.上海市環(huán)境科學(xué)研究院，上海 200233；3.中國(guó)長(zhǎng)江三峽集團(tuán)有限公司，北京 100038；4.上?？睖y(cè)設(shè)計(jì)研究院有限公司，上海 200335)

水體葉綠素a的濃度是表征水體富營(yíng)養(yǎng)化水平的重要指標(biāo)之一。隨著衛(wèi)星遙感技術(shù)的發(fā)展，基于多光譜或高光譜衛(wèi)星遙感影像開(kāi)展水體葉綠素a濃度的遙感反演已取得顯著進(jìn)展，并發(fā)展出標(biāo)準(zhǔn)藍(lán)綠比值法[1]、紅外比值法[2]、歸一化葉綠素指數(shù)(NDCI)法[3]、三波段模型[4]、四波段模型[5]、正交函數(shù)分析[6]等算法。這些反演模型中，半經(jīng)驗(yàn)方法占據(jù)主導(dǎo)地位，經(jīng)驗(yàn)方法次之，分析法在應(yīng)用中使用得較少[7]。

高光譜遙感技術(shù)的快速發(fā)展使得高光譜遙感影像在水色反演中得到越來(lái)越多的應(yīng)用。研究表明，基于高光譜遙感影像的水色反演精度比基于多光譜遙感影像的高，這是因?yàn)槎喙庾V影像波段少且間距寬，難以準(zhǔn)確選擇敏感波段。開(kāi)展高光譜遙感水色反演則可實(shí)現(xiàn)精確選用對(duì)反射率敏感的特定波段的算法，如特征峰[8]、生物光學(xué)模型[9]等，以及基于連續(xù)光譜的一階導(dǎo)數(shù)光譜[10]等算法，還可通過(guò)迭代優(yōu)化[11]等方法進(jìn)行處理。

Kratzer等[12]研究表明，使用衛(wèi)星遙感進(jìn)行葉綠素a濃度的反演對(duì)傳感器的空間分辨率有著明確的要求。高分五號(hào)衛(wèi)星是我國(guó)于2018年5月9日發(fā)射的首顆高光譜綜合觀測(cè)衛(wèi)星，運(yùn)行于太陽(yáng)同步軌道。其光譜分辨率既高又全，在0.4～2.5 μm內(nèi)具有330個(gè)波段，光譜分辨率最高可達(dá)0.03 cm-1，光譜定標(biāo)精度最高可達(dá)0.008 cm-1，還具有30 m的空間分辨率與60 km的幅寬[13]。這些性能足以滿足對(duì)葉綠素a濃度進(jìn)行反演的基本要求，并且韋安娜等[14]已經(jīng)應(yīng)用高分五號(hào)衛(wèi)星的遙感數(shù)據(jù)在內(nèi)陸水體中對(duì)葉綠素a濃度進(jìn)行反演。

長(zhǎng)江河口及其鄰近海域是我國(guó)富營(yíng)養(yǎng)化水平最高的海域，其中長(zhǎng)江口的富營(yíng)養(yǎng)化指數(shù)為5～6[15]，呈現(xiàn)出中重度的富營(yíng)養(yǎng)化[16]。這是因?yàn)殚L(zhǎng)三角城市群的興起使得人口劇增，導(dǎo)致城市生活污水、工農(nóng)業(yè)污染的排放量不斷增加，大量污染物被排放到長(zhǎng)江中，致使生態(tài)環(huán)境惡化[17]?；诟叻治逄?hào)衛(wèi)星高光譜影像提取長(zhǎng)江河口水體的表面反射率，結(jié)合實(shí)際采樣數(shù)據(jù)進(jìn)行數(shù)學(xué)建模，遙感監(jiān)測(cè)研究長(zhǎng)江河口的水質(zhì)，旨在為高分五號(hào)衛(wèi)星多光譜遙感的應(yīng)用找到方法，為高分系列衛(wèi)星的高光譜遙感在更大時(shí)空尺度中的推廣應(yīng)用提供案例，以及為長(zhǎng)江口環(huán)境監(jiān)測(cè)與保護(hù)提供技術(shù)支撐。

1 實(shí)測(cè)與反演研究方法

1.1 研究區(qū)域點(diǎn)位選取與葉綠素a濃度實(shí)測(cè)

長(zhǎng)江河口位于上海市北部，起自江蘇江陰鵝鼻嘴，止于入?？诘碾u骨礁，長(zhǎng)約232 km。其經(jīng)緯度為東經(jīng)120.2°～122.3°、北緯30.8°～32.2°，包含崇明島、長(zhǎng)興島、橫沙島等島嶼與多個(gè)沙洲。區(qū)域內(nèi)徑流大、輸沙量多、余流速度較大，物理化學(xué)環(huán)境受潮汐作用及上游來(lái)水等水文條件交匯混合的復(fù)雜影響。

由于長(zhǎng)江口水動(dòng)力較強(qiáng)，水色變化相對(duì)較快，為保證遙感數(shù)據(jù)與現(xiàn)場(chǎng)數(shù)據(jù)的一致性，現(xiàn)場(chǎng)監(jiān)測(cè)在衛(wèi)星經(jīng)過(guò)時(shí)同步開(kāi)展。通過(guò)查詢衛(wèi)星軌道數(shù)據(jù)，計(jì)算得出衛(wèi)星過(guò)境時(shí)間為2019年5月24日13時(shí)，衛(wèi)星過(guò)境時(shí)刻及其前后5 min各進(jìn)行一次現(xiàn)場(chǎng)采樣，根據(jù)HJ 897—2017《水質(zhì) 葉綠素a 的測(cè)定 分光光度法》進(jìn)行檢測(cè)，并將3組數(shù)據(jù)的平均值作為測(cè)量值。由于長(zhǎng)江河口的水域可達(dá)性較差，為獲得較多的點(diǎn)位數(shù)量，同時(shí)選取部分河流、湖泊的點(diǎn)位進(jìn)行采樣。共獲得9個(gè)點(diǎn)位，各點(diǎn)位的經(jīng)緯度及實(shí)測(cè)葉綠素a濃度cchl-a如表1所示。

表1 采樣點(diǎn)經(jīng)緯度及實(shí)測(cè)葉綠素a濃度

1.2 衛(wèi)星影像預(yù)處理及水體光譜數(shù)據(jù)提取

使用的高分五號(hào)衛(wèi)星影像的基本情況如表2所示。其中，影像L1_60372_GF5_2874264同時(shí)用于樣點(diǎn)數(shù)據(jù)選取及后續(xù)反演，影像L1_60371_GF5_2874334及L1_60373_GF5_2874323僅用于選取樣點(diǎn)數(shù)據(jù)。

表2 高分五號(hào)衛(wèi)星影像的基本信息Table 2 Basic information of GF-5 satellite images

影像數(shù)據(jù)預(yù)處理包括地理校正、輻射定標(biāo)、大氣校正及半手工水體提取。選取5～10個(gè)地面控制點(diǎn)(ground control point, GCP)，采用亞米級(jí)差分GPS現(xiàn)場(chǎng)測(cè)量地理坐標(biāo)，選用polynomial二階校正模型進(jìn)行幾何精校正。使用ENVI 5.3軟件中Radiometric Correction的子工具Radiometric Calibration將高分五號(hào)衛(wèi)星影像數(shù)據(jù)包中的輻射定標(biāo)參數(shù)代入影像后即完成輻射定標(biāo)。采用FLAASH進(jìn)行大氣校正。圖1為處理前和大氣校正后的長(zhǎng)江河口影像。在提取可見(jiàn)-近紅外波段差值的閾值基礎(chǔ)上，手動(dòng)操作去除易與人工建筑光譜混淆的河流、湖泊，以實(shí)現(xiàn)水體的提取。根據(jù)實(shí)際采樣點(diǎn)的經(jīng)緯度，在ENVI 5.3軟件中提取各點(diǎn)位的水體光譜數(shù)據(jù)(見(jiàn)圖2)作為反演模型的輸入?yún)?shù)。

圖1 預(yù)處理前后的遙感影像Fig.1 Remote sensing images before and after preprocessing

圖2 各點(diǎn)位的光譜數(shù)據(jù)Fig.2 Radiation data of each point

1.3 遙感反演建模方法

基于常規(guī)的數(shù)學(xué)反演方法，采用數(shù)種代數(shù)模型通過(guò)計(jì)算機(jī)篩選不同波段組合方式中決定系數(shù)R2最大的波段組合。根據(jù)實(shí)測(cè)的葉綠素a濃度分布情況，選用指數(shù)形式作為代數(shù)模型的基本形式，通過(guò)改變中間參數(shù)X即波段組合方式實(shí)現(xiàn)模型的差異化處理。

cchl-a=CeAX+B

式中：X為中間參數(shù)，是數(shù)個(gè)波段經(jīng)一定方式的組合；A、B、C為待定系數(shù)。利用計(jì)算機(jī)根據(jù)窮舉法初步列舉各種波段組合中各項(xiàng)波段組合參數(shù)與葉綠素a濃度的相關(guān)性，選取各模型中相關(guān)性最強(qiáng)的波段組合參數(shù)及實(shí)測(cè)cchl-a推算待定系數(shù)A、B、C。

1.3.1 簡(jiǎn)單的波段組合參數(shù)

所使用的波段組合方式包括任意兩個(gè)波段的差值與比值的經(jīng)驗(yàn)參數(shù)，如式(1)和(2)所示。

X=R1-R2

(1)

X=R1/R2

(2)

為提升反演精度并降低水體其他參數(shù)的干擾，在式(1)和(2)的基礎(chǔ)上組合得到如式(3)和(4)所示的組合參數(shù)。

(3)

X=(R1-R2)+(R3-R4)

(4)

式中：Rn為在該組合方式中使用的第n個(gè)波段的反射率。

1.3.2 基于生物光學(xué)模型的波段組合參數(shù)

基于生物光學(xué)模型構(gòu)建算法屬于分析或半分析法。相關(guān)算法的波段組合參數(shù)如下：

Mishra等[3]提出的NDCI法的中間參數(shù)如式(5)所示，該算法也可用于監(jiān)測(cè)水體的葉綠素a濃度。

X=(R1-R2)/(R1+R2)

(5)

邢前國(guó)等[18]使用的反射谷算法的中間參數(shù)如式(6)所示。該算法可在降低其他水體參數(shù)的干擾下研究葉綠素a特征吸收峰對(duì)整體光譜的影響。

(6)

式中：λn為在該組合方式中使用的第n個(gè)波段的波長(zhǎng)。

Dall′Olmo等[19]在Gitelson等[4]的基礎(chǔ)上改進(jìn)得到三波段模型，該模型通過(guò)一定假設(shè)消去了除色素顆粒物吸收以外的光學(xué)參數(shù)，得到色素顆粒物對(duì)光的吸收率與3個(gè)波段反射率的關(guān)系，將該關(guān)系作為中間參數(shù)，如式(7)所示。

(7)

Le等[5]提出四波段模型，其將各種物質(zhì)的吸收系數(shù)和反射系數(shù)代入遙感反射率的光學(xué)方程得到組合參數(shù)，如式(8)所示。陳軍等[20]研究指出，4個(gè)波段應(yīng)在波長(zhǎng)660～750 nm內(nèi)選取。

(8)

1.4 遙感反演評(píng)價(jià)體系

反演評(píng)價(jià)體系包括反演模型的統(tǒng)計(jì)學(xué)參數(shù)分析、反演結(jié)果影像分析、誤差穩(wěn)定性分析等3方面，其中，統(tǒng)計(jì)學(xué)參數(shù)分析作定量分析，后兩者作定性分析。使用SPSS 20.0軟件計(jì)算統(tǒng)計(jì)學(xué)參數(shù)。用決定系數(shù)R2表征自變量對(duì)目標(biāo)變量的解釋程度；用平均相對(duì)誤差與均方根誤差ERMSE表征模型的準(zhǔn)確性；用F統(tǒng)計(jì)值及其檢驗(yàn)概率P判斷模型的顯著性，采用的閾值α=0.005。反演結(jié)果影像分析包括圖像判讀和云陰影分析。云陰影分析利用部分影像中的陰影區(qū)域整體光照強(qiáng)度弱于其他區(qū)域，且難以通過(guò)其他已知參數(shù)進(jìn)行穩(wěn)定計(jì)算的特征，對(duì)西北側(cè)距離陳行水庫(kù)較近的一處陰影區(qū)域(見(jiàn)圖3)進(jìn)行判讀。若某一模型在該區(qū)域內(nèi)的反演結(jié)果與其他區(qū)域存在明顯差異，表明該模型受光照條件的影響較大，反之，受不確定因素的影響較小。

圖3 云陰影位置Fig.3 Position of cloud shadow

誤差穩(wěn)定性分析中，為輸入?yún)?shù)的反射率添加隨機(jī)誤差可近似模擬微小隨機(jī)擾動(dòng)對(duì)反演模型的影響。隨機(jī)生成在[0.95,1.05]內(nèi)均勻分布的隨機(jī)變量作為倍率與實(shí)際反射率相乘，得到一個(gè)與實(shí)際反射率有微小差距的模擬反射率，基于模擬反射率使用各種代數(shù)模型進(jìn)行反演得到模擬反演結(jié)果。將全部研究區(qū)域的模擬反演結(jié)果與實(shí)際反演結(jié)果的比值輸出為圖像。

2 結(jié)果與討論

2.1 遙感反演結(jié)果

利用8種中間參數(shù)構(gòu)建的代數(shù)模型如表3所示。

表3 基于8種中間參數(shù)構(gòu)建的代數(shù)模型Table 3 Models built by eight intermediate parameters

8種代數(shù)模型的反演結(jié)果與實(shí)際結(jié)果如圖4所示。由圖4可知：模型f、g、h的葉綠素a濃度反演結(jié)果最接近實(shí)測(cè)結(jié)果；模型a～e中明顯有數(shù)個(gè)點(diǎn)與直線的偏離程度較大，其中模型a在低葉綠素a濃度下擬合精度較好，但在葉綠素a濃度較高時(shí)誤差極大，說(shuō)明模型a難以應(yīng)用于葉綠素a濃度較高的環(huán)境。

圖4 8種代數(shù)模型的葉綠素a濃度反演數(shù)值與實(shí)測(cè)數(shù)值對(duì)比Fig.4 Comparison of inverse values of chlorophyll-a concentration by eight algebraic models and the test results

遙感反演的葉綠素a濃度的整體分布如圖5所示。由圖5可知：模型b、c、e、g的葉綠素a濃度分布總體處于正常水平，且反演濃度大致為模型b>模型e>模型g>模型c。模型a、d、f在部分水域反演得到的葉綠素a濃度顯著高于其他模型，且與實(shí)際認(rèn)知不符，說(shuō)明這4種模型還需作進(jìn)一步檢驗(yàn)與調(diào)整。模型h在九段沙附近反演得到的葉綠素a濃度明顯比其他模型高，說(shuō)明模型h的反演結(jié)果在能觀察到底泥的淺水中會(huì)受底泥的干擾。模型b、e的反演結(jié)果中有大量較為明顯的噪點(diǎn)，模型h也在一定程度上出現(xiàn)噪點(diǎn)。因這些模型主要是基于比值得到的，可見(jiàn)以比值為基礎(chǔ)的代數(shù)模型在應(yīng)用中還需作進(jìn)一步調(diào)整。但同樣基于比值的模型c未出現(xiàn)噪點(diǎn)。

圖5 利用8種代數(shù)模型反演的葉綠素a濃度分布Fig.5 Distribution of chlorophyll-a concentration inversed by eight algebraic models

2.2 模型的統(tǒng)計(jì)學(xué)參數(shù)

8種代數(shù)模型的統(tǒng)計(jì)學(xué)參數(shù)如表4所示。由表4可知：隨著使用的波段數(shù)的增多，模型整體精度逐漸提高。模型a除平均相對(duì)誤差外的各項(xiàng)參數(shù)均最低，但進(jìn)一步計(jì)算發(fā)現(xiàn)，模型a在葉綠素濃度低于20 μg/L時(shí)的反演精度為R2=0.975 7，F(xiàn)=201.111，P=3.13×10-5，說(shuō)明模型a的應(yīng)用仍有一定可行性。模型b的平均相對(duì)誤差明顯較高，這是由該模型在低濃度下有一個(gè)點(diǎn)的反演結(jié)果相對(duì)誤差極高所導(dǎo)致的，表明該模型可能難以在低葉綠素a濃度水體中應(yīng)用。模型h具有最高的R2與F統(tǒng)計(jì)值以及最低的均方根誤差，模型d、f的這兩項(xiàng)指標(biāo)同樣較高；但模型d、f、h的反演結(jié)果圖存在一定缺陷(見(jiàn)第2.1節(jié))，這可能是由影像信號(hào)的極小誤差以及過(guò)擬合所導(dǎo)致的；模型c、e的精度統(tǒng)計(jì)結(jié)果相對(duì)較低，同時(shí)反演結(jié)果圖的整體分布則更符合實(shí)際。因此c、d、e、f、h模型的應(yīng)用可行性還需進(jìn)一步判定。模型g的平均相對(duì)誤差最低，其他精度評(píng)價(jià)指標(biāo)在所有模型中表現(xiàn)良好，整體反演結(jié)果也較為理想。

表4 各模型反演結(jié)果的統(tǒng)計(jì)學(xué)參數(shù)(n=9)Table 4 Statistical parameters of each model(n=9)

將研究結(jié)果與文獻(xiàn)[21-24]進(jìn)行對(duì)比，結(jié)果如表5所示。由表5可知，與使用其他數(shù)據(jù)源的同種模型相比，使用高分五號(hào)衛(wèi)星遙感數(shù)據(jù)構(gòu)建的模型具有更高的R2和更低的平均相對(duì)誤差(波段比值模型對(duì)比結(jié)果除外)，說(shuō)明高分五號(hào)衛(wèi)星數(shù)據(jù)相較于其他數(shù)據(jù)來(lái)源具有明顯優(yōu)勢(shì)。

表5 基于不同數(shù)據(jù)源的反演精度比較

2.3 環(huán)境光照對(duì)反演結(jié)果的影響

圖6為各個(gè)影像截取的中西北側(cè)部分受云陰影影響的區(qū)域。由圖6可知：模型a、b、g在陰影區(qū)域內(nèi)的反演數(shù)值明顯比其他區(qū)域高，即空間光照對(duì)這3種模型的影響較為明顯,表明使用模型a、b、g時(shí)，每次反演都要重新確定光照條件，且當(dāng)同一影像包含光照條件不同的區(qū)域時(shí)也需手動(dòng)修正，這對(duì)實(shí)際操作提出更高的要求。模型f在該區(qū)域出現(xiàn)反演數(shù)值疑似偏低的結(jié)果，說(shuō)明模型f還需進(jìn)一步優(yōu)化。模型c、d、e、h受云陰影影響時(shí)反演結(jié)果未發(fā)生顯著變化，說(shuō)明環(huán)境光照變化不大時(shí)，這4種模型受環(huán)境光照的影響較小，有較高的泛化潛力。但Hu[25]研究發(fā)現(xiàn)，模型c也會(huì)因受到云量以及太陽(yáng)耀斑等因素的影響而出現(xiàn)誤差。

圖6 云陰影對(duì)8種代數(shù)模型反演結(jié)果的影響Fig.6 Effect of cloud shadow on inversion results of eight algebraic models

2.4 隨機(jī)誤差對(duì)反演結(jié)果的影響

各模型添加隨機(jī)數(shù)前后的反演結(jié)果比值分布如圖7所示,其中，紅色表示隨機(jī)噪聲對(duì)反演的影響較大，藍(lán)色表示隨機(jī)噪聲對(duì)反演結(jié)果的影響較小。由圖7可知：在大多數(shù)水域中，小幅度隨機(jī)誤差不會(huì)明顯影響模型a的反演結(jié)果。模型c、e、f、g受小幅度隨機(jī)誤差的影響逐漸變大，且模型f在部分區(qū)域反演結(jié)果受誤差干擾極不穩(wěn)定。模型b、d、h在反射率有較小的變化時(shí)，其反演結(jié)果會(huì)產(chǎn)生較大幅度的變化，特別是模型b在小幅度隨機(jī)誤差影響下的反演結(jié)果受誤差干擾極不穩(wěn)定。說(shuō)明模型b、d、f、h可能難以應(yīng)用于復(fù)雜的實(shí)際環(huán)境?？傮w而言，形式較簡(jiǎn)單的反演模型在隨機(jī)誤差的影響下較為穩(wěn)定，而形式復(fù)雜的模型容易在結(jié)果中放大數(shù)據(jù)源中的隨機(jī)誤差。這可能是因?yàn)樾问綇?fù)雜的模型更注重反射率在精細(xì)條件下的變化。

圖7 各代數(shù)模型的抗噪性能Fig.7 Noise-robust of each algebraic model

3 結(jié) 語(yǔ)

基于實(shí)測(cè)的水體葉綠素a濃度數(shù)據(jù)以及高分五號(hào)衛(wèi)星的高光譜影像，構(gòu)建8種反演代數(shù)模型并將其應(yīng)用于長(zhǎng)江河口的水質(zhì)反演。結(jié)果表明：

(1)使用高分五號(hào)衛(wèi)星在長(zhǎng)江河口進(jìn)行葉綠素a濃度的反演具有可行性，且相比其他衛(wèi)星數(shù)據(jù)源具有明顯優(yōu)勢(shì)。其中模型g的反演精度最高，其在環(huán)境光照條件基本穩(wěn)定時(shí)反演精度最好，僅在云陰影下反演數(shù)值偏高。

(2)光照條件不明時(shí)，模型c、e的反演精度較高，但平均相對(duì)誤差仍較大。

(3)復(fù)雜模型易受環(huán)境中的細(xì)微條件干擾，簡(jiǎn)單模型又因參數(shù)簡(jiǎn)單而擬合結(jié)果較差。故反演模型的復(fù)雜度應(yīng)適中，且應(yīng)盡量減少比值算法的使用。

受自然條件限制，本文采樣點(diǎn)數(shù)量仍有不足，時(shí)間空間范圍有限，未驗(yàn)證更廣的時(shí)空尺度下模型的準(zhǔn)確性。高分五號(hào)衛(wèi)星在更廣時(shí)空尺度下的應(yīng)用，如模型在同一空間不同季節(jié)下的應(yīng)用及同一時(shí)間更廣空間的應(yīng)用，尚待進(jìn)一步研究。