基于Landsat時序數(shù)據(jù)的森林干擾監(jiān)測*

鐘 莉 陳蕓芝 汪小欽

(福州大學空間數(shù)據(jù)挖掘與信息共享教育部重點實驗室 衛(wèi)星空間信息技術(shù)綜合應用國家地方聯(lián)合工程研究中心數(shù)字中國研究院(福建) 福州 350108)

森林是陸地生態(tài)系統(tǒng)的重要組成部分，是地球上最大的碳儲存庫，其動態(tài)變化直接影響全球陸地森林生態(tài)系統(tǒng)的碳循環(huán)和碳蓄積(葛全勝等， 2008； 趙金龍等， 2013)。森林動態(tài)變化包括森林突發(fā)性變化和植被時序趨勢變化。森林干擾是森林突發(fā)性變化，根據(jù)干擾的起因，通?？煞譃槿藶楦蓴_和自然干擾，其形式主要表現(xiàn)為森林火災、森林采伐以及因受雨雪和病蟲災害等所引起的森林改變(Gowardetal., 2013)。近年來，隨著經(jīng)濟發(fā)展，全球各地森林植被頻繁遭到嚴重破壞，人為干擾(包括森林采伐、城市化、改林為耕)也逐漸增加(李維長， 2000)。

遙感技術(shù)的不斷進步，衛(wèi)星遙感數(shù)據(jù)作為時間尺度變化監(jiān)測的重要手段越來越受到研究者重視(Coppinetal.， 1996)。森林干擾研究的數(shù)據(jù)源通常分為基于區(qū)域和全球尺度遙感數(shù)據(jù)2類，其中，AVHRR、MODIS和SPOT VEG等數(shù)據(jù)空間分辨率較低，多用于大尺度的森林變化監(jiān)測(宋富強等， 2011； 王強等， 2012)； 相比而言，Landsat系列影像數(shù)據(jù)具有較高精度，多用于區(qū)域尺度的森林干擾監(jiān)測及其造成的森林碳通量變化研究(楊辰等， 2015)。美國地質(zhì)調(diào)查局于2008年免費開放了Landsat堆棧數(shù)據(jù)(Woodcocketal.， 2008)，該數(shù)據(jù)具有良好的穩(wěn)定性和持續(xù)性，基于Landsat時序數(shù)據(jù)的研究取得了長足發(fā)展(Townshendetal.， 2012； Wulderetal.， 2015)。以往研究中，森林變化監(jiān)測通常是對影像分類后進行兩兩比較(Cohenetal.， 2010)，不僅費時費力導致效率低下，而且針對長時序分析精度不高，不能滿足應用需求(楊辰等， 2013)；影像差異法利用差值、比值、典型相關(guān)性變化等方法構(gòu)建差異影像，雖然消除了分類后帶來的誤差積累，但是受輻射差異影響較大，除進行嚴格的輻射歸一化外，只能通過選擇同一傳感器和同一季相數(shù)據(jù)以減少噪聲(Deroseetal.， 2011； Maseketal.， 2008)。干擾監(jiān)測中，運用幾期影像很難區(qū)分突發(fā)干擾和背景噪聲，而時間序列分析中，可利用背景噪聲和突發(fā)變化的時間序列軌跡不同將其分開(Hostertetal.， 2003)。國外學者利用時間序列影像聯(lián)合分析，通過提取時序特征，進而識別森林干擾，大大提高了工作效率和監(jiān)測精度(Huangetal.， 2009； Jamalietal.， 2015； Kennedyetal.， 2010； Sulla-Menasheetal.， 2014； Verbesseltetal.， 2010)。

福建是我國南方重點林區(qū)之一，森林覆蓋率約66%，居全國首位。森林在全球碳循環(huán)和碳管理中占有重要位置(Faheyetal.， 2010)。本研究以福建省長汀縣為研究區(qū)，基于Landsat時序數(shù)據(jù)，采用時間序列軌跡分析(LandTrendr)方法進行森林干擾監(jiān)測，以期為陸地森林生態(tài)系統(tǒng)碳循環(huán)和碳蓄積以及氣候變化研究提供參考。

1 材料與方法

1.1 研究區(qū)概況

長汀縣地處福建省西部，25°18′40″—26°02′05″N，116°00′45″—116°39′20″E，南與廣東省近鄰，西與江西省接壤，呈四周凸起、中間下凹的盆地特征，地形破碎，嶺谷相間。該縣處中亞熱帶季風氣候區(qū)，平均氣溫18.3 ℃，無霜期年均260天，四季分明。花崗巖發(fā)育的紅壤抗蝕能力弱，原生植被多遭破壞，現(xiàn)有植被主要包括馬尾松(Pinusmassoniana)、灌叢和荒草坡在內(nèi)的次生植被和人工植被。研究區(qū)地理位置見圖1。

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

選用長汀縣2000—2016年15期Landsat影像數(shù)據(jù)(表1)，條帶號121/42，數(shù)據(jù)來源于美國地質(zhì)調(diào)查局(USGS)。所選影像盡量保持同一季相，以減少植被物候差異對光譜識別的影響。長汀縣夏季雨水豐富，為獲取高質(zhì)量影像，將時間集中在9、10月，個別年份由于云覆蓋等因素不滿足要求，將時間間隔延長至2年。數(shù)字高程資料來源于地理空間數(shù)據(jù)云平臺提供的全球30 m分辨率GDEM數(shù)字高程產(chǎn)品。

對每年影像進行輻射定標和集成于ENVI的基于MODTRAN4+輻射傳輸模型的FLAASH大氣校正，得到長汀縣2000—2016年時序地表真實反射率影像。為消除云和云陰影對結(jié)果的干擾，利用Fmask算法(Zhuetal.， 2012)逐一去除云像元。對掩膜后留下的空值采用最鄰近年份的無云像元時間內(nèi)插方法進行填充，對于時間序列影像的第1年和最后1年，使用最鄰近年份的無云像元替代：

xi=xp+(i-p)(xn-xp)/(n-p)。

(1)

式中：xi為第i年的填充值；xp和xn分別為第p年和第n年的指數(shù)值。

圖1 研究區(qū)地理位置

表1 遙感影像數(shù)據(jù)介紹

1.3 光譜指數(shù)選取

選用歸一化燃燒率(normalized burn ratio，NBR)構(gòu)建時序數(shù)據(jù)進行干擾監(jiān)測。國外學者(Kennedyetal.， 2010； Cohenetal.， 2010)基于大量數(shù)據(jù)研究表明，NBR對多種類型干擾事件具有最大的敏感性。NBR由對植被葉綠素含量較敏感的近紅外(NIR)和對水汽含量較敏感的短波紅外(SWIR)組合構(gòu)成，健康植被因高NIR值和低SWIR值，使其具有較高的NBR值，而受干擾的像元由于呈現(xiàn)土壤特征，SWIR值增大，使NBR值降低。據(jù)此，可根據(jù)時序軌跡中NBR值的突變識別森林變化信息：

(2)

1.4 時間序列軌跡擬合算法

LandTrendr算法(Kennedyetal.， 2010)的核心是將復雜的時序軌跡簡化為一組相連的線段，并對分割后的分段進行線性擬合，旨在消除噪聲，突出重要的信息。算法主要參數(shù)設(shè)置如表2所示，關(guān)鍵流程如下：

1) 時序軌跡提取及噪聲消除 逐像元提取時序軌跡(圖2a)，并通過設(shè)置噪聲值進一步去除預處理過程中殘留的細小噪聲(圖2b)。

2) 時序軌跡分割 對整個時間序列，首先將時序中的第1年和最后1年作為第1次分割的起始點，得到觀測年份與指數(shù)值的一階回歸方程【圖2c(1)】，并計算每個觀測年指數(shù)值與預測值的絕對差，選取絕對差最大的一年作為下次分段的分段點，據(jù)此將時間序列分為2部分【圖2c(2)】。然后對2部分分別計算回歸方程和每段的均方誤差(MSE)，選取MSE較大的一段參與下次分割，以此類推，直到分段數(shù)達到設(shè)定的初次分段數(shù)【圖2c(3)～(5)】。為防止分段數(shù)過多造成過擬合現(xiàn)象，選用角度閾值進一步判斷，使分段數(shù)達到最大分段數(shù)。通過計算對比原始時序軌跡中段與段的夾角，去除最“淺”角對應的分段【圖2c(6)】。

3) 時序軌跡擬合 聯(lián)合基于點與點連線及基于回歸連線2種方法進行時序軌跡擬合。對于第1個分段，分別計算2種方法的MSE，選取MSE較小的方法作為分段的連接線。從下一分段開始，每一個連接線必須與前一個相連，同樣計算2種方法的MSE，選取MSE較小的方法連接，最后得到一條完整的擬合軌跡。

4) 模型簡化 采用恢復率閾值(1/n)剔除不合理的分段點，對上一過程擬合的時序軌跡通過迭代方式簡化模型并重新計算擬合。在研究區(qū)內(nèi)，絕大多數(shù)森林干擾為突發(fā)性事件，但植被恢復是個緩慢的過程,通過設(shè)置恢復率閾值(1/n)，將植被恢復時間小于n年的森林干擾視為噪聲。根據(jù)長汀縣植被恢復狀況，將恢復率閾值設(shè)為0.5，恢復時間小于2年，即干擾發(fā)生1年后，植被恢復到干擾發(fā)生前狀態(tài)的區(qū)域視為噪聲。

5) 光譜濾波 通過回歸模型將NBR值轉(zhuǎn)化為植被覆蓋率進而生成干擾量，以植被變化強度閾值(表2濾波過程)為濾波條件，剔除原始分割算法中小的光譜分段點，以減少光譜異常信號或植被物候差異造成的時序軌跡過度擬合。

基于以上過程，可識別時間序列中NBR下降、上升和保持不變3種特性。經(jīng)過噪聲去除及對由光譜異常信號或植被物候差異造成的偽變化進行濾波后，NBR值下降即由森林干擾導致。

圖2 LandTrendr分割流程

表2 LandTrend參數(shù)設(shè)置

1.5 森林干擾制圖

制圖前需先排除水體、裸地等非森林地物?；贜BR指數(shù)的時序分割算法獲取頂點的濕度、亮度和綠度影像，計算每年每個像元的纓帽變換角(TCA)。TCA值越大，對應的森林覆蓋率越高，選擇平均值最高的TCA值，其對應的像元最像森林像元(Huangetal.， 2002)。通過TCA閾值提取森林像元，提取的像元在整個時序中被識別為持續(xù)森林。由于研究所用時間序列較短，可監(jiān)測到的二次干擾區(qū)域很少，故只提取最大干擾作為監(jiān)測結(jié)果。TCA計算公式為：

式中： TCG和TCB分為纓帽變換綠度和濕度分量。

2 結(jié)果與分析

2.1 森林干擾監(jiān)測結(jié)果

圖3為長汀縣2000—2016年森林干擾分布，圖中非森林和持續(xù)森林已被腌膜，以白色背景表示，其余不同顏色代表不同年份的干擾分布。

2.2 森林干擾時空分布特征

圖4為長汀縣森林干擾面積隨時間的變化特征。經(jīng)統(tǒng)計分析， 2000—2016年長汀縣森林干擾總面積192.49 km2，平均每年受擾動森林12.83 km2； 2001年干擾量最小，不足1 km2； 2004、2008和2009年森林受干擾較為嚴重，均在30 km2以上，約占當年森林面積的1.3%，共占干擾總面積的50%，其中2004年干擾面積高達32.85 km2； 2003、2006、2007和2010—2011年干擾面積略大于10 km2，其余年份干擾面積小于或遠小于10 km2。

圖4 長汀縣森林干擾面積隨時間的變化特征

從長汀縣氣象局獲取的數(shù)據(jù)得知， 2004年是長汀縣近20年來森林火災次數(shù)最多的一年， 2008和2009年也爆發(fā)了大面積森林火災(江帆等， 2015)。2008年受“拉尼娜”特大雨雪冰凍災害侵襲，南方森林資源損失慘重(曹坤芳等， 2010)。2002—2004年、2006—2008年干擾面積呈持續(xù)上升趨勢， 2001和2005年干擾降低，特別是從2009年起干擾面積持續(xù)減少直到2016年近乎最低。2000—2016年森林干擾面積在個別年份波動較大，但總體上隨時間呈下降趨勢。

對每年影像進行森林提取，獲得長汀縣歷年森林面積及所占土地總面積比例。由圖5可知， 2004、2008和2009年森林面積下跌較大， 2009年起呈緩慢增長趨勢，與擾動情況相呼應。對比圖4和圖5，2004、2008和2009年的森林干擾面積占當年森林面積的1.3%左右，其余各年均小于當年森林面積的0.6%。

圖5 長汀縣歷年森林面積及所占總面積比例

由圖6可知，長汀縣森林干擾面積隨持續(xù)時間呈極速下降趨勢，后保持低值平緩，說明干擾持續(xù)時間保持較短的年數(shù)，主要為1～2年，少數(shù)為3年，其中發(fā)生在1年的干擾面積比例最大，達82%； 發(fā)生在2和3年的干擾面積分別占13%和4%。結(jié)合Google影像目視解譯，長汀縣森林干擾主要是由森林火災和人工砍伐造成的急劇干擾事件。

圖6 長汀縣森林干擾面積隨持續(xù)時間變化分布

由圖3各年份干擾分布可知，長汀縣森林干擾面積具有較強的空間分布特征。2004、2008和2009年森林干擾斑塊面積較大，主要集中在長汀縣中東部非森林區(qū)域附近，受人為干擾影響。2008年除分布在非森林區(qū)域附近的干擾外，有大量因受雨雪災害影響的細小斑塊分布在長汀縣西部。其他年份亦如此，主要分布在長汀縣非森林區(qū)域附近，四周分布零散細小的斑塊。圖7清晰展現(xiàn)了長汀縣森林干擾面積在各海拔段的分布情況，隨海拔升高，干擾呈明顯下降趨勢，其中海拔380～460 m之間的干擾面積最大，超過60%的干擾發(fā)生在中低海拔地區(qū)，海拔大于940 m的干擾僅占0.45%。

圖7 長汀縣森林干擾面積隨海拔變化分布

2.3 基于地面調(diào)查的精度評估

圖8 長汀縣森林干擾實地調(diào)查

結(jié)合森林干擾監(jiān)測結(jié)果和地面調(diào)查資料對4個典型樣區(qū)(樣區(qū)分布見圖3)進行驗證，獲取干擾發(fā)生的時間和位置。圖8a-d從左到右分別為樣區(qū)干擾年TM影像、監(jiān)測結(jié)果、實地調(diào)查照片和NBR指數(shù)時序軌跡，實地調(diào)查照片均拍攝于2013年7月。當干擾發(fā)生時，NBR指數(shù)大幅度降低，之后逐漸增加，植被慢慢恢復。經(jīng)目視解譯并對照提取結(jié)果可知，干擾斑塊可被完整提出，且邊界準確清晰，細小干擾也能逐一識別。圖8a、b干擾分別發(fā)生在2010和2008年，到2013年植被有所恢復但不是特別理想，從實地調(diào)查照片可知植被主要以灌木為主； 從圖8c、d可以看出，植被恢復效果良好。實地調(diào)查結(jié)果與研究得到的干擾監(jiān)測結(jié)果一致，表明該方法適用于長汀縣森林干擾監(jiān)測，且具有較好的監(jiān)測能力。

2.4 基于像元尺度的精度驗證

采用分層隨機抽樣方法，以像元為單位提取森林干擾樣本進行精度驗證。選用每年作為一層(包括15期森林干擾區(qū)域，森林和非森林各為1層)，共17層。以最少像素層30為基準，對森林干擾層進行等比例隨機抽樣，15期干擾樣本數(shù)量共9 075個，森林和非森林像元數(shù)偏多，適當減少，分別為3 200和1 000個，即共13 275個像元樣本參與驗證。驗證結(jié)果統(tǒng)計為混淆矩陣，見表3。

由表3可知，時間序列軌跡算法可以很好提取長汀縣2000—2016年森林干擾，總體精度達96.26%，Kappa系數(shù)為0.92，具有較高的監(jiān)測精度。各年份用戶精度均在80%以上， 2014和2016年相對其他年份出現(xiàn)一定漏提。2016年處于時序軌跡的末端且前1年影像缺失，加之目視解譯發(fā)現(xiàn)該年干擾范圍和斑塊面積較小，分布零星，導致2016年出現(xiàn)45%左右的漏提。

3 討論

LandTrendr算法根據(jù)數(shù)據(jù)本身決定變化軌跡，利用干擾與背景噪聲不同，通過參數(shù)控制進行分割擬合，消除時序中的噪聲和小的光譜變化點，以突出重要的信息，最終得以識別干擾。該方法不僅可提取多種類型的干擾，還可監(jiān)測同一區(qū)域多次擾動變化，同時獲取植被干擾和恢復相關(guān)信息。

由于云和云陰影的影響，不能保證每年都有滿足要求的高質(zhì)量影像。某些年份影像缺失，可能導致干擾監(jiān)測時間延誤1年，或某些區(qū)域出現(xiàn)植被恢復過快現(xiàn)象導致輕微漏提。后續(xù)工作可重構(gòu)高時間分辨率影像，降低因影像缺失造成的漏提、誤提現(xiàn)象。經(jīng)目視解譯可知，長汀縣森林干擾類型主要為森林火燒和人工采伐，本研究所用方法雖可提取多種類型的干擾，但未能劃分出具體的干擾類型，后續(xù)工作可結(jié)合其他技術(shù)方法分別提取出火燒跡地和采伐跡地。

4 結(jié)論

針對我國森林干擾頻繁、干擾類型復雜多樣的特點，采用LandTrendr時序軌跡分析方法研究福建省長汀縣森林干擾狀況，并進行干擾時空特征分析。結(jié)論如下：

1) 2000—2016年，長汀縣森林干擾面積在個別年份波動較大，但總體上隨時間呈下降趨勢。

2) 長汀縣森林干擾主要是由森林火災和人工砍伐造成的急劇干擾事件，且主要發(fā)生在非森林區(qū)域附近的低海拔地區(qū)。

3) 基于目視解譯和實地調(diào)查結(jié)果與研究得到的干擾監(jiān)測結(jié)果一致，干擾斑塊可被完整提出，且邊界準確清晰，細小干擾也能逐一識別； 基于像元尺度精度驗證的總體精度達96.26%，Kappa系數(shù)為0.92，各年份用戶精度均在80%以上，除個別年份外，生產(chǎn)者精度均在75%以上，具有較高的監(jiān)測精度。