李勇永，郭 濤，張曉萍，王學(xué)恭

（1. 洛陽(yáng)師范學(xué)院 國(guó)土與旅游學(xué)院，河南 洛陽(yáng) 471934；2. 四川省農(nóng)業(yè)科學(xué)院 遙感與數(shù)字農(nóng)業(yè)研究所，四川 成都 610066；3. 天水師范學(xué)院 資源與環(huán)境工程學(xué)院，甘肅 天水 741001）

Google Earth Engine（GEE）是一個(gè)全球尺度地理空間數(shù)據(jù)云服務(wù)平臺(tái)，它繼承了多種地理空間數(shù)據(jù)并提供了超過(guò)1 000種的數(shù)據(jù)類型和算子，用戶無(wú)需下載數(shù)據(jù)即可直接在線分析。近年來(lái)，眾多學(xué)者借助遙感云計(jì)算平臺(tái)GEE 來(lái)獲取水體[1-10]和生態(tài)環(huán)境信息[11-13]。目前，基于GEE平臺(tái)評(píng)價(jià)洪災(zāi)對(duì)鄭州生態(tài)影響的研究較少，洪災(zāi)對(duì)人類造成了巨大傷害，對(duì)不同環(huán)境下植被生長(zhǎng)的影響如何？在GEE平臺(tái)高效的洪災(zāi)動(dòng)態(tài)變化監(jiān)測(cè)基礎(chǔ)上開(kāi)展災(zāi)后生態(tài)影響評(píng)估顯得尤為重要。

本文綜合多源、多時(shí)相遙感數(shù)據(jù)，從不同角度提取鄭州市7·20洪災(zāi)水體信息，并對(duì)洪災(zāi)影響下的城市生態(tài)環(huán)境變化進(jìn)行監(jiān)測(cè)與評(píng)估，進(jìn)而表明多光譜和雷達(dá)影像數(shù)據(jù)都能在洪災(zāi)監(jiān)測(cè)和生態(tài)破壞評(píng)價(jià)中發(fā)揮重要作用，為決策者提供支持，讓政府或公眾可以直觀了解洪災(zāi)的嚴(yán)重程度及其對(duì)生態(tài)的影響程度。

1 研究區(qū)概況與數(shù)據(jù)集

1.1 研究區(qū)概況

鄭州市地處112°42′～114°14′E，34°16′～34°58′N，位于河南省中部偏北地區(qū)，地勢(shì)西南高東北低（圖1）。全市總面積約7 446 km2。鄭州境內(nèi)以黃河、淮河水系為主，旱澇災(zāi)害頻繁，氣候變異性明顯。2021年7 月20 日—7 月30 日，河南省遭遇長(zhǎng)時(shí)間大規(guī)模降雨，鄭州“7·20”特大暴雨強(qiáng)度和范圍突破歷史記錄，全市城鄉(xiāng)大面積受淹，城鎮(zhèn)街道洼地積澇嚴(yán)重、河流水庫(kù)洪水短時(shí)猛漲、山丘區(qū)溪流溝道大量壅水，形成重大自然災(zāi)害[14]。

圖1 研究區(qū)域地形圖

1.2 數(shù)據(jù)集介紹

本文選用Senitnel-1、Sentinel-2 以及MODIS 影像數(shù)據(jù)提取鄭州市洪水災(zāi)害區(qū)域并進(jìn)行災(zāi)后生態(tài)評(píng)價(jià)。選用歐盟委員會(huì)聯(lián)合研究中心（JRC）制作的全球地表水月尺度產(chǎn)品[15]與洪水災(zāi)發(fā)生前、中、后的洪水提取面積進(jìn)行空間分布對(duì)比。選用武漢大學(xué)黃昕教授團(tuán)隊(duì)制作的近30 a中國(guó)土地覆蓋數(shù)據(jù)產(chǎn)品（CLCD）[16]分析不同土地覆蓋條件下洪水災(zāi)害對(duì)區(qū)域生態(tài)環(huán)境的影響。為便于研究，將2021 年6 月6 日—2021 年7 月20日、2021年7月20日—2021年7月30日、2021年7月30日—2021年9月30日界定為洪水發(fā)生前、中和后期3個(gè)階段。文中涉及數(shù)據(jù)集如表1所示。

表1 遙感數(shù)據(jù)集

2 研究方法

2.1 水體提取

洪水淹沒(méi)導(dǎo)致地表覆蓋類型發(fā)生突變，改進(jìn)的歸一化差異水體指數(shù)(MNDWI)[17]可用于目標(biāo)時(shí)段的水體指數(shù)計(jì)算，公式如下：

MNDWI=（GREEN-MIR）/（GREEN+MIR） （1）式中，GREEN為綠光波段；MIR為中紅外波段；分別對(duì)應(yīng)于哨兵2 號(hào)多光譜數(shù)據(jù)的波段3 和波段12。該方法對(duì)于植被、建筑物及裸露地表均具有較好的抑制作用[18]，但是由于洪水降雨期間會(huì)出現(xiàn)云霧遮擋等情況，多光譜成像易受影響。

洪水發(fā)生時(shí)，SAR影像中地表面散射體的雷達(dá)后向散射系數(shù)在時(shí)序上會(huì)呈現(xiàn)異常波動(dòng)，影像中亮度代表回波強(qiáng)度，水體回波強(qiáng)度較小，陸地回波強(qiáng)度較大，水體在影像中呈現(xiàn)暗色或黑色，陸地呈現(xiàn)灰白色或黑灰色[19]。SAR 影像雖然具有較強(qiáng)的穿透性，能克服云霧干擾，但受地形因素影響較大，往往會(huì)把地形陰影識(shí)別為水體。

水體與陸地其他地物之間的差異在MNDWI 指數(shù)和雷達(dá)后向散射值的灰度圖像中較為突出，但還需使用閾值對(duì)圖像進(jìn)行分割?；诮?jīng)預(yù)處理的影像灰度直方圖，以目標(biāo)和背景間方差最大為基礎(chǔ)假設(shè)，采用最大類間方差法動(dòng)態(tài)確定圖像的全局最優(yōu)閾值，實(shí)現(xiàn)水體提取[7]。

2.2 遙感生態(tài)指數(shù)評(píng)價(jià)

從綠度、濕度、熱度及干度四個(gè)方面建立洪水災(zāi)區(qū)的遙感生態(tài)質(zhì)量評(píng)估體系，并采用主成分分析法構(gòu)建綜合遙感生態(tài)指數(shù)[20-22]，即：

式中，RESI 為遙感生態(tài)指數(shù)；NDVI 為植被指數(shù)；WET 為濕度；LST 為地表溫度；NDBSI 為土壤指數(shù)。文中采用的主成分分析是一種采取依次垂直旋轉(zhuǎn)坐標(biāo)軸的方法，將多個(gè)變量的信息通過(guò)線性變換集中到少數(shù)幾個(gè)特征分量的多維數(shù)據(jù)壓縮技術(shù)[23]。該方法的最大優(yōu)點(diǎn)是集成過(guò)程中，各分量指標(biāo)的權(quán)重主要根據(jù)數(shù)據(jù)本身的性質(zhì)以及各指標(biāo)的貢獻(xiàn)度客觀確定，避免人為確定，造成結(jié)果偏差[24]。

2.3 植被覆蓋率評(píng)價(jià)

植被覆蓋率（FVC）指包括地上葉、莖、枝在內(nèi)的植被垂直投影面積占統(tǒng)計(jì)區(qū)域總面積的百分比，是衡量植被生長(zhǎng)狀況的重要指標(biāo)[25-26]。獲取區(qū)域植被覆蓋率及其變化，對(duì)揭示地表生態(tài)環(huán)境空間變化規(guī)律、探索變化驅(qū)動(dòng)因素及評(píng)價(jià)區(qū)域生態(tài)環(huán)境具有重要意義。文中基于NDVI，通過(guò)像元二分模型來(lái)估算不同時(shí)期的植被覆蓋率。

式中，NDVIsoil為純土壤像元的植被指數(shù)；NDVIvege為純植被像元的植被指數(shù)。

2.4 技術(shù)路線

基于GEE平臺(tái)對(duì)洪水災(zāi)害發(fā)生前、中、后的水體進(jìn)行提取，并采用遙感生態(tài)指數(shù)和植被覆蓋率評(píng)價(jià)洪水對(duì)區(qū)域生態(tài)環(huán)境的影響，關(guān)鍵過(guò)程如下。

1）在GEE 平臺(tái)中篩選出研究區(qū)洪水災(zāi)害發(fā)生前（63景）、中（14景）、后（89景）的Sentinel-2 MSI數(shù)據(jù)，然后以無(wú)云像元重構(gòu)合成3 個(gè)目標(biāo)時(shí)段各自的最小云量影像集合，并求取集合中每幅合成影像的MNDWI 值，最后基于MNDWI 指數(shù)集合求取MNDWI均值影像，以克服云霧遮擋。

2）在GEE 平臺(tái)中篩選出研究區(qū)洪水災(zāi)害發(fā)生前（13 景）、中（3 景）、后（20 景）的Sentinel-1 SAR GRD 數(shù)據(jù)，然后選取VV 極化波段并掩膜掉邊緣無(wú)效值形成影像集合，最后基于經(jīng)過(guò)掩膜的VV 波段集合求取SAR均值影像，以克服地形影響。

3）采用最大類間方差法確定MNDWI均值影像和SAR均值影像的水陸分割閾值，分別實(shí)現(xiàn)水體信息粗提取，然后聯(lián)合SAR影像和多光譜影像處理結(jié)果，當(dāng)兩個(gè)提取結(jié)果在同一地理位置都是水體時(shí)，則判定為水體，否則都為非水體，以獲取水體信息的精細(xì)提取結(jié)果。

4）采用吉林一號(hào)高分影像，基于人工標(biāo)注真值進(jìn)行定量驗(yàn)證；采用JRC水體數(shù)據(jù)集，基于空間分布進(jìn)行定性驗(yàn)證；采用GIS 空間統(tǒng)計(jì)技術(shù)分析洪水淹沒(méi)情況和不同土地覆蓋類型的淹沒(méi)程度。

5）在GEE 平臺(tái)中篩選出研究區(qū)2019—2022 每年8 月份的MODIS 數(shù)據(jù)產(chǎn)品，分別用于計(jì)算熱度（MOD11A2 每年各4 景），綠度（MOD13A1 每年各2景）、 濕 度（MOD09A1 每 年 各4 景） 和 干 度（MOD09A1每年各4景），然后通過(guò)像元篩選合成最小云量影像以克服云霧遮擋，通過(guò)JRC水體數(shù)據(jù)進(jìn)行掩膜以消除大量水體對(duì)遙感指數(shù)主成分計(jì)算的影響，最后基于研究區(qū)無(wú)云合成影像計(jì)算綠度、濕度、熱度和干度，并用主成分分析計(jì)算獲取研究區(qū)2019—2022年各年相同時(shí)段的遙感生態(tài)指數(shù)，將2019、2020 和2022三年的均值作為非受災(zāi)年遙感生態(tài)指數(shù)。

6）在GEE 平臺(tái)中篩選出研究區(qū)2019—2022 每年8—9 月份的Sentinel-2 MSI 數(shù)據(jù)（分別有89、84、84、76 景），然后通過(guò)像元篩選合成每年的最小云量影像以克服云霧遮擋，并基于合成的最小云量影像計(jì)算NDVI，最后通過(guò)像元二分模型來(lái)估算2019—2022年各年相同時(shí)段的植被覆蓋率，并將2019、2020 和2022三年的均值作為非受災(zāi)年植被覆蓋率。

7）利用同期對(duì)比法評(píng)價(jià)災(zāi)后不同土地覆蓋類型的生態(tài)狀況。

3 結(jié)果與分析

3.1 洪水淹沒(méi)情況

3.1.1 洪水監(jiān)測(cè)結(jié)果驗(yàn)證

利用吉林一號(hào)高分影像的人工標(biāo)注真值進(jìn)行精度定量評(píng)價(jià)，其中水體和非水體標(biāo)注各120 個(gè)（圖2）。通過(guò)建立混淆矩陣計(jì)算出評(píng)價(jià)洪水監(jiān)測(cè)結(jié)果好壞的總體精度及表征水體提取一致性程度的Kappa 系數(shù)，對(duì)研究結(jié)果進(jìn)行驗(yàn)證評(píng)價(jià)（表2）。

圖2 吉林一號(hào)衛(wèi)星影像與人工標(biāo)注真值（Copyright，2021年7月24日由長(zhǎng)光衛(wèi)星技術(shù)有限公司提供）

由表2 可知，提取結(jié)果達(dá)到高度一致的水平，存在個(gè)別錯(cuò)分的原因有：①洪水期間地表覆蓋的形態(tài)發(fā)生劇烈變化，影像數(shù)據(jù)集無(wú)法精確識(shí)別；②水體提取是以洪水發(fā)生期間，無(wú)云像元中值合成目標(biāo)時(shí)段的最小云量影像為基礎(chǔ)，所得結(jié)果表征的是整個(gè)洪水期間地表水體覆蓋基本情況，而以拍攝于7月24日的高分影像進(jìn)行驗(yàn)證會(huì)客觀地導(dǎo)致部分驗(yàn)證誤差。

表2 監(jiān)測(cè)結(jié)果驗(yàn)證情況

將JRC的地表水空間分布與洪水發(fā)生前中后各階段的地表水空間分布進(jìn)行提取精度定性評(píng)價(jià)（圖3）。經(jīng)空間統(tǒng)計(jì)可知，JRC 的地表水面積為196.27 km2，洪水發(fā)生前為148.64 km2，洪水發(fā)生時(shí)為976.66 km2，洪水發(fā)生后為148.3 km2。JRC 計(jì)算的水域面積為全年總的水流量與季節(jié)性水的總和，文中計(jì)算的水域面積僅為6 月1 日—7 月20 日日的地表水面積，其總量小于JRC 的面積。通過(guò)自西向東統(tǒng)計(jì)不同經(jīng)度范圍內(nèi)的面積（圖3）可知，采用方法與JRC 保持了較高的空間一致性。

圖3 洪災(zāi)前、中、后與JRC的水域面積對(duì)比

3.1.2 洪水淹沒(méi)情況分析

為準(zhǔn)確了解不同土地覆蓋的受災(zāi)狀況，從洪災(zāi)期間水域部分去除JRC水域，得到實(shí)際的洪水覆蓋區(qū)域，并分析該部分水域內(nèi)的土地覆蓋狀況（圖4和表3）。結(jié)果表明，主要受災(zāi)的土地覆蓋類型為耕地和不透水面，其中，耕地淹沒(méi)面積為386 km2，占區(qū)域總面積的比例最大，不透水面淹沒(méi)面積達(dá)274 km2，受影響的比例最大。

表3 主要土地覆蓋受災(zāi)情況統(tǒng)計(jì)

圖4 洪災(zāi)前、中、后與JRC的水域面積對(duì)比

3.2 災(zāi)區(qū)生態(tài)狀況變化

3.2.1 遙感生態(tài)指數(shù)變化

基于空間化的綠度、濕度、熱度和干度等遙感生態(tài)指數(shù)（圖5），通過(guò)主成分分析獲取2019—2022 年每年8 月的RESI，然后將非受災(zāi)年的RESI 均值與受災(zāi)年2021 年進(jìn)行對(duì)比，以綜合反映研究區(qū)在受災(zāi)后的生態(tài)狀況，將RSEI 以0.2 為間隔分為I（0～0.2）、Ⅱ（0.2～0.4）、Ⅲ（0.4～0.6）、Ⅳ（0.6～0.8）、V（0.8～1.0），分別代表生態(tài)狀況的差、較差、中等、良、優(yōu)5個(gè)等級(jí)（圖7）。由非受災(zāi)年至2021年研究區(qū)同期遙感生態(tài)指數(shù)的等級(jí)演變統(tǒng)計(jì)結(jié)果（表4）可知，2021年8月，生態(tài)狀況等級(jí)提高的面積為2 060 km2，占比27.5%；等級(jí)退化的面積為422 km2，占比5.6%；等級(jí)不變的面積為4 616 km2，占比60.9%。從遙感生態(tài)指數(shù)的角度看，研究區(qū)生態(tài)狀況總體較為穩(wěn)定，呈轉(zhuǎn)好趨勢(shì)。以耕地為例，研究區(qū)耕地占比60%以上，大部分位于東部地區(qū)，但淹沒(méi)比例較小，一方面降雨過(guò)后耕地濕度變大，另一方面轉(zhuǎn)好的墑情促使以夏玉米為代表的作物長(zhǎng)勢(shì)向好，耕作區(qū)植被覆蓋率變大。以上兩方面對(duì)生態(tài)環(huán)境質(zhì)量起正反饋?zhàn)饔?，能夠促進(jìn)生態(tài)環(huán)境的好轉(zhuǎn)[27]。

表4 非受災(zāi)年至2021年洪水災(zāi)區(qū)遙感生態(tài)指數(shù)等級(jí)轉(zhuǎn)移矩陣/km2

圖5 主要土地覆蓋受災(zāi)分布圖

圖6 非受災(zāi)年與2021年的遙感生態(tài)指數(shù)各項(xiàng)指標(biāo)變化

圖7 洪水災(zāi)區(qū)遙感生態(tài)指數(shù)變化

3.2.2 植被覆蓋率變化

對(duì)時(shí)空范圍內(nèi)經(jīng)過(guò)預(yù)處理的哨兵2 號(hào)多光譜影像進(jìn)行植被覆蓋率（FVC）提取，其中包含時(shí)間篩選、空間篩選、去云和調(diào)用NDVI 函數(shù)，最后以各時(shí)段內(nèi)的NDVI 均值合成結(jié)果為基礎(chǔ)，評(píng)價(jià)地表植被覆蓋情況。將FVC 以0.2為間隔分為I（0～0.2）、Ⅱ（0.2～0.4）、Ⅲ（0.4～0.46）、Ⅳ（0.6～0.8）、V（0.8～1.0），分別代表生態(tài)狀況的差、較差、中等、良、優(yōu)5個(gè)等級(jí)（圖8）。由非受災(zāi)年至2021 年研究區(qū)同期植被覆蓋率等級(jí)演變統(tǒng)計(jì)結(jié)果（表5）可知，生態(tài)狀況等級(jí)提高的面積為1 704 km2，占研究區(qū)總面積的23.8%；等級(jí)退化的面積為1 725 km2，占研究區(qū)總面積的24.1%；等級(jí)不變的面積為3 730 km2，占研究區(qū)總面積的52.1%。從植被覆蓋率角度看，災(zāi)后生態(tài)環(huán)境質(zhì)量總體穩(wěn)定。由良轉(zhuǎn)為優(yōu)（700 km2）在等級(jí)提高的類型中貢獻(xiàn)較大，表明洪水過(guò)后植被覆蓋率較高的區(qū)域生態(tài)狀況進(jìn)一步改善，這是因?yàn)楹樗饶転橹脖惶峁┧?，又能將養(yǎng)分帶上來(lái)在一定程度上促進(jìn)植被生長(zhǎng)。

表5 非受災(zāi)年至2021年洪水災(zāi)區(qū)植被覆蓋率等級(jí)轉(zhuǎn)移矩陣/km2

圖8 洪水災(zāi)區(qū)植被覆蓋率變化

3.2.3 洪災(zāi)影響下不同土地覆蓋的生態(tài)變化

不同土地覆蓋對(duì)洪水災(zāi)害的敏感性差異較大[28]，為獲取各覆蓋區(qū)域生態(tài)對(duì)洪水的響應(yīng)，首先排除JRC水體區(qū)域，然后統(tǒng)計(jì)淹沒(méi)區(qū)與非淹沒(méi)區(qū)不同土地覆蓋的遙感生態(tài)指數(shù)和植被覆蓋率（表6、7）。

遙感生態(tài)指數(shù)區(qū)域統(tǒng)計(jì)（表6）顯示，淹沒(méi)區(qū)和非淹沒(méi)區(qū)的生態(tài)總體狀況都趨向好轉(zhuǎn)，均值增長(zhǎng)在11%以上。其中，耕地和不透水面的增長(zhǎng)最為顯著，林地和灌木覆蓋條件下的生態(tài)變化不明顯。

表6 不同土地覆蓋的遙感生態(tài)指數(shù)變化/%

植被覆蓋率區(qū)域統(tǒng)計(jì)（表7）顯示，淹沒(méi)區(qū)和非淹沒(méi)區(qū)的生態(tài)狀況總體趨向好轉(zhuǎn)，淹沒(méi)區(qū)均值增長(zhǎng)24%，非淹沒(méi)區(qū)均值增長(zhǎng)12%。耕地和不透水面的生態(tài)狀況顯著轉(zhuǎn)好；裸地區(qū)域沒(méi)有植被覆蓋；林地、灌木和草地的生態(tài)狀況總體下降，尤其是草地類型，在淹沒(méi)區(qū)下降60%，在非淹沒(méi)區(qū)下降22%。需要說(shuō)明的是，表7 中不透水面區(qū)域的植被覆蓋率雖然在洪水過(guò)后都顯著上升，但都在30%以下，根據(jù)不透水面與植被覆蓋率在城市建成區(qū)呈負(fù)相關(guān)關(guān)系，植被覆蓋率均值低于30%可劃分為不透水面類型[29]，證明研究結(jié)果具有較高的可信度。

表7 不同土地覆蓋的植被覆蓋率變化/%

總體來(lái)看，受人類活動(dòng)影響較大的耕地和不透水面覆蓋區(qū)的生態(tài)狀況均有好轉(zhuǎn)，受人類活動(dòng)影響

較小的林地、灌木、草地覆蓋區(qū)生態(tài)狀況輕微變差；基于較高分辨率哨兵2 號(hào)多光譜影像獲取的植被覆蓋率對(duì)不同土地覆蓋下生態(tài)變化的辨別能力要優(yōu)于基于中分辨率MODIS 數(shù)據(jù)獲取的遙感生態(tài)指數(shù)。

4 結(jié)論與討論

基于GEE 云平臺(tái)，使用MNDWI 函數(shù)并參考最大類間方差法對(duì)Senitnel-1和Sentinel-2影像數(shù)據(jù)進(jìn)行處理，提取洪水災(zāi)區(qū)的水體信息，然后采用MODIS 和Sentinel-2 MSI 數(shù)據(jù)分別計(jì)算非受災(zāi)年與2021 年災(zāi)后同期的遙感生態(tài)指數(shù)和植被覆蓋率，實(shí)現(xiàn)生態(tài)狀況快速監(jiān)測(cè)分析。主要結(jié)論如下：

1）采 用GEE 聯(lián) 合Sentinel-1 SAR GRD 和Senti?nel-2 MS數(shù)據(jù)集獲取災(zāi)區(qū)的水體覆蓋精度較高；耕地受災(zāi)面積最大，不透水面受災(zāi)比例最高。

2）相比非受災(zāi)年，研究區(qū)生態(tài)狀況變化整體趨勢(shì)穩(wěn)定，災(zāi)后生態(tài)狀況等級(jí)提高部分的占區(qū)域總面積的20%以上。

3）淹沒(méi)區(qū)和非淹沒(méi)區(qū)的生態(tài)狀況都趨于好轉(zhuǎn)，與人類活動(dòng)密切關(guān)聯(lián)的耕地、不透水面覆蓋區(qū)域的生態(tài)狀況呈顯著好轉(zhuǎn)，受人類活動(dòng)影響較小的林地、灌木及草地的生態(tài)狀況呈輕微惡化。