陜北長城沿線植被恢復(fù)與生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)模擬分析

米朝娟,周自翔,劉 婷,武 佳,孫彥旭

西安科技大學(xué)測繪科學(xué)與技術(shù)學(xué)院,西安 710054

生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)是干旱和半干旱地區(qū)生態(tài)系統(tǒng)最重要的生態(tài)服務(wù)功能之一,實施防風(fēng)固沙工程是抑制和固定風(fēng)沙區(qū)沙塵不被再次移動以實現(xiàn)荒漠化逆轉(zhuǎn),有利于區(qū)域經(jīng)濟和人類福祉的可持續(xù)發(fā)展[1]。生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)通常從風(fēng)蝕角度為防風(fēng)固沙服務(wù)提供研究。學(xué)術(shù)界對風(fēng)蝕的評估始于1940年代,風(fēng)蝕模型是評估土壤風(fēng)蝕狀況的主要技術(shù)手段,包括德克薩斯模型(Texas erosion analysis model,TEAM)[2—3]、風(fēng)蝕預(yù)報系統(tǒng)(WEPS)[4—5]、風(fēng)蝕方程(WEQ)[6]、修正風(fēng)蝕方程模型(Revised Wind Erosion Equation,RWEQ)[7—9]等。由于一些模型主要集中在風(fēng)蝕發(fā)生的可能性,無法評估潛在風(fēng)蝕量。而大多數(shù)學(xué)者考慮到REWQ模型具有良好的適用性,涉及氣候條件、植被覆蓋、土壤易蝕性、土壤結(jié)皮、地表粗糙度等綜合因素來計算植被造成的沙粒滯留量所發(fā)生的風(fēng)力侵蝕[10—12]。部分學(xué)者不斷探索和驗證該模型并通過調(diào)整參數(shù),利用RWEQ模型對寧夏[7—12]、內(nèi)蒙古[13—14]、京津冀[15—17]等中國干旱半干旱風(fēng)沙區(qū)的風(fēng)蝕狀況進行了評估,取得了良好的效果。

植被是影響防風(fēng)固沙生態(tài)功能的關(guān)鍵指標(biāo),也是檢驗防風(fēng)固沙區(qū)生態(tài)保護成效的重要依據(jù)[11]。恢復(fù)植被有利于保持土壤水分、改變地表狀態(tài),從而降低風(fēng)蝕發(fā)生的潛力[12]。風(fēng)蝕的動態(tài)變化是自然氣候變化和人類活動共同作用的結(jié)果,退耕還林使地表覆蓋發(fā)生顯著變化,進而與生態(tài)系統(tǒng)服務(wù)能力的變化相關(guān)度最高[18]。應(yīng)用地理探測器模型分析自然和人為因子對寧夏地區(qū)防風(fēng)固沙服務(wù)功能空間格局形成及演變過程中的貢獻與交互作用[1]發(fā)現(xiàn)土壤類型和植被類型在防風(fēng)固沙服務(wù)功能空間分布格局中發(fā)揮最重要的作用。而以歸一化植被指數(shù)(NDVI)計算植被覆蓋度,結(jié)合氣象要素從不同尺度來評估寧夏靈武白芨灘自然保護區(qū)防風(fēng)固沙功能時[12],植被覆蓋度的增加與研究區(qū)的防風(fēng)固沙服務(wù)成正比。朱趁趁等[14]利用RWEQ模型實現(xiàn)土地利用變化對內(nèi)蒙古荒漠草原的防風(fēng)固沙能力進行評估。結(jié)果表明：土地利用方式以林地恢復(fù)、建設(shè)用地擴張、不同覆蓋度草地轉(zhuǎn)化變化為主對防風(fēng)固沙服務(wù)有一定的增強作用。徐潔等[19]基于RWEQ和HYSPLIT模型模擬了生態(tài)功能區(qū)防風(fēng)固沙服務(wù)和防風(fēng)固沙服務(wù)空間流動路徑,從生態(tài)系統(tǒng)服務(wù)流動的角度建立了生態(tài)功能區(qū)及其防風(fēng)固沙服務(wù)受益區(qū)之間的時空聯(lián)系。而以往的研究大多是大區(qū)域尺度,很少評估土壤風(fēng)蝕對風(fēng)沙區(qū)邊緣地帶的影響,而且采用的氣象數(shù)據(jù)精度不夠。由于目前風(fēng)沙邊緣地區(qū)的氣象環(huán)境變化、植被恢復(fù)與防風(fēng)固沙能力之間的相互作用機制與如何保證防風(fēng)固沙服務(wù)功能的可持續(xù)性還需進一步研究。

基于生態(tài)工程的防風(fēng)固沙、植被恢復(fù)、土壤保育等生態(tài)功能監(jiān)測與評價,可為促進地區(qū)經(jīng)濟和社會的可持續(xù)發(fā)展提供實踐指導(dǎo)[20]。在我國土地沙化依然在擴展的趨勢下,陜北長城沿線處于毛烏素沙漠的邊緣地區(qū),長期受到自然因素和人類活動擾動影響,是最脆弱的生態(tài)系統(tǒng)之一[21]。目前,榆林大地基本上消滅了沙漠化土地現(xiàn)象,在全球沙漠治理的意義上率先實現(xiàn)了荒漠化的逆轉(zhuǎn)。但該地在干旱少雨、生態(tài)脆弱的自然條件下,導(dǎo)致頻發(fā)干旱和風(fēng)蝕現(xiàn)象。因此,從生態(tài)系統(tǒng)服務(wù)的角度探討植被恢復(fù)、氣象環(huán)境變化與防風(fēng)固沙服務(wù)能力之間相互作用機質(zhì)以及如何保證防風(fēng)固沙服務(wù)的可持續(xù)性。本文研究內(nèi)容包括：(1)分析陜北長城沿線2000—2018年植被覆蓋時空變化。(2)利用RWEQ模型模擬長城沿線生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)物理量的時空模式,實現(xiàn)防風(fēng)固沙服務(wù)的量化。(3)利用情景分析評估風(fēng)蝕因子改變的情況下,植被恢復(fù)對防風(fēng)固沙服務(wù)物理量和服務(wù)能力的影響。擬解決的關(guān)鍵問題：(1)以定量化研究的結(jié)果為長城沿線2000—2018年的風(fēng)蝕現(xiàn)象的變化情況提供科學(xué)依據(jù)。(2)探討風(fēng)蝕因子改變,植被恢復(fù)對風(fēng)沙邊緣地區(qū)風(fēng)蝕量的影響及生態(tài)系統(tǒng)防風(fēng)固沙物質(zhì)量和能力的影響。明確氣象環(huán)境改變時,植被恢復(fù)與防風(fēng)固沙服務(wù)的作用機制,并推進土地資源的保護和建設(shè),實現(xiàn)可持續(xù)高質(zhì)量發(fā)展，也為我國陜北地區(qū)防風(fēng)固沙治理工作進程提供一定的科學(xué)參考。

1 研究區(qū)及數(shù)據(jù)源

1.1 研究區(qū)概況

陜北長城沿線(圖1)位于陜西省最北端,包括定邊、靖邊、橫山、榆陽、神木、府谷六個市縣(區(qū))[21—22],北鄰鄂爾多斯草原,西鄰寧夏自治區(qū)。位于36°45′—39°58′N,107°35′—111°29′E之間。作為陜北黃土高原丘陵溝壑區(qū),地處干旱與半干旱地帶,地勢東低西高、年降水量偏少、土壤類型主要以沙土為主,為風(fēng)力侵蝕提供了“有利”的發(fā)生條件。作為干旱頻發(fā)區(qū),面臨著植被退化、水土流失等生態(tài)問題,生態(tài)系統(tǒng)脆弱易受人類活動干擾和氣候變化的影響。

圖1 研究區(qū)地理位置Fig.1 Geographical location of the study area

1.2 數(shù)據(jù)來源

本文使用的數(shù)據(jù)包括NDVI數(shù)據(jù)、數(shù)字高程模型(Digital Elevation Model,DEM)、土壤數(shù)據(jù)、雪蓋數(shù)據(jù)和氣象數(shù)據(jù)。(1)NDVI數(shù)據(jù)來源于美國地質(zhì)勘探局(http://glovis.usgs.gov/)數(shù)據(jù)產(chǎn)品,選用2000—2018年時間段,空間分辨率為1km,利用MRT(MODIS Reprojection Tools)最大合成年數(shù)據(jù)。(2)DEM數(shù)據(jù)來源于中國科學(xué)院資源與環(huán)境數(shù)據(jù)云平臺(http://www.resdc.cn)空間分辨率為1km。(3)土壤數(shù)據(jù)和積雪覆蓋數(shù)據(jù)來源于中國科學(xué)院旱區(qū)寒區(qū)科學(xué)數(shù)據(jù)中心(http://westdc.westgis.ac.cn),土壤數(shù)據(jù)為世界土壤數(shù)據(jù)庫(Harmonized World Soil Database,HWSD)所提供的中國土壤數(shù)據(jù)為1∶100萬土壤數(shù)據(jù)庫。中國雪深長時間序列數(shù)據(jù)集提供1979年1月1日到2020年12月31日逐日的中國范圍的積雪厚度分布數(shù)據(jù),空間分辨率為0.25°。(4)氣象數(shù)據(jù)(風(fēng)速、氣溫、降水)來源于時空三級環(huán)境大數(shù)據(jù)平臺(https://data.tpdc.ac.cn/zh-hans/)提供的中國區(qū)域地面要素數(shù)據(jù)集(China meteorological forcing dataset,CMFD)[23—25],該數(shù)據(jù)時間范圍為1979—2018年,時間分辨率為3小時,空間分辨率為0.1°,數(shù)據(jù)格式為NETCDF格式,運用Python語言提取出研究區(qū)內(nèi)349個規(guī)則點氣象要素,大于實際研究區(qū)4個國家氣象站點,足以從柵格尺度上真實的反映出研究區(qū)內(nèi)風(fēng)速大于5m/s天數(shù)；(5)土壤濕度數(shù)據(jù)來源于世界氣候中心(https://www.worldweather.cn/),空間分辨率為0.1°。最終將所有數(shù)據(jù)統(tǒng)一為CGCS_2000_Albers坐標(biāo)系,利用ArcGIS平臺重采樣為1km分辨率。

2 研究方法

2.1 Theil-Sen斜率估算與M-K顯著性檢驗結(jié)合

為了進一步體現(xiàn)出研究區(qū)植被變化的顯著性,采用斜率估算和顯著性檢驗相結(jié)合的方法。而Theil-Sen斜率估算和Mann-Kendall顯著性檢驗是兩種非參數(shù)檢驗方法,目前已經(jīng)得到廣泛的應(yīng)用[26]。因其計算效率高,能夠有效地測定時間變化趨勢的起始位置,并具有檢測范圍廣、定量化程度高等優(yōu)點,比一元線性回歸趨勢檢驗空間變化更為敏感。Theil-Sen斜率估算通過計算序列中連續(xù)兩期時間序列數(shù)據(jù)之間的斜率,將所有像元的數(shù)值對斜率的中值做時間序列的總變化趨勢。公式如下：

(1)

式中,β是所有數(shù)據(jù)對斜率的中值。其大小表示平均變化率,正負表示時間序列的變化趨勢。Kendall()為取中值函數(shù)。xj和xi分別表示時間序列中的第j項和第i項的值。但由于該斜率估算方法不能實現(xiàn)時間序列的趨勢顯著性判斷。通常使用Mann-Kendall顯著性檢驗來判斷時間序列的趨勢顯著性。

Mann-Kendall方法是由彭曼首次提出并使用,并由Kendall進行改進,對于序列X=(X1,X2,X3...Xn)先確定所有對偶值,之后比較與檢驗統(tǒng)計量S的大小關(guān)系。該方法的原始假設(shè)是：H0：時間序列中的數(shù)據(jù)隨機排列,即不存在顯著趨勢；H1：時間序列數(shù)據(jù)存在上升或下降的單調(diào)趨勢。公式如下：

(2)

(3)

由于時間序列維度大小不同,顯著性檢驗統(tǒng)計量的選取有所不同。當(dāng)時間序列小于10時,直接使用統(tǒng)計量S進行雙邊趨勢檢驗。在給定顯著性水平α下,如果|S|≥Sα/2則拒絕原假設(shè),則存在顯著趨勢。反之則不存在顯著性趨勢。當(dāng)時間序列≥10時,統(tǒng)計量S近似服從標(biāo)準(zhǔn)正態(tài)分布,可將S標(biāo)準(zhǔn)化后得到檢驗統(tǒng)計量Z,Z值計算公式如下：

(4)

(5)

式中,n是時間序列中數(shù)據(jù)個數(shù)；m是序列中重復(fù)出現(xiàn)的數(shù)據(jù)組個數(shù)；ti是第i組重復(fù)數(shù)據(jù)組中的重復(fù)數(shù)據(jù)個數(shù)。同樣,采用雙邊趨勢檢驗,在給定顯著性水平α下,在正態(tài)分布表中查得臨界值Z1-α/2。如果|Z|≤Z1-α/2時,則接受原假設(shè)即趨勢不顯著；反之,則相反。

2.2 生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)物質(zhì)量評價

植被所引起的風(fēng)蝕減少量可視為防風(fēng)固沙服務(wù)的物理量,其表示在裸土條件下潛在風(fēng)蝕量與植被覆蓋條件下實際風(fēng)蝕量之間的差[27]。在充分考慮氣候條件、植被覆蓋、土壤等要素情況下,可根據(jù)RWEQ模型定量評估潛在風(fēng)蝕量和實際風(fēng)蝕量。利用兩者的差值進行生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)物理量的評價[28]。其計算公式如下：

(6)

Qrmax=109.8×(WF×EF×SCF×K′)

(7)

Sr=150.71×(WF×EF×SCF×K′)-0.3711

(8)

(9)

Qmax=109.8×(WF×EF×SCF×K′×C)

(10)

S=150.71×(WF×EF×SCF×K′×C)-0.3711

(11)

G=SLr-SL

(12)

式中,SLr表示單位面積年潛在風(fēng)蝕量(kg/m2)；Qrmax為潛在風(fēng)力最大輸沙能力(kg/m)；Sr為潛在關(guān)鍵地塊長度(m)；SL表示單位面積年實際風(fēng)蝕量(kg/m2)；Qmax為風(fēng)力最大輸沙能力(kg/m);S為關(guān)鍵地塊長度(m)；G表示單位面積防風(fēng)固沙量(kg/m2);z表示下風(fēng)向距離(m),本文計算取50m；WF表示風(fēng)蝕因子(kg/m)；EF表示土壤可蝕性因子(無量綱)；SCF表示為土壤結(jié)皮因子(無量綱)；K′表示為土壤粗糙度因子(無量綱)；C表示植被因子(無量綱)。具體表達式如下：

風(fēng)蝕因子WF在考慮風(fēng)速作用的前提下,還將積雪覆蓋對地表的保護和土壤濕度納入考量,計算公式為：

WF=wf×(ρ/g)×SW×SD

(13)

wf=u2(u2-u1)2×Nd

(14)

式中,wf表示風(fēng)力系數(shù)(m/s3)；g表示重力加速度(m/s2),本文取9.8 m/s2；ρ表示空氣密度(kg/m3),本文中采用20℃時的空氣密度值,即1.205kg/m3；SW表示土壤濕度因子(無量綱)；SD表示積雪覆蓋因子(無量綱),即無積雪覆蓋的天數(shù)與研究總天數(shù)之比,積雪深度小于25.4mm時,表示無積雪；u1為最小起沙風(fēng)速,模型中默認為5m/s;u2表示2m高度處的風(fēng)速(m/s)；Nd為每月風(fēng)速超過閾值風(fēng)速的天數(shù)。

土壤可蝕性因子EF指在一定土壤理化條件下土壤受風(fēng)蝕影響大小。土壤結(jié)皮因子SCF指在一定理化條件下土壤結(jié)皮抵抗風(fēng)蝕能力的大小。其表達式為：

(15)

(16)

式中,sa為土壤砂粒含量(%);si為土壤粉粒含量(%)；cl為土壤粘粒含量(%)；OM為土壤有機質(zhì)含量(%)；CaCO3為碳酸鈣含量(%)。假定土壤可蝕性因子和土壤結(jié)皮因子隨時間保持不變,在應(yīng)用世界土壤數(shù)據(jù)庫的基礎(chǔ)上,結(jié)合陜北長城沿線風(fēng)沙邊緣區(qū)土壤質(zhì)地特征,參考RWEQ模型內(nèi)嵌土壤資料表,計算出土壤可蝕性因子和土壤結(jié)皮因子。

地表粗糙度因子K′是由地形引起的土地表面粗糙程度對土壤風(fēng)蝕的影響,公式如下：

K′=cosα

(17)

式中,α表示為坡度。

植被因子C表示一定的植被覆蓋對土壤風(fēng)蝕的抑制程度,公式如下：

C=e-0.0483(FVC)

(18)

(19)

式中,FVC代表植被覆蓋度(%),NDVIveg和NDVIsoil分別表示植被覆蓋的最大和最小值。本文利用植被指數(shù)法,基于像元二分模型計算得出95%和5%的累積頻率對應(yīng)下的植被覆蓋最大和最小值[29—30]。

生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)功能可以表明植被的實際固沙的能力[31]。利用研究區(qū)面積,定量化出防風(fēng)固沙物理量,以此來評估區(qū)域生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)能力的大小。計算公式為：

QSR=10×G×A

(20)

ASR=QSR/A

(21)

式中,QSR為研究區(qū)防風(fēng)固沙量(t/a)；A為研究區(qū)面積(hm2)；ASR為生態(tài)系統(tǒng)防風(fēng)固沙能力(t/hm2)。

3 結(jié)果與分析

3.1 Theil-Sen斜率估算與M-K檢驗

采用Theil-Sen斜率估算和M-K顯著性檢驗相結(jié)合的方法,利用NDVI數(shù)據(jù)得到的陜北長城沿線風(fēng)沙邊緣區(qū)植被變化情況(圖2)。結(jié)果表明：2000—2018年研究區(qū)內(nèi)98.74%的區(qū)域植被變化表現(xiàn)為改善。其中,約有90.34%植被覆蓋區(qū)域的植被狀況呈現(xiàn)出顯著改善。結(jié)合哨兵二號影像數(shù)據(jù)對比來看,榆陽區(qū)、神木市等因城市擴展,導(dǎo)致植被出現(xiàn)退化。不顯著改善集中表現(xiàn)在定邊縣占比較大,說明風(fēng)沙邊緣地區(qū)的植被還需進一步保護。

雖然Theil-Sen斜率估算和M-K檢驗方法有計算效率高、檢測范圍廣、定量化程度高等優(yōu)點,比一元線性回歸趨勢檢驗空間變化更為敏感。但本文檢測到的陜北長城沿線風(fēng)沙邊緣地區(qū)植被變化空間分布模式與基于一元線性回歸分析方法檢測到的植被空間變化模式大體上一致。

3.2 風(fēng)蝕模型因子時空格局

3.2.1風(fēng)蝕因子結(jié)果

從表1可以看出,2000—2018年風(fēng)蝕因子呈現(xiàn)出減少-增多-減少的趨勢,且風(fēng)蝕因子均值從24.02m/s波動式降至12.26m/s,整體上相對減小。在2010年風(fēng)蝕因子均值達到最大,為74.13m/s,2005年最低,為6.32m/s。從各點大于5m/s的數(shù)值和天數(shù)來看,2010年相對多于其他年份,致使2010年的風(fēng)蝕因子出現(xiàn)最大值。

表1 2000—2018年風(fēng)蝕因子變化

在空間分布來看(圖3),風(fēng)蝕因子分布差異明顯,且風(fēng)蝕因子高值分布在榆陽區(qū)和神木縣周邊。對比分析,2010年較其他年份差異性更為顯著。在2005年定邊縣風(fēng)蝕因子高于其他五個縣區(qū),統(tǒng)計分析這是由于當(dāng)年定邊縣周邊氣象要素點大于起沙風(fēng)速閾值的天數(shù)和個數(shù)偏多。

圖2 2000—2018年植被變化趨勢及顯著性檢驗Fig.2 Vegetation change trend and significance test from 2000 to 2018

圖3 2000—2018年風(fēng)蝕因子空間分布Fig.3 Spatial distribution of wind erosion factors from 2000 to 2018

3.2.2土壤及地形因子結(jié)果

土壤可蝕性因子表征土壤抵抗侵蝕能力的大小,土壤對風(fēng)蝕作用敏感程度越大,可蝕性越高。陜北長城沿線風(fēng)沙邊緣區(qū)EF取值(圖4)在0.26—0.61,平均值為0.48。在空間分布上,榆陽地區(qū)土壤對風(fēng)蝕作用較敏感,EF取值較高,而定邊縣和府谷縣EF取值較低。從土壤理化性質(zhì)出發(fā),紅砂土的EF值最高,這是因為其沙粒含量達88%以上,容易發(fā)生風(fēng)蝕現(xiàn)象。潛育土EF值最低且為0.266,因其粘粒和有機質(zhì)含量達到43.7%相比于其他土壤類型不易起沙。

土壤結(jié)皮因子與土壤易受風(fēng)力侵蝕威脅成反比,即土壤表層結(jié)皮越堅固,土壤結(jié)皮因子SCF值越大,土壤越不易受到風(fēng)蝕威脅。整個研究區(qū)SCF的取值(圖4)在0.06—1.00,平均值為0.56。SCF越高抵抗風(fēng)蝕能力越差。從空間分布來看,SCF高值分布在長城沿線分割線以北,與內(nèi)蒙古毛烏素沙漠相鄰且抵抗風(fēng)蝕能力較差。從土壤理化性質(zhì)來看,紅沙土粘粒含量和有機質(zhì)含量最低,導(dǎo)致SCF最高,為0.837,抗風(fēng)蝕能力最差。始成土的粘粒和有機質(zhì)含量相對較高,SCF最低,說明其抗風(fēng)蝕能力最強,不易發(fā)生風(fēng)蝕。

地表粗糙度因子K′ 反映了地形對風(fēng)蝕強度大小的影響。長城沿線K′(圖4)取值范圍為0.004—1.00,均值為0.65。空間分布上,長城沿線分割線以北地表粗糙度占大部分,主要分布在定邊縣北部、靖邊縣、榆陽區(qū)、神木縣。根據(jù)海拔高度和地形來看,定邊縣和靖邊縣海拔北低南高,榆陽區(qū)、神木縣西北高東南低,且兩兩縣區(qū)高值區(qū)和低值區(qū)相近,造成地表粗糙度高,更易于侵蝕的發(fā)生。且定邊縣和靖邊縣南部山地海拔較高和起伏較大,不利于侵蝕的發(fā)生。

圖4 土壤因子分布Fig.4 Distribution of soil factors

3.2.3植被因子結(jié)果

植被因子表征植被對風(fēng)蝕抑制程度的大小,植被生長越好對風(fēng)蝕的抑制程度越大進而影響生態(tài)系統(tǒng)防風(fēng)固沙能力。研究區(qū)C因子取值范圍為0.007—1.00,植被整體呈現(xiàn)好轉(zhuǎn)趨勢但均值呈輕微下降趨勢,局部存在波動?？臻g分布上(圖5),C因子較高值主要分布在神木市和榆陽區(qū)等地,表現(xiàn)出易被風(fēng)力侵蝕,C因子較低值則相反。

圖5 2000—2018植被因子Fig.5 The vegetation factor from 2000 to 2018

3.3 防風(fēng)固沙物質(zhì)量及時空格局

3.3.1潛在與實際風(fēng)蝕時空變化

潛在風(fēng)蝕能夠反映出土壤可能的最大風(fēng)蝕量,即土壤風(fēng)蝕危害的“危險程度”。從表2可看出,2000—2018年潛在風(fēng)蝕量在4569.18—64164.44萬t,2010年達到潛在風(fēng)蝕量峰值,2005年最少。與2000年相比,2018年減少了13160.93萬t?？梢婈儽遍L城沿線土壤風(fēng)蝕危害的“危險程度”在呈波動式降低。與2000年相比,2018年平均單位面積潛在風(fēng)蝕量下降了3.92kg/m2,總體呈“減少-增加-減少”的波動式變化。2000—2018年實際風(fēng)蝕量在1056.31—24211.89萬t,在2010年達到實際風(fēng)蝕量峰值(24211.89萬t),2005年最少。與2000年相比,2018年增加了894.94萬t。與2010年相比,2018年減少了20161.58萬t,平均單位面積風(fēng)蝕量減少了6.0kg/m2,減少幅度較大,整體呈“減少-增加-減少”的變化特點。

2000—2018年潛在風(fēng)蝕出現(xiàn)顯著的空間差異性。根據(jù)各個因子空間分布,在風(fēng)蝕因子分布較大、地表較粗糙、土壤結(jié)皮較大的地區(qū)潛在風(fēng)蝕強度最為劇烈且變化幅度最為明顯,如榆陽區(qū)和神木市。而定邊縣、靖邊縣和府谷縣潛在風(fēng)蝕強度較輕。實際風(fēng)蝕年際空間變化上,也存在明顯的物質(zhì)量變化。依據(jù)水利部《土壤風(fēng)蝕分類標(biāo)準(zhǔn)(SL190—2007)》將土壤風(fēng)蝕分為微度(<200t/a)、輕度(200—2500t/a)、中度(2500—5000t/a)、強烈(5000—8000t/a)、極強烈(8000—15000t/a)、劇烈(>15000t/a)[32—33]。根據(jù)分類統(tǒng)計,微度土壤風(fēng)蝕面積增加了16.16%,輕度風(fēng)蝕面積減少了15.94%,2010年強烈、極強烈、劇烈風(fēng)蝕面積占比最大,與此對應(yīng)該年的風(fēng)蝕因子最大,導(dǎo)致風(fēng)蝕現(xiàn)象嚴重。2010到2018年,風(fēng)蝕因子逐漸減小,土壤風(fēng)蝕強度也逐漸有了明顯的改善。定邊縣、靖邊縣和府谷縣微度風(fēng)蝕強度占比較多；榆陽區(qū)和神木市潛在風(fēng)蝕最為劇烈且變化幅度最為明顯,也是單位面積實際風(fēng)蝕量較嚴重的區(qū)域(圖6)。從時間統(tǒng)計和空間分布比較來看,每年的潛在風(fēng)蝕量>防風(fēng)固沙量>實際風(fēng)蝕量,說明在最大危害程度內(nèi),防風(fēng)固沙服務(wù)發(fā)揮明顯作用。而從整體情況來看,在2000—2018年陜北長城沿線防風(fēng)固沙治理過程中,土壤風(fēng)蝕現(xiàn)象已得到顯著遏制。相比于潛在風(fēng)蝕,實際風(fēng)蝕多考慮植被因子,由于其他因定不變,主要是風(fēng)蝕和植被因子對風(fēng)蝕的影響最密切。

表2 2000—2018年防風(fēng)固沙物質(zhì)量統(tǒng)計/(kg/m2)

圖6 2000—2018年風(fēng)蝕強度空間分布Fig.6 Spatial distribution of wind erosion intensity from 2000 to 2018

3.3.2防風(fēng)固沙服務(wù)物質(zhì)量時空變化

潛在與實際風(fēng)蝕差值表征土壤風(fēng)蝕模數(shù)(單位面積防風(fēng)固沙物理量)。研究表明：防風(fēng)固沙物質(zhì)量從2000年的21944.56萬t減少到2018年的7987.44萬t。2010年變化波動增長較大且達到最大固沙量。依據(jù)張彪等[31]將防風(fēng)固沙服務(wù)能力ASR分為：低值區(qū)(ASR≤5t/hm2)、較低區(qū)(5≤ASR≤20t/hm2)、一般區(qū)(20≤ASR≤60t/hm2)、較高區(qū)(60≤ASR≤120t/hm2)和高值區(qū)(ASR≥120t/hm2)。從防風(fēng)固沙服務(wù)能力分區(qū)統(tǒng)計來看,如表3所示,2000—2018年生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)能力高值區(qū)和較高值區(qū)出現(xiàn)減小-增多-減小的變化特點,且整體上是減小的趨勢,且對比觀察不同年份的防風(fēng)固沙強度的分布規(guī)律即存在相似又存在差異。長城沿線生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)能力表現(xiàn)出明顯的空間異質(zhì)性,防風(fēng)固沙服務(wù)能力高值區(qū)主要集中在神木市和榆陽區(qū),如圖7所示。根據(jù)模型分析,如果風(fēng)蝕因子偏大且植被因子降低,即風(fēng)場強度增加導(dǎo)致防風(fēng)固沙量顯著減少。風(fēng)速的不穩(wěn)定性會導(dǎo)致風(fēng)蝕因子結(jié)果出現(xiàn)較大偏差,進而對防風(fēng)固沙服務(wù)的空間分布變化產(chǎn)生較大的影響。在抑制土壤風(fēng)蝕現(xiàn)象發(fā)生的基礎(chǔ)上,固沙能力也在減小,后期需要加以重點保護植被。

表3 2000—2018年防風(fēng)固沙能力空間分布/km2

圖7 2000—2018年防風(fēng)固沙能力空間分布Fig.7 Spatial distribution of wind prevention and sand fixation capacity from 2000 to 2018

3.4 植被恢復(fù)對土壤風(fēng)蝕及防風(fēng)固沙服務(wù)影響

為了進一步探究風(fēng)蝕因子變化,植被恢復(fù)對風(fēng)蝕現(xiàn)象及生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)的影響,利用情景模擬法,以最小(2005年)、最大風(fēng)蝕(2010年)因子為基準(zhǔn),分設(shè)兩種情景(即情景1和2)。由于計算潛在風(fēng)蝕時沒有考慮到植被覆蓋,其他因子固定不變,即：2005年的潛在風(fēng)蝕為情景1條件下的潛在風(fēng)蝕,2010年的潛在風(fēng)蝕為情景2條件下的潛在風(fēng)蝕。從物質(zhì)定量化分析角度出發(fā)(表4),橫向?qū)Ρ惹榫?得出：最小風(fēng)蝕因子情況下,植被長時間序列恢復(fù),植被因子降低,抵抗風(fēng)蝕能力增強,實際風(fēng)蝕量從1066.56萬t減少到825.38萬t,相比減少了241.18萬t,防風(fēng)固沙量在逐年提升的同時防風(fēng)固沙服務(wù)能力也顯著提升。與情景1相比,情景2中風(fēng)蝕因子的改變使得土壤風(fēng)蝕量從2000年的11764.29萬t波浪式增加到2018年的22649.44萬t,說明植被不再是影響風(fēng)蝕的主要因素,計算風(fēng)蝕因子的土壤濕度、風(fēng)力系數(shù)等因子的改變使得地表更不易抵抗風(fēng)蝕。但縱向?qū)Ρ惹榫?和2：隨著風(fēng)蝕因子的顯著變化,土壤風(fēng)蝕危害的“危險程度”達到最大(潛在風(fēng)蝕),以2000年為例,實際風(fēng)蝕量從1066.56萬t顯著增長到11764.29萬t,防風(fēng)固沙量從3449.6萬t變到52175.12}。由此可見：其他因素不變,氣象條件惡劣的情況下,極易發(fā)生土壤風(fēng)蝕。但兩種情境下的每年防風(fēng)固沙量>實際風(fēng)蝕量,利用植被固沙的效果也越來越明顯,生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)功能也在增加。

表4 2000—2018年分情景下防風(fēng)固沙物質(zhì)量統(tǒng)計/(kg/m2)

空間分布上(圖8),在風(fēng)蝕因子條件下(情景1),植被恢復(fù),2000—2018年風(fēng)蝕現(xiàn)象逐漸退化。在最大風(fēng)蝕因子下(情景2),造成風(fēng)蝕現(xiàn)象嚴重,且主要集中在榆陽區(qū)和神木縣,導(dǎo)致定邊縣和靖邊縣受到殃及。這四個縣區(qū)的防風(fēng)固沙量也出現(xiàn)明顯差異。可見風(fēng)蝕因子是主導(dǎo)土壤風(fēng)蝕的主導(dǎo)因素。整體來看：采取“鎖綠”來防沙治沙是最為有效的,且隨著風(fēng)蝕因子的改變水土流失治理成效越為明顯,生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)能力也越強。植被恢復(fù)對長城沿線的土壤風(fēng)蝕起到有效的抑制,對生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)起到了一定的增強作用。但風(fēng)蝕因子中風(fēng)速屬于不可防控的氣象因素,只能依靠恢復(fù)植被來增強生態(tài)系統(tǒng)防風(fēng)固沙能力。

圖8 分情境下風(fēng)蝕空間分布Fig.8 The amount of wind erosion under different scenarios

4 結(jié)論

本文基于RWEQ模型分析了陜北長城沿線防風(fēng)固沙服務(wù)物理量的時空模式,實現(xiàn)風(fēng)沙邊緣地區(qū)生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)的量化,并評估風(fēng)蝕和植被因子變化下對生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)物理量和服務(wù)能力的影響。得到以下結(jié)論：(1)通過分析2000—2018年陜北長城沿線六個縣區(qū)年際NDVI的時空變化,風(fēng)沙區(qū)內(nèi)植被呈明顯增加趨勢；(2)2000—2018年整個風(fēng)沙區(qū)實現(xiàn)了從“沙進人退”到“綠進沙退”的轉(zhuǎn)變,土壤風(fēng)蝕危害的“危險程度”也在逐步降低,實際風(fēng)蝕和生態(tài)系統(tǒng)防風(fēng)固沙總量2005年和2010年分別為4569.18萬t和64164.44萬t,且潛在風(fēng)蝕量>防風(fēng)固沙量>實際風(fēng)蝕量,生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)發(fā)揮顯著作用。(3)利用情景變化定量分析,探索植被恢復(fù)可以一定程度上抑制風(fēng)蝕現(xiàn)象。在最小風(fēng)蝕因子情況下,不受外力條件變化,植被防風(fēng)固沙服務(wù)的作用發(fā)揮不明顯,且一定程度抑制風(fēng)蝕現(xiàn)象的發(fā)生。風(fēng)蝕因子改變加重研究區(qū)的土壤風(fēng)蝕,植被所發(fā)揮的生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)越為明顯。而兩因子相結(jié)合(即實際情況下的風(fēng)蝕量),風(fēng)蝕因子占主導(dǎo)地位,對土壤風(fēng)蝕和防風(fēng)固沙量影響作用較大。不論是真實情況還是分情景下的風(fēng)蝕,生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)空間差異顯著,榆陽區(qū)和神木市的土壤風(fēng)蝕現(xiàn)象最為嚴重,也是生態(tài)系統(tǒng)防風(fēng)固沙服物質(zhì)量和能力高值集中區(qū)?？傊?恢復(fù)植被是實現(xiàn)區(qū)域性荒漠化逆轉(zhuǎn)和改善生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)與能力的必要條件,恢復(fù)植被不僅可以減少土壤風(fēng)蝕,也可以增強生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)能力。

5 討論

(1)當(dāng)前土壤風(fēng)力侵蝕的定量化驗證還存在一定的困難。實測的土壤風(fēng)蝕量是開展研究結(jié)果驗證的基礎(chǔ),由于沒有實測的土壤風(fēng)蝕量數(shù)據(jù),無法對利用模型模擬的實際風(fēng)蝕量結(jié)果進行驗證,后續(xù)可進一步采用沙塵暴發(fā)生頻率來驗證本文的可靠性。鑒于此,文中利用情景模擬法可以間接進行驗證,對比分析植被恢復(fù)對風(fēng)蝕現(xiàn)象有一定的遏制作用,與劉碩等[16]的研究結(jié)果一致。在改變風(fēng)力因子時,風(fēng)蝕現(xiàn)象依舊明顯,植被發(fā)揮的生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)也顯著提高。

(2)提高模型的優(yōu)化性。在實現(xiàn)陜北長城沿線2000—2018年生態(tài)系統(tǒng)防風(fēng)固沙量以及服務(wù)能力,由于實測氣象站點少,而采用中國區(qū)域地面要素數(shù)據(jù)集,但該數(shù)據(jù)集時間范圍限制,未能及時做到最新年份；考慮到該數(shù)據(jù)集時間分辨率較大,且有349個規(guī)則分布的點,在選用空間分辨率時采用1km空間分辨率,但對于土壤風(fēng)蝕研究而言,1km分辨率足以反映其空間差異。在實現(xiàn)雪蓋因子計算的過程中,計算的是年積雪覆蓋天數(shù)的概率,得到年尺度風(fēng)蝕量。在后期研究細化時間尺度,從月尺度考慮生態(tài)系統(tǒng)防風(fēng)固沙物質(zhì)量和服務(wù)能力。

(3)本文的研究結(jié)論僅從生態(tài)系統(tǒng)防風(fēng)固沙量和服務(wù)能力角度出發(fā),在情景1條件下植被的固沙作用效果不明顯,但在最大風(fēng)蝕因子條件下,植被對抑制風(fēng)沙吹動,防風(fēng)固沙效果越來越明顯。三種情景每年的潛在風(fēng)蝕量>防風(fēng)固沙量>實際風(fēng)蝕量。風(fēng)蝕因子增大,潛在風(fēng)蝕增加、實際風(fēng)蝕增加、防風(fēng)固沙量增加,生態(tài)系統(tǒng)防風(fēng)固沙服務(wù)功能也增強。這與治沙研究者提出防沙治沙、植樹造林工程以來的預(yù)想一致。從“沙進人退”到“綠進沙退”的轉(zhuǎn)變是數(shù)十年治沙的功勞,“榆林綠”鎖住了毛烏素沙地。作為自然地貌,毛烏素沙地不是真的消失,只是流動沙丘全部得到固定。但這些治理成果還屬人工干預(yù),尚需加快沙區(qū)生態(tài)系統(tǒng)的正向演替,實現(xiàn)自我循環(huán)發(fā)展,從根本上改變沙區(qū)生態(tài)面貌。

本文研究特色在于：采用中國區(qū)域地面要素數(shù)據(jù)集,更好的實現(xiàn)長城沿線風(fēng)蝕因子的空間分布。但還存在不足之處：在考慮土壤理化性質(zhì)隨時間不發(fā)生變化的情況下,但實際情況中,由于陜北長城沿線六個縣區(qū)耕地面積占比較大,長期使用化肥致使土壤理化性質(zhì)發(fā)生輕微改變,導(dǎo)致文章實驗結(jié)果的可靠性存在一定的誤差,土壤理化性質(zhì)的改變對土壤風(fēng)蝕具有一定的影響,對風(fēng)沙區(qū)土壤理化性質(zhì)的測定研究也是今后的研究的一個重要科學(xué)問題。在可控的人為范圍內(nèi),做好治沙工程,堅持不懈地防沙治沙、植樹造林使得生態(tài)環(huán)境持續(xù)改善,為榆林經(jīng)濟社會的快速發(fā)展提供了強有力的基礎(chǔ)保障。