2021十大孔兑区域土壤侵蚀空间分布特征

于海云，覃湘栋，庞治国，3，王志军，付俊娥，3

(1.内蒙古机电职业技术学院，内蒙古 呼和浩特 010070；2.中国水利水电科学研究院，北京 100038；3.水利部防洪抗旱减灾工程技术研究中心，北京 100038；4.内蒙古自治区水利事业发展中心，内蒙古 呼和浩特 010020)

1 研究背景

十大孔兑是指黄河内蒙古段由南向北并列流入黄河的十条一级支流所构成的流域，其不仅是晋陕蒙接壤地区国家能源重化工基地的组成部分，还分布有包钢生产生活水源地和鄂尔多斯市现代农牧业产业带和新兴工业园，整体区位十分重要[1]。十大孔兑位于黄河中上游多沙粗沙区，流域内暴雨频发，加之流域中上游属于极易发生土壤侵蚀的地貌类型，造成十大孔兑区域较为严重的土壤侵蚀，这不仅会破坏流域上游的土地资源和生态环境，同时由于土壤侵蚀，暴雨引发的洪水会携带大量泥沙泄入流域下游与黄河，造成下游农田与工业生产设施冲毁、黄河主河道淤堵等灾害。整体而言，土壤侵蚀对区域的社会经济发展与自然环境保护都极为不利[2-4]，已成为十大孔兑区域最为严重的环境问题之一，开展土壤侵蚀评估与监测研究十分重要。

土壤侵蚀的评估与监测早是通过在监测站点建立标准径流小区[5]和集沙盘[6-7]，测算土壤流失量以评价研究区整体的侵蚀状态，该方法监测面积小、成本高，很难表征土壤侵蚀的空间变化趋势。虽然站点监测有缺陷，但可通过监测数据建立多因子与土壤侵蚀量之间的关系模型[8-11]，这些因子(植被覆盖、坡度、风速、土地利用等)多数是可以通过遥感技术获取，这使得基于遥感手段监测评估土壤侵蚀成为可能。国内土壤侵蚀遥感监测评估研究较多，如童珊等结合RUSLE模型、CA-Markov模型及LMDI模型对祁连山南坡2000到2019的土壤侵蚀进行了时空动态监测，分析了不同地形条件下的土壤侵蚀变化特征，并对土壤侵蚀影响因素进行了定量分析研究[12]；葛德祥等基于USLE模型开展二滩库区蓄水前后的土壤侵蚀变化评估[13]；高泗强等基于修正后的RUSLE模型进行土壤侵蚀模数计算，后借助地理信息系统进行空间数据的统计分析，总结出重庆市县级行政区尺度的土壤侵蚀特征[14]；邹秦英等基于2018的土地利用数据，借助RUSLE模型计算了汾河上游的土壤侵蚀模数，并对汾河上游的土壤侵蚀空间分布特征进行了总结[15]。

十大孔兑作为内蒙古自治区土壤侵蚀的重点关注区域，该区域1951—2010向黄河输送的均泥沙量为1939万t[16]，超过该流域入黄泥沙总量的十分之一，由于严重的土壤侵蚀引发的各类环境问题，我国在十大孔兑区域开展了大量的生态监测与治理研究。如刘殿君等[17]以淤积概化法监测了十大孔兑区域354座淤地坝多个时的减沙量，其认为淤地坝仍是该区域最有效的水土保持措施；滑永春等[18]基于遥感生态指数评价模型分析了十大孔兑区域近20生态环境健康度的分布格局及变化趋势，认为该区域的生态环境仍有待改善。所以为十大孔兑区域生态环境的良性发展，土壤侵蚀的定评估与空间特征分析十分必要，本文基于多源数据结合中国水土流失方程(Chinese Soil Loss Equation)与多地类的风力侵蚀模型对2021的十大孔兑区域的土壤侵蚀进行评估，而后基于评估结果总结十大孔兑土壤侵蚀的空间分布特征，并综合多个因素开展土壤侵蚀特征分析，为十大孔兑的土壤侵蚀治理与生态恢复提供参考。

2 研究区及数据

2.1 研究区十大孔兑位于内蒙古自治区西南部鄂尔多斯市境内，流域总土地面积为1.08万km2，整体地势南高北低，流域由南向北包含砒砂岩丘陵沟壑、库布齐沙漠、冲洪积平原三个侵蚀地貌类型。十大孔兑区域气候上属于典型大陆性气候，平均气温7℃左右，平均降雨在200～400 mm，蒸发量约为2200 mm。十大孔兑的主要土地利用类型是林地，面积为3830.09 km2，其次是草地和耕地，面积分别为3019.5 km2和2013.93 km2，具体分布如图1所示。

图1 十大孔兑区域土地利用空间分布图Fig.1 Sketch map of land use in Ten Kongduis

2.2 基础数据土壤侵蚀模数的计算需要多种数据支撑，包括地形、土地利用、植被覆盖、降雨、风速、表土湿度等，下面介绍数据来源。植被数据采用基于Landsat多光谱数据对MODIS归一化植被指数产品降尺度至10 m分辨率的植被指数数据，每半月一，1共计24数据，每进行更新。土地利用数据采用的是鄂尔多斯市2021土地利用目视解译成果，目视解译使用的遥感影像是基于高分1、2号卫星遥感数据进行数据融合后形成的2 m分辨率数字正射影像，解译成果经过多次修正与实地验证，整体精度较高，每进行更新。降雨、风速、土壤侵蚀因子、表土湿度、坡度坡长因子等数据均采用水利部2021水土流失动态监测下发的官方数据，其中：降雨、风速和表土湿度数据是由长时序监测资料或遥感产品综合计算得出，基本可以反映十大孔兑区域的一间的气象变化规律；土壤侵蚀因子数据是综合水利普查数据和区域内径流小区监测资料计算得出；坡度坡长因子是由30 m数字高程模型计算得出，这些数据水利部每五更新一次，使数据的科学可靠性有所保证。

3 方法

3.1 水蚀模数计算土壤的水力侵蚀计算模型较多，包括通用水土流失方程、改进通用水土流失方程等，其与区域降雨、土壤质地、坡度、植被等多个参数相关。本文参考2022度水土流失动态监测技术指南，采用中国土壤流失方程CSLE[19]计算十大孔兑区域的土壤水力侵蚀模数。方程基本形式如式(1)：

Mw=RKLSBET

(1)

式中：Mw为土壤水力侵蚀模数，t·hm-2·a-1；R为多平均降雨侵蚀力因子；K为土壤可蚀性因子，t·hm2·h·hm-2·MJ-1·mm-1；B为植被覆盖和生物措施因子，利用MODIS归一化植被指数产品和Landsat多光谱影像采用参数修订方法计算得出；L、S、E、T分别为是坡长因子、坡度因子、水保工程措施因子和土地耕作措施因子，均无量纲，各因子具体计算参考2022度水土流失动态监测技术指南[20]。

3.2 风蚀模数计算土壤风力侵蚀主要发生在沙地、裸土地、耕地、草地上，结合十大孔兑的土地利用现状，参考2022度水土流失动态监测技术指南，本文选用了三种风力侵蚀模型[20-21]，分别针对耕地、园林草地和沙地。

耕地风力侵蚀模型的基本形式如式(2)所示：

(2)

式中：Mfa为每半月耕地上的风力侵蚀模数，t·hm-2·a-1；W为半月内的表土湿度因子，其值在0到1之间；Hj为半月内各风速等级累计发生的时间，min，具体计算参考2022度水土流失动态监测技术指南；Z0为土地表面粗糙程度，cm；j为风速等级的序号，风速在5～40 m/s内按照1 m/s为间隔划分成35个等级，取值为1到35的整数；Uj为第j个等级下的平均风速，m/s，例如风速等级是5～6 m/s，Uj的值就是5.5 m/s。

园林草地上的风力侵蚀模型计算公式如式(3)所示：

(3)

式中：Mfg是每半月内的园林草地上的风力侵蚀模数，t·hm-2·a-1；V是植被覆盖度，根据归一化植被指数根据经验公式计算得出。

沙地上的风力侵蚀模数采用式(4)进行计算。

(4)

式中Mfs是每半月内沙地上的风力侵蚀模数，t·hm-2·a-1。

上述风力侵蚀计算都是半月时间尺度，为得到区域上全的风力侵蚀模数，需要进行求和运算，计算公式如式(5)所示：

(5)

式中：Mtotal是全风力侵蚀模数；Mhaf是各种地类上半月内的侵蚀模数。

3.3 侵蚀分级侵蚀模数是根据气象、植被等数据计算的土壤侵蚀估计量，根据侵蚀的估计量可对区域的土壤侵蚀强度进行分级，本文主要参考了土壤侵蚀分类分级标准(SL190-2007)[22]，将风力与水力侵蚀均划分为微度、轻度、中度、强烈、极强烈和剧烈六个等级，具体分级标准见表1。

表1 土壤侵蚀强度分级标准(单位：t·hm-2·a-1)

3.4 侵蚀综合分析模数分级后，同一区域即会有一个风力侵蚀强度，也会有水力侵蚀强度，需要依据侵蚀等级判断该区域发生土壤侵蚀的主要方式。参考十大孔兑往的土壤侵蚀监测成果，十大孔兑区域的水力侵蚀面积大于风力侵蚀面积，所以采取以下分析策略：若该像元风力侵蚀强度大于水力侵蚀强度，则将该像元设定为风力侵蚀区并采用其风力侵蚀强度作为最终土壤侵蚀强度；若该像元的水力侵蚀强度大于等于风力侵蚀强度时，则将其认定为水力侵蚀区并将其水力侵蚀强度作为最终土壤侵蚀强度。

4 十大孔兑2021度土壤侵蚀评估结果

经上述方法计算可以得到十大孔兑区域2021土壤侵蚀的空间分布(如图2(a))。轻度及以上等级被判定为该区域属于土壤侵蚀区，则十大孔兑土壤侵蚀区面积为4398.85 km2，占十大孔兑总土地面积的40.86，其中水力侵蚀面积为2461.56 km2，占总侵蚀面积的55.96，风力侵蚀面积为1937.29 km2，占总侵蚀面积的43.04。本研究也收集了2020的基础数据并按照上述方法评估了该十大孔兑区域的土壤侵蚀(图2b)。相较于2020，2021十大孔兑区域的土壤侵蚀面积呈现减少的趋势，侵蚀总面积减少了55.16 km2，其中水力、风力侵蚀面积分别减少了42.95 km2和12.21 km2。综合两土壤侵蚀的空间分布，水力侵蚀主要分布在十大孔兑上游，且地块较为破碎，整体分布离散；风力侵蚀主要分布于十大孔兑中下游，风力侵蚀地块整体相对完整，分布较为集中。

图2 2021和2020十大孔兑区域土壤侵蚀空间分布图Fig.2 Spatial distribution diagram of soil erosion in Ten Kongduis in 2021 and 2020

2021十大孔兑区域各强度等级的土壤侵蚀面积如图3，轻度等级、中度等级、强烈等级、极强烈等级和剧烈等级土壤侵蚀面积为2614.38 km2、1328.74 km2、276.29 km2、194.45 km2和40.15 km2，分别占总侵蚀面积的59.62、29.68、5.92、4.01和0.77。十大孔兑区域土壤侵蚀面积中：轻度和中度土壤侵蚀强度面积占比较大，两强度面积占总侵蚀面积的89.30，其中风力和水力侵蚀面积基本各占一半，整体别较小；强烈及以上侵蚀面积在总侵蚀面积中占比较小，仅为10.70，且该部分土壤侵蚀类型主要是水力侵蚀，风力侵蚀占比极小。

图3 2021十大孔兑区域土壤侵蚀强度统计柱状图Fig.3 Statistical histogram of soil erosion intensity of Ten Kongduis in 2021

5 分析与讨论

2021鄂尔多斯市水土保持公报[23]显示十大孔兑风力侵蚀面积为1113.77 km2、水力侵蚀面积为2401.02 km2。对比本研究的评估结果，发现水力侵蚀面积较为接近，但风力侵蚀面积存在较大的异。进一步对比两者的空间分布，发现公报数据与本研究评估结果在空间分布上基本一致，但公报的统计范围与本研究有异，公报仅统计十大孔兑上游(区域南部)和部分下游区域(区域北部十条河流、沟道的沿岸)的土壤侵蚀，结合十大孔兑上游水蚀、中下游风蚀的土壤侵蚀空间分布特征，可知本研究相较于公报数据多出的风蚀面积主要源于区域下游未参与公报统计范围的风力侵蚀区。以此而言，本研究土壤侵蚀评估结果具备一定的准确性与合理性的，所以可以围绕本研究评估结果开展下述的讨论分析。

5.1 基于土地利用的土壤侵蚀空间分析将2021与2020两的土壤侵蚀结果按照两土地利用类型进行分区面积统计，得到不同土地利用类型各强度等级土壤侵蚀面积表(表2)，由表2可知，十大孔兑区域土壤侵蚀面积在不同土地利用类型上别较大。近些，草地、林地两种土地利用类型由于面积基数较大，其土壤侵蚀面积分别位居第一和第二位，同时二者在侵蚀强度分布上也极为相似，相较于其他地类侵蚀强度分布较为分散，至于侵蚀面积占该地类面积的比例，可以看出林地要明显低于草地，说明草地相较于林地更容易发生土壤侵蚀，这也是符合侵蚀规律的；耕地上的土壤侵蚀面积位居第三，在侵蚀强度上分布较为集中，轻度侵蚀占该地类总侵蚀面积的96；其他土地(包括盐碱地、沙地和裸土地)上的土壤侵蚀面积位居第四，侵蚀强度上以中度侵蚀为主，且侵蚀面积占该地类面积比例超过90，说明这些地类上极易发生土壤侵蚀；草地、林地、耕地和其他土地是十大孔兑区域土壤侵蚀的主要来源地，其余地类上的侵蚀面积占比较小。

表2 不同土地利用类型(一级类)[24]土壤侵蚀面积及比例

十大孔兑区域不同土地利用类型上的风力、水力侵蚀面积如图4所示，结合图2和图4进行分析。

图4 十大孔兑区域不同土地利用类型侵蚀面积堆积柱状图Fig.4 Schematic diagram of wind and water erosion area of Ten Kongduis in different land use

(1)十大孔兑区域发生水力侵蚀的主要土地利用类型是草地与林地，而这两种土地利用类型集中分布在十大孔兑上游，这是水力侵蚀主要分布在孔兑上游的部分原因。为了更详尽的探究其原因，需要分别对水力侵蚀模型的7个因子进行统计分析，本研究统计了不同侵蚀强度下的各因子平均值(表3)，从中可以发现，降雨侵蚀因子、坡度因子、坡长因子与土壤可蚀性因子的平均值与水力侵蚀强度呈显很强的一致性，说明降雨、地形和地貌对水力侵蚀的空间分布都有贡献，但对比各因子的变化幅度，坡度因子和坡长因子的变化幅度最大，所以推测地形坡度可能是影响水力侵蚀空间分布的重要因素，结合十大孔兑上游特殊的地貌进行分析，区域上游属于砒砂岩丘陵沟壑，大坡度地形出现频繁，且土壤抗冲蚀力不足，最终导致了南部区域严重的水力侵蚀。

表3 不同水力侵蚀强度下的因子平均值

(2)十大孔兑区域的风力侵蚀主要集中在耕地和其他土地(包括盐碱地、沙地和裸土地)上，风力侵蚀可以分为两部分：其一是由西向中部延伸的条形风蚀区，其位于库布齐沙漠，主要的土地利用类型是沙地，该区域地势平坦、风速较大且植被覆盖较，最终导致区域风力侵蚀模数较大，在1000～2500 t·hm-2·a-1之间，侵蚀强度为中度侵蚀；其二是区域下游的轻度风蚀区，其位于冲洪积平原，易发生风力侵蚀，但由于该区域土地利用类型主要是耕地，相较于沙地，耕地拥有作物植被覆盖且地表相对粗糙，具备一定的风蚀抵抗能力，所以在侵蚀模数上相对沙地较小，侵蚀强度整体集中在轻度侵蚀。

5.2 基于植被覆盖的土壤侵蚀空间分析本研究获取了202124(每半月一)的十大孔兑区域的植被覆盖度数据，考虑到季节影响选取植被最旺盛的第14数据分析植被覆盖情况。由图5可知，十大孔兑区域除去中游库布齐沙漠，植被覆盖度整体呈现出北高南低的趋势。十大孔兑下游的植被覆盖度普遍在0.7以上，高植被覆盖说明该区域抵抗水力侵蚀能力较强，模数计算上体现在B因子(植被覆盖和生物措施因子)较小，加之十大孔兑北部多为冲击平原，地势平坦小，坡度坡长因子较小，两者共同作用影响十大孔兑下游的水蚀模数，使得孔兑下游水力侵蚀面积较小；十大孔兑中游库布齐沙漠的植被覆盖度极低，基本在0.1以下，由于植被状态较，导致该区域抵抗风力侵蚀的能力较，造成风力侵蚀模数较大最终形成大面积的中度风蚀区；十大孔兑上游的植被覆盖度普遍在0.35以下，相较于一般的林地和草地，其植被覆盖度偏低，再者考虑上游是砒砂岩丘陵沟壑，部分区域由于坡度较大，致使植被无法有效的固持土壤，抵抗降雨造成的水力侵蚀，最终导致了十大孔兑上游严重的水力侵蚀，同时也阐明了上游区域水力侵蚀地块较为破碎的原因。

图5 2021十大孔兑第14植被覆盖度空间分布图Fig.5 The spatial distribution of vegetation coverage in the 14th phase of Ten Kongduis in 2021

5.3 县区尺度的土壤侵蚀空间分析十大孔兑区域覆盖鄂尔多斯市东胜区、达拉特旗、准格尔旗、杭锦旗4旗县(区)的部分区域，各旗县(区)在十大孔兑区域的占地面积分别为939.54 km2、7919.75 km2、400.91 km2、1504.51 km2，该部分的土壤侵蚀情况如图6、图7所示。分析可以发现：(1)十大孔兑区域的土壤侵蚀主要发生在达拉特旗、杭锦旗和鄂尔多斯东胜区境内，三者侵蚀面积之和占十大孔兑总侵蚀面及的97.47；(2)不同旗县境内部分的风力和水力侵蚀面积分配有明显异，杭锦旗境内部分以风力侵蚀为主，鄂尔多斯市东胜区境内部分则以水力侵蚀为主，准格尔旗与达拉特旗境内则是风力、水力侵蚀面积各占一半，这种侵蚀特征可能与县区所覆盖的侵蚀地貌息息相关，杭锦旗境内部分靠近中游库布齐沙漠，所以侵蚀基本以风蚀为主，鄂尔多斯市东胜区完全由砒砂岩丘陵沟壑覆盖，故以水力侵蚀为主，准格尔旗与达拉特旗则包含两个以上的侵蚀地貌，所以风力、水力侵蚀面积别较小；(3)强烈及以上土壤侵蚀主要属于水力侵蚀，基本发生在达拉特旗和鄂尔多斯东胜区境内。

图6 十大孔兑各县区不同土壤侵蚀强度面积堆积柱状图Fig.6 Accumulation column chart of different soil erosion intensities of Ten Kongduis in each county

图7 十大孔兑各县区不同侵蚀强度面积堆积柱状图Fig.7 Columnar chart of accumulation area of different erosion intensities of Ten Kongduis in each county

6 结论

本研究基于多源数据结合多种土壤侵蚀模型对2021十大孔兑区域的土壤侵蚀做出评估，结果表明2021度十大孔兑区域的土壤侵蚀面积为4398.85 km2，占区域总面积的40.86，其中水力侵蚀和风力侵蚀面积分别为2461.56 km2和1937.29 km2。侵蚀强度上，轻度、中度、强烈、极强烈和剧烈土壤侵蚀面积分别为2614.38 km2、1328.74 km2、276.29 km2、194.45 km2和40.15 km2，整体以中轻度侵蚀为主，研究进而从多角度分析十大孔兑区域的土壤侵蚀特征，最终得到如下结论：

(1)十大孔兑区域的水力侵蚀离散分布于孔兑上游，而上游属于砒砂岩丘陵沟壑区，虽然土地利用类型多为林地和草地，但该区域的坡度变化大、土壤可蚀性较高且植被覆盖度相对较低，使得该区域土壤抵抗降雨侵蚀的能力不足，最终导致了孔兑上游严重的水力侵蚀。所以十大孔兑水力侵蚀的分布是由区域地形地貌、植被覆盖等多因子共合同作用导致的，其中地形地貌对侵蚀分布的影响相对较大。

(2)十大孔兑区域的风力侵蚀聚集分布于孔兑中下游，且集中在沙地和耕地上，孔兑中游的库布齐沙漠区域，植被覆盖度极低，土壤抵抗风蚀能力较，是孔兑中风蚀最为严重的区域，其次是下游的冲洪积平原，但由于耕地的植被覆盖度较高，拥有一定抵抗风蚀的能力，所以侵蚀强度上基本轻度侵蚀。十大孔兑区域风力侵蚀的分布主要受到区域土地利用和植被覆盖度的影响。

(3)十大孔兑区域的土壤侵蚀集中在达拉特旗、东胜区和杭锦旗境内，且由于各县区部分所覆盖的侵蚀地貌不同，各县区在风力和水力侵蚀面积比例上有很大异，东胜区境内的部分被砒砂岩丘陵沟壑区覆盖，土壤侵蚀以水力侵蚀为主，杭锦旗境内的部分靠近中游的库布齐沙漠，故以风力侵蚀为主，达拉特旗和准格尔旗境内的部分包含多种地貌类型，所以风力、水力侵蚀面积相不大。