基于土地利用的淮海经济区碳储量估算与预测

解天琪 李 龙 陈 鑫 周 翼 张国梁 陈龙乾

(1.中国农业大学 土地科学与技术学院，北京 100193；2.中国矿业大学 环境与测绘学院，江苏 徐州 221116；3.荷语布鲁塞尔自由大学 地理系，比利时 布鲁塞尔 1050)

土地利用变化能够改变生物地理群落的结构和功能，从而影响碳循环过程，因此成为造成陆地生态系统碳储量变化关键要素之一[1-2]。除了燃烧化石燃料之外，土地利用变化对大气中CO2浓度增加的影响最大[3-4]。由于人类活动对土地利用的干扰加剧，我国土地利用模式发生了较大变动，使得陆地生态系统的结构和功能发生变化，进而影响碳储量，引发了重大环境效应。因此，从土地利用的角度研究碳排放，不仅可促进生态保护，还有利于耕地合理利用[5]。淮海经济区是我国主要农副产品基地与粮仓，近年来由于该地区人口不断增长，经济快速发展，对土地资源尤其是建设用地资源的需求不断增加[6-7]，导致了过度的土地非农化、城市“摊大饼”蔓延等问题，影响到区域碳收支的平衡，而碳收支平衡对一个区域的生态水平起着至关重要的作用[8]。因此，进行淮海经济区土地利用时空变化特征分析及其碳储量测算研究，对该区域土地的合理利用十分必要[9]。

早在20世纪50年代国际上就有研究估算了全球土壤有机碳库的储量[10]。Rubey等[11]、Bohn等[12-13]、Post等[14]随后分别计算了全球土壤有机碳库储量。90年代以来，不同研究结果差别较小[15-18]。国内对土壤碳储量的研究相对较晚，从20世纪80年代开始，对土壤有机碳库的研究才日益增多，基于不同空间尺度对土壤碳储量进行了系列的估算与分析[19-29]。关于植被碳储量的计算，初期主要通过野外实地调查生物量和面积统计资料，但在野外调查时易受调查人员主观判断的影响，从而导致较大误差[30]；近期的研究普遍采用建立生物量和蓄积量间关系的方法来估算植被碳储量[31]。

已有研究侧重于历年碳汇数据预测碳储量，往往忽略了外部因素的变化对碳储量的影响[32]。因此，为掌握淮海经济区碳时空特征，本研究拟从土地利用的观点出发，基于遥感技术的应用和土地利用数据，分析2006、2011和2017年3 个时期土地利用变化，并依据CA-Markov模型模拟预测淮海经济区2025年的土地利用分类情况，进而对因其变化引起的碳储量变化进行分析，以期预测未来土地利用变化造成的碳收支变化，也为政策制定提供定量依据。

1 研究区概况与数据来源

1.1 研究区概况

淮海经济区位于32°22′43″～36°33′42″ N，113°51′42″～120°55′44″ E，包括苏豫鲁皖4 个省20 个地级市147个县级单位(2019年1月，国务院批复同意撤销地级莱芜市，将其所辖区域划归济南市管辖，但由于本研究所涉及的数据未能根据国家相关规定实时更新，故仍将莱芜市作为地级市进行研究)，总面积达1.781×105km2，占全国土地总面积的1.86%(图1)。研究区气候以半湿润季风气候为主，温和湿润，雨热同期，年降雨量为750～1 050 mm，年均气温在14 ℃左右，光热条件好；植被以落叶阔叶林为主，土壤肥沃；主要作物包括小麦、稻谷、棉花、花生等，作物熟制为一年两熟，轮作方式主要包括小麦-水稻、小麦-玉米等。

1.2 数据来源

根据研究目的及研究区范围，本研究中的遥感数据选择2006、2011和2017年3 个时期MCD12Q1数据，下载于美国地质调查局网站(https:∥earthexplorer.usgs.gov)。

社会经济数据来源于2007、2012、2018年苏豫鲁皖统计年鉴；DEM数字高程数据从Google Earth网页提取；碳密度数据来源于国家生态系统观测研究网络(http:www.cnern.org.cn)和参考文献[23,33-37]的研究结果；道路矢量数据使用水经注万能地图下载器获取。

2 研究方法

2.1 土地利用分类确定方法

本研究结合《土地利用现状分类》(GB/T 21010—2007)和淮海经济区实际情况，将研究区分为水域、林地、草地、耕地、建设用地、园地、裸地和其他用地等8 类(表1)。

图1 淮海经济区区位图Fig.1 Location map of the Huaihai Economic Zone

表1 土地利用分类体系Table 1 Land use classification system

2.2 预测模型创建及精度模拟验证

20世纪50年代初，现代计算机创始人Neumann提出标准元胞自动机(cellular automaton，CA)的雏形，其模型原理：元胞自动机的空间、时间和状态均是离散的，每个元胞也都处于离散状态；各个元胞状态按照相同的转变规则进行同步更新，每个变量状态是有限的，且仅在局部遵循转变规则[38]。Markov模型是基于概率创建的一种随机型的时序模型，根据某一变量的目前形态和趋向预测其未来状态的方法，因此，Markov模型可根据不同时期土地利用分类情况计算得出的转置机率来预测未来某一时期的土地利用情况[30]。然而，Markov模型只能预测土地覆盖分类数量方面的变化，不能预测其在空间方面的变化，即无法预测土地利用类型空间上的分布，元胞自动机恰好具备模拟空间演变这一功能，因此将二者结合组成CA-Markov模型，能实现在数量和空间上对土地利用变化进行模拟预测。

根据CA-Markov模型，结合前期研究区不同时期遥感图像分类结果,计算得到未来某时期土地利用分类情况。本研究通过计算2017年土地利用现状栅格数量和模拟预测栅格数量的误差绝对值来验证模拟精度，计算公式如下[39]：

(1)

式中：E为第i种土地利用类型的数量的差值；Mix、Miy分别为第i种土地利用类型的实际栅格数量和模拟预测栅格数量，个。需要注意的是运用此方法比较的是E的绝对值，其绝对值越小表示模拟精度越高。

本研究所用的CA-Markov预测流程如图2。

图2 CA-Markov预测流程图Fig.2 CA-Markov prediction flow chart

2.3 碳密度数据选定及碳储量估算

考虑到碳密度计算所需的数据的可得性，本研究中的碳密度数据参考徐钊等[33]、贾松伟[34-35]、李伟等[36]、揣小伟等[37]已发表的研究结果；土壤碳密度计算，直接参考揣小伟等[37]研究结果，该研究有区域有5 个在淮海经济区；植被碳密度的计算根据生物量碳密度估算法，即基于生物量推算出不同植被碳密度，通过对淮海经济区所涉及的苏、豫、鲁、皖四省的主要作物面积及产量进行统计分析发现其结构相似，因此，河南、山东、安徽各土地利用类型碳密度可参考江苏省的研究结果[37]，各地类碳密度数据见表2。

利用碳密度计算不同地类的土壤碳储量，计算公式如下[40]:

TOC=SOC+ZOC

(2)

(3)

(4)

式中：TOC为淮海经济区总碳储量，Tg；SOC为土壤总碳储量，Tg；Xi为第i种土地利用类型的面积，km2；SOCDi为第i种土地利用类型的土壤碳密度，kg/m2；ZOC为植被总碳储量，Tg；ZOCDi为第i种土地利用类型的植被碳密度，kg/m2。

表2 淮海经济区不同地类碳密度Table 2 Carbon density of different land types in the Huaihai Economic Zone kg/m2

3 结果与分析

3.1 土地利用变化分析

研究发现淮海经济区土地利用类型以耕地为主，在2006、2011、2017年所占比例分别为87.51%、87.10%和86.21%，呈现出缓慢减少的趋势；建设用地面积持续增加，由2006年的9 697.94 km2增加到2017年的10 962.10 km2；草地不断扩张，占地面积由5 786.50 km2扩张到6 177.75 km2；水域、林地和园地面积所占比例较少，呈缓慢增加趋势；在此期间，其他土地利用类型面积不断减少，2017年总面积仅为142.94 km2，所占比例不足1%。

为显示不同土地类型之间的转化的细节，利用Origin2020制作桑基图(Sankey diagram)(图3)。由图3可知：研究区内主要是耕地、草地、建设用地和园地之间的转化；草地的转出面积与转入面积都是最大的，转入面积大于转出面积导致草地面积总体上有所增加，说明草地主要是空间分布的变化；从转出方向看，草地转为耕地的面积最大，占草地的31.74%，水域次之；从转入方向看，耕地和水域是草地面积增加的主要原因，其中9.86%的耕地转换为草地；建设用地转出方向较多且面积相当，但由于2.49%耕地的转入，建设用地面积在总量上有所增加，表明了淮海经济区正处于经济发展阶段，同时也映射出在人口持续增长的压力下对建设用地的需求不断增加；裸地面积有所增加，原因在于研究区矿产资源丰富，煤、石油、天热气等储量丰富，开采量不断增加，造成区域内裸地面积不断增加；水域呈现持续扩张的趋势，其原因是淮海经济区对生态环境保护的日益重视。

图3 2006—2017年土地利用面积变化转移Fig.3 Change and transfer of land use area from 2006 to 2017

3.2 土地利用模拟

3.2.1CA-Markov模型模拟精度验证

首先对2017年土地利用变化进行模拟，并将其与2017年土地利用现状进行比较以验证模拟精度(表3)。

由表3可知，除了其他用地之外，各个地类面积预测误差均小于5%，表明该模拟的准确程度较高，能够较为精确的预测未来某一时间的土地利用情况。其他用地误差较大的原因可能是其所占面积较小，2017年底其他用地面积仅占研究区域总面积的0.08%，基线值较低，且呈递减趋势，在后续2025年模拟预测中可忽略不计。

表3 2017年土地利用模拟精度Table 3 Land use simulation accuracy in 2017

3.2.22025年土地利用模拟

利用CA-Markov模块制作淮海经济区2025年土地利用类型的适宜性图集从而得到2025年土地利用预测数据(图4)。

对预测结果进行统计，2025年淮海经济区水域、林地、草地、耕地、建设用地、园地、裸地、其他用地面积分别为2 019.30、195.09、6 209.72、151 958.73、11 863.31、1 566.64、4 804.71、87.04 km2，分别占研究区总面积的1.13%、0.11%、3.47%、85.02%、6.66%、0.88%、2.69%、0.05%。

3.3 碳储量时空变化分析

3.3.1碳储量测算

按照2.3土地利用数据和碳储量计算公式，对淮海经济区2006、2011、2017年碳储量计算结果如表4所示，对2025年碳储量预测结果如表5所示。

3.3.2碳储量时空演变特征

1) 碳密度和碳储量时间变化特征

利用Origin绘图功能制作2006—2025年研究区各土地利用类型碳密度和碳储量随着时间的变化情况，如图5所示。

图4 2025年土地利用预测数据Fig.4 Land use forecast data for 2025

表4 各土地利用类型碳储量Table 4 Carbon density and carbon storage of land use types Tg

表5 2025年各土地利用类型预测碳储量Table 5 Carbon density and carbon storage of land use types in 2025 Tg

由图5可见：除耕地和其他用地碳储量在研究期呈递减趋势外，其余地类呈不同程度递增趋势。2006—2025年，林地碳密度最大，变化也最为明显，总体上呈匀速增加趋势，结合这一期间林地面积数据得出林地碳密度变化主要原因在于林地面积的不断扩张；水域、草地、耕地、建设用地碳密度几乎没有变化；园地仅次于林地的碳密度，总体上呈缓慢增加趋势。除耕地外，建设用地碳储量最大且在整个研究时期内总体呈现不断递增趋势；草地在2006—2011年内增长幅度较为明显，在2011—2025年期间内增长速度有所减缓；耕地碳储量总体呈现递减趋势，但减少速度极慢，是因为耕地面积在一定范围内有缩减；水域、林地、园地碳储量变化较小；裸地在2006—2011年碳储量增长幅度较大，在2011—2025年内仍然呈现递增趋势，但增长速度大大减小；其他用地碳储量明显减少，原因在于其面积不断被其他土地利用类型侵占。

2) 碳密度和碳储量空间变化特征

基于ArcGIS绘制整个研究期间淮海经济区碳密度、碳储量空间分布情况(图6)。由图(a)可知：在碳密度方面，2006—2017年碳密度增加区域在研究区内零散分布，其变化范围大多在0.99～3.97 kg/m2，原因在于草地向林地的转换；在2017—2025年，碳密度变化总体上呈现出“多零散，少集聚”的分布，究其原因是济宁市2018年开展实施建设国家森林城市，林地面积大幅度增加，和枣庄市、徐州市接壤处碳密度变化较大并且集中。由图(b)可知：在碳储量方面，2006—2011年，碳储量总体上变化不明显，临沂南部增加0～9.03 Tg，济宁市东南部、宿迁市和淮安市交界处碳储量有所减少；2011—2017年碳储量总体上减少，变化值处在9.08～107.26 Tg，这是由于城市化进程中建设用地需求不断增加，侵占其他地类，导致该区域内碳密度降低从而碳储量减少；2017—2025年，碳储量延续前期减少趋势，但济宁市与枣庄市、徐州市交界处碳储量基本不变，原因在于该区域碳密度有所增加。

在上述结果的基础上，利用ArcGIS的空间分析功能的标准差椭圆对2006—2025年碳储量方向分布变化进行分析，生成椭圆面积表示范围，椭圆扁率表示方向性的明显程度。结果表示：整个研究期间，生成的椭圆集中在淮海经济区东部，表示碳储量变化多集中在研究区中心城市徐州及其周边地区；椭圆面积呈现出先增后减的趋势，究其原因在于2011—2017年淮海经济区经济高速发展，城镇化速度不断加快，地类之间的相互转换涉及范围较广，导致碳储量变化区域增大；2006—2025年椭圆扁率不断减小，意味着碳储量变化方向性越来越模糊，呈现出由枣庄、徐州、宿迁和淮安一带向四周发散的趋势。

4 结论与讨论

本研究以土地利用为出发点，模拟计算了淮海经济区历史时期的碳储量变化，精度较高，并以此预测了未来一段时期的碳储量变化，对其变化情况进行分析，得到的主要研究结论如下：

1)分析土地利用历史时空变化可知，耕地为淮海经济区主要的土地利用类型，应对耕地进行二次划分以提高研究精度。研究发现耕地非农化趋势明显，究其原因是为了追求最大利润，发展乡镇企业与“开发区热”，耕地不断被占用导致[41]；建设用地不断扩张，原因在于淮海经济区各市经济水平不断提高，城镇化进程不断发展，大规模开发，导致了建设用地面积不断增加；水域、林地、草地、园地、裸地面积所占比例较少，总体上呈递增趋势；其他用地面积不断减少。

图5 各土地利用类型碳密度变化和碳储量时间变化Fig.5 Time variation of carbon density and carbon storage by land use type

图6 各土地利用类型碳密度变化和碳储量空间变化Fig.6 Spatial changes in carbon density and carbon storage by land use type

2)利用CA-Markov模型模拟2025年淮海经济区土地利用取得了较好的效果；2017—2025年期间土地利用变化趋势与2006—2017年保持一致：耕地持续非农化；建设用地持续扩张；水域、林地、草地、园地、裸地总体上仍呈递增趋势；其他用地面积持续减少。

3)土地利用碳储量测算结果显示：时序上，2006、2011和2017年碳储量分别为1 720.018、1 720.020、1 716.531 Tg，预计到2025年碳储量会达到1 720.726 Tg，较2017年增加4.195 Tg。由于耕地有利于减缓碳密度和碳储量的下降，虽然研究期间耕地面积不断减少但其仍占很大比例；空间上，碳储量变化总体上呈现集聚特征，其减少区域成片分布但其减少量较小，增加区域分布零散但增加量较大。

因此，淮海经济区土地利用方式的转变导致碳储量的变化，耕地在一定程度上减缓了碳储量的减少，而且研究区碳储量随着林地面积的增加而增加。因此，应在切实保护耕地的基础上，积极实施退耕还林，保护林地资源，既保证了粮食安全又有利于提高研究区生态环境质量。