西南典型喀斯特地區石漠化時空演變特征
——以貴州省普定縣為例

習慧鵬,王世杰,白曉永,*,唐 紅,吳路華,陳 飛,2,肖建勇,2,王明明,李匯文,操 玥,陳 歡,2,冉 晨,2,羅旭玲,2

1 中國科學院地球化學研究所,貴陽 550081 2 貴州師范大學地理與環境科學學院,貴陽 550025 3 中國科學院大學,北京 100049

中國南方喀斯特地區存在著非常嚴重的石漠化問題[1- 3],在此地區的自然災害和貧窮發生率極高,嚴重制約著該地區的發展。目前,研究其空間演變規律已成為治理石漠化的重要前提,以往學者利用不同方法分析了石漠化的空間演變規律：熊康寧等基于2000年Landsat TM影像和實際調查數據,將貴州省石漠化分為6級并探討了貴州省石漠化的現狀及其空間演變趨勢[4]；胡寶清等基于2000年TM影像提取了廣西省石漠化信息,并利用數理統計等方法探討了石漠化的分布演變規律,得出了人類活動是石漠化產生的主要原因[5]；白曉永等基于1986、1995、2000年3期貴州省石漠化數據,利用轉移矩陣等方法探討了貴州省石漠化的時空演變規律[6]。可以看出,以往的研究多基于實際數據和數理統計的方法綜合分析,所取得的成果對早期的石漠化治理提供了巨大的科技支撐,然而研究方向多集中于省域級大尺度低空間分辨率上,且研究年限多在20世紀80年代—21世紀初,所取得的研究成果和經驗方法對已治理了20余年的石漠化指導意義不足。

此外,如何快速高效地獲取石漠化各因子的貢獻率一直是眾多學者的研究重點,在以往研究中學者們通過數學分析的手段總結出了很多經典的方法。藍安軍等基于主成分分析、相關性分析,并建立了有關石漠化的動力指數公式,對貴州省的石漠化眾多影響因子分析得出了包括植被覆蓋度、喀斯特面積、河網密度、耕地、未利用地、草地、海拔、土地開墾率和25°以上的土地面積等9個因子是導致石漠化產生的驅動因素[7]；胡寶清等基于主成分分析,主要分析了巖性、地貌、水熱條件、土壤和NDVI等5個指標,計算得到石漠化地質背景區的環境生態指數回歸模型[8]；皇甫江云等基于專家打分的方法,對各等級石漠化對應的影響因子歸類打分,得出草地石漠化中各因子的貢獻率[9]。可以看出,以往的經典方法受人為影響較大,且由于喀斯特環境的復雜性,石漠化影響因子回歸模型無法推廣。

因此,針對以往研究中存在的兩大問題,本文選取貴州省中部典型喀斯特地區普定縣為研究區,基于Landsat、DEM 等30 m分辨率數據,提取了植被覆蓋度、巖石裸露率和坡度等3個最具代表性的指標來精確反演石漠化[10- 12],旨在解決以下3個問題：(1)尋找更加高效快速的石漠化模擬反演方法。(2)引入地理探測器,采用更科學的方法確定各因子的貢獻率。(3)研究各等級石漠化分布規律,分析新形勢下貴州典型喀斯特地區的石漠化演變速率,揭示石漠化的時空演變特征(技術路線見圖4)。研究結果中獲得的各因子貢獻率為反演石漠化提供了參考,發現的石漠化分布演變規律為石漠化治理提供了方向。

1 研究區與數據源

1.1 研究區概況

位于貴州省中西部地區的普定縣(圖1),全縣國土面積約1081.8 km2,喀斯特面積為869.64 km2,年均降水量約1200 mm,年均溫15.3℃,平均海拔為1100 m,地勢南北高,中間低。截止2013年,普定縣人口密度達到440人/km2,遠遠高于我國西部其他地區。由于巖溶地貌的廣泛發育,石灰土比重占到了63%,連年的高強度降水和人類活動,使得普定縣石漠化極為嚴重,2008年,普定被確定為全國100個、全省55個石漠化綜合治理試點縣,眾多科研單位在此地區開展石漠化治理工作,普定縣成為了石漠化防治研究的典型地區。

圖1 研究區位置圖,左圖為研究區在貴州省的位置,右圖為生成的普定縣TIN(triangular irregular network)Fig.1 Location map of the study area

1.2 數據源

本研究使用的實驗數據主要包括：普定縣巖性圖、普定縣石漠化分布圖、Landsat TM/ETM+、DEM等30 m空間分辨率數據,具體見表1。

表1 主要數據源

2 研究方法

2.1 石漠化指標

2.1.1 巖石裸露率

巖石裸露作為石漠化地區最顯著的特征,在石漠化信息的提取中有著至關重要的作用,本文擬采用NDRI來獲取研究區巖石裸露率,NDRI(Normalized Difference Rock Index)是指歸一化巖石指數,基于喀斯特碳酸鹽巖的光譜信息在短波紅外波段具有較大差異[13-14],基于NDVI思想,在ENVI 5.2軟件下進行波段運算,根據公式(1)得到研究區的NDRI：

(1)

根據像元二分模型原理,假設像元由巖石裸露和非裸露兩部分組成,根據NDRI,計算得到基巖裸露率(Fr)：

(2)

式中, NDRIr為全由巖石裸露得到的像元值,NDRIo為全無巖石裸露得到的像元值。在統計得到的所有像元值中,取頻率1%的像元值為NDRIo,取累計頻率為99%的像元值為NDRIr。

2.1.2 FVC(Fractional Vegetation Cover)

FVC植被覆蓋度指的是植被冠層在地面的垂直投影面積與所處地總面積之比[15-16],同樣根據像元二分模型原理,本文擬采用NDVI來獲取研究區植被覆蓋度,公式如下：

(3)

在計算得到像元的NDVI值后,統計其百分比,取1%的像元值為NDVIo,作為沒有植被存在的像元值,取99%的像元值為NDVIr,作為完全由植被覆蓋的像元值。

2.1.3 Slope

本文選擇的DEM數據為GDEMDEM 30 m分辨率數字高程數據,數據標識為ASTGTM_N26E105,在ArcGIS下對其進行投影轉換為ALBERS等面積投影,坐標系為Krasovsky_1940_Albers,掩膜提取后制作得到坡度圖。

2.2 貢獻率獲取

地理探測器是探測和利用空間分異性的工具[17],在空間一致的情況下,通過采集樣本點上各自變量和因變量的值,可以利用公式計算自變量對因變量的解釋能力大小,公式表達式為：

(4)

圖2 2000年石漠化及其布點圖Fig.2 Karst rocky desertification and its site map in 2000

圖3 各因子十級分布圖Fig.3 Ten level distribution map of each factor

2.3 石漠化分類定級標準

基于前人工作建立石漠化分級分類標準[10, 18- 21](表2),通過地理探測器得到各因子貢獻率之后,在ArcGIS 10.2軟件平臺下按照下表進行重分類,并疊置分析,基于普定縣巖性圖和土地分類結果剔除非喀斯特地區和水體、城鄉建筑等不可能發生石漠化的地方,得到石漠化空間分布圖。

表2 石漠化分類定級標準

2.4 石漠化空間動態變化分析

石漠化動態變化無論是在時間上還是在空間上都是一個復雜的過程,傳統的數字分析表達能力欠佳,因此本文在探討石漠化時空演變時將數學模型與制圖綜合在一起,旨在探討以下科學問題：第一是明確石漠化因子的變化特征,以探討因變量對自變量的影響大小；第二是明確石漠化的轉移方向；第三是明確各等級石漠化的轉移快慢。因此,有必要對石漠化的因子變化和時空轉變分別構建模型討論。

2.4.1 石漠化因子變化

本文中石漠化的影響因子主要由植被覆蓋度、巖石裸露率、坡度等組成,而坡度由于變化較小無法捕捉,故不在此討論。在石漠化演變過程中,像元中的植被覆蓋度和巖石裸露率等級必然也在發生著變化,因此,兩者的變化必然與石漠化演變具有密切關系,基于此,構建下述公式：

Q=q1-q2

(5)

式中,Q指兩個不同年份的因子在柵格計算器下相減后的值,q1指一段時間內該像元變化前的等級,q2指變化后的等級。

2.4.2 石漠化轉移矩陣

在某一時間段內,石漠化的演變方向與演變規模大多較復雜,且時常因地而異,因此需要定量化地討論以明確其大致的演變方向和趨勢,構建轉移變化矩陣。

(6)

式中,S為面積,單位為km2；m,n為一段時期內前期和末期石漠化的等級；Sk,k為在該時期內石漠化的轉移面積,單位為km2。

2.4.3 石漠化轉移變化速率

僅僅探討石漠化的轉移方向和面積對于治理石漠化的指導意義是遠遠不夠的,石漠化的演變快慢反映了一種類型石漠化在該地區的變化發展模式,在橫向比較上更具有科學價值,也更能說明石漠化的演變情況。因此需要建立模型以探討石漠化的演變速率。

(7)

式中,Pi為某一類型的石漠化演變速率,單位為km2/a；S2為一定時間段內末期的面積,S1為前期的面積,單位為km；T為該時間段,單位為a。

2.5 技術路線

圖4 技術路線圖Fig.4 Technical roadmap

3 結果與分析

3.1 權重分析

在地理探測器中根據石漠化和各因子等級計算得到各自所占的權重,結果表明,三個指標中植被覆蓋度和巖石裸露率主要控制了石漠化的分布：FVC(0.47)>NDRI(0.45)>Slope(0.15),由于三者相加大于1,故各因子權重所占的百分比即為本文中得出的石漠化貢獻率：FVC(44%)、NDRI(42%)、Slope(14%)(表3)。

表3 石漠化各因子貢獻率

Q：百分比權重,quantile；P：置信度,power

3.2 石漠化分布及各因子等級圖

按照表2已劃分的石漠化等級標準,在ENVI 5.2、ArcGIS 10.2軟件下制圖分類,得到普定縣1990—2015年FVC、NDRI、Slope等級圖,并在地理探測器中獲得其各自的貢獻率,在ArcGIS下疊置分析,并結合Landsat影像,在ENVI下得到土地利用分類圖,剔除了水體、城鄉建筑地等不可能發生石漠化的地區,得到各年份石漠化空間分布圖。

3.2.1 NDVI提取植被覆蓋度

首先在ENVI 5.2下波段運算得到的1990—2015不同年份的NDVI,然后再根據公式(3)計算得到植被覆蓋度,并在ArcGIS10.2軟件下利用Reclassify工具進行重分類,按照表2將其分為>71%、51%—71%、31%—50%、21%—30%、11%—20%、<10%共六類并分別賦值為1、2、3、4、5、6,得到1990—2015年植被覆蓋度分布圖(圖5)。由圖中可以看出,越靠近河流和城鎮地區,植被覆蓋度的等級越低。

圖5 植被覆蓋度空間分布Fig.5 The spatial distribution of fractional vegetation cover

3.2.2 NDRI提取巖石裸露率

首先根據公式(1),在ENVI 5.2下波段運算得到的1990—2015不同年份的NDRI, 再利用公式(2)計算得到巖石裸露率,并在ArcGIS下重分類,將其按照表2分為<20%、21%—30%、31%—50%、51%—70%、71%—90%、>90%共六類并分別賦值為1、2、3、4、5、6,得到1990—2015年巖石裸露率分布圖(圖6)。由圖中可以看出,越靠近河流和城鎮地區,巖石裸露率的等級越高,且1997年的巖石裸露面積明顯最小。

圖6 巖石裸露率空間分布Fig.6 The spatial distribution of rock exposedness

3.2.3 DEM提取坡度

在投影變換后,利用ArcGIS軟件中的Slope工具根據DEM數據提取坡度,并按照表2的標準,重分類得到普定縣坡度空間分布圖(圖7),由圖中可以看出,坡度較陡的區域多沿著普定縣SW—NE方向分布,總體呈現出南北坡度大,中部較平坦。

圖7 坡度分布圖Fig.7 The spatial distribution of Slope

3.2.4 石漠化分布圖

在得到各年段的植被覆蓋度、巖石裸露率和坡度分布后,利用ArcGIS 10.2中的疊加分析工具,按照地理探測器計算得出的各因子貢獻率疊加制圖,并剔除非喀斯特地區和水體、城鄉建筑等不可能發生石漠化的區域, 得到1990、1997、2006和2015年的石漠化空間分布圖(圖8)。

圖8 石漠化空間分布圖Fig.8 The spatial distribution of karst rocky desertification

3.2.5 精度驗證

為檢驗反演效果,本文選取了喀斯特科學數據中心提供的2005年普定縣石漠化實際分布圖作為標準,將制作得到的2006年數據進行了對比,由圖9和表4可以看出,石漠化在總體上反演精度較高,達到了80%,該石漠化反演方法可信度較高。

3.3 時空演變分析

3.3.1 石漠化面積統計分析

根據石漠化空間分布圖,利用柵格計算器計算得到普定縣不同年份的各等級石漠化面積,統計如下表5,可以看出,無石漠化面積占比最大最穩定,而極重度石漠化由于面積小,變化面積不明顯,而潛在石漠化、輕度石漠化和中度以上石漠化的面積在不同年段發生著劇烈地變化。1990—1997年間,石漠化面積減少了70.95 km2；1997—2006年間,石漠化面積增加了216.45 km2；2006—2015年間,石漠化面積又減少了111.05 km2；1990—2015年間,石漠化面積增加了34.50 km2。 綜上,普定縣輕度以下石漠化分布面積在喀斯特地區均有大量分布,而中度以上多分布在三岔河流域沿岸附近；在1990—2015年間,普定縣石漠化演變復雜劇烈,經歷了由好轉到惡化再到好轉的階段,演變總體趨勢在好轉,但局部在惡化,說明治理石漠化不可掉以輕心,還需要加大治理力度。

圖9 實際與模擬石漠化對比Fig.9 Comparison between actual and simulated karst rocky desertification

年份Year無石漠化NKRD潛在石漠化PKRD輕度石漠化LKRD中度石漠化MKRD重度石漠化SKRD已石漠化KRD2005年/km2379.10164.79225.5792.527.60325.692006年/km2306.99171.23271.96100.6518.80391.36精度 Precision/%80.9896.2482.9491.9240.5383

表中石漠化類型均為英文簡寫,KRD:石漠化,Karst rocky desertification；NKRD:無石漠化,No karst rocky desertification；PKRD:潛在石漠化,Potential karst rocky desertification; LKRD:輕度石漠化,Light karst rocky desertification；MKRD:中度石漠化,Moderate karst rocky desertification；SKRD:重度石漠化,Severe karst rocky desertification

表5 各級石漠化在喀斯特區域的面積和比例

ESKRD：極重度石漠化,Extremely severe karst rocky desertification

3.3.2 石漠化因子時空演變分析

根據圖5以及公式(5),在ArcGIS 10.2軟件下運用柵格計算器進行減法運算,將各年份植被覆蓋度相減,得到1990—1997、1997—2006、2006—2015、1990—2015年間的植被覆蓋度變化圖(圖10),其中,正值表示該地區的植被覆蓋度等級在升高,負值則表示在下降。可以看出,1990—1997年植被覆蓋度大都在減少,而在2006—2015年卻大都在增加,在1997—2006年間,植被覆蓋度的變化則具有明顯的地理位置特征。

圖10 植被覆蓋度時空演變Fig.10 Fractional vegetation cover changes in each period

再次運用柵格計算器工具,將圖6中的各年份巖石裸露率相減,得到1990—1997、1997—2006、2006—2015、1990—2015年間的變化圖(圖11),同樣,正值表示該地區的巖石裸露率等級在升高,負值則表示在下降。可以看出,1997—2006年巖石裸露率等級大都在升高,而在其他年段則大都在降低,結合普定縣行政區分析,升高的地區多位于普定縣中東部城關鎮一帶。

圖11 巖石裸露率時空演變Fig.11 Rocky exposedness changes in each period

3.3.3 石漠化轉移變化時空演變分析

根據公式(6),在柵格計算器下通過柵格疊加得到某一年份對應另一年份的石漠化面積,統計結果見圖12(NKRD—無石漠化；PKRD—潛在石漠化；LKRD—輕度石漠化；MKRD—中度石漠化；SKRD—重度石漠化；ESKRD—極重度石漠化；Non—karst：非喀斯特區域)和圖13,可以看出：1990—1997年,石漠化好轉的比例(25.06%)>惡化(16.76%)；1997—2006年,石漠化好轉(9.73%)<惡化(37.83%)；2006—2015年,石漠化好轉(28.26%)>惡化(14.20%)；1990—2015年,石漠化好轉(24.65%)<惡化(25.47%)；此外,在圖中可以看出由無、潛在石漠化轉變為重度以上石漠化的面積雖然不大,但占比極高,結合影像分析(圖14),發現此種演變方式集中于三岔河流域沿岸附近,由于河道拓寬,導致沿岸一帶植被減少、基巖裸露率增大,且此地區坡度均較大,是石漠化發生急劇變化的主要原因。而重度以上石漠化轉變為無、潛在石漠化面積較少,結合影像分析(圖14),發現此種演變方式多發生于城關鎮一帶,此處靠近普定縣城,且海拔較低,有利于土壤積累,極少部分較為嚴重的石漠化地區突變為了無石漠化或潛在石漠化,筆者認為是居民的開墾荒地導致附近的石漠化土地演變為了耕地或林地,因此可以看到人類活動的兩面性。

綜上,在石漠化演變過程中,輕度以下石漠化演變最為劇烈,尤其是潛在石漠化和輕度石漠化,演變面積大,方向難測；石漠化的演變方式因地而異,在靠近河流沿岸地區,石漠化常發生突變而惡化,而在城鄉居民地附近,常發生重度石漠化突變為無石漠化。

圖12 石漠化轉移變化分布圖Fig.12 Change distribution map of karst rocky desertification圖中石漠化類型均為英文簡寫,KRD：石漠化,Karst rocky desertification；NKRD：無石漠化,No karst rocky desertification；PKRD：潛在石漠化,Potential karst rocky desertification；LKRD：輕度石漠化,Light karst rocky desertification；MKRD：中度石漠化,Moderate karst rocky desertification；SKRD：重度石漠化,Severe karst rocky desertification；ESKRD：極重度石漠化,Extremely severe karst rocky desertification

圖13 1990—2015年石漠化轉移矩陣Fig.13 Karst rocky desertification transfer matrix in 1990—2015

圖14 石漠化解譯與實際影像和衛星圖對比Fig.14 Comparison of karst rocky desertification with actual image and satellite map

3.3.4 石漠化轉移變化速率

再次運用柵格計算器計算得到普定縣石漠化轉移速率變化圖(圖15),并統計其轉移矩陣如圖16,可以看出重度以上石漠化演變速率極低,幾乎等于0,這也驗證了前人的工作基礎[6, 22]。在1990—1997年間,輕度-潛在石漠化變化速率最大,為11.2 km2/a；1997—2006年間,潛在-輕度石漠化變化速率最大,為16.5 km2/a,2006—2015年間,輕度-潛在石漠化變化速率最大,為10.4 km2/a,1990—2015年間,潛在-輕度石漠化變化速率最大,為3.22 km2/a。計算各等級石漠化平均變化速率得出(圖17),潛在石漠化惡化的速率最大,為2.75 km2/a,其中,在1997—2006年間中速率達到了16.5 km2/a,演變十分迅速。此外,在各年段上,中度以上石漠化整體均在好轉,而輕度以下石漠化則反復無常,尤其是潛在石漠化,好轉的速率極小,而惡化則十分迅速,因此需要加強對輕度以下石漠化的監察與防護。

圖15 石漠化轉移速率分布圖Fig.15 Velocity distribution map of karst rocky desertificationX—X為石漠化等級無變化區域,即石漠化演變速率為0；NK為非喀斯特地區,非喀斯特地區不涉石漠化問題

圖16 1990—2015年石漠化轉移速率矩陣Fig.16 Transfer rate matrix of karst rocky desertification in 1990—2015

圖17 各級石漠化平均變化速率 Fig.17 Average change rate of karst rocky desertification at all grades

4 討論

4.1 石漠化時空演變分析

1)1990年至1997年間,石漠化面積減少了71.0 km2,實際上是在好轉,而在1997年至2006年間,石漠化面積又快速增加了217.0 km2,尤其是潛在向輕度石漠化的轉變面積最大,速度最快,這段時間正處于經濟發展的黃金期,人類不合理活動強度的增大是導致惡化的主要原因,直到2005年后,政府開始實施退耕還林政策,大量石漠化土地被禁止開發,石漠化面積才又開始逐漸減少[23-24],從2006年至2015年,石漠化面積又減少了110.0 km2。整體來看,從1990年到2015年,石漠化面積增加了34.5 km2,普定縣石漠化經歷了由好轉到惡化再到好轉的演變,總體趨勢雖然在好轉,但局部地區仍在惡化,石漠化的治理不可掉以輕心。

2)石漠化轉移變化平均速率最快的是潛在石漠化,平均速率為2.75 km2/a,其中在1997—2006年間,速率達到了16.5 km2/a；變化速率最慢的是極重度石漠化,平均速率為0.0002 km2/a,這也驗證了前人的工作成果,石漠化變化速率與等級高低成反比。此外,好轉速率最快的是中度石漠化,平均速率為0.80 km2/a,惡化速率最快的是無石漠化,平均速率為3.10 km2/a,由圖17看出,輕度以下石漠化變化方向不易確定,而中度以上石漠化則幾乎都在好轉,因此,在治理中度以上石漠化的同時,更要注重輕度以下尤其是潛在石漠化的防治工作。

3)結合石漠化因子轉移變化圖與石漠化空間變化圖可以看出,植被覆蓋度與巖石裸露率的變化成相反關系,在植被覆蓋度增加高的區域,石漠化演變多以好轉為主,多集中在普定縣中部、東南一帶；而在巖石裸露率增加高的區域,石漠化演變多以惡化為主,多集中在三岔河一帶。

4)由圖8分析得出,在石漠化空間分布上,本研究發現普定縣中度以上石漠化多集中分布在三岔河沿岸附近,而輕度以下石漠化則分布較散亂,無太大規律。

5)由圖13、圖14看出,在三岔河附近,常發生著無或潛在石漠化突變為重度石漠化的現象,這往往是由河道拓寬,雨水沖刷,導致無或潛在石漠化地區植被劇減,巖石裸露率暴增導致。此外,在城關鎮靠近縣城一帶,也會出現重度石漠化突變為無石漠化的現象,筆者認為由于此處地勢較低,易于土層堆積,人類的開墾荒地等城鎮化活動使得重度或極重度石漠化也變為了無或潛在石漠化,從而這也體現了人類活動對石漠化演變的積極一面。

4.2 展望與不足

本文使用的數據由于依賴于遙感影像,因此難免受到衛星過境時間和天氣云量等的影響,加之獲得的植被覆蓋度、巖石裸露率等常受到季節變化的干擾,很難去計算全年的平均值,此為最大的局限性。本文采取的方法在總體上反演石漠化精度較高,因子權重較為合理,在喀斯特中小尺度地區內反演石漠化具有一定的優勢。

5 結論

本文基于Landsat、DEM等數據提取得到的普定縣1990—2015年4期石漠化影像展開分析,重點討論了普定縣的石漠化空間演變規律,主要認識如下：1)在1990—2015年間,普定縣石漠化演變復雜劇烈,經歷了由好轉到惡化再到好轉的階段,總體演變趨勢雖然好轉,但局部地區仍舊在惡化,石漠化治理力度需要加大。2)普定縣中度以上石漠化集中分布在三岔河沿岸附近,輕度以下石漠化則無明顯地理分布特征。3)潛在和輕度石漠化演變方向復雜,潛在石漠化變化速率最大(2.75 km2/a),最大時在1997—2006年間達到了16.5 km2/a。4)植被覆蓋度與巖石裸露率的變化呈負相關,前者主導了石漠化的好轉,后者控制了石漠化的惡化。5)石漠化演變存在突變,在三岔河沿岸附近,存在無/潛在石漠化有突變為重度石漠化的現象；在靠近城鄉居民地附近,也存在重度以上石漠化突變為無石漠化的現象。6)植被覆蓋度、巖石裸露率、坡度等能較準確反演石漠化,貢獻率分別為44%、42%、14%。