徐崇敏陈 瑾张露丹林 森邱荣祖胡喜生
(福建农林大学 交通与土木工程学院,福建 福州350108)
为缓解城市群社会经济利益和生态环境保护的国土空间冲突,统筹国土空间经济高质量发展,实现区域协调发展,分区协调治理方式势在必行。2010年国务院发布了《全国主体功能区规划》,“十四五”期间将继续落实主体功能区制度,构建高质量国土空间格局[1]。当前,建设主体功能区是生态文明理念下提出的国土开发和保护重大战略,各区域需持续推进主体功能区建设,按照不同主体功能区的任务和重点合理规划国土空间发展方向。
国土空间格局模拟研究已经成为热点,方法上主要分为数量变化和空间变化的模拟模型,数量预测模型常用系统动力学模型(SD)[2]、马尔科夫链模型(Markov chain)[3],空间拟模型常用CLUE-S模型[4]、元胞自动机模型(CA)[5]、SLEUTH[6]、多智能体系统(MAS)[7]、未来土地利用模拟模型(future land use simulation,FLUS)[8]。不同的模型有其相应的功能机制,Markov模型具有土地利用/覆盖数量预测功能[3],但难以模拟空间格局演变,而FLUS模型具有空间演变分布预测功能且模拟精度高[2],并运用在城市群[9-10]、特殊地域[11]、省[12-13]、市[14]、县域等[15]多尺度区域研究中。本文将FLUS模型结合Markov模型数量预测优势展开研究。目前,国内对国土空间结构预测主要是从三生空间的角度,将整片区域根据不同规划目标或不同发展情景设置主导目标和约束条件进行研究[16-18]。但是,实际上大尺度区域范围存在发展不均衡性和地理异质性等问题,如不同主体功能区的发展目标并不相同,因此用地的发展需求也不同,若将研究区域作为整体进行空间预测,则忽略了各区域内部差异性,将导致预测结果的偏差。基于此,从不同主体功能区的角度出发,考虑不同区域的地理特点、发展现状和规划目标等因素设置不同约束条件,模拟结果将更贴合未来国土空间规划发展方向,可以为主体功能区有序发展提供科学依据。
福州都市圈作为福建高质量发展的重要增长极,为深化两岸交流、提升国际化水平提供有力支撑[19]。然而目前福州都市圈仍处于发育阶段,区域内深层次合作发展机制尚不成熟,生态环境和资源承载能力并不均衡[20]。当前,对于福州都市圈国土空间格局模拟的相关研究相对缺乏。鉴于此,本研究以福州都市圈为研究范围,从不同主体功能区视角出发,分析区域土地利用/覆盖类型时空变化,提出基于历史情景和主体功能区规划相结合的邻域因子权重和转换成本矩阵参数设计方法,利用FLUS-Markov模型对比基于整体模拟和主体功能区划的分区模拟的2030年福州都市圈土地利用/覆盖空间分布差异,预测研究区国土空间格局未来发展趋势,促进区域协调发展。
福州都市圈处于中国东南沿海、福建省东北部,面对台湾省,是一个具有山海生态的城市群,地理坐标介于东经117°0′—120°0′,北纬25°0′—28°0′。研究区域包括四市一区,陆域面积达2.60×104km2,约占福建省的21.5%,2020年常住人口约1 300万人,占福建省的33.5%,地区生产总值约1.5万亿元,占福建省生产总值的34.5%。
福州都市圈作为中国批复的第2个国家级都市圈,在海峡西岸地区经济协调发展中承担着重要的战略作用。区域人口聚集,交通发达,生态良好,与台湾省联系密切,其作为21世纪海上丝绸之路的核心区与国际接轨,综合发展潜力大。但都市圈中城市快速扩张区域多分布于东部沿海,生态重点保护区域多分布于西部,存在发展不均衡问题。为保证该区域可持续发展,因地制宜地分区规划对福州都市圈的经济与生态并行发展具有重要现实意义。国家依据区域的环境承载能力、现有开发密度和未来发展潜力,将区域划分成不同主体功能区,进行统筹规划,形成生态与经济协调发展的空间格局。在《全国主体功能区规划》政策文件中,根据开发方式将中国国土空间分为以下主体功能区:优化开发区域、重点开发区域、限制开发区域和禁止开发区域。
福建省主体功能区主要以县级行政区为划分单元,福州都市圈范围(图1)内不同主体功能区的分类情况详见表1。由于限制开发区域分为农产品主产区和重点生态功能区,而重点生态功能区的功能定位以提供生态服务为主,这与禁止开发区域相似,故将重点生态功能区纳入禁止开发区域范围进行模拟。
表1 福州都市圈不同主体功能区分类情况
图1 福州都市圈主体功能区分布
数据主要包括福建省30 m分辨率土地利用/覆盖遥感监测数据(2000,2010和2020年),来源于中国科学院资源环境科学与数据中心(http:∥www.resdc.cn),土地利用/覆盖数据集已经过野外调查点随机抽样核查,总体精度达88.95%[21],并将区域内土地利用/覆盖遥感监测数据采用一级分类系统,根据土地资源及其利用属性,重分类为耕地、林地、草地、水域、建设用地和未利用土地6类;福建省30 m分辨率数字高程数据来源于中国科学院计算机网络信息中心地理空间数据云平台(http:∥www.gscloud.cn);中国河流、铁路站点、高速公路、一级公路、二级公路、三级公路、主干道数据(2014和2020年)来源于Open Street Map(www.openstreetmap.org),经ArcGIS软件欧氏距离工具处理为距离栅格数据[10],数据分辨率为300 m;中国1 km分辨率月平均气温数据集(2010和2020年)、中国1 km分辨率年降水量数据集(2010和2020年)、全国500 m分辨率类NPP-VIIRS夜间灯光数据集(2010和2020年)[22]来源于国家科技基础条件平台—国家地球系统科学数据中心(http:∥www.geodata.cn)。将所有数据格式统一为tif格式,对数据进行重采样,分辨率为30 m,均采用krasovsky_1940_Albers投影坐标系,然后将所有数据图层按掩膜提取至各主体功能区范围,并以2020年福建省土地利用/覆盖栅格数据为基准捕捉栅格,保证各主体功能区所有图层行列数一致。
1.3.1 模拟流程及参数设置 FLUS模型主要用于模拟大尺度范围各种土地利用/覆盖类型的未来空间分布情况,也是目前较为成熟且被广泛应用的土地利用预测模型[10]。该模型由人工神经网络算法(artificial neural network,ANN)和元胞自动机(CA)模块两部分组成,ANN模型可以协同整合自然环境、人类社会、经济发展等多种驱动因子的复杂影响,从而建立起不同土地利用/覆盖类型同驱动因子之间的关联,得到现状下各土地利用/覆盖类型的适宜性概率[2]。CA模型引进了基于轮盘赌选择的自适应惯性竞争机制,用于处理多种土地利用/覆盖类型在驱动因子下的复杂变化,从而实现较高精度的土地利用/覆盖空间分布模拟[8]。从基于不同主体功能区的分区模拟和基于统一约束条件的整体模拟两种模拟方式出发,在关键参数设置的基础上,利用FLUSMarkov模型(图2)分别对2030年国土空间分布格局进行模拟预测,通过单一土地利用动态度、综合土地利用动态度和土地利用程度综合指数等指标定量分析两种方式下模拟结果,以对比两种方式模拟结果与区域政策规划中空间结构的贴合度,以验证分区模拟的科学性。
图2 土地利用/覆盖空间分布格局模拟流程
采用FLUS-Markov模型对土地利用/覆盖空间分布模拟时,模型参数设置为关键环节。现土地利用/覆盖模拟研究中多采用统一参数通过整体模拟方式进行预测,缺乏考虑土地利用/覆盖类型和驱动因子在时间和空间上的差异。故针对关键参数的设置,考虑了以下因素并提出解决方案:①参考相关研究和研究区现状,考虑时空异质的驱动因子选择和分期适宜性概率分布图层计算;②结合区域内土地利用/覆盖空间分布特征,基于历史情景的土地利用/覆盖类型邻域因子权重计算;③依据政府规划文件和国土空间结构,针对政策导向的不同主体功能区内转换成本矩阵设置。
(1)驱动因子选择。土地利用/覆盖变化是由自然环境和人类活动等多种因素共同影响的结果,除高程、降水和气温等自然因素对土地利用/覆盖变化有所影响,区域交通设施情况和经济发展状态会直接或间接影响周围土地利用/覆盖变化[23]。综合对土地利用/覆盖变化驱动因子的分析及探究,最终选取自然(高程、坡向、坡度、气温、降水)、社会(河流、铁路、高速、一级、二级、三级、主干道)、经济(夜间灯光)等方面的13项驱动因子(表2)。
表2 土地利用/覆盖变化驱动因子数据说明
(2)适宜性概率分布计算。各土地利用/覆盖类型的适宜性概率分布图层由ANN模型计算得到,模型中包含输入层、隐藏层和输出层神经元网络,每一个神经元分别代表一个驱动因子,可拟合土地利用/覆盖类型分布概率与多项驱动因子的空间作用关系,用于模拟不同土地利用/覆盖类型发展概率和空间分布[24]。
由于10 a间驱动因子存在时空差异性,导致2010年与2020年的驱动因子有不同作用强度和方式,所以各土地利用/覆盖类型适宜性概率分布可能存在差异。因此,为提高模型模拟精度,在ANN模型中,分别以2010年和2020年土地利用/覆盖栅格数据为基础,输入两期相应驱动因子数据来计算2010年和2020年神经网络适宜性概率空间分布,分别用于CA模型预测2020年和2030年土地利用/覆盖结构。
(3)邻域因子权重参数设置。邻域因子权重表示某种土地利用/覆盖类型的扩张强度,取值范围为0~1,取值越靠近1表示该土地利用/覆盖类型扩张能力越强。针对邻域因子权重的设置,有学者使用各土地利用/覆盖类型历史面积变化量的绝对无量纲值来确定邻域权重参数[10],但由于福州都市圈各土地利用/覆盖类型面积存在较大差异,本文提出一种相对无量纲值计算方式。
由于主体功能区战略是2011年开始实施,计划于2020年实现,规划任务至未来,故2010—2020年不同主体功能区内各土地利用/覆盖类型的面积变化量可以很好地表示分区规划政策下各土地利用/覆盖类型的扩张程度。因此利用ArcGIS软件统计2010年和2020年各土地利用/覆盖类型的面积,计算得到10 a间各土地利用/覆盖类型面积的变化量,将面积变化量的绝对无量纲值和相对无量纲值进行比较,得出更满足数据结构需求的邻域因子权重参数。
绝对无量纲值是将各土地利用/覆盖类型的面积变化量进行归一化处理,使其值处于0~1范围内,计算公式如下:
式中:N代表土地利用/覆盖类型历史面积变化量的绝对无量纲值;X代表某土地利用/覆盖类型时间段内的变化量;Xmin代表所有土地利用/覆盖类型面积变化量中最小值;Xmax代表面积变化量中最大值。
相对无量纲值是将各土地利用/覆盖类型的面积变化量除以初始面积得到各土地利用/覆盖类型的面积变化率,再将其进行归一化处理,使其值处于0~1范围内,计算公式如下:
式中:N*代表土地利用/覆盖类型历史面积变化量的相对无量纲值;R代表某土地利用/覆盖类型时间段内面积变化率;X代表某土地利用/覆盖类型时间段内的面积变化量;X*代表某土地利用/覆盖类型初始面积;Rmin代表时间段内所有土地利用/覆盖类型面积变化率的最小值;Rmax代表时间段内面积变化率的最大值。
2010和2020年福州都市圈整体范围和各主体功能区范围各土地利用/覆盖类型面积及变化量详见表3。2010—2020年各区域范围主要以建设用地面积大幅增加和耕地林地草地面积减少为主要表现,各区域中面积减少占比最大都是耕地,其中重点、优化开发区域中建设用地的扩张还伴随着水域面积的减少,表明各主体功能区以发展为主要目的,经济发展期间建设用地的扩张多占用农业用地,同时伴随着林地和草地等生态用地的调整。与绝对无量纲值方法相比。根据各土地利用/覆盖类型历史面积计算得到的相对无量纲值可以避免由于各地类初始面积占比过大或过小的差异伴随着面积变化剧烈而导致绝对无量纲值极端的问题。因此采用相对无量纲值方法计算结果(表4)作为邻域因子权重值,可以更好地表征区域内土地利用/覆盖类型的扩张强度规律。
表3 福州都市圈2010和2020年各功能区土地利用/覆盖类型面积及变化量 km2
表4 各功能区土地利用/覆盖类型历史变化量相对无量纲值
(4)转换成本矩阵设置。转换成本矩阵是指各土地利用/覆盖类型间的转变规则,0表示各土地利用/覆盖类型间不可以转变,1表示可以转变。理论上各种转变都是被允许的,但是针对不同主体功能区受规划政策影响,在经济和生态等方面发展方向不同,故在不同主体功能区中各土地利用/覆盖类型间转换成本并不相同[13]。
分区模拟中,针对不同主体功能区的发展现状和政策规划,设置了4种不同的转换成本矩阵(表5)。重点开发区域以经济发展为主,主要分布于东部沿海区域,地势相对平坦,适合人类生活和经济发展,区域内应加快推进城镇化,除建设用地不能转换成其他用地,其他土地利用/覆盖类型间可以相互转换;优化开发区域以高质量发展为主,地跨闽江流域,城市发展规模已经较为完善,区域内应控制开发强度,对水域保护进一步加强;限制开发区域以保障农业发展为主,集中分布于研究区的中部且地势较高,区域内包含农产品主产区,应保证耕地发展,控制建设用地的侵占能力;禁止开发区域以保护生态安全为主,森林覆被率高,2020年该主体功能区内林地面积占比达75.35%,包含保障全省生态安全的重要区域,在该区域内各类用地保护等级由高到低排序为:林地、耕地、水域、草地、未利用土地、建设用地。整体模拟中,将区域整体设置统一转换成本矩阵和邻域因子权重,对土地利用/覆盖结构进行模拟预测。在区域中兼顾人地和谐及绿色发展理念,按照综合发展的需求,将各地类转换等级由高到低排序为:建设、林地、耕地、水域、草地及未利用土地,高等级不能向低等级转变,以此设置转换成本矩阵[16]。
表5 各功能区转换成本矩阵
(5)未来像元参数设置。Markov模型基于初始年份各土地利用/覆盖类型的像元数来预测未来土地利用/覆盖类型目标像元数,具有无后效型。为提高模拟精度,本文以2000,2010和2020年土地利用/覆盖数据为基础,分别预测2020和2030年土地利用/覆盖类型目标像元数。计算公式为[13]:
式中:St,St+1分别表示t时刻和t+1时刻的土地利用/覆盖状态;Pij表示在t时刻土地利用/覆盖类型i转变为j的概率。
(6)土地利用/覆盖模拟模型。将适宜性概率分布图层、各地类邻域因子权重、转换成本矩阵和目标像元数等参数输入CA模型,通过不断地循环迭代使模拟不断逼近目标值,模拟预测目标年份的土地利用/覆盖空间分布[8,24]。
1.3.2 土地利用强度计算
(1)单一土地利用类型动态度。单一土地利用类型动态度是衡量单一种土地利用/覆盖类型动态变化指标,体现研究区域在一定时间范围内某种土地利用/覆盖类型的数量变化情况[25],表达式为:
式中:V为研究时段内单一土地利用类型动态度;Aa为研究初期某种土地利用/覆盖类型的面积;Ab为研究末期某种土地利用/覆盖类型的面积;T为研究期时段长。
(2)综合土地利用类型动态度。综合土地利用类型动态度是反映研究区域中时间段内综合土地利用/覆盖类型数量变化程度指标[26],其表达式为:
式中:R表示研究区域土地利用类型综合动态度;ΔAij为研究时段内第i类土地利用/覆盖类型转化为第j类土地利用/覆盖类型面积的绝对值;Ai为研究时段初期第i类土地利用/覆盖类型的面积。
(3)土地利用程度综合指数。反映了人类对土地资源开发利用的广度和深度。根据庄大方[27]所提出的综合分析法,可以计算土地利用程度,掌握区域土地开发利用的综合水平。其计算公式为:
式中:La为土地利用程度综合指数;n为土地利用程度分级数;Ai指第i级土地利用程度分级指数,参照前人研究得到土地利用程度分级赋值[26],耕地、林地、草地、水域、建设用地和未利用土地的分级指数分别为3,2,2,2,4,1;Ci是第i级土地利用程度面积比例。
利用ANN模型生成2010和2020年福州都市圈各土地利用/覆盖类型的适宜性概率分布图层(图3—4),从而获得各地类的空间分布格局,适宜性高的区域多分布于该土地利用/覆盖类型空间周围,耕地的适宜性区域主要分布在地势平缓和沿海人口聚集区域;林地和草地的适宜性区域大面积重合,连片分布于非沿海区域;水域主要分布在闽江流域以及现有水系和沿海的周边区域;建设用地的高适宜性区域主要以福州和沿海城市群为中心,向周围逐步递减,沿海岸线城市呈现出联合发展的可能。
图3 福州都市圈2010年各土地利用/覆盖类型适宜性概率分布
为检验FLUS-Markov模型分区模拟精度,对2020年各主体功能区土地利用/覆盖空间分布进行模拟,将模拟结果与实际情况进行对比(图5),并用kappa系数进行检验。重点、优化、限制和禁止开发区域在20%随机采样时的kappa系数分别为0.94,0.91,0.85,0.98,均大于0.75,总体精度分别为96.27%,94.04%,92.45%和99.30%,均大于90%,表明在以上参数设置的基础上该模型具有较高的可信度,FLUSMarkov模型能反映区域土地利用/覆盖空间变化规律,可以对2030年主体功能区土地利用/覆盖空间分布进行模拟。
图5 福州都市圈2020年土地利用/覆盖现状及模拟结果
将分区模拟和整体模拟两种模拟方式下模拟得到2030年土地利用/覆盖空间分布结果裁剪至不同主体功能区范围,从土地利用动态度和土地利用程度等指标来量化对比(表6),探究两种方式模拟所得土地利用/覆盖空间结构与都市圈各功能区发展目标贴合度。对比单一土地利用动态度指数发现,整体模拟结果中限制开发区域耕地的扩张趋势符合区域规划要求,动态度为0.14%,但所有区域林地都在缩减,重点和优化开发区域的建设用地发展强度不够,比分区模拟结果低0.31%和0.02%,限制和禁止区域的建设用地分别有0.71%和0.78%的扩张。而分区模拟结果中随着区域生态保护力度逐渐加大,水域的单一动态度有递增的梯度变化,随着区域经济发展力度的减弱,建设用地的单一动态度逐渐减小,这与区域发展目标更贴合。对比综合土地利用动态度和土地利用程度指数发现,整体模拟的重点和优化开发区域土地利用/覆盖变化相对更平缓,土地利用程度相对较低,低效闲置土地相对更多,而限制和禁止开发区域土地利用动态度高,土地利用程度相对更高,更容易导致生态脆弱。综上所述,与整体模拟相比,分区模拟的未来土地利用/覆盖模拟结果与主体功能区未来发展目标更贴合,可以根据福州都市圈各主体功能区的发展目标提供土地利用/覆盖发展方向的科学预测。
表6 福州都市圈2020—2030年分区模拟和整体模拟各功能区土地利用动态度对比
图4 福州都市圈2020年各土地利用/覆盖类型适宜性概率分布
绘制福州都市圈各主体功能区2030年土地利用/覆盖分区模拟结果(如图6所示),并统计2020—2030年各主体功能区内土地利用/覆盖类型面积及比例变化量(表8)。重点开发区域中,建设用地作为面积变化量最大的土地利用/覆盖类型,增加了320.24 km2,比例增长1.48%。其中233.69 km2主要来源于耕地,67.20 km2来源于水域。利用桑基图对区域土地利用/覆盖转移方向进行可视化处理(图7)。同时,计算各县区内建设用地单一动态度表示建设用地扩张强度(表7),发现2030年用地空间增长多分布于沿海城市[27],以莆田市、福州市(长乐区、闽侯县、福清市)和平潭综合实验区等城市为中心向外扩张的形势(图8),建设用地联合发展态势进一步加强;该区域内耕地面积减少幅度最大(达153.23 km2),比例减少0.71%,有78.48 km2林地和20.59 km2水域转入,但仍小于耕地转出面积;而林地和水域面积减少62.63和88.79 km2,大部分转入了耕地和建设用地。从空间变化来看,沿海区域与山区发展不平衡问题较为突出,随着沿海城市建设速度的加快,城市扩张和工业发展会占用大量的耕地等生产用地,并割裂生态用地,该模拟结果符合区域对土地利用/覆盖发展方向的基本判断。
表8 福州都市圈2020—2030年各主体功能区土地利用/覆盖类型面积及比例变化量
图6 福州都市圈2030年土地利用/覆盖模拟结果
图7 福州都市圈重点开发区域2020—2030年土地利用/覆盖转移情况(km2)
图8 福州都市圈2020—2030年建设用地扩张空间分布
表7 福州都市圈各县区建设用地扩张强度
优化开发区域中建设用地扩张面积最大,达17.75 km2,所占比例增长1.76%,多集中于福州市仓山区周围;耕地约有16.71 km2减少,所占比例减少1.66%,其中15.60 km2转入建设用地;林地有耕地和草地的转入,面积增加2.56 km2;而草地和水域有不大于3 km2的小幅减少。限制开发区域中各地类面积变化都小于3 km2,且变化分布范围较为零散。建设用地的发展速度有所控制,耕地和草地仍有小面积的减少,且多转向林地,林地增长2.98 km2。禁止开发区域中各地类面积变化幅度较小,均小于1.5 km2。其中建设用地向草地、水域等其他地类转化,面积有小幅缩减;林地和水域有草地的转入,各自有1 km2左右的面积增加,资源承载能力增强。
FLUS-Markov模型中参数设置是提高模型模拟精度的关键。本文考虑区域发展差异性,提出一种基于不同主体功能区的分区模拟预测方案,从而实现研究区未来国土空间格局模拟,增加了模型对大尺度区域模拟的适用性,但在模拟过程中仍然存在一些不足。
(1)在邻域因子权重的确定方面,有学者采用经验赋值法或基于历史面积变化量确定模型参数等方法[10]。本文考虑研究区域各地类初始面积差异性的影响,提出一种客观判断的参数设置方法,将各土地利用/覆盖类型历史面积变化量的相对无量纲值表示地类的扩张强度,并运用于FLUS-Markov模型中模拟2020年土地利用/覆盖空间分布,模拟精度kappa系数均在0.85以上,说明该方法适用于该研究区域未来土地利用/覆盖模拟,但是否具有普适性仍需要更多区域的实践证明。
(2)在驱动因子的设置方面,土地利用/覆盖变化受自然环境、人类社会、政策规划、地域文化等多种驱动因子的复杂影响。但由于部分数据难以量化获取,导致驱动因子的设置并不全面。因此,在未来研究中应全面、综合考虑驱动因子复杂性和代表性,可以更精准地预测土地利用/覆盖空间分布格局。
(3)在限制区域的设置方面,区域内包含大面积耕地保护区和自然保护区等关系到区域粮食安全、生态安全等区域。若以永久基本农田和各级自然保护区为刚性底线融入国土空间模拟预测中,可以得到更切合的模拟结果。
(1)在研究区域内基于历史情景的各土地利用/覆盖类型历史面积变化量的相对无量纲值能很好代表各地类的扩张强度,可以作为邻域因子权重参数设置。因此,将该参数结合考虑空间异质性特征的适宜性概率分布计算和转换成本矩阵等参数的设置,经精度检验,各主体功能区模拟结果总体精度均在90%以上。这表明在这种参数组合下,FLUS-Markov模型适用于对未来土地利用/覆盖进行模拟。
(2)利用土地利用动态度和土地利用程度等指标对分区和整体模拟结果进行量化对比,发现分区模拟结果与主体功能区划的发展目标更加吻合。分区模拟结果表示,2030年各主体功能区均有耕地、草地的减少;重点、优化发展区域面积变化量最大的为建设用地,分别增加了320.24和17.75 km2,多分布于沿海城市,且大部分来源于耕地转入;限制、禁止开发区域各土地利用/覆盖类型面积变化幅度不大,均小于3 km2。
(3)福州都市圈内分区模拟结果基本符合主体功能区划土地利用/覆盖未来发展方向,模拟结果可为研究区未来国土空间规划及生态空间管控提供决策参考。