基于ArcGIS的望谟河流域地势起伏度分析

高蜻,唐丽霞,谷晓平,于飞,成星霖,黄松

(1.贵州大学林学院,550025,贵阳；2.贵州省金沙县禹谟镇政府,551804,贵州金沙；3.贵州省山地环境气候研究所,50000,贵阳)

基于ArcGIS的望谟河流域地势起伏度分析

高蜻1,2,唐丽霞1†,谷晓平3,于飞3,成星霖1,黄松1

(1.贵州大学林学院,550025,贵阳；2.贵州省金沙县禹谟镇政府,551804,贵州金沙；3.贵州省山地环境气候研究所,50000,贵阳)

地势起伏度的研究对中小尺度流域的滑坡、泥石流等地质灾害的敏感性分析和危险性评价具有重要的理论意义和实际应用价值。以山洪地质灾害频发的望谟河流域为研究区,基于ArcGIS平台以望谟河流域1∶1万的DEM为基础数据,采用窗口分析法和均值变点分析法,利用Python模块编程自动提取和计算望谟河流域地势起伏度的最佳统计单元,并对望谟河流域的地势起伏度进行分析。结果表明:望谟河流域地势起伏度最佳统计单元的面积是1.49 km2,地势起伏度在0～648.23 m之间；中部河网密集区起伏度大,流域边缘起伏度小；流域地貌属于小起伏山地,以切割高地和切割山地为主,土壤侵蚀较严重。

地势起伏度；望谟河流域；地质灾害；ArcGIS

地势起伏度(relief amplitude)是以DEM为手段,分析地形起伏特征的宏观指标,近年来被广泛应用于区域水土流失分析、土壤侵蚀敏感性评价、土地利用评价、地质环境评价、地形分析、生态环境评价等多方面[1-6]。此外,开展区域地势起伏度研究可为区域滑坡、泥石流敏感性和危险性评价提供基础数据,其提取精度直接影响区域滑坡、泥石流危险性评价和风险评估的正确性[7]。许多学者围绕着提取和确定地势起伏度最佳统计单元的方法进行了研究；但多数研究主要集中于全国范围和地形地貌复杂的特定大区域上,比例尺多小于1∶5万[1-6,8-13],对于以小流域比例尺为1∶1万DEM为基础数据研究的较少。然而,小流域作为地貌中最基本的自然地理单元,是最能反映地貌宏观整体特征的区域。采用大比例尺、高分辨率的DEM分析地势起伏度,对分析小流域的滑坡、泥石流等地质灾害的发生、发展及小流域的综合治理尤为重要[8-10]。

望谟县是贵州省山洪灾害一级重点防治县,县内河流深切,沟谷发育,地形破碎。近年来更是频发山洪、滑坡、泥石流等地质灾害。在“20110606”特大山洪泥石流灾害中,境内2 h降雨量达250 mm,产生1 000万m3的泥石流,震惊全国。笔者以望谟县的重灾流域——望谟河流域为研究区,以ArcGIS为平台,以望谟河流域1∶1万的DEM为基础数据,采用窗口分析法和均值变点分析法,利用Python模块编程,自动提取和计算望谟河流域地势起伏度的最佳统计单元,并对望谟河流域地势起伏度的情况进行分析,以期为当地灾害危险度区划,当地地质灾害频发的原因分析提供依据。

1 研究区概况

望谟河流域地处贵州省望谟县中部,云贵高原向广西丘陵过渡的斜坡地带,地理坐标位于E 106° 2′～106°12′,N 25°9′～25°23′。总面积194.31 km2,最高海拔1 675 m,最低海拔545 m(图1)。该流域内河流属季节性(雨源性)山区河流,年均降水量1 240.3 mm,汛期多在5—8月,降水集中,降水量占全年的87.3%。自1979年至今,已发生19起大型洪水。近年来崩塌、滑坡和泥石流等各种地质灾害频繁发生,严重威胁当地人民生命财产安全。

图1 研究区地理位置Fig.1 Location of the study region

2 数据来源与研究方法

2.1 数据来源与预处理

本研究采用的基础数据为望谟河流域1∶1万比例的流域地形图矢量化制作的DEM,栅格分辨率采用5 m×5 m。以望谟河流域水文监测站为流域出口,ArcSWAT为平台,划定流域范围和子流域,计算流域面积。根据ArcSWAT软件要求以及各投影特性,选用保留面积性质的等积投影(Albers投影)作为投影方法,第1标准纬线定为N 25°,第2标准纬线定为N 47°,中央经线定为E 105°,基准面选用krasovsay1940。

2.2 地势起伏度的提取

2.2.1 提取原理 按照地貌发育的基本理论,一种地貌类型存在一个使最大高差达到相对稳定的最佳统计窗口,这个最佳统计窗口即为地势起伏度的最佳统计单元[8,11]；因此,地势起伏度指在所指定的分析区域内所有栅格中最大高程与最小高程的差值。

式中:Ri为该领域范围内的高差,m；Hmax为领域内所有像元高程值的最大值,m；Hmin为领域内所有像元高程值的最小值,m。本研究采用窗口分析法提取地势起伏度,窗口分析法的基本原理是通过对栅格数据系统中的1个、多个栅格点或全部数据,开辟一个有固定分析半径的窗口,并在该窗口内进行诸如极值、均值、标准差等一系列统计计算,或者进行差分及与其他层面的信息进行必要的复合分析,从而实现栅格数据有效的在水平方向扩展分析[12-13]。

2.2.2 地势起伏度的提取步骤 第1步,以ArcGIS为平台,制作望谟河流域的DEM,进行洼地填充等预处理,并将其存为.tif格式；第2步,根据ArcGIS窗口分析法提取地势起伏度的基本原理和提取步骤编写代码；第3步,应用ArcGIS的Python模块自动提取网格2×2、3×3、4×4、5×5、…、1 000× 1 000下的地势起伏度；第4步,导出数据,整理成表格。

2.3 最佳统计单元计算

2.3.1 计算原理 确定一个最佳分析区域是地势起伏度提取算法中的核心步骤,也是决定区域地势起伏度提取效果和有效性的关键。对于相同区域的DEM,当分析区域由小变大:区域内最高点与最低点的高差,无论何处都是从小变大,一般情况下,高差从开始时以较快的速度增加,以后增加的速度变缓,在分析区域面积达到某一个阈值之后,高差变化基本稳定在一个数值上。使高差变化率由快变缓的节点即为拐点,这个“拐点”所对应的面积即为最终的统计单元面积,而“拐点”所对应的高差值即为研究区的地势起伏度[14-16]。最佳统计单元的确定属于典型的变点分析问题,均值变点分析法对有一个变点的检验最为有效,在提取最佳统计单元时更为合适[6-8]。均值变点分析法主要方法[17]如下:

设{xt,t=1,2,…,N}为非线性系统的输出,其系统模型或输出序列在某未知时刻起了突然变化,该时刻即为变点,变点的存在会使原样本的统计量S与样本分段后的统计量Si之间的差距增大。

1)令i=2,…,N,对每个i将样本分为2段: x1,x2,…,xi-1和xi,xi+1,…,xN,计算每段样本的算术平均值和,以及统计量

计算期望值

2.3.2 S和Si值的自动计算步骤 根据Python模块的基本功能,利用Python编程来计算最佳统计单元的值。第1步,对平均地势起伏度数据Ri进行处理得到单位面积数据T,再对该数据取对数lnT即得到数据序列{xi}；第2步,将数据序列{xi}以.txt格式保存；第3步,根据均值变点分析法的计算原理和公式编写代码,应用ArcGIS的Python模块对S和Si值进行批量自动计算；第4步,S和Si值的整理与分析。

3 结果与分析

3.1 地势起伏度

以ArcGIS为平台,采用窗口分析法,用Python模块编程自动提取不同网格单元对应下的地势起伏度,导出数据并整理,以网格单元面积为横坐标,相应的平均地势起伏度计算结果为纵坐标,绘制平均地势起伏度与网格单元面积对应散点分析图(图2)。从变化趋势来看,望谟河流域平均地势起伏度随网格单元面积的增大呈增大的趋势,对应散点分析图呈S型分布。在统计单元面积为200万m2之前,研究区内平均地势起伏度随着统计单元面积的增大而增大,且变化幅度较大；在1～200万m2之间,曲线上出现了平均地势起伏度增加由陡变缓的那一点,在增速变缓的过程中存在一处明显拐点使增大趋势陡变缓,在200万m2之后平均地势起伏度增加趋势变缓,并近似趋于平稳。平均地势起伏度与网格单元面积之间存在显著的对数变化关系,由陡变缓的这一点所对应的面积即最佳统计单元的面积,为验证是否可以使用均值变点分析法寻找拐点,对该散点图以Logarithmic对数方程y=B0lnx+B1进行拟合。得出拟合曲线方程y=85.001lnx+209.53时,R2=0.974 5,经拟合优度检验,拟合度好,达到显著水平,通过统计学检验,即该流域存在地势起伏度最佳统计单元,并可以采用均值变点分析法寻找拐点。

图2 平均地势起伏度与网格单元面积对应关系拟合曲线Fig.2 Fitting curve of the relationship between gird units and average relief amplitude

3.2 最佳统计单元

均值变点分析法可以准确计算拐点的位置,但是计算量大,工作复杂,容易出错。采用Python模块编程可以自动计算Si和S值,准确而且高效。通过计算整理可以得到的统计量S=2 619.636 466 95,整理可得Si值与对应点和对应网格之间的关系可以得到均值变点分析的统计结果,见表1。

Si值与对应点和对应网格之间的关系可知存在一个极值点使得Si值先随网格单元面积面积的增大而减小,之后Si值随网格单元面积的增大而增大。S与Si差值的变化拟合曲线见图3,明显看出S与Si差值的变化拟合曲线是单峰曲线,即存在一个点使得S与Si的差值先增大后减小,并逐渐趋近于0,这个点即为拐点,位于i=200～300之间。通过计算,在i=244即第243个点时,S与Si的差值达到最大,而第243个点对应的分析窗口为244×244,面积为148.84万m2；因此,望谟河流域提取地势起伏度的拐点位于网格244×244处,本次实验提取地势起伏度的最佳统计单元面积为148.84万m2。

表1 均值变点分析的统计结果Tab.1 Statistical result of mean change point analysis

图3 S和Si差值的变化曲线Fig.3 Change curve of the difference between S and Si

3.3 望谟河流域地势起伏度特征

根据计算分析选取网格244×244,面积1.49 km2分别为望谟河流域地势起伏度最佳分析窗口、最佳统计单元面积,通过对研究区的DEM处理得到研究区地势起伏度分布范围在0～648.23 m之间,地势起伏度在 200 m以上的地区占总面积的90.24%。中部河网密集区起伏度大,多在200 m以上；起伏度小于200 m的地区集中在流域边缘较为平坦的望谟县主城区和流域出口附近。为了探讨本文所用方法提取的地势起伏度的准确性,进一步分析望谟河流域的地势起伏度特征,分别采用我国数字地貌制图规划标准和刘新华等[4]在对比中国地势起伏度和中国地貌、中国土壤侵蚀特征时采用的分级标准,对望谟河流域地势起伏度进行分级并作特征分析。

首先根据我国数字地貌制图规划标准[18]:平原(＜30 m)、台地(30～70 m)、丘陵(70～200 m)、小起伏山地(200～500 m)、中起伏山地(500～1 000 m)、大起伏山地(1 000～2 500 m)、极大起伏山地(＞2 500 m)对流域地势起伏度进行分级。根据实际情况,将望谟河流域地势起伏度分为＜200 m、200～500 m、500～1 000 m 3级,分级结果见图4 (a)。其中地势起伏度在200 m以下的地区所占面积最小,只有1.33%；研究区主要地势起伏度分布在200～500 m之间,占总面积的88.88%；地势起伏度大于500 m以上的面积占9.7%。流域属于小起伏山地地貌,与实际情况相向吻合。

采用刘新华等[4]在对比中国地势起伏度和中国地貌、中国土壤侵蚀特征时采用的分级标准将流域地势起伏度划分为:＜200 m、200～400 m、400～ 600 m、＞600 m 4个等级。分级结果见图4(b),研究区在200～400 m和400～600 m起伏度上的分布最广,分别占总面积的45.78%和52.43%。对比刘新华等[4]在对比中国地势起伏度和中国地貌、中国土壤侵蚀时的结论:地势起伏度分布在200～400 m地区地貌为切割高地；400～600 m的地区地貌为切割山地,土壤侵蚀均较严重。研究区的地区地貌主要属于切割高地和切割山地,土壤侵蚀较严重,与当地实际情况一致,也符合刘新华等[4]将中国地势起伏度应用于水土流失评价中的结论。

图4 望谟河流域地势起伏度分级图Fig.4 Classification of relief amplitude in the Wangmo River Watershed

为进一步分析地势起伏度和流域地貌的关系,采用ArcSWAT模型的子流域划分功能,以最佳统计单元面积1.49 km2为阈值将望谟河流域划分55个子流域,统计分析地势起伏度与平均坡度、河道比降等因子的关系。流域中部地势起伏度大河网密集,河道主要分布在地势起伏度为200～500 m的区间内,占总河道的90.66%。其中,分布在地势起伏为300～400 m区间内的河道最多,占总河道的42.40%,分布在400～500 m区间上的河道占总河道的33.63%。流域内最大河道比降为14.5%,河道比降随地势起伏度的增加而增加,与地势起伏度之间有显著的正相关关系(P=0.000 1),子流域平均坡度也随地势起伏度的增加而增加与地势起伏度之间呈显著的正相关关系(P=0.000 1)。

海拔落差大、地势起伏度大、河道短等综合原因造成了望谟河河道比降大,水流势能强冲刷力强等特点；望谟河借助流域内的地势起伏度和海拔落差带来的高势能、高动能和极强的冲刷力一路南下,深切河道,河流强烈的侵蚀作用和频发的山洪地质灾害也影响着研究区的地势起伏度。在地势起伏度、地质构造、强烈的侵蚀作用以及频发的山洪地质灾害等原因的共同作用下最终形成了境内北高南低,沟谷纵横,山高谷深、地形破碎的独特的地貌。

4 结论

1)以望谟河流域1∶1万的地形图制作的分辨率为5 m的DEM为基础数据,采用窗口分析法,自动提取网格2×2、3×3、4×4、5×5、…、1 000×1 000下的平均地势起伏度,平均地势起伏度与网格单元面积之间的对应关系呈显著的对数关系。

2)以ArcGIS的Python模块为平台,采用均值变点分析法,利用编程实现了对望谟河流域地势起伏度最佳统计单元的自动计算,望谟河流域地势起伏度最佳统计单元面积为148.84万m2。

3)研究区内海拔落差大,地势起伏度在0～648.23 m之间,地势起伏度在200 m以上的面积占总面积的90.24%；起伏度小于200 m的地区集中在流域边缘较为平坦的望谟县主城区和流域出口附近。地貌属于小起伏山地地貌,以切割高地和切割山地为主,土壤侵蚀较严重。

4)研究区中部地势起伏度大河网密集,河道主要分布在地势起伏度为200～500 m的区间内,占总河道的90.66%；河道比降随地势起伏度的增加而增加,与地势起伏度之间有显著正相关关系(P= 0.0001),子流域平均坡度也随地势起伏度的增加而增加与地势起伏度之间呈显著正相关关系(P= 0.0001)；地势起伏度、地质构造、强烈的侵蚀作用以及频发的山洪地质灾害等原因形成了北高南低,沟谷纵横,山高谷深、地形破碎的独特地貌,又因历史、地理、人文等因数,流域出口和下游等地势较为平坦的地区成为了人口密集的主城区,所以当地频发的山洪地质灾害时刻威胁着人民的生活和安全、经济和社会的发展。从地形因子角度研究望谟河流域地质灾害频发的原因,对当地的灾害危险度进行区划、对地质灾害进行科学评价将是下一步的研究重点。

[1] 范娟,杨武年,南聪强,等.基于GIS的彭州市崩塌滑坡地质灾害地形因素分析[J].地理空间信息,2012, 10(6):19- 21

[2] 毕晓玲,李小娟,胡卓玮,等.地形起伏度提取及其在区域滑坡灾害评价中的应用:以四川省为例[J].安徽农业科学,2011,39(4):2413- 2416

[3] 郭芳芳,杨农,孟晖,等.地形起伏度和坡度分析在区域滑坡灾害评价中的应用[J].中国地质,2008,35 (1):131- 143

[4] 刘新华,杨勤科,汤国安.中国地形起伏度的提取及在水土流失定量评价中的应用[J].水土保持通报, 2001,21(1):57- 59

[5] 赵卫权,苏维词,袁俊.基于地形起伏度的贵州省景观空间格局分异特征[J].水土保持研究,2010,17(2): 105- 110

[6] 张磊.基于地形起伏度的地貌形态划分研究[D].石家庄:河北师范大学,2009

[7] 王志恒,胡卓玮,赵文吉,等.应用累积和分析算法的地形起伏度最佳统计单元确定[J].测绘科学,2014, 39(6):59- 64

[8] 涂汉明,刘振东.中国地势起伏度最佳统计单元的求证[J].湖北大学学报(自然科学版),1990,12(3):266-271

[9] 王康,何俊仕,于德浩,等.采用ArcGIS平台的地势起伏度自动提取技术研究[J].沈阳理工大学学报, 2013,32(2):63- 67

[10]朱红春,陈楠,刘海英,等.自1:10000比例尺DEM提取地形起伏度:以陕北黄土高原的实验为例[J].测绘科学,2005,30(4):86- 88

[11]刘振东,孙玉柱,涂汉明.利用DTM编制小比例尺地势起伏度图的初步研究[J].测绘学报,1990,19(1):57-62

[12]郎玲玲,程维明,朱启疆,等.多尺度DEM提取地势起伏度的对比分析:以福建低山丘陵区为例[J].地球信息科学,2007,9(6):1- 6

[13]张军,李晓东,陈春艳,等.新疆地势起伏度的分析研究[J].兰州大学学报(自然科学版),2008,44(增刊1):10- 13

[14]韩海辉,高婷,易欢,等.基于变点分析法提取地势起伏度:以青藏高原为例[J].地理科学,2012,32(1): 101- 104

[15]王玲,同小娟.基于变点分析的地形起伏度研究[J].地理与地理信息科学,2007,23(6):65- 67

[16]张锦明,游雄.地形起伏度最佳分析区域研究[J].测绘科学技术学报,2011,28(5):369- 373.

[17]项静恬,史久恩.非线性系统中数据处理的统计方法[M].北京:科学出版社,1997:3- 6

[18]中国科学院地貌图集编辑委员会.中华人民共和国地貌图集(1∶100万)[M].北京:科学出版社,2009

(责任编辑:程 云)

Analysis of ArcGIS-based relief amplitude of the Wangmo River Watershed in Guizhou

Gao Qing1,2,Tang Lixia1,Gu Xiaoping3,Yu Fei3,Cheng Xinglin1,Huang Song1
(1.College of Forestry,Guizhou University,550002,Guiyang,China；2.Yumo Government of Bijie City,551804,Jinsha, Guizhou,China；3.Guizhou Institute of Mountainous Environment and Climate,550002,Guiyang,China)

Studies of relief amplitude are of significance and practical application values not only to the sensitivity analysis of geological hazards like landslide and debris,but also to evaluation of the risk of geological hazards.However,little is known about relief amplitude in medium and small scale watersheds.Our objectives were to characterize the best relief amplitude of optimum statistical unit area of medium and small scale watersheds,learn geomorphic features of watersheds of such scales and figure out the reason of geological hazards.We chose the Wangmo River Watershed located in southwestern area of Guizhou Province and famous for its frequent landslides as our study area.We used a geographic information system(ArcGIS)and the DEM of 1:10000 scale to compute the optimum statistical relief amplitude unit of the watershed.Firstly,we adopted window analysis method and python module programming based on ArcGIS to automatically extract the different areas corresponding to the relief；secondly,we employed mean change point analysis method and python module programming to automatically calculate the optimum statistical relief amplitude unit of the watershed.Finally,we analyzed the relief amplitude characteristics of the watershed and the relationships between reliefamplitude,average slope and stream gradient factors.The results revealed that:1)There was a significant logarithmic relationship between relief amplitude and average cell area of the watershed.2) Based on python module programming we can automatically calculate the statistically optimum relief amplitude unit,that is 1.49 km2,and the relief amplitude ranged between 0-648.23 m.3)The relief amplitude is higher in the middle of the areas where river networks are concentrated,while it is lower at the edge of the plateau,and the rivers are mainly distributed in the range of 200-500 m of the terrain. The landscape belongs to small relief mountain.4)There is a significantly positive correlation between stream gradient and relief amplitude,that is,the stream gradient increased with the increase of the relief amplitude.Theseresultsshow thattheWangmoRiverWatershed belongstosmallmountain geomorphology,with cutting highlands and cutting mountains as well as serious soil erosion.

relief amplitude；Wangmo River Watershed；geological hazards；ArcGIS

P951

A

1672-3007(2015)04-0009-06

2014- 10- 30

2015- 04- 27

项目名称:贵州省科学技术基金“喀斯特流域气候/土地利用变化下的水文特征研究”(黔科合J字[2010]2240)；贵州大学引进人才科研项目“基于SWAT模型分析典型喀斯特流域的水文特征”(贵大人基合字(2009)006)；中国气象局第三次青藏高原大气科学试验

高蜻(1989—),女,硕士研究生。主要研究方向:流域综合治理。E-mail:xiuceqing@sina.com

†通信作者简介:唐丽霞(1976—),女,副教授。主要研究方向:地理信息系统在流域水土流失治理中的应用。E-mail:xialitang123@163.com