0 引言
【研究意义】生态系统固碳服务（以下简称“固碳服务”）通过捕获大气中的碳并把捕获到的碳固定起来,它可以抵消人类向大气中排放的一部分二氧化碳,从而起到调节气候的作用,是生态系统服务中调节服务的重要组成部分[1]。【前人研究进展】关于生态系统服务功能的受益人群、生态补偿等研究一直都广受研究者关注[2-3],而生态服务供给与人类需求的相互关系及其空间平衡问题近年来逐渐引起国内外专家的重视[4-8]。服务供给与需求的空间匹配问题是梳理生态系统服务空间流动过程研究的主要障碍之一[9-10]。目前得到较多****认同的空间匹配模型是服务路径属性网络（SPANs）模型[11-13],但其研究主要还处于概念阶段[14-15]。从空间流动的角度来看,固碳服务供需的研究主要可分为典型生态系统与大气之间的固碳服务供给（以下简称“碳供给”）研究[16-17]和区域之间的固碳服务需求（以下简称“碳需求”）研究[18]。这些研究促进了人们对于碳供需的理解,但对于以人类活动和自然环境共同组成的复杂系统而言,单纯的生态碳循环研究或经济学碳供需研究不能全面反映区域中固碳服务的供需平衡[19-20]。近年来该研究领域逐渐突破单一视角的局限,开始从人类福祉的角度出发,开展区域固碳服务供需平衡研究[21],一定程度上填补了区域碳平衡研究的空白,但对碳管理实践而言仍有一定局限性。【本研究切入点】作为影响地表植被覆盖和人们能源消费水平的重要因素,科学合理的土地利用格局配置无疑是实现有效的区域碳管理、引导区域低碳发展的重要工具[22]。国内外****已经在土地利用的碳排放效应、城市化进程中的碳收支等领域做了大量的实践研究,但多侧重于土地类型面积的线性组合和土地利用预测,对区域的空间结构优化研究仍需继续探索。【拟解决的关键问题】本文以关中-天水经济区为研究区域,估算研究区碳需求与供给,量化固碳服务供需平衡状况,模拟区域内固碳服务的流动状况,分析环境变量对碳供给的影响,探讨区域中碳供给的不确定性,并给出空间布局优化策略,促进对固碳服务空间流动的客观认识,为区域碳管理提供科学直观的理论依据和数据支撑。
1 材料与方法
1.1 研究区域
关中-天水经济区（简称“关天经济区”）（图1）的经度范围是104°34′47″E—110°48′38″E,纬度范围是33°21′37″N—35°51′15″N,面积约8.01×104 km2,年均气温6—13℃,年降水量500—800 mm,海拔200—3 700 m。根据地形地貌和气候差异可将研究区分为三大部分：秦岭山脉、关中平原和北山山系。
显示原图|下载原图ZIP|生成PPT图1研究区行政区划图
-->Fig. 1Administrative map of study area
-->
作为陕西省重点生态保护区和南水北调重要水源涵养基地,秦岭的大部分区域都属于限制开发区,生态基础优异,是关中地区的主要林地和天然碳库。其土壤类型以棕壤、黄棕壤、褐土为主,其间散布山地棕壤、暗棕壤、粗骨土和
1.2 数据来源
研究选用的数据源主要包括遥感数据、地理基础数据、气象数据、土地利用专题地图和Harmonized World Soil Database（HWSD）。其中,遥感数据为2014年250 m空间分辨率的MODIS数据产品。地理基础数据包括DEM、坡度、坡向、土壤类型等。气象数据来源于中国气象科学数据共享服务网,包括年平均降水、年平均气温、7月平均气温、地表太阳辐射等。HWSD为本文提供中国1∶100万土壤有机碳含量分布图。土壤有机碳实测数据为2015年和2016年8月分别于研究区内野外采样并通过重铬酸钾法测定所得。1.3 固碳服务供需模型
1.3.1 碳需求模型 碳排放的测算方法主要有实测法、物料衡算法和排放系数法,其中排放系数法（Emission-Factor Approach）是IPCC提出的、适用于多尺度的碳排放估算方法,但单纯的碳排放测算难以体现碳排放的空间差异。为了更加直观地表现区域固碳服务需求的空间分布,本文基于能源统计中人均碳排放数据和人口密度分布计算区域碳源的空间分布。\[C_e=\sum
X_{x=1}\rho(x)\times \varphi(x) \ \ (1)\]
式中,Ce 为人类社会经济活动中的碳排放,即固碳服务需求;ρ(x)为像元 x的人口空间密度;φ(x)为像元x的人均碳排放; X 表示研究区的像元总数。
应用较广的人口数据空间化的主要模型有负指数模型及其改进模型、核密度估计模型、分区密度模型、多元回归模型、多因素融合模型和智能化模型[23]。为了综合考虑自然、经济因素对人口分布的影响程度,本文选用多因素融合模型,基于遥感影像解译结果和居民点分布,将区域分为城镇区域和农村区域;基于多种基础数据分别计算城镇人口密度系数和农村人口密度系数;结合市县人口统计数据和人口密度系数计算出区域人口密度分布状况。基于城镇面积的城镇人口密度系数模型为[24]：
Vij=Aj╳ln Ai╳exp\(\bigg[-(\frac{r_j}{\sqrt{A_i/\pi}})
{sigma}\bigg]\)=Aj╳ln A_i[-1.9874rj1.2Ai-0.6] （2）
式中,Vij为第i个城市中第j个栅格的人口密度系数,Aj为j栅格中的城市用地面积,Ai为第i个城市的面积,rj为第j个栅格中心至城市中心的距离,参数σ反映城市的不同发展阶段。而农村人口密度系数则是根据调研所获得的数据和公式（3）求Drural与Drest比值后,经单要素加权融合、多元逐步回归分析、归一化处理等后得出[25]。
Ptotal-Ptown=Arural╳Drural+Arest+Drest （3）
式中,Ptotal、Ptown、Arural、Drural、Arest、Drest分别表示各县总人口、城镇人口、农村居民点的总面积、农村居民点平均人口密度、远离居民点地区的总面积和平均人口密度。
1.3.2 碳供给模型 碳供给部分包括生物固碳和土壤固碳两部分。考虑到模型复杂度、估算精度、所需数据的可获得性[26]及研究区域实际情况,本研究选用CASA模型计算生物部分的碳库存量,选用由周涛等[27]改进的碳循环过程模型来反演土壤基础呼吸进而估算土壤部分的碳库存量。
CASA模型[28]通过光合有效辐射（APAR）和实际光能利用率ε(x,t)的乘积来表示每个月的植被初级生
产力的变化,计算公式如下：
NPP=APAR(x,t)╳ε(x,t) （4）
式中,NPP为植物净初级生产力,APAR(x,t)表示像元x在t月份吸收的光合有效辐射（MJ·m-2）,ε(x,t)表示像元x在t月份的实际光能利用率（gc·MJ-1）。根据不同植被地上生物量与地下生物量的经验比例值,最终计算区域生物碳库的碳供给量。
通过碳循环过程模型来反演土壤基础呼吸,建立土壤中CO2排放与土壤有机碳（SOC）之间的关系模型。该模型在原碳循环过程的基础上增加了水分因子,用年降水量和年潜在蒸散量结合温度敏感因子描述土壤水分对土壤呼吸的影响。其公式如下：
\[A_ij=\frac{NPP}{exp(b\times T)\times y}\ \ (5)\]
式中,Aij为土壤基础呼吸,NPP为净初级生产力,T为年平均温度,b为温度敏感常数因子,y为水分对土壤呼吸的限制因子。
由于HWSD中的土壤数据来自全国第二次土壤普查,与现有的实际数据可能存在一定的偏差。因此本文先计算了2000年的土壤基础呼吸,将其与HWSD土壤有机碳数据做回归分析,得出土壤基础呼吸与土壤有机碳含量的回归模型;将2014年的土壤基础呼吸代入回归模型,得出2014年的土壤有机碳含量。使用野外调研实测数据验证模拟精度,标准误差为4.9726,研究结果具有较好可信度。
1.4 固碳服务的空间流动模型
从生态学的角度,生态学家们关于“流”的概念最早可以追溯到将生命的代谢过程看作是能量、物质与周围环境不断发生交换的过程,即能量流和物质流[29-30]。生态系统中的“能量流”通常用“林德曼定律”来表示[31],但关于碳的空间流动目前尚没有成熟的理论和方法。“场”在物理学中是指物体在空间中的分布情况,可以用来描述某种物理量在空间的分布和变化规律,用空间位置函数来表征。若物理量是标量,那么空间上每一点都对应着该物理量的一个确定数值,则称此空间为标量场（如电势场、温度场）;若物理量是矢量,那么空间每一点都存在着它的大小和方向,则此空间为矢量场。物理学中常见的有电场、速度场等,地理学中常见的有风场等。受此启发,固碳服务的空间分布也可以理解成一个场：由于环境变量空间分布的差异,区域的固碳服务势能也有所差异,在这一势能的驱使下,区域内部以及区域与环境之间会产生固碳服务在空间上的流动。区域的碳源与碳汇受到来自自然生态和社会经济等多方面的影响,二者之间的平衡数量关系也比较复杂。在此引入流动比率（Ri）的概念,通过计算碳源和碳汇的比值,以表征同一区域中碳源与碳汇平衡的数量关系[19]。其概念模型为：
\[R_i=\frac{O_{p_i}}{I_{p_i}}\ \ （6）\]
式中,Ri为像元i的源汇流动比率,Opi为像元i某时段的碳排放量,即通过碳源模型计算得到的碳排放量;Ipi为同一时段像元i的碳吸收量,即生物固碳量。Ri 能够表达区域固碳服务的自给率,Ri 越小,区域的固碳服务自给率越强,反之则越弱。当Ri<1,表明区域碳供给倾向大于碳需求倾向,即为供给型区域;反之,当Ri>1,表明区域碳需求倾向大于碳供给倾向,即需求型区域。运用地理信息系统手段,对固碳服务的空间格局及其流动比率做空间表面分析,可绘制区域固碳服务空间流动的方向和趋势。
1.5 条件概率与因子状态可视化
生态系统固碳过程涉及大量环境因子,其格局优化需以各因子对固碳功能的影响程度为前提。本文基于贝叶斯信念网络（Bayesian Belief Networks,BBN）,以概率统计为理论基础,借助python编程,模拟和计算各固碳因子的不确定性,并将其分布格局空间可视化（图2）。首先,参考自然分级结果,综合考虑研究区实际情况与坡度等相关等级划分标准,将备选因子的栅格数据按照离散化标准（表1）做离散化处理。考虑到生物碳库和土壤碳库的生物循环和代谢过程存在较大区别,因此本文中二者的备选因子也有所差异。其中,生物碳库的备选因子有高程、7月平均气温、潜在蒸散发、降水、坡度、土壤类型、太阳辐射,土壤碳库的备选因子则包括高程、净初级生产力、降水、坡度、土壤类型、太阳辐射、年平均气温。
其次,基于贝叶斯网络原理,计算碳库中各变量状态的条件概率。设有栅格变量,记作A={A1,…,Ai},其中Ai为栅格变量A的不同状态。各变量不同状态的栅格数量记为SUM Ai,研究区栅格总数记为SUM,那么栅格变量A的不同状态概率PAi为:
\[P_{Ai}=\frac{SUM_{Ai}}{SUM}\ \ （7）\]
又有事件B={B1,…,Bj},事件Ai和Bj同时发生记作Ai Bj,满足此条件的栅格总数记为SUM Ai Bj。那么事件Ai和Bj同时发生的联合概率P（Ai Bj）为：
\[P(A_iB_j)=\frac{SUM_{A_iB_j}}{SUM}\ \ （8）\]
那么事件Bj发生的条件下,事件Ai发生的条件概率P（Ai |Bj）为：
\[P(A_i|B_j)=\frac{P(A_iB_j)}{P(B_j)}\ \ （9）\]
通过Python编程统计各因子对应的栅格数,计算备选变量两两配对的条件概率,生成条件概率表,并利用Python将其可视化为像素图。根据像素图,选取关键变量状态子集,并将关键变量状态组合空间可视化。
然后,筛选关键变量最优状态子集。其中,最优状态是指当固碳状态等级达到最优（状态等级=3）时所对应的最高概率备选变量状态,即最适宜于固碳的状态,可通过分析条件概率像素图得出;关键变量是指备选因子中对固碳结果影响最大的几个变量,本文根据实际情况将关键变量数量设置为2个,可通过计算熵减度筛选出。
Table 1
表1
表1碳库备选因子及其数据状态离散化标准
Table 1Carbon pool alternative factors and their discrete states
| 因子 Factors | 名称 Name | 状态等级 States | 实际值范围 Actual values |
|---|---|---|---|
| 高程 DEM | DEM | 1 | 208m<DEM≤700m |
| 2 | 700m<DEM≤1300m | ||
| 3 | 1300m<DEM≤3754m | ||
| 7月平均气温 Average temperature in July | T_July | 1 | 16℃<July_T≤22℃ |
| 2 | 22℃<July_T≤25℃ | ||
| 3 | 25℃<July_T≤28℃ | ||
| 潜在蒸散发 Potential evapotranspiration | PET | 1 | 839mm·a-1< PET≤1400mm·a-1 |
| 2 | 1400mm·a-1< PET≤1650mm·a-1 | ||
| 3 | 1650mm·a-1< PET≤2114mm·a-1 | ||
| 降水 Precipitation | Prec | 1 | 440mm·a-1≤Prec≤550mm·a-1 |
| 2 | 550mm·a-1< Prec ≤650mm·a-1 | ||
| 3 | 650mm·a-1< Prec ≤740mm·a-1 | ||
| 坡度 Slope | Slope | 1 | 0°≤Slope <15° |
| 2 | 15°≤Slope <30° | ||
| 3 | 30°≤Slope ≤77° | ||
| 土壤类型 Soil types | Soil | 1 |  土,新积土,潮土,紫色土,红壤,风沙土,黄褐土,黄绵土Tier soil, Alluvial soil, Fluvo-aquic soil, Purple soil, Red soil, Aeolian sandy soil, Yellow-cinnamon soil, Loessal soil |
| 2 | 山地草原土,水稻土,石灰岩土,粗骨土,褐土,黄绵土,黄棕壤 Mountain steppe soil, Paddy soil, Limestone soil, Skeletal soil, Cinnamon soil, Loessal soil, Yellow-brown soil | ||
| 3 | 亚高山草甸土,山地棕壤,暗棕壤,棕壤,沼泽土 Subalpine meadow soil, Mountain brown soil, Dark brow forest soil, Brown forest soil, Bog soil | ||
| 太阳辐射 Solar radiation | Sol | 1 | 387 MJ·m-2<sol≤410 MJ·m-2 |
| 2 | 410 MJ·m-2<sol≤430 MJ·m-2 | ||
| 3 | 430 MJ·m-2<sol≤450 MJ·m-2 | ||
| 净初级生产力 Net primary productivity | NPP | 1 | 0< NPP ≤650 g·m-2·a-1 |
| 2 | 650 g·m-2·a-1< NPP ≤1000 g·m-2·a-1 | ||
| 3 | 1000 g·m-2·a-1< NPP ≤1700 g·m-2·a-1 | ||
| 温度 Temperature | T | 1 | 5℃≤ T ≤10℃ |
| 2 | 10℃< T ≤13℃ | ||
| 3 | 13℃< T ≤16℃ | ||
| 生物碳 Biological carbon | C_plant | 1 | 0≤ C_above≤600 g·m-2 |
| 2 | 600 g·m-2< C_above≤1200 g·m-2 | ||
| 3 | 1200 g·m-2< C_above≤2700 g·m-2 | ||
| 土壤有机碳密度 Soil organic carbon density | SOCD | 1 | 0≤ SOCD ≤160000 g·m-2 |
| 2 | 160000 g·m-2< SOCD ≤260000 g·m-2 | ||
| 3 | 260000 g·m-2< SOCD ≤700000 g·m-2 |
新窗口打开
显示原图|下载原图ZIP|生成PPT图2生态系统固碳服务功能关键变量筛选及其可视化
-->Fig. 2Filtering and visualization of key variables in carbon sequestration
-->
通过熵差法计算每个变量对固碳结果的熵减度：
\[I=E(Q)-E(QF)=\sum_q\sum_f\frac{P(q,f)log_2[P(q,f)]}{P(q)P(f)}\ \ (10)\]
其中,Q为生物（或土壤）固碳节点,F为生物（土壤）碳库的备选因子;E(Q)为生物（或土壤）固碳的熵,E(QF)为生物（或土壤）固碳与备选因子联合发生的熵;P(q)为生物（或土壤）固碳状态q发生的概率,P(f)为备选因子状态f发生的概率,P(q,f)为其联合概率。I为熵减度,取值范围在0—1,I越大则该因子对固碳结果的影响越大,反之I越小则该因子对固碳结果的影响越小。
根据上述方法,可以挑选出关键变量的最优状态子集,在关键变量的最优状态能够同时发生的区域,获得高固碳效益的概率最高。
2 结果
2.1 固碳服务供给与需求
土地利用方式是人类社会活动作用于自然环境和气候变化的重要介质,不同的土地利用方式具有不同的生物碳固定能力和土壤碳固定能力。在关天经济区,草地、林地和耕地的生物固碳总量分别占全区的41.90%、28.54%和28.18%,其土壤固碳总量分别占全区的27.59%、36.41%和34.73%,是研究区的主要碳汇区域。而城市用地的单位面积固碳能力为269.47 g·m-2·a-1,单位面积平均碳排放能力却高达128.59 g·m-2·a-1,总碳排放占全区总量的49.35%,是研究区最主要的碳源区域。耕地的单位面积碳排放能力虽不强（4.81 g·m-2·a-1）,但因面积广大而导致其总碳排放能力占全区的33.94%,成为研究区的主要碳源之一。由此可见,土地利用格局对于区域碳源汇平衡具有重要意义。研究区固碳服务在整体上供大于求,但各区县的供需状况差异明显（图3）。分布在秦岭山脉上的凤县、太白县等,依靠秦岭良好的生态基础,拥有较高的固碳服务供给量（1.10×1010—1.33×1010 gC）,而其经济发展水平有限,所需要的固碳服务量相对较低（0.41×105—0.89×105 gC）。天水市麦积区、周至县、宝鸡市陈仓区和洛南县虽然固碳服务需求量较高（8.23×105 —10.80×105gC）,但由于供给量（0.90×1010—1.42×1010 gC）远远大于需求,因此固碳平衡仍为供大于求。而西安多个市区、咸阳市各区、宝鸡市金台区等的固碳服务供给量（≤5.57×108 gC）明显低于其他市区,而需求量却较高（≥6.84×105 gC）。
显示原图|下载原图ZIP|生成PPT图3各区县生态系统固碳服务供需统计图
1：白水县;2：宝鸡市陈仓区;3：宝鸡市金台区;4：宝鸡市渭滨区;5：彬县;6：长武县;7：澄城县;8：淳化县;9：大荔县;10：丹凤县;11：凤县;12：凤翔县;13：扶风县;14：富平县;15：甘谷县;16：高陵县;17：韩城市;18：合阳县;19：华县;20：华阴市;21：户县;22：泾阳县;23：蓝田县;24：麟游县;25：礼泉县;26：陇县;27：洛南县;28：眉县;29：蒲城县;30：乾县;31：千阳县;32：秦安县;33：清水县;34：岐山县;35：三原县;36：商洛市商州区;37：太白县;38：天水市麦积区;39：天水市秦州区;40：铜川市王益区;41：铜川市耀州区;42：铜川市印台区;43：潼关县;44：渭南市临渭区;45：武功县;46：武山县;47：西安市灞桥区;48：西安市碑林区;49：西安市长安区;50：西安市莲湖区;51：西安市临潼区;52：西安市未央区;53：西安市新城区;54：西安市阎良区;55：西安市雁塔区;56：咸阳市秦都区;57：咸阳市渭城区;58：兴平市;59：旬邑县;60：杨陵区;61：宜君县;62：永寿县;63：张家川自治州;64：柞水县;65：周至县
-->Fig. 3Regional ecosystem carbon sequestration service demand and supply statistics
1: Baishui; 2: Chencang; 3: Jintai; 4: Weibing; 5: Binxian; 6: Changwu; 7: Chengcheng; 8: Chunhua; 9: Dali; 10: Danfeng; 11: Fengxian; 12: Fengxiang; 13: Fufeng; 14: Fuping; 15: Gangu; 16: Gaoling; 17: Hancheng; 18: Heyang; 19: Huaxian; 20: Huaying; 21: Huxian; 22: Jingyang; 23: Lantian; 24: Lingyou; 25: Liquan; 26: Longxian; 27: Luonan; 28: Meixian; 29: Pucheng; 30: Qian; 31: Qianyang; 32: Qinan; 33: Qingshui; 34: Qishan; 35: Sanyuan; 36: Shangzhou; 37: Taibai; 38: Maiji; 39: Qinzhou; 40: Wangyi; 41: Yaozhou; 42: Yintai; 43: Tongguan; 44: Linwei; 45: Wugong; 46: Wushan; 47: Baqiao; 48: Beilin; 49: Changan; 50: Lianhu; 51: Lingtong; 52: Weiyang; 53: Xincheng; 54: Yanliang; 55: Yanta; 56: Qindu; 57: Weicheng; 58: Xingping; 59: Xunyi; 60: Yangling; 61: Yijun; 62: Yongshou; 63: Zhangjiachuan; 64: Zhashui; 65: Zhouzhi
-->
2.2 固碳服务供需关系的空间格局
通过固碳服务流动比率（Ri ）空间分布,可以更直观地了解区域的碳平衡状况（图4）。关中平原各市县的城区普遍表现为Ri ≥1,分布特征为Ri 以西安市、咸阳市和渭南市为高值中心向四周递减;而研究区大部分区域Ri <1,属于供给型区域,秦岭山脉和北山山系的Ri 则明显偏低,部分森林核心区域Ri 几乎为零。根据研究区固碳服务供需的格局分区（图5）可知,研究区可分为三个主要的碳源集聚中心,分别是以天水市（Ri>0.06）和彬县（Ri >0.06）为中心的两个碳源集聚中心以及以西安市（Ri >0.20）为中心的多级碳源集聚中心。
显示原图|下载原图ZIP|生成PPT图4研究区生态系统固碳服务流动比率空间分布
-->Fig. 4Current ratio spatial distribution of carbon sequestration
-->
显示原图|下载原图ZIP|生成PPT图5研究区生态系统固碳服务平均流动比率空间分区
-->Fig. 5Spatial partition of average carbon sequestration current ratio in study area
-->
显示原图|下载原图ZIP|生成PPT图6研究区生态系统固碳服务流动比率总和空间分区
-->Fig. 6Spatial partition of total carbon sequestration current ratio in study area
-->
其中,西安市（碳源中心A,Ri >0.20）为一级碳源聚集区,杨陵区（碳源中心B,Ri >0.14）、咸阳市和高陵县（碳源中心C,Ri >0.10）为二级碳源聚集区,渭南市、武功县等（碳源中心D,Ri >0.04）为三级碳源聚集区。研究区各县级以上行政区的碳源汇流动比率总和分布情况主要表现为以关中平原为中心轴,南部流动比率的变化率略大于北部,尤其是眉县、周至县、武功县和杨陵区之间,流动比率呈现出骤降趋势（图6）。供给特征最强的区域（∑Ri≤200）主要为秦岭山脉和北山山系的各县,供给特征较明显的区县（400<∑Ri≤800）主要分布在研究区西侧的宝鸡市和天水市,而关中平原则分布着供给特征较弱的区县（∑Ri >800）。
固碳服务空间分布的不均衡致使需求型区域与供给型区域的流动比率相差上百倍,根据物质循环原理,需求型区域所消耗的碳将由其周边的供给型区域提供,这在一定程度上会引起固碳服务由供给型区到需求区的空间流动。研究区的固碳服务空间流动整体上表现为：在秦岭中段和东段,大量的固碳服务向以西安市为中心的关中城市群涌进,关中平原北侧也有固碳服务流入;在秦岭西段,固碳服务则主要流向天水市,且流量较小;北山的彬县呈现出一定的弱供给特征,麟游县、旬邑县等有小支固碳服务向彬县流入补给;东北侧蒲城县表现出较明显的弱供给特征,其固碳服务的补给源则是铜川市、澄城县、华县等（图7）。
显示原图|下载原图ZIP|生成PPT图7研究区生态系统固碳服务的空间流动示意图
-->Fig. 7Diagram of carbon sequestration spatial flow
-->
2.3 固碳的条件概率与格局分析
2.3.1 条件状态子集筛选 将生物碳库和土壤碳库的固碳能力与其备选因子的条件概率可视化表达为像素图（图8、图9）,图中将每个因子离散化为3个等级,横纵坐标是每个因子的不同状态,纵坐标是事件A,横坐标是事件B,两两状态所对应的色块表示B条件下A发生的条件概率。颜色由蓝到红表示概率从0到1。图中最右一列九宫格表示当某一地理位置的固碳功能处于不同状态时,对应位置的各种变量状态的组合情况。最下方一行九宫格则可以表示当某一地理位置的变量在不同状态时,对应位置的固碳功能的条件概率;其中当固碳能力状态为3时,其所对应的高概率备选变量状态集合为最优固碳状态的高概率条件状态子集（简称“最优状态子集”）。令生物固碳功能处于高中低三种状态时,其对应的备选变量状态组合称为“备选变量状态子集Ⅰ、Ⅱ、Ⅲ”。那么,状态子集Ⅰ= {DEM=1, T_July=2, PET=2, Prec=2, Slope=1, Soil=2, Sol=2}。状态子集Ⅱ= {DEM=2, T_July=2, PET=2, Prec=2, Slope=1, Soil=2, Sol=2}。状态子集Ⅲ= {DEM=3, T_July=2, PET=2, Slope=2, Soil=2, Sol=2}。生物碳库最优状态子集= {DEM=3, T_July=3, PET=1, Prec=1, Slope=3, Soil=3, Sol=2}。其中,P(DEM=3)=0.4513, P(T_July=3)=0.3474, P(PET=1)=0.5169, P(Prec=1)=0.3475, P(Slope=3)= 0.5254, P(Soil=3)= 0.4939, P(Sol=2)=0.3362。
显示原图|下载原图ZIP|生成PPT图8生物碳库变量间的条件概率像素图
-->Fig. 8Conditional probability pixel map of each two variables in biological carbon pool variable set
-->
显示原图|下载原图ZIP|生成PPT图9土壤碳库变量间的条件概率像素图
-->Fig. 9Conditional probability pixel map of each two variables in soil carbon pool variable set
-->
同样地,令土壤固碳功能处于高中低三种状态时,其对应的备选变量状态组合称为“关键变量状态子集Ⅳ、Ⅴ、Ⅵ”。那么子集Ⅳ= {DEM=1, NPP=1, Prec=2, Slope=1, Soil=2, Sol=2, T=2},子集Ⅴ= {NPP=2, Prec=2, Slope=1, Soil=2, Sol=2, T=2},子集Ⅵ= {DEM=3, NPP=3, Slope=2, Soil=2, Sol=2, T=2}。土壤碳库最优状态子集= {DEM=3, NPP=3, Prec=1, Slope=3, Soil=3, Sol=2, T=3}。其中,P(DEM=3)= 0.4089,P(NPP=3)= 0.9013,P(Prec=1)=0.2736, P(Slope=3)= 0.4397, P(Soil=3)= 0.4118, P(Sol=2)=0.0994,P(T=3)= 0.3387。
2.3.2 状态子集的空间格局 将各变量状态子集分别做空间可视化表达可知,生物碳库子集Ⅰ主要分布在富平县东部韩城县中部以及三原县北部和宝鸡市渭滨区北部;子集Ⅱ主要分布在研究区北部及秦岭北向分支;子集Ⅲ主要分布在研究区北部本山山系。生物碳库最优状态子集则主要分布在研究区西南侧秦岭山脉沿线及秦岭北向支线的部分区域（图10）。
土壤碳库的各组变量子集同样也主要分布在秦岭山脉沿线及本山山系（图11）。子集Ⅳ分布区域由研究区的东北角经关中平原延伸至宝鸡市辖区。子集Ⅴ该子集的分布较前一子集而言更为零散,主要分布在研究区北部,秦岭的北向分支有零星分布。子集Ⅵ主要分布在旬邑县的北部及其与耀州区、淳化县交界处和凤县的北部。土壤碳库最优状态子集则主要分布在研究区西南侧秦岭山脉沿线及秦岭北向支线的部分区域。
显示原图|下载原图ZIP|生成PPT图10生物固碳功能关键变量状态子集空间分布图
-->Fig. 10Spatial distribution map of key variable state subset of biological carbon sequestration
-->
显示原图|下载原图ZIP|生成PPT图11土壤碳固定功能关键变量状态子集空间分布图
-->Fig. 11Spatial distribution map of key variable state subset of soil carbon storage
-->
2.4 固碳服务空间格局优化
2.4.1 关键变量的最优状态子集 经过熵减算法得出生物、土壤固碳与其备选因子的熵减度（表2）,熵减度越高表示该因子对固碳能力的影响越大。其中高程（DEM）对生物固碳影响最大,熵减度为0.0160,其次为潜在蒸散发（PET）,熵减度为0.0082。坡度和降水对生物固碳的影响程度相当,熵减度同为0.0032,太阳辐射和7月平均温度的熵减度紧随其后,土壤类型对生物固碳的影响程度最不显著。而净初级生产力（NPP）对土壤固碳影响最大,熵减度为0.0795,其次为高程（DEM）,熵减度为0.0106。气温（T）和太阳辐射（Sol）对土壤固碳的影响程度相当,熵减度同为0.0023,降水（Prec）和坡度（Slope）的熵减度紧随其后,土壤类型（Soil）对土壤固碳的影响程度最不显著。
Table 2
表2
表2备选因子对固碳的敏感性（按熵减度从大到小排序）
Table 2Sensitivity of alternative variable to carbon sequestration (variables are listed in order of influence on sequestration from most to least influential)
| 生物碳库 Biological carbon pool | 土壤碳库 Soil carbon pool | ||
|---|---|---|---|
| 变量 Node | 熵减度 Entropy reduction | 变量 Node | 熵减度 Entropy reduction |
| DEM | 0.0160 | NPP | 0.0795 |
| PET | 0.0082 | DEM | 0.0106 |
| Slope | 0.0032 | T | 0.0023 |
| Prec | 0.0032 | Sol | 0.0023 |
| Sol | 0.0025 | Prec | 0.0017 |
| T_July | 0.0021 | Slope | 0.0011 |
| Soil | 0.0013 | Soil | 0.0004 |
新窗口打开
取熵减度前二的因子作为碳库的关键因子,即生物碳库的关键因子为DEM和潜在蒸散发,土壤碳库的关键因子为NPP和DEM。
结合最优状态分析可知：生物碳库的关键变量最优状态子集={DEM=3, PET=1},其中P(DEM=3)= 0.4513, P(PET =1)= 0.3475;土壤碳库的关键变量最优状态子集={NPP=3, DEM=3},其中P(NPP=3)= 0.9013, P(DEM=3)= 0.4089。如果能够同时满足某一碳库的两个关键变量最优状态条件,那么该区域是该碳库关键变量的最优状态分布区;如果只满足其中一个,那么该区域是该碳库关键变量的次优状态分布区。
生物碳库关键变量的最优状态子集（图12）主要分布在宝鸡市南部秦岭山脉和天水市西南角,该区域中生物固碳状态为最优的概率可以达到54.36%。而生物碳库关键变量的次优状态子集则主要在天水市、秦岭山脉和北山山系,该区域中生物固碳状态为最优的概率为42.49%。其中分布于宝鸡的关键变量最优子集大多已是省级或国家级森林保护单位,受保护程度较高,但天水市西南侧的区域则有待加强。如果要进一步发展生物固碳服务的储备区域,则可以从关键变量次优状态子集的广大区域着手。
土壤碳库关键变量的最优状态子集（图13）主要分布在秦岭山脉沿线、天水市西南角和咸阳市东北角,这些区域中土壤固碳状态为最优的概率高达92.84%。而土壤碳库关键变量的次优状态子集则主要分布在天水市境内,该区域中生物固碳状态为最优的概率为19.32%。土壤固碳能力的强度很大程度上依赖于植被对土壤的有机质输入,因此保护土壤固碳能力的关键在于保护植被。
显示原图|下载原图ZIP|生成PPT图12研究区生物碳库关键变量最优与次优状态子集分布
-->Fig. 12Suitable area and core area of biological carbon pool in study area
-->
显示原图|下载原图ZIP|生成PPT图13研究区土壤碳库关键变量最优与次优状态子集分布
-->Fig. 13Suitable area and core area of soil carbon pool in study area
-->
2.4.2 固碳格局优化建议 低碳目标导向的固碳格局应以生物碳库和土壤碳库的关键变量最优状态分布区域作为维持该区域基本碳汇功能的碳汇核心保护区。在此基础上,以生物碳库和土壤碳库的关键变量次优状态分布区域作为继续提升研究区碳汇功能的优选区域,充分利用固碳关键变量最优状态分布区域,能够在维持现有固碳服务供需格局的基础上,进一步增强研究区的碳汇功能,使研究区的低碳经济发展更加有保障。
对比研究区2014年固碳估算结果与生物、土壤碳库关键变量最优状态子集的空间分布情况,得出研究区生物固碳和土壤固碳的适宜优化区域（图14）。这些区域满足固碳关键变量最优状态组合条件,但固碳能力未达到最优状态,可作为研究区低碳目标导向下的固碳格局优化的首选区域。其中,生物固碳的适宜优化区域明显大于土壤固碳适宜优化区。生物碳库的优化应将天水市武山县和秦州区划为重点碳汇格局优化区,加强绿化工程的建设。此外,麦积区、清水县、张家川回族自治县以及凤县和太白县需根据实际情况开展不同程度的生态修复工程。而土壤固碳格局的优化则需要从秦岭中段着手。由于土壤固碳能力在很大程度上依赖于土壤有机质的输入量,因此需要秦岭中段各市县在保护已有森林的同时,改善适宜优化区的植被覆盖,增加其土壤的有机质输入,从而提高适宜优化区的土壤固碳量。
显示原图|下载原图ZIP|生成PPT图14研究区生物固碳与土壤固碳的适宜优化区域
-->Fig. 14The optimized area for carbon sequestration
-->
3 讨论
3.1 区域生态系统固碳服务的空间流动
目前关于生态系统服务评估的研究大多数依赖于静态指标,例如BRYAN等基于3-PG-2模型预测全球气候变化和土地利用变化背景下的澳大利亚碳汇供给[32]和未来国土部门对于温室气体减排需求[33]。对于固碳服务在供给区与受益区之间的空间流动研究主要还在探索阶段。SERNA-CHAVEZ等[34]提出一个全球尺度的多种生态系统服务流分析框架,在气候调节服务部分中分析了碳封存与碳固定的供给热点（provisioning hotspots,p）、受益区（benefiting areas,b）以及流动区（flow area,f）及其相互关系,认为与b、p不重合的区域都是流动区,固碳服务的空间流动可能具有空间依赖性。但关于固碳服务空间流动的实质性研究成果还较少[14]。本文受到“场”概念和“资源流动”相关研究的启发,以碳源汇流动比率的等值线为基础,借助“等压梯度力”的原理给出研究区碳源汇空间流动示意图,在碳源汇的空间流动上做了可取的尝试,为生态安全、生态补偿等研究领域提供可参考的思路。但是根据SPANs的概念模型,生态系统服务可以通过分布式空间模拟详细地描绘出生态服务的流动路径、流量大小和流动过程等[11-13],这为本文的下一步研究提供了重要的思路。精确地可视化服务流还需要更系统地考虑固碳服务产生和被消费的过程机制,对服务流的载体、路径、流量、衰减和转化等做具体的描述和研究。以此为基础准确界定固碳服务的供给区和受益区,结合流动路径、流量等将生态系统服务流转化为经济价值,将为生态补偿的应用提供有力的科学依据。
3.2 基于贝叶斯与熵减度的空间格局优化
关于土地利用结构的低碳优化模式与策略的相关研究近年来有所发展,例如赵荣钦等[22]基于线性规划评估了南京市土地利用碳减排潜力,并给出低碳土地利用优化建议。但此类研究多注重各土地利用类型面积的分配上,关于土地利用空间布局上的优化讨论比较欠缺。本文基于贝叶斯原理,借助GIS手段计算关键变量对碳固定影响的不确定性,并将不确定性和关键因子分布做可视化分析,研究结果对于区域的空间格局优化具有一定的参考价值。当然,该模型也存在一定的缺陷和不足之处：备选因子的选取、判断最优关键变量子集的阈值以及模型本身的结构优化问题依然有待探讨和完善。在未来的相关研究中,应综合考虑更多的因子（例如人文、社会经济发展等因素）以使研究结果更加贴合实际情况。关于空间格局优化,则可以结合未来气象预测、经济发展等对未来土地利用进行更加细致的预测。另外,受方法所限,文章关于不确定性分析还不够深刻,未来研究中应结合多情景分析探讨变化中的多个因子之间的不确定性关系。
3.3 不确定性分析
本文耦合了多个模型以模拟研究区生态系统固碳能力,从空间流动和不确定性的角度较客观地分析了固碳空间格局,但多个模型的耦合也给结果带来了一定的不确定性。目前人口数据的空间化模型和方法趋于多样化,本文选用的多源数据融合的方法,能够从自然、人文等多方面反映人口的空间分布,对于后续的碳源空间布局分析具有较好的参考价值。但该模型也存在一定的不足,例如参数设置复杂,融合权重的确定较为主观等。目前应用较广的还有基于夜间灯光数据和NDVI及土地利用类型数据的人口数据空间化模型,但该方法更适合于中大尺度的研究,且存在灯光像元过饱和及溢出等问题[23]。随着移动基站等新型数据源的发展,基于“大数据”的人口数据空间化模型也日益引起相关领域的关注[35-36],这可能是未来人口数据空间化的一大趋势。
在固碳供给估算部分,研究通过野外采样和实验,模拟和验证了土壤表层（30 cm）土壤有机碳含量,提高了土壤固碳估算数据的可靠性。然而受客观条件所限,研究缺乏对NPP估算的验证,这给估算过程带来了一定程度的不确定性。在未来的研究中,将更加注重研究的验证部分,增加净初级生产力的野外考察和验证实验。
此外,受先验知识和研究水平所限,本研究只考虑了关天经济区内部各区县的固碳服务流动,而没有讨论研究区边缘区县与研究区外部的固碳服务交流。在未来研究中,我们将把区域与外界的固碳服务交流纳入研究范畴,力争更加全面和系统地分析区域的固碳服务空间流动。同时,生态系统服务流的研究不应局限于单种生态系统服务的研究和探讨,多种生态系统服务流之间的时空耦合与权衡/协同关系是该领域研究未来发展的重要方向。
4 结论
本文结合多源数据融合模型、CASA模型、遥感碳循环过程模型,量化关天经济区固碳服务的供需平衡,引入流动比率,借助贝叶斯原理和熵模型等多种模型和方法,模拟固碳服务空间流动,探讨研究区固碳服务的格局优化对策。结果表明：研究区固碳服务需求小于供给,具有较明显的空间流动,生物碳库的适宜优化主要分布在天水市武山县和秦州区,土壤碳库的适宜优化区主要分布在秦岭中段各区县。本研究以空间可视化的方式梳理生态系统固碳服务从供给区到受益区的流动过程,揭示区域生态系统固碳服务功能转变的空间规律;使用贝叶斯原理和熵减模型,从概率的角度探讨研究区固碳格局优化的不确定性问题,为生态系统服务流和区域固碳格局优化提供了一种新的思路。The authors have declared that no competing interests exist.
