[0001] 本发明涉及数据处理技术领域，尤其涉及一种土地利用变化模拟方法。

[0002] 土地利用变化模拟是土地变化科学中重要的组成部分，能够根据当前或特定时间的发展规律推断未来土地空间格局。通过对未来不同发展情景假设的分析，可利用土地利用变化模拟判断不同发展情景假设、规划、政策等对土地空间格局的影响。

[0003] 现有技术中，通常使用CLUMondo模型进行土地利用变化模拟。但是CLUMondo模型使用Logit回归法来建立适宜性参数与生物物理、社会经济因子的线性关系，不够贴合实际，影响模拟结果的准确性。因此，亟需一种有效的方法以解决上述问题。

[0020] 为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

[0021] 为了便于更加清晰地理解本发明各实施例，首先对一些相关的背景知识进行如下介绍。

[0022] 人类对土地资源的开发和利用程度不断增强，已成为地球表层尤其是陆地表层变化的主要原因。然而人类对土地的改造所带来的环境风险和资源短缺问题已经威胁到人类自身的生存和发展，如何合理利用土地资源已成为当今时代的重要议题。土地利用变化模拟是土地变化科学中重要的组成部分，能够根据当前或特定时间的发展规律推断未来土地空间格局。通过对未来不同发展情景假设的分析，可利用土地利用变化模拟判断不同发展情景假设、规划、政策等对土地空间格局的影响。土地利用变化模拟也是其他研究领域的重要支撑，例如全球环境变化模型的重要输入和环境管理、生物多样性、碳中和等研究领域的核心关注点。

[0023] 在过去的几十年中，科学家们开发了许多不同类型的土地利用变化模拟模型，以实现更加科学和符合预期需求的未来土地利用变化模拟。土地利用变化模拟模型可以分为四种类型：经验统计模型、随机模型、优化模型和基于土地动力学机制的模型，其中，经验统计模型使用回归方法建立土地与可能的外部驱动因素之间的定量关系，被广泛应用于土地变化模拟研究中。

[0024] 在经验统计模型中，CLUMondo模型相比于大多数的其他经验统计模型更具优势，可以支持将土地所承载的功能或服务作为模拟的终止条件。另外，CLUMondo模型支持在土地类型及其服务之间建立多对多的对应关系，通过迭代实现供需平衡，这使得服务可以由提供它们的多种土地类型面积满足，而不是被限制于特定某种土地类型。

[0025] 但是CLUMondo模型使用Logit回归法来建立适宜性参数与生物物理、社会经济因子的线性关系，即CLUMondo模型在建立适宜性参数与变化驱动因子的关系时使用了线性模型，内部运行机制也存在不足，不够贴合实际，影响模拟结果的准确性。并且CLUMondo模型的原始内部运行机制使得每次迭代都会引起很大的变化，这样的框架在早期迭代中有利于需求和供给的快速收敛，但当需求和供给之间的差异很小时则可能不利于模型收敛。

[0026] 因此，本发明提供了一种土地利用变化模拟方法，通过获取目标区域中各子区域的土地系统类型、所述目标区域对应的各变化驱动因子和各土地服务需求量；根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，可以建立土地适宜性和驱动因子之间的复杂和非线性关系，计算各所述土地系统类型的土地服务能力；根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势；根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响；并根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力；将所述土地利用现状数据、各土地转换潜力、各土地服务能力和各土地服务需求量，输入至改进的土地利用模型进行土地利用变化模拟，得到模拟结果，能够基于逐步迭代的方法促进土地利用模型内部迭代实现供需平衡，进而提高土地利用变化模拟的效率和可靠性。

[0027] 下面结合图1对本发明提供的土地利用变化模拟方法及装置进行说明。

[0028] 图1是本发明提供的土地利用变化模拟方法的流程示意图，参见图1所示，包括步骤101‑步骤105，其中：步骤101：获取目标区域中各子区域的土地系统类型、所述目标区域对应的各变化驱动因子和各土地服务需求量。

[0029] 首先需要说明的是，本发明的执行主体可以是模拟土地利用变化的任何电子设备，例如可以为智能手机、智能手表、台式电脑、手提电脑等任何一种。

[0030] 具体地，目标区域是指需要进行土地利用变化模拟的土地，也即研究区。子区域是指对目标区域进行划分的栅格。土地系统类型，也即土地系统数据，是指子区域在土地利用以及密度两个维度上的类型。变化驱动因子包括土地利用变化驱动因子及其对应的数值，是社会经济和生物物理数据对土地利用变化的驱动，变化驱动因子包括八类，分别是社会经济因子、土壤因子、可达性因子、农业和植被因子、地形因子、气候因子、牲畜因子和占比因子。土地服务是指土地所提供的使用服务，包括森林、草地/牧场、灌木地和耕地四类。各土地服务需求量是指对每种土地服务的供给量，如需求面积。

[0031] 实际应用中，可以先确定目标区域，然后将目标区域划分成多个子区域，进而获取各子区域的土地系统类型、目标区域对应的各变化驱动因子和各土地服务需求量。

[0032] 获取各子区域的土地系统类型、目标区域对应的各变化驱动因子和各土地服务需求量的方法有多种，例如，用户通过上传页面上传各子区域的土地系统类型、目标区域对应的各变化驱动因子和各土地服务需求量，相应地，执行主体获取到各子区域的土地系统类型、目标区域对应的各变化驱动因子和各土地服务需求量；又如，执行主体接收到土地利用变化模拟指令或者获取指令，相应地，执行主体从土地利用变化模拟指令或者获取指令所指向的存储区获取各子区域的土地系统类型、目标区域对应的各变化驱动因子和各土地服务需求量。本发明对此不作限定。

[0033] 此外，还可以通过全球变化分析模型获取各土地服务需求量，即从全球变化分析模型的全球模拟结果中获取未来情景下的各土地服务需求量。

[0034] 示例性地，从全球变化分析模型v5.3的模拟结果中，获取未来情景下的土地利用情况。全球变化分析模型v5.3输出70种土地类型，将70种土地类型划分为林地、草地/牧场、灌木地和耕地四类土地服务，得到未来情境下每个流域对这四种土地服务的需求量，即各土地服务需求量。

[0035] 步骤102：根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，并计算各所述土地系统类型的土地服务能力。

[0036] 具体地，土地适应性是指当前土地系统类型被分配到当前子区域的适宜性。土地服务能力是指当前土地系统类型能提供特定土地服务的能力，特征土地服务为任一土地服务。

[0037] 实际应用中，在获取了各土地系统类型和各变化驱动因子的基础上，进一步地，根据各土地系统类型和各变化驱动因子，使用随机森林算法，计算各土地系统类型被分配到某个子区域的适宜性，获取土地适宜性。同时，基于目标区域，计算各土地系统类型的土地服务能力。

[0038] 需要说明的是，为了便于排错，可以将所有变化驱动因子的原有格式转换为ASCII格式。

[0039] 步骤103：根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势。

[0040] 具体地，竞争优势是指当前子区域相对于其他子区域拥有的提供土地服务的优势。

[0041] 实际应用中，得到了各土地系统类型、各土地服务能力和各土地服务需求量的基础上，进一步地，可以将各土地系统类型、各土地服务能力和各土地服务需求量按照设定的竞争优势算法进行计算，得到各子区域的竞争优势。

[0042] 步骤104：根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响；并根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力。

[0043] 具体地，土地利用现状数据表征当前目标区域的土地使用情况，土地利用类型是指子区域在土地利用维度上的类型，土地利用类型包括十类，分别是耕地、林地、草地、灌木地、湿地、水体、苔原、人造地表、裸地和雪/冰。转换阻力通过土地系统数据判断每种土地系统类型能否转换为其他土地系统类型以及转换的难易程度。邻域影响考虑到子区域受到周围子区域（相邻子区域）的影响，所谓影响是指周围子区域的土地系统类型影响到该子区域的土地系统类型的判断。

[0044] 实际应用中，可以按照设定的转换阻力算法，以土地利用现状数据为输入，计算各土地系统类型的转换阻力，并按照设定的邻域影响算法，以土地利用现状数据为输入，计算各土地系统类型的邻域影响。

[0045] 在得到各所述土地适宜性、各转换阻力、各邻域影响、各竞争优势的基础上，进一步地，可以将各土地适宜性、各转换阻力、各邻域影和各竞争优势输入至改进的CLUMondo模型，计算得到土地转换潜力；也可以按照设定的计算规则，以各土地适宜性、各转换阻力、各邻域影和各竞争优势为输入，计算各子区域的土地转换潜力。

[0046] 步骤105：将所述土地利用现状数据、各所述土地转换潜力、各所述土地服务能力和各所述土地服务需求量，输入至改进的土地利用模型进行土地利用变化模拟，得到模拟结果。

[0047] 具体地，土地利用模型也即土地利用变化模拟模型，优选地，改进的土地利用模型为改进的CLUMondo模型。

[0048] 实际应用中，在得到各土地转换潜力、各土地服务能力和各土地服务需求量的基础上，进一步地，将各土地系统数据、各所述土地服务能力、各所述土地服务需求量和土地转换潜力输入改进的CLUMondo模型进行土地利用变化模拟，由改进的土地利用模型对土地系统数据、各所述土地服务能力、各所述土地服务需求量和土地转换潜力进行嵌套迭代后，得到土地利用变化模拟的模拟结果。

[0049] 本发明提供的土地利用变化模拟方法，通过获取目标区域中各子区域的土地系统类型、所述目标区域对应的各变化驱动因子和各土地服务需求量；根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，可以建立土地适宜性和驱动因子之间的复杂和非线性关系，计算各所述土地系统类型的土地服务能力；根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势；根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响；并根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力；将土地系统数据、各土地转换潜力、各土地服务能力和各土地服务需求量，输入至改进的土地利用模型进行土地利用变化模拟，得到模拟结果，能够基于逐步迭代的方法促进土地利用模型内部迭代实现供需平衡，进而提高土地利用变化模拟的效率和可靠性。

[0050] 在本发明一个或更多个可选的实施例中，所述获取目标区域中各子区域的土地系统类型，包括：根据所述土地利用现状数据，确定所述目标区域中各子区域的土地利用类型；
针对每个所述土地利用类型，确定所述土地利用类型对应的各子区域的面积占
比，所述面积占比表征目标面积与所述子区域的总面积的比值，所述目标面积表征所述子区域中属于所述土地利用类型的土地面积；
根据各所述土地利用类型和各所述面积占比，确定所述目标区域中各子区域的土地系统类型。

[0051] 实际应用中，先获取目标区域的土地利用现状数据，然后针对每个子区域，从土地利用现状数据中识别该子区域对应的子土地利用现状数据，并根据该子土地利用现状数据识别该子区域的土地利用类型：根据子土地利用现状数据，确定该子区域中各初始土地利用类型对应的土地面积，将最大土地面积对应的初始土地利用类型确定为该子区域的目标土地利用类型，即子区域的目标土地利用类型为子区域范围内面积占比最高或面积最大的初始土地利用类型。

[0052] 进一步地，针对每个土地利用类型，从密度维度上将当前土地利用类型进行细分，进而得到目标区域中各子区域的土地系统类型：针对每个土地利用类型，确定属于该土地利用类型的目标子区域，然后针对每个目标子区域，确定该目标子区域中该土地利用类型对应的目标面积与该目标子区域的比值，即面积占比，遍历所有目标子区域，得到该土地利用类型对应的各面积占比，然后将面积占比进行密度划分，将该土地利用类型划分为多个土地系统类型，进而得到各目标子区域的土地系统类型。遍历所有土地利用类型，得到目标区域中各目标子区域的土地系统类型。如此，可以基于面积占比确定土地系统类型，可以提高土地系统类型的精准度和可靠性。

[0053] 示例性地，以目标区域为A地区进行说明，利用数据中心发布的最新的Globeland30土地利用数据和A地区的边界数据裁剪下A地区的数据，即土地利用现状数据，其中，Globeland30为30米空间分辨率全球地表覆盖数据。将Globeland30土地利用数据从
30米分辨率聚合到990米分辨率。对于每个990*990米的大栅格，即子区域，定义其土地利用类型为范围内面积占比最高的土地利用类型。针对每种土地利用类型，统计目标区域中所有类型为该土地利用类型的大栅格中，该土地利用类型实际的面积比例，即面积占比。如森林类型的大栅格，其中有森林的面积为900*900，则占比为900*900/（990*990）=100/121。在确定该土地利用类型对应的各子区域的面积占比，根据将面积占比进行密度划分，将该土地利用类型划分为多个土地系统类型，进而得到各目标子区域的土地系统类型，即将土地利用类型进一步细分为低、中、高三个密度级别，得到三个土地系统类型，例如土地利用类型为耕地，则对应的三个土地系统类型分别为低密度耕地、中密度耕地和高密度耕地，进一步得到各目标子区域的土地系统类型。

[0054] 在本发明一个或更多个可选的实施例中，所述根据各所述土地利用类型和各所述面积占比，确定所述目标区域中各子区域的土地系统类型，具体实现过程可以如下：针对每个所述土地利用类型，对所述土地利用类型对应的各子区域的面积占比，调用自然断点算法，确定第一阈值和第二阈值；
根据所述第一阈值和第二阈值，将所述土地利用类型划分为三个土地系统类型；
将所述土地利用类型对应的各子区域的面积占比分别与所述三个土地系统类型
进行比较，确定所述土地利用类型对应的各子区域的土地系统类型。

[0055] 具体地，自然断点算法可以将对分类间隔加以识别，可对相似值进行最恰当地分组，并可使各个类之间的差异最大化；要素将被划分为多个类，对于这些类，会在数据值的差异相对较大的位置处设置其边界，即阈值。

[0056] 实际应用中，针对每个土地利用类型，使用自然断点算法将该土地利用类型对应所有面积比值分为三组，即确定第一阈值和第二阈值。

[0057] 在第一阈值大于第二阈值的情况下：将面积占比大于第一阈值的子区域的土地利用类型的密度级别确定为高密度级别，进而确定该子区域的土地系统类型，例如土地利用类型为湿地，则该子区域的土地系统类型为高密度湿地；将面积占比大于第二阈值且小于等于第一阈值的子区域的土地利用类型的密度级别确定为中密度级别，进而确定该子区域的土地系统类型，例如土地利用类型为草地，则该子区域的土地系统类型为中密度草地；将面积占比小于等于第二阈值的子区域的土地利用类型的密度级别确定为低密度级别，进而确定该子区域的土地系统类型，例如土地利用类型为水体，则该子区域的土地系统类型为低密度水体。

[0058] 在第一阈值小于第二阈值的情况下：将面积占比大于第二阈值的子区域的土地利用类型的密度级别确定为高密度级别，进而确定该子区域的土地系统类型；将面积占比大于第一阈值且小于等于第二阈值的子区域的土地利用类型的密度级别确定为中密度级别，进而确定该子区域的土地系统类型；将面积占比小于等于第一阈值的子区域的土地利用类型的密度级别确定为低密度级别，进而确定该子区域的土地系统类型。

[0059] 具体地，针对每个所述土地利用类型，对所述土地利用类型对应的各子区域的面积占比，调用自然断点算法，确定第一阈值和第二阈值，具体实现过程如下：计算各面积占比的第一方差，并确定将各面积占比分为三组的所有分组结果；
针对每个分组结果，计算所述分组结果中每个分组的第二方差，并将各所述第二方差相加得到所述分组结果的方差和；
根据各所述分组结果的方差和，计算各所述分组结果的方差拟合优度；
根据最大方差拟合优度对应分组，确定第一阈值和第二阈值。

[0060] 具体地，方差是指平均值的偏差平方和。

[0061] 实际应用中，先按照 计算当前土地利用类型对应的各面积占比的第一方差，其中，SDAM为第一方差xi是第i个面积占比，即当前土地利用类型对应的第i个子区域中该土地利用类型的实际面积比例值， 为该土地利用类型对应的所有面积占比的平均值。

[0062] 然后，对于每个可能的分组组合，即分组结果（每个分组结果中包含三个分组），按照 计算每个分组结果的方差和，并找到方差和最小的分组结果，其中，SDCM为方差和，xij是第j组内中第i个面积占比， 是第j组所有面积占比的平均值，为第二方差。由于将土地利用类型分为三个密度级别，因此j为小于等于3的正
整数。

[0063] 进一步地，按照GVF=（SDAM‑SDCM）/SDAM计算各所述分组结果的方差拟合优度，其中，GVF为方差拟合优度。GVF的取值范围为1（分组最好）到0（分组最差）。找到GVF最高的一种分组结果，确定密度级别的两个阈值，即第一阈值和第二阈值。

[0064] 如此，采用自然断点算法，可以提高第一阈值和第二阈值的可靠性和效率，进而提高确定土地系统类型的可靠性和效率。

[0065] 需要说明的是，在确定每一个子区域的土地系统类型后，需要将空间分辨率进行调整，以匹配大部分变化驱动因子的空间分辨率。沿用上例，在确定每一个子区域的土地系统类型后，将空间分辨率由990m*990m重采样至1km*1km，重采样使用最邻近重采样方法，每一个1km*1km的栅格属性值都是距离该栅格中心最近的990m*990m数据样本中心所在的栅格属性值。

[0066] 在本发明一个或更多个可选的实施例中，所述根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，具体实现过程可以如下：针对每个所述子区域，令c=1，调用随机森林算法，基于各所述变化驱动因子，确定第c个所述土地系统类型被分配到所述子区域的土地适宜性；
在c小于N的情况下，令c=c+1，继续执行所述调用随机森林算法，基于各所述变化驱动因子，确定第c个所述土地系统类型被分配到所述子区域的土地适宜性的步骤，N为所述土地系统类型的数量。

[0067] 实际应用中，使用随机森林方法计算土地适宜性：P_locc，j，该参数表示第j个土地系统类型被分配到第c个子区域的适宜性。计算的方法如下：首先准备子区域c上的变化驱动因子（数值）：（X1,c,X2,c,…,Xm,c），其中，m表示变化驱动因子的数量，c表示第c个子区域。对每种土地系统类型j，建立一个随机森林模型，训练随机森林模型时使用的响应变量是0,
1二值的，当子区域c当前的土地系统类型是第j个土地系统类型，则子区域c对应的响应变量的值确定为1，否则为0。所有子区域上的变化驱动因子和响应变量共同构成了样本集。通过训练随机森林模型，建立驱动因子与土地适宜性：P_locc，j之间的关系。进一步地，遍历所有子区域和所有土地系统类型。如此，通过耦合随机森林算法，建立土地适宜性和变化驱动因子之间的复杂和非线性关系，能够更准确地反映变化驱动因子和土地适宜性之间的关系，有利于提高土地利用变化模拟的可靠性。

[0068] 可选地，为提高土地适宜性的准确度，所有的随机森林模型都设置有200棵树。

[0069] 在本发明一个或更多个可选的实施例中，所述计算各所述土地系统类型的土地服务能力，包括：基于所述目标区域中各子区域的土地系统类型，确定所述目标区域的土地系统
图，并基于所述目标区域的土地利用现状数据，确定所述目标区域的土地现状图，所述土地系统图的分辨率高于所述土地现状图的分辨率；
根据所述土地系统图和所述土地现状图，计算各所述土地系统类型的土地服务能力。

[0070] 实际应用中，先基于目标区域中各子区域的土地系统类型得到目标区域的土地系统图，并基于目标区域的土地利用现状数据得到目标区域的土地现状图。然后将土地系统图和土地现状图进行叠加，按照CAj,d=Sj,d/Sj计算出各土地系统类型的土地服务能力，其中，CAj,d是第j种土地系统类型相对于第d种土地服务（森林、草地/牧场、灌木地、耕地中的一种）的土地服务能力，Sj是目标区域内第j种土地系统类型的总面积，Sj,d是在第j种土地系统类型上第d种土地服务的供给量。如此，可以快速计算土地系统类型的土地服务能力，提高土地服务能力的准确度，有利于提高土地利用变化模拟的准确度。

[0071] 示例性地，通过叠加A地区的土地系统图（分辨率为990m）和Globeland30地图（分辨率为30m，即土地现状图），分别计算目标区域的30种土地系统类型（10种土地利用类型分别被划分三种密度级别后，得到30种土地系统类型）的土地服务能力CAj,d。

[0072] 此外，Sj,d的计算方法符合是对下述条件的土地利用类型面积求和：（1）30m*30m的子区域与土地系统图上的第j种土地系统类型重叠；（2）30m*30m的土地利用类型能够提供第d种土地服务。

[0073] 在本发明一个或更多个可选的实施例中，所述根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势，具体实现过程可以如下：针对每个所述子区域，将所述子区域对应的各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量输入至竞争优势函数中进行计算，得到所述子区域的竞争优势；
所述竞争优势函数通过以下公式（1）表示：
（1）
其中，P_cmpc,i,j表示第i次内部迭代中第c个子区域中第j个土地系统类型的竞争优势；CAj,d表示第j个土地系统类型提供第d种土地服务的土地服务能力；CAu,d表示第u个土地系统类型提供第d种土地服务的土地服务能力；第u个土地系统类型是第c个子区域的土地系统类型，也即土地系统现状类型；参数inertiad,i表示第（i‑1）次内部迭代结束时，第d种土地服务和各所述土地系统类型提供的服务之间的累积差异，是一种自适应的惯性机制。

[0074] 实际应用中，在得到根据各土地系统类型、各土地服务能力和各土地服务需求量的基础上，进一步，将各土地系统类型、各土地服务能力和各土地服务需求量代入竞争优势函数进行计算，得到子区域的竞争优势。如此，可以快速准确地确定竞争优势，进而提高土地利用变化模拟的效率和可靠性。

[0075] 其中，惯性机制inertiad,i通过以下公式（2）表示：（2）
其中，Demandd表示对第d种土地服务的需求量，即第d个土地服务需求量；
Supplyd,i‑1表示第（i‑1）次内部迭代结束时的土地系统类型对第d种土地服务的供给量。是速度参数，通过以下公式（3）表示：
（3）
公式（2）和公式（3）中，inertiad,i为第（i‑1）次内部迭代结束时，预设的第d种土地服务和土地系统类型提供的服务之间的累积差异；它受到速度参数 的控制；
由一个随机种子触发，并随着内部迭代次数的增加而增加，优选地，参数Seed和Step的初始值分别为1和0.001。

[0076] 在本发明一个或更多个可选的实施例中，所述根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响，具体实现过程可以如下：基于所述目标区域的土地利用现状数据，评估各所述土地系统类型转换为其他土地系统类型的难易程度，获取各所述土地系统类型的转换阻力；
基于所述土地利用现状数据和各所述子区域的相邻子域的土地系统类型，计算各所述子区域的土地系统类型转化分别为其他土地系统类型的邻域影响。

[0077] 实际应用中，可以基于所述土地利用现状数据，按照设定评估策略，评估第j个土地系统类型转换为其他土地系统类型的难易程度，获取第j个土地系统类型的转换阻力，其中，j为不大于30的正整数。转换阻力 的范围为[0,1]。土地系统类型的转换阻力数值为0表示其极容易转换，数值1则表示禁止转换。

[0078] 此外，需要基于土地系统数据，并结合第c个子区域的各邻域（相邻子域）的土地系统类型，计算第c个子区域的土地系统类型转化为第j个土地系统类型的可能性，也即转换阻力 。 分值越大，说明其被同质化为第j种土地系统类型的可能性越大。

[0079] 在本发明一个或更多个可选的实施例中，所述根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力，具体实现过程可以如下：针对每个子区域，若所述子区域的土地系统类型与待转换土地系统类型相同，则将所述子区域对应的所述土地适宜性、所述转换阻力、所述邻域影响与所述竞争优势之和作为土地转换潜力；若所述子区域的土地系统类型与待转换土地系统类型不同，则将所述子区域对应的所述土地适宜性、所述邻域影响与所述竞争优势之和作为土地转换潜力。

[0080] 具体地，待转换土地系统类型是指将当前的土地系统类型进行转换后的土地系统类型。

[0081] 实际应用中，针对第c个子区域，且待转换土地系统类型为第j个土地系统类型，若第c个子区域的土地系统类型为第j个土地系统类型，则土地转换潜力 等于土地适宜性、转换阻力、邻域影响与竞争优势之和；若第c个子区域的土地现状类型不为第j个土地系统类型，则土地转换潜力 等于土地适宜性、邻域影响与竞争优势之和。

[0082] 实际应用中，参见图2，图2是本发明提供的改进的CLUMondo模型迭代的流程示意图：将各土地适宜性、未来的土地需求（各土地服务需求量）、土地服务能力和竞争优势等输入改进的CLUMondo模型，用于进行基于土地供需服务平衡的土地变化模拟。改进的CLUMondo模型的迭代过程中，引入更精细的迭代，嵌套在每个外部迭代中。如果满足以下条件，更精细的迭代将被触发：一个外部迭代的所有内部迭代被执行，且改进的CLUMondo模型迭代出的土地系统不能满足未来的土地服务需求，即模型未收敛。更精细的迭代首先在所有可以改变的子区域中以及在所有可能的土地系统类型中确定最高的土地转换潜力，即 ，也即 ，其中 表示任意子区域c转化为任意土地系统类型k的土地转换潜力；q是计数变量，取值范围是1至n（与子区域数一致），用来在一次迭代中遍历所有的子区域；k是第k种土地系统类型。然后，它将第 种土地系统类型分配给第 个子区域，接着检查改进的CLUMondo模型是否收敛。如果没有达到收敛，第 个子区域将从进一步考虑中删除，并重复之前的步骤。之后迭代计算，直到改进的CLUMondo模型使模拟的土地数量符合土地需求，或所有的子区域都在更精细的迭代中尝试过土地系统转换。
如此，通过改进的CLUMondo模型内部运行机制，使用改进竞争优势参数和更精细的迭代，让土地利用模型能够精细地调整，可以促进土地利用模型迭代达到供需平衡，提高了土地利用变化模拟的效率和可靠性。

[0083] 更具体地，参见图2：改进的CLUMondo模型通过两个嵌套的迭代来模拟土地的变化：外部迭代t可以涵盖全部的模拟时间段；当需要输出中间某时刻的土地变化模拟结果时，外部迭代t也可以仅涵盖部分的模拟时间段。其中，t的范围为1至T。每一个外部迭代都包含大量的内部迭代，优选地，一个外部迭代中包含20000次内部迭代i，即i的范围为1至20000。

[0084] 在每一次内部迭代中，改进的CLUMondo模型按照公式（4）决定每个子区域的土地系统类型是否发生改变，改变为何种土地系统类型：（4）
公式（4）中，c表示第c个子区域，t表示第t个外部迭代，i表示第t个外部迭代中的第i次内部迭代。 表示第c个子区域在第t个外部迭代中的第i次内部迭代结束后的土地系统类型。 是土地转换潜力，代表第c个子区域的土地系统类型转换（或继续保持）为第j种（1≤j≤n，n为土地系统类型的总数）土地系统类型的概率。φ是研究区（目标区域）中禁止发生土地系统变化的区域（例如生态保护区）。(T(c, t, 0))用于计数，记录第c个子区域的土地系统类型保持与土地系统图中第c个子区域的土地系统类型，即 一致的外部迭代次数。(T(c,t,0))表示第c个子区域的土地系统类型应保持不变的最小外部迭代次数，是一个非负整数。 控制是否允许土地系统类型 转换
为土地系统类型k或第k种土地系统类型。

[0085] 也即，要计算“ ”，先判断c∈φ是否成立，若是，则 = ，若否，继续判断 (T(c, t, 0))＜ (T(c,t,0))是否成立，如果是，则 = ，
如果否，判断 是否成立，若否，则 = ，若是，则
=  。进而判断改进的CLUMondo模型是否收敛，若
是，则中断或停止；若否，则令i=i+1，继续内部迭代。如果i=20000后，改进的CLUMondo模型还是没有收敛，则启动更精细的迭代，针对i（从20001至20001+n），计算
，令 ，其中 和 是使得土地转换潜力 取最大值的参数组合，即第 个
子区域转换为 土地系统类型的土地转换潜力最大。再判断改进的CLUMondo模型是否收敛，若是，则中断或停止；若否，则令 ，继续更精细的迭代。

[0086] 决定 的关键是土地转换潜力 ，其由模型中四个影响土地变化的参数相加得到，四个参数分别为土地适宜性 ，转换阻力 ，邻域影响和竞争优势 。土地转换潜力的计算公式以下公式（5）表示：
（5）
其中转换阻力通过历史土地数据判断过往每种土地系统类型能否转换为其他土
地系统类型以及转换的难易程度，每种土地系统类型对应1个弹性系数，范围为[0,1]。一种土地系统类型的转换阻力数值为0表示其极容易转换，数值1则表示禁止转换。

[0087] 邻域影响则不仅考虑了每个地块自身，而且考虑到该地块周围地块的影响，所谓影响是指周围地块类型影响到中心地块类型转换的判断， 分值越大，说明其被同质化为土地系统类型j的可能性越大。

[0088] 将所述土地适宜性、邻域影响、转换阻力和改进竞争优势计算结果输入所述改进的CLUMondo模型，计算土地转换潜力。输入起始年份的土地系统后，模型基于土地转换潜力配置土地系统，根据每次迭代后每种土地系统类型的面积和土地服务供给能力，计算出土地服务总供给量。当总供给量与预测年份总需求量的差异在阈值内，模拟结束，否则模型根据土地服务需求的供给情况和土地系统类型的竞争能力调整竞争优势，并重新计算土地转换潜力，依据新的土地转换潜力配置新一次迭代的土地系统，再次判断本次迭代后的土地是否满足预测年份的土地服务需求，直到迭代出的土地系统满足未来的土地服务需求。

[0089] 本发明提供的土地利用变化模拟方法，使用所述随机森林算法代替二元Logit回归算法，建立土地利用变化的变化驱动因子与土地系统现状之间的关系，计算土地适宜性，在土地利用模型模拟土地变化时根据建立的随机森林算法给出每种土地系统出现概率的参考，能够更准确地反映变化驱动因子和土地适宜性之间的复杂非线性关系。使用改进竞争优势参数和更精细的迭代，让模型能够精细地调整，实现供需平衡。

[0090] 下面对本发明提供的土地利用变化模拟装置进行描述，下文描述的土地利用变化模拟装置与上文描述的土地利用变化模拟方法可相互对应参照。

[0091] 图3是本发明提供的土地利用变化模拟装置的结构示意图，如图3所示，该土地利用变化模拟装置300包括：获取模块301、第一确定模块302、第二确定模块303、计算模块304和模拟模块305，其中：获取模块301，被配置为获取目标区域中各子区域的土地系统类型、所述目标区域对应的各变化驱动因子和各土地服务需求量；
第一确定模块302，被配置为根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，并计算各所述土地系统类型的土地服务能力；
第二确定模块303，被配置为根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势；
计算模块304，被配置为根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响；并根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力；
模拟模块305，被配置为将所述土地利用现状数据、各所述土地转换潜力、各所述土地服务能力和各所述土地服务需求量，输入至改进的土地利用模型进行土地利用变化模拟，得到模拟结果。

[0092] 本发明提供的土地利用变化模拟装置，通过获取目标区域中各子区域的土地系统类型、所述目标区域对应的各变化驱动因子和各土地服务需求量；根据各所述土地系统类型和各所述变化驱动因子，确定所述各子区域的土地适宜性，可以建立土地适宜性和驱动因子之间的复杂和非线性关系，计算各所述土地系统类型的土地服务能力；根据各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量，确定所述各子区域的竞争优势；根据所述目标区域的土地利用现状数据，计算各土地系统类型的转换阻力和邻域影响；并根据各所述土地适宜性、各所述转换阻力、各所述邻域影响和各所述竞争优势，计算各子区域的土地转换潜力；将所述土地利用现状数据、各所述土地转换潜力、各所述土地服务能力和各所述土地服务需求量，输入至改进的土地利用模型进行土地利用变化模拟，得到模拟结果，能够基于逐步迭代的方法促进土地利用模型内部迭代实现供需平衡，进而提高土地利用变化模拟的效率和可靠性。

[0093] 可选地，所述获取模块301，进一步被配置为：根据所述土地利用现状数据，确定所述目标区域中各子区域的土地利用类型；
针对每个所述土地利用类型，确定所述土地利用类型对应的各子区域的面积占
比，所述面积占比表征目标面积与所述子区域的总面积的比值，所述目标面积表征所述子区域中属于所述土地利用类型的土地面积；
根据各所述土地利用类型和各所述面积占比，确定所述目标区域中各子区域的土地系统类型。

[0094] 可选地，所述获取模块301，进一步被配置为：针对每个所述土地利用类型，对所述土地利用类型对应的各子区域的面积占比，调用自然断点算法，确定第一阈值和第二阈值；
根据所述第一阈值和第二阈值，将所述土地利用类型划分为三个土地系统类型；
将所述土地利用类型对应的各子区域的面积占比分别与所述三个土地系统类型
进行比较，确定所述土地利用类型对应的各子区域的土地系统类型。

[0095] 可选地，所述第一确定模块302，进一步被配置为：针对每个所述子区域，令c=1，调用随机森林算法，基于各所述变化驱动因子，确定第c个所述土地系统类型被分配到所述子区域的土地适宜性；
在c小于N的情况下，令c=c+1，继续执行所述调用随机森林算法，基于各所述变化驱动因子，确定第c个所述土地系统类型被分配到所述子区域的土地适宜性的步骤，N为所述土地系统类型的数量。

[0096] 可选地，所述第一确定模块302，进一步被配置为：基于所述目标区域中各子区域的土地系统类型，确定所述目标区域的土地系统
图，并基于所述目标区域的土地利用现状数据，确定所述目标区域的土地现状图，所述土地系统图的分辨率高于所述土地现状图的分辨率；
根据所述土地系统图和所述土地现状图，计算各所述土地系统类型的土地服务能力。

[0097] 可选地，所述第二确定模块303，进一步被配置为：针对每个所述子区域，将所述子区域对应的各所述土地系统类型、各所述土地服务能力和各所述土地服务需求量输入至竞争优势函数中进行计算，得到所述子区域的竞争优势；
所述竞争优势函数通过以下公式（1）表示：
（1）
其中，P_cmpc,i,j表示第i次内部迭代中第c个子区域中第j个土地系统类型的竞争优势；CAj,d表示第j个土地系统类型提供第种土地服务的土地服务能力；CAu,d表示第u个土地系统类型提供第d种土地服务的土地服务能力；第u个土地系统类型是第c个子区域的土地系统类型；参数inertiad,i表示第（i‑1）次内部迭代结束时，第d种土地服务和各所述土地系统类型提供的服务之间的累积差异。

[0098] 可选地，所述计算模块304，进一步被配置为：基于所述目标区域的土地利用现状数据，评估各所述土地系统类型转换为其他土地系统类型的难易程度，获取各所述土地系统类型的转换阻力；
基于所述土地利用现状数据和各所述子区域的相邻子域的土地系统类型，计算各所述子区域的土地系统类型转化分别为其他土地系统类型的邻域影响。

[0099] 可选地，所述计算模块304，进一步被配置为：针对每个子区域，若所述子区域的土地系统类型与待转换土地系统类型相同，则将所述子区域对应的所述土地适宜性、所述转换阻力、所述邻域影响与所述竞争优势之和作为土地转换潜力；若所述子区域的土地系统类型与待转换土地系统类型不同，则将所述子区域对应的所述土地适宜性、所述邻域影响与所述竞争优势之和作为土地转换潜力。

[0100] 图4示例了一种电子设备的实体结构示意图，如图4所示，该电子设备可以包括：处理器（processor）410、通信接口（Communications Interface）420、存储器（memory）430和通信总线440，其中，处理器410，通信接口420，存储器430通过通信总线440完成相互间的通信。处理器410可以调用存储器430中的逻辑指令，以执行土地利用变化模拟方法。

[0101] 此外，上述的存储器430中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器（ROM，Read‑Only Memory）、随机存取存储器（RAM，Random Access Memory）、磁碟或者光盘等各种可以存储程序代码的介质。

[0102] 另一方面，本发明还提供一种计算机程序产品，所述计算机程序产品包括计算机程序，计算机程序可存储在非暂态计算机可读存储介质上，所述计算机程序被处理器执行时，计算机能够执行上述各方法所提供的土地利用变化模拟方法。

[0103] 又一方面，本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各方法提供的土地利用变化模拟方法。

[0104] 以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下，即可以理解并实施。

[0105] 通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到各实施方式可借助软件加必需的通用硬件平台的方式来实现，当然也可以通过硬件。基于这样的理解，上述技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来，该计算机软件产品可以存储在计算机可读存储介质中，如ROM/RAM、磁碟、光盘等，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行各个实施例或者实施例的某些部分所述的方法。

[0106] 最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。