Online Predictive Monitoring Using Dynamic Imaging of Furnaces with

Feb 2, 2011 - The presented MPCA−HMM model is able to extract the sequence of hidden states and generate the discriminatory power of different fault...
0 downloads 0 Views 4MB Size
ARTICLE pubs.acs.org/IECR

Online Predictive Monitoring Using Dynamic Imaging of Furnaces with the Combinational Method of Multiway Principal Component Analysis and Hidden Markov Model Junghui Chen,*,† Tong-Yang Hsu,† Chih-Chien Chen,‡ and Yi-Cheng Cheng‡ †

R&D Center for Membrane Technology, Department of Chemical Engineering, Chung-Yuan Christian University, Chung-Li, Taiwan 320, Republic of China ‡ Green Energy & Environment Research Laboratories, Industrial Technology Research Institute, Hsinchu, Taiwan 310, Republic of China ABSTRACT: Furnace processes play a very important role in modern manufacturing industries. They are complicated transient processes and almost a “black box” for operators. It is rather difficult to diagnose those using classical methods, such as statistical classifications. In this Article, novel predictive video monitoring that utilizes prediction from the hidden Markov model (HMM) and multiway principal component analysis (MPCA) is proposed. MPCA is used to extract the cross-correlation among spatial relationships in the low dimensional space, while HMM constructs the temporal behavior of the sequence of the spatial features. Also, HMM can provide state-based segments, which allow predictive models to monitor signals at different time points. With the future predictions, the progress of the current operation can be tracked under a simple probability monitoring chart that shows the occurrence of the observable upsets in the future. A real furnace system is used to verify the effectiveness of the proposed method.

1. INTRODUCTION Furnaces are widely used in many industrial processes to produce electrical power and manufacture high-volume products, like steel, glass, and cement. One of the important reasons for developing a monitoring and diagnosis technique in the combustor system is the concern of energy conservation and the environmental impact. It is very important to ensure that the combustors are always in the stable and optimal condition for achieving higher efficiency, better economics, and even reducing the pollutant emissions. However, the conventional monitoring technique can only detect the abnormal situation of the operation process. By the time the faults occur, it is too late to fix the process, and the incident will have impacted on the plant efficiency and resulted in environmental sanctions. Thus, a technique based on the observed system condition for signaling future system degradation is required. The abnormal limit violations can be found at the early time before the serious matter takes place so that the follow-up effective trouble shooting and continuous quality improvement can be made.1-4 An automated monitoring method can employ a flame detector to decide the automatic ignition and extinction of combustors.5,6 Nevertheless, it cannot provide phenomena information for combustion diagnosis. Some flame detectors based on ultraviolet and infrared sensing are used to diagnose the combustion state, but they cover only part of the flame and thus provide insufficient flame information. Some laser-based techniques have been developed to diagnose combustion systems with quantitative measurements of combustion phenomena parameters, such as concentration, temperature, velocity, rate of heat release, and pollutant emission, but the complexity and higher price of laser-based systems limit their deployment in the industry.2,3 Regarding the diagnosis methodology, the heat release has been measured to detect the instability r 2011 American Chemical Society

of the dynamic combustion process by EWMA control chart.7 However, once the residuals data violate the assumption of independent and identical normal distribution, this kind of control chart should not be used. Flame imaging has already been investigated in the past. The simple way to diagnose combustors is viewing the flame through an industrial television set mounted on a wall opposite to a burner nozzle. Another similar way is to visually inspect the flame through a peep window on a furnace wall.5 However, even if the color of flame images provides a lot of information about combustion, it is hard to quantify the combustion performance and phenomena just by human experience. With the recent progress of optical sensing and digital image processing techniques, image-based visualization methods have been used to study combustion flame in laboratories and in the industry. Some researchers employed sensitivity methods based on images to quantify the combustion process, and then the results were used for process diagnosis. Yan et al. (2002) investigated the relationship between image features and two combustion process parameters,8 that is, furnace load and mass air flow rate, to improve understanding of complex combustion phenomena. If the features of the flame image become abnormal, the cause can be determined from the above relationship information. Lu et al. (2005) demonstrated that relevance between image features and O2 in flue gas could be used as the reference of diagnosis.9 Huang et al. (1999) also investigated the relationship between flickers of a diffusion flame and the air flow rate. In fact, there is important and abundant information existing in flame images to represent Received: March 18, 2010 Accepted: December 7, 2010 Revised: October 9, 2010 Published: February 2, 2011 2946

dx.doi.org/10.1021/ie100671j | Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research the quality of combustion.10 Most analyses provide basic relationship information for diagnosis. However, it is well-known that the turbulent flame flicker frequency is not constant. Therefore, it is hard to get the good prediction model of combustion parameters to describe the combustion phenomena exactly. In the past, multivariate techniques have been employed to monitor the combustion condition. Multiway principal component analysis (MPCA) was performed on a digital multivariate image for softwood lumber grading.11 The Hotelling’s T2 and DModX statistic were used to detect caking of HF furnaces,12 and only the spatial characteristics were considered. However, turbulent flame varies with time, and the flame is expected to be serially correlated. Because the variations in flame pixels can be considered as random events, hidden Markov models (HMM) are popular in statistical modeling. HMMs are statistical models of sequential data that have proven successful in speech and language modeling.13,14 Most of these processes have a common characteristic such that there is an underlying but unobservable dynamic system that operates with uncertain dynamics, giving rise to some dependent observations. The chain structure of the parametric model characterized by the state transition probability can be adaptively trained using the observed sequential data. Literature does suggest HMM be applicable to carrying out diagnostics.15-18 Thus, HMM can be a good candidate for modeling fault diagnostics of combustions as an uncertain process that gives rise to dependent imaging, the fault features. This diagnosis technique is the subsequent use of the classified type of fault from the image data. In this Article, HMM is used to perform fault diagnosis. The main motivation here is to handle the fault in trend representations. However, the spatial images based on HMM suffer from the requirement of huge storage and computational load, hindering HMM’s practical applications. Furthermore, literature on predictive monitoring is extremely limited.15 It is still in its infancy even though it has been gaining importance in recent years. Recovery from failures due to changes in the combustor can be very time-consuming and expensive. Early detection of developing abnormalities is especially important. Thus, there are significant incentives to develop predictive techniques that are capable of signaling future behaviors early enough. The effective corrective action can then be taken to prevent the incident from actually occurring. In this Article, a novel approach that combines MPCA with HMM as a feature extraction method is presented to prognose online operating combustors. MPCA is used to extract the spatial features of the trajectory in a smaller space, while HMM is used to train temporal behavior using a dynamic probability model. On the basis of the trained HMM, a trend selection scheme that facilitates automatic extraction of temporal features from observation data in terms of the trend model templates is proposed. The future predictions are made after a fault or an abnormal situation has been detected. The predictions could then be calculated for determining when the control limit violation occurs for this sustained disturbance fault. The proposed method can interpret the future predictions and provide precise information about specific faults. This Article is structured as follows. The next section prepares the ground for the proposed method by presenting a brief overview (1) of HMM as a tool for building up finite stochastic models, and (2) of MPCA for decomposing the process variables and extracting the spatial features. The proposed image-based approaches for monitoring combustors are developed in section 3. Their monitoring statistics are

ARTICLE

Figure 1. Basic structure of a hidden Markov model.

also established to online diagnose the combustor operations. The framework for carrying out predictive monitoring using MPCAHMM is developed in section 4. In section 5, the effectiveness of the proposed method is demonstrated through a real combustion experiment. The example investigates the performance of the proposed method, and a discussion is made regarding the suitability for diagnostic and predictive activities. Finally, concluding remarks are made.

2. BACKGROUND TECHNIQUES: HMM AND MPCA 2.1. HMM. HMM is one of the finite stochastic models. It facilitates the study of the sequential time series data, and it has proven successful in speech modeling and character recognition.14,19 In HMM, the distribution that generates a sequential observation depends on the state of an underlying but unobserved Markov process. Let {ok,k = 0,1,...} denote a sequence of observation and {sk,k = 0,1,...} denote a sequence of hidden state. The basic structure of HMM is illustrated in Figure 1. Thus, an HMM is a combination of two processes, a Markov chain, which determines the state at time k, sk = j, and 1 e j e S (S is number of states), and a state-dependent process, which generates the observation ok depending on the current state sk = j. In reality, only the statedependent process {ok,k = 0,1,...} is observed, while the hidden state process {sk,k = 0,1,...} remains unkown. To capture the interaction in the sequence data, O ={o1o2... od} is collected. HMM is trained to fit the distribution probability P(O) of the output O. In Figure 1, there is a set of hidden variables, s = {s1,...,sd}, and a fixed number of possible states (S). The probability P(O) can be written as a sum of terms P(O|s) over all the possible state sequences: X PðOÞ ¼ PðOjsÞ

¼

X s

¼

X s

K

all s

PðOjsÞ 3 PðsÞ

Π Pðok jsk ÞPðsk jsk - 1 Þ

k¼1

ð1Þ

HMM has three sets of parameters. (1) The initial state distribution, Π = {πj}: It controls the prior distribution of the hidden state at the starting point k = 1, πj = P(s1 = j). 2947

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

(2) The transition distribution, P(sk|sk-1): It determines the probability of transitioning to state sk = j1 given sk-1 = j2. The transition distribution (P(sk|sk-1)) is defined as Pðsk ¼ j1 jsk - 1 ¼ j2 Þ ¼ aj1 j2 , 1 e j1 , j2 e L

ð2Þ

The probability of the transition, A = {aj1j2}, is an S  S state transition probability matrix. (3) The observation distribution, P(ok|sk): The form of the observation distribution depends on the observation as well as the model assumptions. The distribution model is often characterized by a mixture of G Gaussian densities and separate parameter conditions on each setting of the hidden state: Pðok jsk ¼ jÞ ¼

G X g ¼1

G X

cjg Nðμjg , Σjg Þ and

cjg ¼ 1, cjg g 0, j ¼ 1, 2, :::, S

where G is the number of mixtures. c is the weight for each mixture component, but it should sum to 1 for each state. Each Gaussian probability density function is an elliptically symmetric density with the mean vector μjg and the covariance matrix Σjg. The observation probability parameters are B = {cjg,μjg,Σjg}. For a detailed mixture of Gaussian distribution, please refer to ref 20. Thus, the HMM model parameters are λ = [A B Π]. Given a set of sequence observations, the training HMM is used to determine the setting of the parameters (λ) that maximizes the probability of the set of training sequences in the model. The parameter estimation of HMM can be approached by maximum likelihood estimation (MLE). Thus, the logarithmic likelihood function is I X

log PðOi Þ

ð4Þ

i¼1

where D = {O1O2...OI} is a set of sequences, and I is the number of sets of observation sequences. MLE can be used to obtain the probability efficiently: λ ¼ arg max Lðλ; DÞ

with P samples and N measured variables. If all the finite dimensional distribution of X is at the normal operating condition, the distribution of X would follow the Gaussian distribution and the covariance matrix Σ. The data matrix X can be decomposed into R linear principal components with R e N if the multivariate data with collinearities exist: R X tr pr þ E ¼ TPT þ E ð7Þ X ¼ r¼1

ð3Þ

g¼1

Lðλ ; DÞ ¼

sets. Through the combinational variables, it yields a better discriminative power. Assume that the measurement data at the normal operation are collected and put in a two-dimensional matrix X(P  N): 2 3 x1, 1 x1, 2 3 3 3 x1, N 6 6 x2, 1 x2, 2 3 3 3 x2, N 7 7 7 ð6Þ X ¼ 6 6 7 4 5 333 xP , 1 xP , 2 3 3 3 xP , N

ð5Þ

λ

The training is a recursive algorithm consisting of two steps: finding the HMM parameters λ that maximize the probability, and computing the joint probability. Each iteration of the MLE algorithm increases L(λ;D), and, generally, the sequence of reestimated parameters λ converges to a local maximum of L(λ;D). The details of the training algorithm can be found in ref 21. Although the HMM model approach performs well in many different fields, it still has problems in the case of a large number of parameters to be trained, particularly in the case of the image problem when we have a moderately large number of observed variables in the images. With the prohibitively large number of parameters, it is impossible to train the model. 2.2. PCA. PCA, one of the multivariable statistical analysis methods, is widely used in the communication theory and image processing. It combines variables linearly and identifies their properties. Thus, it can achieve dimensionality reduction on the data set and compress noisy and correlated measurements into a simpler and smaller informative subspace for measurement data

where P and T are defined as the matrix of principal component loading and the matrix of principal component scores respectively. Now the measurement data can be represented through a set of linear functions, which models the combinational relationship between measurement variables (X) and latent variables (P). E, an error matrix, captures the variations in the error space spanned by the loading vectors associated with the N-R smallest singular values. The subspace spanned by X̂ = TPT and E is referred to as the score space and the error space. Details of PCA for monitoring processes are readily obtained in many books on multivariate statistics and data analysis.22,23 Note that the PCA model assumes that the statistical data follow Gaussian distribution, constructing the joint probability distribution of all the measurements. Although it extracts the cross-correlation among variable relationships in the low-dimensional space, the temporal behavior of the sequence of the feature variables is not considered.

3. MPCA-BASED HMM MODELS FOR FLAME IMAGES FROM COMBUSTORS The combustors in industrial processes are usually complicated. The turbulent flame flicker frequency is not constant, and it varies with time. This causes flame images to exhibit spatially homogeneous repetitions and temporally dynamic patterns. Because of autocorrelation among image observations, if PCA is applied to the observed images, the statistic distributions of the extracted features may not conform to Gaussian distribution. On the other hand, it is not possible to train the complete image variables directly applied to HMM. For example, the number of observed variables in a single image with (658  492) is more than 971 208 at each time point. The number of parameters is 5 827 260 when HMM has three hidden states. With a prohibitively large number of parameters, it is impossible to train the HMM model. It is desirable and computationally economical to reduce the measurement of the image variables in a lower dimensional space before applying a process monitoring and diagnosis scheme. MPCA-HMM-based process monitoring and diagnosis is achieved by extracting spatial info using PCA and connecting the temporal dynamic pattern using HMM for better fault interpretation and differentiation. 3.1. Extraction of Spatial Features by MPCA. The data set of a sequence of τ images can be arranged as a two-dimensional 2948

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 2. Unfolded structures of the dynamic textures. The original image sequence is unfolded into the trajectory as time evolves.

array as shown in Figure 2. Assume each color image is of the size N1  N2 and the original color encoding is RGB. It can be reshaped into a sequence of column vectors of the length M = N1  N2  3. With a sequence of τ color images, the column vectors (xk ∈ RM, k = 1,...,τ) are ordered from top to the bottom in a two-dimensional array X = [x1 3 3 3 xτ]T∈RτM, where lower indices indicate that the matrix collects frames from 1 to τ. The temporal mean of each pixel is computed and subtracted from X; then PCA is applied to X to extract the score matrix. 3 2 3 2 t11 t12 3 3 3 t1R tT1 7 6 7 6 7 6 7 6 7 6 T7 6 7 7 6 t t t ¼ ð8Þ t T ¼ 6 k1 k2 kR 7 6 k7 6 7 6 7 6 5 4 5 4 tTτ tτ1 tτ2 3 3 3 tτR tTk

The score vector fixes the weights of the linear combination, so the change of the features of the image frames with time is extracted. 3.2. Extraction of Temporal Features by HMM. When HMM is applied to the flame images from the combustor system, the method for modeling sequential features of image data from PCA compression is a chain structure that captures interaction between a sequence of inputs: ok ¼ ½ o1, k ¼ ½ t1, k

o2, k t2, k

::: :::

or, k

::: oR, k

o R þ 1, k T

tr , k

::: tR, k

vk T

Figure 3. Dynamic training pool of the moving window.

variation at time point k, which is defined as vk ¼ eTk ek

To capture the dynamic behavior, a training pool is constructed by dynamically updating a segment of data throughout the training. Figure 3 shows a graphical illustration of the dynamic training pool. The serial correlations are taken by augmenting each observation vector with the previous d observation and stacking the training data in the following manner: D ¼ fOi ; Oi ¼ foi :::oi þ d - 1 g; i ¼ 1, 2, :::, Lg

ð9Þ

where ok consists of extracted tr,k, r = 1,2,...,R over the first R components of T at time point k and the corresponding residual

ð10Þ

ð11Þ

Assume that I sets of observation sequences are independently sampled from the identical distribution and λ = [A B Π] be the set of all the parameters of the model. The probability P(O) of 2949

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 4. Architecture of the MPCA-HMM-based classification system.

HMM can be written as a sum of terms P(O|s) over all the possible state sequences. 3.3. Online Classification of the Process Status. The main reason for developing the MPCA-HMM model is the subsequent use of the classified trend in fault diagnosis. The automated image defect recognition depends on the use of event signatures for fault diagnosis. In the past, the distribution of defects was estimated by the statistical methods to analyze the image.24,25 However, the assumption of the distribution was restricted. It could not advise process engineers of any specific pattern type to correct the defect combustors. Thus, the image processing is combined with MPCA-HMM to characterize spatial features on combustion images. The online classification of the current process status is shown in Figure 4. During training, the image data from L different event conditions are separately collected to train MPCA-HMM models and determine the operation state opt sequence sopt = {sopt 1 ,...,sd }, where L system conditions include the normal process operation and the operation under L-1 different fault events. Each event will produce a distinct λ. Hence, there are L MPCA-HMM models, λl = [Al Bl Πl], l = 1,2,...,L. Each model corresponds to a system condition. The computation of assessing the current status is rather simple. The architecture of the fault detection system is shown in Figure 4. The collected image series are fed into each MPCAHMM model, which gives a probability (P(O|λl)). The maximum one will correspond to the most probable system condition. If the operation at the current time point is normal, keep monitoring

the next new time point; otherwise, further analyze when the suspected fault condition would be out of the normal operation condition.

4. IMAGE-BASED PREDICTIVE MONITORING SCHEME FOR COMBUSTORS Undoubtedly, the process trend for the abnormal operation carries an intuitive meaning about how the abnormal process behavior changes over time. The objective of predictive monitoring is to predict the progression of the fault condition over the safety limit. This limit is used to define emergency condition, which corresponds to physical process limitations when the operation must be shut down. In the fault diagnosis, the faulty states of the combustor are not observable directly. Given a sequence of the observations, the most likely hidden state sequence for the corresponding faulty condition is chosen. There is a need to make the technique cognizant of the process trend so that such a scenario is very likely to be accommodated. To quantify the process trend of each fault, the data from the faulty condition are applied to the MPCA-HMM model. The most likely state sequence (sopt,l) for each distinct condition (l = 1,2,..., L-1) can be found. Note that the normal condition is not applied here because the state sequences of L-1 fault events are used for our prediction monitoring. It results in the observed sequence, which is equivalent to maximizing sopt, l ¼ arg max PðOjs, λl Þ, l ¼ 1, 2, :::; L - 1 ð12Þ s

2950

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 5. Illustration of the process trend approximated as a sequence of piecewise AR models.

The simple way to find the most likely hidden states of the faulty is a dynamic programming technique termed the Veterbi algorithm. The Veterbi algorithm is a recursive optimal solution to the problem of estimating the state sequence of a finite state Markov chain. To classify the trend curve of each condition (l) for the entire sequence of the measured score, from the optimum sequence of state sopt,l obtained from max P(s|Oλl), the optimum subsequence of states for window d is obtained: opt, l opt, l opt, l ð13Þ sopt, l ¼ ½s1 s2 :::sd  To make the notation simple, the condition index is removed in the following derivations. As the state is retained in each segment, the adjacent states with the same state are collected. The optimal state sequence for the time series scores is used to find the break points of the segment lines. Thus, the entire sequence of states opt opt [sopt 1 s2 ... sd ] can be reduced into the sequences of the segments {h1,h2,...,hH}: h1  fs1 , s2 , :::, sk1 g h2  fsk1 þ 1 , sk1 þ 2 , :::, sk1 þ k2 g

analysis of AR models is of interest in particular for highdimensional systems, estimating coefficient matrices from highdimensional time series data is not an easy task. Because the elements of the observation ok are obtained from MPCA decomposition, the coupling effect in the AR system can be overcome. Now, Rþ1 univariate AR models are separately constructed: j

j

ARj, r : or, k þ 1 ¼ a0 wj ðok Þor, k þ a1 wj ðok - 1 Þor, k - 1 þ ::: j

þ ak - dr wj ðok - vj Þor, k - vj when k∈½kj kj þ kj þ 1 

hH  fsk1 þ ::: þ kH - 1 þ 1 , sk1 þ ::: þ kH - 1 þ 2 , :::, sk1 þ k2 þ ::: þ kH g ð14Þ where hj, j = 1,2,...,H, is the sequence of the same states shown in Figure 5. H is the number of segments. Note that all the window data may not follow the same best path. Different selected states for different windows account for the fact that the observed values of the process variables do not perfectly conform to the deterministic model, but the operating process varies in a stable manner around a fixed state. This natural variability is resulted from the uncertain variations and disturbances among the lurking variables that affect the system. To select the most representative states at each segment, the trained model is used to find the highest frequency of the occurrence state. With the trend extraction procedure based on state durations for each window record {o1o2...okj}, the behavior of the system is decomposed into H disjoint regions {hj}j=1,2,...,H, and each region consists of the same hidden states. The trend in each region (hj) over [kjkjþkjþ1] can be approximated by a local autogressive (AR) model defined by ARj : ok ¼

Figure 6. Local segment validity.

ð16Þ

where the order of each segment (vj) is obtained using the Akaike’s information criterion.26 To characterize the importance of local transitions at different segments, the validity observation, wj(ok) g 0, is used to indicate the relative validity of the local segment with the state hopt i : wj ðok Þ ¼

bj ðok Þ H P j¼1

ð17Þ

bj ðok Þ

where bj(ok) is the probability of observation when the observations are ok under the state with hj. Like the conventional HSMM, it is a Gaussian distribution characterized by the mean vector (μj) and the diagonal covariance matrix (Σj) after the MPCA-HSMM model is trained. The relative validity (wj(ok)) is shown in Figure 6. Based on the observation probability, the relative validity of the local segment is chosen to give a reasonable amount of overlap, because points (which are marked in red shown in Figure 6) in the area of overlap between two adjacent segments will have a well-interpolated overall output between the outputs of overlapping segments. To ensure completeness of the model, a global predicted model can be formed: H X ARj, r or, k þ 1 ¼ j¼1

j A 1 ok - 1

j þ A 2 ok - 2

þ ::: þ A jvj ok - vj

when k∈½kj kj þ kj þ 1 

¼

ð15Þ

H X j¼1

where Ajm ∈ R(Rþ1)(Rþ1) are the coefficient matrices, and vj is the length of moving window of the jth AR model. Because the

j

j

a0 wj ðok Þor, k þ a1 wj ðok - 1 Þor, k - 1 þ ::: j

þ ak - dr wj ðok - vj Þor, k - vj 2951

ð18Þ

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

It requires a higher degree of computation complexity, but it is advantageous in that points in the areas of overlap between two segments. Therefore, the output may lead to smoother changes. The identification determines parameters ajs,j = 1,2,...,H and s = 1,2,...,vj, based on a standard least-squares method without much difficulty. With the measured vector (ok), the future values at the next time points can be predicted using eq 18. The corresponding statistics is then calculated for checking when the future behavior of the combustion furnace violates the confidence bound (i.e., warning limit). An early warning of the abnormal situation occurs. The bound of the warning condition can be calculated on the basis of the developed probability distribution under the normal operation. Past historical image data (Oemg) with the emergency conditions, for example, an operation that runs away and an abnormal pressure relief valve, etc., are collected to determine the shutdown of the operation process. The warning limit (Pth) can be calculated by Pth ¼ PðOemg jsopt , λÞ

ð19Þ

for detecting any departure of the process from its standard behavior. The limit Pth indicates the desired value of the probability density function that will cause the emergency condition. If the prediction observations from eq 18 violate the warning limit, an early warning of the operation situation would occur at this time point; otherwise, the above similar procedure would be repeated for the future prediction and monitoring at the next time point. To sum, in combustion predictive monitoring, MPCA-HMM is used in two phases. The first is the modeling phase. It is operated in the learning mode to establish the model parameters for each operation condition, including the normal operation. With the built MPCAHMM that reveals the most likely sequence of the hidden states, the predicted AR models can then be constructed. The second is the diagnosis phase. In this phase, given a new set of observations that describe the status of the current operation, the log-likelihood of each trained MPCA-HMM (P(O|sopt,λl),l = 1,...,L) is calculated. The model with the maximum probability will correspond to the most probable system condition. This means that the observation is a good match to one of the fault models (l). AR models are then used to

predict the future unknown operation behavior and to check when the P(O|sopt,λ) > Pth statistics at the future time points are out of the control limits in the operation model (λ). If the warning period given by the prediction analysis is long enough, the root cause of the abnormal situation should be analyzed, and a proper corrective action should be made.

5. ILLUSTRATIVE EXAMPLE As for industrial inspection, the application of performance monitoring of the data-driven method to the combustion system using the Markov-based model is rarely reported. In this section, the validity and practicability of the proposed method of MPCAbased HMM will be verified in a flame system. A schematic diagram of the experimental furnace with a flame-monitoring system is shown in Figure 7. The burner of the experimental furnace uses industrial heavy oil as fuels. The air source of combustion comes directly from a direct-driven air compressor. A digital color camera is installed to capture the images of flames in the furnace. The camera with the specifications of 658  492 pixels and a resolution of 24 bits per pixel is capable of capturing flame images at 73 frames per second. The output signal of the camera is sent to a computer in IEEE-1394a interface, and the imaging sample rate is set as one frame per 4 s considering the computer loading. To protect the CCD camera from damage caused by high temperature and keep image qualities from dust, an air-cooling device is used. A set of optical filters is connected to

Figure 7. Scheme of the experimental system.

Figure 8. Time-sequential images where the interval between the images is 3 min: (a) normal; (b) abnormal with a linear decrease in the feed fuel rate; and (c) abnormal with a linear decrease in the feed air rate. 2952

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 9. The testing case 1 monitored by the MPCA-HMM control charts with 95% (dash line) and 99% (solid line) control limits and the warning limit (dotted line) at (a) 590, (b) 650, and (c) 700. The predictions are shown as a red solid dashed line. The corresponding selected hidden state at each time point is shown in each subplot below.

Figure 10. The testing case 2 monitored by the MPCA-HMM control charts with 95% (dash line) and 99% (solid line) control limits and the warning limit (dotted line) at (a) 1570, (b) 1600, and (c) 1650. The predictions are shown as a red solid dashed line. The corresponding selected hidden state at each time point is shown in each subplot below.

the front of the camera to gather more robust flame images and prevent the camera from saturation.

The proposed method is implemented real-time on a laptop and tested for a large variety of conditions in comparison with the 2953

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 11. The original data (dashed line in blue) and the predicted results (solid line in red) of (a) the abnormal condition 1 and (b) the abnormal condition 2.

method utilizing only the color variation information. In this study, three operating conditions of the combustion system are conducted: (1) a normal operation, (2) a linear decrease in the feed fuel rate (fault 1), and (3) a linear decrease in the feed air rate (fault 2). All the images under these operations are recorded and retained in a data storage system, providing a wealth of data

for analysis. Note that the shorter sampling intervals can provide early detection of the process behavior, but the dynamic characteristics would not be negligible because the observation series is dependent on each other. The purpose of the monitoring system is to detect the abnormal operation as early as possible via the image information. 2954

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 12. The sequence of the collected images from the combustion for testing.

Figure 13. Diagnosis result for the first case of the combustion system.

For easy visualization and comparisons, time-sequential images in one normal condition and two abnormal ones are shown in Figure 8. For the abnormal conditions, the first three sequential images are given for the initiating fault, and the last three sequential images are given for the fault occurring in a certain period of time. The interval between the images is 3 min. The differences between the normal (Figure 8a) and the fault behaviors cannot be distinguished from the image data directly at the beginning (Figure 8b and c), but the colors of the image sequences for the normal condition (Figure 8a) are darker than those of the abnormal condition (Figure 8b and c) when the fault

exists for a while. However, the spatiotemporal correlations in the image are not apparent because flame flicker frequency, which is not constant, incurs the variations in flame pixels. For the normal condition, a total of 500 image sequence data is collected when no deterministic disturbances are identified in the operating log. It is used to build the normal MPCA-HMM model. For the other two fault events, based on the identified disturbances in the operating log, the time-sequential image data starting from the normal until the new steady state are separately collected to build up MPCA-HMM models. Each MPCA-HMM corresponds to a training pattern. In each MPCA-HMM, the time-sequential 2955

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

Figure 14. Diagnosis result for the second case of the combustion system.

image data set is unfolded in a two-dimensional array (shown in Figure 1). This matrix is then decomposed by MPCA, and the temporal sequential data are constructed in a two-dimensional array (shown in Figure 3). This matrix is then decomposed by MPCA, and via cross-validation, the number of principal components for normal, fault 1, and fault 2 data are 7, 5, and 6, respectively. The three models separately capture 98.6%, 96.1%, and 97.4% of the variation in the process data set. The features are used to train HMM with a proper window. The control charts for the warning condition are shown in Figures 9 and 10, where the warning limit of log P(O|λ) is calculated from the past emergency data. Also, the control limits representing 95% and 99% confidence levels are shown in Figures 9 and 10. The probably values shown in Figures 9 and 10 are for the testing data; they are different from the ones used in training PCA-HMM models and will be explained later. To characterize the dynamic trend models of the two fault conditions, the AR models for each condition should be constructed. On the basis of the analysis of the state sequence, two segments are needed for both conditions. Figure 11 shows the one-step ahead prediction of the scores of two abnormal cases obtained by AR regression models. The predicted outputs by the trend analysis can significantly represent the essential process behavior change of the abnormal operation. To illustrate the power of MPCA-HSMM, an operation sequence is recorded and shown in Figure 12, while the combustion goes from the normal condition to a fault condition, back to the normal, and finally switches to another fault condition. These two fault conditions, which are similar to two assigned faults of the reference patterns, are generated for testing. Two faults occur at 501 and 1501, respectively. In testing the first case of data, three models (including a normal and two abnormal ones) derived by the above training method are used to monitor the operation condition change. In each condition, the time-sequential

image data are fed into each MPCA-HMM model, which gives a probability; the maximum one will correspond to the most probable system condition. The result is shown in Figure 13. The graph shows the changes of log-likelihood with time in each model. As expected, the result clearly shows that each model has the highest log-likelihood at each corresponding time. Before time point 590, the log-likelihood value of the normal model is the largest; however, the log-likelihood model value for the abnormal model 1 is pretty close to that of the normal model but an increasing trend. Thus, at this time point, the results indicate the first testing case is accurately identified as fault condition 1. For the prediction made at k = 590, the corresponding statistic measure for the future time points can be estimated. Because of the space limitation, Figure 9 only shows the online prediction charts at time k = 590, k = 650, and k = 700. The prediction length for all time is 400. The dashed thick lines represent the predicted future value. The dotted lines represent the warning limits. Because the actual violation occurs at k = 990, the predictions for out-of-control limits are made at k = 984, k = 1000, and k = 999 for the above three starting time points, respectively. The early warning time lengths are Δk = 394, Δk = 350, and Δk = 299. Hopefully, this warning period would be sufficient to allow a corrective action and prevent the occurrence of the significant fault. In testing the second case of data, Figure 14 shows the online diagnosis charts using the proposed method. The probability is within the normal condition at the beginning, but after k = 1570, P(Ok|sk,λ1) of fault model 2 has the highest value. The results indicate the event is identified as fault 2. Figure 10 shows the online prediction charts at time k = 1570, k = 1600, and k = 1650. The proposed model can make the predictions at the above three time points successfully. The violation of the control limit occurs at k = 1737, k = 1746, and k = 1723, respectively. The actual violation occurs at k = 1730. It indicates that the future predictions 2956

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research can be made on the basis of the prediction made earlier. Thus, a limit violation can be predicted. Knowledge about the relationships between measurements and the causes of the faults in many instances can be obtained by reviewing the daily monitoring reports for a combustor. The goal of diagnosis is to focus particularly on those reports detailing when faults occurred and the investigation about the cause of the fault. For this work, only two different faults are picked for demonstration. Although only two different faults with fixed levels are tested here, based on the new daily monitoring record, the faults and the levels of the faults correspond to the degree of deterioration can be designated in the same manner. However, to alleviate training each level of deterioration or the insufficient data of a fault, different levels of the faults can be put in the sample training pool to enhance the discrimination among the cause faults. The investigation of this part will be further explored in our future work.

6. CONCLUSIONS A novel predictive monitoring method for the abnormal situation management of the combustors using image signals is developed to determine when the normal operation will drift away in the future and create the initial fault. The goal is to reduce unnecessary maintenance, while equipment does not cause unexpected work stoppages. In this Article, HMM-based predictive monitoring in videos is proposed. The sequences of images from combustors can be regarded as simply two-dimensional signals with the independent variables repeated in the time frame. The spatial images based on HMM suffer from annoying requirements of huge storage and computation load. A predictive model based on the two-stage feature extraction is developed to increase the robustness and computation efficiency of the spatiotemporal pattern by reducing the dimensionality of the data and retains most of the essential features. MPCA is used as a dimensional reduction technique to project the single frame onto a lower dimensional space. HMM is used to represent the evolution of the extracted feature in the lower dimensional space. The merits of the proposed model can be drawn: • MPCA-HMM provides a stochastic model for combustions in the sequential images. It combines fault classification and the stochastic model to detect variations on the basis of the sequence of images. It can be used to deal with complicated, nonstationary, and stochastic flame images with high flicker frequency. • The presented MPCA-HMM model is able to extract the sequence of hidden states and generate the discriminatory power of different fault events of the process dynamics. The dynamic predictive models are successfully built on the basis of the state segments. • The predictive model is able to provide a more rapid detection of the process fault than other conventional monitoring methods. It allows maintenance time to fix the fault, which may cause serious problems. On the basis of the experimental testing on a combustion operation, the presented MPCA-HMM model has good performance, and it indicates the new method is indeed effective. ’ AUTHOR INFORMATION Corresponding Author

*Tel.: þ886-3-2654107. Fax: þ886-3-2654199. E-mail: jason@ wavenet.cycu.edu.tw.

ARTICLE

’ ACKNOWLEDGMENT We wish to express our sincere gratitude to the Centerof-Excellence (COE) Program on Membrane Technology from the Ministry of Education (MOE), R.O.C., and to the Department of Industrial Technology, Ministry of Economics Affairs, Taiwan, for their financial support. ’ REFERENCES (1) Yan, Y.; Lu, G.; Colechin, M. Monitoring and Characterization of Pulverized Coal Flames Using Digital Imaging Techniques. Fuel 2002, 81, 647. (2) Kohse-Hoinghaus, K.; Barlow, R. S.; Alden, M.; Wolfrum, J. Combustion at the Focus: Laser Diagnostic and Control. Proc. Combust. Inst. 2005, 30, 89. (3) Gang, L.; Yong, Y.; Mike, C. A Digital Imaging Based Multifunctional Flame Monitoring System. IEEE Trans. Instrum. Meas. 2004, 53, 1152. (4) Beer, J. M. Combustion Technology Developments in Power Generation in Response to Environmental Challenges. Prog. Energy Combust. Sci. 2000, 26, 301. (5) Nitsuyo, M.; Kurihara, N. Combustion State Diagnostic Method. Unite States Patent 4555800, 1985. (6) Lee, M. H.; Harashima, F. Flame Detection for the Steam Boiler Using Neural Networks and Image Information in the Ulsan Steam Power Generation Plant. IEEE Trans. Ind. Electron. 2006, 53, 338. (7) Fichera, A.; Pagano, A. Monitoring Combustion Unstable Dynamics by Means of Control Charts. Appl. Energy 2009, 86, 1574. (8) Yan, Y.; Lu, G.; Colechin, M. Monitoring and Characterisation of Pulverized Coal Flames Using Digital Imaging Techniques. Fuel 2002, 81, 647. (9) Lu, G.; Gilabert, G.; Yan, Y. Vision Based Monitoring and Characterization of Combustion Flames. J. Phys.: Con. Ser. 2005, 15, 194. (10) Huang, Y.; Yan, Y.; Lu, G.; Reed, A. On-Line Flicker Measurement of Gaseous Flames by Image Processing and Spectral Analysis. Meas. Sci. Technol. 1999, 10, 726. (11) Bharati, M. H.; MacGregor, J. F.; Tropper, W. Softwood Lumber Grading through On-Line Multivariate Image Analysis Techniques. Ind. Eng. Chem. Res. 2003, 42, 5345. (12) Yoon, S.; Landry, J.; Kettaneh, N.; Pepe, W.; Wold, S. Multivariate Process Monitoring and Early Fault Detection (MSPC) Using PCA and PLS. NPRA Plant Automation and Decision Support Conference, 2003. (13) Bengio, Y. Markovianmodels for Sequential Data. Neural Comput. Surveys 1999, 2, 129. (14) Rabiner, L. R. A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition. Proc. IEEE 1989, 77, 257. (15) Baruah, P.; Chinnam, R. B. HMMs for Diagnostics and Prognostics in Machining Processes. Int. J. Prod. Res. 2005, 43, 1275. (16) Ertunc, H. M.; Loparo, K. A.; Ocak, H. Tool Wear Condition Monitoring in Drilling Operations Using Hidden Markov Models (HMMs). Int. J. Mach. Tools Manuf. 2000, 41, 1363. (17) Zhou, S. Y.; Wang, S. Q. On-line Fault Detection and Diagnosis in Industrial Processes Using Hidden Markov Model. Dev. Chem. Eng. Miner. Process. 2008, 13, 397. (18) Wong, J. C.; McDonald, K. A.; Palazoglu, A. Classification of Abnormal Plant Operation Using Multiple Process Variable Trends. J. Process Control 2001, 11, 409. (19) Ge, M.; Du, R.; Xu, Y. Hidden Markov Model Based Fault Diagnosis for Stamping Processes. Mech. Syst. Signal Process. 2004, 18, 391. (20) Bengio, Y.; Lauzon, V. P.; Ducharme, R. Experiments on the Application of IOHMMs to Model Financial Returns Series. IEEE Trans. Neural Networks 2001, 12, 113. (21) Lee, J. M.; Kim, S. J.; Hwang, Y.; Song, C. S. Diagnosis of Mechanical Fault Signal Using Continuous Hidden Markov Model. J. Sound Vib. 2004, 276, 1065. 2957

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958

Industrial & Engineering Chemistry Research

ARTICLE

(22) Chiang, L. H.; Russell, E. L.; Braatz, R. D. Fault Detection and Diagnosis in Industrial Systems; Springer: London, 2001. (23) Jackson, J. E. A User’s Guide to Principal Components; John Wiley & Sons: New York, NY, 1991. (24) Melzner, H. Statistical Modeling and Analysis of Wafer Test Fail Counts. Advanced Semiconductor Manufacturing 2002 IEEE/SEMI Conference and Workshop, 2002; p 266. (25) Hess, C.; Weiland, L. Extraction of Wafer Level Defect Density Distribution to Improve Yield Prediction. IEEE Trans. Semicond. Manuf. 1999, 11, 175. (26) Kashyap, R. Bayesian Comparison of Different Classes of Dynamic Models Using Empirical Data. IEEE Trans. Autom. Control 1977, 22, 715.

2958

dx.doi.org/10.1021/ie100671j |Ind. Eng. Chem. Res. 2011, 50, 2946–2958