基于ca和sar的海岸带土地利用变化模拟方法
【专利摘要】本发明涉及一种基于CA和SAR的海岸带土地利用变化模拟方法，包括以下步骤：1)对土地利用初始年份和结束年份进行遥感分类，基于矢量地图数据和卫星遥感影像获得影响土地利用转化的各空间变量因子数据；2)获取土地利用初始年份状态值、结束年份状态值以及各因子样本数据；3)对有效样本点建立CA规则转换，产生空间变量影响下土地的转化概率图；4)建立基于SAR的地理SARCA模型；5)利用已建立的SARCA模型进行土地利用变化的模拟和预测，并输出土地利用变化的模拟结果；6)对SARCA分别从规则拟合精度及模拟结果两方面进行模拟精度的评定；7)输出并保存模拟结果。与现有技术相比，本发明具有模拟效果好等优点。
【专利说明】
基于CA和SAR的海岸带±地利用变化模拟方法
技术领域
[0001] 本发明设及一种海岸带±地利用变化模拟方法，尤其是设及一种基于CA(元胞自 动机)与SAR(空间自回归）的海岸带±地利用变化模拟方法。
【背景技术】
[0002] 对±地利用变化模拟的方法一般如下：
[0003] 元胞自动机(CA)是一种能够"自下而上"自动演化的模型，常用于模拟地理复杂系 统及其动态演变过程，在城市扩展和±地利用变化领域尤其受到重视。国内外文献从转换 规则、元胞邻域、尺度影响、模型评价、案例应用、精度评定等多角度进行研究。元胞下一时 刻的状态由转换规则决定的，因此转换规则是CA建模的核屯、和难点、地理系统模拟科学问 题之所在，该规则直接决定了 CA模型的特征并影响其模拟性能。
[0004] 目前有很多方法可W用来获取转换规则，简单地归纳为两类：显性规则和隐性规 贝1J。显性规则指的是CA参数可W类比为空间变量的权重或规则具有明确的物理意义，通过 运些权重能够清晰了解每一个变量对±地利用变化的贡献，Logistic回归是其典型代表； 隐性规则指的是构建CA规则的参数没有明确的物理意义，神经网络是其典型代表。除了 Logi St i C回归，显性规则还包括主成分分析、偏最小二乘回归、核方法、粗集、模糊集和案例 推理等，也包括全局最优化方法，如遗传算法、粒子群算法和模拟退火算法等。除了神经网 络，隐性规则还包括蚁群智能、免疫智能、支持向量机、系统动力学、CA-Markov等。国内外文 献中设及的CA模型远不止上述所列，运些模型各有其优势与适用性，成功地应用在国内外 不同区域，实质性地推动了 CA理论与应用的发展。
[0005] 空间实体之间存在空间自相关特性，±地利用与城市格局也不例外。地理学第一 定律表明，相邻的地理实体有较高的相似性，距离较远的实体相似性较小。运意味着±地利 用在空间上可能存在聚类簇，在运些聚类簇内±地利用与空间变量之间的关系是高度相似 的。因此，构建CA模型时如果考虑空间自相关因素，所建立的模型则能很好地体现相关性及 距离因素的影响。空间自回归(Spatial AutoRegression,SAR)是一种考虑空间自相关特性 的统计学方法，在±地利用变化分析中已经得到应用，将其用于CA规则的获取能够使建立 的模型融合空间自相关特性。
【发明内容】
[0006] 本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种模拟效果好的 基于CA和SAR的海岸带±地利用变化模拟方法，采用的影响因素主要为空间距离变量，且采 用的空间自回归的地理CA模型SARCA的模拟效果优于传统的LogCA模型，且能较好的完成对 海岸带±地利用变化的模拟。
[0007] 本发明的目的可W通过W下技术方案来实现：
[000引一种基于CA和SAR的海岸带±地利用变化模拟方法，包括W下步骤：
[0009] 1)对±地利用初始年份和结束年份进行遥感分类，基于矢量地图数据和卫星遥感 影像获得影响±地利用转化的各空间变量因子数据；
[0010] 2)通过采样获取研究区域有效样本点，获取±地利用初始年份状态值、结束年份 状态值W及各因子样本数据；
[0011] 3)根据步骤1)和2)的数据对有效样本点建立CA规则转换，产生空间变量影响下± 地的转化概率图；
[001 ^ 4)利用已建立的转换规则，建立基于SAR的地理SARCA模型；
[0013] 5)利用已建立的SARCA模型进行±地利用变化的模拟和预测，并输出±地利用变 化的模拟结果；
[0014] 6)对SARCA分别从规则拟合精度及模拟结果两方面进行模拟精度的评定；
[0015] 7)输出并保存模拟结果。
[0016] 所述的步骤1)具体为：
[0017] 11)利用两期遥感图像，基于马氏距离监督分类方法获取初始年份和结束年份的 遥感分类±地利用图；
[0018] 12)选取影响±地利用变化的空间变量因子，包括到市中屯、、区中屯、、主干道、农用 地、中环线、外环线的距离；
[0019] 13)通过遥感图像数据、行政区划图、道路交通图在ArcGIS中利用欧几里德距离获 取到市中屯、、区中屯、、主干道、农用地、中环线、外环线距离的变量。
[0020] 所述的步骤2)具体为：
[0021] 14)设定采样分辨率，对±地利用图、因子图层进行采样，通过ArcGIS生成Fishnet 采样点，获取上地利用初始年份状态值、结束年份状态值W及各因子样本数据，为CA规则转 换提供可靠数据。
[0022] 所述的步骤3)具体为；
[0023] 15)利用空间自回归方法SAR，对获取的有效采样点数据W及空间变量因子数据训 练CA转换规则；
[0024] 16)利用已建立的CA转换规则，在30m空间分辨率下获取空间变量影响下的±地转 化概率Pdi分布。
[0025] 所述的获取±地转化概率数据方法如下：
[00%]假定y表示元胞状态是否转变，从时间t到t+1元胞状态从非城市Non-urban转变为 城市Urban, y记为1;从时间巧ljt+1元胞状态未发生改变，则y记为0。利用获取的空间变量数 据进行±地转化概率的测算。
[0027]所述的步骤4)具体为；
[002引17)基于步骤16)的数据建立基于SAR的地理CA模型：
[0029] 任意元胞i在时刻t+1的状态由其在时刻t的状态和元胞邻域决定，则非城市Non- urban 元胞转化为城市化 ban 元胞的全局概率表达为Pgi, t ， 可表达为：
[0030] Pgi't = PdiXPw'tXConXR (1)
[0031] 式中，Pdi是由空间变量决定的±地利用转化概率，PNi,t是半径为N的邻域影响，Con 是限制因素，R是随机因素；
[0032] 对于邻域影响的评估，CA模型采用正方形或圆形的规则邻域，mXm的Moore邻域表 达为：
[0033]
(2)
[0034] 式中，中屯、元胞巧参与计算，Sfxm 口,',t=化帕n)表示mxm邻域范围内城市元胞的 总数，选取Moore 5 X 5作为元胞邻域；
[0035] 限制因素Con表示元胞受到某种限制而不能发展并转化为城市元胞，运种限制包 括大面积水体、基本农田、生态保护区、W及公园和绿地，Con表达为：
[0036] Con = Bin(celli(t)~available) (3)
[0037] 式中，Con取值0或1，0表示该元胞不可W发展为城市元胞，1表示可W发展为城市 元胞；
[0038] 随机因素R用于模拟不确定因素导致的元胞状态转变，某元胞在附近没有城市元 胞的状态下通过随机因素提高发展概率，从非城市状态转化为城市状态；随机因素R表达 为：
[0039] R=i+(lnr)a (4)
[0040] 式中，r是0到1之间的随机数，a是随机因素R的控制参数，取值是0和10之间的整 数；
[0041] 由空间变量决定的转化概率Pdt是转换规则的核屯、部分，代表了运些因素对±地利 用的影响，并通过概率的方式影响到下一时刻的元胞状态，概率Pdt表达为：
[00创
巧)
[0043] 式中，Zi是空间变量对±地利用转变的影响，Z述过logistic回归获取，Z康达为：
[0044]
巧）
[0045] 式中，曰日是常数，aj是CA转换规则的参数，亦即空间变量Xj的权重，其中j = 1， 2，...，k;
[0046] 如果利用空间自回归(SAR)去获取CA参数，则Zi表达为：
[0047]
(7)
[004引式中，W是Zi的标准化空间权重矩阵，P是空间滞后Wzi的参数，X= (Xi，..，Xk)是向 量形式的空间变量，e是空间变量X的参数向量，即空间变量的权重向量，其意义等同于式 (6)中的a八哺In分别是残差项e的方差和空间自相关系数;基于空间邻接关系定义空间权 重W，具体地采用1阶Q邻接(如een Contiguity)权重矩阵；
[0049] 通过式(7)得到式（1)等式右侧的4项，即可测算元胞i的全局转化概率Pgi,t;在实 际计算中，logistic回归用R语言(RGu i)实现，SAR则用GeoDa实现。
[0050] 计算得到元胞全局转换概率Pgi,之后，与设定的阔值Pthd(其值在0~1之间）进行 比较，W确定元胞i在下一时刻的状态是否发生转变，当元胞i的转换概率Pgi,t大于设定阔 值Pthd时转化为城市类型，否则元胞保持状态不变：
[0化 1 ]
(8)
[0化2]
[0053] 所述的步骤5)具体为；
[0054] 18)采用GeoDa实现SAR，选用某年的±地利用格局为初始状态利用SARCA模型运行 M次(初始与结束的年份差)得到±地利用变化的模拟及预测结果；
[0055] 19)将基于两种模型模拟W及预测的±地利用变化结果输出。
[0056] 所述的步骤6)具体为；
[0057] 20)对SARCA模型模拟的±地利用结果分别从规则拟合精度及模拟结果两方面进 行模拟精度的评定；
[0058] 首先通过与遥感分类的±地利用格局进行比较，对模拟结果进行精度计算，主要 指标有总体精度（Overall Accuracy)、Ka卵a系数和景观指数(Xandscape Mehics)等。将 总体精度分解为城市化it)和非城市（Correct Rejection)两类、将错误分解为忽略性 (Mi SS或Omi ssion)和替代性(化Ise或Commi SS ion)两类，利用化ntius方法将错误同时分解 为数量错误(如antity化ror)和位置错误(Allocation Error)。忽略性错误指的是实际为 城市但模拟为非城市，即CA模型未能捕捉到的城市元胞;替代性错误指的是实际为非城市 但模拟为城市，即CA模型错误地增加的城市元胞；
[0059] 将模拟结果与遥感分类结果进行叠加，叠加结果包含5类:实际和模拟均为城市 化it)、实际为非城市模拟为城市(False)、实际为城市模拟为非城市(Miss)、实际和模拟均 为非城市(CR)、W及水体(Water)。根据目视判别显示，SARCA模型的模拟结果与实际分类结 果的差异。
[0060] 与现有技术相比，本发明能够利用空间变量因子数据W及CA规则转换建立SARCA 模型对±地利用变化进行模拟结果预测，同时为检测SARCA模型的性能，利用Logistic回归 建立CA模型化OgCA)并模拟同时期±地利用变化，在规则拟合精度及模拟结果两方面与 SARCA模型进行比较。CA规则的拟合结果表明，SAR在AIC、残差的描述性统计量和空间自相 关指标等方面均优于Logistic回归。因此SARCA模型在模拟W及预测±地利用变化结果方 面具有较高的合理性。
【附图说明】
[0061] 图1为本发明的流程图；
[0062] 图2为实施例案例研究区域图；
[0063] 图3为空间变量决定的±地利用转化概率Pdi示意图；
[0064] 图4为LogCA和SARCA两种模型下2015年的预测结果示意图。
【具体实施方式】
[0065] 下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完 整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部实施例。基于本发 明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实 施例，都应属于本发明保护的范围。
[0066] W2005-2015年间上海城市±地利用为案例区域，本案例区域位置如图2所示。为 验证SARCA模型在±地利用变化模拟上的有效性，案例中将基于Logi St ic回归的CA模型 化OgCA)作为比较对象，也模拟了上海市同期城市±地利用变化过程，结果表明SARCA的模 拟效果优于传统LogCA模型。基于元胞自动机与空间自回归的海岸带±地利用变化模拟模 型方法，包括W下步骤：
[0067] 1)首先选择上海市2005、2015年的遥感图像数据，W及行政区划图、道路交通图， 作为训练CA规则转换W及获取±地转变概率的基础数据；
[0068] 2)根据遥感图像数据利用系统抽样方法对每Ikm2网络抽取一个点，获取各空间变 量的值、±地利用初始年份和结束年份状态值，上海全域共获有6701个有效样本点；
[0069] 3)利用各年遥感图像W及行政区划图层、道路交通图层、然后利用空间分析工具 中的欧几里德距离化UClidean Distance)计算出各区域到市中屯、、到区中屯、、到主干道路、 到农用地、到中环线W及到外环线的距离数据；
[0070] 4)在ENVI软件中利用监督分类中的最大似然分类法，对上海遥感图像进行监督分 类从而对±地利用格局进行解译；
[0071] 5)利用系统抽样方法得到的有效样本点W及各空间变量值，在实际计算中，用R语 言进行Logistic回归的实现，用GeoDa来实现SAR分别获取±地转变概率。表1显示了两种模 型的拟合优度，表明SARCA模型更好，表2显示了两种模型下的转换规则参数，图3显示了两 种模型下的±地转变；
[0072] 6)利用获取的±地转变概率W及CA转换规则建立基于SAR的地理CA模型和LogCA 模型；
[0073] 7)^2005年状态为初始值分别利用SARCA模型W及LogCA模型运行10次进行从而 预测2015年±地利用变化；
[0074] 8)将模拟预测后的±地利用变化值在ArcGIS中加载，得到2015年图4的±地利用 变化图，同时模拟结果的精度显示SARCA显著优于传统LogCA模型；
[0075] 9)输出并保存可视化的结果。
[0076] 表1 CA转换规则获取的拟合效果
[0080] 表3 2015年两种模型的预测精度
[0077]
[007引
[0079]
[0081]
[0082] W上所述，仅为本发明的【具体实施方式】，但本发明的保护范围并不局限于此，任何 熟悉本技术领域的技术人员在本发明掲露的技术范围内，可轻易想到各种等效的修改或替 换，运些修改或替换都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应W权利 要求的保护范围为准。
【主权项】
1. 一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在于，包括以下步骤： 1) 对土地利用初始年份和结束年份进行遥感分类，基于矢量地图数据和卫星遥感影像 获得影响土地利用转化的各空间变量因子数据； 2) 通过采样获取研究区域有效样本点，获取土地利用初始年份状态值、结束年份状态 值以及各因子样本数据； 3) 根据步骤1)和2)的数据对有效样本点建立CA规则转换，产生空间变量影响下土地的 转化概率图； 4) 利用已建立的转换规则，建立基于SAR的地理SARCA模型； 5) 利用已建立的SARCA模型进行土地利用变化的模拟和预测，并输出土地利用变化的 模拟结果； 6) 对SARCA分别从规则拟合精度及模拟结果两方面进行模拟精度的评定； 7) 输出并保存模拟结果。2. 根据权利要求1所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤1)具体为： 11) 利用两期遥感图像，基于马氏距离监督分类方法获取初始年份和结束年份的遥感 分类土地利用图； 12) 选取影响土地利用变化的空间变量因子，包括到市中心、区中心、主干道、农用地、 中环线、外环线的距离； 13) 通过遥感图像数据、行政区划图、道路交通图在ArcGIS中利用欧几里德距离获取到 市中心、区中心、主干道、农用地、中环线、外环线距离的变量。3. 根据权利要求2所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤2)具体为： 14) 设定采样分辨率，对土地利用图、因子图层进行采样，获取土地利用初始年份状态 值、结束年份状态值以及各因子样本数据，为CA规则转换提供可靠数据。4. 根据权利要求3所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤3)具体为； 15) 利用空间自回归方法SAR，对获取的有效采样点数据以及空间变量因子数据训练CA 转换规则； 16) 利用已建立的CA转换规则，在30m空间分辨率下获取空间变量影响下的土地转化概 率Pdi分布。5. 根据权利要求4所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的获取土地转化概率数据方法如下： 假定y表示元胞状态是否转变，从时间t到t+Ι元胞状态从非城市Non-urban转变为城市 Urban，y记为1;从时间t到t+Ι元胞状态未发生改变，则y记为0。利用获取的空间变量数据进 行土地转化概率的测算。6. 根据权利要求4所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤4)具体为； 17) 基于步骤16)的数据建立基于SAR的地理CA模型： 任意元胞i在时刻t+Ι的状态由其在时刻t的状态和元胞邻域决定，则非城市Non-urban 元胞转化为城市Urban元胞的全局概率表达为Pgl, t，可表达为： Pgi, t = Pdi X Pni, t X Con X R (1) 式中，Pd1是由空间变量决定的土地利用转化概率，PNl,t是半径为N的邻域影响，Con是限 制因素，R是随机因素； 对于邻域影响的评估，C A模型采用正方形或圆形的规则邻域，m X m的M ο 〇 r e邻域表达 为：(2)式中，中心元胞i不参与计算， I示mXm邻域范围内城市元胞的总 数，选取Moore 5 X 5作为元胞邻域； 限制因素 Con表示元胞受到某种限制而不能发展并转化为城市元胞，这种限制包括大 面积水体、基本农田、生态保护区、以及公园和绿地，Con表达为： Con = Bin(celli⑴~available) (3) 式中，Con取值0或1，0表示该元胞不可以发展为城市元胞，1表示可以发展为城市元胞； 随机因素 R用于模拟不确定因素导致的元胞状态转变，某元胞在附近没有城市元胞的 状态下通过随机因素提高发展概率，从非城市状态转化为城市状态;随机因素 R表达为： R=l+(lnr)a (4) 式中，r是0到1之间的随机数，α是随机因素 R的控制参数，取值是0和10之间的整数； 由空间变量决定的转化概率Pdt是转换规则的核心部分，代表了这些因素对土地利用的 影响，并通过概率的方式影响到下一时刻的元胞状态，概率Pdt表达为：CS) 式中，Z i是空间变量对十地利用转变的影响，Z i通过logistic回归获取，Z i表达为：(6) 式中，aQ是常数，aj是CA转换规则的参数，亦即空间变量Xj的权重，其中j = l，2, ...，k; 如果利用空间自回归（SAR)去获取CA参数，则zi表达为：(7) 式中，W是Zi的标准化空间权重矩阵，P是空间滞后Wzi的参数，X=(xi，...，xk)是向量形 式的空间变量，β是空间变量X的参数向量，即空间变量的权重向量，其意义等同于式(6)中 的^，〇 2和In分别是残差项ε的方差和空间自相关系数;基于空间邻接关系定义空间权重W; 通过式(7)得到式(1)等式右侧的4项，即可测算元胞i的全局转化概率Pgl,t; 计算得到元胞全局转换概率Pgl, t之后，与设定的阈值Pthd进行比较，以确定元胞i在下一 时刻的状态是否发生转变，当元胞i的转换概率Pgl,t大于设定阈值Pthd时转化为城市类型， 否则元胞保持状态不变： (8) ,7. 根据权利要求6所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤5)具体为； 18) 采用GeoDa实现SAR，选用某年的土地利用格局为初始状态利用SARCA模型运行M次 得到土地利用变化的模拟及预测结果； 19) 将基于两种模型模拟以及预测的土地利用变化结果输出。8. 根据权利要求7所述的一种基于CA和SAR的海岸带土地利用变化模拟方法，其特征在 于，所述的步骤6)具体为； 20) 对SARCA模型模拟的土地利用结果分别从规则拟合精度及模拟结果两方面进行模 拟精度的评定； 首先通过与遥感分类的土地利用格局进行比较，对模拟结果进行精度计算，指标包括 总体精度Overall Accuracy、Kappa系数和景观指数Landscape Metrics; 将总体精度分解为城市Hit和非城市Correct Rejection两类、将错误分解为忽略性和 替代性两类，利用Pontius方法将错误同时分解为数量错误和位置错误;其中忽略性错误指 的是实际为城市但模拟为非城市，即CA模型未能捕捉到的城市元胞;替代性错误指的是实 际为非城市但模拟为城市，即CA模型错误地增加的城市元胞； 将模拟结果与遥感分类结果进行叠加，叠加结果包含5类:实际和模拟均为城市、实际 为非城市模拟为城市、实际为城市模拟为非城市、实际和模拟均为非城市、以及水体;根据 目视判别显示，SARCA模型的模拟结果与实际分类结果的差异。
【文档编号】G06F17/50GK106021751SQ201610356499
【公开日】2016年10月12日
【申请日】2016年5月26日
【发明人】冯永玖, 杨倩倩
【申请人】上海海洋大学