珠江三角洲土地利用变化模拟与景观格局脆弱性分析

［摘 要］ 随着经济的快速发展，珠江三角洲地区土地利用出现了剧烈变化，生态环境受到的胁迫明显增加，地区可持续发展进程受到严峻挑战。基于FLUS模型，以2020年珠江三角洲土地利用数据为基础，模拟2035年3种情景下研究区土地利用变化情况，并在此基础上分析研究区的景观脆弱性。结果表明，3种情景下耕地面积均呈减少趋势，建设用地面积均呈增长趋势，但在协调发展情景下研究区2035年土地利用格局较为合理；3种情景下研究区中部高景观脆弱性范围均显著增加，中部外围低景观脆弱性范围均有所减少，但经济扩张情景下研究区全域景观脆弱性水平最高，协调发展情景下全域景观脆弱性水平最低，土地利用布局最合理。

［关键词］ 景观脆弱性；FLUS模型；土地利用；多情景模拟；空间自相关

［中图分类号］ F301.24；P901 ［文献标志码］ A ［文章编号］ 1674-7909（2022）19--6

0 引言

景观格局指景观的空间结构特征，是景观组成单元的类型、数目及空间分布与配置，是景观异质性在空间上的综合表现［1］。一般利用有生态学意义的景观格局指数描述景观格局的空间配置特征，景观格局的时空动态变化研究能够揭示景观格局与生态过程的相互作用关系［2］。因此，基于景观格局的时空演化过程研究，不仅能反映区域生态环境变化状况，也能有效预测生态过程的特征，进而有效指引区域景观格局优化和管理［3］。

随着社会经济的迅速发展，人类对土地资源的开发不断加重，生态环境受到的压力和胁迫不断增加。这加剧了景观结构和组成的变化，由此产生的生态问题屡见不鲜［4］。在这一背景下开展景观格局脆弱性研究尤其必要。景观格局脆弱性是指景观格局在受到外界扰动时其结构、功能和特性容易发生改变的属性［5］，体现了景观系统对外界干扰的敏感程度及适应能力［6］。任志远等［7］和梁佳欣等［8］分别利用景观脆弱度指数对银川盆地和山东省南四湖湿地景观格局脆弱度进行时空分异特征研究，认为景观脆弱度指数能够直观反映研究区脆弱度状况，揭示人为活动与景观脆弱性的关系。

当前，学者对景观格局脆弱性的研究主要包括景观格局脆弱性时空分异研究［9-10］、景观格局脆弱性的时空分异成因研究［11-12］等方面，针对景观格局脆弱性的预测性研究还有所欠缺。曾永年等［13］和杨阳等［14］分别利用CLUE-S模型，以青海高原东部和福建省平潭县海坛岛为例进行土地利用变化多情景模拟，并利用景观生态风险指数进行区域生态风险评价，研究结果有助于决策者判断不同情境的优劣程度，从而支撑土地利用科学规划。然而，有研究指出CLUE-S模型存在一定的不足［15］。而FLUS模型作为一种新型土地利用变化模拟模型能够克服CLUE-S模型存在的缺陷，获得更好的土地利用模拟效果［16-17］。因此，有必要在利用FLUS模型模拟土地利用多情景变化的基础上，对模拟结果进行景观格局脆弱性评价。所以，笔者以珠江三角洲地区为例，使用系统动力学（System Dynamics，SD）模型，依据对研究区未来发展的不同预计构建3种不同的土地利用变化路径，在此基础上利用FLUS模型模拟2035年土地利用变化状况（以2020年为基期），再使用景观格局脆弱性指数分析研究区不同情景下面临的景观风险，从而为研究区土地利用规划提供决策依据。

1 研究区概况与数据来源

1.1 研究区概况

珠江三角洲地区位于广东省中南部（21°40′～111°20′～115°30′E），是我国重要的经济中心，当地海运便利、交通通达，总面积5.54万km2（见图1）。2020年，珠江三角洲地区人口为7 801.43万人，生产总值为89 523.93亿元，城市化率为86.28%。随着城镇化的快速发展，珠江三角洲地区土地利用出现了剧烈变化，耕地、水体面积明显减小，建设用地面积显著增加，区域生态环境受到的胁迫明显增加，可持续发展进程受到严峻挑战。

1.2 数据来源及处理过程

笔者所使用数据包括土地利用数据和土地利用驱动因子数据，其中驱动因子数据包括基础地理数据、自然环境数据、社会经济数据。土地利用数据来自全球30 m地表覆盖数据集（GlobeLand30），从中提取研究区2020年土地利用数据，并根据研究区实际情况，将土地利用数据划分为耕地、林地、草地、水体、建设用地和未利用地六类。参考以往研究［18］，选取13个影响土地利用变化的驱动因子（见表1），再根据先前研究结论将所有数据重采样为100 m［17］，同时将所有数据统一为UTM投影坐标系。表1中，共有13个土地利用变化驱动因子，余下为指标因子，在系统动力学模型中使用。

2 研究方法

2.1 FLUS模型

FLUS模型由基于SD模型的土地利用数量预测模块和基于元胞自动机（Cellular Automaton，CA）的地类分配模块两部分组成。首先，利用SD模型根据研究区社会经济发展、开发政策导向的不同预期推演研究区未来土地利用数量变化。其次，使用人工神经网络（Artificial Neural Network，ANN）基于土地利用变化驱动因子估计不同地类转化概率，接着将地类转化概率、元胞间相互作用、地类变化趋势等因素结合起来计算元胞转化总概率。最后，利用轮盘赌来确定元胞是否发生类型转化，从而实现土地利用变化模拟。FLUS模型总体结构详见张经度等和Liu等研究［17-18］。

2.2 SD模型

笔者使用Vensim软件构建研究区未来土地利用需求的SD模型（见图2）。模型包括社会经济和土地利用两个子系统，前者以人口增长率和地区生产总值增长率为主导因素，主要模拟人口增长对粮食需求的影响及地区生产总值增长对固定资产投资的影响；后者包括6种土地利用类型，重点模拟社会经济因素及地类间相互作用对土地利用的影响，另选用城镇化率、技术进步（以粮食单位面积产量衡量）及退耕还林政策作为辅助因素。

2.3 景观脆弱度指数

景观脆弱度指数（Landscape Vulnerability Index，LVI）是生态系统对外界干扰的敏感性及受到干扰后自我恢复能力的表征［19］，由景观敏感度指数（Landscape Sensitivity Index，LSI）和景观适应度指数（Landscape Adaptation Index，LAI）组成。其计算公式为

LVI=LSI×（1-LAI）（1）

景观敏感度指数体现了景观受到外界干扰时自身的反映程度，取决于外界干扰因素的强弱及景观变化的方向。其计算公式为

（2）

式（2）中：n为景观类型数量；i为用地类型；Vi为景观易损度指数，指景观在外界干扰下的损害程度［6］，参考前人研究和研究区实际，设定耕地、林地、草地、水体、建设用地和未利用土地的景观类型易损度分别为0.14、0.23、0.23、0.04、0.04和0.32，实际计算时若某一单元格网内不包含某种地类，则其余地类的易损度值按比例扩大；Ui为景观干扰度指数，指景观受到外界干扰时，景观内部由单一、规则、均值和连续的整体向破碎、零散、异质和不连续的破碎斑块变化［5］，其计算公式为

Ui=FNi×FDi×DOi（3）

式（3）中：FNi为景观破碎度指数，值越低表示景观破碎程度越低，反之则越高；FDi为分形维数，值越接近1说明景观形状越简单，反之则越复杂，自然度越强；DOi为景观类型优势度，该指数表示景观内部某种景观的优势程度，值越大则景观优势程度越大，值越小则景观优势程度越小。

景观适应度指数是景观在外界干扰下的适应和恢复能力，与景观的结构、功能、多样性和分布均匀程度密切相关，其计算公式为

LAI=PRD×SHDI×SHEI（4）

式（4）中：PRD为斑块丰度密度指数，表示斑块内物种的丰度；SHDI为香农多样性指数，表示斑块内物种多样性；SHEI为香农均匀度指数，表示不同生态系统的分配均匀度，值越高表明景观格局越均匀，系统越稳定。

2.4 空间变异理论

空间变异理论研究对象为区域化变量，即在空间上分布的变量。通常，空间上分布的变量具有2个性质：①随机性，即在局部某一点，空间变量的取值是随机的；②结构性，即对整个区域而言，存在一个总体或平均的结构，相邻变量的取值具有相关关系［20］。由于土地利用数据为空间变量，针对土地利用数据中存在的随机性和结构性特征，可以利用半变异函数描述其性质［21］。笔者使用GS+7.0软件对研究区景观脆弱性分布进行半变异函数拟合。半变异函数公式为

（5）

式（5）中：γ（h）为间距为h的半方差，N（h）是以h为间距的所有观测点的成对数目，Z（xi）和

Z（xi+h）为不同距离的观测值。

3 实证研究

3.1 2020—2035年土地利用情景模拟与分析

根据对研究区未来发展的不同预计，通过设置不同模型参数及约束条件，构建基准情景、经济扩张情景、协调发展情景，利用FLUS模型分别模拟各情景下2035年研究区土地利用分布格局。

基准情景表示研究区人口、经济发展与技术进步速度延续现有趋势，土地利用需求由2010—2020年土地利用变化量线性拟合得到。该情景体现了研究区土地利用无约束发展情况。系统动力学模型预测得到该情景下研究区2035年耕地、林地、草地、水体、建设用地和未利用土地的面积分别为8 026.41、24 524.79、2 470.72、4 341.90、13 031.51 km2和26.07 km2。

经济扩张情景表示研究区实行经济主导的发展政策，当地政府为经济发展提供充足的土地供应。该情景展现了在较快经济增速的背景下，满足建设用地扩张需求的研究区土地利用变化情况。研究预测得到该情景下研究区2035年耕地、林地、草地、水体、建设用地和未利用土地的面积分别为7 529.50、24 432.77、1 576.30、3 960.79、14 879.45 km2和42.60 km2。

《广东省国土空间规划（2020—2035年）》要求严格落实耕地和永久基本农田保护任务［22］，据此在协调发展情景中引入耕地面积反馈调节机制，使研究区耕地面积维持在基准水平之上，同时以《广东省国土空间规划（2020—2035年）》中珠江三角洲陆域生态屏障为土地利用变化限制区域。该情景体现了在保证研究区生态屏障完整及控制耕地面积前提下的研究区土地利用变化情况。研究预测得到该情景下研究区2035年耕地、林地、草地、水体、建设用地和未利用土地的面积分别为8 778.65、24 783.62、2 671.74、3 803.93、12 357.27 km2和26.19 km2。

基于不同情景预测得到的地类数量，模拟得到研究区3种情景下土地利用分布格局（见图3）。在3种情景下，地类转化的主导过程均为耕地向建设用地转移，其中经济扩张情景下建设用地扩张最为剧烈，侵占了城郊大片耕地。从各土地利用类型变化情况来看，3种情景中耕地均呈现减少趋势且减少的空间分布接近。在基准情景中，佛山市耕地面积减少近600 km2，是研究区耕地减少幅度最大的地级市，东莞市、广州市、中山市耕地也出现明显减少；在经济扩张情景中，佛山市耕地面积流失近700 km2，而广州市、东莞市、深圳市耕地面积大幅减少，减少幅度超过基准情景；在协调发展情景中，佛山市、广州市、东莞市依然是耕地面积减少幅度最大的城市，但减少幅度小于前两种情景。3种情景的林地变化幅度均较小，但变化的空间差异较大。3种情景建设用地均出现增长且空间分布类似，都呈现中心剧烈、外围缓和的增长趋势，但经济扩张情景下建设用地增长幅度明显高于其他两种情景，其中佛山市建设用地增长量达762 km2，其他两种情景建设用地增幅较接近，但协调发展情景下建设用地增幅略小。