伊犁地區植被生產力的時空分異及其與地形因子關系研究

張振東，昝 梅

(1.新疆師范大學地理科學與旅游學院，烏魯木齊 830054；2.山東師范大學地理與環境學院，濟南 250358)

陸地生態系統是碳元素重要的源和匯，在全球碳循環中扮演了至關重要的角色，而植被作為陸地生態系統的主體在保證全球陸地生態系統中的物質循環和能量流動的動態平衡中發揮了重要的作用[1].植被凈初級生產力(net primary productivity，NPP)是指綠色植物通過光合作用產生的全部有機物量扣除自身生長所必須的有機量后的剩余部分[2].由于NPP能夠以統一的標準反應生態系統生產力，可以很好的反映植被對自然資源的利用能力[3-4].因此，從19世紀90年代開始植被凈初級生產力的研究開始得到各國的廣泛關注，遙感技術也開始迅速發展.遙感技術在獲取地面植被動態變化方面具有時間分辨率高、精確度高的優勢，很快成為研究植被凈初級生產力的主流方式.CASA(carnegie ames stanford approach)模型也隨著遙感技術的進步而興起，CASA模型的優勢是利用遙感數據進行估算，參數少，穩定性強，能夠實現動態更新植被凈初級生產力[5].

伊犁地區是絲綢之路經濟帶的重要通道，地域優勢明顯.但隨著氣候變化和人類經濟活動影響的增加，綠洲冷島效應減弱、冰川面積退縮、草場退化和荒漠化等生態環境問題日趨嚴重，威脅著綠洲生態安全和生產的可持續發展[1，6].研究伊犁地區的植被凈初級生產力有助于更好的掌握伊犁地區植被和綠洲生態系統的變化情況，同時為以NPP為基礎的其他有關伊犁地區的科學研究提供理論基礎[3].目前已經有很多學者利用不同模型對中國各類生態系統的NPP進行反演.如朱文泉等采用改進的最小二乘法對中國典型植被的最大光利用率進行了系統的模擬，并針對不同植被分類精度可能帶來的誤差對最大光利用率進行了敏感性分析[7]；艾則孜提約麥爾·麥麥提利用MODIS數據采用GWR(地理加權回歸)建模法對2000年—2014年博斯騰湖流域NPP時空變化特征及影響因子進行了分析[8]；另外有多個地區和省份利用CASA模型求取植被的凈初級生產力，如東北三省、四川、河南、福建、和新疆等[1，4，9-10].但是，某些在社會、經濟和生態方面具有重要地位的熱點區域的植被生產力及影響因素的研究卻未被充分重視，如伊犁地區.因此，本文利用CASA模型，以伊犁地區為研究區，利用MODIS數據、日照數據和DEM數據運用CASA模型對伊犁地區NPP進行估算，進而對伊犁地區2012年—2014年植被NPP時空格局變化特征進行探究；分析伊犁地區不同植被類型NPP的特征；并闡明不同地形條件下伊犁地區NPP變化特征；最后通過對比驗證分析CASA模型在研究區的可用性和估算結果的可信程度[11].

1 材料與方法

1.1 研究區概況

伊犁哈薩克自治州地處新疆維吾爾自治區天山北部的伊犁河谷內(如圖1).位于80°09′42″～84°56′50″E，42°14′16″～44°53′30N″.地形三山環繞，地貌類型復雜多樣，主要劃分為山地、丘陵、平原和山間谷地四類；海拔高程為481～6 308 m .氣候屬于寒溫帶半干旱大陸性氣候；植被類型主要以草地、農作物、稀疏植被為主，分別占研究區總面積的60.02%、27.01%、7.37%，其他植被類型僅占5.62%.植被覆蓋度變幅大，受降水和氣溫的影響，山地區域植被覆蓋度呈上升趨勢，伊犁地區是新疆植被覆蓋度最高的地區，平均覆蓋度約90%.除雪線以上的常年冰雪帶和雪線以下的寒凍裸巖帶外，其余地面均有植被.低山帶為草地，中山帶為林地，高山帶為夏季優良草場[12].

圖1 研究區概況圖Fig.1 Overview of the study area

1.2 數據來源

1.2.1 MODIS的NDVI數據 本文選用的MODIS數據來源于美國國家航空航天局網站(http://ladsweb.nascom.nasa.gov/).下載了2012年—2014年的MOD13A1產品，共144幅圖像.空間分辨率為500 m，時間分辨率為16 d.采用最大值合成法(MVC)得到以月為周期的伊犁地區NDVI數據.

1.2.2 MODIS的土地覆被數據 本文選用MODIS土地覆蓋類型產品(land cover data)MCD12Q1，空間分辨率為500 m，時間分辨率為1 a.根據伊犁地區植被類型實際分布狀況，最終選擇植物功能型方案作為研究區土地覆被的分類方案將伊犁地區的植被類型分為針葉林(常綠針葉林、落葉針葉林)、闊葉林(落葉闊葉林)、草地、農作物(谷類作物、闊葉作物)、稀疏植被五種類型.通過以上方法我們分別獲得了研究區2012年、2013年和2014年的植被類型空間分布圖(如圖2).

1.2.3 氣象數據 從中國氣象數據網(http://data.cma.cn/)下載得到伊犁地區周圍11個氣象臺站2012年—2014年月平均氣溫、月降水量數據和日照時數的日值數據.通過克里格(Kriging)插值法生成500 m空間分辨率的月均氣溫和月降水量的序列柵格數據.將日值日照數據累加，計算得到每個站點逐月日照百分率后采用Kriging插值法獲得研究區逐月日照百分率的空間時序柵格數據，再利用太陽輻射能計算公式得到2012年—2014年逐月的太陽總輻射能序列柵格數據[13].

1.3 CASA模型的構建

本文選用的CASA模型是由Potter等基于Monteith的理論提出的.該模型是基于植被的機理過程建立的，利用APAR和ε來模擬植被的NPP[14-16]，所需參數相對較少且易于獲取.有效的避免了由于參數缺乏而人為簡化或估算產生的誤差[19].相對于其他模型的優勢是利用遙感數據進行估算，覆蓋的范圍廣；數據更新快；時間分辨率高；能夠實現動態更新植被凈初級生產力[5]，該模型的具體公式為：

NPP=APAR(x，t)×ε(x，t)，

(1)

式中，NPP表示植被凈初級生產力，APAR(x，t)表示植被吸收的光合有效輻射，ε(x，t)表示光能利用率.x表示像元的空間位置，t表示時間(月).

APAR主要由SOL和FPAR兩部分決定[5].APAR(x，t)的估算公式為(2)：

APAR(x，t)=SOL(x，t)×FPAR(x，t)×0.5，

(2)

式中，SOL(x，t)表示像元x在t月份的太陽總輻射能；FPAP(x，t)表示植被層對入射光合有效輻射的吸收比例；其中0.5是常數，它表示植被能夠吸收利用的太陽有效輻射量.

圖2 2012年—2014年伊犁地區植被類型空間分布圖(a為2012年、b為2013年、c為2014年)Fig.2 Spatial distribution map of vegetation types in Yili area from 2012 to 2014 (a is 2012，b is 2013，and c is 2014)

相對于地球輻射來說，SOL的波長較短，故又稱為短波輻射.SOL是利用日照時數采用經驗公式進行的估算.公式如下[13]：

(3)

式中a，b是計算太陽總輻射量的參數，其中，a=0.353，b=0.543[2].n/N為日照百分率，Ra為天文輻射.Ra由日地距離d，緯度φ，太陽日落角ω，太陽常數(0.082 MJ·m-2·min-1)，太陽赤緯σ和天數t通過公式(4)計算得到：

(ω·sinφ·sinσ+cosφ·cosσ·sinω)，

(4)

其中

ω=arccos(-tanφ·tanσ)，

根據上述公式可以計算出2012年—2014伊犁地區11個臺站逐月天文輻射.再利用ArcGIS軟件得到空間插值后的研究區天文輻射逐月柵格數據集.

根據文獻資料可以得到FPAR與植被類型和NDVI存在一定的關系.通過Kumar和Monteith等人的驗證，提出利用FPAR和NDVI建立的線性關系來估算FPAR，估算結果具有較高的精度[14，17]，具體公式為：

(FPARmax-FPARmin)+FPARmin，

(5)

式中，NDVI(x，t)是指像元x在t月份的NDVI值.NDVI(i，max)是指i種植被類型NDVI的最大取值.NDVI(i，min)是指i種植被類型NDVI的最小值.FPARmax和FPARmin是指FPAR的最大值和最小值，一般為常數，即FPARmax=0.950，FPARmin=0.001[13].其中，NDVI(i，max)和NDVI(i，min)與植被的類型有關，是把不同植被類型每月NDVI的值從小到大排序累積后選取5%和95%處的值分別作為該月的NDVI最小和最大值，詳見表1.

表1 各年份每種植被類型NDVI最大值最小值統計Tab.1 Statistics of maximum and minimum values of vegetation

光能利用率指植被層吸收入射光合有效輻射并將其轉化為有機物的能力，與溫度和水分等因子有關[15-16]．

ε(x，t)=T1(x，t)×T2(x，t)×W(x，t)×εmax，

(6)

式中，T1(x，t)、T2(x，t)是溫度脅迫因子，W(x，t)是水分脅迫因子.溫度和水分脅迫因子可以參考文獻研究結果得到[18-20].εmax是最大光能轉化率，即在理想狀態下植被吸收光合有效輻射并轉化為有機碳的最大光能轉化率[21]，一般與植被類型有關，本文主要參考朱文泉和馮益明等模擬的εmax結果求取平均值后獲得研究區主要植被類型的εmax[7，22-25]，詳見表2.

表2 各植被類型最大光能利用率Tab.2 Maximum light energy utilization rate of each vegetation type

本文利用MODIS數據(MOD13A1，MCD12Q1)、DEM數據和氣象數據，依照上述模型構建方法對數據進行處理，生成2012年—2014年月尺度NDVI數據、太陽輻射數據以及年尺度植被分類數據，然后進行CASA模型APAR、FPAR、SOL、ε等參數的估算，以此來研究伊犁地區植被生產力的時空分異及其與地形因子的關系，模擬流程如圖3所示.

圖3 植被凈初級生產力模擬流程Fig.3 Simulation process of vegetation net primary productivity

2 結果與分析

2.1 伊犁地區2012年—2014年NPP空間變化特征

模擬結果表明，研究區NPP的空間分布具有較大的差異.由圖4可以看出，2012年—2014年伊犁地區的植被NPP的年平均最大值為868.5 gC·m-2·a-1，年平均最小值為10 gC·m-2·a-1，總體平均值為200 gC·m-2·a-1.研究區2012年—2014年植被NPP的分布總體呈現出東部大于西部的特征.NPP的高值區(>500 gC·m-2·a-1)主要在伊犁東部及北部和南部的山地邊緣地區呈半環狀分布.中值區(200～500 gC·m-2·a-1)主要集中在中部河流流經的水資源較為豐富的地區.低值區(<200 gC·m-2·a-1)主要集中在伊犁中部植被覆蓋率低、城鎮較為集中的地區以及海拔較高的山區.

此外，結合2012年—2014年伊犁地區植被類型圖(圖2)可以看出NPP的高值區在針葉林、闊葉林和草地分布的地區較為集中，中值區在農作物分布的地區較為集中，低值區在稀疏植被分布的地區較為集中.以上分布特征也說明CASA模型估算得到的伊犁地區NPP的結果較為合理.

圖4 2012-2014年NPP空間分布(a為2012年、b為2013年、c為2014年)Fig.4 NPP spatial distribution from 2012 to 2014(a is 2012，b is 2013，and c is 2014)

2.2 伊犁地區2012年—2014年不同植被類型NPP月平均變化特征

由2012年—2014年不同植被類型逐月平均NPP變化趨勢(圖5)可以看出，研究區植被NPP隨著季節的變化呈現出周期性變化特征.不同植被類型的NPP隨時間變化的周期和大小不同，具有較好的區分度.

不同植被類型的逐月平均NPP呈單峰型曲線特征.各種植被類型的生長季均集中在每年的3～11月.各種植被類型的平均NPP基本都在2月達到最小值，7月達到最大值.另外，不同植被類型NPP的月變化趨勢基本相同，即3～7月NPP處于增長期，其中3～6月增長速率較快，6～7月增長速率較慢.各種植被類型NPP在7～12月為下降期，7～11月下降速率較快.

通過比較發現2012年—2014年不同植被類型逐月NPP平均值的變化幅度存在較大的差異.闊葉林月均NPP的變化速率和變換幅度最大.農作物、針葉林和草地月均NPP的變化速率和變化幅度差距不大，處于中等水平.稀疏植被變化速率和變化幅度最小.一般而言林地的固碳能力在同等水熱條件下大于草地和農作物，所以導致林地(闊葉林、針葉林)的月均NPP較其他植被類型大.而農作物和草地的NPP分別受到人為灌溉和物候、氣候環境的影響更大.由圖2研究區植被類型空間分布圖可知，稀疏植被主要分布在伊犁地區的南部高海拔山區，受水熱限制導致該種植被類型的月均NPP普遍較小.

圖5 2012年—2014不同植被類型逐月平均NPP變化趨勢Fig.5 Trends of monthly average NPP of different vegetation types from 2012 to 2014

2.3 不同植被類型對NPP的貢獻率

結合研究區植被類型空間分布圖(圖2)可知，伊犁地區針葉林和闊葉林主要分布在伊犁地區的南部山地邊緣，約占植被分布總面積的5.62%.草地面積占植被分布總面積的60.02%.農作物則主要集中在伊犁河和喀什河流經的水資源較為豐富的地區，占植被分布總面積的27.01%.草地和農作物所占面積超過植被分布總面積的87.02%.本文利用伊犁地區不同植被類型年平均NPP和不同植被類型分布面積的乘積在整個伊犁地區所有植被類型年平均NPP總和中所占的比例作為不同植被類型對該地區NPP的貢獻率.由表3可以看出，研究區草地NPP的貢獻率最大為58.59%；第二是農作物，對伊犁地區NPP的貢獻率為31.45%；針葉林和闊葉林對伊犁地區NPP的貢獻率基本相同，分別為3.75%、3.60%.最小的是稀疏植被，對伊犁地區NPP的貢獻率僅為2.16%.

表3 不同植被類型對NPP值的貢獻率Tab.3 Contribution rates of different vegetation types to NPP values

2.4 地形因子與伊犁地區植被NPP的關系

2.4.1 海拔高度與伊犁地區植被NPP的關系 利用伊犁地區DEM數據對伊犁地區NPP的空間分布進行分析.由圖6可以看出，伊犁地區年均NPP總體上呈現隨海拔增加而先增加后減少的特征，年均NPP的中值區分布在低海拔地區(<1 500 m)；高值區分布在中海拔地區(1 500～3 000 m);低值區分布在高海拔地區(>3 000 m).具體來說：一方面2012年—2014年均NPP隨海拔變化的梯度性明顯，年均NPP為0～100 gC·m-2·a-1主要集中在3 500～5 500 m；年均NPP為100～200 gC·m-2·a-1主要集中在海拔504～1 000 m和3 500～4 000 m；年均NPP為200～400 gC·m-2·a-1主要集中在海拔1 000～1 500 m和3 000～3 500 m；年均NPP為500～600 gC·m-2·a-1主要集中在海拔1 500～3 000 m.另一方面：年均NPP為400～500 gC·m-2·a-1在海拔504～3 500 m范圍內分布較均勻，在不同海拔梯度內所占的比例差距較小，分別為27.15%(500～1 000 m)、20.74%(1 000～1 500 m)、28.31%(1 500～2 000 m)、23.77%(2 500～3 000 m)、20.93%(3 000～3 500 m)；年平均NPP為0～100 gC·m-2·a-1在海拔高于3 500 m時占不同海拔梯度內的比例均為最高，分別為52.18%(3 500～4 000 m)、93.43%(4 000～4 500 m)、99.47%(4 500～5 000 m)、100%(5 000～5 500 m)，這與研究區年均NPP隨海拔的變化情況基本一致.另外，由圖6可以看出，多個年均NPP梯度(500～600 gC·m-2·a-1、600～700 gC·m-2·a-1、700～800 gC·m-2·a-1)在海拔2 000～2 500 m處均呈現出增長的趨勢.出現這種現象的原因主要是由于海拔較低的地區人口稠密，天然植被生長受人類影響強烈，破壞嚴重，但農作物大部分生長于此，光熱充足并且膜下滴管提供充足水分，年平均NPP處于中值水平.隨著海拔的升高，在中山帶形成降雨帶，天然植被生長茂盛，年平均NPP最高.但隨著海拔的持續升高降水減少，溫度降低，土層稀薄，土壤肥力差，年平均NPP降到最低，基本保持在0～100 gC·m-2·a-1范圍內.

圖6 2012年—2014年平均NPP隨高程的變化Fig.6 Average NPP changes with elevation in 2012-2014

2.4.2 坡度與伊犁地區植被NPP的關系 圖7為2012年—2014年平均NPP隨坡度變化圖，可以看出年平均NPP呈現隨坡度增加而先增加后減小的特征.當坡度小于35°時，研究區植被年均NPP 400～500 gC·m-2·a-1主要集中在坡度0～5°；年均NPP為600～700 gC·m-2·a-1的植被主要集中在5～25°坡度范圍內；年均NPP為0～300 gC·m-2·a-1的植被則主要分布在坡度大于35°范圍內，另外，年均NPP 500～600 gC·m-2·a-1在0～35°坡度范圍內分布較均勻，在不同坡度梯度內所占的比例差距較小，分別為29.84%(0～5°)、29.82%(5～15°)、28.53%(15～25°)、24.27%(25～35°).由于坡度較低的地區地勢平坦，水熱條件好，適宜農作物的生長，導致植被年均NPP保持在中值水平；隨著坡度的增加，人類干擾減少，林地面積增大，年均NPP增大；當坡度大于35°后，土壤侵蝕加劇，土層薄，土壤肥力差，不利于植被的生長，導致植被年均NPP降低.

圖7 2012年—2014年平均NPP隨坡度的變化Fig.7 Average NPP changes with slope in 2012-2014

3 討論

3.1 結果驗證

通過查閱文獻可知MODIS的NPP產品已經被許多學者驗證過結果較好，而且部分地區結果精度較高[16].用MOD17A3 NPP數據對研究區利用CASA模型的估算的植被NPP結果進行對比驗證，結果如圖8所示.選用R2作為趨勢線擬合程度的指標，一般在0～1范圍內取值，值越大擬合程度越高.通過對2012年—2014年CASA模型NPP估算值和MOD17A3模擬數據進行比較發現2012年—2014年的R2值分別為0.87、0.86、0.95.說明利用CASA模型估算伊犁地區的NPP與MOD17A3的結果一致性較好，也說明利用CASA模型對伊犁地區植被NPP估算實用性較強.

圖8 2012年—2014年CASA模型估算值與MODIS產品NPP值比較Fig.8 Comparison of 2012-2014 CASA model estimates with NPP values of MODIS products

3.2 植被NPP的時空格局

植被NPP的空間分布受地形、土壤、水分和植被類型等因素的影響，具有很強的空間異質性.通過CASA模型模擬結果表明，伊犁地區的植被NPP分布總體呈現出東部大于西部的特征.綜合伊犁地區地形及植被類型分布情況分析發現(圖2)，伊犁東部多為農作物、林地、草地，該區域位于伊犁河下游，區域水熱條件好，草地質量好，面積廣闊，并且地形較平坦，利于人為灌溉農田系統的植被生長，加上伊犁北部和南部山地的邊緣地區林地面積廣闊，同等水熱條件下林地的固碳能力要大于其他植被類型；從而導致上述地區的植被NPP較大；中部地區植被類型主要為農作物，植被NPP受水分條件影響大，該區域河流廣布，灌溉水源豐富，水土條件優良，因而植被NPP處于中值水平；而南部地區受地形、土壤類型的限制，以山地和荒漠為主，植被稀疏，導致該區域植被NPP較低；伊犁西部地勢相對較低，水熱條件相對不足，導致該區域植被類型單一的草地NPP明顯很低.研究區不同植被類型的年內月均NPP呈單峰型曲線主要是受該區域植被物候周期影響所致.

3.3 植被NPP與地形因子的關系

海拔和坡度是主要的地形指標，伊犁地區年均NPP呈現隨海拔的增加而先增加后減少的特征.年均NPP的中值區分布在低海拔地區(<1 500 m)；高值區分布在中海拔地區(1 500～3 000 m);低值區分布在高海拔地區(>3 000 m).研究區受山盆結構分異造成的水熱條件限制明顯，植被年均NPP也呈現出隨坡度增加先增加后減小的特征.植被分布集中在高程小于3 000 m且坡度在0～35°的區域內，此區域水熱條件也適合植被的生長，因此在此區域植被NPP年均值最高.而在海拔大于3 000 m并且坡度大于35°的區域內土壤侵蝕加劇，土層薄，土壤肥力差，加之水熱條件不足，不利于植被的生長，因此植被年均NPP在此區域出現最低值.

4 結論

本文利用美國國家航空航天局網站上以月為周期的NDVI數據集MODIS13A1、以年為周期的植被類型數據集MCD12Q幅以及各種氣象數據，利用CASA模型估算得到了以月為周期，空間分辨率為500 m的伊犁地區植被NPP.進而分析了研究區2012年—2014年植被NPP的時空分布特征并探討了不同植被類型對NPP的貢獻率.主要結論如下.

1) 不同植被類型的月均NPP呈現3～7月增加，7月到11月下降的趨勢.2月各種植被類型的NPP達到最小值，7月達到最大值；但不同植被類型月均NPP的變化速率和變換幅度存在較大差異.闊葉林月均NPP的變化速率和變換幅度均為最大.農作物、針葉林及草地月均NPP的變化速率和變化幅度差距不大，在所有植被類型中變化速率和變化幅度居中.稀疏植被月均NPP變化速率和變化幅度均最小.

2) 不同植被類型對伊犁地區NPP的貢獻率不同.其中草地的貢獻率最大為58.59%；其次是農作物，貢獻率都為31.45%；針葉林和闊葉林對伊犁地區NPP的貢獻率基本相同，分別為3.75%、3.60%.最小的是稀疏植被，對伊犁地區NPP的貢獻率僅為2.16%.

3) 伊犁地區年均NPP隨海拔和坡度的增加均呈現出先增加后減少的特征；通過驗證說明，利用CASA模型對伊犁地區植被NPP估算的結果與MOD17A3產品一致性較好.

氣候尤其是水熱因子對植被凈初級生產力的影響極為顯著，植被本身的物理特性對植被NPP的影響同樣不可忽略，但本研究僅僅針對地形因子與植被凈初級生產力的關系進行研究，具有一定的片面性；另外，本研究缺乏實地植被NPP的采樣驗證.這些問題都有待于進一步深入探討.