Dynamic simulation of landscape ecological risk in mountain towns based on PLUS model

安宁市位于云南省昆明市西南(24°31′~25°06′N，102°10′~102°37′E)，辖连然、金方、太平新城等9个街道办事处，总面积1 301.81 km²；以山地为主，平均海拔1 800 m。2019年，全市地区国内生产总值(GDP)575.14亿元，常住人口38.9万人，城镇化率78.3%，是云南省最大的冶金、盐、磷、石化工基地。

作为云南省转型升级样板区和滇中绿色发展示范区以及“一带一路”建设先行先试区，安宁市承担了滇中城市群新型城镇化空间拓展和长江上游生态安全屏障的双重使命。协调好城镇快速发展与生态环境保护的关系，成为了安宁市高质量发展的重要命题。

土地覆被数据来源于2000、2010和2020年Globeland 30全球地表覆被数据库(http://www.globallandcover.com)，空间分辨率为30 m×30 m，区内土地覆被类型包括耕地、林地、草地、灌木地、水体、人造地表等6种。高程数据来自地理空间数据云(http://www.gscloud.cn)；国内生产总值、人口空间分布数据来源于中国科学院资源环境科学与数据中心(http://www.resdc.cn)；夜光遥感数据来源于Luojia1-01数据(http://59.175.109.173:8888/index.html)，并采用ENVI 5.3校正；道路数据来源于OSM数据集(https://www.openstreetmap.org)；河流、湖泊数据来源于全国地理信息资源目录服务系统(https://www.webmap.cn)。

①景观生态风险评价单元划分。利用ArcGIS 10.8对土地利用数据进行1 km×1 km网格化，将安宁市划分为1 310个评价单元；利用Fragstats 4.2软件分别计算每一评价单元的景观生态风险指数值，以此作为样地中心点的生态风险值。②景观生态风险指数构建^[20]。从景观格局的角度出发，借鉴已有研究成果^[21]，构建景观生态风险指数。计算公式为：

其中：I_ERk为第k个样本单元内的景观生态风险指数；A_ki为第k个样本单元内景观类型i的面积；A_k为第k个样本单元的面积；R_i为景观类型i的景观损失度指数。

其中：E_i为景观干扰度指数，通过景观破碎度C_i、景观分离度S_i和景观优势度D_i构建，a、b、c为相应各景观指数的权重，且a+b+c＝1，根据已有研究成果及实际情况^[22-23]，分别赋值为0.5、0.3、0.2；F_i为景观脆弱度指数，结合已有研究成果^[24]，通过专家打分法对6类景观类型赋值，将水域赋为6，耕地赋为5，草地赋为4，灌木地赋为3，林地赋为2，人造地表赋为1，最后进行归一化处理得到各类景观的脆弱度指数(F_i)值，分别为0.29、0.24、0.19、0.14、0.10、0.05。③探索性空间数据分析。运用ArcGIS地统计分析模块对样地中心点数据进行半变异函数拟合，得到最佳拟合模型；选取指数克里金插值法分别对2000、2010和2020年3期数据进行空间插值^[25]，采用自然间断点法划分景观生态风险等级。通过Geoda软件计算全局空间自相关莫兰指数(I)，来表征景观生态风险指数在空间分布模式及其集聚效应^[26]，I值为−1~1。I＞0表示正相关，I＜0表示负相关，I＝0表示不相关。

PLUS模型集成了基于土地扩张分析的规则挖掘方法和基于多类型随机种子机制的元胞自动机(CA)模型2个模块，对各类土地利用变化的影响因素解释性更好，模拟结果精度更高^[18]。模型通过提取2期土地变化中各类用地扩展部分，运用随机森林算法获取各类用地的发展概率，再利用基于多类随机斑块种子的CA模型对未来景观格局进行模拟预测。首先，根据安宁市实际情况及数据的可获取性，从自然因素、社会经济因素和可达性3个方面选取高程、坡度、地形起伏度、公里网格GDP、人口密度、夜光遥感、距道路距离、距城镇距离、距湖泊距离和距河流距离等10个影响因子，栅格化后统一成与土地覆被数据相同的投影坐标系及空间分辨率。其次，利用用地扩张分析策略(LEAS)模块运算得到研究区各景观类型的发展概率。最后，结合未来各类用地的目标像元数、转移成本矩阵、随机斑块种子的概率及邻域因子等相关参数，基于多类随机斑块种子的CA模型实现研究区景观类型变化模拟。本研究运用Markov模型进行未来景观类型的需求预测，在已有基础上^{[16, 27-28]}结合研究区景观类型转移情况，反复调试各参数，确定景观类型转移成本矩阵，将随机斑块种子的概率设定为0.01(参数范围0~1，越接近1表示越容易产生新的斑块)，耕地、林地、草地、灌木地、水域及人造地表邻域因子参数设定为0.7、0.4、0.3、0.3、0.2、0.9(参数范围0~1，越接近1表示景观类型的扩张能力越强)。选择Kappa系数和FOM(Figure of Merit)系数进行模拟结果的精度评估，其中Kappa系数为0~1，大于0.7表示模拟结果一致性较高，精度较高^[29]；FOM系数由实际土地变化与预测的变化的交集与两者并集的比值确定，其值为0~1，值越高表示模拟结果的精度越高^[30]。

2000—2020年，安宁市景观类型以林地、耕地和草地为主(表1)。其中，林地是最主要的景观类型，占研究区总面积的43.46%；其次为耕地，占总面积的24.94%，与社会经济关系最为紧密的人造地表占比9.63%。全区的景观类型变化特征表现为人造地表面积持续增加，耕地、灌木地面积持续减少，林地、草地面积先增加后减少，水域面积先减少后增加。人造地表面积呈上升趋势，耕地、灌木地、草地、水域面积总体呈下降趋势，林地面积较为稳定。其中，人造地表面积增加幅度最大，共增加8 464.36 hm²，占总变化量的49.92%；耕地面积减少幅度最大，共减少5 274.19 hm²，占总变化量的31.10%。2000—2010年，各景观类型面积变化最明显的是林地，增加1 358.30 hm²，占总变化量的33.06%；2010—2020年，各景观类型面积变化最明显的是人造地表，增加8 175.00 hm²，占总变化量的48.89%。

2000—2020年，人造地表扩张速度变化最快，主要由耕地、草地、灌木地转入；2010—2020年土地动态度明显高于2000—2010年，其中2000—2010年减少速度较快的景观类型为灌木地和水域，灌木地主要转移为林地；水域主要转移为耕地、林地和草地；2010—2020年人造地表增加速度明显，城镇的快速扩展占用了大量的耕地、林地和草地。城镇化进程中以人造地表景观扩张为主，并由此产生山水林田湖草等生态景观面积减少，景观结构发生剧烈变化，生态系统稳定性减弱。

由图1可知：安宁市2000、2010和2020年景观生态风险的莫兰指数(I)分别为0.525 4、0.552 1和0.492 9，均表现出显著的空间正自相关性(P＜0.01)，说明景观生态风险值在空间上表现为聚集分布，且相互之间存在影响。2000—2010年，景观生态风险的I值上升，2010—2020年下降，整体呈下降趋势，表明安宁市景观生态风险值在空间分布上的依赖性减弱，空间趋同性逐渐降低，人类活动影响着生态风险指数的空间格局，并导致其空间分布由聚集向均匀过渡。

Figure 1.  Landscape ecological risk Moran’s I scatter diagram in Anning, 2000−2020

运用指数克里金插值法得到安宁市各年生态风险值，基于自然断点法将2020年研究区生态风险值分为低生态风险(I_ER＜0.151 1)、较低生态风险(0.151 1≤I_ER＜0.177 8)、中等生态风险(0.177 8≤I_ER＜0.202 6)、较高生态风险(0.202 6≤I_ER＜0.230 3)、高生态风险(I_ER≥0.230 3)。其他各期数据均采用2020年的分级区间，便于各期数据的比较分析，得到安宁市景观生态风险等级空间分布图(图2)。

Figure 2.  Spatial distribution of landscape ecological risk level in Anning, 2000−2020

由表2可知：2000—2020年安宁市景观生态风险整体以中等风险和较高风险区为主，占研究区景观生态风险等级总面积的52.95%~59.76%。结合图2来看，低生态风险区主要分布在研究区南北两端以及中西部山地地带，景观类型主要为林地、草地，人为干扰较少，景观损失度较低；研究期内面积总体上升，主要为青龙、连然、金方、太平新城等街道组成的中心城区的低生态风险面积增加。这是由于青龙街道在此期间建成国家园林城镇，林地等景观面积增加，生态环境改善；随着城镇化的推进，中心城区城镇拓展有序，人造地表景观集中成片，稳定性强，外界干预影响较少。此外，县街和草铺街道作为安宁市重要的工业基地，为打造绿色工业城市的景观风貌，建设了麒麟现代工业景观分区和职业教育景观风貌分区，景观斑块的连接度增强，景观生态风险得到有效防控，低生态风险区面积增加。较低生态风险区面积持续上升，主要分布在低生态风险区的周围。中生态风险区主要分布于八街和县街街道，以及其余各街道耕地景观周边；该区以耕地景观为主，多分布于山地丘陵地带，景观形态复杂，水网密布，林地、草地、灌木地和人造地表景观等相互交错，破碎化程度较高，生态风险加剧，但研究期内该区整体面积变化相对稳定。较高生态风险区分布较为广泛，主要分布于高生态风险区向中等生态风险区过渡区域，研究期内面积占比降低，共降低6.94%。高生态风险区主要分布于各街道城镇扩展的边缘地带，研究期内总面积呈下降趋势，各时期占比分别为15.28%、14.47%、10.64%；说明城镇快速扩张对各类景观类型干扰度较大，尤其在城市边缘地带建设用地较为分散，人文景观与自然景观交替频繁，打断了生态系统的整体，使得景观破碎度高、分离度明显上升，加大了生态风险；其中禄脿街道整体处于高生态风险等级，作为集磷矿工业、特色农业、交通枢纽、特色农业等多种功能的混合街道，景观破碎度和分离度大，高生态风险区面积较多且呈现长期稳定的高风险状态。应用转移矩阵进一步分析研究区景观生态风险各等级间的转移特征，根据表3可以看出：2000—2020年，较低生态风险区主要转移为低生态风险区，中等风险区主要转移为较低、低生态风险区，较高生态风险区主要转移为中等生态风险区，高生态风险区主要转移为较高生态风险区。各等级生态风险区均主要表现为向次级生态风险区转移的特征。由此可以看出，研究区整体生态安全保障水平有所提高，城市经济发展与环境保护逐渐向着良性发展的态势转变。

为验证PLUS模型对景观类型模拟结果的可靠性，本研究以2010年土地覆被数据为基础，结合各景观类型发展概率，并设置相关参数，得到2020年景观类型空间格局模拟结果，并与2020年景观类型现状进行对比，精度验证计算求得总体精度为89.84%，Kappa系数为0.856，FOM值为0.247，表明模型具有较高的可信度，故基于PLUS模型对安宁市未来景观类型变化进行模拟预测。依据安宁市历史年份的景观类型变化特征，利用Markov模型构建3种不同情景下安宁市2030年景观类型。其中自然发展情景以2020年景观类型现状数据为基础，预测2030年各景观类型的面积及空间分布情况；考虑安宁市正处于快速发展的重大战略机遇期，云南石化治理、滇中新区建设等相关政策进一步推进，势必带来建设用地的快速增长，因此城镇发展情景将耕地、林地、草地、灌木地向人造地表的转移概率增加20%，人造地表向除耕地外的其他景观类型转移概率降低30%；考虑安宁市是长江上游生态安全格局的重要组成部分，为了进一步打造全国“工业城市中生态最好，生态城市中工业最强”城市，未来的景观类型变化中要进一步加强对生态用地的保护，在生态保护情景设置中，将耕地、林地向人造地表的转移概率降低30%，草地、灌木地、水域向人造地表的转换概率降低20%，人造地表向林地的转移概率增加10%，并严禁生境退化，将区域各个水库、湖泊作为约束条件，限制其任意转换。

综合2020年景观类型现状和各影响因子数据，运用PLUS模型模拟预测得到安宁市2030年自然发展情景、城镇发展情景和生态保护情景下景观类型空间分布格局。由图3可知：与2020年景观类型现状数据相比，安宁市2030年3种情景下各景观类型面积呈现不同程度的变化(表4)，人造地表、水域的面积持续增加，耕地、林地、草地、灌木地面积持续减少。具体来看，自然发展情景下，安宁市各景观类型发展趋势发生转变，各街道城镇中心的人造地表均有一定扩张。城镇发展情景下，人造地表面积进一步增加，与自然发展情景相比，八街(农业中心)、县街(城市综合服务中心)等街道明显扩张。生态保护情景下，人造地表扩张幅度明显缩小，相比自然发展情景和城镇发展情景，人造地表分别下降2 076.06和3 704.97 hm²，耕地、林地、草地、灌木地面积减少幅度有所缓和，生态用地得到保护。

Figure 3.  Spatial distribution of the landscape types in 2030 under different scenarios in Anning

由图4可知：2030年3种情景下安宁市景观生态风险空间分布与2020年相似，整体上仍以中等生态风险区、较高生态风险区为主，但高生态风险区和低生态风险区面积均有上升。与2020年相比，自然发展情景下的较低、中等及较高生态风险区面积有所下降，而高生态风险区的面积增长幅度较大(表5)，增加了2 436.66 hm²，低生态风险区也有一定的增加，扩大了1 227.79 hm²。城镇发展情景下，较低生态风险区、较高生态风险区面积下降，较高生态风险区和高生态风险区面积均得到扩张，分别增大2 189.08和2 779.71 hm²；相比其他情景，低生态风险区面积最小，较高、高生态风险区面积最大。这是由于按照现有发展趋势，建设用地面积将会进一步提高，各景观类型均发生变化，人造地表持续外延式扩张，侵占周边耕地、林地等生态景观，各景观的连通性、稳定性下降，导致较高、高生态风险面积增加。生态保护情景下，高生态风险区面积增幅明显下降，低生态风险区显著增加，与其他情景相比，高生态风险区面积最小，景观生态风险等级面积差异显著；主要原因一方面是由于生态保护情景约束了人造地表扩张，减缓了城镇边缘地区景观破碎度、分离度的升高；另一方面该情景设置了水域限制区，保护了水域的稳定性，约束了城镇扩展，增强了新增人造地表空间结构的整合性，进而增强了景观稳定性，使得生态保护区的低生态风险区面积明显上升。

Figure 4.  Spatial distribution of landscape ecological risk level in 2030 under different scenarios in Anning

利用ArcGIS空间分析中的栅格计算器进行图层运算，将生态风险等级降低的区域命名为改善区，生态风险等级提升的区域命名为恶化区，生态风险等级未发生改变的区域命名为稳定区，得到安宁市景观生态风险变化图(图5)。3种情景下，景观生态风险改善区空间分布类似，主要分布在草铺、太平新城、县街、八街街道，零散分散于其他街道。不同情景下恶化区空间分布差异较为显著，城镇发展情景下恶化区面积明显高于其他情景，主要分布于连然、温泉、草铺、县街及八街街道；生态保护情景下恶化区有明显约束，八街街道的恶化区明显少于城镇发展情景，县街街道的恶化区明显少于城镇发展情景和自然发展情景。综合来看，人类活动显著影响未来景观生态风险的变化，在城镇合理发展的同时，对城镇建设用地发展进行管控，对生态保护区有效保护，能显著改善景观生态风险。生态保护情景充分考虑了城镇经济建设与生态环境保护协同发展，有利于形成有序的国土资源空间管控局面，更符合安宁市未来生态安全格局的发展。

Figure 5.  Change of landscape ecological risks under different scenarios in Anning, 2020−2030

本研究通过景观生态风险评价模型分析安宁市2000—2020年景观生态风险时空演变特征，并采用PLUS模型对2030年不同情境下安宁市景观生态风险空间分布特征和变化趋势进行模拟预测，结果表明：①安宁市景观类型以林地、耕地和草地为主，人造地表面积呈上升趋势，耕地、灌木地、草地、水域面积总体呈下降趋势，林地面积相对较为稳定。其中，人造地表面积变化最为显著，主要以耕地、灌木地、草地转移入为主。城镇化进程中以人造地表景观扩张为主，并由此产生山水林田湖草等生态景观面积减少，生态系统稳定性减弱。②2000—2020年，安宁市景观生态风险值空间集聚状态明显，但聚集程度呈下降趋势。景观生态风险等级以中等生态风险区、较高生态风险区为主，低、较低生态风险区面积持续增加，高、较高生态风险区面积持续减少，生态环境状况整体向好；高生态风险区主要分布于各城镇扩展的边缘地带；生态风险的空间分布规律与安宁市人类活动干预强度的空间分布规律密切相关，应重点关注城镇用地扩张所导致的景观类型结构及其生态风险的动态变化。③2030年安宁市3种情景下各景观类型呈现不同程度的变化，景观生态风险空间分布特征与2020年相似，但高生态风险区和低生态风险区面积均有上升。3种情景中，生态保护情景下的高生态风险区面积最小、低生态风险区面积最大、恶化区最少。该情景约束了人造地表的扩张，保护了水域等生态用地，减缓了城镇边缘地区景观破碎度、分离度的升高；使新增人造地表空间结构的整合性增强，提高了景观的稳定性。

安宁市是云南省转型升级样板区、滇中绿色发展示范区。研究期间安宁市景观生态风险等级整体呈下降趋势，表明全国长期坚持生态文明建设和云南省打造世界一流的“绿色能源”“绿色食品”“健康生活目的地”3张牌取得了一定的成效，为城镇生态安全格局的改善提供了良好的条件。中小城镇发展转型过程中，城镇的快速扩张显著影响景观生态风险的演变。不同情境下2030年安宁市景观生态风险空间分布和变化趋势预测模拟中高生态风险区均有明显上升，恶化区多为城镇扩张的边缘区域，控制人造地表对耕地、林地、灌木地的侵占是防控景观生态风险的重要环节。结合研究区现状和未来不同情境下模拟的景观生态风险评价结果，重视城镇在经济转型快速发展过程中不同土地利用格局结构和功能对生态风险区变化的影响，安宁市生态风险防范和生态环境治理提出以下建议：①针对高和较高生态风险区，一方面要做好城镇国土空间规划，严格管控城镇开发边界、永久基本农田保护和生态保护红线，科学布局城镇空间、农业空间和生态空间，按照集约适度、绿色发展的要求有序扩张城镇用地，加强对城镇扩张而导致的边缘区域耕地、林地、灌木地的破碎化的整体性保护与修复。另一方面，协调好工业发展、城镇发展与生态环境之间的关系， 坚持安宁市工业总体规划中工业“点—轴—组团”的空间结构，集约利用工业用地；同时开展各工业区、城市区之间的生态隔离区、生态协调区和生态保护区建设，做好各分区景观生态系统综合整治和自然修复。②针对中等生态风险区，加强区域土地整理，增加林草覆被，减少景观破碎程度，加强生态系统稳定性。特别是农业用地的土地整理，协调耕地与生态景观和城镇景观之间的关系，构建“山水林田湖草”一体化生态保护修复新格局，加快“三产”融合，持续推进高原特色农业发展，力争使中等生态风险区转变为低、较低生态风险区。③针对低和较低生态风险区，避免在城镇开发建设中对林地、人造地表景观整体性的破坏，通过合理调整城镇总体规划，在保持现有良好自然环境的基础上，打造由生态背景、开敞空间体系和人文景观构成的多功能景观系统。以林地为主的低、较低生态风险区不仅应该严格管控林地的开发，改善林地的树种结构和林地质量，而且需严防森林火灾，防控由森林火灾引起的景观生态系统的扰动。④未来城镇发展情境选择一方面应考虑城镇发展阶段和未来城镇管控的重点布局，另一方面应综合分析多种情景模式下研究区生态安全格局优化，协调生态保护与经济发展之间的关系，才能在生态环境保护的基础上实现高质量的发展。不同情境对比下，基于生态保护情景的安宁市景观生态风险空间分布特征和变化趋势更接近生态安全格局优化方案，更符合安宁市转型升级样板区、滇中绿色发展示范区定位及绿色可持续发展理念。