殷崎棟,柳彩霞,田 野
1 陜西生態(tài)環(huán)境規(guī)劃設(shè)計(jì)院有限公司, 西安 710075 2 中國(guó)科學(xué)院空天信息創(chuàng)新研究院遙感科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室, 北京 100101 3 中國(guó)科學(xué)院城市環(huán)境研究所城市環(huán)境與健康重點(diǎn)實(shí)驗(yàn)室, 廈門 361021
加強(qiáng)和發(fā)展自然保護(hù)區(qū)事業(yè)是當(dāng)前我國(guó)生態(tài)環(huán)境保護(hù)工作的一項(xiàng)重要而緊迫的任務(wù)。但隨著我國(guó)經(jīng)濟(jì)的快速發(fā)展和對(duì)自然資源的需求日益增加,保護(hù)與開發(fā)的矛盾日益增加,“人工化”、“生境退化”的趨勢(shì)日益嚴(yán)重,自然保護(hù)區(qū)面臨著前所未有的壓力[1]。
遙感技術(shù)的發(fā)展為自然保護(hù)區(qū)保護(hù)提供了切實(shí)可行的監(jiān)測(cè)和保護(hù)手段,特別是陸地衛(wèi)星Landsat提供30m空間分辨率的定期圖像采集,平均重訪周期為16d,可免費(fèi)獲得并追溯到1980年代早期,可支持植被和土地覆蓋變化的長(zhǎng)期觀測(cè)和監(jiān)測(cè)[2]。利用經(jīng)過大氣校正和地形校正的地表反射率產(chǎn)品,確保多個(gè)時(shí)間序列的圖像具有可比性[3],使Landsat圖像在長(zhǎng)時(shí)間序列地表覆蓋變化分析方面具有優(yōu)越性。在長(zhǎng)時(shí)間序列遙感數(shù)據(jù)的支持下,可以在不同的應(yīng)用領(lǐng)域進(jìn)行生態(tài)環(huán)境和土地變化的評(píng)估、追溯和歷史重構(gòu),如植被監(jiān)測(cè)[4-5]、土地利用變化[6-7]、環(huán)境擾動(dòng)[8- 10]等。在國(guó)內(nèi)自然保護(hù)區(qū)遙感監(jiān)測(cè)方面,學(xué)者們已就嶗山[11]、可可西里[12]、十八里峽[13]、佛坪[14]等自然保護(hù)區(qū)的人類活動(dòng)進(jìn)行了遙感監(jiān)測(cè)。然而,先前的研究大部分集中在時(shí)間序列圖像的離散信息提取,年際變化研究存在年份的空白,無(wú)法捕捉連續(xù)的地表變化特征,從而限制了數(shù)據(jù)的長(zhǎng)期監(jiān)測(cè)和生態(tài)環(huán)境穩(wěn)定性評(píng)估能力。
Landsat時(shí)間序列數(shù)據(jù)可以監(jiān)測(cè)森林?jǐn)_動(dòng),但是需要合理處理因大氣和幾何校正、植被物候、太陽(yáng)高度變化和傳感器退化引起的數(shù)據(jù)噪聲。Kennedy等[15]發(fā)現(xiàn),在美國(guó)俄勒岡州西部,不同的光譜指數(shù)在準(zhǔn)確檢測(cè)細(xì)微干擾信號(hào)方面的能力各不相同。歸一化植被指數(shù)(NDVI)的效率低于歸一化燃燒比(NBR)或纓帽變換濕度(TCW),在森林監(jiān)測(cè)領(lǐng)域,很多算法使用有限的光譜波段或指數(shù)來(lái)檢測(cè)干擾,例如NBR[15]、森林指數(shù)(FI)[16]、NDVI和SWIR-NIR比值[17]或纓帽變換角度[18],因此,在檢測(cè)森林?jǐn)_動(dòng)時(shí)仔細(xì)考慮比較光譜指標(biāo)是非常重要的步驟。
目前,已經(jīng)有多種基于Landsat時(shí)間序列的變化檢測(cè)技術(shù)。例如LandTrendr[15]、BFAST[19-20]、CCDC[21-22]。由于時(shí)間序列算法的數(shù)據(jù)預(yù)處理和算法本身都非常耗時(shí),為了達(dá)到更廣泛的研究應(yīng)用群體,LandTrendr算法得以在Google Earth Engine (GEE)平臺(tái)上實(shí)施[23],從而提高了算法計(jì)算速度。與以前在本地服務(wù)器上運(yùn)行變化檢測(cè)算法相比,GEE中的LandTrendr利用云計(jì)算的優(yōu)勢(shì)更有效地檢測(cè)森林的變化,為森林自然保護(hù)區(qū)的保護(hù)提供了大規(guī)模運(yùn)用的機(jī)會(huì)。
了解森林保護(hù)區(qū)生態(tài)系統(tǒng)的組成、結(jié)構(gòu)和功能的現(xiàn)狀和長(zhǎng)期趨勢(shì)發(fā)展,可以為保護(hù)區(qū)管理提供科學(xué)可信的本底資料。陜西省地處黃河中游、長(zhǎng)江上中游地區(qū),屬于國(guó)家生態(tài)環(huán)境建設(shè)重點(diǎn)區(qū)域。全省共有45個(gè)省級(jí)以上森林自然保護(hù)區(qū)和濕地自然保護(hù)區(qū),其中國(guó)家級(jí)19個(gè),省級(jí)26個(gè)。為了檢驗(yàn)LandTrendr算法對(duì)不同生態(tài)氣候分區(qū)的森林植被類型的敏感性是否存在差異,我們選取陜西省柴松和太白山自然保護(hù)區(qū),利用GEE云計(jì)算平臺(tái)中的LandTrendr算法進(jìn)行森林變化檢測(cè)分析,追蹤近40年來(lái)陜西省柴松和太白山自然保護(hù)區(qū)的森林?jǐn)_動(dòng)。主要回答以下問題:(1)柴松省級(jí)自然保護(hù)區(qū)和太白山國(guó)家級(jí)自然保護(hù)區(qū)自保護(hù)區(qū)成立以來(lái)植被有明顯的變化趨勢(shì)嗎?(2)如果保護(hù)區(qū)內(nèi)植被發(fā)生了明顯變化,變化的空間分布如何?(3)發(fā)生變化的原因是什么?
陜西延安柴松省級(jí)自然保護(hù)區(qū)(東經(jīng)108°39′—108°50′,北緯35°56′—36°07′)位于延安市富縣西部(圖1)。保護(hù)區(qū)成立于2004年。境內(nèi)南北長(zhǎng)約18 km,東西寬約16.5 km,海拔1000—1582 m之間,總面積17640 hm2(表1)。保護(hù)區(qū)植被屬華北暖溫帶半濕潤(rùn)落葉闊葉林地帶的北部落葉闊葉地帶的西緣。森林植被均為天然次生林,成為黃土高原上的植物種類獨(dú)特、多樣、豐富的天然綠色“基因庫(kù)”,其中價(jià)值較高的有柴松、油松、側(cè)柏、櫟類和山楊等86種。主要保護(hù)對(duì)象為柴松林及其森林生態(tài)系統(tǒng)。
陜西太白山國(guó)家級(jí)自然保護(hù)區(qū)(東經(jīng)107°22′25″—107°51′30″和北緯33°49′30″—34°05′35″)位于陜西省太白、眉縣、周至三縣交界處(圖1),保護(hù)區(qū)成立于1965年,1986年7月晉升為國(guó)家級(jí)自然保護(hù)區(qū),占地面積56325 hm2(表1),主要保護(hù)對(duì)象為森林生態(tài)系統(tǒng)和自然歷史遺跡。保護(hù)區(qū)地處秦嶺山脈中段,太白山頂海拔3767 m,為秦嶺的主峰。秦嶺山脈為中國(guó)南北氣候和植物的分界線,太白山則處于華北、華中和青藏高原三區(qū)生物交匯過渡地帶,區(qū)內(nèi)動(dòng)植物資源豐富,植被垂直分帶明顯。高等植物有2000余種,國(guó)家重點(diǎn)保護(hù)植物有連香樹、水青樹、星葉草、太白紅杉等21種。高等動(dòng)物有270多種,國(guó)家保護(hù)動(dòng)物有大熊貓、羚牛、豹等20多種。在太白山海拔3000 m以上的高山區(qū),還保存著比較完整的第四紀(jì)冰川遺跡。該保護(hù)區(qū)是開展科學(xué)研究難得的“天然實(shí)驗(yàn)室”。
圖1 研究區(qū)位置意圖Fig.1 Location map of the study area
表1 研究區(qū)概況
利用GEE平臺(tái)自1984—2018年期間的全部存檔Landsat TM/ETM/OLI的地表反射率數(shù)據(jù)構(gòu)建時(shí)間序列數(shù)據(jù)集。由于藍(lán)光波段對(duì)大氣狀況敏感,只選取了5個(gè)波段:綠(green)、紅(red)、近紅外(nir)和2個(gè)短波紅外波段(swir1和swir2)。為了盡量減少因物候或太陽(yáng)高度角變化引起的反射率變化,選取每年6月10日到8月30日的影像。由于Landsat OLI比TM/ETM+具有更高的輻射分辨率(12比特),我們利用了這兩個(gè)傳感器光譜值之間的統(tǒng)計(jì)協(xié)調(diào)函數(shù)對(duì)反射率進(jìn)行歸一化[24],接著進(jìn)行提取年集合的中位反射率值,構(gòu)建了具有最小云量圖像的35年合成圖,進(jìn)而計(jì)算歸一化植被指數(shù)NDVI[25](公式(1))和歸一化燃燒比值NBR[26](公式(2))。通過比較不同比值指數(shù)和波段對(duì)森林?jǐn)_動(dòng)的敏感性(圖2),NDVI和NBR對(duì)柴松保護(hù)區(qū)的森林?jǐn)_動(dòng)識(shí)別效果更好,單一波段受數(shù)據(jù)噪聲的影響較大。綜合考慮以往研究經(jīng)驗(yàn),本文利用NBR指數(shù)探測(cè)森林?jǐn)_動(dòng)。
(1)
(2)
圖2 不同植被指數(shù)和波段對(duì)森林?jǐn)_動(dòng)的敏感性Fig.2 Sensitive of vegetation index and bands for forest disturbance NBR:歸一化燃燒比值, Normalized Burn Ratio; NDVI:歸一化植被指數(shù), Normalized Difference Vegetation Index; swir2:短波紅外波段2;swir1:短波紅外波段1
利用Kennedy, Yang[15]開發(fā)的LandTrendr算法來(lái)檢測(cè)保護(hù)區(qū)內(nèi)森林的變化。LandTrendr算法的核心是一種時(shí)間分割算法,用于捕捉時(shí)間序列中長(zhǎng)期、漸進(jìn)或者短期的劇烈變化,可以分析每個(gè)像元的時(shí)間序列譜線來(lái)監(jiān)測(cè)其是否發(fā)生變化。每個(gè)像素的輸入是一個(gè)光譜指數(shù)或者波段的時(shí)間序列。通過設(shè)置去除噪聲引起的尖峰(spikeThreshold)、識(shí)別潛在的斷點(diǎn)(vertexCountOvershoot)、擬合軌跡和最佳分段數(shù)(pvalThreshold),尋找最佳模型(bestModelProportion)。因此,LandTrendr要求設(shè)定控制參數(shù)以保證變化檢測(cè)的質(zhì)量,通過實(shí)驗(yàn)分析,本研究對(duì)LandTrendr的參數(shù)如表2所示。在這里需要說明的是,基于保護(hù)區(qū)建立時(shí)間的考慮,對(duì)于柴松保護(hù)區(qū),研究集中在2004—2018年,總計(jì)15a。對(duì)于太白山保護(hù)區(qū),研究階段為1984—2014年(總計(jì)31a)。因此對(duì)柴松保護(hù)區(qū)而言,研究時(shí)段相對(duì)較短,在LandTrendr的參數(shù)設(shè)置中,分割單元的數(shù)目可以稍少,maxSegments設(shè)置為4。根據(jù)以往研究經(jīng)驗(yàn),2個(gè)保護(hù)區(qū)的其他參數(shù)設(shè)置一致。
表2 LandTrendr參數(shù)
對(duì)森林?jǐn)_動(dòng)發(fā)生的空間和時(shí)間分布檢測(cè)對(duì)森林保護(hù)和政策管理非常重要。相比傳統(tǒng)的土地覆蓋分類和變化檢測(cè)技術(shù),LandTrendr可以成功捕獲森林?jǐn)_動(dòng)信號(hào)進(jìn)行空間和時(shí)間制圖。Landsat系列衛(wèi)星能提供長(zhǎng)達(dá)近40年的衛(wèi)星觀測(cè)數(shù)據(jù),使對(duì)森林保護(hù)區(qū)的長(zhǎng)期觀測(cè)和監(jiān)測(cè)成為可能。盡管已有很多利用LandTrendr或者相似時(shí)間序列檢測(cè)算法在本地服務(wù)器的使用案例[15, 27- 29],算法實(shí)施帶來(lái)的數(shù)據(jù)準(zhǔn)備和運(yùn)行帶來(lái)巨大的時(shí)間成本和硬件成本。但是,在GEE支持下的LandTrendr算法極大的提高了算法效率。用戶不再需要從美國(guó)地質(zhì)調(diào)查局或者其他衛(wèi)星地面站下載大量的原始有云衛(wèi)星圖像,而是基于Google云平臺(tái)合成用戶感興趣區(qū)域的無(wú)云或者少云高質(zhì)量拼接圖像,進(jìn)而利用LandTrendr的主要作者和GEE合作開發(fā)的算法進(jìn)行森林?jǐn)_動(dòng)探測(cè)。用戶可以通過Java script編程加載感興趣結(jié)果圖層和時(shí)間序列衛(wèi)星圖像,在GEE平臺(tái)中對(duì)算法結(jié)果進(jìn)一步目視判讀和參數(shù)調(diào)整,最后只需要導(dǎo)出結(jié)果圖層,如擾動(dòng)年份等。
為了準(zhǔn)確評(píng)估森林?jǐn)_動(dòng)變化的準(zhǔn)確性,在研究時(shí)間內(nèi),在柴松保護(hù)區(qū)內(nèi)選擇了100個(gè)像元表示變化,另外100個(gè)像元代表不變區(qū)域。對(duì)每個(gè)像元,通過人工目視判讀合成的15幅圖像,并使用高空間分辨率的Google Earth(GE)圖像對(duì)土地覆蓋進(jìn)一步確認(rèn),記錄該像元是否發(fā)生了變化,何時(shí)(發(fā)生年份)發(fā)生了變化,從而組成200個(gè)像元的驗(yàn)證樣本集,生成混淆矩陣,計(jì)算用戶精度、生產(chǎn)者精度和總體精度。
LandTrendr的時(shí)間分割和擬合可以捕捉森林向非森林轉(zhuǎn)變的光譜特征。表3所示的混淆矩陣說明,生產(chǎn)者精度和用戶精度均達(dá)到85%以上,總體精度達(dá)到93%。對(duì)于檢測(cè)出森林變化的像元,進(jìn)一步評(píng)估了變化年份。表4給出了變化年份的精度,2007年之后的變化年份生產(chǎn)精度都達(dá)到了80%以上。用戶精度除2007年和2005年稍低以外(分別是70%和75%, 其他年份的用戶精度達(dá)100%,變化年份檢測(cè)的總體精度為89%。
表4 森林?jǐn)_動(dòng)年份精度評(píng)估
圖3所示的擾動(dòng)圖顯示了2004—2018年間柴松保護(hù)區(qū)森林?jǐn)_動(dòng)的空間分布和變化年份。森林?jǐn)_動(dòng)發(fā)生的區(qū)域主要集中在柴松保護(hù)區(qū)的實(shí)驗(yàn)區(qū)和緩沖區(qū),但是核心區(qū)也存在部分森林?jǐn)_動(dòng)。擾動(dòng)類型主要為人類活動(dòng)所致。總體上,2004—2018年間,共有42.74 hm2的森林發(fā)生了擾動(dòng)(表5),轉(zhuǎn)化為非森林,共計(jì)251個(gè)斑塊,年變化最大面積值和斑塊數(shù)均發(fā)生在2008年,2007、2015年和2005年次之,森林向非森林轉(zhuǎn)化發(fā)生的斑塊數(shù)依次為57、24和20塊(圖3)。
圖3 柴松保護(hù)區(qū)擾動(dòng)空間和歷史分布Fig.3 Disturbance history of Chaisong Reserve(a. 擾動(dòng)年分布; b. 實(shí)驗(yàn)區(qū)擾動(dòng); c. 緩沖區(qū)和核心區(qū)擾動(dòng); d-e. Landsat圖像(NIR/R/G); f-h. GE高清衛(wèi)星影像;P1位置[108.723276, 36.082742]; P2位置[108.697087, 36.024953])
表5 柴松保護(hù)區(qū)森林?jǐn)_動(dòng)圖斑統(tǒng)計(jì)
圖4所示的擾動(dòng)圖顯示了1984—2018年間太白山保護(hù)區(qū)森林?jǐn)_動(dòng)的空間分布和變化年份。森林?jǐn)_動(dòng)發(fā)生的區(qū)域主要集中在太白山保護(hù)區(qū)核心區(qū)。擾動(dòng)主要由自然因素引起,如第四紀(jì)冰川融化形成的侵蝕溝和其他自然因素造成的森林退化。1984—2018年間,共有23.68 hm2的森林發(fā)生了擾動(dòng)(表6),轉(zhuǎn)化為非森林,共計(jì)57個(gè)斑塊數(shù),年變化最大面積值和斑塊數(shù)均發(fā)生在2013年,1994年次之(7塊,1.87 hm2)(圖4)。
本文結(jié)果表明,LandTrendr可以捕獲不同氣候區(qū)森林?jǐn)_動(dòng)的空間分布特征。在降雨量偏少的柴松保護(hù)區(qū),擾動(dòng)最為頻繁的年份集中在2007—2008年間,這是由于保護(hù)區(qū)內(nèi)違法建設(shè)了一批以石油開采為目的的鉆井平臺(tái)及附屬設(shè)施(圖3),從而造成了森林保護(hù)區(qū)森林覆蓋的破碎化。自2016年之后森林?jǐn)_動(dòng)有減少趨勢(shì)。以環(huán)保部為執(zhí)行主體的自然保護(hù)區(qū)“綠盾”行動(dòng)自2017年在國(guó)家級(jí)保護(hù)區(qū)實(shí)施以來(lái),對(duì)省市級(jí)自然保護(hù)區(qū)的保護(hù)和監(jiān)管起到了督促作用,在保護(hù)區(qū)內(nèi)的部分違法違規(guī)設(shè)施正在逐步清退。
圖4 太白山保護(hù)區(qū)擾動(dòng)空間和歷史分布Fig.4 Disturbance history of Taibaishan Reserve Disturbance yeara. 擾動(dòng)年分布; b. 冰川消融擾動(dòng); c. 其他自然擾動(dòng); d-f. Landsat圖像(NIR/R/G); g-i. GE高清衛(wèi)星影像;P1位置[107.79248, 33.94858]; P2位置[107.687327, 33.966819]
表6 太白山保護(hù)區(qū)森林?jǐn)_動(dòng)圖斑統(tǒng)計(jì)
在降雨量偏多的太白山保護(hù)區(qū),森林?jǐn)_動(dòng)活動(dòng)遠(yuǎn)遠(yuǎn)少于柴松保護(hù)區(qū),擾動(dòng)最為頻繁的年份集中在2013年,且擾動(dòng)活動(dòng)與氣候變化和冰川融化造成的地形地貌變化相關(guān)(圖4)。長(zhǎng)遠(yuǎn)來(lái)看,冰川活動(dòng)的變化,特別是冰川融化導(dǎo)致的裸露冰磧石移動(dòng),會(huì)影響到冰川末端的高山森林帶,可能會(huì)改變高寒地區(qū)森林的分布[30-31]。因此,氣候變化對(duì)自然保護(hù)區(qū)的生態(tài)環(huán)境保護(hù)提出了更多挑戰(zhàn)[32],進(jìn)行長(zhǎng)時(shí)間序列的遙感監(jiān)測(cè)能為處于高寒地帶人類難以到達(dá)的森林生態(tài)系統(tǒng)類保護(hù)區(qū)提供可參考依據(jù)。
進(jìn)一步確定自然保護(hù)區(qū)內(nèi)森林退化的誘導(dǎo)因素,哪些是人類活動(dòng),特別是非法人類活動(dòng)造成的,如違章建筑、違章工礦設(shè)施。哪些是自然因素造成的,如森林病蟲害、氣候變化等。本研究涉及的柴松自然保護(hù)區(qū)的關(guān)鍵擾動(dòng)因素是因發(fā)展經(jīng)濟(jì)驅(qū)動(dòng)的人類活動(dòng),有少量活動(dòng)甚至是在保護(hù)區(qū)的核心區(qū)。雖然造成森林?jǐn)_動(dòng)的土地面積與整個(gè)保護(hù)區(qū)的面積相比,占比并不大,但是森林?jǐn)_動(dòng)破壞了保護(hù)區(qū)動(dòng)物棲息地的連通性,造成了自然景觀的破碎化,開采活動(dòng)造成的噪聲、粉塵及其他行為對(duì)動(dòng)植物的影響還需進(jìn)一步評(píng)估。根據(jù)我國(guó)《自然保護(hù)區(qū)條例》規(guī)定:自然保護(hù)區(qū)內(nèi)保存完好的天然狀態(tài)的生態(tài)系統(tǒng)以及珍稀、瀕危動(dòng)植物的集中分布地,應(yīng)當(dāng)劃為核心區(qū),禁止任何單位和個(gè)人進(jìn)入;除依照本條例第二十七條的規(guī)定經(jīng)批準(zhǔn)外,也不允許進(jìn)入從事科學(xué)研究活動(dòng);核心區(qū)外圍可以劃定一定面積的緩沖區(qū),只準(zhǔn)進(jìn)入從事科學(xué)研究觀測(cè)活動(dòng)。緩沖區(qū)外圍劃為實(shí)驗(yàn)區(qū),可以進(jìn)入從事科學(xué)試驗(yàn)、教學(xué)實(shí)習(xí)、參觀考察、旅游以及馴化、繁殖珍稀、瀕危野生動(dòng)植物等活動(dòng)。未來(lái)我們需要確定核心區(qū)、緩沖區(qū)、實(shí)驗(yàn)區(qū)的擾動(dòng)范圍和時(shí)間,結(jié)合動(dòng)植物多樣性觀測(cè)數(shù)據(jù)(如紅外相機(jī)捕捉的野生動(dòng)物活動(dòng)區(qū)域),開展保護(hù)區(qū)森林?jǐn)_動(dòng)對(duì)動(dòng)物棲息地的影響研究[33-34]。
本文基于光譜軌跡變化檢測(cè)算法對(duì)森林保護(hù)區(qū)森林?jǐn)_動(dòng)進(jìn)行了評(píng)估。利用GEE云計(jì)算平臺(tái)的LandTrendr算法和Landsat影像,對(duì)陜西省柴松自然保護(hù)區(qū)和太白山自然保護(hù)區(qū)的森林向非森林轉(zhuǎn)化進(jìn)行了檢測(cè)。對(duì)變化像元和穩(wěn)定像元的總體檢測(cè)精度達(dá)到93%,對(duì)擾動(dòng)年份的總體檢測(cè)精度達(dá)到89%。柴松保護(hù)區(qū)的擾動(dòng)年份主要發(fā)生在2008年,由人類活動(dòng)所致。太白山自然保護(hù)區(qū)的擾動(dòng)年份主要發(fā)生在2013年,由自然因素所致。柴松保護(hù)區(qū)擾動(dòng)面積為42.74 hm2,太白山保護(hù)區(qū)的擾動(dòng)面積為 23.68 hm2,柴松保護(hù)區(qū)的擾動(dòng)斑塊數(shù)目也明顯大于太白山保護(hù)區(qū)。本文研究結(jié)果可以幫助科研人員和管理人員了解當(dāng)前兩個(gè)保護(hù)區(qū)的森林生長(zhǎng)狀態(tài),同時(shí)為自然保護(hù)區(qū)生態(tài)環(huán)境保護(hù)評(píng)估提供基線信息。