王 凡,吳一平,*,李匯文,張廣創(chuàng),連炎清,王文科

1 西安交通大學(xué),人居環(huán)境與建筑工程學(xué)院,西安 710049

2 中國(guó)科學(xué)院地球環(huán)境研究所,西安 710061

3 長(zhǎng)安大學(xué),水利與環(huán)境學(xué)院,西安 710054

生態(tài)系統(tǒng)彈性力描述了生態(tài)系統(tǒng)受到環(huán)境沖擊之后保持其原有結(jié)構(gòu)和功能的能力,它對(duì)于區(qū)域生態(tài)系統(tǒng)安全評(píng)估具有重要意義[1-2]。生態(tài)系統(tǒng)彈性力被定義為在維持系統(tǒng)結(jié)構(gòu)、功能反饋等不變的前提下,通過(guò)調(diào)整系統(tǒng)狀態(tài)變量和驅(qū)動(dòng)變量等參數(shù),系統(tǒng)能吸收的擾動(dòng)量,用以表示生態(tài)系統(tǒng)在面對(duì)外部干擾時(shí)保持生態(tài)系統(tǒng)功能的能力[3- 5]。 Holling等[3]將這種能力量化為一個(gè)生態(tài)系統(tǒng)能夠承受的變化幅度。通常在受到較小的環(huán)境干擾時(shí),生態(tài)系統(tǒng)更容易恢復(fù),當(dāng)生態(tài)系統(tǒng)遭受到嚴(yán)重的擾動(dòng)時(shí)會(huì)趨向遠(yuǎn)離其均衡穩(wěn)態(tài),進(jìn)而轉(zhuǎn)變?yōu)榱硪环N狀態(tài),最終導(dǎo)致系統(tǒng)受損而退化[5]。在沒(méi)有人為干擾的自然生態(tài)系統(tǒng)中,生態(tài)彈性力也可用來(lái)反映自然生態(tài)系統(tǒng)承載力,是區(qū)域生態(tài)承載力的研究基礎(chǔ)和描述脆弱生態(tài)區(qū)狀態(tài)的綜合變量[6- 8]。

量化生態(tài)系統(tǒng)彈性力是預(yù)測(cè)生態(tài)系統(tǒng)應(yīng)對(duì)全球變化的關(guān)鍵[9]。國(guó)內(nèi)外學(xué)者采用不同的方法對(duì)生態(tài)系統(tǒng)彈性力展開(kāi)了大量的研究,Ashutosh等[10]采用水分效率模型來(lái)評(píng)估印度的生態(tài)系統(tǒng)彈性力,發(fā)現(xiàn)生態(tài)系統(tǒng)彈性力可能受土地覆蓋、氣候、水文要素的影響；Stephen等[11]利用全局敏感性和不確定性分析法評(píng)估了生態(tài)系統(tǒng)彈性力,為生態(tài)系統(tǒng)彈性力定量評(píng)價(jià)提供了有效方案；Luciano等[12]基于氣候生態(tài)位理論提出了測(cè)量生態(tài)系統(tǒng)彈性力的方法,發(fā)現(xiàn)持續(xù)的氣候變化可通過(guò)促進(jìn)森林生物多樣性的侵蝕并導(dǎo)致植被覆蓋度降低而加速生態(tài)系統(tǒng)彈性力的喪失；Wu等[9]用葉面積指數(shù)來(lái)評(píng)估陸地生態(tài)系統(tǒng)的彈性力,并將其動(dòng)態(tài)變化歸因于氣候與環(huán)境因素；生態(tài)彈性力模型[13-14]結(jié)合主成分分析、模糊數(shù)學(xué)的綜合指標(biāo)評(píng)估[15-16]方案目前也被廣泛應(yīng)用。對(duì)于生態(tài)系統(tǒng)彈性力的影響因素,已有研究表明生態(tài)系統(tǒng)彈性力與植被生產(chǎn)力、生物多樣性呈正相關(guān)關(guān)系[17]。此外,氣候變化和人類(lèi)活動(dòng)的累計(jì)效應(yīng)影響著生態(tài)系統(tǒng)的結(jié)構(gòu)和動(dòng)態(tài),進(jìn)而顯著影響生態(tài)系統(tǒng)的恢復(fù)力[18]。雖然大量學(xué)者對(duì)生態(tài)系統(tǒng)彈性力的影響因素進(jìn)行了探究,但目前仍缺少一致的觀點(diǎn)和表征生態(tài)系統(tǒng)彈性力的指標(biāo)體系[19]。總的來(lái)說(shuō),現(xiàn)有研究方法多應(yīng)用于區(qū)域尺度上,且考慮的因素還不夠全面,針對(duì)工礦活動(dòng)下的流域生態(tài)系統(tǒng)彈性力的研究還鮮有報(bào)道。

陜北榆神礦區(qū)作為我國(guó)重要的煤炭生產(chǎn)基地之一[20],煤炭資源豐富,人類(lèi)活動(dòng)劇烈,持續(xù)的工礦活動(dòng)引發(fā)了礦區(qū)局部水土流失[21]、植被退化[22]以及地下水位下降[23-24]等一系列生態(tài)環(huán)境問(wèn)題,導(dǎo)致該地區(qū)生態(tài)環(huán)境更加脆弱[25]。自1999年起,該地區(qū)退耕還林還草、天然林保護(hù)和封山禁牧等生態(tài)修復(fù)工程的實(shí)施雖然整體上改善了區(qū)域的生態(tài)環(huán)境[26-27],但工礦活動(dòng)的直接影響和間接影響區(qū)內(nèi)的生態(tài)環(huán)境變化特征尚不清晰,同時(shí)對(duì)近年來(lái)工礦影響區(qū)的生態(tài)系統(tǒng)穩(wěn)定性狀況缺乏系統(tǒng)的了解。煤炭開(kāi)發(fā)區(qū)生態(tài)環(huán)境復(fù)雜,影響其生態(tài)系統(tǒng)穩(wěn)定性因素較多,如何找準(zhǔn)生態(tài)影響因子,客觀合理地評(píng)價(jià)礦區(qū)生態(tài)現(xiàn)狀及其恢復(fù)力顯得尤為重要。本文以榆神礦區(qū)所在流域?yàn)檠芯繉?duì)象,在了解研究區(qū)整體土地利用和煤炭礦區(qū)分布格局的基礎(chǔ)上,依據(jù)研究區(qū)生態(tài)環(huán)境背景,選取水文、土壤、植被3個(gè)方面5個(gè)指標(biāo),運(yùn)用綜合指標(biāo)法和主成分分析法對(duì)研究區(qū)生態(tài)彈性力進(jìn)行評(píng)估,探討研究區(qū)2009—2018年生態(tài)系統(tǒng)彈性力時(shí)空演變特征,分析煤炭開(kāi)發(fā)對(duì)研究區(qū)生態(tài)系統(tǒng)彈性力的影響。本研究通過(guò)對(duì)礦區(qū)生態(tài)系統(tǒng)彈性力的評(píng)估來(lái)量化礦區(qū)生態(tài)系統(tǒng)穩(wěn)定性和可調(diào)節(jié)能力,以期為優(yōu)化礦區(qū)生態(tài)格局和進(jìn)一步完善礦區(qū)可持續(xù)發(fā)展策略提供科學(xué)依據(jù)。

1 研究區(qū)概況

榆神礦區(qū)位于陜晉蒙接壤地帶,是陜北煤炭基地的主要礦區(qū)之一,礦區(qū)面積約5265 km2,2013年,榆神礦區(qū)煤炭產(chǎn)量已達(dá)到2.78億t[24]。礦區(qū)總體被劃分為四個(gè)規(guī)劃區(qū),根據(jù)不同的資源勘探程度,開(kāi)啟了前三期規(guī)劃區(qū)的開(kāi)發(fā),第四期規(guī)劃區(qū)正在準(zhǔn)備階段[23]。本文以榆神礦區(qū)所在的3個(gè)流域?yàn)檠芯繀^(qū)(圖1),包含禿尾河流域(占流域總面積44.7%)、榆溪河流域(占45%)和佳蘆河流域(占10.3%),研究區(qū)面積約為1.07萬(wàn)km2,處于東經(jīng)109°—110.52°,北緯37.98°—39.34°之間。研究區(qū)東南部為黃土丘陵溝壑區(qū),西北部為風(fēng)沙草灘區(qū),地勢(shì)由東南向西北遞增,主要植被有沙柳、檸條、沙蒿、錦雞兒、苜蓿等天然植被和柳樹(shù)、楊樹(shù)、棗樹(shù)等人工植被,主要土壤類(lèi)型為風(fēng)沙土、黃綿土和沖積土[23]。該區(qū)地處西北內(nèi)陸干旱區(qū),降水年際變化大,多年年均降水約400 mm,水資源量短缺,多年年均潛在蒸發(fā)量約1700 mm；區(qū)域內(nèi)年均氣溫在6—14℃ 之間,多年均溫約10℃[28]。由于研究區(qū)氣候干燥、植被稀疏、地形坡度大等自然地理特征使其生態(tài)環(huán)境脆弱,煤炭資源的開(kāi)采,導(dǎo)致研究區(qū)地表植被破壞、水土流失、土地荒漠化程度加劇、地下水位下降、地面塌陷等生態(tài)環(huán)境問(wèn)題的產(chǎn)生和發(fā)展[4],嚴(yán)重影響該區(qū)域生態(tài)系統(tǒng)的穩(wěn)定性[29]。

圖1 研究區(qū)高程及地理位置Fig.1 The elevation and geographical location of the study area

2 數(shù)據(jù)與方法

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

本研究選取2009、2012、2015和2018年4期Landsat TM 30 m分辨率的遙感影像數(shù)據(jù)(行列號(hào)127/33、34),數(shù)據(jù)都集中于6月份,且云量小于5%。對(duì)遙感影像進(jìn)行輻射定標(biāo)、大氣校正、波段融合、影像鑲嵌和裁剪等預(yù)處理工作,依據(jù)全國(guó)生態(tài)系統(tǒng)分類(lèi)體系[30]和研究區(qū)土地利用類(lèi)型的實(shí)際情況,將研究區(qū)的土地利用類(lèi)型劃分為耕地、草地、林地、建設(shè)用地、水域、未利用地和工礦用地7種。通過(guò)人機(jī)交互式解譯獲得研究區(qū)土地利用類(lèi)型空間分布圖,最后利用目視解譯結(jié)果與整體分類(lèi)結(jié)果建立混淆矩陣,計(jì)算Kappa系數(shù)來(lái)驗(yàn)證分類(lèi)精度,四期土地利用分類(lèi)精度達(dá)80%。

為了評(píng)價(jià)研究區(qū)的生態(tài)系統(tǒng)彈性力,本文收集了研究區(qū)2009—2018年生態(tài)環(huán)境因子時(shí)空數(shù)據(jù)。其中地表徑流和地下徑流數(shù)據(jù)基于SWAT(Soil and Water Assessment Tool)模型在水文響應(yīng)單元上的模擬結(jié)果,利用地理信息系統(tǒng)進(jìn)行數(shù)據(jù)格式轉(zhuǎn)換,得到2009—2018年分辨率1 km的地表、地下徑流空間分布圖。植被覆蓋度是SPOT5-VGT2衛(wèi)星采集的旬值數(shù)據(jù),空間分辨率為1 km； 植被生產(chǎn)力(NPP)數(shù)據(jù)來(lái)源于美國(guó)航空航天局(NASA)(https://earthdata.nasa.gov/)提供的MOD17A3數(shù)據(jù)產(chǎn)品,空間分辨率為1 km。土壤保持量基于通用土壤流失方程(revised universal soil loss equation,RUSLE)計(jì)算,利用地形、土壤數(shù)據(jù)并結(jié)合月降雨量、植被覆蓋度計(jì)算得到,其中地形數(shù)據(jù)采用地理空間數(shù)據(jù)云SRTM 90 m 數(shù)字高程模型(DEM)產(chǎn)品(http://www.gscloud.cn/),土壤數(shù)據(jù)來(lái)源于寒區(qū)旱區(qū)科學(xué)數(shù)據(jù)中心(http://westdc.westgis.ac.cn)提供的世界土壤數(shù)據(jù)庫(kù),空間分辨率為1 km,所需數(shù)據(jù)包括土壤機(jī)械組成和土壤有機(jī)碳數(shù)據(jù)。

2.2 評(píng)價(jià)方法

2.2.1生態(tài)彈性力指標(biāo)體系構(gòu)建

煤炭礦區(qū)生態(tài)彈性力反映了生態(tài)系統(tǒng)對(duì)礦區(qū)社會(huì)經(jīng)濟(jì)活動(dòng)壓力的恢復(fù)能力[31],根據(jù)礦區(qū)實(shí)際生態(tài)環(huán)境狀況,選取水文、土壤和植被三個(gè)方面的指標(biāo)(表1),通過(guò)構(gòu)建生態(tài)彈性力評(píng)價(jià)指標(biāo)體系,結(jié)合空間主成分分析法量化研究區(qū)的生態(tài)彈性力。根據(jù)研究目標(biāo),按照目標(biāo)分層法建立指標(biāo)體系,其中生態(tài)彈性力作為目標(biāo)層,水文、土壤和植被三個(gè)要素作為準(zhǔn)則層,基于準(zhǔn)則層選取相關(guān)指標(biāo)作為指標(biāo)層。具體步驟如下：

表1 生態(tài)彈性力評(píng)價(jià)指標(biāo)體系Table 1 Ecological resilience evaluation indicator system

①為解決各個(gè)指標(biāo)量綱不同的問(wèn)題,首先需要對(duì)原始數(shù)據(jù)進(jìn)行歸一化處理；

(1)

式中,Yi為指標(biāo)標(biāo)準(zhǔn)值；Xi為第i項(xiàng)指標(biāo)的值；max(Xi)與min(Xi)為第i項(xiàng)指標(biāo)的最大值和最小值,i的范圍為1—n,n表示柵格個(gè)數(shù)。

②根據(jù)歸一化數(shù)據(jù)判定指標(biāo)之間的相關(guān)性,建立相關(guān)系數(shù)矩陣,指標(biāo)之間的相關(guān)系數(shù)越大表明有必要對(duì)數(shù)據(jù)進(jìn)行主成分分析；

③計(jì)算相關(guān)系數(shù)矩陣的特征值、主成分貢獻(xiàn)率和累計(jì)貢獻(xiàn)率,確定各主成分個(gè)數(shù),并計(jì)算主成分權(quán)重；

(2)

式中,Gj為主成分j的權(quán)重,λ為特征值的方差貢獻(xiàn)率,m為主成分個(gè)數(shù)(一般特征值λm的累計(jì)貢獻(xiàn)率達(dá)80%—95%)；

④計(jì)算生態(tài)彈性力指數(shù)；

(3)

式中,Z為生態(tài)彈性力指數(shù),Fj為主成分j的綜合得分。

2.2.2水文要素計(jì)算

SWAT模型是對(duì)流域尺度水文過(guò)程進(jìn)行連續(xù)時(shí)段模擬的物理模型,能夠模擬多種土地利用、土壤和管理?xiàng)l件下的水文過(guò)程[32]。本文采用該模型對(duì)研究區(qū)的降雨徑流進(jìn)行模擬,使用逐月徑流對(duì)模型進(jìn)行率定和驗(yàn)證。選用決定系數(shù)R2、百分?jǐn)?shù)偏差PBIAS評(píng)價(jià) SWAT 模型在率定期和驗(yàn)證期的表現(xiàn),具體的計(jì)算公式參照已有文獻(xiàn)[33]。 基于 Moriasi 等[34]人 2007 年提出的標(biāo)準(zhǔn),兩個(gè)水文站點(diǎn)的模擬結(jié)果在率定期和驗(yàn)證期憑借R2>0.50, |PBIAS|<25%被判定為良好(表2)。

表2 高家川和申家灣水文站點(diǎn)徑流模擬結(jié)果評(píng)價(jià)Table 2 Evaluation of runoff simulation results at Gaojiachuan and Shenjiawan hydrological stations

2.2.3土壤保持量計(jì)算

土壤保持是指生態(tài)系統(tǒng)防止水土流失的侵蝕調(diào)控及對(duì)泥沙的儲(chǔ)積保持能力,可以用一個(gè)地區(qū)的土壤保持量來(lái)衡量土壤保持能力[35]。目前往往采用修正后的通用土壤流失方程評(píng)估土壤保持能力[36],認(rèn)為潛在土壤侵蝕量與實(shí)際侵蝕量的差值是在植被和人為管理作用下的土壤保持量,其計(jì)算公式為：

SR=APn-ARn

(4)

APn=R×K×LS

(5)

ARn=R×K×LS×C×P

(6)

式中,SR為土壤保持量(t hm-2a-1)；APn為潛在侵蝕量；ARn為在植被覆蓋和水土保持措施下發(fā)生的實(shí)際侵蝕量；R為降雨侵蝕力因子(MJ mm h-1hm-2a-1)；K為土壤可蝕性因子(t hm2h hm-2MJ-1mm-1)；LS為坡長(zhǎng)坡度因子；C為植被覆蓋與管理因子(0—1)；P為水土保持措施因子(0—1)。

降雨侵蝕力因子(R)可以反映降雨引起土壤侵蝕的潛在能力[37]。本研究采用年降雨侵蝕力模型計(jì)算降雨侵蝕力因子[38],其計(jì)算公式為：

(7)

式中,Prei為月降雨量(mm),Pre為年總降雨量(mm)。

土壤可侵蝕性因子(K)反映了土壤的抗侵蝕能力,其大小與土壤結(jié)構(gòu)和有機(jī)碳含量密切相關(guān),具體計(jì)算公式[37]如下：

(8)

式中,sa為沙粒含量(%),si為粉粒含量(%),cl為粘粒含量(%),c為有機(jī)碳含量(%),sn=1-sa/100。

植被覆蓋因子(C)反映的是植被覆蓋及相關(guān)管理措施對(duì)土壤的綜合作用,主要受植被覆蓋度和土地利用類(lèi)型等因素的影響,C值越小,表明地面植被覆蓋度越高。其計(jì)算公式[39]為：

(9)

式中,fc為植被覆蓋度。

坡長(zhǎng)坡度因子(LS)是影響坡面土壤侵蝕的重要因素,其計(jì)算公式[39]如下：

將上面的基于單處理器的算法思想拓展為支持多處理器的算法，即可得到求解原問(wèn)題的優(yōu)先關(guān)系調(diào)度算法(PR)：

(10)

(11)

(12)

式中,λ為坡長(zhǎng)(m),θ為坡度,n為坡長(zhǎng)指數(shù),S為坡度指數(shù)。

水土保持措施因子(P)是特定的水土保持措施下土壤流失量與未實(shí)施水土保持措施時(shí)的土壤流失量比值,研究發(fā)現(xiàn)水土保持措施因子與坡度存在線性關(guān)系[40],其計(jì)算公式為：

P=0.2+0.03θ

(13)

3 結(jié)果與分析

3.1 土地利用分布特征

研究區(qū)土地利用類(lèi)型以林地、草地、未利用地為主(圖2、表3),耕地、林地和草地三者面積之和約占整個(gè)研究區(qū)的70%,未利用地面積約占研究區(qū)的30%。研究期間,未利用地面積不斷減少,減少了780.2 km2,主要表現(xiàn)在研究區(qū)的中部和西北部；其他土地利用類(lèi)型都表現(xiàn)出了不同程度的增加,其中草地、林地增加幅度最大,分別增加了351.2、 217.7 km2；建設(shè)用地增加了59.1%,主要體現(xiàn)在榆陽(yáng)區(qū)和神木縣市區(qū)。對(duì)比2009年和2018年土地利用變化,土地利用轉(zhuǎn)移主要發(fā)生耕地、林地、草地和未利用地之間,同時(shí),89 km2的未利用地和112 km2的草地轉(zhuǎn)換成建設(shè)用地,3.4 km2的未利用地和4.4 km2的草地轉(zhuǎn)換成工礦用地?？傮w來(lái)看,隨著生態(tài)恢復(fù)工程的實(shí)施和社會(huì)經(jīng)濟(jì)的快速發(fā)展,研究區(qū)土地利用空間格局發(fā)生了很大的變化,生態(tài)環(huán)境明顯改善,城市發(fā)生了進(jìn)一步的擴(kuò)張[41]。

表3 研究區(qū)2009—2018年土地利用面積轉(zhuǎn)移矩陣/km2Table 3 Land use area transfer matrix of study area from 2009 to 2018

研究期間,工礦用地面積明顯增加(圖3),從8.4 km2增加到14.2 km2。為了解煤炭開(kāi)發(fā)對(duì)礦區(qū)生態(tài)環(huán)境的影響,以工礦用地為中心分別設(shè)置半徑1 km和3 km的緩沖區(qū),1 km緩沖區(qū)作為煤炭開(kāi)發(fā)的直接影響區(qū),1—3 km緩沖區(qū)作為間接影響區(qū)。2009年、2012年、2015年、2018年直接影響區(qū)分別占整個(gè)區(qū)域面積的14.8%、14.7%、19.6%、24%,間接影響區(qū)面積也逐漸增加,約占研究區(qū)面積的45%。

圖3 煤炭開(kāi)發(fā)區(qū)及其影響區(qū)域分布Fig.3 Distribution of coal development zones and influenced areas

3.2 生態(tài)彈性力評(píng)價(jià)

3.2.1生態(tài)彈性力評(píng)價(jià)指標(biāo)分析

經(jīng)過(guò)分析,得到了研究區(qū)2009—2018年生態(tài)彈性力評(píng)價(jià)指標(biāo)時(shí)間趨勢(shì)(圖4)和空間分布(圖5)。

圖4 生態(tài)彈性力評(píng)價(jià)指標(biāo)時(shí)間變化Fig.4 Time variation of evaluation indicators of ecological resilience

圖5 生態(tài)彈性力評(píng)價(jià)指標(biāo)空間分布Fig.5 Spatial distribution of evaluation indicators of ecological resilience

2009—2018年,研究區(qū)地表、地下徑流量均呈現(xiàn)顯著增加的趨勢(shì),增加速率分別為6.27、8.36 mm/a。對(duì)于長(zhǎng)期來(lái)說(shuō),隨著植被的增加研究區(qū)徑流量呈下降趨勢(shì),但是近年來(lái),由于降雨量的顯著增加導(dǎo)致徑流量增加,這與實(shí)測(cè)資料分析結(jié)果一致[42]；此外,在植被增加的同時(shí)研究區(qū)建設(shè)用地面積也增加,通過(guò)統(tǒng)計(jì)不同土地利用對(duì)應(yīng)的地表徑流量(表4),發(fā)現(xiàn)林、草地僅是建設(shè)用地地表徑流量的50%左右。從徑流量在水文響應(yīng)單元尺度上的空間分布特征來(lái)看,研究區(qū)地表徑流量主要處于6—30 mm之間,在研究區(qū)北部紅堿淖湖周?chē)乇韽搅髁孔钚?<6.5 mm),總體上呈現(xiàn)北部低西南部高的空間分布特征。研究區(qū)地下徑流量較地表徑流量大,主要由于研究區(qū)大面積為風(fēng)沙區(qū),地表下滲量大；地下徑流量在空間上表現(xiàn)出東北部高、西部偏低的特征,東南部佳蘆河流域地下徑流量最小,主要處于29—70 mm之間。

表4 不同土地利用地表徑流量/mmTable 4 Surface runoff of different land uses

研究期間,植被覆蓋情況顯著改善,增加速率為0.63%/a；植被生產(chǎn)力變化趨勢(shì)基本與植被覆蓋度一致,其增加速率為6.85 gC/m2?？臻g上,研究區(qū)植被覆蓋度處在0—0.29之間,呈現(xiàn)出東南部明顯高于西北部的特征。研究區(qū)植被生產(chǎn)力多年均值處在8—355 gC/m2之間,東南部和北部林草覆蓋較多的區(qū)域植被生產(chǎn)力較高,且沿著河道的耕地區(qū)植被生產(chǎn)力最高(> 200 gC/m2),由于研究區(qū)的西北部主要為未利用地, 植被生產(chǎn)力也處于較低的水平,主要集中在28—140 gC/m2之間。

3.2.2生態(tài)系統(tǒng)彈性力分析

根據(jù)生態(tài)彈性力評(píng)價(jià)指標(biāo)間存在相關(guān)系數(shù)大于0.5的結(jié)果得知,指標(biāo)之間存在一定的相關(guān)性,適合采用主成分分析法確定權(quán)重。依據(jù)選取原則確定了兩個(gè)主成分,其權(quán)重如表5所示。根據(jù)主成分荷載,主成分一主要支配指標(biāo)為地下徑流量、植被覆蓋度,主成分二主要支配指標(biāo)為地表徑流量、植被生產(chǎn)力。

表5 空間主成分分析的各主成分的特征值、比值及權(quán)重Table 5 Percent and accumulative eigenvalues

研究區(qū)2009—2018年生態(tài)彈性力的時(shí)間變化趨勢(shì)如圖6所示,生態(tài)彈性力指數(shù)總體呈現(xiàn)增加趨勢(shì),從2009年的0.53增加到2018年的0.60,但中間過(guò)程波動(dòng)幅度較大,主要受徑流量的影響,這是由于近年來(lái)陜北地區(qū)降雨量顯著增加且年際變化大。2009—2017年,生態(tài)彈性力呈現(xiàn)“W”型變化特征,其中2009—2011年,彈性力由0.53增加至0.56,基本呈現(xiàn)平穩(wěn)趨勢(shì),2012年后彈性力呈現(xiàn)先增加后減小的變化趨勢(shì),2015年降低至研究期間的最低值,2015—2017年彈性力指數(shù)顯著增加,2017年生態(tài)彈性力達(dá)到最大值0.77,2017—2018年彈性力又表現(xiàn)出明顯的減小?？傮w上,研究區(qū)生態(tài)系統(tǒng)抗干擾和自我調(diào)節(jié)能力不斷提高但不穩(wěn)定,與鄭欣等[16]在鄂爾多斯的研究結(jié)果相似。

從各個(gè)主成分來(lái)看,第一主成分與植被覆蓋度和地下徑流量都存在較高的正相關(guān)關(guān)系。由于第一主成分的貢獻(xiàn)率較大,其變化趨勢(shì)與生態(tài)彈性力基本一致。因此,第一主成分對(duì)生態(tài)彈性力起決定性作用,而植被覆蓋度和地下徑流量主要決定了研究區(qū)生態(tài)彈性力的變化趨勢(shì)。從圖6中可以看出,主成分二減緩了生態(tài)彈性力的波動(dòng)幅度,而主成分二的變化主要取決于地表徑流量和植被生產(chǎn)力。因此,地表徑流量和植被生產(chǎn)力對(duì)減緩研究區(qū)生態(tài)系統(tǒng)穩(wěn)定性的變化起到正向作用。分析表明,本研究區(qū)生態(tài)彈性力對(duì)水文、植被條件變化比較敏感,主要受植被覆蓋度和地下徑流量的影響。

由于考慮了多個(gè)空間指標(biāo),使得生態(tài)彈性力空間分布差異較大(圖6),總體呈現(xiàn)東南部偏高西北部低的特征,表明研究區(qū)東南部的生態(tài)系統(tǒng)較西北部穩(wěn)定。由于研究區(qū)東南部植被覆蓋度大,又屬于流域的下游,水資源相對(duì)豐富,生態(tài)彈性力較高,而西北部受氣候和土壤條件影響,裸地分布較廣,植被覆蓋度和植被生產(chǎn)力都較低,因此生態(tài)彈性力較低。局部來(lái)看,在榆陽(yáng)區(qū)和神木市縣區(qū)周?chē)把鼐€生態(tài)彈性力較低,主要?dú)w因于城市建成區(qū)人類(lèi)活動(dòng)頻繁,生態(tài)環(huán)境較差。

圖6 2009—2018年生態(tài)彈性力時(shí)空分布特征Fig.6 Spatial and temporal distribution characteristics of ecological resilience from 2009 to 2018

3.2.3榆神礦區(qū)生態(tài)彈性力變化

統(tǒng)計(jì)不同土地利用對(duì)應(yīng)的生態(tài)彈性力(圖7),得到不同土地利用生態(tài)彈性力排序：林地>草地>耕地>建設(shè)用地>工礦用地,由于林地的植被生產(chǎn)力、植被覆蓋度高且水土保持能力強(qiáng),因此生態(tài)系統(tǒng)穩(wěn)定性較高。2009—2015年,生態(tài)彈性力呈現(xiàn)下降趨勢(shì),林草地的彈性力下降最為明顯。雖然近年來(lái)隨著生態(tài)修復(fù)措施的實(shí)施,研究區(qū)植被覆蓋度和植被生產(chǎn)力增加,但生態(tài)彈性力受水文要素的影響很大,由于降雨量的降低使得徑流量減少,導(dǎo)致生態(tài)彈性力降低。

圖7 不同土地利用生態(tài)系統(tǒng)彈性力Fig.7 Ecological resilience of different land use

煤炭礦區(qū)直接影響區(qū)是煤炭開(kāi)發(fā)的內(nèi)圍區(qū)域,生態(tài)系統(tǒng)彈性力較間接影響區(qū)低,平均低9.23%。其中2009年兩個(gè)影響區(qū)生態(tài)系統(tǒng)彈性力差距最大,間接影響區(qū)比直接影響區(qū)彈性力高0.07,2015年差距最小,整個(gè)區(qū)域的彈性力介于直接影響區(qū)和間接影響區(qū)之間(圖8)?？傮w而言,由于耕地的植被覆蓋度和植被生產(chǎn)力較高,生態(tài)系統(tǒng)彈性力也較高,而直接影響區(qū)中耕地面積比例較大,導(dǎo)致直接影響區(qū)的生態(tài)彈性力偏高。由于間接影響區(qū)中有很大比例是裸地,裸地對(duì)應(yīng)的生態(tài)彈性力低。因此,間接影響區(qū)相比直接影響區(qū)的彈性力雖有增加,但差距不明顯。

圖8 煤炭開(kāi)發(fā)影響區(qū)生態(tài)系統(tǒng)彈性力Fig.8 Ecological resilience in coal exploitation affected areas

4 結(jié)論與討論

4.1 討論

4.1.1煤礦開(kāi)采對(duì)礦區(qū)生態(tài)系統(tǒng)的影響

隨著礦區(qū)煤炭基地的建成,工礦用地侵占草地和灌叢面積,加上生態(tài)環(huán)保政策落實(shí)不到位,加劇了礦區(qū)植被退化。由于工礦用地人類(lèi)活動(dòng)頻繁,煤礦開(kāi)采直接影響區(qū)生態(tài)系統(tǒng)穩(wěn)定性較差。大規(guī)模開(kāi)采活動(dòng)可能造成礦區(qū)地面塌陷、地裂縫,導(dǎo)致植物根系直接拉斷,并且土壤物理結(jié)構(gòu)的破壞和水分養(yǎng)分的流失也影響植物生長(zhǎng),造成開(kāi)采區(qū)植被退化和生產(chǎn)力降低。隨著采煤面積的增大可能會(huì)引起地裂隙加大、加深,地表水通過(guò)裂縫和破碎區(qū)連通的位置滲入地下,造成地表水流失；同時(shí)裂隙可能穿透上部隔水層,從而造成地下水滲漏,地下水位下降。在煤礦開(kāi)采過(guò)程涉及和影響的土地范圍內(nèi),原生森林植被生態(tài)系統(tǒng)、土壤結(jié)構(gòu)特征以及其天然的水土保持功能被嚴(yán)重削弱,表層土壤和植被的水土保持能力大幅下降。因此,在煤炭開(kāi)采區(qū)應(yīng)及時(shí)采取生態(tài)修復(fù)治理措施,降低采礦活動(dòng)引起的地質(zhì)環(huán)境問(wèn)題,從而提高礦區(qū)生態(tài)系統(tǒng)彈性力。

4.1.2生態(tài)彈性力指標(biāo)選取合理性分析

本文從水文、土壤、植被3個(gè)方面選取指標(biāo)綜合評(píng)估榆神礦區(qū)生態(tài)系統(tǒng)彈性力,而氣候條件也是決定生態(tài)彈性力的主要因素。通常在區(qū)域尺度較大的情況下評(píng)估彈性力需要考慮氣候要素,而本研究區(qū)氣象要素空間差異較小,因此沒(méi)有直接選取氣象指標(biāo)。由于礦區(qū)本身屬于水資源缺乏地區(qū),水資源的供給能力是生態(tài)彈性力提高的主要限制因素,選用了地表、地下徑流兩個(gè)水文指標(biāo)?？紤]到主成分分析時(shí)各指標(biāo)間需要存在一定的相關(guān)性,選用了土壤保持量、植被覆蓋度和生產(chǎn)力等相關(guān)指標(biāo)來(lái)反映研究區(qū)生態(tài)環(huán)境變化。植被覆蓋度和生產(chǎn)力均可表征生態(tài)系統(tǒng)的自我恢復(fù)能力。植被覆蓋度對(duì)區(qū)域生態(tài)環(huán)境變化有重要指示作用,能夠反映區(qū)域的水土保持、氣候調(diào)節(jié)以及生態(tài)系統(tǒng)的穩(wěn)定性強(qiáng)弱[43],多項(xiàng)研究直接從植被覆蓋度、葉面積指數(shù)等遙感數(shù)據(jù)中檢測(cè)植被異常以測(cè)量生態(tài)系統(tǒng)的恢復(fù)力[9,44]。而植被生產(chǎn)力可表征一個(gè)區(qū)域地形地貌、水分熱力狀況以及植被的生產(chǎn)適宜性,可通過(guò)植被生產(chǎn)力來(lái)客觀反應(yīng)生態(tài)系統(tǒng)各要素的健康程度[43,45],Frazier等利用遙感反演的NPP作為生態(tài)系統(tǒng)健康指標(biāo)監(jiān)測(cè)生態(tài)系統(tǒng)恢復(fù)情況[46]。植被灌層通過(guò)攔截降雨、減少雨滴擊濺、阻擋泥沙運(yùn)輸來(lái)減少水土流失,不同的植被和類(lèi)型都會(huì)影響一個(gè)區(qū)域的水土保持能力[47],而對(duì)于陜北黃土溝壑地區(qū),有研究表明降雨是影響水土流失的主要因素,植被覆蓋與土壤侵蝕的關(guān)系不明顯[48]。同時(shí)土壤組成也是土壤保持量的決定性因素,研究發(fā)現(xiàn)不同性質(zhì)的土壤生態(tài)彈性力存在很大差異[49]。

4.1.3優(yōu)勢(shì)與不足

利用空間二維柵格數(shù)據(jù)進(jìn)行空間主成分分析來(lái)評(píng)估生態(tài)彈性力,高度保留了各指標(biāo)的空間信息,從而可以反映區(qū)域生態(tài)彈性力的空間差異,便于識(shí)別生態(tài)問(wèn)題嚴(yán)重的區(qū)域,可為生態(tài)環(huán)境修復(fù)提供參考依據(jù)。以往的研究大多考慮了生態(tài)彈性力的時(shí)間趨勢(shì),反映的是區(qū)域整體生態(tài)彈性力相對(duì)變化情況。本研究方法更能反映生態(tài)系統(tǒng)彈性力區(qū)域空間格局情況。

由于數(shù)據(jù)分辨率以及指標(biāo)獲取等限制,本研究仍存在一定局限性：數(shù)據(jù)的分辨率較低導(dǎo)致綜合評(píng)估得到的生態(tài)彈性力空間細(xì)節(jié)描述較弱；水文模型校準(zhǔn)時(shí)采用的觀測(cè)數(shù)據(jù)年份較早,且在研究區(qū)的模擬效果只達(dá)到了滿意,得到的徑流量存在不確定性；用于表征生態(tài)環(huán)境變化的因素還有很多,本研究中由于數(shù)據(jù)可獲取性問(wèn)題,所選取的指標(biāo)不夠全面,我們認(rèn)為盡可能全面的評(píng)估對(duì)礦區(qū)開(kāi)發(fā)規(guī)劃更具有指導(dǎo)意義,但本研究在區(qū)域生態(tài)環(huán)境健康評(píng)價(jià)方面提供了一套有效的思路和方法。

4.2 結(jié)論

本研究通過(guò)遙感解譯獲得了榆神礦區(qū)土地利用分布格局,在充分明確礦區(qū)生態(tài)環(huán)境情況的基礎(chǔ)上評(píng)估了榆神礦區(qū)生態(tài)系統(tǒng)穩(wěn)定性,分析了煤炭開(kāi)發(fā)對(duì)礦區(qū)生態(tài)環(huán)境質(zhì)量的影響程度。主要結(jié)論：

近10年來(lái),隨著生態(tài)修復(fù)措施的實(shí)施和社會(huì)經(jīng)濟(jì)的快速發(fā)展,研究區(qū)生態(tài)環(huán)境明顯改善,林地草地面積明顯增加,但是煤炭開(kāi)發(fā)用地也大幅增加,從8.4 km2增加至14.2 km2。

研究期間,榆神礦區(qū)生態(tài)系統(tǒng)抗干擾能力和可調(diào)節(jié)能力總體不斷提高但伴有波動(dòng),生態(tài)系統(tǒng)彈性力從2009 的0.53增加至2018年的0.60,空間上呈現(xiàn)東南部高西北部低的特征,該區(qū)生態(tài)系統(tǒng)彈性力主要受植被覆蓋度、地下徑流量的影響。

煤炭開(kāi)采活動(dòng)引發(fā)的地面塌陷、地裂縫等地質(zhì)環(huán)境問(wèn)題致使煤礦開(kāi)采區(qū)植被退化、地下徑流、地表徑流量降低以及水土保持能力下降,導(dǎo)致煤炭開(kāi)采直接影響區(qū)的生態(tài)彈性力低于間接影響區(qū),平均低9.23%。因此,及時(shí)采取修復(fù)治理措施降低煤炭開(kāi)發(fā)引起的地質(zhì)環(huán)境問(wèn)題是提高礦區(qū)生態(tài)彈性力有效途徑。