A method for detecting slowly changing faults in chemical processes based on typical variable heterogeneity analysis

Abstract

The present invention relates to a chemical process slow-changing fault detection method based on weighted probability typical variable dissimilarity analysis, the steps of which are: standardizing the training data and constructing a historical data matrix and a future data matrix, establishing a CVDA model and extracting the CVDA features of the training data, and calculating the corresponding training WD features, and further calculating the statistics T corresponding to the training WD features. 2 , statistics Q and statistics D and determine the corresponding control limits; collect test data, construct historical data vectors and future data vectors of test data, use the established CVDA model to extract the corresponding online CVDA features, calculate the online WD features, further calculate the online weighted WD features, and calculate the statistic T based on the online weighted WD features 2 , statistics Q and statistics D, and monitor using control limits. The fault detection method of the present invention introduces Wasserstein distance and online weighting strategy, which improves the detection performance of slowly varying faults in chemical process data.

Description

Technical Field

[0001] The present invention belongs to the technical field of chemical process fault detection, and relates to a chemical process slow-changing fault detection method based on canonical variable dissimilarity analysis (CVDA for short). Background Art

[0002] With the continuous development of science and technology, the scale and complexity of modern chemical systems are constantly increasing, and at the same time, higher requirements are put forward for the system performance and production safety of chemical processes. It is precisely because of the large-scale and high complexity of modern chemical systems that fault events inevitably occur frequently. Slow-varying faults have the characteristics of slow change and small initial amplitude. If slow-varying faults can be detected in time within the controllable range of the chemical process and corresponding measures can be taken, it will effectively reduce productivity losses and ensure the safety of workers. It has important practical significance for the safety monitoring of chemical processes. With the development of computer technology, a large amount of data of chemical processes has been stored, which has prompted the widespread application of data-driven fault diagnosis technology.

[0003] Canonical variable dissimilarity analysis is a new algorithm in the field of chemical process monitoring, first proposed by Pilario and Cao in 2018. However, in process monitoring, the traditional CVDA method only uses the original measurement data to construct statistics and all CVDA features are modeled and monitored on an equal footing, without considering the differences in the probability information contained in the data and the fault information carried by the features, resulting in poor detection of slow-varying faults. Therefore, how to consider the probability information and feature differences contained in the process data to improve the detection performance of slow-varying faults is a key issue that needs to be solved in chemical process monitoring. Summary of the invention

[0004] The present invention aims at the problem that the traditional CVDA method does not consider the probability information included in the process data and the difference of the fault information carried by the features, and provides a chemical process slow-changing fault detection method based on typical variable dissimilarity analysis. The method introduces Wasserstein distance (WD) into the CVDA method, uses WD to mine the probability information contained in the slow-changing fault data, and further considers the feature difference, introduces the weight coefficient, and performs online weighting on the features according to the size of the fault information carried by the features, so as to achieve the purpose of improving the performance of slow-changing fault detection.

[0005] In order to achieve the above object, the present invention provides a method for detecting slow-changing faults in a chemical process based on typical variable heterogeneity analysis, comprising the following steps:

[0006] (I) Collect process data under normal operating conditions in the historical database as training data X 0 , and using the training data X 0 The mean and standard deviation For the training data X 0 Perform standardization to obtain standardized training data X; construct training history data matrix P and future data matrix F from training data X;

[0007] (ii) performing canonical variable dissimilarity analysis on the historical and future data matrices P and F to obtain a CVDA model, extracting the state vector s, residual vector e and dissimilarity feature d of the historical and future data matrices P and F using the CVDA model, and calculating the training WD feature s corresponding to the state vector s, residual vector e and dissimilarity feature d using the Wasserstein distance. WD 、e WD and d WD ;

[0008] (III) Based on the training WD feature s WD 、e WD and d WD Calculate the corresponding statistic T 2 , statistics Q and statistics D, given the confidence level α, calculate the statistic T by the kernel density estimation method 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim ;

[0009] (IV) Collect test data x new , using the training data X 0 The mean and standard deviation For the test data x new Perform standardization to obtain the standardized test data x t , from the test data x t Construct test history data vector p t and the future data vector f t ;

[0010] (V) Using the CVDA model obtained in step (II), obtain the test history data vector p t and the future data vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t , and use Wasserstein distance to calculate the state vector s t, residual vector e t and heterogeneity characteristics d t Corresponding online WD features and

[0011] (VI) Based on the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , obtain online weighted WD features and

[0012] (VII) The online weighted WD feature and Calculate the corresponding new statistic T t 2 , statistic Q t and statistic D t , and according to the statistic T t 2 , statistic Q t and statistic D t Whether it exceeds the corresponding control limit, judge the test data x t Whether a fault occurs.

[0013] Furthermore, in step (a), the training data X 0 The specific process of standardization and establishing the training history data matrix P and future data matrix F is as follows:

[0014] First, use the training data X 0 The mean and standard deviation The training data X is calculated by formula (1) 0 After standardization, the expression of formula (1) is:

[0015]

[0016] Training data X 0 After the normalization process according to the above formula (1), the normalized training data X can be obtained.

[0017] The training data X=[x 1 ,x 2 ,...,x i ,…,x n ] T ∈R n×m Medium i The corresponding training history data vector p i and the future data vector f i , where i represents the sampling time, x i represents the process variable collected at the i-th moment, n represents the number of samples, and m represents the number of variables. Formula (2) and formula (3) are expressed as follows:

[0018]

[0019]

[0020] In the formula, h represents the time lag number,

[0021] The historical data matrix P and the future data matrix F corresponding to the training historical data vector and the future data vector are calculated by formulas (4) and (5), which are expressed as follows:

[0022]

[0023]

[0024] Where, M = n-2h + 1, p h+1 represents the historical data vector corresponding to the h+1th moment, f h+1 represents the future data vector corresponding to the h+1th moment;

[0025] Furthermore, in the step (ii), the state vector s, residual vector e and dissimilarity feature d of the historical and future data matrices P and F are extracted for establishing the CVDA model, and the training WD feature s corresponding to the state vector s, residual vector e and dissimilarity feature d is calculated. WD 、e WD and d WD The specific steps are:

[0026] First, the covariance matrix Σ corresponding to the historical and future data matrices P and F is calculated by formulas (6)-(8): pp and Σ ff and the cross-covariance matrix Σ pf , formulas (6)-(8) are expressed as:

[0027]

[0028]

[0029]

[0030] For the covariance matrix Σ pp and Σ ff and the cross-covariance matrix Σ pf Carry out the generalized eigenvalue decomposition shown in formula (9), the expression of formula (9) is:

[0031]

[0032] Where U represents the left singular matrix, V represents the right singular matrix, The eigenvalue The diagonal matrix formed by descending order;

[0033] By solving formula (9), we can carry out formulas (10) and (11) to calculate the projection matrix and Formulas (10) and (11) are expressed as follows:

[0034]

[0035]

[0036] The projection vector pair {a i ,b i} corresponds to the i-th correlation coefficient λ i , according to λ i The size of the value, the projection vector pair is divided into a typical variable subspace and a residual subspace, the first k eigenvalues ​​λ 1 >λ 2 >…>λ k Represents the main correlation coefficient corresponding to the typical variable subspace projection matrix A k =[a 1 ,a 2 ,…,a k ], the remaining Eigenvalues represents the residual subspace projection matrix corresponding to the weak correlation coefficient, and the value of k is determined by the average method;

[0037] For the training history data vector p corresponding to the i-th sampling time in the training data X i and the future data vector f i , extract the corresponding state vector s through formulas (12)-(13) i , residual vector e i and heterogeneity characteristics d i, the expressions of formulas (12)-(13) are:

[0038]

[0039]

[0040]

[0041] Where U k represents the first k columns of the matrix U, and the dissimilarity feature measure is the future data vector f i The predicted typical variables and historical data vector p i The difference in canonical variates obtained.

[0042] Calculate the state vector s = [s 1 ,s 2 ,…,s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ;

[0043] Use sliding window to obtain state vector s i,j The local mean and variance of the state vector s is calculated by formula (15): i,j Corresponding WD characteristics The expression of formula (15) is:

[0044]

[0045] In the formula, Represents the state vector s required using the sliding window i,j The local mean of Represents the state vector s required using the sliding window i,j The local variance of .

[0046] Train WD features by computing Similarly, the residual vector e is calculated i and heterogeneity characteristics d i The corresponding training WD features and The residual vector e i,j Corresponding WD characteristics Dissimilarity feature d i,j Corresponding WD characteristics

[0047] Furthermore, in step (iii), the training WD feature s is calculated WD 、e WD and dWD The corresponding statistic T 2 The specific steps of calculating statistics Q and statistics D are:

[0048] Based on WD characteristics WD 、e WD and d WD , the statistics T are constructed by formulas (16)-(18) respectively 2 , statistics Q and statistics D, the expressions of formulas (16)-(18) are:

[0049]

[0050]

[0051]

[0052] In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d.

[0053] Given a confidence level α, the statistic T is calculated using the kernel density estimation method 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim .

[0054] Furthermore, in step (iv), using the training data X 0 The mean and standard deviation The test data x is calculated by formula (19) new After standardization, the expression of formula (19) is:

[0055]

[0056] The test data x collected at sampling time t new After the above formula (19) is used for standardization, the standardized test data x can be obtained. t .

[0057] The test data x t Corresponding test history data vector p t and the future data vector f t Calculated by formula (20) and (21), the expressions of formula (20) and (21) are:

[0058]

[0059]

[0060] Furthermore, in step (V), the test history vector p t and the future vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t Calculated by formula (22)-(24), the expressions of formula (22)-(24) are:

[0061]

[0062]

[0063]

[0064] Where U k To solve formula (9), we obtain the first k columns of the left singular matrix U, A k and B k To solve equations (10) and (11) to obtain the first k pairs of projection vectors;

[0065] The state vector s=[s 1 ,s 2 ,…,s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ;

[0066] Use sliding window to obtain state vector s t,j The local mean and variance of the state vector s is calculated by formula (25) t,j Corresponding WD characteristics The expression of formula (25) is:

[0067]

[0068] In the formula, Represents the state vector s required using the sliding window t,j The local mean of Represents the state vector s required using the sliding window t,j The local variance of .

[0069] Testing WD characteristics by calculation Similarly, the residual vector e is calculated t and heterogeneity characteristics d t The corresponding training WD features and The residual vector e t,j Corresponding WD characteristics Dissimilarity feature d t,j Corresponding WD characteristics

[0070] Furthermore, in step (six), the online weighted WD feature is obtained and The specific steps are:

[0071] First, according to the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , calculated by formula (26)-(28), the expressions of formula (26)-(28) are:

[0072]

[0073]

[0074]

[0075] In the formula, the weight w s,t =[w s,t (1),…,w s,t (j),…,w s,t (k)] T w in s,t (j) represents the jth online WD feature The weighting coefficient, weight w in e,t (j) represents the jth online WD feature The weighting coefficient, weight w d,t =[w d,t (1),…,w d,t (j),…,w d,t (k)] T w in d,t (j) represents the jth online WD feature The weighting coefficient of

[0076] The weighted threshold and Generally choose to train WD features and 1 to 2 times the maximum value among the elements.

[0077] Get online weighted WD features and Calculated by formula (29)-(31), the expressions of formula (29)-(31) are:

[0078]

[0079]

[0080]

[0081] Furthermore, in step (VII), the online weighted WD feature and Calculate the corresponding new statistic T t 2 , statistic Q t and statistic D t The specific steps are:

[0082] Based on online weighted WD features and The statistics T are constructed by formulas (32)-(34) respectively 2 , statistics Q and statistics D, the expressions of formulas (32)-(34) are:

[0083]

[0084]

[0085]

[0086] In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d.

[0087] Determine the test data x t The steps to check if a fault occurs are: Q t ≤Q lim And D t ≤D lim When the chemical process is in normal working condition, it is considered that the chemical process is in normal working condition; otherwise, it is considered that the chemical process is in failure.

[0088] Compared with the prior art, the present invention has the following beneficial effects:

[0089] The chemical process slow-varying fault detection method based on typical variable dissimilarity analysis provided by the present invention utilizes the Wasserstein distance to measure the probability distribution changes of state vectors, residuals and dissimilarity features, thereby realizing the mining of process data probability information. Meanwhile, the feature differences are further considered, and the features are weighted by weight coefficients to increase the weights of features carrying more fault information, so that the statistics can more obviously reflect the slow-varying fault information in the chemical process, thereby improving the slow-varying fault detection results and increasing the slow-varying fault detection rate. BRIEF DESCRIPTION OF THE DRAWINGS

[0090] Figure 1 It is a flow chart of the chemical process slow-changing fault detection method based on typical variable heterogeneity analysis of the present invention;

[0091] Figure 2 A structural diagram of the Tennessee-Eastman (TE) process according to an embodiment of the present invention;

[0092] Figure 3a This is a schematic diagram of monitoring results of a TE process fault 21 using a traditional CVDA method according to an embodiment of the present invention;

[0093] Figure 3b This is a schematic diagram of monitoring results of a TE process fault 21 using the chemical process slow-changing fault detection method based on typical variable heterogeneity analysis according to an embodiment of the present invention; DETAILED DESCRIPTION

[0094] The present invention is described in detail below by way of exemplary embodiments. However, it should be understood that elements, structures, and features in one embodiment may also be beneficially combined in other embodiments without further description.

[0095] See also Figure 1 The present invention discloses a method for detecting slow-changing faults in a chemical process based on typical variable heterogeneity analysis, which comprises the following steps:

[0096] (I) Collect process data under normal operating conditions in the historical database as training data X 0 , and using the training data X 0 The mean and standard deviation For the training data X 0 Perform standardization to obtain standardized training data X; construct training history data matrix P and future data matrix F from training data X. The specific steps are as follows:

[0097] Using training data X 0 The mean and standard deviation The training data X is calculated by formula (1)0 After standardization, the expression of formula (1) is:

[0098]

[0099] The training data X=[x 1 ,x 2 ,...,x i ,…,x n ] T ∈R n×m Medium i The corresponding training history data vector p i and the future data vector f i , where i represents the sampling time, x i represents the process variable collected at the i-th moment, n represents the number of samples, and m represents the number of variables. Formula (2) and formula (3) are expressed as follows:

[0100]

[0101]

[0102] In the formula, x i represents the process variable collected at the i-th moment, h represents the time lag number,

[0103] The historical data matrix P and the future data matrix F corresponding to the training historical data vector and the future data vector are calculated by formulas (4) and (5), which are expressed as follows:

[0104]

[0105]

[0106] Where, M = n-2h + 1, p h+1 represents the historical data vector corresponding to the h+1th moment, f h+1 represents the future data vector corresponding to the h+1th moment;

[0107] (ii) performing canonical variable dissimilarity analysis on the historical and future data matrices P and F to obtain a CVDA model, extracting the state vector s, residual vector e and dissimilarity feature d of the historical and future data matrices P and F using the CVDA model, and calculating the training WD feature s corresponding to the state vector s, residual vector e and dissimilarity feature d using the Wasserstein distance. WD 、e WD and d WD , the specific steps are:

[0108] First, the covariance matrix Σ corresponding to the historical and future data matrices P and F is calculated by formulas (6)-(8): pp and Σ ff and the cross-covariance matrix Σ pf , formulas (6)-(8) are expressed as:

[0109]

[0110]

[0111]

[0112] For the covariance matrix Σ pp and Σ ff and the cross-covariance matrix Σ pf Carry out the generalized eigenvalue decomposition shown in formula (9), the expression of formula (9) is:

[0113]

[0114] Where U represents the left singular matrix, V represents the right singular matrix, The eigenvalue The diagonal matrix formed by descending order;

[0115] By solving formula (9), we can carry out formulas (10) and (11) to calculate the projection matrix and Formulas (10) and (11) are expressed as follows:

[0116]

[0117]

[0118] The projection vector pair {a i ,b i} corresponds to the i-th correlation coefficient λ i , according to λ i The size of the value, the projection vector pair is divided into a typical variable subspace and a residual subspace, the first k eigenvalues ​​λ 1 >λ 2 >…>λ k Represents the main correlation coefficient corresponding to the typical variable subspace projection matrix A k =[a 1 ,a 2 ,…,a k ], the remaining Eigenvalues represents the residual subspace projection matrix corresponding to the weak correlation coefficient, and the value of k is determined by the average method;

[0119] For the training history data vector p corresponding to the i-th sampling time in the training data X i and the future data vector f i , extract the corresponding state vector s through formulas (12)-(13) i , residual vector e i and heterogeneity characteristics d i , the expressions of formulas (12)-(13) are:

[0120]

[0121]

[0122]

[0123] Where U k represents the first k columns of the matrix U, and the dissimilarity feature measure is the future data vector f i The predicted typical variables and historical data vector p i The difference in canonical variates obtained.

[0124] Calculate the state vector s = [s 1 ,s 2 ,…,s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ;

[0125] Use sliding window to obtain state vector s i,j The local mean and variance of the state vector s is calculated by formula (15): i,j Corresponding WD characteristics The expression of formula (15) is:

[0126]

[0127] In the formula, Represents the state vector s required using the sliding window i,j The local mean of Represents the state vector s required using the sliding window i,j The local variance of .

[0128] Train WD features by computing Similarly, the residual vector e is calculated i and heterogeneity characteristics d i The corresponding training WD features and The residual vector ei,j Corresponding WD characteristics Dissimilarity feature d i,j Corresponding WD characteristics

[0129] (III) Calculate and train WD features s WD 、e WD and d WD The corresponding statistic T 2 , statistics Q and statistics D and calculate statistics T 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim , the specific steps are:

[0130] Based on WD characteristics WD 、e WD and d WD , the statistics T are constructed by formulas (16)-(18) respectively 2 , statistics Q and statistics D, the expressions of formulas (16)-(18) are:

[0131]

[0132]

[0133]

[0134] In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d.

[0135] Given a confidence level α, the statistic T is calculated using the kernel density estimation method 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim .

[0136] (IV) Collect test data x new , using the training data X 0 The mean and standard deviation For the test data x new Perform standardization to obtain the standardized test data x t , from the test data x t Construct test history vector p t and the future vector f t , the specific steps are:

[0137] Using training data X 0 The mean and standard deviation The test data x is calculated by formula (19) new After standardization, the expression of formula (19) is:

[0138]

[0139] The test data x collected at sampling time t new After the above formula (19) is used for standardization, the standardized test data x can be obtained. t .

[0140] The test data x t Corresponding test history data vector p t and the future data vector f t Calculated by formula (20) and (21), the expressions of formula (20) and (21) are:

[0141]

[0142]

[0143] (V) Using the CVDA model obtained in step (II), obtain the test history data vector p t and the future data vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t , and use Wasserstein distance to calculate the state vector s t , residual vector e t and heterogeneity characteristics d t Corresponding online WD features and The specific steps are:

[0144] First, the test history data vector p t and the future data vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t Calculated by formula (22)-(24), the expressions of formula (22)-(24) are:

[0145]

[0146]

[0147]

[0148] Where Uk To solve formula (9), we obtain the first k columns of the left singular matrix U, A k and B k To solve equations (10) and (11) to obtain the first k pairs of projection vectors;

[0149] The state vector s=[s 1 ,s 2 ,…,s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ;

[0150] Use sliding window to obtain state vector s t,j The local mean and variance of the state vector s is calculated by formula (25) t,j Corresponding WD characteristics The expression of formula (25) is:

[0151]

[0152] In the formula, Represents the state vector s required using the sliding window t,j The local mean of Represents the state vector s required using the sliding window t,j The local variance of .

[0153] Testing WD characteristics by calculation Similarly, the residual vector e is calculated t and heterogeneity characteristics d t The corresponding training WD features and The residual vector e t,j Corresponding WD characteristics Dissimilarity feature d t,j Corresponding WD characteristics

[0154] (VI) Based on the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , obtain online weighted WD features and The specific steps are:

[0155] First, according to the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , calculated by formula (26)-(28), the expressions of formula (26)-(28) are:

[0156]

[0157]

[0158]

[0159] In the formula, the weight w s,t =[w s,t (1),…,w s,t (j),…,w s,t (k)] T w in s,t (j) represents the jth online WD feature The weighting coefficient, weight w in e,t (j) represents the jth online WD feature The weighting coefficient, weight w d,t =[w d,t (1),…,w d,t (j),…,w d,t (k)] T w in d,t (j) represents the jth online WD feature The weighting coefficient of

[0160] The weighted threshold and Generally choose to train WD features and 1 to 2 times the maximum value among the elements.

[0161] Get online weighted WD features and Calculated by formula (29)-(31), the expressions of formula (29)-(31) are:

[0162]

[0163]

[0164]

[0165] (VII) The online weighted WD feature and Calculate the corresponding new statistic T t 2 , statistic Q t and statistic D t , and according to the statistic T t 2 , statistic Q t and statistic D t Whether it exceeds the corresponding control limit, judge the test data x t Whether a fault occurs, the specific steps are as follows:

[0166] Based on online weighted WD features and The statistics T are constructed by formulas (32)-(34) respectively 2 , statistics Q and statistics D, the expressions of formulas (32)-(34) are:

[0167]

[0168]

[0169]

[0170] In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d.

[0171] Determine the test data x t The steps to determine whether a fault has occurred are: when T t 2 ≤T li 2 m , Q t ≤Q lim And D t ≤D lim When the chemical process is in normal working condition, it is considered that the chemical process is in normal working condition; otherwise, it is considered that the chemical process is in failure.

[0172] In the above method, steps (1) to (3) are the offline modeling stage, and steps (4) to (7) are the online monitoring stage.

[0173] The fault detection method of the present invention standardizes the training data, uses the standardized training data to establish historical and future data matrices, and then establishes a CVDA model based on the historical and future data matrices, extracts CVDA features from the training data, and calculates the corresponding training WD features, and further calculates the statistics T corresponding to the training WD features. 2 , statistics Q and statistics D and determine the corresponding control limits; collect test data, construct historical data vectors and future data vectors of test data, use the established CVDA model to extract the corresponding online CVDA features, calculate the online WD features, further calculate the online weighted WD features, and calculate the statistic T based on the online weighted WD features 2 , statistics Q and statistics D, and monitor using control limits. The fault detection method of the present invention introduces Wasserstein distance and online weighting strategy, which improves the detection performance of slowly varying faults in chemical process data.

[0174] In order to more clearly illustrate the beneficial effects of the above-mentioned fault detection method of the present invention, the above-mentioned fault detection method of the present invention is further described below in conjunction with embodiments.

[0175] Example:

[0176] The Tennessee-Eastman (TE) process is an experimental platform established by Downs and Vogel of Eastman Chemical Company in the United States based on an actual chemical process. It is now widely used to verify the advantages and disadvantages of control algorithms and process monitoring methods. Figure 2 The TE process is mainly composed of five units, including a reactor, a product condenser, a gas-liquid separator, a circulating compressor and a stripping tower. The TE process has a total of 53 variables, including 22 continuously measured variables, 19 component variables and 12 operating variables. In this embodiment, see Table 1, select 33 variables in the TE process; see Table 2, there are a total of 4 faults.

[0177] Table 1

[0178] Variable Label Variable Description Variable Label Variable Description 1 A feed (stream 1) 18 Stripper temperature 2 D Feed (Stream 2) 19 Stripper flow 3 E feed (stream 3) 20 Compressor power 4 A and C feed (stream 4) 21 Reactor cooling water outlet temperature 5 Recirculation flow (flow 8) 22 Separator cooling water outlet temperature 6 Reactor feed rate (stream 6) 23 D Feed amount (flow 2) 7 Reactor pressure 24 E feed amount (stream 3) 8 Reactor level 25 A feed amount (flow 1) 9 Reactor temperature 26 A and C material quantity (flow 4) 10 Venting rate 27 Compressor recirculation valve 11 Product separator temperature 28 Vent valve 12 Product separator level 29 Separator tank liquid flow (flow 10) 13 Product separator pressure 30 Stripper liquid product flow (stream 11) 14 Product separator bottom flow (stream 10) 31 Stripper water flow valve 15 Stripper Grade 32 Reactor cooling water flow 16 Stripper pressure 33 Condenser cold water flow 17 Stripper bottom flow (flow 11)

[0179] Table 2

[0180] Fault name describe type Fll Reactor cooling water inlet temperature random F15 Condenser cooling water valve Sticky F19 unknown unknown F21 Stream 4 valve sticks Constant

[0181] In this embodiment, the traditional CVDA method and the WPCVDA method of the present invention are used for simulation comparison. Three performance indicators, namely the fault detection rate FDR, the fault false alarm rate FAR and the fault detection time FDT, are used to evaluate the fault detection performance of different fault detection methods. Specifically, the fault detection rate FDR is defined as the percentage of fault samples exceeding the control limit to the actual total fault samples; the fault false alarm rate FAR is defined as the percentage of normal samples exceeding the control limit to the actual total normal samples; the fault detection time FDT is defined as the first fault sample when six consecutive fault samples exceed the control limit. A 99% confidence level is set to obtain the statistical confidence limit. In the WPCVDA method of the present invention, the sliding window width is 30, the feature weighted threshold is set to 1.5 times, the first 160 samples of the test data are normal samples, and the subsequent samples are fault samples.

[0182] Take fault 21 as an example to illustrate the fault detection effect. Fault 21 is caused by the sticking of the valve of flow 4. The conventional CVDA method and the WPCVDA method of the present invention are used for simulation comparison. For the monitoring diagram of fault 21, see Figure 3a and Figure 3b See also Figure 3a ,Due to the slow change of faults, the traditional CVDA method has a certain detection effect but the detection effect is not good. Its statistic D has a detection rate of 74.75% and statistic T 2 The detection rate is 46.88%, the detection rate of statistic Q is 38.38%, and the earliest fault is detected in the 360th sample. Based on the CVDA method, the WPCVDA provided by the present invention considers the probability information and feature differences of process data, and the detection rate of statistic D is increased to 79.13%, and the detection rate of statistic T is increased to 79.13%. 2 The detection rate was increased to 82.75%, the detection rate of the statistic Q was increased to 67.13%, the occurrence of the fault was detected at 190 sample points, the fault detection time was advanced by 170 sample points, and it had the highest fault detection rate compared with the CVDA method. Therefore, the WPCVDA method provided by the present invention significantly improves the detection performance of the TE process fault 21.

[0183] Table 3 and Table 4 respectively give the fault detection rates and fault detection times of the conventional CVDA method and the WPCVDA method of the present invention for TE process faults 11, 15, 19 and 21.

[0184] Table 3

[0185]

[0186] Table 4

[0187]

[0188] As can be seen from Tables 3 and 4, the WPCVDA method provided by the present invention has achieved the best monitoring results overall and has the best fault detection performance. In addition, the monitoring effect of faults 11, 15, 19 and 21 has been improved significantly.

[0189] Based on the above analysis, the PWCVDA method provided by the present invention uses the Wasserstein distance to mine the probability information contained in the chemical process data and considers the differences in the fault information of the features to perform online weighting on the features. Its fault detection effect is better than that of the traditional CVDA method.

[0190] The above embodiments are only used to illustrate the present invention and are not intended to limit the scope of protection of the present invention. Within the scope of the technical solution described in the present invention, various simple deformations and modifications made by technicians in the relevant technical field should be included in the above patent application scope.

Claims

1. A method for detecting slow-changing faults in chemical processes based on typical variable heterogeneity analysis, characterized in that: Contains the following steps: ( 1) Collect process data under normal operating conditions in the historical database as training data X0, and use the mean of training data X0 and standard deviation The training data X0 is standardized to obtain the standardized training data X; the training history data matrix P and the future data matrix F are constructed from the training data X; (ii) performing canonical variable dissimilarity analysis on the historical and future data matrices P and F to obtain a CVDA model, extracting the state vector s, residual vector e and dissimilarity feature d of the historical and future data matrices P and F using the CVDA model, and calculating the training WD feature s corresponding to the state vector s, residual vector e and dissimilarity feature d using the Wasserstein distance. WD 、e WD and d WD ; (III) Based on the training WD feature s WD 、e WD and d WD Calculate the corresponding statistic T 2 , statistics Q and statistics D, given the confidence level α, calculate the statistic T by the kernel density estimation method 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim ; (IV) Collect test data x new , using the mean of the training data X0 and standard deviation For the test data x new Perform standardization to obtain the standardized test data x t , from the test data x t Construct test history data vector p t and the future data vector f t ; (V) Using the CVDA model obtained in step (II), obtain the test history data vector p t and the future data vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t , and use Wasserstein distance to calculate the state vector s t , residual vector e t and heterogeneity characteristics d t Corresponding online WD features and (VI) Based on the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , obtain online weighted WD features and (VII) The online weighted WD feature and Calculate the corresponding new statistic T t 2 , statistic Q t and statistic D t , and according to the statistic T t 2 , statistic Q t and statistic D t Whether it exceeds the corresponding control limit, judge the test data x t Whether a fault occurs.

2. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 1, characterized in that: In step (a), firstly, the mean of the training data X0 is used and standard deviation The training data X0 is standardized by formula (1), and the expression of formula (1) is: After the training data X0 is normalized by the above formula (1), the standardized training data X can be obtained. The training data X=[x1, x2, ..., x i ,…,x n ] T ∈R n×m Medium i The corresponding training history data vector p i and the future data vector f i , where i represents the sampling time, x i represents the process variable collected at the i-th moment, n represents the number of samples, and m represents the number of variables. Formula (2) and formula (3) are expressed as follows: In the formula, h represents the time lag number, The historical data matrix P and the future data matrix F corresponding to the training historical data vector and the future data vector are calculated by formulas (4) and (5), which are expressed as follows: Where, M = n-2h + 1, p h+1 represents the historical data vector corresponding to the h+1th moment, f h+1 Represents the future data vector corresponding to the h+1th moment.

3. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 2, characterized in that: In the step (ii), the state vector s, residual vector e and dissimilarity feature d of the historical and future data matrices P and F are extracted for establishing the CVDA model, and the training WD feature s corresponding to the state vector s, residual vector e and dissimilarity feature d is calculated. WD 、e WD and d WD The specific steps are: First, the covariance matrix Σ corresponding to the historical and future data matrices P and F is calculated by formulas (6)-(8): pp and Σ ff and the cross-covariance matrix Σ pf , formulas (6)-(8) are expressed as: For the covariance matrix Σ pp and Σ ff and the cross-covariance matrix Σ pf Carry out the generalized eigenvalue decomposition shown in formula (9), the expression of formula (9) is: Where U represents the left singular matrix, V represents the right singular matrix, The eigenvalue The diagonal matrix formed by descending order; By solving formula (9), we can carry out formulas (10) and (11) to calculate the projection matrix and Formulas (10) and (11) are expressed as follows: The projection vector pair {a i ,b i } corresponds to the i-th correlation coefficient λ i , according to λ i The size of the value, the projection vector pair is divided into a typical variable subspace and a residual subspace, the first k eigenvalues ​​λ1>λ2>…>λ k Represents the main correlation coefficient corresponding to the typical variable subspace projection matrix A k =[a1,a2,…,a k ], the remaining Eigenvalues represents the residual subspace projection matrix corresponding to the weak correlation coefficient, and the value of k is determined by the average method; For the training history data vector p corresponding to the i-th sampling time in the training data X i and the future data vector f i , extract the corresponding state vector s through formulas (12)-(13) i , residual vector e i and heterogeneity characteristics d i , the expressions of formulas (12) - (13) are: Where U k represents the first k columns of the matrix U, and the dissimilarity feature measure is the future data vector f i The predicted typical variables and historical data vector p i The difference in canonical variates obtained; Calculate the state vector s = [s1, s2, ..., s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ; Use the sliding window to obtain the state vector s i,j The local mean and variance of the state vector s is calculated by formula (15): i,j Corresponding WD characteristics The expression of formula (15) is: In the formula, Represents the state vector s required using the sliding window i,j The local mean of Represents the state vector s required using the sliding window i,j The local variance of . Train WD features by computing Similarly, the residual vector e is calculated i and heterogeneity characteristics d i The corresponding training WD features and The residual vector e i,j Corresponding WD characteristics Dissimilarity feature i,j Corresponding WD characteristics 4. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 3, characterized in that: In step (iii), the training WD feature s is calculated WD 、e WD and d WD The corresponding statistic T 2 The specific steps of calculating statistics Q and statistics D are: Based on WD characteristics WD 、e WD and d WD , the statistics T are constructed by formulas (16)-(18) respectively 2 , statistics Q and statistics D, the expressions of formulas (16)-(18) are: In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d; Given a confidence level α, the statistic T is calculated using the kernel density estimation method 2 The corresponding control limits The control limit Q corresponding to the statistic Q lim and the control limit D corresponding to the statistic D lim .

5. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 4, characterized in that: In step (iv), the mean of the training data X0 is used and standard deviation The test data x is calculated by formula (19) new After normalization, the expression of formula (19) is: The test data x collected at the tth sampling time new After the above formula (19) is used for standardization, the standardized test data x can be obtained. t； The test data x t Corresponding test history data vector p t and the future data vector f t Calculated by formula (20) and (21), the expressions of formula (20) and (21) are: 。 6. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 5, characterized in that: In the step (V), the test history data vector p t and the future data vector f t The corresponding state vector s t , residual vector e t and heterogeneity characteristics d t Calculated by formula (22)-(24), the expressions of formula (22)-(24) are: Where U k To solve formula (9), we obtain the first k columns of the left singular matrix U, A k and B k To solve equations (10) and (11) to obtain the first k pairs of projection vectors; The state vector s = [s1, s2, ..., s j ,…s k ] is the jth state vector s j The corresponding mean μ j With variance σ j ; Use sliding window to obtain state vector s t,j The local mean and variance of the state vector s is calculated by formula (25) t,j Corresponding WD characteristics The expression of formula (25) is: In the formula, Represents the state vector s required using the sliding window t,j The local mean of Represents the state vector s required using the sliding window t,j The local variance of Testing WD characteristics by calculation Similarly, the residual vector e is calculated t and heterogeneity characteristics d t The corresponding training WD features and The residual vector e t,j Corresponding WD characteristics Dissimilarity feature t,j Corresponding WD characteristics 7. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 6, characterized in that: In the step (six), the online weighted WD feature and The specific steps are: First, according to the WD characteristics Relative to the training WD features s WD Changes in WD characteristics Assign weight w s,t , WD characteristics Relative to the training WD feature e WD Changes to the Online WD Features Assign weight w e,t , WD characteristics Relative to the training WD feature d WD Changes to the Online WD Features Assign weight w d,t , calculated by formula (26)-(28), the expressions of formula (26)-(28) are: In the formula, the weight w s,t =[w s,t (1),…,w s,t (j),…,w s,t (k)] T w in s,t (j) represents the jth online WD feature The weighting coefficient, weight w in e,t (j) represents the jth online WD feature The weighting coefficient, weight w d,t =[w d,t (1),…,w d,t (j),…,w d,t (k)] T w in d,t (j) represents the jth online WD feature The weighting coefficient of The weighted threshold and Generally choose to train WD features and 1 to 2 times the maximum value of the element; Get online weighted WD features and Calculated by formula (29)-(31), the expressions of formula (29)-(31) are: 。 8. A chemical process slow-changing fault detection method based on typical variable heterogeneity analysis as claimed in claim 7, characterized in that: In step (VII), the online weighted WD feature and Calculate the corresponding new statistic T t 2 , statistic Q t and statistic D t The specific steps are: Based on online weighted WD features and The statistics T are constructed by formulas (32)-(34) respectively 2 , statistics Q and statistics D, the expressions of formulas (32)-(34) are: In the formula, Λ d Represents the covariance matrix corresponding to the dissimilarity feature d; Determine the test data x t The steps to check whether a fault occurs are: Q t ≤Q lim And D t ≤D lim When the chemical process is in normal working condition, it is considered that the chemical process is in normal working condition; otherwise, it is considered that the chemical process is in failure.

Images