黄土高原是世界上面积最大的黄土堆积区，黄土厚度大多在100~300 m左右^[1]。严重的水土流失和由此引发的各种生态危机，对构建和谐社会造成了严重障碍^[2]，严重威胁中国华北地区的防洪安全和生态安全，而且导致黄土高原成为中国有名的经济落后地区^[3]。宁夏固原市位于黄土高原的西北边际，生态环境脆弱，自然灾害频发，水土流失严重^[4]。中国经济社会发展进入新时代，水土保持面临新的更大的战略地位与责任^[5]，在新时代背景下黄土高原水土流失治理也进入了新的时期，面临新的问题^[6]。在固原市水土流失面积大、分布广且人力、财力和物力投入有限的情况下，确保水土流失治理投入收益最大化、水土保持生态效益充分发挥和深入剖析水土流失与人居环境的内在关系成为了地区水土流失治理和生态环境修复中亟待解决的关键问题之一。

以往研究在水土流失风险评价方面结合USLE^[7]，RUSLE^[8]，CSLE等^[9]土壤流失方程，利用遥感、GIS等技术对区域水土流失风险进行评估，定量分析水土流失风险时空分异特征。对于生态风险与人居环境耦合关系的研究主要聚焦在生态环境与社会经济系统^{[8, 10]}、生态风险与城镇化发展等^[11]方面，但较少针对水土流失风险等具体因素与人居环境的关系来探究生态环境与风险对人类自然生存条件的关系。对于水土流失风险与人居环境耦合关系的研究主要是基于风险扰动的适应性分析框架，采用风险—适应能力指数评估区域人地系统适应性，分析水土流失风险以及人地系统适应能力时空演化特征^[12]，但由于人地系统是经济—社会—自然组成的复杂综合系统，研究中对人地系统适应性评估进行了简化，且分析仅由5个指标计算得到的水土流失风险与经济—社会—自然组成的复杂综合系统的关系，其分析结果的准确性和精度有待提高。

本研究基于GIS技术，结合PSR模型、加权叠加及层次分析等方法，对固原市小流域水土流失风险和人居环境自然适宜性进行评价，分析两者相互作用及产生的影响，计算两者的相关性、耦合度、耦合协调度，在此基础上充分探究各小流域水土流失风险与人居环境自然适宜性的耦合协调类型，即人、地系统间的相互促进或相互拮抗程度和类型，分析人地系统适应性。以精准治理水土流失为突破点，以提高水土资源条件、人居环境自然条件为目标，以期为构建固原市更深层次、更高质量水土流失综合治理体系夯实高质量发展基础，为同类型地区水土保持精准治理提供借鉴和参考。

1 研究区概况与数据来源 1.1 研究区概况

固原市位于中国西北部黄土高原西部边缘、宁夏回族自治区南部，具有明显的黄土高原特点，生态环境脆弱，自然灾害频发，水土流失严重^[4]，是典型的暖温半干旱性气候。四季分别明显，冬季漫长且寒冷，夏季短暂且干燥，春季与秋季的温度变化快，降雨多集中在夏季，春、秋、冬季降水较少，不同区域的降水差异比较大。年平均气温略微偏低，为6.1 ℃，年平均日照时数超过2 500 h，全年昼夜温差大，年平均降水量492. 2 mm。固原市中部的六盘山由南至北贯穿全市，总体海拔由南至北逐渐降低，海拔在1 300~3 000 m之间。下辖原州区、西吉县、彭阳县、隆德县与泾源县。

1.2 数据来源

数字高程模型(分辨率30 m×30 m)来源于地理空间数据云网站(http://www.gscloud.cn/)；土壤类型数据(分辨率1 km×1 km)来源于联合国粮农组织(FAO)和维也纳国际应用系统研究所(ⅡASA)所构建的世界土壤数据库(harmonized world soil database version 1.1, HWSD)；降雨数据、风速数据、地表粗糙度、植被覆盖度、土地利用类型数据来源于宁夏回族自治区2021年度动态监测管理系统；其中降雨数据、风速数据的时间跨度为1991年1月1日至2021年12月31日，精度为日尺度；人口密度数据来源于由佛罗里达大学地理系和新型病原研究所主持的项目WorldPop(https://www.worldpop.org/)。

小流域划分数据来源于宁夏回族自治区水利厅印发的《宁夏回族自治区小流域划分名录》；水土保持4级区划分数据来源于《宁夏回族自治区水土保持规划(2016—2030年)》。

2 研究方法 2.1 水土流失风险评价 2.1.1 评价指标分级与赋值

本研究根据PSR模型的原理将评价体系分为目标层(G)、系统层(S)、指标层(I)3个层次。目标层(G)表示水土流失风险水平；系统层(S)包含压力指标、状态指标和响应指标3个层次；指标层(I)包含水土流失风险评价的具体指标。PSR(pressure-state-response)模型在20世纪80年代由经济合作和发展组织(OECD)与联合国环境规划署(UNEP)确立为研究环境问题的概念模型^[13]。通过对固原市水土流失现状和特点的分析，主要根据各指标的科学性、实用性、逻辑性、独立性、可获取性和可量化性6个原则并参考RUSLE土壤侵蚀模型因子来选择本研究水土流失风险评价的评价指标^{[9, 14]}。将11个评价指标分为风险源指标和受体重要性指标两类。

为便于计算，将该研究中评价指标进行分级并赋值(见表 1)。

年降雨量、植被覆盖度指标采用相等间隔分类法划分。年均风速、地表粗糙度、受影响人数指标采用自然分类法划分。坡度指标以《土壤侵蚀分类分级标准》^[15]中的面蚀(片蚀)分级指标为主要参考依据作出分级。坡向指标按阴坡到阳坡划分为5个等级^[16]。土壤类型指标划分依据不同土壤的水土流失敏感性进行划分^{[14, 17]}。人类活动和潜在经济损失指标根据土地利用类型和《土壤侵蚀分类分级标准》中的滑坡、泥石流危险度分级进行划分^[14]。

2.1.2 设定指标权重

按照《宁夏回族自治区水土保持规划(2016—2030年)》，固原市共划分为4个水土保持4级区：Ⅰ—六盘山水蚀区，Ⅱ—黄土丘陵沟壑残塬水蚀区，Ⅲ—黄土丘陵沟壑水蚀区，Ⅳ—黄土丘陵沟壑水风蚀交错区。本研究采用层次分析法分别确定4个水土保持区的水土流失风险指数权重。指标权重的大小反映该指标对水土流失风险的影响程度大小。此方法具有实用、简单、方便快捷等特点，其计算结果的准确性能够得到保障^[18-20]。

(1) 建立层次结构模型。递阶层次结构如图 1所示。

(2) 构造判断矩阵，进行层次单排序并进行一致性检验。邀请多位权威专家对同一层次因素进行比较，按照1—9标度法对比较的重要性程度进行赋值，构造出判断矩阵并进行一致性检验。以危险源指标P₁的判断矩阵为例：

① P₁与X构造矩阵并计算权重。

② X₁与X₁₁，X₁₂构造矩阵并计算权重。

③ X₂与X₂₁，X₂₂，X₂₃，X₂₄构造矩阵并计算权重。若各判断矩阵的CR值均小于0.10，可以进行后续分析。

④ 层次总排序与一致性检验。前3个步骤得出的是各组指标相对上一层指标的权重值，最终需要计算措施层因素针对目标层的相对权重，即进行层次总排序。

若CR值小于0.10，则认为该专家的判断具有满意的一致性，在一致性比例均符合要求的前提下，取平均值即得到最终权重，结果见表 2。

2.1.3 小流域水土流失风险计算及等级划分

确定各指标的权重后，在GIS中进行加权求和计算，不同尺度栅格数据加权求和后的栅格大小以最大尺度栅格为准，即利用各评价指标分级后精度为1 km×1 km的栅格数据按照对应权重加权后加和，再对危险源指标与受体重要性指标进行加权求和，得到固原市栅格单元的水土流失风险值。在GIS中使用小流域划分矢量数据对上述评价结果进行区域分析，将各小流域栅格单元水土流失风险平均值确定为该小流域的风险值。最后将得到的以小流域为单元的水土流失风险指数C进行归一化处理。小流域水土流失风险指数C分为5级，各水土流失风险等级特点见表 3。

2.2 人居环境自然适宜性评价 2.2.1 评价指标选取与技术流程

影响人居环境的自然因素众多，但最为根本且决定着其他自然因素、对人居环境自然适宜性起主导作用的，主要包括地形条件、水热气候条件和水文状况，以及综合反映区域自然条件的土地利用/覆盖特征^[21]。为系统评价固原市人居环境自然适宜性，首先建立包括地形、气候、水文和地被条件的4大单要素评价模型，构建人居环境自然适宜性综合评价模型——人居环境指数模型，计算栅格尺度固原市不同地区的人居环境指数，在GIS中使用小流域划分矢量数据对上述评价结果进行区域分析，定量评价固原市小流域人居环境自然适宜性^[22-24]。

2.2.2 栅格单元人居环境自然适宜性评价指标计算

(1) 地形起伏度的提取依据地形起伏度的提取方法^[22-23]，参考固原市地形特征、海拔高度、空间分布，将地形起伏度定义为：

$ \begin{aligned} \mathrm{RDLS}= & \operatorname{ALT} / 1000+\{[\max (H)-\min (H)] \times \\ & {[1-P(A) / A]\} / 500 } \end{aligned} $

(1)

式中：RDLS表示地形起伏度指数; ALT表示在一定区域内的平均海拔, 在本研究中定义为3×3窗口内的9个栅格的平均海拔; max(H)和min(H)分别表示3×3窗口内的最高海拔与最低海拔, 以上参数均由DEM数据借助邻域分析工具计算得到, 将研究区域内坡度≤5°的栅格定义为平地, 并定义3×3的矩形窗口大小, 则P(A)表示窗口内平地的栅格数; A表示窗口内总栅格数, 即为9, 通过领域分析工具统计窗口内平地栅格占总栅格数占比即可得P(A)/A。

(2) 气候条件是影响人类活动和人居环境的重要因子，采用温湿指数模型作为人居环境气候适宜性评价的指标。其计算公式为：

$ \mathrm{THI}=1.8 t+32-0.55(1-f)(1.8 t-26) $

(2)

式中：t是年平均摄氏温度(℃); f是年平均空气相对湿度(%)^[25]。

(3) 依据国家生态环境部生态环境评价标准，采用降水量和水域面积比重来表征区域水资源的丰缺程度。前者表示天然状态下区域自然给水能力的大小，后者表征区域集水与汇水能力的强弱。水文指数的具体计算公式为^[26]：

$ \mathrm{WRI}=\alpha P+\beta W a $

(3)

式中：WRI表示水文指数; P表示归一化的降水量, 通过对气象站点降水量数据进行克里金插值得到; Wa表示归一化的水域面积, 利用固原市NDVI图像, 通过栅格计算器选择-1≤NDVI < 0的栅格得到水域栅格图, 再结合小流域边界图, 计算每个小流域内水域面积占比即可得到归一化水域面积; α和β表示降水量和水域面积的比例权重, 通过专家打分法将α取值0.8, β取值0.2^[22]。

(4) 依据生态环境状况评价技术规范中的植被覆盖指数，构建地被指数来表征地表植被的覆盖状况。其计算公式为：

$ \mathrm{LCI}=\mathrm{NDVI} \times \mathrm{LT}_i $

(4)

式中：LCI为地被指数; NDVI为该单元格的归一化植被指数; LT_i为各土地利用类型的权重, i=1, 2…20, 分别代表耕地、林地、草地、水域、建设用地与未利用地等5大用地类型中的水田、旱地等20类二级土地利用类型, 各用地类型的权重采用专家打分法将土地利用类型划分为5大类, 即耕地、林地、水域、建设用地、未利用土地^[22]。

2.2.3 小流域人居环境自然适宜性评价

人居环境指数由地形、气候、水文及地被4大自然因子构成，表征区域人居环境的自然适宜性，为增进各指标之间的横向可比性，对各单项指标进行了标准化处理。

人居环境指数计算公式为：

$ \begin{aligned} \mathrm{HEI}= & \alpha \times \mathrm{NRDLS}+\beta \times \mathrm{NTHI}+ \\ & \chi \times \mathrm{NWRI}+\delta \times \mathrm{NLCI} \end{aligned} $

(5)

式中：HEI为人居环境指数; NRDLS为地形起伏度; NTHI为温湿指数; NWRI为水文指数; NLCI为地被指数; α，β，χ和δ分别为各自然地理区地形起伏度、温湿指数、水文指数和地被指数对应的权重值分别为0.30, 0.31, 0.14, 0.25^{[22, 27]}。

在对固原市人居环境自然适宜性进行归一化处理的基础上，结合小流域划分数据，利用GIS的分区统计功能计算出各小流域中所有栅格单元的平均值，确定小流域人居环境自然适宜性等级。小流域人居环境自然适宜性等级划分及等级特点见表 4。

2.3 小流域水土流失风险与人居环境自然适宜性耦合

为分析水土流失风险与人居环境自然适宜性的相互作用及产生的影响，需要计算两者的耦合关系。耦合指多个体系或运动形式之间通过各种相互作用而彼此影响以至协同的现象^[28]，耦合协调度反映不同系统间的相互促进或相互拮抗程度^[29]。本文采用标准化后数据计算两者同步性及整体协调发展水平，将两组数据均转化为正向指标，即数值越大水土流失风险越低、人居环境自然适宜性越高。计算公式为：

$ C=f(x)^k \times g(y)^k /[\alpha f(x)+\beta g(y)]^{2 k} $

(6)

$ T=\alpha f(x)+\beta g(y) $

(7)

$ D=\sqrt{(C \times T)} $

(8)

式中：C为系统耦合度; f(x)为研究区水土流失风险归一化值; g(y)为研究区人居自然环境适宜性归一化值; k为调节系数, 可调节区分度, 本研究中k=3; α和β为待定系数, 表示水土流失风险和人居环境自然适宜性的重要性关系, 两者均取0.5; T为水土流失风险与人居环境自然适宜性综合性指数; D为水土流失风险与人居环境自然适宜性的耦合协调度, D值越大, 表示水土流失风险与人居环境自然适宜性协调程度越高。根据研究区实际情况及研究需求, 并借鉴相关学者的研究成果^[30-31], 耦合协调度等级划分标准见表 5。

3 结果与分析 3.1 小流域水土流失风险

如图 2所示，固原市小流域水土流失风险等级为安全和较安全的小流域共139条，所占比例超过45%；中度危险与高度危险的小流域共87条，占28%。总体分布呈现出南部高于北部，西部高于东部，南部中间低东西两侧高的特征。中高风险等级的小流域主要分布在Ⅲ(黄土丘陵沟壑水蚀区中部)、Ⅰ(六盘山水蚀区)和Ⅱ(黄土丘陵沟壑残塬水蚀区)南部交界区。各水土流失风险等级的小流域数量、比例及分布特征见表 6。

3.2 小流域人居环境自然适宜性

如图 3所示，基于地形、植被、气候、水文条件的人居环境自然适宜性评价研究表明，固原市人居环境自然适宜性由高到低依次是原州区、泾源县、彭阳县、隆德县、西吉县。总体分布特征为中部、东北部、西北部较高，西南部、东南部较低，形成了“两低三高”的总体分布格局，即2个低适宜性区和3个高适宜性区。

人居环境自然适宜性等级为不适宜的小流域主要集中在西吉县南部，东至唐家河小流域，西至芦子滩小流域，北至白台小流域，南至渝河中游小流域，共计27条不适宜小流域，临界适宜小流域分布在其周边，形成了西吉县西南部低人居环境自然适宜性小流域集中连片区。该区域位于葫芦河西吉段下游河道周边，主要受植被覆盖指数影响，人居环境自然适宜性较低。因此，西吉县是固原市人居环境自然适宜性等级最低的县区，除北部葫芦河西吉段上游、西吉县城区及周边区域外，全县超50%的小流域处于不适宜区和临界适宜区，面积占全市总面积的18%。固原市另一低人居环境自然适宜性小流域集中连片区位于彭阳县西南部、泾源县东北部和原州区东南部区域，包括沙塘小流域等6条不适宜小流域，约占固原市总面积10%。固原市面积最大的中高适宜性小流域集中连片区位于原州区北部与彭阳县东北部。其他中高适宜区分布在南部六盘山山区和西吉县葫芦河流域。另外，一般适宜区位于以上两个低适宜性小流域群和3个高适宜性小流域群的中间过渡区域。

3.3 水土流失风险—人居环境自然适宜相关性及耦合度

如图 4所示，固原市小流域水土流失风险与人居环境自然适宜性线性回归方程为，其中相关系数R=0.61，表征两组变量呈现中度相关关系。10周期移动平均趋势线和线性回归方程显示，水土流失风险越高，人居环境自然适宜性等级越低，两组变量呈现负相关。

耦合度表征双方相互作用的强弱程度，不分相互促进或抑制。将两组数据均进行正向化处理，即数值越大自然状况越好，再进行耦合度计算。如图 5耦合度统计结果所示，74%的小流域水土流失风险与人居环境自然适宜性耦合度为高水平耦合与良性共振耦合，仅有13%的小流域处于拮抗阶段和分离阶段。由此可见，固原市小流域水土流失风险与人居环境自然适宜性呈现较高的耦合性。

3.4 耦合协调度与耦合协调类型

如图 6所示，固原市64%的小流域耦合协调度等级为协调发展和优质协调，在全市各个区县均有分布，主要集中在北部绝大部分区域及中南部六盘山区且两区域相互连接，其中优质协调小流域集中区位于研究区东北、西北及中南部，小流域数量占全市23%。现有两个耦合协调度等级为濒临失调及以下的低耦合度区域，分别位于研究区西南部与东南部，两区由六盘山阻断，小流域数量占比为36%，其中耦合协调度等级为极度失调的小流域集中区位于西南部低耦合协调度区的中心位置，小流域数量占比为5%。综合考虑各县区不同耦合协调度小流域数量及占比，全市四县一区耦合协调度由高到低分别为原州区、彭阳县、泾源县、隆德县和西吉县。

如图 7所示，固原市共有两个低耦合协调度区，分别位于西南部和东南部，其中西南部低协调度区小流域主要类型为失调共损型和失调环境适宜性滞后型，东南部低协调度区小流域主要类型为失调水土流失风险滞后型。协调同步型和协调环境适宜性滞后型小流域在全市分布范围大且无明显聚集区，协调水土流失风险滞后型小流域主要集中分布在东北地区，少量分布在南部地区。

耦合协调类型在固原市各个县区的数量比例如图 8所示，协调同步型与失调共损型分别在协调区和失调区比例最高，协调区中人居环境自然适宜性相对滞后，失调区中水土流失风险相对滞后；协调同步型小流域主要分布在彭阳县(28条)、原州区(27条)，失调共损型小流域主要分布在西吉县(16条)、隆德县(10条)。协调水土流失风险滞后型分布在原州区(11条)、彭阳县(24条)、泾源县(8条)，隆德县和西吉县没有协调水土流失风险滞后型小流域，失调水土流失风险滞后型小流域在彭阳县(12条)、隆德县(9条)和泾源县(9条)分布较多。

失调型小流域主要分布在西吉县南部与隆德县西北部交界处、泾源县东北部和彭阳县西南部交界处，形成了固原市西南部和东南部两大失调型小流域集中区。黄河支流渭河的重要支流葫芦河由北至南贯穿固原市西吉县，由于西吉县城区位于葫芦河西吉段上游段，西吉县北部一直以来都是水土流失防治和生态环境保护与修复的重点区域，流域水土流失得到有效控制，生态环境状况良好。但葫芦河西吉段下游河道周边主要受植被覆盖指数影响，人居环境自然适宜性较低；由于该区域坡度较小且受影响人数较少，所以水土流失风险相对不高，因此沿葫芦河西吉段下游河道周边形成了失调环境适宜性滞后型小流域集中区，在该区域外围，水土流失风险受坡度影响而升高，外围失调共损型小流域呈环状分布，因此形成了西南部失调型小流域集中区。东南部失调型小流域集中区位于六盘山宁夏段东部，该区域虽植被覆盖度高、水文状况良好，但地表起伏较大，导致人居环境自然适宜性等级为一般适宜及以下；由于该区域坡度较大且生态服务价值较高，导致水土流失风险与人居环境自然适宜性相比更加滞后，因此形成了六盘山宁夏段东部失调水土流失风险滞后型小流域集中区。

4 讨论与结论 4.1 讨论

本研究详细阐述了水土流失风险评价中指标选取、指标分级、权重赋值等过程，形成了以PSR模型、加权叠加法及层次分析法为核心的水土流失风险评价方法，将生态风险评价细化到了水土流失风险评价，是对黄土高原水土流失风险评价的补充。将由地形、气候、水文和地被条件计算得到的人居环境自然适宜性作为水土流失风险与人地系统关系探索中的要素。在上述评价结果的基础上以小流域为单元对水土流失风险与人居环境自然适宜性进行耦合分析，探究两者的耦合度、耦合协调度和耦合协调类型。水土流失与环境适宜性耦合关系探索是对西北黄土高原地区进行人地关系研究中水土流失与人地系统适应关系领域的补充，该研究方法适用于大部分地区。其中人居环境自然适宜性评价结果与陈红翔等^[32]对宁夏南部山区生态脆弱性评价结果、朱志玲等^[33]对宁夏生态环境敏感性评价结果具有较高的一致性；水土流失风险评价结果与朱志玲等^[33]、郑雪慧等^[34]对宁夏土壤侵蚀敏感性评价结果以及水利部发布的2021年宁夏回族自治区土壤侵蚀动态监测结果均具有较高的一致性。

本研究的创新点在于将水土流失风险与人居环境自然适宜性评价的单元由行政区域、河流流域等宏观尺度单元或栅格等小尺度单元扩展到以小流域为评价单元。对黄土高原人口稀少地区水土流失风险与人居环境自然适宜性进行耦合关系分析以表达人地关系也是本研究的重要创新点和对该研究领域的补充。此外，因以往研究主要探讨的是城镇化水平、社会经济系统等与生态环境的关系，故研究区往往是城镇或人口密集区域，对农村地区及人口稀少地区的研究较少，本文对西北黄土高原地区进行人地关系研究是对水土流失与人地系统适应关系研究领域的补充。

水土流失风险和人居环境自然适宜性评价是基于规律性的评价，其本身存在较强的推测性质，且本研究未对评价结果进行实地调查验证，可能会小幅度影响评价结果及耦合关系分析的精确性。由于风险指数(C)、人居环境指数(HEI)和耦合协调度(D)等数据均进行标准化处理，评价结果仅具有相对性，无法显示其水土流失风险和环境适宜性的绝对大小，即本研究中水土流失风险和人居环境自然适宜性只能反映小流域在固原市所有小流域中的排名或相对高低，导致评价结果与更大尺度的评价结果存在一定差异。此外，在水土流失风险和人居环境自然适宜性评价中对各因子的权重赋值及影响因子的选择仍然需要进一步研究。

4.2 结论

(1) 固原市水土流失需要优先治理和重点防护的小流域包括杨白虎等87条小流域，即水土流失风险为中度危险与高度危险小流域；人居环境自然适宜性等级为不适宜的小流域主要分布在西吉县南部和隆德县东北部，包括聂家河等36条小流域。

(2) 固原市小流域水土流失风险与人居环境自然适宜性具有较高的相关性和耦合度，即两个体系之间通过各种相互作用而彼此影响以至协同。耦合度分布特征为东北至西南逐渐降低。

(3) 固原市小流域水土流失风险与人居环境自然适宜性具有明显的相互促进作用，但在13%的小流域中两系统存在相互拮抗现象，低协调度小流域包括玉桥小流域等39条，主要集中分布在固原市西南部与东南部。

(4) 潘沟等79条协调、失调水土流失风险滞后型小流域需要进行针对性水土流失治理以免造成水土流失对人居环境的破坏。水洛河等36条失调共损型小流域是固原市自然环境最恶劣的区域，应结合地区实际需要对该区域自然环境进行重点治理和改造，从而降低水土流失风险，提高人居环境自然适宜性。