基于地理探測器的江蘇省沿江地區(qū)土壤侵蝕定量分析

吳鎮(zhèn)宇, 羅夢琦, 郭紅麗, 吳 芳, 黃昱楠, 戴玉婷, 劉 霞

(1.南京林業(yè)大學(xué) 南方現(xiàn)代林業(yè)協(xié)同創(chuàng)新中心 江蘇省水土保持與生態(tài)修復(fù)重點實驗室， 江蘇 南京210037; 2.江蘇省水文水資源勘測， 江蘇 南京 210029)

長江經(jīng)濟帶發(fā)展是國家《促進中部地區(qū)崛起“十三五”規(guī)劃》中“三大戰(zhàn)略”之一。水土資源是制約長江經(jīng)濟帶經(jīng)濟社會發(fā)展的主要因素之一，根據(jù)水利部發(fā)布的2019年中國水土保持公報[1]，長江經(jīng)濟帶水土流失面積3.95×105km2，占區(qū)域土地總面積19.18%。在“以生態(tài)優(yōu)先、綠色發(fā)展為引領(lǐng)，依托長江黃金水道，推動長江上中下游地區(qū)協(xié)調(diào)發(fā)展和沿江地區(qū)高質(zhì)量發(fā)展”總體要求下[2]，區(qū)域資源環(huán)境承載力、水土流失綜合整治以及生態(tài)環(huán)境保護任務(wù)更加緊迫。研究區(qū)域土壤侵蝕空間分布特征、識別土壤侵蝕的主要影響因素對指導(dǎo)區(qū)域水土流失綜合防治具有重要的意義。任樂萌、王慧琴[3-4]分別采用地理加權(quán)回歸模型、灰色關(guān)聯(lián)度判定、相關(guān)性分析等方法開展丹江口庫區(qū)、華北北部區(qū)域土壤侵蝕驅(qū)動因子識別及影響程度分析。何琪琳等[5]采用聚類分析與回歸分析進行東北黑土區(qū)坡耕地土壤侵蝕影響因素識別。傳統(tǒng)的回歸分析、相關(guān)系數(shù)分析等方法可以進行土壤侵蝕驅(qū)動力因子識別及影響評價，但難以定量分析多個變量間的耦合作用及對土壤侵蝕空間分異的影響[6-7]。地理探測器是一種可探測變量空間分異性及其驅(qū)動力的統(tǒng)計學(xué)方法，既能探測單個自變量對因變量的影響程度，也能判斷兩兩因子交互作用的強弱與類型，并進行風(fēng)險區(qū)域識別[8]。目前，地理探測器在喀斯特地貌類型區(qū)[9-11]、陜北黃土丘陵溝壑區(qū)已有探索[8]，但在江蘇沿江地區(qū)未見研究與實踐。江蘇沿江地區(qū)是長江經(jīng)濟帶的重要組成部分，區(qū)域經(jīng)濟發(fā)達，開發(fā)利用程度高，人口密度大，人為活動強烈，經(jīng)濟社會發(fā)展需求與水土資源承載能力的矛盾比較突出。為此，本文以江蘇沿江地區(qū)為研究區(qū)域，基于地理探測器對區(qū)域土壤侵蝕狀況與驅(qū)動力的多因子耦合影響進行定量歸因，為區(qū)域生態(tài)環(huán)境保護及社會經(jīng)濟可持續(xù)發(fā)展提供理論依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

江蘇沿江地區(qū)位于江蘇省南部(118°12′—122°0′E，30°45′—33°14′N)，包括南京市、鎮(zhèn)江市、常州市、無錫市、蘇州市、揚州市、泰州市和南通市8市59個縣市區(qū)，土地面積5.10×104km2，占江蘇省土地總面積的47.57%。研究區(qū)地勢平緩，主要地貌類型為平原、丘陵山區(qū)、河湖水域，高程一般在200 m以下，分布有老山山脈、寧鎮(zhèn)山脈、茅山山脈和宜溧山脈，主要河流水系有長江與太湖水系；屬亞熱帶季風(fēng)氣候，多年平均降雨量1 147.19 mm；以棕壤土、潮土、水稻土、褐土等為主；地帶性植被為亞熱帶常綠闊葉林，主要樹種香樟(Cinnamomumbodinieri)、馬尾松(Pinusmassoniana)、杉木(Cunninghamialanceolata)、水杉(Metasequoiaglyptostroboides)、毛竹(Phyllostachysheterocycla)等。研究區(qū)人口密度983人/km2，區(qū)域生產(chǎn)總值7.30×1012元，占江蘇省的78.83%。

1.2 研究方法

1.2.1 數(shù)據(jù)源與處理 基于2 m分辨率GF-1與ZY-3遙感影像(2017年3—11月，2018年1—6月)，結(jié)合野外調(diào)查，通過人機交互解譯獲取土地利用數(shù)據(jù)。采用250 m空間分辨率MODIS植被指數(shù)產(chǎn)品計算標準植被指數(shù)(NDVI)與植被覆蓋度(FVC)。采用5 m分辨率數(shù)字高程模型(DEM)提取坡度、坡長。采用1981—2017年研究區(qū)188個站點逐日降雨量數(shù)據(jù)計算多年平均降雨量及降雨侵蝕力。土壤可蝕性因子數(shù)據(jù)來源于第一次全國水利普查成果[9]。

1.2.2 中國土壤流失方程 采用中國土壤流失方程(CSLE)計算土壤侵蝕模數(shù)，計算公式為：

A=RKLSBET

(1)

式中：A為土壤流失量〔t/(hm2·a)〕;R為降雨侵蝕力因子〔(MJ·mm)/(hm2·h·a)〕;K為土壤可蝕性因子〔(t·hm2·h)/(hm2·MJ·mm)〕;L為坡長因子(無量綱);S為坡度因子(無量綱);B為植被覆蓋與生物管理措施因子(無量綱);E為工程措施因子(無量綱)；T為耕作措施因子(無量綱)。模型各因子計算方法參照文獻《中國水力侵蝕抽樣調(diào)查》[10]。

1.2.3 地理探測器原理及應(yīng)用 地理探測器法是通過探測要素空間異質(zhì)性來揭示其驅(qū)動力的空間統(tǒng)計學(xué)方法，包括4個模塊[7]，采用解釋力(q值)度量自變量對因變量空間分異性的影響程度。本文以地貌形態(tài)類型、坡度、植被覆蓋度、多年平均降雨量、土地利用類型為自變量，土壤侵蝕模數(shù)為因變量，采用因子探測器、交互作用探測器對土壤侵蝕影響因素的定量歸因進行探索，采用風(fēng)險探測器識別土壤侵蝕風(fēng)險區(qū)域。

根據(jù)地理探測器要求，需要對土壤侵蝕影響因子數(shù)據(jù)進行離散化處理[7]。結(jié)合研究區(qū)特征，將地貌形態(tài)類型劃分為平原區(qū)、丘陵山區(qū)與水域3類；坡度分為<2°，2°～5°，5°～10°，10°～15°，15°～20°，20°～25°，25°～30°，30°～35°，>35°共9個等級；植被覆蓋度以0.1為步長，將區(qū)間范圍0～1按等距法分為10個等級；多年平均降雨量以50 mm為步長，將區(qū)間范圍993.04～1 452.55 mm按等距法分為10個等級；土地利用類型劃分水田、旱地、果園、茶園、林地、草地、城鎮(zhèn)村建設(shè)用地、生產(chǎn)建設(shè)活動擾動、交通運輸用地、水域及水利設(shè)施用地10個類型。通過對數(shù)據(jù)重采樣及坐標系轉(zhuǎn)換，在研究區(qū)范圍內(nèi)采用行列數(shù)均為300的漁網(wǎng)均勻布設(shè)17 492個樣點，提取各影響因子基礎(chǔ)信息，作為地理探測器運行數(shù)據(jù)。各因子如附圖1—5所示。

2 結(jié)果與分析

2.1 區(qū)域水土流失狀況

2.1.1 水土流失空間分異性 研究區(qū)水土流失面積823.93 km2，占土地總面積1.62%，空間分異性明顯，主要分布在區(qū)域西南部，以南京、鎮(zhèn)江、常州、無錫為主(附圖6)。侵蝕強度以輕度為主，占水土流失面積的77.1%；中度及其以上侵蝕強度面積占比較小。

2.1.2 不同地貌形態(tài)類型下水土流失特征 研究區(qū)平原區(qū)占總面積的70.2%；丘陵山區(qū)占17.96%，主要分布于區(qū)域西南部；水域占11.84%。水土流失面積丘陵山區(qū)占91.73%；平原區(qū)占8.27%。

2.1.3 不同土地利用類型下水土流失特征 研究區(qū)水田占總面積的33.77%，水域及水利設(shè)施用地占27.53%，城鎮(zhèn)村建設(shè)用地占18.43%，林地與交通運輸用地分別占8.76%，7.06%，果園、茶園、旱地等僅占4.45%。水土流失主要分布于林地、生產(chǎn)建設(shè)活動擾動、果園等土地利用類型，分別占水土流失面積的51.1%，17.28%和13.33%，均以輕度侵蝕為主；旱地、茶園、草地、城鎮(zhèn)村建設(shè)用地與交通運輸用地的水土流失面積較少，共占18.29%；中度及其以上強度侵蝕集中分布于生產(chǎn)建設(shè)活動擾動。但從水土流失與土地利用類型面積占比看，生產(chǎn)建設(shè)活動擾動、果園、茶園水土流失占該土地利用類型面積的比例均在10%以上，說明這3類是區(qū)域內(nèi)容易發(fā)生水土流失的土地利用類型。

2.2 土壤侵蝕影響因子顯著性研究

以地貌形態(tài)類型與土地利用類型為切入點，運行因子探測器對區(qū)域土壤侵蝕狀況進行分類分區(qū)定量歸因?？紤]水域地貌形態(tài)類型無侵蝕、水田與水域及水利設(shè)施用地土地利用類型無侵蝕，本文僅圍繞平原區(qū)與丘陵山區(qū)，對8個土地利用類型開展分析。

2.2.1 不同地貌形態(tài)類型下土壤侵蝕影響因子顯著性 不同地貌形態(tài)類型下研究區(qū)土壤侵蝕的主導(dǎo)影響因子均為土地利用類型和坡度，但解釋力q值有差異(表1)。平原區(qū)以坡度為主導(dǎo)因子(0.055 1)，其次是土地利用類型(0.013 0)，由于平原區(qū)土壤侵蝕的空間分異性較小，各影響因子的解釋力均較??；丘陵山區(qū)則以土地利用類型為主導(dǎo)因子(0.112 6)，其次是坡度(0.036 8)。在丘陵山區(qū)，除坡度因子外，其余因子解釋力均高于平原區(qū)；可能因為丘陵山區(qū)地形較為復(fù)雜，導(dǎo)致植被覆蓋、土地利用的空間差異性增大，多重因素組合對土壤侵蝕產(chǎn)生影響，增加了土地利用、植被覆蓋、降雨量同土壤侵蝕模數(shù)空間變異的相似性。

表1 不同地貌形態(tài)類型影響因子顯著性

2.2.2 不同土地利用類型下土壤侵蝕影響因子顯著性 不同土地利用類型下研究區(qū)土壤侵蝕的主導(dǎo)影響因子及解釋力q值有顯著差異(表2)。茶園以植被覆蓋度為主導(dǎo)因子，除茶園外的其他類型土地均以坡度為主導(dǎo)因子。從不同影響因子來看，坡度對果園、生產(chǎn)建設(shè)活動擾動的影響最大，解釋力分別為0.607 1和0.590 4，其次為草地和城鎮(zhèn)村建設(shè)用地，茶園和林地解釋力在0.1以下；地貌形態(tài)類型對草地、生產(chǎn)建設(shè)活動擾動的影響最大，解釋力分別為0.131 2和0.109 2，其余土地利用類型解釋力均在0.1以下；降雨量對果園的影響最大，為0.167，其余土地利用類型解釋力均在0.1以下。林地各影響因子解釋力均在0.1以下，與各因子的空間相似性都比較低有關(guān)，可能是因為其層次結(jié)構(gòu)比較復(fù)雜，受林下灌草和枯落物影響較大。

2.3 土壤侵蝕影響因子交互作用研究

2.3.1 不同地貌形態(tài)類型下土壤侵蝕影響因子交互作用 采用交互探測器可知，不同地貌形態(tài)類型下，土壤侵蝕影響因子兩兩交互作用均呈現(xiàn)非線性增強的趨勢，以坡度與土地利用類型的協(xié)同作用為首位，遠高于其他交互作用(表3)。由表3可知，平原區(qū)主導(dǎo)交互作用解釋力(0.700 6)明顯高于丘陵山區(qū)主導(dǎo)交互作用(0.484 9)，說明平原區(qū)更易受坡度和土地利用類型的影響，與平原區(qū)水土流失主要位于生產(chǎn)建設(shè)活動擾動有關(guān)。平原區(qū)土壤侵蝕的交互作用前3位以坡度為核心疊加其他因子，表明在地勢相對平坦的區(qū)域，坡度是土壤侵蝕的主導(dǎo)影響因子。丘陵山區(qū)土壤侵蝕的交互作用前3三位以土地利用類型為核心疊加其他因子，可能因為丘陵山區(qū)地形復(fù)雜，土地利用類型的空間差異性增大，與坡度、降雨量、植被覆蓋因子疊加增加了與土壤侵蝕空間變異的相似性。

表2 不同土地利用類型影響因子顯著性

表3 不同地貌形態(tài)類型土壤侵蝕影響因子交互作用及其解釋力(q)

2.3.2 不同土地利用類型下土壤侵蝕影響因子交互作用 不同土地利用類型下，土壤侵蝕影響因子交互作用呈現(xiàn)非線性增強和雙因子增強趨勢，排在前三位的交互作用詳見表4。由表4可知，城鎮(zhèn)村建設(shè)用地、果園、生產(chǎn)建設(shè)活動擾動、交通運輸用地、旱地、林地的交互作用排在第一位的均為坡度∩降雨量，城鎮(zhèn)村建設(shè)用地與果園解釋力最高，分別為0.959 4和0.766 8；茶園與草地的交互作用排在第一位的均為坡度∩植被覆蓋度，解釋力分別為0.919 1和0.725 8。除果園、茶園、林地與草地外，其他土地利用類型排在第二位的交互作用均為坡度∩地貌形態(tài)類型。果園、茶園排在前3位的交互作用解釋力較強，均在0.6以上；林地的交互作用解釋力均低于0.3，與林地各影響因子空間差異性較小有關(guān)。

2.4 土壤侵蝕高風(fēng)險區(qū)域識別

采用風(fēng)險探測器對研究區(qū)水土流失集中分布的果園、茶園、林地與生產(chǎn)建設(shè)活動擾動等土地利用類型進行高風(fēng)險區(qū)域識別。果園、茶園、林地的土壤侵蝕高風(fēng)險區(qū)域均位于坡度>20°，降雨量1 300～1 400 mm、植被高覆蓋度的丘陵山區(qū)。其原因如下： ①由于這些區(qū)域植被覆蓋度雖高，但林下蓋度較低，林分結(jié)構(gòu)相對單一，林下灌草與枯落物少； ②由于果園、林地土壤侵蝕為坡度主導(dǎo)，茶園為坡度∩植被覆蓋度主導(dǎo)，丘陵山區(qū)坡度較大導(dǎo)致土壤侵蝕風(fēng)險較高。生產(chǎn)建設(shè)活動擾動的土壤侵蝕高風(fēng)險區(qū)域位于坡度25°～30°，降雨量1 300～1 400 mm的丘陵山區(qū)。

3 結(jié) 論

(1) 根據(jù)因子探測器結(jié)果，平原區(qū)土壤侵蝕以坡度為主導(dǎo)因子，但解釋力不高(0.055 1)；丘陵山區(qū)土壤侵蝕以土地利用類型為主導(dǎo)因子，解釋力為0.112 6。不同土地利用類型土壤侵蝕影響因子的層間差異顯著。茶園土壤侵蝕以植被覆蓋度為主導(dǎo)因子，解釋力高達0.867 1；除茶園外的其他類型土地的土壤侵蝕以坡度為主導(dǎo)因子，解釋力范圍為0.037 7～0.607 1。

(2) 兩兩因子交互作用均能增強對土壤侵蝕的解釋力。平原區(qū)與丘陵山區(qū)主導(dǎo)交互作用的因子為坡度疊加土地利用類型，解釋力分別為0.700 6和0.484 9，且顯著高于其他交互作用。不同土地利用類型的主導(dǎo)交互作用的因子主要為坡度疊加降雨量，對城鎮(zhèn)村建設(shè)用地與果園影響較大，解釋力分別為0.959 4和0.766 8；其次為坡度疊加植被覆蓋度，對茶園、草地影響較大，解釋力分別為0.919 1和0.725 8。

(3) 不同土地利用類型土壤侵蝕高風(fēng)險區(qū)域均位于坡度>20°，降雨量1 300～1 400 mm，植被高覆蓋度的丘陵山區(qū)。

4 討 論

對土壤侵蝕影響因子進行定量歸因，有助于開展區(qū)域水土流失綜合防治工作。但土壤侵蝕是由多因子耦合影響的。目前研究土壤侵蝕影響因素主要采用多元回歸分析、相關(guān)分析、主成分分析、灰色關(guān)聯(lián)度等方法。上述方法對數(shù)據(jù)有限制，且僅能分析自變量對因變量的影響，難以反映其耦合影響和空間分異情況。本文使用的地理探測器對自變量與因變量無需假設(shè)與限定，可探測數(shù)據(jù)型數(shù)據(jù)與定性數(shù)據(jù)。該方法操作簡便，結(jié)果包括單因子影響、兩因子交互作用及風(fēng)險區(qū)域分析。目前已在西部喀斯特地貌類型區(qū)、黃土丘陵溝壑區(qū)有過嘗試，在東部江蘇沿江地區(qū)進行探索十分必要。本文僅從區(qū)域?qū)用娣治鐾寥狼治g影響因子，但江蘇沿江地區(qū)的平原沙土區(qū)地貌以及生產(chǎn)建設(shè)項目密集的特征，同時也是該區(qū)域土壤侵蝕關(guān)注的重點，可作為下一步的研究方向，結(jié)合地理探測器進行研究和探索。