基于GEE 平臺的黃河流域水體指數(shù)研究

孫根云，邵寶婕，丁孫金衍，姚延娟，付 航，王霓妮

（1.中國石油大學(xué)（華東） 海洋與空間信息學(xué)院，山東 青島 266580；2.青島海洋科學(xué)與技術(shù)試點國家實驗室 海洋礦產(chǎn)資源評價與探測技術(shù)功能實驗室，山東 青島 266071；3.生態(tài)環(huán)境部 衛(wèi)星環(huán)境應(yīng)用中心 國家環(huán)境保護衛(wèi)星遙感重點實驗室，北京 100094）

1 引 言

黃河流域是我國重要的生態(tài)安全屏障［1］和經(jīng)濟發(fā)展的重要區(qū)域，快速、準(zhǔn)確地提取流域水體對流域水資源調(diào)查、綜合治理等具有重要意義。當(dāng)前，基于遙感技術(shù)的水體提取方法已被廣泛應(yīng)用于水資源調(diào)查、海岸線動態(tài)變化監(jiān)測等領(lǐng)域［2］，水體提取方法主要分為水體指數(shù)方法和機器學(xué)習(xí)方法［3－5］。其中機器學(xué)習(xí)方法基于人工提取特征，使用分類器來實現(xiàn)水體和背景地物的區(qū)分，通常需要大量樣本訓(xùn)練分類器，樣本質(zhì)量會直接影響最終的水體提取精度，并且對于大范圍的水體，難以獲得優(yōu)質(zhì)、充足的訓(xùn)練樣本［6－10］。相較而言，水體指數(shù)方法利用波段間的代數(shù)運算來凸顯水體與背景地物的差異，此類方法形式簡單、處理數(shù)據(jù)快速有效、分類精度較高，在大尺度水體提取中具有很大優(yōu)勢［11］。近年來，相關(guān)研究人員提出了一系列水體指數(shù)方法，如：MCFEETERS［12］利用水體在綠色波段反射率高、在近紅外波段反射率低的特點，提出了歸一化差異水體指數(shù)（NDWI），NDWI方法可以有效提取水體，但提取過程中容易將建筑物錯分為水體。針對這一問題，徐涵秋［2］采用中紅外波段代替近紅外波段，提出了改進的歸一化差異水體指數(shù)（MNDWI），但該方法容易受山體陰影（簡稱山影）的影響。Jiang 等［13］基于Sentinel－2 影像提出了陰影水體指數(shù)（SWI），其對于純水、渾濁水、鹽水和浮冰有較好的提取效果。盡管上述水體指數(shù)方法取得了不錯的提取效果，但對于大范圍復(fù)雜水體的提取效果仍有待驗證。

對于大尺度流域水體的提取，通常需要大量的遙感影像才可實現(xiàn)區(qū)域全覆蓋，而傳統(tǒng)的遙感處理平臺難以實現(xiàn)在某一特定時期、較大范圍內(nèi)對大量衛(wèi)星影像的快速調(diào)用［14］。GEE 是專門用于處理衛(wèi)星圖像和其他地理數(shù)據(jù)的平臺，現(xiàn)階段已被廣泛用于地理數(shù)據(jù)大規(guī)模分析［11，15］。本文基于GEE 平臺，在深入分析黃河流域不同背景地物、不同狀態(tài)水體光譜特征的基礎(chǔ)上，提出一種黃河流域水體指數(shù)（YRBWI），并與已有的水體指數(shù)NDWI、MNDWI和SWI進行對比，驗證其去除遙感影像固有山影的有效性以及提取精度。

2 研究區(qū)概況與數(shù)據(jù)選取

2.1 研究區(qū)概況

黃河流經(jīng)黃土高原水土流失區(qū)，是我國西北、華北地區(qū)的重要生態(tài)屏障，為我國北方生態(tài)文明建設(shè)提供了載體，源區(qū)被譽為“中華水塔”［16］。流域面積為79.5萬km2，從西到東橫跨青藏高原、內(nèi)蒙古高原、黃土高原和黃淮海平原4 個地貌單元，是我國第二大流域［17］。黃河源頭常年積雪，上游河水清澈，大面積沼澤地和湖泊較多；中游流經(jīng)黃土高原地區(qū)，水土流失嚴(yán)重；下游泥沙量大，河道泥沙淤積現(xiàn)象嚴(yán)重［14，16］。流域獨特的地貌特征和水文現(xiàn)象決定了常規(guī)的水體提取方法無法滿足全流域水體自動提取的要求。流域水體提取主要存在以下難點：流域地物復(fù)雜、地形多變，山影嚴(yán)重干擾水體提取；流域水體復(fù)雜，高含沙量渾濁水體與湖泊、水庫等清澈水體的光譜特征差別較大，傳統(tǒng)水體指數(shù)大多只對清澈水體有較好的提取效果，無法滿足流域內(nèi)復(fù)雜水體的高精度提取。

2.2 數(shù)據(jù)選取

在眾多光學(xué)數(shù)據(jù)中，Sentinel－2 影像數(shù)據(jù)是唯一在紅邊范圍內(nèi)含有3 個波段（B5、B6 和B7）的數(shù)據(jù)。選用GEE 平臺中2018 年所有可用的黃河流域Sentinel－2 影像，首先通過QA60 波段對云陰影等觀測質(zhì)量較差的像素進行掩膜，再對所有影像進行中值合成，最終將其用于2018 年黃河流域水體提取。結(jié)合Sentinel－2 影像和Google Earth 高分辨率遙感影像，選取2 種典型水體（清澈水體、渾濁水體）和6 種典型背景地物（農(nóng)田、沙漠、山影、植被、建筑物、道路），對每種水體和背景地物分別選取100 個樣本點用于繪制光譜曲線。通過分層隨機抽樣，結(jié)合Sentinel－2 和Google Earth 影像進行目視解譯生成940 個水體樣本和1 073 個非水體樣本，用于水體提取結(jié)果的精度評價。Sentinel－2 影像光譜覆蓋了430～2 300 nm 范圍內(nèi)的13個波段，可見光、近紅外波段的空間分辨率為10 m，紅邊（Red＿edge）、短波紅外（SWIR）波段的空間分辨率為20 m，衛(wèi)星重訪周期為5 d。

3 技術(shù)與方法

3.1 YRBWI 的構(gòu)建

黃河流域地表水體的反射特性經(jīng)常受到周圍環(huán)境（山影、泥沙等）的影響。為研究不同水體及背景地物的光譜特征，分別繪制其光譜曲線（見圖1），參考已有研究［13，18］，選擇Sentinel－2 影像的B5 波段作為紅邊波段。可以看出，從藍色波段到紅邊波段，清澈水體和渾濁水體的反射率變化呈現(xiàn)相反趨勢，且二者與背景地物的反射率相近，無法同時實現(xiàn)兩種水體的提取。從紅邊波段到近紅外（NIR）波段，背景地物的反射率均表現(xiàn)為上升趨勢，而清澈水體和渾濁水體的反射率均表現(xiàn)為下降趨勢，水體與背景地物的反射率差異明顯，容易進行區(qū)分。從近紅外波段到SWIR2 波段，水體與部分背景地物的反射率相近且均表現(xiàn)為下降趨勢，無法進行區(qū)分。

圖1 黃河流域不同水體和背景地物的光譜曲線

對圖1 分析可知，從紅邊波段到近紅外波段，水體與背景地物的反射率差異明顯、變化趨勢相反，因此利用這兩個波段構(gòu)造水體指數(shù)，以達到同時提取清澈水體、渾濁水體并消除背景地物影響的效果。據(jù)此構(gòu)建的初步水體指數(shù)（PWI）如下：

式中：ρRed＿edge、ρNIR分別為影像在紅邊和近紅外波段的反射率。

水體與背景地物相反的反射率變化趨勢是區(qū)分水體與背景地物的關(guān)鍵。然而，根據(jù)圖1 可知，清澈水體和渾濁水體的反射率在紅邊波段和近紅外波段的下降幅度較小，這可能導(dǎo)致水體提取效果不明顯。為凸顯這種變化趨勢，進一步對PWI進行差異增強處理。對于式（1）中PWI的分子項，使用斜率大于1 的單調(diào)遞增函數(shù)可以增強水體反射率的下降趨勢?？紤]到各種背景地物在兩個波段的反射率都在0～1 之間，選取常用的y＝sinx、y＝tanx、y＝ex、y＝lnx這4 種單調(diào)遞增函數(shù)作為增強函數(shù)。由于斜率為1 的y＝x函數(shù)能夠反映出背景地物在兩個波段原本的反射率差值，因此以y＝x作為對照函數(shù)，5 種函數(shù)在反射率為0～1 之間的函數(shù)曲線見圖2。以渾濁水體為例，5 種函數(shù)處理后渾濁水體在紅邊和近紅外波段間的反射率差值見圖3。

圖2 5 種函數(shù)在反射率為0～1 之間的函數(shù)曲線

圖3 5種函數(shù)處理后渾濁水體在紅邊和近紅外波段間的反射率差值

分析圖2 可知，反射率在0～1 之間時y＝lnx函數(shù)的斜率（增長速率）大于其他函數(shù)，且反射率小于0.3時y＝lnx函數(shù)的增長趨勢更加明顯，更加有利于加大水體在紅邊和近紅外波段之間的反射率差異。對比圖3中不同函數(shù)的反射率拉伸效果可知，y＝sinx、y＝tanx函數(shù)幾乎沒有拉伸效果，y＝ex函數(shù)的效果稍好，y＝lnx函數(shù)的拉伸效果最好。因此，選取y＝lnx作為拉伸函數(shù)，對水體在紅邊和近紅外波段的反射率差異進行增強處理。對于式（1）中PWI的分母項，考慮到PWI分母變大會導(dǎo)致其整體減小，從而造成反射率差異效果的抵消，因此對分母項保持原有形式，不做拉伸處理，最終得到黃河流域水體指數(shù)（YRBWI）：

3.2 對比水體指數(shù)

為驗證YRBWI水體指數(shù)的有效性，選取NDWI、MNDWI和SWI水體指數(shù)進行對比，這3 種水體指數(shù)公式分別為

式中：ρGreen、ρSWIR1分別為影像在綠色、SWIR1 波段的反射率。

3.3 精度評價指標(biāo)的確定

采用總體精度（OA）和Kappa系數(shù)評定各水體指數(shù)提取水體的精度?？傮w精度指被正確分類的類別像元數(shù)（樣本數(shù)）與總類別像元數(shù)的比值；Kappa系數(shù)代表正確分類與完全隨機分類產(chǎn)生錯誤減少的比例。OA和Kappa的計算公式如下：

式中：pkk為被正確分為k類的樣本數(shù)；p為樣本總數(shù)；pii為被正確分類為i類的樣本數(shù)；pi＋為分類所得到的第i類的樣本總數(shù)；p＋i為目視解譯得到的第i類的樣本總數(shù)。

4 結(jié)果分析

4.1 精度評價

各水體指數(shù)的最佳閾值以人工方式確定，將0 作為基本閾值，在此基礎(chǔ)上多次調(diào)整，選取水體提取精度最高時對應(yīng)的閾值作為該指數(shù)的最佳閾值（T）。經(jīng)計算，YRBWI、NDWI、MNDWI、SWI的最佳閾值分別為－0.05、－0.10、0.25、0.15，4 種水體指數(shù)在最佳閾值下的OA值分別為96.87%、93.64%、95.78%、96.08%，Kappa系數(shù)分別為0.937 0、0.872 4、0.914 9、0.921 1，因此YRBWI的OA和Kappa值最大，精度最高。

4 種水體指數(shù)在最佳閾值下對渾濁水體、清澈水體和山影的正確分類率見表1，可以看出，與NDWI、MNDWI和SWI相比，YRBWI對渾濁水體和山影有最高的正確分類率，盡管對于清澈水體來說，YRBWI的正確分類率略低于NDWI，但仍取得了99%以上的效果。說明YRBWI在準(zhǔn)確提取黃河流域高含沙量水體的同時，還能有效去除流域內(nèi)大量的山影噪聲，YRBWI能夠很好地適應(yīng)該流域的廣闊范圍和特殊的地理環(huán)境。

表1 4 種水體指數(shù)對渾濁水體、清澈水體和山影的正確分類率

4.2 水體指數(shù)增強性能

為比較不同水體指數(shù)對水體與背景地物之間差異的增強性能，不同水體和背景地物的各水體指數(shù)均值見圖4（a），圖4（b）是對圖4（a）中NDWI、MNDWI和SWI均值的放大圖。與其他3 種指數(shù)相比，不同水體和背景地物的YRBWI均值的差異最大，說明YRBWI有最好的增強性能。另外，清澈水體、渾濁水體的YR?BWI均值均大于0，背景地物均值均小于0，說明YRBWI的分類效果較好。分析圖4（b）可知，NDWI、MNDWI和SWI對于水體與背景地物之間的差異有一定的增強效果，但渾濁水體的NDWI均值小于0，說明NDWI存在把渾濁水體錯分為背景地物的情況；山影的MNDWI和SWI均值大于0，說明MNDWI和SWI存在把山影錯分為水體的情況。由此可知，與其他3 種指數(shù)相比，YRBWI對水體與背景地物之間差異的增強效果明顯，具有較好的分類效果。

圖4 不同地物類型的各水體指數(shù)平均值

4.3 不同水體指數(shù)的提取效果對比

為驗證YRBWI在全流域內(nèi)不同區(qū)域的提取效果，選取3 個代表性區(qū)域（見圖5），不同區(qū)域的各水體指數(shù)提取結(jié)果見圖6～圖8。圖6 中細長的線狀部分是水體，YRBWI提取的水體清晰明亮且與周圍背景地物差異明顯；MNDWI和SWI提取的水體較為細小，存在部分水體的缺失；NDWI提取的水體雖然較為完整，但是周圍山影噪聲較多，水體邊界不清晰。圖6 中紅框部分是山影，YRBWI提取的水體幾乎不存在山影噪聲，而NDWI、MNDWI和SWI提取的水體均存在大量山影噪聲，錯分情況嚴(yán)重。圖7 中綠色橢圓部分是清澈的湖泊，藍色橢圓部分是渾濁的黃河水體，YRBWI、NDWI和SWI提取的清澈水體、渾濁水體都較為完整，相比之下，MNDWI提取的渾濁水體略有缺失；圖7 中的紅框部分是城市區(qū)域，YRBWI、MNDWI和SWI能夠在錯分極少的情況下很好地提取城市中較為細碎的水體，而NDWI存在嚴(yán)重的背景地物錯分現(xiàn)象。圖8 中綠色橢圓部分是黃河支流，YRBWI和NDWI提取的黃河支流水體清晰明亮，MNDWI和SWI能夠識別黃河支流但是水體缺失嚴(yán)重；圖8 中的紅框部分是大量積雪和部分山影，YRBWI對于積雪和山影噪聲的去除效果很好，NDWI提取的水體存在少量積雪，而MNDWI和SWI對于積雪的去除效果極差，積雪的明亮程度甚至超過綠色橢圓部分的水體，說明積雪錯分情況嚴(yán)重。綜上來看，與其他3 種水體指數(shù)相比，YRBWI對黃河流域的水體提取效果最好。

圖5 黃河流域內(nèi)選取的3 個代表性區(qū)域

圖6 A 區(qū)域的各水體指數(shù)提取結(jié)果

圖7 B 區(qū)域的各水體指數(shù)提取結(jié)果

圖8 C 區(qū)域的各水體指數(shù)提取結(jié)果

4.4 閾值敏感性分析

為比較各水體指數(shù)對于閾值的敏感性，計算4 種水體指數(shù)在不同閾值下的精度值。根據(jù)上文確定的YRBWI、NDWI、MNDWI、SWI的最佳閾值分別為－0.05、－0.10、0.25、0.15，以0.05 為間隔選取最佳閾值附近的8 個閾值。其中：YRBWI的閾值范圍是－0.10～0.25，NDWI的閾值范圍是－0.20～0.15，MNDWI的閾值范圍是－0.05～0.30，SWI的閾值范圍是－0.10～0.25。各水體指數(shù)在不同閾值下的OA和Kappa值的變化見圖9?？梢钥闯觯S閾值的增大，4 種水體指數(shù)的OA和Kappa值都呈現(xiàn)先增大后減小的趨勢，但精度值的變化范圍有所不同。以O(shè)A為例，YRBWI的OA值在最佳閾值后雖有所下降，但總體上變化較為平緩，在－0.10～0.25閾值范圍內(nèi)其OA的最大值和最小值相差僅5%左右。而NDWI在－0.20～0.15 閾值范圍內(nèi)其OA的最大值和最小值相差約20%，MNDWI和SWI在各閾值范圍內(nèi)其OA的最大值和最小值均相差15%以上。由此可見，YRBWI的水體提取精度隨閾值變化的波動較小，對閾值變化不敏感，抗干擾性強。而其他3種水體指數(shù)尤其是NDWI的水體提取精度隨閾值變化波動較大，對閾值變化較敏感，抗干擾性弱。

圖9 4 種水體指數(shù)在不同閾值下的精度變化

此外，在－0.10～0.10 閾值范圍內(nèi)，YRBWI的精度最高，而閾值分別大于0.10 和0.15 時，其精度分別低于SWI和MNDWI的精度。以SWI為例，結(jié)合圖4（b）可知，山影的SWI均值大于0，說明大部分山影的SWI值大于0，因此在閾值大于0 的基礎(chǔ)上適度增大閾值可以去除SWI提取結(jié)果中山影噪聲的影響，從而提高SWI的提取精度。而對于YRBWI，結(jié)合圖4（a）可知，水體的YRBWI均值大于0，背景地物的YRBWI均值小于0，因此在大于0 的基礎(chǔ)上增大閾值會使YRBWI存在水體漏分現(xiàn)象，從而導(dǎo)致精度降低。

5 結(jié) 論

針對黃河流域的復(fù)雜水體信息，提出一種新的水體指數(shù)YRBWI用于黃河流域的大尺度水體提取，結(jié)論如下。

（1）使用紅邊和近紅外波段構(gòu)造PWI，在對比多個增長函數(shù)后選擇y＝lnx函數(shù)對PWI進行增強處理得到Y(jié)RBWI。與NDWI、MNDWI和SWI這3 種水體指數(shù)相比，YRBWI在黃河流域水體提取中表現(xiàn)優(yōu)異，既能有效增強渾濁水體和清澈水體與背景地物的光譜差異，準(zhǔn)確提取高含沙量渾濁水體，也能解決大部分指數(shù)存在的山影問題，可實現(xiàn)快速、自動、高精度的水體提取。

（2）相比于其他指數(shù)，YRBWI對閾值變化不敏感，具有良好的抗干擾性。在未來研究中，可考慮結(jié)合水體的時間序列信息來進行黃河流域水體提取和變化研究。