基于物理信息机器学习方法，实现微环境中实时、高空间分辨率的空气污染暴露估算 约翰·G·巴尔齐斯（John G. Bartzis） 伊奥annis·A·萨凯拉里斯（Ioannis A. Sakellaris） 斯皮罗斯·安德罗诺普洛斯（Spyros Andronopoulos） 亚历山德罗斯·韦内塔诺斯（Alexandros Venetsanos） 费尔南多·马丁-洛伦特（Fernando Martín-Llorente） 斯泰因·詹森（Stijn Janssen）

《Environments》：Towards Real-Time, High-Spatial-Resolution Air Pollution Exposure Estimation in Microenvironments Supported by Physics-Informed Machine Learning Approaches John G. Bartzis, Ioannis A. Sakellaris, Spyros Andronopoulos, Alexandros Venetsanos, Fernando Martín-Llorente and Stijn Janssen

　　摘要：在复杂的城市环境中，高空间和时间分辨率下可靠且及时地估算空气污染暴露量仍然具有挑战性，因为污染物浓度受到交通排放、城市形态和气象条件的影响而变化。本研究提出了一种基于物理知识的机器学习框架，用于近乎实时的NO2浓度估算，能够处理精细的空间尺度。该方法结合了有限的稳态计算流体动力学（CFD）模拟与操作性气象和空气质量数据。首先使用特定风向下的CFD模拟来表征特定地点的扩散模式，然后利用每小时的气象观测数据对这些输出进行缩放，以生成基于物理的浓度描述符。通过随机森林和极端梯度提升（Extreme Gradient Boosting）实现的机器学习预测器被训练来通过整合额外的环境和观测特征来改进这些估算结果。该方法应用于比利时安特卫普的一个1公里×1公里的城市区域，并在FAIRMODE跨比较框架内进行了验证。与一个月内收集的105个被动样本的测量结果相比，显示出显著的改进，决定系数高达R2 = 0.965，并且各位置的偏差减少。这些发现表明，将物理建模与机器学习相结合可以在城市环境中实现准确且计算效率高的高分辨率暴露评估。

1. 引言
城市空气污染仍然是一个重大的环境和公共卫生挑战，特别是在人口密集的城市环境中，如二氧化氮（NO2）、颗粒物（PM）和臭氧（O3）等优先污染物的暴露在精细尺度上表现出强烈的空间和时间变异性[1,2,3]。在以交通为主导的城市微环境中，污染物浓度受高度局部化的排放源、复杂的街道尺度流动和扩散过程以及快速变化的气象条件的影响，导致短距离内的浓度梯度明显[4,5]。因此，准确描述这些微环境中的群体暴露对于有效的空气质量管理、流行病学评估和有针对性的缓解策略设计至关重要。

为了应对这些挑战，欧盟通过《环境空气质量指令》（2008/50/EC）建立了全面的监管框架，该指令为关键污染物设定了具有法律约束力的限值和目标阈值[6]。这些目标要求成员国不仅监测空气质量，还要评估合规情况，及时向公众通报，并在标准未达到时实施有效的空气质量计划。实现这些目标需要关于污染物浓度的准确、空间分辨和时间响应的信息，特别是在暴露变异性最高的城市微环境中。

在受监管的城市空气污染物中，NO2特别重要，因为它是一个主要的交通相关排放和道路附近暴露的指标[6,7,8]。由于其与燃烧过程的强关联以及众所周知的健康影响，NO2被广泛用作城市环境中交通引起的空气污染的追踪器[9]。因此，NO2是开发和评估高分辨率暴露估算方法的首选对象，同时也代表了与城市空气质量评估相关的更广泛挑战。

监管的空气质量监测网络为包括NO2、PM和O3在内的主要优先污染物提供连续的在线测量，具有高精度和数据质量。然而，由于经济、基础设施和运营限制，这些网络的空间覆盖范围通常仅限于相对较少的固定监测站[7,8,9,10]。因此，城市区域的大部分——特别是在街道和微环境尺度上——仍然没有实时监测。因此，未测量位置的污染物浓度估算依赖于基于模型的空间外推方法，包括简化的扩散模型、土地利用回归技术和城市尺度的化学传输模型[7,11,12]。在复杂的城市环境中，这些方法在近乎实时的操作环境中应用时可能会伴随较大的不确定性[13]。

高分辨率建模框架，如计算流体动力学（CFD）和先进的城市化学传输模型，可以更好地表示流动、扩散和化学过程，从而能够解析多种污染物的精细尺度浓度模式[14,15,16]。然而，它们的高计算要求通常限制了它们在离线分析或准操作场景中的应用，限制了它们在城市微环境中直接用于实时暴露估算的适用性。

最近在人工智能（AI）和机器学习（ML）方面的进展为解决这些限制开辟了新的途径，使得能够将观测数据与大气过程的物理知识相结合[17,18,19]。机器学习技术越来越多地应用于城市空气污染物浓度（包括NO2、PM和O3）的估算和预测，显示出在捕捉排放、气象和浓度之间的非线性关系方面的强大能力[20,21,22,23,24,25]。然而，纯数据驱动的方法在应用于复杂的扩散机制时可能缺乏物理一致性和鲁棒性。基于物理知识或物理知识的机器学习方法将物理约束或模型推导的信息嵌入到数据驱动的框架中，为数据稀疏的城市环境中提供可靠、高分辨率、近乎实时的空气污染暴露估算提供了有希望的途径[26,27]。

这里的挑战是：从现有的基于物理知识的机器学习方法的经验出发，我们如何进一步发展，以便能够解决复杂的空气质量问题（例如城市环境），甚至实时为用户提供答案。Bartzis等人的最新工作[28]首次采用了这种方法，针对城市微环境中的脉冲排放量化问题进行了研究。

本研究旨在进一步发展一种新的基于物理知识、由AI驱动的方法，用于近乎实时的空气污染估算，能够处理任何复杂程度的微环境，为选定的位置提供暴露估算。所提出的框架适用于主要优先污染物（NO2、PM、O3），本研究主要关注与交通排放和道路附近暴露密切相关的NO2。该方法整合了（a）历史观测数据与流动和扩散的高级建模输出，以及（b）通常仅在少数监测位置可用的同时在线空气质量观测数据。通过结合对物理过程的理解和机器学习技术，该方法旨在为城市内的任何感兴趣的位置提供可靠的暴露估算。重要的是，其目标是通过改善污染热点的识别、增强暴露评估以及为城市空气质量管理提供实时决策支持来支持欧盟的关键空气质量目标。该方法在比利时安特卫普的一个复杂城市区域进行了演示，并在FAIRMODE跨比较框架内进行了验证，重点验证其提供准确、计算效率高且与政策相关的空气质量信息的能力。最终，这项工作有助于弥合高分辨率科学建模与有效实施欧盟空气质量政策所需的操作空气质量管理系统之间的差距。应当注意的是，这个案例研究来自FAIRMODE跨比较练习[29]。该练习由FAIRMODE WG4开发，该工作组负责欧盟空气质量指令（AAQD）背景下的城市微尺度空气质量建模。其中一个主要目标是检验此类模型在AAQD框架内用于空气质量评估和规划的适用程度。

2. 方法论
2.1. 研究地点
所选的地点是比利时安特卫普市，用于FAIRMODE联合跨比较练习，主要目的是评估现有建模方法估算城市热点长期平均空气污染物浓度图的适用性[29]。所选地点是安特卫普市的一个大约1公里×1公里的城市建成区，具有典型的西北欧城市特征，由街道峡谷和开放区域组成。图1展示了该地点的示意图。地点的详细数据可以在Martín等人（2024）[29]中找到。数据总结见表1。

图1. 地点概览。红点表示NO2监测位置，绿点表示被动采样位置。
表1. 地点主要数据。

如图1所示，两个监测站（42R801和42R802）提供每小时平均NO2浓度。此外，通过在该地区安装的105个被动采样网络提供了额外的实验NO2数据，这些样本在2016年4月30日至2016年5月28日期间收集NO2，覆盖了几乎一个月的时间。每个采样位置产生的数据是整个采样期间的NO2浓度平均值。

2.2. 本研究方法
基本思想是利用定期监测数据来构建一种适当的方法，以便在周围区域的任何感兴趣的位置提供近乎实时的暴露估算。本研究方法采用了一种机器学习回归器（MLR）来完成这项任务。所选的MLR使用可用的每小时数据进行训练和测试。挑战在于建立一个可靠、近乎实时的操作方法。该方法将在进行了被动采样的105个位置进行验证（表1）。因此，所提出的方法必须（a）相对简单，最大限度地利用现有专业知识；（b）操作性强，且计算时间和资源消耗最小。具体步骤如下：

2.2.1. 参考模拟
在参考稳态流动和扩散条件下以及规定的排放空间剖面下进行了有限数量的CFD模拟。为此目的需要一个合适的流动和扩散模型。在本研究中，使用了ADREA-HF CFD RANS代码，该代码在过去已被广泛用于此类研究[33]。对于空气湍流闭合，采用了标准的k-ε湍流方案。ADREA-HF的一个特点是它将真实的城市地形引入笛卡尔域，不仅作为体积和表面孔隙率，还作为相关的建筑物/地面/排放表面，通过它们的实际面积和方向。值得一提的是，每条道路都被视为一个近地面排放表面。

研究的计算域尺寸为1020米×1020米×300米，分别沿西-东、南-北方向和垂直向上方向。该域被划分为204×204×60的矩形单元，水平分辨率为5米×5米，垂直网格从地面到顶部不均匀，范围从2米到10米。

根据本研究方法使用的CFD模拟输出包括特定位置的稳态NO2浓度，包括（a）两个监测位置42R801和42R802，以及（b）上述提到的105个被动采样位置。

在以下条件下进行了有限数量的CFD流动和扩散模拟：
（a）所有模拟都在域顶部定义的单一真实参考风速（VREF = 10米/秒）下进行。（即300米）。
（b）每个模拟对应一个风向。选择了覆盖所有360度的等距风向，以确保在需要时能够进行可靠的空间插值。在FAIRMODE跨比较练习中，22.5到45度的增量（即8-16个扇区）似乎可以给出可接受的结果。出于保守考虑，本研究选择了32个方向（NSIM = 32），相当于风向增量11.25度。因此，选定的风向如下：
（c）通过使用ADREA-HF代码模拟相应的1D中性流，并设置上述VREF的速度来计算进气边界条件。
（d）在每个道路表面上，引入了相关的年度交通排放率，以排放通量（μg/m3/秒）的形式表示。

进行了32次CFD模拟，并计算了上述107个位置（即105个被动采样点和两个污染物站点）的“参考”NO2浓度（ns = 1, …, NSIM），以便在需要时在线使用。

2.2.2. 在线流动和扩散模型应用
下一步是避免直接进行CFD模拟，估算2016年的NO2每小时浓度。由于（a）高计算能力和时间要求，以及（b）难以以所需精度获得额外进气边界条件，这样的CFD模拟可能不切实际。一个有吸引力的替代方案是依赖于第一步的参考CFD模拟，并利用迄今为止获得的专业知识，通过适当的缩放/外推近似来获得结果。在这一方向上，近年来人们做出了重大努力，特别是在高时间分辨率下，关注污染物的浓度和扩散行为，包括在复杂地形中的短时间释放[28]。所提出的方法已成功应用于当前的安特卫普数据，具体步骤如下：根据Bartzis等人的研究[28]，考虑了每小时观测到的风速和风向对[(Vn,θn)，n = 1, …, NWIND]。共有8784对风速数据，对应于2016年闰年的总小时数。对于每个传感器和每个小时，使用观测到的风向θn来估计该特定方向的“参考”浓度，这些“参考”浓度是通过线性插值从上述计算结果中得出的。然后，应用1/V缩放规则[28]来估计当前的基于CFD的模型浓度CSMn，具体公式如下：(1) 其中 是由上述安特卫普数据给出的排放时间因子。需要注意的是，如果将年排放量乘以该因子，得到的值就是针对特定n小时的排放量。需要注意的是，CSMn仅由当地排放引起。为了获得最终的模型浓度，我们必须加上安特卫普数据提供的背景浓度(CB)（见表1）。(2) 其中 是反映域外NO2传输贡献的给定每小时背景浓度。值得注意的是，这些数据也被提供用于FAIRMODE模型对比研究[29]。然而，在这个阶段需要做出以下评论：(a) 关于公式(1)的有效性，高斯模型[34]已经成功地应用了这种关系。在复杂环境中，它通常能给出良好的结果。此外，上述缩放方法也已成功应用于当前的安特卫普数据[35]。另一方面，有迹象表明其有效性在很大程度上受到CFD RANS方法本身的有效性限制[36]。在局部流动不稳定且固有的亚小时浓度波动相对较大的特定位置，预计会出现恶化[36]。(b) 结果在很大程度上取决于背景浓度，而这些浓度并不总是容易和直接估计的。(c) 所有的当前计算都是在中性大气稳定性条件下进行的。在城市地区，中性稳定性条件很常见，因为城市树冠的复杂性通常会产生显著的局部风速梯度，增加局部湍流的产生并增强中性稳定性。然而，不稳定的条件也是可能的，不能预先忽略。初步数据分析表明，在这个特定的案例研究中，大部分时间的环境大气稳定性条件确实是中性的[35]。(d) NO2是一种化学活性较强的物质，这可能会影响NO2的预算，尽管对于1公里或更短的距离，通常不会产生显著影响。另一方面，FAIRMODE模型对比研究的结果表明，假设中性稳定性和非化学反应性的CFD模型与NO2观测值相比，能够得到良好的结果。(e) 为了完整性，不应忘记来自CFD建模假设的额外误差，包括流入边界条件。预计通过AI辅助方法来改进结果，可以直接或间接地消除上述一些缺点。2.2.3 机器学习回归器（MLR）目前的重点是合适的监督机器学习方法，能够处理表格数据，如上所述。一个重要的特性是我们正在寻找能够成功满足时间约束的方法。在空气污染应用方面的过往经验是选择方法时的另一个标准。在成熟的回归方法中，随机森林回归（RFR）和极端梯度提升回归（XGBR）被选为起始方法，因为它们在气象学和空气质量问题中得到了广泛的应用。RFR由于其鲁棒性和强大的泛化能力，经常被用于污染物浓度预测和环境回归任务[37,38]。同样，XGBR在空气质量估计和相关气象建模中也展示了高预测准确性[39,40,41]。所选模型的一个重要特性是它们能够在计算时间约束下高效运行。基于集成树的方法，如随机森林（RF）和极端梯度提升（XGBoost），以其可并行化和可扩展性而闻名。RF独立构建树，允许高效并行计算[41]。XGBoost通过优化梯度提升、并行树构建、缓存感知的内存访问和稀疏感知学习进一步提高了计算效率[39]。这些特性使得这两种方法特别适合需要快速模型训练和预测的应用，例如环境监测和操作性空气质量预测系统。MLR特征开发特定MLR算法的关键因素是选择适当的特征变量。需要注意的是，所考虑的站点每小时都进行常规监测，包括理论上影响NO2扩散的参数，如环境风速和方向、温度、湿度和太阳辐射。观察到NO2浓度时间序列显示出大约24小时的日变化行为。因此，一天中的小时也是一个需要考虑的额外特征。另一方面，在本研究中，如表1所示，NO2的每小时浓度在两个地点测量：42R801（背景）和42R802（交通）。逻辑上，应将局部背景（42R801）浓度作为输入特征。因此，局部交通（22R802）浓度将作为目标变量的一部分，如下所述。影响流动和扩散的一个关键参数是地面形态。对于如此高的空间分辨率，这不能仅用几何术语来描述。根据Bartzis等人[28]引入的概念框架，并调整到当前问题，模型浓度CSM（由公式(1)给出）非常重要。为了避免通常有问题的外部传输浓度的影响，提出了以下输入变量(XMODEL)：(3) 其中COBS表示两个监测站测量的NO2每小时浓度。从概念上讲，XMODEL数组以一种相当简单的方式表达了两个重要的决定因素：(1) 通过CSM将详细的地面形态投影到流动和扩散阶段；(2) 通过单个站点（42R801）观测到的浓度影响污染物向局部域的传输，从而抵消了通常有问题的“背景”浓度(CB)的影响。关于用于训练/测试所选MLR的目标变量(TARGET)，选择相当直接：使用第二个站点42R802的观测到的每小时浓度，即(4) 总之，(a) MLR输入变量包括(i) 观测到的气象站风速、风向、温度、相对湿度和太阳辐射的每小时数组，(ii) 日间小时（1–24小时）以及(iii) 由公式(3)给出的XMODEL参数。(b) MLR输出变量TARGET是指第二个站点42R802的观测到的NO2每小时浓度数组，如公式(4)所示。MLR训练和测试RFR和XGBR都是使用Python 3.9和scikit-learn创建的。训练和测试的特征矩阵始终相同，包括“MLR特征”部分描述的变量。关于行数，选择了所有没有缺失数据的小时。因此，考虑了8784个2016年的数据点中的8342个每小时数据点。其中80%用于训练，20%用于测试。两种方法获得的结果都给出了良好且可比的确定系数；即R2参数分别为：对于测试数据，R2 = 0.936（RFR）和R2 = 0.965（XGBR）；对于训练数据，R2 = 0.960（RFR）和R2 = 0.989（XGBR）。各种输入变量的重要性（即杂质减少量和排列重要性[29]）在图2中进行了说明。所得值对于两种MLR来说似乎是可比的。主要输入变量似乎是XMODEL，对于RFR和XGBR都是如此。这是合理的，因为这些参数结合了形态和关键的物理因素。风参数的重要性较低可以解释为风效应在很大程度上已经包含在XMODEL中。与热相关的参数（温度、湿度、太阳辐射）的作用似乎有限，这支持了在中性大气稳定性近似在这种特定应用中相当有效的假设。RFR表明XMODEL的重要性相对较高，而其他参数的重要性相对较低。图2. 当前MLR输入特征的重要性。当前MLR在测试数据上的表现如图3a,b所示，显示了预测值与观测到的NO2浓度。结果也与上述预期相符。对于两种方法来说，所得结果似乎都是合理的，支持了当前概念的有效性。XGBR的预测似乎略好一些。图3. 当前的MLR，即RFR（a）和XGBR（b）。在指定用于测试的小时，对42R802站点观测到的NO2浓度的预测。这里很清楚，当前方法实际上是对通过CFD建模得到的给定小时污染物浓度空间分布的校正器，通过考虑额外的定义数据，如(a)气象站观测数据（即温度、相对湿度和大气稳定性）和(b)仅来自一个本地站点的污染物每小时浓度数据。3. 验证练习验证练习涉及2016年4月30日至5月28日期间通过被动采样测量的105个地点的NO2平均浓度。这些地点如图1所示。就地理分布而言，它们在水平和垂直方向上都覆盖了相当的范围。更具体地说，相对于42R801参考站，X轴的范围是-491米到491米，Y轴的范围是-506米到457米。关于Z轴，样本在不同的高度分布如下：(a) 地面附近的四个样本（高度2.5米到3.5米）；(b) 一楼的85个样本（高度3.23–7.6米）；(c) 二楼的11个样本（高度4.34到10.85米）；(d) 三楼的五个样本（高度10.6米）；(e) 六楼的一个样本（高度19米）。每个MLR用于在上述采样期间（即696小时）提供105个地点的NO2每小时浓度，然后计算了整个活动期间的平均值。这些平均值与观测值进行了比较，以寻找一致性。观察这段时间的数据，发现来自选定参考站点42R801的几个NO2每小时浓度数据点缺失。实际上，696个数据中有31个缺失。因此需要应用相应的插补过程。为了简单和实用，使用RFR方法进行了插补处理：现在42R801站点的NO2每小时浓度成为目标变量。关于输入变量，还包括上述气象站观测数据以及日间小时。额外的输入参数是给定的NO2背景每小时浓度（表1），这是NO2浓度时间分布的度量。使用这些数据构建了一个新的RFR预测器。遵循相同的训练/测试策略（即80/20）。获得的结果给出了良好的R2统计值：测试数据的R2 = 0.912，训练数据的R2 = 0.945。将预测的缺失浓度添加到42R801的NO2浓度中，以构建一个完整的数组，用于当前的验证练习。值得注意的是，由于计划使用相同的数据集进行验证阶段的两种方法，因此只应用了一种方法（即RF方法）。核心验证工作关注所选两种机器学习回归（MLR）模型（RFR和XGBR）在预测105个采样点平均NO2浓度方面的表现；每个样本的结果如图4所示。样本按浓度大小顺序排列。有趣的是，两种MLR模型的结果几乎相同。显然，MLR模型取得了显著的改进。更具体地说，在中等浓度位置的预测与观测值之间的一致性非常显著。对于高浓度区域，预测结果与观测趋势相符。考虑到这些位置经常受到局部性强且不规则排放的影响，这里的比较似乎也是合理的。这种改进在全球范围内也有体现，如图5所示，该图展示了所有105个样本的平均浓度对比。另一方面，MLR模型显著修正了计算流体动力学（CFD）模型的预测结果。需要进一步研究CFD模型的低估程度及其受特定CB阵列输入的影响。此外，在高暴露区域，预测峰值与实际趋势一致。关于差异，预测结果符合流动和扩散模型的预期，因为MLR模型是依据该模型进行训练的。图4显示了每个被动采样点测得的平均NO2浓度的MLR模型（RFR和XGBR）预测结果。样本编号是根据测量得到的平均浓度值从低到高进行的。图5显示了所有样本平均NO2浓度的MLR模型（RFR和XGBR）预测结果。值得注意的是，采样位置的高度对预测性能有显著影响。表2显示了各楼层的平均浓度。显然，MLR模型的改进效果显著。第六层测量的差异最为明显，这可能是由于非局部传输造成的，而在更高海拔处这种影响更为显著；然而，由于该高度只有一个样本，因此不应过分重视这一差异。表2列出了不同垂直位置的样本平均NO2浓度（μg/m3）。

4. 结论性意见
本研究提出了一种基于物理原理、由人工智能驱动的方法，用于近实时估计复杂城市微环境中的NO2浓度。该方法结合了有限但物理上一致的CFD参考模拟结果、操作气象数据和监测数据，并使用MLR模型作为校正和适应组件。关键创新在于将基于形态和物理的描述符（XMODEL）嵌入机器学习预测器的特征空间中，从而在保留基本扩散行为的同时，允许根据数据进行残差误差的校正。
在安特卫普FAIRMODE案例研究中的应用表明，随机森林回归（Random Forest Regression）和极端梯度提升（Extreme Gradient Boosting）模型显著提高了浓度估计的准确性，相较于纯扩散模型有明显改进。这种改进在大多数采样高度和空间位置上都是一致的，仅在较高海拔处由于非局部传输效应的增强而出现较大差异。结果证实，所提出的框架有效地起到了物理引导的模型校正器作用，能够捕捉大气稳定性和二次效应等额外因素，而无需进行计算密集型的实时CFD模拟。
该方法满足了两个主要操作要求：（a）适合近实时应用的计算效率；（b）充分利用现有的物理建模知识。尽管本研究关注的是与交通相关的NO2污染物，但其概念方法也可用于解决其他优先污染物（如PM和O3）的相关问题，前提是加入适当的背景处理和化学考虑因素。
还需要强调的是，基于CFD的建模被大幅简化了，包括：（a）假设大气条件为中性，以及忽略了NO2的化学反应；（b）即使在局部流动不稳定或亚小时浓度波动较大的情况下，也应用了RANS方法和风速缩放（方程（1））。需要进一步研究如何有效应对这些问题，这可能涉及CFD RANS方法本身的改进。需要强调的是，当前的基于人工智能的方法有望以校正的方式解决上述不足。未来的工作还应包括系统的多站点验证、引入依赖于稳定性的校正机制以及优化背景浓度处理方法。
将该方法整合到运行中的空气质量管理系统和智慧城市平台中是自然的下一步。应用该方法可以为管理者/操作员提供近实时的响应，而无需依赖复杂的计算资源，从而实现及时干预。总体而言，本研究为基于物理原理的机器学习提供了可靠、高分辨率的实时城市暴露评估的坚实方法论基础。