CN110909924A - 基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法 - Google Patents

Abstract

本发明涉及基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，包括：1)：对卫星遥感影像进行监督分类获得土地利用分类图，并建立空间变量因子数据；2)：基于空间变量因子数据获取研究区域内的有效样本点；3)：建立CA城市扩张模拟原型模型，并基于有效样本点数据，对CA参数进行获取；4)：建立优化CA参数的相关目标函数，并利用交叉熵优化器优化CA参数；5)：建立CA转换规则，并获取转化概率图；6)：建立城市扩张模拟CA_CEO模型，并模拟预测城市扩张动态和未来可能情景；7)：对CA_CEO模型及其模拟预测结果，进行精度评定，并输出保存模拟结果。与现有技术相比，本发明有效优化CA模型并通过客观定权实现多目标城市扩张情景预测。

Description

基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法

技术领域

本发明涉及一种城市扩张多情景模拟方法，尤其是涉及一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法。

背景技术

城市发展是土地利用方式变化的结果，而土地利用方式变化又受到人类活动的强烈影响，并导致一系列社会和环境问题。当前，全球城市人口的快速增长使得对城市土地利用的需求不断增长，并且这种需求预计将持续数十年。通过构建模型来提高土地利用模拟的准确性，它可以作为未来区域土地利用规划和城市发展决策的重要参考价值。CA是一种时空动态模型，具有明显的时空耦合特征，特别适用于复杂土地利用变化研究的动态模拟。CA模型可以采用迭代自下而上的计算，结合邻域配置、约束效应和转移规则而实现。结合GIS的空间分析功能和CA技术可以有效模拟和预测城市空间扩张的时空过程。

转换规则是CA模型的核心，它决定了模拟的空间过程和逻辑关系。传统的CA模型主要使用统计学方法来定义转换规则。大多数这些转换规则由数学公式表示，并且公式中参数的确定非常困难。近年来，先后发展了获取CA模型转换规则的智能方法，包括遗传算法GA、蚁群算法ABC、人工神经网络ANN和粒子群优化PSO。

交叉熵优化CEO算法是近年出现的一种新型随机优化算法，它只依赖于适应度函数，不需要梯度等信息。CEO可用于衡量特征工程中两个随机变量之间的相似性。语言模型是通过训练集获得的，交叉熵优化是测量模型在测试集上的正确率；在逻辑回归中，CEO测量逻辑回归的预测结果与实际结果之间的差异程度，称为交叉熵损失函数。通过CEO优化器可以克服传统优化器容易陷入局部最优，或获取的参数没有明确的物理意义的弊端，同时能够实现多重条件下城市空间动态和土地利用变化多重情景的模拟预测。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，采用交叉熵优化元胞自动机方法能够实现logit参数的优化，同时能够为未来情景多重目标预测提供客观定权，实现城市扩张模拟与多情景预测。

本发明的目的可以通过以下技术方案来实现：

一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，包括以下步骤：

步骤1：对卫星遥感影像进行监督分类获得初始和结束年份的土地利用分类图，并建立影响城市扩张的空间变量因子数据；

步骤2：基于空间变量因子数据通过分层随机抽样方法获取研究区域内的有效样本点；

步骤3：建立CA城市扩张模拟原型模型，并基于有效样本点数据，利用logit方法对CA城市扩张模拟原型模型的CA参数进行获取；

步骤4：建立优化CA参数的相关目标函数，并利用交叉熵优化器在目标函数的引导下优化利用logit方法获取的CA参数；

步骤5：基于交叉熵优化器和优化完毕的CA参数，建立CA转换规则，并利用已建立的CA转换规则获取空间变量影响下土地的转化概率图；

步骤6：利用转换规则和转化概率图建立基于CEO的城市扩张模拟CA_CEO模型，并利用该模型模拟预测城市扩张动态和未来可能情景；

步骤7：对CA_CEO模型及其模拟预测结果，分别从规则拟合和模拟结果两方面进行精度评定，并输出保存模拟结果。

进一步地，所述的步骤1包括以下分步骤：

步骤11：利用两期陆地卫星遥感图像，基于马氏距离监督分类方法获取初始年份和结束年份的城市土地利用分类图；

步骤12：获取地形、经济和基础设施矢量和栅格数据，提取影响城市扩张的空间驱动因子数据；

步骤13：综合遥感图像数据、行政区划图、道路交通图等，在GIS环境中利用欧几里德距离方法获取到现有城市、区中心、主要道路、海岸线的距离数据。

进一步地，所述的步骤2具体包括：利用分层随机抽样方法对城市土地利用分类图和驱动因子图层进行抽样，获取研究区域内的有效样本点，为建立CA转换规则提供训练样本。

进一步地，所述的步骤3中的CA城市扩张模拟原型模型，其描述公式为：

式中，P_all表示土地转换概率，P_th表示决策规则的阈值，UrbanCATranRule表示UrbanCA的转移函数，C.s_n表示第n个单元状态，N.e表示领域效应，T.p表示驱动因子，OS.c表示总体和空间约束。

进一步地，所述的步骤3中的CA城市扩张模拟原型模型的综合的土地转换概率，其描述公式为：

式中，P_T.p表示由驱动因子决定的转移概率，P_n.e,t表示领域窗口中的土地转换率，PS表示概率缩放范围，NS表示领域缩放范围。

进一步地，所述的步骤4包括以下分步骤：

步骤41：建立用于交叉熵优化器优化的目标函数，以转换规则拟合的相对误差为基准，当相对误差目标函数越小，则其对应的CA参数则越优；

步骤42：利用交叉熵优化器在目标函数的引导下优化利用logit方法获取的CA参数。

进一步地，所述的步骤42中的交叉熵优化器中的信息的实际分布与拟合分布的交叉熵的计算公式为：

H(p,q)＝H(p)+D_KL(p||q)

式中，H(p,q)表示信息的实际分布p与拟合分布q的交叉熵，H(p)表示信息的实际分布p的熵，D_KL(p||q)使用Kullback-Leibler散度定义的相对熵。

进一步地，所述的步骤5包括以下分步骤：

步骤51：基于交叉熵优化器和优化完毕的CA参数，对获取的有效样本点数据以及空间变量因子数据训练CA转换规则；

步骤52：调整CEO的上下界和运行参数，得到不同条件下的CA参数组合，即多重CA参数集合；

步骤53：利用已建立的CA转换规则，在设定空间分辨率下获取空间变量影响下土地的转化概率图。

进一步地，所述的步骤6具体包括：利用转换规则和转化概率图建立基于CEO的城市扩张模拟CA_CEO模型，在GIS建模和模拟环境下，选用某年的城市土地利用格局为初始状态，利用模型运行M次，其中M表示初始与结束的年份差，得到城市土地利用变化的模拟及预测结果，并将该模型模拟和预测的土地利用变化结果输出。

进一步地，所述的步骤7包括以下分步骤：

步骤71：通过与遥感分类的土地利用格局进行比较，对CA_CEO模型模拟结果进行精度计算并评定，所述精度计算的指标包括：准确性Accuracy、莫兰指数Global Moran’s I、高/低聚类指数General G和景观指数；

步骤72：将CA_CEO模型模拟结果与遥感分类结果进行叠加并评定，所叠加的结果包括：初始城市Initial urban、实际和模拟均为城市Hit、实际为非城市模拟为城市False、实际为城市模拟为非城市Miss；

步骤73：在元胞自动机建模环境下输出并保存模拟结果。

与现有技术相比，本发明具有以下优点：

(1)本发明基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，采用交叉熵优化元胞自动机方法能够实现logit参数的优化，同时能够为未来情景多重目标预测提供客观定权，实现城市扩张模拟与多情景预测。

附图说明

图1为本发明的流程图；

图2为本发明实施例案例区域图；

图3为1995、2005和2015年的真实城市空间格局示意图；

图4为城市空间扩展的驱动因子图；

图5为CA_CEO模型在不同情景下的城市扩张概率图；

图6为CA_CEO模型在2005和2015年的模拟结果示意图；

图7为CA_CEO模型在2005和2015年的模拟评价示意图；

图8为CA_CEO模型在2025年的城市扩张情景预测示意图；

图9为CA_CEO模型在2035年的城市扩张情景预测示意图。

具体实施方式

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

本发明通过以下技术方案(如图1所示)来实现：

1)基于卫星遥感影像马氏距离监督分类，对初始年份和结束年份土地利用进行制图；基于矢量地图数据和城市土地利用变化图，建立影响城市扩张的空间变量因子数据；

2)通过分层随机抽样方法获取研究区域有效样本点，获取土地利用初始年份状态值、结束年份状态值、以及各因子样本数据；

3)首先建立CA城市扩张模拟原型模型；其次，基于上述样本点数据，利用logit方法对CA参数进行获取；

4)建立优化CA参数的目标函数-以相对误差作为目标函数，利用交叉熵优化器(CEO)在目标函数的引导下优化logit获取的CA参数，建立CA规则转换并产生空间变量影响下土地的转换概率图；

5)利用已建立的转换规则、结合邻域、限制、随机等模型元素，建立基于CEO的城市扩张模拟CA_CEO模型；利用已建立的CA_CEO模型模拟预测城市扩张动态和未来可能情景；

6)对CA_CEO模型及其模拟结果，分别从规则拟合和模拟结果两方面进行精度评定；

7)输出并保存模拟结果。

所述步骤1)具体为：

8)首先，利用两期陆地卫星遥感图像，基于马氏距离监督分类方法获取初始年份和结束年份的城市土地利用图；

9)其次，获取地形、经济和基础设施等矢量和栅格数据，提取影响城市扩张的空间驱动因子；

10)综合遥感图像数据、行政区划图、道路交通图等，在GIS环境中利用欧几里德距离方法获取到现有城市、区中心、主要道路、海岸线等关键要素的距离。

所述步骤2)具体为：

11)利用分层随机抽样方法对城市土地利用(城市空间格局)图、驱动因子图层进行抽样，要求至少5000个有效样本点，为建立CA规则转换提供训练样本；

所述步骤3)具体为：

12)CA原型模型包括四个基本部分：单元状态，邻域效应，由驱动因子产生的转换概率以及总体和空间约束。因此，可以通过以下方式定义UrbanCA中的CA模型：

式中，P_all表示土地转换概率，P_th表示决策规则的阈值，UrbanCATranRule表示UrbanCA的转移函数，C.s_n表示第n个单元状态，N.e表示领域效应，T.p表示驱动因子，OS.c表示总体和空间约束。

邻域是指以特定形状为边界的单元的组合，它确定下一次迭代的单元状态。对于一个单元，其邻域单元的状态在整个模拟过程中可能会连续变化，因此该单元邻域效应的概率也在变化。邻域效应通常由m×m个单元定义，并且邻域窗口中的土地转换率可以表达为：

其中Binaryization()是一个条件函数，如果单元格j的单元格状态C.s_j为城市，则返回1；否则，返回0。本发明应用5×5正方形邻域配置。

在本发明中，总体约束是可用于开发的全部城市单元，而空间约束是特定位置(例如水体)的城市发展约束。综合土地转换概率P_all可以定义为：

其中，PS是介于0.5-1.0之间的概率缩放范围，用来修改概率效应；而NS是介于0.0-0.1之间的邻域缩放范围，以修改邻域效应。

13)其次，定义由驱动因子决定的转移概率P_T.p，即单元格从非城市状态向城市状态的局部转移概率，代表驱动因素对土地利用变化的影响。假定y表示元胞状态是否转变，从时间t到t+1元胞状态从非城市(Non-urban)转变为城市(Urban)，y记为1；从时间t到t+1元胞状态未发生改变，则y记为0。当采用logit时可以表示为：

其中，a₀是常数，a_j是获取的CA转换规则的参数，x_j是自变量。

所述步骤4)具体为：

14)建立用于交叉熵优化器(CEO)优化的目标函数，以转换规则拟合的相对误差为基准，当相对误差目标函数越小，则其对应的CA参数则越优。

CEO的主要思想是将组合优化问题转化为关联随机优化问题,通过构造随机序列,使其以一定的概率收敛到最优或次优的结果。首先根据随机性或概率分布构造一个随机序列样本{X₁,X₂,…X_n},并计算每个样本对应的目标函数值S(1),S(2),…,S(n)。目标函数值的(1-α)分位点可以表达为：

当迭代次数t≥a且

时，可中止迭代，此时生成的是最优的随机序列样本；否则令t＝t+1，进行重新迭代。

15)CEO优化方法可通过下列过程实现。在信息论中，熵是由一组随机数据生成的平均信息含量。事件的概率分布和每个事件的信息内容构成一个随机变量，所有随机变量的均值(即期望值)是该分布所产生的信息内容(即熵)的平均值。熵是对随机变量不确定性的一种度量，是信息内容的期望值，离散分布的期望值可以通过以下公式得出：

其中，P(x_i)为离散随机变量X的概率分布，I(x_i)为离散随机变量X的信息量。CEO是深度学习中的常用概念，通常用于发现目标值和预测值之间的差异。给定集合上的分布p和q的CEO定义如下：

H(p,q)＝E_p[-log q]

其中，p是信息的实际分布，q是拟合分布。由于实分布的平均分布长度由拟合分布表示，因此计算p分布的平均值E_p[-log q]。

可以使用Kullback-Leibler散度来定义：

H(p,q)＝H(p)+D_KL(p||q)

式中，H(p,q)表示信息的实际分布p与拟合分布q的交叉熵，H(p)表示信息的实际分布p的熵，D_KL(p||q)使用Kullback-Leibler散度定义的相对熵。

根据吉布斯不等式，实分布的信息熵H(p)小于交叉熵H(p,q)，较小的值是相对熵(D_KL(p||q))。

16)利用步骤16)所述CEO方法和已建立的转换规则(步骤15)，对获取的有效采样点数据以及空间变量因子数据训练CA转换规则；

17)调整CEO的上下界和运行参数，得到不同条件下的CA参数组合，即多重CA参数集合；

18)利用已建立的CA转换规则，在30m空间分辨率下获取空间变量影响下的土地转换概率分布。

所述步骤5)具体为：

19)在GIS建模和模拟环境下，选用某年的城市土地利用格局为初始状态，利用模型运行M次(初始与结束的年份差)，得到城市土地利用模拟和预测结果。

20)将该模型模拟和预测的土地利用变化结果输出。

所述步骤6)具体为：

21)对模型模拟的城市土地空间线形态，分别从规则拟合和模拟结果两方面进行精度评定。评价方法如下：

通过与遥感分类的土地利用格局进行比较，对模拟结果进行精度计算，主要指标有准确性(Accuracy)、莫兰指数(Global Moran's I),高/低聚类指数(General G)，和景观指数等。将模拟结果与遥感分类结果进行叠加，结果包含：初始城市(Initial urban)、实际和模拟均为城市(Hit)、实际为非城市模拟为城市(False)、实际为城市模拟为非城市(Miss)。在评估准确性方面，应用总体精度(OA)和图形优度(FOM)两个指标，总体精度是使用最广泛的度量指标，而FOM注重评估元胞状态变化正确性的指标。使用Global Moran's I来测量空间自相关，而General G指数用于测量聚类程度。

所述步骤7)具体为：

22)在元胞自动机建模环境下输出并保存模拟结果。

本发明方法的实际实施例：

以1995-2015年间温州市历史城市空间扩张为案例，本案例区域位置如图2所示。为验证CA_CEO模型在城市空间扩张模拟上的有效性，模拟了温州市该时段的城市空间扩张过程；并利用CA_CEO模型面向四个不同的城市情景建立了四个不同的转换规则，预测并比较了温州市2025和2035年四种情景下的城市空间扩张。

基于交叉熵优化器的城市扩张多情景模拟元胞自动机优化方法，包括以下步骤：

1)首先，选择温州市1995、2005和2015年的陆地遥感图像，以及行政区划、道路交通、人口数据等，作为训练CA规则转换和获取土地转变概率的基础数据；

2)利用马氏距离监督分类方法对温州市遥感图像进行分类，产生温州市1995、2005和2015年的真实城市空间格局图(图3)；

3)基于各年遥感图像、行政区划图层、道路交通图层，在GIS环境下利用空间分析工具中的欧几里德距离计算出各元胞到已有城市、区中心、主干道路、海岸线等的距离，形成城市空间扩张的驱动因子图(图4)；

4)利用分层随机抽样方法对各空间变量的值、初始和结束年份的城市土地利用图进行采样，得到5691个有效样本点；

5)根据抽样得到的有效样本点，在R语言环境下利用CEO获取CA参数(表1)，从而产生城市土地转变概率图(图5)，进一步建立基于CEO的地理CA模型CA_CEO。图5为四种不同情景下的城市土地转换概率，表2为概率的统计信息。

6)以1995年真实城市状态为起点，将CA_CEO模型运行10次获取2005年模拟城市格局；以2005年真实城市状态为起点，将CA_CEO模型运行10次获取2015年模拟城市格局(图6)。对上述模拟结果进行图形化评估(图7)。

7)以2015年状态为初始状态，利用马尔可夫链方法对2025和2035年的城市像元进行预测，从而确定2025和2035年的城市土地总量。

8)利用已建立的CA_CEO模型预测2025和2035年城市情景，图8、图9分别显示2025年和2035年四种情景下城市扩张预测；表3和表4是对四种不同情景预测格局进行评估，结果表明本发明的方法能够为未来情景多重目标预测提供客观定权，实现城市扩张模拟与多情景预测。

表1.四种情景下的CA转换规则参数

表2.四种情景下的城市扩张转换概率统计

	最小值	最大值	值域	平均数	标准差
						情景1-BAU	0.0	1.0	1.0	0.07	0.11
情景2-区中心	0.0	1.0	1.0	0.01	0.04
						情景3-主要道路	0.0	1.0	1.0	0.03	0.07
情景4-海岸线	0.0	1.0	1.0	0.02	0.07

表3. 2025年城市扩展情景空间格局分析

表4. 2035年城市扩展情景空间格局分析

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到各种等效的修改或替换，这些修改或替换都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。

Claims (10)

1.一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，包括以下步骤：

步骤1：对卫星遥感影像进行监督分类获得初始和结束年份的土地利用分类图，并建立影响城市扩张的空间变量因子数据；

步骤2：基于空间变量因子数据通过分层随机抽样方法获取研究区域内的有效样本点；

步骤3：建立CA城市扩张模拟原型模型，并基于有效样本点数据，利用logit方法对CA城市扩张模拟原型模型的CA参数进行获取；

步骤4：建立优化CA参数的相关目标函数，并利用交叉熵优化器在目标函数的引导下优化利用logit方法获取的CA参数；

步骤5：基于交叉熵优化器和优化完毕的CA参数，建立CA转换规则，并利用已建立的CA转换规则获取空间变量影响下土地的转化概率图；

步骤6：利用转换规则和转化概率图建立基于CEO的城市扩张模拟CA_CEO模型，并利用该模型模拟预测城市扩张动态和未来可能情景；

步骤7：对CA_CEO模型及其模拟预测结果，分别从规则拟合和模拟结果两方面进行精度评定，并输出保存模拟结果。

2.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤1包括以下分步骤：

步骤11：利用两期陆地卫星遥感图像，基于马氏距离监督分类方法获取初始年份和结束年份的城市土地利用分类图；

步骤12：获取地形、经济和基础设施矢量和栅格数据，提取影响城市扩张的空间驱动因子数据；

步骤13：综合遥感图像数据、行政区划图、道路交通图等，在GIS环境中利用欧几里德距离方法获取到现有城市、区中心、主要道路、海岸线的距离数据。

3.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤2具体包括：利用分层随机抽样方法对城市土地利用分类图和驱动因子图层进行抽样，获取研究区域内的有效样本点，为建立CA转换规则提供训练样本。

4.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤3中的CA城市扩张模拟原型模型，其描述公式为：

式中，P_all表示土地转换概率，P_th表示决策规则的阈值，UrbanCATranRule表示UrbanCA的转移函数，C.s_n表示第n个单元状态，N.e表示领域效应，T.p表示驱动因子，OS.c表示总体和空间约束。

5.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤3中的CA城市扩张模拟原型模型的综合的土地转换概率，其描述公式为：

式中，P_T.p表示由驱动因子决定的转移概率，P_n.e，t表示领域窗口中的土地转换率，PS表示概率缩放范围，NS表示领域缩放范围。

6.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤4包括以下分步骤：

步骤41：建立用于交叉熵优化器优化的目标函数，以转换规则拟合的相对误差为基准，当相对误差目标函数越小，则其对应的CA参数则越优；

步骤42：利用交叉熵优化器在目标函数的引导下优化利用logit方法获取的CA参数。

7.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤42中的交叉熵优化器中的信息的实际分布与拟合分布的交叉熵的计算公式为：

H(p，q)＝H(p)+D_KL(p||q)

式中，H(p，q)表示信息的实际分布p与拟合分布q的交叉熵，H(p)表示信息的实际分布p的熵，D_KL(p||q)使用Kullback-Leibler散度定义的相对熵。

8.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤5包括以下分步骤：

步骤51：基于交叉熵优化器和优化完毕的CA参数，对获取的有效样本点数据以及空间变量因子数据训练CA转换规则；

步骤52：调整CEO的上下界和运行参数，得到不同条件下的CA参数组合，即多重CA参数集合；

步骤53：利用已建立的CA转换规则，在设定空间分辨率下获取空间变量影响下土地的转化概率图。

9.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤6具体包括：利用转换规则和转化概率图建立基于CEO的城市扩张模拟CA_CEO模型，在GIS建模和模拟环境下，选用某年的城市土地利用格局为初始状态，利用模型运行M次，其中M表示初始与结束的年份差，得到城市土地利用变化的模拟及预测结果，并将该模型模拟和预测的土地利用变化结果输出。

10.根据权利要求1所述的一种基于交叉熵优化器的城市扩张多情景模拟元胞自动机方法，其特征在于，所述的步骤7包括以下分步骤：

步骤71：通过与遥感分类的土地利用格局进行比较，对CA_CEO模型模拟结果进行精度计算并评定，所述精度计算的指标包括：准确性Accuracy、莫兰指数Global Moran’s I、高/低聚类指数General G和景观指数；

步骤72：将CA_CEO模型模拟结果与遥感分类结果进行叠加并评定，所叠加的结果包括：初始城市Initial urban、实际和模拟均为城市Hit、实际为非城市模拟为城市False、实际为城市模拟为非城市Miss；

步骤73：在元胞自动机建模环境下输出并保存模拟结果。

Images