Fault Diagnosis for a Multiblock Batch Process Based on

Therefore, in this paper, an intermediate block dependency analysis (IBDA) reconstruction algorithm is proposed, based on which normal information is ...
1 downloads 11 Views 2MB Size
Article pubs.acs.org/IECR

Fault Diagnosis for a Multiblock Batch Process Based on Intermediate Block Dependency Analysis Reconstruction Rongrong Sun† and Yingwei Zhang*,‡ †

College of Information Science and Engineering and ‡State Key Laboratory of Synthetical Automation for Process Industries, Northeastern University, Liaoning 100819, P. R. China ABSTRACT: For multiblock batch process, the relationship of variables between two adjacent blocks was not fully considered by the conventional fault diagnosis method. Therefore, in this paper, an intermediate block dependency analysis (IBDA) reconstruction algorithm is proposed, based on which normal information is separated, and online fault detection along the fault direction is achieved. The main contributions are as follows: (1) IBDA is explored to analyze the relative change of variables between two adjacent blocks, and then the increased subspace, normal subspace, and decreased subspace in each block space are separated; (2) to improve the accuracy of fault correction, the fault directions are accurately selected by calculating the fault amplitude of each direction, as a consequence of which both the responsible variables and the uninformative variables are selected; (3) fault is reconstructed by calculating fault information along fault directions, and a reconstruction-based detection model is given, which is more effective than the conventional method. The advantage and effectiveness of the proposed method are illustrated with penicillin fermentation process.

1. INTRODUCTION Because of the wide application of batch production in modern industrial process, fault detection and diagnosis are becoming an extremely important part of plant safety, plant efficiency, product quality, etc.1−3 Many scholars devoted themselves to the research and already made some achievements.4,5 Multiway statistical process monitoring techniques such as multiway principal component analysis (MPCA)6 and multiway partial least-squares (MPLS)7 have been widely applied to process monitoring and received great success. All the data of a batch are preserved as a statistical sample to establish the MPCA or MPLS model in initial research. However, the different features between the two adjacent blocks are not considered, and a high error rate exists as the future traces of variables need predicting. As a result, Kosanovich et al.8 and Dong et al.9 proposed two MPCA/ nonlinear MPCA models to analyze the block-specific nature of a two-block jacketed exothermic batch chemical reactor. And Lu and others10,11 developed a sub-PCA method. Though the detection effect was power, the dependency of variables between two adjacent blocks was not considered. On the basis of the dependency of variables between two adjacent blocks, the specific changes (increase, decrease and normal) of variables characteristic from one block to an adjacent block can be obtained. Therefore, different detection modeling and statistics can be given on the basis of different features of each subspace. Compared to the case for sub-PCA modeling of different blocks, the dependency of variables between two adjacent blocks can more accurately describe the characteristics © XXXX American Chemical Society

of the different variables. Furthermore, some important information shows that the increased part in residual subspace (RS) plays a role in alarm. However, RS is not considered by sub-PCA method. As a result, the accuracy of detection modeling is decreasing. Again, instead of affirming the cause of fault, the models can only find fault. By assuming that the faulty variables make a high contribution to the output, the method of contribution plot12−14 has been widely used to isolate faulty variables. However, the contribution of faulty variables may be inaccurate for the correlation between variables. To resolve the problem, Dunia15 proposed the variance of reconstruction error (VRE) method to reconstruct fault. Yue and Qin16 generalized the reconstruction method to the fault isolation. Tsai and others17,18 applied the reconstruction method to the batch process. Although the reconstruction method isolates part of the fault information, it may not effectively select fault directions that are critical to the reconstruction of fault data. As a result, some important information (the increased variables in RS) may be lost. The variables in RS are divided into three subspaces (increased subspace, decreased subspace, and normal subspace) based on the relative change between two adjacent blocks. The increased variables in RS may play a role in out-of-control T2. If they persist in RS, the fault Received: May 18, 2016 Revised: October 16, 2016 Accepted: October 21, 2016

A

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research directions may not be accurately extracted. And then, the fault may not be diagnosed.19 With the above considerations, fault reconstruction based on intermediate block dependency analysis (IBDA) is proposed to complete fault diagnosis for a multiblock batch process. First, considering the non-Gauss and nonlinear feature of the batch process, multiway kernel independent component analysis (MKICA)20−23 is presented to extract the independent components (ICs). Compared to principal components (PCs), ICs have higher order and independent statistics, which can increase the isolability of the variables. After that, the relative change of the intermediate block is analyzed and it reveals how the variables change from one block to another. And then, a comprehensive subspace decomposition is completed. Lastly, fault directions and fault amplitude are extracted. As a result, fault data are reconstructed along fault directions to isolate normal information from fault data. The remaining sections of this paper are organized as follows. Section 2 describes the IBDA in detail. Section 3 is developed to demonstrate the reconstruction method of a multiblock batch. An online detection strategy is also given in this section. The experimental illustration is given in section 4. At last, conclusions are drawn in Section 5.

Figure 2. Illustration of intermediate block analysis.

2.1. SSPP Algorithm. For the time-slice data Xk(I×J), the KICA algorithm is applied to obtain the initial time-slice models (KICA model). Then we calculate the SPE statistic at each time and determine the confidence limit (Ctrk). From the beginning of batch processes, we add the next time slice one by one to the existing ones and variable-unfold them within the current time region, Xv,k(IK×J). The KICA algorithm is applied to obtain the KICA model. Similarly, the SPE statistic is calculated and the confidence limit (Ctrv,k) is given. On the basis of Ctrv,k > αCtrk, the time k*, which represents the block division time, is obtained. And α is a constant attached to Ctrk. The time slices before k* are denoted as one sub-block. By recursively repeating the above step of the updated beginning of the batch process, we find the following segments. 2.2. Subspace Separation. We collect the normal process data X ∈ RI×K×J and unfold variable-wise as X ∈ RIK×J. And then, SSPP is applied to obtain Xc,i ∈ RIKc,i×J (c = 1, 2, ..., C; i = 1, 2), where Kc,i is the number of time samples belonging to each half (i = 1 indicates the early half and i = 2 indicates the latter half within each block). We select two data sets (Xc,2 ∈ RIKc,2×J and Xc+1,1 ∈ RIKc+1,1×J) in two adjacent halves of two two adjacent blocks, which contain the same number of variables and maybe different numbers of samples. For simplicity, the block data by variable-wise unfolding normalized time slices are described by Xc and Xc+1, respectively. First, KICA is applied to the observed data to obtain the independent component, Sc ∈ Ro×IKc and Sc+1 ∈ Ro×IKc+1. And then, they are used to establish the dependency between Xc and Xc+1 to get the relative change between them. At last, the more accurate detection model is obtained. Because the useful information on small variance data may be covered by the big variance data, Sc is divided into principal component subspace (PCS) and RS:

2. IBDA USING MKICA Extraction of main information is a critical element to accurately divide subspace. Therefore, in this section, MKICA is applied to extract the ICs. For a multiblock batch process, whenever sample data, X(I×K×J) (I, J, and K represent batch, variable and time respectively), are available, it is first preprocessed using batch unfolding or variable unfolding. For the sample data needing preprocessing, accurate characterization of samples is a critical problem for fault diagnosis; however, errors can exist with batch unfolding in the error estimation of the measured value from the start to the end. Therefore, variable unfolding is used in this paper. As is shown in Figure 1, the three-dimensional array X(I×K×J) is the same

⎧ ST = STL LT + STL LT = S ̂ + E ⎪ c c c c c c ,e c ,e c c ⎨ ⎪ G = S TL ⎩ c c c

(1)

where Sĉ and Ec represent principal and residual, Gc ∈ R and Lc ∈ Ro×Rc are the score and loading in PS, Lc,e ∈ Ro×Rc,e is the loading in RCS, and Rc,e = o − Rc. 2.2.1. Subspace Separation in PCS. To complete the subspace separation, the relationship between c and c + 1 block should be established. So, Sc+1 is projected onto Lc to get the correlative score Gc+1 = Sc+1TLc. And then we define the variance ratio of the score between two blocks as follows: IKc×Rc

Figure 1. Variable unfolding for multiblock batch data.

as an unfolded two-dimensional array Xm(IK×J) in the variable direction. And then, block division is completed by the stepwise sequential phase partition (SSPP) algorithm.24 As a result, the whole process is divided into many blocks and every block is divided into two halves, as shown in Figure 2, where the latter half of the previous block and the early half of the latter block is called the intermediate block. Then, KICA is applied to each block to obtain ICs of each block. Lastly, the dependency between two adjacent blocks is analyzed to obtain more accurate subspace.

RTc , i =

var(Gc + 1(:,i)) var(Gc(:,i))

(i = 1, 2, ..., Rc)

(2)

where var(•) denotes the variance of scores and (:,i) denotes the ith column of Gc+1 and Gc, RTc is an Rc-dimensional vector composed of RTc,i. The larger the RTc,i index is, the more dramatic the effect in T2 as it is calculated by score. Here, a confidence region (an B

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

threshold E CtrU and a lower threshold E CtrL) is defined to get the relative change from c block to c + 1 block. Similar to the detection directions partition method in PCS, the detection direction in RS is divided to three parts which are Lc,e,l(o×Rc,e,l), Lc,e,s(o×Rc,e,s), and Lc,e,n(o×Rc,e,n). The variables along Lc,e,l may play a role in out-of-control T2. And the other two variables along Lc,e,s and Lc,e,n are still residual.19 So, the increased variations in RS reconstructed by Lc,e,l are Sc+1,e,l = STc+1Lc,e,lLTc,e,l. Then, PCA is applied to Sc+1,e,l to obtain the principal systematic variations in order with Rc+1,e,l:

upper threshold TCtrU and a lower threshold TCtrL) is defined to get the relative change from the c block to the c + 1 block.19 (a) When the RTc,i index is within the confidence region, the variations in block c + 1 are similar to those in block c. And they can be detected by the same detection model. We keep Gc(:,i)(IKc×Rc,n) by the RTi index, which is within TCtrU and TCtrL. Correspondingly, the concerned detection directions Lc(:,i) are obtained and they are recorded as Lc,n(o×Rc,n). (b) If the RTc,i index is larger than TCtrU, it means that the variations in block c + 1 cause T2 to increase dramatically. We keep Gc(:,i) (IKc×Rc,l) by the RTi index, which is more than TCtrU. Correspondingly, the concerned detection directions Lc (:,i) are obtained and they are recorded as Lc,l(o×Rc,n). (c) If the RTc,i index is smaller than TCtrL, it means that the variations in block c + 1 cause T2 to decrease dramatically. And the detection model in block c will not cause alarm. We keep Gc(:,i) by the RTi index, which is smaller than TCtrL. Correspondingly, the concerned detection directions Lc(:,i) are obtained and they are recorded as Lc,s(o×Rc,s). On the basis of the above analysis, the three detection directions (Lc,n, Lc,l, Lc,s) are obtained and they are orthogonal with each other. They are used to reconstruct the process variations respectively: ⎧S = ScT+ 1Lc , nLcT, n ⎪ c + 1, n ⎪ ⎨ Sc + 1, l = ScT+ 1Lc , lLcT, l ⎪ T T ⎪S ⎩ c + 1, s = Sc + 1Lc , sLc , s

Sc + 1, e = Sc + 1, e , l Lc + 1, eLcT+ 1, e = Gc + 1, eLcT+ 1, e

where Rc+1,e,l is determined by the cumulative explained variance T T rate and the remaining part Lc,e,lLc,e,l − Lc+1,eLc+1,e is also residual. At last, the residuals in c + 1 block are calculated as

Ec + 1, s = ScT+ 1Lc , e , sLcT, e , s

(8)

Intermediate block dependency and how variables change from one block to another are obtained according to the above analysis. However, how variables change along the time direction in different subspaces and what happened during the operation process are not. 2.3. Online Detection Strategy. For the above problems, detection statistics are described as follows. Because the calculation method of block c is similar to that for block c + 1, only detection statistics of block c + 1 are given. For time slices within the early half: ⎧T 2 = gT −1 gc + 1, n ,new Λ c + 1, n ,new c + 1, n ⎪ c + 1, n −1 T T ⎪ = s new Lc + 1, nΛc + 1, n Lc + 1, nsnew ⎪ ⎪T 2 = gT Λ −1gc + 1, l ,new c + 1, l ,new c + 1, l ⎪ c + 1, l ⎪ = s Tnew Lc + 1, lΛc + 1, l −1LcT+ 1, lsnew ⎪ ⎪T 2 = gT −1 gc + 1, s ,new Λ c + 1, s ,new c + 1, s ⎪ c + 1, s ⎪ −1 T T = s new Lc + 1, sΛc + 1, s Lc + 1, ssnew ⎨ ⎪ −1 2 ⎪ Tc + 1, e = g T gc + 1, e ,new Λ c + 1, e ,new c + 1, e ⎪ −1 T T = s new Lc + 1, eΛc + 1, e Lc + 1, esnew ⎪ ⎪ ⎪ SPEc + 1, s = ecT, s ,new ec , s , new = || s Tnew Lc , e , sLcT, e , s ||2 ⎪ ⎪ SPE = e T c+1 c + 1,new ec + 1,new ⎪ 2 T T T ⎪ = || s T L LT − s T L new c , e , l c , e , l new c + 1, eLc + 1, e + s new Lc , e , nLc , e , n || ⎩

In this way, the data in block c + 1 with increased variations, decreased variations, and unchanged variations are obtained. And then, eigenvalue decomposition is used to the covariance matrix of Sc+1,n, Sc+1,l, and Sc+1,s to obtain the loading, which explain their major systematic variations in order with Rc+1,n, Rc+1,i, and Rc+1,s components, respectively:

(4)

where Rc+1,n = rank(Sc+1,nLc,n), Rc+1,l = rank(Sc+1,lLc,l), and Rc+1,s = rank(Sc+1,sLc,s). 2.2.2. Subspace Separation in RS. The residual model of block c is used to check block c + 1. Sc+1 is projected to the residual of block c to get the correlative residual Ec+1 = STc+1Lc,eLTc,e = ∑Ri=1c,e (STc+1lc,e,ilTc,e,i). And then we define the variation difference between two blocks as follows: 2

(9)

where snew and enew are independent component and residual of new data respectively, Λc+1,n, Λc+1,l, and Λc+1,s are diagonal matrices including the eigenvalue of the covariance matrices Sc+1,n, Sc+1,l, and Sc+1,s, respectively. Sc+1,n = STc+1Lc,nLTc,n, Sc+1,l = STc+1Lc,lLTc,l, and Sc+1,s = STc+1Lc,sLTc,s. Statistics of the latter half of block c are similar to those of the early half of block c + 1. If all the detection statistics are below the confidence limits, the process is operating under the normal condition. Otherwise, the process is operating under a fault status.

2

Δc , e , i = || ScT+ 1lc , e , ilcT, e , i || − || ScTlc , e , ilcT, e , i || (i = 1, 2, ..., Rc , e)

(7)

Ec + 1 = ScT+ 1Lc , e , lLcT, e , l − Sc + 1, e , l Lc + 1, eLcT+ 1, e + ScT+ 1Lc , e , nLcT, e , n

(3)

⎧S = Sc + 1, nLc + 1, nLcT+ 1, n = Gc + 1, nLcT+ 1, n ⎪ c + 1, n ⎪ ⎨ Sc + 1, l = Sc + 1, l Lc + 1, lLcT+ 1, l = Gc + 1, lLcT+ 1, l ⎪ T T ⎪S ⎩ c + 1, s = Sc + 1, sLc + 1, sLc + 1, s = Gc + 1, sLc + 1, s

(6)

(5)

3. IBDA-BASED FAULT RECONSTRUCTION In this section, the fault reconstruction algorithm based on IBDA is given for fault diagnosis. Instead of hypothesizing fault direction, the fault amplitude of each direction is calculated to determine fault direction. And then, the fault reconstruction

where ∥•∥ denotes the calculation of Euclidean length and Δc,e is an Rc,e-dimensional vector composed of Δc,e,i. The larger the Δc,e,i index is, the more it contributes to the out-of-control SPE as the SPE statistic is calculated by the square of residual. Here, a confidence region (an upper C

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research ⎧S T ⎪ c + 1, n ,f = Sc + 1, n ,f Lc + 1, n ,f Lc + 1, n ,f ⎪ = Gc + 1, n ,f LcT+ 1, n ,f ⎪ ⎪ Sc + 1, l ,f = Sc + 1, l ,f Lc + 1, l ,f LcT+ 1, l ,f = Gc + 1, l ,f LcT+ 1, l ,f ⎪ T T ⎪S ⎪ c + 1, s ,f = Sc + 1, s ,f Lc + 1, s ,f Lc + 1, s ,f = Gc + 1, s ,f Lc + 1, s ,f ⎨ ⎪ Sc + 1, e ,f = Sc + 1, e ,f Lc + 1, e ,f LcT+ 1, e ,f = Gc + 1, e ,f LcT+ 1, e ,f ⎪ ⎪E = ScT+ 1,f Lc , e , s ,f LcT, e , s ,f ⎪ c + 1, s ,f ⎪ T T ⎪ Ec + 1,f = Sc + 1,f Lc , e , l ,f Lc , e , l ,f − Sc + 1, e , l , f Lc + 1, e ,f ⎪ LcT+ 1, e ,f + ScT+ 1,f Lc , e , n ,f LcT, e , n ,f ⎩

model of each subspace based on fault direction is determined. At last, a fault diagnosis strategy is given. 3.1. Selection of the Fault Principal Component. The main idea of fault reconstruction can be interpreted by Sf = S + Σξ

(10)

where Sf, S, Σ, and ξ represent independent components of fault data, independent components have been reconstructed (fault direction and fault amplitude, respectively), and subscript f denotes fault. Fault information is isolated by reconstruction, and thus fault data are divided to normal information (S) and fault information (Σξ). For eq 10, we multiply both sides by load (L), which contains the loads of principal subspace (Lp) and residual subspace (Le): ⎧ Sf L = SL + Σξ L ⎨ ⎩Gf = G + Σf ξf

where Sc+1,n,f, Sc+1,l,f, Sc+1,s,f. and Sc+1,e,f are the principal subspaces, Ec+1,s,f and Ec+1,f are the residual subspaces, Gc+1,n,f, Gc+1,l,f, Gc+1,s,f, and Gc+1,e,f are the scores of Sc+1,n,f, Sc+1,l,f, Sc+1,s,f, and Sc+1,e,f, respectively, and Lc+1,n,f, Lc+1,l,f, Lc+1,s,f, and Lc+1,e,f are the loads of Sc+1,n,f, Sc+1,l,f, Sc+1,s,f, and Sc+1,e,f, respectively. 3.2.1. Reconstruction of the Principal Subspace. 3.2.1.1. Reconstruction of Sc+1,n,f. Gc+1,n,f is projected onto Lc+1,n: Sc+1,n,f,r = Gc+1,n,fLTc+1,n. On the basis of section 3.1, the fault principal component of Sc+1,n,f,r can be obtained, which is recorded as Gc+1,n,f,r * (IKc+1×Pn). And the corresponding load is recorded as L*c+1,n,f,r(o×Pn). The fault information on Sc+1,n,f T causing alarm can be written as S*c+1,n,f,r = G*c+1,n,f,rL*c+1,n,f,r . Considering that some data in Lc+1,n,f,r * are useless for alarm, PCA is applied to Sc+1,n,f,r * to obtain the fault direction Lc+1,n,f,r,p(o×qn). So, the reconstruction model of Sc+1,n,f can be written





(11)

where Σf represents fault directions, ξf represents fault amplitude, Gf and G represent principal components of fault data (fault principal component), and principal components have been reconstructed, respectively. As worthy of mention, there are two parts in Gf, Gf,p and Gf,e, where subscript p denotes the principal subspace and e denotes the residual subspace. The reconstructed model (the fault principal component is reconstructed along fault directions) is illustrated by the above equation. Here, we note that there is fault on each principal component and thus Σf denotes the unit matrix. Therefore, ξf can be expressed as ξf, j = || Gf (j) − G(j)||2 = || ΔG(j)||2

(13)

Sc + 1, n ,f = S*c + 1, n ,f + Sc + 1, n ,f, o = Sc + 1, n ,f Lc + 1, n ,f, r , pLcT+ 1, n ,f, r , p + Sc + 1, n ,f, o

j = 1, 2, ..., Z

= Gc + 1, n ,f, r , pLcT+ 1, n ,f, r , p + Sc + 1, n ,f, o

(12)

(14)

where S*c+1,n,f and Sc+1,n,f,o are the fault information and normal information on S c+1,n,f , respectively, and G c+1,n,f,r,p = Sc+1,n,fLc+1,n,f,r,p. 3.2.1.2. Reconstruction of Sc+1,l,f . Gc+1,l,f is projected onto Lc+1,l: Sc+1,f,r = Gc+1,l,fLTc+1,l. On the basis of section 3.1, the fault principal component of Sc+1,l,f,r can be obtained, which is recorded as Gc+1,l,f,r * (IKc+1×Pl). And the corresponding load is recorded as Lc+1,l,f,r * (o×Pl). The fault information on Sc+1,l,f T causing alarm can be written as S*c+1,l,f,r = G*c+1,l,f,rL*c+1,l,f,r . Considering that some data in Lc+1,l,f,r * are useless for alarm, PCA is applied to Sc+1,l,f,r * to obtain the fault direction Lc+1,l,f,r,p(o×ql). So, the reconstruction model of Sc+1,l,f can be written

where ξf,j is the j element of ξf, Z is the magnitude of principal component, ΔG(j) = Gf(j) − G(j), and ΔG(j) is the j column of ΔG. ΔGp(j) = Gf,p(j) − Gp(j) and ΔGe(j) = Gf,e(j) − Ge(j) are the fault principal components of the principal subspace and the residual subspace, respectively. The larger the ξf,j index is, the greater the influence on fault information. And then, the fault principal component is determined on the basis of the accumulation contribution rate. Considering that MKICA is used to complete subspace separation, in a simulation process, to avoid sources with big variance overlap the changes in sources with small variance, the fault principal component is extracted by keeping 90%. 3.2. Fault Reconstruction Model. We collect the fault data Xf ∈ RI×K×J, complete block partition (Xc,f ∈ RIKc×J and Xc+1,f ∈ RIKc+1×J) based on SSPP, and get independent components (Sc,f ∈ RO×IKc and Sc+1,f ∈ RO×IKc+1). The relationship between the fault data (Sc,f and Sc+1,f) and normal data (Sc and Sc+1) needs to be established because alarm is raised by fault information, which is used to reconstruct fault data. Considering that the reconstruction algorithms of c block and c + 1 block are similar, the specific steps of reconstruction of the c + 1 block are elaborated below. On the basis of section 2.1, the fault model of the c + 1 block in principal subspace and residual subspace can be given:

Sc + 1, l ,f = S*c + 1, l ,f + Sc + 1, l ,f, o = Sc + 1, l ,f Lc + 1, l ,f, r , pLcT+ 1, l ,f, r , p + Sc + 1, l ,f, o = Gc + 1, l ,f, r , pLcT+ 1, l ,f, r , p + Sc + 1, l ,f, o

(15)

where S*c+1,l,f and Sc+1,l,f,o are the fault information and normal information on Sc+1,l,f, respectively, Gc+1,l,f,r,p = Sc+1,l,fLc+1,l,f,r,p. 3.2.1.3. Reconstruction of Sc+1,s,f. Gc+1,l,f is projected onto Lc+1,s: Sc+1,s,f,r = Gc+1,s,fLTc+1,s. On the basis of section 3.1, the fault principal component of Sc+1,s,f,r can be obtained, which is recorded as Gc+1,s,f,r * (IKc+1×Ps). And the corresponding load is recorded as Lc+1,l,f,r * (o×Ps). The fault information on Sc+1,s,f T causing alarm can be written as S*c+1,s,f,r = G*c+1,s,f,rL*c+1,s,f,r . D

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

where Gc+1(IKc+1×o) and Gc+1,f(IKc+1×o) are score matrix of Ec+1 and Ec+1,f, respectively, Lc+1(o×o) and Lc+1,f(o×o) are load matrix of Ec+1 and Ec+1,f, respectively. Gc+1,f is projected to Lc+1: Ec+1,f,r = Gc+1,fLTc+1. On the basis of section 3.1, the fault principal component of Ec+1,f,r can be obtained, which is recorded as Gc+1,d,f,r * (IKc+1×Pd). And the corresponding load is recorded as Lc+1,d,f,r * (o×Pd). The fault information on Ec+1,f causing alarm can be written as E*c+1,f,r = T G*c+1,d,f,rL*c+1,d,f,r . Considering that some data in Lc+1,d,f,r * are useless for alarm, PCA is applied to Ec+1,f,r * to obtain the fault direction Lc+1,d,f,r,p(o×qd). So, the reconstruction model of Ec+1,f can be written

Considering that some data in L*c+1,s,f,r are useless for alarm, PCA is applied to Sc+1,s,f,r * to obtain the fault direction Lc+1,s,f,r,p(o×qs). So, the reconstruction model of Sc+1,s,f can be written Sc + 1, s ,f = S*c + 1, s ,f + Sc + 1, s ,f, o = Sc + 1, s ,f Lc + 1, s ,f, r , pLcT+ 1, s ,f, r , p + Sc + 1, s ,f, o = Gc + 1, s ,f, r , pLcT+ 1, s ,f, r , p + Sc + 1, s ,f, o

(16)

where S*c+1,s,f and Sc+1,s,f,o are the fault information and normal information on Sc+1,s,f respectively, Gc+1,s,f,r,p = Sc+1,s,fLc+1,s,f,r,p. 3.2.1.4. Reconstruction of Sc+1,e,f. Gc+1,e,f is projected onto Lc+1,e: Sc+1,e,f,r = Gc+1,e,fLTc+1,e. On the basis of section 3.1, the fault principal component of Sc+1,e,f,r can be obtained, which is recorded as Gc+1,e,f,r * (IKc+1×Pe). And the corresponding load is recorded as Lc+1,e,f,r * (o×Pe). The fault information on Sc+1,e,f T causing alarm can be written as S*c+1,e,f,r = G*c+1,e,f,rL*c+1,e,f,r . Considering that some data in L*c+1,e,f,r is useless for alarm, PCA is applied to Sc+1,e,f,r * to obtain the fault direction Lc+1,e,f,r,p(o×qe). So, the reconstruction model of Sc+1,e,f can be written

Ec + 1,f = E*c + 1,f + Ec + 1,f, o = Ec + 1,f Lc + 1, d ,f, r , pLcT+ 1, d ,f, r , p + Ec + 1,f, o = Gc + 1, d ,f, r , pLcT+ 1, d ,f, r , p + Ec + 1,f, o

where Ec+1,f * and Ec+1,f,o are the fault information and normal information on Ec+1,f, respectively, and Gc+1,f,r,p = Ec+1,fLc+1,d,f,r,p. 3.3. Fault Diagnosis Based on IBDA Reconstruction. On the basis of the analysis of IBDA reconstruction, the reconstruction model of different subspaces is given by calculating fault directions. In this way, IBDA reconstruction algorithm is understood and made good use of for fault diagnosis to judge fault cause. Fault is reconstructed on the premise that the fault has been detected by IPDH. And then, different reconstruction models will be tried to find which one can well correct the alarms of detection statistics where the corresponding faulty variables are marked and corrected for each fault candidate. Whenever new data at the kth time interval xnew ∈ RJ×1 are found, we normal them and calculate them as independent Snew and residual enew. And then, the detection statistics for different subspaces are calculated. The systematic score in four principal subspaces is calculated as

Sc + 1, e ,f = S*c + 1, e ,f + Sc + 1, e ,f, o = Sc + 1, e ,f Lc + 1, e ,f, r , pLcT+ 1, e ,f, r , p + Sc + 1, e ,f, o = Gc + 1, e ,f, r , pLcT+ 1, e ,f, r , p + Sc + 1, e ,f, o

(17)

where S*c+1,e,f and Sc+1,e,f,o are the fault information and normal information on Sc+1,e,f, respectively, Gc+1,e,f,r,p = Sc+1,e,fLc+1,e,f,r,p. 3.2.2. Reconstruction of Residual Subspace. 3.2.2.1. Reconstruction of Ec+1,s,f. Gc+1,s,f is projected onto Lc,e,s: Ec+1,s,f,r = T Gc+1,s,fLc,e,s = STc+1,fLc,e,s,fLTc,e,s, where Gc+1,s,f = STc+1,fLc,e,s,f and it is the score of Ec+1,s,f. On the basis of section 3.1, the fault principal component of Ec+1,s,f, r can be obtained, which is recorded as Gc+1,d,s,f,r * (IKc+1×Pd,s). And the corresponding load is recorded as L*c+1,d,s,f,r(o×Pd,s). The fault information on Ec+1,s,f causing alarm T can be written as E*c+1,s,f,r = G*c+1,d,s,f,rL*c+1,d,s,f,r . Considering that some data in Lc+1,d,s,f,r * are useless for alarm, PCA is applied to Ec+1,d,s,f,r * to obtain the fault direction Lc+1,d,s,f,r,p(o×qd,s). So, the reconstruction model of Ec+1,s,f can be written

gc + 1, n ,f, o ,new = sc + 1, n ,f, o ,new Lc + 1, n = (sc + 1, n ,f,new − s*c + 1, n ,f,new )Lc + 1, n = (sc + 1, n ,f,new − sc + 1, n ,f,new Lc + 1, n , f , r , pLcT+ 1, n ,f, r , p)Lc + 1, n

Ec + 1, s ,f = E*c + 1, s ,f + Ec + 1, s ,f, o

(22)

= Ec + 1, s ,f Lc + 1, d , s ,f, r , pLcT+ 1, d , s ,f, r , p + Ec + 1, s ,f, o =

Gc + 1, d , s ,f, r , pLcT+ 1, d , s ,f, r , p

(21)

+ Ec + 1, s ,f, o

gc + 1, l ,f, o ,new = sc + 1, l ,f, o ,new Lc + 1, l = (sc + 1, l ,f,new − s*c + 1, l ,f,new )Lc + 1, l

(18)

= (sc + 1, l ,f,new − sc + 1, l ,f,new Lc + 1, l ,f, r , pLcT+ 1, l ,f, r , p)Lc + 1, l

where Ec+1,s,f * and Ec+1,s,f,o are the fault information and normal information on Ec+1,s,f, respectively, Gc+1,d,s,f,r,p = Ec+1,s,fLc+1,d,s,f,r,p. 3.2.2.2. Reconstruction of Ec+1,f. Singular value decomposition (SVD) is used to the covariance matrix Ec+1 and Ec+1,f to get the score matrix: ⎧1 T ⎪ Ec + 1Ec + 1 = LcT+ 1Λc + 1LcT+ 1 ⎨I ⎪G = E L ⎩ c+1 c+1 c+1 ⎧1 T T T ⎪ Ec + 1,f Ec + 1,f = Lc + 1,f Λc + 1,f Lc + 1,f I ⎨ ⎪G ⎩ c + 1,f = Ec + 1,f Lc + 1,f

(23)

gc + 1, s ,f, o ,new = sc + 1, s ,f, o ,new Lc + 1, s = (sc + 1, s ,f,new − s*c + 1, s ,f,new )Lc + 1, s = (sc + 1, s ,f,new − sc + 1, s ,f,new Lc + 1, s ,f, r , pLcT+ 1, s ,f, r , p)Lc + 1, s (24)

(19)

gc + 1, e ,f, o , new new = sc + 1, e ,f, o ,new Lc + 1, e = (sc + 1, e ,f,new − s*c + 1, e ,f,new )Lc + 1, e = (sc + 1, e ,f,new − sc + 1, e ,f,new Lc + 1, e ,f, r , pLcT+ 1, e ,f, r , p)Lc + 1, e

(20)

(25) E

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

research and the project application. A typical penicillin fermentation process consists of three major operation blocks: bacteria growth block, penicillin secretory block, and somatic cell block. Again, penicillin in secretory block accounts for about 70−80% of the total. Therefore, shortening the bacteria growth block and somatic cell block, extending the penicillin secretory block, and keeping the biggest growth rate of penicillin production are the key to improving penicillin production. It is a typical multiblock batch process and has been widely used for fault diagnosis, process detection, optimization control, etc. It can be readily implemented for experiments to study the effectiveness and advantage of IBDAbased reconstruction. The flow diagram of the penicillin fermentation process is shown in Figure 4.

The residual in two residual subspaces is calculated as ec + 1, s ,f, o ,new = ec + 1, s ,f,new − e*c + 1, s ,f,new = ec + 1, s ,f,new − ec + 1, s ,f,new Lc + 1, d , s ,f, r , pLcT+ 1, d , s ,f, r , p (26)

ec + 1,f, o ,new = ec + 1,f,new − e*c + 1,f,new = ec + 1,f,new − ec + 1,f,new Lc + 1, d ,f, r , pLcT+ 1, d ,f, r , p

(27)

The detection statistics can be calculated as 2 ⎧T = gc + 1, n ,f, o ,newΛc + 1, n−1g cT+ 1, n ,f, o ,new ⎪ c + 1, n ,f, o ⎪ −1 T 2 ⎪Tc + 1, l ,f, o = gc + 1, l ,f, o ,newΛc + 1, l g c + 1, l ,f, o ,new ⎪ ⎪Tc + 1, s ,f, o 2 = gc + 1, s ,f, o ,newΛc + 1, s−1g cT+ 1, s ,f, o ,new ⎨ 2 ⎪T = gc + 1, e ,f, o ,newΛc + 1, e−1g cT+ 1, e ,f, o ,new ⎪ c + 1, e ,f, o ⎪ 2 ⎪ SPEc + 1,f, o = || ec + 1,fo,new || ⎪ 2 ⎪ SPE c + 1, s ,f, o = || ec + 1,s,fo,new || ⎩

(28)

where Λc+1,n−1, Λc+1,l−1, Λc+1,s−1, and Λc+1,e−1 are a diagonal matrix including the eigenvalue of the covariance matrix of Sc+1,n,f, Sc+1,l,f, Sc+1,s,f, and Sc+1,e,f, respectively. In comparison with statistics detected by IBDA, if all statistics stayed within the control limit,25−27 the fault cause is determined and it is equal to the fault of reconstruction model that brings the out-of-control signals back to normal. Otherwise, the fault is a new fault or combination of the intrinsic fault. The online fault diagnosis procedure is implemented as is shown in Figure 3.

Figure 4. Penicillin fermentation process diagram.

The data used in the experiments are generated by Pensim V2.0. Thirty normal batches are produced as the initial training sample data in the simulation. Each batch of penicillin fermentation lasts about 400 h and time intervals are set as 1 h. Eight process variables as is shown in Table 1 are selected,

4. SIMULATION 4.1. Penicillin Fermentation Process. The penicillin fermentation process, a typical nonlinear, non-Gaussian, and dynamic production process, has important value in academic

Table 1. Process Variables Used in the Modeling no.

process variable

unit

1 2 3 4 5 6 7 8

aeration rate substrate feed rate substrate concentration biomass concentration penicillin concentration culture volume CO2 concentration generate heat

L/h L/h g/L g/L g/L L mmol/L Cal

which, thus, results in the offline modeling data array X(30×8×400). From the trajectories shown in Figure 5, it can be seen that only the aeration rate nearly obeys Gaussian distribution, and other variables obey non-Gaussian distribution. MKICA requires no more than one variable to obey Gaussian distribution. So, the eight variables are well suitable for modeling and detecting. Considering the abnormality starts from the very beginning of each batch cycle, the detection results during the early half of block 1 are presented. The duration of each block of the penicillin fermentation process is 1−50 h (bacteria growth

Figure 3. Fault diagnosis scheme. F

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 5. Trajectories of the sampled variables.

For fault 2, the fault amplitude is shown in Table 3. Similarly with the judgment method of fault 1, the third principal component and fourth principal component are considered the fault principal components of Sc+1,e,f and Ec+1,f; the third principal component is considered the fault principal component of Sc+1,n,f and Sc+1,s,f; and the fourth principal component is considered the fault principal component of Ec+1,s,f. 4.3. Fault Diagnosis Based on Fault Principal Component. From section 4.2, the fault principal components of each subspace for two different faults are given. And then, for each fault, on the basis of eqs 23−27, the principal component and residual of each subspace are obtained. As a result, the detection statistics can be obtained on the bsais of eq 27. To illustrate the effectiveness of the proposed method, the reconstruction model of fault 1 is applied to diagnose fault 1 and fault 2, respectively. Similarly, the reconstruction model of fault 2 is applied to diagnose fault 1 and fault 2, respectively. 4.4. Analysis of Simulation Results. The first normal time (FNT), fault detection rate (FDR), and alarm rate (AR) are checked as the fault detection performance. FNT is defined as the time when three or more continuous samples stay with the control region. FDR is defined as Q1/Q2, where Q1 is the number of faults detected in the testing result and Q2 is the number of all faults. It represents the testing result. AR is defined as Q3/Q4, where Q3 is the number of faults detected in the reconstruction result and Q4 is the number of all sample faults. It represents the reconstruction result. The IBDA-based fault reconstruction algorithm, sub-PCA reconstruction algorithm, VRE algorithm, and artificial immune systems (AIS)28 are applied to the experiment. For the proposed method, as is shown in Figure 6, it can be clearly seen that Te2 and SPEs represent the detection results in PCS and RS, respectively, as they have rapid detection and low alarm rates. From the testing result, all of the detection statistics stayed in the upper the control limit. This demonstrates that the proposed method can effectively monitor fault. However, using the sub-PCA method, the T2 statistic detects the fault with three-sampled time delay, and many samples stay within the control limit in both T2 and SPE statistics, which is shown in Figure 7. Again, the alarm rate of T2 and SPE is higher than that of Te2 and SPEs. For the two faults, fault detection results during the first block are summarized in Table 4 as evaluated by FNT and FDR. The results are calculated by mean and mean absolute deviation (MAD) based on 30 batches for each fault case. From

Table 2. Fault Amplitude of Fault 1 along the Principal Component Direction principal direction

Sc+1,n,f

1 2 3 4

0.035 5.686 4.964 0.946

Sc+1,l,f

Sc+1,s,f

Sc+1,e,f

Sc+1,s,f

Sc+1,f

9.236 158.476 134.564 0.015

0.028 3.254 2.146 0.781

0.564 1.326 59.781 6.158

115.489 5948.453 264.516 0.254

Table 3. Fault Amplitude of Fault 2 along the Principal Component Direction principal direction

Sc+1,n,f

1 2 3 4

23.564 11.235 356.149 0.0123

Sc+1,l,f

Sc+1,s,f

Sc+1,e,f

Ec+1,s,f

Ec+1,f

0.235 0.281 3.159

3.914 0.791 39.481 42.819

0.394 4.186 0.364 348.729

19.641 35.716 6948.261 5694.871

block), 51−290 h (penicillin secretory block), and 291−400 h (somatic cell block), respectively, by the SSPP algorithm. The results are agree with the actual situation of the penicillin fermentation process. The early half of the bacteria growth block lasts 21 h. 4.2. Selection of Fault Principal Component. To prove the validity of the presented method, two kinds of fault data are collected: fault 1 is implemented by increasing the aeration rate sharply from start to 50 h; fault 2 is implemented by increasing the substrate feed rate gradually from start to 50 h. For fault 1 and fault 2, the selection of the fault principal component is a key factor to isolate normal information from fault data. The fault amplitude of fault 1 along the principal component direction is shown in Table 2, where the fault amplitudes of the second principal direction and third principal direction are much more than other amplitudes in Sc+1,n,f, Sc+1,s,f, and Sc+1,e,f. So, the second principal component and third principal component are considered their fault principal components. The fault amplitude of the second principal component is much more than the amplitude of the others in Ec+1,f, the fault amplitude of third principal component is much more than the amplitude of the others in Ec+1,s,f, and as a result, the second principal component and third principal component are considered the fault principal components of Ec+1,f and Ec+1,s,f, respectively. In addition, it is noted that there are no variables in Sc+1,l,f. G

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 6. Fault detection results for one batch from fault 1 by the proposed method.

Table 4, it can be clearly observed that the FNT and FDR of PCS for fault 1 resulting from the proposed method are similar to those resulted from the sub-PCA method, indicating that the proposed method and sub-PCA method have similar effect in fault detection. However, the fault in PCS is not detected by the sub-PCA method for fault 2. In PCS, Table 4 shows FNT by the proposed method is less than those by the sub-PCA method, and FDR by the proposed method is more than those by the sub-PCA method, which means that the proposed method is better to detect fault than the sub-PCA method. For fault diagnosis, Figure 6 shows that nearly all of the detection statistics stay within the control limit using the reconstruction model of fault 1. However, all of the detection statistics exceed the control limit and are even higher than the primary fault when the reconstruction model of fault 2 is applied to detect fault 1. This demonstrates that the proposed

method can effectively diagnose fault. In comparison, with the VRE method, four samples stay in the upper control limit in the T2 statistic, and the SPE statistic diagnoses the fault with a three-sampled time delay, which is shown in Figure 8. Similarly, fault diagnosis results are summarized in Table 5 as evaluated by FNT and AR. The calculation method is similar to fault detection. From Table 5, it can be clearly observed that both FNT and AR by the proposed method are less than those by VRE method. This demonstrates the advantage of the proposed method in fault diagnosis. Fault diagnosis results of the proposed method and AIS are summarized in Table 6 as evaluated by FNT and AR. It can be clearly seen that nearly all of FNT of AIS are smaller than the proposed method. It means that AIS can more quickly diagnose fault. However, from AR, we can see that the proposed method can diagnose the fault more accurately. H

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 7. Fault detection results for one batch from fault 1 by the sub-PCA method.

Table 4. Fault Detection Results for Two Faults Evaluated by FNT and FDR (Mean ± MAD) PCS

RS

proposed algorithm

sub-PCA algorithm

proposed algorithm

sub-PCA algorithm

fault no.

FNT

FDR

FNT

FDR

FNT

FDR

FNT

FDR

1 2

3.0 ± 0.2 5.1 ± 0.4

99.5 ± 1.4% 99.3 ± 3.5%

3.5 ± 0.6

98.6 ± 2.5%

3.0 ± 0.0 3.1 ± 0.2

98.7 ± 1.2% 99.3 ± 1.8%

4.6 ± 1.9 6.3 ± 0.7

97.7 ± 4.1% 42.9 ± 1.8%

5. CONCLUSION

On the basis of the above analysis, clearly, the proposed method can detect and diagnose the fault more accurately and

In the present work, a fault reconstruction algorithm on the basis of IBDA is proposed to complete the fault diagnosis for a multiblock batch process. First, subspace separation is completed on the basis of IBDA where the relative change between two adjacent blocks is analyzed by MKICA. And, then, fault directions of each subspace are extracted to get a

timely in comparison with the conventional sub-PCA method and the VRE method. The comparison result between the proposed method and AIS can illustrate that the proposed method can diagnose the fault more accurately. I

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 8. Fault detection results for one batch from fault 1 by the VRE method.

Table 5. Fault Diagnosis Results for Two Faults Evaluated by FNT and AR (Mean ± MAD) PCS proposed algorithm

RS VRE algorithm

proposed algorithm

VRE algorithm

fault no.

FNT

AR

FNT

AR

FNT

AR

FNT

AR

1 2

3.0 ± 0.2 5.1 ± 0.4

1.7 ± 1.3% 2.1 ± 0.4%

3.4 ± 0.5 6.5 ± 0.8

7.9 ± 4.1% 4.1 ± 1.6%

3.0 ± 0.0 3.1 ± 0.2

0.5 ± 0.3% 0.2 ± 0.1%

5.3 ± 1.4 4.2 ± 0.5

1.4 ± 0.7% 0.6 ± 0.3%

Table 6. Fault Diagnosis Results for Two Faults Evaluated by FNT and AR (Mean ± MAD) PCS

RS

proposed algorithm

AIS algorithm

proposed algorithm

AIS algorithm

fault no.

FNT

AR

FNT

AR

FNT

AR

FNT

AR

1 2

3.0 ± 0.2 5.1 ± 0.4

1.7 ± 1.3% 2.1 ± 0.4%

3.0 ± 0.1 3.0 ± 0.3

1.9 ± 1.2% 3.2 ± 0.8%

3.0 ± 0.0 3.1 ± 0.2

0.5 ± 0.3% 0.2 ± 0.1%

3.0 ± 0.1 3.0 ± 0.1

0.9 ± 0.2% 0.6 ± 0.3%

J

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

(15) Qin, S. J.; Dunia, R. Determining the number of principal components for best reconstruction. J. Process Control 2000, 10, 245− 250. (16) Yue, H.; Qin, S. J. Reconstruction based fault identification using a combined index. Ind. Eng. Chem. Res. 2001, 40, 4403−4414. (17) Tsai, C. S.; Chang, C. T.; Chen, C. S. Fault detection and diagnosis in batch and semi-batch processes using artificial neural networks. Chem. Eng. Commun. 1996, 143, 39−71. (18) Tarifa, E. E.; Scenna, N. J. A methodology for fault diagnosis in batch reactors. Comput. Chem. Eng. 1999, 23, S223−S226. (19) Zhang, Y. W.; Qin, S. J. Fault detection of nonlinear processes using multiway kernel independent component analysis. Ind. Eng. Chem. Res. 2007, 46, 7780−7787. (20) Zhang, J.; Xu, Z. D.; Zhou, X. C.; Zhang, H. B.; Yang, N.; Wu, Y. S.; Chen, Y. B.; Yang, G. S.; Ren, T. T. Loss of expression of MHC class I-related chain A (MICA) is a frequent event and predicts poor survival in patients with hepatocellular carcinoma. Int. J. Clin. Exp. Pathol. 2014, 7, 3123−3131. (21) Ge, Z. Q.; Song, Z. H. Performance-driven ensemble learning ICA model for improved non-Gaussian process monitoring. Chemom. Intell. Lab. Syst. 2014, 130, 115−122. (22) Guo, H.; Li, H. G. On-line Batch Process Monitoring with Improved Multi-way Independent Component Analysis. Chin. J. Chem. Eng. 2013, 21, 263−270. (23) Zhao, C. H.; Sun, Y. X. Step-wise sequential block partition (SSPP) algorithm based statistical modeling and online process monitoring. Chemom. Intell. Lab. Syst. 2013, 125, 109−120. (24) Zhao, C. H.; Gao, F. R. Statistical modeling and online fault detection for multiphase batch processes with analysis of betweenphase relative changes. Chemom. Intell. Lab. Syst. 2014, 130, 58−67. (25) Alcala, C. F.; Qin, S. J. Reconstruction-based contribution for process monitoring with kernel principal feature analysis. Ind. Eng. Chem. Res. 2010, 49, 7849−7857. (26) Qin, S. J.; Dunia, R. Determining the number of principal features for best reconstruction. J. Process Control 2000, 10, 245−250. (27) Li, G.; Qin, S. J.; Zhou, D. Output Relevant Fault Reconstruction and Fault Subspace Extraction in Total Projection to Latent Structures Models. Ind. Eng. Chem. Res. 2010, 49, 9175−9183. (28) Shu, Y. D.; Zhao, J. S. Fault Diagnosis of Chemical Processes Using Artificial Immune System with Vaccine Transplant. Ind. Eng. Chem. Res. 2016, 55, 3360−3371.

reconstruction model for fault diagnosis. Finally, normal information is extracted from the fault data by the reconstruction model, which are established by the relationship between fault data and normal data and fault directions. As worthy of mention, for the online detection fault, it is successively compared with historical fault data in the fault library built by reconstructing different fault directions of intrinsic fault. In the penicillin fermentation process example, two different faults are simulated to illustrate that the proposed methods can effectively distinguish different failures. On the basis of the IBDA, how to improve the accuracy of detection may need further exploration in the future.



AUTHOR INFORMATION

Corresponding Author

*Y. Zhang. E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The work is supported by China’s National 973 program (2009CB320600) and National Natural Science Foundation of China (No. 61325015 and No. 61273163).



REFERENCES

(1) Yu, J.; Qin, S. J. Multiway Gaussian Mixture Model Based Multiphase Batch Process Monitoring. Ind. Eng. Chem. Res. 2009, 48, 8585−8594. (2) Hong, J. J.; Zhang, J.; Morris, J. Progressive multi-block modelling for enhanced fault isolation in batch processes. J. Process Control 2014, 24, 13−26. (3) 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. (4) Kosanovich, K. A.; Dahl, K. S.; Piovoso, M. J. Improved process understanding using multiway principal component analysis. Ind. Eng. Chem. Res. 1996, 35, 138−146. (5) Shu, Y. D.; Ming, L.; Cheng, F. F.; Zhang, Z. P.; Zhao, J. S. Abnormal situation management: Challenges and opportunities in the big data era. Comput. Chem. Eng. 2016, 91, 104−113. (6) Chen, J. H.; Chen, H. H. On-line batch process monitoring using MHMT-based MPCA. Chem. Eng. Sci. 2006, 61, 3223−3239. (7) Liu, Q.; Qin, S. J.; Chai, T. Y. Multiblock Concurrent PLS for Decentralized Monitoring of Continuous Annealing Processes. IEEE Trans. Ind. Electron. 2014, 61, 6429−6437. (8) Kosanovich, K. A.; Piovos, M. J.; Dahl, K. S. Multi-way PCA applied to an industrial batch process. Proc. Am. Control Conf. 1994, 1294−1298. (9) Dong, D.; McAvoy, T. J. Multi-stage batch process monitoring. Proc. Am. Control Conf. 1995, 1857−1861. (10) Lu, N. Y.; Gao, F. R.; Wang, F. L. Sub-PCA modeling and online monitoring strategy for batch process. AIChE J. 2004, 50, 255− 259. (11) Wang, J.; Wei, H. T.; Cao, L. L.; Jin, Q. B. Soft-Transition SubPCA Fault Monitoring of Batch Processes. Ind. Eng. Chem. Res. 2013, 52, 9879−9888. (12) Gertler, J.; Li, W.; Huang, Y.; McAvoy, T. Isolation enhanced principal component analysis. AIChE J. 1999, 45, 323−334. (13) Miller, P.; Swanson, R. E.; Heckler, C. E. Contribution plots: a missing link in multivariate quality control. Appl. Math. And Comp. Sci. 1998, 8, 775−792. (14) Westerhuis, J. A.; Gurden, S. P.; Smilde, A. K. Generalized contribution plots in multivariate statistical process monitoring. Chemom. Intell. Lab. Syst. 2000, 51, 95−114. K

DOI: 10.1021/acs.iecr.6b01923 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX