北京市燃气工业锅炉CO2排放因子及不确定性分析

陈操操, 邱大庆, 于凤菊, 孙 粉, 胡永锋, 王立波, 孙大利

北京市应对气候变化研究中心，北京 100161

气候变化不仅是环境问题，还是发展问题，已经成为影响区域经济和社会可持续发展的重要因素[1]. 明确温室气体排放清单是推动碳减排控制目标的一项重要基础性工作[2]，降低碳排放评估的不确定性对于碳减排政策选择和成本效益分析十分重要[3-4]. 我国对温室气体排放评估工作大部分始于20世纪90年代末，省市级和区域碳排放核算大多采用IPCC国家温室气体清单指南，主要是基于源清单排放因子方法[5-7]，也是基于生产侧碳排放核算方法. 但是问题在于，IPCC国家温室气体清单指南默认的排放因子库主要获取于发达国家，与我国的能源和计量体系存在一定差异，不一定能反映国内实际情况，直接使用IPCC国家温室气体清单指南中提供的缺省数据，可能会给国内能源消费温室气体排放评估带来较大的不确定性.

在国内城市和区域尺度上，采用具有本地化特征的排放因子较为少见[8]，缺乏在测试参数层面进行的系统性分析[9-11]，或者不同研究者对同一问题的研究结论存在争议，并且煤炭等固体燃料研究相对较多，而其他燃料品种较为罕见，有必要在条件成熟的区域或城市建立标准的地方CO2排放核算方法和数据库[12]. 例如，LIU等[13]在Nature期刊撰文，基于中国602个煤样数据，估计了中国化石燃料CO2排放因子，认为IPCC提供的煤炭缺省因子要比中国煤炭排放因子高40%，中国能源燃烧和水泥生产过程CO2排放估计具有较大的不确定性. 腾飞等[14]研究认为文献[13]误将中国不同煤种平均排放因子与IPCC焦煤缺省排放因子进行比较，但二者煤种排放数据集并不可比. 朱德臣[15]评估了国内燃煤发电企业元素碳含量，检测得到单位热值含碳量与省级排放清单推荐缺省值的差异及其对温室气体排放的影响. 其他典型研究还包括，国家发展和改革委员会应对气候变化司组织开展中国2005年温室气体排放清单研究，对不同品种的燃煤进行测试，认为明确电煤和原煤的热值对于降低排放清单的不确定性具有重大意义[16]. 蔡博峰[17]综述了当前国际煤炭燃烧CO2排放因子的主要进展和热点问题，提出目前研究焦点主要是对IPCC国家温室气体清单指南的排放因子进行区域级的修正和对其分类体系改进. 高欣等[18]基于中国典型动力煤质数据，建立了煤质和煤炭排放因子的经验关系式. 总体上，国内基于测试数据开展的温室气体排放因子评估及参数细化研究还较少[10-11,19]，尤其是通过对设备运行情况测量得到的排放因子数据更少，现有部分研究结论也存在争议.

当前北京市能源结构中天然气已经取代燃煤，成为能源消耗的重要构成部分，天然气锅炉也成为城市主要能耗设备. 为压减燃煤、减少污染物排放，北京市大规模开展天然气锅炉替代燃煤锅炉，2013—2017年北京市99.8%的燃煤锅炉已被淘汰. 到2017年底，北京全市共淘汰燃煤锅炉 4 453 台、13 259 蒸吨，完成燃气锅炉低氮改造约 7 000 台、23 000 蒸吨[20]. 2017年北京市天然气消费折标煤 2 268×104t[21]，占全市能源消费总量的31.8%. 因此对北京市天然气锅炉CO2排放因子的细化评估能够有效提高城市温室气体清单的准确性，明确北京实测的CO2排放因子与IPCC国家温室气体清单指南相应缺省值是否存在差异，对天然气锅炉碳排放的评估结果有无显著区别，了解不同天然气锅炉出力是否会对排放因子带来影响等问题，具有重要的科学与现实意义.

鉴于此，该研究基于燃烧设备实测样本和文献调研方法，通过分析影响排放因子的多个关键参数特征，对北京市天然气锅炉CO2排放因子进行评估，给出排放因子推荐值，运用概率分析方法识别和量化排放因子分析的不确定性，得到统计分布区间和不确定性，进而提高本地化排放因子的评估精度.

1 研究方法

1.1 排放因子计算方法

IPCC国家温室气体清单指南和中国省级温室气体清单指南都给出了不同燃料燃烧CO2排放因子的计算方法，可根据设施燃烧燃料品种平均低位发热值、单位热值含碳量和碳氧化率来量化燃烧燃料产生的CO2排放因子. 公式如下：

(1)

式中：EF为单位化石燃料燃烧CO2的排放量，t/t或t/m3；NCV为化石燃料的低位热值；MJ/kg或MJ/m3；CC为化石燃料的潜在单位热值含碳量，t/TJ；OF为燃料的碳氧化率，%；44/12为CO2与C的分子量之比. 式(1)可以用于计算天然气锅炉排放因子，其中天然气的潜在单位热量含碳量计算方法见式(2).

(2)

式中：CCi为天然气潜在单位热值含碳量，t/TJ；bi为天然气各气体组分中分子式的碳原子个数；ai为天然气中各气体组分的摩尔百分比，%；NCV为天然气理想发热值，MJ/m3. 由于企业测试报告仅给出了天然气各组分的体积百分比，因此该研究参考GB/T 11062—2014《天然气发热量、密度、相对密度和沃泊指数的计算方法》[22]，查表获取天然气中甲烷、乙烷、丙烷等各组分的压缩系数，将体积百分比转换为摩尔百分比，通过摩尔发热量和压缩因子进一步计算理想混合气体的理想低位发热值.

传统工业锅炉的能源计量是对工业锅炉效率的测试和研究，工业锅炉能源效率等于燃烧效率和传热效率的乘积. 燃气锅炉碳氧化率与燃烧效率直接相关，与天然气燃料有很大的联系. 《2005中国温室气体清单研究》[16]中将天然气缺省氧化率取值为99%. 有研究者曾通过正平衡、反平衡等3种方法评估工业锅炉的碳氧化率，得到燃气锅炉平均碳氧化率在90.91%～92.16%之间，但该数值明显偏低[23]. IPCC国家温室气体清单指南没有推荐不同燃烧设备的碳氧化率因子，而是默认为100%. 李青青等[9]对比了IEA、CDIAC、EDGAR和EIA国际四大典型碳数据库，并对中国碳排放核算进行参数级别的对比，发现气体燃料碳氧化率取值范围为98%～99.5%. 综上，该研究选择燃气锅炉碳氧化率取值为99%，与《2005中国温室气体清单研究》[16]的口径保持一致.

1.2 蒙特卡洛分析

蒙特卡洛分析是一种通过随机变量的统计试验分析及随机模拟获得近似结果的方法. 工业天然气锅炉二氧化碳排放因子评估的不确定性来源于两个方面：①计算方法的不确定性；②活动水平数据和模型各参数的不确定性. 该研究中不确定性主要考虑活动水平数据观测随机误差分布以及模型中参数优化方法导致的不确定性，没有考虑计算方法引起的不确定性.

该研究采用蒙特卡洛(Monte Carlo)方法研究各参数对天然气锅炉CO2排放因子的不确定性影响. 这是一种采用统计抽样理论近似求解数学问题或物理问题的方法，其基本步骤[24]：①选择需进行不确定性分析的输入参数；②确定输入参数的概率分布与范围(probability distribution functions, PDF)；③产生模型输出值的概率分布，包括输入参数随机变量抽样和参数不确定性在模型中增值传播的模拟. 该研究参照以上步骤，其中，输入为北京市天然气工业锅炉排放因子涉及参数，可利用监测数据根据经验数据统计分析或专家判断法获取，数据分布代表了参数的变化范围，用概率分布模型进行描述. 常用的概率分布模型包括正态分布(Normal)、对数正态分布(Binomial)、伽马分布(Gamma)、贝塔分布(Beta)、韦伯分布(Weibull)、均匀分布(uniform)和三角分布(triangle)等. Monte Carlo方法最终利用所得参数概率密度分布结果作为输入变量，从模型输入中随机抽取值并计算相应模型输出大小，经过n次重复抽样，由n次模型输出构成的概率分布代表了模型输出的不确定性.

同时，该研究采用敏感性分析(sensitivity analysis)方法识别排放因子输入参数变化对蒙特卡洛模型输出变化的影响排序，通过皮尔逊(Pearson)相关系数法评估模型输入和输出的线性相关程度，将所有相关系数归一化到100%，从而解释模型输入参数对天然气工业锅炉排放因子不确定性的贡献.

1.3 数据来源与分析

为了保证数据的代表性和权威性，通过北京市企业(单位)碳排放报送和填报系统，共获取北京市227台天然气工业锅炉能效的检测报告. 这些检测报告全部具备中国计量认证(CMA)标志. 所测试的工业锅炉主要分为天然气蒸汽锅炉和热水锅炉，并不包括电厂锅炉，设备出力为0.35～165蒸吨范围不等，能够代表北京市的天然气工业锅炉器具. 涉及锅炉检测参数的数据主要包括：甲烷等12种天然气气体的体积成分、低位发热量、设计出力、实际出力、运行效率等. 为了识别和剔除统计异常值，采用Grubbs准则检验在0.05显著性水平下有无离群值，使测试数据样本更符合客观实际. 主要测试指标如表1所示.

表1 天然气锅炉主要测试指标

2 排放因子评估结果

2.1 北京市天然气锅炉排放因子计算结果

该研究共获得227个北京市天然气锅炉的有效锅炉实测检测数据，约占全市全部燃气锅炉数量的3.2%，通过式(1)～(3)得到单位热值含碳量、低位热值等表征关键参数，进而计算得到天然气锅炉实测排放因子. 按实物量口径，北京市天然气锅炉CO2排放因子在18.416～20.702 t/(104m3)之间，平均值和中位数分别为20.300和20.517 t/(104m3)，在0.05显著性水平下Shapiro-Wilk检验并不显著服从正态分布(见图1). IPCC1996国家温室气体清单指南和IPCC2006温室气体清单优良做法指南(分别简称“IPCC1996”和“IPCC2006”)中能源燃烧天然气碳排放因子单位为t/TJ，相比于IPCC1996，IPCC2006进一步给出了碳排放因子的高低值区间. 这些排放因子并不特别针对中国的燃料，且在实际计算中需要利用CO2和C的分子量比值进一步转化. 为了方便比较，该研究将基于体积的CO2排放因子转换到同口径基于低位热值的CO2排放因子(见表2). 需要指出，IPCC1996和IPCC2006中天然气缺省的氧化因子均为1，而笔者参考《2005中国温室气体清单研究》[16]，选择天然气锅炉碳氧化率取值为99%.

当数据为正态分布时，算术平均值为描述数据集中趋势的最优选择，但算术平均值易受极端数据的影响；对于不满足正态分布、比较分散的数据，常用中位数表征其平均分布特征[25]. 由于北京市天然气锅炉CO2排放因子并不服从于正态分布且有偏，笔者推荐取中位数值20.517 t/(104m3)或15.378 t/TJ. 由表2可见，北京市地方特征的天然气锅炉CO2排放因子稍高于IPCC1996和IPCC2006缺省值，且前者取值范围要比后者更集中.

图1 北京市天然气锅炉CO2排放因子分布特征Fig.1 Distribution chart of carbon dioxide emission factor of natural gas boilers in Beijing

表2 北京市天然气锅炉实测排放因子

2.2本地化的关键参数特征

图2 北京市天然气单位热值含碳量频数分布Fig.2 Frequency distribution of carbon content per unit low calorific value of natural gas in Beijing

根据式(1)，该研究中影响本地化特征排放因子的关键参数主要有3个——天然气单位热值含碳量(CC)、天然气单位低热值(NCF)和碳氧化率(OF). 从前两个参数频率分布和基本统计结果来看，在0.05显著性水平下北京市天然气单位热值含碳量不符合正态分布，其平均值和中位数分别为15.367和15.378 kg/GJ，根据数据分布形态建议取值为中位数(15.378 kg/GJ)，北京市天然气单位热值含碳量稍高于IPCC1996和IPCC2006默认值(15.3 kg/GJ). 由图2可见，在0.05显著性水平下北京市天然气收到基低位热值不符合正态分布，平均值和中位数分别为 36 391 和 36 750 kJ/m3，分布建议取值为中位数(36 750 kJ/m3)，与IPCC1996和IPCC2006(38 931 kJ/m3)相比，北京市天然气低位热值稍小. 碳氧化率取值范围在0.99～1之间，检测报告并没有提供碳氧化率检测参数，因此该值参考《2005中国温室气体清单研究》[16]，确定该研究中碳氧化率取值为0.99.

2.3 锅炉吨位和效率对排放因子影响

为了检验天然气锅炉吨位是否对排放因子产生影响，该研究依据227个锅炉样本(其中213个有效样)，按设计吨位(M)分成4组(M≤2 t/h，2 t/h10 t/h)，未对原始数据进行转换，采用经典的单因素方差ANOVA分析方法，基于Student-Newman-Keuls (SNK)和Ducan比较四组平均值差异. 经检验，在0.05显著性水平下北京市4个分组天然气锅炉排放因子未发现各设计吨位分组间存在显著性差异(见表3). 这表明天然气作为气态燃料，燃烧较为充分，锅炉吨位对排放因子的影响不大.

表3 基于不同吨位北京市天然气锅炉CO2排放因子单因素方差分析

从天然气锅炉样本的运行效率看，运行效率平均值在90.0%以上，中位数为90.5%，经过Anderson-Darling非参数检验方法(简称“A-D test”)检验，分布区间符合Weibull分布. 需要指出，该研究的调查样本来自工业锅炉的热效率检测，而锅炉热效率是燃料燃烧效率或氧化率和传热效率的乘积[14]. 由于缺乏计算传热效率的烟气排放损失热量、工业锅炉炉壁散热损失热量等有关检测参数，该研究无法直接建立锅炉运行效率和天然气燃烧氧化率之间的定量关系. 从锅炉出力不同水平看，单因素方差分析显示天然气锅炉运行效率之间无显著差异，间接说明燃料氧化率对天然气锅炉排放因子的影响较小.

2.4 基于蒙特卡洛模型的不确定性分析

界定参数分布及取值范围是蒙特卡洛模型分析的关键步骤. 该研究中参数界定有两种途径：①基于样本数据利用非参数检验方法(A-D test)，尝试可能的分布函数并排序，评估确定参数概率分布并给出最优估值；②依靠IPCC、中国国家温室气体排放清单推荐参数取值范围，综合得到参数概率分布和估值.

输入参数选取如下：北京市天然气单位热值(NCF)样本符合最小极值分布，在90%概率区间(5%～95%)内取值为 35 402～37 234 kJ/m3，推荐取值为中位数 36 750 kJ/m3. 天然气单位热值含碳量(CC)样本符合最小极值分布，在90%概率区间内取值为15.28～15.42 kg/GJ，推荐取值为中位数15.38 kg/GJ. 氧化率(OF)样本符合三角形分布，在90%概率区间内取值为0.99～1，推荐取值为0.99. 模型设定各参数按照概率分布模拟运行 10 000 次试验.

图3为北京市天然气锅炉CO2排放因子的直方图，直方图代表了排放因子的概率分布，即不确定性分布. 由图3可见，北京市天然气锅炉CO2排放因子90%的概率分布范围为19.78～20.84 t/(104m3)，与推荐值〔20.517 t/(104m3)〕相比，排放因子的不确定范围为-3.59%～1.57%. 模型输出不确定范围并不沿平均值对称分布，A-D test检验显示，天然气锅炉CO2排放因子符合Weibull分布. 这是由模型输入北京市天然气单位热值(NCF)、天然气单位热值含碳量(CC)和氧化率(OF)3个参数不确定性范围的分布特点决定的. 为了解各因素变化对CO2排放因子的影响程度，通过方差贡献来验证天然气低位热值、单位热量含碳量和氧化率这3个参数对排放因子的敏感度. 由表4可见，天然气低位热值对排放因子的影响最大，其方差贡献率达到94%，而单位热值含碳量和氧化率二者的方差贡献率仅占6%. 对天然气低位热值进行准确评估，改进活动水平数据收集工作，有助于指导并降低天然气锅炉CO2排放因子不确定性.

注：灰色部分代表CO2排放因子90%的概率区间.图3 北京市天然气锅炉CO2排放因子蒙特卡洛模拟直方图Fig.3 Monte Carlo simulated histogram of carbon dioxide emission factor of Beijing natural gas boiler

2.5 与已有类似研究的对比

为了更好地理解北京市天然气锅炉排放因子所处水平，该研究将北京市天然气锅炉CO2排放因子

表4 蒙特卡洛模型参数敏感度分析

在单位实物量和热值不同单位情况下的数值与典型国家清单、城市清单的数据进行横向对比. 由表5可见：北京市天然气锅炉CO2排放因子与中国、美国、英国国家清单评估数值接近，稍低于IPCC2006推荐缺省值，位于合理的值域区间. 北京市天然气锅炉CO2排放因子数值比纽约和加州同类数值高出9.98%. 查阅发现，纽约市温室气体排放清单并没有给出详细的排放因子数据来源; 加州天然气排放因子与纽约市基本一致，加州温室气体清单排放因子评估方法综合了基于设备设施连续监测(CEMS)及根据化石燃料成分检测确定缺省因子方式. 值得注意的是，美国国家、州和城市清单燃料热值均采用高位热值(HHV)，而非IPCC、中国国家温室气体清单或省级温室气体清单惯常使用的低位热值(LHV)的数值，高位热值和低位热值一般相差10%，这可能是造成北京与美国加州及纽约天然气锅炉CO2排放因子差异的主要原因. 在今后进行中外尤其是中美城市之间对比时，排放因子的差异值得关注.

表5 典型国家或城市天然气固定设施排放因子对比

4 结论

a) 该研究使用北京市天然气锅炉检测样本数据，针对关键参数的分布特征，根据专家经验和国内外相似研究确定参数合理范围，将蒙特卡洛模拟应用于本地天然气锅炉排放因子评估，给出推荐排放因子及具有一定概率分布的范围，提高了排放因子评估精度.

b) 据测算，按实物量口径，北京市天然气锅炉CO2排放因子推荐值为2.052 kg/m3，90%概率分布范围为1.982～2.086 kg/m3；基于热值CO2排放因子推荐值为55.829 kg/GJ，90%概率分布范围平均值为55.788～55.908 kg/GJ，排放因子的不确定范围为-3.59%～1.57%. 经统计学检验，北京市天然气锅炉CO2排放因子符合Weibull分布.

c) 北京市天然气锅炉CO2排放因子稍低于IPCC2006推荐缺省值(0.5%)，与中国、美国、英国国家清单推荐取值相近，但是比纽约和加州同类数值高出9.98%.

d) 从参数敏感度贡献看，天然气低位热值对排放因子的影响最大，其方差贡献率达94%，而单位热值含碳量和氧化率二者的方差贡献率仅占6%. 对天然气低位热值进行准确评估，改进北京市活动水平数据收集工作，有助于指导并降低天然气锅炉CO2排放因子的不确定性.

e) 该研究仅开展了北京市燃气锅炉CO2排放因子的不确定性评估，建议后续同类研究继续加强对国内排放因子的特征参数分析，扩大不同行业、不同温室气体的活动水平数据及排放因子不确定性的研讨，形成分区域的基于实测的排放因子库，增加核算的准确度.