基于PLUS-SD耦合模型的黑河流域中游甘临高地区土地利用研究

蒋小芳， 段翰晨， 廖 杰， 宋 翔， 薛 娴

（1.中国科学院西北生态环境资源研究院沙漠与沙漠化重点实验室，甘肃 兰州 730000；
2.中国科学院大学，北京 100049；
3.中国科学院西北生态环境资源研究院干旱区盐渍化研究站，甘肃 兰州 730000）

20 世纪90 年代土地利用结构变化问题开始受到关注［1］，现已成为LUCC 研究的热点。LUCC 模型是研究土地利用动态变化的重要手段，根据各种模型的优势可以分为数量模拟模型、空间模拟模型及耦合模型［2-5］。数量模拟模型主要包括SD（System Dynamics）、GF（Gray Forecast）和MC（Markov Chain）等模型，这类模型能有效地模拟和预测土地利用类型的数量结构，但难以实现空间结构模拟［6-8］。空间模拟模型有CA（Cellular Automaton）、CLUE-S（Conversion of Land Use and Its Effects at Small Regional Extent）、FLUS（Future Land Use Simulation Model）、ABM（Agent-Based Model）和PLUS（Patch-Generating Land Use Simulation Model）等模型，它们均可用于模拟地类的空间布局，但在数量模拟方面存在不足［9-11］。耦合模型旨在集合多个模型的优势，提高了单一模型的模拟精度，为土地利用研究提供有力保障［12-14］。

目前，用于土地利用/土地覆盖变化模拟的几种耦合模型主要是通过数量模型得到各地类的数量结构，然后用空间模型得到未来不同情景下的空间布局。He 等［15］以成都市为研究区，运用SD 模型在宏观层面上对各地类进行了数值上的模拟，采用CLUE-S 模型在微观层面上对该数值进行了空间上的分配，成功地探索了不同情景下的未来城市增长模式。Liang 等［16］基于SD-CLUE-S 和InVEST（Intergrated Valuation of Environmental Services and Tradeoffs）模型模拟和预测了2000—2018年张掖绿洲土地利用变化对像元尺度和区域尺度碳储量的影响。Ding 等［17］建立了一个基于FLUS 与InVEST 模型的耦合模型，有效地模拟了东营市2030年不同情景的人居环境质量变化。

研究表明，在耦合模型中，空间模拟模型的性能尤为重要，而不同模型的适用区域存在差异。PLUS模型通过随机森林算法挖掘地类转化概率，善于处理高维数据［18］。CLUE-S 模型来源于CLUE 模型，解决了CLUE 模型不适用于小区域模拟的局限性［4,19］。FLUS模型采用神经网络算法获得适宜性概率，有利于保证模拟精度的提升［20］。Peng等［21-22，4］将随机森林算法与CLUE-S 模型结合，充分利用随机森林算法在处理高维共线性问题方面的优势。林丽等［23］选择LCM（Land Change Modeler）、FLUS、CAMarkov和CLUE-S模型对地形复杂的云南山区进行研究，最终发现不同土地利用类型的适用模型存在差异，FLUS模型更适于模拟建设用地，CLUE-S模型在园地存在过拟合现象。Liang 等［18］对武汉市进行研究，发现PLUS 模型在该区域的模拟性能优于FLUS模型。由此可见，不同空间模拟模型在不同区域表现出不同的模拟结果。

位于中国西北干旱区的张掖市甘州区、临泽县和高台县（本文简称甘临高），处于河西走廊中部的黑河流域中游、青藏高原和蒙古高原的连接处，具有典型的区域代表性，加之该区生态环境和社会经济发展对土地利用变化较为敏感。因此，研究该区的土地利用变化具有重要意义。为深入探究PLUS、CLUE-S 和FLUS 空间模拟模型在该区的适用性，更好地预测土地利用变化，笔者选取甘临高地区为研究区，通过对比PLUS、FLUS 和CLUE-S 模型模拟结果，筛选最适于该区的空间模拟模型；
然后将其与SD 数量模拟模型耦合，预测2030 年不同情景下研究区的土地利用格局，为区域社会经济可持续发展提供参考。

1.1 研究区概况

黑河是中国第二大内陆河，黑河流域中游的张掖市甘州区、临泽县和高台县的是本文的研究区。该区地理位置为98°57′～100°52′E、38°39′～39°59′N（图1）。研究区的气候类型为温带大陆性气候，年均降水量约为100～150 mm，年均蒸发量高达2047.9 mm，年均气温在8～9 ℃区间波动。研究区地势南北高、中间低，主要的地表覆盖类型为耕地、草地、荒漠等［24］。该区位于干旱区，生态环境脆弱，草地的萎缩和耕地的不合理扩张加剧了荒漠化。

图1 研究区在甘肃省的地理位置（a）、高程（b）、地貌类型（c）和2015年土地利用空间分布（d）Fig.1 Geographical location of the study area in Gansu Province(a),elevation(b),landform type(c)and the spatial distribution map of land use in 2015(d)

1.2 数据来源及处理

土地利用数据、自然环境数据和社会经济数据均来源于中国科学院资源环境数据共享中心（http://www.resdc.cn）。除气温、降水、人口和GDP 数据的空间分辨率为1 km，其余数据的空间分辨率为30 m（表1）。本研究将所有数据的空间分辨率重采样为30 m。不同像元与道路、河流、居民点的距离数据均是ArcGIS 10.2 软件中基于相应矢量数据计算欧氏距离的结果。

表1 研究数据源信息Tab.1 Data source information of the research

1.3 模型与方法

1.3.1 土地利用空间格局模拟及验证

（1）3种空间模型介绍

PLUS 模型包括基于土地扩张分析策略的转化规则挖掘框架和基于多类型随机斑块种子机制的CA 模型两大模块［18］（图2）。CLUE-S 模型包括非空间土地需求模块和空间分配模块两大部分［4,19］。本研究采用线性插值法计算2000—2015 年间不同地类的需求面积，不同地类与驱动因素之间的关系数据来源于Logistics 回归方程。研究采用ROC（Relative Operating Characteristics）方法检验Logistics 回归的结果，当ROC＞0.7，说明回归方程能够较好地解释地类结构［25］，然后采用空间分配模块进行空间结构模拟。FLUS模型由2个模块构成［20］，第一个模块采用神经网络算法基于初始年份的土地利用数据与驱动因素栅格数据获得适宜性概率，第二个模块基于轮盘赌法则的自适应惯性竞争机制研究不同驱动因素影响下各地类的转化。

图2 本研究流程Fig.2 Flow chart of this study

（2）模型精度验证指标

主要采用3种模型精度验证方法。第一种方法是逐像元精度评价方法，评价指标为FoM［26-27］。第二种方法是空间格局比较法，评价指标为景观格局相似度，本研究采用Fragstats 4.2.1软件计算模拟与实际的土地利用结构图的景观指数，进而获得景观格局相似度Sim［28］，公式如下：

式中：n表示指标总数；
li,s表示模拟土地利用结构图的景观指数；
li,o表示实际土地利用结构图的景观指数。第三种方法是Kappa 系数，该系数主要被用于验证空间模拟精度。本文使用PONTIUS 发展的Kappa系数量化位置错误和数量错误［29］，采用IDRISI 17.0 软件计算数量Kappa 系数（Kno）、位置Kappa系数（Klocation）和标准Kappa系数（Kstandard）。

1.3.2 基于SD模型的未来情景预测

（1）SD模型介绍

SD 模型于1956 年建立，该模型是生态环境和经济发展复合系统的仿真模拟实验室，在数据模拟方面功能强大［10,30］。本研究使用的SD 模型包括人口、经济、生产力和气候4 个子系统（图3）。甘临高地区的各变量数据参考张掖地区，各系统变量的数值来源于《甘肃发展年鉴》《张掖统计年鉴》和《中国城市统计年鉴》等文献［31-34］。使用SPSS 23软件进行统计分析，获得不同变量间的数量关系，SD 模型模拟的时间范围为2000—2030年，时间步长为1 a。

图3 SD模型结构框架Fig.3 Flow chart of SD model

（2）SD模型检验

SD模型历史仿真检验的相对误差计算公式为：

为验证模型仿真效果，本文确定历史检验时间为2015年，历史检验的相对误差均低于5%，这说明模型的模拟精度较高，能够正确预测土地系统结构（表2）。

表2 SD模型模拟结果的历史性检验Tab.2 Historical test of SD model’s simulation results /hm2

（3）情景选择

“十二五”规划期间张掖市GDP 年均增速为10.5%，人口年均增长率为5.2‰，“十三五”规划期间GDP 年均增长率为8.5%，综合考虑将HD、ED 和EP 情景下的GDP 增长率分别设置为12%、14%和9%（表3）；
由于二孩政策的实施，本研究将HD、ED和EP情景下的人口增长率分别设置为6‰、10‰和4‰。甘肃省城镇体系规划（2013—2030年）中要求规划期末张掖绿洲生态城乡统筹发展区城镇化率应达到61.8%，据此本研究将HD、ED 和EP 情景下的城镇化率分别设置为61.8%、65%和60%。IPCC报告指出2016—2035 年全球气温和降水会逐年增加［35］。参考前人的研究［36］，HD、ED 和EP 情景下的年均气温与降水变化量都设置为正值。根据统计年鉴中的技术研发投入数据确定技术参数。

表3 不同情景的参数设定Tab.3 Parameter setting of different scenarios

2.1 PLUS模型与其他模型的模拟结果

为对比分析CLUE-S、PLUS和FLUS模型的模拟效果，通过重采样将研究像元尺度统一为50 m。当模拟尺度为50 m时，耕地、林地、农村建设用地和未利用地的ROC最高（表4）。由此可知，CLUE-S模型在本研究区的最佳研究尺度为50 m。PLUS 模型在30 m、50 m、100 m 和200 m 栅格尺度下的FoM 分别为0.234、0.184、0.145 和0.152，相同像元尺度下FLUS 模型的模拟精度分别为0.147、0.133、0.128 和0.117。据此可得，PLUS和FLUS模型的最佳拟合尺度均为30 m。

表4 CLUE-S模型的ROC指数Tab.4 ROC index of CLUE-S model

在50 m 研究尺度下PLUS、FLUS 和CLUE-S 模型的FoM 值分别为0.184、0.133 和0.133，PLUS 模型在模拟精度上优于FLUS 和CLUE-S 模型。PLUS、CLUE-S 和FLUS 模型的Kno分别为0.9625、0.9444 和0.9188，Klocation分别为0.9456、0.9193 和0.8818，Kstandard分别为0.9455、0.9191 和0.8818。3 个模型的Kappa系数值均大于0.8，且Kno高于Klocation，说明3个模型的数量预测能力均优于位置预测能力。PLUS 模型的3 种Kappa 系数最高，且均大于0.95，说明模型模拟结果与实际情况基本一致。

本文采用15 个景观指数比较3 个模型在挖掘景观变化方面的能力。Fragstats 4.2.1软件计算的结果表明（表5），PLUS 模型有9 个景观指数更贴近实际的景观指数：NP、PARA_MN、PARA_MD、PARA_RA、PARA_CV、ENN_MN、ENN_RA、ENN_MD和ENN_SD；
CLUE-S 模型有4个景观指数与实际的景观指数相符：PARA_AM、PLADJ、ENN_MD 和ENN_CV；
FLUS模型的LPI、PARA_SD和ENN_AM 3个景观指数值与实际的景观指数相近。PLUS、FLUS 和CLUE-S 模型的景观相似度Sim分别为0.64、0.21 和0.46，这说明PLUS 模型在景观格局模拟方面优势明显。

表5 不同模型2015年模拟结果的景观指数Tab.5 Landscape data of simulation results of different models in 2015

总体上，PLUS 模型的整体模拟效果明显优于CLUE-S和FLUS模型，CLUE-S模型在地势平坦地区模拟准确度高于FLUS模型，而FLUS模型在地势起伏较大地区的模拟效果优于CLUE-S 模型（图4）。（1）在临泽县、高台县和甘州区，CLUE-S 模型模拟的地形平坦处的耕地范围略大于实际范围，但是准确度较高。在高台县南部海拔高处，CLUE-S模型模拟的草地范围大于实际范围；
在临泽县和甘州区的北部海拔高处，CLUE-S模型模拟的草地范围小于实际范围；
在地形平坦的中部地区，CLUE-S 模型模拟的建设用地范围略小于实际范围，在海拔高处，该模型模拟的林地和建设用地范围大于实际范围。（2）在临泽县和甘州区，FLUS模型模拟的耕地范围显著大于实际范围，在高台县FLUS 模型模拟的耕地范围小于实际范围；
FLUS模型模拟的海拔高处的草地范围与实际情况较为接近，模拟精度较CLUE-S模型高，但略低于PLUS模型；
FLUS模型模拟的林地和建设用地范围较实际范围略大。（3）PLUS模型的模拟准确度显著高于CLUE-S和FLUS模型，其模拟结果基本与实际范围相符。与FLUS 模型相比，CLUE-S和PLUS模型模拟的农村建设用地范围与实际更相符。甘州区、临泽县和高台县的城市建设用地发展过程中均出现了“飞地式”扩张现象，3 个模型在模拟这类建设用地时均失效，这或许是驱动因素与地类变化之间时间滞后的结果。

图4 不同模型2015年的模拟结果对比Fig.4 Comparison of simulation results of different models in 2015

2.2 PLUS-SD耦合模型的不同情景预测结果

甘州区、临泽县和高台县各地类变化的区域多位于地类边缘，以边缘式扩张为主（图5，表6～表7）。（1）基于30 m分辨率土地利用数据模拟的结果表明，EP和HD情景中2015—2030年期间耕地面积增加，其来源均有草地。EP 情景下耕地面积减少，主要转化为林地和草地，说明政府采取退耕还林还草措施保护生态环境。（2）ED情景中城市化水平显著提高，技术发展迅速，2030年林地面积逐渐减少，主要转化为耕地。EP和HD情景林地面积增加，HD情景中部分未利用地转化为林地，EP情景增加的林地主要来自耕地。（3）2030年ED和HD情景的草地面积低于2015 年，EP 情景则相反。（4）城市发展导致热岛效应和雨岛效应，全球变暖促使降水增加，模拟结果显示3 种情景的水域面积均增加。（5）2030 年3 种情景的城市建设用地面积均增加，主要来自于未利用地。ED、EP和HD情景的农村建设用地均呈不同程度下降，可能是大量农村人口向城市人口转变。（6）3种情景的未利用地主要向耕地和城市建设用地转化，ED情景转化面积最大，HD情景次之，EP情景转化面积最小。

图5 不同情景中的2030年土地利用结构Fig.5 Land use structure in different scenarios in 2030

表6 不同情景的未来各土地利用类型面积预测结果Tab.6 Future prediction results of land use types under different scenarios /hm2

表7 不同情景的未来各土地利用转移矩阵Tab.7 Future land use transition matrix under different scenarios /hm2

3.1 PLUS-SD模型在甘临高地区的适用性分析

（1）通过PLUS、FLUS 和CLUE-S 模型在甘临高地区的模拟结果，可以发现PLUS模型在位置预测、数量预测和景观格局相似度方面均优于FLUS 和CLUE-S 模型，CLUE-S 模型在地形较为平坦地区模拟精度较高，而FLUS 模型在海拔高处模拟效果略优于CLUE-S 模型。前人的研究得出了类似结论，Liang 等［18］的研究表明，PLUS 模型在模拟历史土地利用变化过程方面优于FLUS 模型；
济南市区土地利用变化模拟结果显示FLUS模型模拟结果的Kappa 系数略低于CLUE-S 模型［37］。（2）本研究表明，FLUS 与PLUS 模型的研究格网愈密集，模拟结果愈优，而CLUE-S 模型并未随着格网密度增加而得到更优的模拟结果。张永民等［38］基于Logistic 回归模型采用100 m、250 m、500 m、750 m和1000 m这5种栅格尺度研究土地利用与不同影响因子的关系，结果显示500 m 尺度下二者的相关性最高。CLUE-S模型基于传统的二值逻辑斯蒂回归确定地类空间结构与影响因素之间的相关关系，其研究尺度不宜过度密集；
而PLUS 模型在模拟细尺度元胞变化方面能力较强；
FLUS模型通过增强空间分配和需求模拟之间的相互映射以达到“紧耦合的目标”。（3）SD模型适用于模拟土地结构中不同影响因素之间的数量关系，通过情景分析模拟政策变化的效果［39］。本研究将SD 与PLUS 模型耦合，获得了不同情景下的土地利用结构，有利于生态和经济的可持续发展。不少研究对未来不同情景下的土地利用格局进行预测［40-41］。Lu等［42］将SD和CCDM（Coupling Coordination Degree Model）结合，对武汉市经济-资源-环境系统的模拟结果进行了分析和评价，清晰地揭示了经济-资源-环境系统中各因素的耦合关系。Wu 等［43］建立SD-CLUE-S 耦合模型，发现综合性模型能够更好地模拟上海宝山地区景观生态系统价值的动态变化。

3.2 不同情景模拟结果的对比分析

甘临高地区的地表环境较恶劣，沙漠和戈壁面积较大。本研究的ED情景和EP情景偏重于经济发展或环境保护，不利于生态和经济协调发展;HD 情景兼具生态友好和经济建设的优点。前人研究表明，甘临高地区存在草地荒漠化现象，草地生态环境脆弱，耕地对水分和养分的需求量大，不合理地开垦耕地将加剧土地荒漠化进程［44］。2030年的HD情景中有少量草地和农村建设用地转化为耕地，耕地主要来自于未利用地；
退耕还林的加大，对未利用地的开发力度促使林地面积增加，有利于发挥国土空间的生态功能。HD情景中城市扩张的主要方向也是未利用地，随着技术的发展，未利用地的利用率和转化率显著提升。近年来高台等地的绿洲扩张明显，且绿洲波动带多位于绿洲和荒漠交界处［45］。随着城镇化进程的加快会导致建设用地扩展，2033年建设用地面积将高达4948.90×104hm2［46］。在HD情景中，政府需严格管控建设用地，防止建设用地过度占用耕地，保护永久性基本农田，要集约且合理地安排建设用地，形成科学有序的国土空间开发格局［47］。

本研究发现PLUS 模型在景观格局、位置和数量模拟方面的精度均优于CLUE-S 和FLUS 模型，CLUE-S模型在研究区地势低平地区模拟精度较高，而FLUS模型在海拔高处模拟效果较CLUE-S优；
对比分析PLUS-SD 耦合模型模拟的甘临高地区2030年ED、EP 和HD 情景下的土地利用情况，发现3 种情景在土地利用类型变化速率方面存在差异；
ED情景中土地利用类型变化速率最快，经济迅猛发展，但是环境破坏严重，林地大面积减少，EP 情景中土地利用类型变化速率偏慢，存在明显的退耕还林还草现象，但经济发展滞缓，HD情景中建设用地占用耕地的现象较少，经济平稳运行，有利于社会可持续发展。