0 引言

科学研究已经明确证明，人类通过燃烧化石燃料排放的温室气体造成了史无前例的全球变暖问题。想要遏制全球变暖，就必须对温室气体排放进行控制，因此全球将全面进入“碳中和”时代。实现碳中和的第一步，就是建立一个完整、精准的碳监测体系。中国环境监测总站于2021年9月发布了《碳监测评估试点工作方案》，指出聚焦区域、城市和重点行业开展碳监测评估试点，目标到2022年，探索建立碳监测评估的技术方法体系，发挥示范效应。建立二氧化碳(CO₂)等温室气体碳汇监测技术体系与核算支持业务体系，为我国碳源汇数据的可测量、可报告、可核查提供技术支撑。各个城市需要遵守特定的标准机制来监测城市的碳排放，并依靠根据燃料使用、能源使用等数据编制的清单来估计总碳排放量和随时间的变化，从而确定各种减排政策的效果。

目前，美国多个城市正在建立碳监测网络，通常是在屋顶或塔顶上对CO₂和甲烷(CH₄)进行原位测量^[1-3]，或使用遥感等方法进行碳监测^[4-7]。美国国家标准与技术研究所(NIST)曾选择洛杉矶和东北走廊城市群组建了高精度CO₂和CH₄塔基原位监测网络^[8-9]，从而开发和改进了城市温室气体的监测技术，并掌握了不同时空尺度上排放估算的不确定性。其中，在东北走廊城市群的碳监测点位选址研究中，将城市现有高塔作为备选点位，基于大气扩散模型和聚类分析的迭代算法选取出最优的建站位置，这些点位对城市区域的温室气体排放具有较强的敏感性。目前国内缺少城市尺度、高时间分辨率的温室气体监测机制，无法很好地掌握城市碳排放现状，进而做出精准的治理、管控措施。

现根据《城市大气温室气体监测点位布设技术指南(第一版)》的要求，利用拉格朗日粒子扩散模型(LPDM)和K-均值(K-Means)聚类算法，对南京所有备选站点进行迭代聚类筛选，根据需求并结合实际情况，获取理论上最优的建站位置，最终筛选出适合的温室气体监测站位置，为科学、准确、完备地掌握南京市温室气体排放特征提供科学支撑。

1 研究方法 1.1 大气传输模型

LPDM模型^[10]不同于通常计算单点的轨迹模式，而是通过计算气块群的运动轨迹，进而实现对大气物质的输送和扩散过程的模拟。该模型可以通过计算大量的气块轨迹去描述空气的输送和扩散，这里的气块并不是真实的气块，而是指无穷小的空气块。作为新一代的拉格朗日模型工具，它可以用来模拟点源、线源、面源和三维大气源排放的大气示踪物质的长期和中尺度输送、扩散、干湿沉降及其辐射衰减过程，其时间上的前向轨迹模拟可以确定排放源的轨迹-扩散过程，而后向轨迹模拟可以确定排放源的影响区域，获取城市区域大气对站点的大气影响函数。

1.2 WRF模式

LPDM模型需要输入气象场数据来进行驱动，使用天气研究和预报模式(WRF)模拟南京市的气象场。该模式是为大气研究和业务预报应用而设计的次世代中尺度数值天气预报系统。它具有2个动态核心、1个数据同化系统和1个支持并行计算和系统可扩展性的软件架构，适用于从几十米到几千公里的各种气象应用。对于研究人员而言，WRF模式可以根据实际大气条件(即通过观测和分析)或理想化条件进行模拟。该模式先进可靠，是获取气象场数据的优良选择，目前在美国国家环境预报中心(NCEP)和其他国家气象中心以及实验室等机构的实时预报配置中被广泛使用。WRF参数设置见表 1。

1.3 K-Means聚类算法

将物理或抽象对象的集合分成由类似的对象组成的多个类的过程被称为聚类。K-Means聚类算法是目前应用最广泛的聚类算法。标准K-Means算法的思路是对于给定的样本集，按照样本之间的距离大小，将输入的样本集D ={x₁，x₂，…，x_m}，划分为k个簇，让簇内的点尽量紧密地连在一起，而让簇间的距离尽量大。假设输出的簇划分为(C₁，C₂，…，C_k)，则目标是最小化平方误差(E)，表达式见式(1)。

$ E=\sum\limits_{i=1}^k \sum_{x \in C_i}\left\|x-\mu_i\right\|_2^2 $

(1)

式中：x——样本；μ_i——簇C_i的均值向量，也称为质心，其表达式见式(2)。

$ \mu_i=\frac{1}{\left|C_i\right|} \sum_{x \in C_i} x $

(2)

k个初始化的质心位置选择对最后的聚类结果和运行时间均有很大的影响，因此需要选择合适的k个质心。如果仅仅是完全随机的选择，有可能导致算法收敛很慢。K-Means++算法就是对K-Means随机初始化质心的方法优化。

K-Means++对于初始化质心的优化策略如下：

(1) 从输入的数据点集合中随机选择1个点作为第1个聚类质心μ₁；

(2) 对于数据集中的每个点(x_i)，计算它与已选择的最近聚类质心的距离D(x_i)，计算公式见式(3)；

$ D\left(x_i\right)=\operatorname{argmin}\left\|x_i-\mu_r\right\|_2^2\left(r=1, 2, \cdots, k_{\text {selected }}\right) $

(3)

式中：k_selected——聚类质心数量；μ——聚类质心。

(3) 选择1个新的数据点作为新的聚类质心，选择的原则是：D(x)较大的点，被选取作为聚类质心的概率较大；

(4) 重复步骤(2)和(3)，直到选择出k个聚类质心，并利用这k个质心作为初始化质心去运行标准的K-Means算法。

本研究使用了改进后的K-Means++算法，能够显著改善分类结果的最终误差，并降低算法实际上需要的计算时间。通过聚类，可以使同类型站点的探测能力相似性最大，而不同类型站点之间的差异性最大。

1.4 数据来源

(1) 采用上海市环境科学研究院编制的2017年长三角排放清单^[11]，空间分辨率为4 km。物种包含细颗粒物(PM_2.5)、可吸入颗粒物(PM₁₀)、一氧化碳(CO)、二氧化硫(SO₂)、氮氧化物(NO_X)、挥发性有机物(VOCs)、氨气(NH₃)和CO₂。

(2) 南京市土地利用类型数据来自30 m全球地表覆盖数据GlobeLand 30(http://www.globallandcover.com/defaults.html?type=data&src=/Scripts/map/defaults/browse.html&head=browse-&type=data)^[12]，它是中国研制的30 m空间分辨率全球地表覆盖数据，总体精度为85.72%，卡帕(Kappa)系数为0.82，数据访问日期为2022年9月14日。GlobeLand 30数据共包括10个一级类型，分别是：耕地、林地、草地、灌木地、湿地、水体、苔原、人造地表、裸地、冰川和永久积雪。

(3) WRF模式的气象驱动数据来自NCEP的全球再分析资料(FNL)(http://rda.ucar.edu/datasets/ds083.2/)，时间分辨率为6 h，空间分辨率为1° × 1°。数据时间为2021年1和7月。用于WRF准确性评估的气象数据来自中国气象数据网(http://data.cma.cn/site/index.html)。

1.5 技术路线

技术路线示意见图 1。

(1) CO₂监测站需要建设在50 m以上的高空。根据城市现有高塔信息，选择50~150 m的高塔位置作为初始点位。根据建设需求和城区分布对高塔进行初筛，初筛后的点位作为备选站点。

(2) 基于大气传输模型(本研究采用LPDM)，获得城市排放通量对所有备选站点的大气影响函数，并获得网格化的大气影响函数结果，对备选站点进行多次K-Means聚类分析。

(3) 使用K-Means聚类分析的迭代算法，对点位进行筛选，剔除掉每次聚类结果中对城市贡献度最小的站点，并将剩余站点进入下次的K-Means迭代计算中。经过多次迭代后，最终留下的站点即为理论上最优的建站位置。

2 结果与讨论 2.1 WRF数据模拟验证

本研究的WRF模拟评估区域拟采用3层嵌套模式，空间分辨率为27，9，3 km。最外层嵌套区域覆盖中国大部分地区，第2层区域覆盖江苏及周边省份，最内层区域覆盖南京市及周边地区。

气象因素在大气污染物的运输、转化和沉积中起着重要的作用，气象数据模拟的准确与否对污染物的模拟至关重要。本研究中，利用南京气象观测站2021年1和7月每日2 m平均气温(T2)、2 m相对湿度(RH2)以及10 m风速(WS10)的观测数据，评估了模式模拟值的准确性。计算其模拟均值，观测均值、偏差、均方根误差和相关系数。气象要素模拟效果评价见表 2。

由表 2可见，1和7月T2，RH2和WS10的均方根误差的平均值为2.13，10.02和1.97，表明模式基本再现了气温、湿度和风速的变化趋势。T2，RH2和WS10的模拟均值与观测均值的相关系数为0.88，0.79，0.44(1月)和0.59，0.68，0.25(7月)，说明观测数据与模拟结果具有较好的一致性。

2.2 CO₂监测站可建设位置评估

南京市土地利用类型面积占比及对应的CO₂排放量占比见表 3。由表 3可见，耕地为南京市面积占比最大的土地利用类型，占比为58.06%，其次为城市，占比为23.37%，此外林地和水体分别占10.22%和7.97%，其他土地利用类型面积占比均＜1%。从CO₂排放量来看，面积占比第2的城市，其CO₂排放量占比最高，为46%，排放主要来自电、热生产活动，制造产业，建筑业和交通运输业，这些行业是碳排放最主要的来源；CO₂排放量占比第2的是耕地，为35%，碳排放可能来自农业生产过程中的土壤管理、秸秆燃烧和化肥使用^[13]。

2020年南京市土地利用类型及高塔分布情况见图 2。由图 2可见，南京市城区一共有210个可建设监测站点的高塔，需要根据建设需求和城区分布对高塔进行初筛：(1)剔除核心城区之外的高塔，(2)仅保留分布在人工地表之上的高塔。经筛选，得到91个备选站点(城市高塔)。

2.3 城市排放通量对备选站点的影响评估

结合WRF模拟结果，使用LPDM模型针对每个备选站点模拟气团的传输过程，分别计算气团来源，进而获得每个网格通量对91个备选站点的大气影响函数。每个备选站点的大气影响函数是1个包含网格经纬度坐标和气团足迹数值的函数序列，可以直观地看出周边区域对站点影响程度的分布。同时因为数据的位置是乱序，而且不计0值，因此需要把数据先填入1个经纬度坐标的二维网格中，再将这个二维网格展开成一维数据来获得1个有序固定的数据。LPDM模型得到的大气影响函数数据属于每个站点的特征数据，K-Means算法正是利用站点的特征数据来对各备选站点进行聚类。

随着季节的变化，气象因素也会产生明显的变化，LPDM模型计算出的大气影响函数和站点探测能力也会发生变化，所以需要分季节进行LPDM模型运算以及后续的聚类计算。因而选择1和7月分别代表冬、夏季开展南京市温室气体监测点位的研究。假如算力、时间充足，还可以对春、秋两季进行同样的计算，也可以把计算的时间尺度从1个月延长至2~3个月。

2.4 K-Means迭代聚类过程

基于2.3中利用LPDM模型计算出的网格化大气影响函数数据，对91个备选站点进行多次K-Means聚类。城市聚类结果和CO₂贡献度评估筛选过程见图 3(a)—(d)，以夏季(2021年7月)为例。

由图 3可见，第1次聚类结束后，91个备选站点被分为8组，标记为8种颜色。值得注意的是，聚类结果明显呈现出地理位置上的聚团效应，这是因为聚类结果中归为同组的相邻站点往往具有相似的气象特征，且站点影响函数结果也较为一致。

对于同组的备选站点，根据城市CO₂贡献度的大小进行剔除。城市CO₂贡献度代表了该站点对城市区域网格的监测覆盖比例，由于温室气体主要在城市区域网格产生，城市贡献度越高，该站点监测效果越好。基于城市土地利用类型和站点的大气影响函数，计算站点的城市CO₂贡献度。第1次聚类结束后，剔除每组聚类结果中城市贡献度最小的站点，然后再次进行聚类，直至只剩最优站点。站点的CO₂城市贡献度计算公式见式(4)。

$ 站点的{\rm{CO}}_{2}城市贡献度=\frac{站点的大气影响函数(城区网格之和)}{站点的大气影响函数(区域所有网格之和)} $

(4)

式中：分子、分母对应的分别是城区网格和区域所有网格。城区网格是基于城市土地利用类型(图 2)中的“城市”这一类型来统计。

2.5 选址位置筛选结果

通过使用冬、夏季数据进行聚类筛选，最终分别得到了南京市冬、夏季的8个最优选址，见图 4(a)(b)(其中编号0—4为2个季节位置相同的点，5—7仅冬季，8—10仅夏季，共计11个)。将最优选址的大气影响函数叠加绘制，可以看出，在最优选址处建设监测站，能对南京市大部分城区的温室气体排放具有相对最大的覆盖范围和较高的敏感性。

南京市CO₂排放量分布情况见图 5，图中对选址站点CO₂排放量进行了标注。由图 5可见，冬、夏季的11个选址站点周边的CO₂排放量由低至高均匀分布，因此在这些点位布设监测站点可以有效观测南京市CO₂排放量的梯度变化，合理反映南京市不同区域内的温室气体排放量。由于冬、夏季气象条件不同，8个最优选址也有所不同，有5个最优选址位置重合(站点编号0—4)，表明这5个重合站点在不同的季节对城区温室气体排放均有较强的敏感性，因而建议作为首选站点选址。而冬、夏季剩余6个站点，可以基于南京市地理位置、气候条件、城市布局及环境管理需求，主要考虑布设城市点、工业区点及城市背景点(边界点)，选取符合要求的监测点位。

3 结语

利用LPDM模型和K-Means聚类算法，通过对南京市现有的210个高塔点位进行迭代聚类筛选，分别得到冬季和夏季各8个温室气体监测建议点位。同时，结合南京市CO₂排放量能够看出，筛选出来的建议点位周边的CO₂排放量由低至高呈均匀分布，说明这些点位可以有效体现出南京市CO₂排放量的梯度变化，合理反映南京市不同区域内的温室气体排放量，因此在最优选址处建设监测站，能对南京市大部分城区的温室气体排放具有相对最大的覆盖范围和较高的敏感性。

在未来的工作中，可以结合城市布局及环境管理需求，从剩余点位中选定适合的监测点位，逐步建立健全温室气体监测网络，从而对南京市不同时空尺度上的温室气体排放达到相对最好的监测效果。