Dynamic Monitoring Method for Multiscale Fault Detection and

Jul 24, 2002 - called the D statistic, is used to detect faults and events of the current operating condition. Additionally, monitoring of individual ...
0 downloads 0 Views 527KB Size
Ind. Eng. Chem. Res. 2002, 41, 4303-4317

4303

Dynamic Monitoring Method for Multiscale Fault Detection and Diagnosis in MSPC Chang Kyoo Yoo,† Sang Wook Choi,‡ and In-Beum Lee*,† Department of Chemical Engineering and School of Environmental Engineering, Pohang University of Science and Technology, San 31 Hyoja Dong, Pohang, 790-784, Korea

A dynamic monitoring method for multiscale fault detection and diagnosis (MFDD) in the wastewater treatment process (WWTP) is proposed. This method is based on dynamic principal component analysis (DPCA) and the D statistic and the monitoring of individual eigenvalues of generic dissimilarity measure (GDM). DPCA method, where the original measurements are lagged, describes a cause and effect of the process. Then the GDM of the DPCA score value, called the D statistic, is used to detect faults and events of the current operating condition. Additionally, monitoring of individual eigenvalues enables the diagnosis of diverse kinds of fault and disturbance sources. The DPCA and the dynamic GDM of the MFDD can be used as a remedy to the problem of analyzing nonstationary processes. The proposed method was applied to fault monitoring and isolation in simulation benchmark data and real plant data. The simulation results clearly show that the method effectively detects faults in a dynamic, multivariate, and multiscale process. Moreover, the proposed approach not only detects faults but also isolates the sources of them to some degree. 1. Introduction Multivariate statistical approaches to process monitoring, fault detection, and diagnosis have rapidly developed over the past 15 years. Advances in sensor technology and data gathering equipment have enabled the collection of data more frequently and variously from almost any process. These improvements in data collection have led to advances in the type and quantity of information that can be gleaned about processes. Many approaches have been proposed in various fields to the problem of extracting process information from data and then interpreting this information. Until about 15 years ago, univariate control and monitoring systems such as Shewhart charts, cumulative sum (CUSUM) charts, and exponentially weighted moving average (EWMA) charts were the methods of choice, and they were occasionally extended to multivariate processes. However, univariate approaches became inadequate as the number of variables increased. To address this problem, univariate analyses have been recently extended to multivariate analyses. In particular, multivariate statistical projection methods such as principal component analysis (PCA) and partial least-squares or projection to latent structures (PLS) have been widely used in modeling and analysis of data sets.1,2 In addition, multivariate monitoring methods using T2 charts, SPE charts and contribution plots have been used.3-5 There are many review papers6-9 on the chemometric techniques that focus on multivariate process control and monitoring. Monitoring of biological treatment processes is very important because recovery from failures is time* To whom correspondence should be addressed. Telephone: 82-54-279-2274. Fax: 82-54-279-3499. E-mail: iblee@ postech.edu. † Department of Chemical Engineering, Pohang University of Science and Technology. ‡ School of Environmental Engineering, Pohang University of Science and Technology.

consuming and expensive. The problems associated with process recovery arise because most of the changes in biological treatment processes are very slow when the process is recovering from a “bad” state to a “normal” state and from a “bulking” state to a “normal” state. Early fault detection and isolation in the biological process are very efficient since they make it possible to execute corrective action well before the onset of a dangerous situation. Applications of multivariate statistical process control (MSPC) to biological processes have recently drawn great interest from a number of researchers. Krofta et al.10 applied MSPC techniques to dissolved air flotation. Rosen11 adapted multivariate statistics to wastewater treatment monitoring using simulated and real process data, while Teppola12 used the combined approach of multivariate techniques, fuzzy and possibilistic clustering, and multiresolution analysis for wastewater treatment monitoring. However, the biological treatment process has several peculiar features that distinguish it from the processes encountered in chemical and industrial engineering. Above all, it is “nonstationary”, which means that the process changes gradually over time. For example, seasonal variations impose a dynamic patter on the process normal condition. In addition, biological treatment involves many underlying phenomena that take place simultaneously, and it may be difficult to separate the effects of specific phenomena. This feature, namely that multiple phenomena simultaneously affect the data at different time or frequency scales, leads to the description of biological processes as having “multiscale” characteristics. If these characteristics interfere with or mask other time or frequency variations, called the “masking effect” , the situation becomes troublesome because the multiscale variations are enlarged up the confidence limit. This is unfavorable since it can lead to situations in which the monitoring algorithm does not detect events that are already underway in the plant.

10.1021/ie0105730 CCC: $22.00 © 2002 American Chemical Society Published on Web 07/24/2002

4304

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Several methods have been recently suggested to solve these problems; these methods utilize adaptive, dynamic PCA (DPCA) and multiscale analysis based on the DPCA13-15 and multiresolution analysis with wavelets.16-24 Ku et al.13 proposed a formulation of DPCA that extracts time-dependent relationships in the measurements by augmenting the data matrix with time-lagged variables. Bakshi16 suggested the multiscale PCA (MSPCA) technique, which combines the ability of PCA to decorrelate variables by extracting a linear relationship with the ability of wavelet analysis to extract deterministic features and approximately decorrelate autocorrelated measurements. Approximate decorrelation of wavelet coefficients makes MSPCA effective for monitoring the autocorrelated measurements without matrix augmentation or time-series modeling. Luo et al.14,15 suggested a multiscale analysis for the DPCA in order to detect faults and to isolate the effects of noise and process changes from the effects of physical changes in the sensor itself. Kano et al.17 proposed a fault detection method which monitors the degree of dissimilarity based on the distribution of timeseries process data. Rosen and Lennox18 applied and developed the adaptive PCA and the multiresolution analysis of wavelets. Teppola and Minkkinen19,20 suggested several types of multiresolution analysis that used a wavelet-PLS regression model to interpret and scrutinize a multivariate model. Choi et al.21 proposed a generic monitoring algorithm that utilized a modified dissimilarity index in the benchmark simulation. Ying and Joseph22 evaluated the feasibility of sensor fault detection using multifrequency signal analysis of noise. Yoo et al.23 extended the dissimilarity concept monitoring method and applied it in a real plant. In addition, Yoo et al.24 suggested a method for simultaneous PLS modeling and fault isolation in a real wastewater treatment plant using the PLS and generic dissimilarity measure (GDM) method. One notable approach is that of Negiz and Cinar,25 who proposed a MSPC method that utilizes a state space identification technique based on canonical variate analysis (CVA) to solve the dynamic problem. This method takes serial correlations into account during the dimension reduction step, like DPCA, and uses the state variables for computing the monitoring statistic in order to remove the serial correlation. This paper concentrates on formulating a solution to the problem of “dynamic multiscale” characteristics and fault isolation in biological processes. The method developed here uses PCA or DPCA for process analysis and solving the multivariate problem and the D statistic for fault detection and diagnosis. These approaches are organized by combining the DPCA model and dynamic multiscale analysis. In the first approach, the DPCA model is used to remove variable autocorrelation and cross-correlation, and for dimension reduction. In the second approach, dynamic multiscale fault detection and diagnosis (MFDD) using a new D statistic and multiscale analysis for the PCA score matrix is proposed. This approach additionally provides a solution to the nonstationary problem with a dynamic PCA and dynamic GDM. The statistical confidence limit of fault detection and isolation is suggested, and its utility is verified by application to data obtained from a simulation benchmark and from a full-scale biological treatment plant. This paper is structured as follows. In the next section, we explain the proposed method. The first

subsection briefly introduces the basic principles of PCA and the conventional monitoring method. In the second subsection, the GDM and dynamic GDM are described and the MFDD monitoring algorithm is outlined. Simulation results are then presented and discussed. 2. Principal Component Analysis (PCA) and Monitoring Statistics PCA is an optimal dimensionality reduction technique in terms of capturing the variance of the data. It decomposes the data matrix X into the sum of the outer product of vectors ti and pi and the residual matrix, E: l

X ) TPT + E )

tipTi + E ∑ i)1

(1)

where ti is a score vector that contains information about the relation between the samples and pi is a loading vector which contains information about the relation between variables. A score vector is orthogonal and a loading vector is orthonormal. PCA can be performed by several methods such as singular value decomposition (SVD) and nonlinear iterative partial least squares (NIPALS). Latent projection into principal component space reduces the original set of variables to l latent variables (LVs). Usually, not all of the score and loading vectors are required to explain important information in the data. In practice, just a few LVs are often sufficient to explain most of the variations in the data. The proper number of LVs can be determined by calculating the proportion of the total variation explained by selected eivenvalues or by cross validation.8,9 A measure of the variation within the PCA model is given by Hotelling’s T2 statistic. T2 at sample k is the sum of the normalized squared scores and is defined as

T2(k) ) t(k)Λ-1t(k)T

(2)

where Λ-1 is the diagonal matrix of the inverse of the eigenvalues associated with the retained principal components. The confidence limit for T2 is obtained using the F distribution.

Tm,n,R2 )

m(n - 1) F n - m m,n-m,R

(3)

where n is the number of samples in the model, m is the number of principal components, and R is an appropriate level of significance for performing the test which typically takes the value of 0.05 or 0.01 for the warning and action limits, respectively. The portion of the measurement space corresponding to the lowest d - m singular values can be monitored by using the squared prediction error (SPE), also called the Q statistic. The SPE is defined as the sum of squares of each row (sample) of E; for example, for the kth sample in X, x(k)

SPE(k) ) e(k)e(k)T ) x(k)(I - PmPmT)x(k)T (4) where e(k) is the kth row of E, Pm is the matrix of the first m loading vectors retained in the PCA model, and I is the identity matrix. The confidence limit for the SPE can be computed from its approximate distribution

[

]

cRx2Θ2h02 Θ2h0(h0 - 1) SPER ) Θ1 +1+ Θ1 Θ2 1

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4305

1/h0

(5)

where cR is the standard normal deviation corresponding to the upper (1 - R) percentile, λj is the eigenvalue associated with the jth loading vector, d

Θi )

λji ∑ j)m+1

for i ) 1, 2, and 3 and

h0 ) 1 -

2Θ1Θ2 3Θ22

The principal component loadings and detection limits for multivariate statistical process control (MSPC) are computed from data representative of normal operating conditions. The process data obtained under normal operating conditions contain only common cause variations, i.e., a variation in the process that is not due to a fault or disturbance. Confidence limits can be developed around the common cause variation for both the systematic component of the variation and the residual component. The systematic component of the process data, which is described by the process model, is monitored using the D statistic chart. The D statistic, which is similar to the Hotelling’s T2 statistic, is a Mahalanobis distance between new data and the center of the normal operating condition data. The pattern of the residuals is monitored using the Q statistic, which is a summation of the squared residuals of a specific time. The D statistic monitors systematic variations in the latent variable space, while the Q statistic represents variations not explained by the retained LVs. That is, faults in the process that violate the normal correlation of variables are detected in the principle component (PC) subspace by the D statistic, whereas faults that violate the PCA models are detected in the residual space by the Q statistic.8,9,11,12,16 Once a fault is detected by MSPC, multivariate charts greatly enhance the capability for diagnosing assignable causes. The basic tool used for fault isolation in MSPC is the contribution plot. Contribution plots are graphical representations of the contribution of each variable to the deviation of the current operating point from that defined by the PCA model for normal in-control operation. By interrogating the underlying PCA model at the point where an event has been detected, one can extract diagnostic or contribution plots which reveal the group of process variables making the greatest contributions to the deviations in the score and SPE. Although these plots do not enable an unequivocal diagnosis of the cause, they do provide considerable insight into the possible causes of an event and thus greatly narrow the search. These plots are usually very powerful at isolating simple or actuator faults, but are much less successful at clearly isolating complex faults. Given that the only information used in this approach is normal process operating data, the fault isolation capabilities outlined above are the best that can be expected. To improve these capabilities would require additional information, and several approaches using prior fault histories have been proposed. These include integrated solution knowledge based supervisory systems (KBS)

exploiting both data and expert knowledge, contribution plots combined with soft computing techniques (e.g., neural networks and fuzzy and genetic algorithms) and causal model approaches (e.g., parity relation and structured residual generation).27 These approaches all isolate faults through the use some form of fault signature developed from analysis of prior faults. Conventional PCA implicitly assumes that the observations at one time are statistically independent of observations at any past time. That is, it implicitly assumes that the measured variable at one time instance has not only serial independence within each variable series at past time instances but also statistical interindependence between the different measured variable series at past time instances. However, the dynamics of a typical chemical or biological process cause the measurements to be time dependent, which means that the data may have both cross-correlation and autocorrelation. PCA methods can be extended to the modeling and monitoring of dynamic systems by augmenting each observation vector with the previous l observations and stacking the data matrix in the following manner.

[

xkT x T X(l) ) · k-1 ·· xk+l-nT

xk-1T xk-2T ·· · xk+l-n-1T

··· ··· ·· · ···

xk-lT xk-l-1T ·· · xk-nT

]

(6)

where xkT is the m-dimensional observation vector in the training set at time instance k and l is the number of lagged measurements. When PCA is performed on the data matrix in eq 6, a dynamic PCA (DPCA) model is extracted directly from the data.13,15 Note that a statistically justified method such as Akaike’s information criterion (AIC) and subspace identification method can be used for selecting the number of lags l to include in the data. The method for automatically determining the number of lags described by Ku et al.13 was not used in the present work. Experience indicates that a value of l ) 1 or 2 is usually appropriate when DPCA is used for process monitoring.27 The D and Q statistic methods provide reliable tools for detecting that a multivariable process has gone outof-control. However, the conventional MSPM method utilizing T2 and Q statistics does not always function well, because it cannot detect changes of correlation among process variables if the T2 and Q statistics are inside the confidence limits. For example, Bakshi16 tested the PCA monitoring of autocorrelated measurements by introducing a disturbance in the form of a mean shift of small magnitude in the input variables. Although such a disturbance will affect the T2 values, it could not be detected with reasonable confidence from the T2 chart. Bakshi showed that conventional PCA monitoring resulted in too many false alarms. Similarly, Lin et al.28 uncovered a fundamental weak point in conventional T2 statistics when they showed that it was unable to detect a linear shift type fault in a multivariate autoregressive process, where the measured data were corrupted by a linear shift of 10% and 20% slope. Lin et al. also assessed the T2 and Q charts for dynamic PCA monitoring of disturbances involving a linear shift. They found that the Q chart successfully detected the linear shift at the fault and all measurements were well above the confidence limit, whereas the shift could not be detected in the T2 chart. A PCA model, if built

4306

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Figure 1. Motivation example for monitoring concept of a distribution of multivariate time series: (O) old operating condition; (b) new operating condition.

properly, extracts the correlation among variables, and a significant increase in the Q statistic indicates a shift from the normal condition. Since a linear shift in the measurements does not violate the PCA model, it cannot be detected at times by the conventional T2 statistic.16,28 A new algorithm for fault detection and diagnosis is required to cope with these problems. 3. Generic Dissimilarity Measure (GDM) 3.1. Theory. Recently, several indices for measuring the dissimilarity between the time series distributions of two data sets have emerged.17,21,23,24 These indices are based on the idea that a change of process operating condition can be detected by comparing the distribution of successive data sets because the data distribution reflects the corresponding process operating condition. The original dissimilarity measure was based on Karhunen-Loeve (KL) expansion and is identical to PCA.17 It compares the covariance structures of two data sets and represents the degree of dissimilarity between them. In the computational procedure, the variance of a transformed data vector is normalized by its corresponding eigenvalue. The dissimilarity measure therefore considers not the absolute magnitude but the relative magnitude of the variance change and consequently neglects the importance of each transformed variable. To consider the importance of each transformed variable, Choi et al.21 and Yoo et al.24 proposed the generic dissimilarity measure (GDM). This measure compares the covariance of two successive data sets and then represents the degree of dissimilarity between them by considering the importance of each transformed variable. Figure 1 shows a motivation example for monitoring concept of a distribution of multivariate time series. The conventional MSPC methods cannot detect the change of operating condition which occurs within the confidence limit of the conventional MSPC methods. The change can be detected by monitoring the distribution of process data or the sequence of time series data. As shown in Figure 1, faults and outliers outside the confidence limit can easily be detected by most MSPC methods. But in order to detect the process change within the normal operation, monitoring techniques for monitoring the distribution of multivariate data have been suggested.21,23,24 In this section, we call a generic dissimilarity measure (GDM) algorithm to mind. This method remains valid for nonstationary processes and does not depend on the normal distribution. It is used for detecting the existence of disturbances as well as for the isolation of the disturbance type through eigenvalue monitoring.

Figure 2. Moving windows between two successive data sets.

The GDM is divided into two major steps: the training of historical data sets under normal conditions and on-line monitoring of a new data set that may contain abnormal events. During the training phase using a normal data set, the intervals and limits of characteristic values are defined. On the basis of the fact that the covariance matrix of the pooled matrix of two data sets can be decomposed using SVD, the training algorithm of the GDM is as follows. First, build two successive data sets (X1 and X2) with a moving window and normalize them with the sample mean and sample variance. Figure 2 illustrates the concept of window and step size using a moving window of each data matrix Xi. The window size refers to the number of samples in each data set and the step size refers to the monitoring interval. Then, find the sample covariance matrix and apply SVD to it.

1 X TX , i ) 1,2 Ni - 1 i i

Si ) S)

(7)

[ ][ ]

1 X1 N - 1 X2

T

N1 - 1 N2 - 1 X1 X2 ) N - 1 S1 + N - 1 S2 (8)

where Si is the sample covariance of data set i, S is the pooled sample covariance, N1 and N2 are the numbers of samples in data sets 1 and 2, respectively, and N is the total sample number (i.e., N1 + N2). Second, apply SVD to S and transform the data matrix (Xi) to orthogonal variables Yi.

SP ) PΛ Yi )

x

Ni - 1 XP) N-1 i

(9)

x

Ni - 1 T N-1 i

(10)

where P is the loading matrix, Λ is the diagonal matrix, and Yi is the transformed variable. Third, find the sample covariance matrices R1 and R2 of Y1 and Y2 and apply SVD to Ri.

R1 + R2 ) Λ

(11)

Riqji ) λjiqji, i ) 1, 2 and j ) 1, ..., r

(12)

where qji is the loading vector, λji is the eigenvalue, and r is the number of the retained principle components. By combining eqs 11 and 12,

R1qj1 ) λj1qj1 and R2qj1 ) (Λjj - λj1)qj1 (13)

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4307

After these linear transformations, two sample covariances of the transformed matrices share the eigenvectors, and then the eigenvalues satisfy the following equation

λj1 + λj2 ) Λjj

(14)

where λji is the jth eigenvalue in the ith data set and Λjj is the jth eigenvalue in the total data set. The greater the similarity between two data sets, the closer their eigenvalues are to 0.5Λjj. As j increases, Λjj sharply decreases. In general, the first few principal components explain most of the variation of the data sets. Then, the following GDM (D) is defined for measuring the dissimilarity of two data sets. r

4 D)

(

λj ∑ j)1 r

)

Λjj 2

2

(15)

∑ Λjj2

j)1

D has a value between 0 and 1. The more similar two data points are, the closer D is to 0. The more dissimilar two data points are, the closer D is to 1. Finally, find the (1 - R)100% confidence interval of each eigenvalue. On the basis of the data obtained under normal operating conditions, the interval (1-R)100% of the eigenvalues is expressed as

- t(1 - R/2; N - 1)s{λji} + λhji e λji e t(1 - R/2; N - 1)s{λji} + λhji (16) where λhji is the sample mean and s{λji} is the sample variance. Statistically, the percentage of λji values below the limit value is (1 - R)100% and the remainder are above the limit.9 Typically the value of R is 0.05 for the warning limit and 0.01 for the action limit. In the present work, we used these typical values of R, which correspond to a 95% and 99% confidence limit, respectively. The control limits of each index, R, are determined so that the number of samples outside the control limit is 5% and 1% of the total samples. 3.2. Monitoring Algorithm. For the monitoring phase of the GDM concept, the confidence limits of the normal data set are calculated through the step outlined above. The sample representing the current operating condition is scaled by the sample mean and sample variance obtained from the previous steps. Then the GDM and the corresponding eigenvalues are calculated using the method outlined above. The GDM quantitatively evaluates the difference between two successive data sets; it can monitor a distribution of time series data and detect a change of the process operating conditions. During normal operation the GDM is zero; an abrupt increase n the GDM indicated a change in the operating condition and the existence of a disturbance. Having detected a possible disturbance, we focus on the individual variations of eigenvalues representative of different scales. As is widely known, the eigenvalues are measured by the amount of variance described by the SVD decomposition. In general, the variance is thought of as information. Because SVD is in descending order of eigenvalues, the first eigenvalue captures the largest amont of information in the decomposition. That is, the

first eigenvalue captures the greatest amount of variation in the data, and the second eigenvalue captures the second information in the data. In this analysis, it is assumed that eigenvalues with large magnitudes represent the effect of low frequency changes such as a large change in the process or the occurrence of a large and long disturbance. Eigenvalues of intermediate magnitude represent small changes in the operating condition such as a short external disturbance, while eigenvalues of small magnitude represent high-frequency phenomena such as sensor faults and measurement noise. The fault isolation approach therefore provides intelligence on the scale at which a disturbance occurs, and can be used to analyze and interpret the physical cause and effect of disturbances. Hence, the eigenvalue at each scale can discern the dominant dynamics and detect the scale on which a fault or disturbance occurs. Most of the variation is captured by the first few eigenvectors, and so only a small number of eigenvalues need to be considered as a monitoring index. If any eigenvalue exceeds its confidence limit, we know that the process operating condition at that scale has changed and therefore that some kind of event or disturbance has occurred. An important consideration in the application of the proposed method is the determination of appropriate window and step sizes. These quantities should be carefully selected by taking into consideration the process characteristics, where the window size should include the process dynamics and the step size should preserve the adequate monitoring performance. We suggest that the window size should be large in comparison to the time constant of the process and the step size should be small in comparison to the sampling time.23 The present work focuses on dynamic monitoring methods, and for dynamic monitoring we suggest the use of a dynamic GDM method. Dynamic GDM is achieved by adding the lagged transformed value to the pooled matrix of the GDM. The augmented pooled data matrix permits implicit time-series modeling to enchase it with cross-correlation within the variables. The data matrix of the dynamic GDM is constructed as shown below

[

Yi(k) ) ti(k - 1) · · · ti(k - l) ti(k - 2) · · · ti(k - 1 - l) ·· · ti(k - w + 1) ti(k - w) · · · ti(k - w - l + 1)

ti(k) ti(k - 1)

]

(17)

where ti(k) is the p-dimensional score vector of the PCA model at time instance k, l is the number of lagged scores, and w is the window size. As mentioned above, the performance of dynamic monitoring is connected to the number of added lags, which should be determined by considering the tradeoff between performance and complexity. In addition, it is known that the added lag variables cause the confidence limit of the monitoring method to widen. If the dynamic PCA does not tackle the dynamic behavior in the PCA model sufficiently, the dynamic GDM which is constructed by the lagged scores of the PCA model can be used as a supplement. 4. Multiscale Fault Detection and Diagnosis (MFDD) The basic idea of the D statistic is that different faults or events may result in different process measurement

4308

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Figure 3. Methodology for dynamic multiscale fault detection and diagnosis.

values, which could be changed in the subspace spanned by selected PCs. Moreover, it manifestes into different areas of the score space of the PCA and then be projected into the change of data distribution of score value. The D statistic provides indispensable insight into process characteristics such as changes in the operating condition, sensor failures, or process upsets. Also, it is well-known that score values of PCA model are distributed normally than the original variables themselves as a consequence of the central limit theorem (CLT). Thus, we would expect the scores, which are a weighted-sum-like mean, to be distributed approximately normally.12 The proposed MFDD method combines the ability of PCA and dynamic PCA to remove cross- and autocorrelation with the ability of the GDM method to both detect faults and diagnose events. The methodology of the MFDD algorithm is shown in Figure 3. The first step in this methodology is the construction of the dynamic PCA model using normal historical data in order to construct the common-causal model. Then, fault detection for the PCA score values is executed by the proposed D statistic and multiscale analysis of each principal eigenvalue. Figure 4 shows the flow diagram of the proposed MFDD algorithm. As in conventional MSPC, if the D statistic exceeds the confidence limit, a process change has occurred. For each eigenvalue plot, if a specific eigenvalue goes outside the normal region and stays there, the eigenvalue will be statistically significant only when the change first occurs. The change is detected first in the scale of eigenvalue that covers the specific feature representing the abnormal operation. Additionally, if the Q statistic of the residual space exceeds the confidence interval, it indicates that changes that violate the PCA model are detected in the residual space. The monitoring scheme is as follows. The change of the PCA scores is detected using the D statistic; that is, the GDM of the DPCA score obtained using a moving window reflects the corresponding operating condition. Faults are diagnosed by monitoring individual eigenvalues to identify the type of event that has occurred within the PCA model. In addition, faults and disturbances that violate the PCA model are detected in residual space using the Q statistic. After the fault detection, fault isolation using eigenvalue monitoring at each scale is used to determine which scale fault has occurred. The multiscale technique therefore gives information on the scale at which process changes or faults occur and facilitates the analysis of the physical or biological phenomena behind the changes. If each

Figure 4. Flow diagram of the proposed MFDD algorithm.

eigenvalue exceeds its corresponding confidence limit, the current process at that scale is changing and an event is occurring. By monitoring at each scale, we can diagnose diverse process variations and events. Thus, it is possible to distinguish between slow variations (seasonal fluctuations or other long-term dynamics), middle-scale variations (internal disturbances or changes in process operating condition), and instantaneous variations (input disturbances, faults, or sensor noises). The approach proposed here therefore gives us the capability to diagnose and interpret the sources of events and faults. Additionally, the variable contribution plots of conventional MSPC can be used as an independent means of interpreting and isolating faults. These plots have been found to be very useful in many applications, because they provide a fault signature that focuses the operator’s attention on a small subset of the large number of process variables, thereby greatly reducing the number of fault possibilities that the engineer need consider. Although these plots do not provide direct fault isolation, they show which group of variables are highly correlated with the fault; it is then up to the process operator to use his or her insight to provide feasible interpretations.26 The proposed method extracts time-dependent relationships in the measurements with the dynamic PCA, augmenting the data matrix by including lagged original variables. If the dynamic PCA does not tackle the dynamic behavior in the PCA model sufficiently, dynamic GDM, which adds the lagged scores of the PCA model in the moving window, is used to remove the remaining dynamic characteristics of the system. The proposed method is adequate for on-line and dynamic monitoring because it is based on a moving window formulation of the GDM. Moreover, it eliminates false alarms after a process returns to normal operation.

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4309

Figure 5. Process layout of the simulation benchmark. Table 1. Various Disturbances in the Simulation Benchmark disturbance type

disturbances

external internal

storm events decreasing nitrification

internal

decreasing settling velocity

sensor faults

nitrate sensor failure

operation change

setpoint change of DO controller

simulation conditions abrupt change of influent flow rate around samples 850 and 1050 specific growth rate for autotrophs: from 0.5 to 0.3 day-1 in a linear fashion during samples 300-500 settling velocity in a secondary settler: from 250 to 150 mday-1 in a linear fashion during samples 300-500 nitrate sensor noises in the second anoxic tank: noise mean changed from 0 to 1 mg N/l during sample 300-400 DO controller setpoint: changed from 2 to 1 mg/L at sample 300

5. Results and Discussion 5.1. Simulation of the Benchmark Plant. The proposed method was applied to the detection of various faults in data obtained from a benchmark simulation.29 The layout of the simulated plant is shown in Figure 5. The plant consists of a five-compartment bioreactor (6000 m3) and a secondary settler (6000 m3). It combines nitrification with predenitrification, which is most commonly used for nitrogen removal. The first two compartments of the bioreactor are not aerated whereas the others are aerated. The International Association of Water Quality (IAWQ) model No. 1 and a 10-layer onedimensional settler model were used to simulate the biological reactions and the settling process, respectively. The influent data used in the simulation were those developed by COST 624, a working group on the benchmarking of wastewater treatment plants.29 The return sludge flow rate (Qr) was set to 100% of the influent flow rate and internal recirculation (Qa) was controlled using a set point (SNO,ref) of 1.0 mg N/L for the nitrate concentration in the second aerator. The aeration (KLa) in aerators 3 and 4 was set to a constant value of 240 day-1. The dissolved oxygen (DO) concentration in aerator 5 was controlled using a set point of 2.0 mg/L. A basic control strategy was proposed to test the benchmark. The aim of this strategy was to control the DO concentration in the final compartment of the reactor by manipulation of the oxygen transfer coefficient and to control the nitrate level in the last anoxic compartment by manipulation of the internal recycle flow rate. The influent characteristics were available from three data sets based on different weather conditions. Each data set started with 1 week of dry weather, which was followed by 1 week of either dry weather, storm events, or prolonged rain. The weekend effect on influent flow and composition was taken into account. Every control strategy was tested using each of these weather data sets. Most WWTPs are subject to large

diurnal fluctuations in the flow rate and composition of the feed stream. Hence, these biological processes exhibit periodic characteristics. The proposed MFDD method was applied to the detection of various disturbances in process data obtained from the benchmark simulation (see Table 1).21 Two types of disturbance were tested using the proposed method, an external disturbance and an internal disturbance. Two short storm events were simulated as the external disturbance. The internal disturbance was simulated by including a linear decrease in the nitrification rate, a linear decrease in the settling velocity in a clarifier, failures of the nitrate sensor in the second aerator, or a set point change of the DO controller in the final aeration basin. For on-line monitoring, we used the following condition. The variables used to build the X-block in the fault detection were the influent ammonia concentration (SNH,in), influent flow rate (Qin), total suspended solid in aerator 4 (TSS4), DO concentration in aerators 3 and 4 (SO3, SO4), and oxygen transfer coefficient in aerator 5 (KLa5). The window size was 20 samples (5 h), and the step size was 5 samples (1.25 h). The mean and variance values for scaling were calculated from the normal data. All of the examples below used a MFDD monitoring algorithm with a 99% confidence limit. The PCA model is made up of three principal components. After the PCA modeling, two lagged values of each PCA score value are added to the dynamic GDM algorithm of the MFDD method. 5.1.1. External Process Disturbances. Hydraulic disturbances are common in most wastewater treatment plants. We simulate a storm event comprising two storms which suddenly occur after a long period of dry weather. The patterns of the measurement variables during the storm weeks are presented in Figure 6. This example shows how external disturbances appear within the proposed method. Figure 7 shows the monitoring results of the MFDD method during the storm weeks. At around samples 850 and 1050, the times of the first and second storm

4310

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Figure 6. X-block variables during the storm weeks.

Figure 7. MFDD monitoring of an external disturbance: (a) D statistic; (b) Q statistic; (c-f) individual eigenvalue plots.

respectively, the D and Q statistics sharply increase and then decrease. From Figure 7, parts a and b, we know that the external disturbance has occurred within the PCA model and can be explained by the PCA model. The plots of the D statistic and each eigenvalue show repeated spikes, which reflect the daily influent variations. Spikes corresponding to variations in the influent were not found in our previous work using GDM method,21,23 which indicates that the dynamic PCA and

GDM make the proposed method sensitive to nonstationary periodic variations. The two storm events are exactly detected by the change of the first eigenvalue, as depicted in Figure 7c. This indicates that external disturbances such as a high flow rate or load are reflected in the first eigenvalue. Meanwhile, the Q statistic has high peak values that coincide with the two events, indicating that an external disturbance also takes place within the PCA model. This result shows

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4311

Figure 8. MFDD monitoring of the deteriorating nitrification case in the internal disturbance: (a) D statistic; (b) Q statistic; (c-f) individual eigenvalue plots.

that MFDD with the D statistic and confidence limit of each eigenvalue detects the fault and diagnoses the fault scale using the first four eigenvalues. 5.1.2. Internal Process Disturbances. In this section, four internal disturbances are considered in simulations in which dry weather is the normal condition.29 First, we simulate a linearly decreasing nitrification rate in the biological reactor, after which we simulate bulking phenomena in the settling tank. Third, we impose a failure of the nitrate sensor of the second aerator. Finally, a set point change of the DO controller in the final aeration tank is included. Table 1 presents a detailed description of the disturbances. The decreasing nitrification rate in the biological reactor was imposed through a linear decrease in the specific growth rate for the autotrophs, µA. The autotrophic growth rate was linearly decreased from 0.5 to 0.3 day-1 between samples 300 and 500, a type of event that is difficult to detect. The X-block used was the same as that used for the external disturbance. The T2 plot from the PCA monitoring has a delay of about 200 samples in detecting this type of event. As shown in Figure 8a, the deterioration of nitrification is detected in the D statistic at around sample 420, which is 110 samples after the event occurred. The D statistic shows periodic peaks which originate from changes in the diurnal influent load, where strong diurnal changes are observed in the flow rate and composition of the feed waste stream. The D statistic shows a relatively large peak value at sample 420. The decreased nitrification rate changed a distribution of multivariate time series data from sample 300, which shifted to the other operating condition with a different distribution. This internal disturbance within the process is explained

mainly by the first eigenvalue in Figure 8a, which gradually increases in the middle of large process change. On the other hand, the Q statistic (Figure 8b) increases continuously and goes out of bounds, which shows that this event occurs outside of the PCA model. As the second internal disturbance, the settling velocity in the secondary clarifier was linearly decreased between samples 300 and 500. This pattern is similar to the case of the deterioration of nitrification. For this analysis it was necessary to add the measurement of effluent total suspended solid (TSSe) to the X-block used for the external disturbance. Figure 9a shows that the D statistic increases after sample 320. As in the case of the deterioration of nitrification, the D statistic remains constant for 20 samples after the settling velocity begins to decrease. The first eigenvalue gradually increases at about sample 400, indicating that it contributes greatly to the fundamental operating condition change. On the other hand, the Q statistic increases continuously, which shows that this event occurs outside the PCA model. To test the ability of MFDD to detect sensor faults in the simulation, we corrupted the nitrate sensor in the secondary anoxic tank. For this analysis it was necessary to add the measurement of nitrate concentration (SNO,2) to the X-block of the external disturbance. The monitoring results are presented in Figure 10. The fault occurred during the period between samples 300 and 400. Prior to that period, the sensor operated properly except for some sensor noise at around sample 140. The sensor fault is first detected by the fourth eigenvalue and then moves to the third eigenvalue, as shown in Figure 10d-f. Thus, the sensor fault initially resulted in a change along the fourth eigenvalue contribution and was subsequently transferred to the third eigenvalue

4312

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Figure 9. MFDD monitoring of the settler bulking case in the internal disturbance: (a) D statistic; (b) Q statistic; (c-f) individual eigenvalue plots.

Figure 10. MFDD monitoring of sensor faults: (a) D statistic; (b) Q statistic; (c-e) individual eigenvalue plots.

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4313

Figure 11. MFDD monitoring performances of a set point change: (a) D statistic; (b) Q statistic; (c-f) individual eigenvalue plots.

contribution, which was mostly dominant thereafter. In addition, the Q statistic (Figure 10b) sharply increases and then decreases between samples 300 and 350, which shows that this type of event occurs within the PCA model. The final internal disturbance considered was a set point change in which the DO controller set point in the fifth biological reactor was suddenly changed from 2 to 1 mg/L at sample 300. The results of these simulations are shown in Figure 11. The X-block used was the same as that used for the external disturbance. The set point change is detected in the D statistic at around sample 330, which is 30 samples after the event occurred. This event is first detected in the third eigenvalue, and then the first and second eigenvalues slowly decrease in absolute magnitude. On the other hand, the Q statistic exactly catches the moment of the set point change outside the PCA model. The eigenvalue plots of Figure 11 show periodic and nonstationary characteristics. This periodicity in the eigenvalue plots originates from the diurnal influent load, where strong diurnal changes are observed in the flow rate and composition of the feed waste stream. These variables tend to fluctuate considerably over a period; hence, it may not be reasonable to assume the process mean and covariance remain constant over time. Because of this characteristic, conventional MSPC methods can lead to frequent false alarms and missed faults when applied under the assumption that the underlying process is stationary. This nonstationary (periodic) behavior of T2 score values in the PCA model is the cause of false alarms and missed faults due to widened confidence limit. 5.2. Full Scale Wastewater Treatment Plant. Process data were collected from a biological wastewater

Figure 12. Plant layout of full scale coke wastewater treatment process.

treatment plant which treats coke plant wastewater from an iron and steel making plant in Korea. The data contained daily average values measured between January 1, 1998, and November 9, 2000, with a total of 1034 samples. This treatment plant uses a general activated sludge process that has five aeration basins (each of size 900 m3) and a secondary clarifier (1200 m3). The plant layout is shown in Figure 12. The treatment plant has two sources: wastewater arrives either directly from a coke making plant (called BET3) or as pretreated wastewater from an upstream WWTP at another coke making plant (called BET2). The coke-oven plant wastewater is produced during the conversion of coal to coke. This type of wastewater is extremely difficult to treat because it is highly polluted and most of the chemical oxygen demand (COD) originates from large quantities of toxic, inhibitory compounds and coal-derived liquors (e.g., phenolics, thiocyanate, cyanides, polyhydrocarbons, and ammonium). In particular, the cyanide (CN) concentration of coke wastewater is a very important variable among the influent loads.

4314

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Table 2. Process Input/Output Variables in the Real WWTP no.

variable

description

X1 X2 X3 X4 X5 X6 X7 X8 X9

Q2 Q3 CN2 CN3 COD2 CNPS COD3 CODPS MLSS•%E

X10 X11 X12 X13 X14 X15

MLSS•R DOaerator Tinfluent Taerator pHAT SVIAT

flowrate from BET2 flowrate from BET2 cyanide from BET2 cyanide from BET3 COD from BET2 effluent cyanide COD from BET3 effluent COD MLVSS at final aeration basin MLSS in recycle DO at final aeration basin influent temp temp at final aerator pH at final aeration basin solid vol index (SVI) at settler solid vol index at settler

X16 SVIsettler

unit

mean

m3/h 179.4 m3/h 85.53 mg/L 2.455 mg/L 14.35 mg/L 156.4 mg/L 1.49 mg/L 2083 mg/L 143.28 mg/L 1605

SD 15.98 8.726 0.3764 2.491 20.28 0.25 295.5 17.66 409.3

Table 3. Variance Captured by the Dynamic PCA Model PC no.

variance captured (cumulative)

variance explained in cross-validation (cumulative)

1 2 3 4 5 6 7

0.167 0.310 0.393 0.465 0.529 0.582 0.634

0.137 0.262 0.325 0.383 0.432 0.469 0.513

mg/L 7194 3444 mg/L 2.064 0.9979 °C 37.6 2.513 °C 30.74 2.379 7.24 0.22 mg/L 11.42 3.15 mg/L

63.31

21.73

5.2.1. The Dynamic PCA Model. Of the more than 30 variables collected from the WWTP, 16 were selected after considering data redundancy and validity. Because the rest variables were determined from an irregular sampling of 3-7 days, they were excluded from the measurement data set. The values in the data set were mean centered and autoscaled to unit variance. Table 2 presents the variables in the data block along with their means and standard deviations (SDs). Of the variables considered, the mixed liquor suspended solid (MLSS) recycled shows a particularly large variation with respect to its mean. The data were divided into two parts. The first 720 observations were used for the development of the DPCA model. In this period, the WWTP operated in an almost normal state. Several observations that were indicative of an abnormal situation were omitted in order to ensure that the training data represented the normal condition, which was confirmed by the expert knowledge of two operators and a researcher. There are several measures of a treatment management index such as the removal efficiency of COD and CN and food to microorganism (F/M) ratio, but the final judge of the normal operation was dependent on operation experts’ knowledge. The 95% confidence limits of the dynamic PCA and MFDD model were decided on the basis of the normal operating data set spanning 2 years. The remaining 314 observations were used as a test data set in order to verify the proposed method. Two lagged variables of each measurement for the DPCA were added in the training data matrix, where 2-day lags correspond to the average hydraulic retention time (HRT) of the system. The method for automatically determining the number of lags described by Ku et al.13 was not used in the present work. For a design of great simplicity of DPCA modeling, we used an useful system information for determining the lag number in this research, which is a hydraulic retention time (HRT) of aerators and settler of the treatment plant. The DPCA model managed to capture about 63% of the X-block variance by projecting the variables from 48 dimensions to 7 dimensions determined by crossvalidation. The variances captured by the DPCA are listed in Table 3. Figure 13 shows the loading plot of the DPCA in the first two principal component dimensions. Measurement variables in the loading plot are composed of one current variable and two lagged variables. The data in Figure 13 confirm that PCA distinguishes the chemical and biological variables, which occupy the different regions

Figure 13. Loading plot of dynamic PCA.

of the plot and exhibit a well-defined pattern. The group in the top left region in Figure 13 contains the variables COD2, COD3, COD•PS, and T•AT, demonstrating that the COD removal rate is strongly correlated with the COD load and the temperature in the aerators. This is exemplified by the heterotrophic biomass activity effects of the temperature in a biological treatment for the carbonaceous nutrients. These variables are uncontrolled or partially controlled throughout the process and therefore exhibit greater variations. The second group, located in the lower region of the loading plot, contains CN2, CN3, CN•PS, Tinput, Q2, and Q3. This grouping indicates that cyanide removal is affected by the cyanide load, influent flow rate and influent temperature. In particular, microorganisms related to cyanide are counter-correlated with the heterotrophic organisms because cyanide compounds are toxic and inhibitory to the growth of heterotrophs of MLSS. This relationship manifests in the opposing directions of the first and second clusters in the loading plot. Hence, shock loading of cyanides in the wastewater influent causes a deterioration of the biological treatment process. The adverse effects of cyanides have been well established in previous experimental results.23,24,30 The third group, located in the upper right region of the loading plot, is made up of MLSS•%E, MLSS•RS, solid volume indexes (SVI) SVI•%E and SVI•RS, and DO in the aerator. This grouping indicates that the settleability of biomass is strongly related to the amount of microorganisms in the total system (aerator and settler). Additionally, the settleability of SVI is related to the DO concentration in the final aerator, which is verified by previous reports that low DO concentration may invoke pinflock bulking phenomena, and excessively high DO concentrations may lead to a deterioration in

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4315

Figure 14. T2 and SPE charts of dynamic PCA in the full scale WWTP.

sludge quality.30 This suggests that the DO concentration in the aeration tank should be controlled. Figure 14 describes the Hotelling’s T2 and SPE charts based on the dynamic PCA model. The 95% confidence limits of the T2 and SPE statistics of the dynamic PCA model with the normal operating data set spanning 2 years were used. Between samples 110 and 125, the T2 values exceed the confidence limit, indicating that a process change may have occurred in this interval. Such deviations of the T2 value can be used to detect a fault more rapidly than static PCA.23 To identify the most likely cause of the deviation, the contributions from every measurement variable were considered, as shown in Figure 15. In this period, the wastewater treatment plant received influent with a high cyanide and COD load and a small flow rate, i.e., a highly concentrated load. The SPE chart increases at sample 100 and then exceeds the confidence limit. The chart subsequently decreases below the confidence limit at sample 150. After sample 200, it fluctuates above and below the confidence limit. We believe that this fluctuation of the SPE chart originates from the new operating condition, in which the microorganism has adapted to the disturbance and the WWTP is operating normally under the new condition. That is, because of the microorganism’s

Figure 15. Contribution plot of T2 chart during samples 110-125.

adaptation ability and control actions existing in the WWTP that would bring the system to a new steady state after process changes or disturbances occurred, one could assume that the relation with the similar variance and different mean between the old and new steady state approximately hold. Consequently, the PCA detection result would be scenario II (Both T2 and SPE charts exceed the control limit) at the beginning of disturbance and would be scenario I (T2 exceeds its control limit, but SPE does not).27,31 5.2.2. MFDD Monitoring. After the construction of the dynamic PCA model, MFDD was applied to the score matrix (T) of the PCA model. To properly monitor the process change or fault and event, the window size for the GDM was set to 20 samples considering the sludge retention time (SRT) and the step size was set to six samples considering the addition of the two lagged values of each variable for the dynamic PCA and the comparison with the conventional T2 chart. The MFDD method uses the statistic of 95% confidence limit. MFDD of the dynamic PCA score values of the test data set is shown in Figure 16. As shown in Figure 16a, the D statistic has large peaks at samples 90 and 170 and small peaks at samples 230 and 280, indicating that process changes happened at these times. Thus, the D statistic shows more rapid and critical detection ability than the conventional T2 chart. Seven eigenvalues are depicted in Figure 16b-h, each of which corresponds to a specific scale of disturbance. The 95% confidence limits of the seven eigenvalues were decided on the basis of the normal operating data set spanning 2 years. The process change is identified by an increase in the sixth eigenvalue and subsequent gradual variations in the fifth, fourth, third and second eigenvalues. From this pattern, we can see that the second and third eigenvalues are similar to the D statistic. Because the second and third eigenvalues are representative of a middlescale fault, we know that there is gradual change in the process operating condition. In detail, the wastewater treatment plant received influent with a high cyanide and COD load and a small influent flow rate, i.e., a highly concentrated load. This influent reduced the activity of the microorganisms and diminished the settling performance, resulting in a SVI increase in the secondary settler. These variations of microorganism

4316

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002

Figure 16. MFDD monitoring performances in the full scale WWTP. (a) D statistic; (b-h) individual eigenvalue plots.

and sludge floc characteristics which were caused by influent load made the biological process the gradual operating change. In addition, we found that an external disturbance, the large influent load, was transformed into an internal disturbance that changed the process operating region in the activated sludge process. The D statistic shows another deviation from sample 230 to the end of the test data set (August 16 through November 9, 2000). At this time, the wastewater treatment plant was modified, with the addition of facilities and treatment equipment. These changes made it feasible for operators to adjust the operation policy, which increased the MLSS concentration and maintained the DO concentration at a high value. The modifications to the plant caused substantial changes to the process and affected all the process conditions after this time, which manifests as a gradual increase of all of the eigenvalues in Figure 16b-h. This result confirms that the proposed MFDD method is distinctly better than conventional methods for analyzing a multiscale process change in a nonstationary signal of unknown characteristics. The proposed MFDD method is therefore an effective technique for extracting information related to changes in a process operating condition and can be used to localize process faults and events of similar scale. 6. Conclusion In this paper, we have presented a new approach to process monitoring based on the GDM of the PCA score value and multiscale analysis of individual eigenvalues. The proposed strategy is able to detect and isolate the effect of multiscale faults in the chemical and biological process. This is achieved by combining dynamic PCA, the dynamic GDM for the PCA score values, and the monitoring of individual eigenvalues. Application of PCA to the lagged matrix (DPCA) has the effect of removing the major dynamics from the process, resulting in residual, that are often uncorrelated. In addition, the present work shows that D statistics with the GDM algorithm and eigenspace utilization can be used to monitor the current operation condition, detect faults and isolate similar kinds of faults. The simulation

results show that the proposed method efficiently detects faults and process changes in a dynamic, multivariate, and multiscale WWTP. Although the approach described here cannot exactly classify the signals at a particular scale, the tracking of separate eigenvalues allows the approximate isolation of various disturbances. Acknowledgment This work was supported by the Brain Korea 21 project. Nomenclature COD ) chemical oxygen demand D ) D statistic E ) residual error vector EVi ) ith eigenvalue GDM ) generic dissimilarity measure l ) time lags Ni, N ) number of samples P ) loading matrix consisting of eigenvectors of S qji ) loading vector of Ri Qr ) return sludge flow rate Qint ) internal flow rate R ) sample covariance matrix of two transformed data sets (Y1 and Y2) Ri ) sample covariance matrix of the transformed variables of data set Yi r ) dimension of compressed data through PCA S ) sample covariance matrix of total data set Si ) sample covariance matrix of data set i SND ) soluble biodegradable organic nitrogen concentration SNH ) NH4+ + NH3 nitrogen concentration SNO ) nitrate and nitrite nitrogen concentration SO ) dissolved oxygen concentration SRT ) sludge retention time s2(i) ) sample variance of xk s{λji} ) estimated standard deviation of λji T2 ) Hotelling’s T2 statistic Ti ) score matrix TNk ) total sample number TSS ) total suspended solids W ) windows size

Ind. Eng. Chem. Res., Vol. 41, No. 17, 2002 4317 Xc ) newly obtained data set X ) total data set Xi ) data set i xjk ) sample mean of xk XB,H ) active heterotrophic biomass concentration Greek Symbols Λ ) diagonal matrix whose diagonal elements are eigenvalues of S Λjj ) eigenvalue of total data set λji ) eigenvalue of Ri λhji ) sample mean of λji µA ) growth rate of autotrophs

Literature Cited (1) Wold, S; Esbensen, K.; Geladi, P. Principal component analysis. Chemom. Intell. Lab. Syst. 1987, 2, 37. (2) Geladi, P.; Kowalski, B. R. Partial least-squares regression: a tutorial. Anal. Chim. Acta 1986, 185, 1. (3) Kresta, J.; MacGregor, J. F.; Marlin, T. E. Multivariate statistical monitoring of process operating performance. Can. J. Chem. Eng. 1991, 69, 35. (4) Wise, B. M.; Ricker, N. L. Recent advances in multivariate statistical process control: improving robustness and sensitivity. IFAC International Symposium, ADCHEM ‘91, Proceedings; Paul Sabatier University: Toulouse, France, 1991, 125. (5) Wise, B. M.; Richer, N. L.; Veltkamp, D. F.; Kowalski, B. R. A theoretical basis for the use of principal component models for monitoring multivariate processes. Process Control Qual. 1990, 1, 41. (6) MacGregor, J. F.; Kourti, T. Statistical Process Control of Multivariate Processes. Control Eng. Practice 1995, 3, 403. (7) Wise, B. M.; Gallagher, N. B. The Process Chemometrics Approach to Process monitoring and Fault detection. J. Process Contol 1996, 6, 329. (8) Jackson, J. E. A users guide to principal component analysis; Wiley: New York, 1991. (9) Johnson, R. A.; Wichern, D. W. Applied Multivariate Statistical Analysis, 3rd ed.; Prentice Hall: Englewood Cliffs, NJ, 1992. (10) Krofta, M.; Herath, B.; Burgess, D.; Lampman, L. An Attempt to understand Dissolved Air Flotation using Multivariate Analysis. Water Sci. Technol. 1995, 31, 191. (11) Rosen, C. Monitoring wastewater treatment system. MS Thesis, Lund, Sweden, 1998. (12) Teppola, P. Multivariate Process Monitoring of Sequential Process DatasA Chemometric Approach. Ph.D. Thesis, Lappeenranta University, Finland, 1999. (13) Ku, W.; Storer, R. H.; Georgakis, C. Disturbance detection and isolation by dynamic principal component analysis. Chem. Intell. Lab. 1995, 30, 179. (14) Luo, R.; Misra, M.; Qin, S. J.; Barton, R.; Himmelblau, D. M. Sensor Fault Detection via Multiscale Analysis and Nonparametric Statistical Inference. Ind. Eng. Chem. Res. 1998, 37, 1024.

(15) Luo, R.; Misra, M.; Himmelblau, D. M. Sensor Fault Detection via Multiscale Analysis and Dynamic PCA, Ind. Eng. Chem. Res. 1999, 38, 1489. (16) Bakshi, B. R. Multiscale PCA with Application to Multivariate Statistical Process Monitoring. AIChE J. 1998, 44, 1596. (17) Kano, M.; Nagao, K.; Ohno, H.; Hasebe, S.; Hashimoto, I. Dissimilarity of process data for statistical process monitoring. Presented at ADCHEM, Pisa, Italy, 2000. (18) Rosen, C.; Lennox, J. A. Adaptive and Multiscale Monitoring of Wastewater Treatment Operation. Water Res. 2000, in press. (19) Teppola, P.; Minkkinen, P. Wavelet-PLS regression models for both exploratory data analysis and process monitoring. J. Chemometrics 2000, 14, 383. (20) Teppola, P.; Minkkinen, P. Wavelets for scrutinizing multivariate exploratory models through multiresolution analysis. J. Chemometrics 2000, 15, 1. (21) Choi, S. W.; Yoo, C. K.; Lee, I. Generic monitoring system in the biological wastewater treatment process. J. Chem. Eng. Jpn. 2001, in press. (22) Ying, C.-M.; Joseph, B. Sensor Fault Detection Using Noise Analysis. Ind. Eng. Chem. Res. 2000, 39, 396. (23) Yoo, C. K.; Choi, S. W.; Lee, I. Disturbance Detection and Isolation in the Activated Sludge Process. Water Sci. Technol. 2002, 45 (4-5), 217-226. (24) Yoo, C. K.; Lee, I. Modeling and Multiresolution Analysis in Wastewater Treatment Process. J. Environ. Eng. 2002 (in revision). (25) Negiz, A.; Cinar, A. Statistical Monitoring of Multivariable Continuous Process with State Space Models. AIChE J.. 1997, 43, 2002. (26) Yoon, S.; MacGregor, J. F. Statistical and Causal ModelBased Approaches to Fault Detection and Isolation. AIChE J. 2000, 46, 1813. (27) Chiang, L. H.; Russel, E. L.; Braatz, R. D. Fault Detection and Diagnosis in Industrial Systems; Springer-Verlag: London, 2001. (28) Lin, W.; Qian, Y.; Li, X. Nonlinear dynamic principal component analysis for on-line process monitoring and diagnosis. Comput. Chem. Eng. 2000, 24, 423. (29) Pons, M. N.; Spanjers, H.; Jeppsson, U. Towards a benchmark for evaluating control strategies in wastewater treatment plants by simulation. Escape 9; Budapest, 1999. (30) Lee, D. S.; Neural Network Modeling of Biological Wastewater Treatment Processes. Ph.D. Dissertation, POSTECH, Korea, 2000. (31) Wang, H.; Song, Z.; Li, P. Fault Detection Behavior and Performance Analysis of Principal Component Analysis Based Process Monitoring Methods. Ind. Eng. Chem. Res. 2002, in press.

Received for review July 5, 2001 Accepted May 24, 2002 IE0105730