本文考虑地块内部异质性，提出一个结合哈希函数和图注意力网络（GAT）的矢量元胞自动机（VCA）方法，用于研究城市土地利用变化；并将该模型应用于模拟深圳市2009年至2012年的城市土地利用变化，结果表明，HashGAT-VCA模型的模拟准确性显著优于其他VCA模型。
【论文题目】
HashGAT-VCA: A vector cellular automata model with hash function and graph attention network for urban land-use change simulation
【题目翻译】
HashGAT-VCA：一种结合哈希函数和图注意力网络的矢量元胞自动机模型，用于城市土地利用变化模拟
【期刊信息】
Landscape and Urban Planning；Volume 250, October 2024, 105145
【作者信息】
Qingfeng Guan, 中国地质大学地理与信息工程学院，湖北省武汉市430078，中国地质大学地理信息系统国家工程研究中心，湖北省武汉市430078，[email protected]
Jianfeng Li, 中国地质大学地理与信息工程学院，湖北省武汉市430078，中国地质大学地理信息系统国家工程研究中心，湖北省武汉市430078，[email protected]
Yaqian Zhai, 河南财经政法大学资源与环境学院，中国河南省郑州市金水东路180号，450016，[email protected]
Xun Liang, 中国地质大学地理与信息工程学院，湖北省武汉市430078，中国地质大学地理信息系统国家工程研究中心，湖北省武汉市430078，[email protected]
Yao Yao，中国地质大学地理与信息工程学院，湖北省武汉市430078，中国地质大学地理信息系统国家工程研究中心，湖北省武汉市430078，粤港澳智慧城市联合实验室，中国深圳，[email protected]
【论文链接】
https://doi.org/10.1016/j.landurbplan.2024.105145
【关键词】
矢量元胞自动机、图注意力网络、哈希函数、土地利用变化模拟、空间交互
【本文亮点】
-
提出了一个考虑地块内部异质性的框架。
-
使用哈希函数处理地块内部驱动因素的异质性。
-
使用图注意力网络（GAT）构建不同的邻域并获取交互关系。
-
不同的哈希函数和邻域对结果具有显著影响。
【摘要】
矢量元胞自动机（VCA）模型在表示不规则形状地块的时空动态方面表现出色，已广泛应用于土地利用变化模拟。然而，当前研究面临以下问题：(1) 大多数VCA模型在评估环境驱动效应时忽略了每个地块内部驱动因素的空间异质性；(2) 在计算邻域效应时，通常采用邻近地块的土地利用类型的简单统计，忽视了邻近地块内驱动因素的影响；(3) 对地块间交互作用的探索能力有限。为解决上述问题，本研究提出了一种用于研究城市土地利用变化的HashGAT-VCA模型。该模型利用哈希函数将每个不规则地块内驱动因素的非均匀分布编码为固定长度的向量，并基于地块间的空间拓扑关系构建图结构。通过采用图注意力网络（GAT），模型探索了环境驱动效应和地块间交互作用机制，以计算每个地块的土地利用变化概率。该HashGAT-VCA模型被应用于模拟中国深圳市2009年至2012年的城市土地利用变化。与其他VCA模型相比，HashGAT-VCA展示了更高的模拟精度。结果表明，HashGAT-VCA能够有效捕捉异质分布的驱动因素和地块间交互作用对土地利用变化的影响。此外，本研究还模拟了在生态控制策略下2025年和2030年的土地利用模式，为城市规划提供了决策支持。
【引言】
在中国，随着城市人口的增加，城市化率从1978年的17.92%迅速上升到2012年的52.57%，并在2022年达到了65.22%。这一持续的城市化进程具有重大的经济、社会和文化意义。然而，它也带来了许多挑战和问题（Deng et al., 2015）。城市化是一个高度复杂的社会经济过程（Deng et al., 2009）。随着城市化的加速，城市变化的模式变得越来越复杂（Tao and Ye, 2022）。实际上，有许多因素影响城市化，使得有效探索城市土地利用变化机制成为一个关键问题。
学者们已经开发了各种地理方法和模型来探索城市扩展机制并模拟城市土地利用变化（Li et al., 2017; Noszczyk, 2019）。其中，最广泛使用的是元胞自动机（CA）模型，它通过自下而上的方法来建模复杂的地理现象（Batty and Xie, 1994; Li et al., 2022; Liao et al., 2019）。
基于栅格的元胞自动机（RCA）模型是CA模型的一种典型形式，通常使用规则形状的单元格作为基本单元来研究复杂的地理动态（Chen et al., 2020; Li and Yeh, 2000; Liang et al., 2021a,b）。然而，由于城市实体的形状不规则，基于栅格的CA模型在进行细粒度模拟时面临挑战（Liao et al., 2016）。为了解决这一问题，一些研究者采用了基于斑块的CA模型，将多个栅格单元合并为斑块，用于城市扩展和土地利用变化模拟（Chen et al., 2014; Chen et al., 2016; Liang et al., 2021a,b; Meentemeyer et al., 2013）。由于土地分块通常由不规则形状的单元表示，使用栅格单元研究这些地块不利于准确捕捉地块的形状特征。
为了解决上述限制，已经开发了矢量元胞自动机（VCA）模型，该模型使用矢量多边形来表示土地分块实体。这种方法成功地解决了任意形状空间单元的细粒度划分问题，并在模拟城市扩展和土地利用变化方面显示出良好的潜力（O’Sullivan, 2001）。VCA将单元空间定义为不规则地理对象的集合，其中每个单元表示一个不规则地理实体。一个单元的邻域定义为一组影响目标单元状态变化的邻近单元（Moreno et al., 2008）。对于地块逐渐从一种类型过渡到另一种类型的情况，提出了DLPS-VCA模型（Yao et al., 2017）。DLPS-VCA使用动态地块分割方法，将土地分块划分为更小的基本单元，旨在更真实地模拟城市土地利用变化的渐变和碎片化过程。一些学者采用聚类技术对城市空间结构进行分类，并结合局部邻域效应和环境影响，探索了不同类型区域的土地利用过渡模式，构建了基于分区的VCA模型，取得了良好的结果（Yang et al., 2023）。
在VCA模型中，地块的土地利用转移概率通常由四个组件组成：环境驱动效应、邻域效应、约束和随机因素（Li and Yeh, 2002）。环境驱动效应指的是环境驱动因素（如地形特征、位置条件和社会经济条件）对地块内部土地利用转换的影响。邻域效应指的是周围地块对中心地块土地利用变化的影响。约束表示对特定地块某些土地利用转换的禁止。随机因素用于反映复杂城市系统发展过程中的随机性（Li et al., 2008）。
研究人员已经对上述组件进行了研究。研究评估并比较了不同邻域设置下的模拟精度（Stevens and Dragicevic, 2007）。为解决邻域的敏感性问题，研究人员开发了30种不同的矢量元胞邻域结构，基于类型和大小（Dahal and Chow, 2015）。Liang等（2018）基于各种类型单元的比例、邻域丰富度和景观模式指标计算了邻域效应。针对约束问题，Yu等（2011）将社会经济、人口和土地管理因素作为约束变量。Lauf等（2012）将家庭和住房因素作为约束变量。Yao等（2017）引入了一些自然元素作为约束。
环境驱动效应的计算依赖于探索驱动因素与土地利用类型/变化之间的复杂非线性关系。人工智能（AI）算法具有适应性、自学习和自组织等特性，使其特别适合处理非线性和复杂问题。许多学者在CA模型中应用了各种AI算法，以揭示驱动因素对土地利用变化的影响机制。典型的例子包括随机森林（RF）（Kamusoko and Gamba, 2015; Liang et al., 2021a,b; Lv et al., 2021）、支持向量机（SVM）（Feng et al., 2016; Ren et al., 2019; Yang et al., 2008）、人工神经网络（ANN）（Li and Yeh, 2002; Pijanowski et al., 2002; Xu et al., 2019）、卷积神经网络（CNN）（Xing et al., 2020; Zhai et al., 2020）以及递归神经网络（RNN）。
然而，传统的VCA模型大多假设驱动因素在每个土地分块内均匀分布，使用平均值来表示土地分块特征。由于模型中使用的驱动因素数据基于不均匀分布的离散单元，使用平均值忽略了每个土地分块内部驱动因素的空间异质性。CNN-VCA（Zhai et al., 2020）通过利用空间连续的驱动因素栅格数据来解决这一限制。该模型使用卷积神经网络（CNN）提取每个单元中心周围一定面积内的驱动因素的空间分布特征及其与多种土地利用类型适宜性的关系，从而提高模拟精度。然而，CNN的感受野受限于卷积核的规则形状和大小，使其在提取VCA中不规则形状和大小单元内驱动因素的分布特征时效果不佳。
在邻域效应的计算中，大多数VCA模型采用邻近单元中土地利用类型的简单统计数据，如不同土地利用类型的数量和多样性（Guan et al., 2023; Kamusoko and Gamba, 2015），来估计邻近单元对中心单元土地利用变化的影响。然而，这种简单统计方法在全面捕捉细胞间相互作用机制方面存在局限，并且忽略了邻近单元内环境驱动因素对中心单元土地利用变化的影响。近年来，图神经网络（GNNs）因其能够捕捉多个实体（即图中的节点）之间的交互信息而被广泛应用于交互建模。
与传统机器学习方法相比，GNN在消息传递过程中引入了邻接矩阵，这描述了实体之间的空间关系（Gori et al., 2005; Scarselli et al., 2009）。在GNN中，每个节点根据其邻近节点的信息进行更新，以获得新的隐藏表示。消息传递过程通常包括两个步骤：消息聚合和消息更新。然而，早期的GNN模型存在计算效率低、迭代过程中参数共享以及学习边缘隐藏状态的困难等问题。为了解决这些挑战，图卷积网络（GCNs）（Kipf and Welling, 2016）将卷积操作扩展到图结构数据中，并取得了显著成功。为了为每个邻近节点分配不同的注意力分数，并识别更重要的邻居，图注意力网络（GAT）（Velikovi et al., 2017）在消息传递过程中引入了注意力机制，其中节点的隐藏状态通过计算其邻居的注意力分数来进行计算。
在HGAT-VCA（Guan et al., 2023）中，基于土地分块的邻接拓扑关系构建了图模型。该模型利用图注意力网络（GAT）揭示了土地分块内部的环境驱动机制以及邻近分块之间的相互作用机制。然而，与大多数先前的VCA模型类似，HGAT-VCA忽略了每个土地分块内部驱动因素的异质分布模式，假设每个分块内的驱动因素均匀分布。
基于上述现有研究的回顾，可以总结出当前研究中的几个问题：（1）大多数VCA模型在评估环境驱动效应时忽略了每个土地分块内部驱动因素的空间异质性；（2）在计算邻域效应时，常常使用邻近分块中的土地利用类型的简单统计数据，忽略了邻近分块内驱动因素对目标单元的影响；（3）探索土地分块之间相互作用的能力通常有限。
因此，需要一种方法，能够充分探索并利用驱动因素在土地分块内的不均匀分布对土地利用变化的影响，考虑其不同的形状和大小。此外，这种方法还应能更好地捕捉土地分块之间的相互作用机制，不仅包括邻近土地分块的土地利用类型对目标单元的影响，还包括它们的驱动因素对目标单元土地利用变化的影响，以提高VCA的模拟和预测能力。
【方法】
这项研究提出了一种 HashGAT-VCA 模型，用于模拟城市区域中不规则形状土地单元的多类别土地利用变化。如图 1 所示，该模型由四个主要模块组成：（1）信息提取和编码：提取每个土地单元中每个驱动因素的栅格数据。然后应用 SHA512 算法，将驱动因素的栅格数据编码为等长度向量，这些向量代表每个土地单元的驱动因素的空间分布特征。（2）图构建和训练：构建土地单元图。每个土地单元被视为一个节点，将编码后的等长度驱动因素向量与其他特征（如土地利用类型）结合，形成土地单元的特征。根据土地单元之间的邻接关系构建边，形成土地单元的图网络。通过训练图注意力网络（GAT），同时探索土地单元内部的环境驱动效应和土地单元之间的邻域效应。（3）VCA模拟：利用训练好的 GAT 计算每个土地单元在环境因素和土地单元间交互作用的综合影响下的土地利用转换适宜性。通过结合约束条件和随机因素，得到每个土地单元的土地利用转换概率。使用历史统计数据进行参数调整，并经过多次迭代，通过轮盘选择过程确定每个土地单元的最终土地利用类型。（4）准确性评估和未来情景预测：训练好的模型可以模拟历史土地利用变化，从而评估其准确性和可靠性。该模型可用于生成不同情景下的未来土地利用空间模式，为政府部门的决策和规划提供支持。
图1. 基于提出的HashGAT-VCA 模型的城市土地利用变化模拟流程图。
2.1. 信息提取和编码
如图2(a)所示，在原始数据集中，地块数据由描述地块位置和形状的矢量多边形组成。特定环境驱动因子的数据显示为具有规则形状和大小的空间单元（例如30米× 30米的像素）的连续栅格表面，表示研究区域内驱动因子值的空间分布（Zhai et al., 2020）。
每个驱动因子的栅格数据首先通过地块的矢量多边形进行裁剪，生成所有地块的一组栅格补丁，每个补丁表示特定地块内驱动因子的空间分布。当栅格单元的质心位于矢量地块的边界内时，它被认为代表地块的一部分（如图2(b)所示）。所有驱动因子的栅格补丁具有不规则形状和不一致的像素数量，与其对应的地块一致。接下来，为每个栅格补丁创建一个最小包围矩形（MBR），并进行填充处理，将MBR中补丁外的像素填充为特定的NoData值（例如，−1或−9999）。最终得到的MBR保留了地块的形状和大小信息，以及地块内驱动因子的空间分布（如图2(c)所示）。
如前所述，图神经网络模型要求所有节点的属性特征维度一致。因此，在HashGAT-VCA模型中，使用Hash函数将各种大小的驱动因子MBR编码为统一长度的向量（如图2(d)所示）。一旦构建了地块图，每个驱动因子向量将作为图中相应节点（即地块）的属性特征的一部分（如图2(e)所示）。其他属性（如地块的土地利用类型）也将添加到相应节点的属性特征中。
图2. 提取和处理过程
哈希函数是广泛使用的加密函数，可以处理几乎任何长度的信息，并输出固定长度的哈希向量（Brickell et al., 2000）。这种特性适合于处理大小不同的驱动因素数据块。哈希函数的主要安全要求是其抗碰撞性，因为减少碰撞可以提高信息编码的效率（Kelarev, 2006）。典型的哈希函数包括 MD4、MD5、SHA-0、SHA-1和 SHA-2。其中，SHA-2具有四种压缩算法：SHA-224、SHA-256、SHA-384和 SHA-512，每种算法具有不同的哈希长度（Dinur et al., 2012; Sklavos and Koufopavlou, 2003）。
HashGAT-VCA 模型使用 SHA-512 算法，该算法包括四个主要步骤：填充、分块、初始化和迭代加密（Aoki et al., 2009; Lien et al., 2004）。其目的是以最小的信息丢失对每个土地分 parcel 内的MBR 进行编码，并生成一个具有统一长度的特征向量，该特征向量将作为图中相应节点特征的一部分。与其他哈希算法相比，SHA-512 具有更长的长度，可以保留更多信息。一些词嵌入方法也可以将不同长度的向量映射为统一长度的向量，例如 Word2Vec 和Doc2Vec，但这些方法需要上下文信息并用于处理离散数据，因此不符合我们模型的要求。SHA-512 的详细过程如下：
（1）填充: SHA-512 的信息长度为 1024 位。假设初始信息的长度为 H，则所需的填充长度为 h，应满足公式（1）。填充信息以 1 开头，后面跟着 h-1 个零。原始信息长度计算后作为 128 位的信息字段，随后附加到末尾。
（2）分块:信息被分为多个块，每个块的长度为 1024 位。每个1024 位块进一步分为 16 个子块，每个子块包含64 位。
（3）初始化: 需要初始化 8 个变量，记作 a-h，以及 80 个常量。8 个链接变量的值通过取前 8 个质数平方根的小数部分，并保留前 64 位来获得。80个常量的值通过取前 80 个质数立方根的小数部分，并保留前64 位来获得。
（4）迭代加密: 迭代加密过程包括 80 轮。迭代加密的算法步骤如图 3 所示。相关的计算方法将在以下文本中解释。在图3 中，a-h 代表初始化变量。Maj、Ch、BISG0和 BISG1 的计算公式由方程（2-5）给出。W [t] 代表从分块过程获得的信息块，其中前16 个词对应于原始信息子块，其余 64 个词使用公式（6-8）计算，其中 K [t] 代表常量。通过按照图2 所示的过程进行多次迭代，得到 a-h 的最终值，总计512 位。为了降低计算成本，这些值可以转换为一个 128 位的十六进制向量，作为网络中节点的特征向量。
经过80 次迭代，驱动因素的栅格 MBR（最小边界矩形）被处理为相同长度。这个处理后的向量被整合到图结构中的节点特征中，为后续工作提供数据基础。除了驱动因素向量，每个地块的其他特征也被整合到对应节点的特征中，包括地块的周长和面积、MBR 的长度和宽度，以及原始地块土地使用类型的独热编码（one-hot encoding）。
图 3. SHA-512 算法的迭代加密过程。
2.2. 图构建与训练
如图1(b)所示，图注意力网络（GAT）需要两个输入：图结构和节点特征，通常表示为G = (V, E, X, A) 。其中，V = {v1, v2, ..., vn}表示节点集合，其中n是节点的数量。E表示边的集合，A是大小为n*n的邻接矩阵。假设d是节点特征向量的长度，X是大小为n*d的特征矩阵。为了使用GAT探索环境驱动效应和邻域交互效应，过程可以分为两个主要部分：图构建和GAT训练。
2.2.1. 图构建
图中的地块图构建如下：每个节点表示一个地块，如果两个地块被认为是邻近的（或相连的），则在它们之间连接一条边。图的邻接矩阵 \( A \) 使用公式 (9) 构建，其中如果两个地块被认为是邻近的，则邻接矩阵的值为 1，否则为 0。尽管这可能忽略了地块之间连接的强度，但 GAT 通过注意力机制学习地块之间交互的强度。
在VCA 模型中，某一地块的邻域效应在很大程度上取决于邻近地块的选择方法，而不同的邻近选择方法会导致不同的邻域配置，从而产生不同的邻域效应。通常有两种邻近选择方法：基于缓冲区的方法和基于拓扑的方法（Dahal 和Chow，2015）。基于缓冲区的方法可以进一步分为两种：边缘基准和质心基准（Barreira-González等，2015）。基于缓冲区的方法涉及在感兴趣地块的质心或边缘周围创建不同大小的缓冲区，如图 4(a)(b) 所示。与缓冲区相交的周围地块被视为与中心地块连接。由于地块形状复杂，地块的质心可能不在地块内部。此外，当地块较大时，缓冲区可能被包含在地块内部，如图 4(c)(d) 所示。因此，基于边缘的建模方法在理论上会产生更好的结果。然而，无论使用哪种方法，都需要通过反复实验来确定适当的缓冲距离，以提高模拟精度。
基于拓扑的方法则将地块视为分析单位。与中心地块连接的地块称为一阶地块，与一阶地块连接的地块称为二阶地块，如图 4(e) 所示。在图论中，通过对邻接矩阵进行相关操作来实现中心地块与高阶地块之间的连接。高阶邻接矩阵 \( \tilde{A}_K \) 的计算使用公式 (10)，其中 \( K \) 为阶数。
图4展示了图结构的构建过程。
2.2.2. 图注意力网络（GAT）训练
在构建土地分 parcel 图之后，使用注意力机制来学习连接节点之间的注意力系数，这些系数表示土地分 parcel 之间的互动强度。中心土地分parcel 的土地使用变化适宜性受其自身驱动因素以及周围土地分 parcel 的影响。
HashGAT-VCA 中使用的 GAT 具有两层。图5展示了具有多头注意力机制的 GAT 结构。其中，h1̅→(0)代表中心节点，其他节点代表邻近节点，0 代表节点的初始输入特征。在连接中心和周围 parcel 的特征后，使用单层前馈神经网络作为可学习的向量执行注意力操作，通过 LeakyReLU 激活函数（Velikovi 等，2017）获得节点之间的注意力值。在对所有节点执行上述操作后，使用 Softmax 得到图的注意力矩阵，并计算第一层的输出结果。此时得到的节点特征由第一层网络的计算输出表示为 h1̅→(1)。随后，通过第二层网络计算第二层的输出，表示多个土地使用过渡的适宜性值，从而得到结果 h1̅→(2)。最终结果（即过渡类型）通过 Softmax 分类器获得。这种GAT 使用历史土地使用变化数据进行训练，以获得图中连接节点之间的注意力系数。
如果周围 parcel j 对中心 parcel i 的影响显著，则对应的注意力系数 eij 将较高，如公式（11）所示。
为了稳定注意力的学习过程，GAT 引入了多头注意力机制，这类似于（Vaswani等，2017）的工作。在本研究中，多头注意力由公式（14）表示。
其中，“||”表示拼接，αkij 代表第k 个头的归一化注意力系数。Wk 是对应输入线性变换的权重矩阵。输出特征hi→ʹ代表多头注意力机制的结果。
图 5. GAT 的实现细节。
2.3. 基于 GAT 的 VCA 模拟
在许多CA（细胞自动机）模型中，土地分 parcel 最终的转化概率由四个组成部分构成：环境驱动效应Pg、邻域效应Ω、约束条件 Pc和随机因素RA。因此，土地 Pi在时间t发生第k种土地利用变化的概率Pik,t可以表示为公式 (15)。在本研究中，使用图网络方法来汇聚周围土地 parcel 的信息，使得中央土地 parcel 能够综合考虑环境驱动效应和邻域效应，这些效应被表示为转化适宜性Ph。因此，在HashGAT-VCA 中，最终的转化概率 \(OP\) 可以通过公式 (16) 给出。
在公式(16) 中，Phik,t代表土地 Pi在时间t发生第k种土地利用变化的适宜性，这表示考虑了环境因素和邻域 parcel 影响的土地利用变化概率。在本研究中，GAT模块用于提取每个土地 parcel 在各种类型土地利用变化下的适宜性。
Pcik，t指的是土地Pi是否允许被开发为城市用地（即包括公共服务、商业、住宅或工业用地），进行第k种土地利用变化的情况。在本研究中，三种土地利用类型（特殊土地、水体和保护区）被指定为限制开发区域。这意味着，如果土地Pi属于限制区域，其约束系数Pck,t设置为0，表示不允许开发。相反，如果土地Pi不属于限制区域，其约束系数Pck,t设置为1，表示允许开发。
RA代表土地Pi在时间t发生第k种土地利用变化的随机因素。它反映了城市区域土地利用动态中的随机效应。在本研究中，RA的值通过一个随机函数计算，范围从 0 到1。具体计算过程为 A = 1 + (−ln y)α，其中α取值范围在 1 到 10 之间，y是一个介于 0 和 1 之间的随机值。
在获得OP后，使用轮盘选择机制来确定最终的土地 parcel 转换（Liang et al., 2021a,b; Liu et al., 2017）。对于每个土地 parcel i 在每个时间步骤t上，假设一个随机数r∈[0,∑k=05 OPik,t)。如果∑k=05 OPik,t的值等于OPi0,t的值，表示该 parcel 将不会发生变化。然而，如果它们不相等，则通过以下公式 (17) 确定最终的土地parcel 转换类型。最终，获得了土地 parcel i在时间t上各种土地利用变化类型的概率，并使用轮盘选择确定最终的转换 parcel。该模型以时间步长表示一年，进行迭代模拟过程。
在模拟过程中，可以根据研究区域的历史土地利用数据设置和校准相关参数，如土地利用变化的总面积、每种土地利用变化类型的最大土地 parcel 面积以及每种土地利用变化的土地 parcel 总数(Yao et al., 2017)。
2.4. 准确性评估
为了定量评估土地利用变化模拟结果的准确性，使用了优值（Figure of Merit, FoM）、生产者准确率（Producer’s Accuracy, PA）和用户准确率（User’s Accuracy, UA）。这些指标主要关注模拟过程中变化的数量（Pontius et al., 2008）。具体计算公式如下（方程18–20）：
A：代表模拟结果保持不变时所导致的错误，而实际上存在真实的土地利用变化。 B：代表模拟出的土地利用变化正确地与现实中的真实变化相匹配。 C：代表模拟结果捕捉到了土地利用变化，但模拟的变化类型与实际变化不匹配。 D：代表模拟结果显示了土地利用变化，而实际上没有真实的变化。
此外，使用一组景观指数（Landscape Indices, LIs）来评估实际与模拟土地利用之间的景观模式相似性（Chen et al., 2014; Hagen-Zanker and Lajoie, 2008）。通过Fragstats 4.2计算了若干个LIs（包括形状指数的范围数（SHAPE）、周长-面积比的范围数（PARA）、相关外接圆的范围数（CIRCLE）、邻接指数的范围数（CONTIG））。相似性通过上述LIs的平均差异来估算。公式如下：
其中Ii,s和Ii,o分别表示第i个景观指数在模拟和实际土地利用中的值，δi表示第i个景观指数的差异。 在公式（22）中，αI是实际和模拟土地利用之间的景观模式相似性，n代表景观指数的数量。
【研究区域和数据】
深圳市位于中国广东省南部珠江口的东海岸，总面积1997.47平方公里，人口约为1252.8万（https://www.sz.gov.cn/cn/zjsz/gl/）。截至2021年，深圳的GDP达到了30664.85亿元人民币，位列中国第三（https://tjj.sz.gov.cn/）。深圳被公认为中国的经济、金融和科技创新中心，拥有文化、创意、高科技、现代物流和金融四大支柱产业。这些产业使深圳成为一个著名的创新和创业中心。
深圳市下辖10个行政区，其中福田区和罗湖区是最早开发和经济最发达的中心城区，是深圳的行政、金融、商业和外交中心（Yiyong and Zacharias, 2015）。南山区是深圳的科技研究、教育和体育中心。大鹏新区有大量自然风景区，但经济发展较慢。宝安区人口众多，工业发展历史较长（http://tjj.sz.gov.cn/）。行政区的分布见图6。
图 6. 研究区域和生态控制区域。
图 7. 深圳市2009年、2012年和2014年的城市土地利用数据。
表 1 2009年至2014年不同土地利用变化类型的总面积（单位：km²）。
为了保护和维护城市生态环境，确保可持续城市发展，防止城市建设用地的无序扩张，深圳在中国生态保护红线的划定和管理方面走在前列。该措施的实施旨在强调城市生态环境的保护和重要性，限制城市建设用地的过度扩张，确保生态系统的健康和稳定。图6的左侧显示了研究区域的情况，右侧显示了生态控制区域的情况。
本研究使用了深圳市土地资源局提供的2009年、2012年和2014年的土地利用数据。土地利用数据分为五类：非城市用地、公共服务、商业、住宅和工业，如图7所示。左侧列展示了2009年至2014年深圳的土地利用情况，右侧两列的红框代表了土地利用的变化。经过统计分析，2009年地块的平均周长为530.31米，中位数为292.74米；平均面积为17178.38平方米，中位数为3062.83平方米。
根据土地利用数据统计，2009年、2012年和2014年地块的数量分别为104,608、114,541和123,325个。其中，2009年非城市用地的比例为70.52%，2012年为69.24%，2014年为68.73%。这表明在2009年至2014年期间，部分非城市用地已转变为城市用地。
为了消除道路对构建地图的影响，相关的道路地块被移除。这并不意味着在模拟中没有考虑交通的影响；只是这些道路地块在图的构建和模拟中不予考虑。这是因为在构建图时，一些地块可能通过道路与远处的地块连接。在后续引入驱动因素时，交通对城市土地利用变化的影响将被考虑。删除道路后的2009年地块数量为103,844。移除道路后的土地转换具体情况见表1，纵轴表示2009年，横轴表示2014年。最后一列表示某种类型地块的总变化面积，最后一行表示已转换为某种类型地块的总面积。
从上述表格中可以看出，2009年至2014年间，深圳市土地面积中有2.73%经历了土地利用变化。这些变化大多数涉及非城市用地转变为四种城市用地类型（公共服务、商业、住宅和工业），而城市用地类型转换为非城市用地的情况较少。城市用地类型之间几乎没有内部转换。因此，本研究仅关注非城市用地转变为城市用地的情况，而未考虑其他类型的转换。非城市用地转变为不同类型城市用地的情况见图8。
城市土地利用变化通常受位置、地形、社会经济活动、几何特征和政策等因素的影响（Santé et al., 2010; Wu et al., 2021; Yao et al., 2017）。然而，由于缺乏详细的人口和经济数据，并且难以将政策等因素直接纳入模型，本研究目前未考虑这些因素。本研究中的位置和交通因素包括距离城市中心、距离高速公路、距离铁路和距离道路。地形因素包括数字高程模型（DEM）和坡度。POI因素包括医院、公交车站、餐馆、娱乐场所、公园、超市、购物中心和工厂。几何特征包括地块的周长和面积。图9显示了所有驱动因素的30米分辨率栅格数据集。POI的密度可以视为人类经济活动的反映。
图8. 2009-2014年非城市土地转变为不同类型城市土地的情况
图9. 辅助地理空间数据集：(a)坡度，(b)数字高程模型(DEM)，(c)公园密度，(d)公交车站密度，(e)超市密度，(f)购物中心密度，(g)工厂密度，(h)餐馆密度，(i)医院密度，(j)距高速公路的距离，(k)距铁路的距离，(l)距道路的距离，(m)距区中心的距离，(n)娱乐场所密度。
【结果】
所提出的HashGAT-VCA模型由多个模块实现。SHA512算法使用Python及多个开源库（如PyTorch、NumPy、SciPy和Hashlib）实现。VCA模型使用C++语言以及一些开源库（包括GDAL和QtCore）实现。实验在配备有Intel Core i9-10900X CPU @ 3.70 GHz、32 GB内存和NVIDIA RTX 4090 GPU的计算机上进行。
4.1. 模型应用及结果
使用2009年的土地使用数据和驱动因素数据集，构建并训练了HashGAT-VCA模型。利用校准后的VCA模型，研究模拟了2009至2012年的城市土地使用变化，并评估了模拟结果的准确性。此外，还使用了四种其他的VCA模型进行比较，即CNN-VCA（Zhai et al., 2020）、RF-VCA、ANN-VCA和LR-VCA。研究还讨论了不同哈希函数对结果的影响，并展示了SHA-512的有效性。不同的图构建方法被用来建立图网络，并比较了三种模型构建方法的准确性，同时对GCN和GAT之间的准确性进行了比较。
在实验中，为避免小土地单元的尺度效应，为每个土地单元选择了一个45米的缓冲区，以提取驱动因素数据并构建土地单元的驱动因素数据集。这确保了即使是较小的土地单元也能够构建驱动因素的数据集，否则可能导致一些小土地单元无法有效表达其自身的驱动因素。在GAT的训练过程中，随机选择了2009年80%的样本数据作为训练数据，剩余的20%用于验证数据。经过多次迭代，选择了训练性能最佳且损失最小的GAT模型。图10展示了五种模型在2012年的实际土地使用情况和模拟结果。
图10. 2012年深圳实际与模拟的土地使用情况。(a) 实际土地使用变化，(b) HashGAT-VCA，(c) CNN-VCA，(d) RF-VCA，(e) ANN-VCA，(f) LR-VCA。
4.2. 性能细节
准确性评估见表 2，HashGAT-VCA 模型展示了最佳的模拟结果。与 CNN-VCA 相比，HashGAT-VCA改进了 5.81%。在三年期间，只有少量的地块发生了土地使用变化，这使得揭示土地使用变化的潜在驱动机制变得困难。实验结果表明，HashGAT-VCA 能够有效探索土地使用变化的驱动机制。在实验的五个模型中，HashGAT-VCA 模型表现最好，说明使用 HashGAT 模型可以更好地提取环境驱动效应和邻里效应。将地块表示为图形可以表达空间信息，这有助于分类和预测任务。
表2 不同模型模拟结果的整体FoM。
此外，为了验证模型的适用性，HashGAT-VCA 模型被应用于模拟和评估2014 年整个区域的地块变化。结果显示，FoM 为0.287，UA 和 PA 值分别为 0.513 和0.333。与 CNN-VCA（FoM = 0.286）和 RF-VCA（FoM = 0.254）的结果相比略有改进。这表明城市土地变化是一个高度复杂的过程，具有随机性和不确定性（Yang et al., 2023; Ye et al., 2015）。尽管 2014 年的结果略逊于 2012 年，FoM下降了 31.41%，但仍表现出良好的性能，表明HashGAT-VCA 模型能够有效模拟地块的土地使用变化。
本研究结合了四个景观指数，包括 SHAPE（形状指数）的范围，PARA（周长面积比）的范围，CIRCLE（相关外接圆）的范围和CONTIG（连通性指数）的范围，以评估模拟和实际土地使用的景观相似性。表 3 显示了2012 年五个模型的景观相似性计算结果。与 CNN-VCA、LR-VCA、ANN-VCA和 RF-VCA 相比，HashGAT-VCA模型与实际土地使用的景观相似性最高（0.999），大约比其他四个VCA 模型高 0.5% 到2.5%。
表3 基于不同模型的模拟结果的景观指数。
为了进一步比较不同模型的模拟土地使用变化，图 11 显示了2012 年五个模型的实际土地使用变化和模拟土地使用变化。此外，图 12 选择了几个代表性区域以说明实际变化与模拟变化之间的差异。这些区域中的主要地块发生变化的是大地块和一些簇。图 13 展示了2014 年的模拟变化结果及一些细节。这些细节揭示了 HashGAT-VCA 和CNN-VCA 模型在模拟过程中取得了良好的结果，而其他三个模型则难以逼近真实情况。FoM 值也证实了这一点（FoM = 0.382）。与CNN-VCA 相比，本研究提出的 HashGAT-VCA 模型的模拟结果更接近地块的实际转变模式。HashGAT-VCA 能够有效识别大地块和一些簇的转变。这是因为，与其他模型相比，HashGAT-VCA 模型中的图网络允许提取环境驱动的适宜性，同时聚合邻里效应，从而提高了准确性。
图11. 深圳市2009年至2012年的实际和模拟土地使用变化。(a) 实际土地使用变化，(b) HashGAT-VCA，(c) CNN-VCA，(d) RF-VCA，(e) ANN-VCA，(f) LR-VCA。
图12. 2009年至2012年的实际和模拟土地使用变化细节。(1) 龙岗区；(2) 南山区；(3)坪山区；(4) 盐田区。
图13. 显示了2014年的模拟结果和细节。
从2009 年到 2012 年，一部分地块经历了非城市土地转变为城市土地的过程。通过表 1 的统计分析发现，在四种转变类型中，转变为工业用地的比例最大，而转变为商业用地的比例最小。
对每种土地转变类型进行分析后发现，由于工业用地的特点，非城市到工业用地的转变集中在一些大区域，主要位于城市边缘。这是由于政策的影响，工业用地通常以工业集群的形式存在，并由于潜在的污染问题而位于城市的外围。另一方面，转变为公共设施和住宅用地的区域相对分散，没有集中区域，并发生在较小的局部区域。这是由于人口密度和经济发展水平等因素可能在城市中分布不均，导致住宅和公共服务用地的分布相对分散。商业用地的转变较少且高度分散。商业用地通常受到市场需求和竞争的强烈影响，商业活动的发展往往更灵活和多样化，因此商业用地的分布可能更加分散。
模拟结果显示，在模拟商业用地（FoM = 0.487）、住宅用地（FoM = 0.389）和工业用地（FoM = 0.386）方面表现良好，而公共服务用地的表现较差（FoM = 0.285）。这是因为商业和住宅用地通常在附近配有相应的支持设施，而工业用地通常以工业区的形式存在。因此，使用基于图的方法可以利用周围的环境信息。公共服务用地的发展通常受政府规划的影响，导致准确性较低。
FoM 的整体表现较好，高 PA 值表明准确预测发生实际变化的地块的能力。然而，较低的 UA 值表明HashGAT-VCA 的误差主要来源于实际没有发生变化但被模拟的情况。这种对变化的敏感性是这些误差的一个因素。
根据表4，深圳不同行政区的模拟准确性情况如下。光明区、福田区、龙岗区、盐田区和坪山区的 FoM 高于深圳的平均水平。另一方面，罗湖区、大鹏区和龙华区的模拟表现较差。这可以归因于这三个区在 2009 年到2012 年期间是新兴开发区，受到更多政策干预，增加了模拟的难度。
表4 不同区县的HashGAT-VCA模型模拟结果的FoM值。
4.3. 参数敏感性分析
编码函数的选择直接影响土地地块特征向量的质量，并对最终的模拟结果产生一定影响。为了探索转换和邻里效应的适用性，使用了能够学习地块之间直接注意力的图注意力网络。因此，需要调查模型组件选择对结果的影响。此外，不同的图构建方法可以用来研究不同邻里影响对土地地块的效果，因此还需要讨论在第2.2.1节中提到的三种图构建方法。
4.3.1. 哈希函数和GAT的影响
首先，为了评估基于哈希的编码效果，使用土地块内部驱动因素的平均值作为相应节点的特征，构建了比较模型（即 Average-GAT-VCA）。内部驱动因素的平均值在大多数VCA方法中也被作为常见特征，例如HGAT-VCA和RF-VCA。为了评估哈希函数选择的影响，使用了三种哈希函数（即SHA-512、SHA-256和MD5）来构建比较模型。此外，还使用了图卷积网络（GCN）来构建模型进行比较，以评估GAT的有效性。表5显示了所有这些模型的准确性。表5 不同编码算法和图神经网络模型的模拟准确率。
结果表明，使用驱动因素的平均值作为节点（即土地块）的特征会导致土地块驱动因素的异质性丧失，从而表现较差。为了进一步调查这些差异，图14展示了使用哈希函数与使用平均值作为节点特征之间注意力系数的差异。图中的中央土地块标记为1。注意力系数表示邻近土地块对中央土地块变化的影响强度。图(c)和(d)中的线条厚度和颜色强度表示注意力系数的强度，指示周围土地块对中央土地块转化的影响。在所有这些情况下，基于哈希的模型成功预测了中央土地块的土地利用变化，而基于平均值的模型则未能做到。这些案例表明，使用驱动因素的平均值作为节点特征忽略了每个土地块内驱动因素的空间异质性，从而导致注意力焦点的变化，导致准确性较低。
图 14. 不同注意力系数的比较。(a) 2009 年土地利用类型。(b) 2012 年土地利用类型。(c) 使用哈希函数的模型的图注意力系数。(d) 使用驱动因素平均值作为节点特征的模型的图注意力系数。
SHA-512由于其较长的编码长度，相比SHA-256和MD5保留了更多的驱动因素信息。由于MD5具有最短的长度（128位），其模拟性能在这三种哈希函数中最差。此外，实验还表明，具有注意力机制的GAT能够更好地探索土地块之间的影响效应，并且相比GCN实现了更高的模拟准确性。
4.3.2. 邻域选择的影响
为了找到合适的邻域选择方法，对图 4 中的三种不同图构建方法进行了比较。在高阶邻域选择方法中，进行了四次实验，最高阶数设置为 4。对于不同的邻域缓冲区，进行了八次实验，缓冲距离范围从 50 米到200 米，间隔为 50 米。实验的一半使用了基于质心的缓冲区，另一半则使用了基于边界的缓冲区。实验结果见表 6。
表 6 不同领域配置的模拟结果。
实验结果表明，当 K = 1 时，基于拓扑的选择方法表现最佳。这是因为深圳的土地使用情况复杂，许多大地块与许多其他地块连接，这会影响结果。一些地块可能有高阶连接，通过这些大地块延伸到远离的位置，这也会影响模拟结果。例如，一些地块代表了山脉、湖泊或其他长条形地块。
基于边界的缓冲方法表现良好，当缓冲距离设置为 50 米时，取得了最佳性能。考虑到边缘数量的增加可能会影响训练效果（Chiang et al., 2019; Z. et al., 2021），HashGAT-VCA中深圳的最佳邻域影响范围是 50 米的边缘缓冲区。随着缓冲距离的增加，边缘数量也增加，导致准确率逐渐下降。
基于质心的缓冲方法在 100 米的缓冲距离下表现最佳，而在 50 米时表现较差。这是因为 50 米的范围相对较小，可能导致缓冲区仍在地块内部，或者地块形状特别不规则，导致质心位于地块外部，如图 4 所示。这些因素都可能影响训练效果。尽管在100 米时取得了最佳值，但总体而言，基于边界的缓冲方法最为合适。然而，确定特定研究区域的最佳影响范围仍需进行多次实验。
4.4. 未来情景模拟
在中国，城市的发展与政府决策密切相关。然而，政策制定通常涉及较大的不确定性和不可预测性，这给城市土地使用变化的模拟带来了巨大的挑战。生态控制策略是中国政府为保护自然生态系统而实施的一项政策。在这一策略下，生态保护红线内的土地不允许开发为城市用地。提出的HashGAT-VCA模型可以通过采用受限的马尔可夫链方法和历史数据来进行未来情景预测（Fan et al., 2008; Yao et al., 2021）。首先，使用受限的马尔可夫链方法基于历史数据预测2014年到2025年和2030年的变化面积。然后，基于HashGAT-VCA模型进行多次迭代。如果达到预定的变化面积，则停止迭代。实验分析在考虑了生态保护红线的情况下进行。图15展示了模拟结果。
图 15 在生态控制下，2025年和2030年的模拟结果。
根据表7，从2014年到2025年，非城市用地转变为城市用地的面积为97.101平方公里。从2025年到2030年，转换面积为45.472平方公里。最右侧列表示了2025年和2030年的模拟土地使用情况，其中用红色框突出显示了这些年份发生变化的单个地块。图16展示了从2014年到2025年和从2025年到2030年的详细模拟变化。这些区域都是人口密集和发展迅速的地区。在2025年，宝安区和坪山区将有显著的工业发展。此外，龙华区和宝安区南部将有大量非城市用地转变为住宅用地。另一方面，大鹏新区的发展有限，仅有少数地块发生了转变。
表 7 生态控制下，各类土地使用变化的总模拟面积（单位：km²）。
图 16 生态控制下，2014–2025年和2025–2030年的模拟变化结果详情。(a1-a4): 2014–2025年模拟变化结果详情。(b1-b4): 2025–2030年模拟变化结果详情。(1) 宝安区；(2)龙华区和南山区；(3) 大鹏区；(4) 龙岗区和坪山区。
到2030年，宝安区、龙华区、南山区和龙岗区将出现大量住宅用地开发。这是由于工业区的扩张、人口增长以及需要平衡人口与住宅用地开发（Peng, 1997）。坪山区将继续开发少量工业用地，而大鹏新区只有少数地块发生变化。这是因为大鹏新区的城市地块较少，发展缓慢，大多数地块位于生态控制区内。
模拟结果表明，深圳市未来将经历显著的工业发展，特别是在宝安区、光明区和龙华区等西部地区，以及坪山区东部的工业集群。随着住宅用地的增加，人口也在不断增长。然而，商业用地的变化相对较少。平衡就业、住房和娱乐是一个关键问题。
预计宝安区和龙华区将实现均衡发展，商业、工业和住宅区都将增长。然而，光明区和坪山区作为重要的工业和人口中心，预计商业发展相对较弱。必须有效规划商业区域，以满足不断增长的人口的就业和生活需求。同时，改善相应的公共基础设施也是支持这些发展的关键。
【讨论和结论】
为了在评估环境驱动效应时考虑土地地块内驱动因素的空间异质性，并更好地捕捉土地利用变化的邻域效应，本研究提出了一种结合哈希函数和图注意力网络（GAT）的VCA方法（HashGAT-VCA）。哈希函数（即SHA-512）用于将不规则形状和大小的土地地块内的驱动因素值编码为固定长度的向量，并与土地地块的其他特征（如土地利用类型）结合，形成土地地块的特征向量。通过构建土地地块的图来表示地块之间的空间拓扑关系。GAT同时提取每个土地地块内及其邻近地块的驱动因素影响和地块间的相互关系，从而评估土地地块的土地利用转型适宜性。
所提出的HashGAT-VCA模型被应用于模拟深圳市2009年至2012年的城市土地利用变化。结果表明，与其他VCA模型（CNN-VCA、RF-VCA、ANN-VCA和LR-VCA）相比，HashGAT-VCA模型的模拟准确性最佳，相比CNN-VCA模型提高了5.81%。HashGAT-VCA模型模拟的土地利用细节与实际情况大致一致，特别是在实际发生变化的地块上。在形态特征上，许多都相同。与CNN-VCA模型相比，HashGAT-VCA模型对变化地块的敏感性更强，有效捕捉了哪些地块发生了变化及变化类型。这些结果表明，基于GAT的方法能够高效提取每个地块邻域内驱动因素的高级特征，从而准确识别城市土地利用变化模式，提高模拟准确性。
为了验证哈希函数的有效性和GAT的性能，使用了土地地块内驱动因素的平均值作为对应特征。结果表明，使用平均值可能忽略地块内的空间异质性，导致注意力的差异和模拟准确性的降低。此外，使用不同的哈希函数进行实验，并将GCN作为GAT的替代方案。实验结果表明，使用平均值作为节点特征会丧失地块内驱动因素的异质性，从而导致准确性的下降。SHA-512由于编码长度较长，可以保留更多的驱动因素信息，因此能够实现更好的模拟准确性，相比于其他哈希函数如SHA-256和MD5。与GCN相比，GAT在提取地块之间的邻域关系方面表现更佳，从而提高了准确性。
通过使用约束Markov链方法和所提出的HashGAT-VCA模型，预测了深圳市生态控制战略下2025年和2030年的土地利用。研究结果表明，深圳市的工业集群发展显著，而商业区域的发展相对有限。未来应重点推动光明区和坪山区的商业发展，同时平衡人口、就业和住宅区域之间的关系。研究结果表明，所提出的模型可以作为预测城市土地利用变化和支持未来土地利用政策决策的工具。
作为一种复杂的社会经济现象，城市土地利用变化受到多种因素的影响，这些因素具有不确定性和随机性（Yang et al., 2023）。准确模拟土地利用变化是具有挑战性的，特别是在城市发展、政府政策和城市规划之间存在紧密关系的情况下（Tian and Qiao, 2014）。在大面积区域中土地转换不频繁的情况下，准确匹配实际变化的难度更大。所提出的HashGAT-VCA模型在模拟性能上表现出色，代表了VCA模型的一项进展。
本研究也存在一些需要在未来工作中解决的局限性。在快速发展的城市中，城市土地利用受到政策和土地规划的巨大影响。未来的研究需要将政策的影响纳入模型中。由于获取详细的人口和经济数据存在困难，考虑如何将这些影响因素纳入模型也很重要。在GAT的训练过程中，随着缓冲区的增加，连接边数也增加，这可能增加训练难度。
资金来源
这项工作得到了以下资助支持：国家重点研发计划（中国）[2023YFB3906803]，国家自然科学基金（中国）[42171466]；中国地质大学（武汉）“CUG Scholar”科研基金 [2022034]，以及粤港澳联合实验室计划 [2020B1212030009]。
作者贡献申明
Qingfeng Guan: 撰写—审稿与编辑，撰写—初稿，项目管理，调查。Jianfeng Li: 撰写—初稿，方法论，数据整理。Yaqian Zhai: 验证，方法论，数据整理。Xun Liang: 撰写—审稿与编辑，验证，方法论。Yao Yao: 撰写—审稿与编辑，撰写—初稿，监督，项目管理，方法论，资金获取，概念构建。
