Process Monitoring Based on Independent Component Analysis

This paper proposes a new monitoring method based on independent component analysis−principal component analysis (ICA−PCA). The Gaussian and ...
0 downloads 0 Views 147KB Size
Download PDF
2054

Ind. Eng. Chem. Res. 2007, 46, 2054-2063

Process Monitoring Based on Independent Component Analysis-Principal Component Analysis (ICA-PCA) and Similarity Factors Zhiqiang Ge and Zhihuan Song* National Laboratory of Industrial Control Technology, Institute of Industrial Process Control, Zhejiang UniVersity, Hangzhou 310027, Zhejiang, People’s Republic of China

Many of the current multivariate statistical process monitoring techniques (such as principal component analysis (PCA) or partial least squares (PLS)) do not utilize the non-Gaussian information of process data. This paper proposes a new monitoring method based on independent component analysis-principal component analysis (ICA-PCA). The Gaussian and non-Gaussian information can be extracted for fault detection and diagnosis. Moreover, a new mixed similarity factor is proposed. This similarity factor is used to identify the fault mode. Because of the non-orthogonal nature of the extracted independent components, a “main angle” is proposed to calculate the ICA-based similarity factor. To handle the cases where two datasets have similar principal components or independent components but the numerical values of the process variables are different, two distance similarity factors are used for complement. A case study of the Tennessee Eastman (TE) benchmark process indicates that the proposed fault detection and mode identification methods are more efficient, compared to the alternative methods. 1. Introduction Because of increasing interest in regard to safety and consistent product quality, the real-time monitoring of the process performance has become a key issue for both productivity and safety enhancement. For the successful operation of any process, it is important to detect process upset, equipment malfunctions, or other special events as early as possible, and then to identify and remove the factors causing those events. In chemical processes, data-based process monitoring methods, which are called statistical process control (SPC), have been widely used.1,2 To monitor multivariable processes, multivariate statistical process control (MSPC) has been developed. Kresta et al.3 demonstrated the usefulness of principal component analysis (PCA) and partial least squares (PLS) as tools of MSPC with their applications to simulated data collected from a fluidized-bed reactor and an extractive distillation column. In the past decade, various extensions of MSPC have been proposed.4-8 To monitor time-varying batch processes, the MSPC functioning was extended to multiway PCA and PLS;9-13 with multi-block PCA and PLS,14 the MSPC functioning was extended to monitor large-scale chemical processes, dynamic PCA (to include process dynamics in a PCA model),15 multiscale PCA based on wavelet analysis for monitoring signals at several different frequency ranges,16 model-based PCA for integrating a process model with PCA,17-19 and a method based on constrained PCA for incorporating external information into a PCA model.20 In addition, the integration of several advanced monitoring methods has been discussed recently.21,22 Many successful applications have shown the practicability of MSPC. However, the previously described MSPC methods do not always function very well. Their achievable performance is limited, because of the assumption that the monitored variables have a Gaussian distribution. In fact, most of the process variables do not strictly form a Gaussian distribution. On the other hand, it is well-known that many of the variables monitored in process systems are not independent; chemical * To whom correspondence should be addressed. Tel.: +86-57187951442-808.Fax: +86-571-87951921.E-mail: [email protected].

processes are usually driven by fewer essential variables, which may not be measured, and the measured process variables may be combinations of these independent latent variables. Independent component analysis (ICA) is an emerging technique for determining several independent variables as linear combinations of measured variables. What distinguishes ICA from PCA is that it searches for components that are both statistically independent and non-Gaussian. PCA can only impose independence up to second-order statistics information, and the direction vectors are constrained to be orthogonal, whereas ICA has no orthogonal constraint and also involves higher-order statistics. ICA is interested in the non-Gaussian information, whereas PCA is looking for the Gaussian information. In a word, ICA may reveal more meaningful information in the non-Gaussian data than PCA. Many applications of ICA have been reported in speech processing, biomedical signal processing, machine vibration analysis, nuclear magnetic resonance spectroscopy, infrared optical source separation, radio communications, etc.23-26 Kano et al.27 developed a unified framework for MSPC, which combined PCA-based SPC and ICA-based SPC. The proposed combined MSPC (CMSPC) could monitor both of the Gaussian and non-Gaussian information of the process; good performance was shown in a multivariate system and a continuously stirred tank reactor (CSTR) process. However, when the non-Gaussian information had been extracted, all the Gaussian variables were chosen for principal components; therefore, noise might be incorporated into the model. Besides, the independent components were monitored separately, so when a large number of independent components are extracted, it is a tedious task for the operator. It is would be difficult to interpret the meanings of independent components and relate them with process operation; therefore, it may be better to incorporate the independent components together for fault detection, as PCA does for the principal components. To our knowledge, the MSPC methods mainly focus on fault detection, but once a fault occurs, we should identify it, isolate it, and then repair it to get the process back to the normal state. There are always many faults in the process, and how these fault modes are identified is the essential step that we should do. Typically, the contribution plots are used to identify the

10.1021/ie061083g CCC: $37.00 © 2007 American Chemical Society Published on Web 03/02/2007

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007 2055

root cause of the abnormal operation. However, they can only identify some simple faults. Dunia et al.28 analyzed the fundamental issues of detectability, reconstructability, and isolatability for multidimensional faults; however, it is based on the assumption that the faults do not transfer to other parts of the process. Zhang et al.29 proposed the characteristic direction method to identify faults; they took the first loading vector as the characteristic direction. Although the first principal component extracted most of the fault information, much of the other information has been omitted. Fortunately, Krzanowski30 developed a method for measuring the similarity of two datasets using a PCA similarity factor, SPCA. Johannesmeyer et al.31 used that PCA similarity factor to identify the similarity between different operation modes. Within that method, all the information that the PCA model had covered was used, so the power of fault identification was greatly improved. However, two datasets may have the same spatial orientation, and they may have similar principal components, but that does not mean that they have same numerical values of the process variables; in a word, they may be located far apart. Fortunately, Singhal and Seborg32-34 have proposed a distance similarity factor to handle this problem; however, the method is based on the assumption that the data formed a Gaussian distribution. In this paper, a new distance similarity factor is proposed, which is a complement to the PCA similarity factor. Besides, the ICA similarity factor and the corresponding distance similarity factor are introduced to gain better results. The independent components have no orthogonal constraint, and these extracted independent components (ICs) have equal roles in the ICA model; therefore, a “main angle” is proposed to calculate similarity between the two ICA models. The motivation of the present work is that the monitoring performance can be improved by focusing on essential variables that drive the process, which are called independent components, and monitoring the non-Gaussian and Gaussian parts separately. By introducing a distance similarity factor and an ICA similarity factor, the fault identification performance could be greatly improved. In section 2, some preliminaries such as ICA and PCA are discussed. In the next section, we propose a new online fault detection method that is based on ICA-PCA. A new fault identification method is proposed in section 4. To investigate the feasibility of the proposed methods, the performance is evaluated and compared with that of the conventional methods on the Tennessee Eastman (TE )benchmark process in section 5. Finally, conclusions are drawn in section 6. 2. Preliminaries 2.1. Independent Component Analysis (ICA). Independent component analysis (ICA) is a statistical and computational technique for revealing hidden factors that underlie sets of random variables, measurements, or signals.23 ICA was originally proposed to solve the blind source separation problem. To introduce the ICA algorithm, it is assumed that l measured variables, x(k) ) [x1(k), x2(k), ..., xl(k)], at sample k can be expressed as linear combinations of r unknown independent components [s1, s2, ..., sr]T (where r e l); the relationship between them is given by

X ) A‚S + E

(1)

where n is the number of measurements, X ) [x1, x2, ..., xn] ∈Rl×n is the data matrix, A ) [a1, a2,, ..., ar] ∈ Rl×r is the mixing matrix, S ) [s1, s2, ..., sr] ∈ Rr×n is the independent component matrix, and E ∈ Rl×n is the residual matrix. The basic problem

of ICA is to estimate the original component S and the mixing matrix A from X; therefore, the objective of ICA is to calculate a separating matrix W so that the components of the reconstructed data matrix S become as independent of each other as possible, given as

Sˆ ) W‚X

(2)

Before applying the ICA algorithm, the data matrix X should be whitened, to eliminate the cross-correlations among random variables. One popular method for whitening is to use the eigenvalue decomposition, considering x(k) with its covariance Rx ) E{x(k)‚x(k)T}, the eigenvalue decomposition of Rx is given by

Rx ) UΛUT

(3)

The whitening transformation is expressed as

z(k) ) Q‚x(k)

(4)

where Q ) Λ-1/2‚UT. One can easily verify that Rz ) E{z(k)‚ z(k)T} is the identity matrix under this transformation. After the whitening transformation, we have

z(k) ) Q‚x(k) ) Q‚A‚s(k) ) B‚s(k)

(5)

where B is an orthogonal matrix, giving the following verification:

E{z(k)‚z(k)T} ) B‚E{s(k)‚s(k)T}‚BT ) B‚BT ) I (6) Therefore, we have reduced the problem of finding an arbitrary full-rank matrix A to the simple problem of finding an orthogonal matrix B. We then can make the following estimate:

sˆ(k) ) BT‚z(k) ) BT‚Q‚x(k)

(7)

From eqs 2 and 7, we can get the relationship between W and B:

W ) BT‚Q

(8)

There are two classic measures of non-Gaussianity: kurtosis and negentropy.23 Kurtosis is the fourth-order cumulant of a random variable; unfortunately, it is sensitive to outliers. On the other hand, negentropy is based on the information-theoretic quantity of differential entropy. We choose negentropy to measure non-Gaussianity in our present study. When we have received the approximate forms for negentropy,23 Hyva¨rinen35 introduced a very simple and efficient fixed-point algorithm for ICA. It is calculated over a sphered zero-mean vector z. This algorithm calculates the independent components one by one, through iterative steps. After all vectors bi (i ) 1, ..., m) were calculated, they made the orthogonal mixing matrix B; we then can obtain sˆ(k) and demixing matrix W form eqs 7 and 8, respectively. Dimension reduction of ICA is based on the idea that these measured variables are the mixture of some IC variables. The performance and interpretation of ICA monitoring each is dependent on the correct choice of the ordering and dimension of the ICA model. Because negentropy is used to measure nonGaussianity, we can select a suitable number of ICs by checking their negentropy values: if the negentropy value of the current IC is zero or approximately zero, it indicates that the nonGaussian information of the process has already been extracted,

2056

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007

and the rest of the process is Gaussian, which can be analyzed by conventional MSPC methods, such as PCA. 2.2. Principal Component Analysis (PCA). Consider an n × m data matrix X, for m process variables with n measurements of each variable. Here, we assume that the data for each variable have been scaled to zero mean and unit variance. The covariance matrix of X is defined as

∑

)

XT‚X n-1

(9)

Let λi (i ) 1, 2, ..., m) be the eigenvalues of the matrix of ∑, which are arranged in descending order to determine the principal components (PCs), and with the principal component loadings pi to be their corresponding eigenvectors. The first k PCs then are selected to build the PCA model, so the data matrix X can be expanded using PC loadings pi, score vectors ti, and a residual matrix E:3 k

X)

ti‚pTi + E ∑ i)1

(10)

For a new observation xnew, the prediction of the PCA model is given by

xˆ new ) tnew‚PT

(11)

where tnew is the score vector (tnew ) xnew‚PT). The resulting residual is defined as

e ) xnew - xˆnew

(12)

For monitoring purposes, two statistical variablessT2 and SPEs can be calculated by t and e, respectively:

T2 ) tnew‚Λ‚tTnew SPE ) enew‚eTnew

k

E)

(14)

3. ICA-PCA Based Fault Detection Method Given the data matrix X, suppose that r independent components are extracted, S ) [s1, s2, ..., sr] ∈ Rr×n. To monitor the non-Gaussian part of the process, a new statistic variable was defined:24,25

(15)

ti‚pTi + F ∑ i)1

(16)

where F is the residual resulting from the PCA model. Here, we define the limits of T2 and SPE statistics to monitor the remaining Gaussian part of the process: k

T ) 2

ti2

∑ i)1 λ

[

]

k(n - 1)

e i

n-k

Fk,(n-k),R

SPE ) f‚fT ) e(I - PPT)eT e SPER

[

(17)

]

cRx2θ2h02 θ2h0(h0 - 1) + SPER ) θ1‚ 1 + θ1 θ2 1

(18)

1/h0

(19)

m where k is the number of PCs, θi ) ∑j)k+1 λij (for i ) 1, 2, 3), 2 h0 ) 1 - (2θ1θ3/3θ2 ), R is the significance level, and cR is the normal deviate corresponding to the upper 1 - R percentile. In PCA monitoring, the confidence limits are based on a specified distribution (shown in eqs 17-19) that is based on the assumption that the latent variables follow a Gaussian distribution. However, in ICA monitoring, the ICs do not conform to a specific distribution; hence, the confidence limit of the I2 statistic cannot be determined directly from a particular approximate distribution. An alternative approach to defining the nominal operating region of the I2 statistic is to use kernel density estimation (KDE).36,37 Here, we only need a univariate kernel estimator, which is defined by

(13)

PCA has been widely reported; therefore, only a brief summary is presented here. 2.3. Problem Statement. Although many industrial applications have shown success in regard to the use of conventional MSPC methods, they are performed under the assumption that the monitored variables have a Gaussian distribution. However, most of the process variables are non-Gaussian, and they may be combinations of fewer essential variables that are not measured. ICA can extract these essential variables efficiently; however, after the non-Gaussian information is extracted from the process, the remaining Gaussian part should also be analyzed. In this paper, PCA is chosen to analyze the Gaussian part of the process, whereas other conventional MSPC methods also can be used. When a fault is detected, a new fault identification method should be developed, which can identify non-Gaussian information between two operation modes. In the next two sections, a new fault detection method and a novel mode identification strategy are proposed.

I2 ) ST‚S

After the non-Gaussian information has been extracted, the residual matrix E is obtained. We can use PCA to analyze it, expanding E as shown below:

ˆf (x,H) )

1

n

K[H-1/2(x - xi)] ∑ n i)1

(20)

where x is the data point under consideration, xi an observation value from the data set, H the bandwidth matrix, and K a kernel function, which satisfies the following conditions:

K(x) g 0

∫R

p

K(x) dx ) 1

(21)

There are several kernel functions. However, in practice, the Gaussian kernel function is the most commonly used. Many measures have been proposed for the estimation of H, the window width, or the smoothing parameter. The problem of choosing how much to smooth is crucially important in density estimation. Several methods exist that automatically choose an optimal value of H: one efficient method is MISE. There are three choices for H: (1) For a problem with p variables, a full symmetrical bandwidth matrix

[

h112 h 2 H ) 21 l hp12

h122 ‚‚‚ h1p2 h222 ‚‚‚ h2p2 l l hp22 ‚‚‚ hpp2

]

which is a positive-definite matrix with p(p + 1)/2 parameters, in which h2ik ) h2ki;

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007 2057

(2) A diagonal matrix with only p parameters,

H ) diag(h21, h22, ..., h2p) (3) A diagonal matrix with one parameter,

H ) h2I where I is an identity matrix. Approach 1 is the most accurate but is unrealistic, in terms of computational load and time. Approach 3 is the simplest, and it is the method used in the case of univariate data. It can be adopted with slight modification; however, it can cause problems in some situations, because of the loss of accuracy by forcing the bandwidths in all dimensions to be the same. This can introduce significant error in the shape of the density function. Hence, approach 2 seems to be the appropriate choice, as a compromise; however, the computational load is still unacceptable if the dimensionality of the problem is high and the size of the sample is large. Here, we choose the Gaussian kernel function for K, and H ) diag(h21, h22, ..., h2p); thus, we can decide the 99% or 95% confidence limit for I2 monitoring. Although the confidence limits of T2 and SPE can be decided by eqs 17-19, they can also be decided by KDE. In some sense, it may be more efficient, but also costly for more computation. The implementation of the monitoring method as mentioned previously consists of two procedures: offline modeling and online monitoring. In the offline modeling procedure, a ICAPCA monitoring model is developed under the normal operating condition (NOC). The fault detection and isolation are executed using this monitoring model in the online monitoring procedure. The detailed algorithm flow of the monitoring procedure is summarized below: (1) Offline modeling procedure Step 1: Acquire an operating data set X during the normal process; Step 2: Normalize and whiten the data matrix X; Step 3: Perform the ICA algorithm to obtain the matrices W, B, and Sˆ , so that Sˆ has great non-Gaussianity. We can also obtain A ) (QTQ)-1QTB and W ) BTQ; Step 4: For each sample, calculate its I2 value,

I2(k) ) sˆ(k)T‚sˆ(k)

(k ) 1, 2, ..., n)

using the KDE to decide its 99% or 95% confidence limit (Ilim2); Step 5: Perform the PCA algorithm to obtain the score matrix T and the load matrix P; then, decide the confidence limits Tlim2 and SPElim for T2 and SPE. (2) Online monitoring procedure Step 1: For new sample data xnew, use the same scaling method as that used in the modeling steps; Step 2: Calculate the ICs of the new sample data (xnew: sˆnew ) W‚xnew); Step 3: Calculate the I2 statistic value for this new sample data:

Inew2 ) sˆTnew‚sˆnew

Figure 1. Fault detection strategy for the proposed independent component analysis-principal component analysis (ICA-PCA) method.

Step 5: Calculate Tnew2 and SPEnew, using eqs 17 and 18; Step 6: If the statistics reject their corresponding limits, some type of fault has been detected; otherwise, go to step 1 and continue monitoring. The online process fault detection scheme of the ICA-PCAbased method is shown in Figure 1. 4. Similarity Factors for Mode Identification Krzanowski30 developed a method for measuring the similarity of two datasets using a PCA similarity factor, SPCA. Let us consider two datasets, D1 and D2; they have the same variables but not necessarily the same number of measurements. Suppose we have chosen k PCs for the two PCA models; those are their corresponding respective subspaces. The PCA similarity factor is a measure of the similarity between datasets D1 and D2, and, to be exact, it lies between the two reduced subspaces. If we let L and M represent the two subspaces, the PCA similarity factor can be calculated from the angles between the PCs:

SPCA )

1 k

k

k

cos2 θij ∑ ∑ i)1 j)1

(22)

where θij is the angle between the ith PC of dataset D1 and the jth PC of dataset D2, which is defined as

cos θij )

LTi Mj ||Li||2||Mj||2

Here, ||Li||2 ) ||Mj||2 ) 1 and 0 eθ eπ/2. Besides, SPCA can also be expressed by subspaces L and M:

SPCA )

trace(LTMMTL) k

(23)

Step 4: The remaining portion is given by

enew ) xnew - A‚sˆnew Then calculate the score vector, tnew ) enew‚P, and we also obtain the residual, fnew ) enew - eˆ new;

If the value of SPCA exceeds a specified cutoff value, the two datasets D1 and D2 are considered to be similar. However, the PCA similarity factor is just a measure of the similarity between the Gaussian part of datasets D1 and D2. For the similarity measure of the non-Gaussian part, we establish

2058

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007

the ICA similarity factor SICA for the corresponding subspaces A and B. What distinguishes SICA from SPCA is the fact that the ICs have no orthogonal constraint, and these extracted ICs have equal roles in the ICA model. Therefore, we propose new angles between the two subspaces for ICA, which we call “main angles” in this paper, and these new “main angles” are defined as follows:

(

)

aT1 b1 aTb cos θ1 ) max max ) (24) ||a1||2||b1||2 a∈A b∈B ||a||2||b||2 where θ1 is defined as the first main angle; the vectors a1 and b1 then are excluded from the subspaces A and B, and, thus, a1 is removed from subspace A and b1 is removed from subspace B. The new subspaces as are denoted as A1 and B1. Similarly, the second main angle θ2 is defined as

(

)

aT2 b2 aTb ) (25) cos θ2 ) max max ||a2||2||b2||2 a∈A1 b∈B1 ||a||2||b||2 Therefore, the ith main angle θi is defined as

cos θ1 ) max max a∈Ai-1 b∈Bi-1

(

)

T

a i bi aTb ) ||a||2||b||2 ||a2||i||bi||2 (26)

After all the main angles have been calculated, the vectors in the two subspaces are rearranged as A′ and B′. In these subspaces, the first vector of A′ and B′ corresponds to the first main angle, the second pair of vectors correspond to the second main angle, and so on. We assume to choose r ICs for the two ICA models. Therefore, the ICA similarity factor can be calculated from the main angles between the two subspaces:

SICA )

1

r

∑ cos2 θi r i)1

(27)

where θi is the ith main angle between the two subspaces A and B. However, there exist such cases that the two datasets have similar spatial orientation but are located apart. In other words, the two datasets have similar PCs or ICs, but the numerical values of the process variables are very different. To handle theses cases, two distance similarity factors, SPC_dist and SIC_dist, are proposed: n

SPC_dist )

xe-d (T ∑ n i)1

SIC_dist )

xe-d (I ∑ n i)1

1

1

2

2(i),T 2(i)) M

L

(28)

n

2

2(i),I 2(i)) B

A

(29)

where n is the total measurements, TL2(i) and TM2(i) are then T2 values of the ith measurement for the corresponding subspaces, and d(‚‚‚) represents their distance. Similarly, IA2(i) and IB2(i) are the I2 values of the ith measurement for the corresponding subspaces, and d(‚‚‚) represents the distance. For simplicity, we can use the mean value of T2 and I2 to calculate the two distance similarity factors:

SPC_dist ) xe-d (TL (mean),TM (mean))

(30)

SIC_dist ) xe-d (IA (mean),IB (mean))

(31)

2

2

2

2

2

2

In summary, there are four similarity factors, and they can be used independently. If more than one similarity factor is used, a weighted average of the similarity factors can be calculated, as

Smix ) a‚SPCA + b‚SICA + c‚SPC_dist + d‚SIC_dist

(32)

where a, b, c, and d are weighting factors (0 e a, b, c, d e 1 and a + b + c + d ) 1). In our experience, selecting a and b values that are bigger than c and d is more suitable for mode identification, because SPCA and SICA are the main reasons why the two modes are different from each other, and SPC_dist and SIC_dist are used as complements. However, when different magnitudes of faults occur in one operation mode, SPC_dist and SIC_dist become important, and SPCA and SICA may be used as complements. Besides, the method of how to choose the values of the weighting factors also is dependent on the features of the process. When the process is extremely non-Gaussian, the weighting factors of ICA-based similarity factors may be chosen to be larger than those of PCA-based similarity factors. On the other hand, when the process is extremely Gaussian, the weighting factors of PCA-based similarity factors may be chosen to be larger than those of ICA-based similarity factors. Before identifying a detected fault, a mode information storeroom should be constructed. This storeroom contains all the operation mode (normal and fault) information of the process, including parameters and matrices for mode identification, such as A, P, I2, and T2. The steps that are required to build the mode information storeroom are as follows: Step 1: Acquire data set X for all operation modes of the process; Step 2: For each operation mode, calculate its mixing matrix A, load matrix P, and I2 and T2 values for every sample point; Step 3: Calculate the similarity factor values for every pair of operation modes, including SPCA, SICA, SPC_dist, SIC_dist, and Smix; Step 4: Determine cut-off values for these similarity factors; and Step 5: Preserve all the parameters and matrices in the mode information storeroom. The flow of building the modes information storeroom is shown in Figure 2. When a fault is detected, mode identification is performed as follows: Step 1: Acquire data set Xfault for the fault operation mode of the process; Step 2: Calculate its mixing matrix Afault, load matrix Pfault, Ifault2, and Tfault2 values for every sample point; Step 3: Calculate the similarity factor values, which lie between the fault operation mode and every mode in the storeroom; Step 4: If the similarity factor value exceeds the cut-off value, the detected fault mode is considered to be the corresponding operation mode in the storeroom; Step 5: If none of the similarity factor values exceeds the cut-off value, the detected fault mode is considered to be a new operation mode, and the storeroom should be updated. The mode identification scheme of the similarity factor based method is shown in Figure 3.

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007 2059

Figure 2. Development of the mode information storeroom.

Figure 3. Mode identification strategy for the proposed method.

Figure 4. Control system of the Tennessee Eastman (TE) process. Table 1. Monitoring Variables in the TE Process No.1

measured variable

No.1

measured variable

No.1

manipulated variable

1 2 3 4 5 6 7 8 9 10 11

A feed D feed E feed total feed recycle flow reactor feed rate reactor pressure reactor level reactor temperature purge rate product separator temperature

12 13 14 15 16 17 18 19 20 21 22

product separator level product separator pressure product separator underflow stripper level stripper pressure stripper underflow stripper temperature stripper steam flow compressor work reactor cooling water outlet temperature separator cooling water outlet temperature

23 24 25 26 27 28 29 30 31 32 33

D feed flow valve E feed flow valve A feed flow valve total feed flow valve compressor recycle valve purge valve separator pot liquid flow valve stripper liquid product flow valve stripper steam valve reactor cooling water flow condenser cooling water flow

5. Simulation and Results In this section, the fault detection and mode identification performance of the proposed methods are compared with that of conventional methods through the Tennessee Eastman (TE) process.38 5.1. Tennessee Eastman (TE) Process. As a benchmark simulation, the TE process has been widely used to test the performance of various monitoring approaches.39-41 This process consists of five major unit operations: a reactor, a condenser, a compressor, a separator, and a stripper. The control structure is shown schematically in Figure 4, which is the second structure listed in the work by Lyman and Georgakis.42 The TE process has 41 measured variables (22 continuous process measurements and 19 composition measurements) and 12 manipulated variables; a set of 21 programmed faults are introduced to the process. The details of the process description are explained well in a book authored by Chiang et al.41

In this paper, the variables that have been selected for monitoring are listed in Table 1. There are 33 variables in the table: 22 continuous process measurements and 11 manipulated variables. The agitation speed is not included, because it is not manipulated. Moreover, we also exclude 19 composition measurements, because they are difficult to measure online in real processes. The simulation data that we have collected were separated into two parts: training datasets and testing datasets. They consisted of 960 observations for each mode (1 normal and 21 fault), respectively, and their sampling interval was 3 min. All faults in the testing datasets were introduced in the process at sample 160, which are tabulated in Table 2. Faults 1-7 are step changes of the process variables; faults 8-12 are random changes of the variables; fault 13 is a slow shift of the reaction kinetics; faults 14, 15, and 21 are related to valve sticking; and faults 16-20 are types of unknown faults. Among these faults, some faults are easy to be detected, because they

2060

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007

Table 2. Process Disturbances fault number 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

process variable

type

A/C feed ratio, B composition constant (stream 4) B composition, A/C ratio constant (stream 4) D feed temperature (stream 2) reactor cooling water inlet temperature condenser cooling water inlet temperature A feed loss (stream 1) C header pressure loss-reduced availability (stream 4) A, B, C feed composition (stream 4) D feed temperature (stream 2) C feed temperature (stream 4) reactor cooling water inlet temperature condenser cooling water inlet temperature reaction kinetics reactor cooling water valve condenser cooling water valve unknown unknown unknown unknown unknown valve position constant (stream 4)

step step step step step step step random variation random variation random variation random variation random variation slow drift sticking sticking unknown unknown unknown unknown unknown constant position

Table 3. Monitoring Results of the TE Process PCA

ICA

ICA-PCA

fault mode

T2

SPE

I2

SPE

I2

T2

SPE

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

0.007 0.022 1.000 0.975 0.750 0.010 0.090 0.035 1.000 0.710 0.811 0.031 0.065 0.162 0.999 0.812 0.241 0.112 0.998 0.712 0.725

0.004 0.015 0.998 0.045 0.745 0 0 0.025 0.997 0.651 0.367 0.027 0.050 0 0.985 0.755 0.112 0.100 0.762 0.524 0.621

0.006 0.020 0.995 0.035 0 0 0.065 0.030 0.995 0.215 0.756 0.025 0.025 0.045 0.426 0.213 0.211 0.102 0.351 0.426 0.635

0.002 0.010 0.990 0 0 0 0 0.021 0.996 0.200 0.182 0 0 0 0.328 0.125 0.196 0.085 0.322 0.358 0.528

0.005 0.018 0.996 0.032 0 0 0.050 0.028 0.994 0.212 0.712 0.024 0.016 0.042 0.412 0.185 0.206 0.100 0.346 0.412 0.622

0.003 0.015 0.993 0.020 0 0 0.035 0.022 0.993 0.156 0.526 0.005 0.011 0.013 0.268 0.146 0.107 0.085 0.238 0.247 0.549

0 0.008 0.985 0 0 0 0 0.011 0.991 0.162 0.182 0 0 0 0.196 0.125 0.195 0.098 0.208 0.306 0.526

greatly affect the process and change the relationships between the process variables. However, there are also faults that are difficult to detect (e.g., faults 3, 9, and 15), because they are very small and have little influence on the process. 5.2. Online Fault Detection of the TE Process. So far, we have made 22 datasets, corresponding to the 22 different modes (1 normal and 21 fault). Before the application of PCA, ICA, and ICA-PCA, all the data were auto-scaled. Nine ICs and 15 PCs were selected for ICA and PCA, respectively, by crossvalidation. In the study of the ICA-PCA method, we select the same number of ICs (9) and PCs (15), to compare the monitoring performance to ICA and PCA. The 99% confidence limits of all the statistic variables were determined by KDE. For each statistic, the detection rates for all 21 fault modes were calculated and tabulated in Table 3. The minimum missing detection rate achieved for each mode is marked with a bold number, except the modes of faults 3 and 9 (these latter faults are quite small and have almost no effect on the overall process, so they were excluded in our research). As shown in

Table 3, ICA outperforms PCA for most fault modes, and, particularly, the missing detection rates of ICA for the modes of faults 5, 10, and 16 are much lower than that of PCA, which indicates that ICA can detect small events that are difficult to detect using PCA. However, the missing detection rates of ICAPCA for most modes are even lower than those of ICA, because ICA-PCA not only monitors the non-Gaussian information of the process by ICA, but the Gaussian information also is monitored by PCA. As a result of the simulation, the best monitoring performance is observed in the case of ICA-PCA. In most modes, the missing detection rate of ICA-PCA is the lowest. The monitoring results of fault 5 are shown in Figure 5. The condenser cooling water inlet temperature is step-changed in the mode of fault 5. In that mode, the flow rate of the outlet stream from the condenser to the separator also increases, which causes an increase in temperature in the separator, thus affecting the separator cooling water outlet temperature. Because we have established a control structure for this process, the control loops will act to compensate for the change, and the temperature in the separator will return to its setpoint. A period of ∼10 h is required to attain the steady state again. As shown in Figure 5a, PCA have detected the fault at sample 160 (approximately). However, the fault cannot be detected after approximately sample 350, because most variables returned to their setpoints. In this case, if a process operator judges the status of the process based on PCA, he would probably conclude that a fault entered the process and then corrected itself within ∼10 h. However, as a matter of fact, the condenser cooling water inlet temperature is still higher than that under the normal condition after sample 350, and the condenser cooling water flow rate is continuously manipulated while most variables return to their setpoints. All of this information indicates that a fault still remains in the process. Because PCA is based on only second-order statistics information, the effect of the condenser cooling water flow may be neglected, compared to the effect of other variables; thus, the problem occurs. On the other hand, as shown in Figures 5b and 5c, all the ICA and ICA-PCA statistical variables remained above their confidence limit, which indicates that a fault remains in the process. What distinguishes ICA from PCA is that ICA involves higher-order statistical information, which is nonGaussian and, thus, can correctly reflect the effect of the condenser cooling water flow rate. Such a persistent fault detection statistic will continue to inform the operator that a process abnormality remains in the process, although all the process variables will seem to have returned to their normal values through control loops.42 In the method of ICA-PCA, ICA has extracted essential components that underlie the process; thus, the following PCA procedure can reflect the effect of the condenser cooling water flow rate. Therefore, ICA and ICA-PCA may provide the process operator with more-reliable information. In addition, comparing the result of ICA-PCA with that of ICA, we can conclude that the effect of the fault is magnified more in the former statistical variables. 5.3. Fault Mode Identification. When a fault has been detected, the next procedure that must be performed is its identification. As we know, the TE simulation process has 21 fault modes; some faults are easy to identify (such as faults 1 and 6), and others are difficult to identify (including faults 3, 9, and 15). In this section, the performances of the different mode identification methods for the TE process are compared. These approaches all use some form of fault signatures that are developed from prior faults: when a new fault occurs, its

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007 2061 Table 4. Results of the Three Methods for Mode Identification

Figure 5. Monitoring results for the mode of fault 5 of the TE process: (a) result based on PCA, (b) result based on ICA, and (c) result based on ICA-PCA.

signature is compared with those in the fault bank to identify the most likely case or mode. Most importantly, as mentioned in section 4, we should learn fault signatures from the training datasets and store them in the database. Therefore, the training datasets are used to develop the ICA-PCA models for all 22 modes, and these models are all stored in the database. Moreover, all the cut-off values of the similarity factors should be determined. In the present work, we use 25 test datasets for each mode to calculate the reliability. In the simulation, all 550 test datasets are chosen randomly, and the results of mode identification using the three different methods are tabulated in Table 4, where N is the number of

mode ID

N

0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25 25

average

25

Ncs 25 23 22 0 6 7 24 17 10 1 8 6 9 8 5 0 1 7 11 0 3 22 9.77

ηcs (%) 100 92 88 0 24 28 96 68 40 4 32 24 36 32 20 0 4 28 44 0 12 88 39.08

Npca 25 25 21 2 9 10 25 19 12 12 24 22 11 10 16 3 18 24 19 23 16 23 16.77

ηpca (%)

Nmix

ηmix (%)

100 100 84 8 36 40 100 76 48 48 96 88 44 40 64 12 72 96 76 92 64 92

25 25 23 3 14 14 25 22 21 17 25 25 14 15 23 5 22 21 24 25 23 24

100 100 92 12 56 56 100 88 84 68 100 100 56 60 92 20 88 84 96 100 92 96

67.08

19.68

78.72

test datasets for every operation mode, Ncs is the number of test datasets that are identified by the characteristic subspace method, ηcs is the mode identification efficiency of the characteristic subspace method, the terms Npca and ηpca correspond to SPCA, and the terms Nmix and ηmix correspond to Smix. Here, we choose the following values: a ) b ) 0.4 and c ) d ) 0.1. As shown in Table 4, the results for operating modes 4, 5, 8, and 14 are greatly improved by Smix, whereas they are difficult to identify by SCS and SPCA, and the results for other operating modes also are, more or less, improved by Smix. However, Smix cannot successfully identify operating modes 3 and 15. Chiang et al.41 also determined that the misclassification rates using dynamica PCA and Fisher discriminant analysis for operating modes 3 and 15 were very high. Raich and C¸ inar43 reported zero success in distinguishing operating mode 3 from the other operating modes using SPCA. The characteristic subspace method shows poor identification accuracy for most of the operating modes, because it only uses the first PC, and much important information is omitted. Identification accuracy is improved using SPCA, which uses information from all of the PCs, but it still gives low identification accuracy for several operating modes, such as modes 3, 4, 5, and 15, because they are very difficult to identify. However, Smix not only uses the Gaussian information (extracted by PCA), but also uses the non-Gaussian information that is extracted by ICA. Moreover, the distance factors are also used as complements; however, unfortunately, high misclassification rates are still shown for operating modes 3 and 15. The detailed results of Smix for the misclassification of all operating modes are presented in Table 5. Because of length limitations of the paper, the results of SPCA and SCS are not presented. The numerical values in the table indicate that many operating modes are misclassified as other operating modes. As seen from Table 5, the numerical values in the second column indicate that many of the operating modes are misclassified as the normal operating mode. We can find that almost every fault mode is misclassified as the normal operating mode. The reason for the misclassifications can be obtained by observing the values of similarity factors for different pairs of operating modes, which are tabulated in Table 6, indicating one identification result for all 22 operating modes. The similarity values of Smix between different operating modes are smaller,

2062

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007

Table 5. Detail Mode Identification Result of the Smix Method ID 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

0 25 0 2 12 5 6 0 3 4 5 0 0 4 3 1 4 3 4 1 0 1 0

1 0 25 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

2 0 0 23 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1

3 0 0 0 3 2 1 0 0 0 3 0 0 2 0 0 5 0 0 0 0 0 0

4 0 0 0 0 14 0 0 0 0 0 0 0 0 0 1 2 0 0 0 0 0 0

5 0 0 0 3 0 14 0 0 0 0 0 0 0 0 0 6 0 0 0 0 0 0

6 0 0 0 0 0 0 25 0 0 0 0 0 0 3 0 0 0 0 0 0 0 0

7 0 0 0 0 0 0 0 22 0 0 0 0 0 0 0 0 0 0 0 0 0 0

8 0 0 0 0 0 0 0 0 21 0 0 0 0 0 0 0 0 0 0 0 0 0

9 0 0 0 2 0 0 0 0 0 17 0 0 4 0 0 1 0 0 0 0 0 0

10 0 0 0 0 0 0 0 0 0 0 25 0 0 0 0 0 0 0 0 0 0 0

11 0 0 0 0 0 0 0 0 0 0 0 25 0 0 0 0 0 0 0 0 0 0

12 0 0 0 0 0 0 0 0 0 0 0 0 14 0 0 0 0 0 0 0 0 0

13 0 0 0 0 0 0 0 0 0 0 0 0 0 15 0 2 0 0 0 0 0 0

14 0 0 0 0 0 0 0 0 0 0 0 0 0 0 23 0 0 0 0 0 0 0

15 0 0 0 3 4 4 0 0 0 0 0 0 1 4 0 5 0 0 0 0 0 0

16 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 22 0 0 0 0 0

17 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 21 0 0 0 0

18 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 24 0 1 0

19 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 25 0 0

20 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 23 0

21 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 24

10 0.76 0.65 0.71 0.95 0.57 0.71 0.53 0.60 0.77 0.90 1.00 0.69 0.70 0.71 0.70 0.64 0.68 0.58 0.58 0.66 0.58 0.59

11 0.66 0.52 0.67 0.87 0.60 0.61 0.48 0.49 0.66 0.79 0.69 1.00 0.72 0.81 0.72 0.74 0.79 0.71 0.60 0.74 0.71 0.66

12 0.66 0.69 0.64 0.82 0.61 0.75 0.59 0.63 0.79 0.77 0.70 0.72 1.00 0.77 0.55 0.63 0.63 0.58 0.65 0.60 0.59 0.58

13 0.69 0.76 0.75 0.81 0.65 0.72 0.59 0.68 0.86 0.84 0.71 0.81 0.77 1.00 0.60 0.64 0.65 0.65 0.69 0.58 0.68 0.59

14 0.74 0.53 0.64 0.84 0.64 0.57 0.51 0.48 0.58 0.79 0.70 0.72 0.55 0.60 1.00 0.86 0.83 0.83 0.65 0.85 0.74 0.76

15 0.73 0.59 0.65 0.89 0.63 0.59 0.50 0.50 0.72 0.81 0.64 0.74 0.63 0.64 0.86 1.00 0.86 0.70 0.67 0.79 0.68 0.68

16 0.75 0.56 0.69 0.94 0.68 0.65 0.56 0.53 0.73 0.88 0.68 0.79 0.63 0.65 0.83 0.86 1.00 0.80 0.70 0.74 0.75 0.72

17 0.66 0.54 0.65 0.92 0.74 0.59 0.54 0.41 0.62 0.80 0.58 0.71 0.58 0.65 0.83 0.70 0.80 1.00 0.52 0.68 0.65 0.52

18 0.55 0.56 0.59 0.71 0.56 0.62 0.51 0.55 0.67 0.66 0.58 0.60 0.65 0.69 0.65 0.67 0.70 0.52 1.00 0.51 0.53 0.48

19 0.71 0.52 0.61 0.92 0.69 0.55 0.52 0.42 0.63 0.82 0.66 0.74 0.60 0.58 0.85 0.79 0.74 0.68 0.51 1.00 0.67 0.68

20 0.69 0.64 0.60 0.98 0.64 0.64 0.51 0.52 0.64 0.76 0.58 0.71 0.59 0.68 0.74 0.68 0.75 0.65 0.53 0.67 1.00 0.63

21 0.67 0.55 0.71 0.75 0.61 0.51 0.57 0.46 0.60 0.71 0.59 0.66 0.58 0.59 0.76 0.68 0.72 0.52 0.48 0.68 0.63 1.00

Table 6. Values of the Smix Similarity Factor for All Operation Modes ID 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

0 1.00 0.63 0.62 0.71 0.67 0.68 0.51 0.69 0.68 0.82 0.76 0.66 0.66 0.69 0.74 0.73 0.75 0.66 0.55 0.71 0.69 0.67

1 0.63 1.00 0.63 0.57 0.51 0.70 0.56 0.70 0.67 0.61 0.65 0.52 0.69 0.76 0.53 0.59 0.56 0.54 0.56 0.52 0.64 0.55

2 0.62 0.63 1.00 0.58 0.66 0.69 0.65 0.67 0.74 0.64 0.71 0.67 0.64 0.75 0.64 0.65 0.69 0.65 0.59 0.61 0.60 0.71

3 0.71 0.57 0.58 1.00 0.90 0.82 0.72 0.82 0.83 0.91 0.95 0.87 0.82 0.81 0.84 0.89 0.94 0.92 0.71 0.92 0.98 0.75

4 0.67 0.51 0.66 0.90 1.00 0.74 0.48 0.61 0.58 0.71 0.57 0.60 0.61 0.65 0.64 0.63 0.68 0.74 0.56 0.69 0.64 0.61

5 0.68 0.70 0.69 0.82 0.74 1.00 0.50 0.79 0.68 0.64 0.71 0.61 0.75 0.72 0.57 0.59 0.65 0.59 0.62 0.55 0.64 0.51

6 0.51 0.56 0.65 0.72 0.48 0.50 1.00 0.53 0.55 0.56 0.53 0.48 0.59 0.59 0.51 0.50 0.56 0.54 0.51 0.52 0.51 0.57

7 0.69 0.70 0.67 0.82 0.61 0.79 0.53 1.00 0.59 0.49 0.60 0.49 0.63 0.68 0.48 0.50 0.53 0.41 0.55 0.42 0.52 0.46

8 0.68 0.67 0.74 0.83 0.58 0.68 0.55 0.59 1.00 0.61 0.77 0.66 0.79 0.86 0.58 0.72 0.73 0.62 0.67 0.63 0.64 0.60

9 0.82 0.61 0.64 0.91 0.71 0.64 0.56 0.49 0.61 1.00 0.90 0.79 0.77 0.84 0.79 0.81 0.88 0.80 0.66 0.82 0.76 0.71

compared to those of SPCA and SCS, indicating that the operating modes are different from each other. Successful identification is dependent on the discrimination potential of the disturbance models, because the high similarity between models would limit the potential to distinguish between the corresponding disturbances. Therefore, the low similarity values of Smix show the power of the new method for mode identification. The values of the similarity factor SPCA are much larger than those of Smix, which indicates that modes are more difficult to identify using the PCA method. However, it does classify many operating modes successfully. When we consider the similarity values of SCS, many values are ∼1. Because of the fact that this method only uses the information of the first principal component, the misclassification rates are high for most of the operating modes. When the first principal direction is similar between the two operating modes, the similarity values are equal to 1. Besides, the similarity values between fault modes 2, 6, 21, and other modes are close to zero, which indicates that the first principle directions of modes 2, 6, and 21 are very different from that of the other operation modes; therefore, these three modes are successfully identified, compared to others. In summary, Smix shows the most efficient result for operation mode identification;

however, SCS gives the worst identification result, and the result of SPCA resides between these two results. 6. Conclusion A novel strategy has been developed for fault detection and mode identification, based on locating previous occurrences in large historical databases. The new fault detection method is based on independent component analysis-principal component analysis (ICA-PCA), which can extract statistical independent features from process data with a non-Gaussian distribution. Given that ICA is more efficient than PCA for extracting nonGaussian information, the use of ICA can greatly improve the monitoring performance. Besides, a new similarity factor is proposed for mode identification, which consists of four similarity factors: SICA, SPCA, and their corresponding distance similarity factors. Because the independent components (ICs) have no orthogonal constraint, and these extracted ICs have equal roles in the ICA model, the “main angle” is proposed to calculate the similarity between the two ICA models. To handle the cases where two datasets have similar principal components (PCs) or ICs but the numerical values of the process variables

Ind. Eng. Chem. Res., Vol. 46, No. 7, 2007 2063

are very different, two distance similarity factors are used as complements. The proposed fault detection and mode identification methods have been evaluated in the Tennessee Eastman (TE) process. The results presented in the case study show that the proposed methods provide superior power of fault detection and mode identification, compared to alternative PCA-based methods. However, these proposed methods are limited to capturing the linear nature of the process: how to extend them to nonlinear processes is one of our research directions. On the other hand, batch process monitoring is also included in our research interests; in contrast to continuous processes, batch process monitoring is more difficult, and much future work must be done in this area. Acknowledgment This work was financially supported by the National Natural Science Foundation of China (No. 20576116). Literature Cited (1) Lowry, C.; Woodall, W.; Champ, C.; Rigdon, S. A multivariate exponentially weighted moving average control chart. Technometrics 1992, 34, 46-53. (2) Pignatiello, J.; Runger, G. Comparisons of multivariate CUSUM charts. J. Qual. Technol. 1990, 22, 173-186. (3) Kresta, J.; MacGregor, J.; Marlin, T. Multivariate statistical monitoring of process operating performance. Can. J. Chem. Eng. 1991, 69, 3547. (4) Simoglou, A.; Georgieva, P.; Martin, E. B.; Morris, A. J.; Feyo de Azevedo, S. On-line monitoring of a sugar crystallization process. Comput. Chem. Eng. 2005, 29, 1411-1422. (5) Hwang, D. H.; Han, C. H. Real-time monitoring for a process with multiple operating modes. Control Eng. Pract. 1999, 7, 891-902. (6) Kano, M.; Nagao, K.; Hasebe, H.; Hashimoto, I.; Ohno, H.; Strauss, R.; Bakshi, B. R. Comparison of multivariate statistical process monitoring methods with applications to the Eastman challenge problem. Comput. Chem. Eng. 2002, 26, 161-174. (7) Kano, M.; Maruta, H.; Tanaka, S.; Hasebe, S.; Hashimoto, I.; Ohno, H. Multivariate statistical process control with dynamic external analysis. In Proceedings of the SICE 41st Annual Conference; The Society of Instrument and Control Engineers (SICE): Tokyo, 2002; pp 3081-3086. (8) Zhao, S. J.; Zhang, J.; Xu, Y. M. Monitoring of processes with multiple operation modes through multiple principle component analysis models. Ind. Eng. Chem. Res. 2004, 43, 7025-7035. (9) Lee, J. M.; Yoo, C. K.; Lee, I. B. Enhanced process monitoring of fed-batch penicillin cultivation using time-varying and multivariate statistical analysis. J. Biotechnol. 2004, 110, 119-136. (10) Chen, J. H.; Liu, K. C. On-line batch process monitoring using dynamic PCA and dynamic PLS models. Chem. Eng. Sci. 2002, 57, 6375. (11) Nomikos, P.; MacGregor, J. F. Monitoring batch processes using multiway principal component analysis. AIChE J. 1994, 40, 1361-1375. (12) Nomikos, P.; MacGregor, J. F. Multivariate SPC charts for monitoring batch process. Technometrics 1995, 37, 41-59. (13) Nomikos, P.; MacGregor, J. F. Multi-way partial least square in monitoring batch processes. Chem. Intell. Lab. Syst. 1995, 30, 97-108. (14) Smilde, A. K.; Westerhuis, J. A; Boque, R. Multiway multiblock component and covariates regression models. J. Chemom. 2000, 14, 301331. (15) Ku, W.; Storer, R. H.; Georgakis, C. Disturbance detection and isolation by dynamic principal component analysis. Chem. Intell. Lab. Syst. 1995, 30, 179-196. (16) Bakshi, B. R. Multiscale PCA with applications to multivariate statistical process monitoring. AIChE J. 1998, 44, 1596-1610. (17) Isermann, R. Model-based fault-detection and diagnosis-status and applications. Ann. ReV. Control 2005, 29, 71-85. (18) Yoon, S.; MacGregor, J. F. Statistical and causal model-based approaches to fault detection and isolation. AIChE J. 2000, 46, 1813-1824. (19) Yoon, S.; MacGregor, J. F. Fault diagnosis with multivariate statistical models: I. Using steady state fault signature. J. Process Control 2001, 387-400.

(20) Kano, M.; Tanaka, S.; Hasebe, S.; Hashimoto, I.; Ohno, H. Monitoring independent components for fault detection. AIChE J. 2003, 49, 969-976. (21) Cheng, C.; Chiu, M. S. Nonlinear process monitoring using JITLPCA. Chem. Intell. Lab. Syst. 2005, 76, 1-13. (22) Sarbu, C.; Pop, H. F. Principle component analysis versus fuzzy principle component analysis A case study: the quality of danube water (1985-1996). Talanta 2005, 65, 1215-1220. (23) Hyvarinen, A.; Oja, E. Independent component analysis: algorithms and applications. Neural Network 2000, 13, 411-430. (24) Lee, J. M.; Yoo, C. K.; Lee, I. B. Statistical process monitoring with independent component analysis. J. Process Control 2004, 14, 467485. (25) Lee, J. M.; Yoo, C. K.; Lee, I. B. Statistical monitoring of dynamic processes based on dynamic independent component analysis. Chem. Eng. Sci. 2004, 59, 2995-3006. (26) Kano, M.; Hasebe, S.; Hashimoto, I.; Ohno, H. Evolution of multivariate statistical process control: application of independent component analysis and external analysis. Comput. Chem. Eng. 2004, 28, 11571166. (27) Kano, M.; Tanaka, S.; Hasebe, S.; Hashimoto, I.; Ohno, H. Combined multivariate statistical process control. In Proceedings of the IFAC Symposium on AdVanced Control of Chemical Processes (ADCHEM 2004), January 11-14, 2004; Elsevier: Amsterdam, 2004: pp 303-308. (28) Dunia, R.; Qin, S. J. Subspace approach to multidimensional fault identification and reconstruction. AIChE J. 1998, 44, 1813-1831. (29) Zhang, J.; Martin, E. B.; Morris, A. J. Fault detection and diagnosis using multivariate statistical techniques. Chem. Eng. Res. Des. 1996, 74 (1), 89-96. (30) Krzanowski, W. J. Between-groups comparison of principal components. JASA, J. Am. Stat. Assoc. 1979, 74, 703-707. (31) Johannesmeyer, M. C.; Singhai, A.; Seborg, D. E. Pattern matching in historical data. AIChE J. 2002, 48, 2022-2038. (32) Singhai, A.; Seborg, D. E. Clustering of multivariate time-series data. Proc. Am. Control Conf. 2002, 3931-3936. (33) Singhai, A.; Seborg, D. E. Pattern matching in multivariate time series databases using a window-moving approach. Ind. Eng. Chem. Res. 2002, 41, 3822-3838. (34) Singhai, A.; Seborg, D. E. Data compression issues with pattern matching in historical data. Proc. Am. Control Conf. 2003, 3696-3701. (35) Hyvarinen, A. A fast fixed-point algorithm for independent component analysis. Neural Comput. 1997, 9, 1483-1492. (36) Chen, Q.; Wynne, R. J.; Goulding, P.; Sandoz, D. The application of principal component analysis and kernel density estimation to enhance process monitoring. Control Eng. Pract. 2000, 8, 531-543. (37) Chen, Q.; Kruger, U.; Leung, A. T. Y. Regularised kernel density estimation for clustered process data. Control Eng. Pract. 2004, 12, 267274. (38) Downs, J. J.; Vogel, E. F. A plant-wide industrial process control problem. Comput. Chem. Eng. 1993, 17 (3), 245-255. (39) Raich, A. C.; C¸ inar, A. Multivariate statistical methods for monitoring continuous processes: assessment of discrimination power of disturbance models and diagnosis of multiple disturbances. Chem. Intell. Lab. Syst. 1995, 30, 37-48. (40) Singhai, A.; Seborg, D. E. Evaluation of a pattern matching method for the Tennessee Eastman challenge process. J. Process Control 2006, 601-613. (41) Chiang, L. H.; Russell, E. L.; Braatz, R. D. Fault Detection and Diagnosis in Industrial Systems; Springer-Verlag: London, 2001. (42) Lyman, P. R.; Georgakist, C. Plant-wide Control of the Tennessee Eastman Problem. Comput. Chem. Eng. 1995, 19, 321-331. (43) Raich, A.; Cinar, A. Diagnosis of process disturbance by statistical distance and angle measures. Comput. Chem. Eng. 1997, 21, 661.

ReceiVed for reView August 15, 2006 ReVised manuscript receiVed January 23, 2007 Accepted January 24, 2007 IE061083G