魯甸地震的滑坡物質運移規律與地形特征

陳曉利 劉春國 傳一健 蘭 劍 魏延坤

1)中國地震局地質研究所，北京 100029 2)中國地震臺網中心，北京 100045 3)北京大學數學科學學院，北京 100871

0 引言

滑坡在構造與侵蝕的相互作用過程中起著重要作用(Montgomery，2001；Larsenetal.，2010)。構造運動使山脈隆升、海拔高程及局部高差增大；侵蝕作用通過使山坡物質從高處向低處運移來限制山脈的高度，并調節或平衡了構造隆升和河流下切作用。在構造活動區，滑坡使山坡變陡的速率可能超過河流下切作用的影響，它引起的物質運移對局部地貌演化有著重要影響(Malamudetal.，2004；Fanetal.，2012)。與長時間緩慢的侵蝕過程不同，地震滑坡是瞬間的突發性事件，一次大地震能產生成百上千的滑坡，所造成的剝蝕、堆積等將使震區地貌發生顯著改變。一些規模較大的滑坡會對局部地形產生重要影響，而地形條件是地震滑坡的一個重要控制因素(Hsü，1975；Ashfordetal.，1997；Chigiraetal.，2006；Meunieretal.，2008；Huangetal.，2009)。

地震滑坡涉及從山坡臨界狀態的形成、地震事件的觸發到滑坡體運動、堆積和震后滑坡體演變的過程(Montgomery，2001)。每一次地震滑坡之后，山坡物質的分布狀態就會重新進行調整：一方面，臨界失穩狀態的物質以 “滑坡”的方式運移進入一個相對穩定的狀態；另一方面，物質運移的過程形成了一些新的不穩定區域，如滑坡體的后緣將成為潛在滑坡危險區。總體上，地震誘發的滑坡是對地表物質均衡狀態的調整。

Keefer(2002)及Malamud等(2004)基于大量實例研究了地震滑坡引起的侵蝕速率。根據震級與誘發滑坡數量的關系，通過樣本統計分析獲得一定震級下相應的地震誘發滑坡總體積的計算公式，即地震事件造成的物質剝蝕量。這類關于地震滑坡剝蝕量的計算一般將整體的滑坡體積(運移物質)平均到一定范圍，獲得的是該區域內的平均剝蝕量。例如，汶川地震后，Parker等(2011)依據目前廣泛應用的地震滑坡體積計算公式和震后SAR數據反映出的地形變化，計算了同震滑坡產生的物質剝蝕量和構造隆升形成的物質凈增量，認為汶川地震在地震復發周期內將導致山體物質的減少。實際上，地震誘發滑坡的空間分布受到震級、斷裂、河流、地形等的影響，使滑坡相對集中地分布在一些區域，并造成這些區域的地形變化。例如，汶川地震誘發了100多個規模>50i000m3的大型滑坡，這些滑坡受到發震斷裂的控制，主要分布在地表破裂帶附近5km范圍內(許強等，2010)；由于滑坡規模大，地表破裂帶附近的地貌受到較為明顯的影響(王運生等，2008；祁生文等，2009；張永雙等，2009；吳樹仁等，2010)。除此之外，此次地震還形成了828個堰塞湖，這些堰塞湖規模不等，泄洪時間長短不一，無疑會影響河流附近的地貌形態(Fanetal.，2012)。

從汶川地震誘發的滑坡可以看出，發震斷裂和河流附近是大型滑坡的聚集區，也是地表物質侵蝕的主要場所。中強地震所誘發滑坡的規模和分布范圍相對于高震級地震都會有所降低，但在一定條件下也可能觸發較大規模的滑坡，造成顯著的地貌改觀，如2014年魯甸MS6.5地震觸發的紅石巖滑坡形成了高達120m的堰塞壩(Changetal.，2016)。研究表明，該次地震觸發滑坡的空間分布沒有明顯受震中或斷裂的控制，一些規模較大的滑坡發生在牛欄江及其支流等具有較陡邊坡和較大地形高差的河谷兩岸，顯示出局部地形對滑坡分布的重要控制作用(Chenetal.，2015)。而魯甸震區現有的地形條件是該區長期地貌演化的結果。

2014年魯甸地震的震區位于青藏高原東南部邊緣與華南塊體之間，地貌上屬于云貴高原的一部分。中新世以前，這里形成了較完整統一的準平原。由于經歷了中新世末及以后的多次構造運動的強烈抬升，以及伴生的侵蝕、溶蝕等外動力作用，原來的準平原被雕塑、改造成現代的高山峽谷相間的山間凹地景觀(Changetal.，2016)。河流下切使河岸變陡、岸坡巖體強度降低，為沿岸滑坡的發生創造了條件。魯甸地震觸發的紅石巖大滑坡的數值模擬結果表明，地震動的觸發起到了重要的作用(陳曉利等，2015)。然而，地震作用僅是大型滑坡的一種觸發條件，并不是必要條件。2017年茂縣新磨村滑坡(Suetal.，2017；Huetal.，2018；Wangetal.，2018)、2018年10月發生在金沙江上游地區的滑坡堵江等事件都不是由地震直接觸發的，而是其他多種因素長期效應的累積結果，其中，地形可能是較重要的條件(Yangetal.，2019)。

本文以魯甸地震影響區域為研究區，以網格為基本分析單元對研究區的地貌特點進行定量描述和統計分析，在此基礎上開展地貌參數與魯甸地震滑坡的相關性分析，揭示滑坡分布與地貌特點的耦合關系。與以往從滑坡體本身所處的地質地貌特征分析研究滑坡與影響因素之間的相互關系不同，本文著重從區域的地貌特點出發，通過地貌特征參數自相關性的研究與比較判斷地表物質的分布特征，尋找其中的分異性單元(滑坡易發區)，從而認識滑坡發生的地貌特點。之后，進一步基于研究區地貌分布特點對未來研究區內可能發生物質運移的位置進行了預測。

1 數據與方法

1.1 數據

本文研究區是圍繞魯甸地震震中的一個矩形范圍，幾乎覆蓋了此次地震誘發的全部滑坡(圖1)，整體面積約為263km2。將研究區按照0.01°(經度)×0.1°(緯度)大小進行網格劃分(每個網格的面積為1.104i1km2)，得到238個單元。每個單元具有4個基本屬性，即表示滑坡嚴重程度的滑坡發生率和3個表示地形地貌特征的參數：平均坡度、平均高程和地形高差。各個參數的定義如下：

圖1 魯甸地震誘發的滑坡分布Fig.1 Distribution of landslides triggered by 2014 Ludian MS6.5 earthquake.數字表示單元內的滑坡發生面積比率(LAR)

(1)滑坡面積發生率

計算每個網格內的滑坡面積比(Landslide Area Ratio，LAR)，即滑坡面積占網格單元面積的百分比：LAR=滑坡面積/單元面積(%)，代表滑坡發生率。圖1 為2014年魯甸地震誘發滑坡面積的發生率分布。研究區共有111個網格單元受到影響，這些單元的滑坡發生率變化范圍為0.03%～36.1%，對應的滑坡平面面積為300～398i544m2。其中20個單元的LAR>10% ，是受到魯甸地震誘發滑坡影響較嚴重的地區，這些區域幾乎都沿河谷兩岸分布。如SE向距離震中約8km的龍泉河、沙壩河沿岸有2處的滑坡發生率達到36%；在震中西部震中距為4km處，沿牛欄江河谷也有一處的滑坡發生率達到26.4%。

(2)地形參數

基于研究區震前的SRTM 30m數字高程模型(DEM)(Shuttle Radar Topography Mission 30)計算得到每個網格單元的平均坡度(S)、平均高程(H)與地形高差(R)。這些地形參數所表達的含義如圖2 所示。

圖2 網格單元地形參數釋義Fig.2 Topographic parameters of grid cell.S1—S5為通過DEM數據計算出的坡度；L為單元邊長，R為局部地形高差，θ為期望坡度

每個網格單元的平均坡度(S)定義為該網格所包含的從DEM數據點中獲得坡度的算數平均值，即

(1)

式中，n為網格內的坡度數據點的個數，Si為網格內各數據點的坡度。研究區的平均坡度分布如圖3 所示。

圖3 魯甸地震誘發的滑坡分布與平均坡度Fig.3 Landslides and average slope angles in the Ludian seismic area.數字表示單元內的平均坡度

與平均坡度的計算類似，平均高程表示為

(2)

式中，n為網格內的高程數據點的個數，Hi為網格內各個點的高程。

局部地形高差(R)是單元格內最大高程與最小高程的差值，即

R=max(elevation)-min(elevation)

(3)

1.2 方法

一定區域范圍內，地形條件對形成滑坡具有重要的作用。在對研究區進行網格劃分的基礎上，可以通過統計分析獲得這些子區域地貌特征的相關性認識。如果某些子區域的地貌參數特點與區域整體地貌特點之間存在較大差異，則認為這些子區域的地貌特征與區域特征不符，是需要進行物質調整的單元，即易于發生滑坡的單元。

如何尋找可以衡量某一類地貌參數變化的物理量則需要根據不同的研究目的設定。楊曉平等(2019)通過分析南迦巴瓦地區的地形坡度和高程變異系數對活動斷裂進行了識別。本文提出了一個新的網格單元坡度的度量值：期望坡度(Expected Slope Angle，簡稱ES)，表示數據的差異程度。對于每個單元，根據其地形高差(R)和單元的邊長(L)計算出坡度角θ，

θ=tan-1(R/L)

(4)

該坡度角θ即為期望坡度ES。

從圖2 及其計算方式可以看出，期望坡度(θ)和算數平均坡度(S)之間在物理意義和計算方法上均存在差異。期望坡度是以網格單元為度量尺度的描述地形特點的一種表示方法，僅與單元的局部地形高差和單元邊長密切相關；而算數平均坡度與單元內坡度計算點本身的坡度數值及計算點的數量相關。可以認為，期望坡度是對一個分析單元實際地形的一種平滑。研究區期望坡度的分布見圖4。

圖4 魯甸地震誘發的滑坡分布與期望坡度Fig.4 Landslides and expected slope angles in the Ludian seismic area.數字表示單元內的期望坡度

這里，引入統計學中 “變異系數”的概念來解釋期望坡度和平均坡度的關系。變異系數又稱 “標準差率”，是衡量資料中各觀測值變異程度的另一個統計量。具體到本文所研究的問題，即通過局部地形高差獲得的ES就是一個表示差異程度的量，而平均坡度又是一個表示平均數的量。從這種角度上，可以定義一個用期望坡度/平均坡度的指標，該指標能夠在去除量綱的情況下分析不同地形的差異。

依據其統計學意義，可以用期望坡度和平均坡度這二者的相互關系來刻畫地貌特征。具體而言，依據這二者之間的相對大小，可以將地貌的類型分為3類，如圖5 所示(即后文圖7 中A、B、C 3塊區域的劃分依據)。

圖5 平均坡度和期望坡度的相互關系代表的不同地貌特征示意圖Fig.5 Sketch showing topographic features derived from relationship between average and expected slope angles.θ為期望坡度角，A、B、C曲線為簡化的地貌形態，其含義見正文。S1—S5是網格單元內數據點的分區，為便于理解，按照海拔由高到低的順序排列數據點S1—S5

圖6 坡度、高程、地形高差的關系Fig.6 Relationships between slope angles，elevations and topographic reliefs in the study area.

圖7 研究區網格單元的平均坡度和期望坡度的關系及滑坡發生率Fig.7 Average slope angles(S)versus expected slope angles(ES)of grid cells and landslide occurrence rates denoted by landslide area ratios(LAR)in the study area.橙色點代表2014年魯甸地震滑坡發生率LAR>10% 的單元，數字是其LAR值。A、B、C代表的區域地貌特征見正文解釋

(1)A型地貌：地形單元的期望坡度>平均坡度，單元的局部地形高差較大，大部分地區的坡度較為緩和；

(2)B型地貌：地形單元的期望坡度<平均坡度，單元內的坡度普遍較陡；

(3)C型地貌：地形單元的期望坡度幾乎與平均坡度相等，表明單元內的坡度分布或梯度變化較為均勻。

在以上地貌類型劃分依據基礎上，可進一步分析地貌類型與同震滑坡的關系。

2 結果

(1)總體地形特征

統計分析表明，研究區全部分析單元的平均坡度、地形高差和高程三者兩兩之間存在不同的相關性(圖6)：1)平均坡度與局部地形高差呈正相關，局部地形高差值越大則平均坡度值就越大；2)平均坡度和高程呈負相關，高程越大，則平均坡度值越小；3)地形高差與高程呈負相關，隨著高程值的逐漸增大，局部地形高差則變小。

Liu-Zeng等(2008)在研究青藏高原北部、中部、東南部地貌特征時提出了正地形和負地形的概念：如果高程與地形高差、坡度正相關，則為正地形；如果高程與地形高差、坡度負相關，則為負地形。其研究表明，青藏高原北部和中部以正地形為特征，東南部以負地形為特征，反映了水系侵蝕對于控制高原地貌的重要性。本文的研究區位于青藏高原東南部，高程與坡度和地形高差呈負相關，故屬于負地形地貌，這與Liu-Zeng等(2008)的認識是一致的，說明研究區內河流侵蝕是影響地貌特征的重要因素。此次魯甸地震誘發滑坡主要沿河谷分布，也是河流侵蝕作用的一種間接表現，說明局部地形條件對地震滑坡具有重要控制作用。

(2)滑坡分布特征

研究區各個單元的期望坡度和平均坡度的相關性分析結果表明，二者之間存在很好的正線性相關關系(圖7)，其擬合公式為

y=0.940i5x+5.099i9 (R2=0.822i9)

(5)

R2>0.8，表明具有顯著的回歸效果。

然而，除了大量樣本點沿擬合線密集分布外，還有一些樣本點呈現出較為明顯的離散性，接下來將對這部分樣本點進行分析。為了分析這些離散程度較大單元的特點，首先需要確定哪些單元點是離散點。為此，采用

y=0.940i5x+5.099i9±2.5

(6)

作為上限和下限篩選出與回歸擬合線相差較大的點(即ES與S的差異較大，位于圖7 中的矩形框外)。在這樣的篩選準則下，R2=0.822i9這一擬合結果能夠保證62%的單元位于擬合線±2.5的范圍內(即圖7 中的矩形框內)，這部分樣本點也就是該線性關系的主體數據。之后，再將剩余的38%的離散點根據其位于矩形框的左上側(18%)與右下側(20%)分為2類(圖7 中矩形框外的A區和B區)。

在這樣的劃分原則下，A區共有43個單元，其中有7個滑坡面積 ≥10% 的單元；B區有49個單元，其中也有7個滑坡面積 ≥10% 的單元。因此，在研究區內20個LAR≥10% 的單元中，共有14個單元的期望坡度與平均坡度有較大的差異，位于圖7 中矩形框外的A區和B區內，占總數的70%(14/20=70%)。該結果表明，不同的地貌環境與滑坡空間分布之間具有較為密切的關系，滑坡發生密度較大的單元其地貌特點與研究區整體的地貌特點具有較為明顯的差異。因此，通過分析研究區整體的地貌特點，可以較好地識別出潛在的滑坡發生單元。我們的分析結果也認為魯甸震區的地貌主要以C類型為主，即單元內以坡度分布較為均勻的地貌類型占主導。

3 討論

對于一定范圍的區域，高低起伏的地形是該區長期物質運移、分配的結果。不同地形具有不同的重力勢能。在地貌演化的各階段中，無論局部地形高差大或小、地形坡度陡或緩，每個階段都存在一個相對穩定的狀態，它反映了構造隆升與侵蝕之間的平衡(Montgomery，2002)。如果用地貌參數的相互關系來體現這個穩定狀態，則應存在某種相對穩定的關系，反映出地貌侵蝕的主要機制。因此，在地貌參數數據的分析上，這種相對穩定的關系可表現為一定范圍內各個單元都遵循某種規則，存在一定的自相關性。

本文研究區的高程與坡度和地形高差呈負相關，表明研究區內河流侵蝕是影響地貌特征的重要因素。在這樣的地貌侵蝕機制下，盡管本研究所提出的期望坡度與平均坡度是從不同尺度對地形特點進行描述，但二者之間表現出較高的相關性(R2=0.822i9)，表明2種衡量地形變化程度的參數之間存在相對穩定的關系，代表一種占主導的地貌形態所具有的特征。因此，可利用期望坡度作為衡量標準來剔除一些 “特殊”的單元，這些特殊單元與Bl?the等(2015)提出的過量載荷地形具有一定的相似之處。Bl?the等(2015)認為，若一定區域范圍內達到穩態的地貌形態，坡體應同時保持一個特征性的臨界坡度，超過這個臨界角度的地形就稱為過量地形。過量地形代表了一種不穩定狀態的地形，其在地形特征分析上應該表現出一定的 “特殊性”。本文中與穩定狀態不符的單元，則可能是存在過量地形載荷分布的區域。這些單元將最終通過滑坡等物質運移方式進行地形高差或坡度的調整，從而靠近該區的分布規律(如擬合線)以達到與整體的和諧。

這種地形地貌分布上的差異性與滑坡空間分布及滑坡規模等之間的關系值得進行深入探討。2017年8月8日發生于四川省阿壩州的九寨溝MS7.0地震觸發了大量滑坡崩塌。與魯甸地震誘發的滑坡相比，九寨溝地震所誘發滑坡的規模均較小，以淺表型滑坡為主。這2次地震發生處的地質構造背景和地貌環境也不相同。九寨溝地震發生在青藏高原中部巴顏喀拉次級塊體的東緣地區，南北地震帶中段的東昆侖斷裂帶東端和岷江斷裂帶的交會區域。地貌上，九寨溝地震震區位于青藏高原與四川盆地兩大地貌單元過渡的、具有河流深切割特征的高山峽谷地帶，地勢由西向東快速降低，形成較陡的地形，最高海拔>4i500m，最低海拔約為1km，坡度一般>30°。采用本文提出的期望坡度和平均坡度的分析方法對九寨溝震區的地貌特征進行分析的結果如圖8 所示，與魯甸地震震區的分析結果具有顯著不同。九寨溝地震震區的平均高程和平均坡度之間、平均高程與地形高差之間沒有明顯的相關性，難以用正地形或負地形進行描述；而期望坡度與平均坡度間的關系更為復雜。如果魯甸地震震區的地貌參數特征可以代表以河流侵蝕為主的地貌類型，且河谷是大型滑坡主要的發育場所，那么九寨溝地震震區的地貌參數分析結果所代表的地貌類型與其眾多的淺表性滑坡分布之間還需進行進一步研究。

圖8 九寨溝地震震區平均高程與平均坡度和地形高差的關系(a)、平均坡度與期望坡度的關系(b)Fig.8 The relationship between average elevation and average slope and topographic relief(a)and the relationship between average slope and expected slope(b)in the Jiuzhaigou earthquake affected area.

4 結論

地貌演化具有階段性和動態性的特點，一定區域內，各個地貌單元具有向相對穩定狀態演化的趨勢。滑坡作為山坡物質運移的方式之一，在地貌演化過程中起到了重要的作用。2014年魯甸MS6.5地震誘發滑坡是該區在短時間內發生的地貌物質運移的集中調整過程。

這些地震觸發的滑坡主要沿河流分布，表明河流侵蝕使河岸地形變陡、強度降低，形成發生物質運移的有利條件，也增加了地震滑坡的易發性。本文基于SRTM 30m數字高程模型(DEM)，利用網格方法對魯甸地震區的地形特征進行定量分析，結果表明，研究區的高程與坡度、地形高差呈負相關，反映出顯著的河流侵蝕效應；地形特征在期望坡度與平均坡度2個不同尺度下表現出很高的一致性，可能代表著地貌在穩定演化階段所具有的一種特征。而與該特征不符的地貌單元，則是可能發生滑坡進行物質調整的區域。2014年魯甸地震觸發的大部分滑坡發生在期望坡度與平均坡度差異較大的區域，表明了區域地貌特點對滑坡總體上的控制作用。作為對比的九寨溝地震震區的地貌參數分析結果則表現出不同的特點。這種對地貌特點與滑坡分布及規模之間相互關系的認識有助于提高地震誘發滑坡預測的準確性及對潛在較大規模的滑坡場所進行早期識別。

致謝審稿人的意見和建議使本文得以完善，在此表示衷心感謝!