2004-2015年賀蘭山自然保護區植被NPP時空變化與氣候響應

李婷婷, 馬 超,2, 郭增長

(1.河南理工大學 測繪與國土信息工程學院, 河南 焦作 454003; 2.河南理工大學 自然資源部礦山時空信息與生態修復重點實驗室, 河南 焦作 454003; 3.河南測繪職業學院, 鄭州 451464)

聯合國大會于2015年9月25日一致通過了題為“改變我們的世界：2030年可持續發展議程”的第70/1號決議，將保護、恢復和促進生態系統可持續發展(目標15)列入全球17項可持續發展目標[1]。目標15具體側重于可持續地管理森林、恢復退化的土地和成功防治荒漠化、減少退化的自然生境和結束生物多樣性喪失。基于遙感技術的長時序變化監測和評估研究，可明確生態變化趨勢，提升地方政府實施、監督和評估可持續發展目標的能力，為“中國落實2030年可持續發展議程進展報告”提供數據參考[2]。

在研究生態系統對全球變化的影響、響應和對策中，植被凈初級生產力成為一項不可或缺的指標和核心內容[3]。植被凈初級生產力(net primary productivity,NPP)是綠色植物總初級生產力(在單位時間和單位面積上通過光合作用產生的全部有機物同化量)扣除自養呼吸后的剩余部分[4]。NPP表示植物光合作用產物固定和轉化的效率，直接反映植物在自然環境條件下的生產能力，是評價生態系統結構與功能協調及生物圈人口承載力的狀況指標[5]。祝萍等[6]將植被NPP作為評價指標對中國典型國家級自然保護區生境狀況時空變化特征進行研究。張鐿鋰等[7]將青藏高原自然保護區的植被NPP變化過程數據與其相鄰等面積區域的植被NPP變化差異進行對比分析，評估自然保護區的保護成效。王培娟等[8]針對長白山自然保護區植被NPP主要影響因子的敏感性進行研究。可見在對自然保護區生態系統的變化檢測和評估研究方面，植被NPP已成為一項重要指標。

中國境內山脈眾多，其中不乏名山，但能以大尺度地理界線聞名者寥寥無幾，賀蘭山幸為其一[2]。賀蘭山是中國北方干旱地理帶的一片綠洲，是寧夏平原的氣候屏障和區域生態環境變化的指示器，其生態走向關乎瀕危物種的延續與北方干旱荒漠帶寧夏段的生態平衡。植被作為干旱區敏感的環境因素，植被覆蓋變化的研究已成為揭示干旱區自然環境變化及其規律的重要手段[9]。周夢云等[10]對寧夏賀蘭山國家級自然保護區建立前后區域生態脆弱性時空格局變化進行研究。劉勝濤等[11]對寧夏賀蘭山自然保護區森林生態系統凈化大氣的環境功能進行研究。顧延生等[12]對賀蘭山中段植被類型及其覆蓋變化進行研究。通過文獻調研發現，以往的研究多集中于對賀蘭山植被進行某些特定生理特征的研究，而沒有對賀蘭山地區的植被生物量時空變化同氣候變化相結合，獲得區域生態演化規律。

由于賀蘭山自然保護區的局部山地效應，使對該區氣候與植被的相關研究變得意義非凡。研究賀蘭山自然保護區植被凈初級生產力時空變化及其與氣候因子的響應關系，可為地方政府提供區域植被生理狀況和生產能力的相關信息，對荒漠化生態治理、自然保護區的管理、生態狀況評價和促進社會經濟可持續發展提供參考依據。

1 研究區概況

寧夏賀蘭山為西北為數不多的森林分布區之一，既是荒漠區與荒漠草原區的分界線，又是干旱區與半干旱區的分界線，具有涵養水源和調節氣候的作用[11]。由于其地理位置的特殊性和生態系統的脆弱性，早在20世紀80年代就已被列為國家級生態保護區。(賀蘭山的管理機構始建于1950年。1982年7月1日，寧夏人大劃定賀蘭山為省級自然保護區，1988年5月國務院批準寧夏賀蘭山自然保護區晉升為國家級自然保護區。)

以1 500 m等高線為基準，并向外延伸1 km來確定研究區的邊界范圍(圖1A)。研究區地處寧夏賀蘭山自然保護區核心區，近南北走向，西鄰阿拉善高原，東接銀川平原。地理坐標為北緯38°07′—39°32′，東經105°17′—106°40′，南北長157.63 km，東西寬120.40 km，面積5 769.25 km2，海拔為1 168～3 497 m。

寧夏賀蘭山是我國東西部氣候的分界，屬季風氣候向大陸性氣候的過渡，也是我國200 mm等降水量線所經過的地區[12]。年均氣溫6～8℃，寒暑變化強烈，年均降水176.4～313.0 mm，全年干旱少雨。

賀蘭山因海拔高，地勢東緩西陡，導致自然氣候產生顯著垂直分異，構成獨特的干旱區山地植被垂直帶，孕育出荒漠草原，耐旱喬木、灌木，油松、山楊林和青海云杉等植被類型[13]，利用中國1∶100萬植被類型空間分布圖將研究區的植被類型進行合并，得到寧夏賀蘭山植被類型空間分布(圖1B)。

圖1 寧夏賀蘭山的地形與植被分布

2 數據與方法

2.1 研究數據

2.1.1 數據來源與處理 表達植物在自然環境下的生產能力的植被凈初級生產力NPP數據，源于美國一級大氣存檔和分布式數據發布中心(Level-1 and Atmosphere Archive & Distribution System Distributed Active Archive Center, LAADS DAAC)，研究選取2004—2015年的MODIS (Moderate Resolution Imaging Spectro radiometer)數據產品：MOD17A3H NPP數據(https:∥ladsweb.modaps.eosdis.nasa.gov)，其空間分辨率為500 m，主要用于研究NPP的變化與氣候響應的關系。該NPP數據是基于MODIS/TERRA衛星的遙感參數，通過利用BIOME-BGC模型與光能利用率模型建立的NPP估算模型模擬得到全球陸地植被凈初級生產力年際變化數據集，目前已在全球不同區域的植被生長狀況、生物量估算、環境變化檢測等研究中得到驗證和廣泛應用[14-15]。利用MODIS數據重投影工具MRT (Modis Reprojection Tool)將NPP數據定義為UTM投影，WGS_84橢球的地理數據，用于分析2004—2015年賀蘭山植被NPP變化特征。

高程數據源于美國太空總署(NASA)和國防部國家測繪局(NIMA)聯合測量的航天飛機雷達地形測繪使命數字高程模型(SRTM3 DEM,http:∥srtm.csi.cgiar.org/SELECTION/inputCoord.asp)，水平分辨率為90 m，主要用于以1 500 m等高線來劃定研究區范圍。

氣象數據源于資源環境數據云平臺(http:∥www.resdc.cn/Default.aspx)，分別為2004—2015年年降水量空間插值數據集，用于分析降水量對植被NPP的影響；2004—2015年年平均氣溫空間插值數據集，用于分析氣溫對植被NPP的影響。

植被類型空間分布數據源于資源環境數據云平臺(http:∥www.resdc.cn/Default.aspx)1∶100萬中國植被圖，通過ArcGIS對研究區植被類型空間分布圖進行投影變換，將植被類型合并為闊葉林、針葉林、灌叢、草地和荒漠這5種，為植被NPP閾值分割提供基礎參考。

2.1.2 數據質量驗證 研究區屬國家級自然保護區，大范圍實地調查難以實施，加之研究時間序列長(12 a)，生物量分布的時間和空間異質性大，通過有限的實地調查所獲得的平均生物量來推算整個研究區的生物量也會存在較大誤差。針對類似情況，已有研究者[16]將NPP估算結果與土地利用分類數據產品進行對比，通過局部土地利用變化趨勢來驗證整體精度，具有片面性。本文在研究區均勻選擇120個樣本點，利用遙感數據產品MOD13A1，通過將樣本點12 a的平均NPP與平均NDVI做對比來驗證精度，檢驗結果顯示二者具有非常高的擬合度(r=0.944，p<0.01)，MOD17A3H NPP數據質量可靠。

2.2 研究方法

2.2.1 偏差分析 平均NPP是基于像元對多年NPP數據進行計算求得的平均值，用于得到研究區多年的平均NPP空間分布圖，其計算公式如下：

(1)

偏差反映某時間段內NPP偏離多年NPP均值的程度，即某段時間內NPP的盈虧情況[17]。其計算公式如下：

(2)

2.2.2 相關性分析 相關性分析[18]可以反映要素之間的相關程度和相關方向，本文采用Pearson相關系數法對植被年NPP與氣溫、降水量的響應關系進行研究，其計算公式如下：

(3)

2.2.3 趨勢分析 一元線性回歸法[17]可以對研究數據進行趨勢分析，得到其變化速率，本文采用一元線性回歸分析法對研究區總NPP、氣溫和降水量進行研究，其計算公式如下：

Yi=β0Xi+β1

(4)

式中：Yi為因變量；Xi為自變量；β0為斜率；β1為一隨機變量。

(5)

式中：β0為變化速率；n為年數(時間序列為2004—2015年，即n=12)；k為時間序列；Xk為第k年的自變量。

3 結果與分析

3.1 MODIS NPP分類處理

通過對2004—2015年的NPP數據進行計算，得到賀蘭山12 a的平均NPP空間分布圖(圖2)。依據賀蘭山自然保護區主要植被類型的空間分布情況和NPP數值的對應關系，將研究區劃分荒漠、草原、灌叢、針葉林、闊葉林5個等級，其分別對應NPP閾值為[0,65]，(65,85]，(85,150]，(150,250]，(250,500]，單位為g C/(m2·a)[13]。

12年的NPP均值分割圖顯示出賀蘭山自然保護區生態系統的景觀結構，植被NPP由外向內逐漸遞增，即從低山淺山的荒漠草原帶，到中低山區的灌叢林帶，再到高海拔山區的針闊混交林帶，植物按海拔垂直分布的層次明顯，構成了當地典型的山地生態景觀。其景觀辨識度優于1∶100萬中國植被覆蓋專題圖。

圖2 2004-2015年平均NPP密度分割圖

3.2 MODIS NPP時間變化分析

研究對賀蘭山自然保護區2004—2015年植被初級凈生產力(NPP)進行了統計分析(表1)。

表1 植被凈初級生產力(NPP)統計

NPP區間在[0,65]的荒漠面積占比總體呈波動性減少的趨勢，荒漠面積在2005年、2011年和2015年荒漠面積達3 500 km2以上，占研究區總面積的60%左右。12 a間荒漠面積占比的最大變幅為35.13%；NPP區間在(65,85]的草原面積占比總體呈波動性增加的趨勢，草原面積在2010年達到961 km2，占研究區面積的16.66%。12 a間草原面積占比的最大變幅達9.65%；NPP區間在(85,150]的灌叢面積占比總體呈波動增加的趨勢，灌叢面積在2012年達到最大值1 832.75 km2，占研究區面積的31.77%。12 a間灌叢面積占比的最大變幅達19.99%；NPP區間在(150,250]的針葉林面積占比較為平穩，約768.75 km2，占研究區面積的13.32%。12 a間針葉林面積占比的最大變幅為5.82%；NPP區間在(250,500]的闊葉林面積占比最少，且波動變化較大，整體呈增加趨勢。闊葉林面積達225 km2以上，占研究區總面積的4%左右。12 a間闊葉林面積占比的最大變幅為4.33%。

2004—2015年研究區植被NPP均值為71.22 g C/(m2·a)，同期全國植被年平均NPP值273.5 g C/(m2·a)[19]，研究區植被NPP均值僅為同期全國植被年平均NPP值的26.04%；植被總NPP變化范圍為0.324 4～0.501 5 Tg C/a，波動幅度達154.6%，總體在變化區間0.177 1 Tg C/a以內呈緩慢增加趨勢；植被總NPP增長速率為0.281 6 g C/(m2·a)，遠低于同期全國植被NPP增長率1.141 5 g C/(m2·a)[19]。偏差分析得到總NPP偏差值在12 a中有7 a低于平均水平，出現生態虧缺的年份多于生態盈余。

3.3 MODIS NPP空間變化分析

3.3.1 蘇峪口國家森林公園剖面水平地帶性變化 蘇峪口國家森林公園是賀蘭山自然保護區植被覆蓋較好的區域，區內設有賀蘭山森林生態系統定位觀測研究站，研究其生態剖面對認識整個自然保護區生態演化具有代表性。以蘇峪口國家森林公園核心區中點(東經105°54′60.00″，北緯38°42′53.96″)為原點，設置W—E和N—S兩個剖面，分析核心區兩個地理方位的植被生產力變化，可以形成二維的時空變化暈渲圖(圖3)，圖中黑色實橫線位置為蘇峪口國家森林公園核心區中點位置所在，在W—E剖面上像素點的位置為28 pixel，在N—S剖面上像素點的位置為77 pixel。

圖3 蘇峪口國家森林公園NPP剖面暈渲圖

在W—E剖面時序變化圖中(圖3A)，蘇峪口國家森林公園植被凈初級生產力在水平地帶性分布上東部要優于西部，其中W—E剖面東部NPP總量為0.014 5 Tg C，占比66.55%，西部NPP總量為0.007 3 Tg C，占比33.45%。由于地理分界的原因，西部為大面積的荒漠和草原所占據，東部淺山地帶有部分草原和灌叢區。從時間序列的整體來看，2004年、2007年、2012—2014年(5 a)植被凈初級生產力整體要優于2005—2006，2008—2011年、2015年(7 a)，這與總NPP變化的偏差分析結果一致。在植被NPP較好的5 a里，荒漠面積明顯減少，草原和灌叢面積相應明顯增加，針葉林和闊葉林有小范圍面積的增加。

在N—S剖面時序變化圖中(圖3B)，蘇峪口國家森林公園植被凈初級生產力在水平地帶性分布上北部優于南部，其中N—S剖面北部NPP總量為0.042 5 Tg C，占比57.35%，南部NPP總量為0.031 6 Tg C，占比42.65%。受地形及日照影響，北半球植被NPP的高值部分主要分布在坡面北側(陰坡)[20]，因沿N—S剖面線存在多個山坡，每座山坡的北側植被NPP均高于南側，所以植被NPP的高值部分呈北高南低交替分布且NPP等值線北側梯度大，南側梯度小。2004—2015年N—S剖面線北部除2004年、2012年因降水豐富使得灌叢面積增加近20個像素(約10 km)外，其余年份絕大部分均由荒漠和草原覆蓋，且面積趨于穩定；2004—2015年N—S剖面線南部由荒漠、草原和灌叢交錯分布，2004年、2007年、2012—2014年針葉林邊際線向南部擴張了近40個像素(約20 km)，其余年份幾乎全部由荒漠、草原和灌叢覆蓋；2004—2015年N—S剖面線中部由針葉林和闊葉林交叉分布，且闊葉林林線在2004年和2012年所占范圍居多。

3.3.2 賀蘭山山脊線生態剖面垂直地帶性變化 因研究區呈地質學上的NNE—SSW向狹長分布，為全面表達研究區生態剖面上NPP變化細節，沿賀蘭山自然保護區地形山脊線建立生態剖面線，為研究生態剖面上NPP變化細節與DEM之間的對應關系，對研究區時間序列2004—2015年NPP剖面與DEM進行疊加分析(圖4)。

圖4 山脊線剖面與DEM疊加分析時序變化

由山脊線剖面與DEM疊加分析時序變化(圖4)可知，植被NPP垂直地帶性分布趨勢明顯，整體呈針葉林—闊葉林—灌叢—草原—荒漠的結構變化。由圖4的研究區生態剖面的NPP與DEM疊加分析可知，植被NPP變化情況與DEM關系密切，越接近山脊線中點(越在研究區中心)，DEM值越大，植被NPP值就越大。海拔在1 500 m左右幾乎無植被覆蓋，植被NPP值接近0值；海拔在1 500～2 000 m對應的植被NPP值的變化范圍為0～100 g C/(m2·a)；海拔在2 000～2 500 m對應的NPP值的變化范圍為100～200 g C/(m2·a)，且變化趨勢高度一致；海拔在于2 500 m時，植被NPP值基本維持在300 g C/(m2·a)左右。當海拔高于2 500 m以后，植被NPP值開始呈下降趨勢，植被NPP值在200～300 g C/(m2·a)變化。

3.4 MODIS NPP的氣候響應

3.4.1 氣候變化情況 2004—2015年賀蘭山區局地氣候呈濕冷化(圖5)。年均氣溫呈波動下降趨勢，年均溫度雖存在波動性但變化較為平穩，變化區間在1.5℃以內，最高值出現在2013年，為8.0℃，最低值出現在2011年，為6.6℃，氣溫變率為-0.01℃/a；年均降水量呈波動上升趨勢，研究區年均降水量變化區間在14 mm以內，在2012年出現的最高年均降水量為312.9 mm，在2005年出現的最低年均降水量為176.4 mm，降水變率為2.77 mm/a。

圖5 2014-2015年研究區氣候變化趨勢

3.4.2 氣候與NPP的響應分析 利用SPSS?Pearson相關性分析得到年均氣溫與總NPP間的相關性系數r=0.050 2，置信度p=0.876，則年均氣溫與總NPP不相關。通過圖6A可直觀看出氣溫與NPP的空間相關性，[0.45，0.95]高度正相關區間主要集中在蘇峪口國家森林公園附近區域，面積約為21.75 km2，占研究區面積的0.38%；[0.25，0.45)中度正相關區間主要零星分布在研究區中部，面積約為305.25 km2，占研究區面積的5.29%；93.13%的區域NPP與氣溫的相關系數均在[-0.25，0.25)低度相關區間范圍內；剩余0.12%的區域NPP與氣溫的相關系數均在[-0.45，-0.95)中度負相關區間范圍內，主要分布在于寧夏平原相交的邊界區域。[-0.95,-0.45)高度負相關區間所占面積為0。

圖6 氣候與NPP相關系數空間分布

利用SPSS?Pearson相關性分析得到兩者之年均降水與NPP間的相關性系數r=0.646 8，置信度p=0.023，則年均降水與NPP呈顯著正相關。通過圖6B可直觀表示年級降水與總NPP的空間相關性，[0.45，0.95]高度正相關區間面積約為2 976.50 km2，占研究區面積的51.59%，所占比例最大，主要分布在草原和灌叢等對降水量響應較強的區域；[0.25，0.45)中度正相關區間面積約為787.25 km2，占研究區面積的13.64%，主要分布在針葉林區域；[0.00，0.25)低度正相關區間面積約為1 961.50 km2，占研究區面積的34.00%，主要分布在海拔較高的闊葉林區域；[-0.25，0.00)低度負相關區間面積約為41.25 km2，占研究區面積的0.71%；[-0.45，-0.25)中度負相關區間面積約為3.46 km2，占研究區面積的0.06%；[-0.95，-0.45)高度負相關區間面積約為0。

4 討論與結論

4.1 討 論

全球變化對陸地生態系統影響巨大，預計將成為21世紀大量物種瀕于滅絕的主要因素[9]。近年來，國際上對氣候變化影響生態系統NPP方面開展了較多工作，研究表明生態系統NPP的正遭受氣候變化和人類活動的雙重影響[17-18]。

賀蘭山作為國家級自然保護區，濫砍盜伐得以有效禁止，且地處高寒、遠離城市，幾乎不受人類活動影響。因此，氣候變化是決定研究區植被生態系統的主要因素，直接或間接地影響著植被的凈初級生產力。且中國寧夏賀蘭山自然保護區位于季風氣候向大陸性氣候的過渡區，對氣候變化的響應極為敏感。研究表明，研究區植被NPP在12 a內總體呈增加趨勢，植被NPP與降水的相關性大于與氣溫的相關性。由于2004年9月—2005年11月在寧夏發生歷史性旱災[21]，大部分區域降水僅達到199.9 mm，為歷史上僅次于1982年(192.4 mm)的第二低值年，年降水日數僅為37 d，導致2005年總NPP成為2004—2015年研究期最小值；2006年寧夏地區出現的罕見暴雨天氣扭轉了歷史性干旱的局勢[22]，使得2006年總NPP未受到2004—2005年干旱年的影響，NPP總量得以恢復；2012年寧夏賀蘭氣象站點監測到日降水量創1951年有氣象記錄以來的最高值[23]，使得2012年NPP值顯著增大。說明在干旱與半干旱地區，降水變化是影響植被NPP的主導因素。這一結果與馬安青[24]、劉雪佳[25]、焦偉[26]等的研究結果一致。

研究還發現，在海拔3 000 m左右的區域，植被NPP與氣溫呈中高度正相關(r∈[0.25，0.45])的同時，對降水呈低度正相關(r∈[0.00，0.25])，說明海拔改變了氣候對研究區植被NPP的影響規律。但本文并未將高程、坡度、坡向等地形因素與氣候因素結合起來研究其與植被NPP的響應關系，這應是下一步的工作重點。

4.2 結 論

(1) 時間變化規律。賀蘭山2004—2015年平均NPP為97.91 g C/(m2·a)；總NPP變化范圍為0.324 4～0.501 5 Tg C/a；NPP增長速率為0.281 6 g C/(m2·a)。賀蘭山自然保護區是以山地針葉林和青海云杉為代表的荒漠區與半荒漠區典型森林生態系統，森林面積占比僅為10%～20%，森林生態狀況不容樂觀。

(2) 空間分布規律。以蘇峪口國家森林公園核心區中點做W—E，N—S剖面，得到長時序水平地帶分布規律：東部植被凈初級生產力要優于西部，北部植被凈初級生產力要優于南部；將賀蘭山沿山脊線生態剖面與DEM疊加分析得到垂直地帶分布規律：植被NPP值整體呈現出針葉林—闊葉林—灌叢—草原—荒漠的垂直景觀結構變化，NPP高值集中分布在海拔2 500 m以上的高海拔地區，在海拔1 500 m左右幾乎無植被覆蓋。

(3) 氣候響應規律。賀蘭山區局地氣候呈冷濕化趨勢，氣溫變率為-0.01℃/a，降水變率為2.77 mm/a。相較于氣溫而言，降水量對研究區NPP的影響較大。植被NPP與氣溫不相關(r=0.050 2，p=0.876)；植被NPP與呈降水量顯著正相關(r=0.646 8，p=0.023)，且高度正相關區域主要分布在草原和灌叢區。