1.本发明涉及国土空间规划和管理技术领域，具体地说，涉及一种基于生态保护重要性评价的城市群弹性开发边界划定方法及系统。
背景技术：2.城市蔓延是一种过度的城镇空间增长形式，容易导致城镇脱离合理的发展轨道而出现大量耕地侵蚀、土地利用率低和生态环境破坏等问题。城镇弹性开发边界是指一定时间阶段内允许城镇开发建设和禁止城镇开发建设区域之间的空间界线，该边界划定的本质是为了确定权衡生态保护和经济发展矛盾的城镇开发空间的范围，在守住自然生态安全边界的同时，引导城镇建设用地有序开发，从而有效治理城市蔓延等城镇化衍生问题，为可持续发展预留更多生态空间和农业空间，最终实现更大范围和更高水平的国土空间保护。城市群是城镇化过程发展到高级阶段出现的空间组织形式，也是承载区域发展要素的主要空间形式，科学划定城市群弹性开发边界，对于保障区域生态安全和引导国土空间协调发展具有重要意义。
3.在自然资源约束不断增强背景下，基于生态约束的反规划思维划定弹性开发边界已成为普遍共识。划定弹性开发边界的关键在于两点，一是生态约束的切入视角，二是边界划定过程中采用的土地利用变化模拟模型。一方面，既往研究不管是从生态系统服务等单一生态视角，抑或是从多维度所得综合生态格局等多重生态视角判定生态约束，其本质均是通过间接评价各地域的生态保护等级，来明确赋予到各地的生态约束强度，目前尚缺乏直接评价各地域的生态保护等级，特别是从生态系统的外部功能性和内在稳定性，即基于生态系统服务和生态脆弱性所得的生态保护重要性，判定各地域生态约束强度的技术实践。
4.另一方面，许多数学模型，如ca-markov模型、约束性ca模型、神经网络模型、geosos 模型和flus模型等，已被广泛应用于土地利用变化模拟进而划定弹性开发边界方面的研究。然而上述模型通常仅适用于模拟城市用地的空间变化，对自然用地的动态模拟以及多种土地利用斑块同步演化能力较弱。斑块生成土地利用变化模拟(patch-generating land usesimulation，plus)模型可用于挖掘土地利用转移规则并预测多种土地利用空间结构的斑块级演化，目前仍缺乏该模型在弹性开发边界划定研究中的相应探索。此外，国内外关于弹性开发边界的既往研究在研究尺度和研究对象上主要以中观的市、区县、自治州为主，城市群等宏观尺度的开发边界划定研究则相对较少。
技术实现要素：5.本发明的目的在于解决现有技术忽略在城市群尺度上基于生态保护重要性评价和plus 模型进行弹性开发边界划定的问题，提供一种城市群弹性开发边界划定方法，步骤科学、通用且可复制，适合应用于任何地区特定时间框架下城市群弹性开发边界划定。
6.为了实现上述目的，本发明所采用的技术方案如下：
7.一种城市群弹性开发边界划定方法，包括以下步骤：
8.s1：选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系，开展生态保护重要性评价指标的具体评价，建立生态保护重要性评价指标空间数据库；
9.s2：根据生态保护重要性评价指标空间数据库进行主成分分析，生成主成分指标数据和主成分分析报告，进行kmo和bartlett's检验以及主成分方差贡献率累积值检验，确定主成分指标数据及各主成分指标的准则权重；
10.s3：设置多种决策风险系数情景，根据确定的主成分指标数据，计算主成分指标的位序权重和权衡度，运用owa算子聚合不同情景下的主成分指标的位序权重和准则权重，进行多情景的生态保护重要性评价模拟；
11.s4：根据生态保护重要性评价模拟结果，提取各种情景下的生态保护优先区，对比不同情景的生态保护优先区对于各情景的生态保护重要性评价结果的保护效率和权衡度，筛选出保护效率与权衡度共同最优情景的生态保护优先区；
12.s5：选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系，开展土地利用变化驱动因子的具体评价，建立土地利用变化驱动因子空间数据库；
13.s6：根据历史土地利用数据，获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵，根据历史土地利用数据，提取各土地利用类型的扩张斑块数据，计算生态保护情景的土地利用转移概率矩阵，获得生态保护情景的目标年土地利用数量结构；
14.s7：根据土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，输出各土地利用类型的发展概率数据，根据历史土地利用数据和各土地利用类型的发展概率数据，模拟基期年的土地利用空间结构数据，并进行模拟精度验证；
15.s8：基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟自然发展情景目标年土地利用空间结构数据；基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟生态保护情景的目标年土地利用空间结构数据，并划定弹性开发边界。
16.进一步地，步骤s1具体为：选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系；结合多源数据、arcgis软件和invest模型，开展生态保护重要性评价指标的具体评价；评价结果统一进行归一化处理，以消除不同评价指标之间的量纲影响；借助arcgis平台整合所有归一化后的评价指标，建立研究区生态保护重要性评价指标空间数据库；
17.其中，从属于生态系统服务的供给服务、调节服务、支持服务以及从属于生态脆弱性的暴露度、敏感性和适应力为一级指标层，选取的具体指标为生态保护重要性评价指标体系的二级指标层，根据对目标层的影响，所选具体指标中从属于适应力的指标为负向指标，其他指标均为正向指标；
18.正向指标进行归一化的计算方法为：
[0019][0020]
负向指标进行归一化的计算方法为：
[0021][0022]
式中，x
norm
为归一化后的指标像元值，x
x
为指标像元值，x
max
为指标所有像元值的最大值，x
min
为指标所有像元值的最小值。
[0023]
进一步地，步骤s2具体为：使用qgis软件的对栅格取样工具，将所有评价指标的栅格图层取样至由研究所需分辨率的栅格转换而成的点图层，然后将点图层的属性表数据导入 spss软件进行主成分分析，生成主成分指标数据和主成分分析报告；进行kmo和bartlett's 检验，将特征值大于1作为选取主成分的标准，进行主成分方差贡献率累积值检验；确定基于生态保护重要性评价指标的主成分指标数据，并根据主成分方差贡献率确定和为1的各主成分指标的准则权重；
[0024]
其中，进行kmo和bartlett's检验，具体为：
[0025]
若主成分分析报告中的kmo和bartlett's检验结果显示：kmo取样适切性量数大于0.7，并且巴特利特球形度检验的显著度小于0.001，则研究数据适合进行主成分分析；否则，研究数据不适合进行主成分分析；
[0026]
进行主成分方差贡献率累积值检验，具体为：
[0027]
若主成分分析报告中的总方差解释结果显示：特征值大于1的主成分方差贡献率累积值超过60％，则生成的主成分指标数据较为合理且有效；否则，生成的主成分指标数据不合理且无效。
[0028]
进一步地，步骤s3具体为：设置z(z≥5)种决策风险系数情景，根据确定的主成分指标数据，使用单调规则递增方法计算主成分指标的位序权重，具体计算方法如下：
[0029][0030][0031]
式中，j为位序；vj为位序权重，vj∈[0,1]；n为主成分指标数量；α为决策风险系数，α∈(0, ∞)；wk为主成分指标重要等级；rk为主成分指标的赋值，按照指标的像元值叠加准则权重后的大小对指标的进行赋值，最大值赋1，次大值赋2，最小值赋n；
[0032]
权衡度为不同决策风险系数下主成分指标间的补偿程度，具体计算方法如下：
[0033][0034]
其中，tradeoff为权衡度，且0≤tradeoff≤1；n为主成分指标数量；wk为第k个主成分指标重要等级；
[0035]
运用terrset软件的mce模块的owa算子机制，聚合不同情景下的主成分指标的位序权重和准则权重，从而获得z种情景下的生态保护重要性空间评价结果。
[0036]
进一步地，步骤s4具体为：基于z种情景下的生态保护重要性评价模拟结果，进行降序排列下的累积生态保护重要性像元值与总的生态保护重要性像元值的比值计算，将比值为 50％对应的生态保护重要性像元值作为各情景生态保护重要性阈值；
[0037]
在arcgis软件中使用聚合工具，将大于或等于生态保护重要性阈值的区域中相对聚集或邻近的图斑聚合为相对完整连片图斑，聚合距离为研究所需分辨率，聚合结果作为各情景初步生态保护优先区；
[0038]
为降低生态保护优先区的破碎化程度，利用情景z1的生态保护优先区图斑确定需要剔除的图斑面积阈值；其中，情景z1为z种情景中各主成分指标位序权重相等的情景，统计分析情景z1的生态保护优先区中面积小于阈值的图斑与所有图斑的数量比值，以及其总面积随面积阈值变化的情况，将拐点对应的最小面积阈值作为最终的图斑面积阈值，并根据阈值剔除相应独立分散小图斑，在此基础上合并国家级和省级的自然保护区、森林公园、风景名胜区等现有各类生态保护区，最终获得各情景最终的生态保护优先区；
[0039]
对比不同情景最终的生态保护优先区对于情景z1的生态保护重要性评价结果的保护效率和权衡度，从而筛选出保护效率与权衡度共同最优情景的最终生态保护优先区；
[0040]
其中，保护效率的计算方法如下：
[0041][0042]
式中，z为情景序号，z＝1,2，...z，pz为情景z的保护效率，为情景z的生态保护优先区范围内情景z1的生态保护重要性平均值，为整个研究区范围内情景z1的生态保护重要性平均值。
[0043]
进一步地，步骤s5具体为：选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系；结合多源数据和arcgis软件，开展土地利用变化驱动因子的具体评价；所得评价结果统一进行归一化处理，以消除不同驱动因子之间的量纲影响；借助arcgis平台整合所有归一化后的驱动因子，建立研究区土地利用变化驱动因子空间数据库；
[0044]
其中，从属于自然因素的地形和水域以及从属于人文因素的人口、经济、区位、基础设施为一级指标层，基础设施包括交通、医疗、教育和商服，选取的具体指标为二级指标层的土地利用变化驱动因子体系。
[0045]
进一步地，步骤s6具体为：基于两期历史土地利用数据，运用plus模型的markov chain 模块获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵；
[0046]
基于两期历史土地利用数据，运用plus模型的提取土地利用扩张模块提取各土地利用类型的扩张斑块数据，以此为基础计算城镇建设用地扩张斑块面积中各非建设用地的占比，并将所得占比数据作为非建设用地向城镇建设用地转移概率的负增长率，从而基于自然发展情景的土地利用转移概率矩阵计算生态保护情景的土地利用转移概率矩阵；
[0047]
基于markov chain原理和生态保护情景的土地利用转移概率矩阵计算生态保护情景的目标年土地利用数量结构；
[0048]
其中，历史土地利用数据包括耕地、林地、草地、水域、城镇建设用地、农村居民点、未利用地共七种土地利用类型，其中城镇建设用地和农村居民点属于建设用地，其余土地利用类型属于非建设用地。
[0049]
进一步地，步骤s7具体为：基于土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，运用plus模型的leas模块输出各土地利用类型的发展概率数据；计
算和设置plus模型的cars模块所需参数；以第二期历史土地利用数据的土地利用数量结构为土地利用需求数据，基于第一期历史土地利用数据和各土地利用类型的发展概率数据，运用 plus模型的cars模块模拟基期年的土地利用空间结构；基于两期历史土利用数据和基期年土地利用空间结构数据，运用plus模型confusion matrix and fom模块进行模拟精度验证；
[0050]
其中，计算和设置plus模型的cars模块所需参数，具体为：
[0051]
基于各土地利用类型的扩张斑块数据计算各土地利用类型的扩张斑块面积在总扩张面积中的占比，将所得占比数据作为邻域权重参数；将土地转移矩阵中各土地利用类型的是否允许转移值设为1，不限制任何土地利用类型之间的转移；递减阈值的衰减系数、随机斑块种子的概率、随机斑块种子最大比例、邻域范围、并行线程数量等参数采用默认值；
[0052]
进行模拟精度验证，具体为：
[0053]
若kappa系数大于0.8且fom系数大于0.2，则plus模型的模拟结果和真实结果基本一致，模拟结果可信；否则，需不断调整plus模型的输入数据和相关参数，直至满足上述kappa 系数和fom系数检验要求才能继续进行后续步骤。
[0054]
进一步地，步骤s8具体为：以自然发展情景目标年土地利用数量结构为土地利用需求数据，以及现有的各类生态保护区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟自然发展情景目标年土地利用空间结构数据；
[0055]
以生态保护情景目标年土地利用数量结构为土地利用需求数据，以及生态保护优先区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟生态保护情景目标年土地利用空间结构数据；
[0056]
将自然发展情景和生态保护情景目标年土地利用空间结构数据中城镇建设用地边界，分别划定为研究区自然发展情景和生态保护情景目标年的弹性开发边界。
[0057]
一种城市群弹性开发边界划定系统，包括：
[0058]
评价指标模块，用于建立生态保护重要性评价指标空间数据库；
[0059]
主成分分析模块，用于根据生态保护重要性评价指标空间数据库进行主成分分析；
[0060]
owa算子模块，用于基于owa算子进行多情景的生态保护重要性评价模拟；
[0061]
情景筛选模块，用于提取各情景的生态保护优先区并筛选出最优情景的生态保护优先区；
[0062]
驱动因子模块，用于建立土地利用变化驱动因子空间数据库；
[0063]
需求预测模块，用于预测自然发展情景和生态保护情景的目标年土地利用数量结构；
[0064]
模拟与验证模块，用于模拟基期年的土地利用空间结构并进行模拟精度验证；
[0065]
边界划定模块，用于模拟自然发展情景和生态保护情景的目标年土地利用空间结构并划定弹性开发边界。
[0066]
与现有技术相比，本发明具有以下优点和有益效果：
[0067]
(1)本发明通过耦合生态系统服务指标和生态脆弱性指标的生态保护重要性评
价，直接判定各地域生态约束强度和生态保护等级，然后以生态保护重要性评价结果提取生态保护优先区，接着运用由生态保护优先区所得限制发展区进行土地利用模拟进而划定开发边界，有效避免了现有技术忽略基于生态保护重要性评价进行弹性开发边界划定的问题。基于本发明技术划定的开发边界有助于保护优质生态空间、维持区域生态安全和保障人类福祉。同时创新提出基于保护效率和权衡度从多情景中筛选出最优情景的生态保护优先区的方法，为自然资源管理和生态保护的决策提供了有价值的技术指导。
[0068]
(2)本发明运用plus模型进行土地利用模拟进而划定城市群开发边界，提供了将plus 模型用于城市群开发边界划定的机制和方法，有效避免了现有技术忽略基于plus模型进行城市群弹性开发边界划定的问题，有助于弥补现有开发边界划定技术的不足，以及提高开发边界划定的科学性。同时本发明适用于城市群尺度的开发边界划定，弥补了现有技术在城市群等宏观尺度的开发边界划定技术研究方面的不足。
[0069]
(3)本发明提供了自然发展情景和生态保护情景的开发边界划定方法，满足了不同情景的边界划定方法需求，以便不同情景的边界划定结果对比分析，从而提高边界划定的科学性，有助于构建更加科学合理的国土空间开发保护格局。同时本发明步骤科学、通用且可复制，过程涉及的数据处理与空间分析，通过目前被生态、地理、规划等多领域广泛使用的arcgis 等软件即可实现，适用于任何地区特定时间框架下研究区城市群开发边界划定。
附图说明
[0070]
图1为实施例中城市群弹性开发边界划定方法的流程示意图。
[0071]
图2为实施例中生态保护重要性评价指标空间分布地图。
[0072]
图3为实施例中9种情景下的生态保护重要性空间评价分布地图。
[0073]
图4为实施例中情景5的生态保护优先区中面积小于阈值的图斑与所有图斑的数量比值以及其总面积随面积阈值变化统计分析图。
[0074]
图5为实施例中各情景最终的生态保护优先区空间分布地图。
[0075]
图6为实施例中土地利用变化驱动因子空间分布地图。
[0076]
图7为实施例中各土地利用类型的发展概率空间分布地图。
[0077]
图8为实施例中自然发展情景2050年土地利用空间结构空间分布地图。
[0078]
图9为实施例中生态保护情景2050年土地利用空间结构空间分布地图。
[0079]
图10为实施例中自然发展情景2050年弹性开发边界空间分布地图。
[0080]
图11为实施例中生态保护情景2050年弹性开发边界空间分布地图。
具体实施方式
[0081]
下面结合附图和具体实施例对本发明城市群弹性开发边界划定方法及系统作进一步说明。
[0082]
请参阅图1，本发明公开了一种城市群弹性开发边界划定方法，包括以下步骤：
[0083]
s1：选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系，开展生态保护重要性评价指标的具体评价，建立生态保护重要性评价指标空间数据库。
[0084]
s2：根据生态保护重要性评价指标空间数据库进行主成分分析，生成主成分指标数据和主成分分析报告，进行kmo和bartlett's检验以及主成分方差贡献率累积值检验，确定主成分指标数据及各主成分指标的准则权重。
[0085]
s3：设置多种决策风险系数情景，根据确定的主成分指标数据，计算主成分指标的位序权重和权衡度，运用owa算子聚合不同情景下的主成分指标的位序权重和准则权重，进行多情景的生态保护重要性评价模拟。
[0086]
s4：根据生态保护重要性评价模拟结果，提取各种情景下的生态保护优先区，对比不同情景的生态保护优先区对于各情景的生态保护重要性评价结果的保护效率和权衡度，筛选出保护效率与权衡度共同最优情景的生态保护优先区。
[0087]
s5：选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系，开展土地利用变化驱动因子的具体评价，建立土地利用变化驱动因子空间数据库。
[0088]
s6：根据历史土地利用数据，获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵，根据历史土地利用数据，提取各土地利用类型的扩张斑块数据，计算生态保护情景的土地利用转移概率矩阵，获得生态保护情景的目标年土地利用数量结构。
[0089]
s7：根据土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，输出各土地利用类型的发展概率数据，根据历史土地利用数据和各土地利用类型的发展概率数据，模拟基期年的土地利用空间结构数据，并进行模拟精度验证。
[0090]
s8：基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟自然发展情景目标年土地利用空间结构数据；基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟生态保护情景的目标年土地利用空间结构数据，并划定弹性开发边界。
[0091]
步骤s1中，选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层，生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系。从属于生态系统服务的供给服务、调节服务、支持服务以及从属于生态脆弱性的暴露度、敏感性和适应力为一级指标层，选取的具体指标为生态保护重要性评价指标体系的二级指标层。根据对目标层的影响，将所选具体指标中从属于适应力的指标分为负向指标，其余指标分为正向指标。然后结合多源数据、arcgis软件和invest模型，开展生态保护重要性评价指标的具体评价。所得评价结果统一进行归一化处理，以消除不同评价指标之间的量纲影响。借助arcgis平台整合所有归一化后的评价指标，建立研究区生态保护重要性评价指标空间数据库。
[0092]
对于正向指标和负向指标分别采用公式(1)和公式(2)进行极差法归一化，其中，正向指标进行极差法归一化的计算方法为：
[0093][0094]
负向指标进行极差法归一化的计算方法为：
[0095][0096]
式中，x
norm
为归一化后的指标像元值，x
x
为指标像元值，x
max
为指标所有像元值的最大值，x
min
为指标所有像元值的最小值。
[0097]
步骤s2中，基于生态保护重要性评价指标空间数据库，使用qgis软件的对栅格取
对应的最小面积阈值作为最终的图斑面积阈值，并根据阈值剔除相应独立分散小图斑，在此基础上合并国家级和省级的自然保护区、森林公园、风景名胜区等现有各类生态保护区，最终获得各情景最终的生态保护优先区。对比不同情景最终的生态保护优先区对于情景z1的生态保护重要性评价结果的保护效率和权衡度，从而筛选出保护效率与权衡度共同最优情景(保护效率越高，权衡度越低，相应情景越优)的最终生态保护优先区；
[0110]
保护效率的具体计算方法如公式(6)所示：
[0111][0112]
式中，z为情景序号，z＝1,2，...z，pz为情景z的保护效率，为情景z的生态保护优先区范围内情景z1的生态保护重要性平均值，为整个研究区范围内情景z1的生态保护重要性平均值。
[0113]
步骤s5中，选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系。从属于自然因素的地形和水域以及从属于人文因素的人口、经济、区位、基础设施(包括交通、医疗、教育和商服)为一级指标层，选取的具体指标为土地利用变化驱动因子体系的二级指标层。然后结合多源数据和arcgis软件，开展土地利用变化驱动因子的具体评价。所得评价结果统一运用公式(1)进行归一化处理，以消除不同驱动因子之间的量纲影响。借助arcgis平台整合所有归一化后的驱动因子，建立研究区土地利用变化驱动因子空间数据库。
[0114]
步骤s6中，基于两期历史土地利用数据(包括耕地、林地、草地、水域、城镇建设用地、农村居民点、未利用地共七种土地利用类型，其中城镇建设用地和农村居民点属于建设用地，其余土地利用类型属于非建设用地)，运用plus模型的markov chain模块，获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵。
[0115]
接着基于两期历史土地利用数据，运用plus模型的提取土地利用扩张模块，提取各土地利用类型的扩张斑块数据，以此为基础计算城镇建设用地扩张斑块面积中各非建设用地的占比，并将所得占比数据作为非建设用地向城镇建设用地转移概率的负增长率，从而基于自然发展情景的土地利用转移概率矩阵，计算生态保护情景的土地利用转移概率矩阵。最后，基于markov chain原理和生态保护情景的土地利用转移概率矩阵，计算生态保护情景的目标年土地利用数量结构。
[0116]
步骤s7中，基于土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，运用plus模型的leas模块，输出各土地利用类型的发展概率数据。然后计算和设置plus 模型的cars模块所需参数。接着以第二期历史土地利用数据的土地利用数量结构为土地利用需求数据，基于第一期历史土地利用数据和各土地利用类型的发展概率数据，运用plus 模型的cars模块模拟第二期历史土地利用数据所属年份，即基期年的土地利用空间结构。最后基于两期历史土地利用数据和模拟所得的基期年土地利用空间结构数据，运用plus模型confusion matrix and fom模块进行模拟精度验证。
[0117]
基于各土地利用类型的扩张斑块数据，计算各土地利用类型的扩张斑块面积在总扩张面积中的占比，将所得占比数据作为邻域权重参数。接着，将土地转移矩阵中各土地利用类型的是否允许转移值设为1，不限制任何土地利用类型之间的转移。递减阈值的衰减系
数、随机斑块种子的概率、随机斑块种子最大比例、邻域范围、并行线程数量等参数采用默认值。
[0118]
若kappa系数大于0.8且fom系数大于0.2，则plus模型的模拟结果和真实结果基本一致，模拟结果可信；否则，需不断调整plus模型的输入数据和相关参数，直至满足上述kappa 系数和fom系数检验要求才能继续进行后续步骤。
[0119]
步骤s8中，以自然发展情景目标年土地利用数量结构为土地利用需求数据，以及现有的各类生态保护区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块，模拟自然发展情景目标年土地利用空间结构数据。
[0120]
接着以生态保护情景目标年土地利用数量结构为土地利用需求数据，以及生态保护优先区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus 模型的cars模块，模拟生态保护情景目标年土地利用空间结构数据。
[0121]
最后将自然发展情景和生态保护情景目标年土地利用空间结构数据中城镇建设用地边界分别划定为研究区自然发展情景和生态保护情景目标年的弹性开发边界。
[0122]
由上可知，本发明提供的一种城市群弹性开发边界划定方法，首先选取指标构建生态保护重要性评价指标体系，并开展生态保护重要性评价指标具体评价以建立生态保护重要性评价指标空间数据库。然后根据评价指标空间数据库进行主成分分析，并进行kmo和bartlett's 检验以及主成分方差贡献率累积值检验，在此基础上确定主成分指标数据及其准则权重。设置z(z≥5)种决策风险系数情景，并计算主成分指标的位序权重和权衡度，在此基础上基于 owa算子进行多情景的生态保护重要性评价模拟。根据生态保护重要性评价模拟结果提取各种情景下的生态保护优先区，并计算各情景的生态保护优先区的保护效率，从而筛选出保护效率与权衡度共同最优情景的生态保护优先区。
[0123]
接着选取指标构建土地利用变化驱动因子体系，并建立土地利用变化驱动因子空间数据库。然后运用plus模型预测自然发展情景和生态保护情景的目标年土地利用数量结构，以及运用plus模型进行基期年的土地利用空间结构模拟，并进行模拟精度验证。最后运用 plus模型进行自然发展情景和生态保护情景的目标年土地利用空间结构模拟，将两种情景中城镇建设用地边界作为研究区目标年的弹性开发边界。本发明的城市群弹性开发边界划定方法步骤科学、通用且可复制，适合应用于任何地区特定时间框架下城市群弹性开发边界划定，为国土空间规划的开发边界决策和空间管控提供理论依据和技术支持。
[0124]
下面选择粤港澳大湾区为研究区域，对本发明的城市群弹性开发边界划定方法进行详细说明。
[0125]
步骤一、建立生态保护重要性评价指标空间数据库。
[0126]
基于研究区域生态环境特点以及相关的文献归纳和政策解读，选取17项指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系。从属于生态系统服务的供给服务、调节服务、支持服务以及从属于生态脆弱性的暴露度、敏感性和适应力为一级指标层，选取的具体指标为生态保护重要性评价指标体系的二级指标层。根据对目标层的影响，将所选具体指标中从属于适应力的指标分为负向指标，其余指标分为正向指标。生态保护重要性评价指标体系如表1所示。
[0127]
表1生态保护重要性评价指标体系
[0128][0129]
然后结合多源数据(相关数据类型及来源详见表2所示)、arcgis软件和invest模型，开展生态保护重要性评价指标的具体评价。
[0130]
表2数据说明
[0131]
[0132][0133]
其中，各评价指标具体评价方法如表3所示。
[0134]
表3粤港澳大湾区生态保护重要性评价指标的评价方法说明
[0135]
[0136][0137]
以上所有数据均通过arcgis软件进行重采样或空间插值到100m空间分辨率。
[0138]
所得评价结果统一进行归一化处理，以消除不同评价指标之间的量纲影响。借助arcgis 平台整合所有归一化后的评价指标，建立研究区生态保护重要性评价指标空间数据库。研究区生态保护重要性评价指标空间分布地图如图2所示。
[0139]
正向指标和负向指标分别采用公式(1)和公式(2)进行极差法归一化，其中，正向指标进行极差法归一化的计算方法为：
[0140][0141]
负向指标进行极差法归一化的计算方法为：
[0142][0143]
其中，x
norm
为归一化后的指标像元值，x
x
为指标像元值，x
max
为指标所有像元值的最大值，x
min
为指标所有像元值的最小值。
[0144]
步骤二、根据评价指标空间数据库进行主成分分析。
[0145]
基于生态保护重要性评价指标空间数据库，使用qgis软件的对栅格取样工具，将所有评价指标栅格图层取样至由100m分辨率栅格转换而成的点图层，然后将点图层的属性表数据导入spss软件进行主成分分析，生成主成分指标数据和主成分分析报告。然后进行kmo 和bartlett's检验，然后将特征值大于1作为选取主成分的标准，进行主成分方差贡献率累积值检验。最终确定了基于生态保护重要性评价指标的5个主成分指标数据，并根据主成分方差贡献率确定和为1的各主成分指标的准则权重，详见表4。
[0146]
表4各主成分指标特征值、方差贡献率和准则权重
[0147][0148]
本实施例的主成分分析报告显示：kmo取样适切性量数为0.813(＞0.7)，巴特利特球形度检验的显著度为0.000(＜0.001)，因此本实施例数据通过kmo检验和bartlett's检验，适宜进行主成分提取。同时特征值大于1的5个主成分指标的方差贡献率累积值超过60％，因此本实施例数据通过主成分方差贡献率累积值检验，生成的主成分指标数据较为合理且有效。
[0149]
步骤三、基于owa算子进行多情景的生态保护重要性评价模拟。
[0150]
设置9种决策风险系数情景(α＝0.00001，0.1，0.2，0.5，1，2，5，10，100000)，然后使用单调规则递增方法计算基于5个主成分指标的位序权重，具体计算方法如公式3、4所示：
[0151][0152][0153]
式中，j为位序；vj为位序权重，vj∈[0,1]；n为主成分指标数量；α为决策风险系数，α∈(0,∞)；wk为主成分指标重要等级；rk为主成分指标的赋值，按照指标的像元值叠加准则权重后的大小对指标的进行赋值，最大值赋1，次大值赋2，最小值赋n。
[0154]
权衡度(trade-off)表示为不同决策风险系数下主成分指标间的补偿程度，具体计算方法如公式5所示：
[0155][0156]
其中，tradeoff为权衡度，且0≤tradeoff≤1；n为主成分指标数量；wk为第k个主成分指标重要等级。9种情景下的决策风险系数、权衡度和位序权重如表5所示。
[0157]
表5不同情景下的决策风险系数、权衡度和位序权重
[0158][0159]
运用terrset软件的mce模块的owa算子机制，聚合不同情景下的主成分指标的位序权重和准则权重，从而获得9种情景下的生态保护重要性空间评价分布地图，详见图3。
[0160]
步骤四、提取各情景的生态保护优先区并择优。
[0161]
基于9种情景下的生态保护重要性评价模拟结果，进行降序排列下的累积生态保护重要性像元值与总的生态保护重要性像元值的比值计算，将比值为50％对应的生态保护重要性像元值作为各情景生态保护重要性阈值。接着在arcgis软件中使用聚合工具，将大于或等于生态保护重要性阈值的区域中相对聚集或邻近的图斑聚合为相对完整连片图斑，聚合距离为100 m，聚合结果作为各情景初步生态保护优先区。
[0162]
为降低生态保护优先区的破碎化程度，利用情景5的生态保护优先区图斑确定需要剔除的图斑面积阈值。其中，情景5为9种情景中各指标位序权重相等的情景。以2km2为初始值和2km2为步长，统计分析情景5的生态保护优先区中面积小于阈值的图斑与所有图斑的数量比值，以及其总面积随面积阈值变化的情况，详见图4。本实施例将图4的拐点(最先出现增长率为0的转折点)对应的最小面积阈值26km2作为最终的图斑面积阈值将拐点，并根据阈值剔除相应独立分散小图斑，在此基础上合并国家级和省级的自然保护区、森林公园、风景名胜区等现有各类生态保护区，最终获得各情景最终的生态保护优先区空间分布地图，详见图5。对比不同情景最终的生态保护优先区对于情景5的生态保护重要性评价结果的保护效率和权衡度，从而筛选出保护效率与权衡度共同最优情景(保护效率越高，权衡度越低，相应情景越优)的最终生态保护优先区。保护效率的具体计算方法如公式(6)所示：
[0163][0164]
式中，z为情景序号，z＝1,2，...z，pz为情景z的保护效率，为情景z的生态保护优先区范围内情景z1的生态保护重要性平均值，为整个研究区范围内情景z1的生态保护重要性平均值，约为0.4069。各情景生态保护优先区范围内情景5的生态保护重要性平均值和保护效率详见表6。
[0165]
表6不同情景的生态保护优先区范围内情景5的生态保护重要性平均值以及保护
效率
[0166][0167]
基于上述分析，从保护效率来看，情景4和情景5的保护效率最高，同为1.0553。再结合权衡度，情景4比情景5的权衡度低，即情景4相比于情景5在一定程度上弥补了主成分指标间的权衡影响。因此，本实施例最终选择情景4的生态保护优先区作为最优情景的态保护优先区。
[0168]
步骤五、建立土地利用变化驱动因子空间数据库。
[0169]
基于研究区的自然本底条件和社会经济条件，选取19个驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系，从属于自然因素的地形和水域以及从属于人文因素的人口、经济、区位、基础设施(包括交通、医疗、教育和商服)为一级指标层，选取的具体指标为土地利用变化驱动因子体系的二级指标层。粤港澳大湾区土地利用变化驱动因子体系如表7所示。
[0170]
表7粤港澳大湾区土地利用变化驱动因子体系
[0171][0172]
然后结合多源数据(相关数据类型及来源详见表2)和arcgis软件，开展土地利用变化驱动因子的具体评价，各驱动因子具体评价方法如表8所示。
[0173]
表8粤港澳大湾区土地利用变化驱动因子的评价方法说明
[0174]
[0175][0176]
评价结果统一运用公式进行归一化处理，以消除不同驱动因子之间的量纲影响。借助 arcgis平台整合所有归一化后的驱动因子，建立研究区土地利用变化驱动因子空间数据库。研究区土地利用变化驱动因子空间分布地图如图6所示。
[0177]
步骤六、预测自然发展情景和生态保护情景的目标年土地利用数量结构。
[0178]
基于2010年和2020年历史土地利用数据(包括耕地、林地、草地、水域、城镇建设用地、农村居民点、未利用地共七种土地利用类型，其中城镇建设用地和农村居民点属于建设用地，其余土地利用类型属于非建设用地)，运用plus模型的markov chain模块，获得自然发展情景的2050年土地利用数量结构(如表12中所示)和土地利用转移概率矩阵(如表 9所示)。
[0179]
表9自然发展情景的土地利用转移概率矩阵
[0180]
[0181][0182]
接着基于2010年和2020年历史土地利用空间结构数据，运用plus模型的提取土地利用扩张模块提取各土地利用类型的扩张斑块数据，以此为基础计算城镇建设用地扩张斑块面积中各非建设用地的占比，如表10所示。并将所得占比数据作为非建设用地向城镇建设用地转移概率的负增长率，即耕地、林地、草地、水域、未利用地向城镇建设用地的转移概率分别降低45.64％、27.87％、2.21％、20.19％、0.19％。
[0183]
表10 2010—2020年城镇建设用地扩张面积中各非建设用地类型的占比
[0184][0185]
基于自然发展情景的土地利用转移概率矩阵，计算生态保护情景的土地利用转移概率矩阵，如表11所示。
[0186]
表11生态保护情景的土地利用转移概率矩阵
[0187][0188]
最后基于markov chain原理和生态保护情景的土地利用转移概率矩阵，计算生态保护情景的2050年土地利用数量结构。如表12所示。
[0189]
表12基于markov chain预测的2050土地利用数量结构
[0190]
[0191][0192]
步骤七、模拟基期年的土地利用空间结构并进行模拟精度验证。
[0193]
基于土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，运用plus 模型的leas模块输出各土地利用类型的发展概率空间分布地图详见图7。
[0194]
然后计算和设置plus模型的cars模块所需参数。基于各土地利用类型的扩张斑块数据计算各土地利用类型的扩张斑块面积在总扩张面积中的占比，将所得占比数据作为邻域权重参数(如表13)。接着将土地转移矩阵中各土地利用类型的是否允许转移值设为1，不限制任何土地利用类型之间的转移。递减阈值的衰减系数(0.9)、随机斑块种子的概率(0.1)、随机斑块种子最大比例(0.001)、邻域范围(3)、并行线程数量(1)等参数采用默认值。
[0195]
表13各土地利用类型的邻域权重参数
[0196][0197]
接着以2020年历史土地利用数量结构为土地利用需求数据，基于2010年历史土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟2020年的土地利用空间结构数据。最后基于2010年和2020年历史土地利用地数据和模拟所得的2020年土地利用空间结构数据，运用plus模型confusion matrix and fom模块进行模拟精度验证。 confusion matrix and fom模块输出结果显示：kappa系数为0.9389(》0.8)，因此plus模型的模拟结果和真实结果基本一致，模拟结果可信。
[0198]
步骤八、模拟自然发展情景和生态保护情景的目标年土地利用空间结构并划定弹性开发边界。
[0199]
以自然发展情景2050土地利用数量结构为土地利用需求数据，以及现有的各类生态保护区为限制发展区，基于2020年土地利用空间结构数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟得到自然发展情景2050年土地利用空间结构空间分布地图，详见图8。接着以生态保护情景2050年土地利用数量结构为土地利用需求数据，以及生态保护优先区为限制发展区，基于2020年土地利用空间结构数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟生态保护情景2050年土地利用空间结构空间分布地图，详见图9。最后将自然发展情景和生态保护情景2050年土地利用空间结构空间分布地图中城镇建设用地边界，分别划定为研究区自然发展情景和生态保护情景的2050年的弹性开发边界，详见图10和图11。
[0200]
综上，本发明的城市群弹性开发边界划定方法具有如下优点和有益效果：
[0201]
(1)本发明通过耦合生态系统服务指标和生态脆弱性指标的生态保护重要性评价，直接判定各地域生态约束强度和生态保护等级，然后以生态保护重要性评价结果提取
生态保护优先区，接着运用由生态保护优先区所得限制发展区进行土地利用模拟进而划定开发边界，有效避免了现有技术忽略基于生态保护重要性评价进行弹性开发边界划定的问题。基于本发明技术划定的开发边界有助于保护优质生态空间、维持区域生态安全和保障人类福祉。同时创新提出基于保护效率和权衡度从多情景中筛选出最优情景的生态保护优先区的方法，为自然资源管理和生态保护的决策提供了有价值的技术指导。
[0202]
(2)本发明运用plus模型进行土地利用模拟进而划定城市群开发边界，提供了将plus 模型用于城市群开发边界划定的机制和方法，有效避免了现有技术忽略基于plus模型进行城市群弹性开发边界划定的问题，有助于弥补现有开发边界划定技术的不足，以及提高开发边界划定的科学性。同时本发明适用于城市群尺度的开发边界划定，弥补了现有技术在城市群等宏观尺度的开发边界划定技术研究方面的不足。
[0203]
(3)本发明提供了自然发展情景和生态保护情景的开发边界划定方法，满足了不同情景的边界划定方法需求，以便不同情景的边界划定结果对比分析，从而提高边界划定的科学性，有助于构建更加科学合理的国土空间开发保护格局。同时本发明步骤科学、通用且可复制，过程涉及的数据处理与空间分析，通过目前被生态、地理、规划等多领域广泛使用的arcgis 等软件即可实现，适用于任何地区特定时间框架下研究区城市群开发边界划定。
[0204]
本发明还公开了一种城市群弹性开发边界划定系统，包括：评价指标模块，用于建立生态保护重要性评价指标空间数据库；主成分分析模块，用于根据生态保护重要性评价指标空间数据库进行主成分分析；owa算子模块，用于基于owa算子进行多情景的生态保护重要性评价模拟；情景筛选模块，用于提取各情景的生态保护优先区并筛选出最优情景的生态保护优先区；驱动因子模块，用于建立土地利用变化驱动因子空间数据库；需求预测模块，用于预测自然发展情景和生态保护情景的目标年土地利用数量结构；模拟与验证模块，用于模拟基期年的土地利用空间结构并进行模拟精度验证；边界划定模块，用于模拟自然发展情景和生态保护情景的目标年土地利用空间结构并划定弹性开发边界。本发明的城市群弹性开发边界划定系统，可执行本发明的城市群弹性开发边界划定方法，可执行方法实施例的任意组合实施步骤，具备该方法相应的功能和有益效果。
[0205]
虽然在功能性模块的背景下描述了本发明，但应当理解的是，除非另有说明，功能和/或特征中的一个或多个可以被集成在单个物理装置和/或软件模块中，或者一个或多个功能和/ 或特征可以在单独的物理装置或软件模块中被实现。还可以理解的是，有关每个模块的实际实现的详细讨论对于理解本发明是不必要的。更确切地说，考虑到在本文中公开的系统中各种功能模块的属性、功能和内部关系的情况下，在工程师的常规技术内将会了解该模块的实际实现。因此，本领域技术人员运用普通技术就能够在无需过度试验的情况下实现在权利要求书中所阐明的本发明。还可以理解的是，所公开的特定概念仅仅是说明性的，并不意在限制本发明的范围，本发明的范围由所附权利要求书及其等同方案的全部范围来决定。
[0206]
功能如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计
算机，服务器，或者网络设备等)执行本发明各个实施例方法的全部或部分步骤。而前述的存储介质包括： u盘、移动硬盘、只读存储器(rom，read-only memory)、随机存取存储器(ram，randomaccess memory)、磁碟或者光盘等各种可以存储程序代码的介质。
[0207]
计算机可读介质的更具体的示例(非穷尽性列表)包括以下：具有一个或多个布线的电连接部(电子装置)，便携式计算机盘盒(磁装置)，随机存取存储器(ram)，只读存储器(rom)，可擦除可编辑只读存储器(eprom或闪速存储器)，光纤装置，以及便携式光盘只读存储器(cdrom)。另外，计算机可读介质甚至可以是可在其上打印程序的纸或其他合适的介质，因为可以例如通过对纸或其他介质进行光学扫描，接着进行编辑、解译或必要时以其他合适方式进行处理来以电子方式获得程序，然后将其存储在计算机存储器中。
[0208]
本发明的各部分可以用硬件、软件、固件或它们的组合来实现。在上述实施方式中，多个步骤或方法可以用存储在存储器中且由合适的指令执行系统执行的软件或固件来实现。例如，如果用硬件来实现，和在另一实施方式中一样，可用本领域公知的下列技术中的任一项或他们的组合来实现：具有用于对数据信号实现逻辑功能的逻辑门电路的离散逻辑电路，具有合适的组合逻辑门电路的专用集成电路，可编程门阵列(pga)，现场可编程门阵列(fpga) 等。
[0209]
上述说明是针对本发明较佳可行实施例的详细说明，但实施例并非用以限定本发明的专利申请范围，凡本发明所揭示的技术精神下所完成的同等变化或修饰变更，均应属于本发明所涵盖专利范围。
技术特征：1.一种城市群弹性开发边界划定方法，其特征在于，包括以下步骤：s1：选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系，开展生态保护重要性评价指标的具体评价，建立生态保护重要性评价指标空间数据库；s2：根据生态保护重要性评价指标空间数据库进行主成分分析，生成主成分指标数据和主成分分析报告，进行kmo和bartlett's检验以及主成分方差贡献率累积值检验，确定主成分指标数据及各主成分指标的准则权重；s3：设置多种决策风险系数情景，根据确定的主成分指标数据，计算主成分指标的位序权重和权衡度，运用owa算子聚合不同情景下的主成分指标的位序权重和准则权重，进行多情景的生态保护重要性评价模拟；s4：根据生态保护重要性评价模拟结果，提取各种情景下的生态保护优先区，对比不同情景的生态保护优先区对于各情景的生态保护重要性评价结果的保护效率和权衡度，筛选出保护效率与权衡度共同最优情景的生态保护优先区；s5：选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系，开展土地利用变化驱动因子的具体评价，建立土地利用变化驱动因子空间数据库；s6：根据历史土地利用数据，获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵，根据历史土地利用数据，提取各土地利用类型的扩张斑块数据，计算生态保护情景的土地利用转移概率矩阵，获得生态保护情景的目标年土地利用数量结构；s7：根据土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，输出各土地利用类型的发展概率数据，根据历史土地利用数据和各土地利用类型的发展概率数据，模拟基期年的土地利用空间结构数据，并进行模拟精度验证；s8：基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟自然发展情景目标年土地利用空间结构数据；基于基期年土地利用数据和各土地利用类型的发展概率数据，模拟生态保护情景的目标年土地利用空间结构数据，并划定弹性开发边界。2.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s1具体为：选取生态保护重要性的评价指标，构建以生态保护重要性评价为目标层、生态系统服务和生态脆弱性为准则层的生态保护重要性评价指标体系；结合多源数据、arcgis软件和invest模型，开展生态保护重要性评价指标的具体评价；评价结果统一进行归一化处理，以消除不同评价指标之间的量纲影响；借助arcgis平台整合所有归一化后的评价指标，建立研究区生态保护重要性评价指标空间数据库；其中，从属于生态系统服务的供给服务、调节服务、支持服务以及从属于生态脆弱性的暴露度、敏感性和适应力为一级指标层，选取的具体指标为生态保护重要性评价指标体系的二级指标层，根据对目标层的影响，所选具体指标中从属于适应力的指标为负向指标，其他指标均为正向指标；正向指标进行归一化的计算方法为：负向指标进行归一化的计算方法为：
式中，x
norm
为归一化后的指标像元值，x
x
为指标像元值，x
max
为指标所有像元值的最大值，x
min
为指标所有像元值的最小值。3.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s2具体为：使用qgis软件的对栅格取样工具，将所有评价指标的栅格图层取样至由研究所需分辨率的栅格转换而成的点图层，然后将点图层的属性表数据导入spss软件进行主成分分析，生成主成分指标数据和主成分分析报告；进行kmo和bartlett's检验，将特征值大于1作为选取主成分的标准，进行主成分方差贡献率累积值检验；确定基于生态保护重要性评价指标的主成分指标数据，并根据主成分方差贡献率确定和为1的各主成分指标的准则权重；其中，进行kmo和bartlett's检验，具体为：若主成分分析报告中的kmo和bartlett's检验结果显示：kmo取样适切性量数大于0.7，并且巴特利特球形度检验的显著度小于0.001，则研究数据适合进行主成分分析；否则，研究数据不适合进行主成分分析；进行主成分方差贡献率累积值检验，具体为：若主成分分析报告中的总方差解释结果显示：特征值大于1的主成分方差贡献率累积值超过60％，则生成的主成分指标数据较为合理且有效；否则，生成的主成分指标数据不合理且无效。4.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s3具体为：设置z(z≥5)种决策风险系数情景，根据确定的主成分指标数据，使用单调规则递增方法计算主成分指标的位序权重，具体计算方法如下：具体计算方法如下：式中，j为位序；v
j
为位序权重，v
j
∈[0,1]；n为主成分指标数量；α为决策风险系数，α∈(0,∞)；w
k
为主成分指标重要等级；r
k
为主成分指标的赋值，按照指标的像元值叠加准则权重后的大小对指标的进行赋值，最大值赋1，次大值赋2，最小值赋n；权衡度为不同决策风险系数下主成分指标间的补偿程度，具体计算方法如下：其中，tradeoff为权衡度，且0≤tradeoff≤1；n为主成分指标数量；w
k
为第k个主成分指标重要等级；运用terrset软件的mce模块的owa算子机制，聚合不同情景下的主成分指标的位序权重和准则权重，从而获得z种情景下的生态保护重要性空间评价结果。5.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s4具体为：基于z种情景下的生态保护重要性评价模拟结果，进行降序排列下的累积生态保护重
要性像元值与总的生态保护重要性像元值的比值计算，将比值为50％对应的生态保护重要性像元值作为各情景生态保护重要性阈值；在arcgis软件中使用聚合工具，将大于或等于生态保护重要性阈值的区域中相对聚集或邻近的图斑聚合为相对完整连片图斑，聚合距离为研究所需分辨率，聚合结果作为各情景初步生态保护优先区；为降低生态保护优先区的破碎化程度，利用情景z1的生态保护优先区图斑确定需要剔除的图斑面积阈值；其中，情景z1为z种情景中各主成分指标位序权重相等的情景，统计分析情景z1的生态保护优先区中面积小于阈值的图斑与所有图斑的数量比值，以及其总面积随面积阈值变化的情况，将拐点对应的最小面积阈值作为最终的图斑面积阈值，并根据阈值剔除相应独立分散小图斑，在此基础上合并国家级和省级的自然保护区、森林公园、风景名胜区等现有各类生态保护区，最终获得各情景最终的生态保护优先区；对比不同情景最终的生态保护优先区对于情景z1的生态保护重要性评价结果的保护效率和权衡度，从而筛选出保护效率与权衡度共同最优情景的最终生态保护优先区；其中，保护效率的计算方法如下：式中，z为情景序号，z＝1,2，...z，p
z
为情景z的保护效率，为情景z的生态保护优先区范围内情景z1的生态保护重要性平均值，为整个研究区范围内情景z1的生态保护重要性平均值。6.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s5具体为：选取土地利用变化的驱动因子，构建以自然因素和人文因素为准则层的土地利用变化驱动因子体系；结合多源数据和arcgis软件，开展土地利用变化驱动因子的具体评价；所得评价结果统一进行归一化处理，以消除不同驱动因子之间的量纲影响；借助arcgis平台整合所有归一化后的驱动因子，建立研究区土地利用变化驱动因子空间数据库；其中，从属于自然因素的地形和水域以及从属于人文因素的人口、经济、区位、基础设施为一级指标层，基础设施包括交通、医疗、教育和商服，选取的具体指标为二级指标层的土地利用变化驱动因子体系。7.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s6具体为：基于两期历史土地利用数据，运用plus模型的markov chain模块获得自然发展情景的目标年土地利用数量结构和土地利用转移概率矩阵；基于两期历史土地利用数据，运用plus模型的提取土地利用扩张模块提取各土地利用类型的扩张斑块数据，以此为基础计算城镇建设用地扩张斑块面积中各非建设用地的占比，并将所得占比数据作为非建设用地向城镇建设用地转移概率的负增长率，从而基于自然发展情景的土地利用转移概率矩阵计算生态保护情景的土地利用转移概率矩阵；基于markov chain原理和生态保护情景的土地利用转移概率矩阵计算生态保护情景的目标年土地利用数量结构；其中，历史土地利用数据包括耕地、林地、草地、水域、城镇建设用地、农村居民点、未利用地共七种土地利用类型，其中城镇建设用地和农村居民点属于建设用地，其余土地利用
类型属于非建设用地。8.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s7具体为：基于土地利用变化驱动因子空间数据库和各土地利用类型的扩张斑块数据，运用plus模型的leas模块输出各土地利用类型的发展概率数据；计算和设置plus模型的cars模块所需参数；以第二期历史土地利用数据的土地利用数量结构为土地利用需求数据，基于第一期历史土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟基期年的土地利用空间结构；基于两期历史土利用数据和基期年土地利用空间结构数据，运用plus模型confusion matrix and fom模块进行模拟精度验证；其中，计算和设置plus模型的cars模块所需参数，具体为：基于各土地利用类型的扩张斑块数据计算各土地利用类型的扩张斑块面积在总扩张面积中的占比，将所得占比数据作为邻域权重参数；将土地转移矩阵中各土地利用类型的是否允许转移值设为1，不限制任何土地利用类型之间的转移；递减阈值的衰减系数、随机斑块种子的概率、随机斑块种子最大比例、邻域范围、并行线程数量等参数采用默认值；进行模拟精度验证，具体为：若kappa系数大于0.8且fom系数大于0.2，则plus模型的模拟结果和真实结果基本一致，模拟结果可信；否则，需不断调整plus模型的输入数据和相关参数，直至满足上述kappa系数和fom系数检验要求才能继续进行后续步骤。9.根据权利要求1的城市群弹性开发边界划定方法，其特征在于，步骤s8具体为：以自然发展情景目标年土地利用数量结构为土地利用需求数据，以及现有的各类生态保护区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟自然发展情景目标年土地利用空间结构数据；以生态保护情景目标年土地利用数量结构为土地利用需求数据，以及生态保护优先区为限制发展区，基于基期年土地利用数据和各土地利用类型的发展概率数据，运用plus模型的cars模块模拟生态保护情景目标年土地利用空间结构数据；将自然发展情景和生态保护情景目标年土地利用空间结构数据中城镇建设用地边界，分别划定为研究区自然发展情景和生态保护情景目标年的弹性开发边界。10.一种城市群弹性开发边界划定系统，其特征在于，包括：评价指标模块，用于建立生态保护重要性评价指标空间数据库；主成分分析模块，用于根据生态保护重要性评价指标空间数据库进行主成分分析；owa算子模块，用于基于owa算子进行多情景的生态保护重要性评价模拟；情景筛选模块，用于提取各情景的生态保护优先区并筛选出最优情景的生态保护优先区；驱动因子模块，用于建立土地利用变化驱动因子空间数据库；需求预测模块，用于预测自然发展情景和生态保护情景的目标年土地利用数量结构；模拟与验证模块，用于模拟基期年的土地利用空间结构并进行模拟精度验证；边界划定模块，用于模拟自然发展情景和生态保护情景的目标年土地利用空间结构并划定弹性开发边界。
技术总结本发明公开了一种城市群弹性开发边界划定方法，包括：S1、建立生态保护重要性评价指标空间数据库；S2、根据评价指标空间数据库进行主成分分析；S3、基于OWA算子进行多情景的生态保护重要性评价模拟；S4、提取各情景的生态保护优先区并择优；S5、建立土地利用变化驱动因子空间数据库；S6、预测自然发展情景和生态保护情景的目标年土地利用数量结构；S7、模拟基期年的土地利用空间结构并进行模拟精度验证；S8、模拟自然发展情景和生态保护情景的目标年土地利用空间结构，并划定弹性开发边界。本发明模拟多种生态保护重要性评价情景，基于保护效率与权衡度筛选出最优情景的生态保护优先区，模拟了目标年的土地利用空间结构，划定城市群弹性开发边界。市群弹性开发边界。市群弹性开发边界。
技术研发人员：陈康富 吴隽宇
受保护的技术使用者：华南理工大学
技术研发日：2022.07.04
技术公布日：2022/11/1
