景观生态风险是指在自然因素或者是人类社会因素的干扰下, 不同因素相互交错作用而对景观生态环境可能产生的不利影响^[1]。流域景观生态风险评价作为区域生态风险的重要分支领域, 日益成为"景观格局——生态过程关联"研究的热点区域^[2]。早期关于生态风险评价的研究主要集中于对不透水面增加, 农业污染源的扩散, 土壤侵蚀等特定风险源进行评价, 缺乏全局视角的研究^{[1, 3]}。而流域生态系统作为一个特殊的自然综合体, 不仅受到特定景观格局的影响, 同时自然与社会因素也会对其产生作用, 因此有必要从多维度、多级别以及多尺度对流域景观生态风险进行综合分析^[1]。目前对于流域景观生态风险评价的研究正在向多重受体和多样风险源的的模式发生转变^[4], 但是从"自然—社会—景观格局"3个维度对流域景观生态风险进行综合评价的相关研究仍不足^[5]。因此, 将景观格局指数与流域"自然—社会"多源要素进行整合, 并将其涵盖进流域景观生态风险评价体系中, 可以为后续的景观格局优化研究工作提供科学依据^[6]。近年来, 国内外的专家学者对不同尺度上的景观格局生态风险评价的有关内容开展了大量的研究, 研究对象包括行政区域^[7]、城市用地^[8]、高原^[9-10]、自然保护区^[11]等。研究方法囊括生态模型、数字模型、景观生态模型和数字地形模型等^[12]。这类模型方法大多数以行政区为研究单元, 评价结果难以反映区域景观生态风险的空间差异特征, 更难反映区域景观生态风险的空间邻近关系^[13]。层次分析法^[14]和基于"源—汇"理论的综合指数构建的方法^[15]也被广泛采用, 但是此类方法对于权重值的确定具有很强的主观性, 不利于在众多影响因素中识别主要的风险要素^[16]。而本研究所采用的空间主成分分析法能够客观的计算出每个影响因素的权重, 避免多指标判断权重的不完善性, 且能够将主成分因子的分析结果落实在地理空间所对应的每一个栅格上, 具有良好的空间可视化效果。

景观生态风险评价与景观格局优化两者之间关系密切, 进行生态风险评价能够反映出区域生态安全状况, 可以为区域景观格局优化提供方向；而进行景观格局优化则需要借助人类对生态环境的主动干预来实现各要素之间的优化配置, 减弱生态风险, 最终促进区域生态安全格局的构建^[17]。最小累积阻力模型(Minimum Cumulative Resistance Model MCR)最早是由Knaapen用于景观格局优化^{[6, 18]}。根据景观生态学相关原理, 在具有异质性特征的生物扩散过程中, 物质能量的流动、交换以及生境的维护都需要克服景观阻力^[19]。MCR模型的优势在于它综合考虑了景观单元内部单元之间的水平联系, 而非景观单元内部的垂直过程, 即考虑了距离和阻力两种因素, 在景观格局优化的过程中具有良好的适用性和延伸性^[20-21]。为此, 本文基于MCR模型, 结合网络分析方法, 更好的体现区域各景观之间的有机联系。

流域作为特殊的一类景观生态系统, 对特定区域的水循环、生物多样性以及生态环境有着重要的影响^[15-16]。涪江流域作为长江流域上游重要的水源地以及生态屏障, 也存在一些问题, 如水土流失严重, 自然灾害频发, 植被覆盖率减弱等, 维护其生态环境的健康发展对于整个川渝地区, 乃至整个长江流域的生态安全以及水源补给安全都有着重要的作用。因此, 对此区域进行景观生态风险评价以及景观格局优化的需要显得尤为迫切。

1 研究区概况、数据与方法 1.1 研究区概况

涪江为长江的二级支流, 嘉陵江右岸最大支流, 发源于四川省境内的雪宝顶, 流经四川省的中部, 自西北向东南斜穿过四川盆地, 在重庆市合川区汇入嘉陵江, 如图 1所示, 整个流域介于29.10°—33.04°N, 103.30°—106.30°E之间, 流域总面积为3.92×10⁴ km², 流域涉及了川渝地区的30个县区。区内海拔高差较大, 高达5300 m, 上游经过川西北丘状高原山区, 地势起伏较大；中下游流经四川盆地盆中丘陵区, 地势较为低平。主要的气候类型为川西高原气候和亚热带湿润季风气候, 降水丰沛, 雨热同期, 但时空分异特征较为显著。整体上研究区景观格局空间分异特征明显, 上游主要以林地为主, 植被覆盖率较高, 中下游地区多以耕地和建设用地为主。由于近些年自然环境的整体变化和人类活动的频繁的影响, 流域内滑坡、崩塌、泥石流等自然灾害频发, 导致区域内景观格局愈加的破碎, 流域所面临的景观生态风险问题日益严峻。

1.2 数据来源

研究数据包括地表覆被数据、DEM数据、土壤类型数据以及其它基础地理数据。①地表覆被数据：2020年GlobeLand 30数据, 空间分辨率为30 m, 来源于http://www.globeland30.org, 总体精度为85.72%, Kappa系数0.82, 满足研究需求。②DEM数据：采用GDEMV2数字高程数据, 空间分辨率为30 m, 来源于中国科学院计算机网络信息中心地理空间数据云平台(http://www.gscloud.cn/)。③土壤类型数据：来源于世界土壤数据库(HWSD)的中国土壤数据集(http://www.fao.org/home/en/), 空间分辨率为1000 m, 比例尺为1 ∶ 100万。④其它基础地理数据：包括居民点数据、工矿用地以及水体数据等, 来源于全球地理信息资源目录服务系统(https://www.webmap.cn/), 比例尺为1 ∶ 25万, 为shp格式文件。通过ArcGIS平台的距离计算工具进行栅格化处理, 栅格大小为30 m。

1.3 研究方法 1.3.1 生态风险评价因子选取

根据涪江流域的生态环境现状以及经济发展情况, 从自然、社会、景观格局3个维度选取了高程、坡度、土壤类型、地表覆盖类型、距河流距离、距工矿用地距离、距居民点距离、香农均匀度指数(SHEI)、蔓延度指数(CONTAG)、植被覆盖度指数(NDVI)等10个指标因子构建流域景观生态风险评价指标体系, 植被覆盖度指数(NDVI)在ENVI平台中基于Landsat8 OLI_TIRS_2018年的遥感影像计算得出；香浓均匀度指数(SHEI)和蔓延度指数(CONTAG)通过Fragstats软件中移动窗口法进行了可视化处理^[22]。其中土壤类型的等级划分根据所对应的固碳能力和抗外力侵蚀能力的差异, 能力越低的区域, 景观生态风险等级越高^[23]；地表覆盖类型参考文献[6]进行分类；其余评价因子利用自然断点法进行分级, 1—5级分别表示较低、低度、中度、较高和高度风险, 分级详情见表 1。生态风险评价因子等级空间分布见图 2。

1.3.2 空间主成分分析

用于评价生态风险的方法中, 空间主成分分析法是将输入的多波段数据变换到一个新的空间中, 对原始空间轴进行旋转而成新的多元属性空间^[24]。通过空间主成分分析法(Spatial principal components SPCA)计算得出空间载荷图、各个成分的贡献率以及累计贡献率, 将累计贡献率超过90%的成分确定为有统计学意义的主成分^[25], 再根据数理统计方法计算每个影响因子的权重, 然后进行加权叠加, 最终得到评价结果。具体的公式表达如下^[26]：

(1)

式中, E表示景观生态风险评价综合结果；a_ij表示第i个栅格对应的第j个主成分；F_j表示第j个主成分的特征贡献率。

1.3.3 景观格局优化

根据俞孔坚^[27]、刁菲菲^[28]等对于景观生态安全格局构建的方法, 本研究根据生态源地以及景观生态风险评价的结果构建景观格局累积阻力表面, 利用最小累积阻力模型(MCR)和网络分析等建立生态廊道、识别生态节点, 实现景观格局优化。公式如下：

(2)

式中, MCR表示生态源地j到任意栅格i之间的最小阻力累计值；D_ij表示景观格局阻力上第i个栅格到第j个生态源地之间的距离；W_i表示在景观格局阻力表面上第i个栅格对生态流运行所产生阻碍的阻力值。

(1)"生态源地"的识别。生态源地是指对区域生态环境以及人类社会可持续发展具有重要作用的斑块, 其能够提供关键的生态服务以及防止生态系统退化而引发的各种生态问题^[29]。本研究对于生态源地的选取为初次筛选的面积大于5 km²的林地、湿地以及水体斑块导入Conefor_input, 计算得出景观一致性概率(Landscape coincidence probability LCP)、整体连通性(Integral index of connectivity IIC)、可能连通性(Probability of connectivity PC)、斑块重要性(DI)等4个景观格局指数, 阈值设为2000, 连通概率设为0.5, 对核心区进行连接度评价, 并将核心区域重要值>1的斑块确认为最终生态源地。

(2) 生态阻力面的确定。生态阻力面是生态源地之间进行物质交换、能量传输以及物种迁徙等生态过程中所受到的阻力^[30]。以流域景观生态风险评价的结果和筛选出的生态源地作为景观格局阻力面生成的基础, 利用成本距离(Cost distance)工具计算得到涪江流域景观格局累积阻力表面, 最后用自然断点法将综合阻力表面的大小分为1—5个等级：低阻力、较低阻力、中阻力、较高阻力、高阻力。景观格局阻力表面的分级标准见表 2。

(3) 生态廊道的建立。廊道是指具有通道或屏障功能的线状和带状的景观要素, 是"源地"之间能量流通的重要桥梁和纽带^[31]。对于生态廊道的提取本文主要是基于选取的"生态源地"以及景观格局累积阻力表面利用Model Builder工具进行建模, 计算出每个生态源地到其余生态源地之间的最短耗费路径, 然后根据长度将其分为5个等级。(表 3)

(4) 生态节点的判别。生态节点是指在整个景观格局中需要加以识别和保护相对脆弱的区域, 利用水文分析的相关工具提取景观格局阻力面的"山脊线", 将"山脊线"与生态廊道相交, 提取生态节点^[32-33]。

(5) 景观格局优化网络构建。为了进一步促进生物流在研究区流通, 以构建的生态廊道和识别的生态节点为基础, 利用网络分析法, 再次计算最短路径, 形成了一条贯穿整个流域的主廊道, 模拟和分析能量在生态廊道网络中的流动和分配情况。

(6) 景观格局优化效果评估。景观连接度可以有效的对景观要素在景观功能和过程上的有机联系进行描述^[6]。前人的研究大多选取景观一致性概率(LCP)、种类一致性概率(CCP)、整体连接度指数(IIC)和可能连接度指数(PC)等对景观格局优化前后效果进行评估^{[6, 34-36]}。其中IIC和PC指数能够较好的表现景观格局斑块与重要生态要素之间的连通情况, 因此本研究选取IIC指数和PC指数对涪江景观格局优化前后效果进行定量分析, 表达式如下^[37]：

(3)

(4)

式中, a_i与a_j分别表示任意斑块i和斑块j的面积；A_L表示流域总面积；NL_ij表示斑块i和斑块j之间的连接数。IIC的取值介于0—1之间, 值越大表示景观格局连接性越高。P_ij^*表示斑块i和斑块j之间所有可能的路径的最大乘积率, PC的取值介于0—1之间, 并且随着连通性的提高而增加。

2 结果分析 2.1 流域景观生态风险评价 2.1.1 指标权重的确定

选取海拔、坡度、土壤类型、地表覆盖类型、距河流距离等10个影响因子, 通过空间主成分分析得出前7个主成分的累计贡献率>90%(表 4), 即前7个主成分能够对研究区景观生态风险的信息进行有效的概括。根据前7各主成分的初始特征根以及贡献率计算得出每个指标因子的权重值(表 5)。从计算结果可以得出蔓延度指数(CONTAG)、距水体距离、香农均匀度指数(SHEI)、坡度占比较大, 分别为0.1293、0.1280、0.1197、0.1098, 说明其对涪江流域景观生态风险评价结果的影响较大。因此, 从权重角度出发, 涪江流域景观生态风险评价结果受自然因子和景观格局因子影响较大, 社会因子影响相对较小。

2.1.2 涪江流域景观生态风险评价结果

研究区景观生态风险空间分布特征见图 3：①高度生态风险区, 风险指数为3.45—4.51, 面积为5697.60 km², 占了整个流域面积的14.55%, 主要位于涪江流域的上游部分地区以及中下游地区有工矿用地以及居民点分布的的区域, 其主要原因是上游地区海拔较高、坡度较陡、地表起伏较大, 在雨季节易发生崩塌、滑坡、泥石流等自然次生灾害, 从而导致面临的景观生态风险问题较为严峻；②较高生态风险区, 风险指数为3.12—3.44, 面积为9941.48 km², 所占整个研究区面积25.38%, 主要在上游高生态风险区的外围区域分布、中下游人类活动影响较为强烈的区域有零星分布如城镇建设用地和工矿用地等；③中度生态风险区, 风险指数为2.80—3.11, 分布面积为9957.43 km², 所占整个流域面积的25.43%, 在中下游区域零散分布；④低度生态风险区, 生态风险指数为2.44—2.79, 分布面积为8117.26 km², 占流域总面积的20.72%, 主要分布于中度生态风险区的外围地区, 覆被类型主要以耕地、草地以及灌木为主。⑤较低度生态风险区, 生态风险指数为1.30—2.43, 面积为5453.52 km², 占流域总面积的13.92%, 占比最小, 主要位于植被覆盖度较高区域、河流流经区域等, 主要是由于该类区域的生态系统服务价值较高、生态环境稳定且受人类活动干预较少, 有利于自然生物在此类地区聚集繁衍以及维持生态平衡。

2.2 涪江流域景观格局优化 2.2.1 流域生态源地的提取

经过筛选, 选取出15个斑块作为最终生态源地, 占地面积为11194.28 km², 占总流域生态用地面积的28.58%, 以大面积分布的林地为主, 河流以及湿地占比较小, 主要分布于上游地区的平武、北川、安州和江油, 该地区植被覆盖度较高, 生态系统较为稳定, 生物多样性丰富, 有利于物种的扩散和维持。生态源地分布见图 5。

2.2.2 流域景观格局阻力表面

由图 4可见, 涪江流域景观格局累积阻力表面的空间分布特征为：①较低阻力地区所占整个研究区阻力面的范围最大, 面积为18143.05 km², 占了46.24%, 主要分布于上游区域植被覆盖度较高地区和涪江干流流经区域等。②低阻力区占了21.93%, 面积为8606.351 km², 占地面积较大, 主要位于梓潼、涪城、中江、射洪等区县。③中等阻力区所占面积为7189.816 km², 占了总面积的18.32%, 主要分布在低阻力区域的外围地区。④较高阻力区主要在分布研究区城镇建设用地的外围区域, 占比为10.39%, 面积为4077.067 km²。⑤高阻力区占比为3.11%, 面积为1220.831 km², 主要分布在中下游地区城市建筑较为密集的区域。

2.2.3 流域生态廊道的构建与优化

从图 5和表 6可见, 共识别出研究区的生态廊道41条, 长度共有5229.038 km。其中原有的生态廊道29条, 新添12条潜在生态廊道, 根据廊道经过的生态用地类型, 将其分为绿带型、河流型以及道路型三类^[16]。其中：①一级生态廊道有3条, 总长度为1079.234 km, 所占比例为20.64%, 最长的一条为22号廊道, 从四川省的九寨沟县一直延伸到重庆市的合川区, 往后依次为26号以及20号廊道。三条一级生态廊道大致贯穿了整个涪江流域, 有利于物种在流域内的迁徙和物质能量在流域范围内的交换, 是上下游维持生态平衡的关键, 需要给予重点关注。②二级廊道有5条, 总长为1392.163 km, 占了生态廊道总长度的26.62%。分别为14号、10号13号、6号和18号生态廊道。二级生态廊道贯穿了涪江流域的中下游区域, 呈闭合状, 有效的连接了中下游的生态源地, 促进当地物种的迁移和保护, 提高生物多样性, 增加生态服务价值。③三级生态廊道共有13条, 总长度为1854.142 km, 长度占比为35.46%, 是占比最大的一类生态廊道, 涉及4号、39号、19号、17号、30号、34号、9号、27号、35号、36号、28号、3号和41号廊道。三级生态廊道虽然比较破碎, 但是涉及到了涪江流域的上、中、下游, 有效的连接了流域内分散分布的面积较小的生态源地, 增加研究区内生态网络的连通性, 促进了较为完整的生态网络结构体系的构成, 有利于物种的保护。④四级生态廊道共有14条, 其总长度为830.744 km, 占所有廊道总长比例的15.89%, 其中21号、15号、24号、5号、23号、25号、12号以及7号此8条廊道在平武、北川、安县以及江油境内, 在以上4个区县内形成了较为完整的廊道网络, 满足了生物沟通的需要, 因此, 此区域内自然保护区众多, 生物多样性较为丰富；其余的37号、32号、40号、2号、40号、37号和1号则较为分散的分布于下游区域, 没有形成较为完整的廊道网络体系, 生态连接度较差, 无法满足生物沟通的需要。⑤五级廊道有6条, 总长度为72.755 km, 仅占所有廊道长度的1.39%, 为占比最小的廊道, 分别为8号、11号、16号、29号、38号和33号廊道。

研究区内的生态廊道类型主要是以绿带型和河流型廊道为主, 主要受生态源地分布的影响。对于绿带型廊道可以进一步增加林地型生态源地的面积并且增大其之间的连通性, 减小景观格局阻力。而对于河流型的生态廊道, 为了进一步减少沿江一带的景观格局阻力, 可以根据河流型廊道的宽度建立缓冲区, 在缓冲区边缘设置绿化带, 促进能量流和生物流在此类型廊道中的流通。

2.2.4 流域生态节点的提取

通过提取景观格局累积阻力面"山谷线", 将其与生态廊道进行叠加分析, 提取出53个生态节点(图 5)。根据其分布的具体位置可以看出, 位于河流流经区的生态节点有17个, 位于生态源地中的生态节点有26个, 位于生态源地之间的生态节点有10个。根据生态节点的概念, 将生态节点分为了3类, 分别为河流湿地区生态节点、重要生态区生态节点和生态脆弱区生态节点, 前两者主要位于斑块面积较大的生态源地内和廊道与廊道交接处, 生物流的流量较大, 在形成的廊道中起到了"踏脚石"的作用；而生态脆弱区的生态节点一般距离生态源地较远, 位于前两类生态节点之间, 是研究区形成生态网络的关键, 因此应该加强对生态脆弱区生态节点的保护和建设, 让其对周围较为脆弱的景观产生辐射作用。

2.2.5 涪江流域景观格局优化网络构建

利用网络分析方法, 形成了一条贯穿整个流域范围的主廊道, 长约1607.38 km(图 5)。因此优化过后的景观生态格局是以主生态廊道为"中轴"、识别出的现有的生态廊道以及潜在生态廊道为"辅助"、生态节点为"枢纽", 打通整个廊道, 形成了涪江流域较为完整的生态网络框架。

2.2.6 优化效果评估

对比涪江流域景观格局优化前后的IIC指数和PC指数(表 7), 可以得出：在优化之前两个指数都维持在0.48左右。在优化后, 研究区的IIC指数和PC指数分别为0.083624、0.079675, 变化率分别呈74.12%和67.87%的正增长。根据相关研究表明^{[6, 36]}, 景观连通性指数较高的区域, 其生态景观格局具有更高的稳定性。因此, 进行景观格局优化有利于提高涪江流域整体景观格局的稳定性。

3 讨论

当前学界对于多源景观生态风险评价的关注逐渐增大, 但是从"自然-社会-景观格局"3个维度综合选取评价因子对景观生态风险进行研究的案例还相对较少^{[1, 4, 19]}。景观格局优化是景观生态风险评价的最终目标^[38]。本文基于景观生态风险评价的结果, 进一步明确了涪江流域目前的景观格局现状以及所面临的生态风险问题, 为后续的景观格局优化奠定基础。在景观格局优化的研究中, 生态节点的数量、质量和空间分布状况会对物种迁移的时间、周期和成功率产生重要的影响, 加大对流域生态节点的建设和保护, 能够促进物质、能量流在整个研究区生态系统内更加流畅的循坏运转, 对于整个流域内生态环境质量的提升和生物多样性的保护具有极其重要的意义^[39]。利用网络分析方法构建的网络生态结构提高了整个区域的连通度指数, 景观生态流可以直接渗透到景观格局阻力较高的区域, 促进研究区生态安全水平的提高。因此应做好上游地区如九寨沟、雪宝顶、王朗等国家级自然保护区的天然林保护工作, 推进流域内森林、草原以及河流湿地生态系统建设工作, 有效的保护和建设好生态网络, 进一步推进研究区内景观生态安全格局的构建。与常规的流域景观格局研究方法相比较, 本文创新之处在于：①从景观生态风险源的选取上, 综合的从自然、社会、景观格局三方面选取评价因子, 并且利用SPCA方法计算得出每个评价因子的权重值, 可以减少由于人为经验和判断所赋予的权重值的主观性, 较为客观的反应出不同评价因子之间重要性差异, 增加评价结果的准确性, 进一步揭示景观生态风险评价结果中所蕴含的生态学内涵。②生态源地的选择是进行景观格局优化的重要环节, 在以往的众多研究中, 生态源地仅仅只是以面积大小为决定因素, 选取的林地、水域等具有生态服务功能的景观类型作为生态源地, 具有较强的主观性^{[6-7, 17]}。而在本研究中对于生态源地的选取以4个表示景观连通度以及斑块重要性的景观格局指数作为依据(LCP、LLC、PC、DI), 最终选取出15个具有较强的稳定性和连通性的生态斑块作为生态源地, 提升对生态源地选择的科学性。③本研究在MCR模型构建出的生态廊道的基础上, 采用了网络分析的方法, 弥补了对于潜在廊道分析提取不足的情况, 使得构建的生态廊道具有更加丰富的层次性, 提取的生态节点在空间定位中更加的清晰, 并且通过计算景观格局优化前后得出的LCP、CP景观指数结果可以看到, 优化后的景观连通度得到了极大的提升。

本研究的所采用的研究方法和研究结果可为涪江流域制定景观格局优化措施提供具有参考意义的布局模式, 并且对于整个长江流域未来景观生态安全构建具有一定的理论和实践意义。但是, 本研究只是得出了优化布局的结构, 对于具体的优化措施, 如构建生态廊道的宽度, 生态节点的具体保护措施等方面, 在未来的研究中都将结合研究区的实际情况进行更深一步的研究和探讨。

4 结论

本文基于空间主成分分析法(SPCA)评估研究区的生态风险等级, 利用最小累积阻力模型(MCR)和网络分析等方法对涪江流域的景观生态结构进行优化, 得出以下结论：

(1) 通过分析得出, 涪江流域的景观生态风险状况在"自然—社会—景观格局"三方面因素的共同作用下, 受自然和景观格局因子的影响较大, 社会因子影响较小。

(2) 从整体上来看, 涪江流域的景观生态风险等级在空间分布上大致呈西北部高于东南部, 生态风险等级为中等、较高、高度的区域的面积占研究区总面积的65.35%, 研究区生态环境所面临的生态风险较高。

(3) 生态源地主要以林地、河流以及湿地为主, 主要分布松潘、平武、北川区县, 总面积为11194.28 km², 占流域总面积的25.58%；构建涪生态廊道共41条, 其中现有的生态廊道有29条, 潜在的生态廊道有12条；识别出生态节点53个。

(4) 构建了以主廊道为"中轴", 现有廊道及潜在廊道为"辅助", 生态节点为"枢纽"的较为完整的网络型生态景观格局；对比优化前后的连通度进行比较, 优化后景观格局的连通度有了极大的提升。