Chemical production anomaly detection method based on space-time diagram variational encoder

Abstract

The invention relates to a chemical production anomaly detection method based on a space-time diagram variational encoder. The method comprises the following steps: selecting public chemical normal working condition data D; an anomaly monitoring model ST-GVAE is constructed, a loss function L (.) is constructed, D serves as input, an Adam optimizer is adopted to train ST-GVAE, model parameters are reversely updated through gradient descent, and trained ST-GVAE'is obtained; a to-be-detected target variable a is selected, normal working condition historical data corresponding to a are input into ST-GVAE ', and an anomaly detection score corresponding to a is obtained and serves as an anomaly detection threshold value of the target variable a; inputting the current chemical working condition data of a into ST-GVAE'to obtain an anomaly detection score corresponding to the current state of a, and if the score is greater than a threshold value, triggering an alarm; and then calculating a reconstruction deviation score B of each node in the PID, taking the node corresponding to the B greater than a deviation threshold as an abnormal node, and outputting the abnormal node. By using the method provided by the invention, the variable abnormal condition of the existing chemical process can be accurately detected and positioned.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of industrial monitoring and fault diagnosis, and particularly relates to a chemical production anomaly detection method based on a space-time graph variational encoder. BACKGROUND

[0002] Chemical production is an important pillar industry for national economy and social development, and the safe operation of its process and the guarantee of product quality are the core propositions of enterprises. The chemical process contains complex process flow and control system, and each variable in the production process is related to each other, and is marked by significant continuity and complex reaction mechanism. The safety detection of traditional chemical process mainly relies on distributed control system (DCS) and advanced process control (APC) technology, but the operation under abnormal conditions is still mainly manual operation. Violation of operation by operators and operation errors can easily lead the chemical process to an abnormal state. At the same time, the data recorded by DCS often implies key information of the production status of the chemical process, which lays a solid data foundation for the development of intelligent monitoring system based on data driving.

[0003] In recent years, the modern intelligent monitoring system of chemical process based on data driving has become the mainstream, and the anomaly detection of chemical process as the first line of defense of the process intelligent monitoring system, its essence is to realize the real-time perception of the health status of the production system through multi-modal feature fusion. Its core function architecture contains double goals, one is to establish a process detection threshold to detect whether the production operation state has failed. When an anomaly is detected, an alarm is sent to the engineers and on-site operators in time. The second is to realize fault identification, that is, to locate the fault to the relevant variables, so as to provide a feature basis for the fault diagnosis process, so that the operators can quickly and accurately take measures to eliminate the impact of the fault on production, and keep the production system safe, stable and reliable. With the increasing complexity of industrial systems, anomaly detection of multi-variable time series data faces many challenges: traditional methods are mostly based on single variable or simple combination between multi-variables, which is difficult to effectively model the space-time correlation between variables; sensor anomalies may trigger a chain of abnormalities through the control loop, leading to detection and positioning difficulties; data-driven models lack the fusion of prior knowledge (such as mechanism model, expert experience), which is easy to be out of touch with the actual process; the scarcity of fault labels in industrial scenarios makes unsupervised learning more practical.

[0004] For the complex and changeable chemical production environment, the existing abnormal detection model has unsatisfactory accuracy in abnormal detection, and the performance of single network structure has certain limitations, the basic series network structure lacks sufficient consideration of system complexity, and it is difficult to effectively integrate the multidirectional information in process variables, so that the abnormal detection model has certain gap in actual application. SUMMARY

[0005] In view of the above problems existing in the prior art, the technical problem solved by the present application is how to accurately detect abnormal faults in chemical production process.

[0006] To solve the above technical problems, the present application adopts the following technical solutions:

[0007] A chemical production abnormality detection method based on a space-time graph variational encoder, comprising the following steps:

[0008] S100: selecting public chemical normal working condition data D, D including historical data of all variables under normal working condition and PID graph of all variables;

[0009] S200: constructing an abnormality monitoring model ST-GVAE, the ST-GVAE comprising a data preprocessing module and an abnormality detection module, the data preprocessing module being used for preprocessing the chemical normal working condition data, the preprocessing including topology graph conversion, data feature extraction and data format standardization conversion;

[0010] The abnormality detection module is composed of an LSTM network embedded with GCN (LSTM-GCN), the LSTM-GCN comprising an encoder, an MLP and a decoder, the encoder and the decoder having the same network structure; the abnormality detection module being used for calculating the abnormality detection score of the target variable in the chemical process and the sensor positioning result corresponding to the target variable;

[0011] S300: constructing a loss function L(·), training the ST-GVAE using D, adopting an Adam optimizer, updating the ST-GVAE parameters by gradient descent in reverse, stopping the training when the loss function reaches the minimum or reaches the maximum training iteration number, obtaining the trained ST-GVAE', and the calculation expression of the loss function is as follows:

[0012]

[0013] Wherein, θ represents the decoder network parameter, represents the encoder network parameter, x represents the node, L represents the sampling number of estimated expected value, l represents the lth sampling, D KL represents the KL divergence term, D KL By constraining the hidden variable approximating the prior distribution pθ ;

[0014] S400: input the historical data of all variables under normal working conditions into ST-GVAE', and obtain the abnormality detection scores corresponding to all variables;

[0015] S500: select the target variable a to be detected, input the historical data of a under normal working conditions into ST-GVAE', and obtain the abnormality detection score corresponding to a, which is taken as the abnormality detection threshold of a;

[0016] Then input the current chemical working condition data corresponding to a into ST-GVAE' to obtain the abnormality detection score corresponding to a. When the abnormality detection score is greater than the abnormality detection threshold of a, an alarm is triggered and the next step is performed, otherwise, it indicates that the variable detection is normal and the detection is completed;

[0017] S600: take the sensor position in the PID graph corresponding to the target variable a to be detected as a node, define a deviation threshold, and calculate the reconstruction deviation score B of each node. The node corresponding to B greater than the deviation threshold is taken as an abnormal node, and then the sensor position corresponding to the abnormal node is output.

[0018] As preferred, the step of preprocessing the chemical normal working condition data in S200 is as follows:

[0019] S210: define the chemical priori knowledge as G=(V, E, A, X), wherein V represents a node set composed of sensors, E represents an edge set, A represents an N*N dimensional adjacency matrix, and X represents a node feature matrix of the graph structure;

[0020] S211: convert the PID graph into a topological graph by G using the chemical mechanism, and obtain a topological structure graph;

[0021] S212: select the target variable data for modeling from the historical data of all variables under normal working conditions as variable data features, which are taken as data feature extraction values;

[0022] S213: adopt a filtering method to denoise the topological structure graph and the data feature extraction values, and then perform format standardization conversion on the denoised topological structure graph and the data feature extraction values, and finally obtain the standardized format data required by ST-GVAE.

[0023] As preferred, in S300, the ST-GVAE is trained by using an Adam optimizer, and the steps are as follows:

[0024] S310: process D by using a data preprocessing module to obtain corresponding standardized format data, divide the standardized format data into time windows, construct a spatio-temporal graph data set, and the time window is defined as follows:

[0025] X t-w+1:t = [X t-w+1 ,..., X t-1 , X t ]

[0026] wherein X t represents the time window data at time t, X t-w+1:t represents a time window interval, and w represents the length of the time window;

[0027] S320: send the connection relationship A between each time window and node in the spatio-temporal graph dataset to the encoder, generate the hidden state sequence [h t-w+1 ,..., h t-1 , h t ] step by step through LSTM-GCN, and take the hidden state h t at the final time step as the global spatio-temporal feature representation; wherein the hidden state at each time step in the hidden state sequence is updated as follows:

[0028] h T = LSTM-GCN(X T , A, h T-1 ), T e [t-w+1, t]

[0029] wherein T represents the time step of the next hidden state;

[0030] S330: input h t to the MLP to calculate the hidden variable set Z t , and the calculation expression is as follows:

[0031] μ i = W μ (i) ·h t (i) +b μ (i)

[0032] σ i = W σ (i) ·h t (i) +h σ (i)

[0033]

[0034] Z t (i) = μ i + σ i , σ i ~ N(0, 1)

[0035] Z t = {Z t (1) ,…,Z t (I) ,…Z t (N)}

[0036] where μ i represents the hidden variable generation mean parameter of node i, σ i represents the hidden variable standard deviation parameter of i, W μ (I) , W σ (I) represents the weight parameter, b μ (i) , b σ (i) represents the bias term, q(Z t |X t-w+1：t ) represents the global hidden variable distribution, represents the hidden variable of node i, and N(0, 1) represents the normal distribution.

[0037] S340: input Z t into the decoder to obtain the anomaly detection score of all variables at time t;

[0038] Calculate L(·) and determine whether the training stop condition is reached. If yes, determine the current ST-GVAE as the trained ST-GVAE', otherwise return to S320. The training stop condition is that the loss function reaches the minimum or reaches the maximum training iteration number.

[0039] As a preferred, the calculation content of calculating the reconstruction deviation score B of each node in S600 is as follows:

[0040] Reconstruct the reverse time sequence in the decoder from time step t to t-w+1 in reverse, and then input Z t as input, and map Z t back to the original data space according to the reverse time sequence to obtain the reconstructed data.

[0041] Calculate the square of the absolute error of the normal working condition historical data corresponding to a and the reconstructed data corresponding to a, and take the square of the absolute error as the reconstruction deviation score B of the node in the PID diagram corresponding to a. The calculation formula is as follows:

[0042] B = (x' t -x t ) 2

[0043] where x' t represents the reconstructed data, and xt represent original working condition data.

[0044] Compared with the prior art, the present application has at least the following advantages:

[0045] The method of the present application integrates two powerful technologies, STGNN and VAE, into a unified framework; STGNN is responsible for effectively extracting and modeling the spatio-temporal dependencies between process variables, and VAE is responsible for learning the complex normal patterns formed by these spatio-temporal features and detecting abnormalities accordingly; this integration aims to overcome the limitations of single models in handling the spatio-temporal characteristics of complex chemical process data.

[0046] The present application proposes an unsupervised chemical process anomaly detection model ST-GVAE based on the topology of process mechanism. The prior knowledge is integrated into the process data combined with the topological information between sensors, and the graph convolution network, LSTM and VAE technology are integrated to efficiently capture the spatio-temporal dependencies between sensors; the model uses normal working condition data and constructs an anomaly score for fault data based on reconstruction error, and uses normal fault-free data to construct a detection threshold for the anomaly score of fault data; the performance of the proposed framework is verified using various fault conditions in the Tennessee Eastman process, and compared with LSTM-AE, LSTM-VAE and LSTM-GCN-AE, effectively proving the superior performance of the ST-GVAE model in chemical process anomaly detection and positioning. BRIEF DESCRIPTION OF DRAWINGS

[0047] Figure 1 It is an abnormal detection flowchart based on ST-GVAE of the present application.

[0048] Figure 2 It is a technical route diagram for constructing the unsupervised anomaly detection model ST-GVAE of the present application.

[0049] Figure 3 It is a diagram for constructing a graph structure based on prior knowledge and data of a chemical process of the present application.

[0050] Figure 4 It is the internal structure of LSTM-GCN in the present application.

[0051] Figure 5 It is the encoding-decoding architecture of ST-GVAE in the present application.

[0052] Figure 6 It is the PID graph of the TE industrial process in the present application.

[0053] Figure 7For comparison before and after Savitzky-Golay filter denoising in the experiment; Figure (a) shows the feed flow rate of component A, Figure (b) shows the pressure of the reactor, Figure (c) shows the liquid level of the reactor, and Figure (d) shows the temperature of the reactor.

[0054] Figure 8 For the relationship topology diagram of the TE process sensor in the experiment.

[0055] Figure 9 For the detection results of 20 fault types by the ST-GVAE model.

[0056] Figure 10 For IDV (1): the bias score calculated according to the prediction error of the ST-GVAE model.

[0057] Figure 11 For IDV (8): the bias score calculated according to the prediction error of the ST-GVAE model.

[0058] Figure 12 For IDV (13): the bias score calculated according to the prediction error of the ST-GVAE model.

[0059] Figure 13 For IDV (14): the bias score calculated according to the prediction error of the ST-GVAE model. DETAILED DESCRIPTION

[0060] The present application will be further described in detail below.

[0061] The present application patent reasonably utilizes advanced technologies such as graph convolution network, long short-term memory network, and variational autoencoder, fuses chemical process mechanism, constructs an unsupervised abnormality detection framework that fuses chemical process information, develops multivariate spatiotemporal correlation modeling technology, captures the dynamic coupling characteristics of sensor data, designs an accurate abnormality positioning mechanism, and locates the process abnormality to the relevant variables, thereby improving the accuracy and interpretability of abnormality detection in complex industrial scenarios.

[0062] The offline training part is to train the model by using existing normal working condition data, so as to obtain an optimal model; the online detection is to check the abnormal state of the current variable to be detected; see Figures 1-2 A chemical production abnormality detection method based on a spatiotemporal graph variational encoder, comprising the following steps:

[0063] S100: selecting public chemical normal working condition data D, wherein D includes historical data of all variables under normal working conditions and PID graphs of all variables;

[0064] S200: Construct an anomaly monitoring model ST-GVAE, the ST-GVAE includes a data preprocessing module and an anomaly detection module, the data preprocessing module is used for preprocessing the chemical normal working condition data, the preprocessing includes topology graph transformation, data feature extraction and data format standardization conversion;

[0065] The anomaly detection module is composed of an LSTM network embedded with a GCN (LSTM-GCN), the LSTM-GCN includes an encoder, an MLP and a decoder, the encoder and the decoder have the same network structure; the anomaly detection module is used for calculating and outputting the anomaly detection score of the target variable in the chemical process and the sensor positioning result corresponding to the target variable; the MLP, the GCN network and the LSTM network are all prior art; the graph convolution network is used to learn the correlation between sensor variables and embed the LSTM network to learn the time correlation of process variables, so as to realize the spatio-temporal fusion of process data;

[0066] The topology graph constructed based on the PID graph and the conservation law defines the physical interaction direction between variables and provides prior knowledge constraints for the model. From the perspective of information transmission, different sensor nodes describe the process information carried by the current equipment, and the directed edge describes the direction of information transmission or message transmission. In addition, the graph structure based on the knowledge of the chemical process qualitatively indicates the relationship between different sensors, and also quantitatively learns the correlation strength between variables from time series data through the graph convolution network. Therefore, this process converts multivariate time series into graph data containing physical semantics, so that the model can not only utilize the explicit rules of chemical mechanism, but also mine the dynamic correlation implied in the data, forming a deep extraction of the correlation between process variables;

[0067] The step of preprocessing the chemical normal working condition data in S200 is as follows:

[0068] S210: Define the chemical prior knowledge as G=(V,E,A,X), wherein V represents a node set composed of sensors, E represents an edge set, A represents an N*N dimensional adjacency matrix, and X represents a node feature matrix of the graph structure; let v i ∈V represent a node, e ij =(v i ,v j )∈E represent an edge from v i to v j . If e ij ∈E, then A ij =1, and if , then A ij =0;

[0069] S211: Topology graph transformation of the PID graph through G is performed by using the chemical mechanism, and a topology structure graph is obtained; the chemical mechanism is prior art;

[0070] S212: Select target variable data for modeling from the historical data of all variables under normal working conditions as variable data features, which are data feature extraction values;

[0071] S213: Adopt a filtering method to denoise the topological structure graph and the data feature extraction values, and then perform format standardization conversion on the denoised topological structure graph and the data feature extraction values, to finally obtain standardized format data required by the ST-GVAE, the purpose of data format standardization is to avoid the influence of dimension difference on the performance of the model, and the method used for data format standardization is z-score; wherein, the filtering method and the z-score are prior art;

[0072] The topological graph constructed based on the PID graph and the conservation law defines the physical interaction direction between variables, provides prior knowledge constraints for the model, and from the perspective of information transmission, different sensor nodes describe the process information carried by the current device, and the directed edge describes the direction of information transmission or message transmission; in addition, the graph structure based on the chemical process knowledge qualitatively indicates the relationship between different sensors, and quantitatively learns the correlation strength between variables from time series data through the graph convolution network; therefore, the process converts the multivariate time series into graph data containing physical semantics, so that the model can not only utilize the explicit rules of chemical mechanism, but also mine the dynamic correlation implied in the data, forming a deep extraction of the correlation between process variables;

[0073] S300: Construct a loss function L(·), train the ST-GVAE using D, use the Adam optimizer to train the ST-GVAE, update the ST-GVAE parameters by gradient descent in reverse, stop training when the loss function reaches the minimum or reaches the maximum training iteration number, obtain the trained ST-GVAE', and the calculation expression of the loss function is as follows:

[0074]

[0075] Wherein, θ represents the decoder network parameter, φ represents the encoder network parameter, x represents the node, L represents the number of estimated expected values, l represents the lth sample, D KL represents the KL divergence term, D KL By constraining the latent variable to approximate the prior distribution p θ ; the ST-GVAE model realizes the optimization of the model parameters by gradually reducing the loss function, in the present application, L=1, as long as the small batch is large enough, a single latent variable sample can effectively estimate the expected value, balancing the calculation efficiency and convergence stability;

[0076] In the S300, the Adam optimizer is used to train the ST-GVAE, and the steps are as follows:

[0077] S310: The data preprocessing module is used to process D to obtain corresponding standardized format data, the standardized format data is divided into time windows, and a space-time graph data set is constructed, and the time window is defined as follows:

[0078] X t-w+1:t = [X t-w+1 ,..., X t-1 , X t ]

[0079] Wherein, X t represents the time window data at t time, X t-w+1:t represents a time window interval, and w represents the length of the time window; the data of a time window under normal operating conditions is taken as a sample, and taking the TE process as an example, there is a point every 3 minutes, and the size of the time window is 20, that is, the size of each sample is 1h;

[0080] S320: Send the connection relationship A between each time window and node in the space-time graph data set to the encoder, and generate the hidden state sequence [h t-w+1 ,..., h t-1 , h t ] through LSTM-GCN step by step, and take the hidden state h t at the final time step as the global space-time feature representation; wherein the hidden state of each time step in the hidden state sequence is updated as follows:

[0081] h T = LSTM-GCN (X T , A, h T-1 ), T ∈ [t-w+1, t]

[0082] Wherein, T represents the time step of the next hidden state; the hidden state h t is further transmitted to a group of multilayer perceptron, the sensors in the industrial system often have heterogeneity, and the fault modes of different devices are significantly different. The global shared multilayer perceptron (MLP) will force all nodes to comply with the same parameter generation rule, resulting in the decline of the response ability of the model to local weak anomalies. Independent parameterization makes the hidden variable distribution of each node have clear physical correspondence, which is convenient for engineers to locate the root cause of the anomaly;

[0083] S330: Input h t into the MLP to obtain a set of hidden variables Z t , and the calculation expression is as follows:

[0084] μ i = W μ (i) ·h t(i) +b μ (i)

[0085] σ i =W σ (i) ·h t (i) +b σ (i)

[0086]

[0087] Z t (i) =μ i +σ i ,σ i ~N(0,1)

[0088] Z t ={Z t (1) ,…,Z t (i) ,…Z t (N)}

[0089] where μ i represents the hidden variable generation mean parameter of node i, σ i represents the standard deviation parameter of the hidden variable of i, W μ (i) , W σ (i) represents a weight parameter, b μ (i) , b σ (i) represents a bias term, q(Z t |X t-w+1：t ) represents a global hidden variable distribution, represents the hidden variable of node i, and N(0,1) represents a normal distribution.

[0090] S340: input Z t into the decoder to obtain an anomaly detection score of all variables at time t;

[0091] Calculate L(·) and determine whether a training stop condition is reached. If the training stop condition is reached, the current ST-GVAE is determined as the trained ST-GVAE', otherwise, return to S320. The training stop condition is that the loss function reaches a minimum or a maximum training iteration number is reached.

[0092] Select the working condition data to be detected for abnormal variables, that is, perform feature extraction on related variables; the target of the online detection stage is to monitor process data in real time, trigger abnormal early warning, and locate abnormal variables.

[0093] When the ST-GVAE model trained is used for anomaly detection on the online test set data, the principle is to train the ST-GVAE model on the normal working condition data without faults, and detect the anomaly score for anomaly detection. The threshold is calculated by the anomaly score of the normal working condition. By calculating the difference between the reconstructed data and the test data at each time, the anomaly score can reveal the abnormality of the observation. Based on the reconstruction probability distribution of the input data X t at time t, the anomaly score tempS t can be established. As shown in the following formula, the negative logarithm of the reconstruction probability at time t maps the normal working condition data with high probability to low anomaly score, and the abnormal working condition data with low probability to high anomaly score:

[0094]

[0095] Abnormal observations usually occur continuously, and in practice, it is acceptable to trigger an alarm within a short delay; the anomaly score of an observation X t is adjusted according to some subsequent observations; as a stochastic sequence-to-sequence model, the ST-GVAE model not only learns to reconstruct X t , but also learns to reconstruct the w-1 observations before X t , which facilitates updating the anomaly score of the observation X t , and the final anomaly score S t of the observation X t is calculated as follows:

[0096]

[0097] where L denotes the number of samples, and L = 1 is taken, following the single-sample estimation strategy during training, and D is the detection delay tolerance, allowing the score of the current time step to be updated in combination with the subsequent D-1 time steps; specifically, at time t, the instantaneous score tempS t is calculated. In the subsequent D-1 time windows, tempS t+d is recalculated, and the average of the instantaneous scores S t calculated for all D times is obtained. The anomaly score can reveal the abnormality of the observation, but in practice, a threshold is still needed to trigger an alarm and indicate the operator to take action. This study proposes a threshold selection method based on the anomaly score sequence {S t} extracted from the normal working condition training set, and calculates the mean and standard deviation thereof. By the 3σ criterion, the threshold Γ is set to cover 99.7% of the normal data. As shown in the following formula:

[0098] μ nor = E{S t}

[0099]

[0100] Gamma = mu nor + 3sigma nor

[0101] wherein mu nor represents a normal mean parameter, and sigma nor represents a normal standard deviation;

[0102] S400: input the historical data of all variables under normal working conditions into the ST-GVAE', to obtain the abnormal detection scores corresponding to all variables;

[0103] S500: online detection: select a target variable a to be detected, input the historical data of a corresponding to the normal working conditions into the ST-GVAE', to obtain the abnormal detection score corresponding to a, and take the abnormal detection score as the abnormal detection threshold of a;

[0104] Then input the current data of a corresponding to the chemical working conditions into the ST-GVAE' to obtain the abnormal detection score corresponding to a, when the abnormal detection score is greater than the abnormal detection threshold of a, an alarm is triggered and the next step is performed, otherwise, it indicates that the variable detection is normal, and the detection is ended;

[0105] S600: take the sensor position in the PID diagram corresponding to the target variable a to be detected as a node, define a deviation threshold, and calculate the reconstruction deviation score B of each node, the deviation score is obtained by using MSE calculation, take the node corresponding to B greater than the deviation threshold as an abnormal node, and then output the sensor position corresponding to the abnormal node. The sensor is taken as a node, the position of which corresponds to the sensor position in the PID diagram one by one, the relationship between the sensors is constructed according to the PID diagram and chemical knowledge, and this step only needs to know which sensor has an influence on another, that is, the graph structure is derived according to the influence relationship, and the connection strength between the graph structures is learned according to the D data;

[0106] The model can calculate the reconstructed output as a whole for each group of input data, the overall abnormal score of the model output is smaller for normal data, and the abnormal score of abnormal data is larger; when an abnormality occurs, that is, the abnormal score is larger, the deviation score of each sensor variable between the reconstructed data and the test data is calculated to locate the abnormality caused by which variables; after locating the variables, according to the time sequence of the deviation variables and the chemical common sense, the possible starting point position of the abnormality can be known;

[0107] The calculation content of the reconstruction deviation score B of each node in S600 is as follows:

[0108] The reverse time series is reconstructed in reverse from time step t to t-w+1 in the decoder, and then Z t As input, Z t is mapped back to the original data space according to the reverse time series to obtain reconstructed data.

[0109] The absolute error square of the corresponding normal working condition historical data of a and the reconstructed data of a is calculated, and the absolute error square is taken as the reconstructed deviation score B of the node in the PID graph corresponding to a. The calculation formula is as follows:

[0110] B = (x' t -x t ) 2

[0111] Wherein, x' t represents reconstructed data, and x t represents original working condition data.

[0112] The specific process is as follows: for abnormal working conditions, when the abnormal score exceeds the threshold value, the abnormal variable is determined by calculating the deviation score of the reconstructed data of each sensor variable and the test data; after locating these variables, the possible starting point position of the abnormality can be inferred according to the time sequence of the deviation variable and the chemical process knowledge; based on the small deviation score between normal data and reconstructed data, but the reconstructed data of abnormal data and the abnormal data itself will produce a large deviation score, the fault is located to the related variable, and through the analysis of the information transmission source of the chemical process, the variable range causing the fault is further narrowed and the specific variable is locked.

[0113] Experimental content and results

[0114] The Tennessee Eastman (TE) process was created by the process control group of Eastman Chemical Company in 1992, which accurately simulates the coupling behavior of multiple process units through modular dynamic models. Because the data generated by the process is highly close to the real industrial scene in dynamic characteristics and noise distribution, the TE process is one of the most authoritative evaluation benchmarks in the field of data-driven chemical process anomaly detection and fault diagnosis.

[0115] The flow chart of the TE process is as follows Figure 6As shown, gaseous reactants A, C, D, E are fed into the reactor, where they react to form liquid products. The reaction is catalyzed by a non-volatile catalyst dissolved in the liquid phase, and the reactor has an internal cooling bundle to remove the heat of reaction. The products exit the reactor in vapor form along with unreacted feed, and the catalyst remains in the reactor. The reactor product stream passes through a cooler to produce condensed products, then enters a vapor-liquid separator. The non-condensed component B is recycled back to the reactor feed stream by a centrifugal compressor. The condensed components flow to a product stripping column, which is stripped by feed stream 4 to remove remaining reactants. Products G and H flow from the stripping column base and are separated in a downstream refining section. Inert gas and byproducts are purged from the vapor-liquid separator primarily as vapor.

[0116] The process consists of five main unit operations: reactor, product condenser, vapor-liquid separator, recycle compressor, and product stripping column. The TE process produces two products from four reactants, making a total of eight components: A, B, C, D, E, F, G, and H. Specifically, gaseous reactants A, C, D, and E produce main products G and H and byproduct F. B is an inert gas that does not participate in the reaction process. The chemical reactions are shown in the following equations:

[0117] a(g) + C(g) + D(g) → G(liq) Product 1

[0118] A(g) + C(g) + E(g) → H(liq) Product 2

[0119] A(g) + E(g) → F(liq) Byproduct

[0120] 3D(g) → 2F(liq) Byproduct

[0121] where g represents gas and liq represents liquid. All reactions are irreversible and exothermic. The reaction rate is a function of Arrhenius expression and temperature.

[0122] The TE process for anomaly detection is a dataset containing normal conditions and 20 pre-defined fault types, which is available on a public website; it includes training and testing datasets for normal conditions and process simulation runs for each fault, with different random seeds for each fault and normal operation, a sampling interval of 3 minutes, a training set containing 500 time series samples for each experiment, a testing set containing 960 test samples, and up to 500 groups of experiments.

[0123] The process contains 53 process variables, including 22 continuous measurement variables, 19 sampled measurement variables, and 12 manipulated variables. In this study, the variables related to composition measurement were removed because the low sampling frequency of composition analysis data does not match the requirement of continuous process monitoring, and in actual chemical processes, the detection of composition is mostly intermittent. In addition, during feature screening, it was found that the continuous variable "xmeas_19" and the operation variables "xmv_5, xmv_9 and xmv_12" had no significant fluctuations and could be considered as system constants. Therefore, they were deleted from the variable list, and the final variables included in the analysis are shown in Table 1.

[0124] Table 1 Description of TE process variables involved in the detection model

[0125]

[0126]

[0127] To verify the anomaly detection performance of the ST-GVAE model, a test set was constructed based on the TE process data set containing 20 pre-defined fault types. Only 960 time series samples labeled as test data in the TE process data set were extracted as a group of data for each type of fault sample, and these faults were introduced from the 161st data point, as shown in Table 2. The test data contains 15 known and 5 unknown fault types, in which faults 1-7 are related to step changes in process variables, faults 8-12 are related to random changes in specific process variables, fault 13 is a slow drift of reaction kinetics, faults 14 and 15 are caused by a valve, and faults 16-20 are unknown types. Based on these test data, the generalization ability of the model to untrained abnormal patterns was evaluated.

[0128] Table 2 TE process fault types

[0129]

[0130]

[0131] (1) Data denoising and data smoothing

[0132] Since the ST-GVAE model adopts an unsupervised training strategy, a training set was constructed based on 20 sets of normal operating condition data (each set of normal operating condition data contains 500 data, a total of 10,000 data). The test set selected 960 data for each fault type to evaluate the sensitivity of the model to 20 types of faults. To address the problem of noise interference in time series data, the Savitzky-Golay filter was used to perform data denoising. The Savitzky-Golay filtering algorithm performs least squares fitting of a polynomial within a sliding window to achieve signal smoothing and high-frequency noise suppression, while maintaining the waveform characteristics of the original data signal. This study used a moving time window containing 21 points and a fitting method based on a cubic polynomial to reduce the impact of noise. To demonstrate the effect of data smoothing and denoising, four key variables under IDV(1) were selected: Nos.(1), (7), (8), and (9) for example, as shown in FIG. 6. By comparing the data before and after denoising Figure 7 , it can be seen that the denoised data maintains the integrity of the significant features in the parameters of component A feed flow, reactor pressure, liquid level, and temperature, proving that the Savitzky-Golay filter does not lose data characteristics when dealing with noise, providing high-quality input for the subsequent anomaly detection model. Figure 7

[0133] (2) Cross-scale data standardization

[0134] To eliminate the impact of dimensional differences on model performance, z-score standardization was performed on the denoised data. The standardization parameters were calculated based on the denoised normal operating condition training data and applied to the test set simultaneously to ensure consistency in data distribution.

[0135] (3) Graph structure construction based on process mechanism

[0136] Based on the PID graph of the TE process Figure 6 and the analysis of the mutual influence between the TE process variables based on chemical prior knowledge, a directed topological graph structure Figure 8 between sensor nodes was constructed. The node numbers and color coding involved in the graph structure are consistent with Figure 6 . The direction of the arrow represents the influence relationship between one sensor node and another sensor node, and also represents the direction of information flow in the process

[0137] In this process, the mutual influence relationship between the operating variables and the corresponding variables needs to be considered. For example, for the reactor, the feed flow first affects the pressure and liquid level of the reactor, and then affects the temperature of the reactor. The temperature of the reactor is also affected by the inlet flow of the reactor cooling water. The temperature of the reactor also affects the temperature of the outlet of the reactor cooling water. ​

[0138] As Figure 9 shown, the ST-GVAE model shows significant anomaly detection capability in the 20 fault test sets (including 15 known faults and 5 unknown faults). The threshold based on normal data training (green dotted line)

[0139] The training set data μ nor = 3.7379, σ nor = 1.91, Γ = 9.4679.

[0140] The normal and abnormal states can be effectively distinguished (blue abnormal score curve), and the score rises rapidly after the abnormal injection time (red vertical line).

[0141] Experiments show that for high sensitivity faults, including fault types IDV(1), IDV(2), IDV(4), IDV(5), IDV(6), IDV(7), IDV(8), IDV(10), IDV(11), IDV(12), IDV(13), IDV(14), the average detection delay is less than 3 sampling time (less than 9 min). For unknown fault types IDV(16), IDV(17), IDV(18), IDV(19), IDV(20), the average detection delay time is 60 sampling points (less than 1 h). For low sensitivity fault types IDV(3), IDV(9), IDV(15), the detection accuracy is low. Even so, the model can detect the occurrence of these abnormalities to a certain extent under the condition of setting a lower threshold. Analysis shows that the reason for the difficulty in detecting low sensitivity fault types is that these conditions only cause slight fluctuations in the operating variables and measured variables after the fault is introduced.

[0142] Model detection performance comparison

[0143] In order to verify the advancement of the ST-GVAE model, comparative experiments are carried out with three baseline models, LSTM-AE, LSTM-VAE, and LSTM-GCN-AE. Based on the unified threshold setting strategy, the abnormal detection thresholds of each model are set as shown in Table 3.

[0144] Table 3 Abnormal detection thresholds of different methods

[0145] Method LSTM-AE LSTM-VAE LSTM-GCN-AE ST-GVAE <![CDATA[μ nor ]]> 0.0565 0.9248 1.1765 3.7397 nor ]]> ​ 0.0298 1.8018 0.7826 1.91 Threshold value Γ 0.1459 6.3302 3.5243 9.4697

[0146] Table 4 Comparison of abnormal detection performance of different models

[0147]

[0148]

[0149] According to the set anomaly detection threshold, the anomaly detection performance of the model for each type of fault can be calculated, as shown in Table 4. The evaluation criteria of the model use accuracy (accyracy) and F1 score. The experimental results are analyzed from the dimensions of global performance, typical fault response, etc.

[0150] From the perspective of the global performance of the model, remove the low sensitivity fault classes: IDV(3), IDV(9), IDV(15). From the analysis in Table 4, the average accuracy of the baseline model LSTM-AE is (80.7%), and the average F1 score is (82.8%). After introducing graph structure learning, the accuracy of LSTM-GCN-AE is 79.1%, and the F1 score is 87.0%. The results show that introducing graph structure helps to improve the overall performance of anomaly detection. The average accuracy of the ST-GVAE model is 87.1%, and the average F1 score is 89.8%. The accuracy and F1 score of the model are improved by 8% and 2.8% compared with LSTM-GCN-AE, respectively. This proves the effectiveness of the design of the variational autoencoder.

[0151] For some typical fault scenarios, such as IDV(13) representing slow drift of reaction kinetics, this fault involves a multi-variable coupling process. IDV(14) represents the sticking of the reactor cooling valve, which has concealment and requires high precision for chemical process anomaly detection. From Table 4, the accuracy of IDV(13) is 0.972, and the F1 score is 0.987. The accuracy of IDV(14) is 0.927, and the F1 score is 0.958, which is more than 5% higher than the LSTM-AE baseline model. For small anomalies IDV(3), IDV(9), IDV(15), the detection performance of the ST-GVAE model is lower, but it also has certain improvement compared with the baseline model. For unknown faults, the model has also improved compared with other model-based models. In summary, ST-GVAE balances the precision and interpretability in the TE process fault detection task by fusing chemical mechanism and spatio-temporal graph variational encoding mechanism.

[0152] Anomaly variable positioning

[0153] To locate the detected anomaly to the relevant variable, based on the 20 fault types contained in the TE process fault case, according to the introduction mode of the fault, it can be divided into 5 categories, of which the first 4 categories are known faults. Based on the propagation mechanism of the fault, this study selects four typical faults IDV(1), IDV(8), IDV(13), IDV(14) to verify the anomaly localization ability of ST-GVAE. This study generates variable-level bias score B by calculating the spatiotemporal reconstruction error, as shown in the formula. By analyzing the reconstruction error of each variable and sorting combined with the mechanism of chemical process, the cause of the fault can be located to the relevant variable. By visualizing the bias score B of the original test data and the model reconstruction data, the anomaly of each fault type can be clearly located to the relevant variable.

[0154] 1) IDV(1): Step anomaly of feed ratio

[0155] IDV(1) fault type is caused by the step of A / C feed ratio in feed stream 4, which reduces the material of A in the recycle stream, and finally affects the composition of A in material stream 6. The feedback controller triggers the manipulated feed flow valve variable (24) to increase the A feed flow variable (1). This disturbance propagates throughout the process, affecting various variables until the controller resolves the impact of the disturbance and establishes a new steady state.

[0156] As shown in the IDV(1) figure in Figure 10 , the anomaly detection model detects that the anomaly score breaks through the threshold at 161 sample points. After detecting the anomaly, the model performs anomaly localization to locate the detected fault to the corresponding sensor variable.

[0157] The visualization of the bias score of the reconstruction data and the test data of IDV(1) is shown in Figure 10 , based on the bias score analysis of the variable, the fault can be located to the measured variables (1), (4), (7), (11), (13), (16), (18), (21), and the operating variables (24), (25). Among them, according to the severity of the bias, the abnormal variables are divided into primary abnormal variables: (1), (24) and secondary abnormal variables: (4), (7), (11), (13), (16), (18), (21), (25).

[0158] Further analysis based on chemical mechanism revealed a fault in the TE process, causing the feedback controller to require continuous control of the feed rate of component A in stream 1 via the operational variable (24). This further caused changes in the reactor pressure variable (7), the distillation column pressure variable (16), and the gas-liquid separator pressure, which in turn caused temperature changes. Except for the measured variable (1) and the corresponding operational variable (24), which required continuous adjustment, all other variables returned to normal. Based on the constructed graph structure and the sequence of variable information transmission signals, it can be found that the measured variable (1) and the corresponding operational variable (24) are the starting nodes of the nodes, and the deviation between these two variables is continuously abnormal. Therefore, the cause of the fault is located to the abnormal feed of component A in stream 4, which requires stream 1 to change the feed rate of component A to restore system stability. This is basically consistent with the cause of the fault.

[0159] 2) IDV(8): Random perturbation of feed ratio

[0160] The IDV(8) failure is caused by random variations in the feed ratio of components A, B, and C in stream 4, such as Figure 9 As shown in IDV(8), the model detected abnormal scores exceeding the threshold at 176 time sample points. The detection time delay was 15 sample points, totaling 45 minutes.

[0161] Visualization of the deviation score between the reconstructed data and the test data of IDV(8) is as follows: Figure 11 As shown, based on the deviation score analysis of variables, the fault can be located to the measured variables (1), (3), (4), (7), (8), (10), (11), (13), (16), (18), (19), (20), (21) and the manipulated variables (23), (24), (25), (26), (29). Among them, according to the severity of the deviation, the abnormal variables are divided into primary abnormal variables: (4), (7), (11), (13), (16), (18), (21), (26), (29) and secondary abnormal variables: (1), (3), (8), (10), (19), (23), (24), (25).

[0162] based on Figure 11The visualization results analysis revealed that the abnormality was not continuous but rather randomly changing. Further analysis based on variables and graph structure identified the source of the fault propagation as the feed flow rate of A in stream 1 and the operating variable (24), and the feed flow rate of stream 4 and the corresponding operating variable (25). This triggered corresponding pressure anomalies (7), (13), (16), and power anomalies (19). The pressure affected the temperatures of the corresponding units (11), (18), (20), and (21). After the temperature anomaly, the feedback regulation was to restore the temperature to the normal range, so the operating variable (29) related to temperature control became abnormal. At the same time, this fault source also caused abnormal liquid level (8) in the abnormal reactor and abnormal discharge flow rate (10) of the gas-liquid separator. Thus, the cause of the fault can be located to the random changes of the measured variables and operating variables (24) in stream 1 and the measured variables (4) and the corresponding operating variables (25) in stream 4, which is consistent with the cause of the fault.

[0163] 3) IDV(13): Slow drift in reaction kinetics

[0164] IDV(13) failures are caused by slow drift in reactor reaction kinetics and are a type of multivariable coupled failure. For example... Figure 9 As shown in IDV(13), the model detected that the abnormal score exceeded the threshold at the 210th time sample point. The detection delay was 49 sample points, or 2.45h.

[0165] Visualization of the deviation score between the reconstructed data and the test data of IDV(13) as follows: Figure 12 As shown, based on the deviation score analysis of the variables, the fault can be located to the measured variables (7), (10), (11), (13), (16), (18), (19), (20), (21), (23) and the manipulated variables (26). Among them, according to the severity of the deviation, the abnormal variables are divided into primary abnormal variables: (7), (10), (11), (13), (16), (18), (19), (26) and secondary abnormal variables: (19), (20), (23).

[0166] Based on the graph structure and process mechanism of the TE process, the source of the anomaly can be located to the relevant variables of the reactor. Due to feedback regulation, the generation of the anomaly does not have a large impact on the flow variables input to the reactor. However, the pressure measurement variables (7), (13), (16) of each unit, as well as the compressor power (19) closely related to the pressure, the vent flow (10) and the operating variable (26) of the vent flow, generate anomalies on the basis of constant flow. In turn, the relevant temperature variables (11), (18), (21) are affected. According to the order of information transmission combined with the mode of feedback regulation, the anomaly can be located to the relevant variables of the reactor.

[0167] 4) IDV (14): Sticking fault of cooling valve

[0168] The IDV (14) fault is caused by the sticking of the reactor cooling water valve. As shown in IDV (14) in Table 1, the model detects that the anomaly score exceeds the threshold at the 161st time sample point, and the delay of anomaly detection is 0, which has timely detection performance. Figure 9

[0169] The deviation score visualization of the reconstructed data of IDV (14) and the test data is shown in FIG. 6. Based on the variable deviation score analysis, the fault can be located to the measurement variables (9), (20) and the operating variable (29). The variables are the reactor temperature (9), the outlet temperature of the reaction cooling water (20) and the reactor cooling water flow (29) respectively. It can be observed from FIG. 6 that the reactor temperature (9) and the reactor cooling water flow (29) first have anomalies, which cause the outlet temperature of the reaction cooling water (20) to have an anomaly. Therefore, the root cause of the anomaly is determined to be related to the cooling of the reactor. Figure 13 Figure 13

[0170] The average anomaly detection accuracy of the ST-GVAE model is 87.1%, and the average anomaly detection F1 score is 89.8%. The accuracy and F1 score of the model are improved by 8% and 2.8% compared with the LSTM-GCN-AE respectively. Based on the detected anomaly score, the variable deviation score is constructed to locate the fault to the relevant variables and visualize, and further through the information transmission of the variables and the process mechanism of the chemical process, the variables are further located to a smaller variable range.

[0171] Finally, it should be pointed out that the above embodiments are only used to illustrate the technical solutions of the present application and not to limit it. Although the present application has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present application can be modified or replaced by equivalents without departing from the purpose and scope of the present application, which should be covered in the scope of the claims of the present application.​​​

Claims

1. A method for chemical production anomaly detection based on spatio-temporal graph variational encoder, characterized in that: Comprise the following steps: S100: select the disclosed chemical normal condition data D, D includes all variables history data and all variables PID chart under normal condition; S200: build an abnormal monitoring model ST-GVAE, ST-GVAE includes data preprocessing module and abnormal detection module, the data preprocessing module is used to preprocess the chemical normal condition data, preprocessing includes topological graph transformation, data feature extraction and data format standardization conversion; The abnormal detection module is composed of LSTM network embedded GCN (LSTM-GCN), LSTM-GCN includes encoder, one MLP and decoder, the encoder and decoder have the same network structure;Abnormal detection module is used to calculate the abnormal detection score of target variable in chemical process and the positioning result of target variable corresponding to sensor; S300: build loss function L(·), use D to train ST-GVAE, adopt Adam optimizer, update ST-GVAE parameters by gradient descent in reverse, stop training when the loss function reaches minimum or reaches maximum training iteration times, get trained ST-GVAE', the calculation expression of loss function is as follows: where θ denotes the decoder network parameters, denotes the encoder network parameters, x denotes the nodes, L denotes the number of samples estimating the expectation value, l denotes the l-th sample, D KL denotes the KL divergence term, D KL by constraining the latent variable approximates the prior distribution p θ ; S400: input all variables history data under normal condition into ST-GVAE', get the abnormal detection score corresponding to all variables; S500: select the target variable a to be detected, input the normal condition history data corresponding to a into ST-GVAE', get the abnormal detection score corresponding to a, take the abnormal detection score as the abnormal detection threshold of a; Then input the current chemical condition data corresponding to a into ST-GVAE', get the abnormal detection score of current chemical condition corresponding to a, when the abnormal detection score is greater than the abnormal detection threshold of a, trigger alarm and proceed to the next step, otherwise, it means that the target variable detection is normal, and the detection is finished; S600: take the sensor position in the PID chart corresponding to the target variable a to be detected as node, define deviation threshold, and calculate the reconstruction deviation score B of each node, take the node corresponding to B greater than deviation threshold as abnormal node, and then output the sensor position corresponding to abnormal node.

2. The chemical production anomaly detection method based on spatio-temporal graph variational encoder of claim 1, wherein: The preprocessing step of the chemical normal condition data in S200 is as follows: S210: define the chemical priori knowledge as G=(V,E,A,X), wherein, V represents the node set composed of sensors, E represents the edge set, A represents the N×N dimensional adjacency matrix, and X represents the node feature matrix of graph structure; S211: topological graph transformation is carried out on PID chart by using G through chemical mechanism, and topological structure chart is obtained; S212: select target variable data for modeling from all variables history data under normal condition as variable data feature, which is used as data feature extraction value; S213: adopt filtering method to denoise topological structure chart and data feature extraction value, then carry out format standardization conversion on denoised topological structure chart and data feature extraction value, finally obtain standardized format data required by ST-GVAE.

3. The chemical production anomaly detection method based on spatio-temporal graph variational encoder of claim 2, wherein: In the S300, the ST-GVAE is trained by using an Adam optimizer, and the steps are as follows: In the S310, the data preprocessing module is used to process the D to obtain corresponding standardized format data, the standardized format data is divided into time windows, and a spatio-temporal graph data set is constructed, and the time window is defined as follows: X t-w+1:t = [X t-w+1 ,..., X t-1 , X t ] wherein X t represents the time window data at time t, X t-w+1:t represents a time window interval, and w represents the time window length; S320: send the connection relationship A between each time window and node in the spatio-temporal graph dataset to the encoder, generate the hidden state sequence [h t-w+1 ,..., h t-1 , h t ] step by step through LSTM-GCN, and take the hidden state h t at the final time step as the global spatio-temporal feature representation; wherein the hidden state at each time step in the hidden state sequence is updated as follows: h T = LSTM-GCN(X T ,A,h T-1 ), T e [t-w+1, t] Wherein, T represents the time step of the next hidden state; S330: h t input to the MLP to compute the set of latent variables Z t The computation expression is as follows: mu i = W μ (i) • h t (i) + b μ (i) σ i = W σ (i) • h t (i) + b σ (i) Z t (i) =μ i +s i ,s i ~N(0,1) Z t = {Z t (1) ,…,Z t (i) ,…Z t (N)} where μ i represents the hidden variable generation mean parameter for node i, σ i represents the hidden variable standard deviation parameter for i, W μ (i) , W σ (i) represents the weight parameter, b μ (i) , b σ (i) represents the bias term, q(Z t | X t-w+1：t ) represents the global hidden variable distribution, represents the hidden variable for node i, N(0, 1) represents a normal distribution; S340: obtain Z t The abnormality detection scores of all variables at time t in the input decoder are obtained. L(·) is calculated, and it is judged whether the training stop condition is reached. If yes, the current ST-GVAE is taken as the trained ST-GVAE', otherwise, the S320 is returned. The training stop condition is that the loss function reaches the minimum or reaches the maximum training iteration number.

4. The chemical production anomaly detection method based on spatio-temporal graph variational encoder of claim 3, wherein: In the S600, the calculation of the reconstruction deviation score B of each node is as follows: The reverse temporal sequence is reconstructed in reverse in the decoder from time step t to t-w+1, and then Z t As input, Z t is mapped back to the original data space according to the reverse temporal sequence to obtain reconstructed data; The absolute error square of the corresponding normal working condition historical data a and the reconstruction data a is calculated, and the absolute error square is taken as the reconstruction deviation score B of the node in the PID graph corresponding to a. The calculation formula is as follows: B = (x' - x) / (y - y') t - x t ) 2 where x' = x - x0 t represents the reconstructed data, x t represents the original operating condition data.

Images