A PM based on NGBoost algorithm 2.5 Concentration long time series prediction method

Abstract

The present invention discloses a PM based on NGBoost algorithm 2.5 The method combines the influence of meteorology and other pollutants to propose a PM2.5 concentration long series prediction method based on natural gradient improvement. 2.5 Concentration prediction model to improve long-term PM 2.5 The accuracy of concentration prediction. The Isolation Forest algorithm was used to determine the outliers of the air pollutants and meteorological data collected from the target site, and the Spearman correlation coefficient was used for feature selection and other preprocessing as the total data set; a data time series discretization method was proposed, and the data of the same hour in all days was used to form the model input. Then, data normalization, data set division ratio, threshold and other parameter tuning experiments were carried out to optimize the model. Finally, the NGBoost algorithm was used to obtain the target site’s PM for the next day using the discretized multidimensional data set. 2.5 The present invention provides an accurate hourly prediction value of PM concentration. 2.5 Long-term concentration series prediction provides an open research framework.

Description

Technical Field

[0001] The present invention relates to air pollutant concentration prediction, machine learning-based modeling and other technologies. Specifically, it relates to a method for predicting PM2.5 concentration in a long time series. 2.5 Concentration method. Background Art

[0002] In recent years, with the development of urbanization and industrialization, air pollution has gradually become a serious social problem. People are exposed to high concentrations of PM for a long time. 2.5 The risk of cancer and respiratory diseases is increased. 2.5 The smog formed by staying in the air for too long will reduce visibility, causing traffic congestion, flight delays and other problems, greatly affecting people's normal work and life. 2.5 Concentration prediction technology, as one of the important means of air quality monitoring, can provide early warning so that people can take protective measures in advance. It is an effective way to improve urban traffic and improve the health index of residents.

[0003] Faced with normalized smog pollution and rising PM 2.5 concentration values, we urgently need a method that can predict PM 2.5 Although there are many methods to calculate the concentration of PM 2.5 concentration prediction models, but current models mainly focus on PM 2.5 Single-step prediction of concentration, i.e. prediction of PM 1 hour later 2.5 concentration, is only applicable to short-term decision-making problems. 2.5 Multi-step concentration prediction mainly relies on deep neural network technology, which has strong learning ability and also requires high-performance computing and storage support. Therefore, in order to facilitate people's travel and protect the public's health, we urgently need to implement a long-term PM based on the NGBoost algorithm. 2.5 The concentration prediction method provides key support for management decision-making. Summary of the invention

[0004] The problem solved by the present invention is: a PM based on NGBoost algorithm is proposed 2.5 The long-term concentration series prediction method makes up for the current inability to efficiently and accurately predict the PM2.5 concentration in the next day. 2.5 The present invention comprehensively analyzes the effects of other air pollutants and meteorological factors on PM 2.5 The influence of concentration increases the accuracy of prediction. In addition, the present invention adopts NGBoost algorithm for model training, which reduces the model complexity and training time, and can obtain PM faster. 2.5Concentration prediction value. Most importantly, the present invention proposes a method of data time series discretization to process the training set, avoiding the use of complex neural networks for prediction while improving model performance and accuracy.

[0005] The technical solution of the present invention is: the present invention proposes a PM based on NGBoost algorithm 2.5 The method uses the IsolationForest algorithm to filter out the outliers in the sample and the Spearman correlation coefficient to filter out the outliers related to PM. 2.5 The characteristic information is strongly related to the concentration. On this basis, the time series features of the data are discretized to obtain a two-dimensional data list of the same hour in different days, and the three-dimensional data matrix is ​​communicated with the day data as the model input. Then the data normalization, data set division ratio, threshold and other parameters are adjusted to assist in optimizing the model. Finally, the NGBoost algorithm is used to train the multidimensional feature data to obtain the hourly PM of the target site in the next day. 2.5 Predicted concentration. BRIEF DESCRIPTION OF THE DRAWINGS

[0006] Figure 1 It is the overall step flow chart of the present invention DETAILED DESCRIPTION

[0007] (1) Data preparation:

[0008] (a) Obtain hourly monitoring data from a target air monitoring station in a certain area for N consecutive days to form the PM 2.5 Run data for multi-step prediction of concentration.

[0009] The air pollutant data mainly include the following fields: PM 2.5 Concentration, PM 10 Concentration, NO 2 concentration, CO concentration, O 3 concentration, SO 2 concentration.

[0010] Meteorological data mainly includes the following fields: temperature, air pressure, relative humidity, wind direction, and wind speed.

[0011] The above two parts of data together constitute the total data set Z. Z is a two-dimensional matrix containing 12 columns of data, that is, the 11 feature fields of the air pollutants and meteorological tables arranged in the above order plus the corresponding hourly time fields together constitute the feature set FS = {fs 1 ,fs 2 ,fs 3 …fs 12}, where the value corresponding to the jth feature in the i-th row is z ij, let α be the number of rows in the total data set minus 1, n be the sum of the number of all samples, then 0≤i≤α and 0≤j≤11.

[0012] (2) Data preprocessing, the specific steps include:

[0013] (a) Processing of outliers. This paper mainly uses the Isolation Forest algorithm to remove outliers from the total data set. First, based on real life and common sense of monitoring data, data with air pollutant concentrations less than 0 are converted to null values. Then, the Isolation Forest algorithm is used to randomly divide the samples into several trees. Each tree is then trained by randomly sampling data to ensure that the variance of each tree meets the operator's requirements. Finally, the outlier score matrix SC of all sample points in Z is calculated using the following formula, where SC = {sc 00 ,sc 01 ,sc 02 …sc xy}(x=12,y=α).

[0014]

[0015] Where E[h(z ij )] is z ij The expected value of the path length in multiple trees. c(n) is the average path length of a tree containing n sample data sets, which is used to standardize the record z ij The path length. H(*) is the harmonic number, and ξ is the Euler constant, which is approximately equal to 0.57722.

[0016] Then the outliers in the data set are returned. The threshold of the outlier score is set to ρ. When the outlier score of any data point is greater than the threshold, that is, When the outlier point z in the data set is ij Set to null value. In this way, a new total data set Z' is obtained.

[0017] (b) Eigenvalue selection. This paper uses the Spearman correlation coefficient to check PM 2.5 The correlation between concentration and other features and the correlation between other features. First, use the following formula to calculate any two features f i 、f j The correlation between ij :

[0018]

[0019] Where n is the number of samples and α is the number of the last row of the index starting from 0, and d k Represents the difference in the order of data in two feature data sets, Z' i and Z'j are the numerical columns corresponding to the i-th feature and the j-th feature in the dataset Z'.

[0020] Then check PM 2.5 The minimum threshold of the correlation coefficient between concentration and other features is η, and the maximum threshold of the correlation coefficient between other features and other features is θ. or , then delete the value column Z corresponding to the jth feature j , not used as input data for subsequent models, and other data columns remain unchanged. According to the above operations, the input data set Z of the model is finally obtained. input .

[0021] (3) Discretization of data time series. The specific steps include:

[0022] (a) Separation of training set and test set data. From the total dataset Z after preprocessing input All data from the 1st day to the dth day of the target monitoring station are selected as the training set, and all data from the d+1th day to the Nth day are selected as the test set.

[0023] (b) Time series discretization of training set data. First, arrange the training set data in chronological order. Group them by time, extract the feature data of the same time on different days in the training set and use them to form 24 two-dimensional data matrices. Let the data m of the i-th hour be i ∈M T×F , where T represents the number of days in the training set and T≤d, and L represents the number of features remaining after preprocessing. Thus, the three-dimensional data matrix M H×T×F , where H = 24 represents 24 hours in a day. Then take out PM 2.5 The data corresponding to the concentration is set as Y H×T As the label of the model, let the label data y of the i-th day and j-th hour be ij ∈Y. Then the remaining data matrix M' H×T×(F-1) Set as the three-dimensional feature matrix to be molded.

[0024] (4) Parameter tuning. The specific steps include:

[0025] (a) Data normalization. In order to eliminate the mutual influence of different feature variables, the data is normalized. Assume that the mapping interval is [a, 1], m' ijk ∈M H×T×F , use the following formula to perform maximum-minimum normalization of data:

[0026]

[0027] The experiments were conducted by setting a to -1 and 0, and the mapping interval that enables the best model performance was obtained by comparing the experimental results.

[0028] (b) Dataset partition ratio selection. The ratio of the training set to the test set was set to (6:4), (7:3), (8:2), and (9:1) for comparative experiments, and the optimal partition ratio was selected based on the experimental results.

[0029] (c) Threshold selection. Different thresholds are selected to reduce the inaccuracy of model prediction. The threshold range is set to [0.75, 0.9], and the search is performed with a step size of 0.01. The threshold that can make the model perform best is found through experimental results.

[0030] (d) In the comparative experiments of (a)-(c) above, each time the model parameters were adjusted, ten random samplings were performed in the data set. The model was trained under the same model parameters using the ten sets of sampled data. The prediction performance values ​​obtained from the ten trainings were averaged as the final prediction performance of the model parameters. The present invention uses this method to more accurately extract the optimal parameters.

[0031] (5) Use the NGBoost algorithm to train the model. The specific steps include:

[0032] (a) Select an appropriate scoring rule. The scoring rule takes the predicted probability distribution P and an observed value y as input and needs to satisfy the following formula:

[0033]

[0034] Where S is the scoring result, Q is the true distribution of the result y, and P is another distribution different from Q. Let parameter θ be the input M' T×(L-1) , then P θ represents the probability distribution of the parameter. In this paper, we use the maximum likelihood estimation method L as the scoring rule, and the scoring calculation formula is as follows:

[0035] L(θ,y)=-log P θ (y)#(5)

[0036] (b) Define the generalized natural gradient. First, according to the maximum likelihood estimation method we choose, we need to calculate the Kullback-Leibler divergence. According to the calculation method of the score L in the maximum likelihood estimation method in formula (5), we can get the calculation formula of the Kullback-Leibler divergence as follows:

[0037]

[0038] Then the Fisher information Γ is calculated L , the formula is as follows:

[0039]

[0040] in It represents the gradient of the score L in the parameter θ and the result value y. The gradient definition method is:

[0041]

[0042] Then, the metric in formula (8) is replaced by the divergence calculated by formula (6), and the generalized natural gradient can be obtained Definition method:

[0043]

[0044] Finally, according to the definition of generalized natural gradient in formula (9), the corresponding generalized natural gradient under the maximum likelihood estimation scoring rule is defined, and the results are as follows:

[0045]

[0046] (c) Use the NGBoost algorithm to predict the PM2.5 concentration in the next day. First, the three-dimensional feature matrix M' obtained in step 3 is H×T×(F-1) and label data Y H×T Loop through and select only the features and label data of the i-th hour in all days each time, that is and Then calculate the initialization parameter θ 0 , the calculation formula is as follows:

[0047]

[0048] in Then, U base learners f are constructed, and the current cycle is assumed to be the data of the kth day in the uth learner. The learning results of each learner are based on the parameter gradient of the learner set in the previous stage. This paper uses a shallow decision tree as the base learner, and the parameter gradient The calculation method of is referred to formula (10). Then the parameter gradient is fitted and the scaling factor φ of the learner is obtained using the following formula u :

[0049]

[0050] in Then, traverse again according to the number of days T, set the current cycle to the data of the kth day, update the parameters according to the following formula, and enter the next cycle:

[0051] θ uk =θ (u-1)k-r(φ u ·f u (x k ))#(13)

[0052] Where r is the learning rate, and u ≥ 0 is required to calculate formula (13). After all days and base learners are traversed, the PM of the ith hour in the T+1th day can be obtained. 2.5 After waiting for all the 24-hour data to be traversed, the PM value for the entire day of day T+1 can be obtained. 2.5 Predicted concentration.

[0053] (6) Output the PM of the site in the next day 2.5 Predicted concentration.

Claims

1. A PM based on NGBoost algorithm 2.5 The method for predicting long concentration series is characterized by: The implementation steps of this method are as follows: (1) Obtain air pollutant and meteorological data from target air monitoring stations in a certain area, and store and pre-process the data; (2) discretizing the data obtained in step (1) in time series, extracting the feature data of the same hour in different days to form a feature matrix; (3) Normalize the data, adjust the model data set division ratio and threshold parameters, conduct multiple groups of comparative experiments, and complete parameter tuning; (4) Use the NGBoost algorithm to train and verify the discretized two-dimensional feature matrix; (5) Output the hourly PM of the target site for the next day 2.5 Concentration prediction results; In the step (1), data preprocessing includes: (1) The Isolation Forest algorithm is used to remove outliers from the total data set. First, the data with air pollutant concentration less than 0 is converted into null values. The Isolation Forest algorithm is used to randomly divide the samples into several trees. Each tree is trained by randomly sampling data to ensure that the variance of each tree meets the operator's requirements. Finally, the outlier score matrix SC of all sample points in Z is calculated using the following formula, where SC = {sc 00 ,sc 01 ,sc 02 …sc xy }, x = 12, y = α; Where E[h(z ij )] is z ij The expected value of the path length in multiple trees; c(n) is the average path length of a tree containing n sample data sets, which is used to standardize the record z ij ; H(*) is the harmonic number, ξ is the Euler constant; Then the outliers in the data set are returned. The threshold of the outlier score is set to ρ. When the outlier score of any data point is greater than the threshold, that is, When the outlier point z in the data set is ij Set to null value; get the new total data set Z'; (2) Calculation of PM 10 concentration, NO2 concentration, CO concentration, O3 concentration, SO2 concentration, temperature, air pressure, relative humidity, wind direction, wind speed characteristic sequence and PM 2.5 Spearman correlation coefficient of concentration series, selected with PM 2.5 Factors with strong correlation are used as features of the input model; first, the following formula is used to calculate any two features f i 、f j The correlation between ij : Where n is the number of samples and α is the number of the last row of the index starting from 0, and d k Represents the difference in the order of data in two feature data sets, Z' i and Z' j are the numerical columns corresponding to the i-th feature and the j-th feature in the data set Z' respectively; Then check PM 2.5 The minimum threshold of the correlation coefficient between concentration and other features is η, and the maximum threshold of the correlation coefficient between other features and other features is θ; or , then delete the value column Z corresponding to the jth feature j , not used as input data for subsequent models, other data columns remain unchanged; according to the above operations, the input data set Z of the model is finally obtained input ; (3) constructing a three-dimensional data matrix of the target site based on the characteristic sequences obtained by screening; In the step (2), data discretization includes: First, separate the training set and test set data; from the preprocessed input data set Z input All data from the 1st day to the dth day of the target monitoring station are selected as the training set, and all data from the d+1th day to the Nth day are selected as the test set. Then, the training set data is discretized in time series; first, the training set data is arranged in chronological order; grouped by time, the feature data of the same time on different days in the training set are taken out and used to form 24 two-dimensional data matrices; let the data m of the i-th hour be i ∈M T×F , where T represents the number of days in the training set and T≤d, and L represents the number of features remaining after preprocessing; thus, the three-dimensional data matrix M H×T×F , where H = 24 represents 24 hours in a day; then take out PM 2.5 The data corresponding to the concentration is set as Y H×T As the label of the model, let the label data y of the i-th day and j-th hour be ij ∈Y H×T ; Then the remaining data matrix M' H×T×(F-1) Set as the three-dimensional feature matrix to be molded; In the step (3), parameter tuning includes: (1) Data normalization: In order to eliminate the mutual influence of different feature variables, the data is normalized. Let the mapping interval be [a, 1], m' ijk ∈M H×T×F , use the following formula to perform maximum-minimum normalization of data: Set a to -1 and 0 respectively to conduct experiments, and compare the experimental results to obtain the mapping interval that can make the model perform best; (2) Data set partition ratio selection: set the ratio of training set to test set to (6:4), (7:3), (8:2), (9:1) for comparative experiments, and select the optimal partition ratio based on the experimental results; (3) Threshold selection: Different thresholds are selected to reduce the inaccuracy of model prediction. The threshold range is set to [0.75, 0.9], and the search is performed with a step size of 0.

01. The threshold that can make the model perform best is found through experimental results. (4) In the comparative experiments (1) to (3) above, each time the model parameters were adjusted, ten random samples were taken from the data set; the model was trained using the ten sets of sampled data under the same model parameter conditions; the prediction performance values ​​obtained from the ten training runs were averaged as the final prediction performance of the model parameters.

Images