A Novel Statistical-Based Monitoring Approach for Complex

Apr 14, 2009 - Vine Copula-Based Dependence Description for Multivariate Multimode Process Monitoring. Xiang Ren , Ying Tian , and Shaojun Li. Industr...
0 downloads 0 Views 170KB Size
4892

Ind. Eng. Chem. Res. 2009, 48, 4892–4898

A Novel Statistical-Based Monitoring Approach for Complex Multivariate Processes Zhiqiang Ge,† Lei Xie,‡ and Zhihuan Song*,† State Key Laboratory of Industrial Control Technology, Institute of Industrial Process Control, Zhejiang UniVersity, Hangzhou 310027, Zhejiang, China, and State Key Laboratory of Industrial Control Technology, Institute of Cyber Systems and Control, Zhejiang UniVersity, Hangzhou 310027, Zhejiang, China

Conventional methods are under the assumption that a process is driven by either non-Gaussian or Gaussian essential variables. However, many complex processes may be simultaneously driven by these two types of essential source. This paper proposes a novel independent component analysis and factor analysis (ICA-FA) method to capture the non-Gaussian and Gaussian essential variables. The non-Gaussian part is first extracted by ICA and support vector data description is utilized to obtain tight confidence limit. A probabilistic approach is subsequently incorporated to separate the residual Gaussian part into latent influential factors and unmodeled uncertainty. By retrieving the underlying process data generating structure, ICA-FA facilitates the diagnosis of process faults that occur in different sources. A further contribution of this paper is the definition of a new similarity factor based on the ICA-FA for fault identification. The efficiency of the proposed method is shown by a case study on the TE benchmark process. 1. Introduction Multivariate statistical process control (MSPC) has been intensively researched and applied to chemical plants.1-6 However, the traditional MSPC schemes are formulated on the assumption that both the process variables and latent impelling variables are Gaussian distributed, and thus the process can be described by the mean and covariance of the variables. However, process variables may not follow Gaussian distribution exactly, due to recycles, system delays, and unmeasured disturbances of the process, etc. This implies that the assumption which the conventional MSPC is based upon is violated, which, in turn, causes higher type I and II monitoring errors. As an emerging technology for dealing with non-Gaussian signals,7 independent component analysis (ICA) was introduced into process monitoring area recently to extract the non-Gaussian sources out of the process recordings.8-12 It is also important to note that all reported ICA-based monitoring methods are based on the assumption that the essential variables which drive the process are non-Gaussian and discarded Gaussian parts are utilized to set up a single SPE monitoring chart. Actually, complex multivariate processes are always driven by nonGaussian and Gaussian essential variables simultaneously. Separating the Gaussian information from the unmodeled Gaussian uncertainty will facilitate the diagnosis of faults which occur in different sources. Due to complicated confidence limit determination of the traditional ICA-based monitoring method, which is based on the kernel density estimation (KDE), this paper introduces support vector data description (SVDD) to improve the determination of confidence limits for the non-Gaussian components.13,14 Compared to KDE,10 the SVDD-based method is more computationally efficient, which relies on a quadratic programming cost function.13 To monitor discarded Gaussian parts of the process information, PCA might be the straightforward but not the optimal approach in which the projected * To whom correspondence should be addressed. Tel.: +86-57187951442-808. E-mail: [email protected]. † Institute of Industrial Process Control, Zhejiang University. ‡ Institute of Cyber Systems and Control, Zhejiang University.

directions have the maximum variations irrespective of the underlying information generating structure. More precisely, the process recordings always consist of inevitable measurement noises, while PCA is not capable of characterizing the different influences from common process Gaussian impelling factors and random variations due to measurement equipments.15 A probabilistic generative model (probabilistic PCA, PPCA) was employed to address the above issue via the expectation and maximization (EM) algorithm. However, PPCA assumes the same noise level of all measurement variables, which cannot be satisfied in practice. In the present paper, a general extension of PPCA model, factor analysis (FA) can model different noise levels of measurement variables. To determine the number of essential variables, lots of useful approaches have been reported, including Akaike information criterion (AIC), Bayesian information criterion (BIC), SCREE test, PRESS cross validation test,16 and the more recent fault detection criteria.17-19 Similarly, this paper addresses this issue by selecting the optimal factor number (PFs) which achieves the best fault detection performance. The final contribution is the development of a novel similarity factor based on non-Gaussian and Gaussian systematic subspaces and the noise subspace. The proposed similarity index is used to classify various different process faults. The PCA similarity factor was first proposed by Krzanowski,20 and later improved by several researchers.21-25 Recently, an ICA-based similarity factor was also proposed as an extension of PCA counterpart in the non-Gaussian subspace.26 However, none of the similarity factors took the noise subspace into consideration, which is engendered by the projection method. In the ICA-FA framework proposed in this paper, the noise subspace characterization can be well described, which, in turn, leads to the definition of a new similarity factor. The rest of this paper is structured as follows. First, some preliminaries are given, including ICA, SVDD, and FA algorithms. The new monitoring scheme is described in the next section. Then the proposed monitoring and fault identification methods are demonstrated in section 4, followed by the TE benchmark process simulation study. Finally, some conclusions are made in section 5.

10.1021/ie800935e CCC: $40.75  2009 American Chemical Society Published on Web 04/14/2009

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4893

2. Preliminaries Brief descriptions of the ICA, SVDD, and FA algorithms are given in this section. 2.1. Independent Component Analysis (ICA). ICA was originally proposed to solve the blind source separation problem.7 The process measured variables, z ∈ Rl can be expressed as linear combinations of r(el) unknown independent components s, and the relationship between them is given by27 z ) As (1) where z is assumed to be whitened and A ∈ Rl×r is the mixing matrix. The basic problem of ICA is to calculate a separating matrix W (which is the inverse of A) so that the components of the reconstructed independent components (ICs) sˆ become as independent of each other as possible. Besides, they are nonGaussian distributed. The relation is given as sˆ ) Wz ≈ s (2) There are several methods for non-Gaussian measurement, among which negentropy is used in the present paper.7 2.2. Support Vector Data Description (SVDD). The main idea of SVDD is to map the input vectors to a feature space, and to find hypersphere with the minimum volume which separates the transferred data from the rest of the feature space. Applications have shown a high generalization performance of SVDD if a large reference data set with few abnormal samples is available.14 Suppose a data set containing n samples sˆi ∈ Rr, i ) 1,2,...,n, is given, the mapping from the original space to the feature space Φ:sˆfF can be simply done by a given kernel K(sˆi,sˆj) ) 〈Φ(si),Φ(sˆj)〉, which computes the inner product in the feature space. In the present paper, the most popular Gaussian kernel is used. To construct the minimum volume of the hypersphere, SVDD solves the following optimization problem: n

min R2 + C R,a,ξ

∑ξ

i

s . t . |Φ(sˆi) - a| e R + ξi,

n

Ri

n

n

∑ R K(sˆ , sˆ )-∑ ∑ R R K(sˆ , sˆ ) i

i

j

i)1

i j

i

j

i)1 j)1 n

s . t . 0 e Ri e C,



(4) Ri ) 1

i)1

where Ri is Lagrange multipliers. After the hpyersphere in the feature space has been constructed, the hypothesis that a new sample sˆnew is normal is accepted if the distance d(Φ(sˆnew)) ) |Φ(sˆnew) - a| e R, which is given by n

d(Φ(sˆnew)) ) [K(sˆnew, sˆnew) - 2

∑ R K(sˆ , sˆ i

i

new) +

i)1

n

x ) Pt + e

(6)

where P ) [p1,p2,...,pl] ∈ Rm×l is the loading matrix, just as that in PCA. The variances matrix of measurement noise e are represented by Σ ) diag{λi}1,2,...,m, in which different noise levels of measurement variables are assumed. If all of λ are assumed to have the same value, then FA is equivalent to PPCA. If all of λ are assumed to be zero, then FA becomes PCA. Therefore, FA is the general description of the Gaussian latent model structure. PPCA and PCA are just two special cases of FA. FA and PPCA are achieved in probability space, which is the key difference between probabilistic methods and the projection method (PCA). The parameter set Θ ) {P,Σ} can be estimated by the expectation and maximization (EM) algorithm, which is an iterative likelihood maximization algorithm.15 EM is a widely used algorithm for parameter estimation due to its simplicity and efficiency. Besides, it can also handle incomplete data set, which is also very attractive to the application in process industry. 3. Independent Component Analysis-Factor Analysis (ICA-FA)

ξi g 0

where a is the center of the hypersphere in the feature space. The variable C gives the trade-off between the volume of the sphere and the number of errors (number of target objects rejected). ξi represents the slack variable which allows the probability that some of the training samples can be wrongly classified. The dual form of the optimization problem can be obtained as follows: min

tributions are Gaussian, and the original measurement variables x are treated as a linear combination of t plus small additive white noise e. The aim of FA is to find the most probable parameter set Θ ) {P,Σ} in the model structure, which is described as follows:

(3)

i)1 2

2

Figure 1. Description of ICA-FA versus ICA and FA.

n

∑ ∑ R R K(sˆ , sˆ )]

1/2

i j

i

j

(5)

i)1 j)1

2.3. Factor Analysis (FA). Similar to PPCA described in ref 15, FA concentrates on the latent variables t whose dis-

As demonstrated in the preliminaries, conventional statistical methods are under the assumption that processes are driven by fewer essential variables. These essential variables are either non-Gaussian (ICA) or Gaussian distributed (FA). However, under more general circumstances, processes are driven by nonGaussian and Gaussian essential variables simultaneously, which formulates the key assumption of the proposed mixed essential component analysis (ICA-FA) method. Figure 1 shows the difference between ICA-FA and the conventional methods (ICA and FA). The relationship of measurement variables and mixed essential variables is given as follows: x ) x1 + x2 ) As + Pt + e

(7)

where x1 and x2 represent non-Gaussian and Gaussian part of the process, respectively, A is the mixing matrix of non-Gaussian essential variables s, P is the loading matrix of Gaussian essential variables t, and e is Gaussian unmodeled uncertainty, relating to measurement noises, etc., with zero means and different variance Σ ) diag{λi}1,2,...,m, which permits measurement variables to have different noise levels. As a general form of PPCA, FA performs better than PCA and PPCA. Since process measurements are always corrupted by noise, the ignorance of such information generation structure will make the PCA method leave some essential information in the discarded space and trigger monitoring errors. In contrast, the

4894

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

noise structure is considered in PPCA and FA, and thus the useful information can be correctly extracted and monitored separately. To estimate the new parameter set Θ ) {A,P,Σ} and extract essential variables EA ) {s,t} from measurement variables x, a two-step estimation and extraction strategy can be developed. First, non-Gaussian essential variables can be extracted by the efficient FastICA or our previously proposed algorithm,12 and thus the mixing matrix A is obtained. Then the parameter P and Σ can be estimated by the EM algorithm, which is described as follows. According to the PPCA method,15 the mean and variance of the Gaussian part of measurement variables x2 can be calculated as

produces a unique analytical solution. Given the training data set S (the extracted independent components), the hypersphere is constructed. The center a and the radius R can be determined by14

µx2 ) PΕ[t] + Ε[e] ) 0 Σx2 ) PΕ[ttT]PT + Ε[eeT] ) PPT + Σ (8)

NGSi ) d2(Φ(sˆi)) ) |Φ(sˆi) - a|2 e NGSlim ) R2

Hence, the distribution of x2 is p(x2) ) Ν{0, PPT + Σ}. EM is an iterative algorithm; it consists of two steps:15 E-step and M-step. In E-step, two statistics of the Gaussian essential variable t is estimated ti ) Ε[ti|x2i] ) Qx2i T Ε[titi |x2i] ) I - QP + Qx2ix2TiQT ˆ

(9)

where Q ) PT(PPT + Σ)-1. In the M-step, the parameter P and Σ are estimated so that the value of the likelihood function is maximized, the estimated parameters are calculated as n

P)(

n

∑ x tˆ )(∑ Ε[t t |x ])) i)1 n

Σ)

∑ i)1

T 2i i

T i i

i)1

-1

2i

(10)

diag(x2ix2Ti - Ptˆix2Ti) ⁄ n

where n is the number of samples, and diag(Z) means to make a diagonal matrix using diagonal elements of matrix Z. Equations 9 and 10 are calculated iteratively until both of the parameters (P and Σ) are converged. In summary, when the reference data set X ) {xi}i)1,2,...,n is available, with the ICA-FA model defined by eq 7, process information can be separated into three different parts. The systematic information includes non-Gaussian and Gaussian parts which are driven by their corresponding essential variables, while the noise information and other unmeasured disturbances are represented by e. 4. ICA-FA Based Process Monitoring and Fault Identification Scheme In this section, the new proposed monitoring scheme is demonstrated. The whole monitoring scheme is discussed in subsection 4.1, followed by the essential variable number selection method in the subsection 4.2. Then the fault identification scheme will be proposed. 4.1. Monitoring Scheme. As described in section 3, the entire information of the process is partitioned into three parts: non-Gaussian systematic information, Gaussian systematic information, and noise information. Therefore, three different monitoring statistics are developed for fault detection. First, for monitoring non-Gaussian information, SVDD-based method is used, which is briefly described in the preliminaries. In contrast to the confidence limit determination method KDE, the SVDD reduces the complexity to quadratic programming problem that

n

a)

∑ R Φ(sˆ ) i

i

i)1

n

R ) [1 - 2

∑ i)1

n

RiK(sˆz, sˆj) +

(11)

n

∑∑

RiRjK(sˆi, sˆj)]1/2

i)1 j)1

The non-Gaussian statistic (NGS) can be developed as the distance between the data sample and the center of the hypersphere, thus (12)

where sˆi is the extracted independent component of the new sample, d(Φ(sˆi)) is given in eq 5. If NGSi e R2, the nonGaussian information is regarded as normal, else it will be regarded as an outlier or some fault and disturbance has happened. Since the distribution of Gaussian essential variables t is assumed to be Ν{0,I}, the Mahalanobis norm of t is identical to the Euclidian norm. Different from the projection method, the Gaussian essential variables cannot be calculated directly. However, it can be substituted by their estimation, which is given in eq 9. It is known that the squared Mahalanobis norm of t follows χ2 distribution, the confidence limit of T2 statistic can be determined as follows: T2i ) |tˆi |2 ) x2TiQTQx2i e χ2R,l

(13)

2 where χR,l represent the confidence limit with R confidence and l is the corresponding freedom parameter. If the value of T2 statistic goes beyond the confidence limit, the Gaussian systematic process is regarded as out of control. Similarly, the noise part of the process information can also be monitored by the squared Mahalanobis distance. Since the Gaussian information matrix x2 can be represented as

x2 ) Pt + e

(14)

Note eqs 9 and 14: the noise information e can be estimated as follows: eˆi ) Ε[ei|x2i] ) (I - PQ)x2i

(15)

Then the confidence limit of squared prediction error (SPE) statistic can be determined by 2 SPEi ) |Σ-1⁄2ei |2 ) x2Ti(I - PQ)TΣ-1(I - PQ)x2i e χR,m

(16) where χ2R,m represent the confidence limit with R confidence and m is the corresponding freedom parameter. If the value of SPE statistic goes beyond the confidence limit, the noise part of the process is regarded abnormal. This may due to the change of noise level, some unmeasured disturbances, and other process structure changes. Note that these two Gaussian parts of information both follow χ2 distribution; they can be monitored together. According to eq 8, the entire Gaussian statistic (GS) can be constructed as follows: 2 (17) GSi ) |(PPT + Σ)-1⁄2x2i |2 ) x2Ti(PPT + Σ)-1x2i e χR,m

Therefore, if changes happen in Gaussian part, the monitored entire Gaussian statistic (GS) should go beyond its confidence

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

limit. However, if one wants to isolate the change between the systematic Gaussian part and the noise part, T2 and SPE statistics should be chosen. In summary, four statistics (NGS, T2, SPE, and GS) are developed for process monitoring. By monitoring NGS and GS, the change between non-Gaussian and Gaussian information can be isolated. Furthermore, the change of Gaussian information can be isolated by monitoring T2 and SPE separately. Hence, three types of changes can be monitored separately by their corresponding statistics. Change of any type or types of information will be reflected in related monitoring charts. 4.2. Number of Essential Variables Selection. The number selection of essential variables consists of two parts: ICs number selection and PFs number selection. For simplicity, the ICs number can be selected by the initial PCA step in the ICA calculation. Conventionally, the selection of the number of PCs under the probabilistic model, as discussed in ref 15, can be achieved by the average variance explanation ratio of the PCs or the explanation ratio of the systematic part over the entire process information. In this paper, the similar idea can be used to determine the number of retained PFs. However, when some kind of process knowledge is available, the selection of PFs can be improved. Actually, the conventional selection of PFs/ PCs is very subjective and has not considered the monitoring performance, e.g., fault detection. Based on the history data records and some available process knowledge, a fault set can be defined to include all appreciate important faults in the process.19 Suppose the fault set includes K faults, which is described as {Fi, i ) 1,2,...,K}. Based on the process knowledge, a fault subspace can be defined by the propagation manner of the fault.28 Denote {Θi, i ) 1,2,...,K} as the fault subspace with dimensions {dim(Θi) g 1, i ) 1,2,...,K} for the defined fault set. The key idea of this number selection method is to explore the fault detection performance under different PF numbers. The number of PFs which yields the best performance will be selected. The critical fault magnitude (CFM) analysis technique which is proposed in Wang et al.18 can be used for performance evaluation. First, consider the fault in the principal factor space, when a fault Fi happened, the Gaussian process data can be described as28 x2 ) x*2 + Θifi

(18)

where the quantity x*2 is the normal process data, and the fault magnitude is |fi|. According to eqs 13 and 18, the principal factor and its corresponding T2 statistic can be obtained as: T2 ) |tˆ |2 ) |Qx2 |2 ) |Q(x*2 + Θifi)|2

(20)

2 Since x*2 represents the normal process data, |Qx*2 | e χR,l . 2 Thus, the T statistic will detect Fi if the following condition is satisfied:28

|QΘifi | g 2√χ2R,l

(21)

Similarly, the sufficient condition for the fault Fi to be detected by the SPE statistic in the noise space can be obtained in the same way: 2 |Σ-1⁄2(I - PQ)Θifi | g 2√χR,m

(22)

Noting eq 21, the T2 statistic based sufficient fault magnitude with l principal factors should satisfy

(23)

where σmax(QΘi) is the maximum singular value of QΘi. Similarly, according to eq 22, the SPE statistic based sufficient fault magnitude with l principal factors should satisfy l 2 -1⁄2 | ≈ 2σ-1 (I - PQ)Θi)√χR,m |fi | g |f i,SPE max(Σ

(24)

-1/2

where σmax(Σ (I - PQ)Θi) is the maximum singular value of Σ-1/2(I - PQ)Θi. Integrating the two statistics T2 and SPE into the GS statistic, the sufficient fault magnitude with l principal factors satisfies l 2 T -1⁄2 | ≈ 2σ-1 Θi)√χR,m |fi | g |f i,GS max((PP + Σ)

(25)

-1/2

where σmax((PP + Σ) Θi) is the maximum singular value of (PPT + Σ)-1/2Θi. Consider the fault detection performance: the smaller the magnitude of the fault to be detected, the better monitoring performance is to be yielded. Noting the eqs 23-25, the values of fault magnitude are function of the number of principal factors l. The optimal number of principal factors lopt can be obtained to minimize of the follow optimization function: T

K

lopt ) arg min{FDP} ) arg min{ l

l

K

)arg min{ l

∑ i)1

∑ (w

l i,GS|f i,GS |)}

i)1

(26)

2 T -1⁄2 (2wi,GSσ-1 Θi)√χR,m )} max((PP + Σ)

where FDP represents “fault detection performance”, which resembles OCFM in ref 19. wi,GS is the weight of every fault corresponding to the GS statistics. To emphasize the importance of some faults, different weights could be given to different faults. Minimum the entire fault detection performance function, the optimal number of principal factor l ) lopt will be obtained. However, in some cases, only some of the process fault knowledge is available. In these cases, the available knowledge or some important faults can be used to determine the optimal number of principal factor lopt. 4.3. ICA-FA Based Fault Identification Scheme. In this subsection, the similarity factor is developed to measure the similarity between the faulty data sets and it can be used to identify the origin of the fault happening by comparing it to a historical faulty database. The PCA similarity factor SPCA was first proposed by Krzanowski20 and is calculated as follows: l

SPCA )

(19)

Using the triangle inequality of the norm operator, we obtain |Q(x*2 + Θifi)| g |QΘifi | - |Qx*2 |

-1 2 l |fi | g |f i,T 2 | ≈ 2σmax(QΘi)√χR,m

4895



1 cos2 θi l i)1

(27)

where θi is the angle between the i-th principal component two PCA subspaces. By introducing different weight for each principal component, this similarity factor has been modified by Singhal and Seborg.22 Besides, Zhao et al.25 also modified the similarity factor by introducing the principal angle, in order to mine the most similar information between two data sets. For the similarity measure between two non-Gaussian models (ICA), ICA-based similarity factor SICAwas developed in our previous work26 which is calculated as follows: r

SICA )



1 cos2 θi r i)1

(28)

where θi is the i-th main angle between two ICA subspaces. Compared with the projection method (PCA), FA estimates the loading matrix P in the probability space. The loading matrix Ppca is ranged due to the value of its corresponding eigenvalue,

4896

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

and the matrix is restricted that PTpcaPpca ) I. However, the new estimated loading matrix P has no such properties. Besides, the noise part is also estimated in FA, which has not been considered in PCA. Suppose two FA loading matrices are represented as U and V, according to the concept of main angle,26 the similarity factor of FA can be developed as follows: l

SFA )



1 cos2 θi l i)1

(29)

where θi is the i-th main angle between U and V, which is defined as cos θi ) max max

(

u∈Ui-1 v∈Vi-1

)

uTi vi uTv ) |u|2|v|2 |ui |2|vi |2

(30)

where ui and vi are vectors corresponding to the i-th main angle and Ui-1 and Vi-1 are residual subspaces when the first i - 1 pairs of main vectors are removed. Hence, the most similar information between the two data sets can be extracted by FAbased similarity factor SFA. Next, we consider the similarity between the noise parts. Suppose the estimated noise matrices of the two fault data sets are Σ1 and Σ2. Due to the diagonal property of the noise covariance matrix, the diagonal elements can be extracted thus: λ1 ) diag(Σ1) λ2 ) diag(Σ2)

(31)

where λ1 and λ2 are two vectors that consist of diagonal elements of Σ1 and Σ2. The similarity factor between the noise parts of two fault patterns can be defined as SN ) √e-d (λ1,λ2) d(λ1, λ2) ) |λ1, λ2 |2 2

(32)

In summary, there are three similarity factors SICA, SFA, and SN defined for similarity measurement of two fault data sets. When these three similarity factors are used simultaneously, a weighted average of these similarity factors can be defined as follows: Savg ) R · SICA + β · SFA + γ · SN

(33)

where 0 e R,β,γ e 1 and R + β + γ ) 1 and are weighted factors. How to choose theses weighted factors is an open question. However, in our experience, R and β can be chosen bigger than γ, because the systematic part takes more responsibility that two data sets differentiate from each other. The choices of the three weighted factors can be determined as R r ) β l (34) 1 γ ) min{R, β} 2 where r and l are numbers of non-Gaussian and Gaussian essential variables. Therefore, if r g l, eq 33 becomes 2r 2l l ·S + ·S + ·S 2r + 3l ICA 2r + 3l FA 2r + 3l N else if r < l, eq 33 becomes Savg )

Savg )

2r 2l r ·S + ·S + ·S 3r + 2l ICA 3r + 2l FA 3r + 2l N

(35)

(36)

5. Case Study of TE Benchmark Process In this section, the performance of the proposed method is tested through the TE process.29 As a benchmark simulation,

Figure 2. Scatter plot of TE normal operation condition.

the Tennessee Eastman process has been widely used to test the performance of various monitoring approaches.16,30 The second control structure listed in Lyman and Georgaki is used.31 The details of the process description are well explained in the book by Chiang et al.16 The variables selected for monitoring are the same as that for control structure development in Lyman and Georgakist.31 Simulation data which is acquired in the normal process consists of 960 samples, and their sampling interval was 3 min. In order to evaluate fault detection and identification performance of the new proposed method, simulation data under 21 fault processes are also obtained. During the simulation, each fault is introduced in the 161 sample. To select the optimal number of PFs for fault detection, some import process faults are considered, which are listed in Song et al.32 For monitoring of the TE process, the number of ICs and PFs should be chosen appropriately. As a result, 3 ICs and 6 PFs are selected. Then the ICA-FA model can be developed. The kernel parameter for the Gaussian function is chosen as σ ) 2.0, the variable C is selected as C ) 0.0625. For comparison, the ICA, PCA, and ICA-PCA algorithms are also carried out. The number of ICs for ICA is selected as 3, while the principal component (PC) number of PCA is chosen as 6. All confidence limits are set as 99%. The scatter plot of the first independent component versus the second independent component and two confidence limits of NGS and I2 is given in Figure 2. Under normal operation, both of the confidence limits indicate that 99% samples are inside their corresponding regions. However, as discussed previously, the confidence limit of NGS is much tighter than that of I2, which is shown correctly in the figure. Therefore, when some fault happens, NGS is more sensitive to detect this fault. Since samples lie in the gap between the NGS confidence limit and the I2 confidence limit, only NGS can detect them. The I2 statistic regards it as a normal operation because the confidence limit has not been rejected. For each statistic, the monitoring results (type II error) are tabulated in Table 1. The minimum value achieved for each mode is marked with a bold number. As shown in the table, ICA-FA outperforms ICA, PCA, and ICA-PCA for most of fault modes. Particularly, the missing detection rates of ICA-FA for faults 10, 15, 16, and 19-21 are much lower than that of other three methods. As indicated in the 2nd-4th columns of the table, the best monitoring performance of some faults is gained by the non-Gaussian statistic NGS, some are gained by T2, while others are gained by the SPE statistic, which shows the fault isolation ability of the proposed ICA-FA method. In order to test the fault identification performance of the proposed method, 22 data sets (1 normal and 21 faulty) were obtained from simulation. Signature of each data set is cal-

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

4897

Table 1. Monitoring Results of TE Process (Type II Error) ICA-FA

ICA

PCA

ICA-PCA

modes/statistics

NGS

T2

SPE

I2

SPE

T2

SPE

I2

T2

SPE

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

0.990 0.015 0.018 0.911 0.963 0.764 0.008 0.596 0.050 0.933 0.385 0.886 0.124 0.060 0.061 0.900 0.664 0.261 0.091 0.869 0.449 0.374

0.994 0.003 0.018 0.993 0.986 0.778 0.000 0.664 0.029 0.985 0.336 0.789 0.009 0.074 0.001 0.978 0.853 0.103 0.101 0.981 0.490 0.806

0.994 0.005 0.021 0.995 0.994 0.855 0.003 0.811 0.118 0.995 0.271 0.556 0.155 0.083 0.000 0.993 0.841 0.039 0.106 0.809 0.454 0.755

0.988 0.003 0.023 0.974 0.990 0.790 0.014 0.641 0.089 0.981 0.579 0.961 0.081 0.068 0.100 0.976 0.748 0.243 0.105 0.971 0.594 0.434

0.991 0.005 0.018 0.983 0.989 0.773 0.000 0.670 0.031 0.976 0.596 0.605 0.024 0.074 0.001 0.975 0.886 0.059 0.104 0.814 0.594 0.853

0.986 0.008 0.034 0.984 0.990 0.773 0.001 0.618 0.059 0.980 0.685 0.859 0.029 0.064 0.058 0.968 0.788 0.184 0.108 0.905 0.780 0.689

0.988 0.005 0.015 0.974 0.981 0.814 0.000 0.735 0.039 0.974 0.728 0.579 0.088 0.056 0.000 0.975 0.850 0.046 0.101 0.851 0.684 0.554

0.988 0.003 0.023 0.974 0.990 0.790 0.014 0.641 0.089 0.981 0.579 0.961 0.081 0.068 0.100 0.976 0.748 0.243 0.105 0.971 0.594 0.434

0.991 0.011 0.018 0.990 0.993 0.785 0.007 0.681 0.064 0.988 0.841 0.815 0.059 0.098 0.001 0.984 0.913 0.179 0.109 0.829 0.818 0.921

0.989 0.003 0.016 0.979 0.984 0.868 0.000 0.751 0.069 0.975 0.519 0.538 0.100 0.079 0.000 0.980 0.929 0.041 0.098 0.866 0.524 0.874

Table 2. Average Similarity Factor Value of Each Operation Mode fault mode

Savg

Smix

SPCA

fault mode

Savg

Smix

SPCA

0 1 2 3 4 5 6 7 8 9 10 11

0.544 0.583 0.393 0.734 0.633 0.616 0.377 0.539 0.613 0.567 0.715 0.709

0.695 0.622 0.672 0.830 0.658 0.670 0.567 0.598 0.691 0.751 0.694 0.693

0.850 0.824 0.815 0.874 0.835 0.864 0.671 0.854 0.826 0.873 0.879 0.851

12 13 14 15 16 17 18 19 20 21 average

0.443 0.591 0.520 0.767 0.700 0.421 0.364 0.781 0.611 0.737 0.589

0.680 0.715 0.700 0.705 0.733 0.668 0.616 0.677 0.674 0.638 0.679

0.841 0.821 0.831 0.872 0.876 0.844 0.692 0.833 0.850 0.810 0.831

culated; thus the parameter set {Anew,Pnew,Σnew}new)0,1,...21 is obtained. Then the similarities between the new signatures (operation modes) and the preserved signatures are calculated through eq 36. According to our previous work, the similarity values of Smix (the mixed similarity factor) between different operating modes are smaller compared to that of SPCA and SCS (the similarity factor use the first principal vector), indicating that operating modes are different from each other. Successful identification depends on the discrimination potential of the disturbance modes, since high similarity between modes would limit the potential to distinguish between the corresponding faults. So the low similarity values of Smix show power for mode identification. To compare the values of different similarity factors, the average of each fault identification value is tabulated in Table 2, including similarity factor Savg, Smix, and SPCA. The minimum value achieved for each mode is marked with a bold number. The smaller the similarity factor value achieved, the more powerful discrimination potential the corresponding similarity factor shows. Therefore, the minimum value of the similarity factor means the best fault identification performance of that similarity factor. As indicated in Table 2, most of the minimum values are achieved by the Savg similarity factor, while only 5 operation modes show the superiority of the Smix similarity factor. In comparison, SPCA gives the worst fault identification performance. Besides, the average values of all of the 22 operation modes are also tabulated, in which the best average performance is also achieved by the Savg similarity factor.

6. Conclusions The present paper has studied the monitoring of complex multivariate processes, which are not only driven by Gaussian essential variables but also driven by non-Gaussian essential variables. The novel ideas and improvements of the proposed ICA-FA method are summarized as follows. Compared to the conventional methods, which assumes that the process is driven by non-Gaussian essential variables or Gaussian essential variables, the ICA-FA method can be regarded as a generalization of these methods. By monitoring the non-Gaussian and Gaussian information separately, these two types of changes can be isolated simultaneously with fault detection. One step further, by considering the noise part of the Gaussian information, the probabilistic method FA is introduced to separate the two types of changes in Gaussian information (systematic and noise parts). In contrast, two types of changes of Gaussian information are confused in conventional ICA, while PCA confuses the non-Gaussian and Gaussian changes. These confusions result in high miss detection rate and cause conventional methods to be insensitive. By introducing SVDD to construct the new non-Gaussian monitoring statistic and determine its corresponding confidence limit, better description of non-Gaussian components and tighter confidence limit are obtained, which make the new method more sensitive and produce less number of type II errors. To determine the number of PFs, a performance function was developed, which is related to fault detection performance of some important process faults. Besides, a new similarity factor based on the proposed ICA-FA method was developed for fault identification. Selection criterion of the weight factor for each individual similarity factor to form an average similarity factor was also given. Higher identification rate and discriminating potential were shown in the TE simulation example. Future work based on this monitoring scheme will consist of dynamic and nonlinear process behaviors. Acknowledgment This work was supported by the National Natural Science Foundation of China. (Grant no. 60774067, 60736021).

4898

Ind. Eng. Chem. Res., Vol. 48, No. 10, 2009

Literature Cited (1) Kresta, J.; MacGregor, J. F.; Marlin, T. Multivariate statistical monitoring of process operating performance. Can. J. Chem. Eng. 1991, 69, 35–47. (2) Nomikos, P.; MacGregor, J. F. Monitoring batch processes using multiway principal component analysis. AIChE J. 1994, 44, 1361–1375. (3) Nomikos, P.; MacGregor, J. F. Multivariate SPC charts for monitoring batch process. Technometrics 1995, 37, 41–59. (4) Nomikos, P.; MacGregor, J. F. Multi-way partial least square in monitoring batch processes. Chem. Intell. Lab. Syst. 1995, 30, 97–108. (5) Ku, W.; Storer, R. H.; Georgakis, C. Disturbance detection and isolation by dynamic principal component analysis. Chem. Intell. Lab. Syst. 1995, 30, 179–196. (6) Bakshi, B. R. Multiscale PCA with applications to multivariate statistical process monitoring. AIChE J. 1998, 40, 1596–1610. (7) Hyvarinen, A.; Oja, E. Independent component analysis: algorithms and applications. Neural Network 2000, 13, 411–430. (8) Kano, M.; Tanaka, S.; Hasebe, S.; Hashimoto, I.; Ohno, H. Monitoring independent components for fault detection. AIChE J. 2003, 49, 969–976. (9) 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, 1157–1166. (10) Lee, J. M.; Yoo, C. K.; Lee, I. B. Statistical process monitoring with independent component analysis. J. Process Control 2004, 14, 467–485. (11) Lee, J. M.; Qin, S. J.; Lee, I. B. Fault detection and diagnosis based on modified independent component analysis. AIChE J. 2006, 52, 3501– 3514. (12) Xie, L.; Wu, J. Global optimal ICA and its application in MEG data analysis. Neurocomputing 2006, 69, 2438–2442. (13) Tax, D. M. J.; Duin, R. P. W. Support vector domain description. Pattern Recognit. Lett. 1999, 20, 1191–1199. (14) Tax, D. M. J.; Duin, R. P. W. Support vector domain description. Machine Learning 2004, 54, 45–66. (15) Kim, D.; Lee, I. B. Process monitoring based on probabilistic PCA. Chem. Intell. Lab. Syst. 2003, 67, 109–123. (16) Chiang, L. H.; Russell, E. L.; Braatz, R. D. Fault detection and diagnosis in industrial systems; Springer-Verlag: London, 2001. (17) Valle, S.; Li, W.; Qin, S. J. Selection of the number of principal components: The variance of the reconstruction error criterion with a comparison to other methods. Ind. Eng. Chem. Res. 1999, 38, 4389–4410. (18) Wang, H. Q.; Song, Z. H.; Li, P. Fault detection behavior and performance analysis of PCA-based process monitoring methods. Ind. Eng. Chem. Res. 2002, 41, 2455–2464.

(19) Wang, H. Q.; Zhou, H. L.; Hang, B. L. Number selection of principal components with optimized process monitoring performance. IEEE Conf. Decision Control 2004, 4726–4731. (20) Krzanowski, W. J. Between-groups comparison of principal components. J. Am. Star. Assoc. 1979, 74, 703–707. (21) Johannesmeyer, M. C.; Singhai, A.; Seborg, D. E. Pattern matching in historical data. AIChE J. 2002, 48, 2022–2038. (22) Singhai, A.; Seborg, D. E. Clustering of multivariate time-series data. Proc. Am. Control Conf. 2002, 3931–3936. (23) 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. (24) Singhai, A.; Seborg, D. E. Evaluation of a pattern matching method for the Tennessee Eastman challenge process. J. Process Control 2006, 16, 601–613. (25) 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. (26) Ge, Z. Q.; Song, Z. H. Process monitoring based on independent component analysis-principal component analysis (ICA-PCA) and similarity factors. Ind. Eng. Chem. Res. 2007, 46, 2054–2063. (27) Hyvarinen, A. A fast fixed-point algorithm for independent component analysis. Neural Comput. 1997, 9, 1483–1492. (28) Dunia, R.; Qin, S. J. Subspace approach to multidimensional fault identification and reconstruction. AIChE J. 1998, 44, 1813–1831. (29) Downs, J. J.; Vogel, E. F. A plant-wide industrial process control problem. Comput. Chem. Eng. 1993, 17, 245–255. (30) 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. (31) Lyman, P. R.; Georgakist, C. Plant-wide Control of the Tennessee Eastman Problem. Comput. Chem. Eng. 1995, 19, 321–331. (32) Song, K.; Wang, H. Q.; Li, P. PLS-based optimal quality control model for TE process. Proceedings of the 2004 IEEE International Conference on Systems, Man and Cybernetics, The Hague, The Netherlands; IEEE: Los Alamitos, CA, 2004; pp 1354-1359.

ReceiVed for reView June 14, 2008 ReVised manuscript receiVed March 10, 2009 Accepted March 24, 2009 IE800935E