王雪璐，王瑋，馮琦勝，呂志邦，梁天剛

（蘭州大學草地農(nóng)業(yè)科技學院 草地農(nóng)業(yè)生態(tài)系統(tǒng)國家重點實驗室，甘肅 蘭州730020）

雪蓋的時空動態(tài)分布是全球大氣模擬的關(guān)鍵因子，對全球水循環(huán)和氣候變化的研究起到至關(guān)重要的作用［1］。青海高原地區(qū)是我國主要的積雪分布區(qū)之一，同時也是我國六大牧區(qū)之一，雪災是該地區(qū)冬春季最主要的自然災害，嚴重影響草地畜牧業(yè)的可持續(xù)發(fā)展。而遙感技術(shù)以其視野廣、時間分辨率高、信息量大、客觀真實性強等優(yōu)點，可實時監(jiān)測大范圍雪災的成災過程。因此，提高對該地區(qū)積雪范圍的監(jiān)測精度，對牧區(qū)防災減災具有特別重要的意義［2－6］。

積雪在可見光波段具有較高的反射率，在短波紅外波段有較強的吸收特性，而大多數(shù)云在可見光波段和短段紅外波段均具有較高的反射率?；谶@種特殊的光譜特性，研究發(fā)現(xiàn)歸一化差值雪被指數(shù) （normalized difference snow index，NDSI）不僅可以識別出積雪，而且還是一個理想的云雪識別器，因此可以使用中分辨率成像光譜儀（moderate resolution imaging spectroradiometer，MODIS）第4波段（0.545～0.565μm）和第6波段（1.628～1.652μm）反射率計算的NDSI，作為監(jiān)測積雪的重要指標［7］。國內(nèi)外在MODIS積雪產(chǎn)品精度評價方面，已開展了大量的研究工作。在北美、歐洲和我國北疆等地的研究結(jié)果表明，在晴空條件下MOD10A1產(chǎn)品的積雪判識精度在80%以上。Liang等［8］研究結(jié)果表明北疆地區(qū)的該產(chǎn)品積雪分類精度為98.2%；Keith和Amir［9］在亨伯河流域的研究結(jié)果表明該產(chǎn)品2008年1－3月雪蓋比例在75.08%～93.87%。然而，在我國青藏高原和高山地區(qū)的積雪分類精度相對較低，存在少分積雪而多分陸地的現(xiàn)象。張學通［10］研究表明晴空條件下青海省MOD10A1積雪分類精度僅為43.5%。曹云剛等［11］研究表明由于青藏高原山區(qū)地表坡度較大、地表粗糙度高、地表覆蓋類型復雜等原因，容易出現(xiàn)斑塊狀積雪，因此嚴重影響了MOD10A1的積雪監(jiān)測精度。郝曉華等［12］利用Landsat ETM＋以祁連山中部地區(qū)為例對MOD10A1產(chǎn)品的NDSI閾值進行了檢驗，結(jié)果表明該產(chǎn)品采用的NDSI閾值0.4偏高，造成了山區(qū)積雪面積低估的現(xiàn)象。本研究以青海省為例，利用美國國家雪冰數(shù)據(jù)中心（National Snow and Ice Data Center，NSIDC）的 MODIS地表反射產(chǎn)品 MOD09GA，結(jié)合青海省氣象臺站數(shù)據(jù)，探索研究了適合青海高原積雪監(jiān)測的NDSI閾值，以期為準確監(jiān)測青藏高原地區(qū)積雪分布范圍提供科學依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

青海省位于青藏高原的東北部，北緯31°39′～36°12′，東經(jīng)89°45′～102°23′，東西長度約1 200km，南北寬約800km，總面積72.23萬km2。東部地區(qū)地勢起伏多山，但海拔較低，西部地區(qū)地勢平緩，主要由高原和盆地結(jié)構(gòu)組成。全省高原地區(qū)面積占總面積的80%以上，平均海拔在3 000m以上，海拔3 000m以下區(qū)域僅占全省總面積的26.3%；海拔3 000～5 000m的地區(qū)占地面積最大，約為全省面積的67%；5 000m以上的地區(qū)僅占全省面積的5%。全省依地形可分為祁連山地、柴達木盆地和青南高原三區(qū)。氣候?qū)俚湫偷母咴箨懶詺夂颍昶骄鶜鉁貫椋?.6～8.6℃、年降水量為15～750mm。天然草地以草甸植被為主，可利用草場占全省總國土面積的43.8%，其中高寒草甸和高寒草原總面積達2.448 16×107hm2，占全省草地面積的80.88%，其余部分為少量荒漠草地植被［13－15］。雪災是該省境內(nèi)主要的自然災害之一，主要分布在青南高原地區(qū)，包括玉樹、果洛兩州及海南、黃南兩州的南部，其中玉樹、雜多、曲麻萊、清水河、瑪多、達日高海拔區(qū)域雪災尤為頻繁，且范圍廣，影響大，對草地畜牧業(yè)所造成的損失也最嚴重［16］。

1.2 研究數(shù)據(jù)

通過NSIDC網(wǎng)站，訂購并下載了2007年10月－2008年3月的每日MOD09GA（Surface Reflectance Daily L2GGlobal 1km and 500m）和 MOD10A1（MODIS／Terra Snow Cover Daily L3Global 500mGrid）產(chǎn)品。覆蓋研究區(qū)的圖像軌道號分別為h25v05和h26v05，共183d，共計732幅。MOD09GA和MOD10A1產(chǎn)品均以HDF－EOS分層數(shù)據(jù)格式（hierarchical data format，HDF）保存［17］。同時，為了開展積雪監(jiān)測精度研究，收集了青海省2007年10月－2008年3月共計183d47個臺站的氣象觀測數(shù)據(jù)，包括臺站經(jīng)緯度、海拔、雪深、平均氣溫等資料。在積雪觀測時，雪深以整數(shù)登記，將小于0.5cm的積雪記錄為無雪。共計8 601條記錄。

1.3 積雪圖像的分類方法

積雪分類算法的核心是歸一化雪被指數(shù)的計算，此外還要排除水體和一些暗物體的影響，因此利用MOD09GA資料判別積雪的公式如下：

式中，Band2、Band4和Band6分別代表MOD09GA地表反射產(chǎn)品第2波段、第4波段和第6波段的反射率值，其空間分辨率為500m。

在NSIDC發(fā)布的全球MODIS雪被產(chǎn)品（如MOD10A1等）中，積雪制圖算法主要包括以下3個條件：當NDSI≥0.4時，認為可識別出雪蓋面積占50%以上的單個像元，為積雪判識的主要條件之一。水體在可見光波段和短波近紅外波段的反射率和積雪類似，但可通過其第2波段的反射率小于11%的特性將其排除；同時通過暗物質(zhì)在第4波段反射率小于10%這一約束，可防止一些可見光波段反射率極低的暗色物體被誤判為積雪。

在森林覆蓋區(qū)，由于植被對雪蓋反射率的影響，需要對積雪判別條件進行相應的調(diào)整，將滿足0.1＜NDSI＜0.4且 NDVI＞0.1的像元也識別為有雪［4］。

式中，Band1和Band2分別代表MOD09GA地表反射產(chǎn)品第1波段和第2波段的反射率值，其空間分辨率為500m。

1.4 研究區(qū)NDSI閾值的確定

影像的積雪分類精度是確定NDSI閾值的關(guān)鍵。依據(jù)青海省47個氣象臺站的實際觀測值，比較分析不同NDSI閾值條件下的積雪分類圖像精度，確定適合研究區(qū)積雪分類的NDSI閾值。積雪分類圖像的精度主要通過積雪分類精度R（%）、總精度P（%）、漏測誤差D（%）和多測誤差V（%）這些指標進行綜合評價。其計算公式如下：

式中，K代表影像和地面臺站均有雪的樣本數(shù)，E代表地面臺站有雪的樣本數(shù)，S代表影像和地面臺站均無雪的樣本數(shù)，T代表地面臺站總樣本數(shù)，M代表影像數(shù)據(jù)無雪但地面臺站有雪的樣本數(shù)，N代表影像數(shù)據(jù)有雪但地面臺站無雪的樣本數(shù)，Z代表地面無雪樣本數(shù)。

1.5 自定義積雪分類圖像的合成

對MOD09GA和MOD10A1數(shù)據(jù)進行預處理，分別得到MODIS的1，2，4，6波段反射率數(shù)據(jù)。依據(jù)公式（1）計算研究區(qū)逐日NDSI圖像，采用不同的NDSI閾值分別對積雪進行分類，并結(jié)合NDVI計算林區(qū)雪蓋，使用band2＞0.11和band4＞0.1兩個閾值消除水體和暗物質(zhì)的干擾，利用MOD10A1產(chǎn)品確定云覆蓋范圍，合成不同NDSI閾值的積雪分類圖像（圖1）。合成圖像的編碼方式和MOD10A1積雪分類圖像一致，即25為陸地，50為云，200為積雪。

圖1 用戶自定義積雪分類圖像合成流程圖Fig.1 Flowchart for composition of user－defined snow classification image

2 結(jié)果與分析

2.1 MOD10A1積雪分類精度分析

青海省2007年10月－2008年3月的積雪季共有183個時相的MOD10A1積雪分類產(chǎn)品，結(jié)合47個氣象臺站觀測的8 601條實測雪情記錄數(shù)據(jù)，對該產(chǎn)品進行精度分析的結(jié)果表明，MOD10A1產(chǎn)品在晴空狀況下的積雪分類精度僅為48.97%，總精度為95.56%。當雪深在1～3cm時，積雪漏測誤差為72.65%，受其影響積雪分類精度僅為27.35%，總精度為95.85%；當雪深大于3cm時，積雪分類精度為86.01%（表1），總精度為99.21%。由此可見，研究區(qū)MOD10A1積雪分類精度明顯低于該產(chǎn)品在其他研究區(qū)的精度驗證結(jié)果。由于MOD10A1積雪分類產(chǎn)品是以MODIS L1B數(shù)據(jù)為基礎(chǔ)，通過大氣校正、幾何糾正，結(jié)合云掩膜產(chǎn)品和水體掩膜產(chǎn)品，利用歸一化雪被指數(shù)NDSI的閾值來判識積雪。在NSIDC發(fā)布的MODIS雪蓋產(chǎn)品中，通用的NDSI閾值為0.40。該閾值是Hall和Riggs［18］在美國阿拉斯加、加利福尼亞、蒙大拿、明尼蘇達地區(qū)利用TM雪蓋圖驗證確定的，他指出部分區(qū)域的NDSI閾值是否合理仍然需要試驗。因此，有必要對青藏高原地區(qū)積雪分類圖像的NDSI閾值進行調(diào)整，以提高該地區(qū)的積雪分類精度。

表1 MOD10A1產(chǎn)品的積雪分類精度Table 1 Snow classification accuracies of composite image MOD10A1for different snow depths

2.2 NDSI閾值比較分析

已有研究表明，當雪深介于0.5～3cm時，MODIS傳感器的積雪識別精度較低［8］。黃曉東等［19］對北疆牧區(qū)MODIS積雪制圖MOD10A1精度評價結(jié)果表明，積雪深度1～3cm時積雪分類精度為7.5%。張學通［10］對青海省2002年10月－2008年3月的6個積雪季MOD10A1精度評價研究結(jié)果表明，1～3cm厚度的積雪分類精度僅為29.93%。由于MODIS圖像的空間分辨率較低（500m），部分地區(qū)的地表坡度較大、常有裸露的石塊和少量植被等地物干擾，所以積雪深度小于3cm時，地表難以被積雪均勻地覆蓋，導致MODSI的積雪分類精度偏低。因此，本研究在探討青海省積雪分類精度時，主要利用雪深大于3cm的臺站積雪實測樣本數(shù)據(jù)，分析確定該研究區(qū)的最佳NDSI閾值。

為了更精準地監(jiān)測青海省雪蓋變化情況，提高積雪分類圖像的精度，本研究分別合成了NDSI≥0.50，0.49，……，0.21和0.20時的雪被分類圖像，結(jié)合氣象臺站實測數(shù)據(jù)，分別統(tǒng)計了臺站有雪圖像也有雪、臺站有雪而圖像無雪、臺站無雪但圖像有雪、臺站無雪圖像也無雪的樣本數(shù)，并計算了相關(guān)精度評價指標。研究結(jié)果表明，晴空下當雪深大于3cm時，隨著NDSI閾值的減小，積雪分類精度和多測誤差均逐漸增高，漏測誤差逐漸減小。即NDSI閾值越大，識別出來的有雪像元越少。當NDSI閾值小于0.37時，較多非雪像元被誤判為有雪像元，而當NDSI閾值大于0.44時，許多有雪像元不能被識別出來。因此，該研究區(qū)NDSI的理想閾值應該在0.37～0.44。分析總精度在不同NDSI閾值下的變化曲線，發(fā)現(xiàn)當NDSI閾值分別為0.44，0.43，0.38，0.37時，雪被分類圖像的總精度最高，均為99.51%（表2）。當NDSI閾值為0.38和0.37時，積雪分類精度明顯高于當NDSI閾值為0.44和0.43時的積雪分類精度，而相應的漏測誤差也降低了近2%，多測誤差基本保持不變（只上升了0.02%）。對NDSI≥0.38和NDSI≥0.37的雪被圖像進行精度比較發(fā)現(xiàn)，當NDSI≥0.37時，其積雪分類精度比NDSI≥0.38的積雪分類精度高出0.7%，且漏測誤差低于NDSI≥0.38時的雪被圖像，所以最終確定NDSI≥0.37作為青海省積雪分類的閾值。

為了更全面地評價自定義NDSI閾值的合理性，本研究分析了NDSI≥0.37時雪被圖像的雪深分段精度。結(jié)果表明，用戶自定義NDSI閾值的雪被圖像積雪分類精度也具有隨積雪深度的增大而提高的特點（表3）。當雪深在1～3cm時，雪被圖像的積雪分類精度相對較低，僅為43.56%。這主要是由于研究區(qū)內(nèi)雪深較淺時NDSI雪被算法受其他地物的影響所造成的。當雪深大于3cm時，積雪分類精度可達90.37%，總精度為99.51%。由此可見，研究區(qū)NDSI閾值確定為0.37時，可以有效地提高積雪分類精度和總精度。

2.3 自定義圖像與MOD10A1的精度比較

利用青海省2007年10月—2008年3月的MOD10A1雪被產(chǎn)品及相應的氣象臺站實測數(shù)據(jù)，比較分析了自定義合成雪被圖像同MODIS每日雪被產(chǎn)品的積雪分類精度。研究結(jié)果表明，自定義雪被產(chǎn)品的積雪分類精度為62.31%，總精度為97.18%，而 MOD10A1產(chǎn)品的積雪分類精度為48.97%，總精度為95.56%。同MOD10A1產(chǎn)品相比較，自定義雪被圖像的積雪分類精度和總精度都有不同程度的提高，尤其是積雪分類精度可提高13%，多測誤差從0.36%降為0.22%，漏測誤差從52.12%降至37.69%。因此，改進后的方法無論在漏測與多測誤差上還是分類精度上都好于NSIDC的MODIS雪被產(chǎn)品，更適合于青海高原的積雪分類和動態(tài)監(jiān)測。

對比分析在不同雪深情況下MOD10A1雪被產(chǎn)品和用戶自定義圖像的積雪分類精度表明，當氣象臺站記錄雪深為0cm時，即實測值為陸地的情況下，自定義雪被圖像和MOD10A1雪被產(chǎn)品的多測誤差分別為0.22%和0.36%。這說明NDSI閾值的適度下調(diào)并沒有對多測誤差和總精度造成不良的影響。當雪深為1～3cm時，自定義雪被圖像和MOD10A1產(chǎn)品的積雪漏測樣本數(shù)分別是114和178，漏測誤差分別為56.44%和72.65%。自定義圖像積雪分類精度為43.56%，MOD10A1積雪分類精度為27.35%。說明自定義圖像在對淺雪（雪深為1～3cm）區(qū)域積雪漏測的狀況有明顯改善。隨著積雪深度的增加，自定義雪被圖像和MOD10A1雪被圖像的積雪分類精度都在提高，當雪深在4～6cm時，自定義雪被圖像的積雪分類精度為89.09%，明顯高于MOD10A1的積雪分類精度76.92%。當雪深≥7cm時，其積雪分類精度分別為91.25%和91.36%（表4）。

表2 不同NDSI閾值下的積雪判識精度評價結(jié)果Table 2 Accuracy and errors for snow cover images under different NDSI thresholds %

表3 NDSI≥0.37時合成圖像的積雪分類精度Table 3 Snow classification accuracies of user－defined snow cover images under NDSI≥0.37for different snow depths

表4 不同雪深下自定義產(chǎn)品與MOD10A1產(chǎn)品的積雪分類精度Table 4 Snow classification accuracies of user－defined snow cover and MOD10A1images for different snow depths

3 結(jié)論與討論

本研究通過對青海省范圍內(nèi)MOD10A1雪蓋產(chǎn)品精度分析，發(fā)現(xiàn)其積雪分類精度并不理想，在雪深大于3 cm時的晴空條件下為80.61%；同時結(jié)合青海省氣象臺站的實測雪情數(shù)據(jù)，比較了不同NDSI閾值雪被圖像的精度，確定出晴空下雪深大于3cm時NDSI閾值為0.37的雪被圖像精度最高，其積雪分類精度可達90.37%，總精度99.51%，多測誤差0.22%，漏測誤差為9.63%；最后在比較了MOD10A1產(chǎn)品和自定義NDSI≥0.37的積雪分類圖像精度后發(fā)現(xiàn)，雪深分段統(tǒng)計的結(jié)果表明，當雪深＜7cm時雪被圖像的分類精度較MOD10A1產(chǎn)品有很大的改善作用；當雪深≥7cm時，積雪分類精度無明顯變化。

由此可見，通過下調(diào)NDSI閾值至0.37能夠較好地提高積雪識別的準確率，特別是當雪深＜7cm時，積雪分類精度較MOD10A1產(chǎn)品有明顯的提高。這是由于在薄雪地區(qū)，積雪容易受地表其他地物的影響，很難完整地覆蓋圖像中的整個像元，導致用來判識積雪的反射率特征不夠明顯。因此，下調(diào)NDSI域值可以使更多的淺雪像元識別出來，從而提高淺層積雪和山區(qū)積雪的識別精度；當雪深≥7cm時，自定義雪被產(chǎn)品的積雪分類精度同MODIS每日雪被產(chǎn)品的精度基本相同，都高于90%。說明本研究雪被算法的改進并未對≥7cm的積雪識別起到明顯效果。這可能是由于本研究是基于NDSI閾值的積雪分類算法改進，以調(diào)整NDSI閾值的方法來提高積雪識別精度。而這種方法是通過降低或提升積雪判別條件來實現(xiàn)的，因此對雪蓋反射率特征顯著的深層積雪影響不大。

總體來說，本研究通過改進青海省積雪覆蓋度圖像算法，使自定義積雪分類圖像的積雪精度明顯高于MOD10A1產(chǎn)品，尤其是提高了對淺層積雪的識別率，不僅適用于提取該研究區(qū)內(nèi)的雪被圖像，而且更加符合我國青藏高原地區(qū)的實際情況，對該研究區(qū)的積雪監(jiān)測具有重要意義。

該項研究存在的一些不足和有待進一步解決的問題主要有以下幾點。首先，精度評價準確性受氣象臺站分布位置限制。本研究所利用的積雪實測數(shù)據(jù)全部來自于青海省內(nèi)的47個氣象臺站，其中大部分臺站分布在青海省的東北部，導致實測數(shù)據(jù)在地理上的分布不均，不能夠準確而全面地代表整個研究區(qū)的情況，從而影響了研究的最終結(jié)果。其次，由于MOD09GA數(shù)據(jù)量較大，從訂購下載到保存再到數(shù)據(jù)的處理的過程比較復雜，因此本研究只選擇了2007—2008年這一積雪季對該區(qū)域NDSI閾值進行了初步的探索。在以后的研究中，可擴展MODIS資料的時空分布范圍，如開展多個積雪季的整個青藏高原地區(qū)的研究，檢驗自定義NDSI閾值的合理性，從而提高研究結(jié)果的準確性。

［1］梁繼，張新煥，王建.基于NDVI背景場的雪蓋制圖算法探索［J］.遙感學報，2007，11（1）：85－93.

［2］張娟，馮蜀青，徐維新.MODIS數(shù)據(jù)在青海省積雪監(jiān)測中的應用［J］.青海氣象，2006，1：55－57.

［3］周陸生，李海紅，汪青春.青藏高原東部牧區(qū)大—暴雪過程及雪災分布的基本特征［J］.高原氣象，2000，19（4）：480－458.

［4］王麗紅，付培建，魯安新.遙感技術(shù)在牧區(qū)雪災監(jiān)測研究中的應用［J］.遙感技術(shù)與應用，1998，13（2）：32－36.

［5］劉興元，龍瑞軍，尚占環(huán).草地生態(tài)系統(tǒng)服務功能及其價值評估方法研究［J］.草業(yè)學報，2011，20（1）：167－174.

［6］韓立輝，尚占環(huán)，任國華，等.青藏高原“黑土灘”退化草地植物和土壤對禿斑面積變化的響應［J］.草業(yè)學報，2011，20（1）：1－6.

［7］Hall D K，Riggs G A，Salomonson V V.Algorithm Theoretical Basis Document（ATBD）for the MODIS Snow and Sea Ice－Mapping Algorithms［DB／OL］.［2011－11－25］.http：／／modis.gsfc.nasa.gov／data／atbd／atbd＿mod10.pdf.2001.

［8］Liang G，Huang D，Wu X.An application of MODIS data to snow cover monitoring in a pastoral area：A case study in Northern Xinjiang，China［J］.Remote Sensing of Environment，2008，112：1514－1526.

［9］Keith N A，Amir A K.Investigation of MODIS Snow Mapping in A Boreal Forest Watershed［C］.Baltimore：ASPRS，2009.

［10］張學通.青海省積雪監(jiān)測與青南牧區(qū)雪災預警研究［D］.蘭州：蘭州大學，2010.

［11］曹云剛，楊秀春，徐斌，等.MODIS在青藏高原大范圍積雪制圖中的應用及存在的問題［J］.科技導報，2007，25（21）：51－54.

［12］郝曉華，王建，李弘毅.MODIS雪蓋制圖中NDSI閾值的檢驗——以祁連山中部山區(qū)為例［J］.冰川凍土，2008，30（1）：132－138.

［13］張國勝，伏洋，楊瓊，等.青海省天然草地類型空間分布特征及氣候分區(qū)［J］.草業(yè)科學，2009，26（1）：23－29.

［14］陳鐘，馮琦勝，黃曉東，等.青海省主要類型天然草地EVI動態(tài)變化［J］.草業(yè)科學，2010，27（5）：23－29.

［15］李東，黃耀，吳琴，等.青藏高原高寒草甸生態(tài)系統(tǒng)土壤有機碳動態(tài)模擬研究［J］.草業(yè)學報，2010，19（2）：160－168.

［16］徐裕華.高原牧區(qū)雪凍災害的氣候特征［J］.四川氣象，1998，（2）：28－32.

［17］Riggs G A，Hall D K，Salomonson V V.MODIS Snow Products User Guide Collection 5［DB／OL］.［2011－01－20］.http：／／modis－snow－ice.gsfc.nasa.gov／uploads／sug＿c5.pdf.2006.

［18］Hall D K，Riggs G A.Accuracy assessment of the MODIS snow－cover products［J］.Hydrological Processes，2007，21：1534－1547.

［19］黃曉東，張學通，李霞，等.北疆牧區(qū)MODIS積雪產(chǎn)品MOD10A1和MOD10A2的精度分析與評價［J］.冰川凍土，2007，29（5）：722－729.