一种基于多尺度邻域分析的景观破碎度模型的构建方法与流程

1.本发明涉及生态系统景观指数模型和计算方法领域，尤其涉及一种基于多尺度邻域分析的景观破碎度模型的构建方法。


背景技术：

2.全球城市化进程不断加快，城市化强度逐渐增高，导致土地利用和土地覆盖发生剧烈变化，森林和耕地等逐渐被城市的不透水面替代，景观破碎程度逐渐加剧，城市热岛效应愈加显著，对生态环境和局部气候产生显著影响。生态景观格局的时空动态研究是环境可持续发展和规划中的重要组成部分，建立景观指数，准确衡量景观格局动态变化过程，是简单有效的方法，可帮助决策者理解景观变化，对于未来景观规划也具有重要意义。机器学习在很多领域中都表现了优越的性能，通过机器学习的方法，建立景观破碎变化的变化过程驱动力分析模型，具有很好的可行性，可用于分析将来可能发生景观破碎，为决策制定提供依据。鉴于机器学习结构和驱动力分析结果的“黑箱”特征，可以使用因子敏感性分析方法，探测驱动力因子的具体影响。
3.对地观测技术的发展，给大范围和连续变化的生态景观观测，提供了可行性，空间分析方法为景观格局研究提供了技术条件。其中邻域(邻近区域)分析方法，是空间信息建模和分析的重要依据，可简单表述为相距越近的地理实体，其属性间的关联性越强，利用该理论，可以从局部视角，详细分析景观的破碎程度及其变化相关联的因子。在分析景观格局时，空间尺度效应对于正确认识和理解景观破碎的分布规律和空间异质性具有重要意义。空间尺度效应指的是当考虑的邻域的大小、形状和方向不同时，分析结果也随之变化的现象。邻域分析中常用的方法，是设置一系列连续的邻域窗口大小，分别分析景观指数及其相关联因子的关系强弱，从而理解景观指数的时空变化规则。
4.本技术发明人在实施本发明的过程中，发现现有技术中存在如下技术问题：
5.目前常使用景观指数来反映景观空间动态变化，但是直接或间接用来衡量景观破碎的指数均存在一些不足。例如平均最邻近(mean nearest neighborhood，mnn)距离，通过计算与同类型斑块(或像素)间距离远近反映破碎情况，距离远则表示分布离散，距离近则表示分布聚集，但存在聚集程度不同但mnn指数值相同的情况，且也无法表示出尺度上的差异；邻近度指数(proximity index,prox)和相似性指数(similarity index,simi)表示景观破碎时没有量纲，它们的值缺少解释价值，常作为一种比较性的指标；聚集指数(ai)根据邻接边数目计算景观破碎度，精确表达景观差异上仍然存在不足，例如针对图1所示的三种植被覆盖的景观布局(灰色块表示绿地，即目标类型c)，虽然三种情况下的ai值均相同(均为0)，但实际的景观布局差异却很大。
6.由此可知，现有技术中的方法存在景观聚集或破碎度计算结果的区分度不高、解释不明晰的问题。


技术实现要素：

7.本发明提出一种景观破碎指数(即景观破碎度模型)的构建方法和景观变化分析方法，采用了基于邻域分析的方法，构建了不同尺度生态系统景观破碎度模型，用于解决现有技术中的方法存在的景观破碎度计算结果表达不充分的技术问题。
8.为了解决上述技术问题，本发明第一方面提供了一种基于多尺度邻域分析的景观破碎度模型的构建方法，包括：
9.s1：从遥感(或调查)数据中获取研究区土地覆盖时空数据(数据为栅格格式)；
10.s2：以栅格像素为基本单位，对研究区时空土地覆盖数据中感兴趣的土地覆盖类型c进行二值化；
11.s3：采用邻域分析方法计算每一个像素类型为c的像素与同类型其它像素的最邻近距离；
12.s4：在给定的尺度上，根据最邻近距离，确定能同时表征目标类型c景观分布差异的指标(第一指标，即指标1)以及用以表征景观分布在尺度和分布密度上差异的指标(第二指标，即指标2)，并基于上述两个指标，构建景观破碎度模型。
13.在一种实施方式中，步骤s2包括：
14.将时空土地覆盖数据中覆盖类型为c的像素值置为1，而其它类型置为0。
15.在一种实施方式中，步骤s3包括：
16.循环过滤二值化后的像素中像素值为1(即目标类型c)的像素p
i
，在给定的以p
i
为中心、边长为奇数w个像素的最大预定义正方形窗口w范围(即邻域窗口)内，寻找其与同类型像素p
j
的间距，像素p
i
与同类型像素p
j
距离dis
i,j
的计算公式为：
[0017][0018]
式中，δy表示p
i
与像素p
j
之间在纵向的间距(像素数)，δx表示p
i
与像素pj之间在横向的间距(像素数)；
[0019]
根据计算得到的像素p
i
与同类型像素p
j
距离，得到最邻近距离：
[0020]
nd
i
＝min{dis
i,j
,j＝1,2,
…
,m}
ꢀꢀꢀ
(2)
[0021]
式中，nd
i
表示最邻近距离，为dis
i,j
中的最小值，m为窗口w范围内除p
i
外标记为1(类型c)的像素总数。
[0022]
在一种实施方式中，步骤s4包括：
[0023]
对不同邻域窗口范围的最邻近距离特征统计分析，得到邻域窗口内的平均景观破碎特征，作为表征目标类型c景观分布差异的指标，其中，平均景观破碎特征的计算公式为：
[0024]
avg_nd
s,i
＝∑nd
s,i
/n
n
ꢀꢀꢀ
(3)
[0025]
其中，σnd
s,i
表示对当前以s为尺度、i为中心点的景观破碎总特征的表示，avg_nd
s,i
表示土地覆盖的类型为c时，同类型分布最短距离的平均值(指标1)，其中尺度为窗口w定义的空间大小，nd
k
表示编号为k的像素的最邻近距离(可通过公式2计算得到)，n
n
为当前s尺度下窗口w内像素值标记为1的像素总数，k为窗口内循环计算的像素编号；
[0026]
根据景观破碎总特征得到用以表征尺度和密度差异的指标avg_nd
s
(指标2)：
[0027]
avg_nd
s
＝∑nd
s,i
/n
ꢀꢀꢀ
(4)
[0028]
n表示当前尺度s邻域内所有像素的总数目；
[0029]
根据平均景观破碎特征和用以表征尺度和密度差异的指标，构建景观破碎度模型：
[0030]
frag＝avg_nd
s,i
‑
avg_nd
s
ꢀꢀꢀ
(5)
[0031]
其中，frag表示景观破碎度模型，且具有以像素为单元的距离的量纲。
[0032]
基于同样的发明构思，本发明第二方面提供了一种基于第一方面所构建的景观破碎度模型的景观变化分析方法，包括：
[0033]
选取两个时间点形成一个从t1到t2的时间段，作为感兴趣时间区间；
[0034]
在不同空间尺度下，采用景观破碎度模型对两个时间点的景观破碎度进行计算；
[0035]
根据景观破碎度计算结果，比较从t1到t2时间段景观破碎度的变化量；
[0036]
以自然因素和人为因素作为关联因子，以关联因子作为特征，建立机器学习模型；
[0037]
基于建立的机器学习进行景观破碎度的时空变化与关联的驱动因子关系分析以及尺度影响分析。
[0038]
本技术实施例中的上述一个或多个技术方案，至少具有如下一种或多种技术效果：
[0039]
本发明提出了一种基于多尺度邻域分析的生态系统景观破碎指数、计算方法，首先从遥感数据中获取研究区时空土地覆盖栅格数据；然后以像素为基本单位，对感兴趣的土地覆盖类型c进行二值提取，接着根据土地覆盖数据计算相同要素最邻近距离。由于本发明使用邻域分析的方法，在邻域内计算要素的最邻近平均距离，并确定了用以表征要素的分布差异的指标以及用以表征尺度和密度差异的指标，从而构建了精确的景观破碎度模型，基于该模型可以提高景观破碎结果的区分度，并进一步提高景观破碎度结果表达的精确性。
[0040]
进一步地，本发明还基于所构建的frag模型，提供了一个景观变化分析方法，结合邻域分析和机器学习的方法，建立景观破碎动态变化的分析模型，分析了景观破碎度的时空变化与关联因子之间的关系以及尺度对模型性能的影响。
附图说明
[0041]
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
[0042]
图1为现有技术中三种植被覆盖景观布局示意图；
[0043]
图2为本发明实施例中景观破碎度模型构建的流程图。
[0044]
图3为本发明实施例中基于frag模型进行景观分析的流程图；
[0045]
图4为本发明实施例中不同情景下景观破碎度frag计算结果比较示意图；
[0046]
图5为本发明实施例中武汉市2015与2020年土地覆盖图；
[0047]
图6为本发明实施例中武汉市2015与2020年植被覆盖破碎度frag及其变化图；
[0048]
图7为本发明实施例中尺度对分析结果的影响示意图；
[0049]
图8为本发明实施例中关联因子对景观破碎度变化的影响示意图。
具体实施方式
[0050]
本技术发明人通过大量的研究与实践发现，急需一种景观指数模型来精确地量化景观破碎度，以及相应的景观动态变化分析模型和方法，而如何建立一种精确的景观破碎指标，并基于景观破碎指数建立合适的分析方法反映景观破碎动态变化过程为解决城市生态变化问题的关键。因此，本发明基于已有的技术方法，提供一种精确表示景观破碎程度的指标，根据自然因素(地形，具体包括高程和坡度)和人为因素(土地利用和土地覆盖转换)，结合邻域分析和机器学习的方法，以植被覆盖为例，建立景观破碎动态变化的分析模型，并分析尺度对模型性能的影响。
[0051]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0052]
实施例一
[0053]
本发明实施例一提供了一种基于多尺度邻域分析的景观破碎度模型的构建方法，包括：
[0054]
s1：从遥感数据中获取研究区时空土地覆盖栅格数据；
[0055]
s2：以像素为基本单位，对研究区时空土地覆盖数据中感兴趣的植被覆盖类型c(即植被类型，包括森林、农田和草地，以下实例中c均指植被类型)进行二值化；
[0056]
s3：采用邻域分析方法计算每一个像素类型为植被c的像素与同类型其它像素的最邻近距离；
[0057]
s4：在给定的尺度上，根据最邻近距离，确定用以表征要素的分布差异的指标以及用以表征尺度和密度差异的指标，并该表征要素的分布差异的指标以及用以表征尺度和密度差异的指标，构建景观破碎度模型。
[0058]
具体来说，如图2所示，为本发明实施例中景观破碎度模型的构建流程图。从frag模型(景观破碎度模型)方面，首先从遥感数据中分类出研究区高精度土地覆盖数据，然后以像素为基本单位，对感兴趣的植被类型c，进行二值提取，即研究区内将覆盖类型为植被(c)的像素值置为1，而其它类型置为0；接着计算每一个像素值为1(表示类型为植被)与同类型其它像素的最邻近距离，最后在给定的尺度s上，统计最邻近距离分布，构建景观破碎度模型(公式5)，绘制景观破碎变化的时空分布。
[0059]
步骤s1从遥感数据中获取的研究区时空土地覆盖数据为高精度遥感土地覆盖栅格数据，可以避免数据误差对frag模型和分析结果产生直接影响。
[0060]
本技术发明人通过大量的研究与实践发现，通过距离来衡量要素的聚集和离散分布是最简单有效的方法之一，结合要素的实际分布差异性和要素密度的差异性，建立的景观破碎度模型frag指数是合理的方法。
[0061]
步骤s3是以像素为单位，根据土地覆盖数据计算所有标记为植被类型的像素最邻近距离。植被类型像素最邻近距离的计算包括如下内容，首先选择合适大小的邻域范围，寻找邻域内与邻域中心均为植被类型的像素坐标，并计算每个非中心像素到邻域中心像素的欧式距离(邻边相接距离为1，对角相接距离为以此类推，均以像素为单元)，将得到的距离进行排序取最小距离，作为当前邻域中心要素的最邻近距离。对研究区内所有要素进行
同样的计算，得到研究范围内要素的最邻近距离远近关系。再使用邻域分析的方法，在邻域内计算c类型覆盖的像素的最邻近平均距离。
[0062]
步骤s4是构建景观破碎度模型，寻找邻域内所有c类型像素，对这些像素的最邻近距离求和之后除以邻域内该要素数目，得表征要素的分布差异的指标，用于表示要素的分布差异。结合与邻域大小有关的变量，反映尺度和密度上的差异性。邻域大小相关的变量使用邻域内所有c类型要素的总数目，用最邻近距离和除以该变量得到表征尺度和密度差异的指标，既可以表示不同尺度上的植被覆盖破碎差异，又在一定程度上表示相同尺度下植被覆盖密度差异。最后根据确定的指标构建更为精确的景观破碎度模型。
[0063]
在一种实施方式中，步骤s2包括：
[0064]
将时空土地覆盖数据中覆盖类型为c的像素值置为1，而其它类型置为0。
[0065]
具体来说，将研究区时空土地覆盖数据land的像素，按照其标记类型为植被c和非植被类型，分别赋予1和0值，形成待分析的二值影像数据集landb。
[0066]
在一种实施方式中，步骤s3包括：
[0067]
循环过滤二值化后的像素中像素值为1的像素p
i
，在给定的以p
i
为中心、边长为奇数w个像素的最大预定义正方形窗口w范围内，寻找其与同类型像素p
j
的间距，像素p
i
与同类型像素p
j
距离dis
i,j
的计算公式为：
[0068][0069]
式中，δy表示p
i
与像素p
j
之间在纵向的间距，δx表示p
i
与像素pj之间在横向的间距；
[0070]
根据计算得到的像素p
i
与同类型像素p
j
距离，得到最邻近距离：
[0071]
nd
i
＝min{dis
i,j
,j＝1,2,
…
,m}
ꢀꢀꢀ
(2)
[0072]
式中，nd
i
(nearest distance,nd)表示最邻近距离，为dis
i,j
中的最小值，m为窗口w范围内除p
i
外标记为1的像素总数。
[0073]
具体来说，通过距离来衡量要素的聚集和离散分布，是最简单有效的方法之一，结合要素的实际分布差异性和要素密度的差异性，发现最邻近同类像素的间距nd，建立frag模型用以指示类型c在研究区域的聚集或离散程度。
[0074]
寻找其与同类型像素p
j
的间距，是指在给定的窗口范围内，寻找与像素p
i
类型相同(即像素值为1)的像素p
j
，并计算二者的欧式距离。从纵向(y)和横向(x)两个方向上，得到两个维度上的间距(分布记为δy和δx)。由公式(1)可知，如果p
j
与p
i
共享邻边，则dis
i,j
距离为1，若p
i
存在对角相接像素p
j
，则dis
i,j
距离为以此类推。则p
i
的最邻近距离nd
i
为dis
i,j
的最小值。
[0075]
公式(2)中，当m＝0时，为方便后续计算，最邻近距离nd
i
取(中心点距邻域窗口边界的最大距离)。
[0076]
在一种实施方式中，步骤s4包括：
[0077]
对不同邻域窗口范围的最邻近距离特征统计分析，得到邻域窗口内的植被覆盖类型的平均景观破碎特征，作为表征植被覆盖的分布差异的指标，其中，平均景观破碎特征的计算公式为：
[0078]
avg_nd
s,i
＝∑nd
s,i
/n
n
ꢀꢀꢀ
(3)
[0079]
其中，σnd
s,i
表示对当前以s为尺度、i为中心点的景观破碎总特征的表示，尺度为窗口大小，nd
k
表示编号为k的像素的最邻近距离，n
n
为当前s尺度下窗口内像素值标记为1的像素总数，k为窗口内循环计算的像素编号；
[0080]
根据景观破碎总特征得到用以表征尺度和密度差异的指标avg_nd
s
：
[0081]
avg_nd
s
＝∑nd
s,i
/n
ꢀꢀꢀ
(4)
[0082]
n表示当前尺度s邻域内所有像素的总数目；
[0083]
根据植被覆盖的平均景观破碎特征和用以表征尺度和密度差异的指标，构建景观破碎度模型：
[0084]
frag＝avg_nd
s,i
‑
avg_nd
s
ꢀꢀꢀ
(5)
[0085]
其中，frag表示植被覆盖景观破碎度模型，表示考虑了景观尺度和聚集结构情况下、以像素为基本单元的植被覆盖类型的平均距离。
[0086]
具体来说，首先进行尺度s情景下景观破碎特征计算(即表征要素的分布差异的指标的计算)。基于本发明思想，计算植被覆盖景观破碎度，可顾及多尺度，通过对不同邻域窗口范围的最邻近距离特征统计分析，给出一个具有综合特征的计算指标(表征要素的分布差异的指标)。对于空间尺度s，寻找以当前位置点为中心、以s为邻域范围内所有值为1的像素，对这些位置上nd求和，即：
[0087]
由于仅用表征要素的分布差异指标avg_nd
s,i
无法体现出尺度和密度对植被覆盖景观破碎度的影响，为了反映尺度s和类型c分布密度上的差异，对与邻域大小(尺度)相关的变量，假设当前尺度s邻域内所有像素(标记为1或0)的总数目为n，用景观破碎总特征除以该变量，可进一步反映感兴趣类型c分布密度相对于总邻域大小的尺度差异，即景观破碎总特征分配到每个邻域内的全部像素后的平均值avg_nd
s
。
[0088]
最后，根据前述确定的指标，建立的修正的基于邻域分析的多尺度生态系统景观破碎度frag模型，即公式(5)。若n
n
＝n，则frag＝0，表示邻域内所有像素均为类型c，没有破碎。因此frag取值范围≥0时，值越小表示景观破碎度越小，反之则越大。avg_nd
s,i
反映了邻域内感兴趣类型c间的平均相距距离，而avg_nd
s
是与邻域大小相关的量，在不同的邻域范围内，能够反映出尺度大小的差异。通过该方法建立的景观破碎指数frag，能够更精确地反映要素分布的细微变化，表达效果更能体现实际情况，应用简单便捷。图4比较了用景观聚集(aggregation index，相当于破碎度frag的反向表示方法)和本发明frag计算结果比较，可见相对于ai，frag很好地区分了不同情景下的景观布局，能够区分指示不同的景观特征(包括图1中所示的3种ai均为0而无法区分植被景观特征的情景)。
[0089]
需要说明的是，文献(he et al.,2000.an aggregation index(ai)to quantify spatial patterns of landscapes.landscape ecology,15(7):591
‑
601)提出了景观聚集度ai的计算，也是公认的反映景观聚集(聚集度越大，破碎度越小，反之亦然)模型，其计算方法为：ai＝g/max_e，其中g为邻域内类型为c的像素共边的实际数目，max_e的计算按如下规则：对于一个栅格数据，设当前计算的邻域区域面积为a，如果n为由像素组成但面积不超过a的最大正方形边长，m＝a
‑
n2，则当m＝0时，max_e＝2n(n
‑
1)，当m<n时，max_e＝2n(n
‑
1) 2m
–
1，当m≥0时，max_e＝2n(n
‑
1) 2m
–
2。
[0090]
实施例二
[0091]
基于同样的发明构思，本发明实施例二提供了一种基于实施例一所构建的景观破碎度模型的变化分析方法，包括：
[0092]
选取两个时间点形成一个从t1到t2的时间段，作为感兴趣时间区间；
[0093]
在不同空间尺度下，采用景观破碎度模型对两个时间点的景观破碎度进行计算；
[0094]
根据景观破碎度计算结果，比较从t1到t2时间段景观破碎度的变化量frag
diff
；
[0095]
以自然因素和人为因素作为关联驱动因子，即以关联因子作为输入特征，景观指数frag的变化作为目标因子(frag
diff
)，建立机器学习模型；
[0096]
基于建立的机器学习进行景观破碎度的时空变化与关联因子关系分析以及尺度影响分析。
[0097]
具体来说，本实施例的景观变化分析方法为景观变化驱动因子的分析方法，通过实施例一中的公式(1)计算t1时刻与t2时刻的景观破碎度，得到各时刻景观破碎的空间分布图，作为时空变化分析的基础。图3为本发明实施例中基于模型进行景观分析的流程图。基于frag景观破碎度模型，可以实现重要的景观变化分析，图3中，为了理解不同时间阶段(从t1到t2阶段)景观破碎度指数的时空变化的驱动因子，选取此两个时间点形成一个从t1到t2的时间段(1)，循环分析不同空间尺度(s＝s1,s2,
…
,s
n
)下，两个时间点的景观破碎度指数(2)，并比较从t1到t2时间段景观破碎度的变化量(3)，以自然因素和人为因素作为关联分析因子，以邻域关联因子作为特征，建立机器学习模型(4)，确定该阶段不同尺度下，景观破碎度变化与不同关联因子间的定量关系，获取驱动因子和优化的分析尺度。
[0098]
在具体实施过程中，可以通过下述方案来实现：
[0099]
首先确定t1到t2时段。
[0100]
通过邻域范围大小确定，设置一系列连续变化的邻域范围，根据研究对象的性质，邻域半径采用从1
‑
10(对应邻域窗口大小从3
×
3到21
×
21个像素)的变化，由景观破碎度frag模型，计算出t1时刻与t2时刻的景观破碎度，得到各时刻景观破碎的空间分布图，作为模型计算结果并支持其变化驱动力分析的数据基础。
[0101]
frag定性变化和定量变化运算。采用不同时刻(t1到t2)的景观破碎度空间分布做差值运算(frag
diff
)，当差值frag
diff
大于0表示景观破碎度增加，小于0表示景观破碎度减小，等于0表示景观破碎度没有发生变化；或直接利用差值做定量变化运算。
[0102]
驱动力关联建模。景观破碎的原因一般认为是人为因素和自然因素共同作用的结果，选择人为因素和自然因素关联因子，作为自变量，将frag定性变化和定量变化运算结果(不同时序上的景观破碎做差得到破碎度时空变化分布，破碎增加、减少和不变)作为因变量。通过随机采样的方式在研究区域空间采集训练和测试样本数据，作为机器学习的输入数据，进行关系建模和验证。其中，定量关系分析采用逻辑回归，定性变化分析采用包括支持向量机(svm)在内的多种机器学习方法，可根据需要，支持任意方法的选择。机器学习中超参数的确定。将采样样本分为训练集、测试集和验证集，在测试集上确定最佳超参数，用于训练学习。机器学习的可解释性不强，需要驱动因子的敏感性分析辅助解释。敏感性分析分别对每个关联因子增加或减少相应的量，根据结果的变化程度和方向判断关联因子如何对分析结果产生影响。
[0103]
驱动因子敏感性分析和尺度分析。对于定性变化分析，针对每个可能影响的驱动因子(如城市化、地形等)，在研究区内计算每个因子的均值，采用自变量增加或减少各自
10％均值的策略，计算由于因变量值的微小改变，而引起的frag变化的程度和方向，判断自变量因子如何对frag变化的大小和方向产生影响。分析不同邻域范围内结果的差异，确定景观破碎度指数应用的最适宜尺度，并分析邻域尺度对frag变化的驱动力分析结果的影响。例如在确定尺度对分析结果的影响时，尺度的差异通过邻域范围大小确定，设置一系列连续变化的邻域范围，邻域半径从1
‑
10(对应邻域窗口大小从3
×
3到21
×
21)，分析不同邻域范围内结果的差异，确定景观破碎度指数应用的最适宜尺度和尺度对分析结果的影响。
[0104]
为了更好地说明本发明的技术方案，将通过具体示例对本发明提供的模型构建方法(实施例一)以及景观变化分析方法(实施例二)进行介绍。
[0105]
下面给出本发明在分析武汉市2015(t1)到2020(t2)年城市化过程中人为因素和自然因素导致的城市绿化景观破碎中的应用。武汉市是具有千万人口的大城市，城市绿化景观变化对于城市环境十分重要，也是城市规划需要考虑的重要因素。
[0106]
1、背景
[0107]
城市的扩张是城市环境研究邻域的重点方向，不透水面替代了原始的耕地、林地和草地，导致城市热岛效应的加剧，城市绿化景观(主要是植被)破碎增加，局部气候产生变化和植被生产力的下降，对城市生态环境产生深远影响。城市化中的绿化覆盖及绿化破碎度，也是生态城市建设应重点考虑的因素,一个城市的绿化破碎度如果持续降低，可能会导致生态环境恶化，生物多样性减少，不利于城市的长远发展，通过客观地分析评价绿化破碎度的分布、变化及影响其变化的相关因素，对于城市规划具有十分重要的作用。我国提出建设环境友好型城市建设，绿地景观的分布情况和变化特征，是分析城市发展对绿化环境的重要组成部分，理解它们之间的影响关系，有助于城市的生态良好性发展。本例通过应用本发明中的景观破碎度模型和分析方法展示了武汉市城市发展过程对绿地景观产生的影响。
[0108]
2、数据
[0109]
基础数据包括土地覆盖数据和地形数据。土地覆盖采用中科院空天信息创新研究院发布的30米精细地表覆盖产品(30个地表覆盖类型)，该数据结合遥感影像和地面采样方法，具有较高的数据精度，根据应用需求，重分类为植被、城市和水体三种土地覆盖类型(图5)，t1和t2时刻的frag由此数据计算。地形数据采用全球数字高程模型srtm30 dem数据，坡度使用高程数据派生而来。
[0110]
景观破碎度frag计算
[0111]
1)获取武汉市高精度遥感土地覆盖产品数据，下载中科院空天信息创新研究院发布的30米精细地表覆盖产品，将包含武汉市所有瓦片数据进行拼接，使用武汉市矢量边界提取研究区内的土地覆盖数据。根据应用需求将土地覆盖重分类为植被、城市和水体三种土地覆盖类型用于计算分析，其中植被覆盖为感兴趣的类型c。
[0112]
2)植被最邻近距离的计算。现实植被很少存在孤立的情况，因此选择邻域半径为10(21
×
21的计算窗口，但可根据具体计算进行调整)。在此邻域内，若邻域中心为植被像素，则寻找邻域内所有植被像素，并按式(1)计算它们到邻域中心箱数的欧式距离(邻边相接距离为1，对角相接距离为以此类推)，将得到的距离进行排序取最小距离，作为当前邻域植被像素的最邻近距离；若邻域中心像素为非植被类型，则跳过该像素，不进行最邻近距离计算。使用移动窗口方法，对研究区内所有植被像素进行同样的计算，得到研究范围内植被最邻近距离远近关系。
[0113]
3)再使用邻域分析的方法，在邻域内计算植被的最邻近平均距离。寻找邻域内所有植被像素，对这些植被像素的最邻近距离求和之后除以邻域内植被像素总数目。用于粗略表示植被破碎。
[0114]
4)最后结合与邻域大小有关的变量，反映植被破碎在尺度和密度上的差异性。邻域大小相关的变量使用邻域内所有像素的总数目，用最邻近距离和除以该变量，既可以表示不同尺度上的植被破碎差异，又在一定程度上表示相同尺度下植被的密度差异。
[0115]
由公式(5)计算2015年与2020年的植被破碎度，得到两个年份植被破碎度的空间分布图，作为时空变化分析的基础(图6的(a)部分和(b)部分)。其中，图5表明，2015
‑
2020年期间，武汉市城市化化面积有明显增加，同时造成对绿化面积的挤占和可能造成绿化破碎度的变化。图6的(a)、(b)分别表示按公式(5)计算的2015年与2020年武汉市植被破碎度空间分布情况，数值越大，绿化覆盖破碎度越大；图6的(c)表示从2015年到2020年武汉市植被破碎度变化的分布情况，中度灰色表示植被破碎度增加(占总面积的42.67％)，浅灰色表示植被破碎度减小(19.56％)，深灰色表示植被破碎度未变化(37.76％)。
[0116]
4、景观分析(frag变化的驱动力因素分析)
[0117]
1)根据不同时刻的景观破碎度空间分布图做差(2020
‑
2015)，获得从2015年到2020年内植被破碎度的变化分布图。差值大于0表示植被破碎度增加，小于0表示植被破碎度减小，等于0表示植被破碎度没有发生变化(图6的(c)部分)。
[0118]
2)景观破碎frag变化的原因，是人为因素和自然因素共同作用的结果，选择各像素距新增城市的距离(d2u)和高程、坡度3个因子，作为人为因素和自然因素因子(自变量)；将不同时序上的植被破碎度frag，得到植被破碎度差值作为因变量。通过随机采样的方式获取样本数据，作为机器学习的输入数据，进行破碎变化的分类(增加、减少、不变)。本实例中采用svm计算frag变化与自变量的关系。
[0119]
3)svm中超参数的确定。将采样样本分为训练集、测试集和验证集，在测试集使用网格搜索技术，选择分类结果最好对应的参数作为svm模型超参数，用于训练学习。
[0120]
4)为确定自变量如何对因变量产生影响，进行了敏感性分析，辅助结果的解释。敏感性分析分别对每个关联因子增加或减少均值的10％，根据分类结果的变化程度和方向，判断关联因子如何对分析结果产生影响。分布改变3个自变量因子的值，重新使用训练好的模型进行分类，若分类结果差异越大，则该因子对植被破碎化的影响就越大，反之亦然；若结果中破碎度增加的类型变多，则该因子的增加会促进植被的破碎，反之亦然。
[0121]
5)确定尺度对分析结果的影响。尺度的差异通过邻域大小确定，设置一系列连续变化的邻域范围，邻域半径从1
‑
10(邻域窗口大小从3
×
3到21
×
21)，分析不同邻域范围内svm分类结果的精度差异，确定植被景观破碎度指数在应用时的最适宜尺度和尺度对svm分类结果的影响。
[0122]
5、主要结果和结论
[0123]
包含原始的2015年和2020年土地覆盖类型(图5)，植被破碎度分布及其在时空分布上的差异性(图6)，尺度对于分析结果的影响(图7)，及敏感性分析判断关联因子对植被破碎的影响强度和方向(图8)。表明，基于frag模型可以很好地表达出武汉市绿化破碎度的时空变化，基于svm的frag模型模拟，可以准确挖掘出影响frag变化的尺度效应，在5
×
5的邻域尺度下，svm模型预测frag变化的总体精度，可以达到85％以上(图7)，且能很好地反应
frag变化受不同因子的影响，例如远离城区时(通过d2u城市化因子)，frag的增加可能性会明显下降，表明城市化促进了绿化的破碎度，而坡度和高程对破碎度的变化的影响相对较弱(图8)，因此在中心城区(城市化核心区)，如何结合地形地貌，对绿化布局进行合理规划，是城市生态发展需要考虑的重要内容。
[0124]
其中，图7中，灰色折线表示svm分类结果的总体精度随邻域大小的变化情况；浅色柱状图表示svm分类结果的kappa系数随邻域大小的变化。横轴表示邻域大小(例如3x3表示由长宽均为3个像素的矩形邻域范围)，纵轴分别表示总体精度和kappa系数。可见当邻域大小为5
×
5时，模型总体精度和kappa值都最大，表明该尺度(5
×
5)下，svm对frag变化的解释和预测最好。
[0125]
图8的(a)、(b)部分分别表示对3个解释因子，即高程、d2u(距城区的距离因子)和坡度因子，在原始值增加和减少10％后，frag模拟结果相对于解释因子未变化前的更新情况。横轴为邻域大小，纵轴为frag变化(增加或减小)的百分比。柱状快的高低表示关联因子对结果影响的程度大小，正负表示关联因子对结果的影响方向，正值表示对应因子的改变，促进植被的破碎(frag增加)，负数则表示对应因子的改变减缓植被破碎(frag减小)。例如图8的(a)部分，距城市的距离(d2u)对结果的影响最大，其次是高程，坡度影响最小。正常情况下，d2u与绿化破碎度的增加呈负相关，高程和坡度的增加会加剧植被的破碎。
[0126]
以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。