京津冀及周边城市PM2.5污染空间关联网络及季节演化研究

安海岗，李佳培，张翠芝，董志良

河北地质大学，河北 石家庄 050031

自 2003年起，雾霾天气就经常出现，雾霾污染不仅对经济造成重大损失，而且对人们的身体健康造成重大伤害。雾霾污染并不仅仅是某一地区的环境问题，而且具有复杂的空间关联与溢出效应。它会随着空气流动、交通运输、产业转移等经济活动扩散到其他地区（马丽梅等，2014），不同地区间的雾霾污染存在着复杂的相互依赖和很强的空间交互影响，某一地区的雾霾污染会扩散到周边甚至更远区域（Lüthi et al.，2015）。对于雾霾污染的空间关联性，多数学者运用莫兰指数（刘晓红等，2017）、GARCH模型、二次指派程序及引力模型等对其进行检验（逯苗苗等，2017）。由于城市雾霾污染数据较易获取，国内外学者多对城市雾霾污染进行研究，发现不同城市间的雾霾污染存在空间溢出效应与外源影响（王立平等，2016），地理近邻效应对环境污染溢出关系具有显著的正影响（刘华军等，2018），雾霾污染跨区域输送还存在一些规律（王树强等，2019）。近些年来，京津冀地区的雾霾污染已经严重影响到了人们的生活健康，引起了政府与社会的高度关注，国内外学者也对该区域的环境污染问题进行了大量的研究，取得了丰硕的成果。研究发现京津冀地区雾霾污染的主要因素是工业与建筑粉尘、煤炭消耗与汽车尾气（李云燕等，2016）；在京津冀地区，雾霾污染经常呈现一种明显的空间自相关（王一辰等，2017），不同城市之间存在较强的空间溢出效应。

自从WS小世界模型（Watts et al.，1998）和BA无标度模型（Barabási et al.，1999）的开创性工作发表以后，在多个领域掀起了复杂系统与复杂网络的研究热潮，把复杂系统与复杂网络理论引入社会治理的研究范畴，是近年来管理学领域研究的显著趋势。复杂网络的理论与成果目前广泛应用于多个领域，但是在环境管理工程研究中应用较少。在雾霾污染研究方面，复杂网络可用来对中国城市PM2.5按全年与季节进行区域划分，并揭示不同区域扩散物理过程（薛安等，2015）；可与空间计量学相结合识别雾霾污染的空间关联网络形态并揭示其整体特征和微观模式，进而构建空间网络权重对雾霾污染与经济发展的关系进行研究（刘华军等，2018）；可与QAP方法结合重构引力模型，研究不同子群内部之间的密切关系，不同因素对网络结构的不同影响作用（逯苗苗等，2017）。

由以上研究可知，目前运用复杂网络理论对城市雾霾污染的研究多是从全局视角出发，对不同城市雾霾污染的空间分布与相互影响进行分析，而很少有学者对区域内城市雾霾污染的空间关联及空间关联随季节演化情况进行研究。京津冀及周边城市都处于北方，属于雾霾污染重灾区，这些城市在冬季均有采暖期，在采暖期内，雾霾污染尤其严重。对雾霾污染季节演化进行研究时，国内外学者也多对其在传统四季的演化进行分析，鲜有学者对采暖期与非采暖期内雾霾污染关联对比与演化情况进行分析。京津冀及周边城市之间的雾霾污染是个相互关联、相互影响的复杂系统，发现该区域雾霾污染相互关联与随季节演化的规律，可为区域污染联合治理、产业与能源结构调整提供政策依据（蒋超等，2018）。

PM2.5作为雾霾的主要成分，具有颗粒小、易吸附有毒物质等显著特征，对人体健康危害极大（齐梦溪等，2019）。本文以京津冀及周边共 31个城市为研究对象，从中国空气质量在线监测分析平台搜集了2015年3月16日—2020年3月15日期间，共5年的城市PM2.5日均浓度数据；计算出了不同城市 PM2.5日均浓度数据之间的皮尔逊相关系数，以城市为节点，相关系数与城市最短距离的比值为权重，相关系数的均值作为是否存在连边的阈值，构建了京津冀及周边城市 PM2.5污染空间关联加权复杂网络，并对此网络的整体结构、节点度、聚集系数、加权聚集系数、中心线等进行了分析。对 PM2.5污染季节演化进行研究时，将其划分为采暖季与非采暖季，以2015年3月16日—2020年3月 15日样本数据为例，对其关联性做对比分析；构建 2015—2016年非采暖季城市关联网络图与采暖季城市关联网络图，并对不同城市空间关联网络随采暖季节演化情况进行了分析。

1 数据与方法

1.1 数据

京津冀及周边城市环境污染的主要污染物是PM2.5，因此，本文采用PM2.5来度量雾霾污染的严重程度。以京津冀及周边共31个城市为研究对象，选取2015年3月16日—2020年3月15日之间的PM2.5日均质量浓度数据作为样本，单位为μg·m-3，数据来自于中国空气质量在线监测分析平台（https://www.aqistudy.cn/）。31个城市分别为北京、天津、石家庄、保定、邯郸、邢台、唐山、沧州、廊坊、衡水、秦皇岛、张家口、承德、太原、大同、阳泉、晋中、长治、晋城、朔州、忻州、济南、德州、聊城、郑州、开封、安阳、焦作、濮阳、鹤壁、新乡。由于PM2.5污染的持续性，对于缺失的数据，可用前一天与后一天的平均值代替。

PM2.5污染在不同区域之间的溢出与传输与多种因素有关，包括空间距离、气象条件及地形地貌等，如果只考虑相关系数不能很好地体现城市PM2.5污染之间的空间关联。因此，本文在构建城市PM2.5污染空间关联网络过程中，引入了城市距离这个变量，不同城市之间的直接距离可以通过城市间不同的经纬度计算得出。在本研究中，不同城市的经纬度来源于百度搜索引擎。

1.2 研究方法

1.2.1 PM2.5时空分布特征

为了展现京津冀及周边城市 PM2.5污染在时序上的变化，首先求出不同城市 PM2.5的月度均值，然后把不同城市 PM2.5以折线图形式进行展示；对于京津冀及周边城市 PM2.5污染在空间的分布，以2016年12月—2017年11月期间不同城市PM2.5的月度均值为样本，利用 ArcGis软件对不同城市PM2.5月度均值的空间分布进行展示。

1.2.2 皮尔逊相关系数与城市 PM2.5污染空间关联网络构建

皮尔逊相关系数是由英国统计学家卡尔·皮尔逊提出的，常用来度量两个时序变量之间的相关性。皮尔逊相关系数值介于-1与1之间，这个相关系数也称作“皮尔逊积矩相关系数”，常用小写字母r表示，计算公式为：

式中，X，Y分别为两个变量；Xi，Yi是变量X，Y对应点i的观测值；是X样本的平均数；是Y样本的平均数。

复杂网络是由大量节点与边构成的有向或无向网路图。复杂网络中的图可定义为：G=(V,E)，其中，V是网络中所有节点的集合，E是所有边的集合。网络的邻接矩阵M表示各个节点间的连接关系，对于无向图，矩阵元mij可以表示为：

构建城市 PM2.5污染空间关联网络时，本文以皮尔逊相关系数矩阵除对角线以外的所有数值的均值作为阈值，这种阈值选择方式能够使得网络的拓扑结构较为稳定（吴翎燕等，2013）。如果两个城市之间的相关系数大于该阈值，则存在连边，反之，则没有。由此，可构建出京津冀及周边城市PM2.5污染相关邻接矩阵M。为了体现不同城市之间的空间关系，本文将相关系数与城市间最短距离的比值作为网络结构中边的权重（薛安等，2015）。在实际构建网络过程中，为了消除城市距离与边权的量纲，对城市距离与权重进行了归一化处理。

1.2.3 关联网络拓扑性质与季节演化

1.2.3.1 度的概念

度（Degree）表示网络中一个节点与其他所有节点直接连线数的总和。复杂网络分为有向网络和无向网络，在有向网络中，度分为出度与入度。度是在网络分析中刻画节点重要性的最直接度量指标。一个节点的节点度越大就意味着这个节点的度中心性越高，该节点在网络中也就越重要。度可描述为：

其中，ki为节点i的度数。

1.2.3.2 聚类系数与修正后聚类系数

聚类系数是用来描述网络中节点连接紧密程度的一个指标，即小集团结构的完美程度。在无权网络中，其计算公式可描述为：

由聚类系数的计算公式可知，该指标仅仅考虑了网络结构的完美程度。当网络中某个节点度数较小时，可能聚类系数较大。此时，如果只考虑其聚类系数并不能很好地度量该节点在网络的地位与重要性（马永军等，2019），因此，本文将度的影响纳入到聚类系数中，提出了修正后聚类系数，其计算公式为（马宇博等，2018）：

其中，Cri表示修正后聚类系数；wi为节点i度的权重，为网络中所有节点的平均度。

1.2.3.3 中介中心性与接近中心性

在复杂网络中，中介中心性又称中间中心性，多用来描述该节点控制其他节点之间的交往能力，一般通过经过该节点的最短路径数目来刻画其重要性，其计算公式为（周磊等，2008）：

其中，c(i,j)是节点i与j之间所有最短路径的总数；ck为通过中间节点k的路径数。

接近中心性是用来度量网络中某一节点与其他节点接近程度的一个指标，可用该节点到所有其他节点的最短路径距离累加起来的倒数来表示。换言之，某一节点与其他节点越接近，它的接近中心性越高（Alexander et al.，2019）。其计算公式为（魏兰清等，2019）：

1.2.3.4 季节关联性对比与网络演化分析

为了体现城市 PM2.5指数关联性在采暖季与非采暖季的差异，将2015年3月16日—2020年3月15日共5年的数据，分为采暖季数据与非采暖季数据，计算出部分重要城市与其他城市在不同季节的皮尔逊相关系数，然后以柱状图形式作对比分析；构建2015年非采暖季与采暖季城市PM2.5污染空间关联网络图，并对其平均度、平均聚类系数、图密度及连边数等指标作对比分析。

2 结果与分析

2.1 城市PM2.5时空分布特征

2.1.1 时间动态变化特征

图1为2015年3月16日—2020年3月15日期间京津冀及周边31个城市的PM2.5浓度月平均值变化趋势图。由此图可知，这些城市的PM2.5月均浓度在11月至次年2月期间最高，这个时间段内，北方城市处于冬季取暖期，PM2.5污染最严重。进入3月以后，随着取暖期的结束，PM2.5浓度逐渐下降，空气质量逐渐改善。在3—10月这段时间内，空气质量较好，维持在较为稳定的状态。进入11月以后，PM2.5污染又开始攀升，进入了另一个循环。

由图1可知，在每年的12月，京津冀及周边城市PM2.5值最高。其中，2015年12月，保定、衡水、德州等城市的 PM2.5污染最严重，月度均值达到了 200 μg·m-3以上，属于重度污染。在 2016年 12月，石家庄的 PM2.5月均值最高，达到 250 μg·m-3以上，属于严重污染；邯郸、安阳的污染也比较严重，月度均值达到了200 μg·m-3以上。由于污染的持续严重，政府加大了对雾霾的治理，在2017年冬季，PM2.5污染情况有所改善，与2016年相比，PM2.5月度均值有所降低，但是在2018—2019年冬季，PM2.5污染又出现了反弹。

图1 京津冀及周边城市PM2.5质量浓度月平均值Fig. 1 Monthly average value of PM2.5 mass concentration in Beijing-Tianjin-Hebei and surrounding cities

2.1.2 空间分布特征

为了体现京津冀及周边城市 PM2.5浓度在空间上的分布，选取2016年12月—2017年11月共12个月的数据，计算出每个月份的 PM2.5平均浓度，并运用 ArcGis软件进行可视化展示，结果如图 2所示。由图可知，2016年12月—2017年2月之间，京津冀及周边城市 PM2.5污染较为严重，尤其是在2016年12月。在空间分布上，中间城市带污染最为严重，具有以石家庄、保定、衡水、邢台、邯郸等城市为核心向外围逐渐减轻的分布特征；北部城市受 PM2.5污染侵害较小，以承德、张家口、大同为例，这3个城市一年四季月度平均PM2.5都在良以上。从时间维度上可知，从2016年12月开始，PM2.5污染逐渐减弱。在2017年3月，只有石家庄、保定空气质量为轻度污染，其他城市空气质量均为良，2017年4月，只有邯郸空气质量为轻度污染。2017年5—9月京津冀及周边城市空气质量月度均值均达到良或优的等级。从2017年10月开始，PM2.5又卷土重来，石家庄、邯郸、邢台最为严重，空气质量变为轻度污染。

2.2 污染关联网络图分析

通过计算京津冀及周边31个城市PM2.5日均浓度相关系数，并求平均可得此网络相关系数的均值为0.6032，因此本文以0.6032为阈值构建城市之间的连边。如果两个城市之间的相关系数大于或等于0.6032，则存在连边，否则不建立连边，由此，可构建出一个包含31个节点，220条边的无向加权网络。

2.2.1 整体网络特征与度分布

图3为京津冀及周边31个城市的PM2.5污染空间网络图，此图密度为 0.473，平均度为 14.194，网络图比较紧密。由此图可知，北京与廊坊、安阳与邯郸、安阳与鹤壁、焦作与新乡、晋中与太原等这些城市对之间 PM2.5污染关联比较紧密。对于这些关联比较紧密的城市，它们的 PM2.5污染治理政策应尽量保持一致性。此网络中度数最大的节点为石家庄与邢台，这两个城市虽然与其他城市都存在PM2.5污染空间关联，但是关联并不紧密。

表1为京津冀及周边31城市的网络度数。由表可知，石家庄、邢台、邯郸、鹤壁度值最高，均达到20以上，这些城市与其他城市PM2.5污染空间关联较为紧密，在 PM2.5污染治理中要重点关注。大同、朔州、张家口度值最低，尤其是张家口，度数为 2，它仅仅与大同与承德存在空间关联，且张家口的空气质量一年四季较为优良，因此，对于这类城市，重点要制定政策保持城市环境的良性循环，减弱其他污染严重城市对其环境的空间影响。

图2 2016年12月—2017年11月京津冀及周边城市PM2.5平均浓度空间分布Fig. 2 Spatial distribution of PM2.5 average concentration in Beijing-Tianjin-Hebei and surrounding cities from December 2016 to November 2017

2.2.2 聚类系数与修正后聚类系数

通过统计计算，可求得此网络平均聚类系数为0.745，网络聚集性较强。在此网络中，聚类系数最大的城市是秦皇岛，值为 1，最小为张家口，值为0。图 4为秦皇岛的自我网络图，由此图可知，此小集团结构的完美程度高。但是，秦皇岛这个节点在所有节点中度值排序位列第五，此时，如果只考虑它的聚类系数没有什么现实意义。因此，本文在考虑度数的基础上，提出修正后聚类系数这个概念，即将度除以平均度作为权重，然后再乘以聚类系数。修正后的聚类系数能够更好地反映城市节点在网络中的重要性。

图3 京津冀及周边城市PM2.5污染空间关联网络图Fig. 3 Spatial network of PM2.5 pollution in Beijing-Tianjin-Hebei and surrounding cities

由表2可知，邯郸、聊城、衡水、石家庄、邢台、安阳、鹤壁为修正后聚类系数最高的7个城市，这些城市度值也比较高，其中邯郸的修正后聚类系数最高，石家庄的度值最高。结合图2中2016年12月城市PM2.5平均浓度空间分布图可知，这些城市大都位于中间城市带，污染均比较严重，说明在京津冀及周边区域内部这些城市在 PM2.5污染上有较强的空间溢出效应，在 PM2.5治理时应对它们进行重点防控。图5为石家庄的自我网络图，由图5可知，石家庄与多数城市均存在 PM2.5污染空间关联，仅仅没有与晋州、朔州、大同、张家口相连。

图4 秦皇岛自我网络Fig. 4 Ego network of Qinghuangdao

2.2.3 中心性分析

对京津冀及周边31个城市的接近中心性与中介中心性做统计计算，结果如图 6所示。由图可知，接近中心性值最高的5个城市为石家庄、邢台、邯郸、鹤壁及保定，中介中心性值最高的5个城市为石家庄、邢台、忻州、承德及晋中，石家庄与邢台这两个城市的接近中心性与中介中心性均最高。这说明石家庄与邢台位于网络的中心，与其他城市PM2.5污染关联最接近，其他城市节点通过它们发生关联的概率更高。利用Gephi软件对京津冀及周边城市PM2.5污染空间关联网络图做模块化处理，并运用ForceAtlas 2算法进行布局展示，结果如图7所示。由图7也可知，31个城市可分为3个子群，石家庄与邢台位于不同子群链接的重要位置。石家庄与邢台在此PM2.5污染空间关联网络中的重要性与它们在该区域所处位置存在很大关联，由图 2可知，石家庄与邢台处于京津冀及周边城市的核心区域。

表1 京津冀及31个周边城市度数Table 1 Degree of 31 cities in Beijing-Tianjin-Hebei and surrounding areas

表2 京津冀及周边城市修正后聚类系数Table 2 Modified clustering coefficient (Mcc) of Beijing-Tianjin-Hebei and surrounding cities

图5 石家庄自我网络图Fig. 5 Ego network of Shijiazhuang

2.3 季节关联性对比与季节关联网络演化

京津冀及周边31个城市均处于中国北方城市，每年存在采暖季与非采暖季，为了探究不同城市在不同季节 PM2.5指数相关性的差异，本文选取石家庄（采暖季PM2.5日均浓度为中度污染水平以上）、北京与天津（轻度污染，经济发达）、张家口（空气质量为良，度值较小）4个城市作为样本，对这4个城市在采暖季与非采暖季期间与其他城市PM2.5指数的相关性作了统计计算与对比，结果如图8所示。

图7 京津冀及周边城市模块化与力引导布局展示Fig. 7 Modularization and strength guidance layout display of Beijing-Tianjin-Hebei and surrounding cities

由图8可知，石家庄与其他城市的采暖季相关系数大都比非采暖季相关系数高，部分城市甚至差异较大，例如邯郸、廊坊、济南、聊城、郑州、安阳、焦作及鹤壁，其他城市差异较小；北京、天津与其他城市的采暖季相关系数与非采暖系数大都差异不大，只有个别城市存在较大差异，且采暖季相关系数较高，例如北京与邯郸、北京与安阳、北京与鹤壁、北京与济南、天津与邯郸、天津与济南；对张家口与其他城市的采暖季相关系数与非采暖季相关系数作对比，发现差异较大的城市为大同、阳泉、晋中、长治、朔州、忻州，但是非采暖季相关系数较高。

图6 京津冀及周边城市接近中心性与中介中心性Fig. 6 Closeness centrality and between centrality of Beijing-Tianjin-Hebei and surrounding cities

图8 部分城市采暖季与非采暖季相关系数对比示意图Fig. 8 Comparison diagram of Pearson correlation coefficient between heating season and non-heating season in some cities

为了对采暖季与非采暖季不同城市 PM2.5污染的空间关联性差异作进一步分析，本文选取 2015年3月16日—2016年3月15日之间的数据作为样本，构建了采暖季与非采暖季城市 PM2.5污染空间关联网络图，结果如图9所示。

由图9可知，采暖季关联网络图与非采暖季关联网络图相比，网络密度更高，边数更多。对平均度与平均聚类系数进行统计计算，采暖季平均度为18.39、平均聚类系数为 0.82，非采暖季平均度为7.94、平均聚类系数为 0.66。由以上可知，在采暖季期间，城市 PM2.5污染指数 PM2.5相关性更强，PM2.5污染空间关联网络图连接更紧密，城市节点度数更高。这说明，在采暖季期间，京津冀及周边31个城市PM2.5污染较为严重时，更容易受到其他城市 PM2.5污染的影响，它们之间存在较强的空间溢出效应；在非采暖季，这些城市空气质量关联性较小，空间溢出效应相对较弱。

3 结论与建议

3.1 结论

图9 2015—2016年采暖季与非采暖季城市PM2.5污染空间关联网络图Fig. 9 Spatial correlation network of Urban PM2.5 pollution in heating season and non-heating season from 2015 to 2016

（1）京津冀及周边31个城市PM2.5浓度值随季节变化差异较大，在采暖期内最高，进入3月后，随着采暖期的结束，浓度值逐渐降低；在空间分布上，PM2.5浓度值在中心区域（石家庄、保定、邢台等城市）最严重，具备从核心向四周逐渐降低的特征。

（2）京津冀及周边城市 PM2.5污染空间关联网络较为紧密，度值较高的节点，聚集性较强，中心性较高，它们控制其他城市 PM2.5污染传输的能力最强，与其他城市PM2.5污染关联最紧密。

（3）在采暖季与非采暖季，不同类型城市PM2.5指数相关性对比存在较大差异。中度污染水平以上城市，采暖季相关系数大都比非采暖季相关系数高；轻度污染且经济发达城市采暖季相关系数与非采暖系数大都差异不大；空气质量较好的城市，部分城市相关系数对比差异较大，但是非采暖季相关系数较高。

（4）采暖季与非采暖季相比，城市 PM2.5污染空间关联网络更为紧密，密度更高。

3.2 建议

（1）每年的12月、1月、2月是PM2.5污染最为严重的时间段，在该段时期内，京津冀及周边城市均处于采暖期，燃煤消耗较为严重，环保部门对于污染排放大户应加强监管，制定更为严厉的政策控制污染源的排放。

（2）位于京津冀及周边中间城市带的城市，污染均比较严重，在 PM2.5污染上有较强的空间溢出效应，对于这类城市当地环保机构应在政府部门与生态环境部的领导下，建立重点监控监测机制，并制定措施切断或者减弱这些城市与其地区 PM2.5污染的传输通道。

（3）对于 PM2.5污染关联比较紧密的城市（例如北京与廊坊），两地环保部门应与京津冀及周边地区大气污染防治领导小组进行及时的沟通与协商，建立专门负责两地 PM2.5污染治理的联合监督监管机构，并根据两地的经济、人口、资源与环境等实际情况，制定较为一致的 PM2.5污染协同治理政策。

（4）在制定污染防治与治理政策时，应考虑不同类型城市的实际情况与季节因素，因时因地的制定与之相适应的污染治理政策，尤其是采暖季期间，应制定相应的政策减弱不同城市之间 PM2.5污染的空间溢出效应。

（5）京津冀及周边城市PM2.5污染关联紧密，应当推进京津冀及周边城市 PM2.5污染的协同与联合治理，实现 PM2.5污染预防与治理一体化，构建区域内全联通的无死角的 PM2.5污染防治立体网络。