A system and method for estimating atmospheric pollutants based on a geographically weighted deep forest

Abstract

The application discloses a kind of atmospheric pollutant estimation system and method based on geographic weighted depth forest, it is related to environmental monitoring, atmospheric science and machine learning cross technical field.The application breaks through the shallow application of traditional geographic weighted random forest to geographic weight, first complete multi-source heterogeneous data space-time matching and preprocessing;Geographic weight loss random forest (GLRF) is constructed, generates geographic weight by adaptive kernel function and reconstructs split criterion by embedding loss function in depth, realizes model innovation in combination with adaptive bandwidth and continuous space-time feature coding.Using Monte Carlo resampling to weaken meteorological interference, generate high spatio-temporal resolution meteorological normalized PM2.5 data nationwide;The application solves the problem of insufficient modeling accuracy of non-stationary space at national scale, overcomes the defects of poor adaptation, low fitting and difficult to identify mutation of traditional model, can eliminate meteorological fluctuations, capture human pollution mutation such as heating, and provide support for environmental protection policy evaluation and clean energy transformation.

Description

Technical Field

[0001] This invention relates to the interdisciplinary fields of environmental monitoring, atmospheric science, and machine learning, specifically to a geographically weighted deep forest atmospheric pollutant estimation system and method. Background Technology

[0002] The current mainstream meteorological normalization technology and traditional geographic weighted random forest models in the field of air quality monitoring mainly rely on two paths—numerical simulation and machine learning—to construct environmental particulate matter concentration estimation systems. Traditional chemical transport models achieve physical simulation through high-precision emission inventories, but their real-time policy assessment capabilities are limited by the lag in spatiotemporal dynamic updates and data uncertainty. Ordinary random forest models assume that samples are independent and identically distributed, failing to capture the spatial non-stationary characteristics of pollution. While existing geographic weighted random forests introduce spatial weights, they often treat them as shallow features or use a global fixed bandwidth, making it difficult to adapt to the "dense in the east and sparse in the west" distribution pattern of monitoring stations. In sparse areas, they are prone to sharp drops in accuracy and insufficient generalization ability.

[0003] The core algorithm architecture has inherent bottlenecks. Existing machine learning solutions have a dual limitation in their modeling logic: on the one hand, ordinary random forests ignore the spatial geographic patterns of atmospheric transport, while traditional geographically weighted random forests only use spatial weights as auxiliary features for sample sampling, failing to deeply reconstruct the core logic of decision tree splitting criteria and loss functions, resulting in the inability to adaptively learn the local environmental patterns of different regions from the algorithm's underlying layers; on the other hand, the bandwidth selection mechanism is rigid, relying heavily on empirical values ​​or globally fixed parameters, lacking a dynamic adaptation mechanism for sparse sites, which severely restricts the model's forecast accuracy and stability under complex terrain at the national scale.

[0004] Furthermore, the entire system lacks collaborative integration. Existing technical solutions generally present a modular and fragmented design, with data preprocessing, model training, and result analysis operating independently, lacking end-to-end fusion. Simultaneously, spatiotemporal feature encoding often employs discretized formats, resulting in numerical gaps across years, and no dedicated pollution mutation detection algorithm has been designed, making it impossible to accurately identify concentration mutations driven by anthropogenic emissions such as heating, thus failing to meet the needs for refined and real-time monitoring.

[0005] In summary, a geographically weighted, deep forest atmospheric pollutant estimation system and method are proposed to address the aforementioned problems. (Invention Content)

[0006] 1. The technical problem to be solved by the present invention

[0007] The purpose of this invention is to propose a geographically weighted, depth-based forest atmospheric pollutant estimation system and method to address the following problems existing in the prior art: (1) Algorithm level: Traditional geographic weighted random forests only use geographic weights as shallow features / sampling weights and do not deeply integrate them into the core of the model, so they cannot accurately capture the spatial non-stationarity of air pollution at national-scale stations.

[0008] (2) Adaptability: The bandwidth selection of traditional models depends on empirical values / global fixed parameters, which cannot adapt to the distribution pattern of "dense in the east and sparse in the west" of monitoring stations in my country. The estimation accuracy of sparse areas drops sharply and the generalization ability is insufficient.

[0009] (3) Meteorological interference removal: Existing technologies cannot accurately eliminate short-term meteorological fluctuation interference in large-scale monitoring networks, and it is difficult to obtain anthropogenic emission-related concentration sequences with high spatial representativeness.

[0010] (4) Spatiotemporal coding level: Traditional spatiotemporal feature coding adopts a discretized format, which has a cross-year numerical discontinuity phenomenon, making it impossible for the model to accurately learn the continuous spatiotemporal trend of pollutant emissions.

[0011] (5) Mutation identification level: Existing technologies have not designed a dedicated pollution mutation detection algorithm, which cannot accurately identify concentration mutations driven by human emissions such as heating, and easily confuse meteorological fluctuations with emission mutations.

[0012] (6) System integration level: Existing technical solutions are modular and fragmented, with each link of data preprocessing, model training and result analysis being independent of each other, lacking end-to-end full-link integration, and unable to support efficient processing of large-scale spatiotemporal data at the national scale.

[0013] 2. Technical Solution To achieve the above objectives, the present invention provides the following technical solution: a method for estimating atmospheric pollutants in geographically weighted deep forests, comprising the following steps: S1. Multi-source data acquisition and preprocessing: Acquire multi-source heterogeneous data within the target area, including measured air pollutant data from national monitoring stations (taking PM2.5 as an example), ERA5 reanalysis meteorological driving factor data, and spatiotemporal coding features; perform quality control and outlier removal on the multi-source heterogeneous data, and perform spatiotemporal matching and alignment to construct a training dataset; S2. Geographic weight matrix construction: Calculate the geographic Euclidean distance ||pj-pi|| between any two monitoring stations i and j within the target area to form a distance matrix dist={disti|1≤i≤2028,i∈Z}; Based on the distance matrix and the preset kernel function, construct the geographic weight matrix w', where w' is the normalized standard geographic weight; S3. Geographically Weighted Deep Forest (GLRF) Model Training: Using the meteorological driving factors and continuous spatiotemporal coding features preprocessed in S1 as input features, and measured data of air pollutants (taking PM2.5 as an example) as the target variable, the geographic weight matrix W_Geo is deeply integrated into the splitting criterion of the deep forest decision tree to reconstruct the loss function (different from the traditional geographically weighted random forest application method that only uses geographic weights as sample features / sampling weights); using a parallel computing architecture, the deep forest model after reconstructing the splitting criterion is trained on a high-performance computing cluster. The Lin model is iteratively trained to obtain a trained GLRF model. The construction of the GLRF model includes an adaptive bandwidth selection step to adapt to the distribution of sites nationwide. S4. Meteorological Normalization Processing: Using a trained GLRF model, under fixed time and seasonal parameters, Monte Carlo random resampling is performed on meteorological driving factors to generate a simulated input dataset; the simulated input dataset is input into the GLRF model for prediction, and the expected mean of the predicted values ​​is calculated to obtain the meteorologically normalized atmospheric pollutant (PM2.5) sequence PMdw after eliminating meteorological disturbances; S5. Result Application Analysis: Based on the meteorologically normalized atmospheric pollutant PMdw sequence, the abrupt change characteristics driven by pollution emissions are identified, and the contribution ratio of typical anthropogenic emission sources to the concentration of atmospheric pollutants is estimated; This method achieves meteorological normalized atmospheric pollutant estimation through loss function reconstruction, adaptive bandwidth selection, continuous spatiotemporal feature encoding, sliding window abrupt change detection, and feature combination of full-link system fusion. Existing technologies have not integrated and applied the above features.

[0014] The difference between this technical solution and the existing technology lies in the following: (1) Algorithm; This section presents the core algorithm, which, by constructing a geographic weight matrix and reconstructing the deep forest splitting criterion, achieves precise analysis of the spatial nonstationarity of air pollution. Precise capture is performed to complete GLRF model training. The geographic Euclidean distance ||pj-pi|| between any two monitoring stations i and j within the target area is calculated, forming a distance matrix dist={disti|1≤i≤2028,i∈Z}. Based on the distance matrix and a preset kernel function, a geographic weighting is constructed. The weighted matrix w', where w' is the normalized standard geographic weight; Using preprocessed meteorological driving factors and continuous spatiotemporal coding features as input features, and measured data of air pollutants (taking PM2.5 as an example) as target variables, the geographical weight matrix W_Geo is deeply integrated into the splitting criterion of the deep forest decision tree to reconstruct the loss function (different from the traditional geographically weighted random forest that only uses geographical weights as sample features / sampling weights); using a parallel computing architecture, the deep forest model after reconstructing the splitting criterion is iteratively trained on a high-performance computing cluster to obtain a trained GLRF model. The construction of the GLRF model includes an adaptive bandwidth selection step to adapt to the distribution of stations nationwide; The loss function of the splitting criterion for deep forest decision trees is reconstructed as follows: breaking through the shallow application logic of traditional geographically weighted random forests that only use spatial distance weights as sample features / sampling weights, the spatial distance weights are directly integrated into the core process of node splitting in random forest (RF) decision trees. By minimizing the sum of geographically weighted squared errors, the optimal splitting features and splitting points are determined, so that geographic weights constrain the growth of decision trees from the bottom layer of model training, thereby achieving accurate capture of the spatial nonstationarity of air pollution.

[0015] The reconstructed formula for the geographic weighted splitting criterion is:

[0016] In the formula, , Indicates the left and right subregions after splitting; Represents the normalized geographical weight of sample i; This represents the observed atmospheric pollutant values ​​for sample i; , The weighted representative values ​​of samples in the two sub-regions (preferred) , .

[0017] (2) Adaptability; This section optimizes the calculation of geographic weights, bandwidth selection, and computational architecture to improve model adaptability and processing speed, taking into account the characteristics of the nationwide site distribution and computational efficiency requirements. The geographic weight is ωij, which is calculated using a Gaussian kernel function, and its formula is: w(pj,pi,bt)=exp-||p-pi||² / (2bt²)) In the formula, ||pj-pi|| is the geographic Euclidean distance between samples i and j, bt is the t-th candidate adaptive bandwidth parameter, where t∈{1,2,…,50}, and the calculated original weight w is normalized to obtain the standard geographic weight w'. The construction process of the trained GLRF model also includes an adaptive bandwidth selection step designed for the dense distribution of monitoring stations in the east and sparse distribution in the west across the country. Specifically, this involves: Unlike the globally fixed bandwidth / empirical bandwidth used in traditional geographically weighted random forests, this addresses the problems of poor bandwidth adaptability and low local fitting accuracy in traditional methods. Based on the spatial locations of all monitoring stations nationwide, a distance matrix between stations is constructed, and the median M of this distance matrix is ​​calculated. For all monitoring stations nationwide as a whole, a set of candidate bandwidth parameters b = {bt|1≤t≤50,t∈Z} is set (the candidate bandwidth range is [M / 100,M] of the aforementioned median M, and 50 candidate values ​​are generated at equal intervals). Five-fold cross-validation (5-fold CV) is used to train the model under different candidate bandwidths, and its performance metric, root mean square error (RMSE), is evaluated. The average root mean square error (RMSE) corresponding to each candidate bandwidth is calculated. The candidate bandwidth bt that minimizes the average root mean square error (RMSE) is selected as the optimal bandwidth b_opt for the entire country. Used for final model training. The high-performance computing cluster is configured with no fewer than 256 CPU cores, with each node configured with no fewer than 64GB of DDR5 memory and all-solid-state NVMe storage. The Dask distributed framework is used to synchronously distribute the adaptive bandwidth selection and spatiotemporal cross-validation tasks to each CPU core for execution, thereby shortening the processing time of nationwide scale data.

[0018] (3) Eliminate meteorological interference; This section utilizes Monte Carlo random resampling combined with GLRF model prediction to eliminate meteorological interference and obtain meteorologically normalized pollutant sequences driven solely by anthropogenic emissions. Using a trained GLRF model, under fixed time and seasonal parameters, Monte Carlo random resampling is performed on meteorological driving factors to generate a simulated input dataset. This simulated input dataset is then input into the GLRF model for prediction, and the expected mean of the predicted values ​​is calculated to obtain the meteorologically normalized atmospheric pollutant (PM2.5) sequence PMdw, eliminating meteorological disturbances. The specific steps of Monte Carlo random resampling include: fixing the NSC seasonal coefficient and Unix time parameter in the spatiotemporal coding features to maintain the long-term trend and periodic fluctuations of anthropogenic emissions; performing 1000 independent random samplings on the ERA5 reanalysis meteorological driving factors from a preset historical meteorological distribution to generate 1000 sets of simulated input meteorological factors; inputting these 1000 sets of simulated input datasets into the trained GLRF model to obtain 1000 predicted values ​​of atmospheric pollutants (PM2.5); and calculating the expected mean of these 1000 predicted values ​​to obtain... The final meteorological normalized PMdw sequence.

[0019] (4) Spatiotemporal coding; This section is responsible for acquiring, preprocessing, and constructing spatiotemporal features of multi-source heterogeneous data, providing high-quality, spatiotemporally consistent input data for subsequent model training. It acquires multi-source heterogeneous data within the target area, including measured air pollutant data from national monitoring stations (taking PM2.5 as an example), ERA5 reanalysis meteorological driving factor data, and spatiotemporal coding features. The multi-source heterogeneous data undergoes quality control, outlier removal, and spatiotemporal matching and alignment to construct the training dataset. The spatiotemporal coding features include NSC seasonal coefficients and Unix time parameters to maintain the long-term trend and periodic fluctuations of anthropogenic emissions.

[0020] (5) Mutation identification; This section, based on meteorological normalized sequences, identifies abrupt changes driven by pollution emissions and analyzes the contribution ratio of typical anthropogenic emission sources, providing support for pollution cause assessment. Based on the meteorological normalized air pollutant sequences, such as PMdw sequences, it identifies abrupt changes driven by pollution emissions and estimates the contribution ratio of typical anthropogenic emission sources to air pollutant concentrations. A sliding window abrupt change detection algorithm is used to identify pollution abrupt change dates (AIDs), calculate the growth slope (k), and analyze the trends in the proportion of typical anthropogenic emission sources such as residential heating and industrial emissions.

[0021] (6) System integration; This section integrates the core technology links and functional modules mentioned above into a complete estimation system, achieving a closed loop from data input to result application. This method achieves meteorological normalized air pollutant estimation through loss function reconstruction, adaptive bandwidth selection, continuous spatiotemporal feature encoding, sliding window abrupt change detection, and feature combination of full-link system fusion. Existing technologies have not integrated and applied these features. A meteorological normalized air pollutant estimation system based on geographic weighted deep forest includes: a data input module for acquiring measured air pollutant data from national monitoring stations, ERA5 reanalysis meteorological driving factor data, and spatiotemporal encoding features; a preprocessing module for performing outlier removal, missing value imputation, spatiotemporal matching, and cross-validation partitioning on the multi-source heterogeneous data acquired in the data input module; and a high-performance computing module including a cluster with no fewer than 256 CPU cores and parallel processing. The scheduling engine supports large-scale bandwidth search and model validation tasks; the GLRF model building module calculates the geographic distance matrix, constructs a geographic weight loss function based on adaptive bandwidth selection, and utilizes a high-performance computing module to perform dual reconstruction and iterative training of the loss function and splitting criterion for deep forests, which differs from the model building logic of traditional geographic weighted random forests; the meteorological normalization processing module performs Monte Carlo random resampling and GLRF model prediction under fixed parameters to obtain PMdw sequences; and the results application and analysis module identifies pollution mutation dates (AID) using a sliding window mutation detection algorithm, calculates the growth slope (k), and analyzes the trend of the proportion of typical anthropogenic emission sources such as residential heating and industrial emissions.

[0022] Compared with existing technologies, the geographically weighted depth forest atmospheric pollutant estimation system and method provided by this invention have the following beneficial effects: (1) Algorithm accuracy: The geographically weighted loss function is reconstructed and the spatial weights are deeply integrated into the decision tree node splitting core. By minimizing the geographically weighted squared error and determining the optimal splitting scheme, the model growth is constrained from the bottom layer of the algorithm, accurately capturing the spatial non-stationarity of air pollution. This improves the national PM2.5 estimation R² to 0.90 and reduces RMSE to 12.51 μg / m³, which is significantly better than traditional methods.

[0023] (2) Computational efficiency: For the distribution of stations in the east and west of the country, an adaptive bandwidth selection strategy based on the median of the distance matrix and five-fold cross-validation is designed to match the optimal bandwidth and independent geographical weight matrix for each station, which completely solves the problem of poor adaptability of traditional fixed bandwidth and greatly improves the estimation accuracy of sparse areas such as the northwest and southwest.

[0024] (3) Event Identification: By using continuous spatiotemporal feature coding, a continuous spatiotemporal coding scheme of UNIXTIME, NSC seasonal coefficient, and DAYOFWEEK is designed to replace the traditional discrete coding, solve the problem of cross-year numerical discontinuity, and enable the model to deeply learn the long-term trend, seasonal fluctuation and intra-week cycle characteristics of pollutant emissions, thereby improving the physical continuity of the normalized sequence. The Monte Carlo random resampling method is used to perform 1000 random samplings of meteorological factors after fixing the spatiotemporal coding parameters. The expected mean is obtained by predicting with the GLRF model and calculating the PMdw sequence to eliminate meteorological disturbances. The interference of short-term meteorological fluctuations is accurately removed, and the pure anthropogenic emission-driven concentration data is extracted. Through accurate identification of pollution mutations, a 30-day sliding window mutation detection algorithm is built. The emission mutation date (AID) is locked by fitting the slope k of the concentration change within the window, and the mutation intensity is quantified by calculating indicators such as ΔPM and HID. The mutations of meteorological fluctuations and anthropogenic emission mutations such as heating are accurately distinguished, providing a scientific basis for identifying the contribution of emission sources.

[0025] (4) System Collaboration: Construct an end-to-end, fully integrated system that integrates six modules, including data processing, model building, and result analysis; relying on A high-performance cluster with 256 cores or more uses the Dask distributed framework to achieve parallel computing throughout the entire process, breaking through the computing power bottleneck, significantly shortening the processing cycle of massive spatiotemporal data, and forming a synergistic effect among various innovative features to achieve efficient processing and estimation of data at the national scale. Attached Figure Description

[0026] ] Figure 1 This is a schematic diagram of the overall static architecture of the meteorological normalized atmospheric pollutant (PM2.5) estimation system based on high-performance parallel computing according to the present invention. Figure 2 This is a schematic diagram illustrating the construction principle of the Geographically Weighted Loss Function (GLRF) in an embodiment of the present invention; Figure 3 This is a scatter plot illustrating the comparison of prediction accuracy between the GLRF model of this invention and the conventional random forest model (Baseline). Figure 4 This is a schematic diagram comparing the time series concentrations of air pollutants (taking PM2.5 as an example) before and after meteorological normalization in a typical city in an embodiment of the present invention. Figure 5 This is a schematic diagram illustrating the use of the sliding window algorithm to identify abrupt changes in meteorological normalized concentration (PMdw) in an embodiment of the present invention. Detailed Implementation

[0027] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0028] The present invention will now be described in further detail with reference to the accompanying drawings and embodiments.

[0029] For an example, please refer to [example]. Figures 1 to 5 As shown: To address the problems mentioned in the technical solutions, this application provides a geographically weighted deep forest-based atmospheric pollutant estimation system and method. This embodiment uses PM2.5 data from 2028 national-level air quality monitoring stations across China as the research object. Relying on a high-performance computing cluster with 256 cores or more, it achieves meteorologically normalized PM2.5 concentration estimation and anthropogenic emission-driven pollution abrupt change identification through six core steps: data acquisition and processing, geographic weight matrix construction, GLRF model training, parallel computing implementation, meteorological normalization processing, and emission mutation feature identification. All implementation steps address the technical shortcomings of traditional geographically weighted random forests and integrate five innovative features, including loss function reconstruction and adaptive bandwidth selection. The specific steps are as follows: Step 1: Data acquisition and multidimensional feature parameterization; The shortcomings of existing technologies are: traditional geographically weighted random forests mostly use discrete time feature encoding, which has the problem of numerical discontinuity across years and cannot accurately capture the continuous spatiotemporal trend of pollutant emissions; multi-source data lacks rigorous quality control and the spatiotemporal alignment accuracy is insufficient, which affects the reliability of model input.

[0030] The advantages of this approach are: designing continuous physics-driven temporal feature encoding (UNIXTIME, NSC seasonal coefficient, etc.) to solve the numerical discontinuity problem; constructing multi-source heterogeneous datasets and performing rigorous quality control and spatiotemporal alignment; and improving model training stability and convergence speed through feature standardization, thereby providing high-quality physics-driven input for the GLRF model.

[0031] Summary of steps: This step constructs a high spatiotemporal resolution analysis foundation through multi-source data acquisition, continuous spatiotemporal feature encoding, rigorous data preprocessing, and spatiotemporal alignment, providing reliable and highly continuous input data support for subsequent model training.

[0032] Please refer to Figure 1As shown, this embodiment constructs a high-dimensional spatiotemporal feature matrix to provide physically driven data input for the Geographically Weighted Deep Forest (GLRF) model. This scheme establishes a high spatiotemporal resolution analytical foundation by integrating nationwide monitoring data with reanalysis meteorological data.

[0033] (1) Construction and spatiotemporal alignment of multi-source heterogeneous datasets; The measured pollution dataset covers the entire period from 2013 to 2020. These sites remained relatively stable after 2015, providing a reliable target variable Y for the model to capture long-term emission trends. The adaptive bandwidth selection strategy of this invention is designed specifically for this spatial distribution characteristic, differing from the global bandwidth scheme of traditional geographically weighted random forests.

[0034] Meteorological driving factor set: Simultaneous acquisition of 0.25°×0.25° high-resolution ERA5 reanalysis data released by the European Centre for Medium-Range Weather Forecasts (ECMWF). This invention extracts core physical quantities closely related to the evolution of atmospheric pollutants (taking PM2.5 as an example), specifically including: surface air temperature (T), relative humidity (RH), surface air pressure (MSLP), 10-meter easterly wind component (10U), 10-meter northerly wind component (10V), and planetary boundary... Floor height (BLH).

[0035] Spatial mapping and temporal alignment: Due to the high spatiotemporal continuity of meteorological factors, this invention performs spatiotemporal matching between gridded meteorological data and discrete monitoring stations. In the time dimension, all meteorological factors are uniformly aligned to the daily or hourly resolution of ground observations, ensuring that all input features and measured data of air pollutants (taking PM2.5 as an example) completely correspond in spatiotemporal scale, forming standardized spatiotemporal matching observation pairs.

[0036] (2) Physically driven temporal feature encoding and mathematical conversion; Traditional geographically weighted random forests often employ discrete time feature encoding (such as months and dates), which suffers from numerical discontinuities across years and fails to accurately capture the continuous spatiotemporal trends of pollutant emissions. This invention addresses this deficiency by designing a continuous, physically driven time feature encoding scheme. To enable machine learning models to deeply understand the periodic fluctuations and instantaneous changes in pollutant emissions, this invention does not directly input the original date but reconstructs it into parameterized features with explicit physical meaning: UNIXTIME (Long-Term Trend Feature): Converts the observation time into the cumulative number of seconds since January 1, 1970. As a continuous time variable, this feature allows the model to quantify the nonlinear emission trend of air pollutant (taking PM2.5 as an example) concentrations decreasing year by year due to the implementation of policies such as the Clean Air Action Plan.

[0037] [NSC Seasonal Coefficient (Periodic Fluctuation Characteristics): Addressing the significant seasonality of air pollution, this invention constructs a nonlinear seasonal coefficient (NSC), calculated based on ordinal days of the year, with values ​​set to -1 at the winter solstice and 1 at the summer solstice. This transformation ensures the physical continuity of temporal characteristics, effectively resolving numerical discontinuities at year-ends in traditional date formats, and helping models accurately capture the seasonal patterns of emissions.]

[0038] DAYOFWEEK (Weekly Cycle Feature): Encodes dates as days of the week to characterize the regular differences in human-caused emissions between weekdays and weekends, enhancing the model's ability to identify social activity cycles.

[0039] The continuous spatiotemporal feature encoding designed in this invention differs from the discrete encoding method of traditional geographic weighted random forests. It solves the numerical discontinuity problem, enabling the model to deeply understand the long-term trend, periodic fluctuations and intra-weekly cycle characteristics of pollutant emissions, and providing highly continuous temporal feature inputs for GLRF models.

[0040] (3) Data quality control and standardized preprocessing; To ensure the reliability of the analysis, this invention employs a rigorous data cleaning process: Outlier removal: First, physically invalid values ​​(less than 0 μg / m³ or greater than 1000 μg / m³) are excluded. For the remaining data, outliers are identified within a 15-day sliding window using the ±3 standard deviation rule, and then globally verified using the 0.1% and 99.9% quantiles.

[0041] Missing value interpolation: For identified noise and short-term missing values, linear interpolation is used to fill in the missing values. The maximum interpolation window is set to 15 days. Missing values ​​exceeding this range will not be filled in to prevent the introduction of highly uncertain synthetic data.

[0042] Feature standardization: The model performs standardization preprocessing for meteorological factors with different dimensions. This process ensures that each feature has a consistent influence weight in numerical calculation, which can significantly improve the numerical stability and convergence speed of the subsequent geographic weighted deep forest algorithm during the training phase.

[0043] Step 2: Construction of the geographic weight matrix and characterization of spatial heterogeneity; The shortcomings of existing technologies are: traditional geographically weighted random forests use global fixed bandwidth or empirical bandwidth, which cannot adapt to the distribution characteristics of monitoring stations in my country, which are dense in the east and sparse in the west, resulting in poor local fitting or excessive smoothing; geographic weights are only applied superficially and have not formed a complete bandwidth optimization and parallel acceleration logic, making it difficult to support nationwide-scale calculations.

[0044] Advantages of this solution: It designs an adaptive bandwidth selection strategy based on the median of the distance matrix and five-fold cross-validation to accurately adapt to the site distribution characteristics; it maps geographical distance to non-linear weights through the Gaussian kernel function and combines it with the Dask distributed framework to achieve parallel acceleration of bandwidth optimization for 2028 sites across the country, which greatly improves computational efficiency and spatial fitting accuracy.

[0045] Summary of steps: This step constructs a geographic weight matrix that adapts to the distribution of stations nationwide through geographic Euclidean distance calculation, adaptive bandwidth selection, geographic weight mapping, and parallel bandwidth optimization. This matrix accurately represents the spatial non-stationarity of air pollution and provides a spatial constraint basis for the GLRF model.

[0046] Please refer to Figure 2 As shown, it breaks down the mathematical logic steps of geographic distance calculation, adaptive bandwidth selection, geographic weight mapping, and loss function reconstruction in detail. This construction principle breaks through the shallow application logic of traditional geographic weighted random forests and is a unique design specific to this invention. This embodiment constructs a spatial geographic weight matrix, transforming the spatial location characteristics of each monitoring station into numerical weights that the model can extract, thereby accurately characterizing the spatial non-stationarity of air pollution.

[0047] (1) Precise calculation of the geographic Euclidean distance matrix; This invention first constructs a spatial location set p = {pk|1≤k≤2028,k∈Z}. Let station pi be the target modeling station. The system needs to calculate the geographic Euclidean distance ||pj-pi|| between the surrounding monitoring stations pj and the modeling station pi within the target area, where the station indices i,j∈{1,2,…,2028}. For each target modeling station, the system needs to iterate and calculate its physical distance to all other monitoring stations nationwide. This calculation process ultimately produces a geographic Euclidean distance matrix D with dimensions 2028×2028, where the element in the i-th row and j-th column is ||pj-pi||. Simultaneously, each row of the matrix is ​​extracted as a distance vector disti (disti ∈ D, i-th row) with the i-th station as the modeling center, denoted as dist = {disti | 1 ≤ i ≤ 2028, i ∈ Z}. Here, disti represents the distance vector from the i-th station to all its neighboring stations when the i-th station is the modeling center, serving as the mathematical basis for the derived space weights.

[0048] (2) Adaptive bandwidth selection strategy; To address the shortcomings of traditional geographically weighted random forests, which use a fixed global bandwidth and cannot adapt to the sparse distribution of monitoring stations in the east and west, this invention designs an adaptive bandwidth selection strategy based on the median of the distance matrix and five-fold cross-validation. Traditional geographically weighted random forests often use a fixed global bandwidth or empirical bandwidth, which cannot adapt to the distribution characteristics of monitoring stations in my country, which are densely distributed along the eastern coast and sparsely distributed in the western inland areas, leading to poor local fitting or over-smoothing. This invention addresses this deficiency with a dedicated adaptive bandwidth selection strategy. The bandwidth parameter determines the rate at which geographic weights decay with increasing distance. A larger bandwidth results in a smoother weight decay, and the model considers a wider range of spatial influences; a smaller bandwidth results in a more drastic weight decay, and the model focuses more on local features. Since the distribution of monitoring stations in China is extremely uneven, with a dense distribution along the eastern coast and a sparse distribution in the western inland areas, this invention employs adaptive bandwidth selection logic to apply appropriate local weights. Search space definition: Let M be the median (Median) in the distance matrix disti of the current target site.

[0049] Candidate bandwidth generation: Within the range [M / 100, M], 50 candidate bandwidth parameters are generated using an equal-interval strategy, denoted as b={bt|1≤t≤50,t∈Z}. This range ensures that the candidate bandwidths can capture both very close-range strong spatial dependencies and cover medium-scale spatial correlations.

[0050] (3) Geographic weight mapping based on Gaussian kernel function; This embodiment uses a Gaussian kernel function to map the physical distance between sites into nonlinear weights required for model construction. Under the t-th candidate bandwidth bt, the original geographic weight w(pj,pi,bt) of the surrounding site j relative to the target modeling site i is calculated as follows: w(pj,pi,bt)=exp(-||pj-pi||² / (2bt²)) In the formula, ||pj-pi|| represents the Euclidean distance between two points. The calculated original weights w(pj,pi,bt) are normalized to obtain the standard geographic weights w'ij (abbreviated as w'). These weights ensure that spatially closer data points contribute more to model construction, thus making the geographic context the basis of model building.

[0051] Step 2.4: Bandwidth optimization and parallel acceleration based on five-fold cross-validation; To determine the optimal bandwidth that best represents the spatial non-stationary characteristics of the target site, this embodiment executes the following automated optimization process: Cyclic validation: The 50 generated candidate bandwidths are sequentially substituted into the improved geographically weighted deep forest model (GLRF) for testing.

[0052] Spatiotemporal cross-validation: Five-fold cross-validation (5-fold CV) is used to calculate the average root mean square error (RMSE) of the model on the validation set for each candidate bandwidth.

[0053] Optimal parameter locking: Select the candidate bandwidth bt (denoted as b_opt) and its corresponding standard weight w' that minimizes the average RMSE as the final modeling parameters for the site.

[0054] Because high-frequency bandwidth cycling and cross-validation are required for all 2028 sites nationwide, the total computational load is enormous. This invention relies on a high-performance computing cluster equipped with more than 256 CPU cores (each node configured with at least 64GB of DDR5 memory and all-solid-state NVMe storage), and employs the Dask asynchronous distributed computing framework as its parallel scheduling engine. The system independently distributes bandwidth search tasks for different sites to each processor core, utilizing massive parallel computing power to achieve rapid construction and optimization of adaptive spatial weights at a nationwide scale.

[0055] The adaptive bandwidth optimization logic of "median of distance matrix benchmark + five-fold cross-validation" designed in this invention is a dedicated design to adapt to the distribution of sites at the national scale. It is different from the fixed / empirical bandwidth scheme of traditional geographic weighted random forest. Combined with the parallel acceleration of the Dask distributed framework, it realizes the rapid construction and optimization of adaptive spatial weights at the national scale, and solves the core problem of poor adaptability of traditional methods.

[0056] Step 3: Reconstruct the loss function of Geographically Weighted Deep Forest (GLRF); The shortcomings of existing technologies: Traditional geographically weighted random forests only use geographic weights as sample features or sampling weights, without involving the core construction of decision tree splitting criteria and loss functions. They cannot constrain model growth from the bottom up and are difficult to accurately capture the spatial non-stationarity of air pollution.

[0057] The advantages of this solution are: it reconstructs both the loss function and the splitting criterion of deep forests, deeply integrates geographical weights into the core process of decision tree node splitting, minimizes geographical weighted squared error and selects the optimal splitting point; combined with high-performance parallel integration and model optimization, it enables independent modeling of sites nationwide, significantly improves model fitting accuracy, and especially optimizes the performance of sparse areas of sites.

[0058] Summary of steps: This step uses random forest as the basic architecture, reconstructs the embedded geographic weight constraints through the loss function, and builds multiple geographic weighted decision trees in parallel to form a GLRF model. This solves the core defect of traditional models that are not good at capturing spatial non-stationarity, and improves the model fitting accuracy and computational efficiency.

[0059] Traditional geographically weighted random forests only use geographic weights as input features or sampling weights, without involving the core construction process of decision tree splitting criteria and loss functions. This fails to constrain model growth from the bottom up to capture spatial non-stationarity. This invention addresses this core deficiency by doubly reconstructing the loss function and splitting criteria of deep forests to construct a geographically weighted deep forest model. Please refer to... Figure 2 As shown, this invention uses RandomForest as its basic architecture and reconstructs the loss function of its core component, RegressionTree, to transform spatial location attributes into the driving force for model training, thereby constructing a geographically weighted deep forest (GLRF) model.

[0060] (1) Geographically weighted logic of the loss function in regression decision trees; Commercially available geographically weighted random forests only use geographic weights for sample selection or model output weighting, without involving the core splitting process of regression decision trees. This fails to fundamentally reflect the constraining effect of spatial weights on model training. This invention makes a substantial improvement in this regard. Regression decision trees are the basic building blocks of deep forest models. In this embodiment, by assigning normalized geographic weights w'ij (where w'ij is the normalized geographic weight of sample i relative to modeling station j) to each training data point as determined by the steps described in Chapter 2, the model's growth process is effectively constrained by spatial non-stationarity. Let the input feature vector be X, and the output target variable be Y (i.e., the mass concentration of air pollutants (taking PM2.5 as an example)). Given a training dataset D={(x1,y1),(x2,y2),…,(x_N,y_N)}, where xi is a feature vector containing meteorological factors and spatiotemporal encoding, yi is the corresponding observation instance, and N is the total number of samples. Let w' = {w'ij | 1 ≤ i, j ≤ 2028, i ∈ Z, j ∈ Z} be the standard geographic weight matrix calculated for the current target modeling site. This weight matrix reflects the spatial proximity of each data point relative to the current modeling center. Taking CART (Classification and Regression Trees) regression trees as an example, constructing a geographically weighted loss forest mainly includes the selection of optimal split points, The process includes region division and decision tree generation.

[0061] (2) Optimal split point selection based on geographically weighted sum of squared errors; The system employs a heuristic method to traverse all feature variables j and their corresponding split points s, performing optimal partitioning of the feature space. For a feature variable x(j) and its value s in the training set, the current sample set is divided into two sub-regions R1(j,s)={x|x(j)≤s} and R2(j,s)={x|x(j)>s}. Traditional regression decision trees typically use minimizing the squared error as the splitting criterion. minj,s[min_c1Σ(xi∈R1(j,s))(y-c1)²+min_c2Σ(xi∈R2(j,s))(y-c2)²] The reconstructed geographic weighted loss function and splitting criterion of this invention deeply integrate geographic weights into the core process of decision tree node splitting, which differs from the shallow application logic of traditional geographic weighted random forests. This formula is a unique design of this invention, and there is no publicly available solution in the prior art that embeds geographic weights into this loss function and uses it as the core basis for the splitting criterion. The core improvement of this invention is that when searching for the optimal feature variable j and split point s, geographic weight w'i is introduced into the objective function to reconstruct the geographicallyweighted loss function, specifically expressed as formula (1):

[0062] In the formula, and The weighted representative value (predicted output value) of the samples in the two sub-regions is calculated as the arithmetic mean of yi in the corresponding regions: The system iterates through all feature variables j and scans all possible split points s for each variable, selecting the parameter pair (j,s) that minimizes the aforementioned geographically weighted loss function. This improvement ensures that, during model growth and splitting, neighboring site samples with higher geographical weights have a greater influence on the splitting decision, thus accurately representing local spatial relationships.

[0063] (3) Region division and output value determination (sub-process two); Using the selected optimal parameters (j,s) mentioned above, the feature space is divided into regions, and the corresponding output value for each sub-region is determined. :

[0064] In the formula, N m This represents the number of samples that fall into the subregion Rm.

[0065] (4) Recursive construction and decision tree generation (sub-processes three and four); The system recursively calls the optimal split point selection steps for the two divided sub-regions until a preset stopping condition is met (such as reaching the maximum tree depth, the number of samples within a node being less than a threshold, or the error descent gradient being lower than a preset value). Finally, the high-dimensional input feature space is divided into M non-overlapping regions R1, R2, ..., R_M, and a single geographically weighted decision tree model f(x) is generated.

[0066] In the formula, I is the indicator function, which is used when sample x falls into region R. m The value is 1 if it is true, and 0 otherwise.

[0067] (5) High-performance parallel integration and model optimization; In a high-performance computing cluster environment (256 or more CPU cores, with a single node configured with no less than 64GB DDR5 memory and all-solid-state NVMe storage), the system randomly extracts sample and feature subsets using bootstrapping, and constructs multiple independent geographically weighted decision trees in parallel to form a GLRF model. This embodiment sets the number of trees to 500, using 80% of the dataset for training and 20% for testing. Since each site has an independent geographical weight matrix w' and an adaptive optimal bandwidth b_opt (b_opt is the candidate bandwidth bt corresponding to the minimum RMSE selected through five-fold cross-validation), the system uses the Dask distributed framework to distribute the modeling tasks for 2028 sites to the computing cores. This parallelized loss function calculation and parameter optimization logic enables the model to maintain extremely high computational efficiency and spatial fitting accuracy even with large-scale spatiotemporal data. Selecting the optimal bandwidth through five-fold cross-validation (CV) significantly improves the model's performance in data-sparse regions (such as the Northwest and Southwest regions). This invention configures an independent geographic weight matrix and adaptive optimal bandwidth for each site, and combines it with the Dask distributed framework to achieve parallel modeling. This optimization effect is something that traditional geographic weighted random forests cannot achieve due to the lack of dedicated bandwidth design and full-link parallel architecture. In particular, the model fitting accuracy is significantly improved in sparsely populated areas such as the Northwest and Southwest.

[0068] Step 4: Parallel implementation scheme based on HPC clusters with 256 cores or more; The shortcomings of existing technologies: the parallel computing of traditional geographic weighted random forests is only for a single modeling stage, the modules lack coordination, and there is no dedicated memory optimization and task scheduling logic designed for large-scale spatiotemporal data at the national scale. It cannot support high-frequency iteration and parameter optimization, resulting in low computational efficiency.

[0069] Advantages of this solution: It designs a parallel implementation scheme adapted to the entire GLRF model chain, relies on HPC clusters with 256 cores or more and Dask asynchronous distributed computing framework to distribute site-dimensional modeling tasks in parallel; it adopts a shared memory mechanism to optimize memory management and combines a parallel acceleration library to improve computing efficiency, which can efficiently support the processing of large-scale data spanning 2028 sites across the country and ten years.

[0070] Summary of steps: This step optimizes the distributed computing architecture, parallel task scheduling, and memory management to achieve full-process parallelism in data processing, bandwidth optimization, and model training, providing computing power support for the large-scale and efficient implementation of the GLRF model and verifying the model's superiority.

[0071] Traditional geographic weighted random forests (GLRFs) rely on parallel computation for only a single modeling stage, and the modules lack coordination, making them unsuitable for efficient processing of large-scale spatiotemporal data at a national scale. This invention addresses this deficiency by designing a parallel implementation scheme adapted to the entire GLRF model chain, achieving parallel processing of data processing, bandwidth optimization, and model training. (Refer to...) Figure 1 As shown, this invention designs a parallel implementation scheme based on a high-performance computing cluster (HPC) for a large-scale hourly dataset spanning 2028 sites across the country over a ten-year period, in order to support high-frequency iteration and parameter optimization of the GLRF model.

[0072] Distributed computing architecture and environment configuration. This solution is implemented on a hardware computing array equipped with more than 256 central processing unit (CPU) cores. Each node is configured with no less than 64GB of DDR5 memory and all-solid-state NVMe storage to meet the high bandwidth requirements of geographic distance matrix calculation and high-frequency random forest split search. The software environment uses the scientific computing stack of Python 3.9 or above and integrates the Dask asynchronous distributed computing framework.

[0073] Parallelized task scheduling logic. Since the GLRF model requires independently calculating geographic weights for each site and performing five-fold cross-validation on 50 candidate bandwidths, its total computational load is approximately 2028 × 50 × 5 model evaluations. This invention encapsulates the site-level modeling tasks into independent computational operators using a parallel scheduling engine, and distributes them to idle cores within the cluster for synchronous execution using Dask's distributed scheduler.

[0074] Memory optimization and middleware management. During parallel training, the system employs a shared memory mechanism to manage ERA5 global meteorological field data, effectively avoiding the risk of memory overflow caused by repeated data loading. To address the frequent calculations of the geographically weighted loss function, the system utilizes a parallel acceleration library to optimize matrix multiplication operations, ensuring the real-time performance and accuracy of splitting criterion search with millions of samples.

[0075] ] Figure 3 The scatter plots comparing the prediction accuracy of the GLRF model provided by this invention with that of the conventional random forest model (Baseline) are as follows: (a) GLRF, (b) the model integrating the spatial-temporal features of STRF into the GeoWeightedLossForest framework (GeoSTRF), (c) spatial-temporal random forest (STRF), (d) conventional random forest (Joint-RF) used for joint training, and (e) conventional random forest (Single-RF) used for single-point training. The figures demonstrate the superiority of the model (a) GLRF of this invention in improving the accuracy of meteorological normalization through metrics such as R2 and RMSE.

[0076] The results demonstrate that the GLRF model of this invention significantly outperforms traditional geographically weighted random forests and other conventional random forest models, validating the innovative value of loss function reconstruction and adaptive bandwidth selection. Leveraging its core design of geographically weighted loss function reconstruction and adaptive bandwidth optimization, the GLRF model achieves significantly higher R² improvements compared to all benchmark models such as Joint-RF, Single-RF, STRF, and GeoSTRF in regions with low data density and sparse site distribution. This fully reflects the accurate ability of the geographic weighting mechanism to capture spatial non-stationarity, effectively compensating for the insufficient fitting accuracy of traditional models in sparsely distributed areas, and its performance advantage is particularly prominent in data-scarce scenarios.

[0077] ] Figure 3 As can be seen, thanks to the core design of geographic weighted loss function reconstruction and adaptive bandwidth optimization, the GLRF model achieves a significantly higher R² improvement compared to all benchmark models such as Joint-RF, Single-RF, STRF, and GeoSTRF in areas with low data density and sparse site distribution. This fully demonstrates the accurate ability of the geographic weighting mechanism to capture spatial non-stationarity and effectively makes up for the shortcomings of traditional models in terms of insufficient fitting accuracy in sparse site areas. Its performance advantage is particularly prominent in data-scarce scenarios.

[0078] Step 5: Weather normalization based on Monte Carlo resampling (De-weathering); The limitations of existing technologies: Traditional geographic weighted random forests and existing meteorological normalization methods cannot completely eliminate the short-term interference of meteorological factors on pollutant concentrations, and it is difficult to extract concentration components driven purely by human emissions, thus affecting the accuracy of pollution attribution.

[0079] The advantages of this approach are: by using a trained GLRF model to construct a simulated scenario through fixed time feature parameters, performing 1000 Monte Carlo random resampling operations on meteorological driving factors, and obtaining the expected mean of the predicted values ​​to effectively remove meteorological interference and obtain an accurate anthropogenic emission driving concentration sequence (PMdw).

[0080] Summary of steps: This step uses Monte Carlo resampling combined with GLRF model inference to eliminate short-term interference from meteorological factors and extract pollutant concentration sequences driven purely by anthropogenic emissions, providing a reliable data foundation for subsequent identification of mutation features.

[0081] This invention utilizes a trained GLRF model to eliminate short-term interference from meteorological factors on the concentration of air pollutants (taking PM2.5 as an example) through statistical simulation, thereby extracting concentration components purely driven by emissions.

[0082] Parameter fixing and simulation scenario construction. During normalization, the system first locks in the time characteristic parameters. By keeping UNIXTIME (long-term trend), NSC seasonal coefficient (cyclical fluctuation), and DAYOFWEEK (weekly cycle) constant at the target time, the system represents the anthropogenic emission level and background environment at that specific moment.

[0083] Monte Carlo random resampling logic. Based on the historical meteorological probability distribution of the station location, the system performs 1000 independent random samplings on the ERA5 reanalysis meteorological driving factors, including surface air temperature (T), wind field (10U, 10V), relative humidity (RH), surface air pressure (MSLP), and boundary layer height (BLH). The sampling process covers the joint distribution of each factor, generating one thousand sets of simulated input vectors representing different possible meteorological scenarios.

[0084] The mathematical expectation of normalized concentration (PMdw) is calculated as follows: One thousand sets of simulated scenarios are input into the optimized GLRF model for inference and prediction, yielding one thousand corresponding predicted values ​​for air pollutants (taking PM2.5 as an example). The final meteorological normalized concentration is obtained by calculating the expected mean of the predicted values, using the following formula: PM=(1 / 1000)×Σk=1¹000GLRF(Meteorology,Time) In the formula, PM represents the meteorologically normalized concentration sequence. This method achieves "de-weather-based" processing of the pollution sequence by averaging all possible meteorological fluctuations, enabling it to intuitively reflect the dynamic changes in anthropogenic emission intensity. Figure 4 This paper presents a comparison of time series concentrations of air pollutants (taking PM2.5 as an example) before and after meteorological normalization in typical cities according to embodiments of the present invention. Taking cities such as Harbin and Beijing as examples, it demonstrates the stable emission trend of normalized concentration PMdw after eliminating short-term meteorological fluctuations. These results prove that the meteorological normalization process of the present invention can effectively remove meteorological interference and obtain a pure emission-driven concentration series, overcoming the deficiency of traditional geographically weighted random forests in completely removing meteorological interference.

[0085] Step Six: Identification and Application Analysis of Emission-Driven Mutation Characteristics; The shortcomings of existing technologies: Traditional geographic weighted random forests and existing meteorological normalization methods have not designed dedicated mutation detection algorithms for emission-driven processes. They mostly use fixed threshold methods or simple trend fitting, which are easily affected by meteorological fluctuations, resulting in inaccurate identification of mutation points and inaccurate attribution of emission contributions.

[0086] Advantages of this solution: Based on the meteorological normalized PMdw sequence, a dedicated sliding window mutation detection (AID) algorithm is designed to accurately pinpoint the trigger date of anthropogenic emission mutations, quantify the mutation intensity, and, combined with industry emission correlation analysis, achieve accurate attribution of anthropogenic emission source contributions, providing a basis for decision-making in environmental policy assessment.

[0087] Summary of steps: This step uses a dedicated sliding window mutation detection algorithm to identify the characteristics of pollution mutations driven by human emissions, quantify the mutation intensity, and complete the attribution of emission source contributions. This solves the problems of inaccurate mutation identification and inaccurate attribution in traditional methods, realizing the practical application value of the technical solution.

[0088] Traditional geographic weighted random forests and existing meteorological normalization methods lack dedicated algorithms for detecting emission-driven abrupt changes. They often employ fixed threshold methods or simple trend fitting, making them susceptible to meteorological fluctuations and resulting in inaccurate abrupt change identification. This invention addresses this deficiency by designing a sliding window abrupt change detection algorithm based on meteorologically normalized PMdw sequences. Based on the generated PM sequences, this invention designs an algorithm to accurately capture anthropogenic emission abrupt changes, used to assess the emission response effects of environmental policies or sudden social events.

[0089] The system employs a sliding window algorithm for mutation detection (AID). A 30-day sliding window moves across the PM sequence. Within each window, a linear regression equation is fitted using the least squares method to extract the slope value k of the concentration change over time. PM = k·t + b In the formula, t represents the normalized time dimension coordinate. By searching for the point in the year when the slope k reaches its maximum peak, the system automatically identifies the trigger date for a significant shift in anthropogenic emissions, i.e., the mutation trigger date (AID).

[0090] Emission-driven mutation intensity quantification (ΔPM). After determining AID, the system further quantifies the emission mutation intensity. By calculating the difference in mean PM between the heating season and the summer season (ΔPM), and indicators such as the number of days with high concentration (HID), the interference of weather factors is eliminated, and the mutation intensity driven purely by emissions is defined.

[0091] Attribution of emissions from typical industries (such as heating). Please refer to [link / reference]. Figure 5 As shown in the diagram, this embodiment of the invention utilizes a sliding window algorithm to identify abrupt changes in meteorological normalized concentration (PMdw); c1: Harbin, c2: Beijing, c3: Hangzhou, c4: Guangzhou, etc. Figure 5 (c1) The mutation detection results in Harbin show that the maximum value of the slope k pinpoints the specific identification logic of the emission mutation trigger date (AID). This algorithm is a dedicated design for emission-driven mutation detection, differing from the general detection logic of traditional methods, and can accurately pinpoint anthropogenic emission mutation points such as heating. The blue vertical line represents the actual trigger date of the detected emission mutation. This invention correlates the AID identification results with the winter municipal heating dates in northern China. By comparing the AID triggering order and concentration increase at stations at different latitudes, important information related to anthropogenic emissions during winter heating was determined. This provides the government with a decision-making basis supported by physical mechanisms for assessing the effectiveness of emission reduction policies and optimizing energy structure adjustments.

[0092] The sliding window mutation detection algorithm designed in this invention differs from the general detection logic of traditional methods. It can accurately eliminate the interference of meteorological fluctuations, lock the date of pollution mutation driven by human emissions, and, combined with the emission source proportion analysis, provide an accurate basis for assessing the contribution of human sources such as heating and industrial emissions. This effect cannot be achieved by traditional geographically weighted random forests.

[0093] Please refer to the above work process. Figures 1 to 5 .

[0094] It should be noted that the term "comprising" or any other variation thereof is intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitation, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes said element.

[0095] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and variations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for estimating atmospheric pollutants in geographically weighted, deep forests, characterized in that: Includes the following steps: S1. Multi-source data acquisition and preprocessing: Acquire multi-source heterogeneous data within the target area, including measured air pollutant data from national monitoring stations, ERA5 reanalysis meteorological driving factor data, and spatiotemporal coding features; The multi-source heterogeneous data is subjected to quality control, outlier removal, and spatiotemporal matching and alignment to construct a training dataset. S2. Geographic Weight Matrix Construction: Calculate the geographic Euclidean distance ||p between any two monitoring stations i and j within the target area. j p i ||, forming a distance matrix dist, where dist = {dist i | 1 ≤ i ≤ 2028, i∈ Z}; Based on the distance matrix and the preset kernel function, construct the geographic weight matrix w', where w' is the normalized standard geographic weight; S3. Geographically Weighted Deep Forest (GLRF) Model Training: Using the meteorological driving factors and continuous spatiotemporal coding features preprocessed in S1 as input features, and measured air pollutant data as target variables, the geographic weight matrix W_Geo is deeply integrated into the splitting criterion of the deep forest decision tree to reconstruct the loss function; using a parallel computing architecture, the deep forest model after reconstructing the splitting criterion is iteratively trained on a high-performance computing cluster to obtain the trained GLRF model. The construction of the GLRF model includes an adaptive bandwidth selection step to adapt to the distribution of stations nationwide; S4. Meteorological Normalization Processing: Using a trained GLRF model, under fixed time and seasonal parameters, Monte Carlo random resampling is performed on meteorological driving factors to generate a simulated input dataset. This simulated input dataset is then input into the GLRF model for prediction. The expected mean of the predicted values ​​is calculated to obtain the meteorologically normalized atmospheric pollutant (PM) sequence after eliminating meteorological disturbances. dw ; S5. Results Application and Analysis: Based on the meteorological normalized atmospheric pollutant PM2.5... dw Sequences were analyzed to identify mutational features driven by pollution emissions and to estimate the contribution of typical anthropogenic emission sources to atmospheric pollutant concentrations.

2. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 1, characterized in that, The splitting criterion loss function of the deep forest decision tree described in S3 is reconstructed as follows: breaking through the shallow application logic of traditional geographically weighted random forests that only use spatial distance weights as sample features / sampling weights, the spatial distance weights are directly integrated into the core process of node splitting in the random forest (RF) decision tree. By minimizing the sum of geographically weighted squared errors, the optimal splitting features and splitting points are determined, so that the geographic weights constrain the growth of the decision tree from the bottom layer of model training, thereby achieving accurate capture of the spatial non-stationarity of air pollution.

3. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 2, characterized in that, The reconstructed geographically weighted splitting criterion formula is as follows: minutes j , s [Σ(x i ∈R1(j,s)) w' i (y i 1)² + Σ(x i ∈R2(j,s)) w' i (y i 2)² ] In the formula, R1(j,s) = {x | x (j) ≤ s}、R2(j,s) = {x | x (j) >s} represents the left and right subregions after splitting; w' i y represents the normalized geographic weight of sample i; i This represents the observed values ​​of air pollutants for sample i; 1. 2 represents the weighted representative value of the samples in the two sub-regions, respectively.

4. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 3, characterized in that, The geographic weight mentioned in S2 is ω ij The ω ij The Gaussian kernel function is used for calculation, and its formula is as follows: w (p j , p i , b t ) = exp ( ||p j p i ||² / (2b t ²)) In the formula, ||p j p i || represents the geographic Euclidean distance between samples i and j, b t Let w be the t-th candidate adaptive bandwidth parameter, where t ∈ {1,2,…,50}. The calculated original weight w is normalized to obtain the standard geographic weight w'.

5. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 4, characterized in that, The construction process of the trained GLRF model described in S3 also includes an adaptive bandwidth selection step designed for the distribution characteristics of the national sites, which are denser in the east and sparser in the west. Specifically, this includes: S3.1 Based on the spatial location of all monitoring stations nationwide, construct a distance matrix between stations and calculate the median M of the distance matrix; S3.

2. For all monitoring stations nationwide as a whole, set a set of candidate bandwidth parameters b = {b t | 1 ≤ t ≤ 50, t ∈ Z} (The candidate bandwidth range is [M / 100, M] of the median M above, and 50 candidate values ​​are generated at equal intervals); S3.

3. The five-fold cross-validation method is used to train the model under different candidate bandwidths and evaluate its performance index, root mean square error (RMSE). The average root mean square error (RMSE) corresponding to each candidate bandwidth is calculated. S3.4 Select the candidate bandwidth b that minimizes the root mean square error (RMSE). t b_opt, the nationally unified optimal bandwidth, is used for final model training.

6. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 5, characterized in that, The high-performance computing cluster described in S3 is configured with no less than 256 CPU cores, with each node configured with no less than 64GB of DDR5 memory and all-solid-state NVMe storage. The Dask distributed framework is used to synchronously distribute the adaptive bandwidth selection and spatiotemporal cross-validation tasks to each CPU core for execution, so as to shorten the processing time of nationwide data.

7. The method for estimating atmospheric pollutants in geographically weighted deep forests according to claim 1, characterized in that, The specific steps of Monte Carlo random resampling described in S4 include: S4.

1. Fixed spatiotemporal coding features include NSC seasonal coefficient and Unix time parameter to maintain the long-term trend and periodic fluctuations of anthropogenic emissions; S4.2 From the preset historical meteorological distribution, perform 1000 independent random samplings on the meteorological driving factors of ERA5 reanalysis to generate 1000 sets of simulated input meteorological factors; S4.

3. Input the 1000 sets of simulated input datasets into the trained GLRF model to obtain 1000 air pollutants (PM2.5). 2.5 Predicted value sequence; S4.4 Calculate the expected mean of the 1000 predicted value sequences to obtain the final meteorological normalized PM2.

5. dw sequence.

8. A meteorological normalized atmospheric pollutant estimation system based on geographically weighted deep forest, used to implement the method according to any one of claims 1-7, characterized in that, include: Data input module: used to acquire measured air pollutant data from national monitoring stations, ERA5 reanalysis meteorological driving factor data, and spatiotemporal coding features; Preprocessing module: used to perform outlier removal, missing value imputation, spatiotemporal matching, and cross-validation partitioning on the multi-source heterogeneous data obtained from the data input module; High-performance computing module: including a cluster of no less than 256 CPU cores and a parallel scheduling engine to support large-scale bandwidth search and model validation tasks; The GLRF model building module is used to calculate the geographic distance matrix, construct the geographic weight loss function based on adaptive bandwidth selection, and use the high-performance computing module to perform dual reconstruction and iterative training of the loss function and splitting criterion for deep forests, which is different from the model building logic of traditional geographic weighted random forests. Meteorological normalization processing module: used to perform Monte Carlo random resampling and GLRF model prediction under fixed parameters to obtain PM2.5 concentration. dw sequence; Results Application Analysis Module: Used to identify pollution mutation dates (AID) using a sliding window mutation detection algorithm, calculate the growth slope (k), and analyze the trend of the proportion of typical anthropogenic emission sources such as residential heating and industrial emissions.

9. The meteorological normalized atmospheric pollutant estimation system based on geographically weighted deep forest as described in claim 8, characterized in that, Each module forms a closed-loop collaborative relationship based on data flow and computational dependencies: The output of the data input module is connected to the input of the preprocessing module, transmitting multi-source heterogeneous data to the preprocessing module to complete quality control and spatiotemporal matching; The output of the preprocessing module is connected to both the high-performance computing module and the GLRF model building module. The standardized dataset is distributed to the high-performance computing module for parallel scheduling and transmitted to the GLRF model building module for geographic weight calculation and model training. The scheduling output of the high-performance computing module is connected to the GLRF model building module, providing parallel computing power support for adaptive bandwidth selection and model iteration; The model output of the GLRF model building module is connected to the meteorological normalization processing module, which transmits the trained GLRF model to the meteorological normalization processing module for Monte Carlo resampling prediction. The PMdw sequence output of the meteorological normalization processing module is connected to the result application analysis module, which transmits the meteorological normalized pollutant sequence to the result application analysis module for pollution mutation identification and emission source contribution ratio estimation. Each module works together through data interfaces and computing power to form an end-to-end execution relationship from data input, preprocessing, model building, meteorological normalization to result analysis, ensuring the continuity and accuracy of national-scale air pollutant estimation.

Images