广州市土地利用多情景模拟及其生态风险时空响应

李玮麒,兰泽英,陈德权,郑至键

(广东工业大学 管理学院,广东 广州 510520)

土地利用/覆被变化(LUCC)作为土地系统的重要表现形式，涉及到社会、经济、生态环境等方方面面的因素，体现了人类社会与自然环境最直接的相互作用关系[1-2]。近年来，随着城市化和工业化进程加快推进，对土地资源尤其是城市建设用地需求加大，导致土地资源非农化倾向增强，人地矛盾突出，由此引发的资源匮乏、生态环境恶化等问题正威胁着人类的生存与发展，因而与之息息相关的土地利用/覆被变化已成为当今全球环境变化研究的重要内容[3-4]。与此同时，计算机技术和3S技术的不断发展，为土地利用动态模拟模型的开发应用提供了强大的技术支撑，推进了土地利用变化监测分析与预测研究的展开[5-7]。土地利用变化与生态环境有密切联系，对区域景观格局与生态环境质量产生直接或间接的影响，这些影响在一定程度上会催生各种生态风险，生态风险评价在自然环境管理中应运而生并逐渐兴起[8-10]。从景观生态学的角度构建生态风险评价模型，综合评估土地利用演变对生态环境潜在的影响，对区域生态安全具有重要的预警作用，有利于探求更适宜的土地利用方式，促进区域社会经济与生态环境的协调发展[11-13]。目前，众多学者主要是基于土地利用变化的历史数据展开生态风险相关研究，且研究区域大多是小城镇、流域、湿地、海岸带以及生态脆弱区等，对于城市化水平较高的大城市进行未来景观生态风险预测的研究较少。为此，本文以广州市为研究区域，采用未来土地利用模拟(FLUS)模型和景观生态风险评价模型，在分析研究区景观生态风险时空分异特征的基础上，对模拟得到的两种情景下未来土地利用空间格局进行生态风险响应研究，以期了解景观格局变化特征与其生态效应的相互作用关系，对研究区土地资源可持续利用和维护区域生态安全具有重要理论与实践意义[14-15]。

1 研究区概况

广州市位于广东省中南部，西江、北江、东江三江汇合处，濒临南海，接近珠江流域下游入海口。其范围是东经112°57′—114°03′，北纬22°26′—23°56′，总面积约为7 434.4 km2。广州市地形多样，地势东北高西南低。东北部是以林地为主的中低山区，中部是丘陵盆地，南部是沿海冲积平原。2018年，广州市实现地区生产总值22 859.35亿元，常住人口1 490.44万人，城镇化率为86.38%。作为中国改革开放的前沿，广州市是城镇化进程快速推进的代表地区，区域内土地利用变化显著，城市发展与生态系统的矛盾日益突出，生态风险持续加剧。

2 数据来源与研究方法

2.1 数据来源与预处理

本文选用的2001，2009年和2017年3个时期的Landsat系列遥感影像数据与DEM影像数据来源于地理空间数据云平台(http:∥www.gscloud.cn)，空间分辨率均为30 m。河流水系、交通道路、行政区界和行政中心等基础地理数据均根据广州市土地利用数据提取处理获得。相关规划及政策约束数据来源于广州市规划和自然资源局。社会经济数据来源于《广州市统计年鉴》和《广州市国民经济和社会发展统计公报》等。基于ENVI 5.4软件，经过影像配准、几何精校正、去条带、镶嵌与裁剪等预处理后，采用神经网络法与人工目视解译法相结合的方法对遥感影像进行分类，结合广州市实际土地利用特征，并参考《土地利用现状分类(GB/T21010-2017)》，将土地利用类型划分为林地、草地、耕地、水域、建设用地和未利用地6类，由此得到广州市2001，2009和2017年的土地利用分类图。在此基础上进行精度检验，得到kappa系数均在80%以上，基本满足研究需要。

2.2 研究方法

2.2.1 土地利用动态度 土地利用动态度可以定量分析区域土地利用的动态变化，对预测未来区域土地利用格局变化发挥重要作用[16]。本文采用单一土地利用动态度衡量研究区域在一定时段内某种土地利用类型的变化幅度，其为正值时表示某类型面积增加，为负值时表示该类型面积减少。可表示为：

(1)

式中：K表示为研究时段内某种土地利用类型动态度；Ua，Ub分别为某种土地利用类型在研究期初和期末的数量；T为研究时段(a)。

2.2.2 FLUS模型 FLUS模型是用于模拟人类活动与自然影响下的土地利用变化以及未来土地利用情景的模型[17]，兼顾了“自上而下”的系统动力学(SD)模型和“自下而上”的元胞自动机(CA)模型，增强了对未来土地利用格局的模拟预测能力[18-19]。主要由基于神经网络的适宜性概率计算和基于自适应惯性机制的元胞自动机预测两个部分组成。首先，根据广州市实际情况及数据的可获取性，从自然地理、交通区位和社会经济三个方面选取了高程、坡度、距铁路距离、距高速公路距离、距国道距离、距省道距离、距水系距离、距行政中心距离、人口密度和地均GDP共10个驱动因子。在整合土地利用数据和各驱动因子数据基础上，采用神经网络算法运算得到研究区各用地类型的适宜性概率。其次，结合未来各类用地的目标像元数、转移成本矩阵和邻域因子等相关参数，基于轮盘赌选择的自适宜惯性机制CA模型实现研究区土地利用变化模拟。本研究采用Markov模型得到自然发展趋势下各土地利用类型的需求预测。转换成本矩阵表示研究时段内任意两种土地利用类型相互转换的可能性，根据研究区实际用地情况可知建设用地向其他用地类型转换的可能性较低，设为0；而其他各类型用地间均可相互转换，设为1。而邻域因子参数表示该土地利用类型扩张能力的强弱，在参考相关文献[20-21]的基础上，结合研究区土地利用类型转移情况，预先设置各地类邻域因子参数进行模拟试验，并将模拟结果与实际土地利用情况进行对比，经过模型的反复调试，最终设置研究区林地、耕地、草地、建设用地、水域及未利用地的邻域因子参数分别为：0.8，0.6，0.1，0.9，0.7，0.1(参数范围为0～1，越接近1表示土地利用类型的扩张能力越强)。

2.2.3 景观生态风险分析

(1) 生态风险单元划分。根据研究区面积及各土地利用类型分布的实际情况，采用等间距系统采样法，将广州市土地利用数据进行3 km×3 km的空间网格划分，共得到950个风险单元。在此基础上计算每个格网单元的生态风险指数，并将各个景观生态风险值赋予相应格网单元的中心点。

(2) 景观生态风险评价模型构建。本研究综合考虑区域景观格局特征及各土地利用类型对生态风险的影响程度，耦合不同土地利用类型面积比例、景观损失度指数以及土地利用生态风险强度参数得到景观生态风险指数的评价模型。计算公式为：

(2)

式中：ERIk为第k个风险单元的景观生态风险指数；n为景观类型数量；Aki为第k个风险单元景观类型i的面积；Ak为第k个风险单元的总面积；Ri为景观损失度；Wi为研究区内景观类型i所反映的生态风险强度参数。

①景观损失度指数Ri。表示各景观类型受到自然环境或人类活动等干扰时所造成的生态损失程度，是景观干扰度和景观脆弱度的综合反映。计算公式为：

(3)

式中：Ei为景观干扰度指数，通过景观破碎度C，景观分离度Si和景观优势度Ui构建，各景观指数可由Fragstats 4.2软件计算得到；a，b，c为相应景观指数的影响系数，且总和为1，根据相关研究成果[22]，将3个指数分别赋值为0.5，0.3，0.2；为景观脆弱度指数，结合研究区实际情况，通过专家咨询法对6类景观类型赋值，将未利用地赋为6，水域赋为5，耕地赋为4，草地赋为3，林地赋为2，建设用地赋为1，并进行归一化处理。

②生态风险强度参数Wi。代表各景观类型生态风险的主观权重值，可依据层次分析法原理和方法在AHP软件中确定[23]。通过建立层次结构模型，构造判断矩阵，计算指标的标准化权重，并对总排序进行一致性检验，最终计算一致性比例为0.024(<0.1)，表明判断矩阵具有满意的一致性。从而得到各土地利用类型的生态风险主观权重分别为：林地0.045 4，耕地0.147 1，草地0.111 3，建设用地0.396 5，水域0.077 2，未利用地0.222 5。

(3) 地统计分析法。地统计学是描述空间上某一系列变量的相关关系和格局的统计方法[24]。本文将950个样区的生态风险指数赋给各样区中心点，借助地统计软件GS+7.0实现样本半变异函数理论模型的最优拟合，在此基础上利用普通克里金法对样本点进行空间插值，得出研究区生态风险指数空间分布图。具体计算公式为：

(4)

式中：γ(h)为半变异函数；h为样点空间间隔距离，称为步长；N(h)为间隔距离为h的样点数；z(xi)，z(xi+h)分别为在 两点处的生态风险指数值。

为了直观地描述景观生态风险程度的空间分布特征，本文采用自然间断点分级法将研究区生态风险值划分为5个等级：低生态风险区(ERI≤0.012)、较低生态风险区(0.0120.033)。

3 结果与分析

3.1 土地利用变化分析

3.1.1 土地利用面积变化 根据Landsat系列遥感影像解译获取广州市土地利用变化数据(见表1)。林地、耕地和建设用地是广州市主要的土地利用类型，分布范围较为广泛。2001—2017年，仅建设用地面积持续增加，共增加了662.48 km2；而其他用地类型规模均不断减少，其中耕地的减少最为突出，减少了355.04 km2，主要由于城市化进程的快速推进致使人地矛盾突出，而耕地大多分布于城市或村镇周边，加之地势平缓，水土条件较好，而且征收成本低廉，因此耕地非农化现象频发，耕地面积的锐减为建设用地的扩张提供大量补充用地。

表1 广州市2001-2017年土地利用变化

3.1.2 土地利用动态度变化 根据公式计算得到广州市不同时段内的单一土地利用动态度(见表2)，以此描述广州市土地利用类型变化幅度和趋势。2001—2017年，广州市草地、建设用地和未利用地的变化幅度较大，其中建设用地为4.76%，呈明显增加趋势，是因为城市化进程中以建设用地扩张为主；草地和未利用地分别为-5.22%和-4.04%，迅速减少，则是由于二者面积较小，占土地总面积的比例非常小，且分布较不稳定，减少的量相对于其基数较大导致变化幅度明显。而林地、耕地和水域变化幅度较小，皆呈逐渐减少的趋势，且减少的幅度相对比较稳定。

表2 广州市2001-2017年土地利用动态度

3.2 土地利用变化模拟预测

为了验证FLUS模型对土地利用模拟结果的可靠性，本研究以广州市2009年土地利用数据为基础，结合土地利用适宜性概率分布情况，并设置相关参数，得到2017年土地利用空间格局模拟结果，并将2017年实际土地利用状况与预测结果进行精度验证。最终计算求得总体精度为80.28%，kappa系数为0.763 4。说明此次模拟结果具有较高的准确度，故可基于FLUS模拟对广州市未来土地利用变化进行模拟预测。依据广州市历史年份的土地利用变化特征及未来区域空间发展规划，分别设定自然发展情景和“三线”约束情景两种模式对其2025年土地利用变化进行模拟预测。在自然发展情景(情景Ⅰ)下，不考虑任何规划政策对土地利用变化产生的约束性影响，按照2009—2017年的土地利用类型转移概率矩阵，基于Markov模型预测得到2025年土地利用目标像元数。在“三线”约束情景(情景Ⅱ)下，结合考虑生态、农业、城镇等用地结构，选取《广州市国土空间总体规划(2018—2035年)》中划定的生态保护红线、永久基本农田保护红线和城镇开发边界线作为土地利用限制转化区域，将耕地、林地转换为建设用地的概率减少30%，草地和水域转换为建设用地的概率减少20%，由此得到广州市2025年“三线”约束情景下土地利用目标像元数。在此基础上，综合2017年土地利用数据和各类驱动因子数据，运用神经网络模型生成2025年各类用地的适宜性概率文件，通过设置相关参数，模拟预测得到广州市2025年自然发展情景和“三线”约束情景下的土地利用空间分布格局(见封2，附图1)。与2017年对比，广州市2025年两种情景下各土地利用类型呈现不同程度的变化(见表3)。

表3 广州市2017与2025年两情景各用地类型面积对比 km2

2025年两种土地利用情景下草地和未利用地的变化程度不明显，建设用地面积均有所增加，林地和水域面积均有所减少，而耕地面积则呈现不同的变化。在自然发展情景下，广州市各类用地按照原先的转移概率发生转变，建设用地在原有基础上呈现外延式扩张，逐渐向中心城区周边区域扩展，其面积增加幅度较为明显，相比2017年增加了225.99 km2；其余用地类型面积均有所减少，其中耕地面积减少量最多，减少了100.58 km2。这说明仍有大量耕地、林地、水域等用地转化为建设用地，主要集中在花都区、白云区、番禺区和增城区等。而在“三线”约束情景下，建设用地扩张幅度明显缩小，面积仅增加了92.29 km2，且林地和水域面积减少幅度有所缓和，耕地面积较2017年相比反而有小幅增加，增加了9.21 km2。相较于自然发展情景，各用地类型的转变在不同程度上得到调控，建设用地的扩张受到一定制约，生态用地得到进一步保护。

3.3 基于土地利用变化的景观生态风险评价

通过上述景观生态风险评价模型可得，广州市2001，2009，2017年景观生态风险指数分别为0.019 7，0.020 3和0.020 9，而2025年自然发展情景与“三线”约束情景下的景观生态风险指数则为0.021 2和0.021 0，表明研究区的生态安全状况趋于恶化，整体处于较低风险等级向中风险等级过渡的阶段。相较于2017年和2025年两种情景下的生态风险指数均有所升高，而自然发展情景下土地利用格局的生态风险指数增幅稍为明显。分析广州市景观生态风险空间分布特征(见图1)和各生态风险等级区域面积变化情况(见表4)，可以得到，广州市景观生态风险等级在空间上均呈现圈层状分布，由中心城区向外围递减。2001—2017年广州市景观生态风险等级面积变化显著，等级低的风险区面积减少，等级高的风险区面积增加，说明研究区景观生态风险不断加剧。而2025年两种土地利用模拟情景下各生态风险等级区域的空间分布大致稳定，低生态风险区域面积均最大，高生态风险区面积最小。但与2017年土地利用生态风险空间分布图相比，二者呈现不同的变化趋势。

图1 2001-2025年广州市景观生态风险空间分布

表4 广州市各生态风险等级区域面积变化

3.3.1 广州市2001—2017年景观生态风险变化 处于低、较低生态风险范围的区域占比较大，其中低生态风险区面积先增后减，总体增加了92.47 km2；而较低生态风险区面积的变化最为剧烈，减少了672.09 km2，二者占土地总面积的比例由53.58%减少至45.75%。该区域主要分布在研究区东北部林地密集及南部水域分布的区域。随着城市人口增加及城镇化规模扩张，土地开发利用程度提高，连片分布的林地、耕地等农用地和水域受外界干扰的现象越来越突出，面积有所减少，从而导致景观生态风险增加，较低等级风险区域相应地减少。中生态风险区面积居中，属于低生态风险区和高生态风险区的过渡地带，呈逐渐缩小的趋势，面积减少了121.32 km2，占土地总面积的比例由28.06%减少至26.42%。该区域在研究区内分布较为广泛，主要集中在耕地和水域连片分布的区域，并与其变化范围近似一致。主要是由于农业结构调整以及城市新增建设用地侵占耕地、水域等，使其面积锐减，景观破碎度和分离度上升，生态损失程度也加深，在一定程度上加剧了其景观生态风险程度。处于较高、高生态风险范围的区域面积较小，二者呈不断增大的趋势，面积共增加了700.94 km2，占土地总面积的比例由18.36%增加至27.83%。该区域集中分布在广州市中心城区附近，且在研究时段内不断向外围地形较为平坦的区域延伸，与建设用地扩张范围有较高吻合度。说明广州市处于快速城市化进程，建设用地增加对其景观生态风险程度有较为显著的影响。

3.3.2 广州市2025年两种情景下景观生态风险变化 与2017年相比，2025年自然发展情景下的低、较低及高生态风险区面积有所下降，而较高生态风险区的面积增长幅度较大，增加了303.62 km2。这主要因为在自然发展情景下，土地利用按照现有发展趋势发生变化，建设用地持续外延式扩张，侵占周边的耕地、林地等生态用地，使生态用地的景观破碎度和分离度增高，因而景观损失度及生态风险指数也随之增大，最终导致低、较低风险区域面积减小，中、较高生态风险区面积有所增加。2025年“三线”约束情景较2017年相比，较低、中生态风险区面积趋于下降，低、较高及高生态风险区面积则有所增长。而与自然发展情景相比，二者的生态风险差异较为显著，其中最明显的差别在于土地利用的低、较低及高生态风险区规模更大。主要原因是该情景下在设置土地利用限制转化区域时，考虑了生态保护红线、永久基本农田保护红线和城镇开发边界对土地利用变化所起到的约束作用，一方面限制了林地、耕地和水域等生态用地向其他地类的转化，增强了其稳定性，使其减少规模得以控制，因而生态用地集中分布的低与较低生态风险区面积较大；另一方面城市建设用地的外延式扩张受到制约，新增面积仅集中在局部小范围内，空间结构的整合性不强，故高生态风险区规模相对较大且连续性增强。由此可见，两种情景下的生态风险空间分布是对情景模拟的预测结果的响应。“三线”约束情景下充分考虑了农业、生态和城镇空间协调可持续发展的要求，有利于严格控制国土空间开发强度，促进节约集约用地，提高土地利用效率，更符合广州市未来城市发展和生态建设的需要。

4 结 论

(1) 2001—2017年，林地、耕地和建设用地是广州市主要土地利用类型。研究区建设用地的面积持续增加，而其他用地类型总体减少。此外，建设用地、草地和未利用地的变化幅度较大，建设用地明显增加，而草地和未利用地迅速减少。

(2) 相较于2017年，广州市2025年自然发展情景下建设用地大幅扩张，总体规模增加225.99 km2，侵占耕地、林地和水域等生态用地现象突出；而“三线”约束情景下，各用地类型的转变在不同程度上得到调控，建设用地的扩张得到制约，林地、水域等用地面积减少幅度有所缓和，耕地反而略微增加，生态用地规模保持相对稳定。

(3) 2001—2017年，广州市土地利用的生态风险指数不断上升，生态安全状况趋于恶化。从空间分异来看，广州市景观生态风险等级呈现圈层状分布，由中心城区向外围递减。低、较低生态风险区主要分布在研究区东北部林地密集及南部水域分布的区域，所占面积总体缩小；而较高、高生态风险区则分布在中心城区建设用地密集的区域，二者规模不断扩大。

(4) 2025年两种情景下的景观生态风险指数相较于2017年均有所升高，但自然发展情景下生态风险指数增幅略高。与自然发展情景相比，“三线”约束情景下土地利用的低、较低及高生态风险区面积更大，说明在加强研究区耕地、生态保护的同时，也保障了城市建设用地的发展，增强了其土地利用结构的稳定性，且土地利用空间布局更有利于减缓研究区生态风险加剧的趋势，促进区域生态系统与社会经济协调发展。