基于Hkr值的北天山冰川侵蚀空间分布特征及其主控因素研究

姚 盼，王 杰，林文旺，曾兰华，陈仁容

（1.嘉应学院地理科学与旅游学院，广东梅州 514015；2.兰州大学资源环境学院西部环境教育部重点实验室，甘肃兰州 730000）

0 引言

冰川侵蚀可以直接塑造和影响地形，是理解晚第四纪气候变冷在山区地貌演化中作用的关键环节。然而，当下有关河流侵蚀对地形影响的研究较多［1-7］，而冰川侵蚀对地形影响的研究却相对较少［8-12］，这主要是因为长尺度冰川侵蚀速率实测困难，与之密切相关的侵蚀过程和地形演化的研究难以开展，认识也较为有限，而基于一定的侵蚀指标去分析冰川侵蚀的控制因素则是了解其作用机制以及与构造、气候间相互关系的有效途径。前人的研究已经表明冰川的侵蚀受诸多因素的影响，如构造、冰川性质、气候、地形、冰量和冰川作用时间等［13-18］，但是究竟哪者发挥着主控作用则仍不明确。已有的研究显示降水在全球不同地区冰川侵蚀速率的变化中发挥着重要的作用［16，19］，但在中高纬度地区却是冰川底部热力性质（可用气温作为代用指标）控制着冰川的侵蚀速率［20］，在中低纬度的某些地区，主控因素则为地形或温度所决定的冰川规模或冰量［17，21］。而且即使主控因素一样，其影响冰川侵蚀的过程可能也不同［20-22］。目前，学者在研究冰川侵蚀的主控因素时，仅在某些构造活动略单一的区域开展，致使对构造的影响认识不足［17，20］。虽然有部分学者通过数值模拟的方式探讨过构造对冰川侵蚀的主控作用［21］，但是在实际情况中构造是否是冰川侵蚀的主控因素呢？又是通过怎样的过程施加影响？针对上述问题，在构造及其他影响因素更加多样化的区域探讨冰川侵蚀的主控因素也是非常必要的。

北天山是中亚地区深受第四纪冰川影响的区域，冰川侵蚀地貌较为发育。且受印度—亚欧板块的碰撞挤压，北天山的构造活动性自西向东逐渐减弱［23-25］，使其成为探讨不同构造抬升背景下冰川侵蚀主控因素的理想场所。鉴于此，本文采用Hkr方法，以北天山北麓7个现代冰川规模较小的冰川流域为研究对象，首先获得了北天山冰川侵蚀的分布特征，在此基础上分析造成各流域侵蚀差异的主控因素，以对冰川侵蚀的过程有更深入的认识。

1 研究区概况

本文研究区位于85°～88° E，42°50′～43°30′N的范围，按前人对我国天山的分类［26］，属于北天山东段，呈北西西走向，平均海拔约2 700 m左右，山顶基本在4 000 m以上，其上有现代冰川存在，北坡以现代冰川为源头还发育了一系列流向北东的河流（图1）。本文所选的冰川流域就是玛纳斯河至乌鲁木齐河之间北坡各河流的源头区域，自东向西，依次编号为1～7（图1）。

图1 北天山区域概况图Fig.1 Overview of the northern Tianshan Mountains

新近纪以来，印度与亚洲板块碰撞的远程效应使最早形成于晚古生代的天山造山带重新复活，形成规模宏大的年轻的板内再生造山带［27］。根据前人对北天山构造活动的分析，该区域总体上构造较为活跃，且抬升速率由西向东呈现逐渐降低的趋势［23-25，28］。本文所选取的流域中，图1中编号6和7构造略活跃，1～5构造活动相对弱一些。研究区各流域岩性多样，流域7大部分区域为砂岩、砾岩、灰岩和盐酸盐，小部分区域为花岗岩；流域6以粉砂岩、砂砾岩为主，局部区域为花岗岩或片麻岩；流域5主要为凝灰岩、硅质岩、砂砾岩夹中酸性火山岩；流域4大部分区域为粉砂岩、砂砾岩，局部区域为花岗岩；流域3以片麻岩、片岩和大理岩为主；流域2以粉砂岩、硅质岩、灰岩为主，局部为闪长岩；流域1主要为火山角砾岩、粉砂岩、凝灰岩、片麻岩、和花岗岩［29］。

北天山地区的气候主要受西风带控制，年降水量一般少于250 mm，受地形影响年降水量在中山区可达650 mm，降水主要集中在夏季；因受北方来的寒流影响，年平均气温为1.3～7.8℃［30］。本文所选7个冰川流域，降水量在300～500 mm之间，自西向东有减少的趋势；气温在-6～-9℃之间（图2）。

图2 北天山7条冰川流域多年平均气温和降水（气温和降水数据来源于World Clim 2.0 Beta version 1）Fig.2 Annual mean temperature and precipitation of seven glacier basins in the northern Tianshan Mountains（tempera⁃ture and precipitation from World Clim 2.0 Beta version 1）

受气候和地形的影响，北天山地区共发育了3 584条冰川，平均规模为0.55 km2，多分布在海拔3 700～4 300 m的区域［31］。本文所选的冰川流域中，除了编号7的现代冰川规模较大，其他流域现代冰川规模均较小，有些甚至快完全消退（图1）。虽然现代冰川规模较小，但是该区域却保存着丰富且典型的古冰川遗迹，末次冰期时的冰川末端可延伸至2 900 m左右［30，32］。天山及周边山地测年的结果显示，研究区内很可能经历了包括MIS6、4、3、2及小冰期和新冰期在内的多次冰川作用［33-42］。

2 方法和数据获取

在冰川侵蚀的定量化方法中，幂函数模型是较为成熟的方法，应用也比较广泛［12，43-47］，但是其在槽谷形态分析中具有复杂性和不确定性，即函数的参数既受拟合过程又受槽谷形成过程的各种因素影响［48］。HI指数是面积-高程分布曲线的积分值，经常用来反映地表被侵蚀的程度，其也曾被用于冰川流域，但有学者认为冰川源头的地形、物质平衡线位置，以及流域面积的选取等均可能造成HI值的变化，致使该指标难以完全指示冰川侵蚀的程度［49-50］。因此，本文采用了将面积-高程分布曲线与冰川槽谷的形态结合起来的Hkr方法来反映冰川侵蚀的程度，该方法相比幂函数模型更加便捷且能覆盖包括冰斗在内的整个冰川作用区，相比HI则更加准确地反映冰川侵蚀程度，在新西兰阿尔卑斯山和祁连山的应用均证实其可靠性［51-52］。

2.1 Hkr值的计算方法

冰川槽谷横剖面的演化研究显示，冰川侵蚀在槽谷底部略偏上的谷坡位置最强［53］，这会导致谷坡上该点的侵蚀量大于谷底，进而造成谷底坡度变缓［51］。基于槽谷形态的这种变化，Sternai等［51］将之转化为了面积-高程分布曲线的变化，并从中提取了能反映冰川侵蚀程度的指标Hkr值（图3）。该方法的原理为：河流“V”形谷时，谷地坡度为一定值Sf；当冰川占据谷地并开始侵蚀后，谷底会越来越缓，坡度小于Sf，而在谷坡上坡度等于Sf点的位置会越来越高，相应的面积-高程分布曲线上所对应的高度值也增加，这一高度值就是Hkr值。因此，Hkr值可作为衡量冰川侵蚀的指标，该值越大，冰川侵蚀越强。Hkr值的计算方法如下：

式中：h'是面积-高程分布曲线的标准化高程；S是面积-高程分布曲线上某点的坡度值；Sf则是河流作用区面积-高程分布曲线［图3（b）中直线］的平均坡度。因河谷与周边地形相比属于负地形，通常位于面积-高程分布曲线的下部，且谷地上部易受其他外力作用的影响，因此在计算“V”形谷坡度Sf时选用了面积-高程分布曲线下半部分的地形（即标准化高程在0～0.5之间的曲线）。

图3 冰川槽谷横剖面形态的演化及相应面积-高程分布曲线反映的Hkr值的变化（图改自参考文献［51］，A是区域总面积，a表示某一高程之下的面积；h表示某高程与最小高程之差，H为最大与最小高程之差）Fig.3 Evolution of glacial valley cross-profile（a）and change of Hkr values reflected by area-elevation distribution curves（b）（modified from reference［51］，a represents the total area of the region，a represents the area below a elevation；h is the difference between a elevation and the minimum elevation，and H is the difference between the maximum elevation and the minimum elevation）

在获取Hkr值时需要确定窗口大小，以避免计算面积不一对指标的影响，窗口的选取可通过下列变差函数实现：

式中：γ表示方差，k=1，2，…，n（栅格个数）；c=25（最小的窗口对应的栅格数目）；N为研究区内栅格的总个数。确定窗口大小的具体步骤为：依据变差函数得到γ随栅格数目变化的曲线；当γ值增大到一定程度时，出现一个拐点，此后无论窗口多大γ值都不再变化，从拐点出现时对应的栅格数目，就可以计算出窗口的大小［51］。

2.2 Hkr值的获取

本文所选用的影像为30 m分辨率的ASTER GDEM数据，随后在此基础上提取Hkr值。首先，从河流作用区的面积-高程分布曲线中提取坡度值Sf。因Hkr值需要通过对冰川与河流地貌进行对比得到，所以前提是明确北天山河流作用的界限。本文将末次冰期多条冰川末端平均海拔2 776 m当作河流作用区的上界（选末次冰期末端是因该时期冰川地貌形态相对完整，年代较为确定，且经历过末次冰期冰川作用的区域很有可能经历过更老冰期的冰川作用），山区河流出山口的平均海拔1 219 m作为下界来计算Sf值，得到了北天山河流作用区的面积-高程分布曲线（图4），计算a/A在0～0.5之间的标准化高程与面积之间的比值，得出Sf=0.99。然后，确定窗口大小。根据式（2）中的方法，计算得出Hkr值的方差随窗口大小变化的曲线（图5），得出窗口大小为100 pixel×100 pixel时较为合适，即实际的窗口约为3 km×3 km。最后，将整个研究区用相同大小的窗口进行分割，运用公式（1）计算每个窗口中的Hkr值，并得到每个冰川流域平均的Hkr值（图6）。由于运算过程复杂，数据均通过MATLAB编程来实现。

图5 Hkr值的方差随窗口大小变化的曲线Fig.5 Variation curve of Hkr value variance with window size

图6 北天山Hkr值的空间分布Fig.6 Spatial distribution of Hkr values in the northern Tianshan Mountains

2.3 影响冰川侵蚀程度的定量和定性因素提取

探讨冰川侵蚀的主控因素是理解冰川如何塑造地形的有效途径，目前已知冰川侵蚀的影响因素很多，包括定性（冰川性质、冰川作用时间、岩性）和定量（冰量、气候和地形等）的因素［14，17，20，21，54-56］，在揭示这其中何种因素为每个山地冰川侵蚀主控因素的研究中，则主要是通过各山地所选流域的冰川侵蚀程度与能定量的影响因素进行相关性分析，以及与无法定量的影响因素的定性分析来获取。本文在上述影响因子之外，将构造因素也纳入进来。气候方面选取气温和降水数据，来源于World Clim 2.0 Beta version 1（http：//worldclim.org）；地形方面采用流域的平均坡度、平均地形起伏、顶点高程和平均高程；冰川作用区岩性根据1：250万的地质图提取［29］；冰川规模则用冰川流域的面积代替［57-58］；因所选研究区冰川的性质和冰川作用时间差异不大［33-42，59］，这两个因子未做分析；构造缩短量则是自西向东逐渐减小［23-25，28］。综上，本文共确定了9个（7个定量和2个定性）影响北天山冰川侵蚀的因素，并提取了相关参数（表1）。此外，各流域冰川侵蚀的数据Hkr也呈现在表1中，因构造抬升速率的缺乏，其定性数据难以呈现在表1中。

表1 北天山7个冰川流域相关参数表Table 1 Related parameters of seven glacial basins in the northern Tianshan Mountains

3 结果和讨论

3.1 北天山冰川侵蚀程度的空间分布特征

北天山Hkr值的分布显示，较大的Hkr值通常出现在河流流经路线及其附近，如玛纳斯河的源头部分，符合一个流域中谷地部分侵蚀最强的特征（图6）。除乌鲁木齐河源外，其他流域末次冰期末端之上侵蚀一般较强，而向下则侵蚀变弱，可能是由于末次冰期末端之前的谷地经历了末次及其之前更早的冰期，而末端之下的谷地仅受到末次之前冰川的作用；乌鲁木齐河距离末次冰期末端11 km附近Hkr值偏大且为1，主要是由于该段河谷横剖面呈上凸状，且坡度小于45°，致使其面积-高程分布曲线上的高差与水平距离之比始终小于Sf，最终造成的结果就是整个窗口内的标准化面积-高程分布图中找不到坡度≥Sf的点，这时就会将Hkr赋值为最大的标准化高程1（图6）。此外，7个冰川流域平均Hkr值的比较结果表明：冰川侵蚀程度在7号流域（玛纳斯河源）最强，其后依次为5号流域（塔西河源）、6号流域（古仁河源）、3号流域（三屯河源）、4号流域（呼图壁河源）、2号流域（头屯河源）、1号流域（乌鲁木齐河源），侵蚀程度总体上呈现自西向东减小的趋势。

3.2 冰川侵蚀程度空间分布的影响因素分析

本文在获得冰川侵蚀分布规律的基础上，进一步分析了造成北天山北坡冰川侵蚀差异的主控因素，首先对7个定量因素进行分析，然后分析2个定性因素的影响。定量分析时为了与青藏高原其他地方进行对比，采用了皮尔逊相关性分析方法。

北天山各冰川流域侵蚀程度及其影响因素间的相关分析表明（表2）：冰川侵蚀的差异与平均降水和顶点高程的变化均存在很好的相关性，系数高达0.799和0.891，与平均地形起伏、流域面积之间也存在较好的相关性，相关系数分别为0.558和0.684，而与其他因子之间的相关性则较弱。此外，各因子交互分析的结果显示两因子共同作用后并未增强对Hkr值的解释（表3）。因此，在所有的定量因素中，降水和顶点高程相较于其他因素对北天山冰川侵蚀程度的影响可能更大，而其中顶点高程为主控因素。

表2 北天山7个冰川流域Hkr与其可能控制因素间的皮尔逊相关性（|R|）矩阵Table 2 Pearson correlation（|R|）matrix for all the potential controlling variables and Hkr of seven glacial basins in the northern Tianshan Mountains

表3 北天山各定量因素间的交互作用Table 3 The interaction effect of quantitative factors in the northern Tianshan Mountains

岩性和构造也是影响冰川作用区侵蚀程度的重要因素［13，54，61-63］。岩性方面主要考虑的是基岩的抗侵蚀能力，抗侵蚀能力越强的基岩越不易被侵蚀，反之则强。一般抗侵蚀能力的强弱按岩石类型排序为：片麻岩或花岗岩>片岩或板岩>砂岩>泥岩或黏土［61，64］。北天山3号流域的基岩抗侵蚀能力较强，其余流域基岩抗侵蚀能力稍弱（表1）。Hkr值与岩性对比结果显示：3号流域基岩抗侵蚀能力强，其Hkr小于5、6、7号流域，却大于1、2、4号流域，并未显示出岩性对冰川侵蚀的主控作用；除3号流域之外的其他流域基岩类型相差不显著，然而Hkr值却存在较大差异，尤其是1、2、4号流域的Hkr较小，而5、6、7号流域的Hkr却很大，这可能源于流域山顶高度和冰川规模的差异。因此，岩性并不是北天山冰川侵蚀的主控因素。

已有研究表明，较高的侵蚀速率发生在构造抬升较快的区域［65］。全球不同区域，冰川侵蚀汇总的结果也显示，无论在何种时间尺度下，构造活跃区域的冰川侵蚀速率大于构造不活跃区域，指示了构造对冰川侵蚀存在重要的影响［13］。北天山地区构造抬升速率自西向东逐渐减小，而7个冰川流域的Hkr值也有相应的变化趋势，这种一致性可能反映了构造活动在北天山冰川侵蚀的分布中也起着重要的作用。

综上所述，冰川侵蚀的7个定量影响因素中，顶点高程发挥着主控作用，2个定性因素中岩性不是主控因素，而构造可能发挥了重要的作用。那么，构造活动是通过怎样的过程影响着冰川的侵蚀呢？前人的研究认为北天山山体海拔自西向东的降低体现了构造活动的逐渐减弱［66］；河流地貌的研究也发现，北天山地形主要反映了构造活动的影响［67］。地形起伏、坡度，以及山顶高度等地貌计量指标通常用来作为构造活动的代用指标［67-70］，皮尔逊相关分析的表格中也表明山顶高度、地形起伏度、平均坡度具有一致的变化趋势，指示了相同的构造活动变化规律（表2）。因此，笔者认为北天山顶点高程与构造活动间存在的内在联系，可能显示了构造在冰川侵蚀中发挥作用的过程。诸多的研究结果显示，构造活动影响冰川侵蚀的途径有三种：（1）山体抬升后，顶点高程和冰川物质平衡线间的高差增加、积累区变大［71］、冰川增厚［58］，最终造成冰川侵蚀增强；（2）构造抬升后，高差增大、地形变陡、可能导致冰川滑动速度增加，冰川侵蚀速率也随之变大［63］；（3）构造活动还有可能使基岩破碎，更易被冰川侵蚀［62］。北天山流域面积与顶点高程（0.626）和Hkr（0.684）之间较好的相关性、Hkr与坡度（0.123）之间的弱相关性（表2），皆表明冰川侵蚀响应构造活动的过程可能是第（1）种方式，而不是第（2）种方式，至于第（3）种方式中基岩破碎程度这一因子是衡量基岩抗侵蚀能力强弱的一部分，从前面的岩性分析中可以知道基岩抗侵蚀能力并不是影响侵蚀的主要因素［72］。综上，在北天山构造活动主要通过第（1）种过程影响冰川的侵蚀，即构造活动越强烈的区域，顶点高程越大，冰川作用正差也可能越大，进而导致积累区面积、冰量也越大，侵蚀更强。根据冰川滑动模型、侵蚀模型以及数值模拟的结果可知冰川侵蚀速率是底部滑动速度的方程，而底部滑动速度与冰川厚度呈正相关关系［16，73-75］，其中冰川厚度又是冰量（或冰川规模）的体现。因此，构造在冰川侵蚀中发挥作用的实质还是通过冰川规模实现的，至于其到底发挥多大程度的作用则需要进一步研究。

祁连山的研究中也发现顶点高程和冰川侵蚀之间存在一定的相关性，但是与流域面积间仅存在非常弱的相关关系，反而与物质平衡线之间存在较好的相关性（北天山各流域ELA相差不显著，顶点高程和ELA之间无明显关系），反映的可能是冰川侵蚀对山顶高程的限制作用［17］。可见同一因子与冰川侵蚀间的关系也是存在较大差异的。此外，从降水量与流域面积间较好的相关性（0.799）可知，降水量虽不是主控因素，其也通过影响冰川规模对北天山的冰川侵蚀产生着重要的影响，这一结论与全球统计的降水与侵蚀速率之间较好的正相关性结论一致［16，19］。综合北天山构造活动、地形与气候均通过冰川规模控制冰川侵蚀的结论，笔者最终得出北天山冰川侵蚀的主控因素是冰川规模。

4 结论与展望

本文基于Hkr方法，通过对北天山北坡7个冰川流域的冰川侵蚀程度进行定量化分析，获得了北天山冰川侵蚀程度的定量数据，并在此基础上运用定量和定性分析相结合的方法探讨了北天山冰川侵蚀程度空间分布差异的主要控制因素，主要结论如下：

（1）北天山玛纳斯河段至乌鲁木齐河段冰川侵蚀程度的空间分布，大致表现为西段冰川流域强于东段，即冰川侵蚀程度在玛纳斯河源最强，其后依次为塔西河源、古仁河源、三屯河源、呼图壁河源、头屯河源、乌鲁木齐河源。

（2）北天山冰川侵蚀程度空间分布上的差异，是构造、地形、气候等多种因素共同作用所致，其中顶点高程、降水均发挥着重要的作用，构造活动与冰川侵蚀自西向东减小的趋势一致，可能体现了构造通过影响顶点高程、积累区面积、冰川规模，进而在冰川侵蚀中发挥作用。

（3）北天山冰川侵蚀程度的空间分布特点与顶点高程的相关性最好，顶点高程通过影响冰川规模、底部滑动速度来控制冰川侵蚀；加之该区域降水也是通过影响冰川规模来对冰川侵蚀施加作用。因此，冰川规模才是控制本研究区冰川侵蚀的根本因素。

本文在北天山冰川侵蚀的主控因素、构造活动影响冰川侵蚀的可能过程方面获得了一些初步的认识，但是尚存一些不足之处。第一，因研究区构造抬升速率的定量数据较为缺乏，本文仅定性分析了构造的影响，对于构造在冰川侵蚀中发挥多大程度的作用仍不明确，今后可在获得定量构造抬升速率的基础上去进一步探讨。第二，迄今为止，学者们仅能用一些代用指标或根据古冰川地貌制图来获得对冰川规模的估计，但是无论哪种方法其实都很难准确获得古冰川规模的数据，这也是为什么文中指示冰川规模的流域面积与冰川侵蚀间的相关性并不是最强的原因。第三，因本文是以造活动自西向东的变化作为条件来选取的冰川流域（冰川流域不能太接近，以免难以反映构造的变化），且衡量侵蚀的方法决定了冰川流域现代冰川的规模不能太大，因此选取的冰川流域较为有限，条件允许的话可以选取更多的冰川流域来进行分析。