基于CASA模型長江流域植被NPP時空演變及與地形因子的關系

李俊豪,梁娟珠

(福州大學 數字中國研究院(福建),福建 福州 350002)

全球氣候變暖、溫室氣體排放和極地冰川融化等環境問題已變成全人類共同面臨的不同尋常的挑戰[1]。我國提出力爭2030年前實現碳達峰,2060年前努力實現碳中和的目標,由相對減排到絕對減排,進而實現零排放[2]。植被凈初級生產力(net primary productivity, NPP)是指在單位面積上植物所積攢的有機物的數量,是由光合作用所固定的有機物總量扣除自身呼吸作用消耗量后的剩余部分,同時,植被NPP是定量描述植被生產能力特征的指標,是衡量區域碳匯能力高低和調節生態過程的重要指標[3],也是陸地生態系統中物質與能量轉換研究的基石?？偠灾?植被NPP就是判定植被生長狀況和固定碳及積累能量的能力。

植被NPP的反演方法主要包括站點實測法和模型法[4-5]。實測法主要在小區域上進行,選取站點的生物量推算NPP,雖然精度高但是費時費力,難以在大范圍區域開展實時動態監測。因此,在全球或區域尺度,模型法成為一種重要的主流的研究方法,而站點實測數據更多用于NPP反演結果的精度驗證。目前,反演植被NPP發展了一系列的模型反演方法,主要包括氣候模型、過程模型和光能利用率模型[6]。

氣候模型主要用于植被NPP早期研究,該類模型相對簡單且所需要的氣象數據易獲取,反映的是理想情況下的植被NPP,即一種區域潛在的植被NPP[7]。而過程模型[8]綜合考慮了植物光合作用、呼吸作用、有機物的分解以及合成作用等多種因素,但是模型本身比較復雜,且輸入的數據較難獲得,區域尺度轉換難以實現。CASA(Carnegie-Ames-Stanford approach)模型是一個基于光能利用率的具有代表性的模型,在全球以及區域生產力的估算中有著較大的影響力[9]。國內外學者利用CASA模型對草地、植被等不同研究對象在全國、省、市群等不同空間尺度開展了大量研究。LIU等[10]通過CASA模型反演,結果表明,1982—2008年全球草地NPP總量為23.98～25.73 PgC,年均增速為0.03 PgC/a。PIAO等[11]利用CASA模型反演了我國1982—1999年植被NPP,結果表明,植被NPP整體上呈增加趨勢,增加速率為0.015 PgC/a。周偉等[12]利用CASA模型模擬我國1982—2010年草地NPP總量為988.3 TgC,其中超過60%的草地表現為增加的趨勢。李軍玲[13]、戴靚[14]等采用CASA模型并借助其他輔助數據分別對河南省和長江三角洲地區的NPP進行分析。賀倩等[15]基于改進的CASA模型反演得到2010—2015年三江源地區植被NPP均值降低51.99 gC/(m2·a)。YANG等[16]基于CASA模型模擬2001—2016年安徽省植被NPP,均值約為228.27 gC/(m2·a),不同土地利用類型的植被NPP均值從高到低為耕地、草地、林地、農用地、水體、濕地、低植被覆蓋區。

在全球植被NPP模擬和碳循環研究中,CASA模型在我國和區域植被生產力的模擬中發揮了很大作用[17];但是由于需要的參數多,精確賦值困難,該模型仍缺乏廣泛和細致的應用,尤其是在流域尺度上,比如我國長江流域地區。因此,本研究選用2001—2018年的MODIS數據和中國氣象數據驅動CASA模型反演長江流域植被NPP,探究長江流域不同時空尺度的植被NPP分布現狀和變化趨勢,以期為流域植被恢復及流域生態環境保護提供科學參考。

1 研究區概況、數據來源

1.1 研究區概況

長江流域是世界三大流域之一,總面積約為180萬平方公里,約占中國陸地總面積的1/5,流域內植被類型豐富[17],主要為常綠闊葉林、常綠針葉林、落葉針葉林、混交林和灌叢等。流域地形地貌類型復雜多樣,地勢西高東低,橫跨我國地勢的三大階梯,主要流經青藏高原、云貴高原、四川盆地、長江中下游平原。流域內氣候類型多樣,其中大多數地區處在亞熱帶季風氣候區。流域內氣溫空間呈東高西低、南高北低的分布趨勢,多年平均氣溫為13 ℃,平均年降水量為1 100 mm,雨量充沛但空間分布不均,呈東南向西北遞減的趨勢。如圖1所示,2001—2018年歸一化植被指數(normalized differential vegetation index,NDVI)多年均值范圍為0～0.95。

圖1 長江流域氣象站點分布及二級流域示意圖Fig.1 Schematic diagram of meteorological station distribution and secondary river basin in the Yangtze River Basin

1.2 數據來源及預處理

遙感數據:NDVI數據通過NASA官網發布的MOD13A1數據集(https://ladsweb.modaps.eosdis.nasa.gov/missions-and-measurements/products/),空間分辨率500 m,其行列號分別為h25v05、h26v05、h26v06、h27v05、h27v06、h28v05、h28v06,共7景影像,并對其進行格式轉換、投影變換、鑲嵌、拼接,通過最大值合成法(maximum value composite, MVC)進一步抵消環境因素的干擾,從而獲得逐月NDVI數據。2005年、2010年、2015年、2018年植被覆蓋數據MCD12Q1,空間分辨率為500 m,提取國際地圈-生物圈計劃(International Geosphere-Biosphere Programme,IGBP)全球植被分類圖。 MOD17A3HGF數據是NPP的數據產品,來源于NASA官網。

氣象數據:降水和氣溫數據來源于國家氣象科學數據中心(http://data.cma.cn/),數據包括長江流域2001—2018年175個氣象站點的地面觀測數據(平均氣溫、降水量等),通過克里金插值法獲取月降水和月平均氣溫數據,并使其空間分辨率與NDVI數據一致。太陽輻射使用全球高分辨率地表太陽輻射數據集,其來源于國家青藏高原科學數據中心(http://data.tpdc.ac.cn),時間單位和空間分辨率分別為月和10 000 m。

地形數據:高程、坡度、坡向數據分別來源于地理空間數據云的SRTMDEM 90M、SRTMSLOPE 90M、SRTMASPECT 90M的數據集。

2 研究方法

2.1 NPP的反演

CASA模型充分考慮了環境條件及植物本身特性,將氣象及環境調控因子作為輸入參數,有效地反演陸地植被NPP,且模型參數可通過遙感或空間分析方法獲得。因此,本文采用CASA模型反演長江流域植被NPP,計算公式[18]如下:

Nx,t=Ax,t×εx,t

(1)

式中:Nx,t為像元x處植被在t月的凈初級生產力;Ax,t為像元x處植被在t月內所吸收的光合有效輻射;εx,t為像元x處植被在t月內的光能轉化率。

Ax,t=Mx,t×Fx,t×0.5

(2)

式中:Mx,t為像元x處植被在t月的太陽總輻射量,MJ/m2;Fx,t為像元x處植被在t月的光合有效輻射吸收比;常數0.5為植被所能利用的太陽有效輻射(波長范圍為0.4～0.7 μm)占太陽總輻射的大致比例。

光能轉化率的計算公式如下:

εx,t=T1,x,t×T2,x,t×Wx,t×εmax

(3)

式中:T1,x,t和T2,x,t為像元x處植被在t月內的溫度因子;Wx,t為像元x處植被在t月內的降水因子;εmax為最大光能轉化率。

由于光合有效輻射吸收比(F)與歸一化植被指數(INDV)和比值植被指數(IRV)分別存在線性關系,本文根據朱文泉[19-20]的研究,計算出適合長江流域的每個植被類型所對應的INDV和IRV的最大值,賦予最小值,使用2個參數分別進行計算后,取其平均值作為F的值。結合文獻[21],設置最大光能轉化率(εmax),并從數據獲取層提升模型運用的區域適用度,從而實現對靜態參數的本地化(表1)。

表1 長江流域植被類型及其參數設置Tab.1 Vegetation types and parameter setting in the Yangtze River Basin

2.2 植被NPP趨勢分析

Theil-Sen Median趨勢分析法可用于逐像元分析植被NPP在時間序列上的變化趨勢及空間差異性,該方法無需數據服從特定分布且不受異常值的干擾,對離散數據具有較強的規避能力[22]。其計算公式為

(4)

式中:C為植被NPP變化速率,C>0說明植被NPP呈上升趨勢,C<0 則植被NPP為下降趨勢;Ni為i年植被NPP值;Nj為j年植被NPP值。

由于Theil-Sen Median趨勢分析法無法進行顯著性檢驗,與Mann-Kendall非參數檢驗方法結合,能較好地克服噪聲對分析結果的影響[23]。因此,植被NPP變化趨勢的顯著性通過Mann-Kendall統計檢驗法判斷,相關計算公式如下:

檢驗值:

(5)

符號函數:

(6)

(7)

標準化統計量:

(8)

式中:sgn是符號函數。趨勢檢驗的方法是零假設,H0:R=0;當|Z|>Z1-α/2時,拒絕零假設。其中,Z1-α/2為標準正態方差;α為顯著性檢驗水平。當|Z|分別大于等于1.28、1.96時,表示分別通過了置信度為90%、95%的顯著性檢驗。

2.3 地形因子的分析

地形因子不僅是可以描述地貌形態的參數,也是植被NPP垂直分布及多樣性分布的關鍵因素。根據長江流域復雜的地形地貌特征,選取的地形因子包括海拔、坡度和坡向。借鑒《水土保持綜合治理規劃通則》[24],結合長江流域地形情況,使用Arcgis 10.2軟件對海拔、坡度和坡向重分類并進行等級劃分。地形因子分類定級見表2。將研究區各級地形因子與NPP進行疊加分析,統計各級地形因子對應的NPP均值,分析NPP隨地形因子等級變化的特征。本文定量分析NPP與高程、坡度和坡向之間的關系,在各地形因子等級選取若干采樣點,通過斯皮爾曼等級相關系數研究在各地形因子等級中變量之間的相關關系。

表2 地形因子分類定級Tab.2 Terrain factor classification and grading

3 結果與分析

3.1 植被NPP反演結果驗證

3.1.1反演結果與遙感數據產品對比

由于長江流域植被NPP實測數據較難獲取,基于2001—2018年長江流域的MODIS NPP數據產品(MOD17A3HGF),經過處理得到長江流域2001—2018年多年NPP的平均值,在研究區內隨機選取675個樣點,將CASA模型反演的數據與MOD17A3HGF產品進行擬合(圖2),以驗證本文研究結果。結果表明,兩者之間的擬合程度好,本文反演值要高于MOD17A3HGF產品。其原因是:MOD17A3HGF產品低估了灌叢、草原草甸的生產力,與其他植被類型的擬合程度較好[25]。MOD17A3HGF產品在長江流域的城鎮及周邊地區,如長三角地區、鄱陽湖水系北部、洞庭湖水系北部、宜昌至湖口中部、漢江東部等出現缺值,突出本文CASA模型用于長江流域研究的優勢,并且二者相關性R2為0.717(P<0.01),具有強相關性,說明該模型可以有效地反演長江流域植被NPP,可為之后的研究提供基礎數據。

圖2 反演值與MOD17A3HGF產品的比較Fig.2 Comparison of the inversion values with the MOD17A3HGF product

3.1.2對比驗證

為了進一步驗證長江流域植被NPP數據的可靠性,將本文的反演值與其他模型的數據進行比較。由表3可知:本文不同植被NPP的反演值在實測值范圍內,并且與前人的反演結果一致,均表現為喬木類NPP較高、稀疏植被NPP較低,說明參數本地化后的CASA模型反演的NPP數據基本可靠,可以很好地反映長江流域植被NPP的實際狀態。

表3 植被NPP反演值與實測值的比較 單位:gC/(m2·a)

3.2 長江流域植被NPP時間變化特征

3.2.1植被NPP年際變化特征

由CASA模型反演的NPP數據繪制2001—2018年長江流域植被NPP年際變化圖,如圖3所示。從圖3可以看出:長江流域在這18年間表現出增長的趨勢,但是增幅較小,NPP年均值范圍為511.92～638.25 gC/(m2·a),NPP多年均值為572.72 gC/(m2·a);峰值638.25 gC/(m2·a)出現在2007年,比NPP多年均值高了11.44%;谷值511.92 gC/(m2·a)出現在2001年,比NPP多年均值低了10.62%。由于受到氣候干旱[29]的影響, 2001年NPP均值明顯小于其他年份?？傮w而言,2001—2007年NPP年均值變化明顯,呈波動式上升,2008—2018年變化相對趨于平穩,表明長江流域生態環境逐漸恢復,區域固碳能力有所增強。

圖3 2001—2018年長江流域植被NPP年際變化圖Fig.3 Interannual change of NPP in the Yangtze River Basin from 2001 to 2018

3.2.2植被NPP年內變化特征

對2001—2018年長江流域植被NPP取月平均值合成,并計算月NPP總量占年NPP總量的比例(貢獻度),如圖4所示。由圖4可知:長江流域NPP年內變化層次感比較明顯,表現出倒“U”型變化特征和季節性。夏半年(4—9月)NPP月均值處在較高位,均大于40 gC/(m2·a),且夏半年NPP總量占全年的77.61%;冬半年(10月—次年3月)NPP月均值整體較低,其中11月—次年2月NPP月均值小于20 gC/(m2·a)。夏季6—8月NPP月均值都超過80 gC/(m2·a),NPP總量占全年NPP的44.66%;春季3—5月、秋季9—11月NPP月均值較低,分別占NPP總量的27.22%、20.35%。冬季對全年NPP的貢獻度僅為7.77%,不足夏季的50%;7、8月NPP月均值最高并且貢獻度最大,均值分別為90.85 gC/(m2·a)和83.92 gC/(m2·a),貢獻度分別是15.86 %和14.65%;12、1月NPP月均值最低,均值分別為13.81、13.43 gC/(m2·月),貢獻度分別為2.35 %和2.41 %。植被NPP表現出的年內變化特征,原因在于:長江流域的植被覆蓋度高的地區位于濕潤地區,屬于亞熱帶季風氣候,夏季高溫多雨,冬季寒冷干燥,從而引起植被在夏半年生長迅速,在冬半年生長緩慢。

圖4 長江流域植被NPP月均值變化圖Fig.4 Monthly mean change chart of NPP of vegetation in the Yangtze River Basin

3.3 長江流域植被NPP空間分布及變化特征

3.3.1植被NPP空間分布格局

對2001—2018年長江流域植被NPP取平均值,繪制空間分布圖。如圖5所示:長江流域多年植被NPP均值在0～1 700.13 gC/(m2·a)范圍內,流域內植被NPP空間分異明顯,整體分布格局上呈階梯狀分布,由第一階梯向第三階梯遞增,自西北向東南遞增,植被NPP高值聚集在青藏高原山麓地區和鄱陽湖水系地區,而青藏高原部分地區和長三角地區植被NPP較低。

圖5 長江流域植被NPP的空間分布圖Fig.5 Spatial distribution diagram of the vegetation NPP in the Yangtze River Basin

為了更好地了解長江流域空間分布格局,分區統計各二級流域2001—2018年植被NPP。由圖6可以看出:12個二級流域中,烏江、嘉陵江、漢江、洞庭湖水系、湖口干流、鄱陽湖水系共6個水系植被NPP均值超過600 gC/(m2·a),其中鄱陽湖水系植被NPP均值679.99 gC/(m2·a)最大,得益于鄱陽湖所在的江西省高度的植被覆蓋度;而金沙江石鼓以上水系的NPP年均值僅為274.00 gC/(m2·a),遠遠低于其他流域,其他5個二級流域(金沙江石鼓以下、岷沱江、宜賓至宜昌、宜昌至湖口、太湖水系)植被NPP均值為500～600 gC/(m2·a)。

圖6 長江流域二級流域的植被NPP均值統計圖Fig.6 Statistical graph of NPP mean value in the secondary level basin of the Yangtze River Basin

3.3.2植被NPP空間變化特征

本文采用Theil-Sen Median趨勢分析法和Mann-Kendall統計檢驗法將長江流域植被變化趨勢劃分為顯著減少(C<0和|Z|≥1.96)、輕微減少(C<0和1.28≤|Z|<1.96)、基本不變(|Z|<1.28)、輕微增加(C>0和1.28≤|Z|<1.96)和顯著增加(C>0和|Z|≥1.96)總共5個分類,繪制2001—2018年的變化速率以及變化趨勢圖。從圖7可知:長江流域植被NPP平均變化速率為0.066 gC/(m2·a),有17.76%的面積表現為增加趨勢,主要分布在岷沱江南部、嘉陵江西部、烏江西部和洞庭湖東部地區,其中岷沱江南部、嘉陵江西部、烏江西部增幅超過5 gC/(m2·a),說明這些區域植被生長環境較好,為植被生長提供豐富的養料。長江流域植被NPP減少趨勢的像元占全域的16.16%,顯著減少的像元占比僅5.48%,主要分布于漢江西部與長江流域西南部,這與當地近年來氣溫下降[30]有關,氣溫較低會導致植被生長放緩。長江流域植被NPP基本不變的像元占比66.08%,主要分布在長江中部和西北部。由此看出,長江流域整體植被生長狀況和固碳能力在改善,但是局部地區有惡化的趨勢。

圖7 2001—2018年長江流域植被NPP變化趨勢圖Fig.7 Trend chart of NPP change of vegetation in the Yangtze River Basin from 2001 to 2018

3.4 長江流域植被NPP變化與地形因子的相關性分析

3.4.1植被NPP與高程的關系

為了解植被NPP與高程之間的關系,統計分析2001—2018年各高程帶植被NPP年均值變化特征,并在各高程帶中選取若干個采樣點分析植被NPP與高程的相關性,如圖8所示。由圖8(a)可以看出:2007年,>500～1 000 m和>1 000～2 000 m高程帶的植被NPP均超過>2 000～3 000 m高程帶的植被NPP,原因是該年降水量激增,在政府有效的調度下[31],干旱缺水的地區有了水分補給,促進了這些地區植被的生長。其他時間,長江流域各高程帶植被NPP變化趨勢與長江流域植被NPP年際變化趨勢整體上一致,都呈現波動性。高程在3 000 m以下時,植被NPP年均值隨著各級高程帶的增加而增加;高程在3 000m以上時,隨著高程的增加植被NPP年均值減少。因此,2007年植被生長最佳的是>1 000～2 000 m高程帶,其他時間是>2 000～3 000 m高程帶植被NPP年均值最高。由圖8(b)可以看出:<200、200～500、>500～1 000、>1 000～2 000 m的高程帶中,植被NPP與高程呈正相關,其中>1 000～2 000 m高程帶呈顯著正相關;>2 000～3 000 m高程帶呈顯著不相關;高程>3 000～4 000 m和>4 000 m的區域呈顯著負相關;整個研究區植被NPP與高程的相關系數是-0.487(P<0.01)。說明高程在>1 000～2 000 m地區,植被NPP隨高程的增加顯著增加;高程在>2 000～3 000 m地區植被NPP變化比較穩定,而高程在>3 000 m地區植被NPP隨高程的增加急劇減少。綜上,高程對NPP的解釋較好,其中>2 000～3 000 m高程帶上植被生長狀況最佳。

3.4.2植被NPP與坡度的關系

2001—2018年長江流域地區植被NPP與坡度的關系如圖9所示。從圖9(a)可以看出:2001—2018年長江流域植被NPP隨坡度的增加而增加,且各坡度范圍植被NPP變化趨勢與研究區NPP年際變化趨勢保持一致。圖9(b)中各坡度分級范圍與植被NPP均呈正相關性。<6°坡度的區域易受到人類活動的影響,因此該地區的坡度對NPP的解釋較差。在6°～15°和>15°～25°坡度范圍,NPP與坡度呈顯著正相關性,即植被NPP隨坡度的增加而增加,這些地區是開展退耕還林、小流域治理、生態移民搬遷等一系列生態工程的重點區域,并且人類活動較少,因而該地區植被生長旺盛,NPP增長。坡度>25°地區多為陡峭的山崖,水分難以保持,氣候條件很難控制,植被面臨光照不足、溫度變化等挑戰;相較坡度<6°地區,該地區的坡度對NPP的解釋略好。因此,<6°坡度地區的NPP是各坡度區間最低的,并且該地區坡度對NPP的解釋較差,在6°～15°和>15°～25°坡度對植被NPP的解釋較好。

圖9 2001—2018年長江流域植被NPP與坡度的關系Fig.9 Relationship between NPP and slope in the Yangtze River Basin from 2001 to 2018

3.4.3植被NPP與坡向的關系

2001—2018年長江流域地區植被NPP與坡向的關系如圖10所示。由圖10(a)可知:平坡、陰坡、半陰坡、陽坡、半陽坡植被NPP變化趨勢與NPP年際變化趨勢一致,但是各坡向NPP年均值相差不大。從圖10(b)可以看出:只有半陰坡與NPP多年均值呈現顯著正相關,其他坡向與NPP相關性系數較低,均小于0.2。因此,除半陰坡外,其他坡向對NPP影響不大。

圖10 2001—2018年長江流域植被NPP與坡向的關系Fig.10 Relationship between NPP and slope direction in the Yangtze River Basin from 2001 to 2018

4 結論與討論

本文基于CASA模型,針對NDVI和最大光能利用率εmax進行本地化,反演出以月為時間尺度的長江流域植被NPP,并分析2001—2018年研究區植被NPP時空變化規律及與地形因子的關系,得出以下結論:

1)通過與MODIS NPP數據產品以及前人的相關研究進行對比,CASA模型可以有效地反演長江流域植被NPP,反映長江流域植被生長的實際狀態。長江流域空間格局顯示東高西低,2001—2007年NPP均值呈波動式上升,2008—2018年變化相對趨于平穩;長江流域植被NPP具有年內季節性的表現,夏季NPP對全年NPP總量的貢獻接近一半。

2)長江流域在2001—2018年18年中植被NPP增加的面積大于減少的面積,增加的區域主要分布在岷沱江水系南部、嘉陵江水系西部、烏江水系西部和洞庭湖水系東部地區,因此長江流域植被生長狀況整體來看有所改善,然而局部地區植被的生長存在惡化的現象。

3)關于地形因素的分析發現,相較于坡度、坡向,高程對長江流域植被NPP的解釋較強。長江流域植被NPP與不同高程帶的響應不同,在<3 000 m的地區植被NPP隨著高程的增加而增加,在高程>3 000 m的地區植被NPP顯著減少,其中最適合植被生長的高程帶是>2 000～3 000 m。在6°～15°和>15°～25°坡度以及半陰坡對植被NPP的解釋較好。

本文反演的NPP數據集時間序列長,能夠較好地反應植被在不同時空尺度下的變化。地形因子反應是植被在自然環境條件下的綜合表現,高程、坡度和坡向等地形因子的改變,影響降水量、氣溫高低、土壤肥力等資源在時空中的分配,進而影響植被生產能力。分析地形因子與植被NPP之間的關系能夠有效地揭露出不同生長條件下植被的分布。除了地形因子的影響,氣溫、太陽輻射、降水等因子均會直接影響植被NPP的時空變化,各因子之間相互限制和相互關聯。如何深入分析植被NPP的各因素之間的綜合影響,是研究植物生態環境和評估區域碳收支的基石,今后將進一步研究各自然因子對植被NPP的綜合作用。