基于MODIS數(shù)據(jù)的淤泥質(zhì)海岸水體指數(shù)比較與分析
——以黃河三角洲海岸為例

吳佳平， 張旸， 張杰， 范勝龍， 楊超， 張小芳

(1.福建農(nóng)林大學(xué)資源與環(huán)境學(xué)院，福州 350002; 2.江西省科學(xué)院鄱陽(yáng)湖研究中心,南昌 330096)

0 引言

國(guó)內(nèi)外學(xué)者已應(yīng)用水體指數(shù)提取各類地表水體，如Hui等[10]利用NDWI和MNDWI提取鄱陽(yáng)湖水體，研究水體的時(shí)空變化； Li等[11]使用Landsat TM，ETM+和ALI數(shù)據(jù)提取長(zhǎng)江水體，比較不同數(shù)據(jù)之間提取水體的準(zhǔn)確性； 廖程浩等[12]應(yīng)用MODIS數(shù)據(jù)，分析不同水體指數(shù)的水體識(shí)別效果，發(fā)現(xiàn)CIWI的水體識(shí)別效果最為理想； 車向紅等[13]采用綜合多種水體指數(shù)的方法，提取青藏高原湖泊水體，監(jiān)測(cè)青藏高原湖泊水體的面積變化； 王凈等[14]以太湖流域?yàn)檠芯繀^(qū)，比較5種水體指數(shù)的水體識(shí)別能力，結(jié)果表明CIWI的水體識(shí)別能力最優(yōu)。

受海洋、河流等多種水動(dòng)力因素的影響，淤泥質(zhì)海岸水體變化的復(fù)雜性遠(yuǎn)高于內(nèi)陸，因此大多數(shù)水體指數(shù)研究均避開該區(qū)域。此外，由于空間分辨率的限制，大多數(shù)水體指數(shù)都是基于Landsat數(shù)據(jù)開發(fā)，而基于MODIS數(shù)據(jù)的水體指數(shù)研究較少，但是對(duì)于淤泥質(zhì)海岸水體的高度動(dòng)態(tài)特性而言，MODIS的高時(shí)間分辨率又具有很強(qiáng)的吸引力。因此，在基于MODIS數(shù)據(jù)的淤泥質(zhì)海岸水體指數(shù)及其性能特點(diǎn)方面，亟需開展深入研究。本文以典型的淤泥質(zhì)海岸地區(qū)——黃河三角洲海岸為研究區(qū)，使用MODIS數(shù)據(jù)，分析不同地表覆蓋類型的光譜特征，對(duì)比研究NDWI，MNDWI，AWEInsh，AWEIsh，TCW和WI2015這6種水體指數(shù)的水體提取精度和提取誤差，分析地表覆蓋因素對(duì)水體指數(shù)的水體提取性能的影響，為進(jìn)一步開展淤泥質(zhì)海岸水體變化監(jiān)測(cè)與研究提供參考依據(jù)。

1 研究區(qū)概況及數(shù)據(jù)源

1.1 研究區(qū)概況

本文研究區(qū)為現(xiàn)代黃河三角洲海岸，西北至潮河口，東南到小清河口(圖1)。黃河流經(jīng)黃土高原所攜帶的大量泥沙，由于其下游和河口地段地形變緩，在沉積作用下形成了黃河沖積扇平原[15]。歷史上，黃河河口尾閭經(jīng)過多次改道, 最終發(fā)育形成了現(xiàn)代黃河三角洲[16]。該區(qū)域經(jīng)過了150余a的自然發(fā)育和人類的填海造陸活動(dòng)，成為了我國(guó)最大的新增大陸。在海洋和河流動(dòng)力的雙重影響下，三角洲海岸環(huán)境變化劇烈[17]，并形成獨(dú)特的潮汐特征[18]。區(qū)內(nèi)植被結(jié)構(gòu)簡(jiǎn)單，類型較為單一，人工濕地面積逐年增多[19]。土地利用/覆蓋類型以旱地和未利用地為主，受人類活動(dòng)影響較大，不同類型之間相互轉(zhuǎn)化頻繁[20]。

圖1 研究區(qū)位置(2009年5月3日Landsat TM B5(R)，B4(G), B3(B)合成影像)Fig.1 Location of the study area

1.2 數(shù)據(jù)源及預(yù)處理

MODIS是美國(guó)地球觀測(cè)系統(tǒng)(earth observation system，EOS)用于地球整體觀測(cè)的衛(wèi)星傳感器，MODIS數(shù)據(jù)具有高時(shí)間分辨率、高光譜分辨率和覆蓋范圍廣的特點(diǎn)[21]。MOD09地表反射率數(shù)據(jù)由MODIS L1B數(shù)據(jù)經(jīng)過大氣校正獲得(https: //ladsweb.modaps.eosdis.nasa.gov)。該數(shù)據(jù)集包含了3種空間分辨率的反射率數(shù)據(jù)，即空間分辨率為250 m(B1—B2)，500 m(B1—B7)和1 km(B1—B16)的每日反射率數(shù)據(jù)[22]。本文應(yīng)用MRT swath(MODIS reprojection tool swath)軟件對(duì)MOD09數(shù)據(jù)集進(jìn)行幾何糾正和格式轉(zhuǎn)換，提取空間分辨率為500 m的反射率數(shù)據(jù)。

以Landsat數(shù)據(jù)作為驗(yàn)證數(shù)據(jù)源(https: //earthexplorer.usgs.gov)，其成像時(shí)間與MODIS/Terra最為接近，可以最大限度地保持研究數(shù)據(jù)與驗(yàn)證數(shù)據(jù)之間的一致性。研究選取兩者的成像時(shí)間差小于20 min的數(shù)據(jù)(表1)。

表1 研究使用的遙感數(shù)據(jù)Tab.1 Description of remote sensing data used in the study

2 研究方法

2.1 樣本分類與訓(xùn)練樣本選擇

1)樣本分類。基于黃河三角洲海岸地表覆蓋類型分布情況[19]，將樣本劃分為水體和非水體2大類。淤泥質(zhì)海岸地區(qū)的不同水體受水深、含沙量等因素的影響，各類水體的反射率存在差異[23]。將水體分為海水(海岸線向海延伸15 km范圍內(nèi)的水域)、潮灘水體(潮灘表面水體)和陸地水體(包括水庫(kù)、養(yǎng)殖水面、水塘、河流等)。將非水體像元?jiǎng)澐譃榻ǔ蓞^(qū)和非建成區(qū)； 其中，非建成區(qū)劃分為土壤(NDVI<0.2)、稀疏植被(0.2≤NDVI<0.5)和茂盛植被(NDVI≥0.5)[7]； 土壤進(jìn)一步劃分為潮灘土壤和耕地土壤2類。

2)訓(xùn)練樣本選擇。根據(jù)影像的紋理、形狀和區(qū)域特征等標(biāo)志人工解譯出每一期影像的海岸線，以海岸線15 km范圍作樣本選取區(qū)。使用ArcGIS10.2軟件將MODIS數(shù)據(jù)轉(zhuǎn)化為格網(wǎng)(500 m×500 m)矢量數(shù)據(jù)，以Landsat B5(R)，B4(G)，B3(B)合成影像作為參考影像，解譯出樣本選取區(qū)和城市集中區(qū)的所有地表覆蓋類型。利用ArcGIS軟件中的創(chuàng)建隨機(jī)點(diǎn)工具，分別隨機(jī)抽取4 000個(gè)水體像元和3 000個(gè)非水體像元作為訓(xùn)練樣本(表2)。

表2 各種地表覆蓋類型的樣本數(shù)分布Tab.2 Distribution of sample numbers of various land cover types

2.2 水體指數(shù)的計(jì)算

本文基于MODIS數(shù)據(jù)計(jì)算6種基于Landsat反射率數(shù)據(jù)構(gòu)建的水體指數(shù)，即TCW，NDWI，MNDWI，AWEInsh，AWEIsh和WI2015。這6種水體指數(shù)適用于多種水體(包括海水、河流、湖泊等)不同環(huán)境下(包括地形陰影、建筑物、云量等)的水體提取。根據(jù)Landsat和MODIS影像各波段的波長(zhǎng)信息，得出MODIS數(shù)據(jù)水體指數(shù)的計(jì)算公式(表3)，利用ENVI5.2軟件中的Band Math工具計(jì)算水體指數(shù)。

職稱評(píng)聘是人事管理的一項(xiàng)基礎(chǔ)性工作，也是識(shí)別和發(fā)現(xiàn)人才的重要手段，直接關(guān)系著專業(yè)技術(shù)人員的切身利益，影響到企業(yè)人才隊(duì)伍建設(shè)和事業(yè)長(zhǎng)遠(yuǎn)發(fā)展。職稱評(píng)審要樹立正確導(dǎo)向，做到客觀公正，體現(xiàn)公平競(jìng)爭(zhēng)，不僅依賴于制定的科學(xué)的評(píng)審標(biāo)準(zhǔn)，更取決于對(duì)標(biāo)準(zhǔn)條件的準(zhǔn)確掌握。建立職稱量化評(píng)價(jià)體系，是提升職稱評(píng)審公信度的重要途徑。本文以勝利油田工程技術(shù)人員為研究對(duì)象，探索建立了工程技術(shù)職務(wù)量化評(píng)價(jià)體系，并在實(shí)踐中運(yùn)用，取得了較好的效果。

表3 水體指數(shù)計(jì)算公式Tab.3 Calculation formulas of water indexes

注：B1—B7分別為MODIS B1—B7波段的反射率。

2.3 精度評(píng)估方法

接受者操作特征曲線(receiver operating characteristics curve，ROC)用于衡量二分類問題的分類精度，表示某個(gè)閾值區(qū)間內(nèi)真正率(true positive rate，TPR)和假正率(false positive rate，F(xiàn)PR)的一一對(duì)應(yīng)關(guān)系。當(dāng)TPR達(dá)到最大值且FPR為最小值時(shí)取得最佳閾值。計(jì)算公式為

TPR=TP/(TP+FN)，

(1)

FPR=FP/(TN+FP)，

(2)

式中：TPR為真正率，即水體的正確率；FPR為假正率，即水體指數(shù)將非水體錯(cuò)分為水體的概率；TP為真正類的數(shù)量，即水體的正確分類數(shù)；TN為真負(fù)類的數(shù)量，即非水體的正確分類數(shù)；FP為假正類的數(shù)量，即非水體錯(cuò)分為水體的數(shù)量；FN為假負(fù)類的數(shù)量，即水體錯(cuò)分為非水體的數(shù)量。ROC曲線下的面積(area under the ROC curve，AUC)用于衡量分類器的優(yōu)劣，通常AUC的值介于[0.5，1]，值越大表明分類器越好[24]。

R語(yǔ)言的pROC程序包可繪制ROC曲線，并能自動(dòng)統(tǒng)計(jì)出AUC值和最佳閾值[25]。圖2為2008年、2009年和2015年各水體指數(shù)的ROC曲線。

(a) 2008年10月7日 (b) 2009年5月3日 (c) 2009年5月19日

圖2-1 不同水體指數(shù)的ROC曲線

Fig.2-1ROCcurvesofdifferentwaterindexes

(d) 2009年6月4日 (e) 2015年3月25日

圖2-2 不同水體指數(shù)的ROC曲線

Fig.2-2ROCcurvesofdifferentwaterindexes

從圖2中可見，水體指數(shù)在不同成像時(shí)間具有不同的最佳閾值、敏感性和特異性。2009年各水體指數(shù)之間的差異最明顯，2008年各水體指數(shù)差異較明顯，2015年各水體指數(shù)之間差異不大。通過最佳閾值統(tǒng)計(jì)每個(gè)成像時(shí)間各水體指數(shù)的總體精度(樣本的正確分類數(shù)/樣本總數(shù))、制圖精度(即TPR)、用戶精度(水體的正確分類數(shù)/水體分類總數(shù))、Kappa系數(shù)、錯(cuò)分誤差(非水體樣本錯(cuò)分為水體的數(shù)量/水體分類總數(shù))和漏分誤差(水體樣本漏分的數(shù)量/水體分類總數(shù))，然后統(tǒng)計(jì)各水體指數(shù)的平均提取精度及平均提取誤差。

3 結(jié)果與分析

3.1 潮灘區(qū)光譜特征分析

潮灘區(qū)的主要地表覆蓋類型為水體、土壤和植被，本文通過各類型樣本的反射率分析其光譜特征(圖3(a))。在可見光波段，潮灘水體受水體含沙量的影響[26]，從藍(lán)光到紅光的反射率略有增強(qiáng)； 在近紅外和短波紅外波段的反射率逐漸下降。潮灘土壤的反射率從可見光到短波紅外2波段緩慢增強(qiáng)，短波紅外3波段轉(zhuǎn)為下降； 該分布特征主要與潮灘土壤的含水量狀態(tài)有關(guān)[27]。潮灘植被多為稀疏植被，泥土對(duì)植被的光譜影響較大[28]，因此其在紅光的反射率略高于綠光； 從紅光到近紅外的反射率陡增，形成反射峰； 短波紅外1波段之后逐漸下降。通過反射率箱線圖和t檢驗(yàn)方法分析潮灘水體與土壤、植被的波段差異性，如圖3，其中×表示均值。可見光波段(圖3(b)—(d))中，潮灘水體與土壤在綠光波段的差異不顯著(t=1.4，p=0.17)，在藍(lán)光和紅光波段的差異顯著(t=15.8和14.7，p<0.001)； 潮灘水體與植被的差異顯著且稍大(t∈ [19.7,29.2]，p<0.001)。在近紅外波段(圖3(e))，潮灘水體與土壤、植被的反射率差異顯著且較大(t= 35.1和45.9，p<0.001)。在短波紅外波段(圖3(f)—(h))，潮灘水體與土壤、植被的反射率差異顯著且均大于其他波段(t∈[52.2,87.3]，p<0.001)?？傮w來(lái)看，短波紅外波段對(duì)潮灘水體與土壤、植被的區(qū)分度最佳，近紅外波段次之，藍(lán)光和紅光波段一般，綠光波段最差。

(a) 光譜特征曲線(b) 藍(lán)光波段B3 (c) 綠光波段B4 (d) 紅光波段B1

(e) 近紅外波段B2 (f) 短波紅外1波段B5 (g) 短波紅外2波段B6 (h) 短波紅外3波段B7

圖3 潮灘區(qū)MODIS地表反射率特征(×表示均值)

Fig.3SurfacereflectancecharacteristicsofMODISinthetidalflats

3.2 水體提取精度分析

通過統(tǒng)計(jì)水體指數(shù)的平均水體提取精度(表4)，可以更好地反映各水體指數(shù)的提取效果。

表4 各水體指數(shù)的平均水體提取精度Tab.4 The average water extraction accuraciesof different water indexes

整體上看，6種水體指數(shù)的AUC值均高于98%，Kappa系數(shù)均高于0.9，F(xiàn)PR均低于4%，說(shuō)明水體提取效果均較好。綜合多種指標(biāo)分析，AWEInsh的水體提取效果最佳，AUC值達(dá)99.54%，總體精度達(dá)97.29%，Kappa系數(shù)達(dá)0.95，均高于其他水體指數(shù)； 同時(shí)具有較高的制圖精度(96.84%)、較高的用戶精度(97.69%)和最低的FPR(2.14%)。WI2015的提取效果次之，制圖精度達(dá)97.22%，高于其他水體指數(shù)； Kappa系數(shù)達(dá)0.94，總體精度和用戶精度分別達(dá)96.69%和97.39%。雖然TCW和NDWI具有較高的用戶精度，但總體精度和制圖精度相對(duì)較低，兩者的總體精度分別為96.68%和96.59%，制圖精度分別為96.63%和95.78%，表明整體上提取效果不如前者。水體提取效果最差的是MNDWI，AUC值為98.06%，Kappa系數(shù)為0.9，總體精度為93.48%，制圖精度為93.39%，均低于其他水體指數(shù)。

圖4為各水體指數(shù)對(duì)3類水體提取的平均制圖精度。海水的制圖精度均高于90%，表明6種水體指數(shù)的海水提取效果均較優(yōu)。MNDWI，AWEInsh，TCW和WI2015對(duì)潮灘水體的提取效果較好，制圖精度均高于90%； NDWI對(duì)潮灘水的提取效果較差，制圖精度低于70%。各水體指數(shù)對(duì)陸地水體的提取效果較差，制圖精度均低于80%。

圖4 各水體指數(shù)的3種水體制圖精度Fig.4 Mapping precisions of three water bodieswith different water indexes

3.3 水體提取誤差分析

表5統(tǒng)計(jì)了各水體指數(shù)的平均水體提取誤差分布情況。通過水體的錯(cuò)分誤差分布情況，分析地表覆蓋因素對(duì)水體提取效果的影響程度。各水體指數(shù)將潮灘土壤錯(cuò)分為水體的概率最高，均高于1.2%，說(shuō)明潮灘土壤對(duì)各水體指數(shù)的提取效果影響較大。耕地土壤錯(cuò)分為水體的概率次之，均低于0.3%。稀疏植被、茂盛植被和建成區(qū)錯(cuò)分為水體的概率最小，均低于0.03%，說(shuō)明這3種地表覆蓋類型對(duì)水體指數(shù)的提取效果影響最小。

表5 不同水體指數(shù)的平均水體提取誤差Tab.5 Average water extraction errors of different water indexes (%)

3種水體的漏分誤差表明，陸地水體的漏分誤差整體較高，最高可達(dá)3.1%。MNDWI的海水漏分誤差較高，達(dá)2.88%。對(duì)MNDWI的水體提取結(jié)果(圖5)進(jìn)行分析，發(fā)現(xiàn)2009年的3景影像均存在海水漏分的現(xiàn)象，漏分的海水像元主要集中在孤東油田和小清河口附近的近海區(qū)域。原因是這些區(qū)域的海水在短波紅外波段的反射率異常高于綠波段的反射率，導(dǎo)致像元的MNDWI值小于最佳閾值而被劃分為非水體。其他水體指數(shù)的海水漏分誤差較低，范圍在0.18%～0.32%之間。潮灘水體漏分誤差普遍較低，范圍在0.06%～0.88%之間。

(a) 2009年5月3日 (b) 2009年5月19日 (c) 2009年6月4日

圖5 MNDWI的水體提取結(jié)果

Fig.5WaterextractionresultdiagramsofMNDWI

4 討論

本文使用MODIS數(shù)據(jù)提取黃河三角洲海岸地區(qū)水體，6種水體指數(shù)都達(dá)到較好的水體提取效果。AWEInsh的提取效果最佳，其次是WI2015，提取效果最差的是MNDWI，但整體差異并不明顯。MNDWI的海水提取精度低于其他水體指數(shù)，漏分誤差最高，這可能受海水含沙量的影響。水體在短波紅外波段的反射率隨著含沙量的增加而升高[26-27]，黃河水?dāng)y帶大量的泥沙入海之后，其泥沙擴(kuò)散范圍限于離海岸20 km[29]，本文選取的海水樣本位于離海岸15 km范圍內(nèi)，恰好在泥沙的擴(kuò)散范圍內(nèi)，造成部分海水像元在短波紅外波段的反射率異常高于綠波段，導(dǎo)致MNDWI值小于最佳閾值，從而被劃分為非水體。

各水體指數(shù)對(duì)潮灘水體的提取結(jié)果與Ryu等[30]的研究結(jié)果基本符合。Ryu等選取的研究區(qū)與黃河三角洲具有一定程度的相似性，其研究結(jié)果表明短波紅外和近紅外波段對(duì)于潮灘水體較為敏感。研究結(jié)果與Fisher等[7]的研究存在一定的差異。Fisher等的研究結(jié)果表明AWEIsh的提取效果優(yōu)于AWEInsh，這可能是由于研究區(qū)和研究數(shù)據(jù)存在差異造成的： Fisher等的研究區(qū)地形起伏較為明顯； 此外，使用的是Landsat數(shù)據(jù)，與驗(yàn)證數(shù)據(jù)的成像時(shí)間差為0～5 d。其考慮到潮汐的影響，未選擇海岸帶地區(qū)的樣本。而本研究中，黃河三角洲海岸的地形平緩，地形因素對(duì)水體提取效果的影響可忽略； 其次，MODIS與Landsat數(shù)據(jù)在波長(zhǎng)寬度和空間分辨率方面存在差異，MODIS波譜范圍廣且單波段波長(zhǎng)寬度較窄，造成這2種數(shù)據(jù)對(duì)同一種地表覆蓋類型的反射率出現(xiàn)一定的差異。

5 結(jié)論

本文以淤泥質(zhì)海岸典型地區(qū)——黃河三角洲海岸為研究區(qū)，使用MODIS每日反射率數(shù)據(jù)，以Landsat數(shù)據(jù)為驗(yàn)證數(shù)據(jù)源，對(duì)比分析6種水體指數(shù)在淤泥質(zhì)海岸地區(qū)及MODIS數(shù)據(jù)的適用性。主要結(jié)論如下：

1)綜合分析多種精度指標(biāo)得出，AWEInsh的水體提取效果最優(yōu)； WI2015次之； NDWI，AWEIsh和TCW的水體提取效果較為接近； MNDWI的水體提取效果最差。

2)各水體指數(shù)提取海水的制圖精度較高； 陸地水體的制圖精度較低。NDWI提取潮灘水體的制圖精度較低，MNDWI，AWEInsh，TCW和WI2015提取潮灘水體的制圖精度較高。

3)各水體指數(shù)的水體漏分誤差主要來(lái)源于陸地水體； 潮灘水體的漏分誤差較低； MNDWI的海水漏分誤差高于其他水體指數(shù)。

4)通過分析各水體指數(shù)的水體錯(cuò)分誤差分布得出，潮灘土壤對(duì)各水體指數(shù)的提取效果影響較大； 耕地土壤的影響次之； 茂盛植被、稀疏植被和建成區(qū)的影響最小。

此外，潮灘區(qū)的光譜特征分析可為今后開發(fā)適用于MODIS數(shù)據(jù)和淤泥質(zhì)海岸的水體指數(shù)提供參考依據(jù)。