Distributed early fault detection method based on first-order perturbation theory
By adopting distributed methods and first-order perturbation theory in industrial processes, combined with CVA and local outlier factor algorithms, the problem of early fault detection is solved, and efficient detection and early warning of early faults in complex industrial processes is achieved.
Patent Information
- Application Number
- CN202510306483.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-03-14
- Publication Date
- 2025-06-27
AI Technical Summary
In complex industrial processes, early failures are difficult to effectively detect through a single centralized global model due to their weak characteristics, concealment and randomness, and the direct dimensionality reduction method will lose fault information and fail to effectively pay attention to local information.
A distributed early fault detection method based on first-order perturbation theory is adopted. Through mutual information, the detection perspective is shifted from global to local, dynamic features are extracted using CVA, and recursive feature updates are performed through first-order perturbation, local outlier factor algorithm is used to pay attention to local information, and finally the monitoring results of each subspace are fused through Bayesian inference.
It effectively reduces the impact of system coupling, enhances the sensitivity to early failures, improves detection performance, and can promptly and effectively detect and early warning of early failures in complex industrial processes.
Smart Images

Figure FT_1 
Figure FT_2 
Figure FT_3
Abstract
Description
Technical Field
[0001] The present invention relates to a data-driven fault detection method, and particularly to a distributed early fault detection method based on the first-order perturbation theory. Background Art
[0002] The complexity of modern industrial processes is increasing day by day, showing the characteristics of large scale, multiple operating units, and strong correlation among various systems. With the development of sensor technology, a large amount of process data has been collected and recorded, and multivariate statistical analysis based on this has attracted the attention of researchers. Common algorithms among them include principal component analysis (PCA), independent component analysis (ICA), canonical correlation analysis (CCA), canonical variable analysis (CVA), etc. However, when a fault occurs in a complex whole-plant process, it may only affect a few variables or operating units. Establishing a single centralized fault detection model to monitor all variables from different units simultaneously is very likely to cause the fault information to be submerged in the large-scale process data. For many serious faults, their initial stages can be regarded as early faults. The characteristics of early faults are weak features, evolving at a low rate and frequency, and being concealed and random. Due to the small fault amplitude, the change in process data is not obvious. Therefore, compared with conventional faults, early faults are more difficult to detect in the large-scale data fluctuations of modern industrial processes.
[0003] For early faults in the whole-plant process, using a single centralized global model often results in poor detection effects. And using the method of directly reducing process variables for dimensionality reduction, although it can reduce the masking of fault information by redundant variables to a certain extent, it will also cause losses to the fault information. In addition, this method does not pay attention to the local information that is crucial for early fault detection. Therefore, in order to reduce the complexity of the model in processing industrial processes, shifting the detection perspective from global to local and dividing the set of large-scale process variables into multiple subspaces by a distributed method is effective; in addition, it is also necessary to enhance the weak features of early faults, highlight the deviation degree of abnormal changes, and improve the detection performance of early faults. Based on the above discussion, in order to effectively detect early faults in complex industrial processes, the present invention proposes a distributed early fault detection method based on the first-order perturbation theory on the basis of the distributed method. Summary of the Invention
[0004] The main technical problems to be solved by the present invention are as follows: First, it is the subspace division problem of the distributed method. According to the mutual information, industrial process variables are divided into different subspaces according to the strength of correlation, shifting the detection perspective from the global to the local and weakening the influence of the system coupling effect. Second, it is the problem of establishing a subspace detection model. In each subspace, first use CVA to extract dynamic features, and then use first-order perturbation to recursively update eigenvalues and eigenvectors, enhancing the sensitivity of the model to early faults, and introducing the local outlier factor algorithm to further focus on the local information of the data to obtain the monitoring results of each subspace. Finally, through Bayesian inference fusion, the monitoring results of each subspace are fused to obtain a comprehensive detection result.
[0005] The technical solution adopted by the present invention to solve the above problems is: A distributed early fault detection method based on the first-order perturbation theory, comprising the following steps:
[0006] Step (1): Industrial process data acquisition and preprocessing link, the specific process is as follows:
[0007] Step (1.1): Collect sample data under normal working conditions in the industrial process as the training data set, denoted as
[0008] Step (1.2): Use the z-score normalization method to preprocess the collected sensor data to obtain the normalized data set X = [x1, x2,..., x m ∈ R n×m , where n is the number of data set samples, m is the number of process variables, and R is the set of real numbers. The specific calculation formula is as follows:
[0009]
[0010] Among them, x i is the data of a single sensor variable in the normalized data set, x i ∈ R n×1 , u is the mean of x i , and σ is the mean of x i .
[0011] Step (2): Use mutual information (MI) to process the linear and non-linear relationships between different variables. According to the correlation judgment between variables, use a data-driven method to divide subspaces, shifting the perspective from the global to the local and highlighting the sensitive information helpful for early fault detection. The specific implementation process is as follows:
[0012] Step (2.1): Calculate the variable x i ∈ R n×1(i = 1, 2, …, m) and the mutual information (MI) between all the remaining variables, and summarize these values in matrix form, which is called the mutual information matrix. The definition of the mutual information matrix is as follows:
[0013]
[0014] Among them, I(x i , x j ) represents the mutual information value between two variables, and p(x i , x j ) represents the joint distribution probability between variables x i and x j ;
[0015] Step (2.2): Use the Spectral Clustering algorithm to construct a graph using the similarity matrix (MI matrix) of the data, and achieve clustering through the eigen-decomposition of the Laplacian matrix of the graph. The detailed process is as follows:
[0016] Regard the similarity matrix W as the adjacency matrix of the graph, where W ij represents the similarity between variables x i and x j . Calculate the degree matrix D based on the similarity matrix:
[0017]
[0018] Among them, the non-diagonal elements are 0, and then calculate the normalized Laplacian matrix (symmetric form) L:
[0019] L = I - D -1 / 2 WD -1 / 2 (5)
[0020] Perform eigen-decomposition on the Laplacian matrix L to obtain the eigenvalues λ1, λ2, …, λ n and the corresponding eigenvectors v1, v2, …, v n . Arrange the first k eigenvectors in descending order to form an n×k eigenvector matrix U:
[0021] U = [v1, v2, …, v k (6)
[0022] Normalize each row of the matrix U to obtain the matrix X:
[0023]
[0024] After that, cluster the row vectors of T to obtain the partitioning result, and construct a subspace based on this partitioning process variable;
[0025] Step (3): The industrial process is filled with a large amount of process data related to time series. Reasonably utilizing this characteristic helps to improve the performance of early fault detection. In the process of constructing a subspace early fault detection model, Canonical Variate Analysis (CVA) can determine the linear combination between past variables and future variables, so it can be used to extract the dynamic characteristics of variables. The specific implementation process is as follows:
[0026] Maximize the correlation coefficient between the past variable y p,r and the future variable y f,r to extract the dynamic characteristics. r represents the current sampling time, p represents the past sampling points, and f represents the future sampling points. The past variable y p,r and the future variable y f,r are described as follows:
[0027]
[0028] The past matrix Y p and the future matrix Y f are described as follows:
[0029] Y p =[y p,p+1 y p,p+2 … y p,p+N ∈R mp×N (10)
[0030] Y f =[y f,f+1 y f,f+2 … y f,f+N ∈R mf×N (11)
[0031] where m represents the number of process variables, p represents the window length, and N is the dimension of the Hankel matrix to be set;
[0032] Then, the sample covariance and cross-covariance of past and future observations can be estimated as:
[0033]
[0034] CVA seeks to obtain the best linear combination that maximizes the correlation coefficient, which can be achieved by performing a singular value decomposition on the scaled Hankel matrix H:
[0035]
[0036] where U and V contain the left and right singular column vectors of H, Σ is a positive semi-definite diagonal matrix, and its diagonal elements are the singular values of H, that is, eigenvalues. The canonical variable s r is calculated as follows:
[0037]
[0038] Regular variable s r It consists of two parts, where the larger singular value represents the principal component space, i.e., the eigenvector, and the smaller singular value represents the residual space. At this point, the feature extraction phase of CVA is completed;
[0039] Step (4): First-order perturbation theory (FOP) is used to study the behavior changes of the system when it is subjected to small disturbances. In the process of feature extraction, FOP theory can identify and highlight the parts that have a greater impact on the industrial process to highlight the main changes in process variables, thereby more effectively analyzing abnormal deviations caused by early faults. Therefore, FOP is used to recursively update features and construct new features that are sensitive to early faults. The specific implementation process is as follows:
[0040] The eigenvalues of the typical variables of CVA are updated as follows:
[0041] λ k,i =(1-∈)λ k-1,i +∈|f| 2 (18)
[0042]
[0043] Among them, ∈ is a positive number close to zero, k is the sampling time, λ i is the diagonal element of Σ in formula (15), i.e., the eigenvalue, s k-1,i is the i-th eigenvector, corresponding to the eigenvalue, belonging to the canonical variable s in formula (16) r The principal component space corresponding to the larger singular value is the feature space; then the feature vector is recursively updated:
[0044]
[0045] Step (5): After obtaining the recursive feature space S that is more sensitive to early faults, the subspace detection index is constructed through the local outlier factor (LOF). LOF is used to find outliers in a multidimensional data set. The size of this value indicates the degree of sample outliers. Since the samples in the training data set are normal, the feature space sample values are similar at each sampling time. Considering that the feature space sample values of the fault detection samples can be used as outliers of normal samples, LOF can be used to construct the detection metric. The specific implementation process is as follows:
[0046]
[0047] Among them, s is a sample in S, s f is the fth nearest neighbor sample point, and is defined as follows:
[0048] Rd(s, s f ) = max{d(s, s f ), Fd(s f )} (23)
[0049]
[0050] where Fd(s f ) represents the distance from s f to its farthest point within the requirement of the number of nearest neighbors;
[0051] Step (6): Use Bayesian inference to perform statistical combination on the statistics of all subspaces to obtain the final comprehensive joint statistic. The specific implementation process is as follows:
[0052] Step (6.1): Use Bayesian inference to perform statistical combination on the statistical data of all subspaces. The failure probability calculation formula for a single subspace x b is;
[0053]
[0054] P(X b ) = P(X b |N)P(N) + P(X b |F)P(F) (26)
[0055] where the conditional probabilities P(X b |N) and P(X b |F) are expressed as
[0056]
[0057] In the above formula, N and F represent the normal state and the failure state respectively; P(N) and P(F) represent the prior probabilities under normal and failure conditions respectively; when the confidence level is determined to be ε, P(N) = ε, P(F) = 1 - ε. J b,new is the statistical result of the b-th subspace in the online detection dataset, and J b,lim is the corresponding control threshold;
[0058] Step (6.2): Combine the detection results of all different subspaces through Bayesian inference and calculate the sum of the joint statistics:
[0059]
[0060] Set a confidence level of 0.99 and construct the statistical confidence limit BIC L in the offline stage to complete the construction of the distributed early fault detection model;
[0061] The above steps (1) to (6) are the offline modeling stage of the method of the present invention, and the following steps (7) to (12) are the implementation process of the online analysis of the method of the present invention.
[0062] Step (7): Collect sample data at a new sampling moment Process the obtained data using the mean and standard deviation obtained in step (1.2) to obtain the standardized test sample data XT = R 1×m ;
[0063] Step (8): Divide the process variable XT according to the data-driven method in step (2) to obtain the subspace division result;
[0064] Step (9): For the current moment data of the process variable in each subspace and the sample data after standardization at the previous 2p - 1 moments, construct an early fault detection model in each subspace according to the formula in the manner of step (3), and use CVA for dynamic feature extraction to obtain the feature vector space;
[0065] Step (10): According to the method in step (4), based on the feature vector space obtained in step (9), combined with FOP recursive update, construct new features sensitive to early faults;
[0066] Step (11): In each subspace in the manner of step (5), use the LOF method to construct a detection metric combined with the sensitive new features of the online data to obtain the LOF value of the detection result for each subspace;
[0067] Step (12): Obtain the Bayesian inference fusion result BIC of the test data set in the manner of step (6), and compare it one by one with the statistical confidence limit BIC in the offline stage L If it exceeds the confidence limit, an alarm is given and it is considered that a fault has occurred; if it does not exceed the limit, it is considered that the system is operating normally.
[0068] Compared with the traditional method, the advantages of the method of the present invention are as follows:
[0069] First, the method of the present invention determines the correlation between process variables through MI, realizes an accurate subspace division result, shifts the fault detection perspective from global to local, and highlights the local characteristics of early faults; second, the method of the present invention extracts the dynamic characteristics of data by using the CVA algorithm, considering the time series correlation of process data; third, the method of the present invention uses the FOP theory to recursively update the feature vector space, realizes the enhancement of fault information, and constructs a sensitive statistic combined with the LOF method, which can detect early faults more efficiently and timely with high performance. It can be said that the method of the present invention is a more excellent early fault detection method. Description of the Drawings
[0070] Figure 1 is the flow chart of the method of the present invention;
[0071] Figure 2 is the industrial flow chart of TEP;
[0072] Figure 3 is the simulation experiment diagram of the method of the present invention; (a) Fault 8WMA-PCA and the simulation results of the present method (b) Fault 13WMA-PCA and the simulation results of the present method. Detailed implementation manners
[0073] The method of the present invention will be described in detail below in conjunction with the accompanying drawings and specific implementation cases.
[0074] The present invention discloses a distributed early fault detection method based on the first-order perturbation theory. As Figure 1 shown, the specific flow chart of the method is presented. The following will illustrate the specific implementation process of the present invention and its superiority over the conventional multivariate statistical quality-related fault detection method in combination with a simulation platform case.
[0075] The Tennessee Eastman (TE) process is a reliable chemical plant simulation program based on actual chemical reaction processes and is widely used as a benchmark. Its specific industrial process is as Figure 2 shown, including five main units: a reactor, a condenser, a compressor, a separator, and a stripper. The TE process includes 41 process variables, 11 operating variables, and 21 preset faults of different fault types. The method of the present invention will use 22 continuous variables and 11 manipulated variables for simulation experiments.
[0076] Table 1: 22 continuous variables in the TE process
[0077]
[0078] Table 2: 11 manipulated variables in the TE process
[0079]
[0080] Table 3: Preset faults in the TE process for comparative experiments
[0081]
[0082] The test set samples with faults are obtained under a 48-hour running simulation. The faults are introduced at 8 hours, and a total of 960 observations are collected. Among them, the first 160 observations are normal data. First, the 960 samples collected are used for offline training to establish a fault detection model, including the following steps:
[0083] Step (1): Collect samples under normal conditions to form a training data set Normalize the data using the z-score method to obtain data X ∈ R 960×33 ;
[0084] Step (2): Use MI to process the correlation between different variables. Based on the correlation judgment between variables, calculate the MI between variable x i ∈ R n×1 (i = 1, 2, …, m) and all other variables, and summarize these values in matrix form, called the mutual information matrix:
[0085]
[0086] Among them, I(x i , x j ) represents the mutual information value between two variables, and p(x i , x j ) represents the joint distribution probability between variable x i and x j ;
[0087] Use the spectral clustering algorithm to construct a graph using the MI matrix of the data, and achieve clustering through the eigenvalue decomposition of the Laplacian matrix of the graph. Regard the similarity matrix W as the adjacency matrix of the graph, where W ij represents the similarity between variable x i and x j . Calculate the degree matrix D based on the similarity matrix:
[0088]
[0089] Among them, the off-diagonal elements are 0, and then calculate the normalized Laplacian matrix (symmetric form) L:
[0090] L = I - D -1 / 2 WD -1 / 2 (33)
[0091] Perform eigenvalue decomposition on the Laplacian matrix L to obtain eigenvalues λ1, λ2, …, λ n and corresponding eigenvectors v1, v2, …, v n , and arrange the first k eigenvectors in descending order to form an n × k eigenvector matrix U:
[0092] U = [v1, v2, …, v k (34)
[0093] Normalize each row of the matrix U to obtain the matrix X:
[0094]
[0095] Then cluster the row vectors of T to get the class division result, which is used to divide the process variables and construct the subspace [X1,X2,…,X B ];
[0096] Step (3): In the constructed subspace, CVA is used to extract the dynamic features of the variables. The specific implementation process is as follows:
[0097] The past matrix Y p and the future matrix Y f The description is as follows:
[0098] Y p =[y p,p+1 y p,p+2 … y p,p+N ]∈R mp×N (36)
[0099] Y f =[y f,f+1 y f,f+2 … y f,f+N ]∈R mf×N (37)
[0100] Among them, p and f represent the window length, which is set to 5;
[0101] Then, the sample covariance and cross-covariance of past and future observations are estimated as:
[0102]
[0103] Perform singular value decomposition on the scaled Hankel matrix h:
[0104]
[0105] where u and v contain the left and right singular column vectors of h, and Σ is a semi-positive definite diagonal matrix whose diagonal elements are the singular values, i.e., eigenvalues, of h. The canonical variable s r The calculation is as follows:
[0106]
[0107] Regular variable s r It consists of two parts, where the larger singular value represents the principal component space, i.e., the eigenvector, and the smaller singular value represents the residual space. At this point, the feature extraction phase of CVA is completed;
[0108] Step (4): Use FOP to recursively update features and construct new features that are sensitive to early faults. The specific implementation process is as follows:
[0109] Update of the eigenvalues of the CVA canonical variables:
[0110] λ K,I =(1 - ∈)λ k-1,i +∈|f i | 2 (44)
[0111]
[0112] Subsequently, recursively update the eigenvectors:
[0113]
[0114] Step (5): Use LOF to construct a subspace detection metric, and the specific implementation process is as follows:
[0115]
[0116] where s is a sample in S, and s f is the F - th nearest neighbor sample point, and is otherwise defined as follows:
[0117] Rd(s, s f ) = max{d(s, s f ), Fd(s f )} (49)
[0118]
[0119] where Fd(s f ) represents the distance from s f to its farthest point within the requirement of the number of nearest neighbors;
[0120] Step (6): Use Bayesian inference to perform statistical merging on the statistics of all subspaces to obtain the final comprehensive joint statistic, and the specific implementation process is as follows:
[0121] Step (6.1): Use Bayesian inference to perform statistical merging on the statistical data of all subspaces. The failure probability calculation formula for a single subspace x b is;
[0122]
[0123] P(X b ) = P(X b |N)P(N) + P(X b |F)P(F) (52)
[0124] where the expressions of the conditional probabilities P(X b |N) and P(X b |F) are respectively
[0125]
[0126] Step (6.2): Combine the detection results of all different subspaces through Bayesian inference, and calculate the sum of the joint statistics:
[0127]
[0128] Set a confidence level of 0.99 and construct the statistical confidence limit BIC for the offline stage L , and complete the construction of the distributed early fault detection model;
[0129] Step (7): Collect noisy sample data at the new sampling moment Process the obtained data using the mean and standard deviation obtained in step (1.2) to obtain the standardized test sample data XT = ∈R 1×33 ;
[0130] Step (8): Divide the process variable XT according to the data-driven method in step (2) to obtain the subspace division result;
[0131] Step (9): For the process variable data at the current moment in each subspace and the sample data standardized in the previous 9 moments, construct an early fault detection model in each subspace in the manner of step (3), and use CVA for dynamic feature extraction to obtain the feature vector space;
[0132] Step (10): According to the method in step (4), based on the feature vector space obtained in step (9), combine FOP recursive update to construct new features sensitive to early faults;
[0133] Step (11): In each subspace in the manner of step (5), use the LOF method to combine the sensitive new features of the online data to construct a detection metric, and obtain the LOF value of the detection result for each subspace;
[0134] Step (12): Obtain the Bayesian inference fusion results BIC related to quality and unrelated to quality of the test data set in the manner of step (6), and compare with the statistical confidence limit BIC in the offline stage L If the statistic BIC of the online collected sample > BIC L , then according to the judgment criterion, the system has a fault and an alarm is issued; if BIC < BIC L , then there is no fault in the system.
[0135] The proposed method is compared with the method of single integrated early fault detection using Moving Weighted Average-Principal Component Analysis (MWA-PCA). The detection results of the two methods for the step fault 08 in TE are shown in Table 1;
[0136] Table 1: Detection results for step fault 8
[0137]
[0138] The detection results of the two methods for the slow drift fault 13 in TE are shown in Table 2:
[0139] Table 2: Detection results for slow drift fault 13.
[0140]
[0141] Finally, the simulation result diagrams of the proposed method for the two faults are shown in Figure 3 . It can be found that the detection effect of the proposed method is excellent, with significantly better accuracy and timeliness, and can achieve timely and effective early warning for early faults.
[0142] The above embodiments are only used to explain the specific implementation of the present invention, rather than to limit the present invention. Therefore, any changes made according to the shape and principle of the present invention should be covered within the scope of the present invention.
Claims
1. A distributed early fault detection method based on first-order perturbation theory, characterized in that: The following steps are involved: Step (1): Industrial process data collection and preprocessing. The specific process is as follows: Step (1.1): Collect sample data under normal working conditions in the industrial process as the training data set, denoted as Step (1.2): The collected sensor data is preprocessed using the z-score standardization method to obtain a standardized data set X = [x1, x2, ..., x m ]∈R n×m , where n is the number of data set samples, m is the number of process variables, and R is a real number set. The specific calculation formula is as follows: where x i is the data of a single sensor variable in the standardized dataset, x i ∈R n×1 , u is x i The mean of i The mean of . Step (2): The correlation between different variables uses mutual information (MI) to process linear and nonlinear relationships. Based on the correlation between variables, a data-driven method is used to divide the subspace, shifting the perspective from global to local, and highlighting sensitive information that helps early fault detection. The specific implementation process is as follows: Step (2.1): Calculate the variable x i ∈R n×1 The MI between (i=1,2,…,m) and all other variables, and summarize these values into a matrix form, called the mutual information matrix. The mutual information matrix is defined as follows: Among them, I(x i ,x j ) represents the mutual information value between two variables, p(x i ,x j ) represents the variable x i and x j The joint distribution probability between and; Step (2.2): Use the spectral clustering algorithm to construct a graph using the similarity matrix (MI matrix) of the data, and implement clustering through the eigendecomposition of the Laplace matrix of the graph. The detailed process is as follows: The similarity matrix W is considered as the adjacency matrix of the graph, where W ij Represents variable x i and x j The similarity between them is calculated based on the similarity matrix D: Among them, the non-diagonal elements are 0, and then the normalized Laplace matrix (symmetric form) L is calculated: L=I-D -1 / 2 WD -1 / 2 (5) Perform eigendecomposition on the Laplace matrix L and obtain the eigenvalues λ1,λ2,…,λ n and the corresponding eigenvectors v1,v2,…,v n , arrange the first k eigenvectors in descending order to form an n×k eigenvector matrix U: U=[v1,v2,…,v k ] (6) Normalize each row of matrix U to get matrix T: Then cluster the row vectors of T to get the partitioning result, which is used to divide the process variables and construct the subspace; Step (3): Industrial processes are full of a large amount of time-series-related process data. Reasonable use of this feature can help improve the performance of early fault detection. In the construction of the subspace early fault detection model, canonical variable analysis (CVA) can determine the linear combination between past variables and future variables, so it can be used to extract the dynamic characteristics of variables. The specific implementation process is as follows: Maximize the past variable y p,r and the future variable y f,r The correlation coefficient between them is used to extract dynamic features. r represents the current sampling time, p represents the past sampling point, f represents the future sampling point, and the past variable y p,r and the future variable y f,r The description is as follows: The past matrix Y p and the future matrix Y f The description is as follows: AND p =[and p,p+1 and p,p+2 … and p,p+N ]∈R mp×N (10) AND f =[and f,f+1 and f,f+2 … and f,f+N ]∈R mf×N (11) Where m represents the number of process variables, p represents the window length, and N is the dimension of the Hankel matrix to be set; Then, the sample covariance and cross-covariance of past and future observations can be estimated as: CVA seeks to obtain the best linear combination that maximizes the correlation coefficient, which can be achieved by performing a singular value decomposition on the scaled Hankel matrix H: where U and V contain the left and right singular column vectors of H, and Σ is a semi-positive definite diagonal matrix whose diagonal elements are the singular values, i.e., eigenvalues, of H. The canonical variable s r The calculation is as follows: Regular variable s r It consists of two parts, where the larger singular value represents the principal component space, i.e., the eigenvector, and the smaller singular value represents the residual space. At this point, the feature extraction phase of CVA is completed; Step (4): First-order perturbation theory (FOP) is used to study the behavior changes of the system when it is subjected to small disturbances. In the process of feature extraction, FOP theory can identify and highlight the parts that have a greater impact on the industrial process to highlight the main changes in process variables, thereby more effectively analyzing abnormal deviations caused by early faults. Therefore, FOP is used to recursively update features and construct new features that are sensitive to early faults. The specific implementation process is as follows: The eigenvalues of the typical variables of CVA are updated as follows: l k,i =(1-∈)λ k-1,i +∈∣f i ∣ 2 (18) Among them, ∈ is a positive number close to zero, k is the sampling time, λ i is the diagonal element of Σ in formula (15), i.e., the eigenvalue, s k-1,i is the i-th eigenvector, corresponding to the eigenvalue, belonging to the canonical variable s in formula (16) r The principal component space corresponding to the larger singular value is the feature space; then the feature vector is recursively updated: Step (5): After obtaining the recursive feature space S that is more sensitive to early faults, the subspace detection index is constructed through the local outlier factor (LOF). LOF is used to find outliers in a multidimensional data set. The size of this value indicates the degree of sample outliers. Since the samples in the training data set are normal, the feature space sample values are similar at each sampling time. Considering that the feature space sample values of the fault detection samples can be used as outliers of normal samples, LOF can be used to construct the detection metric. The specific implementation process is as follows: Among them, s is a sample in S, s f is the fth nearest neighbor sample point, and is defined as follows: Rd(s,s f )=max{d(s,s f ),Fd(s f )} (23) Among them, Fd(s f ) means s f The distance to the farthest point within the required number of nearest neighbors; Step (6): Use Bayesian inference to statistically combine the statistics of all subspaces to obtain the final comprehensive joint statistics. The specific implementation process is as follows: Step (6.1): Use Bayesian inference to statistically combine the statistics of all subspaces, and a single subspace x b The failure probability calculation formula is: P(X b )=P(X b |N)P(N)+P(X b |F)P(F) (26) Among them, the conditional probability P(X b |N) and P(X b |F) are respectively In the above formula, N and F represent the normal state and fault state respectively; P(N) and P(F) represent the prior probability under normal and fault conditions respectively; when the confidence level is determined to be ε, P(N) = ε, P(F) = 1-ε. b,new is the statistical result of the b-th subspace in the online detection data set, J b,lim is the corresponding control threshold; Step (6.2): Combine the detection results of all different subspaces through Bayesian reasoning and calculate the sum of the joint statistics: Set the confidence level to 0.99 and construct the statistical confidence limit BIC for the offline phase L , complete the construction of distributed early fault detection model; The above steps (1) to (6) are the offline modeling stage of the method of the present invention, and the following steps (7) to (12) are the implementation process of the online analysis of the method of the present invention. Step (7): Collect sample data at the new sampling time The obtained data is processed using the mean and standard deviation obtained in step (1.2) to obtain the standardized test sample data XT = R 1×m ; Step (8): Divide the data into process variables XT according to the data-driven method in step (2) to obtain the subspace partitioning result; Step (9): The current moment data of the process variables in each subspace and the sample data of the previous 2p-1 moments after standardization are combined, and an early fault detection model is constructed in each subspace according to the formula in the manner of step (3), and CVA is used to perform dynamic feature extraction to obtain a feature vector space; Step (10): According to the method in step (4), based on the feature vector space obtained in step (9), combined with FOP recursive update, a new feature sensitive to early faults is constructed; Step (11): In each subspace, the LOF method is used in combination with the sensitive new features of the online data to construct a detection metric in the manner of step (5), and the LOF value of the detection result of each subspace is obtained; Step (12): Obtain the Bayesian inference fusion result BIC of the test data set in the same way as in step (6), and compare it with the statistical confidence limit BIC in the offline stage. L A one-by-one comparison is performed. If the confidence limit is exceeded, an alarm is triggered and a fault is considered to have occurred. If the limit is not exceeded, the system is considered to be operating normally.