科爾沁左翼中旗土地沙化時空動態(tài)變化及驅(qū)動力分析

韓阿茹汗 包斯琴 馬雅婷 那日嘎拉 張佳華

(內(nèi)蒙古農(nóng)業(yè)大學(xué)沙漠治理學(xué)院，內(nèi)蒙古 呼和浩特 010018)

土地沙化是指因風(fēng)力、重力等外營力作用，導(dǎo)致地表組成物質(zhì)中的細(xì)顆粒損失，或者原有地表因為外來沙礫覆蓋等原因而導(dǎo)致土地退化的過程[1，2]。科左中旗境內(nèi)有著中國4大沙地之一的科爾沁沙地，是防風(fēng)治沙的重要區(qū)域，曾是水草豐沛的科爾沁大草原，進(jìn)入21世紀(jì)以來，受人口劇增、濫墾濫牧等多種因素的影響，導(dǎo)致土地退化嚴(yán)重、水土流失加劇，生態(tài)環(huán)境急劇惡化[3]。科爾沁沙地殲滅戰(zhàn)是習(xí)近平總書記親自部署的“三北”攻堅戰(zhàn)3大標(biāo)志性戰(zhàn)役之一，2023年7月31日，內(nèi)蒙古自治區(qū)建設(shè)我國北方重要生態(tài)安全屏障促進(jìn)條例正式發(fā)布，明確表示了土地沙化荒漠化綜合防治應(yīng)因地制宜，突出重點(diǎn)治理，全面提升生態(tài)系統(tǒng)的質(zhì)量和穩(wěn)定性，因此分析科左中旗土地沙化時空動態(tài)變化，了解土地沙化的空間分布具有現(xiàn)實(shí)意義。

隨著遙感技術(shù)的快速發(fā)展，衛(wèi)星圖像成為監(jiān)測和分析土地利用時空變化的重要數(shù)據(jù)支撐，谷歌地球引擎遙感大數(shù)據(jù)平臺因其獲取數(shù)據(jù)便捷，處理計算數(shù)據(jù)準(zhǔn)確，已廣泛成為遙感監(jiān)測的重要工具[4]。目前對于土地沙化的監(jiān)測方法，野外實(shí)地調(diào)查與遙感目視解譯結(jié)合，這種方法耗時耗力，不適用于長時間序列的監(jiān)測；多指標(biāo)綜合監(jiān)測法，單一的地表參量指標(biāo)不能夠全面反映區(qū)域沙化，因此想要監(jiān)測結(jié)果更加合理，應(yīng)選取多個地表參量進(jìn)行共同監(jiān)測。如，曾永年等[5]利用歸一化植被指數(shù)(Normalized Difference Vegetation Index，NDVI)和地表反照率指數(shù)(Surface Albedo)構(gòu)建了特征空間，選取沙漠化遙感監(jiān)測差值指數(shù)(Desertification Difference Index，DDI)作為監(jiān)測指標(biāo)。劉曉茜等[6]通過建立Albedo-NDVI特征空間沙化監(jiān)測模型法對騰格里沙漠東南緣格林灘綠洲20年間沙化土地進(jìn)行動態(tài)監(jiān)測及分析。該模型充分利用了多維遙感信息，可以自動提取土地沙化信息，可以避免單一地表參量的誤差，有利于土地沙化的監(jiān)測與分析。

科爾沁左翼中旗作為科爾沁沙地的重要部分，亟需對土地沙化監(jiān)測與驅(qū)動力分析。基于此，本文以科爾沁左翼中旗為研究區(qū)，通過提取歸一化植被指數(shù)與地表反照率指數(shù)分析2000—2021年土地沙化的時空演變特征及轉(zhuǎn)移特征，并利用地理空間探測器模型探究自然驅(qū)動因子和社會經(jīng)濟(jì)驅(qū)動因子對土地沙化的影響，進(jìn)而對土地沙化提供新的治理思路，為防沙治沙提供科學(xué)依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

科爾沁左翼中旗位于內(nèi)蒙古自治區(qū)通遼市(E121°08′～123°32′，N43°32′～44°32′)，地處大興安嶺東南邊緣，西遼河北岸，地勢西北高東南低。旗內(nèi)氣候?qū)僦袦貛Т箨懶约撅L(fēng)氣候，氣候干燥多風(fēng)沙。土壤類型主要包括草甸土、風(fēng)沙土、栗鈣土、鹽土、沼澤土，植被類型有沙生植被、草甸植被、鹽生植被等。

圖1 研究區(qū)地理位置

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

1.2.1 遙感數(shù)據(jù)

科爾沁左翼中旗2000年、2005年、2010年、2015年、2021年的Landsat影像數(shù)據(jù)均下載自GEE(Google Earth Engine)云平臺。其中，2000年、2005年、2010年為TM影像，2015年和2021年為OLI影像，分辨率均為30m。

1.2.2 驅(qū)動因子數(shù)據(jù)

本文選擇自然驅(qū)動因素和社會經(jīng)濟(jì)驅(qū)動因素的10個驅(qū)動因子?xùn)鸥駭?shù)據(jù)進(jìn)行土地沙化驅(qū)動研究。其中，年平均氣溫和年降水量數(shù)據(jù)通過國家氣象科學(xué)中心(http：//www.geodata.cn)途徑獲得，人口空間分布、GDP、土壤類型、土壤侵蝕、平均風(fēng)速、日照時數(shù)均來源于中國科學(xué)院資源環(huán)境科學(xué)與數(shù)據(jù)中心(https：//www.resdc.cn)，坡度數(shù)據(jù)使用Arcgis 10.8坡度功能計算，DEM數(shù)據(jù)來源于地理空間數(shù)據(jù)云(http：//www.gscloud.cn)。

1.3 研究方法

1.3.1 Albedo-NDVI特征空間

地表反照率(Albedo)是重要的遙感反演參數(shù)，反映了地球表面對太陽短波輻射反射特效。Albedo受到地表粗糙度、地表含水量因素影響而發(fā)生變化。參照梁順林[7]建立的反演模型來估算地表反照率，計算公式：

Albedo=0.356B+0.130R+0.373NIR+0.085SWIR1+0.072SWIR2-0.0018

(1)

式中，B為藍(lán)光波段的反射值；R為紅光波段的反射值；SWIR為短波紅外反射值，利用Arcgis柵格計算器計算出各年份的地表反照率。

歸一化植被指數(shù)(NDVI)是一種使用較為廣泛的植被指數(shù)，是現(xiàn)階段充分體現(xiàn)植被具體生長情況的一種指數(shù)，依據(jù)該指數(shù)可知植被的覆蓋量、植物量以及生長狀況等數(shù)據(jù)[8]。可通過遙感影像近紅外波段和紅光波段的反射值計算，計算公式：

NDVI=(NIR-R)/(NIR+R)

(2)

式中，NIR為近紅外波段的反射值；R為紅光波段的反射值。

為了后續(xù)數(shù)據(jù)處理的方便，對二者進(jìn)行歸一化處理：

N=[(NDVI-NDVImin)/(NDVImax-NDVImin)]×100%

(3)

A=[(Albedo-Albedomin)/(Albedomax-Albedomin)]×100%

(4)

式中，N為歸一化后的NDVI；A為歸一化后的Albedo；NDVImin、Albedomin分別為NDVI和Albedo的最小值；NDVImax、Albedomax分別為NDVI和Albedo的最大值。

對歸一化后的Albedo和NDVI數(shù)據(jù)，利用Origin軟件進(jìn)行線性回歸分析并繪制特征空間散點(diǎn)圖，可以得出，在Albedo-NDVI特征空間中，Albedo和NDVI有著明顯的負(fù)線性關(guān)系，計算得到擬合方程Albedo=a×NDVI+b，a為方程斜率，b為方程截距，即隨著土地沙化程度的加劇，NDVI的值不斷變小，而Albedo的值逐漸變大。高中低等級的植被覆蓋區(qū)域可以明顯的劃分，反映了沙化土地的地表信息和水熱組合狀況。其線性擬合方程：

Albedo2000=-0.6243×NDVI+0.8697

(5)

Albedo2005=-0.6490×NDVI+0.8552

(6)

Albedo2010=-0.7785×NDVI+0.9322

(7)

Albedo2015=-0.6413×NDVI+0.8742

(8)

Albedo2021=-0.4038×NDVI+0.6445

(9)

圖2 Albedo-NDVI特征空間散點(diǎn)圖

1.3.2 沙漠化遙感監(jiān)測差值指數(shù)計算

在以上構(gòu)建的Albedo-NDVI特征空間中不難發(fā)現(xiàn)，兩者的關(guān)系為線性負(fù)相關(guān)，若將在代表土地沙化變化趨勢的垂直方向劃分，即在線性擬合的趨勢線上劃分二維空間，由此反演出沙漠化差值指數(shù)模型(DDI)[9]，通過此模型可以區(qū)分不同等級的沙化土地，進(jìn)而可以對不同程度的沙化土地進(jìn)行監(jiān)測和分析，計算公式：

DDI=(-1/a)×NDVI-Albedo

(10)

式中，DDI為沙漠化差值指數(shù)；a為特征空間線性擬合的斜率。計算結(jié)果表達(dá)式：

DDI2000=1.6015×NDVI-Albedo

(11)

DDI2005=1.5408×NDVI-Albedo

(12)

DDI2010=1.2844×NDVI-Albedo

(13)

DDI2015=1.5593×NDVI-Albedo

(14)

DDI2021=2.4765×NDVI-Albedo

(15)

結(jié)合《第五次全國土地荒漠化和沙化監(jiān)測技術(shù)規(guī)定》和Jenk自然間斷點(diǎn)法[10]，將DDI分為極重度沙化、重度沙化、中度沙化、輕度沙化和非沙化。自然間斷點(diǎn)法具有分級的良好界限，可以把同一類別中的數(shù)據(jù)值差異值差異最小化，把類別之間的數(shù)據(jù)值差異最大化。

1.3.3 Theil-Sen Median趨勢分析和Mann-Kendall檢驗分析法

Theil-Sen median斜率估計被廣泛運(yùn)用在分析長時間序列數(shù)據(jù)的趨勢分析中，是一種使用較普遍計算非參數(shù)統(tǒng)計的方法。作為一種檢驗非參數(shù)統(tǒng)計方法，Mann-Kendall優(yōu)點(diǎn)在于對趨勢是否線性不做要求，也無需測量值服從正態(tài)分布[11]。Sen趨勢分析和M-K顯著性檢驗已成為長時間序列數(shù)據(jù)趨勢判斷的重要方法，與線性回歸分析相比，具有很強(qiáng)的抵抗數(shù)據(jù)誤差的能力，不易受異常值干擾，在一定程度上提高了準(zhǔn)確性。對于土地沙化的變化趨勢類型，根據(jù)實(shí)際情況并參考郭曉娜研究方法[12]劃分為嚴(yán)重退化、輕微退化、穩(wěn)定不變、輕微改善、明顯改善。

1.3.4 土地沙化轉(zhuǎn)移矩陣

土地沙化轉(zhuǎn)移矩陣可知其在不同階段，存在不同程度沙化類型相互轉(zhuǎn)化的特征，反映不同程度沙化之間的轉(zhuǎn)移速率和土地沙化的時空格局變化，其表達(dá)式：

(16)

式中，i與j分別為研究始末的土地沙化類型；Pij為土地沙化類型i轉(zhuǎn)為土地沙化類型j的數(shù)量。

1.3.5 地理探測器

地理探測器是王勁峰等[13]通過研究得出一種全新的統(tǒng)計學(xué)方法，用以詮釋背后驅(qū)動因子，其核心理念是探討若某個自變量影響某個因變量，在空間分布上面，自變量與因變量存在相似性特征，可基于空間異質(zhì)性原理分析各驅(qū)動因素對土地沙化的影響。根據(jù)計算得到的q值判斷各影響因子的影響力大小，依據(jù)q值可知，自變量影響因變量的程度，兩者呈現(xiàn)正相關(guān)關(guān)系；p值則代表統(tǒng)計學(xué)中的置信度檢驗標(biāo)準(zhǔn)，p值<0.05則表示為差異顯著。地理探測器由因子、交互、風(fēng)險以及生態(tài)等4種探測手段組成，本文使用因子探測進(jìn)行分析，其計算公式：

(17)

2 結(jié)果與分析

2.1 科爾沁左翼中旗土地沙化空間分布特征

利用Arcgis 10.8軟件，通過沙化差值指數(shù)法計算得到研究區(qū)2000年、2005年、2010年、2015年和2021年5期土地沙化遙感影像，空間分布圖如圖3所示。

圖3 科爾沁左翼中旗土地沙化空間分布圖

由圖4可知，2000—2021年土地分布情況出現(xiàn)較大變化，有一部分土地沙化程度好轉(zhuǎn)，但也有一部分土地沙化程度出現(xiàn)惡化。在西部和東北部出現(xiàn)極重度和重度土地沙化；西北部和東部存在中度沙化土地，其分布較為分散；中部和西南部主要分布的是輕度沙化土地；西南部則主要分布的是非沙化土地。總體而言，研究區(qū)土地沙化程度類型主要以中度沙化和輕度沙化土地分布為主。

圖4 2000—2021年科爾沁左翼中旗土地沙化等級轉(zhuǎn)移矩陣

2.2 科左中旗土地沙化面積變化特征

由表1可以看出，在研究區(qū)中，不同等級沙化土地面積均出現(xiàn)不同變化。2000—2021年，研究區(qū)極重度沙化土地面積出現(xiàn)先減后增加的情況，相較于2000年，2021年極重度沙化土地面積減少了0.02%；重度沙化土地面積整體上呈增長趨勢，在2015年下降至19.72%，隨后在2021年又上漲到21.38%；中度沙化土地面積整體呈減少的趨勢，2021年占比26.15%；輕度沙化土地面積整體呈下降趨勢，在2015年增加到27.94%，隨后在2021年又下降到20.70%；未沙化土地面積在2021年增加到23.53%。

表1 2000—2021年科爾沁左翼中旗沙化土地面積及比例

2.3 土地沙化程度轉(zhuǎn)移分析

由圖5可知，2000—2005年，中度轉(zhuǎn)重度和輕度轉(zhuǎn)中度較為顯著，轉(zhuǎn)化面積分別為901.71km2和813.19km2，同時也有699.88km2的中度轉(zhuǎn)為輕度；2005—2010年，輕度轉(zhuǎn)中度和重度轉(zhuǎn)中度較為明顯，轉(zhuǎn)化面積分別為804.82km2和761.56km2；2010—2015年，面積轉(zhuǎn)化較大的是中度轉(zhuǎn)輕度和重度轉(zhuǎn)中度，轉(zhuǎn)化面積分別為888.78km2和751.10km2，說明此階段有較多面積的沙化土地好轉(zhuǎn)；2015—2021年主要是輕度和中度之間的轉(zhuǎn)化；2000—2021年，整體上也是輕度和中度之間的轉(zhuǎn)化，有704.61km2的中度轉(zhuǎn)成了輕度，說明總體上土地沙化有好轉(zhuǎn)趨勢。

圖5 2000—2021年科爾沁左翼中旗土地沙化趨勢空間分布圖

2.4 土地沙化變化趨勢分析

通過Theil-Sen Median趨勢分析和Mann-Kendall檢驗結(jié)果對2000—2021年研究區(qū)土地沙化22年的變化趨勢進(jìn)行分析，從表2可以看出，變化趨勢以穩(wěn)定不變?yōu)橹鳎佳芯繀^(qū)總面積的46.71%。呈改善趨勢的土地占比為30.75%，主要分布于研究區(qū)中部；但也有22.52%的土地呈退化趨勢，主要分布于研究區(qū)西部和東北部。其中，明顯改善土地和嚴(yán)重退化土地均占0.01%。總的來看，沙化狀況趨于穩(wěn)定發(fā)展，說明應(yīng)在加大治理沙化土地的同時，堅持邊治理邊鞏固治理成果的原則。

表2 2000—2021年科爾沁左翼中旗土地沙化趨勢分析

表3 科爾沁左翼中旗土地沙化驅(qū)動因素解釋力統(tǒng)計

2.5 土地沙化驅(qū)動因素分析

利用地理探測器模型對選取的驅(qū)動因子進(jìn)行離散化處理和分類發(fā)現(xiàn)，研究區(qū)土地沙化受到多種驅(qū)動因素的影響。由圖3可知，各個驅(qū)動因子顯著性影響土地沙化的程度，且得出其影響程度的順序為土壤類型>DEM>平均氣溫>平均風(fēng)速>日照時數(shù)>年降水量>GDP>人口空間分布>土壤侵蝕>坡度。其中，土壤類型、DEM和年平均氣溫3種因子影響程度均超過10%，是主要影響因子。社會經(jīng)濟(jì)因子中GDP和人口空間分布q值排序靠后，這說明自然因素比社會經(jīng)濟(jì)因素對土地沙化影響程度更大。

3 討論

本文基于GEE大數(shù)據(jù)平臺，通過Landsat遙感影像提取了研究區(qū)NDVI和Albedo數(shù)據(jù)，從土地沙化時空演變、轉(zhuǎn)移矩陣、地理探測器模型等方面對科爾沁左翼中旗土地沙化進(jìn)行監(jiān)測與驅(qū)動因素分析。綜合分析，2000—2021年研究區(qū)土地沙化主要分布在中部和東北部，整體呈穩(wěn)定趨勢發(fā)展，改善土地面積大于退化土地面積，這與眾多研究科爾沁沙地沙化的研究結(jié)果較相似[14]。土壤類型、DEM和年平均氣溫3種因子是土地沙化的主要影響因子。科爾沁左翼中旗處于我國北方農(nóng)牧交錯帶，生態(tài)環(huán)境脆弱。從研究方法來看，選用歸一化植被指數(shù)和地表反照率2種地表參量相結(jié)合的方法，綜合了對植被敏感的可見光和近紅外波段的反射光譜信息，二者擬合的線性函數(shù)具有明顯的負(fù)相關(guān)關(guān)系，此研究結(jié)果與曾永年等[5]研究結(jié)果一致，曾永年等指出，隨著植被指數(shù)的增加，地表反照率減小，土地沙化過程得到了明顯的反映。李嘉豪等[15]利用該指數(shù)模型研究得到小灤河流域土地沙化動態(tài)變化，本研究也認(rèn)為，沙漠化差值指數(shù)可作為土地沙化的監(jiān)測指標(biāo)。地理探測器模型現(xiàn)已成為主流地理分析模型，利用此模型可定量化地揭示不同驅(qū)動因子對土地沙化的影響大小。本文選取了地理探測器模型中的因子探測模塊，后續(xù)應(yīng)結(jié)合其他3種模塊做進(jìn)一步分析。土地沙化是我國當(dāng)前最為嚴(yán)重的生態(tài)環(huán)境問題之一，近年來，隨著退耕還林還草、實(shí)施“山水林田湖草沙一體化保護(hù)和修復(fù)工程”等一系列治理措施，科爾沁沙地治理成效顯著，但仍有小部分區(qū)域改善不明顯。因此，應(yīng)加大治沙防沙治理力度，因地制宜，阻止土地沙化進(jìn)一步發(fā)展，優(yōu)化生態(tài)環(huán)境。

4 結(jié)論

2000—2021年，研究區(qū)沙化土地主要是輕中度沙化且集中在東北部和中部。對比2015年和2021年數(shù)據(jù)可知，雖然增加了2.80%的重度和極重度的沙化土地，但同時增加了8.52%的非沙化土地，說明沙化土地有好轉(zhuǎn)的跡象；2000—2021年，研究區(qū)沙化土地變化趨勢以穩(wěn)定不變?yōu)橹鳎娣e占比46.71%，輕微改善面積占比30.75%，輕微退化面積占比22.52%；通過地理探測器模型探測，在自然因素和社會經(jīng)濟(jì)因素中，研究區(qū)土地沙化主要受到自然因素影響，其中土壤類型是最主要的驅(qū)動因子，q值為0.1279，其次是DEM和年平均氣溫，q值分別為0.1244、0.1055。