黄土高原典型流域水蚀因子多重分形特征

李莹莹,魏舟,,童亚文,张青峰

(1.西北农林科技大学资源环境学院,陕西杨凌 712100;2.陕西地矿区研院有限公司,陕西咸阳 712099)

水土流失是一个多重尺度地理界面过程。水蚀因子既是研究水土流失的基本内容,也是水土流失基本影响因子和水土流失评价与预报模型的重要参数。综合来看,水蚀因子大致可以分成几何参数(如高程、坡度、坡向、汇流量、河网数量、土壤孔隙度)和形态参数(如土地利用、植被覆盖)两大类,且相互作用和影响,表现出多尺度性、自相似性、时空耦合性等特征。其中,以高程、坡度、土地利用、植被覆盖为代表的水蚀因子尤为突出。深入量化水蚀因子对于揭示水土界面过程和水蚀机理具有重要意义,而如何定量描述这一非线性动态过程成为重要的科学问题。

分形理论和方法是探究不规则空间的有效手段[1-3]。大量研究表明,简单分形维数只能描述与表征水蚀因子的整体性[4-5],而多重分形作为描述非线性和复杂性事物的有力工具恰好可以探明水蚀因子空间局部特征[6-8],且主要有盒维数法、固定质量法、固定半径法、推广G-P法4种计算方法。近年来,诸多学者运用多重分形模型对地表剖线[9]、高程[10-11]、沟沿线[12]、坡度[13]、沟谷网络等[14]表征地貌形态以及土地利用[15]、植被覆盖等[16]地表覆被的单一因子展开了定量化研究,均侧重于地表形态特征,而对侵蚀发育过程中各参数,特别是针对黄土高原不同地貌类型区小流域的水蚀因子多重分形特征的综合研究依然非常薄弱。深入研究黄土高原典型流域的不同水蚀因子多重分形特征在一定程度上可以让人们深入认识和理解多尺度水蚀现象,进而促进多尺度水蚀模型的构建。

因此,本文分别以黄土高原沟壑区(王东沟流域与泥河沟流域)和黄土丘陵沟壑区(泉家沟流域与纸坊沟流域)的两个典型小流域为例,选取最具代表性的高程、坡度、土地利用、植被覆盖度4个水蚀因子构建多重分形谱,探讨黄土高原典型小流域水蚀因子多重分形特征,以期为黄土高原区域水文过程模拟与水土流失预测模型提供科学依据。

1 研究区概况

黄土高原地貌组合可分为高原沟壑区、丘陵沟壑区两大类型,前者由黄土塬和沟谷组成,后者主要由梁、峁和沟谷组成。高原沟壑区是黄土地貌中相对完整的地貌形态,坡度小,地势较为平坦,水土流失轻微,为重要的农业基地;丘陵沟壑区沟壑纵横,地形十分破碎,15°以上的坡面面积占50%～70%,是中国水土流失最严重的区域。两个地貌类型区的小流域基本情况见表1。

表1 研究区基本概况Table 1 Basic overview of the study area

2 材料与方法

2.1 数据来源及预处理

空间分辨率5 m 的DEM 数据(2018 年),数据来自INTERMap公司的NEXTMap数据集(https:∥www.intermap.com/nextmap);空间分辨率10 m的土地利用数据(2017年),来自FROM_GLC 数据,该数据集的总体精度为72.76%(http:∥data.ess.tsinghua.edu.cn/fromglc10_2017v01.html);空间分辨率10 m 的植被覆盖度数据(2020 年)运用GEE(Google Earth Engine)通过像元二分法从Sentinel-2数据中提取[17]。

所有数据投影坐标均设置成WGS_1984_UTM_Zone_48N,分别将数据重采样成10 m×10 m的栅格数据。

2.2 多重分形参数及其计算

盒维数法常用于计算地形地貌的多重分形,如图1所示。本文采用盒维数法对4个水蚀因子分别进行多重分形计算[18]。

图1 盒维数法示意图Fig.1 Schematic diagram of box dimension method

(1)首先,将研究对象分成一系列的小盒子,小盒子的边长记为ε,定义第i个盒子的密度分布函数Pi为:

式中:ε为一个测量单位,～表示两者存在线性关系。非整数α一般称为奇异指数,指水蚀因子分布概率Pi随ε变化而呈现一致标度关系的各个栅格点子集的性质,其取值和盒子的位置有关。α越大,则该子集的分布概率越小,反之亦然。

(2)其次,将分形上具有相同α值的小盒子数目记为Nα(ε),它是与ε大小有关的,并可写成:

式中:f(α)为具有相同α的值的子集的分形维数,反映不同α下的分形特征。

(3)运用矩表示方法,给定不同的q值对各个盒子的概率进行q次加权求和,即配分函数Xq(ε):

式中:τ(q)为质量指数,多重分形谱成立的前提条件是,在无标度区间内,lnXq(ε)～lnε曲线拟合结果呈线性关系;q∈=(-∞,+∞)。

(4)定义广义分形维数(Dq)。

(5)对τ(q)～q函数进行勒让德变换(Leg-endre transformation),最终得到f(α)～α曲线,即多重分形谱:

多重分形谱中特征参数的含义:f(α)max为水蚀因子的简单维数D,它是水蚀因子空间特征的总体表征,值愈大,说明水蚀因子整体空间结构愈不均匀。αmax和αmin分别为水蚀因子概率子集奇异指数的最大值和最小值;Δα为奇异指数变化范围,定量地反映水蚀因子分布概率的复杂程度,其数值等于αmin与αmax的差值,数值的大小与水蚀因子空间内部的复杂度成正比,值愈大,水蚀因子空间结构愈复杂,内部起伏程度大;f(αmax)和f(αmin)分别表示水蚀因子分布最大、最小概率子集的分形维数,分别表明水蚀因子最大、最小单元子集在流域中所占比例大小;Δf(α)表示最大水蚀因子分布概率与最小水蚀因子分布概率分形维数的差值,间接反映了流域水蚀因子最大奇异指数与最小奇异指数所代表子集数目的比例,其值越接近0,多重分形谱越对称,Δf(α)＜0,在多重分形谱曲线上表现为右钩状,反之,则表现为左钩状。

使用ArcGIS中等间距分级法,将研究区的高程间隔100 m 共分成8类,坡度间隔10°共分成8类(表2),此外,研究区的土地利用主要有7大类别,植被覆盖度参考论文划分成无覆盖区(FVC=0)、低覆盖区(0＜FVC≤0.1)、较低覆盖区(0.1＜FCV≤0.3)、中等覆盖区(0.3＜FVC≤0.5)、较高覆盖区(0.5＜FVC≤0.7)、高覆盖区(FVC＞0.7)6类等级[19]。

表2 水蚀因子的分级Table 2 Classification of water erosion factors

3 结果与分析

3.1 水蚀因子概率分布特征

从图2可以看出,王东沟的高程栅格集中分布在1 000～1 233 m范围内,占比最高的是1 200～1 233 m,占比36.914%。而泥河沟的高程栅格大部分分布在800～1 100 m,占比最高的介于900～1 000 m,占比达32.275%,两个流域的高程各类别栅格占比相差较大。主要是因为王东沟存在部分的梁状丘陵,是山地与黄土塬之间的过渡陆地带,海拔较高。而在坡度因子中,两个小流域的坡度栅格占比特征具有高度的相似性,王东沟和泥河沟的坡度栅格大部分处于0°～40°中小坡度区间,栅格占比＞95%,其中栅格占比最高的均是0°～10°区间。在土地利用因子中,王东沟流域缺少水体用地类型,因而泥河沟的土地利用类型更为丰富。王东沟的土地利用类型栅格集中分布在耕地和林地,而泥河沟的土地利用类型栅格集中分布在耕地、林地和草地。两个流域的土地利用类型栅格占比最高的均是耕地,但泥河沟的草地栅格占比远超王东沟。在植被覆盖因子中,王东沟和泥河沟植被覆盖度的栅格均主要集中在0.1＜FVC≤1 范围内,其中王东沟植被覆盖度栅格占比最高的是0.7＜FVC≤1区间,而泥河沟栅格占比最高的区间是0.3＜FVC≤0.5。整体上,王东沟的植被覆盖比泥河沟好,土壤保持服务更强。

图2 水蚀因子分级统计Fig.2 Classification statistics of water erosion factors

在黄土丘陵沟壑地貌类型区中两个小流域的4个水蚀因子分布特点如图2所示,在高程因子中,泉家沟的高程栅格占比最高为75.367%,位于1 000～1 100 m 范围内。而纸坊沟的高程栅格主要集中在1 100～1 300 m,占比最高达42.214%,介于1 200～1 300 m。在坡度因子中,泉家沟和纸坊沟的坡度栅格均集中分布在10°～30°范围内,栅格占比＞95%,其中栅格占比最高的均是20°～30°区间,且泉家沟在10°～30°范围内栅格占比均比纸坊沟大。在土地利用因子中,泉家沟的土地利用类型栅格主要是裸地,占比总面积55.759%。而纸坊沟的土地利用类型栅格集中分布在耕地与林地,分别占比40.949%,44.337%。在植被覆盖因子中,泉家沟和纸坊沟植被覆盖度栅格均主要集中在0.1＜FVC≤1 范围内,占比面积超过90%。其中,泉家沟植被覆盖度栅格占比最高的在0.7＜FVC≤1区间,而纸坊沟栅格占比最高的是0.3＜FVC≤0.5区间。总体来看,泉家沟与纸坊沟的植被生态环境差异较小。

对比4个研究区的水蚀因子栅格统计结果:高程分布差异性较大,而植被覆盖分布差异性较小,整体植被生态环境良好;除泉家沟以外,研究区的土地利用空间格局具有相似性,均是耕地、林地、草地占区域的主要面积;高原沟壑区与丘陵沟壑区两种地貌类型的流域坡度分布差异性明显,前者地形特征较为均一,后者地形起伏较大,但在同一地貌区内,区域间坡度分布特征相似性极高。

3.2 水蚀因子多重分形谱及其参数

3.2.1 多重分形特征检验 图3仅展示了王东沟流域高程因子的配分函数与盒子大小关系曲线,即高程配分函数Xq(ε)～ε在双对数坐标下的关系曲线,当q取不同值时,配分函数与盒子尺度均呈线性关系,各水蚀因子均符合多重分形计算的条件。

图3 王东沟流域高程X q(ε)～ε 关系曲线Fig.3 X q(ε)～εcurve of elevation relation of Wangdonggou Watershed

如图4所示,王东沟、泥河沟、泉家沟和纸坊沟4个流域不同水蚀因子的多重分形谱均为凸函数,多重分形谱的函数值都是随着奇异指数的增加均呈现出先增加再减小的变化特征。虽然王东沟与泥河沟均位于黄土高原沟壑区,但是在同一水蚀因子的多重分形谱对比中,两个流域在多重分形谱对称性、图像分布范围均有差异。此外,同样位于黄土丘陵沟壑区的泉家沟和纸坊沟流域的水蚀因子多重分形谱亦是如此。多重分形谱参数信息统计结果见表3。

图4 水蚀因子多重分形谱Fig.4 Multifractal spectrum of water erosion factors

表3 水蚀因子多重分形谱的重要参数Table 3 Important parameters of multifractal spectrum of water erosion factor

3.2.2 多重分形的特征 (1)黄土高原沟壑区多重分形特征。从图4与表3中可以看出,黄土高原沟壑区水蚀因子的差异性:王东沟的高程因子和坡度因子的简单维数f(α)max均大于泥河沟的,整体而言泥河沟地貌形态较王东沟更复杂,此时的简单维数不能够完全刻画高程、坡度内部的精细结构,需参考多重分形的其他参数。王东沟的高程与坡度的奇异指数的跨度Δα更大,其次,王东沟的高程最大奇异指数比泥河沟的大,而最小奇异指数比泥河沟的小,说明在这两个研究区中,泥河沟最小高程的栅格数目所占比例更大,而王东沟最大高程栅格概率占比更大,王东沟的高程与坡度内部分布更不规则;泥河沟高程和坡度的Δf(α)＞0,而王东沟高程的Δf(α)＜0,其坡度的Δf(α)＞0,表明泥河沟高程与坡度的多重分形谱均呈右钩状,则泥河沟较大高程与坡度区域对空间变异性的贡献率较大,而王东沟坡度的多重分形谱为左钩状,其空间变异性更多地依赖于高程较低、坡度较大的栅格值,即王东沟地貌形态相对起伏度较大。在土地利用、植被覆盖因子中,泥河沟的简单维数皆大于王东沟,说明从整体上看相比于王东沟,泥河沟的土地利用与植被覆盖的空间结构更复杂。王东沟比泥河沟流域的奇异指数跨度大,表明王东沟地块结构间破碎化更严重;王东沟土地利用的Δf(α)＜0且接近于0,其植被覆盖的Δf(α)＞0,王东沟这两个因子的Δf(α)均小于泥河沟的Δf(α),说明王东沟的土地利用多重分形谱对称性更好,王东沟的植被覆盖的多重分形谱呈左钩状,与泥河沟的两个因子多重分形谱的形状一致,泥河沟裸地、高植被覆盖度概率较大占主导地位,王东沟的耕地与高植被数据贡献率更大,表明泥河沟的高植被覆盖分布区域较均匀,其土地利用类型以及植被生态环境空间异质性较小。

(2)黄土丘陵沟壑区多重分形特征。由图4和表3可知,黄土丘陵沟壑区水蚀因子的差异性:在高程因子和坡度因子中,泉家沟与纸坊沟流域高程的简单维数f(α)max分别为1.940,1.917,从整体上来看,显然泉家沟地貌形态较纸坊沟更复杂。泉家沟与纸坊沟坡度的简单维数分别是1.984,1.982,两者几乎相等,仅仅利用简单维数表征两个流域的地貌特征十分接近,需寻求更为敏感的特征指标对流域进行描述。高程奇异指数的变化范围更大的是泉家沟,纸坊沟的坡度奇异指数的跨度较大,此外,泉家沟与纸坊沟高程的Δf(α)＜0,且参数值相差甚小,这两流域坡度的Δf(α)＞0,纸坊沟的此参数值大于泉家沟,这一系列表明纸坊沟的高程内部差异较大,空间结构不均匀程度更高,纸坊沟地貌表面峰谷数目的比例更大,即地貌相对起伏度更大。在土地利用因子、植被覆盖因子中,土地利用的简单维数较大值是泉家沟,而植被覆盖的简单维数较大值是纸坊沟,从整体上看相比于纸坊沟,泉家沟的土地利用类型空间结构更为不均,但其植被环境内部较为均一。从其他参数对局部特征进行对比分析,泉家沟的最大奇异指数以及奇异指数的跨度Δf(α)均比纸坊沟流域大,同时,泉家沟的Δf(α)均大于纸坊沟,多重分形谱都是呈左钩状,反映了泉家沟的裸地占主导地位,概率分布较大,植被生态环境区域差异大,地表覆被的破碎化程度较严重。

总体而言,在两种地貌类型区中,黄土高原沟壑区流域高程因子的奇异指数跨度Δα均小于黄土丘陵沟壑区的,且丘陵沟壑区的流域Δf(α)＜0,表明丘陵沟壑地貌区域高程分布复杂度更大,峰谷数目所占比例较大,高程起伏变化大,其地貌形态更尖锐;坡度中高原沟壑区Δα与Δf(α)皆大于丘陵沟壑区,表征高原沟壑区的坡度空间分布状况更不规则;土地利用中,高原沟壑区内部两个流域Δα与Δf(α)的差异较大,其土地利用空间格局更不均匀;植被覆盖中,丘陵沟壑区植被覆盖空间变异程度更高。

4 讨论

坡谱信息熵是地貌变异分析的一个重要指标[20],从另一个层面揭示地面的破碎程度,进一步阐述坡度空间变异情况,对地形特征进行补充说明。其次,数据的分辨率对多重分形谱的构建结果具有很大的影响[21],探索多重分形参数随着数据分辨率变化而变化的规律对准确构建多重分形谱具有指导意义。

4.1 黄土沟壑区坡度的地面坡谱特征

张勇[22]提出地面坡谱的理论概念并阐述了不同地貌类型区具有不同的坡谱曲线特征。本研究在坡度概率分布统计过程中,王东沟和泥河沟流域的坡度概率特征曲线即坡谱与王春[23]提取的黄土高原沟壑区坡谱曲线基本吻合,变化趋势是先减少再增加最后又减少;泉家沟和纸坊沟的坡谱同赵牡丹等[24]的研究结果一致,各级面积占比均先增大后减小。李发源[20]、刘双琳[25]等认为坡度信息熵在一定程度上反映了地形起伏的强弱和地表复杂度的大小,进而利用坡度信息熵对研究区坡度概率曲线特征进行补充说明(图5),王东沟、泥河沟、泉家沟与纸坊沟的总信息熵分别是2.544,2.339,2.517,2.553,丘陵沟壑区的两个流域坡度信息熵大于高原沟壑区的熵值,与李发源等[20]研究结果相吻合,说明丘陵沟壑区的两个小流域地形起伏更大,且地形起伏最大的是纸坊沟,此外,在高原沟壑区中,王东沟地貌形态较复杂。

图5 坡谱信息量Fig.5 Slope spectrum information graph

4.2 不同分辨率下的水蚀因子特征

陈旺等[11]计算了不同分辨率DEM 的多重分形谱,揭示随着DEM 分辨率的降低,喀斯特地貌等微地形的地貌形态进一步被概括。利用ArcGIS 将水蚀因子的栅格数据重采样成10,15,30 m,随着栅格数据分辨率的降低,水蚀因子的多重分形参数均发生了不同程度的变化[26]。从表4可知,4个流域的简单维数f(α)max均在逐渐减小,表明整体上地貌形态和地表覆盖状况信息受到概括。4个研究区高程的Δα增加,说明高程信息概括随着DEM 分辨率的降低而越受到限制;王东沟和纸坊沟坡度与土地利用的Δα增加,表明坡度和土地利用信息概括随着DEM 分辨率的降低而越受到限制,而泥河沟和纸坊沟坡度与土地利用的Δα先减小后增加,则坡度和土地利用信息概括先细致而后受到限制;王东沟与泥河沟植被覆盖的Δα先减少后增大,则植被覆盖信息概括随着DEM 分辨率的降低先细致而后受到限制,而泥河沟与纸坊沟植被覆盖的Δα一直在减小,说明植被覆盖信息概括导致复杂度升高。多重分形谱的特征与数据的分辨率息息相关,选择恰当的分辨率数据对准确获取对应因子的多重分形谱具有重要意义。此外,盒子尺寸、q值的选取最终都会影响多重分形结果。当栅格数据分辨率发生变化后,盒子尺寸以及q值是否也需调整变动还有待进一步研究。

表4 不同分辨率下水蚀因子多重分形谱的参数Table 4 Parameters of multifractal spectrum of water erosion factors at different resolutions

5 结论

(1)在不同地貌区与同一地貌区内的流域,高程分布差异性均较大;不同地貌类型区的流域坡度空间分异显著,黄土高原沟壑区流域的坡度主要分布较小坡度区域,而中等坡度范围在黄土丘陵沟壑区坡度占比最多;同一地貌类型区中,流域的坡度分布特征均具高度相似性,即栅格占比特征趋势线相似。除泉家沟以外,其余3个研究区的土地利用类型占比均是耕地、林地、草地占区域的主要面积;研究区之间植被覆盖分布差异性较小,整体植被生态环境良好。

(2)地貌相对起伏最大和土地利用空间分布最不均匀的流域均是王东沟流域,其坡度与土地利用的Δα值分别是0.817,0.310,植被覆盖区域差异最大的是泉家沟流域,其植被覆盖的Δα=1.309。

(3)黄土高原沟壑区较丘陵沟壑区地貌表面更为圆润、坡度空间分布不规则程度更高;黄土丘陵沟壑区的高程较黄土高原沟壑区的相对起伏更大,土地利用破碎化更严重,丘陵沟壑区的地表植被覆盖区域之间复杂程度较高,内部植被生态环境波动较大。

本研究结果与实际情况有着较好的对应,表明多重分形理论与方法对区域开展水蚀因子量化研究是可行的,对多尺度水蚀模型的构建、水文过程模拟以及水土流失预报等具有积极意义。