山東省植被覆蓋的時空變化分析

段晨陽，王曉艷

（1.西安科技大學 測繪科學與技術(shù)學院，陜西 西安 710054）

植被作為生態(tài)系統(tǒng)的平衡調(diào)節(jié)器和生態(tài)變化的指示器，能夠為生態(tài)變化、環(huán)境演化研究提供科學依據(jù)。近年來基于時間序列的衛(wèi)星遙感植被覆蓋度變化分析逐漸增多，為植被監(jiān)測提供了數(shù)據(jù)支撐，歸一化植被指數(shù)（NDVI）對植物檢測靈敏度較高，能表示地表綠度和植被覆蓋特征，可廣泛應用于大尺度地表植被活動監(jiān)測與評估[1-3]。國內(nèi)外學者對植被覆蓋變化采用不同方法進行研究，多是關(guān)于植被變化及產(chǎn)生因素、植被動態(tài)變化監(jiān)測、物種分布等[4-13]?？臻g自相關(guān)是地理空間變量分布與相鄰位置間相關(guān)性的研究，通過監(jiān)測分析區(qū)域的變異是否依賴于臨近位置上的變異來判別是否存在空間相關(guān)性[8]，目前主要應用于經(jīng)濟[9]、人口[10]、土地利用[11]等方面，也有一些學者利用空間自相關(guān)分析植被覆蓋在空間上的相關(guān)性[12-13]。MODIS NDVI在區(qū)域小比例尺范圍內(nèi)植被覆蓋動態(tài)變化分析上具有較好的應用。本文以山東省為研究區(qū)，進行植被覆蓋的動態(tài)變化以及空間自相關(guān)分析，首先對山東省植被覆蓋進行等級劃分，將時間序列變化的MODIS影像與空間自相關(guān)分析相結(jié)合，利用一元線性回歸方法，擬合每一個像元的變化趨勢，并結(jié)合植被覆蓋重心遷移路徑分析，以揭示植被覆蓋在時間和空間上的變化規(guī)律，更好地表示植被覆蓋的動態(tài)演變規(guī)律以及植被覆蓋的空間聚集程度，為山東省生態(tài)環(huán)境保護工作提供科學依據(jù)[14]。

1 研究區(qū)概況

山東省地處中國東部沿海地區(qū)，位于34°22′～38°23′N、114°19′～122°43′E之間，地域遼闊，受溫帶季風氣候影響。圖1為山東省DEM圖。

圖1 研究區(qū)DEM圖

2 數(shù)據(jù)源與研究方法

2.1 數(shù)據(jù)來源及預處理

本研究采用美國地球資源觀測系統(tǒng)數(shù)據(jù)中心提供的植被指數(shù)產(chǎn)品（MOD13Q1），空間分辨率為250 m、時間分辨率為16 d。選取2008—2020年植被旺盛季節(jié)[15]即每年6—9月份影像，合計85景。利用MRT軟件對圖像進行NDVI提取，投影到Albers投影坐標系統(tǒng)下。為了去除云、大氣以及太陽高度角等的影響，采用國際通用的最大值合成法（MVC）對數(shù)據(jù)進行時間尺度擴展計算[16]，利用山東省矢量邊界對圖像進行裁剪。

利用ArcGIS地圖代數(shù)工具對影像進行異常剔除以及歸一化處理[0，1]，利用焦點統(tǒng)計算法，選擇用于計算統(tǒng)計信息的每個單元格周圍區(qū)域的形狀為矩形，統(tǒng)計類型為平均值。隨機生成點文件，將焦點統(tǒng)計后的影像進行重分類，通過提取值到點工具，將柵格數(shù)據(jù)的NDVI值批量提取到對應點的屬性表中，以便后期空間相關(guān)性分析。

采用Anselin等開發(fā)的Geoda軟件進行空間自相關(guān)分析，得到山東省2008—2018年NDVI全局Moran’s I指數(shù)，局部自相關(guān)分析得到空間自相關(guān)的LISA（LISA）聚集圖。

2.2 研究方法

采用MVC法進行NDVI時間序列計算。

式中，NDVIi為第i年NDVI最大值；NDVIt為各像元植被16 d合成NDVI最大值；d為第i年16 d合成NDVI影像總數(shù)目。

采用一元線性回歸方法對每個像元NDVI值進行時間序列變化分析，利用最小二乘法擬合像元不同年份的斜率來反映區(qū)域植被覆蓋變化的時空特征[15]。

式中，θslop為NDVI值變化斜率；n為時間序列值；i=1，2，3…，1為年號；Ni為第i年NDVI值。當斜率為正值時，NDVI增加，表示植被覆蓋度增加，反之減少。

重心可以表示地理要素的時空分布特征，移動方向指向高密度區(qū)域，移動距離反映要素變化的空間差異程度。將研究區(qū)域劃分為若干子區(qū)域，每個子區(qū)域的重心坐標為（Xi，Yi），Mi為該子區(qū)域某種意義下的重量，(Xˉ，Yˉ)為子區(qū)域重心坐標在重量影響下的加權(quán)平均值。本文以植被NDVI作為研究區(qū)的重量，探索植被覆蓋在時空變化下的規(guī)律。

應用空間自相關(guān)方法對山東省植被覆蓋空間分布進行分析研究。通常采用全局和局部兩種指標來衡量，全局指標是用于探測某現(xiàn)象在整個研究區(qū)域的空間分布模式，分析其是否有聚集特性存在；局部指標是反映整個大區(qū)域中，局部小區(qū)域單元上的屬性值與相鄰局部小區(qū)域單元上同一屬性值的相關(guān)程度[17]。本文運用應用較為廣泛的Moran’s I指數(shù)進行研究。

首先確定空間權(quán)重矩陣，通常用一個二元對稱空間權(quán)重矩陣W來表示n個區(qū)域位置的鄰近關(guān)系[18]，如公式（4）：

式中，n為空間單元個數(shù)；wij為區(qū)域i和的鄰接關(guān)系。全局Moran’s I指數(shù)為（5）：

式中，xi、xj為空間位置i和j的觀察值，全局Moran’s I的取值范圍為[-1，1]，樣本值越趨近于1空間差異性越小，反之空間差異性越大。正值表示樣本值趨于聚集，負值表示樣本值在空間上的聚集性逐漸減小或者消失，0表示樣本值在空間上相互獨立，沒有相關(guān)性。

對于全局Moran’s I指數(shù)，可以用標準化統(tǒng)計量Z來檢驗，Moran’s I指數(shù)檢驗的標準化統(tǒng)計量為：

式中，以正態(tài)分布95%置信區(qū)間雙側(cè)檢驗閾值1.96為界限[6]，E(Ii)和Z(Ii)分別為理論期望和理論方差。當Z(Ii)＞1.96時，空間關(guān)系為正相關(guān)；當Z(Ii)＜-1.96時，空間關(guān)系為負相關(guān)；當Z(Ii)介于-1.96～1.96之間時，觀測值呈隨機分布。

全局空間自相關(guān)假定空間同質(zhì)，但區(qū)域要素之間普遍存在異質(zhì)性，相同變量在不同空間尺度上的空間自相關(guān)程度還是有差別的，區(qū)域空間自相關(guān)可以更好地解釋植被空間異質(zhì)特性。

局部Moran’s I公式為：

3 結(jié)果與分析

3.1 NDVI的空間分布

NDVI是植物生長狀態(tài)以及空間分布密度的最佳指示因子，與植被分布密度呈線性關(guān)系。將NDVI指數(shù)劃分為5個級別[2]：小于0.1的區(qū)域為非植被區(qū)；0.1～0.3為植被覆蓋匱乏地區(qū)；0.3～0.6為植被覆蓋較匱乏地區(qū)；0.6～0.8為植被覆蓋良好地區(qū)；大于0.8為植被覆蓋優(yōu)秀地區(qū)。山東省植被類型主要以農(nóng)業(yè)為主，主要分布在魯西北平原。落葉闊葉林、針葉林、灌叢、草地等主要分布在魯中南山地區(qū)陵、膠東丘陵地區(qū)。圖2為山東省2008—2020年的NDVI空間分布圖，近13 a來NDVI較高的地區(qū)在魯西南及魯西北等農(nóng)業(yè)區(qū)，NDVI均在0.6以上。中部山地地區(qū)以及半島丘陵地區(qū)NDVI指數(shù)處適中狀態(tài)。NDVI指數(shù)較低區(qū)域為魯西南地區(qū)的南四湖及東平湖、北部的黃河三角洲以及萊州灣地區(qū)以及城鎮(zhèn)建設(shè)地區(qū)。由此可見，植被覆蓋的空間分布受到土地利用類型的影響。

圖2 NDVI空間分布

3.2 NDVI時間變化

時間序列的遙感影像圖能夠反映植被覆蓋的變化情況，采用一元線性回歸分析對每一個柵格NDVI的變化趨勢進行模擬，斜率代表了該柵格的NDVI在2008—2018年間的變化趨勢。若斜率為正值，表示植被呈增長趨勢；反之呈減少趨勢。

利用MATLAB處理計算，斜率最大值為0.102，最小值為-0.104。依據(jù)斜率的大小情況，劃分為5個等級[19]，-0.104～-0.02為顯著減少；-0.02～-0.005為輕微減少；-0.005～0.005為保持不變；0.005～0.02為輕微增加，0.02～0.102為顯著增加。由表1可知植被輕微增加的像元數(shù)最多占62.5%，相對來說輕微減少的像元所占比重要最小，僅為2.1%，得出近年來山東省植被呈現(xiàn)出穩(wěn)定和增長趨勢。由圖3可見，植被變化趨勢有明顯的空間分布差異，大部分地區(qū)呈穩(wěn)定增長趨勢。其中植被顯著增長地區(qū)主要集中在聊城、德州平原地區(qū)以及濰坊等農(nóng)業(yè)區(qū)，由于實施鄉(xiāng)村振興戰(zhàn)略，農(nóng)業(yè)得到很大發(fā)展，糧食產(chǎn)量創(chuàng)新高，造林綠化使得植被覆蓋增加。而北部的黃河三角洲及萊州灣、湖泊濕地及沿海部分地區(qū)有輕微導致植被覆蓋減少，而湖泊濕地周圍植被極不穩(wěn)定，受環(huán)境人為因素干擾大，出現(xiàn)減少趨勢。

表1 各變化等級統(tǒng)計

圖3 植被覆蓋變化趨勢

3.3 重心演化分析

為揭示山東省植被覆蓋在時間以及空間上的動態(tài)變化，計算2008—2020年時間序列變化的植被空間重心坐標。如圖4所示，取間隔2 a的NDVI重心坐標生成散點圖，箭頭的指向即為重心遷移變化方向。由圖可知，2008—2020年間植被覆蓋的重心發(fā)展方向經(jīng)歷幾個過程，即東北方向、西南方向、西北方向，再到西南方向。結(jié)合圖3植被變化趨勢作對比分析，2008—2012年植被覆蓋重心坐標往東北偏移，濰坊地區(qū)大力發(fā)展農(nóng)業(yè)，打造新型蔬菜大棚，著力于蔬菜出口貿(mào)易，使得植被覆蓋持續(xù)增加。2012—2014年植被覆蓋重心坐標往西南偏移，菏澤、濟寧地區(qū)位于平原地區(qū)，農(nóng)業(yè)大力發(fā)展，使得植被覆蓋增加。2014—2018年植被覆蓋重心坐標往西北偏移，聊城、德州地區(qū)農(nóng)業(yè)生產(chǎn)得到迅速發(fā)展，糧食產(chǎn)量大大提高，黃河三角洲地區(qū)實施生態(tài)保護政策，植被覆蓋度大大提高。2018—2020年往西南方向發(fā)展。植被覆蓋重心的遷移變化受到很多因素的影響，可能還受溫度降水等綜合影響。

圖4 NDVI重心遷移路徑

3.4 空間自相關(guān)分析

利用Geoda軟件進行全局自相關(guān)分析，用Weight工具創(chuàng)建研究區(qū)空間權(quán)重矩陣文件，依照Threshold Distance規(guī)則創(chuàng)建空間權(quán)重矩陣，得到全局Moran’s I為0.898，通過了P＜0.05的顯著性檢驗。Moran’s I系數(shù)反映了空間鄰近點NDVI的相似程度[20]，隨機點集中分布在第三象限，表明山東省表現(xiàn)為顯著的正空間自相關(guān)，植被覆蓋呈現(xiàn)出聚集狀態(tài)。山東省大部分地區(qū)植被類型為1年一熟或者2年三熟的農(nóng)作物以及草地林地過渡區(qū)域，植被覆蓋季節(jié)性變動態(tài)化小，植被覆蓋穩(wěn)定性強[21]，這些因素對于植被覆蓋在空間上的顯著相關(guān)性產(chǎn)生一定影響。

全局空間自相關(guān)分析是對整個研究區(qū)內(nèi)的地理空間要素整體綜合分析，來探求整體表現(xiàn)出來的規(guī)律，是一個整體的反映。但是對于內(nèi)部來說，各個區(qū)域的局部空間自相關(guān)性并不是完全一樣的，經(jīng)常表現(xiàn)出不同程度的空間異質(zhì)性，因此對研究區(qū)進行局部空間自相關(guān)分析，能夠發(fā)現(xiàn)全局分析沒有發(fā)現(xiàn)的規(guī)律，從而更好的揭示空間聚集的一般規(guī)律。由Moran’s I的數(shù)值看出NDVI在空間上表現(xiàn)出較強的正相關(guān)，但是沒有清晰地表明相關(guān)性在空間中的分布情況，而LISA圖能夠很好地彌補這一不足。

圖5為通過P＜0.05顯著性檢驗的NDVI LISA聚集圖，除北部的黃河三角洲和萊州灣、南四湖、沿海以及零星地區(qū)NDVI表現(xiàn)為低-低自相關(guān)以外，山東省大部分地區(qū)呈高-高自相關(guān)趨勢。山東省為農(nóng)業(yè)大省，植被類型大多以農(nóng)用地為主[22]，平原面積占全省面積的55.56%，主要分布在魯西北地區(qū)和魯西南局部地區(qū)，所以這些地區(qū)植被覆蓋總體較好，高-高聚集現(xiàn)象明顯；魯中南為山地區(qū)陵地區(qū)，主要為森林、灌叢和草地等，植被覆蓋總體良好，聚集現(xiàn)象不那么強；魯東沿海地區(qū)，旅游業(yè)發(fā)展迅速，開發(fā)建設(shè)導致植被覆蓋整體較差，低-低聚集明顯；黃河三角洲以及萊州灣地區(qū)一直為NDVI低值聚集區(qū)，距離黃河河道越遠，渤海海岸線越近，NDVI值越小。近海地區(qū)常年受到黃河水沙及人類活動的影響，土壤含鹽量較高，生態(tài)脆弱，植物種類匱乏，覆蓋率低，低-低聚集現(xiàn)象明顯。

圖5 山東省NDVILISA圖

4 結(jié)語

本文以2008—2020年的MODIS影像為研究基礎(chǔ)，對山東省植被覆蓋的時間序列的變化情況，以及植被覆蓋在空間上的關(guān)聯(lián)情況做出研究，得到了以下結(jié)論：

1）近年來，山東植被覆蓋整體表現(xiàn)為穩(wěn)定上升的趨勢，75%的區(qū)域保持增長趨勢，呈增加趨勢的區(qū)域遠遠大于呈減退趨勢的區(qū)域。

2）植被變化趨勢在空間上具有一定差異性，西部地區(qū)為增長型，東部地區(qū)為穩(wěn)定性。植被覆蓋重心往魯西北、魯西南農(nóng)業(yè)區(qū)偏移。

3）全局Moran’s I均值為0.898，顯著性檢驗P＜0.05，置信度較高，該結(jié)果表明山東省植被覆蓋有較強的空間聚集性。局部空間自相關(guān)分析來看，山東省大部分地區(qū)均表現(xiàn)為高-高自相關(guān)，植被覆蓋度較好；魯東沿海地區(qū)、黃河三角洲和萊州灣地區(qū)以及各地零星地區(qū)表現(xiàn)為低-低自相關(guān)，植被覆蓋度較差。

基于時間序列的遙感影像圖結(jié)合空間自相關(guān)方法，對山東省植被覆蓋進行研究，能夠客觀地揭示植被覆蓋動態(tài)變化規(guī)律以及空間聚集情況。統(tǒng)計學方法的數(shù)據(jù)指標能夠反映一定的客觀規(guī)律，但是植被覆蓋的變化與諸多因素有關(guān)。在以后的研究工作中，應綜合多種因素更深入的進行探究。