Incipient Fault Detection for Complex Industrial Processes with

Subscriber access provided by UNIV OF DURHAM

Article

Incipient Fault Detection for Complex Industrial Processes with Stationary and Nonstationary Hybrid Characteristics Chunhui Zhao, and Biao Huang Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b00233 • Publication Date (Web): 20 Mar 2018 Downloaded from http://pubs.acs.org on March 26, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Incipient Fault Detection for Complex Industrial Processes with Stationary and Nonstationary Hybrid Characteristics Chunhui Zhao1,2*, Biao Huang2* 1 State Key Laboratory of Industrial Control Technology, College of Control Science and Engineering, Zhejiang University, Hangzhou, 310027, China 2 Department of Chemical and Materials Engineering, University of Alberta, Edmonton, AB, T6G-2V4, Canada Abstract For nonstationary process which has a time-varying mean or a time-varying variance or both, it can be difficult to detect incipient disturbances which may be hidden by the time-varying process variations. Besides, stationary and nonstationary characteristics may co-exist in complex industrial processes which, however, have not been studied for process monitoring. In the present work, a triple subspace decomposition based dissimilarity analysis algorithm is developed to detect incipient abnormal behaviors in complex industrial processes with both stationary and nonstationary hybrid characteristics. The novelty is how to comprehensively separate the stationary and nonstationary process characteristics and describe them respectively. First, a stationarity evaluation and separation strategy is proposed to decompose the data space into three subspaces, revealing the linear stationary process characteristics, the nonlinear stationary process characteristics and the final nonstationary process characteristics. Then, a triple subspace distribution monitoring strategy is proposed to quantitatively evaluate the changes of linear and nonlinear stationary and nonstationary distribution structures. The paper demonstrates that the new method has better performance in detection of incipient abnormal behaviors that

*

Corresponding author: E-mail address: [email protected], [email protected]

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

are responsible for distortion of the underlying covariance structure in industrial processes with both stationary and nonstationary hybrid process characteristics. Its feasibility and performance are illustrated with two case studies of real industrial processes. Keywords: Fault detection; stationarity and nonstationarity; subspace decomposition; dissimilarity analysis; incipient fault. 1. Introduction With increasing demands in plant safety and operation efficiency, process monitoring has been drawing increasing attentions in recent decades for industrial applications1-8 to detect process upsets, equipment malfunctions, or other special events as early as possible so as to identify and remove the factors causing those events. Nowadays, the plant-wide processes are usually characterized by large-scale, multiple operation units and the resulting complex correlations among variables, and monitoring of such processes has become an important issue. Meanwhile, with the advent of advanced instrumentation and automation systems, process data has become abundant and a large number of highly correlated variables are archived. There have been extensive theoretical and practical studies in which data-driven methods, in particular multivariate statistical process monitoring(MSPM) methods, progress rapidly. During the past thirty years, principal component analysis (PCA)9,10 and partial least squares (PLS)11-12 has been the typical methods that have drawn increasing attentions because they rarely need a priori process knowledge. Process monitoring using PCA or PLS is based on a predefined model built from normal process data. Both methods are able to reduce the dimensionality of the monitoring space by projecting measurement data onto a lower-dimensional latent space. The processes are then monitored in these subspaces by utilizing distance-based statistics


Page 2 of 49

Page 3 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


to timely detect any deviation from normal or “in-control” region. Faults in a plant are unavoidable, revealing the unwarranted change with respect to the nominal fault-free behavior. It is important to detect the early occurrence of a small process deviation before the serious failure of the overall process so that it can be

handled

timely

to

prevent

undesired

consequences.

Many

successful

applications13-20 have reported the practicability of multivariate statistical analysis methods for fault detection. However, most of them consider the detection of significant faults with a certain fault magnitude that cannot be buried by the variability used in calculation of monitoring statistics because these methods evaluate the process deviations by calculating a distance metric. So these methods may not be sensitive enough to detect incipient changes that may only change the variable covariance structure or the underlying distribution. This problem has caused concern and some previous work has been presented to address this issue. A new monitoring index was proposed by Kano et al.17 to check the process changes along each distribution direction derived by PCA, which, however, did not evaluate changes of the process variance. Then a dissimilarity index was proposed by Kano et al.20, called DISSIM method, which can quantitatively calculate the distribution deviation from the normal condition. Recognizing that a change of operating condition can be reflected in both distribution directions and variances, DISSIM method decomposed the distribution of time-series data to reveal the changes of underlying operating condition. They used Karhunen-Loeve (KL) expansion21 to transform the measurement data, resulting in the same eigenvectors, i.e., the same directions. Thus the difference between eigenvalues can reveal the difference of the process variances. Since then, further development of the DISSIM method have been reported for applications to batch processes22-23 and nonlinear processes.



However, although the above mentioned methods have been widely used for fault detection which evaluate the changes of distribution structure and are thus more sensitive to incipient changes, they, however, assume that the process is stationary which can not be well satisfied in practice. In particular, complex industrial processes in general show nonstationary process characteristics besides stationary ones, revealing a time-varying mean or a time-varying variance or both24. Such behavior may be caused by seasonal changes, processes that involve emptying and filling cycles, throughput changes, the presence of unmeasured disturbances, and operator interventions, etc25. As the statistical properties of nonstationary variables vary along time, the traditional multivariate statistical process monitoring techniques are inadequate to deal with the modeling and monitoring of nonstationary processes. Nonstationary process monitoring is a difficult task which has only been sporadically addressed. One typical method is adaptive strategy26 that may consecutively update the monitoring model to capture the nonstationary process behaviors. This approach, however, may undesirably accommodate and adapt to slow-varying fault conditions. Previous work27,28 have mentioned that an Autoregressive Integrated Moving Average (ARIMA) time series model can be applied in some cases where data show evidence of non-stationarity, where an initial differencing step can be applied one or more times to eliminate the non-stationarity. Application of ARIMA for process monitoring has been discussed and described with respect to both single-input single-output systems (SISO) and multiple-input multiple-output (MIMO) processes29. An integrated variable is defined to be a nonstationary variable with the difference stationary. Engle and Granger30 introduced the concept of cointegration of a set of integrated nonstationary variables to describe the relation between the variables. More precisely, they showed that it is possible for a linear combination of a set of integrated


Page 4 of 49

Page 5 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


nonstationary variables to be stationary if the nonstationary variables are integrated of the same order and share common trends. Zhao et al.25 proposed the use of the cointegration analysis (CA) method to estimate both static and dynamic long-run equilibrium relationships between nonstationary time series that produce stationary residual sequences which were further exploited to address the process monitoring problem for nonstationary dynamic processes. Sun et al.31 developed an online fault diagnosis strategy for nonstationary processes with no priori fault information which automatically and in real-time isolated multiple faulty variables that are responsible for the abnormal nonstationary operation by borrowing the least absolute shrinkage and selection operator (LASSO) trick. For the CA method, which is designed to analyze nonstationary variables being integrated of the same order, its scores are not guaranteed to be orthogonal to each other that may cause inconvenience for process monitoring. Besides, it is limited to linear nonstationary relationship, which can not explore the nonlinear stationary components. It is common that the complex industrial processes possess both stationary and nonstationary characteristics, which can be linear or nonlinear. For each single variable, it may contain both stationary and nonstationary components and considering the complex correlations among variables it may not be proper to simply classify the observed single variable to be stationary or nonstationary. However, the above mentioned methods only focused on the analysis and monitoring of nonstationary variables by simply distinguishing them from the stationary variables using nonstationarity test. Besides, they borrowed the traditional distance-based monitoring index of MSPM which may not be sensitive to incipient abnormality. To the best of the authors’ knowledge, a reliable framework of incipient fault detection for processes with both stationary and nonstationary hybrid characteristics has not



been reported in the literature. Here the term of hybrid process is used to define the complex process in which the process variables possess both stationary and nonstationary characteristics. To solve the above mentioned problem, two key points should be addressed: (1) how to separate the stationary and nonstationary process characteristics comprehensively; (2) how to monitor incipient changes relevant to the discussed different kinds of process variations. This paper proposes a triple subspace decomposition based dissimilarity analysis algorithm for online incipient fault detection of complex industrial processes with stationary and nonstationary hybrid characteristics. Here, the underlying characteristics are represented by latent components which can be obtained by linear or nonlinear combination of process variables. First, it gives rise to a total subspace decomposition of data variation which can separate the linear stationary components from the nonlinear stationary components and the final nonstationary ones so that three different process characteristics can be well distinguished from each other. Second, it analyzes the changes of distribution structures for fault detection in all three subspaces to quantitatively reveal incipient changes from different perspectives. The paper demonstrates that the new method has better performance in detecting incipient faults associated with the underlying distribution structure of different process characteristics. The main contribution is summarized as below: (1) It is the first time that stationary and nonstationary hybrid characteristics are noticed and explored for data-driven process monitoring. (2) The stationary and nonstationary process characteristics can be separated from each other both linearly and nonlinearly, which enhances process understanding. (3) Incipient changes can be detected by comprehensively evaluating the changes


Page 6 of 49

Page 7 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of distribution structures for both stationary and nonstationary characteristics. The remainder of the paper is presented as follows. First, the motivation for the current work is discussed. Subsequently, the proposed algorithm is developed for fault detection. Third, we present the application results of the proposed method to real industrial processes. Comparison is made with previous methods. The conclusion is drawn in the last section. 2. The proposed methodology 2.1 Problem statement For complex industrial processes subject to different factors, two recognitions are first given: (1) The process variables may present typical nonstationary characteristics with the linear or nonlinear variable correlations; (2) Although some of the underlying process characteristics are nonstationary, the others may stay stationary. Based on the above recognitions, different underlying process characteristics should be distinguished from each other. It includes two important components, i.e., how to separate the stationary and nonstationary process characteristics more comprehensively and how to capture the changes of their different distribution structures so that incipient changes can be detected in time. Instead of stationarity test for each single variable and completely isolating each into the stationary or nonstationary one, the stationary and nonstationary latent components are used as the analysis subjects which should be separable from each other. These latent components can come from the linear or nonlinear combinations of multiple variables. Considering the problem of nonlinear analysis methods, such as computation complexity, prone to overfitting, misspecification, the analysis of linear relationship should be given the priority. That is, if some variables can be linearly combined to get the stationary features, they should be considered for linear analysis first. Based on



Page 8 of 49

the above consideration, the whole measurement space can be separated into three subspaces as shown in Figure 1 by both linear and nonlinear analyses. They are linear stationary subspace (LSS) and nonlinear stationary subspace (NSS) which are combined to describe the stationary process characteristics. The left subspace only covers nonstationary characteristics, termed nonstationary subspace (NSS). In the following subsections, triple subspace decomposition is conducted and subsequently dissimilarity analysis is implanted in different subspaces for incipient fault detection of complex industrial processes with both stationary and nonstationary characteristics. 2.2 Linear stationary subspace (LSS) decomposition Assume that J process variables are measured online at k=1,2,…, K time instances. It forms the regular data analysis unit, denoted as X ( K × J ) . The variables are normalized to have zero mean and unit standard deviation. In this way, the normalized data set X is available. Here for simplicity, the normalized data sets are denoted by the same symbols. Apply PCA algorithm to get the variations along each variable correlation direction, T = XP  ˆ = X − TPP Τ E = X − X

(1)

where T is PCA scores and P denotes the loadings to calculate scores from X as well ˆ represents the process variation information as the reconstruction relationships. X

that can be recovered by the PCA scores. Here, all PCA scores with nonzero variances are retained which thus can recover all process variations. Apply the Augmented Dickey-Fuller (ADF) test32 to determine the stationarity of the PCA scores that can separate the stationary scores from the nonstationary scores, which are denoted as Ts ( K × Rs ) and Tns ( K × Rns ) , respectively. Rs denotes the number of the identified stationary scores and Rns denotes the number of identified


Page 9 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


nonstationary scores. Note it is not necessary that the PCA scores with the largest variances are the nonstationary ones by ADF test. Correspondingly, the original measurement data space can be separated into two subspaces,

% = T P Τ = XP P Τ  X s s s s s % Τ Τ  X ns = Tns Pns = XPns Pns

(2)

where, Ps denotes the part of loadings P corresponding to the stationary scores while

Pns for nonstationary scores. Since the PCA scores are orthogonal with each other,

% denotes the their linear combinations do not change the stationarity. Therefore, X s stationary subspace resulting from the linear combination of stationary scores and

% denotes the nonstationary subspace from the linear combination of nonstationary X ns scores. The subscripts s and ns denote stationary and nonstationary respectively. In this way, one stationary subspace is first separated from the variables by linear analysis, termed LSS. 2.3 Nonlinear stationary subspace (NSS) decomposition % , since they are linear combination of some orthogonal nonstationary For X ns

components from Eq. (2), all variables are nonstationary from which no stationary component can be derived by linear combination. However, the above analysis is only from the linear perspective and these variables may present nonlinear relationships. In the following, nonlinear analysis can be conducted for further extraction of stationary components. Here, a nonlinearity evaluation analysis strategy5 is used to check

% . If there whether there is a nonlinear relationship among different variables of X ns are nonlinear variable correlations, the following analysis should be conducted. Otherwise, stationary component can only be extracted by linear combination. That is,



Page 10 of 49

no nonlinear stationary subspace can be separated. For extraction of nonlinear stationary components, kernel PCA algorithm33,34 is % . Considering a nonlinear mapping Φ : x ∈ R J → z ∈ R h , used to further analyze X ns

% is transformed into the high-dimensional feature space, the nonstationary data X ns % ) , where the dimension of the feature space, h , can be arbitrarily large or even Φ(X ns infinite. The basic idea is to perform the linear PCA analysis in the high-dimensional feature space which can be regarded to be nonlinear PCA analysis in the original input space. That is, a linear data structure is more likely to be available after high-dimensional nonlinear kernel mapping35 in which the higher-dimensional linear space is regarded as the feature space (F). The scores in the feature space can be obtained using the following transformation: −1 2 % )V = Φ% ( X % )Φ% ( X % )T HΛ −1 2 = KHΛ % TΦ = Φ% ( X ns ns ns

(3)

% ) is the mean centered form of Φ% ( X % ). where, TΦ are the KPCA scores. Φ% ( X ns ns V

is

the

eigenvector

Λ = diag(ξ1 , ξ 2 ,..., ξ K )

matrix

from

% ) . Φ% ( X ns

H = [α1 , α 2 ,..., α K ]

and

are respectively the eigenvector matrix and nonzero

% which are all known. In some cases, the last few eigenvalues in eigenvalues of K Λ are very closely to zero which can be excluded to avoid the singular problem. A

% )Φ ( X % ) T is defined and its visual expression can K × K kernel matrix K = Φ ( X ns ns be readily obtained based on the choice of specific kernel function36 where the Gaussian kernel function is the most commonly used in practice37 and the

% = Φ% ( X % )Φ% ( X % ) T can also be explicitly derived. mean-centered kernel matrix K ns ns Therefore, the KPCA scores ( TΦ ) with nonzero variances are retained and ordered according to the variances denoted by Λ , from which, the correlations are removed


Page 11 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


% ) can be recovered. and all process variations of Φ% ( X ns % ) , have been transformed to In this way, all invisible mapped data matrices, Φ% ( X ns the visible scores in the kernel space. Apply the Augmented Dickey-Fuller (ADF) test37 to determine the stationarity of the kernel PCA scores which can separate the stationary scores from the nonstationary scores, which are denoted as TΦ , s ( K × RΦ , s ) and TΦ ,ns ( K × RΦ ,ns ) , respectively. It is noted it is not necessary that the KPCA scores with the largest variances are the nonstationary ones by ADF test. RΦ , s denotes the number of stationary scores and RΦ ,ns denotes the number of nonstationary scores. Correspondingly, the KPCA loadings are separated into the stationary and nonstationary parts, Vns and Vs , as well as the variances, Λ ns and Λ s . The

% ) can be nonlinearly separated from the others, stationary part of Φ% s ( X ns

% )=T VT Φ% s (X ns Φ,s s % % )V V T = Φ (X ns

s

s

% )Φ% (X % )T H Λ −1H TΦ% (X % ) = Φ% (X ns ns s s s ns T % % −1 % =KH Λ H Φ (X ) s

s

s

ns

% )=T V T Φ% ns (X ns Φ , ns ns % % )V V T = Φ(X ns ns ns % )Φ% (X % )T H Λ −1H TΦ% (X % ) = Φ% (X ns ns ns ns ns ns T % % −1 % =KH Λ H Φ (X ) ns

ns

(4)

ns

(5)

ns

In this way, the second stationary subspace is separated from the variables which is achieved by nonlinear analysis, termed NLSS. Although it is not visually known resulting from the invisible mapped data matrices, the corresponding scores are visible which are used for monitoring. After the separation of the second stationary subspace, the left subspace is deemed to only cover nonstationary characteristics. It is possible that the hybrid characteristics may change with time. In that way, the variable



correlations may change which requires that the models should be updated to extract stationary variations. Therefore, some updating strategy may be needed.

2.4 Property analysis and discussion The triple subspace decomposition algorithm realizes the comprehensive stationary and nonstationary decomposition of X space supervised by linearity and nonlinearity, covering the linear stationary subspace (LSS), nonlinear stationary subspace (NLSS) and the final nonstationary subspace (NSS), which present different process characteristics. The meanings of different subspaces are summarized in Table 1. Some properties are analyzed and discussed as below. Like the PCA algorithm, the orthogonality among subspaces holds, which is shown in the following:

% TX % = 0 , Φ% ( X % )Φ% ( X % )T = 0 . (1) X s ns s ns ns ns

(2)

∀t i , t j ∈ {Ts , Tns } : t i t j = 0, i ≠ j T

and

∀t i , t j ∈ {TΦ , s , TΦ ,ns } : t i T t j = 0, i ≠ j

.

For both, they can be readily proved considering the orthogonality of PCA and KPCA scores.

% (3) No stationary component can be derived from X ns by linear projection.

% ) by nonlinear mapping. (4) No stationary component can be derived from Φ% ns ( X ns % θ and s is stationary resulting from linear Proof: Assume there is s = X ns

% . According to Eq. (2), we know that X % is linear combination combination of X ns ns of some orthogonal nonstationary components so that s = Tns Pns Τθ = Tns θ* . That is, s is in fact a linear combination of some orthogonal nonstationary components which thus should be nonstationary. Therefore, the assumption is not satisfied. The statement in (3) is thus proved. Similarly, the statement in (4) can be proved according to Eq.


Page 12 of 49

Page 13 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(5). This algorithm performs triple subspace decomposition by both linear and

% are nonstationary. Besides, nonlinear analyses. For Eq. (2), all the vectors in X ns from these nonstationary variables, no stationary component can be derived by linear combination. However, nonlinear analysis can be conducted to check whether the

% ) , since they are nonlinear stationary components can be further derived. For Φ% ns ( X ns combination of some orthogonal nonstationary components in the higher-dimensional kernel space from Eq. (5), all the vectors are nonstationary. Besides, since the measurement data have been explored by both linear and nonlinear analysis for extraction of stationary components, the final subspace can not be further decomposed for extraction of stationary components, only revealing nonstationary information. It is noted that there may be two extreme issues in practice, i.e., no linear stationary correlations and no nonlinear stationary correlations. For the first case, no linear stationary scores can be extracted and only KPCA can be used for extraction of nonlinear stationary variations. For the second case, only linear PCA is needed to extract linear stationary scores. In the present work, the linearity and nonlinearity of variable correlations are considered for subspace separation which are two major kinds of variable characteristics. Of course, other criterions can be used for separation, such as consideration of different machines where the variables are collected or date when the production is conducted. Although the stationary feature extraction is conducted based on PCA and KPCA methods, it can be readily extended to other methods. Many linear statistical methods can be readily extended to their nonlinear versions38-41 which can be used for nonlinear analysis of stationary variations.

2.5 Triple subspace based DISSIM analysis



Page 14 of 49

Based on the subspace decomposition mentioned before, thee subspaces are available, including two stationary subspaces, either linear or nonlinear, and one nonstationary subspace. In different subspaces, process monitoring should be performed to check the changes of different process characteristics.

(1). DISSIM for linear stationary subspace % , one reference data set X % ( N × J ) is For the linear stationary subspace X s s ,r s ,r chosen (where N s ,r is the number of reference samples) and then a moving window

% ( N × J ) ) with length of N (X s , w is used for dissimilarity analysis to evaluate the s,w s,w changes of distribution structure along time direction. The subscript r is used for reference distribution and subscript w for moving window data that are consecutively updated along time direction. They cover the same number of variables (J) but may have different number of observations, denoted by N r and N s , w respectively. DISSIM algorithm is performed to explore the distribution difference between two data sets which can be evaluated by defining the following index Ds, A % ,X % ) = 4 ∑ ( λ j − 0.5 )2 Ds = diss ( X s ,r s,w s ,r J j =1

(6)

where, λsj,r denotes the eigenvalues derived from the reference data set and A is

% % limited by the rank of X s , r and X s , w . (2). DISSIM for nonlinear stationary subspace

% ) which cover the nonlinear stationary characteristics, nonlinear For Φ% s (X ns DISSIM algorithm23 is implemented in the high-dimensional feature space to analyze the changes of nonlinear distribution structure.

% ) and a moving window ( Φ% (X % ) ) is The reference data set is set to be Φ% s ,r (X ns s ,w ns used for dissimilarity analysis to evaluate the changes of distribution structure along


Page 15 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


time direction. Here, for expression simplicity, the number of samples in reference data set and moving window is set the same as that for linear stationary subspace. First, the mixed covariance is calculated for eigenvalue decomposition,

% )  T Φ% (X % ) Φ% s ,r (X ns s ,r ns R =     % )  Φ% (X % ) N s ,r Φ% s , w (X ns   s , w ns  % )TΦ% (X % ) Φ% (X % )TΦ% (X % ) Φ% (X s ,r ns s ,w ns = s ,r ns + s ,w ns N s ,r + N s , w N s , r + N s ,w Φ

1 + N s,w

(7)

R Φ PΦ = PΦ ΛΦ

(8)

where, orthogonal matrix P Φ denotes the principal projection directions, and Λ Φ is a diagonal matrix whose diagonal elements are eigenvalues of R Φ .

% ) and Φ% (X % ) are then transformed by using the The data matrices Φ% s ,r (X ns s ,w ns directions with nonzero variances,

YsΦ,r = YsΦ, w =

1

N s ,r N s,r + N s,w N s ,w N s ,r + N s , w

% )P Φ Λ Φ − 2 Φ% s , r (X ns % )P Φ Λ Φ% s , w (X ns

(9)

1 Φ 2 −

The covariance matrices of the transformed data matrix clearly satisfy S Φs ,r + S Φs ,w =

1 1 YsΦ, r T YsΦ,r + YsΦ, w T YsΦ,w N s ,r N s,w −

=

1

% )TΦ% (X % )+Φ% (X % )TΦ% (X % ) ) PΦ ΛΦ Λ Φ 2 P Φ T (Φ% s ,r (X ns s ,r ns s ,w ns s,w ns N s,r + N s,w −

1

= ΛΦ 2 PΦ T R Φ PΦ ΛΦ

−

−

1 2

(10)

1 2

=I where I is R × R identity matrix. The monitoring index is constructed by performing eigenvalue decomposition on the covariance matrix, S Φs ,r or S Φs , w . However, S Φs ,r seems to be invisible to have



Page 16 of 49

% ) , P Φ and Λ Φ . Using the same kernel trick as that unknown components: Φ% s ,r (X ns used in previous subsection, YsΦ, r can be clearly obtained and known,

N s ,r

YsΦ,r =

N s,r + N s,w

1

% )P Φ Λ Φ − 2 Φ% s , r (X ns

% ) T Φ% s , r (X ns % ) Φ% s , r (X  H Φ Π Φ −1 2 Λ Φ −1 2 ns % % N s,r + N s,w Φs , w (X ns )  N s ,r

=

= N s ,r  K s ,rr Φ

(11)

K s ,rw Φ  H Φ Π Φ −1

where, H Φ and Π Φ are respectively the eigenvector matrix and nonzero eigenvalues

Ks

Φ

decomposed

from

K sΦ

which

are

all

known,

in

which,

% ) T % )  Φ% (X % )  T K Φ K Φ  Φ% s ,r (X Φ% s ,r (X ns s ,r ns s , rr s , rw ns Φ Φ Φ −1 2 and =  H Π   =  , P = Φ Φ % % % % % % Φs , w (X ns )  Φs , w (X ns )  Φs , w (X ns )   K s , wr K s , ww 

ΛΦ =

ΠΦ N s ,r + N s , w

.

From

Eq.

(4),

it

is

clear

that

% )Φ% (X % )T K s ,rr Φ =Φ% s ,r (X ns s ,r ns =TΦ , s ,r Vs T Vs TΦ , s , r T . % )Φ% (X % )T H Λ −1/ 2 T T = TΦ ,s ,r Λ s −1/2 H s TΦ% (X Φ , s ,r ns ns s s % Λ −1/2 T T = T Λ −1/2 H T KH Φ ,s ,r

s

s

s

s

Φ , s ,r

Thus, S Φs ,r which is the covariance of YsΦ, r can be obtained for eigenvalue decomposition and the monitoring index DΦ ,s is obtained to evaluate the dissimilarity of nonlinear correlations underlying the distribution structures of two data sets,

% ),Φ% (X % )) = 4 DΦ ,s = diss(Φ% s ,r (X ns s,w ns AΦ , s

AΦ ,s

∑ (λ

Φ, j s ,r

− 0.5) 2

(12)

j =1

where, λsΦ, r, j is the eigenvalue obtained from S Φs ,r and AΦ ,s is limited by the rank

% ) and Φ% (X % ). of Φ% s ,r (X ns s,w ns


Page 17 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(3). DISSIM for the final nonstationary subspace

% ) in the high-dimensional feature For the final nonstationary information Φ% ns ( X ns space, it covers all the left information after extraction of the above mentioned two subspaces. Despite of its nonstationary nature, the differencing information can also be used for monitoring purpose. They are along the remaining directions Vns which are orthogonal to those spanned by the columns of directions Vs . Here, the difference value is calculated as,

% ) = Φ% ( X % ) − Φ% % ∆Φ% ns ( X ns ns , k ns ns , k +1 ( X ns )

(13)

% ) and Φ% % where, Φ% ns ,k ( X ns ns , k +1 ( X ns ) denote the two data sets at the neighboring time. Then DISSIM algorithm is performed to check the changes of distribution in this subspace which is similar to that for nonlinear stationary subspace. The only

% ) are used to replace Φ% ( X % ). difference is that the difference values ∆Φ% ns ( X ns ns , k ns Although the difference values are not visually known, the kernel trick can be used to solve this problem. For example, the Gram matrix of differences can be derived as,

% )∆Φ% ( X % )T ∆Φ% ns ( X ns ns ns T % ) − Φ% % % % % % = (Φ% ns ,k ( X ns ns , k +1 ( X ns ) )( Φns , k ( X ns ) − Φns , k +1 ( X ns ) )

(14)

= K ns ,k ,k Φ − K ns ,k ,k +1Φ − K ns , k +1,k Φ + K ns ,k +1,k +1Φ Thus, the monitoring index DΦ ,ns is obtained to evaluate the dissimilarity of nonlinear correlations underlying the distribution structures of two data sets,

DΦ ,ns

% ),∆Φ% (X % )) = = diss(∆Φ% ns ,r (X ns ns , w ns

4

AΦ ,ns

AΦ ,ns

∑ (λ

Φ, j ns , r

− 0.5)2

(15)

j =1

where, λnsΦ,,rj is the eigenvalue obtained by DISSIM in the high-dimensional feature

% ) and ∆Φ% (X % ) which are space and AΦ ,ns is limited by the rank of ∆Φ% ns ,r (X ns ns , w ns the reference data set and the moving window used for DISSIM in this feature space.



For each of the three subspaces, distribution structure is evaluated to check the changes of operation status. When the concerned process correlations denoted by the two data sets are quite similar to each other, the eigenvalues must be near 0.5, and then monitoring statistics will be near zero. On the other hand, when the two data sets nonlinearly distribute differently from each other, the largest and the smallest eigenvalues should be near one and zero respectively and thus the monitoring statistics will be near one. Along time direction in each subspace, the moving window is updated to get more monitoring statistics to evaluate the distribution difference between the moving windows and the reference one. The control limits are thus determined by using non-parametric kernel density estimation37 for different process characteristics. An algorithmic flowchart has been provided to illustrate the offline modelling for readability as shown in Figure 2. 2.6 Online application For on-line fault detection, the DISSIM index will be evaluated in the above three different subspaces. First, the new sample is normalized using the mean and the variance information obtained from training data. Then, it is decomposed to get three

% different subspaces, the linear stationary subspace ( X new , s ), the nonlinear stationary % % % subspace ( Φ% s (X new , ns ) ) and the final nonstationary subspace ( Φns (X new , ns ) ) along with the difference value in the final nonstationary subspace,


Page 18 of 49

Page 19 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


x% new,s Τ = x new Τ Ps Ps Τ x% new,ns Τ = x new Τ Pns Pns Τ Φ% s (x% new,ns )T = Φ% (x% new,ns )T Vs Vs T % )T H Λ −1H TΦ% (X % ) = Φ% (x% new,ns )TΦ% (X ns s s s ns T % % −1 =k H Λ H Φ (X ) new

s

s

s

ns

(16)

Φ% ns (x% new,ns ) = Φ% (x% new,ns )T Vns Vns T % )T H Λ −1H TΦ% (X % ) = Φ% (x% new,ns )TΦ% (X ns ns ns ns ns −1 T % % =k H Λ H Φ (X ) new

ns

ns

ns

ns

∆Φ% ns (x% new,ns )=Φ% ns (x% new,ns ) − Φ% ns (x% new−1,ns ) where, the related model coefficients are all known which are obtained during the stage of model development. In each subspace, the current new data window that represents the actual operating condition is organized by the new sample and its neighboring historical samples with the predefined window length. The new dissimilarity index D is then calculated in three different subspaces to evaluate the distribution difference between the actual and the reference data set which has been defined during model development. Although the two nonlinear subspaces are not clearly known, the monitoring statistics can be calculated specifically by using kernel trick. If any index is outside the control limit, the current operation condition is judged to present a different covariance structure from the reference one. Besides, it can check which subspace is disturbed by the abnormality. The current data matrix will be updated continuously by moving the time-window forward step-wise and three monitoring statistics are calculated consecutively along time to evaluate the operation status. 2.7 Outline of the algorithm In summary, the procedure of modeling and monitoring for complex industrial processes with stationary and nonstationary characteristics is listed as below. Offline modelling:



(1) For data set X ( K × J ) , data normalization is conducted to make it have zero mean and unit variance. (2) Perform PCA on the normalized data and apply the ADF test to determine the stationarity of the PCA scores and separate the stationary scores from the nonstationary scores which reconstruct the linear stationary subspace. (3) Perform KPCA on the linear nonstationary subspace and apply the ADF test on the KPCA scores to separate the nonlinear stationary scores from the others, which reconstruct the nonlinear stationary subspace. (4) Conduct triple subspace DISSIM to evaluate the changes of process distribution in each subspace in which the reference data as well as the control limit is defined and the difference values are calculated in the final nonstationary subspace. It includes the following substeps: Step (4.1) Acquire normal operating data and determine the size of time window. Then series of data sets can be generated by moving the time window step by step along time. Step (4.2) Randomly choose one from the moving windows and set it as the initial reference model. The other moving windows work as training data sets, which should be normalized by the mean and variance obtained from the reference model. Step (4.3) Calculate the dissimilarity index between the training data sets and the reference model. Then determine the formal control limit using the non-parametric empirical density estimate method. Online monitoring: (1) For the new sample, it is normalized using the information from the training data.


Page 20 of 49

Page 21 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(2) The normalized sample is separated into three different subspaces. (3) In each subspace, the new moving window is obtained. (4) Perform dissimilarity analysis in different subspaces. If they report out-of-control index D, it means the current operation condition is judged to present different covariance structure from the reference one, revealing abnormal changes of the corresponding process characteristics. Otherwise, the process is operating normally in the concerned subspace. 3. Case Studies In this section, the proposed method is applied to two real industrial processes with hybrid stationary and nonstationary characteristics. The data are from different industrial fields. One is for the thermal power process which covers 159 variables and the other is for cigarette production process which has 23 variables. 3.1 Large-scale thermal power process In this subsection, the proposed fault detection method is applied to a real industrial process, 1000MW ultra supercritical unit with large capacity, high parameter and low energy consumption. In a coal-fired power plant, the prime mover is steam driven. By heating, water is transformed into steam, and it is then condensed to push a turbine of power generators to produce electricity. Thermal system of a power plant is the main place where the energy converts and mainly includes two subsystems, steam turbine system and boiler system. The boiler system heats the water to produce steam that has high pressure and temperature, and then transports the steam to the steam turbine system to drive the electrical generator. The chemical energy of coal is transformed by the boiler system to the heat energy. The electrical generator is driven by the steam turbine system which can transform the heat energy to the mechanical energy. Finally, the mechanical energy is transformed to the electrical energy by the electrical



generator. Figure 3 shows a generic diagram of a 1000MW ultra-supercritical thermal power unit. Because of the complexity of the thermal system, the process may have nonlinear, time-varying, multivariable and correlated characteristics, which may cover both stationary and nonstationary characteristics, thus providing a good example for illustrating the proposed algorithms. The data are collected from a 1000MW ultra supercritical unit in one thermal power plant with sampling time of 3min. In the present work, 159 variables are used, which involve pressure and temperature of primary steam and related steam, water level, flow etc. Due to variable load conditions, it is a typical stationary and nonstationary hybrid process. Considering the changes of process correlations under wide variable load conditions, the process data for operation with load varying within a certain region [500 1000] MW are collected which are deemed to have the similar variable correlations around the load condition of 750 MW and nonstationary process characteristics. 960 normal samples are used for model development and another 960 samples are used to test the proposed method, which cover two fault cases. First, the triple subspace decomposition is conducted to separate the original measurement space into three subspaces, including two stationary subspaces from both linear and nonlinear perspectives and one nonstationary subspace. As shown in Figure 4 (a), the fourth variable (Condenser temperature (℃)) is separated into linear stationary and linear nonstationary parts, which reveals that after the separation of stationary part, the nonstationarity is more obvious. In the kernel space, the linear nonstationary part is further separated into the nonlinear stationary part and the final nonstationary part, which are comparatively shown in Figure 4 (b) using kernel scores since the measurement sample is not visually known in the kernel space. From both subplots, the stationary information can be effectively extracted which makes the


Page 22 of 49

Page 23 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


remaining more nonstationary. The decomposition result is summarized in Table 2, where the process variations are quantitatively evaluated for each subspace with the corresponding model order. By comparing these variations, it is clear that only 19.89% variations of the original measurement space are linear stationary while the others are linear nonstationary. Then these linear nonstationary variations are further analyzed by nonlinear analysis. It is found that the first two KPCA scores are judged to represent nonstationary variations, which cover only 74.51% of the total process variations, while most of the nonlinear variations (5.60%) can be transformed into stationary components by nonlinear analysis. Then, the triple subspace monitoring strategy is developed for incipient fault detection where the time window is set to cover 400 samples to represent the actual operating condition which is organized by the new sample and its neighboring historical samples within the predefined window length. The monitoring statistics are calculated for both stationary and nonstationary parts in different subspaces and comparatively shown in Figure 5. In the current work, without special statement, the x-axis label denotes the current time although a moving window with its neighboring historical samples may needed for the proposed method. It is noted that in Figure 5(b), the statistic is calculated for the nonlinear components without differencing. It is clear that the D2 statistic shows an obvious increase trend for the nonstationary parts in both linear and nonlinear subspaces, which can not be enclosed by a proper confidence region and is thus not sensitive to faults. In contrast, the stationary part presents more stable values of monitoring statistic, which, thus, can better enclose the normal region. Here, 99% confidence level is used for definition of control limit. For Fault #1, the press of cooling water slowly increases from the 46th sample in the condenser. The fault detection result based on the proposed method is shown in



Figure 6. From Figure 6, the fault can be correctly detected from the very beginning due to the use of moving window. Also, it is found that the fault can be detected in all subspaces. In contrast, using the conventional PCA method and the conventional DISSIM method without separating stationary and nonstationary information, the detection time shows some delay (50 samples delay for PCA and 37 samples delay for the conventional DISSIM) which are not shown here for simplicity. Therefore, without considering the mixed effects of stationary and nonstationary characteristics, the fault detection sensitivity is jeopardized. For Fault #2, the outlet press of the circulating water pump starts slowly increases from the 410th sample. Using the proposed method as shown in Figure 7, the small disturbances can be first detected from the 143rd time window, that is, the 542nd sample, in the nonlinear stationary subspace. In contrast, using the conventional PCA method, the fault can not be detected as shown in Figure 8(a). Using the conventional DISSIM method without considering the effects of nonstationary characteristics, the fault is not noticed until the 368th time window as shown in Figure 8(b). Comparing the results shown in Figure 7 and Figure 8(b), the fault can be detected using the proposed method from the 366th time window in the linear stationary subspace besides the earlier detection in nonlinear stationary subspace, revealing the improvement of fault detection performance by separation of different process characteristics. For detecting both faults, the improved performance of the proposed method has been verified, revealing its superiority for incipient fault detection. Comparatively, incipient fault may not be effectively and timely detected by the other methods since the incipient disturbance may be hidden by the nonstationary variations which thus can not be highlighted.


Page 24 of 49

Page 25 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3.2 Cigarette production process In this subsection, the proposed fault diagnosis method is tested for cigarette production process. Cut-made process is one of the most important procedures in cigarette production process. It determines the flavor and style characteristic of cigarette which is critical to the final cigarette quality. Two major operation machines, SIROX warming and humidification machine and cylinder drier (KLD), are included for this process as illustrated in Figure 9. The leaf-silk is inflated and dampened in the first machine. Then the leaf-silk is dried in the second machine using the barrel heated with saturated steam to reduce the moisture content from 20% to 12%. Due to throughput changes and operator interventions, the process presents typical stationary and nonstationary mixture characteristics. The data used here are collected from one cigarette factory that are different from those in Ref. 42 which have been consciously selected to make sure it was stationary process. Twenty-three measurement variables are used for modeling as described in Table 3 with sampling interval of 10 seconds. First, during offline modelling stage, 726 normal samples are used and the length of moving window is set to be 100 samples from which the control limit is properly set with 0.99 as the significance level. The decomposition result is summarized in Table 4, where the process variations are quantitatively evaluated for each subspace with the corresponding model order. It is found that 70.65% are linear stationary by the second score through the fourth score which will be further analyzed by nonlinear analysis. It is found that the first two KPCA scores are deemed to represent nonstationary variations, which cover 14.80% of the total process variations, while 14.55% of the variations can be transformed into stationary components by nonlinear analysis. Here to prove the superiority of the proposed distribution-based diagnosis method, one normal case and two small disturbances are considered. The monitoring results



are summarized in Table 5 for three methods evaluated by false alarm ratio (FAR) for normal case and detection time delay (DTD) for fault cases. For the normal case, fewer false alarms are observed for the proposed method in comparison with the other two methods, the conventional PCA and DISSIM methods. For Fault #1, the Region-1 vapour pressure is slowly increased which are imposed from the 110th sample by the operator which will result in changes of other variables because of their close correlations. It reveals different variable correlations from those under normal conditions. Despite the change of Region-1 vapour pressure, the related variables stay within the normal region for a certain time because of their slow changes. For Fault #2, the initial moisture content (Variable #2) is slightly larger than the usual by 1%. It can not be noticed until the outputs (Variables #20 and #22) are measured which slowly increase. Due to the changes of KLD dried moisture content and output cooling moisture content, both the set points of Region-1 and Region-2 wall temperatures increase resulting from the use of close-loop control to remove more moisture to make the output normal. For both faults, the proposed method issues alarms with a smaller time delay in comparison with the other two methods owing to the separation of different process characteristics as well as their close description and supervision. 4. Conclusions An incipient fault detection method has been proposed for complex industrial processes with both stationary and nonstationary characteristics. Different from the other fault detection methods, it can probe into the underlying stationary and nonstationary process characteristics more comprehensively. To achieve this purpose, a triple subspace decomposition strategy is developed to isolate the linear stationary process characteristics, the nonlinear stationary process characteristics and the final


Page 26 of 49

Page 27 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


nonstationary process characteristics. The triple subspace distribution monitoring is conducted to quantitatively evaluate the changes of linear and nonlinear stationary and nonstationary distribution structures respectively. By comprehensive subspace decomposition, it can more closely describe and supervise the distribution changes from different aspects. Both applications to real industrial process data demonstrate the feasibility and effectiveness of the proposed method. Therefore, a fault detection system by comprehensive subspace decomposition will be promising for identification of incipient disturbances in complex industrial processes with stationary and nonstationary hybrid characteristics. Acknowledgments This work was supported by the NSFC-Zhejiang Joint Fund for the Integration of Industrialization and Informatization under Grant U1709211 and the National Natural Science Foundation of China under Grant 61433005.



References [1] Zhao, C.H.; Gao, F.R. Online fault prognosis with relative deviation analysis and vector autoregressive modeling. Chem. Eng. Sci. 2015, 138, 531-543. [2] .Cao, Z.X.; Lu, J.Y.; Zhang, R.D.; Gao, F.R. Online average-based system modelling method for batch process. Comput. Chem. Eng. 2018, 108, 128-138. [3] Cao, Z.X.; Yang, Y.; Lu, J.Y.; Gao, F.R. Constrained two dimensional recursive least squares model identification for batch processes. J. Process Control 2014, 24, 871-879. [4] Zhao, C.H.; Sun, Y.X. Multispace total projection to latent structures and its application to online process monitoring. IEEE Trans. Control Syst. Technol. 2014, 22, 868-883. [5] Li, W.Q.; Zhao, C.H.; Gao, F.R. Linearity evaluation and variable subset partition based hierarchical process modeling and monitoring. IEEE Trans. Ind. Electron. 2018, 65, 2683-2692. [6] Cao, Z.X.; Yang, Y.; Yi, H.; Gao, F.R. Priori knowledge-based online batch-to-batch identification in a closed loop and an application to injection molding. Ind. Eng. Chem. Res. 2016, 55, 8818-8829. [7] Zhao, C.H.; Gao, F.R. Fault-relevant principal component analysis (FPCA) method for multivariate statistical modeling and process monitoring. Chemom. Intell. Lab. Syst. 2014, 133, 1-16. [8] Zhao, C.H.; Gao, F.R. Critical-to-fault-degradation variable analysis and direction extraction for online fault prognostic. IEEE Trans. Control Syst. Technol. 2017, 25, 842-854. [9] Wold, S.; Esbensen, K.; Geladi, P. Principal component analysis. Chemom. Intell. Lab. Syst. 1987, 2, 37-52.


Page 28 of 49

Page 29 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


[10] Jackson, J.E. A user's guide to principal components; John Wiley & Sons: Canada, 2005. [11] Dayal, B.S.; MacGregor, J.F. Improved PLS algorithms. J. Chemom. 1997, 11, 73-85. [12] Burnham, A.J.; Viveros, R.; MacGregor, J.F. Frameworks for latent variable multivariate regression. J. Chemom. 1996, 10, 31-45. [13] Zhao, C.H.; Fu, Y.J. Statistical analysis based online sensor failure detection for continuous glucose monitoring in Type I diabetes. Chemom. Intell. Lab. Syst. 2015, 144, 128-137. [14] Rato, T.J.; Reis, M.S. Defining the structure of DPCA models and its impact on process monitoring and prediction activities. Chemom. Intell. Lab. Syst. 2013, 125, 74-86. [15] Zhao, C.H.; Gao, F.R. A sparse dissimilarity analysis algorithm for incipient fault isolation with no priori fault information. Control Eng. Pract. 2017, 65, 70-82. [16] Zhao, C.H.; Gao, F.R. A nested-loop fisher discriminant analysis algorithm. Chemom. Intell. Lab. Syst. 2015, 146, 396-406. [17] Kano, M.; Hasebe, S.; Hashimoto, I.; Ohno, H. A new multivariate statistical process monitoring method using principal component analysis. Comput. Chem. Eng. 2001, 25, 1103-1113. [18] Zhao, C.H.; Gao, F.R. Fault subspace selection approach combined with analysis of relative changes for reconstruction modeling and multifault diagnosis. IEEE Trans. Control Syst. Technol. 2016, 24, 928-939 [19] Basseville, M. On-board component fault detection and isolation using the statistical local approach. Automatica, 1998, 34, 1391-1415.



[20] Kano, M.; Hasebe, S.; Hashimoto, I.; Ohno, H. Statistical process monitoring based on dissimilarity of process data. AIChE J. 2002, 48, 1231-1240. [21] Fukunaga, K.; Koontz, W.L. Application of the Karhunen-Loeve expansion to feature selection and ordering. IEEE Trans. Comput. 1970, 19, 311-318. [22] Zhao, C.H.; Wang, F.L.; Jia, M.X. Dissimilarity analysis based batch process monitoring using moving windows. AIChE J. 2007, 53, 1267-1277. [23] Zhao, C.H.; Wang, F.L.; Zhang, Y.W. Nonlinear process monitoring based on kernel dissimilarity analysis. Control Eng. Pract. 2009, 17, 221-230. [24] Brockwell, P.J.; Davis, R.A. Time series: theory and methods; Springer Science Business Media: New York, 2006. [25] Zhao, C.H.; Huang, B. A full-condition monitoring method for nonstationary dynamic chemical processes with cointegration and slow feature analysis. AIChE J. DOI: 10.1002/aic.16048. [26] Xian, X.C.; Wang, A.D.; Liu, K.B. A nonparametric adaptive sampling strategy for online monitoring of big data streams. Technometrics 2018, 60, 14-25. [27] Box, G.E.P.; Jenkins, G.M.; Reinsel, G.C. Time series analysis-forecasting and control; Prentice-Hall: Englewood Cliffs, 1994. [28] Asteriou, D.; Hall, S.G. Applied Econometrics; Palgrave MacMillan: New York, 2011. [29] Harrisa, T.J.; Seppalaa, C.T.; Desborough, L.D. A review of performance monitoring and assessment techniques for univariate and multivariate control systems. J. Process Control. 1999, 9, 1-17. [30] Engle, R.F.; Granger, C.W.J. Cointegration and error-correction: representation, estimation and testing. Econometrica 1987, 55, 251-276.


Page 30 of 49

Page 31 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


[31] Sun, H.; Zhang, S.M.; Zhao, C.H.; Gao, F.R. A sparse reconstruction strategy for online fault diagnosis in nonstationary processes with no priori fault information. Ind. Eng. Chem. Res. 2017, 56, 6993-7008. [32] Dickey, D.A.; Fuller, W.A. Distributions of the estimators for auto-regressive time series with a unit root. J. Am. Stat. Assoc. 1979, 74, 427-431. [33] Lee, J.M.; Yoo, C.K.; Choi, S.W.; Vanrolleghem, P.A.; Lee, I.B. Nonlinear process monitoring using kernel principal component analysis. Chem. Eng. Sci. 2004, 59, 223-234. [34] Choi, S.W.; Lee, C.; Lee, J.M.; Park, J.H.; Lee, I.B. Fault detection and identification of nonlinear processes based on kernel PCA. Chemom. Intell. Lab. Syst. 2005, 75, 55-67. [35] Haykin, S. Neural Networks; Prentice-Hall: Englewood Cliffs, 1999. [36] Christianini, N.; Shawe-Taylor, J. An introduction to support vector machines and other kernel-based learning methods; Cambridge University Press: UK, 2000. [37] Silverman, B.W. Density estimation for statistics and data analysis. Chapman & Hall: UK, 1986. [38] Zhao, C.H.; Gao, F.R.; Wang, F.L. Nonlinear batch process monitoring using phase-based kernel independent component analysis-principal component analysis (KICA-PCA). Ind. Eng. Chem. Res. 2009, 48, 9163-9174. [39] Zhang, S.M.; Zhao, C.H. Stationarity test and Bayesian monitoring strategy for fault detection in nonlinear multimode processes. Chemom. Intell. Lab. Syst. 2017, 168, 45-61. [40] Yu, W.K.; Zhao, C.H. Sparse exponential discriminant analysis and its application to fault diagnosis. IEEE Trans. Ind. Electron. 2018, 65, 5931-5940.



[41] Yu, W.K.; Zhao, C.H. Online fault diagnosis in industrial processes using multi-model exponential discriminant analysis algorithm. IEEE Trans. Control Syst. Technol. DOI: 10.1109/TCST.2017.2789188. [42] Zhao, C.H.; Wang, W.; Gao, F.R. Probabilistic fault diagnosis based on Monte Carlo and NeLFDA for industrial processes. Ind. Eng. Chem. Res. 2016, 55, 12896-12908.


Page 32 of 49

Page 33 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


List of Figure Captions Figure 1 The illustration of three subspaces separated from the process variables Figure 2 An algorithmic flowchart of the proposed modelling method Figure 3 Systematic configuration for 1000 MW ultra-supercritical thermal power unit Figure 4 Illustration of the separated stationary part and nonstationary part in (a) linear subspace and (b) nonlinear subspace for normal case Figure 5 DISSIM monitoring statistics for the separated stationary part and nonstationary part in (a) linear subspace and (b) nonlinear subspace for normal case Figure 6 Monitoring results for Fault #1 using the proposed method in three different subspaces Figure 7 Monitoring results for Fault #2 using the proposed method in three different subspaces Figure 8 Monitoring results for Fault #2 using (a) PCA method and (b) the conventional DISSIM method Figure 9 Diagram of (a) SIROX and (b) KLD which are two major operation machines in cigarette production process



LSS NSS NLSS

Figure 1 The illustration of three subspaces separated from the process variables (LSS: linear stationary subspace; NLSS: nonlinear stationary subspace; NSS: nonstationary subspace)


Page 34 of 49

Page 35 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Normal measurement data after data normalization

Apply PCA to get the PCA scores T along variation directions

Apply the ADF test to determine the stationarity of the PCA scores

Separate the stationary scores Ts from the nonstationary scores Tns

Separate the original measurement data space into % and nonstationary subspace X % stationary subspace X s ns

% Perform KPCA on X ns to extract nonlinear PCA scores TΦ

Apply the ADF test to determine the stationarity of the KPCA scores

Separate the stationary scores TΦ ,s from the nonstationary scores TΦ ,ns

Separate the original measurement data space % ) and into stationary subspace Φ% s ( X ns % ) nonstationary subspace Φ% ns ( X ns % , Φ% ( X % ) and Perform DISSIM in X s s ns % ) respectively for monitoring model Φ% ( X ns

ns

development

Figure 2 An algorithmic flowchart of the proposed modelling method


Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47

Figure 3 Systematic configuration for 1000 MW ultra-supercritical thermal power unit


Page 36 of 49

Page 37 of 49

29

28

Condensor temperature (°C)

27

26

25

24

23 Linear stationary part Linear nonstationary part

22

21

20

0

100

200

300

400

500 Samples

600

700

800

900

1000

(a)

0.8

Nonlinear stationary part Nonlinear nonstationary part

0.6

KPCA scores

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.4

0.2

0

-0.2

0

100

200

300

400

500 Samples

600

700

800

900

1000

(b) Figure 4 Illustration of the separated stationary part and nonstationary part in (a) linear subspace and (b) nonlinear subspace for normal case



0.25 for linear stationary for linear nonstationary

D2 in linear subspace

0.2

0.15

0.1

0.05

0

0

100

200

300 Samples

400

500

600

500

600

(a)

0.75 for nonlinear stationary for nonlinear nonstationary

0.7 0.65 D2 in nonlinear subspace

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 49

0.6 0.55 0.5 0.45 0.4 0.35

0

100

200

300 Samples

400

(b) Figure 5 DISSIM monitoring statistics for the separated stationary part and nonstationary part in (a) linear subspace and (b) nonlinear subspace for normal case


0.35

0.3 Monitoring statistics Confidence limit

0.25

0.2

D2 in nonlinear stationary subspace

0

100

200

300 Samples

400

500

600

500

600

500

600

0.65 0.6 0.55

Monitoring statistics Confidence limit

0.5 0.45 0.4

0

100

200

300 Samples

400

0.9 D2 in nonstationary subspace

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


D2 in linear stationary subspace

Page 39 of 49


0.85 0.8 0.75 0.7 0.65

0

100

200

300 Samples

400

Figure 6 Monitoring results for Fault #1 using the proposed method in three different subspaces


D2 in nonlinear stationary subspace

Page 40 of 49


0.2 0.19 0.18 0.17 0.16

0

100

200

300 Samples

400

500

600


0.46 0.44 0.42 0.4 0.38

0

100

200

300 Samples

400

500

600

300 Samples

400

500

600

0.8 D2 in nonstationary subspace

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

D2 in linear stationary subspace


0.6

0.4

0.2

0

Monitoring statistics Confidence limit 0

100

200

Figure 7 Monitoring results for Fault #2 using the proposed method in three different subspaces


Page 41 of 49

300 250


T2

200 150 100 50 0

0

100

200

300

400

500 600 Samples

700

800

900

1000

700

800

900

1000

200 Monitoring statistics Confidence limit

150

SPE

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


100

50

0

0

100

200

300

400

500 600 Samples

(a)



0.225 0.22 Monitoring statistics Confidence limit

0.215 0.21 0.205 D2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 49

0.2 0.195 0.19 0.185 0.18 0.175

0

100

200

300 Samples

400

500

600

(b) Figure 8 Monitoring results for Fault #2 using (a) PCA method and (b) the conventional DISSIM method


Page 43 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(a)

(b) Figure 9 Diagram of (a) SIROX and (b) KLD which are two major operation machines in cigarette production process



Table 1 Definitions of the different monitoring subspaces Subspace Dimension Description Monitoring statistics The subspace

% X s

separated from

% ,X % ) Ds = diss ( X s ,r s,w

Rs X-space that is linear stationary (LSS) The subspace separated from

% ) Φ% s (X ns

RΦ ,s

% ) -space that is Φ( X ns

% ),Φ% (X % )) DΦ , s = diss(Φ% s ,r (X ns s,w ns

nonlinear stationary (NLSS)

The final subspace

% ) Φ% ns (X ns

RΦ ,ns

that is nonstationary

% ),∆Φ% (X % )) DΦ ,ns = diss(∆Φ% ns ,r (X ns ns , w ns

(NSS)


Page 44 of 49

Page 45 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 2 Subspace decomposition result for thermal power plant Subspace

Model order

Linear stationary subspace (LSS)

Nonlinear stationary subspace (NLSS)

Nonstationary subspace (NSS)

155

724

2

19.89

74.51

5.60

Modeling variations (%)



Page 46 of 49

Table 3 Variables in cut-made process of cigarette used for modeling Variable Variable Variable Variable Variable Variable No.

Description

No.

Initial cut 1

tobacco flow

3

Initial moisture content (%) SIROX vapour pressure (bar)

9

temperature (℃)

10

6

volume flow

11

KLD vapour pressure (bar)

vapour pressure

12

Region-1 wall temperature

vapour pressure

KLD hot wind speed (m/s) KLD water

19

removal mass (l/h) KLD dried

20

moisture content (%) KLD dried

21

temperature

rate (m3/h)

(bar)

(℃ )

SIROX vapour

Region-2 wall

Cooling

mass flow rate

14

SIROX vapour diaphragm valve

15

temperature

Region-1

Cooling

condensed

temperature

water

Region-2

exhaust negative

condensed

pressure (ubar)

16

moisture content (%)

(℃ )

KLD moisture

22

(℃ )

temperature

opening (%)

8

18

Region-2 13

temperature (℃ )

(bar)

(kg/h)

7

17

Region-1

SIROX vapour 5

exhaust damper

Description KLD hot wind

opening (%)

SIROX 4

No.

KLD moisture

rate (kg/h) 2

Description

water temperature (℃ )


23

(℃ )

Page 47 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 4 Subspace decomposition result for cut-made process of cigarette Subspace

Model order

Linear stationary subspace (LSS)

Nonlinear stationary subspace (NLSS)

Nonstationary subspace (NSS)

20

724

2

70.65

14.55

14.80

Modeling variations (%)



Table 5 Online monitoring results for cigarette production process and three methods evaluated by false alarm ratio (FAR) for normal case and detection time delay (DTD) The proposed method DISSIM method PCA method FAR(%) DTD(sample) FAR(%) DTD(sample) FAR(%) DTD(sample) Normal 5.74 9.57 9.89 Fault 10 41 57 #1 Fault 20 73 80 #2 Case


Page 48 of 49

Page 49 of 49 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC graphic page


Incipient Fault Detection for Complex Industrial Processes with

Recommend Documents