由各类建设用地增加引起的土地利用变化是城市化最为显著的特征之一。可以说, 城市化实质是土地可利用功能的变性过程, 无论是人口集聚、产业结构调整, 还是基础设施建设, 都必须通过土地的重新配置来实现^[1-2]。京津冀城市群作为我国经济发展的核心区域, 在促进全国区域协调发展与保障生态安全格局具有多重战略意义^[3-4], 然而, 由于快速的城市扩张引发的土地利用显著变化, 以及一系列的生态环境问题, 已经严重制约着京津冀的可持续发展。区域尺度人口的大规模聚集增长与土地资源的有限性和可利用性, 以及社会经济需求的增长与资源禀赋之间的矛盾日益加重^[5]。因此, 区域城市发展布局与生态安全多目标的权衡与优化是保障京津冀可持续发展的关键问题。

土地利用结构优化的概念, 最初来源于Plummer于1993年给出的具体表述, 即为达到一定的社会、经济及生态保护等目标, 依据土地适宜性评价和土地资源的自身特性, 对区域内土地资源的各种利用方式进行更加合理的数量安排和空间布局, 从而提高土地利用效率及效益, 维持土地生态系统的相对平衡, 实现土地资源的可持续利用^[6]。简言之, 土地利用结构优化是通过合理的数量安排和空间布局实现土地资源的可持续利用^[7]。

目前, 国内外土地利用结构优化规划中, 优化方法主要有土宜法、综合法和模型法^[8]。常用的数学模型有综合平衡法、线性规划法、系统动力学模型、灰色线性规划法与模糊线性规划法等。多种模型方法各有优劣:例如, 线性规划模型在土地利用结构优化中广泛应用, 较为简单, 但其模型目标单一, 缺乏空间性和时间动态性, 应用范围有限^[9]; 灰色线性模型是研究不确定性决策问题的有效方法, 但仍为单目标的, 缺乏空间分异性, 应用尚不成熟^[9]; 系统动力学模型是一种系统机理分析模型, 要求对模拟的系统内部机制充分了解, 在因果关系不明确的情况下难以与复杂的空间信息相融合^[10]。

近年来, 适应土地可持续利用规划的多目标优化方法由于其既具有系统动力学模型的正负反馈特征, 但比系统动力学方法应用更简单, 同时又克服了线性规划模型的目标函数一维性和可行解区局限性的弊端, 而逐渐为人们所应用。王红瑞等通过建立北京市丰台区土地利用多目标优化模型, 结合相关政策分析设计了土地利用结构调整方案^[11]; 向芸芸和蒙吉军以武汉城市圈为研究区, 以经济效益和生态效益的最大化为目标, 获得了符合城市圈发展需求的土地利用结构优化方案, 并提出了土地资源可持续利用的建议^[12]。同时, 遗传算法已被认为可能是最适合于多目标优化的方法, 而近年来对遗传算法求解多目标优化问题的研究也逐渐成为热点^[11-12]。遗传算法作用于整个种群, 同时强调个体的整合, 并具有处理目标函数的间断性和多峰性等复杂问题的能力, 能有效应用于多目标搜索, 适用于解决复杂的非线性和多维空间寻优问题^[13-15]。因而, 一系列的研究表明, 利用遗传算法的内在并行机制和全局优化特性, 求解多目标土地利用结构优化模型, 从而可以较好地解决土地利用结构的优化配置问题^[16-17]。

然而, 若仅考虑构建多目标优化模型, 并依靠遗传算法求解, 只能从数量结构上对土地利用进行优化配置, 很难从空间分配上进行优化, 从而有其局限性。随着GIS技术的不断发展, 一系列土地利用空间分析模型逐渐进入大众视野, 其中CLUE-S模型可综合考虑空间限制性区域、土地适宜性等特征合理地进行土地利用的空间分配, 应用广泛^[18-21]。将其与多目标土地利用结构优化模型相结合, 可以很好地弥补多目标优化模型空间优化不足的缺陷, 从数量安排和空间布局两个维度对土地利用进行优化配置^[22-30]。

综上所述, 本研究将针对城市群这一重要的区域城市化模式, 以京津冀城市群为例, 一方面, 综合考虑经济、生态以及社会效益, 建立权衡城市扩张、耕地保护与生态效益的多目标优化模型; 另一方面, 结合京津冀城市群景观格局的动态演变过程与驱动机制, 构建基于重要生态保护区发展限制的CLUE-S模型。两者相结合, 提出多目标权衡下的2020年京津冀城市群土地扩张优化情景预案。本研究的方法及结果对京津冀城市群可持续发展的规划与决策制定具有科学的参考价值。

1 研究方法 1.1 研究区概况

京津冀城市群地处环渤海核心地带, 面积约21.4万km², 地势西北高、东南低, 区位优势明显。京津冀城市群包括北京、天津2大直辖市以及河北省的石家庄、唐山、保定、廊坊、张家口、承德、秦皇岛、沧州、衡水、邢台与邯郸11个地级市(图 1)。截至2014年, 作为我国最具影响力的城市群, 京津冀城市群土地面积仅占全国2.3%, 而生产总值为6.65万亿元, 占全国10.5%, 总人口占全国8.1%, 其中, 北京和天津等中心城市的城镇化率分别高达86.4%和82.3%, 高于全国平均水平54.8%(数据来源:北京市统计局)。在社会经济迅速发展的同时, 高强度的城市扩张也对京津冀原本脆弱的生态环境造成了巨大压力, 例如, 水资源短缺、耕地减少、严重的雾霾等问题已成为制约区域可持续发展的重要限制因素。一些地区的生态服务功能明显下降, 例如, 西部太行山、燕山土壤侵蚀和坝上高原荒漠化加剧, 平原洼淀萎缩及消失, 以及地面沉降加剧^[31]。2016年2月, 全国第一个跨省市的区域“十三五”规划印发实施, 明确了京津冀地区未来五年的发展目标, 提出到2020年, 京津冀地区的整体实力将进一步提升, 应在经济保持中高速发展的同时, 使生态环境质量得到明显改善。

1.2 多目标土地利用优化模型 1.2.1 目标函数的构建

目标函数的选取方面, 考虑到京津冀城市群当前土地资源特征与区域发展要求, 本研究以权衡社会、经济与生态三大效益为目标, 主要选取了经济效益和生态效益两个目标, 而社会效益方面, 由于社会效益主要指社会各部门对土地的需求程度, 包括人均建设用地、人均耕地、城镇化水平、人均绿地面积等, 涉及面广, 很难用单一的最大或最小目标函数进行量化, 且社会效益主要体现在相关政策文件对各种土地利用类型的限制或保护, 故而本研究将社会效益转化为约束条件进行处理。

(1) 经济效益目标

(1)

式中, F₁为研究区内的国内生产总值(元), 用以表征研究区经济效益总值; B_i为经济效益系数, 即土地利用类型i所产出的单位面积GDP(元/hm²); x_i为土地利用类型i的规划面积(hm²)。

式中B_i的计算是结合2005-2010年《城市统计年鉴》中的经济数据(林业产值、牧草业产值、渔业产值、农业产值等), 通过表 1中的公式换算得到逐年的经济效益系数, 再利用GM(1, 1)灰色模型预测求得的2020年经济效益系数(其他地类设置为0.01)。灰色预测GM(1, 1)模型是通过部分已知的数据, 将系统信息抽象的概念量化, 量化的概念模型化, 最终进行模型优化来预测部分未知的内容^[32]。经济效益系数的具体含义如表 1。

将经济效益系数代入式(1)可得经济效益目标函数如式(2):

(2)

(2) 生态效益目标

1997年, Costanza等的研究成果使得生态系统价值评估的原理和方法从科学意义上明确, 谢高地等在参照Costanza全球生态系统服务功能评价模型研究的基础上, 修正后制定了中国陆地生态系统单位面积生态服务价值系数表^[33-35], 其中主要考虑的生态系统服务为气候调节、气体调节、水源涵养、土壤形成与保护、废物处理、生物多样性保护、原材料生产和娱乐文化价值。由于目前国内外生态系统服务价值的评估都难以得出让公众和学术界普遍接受并认可的结果, 中国600多位具有生态学背景的专家通过5年调查形成了一个基于专家知识的生态系统服务价值单价体系, 在没有更好更有效的方法时, 可将其用作适合中国生态系统状况和经济发展水平的生态系统服务价值评估工具^[36]。

生态价值目标函数如下

(3)

式中, F₂为研究区内生态系统服务功能价值总量(元), C_i为土地利用类型i相对应的生态价值系数(元/hm²), x_i为土地利用类型i的规划面积(hm²)。

生态系统服务价值系数指生态系统产生的生态服务的相对贡献大小的潜在能力, 定义为1 hm²全国平均产量的农田每年自然粮食产量的经济价值, 以此可将权重因子表转换为当年生态系统服务单价表^[37]。本文采用了谢高地等制定的中国陆地生态系统单位面积生态服务价值系数表, 并通过修正系数^[29]核算得到京津冀各土地利用类型的单位面积生态系统服务价值。系数校正方法为:以2001-2009年京津冀区域平均粮食单产428750 kg/km²为基准, 粮食单价按1.4元/kg(2002-2010年国家退耕还林补助标准)为基准, 确认一个没有人力投入的自然生态系统提供的经济价值量为全国平均粮食单产市场价值的1/7^[37]。校正后的生态价值系数如表 2所示。

将生态效益系数代入式(3)可得生态效益目标函数如式(4):

(4)

1.2.2 约束条件

大量的案例研究都证明了政策对土地利用变化的重要作用, 例如《京津冀协同发展土地利用总体规划》的发布等。本研究主要考虑了土地总面积约束、人口总量约束、土地开发强度约束、湿地面积约束、耕地保有量约束、其他用地开发强度约束、林地保有量约束、草地保有量约束和数学模型约束, 选取的约束指标总体上反映社会效益要求, 以求实现基于政策的对建设用地扩张加以限制以及对林地、湿地、耕地等生态用地数量加以保护。数据来源主要是《北京市土地利用总体规划(2006-2020年)》、《天津市土地利用总体规划(2006-2020年)》、《河北省土地利用总体规划(2006-2020年)》等, 具体设定如下。

(1) 土地总面积约束。各土地利用类型规划面积的总和应等于研究区总面积, 即:

(5)

式中, A为研究区总面积(hm²), 将京津冀城市群总面积代入上式求得, 见式(6)。

(6)

(2) 人口总量约束。按照城镇用地和生态用地平均人口密度预测, 区域内的人口总量应控制在预期年(2020年)规划人口数之内, 从而满足区域内的人口承载力要求, 即:

(7)

式中, D₁、D₂分别为基于规划文件的农用地及城镇用地的2020年预测人口密度(人/hm²); x_i、x_j分别为生态用地和城镇用地面积(hm²), i=1, 2, 3, 4;j=5(其中x₁为林地; x₂为草地; x₃为湿地; x₄为耕地; x₅为人工表面); P为2020年京津冀规划总人口(人), 由《河北省土地利用总体规划(2006-2020年)》、《天津市土地利用总体规划(2006-2020年)》和《京津冀协同发展规划纲要》中获取数据, 可得2020年河北省、天津市、北京市总人口上限分别为7480万、1500万和2300万。将京津冀城市群2020年规划总人口上限赋为P值, 代入数据得式(8)。

(8)

(3) 土地开发强度约束。为防止城市建设用地扩张过快, 对生态环境胁迫加重, 故而利用预期年土地利用规划的建设用地面积加以限制, 设约束条件如下:

(9)

式中, A_{5, 0}为人工表面的现状面积(hm²); A_{5, p}为人工表面2020年规划设定的人工表面面积上限(hm²)。现状面积来源于2010年土地利用分布图, 规划面积来源于《河北省土地利用总体规划(2006-2020年)》、《天津市土地利用总体规划(2006-2020年)》和《北京市土地利用总体规划(2006-2020年)》, 代入数据得式(10):

(10)

(4) 湿地约束。近年来, 京津冀城市群水资源承载力不断下降, 水资源破坏及污染严重, 为保护区域水环境, 设定湿地约束条件为预期年湿地面积不得少于基于历史发展预测得到的预测值, 如下表示:

(11)

式中, A_{3, p}为基于历史演变的湿地预测面积(hm²)(具体方法见1.3.2), 代入数据得式(12):

(12)

(5) 耕地保有量约束。由河北、天津、北京土地利用总体规划(2005-2020)确定耕地保护红线, 结合基于历史发展模拟得到的耕地预测值, 设置约束条件如下:

(13)

式中, A_{4, p}为《河北省土地利用总体规划(2006-2020年)》、《天津市土地利用总体规划(2006-2020年)》和《北京市土地利用总体规划(2006-2020年)》中2020年耕地保有量总和(hm²)。A_{4, p}′为基于历史演变的耕地预测面积(hm²)(具体方法见1.3.2), 代入数据得式(14):

(14)

(6) 其他用地开发强度约束, 如下式:

(15)

式中, r为土地开发利用率, 有关京津冀的土地开发利用率暂未找到标准, 参考其他城市设定为92%, 代入数据得式(16):

(16)

(7) 林地保有量约束, 如下式:

(17)

式中, t为森林覆盖率。日前国家林业局分别与北京、天津、河北三省市人民政府签署《共同推进京津冀协同发展林业生态率先突破框架协议》中提出, 2020年京津冀城市群森林覆盖率应达35%。即t=0.35。代入数据得式(18):

(18)

(8) 草地保有量约束。考虑到京津冀城市群生态环境质量, 设定预期年草地面积不低于现值, 如下:

(19)

式中, A_{2, 0}为草地的现状面积(hm²)(具体见1.3.1), 代入数据得式(20):

(20)

(9) 数学模型约束。各变量非负, 即:

(21)

1.2.3 遗传算法

遗传算法从一组随机产生的初始解(种群)开始搜索过程, 种群中每个个体, 即一个随机变量, 称为染色体。这些染色体在迭代中通过交叉、变异、选择而不断进化, 产生下一代染色体, 即后代。后代中染色体是否可以作为下一代种群参与进化通过适应度来度量。通过若干迭代后, 算法大概率收敛于最优染色体, 即最优解或非劣解^[13]。

由于多目标优化问题往往需要的是整体最优, 而非单个最优, 且多个目标之间可能是相互竞争的关系, 所以需要求出一组解决方案的集合, 称为Pareto集。本研究采用MATLAB R2012b软件中的多目标遗传算法(MOGA)工具箱对土地利用结构多目标优化模型进行求解。

1.3 土地利用空间优化配置方法 1.3.1 数据来源

京津冀城市群的土地利用现状数据采用的是“全国生态环境十年变化(2000-2010年)遥感调查评估”专项的全国土地利用分类体系的一级分类结果, 主要包括林地、草地、湿地、耕地、人工表面与其他6大类地物^[38]。土地利用特征的提取主要基于Landsat-5 TM遥感影像, 利用面向对象和回溯相结合的土地利用/覆被信息提取方法, 获取得到2000年、2010年的土地利用分布图。并结合高分辨率遥感数据和野外实地调查数据对遥感解译数据进行修正及验证, 以保证分类精度^[38-40]。

1.3.2 CLUE-S模型

已有研究表明, 模型是模拟与预测城市时空发展格局的有效方法, 国内外很多学者通过构建模型的方法对不同类型的、单一的城市扩张演变过程进行了模拟与预测。其中, CLUE-S(the Conversion of Land Use and its Effects at Small regional extent)模型是荷兰瓦赫宁根大学“土地利用变化和影响”研究小组在CLUE模型的基础上开发的^[24]土地利用模拟模型, 目前已在国内外得到广泛应用^[21-23]。CLUE-S模型有两个基本假设, 即(1)某地区的土地利用变化是受该地区的土地利用需求所驱动的, 并且(2)一个地区的土地利用分布格局总是与土地需求及该地区的自然环境和社会经济状况处于动态平衡状态^[28]。以这两个假设为前提, CLUE-S模型的核心原理是基于Logistic回归的区位适宜性, 并运用系统论的方法模拟不同土地覆盖类型随时空变化过程中的竞争关系和相互作用的经验性分析。

本研究利用CLUE-S模型将多目标优化结果分配至空间上, 从而更直观地体现土地利用优化配置在京津冀城市群的表现形式^[21]。

1.3.3 综合生态系统服务功能的保护

已有研究表明, 生物多样性、防风固沙、土壤保持和水源涵养这4种生态服务功能对区域的生态安全具有重要影响。研究通过建立生态服务综合指标来量化空间上每个栅格单元对维持区域生态安全的相对重要性。方法详见环境保护部与中国科学院“全国生态环境十年变化(2000-2010年)调查评估”中应用的区域生态系统服务的定量评价方法与结果^[31-33]。

由于京津冀城市群基于历史条件演变的未来城市扩张将对生态系统服务功能产生严重的侵占, 不利于区域生态安全^[41], 故而本研究将综合生态系统服务功能作为限制发展区域输入CLUE-S模型, 以实现未来优化配置中从空间层面对生态环境的保护。

2 研究结果与讨论 2.1 土地利用结构优化方案解集

一般情况下, 多目标优化问题中的多个目标不可能同时达到最优, 甚至可能是相互冲突的, 而由于决策者自身经验及偏好的不同, 可能会得到不同的最优解, 而在多目标优化中这些可能的最优解称为Pareto非劣解, 利用多目标遗传算法可以同时得到多个Pareto非劣解, 从而提供更多的选择空间。

本研究通过MATLAB中的的MUGA模块进行迭代运算, 直到算法运行满足停止准则, 本研究的终止准则为Pareto解的展开速率低于默认值。图 2为MOGA搜索到的Pareto非劣解前端, 非劣解前端大体呈曲线分布, 说明大多数非劣解都可以搜索得到。

研究对表中的30个Pareto非劣解方案进行比较分析, 可以发现各个方案均能够得到较好的综合效益, 且相差并不大。本研究从中选择2种典型方案(表 3)进行下一步分析研究, 其中, 方案1的湿地面积是所有方案中最大的, 这主要也是由于湿地在生态效益目标函数中所占比重较大的原因; 方案2中耕地和人工表面面积都相对较大, 既可以保证一定的区域经济发展, 又可以在一定程度上实现耕地保护, 是所有方案中相对最优的方案, 故而将方案2作为土地利用数量结构优化方案, 尝试优化。

2.2 基于多目标的土地利用优化配置情景模拟

将多目标土地利用结构优化结果(方案2)输入CLUE-S模型, 同时考虑综合生态系统服务功能的限制开发, 结合个土地利用类型的历史转移规则, 预测模拟得到京津冀城市群2020年土地利用格局优化情景预案, 与对应的基于历史条件演变的未做优化的土地利用格局情景(图 3)相比^[21], 优化后的土地利用格局变化情景(图 3)存在明显差异。

总体来看, 相较于基于历史条件下的2020年土地利用格局变化情景, 由于2020年《共同推进京津冀协同发展林业生态率先突破框架协议》中提出的2020年京津冀城市群森林覆盖率35%的要求, 故而优化后草地增加速率相对减缓, 同时林地面积大幅增加, 未优化情景中林地面积相较于2010年增加2.37%, 而优化后林地面积增加4.12%, 提高了区域整体的森林覆盖率。湿地面积相较于未优化情景衰退速率有小幅降低, 但整体上变化不大。由于本优化方案所考虑的主要因素为耕地保护和城市扩张相协调, 故而耕地面积相较原方案降幅减小, 基于历史条件演变的原方案中耕地面积较2010年下降5.56%, 而优化方案中下降4.40%, 符合耕地保有量要求。人工表面变化方面, 由于耕地保护和林地保有量增加的双方要求, 故而人工表面增加比例大幅下降, 未优化情景中, 人工表面较2010年增加17.02%, 而优化后的人工表面增加比例仅为7.71%, 由此可以看出, 优化方案中区域经济发展为生态环境保护做出一定让步, 同时, 由于经济集约化发展以及产业结构调整, 增速放缓后的人工表面面积依然有可能符合京津冀城市群稳定、可持续发展的要求(表 4)。

3 结论

本文通过构建权衡京津冀城市群社会、经济发展与生态效益的多目标土地利用结构优化模型, 并利用遗传算法求解, 得到土地利用结构优化配置的情景预案, 从而对该区域内的土地利用结构进行优化, 最终结合空间显性CLUE-S模型, 模拟得到未来(2020年)的京津冀土地利用结构的空间优化配置情景预案。主要结论如下:

(1) 利用多目标遗传算法(MOGA)对满足社会、经济与生态效益等不同目标的土地利用结构优化配置问题进行求解, 可以得到一组使综合效益达到较优水平的Pareto非劣解方案, 从一定程度上为决策人预案提供多种选择方案, 并通过与空间显性模型CLUE-S相结合, 从而实现对权衡不同发展目标的土地利用结构及格局分布的情景模拟。

(2) 优化方案与基于历史条件演变的原方案相比, 林地的增加可以相对满足政策要求, 增幅较大, 但同时草地面积增长放缓。湿地面积与原方案相差不大, 这也是由该方案的耕地保护倾向性决定的, 但其在空间分配上相对更为合理。耕地较原方案减幅放缓, 符合耕地保有量的要求, 同时, 人工表面为区域生态环境保护做出一定让步, 但结合集约化经济的推行以及产业结构的调整, 仍能满足社会经济发展需求。

(3) 由于城市群的范围较大, 影响其发展的因素很多, 不同因素的数据可获得性、定量化、空间化, 以及匹配度存在差异, 对土地覆盖/利用空间配置的模拟结果会有明显影响, 也是实施本研究所采用方法的局限性。如何更加全面、精细地定量化与空间化影响城市群发展的关键因素, 有待进一步补充完善。