1870
Ind. Eng. Chem. Res. 2003, 42, 1870-1878
PROCESS DESIGN AND CONTROL Improved Threshold for the Local Approach in Detecting Faults Lechang Cheng,† K. Ezra Kwok,*,† Biao Huang,‡ and Ping Li§ Department of Chemical & Biological Engineering, University of British Columbia, 2216 Main Mall, Vancouver, British Columbia, Canada V6T 1Z4, Department of Chemical & Materials Engineering, University of Alberta, Edmonton, Alberta, Canada T6G 2G6, and Institute of Industrial Process Control, Zhejiang University, 310027 Hangzhou, People’s Republic of China
Fault detection and isolation (FDI) has become a critical issue for increasing availability, reliability, and production safety in industrial process monitoring. While model-based FDI methods rely on the mathematical model and input-output data of a process to perform detection, the local approach is a new model-based FDI method that aims to detect slight changes of parametric properties of a system. Robustness with respect to model uncertainties is an important issue for the local approach. To reduce false alarms caused by the estimation error of process parameters, a new algorithm is proposed to recalculate the threshold based on the original threshold and covariance matrix of the estimated parameters. A similar algorithm is also provided to recalculate the threshold to reduce false alarms due to regular parameter fluctuations. The effectiveness of the proposed algorithms is shown by simulation studies on a papermaking crossdirection control system. 1. Introduction With increasing productivity requirements and performance specifications, chemical processes and automatic control systems are becoming more and more complicated. Such conditions increase the possibility of system failures, which can be characterized by unpredictable changes in the system dynamics. Consequently, there is a growing demand for fault tolerance for the purpose of availability, reliability, and production safety. Fault tolerance can be achieved not only by providing hardware redundancy and improving the reliability of individual functional units in the system but also by implementing plant health monitoring systems and online fault detection. Model-based fault detection and isolation (FDI) have drawn growing attention and have become a critical issue for industrial system monitoring in the past several decades.1-4 Among various model-based fault detection methods, a local approach based on system identification and local testing has been developed by Benveniste et al. (1987), Zhang et al. (1994), and Basseville (1998).5-7 The main idea of the local approach is to transfer fault detection problems concerning a parametrized process into a universal problem of monitoring the mean of a Gaussian vector. Benveniste et al. (1987)5 first proposed the basic idea of the local approach in change detection and model validation. Building on this, Huang applied it to the problem of process and control loop performance monitoring.8 * To whom all correspondence should be addressed. Email:
[email protected]. † University of British Columbia. ‡ University of Alberta. § Zhejiang University.
All model-based FDI methods rely on the mathematical model of the monitored system to detect and diagnose. If the system model is accurate, fault detection may be very straightforward; however, noise and model uncertainties are inevitable in a practical system. The consequence is that, even in a fault-free case, the nominal relationship between the system input and output will not be satisfied. This may lead to erroneous decisions on the presence of faults in the system. To improve the reliability and performance of fault detection algorithms, it is of great importance to consider the effects of noise and model uncertainties on the detection algorithms, i.e., the robustness of FDI algorithms. Because a number of methods have been proposed in the stage of residual generation to improve the insensitivity of the residual to model uncertainties, the robustness of the detection algorithms can also be improved by appropriate selection of thresholds to reduce the false alarm rate. From the application point of view, the selection of thresholds is as important as the generation of the residual. In this paper, a new algorithm is proposed to revise the threshold in order that the false alarm rate be reduced to the required level while maintaining the sensitivity of fault detection algorithms toward faults. The result is also extended to the situation of regular system parameter fluctuations. A similar algorithm is provided to calculate the threshold to be used in the residual evaluation. With the revised threshold, the detection algorithm will not generate an alarm under acceptable parameter changes. To test the applicability of the local approach and the proposed algorithm of revising the threshold, we applied them to a papermaking cross-direction (CD) control system. Simulation results show the improvement of the revised threshold over the original threshold.
10.1021/ie010797i CCC: $25.00 © 2003 American Chemical Society Published on Web 04/05/2003
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003 1871
2. Principle of the Local Approach In fault detection, a linear system characterized by the following parametric model is considered:
yk ) φkTθ + ek
(1)
where yk, φk, and ek are the system output, regression variable, and a white noise series at time k, respectively, and θ is the vector of the system’s parameters. A function H(θ0,yk,φk) can be called a valid primary residual if
EθH(θ0,yk,φk) ) 0
if θ ) θ0
(2)
EθH(θ0,yk,φk) * 0
if θ * θ0
(3)
where Eθ means the expectation when the value of system parameter is θ. For a linear system given by eq 1, a valid primary residual can be calculated as
H(θ0,yk,φk) ) φk(yk - φkTθ)
1 H(θ0,yk,φk) xN
(5)
The sensitivity matrix and covariance matrix are given respectively as9
∂ H(θ,yk,φk)|θ)θ0 ∂θ
(6)
Σ(θ0) ) lim Eθ0ζN(θ0) ζN(θ0)T
(7)
M(θ0) ) Eθ0 Nf∞
t ) ζNT(θ0) Σ-1ζN(θ0)
Usually, the nominal parameter θ0 of a system is unknown and is estimated through system identification. Because of noise and disturbance, the estimated parameter θˆ 0 is not equal to θ0. Consequently, the primary residual H(θˆ 0,yk,φk) and normalized residual ζN(θˆ 0) calculated with the estimated parameter θˆ 0 have a nonzero mean even in a fault-free case. With the estimated parameter θˆ 0, the χ2 value of the normalized residual of ζN(θˆ 0) becomes
ht ) ζNT(θˆ 0) Σ-1ζN(θˆ 0)
(8)
For the detection of a small deviation in the system parameter, the local approach is applied to decide between two hypotheses:
H0: θ ) θ0 H1: θ ) θ0 + λ/xN
(9)
If the sensitivity matrix M is of full column rank and the covariance matrix Σ is positive definite, when N f ∞,9 the normalized residual will become
ζN(θ0) ∼ N(0,Σ) ζN(θ0) ∼ N(-Mλ,Σ)
under H0 under H1
(14)
ht is a noncentral χ2-distributed random variable, and the possibility of ht to exceed χR2 is larger than R. Therefore, the threshold needs to be revised so that the false alarm rate can be kept below the specified level. However, the problem remains that the larger the threshold is, the less sensitive the detection algorithm is to faults. To minimize the loss of detection performance, the threshold needs to be set as small as possible. The problem of revising the threshold can be summarized as follows: (1) Use the original scheme of the local approach and the estimated parameter θˆ 0 to calculate the primary residual and the normalized residual and to do a χ2 test. To keep the false alarm rate in the specified level, a new threshold χˆ R2 should be found such that
For the sake of simplicity, we will use the following notations:
M ) M(θ0), Σ ) Σ(θ0)
(13)
3. Revision of the Threshold
(4)
Given a primary residual H(θ0,yk,φk) and an N-size sample, the normalized residual is defined as
ζN(θ0) )
from a standard χ2 table. The χ2 value of the normalized residual t is compared with χR2, and a fault is detected when t is larger than χR2. When the sensitivity matrix M is square and invertible, t can be simplified as
P(th>χˆ R2|H0) e R
(15)
where ht is calculated using eq 14 and H0 means a faultfree case. (2) To keep the sensitivity of the detection algorithm, the threshold χˆ R2 should be set as small as possible. 3.1. Distribution of the Normalized Residual. The first step of finding a suitable threshold is to find the distribution of the normalized residual with the estimated parameter θˆ 0. Theorem 1. For the linear system given by eq 1 and the primary residual calculated with eq 4, the distribution of the normalized residual calculated with the estimated parameters θˆ 0 is as follows:
ζN(θˆ 0) ∼ N(Ma,Σ)
(10)
under H0
ζN(θˆ 0) ∼ N(Ma - Mλ,Σ)
under H1
(16) (17)
(11)
and the generalized likelihood ratio of the normalized residual ζN(θ0)
t ) ζNT(θ0) Σ-1M(MTΣ-1M)-1MTΣ-1ζN(θ0) (12) is asymptotically χ2 distributed, with the number of degrees of freedom equal to the dimension of θ. It is centered under H0 and noncentrality under H1. With a preassigned false alarm rate R and the dimension of the normalized residual, a threshold χR2 can be obtained
where
a ) xN(θˆ 0 - θ0)
(18)
(The proof can be found in appendix 1.) Comparing eqs 16 and 17 with eqs 10 and 11, one can conclude that the estimated error (θˆ 0 - θ0) does not affect the covariance matrix of the normalized residual. Rather, it only contributes a constant bias Ma to the mean of the normalized residual under both H0 and H1. However, the real value of the estimation error θ0 - θˆ 0
1872
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003
is unknown; hence, it is necessary to use the original scheme to perform a χ2 test on ζN(θˆ 0) as follows:
Because ht is a noncentral, χ2-distributed, random variable in a fault-free case, the probability of ht exceeding χR2 in such a case (i.e., the false alarm rate) can be larger than the preassigned value R. Therefore, the threshold χR2 needs to be revised as χˆ R2 so that
eters are caused by the change in room temperature during a day. Because these parameter changes are acceptable, we need to revise the threshold to reduce false alarms due to regular fluctuations of system parameters. Supposing that A is the regular range of the parameter fluctuations, the algorithm to revise the threshold is as follows: At one moment, suppose that the real parameter value is θ and it is in the range of regular fluctuation. The normalized residual is
P(th>χˆ R2|H0) e R
ζN(θ0) ∼ N(Ma,Σ)
(26)
a ) xN(θ - θ0)
(27)
ht ) ζNT(θˆ 0)Σ-1ζN(θˆ 0)
(19)
3.2. Revision of the Threshold. Our objective is to find a new threshold χˆ R2 so that the false alarm rate P(t>χˆ R2|H0) is smaller than R. If, for a threshold χˆ R2, t e χR2, then the following inequality is always true: 2
ht e χˆ R
Hence, P(teχˆ R2|H0) > P(theχR2|H0) and - P(teχˆ R2|H0) e1 - P(theχR2|H0) ) R.
The χ2 value of ζN(θ0) is
(20) 2|H
P(t>χˆ R
0)
k χ 2 + kaTMTΣ-1Ma k-1 R
t ) ζN(θ0) Σ-1ζN(θ0)
)1
Therefore, the problem is transferred to finding a threshold χˆ R2 such that ht e χˆ R2 is always true when t e χR2. Besides, χˆ R2 should be set as small as possible to maintain quality detection performance. Theorem 2. For every k > 1, the following threshold χ` R2 guarantees that when t e χR2, we always have ht e χˆ R2:
χˆ R2 )
where
(21)
Because the mean of ζN(θ0) is not zero, then t is a noncentral χ2-distributed random variable. If the original threshold is used, this will lead to increased false alarms. To reduce the false alarm rate to the preassigned value R, a new threshold should be found such that the false alarm rate
where 2
∆
T
T -1
χβ ) a M Σ Ma (The proof can be seen in appendix 2.) In eq 26, χβ2 remains to be determined. Usually the value of nominal parameter θˆ 0 is obtained through system identification and it is assumed that
θˆ 0 ∼ N(θ0,Σ1)
(23)
a ) (θˆ 0 - θ0)xN ∼ N(0,NΣ1)
(24)
Then,
and χβ2 can be calculated as 2
T
T -1
χβ ) E(a M Σ Ma)
(25)
When there is no model uncertainty, θˆ 0 ) θ0.
a ) (θˆ 0 - θ0)xN ) 0
P(t>χˆ R2|θ∈A) e R
(29)
ht ) (ζN(θ0) - Ma)TΣ-1(ζN(θ0) - Ma)
(30)
Let
and the minumum threshold is
χˆ R,min2 ) (xχR2 + xχβ2)2 ) χR2 + χβ2 + 2xχR2χβ2 (22)
(28)
Because
ζN(θ0) - Ma ∼ N(0,Σ) ht is a central χ2 random variable with the degree of freedom equal to the dimension of the normalized residual. Therefore, P(th>χR|θ∈A) ) R. If a new threshold can be found that when t e χR is always true when ht < χR, then
P(t P(thχˆ R2|θ∈A) ) 1 - P(teχˆ R2|θ∈A) e 1 P(theχR2|θ∈A) ) R Using the result in section 3.3, the revised threshold can be calculated as
χˆ R2 ) (xχR2 + xχβ2)2
(31)
where χβ2 should be computed as
χβ2 ) max(N(θ - θ0)TMTΣ-1M(θ - θ0)), θ ∈ A (32)
and χβ2 ) 0.
χˆ R,min2 ) (xχR2 + xχβ2)2 ) χR2 3.3. Extension to Parameter Fluctuation. Sometimes the system undergoes acceptable parameter fluctuations. For example, changes of the system param-
4. Numerical Example In the following simulation studies, a first-order single input-single output system is used:
y(k) + cy(k - 1) ) bu(k - 1) + e(k)
(33)
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003 1873
Figure 1. χ2 estimation of the normalized residual in a fault-free case.
Figure 2. Faulty case due to a change in the parameter c.
where u(k) and y(k) are the input and output of the system at time k, respectively, and e(k) is a white noise. The nominal values of the system parameters are defined as c0 ) -0.6 and b0 ) 0.4. The primary residual is calculated as
H(u(k),y(k)) )
[
]
-y(k - 1) (y(k) + cy(k - 1) u(k - 1) bu(k - 1)) (34)
The system parameters c and b were estimated using the least-squares method as
cˆ 0 ) -0.5899, bˆ 0 ) 0.4691 The sensitivity matrix M and the covariance matrix Σ were estimated using training data. The false alarm ratio R was set at 0.01. The original threshold was obtained from the standard χ2 table as χR2 ) 9.21. With the covariance matrix of the estimated parameters cˆ 0 and bˆ 0, the sensitivity matrix M, and the covariance matrix Σ, the new threshold could be calculated as χˆ R2 ) 13.57. To show the improvement of the new threshold, it was tested under both fault-free and faulty cases.
4.1. Results in a Fault-Free Case. The process model (33) was excited with a white noise signal with a variance of 1 for 10 000 time units. A χ2 value was calculated at every 10 time units in order to reduce the computation load. Every χ2 value was based on the past 100 time units (i.e., N ) 100). A plot of χ2 value versus time units is shown in Figure 1. It indicates a false alarm rate of 2.3% at the original threshold and 0.7% at the revised threshold. Indeed, the revised threshold was more suitable because it reduced the false alarm rate to approximately the preset value of 1%. 4.2. Results in Faulty Cases. The following simulations were aimed at testing the effect of the revised threshold when faults occur. The process model was excited with a white noise signal with a variance of 1 for 5000 time units. A χ2 value was calculated at every 10 time units. Every χ2 value was based on the past 100 time units. In Figure 2, the parameter c was changed from -0.6 to -0.8 at the 2500th time unit. In Figure 3, the parameter b was changed from 0.4 to 1.3 at the 2500th time unit. The plots of χ2 value versus time units are shown in Figures 2 and 3. They show that when faults occurred, the effect of the revised threshold on the fault detection is not significant. The local approach
1874
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003
Figure 3. Faulty case due to a change in the parameter b.
Figure 4. Scheme of a typical Fourdrinier paper machine.
could detect changes in the system parameters c and b effectively. The simulation results proved the effectiveness of the revised threshold because it reduced the false alarm rate to below the expected value while maintaining the detection performance. 5. Application of Fault Detection in Papermaking Systems To test the applicability of the local approach and the revised threshold, they were applied to an industrial papermaking CD control system. Figure 4 shows a typical paper machine. In making paper, a mixture of fibers, water, and chemicals is pumped into the headbox, which is designed to evenly distribute the mixture onto a moving wire. After the press and dryer sections, the paper is calendared and packaged.11 In papermaking systems, there are three important properties: basis weight, moisture, and caliper. To
obtain a high-quality paper product and a high product rate, automatic control systems have been designed to control each of the above three properties. With the development of measurement devices and actuators, CD control has been rapidly gaining prominence since the 1980s.12 The adjustment of CD properties is done through a series of actuators mounted across the machine. The performance of a CD control system depends not only on the control algorithm but also on the reliability of the actuators. With the growing demand for fault tolerance, FDI has become an important issue for CD control algorithms. It is helpful to monitor the system online using system input and output data and to detect the actuator faults. By detecting the actuator faults and reporting them to the supervision system, one can repair the faulty actuators, improve the availability of the system, and maintain the overall control performance. In the following simulation,
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003 1875
the local approach will be applied to the problem of actuator fault detection for CD control systems of papermaking systems. 5.1. System Model. An open-loop CD system can be described by the ARX model below. To simplify the problem, it is assumed that the dynamics for all outputs are the same. Under the above assumption, the dynamic matrix can be written as follows:
y(k) + cy(k - 1) ) Bu(k - d) + e(k)
(35)
where (i) u(k) and y(k) are vectors of the system inputs and outputs at time k and e(k) is a vector of white noise and (ii) B is the interaction matrix, which contains the coupling information, d is the time delay, and c is the dynamic coefficient of the system. It is assumed that the faults only cause changes in interaction matrix B. Therefore, the actuator faults can be detected by detecting changes of elements in the interaction matrix B. Because of its high dimension, it is very difficult to detect the change of one single element in the interaction matrix B. Because the fault with an actuator will cause changes in all of the elements in the column corresponding to this actuator, it is assumed that the fault affects the elements at the same rate. This means that when a fault happens with the pth actuator, the interaction matrix B becomes
[
]
b11 b12 ‚‚‚ rb1p ‚‚‚ b1n b b ‚‚‚ rb2p ‚‚‚ b2n ) Bf B h ) 21 22 l l l l l l bm1 bm2 ‚‚‚ rbmp ‚‚‚ bmn where
yi(k) + cyi(k - 1) ) Biui(k - d)fi + ei(k), for i ) 1, ..., l (37) For each of those systems, one can detect the faults of those actuators that have an effect on the outputs in the particular section. Primary residuals can be calculated as
Hi(k) ) (Biui(k - d))T[yi(k) + cyi(k - 1) Biu(k - d)fi0], for i ) 1, ..., l (38) 5.3. Simulation Results. The best way to test the applicability of the proposed fault detection algorithm is to intentionally disable one or several actuators in a CD control system and see whether the algorithms can detect the faults or not. However, it is not feasible to do this because it may deteriorate the control performance. An alternative way is to collect data from a CD control system and simulate actuator faults by modifying the input and output data. If the proposed algorithms can successfully detect the simulated faults, one can anticipate that it will also detect real actuator faults in an actual CD control system. The data used in the simulation studies were collected from the moisture profile CD control system with a steam box application in a Canadian paper mill producing newsprint grade paper. Steam boxes are quite slow systems because the dynamics of the systems depend on the thermal momentum of the paper machine. In the moisture CD control system, there are 50 actuators and 1500 measurements for each scan. The process model was obtained through a bump test, and the parameters were determined through the Levenburg-Marquhardt algorithm.11 The model of the CD control can be described as
f ) [1 ‚‚‚ r ‚‚‚ 1 ]T
∆y(k) ) c∆y(k - 1) + B∆u(k - d)
Therefore, the system model (34) can be rewritten as
where y(k) is the CD profile, u(k) is the vector of the setpoints for the actuators, and ∆ ) 1 - q-1. In this system, c ) 0.9432 and d ) 2. The move of one actuator will cause the change of approximately 10 CD measurements around the actuator. In the following simulation study, closed-loop data obtained from an industrial CD control system were used to assess the performance of the proposed algorithms. An artificial fault was introduced by modifying the output data. Suppose that u(k) is the vector of setpoints for the actuators and uˆ (k) is the actual position of the actuators. In the fault-free case, u(k) should be equal to u ˆ (k); however, actuator faults should cause a difference between u(k) and u ˆ (k) and further lead to the change of system output. Let
y(k) + cy(k - 1) ) Bu(k - d)f + e(k)
(36)
where f is an n × 1 vector. Under nominal conditions, all elements of f are 1, i.e., f 0 ) [1, ..., 1]T. f can be called the index of fault because changes of its elements directly indicate faults of the corresponding actuators. The problem of detecting actuator faults is then transferred to the problem of detecting changes of f with the system input u and output y. 5.2. Actuator Fault Detection. With model uncertainty and noise present in the system model, the change of one actuator may not have enough effect on the system outputs. This may lead to a fault not being detected by the algorithm. In fact, the interaction matrix is a sparse matrix. In other words, based on the sparse matrix, any individual output is only affected by the move of several actuators and the impact of one actuator will only appear at several system outputs. Therefore, one can divide the system actuators into several sections and use all of the outputs that are affected by those actuators to detect the faults. Because the dimension of the output within a section is much less than that of the total output, it is anticipated that, by appropriate division of the system output into sections, the detection algorithms can be more sensitive to the fault of a singleactuator. Suppose the system actuators are divided into l sections and each of the subsystems can be rewritten as
ˆ (k) ue(k) ) u(k) - u
(39)
(40)
With the system model (39), the changes of the system output due to actuator faults can be calculated as
ye(k) ) cye(k-1) + Bue(k - d)
(41)
Then the actual measurements should be
yˆ (k) ) ye(k) + y(k)
(42)
In practice, only the faulty output yˆ (k) and the control signal u(k) can be obtained. These will be used to detect the simulated actuator faults.
1876
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003
proposed to recalculate the threshold value to reduce the false alarm rate to the preset value. A similar algorithm was also proposed to calculate the threshold to accommodate the regular parameter fluctuations. Simulation results show that the revised threshold can reduce the false alarm rate to a value close to the preassigned level and maintain the ability of detecting faults in the case of both parameter uncertainties and parameter fluctuations. The local approach can be used to detect actuator faults in papermaking CD control systems, and the proposed algorithm can be used to reduce false alarms due to noise and model uncertainties. Appendix 1: Proof of Theorem 1 The primary residual calculated with the estimated parameter θˆ 0 is Figure 5. Fault-free case for the 25th actuator.
In the following simulation study, a set of raw data of 400 scans was obtained from a CD control system. For each of the sections, the sensitivity matrix and covariance matrix were estimated with the extracted CD component and input data. A threshold was also calculated using the algorithms proposed in section 3. Because of the limited number of data, the moving window of eight scans was used. A χ2 value was calculated for each scan. Figures 5 and 6 show the results of χ2 tests in a fault-free case and a faulty case, respectively. In the fault-free case, the revised threshold greatly reduced false alarms. In the faulty case, the local approach successfully detected the abrupt fault of a single actuator with the presence of model uncertainty and disturbance. 6. Conclusions The robustness of the local approach with respect to model uncertainties was investigated. An algorithm was
Figure 6. Detection of an abrupt fault at the 25th actuator.
H(θˆ 0,yk,φk) ) φk(yk - φkTθˆ 0) ) φk(yk - φkTθ0) + φkφkT(θˆ 0 - θ0) ) H(θˆ 0,yk,φk) + φkφkT(θˆ 0 - θ0)
(A1)
Therefore,
ζN(θˆ 0) )
1
N
∑ H(θˆ 0,yk,φk) )
xNk)1
1
N
∑ H(θ0,yk,φk) +
xNk)1
1
) ζN(θ0) +
[∑ ] 1
N
Nk)1
xN
N
φkφkT(θˆ 0 - θ0) ∑ k)1
φkφkT (θˆ 0 - θ0)xN
When N approaches infinity
(A2)
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003 1877
1
N
∑ φkφkT f cov(φk) Nk)1
(A3)
The primary residual of linear system (1) is given by eq 4; then
∂ H(θ,yk,φk) ) Eθ0(φkφkT) M(θ0) ) Eθ0 ∂θ
(A4)
(ζN(θˆ 0) - Ma)TΣ-1(ζN(θˆ 0) - Ma) < χR2 S 1 ζNT(θˆ 0) Σ-1ζN(θˆ 0) < χR2 + ζNT(θˆ 0) Σ-1ζN(θˆ 0) + k kaTMTΣ-1Ma - aTMTΣ-1Ma
(
S 1-
1 T ζ (θˆ ) Σ-1ζN(θˆ 0) < χR2 + k N 0 (k - 1)aTMTΣ-1Ma
)
Therefore, when k > 1
Therefore:
1
N
∑ φkφkT f M(θ0)
(A5)
Nk)1
(ζN(θˆ 0) - Ma)TΣ-1(ζN(θˆ 0) - Ma) < χR2 S k ζNT(θˆ 0) Σ-1ζN(θˆ 0) < χ 2 + kaTMTΣ-1Ma (A14) k-1 R Therefore, if the threshold χˆ R2 is set as
because
ζN(θ0) ∼ N(0,Σ)
under H0
ζN(θ0) ∼ N(-Mλ,Σ)
(A6)
under H1
(A7)
With the notation
M ) M(θ0), a ) (θˆ 0 - θ0)xN
χˆ R2 )
k
ζN(θˆ 0) ∼ N(Ma,Σ)
under H0
ζN(θˆ 0) ∼ N(Ma - Mλ,Σ)
(A8)
under H1
(A9)
Appendix 2: Proof of Theorem 2 We are looking for a threshold χˆ R2 so that when t < χR we always have ht < χˆ R2. 2,
t ) ζN(θ0)TΣ-1ζN(θ0) ) (ζN(θˆ 0) - Ma)TΣ-1(ζN(θˆ 0) Ma) (A10) as
(ζN(θˆ 0) - Ma)TΣ-1(ζN(θˆ 0) - Ma) < χR2 S T -1
a M Σ Ma (A11) where S means equivalent to. It is easy to prove that for any n-dimensional vectors X and Y
1 2XTY e XTX + kYTY k
(A12)
where k is any positive real number. Therefore, for positive number k
1 2ζN (θˆ 0) Ma < ζNT(θˆ 0) Σ-1ζN(θˆ 0) + kaTMTΣ-1Ma k (A13) T
Substituting eq A13 into eq A11,
(k -k 1χ
2
R
)
+ kaTMTΣ-1Ma ) (xχR + xχβ2)2 (A16) 2
Nomenclature a ) difference between the estimated parameter and the nominal parameter B ) interaction matrix of the paper machine CD control system b, c ) parameters of a linear system E ) expectation e(k) ) white noise at time k f ) index vector of faults in CD control systems H(θ0,yk,φk) ) primary residual at time k k ) time M(θ0) ) sensitivity matrix of the normalized residual N ) number of samples t ) generalized likelihood ratio u(k) ) system input at time k ue(k) ) change of the system input due to actuator faults y(k), yk ) system output at time k Greek Letters
ζNT(θˆ 0)Σ-1ζN(θˆ 0) < χR2 + 2ζNT(θˆ 0) Σ-1Ma T
(A15)
ζNT(θˆ 0) Σ-1ζN(θˆ 0) e χˆ R2 will always be true when (ζN(θˆ 0) - Ma)TΣ-1(ζN(θˆ 0) - Ma) < χR2 is true. Optimizing χˆ R2 with respect to k gives
χˆ R,min2 ) min
the distribution of the normalized residual calculated ζN(θˆ 0) with the estimated parameter θˆ 0 can then be summarized as
k χ 2 + kaTMTΣ-1Ma k-1 R
θ ) system parameter vector θ0 ) nominal value of the system parameter θˆ 0 ) estimated value of the system parameter ζN(θ0) ) normalized residual λ ) change of the system open-loop vector Σ ) covariance matrix of the normalized residual χR2 ) threshold of the generalized likelihood test χˆ R2 ) revised threshold χβ2 ) threshold due to model uncertainty R ) false alarm ratio of the generalized likelihood test
Literature Cited (1) Frank, P. M. Fault Diagnosis in Dynamic Systems Using Analytical and Knowledge-based RedundancysA Survey and Some New Results. Automatic 1990, 26 (3), 459-474. (2) Gertler, J. J. Survey of model-based failure detection and isolation in complex plants. IEEE Control Syst. Mag. 1988, 3-11. (3) Basseville, M. Detecting changes in signals and systemss a survey. Automatica 1988, 24 (3), 309-326.
1878
Ind. Eng. Chem. Res., Vol. 42, No. 9, 2003
(4) Isermann, R. Fault diagnosis of machines via parameter estimation and knowledge processingstutorial paper. Automatic 1993, 29 (4), 815-835. (5) Benveniste, A.; Basseville, M.; Moustakides, G. V. The asymptotic local approach to change detection and model validation. IEEE Trans. Autom. Control 1987, 32 (7), 583-592. (6) Basseville, M. On-board component fault detection and isolation using the statistical local approach. Automatica 1998, 34, 1391-1415. (7) Zhang, Q.; Basseville, M.; Benveniste, A. Early warning of slight changes in systems. Automatica 1994, 30, 95-113. (8) Huang, B. Process and control loop performance monitoring through detection of abrupt parameter changes. Proceedings of the 1999 IEEE Canadian Conference on Electrical and Computer Engineering, Edmonton, Alberta, Canada, 1999. (9) Zhang, Q.; Basseville, M. Monitoring nonlinear dynamical systems: a combined observer-based and local approach. Proceed-
ings of the 37th IEEE Conference on Decision & Control, Tampa, FL, 1998. (10) Zhang, Q.; Basseville, M.; Benveniste, A. Fault detection and isolation in nonlinear dynamic system: a combined inputoutput and local approach. Automatic 1998, 34 (11), 1359-1373. (11) Ghofraniha, J. Cross directional response modeling, identification and control of dry weight profile on paper machines. Ph.D. Thesis, Department of Electrical Engineering, University of British Columbia, Vancouver, British Columbia, Canada, 1998. (12) Kwok, E.; Li, P. MD-CD Interaction Modeling and Control of Paper Machines. Can. J. Chem. Eng. 1997, 75 (1), 143-151.
Resubmitted for review October 28, 2002 Revised manuscript received December 23, 2002 Accepted February 5, 2003 IE010797I
