王凱利,王志慧,肖培青,*,王鐵生,張 攀

1 華北水利水電大學(xué)測(cè)繪與地理信息學(xué)院,鄭州 450045 2 黃河水利科學(xué)研究院,水利部黃土高原水土保持重點(diǎn)實(shí)驗(yàn)室,鄭州 450003

黃土高原是全球水土流失最為嚴(yán)重的區(qū)域,生態(tài)系統(tǒng)極為脆弱[1],隨著國(guó)家相繼實(shí)施退耕還林還草、封山禁牧等生態(tài)修復(fù)工程以來(lái),植被覆蓋度已由1999年的31.6%提高到2018年的63%[2]；與此同時(shí)黃河中游徑流泥沙銳減[3—4],生態(tài)環(huán)境明顯改善[5—6]。已有研究表明,植被恢復(fù)是黃土高原入黃泥沙減少的重要原因[7—11]。但是,大規(guī)模的植被恢復(fù)使更多的降水被植被冠層截留蒸發(fā)[12],同時(shí)也導(dǎo)致了植被耗水量的劇烈增加[13—14],陸地長(zhǎng)期蓄水變化量減少[15],部分造林密度較大區(qū)的深層土壤已出現(xiàn)干層現(xiàn)象[16—17]。因此,在氣候變化和人類活動(dòng)的雙重驅(qū)動(dòng)影響下,植被生長(zhǎng)耗水需求與天然降雨補(bǔ)給之間的矛盾與日俱增。在土壤水分長(zhǎng)期虧損狀態(tài)下植被生態(tài)系統(tǒng)會(huì)發(fā)生嚴(yán)重退化,植被減沙功能則會(huì)顯著降低,這將對(duì)黃土高原流域未來(lái)入黃泥沙量產(chǎn)生重大影響。因此,提出能夠維持黃土高原植被可持續(xù)性演替的植被恢復(fù)閾值,對(duì)于黃土高原生態(tài)建設(shè)、水土保持措施科學(xué)布局以及未來(lái)黃河流域水沙變化趨勢(shì)預(yù)測(cè)都具有非常重要的科學(xué)與現(xiàn)實(shí)意義。

由于植被恢復(fù)潛力被認(rèn)為是植被恢復(fù)閾值與現(xiàn)狀植被的差距,因而計(jì)算植被恢復(fù)閾值對(duì)于評(píng)估植被恢復(fù)潛力至關(guān)重要[18]。王延平等[19]利用水量平衡模型估算陜北米脂紫花苜蓿的植被生物量閾值為2600—3500 kg/hm2；Fu等[20]通過(guò)SWCCV模型確定六道溝流域錦雞兒和沙柳的LAI閾值為1.27和0.7；劉丙霞等[21]運(yùn)用SHAW模型模擬發(fā)現(xiàn)神木六道溝流域的苜蓿和檸條生物量閾值分別為1980 kg/hm2和5050 kg/hm2。這些研究主要聚焦于黃土高原地區(qū)典型灌草物種在樣地或小流域尺度的恢復(fù)閾值,難以反映整個(gè)黃河流域植被恢復(fù)閾值的空間分異特征。遙感技術(shù)的快速發(fā)展為區(qū)域尺度植被恢復(fù)閾值的確定提供了新的技術(shù)方法。趙廣舉[2]、高海東等[18]基于“相似生境原則”,利用直方圖統(tǒng)計(jì)法提取了不同生境斑塊的植被恢復(fù)閾值。為了解決生境條件劃分標(biāo)準(zhǔn)的空間異質(zhì)性問(wèn)題,Zhang等[22]提出了“基于滑動(dòng)窗口的相似生境植被恢復(fù)閾值”改進(jìn)計(jì)算方法。這類方法雖簡(jiǎn)單易行,但在劃分具有“相似生境”斑塊時(shí)卻難以充分考慮微地形變化對(duì)生境條件的影響,且提取的閾值僅代表研究時(shí)段斑塊中植被參數(shù)最大值,并不代表維持植被耗水與降雨供給平衡的植被閾值。Feng等[23]構(gòu)建了黃土高原生態(tài)系統(tǒng)蒸散發(fā)(ET)與總初級(jí)生產(chǎn)力(GPP)的定量響應(yīng)關(guān)系,并模擬了不同生態(tài)系統(tǒng)可供水量情景下的GPP閾值。Liang等[24]在此基礎(chǔ)上進(jìn)一步考慮了立地條件異質(zhì)性對(duì)ET與GPP定量響應(yīng)關(guān)系的影響。但GPP表征植被光合作用固碳能力(gC m-2a-1),該指標(biāo)雖與ET存在緊密生物物理關(guān)聯(lián),卻難以直接應(yīng)用于生態(tài)建設(shè)工程實(shí)踐。Zhang等[25]利用Eagleson生態(tài)水文模型與生態(tài)最優(yōu)理論模擬了不同氣候情景下黃土高原植被覆蓋度恢復(fù)閾值及其空間分布。由于生態(tài)水文模型計(jì)算過(guò)程復(fù)雜,且模型中土壤和植被參數(shù)難以準(zhǔn)確獲取,對(duì)模擬精度帶來(lái)較大不確定性。

為此,本研究利用已有ET和植被參數(shù)遙感產(chǎn)品,在考慮ET對(duì)植被的線性和非線性響應(yīng)關(guān)系以及植被與氣象要素相互影響作用的基礎(chǔ)上,利用逐步回歸方法逐像元構(gòu)建ET與氣象要素和植被指標(biāo)的最優(yōu)響應(yīng)關(guān)系,基于水量平衡原理計(jì)算不同氣候條件下林草植被覆蓋度恢復(fù)閾值。以期為黃土高原生態(tài)建設(shè)提供重要參考與科學(xué)技術(shù)支撐。

1 研究區(qū)與方法

1.1 研究區(qū)概況

黃土高原地勢(shì)東南低西北高,地形破碎,溝壑縱橫,屬半干旱大陸性季風(fēng)氣候,年降水量為443.93 mm,年均氣溫為9.57℃(圖1),年均潛在蒸散發(fā)高達(dá)1319.57 mm。植被類型呈現(xiàn)地帶性規(guī)律分布,從東南到西北依次為森林帶、森林草原帶、典型草原帶、荒漠草原帶和草原化荒漠帶。受土壤、降雨、地形地貌、人類活動(dòng)等的影響,植被覆蓋度較低,水土流失嚴(yán)重。為改善區(qū)域生態(tài)環(huán)境,自1999年開(kāi)始經(jīng)過(guò)多年的生態(tài)修復(fù)治理,黃土高原植被恢復(fù)效果顯著,水土流失得到有效抑制[26]。

圖1 黃土高原年均降雨量與年均溫度Fig.1 Annual mean precipitation and annual mean temperature on the Loess Plateau

1.2 數(shù)據(jù)來(lái)源與處理

氣象數(shù)據(jù)采用中國(guó)氣象數(shù)據(jù)網(wǎng)(https://data.cma.cn/)共享發(fā)布的1960—2018年黃土高原295個(gè)氣象站點(diǎn)日尺度氣象數(shù)據(jù),包括降雨、溫度、日照時(shí)數(shù)、相對(duì)濕度、風(fēng)速。利用AUSPLINE氣象專業(yè)插值軟件將所有氣象站點(diǎn)數(shù)據(jù)進(jìn)行空間插值,并利用FAO Penman-Monteith公式[27—28]和Yuan等[29]的方法計(jì)算潛在蒸散發(fā)(ET0)和飽和水汽壓差(VPD)。另外,本研究獲取了黃土高原及鄰近區(qū)域的4個(gè)Flux通量站(Du2、Ha2、HaM、Sw2)ET觀測(cè)數(shù)據(jù),黃河干流5個(gè)主要控制水文站(唐乃亥、青銅峽、頭道拐、龍門、花園口)年徑流量實(shí)測(cè)數(shù)據(jù)(圖2),數(shù)據(jù)測(cè)量準(zhǔn)確度高。其中,唐乃亥水文站以上為黃河源區(qū),青銅峽水文站以上為黃河流域主要產(chǎn)水區(qū),頭道拐水文站為黃河上游與中游的分界,龍門水文站為黃河重點(diǎn)報(bào)汛站,花園口水文站為黃河中游與下游的分界。而“唐乃亥-青銅峽”為黃河上游主要來(lái)水區(qū),“青銅峽-頭道拐”為黃河沙漠寬谷段,“頭道拐-龍門”為黃河中游主要產(chǎn)沙區(qū),“龍門-花園口”為黃河中游主要產(chǎn)水區(qū)；四個(gè)區(qū)間覆蓋了黃土高原90%以上區(qū)域,不同區(qū)間的氣候、地形地貌和植被分布均具有典型特征,區(qū)域蒸散發(fā)可代表整個(gè)黃土高原蒸散發(fā)及其空間異質(zhì)性。干旱指數(shù)r為多年平均潛在蒸散發(fā)與多年平均降雨量的比值,根據(jù)r直方圖分布[30]將其劃分為0.85—2、2—2.5、2.5—3、3—4、4—16共5個(gè)干旱指數(shù)區(qū)間(圖2),以保證不同區(qū)間的像元個(gè)數(shù)大致相同。

ET遙感產(chǎn)品選用PML_V2和MOD16 A2 GF,空間分辨率為500 m。由于增強(qiáng)型植被指數(shù)(Enhanced Vegetation Indices,EVI)已消除土壤背景和大氣對(duì)光譜反射率的影響,且能克服高植被覆蓋區(qū)植被指數(shù)飽和效應(yīng)問(wèn)題,比其它植被指數(shù)更能準(zhǔn)確捕捉植被時(shí)空動(dòng)態(tài)變化[22],選用MODIS13 A1 EVI表征植被參數(shù),空間分辨率為500 m。地表覆蓋數(shù)據(jù)采用GLC_FCS30產(chǎn)品,空間分辨率為30 m,該產(chǎn)品全球范圍驗(yàn)證精度82.50%,kappa系數(shù)為0.78[31]；將其重分類為農(nóng)田、森林、灌草地、其它(圖2),其中森林與灌草總稱為林草植被[32]。所有遙感產(chǎn)品均重采樣到500 m分辨率。利用GRACE產(chǎn)品CSR RL06_mascons提取不同區(qū)間(唐乃亥-青銅峽、青銅峽-頭道拐、頭道拐-龍門、龍門-花園口)的統(tǒng)計(jì)平均值,即陸地水儲(chǔ)量的變化量(TWSC)。

圖2 黃土高原地面站點(diǎn)及自然概況Fig.2 Ground stations and natural survey of the Loess Plateau

1.3 研究方法

1.3.1ET產(chǎn)品精度評(píng)價(jià)方法

提取Flux站坐標(biāo)鄰近3×3像元(1.5 km×1.5 km)ET產(chǎn)品空間平均值,并利用均方根誤差(RMSE)、平均相對(duì)誤差(MRE)和納什效率系數(shù)(Nash)指標(biāo)評(píng)估PML_V2 ET和MOD16 A2 GFET產(chǎn)品精度。另外,在唐乃亥-青銅峽、青銅峽-頭道拐、頭道拐-龍門、龍門-花園口4個(gè)區(qū)間上分別利用水量平衡原理計(jì)算區(qū)域尺度ET,計(jì)算公式為[33]：

ETWBi=Pi-Ri-ΔSi

(1)

式中,ETWBi為第i年區(qū)間蒸散發(fā)量(mm)；Pi為第i年區(qū)間實(shí)測(cè)降雨量(mm)；Ri為第i年區(qū)間實(shí)測(cè)徑流量(mm)；ΔSi為第i年GRACE重力衛(wèi)星反演的區(qū)間TWSC(mm)。

以基于水量平衡原理計(jì)算ET為區(qū)域ET實(shí)測(cè)值,利用RMSE和MRE指標(biāo)評(píng)估PML_V2和MOD16 A2 GF ET產(chǎn)品在區(qū)域尺度上的模擬精度。

1.3.2趨勢(shì)斜率分析法

采用線性斜率法和F檢驗(yàn)計(jì)算黃土高原2000—2018年ET變化趨勢(shì)及其顯著性水平。線性斜率的計(jì)算公式如下[34]：

(2)

式中,Xi為第i年的ET值;n為總年數(shù)。當(dāng)Slope>0時(shí),ET呈增加趨勢(shì)；當(dāng)Slope<0時(shí),ET呈減小趨勢(shì)。當(dāng)ET變化趨勢(shì)的F檢驗(yàn)顯著性水平P<0.05時(shí),ET呈顯著變化。

1.3.3ET與氣象要素和EVI響應(yīng)關(guān)系構(gòu)建

借鑒前人研究成果與經(jīng)驗(yàn),本研究將構(gòu)建ET與降雨、溫度、日照時(shí)數(shù)、飽和水汽壓差和EVI之間定量響應(yīng)關(guān)系,且在構(gòu)建響應(yīng)關(guān)系過(guò)程中充分考慮植被與ET的非線性關(guān)系和氣象要素與植被之間的交互作用,具體公式如下[35]：

ET=β0+βEVI×EVI+βi×Xi+ε

(3)

Xi={EVI2,Cj,Cj×EVI},j=1,2,3,4

(4)

式中,ET為實(shí)際蒸散發(fā);βEVI為ET對(duì)EVI的敏感性系數(shù);Xi為其余變量子集;βi為所選變量子集Xi中各變量的敏感性系數(shù)。Xi中EVI2表征植被的非線性作用;Cj表征4種氣象要素(降雨、溫度、日照時(shí)數(shù)和飽和水汽壓差)的線性作用;Cj×EVI表征各氣象要素與植被之間的交互作用。

本研究利用多元逐步回歸方法[35]篩選出每個(gè)柵格的最佳變量子集Xi,從而在每一個(gè)柵格上建立“最優(yōu)”的ET與氣象因子和EVI定量響應(yīng)關(guān)系。

1.3.4豐、平、枯水年劃分

根據(jù)1960—2018年降水系列資料,運(yùn)用P-Ⅲ型頻率分布曲線分析逐個(gè)像元的降雨頻率。根據(jù)莫淑紅等[36]的分類方法將降雨頻率小于25%的年份劃分為豐水年；降雨頻率在25%—75%之間的為平水年；降雨頻率大于75%的為枯水年。P-Ⅲ型頻率分布曲線概率密度函數(shù)如下[37]：

(5)

式中,Γ(α)是伽瑪函數(shù);α、β、a0分別為P-Ⅲ型頻率分布的形狀參數(shù)、尺度參數(shù)和位置參數(shù),且α>0、β>0。根據(jù)劉煥彬[38]方法計(jì)算公式中所需各參數(shù)。

1.3.5不同情景下植被覆蓋度恢復(fù)潛力估算

根據(jù)水量平衡原理,蒸散發(fā)是由降雨量扣除產(chǎn)流量和陸地水儲(chǔ)量變化量之后的剩余部分。當(dāng)降雨補(bǔ)給能夠滿足生態(tài)系統(tǒng)耗水需求時(shí),植被群落演替具有可持續(xù)性。在維持降雨補(bǔ)給與生態(tài)系統(tǒng)ET的平衡關(guān)系條件下,ET的最大值為降雨量,即不發(fā)生任何產(chǎn)流和蓄水量的變化。具體公式如下[33]：

ET=P-R-ΔS

(6)

式中,ET為蒸散發(fā);P為降雨量;R為徑流量,ΔS為陸地水儲(chǔ)量變化量。

將ETmax=P帶入公式(3)反算EVI值,該值即為能夠維持降雨供給和生態(tài)系統(tǒng)耗水平衡的EVI閾值。

結(jié)合地表覆蓋產(chǎn)品數(shù)據(jù),利用像元二分法將EVI閾值轉(zhuǎn)化為林草植被覆蓋度,公式如下[39]：

(7)

式中,FVCi,j為第j類林草植被類型中的第i個(gè)像元的植被覆蓋度;EVIi,j為第j類林草植被類型中第i個(gè)像元的EVI值;EVImax,j為第j類林草植被類型的EVI年最大統(tǒng)計(jì)直方圖95%分位數(shù)對(duì)應(yīng)值;EVImin為裸地類型的EVI年最小值統(tǒng)計(jì)直方圖5%分位數(shù)對(duì)應(yīng)值。

植被覆蓋度恢復(fù)閾值與當(dāng)前植被覆蓋度現(xiàn)狀的差值即為植被覆蓋度恢復(fù)潛力,計(jì)算公式如下[22]：

FVCpot=FVCthr-FVCcur

(8)

式中,FVCpot為植被覆蓋度恢復(fù)潛力;FVCthr為植被覆蓋度恢復(fù)閾值;FVCcur為植被覆蓋度現(xiàn)狀。

2 結(jié)果與分析

2.1 黃土高原ET與EVI時(shí)空變化特征

圖3分別展示PML_V2產(chǎn)品和MOD16 A2 GF產(chǎn)品在4個(gè)通量站點(diǎn)和4個(gè)流域區(qū)間的驗(yàn)證精度,可看出站點(diǎn)尺度ETPML的RMSE、MRE和Nash系數(shù)分別為4.5 mm/8d、-14.7%和0.74,而區(qū)域尺度ETPML的RMSE和MRE分別為34.3 mm/a和1.57%。站點(diǎn)尺度ETMODIS的RMSE、MRE和Nash系數(shù)分別為6.17 mm/8d、168.4%和0.54,而區(qū)域尺度ETMODIS的RMSE和MRE分別為96.1 mm/a和-15.1%。基于站點(diǎn)和區(qū)域尺度驗(yàn)證結(jié)果均表明,PML_V2 ET產(chǎn)品在黃土高原地區(qū)的模擬精度較高,更能夠準(zhǔn)確捕捉ET時(shí)空動(dòng)態(tài)變化。黃土高原2000—2018年多年平均ET為445.36 mm,多年平均EVI為0.17,均呈現(xiàn)由東南向西北遞減的空間分布格局(圖4),空間分布格局基本一致。

圖3 PML_V2和MOD16 A2 GF ET產(chǎn)品驗(yàn)證結(jié)果Fig.3 Validation results of PML_V2 and MOD16 A2 GF evapotranspiration products Du2：多倫_草地,Duolun_grassland (D01)；Ha2：海北灌叢,Haibei Shrubland；HaM：西藏海北高山遺址,Haibei Alpine Tibet site；Sw2：四子王旗,Siziwang Grazed (SZWG)；RMSE：均方根誤差,root mean squared error；MRE：平均相對(duì)誤差,mean relative error；Nash：納什效率系數(shù),nash-sutcliffe efficiency coefficient；Tang：唐乃亥；Qing：青銅峽；Tou：頭道拐；Long：龍門；Hua：花園口

圖4 多年平均ET和EVI的空間分布Fig.4 Spatial distribution of perennial mean ET and EVI ET：蒸散發(fā),evapotranspiration； EVI：增強(qiáng)型植被指數(shù),enhanced vegetation index

從圖5中可看出,黃土高原ET呈增加趨勢(shì)的區(qū)域面積占比為92.18%,其中62.76%的區(qū)域呈顯著增加趨勢(shì),主要分布于山西省和青海省、陜西省北部等退耕還林草工程重點(diǎn)實(shí)施區(qū)域。黃土高原西北部和東南部平原地帶ET呈顯著下降趨勢(shì),占比7.82%,這主要與城鎮(zhèn)化等劇烈人類活動(dòng)有關(guān)。黃土高原EVI呈增加趨勢(shì)的面積占比為96.04%,其中85.56%的區(qū)域呈現(xiàn)顯著增加趨勢(shì),可見(jiàn)退耕還林還草效果顯著。內(nèi)蒙古和寧夏北部、陜西省南部、青海和山西省等零星地區(qū)的EVI呈現(xiàn)顯著下降趨勢(shì),面積占比僅為1.35%。

黃土高原林草植被ET年際增長(zhǎng)速率為5.08 mm/a,其中不同類型植被生態(tài)系統(tǒng)ET年際增長(zhǎng)速率排序?yàn)椋汗嗖?森林(圖6)。黃土高原林草植被EVI年際增長(zhǎng)速率為0.0026/a,其中不同植被類型的EVI年際增長(zhǎng)速率為：森林>灌草(圖6)。

圖5 2000—2018年ET與EVI時(shí)空分布Fig.5 Temporal and spatial distribution of ET and EVI, 2000—2018P<0.05：F檢驗(yàn)顯著性水平P<0.05時(shí),ET與EVI呈顯著變化

圖6 不同植被類型ET與EVI時(shí)序動(dòng)態(tài)變化Fig.6 Temporal dynamic changes of ET and EVI in different vegetation cover type zone

2.2 ET與EVI響應(yīng)關(guān)系空間分布特征

利用逐步多元回歸方法逐個(gè)像元構(gòu)建ET與氣象要素和EVI響應(yīng)關(guān)系。ET與EVI之間的響應(yīng)關(guān)系形式主要分為：1)僅EVI與ET相關(guān)；2)EVI和EVI2與ET相關(guān)；3)EVI和EVI×Cj與ET相關(guān)；4)EVI、EVI2和EVI×Cj均與ET相關(guān)。不同響應(yīng)關(guān)系形式的空間分布如圖7所示,統(tǒng)計(jì)結(jié)果表明ET與EVI呈線性關(guān)系的區(qū)域面積占比為60.8%,其中響應(yīng)關(guān)系受到EVI與氣象要素相互作用影響的區(qū)域占比為24.3%,主要分布于山西省北部和東部、陜西省中北部、寧夏南部、甘肅省東部和西部等部分區(qū)域。ET與EVI呈非線性關(guān)系的區(qū)域面積占比為39.2%,其中受到EVI與氣象要素相互作用影響的區(qū)域占比為14.9%,主要分布于甘肅省中西部、陜西省西北部、山西省北部以及內(nèi)蒙古東南部。多元回歸模型整體R2大小為0.66,其中,僅EVI與ET相關(guān)、EVI和EVI2與ET相關(guān)、EVI和EVI×Cj與ET相關(guān)、EVI、EVI2和EVI×Cj均與ET相關(guān)四種響應(yīng)關(guān)系模型的R2分別為0.61、0.67、0.67、0.74；73.35%的區(qū)域所構(gòu)建的響應(yīng)關(guān)系模型具有顯著性(圖7),模型效果較好。本文方法所構(gòu)建的響應(yīng)關(guān)系模型的平均均方根誤差(RMSE)為27.3 mm。為了減小ET產(chǎn)品驗(yàn)證的不確定性,站點(diǎn)尺度和區(qū)域尺度PML_V2產(chǎn)品的平均RMSE為35.15 mm,利用誤差傳播定律[40]計(jì)算多元回歸模型ET模擬值的平均RMSE為44.5 mm。

圖7 ET與EVI響應(yīng)關(guān)系及其顯著性Fig.7 Response relationship and significance of ET and EVIR2：可決系數(shù),coefficient of determination

圖8 不同植被類型區(qū)ET與EVI響應(yīng)關(guān)系統(tǒng)計(jì) Fig.8 The response relationship between ET and EVI in different vegetation cover type zone

不同植被ET與氣象因子和EVI的響應(yīng)關(guān)系具有差異性。圖8統(tǒng)計(jì)結(jié)果表明黃土高原林草植被ET與EVI呈線性關(guān)系的比例為62.23%,其中受到EVI與氣象因素相互作用的面積比例為27.19%。ET與EVI呈現(xiàn)非線性關(guān)系面積比例為37.77%,其中受到EVI與氣象要素相互作用的比例為14.54%。森林和灌草ET與EVI響應(yīng)關(guān)系以線性為主,面積占比為59.27%、63.20%；其中受到EVI與氣象要素相互作用影響的面積比分別為19.84%、29.60%。森林ET與EVI呈現(xiàn)非線性響應(yīng)關(guān)系的面積比(40.73%)高于灌草(36.80%),其中受到EVI與氣象要素相互作用影響的面積比分別為14.04%、14.71%。

2.3 豐、平、枯水年情景下林草植被覆蓋度恢復(fù)潛力

利用公式(3)—(8)對(duì)黃土高原地區(qū)豐、平、枯水年不同情景下的林草覆蓋度恢復(fù)閾值與恢復(fù)潛力進(jìn)行定量計(jì)算,結(jié)果如圖9所示；且根據(jù)2.2節(jié)中ET模擬值RMSE可計(jì)算得到林草覆蓋度恢復(fù)閾值與恢復(fù)潛力的平均RMSE為5.7%。

豐水年情景下,黃土高原林草植被覆蓋度平均恢復(fù)閾值為(71.5±37.3)%,平均恢復(fù)潛力為(9.4±30.4)%。其中,林草植被覆蓋度已超過(guò)恢復(fù)閾值的面積占比為33.35%,主要分布于甘肅省西北部、內(nèi)蒙古西南部、青海省南部、寧夏北部和山西省部分區(qū)域。而林草植被覆蓋度仍有恢復(fù)潛力的面積占比為66.65%,主要分布于內(nèi)蒙古和甘肅省東部、陜西省北部、山西省、青海省北部等區(qū)域,其中河龍區(qū)間中北部等區(qū)域林草植被覆蓋度恢復(fù)潛力可達(dá)到30%—60%。

平水年情景下,黃土高原林草植被覆蓋度平均恢復(fù)閾值為(55.6±35.9)%,平均恢復(fù)潛力為(-6.4±23.8)%。其中,林草植被覆蓋度已超過(guò)恢復(fù)閾值的面積占比為65.89%,主要分布于內(nèi)蒙古西南部、寧夏、甘肅、青海和陜西省南部、山西省部分區(qū)域。而林草植被覆蓋度仍有恢復(fù)潛力的面積占比為34.11%,主要分布于陜西省中北部、內(nèi)蒙古東部部分地區(qū)、青海省北部、山西省部分區(qū)域,其中陜西省中北部、青海省北部林草植被覆蓋度恢復(fù)潛力可達(dá)到15%—35%。

枯水年情景下,黃土高原林草植被覆蓋度平均恢復(fù)閾值僅為(22.4±26.0)%,平均恢復(fù)潛力為(-39.7±26.2)%。其中,林草植被覆蓋度已超過(guò)恢復(fù)閾值的面積占比為96.77%,在整個(gè)黃土高原連片分布。而林草植被覆蓋度仍有恢復(fù)潛力的面積占比僅為3.23%,主要分布于內(nèi)蒙古西北部和青海省北部等零星地區(qū),其林草植被覆蓋度恢復(fù)潛力僅為0—10%。

圖9 林草植被覆蓋度恢復(fù)閾值與恢復(fù)潛力Fig.9 Restoration threshold and restoration potential of shrubs-herbs-arbor vegetation coverFVCthr：植被覆蓋度恢復(fù)閾值Restoration threshold of fractional vegetation cover；FVCpot：植被覆蓋度恢復(fù)潛力Restoration potential of fractional vegetation cover

2.4 林草植被覆蓋度恢復(fù)潛力空間分異特征

降水量是限制黃土高原植被生長(zhǎng)的關(guān)鍵性因素,不同植被對(duì)干旱的敏感性不同。黃土高原灌草地多為中、低覆蓋,水分條件較差,植被生長(zhǎng)對(duì)干旱的響應(yīng)較為敏感。而東南部多為森林,降水豐富,較少受水分脅迫,植被生長(zhǎng)受干旱影響較小。

根據(jù)圖10結(jié)果可以看出,不同降雨量情景下的林草植被覆蓋度恢復(fù)潛力存在空間異質(zhì)性。豐水年情景下,灌草平均恢復(fù)潛力(10.2%)大于森林(6.3%)。林草植被覆蓋度平均恢復(fù)潛力隨干旱指數(shù)增加呈先增加后減小的變化特征,存在閾值效應(yīng)[41],峰值(25.3%)出現(xiàn)在干旱指數(shù)3—4區(qū)間。森林平均恢復(fù)潛力最大值(18.2%)出現(xiàn)在干旱指數(shù)2.5—3區(qū)間。平水年情景下,森林和灌草已分別超過(guò)平均恢復(fù)閾值6.5%和6.3%。林草植被覆蓋度超過(guò)恢復(fù)閾值的程度隨干旱指數(shù)增加而加劇?？菟昵榫跋?森林和灌草已分別超過(guò)其平均恢復(fù)閾值41.7%和39.0%。林草植被覆蓋度超過(guò)恢復(fù)閾值的程度隨干旱指數(shù)增加呈先加劇后減輕的變化特征,最大差值出現(xiàn)在干旱指數(shù)2.5—3區(qū)間。通過(guò)對(duì)比可知,不同歷史降雨情景對(duì)干旱指數(shù)2.5—4區(qū)間的林草植被恢復(fù)潛力影響最為劇烈。

圖10 不同植被類型植被恢復(fù)潛力隨干旱指數(shù)變化規(guī)律Fig.10 Changes of vegetation restoration potential with the increase of drought index in different vegetation coverage type zone

3 討論

退耕還林還草工程實(shí)施以來(lái),黃土高原植被恢復(fù)效果顯著,水土流失減弱[42—43],各項(xiàng)生態(tài)系統(tǒng)服務(wù)顯著提升[44]。但是,土壤水分已成為干旱/半干旱區(qū)植被演替的主要限制因素?；跇拥夭煌疃韧寥罎穸葘?shí)測(cè)數(shù)據(jù)的大量研究表明[45—48],人工建造喬木林對(duì)土壤水分的過(guò)度消耗超過(guò)了當(dāng)?shù)亟涤暄a(bǔ)給,已造成深層土壤出現(xiàn)干層現(xiàn)象。在區(qū)域尺度方面,大量前人研究結(jié)果[2,22,25]均表明黃土高原植被恢復(fù)已接近該區(qū)域水資源承載力閾值,這與本研究定性結(jié)論也基本保持一致。

前人研究基于“相似生境原則”[2,18,22],利用直方圖統(tǒng)計(jì)法提取不同生境斑塊的植被恢復(fù)閾值,但“相似生境”斑塊劃分誤差會(huì)嚴(yán)重影響提取結(jié)果,且該方法計(jì)算過(guò)程并未考慮降雨供給與生態(tài)系統(tǒng)ET的平衡關(guān)系,因此其提取結(jié)果可能是在消耗過(guò)量土壤水分條件下的植被最大覆蓋度。本研究從水量平衡原理的角度出發(fā),構(gòu)建了生態(tài)系統(tǒng)ET與氣象要素和植被指數(shù)之間的定量響應(yīng)關(guān)系,定量計(jì)算了能夠維持降雨供給與生態(tài)系統(tǒng)耗水平衡關(guān)系的林草植被恢復(fù)閾值,可用于林草植被可持續(xù)性評(píng)價(jià)。這與Feng[23]、Liang等[24]估算植被恢復(fù)閾值的研究思路類似,但前人研究尺度為縣域,且所用植被表征指標(biāo)GPP(gC m-2a-1)并不便于直接應(yīng)用于生態(tài)建設(shè)工程實(shí)踐,而本研究計(jì)算了每個(gè)500 m×500 m網(wǎng)格內(nèi)的林草植被恢復(fù)閾值,其結(jié)果可直接用于工程實(shí)踐。另外,前人研究主要利用近年遙感觀測(cè)數(shù)據(jù)提取植被閾值,僅代表了近十年氣候條件下林草植被動(dòng)態(tài)變化規(guī)律與恢復(fù)閾值,并未全面考慮未來(lái)不同氣候情景對(duì)植被閾值的影響。而本研究基于水量平衡原理和響應(yīng)關(guān)系模型,預(yù)測(cè)了未來(lái)豐、平、枯水年不同氣候條件下黃土高原林草植被覆蓋度恢復(fù)閾值及空間分布,可為未來(lái)黃土高原植被的可持續(xù)性評(píng)價(jià)與生態(tài)建設(shè)規(guī)劃提供重要參考。

生態(tài)系統(tǒng)ET的反演精度直接決定著本方法估算植被恢復(fù)閾值的準(zhǔn)確性,由于可收集共享數(shù)據(jù)有限,本研究驗(yàn)證精度存在一定不確定性。未來(lái)可收集更多實(shí)地測(cè)量土壤濕度和ET通量數(shù)據(jù)對(duì)多種ET反演模型進(jìn)行驗(yàn)證與改進(jìn)。生態(tài)系統(tǒng)ET與多種因素實(shí)際存在復(fù)雜函數(shù)響應(yīng)關(guān)系,本研究利用簡(jiǎn)單數(shù)學(xué)關(guān)系構(gòu)建ET與氣象因素與植被之間的響應(yīng)關(guān)系,存在一定不確定性,但本研究在構(gòu)建關(guān)系過(guò)程中充分考慮了ET對(duì)植被的非線性響應(yīng)和氣象因素與植被之間的相互影響,可在一定程度上減小該不確定性,未來(lái)可結(jié)合統(tǒng)計(jì)方法、水文模型、陸面模式等手段綜合研判植被與ET的響應(yīng)關(guān)系。另外,在計(jì)算響應(yīng)關(guān)系模型的不確定性過(guò)程中,并沒(méi)有詳細(xì)考慮各氣象要素觀測(cè)誤差對(duì)反演結(jié)果不確定性的影響,不確定性定量評(píng)價(jià)仍有待加強(qiáng)。本研究在計(jì)算ET最大值時(shí)并未考慮經(jīng)濟(jì)社會(huì)用水情況,僅表示天然水循環(huán)情景下生態(tài)系統(tǒng)可利用的最大水量。

未來(lái)生態(tài)建設(shè)工程規(guī)劃中,建議對(duì)于已超過(guò)植被覆蓋度恢復(fù)閾值的區(qū)域,應(yīng)采取有效調(diào)控措施,如人工建造植被多以灌草為主或輕度放牧等。對(duì)于仍有恢復(fù)潛力的區(qū)域,也需要根據(jù)氣候帶和立地條件開(kāi)展差異化生態(tài)建設(shè),做到“分區(qū)分類、因地制宜”,改善優(yōu)化植被結(jié)構(gòu)與格局。采用“植被適度恢復(fù)以實(shí)現(xiàn)高質(zhì)量的生態(tài)系統(tǒng)服務(wù)”的思路對(duì)植被恢復(fù)的進(jìn)程和目標(biāo)進(jìn)行規(guī)劃,降低水資源的過(guò)度消耗,提升黃土高原區(qū)生態(tài)系統(tǒng)服務(wù)功能,實(shí)現(xiàn)黃河流域高質(zhì)量發(fā)展[49]。

4 結(jié)論

本研究基于天然水循環(huán)水量平衡原理,利用多源遙感產(chǎn)品和實(shí)測(cè)地面資料,構(gòu)建了ET與氣象要素和EVI之間的定量響應(yīng)關(guān)系,定量估算了不同降雨情景下黃土高原林草植被覆蓋度恢復(fù)閾值與恢復(fù)潛力。主要研究結(jié)論如下：

(1)PML_V2產(chǎn)品基于通量站點(diǎn)觀測(cè)值與水量平衡公式驗(yàn)證的誤差RMSE分別為4.5 mm/8d、34.3 mm/a,PML_V2產(chǎn)品更能夠準(zhǔn)確捕捉ET時(shí)空動(dòng)態(tài)變化。黃土高原2000—2018年多年平均ET與EVI為445.36 mm、0.17；林草植被ET與EVI增長(zhǎng)速率為5.08 mm/a、0.0026/a。

(2)建立的ET與EVI和氣象要素的響應(yīng)關(guān)系模型的估算誤差RMSE為44.5 mm。不同植被ET與氣象因子和EVI的響應(yīng)關(guān)系具有明顯空間異質(zhì)性。

(3)豐、平、枯水年情景下,黃土高原林草植被覆蓋度平均恢復(fù)閾值分別為(71.5±37.3)%、(55.6±35.9)%和(22.4±26.0)%。與林草植被覆蓋度現(xiàn)狀(62.1%)相比,豐水年情景下林草植被覆蓋度仍有(9.4±30.4)%的增長(zhǎng)空間,平水年和枯水年情景下林草蓋度已超過(guò)恢復(fù)閾值(6.4±23.8)%和(39.7±26.2)%。

(4)黃土高原林草植被覆蓋度恢復(fù)潛力與區(qū)域干旱指數(shù)密切相關(guān)。豐水年情景下,林草植被覆蓋度平均恢復(fù)潛力隨干旱指數(shù)增加呈先增加后減小的變化特征；平水年情景下,林草植被覆蓋度超過(guò)恢復(fù)閾值的程度隨干旱指數(shù)增加而加??；枯水年情景下,林草植被覆蓋度超過(guò)恢復(fù)閾值的程度隨干旱指數(shù)增加呈先加劇后減輕的變化特征。