1 Performance trade-offs for networked jump observer-based fault diagnosis Daniel Dolz, Ignacio Pe˜narrocha, Roberto Sanchis Abstract—In this paper, we address the fault diagnosis problem for discrete-time multi-sensor systems over communication networks with measurement dropouts. We use the measurement outcomes to model the measurement reception scenarios. Based on this, we propose the use of a jump observer to diagnose multiple faults. We model the faults as slow time-varying signals and introduce this dynamic in the observer to estimate the faults and to generate a residual. The fault detection is assured by comparing the residual signal with a prescribed threshold. We design the jump observer, the residual and the threshold to attain disturbance attenuation, fault tracking and detection conditions and a given false alarm rate. The false alarm rate is upper bounded by means of Markov’s inequality. We explore the tradeoffs between the minimum detectable faults, the false alarm rate and the response time to faults of the fault diagnoser. By imposing the disturbances and measurement noises to be Gaussian, we tighten the false alarm rate bound which improves the time needed to detect a fault. A numerical example is provided to illustrate the effectiveness of the theory developed in the paper. Index Terms—Fault diagnosis, false alarm rate, time to detect faults, jump linear system, dropouts. I. INTRODUCTION Networked control systems have been extended to many industrial applications due to the diverse offered advantages, as the reduction on the installation cost or the increase on the flexibility, provided by the communication network [1]. In these kinds of systems, the controller unit, the sensors and the actuator are not collocated and the exchange of information is done through a shared network, leading to some networkinduced issues as time delays and dropouts [2], [3]. Owing to the need for reliability, safety and efficient operation of these networked systems, model-based fault diagnosis methods [4] have been recently introduced to operate over networks [5]. Fault detection over communication networks when using an observer-based fault detection scheme is addressed by the comparison between a residual signal generated with the estimated system outputs and a threshold. The residual is conceived to balance the robustness against network effects and disturbances, and the fault sensitivity [6]–[9]. Assuring a predefined false alarm rate (FAR) is a key problem. In the majority of the networked fault detection Copyright (c) 2015 IEEE. Personal use of this material is permitted. However, permission to use this material for any other purposes must be obtained from the IEEE by sending a request to
[email protected]. Daniel Dolz, Ignacio Pe˜narrocha, Roberto Sanchis are with Department of Industrial System Engineering and Design, University Jaume I of Castell´on , Spain {ddolz,ipenarro,rsanchis}@uji.es. This work has been funded by MICINN project DPI2011-27845-C0202 from the Spanish government and by project P1·1B2013-51 and grant PREDOC/2011/37 from Universitat Jaume I. proposals, the threshold is chosen to reduce the FAR to the minimum [10], [11], but without quantifying it. Some works as [7], [8], [12] characterize the mean and variance of the residual and use Markov’s inequality to impose a desired FAR bound. However, Markov’s inequality is known to be conservative [13]. The main problem to get a proper FAR bound is to obtain the probability distribution of the residual signal. In [14] the residual was computed as a quadratic form of the outputs estimation error by means of the inverse of the outputs estimation error matrix covariance given by a Kalman filter. With that, their residual signal follows a chi-squared distribution and an exact FAR can be fixed. But, to the best of the authors’ knowledge, the extension to observers with predefined gains (which have less implementation cost) for networked systems with dropouts has not been addressed. Regarding the fault estimation problem, the most common approach is to make the residual track the fault or a weighted fault signal by guaranteeing some performances of the fault estimation error under disturbances and the network issues [15]– [19]. Recently, to improve the fault estimation performances, the authors in [20] introduced a dynamic of the fault signal on the fault estimator. Fault detection and estimation can be combined to attain fault diagnosis. According to [21], the performance of a fault detection algorithm is defined by means of the trade-offs between the time to detect a fault and the FAR. This definition can be extended to the fault diagnosis case by considering also the convergence speed of a norm of the fault estimation error. The authors in [22] show that there exists a trade-off between the fault detection rate and the FAR. More recently, the existence of a compromise between the time to detect a fault and the fault sensitivity has been demonstrated in [23]. Nevertheless, none of them explores the compromises between the minimum detectable faults, the FAR and the fault diagnosis (detection and estimation) speed. The dropouts in the fault diagnosis problem over communication networks have been mainly studied in the packetized case [7], [8], [15]. The multi-sensor case was studied in [24] with an invariant observer gain approach, however the use of jump observers that adapt their gains to the network scenario has been proved to enhance the estimation performances [25], [26]. Networked jump observer-based fault estimators have recently started to receive attention [6], [8]. Motivated by the previous analysis, in this paper we face the fault diagnosis problem for multi-sensor systems with dropouts through the combination of fault detection and fault estimation. The faults are characterized as slow-time varying signals and the network dropouts are modeled with the combination of
2 available measurements at the fault diagnoser. We introduce a jump observer to estimate the faults and define the residual signal as a quadratic form of the estimated fault vector. The design of the jump observer and residual is addressed through an iterative linear matrix inequalities (LMIs) procedure that allows obtaining the predefined set of observer gains and the fault detector parameters. The design is carried out to achieve disturbance and measurement noise attenuation, and fault diagnosis performances under a prescribed FAR. We propose two design strategies: the first one consists of fixing the response speed to faults and minimizing the minimum detectable fault, and the second one consists of fixing the minimum detectable fault and minimizing the response time. The trade-offs between the minimum detectable faults, the FAR and the delay between fault occurrence and detection (response time of the fault estimator) are highlighted. Furthermore, we derive two ways of bounding the FAR depending on whether the residual signal probability distribution is unknown (Markov’s inequality approach) or known as a result of assuming Gaussian disturbances and measurement noises (chi-squared approach). Notation : Let Aand Bbe some matrices. A(i, i)defines the i-th diagonal element of A. The maximum and minimum eigenvalues of Aare denoted by λ(A)and λ(A)respectively. ABmeans that matrix A−Bis negative semidefinite. Similar applies to . The direct sum is represented by L, where ALBis a block diagonal matrix with Aand Bon its diagonal. Operator vec(A)generates a vector by stacking the columns of matrix A. Let x[t]∈Rnbe a stochastic process. Expected value and probability are denoted as E{·} and Pr{·}. We write kx[t]k2 2,x[t]Tx[t]for the ℓ2norm of x[t],kxk∞,maxtmaxi|xi[t]|for the ℓ∞norm of xand kxk2 RMS ,limK→∞ PK−1 t=0 1 Kkx[t]k2 2for its RMS norm. II. PROBLEM FORMULATION Let us consider linear time invariant discrete-time systems defined by equations x[t+ 1] = A x[t] + Buu[t] + Bww[t] + Bff[t],(1) where x∈Rnis the state, u∈Rnuis the vector of known inputs, w∈Rnwis the state disturbance assumed as a random signal, uncorrelated in time, with zero mean and known covariance matrix E{w[t]Tw[t]}=Wfor all t, and f∈Rnfis the fault vector. Throughout this work we assume that the known input uis causally available at all times, see Fig. 1. This general model includes as a particular case a system without known inputs, by simply taking Bu= 0. The measurable outputs of the system are modeled by equation y[t] = C x[t](2) where y∈Rnyis the output vector. Different sensors with different characteristics on sampling rate or noise, that may have faults, can be connected to one single measurable output, but at least each measurable output is measured by one sensor, having nm≥nysensors. We define the measurement value as mj[t] = cjx[t] + hjf[t] + vj[t], j = 1,...,nm(3) fs2[t] y1[t] transmission outcome ma 2[t], α2[t] ma 1[t], α1[t] dropouts Network Plant y2[t] me 1[t] me 2[t] r[t] u[t] Fault Diagnoser Unit ˆ f[t] • Sensor 1 Sensor 2 fp[t] fs1[t] Fault Estimator Fault Detector • Residual generator fault? Fig. 1. Networked fault diagnosis problem for two sensors with possible faults in the plant (fpfor actuators and other faulty components) and in the sensors (fs1, fs2). where mj[t]∈Rrepresents the t-th measurement of the j-th sensor and vj[t]∈Rthe j-th sensor noise assumed as a zero mean random signal with known variance E{vj[t]2}=σ2 jfor all t, that is uncorrelated with respect to the time index t. We also consider that viis mutually uncorrelated with vj6=i.cj denotes one row of matrix C(several cjcould be equal and correspond to the same row of C) and hjdenotes each one of the rows of matrix H. In the current work, we model the fault signal as a slow time-varying one (cf. [20], [27]), i.e., f[t+ 1] = f[t] + ∆f[t](4) where ∆f[t]is the variation of the fault from instant tto t+1. Equation (4) allows modeling, for instance, step signals (∆f[t] only takes a nonzero value at the time the fault appears) or ramp signals (∆f[t]takes a constant value), that have been widely used in the literature to analyze the behavior of fault detection algorithms [4]. Along this paper, we consider that w[t],vj[t]for all j= 1,...,nmand ∆f[t]are mutually uncorrelated for all t. We introduce an extended order model to include the fault dynamic as z[t+ 1]= ¯ Az[t]+ ¯ Buu[t]+ ¯ Bww[t] + ¯ Bf∆f[t](5) with z[t] = x[t]Tf[t]TTand ¯ A=A Bf 0I,¯ Bu=Bu 0,¯ Bw=Bw 0,¯ Bf=0 I where z∈R¯nwith ¯n=n+nf. In this work we intend to detect and estimate (diagnose) the possible system faults (represented by vector f[t]) when the measurements are transmitted through a communication network that may induce dropouts. In this case, the system output measurements are not available at every discrete time instant. When the dropout rate is high, the fault estimation problem becomes more difficult and the importance of a fast response to faults and a low FAR becomes more evident.
3 The measured and transmitted value from sensor jat instant tis me j[t] = ¯cjz[t] + vj[t],(6) with ¯cj= [cjhj]and j= 1,...,nm. We assume that the pair (¯ A, ¯ C)is detectable (being ¯ Cthe matrix whose rows are ¯cj). Remark 1. If the pair (¯ A, ¯ C)is not detectable (i.e., nf> nm), only a combination of the faults can be detected. Then, a previous transformation of the system, as proposed in [28], must be done (leading to new ¯nffaults such that ¯nf≤nm) before the proposed technique becomes applicable. A. Network transmissions characterization Each sensor samples its output synchronously with the known input update and sends independently a time-tagged packet with the measurement me j[t]to the fault diagnoser station, through an unreliable communication network with packet dropouts. We define the binary variable αj[t]that indicates the availability of the j-th sensor measurement (j= 1,...,nm) at each instant t, as. αj[t] = 0if me j[t]is not received at t, 1if me j[t]is received at t. (7) Then, the availability matrix α[t] = Lnm j=1 αj[t]is a binary diagonal matrix that can only have ones in its diagonal. Thus, using α[t]we can redefine the available measurements at instant tas ma[t] = α[t]¯ Cz[t] + v[t].(8) Note that a component of vector ma[t]is null when the corresponding measurement is not available. v[t] = [v1[t]···vnm[t]]Tis the measurement noise vector with covariance E{v[t]v[t]T}=V=Lnm j=1 σ2 j(for all t). The possible values of α[t]at each instant tbelong to a known finite set α[t]∈Ξ = {η0, η1,...,ηq}, q = 2nm−1,(9) where ηidenotes each possible combination of the available measurements at the fault diagnoser station (measurement reception scenario). Matrix η0denotes the scenario in which there is no measurement available, i.e., η0= 0. We characterize the network behavior using the total probability of each scenario in Ξ. We denote by pi=Pr{α[t] = ηi} the probability of having the measurement reception scenario ηiat instant t.p0denotes the probability of having no measurements. In the current paper, we assume that the arrival probability from each sensor is governed by an independent and identically distributed process [29]. We denote by βjthe probability of having available the measurement from sensor jat instant t, i.e., βj=Pr{αj[t] = 1}. Then, the probability of having a given combination of available measurements ηi∈Ξis pi=Pr{α[t] = ηi}=Y j∈I(ηi) (1 −βj)Y j6∈I(ηi) βj(10) for all i= 0,...,q where I(ηi),{j|ηi(j, j) = 0}. B. Fault diagnosis method We propose the following fault estimation algorithm for system (5)-(6). At each instant t, the model is run in open loop leading to ˆz[t−] = ¯ Aˆz[t−1] + ¯ Buu[t−1].(11) If no measurement is received, we keep the open loop estimation, i.e., ˆz[t] = ˆz[t−]. If a measurement arrives at instant t=tk, the state is updated as ˆz[tk] = ˆz[t− k] + L[tk] (ma[tk]−α[tk]¯ Cˆz[t− k]),(12) where L[tk]is the updating gain matrix and ma[tk]is defined in (8). Remark 2. While t∈Nrefers to each time instant, tk (with k∈N) enumerates only the instants where some measurements are received. For instance, if we receive some measurements only at instants tk= 8 and tk+1 = 10, but not at t= 9, then instant tk+ 1 = 9 (or tk+1 −1 = 9) refers to instant 9, when no measurement is received. Let us denote z[tk]by zk. Defining the extended state estimation error at updating instants as ˜zk=zk−ˆzk, the estimation error dynamics is given by ˜zk=(I−Lkαk¯ C)¯ ANk˜zk−1−Lkαkvk + Nk X l=1 (I−Lkαk¯ C)¯ Al−1BWW[tk−1+l−1] (13) being BW= [ ¯ Bw¯ Bf]and W[tk−1+l−1] = [w[tk−1+l− 1]T∆f[tk−1+l−1]T]T.Nkdenotes the number of consecutive instants without measurements (which is unbounded), i.e., Nk=tk−tk−1. The fault detection algorithm uses the estimated faults to compute a residual signal at instants t=tkas rk=ˆ fT kF−1ˆ fk,(14) where the common fault detection decision is given by if rk≤rth No fault if rk> rth Fault being rth >0a threshold to be defined. Then, fault isolation is achieved by means of the combination of fault detection and fault estimation, allowing us to identify which is the origin of the fault. Remark 3. According to [4], the minimum detectable fault is a fault that drives the residual to its threshold, provided no other faults, disturbances and measurement noises are present. Then, assuming a zero fault estimation error (i.e. ˆ f=f), each diagonal element of Fin (14) multiplied by rth defines the minimum detectable fault as fmin,l =rthF(l, l)for the corresponding channel (l= 1,...,nf). Considering the fault detection logic, the FAR is defined as the average probability of rising false alarms over an infinitetime window, i.e. Ψ = lim K→∞ K−1 X k=0 Pr{rk> rth |fk= 0}.(15)
4 The aim of this work is to compute the gain matrices Lk, the matrix F, and the threshold rth such that the fault diagnoser attains disturbance and measurement noise attenuation, and fault diagnosis performances for a given FAR. These objectives can be reached with an invariant observer gain (as in the majority of reviewed works), or with a jump one (e.g. [6], [8]). In this work, we relate the gain Lkto the sampling scenario αk, as Lk=L(αk), with the law Lk=Liwhen αk=ηi for αk=η1,...,ηq. Then, the matrices are computed off-line leading to the finite set Lk∈ L ={L1,...,Lq}.(16) III. FAULT DIAGNOSER DESIGN:DROPOUT-FREE Let us first consider the case without measurement dropouts, i.e., βj= 1 for all j= 1,...,nm. In this case, α[t]is always the identity, which implies that each instant tis a measurement instant (tk=t) leading to Lk=Land Nk= 1, for all k. The following theorem presents how to design the observer gain L and the matrix Fthat defines the residual (14) based on the H2norm of system (13). Theorem 1. Consider the estimation algorithm (11)-(12) applied to system (1)-(4) with standard sampling. If there exist symmetric matrices P,F,Γw,Γv,Γf, and full matrix X fulfilling P¯ A0 ¯ ATP¯ Bf 0¯ BT fF 0,(17a) P¯ Bµ ¯ BT µΓµ0, µ ={w, v, f}(17b) with ¯ A= (P−X¯ C)¯ A, ¯ Bw= (P−X¯ C)¯ Bw, ¯ Bv=−X, ¯ Bf= (P−X¯ C)¯ Bf then, defining the observer gain matrices as L=P−1X, the following statements hold: i) In the absence of disturbances, faults, and measurement noises the extended state estimation error (13) converges to zero. ii) Under zero initial conditions (i.e., ˜z0= 0), the fault estimation error is bounded by E{k˜ fk2 RMS} ≤ λ(F)tr(¯ Γ) + nfλ(Γf)∆f2 max,(18) where ¯ Γ = ΓwW+ ΓvVand k∆fk∞≤∆fmax. Proof. See Appendix A. The above theorem states that Fis related to the expected value of the squared RMS norm of the fault estimation error. We can extract from (18) that the fault estimation (and therefore the residual signal) is more sensitive to disturbances and measurement noises when the maximum of the minimum detectable faults (by means of λ(F)) is higher. Furthermore, the lower the value λ(Γf), the lower the effect of the faults on the estimation error. The next theorem extends the results of the previous one to bound the FAR. Theorem 2. For a given threshold rth >0and 0≤φ≤1, and under the premisses of Theorem 1, if tr(ΓwW) + tr(ΓvV) = φ rth,(19) and constraints (17) are fulfilled, then, the following additional statement holds: iii) In the absence of faults and under zero initial conditions, the fault detection algorithm assures a FAR (15) bounded by φ. Proof. See Appendix B. The next theorem extends the previous one showing how the fault estimation error decays at each measurement instant. Theorem 3. For a given threshold rth >0and 0≤φ≤1, and under the premisses of Theorem 2, if Γf−¯ BT fP¯ Bf0,(20) and constraints (17),(19) are fulfilled, then, the following additional statement holds: iv) The fault estimation error given by E{k˜ fkk2 2}decays with ρ= 1 −1 λ(ΓfF).(21) Proof. See Appendix C. The above theorem shows that E{k˜ fkk2 2}decays with ρ, from the initial conditions to the steady state region (see (41)). ρdepends on the maximum eigenvalue of the product ΓfF. If Fis fixed to assure the detection of some given minimum faults, Γfdetermines the response time of the fault estimator (by means of ρ) and therefore the time to detect a fault (as the residual is defined with the estimated faults). Remark 4. Under a step-like fault, the number of instants with measurement reception, denoted by K, until the initial value of the fault estimation error is decreased below ξ(%), characterizes the settling time of the fault estimation vector (time to achieve the 100 −ξof the final value). Kcan be obtained approximately by solving equation ρk+1 =ξ/100, (see (41)) leading to K=log(ξ/100) log(ρ)−1(22) where ⌈·⌉ is the operator that rounds its argument to the nearest integer towards infinity. One of the most used values for ξin system theory is ξ= 2%. Thus, ξ= 2% refers to the number of time instants until reaching the 98% of the fault estimation final value. Remark 5. For a fixed value of F, increasing the FAR by means of φleads to an increase in the values of Γwand Γv, see (19). Higher values on these variables alleviate the constraints over Pin (17b), increasing the solution space in the search for a feasible matrix P. This would allow, for instance, structure constraints over matrix P. Matrix Γfin (20) constrains the last diagonal block on P. Then, increasing φcan enlarge the solution space to find lower values on Γf, which, in turn, lead to lower values of ρ(faster
5 fault diagnosers). These ideas are analyzed in the examples section. We used Markov’s inequality in Theorem 2 to bound the FAR. However, it is well known that the bound yielded by Markov’s inequality may be very conservative (see [13]) because it does not consider the probability distribution of the residual rk. This may result in a real FAR that is some orders of magnitude lower than the desired one, which, as shown in the examples, may lead to a very slow response of the fault diagnoser (characterized by ρin Theorem 3). Most of the works in the literature share this important drawback. In order to overcome this, a more accurate bound on the FAR would be desirable. Assuming that the disturbances wkand the measurement noises vkare Gaussian, we show in the next theorem how to impose an appropriate value to matrix F to force the residual rkfollow a chi-squared distribution, which allows us to tighten the FAR bound. Theorem 4. For the fixed threshold rth =nfand for a given 0≤φ≤1, under the premisses of Theorem 3, if F=φ−1Σf(23) and constraints (17),(19),(20) are fulfilled, with 1 vec(Σf) = I−G¯ A⊗G¯ A−1vec (Y1),(24) Y1=G¯ BwW¯ BT wGT+P−1X V (P−1X)T, G= (I−P−1X¯ C), then, in the absence of faults, under zero initial conditions and Gaussian disturbances and measurement noises, if the fault diagnoser gain is defined as L=P−1X, then the FAR is given by Ψ = 1 −CDFX2 nfrth φ(25) where CDFX2 nf(rth φ) = Pr{rk φ≤rth φ}denotes the cumulative distribution function (CDF) of a chi-squared random variable with nfdegrees of freedom, X2 nf. Proof. See Appendix D. Remark 6. Following the definition of the CDF of a chisquared random variable, the value of φneeded to obtain a desired FAR ψwith the chi-squared approach is always higher (for any value of nf) that the one required with the Markov’s inequality approach. For instance, if nf= 2 and ψ= 10−3using Theorem 2 requires φ= 10−3while Theorem 4 requires φ= 0.145. Following Remark 5, this implies that using the chi-squared approach could lead to fault diagnosers with a faster response to faults than employing the result on Theorem 2. However, Markov’s inequality approach (from Theorem 2) has wider applications because it does not require Gaussian disturbances and noises, as the chi-square approach (from Theorem 4) does. 1The Kronecker product between A∈Rn×mand B∈Rp×qis a block matrix such as A⊗B a11 A· · · a1mB . . ..... . . an1B· · · anmB ∈Rnp×mq . The vectorization of matrix Ais vec(A) = a11 · · · an1a12 · · · anmT. Theorem 4 has shown how to reduce the conservatism of the approach when assuming Gaussian disturbances, but at the cost of including new nonlinear equality constraints that are hard to handle. We will show in the design strategies section how to overcome this issue. IV. APPLICATION TO NETWORKED TRANSMISSION In the previous section we presented how to design the fault diagnoser and to characterize the obtained FAR and response time to faults for measurement transmission without dropouts. In this section we extend the previous results to a more interesting case where measurement information is not always available due to network dropouts. This will stress the need of fast fault detection with a low FAR. The following theorem extends Theorem 3 and shows how to find the set of observer gain matrices (16) and the matrix that defines the residual (14). Theorem 5. For a given threshold rth >0and 0≤φ≤1, consider the estimation algorithm (11)-(12) applied to system (1)-(4). Assume that there can be qdifferent measurement reception scenarios ηi(i= 0,...,q) with a probability pi. If there exist symmetric matrices P,Q,F,Γv,Γw,Γf,, and full matrices Xifulfilling P−M1¯ Bf ¯ BT fF0,(26a) Γw−¯ BT wM2¯ Bw0,(26b) Lq l=1 P M3 MT 3Γv0,(26c) Γf−¯ BT f(M5+M6)¯ Bf0,(26d) Lq l=1 P M4 MT 4Q0,(26e) tr(ΓwW) + tr(ΓvV) = φ rth,(26f) Γf−¯ BT fP¯ Bf0,(26g) with vec(M1) = ϕ(¯ A)−1vec( ¯ ATQ¯ A), M2= (1 −p0)M5+p0 1−p0 M1, M5=1 (1 −p0)2Q, M3= q1 1−p0p1X1η1 . . . q1 1−p0pqXqηq , M4= √p1(P−X1η1¯ C) . . . √pq(P−Xqηq¯ C) , vec(M6) = ϕ(¯ A)−1vec( ¯ ATM5¯ A) + p0 1−p0 vec(M1), and ϕ(¯ A) = I−p0¯ AT⊗¯ AT, then, defining the observer gain matrices as Li=P−1Xi, the following statements hold: i) In the absence of disturbances, faults and measurement noises, (13) converges to zero in average. ii) Under zero initial conditions, the fault estimation error is bounded by E{k˜ fk2 RMS} ≤ λ(F)·φ rth +nfλ(Γf)∆f2 max,(27)
6 where ∆fmax is a constant that depends on the fault magnitude that bounds vector ∆fk∈Rnfas k∆fk∞≤ ∆fmax, being ∆fka vector that fulfills for all kthat ∞ X N=1 NpN−1 0 N−1 X l=0 (⋆)TQ¯ Al¯ Bf∆f[tk+l] | {z } ⋆ = ∆fk T ∞ X N=1 NpN−1 0 N−1 X l=0 ¯ BT f(¯ Al)TQ¯ Al¯ Bf!∆fk. (28) iii) Under zero initial conditions and in the absence of faults, the residual evaluation assures a FAR (15) bounded by φ. iv) The fault estimation error given by E{k˜ fkk2 2}decays with ρ= 1 −1 λ(ΓfF).(29) Proof. See Appendix E. Remark 7. The existence of vector ∆fkdefined in (28) is assured because it represents an equality constrained problem with one equation and nfdegrees of freedom. For instance, under ramp-like faults (∆f[tk+l]is constant), ∆f[tk+l] = ∆fk (for all l= 0,1,...) and ∆fmax =k∆fk∞. Furthermore, the exact value of ∆fmax is not relevant for the analysis. In the aim of reducing the conservativeness introduced by Markov’s inequality to bound the FAR, the next theorem extends Theorem 4 by forcing rkto follow a chi-squared distribution when measurements are subject to dropouts. Theorem 6. If the threshold is set as rth =nfand for a given 0≤φ≤1, under the premisses of Theorem 5, if F=φ−1Σf(30) and constraints (26) are fulfilled for i= 1,...,q, with Σf=¯ BT f(R−p0¯ AR ¯ AT)¯ Bf,vec(R) = Y−1 1vec(Y2), Y1=ϕ(¯ A)−( q X i=1 pi(Gi¯ A)⊗(Gi¯ A)), ϕ(¯ A) = I−p0¯ AT⊗¯ AT, Y2=1 1−p0 q X i=1 piLiηiV ηT iLT i+ q X i=1 piGi(SW)GT i, SW=1 1−p0¯ BwW¯ BT w+p0¯ ASW,∞¯ AT, vec(SW,∞) = ϕ(¯ A)−1vec( ¯ BwW¯ BT w), Li=P−1Xi, Gi=I−Liηi¯ C then, in the absence of faults, under zero initial conditions and Gaussian disturbances and measurement noises, the FAR is given by (25). Proof. See Appendix F. V. FAULT DIAGNOSIS STRATEGIES Based on the derived results on Theorem 5, we propose the following two strategies to address the design of a fault diagnoser depending on the needs of the application. First, let us consider that we desire to detect faults over a certain value, i.e to fix the minimum detectable fault on each channel fmin,l (for l= 1,...,nf), with a guaranteed FAR, and to detect as fast as possible the appearance of faults (i.e., with the lowest ρ). The next optimization problem deals with this design problem. Strategy 1. For a given threshold rth >0, let ψbe the desired FAR, fix φto be φ=ψ, and let Fbe a diagonal matrix such that F=Lnf l=1 f2 min,l/rth. Then, the minimization problem minimize γ subject to X1={(26), F F,ΓfF γI}(31) along variables γ,P,Q,F,Γv,Γw,Γf, and Xi(with i= 1,...,q), leads to the fault diagnoser with the fastest response under faults, able to detect faults over fmin,l (with l= 1,...,nf) with a FAR below ψ. Remark 8. The computational complexity of Strategy 1 can be described as follows. The size of the full involved LMI constraint is (2nm+ 1)(n+nf) + 5nf+nw+nm+ 1. Symmetric matrices as P∈R(n+nf)×(n+nf)have (n+ nf)(1 + n+nf)/2decision variables, while full matrices as Xi∈R(n+nf)×nmhave (n+nf)nm. Furthermore, Strategy 1 is based on semidefinite programming and therefore does not require a high computing capacity. This kind of problems can be solved using MATLAB toolboxes as Yalmip [30] (which can handle large scale problems). Second, let us assume that we desire to impose the response speed under the appearance of faults (by means of ρ) with a guaranteed FAR. Then, the minimum detectable faults can be minimized through the next optimization problem. Strategy 2. For a given threshold rth >0, let ψbe the desired FAR, fix φto be φ=ψ, and let ¯ρbe the given upper bound on how the fault estimation error decays, i.e., ρ≤¯ρ. Then, the minimization problem minimize γ subject to X2=(26),tr(F)≤γ, ΓfF≤(1 −¯ρ)−1I(32) along variables γ,P,Q,F,Γv,Γw,Γfand Xi(with i= 1,...,q), leads to the fault diagnoser with the minimum value of the sum of the squared minimum detectable faults (defined by matrix F) with ρ≤¯ρand a FAR below ψ. Remark 9. Optimization problem (32) is nonlinear because of the bilinear matrix inequality (BMI) that affects the product ΓfF. This can be solved with the following rank constrained problem (1 −¯ρ)−1F F FΛ0,rank ΓfI IΛ≤nf where a new symmetric decision matrix Λhas been added. This problem can be iteratively handled with the well known cone complementarity linearization (CCL) algorithm [31]
7 (which only addresses feasibility by relaxing the rank constraint with a positive semidefinite constraint on the involved matrix) over a bisection algorithm. Solving Strategy 2 is more time consuming than Strategy 1 because of the iterations introduced by the CCL and the bisection algorithm. Nevertheless, it only introduces nf(1 + nf)decision variables (due to Λ) and only increases the full LMI size in 2nf+ 1 over a semidefinite programming problem, and therefore, the computational complexity is not really an issue. Both design strategies are still valid when including nonlinear equality constraints (30) but need more computational effort. The next strategy extends the previous ones to consider the chi-squared approach presented in Theorem 6. Strategy 3. The minimization problem minimize γ subject to Xj,(30), ψ= 1 −CDFX2 nfrth φ(33) along variables γ,P,Q,F,Γv,Γw,Γf, and Xi(with i= 1,...,q) with rth =nf, extends the design made in Strategy 1, if j= 1, or in Strategy 2, if j= 2, to tighten the FAR bound with the chi-squared approach. Remark 10. Optimization problem (33) is nonlinear due to constraint (30). This optimization problem can be solved iteratively with LMI constraints by forcing matrix Fat each step kto be as Fφ−1Σf(Lk−1), until Σf(Lk−1)converges to a constant value, where Σf(Lk−1)is the covariance matrix in (30) evaluated with the observer gains at step k−1. The computational burden of each of the iterations is nearly the same as in Strategy 1 (or Strategy 2), but the total computing time is multiplied at most by the number of iterations. However, we are again dealing with a semidefinite programming problem, therefore the computational load is not a problem. Remark 11. Strategy 3 will lead, in general, to minimum detectable faults under fmin,l (for l= 1,...,nf). If we do not intend to detect faults under fmin,l, we can first solve the optimization problem involved in Strategy 3 and then use rk= ˆ fT kF−1ˆ fkin the real-time implementation (where Fincludes the original prescribed minimum detectable faults, fmin,l). In this case, as we impose in the design that φ−1Σf F, the obtained FAR will be upper-bounded by (43). VI. EXAMPLE Let us consider an industrial continuous-stirred tank reactor process (borrowed from [32]) where the discretized state-space model is A=0.972 −0.001 −0.034 0.863 , Bu=−0.084 0.023 0.076 0.414, Bw=Bu, C =1 0 0 1. We desire to detect faults from the second actuator and the first sensor, i.e. Bf=0.023 0 0.414 0, H =0 1 0 0. The state disturbances and measurement noises are Gaussian with covariance matrices W=0.11 0.03 0.03 0.13, V =0.01 0 0 0.01. We consider that the measurements are independently acquired through a communication network where the probabilities of having available the measurements from each sensor are β= [0.58 0.46]. For ease of analysis, in this example we will only explore the case when we impose that the minimum detectable faults are below some given values and we try to obtain the fastest response to faults of the fault diagnoser, i.e. we will only analyze Strategies 1 and 3. For ease of notation, let us assume that the requirement over the minimum detectable faults is such that FfminI. In the next, we impose the threshold to be rth =nf. First, let us study the compromises between the minimum detectable faults fmin, the desired FAR ψand the speed of the fault diagnoser by means of ρin the design procedure. Fig. 2 illustrates these trade-offs for five different desired FARs with ψ= [10−110−210−310−410−5]and for the two presented approaches to assure them: through Markov’s inequality (left hand side figure, Strategy 1) and through characterizing the probability distribution of the residual signal (right hand side figure, Strategy 3). We note that imposing smaller minimum detectable faults or lower FARs results in a slower response time to faults (ρhigher). We also find that forcing Fto be as defined in (30) (chi-squared approach) results in a faster response under faults (ρsmaller) for the same minimum detectable faults than using Markov’s inequality approach. Furthermore, Fig. 2 shows an asymptotic behavior of ρwith respect to fmin, leading to a minimum achievable value. Second, let us study the behaviour of some fault diagnosers in simulation, where u[t] = 0 for all t. Table I compares the fault diagnosis performances for the case when Fis unconstrained, case C1 (where Markov’s inequality approach is used, Strategy 1) and when Fis constrained to be as in (30), case C2 (where the chi-squared approach is used, Strategy 3). For both cases we impose ψ= 10−3and fmin = 0.6. We also include in Table I a case C3 where we reduced the fmin from case C2 to the half. The matrices Fobtained for the three cases are: FC1 =0.18 0 0 0.18, FC2 =0.161 −0.025 −0.025 0.107 , FC3 =0.022 −0.008 −0.008 0.041 . As illustrated in Table I, for case C3, we can detect smaller faults than in case C2 at the expense of being slower than in case C2. However, we still are much faster than in case C1 where the guaranteed detectable faults were higher. Moreover, as stated in Remark 11 cases C2 and C3 can detect faults below the imposed fmin (fmin,1for the actuator fault and fmin,2for the sensor fault). Concerning the computational burden, obtaining C1 takes 0.4sec (using Yalmip with SeDuMi solver [33] in a i7-3770 processor at 3.40 GHz) while C2
8 ρ= 0.72 fmin ρ, convergence rate of E{k ˜ fkk2 2} Trade-offs with Markov’s inequality approach ρ= 0.75 fmin Trade-offs with chi-squared approach ψ= 10−1 ψ= 10−2 ψ= 10−3 ψ= 10−4 ψ= 10−5 00.4 0.8 1.2 1.6 22.4 2.80 0.4 0.8 1.2 1.6 22.4 2.8 0.7 0.8 0.9 1 Fig. 2. Trade-offs on the observer-based fault diagnoser design. requires 2.5sec (with 10 iterations). Note that as previously stated in Remark 8 and 10, the computational cost is not an issue. After a simulation of 106instants with no faults, we verify that the FAR obtained in simulation (by dividing the number of risen alarms by the total number of simulation time instants) for case C2 and case C3 is the same as forecasted in the design, but for case C1 is much lower (several orders of magnitude) than the imposed bound. This conservativeness of the Markov’s approach results in an extremely slow residual dynamics (as seen in Fig. 3), and a huge time to detect the fault (characterized by 6101 measurement instants, see (22)), that is useless in practice. To alleviate this conservativeness, we add to the analysis a fourth case C4 (with FC4=FC1) where, as a difference from case C1, we impose φ= 0.1 (ψ≤0.1). Then, we obtain a fault diagnoser similar to C2 with a FAR in simulation of 10−4(see Table I), which is under the desired one of 10−3. This shows that we can compensate the conservativeness of the Markov’s approach by increasing the value of φand then verifying in simulation if the prescribed bound is fulfilled, but we cannot guarantee a priori a given tight false alarm rate or minimum detectable faults. Fig. 3 and Fig. 4 show the fault estimation and fault detection performances resulting from simulating the fault diagnosers from Table I under the appearance of two step faults, one for each channel, of an amplitude of 0.7at time t= 100 (disappearing at t= 400) for f1, and at t= 200 (disappearing at t= 500) for f2. The fault diagnosers for case C2 and C4 are the fastest ones to detect the faults and their estimation of the faults have the lowest settling time. However they are the most sensitive under state disturbances and measurement noises (as they have the highest φλ(F)product, see (27)). For case C1, the fault detector cannot detect the faults on time because it has a too slow dynamic due to the conservativeness introduced by the Markov’s inequality. Case C3, is an intermediate case between C1 and C2. Even if for case C3 the estimated faults converge slower to the faults than for cases C2 and C4, the detection mechanism only takes 6 more instants to detect the fault. This is due to the fact that C3 can detect lower faults than C2 and C4 (note that the diagonal of F−1 C3are higher than the ones of F−1 C2and F−1 C4). Finally, note that the settling time at the 98% (ξ= 2%) for the fault estimation, measured in terms of the number of measurement instants, is in the order of K(defined in (22)). For example, for case C3, the settling time is of 60 measurement instants for ˆ f1and of 130 for ˆ f2, while it was characterized by K= 167 from (22). VII. CONCLUSION In the current work, we designed a jump observer-based fault diagnoser to detect and estimate faults under measurements dropouts. We constructed the residual signal using a quadratic form of the estimated faults. A finite set of observer gains is used to estimate the faults and each gain is applied depending on the measurement outcomes. We employed the measurement successful reception probabilities from each sensor to describe the possible measurement reception scenarios. The proposed design method allows finding a trade off between the achievable minimum detectable faults and the response time to faults, while guaranteeing a prescribed false alarm rate. Two design strategies can be used: fixing the minimum detectable faults and then minimizing the response time, or fixing the response time and then minimizing the minumum detectable faults. We developed two ways of imposing a desired false alarm rate depending on the assumed knowledge about the probability distribution of the residual signal. If no information is assumed to be known, the Markov’s inequality leads to a very conservative bound on the false alarm rate. If the disturbances and noise are assumed to be Gaussian, a certain condition imposed on matrix Fleads to a chi-squared residual distribution. In this case a very precise bound on the false alarm rate is attained, improving the fault diagnosis performance. Further research may include extensions to delayed measurements with Markovian models for the missing measurements and analytical characterization of the missing fault rate.
9 TABLE I FAULT DIAGNOSERS COMPARISON. Case Design Simulation fmin fmin,1fmin,2φ ψ ρ KFAR C1 0.6 0.6 0.6 10−310−30.999 6101 0 C2 0.56 0.46 0.52 0.145 10−30.808 18 10−3 C3 0.21 0.29 0.29 0.145 10−30.977 167 10−3 C4 0.6 0.6 0.6 0.1 0.1 0.798 17 10−4 Estimated fault, f1[t] Case C1 and C4 Estimated fault, f2[t] simulation instants, t Case C2 simulation instants, t Case C3 simulation instants, t C4 C1 C1 C4 0200 400 6000 200 400 6000 200 400 600 -0.4 0 0.4 0.8 1.2 -0.4 0 0.4 0.8 1.2 Fig. 3. Fault estimation performances for the analyzed cases on Table I. simulation instants, t Residual signal r[t]and threshold rth Case C1 and C4 fault rth simulation instants, t Case C2 fault rth simulation instants, t Case C3 fault rth detection at t= 106 C1 C4 detection at t= 106 detection at t= 112 0 100 200 300 0 100 200 3000 100 200 300 0 5 10 15 20 25 Fig. 4. Fault detection performances for the analyzed cases on Table I. Fig. 4. Fault detection performances for the analyzed cases on Table I. APPENDIX Let us first introduce the following lemmas. Lemma 1 ( [34]).Let ωbe a stochastic vector with mean µ and a covariance matrix W, and Pa symmetric matrix. Then E{ωTPω}=µTPµ + tr(PW). Lemma 2 ( [35]).Let P be a positive semidefinite matrix, xia vector with appropriate dimensions and µi≥0scalar constants (with i= 1,2,...). If the series concerned is convergent, then we have ∞ X i=1 µixi!T P ∞ X i=1 µixi!≤ ∞ X i=1 µi!∞ X i=1 µixT iPxi. A. Proof of Theorem 1 Let us define the Lyapunov function at instant t=tkas Vk= ˜zT kP˜zk. i) In the absence of disturbances, faults and measurement noises, after taking Schur’s complements on (17a) and premultiplying the result by ˜zT kand postmultipliying by its transpose, we obtain that Vk+1 −Vk≤0that assures that the extended