A non-Gaussian nonlinear process monitoring method based on independent slow feature analysis

Abstract

The present invention discloses a non-Gaussian nonlinear process monitoring method based on independent slow feature analysis, which belongs to the field of industrial process monitoring and fault diagnosis. The method includes the following steps: 1. Collecting observation data under historical normal working conditions, and using ISFA to extract hidden variables that are both independent and slowly changing from the observation data; 2. Establishing relevant statistical indicators, and using the kernel density estimation method to estimate the probability density function of the statistical indicators to determine the control limit; 3. Collecting online process data, calculating the statistical indicators of the online data and comparing them with the control limit, if the control limit is exceeded, it is judged that there is a fault in the process. The present invention can better extract hidden variables that can characterize the nonlinearity of the process in a non-Gaussian, nonlinear process in which the observation data is a nonlinear mixture of hidden variables, thereby obtaining more accurate detection results than traditional data-driven methods, thereby reducing the losses caused by faults and improving the safety of process operation.

Description

Technical Field

[0001] The invention belongs to the field of industrial process monitoring and fault diagnosis, and in particular relates to a non-Gaussian nonlinear process monitoring method based on independent slow features. Background Art

[0002] In industrial processes, since process variables are significantly correlated and information is highly redundant, the main changes in the process are usually dominated by a small number of latent variables, and the dimensions of these dominant latent variables are often much smaller than the dimensions of process variables. Therefore, as long as the dominant variables are monitored for abnormalities, it can be determined whether the industrial process has a fault. Some typical data-driven latent variable models such as principal component analysis (PCA) and independent component analysis (ICA) have received widespread attention and great development and progress in the field of multivariate statistical process monitoring (MSPM). However, these traditional methods have the following problems:

[0003] (1) Traditional PCA assumes that latent variables are statistically uncorrelated and obey Gaussian distribution. In actual industrial processes, such as chemical and biochemical processes, these assumptions are usually not true. These latent variables often obey non-Gaussian distribution. In this case, the accuracy of traditional PCA-based process monitoring methods is very low.

[0004] (2) Although ICA does not require data to follow a Gaussian distribution and further uses higher-order information (skewness and kurtosis) based on PCA, another assumption that ICA can restore latent variables from observed data is that the observed data is a linear combination of latent variables that are statistically independent of each other. In reality, due to the strong nonlinearity of process variables, the assumption that the observed data is a linear mixture of latent variables is difficult to meet, and the effectiveness of ICA for process monitoring will be greatly reduced.

[0005] As today's industrial processes become more and more complex and measurement equipment is constantly upgraded, the measurable process variables are more diverse. These variables often have the characteristics of strong correlation, nonlinearity, high dimension, and non-Gaussianity. Therefore, the dilemma faced by traditional methods is not uncommon. Summary of the invention

[0006] The purpose of the present invention is to provide a non-Gaussian nonlinear process monitoring method based on independent slow feature analysis in view of the deficiencies of the prior art.

[0007] The objective of the present invention is achieved through the following technical solutions:

[0008] A non-Gaussian nonlinear process monitoring method based on independent slow feature analysis comprises the following steps:

[0009] Step 2: Online monitoring: Independent slow features and residuals are extracted from the online data of the process through ISFA, and the statistical indicators of the online data are calculated using the extracted independent slow features and residuals. If the value of any calculated online data statistical indicator exceeds the control limit established in the offline training stage of step 1, it is judged that there is a fault in the process.

[0010] ]The offline training process described in step 1 is as follows:

[0011] 1.1) Assume that the observed data is a set of N-dimensional zero-mean input signals x(t) = [x1(t), ..., x N (t)] T , t is the set of all sampling times under the process running state, and T is the transpose of the matrix. First, the input data x(t) is nonlinearly expanded z(t)=h(x(t)), where z(t) is an L-dimensional signal with a mean of 0, and h is some nonlinear function, such as the binomial expansion h(x(t)) of x(t) has the following form h(x)=[x1,...,x N ,x1x1,x1x2,...,x N x N ] T -h0 T , h0 T is a constant vector whose value corresponds to the mean of each dimension, so that the mean of each dimension of h(x) is 0. The signal z(t) after nonlinear expansion is whitened (spherized) to obtain y(t) = Wz(t), so that the components y of the whitened signal y(t) are i The mean of (t) is 0, the variance is 1, i∈[1,L], representing the dimension of the signal y(t); W is the whitening matrix W=Λ -1 / 2 U T , T is the transpose of the matrix. Where Λ and U are the covariance matrix C of the data z(t) z =E(z(t)z T (t)) Perform singular value decomposition C z =UΛU T We get that T is the transpose of the matrix;

[0012] 1.2) Perform independent slow feature analysis on the whitened signal y(t) to obtain the L×L matrix Q, so that the first R components u1(t),...,u R(t) are as statistically independent as possible and change slowly over time. These R components are called independent slow features. The definitions of y(t), z(t), h and W are the same as in step 1.1. The value of R is much smaller than L. The remaining LR components u R+1 (t),...,u L (t) Compared with the first R components, which change faster and may not be independent of each other, the optimization goal of independent slow feature analysis is written as minimizing the objective function:

[0013]

[0014] where b ICA and b SFA Represent the weights of the independent component analysis (ICA) part and the slow feature analysis (SFA) part, both are set to 1, τ is the delay, increasing the delay can better describe the dynamics of the system, and setting multiple delay τ values ​​can have better robustness. ICA represents the set of desirable time delays τ for the objective function of ICA, represents the weight of the ICA objective function for different values ​​of the delay τ. The ICA part of the ISFA objective function ensures the independence of the signal, and the SFA part ensures the slowness of the signal changing over time. The goal of ISFA is to find an orthogonal matrix Q such that u(t) = Qy(t) is the minimum for the ISFA objective function. The ISFA objective function based on Givens rotation is written as

[0015] The method for obtaining the orthogonal matrix Q is continuous Givens rotation. Givens rotation refers to the rotation transformation of selecting two parameters μ and ν to rotate the vector counterclockwise around the origin in the (μ, ν) plane by φ radians. The rotation matrix has the following form:

[0016]

[0017] where u′ is some instantaneous signal during the rotation optimization step;

[0018] 1.3) Let the unmixing matrix B = QW, where the definitions of W and Q are the same as in steps 1.1 and 1.2. According to the Euclidean norm of each row of matrix B, select the largest d rows and arrange them in descending order of Euclidean norm to form matrix B. d , the remaining rows are arranged in descending order of Euclidean norm to form matrix B e ; Calculate the subspace u composed of independent slow features d (t) = B d h(x(t)) and residual subspace u e (t) = B e h(x(t)) and its statistical index I2 (t) = u d (t) T u d (t), I e 2 (t) = u e (t) T u e (t);

[0019] 1.4) Use kernel density estimation method to estimate statistical indicators I 2 (t) and I e 2 The probability density function of (t) is used to calculate the control limits according to the confidence level.

[0020] The online detection process described in step 2 is as follows:

[0021] 2.1) Select the same nonlinear function h as in the offline training phase in step 1 to deal with the test data x new Do nonlinear expansion and get z new (t) = h(x new (t));

[0022] 2.2) Calculate its independent slow feature u dnew =B d z new and the residual u enew =B e z new ;

[0023] 2.3) Calculate statistical indicators

[0024] 2.4) If and The data to be tested in this group belongs to the data under normal operating conditions, and no fault occurs during the process. and If any of them exceeds its corresponding control limit, the process is judged to have failed.

[0025] Beneficial effects of the present invention:

[0026] The present invention can better extract latent variables that can characterize the nonlinearity of the process in a non-Gaussian, nonlinear process in which the observed data is a nonlinear mixture of latent variables, thereby obtaining more accurate detection results than traditional data-driven methods, thereby reducing the losses caused by faults and improving the safety of process operation. BRIEF DESCRIPTION OF THE DRAWINGS

[0027] Figure 1 It is a flow chart of the present invention.

[0028] Figure 2 This is an overview diagram of the present invention

[0029] Figure 3 TE process control structure schematic diagram

[0030] Figure 4 It is the result of the ICA-based detection of Tennessee-Eastman process fault 11.

[0031] Figure 5 The results of dICA based fault 11 detection on the Tennessee-Eastman process.

[0032] Figure 6 It is based on the ISFA's detection results of Tennessee-Eastman process fault 11. DETAILED DESCRIPTION

[0033] The present invention is further described below in conjunction with the accompanying drawings and embodiments.

[0034] like Figure 1 , Figure 2 As shown, the steps of a non-Gaussian nonlinear process monitoring method based on independent slow feature analysis of the present invention are as follows:

[0035] Step 2: Online monitoring: Independent slow features and residuals are extracted from the online data of the process through ISFA, and the statistical indicators of the online data are calculated using the extracted independent slow features and residuals. If the value of any calculated online data statistical indicator exceeds the control limit established in the offline training stage of step 1, it is judged that there is a fault in the process.

[0036] The offline training process described in step 1 is as follows:

[0037] 1.1) Assume that the observed data is a set of N-dimensional zero-mean input signals x(t) = [x1(t), ..., x N (t)] T , t is the set of all sampling times under the process running state, and T is the transpose of the matrix. First, the input data x(t) is nonlinearly expanded z(t)=h(x(t)), where z(t) is an L-dimensional signal with a mean of 0, and h is some nonlinear function, such as the binomial expansion h(x(t)) of x(t) has the following form h(x)=[x1,...,x N ,x1x1,x1x2,...,x N x N ] T -h0 T , h0 Tis a constant vector whose value corresponds to the mean of each dimension, so that the mean of each dimension of h(x) is 0. The signal z(t) after nonlinear expansion is whitened (spherized) to obtain y(t) = Wz(t), so that the components y of the whitened signal y(t) are i The mean of (t) is 0, the variance is 1, i∈[1,L], representing the dimension of the signal y(t); W is the whitening matrix W=Λ -1 / 2 U T , T is the transpose of the matrix. Where Λ and U are the covariance matrix C of the data z(t) z =E(z(t)z T (t)) Perform singular value decomposition C z =UΛU T We get that T is the transpose of the matrix;

[0038] 1.2) Perform independent slow feature analysis on the whitened signal y(t) to obtain the L×L matrix Q, so that the first R components u1(t),...,u R (t) are as statistically independent as possible and change slowly over time. These R components are called independent slow features. The definitions of y(t), z(t), h and W are the same as in step 1.1. The value of R is much smaller than L. The remaining LR components u R+1 (t),...,u L (t) Compared with the first R components, which change faster and may not be independent of each other, the optimization goal of independent slow feature analysis is written as minimizing the objective function:

[0039]

[0040] where b ICA and b SFA Represent the weights of the independent component analysis (ICA) part and the slow feature analysis (SFA) part, both are set to 1, τ is the delay, increasing the delay can better describe the dynamics of the system, and setting multiple delay τ values ​​can have better robustness. ICA represents the set of desirable time delays τ for the objective function of ICA, represents the weight of the ICA objective function for different values ​​of the delay τ. The ICA part of the ISFA objective function ensures the independence of the signal, and the SFA part ensures the slowness of the signal changing over time. The goal of ISFA is to find an orthogonal matrix Q such that u(t) = Qy(t) is the minimum for the ISFA objective function. The ISFA objective function based on Givens rotation is written as

[0041]

[0042] The method for obtaining the orthogonal matrix Q is continuous Givens rotation. Givens rotation refers to the rotation transformation of selecting two parameters μ and ν to rotate the vector counterclockwise around the origin in the (μ, ν) plane by φ radians. The rotation matrix has the following form:

[0043]

[0044] where u′ is some instantaneous signal during the rotation optimization step;

[0045] 1.3) Let the unmixing matrix B = QW, where the definitions of W and Q are the same as in steps 1.1 and 1.2. According to the Euclidean norm of each row of matrix B, select the largest d rows and arrange them in descending order of Euclidean norm to form matrix B. d , the remaining rows are arranged in descending order of Euclidean norm to form matrix B e ; Calculate the subspace u composed of independent slow features d (t) = B d h(x(t)) and residual subspace u e (t) = B e h(x(t)) and its statistical index I 2 (t) = u d (t) T u d (t), I e 2 (t) = u e (t) T u e (t);

[0046] 1.4) Use kernel density estimation method to estimate statistical indicators I 2 (t) and I e 2 The probability density function of (t) is used to calculate the control limits according to the confidence level.

[0047] The online detection process described in step 2 is as follows:

[0048] 2.1) Select the same nonlinear function h as in the offline training phase in step 1 to deal with the test data x new Do nonlinear expansion and get z new (t) = h(x new (t));

[0049] 2.2) Calculate its independent slow feature u dnew =B d z new and the residual u enew =B e z new ;

[0050] 2.3) Calculate statistical indicators

[0051] 2.4) If and The data to be tested in this group belongs to the data under normal operating conditions, and no fault occurs during the process. and If any of them exceeds its corresponding control limit, the process is judged to have failed.

[0052] Application Examples

[0053] The effectiveness of the method of the present invention is illustrated below by combining a simple multivariable dynamic process and the Tennessee Eastman (TE) process.

[0054] Consider the following numerical model of a multivariable dynamic process, which is a modified version of a numerical model widely used in the field of process monitoring. Next, the implementation steps of the present invention are described in detail in conjunction with this specific process:

[0055]

[0056] y(t)=z(t)+v(t)

[0057] Where u is an autocorrelation input:

[0058]

[0059] in

[0060]

[0061] It is a nonlinear combination of latent variables sig1 and sig2, which has extremely extreme nonlinearity.

[0062]

[0063]

[0064] The output y is equal to z plus a random noise v, where the noise v satisfies the mean of 0 and the variance of 0.1. The input signal u and the output signal y can be measured, while the remaining z and w cannot be measured. The measurable u and y form the input variable x(t) = [y(t) T u(t) T ] TThe training set consists of 200 samples of x(t) under normal working conditions. After data preprocessing, the mean of x(t) is scaled to 0. In this multivariate dynamic process, we mix the latent variables in an extreme nonlinear form, and artificially add faults to the process by changing the latent variables. PCA, ICA and ISFA are used to monitor the process.

[0065] Traditional PCA not only assumes that the sequences between the current moment and the historical moment for the same latent variable are independent, but also assumes that different latent variables are independent of each other. However, such an assumption is difficult to hold in real industrial processes, because the data in industrial processes not only have cross-correlation characteristics but also autocorrelation characteristics. The assumption of traditional ICA to extract latent variables from observed data is that the observed data is composed of a linear combination of latent variables. Such an assumption is often difficult to meet in real complex industrial processes. In this multivariate dynamic process, neither the assumptions of traditional PCA nor the assumptions of traditional ICA are met. The PCA used here is the traditional PCA that everyone is familiar with, and the statistical indicator used is Hotelling's T. 2 and SPE, the number of principal components selected is 4. There are many ICA methods for extracting statistically independent latent variables. In this experiment, we use the most mature and widely used FastICA method. The number of independent components selected is 2, and the statistical index used is I 2 and In order to make a fairer comparison and reflect the superiority of the independent slow features extracted by ISFA over the independent components extracted by ICA, when using ISFA, the number of independent slow features is also selected as 2. The delay in the independent slow feature optimization step is 2, R = 5, and R = 5 is selected to keep the same dimension as the original observation data. The confidence level of the three methods is selected as 95%.

[0066] The number of samples in the fault test set is 500, which is used to detect and diagnose artificially set faults as follows:

[0067] Fault 1: The hidden variable sig1 begins to undergo a step change of 0.4 at the 50th moment. This is a relatively minor fault, and the fault will be hidden in the noise information.

[0068] Fault 2: Add 0.05(t-50) to the value of the latent variable sig1 at time t, and t is 50-149. This makes the latent variable sig1 increase linearly from the 50th moment to the 149th moment.

[0069] Its detection rates for two artificially set faults are shown in Table 1. The highest detection rate for the same fault is marked in bold. When the number of latent variables is small, ISFA can still obtain higher accuracy than the other two methods in two different types of faults, indicating that in this multivariate dynamic process, ISFA has a stronger ability to extract latent variables, and the extracted latent variables can more accurately monitor the changes of the entire process. In summary, in this multivariate dynamic process where the observed data is a nonlinear mixture of latent variables, ISFA is significantly better than PCA and ICA in fault detection for two different fault types. The experimental results strongly support the theoretical analysis.

[0070] Table 1 Fault detection rate of each monitoring method for multivariable dynamic process (%)

[0071]

[0072] The control structure schematic diagram of the TE process is as follows: Figure 3 As shown. The TE process data consists of a total of 53 observed variables, including 12 control variables XMV (1-12), 22 process measurements XMEAS (1-22) and 19 component measurements XMEAS (23-41). The sampling interval of the TE process is 3 minutes. The training set samples are obtained under 25 hours of simulation operation, with a total of 500 observation samples. The test set samples are obtained under 48 hours of simulation operation. The fault is introduced in the 8th hour. The test set collects a total of 960 observation samples. The first 160 observations are fault-free data under normal conditions. The TE process has a total of 21 different faults. Since the 19 component measurement variables are difficult to be detected in real time in the actual chemical process, the control variable XMV (12) stirring speed is not actually controlled. This experiment selects 11 control variables XMV (1-11) and 22 process variables XMEAS (23-41), a total of 33 variables, as input variables.

[0073] Next, the implementation steps of the present invention are described in detail in conjunction with the specific process:

[0074] The training data are all preprocessed by removing the mean. In this experiment, the number of independent components selected by ICA is 9, the number of independent components selected by dynamic ICA is 22, and the delay is T=2. In order to provide a more fair and reasonable comparison result, the proposed ISFA selects 22 independent slow features, the delay T=2 is the same as dynamic ICA, and R=33 maintains the same dimension as the original input variable. It is worth noting that after nonlinear expansion, it is easy to have data dimension explosion and highly redundant information. While increasing the amount of calculation, highly redundant information will also have an adverse effect on the final result. The purpose of nonlinear expansion is to solve nonlinear problems linearly in a higher spatial dimension, but for industrial data such as TE process, the data itself is nonlinear. Strong correlation and other characteristics, so in this part, ISFA omits the data preprocessing step of nonlinear expansion. The confidence level of each statistical indicator of the three methods is selected as 99%.

[0075] The detection rates of 21 faults are shown in Table 1. The highest detection rate for the same fault is marked in bold. From the final results in the table, except for faults 10 and 16, the fault detection accuracy of the proposed method is higher than that of the other two comparison methods. The detection rate of ISFA for fault 10 is higher than that of dynamic ICA, which is only 2.63% lower than ICA. The detection rate of ISFA for fault 16 is higher than that of dynamic ICA, which is only 0.13% lower than ICA. In general, the performance of ISFA in the field of process monitoring is better than that of traditional ICA and dynamic ICA. When the same number of latent variables are selected to construct the statistical indicators, the fault detection rate of ISFA is higher than that of dynamic ICA, which further proves that ISFA has a stronger ability to extract latent variables than the other two methods, and the extracted latent variables can more accurately monitor the changes of the entire process. Figure 4 , Figure 5 and Figure 6 The following is the monitoring effect diagram of fault 11 in the TE process using ICA, dICA and ISFA respectively. The first 160 curves marked with light solid lines represent normal sample data under normal working conditions, and the curves marked with dark solid lines after 800 represent fault sample data. The dotted line represents the control limit. It can be clearly seen from the figure that ISFA is superior to the other two methods.

[0076] Table 2 Fault detection rate of each monitoring method in TE process (%)

[0077]

[0078] The above embodiments are used to illustrate the present invention rather than to limit the present invention. Any modification and change made to the present invention within the spirit of the present invention and the protection scope of the claims shall fall within the protection scope of the present invention.

Claims

1. A non-Gaussian nonlinear process monitoring method based on Independent Slow Feature Analysis (ISFA), It is characterized in that The following steps are involved: Step 1: Offline training, collect data of normal operation of the process, ISFA extracts hidden variables and residuals that are both statistically independent and change slowly over time from the observed data, and uses the subspace composed of the extracted hidden variables and the residual subspace to establish corresponding statistical indicators, and uses the kernel density to estimate the probability density function of the statistical indicators, and calculates the control limits according to the confidence level. The hidden variables extracted by ISFA are called independent slow features; Step 2: Online monitoring, using ISFA to extract independent slow features and residuals from the data of the process online operation, and using the extracted independent slow features and residuals to calculate the statistical indicators of the online data respectively. If the value of any calculated online data statistical indicator exceeds the control limit established in the offline training phase of step 1, it is judged that there is a fault in the process; The offline training process described in step 1 is as follows: 1.1) Assume that the observed data is a set of N-dimensional zero-mean input signals x(t) = [x 1 (t),…,x N (t)] T , t is the set of all sampling times under the process running state, T is the transpose of the matrix, firstly, the input data x(t) is nonlinearly expanded z(t)=h(x(t)), z(t) is an L-dimensional signal with a mean of 0, h is a nonlinear function, if the binomial expansion h(x(t)) of x(t) has the following form h(x)=[x 1 ,…,x N ,x 1 x 1 ,x 1 x 2 ,…,x N x N ] T -h 0 T ,h 0 T is a constant vector, and each component of the vector corresponds to the mean of each dimension of the signal, so that the mean of h(x) is 0; the signal z(t) after nonlinear expansion is whitened to obtain y(t) = Wz(t), so that the components y of the whitened signal y(t) i The mean of (t) is 0, the variance is 1, i∈[1,L], representing the dimension of the signal y(t); W is the whitening matrix W=Λ -1 / 2 U T , where Λ and U are the covariance matrices C of z(t) z =E(z(t)z T (t)) Row singular value decomposition C z =UΛU T get; 1.2) Perform independent slow feature analysis on the whitened signal y(t) to obtain the L×L matrix Q, so that u(t)=Qy(t)=QWz(t)=QWh(x(t)) has the first R components u 1 (t),…,u R (t) are as statistically independent as possible and change slowly over time. These R components are called independent slow features. The definitions of y(t), z(t), h and W are the same as in step 1.1). The value of R is much smaller than L. The remaining LR components u R+1 (t),…,u L (t) Compared with the first R components, which change faster and may not be independent of each other, the optimization goal of independent slow feature analysis is written as minimizing the objective function: where b ICA and b SFA They represent the weights of the independent component analysis (ICA) part and the slow feature analysis (SFA) part, both of which are set to 1. τ is the delay. Increasing the delay can better describe the dynamics of the system. Setting multiple delay τ values ​​can have better robustness. T ICA represents the set of desirable time delays τ for the objective function of ICA, represents the weight of the ICA objective function for different values ​​of the delay τ. The ICA part of the ISFA objective function ensures the independence of the signal, and the SFA part ensures the slowness of the signal changing over time. The goal of ISFA is to find an orthogonal matrix Q such that u(t) = Qy(t) is the minimum for the ISFA objective function. The ISFA objective function based on Givens rotation is written as where u′ is some instantaneous signal during the rotation optimization step; The method for obtaining the orthogonal matrix Q is continuous Givens rotation. Givens rotation refers to the rotation transformation of selecting two parameters μ and ν to rotate the vector counterclockwise around the origin in the (μ, ν) plane by φ radians. The rotation matrix has the following form: 1.3) Let the unmixing matrix B = QW, where the definitions of W and Q are the same as in steps 1.1 and 1.

2. According to the Euclidean norm of each row of matrix B, the rows in B are arranged in descending order according to the size of the Euclidean norm to form a new matrix, and the first d rows are selected as matrix B. d , and the remaining rows form the matrix B e ; Calculate the subspace u composed of independent slow features d (t) = B d h(x(t)) and residual subspace u e (t) = B e h(x(t)) and its statistical index I 2 (t) = u d (t) T u d (t) and I e 2 (t) = u e (t) T u e (t); 1.4) Use kernel density estimation method to estimate statistical indicators I 2 (t) and I e 2 The probability density function of (t) is used to calculate the control limits according to the confidence level. and 2. The method according to claim 1, It is characterized in that The online monitoring process described in step 2 is as follows: 2.1) Select the same nonlinear function h as in the offline training phase in step 1 for the online data x new Do nonlinear expansion and get z new (t) = h(x new (t)); 2.2) Calculate its independent slow feature u dnew =B d z new and the residual u enew =B e z new ; 2.3) Calculate statistical indicators 2.4) If and This group of online data belongs to the data under normal operating conditions, and no fault occurs during the process. and If any of them exceeds its corresponding control limit, the process is judged to have failed.

Images