拜泉县的泥沙连通性研究

Research on Sediment Connectivity in Baiquan County

摘要: 东北黑土地区的黑土地地况决定粮食的生产保障，对中国的粮食安全问题至关重要。在东北黑土区，关于泥沙连通性的研究仍处在起步阶段。以拜泉县为研究区域，结合拜泉县的高程DEM、归一化植被指数NDVI等数据，通过计算泥沙连通性指数(IC)的大小，得出1990年、1995年、2000年、2005年、2010年、2015年和2020年的泥沙连通性指数的时空分布特征，将其运用Sen趋势分析、M-K检验分析得到了1990~2020年拜泉县的IC趋势变化情况，并且将拜泉县不同高程、坡度分别与IC进行叠加分析。结果表明：拜泉县1990~2020年泥沙连通性指数均值整体上呈现减小趋势，其中，1990年、1995年、2000年、2005年、2010年、2015年、2020年IC泥沙连通性指数IC均值分别为−9.27、−10.45、−12.34、−12.50、−12.57、−12.53、−12.40，泥沙连通性强弱与泥沙连通性指数的大小呈正相关，进而拜泉县整体泥沙连通性呈减弱趋势；拜泉县1990~2020年IC趋势特征为：其西北、西南高程较低、坡度较低区域的IC趋势为减少，其东南和中间河流区域的IC趋势为增加，其余区域的IC趋势为不显著变化；拜泉县高程在210~270 m之间、坡度在1˚~4˚之间、植被覆盖少的区域的泥沙连通性指数值较大，泥沙连通性较强，这些区域往往更容易发生黑土侵蚀。

Abstract: The fertile black soil in Northeast China plays a pivotal role in ensuring the nation’s food security, making it a cornerstone of China’s agricultural stability. In the black soil region of Northeast China, research on sediment connectivity is still in its nascent stages. Taking Baiquan County as the research area, by integrating data such as the Digital Elevation Model (DEM) and the Normalized Difference Vegetation Index (NDVI) of Baiquan County and calculating the magnitude of the sediment connectivity index (IC), the spatio-temporal distribution characteristics of the sediment connectivity index from 1990, 1995, 2000, 2005, 2010, 2015 to 2020 were obtained. Subsequently, by applying Sen’s slope estimator and the Mann-Kendall test, the trend changes of the sediment connectivity index (IC) in Baiquan County from 1990 to 2020 were analyzed. Additionally, an overlay analysis was conducted to examine the relationship between the IC and different elevations and slopes in Baiquan County. The results show that the mean value of the sediment connectivity index in Baizhan County from 1990 to 2020 generally shows a decreasing trend. Among them, the mean values of the IC sediment connectivity index IC in 1990, 1995, 2000, 2005, 2010, 2015, and 2020 were −9.27, −10.45, −12.34, −12.50, −12.57, −12.53, and −12.40 respectively. The strength of sediment connectivity is positively correlated with the magnitude of the sediment connectivity index. Therefore, the overall sediment connectivity in Baizhan County shows a weakening trend. The IC trend characteristics of Baiquan County from 1990 to 2020 are as follows: In the areas in the northwest and southwest with lower elevations and lower slopes, the IC trend is decreasing; in the southeastern and central river areas, the IC trend is increasing; and in the remaining areas, the IC trend shows no significant change. In the areas of Baiquan County with elevations ranging from 210 to 270 meters, slopes ranging from 1˚ to 4˚, and less vegetation coverage, the sediment connectivity index values are larger and the sediment connectivity is stronger. These areas are often more prone to black soil erosion.

1. 引言

黑土是东北地区耕地土壤较好的一个土壤类型，它是腐殖质积累与淋溶共同作用的产物。黑土区是中国主要粮食产区和重要的商品粮基地，由于社会系统与生态系统的互动性加强，人类社会经济发展需求的过度干扰，导致黑土区土壤质量和生态功能严重退化，并影响东北黑土区农业的可持续发展[1]。本研究所选择的拜泉县，其具有“全国产粮百强县”的美称，对粮食产量贡献较大，由此可见，该县耕地粮食产量的稳定和增加对区域粮食贡献具有重大作用，然而近年来耕地侵蚀沟不断增加，水土流失严重[2]。泥沙连通性概念有助于理解泥沙在地貌景观单元间的泥沙输移过程，量化泥沙在地貌单元之间的连接强度，可进一步掌握流域内泥沙输移的内部动态，调控径流泥沙输移路径[3]。研究泥沙连通性的方法包括野外调查法、图论法、指标法和模型模拟法，不同方法的理论基础、数据要求、实施过程与结果存在一定差异，连通性指数目前应用最为广泛，但该指数更多强调结构连通性[4]。连通性大的区域，土壤侵蚀越容易发生，反之，则土壤侵蚀越不容易发生。东北黑土侵蚀问题越来越受到大家的重视，黑土层也在逐年变薄，而形成1 cm黑土层需要几百年时间。泥沙连通性是表征泥沙在地貌单元间物理级联关系的重要参数，是探究泥沙来源、识别水土保持重点区域的基础，探讨泥沙连通性具有重要意义[5]。而且，泥沙连通性可以反映泥沙源汇的潜在联系，识别流域水土流失热点区域及泥沙迁移路径。研究泥沙连通性的影响因素有助于更好地理解泥沙连通性的时空变化特征[6]。连通性指数IC由于计算过程简单，所需数据(DEM、土地覆盖数据)较少，机理性强，被广泛应用于各个地区的泥沙连通性的研究中，但是目前与东北黑土地区有关的研究较少。土壤侵蚀与水文泥沙连通性呈显著的正相关，主要体现在坡度、坡向、土地利用和高程上[7]。本研究以拜泉县为研究区域采用泥沙连通性指数IC的大小来表征这个区域的泥沙连通性的强弱，结合近30年的IC值的时空变化特征，来分析近30年的泥沙连通情况，进而分析出该县的哪些区域往往更容易发生黑土侵蚀。

2. 材料与方法

2.1. 研究区概况

拜泉县，位于黑龙江省中西部、齐齐哈尔市东部。乌裕尔河与通肯河之间，松嫩平原北边缘。地理位置介于东经125˚30′~126˚31′，北纬47˚20′~47˚55′之间，东西宽66 km，总面积3599.15 km²。高程为173~339 m，地势由东北向西南倾斜，东北高，西南低。拜泉县属寒温带大陆性季风气候，春季降水少而且风大，夏季短而多雨，秋季早霜，冬季漫长寒冷干燥，四季温差较大，年平均气温1.2℃。

2.2. 数据来源

DEM数据采用中国科学院计算机网络数据中心地理空间数据云平台提供的分辨率为30 m的ASTER GDEM高程数据，归一化植被指数NDVI数据来源于中国科学院资源环境科学与数据中心，分辨率为30 m。

2.3. 研究方法

2.3.1. 连通性指数的计算

利用Borselli等[8]提出的连通性指数(IC，式中记为IC)，量化拜泉县的泥沙连通程度。该指数综合考虑地表植被覆盖和地形状况，将流域的泥沙连通过程分为集水区泥沙向下输移和泥沙被输移至最近泥沙汇入两个过程，以此表征泥沙在源–汇间的连通可能性，其表达式为

$\text{IC}={\log }_{10}\frac{{\text{D}}_{\text{up}}}{{\text{D}}_{\text{dn}}}={\log }_{10}\left(\frac{\overline{\text{W}}\cdot \overline{\text{S}}\cdot \sqrt{\text{A}}}{\sum \frac{{\text{d}}_{\text{i}}}{{\text{w}}_{\text{i}}\cdot {\text{s}}_{\text{i}}}}\right)$ ⑴

式中：D_up为集水区泥沙向下输移的潜在可能性；D_dn为泥沙经过水流路径输移至最近泥沙汇的潜在可能性；$\overline{\text{W}}$、$\overline{\text{S}}$、A分别代表集水区的地表状况和土地利用的平均加权系数、平均坡度、面积；w_i、s_i、d_i分别为第i个栅格单元的加权系数、第i个单元的坡度、第i个单元沿最陡坡度到达交叉汇点的运移距离。IC的范围为(−∞, +∞)，值越大代表该流域的泥沙连通可能性越大。

IC的计算以流入主河道为参照，其模型的搭建与计算分析在ArcGIS Pro3.0.2软件中完成，共得到1990年、1995年、2000年、2005年、2010年、2015年、2020年7个年份IC值。其中，集水区面积A、坡度S和径流路径长度d根据DEM数据计算得到；权重因子W应用修正通用土壤流失方程(RUSLE)中的植被覆盖及管理因子(C因子)表示，本研究为了得到连续年份的C值，采用蔡崇法等[9]建立的C因子与植被覆盖度$\text{fc}$ 之间的关系式计算拜泉县30年来的植被覆盖管理措施因子，方法如下：

$fc=\frac{\text{NDVI}-{\text{NDVI}}_{\text{min}}}{{\text{NDVI}}_{\text{max}}-{\text{NDVI}}_{\text{min}}}$ (2)

$\text{C}=\{\begin{array}{c}1\\ 0.6508-0.3436\lg fc\\ 0\end{array}\begin{array}{c}0\le fc<0.1\%\\ 0.1\text{\%}\le fc<78.3\%\\ 78.3\text{\%}\le \text{fc}\end{array}$ (3)

2.3.2. 统计分析

利用ArcGISPro3.0.2获取拜泉县1990年、1995年、2000年、2005年、2010年、2015年、2020年的IC栅格数据的分布及其均值，通过将得到的栅格数据在Matlab R2024b中运用Theil-Sen Median趋势分析、Mann-Kendall检验分析得到了1990~2020年拜泉县的IC趋势变化情况。

Theil-Sen Median方法又称为Sen斜率估计，是一种稳健的非参数统计的趋势计算方法[10]。该方法计算效率高，对于测量误差和利群数据不敏感，适用于长时间序列数据的趋势分析[11]。其计算公式为：

$\beta =\text{Median}\left(\frac{{\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}}{\text{j}-\text{i}}\right)\forall \text{j}>\text{i}$ (4)

式中：Median ( )代表取中值，若$\beta$ 大于零，则表明IC为增长趋势，反之为下降趋势。Mann-Kendall(MK)检验是一种非参数的时间序列趋势性检验方法，其不需要测量值服从正态分布，不受缺失值和异常值的影响，适用于长时间序列数据的趋势显著检验[12]。其过程如下：对于序列X_t = x₁，x₂，……x_n，先确定所有对偶值(x_i, x_j, j > i)中x_i与x_j的大小关系(设为S)。作如下假设：|H| < 0，序列中的数据随机排列，即无显著趋势；|H| > 1，序列存在上升或下降趋势。检验统计量S计算公式为：(在Z的计算公式中，当S > 0时，分子为S − 1) [13]。

$\text{s}={\displaystyle \sum }_{\text{i}=1}^{\text{n}-1}{\displaystyle \sum }_{\text{j}=\text{i}+1}^{\text{n}}\text{sgn}\left({\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}\right)$ (5)

式中：n是时间序列中数据点总数，x_i和x_j分别是时间点i和j的观测值(j > i)，sgn ()为符号函数，计算公式为：

$\text{sgn}\left({\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}\right)=\{\begin{array}{cc}+1& {\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}>0\\ 0& {\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}=0\\ -1& {\text{x}}_{\text{j}}-{\text{x}}_{\text{i}}<0\end{array}$ (6)

使用检验统计量Z进行趋势检验，Z值的计算方法如下：

$\text{Z}=\{\begin{array}{cc}\frac{\text{S}-1}{\sqrt{\text{Var}\left(\text{S}\right)}}& \left(\text{S}>0\right)\\ 0& \left(\text{S}=0\right)\\ \frac{\text{S}+1}{\sqrt{\text{Var}\left(\text{S}\right)}}& \left(\text{S}<0\right)\end{array}$ (7)

式中：Var的计算公式为：

$\text{Var}\left(\text{S}\right)=\frac{\text{n}\left(\text{n}-1\right)\left(2\text{n}+5\right)-{{\displaystyle \sum }}_{\text{p}=1}^{\text{m}}{\text{t}}_{\text{p}}\left({\text{t}}_{\text{p}}-1\right)\left(2{\text{t}}_{\text{p}}+5\right)}{18}$ (8)

式中：n是数据点总数，m是数据中“结”(具有相同数值的组)的数量，t_p是第p个“结”的大小(即该数值出现的次数)。

同样采用双边趋势检验，在给定显著性水平下，在正态分布表中查得临界值为Z1 − α/2。当|Z| ≤ Z1 − α/2时，接受原假设，即趋势不显著；若|Z| > Z1 − α/2，则拒绝原假设，即认为趋势显著。本文给定显著性水平α = 0.05，则临界值Z1 − α/2 = ± 1.96，当Z的绝对值大于1.65、1.96和2.58时，表示趋势分别通过了信度为90%、95%和99%的显著性检验[14]。趋势显著性的判断方法见表1。

Table 1. Mann-Kendall test trend categories

表1. Mann-Kendall检验趋势类别

β	Z	趋势类别	趋势特征
β > 0	2.58 < Z	4	极显著增加
	1.96 < Z ≤ 2.58	3	显著增加
	1.65 < Z ≤ 1.96	2	微显著增加
	Z ≤ 1.65	1	不显著增加
β = 0	Z	0	无变化
β < 0	Z ≤ 1.65	−1	不显著减少
	1.65 < Z ≤ 1.96	−2	微显著减少
	1.96 < Z ≤ 2.58	−3	显著减少
	2.58 < Z	−4	极显著减少

3. 结果与分析

3.1. 泥沙连通性指数分布与趋势特征

注：该图基于天地图官网下载的审图号为GS (2024) 0650号的标准地图制作，底图无修改。

Figure 1. Spatial distribution maps of sediment connectivity in Baiquan County from 1990~2020

图1. 1990~2020年拜泉县泥沙连通性空间分布图

计算拜泉县内每个像元的IC年均值，得到了1990年、1995年、2000年、2005年、2010年、2015年、2020年泥沙连通性指数(IC)均值分别为−9.27、−10.45、−12.34、−12.50、−12.57、−12.53、−12.40，拜泉县1990~2020年总体IC均值呈现减小的趋势。从1990~2020年的泥沙连通性指数(IC)分布图(见图1)中可以看出IC大小分布特征为拜泉县四周IC值偏大，中部除河流外IC值偏小，大部分区域随着高程的降低，IC值变大，主要因为随着高程的降低，流域河网密度增大，低高程带地区泥沙更容易输移至河道，这与路志强[15]的研究基本一致，随着坡度的增加表现为IC值先增后减的变化走向，拜泉县IC值较大的地方大致分布在河流湖泊、水库附近，而泥沙连通性的强弱往往与IC值的大小呈正相关，因此可以通过拜泉县IC值大小的分布情况来判断拜泉县泥沙连通性强弱的分布情况，拜泉县泥沙连通性较强的分布为其区域四周以及县内河流、水库。从1990年到2010年拜泉县泥沙连通性指数均值减小，而从2010年到2020年拜泉县泥沙连通性指数均值略微增大，但是1990年到2020年拜泉县泥沙连通性指数整体开始减小，全县整体上来看泥沙连通性均有不同程度的减弱。

将拜泉县1990年、1995年、2000年、2005年、2010年、2020年的IC栅格图在Matlab R2024b中运用Theil-Sen Median趋势分析、Mann-Kendall检验分析得到了1990~2020年拜泉县的IC趋势变化情况(见图2)，可以看到1990~2020年拜泉县其西北、西南高程较低、坡度较低区域的IC趋势为减少，其东南和中间河流区域的IC趋势为增加，其余区域的IC趋势为不显著变化，这与计算出的1990~2020年拜泉县IC大小整体呈减小的趋势是一致的。

注：该图基于天地图官网下载的审图号为GS (2024) 0650号的标准地图制作，底图无修改。

Figure 2. Trend characteristic map of IC values from 1990~2020

图2. 1990~2020年IC值趋势特征图

3.2. 不同高程下泥沙连通性特征

不同高程下拜泉县IC值的趋势特征以及泥沙连通性强弱的关系，将拜泉县的坡度按<190 m、190~210 m、210~230 m、230~250 m、250~270 m、270~290 m、290~310 m和>310 m共8个高程区间对应划分为1~8等级，将拜泉县高程栅格图在ArcMap10.8上按这8个等级进行重分类(见图3)。

将拜泉县高程重分类后的栅格图和拜泉县IC值趋势分析栅格图叠加后，得到了不同高程区间的IC值趋势特征，拜泉县1990~2020年的IC值的趋势特征为极显著减少和显著减少区域主要分布在210~230 m、230~250 m和250~270 m高程区间，随着高程的增加，拜泉县1990~2020年的其IC值的趋势特征的变化为先变大后变小，泥沙连通性先变强后变弱。这说明高程更高的地方泥沙从源到汇其输移路径更长，泥沙输移也会越困难。总的来说，随着时间推移，拜泉县在不同高程带泥沙连通性得到了降低，其泥沙输移难度增加。拜泉县泥沙连通性随高程变化分异特征明显，这主要因为拜泉县地势呈东北向西南倾斜，随着高程的降低，流域河网密度增大，低高程地区泥沙到达河道的路径短，泥沙更容易输移至河道，泥沙连通性变大。

3.3. 不同坡度下泥沙连通性特征

不同坡度下拜泉县IC值的趋势特征以及泥沙连通性强弱的关系，将拜泉县的坡度按0˚~1˚、1˚~2˚、2˚~3˚、3˚~4˚、4˚~5˚、5˚~6˚、6˚~7˚和>7˚共8个坡度区间对应划分为1~8等级，将拜泉县高程栅格图在ArcMap10.8上提取坡度后按这8个等级重分类(见图3)。

将拜泉县坡度重分类后的栅格图和拜泉县IC值趋势分析栅格图叠加后，得到了不同坡度区间的IC值趋势特征，拜泉县1990~2020年的IC值的趋势特征为极显著减少和显著减少区域主要分布在1˚~2˚、2˚~3˚和3˚~4˚坡度区间，随着坡度的升高，其IC值的趋势特征的变化为先变大后变小，随着坡度增大，泥沙连通性变强，但大于6˚的坡度之后其IC值的趋势特征变化不明显，表明其泥沙连通性开始减弱，该区间对于泥沙连通性强弱起主导作用的因素可能是坡长、汇水面积、坡面形态等其他地形地貌因素，而不是坡度因素。

注：该图基于天地图官网下载的审图号为GS (2024) 0650号的标准地图制作，底图无修改。

Figure 3. Baiquan County slope classification and elevation classification

图3. 拜泉县坡度重分类和高程重分类

4. 结论与讨论

1990~2020年拜泉县泥沙连通性指数(IC值)整体呈减小趋势，其IC均值也是呈现出整体减小的趋势，泥沙连通性呈减弱趋势。不同高程下拜泉县IC值的趋势特征以及泥沙连通性强弱的关系为：拜泉县1990~2020年的IC值的趋势特征为极显著减少和显著减少区域主要分布在210~230 m、230~250 m和250~270 m高程区间，随着高程的增加，拜泉县1990~2020年的其IC值的趋势特征的变化为先变大后变小，泥沙连通性先变强后变弱；不同坡度下拜泉县IC值的趋势特征以及泥沙连通性强弱的关系为：拜泉县1990~2020年的IC值的趋势特征为极显著减少和显著减少区域主要分布在1˚~2˚、2˚~3˚和3˚~4˚坡度区间，随着坡度的升高，其IC值的趋势特征的变化为先变大后变小，随着坡度增大，泥沙连通性变强，但大于6˚的坡度之后其IC值的趋势特征变化不明显，表明其泥沙连通性开始减弱。但是，本研究使用的研究方法存在一定的局限性，比如泥沙连通性指数(IC)的计算的权重因子为DEM和植被覆盖及管理因子(C因子)，未考虑降雨、人为活动等动态的因素，计算出来的IC值是结构上的，可能会导致在后续的不同高程、坡度分别与IC叠加分析的准确性。后期研究会对泥沙连通性指数的计算考虑引入更多的权重因子，以及进行更细致、全面的分析。

本文通过在ArcGIS Pro3.0.2上搭建泥沙连通性指数计算的模型来计算泥沙连通性指数的大小，进而分析了1990年、1995年、2000年、2005年、2010年、2015年和2020年的泥沙连通性指数的时空分布特征，将这7年的泥沙连通性指数的栅格数据在Matlab R2024b上运用Sen趋势分析、M-K检验分析得到了1990~2020年拜泉县的IC趋势变化情况，并将拜泉县不同高程、坡度分别与IC叠加分析，得出拜泉县1990~2020年泥沙连通性指数大小呈减小趋势、泥沙连通性呈减弱趋势，拜泉县高程在210~270 m之间、坡度在1˚~4˚之间、植被覆盖少的区域的泥沙连通性指数值较大，泥沙连通性较强，这些区域往往更容易发生黑土侵蚀，这对拜泉县的水土流失治理可以提供参考，也为东北黑土侵蚀的分析与治理起到帮助。