鄂爾多斯—榆林地區景觀生態風險評價及其驅動因子分析

馬 駿, 裴燕如, 王慧媛, 于 強, 牛 騰, 岳德鵬

(北京林業大學 精準林業北京市重點實驗室, 北京 100083)

隨著社會經濟發展，人類城鎮化加快。一方面土地利用類型不合理利用，人地矛盾尖銳；另一方面大多數國內的礦產資源不能夠合理地開發利用，環境代價巨大，使絕大多數的自然生態系統都受到人類活動的影響，而這些人類活動影響會造成諸多的生態風險，因此科學管理這些生態風險成為人類與自然和諧發展的重要前提[1]。生態風險評價作為現階段生態環境修復、治理以及生態安全格局構建和生態風險預警等諸多工作的重要參考依據，引起了全球學者的高度重視。

生態風險指生態系統及其組分在受到外界干擾后所承受的風險，其會對生態系統的結構和功能會造成不利影響[1-2]，生態風險起源于小尺度的人體健康評價，后來演變為衡量一種或多種脅迫因素對生態系統造成正在發生或潛在發生的具有負生態效益的可能性的研究手段[2-4]。隨著生態風險評價理論和方法的完善，區域生態風險評價的提出，將空間異質性和風險因子等級關系納入到生態風險評價體系。區域生態風險主要應用在土壤、大氣以及水環境的生態評估。趙曉光等人[5]對礦區不同土地利用類型的土壤重金屬污染程度和綜合潛在生態風險危害進行評價,何瑞東等人[6]分析鄭州市某生活區大氣重金屬元素的污染特征，對鄭州市某生活區進行潛在生態風險和居民健康風險的評估。李捷等學者[7]探究了北方六湖沉積物中重金屬污染源的差異性。為了探究景觀格局對生態風險的定量影響，基于區域生態風險評價提出景觀生態風險評價，景觀生態風險評價更加注重尺度效應、空間異質性以及景觀格局對的生態功能、過程的影響。目前國內外對景觀生態風險的研究已經相對成熟，從研究對象上來看，流域[8-10]、礦區[11-12]、城市[13]、自然保護區[14]等已經成為研究的熱點，從研究尺度來看，綜合考慮研究區面積和景觀斑塊的平均面積來劃分網格是主流的研究手段[15-17]，也有部分學者從縣域、市域等行政單位為研究單元進行生態風險研究[18-19]；從評價方法上來看，大多基于風險源匯和景觀格局兩種方法，有的學者還基于生態系統服務功能進行景觀生態風險評估的研究[20]，當前關于景觀生態風險評價的研究主要將土地利用與生態風險相結合并探究其時空分布，鮮有涉及到不同驅動力影響生態風險演變過程的定量分析，在本研究中使用地理探測器探測鄂爾多斯—榆林地區(以下簡稱“鄂榆地區”)生態風險演變過程中的關鍵驅動因子。

鄂榆地區礦產資源豐富，被譽為中國的“能源走廊”。進入21世紀，西部大開發戰略實施，西部資源開發產業化進程加快，煤炭能源開發成為鄂榆地區的支柱產業[21]。隨著人為開發活動的進行，對整個鄂榆地區的草地、林地、未利用地等景觀要素都造成了不同程度的影響，在很大程度上改變了景觀結構的布局，并對生態安全產生威脅。基于此，本文以鄂榆地區為研究區，基于土地利用變化信息和景觀格局構建景觀生態風險評價模型，測定各評價單元的景觀生態風險值，通過地理探測器對引起生態風險的驅動力因素進行探測，確定生態風險時空演變的主要驅動因子，為生態環境治理、修復以及生態安全格局的構建提供科學依據。

1 數據來源與研究方法

1.1 研究區概況

研究區位于我國第二大沉積盆地——鄂爾多斯盆地腹地，地處內蒙古自治區鄂爾多斯市和陜西省榆林市兩個地級市，并以此為邊界向外擴長5 km,將對鄂榆地區生態環境有重要影響的黃河及黃河周邊的地區也能納入研究區范圍內。

鄂爾多斯市位于內蒙古自治區西南部，介于北緯37°35′24″—40°51′40″，東經106°42′40″—111°27′20″，總面積為86 752 km2。榆林市位于陜西省最北部，介于北緯36°57′—39°35′，東經107°28′—111°15′，總面積為43 578 km2。鄂榆地區的東部、北部、西部三面被黃河環繞，北側與庫布齊沙漠相接，毛烏素沙地則地處陜西省榆林地區和內蒙古自治區鄂爾多斯市之間，屬于北溫帶半干旱大陸性氣候區，地勢由西北向東南傾斜，土地利用類型主要以草地和未利用地為主，冬夏寒暑變化大，多年平均氣溫7.2 ℃，日最高溫度39 ℃，日最低氣溫-31.4 ℃。多年平均降水約357 mm，降水主要集中在7—9月，占全年總降水的70%左右。鄂榆地區人口結構主要以漢族占大多數，根據第7次人口普查數據，鄂爾多斯市常住人口為2 153 638人，榆林市常住人口為3 624 750人。鄂榆地區的產業結構主要以第二產業為主，第二產業的比例大致保持在50%左右。2020年，鄂爾多斯市實現生產總值(GDP)為3 533.66億元，榆林市實現生產總值(GDP)為4 089.66億元。

1.2 數據來源及處理

本文數據主要選取鄂榆地區2000和2010年30 m空間分辨率的Landsat 7 TM/ETM+遙感影像和2020年30 m空間分辨率的Landsat 8OLI遙感影像為數據源，將其在ENVI 5.3中對影像進行輻射校正、大氣校正、幾何校正、圖像增強以及鑲嵌裁剪等預處理，得到鄂榆地區3期的遙感影像數據。利用極大似然法進行監督分類，根據目視解譯并結合研究區實際情況，采取LUCC分類體系將土地利用類型分為耕地、林地、草地、水域、城鄉，工礦，居民用地(以下稱建設用地)及未利用地6類，在解譯結果隨機抽取驗證樣地與Google Earth和實地進行選取的1 217個驗證點進行驗證，總體分類精度均達88%以上，滿足研究需求。驅動因子分析中使用的人為干擾度數據是通過于立忠等[22]提出的人為干擾度賦值表，參考谷東起等[23]建立的人為干擾度計算模型計算并將其空間化所得；NDVI數據為NASA提供的MOD13A1產品(http:∥ladsweb.nascom.nasa.gov/)；年降水量數據、年均氣溫數據來自于國家地球系統科學數據共享服務平臺(http:∥www.geodata.cn/);DEM數據來自地理空間數據云平臺(www.gscloud.cn/)；地表蒸散發數據來源于NASA官方遙感影像數據(https:∥ladsweb.modaps.eosdis.nasa.gov/)。

1.3 研究方法

1.3.1 土地利用變化指標構建

(1) 土地利用動態度。在景觀生態風險的研究中，土地利用類型通常被看作是景觀嵌體的類型。土地利用動態度可以描述景觀更替變換的速度與強度，通常定義為指定時間范圍內區域土地利用類型變化的快慢[24]，以表現區域內土地利用的劇烈程度。通常分為以下兩類：

①單一土地利用類型動態變化度。單一土地利用類型動態變化度指的是某研究區一定時間范圍內，某一土地利用類型的面積變化情況[24]，計算公式為：

(1)

式中：K為研究時段內某一土地利用類型動態度；Ua，Ub分別問研究初、研究末期某一土地利用類型面積；T為研究時段(a)。

②綜合土地利用類型動態變化度。綜合土地利用類型動態變化度表示所有土地利用類型變化的整體情況，計算公式為：

(2)

式中：LC為研究時段內土地利用類型年變化率； LUi為研究初期第i類土地利用類型面積； ΔLUi-j為研究初期至研究末期時段內第i類土地利用類型轉為非i類土地利用類型面積的絕對值。

(2) 土地利用類型轉移矩陣。土地利用作為景觀最直觀的表達形式，當強烈人類活動或自然災害發生時，會對土地利用類型造成變更，景觀要素的空間結構以及生態功能也隨之發生改變，從而使區域的景觀生態風險發生變化。土地利用類型轉移矩陣描述研究時段內土地利用類型的變更過程，能夠充分反映研究時段內土地利用類型的變化方向，計算公式為：

(3)

式中：Aij為土地利用類型面積；n為土地利用類型數量；i，j為轉移前后的土地利用類型。

1.3.2 景觀生態風險指數構建 土地利用速度、強度以及方向的變化會引起景觀格局的改變，不同的景觀格局對生態風險的抵御能力是不相同的。為了將景觀生態風險在空間上以定量化的形式呈現，并改善現狀中不利的環境要素，在未來的設計規劃中提出發展性意見，通常通過計算多種景觀指數，按照不同的權重構建景觀生態風險指數。

(1) 劃分評價單元。為了將生態風險指數空間化，結合蘇海民等研究[25]，評價單元大小的選擇應為斑塊平均面積的2～5倍，利用ArcGIS 10.2漁網構建工具對研究區處理，分割為5 km×5 km大小相同的生態風險單元，共劃分5 864個，作為評價單元。

(2) 景觀生態風險模型構建。景觀格局是在多種驅動力在景觀結構綜合作用形成的，反映了人類活動對景觀結構作用的強度、頻率以及范圍，同時也影響著生態系統的格局和過程。根據相關文獻[17,23]利用景觀干擾度和景觀脆弱度構建景觀生態風險評價模型，計算每個評價單元的景觀生態風險指數，將其作為評價單元的中心值。

景觀干擾度是景觀受到外界因素干擾程度的一種定量表征，將景觀破碎度、景觀分離度、景觀優勢度按一定的權重疊加計算所得[26-30]，結果詳見表1。

表1 景觀指數計算公式及意義

景觀脆弱度(Fi)表示不同類型景觀對外界干擾的敏感程度和易損程度，景觀脆弱度在不同階段的景觀自然演替過程中是不同的[31]。根據前人的研究成果[32]，將研究區的景觀分別進行如下賦值：未利用地6，水域5，耕地4，草地3，林地2，建設用地1，歸一化得到景觀類型脆弱度指數。

通過景觀生態風險指數可以建立起景觀結構與生態風險的關系，能夠定量地反映景觀格局與生態風險的相關性，計算方法如下為：

(4)

式中：Aki表示第k個樣區i類景觀的面積；Ak表示第k類樣區的總面積。

1.3.3 空間自相關分析 景觀生態風險在空間上量化后，通常會在空間上按一定的規律分布，而空間相關性分析是對景觀生態風險空間分布模式描述的一種有效方法。空間相關性是檢測空間中某一要素的屬性值是否與相鄰空間點屬性值相關聯的指標。正相關表示某空間單元的屬性值與其相鄰單元的屬性值有相同的變化趨勢，負相關則相反。常用的指標有Moran’sI指數和LISA指數。

(1) Moran’sI指數(全局自相關)是描述空間要素及屬性值在區域內的空間依賴程度，其計算公式為：

(5)

(2) LISA指數(局部自相關)是描述空間單元與其相鄰單元的相似程度，其計算公式為：

(6)

1.3.4 地理探測器 景觀生態風險是由多種風險因子對區域內景觀的綜合影響所產生，是由多種驅動力，包括自然環境、社會經濟環境以及人類活動等多種條件綜合作用所致。地理探測器是探測空間分異性并揭示其背后驅動因子的一種統計學方法。其核心思想是：如果自變量對某個因變量有很強的關系，那么自變量和因變量在空間分布上具有一定的相似性。其中地理探測器中的主要用來因子探測器用來探測自變量對因變量的解釋力,其大小用q值來衡量[33]。因此本文選取因子探測器來分析影響鄂榆地區景觀生態風險的各個驅動因子,計算公式為：

(7)

(8)

2 結果與分析

2.1 土地利用類型時空變化分析

2.1.1 土地利用結構變化 根據圖1和表2可得，在鄂榆地區，草地、耕地和未利用地所占的面積比例最大，共占90%以上，在2000—2010年，林地、建設用地、未利用地呈現增加的趨勢，其中建設用地的面積增加最多，相比2000年增加了0.29%；耕地、草地及水域分別呈現減少的趨勢，耕地的面積減少最多，相比2000年減少了0.21%，其他土地利用類型的變化幅度比較小，總的來說2000—2010年的趨于相對穩定的狀態；在2010—2020年，耕地、水域和建設用地的面積有較明顯的增長，其中耕地面積增加比例最大，相比2010年增加2.01%，其次為建設用地，相比2010年增加1.32%。相反草地的面積減少比例最大，相比2010年減少2.17%，其他土地利用類型變化幅度較小，處于相對穩定的狀態。

圖1 鄂榆地區2000，2010和2020年土地利用變化

表2 鄂榆地區2000-2020年土地利用類型面積變化

2.1.2 土地利用速度變化 鄂榆地區的單一土地利用動態度和綜合土地利用動態度詳見表3，在2000—2010年和2010—2020年2個時間段內，建設用地的土地利用動態度均為最高，分別為5.75和16.47，表明在2000—2020年建設用地的增加速度始終正向增長，且2010—2020年建設用地的擴張增加速度比2000—2010年要快，增加的建設用地主要分布在烏海市、達拉特旗、鄂爾多斯市、準格爾旗伊金霍洛旗、府谷縣、神木縣、榆陽區、靖邊縣以及定邊縣等地區。在研究期間，2000—2010年與2010—2020年的綜合土地利用動態變化度分別為0.19與0.69，表明在2010—2020年研究區整體的土地利用類型間變化比較活躍。

表3 鄂榆地區2000-2020年單一動態度和綜合動態度

2.1.3 土地利用方向變化 2000—2020年土地利用類型發生了明顯的變化(表4)，其中在2000—2010年，草地是主要的轉出類型也是主要的轉入類型，草地轉出為耕地的面積為987.18 km2，占轉出類型總比重的45.52%，轉出為未利用地的面積為504.17 km2，占轉出類型總比重的23.25%，轉出為建設用地的面積為211.26 km2，占轉出比例的9.74%;轉入草地的來源主要是耕地，轉入面積為1 152.12 km2，占總轉入比例的57.93%，其次為未利用地，轉入面積為441.33 km2，占總轉入比例的22.19%。其他土地利用類型的轉入和轉出面積較少。在2010—2020年，草地仍然為主要的轉出類型和主要的轉入類型，轉出耕地面積為5 154.39 km2，占總轉出比例的51.30%，轉出為未利用地的面積為2 156.52 km2，占總轉出比例的21.46%，轉出為建設用地的面積為1 375.14 km2,占總轉出比例的13.69%；轉入草地的主要來源是未利用地和耕地，其中未利用地轉入面積為3 131.30 km2，占轉入比例的44.71%，耕地轉入面積為2 483.85 km2，占轉入面積的35.47%，其他土地利用類型的轉入轉出變化很小。草地景觀作為一種生態脆弱性比較低的景觀，且在鄂榆地區分布相當廣泛，當受到自然和人為因素等作用容易產生性變，而草地作為鄂榆地區主要的景觀類型，其數量的增減會導致景觀生態風險的直接變化。

表4 鄂榆地區2000-2020年土地利用類型面積轉移矩陣 km2

2.2 景觀生態風險時空變化分析

在ArcGIS 10.2中采用克里金插值分析方法，對研究區生態風險進行空間插值，然后根據自然斷點法，將生態風險區劃分為低生態風險、較低生態風險、中生態風險、較高生態風險和高生態風險5類風險區域。

由圖2可知，總體上鄂榆地區的中生態風險、較高生態風險和高生態風險地區主要分布在北部的庫布齊沙漠和中部的毛烏素沙漠地區，較低生態風險和低生態風險地區主要分布在東部、北部邊緣地區、西南部分地區以及庫布齊沙漠和毛烏素沙漠之間的地區。由表5可知，2000—2010年，鄂榆地區的中生態風險地區、較高生態風險地區以及高生態風險地區面積呈現擴張趨勢，擴大比例分別2.02%，0.22%和0.57%，與較高生態風險地區和高生態風險地區相比，中生態風險地區的面積增幅較大，主要由于在研究期間鄂榆地區實施退耕還林(草)、沙漠農業試驗田以及圍封禁牧等大量生態保護政策，對鄂榆地區庫布齊沙漠和毛烏素沙漠地區沙漠化逆轉起到一定的作用，沙漠地區出現草地、林地和耕地等景觀類型，植被作物從無到有，區域景觀類型多樣性增多，導致斑塊破碎度增大，脆弱性增強，生態風險呈現升高趨勢。

表5 鄂榆地區景觀生態風險等級面積及比例

在2010—2020年，鄂榆地區的中生態風險地區、較高生態風險地區及高生態風險地區的面積呈現大幅降低趨勢，減少比例分別為4.28%，1.09%和2.06%，同時低生態風險地區和較低生態風險地區面積大幅增加，增加比例分別為4.48%和2.95%。由圖2可知，在西北部的庫布齊沙漠地區較高生態風險區和高生態風險區面積減少顯著，主要原因是隨著多年生態保護政策的實施，鄂榆地區生態環境演變具有明顯的生態累積效應，生態功能的可恢復性發揮了作用，草地、耕地景觀類型在空間分布上逐漸呈現均勻趨勢，景觀優勢度增加且破碎度減少，脆弱性降低，抵抗外界干擾的能力增強，最終使得生態風險降低。

圖2 鄂榆地區2000-2020年景觀生態風險空間分布

2.3 景觀生態風險空間聚集特征分析

鄂榆地區2000，2010和2020年的全局Moran’sI指數均在0.5以上，分別為0.823，0.826和0.845，說明研究區景觀生態風險的空間分布呈顯著的正相關關系。散點主要集中在一、三象限，相似值呈現聚集狀態，大多數的風險單元在空間上的分布規律呈現高—高和低—低，高生態風險區周圍分布著相應的高生態風險區，低生態風險地區分布著相應的低生態風險區。全局自相關Moran’sI指數在2000—2020年呈現上升趨勢，表明鄂榆地區風險聚集趨勢呈現上升趨勢。

鄂榆地區2000，2010和2020年的局部自相關LISA指數分布圖如3所示。從空間分布上看，3個時期的高—高生態風險地區主要分布在杭錦旗西北部和東北部、鄂托克前旗的西部、烏審旗及榆陽區，少量高—高生態風險地區分布在靖邊縣和準格爾旗周圍，這些區域往往景觀類型多樣性高，景觀優勢度低且景觀分布破碎化，人類活動極易引起以林草退化為主等生態問題。低—低生態風險地區主要分布在鄂榆地區杭錦旗、達拉特旗北部邊緣地區、準格爾旗、府谷縣、神木縣、佳縣、吳堡縣和清澗縣東部邊緣地區，以及正定縣、烏海市、鄂托克旗中部、杭錦旗東南部、鄂爾多斯市、伊金霍洛旗西北部和神木縣北部地區，這些地區景觀類型以耕地為主，景觀優勢度高，且斑塊的破碎度較低。極少數的高—低和低—高生態風險區主要在分布在高—高和低—低生態風險聚集區周圍，說明研究區的生態風險比較穩定，出現急劇變化的可能性較低。

圖3 鄂榆地區2000-2020年景觀生態風險LISA指數分布

2.4 景觀生態風險驅動力分析

鄂榆地區景觀生態風險的時空演變受多種外界干擾因素綜合影響。例如人為、自然及社會經濟因素等。本文采用地理探測器探測各個驅動因子對2000，2010及2020年景觀生態風險的解釋力(圖4)。本文選取人為干擾度、人造地表作為人為因素，利用NDVI值、降水量、年平均氣溫、高程以及地表蒸散發作為自然因素綜合分析鄂榆地區景觀生態風險影響機制的分析指標。基于ArcGIS平臺對7個因子使用自然斷點法進行分級，使用地理探測器測算每個驅動因子的q值，q值越大說明對景觀生態風險的解釋力越強。

圖4 鄂榆地區2000-2020年鄂榆地區景觀生態風險因子解釋力

由圖4可知，在2000—2020年不同時期影響鄂榆地區景觀生態風險的驅動因子的解釋力變化相對顯著。在人為因素中，人為干擾度的q值均大于0.7，且逐年增加，反映了人類活動強度對景觀格局的擾動呈現增強趨勢，人為干擾度對鄂榆地區的景觀生態風險的解釋力越來越增強，成為影響鄂榆地區景觀生態風險的最主要驅動因素；而人造地表的q值小于0.1，對鄂榆地區的景觀生態風險解釋力并不強；在自然因素中，NDVI的q值均在0.6左右，對景觀生態風險有很強的解釋力，年降水量的q值在2000和2010年大于0.3，2020年大于0.1，年平均氣溫的q值保持在0.1附近，這兩種驅動因子在自然因素中僅次于NDVI，但作為植被生長必不可少的因素，在一定程度上也影響著生態風險的變化。高程和地表蒸散發的q值均在0.1以下，其對生態風險的影響作用微乎其微。

在自然因素中，NDVI、年降水量和年均氣溫是主要的驅動因子，其中景觀生態風險與植被覆蓋、降水和氣溫因子的起伏變化密切相關，結合鄂榆地區土地利用數據以及NDVI數據，當地的植被覆蓋范圍在2000—2020年期間有明顯的擴大，庫布奇沙漠和毛烏素沙地的治理中，以生態保護為主的防沙治沙取得到了明顯的成效，沙區的生態環境從不穩定逐漸趨于穩定，由人工干預的植被群落結構也趨于穩定，景觀生態風險逐漸變低；根據鄂榆地區的氣象資料，2000—2010年，鄂榆地區的氣溫呈現降低的趨勢，降水呈現波動上升的趨勢。2010—2020年，氣溫呈現平穩狀態，降水仍呈現波動上升趨勢。與全球變暖和全球降水呈現極端的趨勢不同，最可能原因就是人類活動的影響，改變了地表景觀類型布局，從而改變了當地的氣候條件，進而導致景觀生態風險也隨之發生變化。

在人為因素中，以生態保護為主的人類活動是造成鄂榆地區生態風險變化的主要驅動因子。為切實推動生態建設發展，從21世紀初西部大開發中央政府政策的實施，眾多生態保護措施接連落實。這些年來，鄂榆地區普遍實施退耕還林、退牧還草，禁牧、休牧、輪牧、生態移民等保護措施的推進。鄂榆地區的植被在恢復的過程中，景觀的結構發生劇烈變化，景觀的數量從單一向綜合發展，景觀的布局也由簡單向復雜轉變，從而導致景觀格局在空間的布局發生顯著變化，景觀生態風險也隨之出現先增減后降低的變化。

3 結 論

(1) 鄂榆地區主要的土地利用類型為草地、耕地和未利用地，約占整個研究區的90%以上，在2010—2020年，土地利用強度和速度最為活躍。在研究期間，草地作為主要的轉入轉出類型，轉出為耕地的比例增加5.68%，轉出為建設用地比例增加3.95%，轉出為未利用地的比例減少1.79%，耕地轉入的比例減少22.46%，未利用地轉入的比例增加22.52%。

(2) 鄂榆地區在2000—2010年，中生態風險、較高生態風險和高生態風險地區面積呈現擴張趨勢，低生態風險、較低生態風險地區面積呈現收縮趨勢，其主要原因是人為生態修復活動的影響，沙漠地區逐漸出現耕地、草地、林地等土地利用類型，景觀類型多樣性增強，破碎度變大，景觀優勢度降低，導致生態風險變高。在2010—2020年，景觀中生態風險、較高生態風險和高生態風險低面積呈現收縮趨勢，而低生態風險和較低生態風險地區呈現擴張趨勢，主要由于多年生態政策的實施，景觀類型逐漸均一化，破碎度降低且優勢度增強，導致生態風險降低。

(3) 鄂榆地區2000，2010，2020年的全局自相關分析Moran’sI指數分別為0.823，0.826和0.845，說明生態風險在空間分布上呈顯著正相關，大多數生態風險單元呈現高—高和低—低分布，少數的生態風險單元在高—高和低—低風險單元周圍分布，說明研究區的生態風險比較穩定，出現急劇變化的可能性較低。

(4) 影響鄂榆地區景觀生態風險的主要驅動因子是人為因素中的人為干擾度，自然因素中的主要驅動是NDVI，年降水量和年均氣溫屬于次要驅動因子，人造地表、高程以及地表蒸散發對鄂榆地區景觀生態風險影響不大。